A compact Keplerian-like disc and H30$\alpha$ emission towards the radio jet in the massive protostar W75N(B)-VLA3

TL;DR

This study uses high-resolution ALMA observations to identify a Keplerian-like disc around a massive protostar and reports the potential first detection of H30α emission from a protostellar radio jet, advancing understanding of massive star formation.

Contribution

It presents the first detailed observation of a compact, Keplerian-like disc around W75N(B)-VLA3 and suggests H30α emission originates from its radio jet, a novel detection in massive protostars.

Findings

Disc radius of ~450 au with Keplerian rotation around a 16 M_sun protostar.

Detection of broad H30α recombination line emission towards the jet.

VLA3 adds to rare systems with confirmed small-scale massive discs and jets.

Abstract

Studying young, high-mass stellar objects is challenging for testing models of massive star formation due to their great distances, often kiloparsecs away. This requires extremely high-angular resolution to resolve features like accretion discs around massive protostars. Consequently, while powerful, collimated outflows are evident in massive protostars, the compact accretion discs anticipated at their centres are still proving difficult to pinpoint. This study presents ALMA continuum and molecular line observations at 1.3 mm of the massive protostar W75N(B)-VLA3. The observations achieve an angular resolution of $\sim$0.12$^{\prime\prime}$ ($\sim$156 au). Dense gas tracers reveal a circumstellar disc of $\sim$450 au in radius surrounding VLA3, with an orientation perpendicular to its associated thermal radio jet. From the millimetre continuum, a total mass of $\approx$0.43-1.74 $M_\odot$ is estimated for the disc. The disc's velocity profile is consistent with Keplerian rotation around a protostar of $\approx$16 $M_\odot$. This adds VLA3 to the small number of massive disc-protostar-jet systems documented in the literature with a centrifugally supported disc with a radius less than 500 au. Additionally, we detected H30$\alpha$ recombination line emission towards the radio jet powered by VLA3. Despite limitations in the spatial and spectral resolution, our data reveal a very broad line, indicative of high-velocity motions, as well as a tentative velocity gradient in the jet's direction, thus favouring the H30$\alpha$ emission to originate from the radio jet. Should this interpretation be confirmed with new observations, W75N(B)-VLA3 could represent the first protostellar radio jet for which a thermal radio recombination line has been detected.

Tables

Table 1. Table 1: Spectral setup of the ALMA observations

Spectral window —	Central Freq.	Bandwidth	Channel width
main targeted species	(MHz)	(MHz)	(km s-1)
spw 0 — continuum	231805.970	1875.00	20.000
spw 1 — 12CO	230538.000	0468.75	00.734
spw 2 — 13CO	220398.684	0234.38	00.384
spw 3 — C18O	219560.358	0234.38	00.385
spw 4 — H2CO	218222.192	0468.75	00.775
spw 5 — SiO	217104.980	0468.75	00.779

Table 2. Table 2: Selection of bright spectral lines detected towards VLA3 (see blue markers in Fig.  2 ), and studied in this work

Species	transition	Freq. (MHz)	${\log }_{10}$ [$A_{\mathrm{ij}}$ (s-1)]	$E_{\mathrm{up}}$ (K)	$n_{\mathrm{crit}}$ (cm-3)	tracer of …
SiO	$5$–$4$	217104.9190	$-3.28429$	0031.26	$>$$2.5\times {10}^6$	shocked gas
H2CO	$3_{0,3}$–$2_{0,2}$	218222.1920	$-3.55007$	0020.96	$>$$3.0\times {10}^6$	dense gas
HC3N	$24$–$23$	218324.7230	$-3.08296$	0130.98	$>$$1.7\times {10}^7$	dense gas
SO2	${22}_{2,20}$–${22}_{1,21}$ ; v2=1	219465.5447	$-3.99721$	1012.74	…	hot, dense gas
HNCO	${10}_{3,8}$–$9_{3,7}$	219656.7695	$-3.92039$	0432.96	…	dense gas
	${10}_{3,7}$–$9_{3,6}$	219656.7708	$-3.92039$	0432.96	…	dense gas
CH3CN	${12}_8$–${11}_8$	220475.8078	$-3.29104$	0525.57	$>$$1.9\times {10}^7$	dense gas
H	$31$–$30$ (30$\alpha$)	231900.9280	…	…	…	ionised gas
H2O	$5_{5,0}$–$6_{4,3}$ ; v2=1	232686.7000	$-5.32146$	3461.91	…	hot, dense gas

Table 3. Table 3: Fitted parameters of the Keplerian-rotation models shown in Fig.  10

	$M_{\ast ,\mathrm{fit}}$	$v_0$	$r_0$	$r_0$
Species	($M_{\odot }$)	(km s-1)	(au)	(arcsec)
H2CO	$20.3\pm 1.2$	0$8.5\pm 0.3$	0$-1\pm 14$	0$-0.001\pm 0.010$
HC3N	$16.9\pm 1.0$	$10.0\pm 0.3$	$-41\pm 11$	0$-0.032\pm 0.008$
SO2	$15.9\pm 1.6$	0$8.9\pm 0.5$	$-29\pm 22$	0$-0.022\pm 0.017$
HNCO	$13.6\pm 1.0$	0$8.4\pm 0.4$	$-12\pm 19$	0$-0.009\pm 0.015$
CH3CN	$10.0\pm 0.9$	0$8.7\pm 0.3$	$-44\pm 15$	0$-0.034\pm 0.012$

Equations

$$M_{\mathrm{dust}}=\frac{S_{\nu }{d}^2}{{\kappa }_{\nu }{B}_{\nu }(T_{\mathrm{dust}})},$$

$$M_{\mathrm{dust}+\mathrm{gas}}=1.4\times 10^{-2}{M}_{\odot }\left[\frac{S_{\nu }}{\mathrm{mJy}}\right]\times {\left[\frac{T_{\mathrm{dust}}}{50\mathrm{K}}\right]}^{-1}.$$

$$v=\mathrm{sign}(r-r_0)\cdot \sqrt{\frac{G{M}_{\ast ,\mathrm{fit}}}{|r-r_0|}}+v_0$$

$$Q=\frac{c_{\mathrm{s}}\kappa }{\pi G\Sigma }$$

$$\frac{S_{\mathrm{L}}}{S_{\mathrm{C}}}={\left[0.25{\left(\frac{{\nu }_{\mathrm{L}}}{\mathrm{GHz}}\right)}^{1.1}{\left(\frac{T_{\mathrm{e}}}{10^4\mathrm{K}}\right)}^{-1.1}{\left(\frac{\Delta v}{{\mathrm{km}\mathrm{s}}^{-1}}\right)}^{-1}+1\right]}^{2/3}-1$$

$${\mathrm{V}}_{\min }={\mathrm{V}}_{\mathrm{j}}\cos ({\theta }_{\mathrm{o}}+{\theta }_{\mathrm{j}}),\mathrm{and}{\mathrm{V}}_{\max }={\mathrm{V}}_{\mathrm{j}}\cos ({\theta }_{\mathrm{o}}-{\theta }_{\mathrm{j}}),$$

$$\Delta \mathrm{V}={\mathrm{V}}_{\max }-{\mathrm{V}}_{\min }=2{\mathrm{V}}_{\mathrm{j}}\sin ({\theta }_{\mathrm{o}})\sin ({\theta }_{\mathrm{j}}),$$

$${\mathrm{V}}_{\mathrm{c}}={\mathrm{V}}_{\mathrm{j}}\cos ({\theta }_{\mathrm{o}})\cos ({\theta }_{\mathrm{j}}).$$

$$\begin{array}{rl}\Delta \mathrm{V}=2{\mathrm{V}}_{\max }& =2{\mathrm{V}}_{\mathrm{j}}\cos ({\theta }_{\mathrm{o}}-{\theta }_{\mathrm{j}})\\ & \simeq 2{\mathrm{V}}_{\mathrm{j}}\sin ({\theta }_{\mathrm{j}}),\mathrm{for}{\theta }_{\mathrm{o}}\simeq \pi /2.\end{array}$$

Full text

††pubyear: 2025††pagerange: A compact Keplerian-like disc and H30$\alpha$ emission towards the radio jet in the massive protostar W75N(B)-VLA3–A

1 Introduction

The high-mass star-forming region (HMSFR) W75N(B) is a well-studied region that contains three massive protostars detected at radio continuum and millimetre wavelengths: VLA1, VLA2, and VLA3 (e.g., Hunter et al., 1994; Torrelles et al., 1997; Gómez et al., 2023). These protostars, with luminosities equivalent to B1-B0.5 type stars (e.g., Torrelles et al., 1997; Shepherd, 2001; Anglada et al., 2018), are located within a region of $\sim$1.5 arcsec ($\sim$1950 au at the distance of 1.3 kpc, as estimated by Rygl et al., 2012), where two filaments seen in NH3 seem to be interacting (Carrasco-González et al., 2010). Previous studies have reported outflow activity observed in radio continuum, maser emission (predominantly associated with VLA1 and VLA2), and thermal molecular lines (e.g., Torrelles et al., 1997; Torrelles et al., 2003; Kim et al., 2013; Surcis et al., 2014, 2023; Carrasco-González et al., 2015; Colom et al., 2018, 2021; Rodríguez-Kamenetzky et al., 2020; Gómez et al., 2023).

These studies have revealed that the three massive protostars, at distinct evolutionary stages, exhibit a variety of characteristics. Specifically, VLA1 drives a thermal radio jet at scales of $\sim$0.1 arcsec ($\sim$130 au) and shows signs of early stage photoioinisation. The outflow associated with VLA2, observed in both water (H2O) maser emission and radio continuum, has undergone a remarkable transformation over several decades. Initially, it was almost isotropic, but has since become collimated at scales of $\sim$0.2 arcsec ($\sim$260 au; see Surcis et al., 2023, for a review of the properties of VLA1 and VLA2). Completing the picture, VLA3 drives a compact thermal radio jet at scales of a few hundred au (see Carrasco-González et al., 2010). This protostar is thought to be responsible for exciting the two pairs of obscured Herbig-Haro (HH) objects (Bc[E]-Bc[W], Bd-VLA4) detected at radio wavelengths $\sim$4 arcsec ($\sim$5200 au) south of VLA3. These HH objects exhibit proper motions moving outwards from the protostar at velocities exceeding 100 km s*-1* (see Rodríguez-Kamenetzky et al., 2020; Gómez et al., 2023, for a detailed description of all these properties). Relatively weak H2O maser emission has also been reported towards VLA3 (60 mJy, VLSR = 10.5 km s*-1*; Torrelles et al., 1997). Since VLA1, VLA2, and VLA3 share a common origin from the gas and dust environment of W75N(B), this region has become a prime candidate for testing high-mass star formation models.

While different theoretical scenarios have been proposed to explain the formation of high-mass stars (e.g., competitive accretion in a turbulence-dominated core, Krumholz et al. 2009; competitive accretion driven by a stellar cluster, Bonnell & Bate 2006; Bondi-Hoyle accretion, Keto 2007), all of them predict the formation of circumstellar discs that enable the transport of mass from the envelope/environment to the forming star (see, e.g., Tan et al., 2014; Krumholz, 2015; Motte et al., 2018; Zhao et al., 2020; Avison et al., 2023; Beuther et al., 2025, for different reviews on high-mass star formation). Therefore, at small scales, one expects to find massive protostars surrounded by discs, accompanied by collimated jets, similar to those observed in low-mass star formation. In their comprehensive review, Beltrán & de Wit (2016) summarise the evidence and properties for discs around intermediate-mass ($\sim$2–7 $M_{\odot}$) and high-mass ($\gtrsim$7 $M_{\odot}$) stars (see also Zhao et al., 2020; Ahmadi et al., 2023). Based on this study, while there is little doubt that discs exist around intermediate-mass stars, discs around the most massive stars still remain elusive. Indeed, even though powerful collimated outflows are identified toward several massive protostars, the expected compact accretion discs at their centres are still hard to find (e.g., Goddi et al., 2020). Observations tend to identify massive, large, rotating structures, usually referred to as ‘toroids’ with sizes $\gtrsim$5000 au (e.g., Torrelles et al., 1983; Beltrán et al., 2005, 2011). A few notable examples of massive protostars associated with compact disc–outflow systems, recently discovered thanks to the improved resolving power of current interferometers, include AFGL 4176 (Johnston et al., 2015, 2020), Cepheus A-HW2 (Patel et al., 2005; Curiel et al., 2006; Jiménez-Serra et al., 2007; Sanna et al., 2025), G11.92$-$0.61 MM1 (Ilee et al., 2018; Bayandina et al., 2025), G17.64+0.16 Maud et al. (2018, 2019), G23.01$-$0.41 (Sanna et al., 2019), G35.20$-$0.74 N (Sánchez-Monge et al., 2013a, 2014; Beltrán et al., 2016), GGD 27-MM1 (Carrasco-González et al., 2012; Girart et al., 2018; Fernández-López et al., 2023), IRAS 16547$-$4247 (Rodríguez et al., 2005; Zapata et al., 2019), IRAS 20126$+$4104 (Cesaroni et al., 1999; Cesaroni et al., 2025; Palau et al., 2017), Orion Src I (Hirota et al., 2017; Ginsburg et al., 2018), and W75N(B)-VLA2 (Torrelles et al., 1997; Carrasco-González et al., 2015; Gómez et al., 2023). For other cases, see Ginsburg et al. (2023) and references therein.

In the recent paper by Gómez et al. (2023), we presented ALMA observations of the continuum, SiO and H2CO emission lines at 1.3 mm wavelengths towards the HMSFR W75N(B). These observations detected the VLA1, VLA2, and VLA3 protostars at millimetre continuum wavelengths at a superb angular resolution and sensitivity. The H2CO emission exhibits a fragmented structure surrounding the three massive protostars, while the SiO emission is strongly concentrated towards VLA2. The elongated structure of the SiO emission, perpendicular to the VLA2 thermal radio jet, together with its kinematics, led the authors to interpret it as arising in strong shocks within a toroid and a wide-angle outflow encircling this protostar.

In the current work, we now shift our focus to VLA3, which is the brightest high-mass protostar at millimetre wavelengths in W75N(B). We undertake a detailed analysis of the spatio-kinematical distribution of different molecular dense gas tracers, together with shocks and ionised gas tracers, and the mm-continuum around VLA3. In Sect. 2, we provide a brief summary of the ALMA observations carried out by Gómez et al. (2023). The main results obtained towards VLA3 are presented in Sect. 3, while Sect. 4 discusses these findings. Finally, Sect. 5 summarises our main conclusions.

2 Observations

We used ALMA to observe the continuum emission at 1.3 mm and multiple spectral lines towards W75N(B)-VLA3. The observations were carried out on April 29, August 18, and October 28, 2021, and on August 2, 2022, as part of project 2019.1.00059.S. A full description of the setup of the ALMA observations as well as the calibration and reduction procedures is given in Gómez et al. (2023). The phase centre of the observations is $\alpha$(J2000) = 20h38m36$.\!^{\mathrm{s}}$486, $\delta$(J2000) = 42*∘*37′34$.\!\!^{\prime\prime}$09. The spectral setup (see Table 1) consists of a broadband (1.875 GHz wide) spectral window (spw 0) with low spectral resolution, intended to detect the continuum emission, and five additional high-spectral resolution windows, centred at selected spectral line transitions, and aimed at studying the gas kinematics. The total effective bandwidth ($\sim$3.75 GHz) includes additional spectral lines (see Sect. 3.2).

For imaging, the task tclean of CASA (version 6.2.1.7; CASA Team et al., 2022) was used with a Brigg’s weighting of visibilities and robust parameter of 0.5, obtaining a synthesized beam size = $0\aas@@fstack{\prime\prime}165\times 0\aas@@fstack{\prime\prime}077$, with position angle (PA) $=-2\aas@@fstack{\circ}2$, and an rms $\simeq$ 0.08 mJy beam*-1* for the continuum image. For the spectral line data the rms is $\simeq$ 1 mJy beam*-1* per velocity channel of 0.78 km s*-1*, with synthesised beams comparable to that of the continuum image. All images presented here have been obtained after self-calibration and corrected by the primary beam response, which has a full-width at half-maximum (FWHM) $\simeq$ 25*′′* at 230 GHz (see Gómez et al., 2023, for more details).

3 Results

3.1 Dust 1.3 mm-continuum emission

Continuum images at a wavelength of 1.3 mm of the HMSFR W75N(B) are shown in Fig. 1 through zoomed-in views from the whole region of $25\arcsec\times 25\arcsec$ ($\simeq 0.16~{\rm pc}\times 0.16~{\rm pc}$), to a region of $5\arcsec\times 5\arcsec$ ($\simeq 0.03~{\rm pc}\times 0.03~{\rm pc}$) centered on VLA2, and down to a region of $0\aas@@fstack{\prime\prime}95\times 0\aas@@fstack{\prime\prime}95$ ($\simeq 1230~{\rm au}\times 1230~{\rm au}$) around VLA3, which is the brightest millimetre continuum source in W75N(B) (see Gómez et al., 2023). VLA3 shows very compact, bright continuum emission, with an observed size at FWHM ($0\aas@@fstack{\prime\prime}18\times 0\aas@@fstack{\prime\prime}09$; PA $=-2^{\circ}$) comparable to the synthesized beam size ($0\aas@@fstack{\prime\prime}17\times 0\aas@@fstack{\prime\prime}08$; PA $=-2^{\circ}$). The maximum intensity of the continuum emission in VLA3 is $I_{\nu}=97$ mJy beam*-1*; slightly different from that reported by Gómez et al. (2023) measured with a smaller beam ($I_{\nu}$ = 87 mJy beam*-1*; beam $=0\aas@@fstack{\prime\prime}13\times 0\aas@@fstack{\prime\prime}06$, PA $=-2^{\circ}$), obtained by excluding baselines shorter than 200 k$\lambda$ to better identify weak compact sources from extended emission. There is also relatively weak extended emission surrounding VLA3 at a level of $\sim$3 mJy beam*-1* (corresponding to about 30$\sigma$), and covering a region of \sim$$0\aas@@fstack{\prime\prime}5, extending preferentially in the northeast-southwest direction. The total flux density of VLA3 is $S_{\nu}\simeq 149\pm 1$ mJy, measured within the 3 mJy beam*-1* contour level shown in Fig. 1-right, which corresponds to an area of $\sim$0.16 arcsec2 and an effective circular radius of \sim$$0\aas@@fstack{\prime\prime}23 ($\sim$300 au).

Based on the flux density measured at a wavelength of 3.6 cm for the free-free emission of the thermal radio jet associated with VLA3 (4 mJy), and assuming a spectral index $\alpha=+0.6$ ($S_{\nu}\propto\nu^{\alpha}$) for the free-free continuum emission (Carrasco-González et al., 2010), the expected contribution of this component can be extrapolated to 1.3 mm, yielding a value of $\sim$30 mJy. Subtracting this contribution from the total emission observed at 1.3 mm infers a residual flux density of $\sim$120 mJy, attributable to thermal dust emission.

Assuming that the dust continuum emission at 1.3 mm is optically thin and the dust is isothermal within the observed region around VLA3, the dust mass content can be estimated from the expression:

[TABLE]

where $d$ is the distance to the source, $\kappa_{\nu}$ the dust opacity, and $B_{\nu}(T_{\rm dust})$ the Planck function at the dust temperature ($T_{\rm dust}$). Then, adopting as in Gómez et al. (2023) a dust opacity $\kappa_{\nu}$(1.3 mm) = 0.899 cm2 g*-1* (see Ossenkopf & Henning, 1994), $d=1.3$ kpc, $T_{\rm dust}=50$ K, and gas-to-dust ratio of 100, the following expression for the total mass of dust+gas is obtained:

[TABLE]

We derive $M_{\rm dust+gas}$ $\approx$ 1.74 $M_{\odot}$ from the measured flux density of 120 mJy for thermal dust emission, assuming $T_{\rm dust}=$ 50 K. This mass reduces to $\approx$ 0.43 $M_{\odot}$ if a higher dust temperature of 200 K is considered. Such a high temperature is expected in the inner regions of the source, based on the detection of vibrationally excited lines of different molecular species (see Sect. 3.2). This temperature range aligns well with the expected values derived from relations established by Tobin et al. (2020) for discs associated with low-mass protostars, following extrapolation to VLA3’s properties. We note that if the dust emission is not optically thin, the value of $M_{\rm dust+gas}$ will likely be higher (see e.g., Añez-López et al., 2020). Despite this, the inner, optically-thick portion of the core would also be expected to be hotter than the temperatures we considered, which would yield lower masses and partially compensate for the opacity effects. Given the current spatial resolution and the fact that detailed modelling of the density and temperature structure of the dust surrounding VLA3 is beyond the scope of this work, we have used the plausible range of temperatures indicated above.

From its continuum emission of $\approx 4$ mJy at 3.6 cm (see Carrasco-González et al., 2010; Rodríguez-Kamenetzky et al., 2020) and the empirical correlation between bolometric luminosity and radio continuum luminosity at centimetre wavelengths (see Anglada et al., 2018, and their Eq. 28), a luminosity of $\sim 4\times 10^{4}$ $L_{\odot}$ is estimated for VLA3, equivalent to a stellar mass of 17–19 $M_{\odot}$ (from the mass-luminosity relation for ZAMS stars, e.g. Salaris & Cassisi, 2005; Gräfener et al., 2011). This stellar mass is about one order of magnitude larger than the dust+gas mass derived for the dense core, but comparable to the dynamical stellar mass of $\sim$16 M*⊙* (see Sect. 4.1).

3.2 Molecular emission

Figure 2 presents ALMA spectra extracted towards VLA3. The flux density of these spectra was obtained by integrating over the region defined by the 3 mJy beam*-1* contour level of the 1.3 mm continuum emission (see Fig. 1-right panel). The figure comprises six panels, each displaying a different spectral window observed with ALMA (see Table 1), arranged in order of increasing frequency. Some of the most prominent spectral lines detected within these windows have been identified, and assigned to their corresponding molecular species. The spectral lines highlighted in blue in this figure are those specifically discussed below, in the context of the study on VLA3, and listed in Table 2.

A comprehensive analysis of the complete molecular line identification and inventory in the W75N(B) HMSFR will be presented in a forthcoming publication (Sánchez-Monge et al., in preparation). This study will focus on a region of $\sim$30 arcseconds ($\sim$39,000 au), encompassing the 40 millimetre continuum sources previously reported by Gómez et al. (2023).

3.2.1 Shock and outflow tracers: SiO (5–4)

The SiO molecule has long been considered a potential shock tracer of outflows (e.g., Schilke et al., 1997; Gusdorf et al., 2008; Sánchez-Monge et al., 2013b). In this study, we aim to investigate the outflow associated with the thermal radio jet, given the previous detection of SiO (1–0) with the VLA (Carrasco-González et al., 2010). Our ALMA observations of the SiO (5–4) line do not reveal any clear structure reminiscent of a typical bipolar outflow associated with VLA3. Instead, as previously reported by Gómez et al. (2023), a dominant absorption feature is detected in SiO (5–4) within the velocity range VLSR $\simeq$ 4–8 km s*-1*, coinciding with the brightest and most compact continuum emission of VLA3 (see the SiO velocity channel and moment maps in Figs. 15 and 16). Interestingly, the peak velocity of this absorption feature (VLSR $\simeq$ 6.3 km s*-1*) is blue-shifted by $\sim$2.5 km s*-1* with respect to the velocity of the ambient dense gas, as estimated by the kinematic analysis presented in Sect. 4. A strong absorption feature is also observed in H2CO within the same velocities as in SiO and coincident with the compact continuum emission (see Sect. 3.2.2).

The SiO velocity field derived from the first-order moment map (see left panel of Fig. 3) reveals a tentative velocity gradient with blue-shifted velocities to the northwest and red-shifted velocities towards the southeast. We note that with the current data, we can not rule out the possibility that the blue-shifted emission might be partially associated with extended SiO emission from VLA2 (see Fig. 16). Additionally, the presence of the strong absorption feature prevents us from conducting a more detailed analysis of the SiO kinematics in the near vicinity of VLA3. Interestingly, the orientation of this tentative SiO velocity gradient is consistent with the elongation of the thermal radio jet (PA $=157\pm 3^{\circ}$) reported by Carrasco-González et al. (2010) at 3.6 cm (see yellow ellipse and dashed line in Fig. 3). The central and right panels of Fig. 3 show the radio continuum emission from the radio jet at 1.3 cm and 7 mm observed with the VLA at higher angular resolution (data from Carrasco-González et al., 2015; Rodríguez-Kamenetzky et al., 2020). These two images have been generated with circular beams ($0\aas@@fstack{\prime\prime}12$ at 1.3 cm, and $0\aas@@fstack{\prime\prime}038$ at 7 mm) resulting in an angular resolution comparable or better than the ALMA observations presented here. The radio continuum emission at both 1.3 cm and 7 mm reveal an elongation in the southeast to northwest direction in agreement with the orientation derived by Carrasco-González et al. (2010) and with the tentative SiO velocity gradient. In summary, the SiO (5–4) observations together with the radio continuum emission maps support the presence of a jet powered by VLA3, at an orientation PA $\approx 157^{\circ}$.

3.2.2 Dense gas tracers: SO2, H2CO, HC3N, HNCO, CH3CN, H2O

Figure 4 displays velocity channel maps showing the spatial distribution of the SO2 ($22_{2,20}$–$22_{1,21}$ ; v2=1) emission. These maps cover the velocity range VLSR = $-1.82$ km s*-1* to $+18.17$ km s*-1*, with a channel spacing of $\approx$1.3 km s*-1*. The location of VLA3 is marked with a yellow cross in each velocity panel for reference. A distinct kinematic trend is evident in the SO2 emission: At velocities of $\simeq$ $-1.82$ to $+6$ km s*-1*, the emission peak in each channel is located southwest of VLA3. As the velocity increases towards higher values, from $\simeq$ $+8$ km s*-1* to $+18$ km s*-1*, the emission peak systematically shifts to the northeast of VLA3. This velocity gradient, perpendicular to the orientation of the thermal radio jet and SiO velocity gradient, suggests organised gas motions probably related to the rotation of dense gas around the massive protostar (see Sect. 4).

The SO2 emission is further investigated through moment maps, presented in Fig. 5. The figure displays the integrated intensity ($0^{\rm th}$ order moment: MOM0), velocity field ($1^{\rm st}$ order moment: MOM1), and velocity dispersion ($2^{\rm nd}$ order moment: MOM2) of the SO2 emission in separate rows. The integrated intensity map of SO2 (MOM0, Fig. 5 top-row) reveals a compact, elongated structure oriented in a northeast to southwest direction (PA $\simeq 62^{\circ}$), spanning \sim$$0\aas@@fstack{\prime\prime}7 (\sim$$900 au at a 3$\sigma$ level, slightly larger than the extent of the mm continuum emission, see Sect. 3.1) and exhibiting two intensity peaks separated by \sim$$0\aas@@fstack{\prime\prime}25 (\sim$$300 au), with the southern one being slightly brighter. A 2D Gaussian fit to the SO2 integrated intensity map results in a deconvolved FWHM of $0\aas@@fstack{\prime\prime}51\pm 0\aas@@fstack{\prime\prime}04\times 0\aas@@fstack{\prime\prime}34\pm 0\aas@@fstack{\prime\prime}03$, with PA $=58\pm 9^{\circ}$. The blue ellipse overlaid on the panel indicates the deconvolved size ($0\aas@@fstack{\prime\prime}21\times 0\aas@@fstack{\prime\prime}07$) and orientation (PA $=157^{\circ}$) of the VLA3 thermal radio jet previously reported by Carrasco-González et al. (2010) at 3.6 cm (see also Fig. 3). The thermal radio jet is situated between these two SO2 emission maxima, and displays an orientation perpendicular to the elongated structure of the gas emission. The velocity field of the SO2 emission, as depicted by the first-order moment (Fig 5 middle-row), shows a velocity gradient extending along the elongated structure. This velocity variation was identified earlier through the velocity channel maps shown in Fig. 4. Conversely, the second-order moment map (Fig. 5 bottom-row) reveals a velocity dispersion of $\sim$3–4 km s*-1*, with no significant variation observed along the molecular structure.

The emission from other dense gas tracer molecules identified towards VLA3 (H2CO, HC3N, HNCO, CH3CN; transitions listed in Table 2) displays a spatial and kinematic distribution which in general is similar to that observed in SO2 (see Appendix A, and Figs. 17 to 24), and therefore is likely associated with intrinsic properties of the circumstellar dense gas around the VLA3 protostar. However, these lines, which have lower excitation energies (see Table 2), exhibit a more extended distribution across the systemic velocity (V${}_{\mathrm{LSR}}\simeq$ 8.9 km s*-1*) range of VLA3, likely tracing the surrounding local environment. Yet, their extreme velocity channels distinctly reveal, as in SO2, a velocity gradient oriented southwest to northeast, suggestive of rotation motions.

In particular, the H2CO emission towards VLA3 is detected across a broad velocity range, from VLSR $\simeq-1.7$ to $+18$ km s*-1* (see Fig. 17) similar to that of the SO2 emission (Fig. 4). The integrated intensity (MOM0) of H2CO reveals an elongated structure oriented northeast to southwest, featuring two intensity peaks and a velocity gradient comparable to that observed in SO2 (compare Figs. 18 and 5). Nevertheless, this velocity gradient is partially obscured by intense extended H2CO emission spanning the entire W75N(B) region within V${}_{\rm LSR}\simeq+5$ to $+13$ km s*-1* (see Fig. 17, and Gómez et al. 2023), and by a strong absorption feature at V${}_{\rm LSR}\simeq+6$ km s*-1* (see Fig. 17), analogous to the SiO absorption feature. The blue-shifted absorption seen in the SiO and H2CO channel maps is reminscent of the classical P-Cygni profiles and may suggest the presence of outflowing or expanding motions in the molecular gas around VLA3. We note that the peak brightness temperature of the 1.3 mm continuum emission is $\approx 100$ K, which is higher than the upper energy levels of the SiO and H2CO transitions (about 20–30 K, cf. Table 2) but lower than the energy levels of the other molecular lines studied in this work (between 130 and 3500 K). Higher spatial resolution observations, capable of resolving the most compact structures, may reveal a much brighter continuum source against which other molecular species might also show self-absorbed profiles.

We also detected the vibrationally excited H2O (55,0–64,3 ; $v_{2}$=1) transition at 232.6867 GHz (see Fig. 2 and Table 2). Vibrationally excited H2O lines have been identified as excellent tracers of the inermost, and hottest regions of discs surrounding (massive) YSOs (see e.g., Ginsburg et al., 2018; Maud et al., 2019). Unfortunately, this transition is covered in our setup within the spectral window spw 0 with a low spectral resolution ($\approx 20$ km s*-1*). Figure 6 shows the H2O emission map through the three velocity channels at which we detect emission from the water line. The emission is compact and spatially coincident with the bright 1.3 mm continuum source, with the emission in the central velocity channel slightly elongated in the north-south direction (deconvolved size $0\aas@@fstack{\prime\prime}89\times 0\aas@@fstack{\prime\prime}29$, with PA $=0^{\circ}$). Interestingly, and despite the coarse spectral resolution, the H2O emission peak shifts from west to east with respect to VLA3 when moving from blue- to red-shifted velocities. This is consistent with the shift observed for the other dense gas tracers (see e.g., Fig. 17 and figures in Appendix A), and supports the idea that this H2O vibrational line may also trace rotational motions around VLA3 as found in other massive protostars. However, higher spatial and spectral resolution observations are necessary to properly characterize the structure and kinematics traced by the H2O line.

3.2.3 Ionised gas tracers: H30$\alpha$ recombination line

The hydrogen radio recombination line (RRL) H30$\alpha$, which indicates the presence of ionised gas, falls within the frequency range of our observations, albeit within the spectral window spw 0 which was optimised for continuum observations, and consequently has a lower spectral resolution ($\approx 20$ km s*-1*) compared to the other spw’s (see Table 1). The spectrum integrated towards VLA3, as shown in Fig. 2, displays a spectral feature at about 231.9 GHz that coincides with the H30$\alpha$ frequency. We have searched for other molecular species with spectral lines near the H30$\alpha$ frequency that might account for the detected line. The only possible contaminant is 33SO2, with a group of transitions (3${}_{9,{\rm K_{c}}}$–2${}_{10,{\rm K_{c}}}$) located between 231.89 and 231.90 GHz. However, using a typical isotopic ratio of 32S/33S of 88 (see Yan et al., 2023), the expected contribution from this contaminant to the observed emission line would be only around 15%, based on the intensity of the SO2 lines detected within the same spectral range. Based on this, we assign the bright spectral feature at 231.9 GHz to the H30$\alpha$ line.

Figure 7 shows a close-up view of the spectrum around the frequency of the H30$\alpha$ line. The top panel shows the spectrum integrated over the 3 mJy beam*-1* continuum level (see Fig. 1) used to determined properties of the VLA3 continuum emisison, while the bottom panel shows the spectrum extracted towards the peak position of the continuum source. In order to determine the intensity and linewidth of the H30$\alpha$ line, we have fitted two Gaussians to account for the spectral line features associated with the H30$\alpha$ line itself and with the neighboring SO2 (143,11–42,12; v2=1) line at 231.980 GHz. The best fit is shown in blue in Fig. 7. The fit to the integrated spectrum delivers for the H30$\alpha$ line a peak intensity of $70\pm 5$ K (or $42\pm 6$ mJy beam*-1*), centered at $20\pm 4$ km s*-1* (with respect to the H30$\alpha$ rest frequency), and with a velocity dispersion111This velocity dispersion, $\sigma$, corresponds to a linewidth or full-width at half-maximum (FWHM) of the line, $\Delta v=\sigma\cdot(8\ln 2)^{1/2}$, equal to $98\pm 10$ km s*-1*. of $42\pm 8$ km s*-1*. The Gaussian fit of the neighboring SO2 feature results in a peak intensity of $87\pm 7$ K, centered at $8\pm 2$ km s*-1* (with respect to the SO2 rest frequency of 231.980 GHz), and with a velocity dispersion of $21\pm 4$ km s*-1* (likely severely affected by the coarse spectral resolution of $\approx 20$ km s*-1*). Similarly, the fit to the spectrum extracted towards the peak results in a H30$\alpha$ line intensity of $27\pm 2$ K, centered at $12\pm 2$ km/s, and with a velocity dispersion of $40\pm 6$ km s*-1*. The fit parameters for the SO2 line are an intensity of $11\pm 2$ K, centered at $4\pm 4$ km s*-1*, with a dispersion of $26\pm 11$ km s*-1*. The spectra shown in Fig. 7 shows a tentative departure from Gaussianity at high velocities, likely tracing high-velocity wings. This high-velocity component is also seen in the channel emission maps shown in Fig. 25. In this figure, the emission is dominated by a compact, bright component associated with the H30$\alpha$ line and with its peak at the velocity channel 8.16 km s*-1* (consistent with the systemic velocity of VLA3). This component is surrounded by an extended structure, which is brightest at the velocity channel 28.35 km s*-1* and is likely related to the blended spectral transitions of 33SO2. The compact component that traces the H30$\alpha$ line emission extends up to high velocities of $\pm 150$ km s*-1*.

In the left panels of Fig. 8, we present the emission of the H30$\alpha$ line, overlaid with the 1.3 mm continuum emission contours and the radio jet location identified by Carrasco-González et al. (2010). The peak of the H30$\alpha$ emission is spatially coincident with the 1.3 mm continuum peak and the radio jet, although there is also weaker emission distributed across the map, particularly towards VLA1. The weaker, extended emission may be related to other molecular species that are blended with the H30$\alpha$ line due to the limited spectral resolution. The right panels of Fig. 8 depict the spectra extracted at the position of the radio jet (top) and at a location shifted southwest by $0\aas@@fstack{\prime\prime}07$ ($\approx 95$ au; bottom). The broad spectral line feature corresponding to H30$\alpha$ is most pronounced at the radio jet’s position, nearly vanishing when slightly offset from it, thereby supporting its spatial association with the VLA3 jet and mm continuum peak. Additionally, the consistent intensities of other spectral lines across both positions (cf. Fig. 8) further support the H30$\alpha$ association with the jet or the central region of VLA3, since these other lines are associated with molecular species that trace the whole extent of the dense core and not only the inner region. It is noteworthy that the lines from other molecular species are considerably weaker than the H30$\alpha$ line at the location of the radio jet, with only the vibrationally excited H2O transition at 232.868 GHz remaining prominent (see Sect. 3.2.2 and Fig. 6). However, and contrary to the H2O line, the H30$\alpha$ channel maps shown in Fig. 25 do not reveal a clear spatial shift of the peak emission in the east-west direction. Interestingly, a tentative north-south spatial shift, coinciding with the radio jet direction, is observed instead. This tentative spatial shift is shown in Fig. 9, where the colored circles depict the H30$\alpha$ peak position for each velocity channel, with the colors corresponding to the channel velocity (ranging from $-150$ to $+150$ km s*-1*). The north-south spatial shift also shows a tentative velocity gradient (with red- and blue-shifted velocities located to the north and south of VLA3, respectively), which, if real, would be opposite to the velocity gradient seen in SiO (see Fig. 3). We note that spatial overlap of blue and red-shifted gas can occur in outflows observed close the plane of the sky (see e.g., Cruz-González et al., 2007; Skretas et al., 2023, for the case of DR21), as might be the case for the jet/outflow in VLA3. However, we note that the tentative spatial shift and velocity gradient observed in the H30$\alpha$ line must be taken cautiously, as the limited angular resolution of the observations, especially in the north-south direction, hinders a more precise confirmation. In Sect. 4.3, we discuss the properties of the H30$\alpha$ emission and its potential connection to the radio jet powered by VLA3.

4 Discussion

4.1 A compact, rotating disc around W75N(B)-VLA3

To comprehensively analyse the spatial and kinematic distribution of the dense gas tracers discussed in Sect. 3.2.2, we constructed the position-velocity (PV) diagrams shown in Fig. 10. The position angle of the PV cut direction is PA $=62^{\circ}$, and is determined by the line that connects the two intensity maxima observed in the integrated intensity (MOM0) images of different dense gas tracers (see e.g., Fig. 5). This orientation is in agreement with the direction of the velocity gradient observed in the velocity field (MOM1) maps (see e.g., Figs. 5 and 18). Furthermore, the PV cut is almost perpendicular to the direction of the thermal radio jet (PA $=157\pm 3^{\circ}$; Carrasco-González et al. 2010), thus favoring the interpretation of the observed velocity gradient tracing rotation motions perpendicular to the jet direction. The center of the PV cut corresponds to the location of VLA3, with coordinates $\alpha(\mathrm{J2000})=20^{\mathrm{h}}38^{\mathrm{m}}36.4815^{\mathrm{s}}$ and $\delta(\mathrm{J2000})=+42^{\circ}37^{\prime}33.355^{\prime\prime}$. The PV diagrams shown in Fig. 10 were generated averaging over nine pixels (of $0\aas@@fstack{\prime\prime}016$ each) across the direction of the cut to enhance the sensitivity.

We have used the approach described in Seifried et al. (2016, see also the implementation by , and the KeplerFit Python-based package222KeplerFit, https://github.com/felixbosco/KeplerFit) to fit a Keplerian disc-like rotating model to the observed PV diagrams. The white, solid lines in Fig. 10 show the best fits of the Keplerian rotation around a central object of mass $M_{*}$. This mass can be determined from the velocity distribution of the material around the central object, since the highest velocity at a given radial distance $r$ to the center of mass is limited by the Kepler orbital velocity $v_{\mathrm{Kepler}}(r)\approx\sqrt{GM_{*}/r}$, where $G$ is the gravitational constant. Based on that, we have estimated the highest velocity at a given radial distance (i.e., the outer envelope velocity) as the velocity channel at which the emission in the PV cut lies just above a threshold of 4 times the rms noise level (i.e., 4$\sigma$). This method was found to be the most robust among those tested by Seifried et al. (2016) on synthetic data where the true mass of the central object was known. The lime and cyan symbols in the PV diagrams of Fig. 10 mark the determined highest velocities ($v$) at the different offsets ($r$). The lime symbols are those ($v,r$) pairs used to fit the Keplerian disc-like rotation model, while the cyan symbols are those pairs that were discarded during the fitting process. These discarded pairs correspond to the northeastern region (beyond 500 au or $0\aas@@fstack{\prime\prime}4$) which is highly contaminated by extended and bright emission likely associated with larger-scale structures around VLA3 (see Sect. 4.1.1), and to the inner region with radius 100 au around the central object, which is not properly resolved with the current angular resolution.

The fitted Keplerian disc-like rotation model is given by

[TABLE]

where $M_{*,\mathrm{fit}}$ is the mass of the central object, which is assumed to be larger than the disc mass (see Bosco et al., 2019). The sign-function is introduced to account for the opposite velocity difference from the blue- and red-shifted emission. The parameter $v_{0}$ provides the velocity shift or systemic velocity of the object, while the parameter $r_{0}$ provides the positional shift with respect to the central position of the PV-cut (i.e., VLA3). We note that masses obtained from this method are uncertain due to the unconstrained disc axis inclination $i$ to the line of sight. The actual central protostellar mass ($M_{*}$) can be obtained from $M_{*,\mathrm{fit}}=M_{*}\sin^{2}{i}$ (see Appendix A in Bosco et al., 2019) if the inclination angle is known. In the following analysis we consider $i=90^{\circ}$ (i.e., edge-on333An inclination angle of 70∘ or 50∘, instead of 90∘ (edge-on), would result in a central protostellar mass about 13% or 70% larger, respectively.). This is based on the obscured (radio bright) HH objects detected up to distances of 5200 au from VLA3 and with high-velocity proper motions exceeding 100 km s*-1* (Rodríguez-Kamenetzky et al., 2020), which suggest that the outflow ejecta are most probably close to the plane of the sky, and therefore the circumstellar disc is almost edge-on. The deconvolved size of the SO2 integrated intensity emission (see Sect. 3.2.2) can also be used to infer the inclination angle of the disc under the assumption that the observed ellipticity is only due to the projection of a circular disc observed at a given inclination. Based on this, we obtain $b\approx a\cos(i)$, where $a$ and $b$ are the major and minor axis of the ellipse, which results in $i=50^{\circ}\pm 15^{\circ}$ for the SO2 structure. We note that any departure from a infinitely thin disc would result in a lower limit to the derived inclination angle, i.e., the disc would in reality be closer to edge-on. By considering an edge-on disc ($i=90^{\circ}$), we set a lower limit to the actual central protostellar mass, and thus we undertake a conservative approach in the following analysis (e.g., Sect. 4.2 and 4.3).

In Table 3 we list the fitted parameters that reproduce the Keplerian curves shown in Fig. 10. The average central mass $M_{*,\mathrm{fit}}$ is $16\pm 4$ $M_{\odot}$, after considering the five fitted spectral lines444An inclination angle of $50^{\circ}$ results in a central mass $M_{*}\approx 27$ $M_{\odot}$.. We note that the PV diagram of the vibrationally-excited SO2 transition stands out as the best of the five diagrams to study the kinematics of the compact rotating structure in VLA3. This is because the H2CO and HC2N transitions (with energy levels $\lesssim 130$ K) are more severely affected by extended emission likely not directly related to the disc, and the HNCO and CH3CN transitions lack sensitivity. Based solely on the SO2 transition, we also derive a central mass of $\approx 16$ $M_{\odot}$. This mass is consistent with the stellar mass of 17–19 M⊙ estimated from the correlation between the bolometric luminosity and the radio continuum luminosity at centimetre wavelengths (see Sect. 3.1). The average reference velocity $v_{0}$ is $8.9\pm 0.6$ km s*-1*, and determines the systemic velocity of VLA3. The average offset $r_{0}$ is $-0\aas@@fstack{\prime\prime}019\pm 0\aas@@fstack{\prime\prime}015$ (or $-25\pm 19$ au) with respect to the center of the PV cut (i.e., the coordinates of VLA3), and thus, it is in agreement with the center of mass being located at the position of VLA3. In Fig. 11, we explore how the derived central mass is affected by the threshold level used when determining the outer envelope velocity in the PV diagrams (see also Ginsburg et al., 2018). For this analysis, we focus on the SO2 PV diagram and vary the threshold level between $3\sigma$ and $8\sigma$. The results suggest an uncertainty of about $\pm 3$ $M_{\odot}$ on the derived central mass due to the chosen threshold level.

In summary, the dense gas tracers around VLA3 show emission consistent with a Keplerian-like rotating disc centered on an object of $\approx 16$ $M_{\odot}$ (assuming an edge-on disc). Based on the 1.3 mm continuum emission (see Sect. 3.1), this disc-like structure has an estimated dust+gas mass (hereafter $M_{\mathrm{gas}}$) of 0.43–1.74 $M_{\odot}$ (i.e., about 10–40 times less than the central object) and extends to an estimated radius between $\approx 300$ au (based on the mm continuum) and $\approx 450$ au (based on the molecular lines).

4.1.1 Morphology and substructures in the W75N(B)-VLA3 disc

As described in Sect. 3.2.2, the circumstellar disc surrounding VLA3 exhibits two peaks in the emission maps of the different dense gas tracers, with the protostar situated centrally between them and with the south-western peak being brighter (see Fig. 5 and figures in Appendix A). The two peaks with their different brightness suggests that the disc is not symmetrical. Such asymmetry has also been observed in other massive protostellar discs, as recently reported for Cepheus A HW2 (Sanna et al., 2025). In that case, the asymmetry was attributed to anisotropic accretion streams of gas and dust flowing from the outer envelope into the inner disc. In the case of W75N(B)-VLA3, these streamers are not well identified in the currently available data. However, an elongated, slightly curved structure to the north-east of VLA3 is visible in multiple species (see, e.g., moment maps of HC3N and H2CO in Figs. 20 and 18, respectively), and appears in the PV-plots shown in Fig. 10 as an excess of emission at an offset $>0\aas@@fstack{\prime\prime}4$ and at velocities $\sim 15$ km s*-1* (i.e., red-shifted with respect to the reference velocity of 8.9 km s*-1*). The elongated structure seen in HC3N (cf. Fig. 20) is reminiscent of similar structures found towards other young stellar objects and related to accretion streamers (e.g., Alves et al., 2019; Akiyama et al., 2019; Pineda et al., 2020; Garufi et al., 2022; Sanna et al., 2025). For VLA3, and based on the red-shifted velocities of the northern elongated structure, together with the location of the brightest peak in the southern, blue-shifted region of the disc, one might imagine the existence of an anisotropic accretion streamer that brings gas from the north-eastern part of VLA3 (from the observer’s perspective) and finally impacts the disc in the southern region, where the molecular emission is brighter. Figure 12 depicts a sketch of this hypothesis. However, a more detailed study of the extended emission around VLA3 together with new molecular observations that better trace the lower-density gas in the expected envelope or accretion streamer would be required to confirm this scenario.

Finally, as noted previously, we observed two peaks in the molecular line emission flanking the central source (cf. Figs. 5 and 18). While this morphology might initially suggest a cavity devoid of molecular gas, such as one created by a toroidal gas distribution, we believe a more plausible explanation might lie in the effect of dust opacity. As shown in the upper-right panel of Fig. 5, the continuum emission fills the region corresponding to this apparent molecular cavity. This emission, located between the two molecular line peaks, indicates the presence of dust555The presence of a real cavity, devoid of both gas and dust, would imply that the observed 1.3 mm continuum emission is not related to dust. Instead, the continuum emission would be dominated by thermal free-free emission from the radiojet. This scenario, however, appears unlikely based on the extrapolation of $\sim 20$% for the free-free contribution to the flux at 1.3 mm. in the central region of the disc and therefore, it is expected to not be completely empty. If the continuum emission from this dust is optically thick ($\tau_{\rm dust}>1$), the molecular line emission from the same region could not be observed, since the continuum already reaches the maximum possible brightness temperature. The line emission would only be observed in the outer regions, where the dust emission is optically thin. The result is an apparent hole in molecular line emission after continuum subtraction, even though the molecular gas is present and can actually be dense. This ‘false’ molecular cavity is naturally obtained in models that consider both line and continuum emission in circumstellar discs (see, e.g., models by Gómez & D’Alessio, 2000), as an effect of high dust opacity. Based on the intensity peak of the 1.3 mm continuum emission towards VLA3 (see Sect. 3.1), we estimate a beam-averaged dust opacity of $\sim$1.2 (see Frau et al., 2010, their Appendix A) consistent with the dust being optically thick. To confirm this interpretation for the apparent cavity, observations of the dust emission around VLA3 at different wavelengths are needed to better assess its opacity. In particular, observations at longer wavelengths are expected to mitigate the effects of dust opacity, thus enabling the study of the dust and gas distribution in the central region of VLA3. In this regard, it is noteworthy that the peak emission of the vibrationally excited H2O line (see Fig. 6) is coincident with the 1.3 mm continuum source. This suggests that the H2O emission arise from very hot gas in the innermost circumstellar regions. Our limited angular precludes us from ascertaining the origin of the emission. However, as mentioned in Sect. 3.2.2, the spatial shift seen in different velocity channels suggests that part of the emission may arise from the innermost region of the circumstellar disc. The warm inner region may be mostly devoid from dust due to sublimation, so that the dust opacity does not block the line emission. Since no inner hole is seen in the continuum map, this putative dust-free area must be unresolved in our images. Moreover, part of the H2O emission may also arise from the outflow, as suggested by the apparent north-south elongation of the line emission close to the direction of the jet (see Fig. 6). Higher spatial and spectral resolution observations are necessary to shed light on the structure and morphology traced by this H2O line around VLA3.

4.2 Comparison with other discs/toroids around massive protostars

In Figs. 13 and 14 we show the properties of the rotating structure found around W75N(B)-VLA3 alongside those of other potential rotating discs or toroids associated with intermediate and high-mass YSOs. For this, we use the review by Beltrán & de Wit (2016), as well as recent observational results obtained at high-angular resolution (see Johnston et al., 2020; Ginsburg et al., 2023, and references therein). Generally, we refer to structures smaller than 1000 au as discs, and those larger than 5000 au as toroids (see also below). Note that the work by Beltrán & de Wit (2016) primarily consists of data collected before the recent ALMA era of high-angular resolution studies. These earlier studies typically used low-energy transitions of dense gas tracers such as CH3CN and CH3OH to characterize the dense, rotating structures around massive YSOs. In contrast, the higher resolution and sensitivity observations of the new studies, enable the use of high-excited transitions of species such as H2O, NaCl, and KCl (e.g., Ginsburg et al., 2023). For the comparison of VLA3 with other sources, we follow the approach used by Beltrán & de Wit (2016) and use the gas mass contained in the rotating structure ($M_{\mathrm{gas}}$) and its radius, the central (proto)stellar mass ($M_{\mathrm{star}}$) derived from luminosity estimates, the dynamical mass ($M_{\mathrm{dyn}}$) derived from PV-plot analysis, and different timescales to evaluate the stability and dynamical status of the rotating object around VLA3.

In panel (a) of Fig. 13, we compare the gas-to-star mass ratio ($M_{\mathrm{gas}}/M_{\mathrm{star}}$) as a function of the radius of the rotating structure. The location of VLA3 in the diagram is in agreement with the correlation found between the distance-independent $M_{\mathrm{gas}}/M_{\mathrm{star}}$ and the radius of the circumstellar structures of other intermediate-mass (IM) and high-mass (HM) YSOs. As stated in Sect. 3.1, the mass of the rotating structure around VLA3 ($M_{\mathrm{gas}}\approx 0.43$–1.74 $M_{\odot}$ based on the mm continuum emission) is about one tenth of the mass of the central object ($M_{\mathrm{star}}\approx 17$–19 $M_{\odot}$ based on the relation between the radio and bolometric luminosities666We note that using the dynamically-derived mass of $\approx 16\pm 4$ $M_{\odot}$ for VLA3 (see Sect. 4.1) does not significantly impact the results discussed here, since both estimates for the stellar masses are in good agreement.). This $M_{\mathrm{gas}}/M_{\mathrm{star}}$ ratio estimated for VLA3 is similar to those recently found in other accretion discs around massive protosars (see filled blue symbols in Fig. 13; Johnston et al., 2020; Sanna et al., 2025). Panel (b) illustrates this mass relationship for VLA3 in comparison to the other rotating structures. The figure shows that the mass of the rotating structure (disc or toroid) is in general proportional to the mass of the central object. It is worth noting that in the earlier studies reported in Beltrán & de Wit (2016), a significant number of rotating structures around IM and HM objects have $M_{\mathrm{gas}}>M_{\mathrm{star}}$, in particular for structures with radii generally $\gtrsim 300$ au. Consequently, these rotating structures are likely to be self-gravitating and are probably not in Keplerian rotation, as the gravitational potential of the system is primarily influenced by the surrounding circumstellar structure rather than the central (proto)stellar mass. Based on different studies (e.g. Shu et al., 1990; Laughlin & Bodenheimer, 1994; Cesaroni et al., 2007) discs can be stable if their mass is $\lesssim 0.3M_{\mathrm{star}}$. For larger disc masses (i.e., $>0.3M_{\mathrm{star}}$), gravitational instabilities might develop and rapidly fragment the disc (e.g., Rosen et al., 2019; Ahmadi et al., 2019). This is expected to happen in the massive and rotating large structures found around massive YSOs, the so-called toroids (see, e.g., Beltrán et al., 2011). These toroids, although not being stable entities, might host stable Keplerian-rotating discs in their interiors, which would be only accessible with high-angular resolution observations (e.g., Sánchez-Monge et al., 2013a; Johnston et al., 2015; Cesaroni et al., 2017). The existence of these stable structures is seen in the recent higher-angular resolution studies (see e.g., Johnston et al., 2020; Ginsburg et al., 2023) that are able to resolve and trace the inner regions of the rotating structures (see filled blue circles in Fig. 13 panel b). These observations show that in most cases these discs have masses less than 0.3 times that of the central (proto)star, as also found towards VLA3, suggesting that they are likely close to Keplerian rotation (with the mass of the central object being the dominant contributor). The high angular resolution observations of W75N(B) conducted with ALMA, have resolved the rotating disc in VLA3, adding it to the list of about a dozen such discs identified around massive protostars so far.

The stability of these structures can be evaluated by estimating the Toomre’s $Q$ parameter (Toomre, 1964, see also Ahmadi et al. 2023), which is defined as

[TABLE]

where $c_{\mathrm{s}}$ is the sound speed777The sound speed of the gas is calculated as $\sqrt{(k_{\mathrm{B}}T)/(\mu m_{\mathrm{H}})}$, where $k_{\mathrm{B}}$ is the Boltzmann constant, $m_{\mathrm{H}}$ is the hydrogen mass, and $\mu$ the mean mass assumed to be 2.3. For the temperature $T$ of the gas we use the values 50 and 200 K as in Sect. 3.1 of the gas, $\Sigma$ the surface density, and $\kappa$ the epicyclic frequency of the disc. For Keplerian rotation, the epicyclic frequency corresponds to the angular velocity $\Omega$ (see Beltrán & de Wit, 2016), which following Cesaroni et al. (2007) can be estimated as $\sqrt{GM_{\mathrm{total}}/R^{3}}$, where $M_{\mathrm{total}}$ is the star plus disc total mass of the system and $R$ is the disc radius. Following Beltrán & de Wit (2016), the surface density is estimated as $\Sigma=M_{\mathrm{gas}}/(\pi R^{2}$). For a stellar mass of 16 $M_{\odot}$ and a disc radius of 450 au, we derive a Toomre’s $Q$ parameter of 0.8 if the disc has a temperature of 50 K and a mass of 1.74 $M_{\odot}$, or a value of 6.2 if we consider a temperature of 200 K and a resulting disc mass of 0.43 $M_{\odot}$. These values are shown in the panel (c) of Fig. 13 and compared to other discs and toroids. Based on this, VLA3 appears stable against gravitational instabilities, with a Toomre’s $Q$ parameter comparable to that of the dozen discs resolved around massive YSOs so far. Panel (c) of Fig. 13 shows how high-resolution observations are able to detect discs that populate the upper-right area of the diagram, which corresponds to stable discs rotating around massive stars ($M_{\mathrm{star}}>10$ $M_{\odot}$).

In Fig. 14, we evaluate the dynamical status of the disc around VLA3. In panel (a), we compare the dynamical mass to gas mass ratio of the discs/toroids against the central (proto)stellar mass. The dynamical mass ($M_{\mathrm{dyn}}$) is calculated assuming equilibrium between centrifugal and gravitational forces as $M_{\mathrm{dyn}}=v_{\mathrm{rot}}^{2}R/G\sin^{2}{i}$, where $v_{\mathrm{rot}}$ is the rotational velocity at a radius $R$, which for VLA3 is estimated to be 6 km s*-1* at 450 au, which results in a mass of $\approx 16$ $M_{\odot}$ as derived from the Keplerian fits (see Sect. 4.1). The dashed, horizontal line in the figure distinguishes two regimes: structures with $M_{\mathrm{dyn}}/M_{\mathrm{gas}}>1$ that are centrifugally supported, and objects with $M_{\mathrm{dyn}}/M_{\mathrm{gas}}<1$ that are gravitationally dominated. The high $M_{\mathrm{dyn}}/M_{\mathrm{gas}}$ value for VLA3 places this disc as a dynamically stable rotating structure comparable to some of the newest discs detected around massive YSOs.

Following Beltrán et al. (2011) and Sánchez-Monge et al. (2014), in Fig. 14 panel (b), we plot the ratio of the free-fall timescale ($t_{\mathrm{ff}}$) to the rotational period ($t_{\mathrm{rot}}$) against $M_{\mathrm{gas}}$ of the rotating structure. The free-fall timescale is proportional to the dynamical timescale needed to refresh the material of the rotating structure, while the rotational period at the outer radius is the time needed by the rotating structure to stabilise after incorporating new accreted material. Note that $t_{\mathrm{ff}}$ could be slightly overestimated since spherical symmetry is assumed (see also Beltrán & de Wit, 2016). The dashed curves in the figure correspond to the theoretical $t_{\mathrm{ff}}/t_{\mathrm{rot}}$ curves for spherical clouds of mass $M_{\mathrm{gas}}$ containing a star of mass $M_{*}$ at the center, in which the gas is rotationally supported against the gravity of the gas plus the star. These curves are expressed as $t_{\mathrm{ff}}/t_{\mathrm{rot}}=[(M_{\mathrm{gas}}+M_{*})/(32\,M_{\mathrm{gas}})]^{1/2}$ (see Beltrán et al., 2011; Sánchez-Monge et al., 2014). The results shown in Fig. 14 confirm that the rotating structure around VLA3 is in agreement with other centrifugally-regulated circumstellar discs, and distinguishable from gravitationally-unstable massive toroids.

Overall, the high-angular resolution study of the disc in W75N(B)-VLA3 contributes to the growing number of recently identified discs surrounding massive YSOs. It also suggests that highly excited transitions of molecules such as SO2, H2O, or NaCl might be more effective than the traditionally used transitions of complex organic molecules like CH3CN and CH3OH for resolving the inner regions of gravitationally-stable, Keplerian-like rotating discs.

4.3 Origin of the RRL emission in W75N(B)-VLA3

In Sect. 3.2.3, we have reported the detection of the RRL H30$\alpha$ emission towards VLA3. The limited spectral resolution, approximately 20 km s*-1*, hampers our ability to fully characterise the shape and velocity structure of the RRL; however, the spatial relationship with both the bright 1.3 mm continuum source and the radio continuum source reported by Carrasco-González et al. (2010) suggests that the H30$\alpha$ emission could be associated with the VLA3 thermal radio jet. Specifically, according to Anglada et al. (2018), the ratio of the RRL and continuum flux densities, $S_{\mathrm{L}}/S_{\mathrm{C}}$, for a biconical thermal jet with constant velocity and ionised fraction is given by

[TABLE]

after considering an ionised helium to ionised hydrogen ratio of 0.1, with $\nu_{\mathrm{L}}$ the frequency of the RRL, $T_{\mathrm{e}}$ the electron temperature of the ionised gas, and $\Delta v$ the full width of the RRL. The radio continuum emission at 2 cm has a flux density of 5.7 mJy (see Carrasco-González et al., 2010) which is extrapolated to a continuum flux density of 25 mJy at 1.3 mm, assuming a spectral index of $+0.6$ (see also Sect. 3.1). Using Eq. (5), one expects an H30$\alpha$ thermal flux density of $\approx 37$ mJy, assuming a temperature of 10 000 K and a line width of 100 km s*-1*. This expected H30$\alpha$ flux density is in good agreement with the observed line intensity (see Sect. 3.2.3), after considering that about 15% of the observed flux comes from a blended 33SO2 line. Removing this contamination from the measured flux results in an observed H30$\alpha$ flux density of $\approx 36$ mJy, consistent with the theoretical expectation.

The observed velocity of the H30$\alpha$ line, which is close to the ambient velocity (see Sect. 3.2.3), could also be consistent with a jet moving almost on the plane of the sky. In addition, the observed RRL line width of $\simeq 100$ km s*-1*, might be explained in an approximate way by considering a conical jet with an aperture semi-angle $\theta_{\mathrm{j}}$, following Rivera-Ortíz, Lizano and Cantó (2025; in preparation). For this, we assume that the flow within the jet is radial with respect to the powering source, and has a constant and uniform velocity $\mathrm{V}_{\mathrm{j}}$. We consider an observer’s line of sight that makes an angle $\theta_{\mathrm{o}}$ with the jet’s symmetry axis. Based on these assumptions, an emission line originating from the jet will be confined to the line-of-sight velocity interval $[\mathrm{V}_{\mathrm{min}},\mathrm{V}_{\mathrm{max}}]$, centered at $\mathrm{V}_{\mathrm{c}}$ = $(\mathrm{V}_{\mathrm{min}}+\mathrm{V}_{\mathrm{max}})/2$. The minimum and maximum velocities are given by

[TABLE]

respectively. The expression for $\mathrm{V}_{\mathrm{max}}$ is valid for $\theta_{\mathrm{o}}>\theta_{\mathrm{j}}$. For $\theta_{\mathrm{o}}<\theta_{\mathrm{j}}$, we have $\mathrm{V}_{\mathrm{max}}=\mathrm{V}_{\mathrm{j}}$. Therefore, the emission line will have a full width at zero intensity of

[TABLE]

centered at

[TABLE]

A counter-jet (i.e., a jet with identical physical conditions but moving in an opposite direction) will produce a ‘mirror’ emission line in the velocity interval $[-\mathrm{V}_{\mathrm{min}},-\mathrm{V}_{\mathrm{max}}]$, centered at $\mathrm{V}_{\mathrm{c}}=-\mathrm{V}_{\mathrm{j}}\cos(\theta_{\mathrm{o}})\cos(\theta_{\mathrm{j}})$, with the same full width as the jet’s line. If the jet is unresolved and its axis is close to the plane of the sky (i.e., $\theta_{\mathrm{o}}\simeq\pi/2$), the two lines (those from the jet and counter jet) will blend in a single line, centered at zero velocity (the velocity of the system) and having a full width at zero intensity of

[TABLE]

Applying this scenario to our observations, if we take a velocity for the jet flow of $\mathrm{V}_{\mathrm{j}}=500$ km s*-1* (similar to those observed in other massive protostars, e.g., Curiel et al., 2006; Anglada et al., 2018), and a line width of $100$ km s*-1* (or $300$ km s*-1*, based on the high-velocity emission seen in Fig. 25), we obtain $\theta_{\mathrm{j}}\simeq 6^{\circ}$ (or $\simeq 17^{\circ}$ after considering the highest observed velocities), which is consistent with the aperture semi-angle derived from the size of the radio jet at a wavelength of 3.6 cm ($\simeq 9^{\circ}$; Carrasco-González et al., 2010; Rodríguez-Kamenetzky et al., 2020). Similarly, considering an inclination angle $\theta_{\mathrm{o}}\approx 50^{\circ}$ (see Sect. 4.1) and the aperture semi-angle $\theta_{\mathrm{j}}\simeq 9^{\circ}$ measured for the VLA3 jet (Carrasco-González et al., 2010), a jet flow of $\mathrm{V}_{\mathrm{j}}\approx 70$–$200$ km s*-1* would suffice to reproduce the observed H30$\alpha$ line width of $100$ and $300$ km s*-1*. A more detailed model, able to draw more robust conclusions, would require H30$\alpha$ line data with higher spectral and angular resolution than our current data provide.

All this suggests that W75N(B)-VLA3 is likely the first radio jet in which a thermal RRL has been detected. Nevertheless, we also considered alternative scenarios. For example, the observed RRL emission could originate from an unresolved, very young ultracompact or hypercompact Hii region (UC/HC Hii) (e.g., Sewilo et al., 2004; Keto et al., 2008). This scenario, however, seems unlikey given the extremely broad line width of \sim$$100 km s*-1* observed for the H30$\alpha$ line. In this case, the expected static (non-expanding) Hii region thermal width of $\sim 20$ km s*-1* cannot alone explain the observed line width, neither can the typical Hii region expansion velocities of $\sim 10$–$15$ km s*-1* (e.g., Faerber et al., 2025). Moreover, although extremely dense, hypercompact Hii regions may have broader lines (e.g., Keto et al., 2008; Rivera-Soto et al., 2020), they only reach widths of $\approx 50$ km s*-1* at cm wavelengths, and $\approx 30$ km s*-1* at mm wavelengths. Another possible scenario is that the RRL emission arises from an inner ionised disc (see, e.g., the case of G17.64$+$0.16 by Maud et al., 2019, see also Tanaka et al. 2016; Jáquez-Domínguez et al. 2025). However, the tentative north-south spatial shift and velocity gradient reported in Fig. 9 for the H30$\alpha$ line, if confirmed with further studies, would favour the radio jet scenario over the ionised disc. Moreover, based on the Keplerian-like profiles shown in Fig. 10, one would expect full widths of $\sim 50$ km s*-1* for RRLs tracing the ionised disc at our current spatial resolution. This is about a factor two to six lower than the observed H30$\alpha$ widths ($\sim 100$ km s*-1*, and $\sim 300$ km s*-1* at zero intensity). All this leaves the jet scenario as the most plausible explanation for the origin of the H30$\alpha$ emission in W75N(B)-VLA3. Should this be confirmed with new observations, VLA3 could represent the first protostellar jet with detected thermal radio recombination lines.

Robust detections of thermal RRLs associated with protostellar thermal radio jets have, to date, remained elusive. This is largely due to the limited sensitivity of current centimetre-wavelength interferometers. While recombination maser lines have been reported towards some disc-outflow systems (e.g., MWC349A: Martin-Pintado et al. 1989; MonR2-IRS2: Jiménez-Serra et al. 2013; G45.47+0.05: Zhang et al. 2019) including the Cepheus A HW2 protostellar jet (Jiménez-Serra et al., 2011), general thermal RRL detections are a significant observational challenge. Detecting multiple thermal RRLs arising within the radio jets at both centimetre (e.g., with the ngVLA, and/or SKA) and millimetre (e.g., with ALMA, SMA, and/or NOEMA) wavelengths will constitute a profound breakthrough. Such detections would, for instance, enable us to analyse the morphology and 3D kinematics of protostellar ionised radio jets by examining their line-of-sight velocities and proper motions, thus providing fundamental information to understand their collimation mechanisms (see, e.g., Anglada et al., 2018; Tanaka et al., 2016; Jáquez-Domínguez et al., 2025, and references therein). Our findings identify W75N(B)-VLA3 as a unique candidate for undertaking this kind of studies.

5 Conclusions

This work presents a study of the high-mass protostar W75N(B)-VLA3 using ALMA data at high angular resolution ($\approx 0\aas@@fstack{\prime\prime}12$, corresponding to $\approx 156$ au). We have used the continuum emission at 1.3 mm, together with spectral line emission from different species (e.g., SiO, H2CO, HC3N, SO2, H2O, H30$\alpha$) to characterize the morphology and properties of the dense gas around the massive protostar. Our principal findings include:

1. –

Keplerian-like disc around VLA3: ALMA continuum and molecular line observations (specifically of SO2, H2CO, HC3N, HNCO, CH3CN, and H2O) reveal a compact, rotating molecular disc with a radius of $\sim$450 au surrounding VLA3. The disc has an estimated gas and dust mass of $\approx 0.43$–1.74 $M_{\odot}$ and its velocity profile is consistent with Keplerian motions around a protostar with mass $\approx 16$ $M_{\odot}$. This dynamical mass aligns well with the stellar mass of $\approx 17$–$19$ M⊙ derived from the correlation between radio continuum and bolometric luminosities. The low disc-to-star mass ratio (approximately one-tenth of the central object’s mass) indicates that the disc is a stable, centrifugally supported structure, rather than a gravitationally unstable toroid. Calculations of the Toomre’s Q parameter further support its stability against gravitational instabilities compared to other massive structures. The properties of the disc in VLA3 are comparable to those recently measured towards discs around massive YSOs when observed at high resolution and sensitivity. The disc is oriented perpendicular to the VLA3 thermal radio jet previously reported at cm wavelengths. 2. –

Disc asymmetry: The disc around VLA3 exhibits two intensity peaks in the dense gas molecular tracers, with the southwestern peak appearing brighter. This difference in brightness suggests an asymmetrical disc structure, potentially caused by anisotropic accretion streams, a phenomenon observed in other young stellar objects. An elongated, slightly curved structure visible to the north-east in some dense gas tracers (e.g., HC3N, H2CO) might be indicative of such an accretion streamer. While the two-peak morphology might suggest a molecular cavity around the central source, it is more plausibly interpreted as an apparent effect due to high dust opacity in the optically thick central region. The vibrationally excited H2O emission, which peaks close to the central source, may arise from the innermost, warmest regions of the disc, where dust would have sublimated, with some potential contribution from a north-south outflow. 3. –

Expanding molecular gas: Blue-shifted absorption features observed in the SiO and H2CO spectral lines, reminiscent of classical P-Cygni profiles, suggest the presence of expanding molecular gas in the near vicinity of VLA3. This expansion motions might originate from the outflow activity of the central source or an expanding photoionised region. 4. –

Detection of H30$\alpha$ emission towards the radio jet: The H30$\alpha$ recombination line is spatially coincident with the thermal radio jet, raising the possibility that W75N(B)-VLA3 could be the first known protostellar radio jet from which thermal radio recombination lines have been detected. The extremely broad H30$\alpha$ full width of $\approx 100$ km s*-1* (or $\approx 300$ km s*-1* at zero intensity), indicative of high velocity motions likely associated with outflows and jets, together with a tentative H30$\alpha$ spatial shift and velocity gradient in the jet’s direction, favour the radio jet origin as the most plausible scenario for this line. New observations at millimetre wavelengths with current instruments such as ALMA, NOEMA, and SMA, together with observations at centimetre wavelengths with ngVLA and SKA, will enable a detailed study and modelling of the RRLs emission and their origin in VLA3.

In summary, these results establish W75N(B)-VLA3 as one of the few well-documented massive protostar-disc-jet systems with a disc radius $\lesssim$500 au, confirming its significance as a key system that offers critical observational constraints for high-mass star formation theories; and identify this source as a unique radio jet object, suitable for conducting an in-depth study of its RRL emission.

Acknowledgements

We thank the anonymous referee for their constructive comments that have helped to improve the content and clarity of the manuscript. A.S-M., J.M.T, and J.M.G. acknowledge support from the PID2023-146675NB grant funded by MCIN/AEI/10.13039/501100011033, and by the programme Unidad de Excelencia María de Maeztu CEX2020-001058-M. A.S-M. also acknowledges support from the RyC2021-032892-I grant funded by MCIN/AEI/10.13039/501100011033 and by the European Union ‘Next GenerationEU’/PRTR. J.F.G., G.A., and G.A.F. acknowledge support through the Severo Ochoa grant CEX2021-001131-S and PID2023-146295NB-I00PID, funded by MCIN/AEI/10.13039/5011000110. G.A.F. also gratefully acknowledges the Deutsche Forschungsgemeinschaft (DFG) for funding through SFB 1601 ’Habitats of massive stars across cosmic time’ (sub-project B1), and the University of Cologne and its Global Faculty programme. S.C. acknowledges the support of DGAPA PAPIIT grants IN104521 and IN107324 as well as CONAHCyT grant CF-2023-I-232. C.G. acknowledges financial support by the European Union NextGenerationEU RRF M4C2 1.1 project n. 2022YAPMJH; FAPESP (Fundação de Amparo á Pesquisa do Estado de São Paulo) under grant 2021/01183-8. C.C.G. acknowledges support from UNAM DGAPA PAPIIT grant IG101224 and from CONAHCyT Ciencia de Frontera project ID 86372. A.R.K. acknowledges support from Consejo Nacional de Investigaciones Científicas y Técnicas (CONICET), and support from the European Research Executive Agency HORIZON-MSCA-2021-SE-01 Research and Innovation programme under the Marie Skłodowska-Curie grant agreement number 101086388 (LACEGAL). J.C. acknowledges support from the IG101125 (DGAPA) grant. This paper makes use of ALMA data ADS/JAO.ALMA#2019.1.00059.S. ALMA is a partnership of ESO (representing its member states), NSF (USA) and NINS (Japan), together with NRC (Canada), MOST and ASIAA (Taiwan), and KASI (Republic of Korea), in cooperation with the Republic of Chile. The Joint ALMA Observatory is operated by ESO, AUI/NRAO and NAOJ. The data processing for this paper has been carried out using the Spanish Prototype of an SRC (SPSRC) service and support, funded by the Spanish Ministry of Science, Innovation and Universities, by the Regional Government of Andalusia, by the European Regional Development Funds and by the European Union NextGenerationEU/PRTR.

Data Availability

The ALMA data used in this work are available at the ALMA Science Archive under project number ADS/JAO.ALMA#2019.1.00059.S

Appendix A Additional figures

Figures 15 to 24 show the channel maps and the $0^{\mathrm{th}}$ and $1^{\mathrm{st}}$-order moment maps of the molecular species studied in this work (i.e., SiO, H2CO, HC3N, HNCO and CH3CN). These figures complement similar figures shown in Figs. 4 and 5 for the SO2 species.