MASSIVE STAR FORMATION TOWARD G28.87+0.07 (IRAS 18411−0338) INVESTIGATED BY MEANS OF MASER KINEMATICS AND RADIO TO INFRARED CONTINUUM OBSERVATIONS

J. J. Li; L. Moscadelli; R. Cesaroni; R. S. Furuya; Y. Xu; T. Usuda; K. M. Menten; M. Pestalozzi; D. Elia; E. Schisano

doi:10.1088/0004-637X/749/1/47

1. INTRODUCTION

Understanding the process of high-mass star formation represents a challenge from both a theoretical and observational point of view. By now, three main competing concepts of massive star formation have been discussed in the recent literature (see Zinnecker & Yorke 2007 and references therein): formation through coalescence of lower-mass stars in extremely dense stellar clusters (Bonnell et al. 1998); monolithic collapse of isolated, turbulence supported cores (McKee & Tan 2003); and competitive accretion where the most massive protostars benefit from the gravitational well created by the lower-mass stars surrounding them, thus boosting the accretion rate (e.g., Bonnell & Bate 2006). The available observational evidence is still inadequate to discern between these different scenarios, mostly because observations of massive young stellar objects (MYSOs) are complicated by the large distances of typically several kpc, implying small angular scales, and the fact that such objects are born still deeply enshrouded in dense dusty envelopes, opaque to optical and near-infrared wavelength radiation.

A possible way to distinguish between different models is to explore the physical conditions and the dynamical processes of the gas close to the protostars. However, present (sub)millimeter interferometers and infrared (IR) telescopes attain angular resolutions ⩾0 farcs 25, corresponding to thousands of astronomical units (AU) at the typical distances of MYSOs of several kpc. This is insufficient to investigate the kinematics of the gas in the immediate surroundings of the forming star, where crucial phenomena such as accretion, rotation, and ejection are likely most prominent. The above-mentioned limitation can be bypassed by very long baseline interferometry (VLBI) observations of maser lines, which attain linear resolutions of a few AU.

A few years ago, we started an observational project to study the high-mass star-forming process by comparing interferometric images of thermal lines of molecular tracers with multi-epoch VLBI studies in three well-known maser species (OH, CH₃OH, and H₂O). At present, three sources (IRAS 20126+4104, G16.59−0.05, and G23.01−0.41) have been already analyzed (Moscadelli et al. 2011; Sanna et al. 2010a, 2010b). The results demonstrate that the synergy between VLBI, multi-species maser observations, and (sub)millimeter interferometric observations of thermal molecular lines is crucial to achieve a multi-scale picture of the environment of newly formed massive (proto)stars. In this article we report on our VLBI maser study of the high-mass star-forming region (HMSFR) G28.87+0.07, and compare our results with available complementary interferometric data of continuum and thermally excited molecular line emission, as well as with subarcsecond resolution mid-IR images obtained with the Cooled Mid-Infrared Camera and Spectrometer (COMICS) on the Subaru Telescope.

In Section 2, we provide a review of previous observations toward this region. Section 3 describes the VLBI observations of the maser transitions, the Very Large Array (VLA) observations of the radio continuum emission at 1.3 cm and 3.6 cm, and the mid-IR observations performed with the Subaru/COMICS Telescope. In Section 4, we present the distributions of the three types of masers, the water maser proper motions, and the properties of the radio continuum. The results are discussed in Section 5 and the conclusions are drawn in Section 6.

2. THE HMSFR G28.87+0.07

The HMSFR G28.87+0.07 (alias IRAS 18411−0338) is located at a kinematic distance of 7.4 kpc⁶ (Codella et al. 1997; Furuya et al. 2008) and is associated with a pc-scale clump with a mass in excess of 10³ M_☉ and a bolometric luminosity of ∼2 × 10⁵ L_☉ (López-Sepulcre et al. 2009; Faúndez et al. 2004). From the observation of NH₃, the systemic velocity of the region (V_sys) with respect to the local standard of rest (LSR) is found to be 101.0 km s⁻¹ (Codella et al. 1997).

On a large scale, mid-IR observations at 8, 10.5, and 18.1 μm have resolved the IRAS source in two objects along the north–south direction (De Buizer et al. 2005; Furuya et al. 2008, their Figure 1). Both sources could be associated with a Midcourse Space Experiment (MSX) point source (Egan et al. 2003) and a 2.2 μm source with strong near-infrared excess, lying 3''–4'' to the southeast of the northern, brighter mid-IR object (Testi et al. 1994). Bolometer observations at 1.2 mm (Faúndez et al. 2004) have detected a peak of dust continuum emission associated with the two mid-IR sources and the near-IR source. Single-dish maps in the CS(5–4) and C¹⁸O(2–1) lines have revealed an arcminutes-scale molecular clump and a monopolar (redshifted) outflow is seen in the ¹³CO(2–1) line (Shirley et al. 2003; López-Sepulcre et al. 2009). Single-dish surveys have detected emission in dense gas tracers as CS(7–6), CS(2–1), and CS(1–0) (Plume et al. 1992; Bronfman et al. 1996; Anglada et al. 1996), and in the H110α radio recombination line (Sewilo et al. 2004; Araya et al. 2007), which suggests the presence of a diffuse H ii region. In addition, OH maser emission has been found in the main (1.665 and 1.667 GHz) hyperfine structure lines and in the 1.720 GHz satellite line, whereas the 1.612 GHz OH satellite has been observed in absorption (Little & Cesarsky 1982; Caswell & Haynes 1983). Weak CH₃OH and strong H₂O maser emission have been reported at 6.7 GHz (Caswell et al. 1995) and 22.2 GHz (Genzel & Downes 1979), respectively.

Interferometric observations at millimeter wavelengths have revealed multiple sources in this region (Furuya et al. 2008, their Figure 1) and show that the dust continuum peak coincides with the position of the northern, stronger mid-IR source. The mass (100 M_☉) from the 3 mm dust continuum emission and the rotational temperature (93 K) inferred from the CH₃CN lines suggest that one is dealing with a hot molecular core (HMC). In this study, we will focus on this HMC, as it contains all the tracers of high-mass star formation that we have investigated with our observations, namely, the masers and the free–free and mid-IR continuum sources. As shown by Codella et al. (1997) and Furuya et al. (2008), the HMC is seen in NH₃, HCO⁺, and CH₃CN line emission. The observations of Codella et al. (1997) were not sensitive enough to detect continuum emission at 1.3 cm, whereas our new images demonstrate the existence of a free–free continuum source inside the HMC (see below). A wide-angle, massive (⩾100 M_☉) bipolar outflow has been detected in the ¹²CO(1–0) line, oriented perpendicular to a CH₃CN velocity gradient detected in the HMC (Furuya et al. 2008). The absolute positions of the 1.665 GHz OH and 22.2 GHz H₂O masers, derived by Forster & Caswell (1999) using the VLA with an accuracy of 0 farcs 5, indicate that these maser emissions are associated with the HMC.

All of these features are consistent with the presence of a deeply embedded MYSO(s) inside the HMC, possibly accompanied by lower-mass objects lying in the neighborhoods of the HMC itself.

3. OBSERVATION AND DATA REDUCTION

3.1. VLA: 1.3 cm and 3.6 cm Continuum

We observed the HMSFR G28.87+0.07 (tracking center: R.A.(J2000) = $18^{\rm h}43^{\rm m}46\mbox{$.\!\!^{\mathrm s}$}24$ and decl.(J2000) = −03°35'30 farcs 4) with the VLA at X (3.6 cm) and K (1.3 cm) bands in both the C-array (in 2008 April) and A-array configurations (in 2008 November). At 3.6 cm in both configurations and at 1.3 cm in the C-array configuration, the continuum mode of the correlator was used, resulting in an effective bandwidth of 172 MHz. At 1.3 cm in the A-array we used mode "4" of the correlator, with a pair of intermediate frequency (IF) bands, each with a width of 3.125 MHz (64 channels) centered on the strongest H₂O maser line and a pair of 25 MHz wide IF bands (8 channels each) sufficiently offset from the maser lines to obtain a measurement of the continuum emission. The two bands centered at the same frequency (detecting both circular polarizations) were averaged to decrease the noise.

At X band, 3C 286 (5.2 Jy) and 3C 48 (3.1 Jy) were used as primary flux calibrators, while 1832−105 (1.4 Jy) was the phase calibrator. For the K-band observations, the primary flux calibrator was 3C 286 (2.5 Jy), the phase calibrator 1832−105 (1.0 Jy), and the bandpass calibrator (for the A-array only) 1733−130 (3.7 Jy).

The data were calibrated with the NRAO⁷ Astronomical Image Processing System (AIPS) software package using standard procedures. Only for the A-array data at 1.3 cm were several cycles of self-calibration applied to the strongest maser channel, and the resulting phase and amplitude corrections were eventually transferred to all the other line channels and to the K-band continuum data (Menten & Reid 1995). This procedure resulted in a significant (at least a factor of two) improvement in the signal-to-noise ratio (S/N).

Natural-weighted maps were made with task IMAGR of AIPS for both the continuum and the line data. Simultaneous observation of the line and continuum emission at K band made it possible to obtain the relative position of the continuum image with respect to the H₂O maser spots with great precision (∼0 farcs 1). Finally, we identified the maser spots that appear in both our Very Long Baseline Array (VLBA) and VLA images and, using the accurate absolute position information obtained from the VLBA observations, we re-centered the VLA spots and continuum image. The correction thus applied turns out to be <20 mas in each coordinate. We conclude that the absolute astrometrical precision of the continuum maps presented in this paper is ∼0 farcs 1.

3.2. VLBI: Maser Lines

We conducted VLBI observations of the H₂O and CH₃OH masers (at four epochs) and of the OH maser (at a single epoch) toward G28.87+0.07 at 22.2 GHz, 6.7 GHz, and 1.665 GHz, respectively. We give more technical details on the observations in the following subsections. In order to determine the masers' absolute positions, we performed phase-referencing observations in fast switching mode between the maser source and the calibrator J1834−0301. This calibrator has an angular offset from the maser source of 2 fdg 4 and belongs to the list of sources of the VLBA Calibrator Survey 3 (VCS3). Its absolute position is known with a precision better than ±1.5 mas, and its flux measured with the VLBA at S and X bands is 154 and 156 mJy beam⁻¹, respectively. Five fringe finders (J1642+3948, J1751+0939, J1800+3848, J2101+0341, J2253+1608) were observed for bandpass, single-band delay, and instrumental phase-offset calibration.

The data were reduced with AIPS following the VLBI spectral line procedures. For a description of the general data calibration and the criteria used to identify individual masing clouds, derive their parameters (position, intensity, flux, and size), and measure their (relative and absolute) proper motions, we refer to the recent paper on VLBI observations of H₂O and CH₃OH masers by Sanna et al. (2010a). For each VLBI epoch, H₂O maser absolute positions have been derived by employing two methods, direct phase-referencing and inverse phase-referencing, and the two procedures always gave consistent results. CH₃OH and OH maser absolute positions have been derived only by using direct phase-referencing. The derived absolute proper motions have been corrected for the apparent proper motion due to Earth's orbit around the Sun (parallax), the solar motion and the differential Galactic rotation between our LSR and that of the maser source. We have adopted a flat Galactic rotation curve (R₀ = 8.3 ± 0.23 kpc, Θ₀ = 239 ± 7 km s⁻¹; Brunthaler et al. 2011), and the solar motion of U = 11.1^+0.69_{− 0.75}, V = 12.24^+0.47_{− 0.47}, and W = 7.25^+0.37_{− 0.36} km s⁻¹ by Schönrich et al. (2010), who recently revised the Hipparcos satellite results.

3.2.1. VLBA: 22.2 GHz H₂O Masers

We observed the HMSFR G28.87+0.07 (tracking center: R.A.(J2000) = $18^{\rm h}43^{\rm m}46\mbox{$.\!\!^{\mathrm s}$}22$ and decl.(J2000) = −03°35'29 farcs 9) with the VLBA in the 6₁₆–5₂₃ H₂O transition (rest frequency 22.235079 GHz). The observations (program code: BM244) consisted of four epochs: 2006 April 23, 2006 June 30, 2006 September 28, and 2007 January 18. During a run of about 6 hr per epoch, we recorded the dual circular polarization with 16 MHz bandwidth centered at an LSR velocity (V_LSR) of 101.0 km s⁻¹. The data were processed with the VLBA FX correlator in Socorro (New Mexico) using an averaging time of 1 s and 1024 spectral channels.

Images were produced with natural weighting, cleaned and restored with an elliptical Gaussian with FWHM size of about 1.3 mas × 0.6 mas at a position angle (P.A.) of −4° (east of north), with small variation from epoch to epoch (shown in Table 2). The interferometer instantaneous field of view was limited to ∼2 farcs 7. In each observing epoch the on-source integration time was ∼2.5 hr resulting in an effective rms noise level of the channel maps (σ) in the range 0.006–0.015 Jy beam⁻¹. The spectral resolution was 0.2 km s⁻¹.

3.2.2. EVN: 6.7 GHz CH₃OH Masers

We observed G28.87+0.07 (tracking center: R.A.(J2000) = $18^{\rm h}43^{\rm m}46\mbox{$.\!\!^{\mathrm s}$}22$ and decl.(J2000) = −03°35'29 farcs 9) with the European VLBI Network (EVN) in the 5₁ − 6₀A⁺ CH₃OH transition (rest frequency 6.668519 GHz). Data were taken in four epochs (program code: EM061 and EM069) separated by ∼1 year: 2006 February 28, 2007 March 18, 2008 March 18, and 2009 March 17. In the first two epochs, the antennae involved in the observations were Cambridge, Jodrell2, Effelsberg, Hartebeesthoek, Medicina, Noto, Torun, and Westerbork. Since we discovered that the longest baselines involving the Hartebeesthoek antenna (e.g., the Ef-Hh baseline is ∼8042 km) heavily resolve the maser emission and do not produce fringe-fit solutions, the Hartebeesthoek antenna was replaced with the Onsala antenna in the third epoch. During a run of ∼6 hr, we recorded the dual circular polarization in two IF bands of width of 2 MHz and 16 MHz, both centered at an LSR velocity of 101.0 km s⁻¹. The 16 MHz band was used to increase the S/N of the weak continuum calibrator. The data were processed with the MKIV correlator at the Joint Institute for VLBI in Europe (JIVE, Dwingeloo, the Netherlands) using an averaging time of 1 s and 1024 spectral channels for each observing bandwidth.

Images were produced with natural weighting, cleaned and restored with an elliptical Gaussian with FWHM size of about 12 mas × 6 mas at a P.A. of 21°, with small variations from epoch to epoch (shown in Table 2). The interferometer instantaneous field of view was limited to ∼9 farcs 2. Using an on-source integration time of ∼2.2 hr, the effective rms noise level of the channel map (σ) varied in the range 0.005–0.01 Jy beam⁻¹. The 2 MHz band spectral resolution was 0.09 km s⁻¹.

3.2.3. VLBA: 1.665 GHz OH Masers

We observed G28.87+0.07 (tracking center: R.A.(J2000) = $18^{\rm h}43^{\rm m}46\mbox{$.\!\!^{\mathrm s}$}34$ and decl.(J2000) = −03°35'29 farcs 9) with the VLBA in the ²Π_3/2J = 3/2 OH transition (rest frequency 1.665401 GHz) on 2007 April 28 (program code: BM244O). During a run of ∼6 hr, we recorded the dual circular polarization with two IF bands of width of 1 MHz and 4 MHz, both centered at V_LSR = 103.0 km s⁻¹. The 4 MHz bandwidth was used to increase the S/N of the weak L-band signal of the continuum calibrator. The data were processed with the VLBA FX correlator in two correlation passes using 1024 and 128 spectral channels for the 1 MHz and 4 MHz bands, respectively. In each correlator pass, the data averaging time was 2 s.

Images were produced with natural weighting, cleaned and restored with an elliptical Gaussian with FWHM size of 17.8 mas × 10.9 mas at a P.A. of 1°. The interferometer instantaneous field of view was limited to ∼18 farcs 5. The on-source integration time was ∼1.9 hr, resulting in an effective rms noise level in each velocity channel of ∼0.02 Jy beam⁻¹. The spectral resolution was 0.2 km s⁻¹.

3.3. Subaru/COMICS: 24.5 μm Mid-infrared Continuum

Using the mid-infrared imaging spectrometer (COMICS; Kataza et al. 2000) at the Cassegrain focus of the 8.2 m Subaru Telescope, we carried out imaging observations of the 24.5 μm emission toward G28.87+0.07 on 2008 July 15. For this purpose, we configured the Q24.5 filter, and employed chopping mode for subtracting the sky-background emission. The camera provides a field of view of ∼42'' × 32'' with a pixel size of 0 farcs 13. Flux calibration was performed toward four standard sources listed in Cohen et al. (1999): HD146051, HD186791, HD198542, and HD198542. We estimated the overall uncertainty in the flux calibrations to be less than 10%. Astrometric calibration was performed by smoothing the 24.5 μm image to the same angular resolution as the 24 μm MIPSGAL (survey of the inner Galactic plane using MIPS; Carey et al. 2009) image, i.e., 6'', and overlaying the two images until a good match was found. Note that the heavy saturation of the MIPSGAL image does not hinder the comparison between it and the smoothed COMICS image. In practice, the emission in MIPSGAL is basically pointlike and the MIPS point spread function (PSF) presents a well-defined circular structure even far from the (saturated) peak. The match between the smoothed COMICS image and the MIPSGAL image can thus be performed successfully on such circular features. Inspection by eye indicates a positional error of less than 1'', comparable to the astrometric accuracy of MIPSGAL.

4. RESULTS AND ANALYSIS

This section reports our main results from the observation of the radio continuum, individual maser species, and IR emission toward the HMSFR G28.87+0.07. In the following, the term "spot" is used to refer to maser emission in a single velocity channel, whereas the term "feature" indicates a collection of spots emitting at a similar position over contiguous channels (i.e., an individual masing cloud). In this work, the parameters of maser emission (position, intensity, proper motion, etc.) are derived for maser features.

4.1. 1.3 cm and 3.6 cm Continuum

Our VLA measurements of the 1.3 and 3.6 cm continuum emission toward G28.87+0.07 (in both the A- and C-array configurations) are summarized in Table 1 and shown in Figure 1. At the angular resolution of the VLA C-array, the emission consists of two components, both detected at 1.3 and 3.6 cm. The northern and southern components, separated by ∼3 farcs 1, are labeled, respectively, "HMC" and "A" in Table 1. The former name is due to the source being associated with the HMC imaged by Furuya et al. (2008). Source "A" is partly resolved with the C-array at 1.3 cm, and is completely resolved out with the A-array at both 1.3 and 3.6 cm.

**Figure 1.** Maps (gray scale) of the VLA continuum emission toward G28.87+0.07 at 1.3 and 3.6 cm. (a) 3.6 cm emission with the VLA C-array. (b) 1.3 cm emission with the VLA C-array. (c) 3.6 cm emission with the VLA A-array. (d) 1.3 cm emission with the VLA A-array. Squares and circles report the positions of the 1.665 GHz OH and the 6.7 GHz CH₃OH masers, respectively, detected in this paper (see Sections 4.2.2 and 4.2.3). In each panel, the contour levels range from 30% to 90% of the peak emission at multiples of 10%. The restoring beam is shown in the lower left or right corner of each panel. A rectangle overlaid on the image indicates the field of view expanded in the adjacent panel. The associated labels indicate the names of the continuum sources identified by us (see Section 4.1 and Table 1).
Download figure:
Standard image High-resolution image

Table 1. Radio Continuum Emission in G28.87+0.07

Label	Array	λ	HPBW	P.A.	Image rms	Peak Position		F_peak	F_int
		(cm)	('' × '')	(°)	(mJy beam⁻¹)	R.A.(J2000)	Decl.(J2000)	(mJy beam⁻¹)	(mJy)
						(^h ^m ^s)	(° ' '')
HMC	VLA−A	1.3	0.122 × 0.093	−17	0.06	18 34 46.23	−03 35 29.8	0.62	0.64
	VLA−C	1.3	1.053 × 0.755	−9	0.05	18 43 46.23	−03 35 29.8	1.03	1.20
	VLA−A	3.6	0.33 × 0.24	−1	0.03	18 43 46.23	−03 35 29.9	0.35	0.42
	VLA−C	3.6	2.62 × 2.11	5	0.02	18 43 46.25	−03 35 29.8	0.53	0.55
A	VLA−A	1.3	0.122 × 0.093	−17	0.06	...	...	<0.18(3σ)	...
	VLA−C	1.3	1.053 × 0.755	−9	0.05	18 43 46.30	−03 35 32.8	0.39	0.60
	VLA−A	3.6	0.33 × 0.24	−1	0.03	...	...	<0.09(3σ)	...
	VLA−C	3.6	2.62 × 2.11	5	0.02	18 43 46.29	−03 35 32.7	0.81	0.98

Notes. Column 1 represents the source label; Columns 2 and 3 report the VLA array configuration and the observing band, respectively; Columns 4 and 5 report the (major and minor) HPBW size and the P.A. of the observing beam, respectively; Column 6 represents the image rms; Columns 7–10 represent the absolute position, the peak intensity, and the integral intensity of the emission peak, respectively. The P.A. of the beam is defined east of north.

Download table as: ASCII Typeset image

In the following, we focus our attention on the centimeter-continuum source HMC because this appears to be spatially associated with the observed 22 GHz H₂O, 6.7 GHz CH₃OH, and 1.665 GHz OH masers in G28.87+0.07 (see Figures 1 and 4). While at 1.3 cm the most extended VLA A-array configuration partly resolves the emission of HMC (missing about half of the flux measured with the C-array), at 3.6 cm the emission is essentially compact in both configurations (yielding comparable integral fluxes). Note, however, that the VLA A-array image at 3.6 cm shows a spur protruding to the south. The peak positions of the 1.3 and 3.6 cm emissions coincide within the uncertainties. The absolute position of the A-array 1.3 cm image has been obtained by aligning positions of persistent water maser spots between the VLA and VLBA observations, and should be accurate to ≈0 farcs 1.

Figure 2 shows the spectral energy distribution (SED) at centimeter wavelengths of HMC. A spectral index (α) of 0.6 is derived from a linear fit of the VLA C-array 1.3 and 3.6 cm and B-array 6 cm flux densities. The latter flux is an upper limit obtained from the GPS 6 cm Epoch 2 images downloadable from the MAGPIS Web site,⁸ which are the most sensitive among the other GPS maps and the relevant map from the CORNISH survey⁹ (see Purcell et al. 2008). We also fit the SED at centimeter wavelengths with a model of a homogeneous H ii region, with a Lyman continuum rate of 7.9 × 10⁴⁵ s⁻¹ and a radius of 0 farcs 04.

farcs — **Figure 2.** Spectral energy distribution of the VLA component HMC toward the HMSFR G28.87+0.07. Dots and error bars report the values and the associated errors (20%) of the fluxes at 1.3 and 3.6 cm measured with the VLA C-array (see Table 1). The triangle indicates the upper limit of the flux at 6 cm obtained with the VLA B-array (Purcell et al. 2008). The linear fit to these measurements, indicated by the dotted line, produces a spectral index α = 0.6. The solid line shows the best fit of the measured fluxes with a model of a homogeneous H ii region, yielding a Lyman continuum rate of 7.9 × 10⁴⁵ s⁻¹ and an H ii region radius of 004.
Download figure:
Standard image High-resolution image

We caution that the shape of the SED could be significantly affected by the different sensitivity to extended structures at different wavelengths. For the same VLA configuration, the shorter the wavelength, the larger the fraction of emission resolved out by the interferometer can be. Therefore, one should consider the possibility that the slope between 3.6 cm and 1.3 cm is steeper than 0.6.

4.2. Masers Results

4.2.1. 22.2 GHz Water Maser

Forster & Caswell (1999) observed the HMSFR G28.87+0.07 with the VLA C-array (half-power beam width (HPBW) of 2 farcs 2 × 1 farcs 6) and identified 13 water maser spots, mainly grouped in two clusters separated by ∼0 farcs 2 in the northeast–southwest direction. The absolute position of the brightest spot (with a flux density of 61.3 Jy at the V_LSR of 106.3 km s⁻¹) differs by ∼0 farcs 1 from that of the corresponding spot in our VLBA observation. The field of view of our VLBA observations covers the whole area within which water maser emission was detected by Forster & Caswell (1999).

Figure 3 (upper panel) shows the total-power spectrum of the 22.2 GHz masers toward G28.87+0.07. This profile was obtained by averaging the total-power spectra of all VLBA antennae, after weighting each spectrum with the antenna system temperature (T_sys). Water maser absolute positions and LSR velocities are plotted in Figure 4. Imaging the range of LSR velocities from 93.2 to 114.3 km s⁻¹, we detected 168 distinct water maser features distributed over a region of about 0 farcs 3 × 0 farcs 5. Most (92%) of them concentrate in a small area of about 0 farcs 3 × 0 farcs 2 centered on the VLA 1.3 cm continuum peak. The line-of-sight (l.o.s.) velocity distribution is remarkably bipolar, with redshifted features located to the west of the radio continuum and blueshifted ones to the east. While the blueshifted features mainly concentrate in two small (size <20 mas) clusters, the distribution of redshifted masers consists of one narrow strip elongated east–west (size of ≈100 mas), together with two diffuse "clouds" of features spread to the north and to the south of the densely populated, narrow strip. It is also interesting to note the presence of two streamlines of water features (enclosed in dashed boxes in the upper panel of Figure 4), oriented north–south and extending over a few tenths of an arcsecond, on either side of the VLA 1.3 cm continuum. We took particular care in filtering out spurious emission owing to poor cleaning of the beam side lobes (elongated mostly along the north–south direction) and can confirm that these water maser streamlines are real.

**Figure 4.** Spatial distribution of 22 GHz H₂O (triangles), 6.7 GHz CH₃OH (dots), and 1.665 GHz OH (squares) masers in G28.87+0.07. Different gray tones (colors in the online version) are used to indicate the maser LSR velocities, according to the gray scale on the right-hand side of the plot. Upper panel: global view of the region where the three maser species distribute. The solid contours give the 3.6 cm continuum emission observed with the VLA A-array. Contour levels range from 20% to 90% of the peak emission (0.35 mJy beam⁻¹) at multiples of 10%. Dashed rectangles enclose the two north–south streamlines of water masers observed on either side of the VLA 3.6 cm continuum. A dotted rectangle overlaid on the image indicates the field of view expanded in the right panel. Lower panel: zoom on the area where most of water masers concentrate. The solid contours give the 1.3 cm continuum emission observed with the VLA A-array. Contour levels range from 40% to 90% of the peak emission (0.62 mJy beam⁻¹) at multiples of 10%.
Download figure:
Standard image High-resolution image

The individual feature properties are presented in Table 3. The position offsets are relative to feature 1. The absolute positions of feature 1 at the four observing epochs are given in Table 2. Feature intensities range from 0.1 to 54 Jy. The spread in LSR velocities ranges from 111.99 km s⁻¹ for the most redshifted feature (No. 107) to 93.33 km s⁻¹ for the most blueshifted one (No. 155).

Table 2. Positions and Brightness

Maser	Feature	Epoch	R.A.(J2000)	Decl.(J2000)	F_peak	NW Beam
			(^h ^m ^s)	(° ' '')	(Jy beam⁻¹)	(mas × mas @ deg)
H₂O	1	2006 Apr 23	18 43 46.22467	−03 35 29.8161	57.4	1.2 × 0.4 @ −17
		2006 Jun 30	18 43 46.22465	−03 35 29.8173	50.3	1.3 × 0.4 @ −16
		2006 Sep 28	18 43 46.22455	−03 35 29.8190	32.0	1.4 × 0.5 @ −3
		2007 Jan 18	18 43 46.22448	−03 35 29.8212	20.7	1.5 × 1.0 @ 20
CH₃OH	1	2006 Feb 28	18 43 46.22079	−03 35 29.8379	0.044	9.9 × 5.6 @ 5
		2007 Mar 18	...	...	...	...
		2008 Mar 18	18 43 46.22031	−03 35 29.8521	0.246	12.7 × 5.2 @ 45
		2009 Mar 17	18 43 46.22011	−03 35 29.8595	0.175	12.6 × 5.9 @ 13
OH	1	2007 Apr 28	18 43 46.26771	−03 35 30.4697	1.94	17.8 × 10.9 @ 1

Notes. Column 1 reports the maser species; Column 2 gives the label number of the feature, as given in Tables 4–6; Column 3 lists the observing date; Columns 4 and 5 report the absolute position; Column 6 gives the peak brightness; Column 7 reports the (major and minor) HPBW size and P.A. of the naturally weighted (NW) beam. The P.A. of the beam is defined East of North.

Download table as: ASCII Typeset image

Table 3. Parameters of VLBA 22.2 GHz H₂O Maser Features in G28.87+0.07

Feature	Epochs of	V_LSR	S_ν	Δα	Δδ	V_α	V_δ
Number	Detection	(km s⁻¹)	(Jy)	(mas)	(mas)	(km s⁻¹)	(km s⁻¹)
0	1,2,3,4	103.17	...	66.83 ± 0.05	0.005 ± 0.05	0	0
1	1,2,3,4	106.01	53.99	0	0	−16.9 ± 1.0	1.4 ± 1.0
2	1,2,3,4	105.07	26.07	37.44 ± 0.05	−3.02 ± 0.05	−23.0 ± 1.2	−26.4 ± 1.2
3	1	104.36	25.38	−2.88 ± 0.05	2.60 ± 0.05	...	...
4	1,2,3	105.27	19.92	−0.85 ± 0.05	0.25 ± 0.05	−2.8 ± 2.2	0.2 ± 2.3
5	2,3,4	104.27	17.62	−3.29 ± 0.05	2.79 ± 0.05	−29.8 ± 1.7	0.2 ± 1.8
...

Notes. Column 1 gives the feature label number; Column 2 lists the observing epochs at which the feature was detected; Columns 3 and 4 report the value of the intensity-weighted LSR velocity and flux density of the strongest spot, averaged over the observing epochs; Columns 5 and 6 give the position offsets (with the associated errors) along the R.A. and decl. axes relative to the feature 1, measured at the first epoch of detection; Columns 7 and 8 give the components of the relative proper motion (with the associated errors) along the R.A. and decl. axes, measured with respect to the reference feature 0.

Only a portion of this table is shown here to demonstrate its form and content. A machine-readable version of the full table is available.

Download table as: Data Typeset image

Only 61 features (36% of the total) persisted over at least three epochs, 34 of which lasting four epochs. We calculated the geometric center (hereafter "center of motion," identified with label 0 in Table 4, and indicated by a star in Figure 5) of features with a stable spatial and spectral structure, persisting over the four observing epochs, and refer our measurement of proper motions to this point. With a time baseline of nine months, relative (sky-projected) velocities of water maser features are derived with a mean uncertainty of ∼23%. The magnitude of relative proper motions ranges from 1.9 ± 2.8 km s⁻¹ to 58.3 ± 2.6 km s⁻¹, with a mean value of 24.0 km s⁻¹. Figure 5 shows the derived water maser proper motions. While the two blueshifted clusters move to the east and northeast, most velocities of the redshifted water masers point either to the south or to the west. In particular, looking at the distribution of the maser velocity vectors along the (east–west) narrow strip outlined by the redshifted features, we note that maser features closer to the "center of motion" move mainly to the south. Getting closer to the western end of the maser strip, the orientation of the maser velocity vectors changes gradually from south to west. The general pattern of water maser relative velocities appears to indicate expansion from a point on the sky close to the "center of motion." Note also that at larger distances from the center of motion, maser velocities seem to be larger and better collimated along the (north)east–(south)west direction. In particular, we identify two H₂O maser clusters located at the northeast and southwest ends of the maser distribution (named "J_b" and "J_r," respectively; see Figure 5), which move with high velocities directed close to the northeast–southwest direction.

**Figure 5.** Proper motions of 22 GHz H₂O masers (triangles) relative to their centers of motion, and spatial distribution of 6.7 GHz CH₃OH masers (dots) in G28.87+0.07. Different gray tones (colors in the online version) are used to indicate the maser LSR velocities, according to the gray scale on the right-hand side of the plot. The cones indicate the 3D velocities of the water maser features relative to the center of motion (marked with a cross). Points without an associated cone have been detected only over one or two epochs and the associated proper motion cannot be computed or is considered unreliable. The cone opening angle gives the 1σ uncertainty on the proper motion direction. The length of the cone is proportional to the velocity, with the amplitude scale indicated in the lower left corner of the figure. Open dotted rectangles mark two clusters of water masers located at the northeast and southwest edge of the distribution (labeled "*J_b*" and "*J_r*," respectively), moving fast and close to the northeast–southwest direction. The VLA 1.3 cm continuum emission is plotted with solid contours. Contour levels range from 40% to 90% of the peak emission (0.62 mJy beam⁻¹) at multiples of 10%. Upper panel: view of the whole distribution of water maser proper motions. Lower panel: zoom on the most crowded area of water maser concentration, to the west of the VLA 1.3 cm continuum.
Download figure:
Standard image High-resolution image

Table 4. Parameters of EVN 6.7 GHz CH₃OH Maser Features in G28.87+0.07

Feature	Epochs of	V_LSR	F_peak	I	Δα	Δδ	v_α	v_δ
	Detection	(km s⁻¹)	(Jy beam⁻¹)	(Jy)	(mas)	(mas)	(km s⁻¹)	(km s⁻¹)
1	1, 3, 4	104.88	0.155	0.190	0	0	1.7 ± 6.7	−15.6 ± 8.1
2	3, 4	105.21	0.268	0.283	−12.41 ± 0.18	−11.24 ± 0.18	...	...

Notes. Column 1 gives the feature label number; Column 2 lists the observing epochs at which the feature was detected; Columns 3 and 4 provide the intensity-weighted LSR velocity and the flux density of the strongest spot, respectively, averaged over the observing epochs; Columns 5 and 6 give the position offsets (with the associated errors) along the R.A. and decl. axes, relative to the feature 1; Columns 7 and 8 give the components of the relative proper motion (with the associated errors) along the R.A. and decl. axes, measured with respect to the reference feature 0 (the center of motion) of water masers.

Download table as: ASCII Typeset image

Figure 6 (upper panel) reports the intrinsic absolute proper motions of the water maser features, derived from the observed absolute proper motion of feature 1 corrected for the apparent motion. The latter has been calculated adding the contribution of the parallax, the solar motion with respect to the LSR (Schönrich et al. 2010), and the differential Galactic rotation, from a flat rotation curve with Galactic constants (R₀ = 8.3 ± 0.23 kpc, Θ₀ = 239 ± 7 km s⁻¹) as recently determined by Brunthaler et al. (2011). The obtained absolute proper motions differ from the relative proper motions presented in Figure 5 by a vector of amplitude ≈40 km s⁻¹ pointing to the southwest, which seems to be the dominant component of all proper motions. We think that this difference could be caused by the large systematic error due to the uncertainties in the source distance and solar and Galactic motion, and/or by a peculiar velocity of the maser source with respect to its LSR reference system. Our VLBA water maser observations were not designed to measure the source parallax. With a positional accuracy of ≈0.1–0.2 mas, and four epochs spanning (only) nine months, the observed absolute proper motion of feature 1, not corrected for the apparent motion, looks similar to a straight line. From the small deviations from linear motion, we derive an upper limit of 0.2 mas for the amplitude of the sinusoidal parallax signature. This corresponds to a lower limit of 5 kpc for the source distance, in agreement with the adopted kinematical distance of 7.4 kpc.

We have checked the reliability of the relative, water maser proper motions by following an alternative approach. As described in Section 4.2.2, we have also derived the absolute proper motion of the strongest 6.7 GHz methanol maser feature. To date, VLBI observations of 6.7 GHz (Sanna et al. 2010a, 2010b) and 12 GHz (Moscadelli et al. 2003, 2010, 2011) methanol masers toward several sources indicate that these masers usually show small proper motions, with typical velocities of only a few km s⁻¹. Assuming that also in G28.87+0.07 methanol masers move significantly slower than the water masers, we can use the absolute proper motion of the methanol masers to correct the absolute proper motion of the water masers, and derive the water maser velocities as approximately seen in a reference system comoving with the star. Figure 6 (lower panel) shows the absolute water maser proper motions after subtracting the absolute methanol maser proper motion (evaluated in Section 4.2.2). Note that, within the measurement errors, the new absolute proper motions are fully consistent with the relative proper motions calculated with respect to the center of motion (Figure 5). This result reinforces our assertion that the center of motion can define a suitable reference system to calculate water maser velocities. In the following section, we use the relative velocities of the water masers to describe the gas kinematics close to the (proto)star.

4.2.2. 6.7 GHz Methanol Maser

Figure 3 (middle panel) shows the total-power spectrum of the 6.7 GHz methanol masers observed on 2008 March 18 toward G28.87+0.07 using the Effelsberg 100 m telescope. The 6.7 GHz methanol maser spectrum, with two narrow features at V_LSR of about 97 and 105 km s⁻¹, resembles that of the water masers (Figure 3, upper panel). Note that Caswell et al. (1995) detected only the emission at V_LSR of 105 km s⁻¹ with the Parkes 64 m telescope.

We imaged the range of LSR velocities from 95 to 107 km s⁻¹ with the EVN and detected two distinct methanol maser features, both with V_LSR close to 105 km s⁻¹. Feature 1 is detected at the first, third, and fourth epoch, whereas feature 2 is detected at the third and fourth epoch only. Due to adverse atmospherical conditions and telescope failures, no 6.7 GHz maser signal was detected at the second epoch. Individual feature properties are presented in Table 4. Positional offsets are relative to feature 1, and the absolute positions of feature 1 at the three epochs of detection are given in Table 2. The two methanol maser features lie at the western border of the distribution of the redshifted water masers, ∼0 farcs 2 to the southwest of the VLA 3.6 cm continuum peak (see Figure 4). Inspecting our images, we could not find any reliable emission at the velocity of ∼97 km s⁻¹, in correspondence of the weakest feature detected in the total-power spectrum. The mean rms noise on the channel maps was 20 mJy beam⁻¹, and we could have missed that emission if falling below the 5σ detection threshold. Because of the narrow-line width (⩽1 km s⁻¹) and the offset in V_LSR from the systemic velocity (101.0 km s⁻¹) of the HMC, we tend to exclude that the weak methanol feature at V_LSR ∼97 km s⁻¹ is thermal emission from the HMC resolved out by the EVN observations.

After correcting for the apparent motion, the absolute proper motion for feature 1 is: V_α = −27.8 ± 1.0, V_δ = −45.7 ± 4.6 km s⁻¹. As discussed in Section 4.2.1, we have used this vector to correct the absolute proper motions of the water masers. The relative motion of the methanol maser feature 1 with respect to the center of motion of the water masers is found to be: v_α = 1.7 ± 6.7, v_δ = −15.6 ± 8.0 km s⁻¹. This is consistent with the results from other VLBI observations which show that 6.7 GHz masers move significantly slower than water masers, and justifies employing the methanol maser absolute motion to correct the water maser absolute motions.

4.2.3. 1.665 GHz OH Maser

Forster & Caswell (1999) observed the HMSFR G28.87+0.07 with the VLA A-B hybrid configuration (HPBW of 5 farcs 5 × 1 farcs 2) in the left-circular polarization (LCP) band and detected eight 1.665 GHz OH maser spots. Their brightest spot had a flux density of 9.85 Jy at the V_LSR of 102.8 km s⁻¹.

Figure 3 (lower panel) shows the T_sys-weighted mean of the VLBA antenna total-power spectra of the 1.665 GHz OH maser emission in the LCP and right-circular polarization (RCP) bands. Three strong, narrow features are visible in the LCP band over the LSR velocity range from 102 to 108 km s⁻¹. Two weak features are detected in the RCP band, one of them at an LSR velocity of ∼108 km s⁻¹. These findings are similar to the results obtained by Caswell & Haynes (1983) with the Parkes 64 m telescope.

We imaged the range of LSR velocities from 102 to 110 km s⁻¹, and identified three distinct OH maser features in the LCP band. In the RCP band, by imaging maser data after direct phase-referencing to the calibrator J1834−0301, no signal was detected above a 5σ level of 0.1 Jy beam⁻¹. Individual feature properties are presented in Table 5. The positional offsets are relative to feature 1, whose absolute position is given in Table 2. The derived absolute position is offset by ∼0 farcs 9 from the reference position of Forster & Caswell (1999). Maser intensities range from 0.7 to 1.9 Jy beam⁻¹. The LSR velocities vary from 102.8 km s⁻¹ for the most blueshifted feature (No. 1) to 105.5 km s⁻¹ for the most redshifted feature (No. 3). The FWHM line width of individual maser features ranges from 0.4 to 1.4 km s⁻¹. Positions and LSR velocities of hydroxy1 maser features are plotted in Figure 4. The three detected features lie ≈0 farcs 6 to the southeast of the VLA 3.6 cm continuum source and are distributed along a line extended ≈0 farcs 3.

Table 5. Parameters of VLBA 1.665 GHz OH Maser Features in G28.87+0.07

Feature	Pol	V_LSR	F_peak	I	Δx	Δy
		(km s⁻¹)	(Jy beam⁻¹)	(Jy)	(mas)	(mas)
1	L	102.9	1.94	5.26	0	0
2	L	106.1	1.44	3.52	−233.57 ± 0.11	143.92 ± 0.13
3	L	105.4	0.67	1.78	−139.46 ± 0.17	46.91 ± 0.21

Notes. Column 1 gives the feature label number; Column 2 indicates the circular polarization of the maser emission; Columns 3–5 report the intensity-weighted LSR velocity, and the intensity and the integral flux of the strongest spot, respectively; Columns 6 and 7 list the position offsets (with the associated errors) along the R.A. and decl. axes, relative to the feature 1.

Download table as: ASCII Typeset image

4.3. 24.5 μm Infrared Emission

Figure 7 shows an image of the 24.5 μm continuum emission obtained with Subaru/COMICS at an angular resolution of 0 farcs 75. One can see that the structure of the emitting region is quite complex. Besides the two components associated with the radio sources HMC and A with a total flux of ∼80 Jy, another fainter source (∼3 Jy) is visible to the south, offset by ∼4'' from A.

The emission from the northern region is clearly resolved and appears to be roughly associated with the two radio objects, although one cannot rule out the existence of additional, fainter unresolved sources contributing to the extended halo around the two main sources. In particular, the 24.5 μm emission from HMC seems to split into two sub-components in the E–W direction, with comparable flux densities: ∼32 Jy for the eastern and ∼43 Jy for the western component. These could trace two distinct objects or diffuse emission surrounding a single source at the origin of the HMC radio continuum. A third, more intriguing possibility is that one is looking at an "hour-glass" shaped nebulosity centered on the embedded source at the origin of the free–free emission.

In the latter scenario, the source lies at the center of a flattened core lying approximately along the l.o.s. This dusty core absorbs most of the IR photons emitted at 24 μm toward the observer, while the IR radiation can escape along the minor axis of the core, lying close to the plane of the sky. We have estimated an optical depth of ∼4 for the HMC at 24 μm, assuming the H₂ column density of 4.2 × 10²³ cm⁻² obtained by Furuya et al. (2008) and the opacity adopted by Whitney et al. (2003a, see their Figure 1). Such a value could be sufficiently large to explain the lack of emission toward the embedded star and determine the "hour-glass" shape of the IR emission. Interestingly, the E–W orientation of this putative bipolar nebula is also roughly consistent with the direction of the H₂O maser jet in Figure 5, suggesting that one could be observing a bipolar flow piercing through a dense core, and revealed on different scales in different tracers. However, such a scenario, albeit appealing, is not consistent with the inclination of the maser jet with respect to the plane of the sky. In fact, the blueshifted emission lies to the E (see Figure 5) and this implies that the eastern lobe should be pointing toward the observer and thus be brighter than the western lobe in the IR. Looking at Figure 7 one can see that this is not the case, being the emission slightly stronger to the west (see above).

In conclusion, we consider both the "hour-glass" nebula and the double source hypotheses too speculative on the basis of the observational evidence collected so far, and in the following we make the conservative assumption that the IR emission is associated with HMC, with a flux of ∼75 Jy.

The southern source A has a 24.5 μm flux density of ∼5 Jy. Also in this case, the association between the IR emission and the radio source is questionable because the two are offset each other by ∼1''. However, this is the astrometrical error of the IR image and we thus believe that such an offset is not sufficient evidence against a common origin for the IR and radio emission toward A. Therefore, we conclude that both emissions can originate from the same object.

4.4. Spectral Energy Distribution

Figure 8 shows the SED of G28.87+0.07 from 3.6 μm to 1.1 mm. Besides our COMICS data, we have used data from the following surveys: Galactic Legacy Infrared Mid-Plane Survey Extraordinaire (GLIMPSE; Benjamin et al. 2003), MIPSGAL (Carey et al. 2009), MSX (Egan et al. 2003), IRAS (Neugebauer et al. 1984), the Herschel infrared Galactic Plane Survey (Hi-GAL; Molinari et al. 2010a, 2010b), the APEX Telescope Large Area Survey of the Galaxy (ATLASGAL; Schuller et al. 2009), and the Bolocam Galactic Plane Survey (BGPS; Drosback et al. 2008). For the Herschel/Hi-GAL data, the most recent release (obtained in 2011 July with HIPE version 7.0.0) was used, resulting from data reduction with the RomaGal software (Traficante et al. 2011), where the most relevant image artifacts have been removed.

Apart from the IRAS and MSX flux densities, which have been taken from the corresponding Point Source Catalogues (sources IRAS 18411−0338 and G028.8621+00.0657), all the other values have been calculated by integrating the emission inside the dashed contour in Figure 9. The latter was chosen in such a way to include the whole emitting region at all wavelengths. In this way we overcome the problem of different angular resolutions at different wavelengths. In the calculation of the flux the background, emission was subtracted assuming that the background flux inside the contour is in all points equal to the typical flux density measured along the contour itself. The flux densities used for the SED are given in Table 6.

**Figure 9.** Each panel shows the pc-scale structure of the G28.87+0.07 star-forming region at a different infrared wavelength (indicated in the upper left corner of the panel). Images from 70 to 500 μm are observed with Herschel, and the map at 870 μm is taken with APEX. The gray scale corresponds to increasing intensity from white to black, as indicated in the wedge to the right. The observing beam is reported in the lower right corner of each panel. The dashed black circle denotes the area of integration to calculate the G28.87+0.07 flux at different infrared wavelengths.
Download figure:
Standard image High-resolution image

Table 6. Flux Densities Measured toward the G28.87+0.07 Region

λ	F_ν	Instrument
(μm)	(Jy)
3.6	0.9	Spitzer/GLIMPSE
4.5	2.1	Spitzer/GLIMPSE
5.8	5.3	Spitzer/GLIMPSE
8.0	9.4	Spitzer/GLIMPSE
8.3	3.3	MSX
12	8.6^a	IRAS
12.1	7.1	MSX
14.7	15.3	MSX
21.3	34.0	MSX
24	>21	Spitzer/MIPSGAL
24.6	87.4	Subaru/COMICS
25	105^a	IRAS
60	1720^a	IRAS
70	1790	Herschel/Hi-GAL
100	<3180	IRAS
160	1350	Herschel/Hi-GAL
250	630	Herschel/Hi-GAL
350	200	Herschel/Hi-GAL
500	53	Herschel/Hi-GAL
870	8.8	APEX/ATLASGAL
1100	2.5	CSO/BGPS

Note. ^aAll IRAS fluxes are assumed to be upper limits in Figure 8.

Download table as: ASCII Typeset image

Figure 9 illustrates the pc-scale structure of the G28.87+0.07 star-forming region, imaged with Herschel from 70 to 500 μm and with APEX at 870 μm, and demonstrates that from the mid-IR to the submillimeter the continuum emission is dominated by a compact source surrounded by a weaker halo.

We have fitted the SED with the radiative transfer model developed by Robitaille et al. (2007) and Whitney et al. (2003a, 2003b; using the SED fitting tool available at http://caravan.astro.wisc.edu/protostars/), which assumes a zero-age main sequence (ZAMS) star with a circumstellar disk, embedded in an infalling envelope, and allows derivation of a number of physical parameters. While we are aware that this model cannot adequately reproduce a complex situation such as that of the multiple MYSOs in G28.87+0.07, we believe that the fit is sufficiently good to guarantee good estimates of integrated quantities such as the luminosity and the mass of the envelope. The best-fit values are, respectively, 2.3 × 10⁵ L_☉ and 3.0 × 10³ M_☉. Note that in the fitting procedure, the IRAS fluxes were considered upper limits due to the limited angular resolution, whereas the MIPSGAL 24 μm flux is in fact a lower limit because the corresponding image is saturated and only the non-saturated pixels have been considered in the calculation. Also, the distance was kept fixed to the nominal value of 7.4 kpc. The best-fit value for the visual foreground extinction is 46 mag.

5. DISCUSSION

5.1. Nature of the Radio Continuum

In this section, we focus on the centimeter-continuum source labeled HMC in Table 1, found in good correspondence with the positions of the 22 GHz H₂O, 6.7 GHz CH₃OH, and 1.665 GHz OH masers in G28.87+0.07 (see Figure 4). The continuum spectrum (shown in Figure 2) is consistent with free–free emission from an ultracompact (UC) H ii region, with the turnover frequency in the range 10–20 GHz, and ionized by a Lyman continuum rate of 7.9 × 10⁴⁵ s⁻¹. This would correspond to a ZAMS star of spectral type B1–B0.5 (Panagia 1973), with a luminosity of L_star ∼ 8 × 10³ L_☉, and a mass of M_star ∼ 10 M_☉. It is worth mentioning that these values would be different if more recent studies than the seminal article by Panagia were adopted (e.g., Crowther 2005). However, the basic result that one is dealing with an early B-type star with a luminosity of ∼10⁴ L_☉ would remain valid.

We note that the fit to the continuum spectrum requires a very small size of the UC H ii region, of ∼0 farcs 04 or ∼300 AU, which appears in contradiction with the fact that our A-array 1.3 cm observations partly resolve the continuum emission (see Table 1).

Since we considered a homogeneous distribution of ionized gas in our fit, one possibility is that the H ii region is not homogeneous. The southern spur identified in the A-array 3.6 cm map of the HMC source (see Figure 4) could indeed suggest the existence of a density gradient toward the south. One may speculate that source A being more extended and perhaps more evolved, the gas in its surroundings is less dense. In this scenario, the north–south streamlines of water masers on either sides of the southern spur might trace shocks at the interface between the ionized gas of HMC (expanding to the south) and the surrounding molecular gas.

We also note that the continuum spectrum of HMC between 1.3 and 6 cm can also be fitted with a power law S_ν∝ν^α, with α = 0.6 (see Figure 2). This spectral index is typical of thermal jets (Reynolds 1986) and provides us with an alternative explanation for the origin of the free–free emission. However, a caveat is in order. As discussed in Section 4.1, the slope of the spectrum between 3.6 cm and 1.3 cm could be significantly affected by the interferometer resolving out part of the emission at the shorter wavelength. If the spectrum is steeper than estimated from our data, the slope could be still consistent with a thermal jet, but also with that of an H ii region with a density gradient (see, e.g., Franco et al. 2000).

We conclude that it is impossible to decide the nature of the free–free emission on the basis of the information collected so far. Notwithstanding this caveat, as a working hypothesis we prefer the jet interpretation, because this permits us to interpret the continuum and maser emission in the same context, as discussed in the following section.

5.2. The Water Maser Flow

In this section we propose that both the continuum and the water maser emission are associated with a thermal jet impinging against the surrounding quiescent molecular gas. This idea is based on two observational findings. First of all, Figure 5 shows that, at larger distances from the center of the water maser distribution, maser (relative) velocities appear to be larger and better collimated along the northeast–southwest direction. Second, the axis of the powerful, molecular outflow mapped at arcsecond scales toward G28.87+0.07 has a northeast–southwest orientation (Furuya et al. 2008, their Figure 5). This it is plausible to assume that the two, fast moving, water maser clusters J_b and J_r, located, respectively, at the northeast and southwest edge of the maser distribution, are excited by the same jet also responsible for driving the large-scale molecular outflow. One could object against the "jet" interpretation that the whole pattern of water maser velocities cannot be explained only in terms of a collimated outflow. In fact, at smaller distances from the center of the maser distribution, water maser velocities point to the north or south directions, and form large angles with the direction of the putative jet. We propose that the scattered pattern of proper motions presented by the redshifted water masers can result from the interaction of the jet with dense circumstellar material, which deflects and increases the opening angle of the flow velocities.

If the free–free continuum source HMC is a jet, the position of the (proto)star powering it may be offset from the free–free continuum peak. If a jet powered by the (proto)star drives the water maser motions, and, in particular, the expansion to the south and west of the narrow strip of redshifted water masers, the (proto)star could be located at an intermediate position between HMC and the redshifted maser strip. This point is not far (within tens of mas) from the position of the center of motion of the water masers (cross in Figure 5). Considering the large uncertainties in deriving the center of water maser expansion, in the following discussion we can safely take the center of motion as a reasonable estimate of the (proto)star location.

In the jet scenario, one can estimate the momentum rate of the jet from the average distance of the H₂O masers from the MYSO and the average maser velocity. Assuming an H₂ pre-shock density of 10⁸ cm⁻³, as predicted by models (Elitzur et al. 1989), the momentum rate in the water maser jet is given by the expression $\dot{P} = 1.5 \times 10^{-3}\,V_{10}^{2}\,R_{100}^{2}(\Omega /4\pi)$ M_☉ yr⁻¹ km s⁻¹. Here V₁₀ is the average maser velocity in units of 10 km s⁻¹, R₁₀₀ is the average distance of water masers from the YSO in units of 100 AU, and Ω is the solid angle of the jet. This expression has been obtained by multiplying the momentum rate per unit surface transferred to the ambient gas ( $n_{\rm H_{2}}\, m_{\rm H_{2}}\, V^{2}$ ), by ΩR², under the assumption that the jet is emitted from a source at a distance R from the masers within a beaming angle Ω. To estimate the average distance and velocity of the water masers tracing the jet in G28.87+0.07, we use all the maser features belonging to the two fast moving clusters J_b and J_r (see Figure 5). The resulting average distance is 122 mas (or 902 AU) and the average speed is 43 km s⁻¹. Using these values, the jet momentum rate is $\dot{P} = 2.2 \, (\Omega /4\pi)$ M_☉ yr⁻¹ km s⁻¹.

For shock-induced ionization in a thermal jet, the momentum rate of the jet can be estimated from the continuum flux assuming that the emission is optically thin (see Anglada et al. 1996). One has: $F_\nu d^{2} = 10^{3.5}(\Omega /4\pi)\dot{P}$ , where F_ν is the measured continuum flux density in mJy, $\dot{P}$ is the jet momentum rate in M_☉ yr⁻¹ km s⁻¹, Ω is the jet solid angle in sterad, and d is the source distance in kpc. Using the flux of 1.2 mJy measured at 1.3 cm with the VLA C-array for source HMC, and the distance of 7.4 kpc, one derives $\dot{P} = 2 \times 10^{-2} \,(\Omega /4\pi)^{-1}$ M_☉ yr⁻¹ km s⁻¹.

Note that the momentum rates estimated from the maser motion and the continuum flux have a different dependency on the jet solid angle, which can thus be estimated by equating the two momentum rate estimates.

We obtain an expression for the jet solid angle: 2.2 (Ω/4π) = 2 × 10⁻² (Ω/4π)⁻¹. We find Ω = 1.2 sterad, i.e., a jet semi-opening angle θ_j = 35°. Correspondingly, the jet momentum rate is $\dot{P} = 0.2 \, M_{\odot }$ yr⁻¹ km s⁻¹. This value is large enough to justify the momentum rate of several 10⁻² M_☉ yr⁻¹ km s⁻¹, estimated by Furuya et al. (2008) for the large-scale molecular outflow. This result supports our hypothesis that the small-scale thermal and H₂O maser jet is sufficiently powerful to feed the large-scale molecular outflow observed by Furuya et al. (2008). We note that the derived parameters for the jet in G28.87+0.07 are comparable to those of the few collimated outflows observed up to now in massive star-forming regions. For comparison, the jets in HH 80–81 (Marti et al. 1998), in IRAS 16547−4247 (Rodríguez et al. 2008), and in IRAS 20126+4104 (Moscadelli et al. 2011) have momentum rates of the order of 10⁻¹ M_☉ yr⁻¹ km s⁻¹ and opening angles of 20°–30°.

5.3. Methanol Maser Environment

The two detected 6.7 GHz methanol maser features are approximately aligned along the same northeast–southwest line defined by the two maser jet clusters J_b and J_r, the peak of the VLA A-array 1.3 cm, and the center of motion (see Figure 4). Figure 3 shows that the total-power spectra of the 6.7 GHz methanol and 22 GHz water masers are similar, both presenting a strong redshifted and a much weaker blueshifted component. The positions of the methanol maser features roughly along the jet axis, and the similarities in the spectral shape of the water and methanol maser emissions, make us speculate that the 6.7 GHz methanol masers could also originate in the (proto)stellar jet, albeit associated with more quiescent gas than that traced by the water masers. Maser VLBI observations toward the sources IRAS 20126+4104 and G23.01−0.41 (Moscadelli et al. 2011; Sanna et al. 2010b) have also indicated that the 6.7 GHz masers might trace entrained gas at the interface with a jet.

Figure 4 clearly shows that the blueshifted lobe of the water maser jet is sampled much more sparsely, as it presents fewer features than the redshifted lobe. This is consistent with the blueshifted line being significantly weaker than the redshifted one in the total-power spectrum (see Figure 3, upper panel). This fact suggests that the circumstellar material is not isotropically distributed around the (proto)star, but it is mostly concentrated on the side of the redshifted lobe. That would also account for the stronger emission line in the total-power spectrum of the CH₃OH masers at the redshifted velocities.

5.4. OH Maser Environment

The spectrum of the 1.665 GHz maser is dominated by three narrow LCP features (see Figure 3, lower panel), and is very similar to the spectrum of the 1.667 GHz maser transition observed by Caswell & Haynes (1983). These authors find that most OH interstellar masers are highly polarized, either RCP or LCP, with the same predominant polarization for both the 1.665 GHz and 1.667 GHz transition. They propose an explanation based on large-scale velocity and/or magnetic field gradients across the masing region.

Figure 1 shows that the 1.665 GHz OH maser features are found significantly offset (≈0 farcs 6) from the compact VLA source HMC and the center of water maser activity, and are closer to the more extended source A. A possible interpretation is that OH masers trace more external, less dense parts of the molecular core surrounding the (proto)star responsible for the excitation of the water and methanol masers. Observations indicate that the density of massive protostellar cores decreases from the center outward, with a power-law index of about 1.6 (n∝r^−1.6) (Beuther et al. 2002a). In G28.87+0.07, the ratio of the average distance (from the center of motion of water masers) of the water (≈0 farcs 1) and OH (≈0 farcs 6) masers, is about 6, implying a decrease in gas density by a factor of about 20. That agrees with maser excitation models, which predict an H₂ pre-shock density of 10⁶–10⁷ cm⁻³ for OH 1.665 GHz (Cesaroni & Walmsley 1991; Cragg et al. 2002), and 10⁸ cm⁻³ for H₂O 22 GHz masers (Elitzur et al. 1989), respectively.

5.5. The Nature of Source HMC

So far, we have explored the properties of source HMC making use of the maser, radio, and IR continuum emission. Now, we want to use our findings to shed light on the nature of the MYSO(s) associated with this source.

We have seen that the radio continuum emission and the water maser features suggest the presence of a jet with a large momentum rate, on the order of ∼10⁻¹ M_☉ km s⁻¹ yr⁻¹. Values that high are consistent with the jet being powered by a high-mass (proto)stars, as indicated, e.g., by the relationship between outflow momentum rate and MYSO luminosity obtained by Beuther et al. (2002b). In particular, from their Figure 4(b), one can conclude that the MYSO must have a luminosity of at least several 10³ L_☉. Alternatively, if the radio emission were originating from an H ii region, we have seen (Section 5.1) that the ionizing star should have a luminosity ≲ 10⁴ L_☉. Clearly, both the thermal jet and the H ii region hypotheses imply similar luminosity for the associated MYSO. The question we want to investigate here is whether such a luminosity is consistent with the results obtained from the IR data.

In Section 4.4 we have found that the total luminosity of the region is on the order of 2 × 10⁵ L_☉. However, this includes all MYSOs in a radius ≳10'' around HMC, whereas we need to estimate the contribution of the sole HMC, to achieve a consistent comparison with the luminosity estimate obtained from the water maser and/or radio continuum. Our 24.5 μm continuum image can be used for this purpose. As shown in Figure 7, the HMC source is dominating the emission at this wavelength, contributing to ∼90% of the total flux density, and a similar conclusion holds also for the GLIMPSE images. Due to the lack of subarcsecond imaging at longer wavelengths, it is difficult to establish whether also the flux density measured close to the peak of the SED is mostly due to source HMC. Nonetheless, the Hi-GAL/Herschel data suggest that the 70 μm emission peaks closer to source HMC than source A. A word of caution is in order here. The angular resolution of the Hi-GAL data at 70 μm is estimated to be about 10'' × 9'' (after reduction with the RomaGal software), while the astrometrical error is ∼3 farcs 5 (the latter results from aligning the PACS 70 μm images to the MIPSGAL 24 μm images by matching suitable point sources detected at both wavelengths). Although large, these values still permit to associate the 70 μm emission with source HMC. In fact, the angular separation between the 3.6 cm continuum source HMC and the 70 μm peak is 3 farcs 7, consistent with the positional error, whereas source A lies 6 farcs 5 from the 70 μm peak. This can be appreciated in Figure 10, where one sees that the 70 μm emission peaks to the NW of HMC. Therefore, in the following we assume that the bolometric luminosity estimated from the SED is likely dominated by the HMC source and conclude that the MYSO(s) embedded in the HMC may be as luminous as ∼2 × 10⁵ L_☉.

**Figure 10.** Overlay of the 3.6 cm continuum map (contours; same as in Figure 7) on the Hi-GAL/Herschel image at 70 μm. Note how the far-IR emission peaks close to the free–free source HMC.
Download figure:
Standard image High-resolution image

Notwithstanding the large uncertainties of all our luminosity estimates, the value obtained from the IR images exceeds by an order of magnitude that derived from the maser and radio continuum data. This inconsistency can be explained if HMC is hosting a multiple stellar system rather than a single star. One can get an idea of how this may affect our estimate of the stellar properties (mass, etc.) by assuming that the luminosity of 2 × 10⁵ L_☉ is emitted by a stellar cluster with a Miller & Scalo (1979) initial mass function. One finds that the luminosity of the most massive star may be as small as ∼3 × 10⁴ L_☉, depending on the statistical distribution of the stellar masses in the cluster. Given the uncertainties, this value is only marginally greater than the luminosity estimate from the radio continuum and maser data.

The same method provides us with the total mass of the stellar cluster, which ranges from 5 × 10² M_☉ to 2 × 10³ M_☉, depending on the statistical distribution of the stars in the cluster. This stellar mass can be divided by the total mass of the associated molecular clump, 3 × 10³ M_☉ (see Section 4.4), to calculate the star formation efficiency, 17%–67%. Albeit large, values like these are not implausible, as we are considering the part of the molecular cloud where star formation concentrates (Lada & Lada 2003; Bonnell et al. 2011).

In conclusion, we believe that the observational evidence collected so far for source HMC suggests that one is observing a young star of 20–30 M_☉, still in the main accretion phase and surrounded by a rich stellar cluster. Such a conclusion holds independently of the scenario adopted to explain the free–free emission (thermal jet or weak H ii region) and raises a question: why is only faint free–free emission detected from this object? A ∼3 × 10⁴ L_☉ star is expected to ionize the surrounding gas emitting a flux density of ∼2 mJy (in the optically thin limit), whereas toward source HMC less than 1 mJy is measured. A possible explanation is that the H ii region is at least partly quenched by accretion from the parental core, as depicted, e.g., in the model by Keto (2002, 2003). Although no evidence of accretion is seen in a typical infall tracer such as HCO⁺ (Furuya et al. 2008), this fact is inconclusive because the majority of the detected HCO⁺(1–0) might be part of the outflowing gas, thus hindering the detection of infall.

Finally, it is also worth mentioning that according to recent theoretical models (Hosokawa et al. 2010), large accretion rates are bound to bloat the protostellar radius and decrease the corresponding temperature and Lyman continuum flux. This could be an alternative explanation of the lack of a (strong) H ii region.

6. SUMMARY

We observed the HMSFR G28.87+0.07 with the VLBI phase-referencing technique in three powerful maser transitions: 22.2 GHz H₂O, 6.7 GHz CH₃OH, and 1.665 GHz OH. In addition, we also performed VLA observations of the radio continuum emission at 1.3 and 3.6 cm with both the A- and C-array, and Subaru/COMICS observations of the mid-infrared continuum emission. We also made use of data from the Hi-GAL/Herschel, ATLASGAL/APEX, as well as other surveys of the Galactic plane. From all these observations, we draw the following conclusions.

1.
The bipolar distribution of l.o.s. velocities and the general pattern of observed relative proper motions indicate that the water masers are tracing expansion with an average velocity of ∼24 km s⁻¹ from a point on the sky close to the maser center of motion. While at larger distance from the center of motion, maser velocities are higher and better collimated along the northeast–southwest direction, at closer distances water masers show smaller speeds and a larger scatter in the velocity orientations. We interpret these facts as the result of the interaction of a (proto)stellar jet with dense circumstellar material.
2.
Two continuum sources (HMC and A), separated by ∼31 along the north–south direction, are observed with the VLA at 3.6 and 1.3 cm. The HMC source is compact and spatially associated with the observed H₂O, CH₃OH, and OH masers. Whether the free–free emission arising from this source is due to an H ii region or a thermal ionized jet cannot be unambiguously established with the present data. However, we believe that the jet hypothesis is more consistent with the direction of the water maser motions, which outline expansion in a bipolar flow. Moreover, the orientation and momentum rate of the putative jet are compatible with those of the molecular outflow detected at arcsecond scales. We hence argue that the jet traced by the water masers is also powering the large-scale outflow.
3.
The combined information obtained from subarcsecond imaging with Subaru/COMICS at 24.5 μm and lower angular resolution images from the Hi-GAL/Herschel and ATLASGAL/APEX surveys, permit to estimate with good accuracy the luminosity (2 × 10⁵ L_☉) and gas mass (3 × 10³ M_☉) of the star-forming region of G28.87+0.07, and establish that ∼90% of the luminosity is coming from the radio source HMC. We conclude that this source must contain multiple stars, with the most massive being at least as luminous as 3 × 10⁴ L_☉. The lack of an associated H ii region indicates that the 20–30 M_☉ star must be still undergoing heavy accretion from the surrounding envelope.

We thank the anonymous referee for carefully reading the text and for constructive criticisms that improved the presentation and analysis of our results. This work was supported by the Chinese NSF through grants NSF 11133008, NSF 11073054, NSF 10733030, and NSF 10703010.

MASSIVE STAR FORMATION TOWARD G28.87+0.07 (IRAS 18411−0338) INVESTIGATED BY MEANS OF MASER KINEMATICS AND RADIO TO INFRARED CONTINUUM OBSERVATIONS

Article metrics

Permissions

Author e-mails

Author affiliations

Dates

ABSTRACT

1. INTRODUCTION

2. THE HMSFR G28.87+0.07