Open Access
Issue
A&A
Volume 692, December 2024
Article Number A56
Number of page(s) 14
Section Planets, planetary systems, and small bodies
DOI https://doi.org/10.1051/0004-6361/202451145
Published online 02 December 2024

© The Authors 2024

Licence Creative CommonsOpen 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

Warped disks, whose orbital plane changes with radius, appear broadly in astrophysics. Theoretical predictions of this phenomenon result from a viscous disk evolving under external torques (Papaloizou & Pringle 1983; Pringle 1992). The ubiquity of warped disks spans diverse astronomical environments, as evidenced by observations ranging from galactic disks (e.g., Sánchez-Saavedra et al. 2003; Chen et al. 2019), X-ray binaries (e.g., Wijers & Pringle 1999; Begelman et al. 2006), active galactic nuclei (e.g., Herrnstein et al. 1996; Greenhill et al. 2003), Class 0 objects (IRAS 04368+2557, Sakai et al. 2018), Class 1 objects (L1489 IRS, Sai et al. 2020), and debris disks (e.g., Mouillet et al. 1997; Currie et al. 2012; Kasper et al. 2015). In protoplanetary disks, warps could be induced by a planet in a misaligned orbit (Nealon et al. 2018), or binary companions (Facchini et al. 2018). In the presence of magnetized stars, in which both the magnetic and rotational axes are tilted relative to the disk’s rotational axis, warped disks are expected to be formed (Romanova et al. 2021). Additionally, hydrodynam-ical simulations suggest that misalignments between the outer and inner regions of a primordial disk may arise after late infall events (Kuffmeier et al. 2021). Warps can also be generated by stellar flybys (Xiang-Gruess 2016; Cuello et al. 2019; Mayama et al. 2020), which produce recognizable observable signatures for these kinds of interacting systems (Cuello et al. 2020).

Observational signatures of warped geometries can be seen in near-infrared scattered light images as azimuthal dips in brightness (Marino et al. 2015). Likewise, nulls in scattered light observations can be produced by polarization self-cancellation (Weber et al. 2023), though this mechanism is only significant in observations with multiple illumination sources. These structures, commonly referred to as shadows, can be categorized into narrow and broad classes. Narrow shadows only extend a few degrees and result from a highly inclined inner disk. This feature was first observed, but not reported, in HD 142527 (Fukagawa et al. 2006), and later identified as intensity nulls (Casassus et al. 2012). Narrow shadows have since been reported in various systems, including HD100453 (Wagner et al. 2015; Benisty et al. 2017), HD 97048 (van der Plas et al. 2017), HD 135344B (Stolker et al. 2017), DoAr44 (Casassus et al. 2018), GG Tau (Keppler et al. 2020), and ZZ Tau 1RS (Hashimoto et al. 2024). In contrast, broad shadows obscure a significant portion of the disk and arise from a slightly inclined inner disk. These characteristics can be found in sources such as TW Hya (Debes et al. 2017) and HD 1430006 (Benisty et al. 2018). A remarkable case is HD 139614, where a shadow that spans more than 180° is attributed to two misaligned components (Muro-Arena et al. 2020). The location of the shadows and outer disk height place constraints on the inner disk orientation (Min et al. 2017). Evidence for disk tearing in a triple star system has been reported in GW Orionis (Kraus et al. 2020). Depending on the inner disk precession rate, these sources can exhibit variability. Multi-epoch observations of J1604 revealed variability of a narrow shadow (Pinilla et al. 2018). TW Hya is a unique case; observations of its shadows (Debes et al. 2017) are also reported in the molecular line integrated intensity, but with an azimuthal offset (Teague et al. 2022). Furthermore, a multi-epoch observation study showed that the single moving shadow evolved into two separate ones, suggesting that there are two misaligned precessing components (Debes et al. 2023).

In molecular line observations, warps appear as a characteristic twist in the velocity maps (Juhász & Facchini 2017). This distinct feature has been observed across various molecular tracers, including 12CO(J = 2−1) (Pérez et al. 2018), 12CO(J = 3−2) (Rosenfeld et al. 2014), HCO+(J = 3−2) (Loomis et al. 2017), HCO+(J = 4−3) (Rosenfeld et al. 2014), and 12CO(J = 6−5) (Casassus et al. 2015a), and it has been reported by multiple studies: TW Hya (Rosenfeld et al. 2012), HD142527 (Casassus et al. 2013; Casassus et al. 2015a; Garg et al. 2021), HD100546 (Walsh et al. 2017), AA Tau (Loomis et al. 2017), MWC 758 (Boehler et al. 2018), J1604 (Mayama et al. 2018), HD143006 (Pérez et al. 2018), HD 100453 AB (van der Plas et al. 2019), GW Orionis (Kraus et al. 2020; Bi et al. 2020), Elias 2–27 (Paneque-Carreño et al. 2021), IRS48 (van der Marel et al. 2021), CG Tau (Wölfer et al. 2021), and SY Cha (Orihara et al. 2023). Unfortunately, radial flows produce a similar signature in the velocity isophotes (Rosenfeld et al. 2014), and therefore it is not clear from molecular line observations which phenomenon is present or dominant within the disk. Distinguishing between these event configurations is a nontrivial task solely based on molecular line observations, and the complexity of interpreting these observational signatures stresses the need for complementary modeling approaches to disentangle the underlying physical processes.

High spatial and spectral resolution observations with the Atacama Large Millimeter/submillimeter Array (ALMA) have enabled mapping the kinematic structures of protoplanetary disks. Kinematic observations, commonly referred to as gas kinematics, have emerged as a valuable tool for studying the disk velocity structure and perturbations from a Keplerian velocity field. The widely used package Extracting Disk Dynamics eddy (Teague 2019) fits a Keplerian rotation model to the velocity maps under the assumption of axial symmetry. However, this approach encounters limitations in its assumptions when attempting to characterize disks with radial and azimuthal variations of the position angle, PA, and inclination, i, in the presence of a warp. An alternative package that is considered for these variations is ConeRot (Casassus 2022). It fits concentric annuli allowing for independent PA and i at each radius, although its implementation does not include a method to account for surface shadowing, which can lead to ambiguities, particularly in regions where the disk is self-obscured, compromising the accuracy of the inferred kinematical features. A more advanced 3D modeling method is DiscMiner (Izquierdo et al. 2021). This tool is specifically designed for fitting intensity channel maps, offering a more detailed analysis. However, the current implementation for the line profile kernel only allows the inclusion of Keplerian components, overlooking the influence of the warped structure on the observed line profiles.

In the context of gas kinematics, the common assumption of an azimuthally symmetric disk during model fitting introduces significant limitations, leading to inaccurate velocity field inferences. Young et al. (2022) suggests that evidence of kinematic warp may have been overlooked due to the conventional analysis methods. To this end, we present an extension of eddy that enables the fitting of a parametric warp. This add-on aims to overcome the limitations of assuming azimuthal symmetry, allowing for a more refined exploration of the kinematic structures within circumstellar disks. In this approach, we account self-shadowing effects in the generation of velocity centroids, overcoming the drawbacks of previous packages. To validate our tool, we constructed three types of radiative transfer models: a warped disk, a disk with a radial flow, and a model combining a warped disk with a radial flow. We explore the residual patterns arising after fitting a warped disk to each radiative model and discuss the particular observational features as a criteria to distinguish between cases.

Our work is presented as follows: in Sect. 2, we explain the prescription used for warping a coordinate system, the adopted disk model, the radiative transfer code setup, the process for creating synthetic observations, and the extension to eddy to account for warp curve fitting. The velocity maps, as well as the resulting residuals are presented in Sect. 3. In Sect. 4, we discuss the validity and limitations of the method, and in Sect. 5, we summarize and conclude.

2 Methods

2.1 Warp model and coordinates – sky description

A disk can be modeled as a series of concentric annuli (Rogstad et al. 1974; Abedi et al. 2014; Kamphuis et al. 2015), where each annulus is characterized by an inclination (i) and a position angle (PA). Establishing a suitable reference frame is essential for developing a model of a warped disk. In this work, we explore two descriptions of warped disks.

In a sky description (e.g., Bohn et al. 2022), we establish a coordinate system where the x-axis points north, the y-axis points east, and the ɀ-axis extends toward the observer at infinity. A warped disk can be constructed based on a parametrization of the normal vector describing each annulus. We define i(r)=asky (r,iin ,iout ,adr,ar0)$i\left( r \right) = {a_{{\rm{sky}}}}\left( {r,{i_{{\rm{in}}}},{i_{{\rm{out}}}},{a_{dr}},{a_{{r_0}}}} \right)$ as a rotation with respect to the x-axis of the sky and PA(r)=asky (r,PAin ,PAout,adr,ar0)${\rm{PA}}(r) = {a_{{\rm{sky}}}}\left( {r,{\rm{P}}{{\rm{A}}_{{\rm{in}}}},{\rm{P}}{{\rm{A}}_{{\rm{out}}}},{a_{dr}},{a_{{r_0}}}} \right)$ as a rotation around the z-axis, representing the observed major axis associated with each annulus. Here, r denotes the midplane radial coordinate. We define a radial profile that smoothly transitions from ain to aout across a distance of ain, given by asky (r,psky )=aout (aout ain )(1+exp(rar00.1adr )),${a_{{\rm{sky}}}}\left( {r,{p_{{\rm{sky}}}}} \right) = {a_{{\rm{out}}}} - {{\left( {{a_{{\rm{out}}}} - {a_{{\rm{in}}}}} \right)} \over {\left( {1 + \exp \left( {{{r - {a_{{r_0}}}} \over {0.1{a_{{\rm{dr}}}}}}} \right)} \right)}},$(1)

where psky ={ ain,aout ,adr,ar0 }${p_{{\rm{sky}}}} = \left\{ {{a_{{\rm{in}}}},{a_{{\rm{out}}}},{a_{dr}},{a_{{r_0}}}} \right\}$ is the parameter vector, ain is a parameter that sets the inner disk PA and i, aout sets the outer disk PA and i, ar0${a_{{r_0}}}$ specifies the inflection point of the curve, representing the midpoint of the transition, and adr controls the width of the curve, adjusting the rate of transition.

The coordinate transformation for these disks is described by the rotation matrices Rɀ(PA(r)) and Rx(i(r)), which account for the radial dependence of these angles. (xw,yw,ɀw)T=Rɀ(PA(r))Rx(i(r))(x,y,ɀ)T,${\left( {x_w^\prime ,y_w^\prime ,z_w^\prime } \right)^T} = {R_z}(PA(r)) \cdot {R_x}(i(r)) \cdot {(x,y,z)^T},$(2)

where · denotes matrix multiplication. The terms Rɀ(PA(r)) and Rx(i(r)) correspond to the rotation matrices around the ɀ- and x-axes, respectively. The transformed Cartesian coordinates in the warped frame are denoted as xw,yw, and ɀw. The rotation matrices are defined by Rɀ(θ)=[ cosθsinθ0sinθcosθ0001 ]${R_}(\theta ) = \left[ {\matrix{ {\cos \theta } & { - \sin \theta } & 0 \cr {\sin \theta } & {\cos \theta } & 0 \cr 0 & 0 & 1 \cr } } \right]$(3) Rx(θ)=[ 1000cosθsinθ0sinθcosθ ].${R_x}(\theta ) = \left[ {\matrix{ 1 & 0 & 0 \cr 0 & {\cos \theta } & { - \sin \theta } \cr 0 & {\sin \theta } & {\cos \theta } \cr } } \right].$(4)

thumbnail Fig. 1

Geometrical representation of the warp models for each prescription. Left: disk-frame representation of the warp. The inner disk orientation is characterized by the angle β, representing the inclination with respect to the outer disk normal vector, and the twist angle γ denoting the rotation around the normal vector. The direction of the observer is defined through a rotation around the outer disk ɀ-axis by ϕ, and inclination i away from it, the camera is rotated by PA around the vector pointing toward the observer. Right: sky-frame representation of the warp. Each ring has an associated i and PA. The inclination is defined with respect to the ɀ-axis. The PA rotation is applied around the line-of-sight represented by the black arrow. On the bottom the respective profiles are plotted. The colored dots represent each ring position, the inflection point is demarcated by the black dashed line, and the rectangular region shows the transition width. We note that the angles shown are representative of the highlighted red ring, the red arrow is the normal vector to the surface defined by that ring.

2.2 Warp model and coordinates – Disk description

Warps can also be described in relation to the outer disk, where the orientation of the inner disk is specified relative to the outer disk. We introduce the tilt angle β(r)=adisk (r,βin,adr,ar0)$\beta \left( r \right) = {a_{{\rm{disk}}}}\left( {r,{\beta _{{\rm{in}}}},{a_{dr}},{a_{{r_0}}}} \right)$, which represents the relative inclination between the inner and outer disk midplanes, and the twist angle γ(r)=adisk (r,γin,adr,ar0)$\gamma (r) = {a_{{\rm{disk}}}}\left( {r,{\gamma _{{\rm{in}}}},{a_{dr}},{a_{{r_0}}}} \right)$, denoting the rotation relative to normal vector of the outer disk. The associated parameter vector is pdisk ={ ain,adr,ar0 }.${p_{{\rm{disk}}}} = \left\{ {{{\rm{a}}_{{\rm{in}}}},{a_{dr}},{a_{{r_0}}}} \right\}.$. The function describing the radial profile in this case is given by adisk (r,pdisk )=ain (1+exp(rar00.1adr)),${a_{{\rm{disk}}}}\left( {r,{p_{{\rm{disk}}}}} \right) = {{{a_{{\rm{in}}}}} \over {\left( {1 + \exp \left( {{{r - {a_{{r_0}}}} \over {0.1{a_{{\rm{dr}}}}}}} \right)} \right)}},$(5)

where r is the radius of the unperturbed disk and ain is a parameter that sets the inner disk tilt or twist. We note that the curves share the inflection point ar0${a_{{r_0}}}$ and the transition width adr.

Two additional transformations, accounting for the outer disk PA and inclination, are applied. In this case, the warped coordinates are defined as follows: (xw,yw,ɀw)T=Rɀ(PA)Rx(inc)Rɀ(γ(r))Rx(β(r))(x,y,ɀ)T.${\left( {x_w^\prime ,y_w^\prime ,_w^\prime } \right)^T} = {R_}({\rm{PA}}) \cdot {R_x}({\rm{inc}}) \cdot {R_}(\gamma (r)) \cdot {R_x}(\beta (r)) \cdot {(x,y,)^T}.$(6)

A visual representation of a warped disk model for each description, along with their respective profiles, is presented in Fig. 1.

We note that these descriptions differ from one another. The sky description establishes orientation relative to the observer direction, while the disk description represents the orientation of the inner disk in terms of the outer disk. Each description is characterized by distinct angles, therefore, these descriptions are distinct. However, they can be represented within the other framework. The selection of the appropriate framework depends on the specific modeling requirements.

2.3 Disk model

To assess the efficacy of the proposed method, we applied it to simulated data that represent the ground truth our model aims to recover. In the following, we explain how we generate three distinct model types: a warped disk, an unperturbed disk with radial flow, and a warped disk with a radial flow. The vertical density structure is described by a Gaussian profile, defined as ρ(xu,yu,ɀu)=Σ(xu,yu)Hp(xu,yu)2πexp(ɀu22Hp(xu,yu)2),$\rho \left( {x_u^\prime ,y_u^\prime ,_u^\prime } \right) = {{\Sigma \left( {x_u^\prime ,y_u^\prime } \right)} \over {{H_{\rm{p}}}\left( {x_u^\prime ,y_u^\prime } \right)\sqrt {2\pi } }}\exp \left( { - {{_u^{\prime 2}} \over {2{H_{\rm{p}}}{{\left( {x_u^\prime ,y_u^\prime } \right)}^2}}}} \right),$(7)

where (xu,yu, ɀu) denote the unwarped coordinates: Σ(xu, yu) is the surface density profile, and Hp (xu, yu) represents the pressure scale height. We assume that these two distributions follow a power-law profile, which can be expressed as (xu,yu)=Σ0(R(xu,yu)r0)γ,$\mathop \sum \nolimits^ \left( {x_u^\prime ,y_u^\prime } \right) = {\Sigma _0}{\left( {{{R\left( {x_u^\prime ,y_u^\prime } \right)} \over {{r_0}}}} \right)^{ - \gamma }},$(8) Hp(xu,yu)=H0(R(xu,yu)r0)ζR(xu,yu),${H_{\rm{p}}}\left( {x_u^\prime ,y_u^\prime } \right) = {H_0}{\left( {{{R\left( {x_u^\prime ,y_u^\prime } \right)} \over {{r_0}}}} \right)^\zeta }R\left( {x_u^\prime ,y_u^\prime } \right),$(9)

where r0 is the reference radius, set to 1 au, while Σ0 and H0 represent the surface density and the aspect ratio at that radius, respectively. The gas surface density exponent is denoted by γ = 1, and the flaring index is ζ = 0.25. We did not incorporate an exponential taper in the surface density profile because it does not affect the regions we are interested in.

All models assume a Keplerian velocity field. For models that include a radial flow, we follow a similar approach as outlined in Rosenfeld et al. (2014), where an inward radial velocity component is introduced as a fraction χ (r) of the local Keplerian velocity. A value χ = 0 corresponds to a perfect circular Keplerian velocity, while χ = 1 indicates adding a radial inward velocity with a magnitude equal to the local Keplerian speed. In this case, the velocity field is vx=vKsinϕχ(r)vKcosϕ${\upsilon _x} = - {\upsilon _{\rm{K}}}\sin \phi - \chi (r){\upsilon _{\rm{K}}}\cos \phi $(10) vy=vKcosϕχ(r)vKsinϕ,${\upsilon _y} = {\upsilon _{\rm{K}}}\cos \phi - \chi (r){\upsilon _{\rm{K}}}\sin \phi ,$(11)

where ϕ is the azimuthal angle in the disk plane. We adopt the same functional form as given in Eq. (5) for the function χ (r), expressed as χ(r)=adisk (r,χ0,χr0,χdr)$\chi (r) = {a_{{\rm{disk}}}}\left( {r,{\chi _0},{\chi _{{r_0}}},{\chi _{dr}}} \right)$.

This function describes a smooth curve that decreases mono-tonically from the scaling factor χ0 to 0, radially outward. The net effect is that it introduces a radial velocity gradually, which generates an incremental twist in the isovelocity curves. We fixed the scaling factor χ0 = 1, following the procedure by Rosenfeld et al. (2014). However, the scaling factor can take any value within the range of 0 (indicating an absence of radial flows) to 2$\sqrt 2 $ (radial infall at the local free-fall velocity). We emphasize that the parameters χr0${\chi _{{r_0}}}$ and χdr are independent of those defining the warp. A radial flow has a similar visual appearance in the velocity maps to warps.

2.4 Radiative transfer simulations

We used the 3D radiative transfer simulation code RADMC3D1 (Dullemond et al. 2012) to model the molecular line emission from our prescribed disk configurations. A central star with a temperature of 5780 K, a radius of 1 R, and a mass of 1 M is placed at the grid center, with the system set at a distance of 140 pc. We used the mctherm task to compute the dust temperatures in the disk, assuming them equal to the gas temperature. Line radiative transfer calculations were conducted under the assumption of local thermodynamic equilibrium (LTE). The molecular data were obtained from the LAMDA database (Schöier et al. 2005). We included microscopic turbulence, which approximately follows αSScs$\sqrt {{\alpha _{{\rm{SS}}}}} {c_s}$, where cs is the sound speed and αSS is the viscosity parameter. We assumed αSS = 10−3. For the 12CO molecular abundance, we assumed the typical 10−4 fraction relative to H2 (Young & Scoville 1991; Lacy et al. 1994). In regions where T < 20 K, we set the CO abundance to 0 to mimic freeze-out. We did not include self-shielding of CO.

The models were implemented using spherical coordinates, where we initialized Nr = 100, Nθ = 180, and Nϕ = 180 grid points in the radial, azimuthal, and polar directions, respectively. The radial dimension of the mesh spans from 0.5 au to 100 au, using a logarithmic distribution. The azimuthal and polar directions are distributed linearly in the ranges [0, 2π] and [0, π], respectively. An additional refinement was included in the models that present a warped disk, adding 50 cells in the radial direction. These cells are linearly distributed within the range covering [ ar0adr,ar0+adr ]$\left[ {{a_{{r_0}}} - {a_{dr}},{a_{{r_0}}} + {a_{dr}}} \right]$, thus the radial grid accounts for the sudden change in orientation from the inner to the outer disk. In the case of a disk exhibiting only a radial flow, these cells are distributed within the range [ χr0χdr,χr0+χdr ]$\left[ {{\chi _{{r_0}}} - {\chi _{dr}},{\chi _{{r_0}}} + {\chi _{dr}}} \right]$.

We computed the channel maps for the 12CO(J = 2 − 1) line transition, centered at 230.538 GHz. The velocity range considered was ±15 km s−1, with a channel spacing of 40 m s−1, comparable to the maximum velocity resolution achieved by ALMA in Band 6.

To assess the observability of these structures, we generated synthetic observations from our models. We employed a two-step post-processing procedure to distinguish between the impacts of spatial resolution and thermal noise. In the first step, a circular Gaussian beam was convolved with each data cube. In the second step, a Gaussian noise image was created, convolved with the same beam, rescaled to the characteristic sensitivity of exoALMA, and finally added to the convolved images. We used a 0.03″ FWHM Gaussian beam for the test cases. This represents an optimistic scenario for the observable kinematic signatures achievable with ALMA. For the models shown in Sect. 3, we incorporated noise with a root mean square (RMS) of 3.0 K over 150 m s−1, similar to the exoALMA setup. Our analysis omits consideration of potential image artifacts arising from sparse uv coverage.

2.5 Warp extension for eddy

The eddy2 package is a tool designed for extracting kinematical information from spatially and spectrally resolved line data by fitting a Keplerian rotation pattern (Teague 2019). This package performs a fitting of the first moment map of an optically thick line by modeling the projected rotation velocity vϕ,proj(r, ϕ) = vϕ(r) cos ϕ sin i to recover the geometrical parameters of the disk and identify localized deviations in the velocity structure. This procedure assumes an axisymmetric disk structure. Hence, it does not account for a radial variation of the geometrical parameters. The package has been largely employed to determine the geometry of disks, as well as to identify localized deviations from the Keplerian field (see e.g., Galloway-Sprietsma et al. 2023; Ribas et al. 2023; Garg et al. 2022).

Our extension introduces the capability to incorporate both parametric warp prescriptions and a radial flow3, as described in Sect. 2.1. To address the shadowing problem that arises from projecting a 3D surface to 2D, we implemented an algorithm based on ray tracing. The algorithm takes two matrices with image dimensions nx and ny, providing the position of these pixels, a set of 2D meshes containing the (x′, y′, ɀ′) positions with corresponding warped surface values, and a 3D matrix containing the quantities for interpolation at pixel positions as an input. To define a surface, a triangulation is applied to the set of points defining the disk. This triangulation defines the vertices of the simplices as adjacent grid points. Then, for each image pixel, a ray is cast from infinity toward -ɀ, to search for an intersection with the triangles, and determine the closest intersected triangle to the image plane. For each selected surface, a barycentric interpolation is employed to approximate the object value at the image pixel position. This process enables the recovery of projected and interpolated values into a 2D image view, addressing the challenges posed by the projection of 3D structures onto 2D observational data. This represents a significant improvement over the shadowed approach already implemented in eddy, providing a more accurate and robust method for deprojecting pixel values describing complex 3D surfaces from observational data.

Table 1

Posterior distribution of best-fit models.

3 Results

This work extends the current capabilities of characterizing astrophysical disks using gas kinematics by incorporating warped geometries into eddy. For the upcoming results, we generated models using the sky prescription, with standard parameters set as PAout = 0° and iout = 25°, defining a prograde disk. We first introduced an inclination-only warp model, where only the inclination changes with radius while the position angle remains fixed. We set iin = 35°, allowing us to focus solely on the effects of inclination variations.

For warped and twisted disks, we adopted PAin = −45°, iin = 35°. For warped disks with radial flows, we used generic values of PAin = 10°, iin = 40°.

All models share common warp parameters, with a warp position set at ar0 = 40.0 au and a transition length of adr = 20.0 au. For disks with radial flow, we chose χr0=40.0${\chi _{{r_0}}} = 40.0$ au and χdr = 20.0 au. These values ensure similarity in the morphology of the velocity maps.

Rosenfeld et al. (2014) demonstrate that employing high angular and spectral resolution observations of molecular lines enables the detection of inflows, which manifest as twisted isophotes proximal to the systemic velocity. However, other mechanisms, such as warps, can produce similar patterns. To illustrate the inherent degeneracy in the signals between radial flows and warps, we show the channel maps of the three distinct types of models described in Sect. 2.3 in Fig. A.1. From now on, we refer to these models as Model W (warped disk model), Model F (disk with radial flows), and Model WF (warped disk model with radial flows). All models exhibit a similar overall physical shape in the first moment map, featuring a distinctive twisted emission pattern typical of sources with non-planar or radial motions.

3.1 Simple models

To generate velocity maps, we utilized the bettermoments4 (Teague & Foreman-Mackey 2018) package, employing a Gaussian fit. To validate whether the inclination-only warp, warp model, and radial flow implemented in eddy accurately recover the input warp and radial flow, we conducted a fit to each synthetic model velocity maps with its corresponding mode. In Fig. 2, we present the velocity map of the synthetic observations, alongside their respective best fits. The free parameters considered are: ΘI={ iin,PAin,adr,ar0 }${{\rm{\Theta }}_I} = \left\{ {{i_{{\rm{in}}}},{\rm{P}}{{\rm{A}}_{{\rm{in}}}},{a_{dr}},{a_{{r_0}}}} \right\}$ for the inclination-only warp, ΘW={ PAin,iin,adr,ar0 }${{\rm{\Theta }}_W} = \left\{ {{\rm{P}}{{\rm{A}}_{{\rm{in}}}},{i_{{\rm{in}}}},{a_{dr}},{a_{{r_0}}}} \right\}$ for the warp, and ΘF={ PAout,iout,χdr,χr0 }${{\rm{\Theta }}_F} = \left\{ {{\rm{P}}{{\rm{A}}_{{\rm{out}}}},{i_{{\rm{out}}}},{\chi _{dr}},{\chi _{{r_0}}}} \right\}$ for the radial flow. We restricted the free parameters to those describing the warp or radial flow, respectively, for computational efficiency, however we ensured that this decision does not affect the posterior distribution. We then sampled the posterior using an affine-invariant Markov Chain Monte Carlo (MCMC) (Goodman & Weare 2010) fitting routine implemented with emcee (Foreman-Mackey et al. 2013). We adopted 32 walkers, 1000 burn-in steps, and 1000 steps. The total number of steps ensures integrated autocorrelation time values between 20 and 90 for the chain of each parameter. Visual inspection of the chains confirms walker convergence. We then take N random draws from the sample, and average them to generate a model rotation map, denoted as vmodel. Across all synthetic models, we observe a twisted pattern in the velocity map with a smooth transition to the outer disk.

For Model I, the posterior distribution reveals an inner inclination of iin=34·327.7E04+7.7E04${i_{{\rm{in}}}} = 34_\cdot^ \circ 32_{ - 7.7{\rm{E}} - 04}^{ + 7.7{\rm{E}} - 04}$ and PAin=0·021.1E03+9.8E04${\rm{P}}{{\rm{A}}_{{\rm{in}}}} = 0_\cdot^ \circ 02_{ - 1.1{\rm{E}} - 03}^{ + 9.8{\rm{E}} - 04}$. The inflection point is found at ar0=0.291.1E05+1.2E05${a_{{r_0}}} = {0.29^{\prime \prime }}_{ - 1.1{\rm{E}} - 05}^{ + 1.2{\rm{E}} - 05}$ (equivalent to 40.96 au), and the transition width is adr=0.209.5E05+9.4E05${a_{dr}} = {0.20^{\prime \prime }}_{ - 9.5{\rm{E}} - 05}^{ + 9.4{\rm{E}} - 05}$ (27.67 au). These results are listed in Table 1. In our posterior distributions, we encounter uncertainties that are notably small, often on the order of 10−3 or smaller. This issue is discussed in Sect. 4.

For Model W, we obtain a posterior distribution spanning an inner position angle PAin=44·751.3E03+1.3E03${\rm{P}}{{\rm{A}}_{{\rm{in}}}} = - 44_\cdot^ \circ 75_{ - 1.3{\rm{E}} - 03}^{ + 1.3{\rm{E}} - 03}$, inner inclination iin =33·986.5E04+6.6E04${i_{{\rm{in}}}} = 33_\cdot^ \circ 98_{ - 6.5{\rm{E}} - 04}^{ + 6.6{\rm{E}} - 04}$, inflection point ar0=0.284.4E06+4.1E06${a_{{r_0}}} = {0.28^{\prime \prime }}_{ - 4.4{\rm{E}} - 06}^{ + 4.1{\rm{E}} - 06}$ (equivalent to 38.87 au), and transition width adr=0.183.0E05+2.9E05${a_{dr}} = {0.18^{\prime \prime }}_{ - 3.0{\rm{E}} - 05}^{ + 2.9{\rm{E}} - 05{\rm{E}}}$ (25.22au).

For the Model F, we obtain PAout =0·065.8E04+5.8E04${\rm{P}}{{\rm{A}}_{{\rm{out}}}} = - 0_\cdot^ \circ 06_{ - 5.8{\rm{E}} - 04}^{ + 5.8{\rm{E}} - 04}$, outer inclination iout =24.722.5E04+2.6E04${i_{{\rm{out}}}} = 24.72_{ - 2.5{\rm{E}} - 04}^{ + 2.6{\rm{E}} - 04}$, radial flow inflection point χr0=0.284.4E06+4.2E06${\chi _{{r_0}}} = {0.28^{\prime \prime }}_{ - 4.4{\rm{E}} - 06}^{ + 4.2{\rm{E}} - 06}$ (equivalent to 39.74 au), and transition width χdr=0.192.9E05+3.0E05${\chi _{dr}} = {0.19^{\prime \prime }}_{ - 2.9{\rm{E}} - 05}^{ + 3.0{\rm{E}} - 05}(26.42{\rm{au}})$ (26.42 au). The findings are summarized in Table 1.

The residuals reveal two main features (Figs. 2f, 2i): a central multipole attributed to high-velocity gas components outside the synthetic cube velocity range, and a double arc feature displaying sub- and super-Keplerian rotation at the transition length that emerges since in the radiative transfer model the emission comes from a emitting volume, while the sampler approximates this region as a surface. Minor scale deviations are attributed to injected thermal noise and the disk surface.

We repeated the analysis using the same synthetic data, but convolved with a 0.1″ FWHM beam, representing the angular resolution achieved by exoALMA. For Model I, we find an inner inclination iin =32·665.1E04+5.3E04${i_{{\rm{in}}}} = 32_\cdot^ \circ 66_{ - 5.1{\rm{E}} - 04}^{ + 5.3{\rm{E}} - 04}$, an inner position angle PAin =0·016.9E04+7.0E04${\rm{P}}{{\rm{A}}_{{\rm{in}}}} = 0_\cdot^ \circ 01_{ - 6.9{\rm{E}} - 04}^{ + 7.0{\rm{E}} - 04}$, an inflection point ar0=0.297.4E06+7.6E06${a_{{r_0}}} = {0.29^{\prime \prime }}_{ - 7.4{\rm{E}} - 06}^{ + 7.6{\rm{E}} - 06}$ (equivalent to 41.07 au), and a transition width adr=0.225.6E05+5.4E05${a_{dr}} = {0.22^{\prime \prime }}_{ - 5.6{\rm{E}} - 05}^{ + 5.4{\rm{E - 05}}}(30.75\,{\rm{au}})$ (30.75 au). For the Model W, we find PAin =43·191.2E03+1.3E03${\rm{P}}{{\rm{A}}_{{\rm{in}}}} = - 43_\cdot^ \circ 19_{ - 1.2{\rm{E}} - 03}^{ + 1.3{\rm{E}} - 03}$, an inner inclination iin =31·824.5E04+4.2E04${i_{{\rm{in}}}} = 31_\cdot^ \circ 82_{ - 4.5{\rm{E}} - 04}^{ + 4.2{\rm{E}} - 04}$, an inflection point ar0=0.283.5E06+3.5E06(38.72 au)${a_{{r_0}}} = {0.28^{\prime \prime }}_{ - 3.5{\rm{E}} - 06}^{ + 3.5{\rm{E}} - 06}\left( {38.72\,{\rm{au}}} \right)$, and transition width adr=0.232.1E05+2.2E05(32.69 au)${a_{dr}} = {0.23^{\prime \prime }}_{ - 2.1{\rm{E}} - 05}^{ + 2.2{\rm{E}} - 05}(32.69\,{\rm{au}})$. For Model F, we find PAout =0012.6E04+2.7E04${\rm{P}}{{\rm{A}}_{{\rm{out}}}} = - 0_ \bullet ^ \circ 01_{ - 2.6{\rm{E}} - 04}^{ + 2.7{\rm{E}} - 04}$, an outer inclination iout =24551.0E04+1.0E04${i_{{\rm{out}}}} = 24_ \bullet ^ \circ 55_{ - 1.0{\rm{E}} - 04}^{ + 1.0{\rm{E}} - 04}$, a radial flow inflection point χr0=0.283.0E06+2.8E06${\chi _{{r_0}}} = {0.28^{\prime \prime }}_{ - 3.0{\rm{E}} - 06}^{ + 2.8{\rm{E}} - 06}$ (equivalent to 38.93 au), and transition width χdr=0.272.1E05+2.0E05(37.68 au)${\chi _{dr}} = {0.27^{\prime \prime }}_{ - 2.1{\rm{E}} - 05}^{ + 2.0{\rm{E}} - 05}\left( {37.68\,{\rm{au}}} \right)$. The recovered posteriors exhibit a greater deviation from the input values compared to the 0.03″ FWHM case from the input, with parameters ar0${a_{{r_0}}}$ and χdr being overestimated by approximately 25% to 50% (see Table 1). Upon examining the residuals (see Fig. B.2), we observe a significant deviation from the background noise, indicating that our ability to characterize the kinematics of the system is limited by angular resolution.

thumbnail Fig. 2

Gallery of v0 velocity maps for Model I, W, and F, indicated in the top right corner. The left column displays the synthetic velocity maps, the middle column shows the respective best-fit, and the right column illustrates the residuals after subtracting a best-fit model. The beam size is shown in the lower left corner of each panel, and the black dashed line marks the radius at which the transition occurs.

thumbnail Fig. 3

Velocity map v0 for Model WF, along with its corresponding best-fit and residual maps. In each panel the beam size is indicated in the lower left corner. The black dashed line represents the radius at which the transition occurs.

3.2 Warp and radial flow differentiation

Both a warp and radial flow can coexist in protoplanetary disks, complicating modeling efforts further. Previous attempts to fit a warp-only disk to twisted kinematical data have not fully explained the observed kinematical features in some systems (Walsh et al. 2017; Loomis et al. 2017; Rivière-Marichalar et al. 2019). To address these limitations, we adopted a compound warp and radial flow model in our eddy analysis. We assessed the performance of eddy in recovering inputs from the combined warp and flow model (Model WF). The free parameters considered for this section are ΘWF={ PAin,iin,adr,ar0,χdr,χr0 }${{\rm{\Theta }}_{WF}} = \left\{ {{\rm{P}}{{\rm{A}}_{{\rm{in}}}},{i_{{\rm{in}}}},{a_{dr}},{a_{{r_0}}},{\chi _{dr}},{\chi _{{r_0}}}} \right\}$.

Fitting the compound model in eddy to Model WF, we found PAin=10·081.2E03+1.2E03${\rm{P}}{{\rm{A}}_{{\rm{in}}}} = 10_^ \circ 08_{ - 1.2{\rm{E}} - 03}^{ + 1.2{\rm{E}} - 03}$, an inner inclination iin =38·899.5E04+1.0E03${i_{{\rm{in}}}} = 38_\cdot^ \circ 89_{ - 9.5{\rm{E}} - 04}^{ + 1.0{\rm{E}} - 03}$, a transition width adr=0.186.1E05+6.4E05${a_{dr}} = {0.18^{\prime \prime }}_{ - 6.1{\rm{E}} - 05}^{ + 6.4{\rm{E}} - 05}$ (equivalent to 25.19 au), an inflection point ar0=0.281.2E03+1.2E03(39.71 au)${a_{{r_0}}} = {0.28^{\prime \prime }}_{ - 1.2{\rm{E}} - 03}^{ + 1.2{\rm{E}} - 03}(39.71\,{\rm{au}})$, a radial flow transition width χr0=0.294.0E06+3.6E06${\chi _{{r_0}}} = {0.29^{\prime \prime }}_{ - 4.0{\rm{E}} - 06}^{ + 3.6{\rm{E}} - 06}$ (equivalent to 40.83 au), and radial flow inflection point χr0=0.193.0E05+2.9E05(27.09 au)${\chi _{{r_0}}} = {0.19^{\prime \prime }}_{ - 3.0{\rm{E}} - 05}^{ + 2.9{\rm{E}} - 05}\left( {27.09\,{\rm{au}}} \right)$. In Fig. 3 we present a summary of these results. A similar multipole appears in the residuals (Fig. 3c). We observe that the residuals in the inner disk exceed the noise levels, with uncorrelated signals exhibiting non-informative kinematics.

Additionally, we evaluated the results when fitting this combined model to a warp-only configuration (Model W) and radial flow-only (Model F). This general approach converges to the specific case of having solely a warp (i.e., χr00${\chi _{{r_0}}}0$) in Model W or a radial flow (i.e., iiniout, and PAin ≈ PAout) in Model F. The sensitivity of the routine highlights its capability to discern between different disk configurations.

Another method to test for the presence of a warp in a system is to calculate the angle that defines the disk’s minor axis, considering only a radial flow. This angle, denoted as θsys = arctan (χ−1 cos i), is derived from Eq. (7) in Rosenfeld et al. (2014)5. The equation shows a dependency on both the radial velocity and the inclination angle χ. When gas moves at the local free-fall velocity, the parameters χ sets a maximum angle at which a radial flow can twist the disk, specifically when χ=2$\chi = \sqrt 2 $.

Any angle that defines the minor axis and exceeds this threshold indicates the presence of a warp (see Fig. 4).

thumbnail Fig. 4

Model velocity map v0 for a disk with a radial flow and a warp. The solid black line represents the minor axis of a disk with a radial flow with the local free-fall velocity, while the dashed gray dashed line represents the actual inner minor axis.

3.3 Diagnostics

Distinguishing between the kinematic signatures of a warp, radial flow, or both is complex due to their similar imprints in the velocity maps. Without prior information about the system, misinterpreting the residuals from a velocity map modeling procedure is possible. To overcome this visual ambiguity, we search for characteristic signatures in the residuals when fitting a warp or radial flow to an intrinsically different model. By examining these residuals, we aim to identify unique patterns or deviations that indicate the underlying system kinematics. This approach aids in interpreting observed residual features as indicative of a warp, radial flow, or both.

A summary of the v0 maps of the synthetic models is shown in the left column of Fig. 5. The right column displays the residual after subtracting a best-fit eddy model, with the type of eddy model indicated in the top right corner. For eddy modes indicating a warp, we consider the following free parameters: ΘW={ PAin,iin,adr,ar0 }${{\rm{\Theta }}_W} = \left\{ {{\rm{P}}{{\rm{A}}_{{\rm{in}}}},{i_{{\rm{in}}}},{a_{dr}},{a_{{r_0}}}} \right\}$. For the mode indicating a radial flow, we considered: ΘF={ χdr,χr0 }${{\rm{\Theta }}_F} = \left\{ {{\chi _{dr}},{\chi _{{r_0}}}} \right\}$.

In Fig. 5b, we observe the residuals after fitting a radial flow to a warped model. A strong quadrupole pattern is evident in the region defining inner disk. This systematic residual has a negative value toward inner disk major axis defining the red-shifted gas and a positive value toward the blue shifted part. Similarly, when fitting a warp to a model with a radial flow (Fig. 5d), we identify a quadrupole with inverted values to those seen in Fig 5b.

The kinematics traced by Model WF are shown in Fig. 5e and Fig. 5g. In Fig. 5f, a large perturbation defining a dipole is observed, matching the inner major axis. If we fit a warp (Fig. 5h), a quadrupole is observed with a similar alignment to what is seen in Fig. 5d.

Alongside these large deviations, two arcs at the boundary marking the transition from the inner to the outer region are identified. All residual patterns seem to be aligned with the inner disk major axis. The amplitude of the residual is significantly higher than the defined noise values.

These features consistently emerge across models with different geometries within the same model type. Despite the convergence of eddy, the posterior distributions characterizing the geometrical parameters are not consistent with the synthetic data, as the model describes another system configuration.

4 Discussion

Overall, a comprehensive analysis of the gas kinematics of a disk suspected of having warps, radial flows, or both, must be done adequately. The extension of the fitting routine enables the recovery of geometrical parameters only when data with sufficient angular resolution is available. If any of the signatures identified in Sect. 3.3 arises, the fitted geometrical parameters should not be used, as they not align with the actual disk parameters. In such cases, the residuals should only be considered as an indicative of the kinematics present in the disk.

The parameters obtained from our fitting procedure are consistent with the inputs of the synthetic radiative transfer models. However, two parameters, adr and χdr, exhibit a posterior distribution that deviates approximately 50% from the input model. This discrepancy likely arises because our model assumes that the emission is coming from an emitting surface rather than a realistic emitting volume, leading to inaccuracies in the region defining the transition from the inner to the outer disk, where the double arc residual feature is present. Nevertheless, this offset does not significantly alter the results, and the model remains valid for locations outside this region. The small uncertainties in the posterior distribution reflect the model underestimation of the data complexity, resulting in overly confident parameter estimates and consequently small uncertainties.

Foreground absorption is another source of error in model fitting. This phenomenon is common in nearby molecular clouds (van Kempen et al. 2009) and can significantly impact the analysis of the source. The presence of these foreground clouds distorts line profiles (van der Marel et al. 2013) and appears as reduced emission in the intensity maps. Such contamination can lead to inaccurate velocity maps. A documented example of this in a warped system is HD 142527 (Casassus et al. 2015b).

Our ray-casting tool facilitates the deprojection of any type of surface. In this work, we focus on warped disks, characterized by misaligned inner and outer regions. By adjusting the inclination and position angle profiles to represent the desired geometry, the tool can be adapted to investigate the gas kinematics of more exotic systems. This could include torn disks, such as GW Ori (Kraus et al. 2020), or systems with multiple misalignments, such as HD 139614 (Muro-Arena et al. 2020).

To evaluate the validity of the proposed method, we examined the effects of noise, spatial resolution, and inclination. In Fig. B.1, we display the corresponding residuals after incorporating a noise level with an RMS of 3.5 K over 350 m s−1, similar to the observational setup of the Molecules with ALMA at Planet-forming Scales (MAPS) large program Öberg et al. (2021). The results indicate that the contribution of the inclination i to the amplitude of projected velocity surpasses the effects of thermal noise, suggesting that the validity of the method is not limited by noise constraints.

In protoplanetary disks, warps have been reported at radii smaller than 20 au, posing challenges for current observational capabilities due to the associated angular sizes. In Fig. B.2, we show the corresponding residuals after subtracting the compatible best-fit model from the Model W and Model F data convolved with 0.1″ FWHM beam. The posterior distributions recovered values comparable to those obtained in Sect. 3.1. However, the residuals are too large to confidently categorize the model type, leading to the potential misinterpretation of the features as fitting the wrong model type to the data (see right column of Fig. 5). We found that for an angular resolution finer than 0.04″, this effect is minimized, while for larger resolutions, accurately discerning and characterizing these features becomes challenging. The available spatial resolution impacts the correct characterization of warps, and despite advancements in highresolution imaging techniques, the small radial scale warps are be susceptible to beam dilution. This limitation increases the complexities involved in observing and characterizing warped disks.

We investigated the range of inclinations for which the method remains applicable. Specifically, for the warped disk scenario, we verified values of iin = (15°, 25°, 50°), while for the radial flow case, we tested for an inclination of iout = (15°, 35°, 50°). All other parameters remained consistent with those outlined in Section 3. A summary of the residual maps is shown in Fig. B.3. We observed an increase in the amplitude of the residuals with inclination. We suggest that this method works optimally for inclinations below 40°.

thumbnail Fig. 5

Displaying v0 galleries for various model types. The left column exhibits synthetic v0 maps, with model types labeled in the top right corner. The right column displays residuals after subtraction of the best-fit model, also labeled in the top right corner. Beam size is denoted in the lower left corner of each peach panel. The black dashed line represents the radius at which the transition occurs.

5 Conclusion

The prevalence of misalignments challenges the assumption of planar disks, requiring methods capable of accounting these more complex geometries. Molecular line observations indicate that misalignments can mimic signals attributed to radial flows. Through comprehensive warped kinematics and radiative transfer simulations, we have investigated whether this degeneracy can be broken. This study advances the characterization of astrophysical disks, particularly those exhibiting warped structures. Applying this tool to actual observational data, presents an opportunity to revisit and reassess warped systems that may have been misinterpreted, as discussed by Young et al. (2022). The tool capability to discern between warps and other kinematic features can lead to a more nuanced interpretation of observed structures.

Our main findings can be summarized as follows:

  1. Gas kinematics in protoplanetary disks are complex, particularly when warps and radial flows are involved, posing challenges for accurate modeling and interpretation. The eddy modeling tool, extended to include warped geometries and radial flows, provides a valuable framework for investigating these types of disks;

  2. Noise levels comparable to previous ALMA observations do not affect the capability to discern between signatures, because the influence of the inclination i and radial flows on the data surpasses the impact of thermal noise. The accurate characterization of warps in protoplanetary disks is limited by angular resolution;

  3. While degeneracies between different disk configurations, such as warped-only and radial flow models, complicate the interpretation of observed data, signatures that emerge in the residuals from fitting an incorrect model to a system that exhibits twisted kinematics can be used as a criteria to distinguish between the warped disk radial flow dilemma;

  4. Future observational efforts, particularly high-resolution ALMA observations, hold promise for refining our understanding of disk kinematics and evolution. Advancements in modeling techniques and observational capabilities are crucial for overcoming the limitations and complexities associated with characterizing gas kinematics in protoplanetary disks.

We also discuss the current limitations that observational studies of warps face. Our methodology establishes a foundation for future investigations and holds promising applications in revisiting and reinterpreting signals from potentially misaligned disks.

Acknowledgements

We sincerely thank the anonymous referee for their valuable suggestions, which significantly enhanced the quality of our paper. T.B. acknowledges funding from the European Union under the European Union’s Horizon 2020 research and innovation program under grant agreement no. 714769 and under the Horizon Europe Research and Innovation Program 101124282 (EARLYBIRD) as well as funding by the Deutsche Forschungsgemeinschaft (DFG, German Research Foundation) under grant 325594231, and Germany’s Excellence Strategy – EXC-2094 – 390783311. Views and opinions expressed are, however, those of the authors only and do not necessarily reflect those of the European Union or the European Research Council. Neither the European Union nor the granting authority can be held responsible for them. Software: NumPy (Harris et al. 2020), emcee (Foreman-Mackey et al. 2013), Matplotlib (Hunter 2007).

Appendix A Channel maps

thumbnail Fig. A.1

Synthetic CO J=2−1 channel maps, imaged at 40 m s−1 spacing for each kind of model. The beam size is shown in the bottom left of the panels. The model type is indicated in the upper left: the top row shows the warped model, middle row shows the unperturbed disk with a radial flow, and the bottom panel shows a disk with a warp and a radial flow. We only show the channels close to the systemic velocity vLSR = 0 k m s−1.

Appendix B Extra figures

thumbnail Fig. B.1

Same as Fig. 2 but with a noise RMS of 3.5 K over 350 m s−1.

thumbnail Fig. B.2

Same as Fig. 2 but the synthetic radiative transfer models are convolved with a 0.1″ FWHM beam

thumbnail Fig. B.3

Residual maps after subtraction of their respective best-fit models. The top row shows the residual for Model W with inclinations iin = (15°, 25°, 50°). Bottom row shows the residual for Model F with inclinations iout = (15°, 35°, 50°). The respective inclination it shown in the bottom right corner, and the beam size is shown in the bottom left corner.

References

  1. Abedi, H., Mateu, C., Aguilar, L. A., Figueras, F., & Romero-Gómez, M. 2014, MNRAS, 442, 3627 [NASA ADS] [CrossRef] [Google Scholar]
  2. Begelman, M. C., King, A. R., & Pringle, J. E. 2006, MNRAS, 370, 399 [Google Scholar]
  3. Benisty, M., Stolker, T., Pohl, A., et al. 2017, A&A, 597, A42 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  4. Benisty, M., Juhász, A., Facchini, S., et al. 2018, A&A, 619, A171 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  5. Bi, J., van der Marel, N., Dong, R., et al. 2020, ApJ, 895, L18 [Google Scholar]
  6. Boehler, Y., Ricci, L., Weaver, E., et al. 2018, ApJ, 853, 162 [Google Scholar]
  7. Bohn, A. J., Benisty, M., Perraut, K., et al. 2022, A&A, 658, A183 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  8. Casassus, S. 2022, ConeRot: Velocity perturbations extractor, Astrophysics Source Code Library [record ascl:2207.027] [Google Scholar]
  9. Casassus, S., Perez M., S., Jordán, A., et al. 2012, ApJ, 754, L31 [NASA ADS] [CrossRef] [Google Scholar]
  10. Casassus, S., van der Plas, G. M., Perez, S., et al. 2013, Nature, 493, 191 [Google Scholar]
  11. Casassus, S., Marino, S., Pérez, S., et al. 2015a, ApJ, 811, 92 [Google Scholar]
  12. Casassus, S., Wright, C. M., Marino, S., et al. 2015b, ApJ, 812, 126 [NASA ADS] [CrossRef] [Google Scholar]
  13. Casassus, S., Avenhaus, H., Pérez, S., et al. 2018, MNRAS, 477, 5104 [NASA ADS] [CrossRef] [Google Scholar]
  14. Chen, X., Wang, S., Deng, L., et al. 2019, Nat. Astron., 3, 320 [NASA ADS] [CrossRef] [Google Scholar]
  15. Cuello, N., Dipierro, G., Mentiplay, D., et al. 2019, MNRAS, 483, 4114 [Google Scholar]
  16. Cuello, N., Louvet, F., Mentiplay, D., et al. 2020, MNRAS, 491, 504 [Google Scholar]
  17. Currie, T., Rodigas, T. J., Debes, J., et al. 2012, ApJ, 757, 28 [NASA ADS] [CrossRef] [Google Scholar]
  18. Debes, J. H., Poteet, C. A., Jang-Condell, H., et al. 2017, ApJ, 835, 205 [NASA ADS] [CrossRef] [Google Scholar]
  19. Debes, J., Nealon, R., Alexander, R., et al. 2023, ApJ, 948, 36 [NASA ADS] [CrossRef] [Google Scholar]
  20. Dullemond, C. P., Juhasz, A., Pohl, A., et al. 2012, RADMC-3D: A multipurpose radiative transfer tool, Astrophysics Source Code Library [record ascl:1202.015] [Google Scholar]
  21. Facchini, S., Juhász, A., & Lodato, G. 2018, MNRAS, 473, 4459 [Google Scholar]
  22. Foreman-Mackey, D., Hogg, D. W., Lang, D., & Goodman, J. 2013, PASP, 125, 306 [Google Scholar]
  23. Fukagawa, M., Tamura, M., Itoh, Y., et al. 2006, ApJ, 636, L153 [Google Scholar]
  24. Galloway-Sprietsma, M., Bae, J., Teague, R., et al. 2023, ApJ, 950, 147 [NASA ADS] [CrossRef] [Google Scholar]
  25. Garg, H., Pinte, C., Christiaens, V., et al. 2021, MNRAS, 504, 782 [Google Scholar]
  26. Garg, H., Pinte, C., Hammond, I., et al. 2022, MNRAS, 517, 5942 [CrossRef] [Google Scholar]
  27. Goodman, J., & Weare, J. 2010, Commun. Appl. Math. Computat. Sci., 5, 65 [Google Scholar]
  28. Greenhill, L. J., Booth, R. S., Ellingsen, S. P., et al. 2003, ApJ, 590, 162 [NASA ADS] [CrossRef] [Google Scholar]
  29. Harris, C. R., Millman, K. J., van der Walt, S. J., et al. 2020, Nature, 585, 357 [NASA ADS] [CrossRef] [Google Scholar]
  30. Hashimoto, J., Dong, R., Muto, T., Liu, H. B., & Terada, Y. 2024, arXiv e-prints [arXiv:2401.02004] [Google Scholar]
  31. Herrnstein, J. R., Greenhill, L. J., & Moran, J. M. 1996, ApJ, 468, L17 [CrossRef] [Google Scholar]
  32. Hunter, J. D. 2007, Comput. Sci. Eng., 9, 90 [NASA ADS] [CrossRef] [Google Scholar]
  33. 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]
  34. Juhász, A., & Facchini, S. 2017, MNRAS, 466, 4053 [NASA ADS] [Google Scholar]
  35. Kamphuis, P., Józsa, G. I. G., Oh, S. H., et al. 2015, MNRAS, 452, 3139 [NASA ADS] [CrossRef] [Google Scholar]
  36. Kasper, M., Apai, D., Wagner, K., & Robberto, M. 2015, ApJ, 812, L33 [NASA ADS] [CrossRef] [Google Scholar]
  37. Keppler, M., Penzlin, A., Benisty, M., et al. 2020, A&A, 639, A62 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  38. Kraus, S., Kreplin, A., Young, A. K., et al. 2020, Science, 369, 1233 [NASA ADS] [CrossRef] [Google Scholar]
  39. Kuffmeier, M., Dullemond, C. P., Reissl, S., & Goicovic, F. G. 2021, A&A, 656, A161 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  40. Lacy, J. H., Knacke, R., Geballe, T. R., & Tokunaga, A. T. 1994, ApJ, 428, L69 [CrossRef] [Google Scholar]
  41. Loomis, R. A., Öberg, K. I., Andrews, S. M., & MacGregor, M. A. 2017, ApJ, 840, 23 [NASA ADS] [CrossRef] [Google Scholar]
  42. Marino, S., Perez, S., & Casassus, S. 2015, ApJ, 798, L44 [Google Scholar]
  43. Mayama, S., Akiyama, E., Panic´, O., et al. 2018, ApJ, 868, L3 [NASA ADS] [CrossRef] [Google Scholar]
  44. Mayama, S., Pérez, S., Kusakabe, N., et al. 2020, AJ, 159, 12 [NASA ADS] [CrossRef] [Google Scholar]
  45. Min, M., Stolker, T., Dominik, C., & Benisty, M. 2017, A&A, 604, L10 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  46. Mouillet, D., Larwood, J. D., Papaloizou, J. C. B., & Lagrange, A. M. 1997, MNRAS, 292, 896 [Google Scholar]
  47. Muro-Arena, G. A., Benisty, M., Ginski, C., et al. 2020, A&A, 635, A121 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  48. Nealon, R., Dipierro, G., Alexander, R., Martin, R. G., & Nixon, C. 2018, MNRAS, 481, 20 [Google Scholar]
  49. Öberg, K. I., Guzmán, V. V., Walsh, C., et al. 2021, ApJS, 257, 1 [CrossRef] [Google Scholar]
  50. Orihara, R., Momose, M., Muto, T., et al. 2023, PASJ, 75, 424 [CrossRef] [Google Scholar]
  51. Paneque-Carreño, T., Pérez, L. M., Benisty, M., et al. 2021, ApJ, 914, 88 [CrossRef] [Google Scholar]
  52. Papaloizou, J. C. B., & Pringle, J. E. 1983, MNRAS, 202, 1181 [NASA ADS] [Google Scholar]
  53. Pérez, L. M., Benisty, M., Andrews, S. M., et al. 2018, ApJ, 869, L50 [CrossRef] [Google Scholar]
  54. Pinilla, P., Benisty, M., de Boer, J., et al. 2018, ApJ, 868, 85 [NASA ADS] [CrossRef] [Google Scholar]
  55. Pringle, J. E. 1992, MNRAS, 258, 811 [NASA ADS] [CrossRef] [Google Scholar]
  56. Ribas, Á., Macías, E., Weber, P., et al. 2023, A&A, 673, A77 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  57. Rivière-Marichalar, P., Fuente, A., Baruteau, C., et al. 2019, ApJ, 879, L14 [CrossRef] [Google Scholar]
  58. Rogstad, D. H., Lockhart, I. A., & Wright, M. C. H. 1974, ApJ, 193, 309 [Google Scholar]
  59. Romanova, M. M., Koldoba, A. V., Ustyugova, G. V., et al. 2021, MNRAS, 506, 372 [NASA ADS] [CrossRef] [Google Scholar]
  60. Rosenfeld, K. A., Qi, C., Andrews, S. M., et al. 2012, ApJ, 757, 129 [NASA ADS] [CrossRef] [Google Scholar]
  61. Rosenfeld, K. A., Chiang, E., & Andrews, S. M. 2014, ApJ, 782, 62 [Google Scholar]
  62. Sai, J., Ohashi, N., Saigo, K., et al. 2020, ApJ, 893, 51 [NASA ADS] [CrossRef] [Google Scholar]
  63. Sakai, N., Hanawa, T., Zhang, Y., et al. 2018, Nature, 565, 206 [Google Scholar]
  64. Sánchez-Saavedra, M. L., Battaner, E., Guijarro, A., López-Corredoira, M., & Castro-Rodríguez, N. 2003, A&A, 399, 457 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  65. Schöier, F. L., van der Tak, F. F. S., van Dishoeck, E. F., & Black, J. H. 2005, A&A, 432, 369 [Google Scholar]
  66. Stolker, T., Sitko, M., Lazareff, B., et al. 2017, ApJ, 849, 143 [Google Scholar]
  67. Teague, R. 2019, J. Open Source Softw., 4, 1220 [Google Scholar]
  68. Teague, R., & Foreman-Mackey, D. 2018, bettermoments: A robust method to measure line centroids, https://zenodo.org/records/1419754 [Google Scholar]
  69. Teague, R., Bae, J., Benisty, M., et al. 2022, ApJ, 930, 144 [NASA ADS] [CrossRef] [Google Scholar]
  70. van der Marel, N., van Dishoeck, E. F., Bruderer, S., et al. 2013, Science, 340, 1199 [Google Scholar]
  71. van der Marel, N., Birnstiel, T., Garufi, A., et al. 2021, AJ, 161, 33 [Google Scholar]
  72. van der Plas, G., Wright, C. M., Ménard, F., et al. 2017, A&A, 597, A32 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  73. van der Plas, G., Ménard, F., Gonzalez, J. F., et al. 2019, A&A, 624, A33 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  74. van Kempen, T. A., van Dishoeck, E. F., Salter, D. M., et al. 2009, A&A, 498, 167 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  75. Wagner, K., Apai, D., Kasper, M., & Robberto, M. 2015, ApJ, 813, L2 [Google Scholar]
  76. Walsh, C., Daley, C., Facchini, S., & Juhász, A. 2017, A&A, 607, A114 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  77. Weber, P., Pérez, S., Guidi, G., et al. 2023, MNRAS, 518, 5620 [Google Scholar]
  78. Wijers, R. A. M. J., & Pringle, J. E. 1999, MNRAS, 308, 207 [NASA ADS] [CrossRef] [Google Scholar]
  79. Wölfer, L., Facchini, S., Kurtovic, N. T., et al. 2021, A&A, 648, A19 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  80. Xiang-Gruess, M. 2016, MNRAS, 455, 3086 [NASA ADS] [CrossRef] [Google Scholar]
  81. Young, J. S., & Scoville, N. Z. 1991, ARA&A, 29, 581 [Google Scholar]
  82. Young, A. K., Alexander, R., Rosotti, G., & Pinte, C. 2022, MNRAS, 513, 487 [NASA ADS] [CrossRef] [Google Scholar]

5

We note that in this derivation, the minor axis points north, rather than the east, as defined in our work.

All Tables

Table 1

Posterior distribution of best-fit models.

All Figures

thumbnail Fig. 1

Geometrical representation of the warp models for each prescription. Left: disk-frame representation of the warp. The inner disk orientation is characterized by the angle β, representing the inclination with respect to the outer disk normal vector, and the twist angle γ denoting the rotation around the normal vector. The direction of the observer is defined through a rotation around the outer disk ɀ-axis by ϕ, and inclination i away from it, the camera is rotated by PA around the vector pointing toward the observer. Right: sky-frame representation of the warp. Each ring has an associated i and PA. The inclination is defined with respect to the ɀ-axis. The PA rotation is applied around the line-of-sight represented by the black arrow. On the bottom the respective profiles are plotted. The colored dots represent each ring position, the inflection point is demarcated by the black dashed line, and the rectangular region shows the transition width. We note that the angles shown are representative of the highlighted red ring, the red arrow is the normal vector to the surface defined by that ring.

In the text
thumbnail Fig. 2

Gallery of v0 velocity maps for Model I, W, and F, indicated in the top right corner. The left column displays the synthetic velocity maps, the middle column shows the respective best-fit, and the right column illustrates the residuals after subtracting a best-fit model. The beam size is shown in the lower left corner of each panel, and the black dashed line marks the radius at which the transition occurs.

In the text
thumbnail Fig. 3

Velocity map v0 for Model WF, along with its corresponding best-fit and residual maps. In each panel the beam size is indicated in the lower left corner. The black dashed line represents the radius at which the transition occurs.

In the text
thumbnail Fig. 4

Model velocity map v0 for a disk with a radial flow and a warp. The solid black line represents the minor axis of a disk with a radial flow with the local free-fall velocity, while the dashed gray dashed line represents the actual inner minor axis.

In the text
thumbnail Fig. 5

Displaying v0 galleries for various model types. The left column exhibits synthetic v0 maps, with model types labeled in the top right corner. The right column displays residuals after subtraction of the best-fit model, also labeled in the top right corner. Beam size is denoted in the lower left corner of each peach panel. The black dashed line represents the radius at which the transition occurs.

In the text
thumbnail Fig. A.1

Synthetic CO J=2−1 channel maps, imaged at 40 m s−1 spacing for each kind of model. The beam size is shown in the bottom left of the panels. The model type is indicated in the upper left: the top row shows the warped model, middle row shows the unperturbed disk with a radial flow, and the bottom panel shows a disk with a warp and a radial flow. We only show the channels close to the systemic velocity vLSR = 0 k m s−1.

In the text
thumbnail Fig. B.1

Same as Fig. 2 but with a noise RMS of 3.5 K over 350 m s−1.

In the text
thumbnail Fig. B.2

Same as Fig. 2 but the synthetic radiative transfer models are convolved with a 0.1″ FWHM beam

In the text
thumbnail Fig. B.3

Residual maps after subtraction of their respective best-fit models. The top row shows the residual for Model W with inclinations iin = (15°, 25°, 50°). Bottom row shows the residual for Model F with inclinations iout = (15°, 35°, 50°). The respective inclination it shown in the bottom right corner, and the beam size is shown in the bottom left corner.

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.