© ESO, 2014

1. Introduction

Sagittarius A West, the so-called mini-spiral, is a prominent feature of the interstellar medium (ISM) in the inner ~2 pc of the Milky Way (Ekers et al. 1983; Lo & Claussen 1983). It is made up of three streamers of ionized gas, termed the Northern, Western, and Eastern arms. Near- and mid-infrared observations of the mini-spiral show that the ionized gas is mixed with warm and hot dust that is in part arranged in complex structures like: extended filaments, bow-shocks and compact, barely resolved knots (Clénet et al. 2004b; Viehmann et al. 2006; Mužić et al. 2007). Studies of the radial velocities and proper motions of the gas in the mini-spiral have shown that their dynamics can generally be described by three bundles of quasi-Keplerian orbits in the combined potential of the central stellar cluster of the Milky Way plus a supermassive black hole (e.g., Vollmer & Duschl 2000; Paumard et al. 2004; Zhao et al. 2009).

It is generally thought that the gas in the mini-spiral is ionized by the radiation of more than one hundred massive, young stars that are located mainly within 0.5 pc of the supermassive black hole, Sagittarius A*, and that were formed in a star formation event about 2–6 Myr ago (Paumard et al. 2006; Bartko et al. 2010; Lu et al. 2013). Some of these young stars appear to orbit the black hole in a disk-like structure that rotates in a clockwise direction (referred to as clockwise system, CWS, in the following). This observation, along with the finding that the young stars are concentrated toward the black hole, gave rise to the currently favored hypothesis for their formation in situ in a formerly existing gas disk around Sagittarius A* (Levin & Beloborodov 2003; Paumard et al. 2006; Lu et al. 2009; Bartko et al. 2009). Considerable open questions in this picture remain, however, in particular, concerning the nature of the young stars that are not part of the CWS (see, e.g., discussion in Genzel et al. 2010): is the second postulated system of a counterclockwise streamer or disk (CCWS) real? Is the CWS strongly warped? Can the current kinematical configuration of the young, massive stars be reached when one assumes that all of them are on the CWS? It is therefore of great interest to obtain as much information as possible about the orbital parameters of all young, massive stars.

Very many bow-shocks show evidence of interactions between stellar winds from the young massive stars and the ISM in the mini-spiral (Clénet et al. 2004a; Tanner et al. 2005; Geballe et al. 2006; Mužić et al. 2007; Buchholz et al. 2011). Some of these bow shocks are among the brightest near- to mid-infrared sources in the central parsec of the Galactic center (GC) and have been known for decades (Rieke & Low 1973; Becklin & Neugebauer 1975) before their nature was finally revealed through radiative transfer modeling and, in particular, high angular resolution imaging (Tanner et al. 2002, 2005).

Tanner et al. (2005) analyzed high angular resolution near-infrared (NIR) observations of the sources IRS 1W, 5, 10W, and 21 and modeled their bow shocks. From the large measured stand-off distances of the bow shocks’ apices the authors concluded that the central sources must be Wolf-Rayet stars and thus were created in the above-mentioned star formation event. Because of the highly embedded nature of the sources, no spectroscopic confirmation of this classification could be obtained so far, with the exception of a tentative identification of IRS 1W as a Be star by Paumard et al. (2006). The lack of suitable spectra means that the line-of-sight velocities of the bow-shock stars are unknown and cannot be used to estimate their orbital parameters.

Nevertheless, bow-shock (BS) sources are of great interest because they allow us, in principle, to infer all six of their phase space parameters, and thus one more (the line-of-sight distance) than in the case of other young stars. Their 2D positions and proper motions can be inferred from imaging. Their line-of-sight (LOS) velocities and distances can be inferred from their interaction with the mini-spiral and from the resulting shape of the bow shocks in combination with models of the gas flow in the mini-spiral (e.g., Vollmer & Duschl 2000; Paumard et al. 2004; Zhao et al. 2009). In this way, Tanner et al. (2005) were able to show that the kinematic properties of IRS 1W and 10W, as determined from their measurements, are consistent with the hypothesis that these stars belong to the CWS. However, they did not have the necessary proper motion measurements to perform this kind of analysis on IRS 5. For IRS 21, Tanner et al. (2005) were unable to resolve its bow-shock and concluded that it was seen face-on (see also Tanner et al. 2002). Buchholz et al. (2011), however, were recently able to resolve the bow-shock structure of IRS 21 through deconvolution of high-quality adaptive optics data and could thus discard the face-on hypothesis.

Here we continue the investigation of the bow-shock sources IRS 1W, 5, 21, and 10W by Tanner et al. (2005). We use our own multiwavelength AO assisted imaging as well as AO-assisted sparse aperture masking (SAM; see, e.g., Tuthill et al. 2006a) to investigate the structure of the sources, which we then fit with fully 3D models. The results from the modeling are then combined with new, more accurate proper motion measurements and with models of the mini-spiral flow to investigate the kinematics of the bow-shock stars.

Section 2 describes the acquisition and reduction of the data. The calibration process and image reconstruction of our SAM observations is described in particular detail because we describe a novel way of self-calibrating SAM data in a dense stellar field. This procedure can be of general interest for SAM observations of targets in dense stellar fields. Section 3 describes the results of the proper motion measurements. Section 4 describes the mathematical model used in determining the shape and emission of the bow shocks and the corresponding fitting procedure. Section 5 describes the local kinematics of the Northern Arm at the position of the bow shocks and the procedure of determining the probability density function of the source orbital planes. Finally, in Sect. 6 we discuss our results and present our conclusions.

2. Observations and data reduction

2.1. AO imaging

We used AO-assisted imaging observations of the central parsec of the GC with VLT/NACO to measure the proper motions of the bow-shock sources and to infer their structure. Table 1 summarizes the observations. The bright supergiant IRS 7, located about 5.5′′ north of Sgr A*, was used to close the loop of the AO. IRS 7 is a variable source (Ott et al. 1999). The most recently published values for its flux are K_s = 7.69±0.06 in April 2006 and K_s = 6.96±0.04 in March 2009 (see Schödel et al. 2012). All imaging data were sky subtracted, flat fielded, corrected for bad pixels, and then combined via a shift-and-add algorithm into final mosaics, as described, for instance, in Schödel et al. (2009). The saturated cores of bright stars were repaired in all images from epochs 2002 to 2008. The images from 2009 to 2011 did not need to be repaired since saturation was not a problem because of the short exposure times used for these observations. The data from 2010 and 2011 were taken with NACO’s cube mode, which registers each individual exposure (see Girard et al. 2010). To improve the image quality of the latter data we applied a lucky imaging method and selected the 30% of the frames with the highest Strehl ratios.

To remove the seeing halos in the AO images, we extracted point spread functions (PSFs) from the H- and K_s-band images that were obtained with the S13 camera of NACO and deconvolved the images with the Lucy-Richardson algorithm (Lucy 1974). The deconvolved images were restored with a Gaussian beam of 2 pixels FWHM to create the final images to be used in our analysis.

2.2. Observations with sparse aperture masking

Our SAM observations of the GC central parsec were performed in the nights of 3/4 and 4/5 June 2010 with VLT/NACO. We used the L27 camera (0.027″/pixel) and chose two combinations of filters and mask: the L′ broad-band filter (λ_central = 3.80 μm, Δλ = 0.62 μm) was combined with theBB_9holes mask and the NB_3.74 narrow-band filter (λ_central = 3.74 μm, Δλ = 0.02 μm) with the 9_holes mask. The latter combination reduces bandwidth smearing, while the first combination leads to higher sensitivity. The observations were recorded using pupil tracking.

The wide field of view (FOV; 27′′ × 27′′) of the L27 camera allowed us to observe all targets simultaneously. All SAM data were pre-processed according to standard imaging data reduction, that is, sky subtraction, flat-fielding, and bad-pixel correction.

Depending on the observational data sets, the sky was obtained by calculating either (i) the median of separate sky observations obtained on a dark cloud about 400′′ north and 700′′ west of Sagittarius A*, or (ii) via the median of the appropriately clipped values at each pixel from the dithered images on-target when no dedicated sky observations close in time were available. We note that sky subtraction did not result in a completely flat background as expected. Instead, patterns remained in the images, probably because of the rapid and variable atmospheric conditions. However, these patterns varied on scales of ~2′′, that is, on scales larger than the size of the SAM PSFs. Their impact on the calibration of the SAM observations was therefore considered not significant. A reduced individual image frame (L′, BB_9holes) is shown in Fig. 1.

The so-called cube mode was used to save each single-detector integration time (DIT) frame. The advantage of this mode is that it allows selecting the best frames. In our case, frame selection was necessary because atmospheric conditions (and the AO correction) were affected by highly variable seeing and the passage of clouds. About 80% of the frames were de-selected (using a simple flux criterion on the star IRS 7). The final set of data consists of a group of 13 cubes composed of 20−30 exposures with 5 s and 20 s DIT each for L′ and NB_3.74 filters, respectively. Finally, the selected frames for each of the individual targets were stored into separate cubes of 128 × 128 pixels.

The standard observing technique for SAM requires closely interwoven observations between calibrator and target. The calibrator must be a nearby point-like source (within ~1°–2°). We obtained observations of the calibrator sources HD 161718 and HD 162907. A dither pattern was used in which the source was located at the center of the four quadrants of the detector so that the sky could be extracted by a simple median superposition of the dithered target observations.

2.3. Extraction and calibration of the interferometric observables from the SAM data

For the subsequent processing of the SAM data we used an IDL pipeline developed at Sydney University (priv. comm., M. Ireland). This code constructs a matched filter for the given instrumental setup (pixel scale, observing wavelength, filter, aperture mask) and uses it to measure the squared visibilities, V², for all interferometric baselines as well as the bispectra, which are the ensemble of the visibility products for all non-redundant baseline triangles. The arguments of these bispectra are the closure phases (CPs). The CPs provide information of the point-symmetry of the source and are free of telescope-dependent phase errors induced by atmosphere and telescope vibrations. Finally, to eliminate the baseline-dependent errors, the calibrated V² were obtained by the ratio of the raw V² of the target and calibrator, while the calibrated CPs were obtained by subtracting the calibrator CPs from the CPs of the target.

We note that the u − v coverage (u − v refers to the coordinates in the interferometric plane, i.e., the coverage of the Fourier transform of the image) of the targets in our observation was not limited by the properties of the interferometric masks used, but was increased through Earth rotation synthesis by combining observations of the targets at different parallactic angles.

To accurately calibrate the SAM observations, it is generally recommended to frequently switch between the target and the calibrator while leaving all AO settings unchanged. But this was not possible for our observations because of the highly variable conditions (seeing, clouds) that made changes to the AO setting necessary, delayed calibrator observations, and led to unstable correction. We thus lacked sufficient, adequate calibrator observations.

One of the possible solutions to overcome this problem was using the IRS 7 star as calibrator because it is contained in all our SAM frames. This source was already tested and used by Tuthill et al. (2006b) as a calibrator for SAM observations. However, they performed their observations at shorter wavelengths than we did. IRS 7 is known to be surrounded by extended emission at MIR wavelengths (Viehmann et al. 2006), which may already show some contribution at the long NIR wavelengths used here. In addition, the images of IRS 7 may have been affected by local systematic variations of the background and do not properly reflect the variance of PSF and background for the entire FOV of our observations.

To overcome these problems we developed a technique that we called the synthetic calibrator. This method takes advantage of the many bright sources in a crowded FOV as is the case for the GC. The procedure is very similar to the iterative extraction of an empirical PSF from a crowded field that is implemented in PSF fitting algorithms such as StarFinder (Diolaiti et al. 2000). We used StarFinder to (i) select the brightest stars in the SAM images (excluding known, obviously extended sources such as IRS 1W) to construct a preliminary PSF by median superposition; (ii) clean the images from secondary sources with the help of the preliminary PSF and subtract large-scale variations of the background emission; and (iii) obtain a final PSF from the median superposition of the cleaned images of the reference stars. The final PSF, which corresponds to the interferogram of a point source, was normalized to a total flux of one. We used the sources IRS 7, IRS 3, IRS 29, and IRS 6E to construct the synthetic calibrator.

This technique provides a reliable way to minimize the contribution of all the potential systematic errors related to each one of the individual sources. Even more, the synthetic calibrator is robust, under certain constraints, to the use of slightly extended sources as inputs (like IRS 3). These restrictions require that the in-going sources have no strong closure phase information (i.e. the source is not asymmetric) and that their emission is dominated by the point-like source. We can observe, at first glance, in Fig. 1 that the four sources used to build the synthetic calibrator fulfill the requested constraints and have a significantly crispier PSF than our science targets.

We used the few valid calibrator observations to verify this technique by comparing the calibrated V² and CP values that were obtained with the synthetic calibrator and with a real calibrator. No systematic differences were found. We also compared the resulting calibrated V² and CP obtained with the synthetic calibrator and with IRS 7 as calibrator. Again, the values agreed within their uncertainties. This also confirms that IRS 7 is unresolved with the chosen SAM setup. Figure 2 shows a comparison of the calibrated V² and CPs using the synthetic calibrator technique and IRS 7 as calibrator for one of our sources. We decided to use the synthetic calibrator technique to calibrate the V² and CP of all our targets. Figure 2 illustrates the information that can be obtained from our SAM measurements: IRS 1W is completely resolved-out for baselines longer than 4 m, which indicates that the source angular size is extended and larger than 120 milliarcsec; the closure phases differ significantly from zero, which excludes a point-symmetric structure for the source.

We note that the synthetic calibrator technique can only be used when several sufficiently bright sources are present within the FOV and when the PSF is not expected to vary significantly over the FOV. It may be ideal for observations of stars in dense clusters, where it can both improve the calibration of SAM data and boost the efficiency of SAM observations because it does not require any change of the instrument setup and avoids overheads caused by calibrator observations.

2.4. SAM image reconstruction

Image reconstruction from the interferometric data was performed with the BSMEM package (see, e.g., Buscher 2013; Lawson et al. 2004). This code uses a maximum-entropy algorithm. Table 2 lists the number of iterations and the reduced χ² values at the end of the deconvolution process for all sources. The size of the beam used in the reconstruction procedure is θ = 60 mas, which corresponds to the nominal angular resolution reached by this technique for the filters we used. The pixel scale used for the reconstruction process was 13 mas/px.

3. Proper motions

To measure the proper motions of the stellar sources at the center of the bow shocks we used VLT/NACO K_s-band AO imaging data obtained between 2002 and 2011 (see Table 1). Schödel et al. (2009) obtained proper motion measurements of our sources using PSF fitting to measure their stellar positions. This method is not suited for extended sources like the bow shocks that we analyze here, however. Instead, we measured the stellar positions by centroiding the emission maxima of the central stars in circular apertures of with a radius of 4, 5, and 6 pixels. The positions and their uncertainties were taken as the mean of the three measurements and its corresponding uncertainty. Note that measuring stellar positions by aperture centroiding is significantly less susceptible to systematic errors caused by (moderate) saturation or nonlinearity than is PSF fitting. Moreover, since the bow shocks are extended targets, PSF fitting to their central stars is very difficult and prone to systematic errors.

Astrometric alignment was achieved with the help of the known positions and proper motions of a sample of about 200 stars contained in all images. The astrometric reference stars were selected from the list published by Schödel et al. (2009). The stars had to be brighter than K_s = 12, with no companion within 0.3′′ and the uncertainty of each component of the proper motions of the reference stars had to be smaller than <30 km s^-1. These criteria ensure accurate aperture measurements of their positions.

We transformed the stellar positions into the astrometric reference frame via a second-order polynomial fit (see Appendix A in Schödel et al. 2009). A Monte Carlo approach was used to estimate the uncertainty of this alignment procedure, with 1000 trials. For each trial the positions of the reference stars were varied randomly, with the offsets drawn from Gaussian distributions with standard deviations corresponding to the positional uncertainty of each star at each epoch.

After obtaining the stellar positions for each epoch, we fitted the data points by linear motions. The reduced χ² was generally significantly larger than 1. There were some sources of uncertainty that were not specifically taken into account, such as errors introduced by camera distortion into the mosaics, the shape of the bow shocks, pixel saturation, or faint secondary stars within the aperture laid over the bow shocks. Hence, the uncertainties of the stellar positions may have been under-estimated. This is no serious problem, however, because we can be certain that linear motions are accurate models to describe the proper motions. Therefore, the uncertainties of the linear fits were all rescaled to a reduced χ² = 1. This procedure avoids underestimating the uncertainties (see Schödel et al. 2009, for a more elaborate discussion). The measured angular velocities were converted into km s^-1 by assuming a Galactic center distance of 8 kpc (e.g., Schödel et al. 2010). Plots of position vs. time along with the linear fits for the proper motions are shown in Appendix A.

Table 3 lists the proper motions of the bow-shock sources and of the stars IRS 16NW and IRS 33N as determined here and by Schödel et al. (2009). IRS 16NW and IRS 33N are two point-like sources that we used as a reference point to compare our method and Schödel et al. (2009). The measurements of IRS 16NW and IRS 33N and of the bow-shock sources all agree within their 1σ uncertainties, with the exception of IRS 1W. We believe that the discrepancy on the latter source is probably due to the PSF fitting technique used by Schödel et al. (2009), which is pooröly suited to strongly resolved/extended sources like IRS 1W. Moreover, IRS 1W is one of the brightest sources in the field and is affected by saturation in many images. Aperture fitting is suited to both extended and (moderately) saturated sources. We therefore believe that the proper motions for the bow-shocks’ stellar sources that we have derived here may be more accurate than those published by Schödel et al. (2009). Tanner et al. (2005) measured the proper motions of IRS 21, IRS 1W, and IRS 10W, which are also listed in Table 3. Their measurements generally agree with ours within the quoted 1σ uncertainties, with the exception of the velocity along the right ascension of IRS 21. Note that the proper motions of Tanner et al. (2005) suffer from very large uncertainties, possibly because of the small FOV that they used for determining the proper motions and because they used images reconstructed from speckle observations that have a significantly lower signal-to-noise ratio than we used (AO imaging).

4. Modeling the bow shocks

4.1. Bow-shock images

In Fig. 3 we show a compilation of our images of the bow-shock sources IRS 1W, IRS 5, IRS 10W, and IRS 21. The bow shock is clearly resolved in almost all sources in at least one wavelength. IRS 21 and IRS 10W are more clearly resolved than in the work of Tanner et al. (2005), which can be attributed to technological progress. While Tanner et al. (2005) used low Strehl shift-and-add reconstructed speckle data, here we use high Strehl and high signal-to-noise AO images and images reconstructed from AO assisted sparse aperture masking, a technique that provides a precise calibration of the PSFs. Note that Buchholz et al. (2011) and Buchholz et al. (2013) presented a resolved image of the bow shock around IRS 21. Hence, this source is not seen face-on as hypothesized by Tanner et al. (2005; see also Tanner et al. 2002).

The extended emission becomes stronger toward longer wavelengths, a consequence of the temperature of the dust in the bow shocks. Tanner et al. (2002) derived a dust temperature of ~1000 K for IRS 21. Viehmann et al. (2006) also showed that the temperatures of the bow-shock sources peak at around 4–10 μm, consistent with a dust temperature of a few hundred degrees. On the other hand, the central star is only detected at the shorter wavelengths (see Fig. 3).

4.2. Bow-shock shape

The shape of the bow shocks results from the interaction between the stellar wind and the surrounding ISM, taking into account the relative motion of the star with respect to the ISM. Assuming a steady-state model, the shape of the bow shock results directly from the stand-off distance, R₀, where the ram pressure from the relative motion through the ISM is balanced by the pressure from the stellar wind, ${\mathit{R}}_{\mathrm{0}}\mathrm{=}\sqrt{\frac{\mathit{ṁ}\mathrm{\ast }{\mathit{\nu }}_{\mathrm{w}}}{\mathrm{4}\mathit{\pi }{\mathit{\rho }}_{\mathrm{ISM}}{\mathit{\nu }}_{\mathrm{\ast }}^{\mathrm{2}}}}\mathit{,}$(1)where ṁ_∗ is the mass-loss rate of the star, ν_w is the velocity of the stellar wind, ρ_ISM is the density of the surrounding ISM, and ν_∗ is the relative velocity between the star and the ISM.

To obtain a steady-state solution for the bow-shock shell, Wilkin (1996) assumed (i) that the bow-shock shape is determined by the momentum of the shocking winds; (ii) that the cooling process is instantaneous; (iii) that the shell is a thin layer of negligible width; and (iv) that the interacting winds are homogeneous and do not change with time. The 2D polar canonical form of a bow shock, under the previous conditions, is then defined by $\mathit{R}\mathrm{\left(}\mathit{\theta }\mathrm{\right)}\mathrm{=}{\mathit{R}}_{\mathrm{0}}\csc \mathit{\theta }\sqrt{\mathrm{3}\mathrm{(}\mathrm{1}\mathrm{-}\mathit{\theta }\cot \mathit{\theta }\mathrm{)}}\mathit{,}$(2)where θ is defined as the angle from the axis of symmetry to any point on the 2D-shell’s surface. To extend the previous relation from 2D to 3D, we let the locus of R(θ) rotate by angle φ as shown in Fig. 4, hence ${\mathit{R}}_{\mathrm{3}\mathrm{D}}\mathrm{\left(}\mathit{\theta ,\phi }\mathrm{\right)}\mathrm{=}\mathit{rot}\mathrm{\left\{}\mathit{R}\mathrm{\right(}\mathit{\theta }\mathrm{\left)}\mathit{,\phi }\mathrm{\right\}}\mathit{,}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}{\mathrm{0}}^{\mathrm{◦}}\mathrm{\le }\mathit{\phi }\mathrm{\le }{\mathrm{180}}^{\mathrm{◦}}\mathit{.}$(3)The cells of the bow shock are defined by sampling R_3D(θ,φ) with increments in θ and φ of 0.5° (see Fig. 4 ).

4.3. Bow-shock emission

4.3.1. Polarization and emission mechanisms

Buchholz et al. (2013) described that the bow-shock sources, particularly IRS 21 and IRS 1W, exhibit strong polarization levels at NIR wavelengths. The polarization increases from shorter to longer wavelengths, highlighting the importance of the extended dust component. Polarization of the extended emission and of the embedded bow-shock sources in the Northern Arm (NA) reaches its maximum in the MIR (see Glasse et al. 2003; Buchholz et al. 2011, 2013). Buchholz et al. (2011) found that the polarization of the bow-shock sources IRS21 and IRS 1W agrees with the orientation of the magnetic field in the NA. These authors also suggested that the ordered magnetic field is responsible for the dust grain alignment in the bow-shock shells. The polarization measurements from NIR to MIR thus support the assumption that the emission in the bow shocks is produced by the dust IR thermal emission, which is reprocessed UV light from the central star.

4.3.2. Bow-shock emission model

Although the previously mentioned findings suggest that the bow-shock emission is mainly produced by thermal emission, we also explored the role of scattered light in the bow-shock emission, particularly at shorter wavelengths. The procedure used to simulate the emission maps was similar to the one developed by Mužić et al. (2010). The dust thermal emission, L_th, was assumed to have a blackbody source function attenuated by the corresponding optical depth, that is, L_th ∝ B(T_d)(1−e^{− τ_abs}). The dust temperature distribution in the bow shock (T_d) is assumed to be created by heating the dust by the stellar photospheric UV photons from the central source. Therefore, under this condition, T_d is computed via the van Buren & McCray (1988) solution, ${\mathit{T}}_{\mathrm{d}}\mathrm{=}\mathrm{27}{\mathit{a}}_{\mathit{\mu }\mathrm{m}}^{\mathrm{-}\mathrm{1}\mathit{/}\mathrm{6}}{\mathit{L}}_{\mathrm{\ast }\mathit{,}\mathrm{38}}^{\mathrm{1}\mathit{/}\mathrm{6}}{\mathit{r}}_{\mathrm{pc}}^{\mathrm{-}\mathrm{1}\mathit{/}\mathrm{3}}\mathrm{K}\mathit{,}$(4)where a_μm is the grain size in microns, L_{∗ ,38} is the star UV luminosity in units of 10³⁸ erg/s and r_pc is the distance to the dust from the star in parsecs. After calculating the contribution from each cell, the total L_th was normalized to one.

The scattering emission was computed using the normalized Henyey & Greenstein (1941) source function modified by the optical depth, ${\mathit{L}}_{\mathrm{sca}}\mathrm{\propto }{\mathit{r}}_{\mathrm{cm}}^{-2}\mathit{P}\mathrm{\left(}{\mathit{\theta }}_{\mathrm{sca}}\mathrm{\right)}{\mathrm{e}}^{\mathrm{-}{\mathit{\tau }}_{\mathrm{sca}}}$, where r_cm is the distance to the dust from the star in centimeters. The scattered light is described by $\mathit{P}\mathrm{\left(}{\mathit{\theta }}_{\mathrm{sca}}\mathrm{\right)}\mathrm{=}\frac{\mathrm{(}\mathrm{1}\mathrm{-}{\mathit{g}}^{\mathrm{2}}\mathrm{)}}{\mathrm{(}\mathrm{1}\mathrm{+}{\mathit{g}}^{\mathrm{2}}\mathrm{-}\mathrm{2}\mathit{g}\cos \mathrm{(}{\mathit{\theta }}_{\mathrm{sca}}\mathrm{)}{\mathrm{)}}^{\mathrm{3}\mathit{/}\mathrm{2}}}\mathit{,}$(5)where g is a dimensionless term that depends on the scattering angle θ_sca, with g = 0 when the scattering is fully isotropic and g < 1 when the scattered light preserves the angle of incidence. After calculating the contribution from each cell, the total L_sca was normalized to one. The previous treatment of the scattering emission does not take into account the polarization of light, only the deflection angle of the incident radiation in the dust. A more detailed model of the scattering emission is beyond the scope of this work.

The total bow-shock emission, L_tot, at every cell in the model was computed adding the previously normalized thermal and scattering emissions, L_tot = ϵ_scaL_sca + ϵ_thL_th. The coefficients ϵ_th and ϵ_sca are the relative efficiencies (ϵ_th + ϵ_sca = 1) of the thermal and scattering contribution. The optical depths, τ_abs(λ) and τ_sca(λ), at every cell of the 3D bow-shock shell are calculated as $\begin{array}{ccc}{\mathit{\tau }}_{\mathrm{abs}}\mathrm{\left(}\mathit{\lambda }\mathrm{\right)}& \mathrm{=}& \\ {\mathit{\tau }}_{\mathrm{sca}}\mathrm{\left(}\mathit{\lambda }\mathrm{\right)}& \mathrm{=}& \end{array}$where C_abs(a,λ) and C_sca(a,λ) are the dust thermal and scattering extinction coefficients. They are calculated as C_abs(a,λ) = πa²(Q_abs) and C_sca(a,λ) = πa²(Q_sca). The values of the thermal and scattering efficiencies, Q_abs and Q_sca, were obtained form the Draine database on the optical properties of interstellar dust grains¹ for a standard mixture of 50% graphite and 50% silicate (see Mathis 1990; Draine & Lee 1984). The grain size distribution used, n_d(a), follows an MRN (Mathis et al. 1977) power-law with the grain size a ranging from a − = 0.005 μm to a + = 0.25 μm.

4.4. Bow-shock fitting

We used the 3D steady-state model and the dust emission processes described in the two preceding sections to create models of the bow-shock shells and fit them to the images. The best-fitting models were found by a χ²-minimization as implemented in the IDL MPFIT package by C. Markwardt (Markwardt 2009). Masks were applied to the bow-shock images to suppress features that are not related to the actual bow shocks. Each masked image was normalized to unity inside the applied mask. The same normalization was applied to the projected and masked models. A central point-source was included in the model (with free parameters for its position and flux). This point-like source concentrates most of the NIR-flux and dominates the model fitting at the shortest wavelengths. We note that the pixel uncertainties are hard (or impossible) to determine in Lucy-Richardson deconvolved images and in the BSMEM reconstructed images. In addition, there are several sources of systematic uncertainties, such as faint point sources or ISM features that may be confused with the bow shocks. Therefore, both uniform and flux-weighting schemes were tested to determine the model with the lowest uncertainties that is least affected by the systematics previously discussed. Our tests show that a flux-weighting scheme produced better fits in H and K_s while a uniform weighting works better in L′. The final uncertainties of the best-fit parameters are obtained from their standard deviation of the fits to the statistically independent images from the three different filters. We also expect this procedure to reliably take care of possible contamination from stars and extended emission because the former are expected to be faint at L′ and the latter at H.

The stand-off distance R₀ (reported in Table  4), which is strongly degenerated with the inclination of the bow shock, was linked in our fitting code to the projected stand-off distance in the plane of the sky, ${\mathit{R}}_{\mathrm{0}}^{{}^{\mathrm{\prime }}}\mathrm{=}{\mathit{R}}_{\mathrm{0}}\mathit{/}\cos \mathrm{(}\mathit{i}\mathrm{-}{\mathrm{90}}^{\mathrm{◦}}\mathrm{)}$, where i is the bow-shock inclination angle between the LOS and the bow-shock apex. When the bow-shock apex is pointing to the observer i = 0° and when the bow-shock shell is projected in the plane of the sky i = 90° (see inset in Fig. 4). ${\mathit{R}}_{\mathrm{0}}^{\mathrm{\prime }}$ was obtained from the mean of a Gaussian fit to the emission obtained by three radial cuts along the bow-shock structure. The central cut was made along the symmetry axis of the bow shock, the two additional cuts were at − 10° and + 10° from the position of the central cut. Figure 5 shows as an example the radial cuts, the fitted Gaussian, and the ${\mathit{R}}_{\mathrm{0}}^{\mathrm{\prime }}$ mean value obtained for IRS5 in the H-filter. A similar procedure was performed for the other sources. The rest of the parameters described in Table 4 were fitted directly to the image described at the beginning of this section.

4.5. Contribution from dust scattering to the bow-shock emission

To quantitatively approximate their relative importance, we investigated the effect of different relative thermal and scattering efficiencies on our models. Since IRS 5 appears to resemble the canonical steady-state bow-shock model most, we used the images of this source in H and K_s to perform this test. We fit different sets of models with fixed efficiencies for thermal and scattering emission at different wavelengths and considered cases of isotropic and anisotropic scattering. For each wavelength, model fitting was performed with increments of 10% in the scattering efficiency from 0 to 70% of the total emission observed. We tested the cases for g = 0 in the scattering source function (isotropic scattering) and for g = 0.5 (partially anisotropic scattering).

The best-fit bow-shock model was found in each trial with the procedure described in Sect. 4.4. Residual images were obtained by subtracting the models from the original image. Figure 6 displays the sum of the squared pixel values of the residual image for each trial in H and K_s for the isotropic case. In both bands, we found that the residuals increase with the scattering contribution, the model with 0% scattering had the lowest residuals. Similar results were found for the anisotropic scattering contribution (g = 0.5). As indicated in Sect. 4.1, we expect the main contribution to the bow-shock emission to be thermal emission from the dust grains. Our tests here are consistent with this picture. Therefore, all the models obtained for this work were fit only for the thermal emission of the dust.

4.6. Properties of the individual sources

Here, we discuss the individual bow-shock sources in detail. Their morphologies can be examined in Fig. 3, while Fig. 10 shows the best-fit bow-shock models overplotted to the images.

IRS 5: this bright bow-shock displays similar morphologies in the H, K_s, and L′ filters. The central star is clearly identified in H and K_s. The projected bow-shock shape is axially symmetrical in the plane of the sky with the vector from the central star to the stagnation point. IRS 5 is the source in our sample that can best be described with the isotropic wind steady-state model. Table 4 describes the best-fit parameters for the different wavelengths. Here, and for the other bow-shock sources, we conservatively estimate the uncertainties of the best-fit parameters from their standard deviation that results from the fits to the images at the different wavelengths (i.e., we did not use the uncertainty of the mean). The small uncertainties support the model assumption of a bow shock in the steady-state phase. Our results suggest that the bow shock of IRS 5 is inclined at ~130°.

From the interaction of IRS 5 with the NA and using the NA model of Paumard et al. (2004), we estimate a position of z = 3.81′′ ± 2.66′′ along the line-of-sight (LOS), where z = 0 corresponds to the plane of the sky. Assuming that IRS 5 is gravitationally bound to Sgr A*, we calculate a 3σ upper limit of ±336 km s^-1 for its LOS velocity.

IRS 1W: as Fig. 3 shows, the target is resolved at the three wavelengths. It has a clear bow-shock structure in K_s and L′, but also deviates from the analytical bow-shock model that we used here: first, a rounded shell appears around the central star that is open on its southern edge and is clearly visible at H and K_s. It is not observed in L′, however. Instead, the tails of the bow shock are better distinguished. The origin of this more rounded shell in the H and K_s images is unclear. Second, the shape near the apex is slightly asymmetrical in L′, with the maximum of the emission shifted to the west of the apex. Possible causes are confusion with a star or some feature in the mini-spiral, or an asymmetrical outflow from IRS 1W. To investigate this point further, the left panel in Fig. 7 shows a close-up of the surroundings of IRS 1W from a NACO L′-band image. The right panel shows a highpass-filtered version of the same image. IRS 1W does not appear to be contaminated by a typical mini-spiral feature such as a thin filament, nor by any bright star (see also Fig. 3). We therefore do not find any obvious or easily testable explanation for the observed asymmetry, but note here that the mini-spiral is generally a complex region with many knots and features in the ISM that may be confused in projection and/or interact with the bow shock. Third, IRS 1W displays a long, tightly closed tail. The tail of IRS 1W is clearly visible in the highpass-filtered image and has a length of approximately 0.8″, or 6400 AU.

IRS 1W shows a remarkable interaction with the ISM of the NA. The proper motion of IRS 1W is 350 ± 25 km s^-1 at an angle 332 ± 5° E of N. Assuming a mass of 4.2 ± 0.3 × 10⁶M_⊙ for Sgr A* (Gillessen et al. 2009; Yelda et al. 2011), the escape velocity at the projected distance, R = 5.30 ± 0.01′′, of IRS 1W is 419 ± 15 km s^-1. The LOS distance of IRS 1W can be estimated from its interaction with the Northern Arm of the mini-spiral. Using the model of Paumard et al. (2004), we estimate a LOS distance of z = 7.37 ± 1.5′′ for IRS 1W and obtain an escape velocity of 320 ± 29 km s^-1 at its 3D position. Although this is lower than its proper motion, both values overlap well within their 1σ uncertainties. Using the 3σ upper limit in the escape velocity and the 3σ lower limit on the proper motion, we can set a firm upper limit of V_z = 300 km s^-1 on the LOS velocity of IRS 1W (assuming that it is bound to Sgr A*, of course). This corresponds to a maximum inclination of its velocity vector of 137°, corresponding to 47° out of the plane of the sky, moving away from the observer. Here, the ambiguity between motion toward or away from the observer is broken by our best-fit bow-shock model. From the model of Paumard et al. (2004) we can calculate the properties of the gas flow. At the position of IRS 1W, the gas in the NA has a proper motion of 226 ± 70 km s^-1 at an angle 207 ± 13° E of N. The orientation of the tail-like feature agrees well with the relative proper motion between star and gas. The angle of the gas flow out of the plane of the sky is a mere 7°, with the gas moving toward the observer, which means that the NA flow is practically within the plane of the sky at the position of IRS 1W. We note that with a relative proper motion between star and gas on the order of 550 km s^-1 the channel-like feature that is marked by the tail of IRS 1W may trace the motion of IRS 1W over the past ~50 years. With such short timescales at work, IRS 1W may be an ideal probe for the conditions of the ISM in the NA, for instance, for the homogeneity of the gas density.

Table 5 shows our best-fit bow-shock model for IRS  1W using our steady-state model. The second row of images in Fig. 10 corresponds to the steady-state isotropic wind-fitted models. As can be clearly seen, our model cannot fit the tail of IRS 1W well (specially in L′) because it is much narrower than what can be produced with this model. Possible explanations of this problem include (i) changes in the density of the ISM at the shocking points; (ii) variations in the mass-loss rate and wind velocity; (iii) anisotropic winds; or (iv) that because of the short interaction times the bow shock in IRS 1W has not yet reached the steady-state bow-shock solution (see Mohamed et al. 2012). The large scatter in the best-fit position angles in Table 5 also suggests that the isotropic wind model is not suitable to reproduce the observed morphology of IRS 1W.

To examine whether there are other models that provide a better fit to the IRS 1W morphology, we decided to perform a model fitting using the narrow solution of a steady-state bow-shock shell described in Zhang & Zheng (1997). We assumed that the central source has a wind bound inside a cone with an opening angle of 40° along the polar axis. Table 6 displays the fitted parameters and the third row of images in Fig. 10 represents the best-fits for the narrow solution. A comparison between the isotropic wind model and the narrow solution (see Fig. 10) demonstrates that the latter fits the observed emission of IRS 1W better, with the main difference being in a smaller R₀ and different inclination. The narrow solution models have considerably lower residuals and result in best-fit parameters with very low scatter between the different observing wavelengths.

Another plausible explanation for the shape of IRS 1W can be changes in the density of the ISM. This condition is not rare in complex environments like the NA where filaments and clumps of gas/dust are observed. For example, Rauch et al. (2013) determined the geometry and emission of the asymmetric IRS 8 bow shock using a shell solution for a quadratic density gradient in the NA ISM. Figure 8 displays our best-fit solution to IRS 1W in L′ for an ISM density gradient. The density gradient function has the form ρ = ρ₀(1 + a₁x + a₂x²), where ρ₀ is the mean density of the ISM and x the direction of the ρ₀ increment. In our model, the density increased from east to west and the values of the coefficients a₁ and a₂ were fixed to 2.0 and 5.0, respectively.

This bow-shock gradient solution produced residuals on the order of those of the narrow solution. Similar values were obtained between the narrow and the gradient solution for the fitted parameters. The best-fit gradient solution produced a R₀ = 533 AU, and inclination of i = 103° and a PA = 9°. Nevertheless, to perform a complete least-squares fitting to bow-shock solutions with free gradient coefficients and/or different gradient functions creates strong degeneracies among the fitted parameters. Hence, we did not investigate other gradient solutions and decided to use the narrow solution best-fit model for the following analysis. Although we cannot fully distinguish between the physical mechanisms that produce the observed shape of IRS 1W, our narrow and gradient models suggest that inhomogeneities in the stellar wind and/or ISM are important in shaping the structure of this particular bow shock.

IRS 10W: the extended emission around this source presents a particularly different morphology in L′, more so than in H and K_s. The central star is detected in the H and K_s filters. Extended emission is observed at the position of the central source in L′. It is unclear whether this emission is related to the central star or not. The projected emission in H and K_s displays a symmetric structure, while in L′ the southward pointing tail is more extended than the northward pointing one. This feature may be related to inhomogeneities in the NA ISM at the location of IRS 10W (see Fig. 9). The radio image mapped at 1.3 cm by Zhao et al. (2009, see their Fig. 10) suggests that IRS 10W is part of the NA and is surrounded by ionized gas in a shell-like structure with a completely different morphology than the one observed at NIR wavelengths, showing no clear evidence of a bow-shock.

Because of the asymmetric shape of IRS 10W in L′, we decided to apply the model fitting only to the H and K_s images. Table 7 shows the best-fit parameters. The fourth row of images in Fig. 10 displays the best-fit models in contours. Although the extended emission around IRS 10W can be interpreted as a bow shock, a more detailed analysis of the IRS 10W proper motion and of the NA motion (at the position of IRS 10W) does not support the idea of a bow shock produced by the relative motion of the star throughout the ISM.

IRS 10W has a proper motion of 217 ± 14 km s^-1 with a PA of 355° E of N (see Fig 10 and Table 3). Paumard et al. (2004) reported a gas flow of 189 ± 84 kms^-1 with a PA of 188° E of N at the position of IRS 10W, while the values according to Zhao et al. (2009) are 229 ± 49 km s^-1 with a PA of 191°. Therefore, the relative proper motion of IRS 10W with respect to the NA is 403 ± 86 km s^-1 at a PA of 1° E of N, or 441 ± 37 km s^-1 at a PA of 3° E of N, that is, at a PA perpendicular to the direction of the possible bow-shock symmetry axis (around ~90°). This suggests that the extended emission around IRS 10W is not produced by the bow-shock model described in Sect. 4.2, which implies that the relative motion and the bow-shock axis must be parallel. Alternative explanations for the observed extended emission are (i) a shock between the stellar wind with a dust knot or filament of the NA; (ii) a shock between the stellar winds of IRS 10W and a secondary source; or (iii) a chance superposition of a some NA feature with IRS 10W. These hypotheses are plausible since IRS 10W is located in a very mixed environment of the extended emission in the NA, some filaments can be observed superimposed along the source position and, additionally, some close-by stars have extended emission (see Fig. 9).

IRS 21: the bow-shock of this source can be clearly seen at all wavelengths and is similar to the deconvolved K_s-band image presented by Buchholz et al. (2011). IRS 21 is best resolved in H band, where the central star is also unambiguously detected. We can discard the hypothesis, made by Tanner et al. (2005), that IRS 21 is seen face-on, since these authors could not resolve its structure in their lower Strehl no-AO observations. In our K_s and L′ images the central source is not detected and the bow-shock morphology is less clearly observed because of the lack of angular resolution. We therefore decided to apply the model fitting only to the H-band image. The steady-state model fits the morphology observed in the H-band very well. The proper motions of the star and of the NA environment are consistent with the relative motion of the bow shock obtained from our model. Table 8 lists the best-fit parameters.

Since we only performed the model fitting to the IRS 21 H-band image, the uncertainties in the fitted parameters were estimated in a different way than for the other bow shocks, for which we have at least two models in two different filters. For IRS 21, the uncertainties reported in Table 8 were estimated by testing three different starting values for the bow-shock inclination angle. We chose this approach because the inclination angle strongly affects the determination of the bow-shock stand-off distance. In fact, it is degenerated for low- and high-inclinations and therefore, our code works best at intermediate values of inclination angles (40° ≤ i ≤ 160°). Very low- or high-inclination angles result in bow-shock shapes that resemble spherical shells. Thereupon, our code cannot compute the stand-off distance since it strongly depends on the inclination angle (see Sect. 4.4). For IRS 21 we tested the code with initial conditions of i = 50°, 90°, 145°. The fitted parameters are reported in Table 8 and their errors are the standard deviation of the three trials. Although the model did not always converge to the same solution, which can be largely explained by degeneracies between the model parameters, the best-fit parameters, resulting from the runs with significantly different starting values for the inclination, have overall similar best-fit parameters with a low scatter between them.

5. ISM velocity field and bow-shock orbital parameters

5.1. Northern Arm dust motion

Our best-fit bow-shock models provide estimates of the local motion of the ISM in the NA. Assuming that our objects are Wolf-Rayet stars (Tanner et al. 2005) with a mass-loss rate on the order of $\dot{{\mathit{m}}_{\mathrm{\ast }}}\mathrm{~}\mathrm{1}$–2 × 10^-5M_⊙/yr and stellar winds with velocities around ν_w ~ 1000 km s^-1, and that the density in the NA has a mean particle number n_H ~ 3 × 10⁴ cm^-3 (see Tanner et al. 2002; Rauch et al. 2013), we first calculated the relative velocity between star and ISM with Eq. (1). Subsequently, we d the measured proper motion of the central star to calculate the proper motion of the ISM at the position of the bow-shocks. Figure 10 shows the projected motion of the ISM, V_ISM(red vector), the relative motion of the bow shock, V_{s − ISM}, according to our best-fit model (blue vector), and the proper motion, V_s, of the central source (green vector) for all our targets. All our ISM motion measurements (magnitude and position angle) agree to within 1σ with the models of Paumard et al. (2004) and/or Zhao et al. (2009) except for IRS 10W, which supports our hypothesis that IRS 10W does not represent a standard bow shock (see Sect. 4.5). Table 9 summarizes the NA projected motions obtained in this work and the previously published models of gas and dust. These results corroborate the hypothesis of Wolf-Rayet stars as the central sources of the bow shocks. For comparison, if we compute the projected NA motions using O-stars, ($\dot{{\mathit{m}}_{\mathrm{\ast }}}\mathrm{~}{\mathrm{10}}^{-7}{\mathit{M}}_{\mathrm{\odot }}\mathit{/}$yr ; ν_w ~ 2000 km s^-1; Prinja et al. 1990) as the central engines of the bow-shocks, the resulting NA motions are a factor of ~5–10 lower than the observationally inferred values. Hence, we discard the idea of O-stars as the central sources of the bright NA bow shocks.

5.2. Orbital motion of the bow-shock sources

The location of the sources that we studied within the NA and the determination of their geometry with the models used here provide valuable information to determine the probability density distribution (PDF) of their orbital parameters, PDF(a,e,i,Ω,ω,P,T₀). We used the methodology described by Lu et al. (2009) and the same definitions for the orbital elements. This method consists of a Monte Carlo simulation of 10⁴ trials to convert the 3D positions and velocities of the sources into their six orbital parameters. The projected positions of the bow shocks are obtained from the astrometric measurements used to compute their proper motions (see Sect. 3). The distances of the bow shocks along the LOS and their respective uncertainties were obtained from the NA model of Paumard et al. (2004). The LOS velocities were constrained in two different ways. On the one hand, they were constrained by assuming that the sources are tied to the gravitational potential of Sgr A* and have a 3D motion lower than their escape velocities. This method provides an upper limit for the magnitude of the LOS velocity but with an ambiguity in its sign. The PDFs of the orbital parameters were therefore constructed using both signs for the LOS velocities. The second method computes the LOS velocity using the relative velocity of the sources through the ISM obtained from our model fitting (see Sects. 4.2, 5.1). With this method we determined an allowed 3σ range of LOS velocities. Furthermore, this approximation helped us eliminate the degeneracy of sign presented in the derivation of the LOS velocity using the first method. The computed escape velocities, the 3D positions, the proper motions, and the 3σ ranges for the LOS velocities are reported in Table 10. The 3D positions and proper motions use a Gaussian distribution as prior, while the LOS velocity uses a uniform distribution. These parameters allowed us to determine the PDF(a,e,i,Ω,ω,P,T₀) of each one of our sources.

After obtaining the PDFs of the orbital parameters, we computed the orientation of the sources’ orbital planes in the sky, as viewed from Sgr A*. The orientation of the orbital planes can be described by a unitary vector () originating at the position of Sgr A* and normal to the source’s orbital plane. To determine the Cartesian components of , the PDFs of the inclination angle, i, and the ascending node, Ω, were used according to the following expression: $\mathit{s̅}\mathrm{=}\left(\begin{array}{c}\\ {\mathit{s}}_{\mathit{x}}\\ {\mathit{s}}_{\mathit{y}}\\ {\mathit{s}}_{\mathit{z}}\end{array}\right)\mathrm{=}\left(\begin{array}{c}\\ \sin \hspace{0.17em}\mathit{i}\hspace{0.17em}\hspace{0.17em}\cos \hspace{0.17em}\mathrm{\Omega }\\ \mathrm{-}\sin \hspace{0.17em}\mathit{i}\hspace{0.17em}\hspace{0.17em}\sin \hspace{0.17em}\mathrm{\Omega }\\ \mathrm{-}\cos \hspace{0.17em}\mathit{i}\end{array}\right)\mathit{.}$(8)Figure 11 displays the pointing vectors obtained with the PDFs of the orientation of the sources’ orbital planes using a Hammer projection. The left column of the maps shows the vectors of the PDF(i,Ω) of the bow shocks using the computed 3σ LOS range. The right column shows the vectors of the PDF(i,Ω) using only the 3σ LOS upper limit determined from the escape velocity. For IRS 10W we only computed its probability distribution using the 3σ LOS upper limit since this source was not clearly identified as a bow shock. Additionally, we included the PDF(i,Ω) of IRS 16CC, a source previously studied by Lu et al. (2009), to check the validity of our code. We obtain similar results to those previously published, which confirms the accuracy of our code. Each map in Fig. 11 displays the 1σ and 2σ confidence levels of the PDFs in red and blue. The confidence levels were computed using only orbital solutions tied to Sgr A*. The ellipses plotted in different colors represent the 1σ contour of the orbital planes of some of the main dynamical structures at the GC. Table 11 shows the peak values of the PDFs for the sources identified as bow shocks.

The plots in the left and right columns of Fig. 11 clearly show that the PDFs based on the modeling of the bow-shocks constrain the orbital planes significantly better than those obtained with the constraints from the escape velocities alone. The PDF(i,Ω) of IRS 5 is the best constrained. The derived LOS velocity range of IRS 5 fully eliminates the ambiguity in the sign of the LOS velocity. Conversely, IRS 1W exhibits very similar PDFs using both methods. This effect can be explained by the fact that we only took into account orbits gravitationally bound to Sgr A*. Furthermore, most of the IRS 1W bound orbits obtained from our MC simulation have LOS velocities within the 3σ range derived from our bow-shock modeling. From this analysis we also observe that IRS 5 and IRS 21 clearly are not part of the CWS or the CCWS. IRS 1W and IRS 10W are relatively close to the position of the CWS but still on the edge of the CWS 1σ limits. Hence, we cannot establish with certainty whether they are CWS members or not. None of the sources are moving in the same plane of the NA.

6. Discussion and conclusions

The AO facilities and the sparse aperture masking mode of NACO at the ESO VLT allowed us to perform a NIR multiwavelength campaign to study with high-angular resolution the structure, circumstellar matter, and kinematics of some of the brightest and most embedded sources in the central parsec of the GC. We used AO-deconvolved and beam-restored images in H and K_s and reconstructed SAM images in L′ to determine the morphology of IRS 5, IRS 1W, IRS 10W, and IRS 21. We presented a new technique for calibrating SAM data of dense star fields that does not require observations of an off-target calibrator source. Our new procedure, called the synthesized calibrator, uses our a priori information of the stars in the field and constructs a median superposition of bona-fide point-like sources to build the final calibrator. This method allowed us to minimize the noise, poor background estimates, and problems in the PSF stability. It is also an alternative way of calibrating SAM observations where traditional calibrators are too dim or located far away from the science target, like in the case of the GC. Additionally, we used AO observations in K_s between 2002 to 2011 to derive the proper motions of our sources.

We interpreted our observations with 3D steady-state bow-shock models and confirmed that three of them (IRS 5, IRS 1W, and IRS 21) are bow shocks associated with the relative motion of Wolf-Rayet stars through the ISM. Although this conclusion was suggested previously, our models give the best constraints on the bow-shock shell and kinematics so far because of our detailed modeling, new proper motion data, and the inclusion of new NA proper motion data from Zhao et al. (2009), which allowed us an additional consistency check. In contrast to Tanner et al. (2005), we also found a good model fit to the IRS 1W. On the other hand, the extended emission around IRS 10W does not appear to trace a bow shock, but may be produced by the shock of two stellar winds, or there may be a chance superposition with a mini-spiral feature, or the star may heat a nearby mini-spiral filament. We also confirmed that IRS 21 is an edge-on bow shock, not a face-on bow shock, as was previously hypothesized by Tanner et al. (2005).

Our results support the idea that the observed NIR bow-shock emission is thermal emission from hot dust grains. Tanner et al. (2002) fit a 1000 K black-body to the spectral energy distribution of IRS 21. Dust can easily reach such temperatures in the local environs of massive (O/WR) stars since they are strong UV radiation emitters. For our bow-shock sources, which have IR luminosities of L ~ 10⁴–10⁵L_⊙ (Tanner et al. 2005), the smallest distance from the central source at which dust can exist without evaporation, r_evap, can be calculated, at first glance, as (Barvainis 1987) ${\mathit{r}}_{\mathrm{evap}}\mathrm{=}\mathrm{1.3}\sqrt{{\mathit{L}}_{\mathrm{UV}\mathit{,}\mathrm{39}}}\hspace{0.17em}{\mathit{T}}_{\mathrm{1500}}^{-2.8}\hspace{0.17em}\mathrm{pc}\mathit{,}$(9)where L_UV,39 is in units of 10³⁹ watts and T₁₅₀₀ in units of 1500 K. Hence, r_evap for our sources is ~8 × 10^-5–25 × 10^-5 pc or 16–52 AU. This lies well inside the measured size of the bow-shock shells of our targets at scales lower than our angular resolution. Taking into account the grain sizes used in our study, we can also compute r_evap using Eq. (4). In this case, the smallest distance at which dust can exist is about 24–75 AU for the given luminosities. Both approximations are consistent, also with the observed bow-shock sizes. From Eq. (4) we obtain a rough estimate of 470 K for the temperature near the apices of the bow shocks. These calculations, together with our scattering and thermal-emission modeling, support the assumption that the radiation in the GC bow-shock sources arises from thermal emission in dust grains. Therefore, the mechanism responsible for the observed polarization is the alignment of asymmetric/elongated dust grains. The alignment mechanism implies a compression of the magnetic field of the Northern Arm in the environs of the bow shocks and the increment of the local field density. These changes produce a rapid magnetic alignment of the grains (Buchholz et al. 2013). From the best-fit bow-shock stand-off distances, we have established that their central sources are massive objects (WR stars) with supersonic stellar wind velocities and strong mass-loss rates. They are therefore short-lived objects and have probably formed during the most recent starburst event at the GC.

One of the main sources of uncertainty in determining the orbital planes of the population of massive stars around the GC is the lack of information of the stellar LOS distances (see e.g., Lu et al. 2009; Yelda et al. 2014). Therefore, their LOS distances are significantly better constrained by the result that the extended emission of our sources are bow-shocks produced by the interaction of the stellar source with the NA. This characteristic, unique in this sources, along with their LOS velocity from our best-fit models, their proper motions and models of the ISM in the NA, allowed us to perform a Monte Carlo analysis of the sources’ orbital planes.

IRS 1W stands out from the other sources because of its unusual morphology (see Sect. 4.6), in particular its long tail, which we expect to trace the relative motion between IRS 1W and NA in the past 50 years. From studying the orbital-plane, we derived that the most probable orbit for this source peaks at an inclination of 132° with an ascending node angle of 154°. Since the NA orbital plane has an inclination of 45° with an ascending node of 15°, we discovered that both the IRS 1W and the NA orbital planes are approximately perpendicular to each other. Therefore, IRS 1W is crossing the NA and its tail traces the thickness of the NA stream. Since the projected tail has an extension of 6400 AU, the NA thickness probably is on the order of ~9000 AU at the position of this source.

The derived orbital planes suggest that our sources are not part of any of the suggested disks of young stars in the GC, but our results are consistent with the most recent detailed work on the orbital distribution of the young stars in the GC by Yelda et al. (2014). These authors suggested that only about 20% of 116 O/WR stars, at a projected distances of 0.8′′ to 12′′ from Sgr A*, are part of the CWS. Hence, if our bow shocks were members of the CWS, the number of disk members would be increased by about 20%. However, our results show that the targeted bow-shock sources probably do not form part of the CWS. Therefore, our results provide new evidence for the complexity of the population of young massive stars in the central parsec of the GC and highlight that the disk members are only a small fraction of the young massive stars around Sgr A*. We conclude that models that intend to explain the formation of the young massive stars in the central parsec of the GC and/or their dynamical evolution need to be adjusted to explain the presence of a considerable number of massive stars with apparently random distribution of their orbital planes, which suggests that either not all stars formed in a disk, or that the disk is rapidly dissolving.

Online material

Appendix A: Proper motion measurements

---

Acknowledgments

The authors thank the referee for his/her useful comments that improved this work. J.S.B., R.S. and A.A. acknowledge support by grants AYA2009-13036-CO2-01, AYA2010-17631 and AYA2012-38491-CO2-02 of the Spanish Ministry of Economy and Competitiveness (MINECO) cofunded with FEDER funds, and by grant P08-TIC-4075 of the Junta de Andalucía. R.S. acknowledges support by the Ramón y Cajal program of the Spanish Ministry of Economy and Competitiveness. J.S.B. acknowledges support by the JAE-PreDoc program of the Spanish Consejo Superior de Investigaciones Científicas (CSIC) and to the ESO studentship program. Part of this work was supported by the COST Action MP0905: Black Holes in a violent Universe. This work is based on observations made with ESO Telescopes at the La Silla Paranal Observatory under programs 071.B-0077(A), 073.B-0084(A), 073.B-0085(E), 077.B-0014(A), 081.B-0648(A), 179.B-0261(X), 183.B-0100(D), 183.B-0100(G), 183.B-0100(V), and 085.D-0214(A).

References

All Tables

All Figures