Issue |
A&A
Volume 674, June 2023
|
|
---|---|---|
Article Number | A113 | |
Number of page(s) | 50 | |
Section | Planets and planetary systems | |
DOI | https://doi.org/10.1051/0004-6361/202245425 | |
Published online | 13 June 2023 |
The Disc Miner
II. Revealing gas substructures and kinematic signatures from planet-disc interaction through line profile analysis
1
European Southern Observatory,
Karl-Schwarzschild-Str. 2,
85748
Garching bei München, Germany
2
Leiden Observatory, Leiden University,
PO Box 9513,
2300 RA
Leiden, The Netherlands
e-mail: izquierdo@strw.leidenuniv.nl
3
Dipartimento di Fisica e Astronomia “Augusto Righi” - Alma Mater Studiorum Università di Bologna,
Via Gobetti 93/2,
40129
Bologna, Italy
4
INAF - Osservatorio Astrofisico di Arcetri,
Largo E. Fermi 5,
50125
Firenze, Italy
5
Dipartimento di Fisica, Università degli Studi di Milano,
Via Celoria, 16,
Milano
20133, Italy
6
School of Physics & Astronomy, University of Leicester,
University Road,
Leicester
LE1 7RH, UK
7
Max-Planck-Institut für extraterrestrische Physik,
Gießenbachstr. 1,
85748
Garching bei München, Germany
Received:
10
November
2022
Accepted:
4
April
2023
Detecting planets in the early stages of formation is key to reconstructing the history and diversity of fully developed planetary systems. The aim of this work is to identify potential signatures from planet-disc interaction in the circumstellar discs around MWC 480, HD 163296, AS 209, IM Lup, and GM Aur, through the study of molecular lines observed as part of the ALMA large program MAPS. Extended and localised perturbations in velocity, line width, and intensity have been analysed jointly using the DISCMINER modelling framework, in three bright CO isotopologues, 12CO, 13CO, and C18O J = 2−1, to provide a comprehensive summary of the kinematic and column density substructures that planets might be actively sculpting in these discs. We find convincing evidence for the presence of four giant planets located at wide orbits in three of the discs in the sample: two around HD 163296, one in MWC 480, and one in AS 209. One of the planet candidates in HD 163296, P94, previously associated with velocity signatures detected in lower velocity resolution 12CO data, is confirmed and linked to localised velocity and line width perturbations in 13CO and C18O too. We highlight that line widths are also powerful tracers of planet-forming sites as they are sensitive to turbulent motions triggered by planet-disc interactions. In MWC 480, we identified non-axisymmetric line width enhancements around the radial separation of candidate planet-driven buoyancy spirals, which we used to narrow the location of the possible planet to an orbital radius of R = 245 au and PA = 193°, at a projected distance of 1.33″ from the star. In the disc of AS 209, we found excess 12CO line widths centred at R = 210 au, PA = 151°, at a projected distance of 1.44″, spanning around the immediate vicinity of a circumplanetary disc candidate proposed previously, which further supports its presence. We report no clear localised or extended kinematic signatures in the discs of IM Lup or GM Aur that could be associated with the presence of planets or gravitational instabilities. On the other hand, we demonstrate that pressure minima exhibit line width minima counterparts in optically thick emission, making them robust tracers of gaps in the gas surface density when analysed together with azimuthal velocity flows. Finally, we show that nine out of eleven millimetre dust continuum rings in the sample are co-located with pressure bumps traced by kinematic modulations, indicating that aerodynamic confinement via pressure traps is a common mechanism for the formation of dust substructures in these discs. Overall, our analysis reveals that all discs in the sample present a remarkable level of substructure in all the line profile observables considered, regardless of the CO isotopologue. However, the magnitude and morphology of the substructures vary between discs and tracers, indicating that the kinematics and thermodynamic properties are likely shaped by different physical mechanisms in each object. We propose that the main kinematic signatures identified in the discs of MWC 480, HD 163296, and AS 209 have a planetary origin, although they do not always manifest as highly localised perturbations, while the discs of IM Lup and GM Aur do not yield clear signatures pointing to the presence of massive planets. Our simultaneous analysis of multiple tracers and observables aims to lay the groundwork for robust studies of molecular line properties focused on the search for young planets in discs.
Key words: planets and satellites: detection / planet-disk interactions / protoplanetary disks
© The Authors 2023
Open Access article, published by EDP Sciences, under the terms of the Creative Commons Attribution License (https://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.
This article is published in open access under the Subscribe to Open model. Subscribe to A&A to support open access publication.
1 Introduction
The formidable precision and sensitivity provided by the Atacama Large Millimeter/submillimeter Array (ALMA) have revolutionised the way we study and understand protoplanetary discs. High angular resolution observations of these planet-forming environments have revealed a myriad of dust substructures possibly shaped by the hydrodynamic interaction of embedded planets with the disc itself, but also due to instabilities in the gas to which the dust grains respond according to their composition and morphology (see Testi et al. 2014; Andrews 2020; Bae et al. 2022a, and references therein).
ALMA has also opened a unique window into the bulk of the gas disc through deep and high-resolution observations of molecular lines, which provide a vast amount of information that has only recently started to be exploited in the context of planetary formation (see reviews by Dutrey et al. 2014; Miotello et al. 2022; Pinte et al. 2022). For instance, a simultaneous analysis of multiple molecular lines gives access to the three-dimensional structure of discs as different chemical species are arranged differently across the radial and vertical extent of the disc (van Zadelhoff et al. 2001; Walsh et al. 2010; van Dishoeck & Bergin 2020; Öberg & Bergin 2021; Paneque-Carreño et al. 2023). Along with thermochemical modelling, this can be used as input information for detailed studies of the formation and destruction sites of molecular species (Booth et al. 2021; Leemker et al. 2021; Zhang et al. 2021), which in addition allows one to gain knowledge on internal and external thermodynamic processes that are key to predict the evolution of discs and planets within (Benz et al. 2014). In parallel, by characterising the morphology and Doppler shift of these lines, it is possible to directly probe the kinematics and temperature of the disc, again in a three-dimensional fashion, but also to perform indirect estimates of other intrinsic properties such as the gas surface density and pressure (as in e.g. Rosenfeld et al. 2013; Teague et al. 2018a,b; Dullemond et al. 2020; Yu et al. 2021), quantify non-thermal motions (as in e.g. Teague et al. 2016; Flaherty et al. 2017, 2020; Liu et al. 2018), or even derive stellar (Czekala et al. 2015, 2017) and disc masses (Veronesi et al. 2021; Lodato et al. 20226) through analyses of the disc background velocity field.
Delving into smaller scales, the study of molecular lines is presently the workhorse of kinematic studies focused on the detection of velocity perturbations to the Keplerian rotation of the gas in planet-forming discs, and with good reason. As several studies have theorised, these kinematic features represent a window into a wide variety of physical processes, including the response of the disc to the presence of embedded planets (Perez et al. 2015; Pérez et al. 2018; Dong et al. 2019; Bae et al. 2021; Bollati et al. 2021; Izquierdo et al. 2021; Rabago & Zhu 2021), stellar companions (Facchini et al. 2018a; Rosotti et al. 2020a) and flybys (see Cuello et al. 2023, for a review), but also to several (magneto-)hydrodynamic instabilities (see Pinte et al. 2022; Bae et al. 2022a, for a summary). In the planetary scenario, these studies predict that an embedded planet can simultaneously trigger localised and extended perturbations, not only in velocity, but also in density and temperature, which translate into characteristic features observable in molecular line emission. Indeed, a select number of works have claimed the discovery of planet-driven signatures based on the detection of ‘kinks’ in intensity channel maps of molecular lines (Pinte et al. 2018b, 2019, 2020; Calcino et al. 2022; Verrios et al. 2022), but also on the analysis of small- and large-scale velocity perturbations to the Keplerian rotation of the gas disc (as in Casassus & Pérez 2019; Teague et al. 2018a, 2021, 2022; Izquierdo et al. 2022).
Recently, the Molecules with ALMA at Planet-forming Scales (MAPS) survey mapped the spatial and spectral distribution of more than 20 chemical species in five protoplanetary discs around the young stars MWC 480, HD 163296, AS 209, IM Lup, and GM Aur (Öberg et al. 2021). Several studies of these data have revealed that the molecular line emission from simple and complex species in these discs is rich in substructure, both radially and vertically (Law et al. 2021a,b; Zhang et al. 2021; Guzmán et al. 2021; Ilee et al. 2021; Le Gal et al. 2021; Paneque-Carreño et al. 2023), providing unparalleled clues about the chemical composition of planets that will eventually configure a planetary system such as ours. More specifically, Teague et al. (2021) present detailed kinematic analyses of the discs around MWC 480 and HD 163296, identifying numerous velocity signatures in 12CO, and tentatively also in 13CO and C18O, possibly connected to the hydrodynamic and gravitational interaction of these discs with embedded planets. Recently, the disc of AS 209 has gained additional interest in the context of the detection of planet-driven features as Bae et al. (2022b) propose the presence of a circumplanetary disc (CPD) after observing a localised signal in the intensity of 13CO channels, within a gap on the outskirts of the circumstellar disc, at a radial distance of ~200 au from the star. Nevertheless, given the amount and wealth of the data provided by the MAPS survey, much remains to be investigated, especially in light of new modelling techniques that allow not only the kinematics of discs to be mined, but also other properties encoded in molecular line profiles with greater precision than previous methods.
In this article, we employ the line profile analysis techniques developed in Izquierdo et al. (2021, or Paper I hereafter), implemented in the DISCMINER package, to search for gas substructures and kinematic signatures that may point to, or refute, the presence of embedded planets in the discs observed by the ALMA large program MAPS. The paper is organised as follows. Section 2 gives an overview of the molecular lines and sources studied. Section 2.1 presents a summary of the modelling strategy, as well as the line profile observables considered in Sects. 3 and 4 to characterise kinematic and gas substructures in the discs. In Sect. 5, the interpretation of these features and their possible relationship with embedded planets is discussed, and we conclude with the main findings of the paper in Sect. 6.
2 Observations and models
In this work we make use of publicly available data from molecular line observations of five protoplanetary discs around the young stars MWC 480, HD 163296, AS 209, IM Lup, and GM Aur, as part of the ALMA large program MAPS (Öberg et al. 2021). We use the JvM-corrected images with a robust parameter of 0.5 (Czekala et al. 2021)1. Although there are multiple molecules and transitions at hand within the MAPS survey, our analysis focuses on the three brightest CO isotopologues, 12CO, 13CO and C18O, in the J = 2−1 rotational transition, observed at an angular resolution of 0.15″, and a velocity channel spacing of 0.2 km s−1. The unprecedented angular and spectral resolution together with the high sensitivity of these observations allow us to apply dedicated state-of-the-art techniques to study the kinematics and morphology of these objects with exquisite precision. Cartoons illustrating the orientation and rotation direction of the discs analysed in this work are presented in Fig. 1.
2.1 Discminer models
The majority of the analysis presented in this paper is performed with the DISCMINER channel-map modelling tool, whose main capabilities were introduced in Paper I and first applied to actual data in Izquierdo et al. (2022) with the aim of looking for planetary signatures in the kinematics of the disc around HD 163296 using 12CO J = 2−1 line observations. In this chapter, we describe the parameters and attributes that make up the channel-map models of the five discs of interest, and present a comparison between selected model and data channels.
2.1.1 Model parameters and attributes
DISCMINER2 aims to reproduce the channel-by-channel intensity of molecular line emission from discs by modelling various types of parameters which can be classified in four groups:
Orientation parameters that control the projected appearance of the disc in the sky, such as inclination, i, position angle, PA, and offset, xc and yc.
Velocity parameters that shape the rotation velocity of the gas disc, such as the stellar mass, M*, and systemic velocity, υsys. Other parameters such as the gas mass or the steepness of the radial pressure profile can also play an important role in setting the disc background velocity (see e.g. Rosenfeld et al. 2013). Although accountable by DISCMINER, modelling these parameters is computationally expensive and out of the scope of this paper. Thus, the model rotation velocity is assumed Keplerian for all discs, with vertical differential rotation sustained by balance of forces in the z-direction.
Surface parameters that determine the elevation of upper and lower emission surfaces above and below the disc midplane.
Line profile parameters that account for the shape of the line profiles, controlling the variation of line peak intensity, line width and line slope across the disc radial and vertical extent.
These parameters are then assembled together in a two-component line profile kernel (one per emission surface), or in a single-component profile (if the lower surface is not resolved), to compute the model intensity, Im, as a function of the disc cylindrical coordinates, (R, z), and the velocity channel, υch. To combine the contribution of the upper and lower emission surfaces to the total line profile, we take the highest intensity between the two components on each pixel and velocity channel. For each component, we adopt the same kernel of Paper I, a generalised bell function of the form, (1)
where Ip is the line peak intensity, and is the Keplerian velocity projected along the line of sight. The line width, Lw, is the half width of the profile at half power. The dimensionless attribute, the line slope Ls, controls how steep the signal drops at the wings and in turn it also determines the spectral extent of the plateau at the top of the profile. As explained in Paper I, this kernel performs better than a Gaussian at reproducing optically thick lines which are flat at the top and decay rapidly towards the wings, as is usually the case for 12CO and 13CO lines (Beckwith & Sargent 1993), but at the same time it is flexible enough to account for optically thinner lines if needed. A summary of the assumed functional forms and parameters involved in each of the model attributes is presented in Table 1.
Finally, to find best-fit model channel maps to the data we couple DISCMINER with a Markov chain Monte Carlo (MCMC) random sampler, EMCEE (Foreman-Mackey et al. 2013), which walks over a vast range of parameters to determine a subset of them that best reproduces the projected intensity of the input datacube on each velocity channel. The best-fit model parameters obtained for each disc and tracer are summarised in Tables A.1–A.5. To explore the influence of disc inclination and different tracers on the retrieved dynamical mass, we also made runs assuming disc inclination fixed at the continuum value for all sources and tracers, and reported the obtained parameters in the same tables. We note that differences in M* range between ~5–10% for variations in inclination within ~1–4°, which is a consequence of the strong correlation between these two parameters (see e.g. Fig. 11 of Izquierdo et al. 2022).
Fig. 1 Illustrating the orientation and rotation direction of the discs observed by the ALMA large program MAPS and analysed here. Dashed arrows mark the origin of each disc coordinate system. For convenience, we assume this reference axis to be parallel to the projected semi-major axis that is closest to the north celestial axis in counterclockwise rotation. Also shown are the signatures reported in this work as potential planet-disc interaction features in the form of both coherent and turbulent velocity fluctuations. |
Attributes adopted in the DISCMINER models for this work.
Fig. 2 Selected CO isotopologue intensity channel maps from the disc around MWC 480, along with best-fit channels obtained with DISCMINER and their corresponding residuals. For comparison, isovelocity contours from the model upper and lower surfaces are overlaid on the data channels as solid and dashed lines, respectively. The ellipses in the background represent millimetre dust rings (solid) and gaps (dashed) as reported by Law et al. (2021a), projected on the upper emission surface of the corresponding model. Kink-like features apparent to the eye are marked with arrows. The radial location of the most prominent kink in 12CO is highlighted by the dotted purple ellipse. The beam size is shown in the bottom left corner of the top left panel for each isotopologue. |
2.1.2 Data versus model channel maps
In Fig. 2, we present a comparison of selected intensity channel maps of the three CO isotopologues studied in this work, for the disc of MWC 480, their corresponding model channels computed with DISCMINER, and residuals obtained by subtracting the intensity of each model channel to that of the data. Data and model channel maps for the discs of HD 163296, AS 209, IM Lup and GM Aur are presented in Figs. B.1–B.4.
We note that there are multiple features that can be easily grasped after quick visual inspection of the data channels. First, it is clear that the five discs yield different radial and vertical extents. Second, non-negligible residuals in channel intensities suggest that a smooth Keplerian model is not the perfect match for these data. However, it is precisely through these deviations from Keplerian that the underlying physics acting on the discs can be uncovered. A massive planet, for example, is expected to trigger perturbations to the Keplerian rotation of the gas in its vicinity, and at the same time launch spiral density waves, with characteristic kinematic and intensity structures, that propagate over large azimuthal and radial portions of the disc (Pérez et al. 2018; Bae et al. 2021; Rabago & Zhu 2021). Indeed, many of the intensity channels of all discs in this sample display wiggles, or kinks, which are usually linked, but not limited, to the presence of velocity perturbations in the gas disc. In Sects. 3 and 4, we characterise variations in the line profile observables to study the origin of these fluctuations at first sight evident in channel maps, but also of those that remain hidden to the eye.
We assume axisymmetric intensity fields for all model channels except for the 12CO disc of AS 209. The reason is that the 12CO J = 2−1 line intensity of this disc is on average twice as weak on (almost all of) the blueshifted side than on the red-shifted side due to foreground cloud absorption (Öberg et al. 2011). To take this effect into account, we used a different intensity normalisation factor for each side of the disc: 0.65I0 for −80° < ϕ < 80°, and 1.25I0 otherwise, where I0 is the intensity normalisation reported in Table A.3 for 12CO, and ϕ is the azimuthal coordinate of the disc reference frame.
2.2 Line profile observables to compare data and models
We employ three observables from line intensity profiles, (i) peak intensity, (ii) line width, and (iii) centroid velocity, to search for gas surface density, temperature and velocity substructures in all discs and tracers considered in this work. The three observables are extracted from Gaussian kernels fitted to both the data and model line profiles, and compared against each other to produce residual maps, as illustrated in Fig. 3 for the 12CO disc around MWC 480. Using information of the disc vertical structure and orientation retrieved by DISCMINER, these line profile observables and their residuals can then be deprojected onto the disc reference frame for quantitative studies of the physical location, extent, and magnitude of the retrieved signatures in the gas disc. Deprojected views of the observable residuals for the same disc and tracer are presented in Fig. 4. Radial profiles of each observable extracted from the data alone, as function of disc and isotopologue, including rotation curves and elevation of model emission surfaces, are summarised in Fig. B.5.
In practice, each line profile observable carries information about different thermodynamic and kinematic properties of the disc. (i) The line peak intensity can be used for studies of the local temperature and density of a molecular species (as in e.g. Facchini et al. 2021; Zhang et al. 2021). Likewise, (ii) the line width depends on the gas temperature and density, but it too responds to local turbulent motions (see e.g. Hacar et al. 2016; Flaherty et al. 2020). Finally, (iii) the shift in velocity of the line profile centroid with respect to the systemic velocity, or simply centroid velocity, accounts for organised gas motions driven by different, usually coupled, physical mechanisms such as the gravitational and hydrodynamic interaction between the disc, star(s) and planet(s), as well as a range of (magneto−) hydrodynamic instabilities that may develop in the gas disc under appropriate conditions (see Pinte et al. 2022, and references therein for a summary). In Fig. B.6, we outline the general recipe for computing channel-map models and extracting observables from a given data cube with DISCMINER, and present a summary of the physical properties of discs and planets that can be learned based on the analysis of such observables.
2.2.1 Upper and lower emission surfaces
One of the strengths of DISCMINER is its ability to model line profile properties on the upper and lower emission surfaces of the disc, simultaneously. As discussed extensively in Paper I, the impact of the lower emission surface can be critical when it comes to kinematical analyses of high resolution observations of circumstellar discs (see also Pinte et al. 2022, and references therein). The lower surface can systematically shift the centroid of the observed intensity profile in an uneven fashion as a function of the disc coordinates, affecting the velocities derived via first moment maps or via parametric fits to the line profile (see e.g. Fig. A.2 of Paper I). In practice, this unphysical shift is source-specific as it depends on the aspect ratio and intensity contrast between the two emitting surfaces.
To overcome this, at the cost of velocity precision, some methods measure velocities on, or around, the peak of the line profile to account primarily for the centroid velocity of the disc upper surface (see e.g. Teague & Foreman-Mackey 2018; Teague et al. 2021). However, when the emission is optically thin, or just marginally optically thick, these methods struggle at distinguishing between upper and lower surfaces because the intensity contrast between both is generally small in these scenarios. The inclusion of the lower surface contribution to the line profiles of our model aims to partly alleviate this effect and enables us to consider the bulk of velocity channels to measure the line centroid shift in velocity with little loss of precision.
2.2.2 Decomposition method to extract azimuthal and meridional velocities
In this section, we present a new pipeline to quantify the vertical and azimuthal components of the disc velocity field, υz and υϕ, by decomposing two types of deprojected velocity residual maps; the standard and the absolute residual map. While the former consists of a direct subtraction of the model centroid velocities to those retrieved from the data, υ0,d – υ0,m (as in e.g. Fig. 5), the latter is the difference between the absolute values referred to the model systemic velocity, |υ0,d – υsys,m| – |υ0,m – υsys,m|. Taking one or another allows to highlight variations in the three-dimensional components of the velocity in different ways. In a standard velocity residual map, an axisymmetric azimuthal perturbation, δυϕ, would appear as two joint semi-rings with the same magnitude but different signs split by the disc minor axis (see e.g. Fig. 5 of Teague et al. 2019b). In the absolute velocity residuals, the same perturbation would look as a full ring, namely one where there is no sign flip involved around the disc minor axis. The reciprocal occurs for an axisymmetric vertical perturbation; a ring-like pattern would be observed in the standard residuals, and two semi-rings with different sign in the absolute residuals.
Now, assuming that the three-dimensional components of the velocity field are axisymmetric, one can obtain the azimuthal and vertical components on each radius by computing azimuthal averages of these absolute and standard velocity residuals, respectively. In the first case, the vertical and radial velocity components cancel out over 2π averages, as illustrated below. In the second scenario, it is the azimuthal and radial components that vanish. Radial profiles of the azimuthal and vertical velocity components computed with this pipeline for all discs and tracers are presented in Sect. 4.
Fig. 3 Extraction of the line profile observables considered in this work, using MWC 480 12CO data as example. Peak intensities (top), line widths (middle), and centroid velocities (bottom) are displayed for both data and best-fit model, as well as residual maps showing the difference between them. Sketches on the rightmost column illustrate how each line profile observable is computed on the data cube and subsequently compared to those derived on the smooth, Keplerian model cube obtained by DISCMINER. The black solid lines are reference contours to ease comparison. Their levels are indicated in the colour bars. Residuals on each of the observables are driven by temperature, density and velocity fluctuations in the gas disc. |
Fig. 4 Line width (top), velocity (middle), and peak intensity (bottom) residuals computed for the 12CO disc of MWC 480, as in Fig. 3, but deprojected on the disc reference frame. The vertical axis in this frame is parallel to the projected minor axis on the sky. The solid lines mark the radial location of millimetre dust rings, and the purple dotted line that of the most prominent kink identified by eye in 12CO channels. Also shown is the north-sky (or PA = 0°) axis for reference. The oblique shape of this axis is due to the deprojected geometry of the disc emission surface. |
Computing υz
A standard velocity residual map is the direct difference between the observed and modelled line-of-sight velocities, υ0,d – υ0,m, which can be decomposed into three orthogonal components of the velocity field as follows,
where the subscripts d and m stand for data and model, respectively. Therefore, an azimuthally averaged profile of this map, which is just a sum of integrals applied to each of the above terms, computed on a given annulus of the disc can be written as,
which naturally results in an expression to compute the vertical component of the velocity field, υz,d, thanks to the azimuthal symmetry of the other two components. We note that since the model velocity field is fully azimuthal it is actually irrelevant for this derivation, and thus an azimuthal average of υ0,d – υsys,m alone would result in the same outcome. Either way, the model-dependent parameters involved in this calculation include systemic velocity, disc inclination, and of course the disc vertical structure which in practice determines the topology of the projected annular paths where azimuthal averages are measured.
Finally, it is straightforward to show that if the azimuthal section ψ, over which azimuthal averages are taken, is symmetric around the disc major and minor axes the same equality holds, (4)
Computing υϕ
If the same procedure is applied on absolute velocity residuals instead, following Appendix C, we obtain that the azimuthal component of the velocity perturbation field, Δυϕ, can be written as, (5)
whose generalisation does depend this time on the total extent of the azimuthal section ψ considered, (6)
These general forms of υz and υϕ are useful in cases where either physical or unphysical contaminating velocities present overwide azimuthal sections of the disc need to be masked out of the analysis in order to remove, or at least minimise, their influence on the computation of velocity profiles, as will be clear later in Sect. 3.2 for the 12CO discs of IM Lup and GM Aur.
As a final note, we highlight that velocity perturbations to the Keplerian rotation are not necessarily axisymmetric as originally imposed by this type of decomposition methods (see also e.g. Teague et al. 2019a). For instance, the presence of spiral features is an immediate proof that this assumption can break. Using a representative sample of (transition) discs, Wölfer et al. (2023) showed that non-axisymmetric substructures may be relatively common features in both intensity and velocity fields. However, despite the fact that such features would impact the amplitude of the velocity components retrieved by this method, other less-affected observable quantities such as the radial location of velocity peaks, troughs, and gradients, can still provide invaluable clues about the presence and direction of azimuthal and vertical flows of material, generally interpreted as due to variations in the disc pressure structure, as explained and illustrated in Sect. 4.
Fig. 5 Centroid velocity residuals traced by 12CO emission for all discs, referred to each disc reference frame. The annotated solid lines mark the radial location of millimetre dust rings (denoted as Bxxx), and the purple dotted line that of the most prominent kink identified by eye in 12CO intensity channels (denoted as Kxxx). Also shown is the north-sky (or PA = 0°) axis for reference. The oblique shape of this axis is due to the deprojected geometry of each disc vertical structure retrieved by the best-fit models. The circle in the bottom left corner of each panel indicates the location of the blue- and redshifted side for each disc. These maps illustrate the deviations from Keplerian rotation of the gas component located at high elevations (0.2 < z/r < 0.4) over the midplane. |
3 Kinematic signatures
In this Chapter, we present a gallery of the kinematic substructures identified in the discs of MWC 480, HD 163296, AS 209, IM Lup, and GM Aur, traced by 12CO, 13CO, and C18O J = 2−1 line emission. More specifically, we report the location and magnitude of significantly localised perturbations, possibly driven by unseen massive planets, and comment on the plausible link between the presence of large-scale kinematic signatures and planet candidates proposed in recent literature. In Table 2, we summarise the relative location of the signatures attributed to the presence of embedded planets reported here and in other works.
As was already manifest in the channel maps presented in Figs. 2 and B.1–B.4, the five discs in the sample display intensity fluctuations at several places, especially in 12CO emission. Qualitatively, it is easy to note that in the disc of MWC 480, for example, there are intensity kinks that span over multiple channels at different radial and azimuthal locations (see arrows in Fig. 2). In the disc of HD 163296, there seems to be fewer, yet very prominent kinks as originally reported by Pinte et al. (2018b). Likewise, the disc of AS 209 appears highly perturbed and void around an annular region at a seemingly constant radial distance of ~200 au. At the same radius, a rather striking kink-like feature in the south of the sky is evident in channels near the systemic velocity, but it is apparent too that such a signature extends over almost the entire azimuth of the disc3. Nevertheless, despite being illustrative and guiding, the nature of these apparent signatures can only be revealed by quantitative analyses that provide a precise estimate of the location, magnitude, and extent of the underlying perturbations in the disc velocity field.
To this end, in the following sections we focus on the study of centroid velocity and line width residuals computed according to Sect. 2.2, and illustrated in Figs. 5, B.7, and B.8, for 12CO, 13CO, and C18O velocities, for all discs, and in Figs. B.9, B.10, and B.11 for line widths. We note that some of the prominent signatures identified by eye in intensity channels are now explicit and quantifiable in these residual maps. For instance, the seemingly extended kinks observed in 12CO channels of the disc around MWC 480 render as ring-like structures in centroid velocity residuals, revealing the presence of radially localised deviations from Keplerian rotation in this system. Also, it is now clear that the strong kink observed in the south of the AS 209 disc appears to be the consequence of a large-scale, filament-like kinematic signature that spans over the entire azimuth of the disc, rather than due to a localised feature. In the following sections we assess how localised these velocity perturbations actually are across all discs and tracers in a quantitative manner.
Location of potential planet-driven signatures in discs reported in this and in previous works based on gas velocity perturbations and substructures, assuming that the planet candidate is in the disc midplane.
3.1 Localised velocity and line width perturbations
Strong and localised velocity perturbations are expected in the vicinity of massive planets, that is, around and along spiral wakes triggered by the gravity of the planet and its hydrodynamic interaction with the host disc (see e.g. Pérez et al. 2018; Rabago & Zhu 2021). In Paper I, we demonstrated that this fact can be readily exploited to pinpoint the azimuthal and radial location of embedded planets with remarkable precision. Our detection technique uses a clustering algorithm to spot those areas of the disc where the variance of the velocity (or line width) perturbations is significantly greater than that of the background velocity (or line width) fluctuations. To do so, we first split the disc radial extent into annuli with a width of a quarter of the beam size. Next, in each annulus, we perform an azimuthal scan to extract the peak residual and take record of its magnitude and azimuthal location. Finally, using a K-means algorithm, we group these peak residuals into clusters along both the radial and azimuthal coordinates of the disc. A cluster is classified as significant when the spectral variance of the peak residuals within the cluster is above a threshold of 3σ set by the variances of the other clusters. The region of the disc where the significant clusters intersect in azimuth and radius is attributed to a strong, localised, and coherent perturbation to the gas Keplerian velocity, possibly related to the presence of an unseen planet. In a quest for new planet candidates, we apply this same methodology to the velocity and line width residual maps of the five discs and CO isotopologues studied in this work.
3.1.1 Detections in the disc of HD 163296: Velocity
In Izquierdo et al. (2022), we applied this technique on DSHARP data (Andrews et al. 2018) of the 12CO disc around HD 163296, finding two significantly localised perturbations possibly driven by massive planets embedded in the disc at wide orbits: one at 94 au, referred to as P94, and another at 261 au, or P261. The latter is potentially related to the kink-like feature reported by Pinte et al. (2018b), attributed in the same work to a planet twice the mass of Jupiter. Here, very similar results are obtained if the same analysis is carried out on the MAPS observations of this disc in the three CO isotopologues of interest.
Figure 6 illustrates the folded velocity maps computed from 12CO, 13CO, and C18O centroid velocities. In the same figure, we present the results of applying the aforementioned clustering algorithm to these maps, as well as azimuthal deprojections of the magnitude of the folded velocity residuals. Our algorithm detects two significantly localised velocity perturbations in different isotopologues, one in 12CO and another in 13CO. The 12CO perturbation, detected at an orbital radius of R = 93 au and an azimuth of ϕ = 51° in the disc reference frame, is in excellent agreement with the P94 perturbation found by Izquierdo et al. (2022). The P94 convention will thus be kept throughout the work. The 13CO perturbation, or P81 hereafter, is detected at R = 81 au and ϕ = 70°, namely interior to the D86 dust gap, and centred at an azimuth closer to the disc minor axis than that of the P94 perturbation. This is illustrated in Fig. 7, where we overlay contours of the velocity perturbations detected in 12CO and 13CO, along with more extended kinematic residuals for reference. Azimuthally averaged line widths are coloured in the background to highlight the location of gas gaps traced by line width minima (see Sect. 4). Even though the centres of the 12CO and 13CO perturbations appear distant from each other, we note that both the blue- and the redshifted parts of these signatures overlap in their azimuthal extent, suggesting that they are likely connected to a common origin.
While the 12CO perturbation yields a folded magnitude of ~0.4 km s−1, the 13CO perturbation peaks at ~0.3 km s−1 velocities. In C18O, the clustering algorithm does not detect any significantly localised signature. However, it is not less interesting as it also displays a Doppler flip pattern similar to that of 13CO, at the same location. The folded magnitude of this perturbation is weaker, around ~0.2 km s−1. Taken together, it is possible that we are witnessing for the first time the structure of a planet-driven perturbation traced at different scale heights in a disc. The varying magnitude of the tentative and the two detected perturbations could be indicating that different velocity components are being probed as one goes from layers in the disc atmosphere to regions closer to the midplane. The morphology of the Doppler shift pattern suggests that the radial component of the perturbation is dominant on layers closer to the disc midplane.
3.1.2 Detections in the disc of HD 163296: Line width
In addition to producing organised velocity perturbations to the Keplerian rotation of the gas disc, massive embedded planets can also trigger turbulent motions in the neighbouring fluid (Disk Dynamics Collaboration 2020; Pinte et al. 2022). These non-thermal motions translate into broadening of the line intensity profiles on top of the local thermal broadening, and are expected to be observable if the driving planet is massive enough (Dong et al. 2019). In this section, we report the detection of localised line width enhancements in the disc of HD 163296. Line width perturbations identified as elongated substructures in the other discs of the sample are discussed in Sect. 3.2.
In Fig. 8, we present folded line width residuals for the disc of HD 163296 traced by 12CO, 13CO and C18O. In all cases, we identify positive non-axisymmetric line width fluctuations of the order of 0.1 km s−1, near the orbit of the planet candidate associated with the P94 velocity perturbation found in 12CO. These line width signatures extend for about 120° between the radial location of the dust ring, B67, and the dust gap, D85. Moreover, our clustering algorithm detects two significantly localised perturbations in the line widths of 13CO and C18O, both associated with strong positive residuals. The centre of the 13CO perturbation is detected at R = 79 au, ϕ = 74°, and that of the C18O perturbation is detected at R = 86 au, ϕ = 75°. We note that the location of these perturbations broadly coincides with the centre of the localised velocity perturbation, P81, detected in the kinematics of 13CO as reported in Sect. 3.1. These localised line widths further reinforce the hypothesis of the presence of a giant planet associated with the P94 velocity perturbation, at the radial and azimuthal location proposed by Izquierdo et al. (2022), and represent the first robust detection of localised line broadening features potentially driven by a massive planet embedded in a disc. Although not significantly localised, we also identify 12CO line width enhancements at larger radial distances, near the planet candidate P261 associated with the kink-like feature, K260, which include a strong ~0.2 km s−1 signature between R = 160–260 au and 45° < ϕ < 90°, and weaker ~0.05km s−1 fluctuations around the proposed planet orbit around R = 260 au, spanning over a much larger azimuthal extent.
3.1.3 Non-detections
Figure 9 presents folded velocity maps for the other discs where our clustering algorithm does not detect any localised perturbation simultaneously in the radial and azimuthal coordinate, in any of the CO isotopologues. We remind the reader that these folded velocity residuals highlight asymmetric substructures in the azimuthal component of the velocity perturbations on one side of the disc with respect to the other, split by the disc minor axis as projected on the sky. These folded residual maps also highlight the presence of radial and/or vertical velocity fields regardless of their symmetry around the disc minor axis.
As illustrated in the second and fourth rows of Fig. 9, the asymmetric, or non-azimuthal, kinematic features in these discs appear to span across large spatial extents. In the discs of MWC 480 and AS 209, these velocity perturbations display as coherent filamentary structures, localised in the radial direction, and with azimuthal extents that comprise almost the entire angular coordinate. Similarly, for the disc of IM Lup we report asymmetries in the form of extended filamentary structures, although weaker in magnitude than in the previous cases. No significantly localised perturbations associated with the presence of massive planets are found in this disc. Finally, in the disc of GM Aur, velocity perturbations with an elongated morphology are also found. In 12CO, a strong asymmetry can be seen in the northern part of the system due to the interaction of the uppermost layers of the disc with material possibly infalling from a remnant envelope. These signatures, classified as extended kinematic substructures, will be discussed further in Sect. 3.2 for each disc.
Fig. 6 Analysis of localised velocity perturbations in the gas disc of HD 163296. Top row: folded velocity residuals computed for 12CO, 13CO, and C18O J = 2−1 lines. Middle row: azimuthal and radial extent of clusters identified by our detection algorithm. Regions in bold yellow highlight the accepted clusters of peak residuals whose variances are significantly greater than those in the background clusters. The weighted average of the 2D location of peak residuals within azimuthal and radial clusters determines the centre of the localised velocity perturbations reported in the text. These are marked as green crosses for 12CO and 13CO. In C18O, only radial clusters were identified as localised by the technique and, therefore, no perturbation was tagged as detected in this tracer. Bottom row: azimuthal profiles of folded velocities (i.e. those in the top row) to highlight the blue- redshifted, or Doppler flip, morphology of the perturbations detected in 12CO and 13CO, and of the tentative perturbation in C18O. |
Fig. 7 Summary of the localised velocity perturbations detected by our clustering algorithm in the disc of HD 163296 in 12CO and 13CO velocities, using DSHARP and MAPS data. Localised red and blue contours highlight the morphology of the perturbations detected in 12CO, P94 and P261, whereas green and purple contours outline the localised signature detected in 13CO, P81, closely overlapping with the P94 signal. The green crosses mark the inferred location of the planet candidates linked to the aforementioned perturbations. Deep blue colours in the background highlight the location of 12CO line width minima which are co-spatial with positive gradients of azimuthal velocity flows probing gas surface density gaps. |
Fig. 8 Folded line width residuals to highlight asymmetric line width enhancements near the orbit of the planet candidate P94, as well as the location of localised line width perturbations detected by our clustering algorithm in 13CO and C18O. |
Fig. 9 Folded velocity residuals for the discs where no localised velocity perturbations are identified by our clustering algorithm. |
3.2 Extended velocity and line width perturbations
Embedded planets can also trigger extended velocity perturbations in the gas disc in the form of axisymmetric and non-axisymmetric substructures. Axisymmetric signatures driven by planets are primarily related to pressure modulations induced by these objects in the host disc. If massive enough, the gravitational interaction of a planet exerts a torque onto the disc fluid around its vicinity, which would in turn carve a density gap whose width and depth depend on the planet mass and on the local viscosity of the fluid. These pressure variations translate into ring-like deviations from Keplerian rotation in the azimuthal component of the disc velocity field (Kanagawa et al. 2015; Teague et al. 2018a). It has also been demonstrated that gas gaps can lead to meridional circulation of material from the atmosphere into the disc midplane and vice versa, which is observable through inspection of the vertical component of the velocity field with respect to the disc reference frame (Morbidelli et al. 2014; Teague et al. 2019a; Yu et al. 2021). The location and magnitude of these axisymmetric velocity perturbations linked to the presence of gas substructures are summarised in Sect. 4.
Non-axisymmetric planetary signatures in the kinematics are mainly associated with spiral-like perturbations of two types. Spiral waves triggered by the gravitational potential of the planet at Lindblad resonances, or simply Lindblad spirals, and spirals driven at buoyancy resonances, or buoyancy spirals, which develop when the vertical temperature gradient of the gas disc is positive and its adiabatic index is larger than one, so long as the buoyancy frequency of the fluid at the planet location matches with the planet orbital frequency (Zhu et al. 2012; Lubow & Zhu 2014). As demonstrated by Bae et al. (2021), buoyancy spirals are more prone to develop in the outer disc at relatively high elevations over the midplane where the thermal relaxation times are low due to infrequent gas-dust collisions. These two substructures differ in that Lindblad spirals carry mainly azimuthal and radial velocity perturbations, while vertical motions predominate in buoyancy spirals. Bae et al. (2021) also demonstrate that the pitch angle of buoyancy spirals is generally smaller than that of Lindblad spirals, especially near the radial location of the perturbing planet. Another distinguishing feature between the two signatures is that the magnitude of the velocity perturbations in Lindblad spirals decreases with increasing height over the disc midplane, while it increases with increasing height in buoyancy spirals.
MWC 480
In Fig. 10, we show polar deprojections of velocity and line width residuals obtained in 12CO and 13CO for the disc around MWC 480. As illustrated in the top panel, the 12CO disc exhibits strong and azimuthally extended ~50 m s−1 velocity perturbations that seem to be concentric and confined to narrow ~30–50 au radial sections. As we subsequently demonstrate in Sect. 4 through the analysis of azimuthal averages of velocity residuals, the velocity perturbations spanning between the dust ring B165 and R ~ 300 au, enclosed by rectangles in Fig. 10, represent contiguous anti-parallel vertical flows of material in the atmosphere of the MWC 480 disc. One of these meridional flows points upwards with respect to the disc midplane, at R = 193 au, and the other downwards, at R = 245 au. The azimuthal coherence of these substructures suggests that vertical motions must indeed be dominant in this region of the disc. However, the fact that they are not fully axisymmetric indicates that radial and/or azimuthal perturbations to the gas velocities are also present.
Another interesting feature worthy of highlighting in this region of the disc is that the strong vertical perturbations observed in 12CO are either vanished or loosely evident in 13CO, as illustrated in the second top panel of Fig. 10. The large difference between the heights traced by 12CO and 13CO may be the reason of this modulation. While 12CO is located at a z/R ~ 0.2, both 13CO and C18O barely reach z/R ~ 0.1 scale-heights (Fig. B.5, bottom row). At these altitudes, it is foreseeable that other velocity components may become dominant, making the magnitude and patterns of the kinematic signatures look different to those in 12CO. This result supports the hypothesis of Teague et al. (2021), who proposed that the meridional velocities observed in 12CO alongside the tightly wound morphology of peak intensity substructures could be the result of buoyancy spirals driven by a planet at R = 245 au in this disc. As mentioned earlier, Bae et al. (2021) demonstrated that the magnitude of the vertical motions across the vertical extent of these type of spirals is expected to decrease as the disc midplane is approached, which is in line with our finding.
Furthermore, as illustrated in the bottom panels of Fig. 10, we detect line width enhancements of the order of ~20–50 ms−1, close to the interfaces between the aforementioned vertical flows, between ~175–325 au. These enhancements appear more strongly on the blueshifted side of the disc, and again mainly at elevated layers traced by 12CO. Interestingly, these substructures are not as extended nor as symmetric in azimuth as their velocity counterparts at similar radii. In fact, they seem to open linearly to wider orbits as a function of azimuth on the blueshifted far side, as if they were part of a spiral structure. The radial coincidence of these line width increments possibly driven by radially localised turbulent motions, and the vertical flows detected in velocity residuals is yet another signature that supports the presence of a massive planet in this region of the disc. Line width enhancements have recently been reported as potential indirect signatures of planets in the discs of HD 163296 (Izquierdo et al. 2022), and TW Hydrae (Teague et al. 2022). Indeed, using 3D hydrodynamics and radiative transfer simulations, Dong et al. (2019) showed that a massive planet can trigger turbulent vertical motions on and around its orbit, which should be observable by quantifying the broadening of line profiles, as long as the planet mass and the angular and spectral resolution of the observations are high enough. Considering this and the fact that the vertical velocity signatures are in good agreement with those expected for a planet-driven buoyancy spiral, we propose a planet candidate in the disc of MWC 480 at an orbital radius of R = 245 au, namely at the radial location proposed by Teague et al. (2021), but at an azimuth of ϕ = −128°, where the peak in 12CO line width residuals is found. The location of this planet candidate is marked as a green cross in Fig. 10, bottom panel. See Fig. 4, top panel, for a Cartesian view of the line width perturbations in this disc.
Additionally, to understand how the retrieved velocity perturbations vary with the assumed vertical structure, we take this source as an example to investigate the impact of using an irregular emission surface on the retrieved velocity residuals. For this experiment, we fix the 12CO model surface to the non-parametric emitting layer extracted by Paneque-Carreño et al. (2023), and re-run DISCMINER to find a new set of parameters for the remaining attributes summarised in Table 1. For better comparison, we also fix the disc inclination in both the smooth and irregular-surface model to the value found by Liu et al. (2019) from mm continuum modelling. As illustrated in Fig. B.15, we find that using an irregular surface leads to higher intensity and velocity residuals overall. In particular, we note that the irregular-surface model displays quadruple patterns in the velocity residuals at large radii (R > 400 au), but also around more specific radial separations (e.g. at R = 300 au and R = 400 au). The former suggests that the surface tapering found by the geometrical, non-parametric method is possibly over-estimated due to the influence of the increasingly strong sub-Keplerian rotation at large radial separations. The latter, on the other hand, may indicate that some of the vertical bumps and troughs retrieved by the geometrical method are instead modulated by the radially confined velocity perturbations present in this disc, and not the consequence of actual variations in the disc surface density. Hence, we only consider smooth parametric surfaces found by DISCMINER for the rest of the analysis.
Fig. 10 Azimuthal deprojection of velocity and line width residuals obtained for the disc of MWC 480 as observed in 12CO and 13CO. The horizontal dashed and solid lines indicate the radial location of dust gaps and rings, respectively. The purple dotted line is the radial location of kink-like features observed in the datacube channels. The vertical lines denote the azimuthal location of the disc main axes as projected on the sky. The black lines overlaid in the line width residuals highlight the location of extended line width enhancements, likely associated to non-thermal motions. The green cross marks the location of the planet candidate reported in Sect. 3.2 based on 12CO line width enhancements. Cartesian versions of these maps can be found in Figs. 4 and B.7. |
HD 163296
In the disc of HD 163296 there are particularly prominent extended perturbations, most of them driven by pressure modulations due to the presence of gas gaps which trigger deviations from Keplerian rotation in the azimuthal component of the velocity field (see e.g. Teague et al. 2018a; Zhang et al. 2021). These azimuthal fluctuations are seen as ring-like structures with a sign flip around the disc minor axis owing to the projection of the velocity field along the line of sight of the observer (see Fig. 5, top middle panel). As illustrated in Fig. 5 of Teague et al. (2019b), super-Keplerian azimuthal velocities would appear redshifted on the redshifted side of the disc, and blueshifted on the blueshifted side, and the other way around for sub-Keplerian velocities. Strong sub-Keplerian azimuthal velocities are also seen at large radii (R > 350 au) induced by a sharp pressure gradient possibly driven by the tapered shape of the surface density profile expected near the edge of the disc (Dullemond et al. 2020; Zhang et al. 2021). Additionally, there is an elongated spiral-like perturbation around an orbital radius of R = 260 au spanning over a large azimuthal extent between ϕ ~ 90° to −60°, which has been recently linked to a massive planet at the same radius (Teague et al. 2021; Izquierdo et al. 2022; Calcino et al. 2022), originally proposed by Pinte et al. (2018a) after observation of an intensity kink in 12CO channels. Interestingly, this elongated signature is more strongly seen at high elevations traced by 12CO emission. In 13CO, the perturbation is just barely seen and hindered by the background noise, whereas in C18O the azimuthal component of the perturbation seems to dominate, suggesting that this spiral is also likely triggered by buoyancy resonances as in the disc of MWC 480, and that the perturbations near the disc midplane could be associated with a shallow density gap carved by the planet candidate.
AS209
The disc of AS 209 displays downward and upward vertical flows around contiguous annuli centred at R ~ 130 au and R ~ 190 au, respectively. As illustrated in Fig. 11, and as opposed to MWC 480, the vertical flows in this disc are clearly seen in 12CO and 13CO simultaneously, meaning that these motions span coherently across large vertical extents and are therefore less compatible with planet-driven buoyancy spirals. Instead, non-planetary mechanisms such as the vertical shear instability, which can trigger vertical motions across a wide range of scale heights (see e.g. Barraza-Alfaro et al. 2021), or magnetically driven winds via ambipolar diffusion, studied by Galloway-Sprietsma et al. (2023) for this specific source, are favoured in this case.
Even though we do not detect any localised velocity perturbation in the disc of AS 209 near the CPD proposed by Bae et al. (2022b) at a radius of R = 203 au, we note that the azimuthally coherent morphology of the large-scale upward velocity flow traced by 12CO, at R ~ 190 au, is disrupted around the azimuth of the planet candidate where it extends to wider orbits forming an arc-like pattern that surrounds the suggested location of the planet. Furthermore, we spot a prominent ~50 ms−1 line width enhancement in 12CO at the same location of the CPD candidate that supports the presence of this object. As illustrated in Fig. 11, bottom panels, this feature is ~50 au wide around the vicinity of the CPD, and becomes narrower at larger azimuthal separations, extending coherently for almost 180° near the CPD orbit. The same azimuthally extended substructure is found in 13CO line widths between R = 140–180 au radii. However, the line width bump observed in 12CO around the immediate vicinity of the CPD no longer remains in this tracer, suggesting that the enhanced velocity dispersion near the planet candidate occurs primarily in the atmosphere of the disc. The peak of the 12CO line width residuals surrounding the CPD candidate is located at an orbital radius of R = 210 au, and an azimuth of ϕ = −113° in the disc reference frame.
We stress that the seemingly localised kink-like feature observed in the 12CO channel maps of AS 209 in the south of the sky does not originate from a localised kinematic perturbation. Instead, it is driven by the upward meridional flow which extends across the entire azimuth of the disc as demonstrated in this section. Figure 12 shows isovelocity contours extracted from the data and model channel maps, and highlights the location of vertical flows to illustrate that the kinks observed in intensity channels are far from probing localised features in the kinematics. Even though kinks have been associated with the presence of planet-driven perturbations (Pinte et al. 2018a), such as Lindblad spirals dominated by radial and azimuthal velocities (Rabago & Zhu 2021; Bollati et al. 2021), we emphasise that any type of velocity perturbation, strong enough along the line of sight, can cause kinks in intensity channels regardless of its physical origin and/or predominant velocity component.
Finally, in Fig. B.16 we illustrate that the azimuthal modulation of the intensity field driven by foreground cloud absorption in this source does not translate into variations in the velocity perturbations that can be quantified. The reason is that the relatively low inclination and flat vertical structure of this disc make the emitted line profiles to be nearly symmetric in frequency, and thus uniform intensity variations do not lead to differences in the retrieved velocities.
IM Lup
The disc of IM Lup exhibits strong velocity perturbations near the projected minor axis, potentially dominated by azimuthal velocities (Fig. 5, bottom left panel). However, due to the strong diffuse emission in-between 12CO layers not accounted by the DISCMINER model, these signatures may not represent physical disturbances and should be interpreted with care. On the other hand, although IM Lup is thought to host a massive gas disc (~0.1 Μ⊙, Lodato et al. 2022), it is still unclear whether it is gravitationally unstable. This disc displays complex, nearly symmetric, spiral structures in the mm continuum, typical of GI, but with varying pitch angle as function of radius which would be in better agreement with the presence of an embedded companion (see discussions by Huang et al. 2018b; Verrios et al. 2022). We do not detect clear spiral signatures in the kinematics of any of the CO isotopologues, which would favour a gravitationally (at least marginally) stable scenario in this source. Instead, we identify arc-like fluctuations both in velocity and peak intensity residuals on the blueshifted side of the outer disc between R ~ 300–700 au.
GMAur
The disc of GMAur displays large-scale velocity perturbations both in the atmosphere and at layers closer to the midplane (Figs. 5 and B.7). It is known that this system is potentially in late interaction with material infalling from the remnant envelope (Huang et al. 2021), which leads to the presence of intensity and velocity spirals and trails, or filamentary structures, that appear most notably at high elevations where the interaction is strongest. The prominent wedge of redshifted velocities between ϕ ~ 45−90° in the disc reference frame (Fig. 5), and of line width enhancements in the same region (Fig. B.9), is in good agreement with that inflow of material. This wedge is also seen as an excess of integrated 12CO intensity in Fig. 17 of Huang et al. (2021) after subtraction of a Keplerian mask. We note that multiple spiral-like features are also present in peak intensity residuals of this system as illustrated in Fig. B.12. However, in 13CO we only tentatively identify spirals in the kinematics, highlighted by solid lines in Fig. B.7. In a non-planetary scenario, these can be the consequence of gravitational instabilities developing in a disc (Kratter & Lodato 2016). However, Schwarz et al. (2021) demonstrate that in spite of the possibly high mass reservoir available in this disc (~0.2 Μ⊙, see also Lodato et al. 2022), the gas temperature is high enough across much of the disc extent to prevent gravitational instabilities to develop. Since we only detect spiral features convincingly in intensity and velocity residuals at high elevations traced by 12CO, we argue that they are more likely the result of infalling material perturbing the gas disc surface density and temperature rather than a consequence of GI.
Fig. 11 As Fig. 10 but for the disc of AS 209. The green circle marks the location of the CPD proposed by Bae et al. (2022b). Cartesian versions of these maps can be found in Figs. 5, B.7 and B.9. |
4 Gas substructures as sites to search for planets
In addition to the localised features and the non-axisymmetric perturbations manifested as spirals and filamentary structures in the disc kinematics (see Sect. 3), embedded planets can also induce axisymmetric variations in the gas surface density and temperature of the host disc in the form of annular gaps and cavities. These substructures translate into pressure gradients capable of producing observable, axisymmetric velocity modulations, both in the rotational component of the velocity field (Kanagawa et al. 2015), and in the vertical direction through the meridional circulation of material from the disc atmosphere onto the midplane and vice-versa (Morbidelli et al. 2014). Moreover, these features are also expected to trigger localised fluctuations in the radial profiles of molecular peak intensities (see e.g. Facchini et al. 2018b) and line widths (Izquierdo et al. 2021), which can be useful to understand whether the underlying pressure modulations are dominated by density and/or temperature variations. How do these observables (anti−)correlate with each other, if at all, and where do they sit radially in this sample of discs? In this section, we focus on the analysis of azimuthally averaged radial profiles of the observable quantities explored in this work to provide a general view of the physical structure and dynamics of the discs in the sample, as well as the connection between annular gas substructures and the planet candidates proposed in Sect. 3 and summarised in Table 2.
Fig. 12 Illustrating isovelocity contours and the presence of coherent velocity substructures in the disc of AS 209 as probed by 12CO and 13CO J = 2−1 line emission. Black solid and dashed lines show isovelocity contours from the data and the best-fit Keplerian model, respectively. For reference, the isovelocity levels correspond to the same velocity channels of the intensity maps displayed around the main panels. Coherent sub- and super Keplerian substructures identified by FILFINDER (Koch & Rosolowsky 2015) in velocity residuals (Figs. 5 and B.7) are overlaid as blue and red thick lines. The green circle marks the location of the CPD candidate proposed by Bae et al. (2022b) in the disc reference frame. We note that the intensity kinks most prominently seen in the south of the sky in 12CO channels, between 3.5 and 5.5 km s−1, are not driven by localised velocity perturbations, but instead are the consequence of a large-scale blue-shifted flow that spans over the entire azimuth of the disc at around ~200 au, and appears at both scale-heights traced by 12CO and 13CO. In Sect. 4.1, we show that such a coherent signature is the result of upward vertical flows. |
4.1 Azimuthal and meridional velocity flows
In Sect. 2.2.2, we demonstrate that it is possible to readily obtain an estimate of the azimuthal and vertical components of the velocity perturbation field by computing azimuthal averages of different kinds of velocity residual maps. These azimuthal and meridional velocity flows are presented in the top panels of Fig. 13, for the disc of MWC 480, and in Figs. 14, B.17, B.18 and B.19, for the discs of HD 163296, AS 209, IM Lup, and GM Aur, respectively. Along both types of velocity profiles, we mark the location of pressure maxima, traced by negative velocity gradients in azimuthal velocities, and the direction of the meridional flow on local maxima and minima in vertical velocities. To assess the effect of continuum subtraction on these velocity profiles, we performed the same analysis on intensity cubes of all three CO isotopologues for MWC 480 without continuum subtraction and, as illustrated in Fig. B.20, we found no differences with respect to the continuum-subtracted cubes that could affect any of the results summarised in this section.
MWC 480
From Fig. 13, we note that in the disc of MWC 480 the positive and negative radial gradients in the azimuthal component of the velocity, retrieved for all CO isotopologues, are in excellent correlation with the radial location of gaps and rings observed in the mm continuum, respectively. This result suggests that the pressure modulations present in the gas component of this system have been there long enough for the mm-sized dust grains to be trapped in annular structures at the currently observed locations. It is also worth highlighting the prominent meridional circulation of material flowing up and down with respect to the disc midplane in this system. In the radial section comprised between the dust ring B165 and a radius of R ~ 300 au, we identify strong vertical motions with amplitudes of up to 50 m s−1 in the disc atmosphere, traced by 12CO at R = 193 au and R = 245 au. This is in line with the vertical velocities found by Teague et al. (2021) in 12CO using a different extraction method (see e.g. Fig. 5 in their paper). The fact that the line-of-sight velocities along this structure are highly coherent in the azimuthal coordinate, as shown in Fig. 10, and that its underlying pitch angle is very low as a function of radius, allows us to estimate with good precision the magnitude of the underlying vertical motions using the decomposition method introduced in Sect. 2.2.2. Furthermore, it is also notable that although all CO isotopologues trace similar variations in the azimuthal velocity component, only 12CO (tracing z/r ~ 0.2 scale heights, Fig. B.5) displays strong fluctuations in the vertical component while 13CO and C18O (tracing z/r < 0.1) remain nearly unperturbed. As discussed in Sect. 3.2, this vertical modulation of vertical velocities strongly supports the presence of buoyancy spirals driven by an embedded planet, originally proposed by Teague et al. (2021) through the observation of tightly wound substructures in 12CO temperatures. In Sect. 3.2, we also use line width information to propose a planet candidate at an orbital radius of R = 245 au and an azimuth of ϕ = −128° in this disc.
HD 163296
In the top panel of Fig. 14 we illustrate that most of the gas pressure bumps present in the disc of this source, traced by negative radial gradients in the azimuthal component of the velocity field, are also co-spatial with the radial location of dust rings observed in the mm continuum. This suggests that the latter are indeed the consequence of pressure traps as confirmed by Rosotti et al. (2020b). In fact, the even better match between the location of these rings and pressure bumps traced by 13CO and C18O, which are ~50% closer than 12CO to the disc midplane (see Fig. B.5, bottom row), further supports this finding. Similarly, the positive velocity gradients identified in the same velocity component at the radial locations of R ~ 50, 85, 135 au, in all CO isotopologues, are in excellent agreement with the gas gaps proposed by Isella et al. (2016) based on radiative transfer models of this disc. Also, from the second top panel in Fig. 14, we note that the meridional circulation of material is very pronounced in this system too. Upward and downward gas flows are particularly prominent around the location of millimetre dust and gas substructures between R ~ 50–160 au, and co-spatial with pressure minima and maxima, respectively. Further out, there is an upward velocity flow centred at R ~ 240 au near the planet candidate associated with the K260 kink. This feature is strongly detected in the disc atmosphere traced by 12CO, but decays rapidly towards lower scale heights probed by 13CO and C18O. As discussed in Sect. 3.1, such a prominent vertical dependence could be interpreted as due to planet-driven buoyancy spirals which, as opposed to Lindblad spirals, are expected to develop more significantly in the disc atmosphere and become weaker as the disc midplane is approached. We note that this radial profile of 12CO vertical velocities is equivalent to that of Teague et al. (2019a, 2021), but with opposite sign. We have cross-checked with the corresponding authors via private communication and agreed that the profile presented in Fig. 14, second top panel, is correct.
Fig. 13 Azimuthally averaged profiles of velocity (top), line width (middle) and peak intensity (bottom) residuals obtained for 12CO, 13CO, and C18O, for the disc around MWC 480. The radial location of millimetre dust gaps and rings is illustrated as dashed and solid lines, respectively. The radial distance of the most prominent kink apparent in 12CO channel maps is shown by the dotted purple line. Strong pressure bumps traced by minimal velocity gradients are marked with minuses in the top panel. Peak meridional flows hinting at gas moving away and towards the disc midplane are highlighted with arrows pointing up and down in the second top panel. The planet marker indicates the orbital radius of the planet candidate proposed at R = 245 au (see Table 2). |
Fig. 14 As Fig. 13 but for the disc around HD 163296. The planet markers indicate the orbital radii of the planet candidates associated with the localised velocity perturbations P94 and P261 (see Table 2). |
AS209
The average residual velocities computed for the disc of AS 209 in the three CO isotopologues are the most dissimilar of the sample (see Fig. B.17). Due to the fact that the model background velocity of all molecules is different – because stellar masses and orientation parameters are allowed to vary freely in all models – a velocity offset with respect to the Keplerian model background may be expected. However, taking systematic offsets aside, the three isotopologue velocity profiles are topologically different. On top of that, only the pressure bumps traced by 12CO coincide, within a beam, with the radial location of the millimetre dust rings, B74 and B121. The fact that the 13CO and C18O velocities behave differently, suggests that the radial pressure gradients in this disc fluctuate strongly as the midplane is approached. This was already suggested by Teague et al. (2018b), who proposed a vertically varying pressure profile to explain the radial offset between scattered light and continuum dust rings observed in AS 209, assuming that both are driven by pressure traps. Nevertheless, it is still puzzling that the pressure bumps in 13CO and C18O are even further from the location of millimetre dust rings in spite of being closer than 12CO to the disc midplane4. A plausible scenario could be that some of the velocity components near the midplane of the disc are far from axisymmetric, in which case the azimuthally averaged velocities would not provide a good description of the actual azimuthal component of the velocity as it would be contaminated by non-axisymmetric radial and/or vertical motions. We identify, however, a strong positive gradient in both the azimuthal velocities of 12CO and 13CO between R = 85–110 au, around the millimetre dust gap D100, which Fedele et al. (2018) modelled by invoking a low-mass planet (<0.1 MJup) at an orbital radius of R = 103 au. On the other hand, this disc displays two fairly prominent vertical flows of material, one approaching the midplane at a radius of R ~ 130 au, and another moving away from the midplane at R ~ 190 au, with a comparable magnitude of the order of 50–100 m s−1 traced by both 12CO and 13CO. The latter signature coincides within a beam size with the radial location where a CPD candidate was proposed earlier by Bae et al. (2022b) at -200 au, and finds a plausible interpretation in the work of Galloway-Sprietsma et al. (2023), who conclude that this coherent upward flow could be the consequence of magnetically driven winds dominated by ambipolar diffusion, triggered by the possibly low densities around the orbit of the proposed planet+CPD system. Further details on the kinematic signatures identified around this planet candidate are summarised in Sect. 3.2.
IM Lup
The millimetre dust gap D116 and ring B133 of this disc co-locate very well with pressure minima and maxima traced by azimuthal velocity flows in all CO isotopologues (Fig. B.18), namely, across a wide range of vertical layers (z/r ~ 0.1–0.3, Fig. B.5). The dust ring B220, however, does not seem to align with any pressure bump in the gas, and therefore a different mechanism other than pressure trapping is likely acting on this radial section of the disc midplane. Additionally, we detect a downward meridional flow of material between R ~ 120–135 au, both in 12CO and 13CO, overlapping with the millimetre dust gap D116 and ring B133. Although no localised perturbation is detected in this region of the disc (see Fig. 9), we note that this meridional gas flow is near the planet candidate proposed by Pinte et al. (2020) and Verrios et al. (2022) at R = 117 au.
GM Aur
As illustrated in Fig. B.19, the strong interaction between this disc and material possibly infalling from a remnant envelope (see e.g. Huang et al. 2021), makes both azimuthal and vertical flows at radii R > 150 au to change substantially with varying scale height. At radii interior to this region, the azimuthal component of the velocity follows similar modulations for all CO isotopologues as opposed to the vertical component, which still exhibits large variations in the z–direction. We note that the millimetre dust ring B86 is highly correlated with a gas pressure bump present at all scale heights, traced by negative gradients in the azimuthal component of the velocities in all three CO isotopologues. Also, the most prominent pressure minima are well aligned with the radial location of millimetre dust gaps D68 and D142. Conversely, no clear pressure bump is found at the location of the B163 ring, suggesting that this dust substructure is either not the result of pressure trapping processes or that the external material is already biasing the derivation of rotation velocities as of this radius.
4.2 Surface density gaps traced by azimuthal velocity modulations and line widths
As demonstrated in Fig. 6 of Paper I, the presence of gaps in the surface density of a disc has an observable impact on the radial distribution of line widths from molecular emission. In particular, due to the reduced optical depth in a density gap, azimuthally averaged line widths from optically thick tracers are expected to be significantly smaller than those in the unperturbed background of the disc. Thus, if the gas temperature is high within the gap, this implies that if a positive radial gradient in azimuthal velocities and a minimum in line widths are detected simultaneously, we can assert that the pressure minimum associated with the gap is strongly dominated by low surface density values and not by its temperature structure. As illustrated by Rab et al. (2020, see e.g. their Fig. 7), if the high temperature of the gap was the dominant component, the azimuthal velocity flows would instead display negative radial gradients at the gap centre.
In Fig. 15, we make explicit the correlation between line width minima (right, top panel) and pressure troughs traced by peak velocity gradients (right, bottom), as a function of the orbital radius for the disc around MWC 480 in the three CO isotopologues. We note that the pressure troughs (in green) co-locate with line width minima (in light blue), and that the pressure bumps coincide with the location of dust rings and less shallow line widths. This analysis suggests the presence of at least two gas gaps centred at a radius of R = 76 au and R = 125 au, dominated by low surface densities, and indicates that the dust rings detected in previous mm continuum observations of this system, B98 and B165, are indeed the consequence of pressure traps. We also note that the same signatures hold in both 13CO and C18O, although less clearly than in 12CO. This is expected because for increasingly optically thinner tracers the dependence of line widths on the gas surface density is weaker (Hacar et al. 2016). Despite the fact that there are peaks in the intensity of all three CO isotopologues, likely associated with elevated temperatures at the radial location of these gaps, the simultaneous detection of line width and pressure troughs indicates that the low surface density in both regions dominates the shape of the pressure gradient responsible for the observed kinematic modulations.
Finally, we note that positive gradients in the azimuthal component of the velocity field agree excellently with the location of line width minima in all discs, especially in 12CO and 13CO emission. In Tables 3, 4, A.6–A.8, we present a summary of the radial locations of maxima and minima in velocity gradients and line widths for all sources. Conversely, in the same tables we report that intensity peaks and troughs show no clear correlation or anti-correlation with respect to the latter two observables. As discussed by Turner et al. (2012) and Rab et al. (2020), peaks and dips in the brightness temperature do not necessarily align with the location of gas density rings and gaps (see Fig. 15, top left vs top right panel). Therefore, we emphasise that a simultaneous analysis of the kinematics and widths of the molecular line profile revealed by emission of (at least marginally) optically thick tracers allows understanding the gas disc substructures and thus the potential sites of planet formation with greater robustness.
4.3 Strong vertical dependence of radial pressure in the disc of HD 163296, around the planet candidate P94
Taking a closer look to the residual velocity profiles of the disc of HD 163296 in Fig. 14, it is easy to note that between R = 86 au and R = 141 au, the rotation velocities traced by 12CO vary substantially in magnitude and extent with respect to those in 13CO and C18O. As these isotopologues emit at different altitudes above the midplane, this is evidence of the vertical dependence of the radial pressure gradient of the gas disc. However, the origin of these pressure variations remains unclear as yet. Do density and temperature gradients contribute similarly to the observed pressure gradients, or does one of them predominate?
The other two average profiles, computed from peak intensity and line width residuals, are helpful to gauge the relative contribution of density and temperature fluctuations to the radial pressure gradient of the disc. Consider the radial section of the disc between R ~ 75 au and R ~ 100 au. In this region, the residual peak brightness temperature of 12CO (Fig. 14, bottom panel) displays a positive gradient, followed by a negative gradient until it reaches a temperature trough at R ~ 120 au. The 12CO line widths (middle panel), on the other hand, exhibit a pronounced positive gradient between the line width minimum near D86 and the line width maximum at 125 au. Assuming that fluctuations in the brightness temperatures and line widths of 12CO are tightly linked to local variations in the gas temperature and density, respectively (see e.g. Hacar et al. 2016; Izquierdo et al. 2022), the coincidence of a temperature trough and a line width peak around 120 au suggests that the pressure maximum, revealed by the disc kinematics, is mainly driven by a gas ring at this radius.
Conversely, 13CO and C18O show the opposite behaviour in both the residual brightness temperatures and line widths within the same radial section. Nevertheless, the same argument made for the pressure bump traced by 12CO at R = 120 au holds for these two molecules, but this time closer-in at R = 100 au. The pressure maximum revealed by the negative velocity gradient of both isotopologues at this radius coexists with a line width peak and a temperature trough, indicating that this pressure bump is most likely driven by gas density enhancements5. The radial shift between this pressure maximum, and the one traced by 12CO, is in agreement with the findings of Rab et al. (2020), who used thermochemical models to show that the relative contribution of temperature and density fluctuations to the pressure gradient may vary substantially at different radii and heights in this disc.
Morphologically, we are possibly witnessing a density substructure that is opening along the vertical direction, like a funnel, around the millimetre dust gap, D85, and a gas gap at the same radial location, possibly carved by the planet candidate, P94, proposed by Izquierdo et al. (2022) and ratified in this work through localised velocity and line width perturbations detected in multiple CO isotopologues. This geometry would explain why the radial separation between the pressure minimum around D86, due to a gas gap, and the pressure maximum at R = 120 au traced by 12CO, is almost twice as large as the separation between pressure minima and maxima in the same radial section at the elevations traced by 13CO and C18O.
Fig. 15 Summary of azimuthally averaged residuals extrapolated onto a 2D grid to aid visualisation of (anti-)correlations between the different line profile observables extracted from the disc of MWC 480. The radial location of dust gaps and rings is illustrated as dashed and solid lines, respectively. The radial distance of the most prominent kink apparent in 12CO channel maps is marked by the dotted purple line. |
Summary of radial substructures derived from azimuthally averaged residuals for the disc of MWC 480.
4.4 Non-axisymmetric variations in the radial profiles of the disc around MWC 480
Not all variations in the line profile observables explored in this study are axisymmetric. There are numerous signatures that appear spatially localised (as presented in Sect. 3.1), and others that span over large azimuthal and radial portions, at times forming arcs along concentric annuli, but also filament-like structures that stretch irregularly across the discs (as in Sect. 3.2).
To quantify and better visualise these asymmetries, it is helpful to compute azimuthal averages over narrower azimuthal ranges instead of considering 2π averages. This is illustrated in Fig. 16 for the disc of MWC 480, where we show azimuthal averages of the three residual observables obtained for 12CO, in the four quadrants delimited by the disc main axes as projected in the sky. It is clear that the average residual profiles, especially those from centroid velocities, vary significantly from one quadrant to another. However, we note that between the outermost millimetre dust ring, B165, and a radius of R ~ 350 au, the average velocities in each quadrant of the disc seem to fluctuate less with respect to each other than at closer-in annuli. It is in this region where the azimuthal coherence of the perturbation reveals the three-dimensional nature of the velocity field, seemingly dominated by vertical motions as illustrated in Fig. 13. Likewise, the radial sign flips are the consequence of contiguous anti-parallel vertical flows, which could be related to the presence of planet-driven buoyancy spirals (Bae et al. 2021; Teague et al. 2021), as discussed earlier in Sect. 3.2.
Furthermore, as demonstrated in line width residuals in the top left panel of Fig. B.9, the disc around MWC 480 displays prominent arc-like enhancements in the 12CO line widths, more strongly seen on the blueshifted side around orbital radii between R ~ 225 and R ~ 325 au, which we associate with the presence of an embedded planet at R = 245 au radius and ϕ = −128° azimuth (see Sect. 3.2). With the by-quadrant analysis presented in Fig. 16, we infer that the line width asymmetries at R ~ 225 au, in the third quadrant of the disc, are ~50 m s−1 higher than those in the mirroring fourth quadrant. Moreover, we note that these line width enhancements are anti-correlated with the temperature peaks found at the same radial and azimuthal location, suggesting that thermal broadening could not be the main driver of these arc-like features in line widths. Instead, turbulent motions should thus be favoured as the cause of such signatures.
Fig. 16 Illustrating the presence of non-axisymmetric features in the disc of MWC 480, as observed in 12CO. Left column: velocity, peak intensity, and line width residuals deprojected on the disc reference frame. Right column: azimuthally averaged residuals extracted per quadrant to highlight azimuthal variations in the different observables. The colour of each radial profile indicates the quadrant of the disc where it was computed according to the colouring code in the top right corner. Quadrants to the right of the disc vertical axis correspond to the redshifted side of the disc. |
5 Discussion
5.1 Morphology of the P94 perturbation in HD 163296
The detection of localised velocity signatures in different tracers around the planet candidate P94 allows us to map the vertical structure of the associated velocities. We note that both the perturbation detected in 12CO, or P94, and the one detected in 13CO, or P81, exhibit a Doppler shift pattern expected from the gravitational and hydrodynamic interaction of a planet with the disc fluid in its immediate vicinity (see e.g. Rabago & Zhu 2021). However, the morphology of these two signatures is clearly different between each other. The redshifted part of the P94 perturbation is in an outer orbit with respect to the centre of the proposed planet at R = 94 au, while the blueshifted part is in an inner orbit (see Fig. 7). The structure of the P81 perturbation, on the contrary, manifests as a redshifted part which coincides azimuthally with that of P94, but it is located interior to the orbital radius of the blueshifted component at its peak. Although not significantly detected by our clustering algorithm, we identify a similar pattern in C18O. This suggests that the perturbation sitting at relatively low scale-heights, traced either by 13CO or C18O, is dominated by the radial component of the velocity, which would thus be in agreement with the velocity structure expected for planetary wakes along Lindblad spirals near the disc midplane (see Fig. 4 of Rabago & Zhu 2021). A similar kinematic feature suggestive of a planet-driven Lindblad spiral was recently reported by Teague et al. (2022) in the disc of TW Hydrae.
5.2 Line broadening as a tracer of massive planets
The simultaneous detection of localised velocity and line width perturbations pointing to the presence of both coherent motions and enhanced velocity dispersions around the planet candidate P94 in the disc of HD 163296 opens the door to the use of yet another observable – non-thermal line broadening – to search for embedded planets in discs, as originally predicted by Dong et al. (2019). Consequently, we employ these two observables to propose the most plausible azimuthal location of a protoplanet candidate possibly driving buoyancy spiral perturbations in the disc of MWC 480. In the disc of AS 209, however, while we identify strong line width enhancements around the vicinity of the CPD candidate reported by Bae et al. (2022b) at R = 203 au, we do not detect velocity signatures that can be associated with localised perturbations linked to the proposed object. Nevertheless, we highlight that this disc is known to host a prominent gas gap revealed by integrated intensity maps at an orbital distance of ~200 au (Guzmán et al. 2018; Law et al. 2021a), leaving the planet candidate with little surrounding material and hence making the detection of coherent small-scale velocity variations a harder challenge. The presence of strong line broadening along with the absence of localised coherent velocities around the CPD candidate suggests that the gas motions in this region are predominantly turbulent.
5.3 Localised velocity perturbations are not so common
The fact that highly localised velocity perturbations expected from the interaction between discs and giant planets are only found in one of the five discs in the sample may indicate that massive planets do not commonly form at large separations from the central star, and if they do, they should migrate rapidly towards lower orbits. Alternatively, it is likely that the influence of other physical mechanisms such as gravitational and hydrodynamic instabilities along with the presence of pressure substructures, as well as the effect of multiple embedded planets and stellar companions, is sufficient to trigger deviations from Keplerian rotation strong enough to hide the localised features driven by single massive planets. A systematic characterisation of composite hydrodynamic simulations including several physical mechanisms at once will be key in the field of disc kinematics to unambiguously distinguish between planetary and non-planetary signatures retrieved from high resolution observations of protoplanetary discs.
5.4 Correlation between pressure bumps and mm dust rings
Using multiple molecular lines from the same sample of discs, Jiang et al. (2022) found no clear (anti−)correlations between line intensities and the location of dust rings and gaps observed in mm continuum data, leaving an open debate on the mechanisms that led to the formation of the millimetre dust substructures present in these discs. However, our analysis of azimuthal velocity flows indicates that at least nine of the eleven millimetre dust rings covered in this study do coexist with maxima in the radial gas pressure profile traced by different CO isotopologues. This finding indicates that there is a clear relationship between the gas pressure variations, more effectively captured by velocity modulations, and the substructures observed in the mm continuum, suggesting that gas pressure traps largely dominate the formation of dust rings in these sources.
6 Conclusions
In this work we apply the DISCMINER modelling and analysis techniques on high angular and velocity resolution archival data of 12CO, 13CO, and C18O J = 2−1 line emission from five discs around the young stars MWC 480, HD 163296, AS 209, IM Lup, and GM Aur, observed as part of the ALMA large program MAPS (Öberg et al. 2021). To study deviations from Keplerian rotation and, more generally, to identify fluctuations in the underlying intensity and velocity fields, we first produce best-fit channel maps for all discs and CO isotopologues, assuming they are smooth in intensity and Keplerian in rotation velocity. Next, we compare morphological properties of the data and the resulting model line profiles, such as amplitude, width, and velocity shift, to search for localised and extended fluctuations hinting at velocity and density substructures in the discs. Finally, we apply a clustering algorithm to search for localised velocity and line width perturbations possibly driven by massive unseen planets, and study the presence of extended kinematic substructures that may also be the consequence of embedded planets. Our main findings are as follows,
6.1 Giant planets around HD 163296, MWC 480, AS 209
Localised velocity features are only detected in the disc of HD 163296, where they are robustly seen in 12CO and 13CO simultaneously, and tentatively in C18O. These signatures are identified and confirmed around the location of the planet candidate P94, proposed by Izquierdo et al. (2022) using 12CO observations with lower spectral resolution. Instead, large-scale kinematic and line width substructures are more common in the other discs. In the discs of MWC 480 and AS 209, we spot numerous ring−, arc− and spiral-like signatures in all line profile observables, spanning over large azimuthal extents, which may also be driven by unseen planets.
Some of the large-scale features observed in the kinematics of the discs around MWC 480 and AS 209 are compatible with anti-parallel vertical flows of material. In MWC 480, these meridional flows take place primarily in the atmosphere of the disc, while in AS 209 they propagate freely from layers in the atmosphere down to regions closer to the disc midplane, as traced by 12CO and 13CO emission. This varying modulation of meridional velocities as function of scale height suggests that such signatures are likely linked to different driving mechanisms in each disc. Planet buoyancy spirals are favoured in MWC 480, whereas MHD winds or vertical shear instabilities provide a better interpretation of the kinematic signatures observed in AS 209.
No localised velocity or intensity signatures are found around the location of the CPD candidate proposed by Bae et al. (2022b) in the disc of AS 209 in any of the CO isotopologues. We emphasise that kinks extracted through visual inspection of intensity channels should not be attributed to localised kinematic features without prior inspection of the underlying velocity field. The presence of vertical flows in this disc, spanning over large azimuthal extents, are the cause of the strong kinks observed in 12CO channels around the orbital radius of the CPD candidate.
We detect both localised and extended line width enhancements in at least two CO isotopologues near the location of the planet candidate P94 in HD 163296, and around the vicinity of the CPD candidate proposed by Bae et al. (2022b) in AS 209. These features would be in agreement with turbulent motions triggered by massive embedded planets as theorised by Dong et al. (2019).
Extended line width enhancements are also seen in the atmosphere of the disc around MWC 480, in the same radial section where anti-parallel vertical flows are detected, namely between ~200–300 au. Likewise, this could be a hint of turbulent motions driven by a young massive planet embedded in this region of the disc. Using this information, we propose a planet candidate at an orbital radius of R = 245 au and an azimuth of ϕ = −128° in this system.
Although not statistically representative, a ratio of four massive planets found at large orbital radii in five protoplanetary discs is noteworthy. If this ratio holds up in future surveys, this could indicate that gravitational instabilities and/or rapid planetary migration are common mechanisms of evolution and planet formation in circumstellar discs.
6.2 Non-detections in the discs of IM Lup and GM Aur
The discs of IM Lup and GM Aur do not exhibit localised nor large-scale velocity perturbations that can be associated with the presence of massive planets.
Kinematic analysis of these discs does not reveal clear indications of spiral arms typical of gravitational instabilities either, despite the possibly high mass reservoir available in both systems (Lodato et al. 2022). Strong spiral features are only detected in the temperature structure of the disc of GM Aur most notably at high elevations traced by12 CO, which are likely the result of interaction between the disc atmosphere and infalling material from the remnant envelope rather than a consequence of gravitational instabilities.
The fact that no clear hints of massive planets are found in the probably most massive discs in the sample is interesting but not statistically conclusive. A systematic study of observable signatures driven by embedded planets in massive environments, extracted from hydrodynamic simulations, would help understand if the disc mass has a significant impact on the detectability of planet-driven perturbations.
6.3 Gas gaps, and relationship between pressure bumps and dust substructures
Azimuthal averages of standard and absolute velocity residuals are sufficient to disentangle the vertical and azimuthal components of the velocity perturbations, respectively. Azimuthal velocity flows are prominent in all discs, with magnitudes of the order of 50–150 m s−1, suggestive of strong pressure variations as function of orbital radius. Meridional flows are robustly detected in all discs, with average amplitudes lower than 50 m s−1. No strong correlation between gas gaps, traced by positive azimuthal velocity gradients, and meridional flows is found.
Most of the pressure bumps and minima detected in all discs appear to be dominated by surface density fluctuations. Also, we identify prominent vertical variations in the radial pressure profile of the disc around HD 163296, in the same radial section where the planet candidate P94 is proposed.
Nine out of eleven millimetre dust continuum rings in the sample are co-located with pressure bumps traced by negative rotation velocity gradients, indicating that pressure traps are a common mechanism for the formation of dust substructures in these discs.
A simultaneous analysis of velocity and line width variations in (at least marginally) optically thick tracers provides a robust assessment of the underlying surface density substructures responsible for the configuration of axisymmetric azimuthal and vertical flows in gas discs.
Our study demonstrates that a comprehensive analysis of line profile observables, including velocity, line width and intensity, from multiple molecular tracers, provides unparalleled access to the three-dimensional physical structure of circumstellar discs and is indispensable towards a robust detection and characterisation of planet formation sites.
Acknowledgements
The authors thank the anonymous referee for their insightful comments and discussions that helped us improve the presentation of results and the clarity of the methodologies introduced in this work. This paper makes use of the following ALMA data: ADS/JAO. ALMA#2018.1.01055.L. 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 National Radio Astronomy Observatory is a facility of the National Science Foundation operated under cooperative agreement by Associated Universities, Inc. This work was partly supported by the Italian Ministero dell Istruzione, Université e Ricerca through the grant Progetti Premiali 2012 - iALMA (CUP C52I13000140001), by the Deutsche Forschungs-gemeinschaft (DFG, German Research Foundation) – Ref. no. FOR 2634/1 TE 1024/1-1, and by the DFG cluster of excellence Origins (www.origins-cluster.de). This project has received funding from the European Union’s Horizon 2020 research and innovation programme under the Marie Sklodowska-Curie grant agreement No 823823 (DUSTBUSTERS) and from the European Research Council (ERC) via the ERC Synergy Grant ECOGAL (grant 855130). SF is funded by the European Union under the European Union’s Horizon Europe Research & Innovation Programme 101076613 (UNVEIL). GR is funded by the European Union (ERC DiscEvol, project number 101039651). GR acknowledges support from the Netherlands Organisation for Scientific Research (NWO, program number 016.Veni.192.233) and from an STFC Ernest Rutherford Fellowship (grant number ST/T003855/1). EvD acknowledges support from the ERC grant agreement No. 101019751 MOLDISK. Views and opinions expressed are however those of the author(s) only and do not necessarily reflect those of the European Union or the European Research Council Executive Agency. Neither the European Union nor the granting authority can be held responsible for them.
Appendix A Supporting tables
Best-fit model parameters obtained for channel maps of the disc of MWC 480.
Best-fit model parameters obtained for channel maps of the disc of HD 163296.
Best-fit model parameters obtained for channel maps of the disc of AS 209.
Best-fit model parameters obtained for channel maps of the disc of IM Lup.
Best-fit model parameters obtained for channel maps of the disc of GM Aur.
Appendix B Supporting figures
Here we provide complementary figures that support some of the results of the work, or figures that were not presented in the main text for reasons of space. In Figures B.1–B.4, we compare data and model channel maps for all CO isotopologues and discs but MWC 480, whose channel maps were already shown in Fig. 2. In Figure B.5, we present rotation curves and azimuthally averaged profiles of peak intensities and line widths, measured on the data alone for all isotopologues, considering the disc vertical structure which is also shown. In Figure B.6, we illustrate the DISCMINER modelling flow and analysis capabilities. We show residual maps of 13CO and C18O centroid velocities in Figures B.7 and B.8. Figures B.9, B.10, and B.11, illustrate 12CO, 13CO, and C18O line width residuals, respectively, and Figures B.12, B.13, and B.14 show peak intensity residuals for the same tracers. In Figure B.15, we show the effect of using parametric and non-parametric emission surfaces on the retrieved intensity and velocity residuals from the 12CO disc of MWC 480. In Figure B.16, we illustrate the impact of assuming an azimuthally symmetric and asymmetric intensity model (accounting for foreground absorption) on the resulting intensity and velocity residuals from the 12CO disc of AS 209. In Figures B.17, B.18, and B.19, we present radial profiles of azimuthal and vertical velocity components, as well as azimuthally averaged line width and peak intensity residuals for the discs of AS 209, IM Lup and GM Aur, respectively, for all isotopologues. A version of the same plot can be found in the main text in Fig. 13 for MWC 480, and in Fig. 14 for HD 163296. Finally, in Figure B.20, we show the effect of continuum subtraction on the velocity profiles extracted for the disc of MWC 480.
Fig. B.5 Additional observables obtained from the DISCMINER analysis of the discs around MWC 480, HD 163296, AS 209, (IM Lup and GM Aur in next page) as traced by 12CO, 13CO, and C18O. Top row: Rotation curves computed from azimuthal averages of deprojected line-of-sight velocities. For reference, the grey dashed line in each panel is the Keplerian rotation curve of the best-fit model found for 12CO. Second row: Azimuthally averaged profiles of line widths, defined as the standard deviation of Gaussians fitted to the data line profiles. For reference, the purple dashed line in each panel represents the thermal broadening computed at the 13CO peak brightness temperature. Third row: Azimuthally averaged profiles of peak intensities converted to brightness temperatures using the Rayleigh-Jeans law. Bottom row: Elevation of upper and lower emission surfaces. Apart from intensity profiles, all panels in each row share the same x and y-axis extent. Vertical dashed, solid, and dotted lines mark the radial location of dust gaps, rings, and kinks observed in the channel maps of 12CO, respectively. |
Fig. B.6 Schematic representation of the DISCMINER modelling flow and capabilities. DISCMINER produces best-fit channel maps of the input data cube by modelling the disc vertical structure and orientation, stellar mass, and attributes that shape the line profile morphology, assuming that the disc emission is smooth and Keplerian. The model and data channels are then compared through three types of moment maps. Large-scale analyses of these maps allow to look into the disc dynamical and thermodynamic structure, but are also useful to reveal potential signatures driven by planetary and non-planetary mechanisms, such as meridional flows or spirallike perturbations. Small-scale analysis of residuals from these maps leads to the retrieval of the potential location of embedded planets through the detection of localised velocity perturbations and sites of enhanced velocity dispersion. |
Fig. B.7 As Fig. 5 but for 13CO. This isotopologue traces vertical layers between 0.1 < z/r < 0.15 in this sample of discs. |
Fig. B.8 As Fig. 5 but for C18O. This isotopologue traces vertical layers between 0.0 < z/r < 0.1 in this sample of discs. |
Fig. B.15 Illustrating the impact of using an irregular (top) and a smooth parametric (bottom) emission surface on the peak intensity and velocity residuals extracted for the 12CO disc of MWC 480. Unphysical quadruple-like patterns observed at R = 300, 400 au and at R > 400 au in the velocity residuals resulting from the non-parametric surface suggest that a smooth (and less tapered) layer provides a better representation of the 12CO emission surface of this disc. |
Fig. B.16 Illustrating the influence of having one (top) or two (bottom) model intensity normalisations on peak intensity and velocity residuals extracted for the 12CO disc of AS 209. Although it is irrelevant for the retrieved kinematic signatures, we adopt the two-normalisation approach which aims to emulate absorption by foreground material occurring on the blueshifted side of the disc, allowing us to better capture fluctuations in the intensity field. |
Fig. B.17 As Fig. 13 but for the disc of AS 209. The planet marker indicates the orbital radius of the planet candidate associated with line width enhancements near the location of the CPD candidate proposed by Bae et al. (2022b) (see Table 2). |
Fig. B.18 As Fig. 13 but for the disc of IM Lup. The calculation of the velocity profiles illustrated in this Figure excludes +30° wedges around the (north and south) disc minor axis to avoid contamination from possibly unphysical velocity residuals (see Sect. 3.2, IM Lup). |
Fig. B.19 As Fig. 14 but for the disc of GM Aur. The calculation of velocity profiles illustrated in this Figure excludes +45° wedges around the (north and south) disc minor axis to minimise contamination from infalling material (see Sect. 3.2, GM Aur). |
Fig. B.20 Comparison of azimuthal and vertical velocity profiles extracted with the decomposition method introduced in Sect. 2.2.2 for the disc of MWC 480 from cubes with and without continuum emission. |
Appendix C Azimuthal average of absolute velocity residuals
In this section, we demonstrate that a radial profile of the azimuthal component of the velocity field, υφ, or Δυφ if referred to an axisymmetric model of the same component, can be obtained through the computation of azimuthal averages of absolute velocities, or absolute velocity residuals, respectively. The latter are defined as follows, (C.1)
Now, it is useful to note that under the reasonable assumption that the observed line-of-sight velocities are dominated by azimuthal velocities6, the absolute value of the data and model velocities can also be written as, (C.2)
for clockwise rotation and positive inclination. Using eq. C.2, an azimuthal average applied to eq. C.1 can thus be split into different independent terms,
from which the azimuthal component of the velocity is followed immediately, (C.3)
Finally, it is possible to generalise the above expression if the azimuthal section considered for the azimuthal average, with a total extent of ψ radians and not necessarily continuous, is symmetric around the disc major and minor axes, (C.4)
References
- Andrews, S. M. 2020, ARA&A, 58, 483 [Google Scholar]
- Andrews, S. M., Huang, J., Pérez, L. M., et al. 2018, ApJ, 869, L41 [NASA ADS] [CrossRef] [Google Scholar]
- Bae, J., Teague, R., & Zhu, Z. 2021, ApJ, 912, 56 [NASA ADS] [CrossRef] [Google Scholar]
- Bae, J., Isella, A., Zhu, Z., et al. 2022a, ArXiv e-prints [arXiv:2210.13314] [Google Scholar]
- Bae, J., Teague, R., Andrews, S. M., et al. 2022b, ApJ, 934, L20 [NASA ADS] [CrossRef] [Google Scholar]
- Barraza-Alfaro, M., Flock, M., Marino, S., & Pérez, S. 2021, A&A, 653, A113 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
- Beckwith, S. V. W., & Sargent, A. I. 1993, ApJ, 402, 280 [NASA ADS] [CrossRef] [Google Scholar]
- Benz, W., Ida, S., Alibert, Y., Lin, D., & Mordasini, C. 2014, in Protostars and Planets VI, eds. H. Beuther, R. S. Klessen, C. P. Dullemond, & T. Henning, 691 [Google Scholar]
- Bollati, F., Lodato, G., Price, D. J., & Pinte, C. 2021, MNRAS, 504, 5444 [Google Scholar]
- Booth, A. S., Walsh, C., Terwisscha van Scheltinga, J., et al. 2021, Nat. Astron., 5, 684 [NASA ADS] [CrossRef] [Google Scholar]
- Calcino, J., Hilder, T., Price, D. J., et al. 2022, ApJ, 929, L25 [NASA ADS] [CrossRef] [Google Scholar]
- Casassus, S., & Pérez, S. 2019, ApJ, 883, L41 [Google Scholar]
- Cuello, N., Ménard, F., & Price, D. J. 2023, EPJ, 138, 11 [NASA ADS] [Google Scholar]
- Czekala, I., Andrews, S. M., Jensen, E. L. N., et al. 2015, ApJ, 806, 154 [Google Scholar]
- Czekala, I., Andrews, S. M., Torres, G., et al. 2017, ApJ, 851, 132 [NASA ADS] [CrossRef] [Google Scholar]
- Czekala, I., Loomis, R. A., Teague, R., et al. 2021, ApJS, 257, 2 [NASA ADS] [CrossRef] [Google Scholar]
- Disk Dynamics Collaboration (Armitage, P. J., et al.) 2020, ArXiv e-prints [arXiv:2009.04345] [Google Scholar]
- Dong, R., Liu, S.-Y., & Fung, J. 2019, ApJ, 870, 72 [NASA ADS] [CrossRef] [Google Scholar]
- Dullemond, C. P., Isella, A., Andrews, S. M., Skobleva, I., & Dzyurkevich, N. 2020, A&A, 633, A137 [EDP Sciences] [Google Scholar]
- Dutrey, A., Semenov, D., Chapillon, E., et al. 2014, in Protostars and Planets VI, eds. H. Beuther, R. S. Klessen, C. P. Dullemond, & T. Henning, 317 [Google Scholar]
- Facchini, S., Juhász, A., & Lodato, G. 2018a, MNRAS, 473, 4459 [Google Scholar]
- Facchini, S., Pinilla, P., van Dishoeck, E. F., & de Juan Ovelar, M. 2018b, A&A, 612, A104 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
- Facchini, S., Teague, R., Bae, J., et al. 2021, AJ, 162, 99 [NASA ADS] [CrossRef] [Google Scholar]
- Fedele, D., Tazzari, M., Booth, R., et al. 2018, A&A, 610, A24 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
- Flaherty, K. M., Hughes, A. M., Rose, S. C., et al. 2017, ApJ, 843, 150 [NASA ADS] [CrossRef] [Google Scholar]
- Flaherty, K., Hughes, A. M., Simon, J. B., et al. 2020, ApJ, 895, 109 [NASA ADS] [CrossRef] [Google Scholar]
- Foreman-Mackey, D., Hogg, D. W., Lang, D., & Goodman, J. 2013, PASP, 125, 306 [Google Scholar]
- Galloway-Sprietsma, M., Bae, J., Teague, R., et al. 2023, ApJ, submitted [arXiv:2304.03665] [Google Scholar]
- Guzmán, V. V., Huang, J., Andrews, S. M., et al. 2018, ApJ, 869, L48 [Google Scholar]
- Guzmán, V. V., Bergner, J. B., Law, C. J., et al. 2021, ApJS, 257, 6 [CrossRef] [Google Scholar]
- Hacar, A., Alves, J., Burkert, A., & Goldsmith, P. 2016, A&A, 591, A104 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
- Huang, J., Andrews, S. M., Dullemond, C. P., et al. 2018a, ApJ, 869, L42 [NASA ADS] [CrossRef] [Google Scholar]
- Huang, J., Andrews, S. M., Pérez, L. M., et al. 2018b, ApJ, 869, L43 [Google Scholar]
- Huang, J., Andrews, S. M., Dullemond, C. P., et al. 2020, ApJ, 891, 48 [Google Scholar]
- Huang, J., Bergin, E. A., Öberg, K. I., et al. 2021, ApJS, 257, 19 [NASA ADS] [CrossRef] [Google Scholar]
- Ilee, J. D., Walsh, C., Booth, A. S., et al. 2021, ApJS, 257, 9 [NASA ADS] [CrossRef] [Google Scholar]
- Isella, A., Guidi, G., Testi, L., et al. 2016, Phys. Rev. Lett., 117, 251101 [NASA ADS] [CrossRef] [Google Scholar]
- Izquierdo, A. F., Testi, L., Facchini, S., Rosotti, G. P., & van Dishoeck, E. F. 2021, A&A, 650, A179 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
- Izquierdo, A. F., Facchini, S., Rosotti, G. P., van Dishoeck, E. F., & Testi, L. 2022, ApJ, 928, 2 [NASA ADS] [CrossRef] [Google Scholar]
- Jiang, H., Zhu, W., & Ormel, C. W. 2022, ApJ, 924, L31 [NASA ADS] [CrossRef] [Google Scholar]
- Kanagawa, K. D., Tanaka, H., Muto, T., Tanigawa, T., & Takeuchi, T. 2015, MNRAS, 448, 994 [NASA ADS] [CrossRef] [Google Scholar]
- Koch, E. W., & Rosolowsky, E. W. 2015, MNRAS, 452, 3435 [Google Scholar]
- Kratter, K., & Lodato, G. 2016, ARA&A, 54, 271 [NASA ADS] [CrossRef] [Google Scholar]
- Law, C. J., Loomis, R. A., Teague, R., et al. 2021a, ApJS, 257, 3 [NASA ADS] [CrossRef] [Google Scholar]
- Law, C. J., Teague, R., Loomis, R. A., et al. 2021b, ApJS, 257, 4 [NASA ADS] [CrossRef] [Google Scholar]
- Le Gal, R., Öberg, K. I., Teague, R., et al. 2021, ApJS, 257, 12 [NASA ADS] [CrossRef] [Google Scholar]
- Leemker, M., van’t Hoff, M. L. R., Trapman, L., et al. 2021, A&A, 646, A3 [EDP Sciences] [Google Scholar]
- Liu, S.-F., Jin, S., Li, S., Isella, A., & Li, H. 2018, ApJ, 857, 87 [NASA ADS] [CrossRef] [Google Scholar]
- Liu, Y., Dipierro, G., Ragusa, E., et al. 2019, A&A, 622, A75 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
- Lodato, G., Rampinelli, L., Viscardi, E., et al. 2022, MNRAS, 518, 4481 [NASA ADS] [CrossRef] [Google Scholar]
- Lubow, S. H., & Zhu, Z. 2014, ApJ, 785, 32 [NASA ADS] [CrossRef] [Google Scholar]
- Miotello, A., Kamp, I., Birnstiel, T., Cleeves, L. I., & Kataoka, A. 2022, ArXiv e-prints [arXiv:2203.09818] [Google Scholar]
- Morbidelli, A., Szulágyi, J., Crida, A., et al. 2014, Icarus, 232, 266 [Google Scholar]
- Öberg, K. I., & Bergin, E. A. 2021, Phys. Rep., 893, 1 [Google Scholar]
- Öberg, K. I., Qi, C., Fogel, J. K. J., et al. 2011, ApJ, 734, 98 [Google Scholar]
- Öberg, K. I., Guzmán, V. V., Walsh, C., et al. 2021, ApJS, 257, 1 [CrossRef] [Google Scholar]
- Paneque-Carreño, T., Miotello, A., van Dishoeck, E. F., et al. 2023, A&A, 669, A126 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
- Pérez, S., Casassus, S., & Benítez-Llambay, P. 2018, MNRAS, 480, L12 [Google Scholar]
- Perez, S., Dunhill, A., Casassus, S., et al. 2015, ApJ, 811, L5 [NASA ADS] [CrossRef] [Google Scholar]
- Pinte, C., Ménard, F., Duchêne, G., et al. 2018a, A&A, 609, A47 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
- Pinte, C., Price, D. J., Ménard, F., et al. 2018b, ApJ, 860, L13 [Google Scholar]
- Pinte, C., van der Plas, G., Ménard, F., et al. 2019, Nat. Astron., 3, 1109 [NASA ADS] [CrossRef] [Google Scholar]
- Pinte, C., Price, D. J., Ménard, F., et al. 2020, ApJ, 890, L9 [CrossRef] [Google Scholar]
- Pinte, C., Teague, R., Flaherty, K., et al. 2022, ArXiv e-prints [arXiv:2203.09528] [Google Scholar]
- Rab, C., Kamp, I., Dominik, C., et al. 2020, A&A, 642, A165 [EDP Sciences] [Google Scholar]
- Rabago, I., & Zhu, Z. 2021, MNRAS, 502, 5325 [NASA ADS] [CrossRef] [Google Scholar]
- Rosenfeld, K. A., Andrews, S. M., Hughes, A. M., Wilner, D. J., & Qi, C. 2013, ApJ, 774, 16 [Google Scholar]
- Rosotti, G. P., Benisty, M., Juhász, A., et al. 2020a, MNRAS, 491, 1335 [Google Scholar]
- Rosotti, G. P., Teague, R., Dullemond, C., Booth, R. A., & Clarke, C. J. 2020b, MNRAS, 495, 173 [Google Scholar]
- Schwarz, K. R., Calahan, J. K., Zhang, K., et al. 2021, ApJS, 257, 20 [NASA ADS] [CrossRef] [Google Scholar]
- Teague, R., & Foreman-Mackey, D. 2018, RNAAS, 2, 173 [NASA ADS] [Google Scholar]
- Teague, R., Guilloteau, S., Semenov, D., et al. 2016, A&A, 592, A49 [CrossRef] [EDP Sciences] [Google Scholar]
- Teague, R., Bae, J., Bergin, E. A., Birnstiel, T., & Foreman-Mackey, D. 2018a, ApJ, 860, L12 [NASA ADS] [CrossRef] [Google Scholar]
- Teague, R., Bae, J., Birnstiel, T., & Bergin, E. A. 2018b, ApJ, 868, 113 [Google Scholar]
- Teague, R., Bae, J., & Bergin, E. A. 2019a, Nature, 574, 378 [NASA ADS] [CrossRef] [Google Scholar]
- Teague, R., Bae, J., Huang, J., & Bergin, E. A. 2019b, ApJ, 884, L56 [Google Scholar]
- Teague, R., Bae, J., Aikawa, Y., et al. 2021, ApJS, 257, 18 [NASA ADS] [CrossRef] [Google Scholar]
- Teague, R., Bae, J., Andrews, S. M., et al. 2022, ApJ, 936, 163 [NASA ADS] [CrossRef] [Google Scholar]
- Testi, L., Birnstiel, T., Ricci, L., et al. 2014, in Protostars and Planets VI, eds. H. Beuther, R. S. Klessen, C. P. Dullemond, & T. Henning, 339 [Google Scholar]
- Turner, N. J., Choukroun, M., Castillo-Rogez, J., & Bryden, G. 2012, ApJ, 748, 92 [NASA ADS] [CrossRef] [Google Scholar]
- van Dishoeck, E. F., & Bergin, E. A. 2020, ArXiv e-prints [arXiv:2012.01472] [Google Scholar]
- van Zadelhoff, G. J., van Dishoeck, E. F., Thi, W. F., & Blake, G. A. 2001, A&A, 377, 566 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
- Veronesi, B., Paneque-Carreño, T., Lodato, G., et al. 2021, ApJ, 914, L27 [NASA ADS] [CrossRef] [Google Scholar]
- Verrios, H. J., Price, D. J., Pinte, C., Hilder, T., & Calcino, J. 2022, ApJ, 934, L11 [NASA ADS] [CrossRef] [Google Scholar]
- Walsh, C., Millar, T. J., & Nomura, H. 2010, ApJ, 722, 1607 [Google Scholar]
- Wölfer, L., Facchini, S., van der Marel, N., et al. 2023, A&A, 670, A154 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
- Yu, H., Teague, R., Bae, J., & Öberg, K. 2021, ApJ, 920, L33 [NASA ADS] [CrossRef] [Google Scholar]
- Zhang, K., Booth, A. S., Law, C. J., et al. 2021, ApJS, 257, 5 [NASA ADS] [CrossRef] [Google Scholar]
- Zhu, Z., Stone, J. M., & Rafikov, R. R. 2012, ApJ, 758, L42 [NASA ADS] [CrossRef] [Google Scholar]
Publicly available on https://alma-maps.info/data.html
Publicly available on https://github.com/andizq/discminer
Due to cloud absorption of the 12CO J = 2–1 emission (Öberg et al. 2011), the blueshifted half of the disc of AS 209, to the west of the sky plane, is dimmer than the redshifted side.
12CO traces layers between z/R ~ 0.1–0.2, whereas 13CO and C18O trace z/R < 0.1 scale heights (see Fig. B.5).
We note that the line width peak of both 13COandC18O at R = 100 au is rather subtle compared to that of 12CO at R = 120 au, and yet the decreasing velocity modulation is not negligible in either case. The reason is that the line width residuals of these optically thinner tracers are not as sensitive to density fluctuations as they are for 12CO (Hacar et al. 2016). Instead, local temperature and turbulence may also play a key role, making it harder to distinguish between the different contributions to the broadening of the line profile, unless detailed radiative transfer modelling is introduced.
Valid for intermediate and high inclinations where the sine is similar or much larger than the cosine of the disc inclination, and for almost all azimuthal angles φ, owing to the fact that azimuthal velocities, υφ,d, are normally between one and two orders of magnitude larger than radial and vertical motions, υr,d and υz,d.
All Tables
Location of potential planet-driven signatures in discs reported in this and in previous works based on gas velocity perturbations and substructures, assuming that the planet candidate is in the disc midplane.
Summary of radial substructures derived from azimuthally averaged residuals for the disc of MWC 480.
All Figures
Fig. 1 Illustrating the orientation and rotation direction of the discs observed by the ALMA large program MAPS and analysed here. Dashed arrows mark the origin of each disc coordinate system. For convenience, we assume this reference axis to be parallel to the projected semi-major axis that is closest to the north celestial axis in counterclockwise rotation. Also shown are the signatures reported in this work as potential planet-disc interaction features in the form of both coherent and turbulent velocity fluctuations. |
|
In the text |
Fig. 2 Selected CO isotopologue intensity channel maps from the disc around MWC 480, along with best-fit channels obtained with DISCMINER and their corresponding residuals. For comparison, isovelocity contours from the model upper and lower surfaces are overlaid on the data channels as solid and dashed lines, respectively. The ellipses in the background represent millimetre dust rings (solid) and gaps (dashed) as reported by Law et al. (2021a), projected on the upper emission surface of the corresponding model. Kink-like features apparent to the eye are marked with arrows. The radial location of the most prominent kink in 12CO is highlighted by the dotted purple ellipse. The beam size is shown in the bottom left corner of the top left panel for each isotopologue. |
|
In the text |
Fig. 3 Extraction of the line profile observables considered in this work, using MWC 480 12CO data as example. Peak intensities (top), line widths (middle), and centroid velocities (bottom) are displayed for both data and best-fit model, as well as residual maps showing the difference between them. Sketches on the rightmost column illustrate how each line profile observable is computed on the data cube and subsequently compared to those derived on the smooth, Keplerian model cube obtained by DISCMINER. The black solid lines are reference contours to ease comparison. Their levels are indicated in the colour bars. Residuals on each of the observables are driven by temperature, density and velocity fluctuations in the gas disc. |
|
In the text |
Fig. 4 Line width (top), velocity (middle), and peak intensity (bottom) residuals computed for the 12CO disc of MWC 480, as in Fig. 3, but deprojected on the disc reference frame. The vertical axis in this frame is parallel to the projected minor axis on the sky. The solid lines mark the radial location of millimetre dust rings, and the purple dotted line that of the most prominent kink identified by eye in 12CO channels. Also shown is the north-sky (or PA = 0°) axis for reference. The oblique shape of this axis is due to the deprojected geometry of the disc emission surface. |
|
In the text |
Fig. 5 Centroid velocity residuals traced by 12CO emission for all discs, referred to each disc reference frame. The annotated solid lines mark the radial location of millimetre dust rings (denoted as Bxxx), and the purple dotted line that of the most prominent kink identified by eye in 12CO intensity channels (denoted as Kxxx). Also shown is the north-sky (or PA = 0°) axis for reference. The oblique shape of this axis is due to the deprojected geometry of each disc vertical structure retrieved by the best-fit models. The circle in the bottom left corner of each panel indicates the location of the blue- and redshifted side for each disc. These maps illustrate the deviations from Keplerian rotation of the gas component located at high elevations (0.2 < z/r < 0.4) over the midplane. |
|
In the text |
Fig. 6 Analysis of localised velocity perturbations in the gas disc of HD 163296. Top row: folded velocity residuals computed for 12CO, 13CO, and C18O J = 2−1 lines. Middle row: azimuthal and radial extent of clusters identified by our detection algorithm. Regions in bold yellow highlight the accepted clusters of peak residuals whose variances are significantly greater than those in the background clusters. The weighted average of the 2D location of peak residuals within azimuthal and radial clusters determines the centre of the localised velocity perturbations reported in the text. These are marked as green crosses for 12CO and 13CO. In C18O, only radial clusters were identified as localised by the technique and, therefore, no perturbation was tagged as detected in this tracer. Bottom row: azimuthal profiles of folded velocities (i.e. those in the top row) to highlight the blue- redshifted, or Doppler flip, morphology of the perturbations detected in 12CO and 13CO, and of the tentative perturbation in C18O. |
|
In the text |
Fig. 7 Summary of the localised velocity perturbations detected by our clustering algorithm in the disc of HD 163296 in 12CO and 13CO velocities, using DSHARP and MAPS data. Localised red and blue contours highlight the morphology of the perturbations detected in 12CO, P94 and P261, whereas green and purple contours outline the localised signature detected in 13CO, P81, closely overlapping with the P94 signal. The green crosses mark the inferred location of the planet candidates linked to the aforementioned perturbations. Deep blue colours in the background highlight the location of 12CO line width minima which are co-spatial with positive gradients of azimuthal velocity flows probing gas surface density gaps. |
|
In the text |
Fig. 8 Folded line width residuals to highlight asymmetric line width enhancements near the orbit of the planet candidate P94, as well as the location of localised line width perturbations detected by our clustering algorithm in 13CO and C18O. |
|
In the text |
Fig. 9 Folded velocity residuals for the discs where no localised velocity perturbations are identified by our clustering algorithm. |
|
In the text |
Fig. 10 Azimuthal deprojection of velocity and line width residuals obtained for the disc of MWC 480 as observed in 12CO and 13CO. The horizontal dashed and solid lines indicate the radial location of dust gaps and rings, respectively. The purple dotted line is the radial location of kink-like features observed in the datacube channels. The vertical lines denote the azimuthal location of the disc main axes as projected on the sky. The black lines overlaid in the line width residuals highlight the location of extended line width enhancements, likely associated to non-thermal motions. The green cross marks the location of the planet candidate reported in Sect. 3.2 based on 12CO line width enhancements. Cartesian versions of these maps can be found in Figs. 4 and B.7. |
|
In the text |
Fig. 11 As Fig. 10 but for the disc of AS 209. The green circle marks the location of the CPD proposed by Bae et al. (2022b). Cartesian versions of these maps can be found in Figs. 5, B.7 and B.9. |
|
In the text |
Fig. 12 Illustrating isovelocity contours and the presence of coherent velocity substructures in the disc of AS 209 as probed by 12CO and 13CO J = 2−1 line emission. Black solid and dashed lines show isovelocity contours from the data and the best-fit Keplerian model, respectively. For reference, the isovelocity levels correspond to the same velocity channels of the intensity maps displayed around the main panels. Coherent sub- and super Keplerian substructures identified by FILFINDER (Koch & Rosolowsky 2015) in velocity residuals (Figs. 5 and B.7) are overlaid as blue and red thick lines. The green circle marks the location of the CPD candidate proposed by Bae et al. (2022b) in the disc reference frame. We note that the intensity kinks most prominently seen in the south of the sky in 12CO channels, between 3.5 and 5.5 km s−1, are not driven by localised velocity perturbations, but instead are the consequence of a large-scale blue-shifted flow that spans over the entire azimuth of the disc at around ~200 au, and appears at both scale-heights traced by 12CO and 13CO. In Sect. 4.1, we show that such a coherent signature is the result of upward vertical flows. |
|
In the text |
Fig. 13 Azimuthally averaged profiles of velocity (top), line width (middle) and peak intensity (bottom) residuals obtained for 12CO, 13CO, and C18O, for the disc around MWC 480. The radial location of millimetre dust gaps and rings is illustrated as dashed and solid lines, respectively. The radial distance of the most prominent kink apparent in 12CO channel maps is shown by the dotted purple line. Strong pressure bumps traced by minimal velocity gradients are marked with minuses in the top panel. Peak meridional flows hinting at gas moving away and towards the disc midplane are highlighted with arrows pointing up and down in the second top panel. The planet marker indicates the orbital radius of the planet candidate proposed at R = 245 au (see Table 2). |
|
In the text |
Fig. 14 As Fig. 13 but for the disc around HD 163296. The planet markers indicate the orbital radii of the planet candidates associated with the localised velocity perturbations P94 and P261 (see Table 2). |
|
In the text |
Fig. 15 Summary of azimuthally averaged residuals extrapolated onto a 2D grid to aid visualisation of (anti-)correlations between the different line profile observables extracted from the disc of MWC 480. The radial location of dust gaps and rings is illustrated as dashed and solid lines, respectively. The radial distance of the most prominent kink apparent in 12CO channel maps is marked by the dotted purple line. |
|
In the text |
Fig. 16 Illustrating the presence of non-axisymmetric features in the disc of MWC 480, as observed in 12CO. Left column: velocity, peak intensity, and line width residuals deprojected on the disc reference frame. Right column: azimuthally averaged residuals extracted per quadrant to highlight azimuthal variations in the different observables. The colour of each radial profile indicates the quadrant of the disc where it was computed according to the colouring code in the top right corner. Quadrants to the right of the disc vertical axis correspond to the redshifted side of the disc. |
|
In the text |
Fig. B.1 As Fig. 2 but for the disc of HD 163296. |
|
In the text |
Fig. B.2 As Fig. 2 but for the disc of AS 209. |
|
In the text |
Fig. B.3 As Fig. 2 but for the disc of IM Lup. |
|
In the text |
Fig. B.4 As Fig. 2 but for the disc of GM Aur. |
|
In the text |
Fig. B.5 Additional observables obtained from the DISCMINER analysis of the discs around MWC 480, HD 163296, AS 209, (IM Lup and GM Aur in next page) as traced by 12CO, 13CO, and C18O. Top row: Rotation curves computed from azimuthal averages of deprojected line-of-sight velocities. For reference, the grey dashed line in each panel is the Keplerian rotation curve of the best-fit model found for 12CO. Second row: Azimuthally averaged profiles of line widths, defined as the standard deviation of Gaussians fitted to the data line profiles. For reference, the purple dashed line in each panel represents the thermal broadening computed at the 13CO peak brightness temperature. Third row: Azimuthally averaged profiles of peak intensities converted to brightness temperatures using the Rayleigh-Jeans law. Bottom row: Elevation of upper and lower emission surfaces. Apart from intensity profiles, all panels in each row share the same x and y-axis extent. Vertical dashed, solid, and dotted lines mark the radial location of dust gaps, rings, and kinks observed in the channel maps of 12CO, respectively. |
|
In the text |
Fig. B.6 Schematic representation of the DISCMINER modelling flow and capabilities. DISCMINER produces best-fit channel maps of the input data cube by modelling the disc vertical structure and orientation, stellar mass, and attributes that shape the line profile morphology, assuming that the disc emission is smooth and Keplerian. The model and data channels are then compared through three types of moment maps. Large-scale analyses of these maps allow to look into the disc dynamical and thermodynamic structure, but are also useful to reveal potential signatures driven by planetary and non-planetary mechanisms, such as meridional flows or spirallike perturbations. Small-scale analysis of residuals from these maps leads to the retrieval of the potential location of embedded planets through the detection of localised velocity perturbations and sites of enhanced velocity dispersion. |
|
In the text |
Fig. B.7 As Fig. 5 but for 13CO. This isotopologue traces vertical layers between 0.1 < z/r < 0.15 in this sample of discs. |
|
In the text |
Fig. B.8 As Fig. 5 but for C18O. This isotopologue traces vertical layers between 0.0 < z/r < 0.1 in this sample of discs. |
|
In the text |
Fig. B.9 As Fig. 5 but for line width residuals. |
|
In the text |
Fig. B.10 As Fig. 5 but for line width residuals from 13CO. |
|
In the text |
Fig. B.11 As Fig. 5 but for line width residuals from C18O. |
|
In the text |
Fig. B.12 As Fig. 5 but for peak intensity residuals. |
|
In the text |
Fig. B.13 As Fig. 5 but for peak intensity residuals from 13CO. |
|
In the text |
Fig. B.14 As Fig. 5 but for peak intensity residuals from C18O. |
|
In the text |
Fig. B.15 Illustrating the impact of using an irregular (top) and a smooth parametric (bottom) emission surface on the peak intensity and velocity residuals extracted for the 12CO disc of MWC 480. Unphysical quadruple-like patterns observed at R = 300, 400 au and at R > 400 au in the velocity residuals resulting from the non-parametric surface suggest that a smooth (and less tapered) layer provides a better representation of the 12CO emission surface of this disc. |
|
In the text |
Fig. B.16 Illustrating the influence of having one (top) or two (bottom) model intensity normalisations on peak intensity and velocity residuals extracted for the 12CO disc of AS 209. Although it is irrelevant for the retrieved kinematic signatures, we adopt the two-normalisation approach which aims to emulate absorption by foreground material occurring on the blueshifted side of the disc, allowing us to better capture fluctuations in the intensity field. |
|
In the text |
Fig. B.17 As Fig. 13 but for the disc of AS 209. The planet marker indicates the orbital radius of the planet candidate associated with line width enhancements near the location of the CPD candidate proposed by Bae et al. (2022b) (see Table 2). |
|
In the text |
Fig. B.18 As Fig. 13 but for the disc of IM Lup. The calculation of the velocity profiles illustrated in this Figure excludes +30° wedges around the (north and south) disc minor axis to avoid contamination from possibly unphysical velocity residuals (see Sect. 3.2, IM Lup). |
|
In the text |
Fig. B.19 As Fig. 14 but for the disc of GM Aur. The calculation of velocity profiles illustrated in this Figure excludes +45° wedges around the (north and south) disc minor axis to minimise contamination from infalling material (see Sect. 3.2, GM Aur). |
|
In the text |
Fig. B.20 Comparison of azimuthal and vertical velocity profiles extracted with the decomposition method introduced in Sect. 2.2.2 for the disc of MWC 480 from cubes with and without continuum emission. |
|
In the text |
Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.
Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.
Initial download of the metrics may take a while.