Open Access
Issue
A&A
Volume 670, February 2023
Article Number L1
Number of page(s) 13
Section Letters to the Editor
DOI https://doi.org/10.1051/0004-6361/202245381
Published online 27 January 2023

© The Authors 2023

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

Planet formation appears to be a robust and efficient process, occurring around both single and multiple stellar systems (Kostov et al. 2016) in protoplanetary disks. The advent of high resolution imaging facilities has demonstrated that nearly all bright and extended disks show substructures, in particular in the small (micron-sized) and large (millimeter-sized) dust tracers seen through scattered and thermal light, respectively (e.g., Andrews et al. 2018; Long et al. 2018; Rich et al. 2022; Benisty et al. 2022; Bae et al. 2022a). Such high resolution studies applied to gas tracers allow the overall physical conditions in the disk to be probed, such as its temperature structure, its surface height (Rich et al. 2021; Law et al. 2021), and pressure variations (Teague et al. 2018a,b; Rosotti et al. 2020). Studies of the disk density and the velocity structure reveal a great complexity, including localized non-Keplerian features that can be attributed to embedded massive protoplanets (e.g., Pinte et al. 2022; Wölfer et al. 2023). Such perturbations from smooth density and velocity distributions can directly constrain planet formation since they are expected to leave clear signatures on the disk structure (e.g., Perez et al. 2015; Yun et al. 2019). For example, the mapping of spiral wakes (Calcino et al. 2022) and the detections of so-called Doppler flips (change of sign in the non-Keplerian feature; e.g., Casassus & Pérez 2019; Norfolk et al. 2022), meridional flows within dust-depleted gaps (Teague et al. 2019a), and a velocity perturbation associated with a circumplanetary disk candidate (Bae et al. 2022b) enable us to take a closer look at planet–disk interaction processes. While most localized kinematical perturbations are analyzed empirically, statistical methods for quantifying their significance have been developed and led to the detection of localized signatures possibly associated with unseen planets (Izquierdo et al. 2021, 2022). Prime targets to search for protoplanets still embedded in their birth environment are the so-called transition disks. As in PDS70 (Keppler et al. 2019) or AB Aur (Tang et al. 2017), these disks host a dust-depleted cavity that has possibly been cleared by massive companions (Zhu et al. 2011).

In this Letter we focus on RXJ1604.3−2130 A (d = 144.6 pc, 1.24 M; Gaia Collaboration 2023; Manara et al. 2020, respectively), hereafter J1604, one of the brightest protoplanetary disks of the Upper Scorpius Association in the millimeter regime (Barenfeld et al. 2016), which exhibits a prominent cavity in the dust continuum and CO line emission (Zhang et al. 2014; Dong et al. 2017; van der Marel et al. 2021). J1604 has a stellar companion located at ∼2300 au, itself a binary with separation 13 au (Köhler et al. 2000). The outer disk of J1604 was resolved with the Atacama Large Millimeter/submillimeter Array (ALMA; Mayama et al. 2018) and the Spectro-Polarimetric High-contrast Exoplanet REsearch (SPHERE) instrument on the Very Large Telescope (VLT; Pinilla et al. 2015), indicating a nearly face-on geometry. Complementary observations are indicative of a misaligned inner disk with respect to the outer disk. Its variable light curve is that of an irregular dipper (Ansdell et al. 2020), infrared scattered-light observations show the presence of two shadows with variable morphology on timescales possibly shorter than a day (Pinilla et al. 2018), and ALMA 12CO (J = 3–2) line observations show deviations from Keplerian rotation in the cavity (Mayama et al. 2018). The position of the scattered-light shadows are suggestive of a large misalignment (∼70–90°). Measurements of the projected rotational velocity (v sin i) indicate that the star is aligned with the inner disk, and thus misaligned with the outer disk (Sicilia-Aguilar et al. 2020).

In this work, we present new ALMA observations of J1604 and focus on the kinematics of the 12CO (J = 2–1) line. Section 2 presents the observations and Sects. 3 and 4 our methodology and results, respectively. Section 5 provides a discussion of the results and Sect. 6, our conclusions.

2. Observations, calibration, and imaging

We present new ALMA Band 6 observations (2018.1.01255.S; PI: Benisty) with five executions spread over two years, obtained on 2019 April 4, 2019 July 30 and 31, 2021 April 29, and 2021 September 27. The spectral setup was designed for continuum detection but includes the 12CO J = 2–1 line. The data were combined with archival data from program 2015.1.00964.S (PI Oberg; see Table A.1). The data calibration and imaging were performed following the procedure of Andrews et al. (2018), with CASA v.5.6.1 (McMullin et al. 2007), and is detailed in Appendix A. The synthesized beam of the 12CO line and dust continuum images is 0.18″ × 0.15″ (102°) and 0.060″ × 0.039″ (−78°), respectively. The rms in a line-free channel was measured to be 1.1 mJy beam−1 (4.3 K) for CO and 10 mJy beam−1 for the dust continuum. Figure 1 shows the dust continuum map (left), which displays a cavity and a bright dust ring peaking at R∼ 0.56″ (∼81 au), and the 12CO peak brightness temperature map (center), which indicates a smaller cavity in gas, with a peak at R ∼ 0.39″ (∼56 au). A selection of channel maps can be found in Fig. A.1.

thumbnail Fig. 1.

ALMA observations of J1604. Panel a: 231 GHz dust continuum. The solid black contours are drawn at [5, 15, 25, 35, 45]σ, and the image is plotted with a power-law scaling of γ = 0.6. Panel b: 12CO peak brightness temperature map computed from I0 using the Planck law with black solid contours drawn at [5, 10, 20, 40, 60, 65, 70] σ. Pixels below 5σ are masked. Panel c: peak intensity residuals after subtracting an azimuthally averaged radial profile from the data. We have adjusted the color scale such that residuals smaller than 1σ are white. The beam sizes are shown in the lower-left corner, and the position of the star is marked by a green cross. In panels b and c, we have overlaid the continuum contours in white and black, respectively.

3. Methodology

3.1. Channel maps model

To model the disk line intensity and kinematics, we used the discminer package of Izquierdo et al. (2021). The code uses parametric prescriptions for the line peak intensity, line width, rotational velocity, and disk emission height to produce channel maps; it uses emcee (Foreman-Mackey et al. 2013) to maximize a χ2 log-likelihood function of the difference between the model and input intensity for each pixel in a channel map. To prescribe the model intensity, we used a generalized bell kernel function of the disk cylindrical coordinates (R, z):

(1)

where Ip is the peak intensity, Lw is half the line width at half power, hereafter “the line width,” and Ls is the line slope. The vch is the channel velocity at which the intensity is computed and v the observed Keplerian line-of-sight velocity. As the disk is nearly face-on, the code is unable to infer an emission height, and we therefore assumed a flat emission surface. We additionally fixed the inclination, i, of the disk to the one inferred from the dust continuum (i = 6.0°; Dong et al. 2017) to break the degeneracy of M ⋅ sin i. The fitting procedure and the Markov chain Monte Carlo (MCMC) search are explained in detail in Appendix B, where the functional form of each model parameter together with its best-fit parameters are summarized in Table B.1. We compare selected channel maps to the best-fit model using these parameters in Fig. A.1.

3.2. Moment maps

The moment maps were computed with bettermoments (Teague & Foreman-Mackey 2018). Since the 12CO line emission is optically thick, we fitted the following line profile to both the data and model channel maps:

(2)

where I0 is the peak intensity of the line and the optical depth, τ (v), varies like a Gaussian; v0 is the line centroid, τ0 the peak optical depth, and ΔV the width of the line – where the full width at half maximum (FWHM) = – as used in Teague et al. (2022). In Fig. 1b we show I0 for 12CO in units of brightness temperature. The corresponding v0 maps for the data and model are displayed in panels a and b of Fig. 2, respectively. The moment maps for ΔV and τ0, as well the error of the line centroid fitting, δv0, can be found in Fig. A.4.

thumbnail Fig. 2.

Line-of-sight velocity maps for v0 data (a) and the discminer model, vmod (b). Panel c: velocity residual map after subtracting vmod from v0. The dust continuum is overlaid in solid contours with equal levels, as in Fig. 1. The innermost region was masked during the fit by one beam size in radius, shown as the gray shaded ellipse. The insets in panels a and c zoom into the innermost region of the disk to highlight the non-Keplerian velocities. Contours are drawn at vsys = (4.62 ± 0.60) km s−1 in steps of 0.1 km s−1 and from −60 to 60 m s−1 in steps of 10 m s−1, respectively. All maps show the synthesized beam for CO (black) and the continuum (white) in the lower-left corner. The beams are masked where the CO peak intensity falls below the 5σ level for panel a and where R > Rout for the rest.

4. Results

4.1. Dust and gas radial and azimuthal brightness profiles

Figure 1 shows the 1.3 mm dust continuum together with the peak brightness temperature map, I0, of the 12CO (J = 2 − 1) line emission. Both dust and gas tracers show a cavity, and the 12CO (J = 2 − 1) line emission extends inward of the dust continuum, as expected if the continuum ring results from dust trapping (e.g., Facchini et al. 2018b, see our Fig. A.2). We note that the 12CO cavity appears asymmetric with respect to the position of the star and that we observe a gap-like feature in the 12CO peak intensity map at R ∼ 1.2″, apparent as a plateau of I0 ≈ 31 mJy beam−1 stretching over ΔR ≈ 0.1″ (Fig. A.2). Interestingly, the disk shows significant azimuthal intensity variations (34% at R = 0.56″; 19% at 0.39″ for continuum and gas, respectively; see also Fig. A.3). Figure 1c shows the residuals obtained after subtracting an azimuthally averaged radial profile from the 12CO peak brightness temperature map. Azimuthal variations are clearly apparent within the dust cavity, with residual values of > 10σ. The fainter regions, distributed along the east-west direction, are broadly aligned with fainter regions seen in the continuum (see the contours of Figs. 1b and A.3) and with the shadows reported in scattered light (Pinilla et al. 2018; Kurtovic et al., in prep.).

4.2. Kinematical features

4.2.1. Localized velocity residuals

The centroid residual map in Fig. 2c shows a prominent, localized non-Keplerian velocity feature of δv ≈ 60 m s−1, between ∼0.35″ and 0.55″ (i.e., 50–80 au), which is at the edge of the dust continuum ring and oriented at PA ≈ (270 ± 15)°. To assess its significance, we followed the variance peak method from Izquierdo et al. (2021). First, the centroid velocity residuals were folded and subtracted along the disk minor axis to remove axisymmetric features. Second, a 2D scan was performed to search for peak velocity residuals and obtain their locations in the folded map. Using these detected points, a K-means clustering algorithm then searched for coherent velocity perturbations within predefined radial and azimuthal bins (MacQueen 1967; Pedregosa et al. 2011). We considered seven radial and ten azimuthal bins, which corresponds to a width of roughly one beam size, to identify clusters. The algorithm then subdivided the input residual points such that the center of each cluster is the closest to all points in the cluster, by iteratively minimizing the sum of squared distances from the data points to the center of the cluster. This led to irregularly spaced bin boundaries since the cluster centers are near the densest accumulations of points.

In Fig. 3 we show the folded velocity residual map together with the detected peak velocity residuals (gray points). The location of the detected peak velocity residuals in azimuth and radius, within identified clusters, can be found in Fig. B.1. Clusters with high significance (those with peak velocity residuals larger than three times the variance in other clusters) are located within one radial and azimuthal bin shown in Fig. 3. Taking the centers of the selected clusters allowed us to identify a localized perturbation at 0.28″ ± 0.07″ (R = 41 ± 10 au) and PA = 280° ±2°. The reported errors are the standard deviation of the peak residual point (R, ϕ) locations within the selected clusters. The detection yields a cluster significance of 5.4σϕ in azimuth and 5.3σR in radius, where σ represents the standard deviation of background cluster variances with a mean of σϕ = 0.034 km2 s−2 and σR = 0.018 km2 s−2 (see black crosses in Fig. B.1). We note that a localized signature is robustly detected regardless of the amount of clusters defined, which we tested using 7–12 azimuthal or 5–9 radial clusters. We reported the clusters associated with the highest significance. Additionally, we note that there are other detections with lower significance at 0.65″ (94 au). This means that the radial extent of the prominent perturbation is roughly 0.40″ (58 au) and the global peak of the folded velocity residuals is at 0.39″ (56 au; as can be seen in the middle panel of Fig. B.1). This analysis confirms the presence of a significant localized non-Keplerian feature, as identified visually in Fig. 2c, within the continuum ring.

thumbnail Fig. 3.

Folded velocity residuals (left) and detected clusters of peak velocities (right) in the disk reference frame. The green shaded annulus and wedge in the right plot mark the significant clusters in radius and azimuth, respectively. The position of the localized velocity perturbation inferred from these clusters is marked with a magenta point with error bars. The gray region (one beam size in radius) indicates the masked area, and the black dashed lines the FWHM of the dust ring.

4.2.2. Spiral feature

Figure 2c also shows an extended arc-like positive residual feature, beyond the dust continuum emission and covering nearly 300° in azimuth, that is more evident in the polar deprojection of the velocity residual map (Fig. 4). To assess if this feature is a coherent structure, we used the FilFinder package (Koch & Rosolowsky 2015) implemented in discminer between 0.30″ and 1.25″ (43–180 au; see Fig. A.5). As indicated by the coherent filaments, the strong localized positive velocity residual discussed in Sect. 4.2.1 seems to be the starting point of a spiral tracing outward up to roughly 1.1″ (159 au). In Fig. 4 we overlay an Archimedean (linear) spiral, prescribed by rspiral = a + bϕspiral, using {a, b}={0.48, 0.12}. Computing the pitch angles tan(β) = − (dr/dϕ)/r, we obtain values ranging from 14° to 6° over the spiral extent.

thumbnail Fig. 4.

Polar projection of the velocity residual map. The solid black line shows a linear spiral trace. The gray region indicates the masked area and the black dashed lines, the FWHM of the dust ring. The y-axis extends farther than 360° to enhance the visibility of the spiral.

4.2.3. A possible warp in the 12CO cavity

An additional feature clearly evident from the velocity maps is the highly perturbed inner disk regions. As seen in the inset of the v0 line centroid map in Fig. 2a, the iso-velocity lines show strong bending in the inner region (∼3 beam sizes in diameter), indicative of non-Keplerian velocities. It likely traces a warp or a misaligned inner disk, as reported in Mayama et al. (2018), Pinilla et al. (2018), Sicilia-Aguilar et al. (2020), and Ansdell et al. (2020) to explain the scattered-light shadows and variable photometry of J1604. Higher angular resolution deep gas observations are, however, needed to assess its morphology and kinematics.

4.3. Deprojected velocity components

To understand the contributions from vϕ, vr, vz, we produced three centroid residual maps for each velocity component, after deprojection (Eq. (C.1)), assuming that “all” the velocities are azimuthal, radial, or vertical (see Fig. C.2; Teague et al. 2022). The localized residual feature at the edge of the dust ring appears to trace variations in the vertical motions, vz, or in the rotational motions, vϕ, or a combination of both. Radial perturbations can be ruled out since the feature is located close to the disk redshifted major axis (PAdisk = 258°), where vr, proj ≈ 0. Assuming purely rotational velocities, this corresponds to perturbations as high as δvϕ ≈ 600 m s−1 (∼0.4 × vkepler) due to the low disk inclination.

As seen in Fig. 2c, the spiral-like velocity residual feature does not change sign around the disk major or minor axes, which would occur for rotational, vϕ, or radial, vr, velocity perturbations, respectively (see Eq. (C.1)). We are likely seeing vertical perturbations, which we are most sensitive to in a nearly face-on disk.

Figure C.1 shows the deprojected and azimuthally averaged radial profiles of each velocity component determined with eddy. For vϕ, we observe super-Keplerian rotation from R ∼ 0.35″–0.70″ (51–101 au), peaking at 0.45″ (65 au), right beyond the dust continuum. The rotational velocities then sharply drop to sub-Keplerian in the inner disk regions. However, we stress that the azimuthally averaged velocities at the radial location of the strong localized perturbation (R ∼ 0.3 − 0.6″,43–87 au) are likely affected by the feature. We tentatively observe radial inflow inward of the CO intensity peak but with very large uncertainties on vr. Finally, we mostly detect downward vertical motion of the disk within R ∼ 1.25″ (181 au).

5. Discussion

5.1. Origin of v0 residuals

In this Letter we report the detection of two main non-Keplerian features, in addition to highly perturbed gas velocities in the gas cavity: (1) a localized positive residual near the edge of the dust ring and (2) an extended spiral-like feature, possibly starting from (1). A variety of velocity residual features were detected in other systems, with a diverse range of inclinations (e.g.,Wölfer et al. 2023). In the case of TW Hya, a similarly face-on disk, the detected perturbations are ∼40 m s−1 (Teague et al. 2022), which is lower than what is derived for (1), and in our case can account for 40% of the local Keplerian velocity assuming that the perturbation is purely due to rotational velocities. This is also larger in the magnitude of deviation than the Doppler flip reported in the HD 100546 transition disk (Casassus et al. 2022). These velocity residuals are often interpreted as tracing planet–disk interactions from massive companions (Pinte et al. 2022) that potentially carve out gaps. It is thus worth noting that our inferred planet location (R = 41 ± 10 au) is close to the gap location in 13CO (J = 2 − 1) at 37 au reported in van der Marel et al. (2021). Comparison with simulations (Rabago & Zhu 2021; Izquierdo et al. 2022) or semi-analytical prescriptions (Bollati et al. 2021) allowed us to estimate a possible planet mass from the velocity deviations. To this end, we considered Eq. (14) from Yun et al. (2019), which relates the change in rotational velocity, δvϕ, to the planet mass, Mp, through 2D hydrodynamic simulations.

Since δV is the difference between the super- and sub-Keplerian peak, we considered the peak of the super-Keplerian motion, (vϕ/vkep)/vkep (see Fig. C.1), as a lower limit for the dimensionless amplitude of the perturbed rotational velocity ( = 0.06), and its double () as an upper limit. Assuming (H/R)p = 0.1 at the planet location (R = 41 au), we estimate the planet mass to roughly be in the range Mp ≈ (1.6 − 2.9)Mjup (α/10−3)0.5.

The extended spiral-like feature appears to be related to the significant localized velocity residual. Due to its low pitch angles and the consistent positive velocity residuals, we speculate that the spiral is caused by buoyancy resonances excited by a massive planet located within the dust ring. Indeed, in contrast with Lindblad spirals, buoyancy spirals are shown to exhibit a tightly wound morphology with predominantly vertical motions (Bae et al. 2021). Such a spiral has also been suggested in TW Hya (Teague et al. 2019b), whose radial extent is similar to the one reported here. Interestingly, Wölfer et al. (2023) report a tentative arc feature in J1604, at R ∼ 1.0″ and ranging from PA ≈ 160 − 200°, probed by the kinematics of the 12CO (J = 3 − 2) line emission, that partly coincides with the spiral-like feature that we detect. Additional observations in optically thin tracers would be very useful to assess its nature.

5.2. Kinematic perturbations due to shadows

The localized residual feature seems to roughly align with the shadows detected in scattered light. Comparing Figs. 1c and 2c, positive and negative v0 residuals broadly align with cold and hot regions in the brightness temperature of 12CO, respectively (see also Fig. A.6). In particular, the orientation of the significant localized velocity perturbation coincides with the western shadow. Such a shadow can cool down the disk material and possibly induce a local drop in pressure support and therefore impact the gas velocity. In this section we estimate whether the detected velocity perturbations could be caused by azimuthal variations in temperature. We related the azimuthal change in temperature, ΔTϕ, to variations in rotational velocity, Δvϕ, by solving for the Navier–Stokes equation in cylindrical coordinates. Following the derivation in Appendix D, we obtain

(3)

with H/R the disk aspect ratio, T the 12CO brightness temperature tracing the gas temperature (as 12CO is optically thick), and vkep the Keplerian velocity. Evaluating this equation for a large H/R = 0.2 along a radially averaged annulus centered at R = (0.39 ± 0.07)″, where we experience the strongest azimuthal changes in the 12CO brightness temperature of up to ΔT ≈ 15 K over T = 60 K, we estimate the change in azimuthal velocity to be a mere δvϕ/vkep ≈ 1%. Hence, the shadows cannot be solely responsible for the localized velocity residual feature. In addition, as the shadows are nearly symmetric, we would expect a similar feature opposite to its east.

Montesinos et al. (2016) investigated the development of spirals due to pressure gradients caused by temperature differences between obscured and illuminated regions. In their simulations, symmetric shadows always form two-armed spirals. However, they only develop for massive (0.25 M) and/or strongly illuminated disks (100 L), which does not seem to be the case for J1604 (∼0.02 M and 0.7 L, Manara et al. 2020), where we also only observe one spiral. Additionally, we would also expect such spiral features to appear in the brightness temperature residuals (Fig. 1c). It is therefore unlikely that shadows are responsible for the extended spiral-like velocity residual feature.

5.3. A warped or misaligned inner disk?

The bending of the iso-velocity curves that we observe in the inset of Fig. 2a is reminiscent of a warped or (broken) misaligned inner disk (Juhász & Facchini 2017; Facchini et al. 2018a). However, as the inner disk is unresolved in our observations, the warp morphology cannot be derived. We note that radial inflows are also degenerate in appearance with warps, as shown by Rosenfeld et al. (2014), and that our observations do not allow us to be conclusive regarding the origin of the disturbed kinematics in the innermost disk. We attempted to infer if the position angle or inclination vary with radius by fitting the innermost disk only (R ≤ 0.5″) with eddy and considering a fixed stellar mass, but we did not find any significant variations relative to our best-fit values. We therefore constrain the warp to be confined within one beam size (∼0.18″, i.e., 26 au) from the center. We note that we obtain a 5% higher dynamical mass of the system when fitting for M while masking the innermost beam size in radius, an effect predicted by hydrodynamical simulations of warps (Young et al. 2022).

While our observations cannot provide a full picture of the system due to a limited angular resolution, the very low mass accretion rate and near-infrared excess (Sicilia-Aguilar et al. 2020), as well as the gas cavity in 12CO with non-Keplerian velocities, suggest the presence of an additional, very massive (possibly stellar) companion within the inner ∼0.25″ (∼35 au). Such a companion would need to be on an inclined (nearly polar) orbit to misalign the inner disk (Zhu 2019), which would then lead to the shadows (Nealon et al. 2019) and variable extinction events seen in the light curves (Ansdell et al. 2020; Sicilia-Aguilar et al. 2020). It would not, however, explain the strong localized velocity residual feature that we report, which we speculate traces a planetary-mass object located at the edge of the dust continuum. Detailed modeling of the system is thus needed to assess the need for an additional companion. An interesting comparison is the CS Cha spectro-binary system (separation of ∼7 au), which shows a similar dust continuum and gas emission at a similarly low inclination but no departure from Keplerian rotation in the 12CO kinematics (Kurtovic et al. 2022).

6. Conclusions

In this Letter we present new ALMA observations of the 1.3 mm dust continuum and the 12CO (J = 2 − 1) line emission from the transition disk around RXJ1604.3–2130 A. The dust continuum shows a large cavity enclosing a smaller 12CO cavity. Azimuthal brightness variations in the 12CO line and dust continuum are broadly aligned with shadows detected in scattered light (Pinilla et al. 2018). Using the discminer package (Izquierdo et al. 2021), we modeled the channel-by-channel line emission and calculated the line-of-sight velocity maps. We report the detection of a coherent, localized non-Keplerian feature at R = (41 ± 10) au and PA = 280° ±2°, that is, within the continuum ring. While broadly aligned with the scattered-light shadows, the localized non-Keplerian feature cannot be due to changes in temperature. Instead, we interpret the kinematical perturbation as tracing the presence of a massive companion of Mp ≈ (1.6 − 2.9) Mjup. We also detect a tightly wound spiral that extends over 300° in azimuth, possibly connected to the localized feature and caused by buoyancy resonances driven by planet–disk interactions. Bending of the iso-velocity contours within the gas cavity indicates a highly perturbed inner region, possibly related to the presence of a misaligned inner disk. However, as the putative planet at ∼41 au cannot explain the gas cavity, the low accretion rate, or the misaligned inner disk, we speculate that another massive companion, likely on an inclined orbit, shapes the inner ∼0.25″(∼35 au).

Acknowledgments

We would like to thank the anonymous referee for the constructive feedback, as well as Clement Baruteau, Kees Dullemond, Guillaume Laibe and Andrew Winter for helpful discussions. This Letter makes use of the following ALMA data: ADS/JAO.ALMA#2017.A.01255.S and ADS/JAO.ALMA#2015.1.00964.S. ALMA is a partnership of ESO (representing its member states), NSF (USA), and NINS (Japan), together with NRC (Canada), NSC 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. This project has received funding from the European Research Council (ERC) under the European Union’s Horizon 2020 research and innovation programme (PROTOPLANETS, grant agreement No. 101002188). L.P. gratefully acknowledges support by the ANID BASAL projects ACE210002 and FB210003, and by ANID, – Millennium Science Initiative Program – NCN19_171. Software: CARTA (Comrie et al. 2021), CASA (McMullin et al. 2007), Discminer (Izquierdo et al. 2021), Eddy (Teague 2019a), FilFinder (Koch & Rosolowsky 2015), GoFish (Teague 2019b), Matplotlib (Hunter 2007), Numpy (van der Walt et al. 2011), Scipy (Virtanen et al. 2020).

References

  1. Andrews, S. M., Huang, J., Pérez, L. M., et al. 2018, ApJ, 869, L41 [NASA ADS] [CrossRef] [Google Scholar]
  2. Ansdell, M., Gaidos, E., Hedges, C., et al. 2020, MNRAS, 492, 572 [Google Scholar]
  3. Bae, J., Teague, R., & Zhu, Z. 2021, ApJ, 912, 56 [NASA ADS] [CrossRef] [Google Scholar]
  4. Bae, J., Isella, A., Zhu, Z., et al. 2022a, ArXiv e-prints [arXiv:2210.13314] [Google Scholar]
  5. Bae, J., Teague, R., Andrews, S. M., et al. 2022b, ApJ, 934, L20 [NASA ADS] [CrossRef] [Google Scholar]
  6. Barenfeld, S. A., Carpenter, J. M., Ricci, L., & Isella, A. 2016, ApJ, 827, 142 [Google Scholar]
  7. Benisty, M., Dominik, C., Follette, K., et al. 2022, ArXiv e-prints [arXiv:2203.09991] [Google Scholar]
  8. Bollati, F., Lodato, G., Price, D. J., & Pinte, C. 2021, MNRAS, 504, 5444 [Google Scholar]
  9. Calcino, J., Hilder, T., Price, D. J., et al. 2022, ApJ, 929, L25 [NASA ADS] [CrossRef] [Google Scholar]
  10. Casassus, S., & Pérez, S. 2019, ApJ, 883, L41 [Google Scholar]
  11. Casassus, S., Cárcamo, M., Hales, A., Weber, P., & Dent, B. 2022, ApJ, 933, L4 [NASA ADS] [CrossRef] [Google Scholar]
  12. Comrie, A., Wang, K. S., Hsu, S. C., et al. 2021, https://doi.org/10.5281/zenodo.3746095 [Google Scholar]
  13. Czekala, I., Loomis, R. A., Teague, R., et al. 2021, ApJS, 257, 2 [NASA ADS] [CrossRef] [Google Scholar]
  14. Dong, R., van der Marel, N., Hashimoto, J., et al. 2017, ApJ, 836, 201 [NASA ADS] [CrossRef] [Google Scholar]
  15. Facchini, S., Juhász, A., & Lodato, G. 2018a, MNRAS, 473, 4459 [Google Scholar]
  16. Facchini, S., Pinilla, P., van Dishoeck, E. F., & de Juan Ovelar, M. 2018b, A&A, 612, A104 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  17. Foreman-Mackey, D., Hogg, D. W., Lang, D., & Goodman, J. 2013, PASP, 125, 306 [Google Scholar]
  18. Gaia Collaboration (Vallenari, A., et al.) 2023, A&A, in press, https://doi.org/10.1051/0004-6361/202243940 [Google Scholar]
  19. Hunter, J. D. 2007, Comput. Sci. Eng., 9, 90 [Google Scholar]
  20. 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]
  21. Izquierdo, A. F., Facchini, S., Rosotti, G. P., van Dishoeck, E. F., & Testi, L. 2022, ApJ, 928, 2 [NASA ADS] [CrossRef] [Google Scholar]
  22. Jorsater, S., & van Moorsel, G. A. 1995, AJ, 110, 2037 [Google Scholar]
  23. Juhász, A., & Facchini, S. 2017, MNRAS, 466, 4053 [NASA ADS] [Google Scholar]
  24. Keppler, M., Teague, R., Bae, J., et al. 2019, A&A, 625, A118 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  25. Koch, E. W., & Rosolowsky, E. W. 2015, MNRAS, 452, 3435 [Google Scholar]
  26. Köhler, R., Kunkel, M., Leinert, C., & Zinnecker, H. 2000, A&A, 356, 541 [NASA ADS] [Google Scholar]
  27. Kostov, V. B., Orosz, J. A., Welsh, W. F., et al. 2016, ApJ, 827, 86 [NASA ADS] [CrossRef] [Google Scholar]
  28. Kurtovic, N. T., Pinilla, P., Penzlin, A. B. T., et al. 2022, A&A, 664, A151 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  29. Law, C. J., Teague, R., Loomis, R. A., et al. 2021, ApJS, 257, 4 [NASA ADS] [CrossRef] [Google Scholar]
  30. Long, F., Pinilla, P., Herczeg, G. J., et al. 2018, ApJ, 869, 17 [Google Scholar]
  31. MacQueen, J. 1967, in 5th Berkeley Symp. Math. Statist. Probability, 281 [Google Scholar]
  32. Manara, C. F., Natta, A., Rosotti, G. P., et al. 2020, A&A, 639, A58 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  33. Mayama, S., Akiyama, E., Panić, O., et al. 2018, ApJ, 868, L3 [NASA ADS] [CrossRef] [Google Scholar]
  34. McMullin, J. P., Waters, B., Schiebel, D., Young, W., & Golap, K. 2007, in Astronomical Data Analysis Software and Systems XVI, eds. R. A. Shaw, F. Hill, & D. J. Bell, ASP Conf. Ser., 376, 127 [Google Scholar]
  35. Montesinos, M., Perez, S., Casassus, S., et al. 2016, ApJ, 823, L8 [Google Scholar]
  36. Nealon, R., Pinte, C., Alexander, R., Mentiplay, D., & Dipierro, G. 2019, MNRAS, 484, 4951 [CrossRef] [Google Scholar]
  37. Norfolk, B. J., Pinte, C., Calcino, J., et al. 2022, ApJ, 936, L4 [NASA ADS] [CrossRef] [Google Scholar]
  38. Pedregosa, F., Varoquaux, G., Gramfort, A., et al. 2011, J. Mach. Learn. Res., 12, 2825 [Google Scholar]
  39. Perez, S., Dunhill, A., Casassus, S., et al. 2015, ApJ, 811, L5 [NASA ADS] [CrossRef] [Google Scholar]
  40. Pinilla, P., de Boer, J., Benisty, M., et al. 2015, A&A, 584, L4 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  41. Pinilla, P., Benisty, M., de Boer, J., et al. 2018, ApJ, 868, 85 [NASA ADS] [CrossRef] [Google Scholar]
  42. Pinte, C., Teague, R., Flaherty, K., et al. 2022, ArXiv e-prints [arXiv:2203.09528] [Google Scholar]
  43. Rabago, I., & Zhu, Z. 2021, MNRAS, 502, 5325 [NASA ADS] [CrossRef] [Google Scholar]
  44. Rich, E. A., Teague, R., Monnier, J. D., et al. 2021, ApJ, 913, 138 [NASA ADS] [CrossRef] [Google Scholar]
  45. Rich, E. A., Monnier, J. D., Aarnio, A., et al. 2022, AJ, 164, 109 [NASA ADS] [CrossRef] [Google Scholar]
  46. Rosenfeld, K. A., Chiang, E., & Andrews, S. M. 2014, ApJ, 782, 62 [Google Scholar]
  47. Rosotti, G. P., Benisty, M., Juhász, A., et al. 2020, MNRAS, 491, 1335 [Google Scholar]
  48. Sicilia-Aguilar, A., Manara, C. F., de Boer, J., et al. 2020, A&A, 633, A37 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  49. Tang, Y.-W., Guilloteau, S., Dutrey, A., et al. 2017, ApJ, 840, 32 [NASA ADS] [CrossRef] [Google Scholar]
  50. Teague, R. 2019a, The Journal of Open Source Software, 4, 1220 [NASA ADS] [CrossRef] [Google Scholar]
  51. Teague, R. 2019b, J. Open Source Softw., 4, 1632 [NASA ADS] [CrossRef] [Google Scholar]
  52. Teague, R., & Foreman-Mackey, D. 2018, https://doi.org/10.5281/zenodo.1419754 [Google Scholar]
  53. Teague, R., Bae, J., Bergin, E. A., Birnstiel, T., & Foreman-Mackey, D. 2018a, ApJ, 860, L12 [NASA ADS] [CrossRef] [Google Scholar]
  54. Teague, R., Bae, J., Birnstiel, T., & Bergin, E. A. 2018b, ApJ, 868, 113 [Google Scholar]
  55. Teague, R., Bae, J., & Bergin, E. A. 2019a, Nature, 574, 378 [NASA ADS] [CrossRef] [Google Scholar]
  56. Teague, R., Bae, J., Huang, J., & Bergin, E. A. 2019b, ApJ, 884, L56 [Google Scholar]
  57. Teague, R., Bae, J., Andrews, S. M., et al. 2022, ApJ, 936, 163 [NASA ADS] [CrossRef] [Google Scholar]
  58. van der Marel, N., Birnstiel, T., Garufi, A., et al. 2021, AJ, 161, 33 [Google Scholar]
  59. van der Walt, S., Colbert, S. C., & Varoquaux, G. 2011, Comput. Sci. Eng., 13, 22 [Google Scholar]
  60. Virtanen, P., Gommers, R., Burovski, E., et al. 2020, scipy/scipy: SciPy 1.5.3 [Google Scholar]
  61. Wölfer, L., Facchini, S., van der Marel, N., et al. 2023, A&A, in press, https://doi.org/10.1051/0004-6361/202243601 [Google Scholar]
  62. Young, A. K., Alexander, R., Rosotti, G., & Pinte, C. 2022, MNRAS, 513, 487 [NASA ADS] [CrossRef] [Google Scholar]
  63. Yun, H.-G., Kim, W.-T., Bae, J., & Han, C. 2019, ApJ, 884, 142 [NASA ADS] [CrossRef] [Google Scholar]
  64. Zhang, K., Isella, A., Carpenter, J. M., & Blake, G. A. 2014, ApJ, 791, 42 [NASA ADS] [CrossRef] [Google Scholar]
  65. Zhu, Z. 2019, MNRAS, 483, 4221 [Google Scholar]
  66. Zhu, Z., Nelson, R. P., Hartmann, L., Espaillat, C., & Calvet, N. 2011, ApJ, 729, 47 [NASA ADS] [CrossRef] [Google Scholar]

Appendix A: Observations, calibration, and imaging

To self-calibrate our observations, we proceeded as follows. We first flagged the channels containing the line to produce a continuum data set. We centered the individual execution blocks (EBs) by fitting the continuum visibilities with a ring model, allowing for a different center and amplitude. This enabled us to recover the phase shift and amplitude re-scaling and apply them to the EBs before combining them. To determine a good initial model for the self-calibration, we used multi-scale cleaning with the tclean task using a threshold of two times the rms noise level of the image. Using the tasks gaincal and applycal, we corrected for phase offsets between spectral windows, and between polarizations, considering a solution interval of the scan length (solint=inf). Executions obtained in 2019 were concatenated and self-calibrated together, and a similar procedure was applied to those obtained in 2021. In addition to the first round of self-calibration, two additional iterations of phase self-calibration were done, with solution intervals of 300s and 180s for the 2019 data; only one round was done for the 2021 data, with a solution interval of 360s. For both data sets, a round of amplitude self-calibration was applied with solint=inf. The solutions were then applied to the gas data. While these two epochs will be analyzed separately for the continuum in a forthcoming paper (Kurtovic et al.), to analyze the gas data, we concatenated them after checking that the data do not show significant variations between the two epochs. We imaged the resulting visibilities with the tclean task using the multi-scale CLEAN algorithm with scales of 0, 1, 3, and 6 times the beam FWHM and an elliptic CLEAN mask encompassing the disk emission. The 12CO (2-1) molecular line observations are imaged with a robust value of 1.0 and a channel width of 0.1 km s−1 and masked by the 4.0 σ threshold. The data were tapered to 0.″1, and we used the “JvM correction” (Jorsater & van Moorsel 1995; Czekala et al. 2021).

In Fig. A.7 we show a comparison of velocity residual maps for additional cubes, imaged using different imaging parameters, to assess the robustness of our detections. We compare residual maps for the same cube as used in the main text, but without JvM correction, and for a cube imaged with a different tapering (0.15″ instead of 0.10″). Best-fit Keplerian models were subtracted from each of the cubes. As is evident from the comparison of the residual maps, the detection of the reported non-Keplerian features are robust irrespective of the imaging parameters. However, the detailed morphology of the velocity residual peaks changes with imaging parameters and, as a consequence, the value of the inferred planet location from the discminer analysis. We find that the best-fit discminer models to the non-JvM-corrected cubes are similar within 3% with respect to the best-fit parameters listed in Appendix B, with the exception of the line slope and some of the peak intensity parameters, which vary by up to 15%. While estimating the systematics due to imaging parameters is beyond the scope of this Letter, Fig. A.7 provides evidence that the detection of non-Keplerian features is robust. We measure the rms in a line-free channel to be 2.6 mJy beam−1 and 2.9 mJy beam−1 for the non-JvM-corrected cubes with a taper of 0.10 and 0.15, respectively.

thumbnail Fig. A.1.

Gallery of selected channel maps. Panels show the 12CO data channel maps (top row) and best-fit model channel maps (middle row), together with intensity residuals in kelvins for each channel (bottom row). In the bottom row the colorbar has been adjusted such that residuals smaller than 1σ are white. The beam size is depicted in the lower-left corner of each channel. For reference, the best-fit systemic velocity was found to be vsys = 4.62 km s−1 and the channel spacing is 100 m s−1.

thumbnail Fig. A.2.

Radial profile of the surface brightness for different tracers. Profiles are normalized to the peak of the emission for the 231 GHz continuum and the CO peak flux (both for the data and discminer model) and for the SPHERE scattered-light observation. Shaded regions show the standard deviation of each azimuthal average. The lines in the lower-right corner show the major beam size (resolution) for each profile in the corresponding color.

thumbnail Fig. A.3.

Azimuthal profiles of the surface brightness, normalized to the peak of the emission. Profiles were extracted at an annulus with a width of approximately one corresponding beam size centered at 0.56″ and 0.39″ for the 231 GHz continuum and the CO peak flux, which both show significant azimuthal intensity variations, of 34% and 19%, respectively. Shaded regions show the standard deviation of each radial average.

thumbnail Fig. A.4.

Additional moment maps of the centroid fitting. Panels show the line width, ΔV (left), the peak optical depth, τ0 (center), and the error of the centroid fitting, δv0 (right). Note that for τ0 < 1 one can assume the line profile to be well presented by a Gaussian, while for τ0 > 5 the line profile has a saturated core, i.e., a very flat top (see Eq. 2). The beam size is depicted in the lower-left corner, and only regions where I0 > 5σ with σ = 1.1 mJy beam−1 are shown.

thumbnail Fig. A.5.

Polar map of the velocity residuals. Same as Fig. 4, but now overlaid with filamentary structures found by FilFinder. The overplotted red and blue lines are the medial axes of the filamentary structures found by the algorithm. To trace the apparent spiral in the residuals, we restricted the algorithm to search for filaments in the radial locations r = [0.3, 1.25]″. For the filamentary detection, we assumed a smoothing size of one synthesized beam size and a minimum size of 500 pixels for a filament to be considered.

thumbnail Fig. A.6.

Polar contour map of the centroid residuals. The upper panel shows residuals inside the cavity, the middle panel residuals between the cavity and the outer edge of dust ring, and the lower panel the residuals in the outer disk. The radial spacing between each contour is ∼1.8 au, and the opacity of the lines increases with radius. We would like to emphasize that the bump between PA≈ − (110 − 70)° in the middle panel includes most of the residual points we detect with the peak variance method of discminer, which can be readily seen in the left panel of Fig. B.1.

thumbnail Fig. A.7.

Comparison of the velocity residual maps for different imaging parameters. The left panel corresponds to the cube used in the main text, while the middle panel corresponds to the same cube without JvM correction. The right panel shows the residual maps for a non-JvM-corrected cube with a different taper (0.15″). In all panels the dust continuum is overlaid in solid contours with equal levels, as in Fig. 1. Best-fit Keplerian models were subtracted from each of the cubes. The detection of the non-Keplerian features is quite robust, irrespective of the imaging procedure.

Table A.1.

Summary of the ALMA Band 6 observations of J1604 presented in this paper.

Appendix B: Model best-fit parameters

For the initial emcee run, we used literature values for the position angle and stellar mass (PA = 260°, M = 1.24 M; Dong et al. 2017; Manara et al. 2020, respectively). The initial values of the other parameters were found by comparing the overall morphology between the data and a prototype model. We performed the MCMC search with 150 walkers that evolved for 2000 steps for an initial burn-in stage. We proceeded in two steps. First, we masked the disk region inward of the dust continuum and only fitted the outer disk (R > 90 au) to get a robust estimate of the stellar mass and avoid the code being confusing by the strongly non-Keplerian velocity features in the inner regions. In this run, we interestingly find a strong offset from the disk center in the x direction of −8.0 au. In a second step, we fixed the stellar mass and now fitted for the whole disk, masking an inner region corresponding to one major beam size (26 au) in radius, where effects of beam smearing are strongest. We ran 150 walkers for 20000 steps until reaching convergence, which presented as a nearly normal distribution of the walkers. The variance and median of the parameter walkers remain effectively unchanged after ∼7000 steps. The best-fit parameters are the median of the posterior distributions, and given errors are the 16 and 84 percentiles in the last 5000 steps of the 20000-step run, summarized in Table B.1.

thumbnail Fig. B.1.

Location of the folded peak velocity residuals. The detected points are shown in azimuth (left) and radius (middle) and were obtained using the peak variance method. Colors indicate the seven different radial clusters specified, where blue peak residual points are within a detected significant radial cluster. The black crosses are the velocity variances of the clusters plotted at the (R, ϕ) location of each cluster center. The centers of the accepted clusters (those with peak velocity residuals larger than three times the variance in other clusters) in radius and azimuth are marked with vertical black lines in both panels. The right-hand plot shows the normal distribution of the peak residual points in a histogram. Note that outliers of the distribution are related to the localized perturbation. The maximum value of all peak folded centroid residuals is at 0.39″ (57 au), its mean value is 39 m/s, and 1σv = 20 m/s (not to be mistaken with the cluster variances).

Table B.1.

Attributes of the discminer model for the 12CO intensity channel maps of the disk around J1604. PA is the position angle of the semimajor axis of the disk on the redshifted side, R the cylindrical radius, and D0 = 100 au a normalization constant for the line properties. The (down-sampled) pixel size of the model is 8.8 au.

Appendix C: Decomposition and deprojection of velocity components

To determine the rotation curves of each velocity component, we used the code eddy (Teague 2019a). We followed the method presented in Teague et al. (2018b), which uses a Gaussian process to determine the azimuthal, vϕ, and radial, vr, velocity components along a given annulus. To this end, we divided the disk into concentric annuli with a radial width of one-fourth of the synthesized beam (∼0.05″), which ranges from 0.18″ to 1.85″, and extracted the velocities over 20 iterations to minimize their standard deviation. To obtain the vertical velocity component, vz, we used the measured azimuthally averaged profiles of vϕ and vr and extended them to produce 2D maps, considering the projection of the following components along the line of sight:

(C.1)

where ϕ is the polar angle in the disk frame (such that ϕ= 0 corresponds to the redshifted major axis) and i the inclination of the disk. In the case of J1604, the disk rotates clockwise, which corresponds to a negative inclination in the above definition. We then subtracted these maps together with the systemic velocity, vsys, from the line-of-sight velocity, v0, map (Fig. 2,a) to obtain a map of the vertical velocity component, vz,  proj:

(C.2)

The radial profile of vz was obtained by deprojecting and azimuthally averaging its 2D velocity map. The radial profiles of the deprojected velocity components can be found in Fig. C.1.

thumbnail Fig. C.1.

Azimuthally averaged and deprojected azimuthal, radial, and vertical velocity components. The radial width of each annulus is one-fourth of the synthesized beam size. The error bars are given by the standard deviation for each velocity component averaged over the 20 iterations used.

thumbnail Fig. C.2.

J1604 deprojected velocity components. It is assumed that all velocities are azimuthal (left column), radial (central column), or vertical (right column). For the azimuthal and radial components, wedges along the minor and major axis have been masked as the observations are insensitive to these components (see Eq. C.1). In each panel the synthesized beam is shown in the lower-left corner.

Appendix D: Derivation of Δvϕ as a function of azimuthal temperature variations, ΔT

To relate the change in the brightness temperature, ΔT, of 12CO to variations in the rotational velocity, vϕ, we solved the Navier-Stokes equation in cylindrical coordinates in the ϕ direction:

(D.1)

where R is the cylindrical radius, ρg the gas density, and μ the mean molecular weight. In a first step, we assumed that the radial variations within the chosen annulus are negligible and inserted the disk gas pressure in the vertically isothermal assumption: , where kB is the Boltzmann constant and mp the proton mass. In a second step, we further assumed the gas density to be constant along the annulus ρg = const. and re-arranged the equation:

(D.2)

In the last step, we inserted the gas midplane density, assuming a disk aspect ratio of H/R = 0.2. Assuming that vϕ ≈ vkep and Σg ≈ 1g/cm2 (see Fig. 3 of Dong et al. 2017) at the location of the annulus at R ∼ 0.4″ (58 au), we can assess the order of magnitude of the second term on the left-hand side of the equation, which is only on the order of 10−6. Therefore, we neglected the second-order derivative and identified the sound speed, cs:

(D.3)

Equation D.3 now connects the fractional azimuthal temperature variation, ΔTϕ/T, to the rotational velocity deviation relative to Keplerian, Δvϕ/vkep.

All Tables

Table A.1.

Summary of the ALMA Band 6 observations of J1604 presented in this paper.

Table B.1.

Attributes of the discminer model for the 12CO intensity channel maps of the disk around J1604. PA is the position angle of the semimajor axis of the disk on the redshifted side, R the cylindrical radius, and D0 = 100 au a normalization constant for the line properties. The (down-sampled) pixel size of the model is 8.8 au.

All Figures

thumbnail Fig. 1.

ALMA observations of J1604. Panel a: 231 GHz dust continuum. The solid black contours are drawn at [5, 15, 25, 35, 45]σ, and the image is plotted with a power-law scaling of γ = 0.6. Panel b: 12CO peak brightness temperature map computed from I0 using the Planck law with black solid contours drawn at [5, 10, 20, 40, 60, 65, 70] σ. Pixels below 5σ are masked. Panel c: peak intensity residuals after subtracting an azimuthally averaged radial profile from the data. We have adjusted the color scale such that residuals smaller than 1σ are white. The beam sizes are shown in the lower-left corner, and the position of the star is marked by a green cross. In panels b and c, we have overlaid the continuum contours in white and black, respectively.

In the text
thumbnail Fig. 2.

Line-of-sight velocity maps for v0 data (a) and the discminer model, vmod (b). Panel c: velocity residual map after subtracting vmod from v0. The dust continuum is overlaid in solid contours with equal levels, as in Fig. 1. The innermost region was masked during the fit by one beam size in radius, shown as the gray shaded ellipse. The insets in panels a and c zoom into the innermost region of the disk to highlight the non-Keplerian velocities. Contours are drawn at vsys = (4.62 ± 0.60) km s−1 in steps of 0.1 km s−1 and from −60 to 60 m s−1 in steps of 10 m s−1, respectively. All maps show the synthesized beam for CO (black) and the continuum (white) in the lower-left corner. The beams are masked where the CO peak intensity falls below the 5σ level for panel a and where R > Rout for the rest.

In the text
thumbnail Fig. 3.

Folded velocity residuals (left) and detected clusters of peak velocities (right) in the disk reference frame. The green shaded annulus and wedge in the right plot mark the significant clusters in radius and azimuth, respectively. The position of the localized velocity perturbation inferred from these clusters is marked with a magenta point with error bars. The gray region (one beam size in radius) indicates the masked area, and the black dashed lines the FWHM of the dust ring.

In the text
thumbnail Fig. 4.

Polar projection of the velocity residual map. The solid black line shows a linear spiral trace. The gray region indicates the masked area and the black dashed lines, the FWHM of the dust ring. The y-axis extends farther than 360° to enhance the visibility of the spiral.

In the text
thumbnail Fig. A.1.

Gallery of selected channel maps. Panels show the 12CO data channel maps (top row) and best-fit model channel maps (middle row), together with intensity residuals in kelvins for each channel (bottom row). In the bottom row the colorbar has been adjusted such that residuals smaller than 1σ are white. The beam size is depicted in the lower-left corner of each channel. For reference, the best-fit systemic velocity was found to be vsys = 4.62 km s−1 and the channel spacing is 100 m s−1.

In the text
thumbnail Fig. A.2.

Radial profile of the surface brightness for different tracers. Profiles are normalized to the peak of the emission for the 231 GHz continuum and the CO peak flux (both for the data and discminer model) and for the SPHERE scattered-light observation. Shaded regions show the standard deviation of each azimuthal average. The lines in the lower-right corner show the major beam size (resolution) for each profile in the corresponding color.

In the text
thumbnail Fig. A.3.

Azimuthal profiles of the surface brightness, normalized to the peak of the emission. Profiles were extracted at an annulus with a width of approximately one corresponding beam size centered at 0.56″ and 0.39″ for the 231 GHz continuum and the CO peak flux, which both show significant azimuthal intensity variations, of 34% and 19%, respectively. Shaded regions show the standard deviation of each radial average.

In the text
thumbnail Fig. A.4.

Additional moment maps of the centroid fitting. Panels show the line width, ΔV (left), the peak optical depth, τ0 (center), and the error of the centroid fitting, δv0 (right). Note that for τ0 < 1 one can assume the line profile to be well presented by a Gaussian, while for τ0 > 5 the line profile has a saturated core, i.e., a very flat top (see Eq. 2). The beam size is depicted in the lower-left corner, and only regions where I0 > 5σ with σ = 1.1 mJy beam−1 are shown.

In the text
thumbnail Fig. A.5.

Polar map of the velocity residuals. Same as Fig. 4, but now overlaid with filamentary structures found by FilFinder. The overplotted red and blue lines are the medial axes of the filamentary structures found by the algorithm. To trace the apparent spiral in the residuals, we restricted the algorithm to search for filaments in the radial locations r = [0.3, 1.25]″. For the filamentary detection, we assumed a smoothing size of one synthesized beam size and a minimum size of 500 pixels for a filament to be considered.

In the text
thumbnail Fig. A.6.

Polar contour map of the centroid residuals. The upper panel shows residuals inside the cavity, the middle panel residuals between the cavity and the outer edge of dust ring, and the lower panel the residuals in the outer disk. The radial spacing between each contour is ∼1.8 au, and the opacity of the lines increases with radius. We would like to emphasize that the bump between PA≈ − (110 − 70)° in the middle panel includes most of the residual points we detect with the peak variance method of discminer, which can be readily seen in the left panel of Fig. B.1.

In the text
thumbnail Fig. A.7.

Comparison of the velocity residual maps for different imaging parameters. The left panel corresponds to the cube used in the main text, while the middle panel corresponds to the same cube without JvM correction. The right panel shows the residual maps for a non-JvM-corrected cube with a different taper (0.15″). In all panels the dust continuum is overlaid in solid contours with equal levels, as in Fig. 1. Best-fit Keplerian models were subtracted from each of the cubes. The detection of the non-Keplerian features is quite robust, irrespective of the imaging procedure.

In the text
thumbnail Fig. B.1.

Location of the folded peak velocity residuals. The detected points are shown in azimuth (left) and radius (middle) and were obtained using the peak variance method. Colors indicate the seven different radial clusters specified, where blue peak residual points are within a detected significant radial cluster. The black crosses are the velocity variances of the clusters plotted at the (R, ϕ) location of each cluster center. The centers of the accepted clusters (those with peak velocity residuals larger than three times the variance in other clusters) in radius and azimuth are marked with vertical black lines in both panels. The right-hand plot shows the normal distribution of the peak residual points in a histogram. Note that outliers of the distribution are related to the localized perturbation. The maximum value of all peak folded centroid residuals is at 0.39″ (57 au), its mean value is 39 m/s, and 1σv = 20 m/s (not to be mistaken with the cluster variances).

In the text
thumbnail Fig. C.1.

Azimuthally averaged and deprojected azimuthal, radial, and vertical velocity components. The radial width of each annulus is one-fourth of the synthesized beam size. The error bars are given by the standard deviation for each velocity component averaged over the 20 iterations used.

In the text
thumbnail Fig. C.2.

J1604 deprojected velocity components. It is assumed that all velocities are azimuthal (left column), radial (central column), or vertical (right column). For the azimuthal and radial components, wedges along the minor and major axis have been masked as the observations are insensitive to these components (see Eq. C.1). In each panel the synthesized beam is shown in the lower-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.