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

Variable Outer Disk Shadowing around the Dipper Star RXJ1604.3–2130*

P. Pinilla<sup>1,9</sup> <sup>①</sup>, M. Benisty<sup>2,3</sup> <sup>①</sup>, J. de Boer<sup>4</sup> <sup>①</sup>, C. F. Manara<sup>5</sup> <sup>①</sup>, J. Bouvier<sup>3</sup>, C. Dominik<sup>6</sup> <sup>①</sup>, C. Ginski<sup>4,6</sup>, R. A. Loomis<sup>7</sup> <sup>①</sup>, and A. Sicilia Aguilar<sup>8</sup>

<sup>1</sup> Department of Astronomy/Steward Observatory, The University of Arizona, 933 North Cherry Avenue, Tucson, AZ 85721, USA; pinilla@email.arizona.edu

<sup>2</sup> Unidad Mixta Internacional Franco-Chilena de Astronomía (CNRS, UMI 3386), Departamento de Astronomía, Universidad de Chile, Camino El Observatorio 1515, Las Condes, Santiago, Chile

<sup>3</sup> Univ. Grenoble Alpes, CNRS, IPAG, F-38000 Grenoble, France

<sup>4</sup> Leiden Observatory, Leiden University, P.O. Box 9513, 2300 RA Leiden, The Netherlands

<sup>5</sup> European Southern Observatory, Karl-Schwarzschild-Str. 2, D-85748 Garching, Germany

<sup>6</sup> Anton Pannekoek Institute for Astronomy, University of Amsterdam, Science Park 904,1098XH Amsterdam, The Netherlands

<sup>7</sup> Harvard-Smithsonian Center for Astrophysics, Cambridge, MA 02138, USA

<sup>8</sup> SUPA, School of Science and Engineering, University of Dundee, Nethergate, Dundee DD1 4HN, UK

Received 2018 July 2; revised 2018 October 10; accepted 2018 October 10; published 2018 November 26

Abstract

Low brightness dips have recently been observed in images of protoplanetary disks, and they are believed to be shadows by the inner disk. We present VLT/SPHERE polarimetric differential imaging of the transition disk around the dipper star RX J1604.3–2130. We gathered 11 epochs that cover a large temporal baseline, to search for variability over timescales of years, months, weeks, and days. Our observations unambiguously reveal two dips along an almost face-on narrow ring (with a width of ~20 au), and the location of the peak of this ring is at ~65 au. The ring lies inside the ring-like structure observed with ALMA, which peaks at ~83 au. This segregation can result from particle trapping in pressure bumps, potentially due to planet(s). We find that the dips are variable, both in morphology and in position. The eastern dip, at a position angle (PA) of ~83°.7 ± 13°.7, has an amplitude that varies between 40% to 90%, and its angular width varies from 10° to 34°. The western dip, at a PA of ~265°.90 ± 13°.0, is more variable, with amplitude and width variations of 31% to 95% and 12° to 53°, respectively. The separation between the dips is 178°.3 ± 14°.5, corresponding to a large misalignment between the inner and outer disks, supporting the classification of J1604 as an aperiodic dipper. The variability indicates that the innermost regions are highly dynamic, possibly due to a massive companion or to a complex magnetic field topology.

**Key words:** accretion, accretion disks – circumstellar matter – planets and satellites: formation – protoplanetary disks – stars: individual ([PZ99] J160421.7–213028)

**Supporting material:** animation

1. Introduction

In the last decade, extreme adaptive optics and coronagraphic observations at optical and near-infrared wavelengths in combination with high angular resolution and sensitivity observations at millimeter wavelengths have opened a new window for our understanding of planet formation. Powerful telescopes such as the Atacama Large (sub-)Millimeter Array (ALMA) and the Spectro-Polarimetric High-contrast Exoplanet REsearch (SPHERE) at the Very Large Telescope (VLT) have provided unprecedented angular resolution down to few astronomical units in protoplanetary disks. This has allowed us to resolve regions in disks where planets can leave an imprint on their nascent environment. These recent observations of protoplanetary disks unveiled an incredible variety of structures, such as cavities, rings, gaps, spiral arms, arcs, and low brightness dips, which are likely related to the diversity of physical processes that rule planet formation and disk evolution.

One of the most exciting sets of protoplanetary disks to study the processes of planet formation are those hosting dust depleted cavities, the so-called transition disks. Most of the disks revealing asymmetric features in scattered light, such as spiral arms (e.g., Muto et al. 2012; Benisty et al. 2015; Stolker et al. 2016)

and shadows (e.g., Avenhaus et al. 2014; Benisty et al. 2017; Canovas et al. 2017) are in fact transition disks. At long wavelengths, these disks usually show, as a main feature, a large dust depleted cavity surrounded by a ring-like (symmetric or asymmetric) structure (e.g., Casassus et al. 2013; Pérez et al. 2014; Zhang et al. 2014; Pinilla et al. 2018; van der Marel et al. 2018). This suggests that the physical processes dominating the evolution of these disks may create different structures in small (micron-sized) and large (millimeter/centimeter-sized) particles, producing a diversity of morphologies at multiwavelength observations.

From these sets of observations, low surface brightness dips in the outer parts of disks (>10 au) observed in scattered light are of particular interest because they are believed to result from shadowing by a misaligned inner disk with respect to the outer disk, when the inner and outer disks have different inclinations and position angles (also called warped or broken disks). These observations allow us to directly connect the resolved outer regions with the unresolved inner disk. Scattered light observations of some transition disks suggest a small to intermediate misalignment between the inner and the outer disks (HD 135344B, LkCa 15, TW Hya, DoAr 44, HD 143006; Debes et al. 2016; Oh et al. 2016; Stolker et al. 2016; Benisty et al. 2018; Casassus et al. 2018), while observations of other transition disks suggest a large misalignment (~70°, HD 142527 and HD100453; Marino et al. 2015; Benisty et al. 2017). In addition, some of

<footnote>* Based on observations performed with SPHERE/VLT under program IDs 099.C-0341(A), 097.C-0902(A) and 095.C-0693(A).</footnote>

<footnote n="9">⁹ Hubble Fellow.</footnote>
these shadows, as in the case of HD 135344B, are variable within timescales shorter than a week, suggesting a very dynamic and asymmetric inner region. Interestingly, warps have also been observed at the very late stages of planet formation, such as in the debris disks AU Mic and β Pic (e.g., Golimowski et al. 2006; Dawson et al. 2011; Boccaletti et al. 2015). Recently, Casassus et al. (2018) identified shadows in another transition disk, DoAr 44, that are not only detected in scattered light, but also in millimeter-observations, suggesting that the shadows are effective in cooling millimeter-dust grains. Molecular line ALMA observations, that show clear deviations from Keplerian motion, also support the idea that some transition disks are warped (e.g., Rosenfeld et al. 2012; Casassus et al. 2015; Loomis et al. 2017; Walsh et al. 2017).

Warps are mainly identified in two regimes depending on how they propagate (Papaloizou & Pringle 1983). On one hand, warps can propagate following the diffusion equation when disks are thin, in particular when $h/r \lesssim \alpha$, where $\alpha$ is the dimensionless viscosity parameter (Shakura & Sunyaev 1973), $h$ is the disk scale height, and $r$ is the distance from the central star. On the other hand, when disks are thick and $h/r \gtrsim \alpha$, warps are expected to propagate as bending waves (Papaloizou & Lin 1995). In both cases, if the sound crossing time is longer than the induced precession time, the disk can warp significantly (e.g., Papaloizou & Terquem 1995; Nixon & King 2012; Facchini et al. 2013; Nealon et al. 2015). For the origin of these warps in protoplanetary disks, different possibilities have been proposed, such as a misaligned stellar magnetic field with respect to the rotation axis of the star or/ and the disk (e.g., Bouvier et al. 2007; Romanova et al. 2013); a misaligned circumbinary disk with respect to the binary orbital plane (e.g., Foucart & Lai 2013); and the interaction with a massive planet or a binary companion on an inclined, and perhaps eccentric, orbit (e.g., Xiang-Gruess & Papaloizou 2013; Lubow & Martin 2016; Martin et al. 2016; Owen & Lai 2017; Facchini et al. 2018). The last scenario is likely the case for the transition disk around HD 142527 where a $0.1 M_{\odot}$ companion in an eccentric orbit with an inclination of $\sim 125^\circ$ within a close to face-on disk (Lacour et al. 2016) may explain most of the observed properties, including a large cavity, shadows, spiral arms, deviation from Keplerian velocity in the CO lines, and the horseshoe-like structure at dust millimeter continuum emission (Price et al. 2018).

In this paper, we present new SPHERE polarimetric differential imaging of the transition disk around the star RX J1604.3-2130$^{10}$ (2MASS J16042165-2130284, hereafter J1604) at different wavelengths. This disk is a member of the Upper Scorpius association and it is located at a distance of $150.2 \pm 1.4$ pc$^{11}$ (Gaia Collaboration et al. 2018). The disk hosts one of the largest dust cavities observed at millimeter wavelengths (Mathews et al. 2012; Zhang et al. 2014; Dong et al. 2017). Observations from ALMA indicate that J1604 is as bright as the most luminous disks in younger regions (e.g., in Taurus and Lupus) in the millimeter, with its dust concentrated in a ring-like structure, and a dust mass of $\sim 40-50 M_{\oplus}$ (Barenfeld et al. 2016; Pinilla et al. 2018). J1604 was previously observed with scattered light with HiCIAO at $1.6 \mu$m (Mayama et al. 2012) and ZIMPOL/SPHERE at $0.626 \mu$m (Pinilla et al. 2015). Both of these observations reveal a ring-like structure in scattered light located at $\sim 0.4''/4$ from the star with a single dip

along the ring, located in the east. Interestingly, in these two observations obtained 3 years apart, the dip was detected at different position angles (PA): $\sim 85^\circ$ (HiCIAO) and $\sim 46^\circ$ (ZIMPOL). If the dip originated from the same shadowing structure, its very fast rotation is inconsistent with the local Keplerian velocity at the ring position (Pinilla et al. 2015). Mayama et al. (2012) also reported a marginal detection of a second dip in the west at a PA of 255$^\circ$.

J1604 was identified as a dipper (Ansdell et al. 2016), which are young stellar objects (YSOs) that exhibit light curves punctuated by recurrent (periodic or aperiodic) dimming events on timescales of a few days, for example, from the Convection, Rotation and Planetary Transits satellite (CoRoT), *Spitzer*, and *Kepler* 2. Dippers are fairly common among YSOs ($\sim 30\%-40\%$, Alencar et al. 2010; Cody et al. 2014; Bodman et al. 2017) and are thought to be systems seen at high inclination, such that the dimming events are due to patches of dusty material that repeatedly occult the star as they cross our line of sight (e.g., Scaringi et al. 2016; Schneider et al. 2018). The prototype dipper, AA Tau, exhibits periodic eclipses every 8.2 days (Bouvier et al. 1999) and a model where the inner edge of the accretion disk is warped by its interaction with the inclined stellar magnetosphere successfully explains the light curve and the spectral variability (Bouvier et al. 2007). Aperiodic dippers have also been interpreted as resulting from clumps of dusty material passing our line of sight to the star, by related or possibly different mechanisms, such as vortices and forming planetesimals (Ansdell et al. 2016). J1604 is very interesting in this context because it shows aperiodic dimming events and it is one of the deepest flux dips among the known K2 dippers (Ansdell et al. 2016). However, its outer disk is close to face-on (Mathews et al. 2012; Zhang et al. 2014), and hence its dipper nature (usually in highly inclined disks) suggests a possible strong misalignment between the inner and outer disk. In addition, it hosts a large ($\sim 83$ au) and highly depleted dust and gas cavity (Dong et al. 2017). In addition, J1604 also evidences variable near- and mid-infrared excess (Dahm & Carpenter 2009).

The paper is organized as follows. In Section 2, we describe the observations and the data reduction from our SPHERE observations at different epochs. In Section 3, we present the data analysis, including the inspection of the radial profile for each observation, the characterization of the dips, and the comparison with recent ALMA observations. In Section 4, we discuss these results in the context of different origins for the shadows and the potential connection of their variability with the dipper nature of J1604. In addition, we also discuss in this section the potential origin of the observed cavity and the evidence of particle trapping from multiwavelength observations. Finally, in Section 5, we provide the conclusions of this work.

## 2. Observations

We obtained multiple epochs observations of J1604 at the Very Large Telescope located at Cerro Paranal, Chile, using the SPHERE instrument (Beuzit et al. 2008), a high-contrast imager with an extreme adaptive optics system (Sauvage et al. 2014). In this paper, we present new polarimetric observations covering nine epochs, obtained between 2016 June and 2017 September, in the near-infrared (J- and H-bands) with the IRDIS instrument (Dohlen et al. 2008). Our new data set is complemented by previously published

<footnote>$^{10}$ SIMBAD name.</footnote>

<footnote>$^{11}$ J1604 is close by and the parallax has a small relative uncertainty, which justifies why the uncertainty is taken as symmetric.</footnote>
visible (R'-band) polarimetric data obtained in 2015 June with ZIMPOL (Pinilla et al. 2015) and near-infrared (H-band) polarimetric observations obtained with Subaru/HiCIAO in 2012 April (Mayama et al. 2012). This allows us to cover a large temporal baseline (2012–2017), and to search variability over timescales of years, months, and days. In all our IRDIS observations, we used a 185 mas diameter coronagraph (N_ALC_YJH_S) to enhance the signal-to-noise ratio on the outer disk regions. The plate scale is 12.26 mas per pixel. We observed J1604 in excellent to good seeing conditions (between 0".5 and 1").

The data reduction procedure is similar to the one reported in de Boer et al. (2016) and is only very briefly described here. In polarimetric differential imaging, the stellar light is split into two orthogonal polarization states, and a half-wave plate (HWP) is set to four positions shifted by 22°.5 to construct a set of linear Stokes images. The data is then reduced following the double difference method, from which one can derive the Stokes parameters Q and U. If we assume single scattering events on the protoplanetary disk surface, the scattered light is linearly polarized in the azimuthal direction; therefore, for convenience, we describe the polarization vector field in polar coordinates with the $Q_φ$, $U_φ$ Stokes images (Schmid et al. 2006). In this framework, the $Q_φ$ image should contain all disk signals, while the $U_φ$ image remains free of it.

Figure 1 shows the Stokes parameter $Q_φ$ of all of our SPHERE observations, including our previously published ZIMPOL image, and the H-band polarized intensity observations with HiCIAO reported by Mayama et al. (2012). In all our new observations, we clearly detect two dark regions, hereafter referred to as dips or shadows, along the ring. These dips appear strongly variable with time, both in location and morphology.

### 3. Data Analysis

#### 3.1. Radial Profile

Figure 2 shows the normalized radial profile of the surface brightness obtained after an azimuthal average. The uncertainty of the data corresponds to the standard deviation in each radial bin of the $Q_φ$ images, divided by the square root of the number of pixels in the bin. We find a narrow bell-shaped curve, with a tail extending at larger radii. We note that this extended tail is not an effect of the instrumental point-spread function (PSF). As a test, we convolved synthetic ring models, of various widths, with the PSF of each epoch and verified that this procedure cannot reproduce an extended tail beyond the peak as observed. Similarly, we also find that the PSF does not affect the dips' morphology, which is discussed in Section 3.2. The details of the PSF shape for each epoch are summarized in Table 3 in the Appendix.

To quantify the location of the peak and the width of the ring, we fitted the radial profile of the normalized surface brightness of each epoch using a Markov chain Monte Carlo (MCMC) method. Considering the typical shape of the radial

profile, we used a Lorentzian prescription, following:

$$A_L \times \left( \frac{\gamma}{(r - r_0)^2 + \gamma^2} \right). \qquad (1)$$

The model has three free parameters ($A_L, r_0, \gamma$), the amplitude, the location of the peak, and the half width at half maximum (HWM), respectively. To perform the fit, we use emcee (Foreman-Mackey et al. 2013), which allows us to efficiently sample the parameter space to derive the maximum likelihood result for each model. The parameter space explored by the Markov chain for the location of the peak and the width are: $r_0 \in [0".2, 0".6]$ and $\gamma \in [0".0, 0".2]$, with uniform prior probability distributions. The Markov chain sampled the parameter space for 1000 steps, with 100 walkers for each epoch. Table 1 summarizes the results from the fits for $r_0$, and $\gamma$, and includes the full width at half maximum (FWHM) too. The uncertainties from the MCMC fit are omitted since they are negligible compared to the mean value, with values of the order of $10^{-5}$ to $10^{-4}$ in all cases.

According to our model fit, the peak of the ring is 0".43 (~65 au) in all epochs. The FWHM varies from 0".12 (~18 au; 2017 August 22) to 0".18 (~27 au; 2017 September 4), in our J-band data. The epoch on 2017 September 4 has the highest uncertainty from the IRDIS observations (Figure 3), and if we neglect this epoch, the width is almost constant and only varies from 0".12 to 0".14. The FWHM of the ZIMPOL and HiCIAO data are 0".19 and 0".20, respectively; but these observations have low signal-to-noise ratios compared to the IRDIS observations.

To better quantify the extended tail of the ring, we performed a power-law fit to the surface brightness beyond the peak (i.e., taking the data from 0".43 outward), such that it is proportional to $r^\xi$. For this fit, we simply perform a nonlinear least squares analysis. The results for $\xi$ are also summarized in Table 1, and the best models are overplotted in Figure 3. The values of $\xi$ vary from -3.52 (ZIMPOL epoch on 2015 June 10) to -5.61 (IRDIS epoch at H-band on 2016 June 30). This power-law index can provide information about the vertical scale height of the disk, and its flaring index.

Assuming a single scattering approximation, that is, the scattering of the starlight happens where the optical depth is unity, the surface brightness of the disk is determined by the disk scale height (especially the scale height index $\beta$, such that $h(r) \propto r^\beta$, assuming hydrostatic equilibrium and vertically isothermal disks), and the radial profile of the surface density of the dust grains. For a disk with $h/r$ constant (i.e., a wedge-shaped disk), the surface brightness of the scattered light is expected to scale with radius as $\propto r^{-3}$, while for flared disks it is expected to scale as $\propto r^{\beta-3}$ (e.g., Whitney & Hartmann 1992; Dullemond et al. 2001; Dong et al. 2012). The range of the slope that we found for J1604 is within the values obtained for several other HAeBe disks (Fukagawa et al. 2010), which contradicts the flaring nature of several of these disks. However, as suggested by Fujiwara et al. (2006) and Grady et al. (2007), self-shadowing from, for example, the ring itself or a vertically thick inner wall, can create a steeper decrease of the surface brightness. Dong (2015) showed that indeed, the radial dependence of the surface brightness changes from a shallower $r^{-1.8}$ to a steeper $r^{-2.5}$ due to self-shadowing.

Another way to explain steep power-law values as the ones we find is to consider a depletion of small grains. Pohl et al. (2017) showed that the radial dependence of the surface
<figure><img src="image_2.png" /><figcaption>**Figure 1.** Scattered light observations of the transition disk around J1604 (they are not scaled by $r^2$). The left upper panel corresponds to *H*-band polarized intensity observations with HiCIAO reported by Mayama et al. (2012). The rest of the panels correspond to the Stokes parameter $Q_\phi$ obtained with VLT/SPHERE. The center upper panel corresponds to observations with ZIMPOL at *R'*-band. The right upper panel corresponds to observations with IRDIS at *H*-band, while the rest of the panels are IRDIS observations at *J*-band. In all the panels, the color scale is linear and in arbitrary units, and the dates are reported when the observations started.</figcaption></figure>

(An animation of this figure is available.)
<figure><img src="image_3.png" /><figcaption>**Figure 2.** Normalized radial profile of the surface brightness obtained after an azimuthal average. The uncertainty of the data correspond to the standard deviation in each radial bin, divided by the square root of the number of pixels in that bin. The shaded area corresponds to the coronagraph radius.</figcaption></figure>

<table><caption>**Table 1**
Results from Fitting the Radial Profile</caption><thead><tr><th rowspan="2">Epoch</th><th colspan="3">Lorentzian Fit</th><th rowspan="2">Power-law Fit<br>ξ</th></tr><tr><th>r<sub>0</sub><br>"</th><th>γ<br>"</th><th>FWHM<br>"</th></tr></thead><tbody><tr><td>2012 Apr 11</td><td>0.43</td><td>0.10</td><td>0.19</td><td>-3.92 ± 0.09</td></tr><tr><td>2015 Jun 10</td><td>0.43</td><td>0.10</td><td>0.20</td><td>-3.52 ± 0.03</td></tr><tr><td>2016 Jun 30</td><td>0.43</td><td>0.06</td><td>0.13</td><td>-5.61 ± 0.04</td></tr><tr><td>2017 Aug 13</td><td>0.43</td><td>0.07</td><td>0.14</td><td>-5.14 ± 0.05</td></tr><tr><td>2017 Aug 14</td><td>0.43</td><td>0.07</td><td>0.13</td><td>-5.25 ± 0.07</td></tr><tr><td>2017 Aug 17</td><td>0.43</td><td>0.06</td><td>0.12</td><td>-5.36 ± 0.06</td></tr><tr><td>2017 Aug 18</td><td>0.43</td><td>0.07</td><td>0.14</td><td>-4.98 ± 0.05</td></tr><tr><td>2017 Aug 22</td><td>0.43</td><td>0.06</td><td>0.12</td><td>-5.47 ± 0.06</td></tr><tr><td>2017 Sep 04</td><td>0.43</td><td>0.09</td><td>0.18</td><td>-4.21 ± 0.06</td></tr><tr><td>2017 Sep 06</td><td>0.43</td><td>0.07</td><td>0.13</td><td>-5.11 ± 0.07</td></tr><tr><td>2017 Sep 16</td><td>0.43</td><td>0.07</td><td>0.14</td><td>-4.63 ± 0.03</td></tr></tbody></table>

**Note.** The uncertainties of the Lorentzian fit are omitted since they are negligible compared to the mean value.

brightness due to a large reduction in the distribution of small grains can change from ~$r^{-2}$ to ~$r^{-5}$. These models considered that grains are trapped in pressure maxima due to a planet carving a gap in a disk. Inside such a pressure bump, large grains accumulate, and small grains are simultaneously generated due to fragmentation. At long times of evolution (~1–5 Myr), most of the grains are located inside these pressure maxima (Pinilla et al. 2012). In the case of J1604, both self-shadowing and a steep depletion of small grains beyond the ring can contribute to the values that we find for $\xi$. However, these values appear to change very

quickly, suggesting rapid variations of the flaring, of the self-shadowing, or/and of the distribution of small grains at large radii.

## 3.2. Characterization of the Dips

Figure 4 shows the radial mapping from 0"25–0"7 of all the images shown in Figure 1 (note that the images are not deprojected). The color scale is linear and the maximum value taken in each case is 80% of the peak. In addition, the azimuthal profile is calculated for each epoch after radially averaging between [0"35 – 0"50]. The uncertainty (shaded areas) are the standard deviation of the data in the radial and azimuthal bins (i.e., $\sqrt{\sigma_{radial}^2 + \sigma_{azimuthal}^2}$, where $\sigma_{radial}$ and $\sigma_{azimuthal}$ are the standard deviation of the data in the radial and azimuthal bins in the $Q_{\phi}$ images, respectively) divided by the square root of the total number of pixels within the ring. The data is normalized to the value at zero degrees in each case, which is our reference value to quantify the amplitude of these dips in the analysis that follows. We also test the following analysis normalizing to the peak of the data, which does not change the results.

### 3.2.1. Gaussian Model

To quantify the morphology of the dips and their variability, we performed a Gaussian fitting to each dip, similarly to the procedure followed with the radial profiles (Section 3.1). This means that we use *emcee*, and fit each dip independently with

$$ -A_{E,W} \times (\exp(-(PA - \min_{E,W})^2 / 2\sigma_{E,W}^2) + c_{E,W}). \quad (2) $$
<figure><img src="image_4.png" /><figcaption>**Figure 3.** Lorentzian MCMC best fit of the radial profile overlap with the uncertainty of the data profile for each epoch as shown in Figure 2. Notice that the radial range changes from Figure 2, to emphasize the ring shape and the fitting. In addition, we also overlay the power-law fit using nonlinear least squares analysis and considering the data points beyond the peak.</figcaption></figure>

The free parameters, $A_{E,W}$, $min_{E,W}$, $\sigma_{E,W}$, and $c_{E,W}$, correspond to the amplitude, location of the minimum, Gaussian standard deviation, and shift of the Gaussian peak, respectively, where E,W indicates east or west for each dip. To fit each dip, we consider the data from 28° to 150° and from 196° to 350° for the eastern and western dip, respectively. For the ZIMPOL epoch, we only perform a fit for the eastern dip, since the western dip is not clearly detected.

The parameter space explored by the Markov chain are: $A_{E,W} \in [0.1, 1.0]$, $min_E \in [28^\circ, 150^\circ]$, $min_W \in [196^\circ, 350^\circ]$, $\sigma_{E,W} \in [5^\circ, 80^\circ]$, and $c_{E,W} \in [0^\circ, 100^\circ]$, with uniform prior probability distributions. The Markov chain sample the parameter space for 2000 steps, with 100 walkers for each epoch. The results of the MCMC fit for the parameters that quantify the dip morphology (i.e., $A_{E,W}$, $min_{E,W}$, $\sigma_{E,W}$) are summarized in Table 2, and shown in Figure 5. We note that in most cases, our fitted $min_{E,W}$ values are close to the minimum value of the radial profile obtained from the images, except in cases in which the dip shape is complex, either with multiple minima (e.g., western dip in epoch 2017 August 18), or if there is a large difference in the levels of the radial profile on each side of the dips (e.g., eastern dip in epoch 2017 September 6). Indeed, for that reason, the fit is very poor for the western dip from the epochs of 2017 September 4 and 6. For these two epochs, we performed another fit by taking a narrower range for the PA from 236 to 300° (instead of 196°–350°), which provides a good match to the dip shape (see red dashed lines in Figure 5).

### 3.2.2. Variable Morphology

Figure 6 shows the variations of the parameters (location, amplitude, and width of each dip). For epochs 2017 September 4 and 6, we use the minimum location and the width of the dip from the fit that assumes a narrow PA range, but the amplitude is taken from the fit with the large PA range, since it averages the amplitude before and after the dip and thus gives a more accurate value for this parameter.

**Dip locations**—The location of the eastern dip varies from ~45° (ZIMPOL epoch) to ~98° (epoch of 2017 August 22), although the location of the minimum of the eastern dip from the ZIMPOL epoch is an outlier in our sample (Figure 6). The ZIMPOL data have low signal-to-noise ratios, but nonetheless we do not discard the existence of this dip at 45° in 2015. The mean value of the minimum location of the eastern dip is $min_E = 83^\circ.7 \pm 13^\circ.7$. If we neglect the ZIMPOL epoch for this calculation, the mean value of the minimum location of the eastern dip is $min_E = 87^\circ.6 \pm 6^\circ.3$, and the dip location varies within ~16° (from 82° to 98°). In the left panel of Figure 6, the dashed lines correspond to the mean values. The position of the western dip varies from ~248° (HiCiAO epoch, 2012 April 11) to ~282° (epoch of 2017 August 18), although the location of the minimum of the western dip from HiCiAO observations is not well constrained (see Table 2 and Figure 5). If we neglect this value, the western dip varies within ~32° (from 250° to 282°). The mean value, including all the epochs in the sample is $min_W = 265^\circ.9 \pm 13^\circ.0$, and the difference in location between the two dips is $min_W - min_E = 178^\circ.3 \pm 14^\circ.5$ (this last calculation neglects the ZIMPOL epoch).
<figure><img src="image_5.png" /></figure>

**Figure 4.** Left panels: radial mapping from 0′′25–0′′7 of all the panels shown in Figure 1. The color scale is linear and the maximum value taken in each case is 80% of the maximum. Right panel: azimuthal profile calculated from the mean values obtained between [0.35–0.50]′′. The shaded areas correspond to the uncertainty of the data and come from the standard deviation in the radial and azimuth divided by the square root of the number of pixels. The data is normalized to the value at zero degrees in each case. The dashed lines correspond to the minimum value obtained from the image between 0° and 150°, and 200° to 350°. For the ZIMPOL data only one minimum is shown.

*Dip amplitudes*—The amplitude of the eastern dip, derived from the near-infrared data sets, varies between 40% (epoch of 2017 September 16) to 90% (epoch of 2017 September 6). At the ZIMPOL epoch, the amplitude is 25%, but this value might be affected by the low

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

signal-to-noise ratio of the data. As for the location, the amplitude of the western dip appears to be more variable than that of the eastern dip, with values ranging from 31% to 95%. In five epochs, the amplitude of the two dips appear to be similar.
<table><caption>**Table 2**
Properties of the Dips</caption><thead><tr><th>Epoch</th><th>min<sub>E</sub> [°]<br>(image)</th><th>min<sub>W</sub> [°]<br>(image)</th><th>min<sub>E</sub> [°]<br>(fit)</th><th>min<sub>W</sub> [°]<br>(fit)</th><th>A<sub>E</sub><br>(fit)</th><th>A<sub>W</sub><br>(fit)</th><th>σ<sub>E</sub> [°]<br>(fit)</th><th>σ<sub>W</sub> [°]<br>(fit)</th></tr></thead><tbody><tr><td>2012 Apr 11</td><td>81.6</td><td>250.3</td><td>92.1<sup>+1.8</sup><sub>-1.8</sub></td><td>247.8<sup>+7.2</sup><sub>-7.9</sub></td><td>0.52<sup>+0.05</sup><sub>-0.06</sub></td><td>0.19<sup>+0.16</sup><sub>-0.06</sub></td><td>22.3<sup>+4.0</sup><sub>-2.7</sub></td><td>52.5<sup>+3.6</sup><sub>-2.7</sub></td></tr><tr><td>2015 Jun 10</td><td>42.2</td><td>...</td><td>44.9<sup>+4.7</sup><sub>-3.3</sub></td><td>...</td><td>0.25<sup>+0.04</sup><sub>-0.04</sub></td><td>...</td><td>28.2<sup>+8.8</sup><sub>-6.4</sub></td><td>...</td></tr><tr><td>2016 Jun 30</td><td>78.8</td><td>258.8</td><td>82.4<sup>+1.9</sup><sub>-2.0</sub></td><td>257.2<sup>+1.4</sup><sub>-1.3</sub></td><td>0.51<sup>+0.07</sup><sub>-0.06</sub></td><td>0.57<sup>+0.04</sup><sub>-0.05</sub></td><td>24.3<sup>+5.8</sup><sub>-3.9</sub></td><td>16.4<sup>+1.6</sup><sub>-1.3</sub></td></tr><tr><td>2017 Aug 13</td><td>92.8</td><td>267.2</td><td>83.1<sup>+2.4</sup><sub>-2.4</sub></td><td>274.4<sup>+1.3</sup><sub>-1.2</sub></td><td>0.50<sup>+0.14</sup><sub>-0.08</sub></td><td>0.67<sup>+0.04</sup><sub>-0.04</sub></td><td>31.4<sup>+10.3</sup><sub>-6.3</sub></td><td>21.9<sup>+2.1</sup><sub>-1.7</sub></td></tr><tr><td>2017 Aug 14</td><td>87.2</td><td>275.6</td><td>89.3<sup>+1.5</sup><sub>-1.6</sub></td><td>281.6<sup>+1.4</sup><sub>-1.3</sub></td><td>0.55<sup>+0.05</sup><sub>-0.05</sub></td><td>0.77<sup>+0.04</sup><sub>-0.04</sub></td><td>19.0<sup>+3.0</sup><sub>-2.2</sub></td><td>26.2<sup>+3.0</sup><sub>-2.4</sub></td></tr><tr><td>2017 Aug 17</td><td>90.0</td><td>270.0</td><td>90.7<sup>+1.6</sup><sub>-1.7</sub></td><td>280.5<sup>+2.8</sup><sub>-3.0</sub></td><td>0.52<sup>+0.05</sup><sub>-0.05</sub></td><td>0.31<sup>+0.04</sup><sub>-0.04</sub></td><td>17.3<sup>+2.7</sup><sub>-2.0</sub></td><td>20.7<sup>+4.1</sup><sub>-3.2</sub></td></tr><tr><td>2017 Aug 18</td><td>92.8</td><td>300.9</td><td>95.3<sup>+1.4</sup><sub>-1.4</sub></td><td>282.3<sup>+2.4</sup><sub>-2.4</sub></td><td>0.59<sup>+0.05</sup><sub>-0.05</sub></td><td>0.64<sup>+0.20</sup><sub>-0.15</sub></td><td>16.6<sup>+2.2</sup><sub>-1.6</sub></td><td>52.8<sup>+12.1</sup><sub>-10.6</sub></td></tr><tr><td>2017 Aug 22</td><td>98.4</td><td>247.5</td><td>98.5<sup>+1.5</sup><sub>-8.4</sub></td><td>249.7<sup>+1.6</sup><sub>-1.5</sub></td><td>0.44<sup>+0.10</sup><sub>-0.07</sub></td><td>0.47<sup>+0.04</sup><sub>-0.04</sub></td><td>9.7<sup>+3.5</sup><sub>-1.8</sub></td><td>15.5<sup>+1.5</sup><sub>-1.3</sub></td></tr><tr><td>2017 Sep 04</td><td>87.2</td><td>267.2</td><td>83.8<sup>+1.5</sup><sub>-1.4</sub></td><td>280.7<sup>+2.0</sup><sub>-2.6</sub></td><td>0.47<sup>+0.05</sup><sub>-0.06</sub></td><td><strong>0.78</strong><sup>+0.10</sup><sub>-0.09</sub></td><td>11.8<sup>+2.0</sup><sub>-1.5</sub></td><td>38.5<sup>+6.6</sup><sub>-6.2</sub></td></tr><tr><td></td><td></td><td></td><td></td><td><strong>268.0</strong><sup>+2.6</sup><sub>-2.6</sub></td><td></td><td></td><td>0.58<sup>+0.15</sup><sub>-0.23</sub></td><td><strong>16.3</strong><sup>+8.3</sup><sub>-6.0</sub></td></tr><tr><td>2017 Sep 06</td><td>92.8</td><td>264.4</td><td>77.2<sup>+1.4</sup><sub>-1.6</sub></td><td>288.7<sup>+2.1</sup><sub>-2.2</sub></td><td>0.90<sup>+0.06</sup><sub>-0.09</sub></td><td><strong>0.95</strong><sup>+0.04</sup><sub>-0.08</sub></td><td>33.8<sup>+3.2</sup><sub>-3.1</sub></td><td>54.9<sup>+3.8</sup><sub>-4.9</sub></td></tr><tr><td></td><td></td><td></td><td></td><td><strong>266.4</strong><sup>+2.2</sup><sub>-2.1</sub></td><td></td><td></td><td>0.62<sup>+0.18</sup><sub>-0.23</sub></td><td><strong>18.4</strong><sup>+6.6</sup><sub>-5.0</sub></td></tr><tr><td>2017 Sep 16</td><td>87.2</td><td>247.5</td><td>92.1<sup>+1.9</sup><sub>-2.1</sub></td><td>250.7<sup>+2.7</sup><sub>-2.2</sub></td><td>0.40<sup>+0.05</sup><sub>-0.05</sub></td><td>0.45<sup>+0.30</sup><sub>-0.07</sub></td><td>15.9<sup>+3.76</sup><sub>-2.1</sub></td><td>12.1<sup>+6.4</sup><sub>-2.2</sub></td></tr></tbody></table>

**Note.** Best parameters from fitting a Lorentzian profile (Equation (1)), which are obtained for a range of PA of 96° to 350°. For the epochs of September 4 and 6, we report a second fit result considering a PA range from 236° to 300°. For these two epochs, the values indicated in bold are the ones considered in the analysis. In addition, we report the angle at which the minimum is obtained from the image (i.e., dashed lines in Figure 4).

**Dip widths**—The width of the eastern dip varies from ~10° (2017 August 22) to ~34° (epoch of 2017 September 6, also the one with the largest amplitude). The width of the western dip varies from ~12° (2017 September 16) to ~53° (2017 August 18). The western dip of the epoch on 2017 August 18 seems to be a composition of two different dips, and this may be the reason why the fit of this dip gives as a result of a very wide dip. There does not seem to be a relationship between the variation of the width of the two dips, and, on average, if we only consider the SPHERE data, the averaged width of the eastern and western dips are ~20° and ~22°, respectively.

In addition, we checked the radial profile outside of the dips, along the ring, from 140° to 200°. In most of the cases, the surface brightness distribution is flat varying within 20% of the reference value (at 0°). In the epoch on 2017 August 17, the surface brightness distribution decreases with PA, from values of 1.1 to 0.88. For the epoch on 2017 September 6, the surface brightness distribution monotonically increases with PA, from values of 1.22 to 1.45. We note that this is the epoch with the largest amplitude for the two dips.

### 3.3. Comparison with ALMA Observations

In Pinilla et al. (2018), we performed an analysis of the dust morphology of several transition disks, including J1604, that were observed with ALMA in the (sub-)millimeter regime. This analysis was done in the visibility plane to characterize the total flux, cavity size, and shape of the ring-like structure. Motivated by models of dust trapping in pressure bumps, we fitted a radially asymmetric Gaussian ring for the millimeter intensity, that is, a Gaussian ring whose inner and outer widths differ. For J1604, based on observations obtained with ~0′′26 ×0′′22 resolution, the inner and outer widths of the Gaussian from the best fit are ~0′′08 and ~0′′14, respectively, while the Gaussian peaks at ~0′′55 ± 0′′01. The cavity size is well resolved while the width of the ring (0′′22) remains unresolved.

Figure 7 shows the radial profile of the surface brightness, which is normalized to the peak of emission for SPHERE

versus ALMA. We randomly chose the epoch on 2017 August 14 as a reference of our IRDIS observations. The ring observed in scattered light resides inside the cavity-observed with ALMA. It is expected that micron-sized particles also exist inside the ring observed at millimeter emission (de Juan Ovelar et al. 2013; Pinilla et al. 2016b). However, shadowing from the ring itself can cause the ring observed in scattered light to be detected fully inside the ALMA ring (see, e.g., Figure 3 in Dullemond & Monnier 2010). This shadowing effect supports the steepness of the surface brightness beyond the peak, as explained in Section 3.1. There is a significant separation between the two peaks (0.43 versus 0′′55, SPHERE versus ALMA). To compare with models of particle trapping by embedded massive planets as discussed in Section 4, we calculate the “wall” of the ring observed in scattered light ($w_{SL}$, defined as the radial location where the flux has increased by half from the minimum in the cavity and the peak of the ring, de Juan Ovelar et al. 2013), we obtain that the wall location for the 2017 August 14 epoch is ~0′′36. This implies that the ratio of $w_{SL}$ and the peak of the millimeter emission is ~0.65.

With a ~0′′25 resolution, Dong et al. (2017) analyzed $^{12}$CO, $^{13}$CO, and $^{C^{18}O}$J = 2–1 line emission from ALMA observations of J1604. They concluded that their gas observations are consistent with a gas cavity that is smaller than the millimeter-dust cavity (with an upper limit for the inner radius at 0′′10). From their thermo-chemical models, they suggested that the gas surface density smoothly increases from 0′′10 to the peak of the millimeter emission and they exclude a sharp transition or double-drop models (i.e., models that assume two localized reductions) for the gas surface density. According to their results, the gas is depleted inside the cavity by 2 to 4 orders of magnitude. Therefore, the ring observed in our scattered light observations lies in between the minimum of the gas surface density (inside 0′′10) and the peak of the millimeter emission (0′′55).
<figure><img src="image_7.png" /><figcaption>**Figure 5.** Left panels: Gaussian MCMC best fit of each dip overlap with the uncertainty of the azimuthal profile for each epoch as shown in Figure 4. The fit shown in dashed red lines of the epochs 2017 September 4 and 6 corresponds to a Gaussian fit, but the PA range is taken from 236° to 300° for fitting the western dip of the epochs, instead of 196°–350° as for the rest of the fits.</figcaption></figure>

## 4. Discussion

In this section, we discuss the potential origin of the observed cavity and the evidence for particle trapping from

multiwavelength observations. In addition, we also discuss different origins for the shadows and their variability, which can be potentially connected to the dipper nature of J1604.

### 4.1. Origin of the Cavity and Evidence of Dust Trapping

As a member of the Upper Scorpius OB association (one of the oldest star-forming regions, 5–11 Myr, that host protoplanetary disks), J1604 is an excellent target to investigate a critical stage when disk dissipation should be almost over (e.g., Williams & Cieza 2011). Its spectroscopic signatures indicate very low accretion rates (e.g., Dahm et al. 2012), recently confirmed with X-Shooter spectra (C. F. Manara et al. 2018, in preparation), with $log_{10} \dot{M}_{\star} [M_{\odot} \text{yr}^{-1}] = -10.54$. This low value of the accretion rate may indicate a strongly gas depleted inner region, which is supported by the analysis of molecular gas (Dong et al. 2017). Nonetheless, J1604 evidences variable NIR and MIR excess (Dahm & Carpenter 2009), which can originate from an optically thick ring located at small (~0.1 au) orbital radii. The potential existence of a dusty inner disk (also suggested from the shadows and the dimming events on the optical light curve), in addition to the low-mass accretion rate, and the comparison between the distribution of the gas and small/large grains (Section 3.3) can give constraints on the origin of the cavity.

The spatial segregation between gas and small/large grains suggest particle trapping in pressure maxima, which naturally explains why this disk can remain quite massive in dust despite its old age. If a single planet is the primary cause of the cavity, the ratio of the wall of the ring observed in scattered light ($w_{SL}$), and the peak of the millimeter emission can hint at the mass of the embedded planet. de Juan Ovelar et al. (2013) combined dust evolution models with hydrodynamical simulations of planet disk interaction to obtain the density distribution of different grain sizes when a massive planet in a circular orbit is embedded in the disk and filter dust grains of different sizes. These results, combined with radiative transfer predictions, allow us to infer the mass of the planet from observations of the ring-like structure of transition disks at different wavelengths (see Figure 8 in de Juan Ovelar et al. 2013). The obtained value of the ratio of the wall of the ring observed in scattered light of 0.65 suggests a planet mass of at least $4 M_{Jup}$.

Observations of CO and its isotopologues suggest that the location of this hypothetical planet or companion should be within ~15 au (assuming a distance of 150.2 pc, Gaia Collaboration et al. 2018), which corresponds to the location of the minimum of the gas surface density ($10^{-2} \text{pc}^{-2}$). This location, however, remains unresolved from ALMA observations. The mass and location of such a hypothetical planet are below the upper limits on close companions derived from Kraus et al. (2008; ~70 $M_{Jup}$ between 10 and 20 au), Ireland et al. (2011; ~83 $M_{Jup}$ within 45 au), and Canovas et al. (2017; 2–3 $M_{Jup}$ from 22 to 115 au). Nevertheless, the mass of this potential planet seems to be too high for the system to maintain an inner disk at the age of Upper Sco. Such a massive planet will block most of the dust (of all sizes) at the outer edge of the planetary gap such that after ~5 Myr of evolution no dusty material would remain in the inner disk (Pinilla et al. 2016b).

Alternatively, the cavity can form due to multiple planets, which leads to wider and shallower cavities (Duffell & Dong 2015), and allow a flow of dust from the outer to the inner disk. This dust is expected to drift to the very inner regions (~1 au) and can pile up at these small orbital radii
<figure><img src="image_8.png" /><figcaption>**Figure 6.** Variations of the east and west minimum properties (i.e., results reported in Table 2). From left to right: location, amplitude, and width of each dip. In the left panel, the dashed lines correspond to the mean values including all epochs (i.e., $min_E = 83^{\circ}7 \pm 13^{\circ}7$, $min_W = 265^{\circ}9 \pm 13^{\circ}0$, and $min_W - min_E = 178^{\circ}3 \pm 14^{\circ}5$). For the western dip and the epochs of 2017 September 4 and 6, we take the location of the minimum and the width from the fit that takes the PA range from 236° to 300°, while for the amplitude, we take a PA range of 196°–350°, as for the rest of the fits.</figcaption></figure>

<figure><img src="image_9.png" /><figcaption>**Figure 7.** Radial profile of the surface brightness, normalized to the peak of emission for SPHERE (black; 2017 August 14) vs. ALMA (red, Pinilla et al. 2018). The horizontal gray line corresponds to the beam major axis from ALMA observations.</figcaption></figure>

because of their low drift velocities near the snow line (Pinilla et al. 2016b). Assuming that the disk effective temperature results from stellar irradiation and accretion, the location of the snow line is expected to be at ~1 au, considering a stellar luminosity of ~0.58 $L_\odot$, a stellar temperature of 4500 K, and a stellar mass of ~1 $M_\odot$ (Dahm & Carpenter 2009; Mathews et al. 2012). However, in this scenario, the gas in the cavity would not be as depleted as suggested by the ALMA observations and the accretion rate would not be that low.

The possibility of a close binary is excluded from high-resolution optical spectroscopy, which does not show signs of a double-line spectroscopic binary (Dahm et al. 2012). However, a massive companion in a wide eccentric orbit (as in HD 142527) still remains as a possibility. In this case, the inner disk may be filled by streamers bridging it to the outer disk. When a companion is massive enough (>5 $M_{Jup}$, depending on the disk viscosity), the disk becomes eccentric and streamers are more efficient at transporting material (Ataiee et al. 2013). As a result, accretion of material onto the planet and flows of material from the outer disk to the inner disk can be enhanced (Kley & Dirksen 2006; Ragusa et al. 2017). In this

case, an inner disk can be maintained for longer times of evolution, and the accretion is expected to be variable with time for both the central star and the companion (Muñoz & Lai 2016). Hence, an eccentric planet could be a viable solution to allow the inner disk to be replenished and at the same time be consistent with the lower limit for the planet mass inferred from the dust segregation, as discussed above.

In the inner disk, where the gas density is depleted, kilometer-sized planetesimals would be completely decoupled from the gas and they would not experience fast inward radial drift. These planetesimals can potentially collide, recreating a belt of micron-sized particles in the inner disk. These micron-sized particles are expected to grow efficiently, unless they reside (or grow and drift) inside the snow line where fragmentation of silicates is efficient, keeping a sufficient amount of small grains close to the star. This scenario would imply that the inner disk extent is very small (within the snow line, i.e., ~0.1–1 au).

Other possibilities for the formation of the cavity and particle trapping include photoevaporation and nonideal magnetohydrodynamical (MHD) effects. On one hand, photoevaporation is consistent with the low accretion rate of J1604, and the detection of the [O I] line (C. F. Manara et al. 2018, in preparation). Ercolano et al. (2018) also demonstrated that X-ray photoevaporation in a disk with a moderate gas depletion of carbon and oxygen can create cavities as big as 100 au with a large range of accretion rates ($M_* \sim 10^{-11} - 10^{-8} M_\odot \text{ yr}^{-1}$). However, as pointed out by Dong et al. (2017), models of photoevaporation predict a sharp cavity edge in gas (e.g., see Figure 1 in Ercolano et al. 2018; see also Alexander et al. 2006), which contradicts the results from the analysis of the CO observations from ALMA.

Nonideal MHD effects, such as dead zones can also create cavities observable at different wavelengths (e.g., Flock et al. 2015; Ruge et al. 2016) because particle trapping can occur at the outer edge of the dead zone where there is a bump in the gas density profile due to the change of accretion from the dead to the active MRI regions. This scenario, however, predicts that cavities at short and long wavelength should be of similar size (Pinilla et al. 2016a). A solution to create spatial segregation, as observed in J1604, is the inclusion of an MHD wind (possibly traced by the [O I] emission line; C. F. Manara et al. 2018, in preparation) to the dead zone models, which can create a large difference in the distribution of gas and small/large particles. In
this case, the gas surface density inside the cavity can be depleted
by several orders of magnitude and increases smoothly with
radius (see Figure 6 in Pinilla et al. 2016a), as suggested for
J1604 by ALMA observations.

## 4.2. *Shadows and Their Variability*

Dips in scattered light images have been interpreted as
shadowing from a misaligned inner disk (Marino et al. 2015;
Benisty et al. 2017). Assuming an outer disk inclination and
position angle of $i = 6^\circ$ and $PA = 80^\circ$ (Pinilla et al. 2018),
respectively, and a disk aspect ratio of 0.1 for the scattering
surface at the ring radius (Dong et al. 2017), we find that the
inner disk should be close to edge-on, leading to a misalign-
ment between inner and outer disks of $\sim 70^\circ-90^\circ$ (see
Equations in Min et al. 2017). A strongly inclined inner disk
is consistent with the dipper activity of this object, by which
dimming events can be caused by patches of dusty circum-
stellar material that repeatedly occult the star as they cross the
line of sight.

The intersection of the planes of the inner and outer disks
defines the PA of the shadows (assuming that they are razor-
thin). Due to the finite scale height of the inner disk, the
shadows appear as broad dark regions, and their widths can in
principle constrain the scale height of the inner disk. If the
relative orientation of the disks is fixed with time, the PA of the
shadows should not vary. However, we find that it varies within
$\sim \pm 14^\circ$, from the estimate of the (local) minimum value of the
surface brightness. It is possible that the PA variations are
related to the variations of the dips' shapes and widths, which
can in turn modify the location of the minimum surface
brightness that we estimate (or PA). Such variations in the
widths and PA of the shadows imply that they are not caused
by a symmetric inner disk with a constant misalignment with
respect to the outer disk. Instead, it is likely that the inner disk
is highly structured and asymmetric, and that its scale height
varies with time, in very short timescales (within a day). In
addition, the fact that the shadow properties and variability are
very different for the eastern and western dips, also supports an
asymmetric morphology of the inner disk.

Both the fast dynamics and asymmetric morphology of the
shadows are likely connected with the aperiodicity of the stellar
dimming events (assumed to be due to an inner disk warp),
which can change as much as $\sim 60\%$ also in timescales of a
few days in dipper objects (Ansdell et al. 2016; Cody &
Hillenbrand 2018). These variations likely originate from
variations of the inner disk scale height of the order of 10% or
more (McGinnis et al. 2015).

It is therefore not surprising that the width of the shadows
observed in J1604 significantly vary on day to week timescales,
as they directly reflect the intrinsic variations of the inner disk
scale height on these timescales. Furthermore, no clear
correlation is expected between the widths of the east and
west shadows as long as the inner disk scale height varies on a
timescale shorter than a few days, which appears to be the case
for most dippers (McGinnis et al. 2015). The maximum
amplitude of the dips, however, depends on the optical
thickness of the inner disk warp, i.e., the dust properties. The
seemingly correlated amplitudes of variability of the east and
west shadows might mean that the dust properties (in
particular, opacity) change on longer timescales than the inner
disk warp shape.

A misaligned stellar magnetic field with respect to the
rotation axis of the star can create a magnetically warped inner
disk edge as in the case of the disk around AA Tau (Bouvier
et al. 1999). Such a warp occults AA Tau periodically, and
accounts well for the spectral variability of this system (e.g.,
Bouvier et al. 1999, 2003, 2007; Ménard et al. 2003). This is
because AA Tau is a fully convective star that hosts a very
strong (2–3 kG) dipolar magnetic field, which maximizes the
star–disk interaction, and produces a fairly stable, though
dynamical, inner disk warp. However, J1604 is a more massive
star with a well-developed inner radiative core. Partly radiative
pre-main-sequence stars tend to exhibit weaker fields that are
mostly octupolar (e.g., V2129 Oph, Donati et al. 2011; Gregory
et al. 2012). In this case, more complex accretion flows are
expected (Alencar et al. 2012), possibly leading to more
unstable and aperiodic star–disk interactions. If the octupole
dominates at the disk level, it is possible to have an asymmetric
disk warp with a complex perturbation of the inner disk scale
height as a function of azimuth. The observed variations of the
shadows PA could therefore result from the varying vertical
shape of the inner disk warp as a function of azimuth. Thus, a
complex magnetic field geometry, coupled with a relatively
weak accretion rate, can produce the strong and irregular
variability seen in the shadows of J1604. Nonetheless, there are
currently no observational estimates of the magnetic field
strength and topology for J1604, and whether a complex
magnetic field topology can be the origin of the observed
variability and shadows still has to be investigated with a
dedicated spectro-polarimetric campaign.

Apart from a misaligned magnetic field, the presence of a
yet-undetected inclined massive companion in the cavity of
J1604 could be responsible for a large misalignment between
inner and outer disks. Indeed, in the case of HD 142527, there
is strong evidence that a $0.1 M_\odot$ companion in an eccentric
orbit with an inclination of $\sim 125^\circ$ (Lacour et al. 2016) is the
cause of the misalignment of the inner disk with respect to the
outer disk, which can also explain most of the observed
properties of this disk (Price et al. 2018). While a single
massive planet or companion in a coplanar and circular orbit
might not explain the spatial segregation of the gas and small/
large grains (see Section 4.1), we cannot exclude the presence
of such an inclined and eccentric companion in the cavity of
J1604.

The variation of the shadows' properties might be to a large
number of dust clumps, orbiting with a large range of
inclinations, maybe due to planetesimal collision in the
depleted inner regions. The effect of light travel time can
create shadows with a large range of morphologies, from arc-
shaped to spiral arms, depending on the disk scale height (flat
versus flared disk) and the disk inclination (Kama et al. 2016).
In the case of J1604, due to the narrow extent of the ring in the
outer disk, the expected shadows by clumps in the inner disk
would not look like spirals but instead, as localized dips
within the ring of emission. In any case, if this mechanism is
responsible for the shadows, these clumps must be very
dynamic (changing position and morphology in day timescales)
to explain the observed variability of the dips, and lead to a
sufficient amount of small dust that would lead to a large radial
optical depth. Detailed modeling to assess if this is possible
is required.

Another possibility to explain the presence of shadows
and the observed variability is to consider the effects of
magnetohydrodynamic instabilities. The Parker instability that
occurs when amplified magnetic fields (by disk dynamo) can
escape from the disk due to magnetic buoyancy (Takasao
et al. 2018). As a result of angular momentum exchange
mediated by magnetic fields, the velocity and density above the
disk increases and magnetic fields can escape due to magnetic
buoyancy. As a consequence, due to the MRI-driven turbulence
and eruptions of the magnetic field, the density near the disk
surface can significantly fluctuate spatially and temporally. In
addition, due to MRI turbulence, the upper disk layers that are
magnetically supported can carry dust grains at high altitude,
and lead to shadows with typical timescales from half to
a tenth of the Keplerian period at the inner disk (Turner
et al. 2010), potentially explaining aperiodic dimming and
shadowing events.

In all scenarios, if the shadows are steady over orbital
timescales at the ring radius, and the cooling timescale is
comparably fast, they could lead to a decrease of the dust and
gas temperature in the ring, and also appear as dips in the
millimeter images. There is marginal evidence for shadows
at millimeter emission from ALMA observations (Dong
et al. 2017), similarly to the case of DoAr 44 (Casassus et al.
2018). This aspect will be further investigated in R. A. Loomis
et al. (2018, in preparation).

**5. Conclusions**

We present new VLT/SPHERE polarimetric differential imaging of the transition disk around the dipper star J1604. We gathered a total of 11 epochs of scattered light imaging that span days, weeks, months, and years (Figure 1). Our findings are:

1. All the scattered light epochs show a ring-like emission that peaks at ~0"/43 from the star (Figure 2). The morphology of the outer tail of this ring rapidly changes with time (Figure 3 and Table 1). The width of this ring remains unresolved in our observations (as in the case of millimeter-observations with ALMA), with a value of <~0"/13. This ring lies inside the cavity resolved at millimeter emission, which also shows a ring-like structure but peaking at 0"/55 (Figure 7). This spatial segregation can be a natural result of particle trapping in pressure bumps, created by, for example, embedded planet(s).

2. In the case of a single massive planet being the origin of cavity, at least a 4 M<sub>Jup</sub> mass planet is required inside the cavity to create the observed segregation of small and large grains. However, such a planet cannot explain the gas surface density structure inferred from ALMA observations and the existence of an inner dust optically thick belt, at the age of Upper Sco, to explain the (variable) NIR excess. Potential alternatives include the possibility of a dead zone and a MHD wind acting together in the evolution of J1604.

3. We detect two clear dips of emission along the ring observed in scattered light (Figures 2 and 4). Both dips are highly variable in amplitude and width (Figure 6). The western dip is in general more variable than the eastern dip. For the eastern dip, the amplitude varies from 40% to 90% and its width varies from 10° to 34°. For the western dip, the amplitude varies from 31% to 95% and its width varies from 12° to 53°. From the 11 epochs, the mean position of the dips are ~83°/7 ± 13°/7 and ~265°/9 ± 13°/0 for the eastern and western dip,

respectively. The averaged separation between the dips
is 178°/3 ± 14°/5.

4. Assuming that these dips are shadowing from a misaligned inner disk, we find that the misalignment between the inner and the outer disk is very large (~polar), and similar to the values found for the other two transition disks: HD 142527 and HD 100453. Current available observations do not provide constraints on what can be the origin of the warp, and it remains an open question if it is due to a companion in an highly inclined orbit (as for the case of HD 142527), a misaligned stellar magnetic field with respect to the rotation axis of the star (as for the case of AA Tau), or other alternatives, such as dusty asymmetric clumps from forming planetesimals. The variability of the morphology of the shadows, along with the rapid variations of the morphology of the ring tail, suggest that the innermost regions are highly dynamic and complex.

5. The misalignment between the inner and the outer disk reconciles the dipper activity of J1604 since dimming events in light curves are mainly observed in highly inclined disks. Future demographic studies are needed to test if close-to-face-on dippers show shadows in scattered light and vice versa.

Future VLTI, spectro-polarimetric campaign, high-resolution ALMA observations and simultaneous optical/IR light curves can help us to better characterize the inner disk causing the dipping events and shadows and the stellar magnetic field. In addition, ALMA line observations can also provide potential variations of the Keplerian motion in the inner regions of J1604 from different molecular lines, which will help us to characterize the warp and the possible presence of massive companion(s) in the disk responsible for the observed cavity.

We acknowledge the referee for a constructive report that
helped improve the paper. The authors are thankful to A. Natta,
S. Facchini, A. Pohl, and M. Ansdell for interesting discussions
on J1604. We thank S. Mayama for sharing the HiCIAO image.
P.P. acknowledges support by NASA through Hubble Fellowship
grant HST-HF2-51380.001-A awarded by the Space Telescope
Science Institute, which is operated by the Association of
Universities for Research in Astronomy, Inc., for NASA, under
contract NAS 5-26555. M.B. acknowledges funding from ANR
of France under contract number ANR-16-CE31-0013 (Planet
Forming disks). C.F.M. acknowledges support through the ESO
Fellowship. J.B. acknowledges funding from the European
Research Council (ERC) under the European Union's Horizon
2020 research and innovation programme (grant agreement No.
742095; SPIDI: Star-Planets-Inner Disk-Interactions). C.D.
acknowledges funding from the Netherlands Organisation for
Scientific Research (NWO) TOP-1 grant, project number
614.001.552. R.A.L. gratefully acknowledged funding from
NRAO Student Observing Support.

*Software:* Numpy (Van Der Walt et al. 2011), SciPy (Jones et al. 2001), Matplotlib (Hunter 2007), emcee (Foreman-Mackey et al. 2013).

**Appendix**
**Observational Notes**

Table 3 is an observing log. For each epoch, it gives the
instrument and its central wavelength. The last column
<table><caption>**Table 3**
Observing Log</caption><thead><tr><th>Epoch</th><th>Instrument</th><th>Central λ[μm]</th><th>FWHM [mas × mas]</th></tr></thead><tbody><tr><td>2012 Apr 11</td><td>HiCiAO</td><td>1.6</td><td>70 × 70</td></tr><tr><td>2015 Jun 10</td><td>SPHERE/ZIMPOL</td><td>0.626</td><td>40 × 33</td></tr><tr><td>2016 Jun 30</td><td>SPHERE/IRDIS</td><td>1.625</td><td>56 × 55</td></tr><tr><td>2017 Aug 13</td><td>SPHERE/IRDIS</td><td>1.245</td><td>40 × 46</td></tr><tr><td>2017 Aug 14</td><td>SPHERE/IRDIS</td><td>1.245</td><td>39 × 45</td></tr><tr><td>2017 Aug 17</td><td>SPHERE/IRDIS</td><td>1.245</td><td>40 × 46</td></tr><tr><td>2017 Aug 18</td><td>SPHERE/IRDIS</td><td>1.245</td><td>43 × 44</td></tr><tr><td>2017 Aug 22</td><td>SPHERE/IRDIS</td><td>1.245</td><td>40 × 40</td></tr><tr><td>2017 Sep 04</td><td>SPHERE/IRDIS</td><td>1.245</td><td>...</td></tr><tr><td>2017 Sep 06</td><td>SPHERE/IRDIS</td><td>1.245</td><td>40 × 46</td></tr><tr><td>2017 Sep 16</td><td>SPHERE/IRDIS</td><td>1.245</td><td>42 × 44</td></tr></tbody></table>

Note. The 2017 September 4 epoch is omitted because of its FLUX image has a too low signal-to-noise ratio to be fitted. The HiCIAO PSF was not fitted, we report the values given by Mayama et al. (2012).

provides the FWHM (in milliarcseconds mas) of a 2D Gaussian fit to the FLUX image of the SPHERE data, and to the core of the PSF of the HiCIAO observations from Mayama et al. (2012).

## ORCID iDs

P. Pinilla https://orcid.org/0000-0001-8764-1780

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

J. de Boer https://orcid.org/0000-0002-1028-7522

C. F. Manara https://orcid.org/0000-0003-3562-262X

C. Dominik https://orcid.org/0000-0002-3393-2459

R. A. Loomis https://orcid.org/0000-0002-8932-1219

## References

Alencar, S. H. P., Bouvier, J., Walter, F. M., et al. 2012, A&A, 541, A116

Alencar, S. H. P., Teixeira, P. S., Guimarães, M. M., et al. 2010, A&A, 519, A88

Alexander, R. D., Clarke, C. J., & Pringle, J. E. 2006, MNRAS, 369, 229

Ansdell, M., Gaidos, E., Williams, J. P., et al. 2016, MNRAS, 462, L101

Ataiee, S., Pinilla, P., Zsom, A., et al. 2013, A&A, 553, L3

Avenhaus, H., Quanz, S. P., Schmid, H. M., et al. 2014, ApJ, 781, 87

Barenfeld, S. A., Carpenter, J. M., Ricci, L., & Isella, A. 2016, ApJ, 827, 142

Benisty, M., Juhasz, A., Boccaletti, A., et al. 2015, A&A, 578, L6

Benisty, M., Juhasz, A., Facchini, S., et al. 2018, A&A, in press (arXiv:1809.01082)

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

Beuzit, J.-L., Feldt, M., Dohlen, K., et al. 2008, Proc. SPIE, 7014, 18

Boccaletti, A., Thalmann, C., Lagrange, A.-M., et al. 2015, Natur, 526, 230

Bodman, E. H. L., Quillen, A. C., Ansdell, M., et al. 2017, MNRAS, 470, 202

Bouvier, J., Alencar, S. H. P., Boutelier, T., et al. 2007, A&A, 463, 1017

Bouvier, J., Chelli, A., Allain, S., et al. 1999, A&A, 349, 619

Bouvier, J., Grankin, K. N., Alencar, S. H. P., et al. 2003, A&A, 409, 169

Canovas, H., Hardy, A., Zurlo, A., et al. 2017, A&A, 598, A43

Casassus, S., Avenhaus, H., Perez, S., et al. 2018, MNRAS, 477, 5104

Casassus, S., Marino, S., Pérez, S., et al. 2015, ApJ, 811, 92

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

Cody, A. M., & Hillenbrand, L. A. 2018, AJ, 156, 71

Cody, A. M., Stauffer, J., Baglin, A., et al. 2014, AJ, 147, 82

Dahm, S. E., & Carpenter, J. M. 2009, AJ, 137, 4024

Dahm, S. E., Slesnick, C. L., & White, R. J. 2012, ApJ, 745, 56

Dawson, R. I., Murray-Clay, R. A., & Fabrycky, D. C. 2011, ApJL, 743, L17

de Boer, J., Salter, G., Benisty, M., et al. 2016, A&A, 595, A114

de Juan Ovelar, M., Min, M., Dominik, C., et al. 2013, A&A, 560, A111

Debes, J. H., Jang-Condell, H., & Schneider, G. 2016, ApJL, 819, L1

Dohlen, K., Langlois, M., Saisse, M., et al. 2008, Proc. SPIE, 7014, 70143L

Donati, J.-F., Bouvier, J., Walter, F. M., et al. 2011, MNRAS, 412, 2454

Dong, R. 2015, ApJ, 810, 6

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

Duffell, P. C., & Dong, R. 2015, ApJ, 802, 42

Dullemond, C. P., Dominik, C., & Natta, A. 2001, ApJ, 560, 957

Dullemond, C. P., & Monnier, J. D. 2010, ARA&A, 48, 205

Ercolano, B., Weber, M. L., & Owen, J. E. 2018, MNRAS, 473, L64

Facchini, S., Juhász, A., & Lodato, G. 2018, MNRAS, 473, 4459

Facchini, S., Lodato, G., & Price, D. J. 2013, MNRAS, 433, 2142

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

Foucart, F., & Lai, D. 2013, ApJ, 764, 106

Fujiwara, H., Honda, M., Kataza, H., et al. 2006, ApJL, 644, L133

Fukagawa, M., Tamura, M., Itoh, Y., et al. 2010, PASJ, 62, 347

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

Golimowski, D. A., Ardila, D. R., Krist, J. E., et al. 2006, AJ, 131, 3109

Grady, C. A., Schneider, G., Hamaguchi, K., et al. 2007, ApJ, 665, 1391

Gregory, S. G., Donati, J.-F., Morin, J., et al. 2012, ApJ, 755, 97

Hunter, J. D. 2007, CSE, 9, 90

Ireland, M. J., Kraus, A., Martinache, F., Law, N., & Hillenbrand, L. A. 2011, ApJ, 726, 113

Jones, E., Oliphant, T., Peterson, P., et al. 2001, SciPy: Open Source Scientific Tools for Python

Kama, M., Pinilla, P., & Heays, A. N. 2016, A&A, 593, L20

Kley, W., & Dirksen, G. 2006, A&A, 447, 369

Kraus, A. L., Ireland, M. J., Martinache, F., & Lloyd, J. P. 2008, ApJ, 679, 762

Lacour, S., Biller, B., Cheetham, A., et al. 2016, A&A, 590, A90

Loomis, R. A., Öberg, K. I., Andrews, S. M., & MacGregor, M. A. 2017, ApJ, 840, 23

Lubow, S. H., & Martin, R. G. 2016, ApJ, 817, 30

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

Martin, R. G., Lubow, S. H., Nixon, C., & Armitage, P. J. 2016, MNRAS, 458, 4345

Mathews, G. S., Williams, J. P., & Ménard, F. 2012, ApJ, 753, 59

Mayama, S., Hashimoto, J., Muto, T., et al. 2012, ApJL, 760, L26

McGinnis, P. T., Alencar, S. H. P., Guimarães, M. M., et al. 2015, A&A, 577, A11

Ménard, F., Bouvier, J., Dougados, C., Mel'nikov, S. Y., & Grankin, K. N. 2003, A&A, 409, 163

Min, M., Stolker, T., Dominik, C., & Benisty, M. 2017, A&A, 604, L10

Muñoz, D. J., & Lai, D. 2016, ApJ, 827, 43

Muto, T., Grady, C. A., Hashimoto, J., et al. 2012, ApJL, 748, L22

Nealon, R., Price, D. J., & Nixon, C. J. 2015, MNRAS, 448, 1526

Nixon, C. J., & King, A. R. 2012, MNRAS, 421, 1201

Oh, D., Hashimoto, J., Tamura, M., et al. 2016, PASJ, 68, L3

Owen, J. E., & Lai, D. 2017, MNRAS, 469, 2834

Papaloizou, J. C. B., & Lin, D. N. C. 1995, ApJ, 438, 841

Papaloizou, J. C. B., & Pringle, J. E. 1983, MNRAS, 202, 1181

Papaloizou, J. C. B., & Terquem, C. 1995, MNRAS, 274, 987

Pérez, L. M., Isella, A., Carpenter, J. M., & Chandler, C. J. 2014, ApJL, 783, L13

Pinilla, P., Benisty, M., & Birnstiel, T. 2012, A&A, 545, A81

Pinilla, P., de Boer, J., Benisty, M., et al. 2015, A&A, 584, L4

Pinilla, P., Flock, M., Ovelar, M. d. J., & Birnstiel, T. 2016a, A&A, 596, A81

Pinilla, P., Klarmann, L., Birnstiel, T., et al. 2016b, A&A, 585, A35

Pinilla, P., Tazzari, M., Pascucci, I., et al. 2018, ApJ, 859, 32

Pohl, A., Benisty, M., Pinilla, P., et al. 2017, ApJ, 850, 52
Price, D. J., Cuello, N., Pinte, C., et al. 2018, MNRAS, 477, 1270

Ragusa, E., Dipierro, G., Lodato, G., Laibe, G., & Price, D. J. 2017, MNRAS, 464, 1449

Romanova, M. M., Ustyugova, G. V., Koldoba, A. V., & Lovelace, R. V. E. 2013, MNRAS, 430, 699

Rosenfeld, K. A., Qi, C., Andrews, S. M., et al. 2012, ApJ, 757, 129

Ruge, J. P., Flock, M., Wolf, S., et al. 2016, A&A, 590, A17

Sauvage, J., Fusco, T., Petit, C., et al. 2014, Proc. SPIE, 9148, 914847

Scaringi, S., Manara, C. F., Barenfeld, S. A., et al. 2016, MNRAS, 463, 2265

Schmid, H. M., Joos, F., & Tschan, D. 2006, A&A, 452, 657

Schneider, P. C., Manara, C. F., Facchini, S., et al. 2018, A&A, 614, A108

Shakura, N. I., & Sunyaev, R. A. 1973, A&A, 24, 337

Stolker, T., Dominik, C., Avenhaus, H., et al. 2016, A&A, 595, A113

Takasao, S., Tomida, K., Iwasaki, K., & Suzuki, T. K. 2018, ApJ, 857, 4

Turner, N. J., Carballido, A., & Sano, T. 2010, ApJ, 708, 188

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

Van Der Walt, S., Colbert, S. C., & Varoquaux, G. 2011, CSE, 13, 22

Walsh, C., Daley, C., Facchini, S., & Juhász, A. 2017, A&A, 607, A114

Whitney, B. A., & Hartmann, L. 1992, ApJ, 395, 529

Williams, J. P., & Cieza, L. A. 2011, ARA&A, 49, 67

Xiang-Gruess, M., & Papaloizou, J. C. B. 2013, MNRAS, 431, 1320

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