<figure><img src="image_1.png" /></figure>

# The Disk Substructures at High Angular Resolution Project (DSHARP). IX. A High-definition Study of the HD 163296 Planet-forming Disk

Andrea Isella<fnref n="1" />, Jane Huang<fnref n="2" />, Sean M. Andrews<fnref n="2" />, Cornelis P. Dullemond<fnref n="3" />, Tilman Birnstiel<fnref n="4" />, Shangjia Zhang<fnref n="5" />, Zhaohuan Zhu<fnref n="5" />, Viviana V. Guzmán<fnref n="6" />,<fnref n="7" />, Laura M. Pérez<fnref n="8" />, Xue-Ning Bai<fnref n="9" />, Myriam Benisty<fnref n="10" />,<fnref n="11" />, John M. Carpenter<fnref n="6" />, Luca Ricci<fnref n="12" />, and David J. Wilner<fnref n="2" />

## Abstract

The Atacama Large Millimeter/submillimeter Array observations of protoplanetary disks acquired by the Disk Substructure at High Angular Resolution Project resolve the dust and gas emission on angular scales as small as 3 astronomical units, offering an unprecedented detailed view of the environment where planets form. In this Letter, we present and discuss observations of the HD 163296 protoplanetary disk that imaged the 1.25 mm dust continuum and ¹²CO J = 2–1 rotational line emission at a spatial resolution of 4 and 10 au, respectively. The continuum observations resolve and allow us to characterize the previously discovered dust rings at radii of 68 and 100. They also reveal new small-scale structures, such as a dark gap at 10 au, a bright ring at 15 au, a dust crescent at a radius of 55 au, and several fainter azimuthal asymmetries. The observations of the CO and dust emission provide information about the vertical structure of the disk and allow us to directly constrain the dust extinction optical depth at the dust rings. Furthermore, the observed asymmetries in the dust continuum emission corroborate the hypothesis that the complex structure of the HD 163296 disk is the result of the gravitational interaction with yet-unseen planets.

**Key words:** planet–disk interactions – protoplanetary disks – techniques: interferometric

## 1. Introduction

In recent years, millimeter-wave interferometers and near-infrared high-contrast cameras have imaged nearby protoplanetary systems at unprecedented angular resolution in both continuum and line emission (e.g., ALMA Partnership et al. 2015; Sallum et al. 2015; Andrews et al. 2016; Isella et al. 2016; Keppler et al. 2018). Although still limited to a few objects, these observations provide information about the processes responsible for the formation of planets.

Several circumstellar disks observed at a spatial resolution better than 50 au reveal ring-like features in the emission of small and large dust particles, as well as of the molecular gas (Andrews et al. 2011a, 2011b; Isella et al. 2013, 2014; Zhang et al. 2014, 2016; van der Marel et al. 2015, 2016, 2018a; Dong et al. 2017, 2018; Fedele et al. 2017, 2018; Tang et al. 2017; Boehler et al. 2018). In particular, the homogeneous survey performed by the Disk Substructure at High Angular Resolution Project (DSHARP) has revealed that multiple-ring systems are ubiquitous among the most massive circumstellar disks (Andrews et al. 2018). The number of rings and their structure vary substantially from object to object, and in some cases, even within the same disk (Huang et al. 2018).

Myriad theoretical models have been proposed to explain the formation of rings in circumstellar disks. These include the

interaction between the disk and yet-unseen giant planets (e.g., Bryden et al. 1999; Zhu et al. 2014; Dong et al. 2018), sharp opacity variations at gas-solid phase transitions (Banzatti et al. 2015; Zhang et al. 2015; Okuzumi et al. 2016), dust accumulations at the edge of low-viscosity regions (Flock et al. 2015; Miranda et al. 2017), and zonal flows via spontaneous concentration of net vertical flux (Bai & Stone 2014; Béthune et al. 2017; Suriano et al. 2018). However, to date, planet–disk interaction models have been the most successful in explaining the observed structures (Jin et al. 2016; Dong et al. 2018; Liu et al. 2018). Taken at face value, these results tend to suggest the existence of a population of young gas giant planets orbiting at several tens of au from the central star, which challenges current planet formation models. Furthermore, the link between planets and rings is supported by the direct detection of possible young planets in the PDS 70 and LkCa 15 systems (Sallum et al. 2015; Keppler et al. 2018), although the nature of some candidates is debated (Thalmann et al. 2016; Mendigutía et al. 2018).

In this Letter, we present new Atacama Large Millimeter/ submillimeter Array (ALMA) observations of the 1.25 mm dust continuum and ¹²CO J = 2–1 line emission of the circumstellar disk around the Herbig Ae star HD 163296. The observations achieve a spatial resolution of 4 au and

<footnote n="1">¹Department of Physics and Astronomy, Rice University, 6100 Main Street, MS-108, Houston, TX 77005, USA; isella@rice.edu</footnote>

<footnote n="2">²Harvard-Smithsonian Center for Astrophysics, 60 Garden Street, Cambridge, MA 02138, USA</footnote>

<footnote n="3">³Zentrum fr Astronomie, Heidelberg University, Albert Ueberle Str. 2, D-69120 Heidelberg, Germany</footnote>

<footnote n="4">⁴University Observatory, Faculty of Physics, Ludwig-Maximilians-Universität München, Scheinerstr. 1, D-81679 Munich, Germany</footnote>

<footnote n="5">⁵Department of Physics and Astronomy, University of Nevada, Las Vegas, 4505 South Maryland Parkway, Las Vegas, NV 89154, USA</footnote>

<footnote n="6">⁶Joint ALMA Observatory, Avenida Alonso de Córdova 3107, Vitacura, Santiago, Chile</footnote>

<footnote n="7">⁷Instituto de Astrofísica, Pontificia Universidad Católica de Chile, Av. Vicuña Mackenna 4860, 7820436 Macul, Santiago, Chile</footnote>

<footnote n="8">⁸Departamento de Astronomía, Universidad de Chile, Camino El Observatorio 1515, Las Condes, Santiago, Chile</footnote>

<footnote n="9">⁹Institute for Advanced Study and Tsinghua Center for Astrophysics, Tsinghua University, Beijing 100084, People's Republic of China</footnote>

<footnote n="10">¹⁰Unidad Mixta Internacional Franco-Chilena de Astronomía, CNRS/INSU UMI 3386, Departamento de Astronomía, Universidad de Chile, Camino El Observatorio 1515, Las Condes, Santiago, Chile</footnote>

<footnote n="11">¹¹Univ. Grenoble Alpes, CNRS, IPAG, F-38000 Grenoble, France</footnote>

<footnote n="12">¹²Department of Physics and Astronomy, California State University Northridge, 18111 Nordhoff Street, Northridge, CA 91130, USA
*Received 2018 October 5; revised 2018 December 5; accepted 2018 December 7; published 2018 December 26*</footnote>
deliver the sharpest images of this source obtained to date. At a distance of 101 pc (Gaia Collaboration et al. 2018), HD 163296 is in several respects the examplar of disks thought to be perturbed by planets. The star is surrounded by a huge (more than 1000 au in diameter) Keplerian disk whose mass has been estimated to range between 0.01 and 0.15 $M_\odot$ (Isella et al. 2007; Tilling et al. 2012; Muro-Arena et al. 2018). Previous ALMA observations that resolved the gas and dust emission on spatial scales of 20 au revealed the presence of three circular gaps in the disk density with radii of 45, 87, and 140 au, and, correspondingly, three dense rings centered at 68, 100, and 160 au (Isella et al. 2016; Zhang et al. 2016). The recent detection of deviations from Keplerian rotation at the location of the gaps and rings by Teague et al. (2018) confirms that these structures correspond to variations in the gas pressure (i.e., variations in the gas density and/or temperature), as opposed to radial changes in the dust opacity. The comparison between the observations and planet-disk interaction models indicates that the density gaps might have been carved by planets with masses between 0.5 and $1 M_J$ orbiting at 48, 86, and 131 au from the central star (Liu et al. 2018; van der Marel et al. 2018b). However, models that assume very low viscosity ($\alpha < 10^{-4}$) suggest that the observed multiple ring structure may also result from the gravitational interaction with a single planet that is less massive than Saturn at 100 au orbital radius (Dong et al. 2018). Furthermore, the presence of an additional planet has been proposed based on the detection of local deviations from Keplerian velocity. Current models suggest that this planet might be located at 260 au from the star and has a mass of about $2 M_J$ (Pinte et al. 2018).

The HD 163296 system has also been the target of a number of optical and infrared high-contrast imaging campaigns aimed at characterizing the morphology of the disk in scattered-light emission and detecting low-mass companions (Grady et al. 2000; Wisniewski et al. 2008; Benisty et al. 2010; Garufi et al. 2014; Muro-Arena et al. 2018). Whereas these observations revealed rings that are similar to those observed at millimeter wavelengths, they did not detect any stellar or planetary mass object at orbital radii larger than 25 au and down to a mass sensitivity of a few $M_J$ (Guidi et al. 2018). Thanks to the unprecedented angular resolution, our new ALMA observations provide detailed information on the morphology of the HD 163296 disk that help to understand the origin of the observed structures.

The structure of this Letter is as follows. The observations and data reduction are discussed in Section 2. The map of the dust continuum emission and its analysis are presented in Section 3, while CO data are discussed in Section 4. The implications of our results in the context of constraining the temperature of the disk and the origin of the observed structures are presented in Section 5. Finally, the main results of our investigation are summarized in Section 6.

## 2. Observations

The observations used in this study were obtained during the DSHARP ALMA Large Program (2016.1.00484.L) and are discussed in detail in Andrews et al. (2018). In brief, the HD 163296 disk was observed in Band 6 ($\lambda \sim 1.3$ mm) in 2017 September in configuration C40-8, which delivered baselines $B$ between 41 m and 5.8 km, and a theoretical angular resolution $\lambda/B_{max} \sim 0.04$. Archival data from project 2013.1.00366.S (Flaherty et al. 2015) and project

2013.1.00601.S (Isella et al. 2016) were used to improve sensitivity and uv-coverage on shorter baselines.

An initial calibration of each data set was produced by the ALMA pipeline. A visual inspection of the data obtained from project 2013.1.00366.S revealed a number of visibilities characterized by very noisy amplitudes that required manual flagging. Short-baseline observations were independently self calibrated and recentered to account for shifts in the position of the disk due to the proper motion of the source. We compared the flux calibration of each track and adjusted the flux scale so that the source visibility amplitude measured on overlapping baselines agrees within 5%. Finally, we combined all of the tracks and performed phase and amplitude self-calibration on both short- and long-baseline data. The self-calibration procedure resulted in an improvement of 43% in the peak signal-to-noise ratio (S/N) in the dust continuum image. A complete description of the DSHARP data calibration procedure is presented in Andrews et al. (2018), while the CASA script used to calibrate and image HD 163296 data, including the manual flagging and flux rescaling, is available online at https://almascience.org/ alma-data/lp/DSHARP. The 1.25 mm dust continuum emission was imaged using the CASA task `tclean` and a robust parameter of $-0.5$ resulting in a synthesized beam FWHM of $0.038 \times 0.048$, corresponding to a spatial resolution of $3.8 \times 4.8$ au at the distance of the source. The rms noise is $0.023$ mJy beam$^{-1}$, and the peak S/N is 185.

The complex gain solutions of the self-calibration of the continuum emission were then applied to the $^{12}$CO data. The line was observed at the same angular resolution of the continuum ($\sim 0.04$), but we imaged it using Briggs robust = 0.5 (no uv-tapering) to achieve higher S/N (see Andrews et al. 2018 for more details about CO imaging). The FWHM of the synthesized beam of the presented CO maps is $0.104 \times 0.095$, which is about half the beam of the CO observations published in Isella et al. (2016), and about seven times smaller than the resolution of the ALMA science verification data (de Gregorio-Monsalvo et al. 2013; Rosenfeld et al. 2013). Channels are spaced in velocity by $0.32$ km s$^{-1}$, but, due to Hanning smoothing, the velocity resolution is $0.64$ km s$^{-1}$. The velocity grid for HD 163296 slightly differs from the fiducial channel spacing of $0.35$ km s$^{-1}$ adopted for the other DSHARP sources. The rms noise per channel is $0.84$ mJy beam$^{-1}$.

## 3. Dust Continuum Emission

The map of the 1.25 mm continuum emission (Figure 1) features two bright elliptical rings previously reported by Isella et al. (2016), as well as three new morphological features: an arc of emission inside the first ring (inset a), an inner dark gap and bright ring at about 10 au and 15 au, respectively, and a central azimuthal asymmetry (inset b). In this section we discuss the morphology of the dust rings, constrain the radial profile of the continuum emission, and investigate the presence of asymmetric structures in the distribution of the solid material.

### 3.1. Morphology of the Dust Rings

As a first step toward characterizing the structure of the dust rings, we fit the crests and troughs of the continuum emission with circles defined by the position of the center ($\Delta x_0$, $\Delta y_0$),
<figure><img src="image_2.png" /><figcaption>**Figure 1.** Map of the HD 163296 disk recorded in the 1.25 mm continuum. The angular resolution of the observations is 0.′′038 × 0.′′048 and is indicated by the white ellipse in the bottom-left corner of the left panel. At the source distance of 101 pc, the spatial resolution is 3.8 au × 4.8 au. Inset (a) and (b) show a zoom-in view of two asymmetric features revealed by the observations. The rms noise is 0.023 mJy beam⁻¹. The horizontal segments indicate a spatial scale of 10 au.</figcaption></figure>

<table><caption>**Table 1**

Results of Elliptical Fitting to the Bright and Dark Rings Observed in the 1.25 mm Continuum Map</caption><thead><tr><th>Feature<br>(1)</th><th>Δx<sub>0</sub> (mas)<br>(2)</th><th>Δy<sub>0</sub> (mas)<br>(3)</th><th>r (′′)<br>(4)</th><th>r (au)<br>(5)</th><th>i (°)<br>(6)</th><th>PA (°)<br>(7)</th></tr></thead><tbody><tr><td>D10</td><td>-6.2 ± 0.8</td><td>10.5 ± 0.7</td><td>0.0987 ± 0.0009</td><td>9.96 ± 0.07</td><td>37.9 ± 1.2</td><td>127.90 ± 1.9</td></tr><tr><td>B14</td><td>-4.8 ± 0.7</td><td>11.4 ± 0.7</td><td>0.1430 ± 0.0010</td><td>14.44 ± 0.07</td><td>47.24 ± 0.64</td><td>131.10 ± 0.88</td></tr><tr><td>D45</td><td>-5.4 ± 2.0</td><td>5.2 ± 1.9</td><td>0.4433 ± 0.0026</td><td>44.77 ± 0.19</td><td>42.22 ± 0.67</td><td>133.67 ± 0.95</td></tr><tr><td>B67</td><td>-5.3 ± 1.0</td><td>7.3 ± 1.1</td><td>0.6633 ± 0.0014</td><td>66.99 ± 0.11</td><td>46.78 ± 0.21</td><td>133.13 ± 0.29</td></tr><tr><td>D86</td><td>6.7 ± 2.3</td><td>-0.3 ± 2.3</td><td>0.8575 ± 0.0028</td><td>86.61 ± 0.22</td><td>47.34 ± 0.32</td><td>132.78 ± 0.45</td></tr><tr><td>B100</td><td>-2.3 ± 0.9</td><td>8.6 ± 0.9</td><td>0.9870 ± 0.0011</td><td>99.69 ± 0.08</td><td>46.59 ± 0.11</td><td>133.46 ± 0.15</td></tr><tr><td>D141</td><td>-0.2 ± 9.9</td><td>6.4 ± 9.5</td><td>1.3923 ± 0.0126</td><td>140.62 ± 0.96</td><td>47.2 ± 0.9</td><td>131.3 ± 1.2</td></tr><tr><td>B159</td><td>-16.2 ± 12.0</td><td>1.1 ± 12.0</td><td>1.572 ± 0.016</td><td>158.7 ± 1.2</td><td>45.7 ± 1.0</td><td>132.0 ± 1.4</td></tr><tr><td>mean<sup>a</sup></td><td>-3.5 ± 0.5</td><td>9.1 ± 0.5</td><td>...</td><td>...</td><td>46.7 ± 0.1</td><td>133.3 ± 0.1</td></tr></tbody></table>

**Note.**

<sup>a</sup> The parameters of D10 and D45 gaps were not used to calculate the mean values.

the radius *r*, the inclination *i* at which they are observed, and the position angle (PA) of their projected major axis. The fit is performed as discussed in Huang et al. (2018). In brief, we locate the points (x, y) corresponding to the radial maxima (or minima) of each ring (or gap), and use the python *emcee* package to calculate the ellipses that best reproduce this set of points. The best-fit parameters are listed in Table 1. Bright rings and dark gaps are labeled with the suffixes B and D, respectively, while the number appearing in the feature name corresponds to their radius in au. The uncertainties on the best-fit values correspond to the 16th and 84th percentile of the marginalized posterior probability distribution.

The outermost features B67, D86, B100, D145, and B155, are well described by concentric circles with an average inclination of (46.7 ± 0.1)° and average PA of (133.3 ± 0.1)°.

Taken at face value, the circular model fitting indicates that D10 and D45 have lower inclinations. We believe that the lower inclination of D45 is due to the presence of the dust crescent centered at 55 au, which partially fills the gap making it appear more circular (i.e., less inclined), while the lower inclination of D10 is compatible with the effect of beam smearing, which is discussed in more details below.

Figure 2 shows the deprojected map and the azimuthally averaged radial profile of the continuum emission obtained using the values for the disk inclination and PA derived from the outermost rings. In the deprojected map, the outer rings appear as vertical stripes, confirming that these structures do indeed have a similar inclination and PA, and that they are intrinsically circular in shape. The intensity ratio between B14 and D10 varies with the azimuthal angle reaching a maximum
<figure><img src="image_3.png" /><figcaption>**Figure 2.** Left panel: polar map of the 1.25 mm dust continuum emission deprojected using a disk inclination of 46°.7. Two different color scales are used to highlight the dust intensity structures observed inside and outside of 40 au. The angle $\theta$ indicates the elongation from the disk apparent minor axis measures positive east of north. The top and bottom white bars indicate the position of the bright rings and dark gaps, respectively. Middle panel: azimuthal average of dust continuum emission expressed in units of brightness temperature $T_b^{pl}$. Right panel: the same as in the middle panel but plotted on a logarithmic temperature scale. The blue curves indicate the dispersion around the mean value, while the red curve indicates the temperature corresponding to 3× the rms noise of the azimuthally averaged intensity calculated, taking into consideration both the rms noise of the continuum map and the number of resolution elements in each radial bin. Vertical dashed lines indicate the position of the dark rings, while vertical solid lines indicate the position of bright rings.</figcaption></figure>

of about 1.2 along the disk major axis ($\theta = \pm 90^{\circ}$) and a minimum of about 1 along the disk minor axis ($\theta = 0^{\circ}$). This variation can be explained by the fact that the effective spatial resolution for an inclined disk is higher along the disk major axis compared to the disk minor axis. In practice, B14 and D10 are spatially resolved only along the disk major axis, while they blend together along the disk minor axis.

Azimuthal variations in the contrast between bright rings and dark gaps are also visible for the outer pairs. For example, along the westward semimajor axis ($\theta = -90^{\circ}$), the intensity reaches a minimum that is consistent with the noise of the observations at the center of D45, and a maximum of 0.91 mJy beam$^{-1}$ at the center of the adjacent bright ring B67. This corresponds to an intensity ratio larger than 37. In the opposite direction, i.e., along the eastward disk semimajor axis ($\theta = 90^{\circ}$), the intensity ratio B67/D45 is 13. As comparison, the intensity ratio along the southward ($\theta = 180^{\circ}$) and northward ($\theta = 0^{\circ}$) direction of the disk minor axis is about 7 and 5, respectively (see Table 2 for a summary of the intensity ratio among rings pairs). The difference between intensity ratios measured along the disk major and minor axes is likely due to beam smearing, while the difference between the intensity ratios measured along the disk major axis suggests the presence of azimuthal asymmetries in the dust emission.

The center and right panels of Figure 2 show the azimuthally averaged profile of the deproject continuum emission expressed in units of brightness temperature $T_b^{pl}$ calculated using the full Planck equation as

$$T_b^{pl} = \frac{h\nu}{k_b} \left[ \ln \left( \frac{2h\nu^3\Omega}{c^2F_\nu} + 1 \right) \right]^{-1}, \quad (1)$$

where the $F_\nu$ is the flux density integrated over a beam defined by the solid angle $\Omega = \pi\theta_{\min}\theta_{\max}/(4 \ln 2)$, where $\theta_{\min}$ and $\theta_{\max}$ are the minimum and maximum FWHM of the synthesized clean beam. The brightness temperature varies between about 50 K (corresponding to an intensity of 3.8 mJy beam$^{-1}$) in the innermost disk regions down to about 3 K (or a flux density of 0.0025 mJy beam$^{-1}$) at about 180 au.

<table><caption>**Table 2**
Intensity Ratios between Bright Rings and Dark Gaps along Different Radial Directions</caption><thead><tr><th>Pair<br>(1)</th><th>+90°<br>(Maj E)<br>(2)</th><th>-90°<br>(Maj W)<br>(3)</th><th>0°<br>(Min N)<br>(4)</th><th>180°<br>(Min S)<br>(5)</th><th>Mean<br>(6)</th></tr></thead><tbody><tr><td>B14/D10</td><td>1.2</td><td>1.1</td><td>1.0</td><td>0.9</td><td>1.1</td></tr><tr><td>B67/D45</td><td>13.3</td><td>&gt;37</td><td>6.6</td><td>4.7</td><td>&gt;17</td></tr><tr><td>B100/D86</td><td>11.7</td><td>7.1</td><td>4.5</td><td>2.5</td><td>7.0</td></tr><tr><td>B159/D141</td><td>2.1</td><td>1.1</td><td>1.0</td><td>0.5</td><td>1.3</td></tr></tbody></table>

**Note.** (1) Name of the pair. (2), (3) Intensity ratio along the apparent disk major axis in the eastward ($\theta = +90^{\circ}$) and westward ($\theta = -90^{\circ}$) directions, respectively. (4), (5) Intensity ratio along the apparent disk minor axis in the northward ($\theta = 0^{\circ}$) and southward ($\theta = -180^{\circ}$) directions, respectively. (6) Intensity ratio of the azimuthally averaged intensity profile.

The pairs B67/D45 and B100/D86 have brightness temperature ratios of 4.3 and 2.4, respectively. Note that the difference between brightness temperature and intensity ratios is due to the nonlinearity of Equation (1).

## 3.2. Width of the Continuum Rings

Dullemond et al. (2018) finds that many of the continuum rings revealed by DSHARP observations, including the rings B67 and B100 of the HD 163296 disk, have a radial width that is narrower than the estimated pressure scale height of the gaseous disk, and conclude that this is strong evidence that these rings have formed by radial migration of dust particles toward local maxima of the gas pressure. Dullemond et al. (2018) defines the width of a ring as the dispersion of a Gaussian function that reproduces the deprojected and azimuthally averaged radial profile of the dust emission across the ring itself. This modeling approach is very fast, but it does not account for the observational noise, nor for the discrete sampling of the uv-plane performed by an interferometer, or for the fact that the synthesized beam is not circular.

Here, we compare the widths of the HD 163296 rings measured by Dullemond et al. (2018) with those estimated by a
<table><caption>**Table 3**

Results of Gaussian Model Fitting of the Dust Continuum Emission</caption><thead><tr><th rowspan="2">Feature<br>(1)</th><th colspan="3">Radial Intensity Profile</th><th colspan="3">Visibilities</th></tr><tr><th>r (au)<br>(2)</th><th>w<sub>d</sub> (au)<br>(3)</th><th>I (Jy/as<sup>2</sup>)<br>(4)</th><th>r (au)<br>(5)</th><th>w<sub>d</sub> (au)<br>(6)</th><th>I (Jy/as<sup>2</sup>)<br>(7)</th></tr></thead><tbody><tr><td>...</td><td>...</td><td>...</td><td>...</td><td>3.8 ± 0.1</td><td>2.3 ± 0.1</td><td>1.7 ± 0.2</td></tr><tr><td>B14</td><td>...</td><td>...</td><td>...</td><td>15.5 ± 0.2</td><td>8.7 ± 0.2</td><td>1.1 ± 0.1</td></tr><tr><td>...</td><td>...</td><td>...</td><td>...</td><td>31.8 ± 0.3</td><td>4.4 ± 0.1</td><td>0.62 ± 0.06</td></tr><tr><td>B67</td><td>67.7</td><td>6.84</td><td>0.38</td><td>67.08 ± 0.03</td><td>6.56 ± 0.05</td><td>0.38 ± 0.02</td></tr><tr><td>B100</td><td>100.0</td><td>4.67</td><td>0.24</td><td>101.16 ± 0.04</td><td>5.8 ± 0.1</td><td>0.20 ± 0.01</td></tr></tbody></table>

**Note.** Columns 2–4 list the best-fit parameters for the radius, width, and intensity, respectively, of the Gaussian components used to model the azimuthally averaged profile of the dust continuum emission as in Dullemond et al. (2018). Columns 5 to 7 list the same quantities obtained by modeling the observations in the Fourier space. Additional parameters not reported in the table are the disk inclination $i = 46^{\circ}7 \pm 0^{\circ}1$, the disk position angle PA = $132^{\circ}8 \pm 0^{\circ}1$, and the offset of the center of the emission relative to the phase center of the observations $\delta x_0 = -8.7 \pm 0.1$ mas, $\delta y_0 = -3.0 \pm 0.1$ mas.

more accurate, but much slower, modeling of the continuum visibilities in the uv domain. We start by assuming that the dust emission is axisymmetric and expressed by a sum of n Gaussian functions

$$I(r) = \sum_{i=1}^{n} I_i e^{-(r-r_i)^2/2w_{d,i}^2}, \quad (2)$$

We generate synthetic 2D images of the continuum emission that are inclined and rotated to simulate all possible disk viewing angles. The center of emission is allowed to vary with respect to the center of the image (i.e., the phase center). Synthetic visibilities ($V^{\text{mod}}$) are then calculated by taking the Fourier transform of the image at the spatial frequencies of the observed visibilities ($V^{\text{obs}}$). The goodness of the model fit is evaluated through a usual $\chi^2$ test, where $\chi^2 = \sum w(V^{\text{obs}} - V^{\text{mod}})^2$, and $w$ is the weight of each visibility measurement provided by the ALMA pipeline. The calculation of $\chi^2$ was performed using the CPU-based version of python package galario (Tazzari et al. 2018). Finally, the $\chi^2$ is used to sample the posterior likelihood using the python package emcee (Foreman-Mackey et al. 2013).

An initial exploration of the parameter space indicates that a minimum of five Gaussian components are required to reproduce the observed continuum emission. In total, the model therefore has 19 free parameters: the coordinates $(x_0, y_0)$ of the center of the emission, the disk inclination $i$ and PA, and the radius $r$, peak intensity $I$, and width $w_d$ of five Gaussian rings. We run *emcee* using 200 walkers initialized around the parameter values obtained by a Gaussian fitting of the radial intensity profile as in Dullemond et al. (2018). Each walker is evolved for an initial burn-in run of $10^4$ steps, after which they all converge toward a stable set of parameters. Following Goodman & Weare (2010), we sample the posterior distribution by letting *emcee* run until convergence, which is established based on the mean autocorrelation time of all the model parameters. We find that convergence is achieved after about $4 \times 10^4$ steps. The optimal model parameters are then calculated as the median of the marginalized probability distribution and are listed in Table 3, while the related uncertainties correspond to the 0.3th and 99.7th percentile of the marginalized posterior probability distribution. This would be equivalent to a $3\sigma$ uncertainty for a Gaussian posterior distribution. A comparison between the image of the residual and the map of the continuum emission is presented in the Appendix.

For B67, we find a width of ~6.6 au, that is very close to the 6.8 au derived by Dullemond et al. (2018). Conversely, the width of B100 is about 25% larger than that measured by fitting the intensity profile. As noted in Figure 3 of Dullemond et al. (2018), the intensity radial profile of B100 has a peaked central region and bright wings that are not well reproduced by a Gaussian function. The value of the width in Dullemond et al. (2018) results from fitting only the innermost part of the ring between 94 and 104 au, while our value resulted from modeling the global profile including the bright wings of the ring. The fact that the two modeling procedures achieve the same result for rings that are intrinsically Gaussian (i.e., B67) indicates that the results of Dullemond et al. (2018) are independent of the exact shape of the beam, or, more generally, from nonlinear effects caused by the discrete uv-sampling and image deconvolution. In practice, the uv-sampling of the observation is sufficiently uniform and the S/N of the continuum map is sufficiently high to allow to perform the data analysis in the image domain without loss of accuracy. This is a great advantage because modeling the dust rings in the image domain takes only a few seconds, compared to the several days, or weeks, required to perform the model fitting in the Fourier domain.

Our analysis indicates that the ring B14 has a radius of about 15.5 au, slightly larger than the value calculated by fitting the crest of the emission as discussed in Section 3.1 (see also Table 1). The best-fit model also includes a small ring of 4 au in radius, therefore suggesting the presence of an inner circular cavity like those observed in transitional disks. However, as discussed in the next section, our azimuthally symmetric best-fit model poorly matches the dust emission from the innermost disk regions, which instead seems to indicate the presence of a dust crescent.

Finally, the last component of the model is a ring with a radius of about 32 au that is required to reproduce the kink in the radial profile of the emission observed at about 30 au from the center. Similar kinks in the radial profile of the dust continuum emission are observed in other DSHARP sources and are discussed in Huang et al. (2018).

### 3.3. Asymmetries of the Continuum Emission

The map obtained by imaging the residual visibilities obtained by subtracting the best-fit Gaussian model from the observations (Figure 3) reveals several features with absolute intensities above 5 times the rms noise and peak intensities as high as 30 times the rms noise. A comparison between the
<figure><img src="image_4.png" /><figcaption>**Figure 3.** Top panels: Cartesian and polar maps of the residual intensity obtained by modeling the continuum image of HD 163296 with a symmetric disk model comprising Gaussian rings as discussed in the text. Dashed and solid contours correspond to ±5× and ±10× the noise level, respectively. Bottom panels: azimuthal profile of the residual intensity along the crescents C4 and C55 (left) and the rings B67 and B100 (right). As a reference, the brown dotted curves indicate the azimuthal profile of the intensity at a radius of 160 au, which is indicative of the noise level of the observations.</figcaption></figure>

image of the residual and the continuum map is shown in the Appendix. The two most prominent structures have the shape of crescents and were already introduced at the beginning of Section 3 (see panels (b) and (c) of Figure 1). The innermost crescent is located at about 4 au from the center and it overlaps with the innermost ring component of the Gaussian model discussed above. Given its azimuthally asymmetric structure, and following the nomenclature adopted for the rings and gaps, we label this feature as C4 (C meaning “crescent”). The peak of C4 has a PA of -44° and an intensity is 0.55 mJy beam⁻¹, or about 20% of the continuum intensity measured at the same position. C4 extends by about 180° in azimuth and by about one resolution element in radius. C4 resembles the crescents observed in other Herbig Ae disks that have been attributed to the presence of anticyclonic vortices (e.g., Casassus et al. 2013; Isella et al. 2013; van der Marel et al. 2013; Boehler et al. 2018) or to an optically thick inner disk warped with respect to the outer disk (e.g., Marino et al. 2015; Benisty et al. 2017). However, higher angular resolution observations are required

to investigate the morphology of HD 163296 innermost disk regions in greater details.

The second crescent, labeled as C55, is centered at an orbital radius of about 55 au and PA of 99°. The peak intensity is 0.64 mJy beam⁻¹, corresponding to a S/N of 28. Modeling C55 intensity as a 2D Gaussian function in the polar plane $I(r, θ) ∝ e^{-(r-r_c)^2/2σ_r^2+(θ-θ_c)^2/2σ_θ^2}$, we find $σ_θ = 17$ au and $σ_r = 3$ au. Accounting for beam convolution (Equation (3)), the intrinsic radial and azimuthal widths of C55 are ~2.2 au and ~16.9 au, respectively.

In addition to C4 and C55, the image of the residuals shows fainter asymmetries along the B67 and B100 dust rings. If the dust emission is optically thin, intensity variations along the rings might probe the clumpiness of the dust distribution, while if the emission is optically thick, they reveal information about variations in the dust temperature. The intensity profiles along B67 and B100 are shown in the bottom-right panel of Figure 3. In the case of B67, the residual intensity varies between ±0.15 mJy beam⁻¹, to be compared with an azimuthally
averaged intensity value measured at the same radius of 1.0 mJy beam⁻¹. The azimuthal intensity profile of B67 shows small amplitude variations characterized by an angular scale of about 4°, and larger variations with angular scales of 30°–50°. Whereas the small-scale variations have the same size of the synthesized beam and are most likely caused by the noise of the observations, the larger-scale structures might trace real asymmetries in the dust emission and, as the emission from B67 is estimated to be optically thin (see Section 5.1), in the dust density.

## 4. ¹²CO J = 2–1 Line Emission

Figure 4 shows continuum-subtracted channels maps of the ¹²CO (2–1) emission (hereafter referred simply as CO emission). Below we describe the morphology and kinematics of the CO emission, and use the observations to constrain disk properties such as the disk temperature and the optical depth of the dust continuum emission.

### 4.1. Morphology of the CO Emission

The CO emission extends far beyond the outer edge of millimeter-wave continuum emission up to an angular distance of 5"/6, or about 560 au, from the central star. As discussed in Isella et al. (2016), the discrepancy between the radial extent of the millimeter wave dust and CO emission indicates a sharp drop in either the dust opacity or dust column density beyond about 200 au.

The CO emission observed at intermediate velocities (e.g., v<sub>lsr</sub> = 3.40–5.00 km s⁻¹ and v<sub>lsr</sub> = 6.28–8.20 km s⁻¹) shows a characteristic “dragonfly” structure tracing both the front and rear warm CO-emitting layers of the disk (Figure 5). The emission coming from the rear of the disk is fainter than that from the front due to a combination of lower temperature of the emitting gas and absorption from the intervening material. As discussed in Rosenfeld et al. (2013), the fact that the front and rear sides of the disk appear as two distinct emitting regions is explained by the large optical depth of the CO line (whose intensity therefore mostly depends on the gas temperature) and a vertical gradient of the disk temperature characterized by warm surfaces and a cold disk midplane. Furthermore, the separation between the front and rear CO emission measures the vertical geometry of the CO-emitting surfaces.

We compare the observed geometry of the CO emission with that of a parametric model where the CO line originates from a geometrically thin layers characterized by a distance from the midplane z<sub>co,0</sub> = ±z<sub>co,0</sub> (r/r<sub>0</sub>)<sup>q</sup>. As a demonstration of the technique, we show in Figure 5 parametric isovelocity curves corresponding to the front side of the disk calculated as presented in the Appendix. The geometry of the isovelocity curves depends on z<sub>co,0</sub>, q, as well as on the stellar mass and disk inclination. Assuming the disk inclination derived from the continuum (i = 46°/7) and a stellar mass of 2.0 M<sub>☉</sub>, we obtain a good match between model and observations for q = 0.5 and z<sub>co</sub>(100 au) = 30 au.

Our observations reveal that the CO emission coming from the back side of the disk is dimmer at the location of the bright dust rings B67 and B100 compared to the emission observed at the location of the dark gaps D45 and D87 (Figure 6). This dimming, which is likely due to absorption of the line emission as it crosses the dusty rings, provides a tool to constrain the optical depth of continuum emission independently from the

assumptions on the dust temperature. The CO emission coming from the back side of the disk can be expressed as I<sub>co</sub><sup>obs</sup> = I<sub>co</sub>e<sup>-τ<sub>d</sub><sup>ext</sup></sup>, where τ<sub>d</sub><sup>ext</sup> is the line-of-sight dust extinction optical depth equal to the sum of the absorption and scattering optical depths (τ<sub>d</sub><sup>ext</sup> = τ<sub>d</sub><sup>abs</sup> + τ<sub>d</sub><sup>sca</sup>), while I<sub>co</sub> is the initial unattenuated CO emission. Defining the parameter ϵ as ϵ ≡ τ<sub>d</sub><sup>abs</sup>/τ<sub>d</sub><sup>ext</sup>, the absorption optical depth of the dust continuum can be written as

$$ \tau_d^{\text{abs}} = \epsilon \ln \frac{I_{\text{co}}}{I_{\text{co}}^{\text{obs}}}. \qquad (3) $$

The observations of the CO emission directly measure I<sub>co</sub><sup>obs</sup> along the bright continuum rings, while the unattenuated emission I<sub>co</sub> can be estimated through the procedure illustrated in Figure 6. In practice, in the channels where CO absorption is observed, we measure the intensity of the CO emission along the isovelocity curves corresponding to the back side emission. We then use the values measured within the dust gaps, where, based on the results of the previous section, τ<sub>d</sub> ≪ 1, to estimate I<sub>co</sub> at the position of the bright rings. We measure the ratio I<sub>co</sub>/I<sub>co</sub><sup>obs</sup> for B67 and B100 along 18 different lines of sight (i.e., 18 different azimuthal angles along the dust rings) and find mean values of 1.9 ± 0.1 and 2.1 ± 0.1. We also perform this measurement for B168 finding I<sub>co</sub>/I<sub>co</sub><sup>obs</sup> = 1.1 ± 0.2. The corresponding dust extinction optical depths are τ<sub>d</sub><sup>ext</sup>(B67) = 0.64 ± 0.05, τ<sub>d</sub><sup>ext</sup>(B100) = 0.74 ± 0.05, τ<sub>d</sub><sup>ext</sup>(B168) = 0.1 ± 0.2. It is worth noting that the ratio I<sub>co</sub>/I<sub>co</sub><sup>obs</sup> could vary with the azimuthal angle if τ<sub>d</sub><sup>ext</sup> is not constant along the dust rings. However, the sensitivity of our observations does not allow us to perform such an analysis, and we must instead rely on the azimuthally averaged value of the CO intensity ratio to investigate the dust extinction properties.

As from Equation (3), the absorption optical depth of the dust rings depends on the dust scattering properties through the parameter ϵ. For dust grains that are much smaller than the wavelength of the observations (in our case, for particle sizes a ≪ λ/2π ~ 0.2 mm), scattering is negligible compared to absorption, and ϵ ≃ 1. For grain sizes that are comparable to the wavelength of the observations (a ~ λ/2π ~ 0.2 mm), scattering, and absorption properties of spherical grains follow Mie theory, and ϵ can vary between about 0.1 and 1, depending on the dust composition and internal structure. Finally, if a ≫ 0.2 mm, scattering and absorption are equal, and ϵ = 0.5 (see, e.g., Bohren & Huffman 1998). Without any prior information about the size and composition of dust grains, the observed absorption of the CO emission constrains the absorption optical depth of B67 and B100 between about 0.07 and 0.7, while the lack of CO absorption at the location of B168 sets an upper limit to the absorption optical depth of this ring at about 0.3. A caveat of this approach is that the emission from the front side of the disk could partially contaminate the emission from the back side, particularly close to the star where the separation between the two surfaces is small. Also, the ratio I<sub>co</sub>/I<sub>co</sub><sup>obs</sup> could be affected by beam smearing, implying that our measurements provide only a lower limit for the true contrast between the initial and absorbed CO emission, and consequently, the corresponding dust optical depth. However, the angular resolution of the CO observations is sufficiently high to believe that this effect is small.

If the absorption optical depth of the dust rings is measured from CO absorption, the corresponding continuum intensity I<sub>ν</sub> can be used to estimate the temperature T<sub>d</sub> of the emitting dust.
<figure><img src="image_5.png" /><figcaption>**Figure 4.** Map of the ¹²CO J = 2–1 line emission recorded toward HD 163296. The FWHM of the synthesized beam is 0.′′104 × 0.′′095, corresponding to a spatial scale of about 10 au at the distance of the source. Each panel has a size of 8″ × 8″, and is labeled with the velocity relative to the local standard of rest (v<sub>lsr</sub>). Channels are spaced by 0.32 km s<sup>-1</sup>. The color scale of the map is linear.</figcaption></figure>
<figure><img src="image_6.png" /><figcaption>**Figure 5.** Zoom-in of four of the CO channels shown in the previous figure (top left: $v_{lsr} = 4.68$ km s$^{-1}$, top right: $v_{lsr} = 6.92$ km s$^{-1}$, bottom left: $v_{lsr} = 3.4$ km s$^{-1}$, bottom right: $v_{lsr} = 8.2$ km s$^{-1}$). Blue dots indicate the isovelocity curves relative to the front CO-emitting layer.</figcaption></figure>

If dust scattering is negligible ($\tau_d^{\text{ext}} = \tau_d^{\text{abs}}$ and $\epsilon = 1$), the conversion between flux density and dust temperature is

$$T_d = \frac{h\nu}{k_b} \left[ \ln \left( \frac{2h\nu^3\Omega(1 - e^{-\tau_d^{\text{abs}}})}{c^2 F_\nu} + 1 \right) \right]^{-1}, \quad (4)$$

or, using Equation (3), as

$$T_d = \frac{h\nu}{k_b} \left[ \ln \left( \frac{2h\nu^3\Omega(1 - I_{\text{co}}^{\text{obs}}/I_{\text{co}})}{c^2 F_\nu} + 1 \right) \right]^{-1}. \quad (5)$$

From the measured continuum intensities we obtain dust temperatures of 24 K and 15 K for B67 and B100, respectively, while the upper limits for the optical depth of B168 translates in a minimum temperature of about 6 K.

If dust scattering is not negligible ($\tau_d^{\text{ext}} > \tau_d^{\text{abs}}$ and $\epsilon < 1$), the dust temperature can be estimated using the formalism

presented in Birnstiel et al. (2018). In this case, the intensity is written as $I_{\nu,d} \sim S(T_d, \tau_{\nu,d}^{\text{ext}}, \epsilon)(1 - e^{-\tau_{\nu,d}^{\text{ext}}})$, where the source function $S_\nu$ depends on the dust temperature, and on the extinction and scattering optical depth along the line of sight. Adopting $\epsilon$ in the range between 0.1 and 1, we obtain dust temperatures between 130 K and 24 K for B67, and between 70 K and 15 K for B100. The fact that the dust temperature depends on the dust-scattering properties hampers the application of this technique as an independent measure of disk temperature. However, if the dust temperature is measured through other means, then the dust-scattering properties can by constrained (see Section 5).

## 4.2. Gas Temperature

At typical densities and temperatures of protoplanetary disks, the peak of the $^{12}$CO J = 2–1 line emission is optically thick and the CO rotational levels are in local thermodynamic
<figure><img src="image_7.png" /><figcaption>**Figure 6.** Illustration of the technique used to measure the extinction optical depth of the dusty rings from observations of the CO emission. Panel (a): blue and red colors correspond to CO emission recorded at $v_{lsr} = 4.36$ and $v_{lsr} = 7.24$ km s$^{-1}$, respectively. Dashed ellipses indicate the position of the rings B67, B100, and B160. Solid lines indicate the isovelocity curves corresponding to the emission coming from the back side of the disk. Panel (b): filled circles with error bars indicate the intensity of the CO line measured along the isovelocity curve. The indices 1–6 on the x-axis correspond to the positions 1–6 along the curve shown in the right panel. The blue dashed line indicates the unattenuated CO emission.</figcaption></figure>

<figure><img src="image_8.png" /><figcaption>**Figure 7.** Left panel: peak intensity map of the $^{12}$CO J = 2–1 line emission. Right panel: radial profile of the CO temperature deprojected accounting for the geometry of the emitting region as discussed in the text. Gray points show the scatter of the temperature measurements. The blue solid curve corresponds to the mean, while the dotted curve indicate the standard deviation from the mean.</figcaption></figure>

equilibrium (see, e.g., Weaver et al. 2018). The intensity at the peak of the line therefore measures the temperature of the emitting gas. Following Weaver et al. (2018), we generate a CO temperature map from non-continuum-subtracted line emission to avoid underestimating the line intensity along lines of sight where the continuum is optically thick (Figure 7). The CO temperature reaches values as high as 120 K in the innermost disk regions and as low as 20 K at the disk outer edge. The peak intensity map is characterized by an evident asymmetry with respect to the apparent disk major axis: the southwest side of the disk is hotter than the northeast side. This is explained in first approximation by the fact that the CO-emitting layer is a conic surface observed from a direction

inclined with respect to the axis of the cone. The orientation of the asymmetry implies that the northern side of the disk is the closest to the observer (see also Rosenfeld et al. 2013).

Assuming, as in the previous section, that the peak of the line is emitted from a geometrically thin layer at a distance from the disk midplane given by $z_{co} = z_{co,0}(r/r_0)^q$, it is possible to deproject the observed CO temperature and estimate its profile as a function of the cylindrical radius $r$. Using the values for $z_{co}$ and $q$ presented above, we derive the CO temperature profile shown in the right panel of Figure 7. Between 30 and 500 au, the radial profile of the CO temperature scales almost linearly with the radius and follows the relation $T_{co}$ (K) $\sim$ 87–0.14 (r/au). Within 30 au, the CO brightness temperature drops. This is
<figure><img src="image_9.png" /><figcaption>**Figure 8.** Left panel: intensity-weighted velocity (moment I) of the ¹²CO J = 2–1 line emission calculated by adopting an intensity threshold of three times the noise level. Right panel: position–velocity diagram obtained by averaging along the whole apparent disk minor axis. The position is the angular distance from the star along the projected disk major axis. The gray dashed line indicates the systemic velocity of 5.8 km s⁻¹. The solid blue curve corresponds to the maximum velocity along the line of sight for a Keplerian disk inclined by 46°7 orbiting a 2.0 M☉ star at the distance of 101 pc.</figcaption></figure>

consistent with the expected drop in brightness temperature as the emitting area becomes smaller than the angular resolution of the observations (beam dilution).

### 4.3. Kinematics of the CO Emission

The overall kinematics of the CO gas is well described by Keplerian rotation around a star with a mass of about 2.0 M☉ (Figure 8). Both the integrated velocity map and the position–velocity diagram are in first approximation symmetric with respect to the apparent disk minor axis. The systemic velocity inferred from CO observations is 5.8 km s⁻¹ and is in agreement with previous results. The position–velocity diagram also reveal diffuse emission at a velocities between about 12.5 and 14 km s⁻¹ that might trace a redshifted stellar wind (Klaassen et al. 2013) or perhaps emission from the leftover of the HD 163296 parent cloud.

The CO map registered at a velocity of 6.92 km s⁻¹ (Figures 4 and 5) is characterized by a kink at δR.A. ~ −1", δdecl. ~ 1"5. This feature was previously reported by Pinte et al. (2018) and has been attributed to perturbations in the gas kinematics caused by a planet with a mass of 2 M,J orbiting at 260 au from the star. Unfortunately, despite the better angular resolution, the observations obtained using the extended ALMA configuration achieve a spectral resolution about 10 times worse than that of those presented by Pinte et al. (2018), and do not allow us to investigate the origin of this feature. Similarly, the coarse velocity resolution of our data might not allow to investigate the small deviations from Keplerian rotation observed within the dust gaps and reported by Teague et al. (2018).

## 5. Discussion

The two main results obtained so far are that (i) the combination of continuum and CO observations directly inform about the temperature and geometry of the HD 163296 disk and (ii) the analysis of the continuum emission constrains radial profile and level of asymmetry of the observed dust rings. In

this section, we first compare the disk temperature inferred from the observations with current models for HD 163296 and, second, we discuss the possible origin of the morphology of the gas and dust emission.

### 5.1. Disk Temperature and Dust Scattering

To date, the HD 163296 disk has been the subject of several studies aimed at constraining the density and thermal structure of the disk through forward modeling of the spatially resolved millimeter-wave molecular line and continuum emission (e.g., Isella et al. 2007, 2016; Tilling et al. 2012; Rosenfeld et al. 2013; Flaherty et al. 2015). This approach relies on a large number of assumptions including the choice of the parametric forms for the gas and dust surface density and temperature (when the latter is not calculated using radiative transfer models), the dust and molecular abundances, the dust opacities, the stellar properties (mass, luminosity, distance), and the gas kinematics (e.g., Keplerian motion, turbulent velocity). Because of the large computational time required to generate synthetic images at high spatial and spectral resolution, the comparison between models and observations is generally performed by either manually adjusting the model parameters until obtaining a satisfactory result (e.g., Isella et al. 2016), or by using automatic algorithms to search for best-fit models among a subset of model parameters. Not surprisingly, forward modeling performed by different investigators lead to different results.

Figure 9 compares the temperature inferred from the 1.25 mm dust continuum intensity to the models discussed in Isella et al. (2016), Rosenfeld et al. (2013), and Flaherty et al. (2015). We choose these three studies because they assume a similar parametric form for the disk temperature and are based on high spectral and angular resolution observations of the ¹²CO line emission. If scattering from dust grains is negligible ($\tau_{d}^{\text{ext}} = \tau_{d}^{\text{abs}}$, $\epsilon = 1$), the brightness temperature $T_{b}^{\text{pl}}$, the physical dust temperature $T_{d}$, and the optical depth $\tau_{d}^{\text{abs}}$ are linked by the relation $B_{\nu}(T_{b}^{\text{pl}}) = B_{\nu}(T_{d})(1 - e^{-\tau_{d}^{\text{abs}}})$. Within 30 au, the model temperatures are bracketed between Rosenfeld
<figure><img src="image_10.png" /><figcaption>**Figure 9.** Left panel: dust temperature as a function of the orbital radius. The black and blue dashed region indicate the continuum brightness temperature and the standard deviation from the mean, respectively. Green, red, and brown curves indicate the disk midplane temperature predicted by Flaherty et al. (2015), Rosenfeld et al. (2013), and Isella et al. (2016), respectively. The vertical lines show the range of dust temperatures compatible with the observed absorption of the CO emission arising from the back side of the disk at the location of the bright rings B67, B100, and B168. The horizontal segments indicate temperatures for different assumptions on the dust scattering. From the bottom up, they correspond to $\epsilon = 1$ (no scattering), 0.9, 0.8, 0.7, 0.6, 0.5, 0.4, and 0.3. For completeness, we also show with teal crosses the dust temperature assumed by Dullemond et al. (2018). Right panel: the black curve and the dashed blue region indicate the CO peak brightness temperature and the standard deviation from the mean, respectively. Green, red, and brown solid curves indicate the temperature of the CO-emitting layer predicted by Flaherty et al. (2015), Rosenfeld et al. (2013), and Isella et al. (2016). Dotted and dashed curves indicate the temperature in the disk midplane and disk surface predicted by the same models, as derived by modeling the dust continuum and $^{12}$CO J = 2–1 emission.</figcaption></figure>

et al. (2013) and Isella et al. (2016), and would imply optical depths of about 1.3 and 0.7, respectively. Beyond 30 au, $T_b^{pl}$ drops well below the model predictions, suggesting that the emission becomes more optically thin. In particular, the model temperatures give $\tau_d^{abs} < 0.02$ and $\tau_d^{abs} \sim 0.04-0.05$ at the center of the dark gaps D45 and D86, respectively, and $\tau_d^{abs} \sim 0.5-0.7$ and $\tau_d^{abs} \sim 0.3-0.5$ at the bright rings B67 and B100.

These measurements rely on the assumption that scattering is negligible, but, as discussed in Section 4.1, this is likely not the case if dust grains have sizes comparable or larger than the wavelength of the observations. The inclusion of scattering significantly complicates the interpretation of the observed continuum emission and a simple relation between $T_b^{pl}$, $T_d$, and the absorption and scattering optical depths, $\tau_d^{abs}$ and $\tau_d^{sca}$ does not exist. Fortunately, however, it is possible to derive an approximated analytic solution to constrain $T_d$, $\tau_d^{abs}$ or $\tau_d^{sca}$ if the other two quantities are known (see Section 5 of Birnstiel et al. 2018).

In the case of HD 163296, the direct measurement of the extinction opacity at the bright rings B67 and B100 discussed in Section 4.1 provide a unique tool to constrain the scattering properties of the circumstellar dust. Using the formalism discussed in Birnstiel et al. (2018), we calculate the physical dust temperature required to reproduce the dust continuum emission measured at B67 (0.77 mJy beam⁻¹) and B100 (0.45 mJy beam⁻¹) for values of the scattering parameter $\epsilon$ varying from 0.1 ($\tau_d^{sca} = 9\tau_d^{abs}$) to 1 (no scattering). The comparison between these temperatures, which are indicated with horizontal segments in the left panel of Figure 9, and those predicted by Flaherty et al. (2015), Isella et al. (2016), and Rosenfeld et al. (2013) constrain $\epsilon$, and therefore the

dust-scattering properties. For B67, we find that $\epsilon \gtrsim 0.6$ leads to dust temperature that is consistent with the models. We note that the midplane temperature of Flaherty et al. (2015) is slightly below the minimum temperature consistent with the observations, which corresponds to $\epsilon = 1$. Interestingly, for B100, the consistency between models and observations requires $\epsilon$ between 0.4 and 0.6, while $\epsilon = 1$ would imply a physical dust temperature of only 15 K. The difference in scattering properties between B67 and B100 is intriguing and, taken at face value, might suggest a variation in the dust properties between B67 and B100. However, the difference might also result from the fact that, due to beam smearing, the measured value of $\tau_d^{ext}(B67)$ might be a lower estimate of the extinction opacity, and consequently, of the scattering opacity. Future ALMA observations expressly designed to image the molecular gas emission from the HD 163296 at both high angular resolution and sensitivity are required to place more stringent constraints on the dust extinction optical depth.

We conclude this section by comparing the peak brightness temperature of the CO emission measured in Section 4.2 to the temperature of the CO-emitting layer predicted by Flaherty et al. (2015), Isella et al. (2016), and Rosenfeld et al. (2013) models (see the right panel of Figure 9). In Section 4.1 we found that the CO-emitting layer corresponds to the surface defined by $z(r) = 30$ au ($r/100$ au)⁰.⁵. Within about 100 au, the measured gas temperature is below the model predictions. This is likely due to the effect of beam dilution as discussed in Weaver et al. (2018). Between 100 and 400 au, the temperature of the CO-emitting layer predicted by Rosenfeld et al. (2013) and Isella et al. (2016) is in good agreement with the observations, while the temperature predicted by Flaherty et al. (2015) is about 30% lower. Also, while the temperature of
<figure><img src="image_11.png" /><figcaption>**Figure 10.** Left panel: image of the 1.25 mm dust continuum emission as in Figure 1. Right panel: synthetic image of the 1.25 dust-continuum emission of a disk perturbed by a 0.15 $M_J$ planet orbiting at about 54 au from the central star. The position of the planet is indicated by the yellow cross. The model is part of the suite of models presented in Zhang et al. (2018). The color scale of the two panels is the same and corresponds to the intensity normalized to the peak value. No attempts have been made to match the model to the observations, beside rotating and inclining the model image to match the orientation of the HD 163296 disk.</figcaption></figure>

the CO-emitting layer predicted by Rosenfeld et al. (2013) and Isella et al. (2016) is similar to that of the stellar irradiated disk surface (dashed line), as expected in the case of a very optically thick line, the Flaherty et al. (2015) temperature is closer to that achieved in the disk midplane (dotted lines).

## 5.2. On the Origin of the Dust Rings

In Section 3, we found that the radial profile of the continuum emission from the HD 163296 disk can be represented by a sum of azimuthally symmetric Gaussian rings. However, significant asymmetric structures appear once the symmetric component of the emission is removed. In this section, we discuss these results in the framework of the planet–disk interaction mechanism, while we refer to Huang et al. (2018) for a discussion of other processes capable of generating dust rings. We also refer the reader to Dullemond et al. (2018) for a discussion of the observed rings properties in the context of dust-trap models.

Massive planets are expected to open gaps in the surface density of protoplanetary disks along their orbits through the strong mutual gravitational interaction with the circumstellar material (see, e.g., Bryden et al. 1999). While the planet minimum mass for this to happen varies as a function of the disk physical properties, most of the simulations agree that for typical disk conditions planets more massive than 20–40 Earth masses should produce gaps observable at the spatial scales probed by ALMA (e.g., Dong et al. 2012; Isella & Turner 2018). In the limit of nearly inviscid disks, simulations have also shown that even planets with masses as low as few Earth masses might generate visible rings (Dong et al. 2018). Given the ubiquity of planets in the universe, it seems therefore reasonable to associate the ring structures observed by ALMA with the presence of planetary systems in the act of forming.

In the case of HD 163296, the planet–disk interaction hypothesis is supported by the detection of gas depletion and

deviation from Keplerian rotation within the dark gaps D87 and D143 (Isella et al. 2016; Teague et al. 2018), which imply that the emission gaps correspond to gaps in the gas density. The comparison with hydrodynamic simulations indicates that D45, D87, and D143 might results from the dust clearing operated by planets with masses between 0.5 and 2 $M_J$ (Liu et al. 2018; Teague et al. 2018). However, the ring morphology of the HD 163296 disk might also be consistent with gravitational perturbations by a single planet with a mass as low as 65 $M_E$ orbiting at about 100 au from the star if the disk viscosity is very low ($\alpha < 10^{-4}$, Dong et al. 2018). Overall, these results are consistent with the planet masses discussed in Zhang et al. (2018), which have been estimated from the azimuthally averaged profile of the continuum emission shown in Figure 2. In addition, Zhang et al. (2018) find that the newly discovered gap B10 might be generated by a planet with a mass between 0.2 and 1.5 $M_J$.

Here we want to point out that the detection of asymmetries in the dust emission discussed in Section 3.3 corroborates the hypothesis that the HD 163296 disk is shaped by the interaction with yet-unseen planets. As shown in Figures 6 and 7 of Zhang et al. (2018), features similar to the crescent observed inside the D45 gap (inset b of Figure 1) are naturally produced by the planet–disk interaction whenever the disk viscosity is low ($\alpha < 10^{-3}$) and the planet mass is above about 0.1 $M_J$. The crescents observed in numerical simulations have two different origins: they either trace dust particles trapped in anticyclonic vortices that form at the edge of the gap opened by the planet (Li et al. 2001), or probe material trapped at the Lagrangian points along the planet orbit. Without any pretense of accuracy, we compare in Figure 10 the map of the 1.25 mm dust continuum emission recorded toward HD 163296 with a synthetic image of a disk model characterized by the presence of a 0.15 $M_J$ planet orbiting at 54 au from the central star. The model is part of the suite of models presented in Zhang et al. (2018) and was generated assuming a viscosity parameter
$\alpha = 10^{-4}$ and a disk pressure scale height $h = 0.05r$. This comparison is purely qualitative, but nevertheless shows the resemblance between the corotation features predicted by models and the circular arc observed in HD 163296.

## 6. Conclusion

The unprecedented resolution of the DSHARP/ALMA data of the HD 163296 circumstellar disk enabled us to dive into a detailed characterization of the morphology of the dust continuum emission and of some of the physical properties of the circumstellar material. The comparison between images of the 1.25 mm continuum and $^{12}$CO line emission and parametric disk models led us to the following conclusions.

1. The new ALMA observations confirm the presence of two bright rings in the millimeter-wave dust continuum emission centered at 67 and 101 au, and reveal an additional ring with a radius of 15 au. The radial profile of the dust emission across the rings is well described by a Gaussian profiles with widths of 8.7, 6.6, and 5.8 au, respectively.

2. The dust continuum emission is characterized by several asymmetric structures. The most prominent consist in crescents centered at about 4 and 55 au. We also find asymmetric structures along the ring B67 characterized by amplitude variations of ±15% relative to the ring mean intensity. We argue that this asymmetries trace local variations in the dust densities and support the hypothesis that the HD 163296 disk is shaped by the gravitational interaction with yet-unseen planets.

3. The observations of the $^{12}$CO line emission probe the temperature and geometrical structure of the emitting gas. We find that the peak of the line arises from a disk layer at a vertical distance from the disk midplane given by $z_{co} = 30$ au($r/100$ au)$^{0.5}$. Between 30 au and 500 au, the temperature of the CO-emitting layer scales linearly from the orbital radius and follows the relation $T_{co}(K) \sim 87 - 0.14(r/au)$.

4. The CO maps show the line emission coming from both the front and rear side of the disk. The latter is attenuated as it passes through the dust rings, providing a tool to directly measure the extinction optical depth $\tau_d^{ext}$ of the dust. We measure $\tau_d^{ext} \sim 0.7$ at the position of the B67 and B100 rings, and $\tau_d^{ext} < 0.3$ for B168. We point out that if the absorption optical depth of the dust rings $\tau_d^{abs}$ can be measured through other means, a direct measurement of $\tau_d^{ext}$ allows to constrain the dust-scattering properties.

5. Finally, we compare the CO temperature inferred from our ALMA observations to theoretical predictions based on previous forward modeling of millimeter wave observations of HD 163296. Overall, we find that models predict temperatures within 30% of the measured values. By adopting the predicted temperatures for the disk midplane, we find that $\tau_d^{abs}$ is between about 0.5 and 0.6 at B67, and between 0.3 and 0.5 at B100. This implies $\tau_d^{sca}$ between 0.1 and 0.2 at B67, and between 0.2 and 0.4 at B100.

We thank the anonymous referee for the comments that helped improving the manuscript. A.I. acknowledges support

from the National Aeronautics and Space Administration under grant No. NNX15AB06G issued through the Origins of Solar Systems program, and from the National Science Foundation under grant No. AST-1715719. S.A. and J.H. acknowledge support from the National Aeronautics and Space Administration under grant No. 17-XRP17_2-0012 issued through the Exoplanets Research Program. J.H. acknowledges support from the National Science Foundation Graduate Research Fellowship under grant No. DGE-1144152. L.P. acknowledges support from CONICYT project Basal AFB-170002 and from FCFM/U. de Chile Fondo de Instalación Académica. C.P.D. acknowledges support by the German Science Foundation (DFG) Research Unit FOR 2634, grants DU 414/22-1 and DU 414/23-1. V.V.G. and J.C. acknowledge support from the National Aeronautics and Space Administration under grant No. 15XRP15_20140 issued through the Exoplanets Research Program. T.B. acknowledges funding from the European Research Council (ERC) under the European Unions Horizon 2020 research and innovation programme under grant agreement No. 714769. M.B. acknowledges funding from ANR of France under contract number ANR-16-CE31-0013 (Planet Forming disks). Z.Z. and S.Z. acknowledge support from the National Aeronautics and Space Administration through the Astrophysics Theory Program with grant No. NNX17AK40G and the Sloan Research Fellowship. L.R. acknowledges support from the ngVLA Community Studies program, coordinated by the National Radio Astronomy Observatory, which is a facility of the National Science Foundation operated under cooperative agreement by Associated Universities, Inc. M.B. acknowledges funding from ANR of France under contract number ANR-16-CE31-0013 (Planet Forming disks).

# Appendix
## Calculation of Isovelocity Lines

The isovelocity lines shown in Figure 5 are calculated as follows. First, we assume that CO molecules rotate around the star at the Keplerian velocity $v_k(R) = \sqrt{GM_*/R}$, where R is the spherical radius. second, we assume that the line emission arises from a geometrically thin layer characterized by the distance from the disk midplane $z_{co}(r) = z_0(r/r_0)^q$, where r is the cylindrical radius, $r = \sqrt{R^2 + z_{co}^2}$. Finally, we assume that the disk midplane is inclined with respect to the line of sight by an angle i and is rotated in the sky by the PA. The component of the Keplerian velocity along the line of sight is $v_l(R, \theta) = v_k(R)\sin i \sin \theta$, where $\theta$ is the azimuthal angle. The isovelocity contours corresponding to the velocity v in a spherical reference frame centered at the center of the disk is calculated by finding the radii that satisfy the relation $v - v_l(R, \theta) = 0$, for $\theta$ between 0 and $2\pi$. Finally, the conversion between the reference frame centered on the disk and the image plane is

$$ \delta R.A. = x' \sin \text{PA} - y' \cos \text{PA} \quad (6) $$

$$ \delta \text{Decl.} = x' \cos \text{PA} + y' \sin \text{PA} \quad (7) $$

where

$$ x' = r \cos \theta \cos i + z_{co}(r) \sin i \quad (8) $$
<figure><img src="image_12.png" /><figcaption>**Figure 11.** Left panel: map of the 1.25 mm continuum emission as in Figure 1. Right panel: image of the residual visibilities resulting from the Gaussian model fitting discussed in Section 3.2. Both panels use the same color scale to facilitate the comparison.</figcaption></figure>

$$y' = r \sin \theta. \qquad (9)$$

A comparison between the image of the residual and the continuum map is shown in Figure 11.

## ORCID iDs

Andrea Isella https://orcid.org/0000-0001-8061-2207

Jane Huang https://orcid.org/0000-0001-6947-6072

Sean M. Andrews https://orcid.org/0000-0003-2253-2270

Cornelis P. Dullemond https://orcid.org/0000-0002-7078-5910

Tilman Birnstiel https://orcid.org/0000-0002-1899-8783

Shangjia Zhang https://orcid.org/0000-0002-8537-9114

Zhaohuan Zhu https://orcid.org/0000-0003-3616-6822

Myriam Benisty https://orcid.org/0000-0002-7695-7605

John M. Carpenter https://orcid.org/0000-0003-2251-0602

David J. Wilner https://orcid.org/0000-0003-1526-7587

## References

ALMA Partnership, Brogan, C. L., Pérez, L. M., et al. 2015, ApJL, 808, L3

Andrews, S. M., Huang, J., Pérez, L. M., et al. 2018, ApJL, in press

Andrews, S. M., Rosenfeld, K. A., Wilner, D. J., & Bremer, M. 2011a, ApJL, 742, L5

Andrews, S. M., Wilner, D. J., Espaillat, C., et al. 2011b, ApJ, 732, 42

Andrews, S. M., Wilner, D. J., Zhu, Z., et al. 2016, ApJL, 820, L40

Bai, X.-N., & Stone, J. M. 2014, ApJ, 796, 31

Banzatti, A., Pinilla, P., Ricci, L., et al. 2015, ApJL, 815, L15

Benisty, M., Natta, A., Isella, A., et al. 2010, A&A, 511, A74

Benisty, M., Stolker, T., Pohl, A., et al. 2017, A&A, 597, A42

Béthune, W., Lesur, G., & Ferreira, J. 2017, A&A, 600, A75

Birnstiel, T., Dullemond, C. P., Zhu, Z., et al. 2018, ApJL, in press

Boehler, Y., Ricci, L., Weaver, E., et al. 2018, ApJ, 853, 162

Bohren, C. F., & Huffman, D. R. 1998, Absorption and Scattering of Light by Small Particles (New York: Wiley-Interscience)

Bryden, G., Chen, X., Lin, D. N. C., Nelson, R. P., & Papaloizou, J. C. B. 1999, ApJ, 514, 344

Casassus, S., van der Plas, G., M. S. P., et al. 2013, Natur, 493, 191

de Gregorio-Monsalvo, I., Ménard, F., Dent, W., et al. 2013, A&A, 557, A133

Dong, R., Liu, S.-y., Eisner, J., et al. 2018, ApJ, 860, 124

Dong, R., Rafikov, R., Zhu, Z., et al. 2012, ApJ, 750, 161

Dong, R., van der Marel, N., Hashimoto, J., et al. 2017, ApJ, 836, 201

Dullemond, C. P., Birnstiel, T., Huang, J., et al. 2018, ApJL, in press

Fedele, D., Carney, M., Hogerheijde, M. R., et al. 2017, A&A, 600, A72

Fedele, D., Tazzari, M., Booth, R., et al. 2018, A&A, 610, A24

Flaherty, K. M., Hughes, A. M., Rosenfeld, K. A., et al. 2015, ApJ, 813, 99

Flock, M., Ruge, J. P., Dzyurkevich, N., et al. 2015, A&A, 574, A68

Foreman-Mackey, D., Hogg, D. W., Lang, D., & Goodman, J. 2013, PASP, 125, 306

Gaia Collaboration, Brown, A. G. A., Vallenari, A., et al. 2018, A&A, 616, 1

Garufi, A., Quanz, S. P., Schmid, H. M., et al. 2014, A&A, 568, A40

Goodman, J., & Weare, J. 2010, Comm. App. Math. Comp. Sci., 5, 65

Grady, C. A., Devine, D., Woodgate, B., et al. 2000, ApJ, 544, 895

Guidi, G., Ruane, G., Williams, J. P., et al. 2018, MNRAS, 479, 1505

Huang, J., Andrews, S., Dullemond, C. P., et al. 2018, ApJL, in press

Isella, A., Chandler, C. J., Carpenter, J. M., Pérez, L. M., & Ricci, L. 2014, ApJ, 788, 129

Isella, A., Guidi, G., Testi, L., et al. 2016, PhRvL, 117, 251101

Isella, A., Pérez, L. M., Carpenter, J. M., et al. 2013, ApJ, 775, 30

Isella, A., Testi, L., Natta, A., et al. 2007, A&A, 469, 213

Isella, A., & Turner, N. J. 2018, ApJ, 860, 27

Jin, S., Li, S., Isella, A., Li, H., & Ji, J. 2016, ApJ, 818, 76

Keppler, M., Benisty, M., Müller, A., et al. 2018, A&A, 617, 44

Klaassen, P. D., Juhasz, A., Mathews, G. S., et al. 2013, A&A, 555, A73

Li, H., Colgate, S. A., Wendroff, B., & Liska, R. 2001, ApJ, 551, 874

Liu, S.-F., Jin, S., Li, S., Isella, A., & Li, H. 2018, ApJ, 857, 87

Marino, S., Perez, S., & Casassus, S. 2015, ApJL, 798, L44

Mendigutía, I., Oudmaijer, R. D., Schneider, P. C., et al. 2018, A&A, 618, L9

Miranda, R., Li, H., Li, S., & Jin, S. 2017, ApJ, 835, 118

Muro-Arena, G. A., Dominik, C., Waters, L. B. F. M., et al. 2018, A&A, 614, A24

Okuzumi, S., Momose, M., Sirono, S.-i., Kobayashi, H., & Tanaka, H. 2016, ApJ, 821, 82

Pinte, C., Ménard, F., Duchêne, G., et al. 2018, A&A, 609, A47

Rosenfeld, K. A., Andrews, S. M., Hughes, A. M., Wilner, D. J., & Qi, C. 2013, ApJ, 774, 16

Sallum, S., Follette, K. B., Eisner, J. A., et al. 2015, Natur, 527, 342

Suriano, S. S., Li, Z.-Y., Krasnopolsky, R., & Shang, H. 2018, MNRAS, 477, 1239

Tang, Y.-W., Guilloteau, S., Dutrey, A., et al. 2017, ApJ, 840, 32

Tazzari, M., Testi, L., Natta, A., et al. 2018, MNRAS, 476, 4527

Teague, R., Bae, J., Bergin, E. A., Birnstiel, T., & Foreman-Mackey, D. 2018, ApJL, 860, L12
Thalmann, C., Janson, M., Garufi, A., et al. 2016, ApJL, 828, L17

Tilling, I., Woitke, P., Meeus, G., et al. 2012, A&A, 538, A20

van der Marel, N., Cazzoletti, P., Pinilla, P., & Garufi, A. 2016, ApJ, 832, 178

van der Marel, N., van Dishoeck, E., Bruderer, S., Perez, L. M., & Isella, A. 2015, A&A, 579, 106

van der Marel, N., van Dishoeck, E. F., Bruderer, S., et al. 2013, Sci, 340, 1199

van der Marel, N., Williams, J. P., Ansdell, M., et al. 2018a, ApJ, 854, 177

van der Marel, N., Williams, J. P., & Bruderer, S. 2018b, ApJL, 867, L14

Weaver, E., Isella, A., & Boehler, Y. 2018, ApJ, 853, 113

Wisniewski, J. P., Clampin, M., Grady, C. A., et al. 2008, ApJ, 682, 548

Zhang, K., Bergin, E. A., Blake, G. A., et al. 2016, ApJL, 818, L16

Zhang, K., Blake, G. A., & Bergin, E. A. 2015, ApJL, 806, L7

Zhang, K., Isella, A., Carpenter, J. M., & Blake, G. A. 2014, ApJ, 791, 42

Zhang, S., Zhu, Z., Huang, J., et al. 2018, ApJL, in press

Zhu, Z., Stone, J. M., Rafikov, R. R., & Bai, X.-n. 2014, ApJ, 785, 122