FIRST SCATTERED-LIGHT IMAGE OF THE DEBRIS DISK AROUND HD 131835 WITH THE GEMINI PLANET IMAGER

Debris disks are the remnant products of planet formation processes; within them, planetesimals collisionally evolve to generate dust disks that are visible in thermal emission and scattered-light. The Gemini Planet Imager (GPI; Macintosh 2014) is one of the first high-contrast instruments equipped with the extreme adaptive optics that is specially designed for direct imaging and spectroscopy of exoplanets and debris disks. By resolving a debris disk, we can characterize the spatial distribution of dust grains, infer its dynamical history, and deduce the presence of unseen planets (see Wyatt 2008 for a review). Unlike most debris disks which are gas-depleted, carbon monoxide is detected in the HD 131835 system (Moór et al. 2015). Being a rare resolved debris disk with detected gas, HD 131835 serves as a unique target for studying the relationship between gas-dust physics and planetary science.

HD 131835 is an A2IV star (Houk 1982) in the Upper Centaurus Lupus (UCL) moving group (Rizzuto et al. 2011), a subgroup of the Sco-Cen association. HD 131835 is ∼15 Myr old (Mamajek et al. 2002; Pecaut et al. 2012), at a distance of 123${}_{-13}^{+16}$ pc (van Leeuwen 2007). Its IR emission was discovered by the Infrared Astronomical Satellite (Moór et al. 2006). Chen et al. (2012) presented MIPS observations and showed that it is one of only four UCL/Lower Centaurus Crux A-type stars with ${L}_{{\rm{IR}}}/{L}_{*}\gt {10}^{-3}$, comparable to β Pic.

Recently, Hung et al. (2015) resolved the debris disk around HD 131835 at 11.7 and 18.3 μm using Gemini/T-ReCS. A three-component dust disk model, composed of an unusually warm continuous disk and two rings, was able to simultaneously explain the spectral energy distribution (SED) and the mid-IR thermal images. Compared to recent disk studies with GPI, such as HD 115600 (Currie et al. 2015) and HD 106906 (Kalas et al. 2015), HD 131835 is less inclined and the disk flux is more radially extended (∼35 to ∼400 AU), allowing us to better study its morphological features. Here we report the first scattered-light detection of dust surrounding HD 131835 in polarized light with GPI.

2.1. Observations

2.2. Data Reduction

2.3. Photometric Calibration

We perform photometric calibration on the polarimetry data by considering the satellite spot to star flux ratio R, the stellar flux ${F}_{\star }$ in physical units, and the average satellite spot flux S in analog-to-digital unit (ADU) per coadd. The calibrated image data D_f can be found using

$$\begin{eqnarray}&&{D}_{f}={D}_{i}\displaystyle \frac{R\;{F}_{\star }}{S},\end{eqnarray} \tag{ 1 }$$

where D_i is the image data in ADU coadd⁻¹. In the GPI H-band, $R\sim 2\times {10}^{-4}$ (Wang et al. 2014). We adopt ${F}_{\star }=965\pm 35$ mJy as the H-band flux of HD 131835 from 2MASS (Cutri et al. 2003). We use an elongated aperture, similar to the shape of a running track, to perform aperture photometry on the satellite spots in polarimetry mode. We obtain a conversion factor of $1\ \mathrm{ADU}\;{\mathrm{coadd}}^{-1}\;{{\rm{s}}}^{-1}=(7.8\pm 1.3)\times {10}^{-4}\;{\rm{mJy}}$, with the uncertainty mostly stemming from measurement of S.

The photometric calibration in the spectral mode (Maire et al. 2014) is done using the Calibrated Datacube Extraction recipe²³ via the GPI DRP based on the same principle. However, instead of using the broadband ${F}_{\star }$ and S fluxes in the above equation, we replace them with the host star spectrum and the average satellite spot spectrum. For the stellar spectrum, we use the IDL Astrolib routine ccm_unred.pro (based on Cardelli et al. 1989) to apply A(V) = 0.187 mag. reddening (Chen et al. 2012) to the Kurucz (1993) stellar atmosphere model with ${T}_{{\rm{eff}}}$ = 8770 K, $\mathrm{log}g$ = 4.0, and solar metallicity.

We resolve the debris disk around HD 131835 through polarimetric differential imaging. Figure 1 shows the calibrated GPI H-band polarized intensity of HD 131835 in radial Stokes ${Q}_{{\rm{r}}}$ (Schmid et al. 2006). In the sign convention adopted here, positive ${Q}_{{\rm{r}}}$ shows the tangentially polarized intensity, while negative ${Q}_{{\rm{r}}}$ represents radially polarized intensity. The disk appears to be inclined, with scattered-light extending from ∼75–120 AU. By fitting the location of the flux peak along the major axis on each wing, we find no significant offset that is larger than 300 mas. If we assume an axisymmetric and radially smooth density structure, the projected eccentricity, e, along the major axis is consistent with zero, and $e\gt 0.2$ is rejected at 1σ. Due to the non-detection on the southeast (SE) quadrant, the eccentricity along the minor axis is left unconstrained. The data are limited by instrumental polarization within the central $\sim 0\buildrel{\prime\prime}\over{.} 3$ and by photon noise at larger radii.

Zoom In Zoom Out Reset image size

Figure 1 shows brightness asymmetries along both the minor and major axes. The northwest (NW) side of the disk is significantly brighter than the SE side, which is undetected. This brightness asymmetry is likely due to light scattered in a preferential direction, as seen in the case of HR 4796A (Perrin et al. 2015). In addition, a weaker brightness asymmetry is present along the major axis. To quantify this brightness asymmetry, we use the best-fit geometric parameters found in Section 5 (except for setting the outer radius to be 180 AU due to the limited FOV) and consider the region exterior to 03 on the NW side of the major axis. Since single scattering by circumstellar dust is expected to produce linearly polarized light only in the ${Q}_{{\rm{r}}}$ polarization states, we measure the noise using Stokes ${U}_{{\rm{r}}}$, which corresponds to the linear polarization 45° from ${Q}_{{\rm{r}}}$. The error at each angular separation in the ${Q}_{{\rm{r}}}$ image was estimated by measuring the standard deviation of the three-pixel wide annulus at the same separation in the Stokes ${U}_{{\rm{r}}}$ image. We find that the northeast (NE) side of the disk is 1.30 ± 0.09 times brighter than the southwest (SW) side. In contrast, the thermal imaging shows that the sides are equally bright, with a 30% brightness asymmetry excluded at $\gt 3\sigma$ (Hung et al. 2015).

4.1. Disk Total Intensity

4.2. Exoplanet Search

We process our spectroscopic mode data to optimize planet detection. We first subtract the PSF for each wavelength channel using the TLOCI code (Marois et al. 2014), which assumes an input spectrum to optimize the subtraction while maximizing the signal-to-noise ratio (S/N) of an exoplanet of that specific spectral type. For our analysis, T8 and L8 spectra are chosen as the priors based on our experience in order to cover a wide range of DUSTY (Chabrier et al. 2000) and COND (Baraffe et al. 2003) exoplanets. The final data cube is then collapsed by a weighted mean, considering the input spectrum and the noise.

We searched for point sources with planet-like spectra but detected none. To estimate our upper limits, we derived the TLOCI contrast curves by measuring the standard deviation of the pixel noise in each annulus of $\lambda /D$ width. These contrast curves are then transformed into exoplanet mass upper limits using the BT-settl models (Allard et al. 2012). In our polarimetry mode observation, we derive the point-source contrast curves by dividing the scatter (due to photon and read noise) at each annulus by the stellar flux, similar to how the contrast curves are derived in the spectroscopic mode. The contrast curves for the total and polarized intensities and mass limits derived from the spectral-mode observations are shown in Figure 2. We reject objects with $M\gtrsim 3.5\ {M}_{{\rm{J}}}$ outside of 05.

Zoom In Zoom Out Reset image size

We take a two-step approach to find a model that fits the SED and the GPI image. First, we use a geometric model to retrieve the structure of the scattered-light disk. Then, fixing the disk geometry, we search for a physical model that is built on the model proposed by Hung et al. (2015) to get an estimate of the main dust properties associated with the polarized scattered-light detection. In both steps, we exclude the central region within 03 due to uncorrected systematic errors from the instrumental polarization and cut out the region beyond 260 AU in the disk plane to reduce the number of pixels without a detected disk signal.

To measure the basic geometric properties of the disk, we adopt a simple two-dimensional continuous disk model. Since the disk is not detected at all azimuths, we cut off 140° symmetrically about the SE semiminor axis in our model (white dashed lines in Figure 3) to exclude the region with S/N $\leqslant 1$. The model extends from the inner radius ${r}_{{\rm{in}}}$ to the outer radius ${r}_{{\rm{out}}}$, with the surface brightness varying only as ${r}^{\alpha }$. Along with the position angle PA of the major axis and the inclination i, we use these five parameters to describe the geometric properties of the disk. We fit the data using the emcee python package (Foreman-Mackey et al. 2013) based on the ensemble MCMC method of Goodman & Weare (2010). After a burn-in period, we let the 100 walkers run for 1200 steps. The best-fit parameters and uncertainties listed in Table 1 are found by taking the median of the marginalized probability density distributions and finding the 1σ confidence intervals. The shapes of the posterior distributions are all single-peaked and approximately normal. Assuming that the disk is optically thin and the polarization fraction and the phase function do not depend on the stellocentric radius, the value of α implies a nearly flat surface density profile of ${\rm{\Sigma }}(r)\propto {r}^{-0.3}$.

Zoom In Zoom Out Reset image size

Table 1.  Best-fit Model Parameters

Parameter	Value	Units
Geometric Parameters
PA	61.4 ± 0.4	°
i	${75.1}_{-0.9}^{+0.8}$	°
${r}_{{\rm{in}}}$	${75}_{-4}^{+2}$	AU
${r}_{{\rm{out}}}$	210 ± 10	AU
α	$-{2.3}_{-0.1}^{+0.2}$	⋯
Hotter MCFOST Model Component
PAa	61.4	°
ia	75.1	°
ra	35–400	AU
${r}_{{\rm{slope}}}{}^{a}$	0.5	⋯
aa	0.03–5	μm
${a}_{{\rm{slope}}}{}^{a}$	−3.5
${M}_{{\rm{dust}}}{}^{a}$	$6.66\times {10}^{-4}$	M⊕
${\rm{Composition}}{}^{a}$	amorphous carbon	⋯
Cooler MCFOST Model Component
PAa	61.4	°
ia	75.1	°
ra	75–210	AU
${r}_{{\rm{slope}}}{}^{a}$	−0.3
${a}_{{\rm{min}}}$	1.56	μm
${a}_{{\rm{max}}}{}^{a}$	1.0	mm
${a}_{{\rm{slope}}}$	−3.46
${M}_{{\rm{dust}}}$	$2.66\times {10}^{-1}$	M⊕
${\rm{Composition}}{}^{a}$	50% carbon + 50% silicate	...

Note.

^aKept fixed.

Download table as:  ASCIITypeset image

Next, we use the radiative transfer code MCFOST (Pinte et al. 2006, 2009) to generate the polarized scattered-light image and the SED. Using Mie theory, MCFOST self-consistently computes the absorption and scattering cross-sections as well as the scattering angle-dependent Mueller matrix, producing model images for all Stokes parameters. We assume a geometrically flat disk and start by modifying the two-component model from Hung et al. (2015) to match the GPI image and the SED. We set the hotter extended component to have amorphous carbon (Li & Greenberg 1997), as its composition is based on the suggested high grain temperatures (Hung et al. 2015). We keep all its parameters (PA, i, disk extent r, power-law index of the surface density ${r}_{{\rm{slope}}}$, grain size a, grain size distribution power-law index ${a}_{{\rm{slope}}}$, and mass of the dust ${M}_{{\rm{dust}}}$) fixed as listed in Table 1. The cooler ring component has its geometry fixed to be the values found in the previous paragraph and its composition set to be 50% amorphous carbon and 50% astro-silicate (Draine & Lee 1984). We have considered arguably simpler compositions (pure silicates and pure amorphous carbon) but they produce worse model fits. Adding water ice to the composition or porosity to the grains also leads to poorer model fits.

We simultaneously fit our MCFOST model to the SED and the scattered-light disk at all azimuthal angles. In the fit, we weight the residuals from each pixel and broadband photometry point equally. In addition to fitting photometry points longward of 10 μm (summarized in Hung et al. 2015), we include six new Herschel points from Moór et al. (2015). With the disk geometry and dust composition fixed, the only free parameters are the size distribution and ${M}_{{\rm{dust}}}$. The SED can only place a lower limit on the maximum grain size so we adopt ${a}_{{\rm{max}}}=1\;{\rm{mm}}$. This only leaves the minimum grain size ${a}_{{\rm{min}}}$, ${a}_{{\rm{slope}}}$, and ${M}_{{\rm{disk}}}$ of the cooler component to vary. To explore the parameter space, we use the genetic algorithm (Mathews et al. 2013), which ensures a fast convergence.

The best-fit MCFOST model reproduces the SED and the observed scattered-light disk with the reduced ${\chi }^{2}$ of 0.97. The fit is shown in Figure 3 and the parameters are listed in Table 1. The best-fit model image fits the disk ansae well. The slight over-subtraction near the inner edge is not significant given the noise in these regions. Although the hotter component does not contribute significantly to the scattered-light image, this component is an important source for thermal emission. Our best-fit MCFOST model also roughly reproduces the extended thermal emission. We can further improve our model with detailed analysis, such as more complex dust composition and grain size distribution, but those are beyond the scope of this letter. Nonetheless, we set strong constraints on the disk properties, reproducing the surface brightness of the scattered-light disk with a model that was initially devised exclusively on thermal emission. We find a minimum grain size that is in reasonable agreement with the expected blowout size of 0.91 μm and the grain size power-law index is only slightly steeper than the canonical ${a}_{{\rm{slope}}}=-3.5$. Overall, the grains detected with GPI seem to follow the intuitive expectations of common disk properties.

The mid-IR (Hung et al. 2015) and scattered-light images show different morphology. Unlike the continuous and extended thermal emission, the disk in scattered light has a cleared region inward of ∼75 AU. In other words, the scattered-light disk starts at a radius that is twice as far away from the star compared to the disk in thermal emission. In polarized scattered light, we detect brightness asymmetries strongly along the minor axis and weakly along the major axis. The brightness asymmetry along the minor axis is likely due to asymmetry in the scattering phase function and is present in our best-fit model. The brightness asymmetry along the major axis could be the result of a dust density enhancement, azimuthal variation of grain compositions, or a projection effect of an eccentric disk if this 3σ feature is real. Since the mid-IR data do not show the asymmetry along the major axis, it suggests that the large grain population is more symmetrically arranged than the small grain population. The simulation done by Wyatt (2006) shows that dust that originates in the break-up of planetesimals trapped in resonance with a planet can have moderate-sized grains (a few μm to a few mm) distributed axisymmetrically but small grains (less than a few μm) exhibit trailing spiral structure that emanates from the resonant clumps. Therefore, the mismatched distributions of large and small grains can identify different forces acting on them and highlight potentially interesting dynamical interactions in the system.

HD 131835 is distinctive compared with the other Sco-Cen debris disks that have recently been imaged in scattered light: HIP 79977, HD 115600, HD 106906, and HD 110058 (Thalmann et al. 2013; Currie et al. 2015; Kalas et al. 2015; Kasper et al. 2015). Among those, HD 131835 has the largest inner radius in the scattered-light-detected component and the most radially extended and nearly flat surface density profile (from 75 to 210 AU with ${\rm{\Delta }}r/r\sim 1$). The other disks either have relatively narrow belts (HD 115600, HD 106906, and HD 110058) or have relatively sharp declines in brightness (HIP 79977). The extended and approximately flat surface density profile suggests that the parent body belt of HD 131835 is likely extremely broad, much more so than any other debris disks imaged to date. This novel feature indicates that the silicate component is not distributed in the form of one or two narrow rings as previously suggested by Hung et al. (2015). In addition, among all the Sco-Cen disks, HD 131835 is the only resolved disk with detected CO gas (Moór et al. 2015), making it a unique and valuable target for studying gas-dust interactions. Besides the potential dynamical influence of undetected exoplanets, interactions between dust and gas could also play a significant role in clearing the dust in the inner disk and creating an eccentric ring (Lyra & Kuchner 2012).

Follow-up observations and detailed modeling are required to characterize the disk in detail. Deeper polarimetry observations are needed to confirm the NE–SW asymmetry. Detection of the disk in total intensity can set a firm constraint on the grain shape and porosity by providing the information on the fractional polarization as a function of scattering angle. Multicolor observations can further constrain the grain composition. Since HD 131835 is located in the southern sky, GPI, SPHERE (Spectro-polarimetric High-contrast Exoplanet Research), and ALMA (Atacama Large Millimeter/submillimeter Array) are powerful enough instruments/facilities for conducting follow-up observations.

This research was supported in part by NASA cooperative agreements NNX15AD95G, NNX11AD21G, and NNX14AJ80G, NSF AST-113718, AST-0909188, AST-1411868, AST-1413718, and DE-AC52-07NA27344, and the U.S. Department of Energy by Lawrence Livermore National Laboratory under contract DE-AC52-07NA27344. Work by L.-W. Hung and A. Greenbaum is supported by the National Science Foundation Graduate Research Fellowships DGE-1144087 and DGE-1232825. We acknowledge the Service Commun de Calcul Intensif de l'Observatoire de Grenoble (SCCI) for computations on the supercomputer funded by ANR (contracts ANR-07-BLAN-0221, ANR-2010-JCJC-0504-01 and ANR-2010-JCJC-0501-01) and the European Commission's 7th Framework Program (contract PERG06- GA-2009-256513). This work is based on observations obtained at the Gemini Observatory, which is operated by the Association of Universities for Research in Astronomy, Inc., under a cooperative agreement with the NSF on behalf of the Gemini partnership: the National Science Foundation (United States), the National Research Council (Canada), CONICYT (Chile), the Australian Research Council (Australia), Ministério da Ciência, Tecnologia e Inovação (Brazil), and Ministerio de Ciencia, Tecnología e Innovación Productiva (Argentina).