# Dynamical mass of the Ophiuchus intermediate-mass stellar system S1 with
DYNAMO–VLBA

Jazmín Ordóñez-Toro Instituto de Radioastronomía y Astrofísica, Universidad
Nacional Autónoma de México, Apartado Postal 72-3, Morelia 58089, México
Departamento de Astronomía, Universidad de Guanajuato, Apartado Postal 144,
36000 Guanajuato, México  Sergio A. Dzib Max-Planck-Institut für
Radioastronomie, Auf dem Hügel 69, D-53121 Bonn, Germany  IRAM, 300 rue de la
piscine, F-38406 Saint Martin d’Hères, France  Laurent Loinard Instituto de
Radioastronomía y Astrofísica, Universidad Nacional Autónoma de México,
Apartado Postal 72-3, Morelia 58089, México  Gisela Ortiz-León Instituto
Nacional de Astrofísica, Óptica y Electrónica, Apartado Postal 51 y 216, 72000
Puebla, México  Max Planck Institut für Radioastronomie, Auf dem Hügel 69,
D-53121 Bonn, Germany  Marina A. Kounkel Department of Astronomy, University
of Michigan, 500 Church Street, Ann Arbor, MI 48105, USA  Josep M. Masqué
Departamento de Astronomía, Universidad de Guanajuato, Apartado Postal 144,
36000 Guanajuato, México  S.-N. X. Medina  Max Planck Institut für
Radioastronomie, Auf dem Hügel 69, D-53121 Bonn, Germany  Phillip A. B. Galli
Núcleo de Astrofísica Teórica, Universidade Cidade de Sao Paulo, R. Galvao
Bueno 868, Liberdade, 01506-000 Sao Paulo, SP, Brazil  Trent J. Dupuy
Institute for Astronomy, University of Edinburgh, Royal Observatory, Blackford
Hill, Edinburgh, EH9 3HJ, UK  Luis F. Rodríguez Instituto de Radioastronomía y
Astrofísica, Universidad Nacional Autónoma de México, Apartado Postal 72-3,
Morelia 58089, México  Luis H. Quiroga-Nuñez  Department of Aerospace, Physics
and Space Sciences, Florida Institute of Technology A: 150 W University Blvd,
Melbourne, 32901, FL, USA

###### Abstract

We report dynamical mass measurements of the individual stars in the most
luminous and massive stellar member of the nearby Ophiuchus star-forming
region, the young tight binary system S1. We combine 28 archival datasets with
seven recent, proprietary VLBA observations obtained as part of the Dynamical
Masses of Young Stellar Multiple Systems with the VLBA project (DYNAMO–VLBA),
to constrain the astrometric and orbital parameters of the system, and recover
high accuracy dynamical masses. The primary component, S1A, is found to have a
mass of 4.11$\pm$0.10 M⊙, significantly less than the typical value, $\sim$ 6
M⊙ previously reported in the literature. We show that the spectral energy
distribution of S1A can be reproduced by a reddened blackbody with a
temperature between roughly 14,000 K and 17,000 K. According to evolutionary
models, this temperature range corresponds to stellar masses between 4 M⊙ and
6 M⊙ so the SED is not a priori inconsistent with the dynamical mass of S1A.
The luminosity of S1 derived from SED-fitting, however, is only consistent
with models for stellar masses above 5 M⊙. Thus, we cannot reconcile the
evolutionary models with the dynamical mass measurement of S1A: the models
consistent with the location of S1A in the HR diagram correspond to masses at
least 25% higher than the dynamical mass. For the secondary component, S1B, a
mass of 0.831 $\pm$ 0.014 M⊙ is determined, consistent with a low-mass young
star. While the radio flux of S1A remains roughly constant throughout the
orbit, the flux of S1B is found to be higher near the apastron.

astrometry — stars:formation — stars:kinematics

††facilities: VLBA (Napier et al., 1994), Gaia (Gaia Collaboration et al.,
2016a)††software: AIPS (Greisen, 2003), Astropy444https://www.astropy.org/
(Astropy Collaboration et al., 2013), NumPy555https://www.numpy.org/ (van der
Walt et al., 2011), SciPy666https://www.scipy.org/ (Jones et al., 2001),
Matplotlib777https://matplotlib.org/ (Hunter, 2007), and APLpy (Robitaille &
Bressert, 2012)

##  1 Introduction

The Ophiuchus molecular cloud complex is one of the nearest and best-studied
sites of active star formation (see review by Wilking et al. 2008a). It has,
indeed, played a crucial role in our overall understanding of star formation.
Most of the ongoing star formation activity occurs in the Ophiuchus core (also
known as Lynds 1688; L 1688), located at 137.3 $\pm$ 1.2 pc (Ortiz-León et
al., 2017, OLK17 hereafter). The most luminous and massive member of L 1688,
and of the entire Ophiuchus region, is the young Herbig-Be star Oph-S1 (also
known as S1). Early radio observations of S1 showed that its radio emission
originates mainly from two components: a thermal extended ($\sim
20^{\prime\prime}$) component associated to a compact H ii region and a non-
thermal compact source embedded within the H ii region (Andre et al., 1988).
Very Long Baseline Interferometry (VLBI) observations later showed that the
non-thermal component is associated with a magnetically active star; S1 was
the first young stellar object directly detected with the VLBI technique
(Andre et al., 1991). These radio observations, as well as early infrared
measurements by Lada & Wilking (1984), suggested that S1 was a young
intermediate-mass star of spectral type B3 to B5, with a mass of about 6 M⊙.

While the properties of the compact radio emission in S1 are consistent with a
magnetically active star, the detection of non-thermal radio emission in an
intermediate-mass star is a priori unexpected. This is because the interiors
of young intermediate-mass stars (3–8 M⊙) are thought to be fully radiative
while they move toward the main-sequence. Under these circumstances, the
dynamo effect, which requires a rotating core surrounded by a convective
envelope, cannot operate (Stelzer et al., 2005, 2006; Hubrig et al., 2009).
Observational evidence, however, has shown that a significant fraction of
Herbig Ae/Be stars and main sequence A/B stars show signs of magnetic activity
(Hubrig et al., 2009; Wade et al., 2011; Hubrig & Schöller, 2021). A plausible
explanation is that the magnetic activity is associated with a low-mass
companion star rather than the intermediate-mass star itself (e.g., Stelzer et
al., 2005, 2006). In the specific case of S1, lunar occultation observations
have indeed revealed the existence of a lower mass companion with a separation
around 20 mas (Richichi et al., 1994; Simon et al., 1995). More recently, Very
Long Baseline Array (VLBA) high angular resolution observations obtained as
part of the Gould’s Belt Distances Survey  (GOBELINS; OLK17, ) detected two
stellar components in S1. From these data, OLK17 obtained preliminary
dynamical masses of $5.78\pm 0.15$ and $1.18\pm 0.10$ M⊙ for the two stars in
S1. While these VLBA observations confirmed the existence of a low-mass
companion, they also show that both the low-mass companion and intermediate-
mass star emit non-thermal radio emission. Thus, a radio emission mechanism
immanent to the intermediate star remains necessary. Possible mechanisms
include the existence of a fossil magnetic field, and the temporary apparition
of a thin convective layer on the surface of the star (Palla & Stahler, 1993;
Tout & Pringle, 1995; Skinner et al., 2004).

Regardless of the exact emission mechanism, the binarity of S1 provides us
with a rare opportunity to dynamically measure the mass of an Herbig Be star.
In this paper, we aim at improving on the mass determination of OLK17. We
present new VLBA observations of S1, obtained as part of the Dynamical Masses
of Young Stellar Multiple Systems with the VLBA (DYNAMO -
VLBA)111https://www3.mpifr-bonn.mpg.de/div/radiobinaries/intro.php project. To
improve the constraints on the mass of the stars in S1, we combine these new
observations with previous VLBA data, particularly those obtained as part of
GOBELINS. In the rest of the paper, we will call the primary of the system S1A
(this is the young Herbig Be star; OLK17 and Loinard et al. 2008, LTM08
hereafter) and the secondary component S1B.

##  2 Observations and data reduction

The DYNAMO-VLBA observations (VLBA project code: BD215) used the same
correlator setup as GOBELINS to record the radio continuum flux at $\lambda=$
6.0 cm ($\nu=$ 5 GHz) using the VLBA of the National Radio Astronomy
Observatory (NRAO). S1 was observed six times, approximately every four
months, from early 2018 to late 2019. During the sixth observation, the
antennas on all sites were not pointing correctly for the first half of the
observation time, and an additional observation was obtained four days later.
Thus, S1 was observed a total of seven times as part of the DYNAMO-VLBA
project. The observational details and the measured position are provided in
Table 1.

The total duration of each observation was three hours. The first and last 25
minutes were dedicated to observing about a dozen quasars distributed over the
entire sky. The rest of the observation was spent on the target source and its
gain calibrator (J1627-2427) in cycles of two minutes on-target and one minute
on the calibrator. This specific gain calibrator was chosen because it is
located in the core of L 1688 at an angular distance of $10^{\prime}$ from S1
(Dzib et al., 2013). For all observations, the julian date reported in Table 1
corresponds to the mid-point of the three-hour observing run.

As mentioned earlier, S1 was observed with the VLBA prior to the DYNAMO-VLBA
project. Multi-epoch, phased-referenced observations were obtained as part of
VLBA projects with codes BL128 (six epochs), BT093 (eight epochs) as well as
the GOBELINS survey with code BL175 (14 epochs). These observations were first
reported in LTM08, OLK17, and Ortiz-León et al. (2018) and we include them in
our analysis. Thus, this paper considers a total of 35 VLBA observations of
S1, distributed between 2005 to 2019.

The initial data calibration followed standard procedures for phase-referenced
VLBI observations and was carried out with the Astronomical Image Processing
System (AIPS) software (Greisen, 2003). The steps have been described in
detail in Loinard et al. (2007), Dzib et al. (2010) and OLK17.

The phase calibrator used in the DYNAMO-VLBA and GOBELINS observations is
located only $10^{\prime}$ from the target source, but its flux density is
only about $50$ mJy. On the other hand, the observations of BL128 and BT093
used the stronger ($\sim 1.0\,$ Jy) quasar J1625–2527 located $\sim 1\hbox
to0.0pt{.\hss}^{\circ}1$ south S1. As a consequence, larger residual
calibration errors are expected to affect the final images obtained from BL128
and BT093. S1A is a compact source with a flux density of about 6.0 mJy (; )
so, given the expected noise levels in our observations, it ought to be
detected on individual baselines with a signal-to-noise ratio of order 5 on
integrations of two minutes. Self-calibration with a solution interval as
short as two minutes can therefore be applied (see Cornwell & Fomalont, 1989,
for further details).

The self-calibration procedure used the following steps. First, we imaged the
externally calibrated visibilities with a few clean components so the image
contained only the strongest peak. Second, phase corrections on time intervals
of 15 minutes are obtained using this image as a model. Then, the first and
second steps are repeated two more times, decreasing the solution time
interval to eight and two minutes, respectively. Final corrections in
amplitude and phase are calculated, with time intervals of fifteen and eight
minutes, using the final image of the previous step as the model. The final
image is obtained using these corrections.

After self-calibration, the data were imaged using a pixel size of 100
$\,\mu$as and a natural weighting scheme (robust = 5 in AIPS). The rms noise
level in the images (25–80 $\,\mu$Jy beam-1, see Table 1) is consistent with
the values expected theoretically. The fluxes and positions of the detected
sources were measured using a two-dimensional Gaussian fitting procedure (task
JMFIT in AIPS); the full list is presented in Table 1. Figure 1 shows the
images of S1 as seen in the seven new epochs from the DYNAMO-VLBA project.

For the observations of projects BL128 and BT093 corresponding to the 14 older
epochs, the measured positions were corrected to account for the most recent
measured position of the gain calibrator used for these observations,
J1625–2527. The positional offset was estimated by OLK17 to be
$\Delta\alpha=-23\,\mu$as and $\Delta\delta=-1.2\,$mas. These values were
added to the measured positions of S1A and S1B in all epochs corresponding to
projects BL128 and BT093, and are already taken into account in the values
listed in Table 1.

Table 1: S1A and S1B measured positions and flux densities from the 35 VLBA observations. |  |  | S1A | S1B |
---|---|---|---|---|---
Project | Date UT | Date |  $\alpha$(J2000.0) |  $\delta$ (J2000.0) | $S_{\nu}$ |  $\alpha$(J2000.0) |  $\delta$ (J2000.0) |  $S_{\nu}$aaUpper-limits are three times the noise level in the image.  | $\sigma_{\rm noise}$
name | (yyyy.mm.dd) | Julian Day |  $16^{h}26^{m}$ [s] |  $-24^{\circ}$23′ [′′] | (mJy) |  $16^{h}26^{m}$ [s] |  $-24^{\circ}$23′ [′′] | (mJy) |  ($\mu$Jy bm-1)
BL128 GA | 2005.06.24 | 2453545.731 | 34.17392792(28) | 28.428288(10) | $11.05\pm 0.14$ | 34.1727518(144) | 28.447731(383) | $~{}~{}0.97\pm 0.27$ | 62
BL128 GB | 2005.09.15 | 2453628.504 | 34.17367442(61) | 28.433862(21) | $~{}~{}5.86\pm 0.16$ |  $\cdots$ |  $\cdots$ | $<0.23$ | 75
BL128 GC | 2005.12.17 | 2453722.248 | 34.17434762(96) | 28.442703(36) | $~{}~{}5.94\pm 0.18$ |  $\cdots$ |  $\cdots$ | $<0.17$ | 57
BL128 GD | 2006.03.15 | 2453810.007 | 34.17463122(73) | 28.452589(29) | $~{}~{}5.92\pm 0.18$ |  $\cdots$ |  $\cdots$ | $<0.20$ | 68
BL128 GE | 2006.06.03 | 2453889.789 | 34.17400238(61) | 28.457124(20) | $~{}~{}5.77\pm 0.17$ | 34.17268902(578) | 28.453578(178) | $~{}~{}1.01\pm 0.21$ | 75
BL128 GF | 2006.08.22 | 2453969.571 | 34.17328928(51) | 28.463940(18) | $~{}~{}6.90\pm 0.15$ |  $\cdots$ |  $\cdots$ | $<0.23$ | 75
BT093 BA | 2007.06.05 | 2454256.784 | 34.17387398(51) | 28.480866(18) | $~{}~{}6.42\pm 0.13$ | 34.17222738(509) | 28.498965(275) | $~{}~{}0.62\pm 0.14$ | 70
BT093 BB | 2007.06.09 | 2454260.773 | 34.17402998(67) | 28.478336(24) | $~{}~{}6.15\pm 0.18$ |  $\cdots$ |  $\cdots$ | $<0.23$ | 75
BT093 BC | 2007.06.13 | 2454264.762 | 34.17380163(30) | 28.481044(12) | $~{}~{}9.39\pm 0.13$ |  $\cdots$ |  $\cdots$ | $<0.20$ | 67
BT093 BD | 2007.06.17 | 2454268.751 | 34.17386498(54) | 28.480123(17) | $~{}~{}7.38\pm 0.15$ |  $\cdots$ |  $\cdots$ | $<0.21$ | 69
BT093 BE | 2007.06.21 | 2454272.741 | 34.17374613(40) | 28.481819(15) | $~{}~{}8.30\pm 0.15$ | 34.17203044(111) | 28.499143(43) | $~{}~{}2.50\pm 0.08$ | 79
BT093 BF | 2007.06.25 | 2454276.730 | 34.17377508(36) | 28.483839(13) | $10.36\pm 0.15$ | 34.17205001(165) | 28.500992(66) | $<0.21$ | 70
BT093 BG | 2007.06.29 | 2454280.718 | 34.17368822(56) | 28.482562(22) | $~{}~{}7.88\pm 0.15$ | 34.17193311(538) | 28.49900(22) | $~{}~{}1.02\pm 0.17$ | 74
BT093 BH | 2007.07.03 | 2454284.708 | 34.17364237(47) | 28.482137(17) | $~{}~{}7.82\pm 0.15$ |  $\cdots$ |  $\cdots$ | $<0.23$ | 75
BL175 BN | 2012.12.16 | 2456278.239 | 34.17329746(93) | 28.629756(29) | $~{}~{}1.80\pm 0.07$ |  $\cdots$ |  $\cdots$ | $<0.08$ | 25
BL175 E1 | 2013.09.02 | 2456537.533 | 34.17210801(57) | 28.648987(19) | $~{}~{}3.38\pm 0.06$ |  $\cdots$ |  $\cdots$ | $<0.08$ | 27
BL175 G0 | 2014.03.01 | 2456718.038 | 34.17336345(48) | 28.661664(17) | $~{}~{}3.13\pm 0.05$ |  $\cdots$ |  $\cdots$ | $<0.08$ | 28
BL175 CR | 2014.10.07 | 2456937.933 | 34.17244123(90) | 28.677675(30) | $~{}~{}3.41\pm 0.09$ |  $\cdots$ |  $\cdots$ | $<0.12$ | 41
BL175 EX | 2015.02.27 | 2457081.043 | 34.17321421(58) | 28.692174(21) | $~{}~{}3.28\pm 0.06$ |  $\cdots$ |  $\cdots$ | $<0.11$ | 35
BL175 GW | 2015.10.04 | 2457300.443 | 34.17207617(62) | 28.702069(22) | $~{}~{}2.79\pm 0.05$ |  $\cdots$ |  $\cdots$ | $<0.08$ | 25
BL175 F5 | 2016.03.26 | 2457473.996 | 34.17309283(23) | 28.718105(8) | $~{}~{}6.44\pm 0.04$ | 34.17127710(467) | 28.733943(184) | $~{}~{}0.27\pm 0.02$ | 30
BL175 F8 | 2016.04.28 | 2457506.878 | 34.17291981(39) | 28.720696(13) | $~{}~{}4.41\pm 0.05$ | 34.17100551(704) | 28.733448(334) | $~{}~{}0.17\pm 0.03$ | 30
BL175 IJ | 2016.09.09 | 2457640.513 | 34.17202055(95) | 28.730246(30) | $~{}~{}4.56\pm 0.08$ |  $\cdots$ |  $\cdots$ | $<0.09$ | 30
BL175 IN | 2016.09.16 | 2457647.520 | 34.17201100(30) | 28.730546(11) | $~{}~{}5.65\pm 0.05$ | 34.17029175(616) | 28.732533(285) | $~{}~{}0.20\pm 0.03$ | 25
BL175 JE | 2017.03.18 | 2457831.019 | 34.17259913(35) | 28.745253(13) | $~{}~{}4.58\pm 0.05$ |  $\cdots$ |  $\cdots$ | $<0.09$ | 28
BL175 JF | 2017.03.25 | 2457838.063 | 34.17259863(25) | 28.745577(9) | $~{}~{}6.95\pm 0.05$ |  $\cdots$ |  $\cdots$ | $<0.08$ | 27
BL175 K8 | 2017.09.03 | 2457999.557 | 34.17175982(45) | 28.753704(17) | $~{}~{}3.95\pm 0.05$ |  $\cdots$ |  $\cdots$ | $<0.20$ | 66
BL175 KB | 2017.09.15 | 2458012.492 | 34.17183287(140) | 28.755294(74) | $~{}~{}5.58\pm 0.10$ |  $\cdots$ |  $\cdots$ | $<0.09$ | 30
BD215 B0 | 2018.02.15 | 2458165.043 | 34.17288885(44) | 28.769136(14) | $~{}~{}6.73\pm 0.05$ | 34.1709178(126) | 28.780844(399) | $~{}~{}0.20\pm 0.03$ | 30
BD215 B1 | 2018.06.27 | 2458296.683 | 34.17202556(52) | 28.779521(16) | $~{}~{}6.45\pm 0.05$ | 34.17036273(380) | 28.779616(22) | $~{}~{}0.73\pm 0.03$ | 26
BD215 B2 | 2018.10.30 | 2458422.339 | 34.17162025(30) | 28.788574(10) | $~{}~{}8.30\pm 0.08$ |  $\cdots$ |  $\cdots$ | $<0.11$ | 38
BD215 B3 | 2019.03.01 | 2458544.026 | 34.17249189(21) | 28.795945(7) | $~{}~{}8.85\pm 0.06$ |  $\cdots$ |  $\cdots$ | $<0.09$ | 30
BD215 B4 | 2019.07.08 | 2458672.653 | 34.17179886(47) | 28.803943(14) | $~{}~{}7.42\pm 0.05$ | 34.1703085(132) | 28.823218(465) | $~{}~{}0.19\pm 0.03$ | 26
BD215 B5 | 2019.11.10 | 2458798.313 | 34.17193330(32) | 28.813973(11) | $~{}~{}5.56\pm 0.05$ | 34.16996475(474) | 28.826151(186) | $~{}~{}0.38\pm 0.05$ | 29
BD215 B6 | 2019.11.14 | 2458802.299 | 34.17197439(31) | 28.814119(11) | $~{}~{}9.79\pm 0.09$ | 34.17000425(112) | 28.825751(39) | $~{}~{}2.75\pm 0.09$ | 48

##  3 Astrometric fitting procedure

The displacements of S1A and S1B on the celestial sphere are described by the
combination of their common trigonometric parallax $\pi$, the uniform proper
motion of their center of mass in right ascension and declination,
$\mu_{\alpha}$ and $\mu_{\delta}$, and their orbital motion. For the primary
component, the complete equations describing the motions are:

$\displaystyle\alpha(t)$ $\displaystyle=\alpha_{0}+\mu_{\alpha}t+\pi
f_{\alpha}(t)+a_{1}Q_{\alpha}(t),$ (1) $\displaystyle\delta(t)$
$\displaystyle=\delta_{0}+\mu_{\delta}t+\pi f_{\delta}(t)+a_{1}Q_{\delta}(t),$
(2)

where $\alpha_{0}$ and $\delta_{0}$ are the coordinates of the center of mass
of the system at the chosen reference epoch, $f_{\alpha}$ and $f_{\delta}$ are
the projections of the parallactic ellipse (Seidelmann, 1992), and $a_{1}$ is
the semi-major axis of the orbit of the primary around the center of mass.
$Q_{\alpha}(t)$ and $Q_{\delta}(t)$ represent the projections of the orbital
motions and are functions of the orbital elements:

$\displaystyle Q_{\alpha}(t)$
$\displaystyle=r\left[\cos{\left(\theta+\omega\right)}\sin{\Omega}-\sin{\left(\theta+\omega\right)}\cos{\Omega}\cos{i}\right]/\cos{\delta},$
(3) $\displaystyle Q_{\delta}(t)$
$\displaystyle=r\left[\sin{\left(\theta+\omega\right)}\sin{\Omega}\cos{i}+\cos{\left(\theta+\omega\right)}\cos{\Omega}\right],$
(4)

where $\theta$ is the true anomaly and $r$ is the radius from the center of
mass along a Keplerian orbit. The inclination angle is $i$, the angle of the
line of nodes $\Omega$, and the angle from the ascending node to periastron
$\omega$. For the secondary component, the semi major axis $a_{2}$ is used
instead of $a_{1}$ in equations (1) and (2) and $\theta$ is rotated
$180^{\circ}$ in equations (3) and (4) (see e.g. Kounkel et al. 2017 for
details).

To solve the previous equations we translated the procedure described by
Kounkel et al. (2017) from IDL to python222The full routine is available at
https://github.com/mkounkel/astrometric_binaries/blob/master/astrometry_binary_python3.ipynb.
This routine simultaneously fits equations (1) to (4) to the data by the
method of least squares using MPFIT. It is able to simultaneously process
absolute positions of both stars from VLBA astrometry, as well as relative
positions of the secondary with respect to the primary from optical data if
they exist. The free orbital parameters in this model are the period $P$, the
semimajor axis of primary $a_{1}$, the eccentricity $e$, the inclination $i$,
the angle of the line of nodes $\Omega$, the time of periastron passage $T$,
the angle from node to periastron $\omega$, and the mass ratio $m_{2}/m_{1}$.
The total mass of the system is derived from Kepler’s third law and the
distance – obtained from the very same fit through the trigonometric parallax.

This routine is sensitive to the initial guesses of several non-linear orbital
parameters, namely $P$, $e$, $\omega$, and $\Omega$; with poorly chosen
initial guesses the fitter may not converge, getting stuck in a local minimum
instead. The remaining parameters are not affected by this, and can always be
scaled. Thus, to fully probe the parameter space, we performed the fit 1000
different times with different initial guesses. The routine then iterates on
these guesses, enabling optimal convergence for several dozens of these sets
of iterations, which are then evaluated using $\chi^{2}$. These sets of
iterations with $\chi^{2}<1.1\chi_{\rm best}^{2}$ are selected (where
$\chi_{\rm best}^{2}$ is the minimum chi squared of the entire set. The
uncertainties are evaluated in two ways: from examining the intrinsic scatter
from all of the valid solutions, and from calculating a weighted average error
of all of the errors returned by MPFIT, weighted by $\chi^{2}$. Typically,
both methods produce comparable distributions.

Figure 1: Final radio images of S1 at 4.5 GHz corresponding to each epoch
observed during the DYNAMO-VLBA project. Intensity background images are
clipped to intensities between –0.1 to 0.6 mJy beam-1. Contour levels are at
–4, 4, 6, 10, 15, 30, and 60 times the noise levels of the images, listed in
Table 1. Images are centered on the measured position of S1A, the primary
component of the S1 binary system as listed in Table 1 and shown in the top-
left corner of each plot. The white ellipse at the bottom-left of each plot
represents the synthesized beam of the corresponding observed epoch.

##  4 Results

The S1 system was detected in all seven new VLBA observations (Figure 1). S1A
was detected in all these epochs, with a mean flux density of 7.59 mJy; S1B
was detected five times with a mean flux density of 0.85 mJy. In the ancillary
observations, thanks to the improvement obtained by applying self-calibration,
we have detected S1B in more epochs than previously reported. In addition, we
found that a previous detection of S1B (BL128 GE from 2006 June 3) was
spurious: the emission peak at that position in fact corresponded to a
prominent side-lobe. We do detect S1B in the self-calibrated image of the same
epoch but at a different position. Therefore, the measured S1B position on
2006 June 03 by OLK17 is discarded in our analysis and we use, instead, the
newly measured position from the self-calibrated image. S1A was detected at
all 35 epochs – a 100% detection rate was also achieved in the publications
reporting on older observations ( ; , Ortiz-León et al. 2018). On the other
hand, the companion, S1B, is detected in 14 epochs, of which only four
detections were reported by OLK17. All the measured positions and flux
densities for S1A and S1B are listed in Table 1.

###  4.1 Astrometry

The final measured positions of S1A and S1B were used as inputs to our
astrometric fitting procedure. The resulting parameters from our best fit are
listed in Table 2 and are plotted in Figures 2 (total sky motion) and 3 (orbit
motion). The total and individual masses are derived from Kepler’s third law,
taking into account the obtained values for the period, semi-major axis, and
the mass ratio from our model fit. As a result, we obtain masses for the
primary and secondary components of 4.11$\pm$0.10 M⊙ and 0.83$\pm 0.01$,
respectively. It is important to note that these mass estimates are somewhat
at odds with those reported by OLK17: 5.78$\pm$0.15 M⊙ for the primary
component and 1.18$\pm$0.10 M⊙ for the secondary component. This discrepancy
arises from two main factors. Firstly, we have significantly increased the
number of S1B detections, yielding more robust and accurate measurements.
Secondly, we have corrected the position of S1B at epoch 2006 June 3rd (JD
2453889.789). As mentioned earlier, the detection of S1B reported by OLK17
turned out to be spurious, and the new position is significantly different
from the previous one. It is worth mentioning that OLK17 based their fits on
only four detections of S1B, so the inclusion of a spurious position had an
especially large detrimental effect. The addition of seven new, high-quality
data, obtained as part of DYNAMO–VLBA was particularly important to constrain
the orbital fit. From our study, we conclude that S1A is indeed an
intermediate-mass star, but less massive than previously suggested.

We also fitted the data with the MCMC Orbitize! package (Blunt et al., 2017)
to assess the accuracy of the orbital parameters obtained from the MPFIT
fitting. This algorithm explores the multidimensional parameter space,
effectively sampling the posterior distribution of the orbital parameters. The
results obtained demonstrate remarkable consistency with the values derived
from the MPFIT fit, confirming the reliability of our orbital parameter
estimates. A detailed description of the MCMC analysis and its results is
presented in Appendix References.

Table 2: Parameters obtained for the best-fit astrometric model. | Parameter | Value | Units
---|---|---|---
Astrometric parameters |  |
| $\alpha_{2016.0,{\rm centre}}$ | 16:26:34.1723678(23) | hh:mm:ss
| $\delta_{2016.0,{\rm centre}}$ | –24:23:28.712919(22) | ${}^{\circ}:\,^{\prime}:\,^{\prime\prime}$
|  $\Delta\alpha_{\rm S1A}$aa$\Delta\alpha$ and $\Delta\delta$ values are the orbital offset relative to the mass center at the epoch 2016.0 of components S1A and S1B.  | $3.69\pm 0.04$ | mas
|  $\Delta\delta_{\rm S1A}$aa$\Delta\alpha$ and $\Delta\delta$ values are the orbital offset relative to the mass center at the epoch 2016.0 of components S1A and S1B.  | $3.25\pm 0.02$ | mas
|  $\Delta\alpha_{\rm S1B}$aa$\Delta\alpha$ and $\Delta\delta$ values are the orbital offset relative to the mass center at the epoch 2016.0 of components S1A and S1B.  | $-17.99\pm 0.04$ | mas
|  $\Delta\delta_{\rm S1B}$aa$\Delta\alpha$ and $\Delta\delta$ values are the orbital offset relative to the mass center at the epoch 2016.0 of components S1A and S1B.  | $-15.87\pm 0.02$ | mas
| $\mu_{\alpha}$ | $-2.536\pm 0.002$ |  mas yr-1
| $\mu_{\delta}$ | $-26.853\pm 0.001$ |  mas yr-1
| $\pi$ | $7.28\pm 0.02$ | mas
| d | $137.2\pm 0.4$ | pc
Orbital parameters |  |
| $a_{1}$ | $0.418\pm 0.002$ | AU
| $a_{2}$ | $2.043\pm 0.018$ | AU
| P | $1.734\pm 0.001$ | years
| $T_{0}$ | $2457163.69\pm 0.65$ | Julian date
| $e$ | $0.65\pm 0.01$ |
| $\Omega$ | $249.6\pm 6.3$ | degrees
| $i$ | $25.4\pm 1.9$ | degrees
| $\omega$ | $168.4\pm 6.4$ | degrees
Dynamical masses |  |
| Total mass | $4.942\pm 0.113$ |  M⊙
| Mass 1 | $4.112\pm 0.099$ |  M⊙
| Mass 2 | $0.831\pm 0.014$ |  M⊙

Figure 2: Measured positions of S1A (red dots) and S1B (blue dots) shown as
offsets from the position of S1A in the first detected epoch used in the
present paper (2005 June 24; see Table 1). The red and blue curves show the
best fits, described in the text, to the positions of S1A and S1B,
respectively.

Figure 3: Orbital solution of S1. Top-panel: Stellar orbits relative to the
central of mass (black cross). The orbits from best fit for S1A (red) and S1B
(blue), are accompanied by a subset of 100 orbit pairs (light red for S1A, and
grey for S1B) to show the range of solutions allowed by the errors in our fit.
The subset of orbits shown corresponds to the best 100 of 5,000 orbits
generated by randomly obtaining orbital parameters from normal distributions
with mean values equal to the values obtained from our best fit and with
standard deviations equal to their errors. Measured positions of S1A (red
squares) and S1B (blue circles) after correction of parallax and proper
motions are also shown. Position errors are on the order of 0.1 mas, barely
seen at the plot scale. Bottom-panel:  Stellar relative positions and orbital
fit model of S1. The blue dots indicate the relative positions of S1B with
respect to S1A, and the errorbars consider the position errors of both
components which are added in quadrature. The dashed black line traces the
line of nodes from the model, and the black cross indicates the position of
the primary. The measured positions of S1B obtained from infrared observations
reported in Richichi et al. (1994) and Cheetham et al. (2015) are also
plotted, as red triangles. These infrared positions are listed in Table 4.

###  4.2 Flux variations

Figure 4 shows the flux of S1A and S1B as a function of the orbital phase,
defined as the fractional part of the equation

$\varphi=\frac{T-T_{0}}{P}$

where $T_{0}$ is the time of periastron (JD 2457163.69), and P is the orbital
period – both taken from the fitting results given in Table 2.

Figure 4: Flux as a function of orbital phase for S1A (top) and S1B (bottom)
at 5 GHz (blue dots) and 8 GHz (red squares). For S1A the bins (gray boxes)
show the mean and standard deviation values of fluxes in widths of 0.2 in the
orbital phase. For S1B the green arrows indicate the flux $3\,\sigma_{\rm
rms}$ upper limits when the source is not detected. The error bars are shown
and sometimes are smaller than the marker size. Figure 5: Relative positions
and orbital fit model. The green stars indicate the expected positions of S1B
in undetected epochs. Orbital phases ($\phi$) are indicated with vertical
black lines and labeled accordingly.

The projects BL128 and BT093 were recorded at a frequency of 8 GHz, while the
GOBELINS project used 5 GHz, except for the first epoch which was observed at
8 GHz. As mentioned in Section 2, the DYNAMO–VLBA observations were also
recorded at 5 GHz. To examine flux variations, we have to account for the fact
that the data have been observed at two different frequencies. For S1A from
the individual measured fluxes, we calculated the mean flux densities at both
frequencies: $5.49\pm 0.45$ mJy at 5 GHz and $6.95\pm 0.57$ mJy at 8 GHz. This
is roughly consistent with a flat spectrum over these frequencies, as also
seen in previous VLA observations (Dzib et al., 2013). Then, the flux values
measured at 8 GHz were scaled to the 5 GHz values by adding the difference
between the mean value at both frequencies.

For S1A (top panel of Figure 4), there are no obvious variations with phase in
the individual flux measurements. To search further for evidence of systematic
flux variations, we averaged the data in bins of width 0.2 (in orbital phase).
The horizontal black lines and vertical gray bars in the top panel of Figure 4
indicate the mean and standard deviation in each bin, respectively. From this
analysis, we do not find any evidence for variation with orbital phase.
Variability at the level of about 50% is present in the data, but it appears
to be stochastic, with no relation to the orbital phase.

For S1B, on the other hand, detections (particularly associated with higher
fluxes) appear to be clustered around an orbital phase value of 0.5 (bottom
panel of Figure 4). More quantitatively, S1B was detected in 9 out of 14
observations (64%) carried out in the range of orbital phases between 0.4 and
0.6, while it was detected in only 5 out of 21 observations (24%) at orbital
phases outside of this range. The lack of detections near periastron is
illustrated further in figure 5 where detections are shown in blue and non-
detections in green. The orbital phases ($\phi$) are indicated by tickmarks
along the elliptical trajectory.

We can examine the significance of the detection rates in each orbital phase
range by modeling the detectability of S1B with a binomial distribution. In
this statistical model, if the events are independent and the probability of a
success in a given experiment is $p$, the probability of $k$ successes in $n$
experiments is given by:

$p_{k,n}={n\choose k}p^{k}(1-p)^{n-k}.$

S1B is detected in 14 out 35 experiments, so the overall probability of
detection (a success) is $p=0.4$. Given this probability and assuming that the
chance of detection is independent of orbital phase, there is only a 4% chance
that 9 detections would be obtained out of 14 observations (the results
obtained near apoastron), while there is only a 6% chance that only 5
detections would be obtained out of 21 observations (the results obtained near
periastron). We conclude that the higher detection fraction near apoastron
than periastron is very unlikely to be the result of chance alone and that the
source is intrinsically more likely to be detectable near apoastron.

##  5 Discussion

In this section, we mostly discuss the implications of our mass measurements
for the properties of the S1 system. It is worthwhile mentioning at the outset
that, to the best of our knowledge, no evidence for cicumstellar (or
circumbinary) disks has ever been reported for this system.

###  5.1 Spectral Energy Distribution of S1

Before dynamical measurements became possible, the mass of S1 had been
estimated from its photometry, and some spectroscopy, to be about 6 M⊙. Let us
briefly review these works. Based on early infrared observations, Grasdalen et
al. (1973) classified S1 as a late B star. It was later classified as B3V by
Elias (1978) due to the surrounding Hii region causing strong CO absorption in
the 2.3 $\mu$m band. Combining spectroscopy and photometry, Cohen & Kuhi
(1979) obtained a B2 classification. Andre et al. (1988) estimated a B3.5
spectral type on the basis of the radio flux of the surrounding Hii region.
From spectroscopic observations, Bouvier & Appenzeller (1992) classified it as
a B4 star since it exhibits strong Hα absorption but no lithium absorption;
this same classification is adopted in Andre et al. (1991). Analysis by Lada &
Wilking (1984) of the spectral energy distribution implied that S1A has a
luminosity of 1,500 $L_{\odot}$ and an effective temperature of 16,000 K
consistent with a main-sequence star of spectral type B3-B5. Finally,
spectroscopy by Wilking et al. (2005a) implied a spectral type earlier than B8
– they adopted a B3 classification based on the effective temperature from
Lada & Wilking (1984). The classification we just described implies a mass
between 5 and 6 M⊙ – this is significantly higher than the value we derive
dynamically.

Recently, Mookerjea et al. (2018) compiled photometric observations of S1A
available through the astronomical catalog service VizieR to plot the SED of
S1A and fit it with a reddened blackbody. They adopted a blackbody temperature
of 17,000 K, which corresponds to the effective temperature of a B3/B4 ZAMS
star, as suggested by the studies mentioned above. They assume a stellar
radius $R=4.35$ R⊙ which, combined with the temperature mentioned above,
results in a luminosity $L=1,420$ L⊙. Their best fit implies a reddening
$E(B-V)=3.8$ mag, a selective extinction ratio $R_{V}\equiv
A_{V}/E(B-V)=3.35$, and a total $V$ band extinction $A_{V}=12.7$ mag.

In order to examine the compatibility between our dynamical mass of S1A and
its observational properties, We performed an analysis similar to that of
Mookerjea et al. (2018). We included the photometric points currently
available in VizieR; which cover the optical and infrared bands from 0.4
$\mu$m to 8 $\mu$m. The complete list of included measurements is provided in
table LABEL:vizier of Appendix References. We adopted a standard selective
extinction ratio $R_{V}=3.1$, and the extinction curves of Fitzpatrick (1999)
for the range between 0.4 $\mu$m and 1.25 $\mu$m and Rieke et al. (1989) for
the range between 1.25 $\mu$m and 8 $\mu$m. Figure 6 shows that the data
points can be reproduced by a range of effective temperatures, from about
14,000 K to 17,000 K. The specific values T${}_{\text{eff}}=$ 14,000 K,
15,700, and 17,000 K shown in Figure 6 were chosen because they correspond to
spectral types B5, B4 and B3 according to Pecaut & Mamajek (2013). Values in
the lower part of this range are consistent with the dynamical mass of S1
determined here. The extinction derived for the models shown in Figure 6 vary
from 10.5 (for $T_{\text{eff}}=14,000$ K) to 11.8 ($T_{\text{eff}}=17,000$ K).
Plots of the residuals between the models and the photometric data are shown
in the bottom panels.

Figure 6: Spectral energy distribution of S1A. The data points (black squares)
are taken from VizieR (see Table LABEL:vizier) while the fits (color lines)
are based on the extinction curves of Rieke et al. (1989) and Fitzpatrick
(1999) [See text]. The figure also includes residual plots for each model
(bottom panels). Figure 7: Evolutionary tracks based on the PISA pre-main
sequence stellar evolution models (Tognelli et al., 2011). These models have a
metalicity $Z=0.02$, helium abundance $Y=0.2880$, mixing-length $\alpha=1.90$
and deuterium abundance $X_{D}=4\times 10^{-5}$. The blue, red, and green data
squares indicate, respectively, stars with $T_{\text{eff}}$ = 17,000 K and $L$
= 2,120 L⊙, $T_{\text{eff}}$ = 15,700 K and $L$ = 1,430 L⊙, and
$T_{\text{eff}}$ = 14,000 K and $L$ = 730 L⊙. These correspond to the three
SED models, of the same colors, in Figure 6. The grey square shows the locus
of a star with $T_{\text{eff}}$ = 14,000 K, $L$ = 300 L⊙, which would be
consistent with a 4 $M_{\odot}$ star shortly before reaching the ZAMS.

###  5.2 Stellar evolution models with no rotation for S1A

Dynamical masses obtained in binary systems are essential to constrain stellar
evolution models. Nowadays, this is particularly important for young stars
because pre-main sequence evolution models are less reliable and less mature
than their main sequence counterparts (e.g., Hillenbrand & White, 2004;
Stassun et al., 2014). Models for intermediate mass young stars are
particularly uncertain given the paucity of existing observational
constraints. It is, therefore, worthwhile to compare the observational
properties of S1A with the predictions of theoretical models for the mass
determined dynamically in the present paper.

In Figure 7, we show evolutionary tracks for masses between 4 and 6 M⊙, based
on the PISA pre-main sequence stellar evolution models presented in Tognelli
et al. (2011). The loci of S1A in the HR diagram for the three SED models
shown in Figure 6 are indicated in Figure 7 as color points. The luminosity
and effective temperature of S1A appear to be consistent with masses between
about 5 and 6 $M_{\odot}$, significantly higher than the measured dynamical
mass of 4.1 $M_{\odot}$. An effective temperature of 14,000 K could be
consistent with a mass of 4 $M_{\odot}$, but only provided that the luminosity
was about 300 $L_{\odot}$ – instead of the 730 $L_{\odot}$ obtained from
modelling SED. This could only be achieved if the radius of the star were of
order 3 $R_{\odot}$ – instead of 4.6 $R_{\odot}$ obtained from modelling SED.

For completeness, we have repeated the same analysis with the $Y^{2}$ (Yonsei-
Yale) stellar models (Yi et al., 2001), the Palla-Stahler models (Palla &
Stahler, 1999), the PARSEC models (Nguyen et al., 2022, see below), and the
YaPSI YALE-POTSDAM stellar models (Spada et al., 2017). We obtained similar
conclusions to those reached with the PISA models: the theoretical predictions
for a mass of 4 $M_{\odot}$ are inconsistent with the locus of S1A in the HR
diagram. Specifically, the luminosity of S1A is at least twice larger than
predicted by the models. Conversely, given the locus of S1A in the HR diagram,
the theoretical models would predict a mass for S1A larger than measured
dynamically, by 20 to 50%. Interestingly, Stassun et al. (2014) also reported
systematic discrepancies between masses predicted by evolutionary models and
measured in binary systems, but for stars with masses below 1 $M_{\odot}$. For
masses above 1 $M_{\odot}$, the model predictions and observed masses were
found to be in agreement at the level of about 10%. We note, also, that
previous dynamical mass measurements of intermediate-mass stars (e.g.,
Johnston et al., 2019) did not find a significant discrepancy with
evolutionary models.

To finish with this section, it is worth emphasizing that the location of S1A
above and to the right of the main sequence on the HR diagram implies that it
has not yet reached the main sequence. Since intermediate mass stars have a
very short pre-main sequence lifetime, this places stringent constraints on
the age of the S1 system. For instance, according to the evolutionary models
of Tognelli et al. (2011) shown in Figure 7, a 5 M⊙ star reaches the main
sequence in just about 1 Myr. The possible location of S1A along the track of
a 5 M⊙ star shown by the green square of Figure 7 would indicate an age of 0.7
Myr according to the same models. Of course, we should proceed cautiously with
obtaining age estimates from the models since we have just shown that the
models cannot reproduce the location of S1A on the HR diagram given its
measured dynamical mass. Nevertheless, it seems reasonable to conclude that
the age of S1A is likely between 0.5 and 1 Myr.

###  5.3 Stellar evolution models with rotation for S1A and magnetic activity

The evolutionary models considered in the preceding section do not include
rotation. To investigate the possibility that rotation might help resolve the
discrepancy between the measured dynamical mass of S1A and the evolutionary
tracks, we considered the PARSEC V2.0 evolutionary models presented by Nguyen
et al. (2022). These models include recent re-calibrations of mixing by
overshooting and rotation. They parametrize rotation via the $\omega$
parameter defined as the ratio between the stellar angular velocity $\Omega$
and the breakup angular velocity $\Omega_{c}$ ($\Omega_{c}$ is the angular
velocity for which the centrifugal force becomes equal to the effective
gravity at the equator).

Figure 8 shows the PARSEC pre-main sequence tracks for a 4 M⊙ star of Solar
metallicity, appropriate for S1A, for $\omega$ values between zero (no
rotation) and 0.99 (near breakup). As can be seen, the effect of rotation is
imperceptible until the star approaches the ZAMS. There, the net effect of
rotation is to make the star redder and slightly less luminous. Thus, rotation
does not help reconcile the models with the dynamical mass of S1A since the
models without rotation produced a luminosity that was already lower than
expected for the dynamical mass of S1A.

The PARSEC models account for rotation but do not include magnetic fields and
this may be an important missing element for the interpretation of S1A since
the very existence of non-thermal radio emission requires the presence of a
substantial superficial magnetic field (102–103 G; Andre et al. 1988). The
origin of this field is unclear. In low-mass stars, superficial magnetic
fields of this magnitude are expected from dynamo theory, but this requires
the presence of deep convective envelopes (Parker, 1955). In intermediate mass
stars such as S1A, however, the energy transport from the core to the surface
is radiative, including during the pre-main sequence phase (Kippenhahn &
Weigert, 1990). Several scenarios have been proposed to account for the
presence of strong superficial magnetic fields in intermediate-mass stars.
Perhaps the most accepted possibility is a fossil origin (Schleicher et al.,
2023). An alternative is the generation of magnetic field through the so-
called Tayler-Spruit dynamo mechanism in a rotating star (Nandal et al.,
2023). Regardless of the origin of the magnetic field, its effect on the early
evolution of intermediate-mass stars remains poorly investigated.

Figure 8: PARSEC tracks (Nguyen et al., 2022) for a 4 M⊙ Solar metallicity
star with different rotation velocities from $\omega=0$ to $\omega=0.99$. See
text for the definition of $\omega$.

###  5.4 The nature of S1B

Our mass measurement and analysis of the SED of S1A implies that S1B is a
young (age of order 0.5 to 1 Myr, assuming co-evality with S1A) Solar type
star ($M\sim 0.8$ M⊙). It is interesting to see if this is consistent with the
few existing resolved infrared observations of S1A. The fluxes used in Figure
6 and listed in Table 4 correspond to a K band apparent magnitude of about 6.3
– this is for the entire system, but we can also take this value for S1A alone
since S1B is much weaker. Richichi et al. (1994) find a magnitude difference
$\Delta K$ = 1.6 between S1A and S1B, while Cheetham et al. (2015) report
$\Delta K$ = 1.2 and 2.3 in their two separate observations. This results in
an apparent K-band magnitude for S1B between 7.5 and 8.6. The distance modulus
for S1 is $\mu=5.7$ and the $A_{V}$ extinction of 12.7 reported from our fit
to the SED of S1A corresponds to $A_{K}=1.4$. Thus, we can estimate an
absolute K band magnitude for S1B between 0.4 and 1.5.

According to the pre-main sequence models of Baraffe et al. (2015), a 0.8 M⊙
star with an age between 0.5 and 1 Myr should have an absolute K band
magnitude between 1.5 and 2. We conclude that, within the errors, the observed
K band magnitude of S1A is consistent with its dynamical mass. It is important
to mention, in addition, that our age and magnitude estimates are quite
uncertain. The measurements of Richichi et al. (1994) and Cheetham et al.
(2015) suggest significant variability (either intrinsic to S1B or caused by
changing foreground absorption). The age of S1B is based on the assumption of
co-evality of S1A and S1B and models for S1A that we have shown to be unable
to reproduce the actual location of S1A in the HR diagram (Section 5.2). The
extinction is assumed to be the same for S1A and S1B, but the presence of
circumstellar material probably results in different absorptions for the two
stars. Future dedicated infrared resolved observations ought to enable a more
detailed comparison between the models and the properties of S1B.

###  5.5 Comparison with Gaia Astrometry

The Gaia satellite is the foremost astrometric mission at optical wavelengths
(Gaia Collaboration et al., 2016a). Its most recent data releases (Gaia DR3;
Babusiaux et al. 2022 and Gaia EDR3; Gaia Collaboration et al. 2021; Lindegren
et al. 2021) were published in 2022 and 2021. The angular resolution of Gaia,
$0\hbox to0.0pt{.\hss}^{\prime\prime}18$ (Lindegren et al., 2021), is
insufficient to resolve the S1 system333It is estimated that the final Gaia
resolution, in the upcoming fifth data release, will be $\sim 0\hbox
to0.0pt{.\hss}^{\prime\prime}1$ (https://www.cosmos.esa.int/web/gaia/science-
performance). It will, therefore, remain insufficient to spatially resolve the
components in this system. since the two components have a maximum separation
of about 25 mas (see Figure 3). The Gaia astrometric solution for S1 is listed
in Table 3 (Gaia source ID: Gaia DR3 6049161961233192576) where, for easy
comparison, we also show the astrometric results derived from our VLBA
observations. The Gaia position does not coincide with the VLBA positions for
S1A and S1B, as shown graphically in Figure 9. It falls near the VLBA position
for S1A, but is slightly displaced from it in the direction of S1B (or,
equivalently, in the direction of the center of mass of the system). This is
not unexpected: S1A is about 5 times more massive than S1B, and therefore
presumably about at least 500 times more luminous. Thus, the photocenter of
the S1 system traced by Gaia should indeed be located near S1A.

Figure 9: VLBA positions at epoch 2016.0 (the reference epoch of the Gaia DR3
release) of S1A (yellow triangle), S1B (purple diamond), and the center of
mass of the S1 system (cyan circle). The pink square indicates the Gaia DR3
position and is taken as the origin of the coordinate system. A black cross is
shown associated with each position; the sizes of these black crosses are ten
times the formal position errors in each case.

It is interesting to consider the consequences of the binarity of S1 on the
Gaia results. The Gaia DR3 astrometric fitting assumes that S1 is a single
star that moves in a linear and uniform fashion. Such a fit does not
adequately describe the motion of S1 given the binary nature of the system.
The Renormalised Unit Weight Error (RUWE) is a parameter provided by Gaia DR3
and used to assess if the single-star model provides a good fit to the
astrometric observations (Lindegren et al., 2018). A RUWE value above 1.4 is
considered an indication that a single-star model does not properly fit the
astrometric data. In the case of S1, the RUWE $=$ 2.75, confirming that a
single-star model does not provide a satisfactory description of the data. The
VLBA, on the other hand, resolves the components of the S1 system, and the
astrometric fitting takes into account the orbital motion of the system. Thus,
for S1, the VLBA astrometric results are certainly superior to those of Gaia.
In addition, the VLBA covers a significantly longer time span (14.39 years
against 2.76 for Gaia), further improving the quality of the VLBA results
(e.g., Dzib et al., 2017). The superiority of the VLBA astrometry is clearly
seen in the accuracy of the parallax measurements (0.02 mas for the VLBA vs.
0.14 mas for Gaia). We also note that the two parallax values are barely
consistent at the 3$\sigma$ level.

Table 3: Comparison of DYNAMO-VLBA project astrometry with Gaia Astrometry Parameter | Gaia DR3. | VLBA |  $\Delta$(VLBA – Gaias EDR3)
---|---|---|---
$\alpha_{2016}$ | 246.642385343${}^{\circ}\pm 0.13\,$mas | 246.642384866${}^{\circ}\pm 0.04\,$mas |  $-1.71\pm 0.14$ mas
$\delta_{2016}$ |  $-24.391308440^{\circ}\pm 0.08\,$mas |  $-24.391309144^{\circ}\pm 0.02\,$mas |  $-2.53\pm 0.08$ mas
VLBA $\alpha$ offset (S1A – S1 center) |  | 3.69 mas |  $\quad 1.97\pm 0.14$ mas
VLBA $\delta$ offset (S1A – S1 center) |  | 3.25 mas |  $\quad 0.71\pm 0.08$ mas
VLBA $\alpha$ offset (S1B – S1 center) |  | –17.99 mas |  $-19.71\pm 0.14$ mas
VLBA $\delta$ offset (S1B – S1 center) |  | –15.87 mas |  $-18.40\pm 0.08$ mas
$\pi$ |  $6.82\,\pm 0.14\,$ mas |  $7.28\,\pm 0.02\,$ mas |
$\mu_{\alpha}$ |  $-3.15\,\pm 0.18\,$ mas yr-1 |  $-2.536\,\pm 0.003\,$ mas yr-1 |
$\mu_{\delta}$ |  $-26.68\,\pm 0.15\,$ mas yr-1 |  $-26.853\,\pm 0.001\,$ mas yr-1 |

###  5.6 Comparison with infrared astrometry

The companion of S1A has also been inferred from infrared observations. The
S1B component was first identified using the Lunar occultation (LO) method by
Richichi et al. (1994). Later, Cheetham et al. (2015) also detected the
companion using Sparse Aperture Masking (SAM). Both IR methods provide
measurements of the total angular separation and the position angle (P.A.)
between the components of the binary. The separations found in the IR
observations are between 20 and 30 mas, consistent with the values we have
measured at the VLBA (Figure 3). From our best fit, we determined the orbital
phase at the epochs of the IR observations, and the corresponding expected
separations and P.A. The results are listed in Table 4 and shown graphically
in Figure 3. We find that the IR observations were made close to apastron;
within the errors, the expected separations and P.A. are in reasonable
agreement with the IR results.

Table 4: Comparison of Infrared and VLBA Astrometry

Infrared | VLBA |
---|---|---
Obs. Date | JD | Method | Separation | P.A. | Separation | P.A. | Orbital
(dd.mm.yyyy) |  |  | (mas) | (∘) | (mas) | (∘) | phase
13.06.1992 | 2448786.5 | LO | $20\pm 1$ | 291a | $22.7\pm 3.2$ | $268.6\pm 7.2$ | 0.77
14.04.2012 | 2456031.5 | SAM | $21.75\pm 0.39$ | $205.11\pm 0.74$ | $21.0\pm 3.1$ | $208.3\pm 9.2$ | 0.21
09.03.2013 | 2456360.5 | SAM | $28.89\pm 3.47$ | $268.54\pm 2.13$ | $24.7\pm 3.4$ | $263.5\pm 7.4$ | 0.73

Notes: a As pointed by Richichi et al. (1994), we have used the correction of
180∘ to the scan angle reported by them.

###  5.7 Flux variations

In at least one young stellar system, V 773 Tau, variations of the non-thermal
radio flux synchronized with the orbital phase of the system have been well
documented (Massi et al., 2002; Torres et al., 2012). In that specific case,
the radio flux of both stars appears to increase substantially near the
periastron of the system. In contrast, the radio emission of S1A shows no
systematic variation with the orbital phase, while the flux of the low-mass
companion S1B increases near the apastron (Section 4.2; Figure 5).

The origin of this behavior of S1B is unclear. In the case of V 773 Tau, it is
plausible that the flux increase near periastron might be caused by increased
flaring due to interactions between the individual magnetospheres when the
stars approach each other (Massi et al., 2002). Since flaring is thought to be
initiated by magnetic reconnection events, it is indeed reasonable to expect
an increase in the number of events when the two magnetospheres interact. To
interpret the reverse effect observed here with the same class of model would
require an unlikely scenario where reconnection events are more frequent when
the magnetospheres are further apart. An alternative would be to invoke the
presence of an optically thick region surrounding S1A with a radius of order
20 mas. Each time S1B plunged into this hypothetical region, its emission
would be absorbed by the optically thick material explaining the lack of
detection when the sources are closer than about 20 mas. The main difficulty
with this scenario is to explain why this region does not absorb the emission
from S1A itself. This requires a somewhat ad hoc morphology, with an inner
hole, for the hypothetical absorbing structure –for instance a toroid seen
close to face-on.

##  6 Conclusions

We have used seven new and 28 archival VLBA observations to study the S1
system in Ophiuchus, and constrain the dynamical masses of its components. We
detected the primary, S1A, in all 35 epochs, and the secondary, S1B, in 14
epochs. The best astrometric fit gives dynamical masses of M${}_{\rm
S1A}=4.11\pm 0.10$ M⊙ and M${}_{\rm S1B}=0.83\pm 0.01$ M⊙. This is the most
accurate mass estimation of S1 to date. The mass of S1A appears to be
significantly smaller than previously suggested (of order 6 $M_{\odot}$).
Reanalyzing the spectral energy distribution of S1, we conclude that it is
compatible with effective temperatures in the range of 14,000 to 17,000 K. The
location of S1 in the HR diagram suggests that it has not yet reached its
zero-age main sequence, but we find that pre-main sequence evolutionary models
are unable to explain its locus in the HR diagram for a mass of 4.1 M⊙. For
the models to pass through the allowed range of luminosity and effective
temperature of S1A, a mass at least 20 to 50% higher would be needed. Thus,
evolutionary models appear to over-predict the mass of S1A by that amount.

A comparison between the VLBA radio astrometry and Gaia DR3 shows that the
Gaia position falls close to, but is not coincident with the position of S1A
measured by the VLBA. This is not unexpected since Gaia surely traces the
motion of the photocenter of the system which should be located close to the
primary. However, the fact that Gaia cannot resolve the system nor account for
its orbital motion results in a parallax measurement that is significantly
less accurate than, and barely compatible with, the VLBA result. Finally, we
analyzed the variations of the flux density of S1A and S1B as a function of
the orbital phase. We find that the flux of S1A remains stable over the entire
orbit while the flux of S1B appears to be higher when the system is close to
apastron. The origin of this behavior is unclear.

J.O., G.O.L., J.M.M. and L.F.R. acknowledge the financial support of CONACyT,
México. L.L. acknowledges the support of DGAPA PAPIIT grants IN108324 and
IN112820 and CONACyT-CF grant 263356. S.A.D. acknowledge the M2FINDERS project
from the European Research Council (ERC) under the European Union’s Horizon
2020 research and innovation programme (grant No 101018682). P.A.B.G.
acknowledges financial support from the São Paulo Research Foundation (FAPESP)
under grants 2020/12518-8 and 2021/11778-9. The Long Baseline Observatory is a
facility of the National Science Foundation operated under cooperative
agreement by Associated Universities, Inc. The National Radio Astronomy
Observatory is a facility of the National Science Foundation operated under
cooperative agreement by Associated Universities, Inc. This work has made use
of data from the European Space Agency (ESA) mission Gaia
(https://www.cosmos.esa.int/gaia), processed by the Gaia Data Processing and
Analysis Consortium (DPAC,
https://www.cosmos.esa.int/web/gaia/dpac/consortium). Funding for the DPAC has
been provided by national institutions, in particular the institutions
participating in the Gaia Multilateral Agreement. This document was prepared
using the collaborative tool Overleaf available at https://www.overleaf.com/.

## References

  * Andre et al. (1988) Andre, P., Montmerle, T., Feigelson, E. D., Stine, P. C., & Klein, K.-L. 1988, ApJ, 335, 940, doi: 10.1086/166979
  * Andre et al. (1991) Andre, P., Phillips, R. B., Lestrade, J.-F., & Klein, K.-L. 1991, ApJ, 376, 630, doi: 10.1086/170311
  * Astropy Collaboration et al. (2013) Astropy Collaboration, Robitaille, T. P., Tollerud, E. J., et al. 2013, A&A, 558, A33, doi: 10.1051/0004-6361/201322068
  * Babusiaux et al. (2022) Babusiaux, C., Fabricius, C., Khanna, S., et al. 2022, arXiv e-prints, arXiv:2206.05989.  https://arxiv.org/abs/2206.05989
  * Baraffe et al. (2015) Baraffe, I., Homeier, D., Allard, F., & Chabrier, G. 2015, A&A, 577, A42, doi: 10.1051/0004-6361/201425481
  * Blunt et al. (2017) Blunt, S., Nielsen, E. L., De Rosa, R. J., et al. 2017, AJ, 153, 229, doi: 10.3847/1538-3881/aa6930
  * Blunt et al. (2020) Blunt, S., Wang, J. J., Angelo, I., et al. 2020, AJ, 159, 89, doi: 10.3847/1538-3881/ab6663
  * Bouvier & Appenzeller (1992) Bouvier, J., & Appenzeller, I. 1992, A&AS, 92, 481
  * Chambers et al. (2016) Chambers, K. C., Magnier, E. A., Metcalfe, N., et al. 2016, arXiv e-prints, arXiv:1612.05560, doi: 10.48550/arXiv.1612.05560
  * Cheetham et al. (2015) Cheetham, A. C., Kraus, A. L., Ireland, M. J., et al. 2015, ApJ, 813, 83, doi: 10.1088/0004-637X/813/2/83
  * Cieza et al. (2007) Cieza, L., Padgett, D. L., Stapelfeldt, K. R., et al. 2007, ApJ, 667, 308, doi: 10.1086/520698
  * Cieza et al. (2009) Cieza, L. A., Padgett, D. L., Allen, L. E., et al. 2009, ApJ, 696, L84, doi: 10.1088/0004-637X/696/1/L84
  * Cohen & Kuhi (1979) Cohen, M., & Kuhi, L. V. 1979, ApJS, 41, 743, doi: 10.1086/190641
  * Connelley & Greene (2010) Connelley, M. S., & Greene, T. P. 2010, AJ, 140, 1214, doi: 10.1088/0004-6256/140/5/1214
  * Cornwell & Fomalont (1989) Cornwell, T., & Fomalont, E. B. 1989, in Astronomical Society of the Pacific Conference Series, Vol. 6, Synthesis Imaging in Radio Astronomy, ed. R. A. Perley, F. R. Schwab, & A. H. Bridle, 185
  * Cutri et al. (2003) Cutri, R. M., Skrutskie, M. F., van Dyk, S., et al. 2003, VizieR Online Data Catalog, II/246
  * Cutri et al. (2012) —. 2012, VizieR Online Data Catalog, II/281
  * Dzib et al. (2010) Dzib, S., Loinard, L., Mioduszewski, A. J., et al. 2010, ApJ, 718, 610, doi: 10.1088/0004-637X/718/2/610
  * Dzib et al. (2013) Dzib, S. A., Loinard, L., Mioduszewski, A. J., et al. 2013, ApJ, 775, 63, doi: 10.1088/0004-637X/775/1/63
  * Dzib et al. (2017) Dzib, S. A., Loinard, L., Rodríguez, L. F., et al. 2017, ApJ, 834, 139, doi: 10.3847/1538-4357/834/2/139
  * Elias (1978) Elias, J. H. 1978, ApJ, 224, 453, doi: 10.1086/156393
  * Evans et al. (2003) Evans, Neal J., I., Allen, L. E., Blake, G. A., et al. 2003, PASP, 115, 965, doi: 10.1086/376697
  * Fitzpatrick (1999) Fitzpatrick, E. L. 1999, PASP, 111, 63, doi: 10.1086/316293
  * Gaia Collaboration (2020) Gaia Collaboration. 2020, VizieR Online Data Catalog, I/350, doi: 10.26093/cds/vizier.1350
  * Gaia Collaboration et al. (2016a) Gaia Collaboration, Prusti, T., de Bruijne, J. H. J., et al. 2016a, A&A, 595, A1, doi: 10.1051/0004-6361/201629272
  * Gaia Collaboration et al. (2016b) —. 2016b, A&A, 595, A1, doi: 10.1051/0004-6361/201629272
  * Gaia Collaboration et al. (2018) Gaia Collaboration, Brown, A. G. A., Vallenari, A., et al. 2018, A&A, 616, A1, doi: 10.1051/0004-6361/201833051
  * Gaia Collaboration et al. (2021) —. 2021, A&A, 649, A1, doi: 10.1051/0004-6361/202039657
  * Girard et al. (2011) Girard, T. M., van Altena, W. F., Zacharias, N., et al. 2011, AJ, 142, 15, doi: 10.1088/0004-6256/142/1/15
  * Grasdalen et al. (1973) Grasdalen, G. L., Strom, K. M., & Strom, S. E. 1973, ApJ, 184, L53, doi: 10.1086/181287
  * Grasser et al. (2021) Grasser, N., Ratzenböck, S., Alves, J., et al. 2021, A&A, 652, A2, doi: 10.1051/0004-6361/202140438
  * Greisen (2003) Greisen, E. W. 2003, in Astrophysics and Space Science Library, Vol. 285, Information Handling in Astronomy - Historical Vistas, ed. A. Heck, 109, doi: 10.1007/0-306-48080-8_7
  * Heinze et al. (2018) Heinze, A. N., Tonry, J. L., Denneau, L., et al. 2018, AJ, 156, 241, doi: 10.3847/1538-3881/aae47f
  * Henden et al. (2015) Henden, A. A., Levine, S., Terrell, D., & Welch, D. L. 2015, in American Astronomical Society Meeting Abstracts, Vol. 225, American Astronomical Society Meeting Abstracts #225, 336.16
  * Hillenbrand & White (2004) Hillenbrand, L. A., & White, R. J. 2004, ApJ, 604, 741, doi: 10.1086/382021
  * Howell et al. (2014) Howell, S. B., Sobeck, C., Haas, M., et al. 2014, PASP, 126, 398, doi: 10.1086/676406
  * Hubrig & Schöller (2021) Hubrig, S., & Schöller, M. 2021, Magnetic Fields in O, B, and A Stars (IOP Publishing, Bristol/Tokyo), doi: 10.1088/2514-3433/abefcc
  * Hubrig et al. (2009) Hubrig, S., Stelzer, B., Schöller, M., et al. 2009, A&A, 502, 283, doi: 10.1051/0004-6361/200811533
  * Hunter (2007) Hunter, J. D. 2007, Computing in Science & Engineering, 9, 90, doi: 10.1109/MCSE.2007.55
  * Johnston et al. (2019) Johnston, C., Pavlovski, K., & Tkachenko, A. 2019, A&A, 628, A25, doi: 10.1051/0004-6361/201935235
  * Jones et al. (2001) Jones, E., Oliphant, T., Peterson, P., et al. 2001, SciPy: Open source scientific tools for Python.  http://www.scipy.org/
  * Kippenhahn & Weigert (1990) Kippenhahn, R., & Weigert, A. 1990, Stellar Structure and Evolution
  * Kounkel et al. (2017) Kounkel, M., Hartmann, L., Loinard, L., et al. 2017, ApJ, 834, 142, doi: 10.3847/1538-4357/834/2/142
  * Lada & Wilking (1984) Lada, C. J., & Wilking, B. A. 1984, ApJ, 287, 610, doi: 10.1086/162719
  * Lasker et al. (2008a) Lasker, B. M., Lattanzi, M. G., McLean, B. J., et al. 2008a, AJ, 136, 735, doi: 10.1088/0004-6256/136/2/735
  * Lasker et al. (2008b) —. 2008b, AJ, 136, 735, doi: 10.1088/0004-6256/136/2/735
  * Lawrence et al. (2007) Lawrence, A., Warren, S. J., Almaini, O., et al. 2007, MNRAS, 379, 1599, doi: 10.1111/j.1365-2966.2007.12040.x
  * Lindegren et al. (2018) Lindegren, L., Hernández, J., Bombrun, A., et al. 2018, A&A, 616, A2, doi: 10.1051/0004-6361/201832727
  * Lindegren et al. (2021) Lindegren, L., Klioner, S. A., Hernández, J., et al. 2021, A&A, 649, A2, doi: 10.1051/0004-6361/202039709
  * Loinard et al. (2008) Loinard, L., Torres, R. M., Mioduszewski, A. J., & Rodríguez, L. F. 2008, ApJ, 675, L29, doi: 10.1086/529548
  * Loinard et al. (2007) Loinard, L., Torres, R. M., Mioduszewski, A. J., et al. 2007, ApJ, 671, 546, doi: 10.1086/522493
  * Massi et al. (2002) Massi, M., Menten, K., & Neidhöfer, J. 2002, A&A, 382, 152, doi: 10.1051/0004-6361:20011585
  * Miret-Roig et al. (2022) Miret-Roig, N., Bouy, H., Raymond, S. N., et al. 2022, Nature Astronomy, 6, 89, doi: 10.1038/s41550-021-01513-x
  * Mookerjea et al. (2018) Mookerjea, B., Sandell, G., Vacca, W., Chambers, E., & Güsten, R. 2018, A&A, 616, A31, doi: 10.1051/0004-6361/201832953
  * Nandal et al. (2023) Nandal, D., Meynet, G., Ekström, S., et al. 2023, arXiv e-prints, arXiv:2312.13340, doi: 10.48550/arXiv.2312.13340
  * Napier et al. (1994) Napier, P. J., Bagri, D. S., Clark, B. G., et al. 1994, IEEE Proceedings, 82, 658, doi: 10.1109/5.284733
  * Nguyen et al. (2022) Nguyen, C. T., Costa, G., Girardi, L., et al. 2022, A&A, 665, A126, doi: 10.1051/0004-6361/202244166
  * Ortiz-León et al. (2017) Ortiz-León, G. N., Loinard, L., Kounkel, M. A., et al. 2017, The Astrophysical Journal, 834, 141, doi: 10.3847/1538-4357/834/2/141
  * Ortiz-León et al. (2018) Ortiz-León, G. N., Loinard, L., Dzib, S. A., et al. 2018, ApJ, 869, L33, doi: 10.3847/2041-8213/aaf6ad
  * Paegert et al. (2021) Paegert, M., Stassun, K. G., Collins, K. A., et al. 2021, arXiv e-prints, arXiv:2108.04778, doi: 10.48550/arXiv.2108.04778
  * Palla & Stahler (1993) Palla, F., & Stahler, S. W. 1993, ApJ, 418, 414, doi: 10.1086/173402
  * Palla & Stahler (1999) —. 1999, ApJ, 525, 772, doi: 10.1086/307928
  * Parker (1955) Parker, E. N. 1955, ApJ, 122, 293, doi: 10.1086/146087
  * Pecaut & Mamajek (2013) Pecaut, M. J., & Mamajek, E. E. 2013, ApJS, 208, 9, doi: 10.1088/0067-0049/208/1/9
  * Richichi et al. (1994) Richichi, A., Leinert, C., Jameson, R., & Zinnecker, H. 1994, A&A, 287, 145
  * Rieke et al. (1989) Rieke, G. H., Rieke, M. J., & Paul, A. E. 1989, ApJ, 336, 752, doi: 10.1086/167047
  * Robitaille & Bressert (2012) Robitaille, T., & Bressert, E. 2012, APLpy: Astronomical Plotting Library in Python.  http://ascl.net/1208.017
  * Schleicher et al. (2023) Schleicher, D. R. G., Hidalgo, J. P., & Galli, D. 2023, A&A, 678, A204, doi: 10.1051/0004-6361/202346809
  * Seidelmann (1992) Seidelmann, P. K. 1992, Explanatory Supplement to the Astronomical Almanac. A revision to the Explanatory Supplement to the Astronomical Ephemeris and the American Ephemeris and Nautical Almanac.
  * Simon et al. (1995) Simon, M., Ghez, A. M., Leinert, C., et al. 1995, ApJ, 443, 625, doi: 10.1086/175554
  * Skinner et al. (2004) Skinner, S. L., Güdel, M., Audard, M., & Smith, K. 2004, ApJ, 614, 221, doi: 10.1086/422708
  * Spada et al. (2017) Spada, F., Demarque, P., Kim, Y. C., Boyajian, T. S., & Brewer, J. M. 2017, ApJ, 838, 161, doi: 10.3847/1538-4357/aa661d
  * Spitzer Science Center & Infrared Science Archive (IRSA)(2021) (SSC) Spitzer Science Center (SSC), & Infrared Science Archive (IRSA). 2021, VizieR Online Data Catalog, II/368
  * Stassun et al. (2014) Stassun, K. G., Feiden, G. A., & Torres, G. 2014, New A Rev., 60, 1, doi: 10.1016/j.newar.2014.06.001
  * Stassun et al. (2019) Stassun, K. G., Oelkers, R. J., Paegert, M., et al. 2019, AJ, 158, 138, doi: 10.3847/1538-3881/ab3467
  * Stelzer et al. (2005) Stelzer, B., Flaccomio, E., Montmerle, T., et al. 2005, ApJS, 160, 557, doi: 10.1086/432375
  * Stelzer et al. (2006) Stelzer, B., Micela, G., Hamaguchi, K., & Schmitt, J. H. M. M. 2006, A&A, 457, 223, doi: 10.1051/0004-6361:20065006
  * Strom et al. (1995) Strom, K. M., Kepner, J., & Strom, S. E. 1995, ApJ, 438, 813, doi: 10.1086/175125
  * Tognelli et al. (2011) Tognelli, E., Prada Moroni, P. G., & DeglÍnnocenti, S. 2011, A&A, 533, A109, doi: 10.1051/0004-6361/200913913
  * Torres et al. (2012) Torres, R. M., Loinard, L., Mioduszewski, A. J., et al. 2012, ApJ, 747, 18, doi: 10.1088/0004-637X/747/1/18
  * Tout & Pringle (1995) Tout, C. A., & Pringle, J. E. 1995, MNRAS, 272, 528, doi: 10.1093/mnras/272.3.528
  * van der Walt et al. (2011) van der Walt, S., Colbert, S. C., & Varoquaux, G. 2011, Computing in Science Engineering, 13, 22, doi: 10.1109/MCSE.2011.37
  * Wade et al. (2011) Wade, G. A., Alecian, E., Grunhut, J., et al. 2011, in Astronomical Society of the Pacific Conference Series, Vol. 449, Astronomical Polarimetry 2008: Science from Small to Large Telescopes, ed. P. Bastien, N. Manset, D. P. Clemens, & N. St-Louis, 262.  https://arxiv.org/abs/0901.0347
  * Wilking et al. (2008a) Wilking, B. A., Gagné, M., & Allen, L. E. 2008a, in Handbook of Star Forming Regions, Volume II, ed. B. Reipurth, Vol. 5, 351
  * Wilking et al. (2008b) Wilking, B. A., Gagne, M., & Allen, L. E. 2008b, VizieR Online Data Catalog, 5, 351
  * Wilking et al. (2005a) Wilking, B. A., Meyer, M. R., Robinson, J. G., & Greene, T. P. 2005a, AJ, 130, 1733, doi: 10.1086/432758
  * Wilking et al. (2005b) —. 2005b, AJ, 130, 1733, doi: 10.1086/432758
  * Yershov (2014) Yershov, V. N. 2014, Ap&SS, 354, 97, doi: 10.1007/s10509-014-1944-5
  * Yi et al. (2001) Yi, S., Demarque, P., Kim, Y.-C., et al. 2001, ApJS, 136, 417, doi: 10.1086/321795

##  Appendix A Orbital Parameters Analysis of S1 using MCMC

As an additional step to rigorously quantify the key astrometric and orbital
parameters, we implemented an additional approach that combines the power of
two techniques: Orbits For The Impatient (OFTI) and Markov Chain Monte Carlo
(MCMC). Both algorithms are integrated into the Orbitize!
library888https://orbitize.info/en/latest/, which is a specialized, open-
source Python-based tool designed for characterizing orbits of binary
astronomical systems and operates within a Bayesian analysis framework (Blunt
et al., 2020). For our analysis, we exclusively utilized 14 observations in
which both components of the binary system were simultaneously detected,
adhering to the necessary criteria for orbit analysis using Orbitize! OFTI
works by generating a multitude of possible orbits based on the observed
astrometric data. It effectively explores the parameter space, laying the
groundwork for subsequent analysis with MCMC (Blunt et al., 2017). Based on
the results of the OFTI fit, we defined priors for the orbital elements,
including the semi-major axis ($a$), eccentricity ($e$), inclination ($i$),
argument of periastron ($\omega$), position angle of nodes ($\Omega$), and
epoch of periastron passage ($\tau$), which represents the fraction of the
orbital period that has elapsed from a given reference epoch to periastron,
within a range of [0, 1], and then proceeded to the MCMC phase. To enhance the
MCMC exploration, we employed a substantial ensemble of 10,000 walkers. These
walkers independently traversed the parameter space, ensuring comprehensive
sampling. To ensure the convergence of the MCMC chains, we executed a burn-in
phase with 500 steps before extracting meaningful results. The results are
presented in Table 5, and the corner plot in Figure 10. This corner plot
visually demonstrates the well-behaved nature of parameter sampling, further
reinforcing the reliability of our results compared to MPFIT. We note, in
particular that the total mass of the system derived from Orbitize!
(5.14$\pm$0.13 M⊙) is fully compatible with our MPIFIT result (4.94$\pm$0.11
M⊙).

Table 5: Results obtained through the Orbitize! method. | Parameter | Value | Units
---|---|---|---
Orbital Parameters |  |
| $a$ |  2.494$\pm$0.022 | AU
| P |  1.738$\pm 0.001$ | years
| $T_{0}$ |  2457164.29$\pm$0.02 | Julian date
| $e$ |  0.64$\pm$0.01 |
| $\Omega$ |  243.6$\pm$8.3 | degrees
| $i$ |  29.4$\pm$3.1 | degrees
| $\omega$ |  169.3$\pm$10.5 | degrees
| Total Mass |  5.137$\pm$0.134 | $M_{\odot}$

Figure 10: Corner plot generated from the MCMC analysis using the Orbitize!
library. The diagonal panels display 1D marginalized posterior distributions
for each of the orbital parameters, including $a$, $e$, $i$, $\omega$,
$\Omega$, $\tau$ and total mass. The off-diagonal panels illustrate 2D
covariances between these parameters.
Figure 11: Allowed orbital configurations for the binary system S1, derived
from the MCMC analysis conducted using the Orbitize! library. The black marker
(star) designates the position of the primary component, and the red points
represent radio observations of the secondary component. The upper right panel
depicts the evolution of the angular separation over time, while the lower
right panel illustrates the variation in the position angle throughout the
orbits.

##  Appendix B Photometric Data Table for SED Analysis

Table 6: VizieR S1 Photometric Data $\lambda$ ($\mu m$) |  Flux (Jy) |  Error flux (Jy) |  SED Filter |  ID Source
---|---|---|---|---
2.19 |  1.99 |  0.06 |  Johnson:K |  Strom et al. (1995)
2.19 |  1.94 |  0.03 |  Johnson:K |  Cutri et al. (2003)
1.63 |  1.3 |  0.04 |  Johnson:H |
1.25 |  0.444 |  0.004 |  Johnson:J |
2.16 |  1.98 |  0.03 |  2MASS:Ks |  Evans et al. (2003)
1.65 |  1.28 |  0.04 |  2MASS:H |
1.24 |  0.456 |  0.008 |  2MASS:J |
7.87 |  0.455 |  0.05 |  Spitzer/IRAC:8.0 |
5.73 |  0.741 |  0.063 |  Spitzer/IRAC:5.8 |
4.49 |  1.04 |  0.06 |  Spitzer/IRAC:4.5 |
3.55 |  1.34 |  0.07 |  Spitzer/IRAC:3.6 |
0.789 |  0.0239 |  $\cdots$ |  Cousins:I |  Wilking et al. (2005b)
1.25 |  0.191 |  0.003 |  UKIDSS:J |  Lawrence et al. (2007)
2.16 |  2 |  0.04 |  2MASS:Ks |  Cieza et al. (2007)
1.65 |  1.31 |  0.04 |  2MASS:H |
1.24 |  0.451 |  0.008 |  2MASS:J |
5.73 |  0.701 |  0.02 |  Spitzer/IRAC:5.8 |
4.49 |  0.969 |  0.018 |  Spitzer/IRAC:4.5 |
3.55 |  1.14 |  0.03 |  Spitzer/IRAC:3.6 |
0.784 |  0.0226 |  0.0083 |  POSS-II:i |  Lasker et al. (2008a)
0.64 |  0.00388 |  0.00146 |  POSS-II:F |
0.468 |  0.000156 |  5.9E-05 |  POSS-II:J |
2.16 |  2.04 |  0.04 |  2MASS:Ks |  Lasker et al. (2008b)
1.65 |  1.31 |  0.04 |  2MASS:H |
1.24 |  0.442 |  0.004 |  2MASS:J |
0.763 |  0.0125 |  0.0002 |  SDSS:i |
0.673 |  0.00668 |  1E-05 |  Gaia:G |
0.673 |  0.0269 |  0.0001 |  Gaia:G |
0.625 |  0.00236 |  3E-05 |  SDSS:r |
0.482 |  0.000245 |  5E-06 |  SDSS:g |
1.65 |  1.31 |  0.04 |  2MASS:H |  Wilking et al. (2008b)
1.24 |  0.451 |  0.008 |  2MASS:J |
5.73 |  0.455 |  0.05 |  Spitzer/IRAC:5.8 |  Cieza et al. (2009)
5.73 |  0.741 |  0.063 |  Spitzer/IRAC:5.8 |
4.49 |  1.04 |  0.06 |  Spitzer/IRAC:4.5 |
3.55 |  1.34 |  0.07 |  Spitzer/IRAC:3.6 |
2.15 |  2.1 |  $\cdots$ |  MKO:Ks |  Connelley & Greene (2010)
0.554 |  0.000582 |  $\cdots$ |  Johnson:V |  Girard et al. (2011)
2.16 |  2.04 |  0.05 |  2MASS:Ks |  Cutri et al. (2012)
1.65 |  1.27 |  0.07 |  2MASS:H |
1.24 |  0.155 |  $\cdots$ |  2MASS:J |
1.24 |  0.206 |  $\cdots$ |  2MASS:J |
1.24 |  0.442 |  0.008 |  2MASS:J |
1.24 |  0.223 |  $\cdots$ |  2MASS:J |
0.547 |  0.000656 |  $\cdots$ |  UVOT:V |  Yershov (2014)
0.439 |  0.000106 |  $\cdots$ |  UVOT:B |
0.763 |  0.0124 |  0.001 |  SDSS:i’ |  Henden et al. (2015)
0.96 |  0.0968 |  $\cdots$ |  PAN-STARRS/PS1:y |  Chambers et al. (2016)
0.865 |  0.0445 |  $\cdots$ |  PAN-STARRS/PS1:z |
0.748 |  0.0125 |  0.0002 |  PAN-STARRS/PS1:i |
0.613 |  0.00235 |  3E-05 |  PAN-STARRS/PS1:r |
0.477 |  0.000244 |  5E-06 |  PAN-STARRS/PS1:g |
0.673 |  0.00779 |  1E-05 |  Gaia:G |  Gaia Collaboration et al. (2016b)
0.748 |  0.00833 |  5E-05 |  PAN-STARRS/PS1:i |  Heinze et al. (2018)
0.613 |  0.00188 |  2E-05 |  PAN-STARRS/PS1:r |
0.772 |  0.0234 |  0.0001 |  GAIA/GAIA2:Grp |  Gaia Collaboration et al. (2018)
0.623 |  0.0073 |  4E-06 |  GAIA/GAIA2:G |
0.505 |  0.000654 |  3E-06 |  GAIA/GAIA2:Gbp |
2.16 |  2.01 |  0.03 |  2MASS:Ks |  Howell et al. (2014)
1.65 |  1.31 |  0.04 |  2MASS:H |
1.24 |  0.451 |  0.008 |  2MASS:J |
2.16 |  2.01 |  0.03 |  2MASS:Ks |  Stassun et al. (2019)
1.63 |  1.3 |  0.04 |  Johnson:H |
1.25 |  0.46 |  0.008 |  Johnson:J |
0.762 |  0.0226 |  0.0001 |  GAIA/GAIA3:Grp |  Gaia Collaboration (2020)
0.582 |  0.00732 |  2.00E-05 |  GAIA/GAIA3:G |
0.504 |  0.000647 |  3E-06 |  GAIA/GAIA3:Gbp |
5.73 |  0.548 |  1.6E-04 |  Spitzer/IRAC:5.8 |  Spitzer Science Center & Infrared Science Archive (IRSA)(2021) (SSC)
2.16 |  2.01 |  0.03 |  2MASS:Ks |  Paegert et al. (2021)
1.63 |  1.3 |  0.04 |  Johnson:H |
1.25 |  0.46 |  0.008 |  Johnson:J |
2.19 |  1.94 |  0.02 |  Johnson:K |  Grasser et al. (2021)
1.63 |  1.3 |  0.04 |  Johnson:H |
1.25 |  0.46 |  0.005 |  Johnson:J |
2.16 |  1.08 |  0.05 |  2MASS:Ks |  Miret-Roig et al. (2022)
1.65 |  0.644 |  0.03 |  2MASS:H |
1.24 |  0.38 |  0.017 |  2MASS:J |


Lead provides voice coordination and safety concerns to Wingman too early when the information is not yet relevant. |  This could be due to the lead having incorrect information or trying to coordinate prematurely. This could be unsafe because it could be distracting to both pilots and might prevent necessary information from getting through. [CCS9, CCS12] |  [RL1.5.41] "All flight personnel shall be briefed on communication protocol before takeoff. If the lead is sharing information too early, they shall stop and notify the other parties, and then repeat the information at the appropriate time. If the wingman tries to act on premature information, an effort shall be made to correct the mistake."
Lead provides voice coordination and safety concerns to Wingman too late when the event has already passed. |  This could be due to the lead being unaware of the event or a delay in communications. |  [RL1.5.42] "All flight personnel shall be briefed on communication protocol before takeoff. If the lead fails to share information on time, the wingman shall make a best effort to maintain safe flight if necessary and a contingency plan shall be followed"