# Constraining the properties of transitional discs in Chamaeleon I with Herschel

Á. Ribas, H. Bouy, B. Merín, G. Duchêne, I. Rebollido, C. Espaillat, C. Pinte

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

## To cite this version:

Á. Ribas, H. Bouy, B. Merín, G. Duchêne, I. Rebollido, et al., Constraining the properties of transitional discs in Chamaeleon I with Herschel. Monthly Notices of the Royal Astronomical Society, 2016, 458, pp.1029-1040. 10.1093/mnras/stw333. insu-03634258

HAL Id: insu-03634258

https://insu.hal.science/insu-03634258

Submitted on 7 Apr 2022

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

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

Distributed under a Creative Commons Attribution 4.0 International License
# Constraining the properties of transitional discs in Chamaeleon I with Herschel*

Á. Ribas,<fnref n="1" />† H. Bouy,<fnref n="2" />† B. Merín,<fnref n="3" /> G. Duchêne,<fnref n="4" />,<fnref n="5" /> I. Rebollido,<fnref n="6" /> C. Espaillat<fnref n="1" /> and C. Pinte<fnref n="5" />,<fnref n="7" />

Accepted 2016 February 9. Received 2016 February 8; in original form 2015 November 11

## ABSTRACT

Transitional discs are protoplanetary discs with opacity gaps/cavities in their dust distribution, a feature that may be linked to planet formation. We perform Bayesian modelling of the three transitional discs SZ Cha, CS Cha, and T25 including photometry from the *Herschel Space Observatory* to quantify the improvements added by these new data. We find disc dust masses between $2 \times 10^{-5}$ and $4 \times 10^{-4}$ M$_{\odot}$ and gap radii in the range of 7–18 au, with uncertainties of ~ one order of magnitude and ~4 au, respectively. Our results show that adding *Herschel* data can significantly improve these estimates with respect to mid-infrared data alone, which have roughly twice as large uncertainties on both disc mass and gap radius. We also find weak evidence for different density profiles with respect to full discs. These results open exciting new possibilities to study the distribution of disc masses for large samples of discs.

**Key words:** planets and satellites: formation–planet–disc interactions – protoplanetary discs – stars: pre-main-sequence – infrared: planetary systems.

## 1 INTRODUCTION

Transitional discs (TDs) are one of the main research topics within the current paradigm of planet formation: these are protoplanetary discs with signatures of cavities and/or gaps in their dust distribution, which could be directly linked to forming planets (see Espaillat et al. 2014, for an updated review on this field). For this reason, a proper characterization of these systems could set strong constrains on the conditions under which planets come to be. However, these holes could also be produced by other mechanisms such as photoevaporation, dust growth and settling towards the disc mid-plane, dynamic clearing by (sub)stellar companions, or a combination of several of these (e.g. Ireland & Kraus 2008; Birnstiel, Andrews & Ercolano 2012; Alexander et al. 2014; Espaillat et al. 2014). Despite the large number of studies of TDs covering the whole wavelength domain with various observing techniques (e.g. photometry, spectroscopy, polarimetry, or interferometry), several questions remain open, such as their real connection with planet formation, whether

every protoplanetary disc goes through a transitional phase, or the main processes behind their evolution.

The Spectral Energy Distributions (SEDs) of full, optically thick circumstellar discs (Class II) have significant infrared (IR) excess with respect to the photospheric level from the near-infrared (~ NIR, 1–5 µm) to millimeter wavelengths (Williams & Cieza 2011). In contrast, SEDs of TDs normally display a very distinctive shape, with small or no NIR excess, but with larger excess at mid-infrared (~ MIR, 5–50 µm) and far-infrared (FIR, 50 µm and longer) wavelengths. This lack of NIR emission is usually attributed to dust depleted inner regions in the disc, and the location and shape of this change in the SED (around ~10 µm) can be used to partially characterize the gap structure (e.g. Espaillat et al. 2010, 2011). Objects with small NIR excess and significant MIR and FIR excesses are known as pre-transitional discs (pre-TDs; Espaillat et al. 2007b) and are thought to be an intermediate stage between full discs and TDs: these have some optically thick dust in their inner region, separated from the outer disc by an optically thin gap (see e.g. Espaillat et al. 2014). Given their possible connection with planet formation, (pre)TDs have been extensively studied and modelled in the past with MIR spectra from the IRS (Infrared Spectrograph) instrument (Houck et al. 2004) on board the *Spitzer Space Telescope*, which yielded thousands of IR spectra of circumstellar discs between 5 and 38 µm and provided detailed information about their dust

<footnote n="1">¹Department of Astronomy, Boston University, 725 Commonwealth Avenue, Boston, MA 02215, USA</footnote>

<footnote n="2">²Centro de Astrobiología, INTA-CSIC, PO Box - Apdo. de correos 78, E-28691 Villanueva de la Cañada, Madrid, Spain</footnote>

<footnote n="3">³European Space Astronomy Centre, ESA, PO Box, 78, E-28691 Villanueva de la Cañada, Madrid, Spain</footnote>

<footnote n="4">⁴Astronomy Department, University of California, Berkeley, CA 94720-3411, USA</footnote>

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

<footnote n="6">⁶Departamento de Física Teórica Módulo 15, Universidad Autónoma de Madrid, E-28049 Madrid, Spain</footnote>

<footnote n="7">⁷UMI-FCA, CNRS/INSU France (UMI 3386), and Departamento de Astronomía, Universidad de Chile, Casilla 36-D Santiago, Chile</footnote>

<footnote>*Herschel is an ESA space observatory with science instruments provided by European-led Principal Investigator consortia and with important participation from NASA.</footnote>

<footnote>†E-mail: aribas@bu.edu (AR); hbouy@cab.inta-csic.es (HB)</footnote>
composition and the morphology of the inner regions (e.g. Kim et al. 2009; Merín et al. 2010; Espaillat et al. 2011; Furlan et al. 2011). However, several parameters such as the disc flaring or mass remain poorly or completely unconstrained with SED modelling and MIR data alone. A more in-depth knowledge of TDs can be achieved with complementary (sub)mm data, where the disc becomes mostly optically thin and the flux can be related to the disc mass (e.g. Andrews & Williams 2005; Najita, Strom & Muzerolle 2007; Andrews et al. 2013). The advent of the Herschel Space Observatory (Pilbratt et al. 2010) produced a large number of FIR observations of protoplanetary and TDs in the 70–500 µm regime, providing us with information at larger wavelengths that can be used to better constrain some parameters of TDs.

In this paper, we add Herschel data to the modelling of three (pre)TDs: SZ Cha (a pre-TDs), CS Cha (a binary system surrounded by a disc with a large cavity), and T25 (a TD with no known companion). These objects belong to the Chamaeleon I region, located at 160–180 pc from the Sun (Whittet et al. 1997, and references therein), and with an age estimate of ~2 Myr (Luhman et al. 2008). Given its proximity, Chamaeleon I has been an usual target for star formation and stellar population studies (e.g. Luhman 2007; Luhman et al. 2008; Belloche et al. 2011), which identified more than 200 Young Stellar Objects (YSOs) in the region. It was also one of the clouds observed by Herschel, and some studies have already provided a Herschel view of its YSO population (Winston et al. 2009), discs around low-mass stars (Olofsson et al. 2013), and TDs (Ribas et al. 2013; Rodgers-Lee et al. 2014). Here, we focus on the impact of Herschel data in different parameters of these discs obtained from SED modelling, by combining radiative transfer modelling with Markov Chain Monte Carlo (MCMC) methods to perform a Bayesian analysis of their properties. Section 2 describes the sample and data used. The modelling procedure can be found in Section 3, and the results of the process are described in Section 4. Finally, we discuss the implications of our analysis in Section 5.

# 2 OBSERVATIONS AND DATA REDUCTION

## 2.1 The sample

A total of 12 sources have been previously classified as (pre)TDs in the Chamaeleon I region. Several of these targets have been modelled in detail, mostly based on their Spitzer MIR spectra (e.g. Kim et al. 2009). Additionally, a number of studies have already explored Herschel data of these discs (e.g. Winston et al. 2012; Ribas et al. 2013; Rodgers-Lee et al. 2014). To analyse the feasibility of estimating disc masses using Herschel data of (pre)TDs, we selected three different objects in Chamaeleon I:

(i) SZ Cha, a pre-TD (Kim et al. 2009),

(ii) CS Cha, a disc with a gap surrounding a binary system (binary separation of ~3.5 au, Nagel et al. 2012), i.e. a circumbinary disc (Guenther et al. 2007; Nagel et al. 2012),

(iii) T25, a TD with no known companion (Kim et al. 2009).

These three objects were selected because they represent the main scenarios in our current understanding of disc evolution: from a pre-TD (SZ Cha), to objects with clean opacity holes either caused by binaries (CS Cha) or other mechanisms (T25). For these targets, we used the stellar parameters in Espaillat et al. (2011), which are listed in Table 1. In the case of CS Cha, the orbital motion of the binary system may change the radiation received by different regions of the disc with time, and hence the emission from the inner disc varies with a similar period. Nagel et al. (2012) used two-star models to show that the variability produced by this effect in this case is only of ~1 per cent at the 10 µm peak, and we approximate the binary system by a single star (our photometric uncertainties and model noise are larger than this value, see Sections 2.3 and 2.4). We therefore adopt the spectral type provided in Luhman (2004) for this target, which has also been used in previous modelling efforts (Espaillat et al. 2007a; Kim et al. 2009; Espaillat et al. 2011; Manoj et al. 2011), allowing for meaningful comparisons.

## 2.2 Herschel data

In Ribas et al. (2013), we presented Herschel aperture photometry measurements of the TDs in Chamaeleon I. Given the inherent difficulties in determining whether a source is detected or not in the presence of conspicuous background, in that previous study we visually inspected the position of the known YSOs in the region (see Ribas et al. 2013, for a complete description of the sample). In this paper, we maintain this criterion, but also expand the analysis of T25 which was previously undetected at 500 µm (see below). Later, Rodgers-Lee et al. (2014) identified a systematic discrepancy between the photometry in Ribas et al. (2013) and the one in Winston et al. (2012), which used the *getsources* algorithm (Men'shchikov et al. 2012). We attribute such discrepancies to the different map-making and source extraction algorithms used, but for the sake of completeness we chose to re-process the corresponding data. Here, we describe the adopted methods.

The Chamaeleon I region was observed by Herschel as part of the Gould Belt Key Programme (André et al. 2010). Parallel mode observations from this programme (OBSIDs: 1342213178, 1342213179) provided PACS (Poglitsch et al. 2010) 70 and 160 µm and SPIRE (Griffin et al. 2010) 250, 350, and 500 µm maps at a scan speed of 60 arcsec s⁻¹. Additional PACS 100 and 160 scan observations are also available in the Gould Belt Key Programme at a scan speed of 20 arcsec s⁻¹ (OBSIDs: 1342224782, 1342224783). Although with smaller coverage, this last set of observations is deeper and has a slower scan speed (hence a smaller point spread function), so we chose to use them for the 160 µm band instead of the parallel mode data. We note that T25 is outside the smaller scan maps, and no 100 µm photometric measurement is available for this target. Its 160 µm flux was therefore obtained from the parallel mode data.

We used the Herschel Interactive Processing Environment (HIPE; Ott 2010) version 12.1 to process the maps of the region. We adopted the standard map-making algorithms used in the Herschel Science

<table><caption>**Table 1.** Coordinates and stellar parameters used in this work for the considered sample of (pre)TDs. Stellar parameters as in Espaillat et al. (2011).</caption><thead><tr><th>Name</th><th>RA<sub>J2000</sub></th><th>Dec.<sub>J2000</sub></th><th>A<sub>V</sub><br>(mag)</th><th>SpT</th><th>T<sub>*</sub><br>(K)</th><th>L<sub>*</sub><br>(L<sub>☉</sub>)</th><th>M<sub>*</sub><br>(M<sub>☉</sub>)</th><th>R<sub>*</sub><br>(R<sub>☉</sub>)</th></tr></thead><tbody><tr><td>SZ Cha</td><td>10:58:16.77</td><td>-77:17:17.1</td><td>1.9</td><td>K0</td><td>5250</td><td>1.9</td><td>1.4</td><td>1.7</td></tr><tr><td>CS Cha</td><td>11:02:24.91</td><td>-77:33:35.7</td><td>0.8</td><td>K6</td><td>4205</td><td>1.5</td><td>0.9</td><td>2.3</td></tr><tr><td>T25</td><td>11:07:19.15</td><td>-76:03:04.9</td><td>1.6</td><td>M3</td><td>3470</td><td>0.3</td><td>0.3</td><td>1.5</td></tr></tbody></table>
<table><caption>**Table 2.** Herschel fluxes of the modelled (pre)TDs in this study. Ellipsis indicate non-detected sources.</caption><thead><tr><th>Name</th><th>F70<br>(Jy)</th><th>F100<br>(Jy)</th><th>F160<br>(Jy)</th><th>F250<br>(Jy)</th><th>F350<br>(Jy)</th><th>F500<br>(Jy)</th></tr></thead><tbody><tr><td>SZ Cha</td><td>4.01 ± 0.80</td><td>3.74 ± 0.75</td><td>3.56 ± 0.71</td><td>2.53 ± 0.51</td><td>1.85 ± 0.37</td><td>1.02 ± 0.20</td></tr><tr><td>CS Cha</td><td>3.20 ± 0.64</td><td>2.88 ± 0.58</td><td>2.27 ± 0.45</td><td>1.31 ± 0.26</td><td>1.04 ± 0.21</td><td>...</td></tr><tr><td>T25</td><td>0.53 ± 0.11</td><td>...</td><td>0.38 ± 0.08</td><td>0.25 ± 0.05</td><td>0.17 ± 0.03</td><td>0.06 ± 0.01</td></tr></tbody></table>

Archive, i.e. *jscanam* for PACS maps (a HIPE adaptation of the *SCANAMORPHOS* software; Roussel 2013), and *destriper* for SPIRE maps. In the case of *jscanam*, we remove turnarounds with speeds 50 per cent lower or higher than the nominal speed value, and we do not use the extended emission gain option for *destriper*, as recommended for point source photometry.

We estimated PACS fluxes with the *AnnularSkyAperturePhotometry* task in HIPE. We adopted aperture radii of 15, 15, and 22 arcsec for PACS 70, 100, and 160 bands, respectively. These apertures were specifically selected after inspecting growth curves of each target in each band. The background was estimated within annulus with inner and outer radii of 25 arcsec and 35 arcsec. We applied the corresponding aperture corrections of 1.206, 1.222, and 1.372 (Balog et al. 2014). Given the negligible colour corrections for PACS for temperatures above 20 K (PACS Photometer - Colour Corrections manual, version 1.0) and the uncertainty in determining the slope of the SED close to the emission peak, we chose not to apply them for PACS bands since they are significantly smaller than the adopted photometric uncertainties (see below). For SPIRE, we used the recommended procedure and fit the sources in the timeline (Pearson et al. 2014). This method does not require aperture corrections. T25 was considered as an upper limit at 500 µm in Ribas et al. (2013), but the procedure in this study successfully detected it in this band. The obtained flux value does not conflict with the previous upper limit in Ribas et al. (2013). Given the better method used here and the flux consistency with the overall shape of the SED, we chose to include it. Conversely, the timeline filter does not detect CS Cha at 500 µm probably due to the strong background, and therefore we do not include this wavelength in its SED. We applied colour correction factors corresponding to black-body emission assuming SPIRE bands trace the Rayleigh-Jeans regime (0.945, 0.948, 0.943 for SPIRE 250, 350, and 500 bands, respectively, see SPIRE Handbook Version 2.5). Finally, we adopted conservative photometric uncertainties of 20 per cent to account for different effects (i.e. absolute flux calibration, background estimation). The obtained *Herschel* photometry is provided in Table 2.

## 2.3 NIR/MIR photometric data of TDs

In Ribas et al. (2013), we compiled photometry from several surveys and catalogues to build well-sampled SEDs of the TDs in the region. However, the aim of this paper is to model these SEDs, and hence we only select non-redundant photometric data. We therefore chose to include the following bands in the NIR/MIR: 2MASS *J*, *H*, and *Ks*, and IRAC1 and IRAC2 bands. This selection provides a nice coverage of the 1–6 µm regime, key to separate (pre)TDs from TDs (Espaillat et al. 2010). It also avoids redundancy (including several data in a small wavelength domain), which could give excessive weight to certain parts of the SED in the model fitting process. Typical photometric uncertainties of these measurements are below 5 per cent, but given the main scope of this paper, possible IR variability of the sources should also be considered to derive proper uncertainties in the physical parameters (Muzerolle et al. 2009). To account for these two effects (photometric uncertainties and

variability), we chose to set uncertainties to be a 10 per cent of the observed fluxes.

Finally, all photometric points were derredened using the corresponding $A_V$ (see Table 1) and the extinction law in Indebetouw et al. (2005).

## 2.4 IRS spectra of TDs

We retrieved low resolution IRS spectra from the Cornell Atlas of *Spitzer/IRS* Sources (CASSIS) data base (Lebouteiller et al. 2011) for the (pre)TDs discs in our study. CASSIS provides optimally extracted IRS spectra, and is well suited for our purposes. For these spectra, we first separated the optimal zones of the IRS spectra (7–14 µm and 20.5–35 µm for the first order, <20.5 µm for second one), and rejected bad pixels (e.g. negative or NaN values). As a compromise between estimating monochromatic fluxes required for model fitting (see Section 3) while reducing the impact of possible artefacts in the spectra, we chose to bin them in 10 equally spaced wavelengths throughout the spectra coverage, and estimate the fluxes for each of them as the mean value of 10 pixels centred around each corresponding wavelength. We checked this procedure to produce nice sampling of the IRS spectra (see Fig. 1 for an example), while being a good compromise for the SED fitting: a whole IRS spectra typically contains 300–400 good pixels, and fitting them all would put most of the weight on the IRS spectra itself. By reducing its contribution to a comparable number to that of photometric data (~10) we ensure that all parts in the SED contribute to the fitting procedure in a similar manner. Additionally, this binning choice allows to encapsulate the basics of the silicate feature (i.e. its presence and strength). As in Section 2.3, we assigned 10 per cent uncertainties to the binned data, a typical variability value for these discs (Espaillat et al. 2011).

<figure><img src="image_3.png" /><figcaption>**Figure 1.** Derredened IRS spectrum of SZ Cha. Black dots show the CASSIS spectrum of this source (after bad pixels rejection) with the corresponding error bars. The binned spectrum and assumed uncertainties are shown as red, larger dots. It properly traces the shape of original data including the silicate feature at 10 µm.</figcaption></figure>
## 3 MODELLING

We aimed at modelling the selected targets and quantifying the impact of adding photometric Herschel data to this process. For this reason, we used two different data sets for each object. The first data set comprises the available data from 2MASS, IRAC1/IRAC2, and the binned IRS spectra. The second data set also includes the Herschel photometry.

We used the MCFOST software (Pinte et al. 2006, 2009) version 2.19 to model these discs. MCFOST is a Monte Carlo-based raytracing code which generates synthetic SEDs and images of circumstellar discs. First, it produces temperature and density maps of the disc using the provided stellar and disc parameters. In this case, we used $10^7$ photons in this step (enough to produce smooth and well-sampled temperature maps of the discs). After this step, a list of wavelengths is provided for MCFOST to calculate the corresponding synthetic monochromatic fluxes. We required 2000 photons to be received for each wavelength, corresponding to noise levels of 2–3 per cent in the flux estimates and well below the assumed observational uncertainties (10–20 per cent).

Our models include seven free parameters: disc dust mass ($M_{\text{dust}}$), inner and outer radii ($R_{\text{in}}$, $R_{\text{out}}$), scaleheight at 100 au ($H_{100}$), flaring index ($h$), surface density exponent ($p$), and the maximum grain size ($a_{\text{max}}$). Given the complex structures of circumstellar discs, there are several degeneracies and dependencies between these parameters, and some may even be totally unconstrained with the available data. We did not attempt to fit the mineralogy of the discs. Instead, we assumed typical astronomical silicate compositions, and fixed the minimum grain size to 0.01 µm. A more in-depth study of the mineralogy of these discs would add an important source of complexity to modelling, and we preferred not to include it in our comparative analysis. We chose a power-law index for the surface density profile. More complex structures such as tapered-edge profiles could also be used (e.g. Lynden-Bell & Pringle 1974), but direct high-resolution observations are required to actually trace the mass distribution in the disc. Following Espaillat et al. (2011), we also fixed the inclination of all discs to 60°. Although this inclination is somewhat arbitrary, none of these objects show signatures of high inclination (e.g. silicate features in absorption or underluminous photospheres). Moreover, for wavelengths >13 µm, the MIR continuum is almost insensitive to this effect unless the disc is very close to edge-on (D’Alessio et al. 2006; Furlan et al. 2006).

Among the objects in the sample, SZ Cha has NIR excess over the photospheric emission. This feature is characteristic of pre-TDs, sources with an optically thick inner disc, separated from the outer disc by a gap in the radial dust distribution (e.g. Espaillat et al. 2007a). Additionally, CS Cha has no NIR excess but a prominent silicate emission feature at 10 µm, indicating the presence of optically thin dust in its inner hole. The inner discs of these objects have already been modelled in detail (e.g. Espaillat et al. 2007b; Kim et al. 2009; Manoj et al. 2011) and we do not attempt to fit them: instead, we adopted the parameter results from these previous studies to reproduce the NIR SED shape. The inner discs remained fixed during the fitting process. This may have an impact in our final results, as discussed later in the paper (see Section 3.2).

### 3.1 Methodology

We adopted a Bayesian approach to properly derive confidence intervals for the outer disc parameters. The usage of Bayesian techniques has increased significantly in Astrophysics during the past years, and we do not intend to explain them in detail. Instead, we

refer the interested reader to introductory works such as Trotta (2008). Also, this technique has already been applied for modelling circumstellar discs with Herschel data (e.g. Cieza et al. 2011; Harvey et al. 2012; Spezzi et al. 2013), mainly via model grids. Here, we describe the adopted fitting procedure.

Bayesian analysis requires that we assign priors to model parameters. While the selection of restrictive priors may have a significant effect on the fitting results, priors are also an important tool to force parameters to remain within certain ranges, avoiding non-physical solutions. We used flat (non-informative) priors for all the parameters, and constrain them to reasonable values for TDs. The prior ranges used were as follows:

(i) $\log(M_{\text{dust}}/M_{\odot})$: from -6 to -2,

(ii) $R_{\text{in}}$: from 1 au to $r_{\text{in-out}}$,

(iii) $R_{\text{out}}$: from $r_{\text{in-out}}$ to 500 au,

(iv) $H_{100}$: from 0.5 to 25 au,

(v) $h$: from 0.8 to 1.3,

(vi) $p$: from -2.5 to 1,

(vii) $\log(a_{\text{max}}/\mu\text{m})$: from -1 to 4,

where $r_{\text{in-out}}$ depends on the target, and is a physically meaningless parameter merely used to avoid the outer disc becoming smaller than the inner one during the evolution of the MCMC. Based on previous results (Kim et al. 2009; Espaillat et al. 2011), we set $r_{\text{in-out}}$ to 30, 40, and 50 au for T25, SZ Cha, and CS Cha, respectively. $M_{\text{dust}}$ and $a_{\text{max}}$ can take values within several orders of magnitude, and hence we chose to explore them in logarithmic scale.

We used a modified version of MCMC methods called ensemble samplers with affine invariance (Goodman & Weare 2010). This method uses several ‘walkers’ or individual chains to explore the posterior distributions of parameters, and is especially useful when these distributions have complex forms. We used a slightly modified version of the implementation by Foreman-Mackey et al. (2013), and set the stretch parameter of the walk to 1.5, getting acceptance ratios between 10 and 50 per cent (a good compromise between a random walk and discarding most of the proposed positions in the chain evolution). In every iteration, the chain comprises 100 walkers. We assume Gaussian uncertainties for our observations, and used the corresponding likelihood function.

To avoid dependencies with the distance to the objects, we normalize every model to the *J* band prior to estimating the likelihood. This should have no impact in our results, since the *J* band traces photospheric emission in TDs and does not depend on disc parameters (i.e. all models obtained for a given object have always the same *J* flux).

When available, we set the initial position of the chains around previous results in the literature (Kim et al. 2009; Espaillat et al. 2011). The posterior from MCMCs are only valid once the system has lost memory of their initial values. This can be quantified using the autocorrelation time of the chains, which gives an estimate of the required number of iterations to draw independent samples. For every case, we computed the autocorrelation time for each walker in each parameter, and took the maximum value for conservative purposes. We then left the system evolve for five autocorrelation times (typically ~500 iterations). At this point, the results are independent of the initial position, and the chain is now sampling the posterior distribution. We then estimated the posterior function with other five autocorrelation times (i.e. 50 000 models, the result of the 500 iterations per 100 walkers used).
## 3.2 Model caveats and limitations

Simple parametric modelling like the one used in this paper offers several advantages (e.g. we can compute synthetic SEDs of complex discs without analytic solution), but it also suffers from some caveats that should be kept in mind. Parametric modelling does not guarantee that a combination of parameters is physically consistent, which we have tried to attenuate using physically meaningful priors. We have not included more complex features/models, such as puffed up inner rims (e.g. Dullemond & Dominik 2004; Isella & Natta 2005), non-axisymmetric inhomogeneities (e.g. Andrews et al. 2011; van der Marel et al. 2013), several radial gaps (as ALMA observations have revealed for HL Tau, ALMA Partnership et al. 2015), or fitting for the inner discs (e.g. Espaillat et al. 2010, 2011) or mineralogies. None the less, some of these caveats are likely to have little or no impact in our final results, if applicable at all. The homogeneous treatment of data and fitting procedure used provides a good understanding of the value of each data set, and an adequate frame for comparing the obtained distributions.

# 4 RESULTS

## 4.1 Fitting results without Herschel data

We first use the data set without Herschel data to explore which parameters can be constrained with NIR/MIR information. The posterior distributions for the three targets are shown in Fig. 2, and the obtained values in Table 3. MCMCs also allow us to study degeneracies between different parameters by plotting the chains in different 2D projections. The degeneracies are very similar in all cases, as revealed by the cornerplots in the appendix (Figs A1 to A3).

The following conclusions can be drawn from the obtained posterior distributions for the model parameters by analysing 5–95 per cent confidence intervals (Fig. 3).

(i) The inner radius ($R_{in}$) can be constrained within 5–20 au for all three targets, corresponding to relative uncertainties between 80 and 120 per cent. This is likely due to the fact that most of the NIR/MIR emission arises in the exposed wall, hence probing the location of this parameter.

(ii) NIR/MIR photometry allows us to calculate the scaleheight ($H_{100}$) with uncertainties within 10 au. This is expected, since different scaleheights modify the amount of stellar flux intercepted by the disc, changing the emission from the inner disc.

(iii) The dust mass ($M_{dust}$) and the rest of the geometrical parameters ($R_{out}$, $h$, and $p$) show little or no constraint at all with NIR/MIR data alone.

(iv) Special attention should be paid to the two-peak distributions obtained for $a_{max}$. Any observation at a given wavelength $\lambda$ is only sensible to emission from grains of size $a \sim \lambda$ (Draine 2006). If we allow $a_{max}$ to take small enough values ($<10$ µm), it could produce substantial changes in the SED and therefore play a role in the modelling, which may explain the double peaked posterior distributions. Although much larger grains are generally expected in discs, this data set is not enough to resolve this effect if we allow $a_{max}$ to take small enough values.

(v) As expected, several degeneracies appear in all cases, the most obvious being inner radius with scaleheight, and scaleheight with flaring index. In some cases, the inner radius also shows a dependence with the maximum grain size, in relation with the previous point.

<figure><img src="image_4.png" /><figcaption>**Figure 2.** Posterior distributions of the free parameters for the considered TDs. Results without Herschel data are shown in blue, those including Herschel in red.</figcaption></figure>
<table><caption>**Table 3.** Results of the MCMC fitting for the seven free parameters considered. We tabulate the obtained median value and the 5 per cent and 95 per cent confidence intervals. For each column, results obtained without *Herschel* data (left) and with them (right) are provided.</caption><thead><tr><th>Object</th><th>log M<sub>dust</sub><br>(log M<sub>☉</sub>)<br>no Herschel | with Herschel</th><th>R<sub>in</sub><br>(au)<br>no Herschel | with Herschel</th><th>R<sub>out</sub><br>(au)<br>no Herschel | with Herschel</th><th>H<sub>100</sub><br>(au)<br>no Herschel | with Herschel</th></tr></thead><tbody><tr><td>SZ Cha</td><td>-3.6<sup>+1.5</sup><sub>-2.0</sub> | -3.4<sup>+1.0</sup><sub>-0.5</sub></td><td>17<sup>+7</sup><sub>-7</sub> | 13<sup>+2</sup><sub>-4</sub></td><td>220<sup>+240</sup><sub>-170</sub> | 200<sup>+240</sup><sub>-120</sub></td><td>8.8<sup>+5.0</sup><sub>-3.1</sub> | 7.6<sup>+2.9</sup><sub>-1.7</sub></td></tr><tr><td>CS Cha</td><td>-3.8<sup>+1.6</sup><sub>-2.0</sub> | -3.8<sup>+1.4</sup><sub>-0.8</sub></td><td>19<sup>+16</sup><sub>-7</sub> | 18<sup>+6</sup><sub>-5</sub></td><td>200<sup>+250</sup><sub>-130</sub> | 210<sup>+250</sup><sub>-140</sub></td><td>5.8<sup>+10.4</sup><sub>-2.7</sub> | 4.6<sup>+4.2</sup><sub>-1.3</sub></td></tr><tr><td>T25</td><td>-4.2<sup>+2.0</sup><sub>-1.6</sub> | -4.7<sup>+0.6</sup><sub>-0.5</sub></td><td>8.0<sup>+4.4</sup><sub>-3.4</sub> | 6.6<sup>+1.6</sup><sub>-2.1</sub></td><td>190<sup>+260</sup><sub>-150</sub> | 220<sup>+250</sup><sub>-160</sub></td><td>6.0<sup>+8.1</sup><sub>-3.0</sub> | 7.6<sup>+9.0</sup><sub>-3.8</sub></td></tr><tr><td>Object</td><td></td><td>h<br>no Herschel | with Herschel</td><td>p<br>no Herschel | with Herschel</td><td>log a<sub>max</sub><br>(log µm)<br>no Herschel | with Herschel</td></tr><tr><td>SZ Cha</td><td></td><td>1.1<sup>+0.2</sup><sub>-0.3</sub> | 1.0<sup>+0.2</sup><sub>-0.2</sub></td><td>-1.0<sup>+1.7</sup><sub>-1.3</sub> | -1.6<sup>+1.9</sup><sub>-0.8</sub></td><td>0.2<sup>+3.5</sup><sub>-0.6</sub> | 2.9<sup>+1.0</sup><sub>-1.1</sub></td></tr><tr><td>CS Cha</td><td></td><td>1.1<sup>+0.2</sup><sub>-0.3</sub> | 1.0<sup>+0.3</sup><sub>-0.2</sub></td><td>-1.1<sup>+1.8</sup><sub>-1.3</sub> | -1.6<sup>+1.9</sup><sub>-0.9</sub></td><td>2.4<sup>+1.5</sup><sub>-2.8</sub> | 2.7<sup>+1.1</sup><sub>-2.1</sub></td></tr><tr><td>T25</td><td></td><td>1.1<sup>+0.1</sup><sub>-0.3</sub> | 1.1<sup>+0.2</sup><sub>-0.2</sub></td><td>-1.0<sup>+1.8</sup><sub>-1.3</sub> | -1.1<sup>+1.1</sup><sub>-1.2</sub></td><td>1.8<sup>+2.0</sup><sub>-2.3</sub> | 2.9<sup>+1.0</sup><sub>-0.9</sub></td></tr></tbody></table>

<figure><img src="image_5.png" /><figcaption>**Figure 3.** Derredened SEDs for the TDs in this study. Photometric data from previous studies are shown as blue solid circles, *Herschel* measurements as orange squares. Uncertainties are plotted, although in several cases are smaller than symbol sizes. We also show 100 randomly selected models from the obtained posterior distributions for each case: blue lines correspond to fitting without *Herschel* data, red lines are the resulting models when including *Herschel* photometry. This gives an idea on the total uncertainties in the SEDs of the modelled discs.</figcaption></figure>

## 4.2 Fitting results with Herschel data

We repeated the modelling procedure including *Herschel* data for the three selected TDs. As in the previous case, the resulting posterior distributions are shown in Fig. 2, and full corner plots in the appendix (Figs A1 to A3). Fig. 3 shows the observed SEDs and modelling results, and Table 3 provides the obtained numerical values.

(i) Compared to the *Spitzer*-only fit, the addition of *Herschel* data makes an important difference for $M_{dust}$ and $R_{in}$. For the latter, the posterior distributions are narrowed down by a factor of 2 with respect to the previous case, with the 5–95 per cent confidence intervals covering 5–10 au, or relative errors of 45–60 per cent. For the dust mass, the improvement is substantial, constraining its value within one order of magnitude for SZ Cha and T25, and a broader distribution (~2 dex) for CS Cha, given that it lacks a detection at 500 µm.

(ii) The scaleheight ($H_{100}$) is better constrained with *Herschel* for SZ Cha and CS Cha, reducing the uncertainties by a factor of ~2. For T25, there is no additional improvement. We note that, despite this, the combination of $R_{in}$, $h$, and $H_{100}$ yield very similar values of the scaleheight at the inner radius, with improvements in uncertainty below of 1 au or less.

(iii) For the rest of disc geometry parameters, there is no real improvement compared to the *Spitzer*-only fit. We do however see marginal evidence of anomalous outer discs in these objects, when

combining the preferred values of $h$ and $p$, specially for SZ Cha and CS Cha. We will discuss this in the following section.

(iv) *Herschel* data break the two-peak degeneracy in the maximum grain size ($a_{max}$). Although they are not enough to provide a real estimate of this value, they inform that this value is very likely larger than 100 µm in all cases.

# 5 DISCUSSION

## 5.1 Masses and inner radii

The estimate of two free parameters in our models have been found to improve significantly when including *Herschel* data: the disc dust mass and its inner radius.

The mass of the disc is one of the most important parameters for planet formation. It determines the available reservoir to build up planets and for accretion on the central star, and can even modify the planet formation mechanism (the disc instability scenario requires high $M_{disc}/M_{*}$ values, e.g. Lodato, Delgado-Donate & Clarke 2005). Although the bulk of mass in protoplanetary discs is in gaseous form, photometric IR data are only sensitive to dust emission, which are efficient radiation absorbers and emitters. Hence, only the dust mass can be (partially) constrained with the presented data. A total disc mass estimate requires either gas mass measurements (e.g. $^{13}$CO, $C^{18}$O, Panić et al. 2008) or assumptions on the
gas-to-dust ratio. Since the former is only available for a few discs, most authors assume a typical gas-to-dust ratio of 100. We adopt this value for comparison with previous studies, but are aware that such an assumption may not hold true in many cases, since the gas-to-dust ratio is likely to change with time and from one source to another (Thi et al. 2010, 2014). Regardless of this, our results show that *Herschel* data can be used to constrain the mass of dust in these (pre)TDs within ~ one order of magnitude for the 5–95 per cent confidence interval, which is a tremendous improvement with respect to MIR photometry only and opens the exciting possibility of studying this parameter for the large number of sources observed with *Herschel* and missing mm observations. As expected and illustrated by the case of CS Cha, the longest wavelengths (SPIRE 500 µm) are the most important ones to constrain the mass, as the disc becomes progressively optically thin at longer wavelengths.

The inner radii of the discs are also significantly better constrained with *Herschel* data, decreasing previous uncertainties by a factor of 2. This improvement arises from PACS data, which narrow down the location of the peak emission, directly related with the illuminated inner wall. Using unbinned IRS spectra could also provide even better estimates for these values. Nevertheless, our results show that at least some information about this parameter is contained in *Herschel* data.

We also compare the obtained values with previous estimates in the literature. There are four main studies which can be used for this purpose: Kim et al. (2009), Espaillat et al. (2011), Ubach et al. (2012), and Rodgers-Lee et al. (2014). Kim et al. (2009) presented detailed modelling of the IRS spectra of TDs in Chamaeleon I with an analytic model, including emission from the optically thin disc and wall and emission from the outer disc treated as a blackbody. In Espaillat et al. (2011), the authors used a more complex irradiated disc model (D'Alessio et al. 2006) including shadowing of the outer disc by the inner disc (Espaillat et al. 2010) to analyse variability in the IRS spectra of several TDs. Ubach et al. (2012) presented 3 and 7 mm interferometric measurements of SZ Cha and CS Cha, providing disc mass estimates. More recently, Rodgers-Lee et al. (2014) performed a multiwavelength study of Chamaeleon I including *Herschel* data from Winston et al. (2012), and analysed TDs in the region with a physical disc model (Beckwith et al. 1990).

(i) Mass values computed with two different methods are available for all the targets: via mm data (Kim et al. 2009; Ubach et al. 2012; Rodgers-Lee et al. 2014) and via the accretion-to-viscosity ratio (Espaillat et al. 2011) following D'Alessio et al. (1998). Our results are in very good agreement with these previous values, and match within a factor of 2 for most cases except for the mass value of CS Cha in Espaillat et al. (2011), SZ Cha in Ubach et al. (2012), and T25 from Kim et al. (2009). In the former case, a value of 0.3 M<sub>☉</sub> is quoted, more than one order of magnitude larger than our estimated median value (0.015 M<sub>☉</sub>) but within the 95 per cent confidence interval range. Therefore, the two results are consistent within uncertainties. Additionally, the value in Espaillat et al. (2011) depends on the disc viscosity, which is usually largely uncertain and could account for this difference. In the case of SZ Cha, our results for the 5–95 per cent interval yields values of 1.3 × 10<sup>-2</sup>–0.4 M<sub>☉</sub>, while Ubach et al. (2012) obtained a total disc mass of 9.4 × 10<sup>-3</sup> M<sub>☉</sub>. Considering that this measurement is also subject to uncertainties (between a factor of 2–10, according to Ubach et al. 2012), then our results match completely within the uncertainty range. For T25, Kim et al. (2009) estimate a 0.007 M<sub>☉</sub> disc mass via 1.3 mm fluxes from Henning et al. (1993). Our study yields a disc mass for T25 of 0.002 M<sub>☉</sub>, with their value lying just at the border of the cor-

responding confidence interval. However, as noted by Rodgers-Lee et al. (2014), the 1.3 mm flux value in Kim et al. (2009) for T25 is an upper limits, and therefore their mass estimates should be considered as such, solving the discrepancy. Rodgers-Lee et al. (2014) also found that *Herschel* data within the 160–500 µm range can be used to estimate disc masses within a factor of 3 without the need of detailed modelling. Our results show larger uncertainties (~ one to two orders of magnitude for the 5–95 per cent confidence interval), stressing the importance of considering other sources of uncertainties (such as disc temperature and composition) to compute realistic confidence intervals of model parameters.

(ii) Disc inner radii estimates are available both in Kim et al. (2009) and Espaillat et al. (2011). Our results using MIR data only are in general good agreement with their values, with the discrepancy of CS Cha. These two works estimated its inner disc radii to be 41 and 38 au, respectively, while we obtain 19<sup>+16</sup>/<sub>-7</sub> au with similar data (i.e. excluding *Herschel*). Their results fall outside the 5–95 per cent confidence intervals derived in this study. Two different effects can explain this apparent discrepancy. First, there is no uncertainty estimation in the quoted studies: if we assume their results to have similar uncertainties to ours, the resulting distributions would overlap significantly and yield consistent values. Additionally, these two works included more complex dust compositions, which can modify the grain emissivity and therefore change the location of the inner radius. We also note that Kim et al. (2009) estimated a 29 au gap for T25, although the improved estimate of 18 au in Espaillat et al. (2011) is completely consistent with ours.

The mass ranges of these TDs are similar to those of Class II discs in other star-forming regions (e.g. Ophiuchus, Taurus, Andrews et al. 2010, 2013), a result already found by Andrews et al. (2011) for 12 TDs observed with submm interferometry. This is somehow intriguing: if TDs are an evolved stage of protoplanetary discs, then we would expect them to have significantly lower masses. In fact, other works found TDs to have masses even higher than those of Class II sources (Najita, Strom & Muzerolle 2007; Najita, Andrews & Muzerolle 2015). If that is the case, TDs (at least classical ones, those with large holes in their dust distribution) could be the evolution of high-mass discs which have formed multiple giant planets (explaining their gaps, Zhu et al. 2011), and not a general evolutionary stage for all protoplanetary discs.

We also compare these values with that of the Minimum Mass Solar Nebula (MMSN; Hayashi 1981), the minimum mass required to form the Solar system. A typical value of this quantity is ~0.02 M<sub>☉</sub> (Davis 2005; Desch 2007). Both SZ Cha and CS Cha are above or close to this value, meaning that despite being in a transitional stage, they still have enough mass to form a significant number of planets (although this does not guarantee that planet formation will take place in the future).

## 5.2 Anomalous outer discs

We find a general trend for flaring indexes (*h*) close to ~1, slightly smaller than those usually found in protoplanetary discs (~1.1–1.3, e.g. Chiang & Goldreich 1997; Olofsson et al. 2013). Additionally, *Herschel* data suggest strongly negative surface density profiles, with no peak at −1, as typically assumed and found in protoplanetary discs (e.g. Andrews et al. 2009). Surprisingly, the obtained values are closer to that of the estimated for the MMSN (i.e. −1.5, −2.2, Hayashi 1981; Desch 2007).

These two results are likely accounting for an observed trend in the SEDs of these three targets: they have a significant amount
of excess in the MIR range up to 70–100 µm (already hinted in Cieza et al. 2011; Ribas et al. 2013), but their slopes between 250 and 500 µm are bluer than those of typical Class II discs. This was found in Ribas et al. (2013) when comparing the SEDs of (pre)TDs in Chamaeleon I with the median SED of the Chamaeleon I and II regions. The obtained steep surface density profiles and flaring indexes reduce the flux at longer wavelengths (SPIRE), and increase it at shorter wavelengths (20–150 µm, IRS, MIPS, and PACS). Low flaring indexes could arise if significant dust settling towards the mid-plane has already occurred in these discs, reducing the disc surface exposed to the stellar radiation specially in the outer regions of the disc. On the other hand, smaller (more negative) surface densities imply that more mass is located close to the star, leaving a fainter outer disc which will emit poorly in the FIR regime. Combined, these results suggest that the modelled (pre)TDs have anomalous outer discs compared to Class II objects. The same fact is found for the T Cha TD using *Herschel* data (Cieza et al. 2011) and in the Lupus region (Bustamante et al. 2015), reinforcing this idea.

We stress that this interpretation is based on weak evidence and a very small sample, and should be considered with caution: the posterior functions of these parameters are broad and do not discard canonical values, but simply make them slightly less probable. The hint of this phenomenon arises from the fact that the three sources under study show this same marginal behaviour. The usage of tapered-edge surface densities profiles (e.g. Lynden-Bell & Pringle 1974; Hartmann et al. 1998) or puffed up inner rims (e.g. Dullemond, Dominik & Natta 2001) may also help explaining the anomalous SED slopes. Further evidence for flattened discs can be obtained by combining accretion and [O I] measurements. Keane et al. (2014) found the *Herschel* [O I] flux of 26 TDs to be ~2 times fainter than those of full discs, suggesting smaller gas-to-dust ratios compared to Class II discs, or smaller flaring indices. If the first scenario is ruled out by detecting significant gas reservoirs (via accretion signatures), then the flatter discs explanation would be favoured. Resolved ALMA observations of larger samples of (pre)TDs and full discs will reveal their real gas content and probe their structure, shedding light on this open issue.

## 6 CONCLUSIONS

We use *Herschel* photometry of three TDs in the Chamaeleon I star-forming region to perform detailed MCMC modelling of their SEDs and study the impact of *Herschel* data in the obtained results. We find that *Herschel* photometry, specially from the SPIRE instrument, can be used to constrain the dust mass in discs within one order of magnitude, as shown by the obtained posterior distributions. *Herschel* data can also help narrowing down the location of the inner radius of the disc. Our results are in good agreement with previous studies.

For the modelled targets, we find disc masses comparable to those of Class II sources in other star-forming regions. Because TDs are likely to represent a more evolved stage of disc evolution, the fact that they do not have significantly lower masses could suggest that the typical transitional class (i.e. discs with large gaps in their dust distributions) is the evolutionary outcome of massive Class II sources, with enough mass to form several giant planets which may have cleared their inner regions. Additionally, we find marginal hints of some dust settling and/or stepper surface density profiles in TDs than in protoplanetary discs. However, this result is tentative and requires further analysis. A larger sample of TDs, combined with gas and accretion measurements as well as resolved

images of these targets could help solving this issue and shed light on the origin of TDs and their real connection with planets.

Given the importance of disc masses for planet formation theories, the results obtained in this study open exciting new options to study this parameter for a large number of targets which lack (sub)mm observations but are present in the *Herschel* Science Archive. Further calibration of these values could also be achieved with more precise disc mass measurements from mm observations. Such a large-scale study could identify underlying relations between the stellar properties, disc masses, and the characteristics of planetary systems.

## ACKNOWLEDGEMENTS

We thank the referee Simon Casassus for his review and constructive comments, which helped improving and focusing the paper. This publication has been possible thanks to funding from the National Science Foundation. The data processing and analysis presented in this paper have made extensive use of the following open source software programs: PYTHON, IPYTHON, NUMPY, SCIPY, and MATPLOTLIB. This work has also made a significant use of TOPCAT (http://www.star.bristol.ac.uk/~mbt/topcat/ Taylor 2005). We are grateful to the developers of these softwares for their contributions. This work is based on observations made with the *Spitzer* Space Telescope, which is operated by the Jet Propulsion Laboratory, California Institute of Technology under a contract with NASA. This publication makes use of data products from the Two Micron All Sky Survey, which is a joint project of the University of Massachusetts and the Infrared Processing and Analysis Center/California Institute of Technology, funded by the National Aeronautics and Space Administration and the National Science Foundation.

## REFERENCES

Alexander R., Pascucci I., Andrews S., Armitage P., Cieza L., 2014, in Beuther H., Klessen R. S., Dullemond C. P., Henning T., eds, Protostars and Planets VI. Univ. Arizona Press, Tucson, AZ, p. 475

ALMA Partnership et al., 2015, ApJ, 808, L3

André P. et al., 2010, A&A, 518, L102

Andrews S. M., Williams J. P., 2005, ApJ, 631, 1134

Andrews S. M., Wilner D. J., Hughes A. M., Qi C., Dullemond C. P., 2009, ApJ, 700, 1502

Andrews S. M., Wilner D. J., Hughes A. M., Qi C., Dullemond C. P., 2010, ApJ, 723, 1241

Andrews S. M., Wilner D. J., Espaillat C., Hughes A. M., Dullemond C. P., McClure M. K., Qi C., Brown J. M., 2011, ApJ, 732, 42

Andrews S. M., Rosenfeld K. A., Kraus A. L., Wilner D. J., 2013, ApJ, 771, 129

Balog Z. et al., 2014, Exp. Astron., 37, 129

Beckwith S. V. W., Sargent A. I., Chini R. S., Guesten R., 1990, AJ, 99, 924

Belloche A. et al., 2011, A&A, 527, A145

Birnstiel T., Andrews S. M., Ercolano B., 2012, A&A, 544, A79

Bustamante I., Merín B., Ribas Á., Bouy H., Prusti T., Pilbratt G. L., Andr P., 2015, A&A, 578

Chiang E. I., Goldreich P., 1997, ApJ, 490, 368

Cieza L. A. et al., 2011, ApJ, 741, L25

D'Alessio P., Cantö J., Calvet N., Lizano S., 1998, ApJ, 500, 411

D'Alessio P., Calvet N., Hartmann L., Franco-Hernández R., Servín H., 2006, ApJ, 638, 314

Davis S. S., 2005, ApJ, 627, L153

Desch S. J., 2007, ApJ, 671, 878

Draine B. T., 2006, ApJ, 636, 1114

Dullemond C. P., Dominik C., 2004, A&A, 417, 159

Dullemond C. P., Dominik C., Natta A., 2001, ApJ, 560, 957
Espaillat C. et al., 2007a, ApJ, 664, L111

Espaillat C., Calvet N., D'Alessio P., Hernández J., Qi C., Hartmann L., Furlan E., Watson D. M., 2007b, ApJ, 670, L135

Espaillat C. et al., 2010, ApJ, 717, 441

Espaillat C., Furlan E., D'Alessio P., Sargent B., Nagel E., Calvet N., Watson D. M., Muzerolle J., 2011, ApJ, 728, 49

Espaillat C. et al., 2014, in Beuther H., Klessen R. S., Dullemond C. P., Henning T., eds, Protostars and Planets VI. Univ. Arizona Press, Tucson, AZ, p. 497

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

Furlan E. et al., 2006, ApJS, 165, 568

Furlan E. et al., 2011, ApJS, 195, 3

Goodman J., Weare J., 2010, Commun. Appl. Math. Comput. Sci., 5, 65

Griffin M. J. et al., 2010, A&A, 518, L3

Guenther E. W., Esposito M., Mundt R., Covino E., Alcalá J. M., Cusano F., Stecklum B., 2007, A&A, 467, 1147

Hartmann L., Calvet N., Gullbring E., D'Alessio P., 1998, ApJ, 495, 385

Harvey P. M. et al., 2012, ApJ, 755, 67

Hayashi C., 1981, Prog. Theor. Phys. Suppl., 70, 35

Henning T., Pfau W., Zinnecker H., Prusti T., 1993, A&A, 276, 129

Houck J. R. et al., 2004, ApJS, 154, 18

Indebetouw R. et al., 2005, ApJ, 619, 931

Ireland M. J., Kraus A. L., 2008, ApJ, 678, L59

Isella A., Natta A., 2005, A&A, 438, 899

Keane J. T. et al., 2014, ApJ, 787, 153

Kim K. H. et al., 2009, ApJ, 700, 1017

Lebouteiller V., Barry D. J., Spoon H. W. W., Bernard-Salas J., Sloan G. C., Houck J. R., Weedman D. W., 2011, ApJS, 196, 8

Lodato G., Delgado-Donate E., Clarke C. J., 2005, MNRAS, 364, L91

Luhman K. L., 2004, ApJ, 602, 816

Luhman K. L., 2007, ApJS, 173, 104

Luhman K. L. et al., 2008, ApJ, 675, 1375

Lynden-Bell D., Pringle J. E., 1974, MNRAS, 168, 603

Manoj P. et al., 2011, ApJS, 193, 11

Men'shchikov A., André P., Didelon P., Motte F., Hennemann M., Schneider N., 2012, A&A, 542, A81

Merín B. et al., 2010, ApJ, 718, 1200

Muzerolle J. et al., 2009, ApJ, 704, L15

Nagel E., Espaillat C., D'Alessio P., Calvet N., 2012, ApJ, 747, 139

Najita J. R., Strom S. E., Muzerolle J., 2007, MNRAS, 378, 369

Najita J. R., Andrews S. M., Muzerolle J., 2015, MNRAS, 450, 3559

Olofsson J., Szűcs L., Henning T., Linz H., Pascucci I., Joergens V., 2013, A&A, 560, A100

Ott S., 2010, in Mizumoto Y., Morita K.-I., Ohishi M., eds, ASP Conf. Ser. Vol. 434, Astronomical Data Analysis Software and Systems XIX. Astron. Soc. Pac., San Francisco, p. 139

Panić O., Hogerheijde M. R., Wilner D., Qi C., 2008, A&A, 491, 219

Pearson C. et al., 2014, Exp. Astron., 37, 175

Pilbratt G. L. et al., 2010, A&A, 518, L1

Pinte C., Ménard F., Duchêne G., Bastien P., 2006, A&A, 459, 797

Pinte C., Harries T. J., Min M., Watson A. M., Dullemond C. P., Woitke P., Ménard F., Durán-Rojas M. C., 2009, A&A, 498, 967

Poglitsch A. et al., 2010, A&A, 518, L2

Ribas Á. et al., 2013, A&A, 552, A115

Rodgers-Lee D., Scholz A., Natta A., Ray T., 2014, MNRAS, 443, 1587

Roussel H., 2013, PASP, 125, 1126

Spezzi L. et al., 2013, A&A, 555, A71

Taylor M. B., 2005, in Shopbell P., Britton M., Ebert R., eds, ASP Conf. Ser. Vol. 347, Astronomical Data Analysis Software and Systems XIV. Astron. Soc. Pac., San Francisco, p. 29

Thi W.-F. et al., 2010, A&A, 518, L125

Thi W.-F. et al., 2014, A&A, 561, A50

Trotta R., 2008, Contemp. Phys., 49, 71

Ubach C., Maddison S. T., Wright C. M., Wilner D. J., Lommen D. J. P., Koribalski B., 2012, MNRAS, 425, 3137

van der Marel N. et al., 2013, Science, 340, 1199

Whittet D. C. B., Prusti T., Franco G. A. P., Gerakines P. A., Kilkenny D., Larson K. A., Wesselius P. R., 1997, A&A, 327, 1194

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

Winston E. et al., 2009, AJ, 137, 4777

Winston E. et al., 2012, A&A, 545, A145

Zhu Z., Nelson R. P., Hartmann L., Espaillat C., Calvet N., 2011, ApJ, 729, 47

## APPENDIX A: CORNERPLOTS FOR THE CONSIDERED TRANSITIONAL DISCS

In this appendix, we provide the cornerplots obtained with the adopted MCMC procedure (Figs A1 to A3).
<figure><img src="image_6.png" /><figcaption>**Figure A1.** Cornerplot for SZ Cha. Histograms show the posterior distribution for each free parameter, scatter plots display the position of each chain in two parameter spaces to trace degeneracies. The results without *Herschel* data are shown in blue, those including *Herschel* in red.</figcaption></figure>
<figure><img src="image_7.png" /><figcaption>Figure A2. Cornerplot for CS Cha. Scheme and colours as in Fig. A1.</figcaption></figure>
<figure><img src="image_8.png" /><figcaption>Figure A3. Cornerplot for T25. Scheme and colours as in Fig. A1.</figcaption></figure>

This paper has been typeset from a TeX/LA TeX file prepared by the author.