Astrometrically registered maps of H2O and SiO masers toward VX Sagittarii

The supergiant VX Sagittarii is a strong emitter of both H2O and SiO masers. However, previous VLBI observations have been performed separately, which makes it difficult to spatially trace the outward transfer of the material consecutively. Here we present the astrometrically registered, simultaneous maps of 22.2 GHz H2O and 43.1/42.8/86.2/129.3 GHz SiO masers toward VX Sagittarii. The H2O masers detected above the dust-forming layers have an asymmetric distribution. The multi-transition SiO masers are nearly circular ring, suggesting spherically symmetric wind within a few stellar radii. These results provide the clear evidence that the asymmetry in the outflow is enhanced after the smaller molecular gas clump transform into the inhomogeneous dust layers. The 129.3 GHz maser arises from the outermost region compared to that of 43.1/42.8/86.2 GHz SiO masers. The ring size of the 129.3 GHz maser is maximized around the optical maximum, suggesting that radiative pumping is dominant.

V X Sagittarii (VX Sgr) is a red supergiant with a semiregular variable period of 732 days 1 . The distance measured from SiO maser proper motion is 1.57 kpc 2 . The photospheric diameter is 8.82 mas 3 (13.85 AU) at 2.0 μm. This star shows a heavy mass loss of about 2.5 × 10 −4 M ⊙ year −1 4 . It is well known that VX Sgr hosts strong OH, H 2 O, and SiO maser emitters, compact enough for very long baseline interferometry (VLBI) observations [5][6][7] . The SiO masers are located at 2-4 stellar radii from the stellar surface, inside the dust-formation layer, while the 22.2 GHz H 2 O maser is located outside the dust layer 8,9 , which is undergoing a radial acceleration. With the precise astrometrical registration of the SiO and H 2 O masers that are observed simultaneously, we can directly compare the properties of these masers on the scales of the individual maser gas clumps for tracing the mass transfer between these layers.
The  3 GHz SiO masers using the non-integer source frequency phase-referencing (SFPR) method [11][12][13] . This results provide observational evidence for a break in spherical symmetry between the SiO and H 2 O maser zone. The 129.3 GHz SiO maser from VX Sgr shows that radiative pumping is dominant, arising from the outermost region compared to the 43.1/ 42.8/86.2 GHz SiO masers. Figure 1 shows the astrometrically registered integrated intensity map of the H 2 O and SiO maser lines observed on March 27, 2016 (optical phase φ = 0.67) toward VX Sgr 14 . The distributions of the SiO masers show a typical ring-like structure, while that of the H 2 O maser shows an asymmetric structure spread slightly in the NW and SE direction in~350 mas and relatively dense distributions in the NE and SW direction in~270 mas. These results are consistent with those of multi-element radio linked interferometer network (MERLIN) and very long baseline array (VLBA) observations 15,16 . In addition, there are few H 2 O maser features in the southern direction in contrast to those of the SiO masers. Figure 2 shows the position-velocity spot maps according to each maser lines. The H 2 O maser features show a wide-angle NE-SW biconical structure extended in slightly the NW and SE direction with respect to the position of central star (the mark "x"). The blue-shifted maser spots are dominant in the eastern part, and the red-shifted maser spots are dominant in the western part. In the case of the SiO maser, the maser components we see are mainly pumped in the tangential direction, so they show a ring-like structure in the line-of-sight. However, their partially clumped structure rather than a fully populated ring and various velocity ranges show an incoherency of the maser spot distribution. The pulsation of the stellar photosphere propagates through the SiO maser region, where dust starts to form. Al 2 O 3 is likely to be the first nuclei created, at least in some lower-mass stars, and the associated SiO maser monitoring has been used to relate their appearance to the kinematics within a few R * 17 . Each SiO maser line is distributed irregularly at any single epoch and, over the pulsation cycle, can show outflow and inward motion in different regions 18 , possibly associated with inhomogeneous dust formation. This would lead to local differences in the efficiency of radial acceleration through radiation pressure on grains and could be the cause of the irregular appearance of the 22.  Fig. 3. We can confirm the previous characteristics of the maser features and distributions according to each maser transition. In order to determine the position of the central star, we performed the ring fitting based on all spot distributions of the 43.1/42.8/ 86.2/129.3 GHz SiO masers. We used all of the SiO maser spots to improve the reliability and position accuracy. The mark "x" indicates the position of the central star. The ring radius obtained with the ring fitting in Fig. 2 is compared with that of the Gaussian fitting for the spot distribution histogram. Their radii obtained from both methods are consistent within the errors. Interestingly, the 86.2/129.3 GHz maser spots are located at outer regions from the central star than the 43.1/42.8 GHz maser spots (Fig. 3, Table 1). This trend has persisted at other epochs as shown in Fig. 4 14 . Therefore, we can confirm that the 42.8 GHz SiO maser is located inside the 43.1 GHz maser, the 86.2 GHz maser is located at the outer region of the 43.1 GHz maser, and the 129.3 GHz maser is located at the outermost region toward VX Sgr.   [19][20][21] , which is consistent with the theoretical excitation conditions of the v = 2 maser 22 . However, the SiO masers of the higher rotational transitions at the same vibrational state v = 1 (SiO v = 1, J = 1-0, J = 2-1, J = 3-2) shown in Fig. 3 and Table 1 arise from the outer region, i.e., the higher rotational transition masers arise further from the central star. In the case of the SiO v = 1, J = 1-0 (43.1 GHz) and J = 2-1 (86.2 GHz) masers, shock-enhanced simulation models 23 accurately predict the larger radii of the 86.2 GHz SiO maser compared with those of the 43.1 GHz SiO maser as a function of the stellar pulsation phase. Their model is also supported by our observational results in which the 86.2 GHz SiO maser intensity is stronger than that of the 43.1, 42.8 GHz SiO masers in VX Sgr. The predominantly collisional pumping model predicted that the 86.2 GHz SiO maser amplification is greater than that of the 43.1 GHz SiO masers for intermediate densities 24 . As an another possibility, the line overlap 25 can be proposed for explaining the larger ring size of 86.2 GHz maser compared to that of 43.1 GHz SiO   (Figs. 3 and 4). The single-dish spectrum (Fig. 2) of the 129.3 GHz maser shows more complex features and a wider velocity range (form −5 to 18 km s −1 ). The variation in line profile and peak velocity as a function of the stellar phase in the 129.3 GHz maser was different from that of the 86.2 and 43.1 GHz masers in single-dish monitoring observations 26 . These results may originate from the difference in their excitation conditions as discussed below.

H 2 O and SiO maser maps and their morphological differences.
Pumping mechanism for the 129.3 GHz SiO maser. We assumed that the masing regions of higher rotational transitions would be closely related to the optical depth of the circumstellar envelopes. The SiO maser population inversion (whether radiatively or collisionally pumped) is caused by self-trapping of photons corresponding to the Δv = 1 ro-vibrational transition 27 . At large opacities, spontaneous de-excitation of the SiO molecule via the Δv = 1 ro-vibrational transition becomes more difficult at higher rotational transitions 25 . This is possible even for radiative pumping in flattened regions such as in thin shells, which can be optically thick tangentially and thin radially. Thus population inversion and amplification of the higher rotational transition SiO masers is favored at progressively higher radii. Our results show a greater increase in radii and other differences between J = 3-2 and J = 2-1 masers, compared with those between J = 2-1 and J = 1-0 (Figs. 3 and 4). Thus the higher rotational transition masers in the SiO v = 1, J = 1-0, J = 2-1, J = 3-2 maser lines seem to arise further from the central star. However, we cannot exclude the radiative models including line overlap 25 for the further distance of the SiO v = 1, J = 3-2 maser from the central star. In addition, the ring radius of 129.3 GHz maser increases near optical maximum as shown in Fig. 4, which displays the ring radius variation of the 129.3 GHz maser at multiple epochs. This fact directly supports the radiative pumping for the 129.3 GHz maser and also suggests that the radiative pumping is more dominant than the collisional pumping in the higher J = 3-2 rotational transition for VX Sgr differently from the J = 2-1 maser.

Discussion
In the viewpoint of the collisional pumping of SiO masers caused by shocks, one can expect to trace the outward shock propagation near the photosphere 28 (Table 1, Supplementary  Tables 2 and 4) masers will allow us to trace the stratification structure of the excitation conditions, and inward/outward motions of maser clumps in different regions possibly associated with local shock and inhomogeneous dust formation.
Finally, our simultaneously astrometrically registered maps of the SiO and H 2 O masers will provide important observational constraints for the local difference in the efficiency of radial acceleration through radiation pressure on grains and the cause of the irregular appearance of the 22.2 GHz H 2 O masers above the dust-formation layers. In addition, the variation of the 129.3 GHz SiO ring radius with the optical light curve suggests that radiative pumping is dominant in the red supergiant VX Sgr. The passed signal is recorded onto the 1 Gbps Mark 5B (MK5B) recorder with 16 base band channels (BBCs) and 16 MHz bandwidth (512 channels) for each. We set the 6 intermediate frequencies of BBCs in the K and Q band, and 2 in the W and D band and the frequencies are arranged randomly in order to avoid high side peaks. Recorded data were correlated by the Distributed FX 31,32 software correlator. The synthesized beam size, position angle, and system noise temperature obtained from the observation for each frequency are listed in Supplementary Table 1.

Methods
We used fringe finder and continuum delay calibrator sources J1743-0350 and J1833-2103, respectively. Observations were performed over 5 h with 2 min scans for the target to continuum delay calibrator. Fringe finder calibrator source was observed every 1 h to exclude an instrumental delay and bandpass calibration. Amplitude calibration was done from the Astronomical Image Processing System (AIPS) task APCAL that measured the system noise temperature and gain variation data. There are no bright continuum delay calibrators stronger than 500 mJy (recommend for observations using 1 Gbps recording mode in KVN) and a separation angle within 4°from VX Sgr. Our delay calibrator J1833-2103 is strong but it has a complicated structure of a gravitational lensing 33 and 6.06°separation angle from the target source.
SFPR data reduction. The data reduction was performed using the AIPS package. Basic data treatment followed the standard procedure for the phase-referencing line imaging method in the K band 13 . First, we used the bright fringe finder to remove the residuals of large group delays and delay rates in the target itself. Second step is to trace the residuals of fringe phases using a nearby continuum delay calibrator, whose solutions were applied to the visibilities of the target. This step enables us to compensate for large and rapidly changing phase errors during the observations. We imaged the continuum delay calibrator using DIFMAP 34 in order to produce a source brightness model composed of CLEAN component. This was used in the AIPS task FRING to find the multi-band group delay and phase residuals.
The method consisted of transferring the multi-band delay solutions of the continuum delay calibrator and the phase rate solution of a strong H 2 O maser channel, which is copied from the K band to the calibration solutions of the SiO masers at high frequency band. The solutions were multiplied by the ratio between low and high frequencies. This method is the basis of non-integer SFPR 13 . The main point of the SFPR and non-integer SFPR method is transferring the low frequency phase solution to the high frequency [11][12][13] . Non-integer SFPR enables to make high frequency VLBI images for maser lines, even though a continuum delay calibrator is weak at the high frequency. We adopted a signal-to-noise threshold of maser spot identification to be 5-10, dependent on the image noise level 20 , in the AIPS task SAD.  Supplementary Fig. 2. The flux of the single dish and VLBI were about 1.5 times stronger at 43.1 and 42.8 GHz, but the morphology on the spot distribution was similar to that of March 27, 2016 (φ = 0.67). The 86.2 GHz maser shows a stronger flux than that of 43.1 and 42.8 GHz masers at both two epochs.
The SiO maser ring fitting results are shown in Supplementary Table 2. The size distribution tendency of each frequency radii was the same within the error range as those of March 27, 2016 (φ = 0.67) (42.8 < 43.1 < 86.2 GHz ring radius). Supplementary Fig. 3 shows the registered map of the SiO masers and the spot distribution histograms from the center of the fitting in each maser line. The 42.8 GHz SiO maser is located at the innermost region, and the 86.2 GHz is at the outermost region. The stellar position difference determined from the ring fitting center using all SiO maser spots is ΔR.A. = 1.23 and ΔDec. = −0.23 mas between these two observations (Table 1 and Supplementary Table 2).
Recovering fluxes. Figure 2 and Supplementary Fig. 2 include the spectra of the VLBI total power and recovered flux (dotted and dashed lines) together with a single-dish (solid line) spectra obtained from the closest date to the VLBI observations. Because the beam size of the VLBI is very narrow compared to that of the single dish, there is a missing flux in the VLBI observations. Also, the different baseline removal method causes flux difference between the VLBI total power and single-dish spectra. In the case of a single-dish observation (position-switching mode), we observe the empty sky to remove the baseline from the target spectrum. On the other hand, the VLBI is calibrated using bandpass data obtained from fringe finder observations. Moreover, the maser clumps with a larger angular size compared to the VLBI beam are resolved out in the cross-correlated spectra (recovered flux). The ratios of the VLBI total power to single-dish flux and recovered values (the ratio of recovered to total power fluxes) are listed in Supplementary Table 3. In the 42.8 GHz SiO maser, the relatively high recovery values seem to originate from a more compact maser spot in the 42.8 GHz maser than that in the 43.1 GHz maser.
Astrometric uncertainty. Our registered maser maps consist of the conventional PR maps of the H 2 O maser, with respect to J1833-2103, and SFPR maps of the SiO masers, with respect to the H 2 O maser positions determined from the PR maps. Therefore, we need to analyze astrometric errors related to these two kinds of registered maps. The propagation of the phase errors in the PR and SFPR results has been thoroughly studied using analytic and numerical methods 11,13,35,36 . Our observational parameters, such as 4 min for the source switching cycle time and 6.06°for the separation angle, yielded a somewhat large amount of errors in the PR maps. Following the formulae in previous work 35 , we can estimate the residual phase errors per baseline in the PR method:~90°phase noise from the dynamic tropospheric terms,~120°phase noise from the static tropospheric terms,~2°p hase noise from the dynamic ionospheric terms, and~15°phase noise from the static ionospheric terms. In addition, the large uncertainty of the antenna positions (cm order) of the KVN yields about 20°of phase noise in our PR results. The tropospheric contributions to the error budget can be removed in the SFPR results of the KVN due to the simultaneous multi-frequency observational capability 11 . Therefore, only the ionospheric terms contribute to the phase errors in the SFPR method, which yield 5°of phase noise from the dynamic ionospheric terms and 30°o f phase noise from the static ionospheric terms. Following the astrometric accuracy of the H 2 O maser positions in KVN PR maps, it can be estimated to be around 2 mas 13  O maser emission using the PR at low frequency is important for the accurate estimation of the position of the central star assumed to be located at the center of the SiO maser distribution, which is determined by the SFPR method. However, the uncertainty of the absolute astrometry will be large due to the large separation between the target and delay calibrator. This could be improved by including geodetic blocks or using GPS data to reduce the residual zenith path-length error.
Ring fitting of the SiO masers. We fitted the SiO maser ring structure of 43.1, 42.8, 86.2, and 129.3 GHz maser lines using the least-squares fitting method adopting the IDL online procedure "mpfitellipse.pro" (http://cow.physics.wisc.edu/ craigm/idl/idl.html). We did not apply the weights function because we assumed that the SiO maser is a perfect circular structure. Smallest chi-square error was selected as the best fitting result. The ring fitting results were expressed as a single center position with the ring radius range. The results of the estimated ring fitting uncertainty are listed in Table 1 and Supplementary Table 2.  Table 4). It is because we could not trace delay rate from non-integer SFPR according to the weather and observation conditions. The system temperatures were 620, 600, 550, and 500 K, the peak fluxes were 15.51, 7.27, 9.27, and 4.54 Jy beam −1 , and the rms values were 60. 34, 25.83, 20.42, and 18.65 Jy beam −1 m s −1 at each epoch, respectively. Compared with the 129.3 GHz map of Fig. 2, the blue-shifted maser features in Region D are the strongest on φ = 1.10. The red-shifted features in South (B') appear in all phases and show the strongest intensities among the other features near the minimum phase (φ = 0.50, 0.54, and 0.67). However, the feature in Region C did not appear at three epochs in 2015. It seems that the feature in Region C is newly generated at φ = 0.67. Region B was the strongest at the φ = 0.11 phase and tends to decrease over time that looks independent of the pulsation cycle. In contrast, Region A is the strongest in the last phase at φ = 1.10. Region A becomes richer and stronger from φ = 0.11 to 1.10. Relatively, maser spots are distributed widely at φ = 1.10 compared to the other epochs, and its ring size is the largest.  Table 1 and  Supplementary Tables 2 and 4.
Data availability. Raw data were generated at the Korea-Japan Correction Center (KJCC) in Daejeon. Derived data supporting the findings of this study are available from the corresponding author on request.