Proper motions of water masers in the star-forming region IRAS 23139+5939^†‡

Abstract

We observed H₂O (6₁₆→5₂₃) maser emission associated with the high-mass star-forming region IRAS 23139+5959 using the KaVA, a combination of very long baseline interferometry (VLBI) arrays between the Korean VLBI Network and the VLBI Exploration of Radio Astrometry (Japan). Through multi-epoch KaVA observations, we detected three groups of maser features, two of which coincide with those previously detected by the Karl G. Jansky Very Large Array. By determining the maser proper motions, we found that the first of maser groups exhibits an expanding motion that traces a wide-angle outflow almost along the line of sight, while the second one seems to be associated with the envelope of an H ii region. We discuss the star formation activity in IRAS 23139+5939, which may be reflected in the high variability of H₂O masers associated with an outflow seen from the front.

1 Introduction

Massive young stellar objects (YSOs) are deeply immersed in a large amount of circumstellar material. Consequently, these cannot be observed directly at optical or even in some cases at near-infrared wavelengths. Fortunately, most embedded YSOs can be detected with radio telescopes at centimeter wavelengths, and this allows us to investigate their nature and study their formation process. However, due to their relatively low flux density in radio continuum emission, it has been difficult to perform a detailed study of the continuum emission except in limited cases. For example, this emission is mainly associated with the signposts of newly born massive stars, such as radio jets and H ii regions. In addition to this, there is no entirely accepted model to explain the formation process of massive stars. Because of this, additional observational techniques are required, such as observations with a high angular resolution, allowing us to study the gas very close to massive YSOs.

Within this context, observations of maser transitions of several molecular species, particularly water (H₂O), in the proximity to massive YSOs provide a powerful diagnostic tool to investigate the first stages of massive star formation. H₂O maser emission is a selective tracer that can occur either in molecular outflows (e.g., Chernin 1995; Torrelles et al. 1998, 2014; Furuya et al. 2000; Patel et al. 2000) or, in rare cases, such as circumstellar disks (e.g., Shepherd & Kurtz 1999; Seth et al. 2002; Rodríguez-Esnard et al. 2014). Maser features (a cluster of spots in individual velocity channels) are like “test particles,” tracing the masing gas’ three-dimensional velocity field very close to YSOs. The combination between line-of-sight velocities and the proper motions of maser features is crucial for determining whether a group of features is rotating, contracting, or expanding, without any assumption of the geometry. The frequent tracking of the masers in temporal scales of from a few weeks to months is necessary to minimize the time-evolutionary effects of the maser features themselves and to be able to calculate their reliable proper motions (e.g., Chibueze et al. 2012).

Maser motions may be useful to study the stability of the continuum jet and the possible entrance of material at the interface between the jet and the surrounding interstellar material. One also expects to see maser features moving along the outflows, tracing circumstellar disks, or delineating the edge of isotropic gas ejections, etc. The new detection of an outflow activity may allow us to identify a new and previously unseen center of massive star formation activity (e.g., Torrelles et al. 2001, 2002). In the earliest stages of massive YSO evolution, some of these centers are likely associated with unexpected phenomena such as “short-lived” episodic ejection events characterized by a few decades’ kinematic age. The powering mechanism of non-collimated outflows associated with such episodic events is still unclear in current star formation theories. These facts open new challenges in the study of the earliest stages of stellar evolution (e.g., Torrelles et al. 2011; Trinidad et al. 2013). In order to see the fates of the short-lived events, we should monitor the isotropic gas ejections for several years.

IRAS 23139+5939 (G111.25−0.77) is a high-mass star-forming region with a luminosity of ∼2 × 10⁴ L_⊙ (Sridharan et al. 2002).¹ It is located near the S157 region at a distance of 3.4 kpc (Choi et al. 2014). This region has been studied in the radio continuum and with different molecular tracers (Tofani et al. 1995; Trinidad et al. 2006; Wu et al. 2010; Obonyo et al. 2019). Methanol maser studies have been performed by several authors (Szymczak et al. 2018; Durjasz et al. 2019) showing high variability with time. On the other hand, infrared (IR) and sub-millimeter (sub-mm) observations are reported by the United Kingdom Infrared Telescope (Varricatt et al. 2010), the Wide-field Infrared Survey Explorer (Wright et al. 2010), the Midcourse Space Experiment (Urquhart et al. 2011), the Bolocam Galactic Plane Survey (Ellsworth-Bowers et al. 2015), and the Submillimetre Common User Bolometer Array (Di Francesco et al. 2008) confirming its nature as a YSO. CO observations reveal intense emission by tracing molecular outflow activity (Wouterloot et al. 1989; Beuther et al. 2002). Using tracers of molecular outflow emission such as SiO (J = 2→1), SiO (5→4), and HCO⁺ (1→0), Sánchez-Monge et al. (2013) found extended high-velocity wings. However, a bipolar morphology has not been revealed in detail. The combination of new results in this paper with those of previous observations of H₂O masers (Tofani et al. 1995; Trinidad et al. 2003, 2006; Goddi et al. 2005; Choi et al. 2014) allows us to clarify such an outflow morphology.

This region has also been studied by Trinidad et al. (2006) using the Very Large Array (VLA) in A-configuration at 1.3 and 3.5 cm. They found two radio continuum sources (I23139 and I23139S) at wavelengths of 3.5 cm and two H₂O maser groups. The latter is spatially associated with the continuum sources, obtaining an accuracy of ∼10 mas in relative positions between I23159 and the maser features. The most interesting result comes from the spatial distribution of the water maser group associated with the strongest continuum source (I23159); the maser distribution exhibits a shell-like structure, very close to the continuum source’s peak emission. These masers have blueshifted velocities and are not gravitationally bound. Therefore, these could trace expansion or contraction. The expansion scenario seems to be the most favorable because the continuum source, apparently hosting maser features, could be a thermal jet or stellar wind. Furthermore, I23159 also seems to be the center of a CO molecular outflow observed in the region.

In addition, another maser group is spread out on a strip of |$\sim {0{^{\prime \prime}_{.}}2}$| and seems to be aligned in the northeast–southwest direction pointing toward the weaker continuum source (I23139S). This group is composed of three maser features; two of them have blueshifted velocities, and the other one has a velocity close to the source’s systemic velocity (−44.1 km s⁻¹; Bronfman et al. 1996). The spatial distribution of these features could also trace an outflow.

In order to build the three-dimensional velocity field of the H₂O maser features in IRAS 23139+5939 and investigate their nature, in this paper we present very long baseline interferometry (VLBI) multi-epoch observations carried out with KaVA, the combined VLBI array of three 21 m telescopes of the KVN (Korean VLBI Network; Lee et al. 2014) and four 20 m telescopes of VERA (VLBI Exploration of Radio Astrometry; Kobayashi et al. 2003). We detected three maser groups, two of them associated with those detected with the VLA: Groups 1 and 2. Based on their calculated proper motions, we found that the water maser features spatially associated with the continuum source I23159 (Group 1) are tracing expansion motions, while the maser features located to the southeast of I23139 (Group 2) seem to be associated with an H ii region.

In section 2, we describe the KaVA observations for maser proper motion measurement, and results and discussion are shown in sections 3 and 4, respectively. Our main conclusions are given in section 5.

2 Observations and data analysis

We observed H₂O masers in IRAS 23139+5939 with seven telescopes of KaVA (proposal ID KaVA18A-05) in three sessions in 2018; on February 6, March 29, and May 22, each one with a duration of 6 hr. As the seven KaVA telescopes form baselines from 300 to 2300 km, KaVA is relatively sensitive to extended radio emission. This has been demonstrated in earlier publications (e.g., Niinuma et al. 2014; Matsumoto et al. 2014). The unique capabilities of high-accuracy astrometry with dual beams and four frequency-band simultaneous observations have also been highlighted in VERA (VERA Collaboration 2020) and KVN (e.g., Yoon et al. 2018), respectively, and these are nowadays applicable to these KaVA observations. However, their applications are out of this paper, because the observations were carried out in the KaVA's early phase.

In order to perform the KaVA observations, continuum sources 3C 454.3, J2312+7241, and J0019+7327 were used as flux density, instrumental delay/phase, and band-pass calibrators, respectively. Parameters of the KaVA observations are given in table 1. Received signals were digitized at a rate of 2 Gbps, filtered into 16 base-band channels (BBCs), each one with a bandwidth of 16 MHz, and recorded at a rate of 1 Gbps. One of the BBCs covered a velocity width of 210 km s⁻¹, including an LSR velocity range from ∼−65 to 0 km s⁻¹ to detect water masers. Recorded signals were processed using the Daejeon Correlator in the Korea–Japan Correlation Center to yield 2048 spectral channels per BBC, corresponding to a velocity spacing of ∼0.11 km s⁻¹.

We reduced the data using the standard procedures of the NRAO AIPS software package, which was handled by the Python/ParselTongue scripts dedicated to the KaVA data.² We calibrated the visibility amplitudes using the corresponding antenna gains and system temperature tables.

Instrumental clock parameters and complex band-pass characteristics were determined and calibrated using the scans in the three calibrators mentioned above. We adopted a rest frequency of the H₂O 6₁₆→5₂₃ maser transition line of 22.23508 GHz. The data were cross-calibrated in phase and amplitude using the velocity channel including the most compact and the brightest maser spots at V_LSR = −47.0, −56.7, and −56.7 km s⁻¹ in the three sessions, respectively, for reference. For sessions 2 and 3, the maser spot used for reference was the same, but it was not detected for session 1, so another was used instead (the most intense at that time), located about 210 mas further away from that of sessions 2 and 3. Calibrated visibilities were naturally weighted to produce the synthesized beam of 2.1 × 1.4 [mas] at a position angle of |$-{26{^{\circ}_{.}}6}$| (for session 2).

These procedures allowed us to generate the synthesis maps of the maser emission with an rms noise ranging from ∼50 to ∼350 mJy, depending on the session and the existence of bright maser spots in the map. The dynamic range of the image cube (the highest ratio of the peak intensity of a maser spot to the noise level in the same velocity channel) achieved 365, while maser spots with a signal-to-noise ratio higher than 8 were identified. In this condition, we obtained a relative positional accuracy of about 0.1 mas, which will give us a reliable estimation of maser proper motions.

Maser spots brighter than 0.4 Jy were found automatically within an area of |${0{^{\prime \prime }_{.}}4} \times {0{^{\prime \prime}_{.}}4}$|⁠. All identified maser spots were fitted with two-dimensional elliptical Gaussians, by determining their positions, flux densities, and radial velocities. Then we further flagged out, carefully by eye, false maser spots raised by bright side lobes around bright maser spots. The remaining spots were grouped into maser features; each one of them corresponds to a physical gas clump composed of maser spots located within ≲1 mas in the sky and ≲2 km s⁻¹ in velocity. We only chose maser features that host at least four maser spots and that are isolated in the sky and velocity. Nevertheless, it was difficult to flag out all possible false features and identify all reliable features, in which some fainter but real features close to brighter features were flagged out while some isolated but false features could contaminate.

We determined the absolute coordinates of the maser spots used as reference mentioned above using the fringe-rate mapping method (Walker 1981). For the reference spot observed in the second session, when it was the brightest in the three sessions, we obtained its absolute coordinates: RA(J2000.0) |$= {23^{\rm h}16^{\rm m}1{_{.}^{\rm s}}3676} \pm {0{^{\rm s}_{.}}016}$|⁠, Dec (J2000.0) |$= + {59^{\circ }55^{\prime }28^{\prime \prime }.43} \pm {0{^{\prime \prime }_{.}}10}$|⁠. The coordinate uncertainty was limited to |$\sim {0{^{\prime \prime}_{.}}1}$| due to the relatively low flux densities of reference spots, and the sparse fringe rate calibration using scans on calibrators observed every ∼30 min. Nevertheless, this is enough to associate them with the radio continuum source detected in the region. All maser feature position offsets were measured relative to the stronger maser spot detected in the second session of KaVA observations.

Since the strongest maser feature, which was used as a reference for the calibration, was not persistent in all three sessions (it was not detected in session 1, thus another maser feature was selected), the strongest one was not used as a positional reference to align the maps; other persistent features were identified. In order to align the maps, two maser features of Group 1, persistent in all three sessions, were selected, the radial velocity and position of which (for session 2) are about −42.6 and −35.2 km s⁻¹, and (17.65, 116.04) and (−5.92, 3.81) mas, respectively. The reference position to align the maps was determined by calculating the mean position (RA and Dec) of both selected maser features, which was (5.87 ± 0.01, 59.93 ± 0.02) mas [RA (J2000.0) |$= {23^{\rm h}16^{\rm m}10{^{\rm s}_{.}}3684} \pm {0{^{\rm s}_{.}}0008}$|⁠, Dec (J2000.0) |$=+{59^{\circ }55^{\prime }28.^{\prime \prime }492} \pm {0{^{\prime \prime}_{.}}060}$|⁠; see table 2].

Given that the calibration of the second and third sessions was performed using the same maser feature, we first aligned the maps of these two sessions, and then the map of the first session was adjusted. In order to do this, we aligned the mean position of the two maser features selected in the third session, assuming that these have been moved at a constant velocity from positions in the second session, i.e., for 54 days (time separation between the observations of sessions 2 and 3). Finally, assuming that these features between sessions 2 and 3 are moved with constant velocity, this velocity was calculated using the mean position and the time between both sessions. Accordingly, considering a time span of 51 days between the first and second session, the estimated velocity between sessions 2 and 3 was also applied to the mean position of the two maser features selected in the first session. Thus, after aligning the maps, all the positions of the maser features from the three observed sessions were determined. Then, some of the features detected in two and/or three epochs were selected to calculate their proper motion.

3 Results

We identified 19, 23, and 13 H₂O maser features toward IRAS 23139+5939 in the three observing sessions, respectively, the velocity range of which spans from V_LSR = −57 to −33 km s⁻¹. Table 2 gives the parameters of the identified H₂O maser features and figure 1 shows the distribution of all the maser features detected during session 2 toward IRAS 23139+5939. We notice that the maser features are distributed in three groups, labeled as 1, 2, and 3. Group 1 contains the most maser features, while the furthest one, Group 3, has few maser features. This behavior was observed practically in the three sessions. Specifically, H₂O maser features in Groups 1 and 2 were detected in all sessions, whereas those of Group 3 were only detected in the first and second sessions.

In general, the peak intensity of the maser features range from a few tenths to a few tens of Jy beam⁻¹. However, a powerful maser feature was detected in the second and third sessions, with a peak intensity of ∼100 and ∼67 Jy beam⁻¹, respectively (table 2). This maser feature, with a velocity of −56.8 km s⁻¹, is blueshifted with respect to the systemic velocity (−44.1 km s⁻¹) of the IRAS 23139+5939 region and was used as reference maser for calibration, but it was not detected in the first session. Figure 1 also shows the spectra of the three maser groups from the second session, in which the maser emission was the most populated. Most of the strongest maser features in this session are blueshifted, except in Group 3. A similar trend is also observed in other sessions. This trend is consistent with the variability of the H₂O maser emission toward IRAS 23139+5939 studied by Trinidad et al. (2003).

Figure 2 shows the H₂O maser features superimposed on the contour map of the 1.3 cm continuum emission from I23139 detected with the VLA (Trinidad et al. 2006), as well as the maser features in Groups 1 and 2 detected in the second session of the KaVA observations. The absolute coordinates of these maser features detected with KaVA do not match exactly with those found with the VLA, which could be explained by their uncertainty (which is in fact worse than |$\sim {0{^{\prime \prime}_{.}}1}$|⁠; see section 2).

In order to carry out a visual comparison between the maser feature distributions in the KaVA observations with those in the VLA, a shift of |$(-{0{^{\prime \prime}_{.}}23}, {0{^{\prime \prime}_{.}}21})$| was applied to the position of the KaVA maser features determined with the fringe-rate mapping method (see figure 2). Comparing both maser distributions, we found that the maser features of the central region of Group 1 (of the KaVA observations) seem to be related to the cluster of the masers (detected with the VLA) spatially associated with the radio continuum source I23139. In addition, we found that the maser features of Group 2 are associated with those H₂O masers found in the southeast from I23139 with the VLA. However, as it was mentioned above, this comparison is only for visual purposes, and because the offset applied to the KaVA data is more significant than its positional uncertainty, the findings from this comparison should be carefully considered.

On the other hand, the maser features in Group 3 were not detected in the previous VLA observations or in the third session of the KaVA observations. We could explain the latter by the fact that H₂O masers observed in star-forming regions show significant change in their global distribution in long and even short time scales (e.g., Trinidad et al. 2003).

We calculated the relative proper motions of the maser features detected in two and/or three sessions by assuming that these have linear motions and by applying a linear fitting to their positions. Since the H₂O maser emission turned out to be highly variable on a time scale of few weeks, it was only possible to measure a limited number of maser proper motions; nine in total in the entire IRAS 23139+5939 region (seven in Group 1 and two in Group 2; see table 3). The determined geometric center (“center of motion”) of the maser features having proper motion in Group 1 is RA (J2000.0) |$= {23^{\rm h}16^{\rm m}10{^{\rm s}_{.}}3676} \pm {0{^{\rm s}_{.}}0003}$|⁠, Dec (J2000.0) |$=+{59^{\circ }55^{\prime }28.^{\prime \prime }432} \pm {0{^{\prime \prime}_{.}}042}$|⁠, corresponding to a proper motion of (0, 0) mas yr⁻¹. Calculated proper motions seem to show that all maser features are being moved northward. However, this result is not consistent with previous studies (e.g., Choi et al. 2014). The mean of the north–south components of the proper motions is ∼0.7 mas yr⁻¹ (1 mas yr⁻¹ =16 km s⁻¹), which is larger than that reported by other authors (e.g., −2.3 mas yr⁻¹: Choi et al. 2014). We suppose that this trend is not real and could be attributed to the reference frame moving together with the reference feature. In order to remove this bias, the mean of the north–south components ∼0.7 mas yr⁻¹ was subtracted from all maser proper motions.

Figure 3 shows proper motions, without the bias to the north, of H₂O maser features. We notice that the maser features in Group 1 are distributed within a region of ∼200 mas along the north–south direction around I23139 (figure 2) and are moved essentially westwards. In contrast to this, the maser features in Group 2 are aligned along the northeast–southwest direction and are moved northeast–southwestwards. This tendency reproduces clearly that scenario previously found by Choi et al. (2014) (see figure 4). Moreover, the maser feature distribution seems to be consistent with its evolution between the epochs of the VLA (Trinidad et al. 2006) and our KaVA observations (figure 2). It exhibits east–west expansion between the maser features in Groups 1 and 2 in the time baseline of ∼17 yr among these observations, which is consistent with that expected from relative maser proper motions derived from this paper (∼16 mas). However, given the positional uncertainty of KaVA observations, this result should be taken meticulously and must be confirmed with new observations.

4 Discussion

In recent years, evidence that massive YSOs undergo an outflow stage governed by same processes as low-mass stars has increased. However, the detection of thermal jets has been very scarce. The aforementioned is particularly relevant since jets play a crucial role in forming low-mass stars, being the driving sources of molecular outflows (e.g., Anglada et al. 2018). Therefore, it is not still clear whether high-mass stars are formed in the same way as low-mass stars. In this way, it is necessary to carry out individual studies of massive YSOs. In particular, IRAS 23139+5939 is an excellent candidate to investigate whether the thermal jet could be driving the bipolar outflow observed on a large scale, where H₂O masers have been detected towards a possible thermal radio jet (Trinidad et al. 2006).

Making a comparison with other massive star-forming regions (e.g., Cepheus A, AFGL 2591; Torrelles et al. 2001, 2011; Trinidad et al. 2013), there is a relatively small number of H₂O maser features detected toward IRAS 23139+5939 (see figure 1), and most of them do not persist over timescales of a few weeks. As a result, it has been challenging to find a regular and systematic maser distribution and measure a significant number of proper motions. Because of this, only a qualitative description is given. However, we discovered that the groups of H₂O maser features have persisted since the previous VLA observations (see figure 2). This allows us to track the evolution of the maser distribution and discuss further the properties of the H₂O masers concerning I23139 and the morphology of the outflow from this YSO in more detail.

The relative positions of the H₂O maser features concerning the 1.3 cm continuum source I23139 have an accuracy of 10 mas in the VLA data (Trinidad et al. 2006). Considering that the maser features in Group 1 observed with KaVA are also spatially associated with those found in the VLA observation, these KaVA maser features are also associated with I23139, that being the pumping source. Trinidad et al. (2006) found that the VLA maser features spatially associated with I23139 are not gravitationally bound with this source but rather tracing expansion motions. Measured maser proper motions prove this suggested expansion.

From table 2, we notice that the strongest H₂O maser feature is found in the central subgroup of Group 1, spatially coinciding with the continuum source I23139. In this subgroup, maser features have the lowest and highest radial velocities, i.e., those with the highest blueshifts and redshifts, respectively. In particular, for the second session of KaVA observations, maser features with “blue” velocities of the central subgroup were more intense (up to tens of Jy) than those with “red” velocities (a few tenths of Jy). The latter could indicate that the “blue” maser features are in front of the continuum source, while the “red” ones are behind it, so we speculate that these look less intense due to the “red” maser emission being masked by a dense structure, for example, a circumstellar disk (see figure 5).

In the same way, for the maser source G353.273+0.641, Motogi et al. (2011, 2013) proposed a disk-masking scenario to explain the nature of the strong blueshift dominance, where an optically-thick disk obscures a redshifted lobe of a compact jet. Under this geometry, we expect an accretion disk/torus around I23139 in the plane of the sky, perpendicular to the line-of-sight (pole-on). Herein, the H₂O maser emission, which is also projected on to the plane of the sky, could be tracing a jet (or, more generally, a conical outflow). In this case, we can also expect the velocities of the H₂O maser features to increase when these are moved away from I23139.

Considering that the spatial distribution of H₂O maser features in Group 1 is extended over a scale of ∼200 mas, we expect an optically thick disk of about 300 au. Figure 5 shows the schematic summary of our consideration in terms of the outflow morphology hosting the H₂O masers in I23139. This scenario seems to be consistent with the spatial distribution, velocity field, and intensities of the H₂O maser features in Group 1 and their expanding motions. However, given the angle of inclination of the outflow in the proposed scenario and the few calculated proper motions, we cannot perform a detailed study.

Despite these limitations, this scenario also seems to be congruent with the observed spatial and kinematic distribution of a wide-angle molecular outflow toward IRAS 23139+5939. This outflow is closely oriented to the line-of-sight with extended high-velocity wings (Wouterloot et al. 1989; Beuther et al. 2002; Sánchez-Monge et al. 2013), where emission from blueshifted wings is more intense than that of redshifted wings.

In addition, Trinidad et al. (2006) suggested that I23139 is a thermal radio jet associated with a massive YSO, while the H₂O maser distribution, which we found, predicts the presence of an accretion disk around I23139. These results show that I23139 is forming a disk–YSO–jet system similar to that observed in low-mass YSOs. Under these conditions, I23139 could be an example of a high-mass protostar with a disk–YSO–jet system, suggesting that I23159 is being formed by accretion.

As previously indicated, the spatial extension of maser features is about 200 mas. Assuming that these features are tracing a wide-angle outflow (poorly collimated) and expanding with a velocity of ∼45 km s⁻¹ (obtained using their radial velocity and proper motions), we found a rough estimation of the dynamical age of 35 yr for the ejection traced by maser features. This dynamical age is greater than that calculated for W75-VLA 2 (13 yr; Torrelles et al. 2003) but comparable with Cepheus A HW2-R5 (33 yr; Torrelles et al. 2001).

In both cases, these values have been interpreted as short-lived episodic spherical winds from a massive object located at the center of the structure, which occurs during the early stages of their evolution. Although the dynamical age obtained for the outflow traced by the maser features is consistent with those obtained in other young protostars, this value must be taken carefully because the H₂O maser outflow seems to be closely oriented to the line of sight, and its spatial extension could be uncertain. Either way, if the estimated age is correct, this result could suggest that the continuum source, I23139 (spatially associated with the maser features of group 1), is in an early evolutionary stage of its formation, which is consistent with it being a radio jet.

In addition, we notice that the flux density of the continuum source I23139 also shows variability with time. With VLA observations in the 22 GHz band (both in the A configuration), it had a flux density of 0.98 mJy in 2001 (Trinidad et al. 2006) but 0.25 mJy 12 years later (Moscadelli et al. 2016). A similar behavior is observed in the 8.4 GHz band (0.92 and 0.53 mJy; see Tofani et al. 1995 and Trinidad et al. 2006, respectively). However, the spectral index estimated for I23139 between 8.4 and 22.2 GHz, with simultaneous observations (Trinidad et al. 2006), is very similar to that estimated between 15 and 22 GHz (Moscadelli et al. 2016), 0.64 ± 0.36 and 0.72 ± 0.30, respectively. The difference in the measured flux density could be attributable to the intrinsic variability of the continuum source, which could contribute to the variability of the H₂O masers detected in the region. Furthermore, we notice that the flux density of the H₂O maser emission is correlated with the intensity of the continuum source at the 22 GHz band. This possible correlation may be further tested by revisiting this maser source in the future.

Regarding the maser features of Group 2, Trinidad et al. (2006) found that these masers trace a large filament along the northeast–southwest direction and seem to be pointing toward the 3.6 cm continuum source I23139S, located |${0{^{\prime \prime}_{.}}35}$| to the southwest from these masers and interpreted as a thermal jet. Trinidad et al. (2006) speculated that the H₂O maser features could be associated with and excited by this continuum source. Although the distribution of H₂O masers observed with the VLA is maintained in KaVA observations in general, only those that are further north were detected. However, the directions of the H₂O maser proper motions do not seem to be consistent with I23139S being the driving source. Based on our observations, we found that maser features are roughly aligned with I23139S, but these point in the opposite direction, which rules out I23139S being the exciting source. In addition, these maser features seem to be tracing an expansion motion driven by an undetected young source. This scenario has also been found in other star-forming regions, for example in AFGL 2591 VLA 3-N (Trinidad et al. 2013), although in this case, the maser features trace an expanding spherical movement rather than a collimated motion.

Moscadelli et al. (2016) found a second weak continuum source in the 22 GHz band, which is close to the midpoint of the Group 2 of maser features detected with the VLA and KaVA. This continuum source could be a candidate for the driving source for this group of H₂O masers. However, this new continuum source has only been detected at 22 GHz, so it is not still possible to calculate its spectral index to research its nature. Nevertheless, using the flux density of the continuum source at 22 GHz and an upper limit for its flux density at 8.4 GHz (3σ; Trinidad et al. 2006), we estimate a rough spectral index ≲−0.17, which could be consistent with an optically thin H ii region. If this continuum source is responsible for exciting H₂O masers, these should trace a spherical expansion and not a collimated motion. Despite the proper motions of the maser features in Group 2 pointing in the opposite direction, their velocities are blueshifted, indicating that these features are in front of the continuum source and are moved on the plane of the sky. Under this scenario, the H₂O maser features may be on the surface of the H ii region and pumped by the central star (see figure 2).

5 Conclusions

Using the data of multi-epoch observations of H₂O masers with KaVA, we detected three groups of masers towards IRAS 23139+5939. Only Groups 1 and 2 coincide with the continuum source I23139. KaVA observations discovered Group 3, which was not detected in the previous VLA observations; it is located quite far from I23139 (approximately |${0{^{\prime \prime}_{.}}9}$|⁠). The H₂O maser emission toward IRAS 23139+5939 is highly variable with time. Because of this, very few maser features persist on time scales of several weeks, making it difficult to determine proper motions in a single set of VLBI observations.

Based on the proper motions of maser features, we can confirm that H₂O masers, which are associated with the continuum source I23139 (Group 1), are tracing expansion motions, possibly a wide-angle short-lived outflow with a rough dynamical age of 35 yr. Furthermore, using a very simplistic model, we suggest that these masers should be tracing an outflow seen from the front, where only masers located at the outflow’s base and on the outermost part of it are detected. Moreover, this model predicts the presence of a circumstellar disk around I23159, suggesting that this massive protostar is in an early evolutionary stage and forming in a similar way to low-mass stars. Finally, we found that maser features in the Group 2 seem to be tracing expansion motions on the surface of an H ii region.

Acknowledgements

VERA is operated by Mizusawa VLBI Observatory, respectively, branches of the National Astronomical Observatory of Japan, National Institutes of Natural Sciences. KVN in KaVA and the Daejeon Correlator in KJCC is operated by Korea Astronomy and Space Science Institute (KASI). MAT acknowledges support from DAIP-UG grant 163/2021. HI was supported by the JSPS KAKENHI Grant Number JP16H02167 and JP21H00047. IT–J acknowledges support from CONACyT, Mexico; grant 754851. EdelaF thanks Kagoshima University for financial support during the visit to the CGE, Institute for Comprehensive Education, Kagoshima University in 2019. He also renders thanks to UdG-PRO-SNI 2021 grant, the staff of the Institute for Cosmic Ray Research (ICRR-UTokyo), Coordinación General Académica y de Innovación (CGAI-UDG), Coordinación de Servicios Académicos (CSA-CUCEI-UDG), PRODEP-SEP through Cuerpo académico UDG-CA-499, Carlos Iván Moreno, Cynthia Ruano, Rosario Cedano, Diana Naylleli Navarro, Luis Fernando González Bolaños, Jorge Alberto Rodríguez Castro, Johnatan Barcenas Alvarez, Patricia Retamoza Vega, Laura Gonzalez Jaime, and Dulce Angélica Valdivia Chávez, Kumiko Sugimoto, and Michiru Ito for financial and administrative supports during Sabbatical year stay at the ICRR, the University of Tokyo in 2021.

Footnotes

References