Principal component analysis tomography in near-infrared integral field spectroscopy of young stellar objects. I. Revisiting the high-mass protostar W33A
Felipe Navarete, Augusto Damineli, João E. Steiner, Robert D. Blum
aa r X i v : . [ a s t r o - ph . S R ] F e b MNRAS , 1–19 (2021) Preprint 5 February 2021 Compiled using MNRAS L A TEX style file v3.0
Principal component analysis tomography in near-infrared integralfield spectroscopy of young stellar objects. I. Revisiting thehigh-mass protostar W33A
F. Navarete ★ , A. Damineli , J. E. Steiner † , R. D. Blum Universidade de São Paulo, Instituto de Astronomia, Geofísica e Ciências Atmosféricas, Rua do Matão 1226,Cidade Universitária São Paulo-SP, 05508-090, Brasil NSF’s Optical–Infrared Astronomy Research Laboratory P.O. Box 26732, Tucson, AZ 85719, USA
Accepted 2021 February 4. Received 2021 February 4; in original form 2020 November 29
ABSTRACT
W33A is a well-known example of a high-mass young stellar object showing evidence of acircumstellar disc. We revisited the 𝐾 -band NIFS/Gemini North observations of the W33Aprotostar using principal components analysis tomography and additional post-processingroutines. Our results indicate the presence of a compact rotating disc based on the kinematicsof the CO absorption features. The position-velocity diagram shows that the disc exhibits arotation curve with velocities that rapidly decrease for radii larger than 0 . ′′ ∼
250 AU) fromthe central source, suggesting a structure about four times more compact than previouslyreported. We derived a dynamical mass of 10.0 + . − . M ⊙ for the “disc+protostar” system, about ∼
33% smaller than previously reported, but still compatible with high-mass protostar status. Arelatively compact H wind was identified at the base of the large-scale outflow of W33A, witha mean visual extinction of ∼
63 mag. By taking advantage of supplementary near-infraredmaps, we identified at least two other point-like objects driving extended structures in thevicinity of W33A, suggesting that multiple active protostars are located within the cloud.The closest object (Source B) was also identified in the NIFS field of view as a faint point-like object at a projected distance of ∼ K -bandcontinuum emission detected in the same field. Another source (Source C) is driving a bipolarH jet aligned perpendicular to the rotation axis of W33A. Key words: methods: statistical – techniques: imaging spectroscopy – stars: protostars – stars:pre-main sequence – ISM: jets and outflows – ISM: individual: W33A
In the last few decades, integral field unit (IFU) observations havebecome crucial to understand the kinematics of very complex ex-tended objects over a variety of scales, from protostellar objects andprotoplanetary discs (e.g. Murakawa et al. 2013; Alves et al. 2019)to galaxies (Kewley et al. 2019). The powerful three-dimensionaldata cubes produced by IFU spectrographs are of increasingly largesize, containing an enormous amount of information and noise thatmust be pulled from the cube by sophisticated means. For example,a typical near-infrared data cube may have millions of pixels, andgiven the complexity of the spatial and spectral information, oftentimes much of the information is hard to extract and may be simplyignored in the analysis in preference of a restricted subset of the data(line maps and their products, spectral analysis of bright regions, ★ E-mail: [email protected] (FN) † In Memoriam etc.). Frequently, a large fraction of the data in these data cubescontains just instrumental or sky noise. Thus, robust procedures toextract information from the large data sets are increasingly requiredin astrophysics.In this context, Steiner et al. (2009) developed a method forthe use of the Principal Component Analysis (PCA) Tomography inoptical and near-infrared IFU data cubes. PCA is a non-parametrictechnique, that is, it does not require any assumption or expecta-tion of the physical and geometrical parameters of the data. PCAcompresses the data expressed as a large set of correlated variables(spatial and spectral pixels) in a small but optimal set of uncorrelatedand orthogonal variables (i.e. the eigenvectors and their spatial pro-jections), ordered by their relative importance to explain the initialdata set.PCA tomography has been recently used in extra-galactic stud-ies (e.g. Ricci et al. 2011; Menezes et al. 2013; Ricci et al. 2014,2015; Dahmer-Hahn et al. 2019), but no application to stellar astro-physics cases has been made yet. Thus, we applied this technique to © Navarete et al. the K -band IFU observations of the well-known high-mass protostarW33A investigated by Davies et al. (2010) to check if the PCA To-mography is able to recover their findings. W33A is a well-knownyoung stellar object (YSO) with indications of an active accretionprocess ongoing. This protostar is associated with the Red MSXSource object G012.9090 − ( 𝐿 bol / 𝐿 ⊙ ) = .
51 (Lumsden et al. 2013) and located at thedistance of 2.4 kpc (Immer et al. 2013).Davies et al. (2010) presented the study of this source throughthree-dimensional spectroscopy in the 𝐾 -band (2.2 𝜇 m) usingNIFS/Gemini North IFU observations. Those authors reported theidentification of the key-structures predicted by the disc-mediatedaccretion scenario (Shu et al. 1987): a circumstellar disc traced bythe CO features in absorption; a jet-like structure traced by the H emission at 2.12 𝜇 m; and the cavity of the jets close to their launch-ing point traced by the spectro-astrometry of the Br 𝛾 emission atmilli-arcsecond scales. Those authors derived the Keplerian massof ∼
15 M ⊙ for the “protostar + disc” system, assuming a distanceof 3.8 kpc.Given the importance of W33A as one of the rare cases of high-mass YSOs with clear signs of active accretion, the K -band IFUobservations of the high-mass protostar W33A offer an excellentopportunity to verify the efficiency of PCA Tomography on repro-ducing the results reported by Davies et al. (2010). We are also in-terested to check how the results are sensitive to the additional post-processing routines suggested by Steiner et al. (2009) (e.g. differen-tial atmospheric refraction correction, high-frequency noise filter-ing, and instrumental fingerprint removal; see also Menezes et al.2014), and if PCA Tomography can provide additional informationbased on the correlation of the spectral features present in the NIFS K -band observations.This work is organised as follows. Section 2 presents the obser-vations, the post-processing routines of the data and the descriptionof the PCA tomography technique. Section 3 presents the resultsand interpretation of the data through PCA tomography, the deriva-tion of the physical parameters of the source and the analysis ofauxiliary near-infrared maps. Section 4 presents the discussion ofthe data and comparison with recent findings in the literature. Ourconclusions are given in Section 5. 𝐾 -band IFU observations The K -band ( 𝜆 ≈ 𝜇 m) IFU observations of W33A were ob-tained with the Near-Infrared Integral Field Spectrometer (NIFS,McGregor et al. 2003). NIFS has a full wavelength coveragein the K -band of about 4200 Å, with a linear dispersion of Δ 𝜆 disp = 2.13 Å/pixel − and a nominal spectral resolving power of 𝜆 / Δ 𝜆 disp = 5160 (Blum & McGregor 2008), or ≈
58 km s − in thevelocity scale. NIFS observes a 3 ′′ × ′′ field-of-view in a 60 × array, delivering a spatial sampling of 0 . ′′
050 pixel − .The observations were carried at the Gemini North 8-m tele-scope on the night of May 25th, 2008. The raw data of W33A weredownloaded from the Gemini Science Archive (project GN-2008A-Q-44). The data include a set of calibration frames required forprocessing the science observations: a set of flat-field images takenwith Quartz lamps on (flat-on) and off (flat-off); flat-field imagestaken with a Ronchi mask (Ronchi flats), used to correct the spatialdistortion along the spectral dimension; a set of dark images, takenwith the same exposure time of the science and telluric sources; and an Argon-Neon lamp spectrum (Ar-Ne) for wavelength calibrationof the data. The A2 V double star HR6798A+B was observed asthe telluric standard star, used for correcting the telluric absorptionfeatures, and for performing the flux-calibration of the final datacubes.The data were processed using the usual reduction procedurebased on IRAF/Gemini package. The telluric correction was ap-plied to the data using a one-dimensional template spectrum ex-tracted from the telluric standard star. The telluric Br 𝛾 feature wasremoved from the template by fitting a Voigt profile. Then, the datacubes were constructed using the Gemini NIFCUBE procedure, byinterpolating the original data into a final spatial pixel (spaxel) sizeof 0 . ′′ . ′′ ± . ′′ Additional 𝐽 ( 𝜆 ≈ 𝜇 m), H (1.64 𝜇 m), H (2.14 𝜇 m) and 𝐾 (2.22 𝜇 m) images of a 90 ′′ × ′′ FOV around the W33A protostarwere obtained from the UKIRT Widefield Infrared Survey for H (UWISH2) (Froebrich et al. 2011). The large scale maps have apixel scale of 0 . ′′ − and a mean seeing of about 0 . ′′
96 at the K -band, offering a complementary view of the stellar and nebularcontent in the vicinity of W33A. The data cubes generated by the standard
IRAF reduction procedurewere further processed using a set of routines written in
IDL (Inter-active Data Language). The purpose of these additional proceduresinclude 𝑖 ) the atmospheric refraction correction and combination ofthe data cubes, 𝑖𝑖 ) the flux-calibration, and 𝑖𝑖𝑖 ) the high-frequencynoise suppression from the data. The individual data cubes were corrected for the effects of thedifferential atmospheric refraction (DAR) using the procedure de-scribed in Menezes et al. (2014). Figure 1 compares the data beforeand after the DAR-correction based on the centroid position of thepoint-like source of W33A as a function of the wavelength for theentire spectral region covered by the K -band data cube. The plotshows the rms of the source position reduces from ∼ http://astro.kent.ac.uk/uwish2/ MNRAS000
IDL (Inter-active Data Language). The purpose of these additional proceduresinclude 𝑖 ) the atmospheric refraction correction and combination ofthe data cubes, 𝑖𝑖 ) the flux-calibration, and 𝑖𝑖𝑖 ) the high-frequencynoise suppression from the data. The individual data cubes were corrected for the effects of thedifferential atmospheric refraction (DAR) using the procedure de-scribed in Menezes et al. (2014). Figure 1 compares the data beforeand after the DAR-correction based on the centroid position of thepoint-like source of W33A as a function of the wavelength for theentire spectral region covered by the K -band data cube. The plotshows the rms of the source position reduces from ∼ http://astro.kent.ac.uk/uwish2/ MNRAS000 , 1–19 (2021)
CA Tomography of the NIFS K-band observations of the high-mass protostar W33A lines). Indeed, Menezes et al. (2014) reported that the total displace-ment introduced by the DAR effects can be larger than the typicalNIFS spatial resolution from 𝐽 - to K -band wavelengths. However,we do not expect a critical enhancement of the results introducedby the DAR correction on the analysis of single spectral lines withwidths of a few tens to hundreds of Angstroms. However, the singleline analysis benefits from the post-processing routines aimed at thesuppression of the noise from the datacube, as further presented inSect. 2.3.3 and 2.4.2.A final data cube was obtained by median combining the sixindividual DAR-corrected data cubes. We compared the residualsbetween the median-combined and an average-combined data cubeto check if the fluxes were preserved. We noted that the largestdeviations were about ∼
1% of the flux in the brightest regions ofthe FOV, closer to the W33A point source, and ∼
5% in low signal-to-noise regions, indicating no significant discrepancy between theusage of the median or the average of the fluxes.
The median-combined data cube was flux-calibrated using a custom
IDL routine based on the
STANDARD , SENSFUNC and
CALIBRATE routines of the
IRAF -NOAO package. The three corrections per-formed by each routine are applied simultaneously by computingthe ‘calibration factor’ function ( 𝐶 𝜆 ), defined as: 𝐶 𝜆 = . (cid:18) 𝑂 𝜆 𝑡 exp 𝐵 𝜆 ( 𝑇 eff , 𝑚 𝐾 𝑠 ) (cid:19) + 𝐴 𝐸 𝜆 (1)where, 𝑂 𝜆 corresponds to the observed telluric counts of the obser-vation (in ADU), 𝑡 exp is the exposure time of the observation (in sec-onds), 𝐵 𝜆 ( 𝑇 eff , 𝑚 𝐾 𝑠 ) is the black-body flux (in erg s − cm − 𝜇 m − )at a given effective temperature 𝑇 eff , scaled by the Ks -band magni-tude of the telluric standard, 𝐴 is the average airmass of the observa-tions and 𝐸 𝜆 is the atmospheric extinction curve (in mag airmass − units).We adopted the Ks -band flux of the 2MASS catalogue(Skrutskie et al. 2006) to derive the flux-calibration of the datacubes. Since the HR6798A+B system is unresolved at the 2MASScatalogue but resolved in the NIFS observations, the constructionof 𝑂 𝜆 was performed by computing the integrated spectrum ofthe resolved components of the telluric standard star over the spa-tial pixels (spaxels) with average flux above a 3- 𝜎 threshold abovethe value of the sky background. Then, the continuum emissionwas modelled by fitting a polynomial function of order 𝑛 ≤
4, us-ing the outlier-resistant IDL algorithm robust_poly_fit . Next, 𝐵 𝜆 ( 𝑇 eff , 𝑚 𝐾 𝑠 ) was constructed by following a two-step procedure.First, the shape of the black-body spectrum was obtained by usingthe effective temperature corresponding to the spectral type of thetelluric standard (A2 V), taken from Tokunaga (2000). Then, theblack-body emission was scaled using the Ks -band flux of the tel-luric star ( Ks = 5.87 mag). The 2MASS Ks -band magnitudes ( 𝑚 Ks )were converted to fluxes ( 𝐹 Ks ) using Eq 2: 𝑚 𝐾 𝑠 = − . (cid:18) 𝐹 𝐾 𝑠 𝐹 𝐾 𝑠, (cid:19) (2)where 𝐹 𝐾 𝑠, = 4.28 · − erg s − cm − 𝜇 m − is the zero-point fluxat the Ks -band (Cohen et al. 2003).Finally, the additive term of Eq. (1) requires a model for theatmospheric extinction as a function of the wavelength. Thus, weadopted the mean extinction through the broad K -band filter forMauna Kea, 𝐸 𝜆 = 0.033 mag/airmass, reported by Tokunaga et al.(2002). Figure 1.
Centroid position of the W33A source as a function of the wave-length (from 2.02 to 2.42 𝜇 m, from blue to red) of the same data cube before(top panel) and after the atmospheric refraction correction (bottom). The(0,0) position is indicated as a × symbol and the error bars corresponds tothe standard deviation of the position across the spectral dispersion. After constructing 𝐶 𝜆 , each spaxel of the data cube (here,represented as 𝐷 𝑥𝑦𝜆 ) was divided by the exposure time of theobservations, and was extinction-corrected and flux-calibrated asindicated by Eq. (3):log ( 𝐷 𝑓 𝑐𝑎𝑙,𝑥𝑦𝜆 ) = . · 𝐶 𝜆 · log (cid:18) 𝐷 𝑥𝑦𝜆 𝑡 exp (cid:19) (3)The results of the flux calibration were checked using the integratedtelluric standard spectrum as follows. First, the total flux of the tel-luric standard star was integrated over the bandwidth of the 2MASS Ks -band filter ( 𝜆 𝑐 = 2.159 𝜇 m, Δ 𝜆 = 0.262 𝜇 m, Cohen et al. 2003).Then, the corresponding Ks -band magnitude was evaluated using MNRAS , 1–19 (2021)
Navarete et al.
Figure 2.
Integrated flux of the K -band NIFS observations of W33A. Theimage is shown in the logarithm scale indicated by the colour bar in theright. The contours are placed at the 1, 10 an 50% of the peak intensity ofthe map. The W33A point source is located around the (0 . ′′
0, 1 . ′′
2) positionin the FOV.
Eq. (2), and compared with the 2MASS Ks -band magnitude of thetelluric star (5.87 mag). Their difference, Δ 𝐾 s = − Ks -band flux ratio of ≈ . ′′
0, 1 . ′′
2) and the associated K -band continuumnebulosity towards the South. NIFS data are provided in a grid of 𝑥 , 𝑦 and 𝜆 coordinates, wherethe astrophysical information is often mixed with a complex noisestructure, including the photon noise, high-frequency noise andsystematic instrumental noise that is often associated with astro-nomical data cubes (see, e.g., Menezes et al. 2014). Although thephoton noise is a well-known and inherent source of noise in astro-nomical observations, the nature of the last two components are ingeneral unknown and may be introduced during the data acquisitionor data processing (e.g. the resampling of the original 0 . ′′ × . ′′ . ′′ × . ′′
050 spatial pixels). Despite the unknown originof the noise in IFU observations, a considerable fraction of it canbe successfully suppressed from the data based on the applicationof a low-pass Butterworth filter (explained below) and the PCATomography (Sect. 2.4) on the data.In the Fourier domain, the spatial information can be expressedin terms of frequencies ( 𝜈 x , 𝜈 y ), allowing us to identify and sup-press high-frequency noise using an appropriate low-pass band fil-ter. Thus, we used a two-dimensional Butterworth filter in the spatialfrequency dimension to suppress the high-frequency noise presentin the data cube, preserving the frequencies containing the astro-physical information. The details on the construction of the Butter-worth filter can be found in Sect. 5 from Menezes et al. (2014). Weused a Butterworth filter with order 𝑛 = 2 and cut-off frequenciesdefined in terms of the Nyquist frequency of the data ( 𝑅 Ny ), corre-sponding to the highest-order sampled by the Fast Fourier Transform(FFT). Due to the different size of the original spaxels in the 𝑥 - and 𝑦 -direction of the NIFS data (0 . ′′
103 and 0 . ′′ Figure 3.
Mean Fourier Transform of the spatial dimensions of the datacube before (top panel) and after the Butterworth filtering (bottom). Thewhite ellipsoid indicates the cut-off frequency used for constructing theButterworth filter to suppress the high-frequency noise of the data. cut-off frequencies must be slightly different in each direction. Thedetermination of the cut-off frequency is a crucial step and mustbe performed carefully: if overestimated, it partially corrects thehigh-frequency noise that should be suppressed; if underestimated,it filters frequencies containing the astrophysical signal, causing anenlargement of the Point Spread Function (PSF).For W33A observations, frequencies greater than 𝑓 cutoff , x = 0.5 𝑅 Ny and 𝑓 cutoff , y = 0.7 𝑅 Ny were filtered. Fig-ure 3 presents the mean spatial FFT of the W33A data cubebefore and after the Butterworth filtering. The cut-off frequenciesdefine an ellipsoid in the Fourier space, indicated by the whitecontour on top of the mean FFT maps. We note that almost allthe high-frequency features shown in the top panel of Fig. 3 weresuppressed by the Butterworth filtering as indicated in the bottompanel of the figure.A simple subtraction of the original cube and the filtered cubecan be used to check if the Butterworth filtering was successfully MNRAS000
Mean Fourier Transform of the spatial dimensions of the datacube before (top panel) and after the Butterworth filtering (bottom). Thewhite ellipsoid indicates the cut-off frequency used for constructing theButterworth filter to suppress the high-frequency noise of the data. cut-off frequencies must be slightly different in each direction. Thedetermination of the cut-off frequency is a crucial step and mustbe performed carefully: if overestimated, it partially corrects thehigh-frequency noise that should be suppressed; if underestimated,it filters frequencies containing the astrophysical signal, causing anenlargement of the Point Spread Function (PSF).For W33A observations, frequencies greater than 𝑓 cutoff , x = 0.5 𝑅 Ny and 𝑓 cutoff , y = 0.7 𝑅 Ny were filtered. Fig-ure 3 presents the mean spatial FFT of the W33A data cubebefore and after the Butterworth filtering. The cut-off frequenciesdefine an ellipsoid in the Fourier space, indicated by the whitecontour on top of the mean FFT maps. We note that almost allthe high-frequency features shown in the top panel of Fig. 3 weresuppressed by the Butterworth filtering as indicated in the bottompanel of the figure.A simple subtraction of the original cube and the filtered cubecan be used to check if the Butterworth filtering was successfully MNRAS000 , 1–19 (2021)
CA Tomography of the NIFS K-band observations of the high-mass protostar W33A Figure 4.
Difference between the NIFS data cube before and after the spatialButterworth filtering. The residuals are shown in units of 𝜎 , in a red-to-bluecolour-scale to highlight the details of the noise structure. The contourscorresponds to the 1, 10 and 50% of the peak flux in the original data cube. applied. Figure 4 exhibits the residuals between the average of thefiltered and the original data cubes. Since the noise is proportional tothe signal-to-noise ratio of the data, the residuals are larger towardsthe point-like source in the FOV. Despite this, the residuals exhibit anoisy pattern over the entire FOV. Even though most of the residualscorrespond to only a fraction of the noise of the datacube, the noisesuppression makes a significant difference for the analysis of weakspectral features (e.g. the CO absorption features at 2.3 𝜇 m). This Section presents a brief description of the Principal Compo-nent Analysis (PCA) Tomography technique and its interpretation.Further details can be found in (Steiner et al. 2009) and referencestherein.The PCA is a statistical procedure that transforms a set ofmulti-dimensional variables which are mutually dependent (i.e. cor-related) into a set of orthogonal and linearly uncorrelated variables.The PCA Tomography is applied to the data cubes in order to quan-tify and extract information, ordered by the fraction of the varianceof these new variables: the first of them explains the highest amountof the data variance, and the higher-order variables explain suc-cessively smaller amounts of the data variance. In the transformedcoordinate system, the spectral pixels of the data cube correspondto the variables, and the spaxels are the observables. The new set ofuncorrelated variables – the eigenspectra – are obtained as a func-tion of the wavelength. The projection of the observables (that is,the spaxels) on the eigenspectra are images, referred as tomograms.The new coordinates of the data are called Principal Compo-nents (PCs), and each PC is orthogonal (i.e. uncorrelated) to theothers. Each principal component (PC) is composed of: 𝑖 ) an eigen-value, 𝐴 𝑖 , which indicates the relevance of the PC in explaining theinformation contained in the original data cube; 𝑖𝑖 ) a tomogram, 𝑇 𝑖 ( 𝑥, 𝑦 ) , corresponding to the projection of the PC on the spatialdimensions of the image; and 𝑖𝑖𝑖 ) an eigenspectrum, 𝐸 𝑖 ( 𝜆 ) , which the routines can be downloaded in Table 1.
Scree test of the PCA tomography for the K -band data cube ofW33A. N 𝑖 Variance Cumulative(%) Variance (%)1 99.80360 99.80362 0.18817 99.99183 0.00253 99.99434 0.00109 99.99545 0.00042 99.99586 0.00025 99.99607 0.00021 99.99628 0.00013 99.99649 0.00008 99.9965. . .260 5.06 · − · − is the projection of the PC on the spectral direction. These threecomponents must be analysed together to properly evaluate the in-formation contained in each PC. The maximum number of PCscorresponds to the number of the spectral elements of a given datacube.The PCA also works as a linear dimensional reduction tech-nique by compressing the data, expressed initially as a large set ofcorrelated variables, into a small but optimal set of uncorrelatedvariables, ordered by their eigenvalues. Hence, most of the infor-mation of the original dataset can be explained with a relativelysmall number of PCs. The analysis of the relative value of 𝐴 𝑖 as afunction of the PC number – the so-called Scree test – allows us todetermine the best number of PCs that explains the initial data (seeSect. 2.4.1).Due to the orthogonality of the PCs, the original data cubecan be reconstructed by considering a smaller set of PCs with in-significant losses of its astrophysical information. This procedureallows us to suppress a significant amount of noise which is oftenassociated with high-order PCs (see Sect. 2.4.2). Each Principal Component is ordered by its eigenvalue (largest tosmallest), 𝐴 𝑖 , which indicates the relevance of the PC in explainthe information contained in the data cube. In other words, theeigenvalue is a direct measure of the variance contained in thePC. The fractional variance of each principal component 𝑖 , 𝑉 𝑖 , isexpressed as the ratio between its eigenvalue 𝐴 𝑖 and the sum of theeigenvalues of all the PCs, given by: 𝑉 𝑖 = 𝐴 𝑖𝑁 Õ 𝑖 = 𝐴 𝑖 (4)Table 1 lists the fractional and cumulative variance for selectedPrincipal Components of the PCA Tomography analysis of W33A.We note that the first PC contains about 99.8% of the total variance ofthe data cube, while higher-order PCs exhibit successively smallerfractional variance values. For W33A, the maximum number ofcomponents is equal to 1782, which is the number of spectral pixelsof the data cube.The analysis of the fractional variance alone does not providea useful way to determine how many PCs contains astrophysical MNRAS , 1–19 (2021)
Navarete et al.
Figure 5.
Scree test showing the variance of each Principal Component (PC)versus the Principal Component order. The variance of each PC is presentedas a percentile of the total variance of the data cube (see Eq. (4)). A dashedline indicates the noise level, based on the linear fit of the data associatedwith PCs of order 10 or higher. information in the datacube. This diagnosis can be assessed byconstructing the Scree test of the PCA, shown in Fig. 5, whichpresents the distribution of the fractional variance as a functionof the order of the PC. The Scree test shows a rapid decrease ofthe variance as a function of the number of the PC, exhibiting twodistinct regions: the signal-dominated
PCs (from PC1 to PC5) fromwhich the variances are significantly larger than the noise level, andthe noise-dominated
PCs (PC6 and onward). The noise level wasestimated by performing a linear fit of 𝑉 𝑖 by 𝑖 for the componentsdominated by noise (PC ≥ The first 5 PCs contain about ≈ The analysis of the large-scale K -band emission of W33A was pre-sented by Davies et al. (2010, hereafter D10), exhibiting the embed-ded protostellar object W33A and its extended K -band nebulosity(see their Fig. 1). We took advantage of the available 𝐽𝐻𝐾 and H maps of the UWISH2 survey (Froebrich et al. 2011) to present amore complete view of the infrared content associated with W33A, including any large-scale H emission at 2.12 𝜇 m, a clear signatureof star-forming activity for a broad range of bolometric luminosityvalues (e.g. Varricatt et al. 2010; Navarete et al. 2015).Figure 7 presents a false-colour 𝐽𝐻𝐾 map of a 1.5 arcmin region around W33A from the UWISH2 data, overlaid by contoursof the continuum-subtracted H emission. The top panel, shows the K -band reflection nebula illuminated by W33A, oriented in the SE-NW direction. The weakness of the NW lobe compared to the SEis then due to the orientation of the system combined with a strongextinction gradient from the (SE) blue- to the (NW) red-shifted lobeof the K -band outflow (in Sect. 3.4.2, we show that the H emissionidentified is blue-shifted in the NIFS FOV, and the presumptivered-shifted lobe is off the NIFS FOV to the NW). The NIFS FOVis indicated by the white box, the position of the K -band compactsource is indicated by the black × symbol, and the position angle(PA = 150 ◦ ) of the disc reported by D10 is shown by the blue arrows,which is aligned with the K -band nebulosity.The white contours on top of the NIR map indicate thecontinuum-subtracted H emission. The H contours indicate noclear extended emission aligned with the main K -band nebulosityof W33A. Instead, they are tracing a set of H knots roughly alignedin the E-W direction, as previously identified by Lee et al. (2013)(see their Fig. 2). At least six bright H knots (labelled from “a” to“f”) were identified as arising from a bipolar jet. The jet has a pro-jected length of ∼ ′′ , corresponding to a linear size of 0.58 parsecsat the distance of W33A (2.4 kpc). The large-scale H jet is offsetby a few arcseconds to the S direction of W33A, suggesting that atleast another active protostar is located in the vicinity of W33A.To better constrain the origin of the H jet, a zoomed in view ofa 24 ′′ × ′′ region around W33A is shown in the bottom panel ofFig. 7. Contours of the 𝐽 -, 𝐻 -, and 𝐾 -band filters (blue, green, andred, respectively) are overlaid to highlight the structures around theW33A source, including two additional point-like sources locateda few arcseconds to the S of W33A (labelled as B and C), and anextended 𝐻 -band emission oriented to the NW direction. Source B(RA = 18:14:39.584, Dec. = − 𝐻 -band emission oriented to the NW direction, overlap-ping the 𝐾 -band nebular emission arising from W33A. Source C(RA = 18:14:39.606, Dec. = − 𝐻 -band, located at the centralposition of the line connecting the H knots “d” and “e”, suggestingthis object is actually launching the large-scale H knots indicatedby the white contours. Despite the fact that Sources B and C do notexhibit unequivocal K -band counterparts in the large-scale maps,any K -band emission associated with Source B and its nebulositymight still be detected in the NIFS FOV.In Fig. 8 we present the NIFS FOV overlaid by contours ofthe 𝐽𝐻𝐾 filters. Thanks to the relatively small plate scale of theNIR maps (0 . ′′ − ), we were able to extract maps with thesame size as the NIFS FOV, allowing us to compare the high-resolution NIFS data (0 . ′′
15) with coarser resolution data taken atshorter wavelengths. The map exhibits the point-source of W33Aat the top central region and the main K -band nebulosity to theSouth. The comparison between the black NIFS K -band contoursand the contours of the 𝐻 - and 𝐾 -band (green and red, respectively)suggests that the main nebulosity may be a superposition of distinctstructures within the NIFS FOV: it is partially arising from W33Aat the upper region of the FOV, while the K -band nebula is bettercorrelated with the extended 𝐻 -band emission (green contours)towards the South.The relatively blue colour of the main nebulosity when com-pared to W33A also favours that the extended emission is brighter MNRAS , 1–19 (2021)
CA Tomography of the NIFS K-band observations of the high-mass protostar W33A Figure 6.
Residuals between the original W33A data cube and the reconstructed data cube using the first 5 PCs. The top panels present the spatial distributionof the residuals of individual spectral channels indicated by the vertical filled red lines of the mean spectrum shown on the bottom panel. at shorter wavelengths, as expected for the extended 𝐻 -band emis-sion detected in the NIFS FOV. At the bottom left (SE directionfrom W33A), a faint K -band point-like source peaks at the brightestand central region of the 𝐻 -band contours, matching the expected K -band counterpart of Source B identified in the 𝐽𝐻𝐾 map (seebottom panel of Fig. 7).A detailed inspection in the vicinity of W33A confirms thatmultiple (proto)stars are located at the line-of-sight of the main K -band nebulosity of W33A. The complementary information fromthe 𝐽𝐻𝐾 maps indicates that both the compact emission arising fromSource B and its extended 𝐻 -band emission might contribute to theextended structures identified in the NIFS FOV, as we cannot be surewhether such emission comes from the scattered light in the outflowcavity or from source B. In addition, the 𝐽 − band counterpart seemsto be a chance alignment of an optically visible star as indicated byDSS images (not shown in the present work). In this section, we briefly analyse the overall properties of the K -band NIFS observations of W33A to guide the reader through thePCA Tomography analysis presented in Sect. 3.3.The left panel of Figure 9 presents a composite RGB map ofthe NIFS FOV based on the integrated emission of a 0.4 𝜇 m band-pass centred at 2.05 (blue), 2.12 (green), and 2.35 𝜇 m (red). Thered and blue contours correspond to the regions used for integratingthe spectrum of the W33A protostar and the extended nebulosity,respectively, shown in the right panel. The spectra shown in the rightpanel were integrated over the regions indicated by the contours: theW33A protostar (red), the main nebulosity (blue) and the spaxelscorresponding to Source B (in green). We excluded the contribution of Source B (the corresponding spaxels are flagged by the greencontour in the left panel) from the nebular emission. The K -bandspectrum of Source B (shown as green) is about ∼
100 times fainterthan the extended nebula (in blue), and does not significantly affectthe integrated spectrum of the brighter region. The spectrum ofSource B will be further discussed in Sect. 3.4.1.The integrated spectra indicate that the compact source shows asteeper increase in flux towards larger wavelengths when comparedto the nebulosity, as expected for enshrouded protostellar objects.In addition, the nebulosity is brighter than the protostar at shortwavelengths, suggesting a relatively strong local reddening effectby the dust towards W33A, that decreases at larger distances fromthe source. This may be due to the viewing geometry. Since thecompact source is likely seen through the obscuring disc/torus,it appears more reddened than the cavity of the outflow (i.e. thenebulosity) we see in scattered light. The low- 𝐽 CO absorptionlines support this picture since they are strong toward the pointsource and missing from the nebular spectrum as expected if densecool gas surrounds the point source.Both the protostar and the nebulosity shows the presence ofthe Br 𝛾 feature in emission at 2.166 𝜇 m, the first four CO bandheadfeatures in emission at 𝜆 𝜇 m, and the R- and P-branches ofthe CO absorption features between 2.32 and 2.37 𝜇 m. A weakemission of the Na i doubled at 2.21 𝜇 m is also observed in bothspectra. At least five ro-vibrational transitions of the H moleculeare only detected towards the extended nebula, which is also roughlyoriented in the direction of an extended H emission driven by W33Aas reported by D10. MNRAS , 1–19 (2021)
Navarete et al.
Figure 7.
Top panel: False-colour
𝐽 𝐻 𝐾 map (blue: 𝐽 , green: 𝐻 , red: K ) around W33A from the UWISH2 data. The position of the NIFS mapis indicated by the white box. W33A is indicated by the black × sym-bol, corresponding to the (0,0) position in the map (RA = 18:14:39.53,Decl. = − H image at the levels 𝜎 𝑛 (where 𝑛 = , , ... ). Knots of H emission are labelled in white, and the whitedashed line connects the knots labelled as “d” and “e”. Bottom panel: Acloser view of the 𝐽 𝐻 𝐾 map in a 24 ′′ × ′′ box around the W33A. The RGBmap is overlaid by contours of each near-infrared filter (blue: 𝐽 , green: 𝐻 ,red: K ) to highlight the details closer to the brightest region of the main neb-ula. The H knots “d” and “e” are indicated by the white contours, togetherwith a dashed white line connecting them. Additional point-like sources areindicated by black × symbols and labelled as B and C, respectively. 𝐾 -band NIFS data cube In this section, we present the analysis and interpretation of thefirst five Principal Components of the PCA tomography applied tothe K -band data cube of W33A. The tomograms and eigenspectraassociated to these PCs are shown in Figs. 10 to 15. Figure 8.
False colour composite map of the NIFS field-of-view (blue:2.05 𝜇 m, green: 2.20 𝜇 m, red: 2.35 𝜇 m) overlaid by the 𝐽 𝐻 𝐾 contours (inblue, green and red, respectively). The
𝐽 𝐻 𝐾 contours are placed at 25, 40,75, and 90% of the peak intensity of each map. For comparison, the blackcontours delineate the K -band emission of the NIFS observations. The mainstructures are labelled on the image. PC1
In general, the overall spatial and spectral information from PC1roughly correspond to projections of the integrated properties ofthe data cube (i.e. its mean in the spatial and spectral dimensions).However, these projections are not given in positive or negativefluxes, but rather in correlations and anti-correlations. Such a strongresemblance to the data is expected since PC1 contains the largestrelative variance of the data cube (explaining roughly 99.8% of thedata; see Table 1.The tomogram of PC1 exhibits a peak at the position of thebright point-like source corresponding to the W33A protostar. Themain nebulosity is also identified in the tomogram, exhibiting largercorrelation at the vicinity of W33A and decreasing towards theSouth. The correlation between the compact and the extended emis-sion (i.e. the fact that both appear as positive variances in the to-mogram) indicates that both structures shares the same spectralinformation shown in the eigenspectrum, shown in the top rightpanel of Fig. 10.The eigenspectrum shows that the variance of the continuumemission increases towards larger wavelengths. Such behaviour (inthe flux spectrum) is often associated with thermal emission of arelatively hot source (dust), as expected for a protostar. Discretespectral features such as the Br 𝛾 , the CO bandhead features, and arelatively weak Na i doublet are also observed with positive varianceand, thus, correlated with the continuum emission. In addition, thelow- 𝐽 CO absorption lines are identified as absorption features.
PC2
The higher-order PCs present the correlation and anti-correlationbetween the spectral features in the eigenspectra, and their associa-tion with the structures probed by the tomograms. Thus, for betterhighlighting the variance gradient between distinct structures, thetomograms are shown in a red-to-blue scale, where correlated (pos-itive variance) and anti-correlated regions (negative variance) areshown in blue and red, respectively.
MNRAS , 1–19 (2021)
CA Tomography of the NIFS K-band observations of the high-mass protostar W33A Figure 9.
Left panel: False-colour RGB image of the NIFS FOV (left, blue: 2.05 𝜇 m, green: 2.12 𝜇 m, red:2.35 𝜇 m) using the flux-calibrated data cube. Thecontours indicate the regions where the integrated K -band spectra (right panel) were extracted: the W33A protostar (in red), the main nebulosity (in blue) andthe compact emission corresponding to Source B (in green).. Right: Integrated K -band spectra of each region indicated in the FOV. The spectrum of Source B(in green) was multiplied by a factor of 50. No nebular contribution was subtracted. The main spectral transitions present in the data are indicated by thevertical dashed lines: the H transitions, the Br 𝛾 feature at 2.16 𝜇 m, the Na i doublet at 2.21 𝜇 m, the CO bandhead features in emission at 𝜆 𝜇 m, and theCO absorption features at 𝜆 ∼ 𝜇 m. Figure 10.
The tomogram (left panel) and the corresponding eigenspectrum (right) of the Principal Component 1 of the PCA Tomography of W33A. Thevariance of the tomogram is indicated in a rainbow colour scale, where black indicates weak correlation (min), and red indicates strong correlations (max).The contour-levels are chosen to highlight the extended emission and the point-like source. The eigenspectrum exhibits the variance ( 𝐸 ) as a function ofwavelength. The position of spectral features guiding the analysis of the data are labelled and indicated by vertical dashed lines. In the top left panel of Fig. 11, the tomogram of PC2 exhibitsthe extended nebulosity with positive variance (shown in blue),anti-correlated with a ring-like structure around W33A exhibitingnegative variance values (in red). The eigenspectrum shows contin-uum and discrete spectral features that are related to the structuresobserved in the tomogram. The continuum emission shows pos-itive variance values for 𝜆 . 𝜇 m, and is negative for longerwavelengths. Combining the spectral and spatial information, weinterpret the continuum emission as dominated by the main nebu-losity at short wavelengths, while the structure around the compactsource is brighter at longer wavelengths. This behaviour can also beinterpreted in terms of the local extinction of the dust: the redden-ing is stronger towards the compact source than that observed in theline-of-sight of the nebulosity. This interpretation is also consistentwith the RGB image shown in Fig. 8.The top right and bottom panels of Fig. 11 also exhibits the H transitions, Br 𝛾 and Na i lines correlated with the blue continuum,indicating that these lines are produced in the main nebulosity andabsent in the compact source. Davies et al. (2010) presented spectraof the individual bright regions shown in the nebulosity of Figure 11(see their Figs. 3 and 10). The three bright regions in the cavityfrom PC2 correspond to the regions defined as B, C and D byD10: B is closest to the point source and shows weaker low- 𝐽 COabsorption strength than the point source, and C and D show yet weaker absorption of the narrow CO transitions. This is consistentwith the spots showing a reflected spectrum of the point source.The narrow CO absorption lines detected at 𝜆 & 𝜇 m in theeigenspectrum appear as anti-correlated features, just as the contin-uum emission at longer wavelengths. The CO absorption featuresare produced in a relatively dense and cold environment, wherethe CO molecules are shielded from the radiation of the protostar.These physical parameters are often observed in the outer regionsof circumstellar discs, consistent with their association with theanti-correlated structure around the W33A protostar. The absenceof the CO bandhead features in the eigenspectrum of PC2 is alsoremarkable. These features are produced at the inner region of thecircumstellar disc, at temperatures of around 2,000 K (Blum et al.2004), unresolved at the spatial resolution of the NIFS observations( ∼ . ′′
15, or 360 AU at 2.4 kpc). Indeed, the PC2 tomogram exhibitsa ring-like structure around the W33A protostar, with variance val-ues dropping to zero towards its central region, which is consistentwith the unresolved origin of the CO bandhead emission.
PC3
The tomogram of PC3 exhibits a spatially resolved structure aroundthe compact source, with coherent changes between positive (blue)and negative variance values (red) at the central position of the point-
MNRAS , 1–19 (2021) Navarete et al.
Figure 11.
The tomogram (top left) and the corresponding eigenspectrum (top right) of the Principal Component 2 of the PCA Tomography of W33A. In theleft panel, the contours of the PC1 tomogram are overlaid as black curves. For orders 𝑛 ≥
2, the tomogram is shown in a red-to-blue colour scale indicated inthe top right, where blue indicates the strongest correlations and red indicates the strongest anti-correlations. In the right panel, the eigenspectrum associatedwith PC1 is overlaid as the blue curve. The horizontal dashed line indicates variance equal to zero. A detailed view of the variances of selected spectral linesof the eigenspectrum are shown in the bottom panels (from left to right: the H line at 2.12 𝜇 m, the Br 𝛾 feature at 2.16 𝜇 m, the CO (2-0) bandhead emissionat 2.29 𝜇 m, the CO absorption features at 2.32 𝜇 m, and the Na i doubled at 2.20 𝜇 m). like object. In extra-galactic objects, this pattern is often identified inrotating structures, such stellar discs or torus around super-massiveblack holes (e.g. Ricci et al. 2014). In the context of protostellarobjects, it is likely arising from a rotating circumstellar disc, asidentified by Davies et al. (2010) when analysing the CO absorptionfeatures (see their Fig. 14). The rotation axis of the structure has aposition angle of about − ◦ (from N to E).The eigenspectrum of PC3 exhibits no variance of the contin-uum, but shows a series of narrow features for 𝜆 & 𝜇 m, likelyarising from the low-J CO lines (see details in the bottom panel ofFig. 12). The blue-shifted components are correlated and so, asso-ciated with the positive variance values in the tomogram (shown asblue) while the red-shifted components are associated with negativevalues (shown as red), indicating the disc is rotating from NE to SWdirection.A similar kinematic pattern is also observed for the Br 𝛾 feature(also indicated in the bottom panel of the figure), suggesting thatthe rotating disc exhibits an ionised component (see discussion inSection 4.4).Figure 13 shows a closer view of the rotating structure identi-fied in the tomogram. The rotation axis of the disc is observed witha position angle of PA = 140 ± ◦ . A radial profile of the variancealong the disc plane is presented in the bottom panel of the figure,showing that the disc kinematic signature is only observed withinthe inner ± . ′′ ∼ PC4
The tomogram of PC4 exhibits a correlated structure roughly per-pendicular to the rotating disc identified in PC3 (see Fig. 12).The eigenspectrum of PC4 shows both continuum and discretefeatures, as observed in PC2 (Fig. 11). The continuum emissionbetween 2.17 and 2.35 𝜇 m is anti-correlated with the emission atshorter and at longer wavelengths. We interpret this result as if thesame region is associated with at least three distinct processes: 𝑖 ) the dust content is reflecting the radiation emitted by the protostar atshort wavelengths; 𝑖𝑖 ) the thermal emission of the dust is observed atlonger wavelengths; 𝑖𝑖𝑖 ) the same region exhibits the CO bandheadfeatures in emission, and the profiles are quite narrow when com-pared to the emission observed in PC1. This may indicate that weobserve this emission moving perpendicularly to the line-of-sight.In addition, the eigenspectrum of PC4 shows a strong Br 𝛾 feature in emission, suggesting the presence of partially-ionisedgas that could also contribute as a source of scattering electronsproducing the reflected light at shorter wavelengths. PC5
The tomogram of PC5 exhibits an extended emission with a differentmorphology to the main nebulosity of W33A (shown by the blackcontours), suggesting that this structure is either produced by adistinct mechanism or exhibits distinct physical conditions.Indeed, the eigenvector exhibits at least six H transitions cor-related with the extended emission shown in the tomogram, indicat-ing that PC5 is either tracing the base of a disc wind, or the cavityof the large-scale outflow shown in Fig. 8, observed with a positionangle (PA) of 145 ◦ (from N to E), measured along the brightest MNRAS000
The tomogram of PC5 exhibits an extended emission with a differentmorphology to the main nebulosity of W33A (shown by the blackcontours), suggesting that this structure is either produced by adistinct mechanism or exhibits distinct physical conditions.Indeed, the eigenvector exhibits at least six H transitions cor-related with the extended emission shown in the tomogram, indicat-ing that PC5 is either tracing the base of a disc wind, or the cavityof the large-scale outflow shown in Fig. 8, observed with a positionangle (PA) of 145 ◦ (from N to E), measured along the brightest MNRAS000 , 1–19 (2021)
CA Tomography of the NIFS K-band observations of the high-mass protostar W33A Figure 12.
The tomogram (top left) and the corresponding eigenspectrum (top right) of the Principal Component 3 of the PCA Tomography of W33A. Fordetails, see Figs. 10 and 11. and more elongated region of the molecular emission. In addition,the spatial distribution of the H emission probed by the tomogramof PC5 is similar to the continuum-subtracted map of the H (1-0) S(1) transition at 2.12 𝜇 m presented by D10 (see their Fig. 15).The most significant difference between the tomogram and the linemap is that the first considers the emission of all the H transitionsavailable in the spectral range (producing a higher signal-to-noisemap) while the second method is based only on the intensity of asingle line (associated with a relatively larger noise).In addition, the Br 𝛾 feature and CO bandhead emission areanti-correlated with the H lines, indicating they arise from thecompact region at the launching point of the jet (indicated as red inthe tomogram). In general, the CO features are produced in a denseregion (e.g. the circumstellar disc) while the H is excited throughshocks in a more diffuse medium. Therefore, PC5 also shows thespatial separation of regions with distinct density regimes. The following results were obtained by analysing the data cubereconstructed using the first five Principal Components of the PCATomography (Sect. 2.4.2). We also compare our findings with thosereported by D10 to check whether the post-processing of the NIFSdata cubes led to a significant improvement on the determination ofthe physical parameters.Figure 16 presents a composite RGB image of the NIFS FOVsimilar to the one presented in Fig. 9. The spectra shown in the rightpanel were integrated over the regions indicated by the contours:the compact K -band object (red), the disc wind/cavity of the out-flow probed by the molecular H emission (green) and the mainnebulosity (blue).The integrated spectra shows that the main nebulosity isbrighter than the protostar at shorter wavelengths, while the em-bedded object gets brighter for 𝜆 & 𝜇 m. This behaviour was also observed in the analysis of PC2 (see Fig. 11), showing that thecontinuum emission is correlated with the nebulosity at short wave-lengths due to the reflection of the radiation from the embeddedprotostar, while it is correlated with the compact source for longerwavelengths likely due to local extinction effects.In addition, the spectra of the main nebulosity (in blue) and theH wind (in green) exhibit at least five ro-vibrational transitions ofthe H molecule. The H emission observed in the spectrum of themain nebulosity is likely due to the overlap between this region andthe H wind. The most energetic transition corresponds to the H 𝜇 m, with upper-level energy of 12,550 K. Theabsence of transitions with larger energy indicates that the physicalconditions of the H gas (i.e. temperature and column density)are not sufficient to produce more energetic transitions. We furtherinvestigate the physical parameters of the molecular H emission inSect. 3.4.2. We performed the photometry of the UKIDSS
𝐽𝐻𝐾 images usingthe IRAF/DAOPhot package to better characterise the protostellarobject candidates in the vicinity of W33A (Sources B and C). Weused the 1.4 ′ × ′ FOV shown in the top panel of Fig. 7, containingenough point-like sources to obtain a good photometric calibrationby cross-matching with the UKIDSS catalogue.First, we obtained an average map of the individual 𝐽 , 𝐻 , and 𝐾 images to identify the positions of the sources that appears in, atleast, one of the three filters. We obtained a list of 298 point-likeobjects with peak intensities above a 3- 𝜎 threshold. We performedthe photometry of the individual images, using an aperture withradii of 𝑟 = 1 ′′ , and a ring with inner and outer radius of 1 . ′′ . ′′ 𝐽, 𝐻 and 𝐾 ) with the UKIDSS catalogue, leading to a list of 132 objects. MNRAS , 1–19 (2021) Navarete et al.
Figure 13.
A detailed view of the disc-like structure probed by the PC3tomogram. Top panel: The normalised variance of the tomogram is shownin a divergent red-to-blue colour-scale. The position of the compact sourceis indicated by the black contours. The red arrow and blue dashed linesindicate the position angle of the rotation axis of the disc and its error,respectively, estimated as PA = 140 ± ◦ . The dashed black line shows thedirection from which the radial variance profile was extracted. Bottom: theradial variance profile extracted along the plane of the disc. The points anderror bars correspond to the mean variance and the 1- 𝜎 error of a 3-pixelwidth region perpendicular to the sampled spatial direction. The dashedvertical lines indicate the region where the absolute value of the variance islarger than zero ( | 𝑟 | ≤ . ′′ A final list containing the photometry of 287 objects, was obtainedafter removing the sources detected only in a single filter, excluding,for example, any bright K -band knot associated with the large-scaleoutflow in the FOV from the initial list.The 𝐽𝐻𝐾 magnitudes of the W33A, Source B and Source Care listed in Table 2, and the full list containing the 287 sourcesis available as on-line material (The first ten rows are listed inTable A1. W33A was not detected in the 𝐽 - and 𝐻 -band and, thus,the reported values correspond to upper limits, as well as the 𝐽 -bandmagnitude of Source C. In addition, the reported 𝐽 -band magnitudeof Source B corresponds to a foreground star in the line-of-sight ofthe object detected in the 𝐻 - and 𝐾 -bands.Fig. 17 presents the 𝐾 versus 𝐻 − 𝐾 colour-magnitude di-agram (top panel) and the 𝐽 − 𝐻 versus 𝐻 − 𝐾 diagram (bot-tom) of the sources within the 1.4 ′ × ′ FOV around W33A.The colour-magnitude diagram indicates that the three sources(W33A, Source B and Source C) are relatively bright K -band ob-jects ( 𝐾 <
13 mag), exhibiting 𝐻 − 𝐾 colour indices consistent withembedded young objects ( 𝐻 − 𝐾 > 𝐻 − 𝐾 = 5.54 ± 𝐻 − 𝐾 = 3.62 ± 𝐻 − 𝐾 = 3.46 ± K -band sources ofa sequence of reddened objects, confirming that they are deeply Table 2.
𝐽 𝐻 𝐾 photometry of the protostellar object candidates in thevicinity of W33A.Source J (mag) H (mag) K (mag)W33A 20.77 ∗ ∗ ± + ± ± ∗ ± ± Notes: ∗ indicates upper limit due to the non-detection of a point-like sourceat the given position in the corresponding filter. + a foreground source lyingvery close to the object exhibits 𝐽 = 15.45 ± embedded within their circumstellar environments as expected forprotostars.The bottom panel of Fig. 17 presents a colour-colour diagrambased on the 𝐽 − 𝐻 and 𝐻 − 𝐾 colour indices. Given that the 𝐽 -bandmagnitude of the three sources corresponds to upper limits of theirreal values, their position in the 𝑦 -axis of the diagram correspondto lower limits. Despite of that, the three sources are consistent withreddened objects, corroborating the young nature of these objects.In Fig. 18, we present the spectrum of W33A (shown as thered curve) and Source B (in blue), integrated over a 0 . ′′
15 radius,each with the extended nebular emission subtracted using an an-nular region around each compact sources, with inner and outerradii of 0 . ′′
20 and 0 . ′′
30, respectively (indicated by the white con-tours in Fig. 18). We note that W33A has a steeper continuum dueto the emission of the dust, while Source B has a blue continuum.For comparison, we also show the integrated spectrum of Source Bextracted using the original data cube (shown as the green curve),before the post-processing of the NIFS data cube and its recon-struction (see details in Sect. 2.3). The comparison between thespectrum of Source B before and after the reconstruction of the datacube shows the effectiveness of the noise suppression due to thepost-processing procedures applied to the data cube. In addition,the K -band spectrum of Source B exhibits relatively weak emissionof the Br 𝛾 and the H transition at 2.12 𝜇 m. These features aretypically observed in YSOs, corroborating the photometric resultsfrom Fig. 17. In Figure 19 we present the velocity and the dispersion map ofthe H (1–0) S(1) transition at 2.12 𝜇 m, the brightest H transi-tion identified in the K -band NIFS data cube. The radial velocitywas corrected by the rest velocity of W33A ( 𝑉 LSR = 36.7 km s − ,Lumsden et al. 2013). The top panel of Fig. 19 shows that the bulkof the H emission is located at a radial velocity of −
56 km s − ,showing no significant velocity gradient in the FOV. In addition, nosubstantial change in the full-width at zero intensity (FWZI) of theline profile is observed towards the emitting region, exhibiting val-ues around 110 km s − (middle panel). The bottom panel of Fig. 19shows the high-velocity of the H emission profile, evaluated as theterminal velocity at the blue side of the H profile (estimated fromthe FWZI parameter). The velocity map exhibits relatively highervelocity values close to the W33A protostar ( ∼ −
130 km s − ) and isroughly constant at large distances ( ∼ −
113 km s − ).We further measured the fluxes of the individual H transitionsin the integrated spectrum of the H cavity (shown as the greenregion in Fig. 16) to construct a ro-vibrational Boltzmann diagramof the H molecule and, thus, to evaluate the physical parametersof the jet. Table 3 summarises the fluxes and the physical propertiesof each H transition. Figure 20 presents the Boltzmann diagram MNRAS000
113 km s − ).We further measured the fluxes of the individual H transitionsin the integrated spectrum of the H cavity (shown as the greenregion in Fig. 16) to construct a ro-vibrational Boltzmann diagramof the H molecule and, thus, to evaluate the physical parametersof the jet. Table 3 summarises the fluxes and the physical propertiesof each H transition. Figure 20 presents the Boltzmann diagram MNRAS000 , 1–19 (2021)
CA Tomography of the NIFS K-band observations of the high-mass protostar W33A Figure 14.
The tomogram (top left) and the corresponding eigenspectrum (top right) of the Principal Component 4 of the PCA Tomography of W33A. Fordetails, Figs. 10 and 11.
Figure 15.
The tomogram (top left) and the corresponding eigenspectrum (top right) of the Principal Component 5 of the PCA Tomography of W33A. Fordetails, Figs. 10 and 11. of the H emission, constructed using the methodology adopted byCaratti o Garatti et al. (2006, 2015) and briefly explained below.The upper-level column density of each H transition ( 𝑁 𝜐,𝐽 ,where 𝜐 and 𝐽 denotes the upper vibrational and rotational levels,respectively) is evaluated as: 𝑁 𝜐,𝐽 = 𝜋 Ω ℎ 𝐹 𝜐,𝐽 𝐴 𝜐,𝐽 𝜈 𝜐,𝐽 𝑔 𝜐,𝐽 (5)where ℎ is the Planck constant (in erg s), Ω is the solid angle of the emitting region (in sr), 𝐹 𝜐,𝐽 is the integrated intensity of thetransition (in erg cm − s − sr − ), 𝐴 𝜐,𝐽 corresponds to the Einsteincoefficient of the transition (in s − ), 𝜈 is the frequency (in Hz), and 𝑔 𝜐,𝐽 is the statistical weight of the transition. The upper-level col-umn density of each H transition is plotted against their upper-levelenergy ( 𝐸 𝜐,𝐽 , in K). Assuming the gas is in Local ThermodynamicEquilibrium (i.e. follows the Boltzman distribution), a linear fit tothe data in a log-normal scale (i.e. log ( 𝑦 ) = 𝛼 𝑥 + 𝛽 ) provides anestimate of the excitation temperature ( 𝑇 exc = − / 𝛼 ), and the total MNRAS , 1–19 (2021) Navarete et al.
Figure 16.
Left panel: False-colour RGB image of the NIFS FOV (left, blue: 2.05 𝜇 m, green: 2.2 𝜇 m, red:2.35 𝜇 m) recreated using PCs 1–5. The contoursindicate the regions where the integrated spectra were extracted. Right: integrated K -band spectra of the W33A point-like source (in red), the main nebulosity(in blue), and the H wind (in green). Table 3. H transitions identified in the spectrum of the cavity of the outflow.Transition Flux 𝜎 𝐹 𝜆 𝑔 𝐸 upper 𝐴 (10 − erg cm − s − ) ( 𝜇 m) (K) (10 − s − )1-0 S(2) 8.51 0.36 2.0338 9 7584 3.981-0 S(1) 30.19 0.40 2.1218 21 6956 3.471-0 S(0) 9.09 0.36 2.2235 5 6471 2.532-1 S(1) 5.16 0.40 2.2477 21 12550 4.981-0 Q(1) 58.49 0.44 2.4066 9 6149 4.29 Notes: the energy levels are from Dabrowski (1984), and the Einsteincoefficients are from Turner et al. (1977). column density [ 𝑁 tot = 𝑄 ( 𝑇 ex ) exp ( 𝛽 ) ] of the H gas. The K -bandextinction ( 𝐴 K ) can also be inferred by choosing the best 𝐴 K valuethat maximises the correlation between the quantities shown in theBoltzmann diagram (for details, see Caratti o Garatti et al. 2015).The best-fit to the Boltzmann diagram shows that the moleculargas is associated with a 𝑇 exc value of (2.3 ± · K, and a totalH column density of log( 𝑁 tot /cm − ) = 18.8 ± K -bandreddening is about 4.20 ± 𝐴 𝑉 of63 ± 𝐴 V / 𝐴 K ratio of 14.95 from Damineli et al.(2016).The origin of the H emission can be inferred using the ratioof H transitions such as the (1–0) S(1)/(2–1) S(1). Based on thefluxes reported in Table 3 and correcting them for reddening, weobtained a (1–0) S(1)/(2–1) S(1) ratio of 9.3 ± ∼
10 for thermal excitationof the H emission through shocks, and significantly larger than thepredicted value of 1.8 for the case of radiative excitation throughUV photons from the protostar (Black & Dalgarno 1976).By taking advantage of the spatial information available in theNIFS data, we implemented a pixel-to-pixel evaluation of the phys-ical parameters of the H gas, as presented in Fig. 21. Spatial mapsof 𝐴 K , 𝑇 exc and 𝑁 tot were obtained by constructing and fitting theBoltzmann diagram of the H transitions for the spaxels exhibitingfluxes above a 3- 𝜎 threshold for all the five H transitions listed inTable 3. In addition, the mass of the H gas ( 𝑀 𝐻 ) was computedbased on the total column density of the H molecules as: 𝑀 H2 ( 𝑥, 𝑦 ) = 𝜇𝑚 𝐻 𝐴 pix 𝑁 tot ( 𝑥, 𝑦 ) (6)where the (x,y) indices correspond to the position of a given spaxelof the data cube, 𝜇 = 1.24 is the mean atomic weight of the H , 𝑚 H is the proton mass (in g) and 𝐴 pix corresponds to the projected areaof the spaxel (in cm ). The top left panel of Fig. 21 presents the distribution of the K -band extinction over the H emission associated with the cav-ity of the outflow. The median K -band extinction of the H emis-sion is 𝐴 K = 4.1 ± 𝐴 V of ∼ ± 𝑇 exc valueis (2.2 ± · K. The total column density of the H gas is pre-sented in the bottom left panel. The distribution of 𝑁 tot is similarto the flux of the H (1–0) S(1) line at 2.12 𝜇 m, indicated by thecontours, suggesting a linear relationship between the brightness ofthis transition and the logarithm of the column density. The totalcolumn density of the gas is about log( 𝑁 tot ) ∼ ± gas, com-puted using Eq. (6). The total mass of the emitting H region is about(5.4 ± · − M ⊙ .A rough estimate of the dynamical timescale of the cavity ofthe outflow, probed by the H emission, can be evaluated as 𝑡 dyn = ℓ H2 | 𝑣 H2 | = cos ( 𝑖 ) sin ( 𝑖 ) ℓ proj | 𝑣 rad | ≈ . ( 𝑖 ) ℓ proj ( AU )| 𝑣 rad | ( km s − ) yr (7)where ℓ proj is the projected length of the emitting region in theplane of the sky ( ℓ proj = ℓ H2 · sin ( 𝑖 ) ), 𝑣 rad is its radial velocity of theemission ( 𝑣 rad = 𝑣 H2 · cos ( 𝑖 ) ), and 𝑖 is the inclination angle of thestructure ( 𝑖 = 30 ◦ ). The H emission extends up to 1 . ′′
65 from thedriving source, corresponding to a projected linear scale of 3960 AUat the distance of 2.4 kpc. The emission peaks at 2.1216 𝜇 m, corre-sponding to a velocity of about − − . Considering the restvelocity of W33A, 𝑣 𝑙𝑠𝑟 = 36.7 km s − (Lumsden et al. 2013), theprojected radial velocity of the H emission is about −
65 km s − ,in agreement with the analysis presented in Fig. 19. Using Eq. (7)and propagating the error on each input parameters , we obtaineda dynamic timescale of about (5.0 ± · years, suggesting theH emission arises from a relatively young structure. We further ob-tained a rough estimate of the mass-loss rate of the ejected gas as ¤ 𝑀 H2 = 𝑀 H2 / 𝑡 dyn ∼ (1.1 ± · − M ⊙ yr − . We considered a typical ± ◦ on the inclination angle, ±
29 km s − on thevelocity (half of the velocity resolution of the NIFS data), and an error of10% on the projected length MNRAS000
29 km s − on thevelocity (half of the velocity resolution of the NIFS data), and an error of10% on the projected length MNRAS000 , 1–19 (2021) CA Tomography of the NIFS K-band observations of the high-mass protostar W33A Figure 17. 𝐾 versus 𝐻 − 𝐾 colour-magnitude diagram (top panel) and the 𝐽 − 𝐻 versus 𝐻 − 𝐾 colour-colour diagram (bottom) of the 1.4 ′ × ′ regionaround W33A. W33A is indicated as the red diamond symbol, Source B asthe blue triangle and Source C as the green squared symbol. The photometricerrors are shown as error bars and the arrows indicates lower limits of theactual values (see text). According to the PCA results, the Br 𝛾 emission is observed towardsthe W33A protostar (PC1) and the extended nebulosity to the south(PC1 and PC2), the rotating disc associated with W33A (PC3), andcompact emitting regions closer to the W33A protostar (PC4 andPC5).Davies et al. (2010) reported that the Br 𝛾 emission was prob-ing the cavity of the large-scale K -band outflow at very small scalesusing the spectro-astrometry of the Br 𝛾 line, a technique that com-bines the spectral and spatial information of the data cube to resolvekinematic structures at sub-pixel scales with impressive accuracy(i.e. one tenth of a milliarcsecond). Thus, we adopted the samemethodology presented by D10 to check if the post-processing of the NIFS data led to a significant improvement of the data at sub-pixel scales.The spectro-astrometry of the Br 𝛾 feature of W33A is pre-sented in Fig. 22. The left panel presents the offset of the centroid inthe right ascension (black curve) and declination axes (red curve) asa function of the wavelength (top panel), together with the Br 𝛾 fluxat the same wavelength scale (bottom). The right panel presents thespatial variation of the centroid position as a function of the velocity,matching the colour scale of the Br 𝛾 profile shown in the left panel.The positional error of 0.12 mas, indicated in the top left corner,was based on the centroid position offset of the nearby continuumaround the Br 𝛾 line (corresponding to the spectral regions shownin red).The blue and red high-velocity components (corresponding tovelocities of −
300 to +300 km s − ) are offset by roughly ∼ ≈ ± ◦ (from N to E),which is compatible with the results reported by D10(156 ± ◦ ). We used the
Penalized Pixel-Fitting algorithm (pPXFCappellari & Emsellem 2004) to derive the velocity map ofthe CO absorption features observed in the spectrum of W33A.The pPXF uses a template spectrum to fit and evaluate the radialvelocity at each spaxel in the data cube. The procedure starts withan initial guess of the velocity and width parameters and solves aGauss-Hermite function based on the template spectrum to fit theobserved one. The residuals between the model and the observedspectra are perturbed using a penalised function and the best modelis chosen based on a 𝜒 minimisation procedure, delivering a set ofparameters (V, 𝜎 , ℎ , ... ℎ 𝑚 ) and their uncertainties for each spaxelof the datacube. For this work, we focused the further analysison the interpretation of the radial velocity measurements only.The brightest spaxel of the continuum emission at ∼ 𝜇 m (atthe central position of the point-source object) was adopted as thetemplate spectrum (corresponding to the rest-velocity at 0 km s − )to derive the kinematics of the CO absorption lines using thespectral region ranging from 2.32 to 2.37 𝜇 m, containing both thelow-J CO features of the R and the P branches. To spatially samplethe rotating structure, the calculations were performed only for thespaxels within a 0 . ′′ ± ◦ (from N to E), compatible with the orientationof the rotating axis observed in the tomogram of PC3 (see Fig. 13).The velocity gradient arises from a relatively compact region of ∼ . ′′
5, corresponding to a projected length of about 1,200 AU at thedistance of 2.4 kpc.We used the velocity map delivered by the pPXF analysis to equivalent to the value reported by the authors, PA = 113 ± ◦ , orientedfrom E to NMNRAS , 1–19 (2021) Navarete et al.
Figure 18.
Left panel: False-colour image of the NIFS FOV integrated from 2.03 to 2.04 𝜇 m. The point-like source is saturated to enhance the contrast betweenW33A and the fainter, K -band counterpart of Source B at the SE edge of the datacube. The black contours indicate the regions where the integrated spectrawere extracted. The nebular emission, evaluated over the regions indicated by the white contours, was subtracted from the spectra. Right: integrated K -bandspectra of the W33A point-like source (in red), and the Source B before (green) and after the post-processing routines applied to the NIFS datacube (in blue).The spectra of Source B were multiplied by a factor of 400. extract a position-velocity (PV) diagram perpendicularly to the rota-tion axis of the disc. The results are presented in the middle panel ofFig. 23, exhibiting the mean value and the standard deviation of thevelocity at each position, evaluated within 3-pixel wide bins acrossthe extraction axis (indicated by the filled black line in the top panelof Fig. 23). The PV diagram indicates a clear separation between theblue- and red-shifted components (with peak velocities around + − − , respectively) separated by roughly 0 . ′′
2. The bestKeplerian model, shown as the solid blue curve, indicates that thetotal mass of the “protostar+disc” system is about 8.7 + . − . M ⊙ . Wenote that the velocity falls more rapidly than any Keplerian modelat larger radii. The mass and its standard deviation were evaluatedusing a Monte-Carlo method (see results in the bottom panel ofFig. 23), varying the input data randomly within their errors. Weused a total of 10,000 simulations to achieve the reported values.Assuming the inclination angle of 𝑖 = 60 ◦ adopted by D10 (estimatedfrom interferometric mid-IR observations from de Wit et al. 2010),the inclination-corrected mass of the “protostar+disc” system inW33A is 10.0 + . − . M ⊙ . This work presents a re-analysis of the 𝐾 -band IFU observationsof the well-known high-mass protostar W33A with an eye to fur-ther extend the analysis presented originally by Davies et al. (D10)through application of a Principal Components Analysis.First, we performed the flux-calibration of the data cube fol-lowed by the differential atmospheric refraction (DAR) correction.The DAR effects correspond to a slight positional offset of thesource as a function of the wavelength (see Fig. 1). The proper cor-rection of such effects is necessary to properly compare the spatialinformation at distinct wavelengths in the data cube, and to performmeasurements at a very small fraction of a pixel as presented in thiswork; see also D10.Then, we performed the PCA tomography of the data cube, andfound that most of the astrophysical information of the data cube iswell-explained by the first five Principal Components (PCs): PC1shows that the embedded compact source associated with W33Adominates the emission in the K -band; PC2 presents the differ-ences between the compact and the extended emission regions,mostly showing an anti-correlation introduced by reddening effects;PC3 exhibits the kinematics of a rotating circumstellar disc around W33A, probed by the low- 𝐽 CO absorption lines and exhibitinga contribution of ionised gas; PC4 probes the cavity of the jets,showing a roughly perpendicular structure aligned to the rotatingaxis of the disc, which is likely reflecting the emission of the innerregion of the disc and the protostar itself; and PC5 exhibits the anti-correlation between the extended H emission and the Br 𝛾 and COemission arising from the compact object.Along with the PCA results, the analysis of the near-infrared 𝐽𝐻𝐾 maps of W33A allowed us to obtain an overall interpretation ofthe environment associated with the protostar and its surroundings,as illustrated in Fig. 24.The data cube was reconstructed using the first five PCs, sup-pressing a considerable fraction of noise associated with PC6 andhigher-orders PCs (see Scree test in Fig. 5). Then, we followed asimilar analysis as that performed by D10 to validate our resultsusing the reconstructed data cube. In general, we found that thepost-processing routines and noise suppression of the data cube leda significant improvement of the results, as discussed below.
The comparison of the large-scale NIR maps with the NIFS FOVin Fig.s 7 and 8 shows that different sources are likely contributingto the observed extended K -band emission detected in the field.In Fig. 11, PC2 shows the anti-correlation between the cir-cumstellar region of W33A (probed mostly by the resolved COabsorption features) and the extended emission mostly probed bythe continuum and Br 𝛾 emission. The tomogram of PC2 exhibitsa clear separation between the extended emission closer to W33Aand the larger structure towards the S.A point-like source is observed at the SE edge of the NIFSFOV, consistent with the position of the Source B (see Fig. 7). Atthe distance of W33A, the linear separation of W33A and Source Bis about ∼ 𝐻 -band compactSource C (see Fig. 7) indicates that multiple sources are locatedwithin this region. The 𝐽𝐻𝐾 photometry and the colour-magnitudeand colour-colour diagrams presented in Fig. 17 allowed us to con-firm that Source B and C are compatible with young, protostellarobjects.
MNRAS000
MNRAS000 , 1–19 (2021)
CA Tomography of the NIFS K-band observations of the high-mass protostar W33A We obtained a robust estimate of the total mass of the“disc+protostar” system and its error based on the pPXF analy-sis of the CO absorption features at 2.35 𝜇 m (see Fig. 23). The bestfit of the position-velocity diagram is consistent with a Keplerianrotation structure around a total mass of Msin ( 𝑖 ) = 8.7 + . − . M ⊙ , as-suming a distance of 2.4 kpc. Assuming the inclination angle of60 ◦ as adopted by D10, the mass is 10.0 + . − . M ⊙ . If considering atypical ± ◦ uncertainty on the inclination angle, the upper-limitof the mass increases by a factor of two, 10.0 + . − . M ⊙ , indicatingthat the inclination uncertainty is a significant source of error on thetotal mass. The inclination angle value of 60 ◦ adopted in D10 and inthis work was derived from the analysis of interferometric mid-IRobservations from de Wit et al. (2010). Those authors argued thatthe derived inclination angle of the system is consistent with thefact that W33A exhibits no strong 𝐻 − band emission arising fromthe central star (as one would expect for a less inclined system), andit is also in agreement with the relatively weak large-scale K -bandredshifted outflow lobe (Fig. 7). Ilee et al. (2013) reported a moreaccurate determination of the inclination angle of the W33A systemas 𝑖 ( ◦ ) = 37 + − based on the modelling of the CO bandhead emis-sion at 2.3 𝜇 m using high spectral resolving power observations(R ∼ K -band emission of the high-mass YSO IRAS 11101 − the jet driven by the IRASsource requires that the disc must be roughly edge-on, with lowinclination angle. They conclude that most of the CO emission isobserved as reflected and scattered light in the outflow cavity walls,and therefore the inclination of the system can be incorrectly in-terpreted based on the modelling of the CO bandhead emission.In addition, the identification of the low- 𝐽 CO features in absorp-tion in the spectrum of W33A (e.g. see Fig. 9) is also indicative ofan inclined disc and, therefore, not compatible with such a smallinclination value as that reported by Ilee et al. (2013).Davies et al. (2010) derived the total mass of W33A as15 + − M ⊙ , considering a distance of 3.8 kpc. Scaling their mass tothe revised distance adopted in our work, 2.4 kpc, the total mass ofW33A is 9.5 + . − . M ⊙ , which is in good agreement with the masswe obtained with the pPXF analysis. The pPXF method allows thedetermination of the kinematics of very faint spectral features, suchas the CO absorption features, and also provides a robust estimateof the errors, an important extension of the D10 result.The rotation curve based on the CO absorption lines reportedby D10 was derived by smoothing the original CO absorption ve-locity centroids using a 5-pixel boxcar filtering. This proceduresuggests this component is a relatively extended structure with anangular size of about 1 ′′ around the source (see their Fig. 14), cor-responding to a radius of ∼ . ′′ ∼
500 AU) from the central position. The velocities deviatefrom Keplerian outside the inner 0 . ′′ ∼
500 AUat the distance of the source). These results suggest the circum-stellar disc associated with W33A is a compact structure, similarto the one observed towards the high-mass YSO IRAS 13481-6124by Kraus et al. (2010). Based on interferometric infrared observa-tions, those authors reported a ∼
130 AU flared disc around a 18 M ⊙ protostellar object lying within the IRAS source. Recently, Maud et al. (2017) presented the analysis of interfer-ometric ALMA observations of the W33A protostar at Band 6 and 7(230 and 345 GHz, respectively) with angular resolutions of ∼ . ′′ . ′′ K -band IFU observations. ALMA observations obtainedat larger baseline lengths (e.g. 16 km) are, therefore, necessary todisentangle the structure and kinematics of the innermost region ofW33A. The near-IR maps from Fig. 7 shows no evidence for a large-scale H jet associated with W33A, as one would expect for an active proto-star (e.g. Varricatt et al. 2010; Navarete et al. 2015). At longer wave-lengths, recent ALMA observations of W33A (Maud et al. 2017)have shown that high-velocity emission of CO and SiO lines areprobing a relatively compact jet-like structure with angular sizes of ∼ . ′′
5. Despite of that, no evidence for a collimated jet-like structurewas observed in the NIFS FOV, but a wide-angle H polar emis-sion driven by the protostar. The nature of the H emission can betraced by the other structures associated with W33A. The orien-tation of the sub-pixel structure probed by the spectro-astrometryof the Br 𝛾 feature (PA = 154 ± ◦ , see Fig. 22) is compatible withthe disc rotation axis (150 ± ◦ , see Fig. 23), suggesting that theionised Hydrogen is tracing the base of the jet/wind in high-massYSOs as also seen in IRAS 13481-6124 through near-IR interfer-ometry (Caratti o Garatti et al. 2016). In addition, the orientation ofthe H emission (PA ∼ ◦ , see Fig. 15) is also compatible withthe ionised gas and the disc axis, favouring the interpretation of awide wind traced by the H emission rather than the emission ofscattered light in the cavity of the large-scale outflow. This inter-pretation is also in agreement with the physical parameters of theH emission derived through the H (1–0) S(1)/(2–1) S(1) line ratioand the ro-vibrational analysis presented in Sect. 3.4.2, confirmingthe thermal origin of the emission.The physical parameters of the H emission were derived byconstructing the Boltzmann diagram of the different H transi-tions of the entire H emitting region (see Fig. 20), and also byevaluating the parameters for each spaxel in the data cube (seeFig. 21). Based on the pixel-to-pixel analysis, the H emission ex-hibits a mean 𝐴 K reddening of 4.3 ± × K, and a total H column density of10 . ± . cm − . The total mass of the H emission was estimatedto be M H = (4.6 ± · − M ⊙ . In general, there is a good agree-ment between the physical parameters evaluated from both methods(i.e. the integrated spectrum and the pixel-to-pixel analysis). Thespatially resolved analysis also shows the clumpy character of theH emission. The 𝐴 K map shown in the top left panel of Fig. 21 in-dicates the extinction towards the H jet is strong and variable, withvalues ranging from ∼ K -band. The 𝐴 K values are larger towards the positions where the H emissionis stronger indicating dense clumps of gas and dust perhaps alongthe cavity walls of the large-scale outflow (recall the emission isconsistent with shocked gas). But the conclusion is that the H isprobing a wind.We further compared the physical parameters of the H emis-sion of W33A using the typical values observed toward large-scaleH jets of high-mass protostars from Caratti o Garatti et al. (2015). MNRAS , 1–19 (2021) Navarete et al.
The total column density of the H wind driven by W33A is con-sistent with the range of column density values (10 to 10 cm − )reported by those authors. In addition, they also derived the meanvisual extinction ( 𝐴 𝑉 ) of the H jets ranging from 1 to 50 mag, witha median of 𝐴 𝑉 = 15 mag. The upper limit of the visual extinctionvalues reported by them is compatible with the median 𝐴 𝑉 valuefor the H wind associated with W33A, 𝐴 𝑉 = 63 ± 𝐴 K towards the W33A proto-star based on the H analysis due to the absence of moleculargas emission closer to the protostar. Alternatively, the 𝐴 K valueof the point-like source can be evaluated based on its K -band ex-cess ( 𝐸 𝐻 − 𝐾 ). Assuming the extinction law from Damineli et al.(2016), the 𝐴 K is obtained as 𝐴 𝐾 = . · 𝐸 𝐻 − 𝐾 . For W33A,the observed 𝐻 − 𝐾𝑠 colour index is about 4.0 mag. Consider-ing the intrinsic 𝐻 − 𝐾𝑠 colour index for a B2 V star ( ∼
10 M ⊙ ), ( 𝐻 − 𝐾𝑠 ) = − 𝐾 𝑠 -band of 𝐴 𝐾 𝑠 ≈ 𝐴 𝑉 ≈
80 mag). This value is in agreement with the upper-limit ofthe 𝐴 K values associated with the H emission, shown in Fig. 21. The Br 𝛾 emission arises from gas at typical electron tempera-tures of 8,000-10,000 K and it is observed in distinct structuresof W33A, suggesting that the hydrogen recombination plays a sig-nificant role in cooling the gas. The Principal Component 1 (Fig. 10)indicates the bulk of the Br 𝛾 emission arises at the central sourceof W33A. The source of ionised gas could be an unresolved ultra-compact H ii region around the protostar, in agreement with theabsence of 5 GHz radio continuum emission towards the protostar(Purcell et al. 2013).In Principal Component 2 (Fig. 11), the Br 𝛾 emission mostlyarises from the the extended emission observed in the FOV, eitheroriginated by Source B or W33A itself.In PC3 (Fig. 12), the Br 𝛾 feature exhibits a similar kinematicpattern as the CO absorption features. Our interpretation is thatthe CO features are probing a compact and rotating circumstellardisc around W33A (see Fig. 23), and the surface of the disc ispartially ionised, producing a disc wind, similar to that observed inHerbig Ae/Be stars (Tambovtseva et al. 2016) but perhaps reachinglarger radii owing to the more massive central object. Although theanalysis of the H emission does not indicate the presence of anyphoto-ionised region in the extended wind, we cannot rule out thepresence of photo-evaporating wind in ∼
100 AU scales closer to thesurface of the underlying cold disc, where the H emission is weakprobably because the molecules are being dissociated. As the windmoves away from the central region, the H molecules are excitedbut not dissociated, producing the observed extended structure asobserved in e.g. Fig. 21. This interpretation is complementary tothe spectro-astrometric result shown in Sect. 3.4.3, which showsthe ionised gas tracing the cavity of the outflow at milli-arcsecondscales, while PC3 shows evidences for the presence of ionised gas at ∼
100 AU scales, in the outer regions of the circumstellar disc. TheNIFS resolution is not sufficient to disentangle the likely ionisedsurface from the underlying cold disc probed by the CO features.Thus, as we pointed out above in Sect. 4.2, follow-up observationsat higher spatial resolutions (e.g. using the very long baseline ofthe ALMA array) focusing on the investigation of the circumstellarenvironment may help provide a better and more accurate view ofthe structure on small scales.The PC4 in Fig. 14 shows that strong Br 𝛾 emission is observedvery close to the W33A protostar, while the spectro-astrometry of the Br 𝛾 feature (Fig. 22) shows the innermost region of the H jetcavity at milli-arcsecond scales ( ∼ ◦ ± ◦ . Our results are in agreement with those reported byD10, suggesting that the spectro-astrometry of the Br 𝛾 feature islikely probing the base of the bipolar wind/jet traced by the H emission. We also observe that the structure identified in Fig. 22 isless affected by noise than that presented in D10, indicating that thepost-processing of the NIFS data cubes and noise suppression playsa significant role on the analysis of sub-pixel scale structures withinthe data. We presented a reanalysis of the 𝐾 -band NIFS/Gemini North ob-servations of the protostar W33A, first done by D10, based on Prin-cipal Component Analysis Tomography. The PCS technique is welladapted to the analysis of complex data sets like the one analysedhere.(i) The PCA tomography was able to recover the spectral andspatial information of structures associated with the three key-ingredients of the disc-mediated accretion scenario of the high-massstar formation process: the cavity of the large-scale outflow probedby the ionised hydrogen, a rotating disc-like structure probed by theCO absorption features, and a wide-angle wind traced by the H emission. All these structures were present in the first five PrincipalComponents of the PCA Tomography. In addition, the techniquealso reveals structures with different kinematics (e.g. the rotationof the disc and the extended wind structure) in a much more effi-cient way than investigating each spectral feature using traditionalmethods.(ii) The ionised gas (Br 𝛾 ) contributes to every structure associ-ated with the W33A protostar: although the bulk of the ionised gasis observed at the central source as expected for a high-luminosityprotostar, it also tracks the inner region of the jet and its cavity, aswell as the surface of the disc.The comparison of the structures identified in the NIFS FOVwith large-scale near-infrared maps around W33A led us to thefollowing conclusions:(iii) The propagation of the H jet identified in the NIFS FOVis consistent with the orientation of the 𝐾 -band reflection nebulaidentified in the 𝐽𝐻𝐾 maps. However, no evidence for a large-scaleH jet driven by W33A was found.(iv) A point-like object offset by ∼ ′′ to the SE of W33A wasidentified in the NIFS FOV. This object (Source B) has a relativelyblue spectrum, which is consistent with the bright 𝐽 -Band sourceidentified in the 𝐽𝐻𝐾 maps. Source B is also associated with a 𝐻 -band extended emission oriented to the NW direction that overlapsthe emission from W33A, contributing to the observed 𝐾 -bandcontinuum extended emission observed in the NIFS FOV.(v) A third point-like source (Source C), offset by 6 ′′ to the Sdirection of W33A, was identified in the large-scale near-infraredmaps. Source C is likely the driving source of a set of H knotsoriented in the E-W direction. The identification of sources B and Crevealed that W33A is located in a region with multiple (proto)stellarobjects.We reconstructed the data cube using only the first five Prin-cipal Components which represents 99.996% of the variance andfrom this data cube we derive the following conclusions: MNRAS , 1–19 (2021)
CA Tomography of the NIFS K-band observations of the high-mass protostar W33A (vi) First, we derived the physical parameters of the H jetand found that it is associated with an excitation tempera-ture of (2.2 ± · K, and a mean H column density oflog( 𝑁 tot /cm − ) = 18.7 ± 𝐾 -band extinction is about4.1 ± 𝐴 V of ∼
61 mag), and the massof the H jet is about log(M 𝐻 /M ⊙ ) = − ± 𝛾 emission,confirming that the ionised gas is tracing the cavity of the jets at sub-pixel scales. The positional offset between the blue- and red-shiftedcomponents is ∼ . ′′
55, about three times smaller than previously re-ported ( ∼ . ′′ ∼ . ′′ + . − . M ⊙ assuming a distance of2.4 kpc. The mass is about ∼
33% smaller than previously reportedby Davies et al. (2010), placing W33A closer to the lower mass-limit for a high-mass protostar.(ix) The above-mentioned results were obtained thanks to theapplication of the PCA Tomography towards a bona fide high-massprotostar, showing the technique is encouraging to be applied toother high-mass young stellar object candidates.
ACKNOWLEDGEMENTS
FN thanks to Fundação de Amparo à Pesquisa do Estado de SãoPaulo - FAPESP for support through proc. 2017/18191-8. AD ac-knowledges FAPESP for support through proc. 2019/02029-2. Theauthors would like to thank the valuable comments and sugges-tions made by the anonymous referee that helped to improve themanuscript.
DATA AVAILABILITY
The data sets were derived from observations available in the pub-lic domain: (under the Gemini program GN-2008A-Q-44).
REFERENCES
Alves F. O., Caselli P., Girart J. M., Segura-Cox D., Franco G. A. P.,Schmiedeke A., Zhao B., 2019, Science, 366, 90Black J. H., Dalgarno A., 1976, ApJ, 203, 132Blum R. D., McGregor P. J., 2008, AJ, 135, 1708Blum R. D., Barbosa C. L., Damineli A., Conti P. S., Ridgway S., 2004,ApJ, 617, 1167Cappellari M., Emsellem E., 2004, PASP, 116, 138Caratti o Garatti A., Giannini T., Nisini B., Lorenzetti D., 2006, A&A,449, 1077Caratti o Garatti A., Stecklum B., Linz H., Garcia Lopez R., Sanna A., 2015,A&A, 573, A82Caratti o Garatti A., et al., 2016, A&A, 589, L4Cohen M., Wheaton W. A., Megeath S. T., 2003, AJ, 126, 1090Dabrowski I., 1984, Canadian Journal of Physics, 62, 1639Dahmer-Hahn L. G., et al., 2019, MNRAS, 489, 5653Damineli A., Almeida L. A., Blum R. D., Damineli D. S. C., Navarete F.,Rubinho M. S., Teodoro M., 2016, MNRAS, 463, 2653Davies B., Lumsden S. L., Hoare M. G., Oudmaijer R. D., de Wit W.-J.,2010, MNRAS, 402, 1504 Fedriani R., et al., 2020, A&A, 633, A128Froebrich D., et al., 2011, MNRAS, 413, 480Ilee J. D., et al., 2013, MNRAS, 429, 2960Immer K., Reid M. J., Menten K. M., Brunthaler A., Dame T. M., 2013,A&A, 553, A117Kewley L. J., Nicholls D. C., Sutherland R. S., 2019, ARA&A, 57, 511Kraus S., et al., 2010, Nature, 466, 339Lee H.-T., et al., 2013, ApJS, 208, 23Lumsden S. L., Hoare M. G., Urquhart J. S., Oudmaijer R. D., Davies B.,Mottram J. C., Cooper H. D. B., Moore T. J. T., 2013, ApJS, 208, 11Maud L. T., Hoare M. G., Galván-Madrid R., Zhang Q., de Wit W. J., KetoE., Johnston K. G., Pineda J. E., 2017, MNRAS, 467, L120McGregor P. J., et al., 2003, in Iye M., Moorwood A. F. M., eds, Societyof Photo-Optical Instrumentation Engineers (SPIE) Conference SeriesVol. 4841, Society of Photo-Optical Instrumentation Engineers (SPIE)Conference Series. pp 1581–1591, doi:10.1117/12.459448Menezes R. B., Steiner J. E., Ricci T. V., 2013, ApJ, 765, L40Menezes R. B., Steiner J. E., Ricci T. V., 2014, MNRAS, 438, 2597Murakawa K., Lumsden S. L., Oudmaijer R. D., Davies B., WheelwrightH. E., Hoare M. G., Ilee J. D., 2013, MNRAS, 436, 511Navarete F., Damineli A., Barbosa C. L., Blum R. D., 2015, MNRAS,450, 4364Purcell C. R., et al., 2013, ApJS, 205, 1Ricci T. V., Steiner J. E., Menezes R. B., 2011, ApJ, 734, L10Ricci T. V., Steiner J. E., Menezes R. B., 2014, MNRAS, 440, 2419Ricci T. V., Steiner J. E., Giansante L., 2015, A&A, 576, A58Shu F. H., Adams F. C., Lizano S., 1987, ARA&A, 25, 23Skrutskie M. F., Cutri R. M., Stiening R., Weinberg M. D., et al., 2006, AJ,131, 1163Steiner J. E., Menezes R. B., Ricci T. V., Oliveira A. S., 2009, MNRAS,395, 64Tambovtseva L. V., Grinin V. P., Weigelt G., 2016, A&A, 590, A97Tokunaga A. T., 2000, Infrared Astronomy. Springer, p. 143Tokunaga A. T., Simons D. A., Vacca W. D., 2002, PASP, 114, 180Turner J., Kirby-Docken K., Dalgarno A., 1977, ApJS, 35, 281Varricatt W. P., Davis C. J., Ramsay S., Todd S. P., 2010, MNRAS, 404, 661Wegner W., 2003, Astronomische Nachrichten, 324, 219de Wit W. J., Hoare M. G., Oudmaijer R. D., Lumsden S. L., 2010, A&A,515, A45+
APPENDIX A: APERTURE PHOTOMETRY OF THEUKIDSS JHK IMAGES
The
𝐽𝐻𝐾 magnitudes of the W33A, Source B and Source C arelisted in Table 2. The full list containing the 287 sources is availableonline as supplementary material and in the CDS.
This paper has been typeset from a TEX/L A TEX file prepared by the author.MNRAS , 1–19 (2021) Navarete et al.
Table A1.
𝐽 𝐻 𝐾
Photometry of the sources (first 10 rows). The full table is available online as supplementary material and in the CDS.ID RA Dec UKIDSS J 𝜎 𝐽 H 𝜎 𝐻 K 𝜎 𝐾 𝐽 − 𝐻 𝜎 𝐽 − 𝐻 𝐽 − 𝐾 𝜎 𝐽 − 𝐾 𝐻 − 𝐾 𝜎 𝐻 − 𝐾 (J2000) (J2000) Designation (mag) (mag) (mag) (mag) (mag) (mag) (mag) (mag) (mag) (mag) (mag) (mag)1 𝑎 ∗ – 15.05 ∗ – 9.51 0.19 5.72 0.64 11.26 0.66 5.54 0.212 𝑏 + 𝑐 ∗ – 16.31 0.11 12.69 0.23 4.46 0.64 8.08 0.68 3.62 0.264 273.66562 -17.86503 J181439.75-175154.1 11.37 0.07 10.99 0.08 10.29 0.20 0.38 0.11 1.08 0.21 0.70 0.215 273.66211 -17.85791 J181438.90-175128.3 13.24 0.08 12.63 0.09 12.35 0.21 0.61 0.12 0.89 0.22 0.28 0.236 273.67337 -17.87352 J181441.61-175224.6 13.35 0.08 13.07 0.09 12.93 0.21 0.28 0.12 0.42 0.23 0.14 0.237 273.67474 -17.86697 J181441.94-175201.1 13.57 0.08 13.18 0.09 12.99 0.21 0.39 0.12 0.59 0.23 0.20 0.238 273.66714 -17.88072 J181440.11-175250.5 14.11 0.08 13.14 0.09 12.63 0.21 0.96 0.12 1.48 0.22 0.51 0.239 273.66238 -17.87348 J181438.97-175224.5 14.10 0.08 13.30 0.09 12.88 0.21 0.81 0.12 1.22 0.23 0.42 0.2310 273.66678 -17.86504 J181440.03-175154.1 14.29 0.08 13.40 0.09 12.95 0.21 0.89 0.12 1.35 0.23 0.46 0.23 Notes: the columns are as follows: (1) ID of the source ( 𝑎 : W33A; 𝑏 : Source B; 𝑐 : Source C); (2) Right ascension (J2000, in degrees); (3) Declination (J2000,in degrees); (4) Name of the source or its designation in the UKIDSS catalogue; (5) 𝐽 -band magnitude; (6) 𝐽 -band magnitude error; (7) 𝐻 -band magnitude;(8) 𝐻 -band magnitude error; (9) 𝐾 -band magnitude; (10) 𝐾 -band magnitude error; (11) 𝐽 − 𝐻 colour index (in mag); (12) 𝐽 − 𝐻 colour index error (inmag); (13) 𝐽 − 𝐾 colour index (in mag); (14) 𝐽 − 𝐾 colour index error (in mag); (15) 𝐻 − 𝐾 colour index (in mag); (16) 𝐻 − 𝐾 colour index error (in mag).An ∗ symbol indicates upper limit due to the non-detection of a point-like source at the given position in the corresponding filter. + indicates the 𝐽 -band fluxof a foreground source lying very close to the targeted object. MNRAS , 1–19 (2021) CA Tomography of the NIFS K-band observations of the high-mass protostar W33A Figure 19.
Velocity map of the H (1–0) S(1) transition at 2.12 𝜇 m. Thetop panel presents the centroid velocity of the emission,the middle panelshows the full-width at zero intensity of the line profile, and the bottompanel indicates the terminal velocity at the blue side of the H emission.The black contours corresponds to 10%, 25%, 50% and 75% of the peakintensity of the H emission. Figure 20.
Boltzmann diagram of the H gas associated with W33A. Thedashed line indicates the best fit to the data. The corresponding temperature,total column density and reddening values are indicated in the upper-rightcorner.MNRAS , 1–19 (2021) Navarete et al.
Figure 21.
Pixel-to-pixel derivation of the H parameters for those pixels exhibiting fluxes above a 3- 𝜎 threshold. From top left to bottom right: the 𝐴 K , 𝑇 exc ,column density and mass of the H gas. The black contours indicates the integrated H (1–0) S(1) flux. The red contours indicate the position of the point-likesource and its extended continuum emission. The range of each parameter was chosen in order to improve the visualisation of the maps. The mean value of the 𝐴 K and 𝑇 exc , and the total column density and mass of the H gas are shown at the bottom of their respective panels. MNRAS , 1–19 (2021) CA Tomography of the NIFS K-band observations of the high-mass protostar W33A Figure 22.
Spectro-astrometry of the Br 𝛾 feature. Top left panel: the centroid position in the right ascension (RA, in black) and declination axis (Decl., inred) as a function of the wavelength. The baseline is indicated by a horizontal grey line. The peak of the Br 𝛾 feature (at 0 km s − ) is indicated by a verticalline. Bottom left: the normalised Br 𝛾 flux. The coloured region indicates the spectral interval shown on the right panel. The spectral regions used to estimatethe positional error of the centroid (about ∼ Br 𝛾 centroid position. Thevelocity with respect to the peak of the emission is indicated by the horizontal colour bar shown in the bottom region of the plot. The positional error is shownin the left top corner of the plot. The dashed line indicates the position-angle (PA) of the red and blue high-velocity components.MNRAS , 1–19 (2021) Navarete et al.
Figure 23.
Top panel: Radial velocity map of the CO absorption featurestowards W33A, derived with the pPXF fitting. The black contours are placedat the 3 𝑛 − 𝜎 level ( 𝑛 = , , , ... ) of the integrated continuum emission at ∼ 𝜇 m. The velocity scale is presented on the right side of the map. Thedashed black line indicates the position angle of 150 ◦ (from N to E) of therotation axis of the disc, following the velocity gradient between the blue-and the red-shifted disc components. The filled black line is perpendicularto the disc axis, indicating the direction where the position-velocity diagramwas extracted. Middle: The position-velocity diagram for the CO absorptionfeatures. The black circles and their errors corresponds the measurementsusing the pPXF analysis. The best Keplerian rotation curve for 𝑀 sin ( 𝑖 ) = . ⊙ is indicated by the solid blue line. The dotted and dashed blue linesindicate the Keplerian curves considering the upper and lower 1- 𝜎 errors onthe mass, respectively. Bottom panel: The cumulative distribution function(CDF) of the 𝑀 sin ( 𝑖 ) value evaluated for 10,000 Monte-Carlo simulationsis shown as the solid blue curve. The relative histogram of 𝑀 sin ( 𝑖 ) valuesis shown as the solid black line. The solid vertical line corresponds to themedian value and the dashed lines are placed at the 16 and 84% confidencevalues, corresponding to the 1- 𝜎 error bars. Figure 24.
Interpretation of the near-infrared
𝐽 𝐻 𝐾 maps and NIFS K -bandobservations of W33A. The illustration roughly follows the same orientationof the observations (North is up and East is to the left). The structures arelabelled and shown in different colours, as indicated in the top left corner.MNRAS000