Issue 
A&A
Volume 636, April 2020



Article Number  A110  
Number of page(s)  23  
Section  Interstellar and circumstellar matter  
DOI  https://doi.org/10.1051/00046361/201936039  
Published online  28 April 2020 
Visible and nearinfrared spectrointerferometric analysis of the edgeon Be star o Aquarii
^{1}
Université Côte d’Azur, Observatoire de la Côte d’Azur, CNRS, Laboratoire Lagrange, France
email: elisson.saldanha@oca.eu
^{2}
INAF – Osservatorio Astronomico di Brera,
Via E. Bianchi 46,
23807
Merate, Italy
^{3}
Univ. Lyon, Univ. Lyon1, ENS de Lyon, CNRS, Centre de Recherche Astrophysique de Lyon UMR5574,
69230
SaintGenisLaval, France
^{4}
Instituto de Astronomia, Geofísica e Ciências Atmosféricas, Universidade de São Paulo,
São Paulo,
Brazil
^{5}
CHARA Array  Georgia State University,
Mount Wilson,
CA, USA
Received:
6
June
2019
Accepted:
6
February
2020
Aims. We present a detailed visible and nearinfrared spectrointerferometric analysis of the Beshell star o Aquarii from quasicontemporaneous CHARA/VEGA and VLTI/AMBER observations.
Methods. We analyzed spectrointerferometric data in the Hα (VEGA) and Brγ (AMBER) lines using models of increasing complexity: simple geometric models, kinematic models, and radiative transfer models computed with the 3D nonLTE code HDUST.
Results. We measured the stellar radius of o Aquarii in the visible with a precision of 8%: 4.0 ± 0.3 R_{⊙}. We constrained the circumstellar disk geometry and kinematics using a kinematic model and a MCMC fitting procedure. The emitting disk sizes in the Hα and Brγ lines were found to be similar, at ~10–12 stellar diameters, which is uncommon since most results for Be stars indicate a larger extension in Hα than in Brγ. We found that the inclination angle i derived from Hα is significantly lower (~15°) than the one derived from Brγ: i ~ 61.2° and 75.9°, respectively. While the two lines originate from a similar region of the disk, the disk kinematics were found to be near to the Keplerian rotation (i.e., β = −0.5) in Brγ (β ~ −0.43), but not in Hα (β ~ −0.30). After analyzing all our data using a grid of HDUST models (BeAtlas), we found a common physical description for the circumstellar disk in both lines: a base disk surface density Σ_{0} = 0.12 g cm^{−2} and a radial density law exponent m = 3.0. The same kind of discrepancy, as with the kinematic model, is found in the determination of i using the BeAtlas grid. The stellar rotational rate was found to be very close (~96%) to the critical value. Despite being derived purely from the fit to interferometric data, our bestfit HDUST model provides a very reasonable match to noninterferometric observables of o Aquarii: the observed spectral energy distribution, Hα and Brγ line profiles, and polarimetric quantities. Finally, our analysis of multiepoch Hα profiles and imaging polarimetry indicates that the disk structure has been (globally) stable for at least 20 yr.
Conclusions. Looking at the visible continuum and Brγ emission line only, o Aquarii fits in the global scheme of Be stars and their circumstellar disk: a (nearly) Keplerian rotating disk well described by the viscous decretion disk (VDD) model. However, the data in the Hα line shows a substantially different picture that cannot fully be understood using the current generation of physical models of Be star disks. The Be star o Aquarii presents a stable disk (close to the steadystate), but, as in previous analyses, the measured m is lower than the standard value in the VDD model for the steadystate regime (m = 3.5). This suggests that some assumptions of this model should be reconsidered. Also, such longterm disk stability could be understood in terms of the high rotational rate that we measured for this star, the rate being a main source for the mass injection in the disk. Our results on the stellar rotation and disk stability are consistent with results in the literature showing that latetype Be stars are more likely to be fast rotators and have stable disks.
Key words: stars: individual: o Aquarii / stars: emissionline, Be / circumstellar matter / techniques: interferometric
© E. S. G. de Almeida et al. 2020
Open Access article, published by EDP Sciences, under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.
1 Introduction
Classical Be stars are mainsequence Btype stars that show (or showed at some time) Balmer lines in emission and infrared excess in their spectral energy distribution. The Be phenomenon is found among the entire spectral range of B stars (e.g., Townsend et al. 2004): M_{⋆} from ~3 M_{⊙} (B9, T_{eff} ~ 12 000 K), up to ~18 M_{⊙} (B0, T_{eff} ~ 30 000 K). These observational characteristics are well explained as arising from a dustfree gaseous disk that is supported by rotation with a slow radial velocity (see, e.g., Rivinius et al. 2013). The most successful theory to explain the evolution of the disk structure is the socalled viscous decretion disk (VDD) model, where its dynamics are driven by viscosity (e.g., Lee et al. 1991; Okazaki 2001; Bjorkman & Carciofi 2005).
It is widely accepted that fast rotation plays an important role in the formation of the Be star disk. However, while interferometric analyses typically provide rotational rates v_{rot} ∕v_{crit} ≳ 0.7 (e.g., Meilland et al. 2012; Cochetti et al. 2019), some statistical studies show rates ranging from ~0.3 up to 1.0 (e.g., Cranmer 2005; Zorec et al. 2016). Moreover, it is still not clear whether the rotational rate is correlated to other stellar parameters such as the effective temperature (e.g., Cochetti et al. 2019). Hence, despite the success of the VDD model, the physical mechanism(s) driving the mass injection remains unclear and a detailed physical characterization for the central star and the disk structure is mandatory to better understand the Be phenomenon. By gaining access to geometry on the milliarcsecond scale and kinematics on a few tens of km s^{−1} scale, spectrointerferometry offers a unique opportunity to probe the circumstellar environment and stellar surfaces of Be stars (see, e.g., Chesneau et al. 2012; Stee & Meilland 2012).
The bright, latetype Be star (type B7IVe) o Aquarii (HD 209409) is known to have a fairly stable disk (Sigut et al. 2015). The stability of the circumstellar disk is evidenced by the quasiconstant equivalent width in the Hα line, doublepeak separation, and the absence of longterm violettored (V∕R) peak variations (e.g., Rivinius et al. 2006; Sigut et al. 2015). This star shows a high value of v sin i ~ 282 km s^{−1} (Frémat et al. 2005) and a shell absorption in Hα, thus indicating a high stellar inclination angle of about 70°, as discussed below.
Meilland et al. (2012) presented the first spectrointerferometric analysis of o Aquarii with the VLTI/AMBER instrument as part of their AMBER survey of eight bright Be stars. Despite the low data quality and very limited number of observations (just one measurement), they were able to significantly constrain the disk geometry and kinematics. They found that the disk emission in the Brγ line, modeled as an elliptical Gaussian distribution,had a FWHM of 14 ± 1 D_{⋆} (with R_{⋆} = 4.4 R_{⊙}), where D_{⋆} and R_{⋆} are, respectively, the stellar diameter and radius. They estimated the inclination angle as i = 70 ± 20° and found a stellar rotational rate of v_{rot} ∕v_{crit} = 0.77 ± 0.21 (Ω/Ω_{c} = 0.93), where v_{crit} and Ω_{crit} are, respectively,the linear and angular critical velocity. New VLTI/AMBER spectrointerferometric measurements of o Aquarii were presented in the Be star survey of Cochetti et al. (2019). Here, they obtained seven goodquality measurements for o Aquarii (i.e., 21 baselines). Using a similar model as in Meilland et al. (2012), they found a Brγ emission FWHM significantly smaller than in Meilland et al. (2012), 8 ± 0.5 D_{⋆} (with R_{⋆} = 4.4 R_{⊙}), and better constrained the object inclination angle (70 ± 5°).
A detailed analysis of o Aquarii using Hα spectroscopy and interferometry was performed by Sigut et al. (2015). These authors combined large band (15 nm) interferometric data centered on Hα, obtained from the Navy Precision Optical Interferometer (NPOI), with Hα spectroscopy from the Lowell Observatory SolarStellar Spectrograph. Using the radiative transfer code BEDISK (Sigut & Jones 2007), they were able to reproduce simultaneously the visibility, Hα line profile, and spectral energy distribution (SED), and showed that the disk is quite stable for up to about ten years. Interestingly, they found a disk extension in Hα (Gaussian FWHM of 12.0 ± 0.5 D_{⋆}) close to the one determined by Meilland et al. (2012) in Brγ (FWHM of 14 ± 1 D_{⋆}). They concluded that this is uncommon since most previously studied Be stars exhibit a larger (up to two times) disk emission region in Hα than in Brγ.
In this paper, we present new CHARA/VEGA spectrointerferometric measurements of o Aquarii centered on the Hα emission line (λ = 0.656 μm). They are analyzed conjointly with the AMBER Brγ line (λ = 2.166 μm) measurements from Meilland et al. (2012) and Cochetti et al. (2019), using models of increasing complexity: simple geometric models, kinematic models, and radiative transfer models. This is the first time the code HDUST has been used to model simultaneously spectrointerferometric data from Hα and Brγ. It is the second time for the kinematic model (i.e., after the δ Scorpii data published in Meilland et al. 2011). This multiwavelength and multiline approach allows us to draw a more complete picture of the stellar surface and circumstellar environment of the Be star o Aquarii.
This paper is organized as follows. In Sect. 2, we present the observations and the data reduction process. Our analysis using geometric models of the VEGA calibrated (absolute) visibility is shown in Sect. 3. In Sect. 4, we fit the VEGA and AMBER differential visibility and phase with a kinematic model using a Markov chain Monte Carlo (MCMC) model fitting method. In Sect. 5, all the interferometric data are analyzed in terms of 3D nonLTE radiative transfer models. Our kinematic and radiative transfer models are discussed in Sect. 6. In Sect. 7, our bestfit models are compared to noninterferometric observables: the spectral energy distribution and line profiles (Hα and Brγ). The comparison with polarimetric data is performed in Sect. 8.4.3 in the context of the disk stability. In Sect. 8, we discuss the morphological, kinematic, and physical descriptions for o Aquarii and its circumstellar disk. Our conclusions are summarized in Sect. 9.
2 Observations
2.1 CHARA/VEGA
The VEGA instrument (Mourard et al. 2009) is one of the two visible beam combiners on the CHARA Array (ten Brummelaar et al. 2005). It can simultaneously combine up to four beams, operating at different wavelengths from 450 to 850 nm. VEGA is equipped with two cameras (blue and red detectors) that can observe in two different spectral domains simultaneously (around the Hβ and the Hα lines). Currently, it is the only instrument at the CHARA Array with a spectral resolution high enough to resolve narrow spectral features such as atomic and molecular lines. It offers 3 spectral modes: R = 1000 (LR), R = 6000 (MR), and R = 30 000 (HR).
o Aquarii was observed 50 times with VEGA between 2012 and 2016 in MR mode centered on the Hα emission line at 0.656 μm. The 2012 and 2016 observations were focused on the disk geometry and kinematics and data were taken with small baselines (up to 105 m) and without stellar calibrators. On the other hand, the 2013 and 2014 campaigns were aimed at constraining not only the Hα emission, but also the Rband continuum geometry. Consequently, observations were carried out with longer baselines (up to 330 m) with a standard calibration plan alternating observations of the science target and few calibrator stars chosen using the SearchCal (Bonneau et al. 2006) tool developed by the JeanMarie Mariotti Center (JMMC)^{1}. Table A.1 shows useful information about the stars used as interferometric calibrators during these campaigns. The complete log of observations is presented in Table A.2 and the corresponding uv plane coverage for the VEGA observations is plotted in Fig. 1.
Data were reduced using the standard VEGA data reduction software^{2} described in Mourard et al. (2012). For all programs, differential visibility and phases were computed from the intercorrelation between a fixed 15 nm window centered on Hα and a sliding smaller window (i.e. 1, 2, or 5 Å, depending on the data quality). For the 2013 and 2014 data, the raw squared visibility was computed for o Aquarii, and its calibrators, using the autocorrelation method on a 15 nm band centered on the Hα emission line (649–664 nm) and another band in the closeby continuum (635–650 nm). Then the transfer function was estimated assuming the diameter of the calibrators recorded before and after the science target observation, and its uncertainty using a weighted standard deviation. Finally, for each measurement, the calibrated squared visibility was derived by dividing o Aquarii’s raw squaredvisibility by the estimated transfer function.
Fig. 1
uv plan coverage obtained around Hα (0.656 μm) with CHARA/VEGA (green) and Brγ (2.166 μm) with VLTI/AMBER (red). 
2.2 VLTI/AMBER
The AMBER instrument (Petrov et al. 2007) was a threebeam combiner (decommissioned in 2018) at the Very Large Telescope Interferometer (VLTI). It operated in the H and Kbands with three spectral resolutions: R = 35 (LR), R = 1500 (MR), and R = 12 000 (HR). It offered the highest spectral resolution at the VLTI, being the most adapted for studying the gaseous environment in emission lines.
o Aquarii was observed with AMBER during two observing surveys of Be stars in 2011 (ESO program 087.D0311) and in 2014 (ESO program 094.D0140). The observations were performed in HR mode in Kband centered on the Brγ emission line at 2.166 μm. The data from 2011 was published in Meilland et al. (2012) and the 2014 data in Cochetti et al. (2019). During this second survey, seven measurements were acquired for o Aquarii with three different triplets. The log of AMBER observations is also presented in Table A.2 and the corresponding uv plane coverage is plotted in Fig. 1.
Calibration was performed using similar methods as the one described for VEGA. However, AMBER measurements were often affected by a highly variable transfer function mainly due to the variable quality of the fringe tracking performed by the FINITO fringe tracker, during the long exposure time needed to perform HR mode observations. As it was the case during our o Aquarii observations, we present in this paper only the analysis of differential measurements obtained using the standard AMBER data reduction software amdlib (Tatulli et al. 2007; Chelli et al. 2009).
3 Geometric modeling: VEGA calibrated visibility
In this section, we fit the Hα and continuum squared visibilities (V^{2}) from the VEGA observations where calibrators were observed. We note that as the AMBER data were not calibrated, such analysis cannot be performed on the Kband continuum and Brγ line.
To determine if we can separate the circumstellar disk and the stellar photosphere emissions and constrain their geometry independently, we fitted our data with geometric models of increasing complexity: onecomponent models (uniform disk, UD, or a uniform ellipse) and twocomponent models (UD plus UD, Gaussian disk, or uniform or Gaussian ellipse).
Here, the first component represents the stellar surface and the second one the circumstellar disk. To perform our fit, we used the LITpro model fitting software (TallonBosc et al. 2008) for optical and infrared interferometric observations developed by the JeanMarie Mariotti Center (JMMC)^{3}.
In Fig. 2, we show the comparison between the visibility curves of our bestfit models to the VEGA data both in the continuum and Hα bands. One sees that the object is partially resolved in the continuum and the Hα line. The lower level of the visibility in the band centered on the Hα line clearly shows that the object is larger in Hα than in the closeby continuum region. Assuming that the emission originates from both the stellar photosphere and a circumstellar disk, the lower visibility in Hα is due to a larger fraction of the Hα flux coming from the disk than from the star. In contrast, the flux contribution from the star is greater than that from the disk in the continuumRband.
Our main results are summarized in Table 1. We only show our results using UD models since there is no improvement in terms of reduced χ^{2} () when considering more complex models, that is, with a higher number of free parameters. For the continuum band, there is no significant improvement in terms of reduced χ^{2} between a simple UD and a twocomponent UD model. The central star is clearly resolved by the longer baselines and its extension is significantly constrained with a UD diameter of θ = 0.28 ± 0.01 mas ( ~ 1.1). This value corresponds to an upper limit to the stellar diameter measurements neglecting the putative contribution of the circumstellar disk in the Rband continuum. Adding a second component to the model only marginally reduces the extension of the first component. The contribution of the second component, representing the circumstellar disk, is small (F_{2} = 0.03 ± 0.03), thus the extension of the disk cannot be constrained.
Unlike the continuum case, the situation is quite different in the band centered on the Hα line. The single uniform disk gives a significantly higher ~ 2.8 for a bestfit model with θ = 0.36 mas. In this case, adding a second component reduces by a factor of two, leading to . Using a model with two uniform disks, we converge to a diameter of the first component similar to the one found from the continuum, that is, 0.26 ± 0.02 mas. The flux contribution of the second component and its extension are significantly constrained. However, the uncertainty remains quite large, that is, F_{2} = 0.15 ± 0.03 and θ_{2} = 6.5 ± 2.1 mas (see Table 1).
Considering that the first component of our model represents the stellar photosphere, our measurement is slightly higher than the value assumed in the work of Sigut et al. (2015) of 0.22 mas. However, their adoption for the stellar angular diameter is based on a spectral typeradius relation for B dwarf stars (Townsend et al. 2004). Moreover, this value of 0.22 mas represents the polar radius. o Aquarii is a fast rotator likely to be significantly flattened, and our measurements are spread over different orientations, so that we end up measuring a mean radius of the star projected on the sky. Assuming a distance of 144 pc (derived from the Gaia DR2 parallaxes, Gaia Collaboration 2018), θ_{1} = 0.26 ± 0.02 mas corresponds to a stellar radius R_{⋆} = 4.0 ± 0.3 R_{⊙}.
Finally, to try to detect any possible stellar or circumstellar disk flattening from the squared visibility measurements, we also computed individual uniform disk equivalent diameter for each V^{2} measurement. This analysis of the uniform disk diameter for o Aquarii, as a function of the VEGA baseline orientation, is shown in Fig. 3. As expected from our analysis (considering uniform elliptical models), we do not find any evidences of flattening from modeling our V^{2} dataset since no clear trends are found in the model residual as varying the baseline position angle.
Fig. 2
VEGA V^{2} measurements in the closeby continuum band (top) and in the Hα band (bottom) are shown in red points. Our bestfit models consisting of one (solid line) and two (dashed line) uniform disks are overplotted in blue. See Table 1 and text for discussion. 
Results from the geometric modeling of the VEGA V^{2} data in the closeby continuum band (635–650 nm) and in the band centered on Hα (649–664 nm).
Fig. 3
Top panels: uniform disk diameter derived from each individual VEGA V^{2} measurements (continuum band in the left and Hα band in the right) plotted as a function of the baseline position angle (PA). The red dotted line represents the bestfit diameter from modeling all the data in each band (θ = 0.28 mas in the continuum and θ = 0.36 mas in the Hα band). Bottom panels: corresponding normalized residuals. 
4 Kinematic modeling: VEGA and AMBER differential data
To constrain the geometry and kinematics of the circumstellar gas in the Hα and Brγ lines, we fit the VEGA and AMBER differential visibility and phase measurements using a simple bidimensional kinematic model for a rotating disk^{4}.
4.1 The kinematic model
This kinematic model was already used in a series of papers about spectrointerferometric modeling of Be stars, including Delaa et al. (2011), Meilland et al. (2012), and Cochetti et al. (2019), and is presented in detail in these references.
In short, the intensity map for the central star is modeled as a uniform disk, and the circumstellar disk as two elliptical Gaussian distributions, one for the flux in continuum, and the other one for the flux in line. The disk is geometrically thin so that the ellipse flattening ratio is set to 1∕cosi, where i is the inclination angle. The disk intensity map in the line is computed taking into account the Doppler effect due to the disk rotational velocity in the considered spectral channels. The parameters of our kinematic model are the following:
 (i)
The simulation parameters: size in pixels (n_{xy}), field of view in stellar diameters (fov), number of wavelength points (n_{λ}), central wavelength of the emission line (λ_{0}), step size in wavelength (δλ), and spectral resolution (Δλ).
 (ii)
The global geometric parameters: stellar radius (R_{⋆}), distance (d), inclination angle (i), and disk majoraxis position angle (PA).
 (iii)
The disk continuum parameters: disk majoraxis FWHM in the continuum (a_{c}), disk continuum flux normalized by the total continuum flux (F_{c}).
 (iv)
The disk emission line parameters: disk majoraxis FWHM in the line (a_{line}) and line equivalent width (EW).
 (v)
The kinematic parameters: rotational velocity (v_{rot}) at 1.5 R_{p} (polar radius) and exponent of the rotational velocity powerlaw (β).
4.2 Model fitting using the MCMC method
To perform our model fitting, we used the code emcee (ForemanMackey et al. 2013). This is an implementation in Python of the MCMC method from Goodman & Weare (2010). Some recent works on stellar interferometry used this code (see., e.g., Monnier et al. 2012; Domiciano de Souza et al. 2014, 2018; SanchezBermudez et al. 2017).
The simulation parameters were set as follows: n_{xy} = 256, fov = 60 D_{⋆}, n_{λ} = 60 (VEGA) and 110 (AMBER), λ_{0} = 6563 Å (VEGA) and 21 661 Å(AMBER), δλ = 2.5 Å (VEGA) and 1.0 Å (AMBER), and Δλ = 5.0 Å(VEGA) and 1.8 Å (AMBER). To reduce the number of free parameters, we set R_{⋆} = 4.0 R_{⊙} and d = 144 pc. We also fixed the disk continuum extension a_{c} and flux F_{c} to 0 for VEGA (i.e., neglecting the disk contribution in the continuum, based on our analysis of the VEGA V^{2} data). In the AMBER analysis, we adopted a_{c} = 3 D_{⋆} and F_{c} = 0.2 from Cochetti et al. (2019). The line equivalent width was set to 19.9 Å in Hα (Sigut et al. 2015). For Brγ, we computed the EW using the AMBER spectra from all observations and found a mean value of 13.6 ± 1.1 Å, which is compatible with the value from Meilland et al. (2012), 12.6 Å, but not with the result from Cochetti et al. (2019) of 18.1 Å. Finally, from the ten parameters of the kinematic model, the fitting of the VEGA and AMBER data were performed with at most five free parameters: i, PA, a_{line}, v_{rot}, and β.
The likelihood function (p_{like}) of the MCMC procedure was chosen as , where is the sum of the χ^{2} computed for the differential visibility and the differential phase. Thus, our attempt to converge to samples of parameters that maximizes the likelihood function means the minimization of the total χ^{2} between our interferometric data and the kinematic model.
We performed three different model fitting tests with different constraints on the value of v_{rot}:
 (i)
Five free parameters: i, PA, a_{line}, v_{rot}, and β. Without theinclusion of any prior probability function in the analysis.
 (ii)
Fourfree parameters: i, PA, a_{line}, and β. The stellar rotational velocity v_{rot} is fixed on the critical value of 391 km s^{−1} (Frémat et al. 2005).
 (iii)
Five free parameters: i, PA, a_{line}, v_{rot}, and β. We take into account a prior probability function p_{prior} on v sin i. Adopting μ = 282 km s^{−1} and σ = 20 km s^{−1}, from the measured v sin i = 282 ± 20 km s^{−1} (Frémat et al. 2005), we have the following expression for p_{prior}: (1)
where v sin i is calculated from the sampled MCMC values for the stellar rotational velocity and inclination angle.
Hence, considering a high weight on p_{prior}, the following quantity for the posterior probability function p_{post} is maxi mized: (2)
Note that this is equivalent to the case of equal weights for p_{prior} and p_{like}, but considering a lower error bar on v sin i, namely, σ = 2 km s^{−1}.
We typically used several hundreds of walkers (~300–900) for the MCMC run. Convergence was obtained for about 50 to 100 iteration steps in each walker, but we used a conservative value of 150 steps in the burnin phase and 50 in the main phase to estimate the parameters values and uncertainties. Overall, we found a mean acceptance fraction of ~0.5–0.6 in our MCMC tests. This is close to the optimal range for this parameter of ~0.2–0.5 (see, e.g., ForemanMackey et al. 2013).
Bestfit kinematic models from test (iii) for our VEGA (Hα) and AMBER (Brγ) differential data.
4.3 Bestfits in Hα and Brγ
We modeleda total of 117 (VEGA) and 24 (AMBER) measurements of differential visibility and phase. The bestfit parameters for the MCMC fit with a prior on v sin i (test iii, described above) are presented in Table 2. The corresponding histograms and the twobytwo parameter correlations from this MCMC run (one for VEGA and other for AMBER) are shown in Fig. 4. The corresponding histograms and correlation plots for the other two fits (tests i and ii) are shown in Figs. B.1 and B.2. One sees that the values of i, PA, and a_{line}, derived from each emission line, differ only marginally in all the fitting tests, showing the robustness of the solution for these parameters.
In Fig. 5, we show examples of VEGA and AMBER data in comparison to our bestfit kinematic models. For later discussion in Sect. 6, the visibility and phase from our bestfit HDUST model is also presented here. Our bestfit kinematic models are able to reproduce both the VEGA and AMBER differential data well. We found a reduced χ^{2} of ~4.0 and 1.6 from fitting, in a separate way, respectively, the VEGA and AMBER datasets.
We derived compatible values for the disk PA (~110°) from fittingthe VEGA and AMBER data with an uncertainty up to ~2°. This result agrees well with previous studies (e.g., Meilland et al. 2012; Touhami et al. 2013; Sigut et al. 2015; Cochetti et al. 2019). On the other hand, the inclination angle determined from the fit to the VEGA data is significantly smaller (i = 61.2 ± 1.8°) in comparison to the one determined from fitting AMBER (i = 75.9 ± 0.4°). This latter value is in good agreement with the results for i found by Meilland et al. (2012) and Cochetti et al. (2019). We also constrain the disk extension with a good precision: a_{line} = 10.5 ± 0.3 D_{⋆} in the Hα line and a_{line} = 11.5 ± 0.1 D_{⋆} in the Brγ line. These values are compatible with the ones determined by Sigut et al. (2015) in Hα and Meilland et al. (2012) in Brγ.
Another aspect concerning the disk extension in Brγ is the significant discrepancy seen in comparison to a_{line} = 8.0 ± 0.5 D_{⋆} from Cochetti et al. (2019). However, these authors used a larger value for the stellar radius of 4.4 R_{⊙} and a closer distance of 134 pc (van Leeuwen 2007), having thus the angular size of the stellar diameter larger in ~19% than the one assumed in our kinematic analysis from our results in Sect. 3. Considering all the other parameters fixed, this results in a smaller disk extension in ~19% than one found from our analysis. Nevertheless, the largest contribution to this discrepancy between our results and the ones from Cochetti et al. (2019) is due to their high value of equivalent width in the Brγ line of 18.1 Å, as discussed in Sect. 4.2, that also implies in a smaller disk extension in this line.
From our various tests, we showed that β and v_{rot} are strongly correlated. To precisely determine their dependence, we computed a grid of kinematic models varying just these two parameters in a regular step size. The values for i, PA, a_{line} are fixed from Table 2. The resulting maps are shown in Fig. 6. As expected, one sees that v_{rot} and β are highly correlated for the VEGA and AMBER data. This high degeneracy can be understood since these two parameters provide the rotational velocity structure in the disk: it is hard to distinguish the effects of each one on the modeling of spectrointerferometric (and spectroscopic) data.
Furthermore, we see that β = −0.5 (Keplerian disk) provides unrealistically high values for the stellar rotational velocity ((≳)400 km s^{−1} ; gray region) of o Aquarii (VEGA analysis). For AMBER, v_{rot} is significantly reduced to about 300–400 km s^{−1}. As shown in Fig. 6, our results from AMBER are consistent with a nearly Keplerian rotating disk (β ~ 0.43). However, it is conspicuous that the β value calculated from the VEGA data (β ~ 0.30) shows such a large departure from the Keplerian case.
Cochetti et al. (2019) derived a stellar rotational velocity of 355 ± 50 km s^{−1} and β = −0.45 ± 0.03. This is in fair agreement with our results for both v_{rot} and β. Considering our MCMC test (ii), where v_{rot} is fixed to the critical value and β is a free parameter, the results for β are shifted to higher values (more positive) with β ~ −0.42 (VEGA) and −0.54 (AMBER).
Therefore, regardless the MCMC fitting considered here, we verify a discrepancy of about 0.1 between the value of β derived from the Hα and Brγ lines. Our results from the AMBER analysis (Brγ) seems to be consistent with a nearly Keplerian rotating disk, but we verified a larger departure from β = −0.5 for the VEGA analysis (Hα).
Fig. 4
Histogram distributions and twobytwo correlations (after the burnin phase) for the free parameters of our bestfit kinematic models using MCMC for the VEGA (left panel) and AMBER (right panel) differential data. The median values are shown in solid red lines and the first and third quartiles in dashed red lines. The median and the first and third quartiles estimated for the parameters of our bestfit models (VEGA and AMBER) are presented inTable 2. In the correlation plots, darker points correspond to models with lower values of χ^{2}. See text for discussion. 
Fig. 5
Comparison between our bestfit kinematic models (dashed red; Table 2) and two different VEGA (top panels) and AMBER (bottom panels) measurements (black line). Our bestfit HDUST model is also shown (dashed blue; Table 5; discussion in Sect. 6). δλ of the kinematic model and AMBER data is increased to 1.8 Å in order to compare them to the HDUST model (δλ fixed to 1.8 Å). 
Fig. 6
maps of 40 000 kinematic models as a function of v_{rot} and β from the fit to VEGA (top panel) and AMBER (bottom panel) differential data. Only these two parameters were varied in a regular step in the intervals shown here. The other parameters are fixed (Table 2). Our results found from the MCMC analysis for v_{rot} and β are indicated with red crosses. In order to highlight the correlation between β and v_{rot}, the gray region corresponds to an arbitrary number of models, encompassing about the 5000 best models in both cases. The value of β = −0.5 (Keplerian disk) and our determination for v_{rot} are marked in dashed black line. Note the strong correlation between the stellar rotational velocity and the disk velocity law exponent in both the cases. Also, note that a Keplerian disk is found from modeling the AMBER data, but not from VEGA. 
5 Radiative transfer modeling
5.1 The code HDUST
We used the 3D nonLTE radiative transfer code HDUST^{5} (Carciofi & Bjorkman 2006, 2008) to perform a deeper physical analysis of o Aquarii. In addition to geometric and kinematic parameters, we seek to derive the density and temperature distributions in the disk, and the spectral energy distribution (SED), none of which was provided by the two simpler models considered in the two previous sections. HDUST uses a Monte Carlo method to solve the radiative transfer, statistical and radiative equilibrium equations for arbitrary density and velocity distributions in gaseous (pure hydrogen) or dusty circumstellar environments.
This code is wellsuited to model the circumstellar environment of Be stars as it implements the VDD model. Thus, the disk velocity law is assumed to be Keplerian (β fixed to −0.5). Many previous studies explored formal solutions of the VDD model in several limiting cases. For example, Bjorkman & Carciofi (2005) investigated the isothermal, steadystate case of a disk formed by a steady mass injection rate over a long time. Effects due to nonisothermal temperature structure were studied by Carciofi & Bjorkman (2008). Haubois et al. (2012) studied the temporal evolution of the disk structure that is subject to variable mass inject rates. Finally, the effects of a binary companion on the disk were studied by Okazaki et al. (2002), Oudmaijer & Parr (2010), Panoglou et al. (2016), and Cyr et al. (2017), among others.
From these studies, the radial density profile in Be star disks is found to be quite complex, for example, depending on the disk age, dynamical state, or presence of a binary companion. Despite this complexity, several studies have shown that the global behavior of this density profile is successfully approximated by a simple radial powerlaw (e.g., Touhami et al. 2009; Vieira et al. 2017). Considering also that the vertical density structure is that of an isothermal disk (hydrostatic assumption in the zaxis), the disk density can be parameterized as follows: (3)
where ρ_{0} is the disk base density, R_{eq} is the equatorial radius, and H(r) is the (isothermal) disk scale height given by: (4)
and H_{0} is the scale height at the disk base, (5)
where M_{⋆} is the stellar mass, G the gravitational constant, and c_{s} the sound speed velocity which depends on the local disk temperature T: (6)
where k_{B} is the Boltzmann constant, μ is the mean molecular weight of the gas, m_{H} is the hydrogen mass, and T is adopted as 0.72 T_{pol}, where T_{pol} is the polar effective temperature (see Correia Mota 2019).
HDUST has been used a few times to model spectrointerferometric observations (e.g., Carciofi et al. 2009; Klement et al. 2015; Faes 2015). From the solution of the radiative transfer problem, we are able to calculate synthetic spectra and intensity maps as a function of the wavelength around specific spectral lines. We estimated the stellar and circumstellar disk parameters from the comparison of our spectrointerferometric observations (visible and nearinfrared) with synthetic observables computed from the Fourier transform of HDUST monochromatic intensity maps.
HDUST parameters in the BeAtlas grid.
5.2 BeAtlas grid
Since a few hours are needed to compute a single HDUST model, it is not possible to perform an iterative model fitting procedure similar to the one described in Sect. 4. To overcome this issue, we used a precomputed grid of HDUST models called BeAtlas (Faes 2015; Correia Mota 2019). The BeAtlas grid is presented and described in detail by these references. It consists of ~14 000 models with images (specific intensity maps), SEDs, and spectra calculated in natural and polarized spectra, over several spectral regions, including the Hα and Brγ lines that are of interest for the analysis of our VEGA and AMBER dataset.
In Table 3, we show the parameter space covered by BeAtlas. Five physical parameters are varying in the grid. The stellar mass M_{⋆}, the inclination angle i, and the stellar oblateness R_{eq}∕R_{p}, fully describe the star. Other stellar parameters such as the stellar polar radius (R_{p}), rotational velocity (v_{rot}) and linear and angular rotational rates (v_{rot}∕v_{crit} and Ω∕Ω_{crit}) can be computed from M_{⋆} and R_{eq}∕R_{p} assuming rigid rotation under the Roche model (see, e.g., Carciofi & Bjorkman 2008). The two last parameters in Table 3describe the circumstellar disk structure and are parameterizations of the VDD model: the base surface density (Σ_{0}) and the radial density exponent (m).
The previously described volume mass density (Eq. (3)) and the surface mass density are related as follows: (7) (8)
From that, to facilitate the comparison to other disk models, we note that the relation between the volume and surface mass densities at the base of the disk is given by: (9)
The range of values for Σ_{0} and m in the grid encompasses somewhat extreme cases in the literature for the circumstellar disk of Be stars. For example, see Fig. 7 of Vieira et al. (2017). The listed values of Σ_{0} correspond to ρ_{0} from ~10^{−12} g cm^{−3} to ~10^{−10} g cm^{−3}. Parametric models with m = 3.5 are equivalent to the steadystate solution of the viscous diffusion equation considering an isothermal disk scale height. Thus, concerning the mass density law exponent m, models with m > 3.5 would represent a disk in an accretion phase, while the ones with m < 3.5 a disk in an ongoing process of dissipation (see, e.g., Haubois et al. 2012; Vieira et al. 2017).
5.3 Results
We performed four different analyses of our data using different subsets. For that, the reduced χ^{2} between the predicted interferometric observables from each HDUST model and the data was calculated as follows:
 (i)
calibrated VEGA V ^{2} in the 642.5 nm band (closeby continuum to Hα).
 (ii)
VEGA differential visibility and phase (Hα line).
 (iii)
AMBER differential visibility and phase (Brγ line).
 (iv)
All the quantities above analyzed together.
Analysis (i) was performed to evaluate the constraint on the stellar mass M_{⋆} and oblateness R_{eq} ∕R_{p}. In Fig. 7, we show the lowest value of for each value of stellar oblateness and mass from fitting the VEGA V^{2} data in the continuum band. The predicted V^{2} from our bestfit BeAtlas model (with M_{⋆} = 4.2 M_{⊙} ; Table 5) is overplotted to the VEGA measurements. For comparison, the predicted visibility curve from the BeAtlas model with the highest stellar mass, M_{⋆} = 14.6 M_{⊙}, is also overplotted to the data. These two models have the same values of i, R_{eq} ∕R_{p}, Σ_{0}, and m. In Sect. 3, we presented a similar analysis, but in terms of simple geometric models. For better visualisation, we show in Fig. 7 the local regression fits of as a function of R_{eq}∕R_{p} and M_{⋆}. Like all such calculations in this paper, all these regression fits of are performed with the LOESS method^{6}.
As in the analysis with geometric models, we cannot constrain the stellar oblateness using VEGA V^{2} data. On the other hand, the mass is better constrained with M_{⋆} ~ 4.8 M_{⊙} (B6 dwarf). From Fig. 7, one sees how the measured V^{2} are mismatched by the HDUST model with M_{⋆} = 14.6 M_{⊙} (unrealistic mass value for o Aquarii) due to the larger polar radius of ~7.4 R_{⊙} in this model. Among all the values for M_{⋆} in the grid, M_{⋆} = 4.2 M_{⊙} corresponds to a B7 dwarf star (Townsend et al. 2004). Since o Aquarii shows luminosity class IIIIV, it could be expected to have a mass somewhat higher than a dwarf of same spectral type, which is compatible with our results.
In Fig. 8, we show the lowest for each value of disk majoraxis position angle PA from the fit to the VEGA and AMBER differential visibilities and phases: analyses (ii) and (iii). Here, the stellar mass is fixed to M_{⋆} = 4.2 M_{⊙} from analysis (i), which also allows a better comparison to other studies of o Aquarii (e.g. Sigut et al. 2015). In both cases, of the models is minimized around PA = 110°, a value that we adopt in the remaining of this section. This is in good agreement to our results found with the kinematic model in Sect. 4.
In Fig. 9, we present our results from modeling the VEGA and AMBER differential visibility and phase in a separate way – analyses (ii) and (iii) – as well as from the simultaneous fit to all the interferometric data (analysis (iv)). The lowest is shown as a function of the following HDUST parameters: the inclination angle, stellar oblateness, base disk surface density, and the radial disk density law exponent. In Table 4, we show the statistics from these parameters calculated from the HDUST models within a certain threshold of , which, in each case, is chosen to match a similar number of models (~15–20 bestmodels). In Table 4, the parameters for the models with the lowest value of are also shown. In Table 5, we show the parameters for the best BeAtlas model to explain simultaneously all our different interferometric datasets.
Since our HDUST analysis is limited to the precomputed BeAtlas grid (limited parameter space and selected parameter values), we stress that the results presented here do not correspond to the real χ^{2} minimum toexplain our datasets in the framework of HDUST. Furthermore, the values for the standard deviation are shown in parenthesis in Table 4 since these are not determinations for the error bars on the parameters. They are just an evaluation for the dispersion on the parameters values of the BeAtlas bestmodels (within in a certain threshold of ). For example, from fitting AMBER, we found that all the BeAtlas models have Σ_{0} = 0.12 g cm^{−2}, and m = 3.0, up to, respectively, the top 207% and top 240% bestmodels. For this reason, it is shown, in this case, null standard deviation in Table 4 for these parameters (top 54% bestmodels).
From the separate analysis of the VEGA and AMBER differential datasets, we are able to describe the stellar and disk parameters, in Hα and Brγ, by the same HDUST model with: R_{eq}∕R_{p} = 1.45, Σ_{0} = 0.12 g cm^{−2}, and m = 3.0. One clear exception is found for the inclination angle. From the Hα analysis, is minimized for i = 56.3°. On the other hand, this is achieved with i = 77.2° in the Brγ line. Such discrepancy of ~20° is in agreement with the one found from our kinematic modeling. As expected, the joint analysis to all the data provides an intermediate mean value of ~65° for the inclination angle, showing a larger dispersion (higher standard deviation) in comparison to the results found from the separate analysis for VEGA and AMBER. One sees that the mean value for stellar oblateness is somewhat decreased, when considering all the datasets. However, in this case, the dispersion is significantly increased (±0.09) when compared to the separate VEGA and AMBER differential fits (±0.05–0.07). This happens due to the inclusion of the calibrated VEGA data in the joint analysis that do not allow us to properly infer this parameter (see, again, Fig. 7).
Fig. 7
Analysis (i): lowest value of reduced χ^{2} () for each value of R_{eq}∕R_{p} (left panel) and stellar mass (middle panel) from the HDUST fit to the VEGA V^{2} data (642.5 nm band). Local regression fits to , as a function of the parameter values, are shown in red line. In the right panel, the predicted visibility from our bestfit HDUST model (red points; Table 5) is compared to the VEGA V^{2} measurements in the continuum band (black points). The predicted visibility from the HDUST model with the highest mass in the BeAtlas grid (14.6 M_{⊙}, highest in the middle panel) is shown in blue points. 
Fig. 8
Lowest value of for each value of disk majoraxis position angle from the HDUST fit to the VEGA (top panel, analysis (ii)) and AMBER (bottom panel, analysis (iii)) differential visibility and phase. Local regression fits of as a function of the disk PA are shown as a red line. 
Fig. 9
Lowest value of for each value of stellar inclination angle, oblateness, base disk surface density, and disk density law exponent from the HDUST fit to the: VEGA differential data (top, analysis (ii)), AMBER differential data (middle, analysis (iii)), and all the interferometric data considered in this section (bottom panel, analysis (iv)). The stellar mass is fixed to 4.2 M_{⊙} and disk PA to 110°. Local regression fits to , as a function of the parameter values, are shown as a red line. The mean parameter values for the sets of best models (Table 4) are marked in dashed black line. Our bestfit BeAtlas model to fit all the interferometric data is shown in Table 5. See text for discussion. 
First three columns: mean and standard deviation values for each HDUST parameter of the BeAtlas grid: from analysis (ii) (19 bestfit HDUST models), analysis (iii) (16 bestfit HDUST models), and analysis (iv) (17 bestfit HDUST models).
Parameters of our bestfit HDUST model in the BeAtlas grid to explain the joint analysis of our interferometric data: VEGA calibrated and differential data and AMBER differential data.
Fig. 10
Intensity maps of our bestfit HDUST and kinematic models at different wavelengths around the Hα line (first two rows) and the Brγ line (last two rows). Flux/pixel is in arbitrary units with the same scale in Hα and Brγ. The image integrated in wavelength around each of these lines (Δλ = 2.7 nm around Hα and 3.9 nm around Brγ) are shown in the second column. 
6 Comparison between kinematic and HDUST bestfit models
In Fig. 5, we compare the synthetic differential visibility and phase from our bestfit kinematic and HDUST models to the actual VEGA and AMBER data for a few baselines. Comparisons to noninterferometric observables (spectral energy distribution and line profiles) are presented in Sect. 7. Our bestfit models are compared to all the AMBER data in Fig. C.1. One sees that our bestfit kinematic models do a better job of reproducing both the VEGA and AMBER data. From the separate kinematic modeling of the VEGA and AMBER differential data, the of the model is lower than with HDUST (BeAtlas grid). Fixing the stellar mass to a reliable value for o Aquarii (4.2 M_{⊙}), our bestfit HDUST model has ~ 6.1 and 4.7 for VEGA and AMBER, respectively. From the kinematic modeling, we found ~ 4.0 and 1.6 to explain these same datasets.
For VEGA, in particular, our bestfit HDUST model adjustment for the measured visibility width is worse than with the kinematic model. This particular issue in modeling the VEGA data can be explained; in HDUST, the disk velocity law exponent is fixed by β = −0.5 (Keplerian disk rotation), while in the kinematic model it is a free parameter. As shown in Sect. 4.3, we find values for β that are higher than −0.5, and this is accentuated from the analysis of the VEGA data (β ~ −0.3).
Apart from this issue regarding the analysis in Hα, we are able to describe well the disk density with the same physical parameters in both the Hα and Brγ lines: Σ_{0} = 0.12 g cm^{−2} and m = 3.0. As will be later discussed, this result found using HDUST is consistent with the ones presented in Sect. 4.3, showing a similar disk extension in these lines.
In Fig. 10, the intensity maps for each model are shown at the closeby continuum region and at different wavelength values in both the Hα and Brγ emission lines. The integrated intensity map (around each of these lines) is also presented. For a more realistic comparison, here we consider our bestfit kinematic model with a small flux contribution of 5% from the disk in the continuum nearby to Hα and a_{c} = 2 D_{⋆}. As shown in Table 2, these parameters were adopted as null in the kinematic analysis for the VEGA data, since we were not able to resolve the disk from our analysis of VEGA V^{2} measurements in the continuum band (Sect. 3). Regarding the continuum region close to Brγ, the disk extension and flux contribution are given in Table 2 for the AMBER analysis.
The major difference between the intensity maps in Hα and Brγ is the disk flattening which is due to the different inclination angle derived from these two regions, i ~57° (Hα) and ~72° (Brγ), from the best models provided in Table 4. Moreover, as seen in the images, the stellar flattening is taken into account in the HDUST modeling, but not in the kinematic model (the star is modeled as a uniform disk). Apart from these departures, we see that our bestfit HDUST model presents a fairly similar distribution to the one computed with the kinematic code: a Gaussian distribution represents the circumstellar disk. This can be better noted considering the full integrated images around the emission lines.
7 Comparison to noninterferometric observables
In this section, we compare our bestfit models, found from the analysis of interferometric observables, to the observed spectral energy distribution (SED) and line profiles (Hα and Brγ) of o Aquarii. With respect to polarimetric data, it is discussed in Sect. 8.4.3 when addressing the disk stability.
Fig. 11
Comparison between the observed o Aquarii and model SEDs from the ultraviolet to the farinfrared region. Flux unit is in erg cm^{−2} s^{−1} Å^{−1} and wavelength is shown in logarithmic scale. IUE/SWP and IUE/LWP spectra are shown in black line and photometric data in black points. Top panel: purely photospheric models (color lines) with variation in the stellar radius (no inclusion of geometrical oblateness): R_{⋆} = 3.2 R_{⊙} (orchid), 4.0 R_{⊙} (red), and 4.4 R_{⊙} (green). Bottom panel: photospheric model with 4.0 R_{⊙} (red) and our bestfit HDUST model from fitting all the interferometric data (dashed blue line; Table 5). Note that the UBVbands are better reproduced with R_{⋆} = 4.0–4.4 R_{⊙}. Our best HDUST model reproduces the observed IR excess due to the circumstellar disk well. 
7.1 Spectral energy distribution
In Fig. 11, we present the spectral energy distribution (SED) of o Aquarii from the ultraviolet (IUE/SWP and IUE/LWP spectra^{7}) to the farinfrared region. References for the photometric data are given as follows: UBVJHKbands (Anderson & Francis 2012), iband (Henden et al. 2016), LMbands (Bourges et al. 2017), and IRAS 12, 25, and 60 μm bands (Abrahamyan et al. 2015).
For the spectral region up to the Vband, we compare the data to the SEDs of purely photospheric atmosphere models with solar metallicity (Castelli & Kurucz 2004). In this region, the circumstellar disk flux level is much lower than the photospheric flux, thus allowing a proper probe of the stellar radius (e.g., Meilland et al. 2009). The surface gravity was fixed at log g = 4.0, this being the closest value in Castelli &Kurucz (2004) to logg = 3.9 that is given by our results of M_{⋆} = 4.2 M_{⊙} and R_{⋆} = 4.0 R_{⊙}. The effective temperature was fixed at 13 000 K, following Cochetti et al. (2019). As in the previous sections, we consider the distance to be 144 pc, from the Gaia DR2 parallax.
These synthetic SEDs were calculated for three different stellar radius values, R_{⋆} : 3.2 R_{⊙} (Sigut et al. 2015), 4.0 R_{⊙}, and 4.4 R_{⊙} (Cochetti et al. 2019). The value of 4.0 R_{⊙} corresponds to the stellar radius determined from the fit to the VEGA V^{2} data using a twocomponent model: 4.0 ± 0.3 R_{⊙}. The effect of interstellar medium extinction is not included in these models since it is negligible for o Aquarii. Assuming a total to selective extinction ratio of R_{V} = 3.1, Touhami et al. (2013) derived a color excess of E(B−V) = 0.015 ± 0.008 for this star from their fit to the SED. This means the observed flux is ~96% of the intrinsic one in the Vband (lower by ~0.02 dex). It is beyond the scope of this paper to estimate the extinction due to the circumstellar disk, however, from the comparison to purely photospheric models, we see in Fig. 11 that the effect of extinction (due to the interstellar and circumstellar matter) is conspicuously weak on the 0.220 μm bump.
From Fig. 11, we see that the UV and visible regions are better reproduced for a stellar radius of about 4.0–4.4 R_{⊙}, when compared to 3.2 R_{⊙}, adopted in Sigut et al. (2015), which corresponds to the expected polar radius for a B7 dwarf. We stress that the radius derived by Cochetti et al. (2019) is closer to our results from the fit to the VEGA V^{2} data (Sect. 3). Their result of R_{⋆} = 4.4 R_{⊙} corresponds to a uniform disk diameter of θ ~ 0.28 mas (d = 144 pc). A better comparison to Cochetti et al. (2019) is hard since they do not provide error bars on R_{⋆} from fitting the SED. Furthermore, they derived R_{⋆} = 4.4 R_{⊙} for o Aquarii using a distance of 134 pc from van Leeuwen (2007), rather than the value of 144 pc adopted here. From Fig. 11, this implies a larger discrepancy between the observed and synthetic SED for R_{⋆} = 4.4 R_{⊙}, overestimating the observed flux.
We also compare the predicted SED of our bestfit HDUST model (Table 5) to the SED of the purely photospheric model with 4.0 R_{⊙}. Despite being able to reproduce the UBVbands well, one sees that a purely photospheric model clearly underestimates the observed flux beyond the nearinfrared due to the flux contribution from the circumstellar disk (e.g., Poeckert & Marlborough 1978; Waters 1986). From Fig. 11, it is evident that the SED is much better reproduced up to the farinfrared region when taking into account the IR excess from the gaseous circumstellar disk present in our bestfit HDUST model.
7.2 Hα and Brγ profiles
Our Hα spectra taken with the VEGA instrument (20 spectra, period from 2012 to 2016) are not analyzed in this work since they are saturated. This is a known effect seen in previous works on Be stars and correlated to the magnitude of the object. We stress that this instrumental saturation effect does not impact the visibilities and phases extracted from the fringes measured with VEGA (see, e.g., Delaa et al. 2011). To overcome this problem we used Hα line profiles from the BeSOS^{8} catalog (Arcos et al. 2018; Vanzi et al. 2012), obtained between 2012 and 2015, and thus covering a similar period to our VEGA observations. The typical spectral resolution of the BeSOS spectra is ~0.1 Å.
In Fig. 12, we compare the Hα and Brγ profiles from our bestfit models to observed profiles, namely, the mean Hα line profiles from BeSOS (7 profiles^{9}) and the mean Brγ line profiles from our AMBER observations (8 profiles). The observed profiles in Fig. 12 were binned in wavelength in order to have a spectral resolution equal to one of the synthetic profiles from the kinematic and HDUST models: 1.3 Å (Hα) and 1.8 Å (Brγ). The mean EW in Hα from the BeSOS data is 19.1 Å. This is in agreement with the mean value of 19.9 Å found in Sigut et al. (2015), based on contemporaneous spectra, and adopted in our analysis with the kinematic code (Sect. 4.3).
First, we note that our bestfit kinematic and HDUST models provide a fairly reasonable match to the observed Hα and Brγ line profiles.The kinematic models correspond to our bestfits obtained from modeling the VEGA and AMBER differential data separately (Sect. 4). On the other hand, our bestfit HDUST model shown in Hα and Brγ is derived from the simultaneous fit to all our interferometric data (Table 5). Moreover, we stress the difficulty found by Sigut et al. (2015), using the radiative transfer code BEDISK, to reproduce the line wings and central absorption in the Hα profile of o Aquarii (see their Fig. 5).
However, it can be seen in Fig. 12 that both our bestfit kinematic and HDUST model are not able to properly reproduce, in particular, the wings of the Hα profile. On the other hand, the wings of the Brγ profile are fairly well reproduced by both of them, especially with HDUST.
Therefore, this inability to reproduce the wings of the Hα profile well is likely due to physical processes in the disk that are not taken into account in our models. It is known that the Hα profile wings of Be stars can be highly affected by noncoherent scattering, thus resulting in nonkinematic linebroadening in this transition (see, e.g., Hummel & Dachs 1992; Delaa et al. 2011). It is beyond the scope of this paper to quantify this possible effect in the Hα line of o Aquarii.
Fig. 12
Comparison between our bestfit kinematic models (dashed red; Table 2) and HDUST model (dashed blue, Table 5) in the Hα and Brγ line profiles. Mean observed line profiles of Hα (BeSOS) and Brγ (AMBER) are shown in black line. Our bestfit kinematic and HDUST models provide reasonable synthetic profiles to the observed ones in both Hα and Bγ. 
8 Discussion
8.1 Disk extension in Hα and Brγ
In Sect. 4.3, we showed that the disk extension is similar in the Hα and Brγ lines. Interestingly, from previous studies, we could expect to find a larger disk extension in Hα than Brγ. For example, Meilland et al. (2011) found that δ Scorpii (B0.3IV), which was also observed with the VEGA and AMBER instruments, shows a circumstellar disk 1.65 times larger in Hα than in Brγ. Furthermore, Gies et al. (2007) derived the angular sizes of four Be stars (γ Cassiopeiae, ϕ Persei, ζ Tauri, and κ Draconis) in the Kband region using interferometric data from the CHARA/CLASSIC instrument. They showed that the disk of these stars was significantly larger (up to ~1.5–2.0 times) in the Hα line than in the Kband. However, Carciofi (2011) investigated theoretically, using the code HDUST, the formation loci of Hα and Brγ, and found them to be quite similar at least in the parameter space explored by the authors (see their Fig. 1). Moreover, Stee & Bittar (2001), using the code SIMECA, found that Be star disks can be larger (up to two times) in Brγ than in Hα.
For a quantitative comparison of the disk extension in Hα and Brγ, we fitted simple Gaussian distributions to the intensity map of our bestfit HDUST model for all the values of inclination angle in BeAtlas. In order to remove the contribution from the star and disk continuum, we removed the image from the continuum before performing the fit and we hide the central part of the image which is affected by the stellar contribution.
In Fig. 13, we show the majoraxis FWHM from our fit as a function of the inclination angle for the Hα and Brγ lines. First, one sees that the disk sizeextension (majoraxis FWHM) varies differently in the Hα and Brγ lines as a function of the inclination angle. The disk extension increases in Brγ with the inclination angle. On the other hand, it decreases significantly in Hα up to i ~ 56° and increases after this value. One sees that the ratio between the extension in these lines decreases from about 1.50 at zero inclination to about 1.05 at 63.5°. Furthermore, we note that the disk extensions in these lines are very close to each other for i ~ 56° (Hα) and i ~ 72° (Brγ): majoraxis FWHM ~ 2.45 mas. Considering d = 144 pc, the disk size is ~10 D_{⋆} (close to our findings from the kinematic modeling).
Therefore, from this simple analysis using HDUST models, we verify our findings using the kinematic model: a similar circumstellar disk extension in Hα and Brγ. This arises since the (equivalent) Gaussian disk to our bestfit HDUST model presents quite different changes on its extension in these lines as a function of the inclination angle. Based on that, we can also explain the difference between δ Scorpii and o Aquarii. The former is seen under a low inclination angle (~30°) and exhibits a high ratio between the Hα and Brγ disk sizes. The latter is seen under a higher inclination angle and shows similar disk sizes in both lines. On the other hand, as discussed above, ϕ Persei and ζ Tau show larger disks in Hα than in the Kband and these stars are seen close to edgeon with i = 78° (Mourard et al. 2015) and 85° (Carciofi et al. 2009), respectively. Thus, this similarity in the disk extensions, found for quite different values of inclination angle, could indicate a more complex physical structure of the circumstellar disk than the one assumed by our bestfit HDUST model (based on a vertically isothermal disk).
Fig. 13
Majoraxis FWHM of Gaussian distribution (fitted from our bestfit HDUST model) as a function of the HDUST inclination angle. All the other HDUST parameters are fixed. Blue points correspond to the fit in Hα and red points in Brγ. The vertical dashed lines mark our values for inclination angle derived from the HDUST analysis, fitting the data in Hα (blue) and Brγ (red). Note that the equivalent Gaussian fits show a similar extension (2.45 mas, marked in horizontal dashed line) for these values of i. 
8.2 Inclination angle and vertical disk structure
From our Hα and Brγ differential data analysis, using the kinematic model, we achieved good precision in the determination of the stellar inclination angle: i ~ 61.2 ± 1.8° (VEGA) and i = 75.9 ± 0.4° (AMBER). Nevertheless, there is a clear discrepancy between the inclination angle found from fitting the VEGA and AMBER datasets. The value determined from VEGA is about 15° lower than the one found in the analysis of the AMBER data. We can show that this issue does not stem from an intrinsic limitation of the kinematic code (2D model) for Be stars seen under high inclination angle (i ≳ 60°). Indeed, by using a sophisticated 3D radiative transfer model (HDUST), not subjected to such a limitation, we verified the same discrepancy on i from the fit to these datasets separately (see, again, in Fig. 9, the trend of , as a function of i).
It may be argued that the difference found in inclination angle is hiding a difference in the disk thickness in these lines. Assuming a nongeometrically thin disk, for an ellipse with major and minor axes denoted, respectively, by a and b, the ratio between a and b, the circumstellar disk flattening, is given by (see, e.g., Meilland et al. 2007b): (10)
where i is the stellar inclination angle and Θ the disk opening angle. Since the i derived from Hα using our physical models is much lower than from Brγ (a reliable value when compared to other results in the literature), this would imply a disk thicker (higher opening angle Θ) in Hα than in Brγ. Considering the values described above, the disk opening angle in Hα would be Θ ~ 37° larger in Hα than in Brγ (assuming a geometrically thin disk in Brγ). Such a high value of opening angle is far beyond what is measured and expected by the VDD model, typically less than ~10° (cf. Rivinius et al. 2013). This might indicate the necessity of more complex physical assumptions in the physical properties of our disk model.
Since the code HDUST provides a pure hydrogen modeling for the photosphere plus disk regions, this disagreement between the VEGA (Hα) and AMBER (Brγ) analyses in the determination of i could be due to an opacity effect. It is wellknown that the inclusion of heavy elements can impact the density and temperature stratifications in the circumstellar disk of Be stars by shielding emission from the star. (see, e.g., Sigut & Jones 2007). Furthermore, we stress that our bestfit HDUST model is a parametric model (based on a vertically isothermal structure). Departures from vertically isothermal disks are wellknown in the literature. For example, using the radiative transfer code BEDISK, Sigut et al. (2009) verified that isothermal and selfconsistent hydrostatic models can present large differences regarding the temperature stratification in the disk of Be stars. Using HDUST, Carciofi & Bjorkman (2008) also found that nonisothermal effects can be significant for denser Be star disks. Thus, further investigation is needed concerning this effect on the determination of i for o Aquarii, but that is beyond the scope of this paper.
Finally, another possibility to explain the difference in apparent inclination angle found in our modeling could be a nonnegligible contribution of a polar wind. Clues of the presence of polar wind, or at least of circumstellar material in the polar regions, have been found by Kervella & Domiciano de Souza (2006) and Meilland et al. (2007a). In our models, we assume that all the circumstellar material is in the thin equatorial disk. If a non negligible fraction of the material is located near the poles, although we would expect it to be quite diluted and optically thin (at least in the continuum), it might affect the line emission with a different magnitude in Hα and in Brγ. If one assumes that the hydrogen level populations favor Hα emission over Brγ, the polar contribution of Hα would be higher, and the environment might look less flattened in this line than in Brγ.
8.3 Stellar and disk rotation
In Sect. 5, our results are presented in terms of the stellar obla teness R_{eq} /R_{p} (denoted by f in Eq. (11)). First, we give the relation between the oblateness and the angular Ω∕Ω_{crit} and linear v_{rot}∕v_{crit} rotational rates as follows: (11)
where R_{eq,crit} and R_{eq} (in units of polar radius) are, respectively, the stellar equatorial radius in the case of critical velocity and the actual one (see, e.g., Frémat et al. 2005; Ekström et al. 2008).
Considering only the uncertainties on v_{rot} (Table 2), with the critical velocity v_{crit} fixed to 391 km s^{−1} (Frémat et al. 2005), we obtain a linear rotational rate of v_{rot} /v_{crit} = 0.83 ± 0.02 (v_{rot} = 325 ± 6 km s^{−1}, VEGA) and 0.775 ± 0.005 (v_{rot} = 303 ± 2 km s^{−1}, AMBER). From the HDUST analysis, we find v_{rot} /v_{crit} = 0.96 (VEGA and AMBER) from our bestfit model (no error bars). This difference between the kinematic and HDUST analysis can be explained since the β exponent (velocity law in the disk) is fixed in the HDUST analysis (Keplerian disk, β = −0.5), while it is a free parameter in the kinematic model. We derived values for β from the kinematic analysis that are significantly higher (more positive) than −0.5 (see Table 2).
Apart from these differences, our analysis is consistent with a high rotational rate for o Aquarii, showing v_{rot} /v_{crit} from ~0.8 up to 1.0, depending on the particular analysis considered. The BeAtlas fits to the VEGA and AMBER differential data are significantly worsened (Fig. 9), when considering R_{eq} /R_{p} = 1.20–1.30 (Ω∕Ω_{crit} = 0.88–0.96). Thus, our HDUST analysis indicates that o Aquarii rotates faster than Ω∕Ω_{crit} = 0.96, disfavouring the lower range of Ω∕Ω_{crit} between 0.86 and 0.93 that is derived by Cochetti et al. (2019).
In Sect. 4.3, we found a strong correlation between the velocity at the base of the disk and the β exponent of the rotation law β. The inferred degeneracy, stronger in the case of the VEGA data, which have a lower spectralresolution with higher uncertainties, prevents us from independently constraining these two parameters with our kinematic model when fittingonly our spectrointerferometric data. However, the addition of an external constraint, the measured v sin i, removed this degeneracy, allowing us to derive a more accurate value of β in comparison to the other MCMC fitting tests (Appendix B).
From our MCMC fit to the AMBER dataset, with a preset v sin i, we derived a β of −0.426 ± 0.003. Thus, the disk appears to be rotating in a nearly Keplerian fashion. Despite this very low error on β, note that the error bars on β change with respect to the presented MCMC tests, up to ± 0.008 (see Fig. B.2). On the other hand, the value derived from the fit of the VEGA data is about 0.1 higher than from AMBER. We stress that this discrepancy cannot be explained by a radial dependent rotational law because both lines roughly stem from the same region in the disk (similar disk extensions in these lines).
Moreover, this apparent higher value of β in Hα was also theorigin of some biases that we found when modeling the VEGA differential data alone using HDUST. Without fixing M_{⋆}, the VEGA analysis with HDUST favours unrealistically high values of stellar mass up to ~11 M_{⊙}. This happened due to the fact the higher mass models also correspond to higher rotational velocity at the base of the disk. As the value of β is fixed to −0.5 in the BeAtlas grid of models, this was the only way to increase the rotational velocity in the disk. In Fig. 14, our bestfit kinematic and HDUST models are compared to a higher mass HDUST model for one VEGA differential measurement and the Hα profile. When compared to our best HDUST model (mass fixed to 4.2 M_{⊙}), the visibilitydrop and the Hα profile are better reproduced with HDUST models with higher value of mass, but also considering a larger value of Σ_{0} (0.28 g cm^{−2}). This happens since the stellar radius is also increased for a higher mass model and the flux contribution from the star is larger in the line. In this case, our BeAtlas model is able to produce more similar synthetic Hα visibility and profile to the ones from our bestfit kinematic model.
One possible explanation for the discrepancy between the value of β determined from Hα and Brγ could be the higher effects of nonkinematic broadening on Hα. This is already evidenced by the larger wings, in terms of Doppler shift, for this emission line. Such effects are known to be due to noncoherent scattering in the circumstellar environment, as explained, for example, in Auer & Mihalas (1968). Global effects on interferometric data were discussed by Stee et al. (2012) in the case of the Be star γ Cassiopeiae observed with VEGA. These authors used a similar kinematic model, but with two additional parameters to quantify the noncoherent scattering and found that about half the flux in the line was affected by such an effect. Never theless, the possible bias on the measurement of β in a line strongly affected by such nonkinematic broadening should be investigated further.
We also note that a possible close companion could influence v_{rot}, as well as the disk structure, as previously mentioned. However, the presence of a close companion with a detectable influence on the measured parameters seems excluded from the observed calibrated V^{2}, and in particularfrom the spectrointerferometric differential observables, which both show signatures well reproduced by a symmetric rotating disk.
Fig. 14
Bias effect of the disk velocity on Hα modeling. One VEGA measurement and observed Hα profile (BeSOS, as in Fig. 12) are shown in black lines. Our bestfit kinematic (dashed red) and HDUST (dashed blue) models are shown in Hα visibility and line profile. They are compared to HDUST models with a higher mass of 10.8 M_{⊙} with: Σ_{0} = 0.12 g cm^{−2} (dashed orchid) and Σ_{0} = 0.28 g cm^{−2} (dashed green). See text for discussion. 
Fig. 15
Top panel: 70 observed Hα profiles of o Aquarii from the BeSS database, covering about 17 yr of observations (2001–2018). Bottom panel: Hα equivalent width as a function of the observation time (modified Julian date). Civil dates are indicated for a part of the measured EW. The mean EW (solid) within the standard deviation (dashed lines) is marked in black. Local regression fit of EW, as a function of time, is shown as a solid red line. The time interval covered by our interferometric observations (VEGA) is indicated with dashed red lines. See text for discussion. 
8.4 Disk variability: a multitechnique analysis
8.4.1 Spectroscopy
The Be star o Aquarii is known to possess a stable Hα line profile for up to several years. For example, Sigut et al. (2015) verified that the EW in the Hα line is stable (within about 5%) up to about nine years (from 2005 to 2014).
To go further in the analysis of the disk stability, we analyzed 70 Hα line profiles, spanning from 2001 to 2018, from the BeSS database^{10} (Neiner et al. 2011). Since these observations are performed with several instruments, the line profiles shown here are interpolated to have spectral resolution of 0.5 Å (lowest resolution in the dataset). From these observations, we calculated the equivalent width (EW) in the Hα line. In Fig. 15, we show the analyzed Hα spectra together to the temporal evolution of the Hα EW.
We found that the disk is fairly stable over this 17yr time span with a mean value of = 18.1 ± 1.2 Å. This valueagrees well to older results in the literature. Slettebak & Reynolds (1978) measured Hα EW = 18.80 ± 0.11 Å in 1975 and 18.58 ± 0.21 Å in 1976. From the Hα profile observed in 1981, Andrillat (1983) measured EW = 17.2 Å. Thus, this supports an even longer global disk stability up to at least 40 yr. However, a slight increasing trend in EW is seen between 2001 and 2012. This could suggest an augmentation in the disk density of o Aquarii in this period. Considering the period of our interferometric observations (from 2012 to 2016), it is hard to observe any trend of Hα EW as a function of time.
Fig. 16
VEGA (top panels) and AMBER (bottom panels) differential visibilities extracted from observations at different epochs (black line). VEGA measurements span four years and the AMBER ones span three years. The observation date and the baseline length (projected onto the sky) are indicated in the top of each panel. Our bestfit kinematic models derived from the fit, in a separate way, to each dataset (VEGA and AMBER) are shown in dashed red line. Note that our bestfit kinematic models match well to the differential visibilities obtained at different epochs. 
8.4.2 Interferometry
These results are consistent with our ability to model, with the same model parameters, simultaneously all our VEGA and AMBER data regardless of the epoch. In Fig. 16, we present a temporal evaluation of our spectrointerferometric data in the Hα and Brγ lines. Since the drop in visibility is expected to change due to possible variations in the disk extension, we only show here the differential visibilities from the VEGA and AMBER observations. These measurements are chosen to cover the whole period of our observations from 2011 to 2016. For a more robust comparison, we chose measurements obtained with different baseline lengths (projected onto the sky), and, thus, covering different levels of spatial resolution.
One sees that, regardless of the period of time of the observations, our final kinematic models provide a very reasonable match to both the VEGA and AMBER data. Thus, considering our interferometric data, we are not able to detect any conspicuous variation of the circumstellar disk extension within a period up to five years (from 2011 to 2016). This is in agreement with previous interferometric studies of o Aquarii by Sigut et al. (2015). Besides that, this analysis supports our approach of fitting each one of the interferometric datasets (VEGA and AMBER) without imposing any discrimination based on the observation time.
Fig. 17
Polarimetric quantities of o Aquarii, as a function of the observation time, spanning about six years. Top panel: observed Vband polarization (44 measurements). Bottom panel: ratio between the observed B and Rbands polarization (40 measurements). The mean values of these quantities are shown in dashed line. See text for discussion. 
8.4.3 Polarimetry
Additional multiepoch polarimetric data also support our findings of a stable disk for o Aquarii, close to the steadystate regime. In Fig 17, we show the temporal evolution of broadband linear polarimetry in the Vband (P_{V}) of o Aquarii, as well as the ratio between the B and Rbands polarization (P_{B}∕P_{R}).
These data were obtained over 43 nights, from June 2010 to August 2016, with the IAGPOL polarimeter (Magalhães et al. 1996), mounted on the 0.6 m Boller & Chivens telescope at Observatório do Pico dos Dias (OPD/LNA). This polarimeter is composedby a rotating halfwave retarder and a Savart Plate used as analyser to provide the modulation of the light polarization, and then the polarimetric quantities. Details of data reduction are found in Magalhães et al. (1984) and Bednarski (2016).
From Fig. 17, the mean value of the observed Vband polarization is = 0.48 ± 0.03%. This value, derived from the mean and standard deviation of the Stokes Q and U parameters, is compatible to the one determined by Yudin (2001), namely, 0.52 ± 0.05%. Since the observations from Yudin (2001) predate our OPD/LNA observations by more than a decade, we conclude that the polarization values of o Aquarii remained very constant for over 20 yr.
In order to determine the intrinsic value of polarization, the interstellar contribution to the observed values quoted above must be removed. For that, we observed four main sequence stars, in the BVRIbands, which are angularly close to o Aquarii. A MCMC method was implemented to process the four BVRI data of each field star, generating a sample of the likelihood function in terms of the interstellar Serkowski parameters P_{max} and λ_{max} (Serkowski et al. 1975; Wilking et al. 1982). The best estimates for these parameters are shown in Table D.1.
There is a good agreement among the PA values of the field stars. Moreover, by using Gaia DR2 distances, we found that P_{max} increases linearly along the line of sight of o Aquarii (seeFig. D.1). In this case, it suggests that the alignment of the grains at the interstellar medium is nearly homogeneous (e.g., McLean & Clarke 1979). Thus, from a simple linear fit to P_{max} vs. distance for the field fields, we determined P_{max} for o Aquarii. The derived interstellar polarization parameters for o Aquarii are shown in Table 6, which are in reasonable agreement with the ones reported in Yudin (2001) of P_{max} = 0.20% and PA_{IS} = 125° (no error bars).
Taking into account our results for the interstellar polarization components, we found for the intrinsic Vband polarization and position angle = 0.49 ± 0.03% and PA ^{int} = 2.5 ± 2.7°, respectively. Yudin (2001) determined = 0.60% with PA^{int} = 6.0°, which is close to our PA^{int} value. Moreover, both estimates for PA^{int} are consistent with our determination for the disk majoraxis position angle (~110°), being almost perpendicular to the polarization vector, as expected.
Furthermore, our bestfit HDUST model (Table 5) predicts a polarization degree of 0.41% in the Vband. This agrees well with our measurement for the average intrinsic polarization of the OPD/LNA data. Therefore, besides the independent checks provided by the SED and spectroscopic data (Sect. 7), our polarimetric data also support our physical model for o Aquarii, which was derived purely from the fit to interferometric data (as discussed in Sect. 5.3).
Lastly, Fig. 17 shows that both the polarization degree in the Vband and the ratio between the B and Rbands are almost constant in time, showing a small scatter around the mean value. In particular, this latter quantity is related to the density scale at the inner portion of the disk (Haubois et al. 2014). From the theoretical investigation of Panoglou et al. (2019), the variation on the polarization degree in the Vband (ΔP_{V}) can reach up to about 0.1% due to asymmetries in the disk density structure, caused by a binary companion. Moreover, Haubois et al. (2014) predicted ΔP_{V} of up to 2% due to temporal changes in the mass decretion rate. The standard deviation of our P_{V} distribution (approximately Gaussian), namely, ~0.03%, is quite a bit lower than the above values. It is well explained in terms of the precision of our polarimetric data, as the typical error bar on P_{V} is ~0.01–0.02% (Fig. 17).
Interstellar parameters derived for o Aquarii: the Serkowski parameters, P_{max} and λ_{max}, with the polarization angle PA_{IS}.
8.4.4 A stable disk
Besides the analysis of the Hα EW and broadband polarimetric quantities, our modeling with the code HDUST indicates that the disk must be close to the steadystate regime: having a radial density law exponent of 3.0 (e.g. Haubois et al. 2012; Vieira et al. 2017). Other studies of o Aquarii are in fair agreement to our findings from HDUST. Using the radiative transfer code BEDISK, Silaj et al. (2010) derived m = 3.5 from the fit to the Hα profile, while Sigut et al. (2015) found m = 2.7 as a representative value from the analysis of all the different observables.
Previous and ongoing studies of Be stars with stable disks found similar results to ours. For example, Klement et al. (2015) found m = 2.9 for the latetype Be star β Canis Minoris (B8Ve). Correia Mota (2019) derived m = for α Arae (B2Vne). The B9Ve star α Columbae shows m = (A. Rubio, priv. comm.). Thus, the radial density exponent is consistently equal or somewhat less than 3.0 for these Be stars with stable disks.Also, from analysing the temporal variation of the disk density, Vieira et al. (2017) identified a slightly extended range of m (between ~3.0 and ~3.5) for the steadystate regime, in comparison to the canonical value of 3.5. As pointed out by these authors, this canonical value is based on simplifications of the standard theory, which assumes, for example, vertically isothermal disks and isolated systems (single stars). One possibility to explain the measured m lower than 3.5 could thus rely on nonisothermal effects in the disk structure (see, e.g., Carciofi & Bjorkman 2008).
Finally, we note that such longterm stability of o Aquarii’s disk is consistent with other results in the literature: latetype Be stars are more likely to have more stable disks than earlier Be stars (e.g., Vieira et al. 2017; LabadieBartz et al. 2018; Rímulo et al. 2018). As discussed in Sect. 8.3, the stellar rotation seems to be very close (~96%) to the critical value (391 ± 27 km s^{−1} from Frémat et al. 2005), in particular regarding the HDUST analysis: v_{rot} = 368 km s^{−1} (Table 5). This is consistent with the results from Cranmer (2005): Be stars with lower effective temperature T_{eff} ≲ 21 000 K – that is, later spectral types such as our target – are more likely to have a rotation rate close to one than the earlier Be stars. Thus, one possibility to explain such a longterm stability of the disk of o Aquarii could rely on its fast rotation, ensuring in this case a nearly constant massinjection rate into the disk.
9 Conclusions
We analyzed VEGA V^{2}, as well as VEGA and AMBER differential visibility and phase of the Beshell star o Aquarii. To date, the spectrointerferometric dataset analyzed in this paper is the largest for a Be star, considering quasicontemporaneous observations in both the Hα (VEGA) and Brγ (AMBER) lines.
For the first time, we measured o Aquarii’s stellar radius (R_{⋆} = 4.0 ± 0.3 R_{⊙}) and determined the disk extension in the Hα and Brγ lines as, respectively, 10.5 ± 0.3 D_{⋆} and 11.5 ± 0.1 D_{⋆}. Using radiative transfer models computed with the code HDUST, we explained the quasiidentical extension of the emission in these lines by an opacity effect found for disks seen under a high inclination angle.
We showed that the inclination angle derived from Hα is about 15° lower than the one determined in Brγ, when analysing each line separately with HDUST. More complex physical models, for example, with nonisothermal vertical scaling of the disk or the addition of heavier elements, could resolve this issue and should be investigated in the future.
Our simple kinematic model highlighted the high correlation between the rotational velocity at the base of the disk and the rotational law exponent β. Assuming external constraints, such as v sin i, we managed to constrain this parameter and showed that the disk rotation is nearly Keplerian (β ~ 0.43) from the analysis in the Brγ emission line. As for the inclination angle, the determination of β, using the Hα line (β ~ 0.30), seemsto be significantly biased. Other studies also verified such a large deviation from the Keplerian rotation for Be stars when analysing interferometric quantities measured in Hα (see, e.g., Delaa et al. 2011). One possible explanation would be the higher effect of noncoherent scattering on the Hα line formation than on Brγ.
Despite being derived purely from the fit to interferometric data, our bestfit HDUST model provides a very reasonable match to noninterferometric observables of o Aquarii: the observedSED, Hα and Brγ line profiles,and polarimetric quantities. Thus, this crosscheck provides an independent validation of our bestfit physical model. We found using HDUST a satisfying common physical description for the circumstellar disk in both Hα and Brγ: a base disk surface density Σ_{0} = 0.12 g cm^{−2} (ρ_{0} = 5.0 × 10^{−12} g cm^{−3}) and a radial density law exponent m = 3.0, that is, close to the steadystate regime according to the VDD model (m = 3.5). This result agrees with recent studies of other Be stars with stable disks, and may indicate the necessity to revise m = 3.5 (steadystate standing for single stars with vertically isothermal disks) that is predicted by the VDD theory. Otherwise, this could indicate nonisothermal effects on the disk vertical structure of o Aquarii. The longterm stability of the o Aquarii’s disk is verified by our analysis of a large sample of Hα profiles and polarimetric data, spanning about 20 and six years, respectively. Combined with older results in the literature, a longer global disk stability is suggested for up to at least 40 yr.
The stellar rotation seems to be very close (~96%) to the critical value (391 km s^{−1}), in particular accordingly to our HDUST analysis: v_{rot} = 368 km s^{−1} from the bestfit HDUST model with fixed M_{⋆} = 4.2 M_{⊙} (cf., Sects. 5.2 and 5.3). One possibility to explain such a longterm stability in the disk of o Aquarii could rely on its own high stellar rotation, being, in this case, a main source for the mass injection from the stellar surface to the disk. Thus, apart from the mass decretion due to other possible mechanisms in Be stars, this would provide a constant rate of mass injection. In short, our results on the stellar rotation and on the disk stability are consistent with the literature results showing that latetype Be stars are more likely to be fast rotators and have stable disks (see Sect. 8.4.4).
Finally, to further investigate these issues, our multiwavelength and multiemission line modeling approach must be performed on a larger sample of Be stars with disks of different densities and seen under different inclination angles. The implementation of a MCMC model fitting procedure with the kinematic model, and the use of our grid of HDUST models (BeAtlas), are very promising for the spectrointerferometric analysis of a large survey of Be stars, providing robust model parameters and associated uncertainties. A future project will attempt this task on a few dozen objects observed with VEGA and AMBER.
Acknowledgements
We thank the anonymous referee for helping to improve this paper. E.S.G.de.A. thanks OCA and the “Ville de Nice” (Nice, France) for the financial support to this work through the “Bourse Doctorale Olivier Chesneau” during the period of 2016–2019. E.S.G.de.A. acknowledges A. Rubio for relevant information about her work on α Columbae. R.L. has received funding from the European Union’s Horizon 2020 research and innovation program under the Marie SkłodowskaCurie grant agreement n. 664931. D.M.F acknowledges FAPESP (grant 2016/168441). A.C.C. acknowledges support from CNPq (grant 307594/20157). This work was supported by the “Programme National de Physique Stellaire” (PNPS) of CNRS/INSU cofunded by CEA and CNES. This work made use of the computing facilities of the Laboratory of Astroinformatics (IAG/USP, NAT/Unicsul), whose purchase was made possible by the Brazilian agency FAPESP (grant 2009/540064) and the INCTA. This work is based upon observationsobtained with the Georgia State University Center for High Angular Resolution Astronomy Array at Mount Wilson Observatory. The CHARA Array is supported by the National Science Foundation under Grants No. AST1211929 and AST1411654. This work used BeSOS Catalogue, operated by the Instituto de Física y Astronomía, Universidad de Valparaíso, Chile: http://besos.ifa.uv.cl and funded by Fondecyt iniciación N 11130702. Based on observations collected at the European Southern Observatory under ESO programmes 087.D0311 and 094.D0140. This work has made use of the BeSS database, operated at LESIA, Observatoire de Meudon, France: http://basebe.obspm.fr. Some of the data presented in this paper were obtained from the Mikulski Archive for Space Telescopes (MAST). STScI is operated by the Association of Universities for Research in Astronomy, Inc., under NASA contract NAS526555. Support for MAST for nonHST data is provided by the NASA Office of Space Science via grant NNX13AC07G and by other grants and contracts. This research has made use of the JeanMarie Mariotti Center (JMMC) services LITpro, SearchCal, and AMHRA, codevelopped by CRAL, IPAG and Lagrange. This work has made use of the SIMBAD and VizieR databases, operated at CDS, Strasbourg, France.
Appendix A Observational logs
Stellar calibrators used for the VEGA observations.
VEGA and AMBER observations.
Appendix B MCMC fitting tests: fits to the VEGA and AMBER data with the kinematic model
Fig. B.1
As in Fig. 4, but for the other MCMC fitting tests (test i in the left and test ii in the right) to fit the VEGA differential data. 
Fig. B.2
As in Fig. 4, but for the other MCMC fitting tests (test i in the left and test ii in the right) to fit the AMBER differential data. 
Appendix C Bestfit kinematic and HDUST models: AMBER
Fig. C.1
Comparison between our bestfit kinematic (dashed red; Table 2) and HDUST (dashed blue; Table 5) models to all the AMBER measurements (black). 
Appendix D Interstellar polarization
Fitted Serkowski parameters, with the polarization angle, for the field stars used to derive the interstellar polarization of o Aquarii.
Fig. D.1
Fitted P_{max} for the field stars (Table D.1, open triangles) as a function of the Gaia DR2 distance. From the linear fit to P_{max} vs. d for the field stars (dotted line), we determined P_{max} for o Aquarii (redcross). 
References
 Abrahamyan, H. V., Mickaelian, A. M., & Knyazyan, A. V. 2015, Astron. Comput., 10, 99 [NASA ADS] [CrossRef] [Google Scholar]
 Anderson, E., & Francis, C. 2012, Astron. Lett., 38, 331 [NASA ADS] [CrossRef] [Google Scholar]
 Andrillat, Y. 1983, A&AS, 53, 319 [NASA ADS] [Google Scholar]
 Arcos, C., Kanaan, S., Chávez, J., et al. 2018, MNRAS, 474, 5287 [Google Scholar]
 Auer, L. H., & Mihalas, D. 1968, ApJ, 153, 245 [NASA ADS] [CrossRef] [Google Scholar]
 Bednarski, D. 2016, Master thesis, Universidade de São Paulo, Brazil [Google Scholar]
 Bjorkman, J. E., & Carciofi, A. C. 2005, ASP Conf. Ser., 337, 75 [NASA ADS] [Google Scholar]
 Bonneau, D., Clausse, J.M., Delfosse, X., et al. 2006, A&A, 456, 789 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Bourges, L., Mella, G., Lafrasse, S., et al. 2017, VizieR Online Data Catalog: II/346 [Google Scholar]
 Carciofi, A. C. 2011, IAU Symp., 272, 325 [Google Scholar]
 Carciofi, A. C., & Bjorkman, J. E. 2006, ApJ, 639, 1081 [NASA ADS] [CrossRef] [Google Scholar]
 Carciofi, A. C., & Bjorkman, J. E. 2008, ApJ, 684, 1374 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Carciofi, A. C., Okazaki, A. T., Le Bouquin, J.B., et al. 2009, A&A, 504, 915 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Castelli, F., & Kurucz, R. L. 2004, ArXiv eprints [arXiv:astroph/0405087] [Google Scholar]
 Chelli, A., Utrera, O. H., & Duvert, G. 2009, A&A, 502, 705 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Chesneau, O., Lagadec, E., OtulakowskaHypka, M., et al. 2012, A&A, 545, A63 [Google Scholar]
 Cochetti, Y. R., Arcos, C., Kanaan, S., et al. 2019, A&A, 621, A123 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Correia Mota, B. 2019, PhD thesis, Universidade de São Paulo, Brazil [Google Scholar]
 Cranmer, S. R. 2005, ApJ, 634, 585 [NASA ADS] [CrossRef] [Google Scholar]
 Cyr, I. H., Jones, C. E., Panoglou, D., Carciofi, A. C., & Okazaki, A. T. 2017, MNRAS, 471, 596 [Google Scholar]
 Delaa, O., Stee, P., Meilland, A., et al. 2011, A&A, 529, A87 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Domiciano de Souza, A., Kervella, P., Moser Faes, D., et al. 2014, A&A, 569, A10 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Domiciano de Souza, A., Bouchaud, K., Rieutord, M., Espinosa Lara, F., & Putigny, B. 2018, A&A, 619, A167 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Ekström, S., Meynet, G., Maeder, A., & Barblan, F. 2008, A&A, 478, 467 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Faes, D. M. 2015, PhD thesis, IAGUniversidade de Sao Paulo, Brazil, LagrangeUniversité de Nice, France [Google Scholar]
 ForemanMackey, D., Hogg, D. W., Lang, D., & Goodman, J. 2013, PASP, 125, 306 [Google Scholar]
 Frémat, Y., Zorec, J., Hubert, A.M., & Floquet, M. 2005, A&A, 440, 305 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Gaia Collaboration (Brown, A. G. A., et al.) 2018, A&A, 616, A1 [Google Scholar]
 Gies, D. R., Bagnuolo, Jr. W. G., Baines, E. K., et al. 2007, ApJ, 654, 527 [Google Scholar]
 Goodman, J., & Weare, J. 2010, Comm. App. Math. Comp. Sci., 5, 65 [Google Scholar]
 Haubois, X., Carciofi, A. C., Rivinius, T., Okazaki, A. T., & Bjorkman, J. E. 2012, ApJ, 756, 156 [NASA ADS] [CrossRef] [Google Scholar]
 Haubois, X., Mota, B. C., Carciofi, A. C., et al. 2014, ApJ, 785, 12 [Google Scholar]
 Henden, A. A., Templeton, M., Terrell, D., et al. 2016, VizieR Online Data Catalog: II/336 [Google Scholar]
 Hummel, W., & Dachs, J. 1992, A&A, 262, L17 [NASA ADS] [Google Scholar]
 Kervella, P., & Domiciano de Souza, A. 2006, A&A, 453, 1059 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Klement, R., Carciofi, A. C., Rivinius, T., et al. 2015, A&A, 584, A85 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 LabadieBartz, J., Chojnowski, S. D., Whelan, D. G., et al. 2018, AJ, 155, 53 [Google Scholar]
 Lee, U., Osaki, Y., & Saio, H. 1991, MNRAS, 250, 432 [NASA ADS] [CrossRef] [Google Scholar]
 Magalhães, A. M., Benedetti, E., & Roland, E. H. 1984, PASP, 96, 383 [NASA ADS] [CrossRef] [Google Scholar]
 Magalhães, A. M., Rodrigues, C. V., Margoniner, V. E., Pereyra, A., & Heathcote, S. 1996, ASP Conf. Ser., 97, 118 [Google Scholar]
 McLean, I. S.,& Clarke, D. 1979, MNRAS, 186, 245 [Google Scholar]
 Meilland, A., Millour, F., Stee, P., et al. 2007a, A&A, 464, 73 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Meilland, A., Stee, P., Vannier, M., et al. 2007b, A&A, 464, 59 [Google Scholar]
 Meilland, A., Stee, P., Chesneau, O., & Jones, C. 2009, A&A, 505, 687 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Meilland, A., Delaa, O., Stee, P., et al. 2011, A&A, 532, A80 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Meilland, A., Millour, F., Kanaan, S., et al. 2012, A&A, 538, A110 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Monnier, J. D., Che, X., Zhao, M., et al. 2012, ApJ, 761, L3 [NASA ADS] [CrossRef] [Google Scholar]
 Mourard, D., Clausse, J. M., Marcotto, A., et al. 2009, A&A, 508, 1073 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Mourard, D., Challouf, M., Ligi, R., et al. 2012, Proc. SPIE, 8445, 84450K [CrossRef] [Google Scholar]
 Mourard, D., Monnier, J. D., Meilland, A., et al. 2015, A&A, 577, A51 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Neiner, C., de Batz, B., Cochard, F., et al. 2011, AJ, 142, 149 [Google Scholar]
 Okazaki, A. T. 2001, PASJ, 53, 119 [NASA ADS] [CrossRef] [Google Scholar]
 Okazaki, A. T., Bate, M. R., Ogilvie, G. I., & Pringle, J. E. 2002, MNRAS, 337, 967 [NASA ADS] [CrossRef] [Google Scholar]
 Oudmaijer, R. D., & Parr, A. M. 2010, MNRAS, 405, 2439 [NASA ADS] [Google Scholar]
 Panoglou, D., Carciofi, A. C., Vieira, R. G., et al. 2016, MNRAS, 461, 2616 [NASA ADS] [CrossRef] [Google Scholar]
 Panoglou, D., Borges Fernandes, M., Baade, D., et al. 2019, MNRAS, 486, 5139 [NASA ADS] [CrossRef] [Google Scholar]
 Petrov, R. G., Malbet, F., Weigelt, G., et al. 2007, A&A, 464, 1 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Poeckert, R., & Marlborough, J. M. 1978, ApJS, 38, 229 [Google Scholar]
 Rímulo, L. R., Carciofi, A. C., Vieira, R. G., et al. 2018, MNRAS, 476, 3555 [NASA ADS] [CrossRef] [Google Scholar]
 Rivinius, T., Štefl, S., & Baade, D. 2006, A&A, 459, 137 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Rivinius, T., Carciofi, A. C., & Martayan, C. 2013, A&ARv, 21, 69 [Google Scholar]
 SanchezBermudez, J., Alberdi, A., Barbá, R., et al. 2017, ApJ, 845, 57 [NASA ADS] [CrossRef] [Google Scholar]
 Serkowski, K., Mathewson, D. S., & Ford, V. L. 1975, ApJ, 196, 261 [NASA ADS] [CrossRef] [Google Scholar]
 Sigut, T. A. A., & Jones, C. E. 2007, ApJ, 668, 481 [Google Scholar]
 Sigut, T. A. A., McGill, M. A., & Jones, C. E. 2009, ApJ, 699, 1973 [NASA ADS] [CrossRef] [Google Scholar]
 Sigut, T. A. A., Tycner, C., Jansen, B., & Zavala, R. T. 2015, ApJ, 814, 159 [Google Scholar]
 Silaj, J., Jones, C. E., Tycner, C., Sigut, T. A. A., & Smith, A. D. 2010, ApJS, 187, 228 [NASA ADS] [CrossRef] [Google Scholar]
 Slettebak, A., & Reynolds, R. C. 1978, ApJS, 38, 205 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Stee, P., & Bittar, J. 2001, A&A, 367, 532 [Google Scholar]
 Stee, P., & Meilland, A. 2012, ASP Conf. Ser., 464, 167 [NASA ADS] [Google Scholar]
 Stee, P., Delaa, O., Monnier, J. D., et al. 2012, A&A, 545, A59 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 TallonBosc, I., Tallon, M., Thiébaut, E., et al. 2008, Proc. SPIE, 7013, 70131J [Google Scholar]
 Tatulli, E., Millour, F., Chelli, A., et al. 2007, A&A, 464, 29 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 ten Brummelaar, T. A., McAlister, H. A., Ridgway, S. T., et al. 2005, ApJ, 628, 453 [NASA ADS] [CrossRef] [Google Scholar]
 Touhami, Y., Gies, D., Coudé du Foresto, V., & Schaefer, G. 2009, AAS Meeting Abstracts, 213, 409.18 [Google Scholar]
 Touhami, Y., Gies, D. R., Schaefer, G. H., et al. 2013, ApJ, 768, 128 [Google Scholar]
 Townsend, R. H. D., Owocki, S. P., & Howarth, I. D. 2004, MNRAS, 350, 189 [Google Scholar]
 van Leeuwen, F. 2007, A&A, 474, 653 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Vanzi, L., Chacon, J., Helminiak, K. G., et al. 2012, MNRAS, 424, 2770 [NASA ADS] [CrossRef] [Google Scholar]
 Vieira, R. G., Carciofi, A. C., Bjorkman, J. E., et al. 2017, MNRAS, 464, 3071 [Google Scholar]
 Waters, L. B. F. M. 1986, A&A, 162, 121 [NASA ADS] [Google Scholar]
 Wilking, B. A., Lebofsky, M. J., & Rieke, G. H. 1982, AJ, 87, 695 [NASA ADS] [CrossRef] [Google Scholar]
 Yudin, R. V. 2001, A&A, 368, 912 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Zorec, J., Frémat, Y., Domiciano de Souza, A., et al. 2016, A&A, 595, A132 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
See VEGA group page at https://lagrange.oca.eu/fr/vega.
LITpro software is available at https://www.jmmc.fr/english/tools/dataanalysis/litpro/
Available at the JMMC service AMHRA: https://amhra.oca.eu/AMHRA/
As implemented in R: https://stat.ethz.ch/Rmanual/Rdevel/library/stats/html/loess.html
Public data available in the Barbara A. Mikulski Archive for Space Telescopes (MAST): https://archive.stsci.edu/iue/
Public data available at: http://besos.ifa.uv.cl
Public data available at: http://basebe.obspm.fr
All Tables
Results from the geometric modeling of the VEGA V^{2} data in the closeby continuum band (635–650 nm) and in the band centered on Hα (649–664 nm).
Bestfit kinematic models from test (iii) for our VEGA (Hα) and AMBER (Brγ) differential data.
First three columns: mean and standard deviation values for each HDUST parameter of the BeAtlas grid: from analysis (ii) (19 bestfit HDUST models), analysis (iii) (16 bestfit HDUST models), and analysis (iv) (17 bestfit HDUST models).
Parameters of our bestfit HDUST model in the BeAtlas grid to explain the joint analysis of our interferometric data: VEGA calibrated and differential data and AMBER differential data.
Interstellar parameters derived for o Aquarii: the Serkowski parameters, P_{max} and λ_{max}, with the polarization angle PA_{IS}.
Fitted Serkowski parameters, with the polarization angle, for the field stars used to derive the interstellar polarization of o Aquarii.
All Figures
Fig. 1
uv plan coverage obtained around Hα (0.656 μm) with CHARA/VEGA (green) and Brγ (2.166 μm) with VLTI/AMBER (red). 

In the text 
Fig. 2
VEGA V^{2} measurements in the closeby continuum band (top) and in the Hα band (bottom) are shown in red points. Our bestfit models consisting of one (solid line) and two (dashed line) uniform disks are overplotted in blue. See Table 1 and text for discussion. 

In the text 
Fig. 3
Top panels: uniform disk diameter derived from each individual VEGA V^{2} measurements (continuum band in the left and Hα band in the right) plotted as a function of the baseline position angle (PA). The red dotted line represents the bestfit diameter from modeling all the data in each band (θ = 0.28 mas in the continuum and θ = 0.36 mas in the Hα band). Bottom panels: corresponding normalized residuals. 

In the text 
Fig. 4
Histogram distributions and twobytwo correlations (after the burnin phase) for the free parameters of our bestfit kinematic models using MCMC for the VEGA (left panel) and AMBER (right panel) differential data. The median values are shown in solid red lines and the first and third quartiles in dashed red lines. The median and the first and third quartiles estimated for the parameters of our bestfit models (VEGA and AMBER) are presented inTable 2. In the correlation plots, darker points correspond to models with lower values of χ^{2}. See text for discussion. 

In the text 
Fig. 5
Comparison between our bestfit kinematic models (dashed red; Table 2) and two different VEGA (top panels) and AMBER (bottom panels) measurements (black line). Our bestfit HDUST model is also shown (dashed blue; Table 5; discussion in Sect. 6). δλ of the kinematic model and AMBER data is increased to 1.8 Å in order to compare them to the HDUST model (δλ fixed to 1.8 Å). 

In the text 
Fig. 6
maps of 40 000 kinematic models as a function of v_{rot} and β from the fit to VEGA (top panel) and AMBER (bottom panel) differential data. Only these two parameters were varied in a regular step in the intervals shown here. The other parameters are fixed (Table 2). Our results found from the MCMC analysis for v_{rot} and β are indicated with red crosses. In order to highlight the correlation between β and v_{rot}, the gray region corresponds to an arbitrary number of models, encompassing about the 5000 best models in both cases. The value of β = −0.5 (Keplerian disk) and our determination for v_{rot} are marked in dashed black line. Note the strong correlation between the stellar rotational velocity and the disk velocity law exponent in both the cases. Also, note that a Keplerian disk is found from modeling the AMBER data, but not from VEGA. 

In the text 
Fig. 7
Analysis (i): lowest value of reduced χ^{2} () for each value of R_{eq}∕R_{p} (left panel) and stellar mass (middle panel) from the HDUST fit to the VEGA V^{2} data (642.5 nm band). Local regression fits to , as a function of the parameter values, are shown in red line. In the right panel, the predicted visibility from our bestfit HDUST model (red points; Table 5) is compared to the VEGA V^{2} measurements in the continuum band (black points). The predicted visibility from the HDUST model with the highest mass in the BeAtlas grid (14.6 M_{⊙}, highest in the middle panel) is shown in blue points. 

In the text 
Fig. 8
Lowest value of for each value of disk majoraxis position angle from the HDUST fit to the VEGA (top panel, analysis (ii)) and AMBER (bottom panel, analysis (iii)) differential visibility and phase. Local regression fits of as a function of the disk PA are shown as a red line. 

In the text 
Fig. 9
Lowest value of for each value of stellar inclination angle, oblateness, base disk surface density, and disk density law exponent from the HDUST fit to the: VEGA differential data (top, analysis (ii)), AMBER differential data (middle, analysis (iii)), and all the interferometric data considered in this section (bottom panel, analysis (iv)). The stellar mass is fixed to 4.2 M_{⊙} and disk PA to 110°. Local regression fits to , as a function of the parameter values, are shown as a red line. The mean parameter values for the sets of best models (Table 4) are marked in dashed black line. Our bestfit BeAtlas model to fit all the interferometric data is shown in Table 5. See text for discussion. 

In the text 
Fig. 10
Intensity maps of our bestfit HDUST and kinematic models at different wavelengths around the Hα line (first two rows) and the Brγ line (last two rows). Flux/pixel is in arbitrary units with the same scale in Hα and Brγ. The image integrated in wavelength around each of these lines (Δλ = 2.7 nm around Hα and 3.9 nm around Brγ) are shown in the second column. 

In the text 
Fig. 11
Comparison between the observed o Aquarii and model SEDs from the ultraviolet to the farinfrared region. Flux unit is in erg cm^{−2} s^{−1} Å^{−1} and wavelength is shown in logarithmic scale. IUE/SWP and IUE/LWP spectra are shown in black line and photometric data in black points. Top panel: purely photospheric models (color lines) with variation in the stellar radius (no inclusion of geometrical oblateness): R_{⋆} = 3.2 R_{⊙} (orchid), 4.0 R_{⊙} (red), and 4.4 R_{⊙} (green). Bottom panel: photospheric model with 4.0 R_{⊙} (red) and our bestfit HDUST model from fitting all the interferometric data (dashed blue line; Table 5). Note that the UBVbands are better reproduced with R_{⋆} = 4.0–4.4 R_{⊙}. Our best HDUST model reproduces the observed IR excess due to the circumstellar disk well. 

In the text 
Fig. 12
Comparison between our bestfit kinematic models (dashed red; Table 2) and HDUST model (dashed blue, Table 5) in the Hα and Brγ line profiles. Mean observed line profiles of Hα (BeSOS) and Brγ (AMBER) are shown in black line. Our bestfit kinematic and HDUST models provide reasonable synthetic profiles to the observed ones in both Hα and Bγ. 

In the text 
Fig. 13
Majoraxis FWHM of Gaussian distribution (fitted from our bestfit HDUST model) as a function of the HDUST inclination angle. All the other HDUST parameters are fixed. Blue points correspond to the fit in Hα and red points in Brγ. The vertical dashed lines mark our values for inclination angle derived from the HDUST analysis, fitting the data in Hα (blue) and Brγ (red). Note that the equivalent Gaussian fits show a similar extension (2.45 mas, marked in horizontal dashed line) for these values of i. 

In the text 
Fig. 14
Bias effect of the disk velocity on Hα modeling. One VEGA measurement and observed Hα profile (BeSOS, as in Fig. 12) are shown in black lines. Our bestfit kinematic (dashed red) and HDUST (dashed blue) models are shown in Hα visibility and line profile. They are compared to HDUST models with a higher mass of 10.8 M_{⊙} with: Σ_{0} = 0.12 g cm^{−2} (dashed orchid) and Σ_{0} = 0.28 g cm^{−2} (dashed green). See text for discussion. 

In the text 
Fig. 15
Top panel: 70 observed Hα profiles of o Aquarii from the BeSS database, covering about 17 yr of observations (2001–2018). Bottom panel: Hα equivalent width as a function of the observation time (modified Julian date). Civil dates are indicated for a part of the measured EW. The mean EW (solid) within the standard deviation (dashed lines) is marked in black. Local regression fit of EW, as a function of time, is shown as a solid red line. The time interval covered by our interferometric observations (VEGA) is indicated with dashed red lines. See text for discussion. 

In the text 
Fig. 16
VEGA (top panels) and AMBER (bottom panels) differential visibilities extracted from observations at different epochs (black line). VEGA measurements span four years and the AMBER ones span three years. The observation date and the baseline length (projected onto the sky) are indicated in the top of each panel. Our bestfit kinematic models derived from the fit, in a separate way, to each dataset (VEGA and AMBER) are shown in dashed red line. Note that our bestfit kinematic models match well to the differential visibilities obtained at different epochs. 

In the text 
Fig. 17
Polarimetric quantities of o Aquarii, as a function of the observation time, spanning about six years. Top panel: observed Vband polarization (44 measurements). Bottom panel: ratio between the observed B and Rbands polarization (40 measurements). The mean values of these quantities are shown in dashed line. See text for discussion. 

In the text 
Fig. B.1
As in Fig. 4, but for the other MCMC fitting tests (test i in the left and test ii in the right) to fit the VEGA differential data. 

In the text 
Fig. B.2
As in Fig. 4, but for the other MCMC fitting tests (test i in the left and test ii in the right) to fit the AMBER differential data. 

In the text 
Fig. C.1
Comparison between our bestfit kinematic (dashed red; Table 2) and HDUST (dashed blue; Table 5) models to all the AMBER measurements (black). 

In the text 
Fig. D.1
Fitted P_{max} for the field stars (Table D.1, open triangles) as a function of the Gaia DR2 distance. From the linear fit to P_{max} vs. d for the field stars (dotted line), we determined P_{max} for o Aquarii (redcross). 

In the text 
Current usage metrics show cumulative count of Article Views (fulltext article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.
Data correspond to usage on the plateform after 2015. The current usage metrics is available 4896 hours after online publication and is updated daily on week days.
Initial download of the metrics may take a while.