© ESO, 2012

1. Introduction

High spatial-resolution imaging polarimetry is a rapidly developing field. It is now recognized that imaging polarimetry is a powerful tool for imaging faint matter around a bright star. Polarimetric differential imaging (PDI) is used to suppress the glare of the central star, which is unpolarized, to be able to see the faint reflected light from the material surrounding the star, which becomes polarized by scattering. This has already been successfully used to image circumstellar disks around young stars (see e.g. Close et al. 1997; Hashimoto et al. 2011) and circumstellar dust shells around evolved stars (e.g. Gledhill et al. 2001), and it has been proposed for imaging the light reflected from exoplanets (Stam et al. 2004). Images in polarized intensity not only contain information on the geometry of the circumstellar material, but also carry information on the particles scattering the light from the central star. In the general case, these two characteristics are intertwined and not easily separated.

Perhaps the most challenging case for interpreting polarimetric images is that of optically thick circumstellar disks around newly born stars. First, being optically thick causes multiple scattering, which alters the degree of polarization. The geometry of the surface of these disks can also cause scattering at unexpected angles, in turn causing certain features in polarization (Perrin et al. 2009). These effects can complicate the interpretation of polarimetric images in terms of planet-disk interactions. Besides these geometrical effects on a large scale, a second complicating factor arises when these images are used to derive general information about the dust particles in the disk. Dust particles are expected to be highly complex in shape and might be quite large in size, since they are a mixture of small grains and large, complex fluffy aggregates (Dominik 2009). Computing the optical properties of such particles is very difficult and computationally demanding, if feasible at all. In addition, the different grain sizes are not expected to be mixed well with the gas, but size-dependent settling of the dust particles causes a size distribution of the grains that changes with height in the disk (Dullemond & Dominik 2004). To conclude this chain of complications, dust properties, i.e. grain size and composition, are expected to vary radially (see e.g. van Boekel et al. 2004). Previous studies have already covered large parts of the parameter space of dusty circumstellar disks in the context of polarimetric imaging. The effects of disk mass and radial extend were studied in e.g. Whitney & Hartmann (1992). More recently, Murakawa (2010) has shown computations of near infrared polarimetric imaging of disks with various grain sizes and parameterized thickness of the disk.

An additional complication is the effect of the finite resolution of real telescopes on the final image. The most powerful diagnostic, the degree of polarization of the scattered light, is highly influenced by the glare from the central star. Also, blurring of the high-intensity, polarized, inner regions of the disk will affect the appearance of the disk in polarized intensity. Furthermore, data reduction techniques that correct for instrumental polarization using the central star as an unpolarized calibration point might suffer from the non-zero innermost disk polarization.

In this paper we present computations of the polarimetric images of protoplanetary disks around young, intermediate-mass stars, so-called Herbig Ae stars. We study how the images are effected by the structure of the dust grains and size-dependent grain settling. We varied the size of the dust particles, the geometry of the disk, and the shape of the dust particles. This was done using state-of-the-art models for the optical properties of the dust particles combined with a highly efficient radiative transfer code. An instrument simulator for a dual-beam imaging polarimeter at a 4 m class telescope was used to study the effects of atmospheric seeing and data reduction artifacts.

2. Model setup and computational approach

In this section we describe the setup of the models in terms of dust grains and disk geometry. In Sect. 2.1 we motivate our choice for the parameters used in the computations. The radiative transfer model used to create the polarimetric images and the instrument simulator are also described.

2.1. Parameter space

The characteristics of dusty, protoplanetary disks are described by a large set of parameters ranging from the microscopic properties of the dust grains to the macroscopic geometrical properties of the disk. Whitney & Hartmann (1992) investigated the effects of disk mass, radial extent, and the amount of disk flaring. They conclude that a flat disk would not be detectable. Also, they conclude that the large extent of many observations available at that time was most likely caused by an extended infalling envelope. The appearance of such an envelope in combination with the effects of outflows was studied in Whitney & Hartmann (1993) and Fischer et al. (1996) for a wide range of envelope and outflow parameters. More recently, Murakawa (2010) has presented computations of near infrared polarimetric imaging of disks. These theoretical studies all provide a good analysis of the dependence of polarimetric imaging on certain parameters of the dusty surroundings of young stars. It is therefore not the aim of this paper to repeat this and present a study covering the entire parameter space.

There are a few aspects that have not received much attention in previous studies, so we focus on these aspects and study their effects on the observed polarimetric images. These aspects are

• grain structure;

• size dependent grain settling;

• atmospheric blurring, and other observational aspects.

Below we explain how we set the model up specifically keeping these three main points in mind.

2.2. Dust grain model

2.2.1. Composition

The composition of protoplanetary dust is far from trivial (for a recent review see Min & Flynn 2010). To arrive at a simple, yet quite accurate model of the composition, we follow the procedure also adopted by Min et al. (2011). First we define the abundances of all available dust-forming elements to be the same as the solar composition derived by Grevesse & Sauval (1998). In meteorites and interstellar dust particles (IDPs) in our solar system all of the available sulfur is in the form of iron sulfide. Therefore, we first put all S in FeS thereby taking away about half of the available Fe atoms. From the remaining Fe, Mg, and Si we produce an amorphous Fe/Mg silicate with a stoichiometry in between olivine and pyroxene. We model this mixture by averaging the optical properties of different types of amorphous silicates. In addition to this iron sulfide/silicate mixture we add a contribution from carbonaceous matter. It is found that in our solar system about half of the carbon atoms are in the solid phase (Geiss 1987), so we add this abundance in the form of amorphous carbon grains. The composition resulting from this contains by mass 58% silicates, 18% iron sulfide, and 24% amorphous carbon.

To compute the opacities of the dust grains, we need the refractive index as a function of wavelength. This has been measured in the laboratory. For the silicates we use the measurements by Dorschner et al. (1995) and Henning & Stognienko (1996), for the iron sulfide we take the laboratory data from Begemann et al. (1994), and for the carboneceous material we take the measurements from Preibisch et al. (1993).

It can be argued that grains in the outer regions of protoplanetary disks contain a coating of ice. However, there is a recent debate about whether ice layers can survive in the upper layers of the disk where it is exposed to UV radiation from the central star (see Grigorieva et al. 2007). We checked how the scattering properties are influenced when adding water ice to the mixture in an abundance of 45% by mass, which is the abundance expected when putting all remaining oxygen in water ice (Min et al. 2011). The scattering properties are affected only slightly. The single scattering albedo goes up by less than 15% in all cases. The extinction coefficient, the maximum degree of polarization, and the asymmetry of the phase function are hardly affected. The resulting polarimetric images are thus very similar to those without ice.

Fig. 1 The scattering properties of the dust particles for different values of the maximum particle size, amax. Upper panels: the normalized phase function, F11, as a function of scattering angle, θ. Lower panels: the degree of linear polarization as a function of scattering angle. The wavelength is λ = 600nm. Solid lines are for the fluffy aggregates, dashed lines for the compact grains.

2.2.2. Grain structure

We distinguish between two extreme types of particle structures using the above composition:

• 1.

  compact, homogeneous grains, and

• 2.

  fluffy, inhomogeneous aggregates.

These two different particle structures require different ways of computing the optical properties. To correctly estimate the scattered light intensity and degree of polarization, the model used for the optical properties of the grains is crucial. Homogeneous spherical grains cause resonances at certain scattering angles and generally behave differently from more natural, irregular shapes. However, for reference and comparison to previous work, we also computed several models using compact, homogeneous spheres.

Compact, homogeneous grains: this model considers a mixture of grains, each composed of a single material. The optical properties of such a mixture can be obtained by simply averaging the optical properties of homogeneous grains using the given abundances of the various materials. We chose here to model the optical properties of the dust grains using the statistical approach. In this approach we take a wide distribution of simple particle shapes to simulate the variety of irregularities occurring in nature. We use the distribution of hollow spheres (DHS; see Min et al. 2005) with an “irregularity parameter” f_max = 0.8. It turns out that this choice of this grain-shape model, together with our choices for the grain sizes and compositions, provides a scattering matrix that is quite similar to cometary dust.

Fluffy, inhomogeneous aggregates: when small dust grains collide, they can form large aggregated structures. The optical properties of such aggregates display characteristics of both the constituents and the size of the aggregate as a whole (see e.g. Volten et al. 2007). We use here the approximate method by Min et al. (2008) to compute the optical properties of large, inhomogeneous fluffy aggregates. In this method we use effective medium theory to mix the dust material with vacuum and thus obtain a highly porous, fluffy structure. The vacuum fraction can be computed from the aggregate size and fractal dimension. Though this method was not primarily designed to simulate scattering at optical wavelengths, we find that many of the characteristics attributed to fluffy aggregates are simulated well. In general we find that the grains computed using these methods with the parameters given by Min et al. (2008) are extremely fluffy and as such give an extremely forward-peaked scattering phase function and a very high degree of polarization at a 90° scattering angle.

2.2.3. Grain size

To address the effects of grain growth on the observable polarimetric signal we use a parameterized grain size distribution. We adopt a power law size distribution with n(a) ∝ a^-3. Here a is the volume equivalent radius of the grains and the distribution runs from a_min to a_max. We fix a_min to 0.03μm and vary a_max from 0.3μm to 3mm. The different grain sizes are treated separately in the radiative transfer. This implies that they can have different temperatures and a different vertical distribution caused by grain settling.

Figure 1 displays the single scattering phase function and the degree of linear polarization as functions of the scattering angle, θ, for different values of a_max at a wavelength of λ = 0.6μm. The most important scattering characteristics are summarized in Table 1. Shown are the mass extinction coefficient, κ_ext, the single scattering albedo, ω, the asymmetry parameter of the phase function, ⟨cosθ⟩, and the maximum degree of polarization from single scattering, P_max.

2.3. Disk model

For the disk model we try to avoid parameterization as much as possible. Therefore we compute the vertical scale height of the disk self-consistently. This means that the vertical density distribution is fully determined by the temperature structure of the disk assuming vertical hydrostatic equilibrium (see Dullemond et al. 2007, for a review of the theoretical picture of the structures of protoplanetary disks). In addition, we compute the radial location and structure of the inner rim using evaporation and sublimation physics as described in Kama et al. (2009). This results in an inner edge of the disk that is curved outwards and slightly higher than the regions behind it. For details we refer to Kama et al. (2009).

For the central star we take a Herbig Ae star. We use two different disk types, a flaring disk, where we assume gas and dust to be perfectly coupled, and a settled disk, where we use the settling prescription from Dubrulle et al. (1995) as described above. These two disk types approximately represent the Group I and II sources as defined by Meeus et al. (2001). All star and disk parameters are listed in Table 2.

As the grains grow to larger sizes they decouple from the gas and settle towards the midplane (Dullemond & Dominik 2004). This is thought to be a possible explanation for the flatter disk structures observed around some stars. A full treatment of settling is beyond the scope of this paper and will be implemented in the future in the MCMax radiative transfer code. Here we parameterize settling using the approximate method by Dubrulle et al. (1995). In this way we do simulate the dependence of settling efficiency on grain size, grain density, gas density, and gas temperature. To distinguish between effects of grain settling and grain optical properties we take for the grain settling parameters always those of the material volume equivalent sphere. In principle, settling is more efficient for homogeneous compact grains and less efficient for fluffy aggregates. This is simply because the larger surface area of the fluffy grains allows them to couple much more easily to the gas. However, since we want to be able to separate disk structure effects as much as possible from effects of grain properties, we have chosen to not let this effect influence the settling efficiency and determine the settling properties of all grains as if they were spherical homogeneous grains with the same mass. The description of grain settling we use has one free parameter, which is the strength of the turbulence, α. This parameter defines the balance between settling and vertical mixing. For the settled disk models we take α = 10^-4, while for the well mixed models we set α = ∞.

Fig. 2 Sketch of the two disk geometries. The left panel shows the flaring, well-mixed disk structure, while the right panel depicts the settled, flat disk.

The radial density distribution of the dust disk was parameterized using a radial surface density (Hughes et al. 2008) (1)for r < R_out. Here R₀ is the turnover point from where an exponential decay of the surface density sets in and p sets the power law in the inner region. We fix this to p = 1, which is a commonly used value (see e.g. Dullemond et al. 2006). Although we take the outer radius of the disk to be 2000 AU, the surface density of the disk is already highly diminished after the turnover radius R₀ = 200AU.

In Fig. 2 we depict the two different disk structures considered. This figure illustrates that the scattering angles of the surface layers of the disk depend strongly on the opening angle. In the case of a very flat disk, the scattering angles observed at the front and the back of the disk are symmetric around 90°. With i the inclination angle of the disk the scattering angles at front and back for a very flat disk are θ = 90° − i and θ = 90° + i respectively. For a strong flaring disk with a large opening angle ψ, the scattering angles at front and back are symmetric around 90° − ψ and become θ = 90° − ψ − i and θ = 90° − ψ + i.

We would like to emphasize that by using the self-consistent treatment described above we limit the amount of free parameters significantly. Besides the parameters we do vary, we have two important parameters that determine the observed polarization signal; R₀ and M_dust. The impact of these two parameters on the polarimetric signal of circumstellar disks was studied in, e.g., Whitney & Hartmann (1992).

2.4. Radiative transfer

The transfer of radiation through the disk is done in full using the 3D axisymmetrical radiative transfer code MCMax (Min et al. 2009). This code uses Monte Carlo radiative transfer to compute the disk temperature structure and from this the vertical density scale-height and the sublimation state of the dust (see also Kama et al. 2009). The scattering of radiation uses the full scattering matrix acting on the Stokes vectors. The grains in the disk are assumed to be in random orientation so symmetry arguments can be used to simplify the scattering matrix.

Often images and other observables are extracted from Monte Carlo radiative transfer methods by simply collecting the escaping photons after they are traced through the disk. Although this is an elegantly simple way, we adopt a different technique to construct the images which allows us to obtain accurate results even in the very faint parts of the image. After the temperature structure of the disk is obtained, images are constructed by integrating the formal solution to the radiative transfer. For this, the local scattering source function is also needed; i.e., at each volume element in the disk we need to know how much radiation is scattered towards the observer. This is again computed using a Monte Carlo approach. However, this time we store the full 3D radiation field at each location in the disk by using the full path traveled by the Monte Carlo photons (as inspired by the path length method for temperature determination by Lucy 1999). In combination with the local thermal source function, which is obtained from the temperature, the formal solution can be easily integrated to obtain smooth images. This is done by integrating from each volume element in the disk the total emitted flux (as obtained from the local temperature) and the total scattered flux (as obtained from the 3D radiation field) diminished by the extinction through the disk along the line of sight towards the observer.

The advantage of the above method over “classical Monte Carlo” photon-collecting methods is that the information accumulated during the Monte Carlo photon-tracing process is used to the fullest to obtain the local 3D radiation field. This way smooth images can also be produced in regions where only a small fraction of the photons arrive. This is especially important for imaging polarimetry where one of the main purposes is to gain a very high dynamic range in the image. The method was benchmarked against other radiative transfer codes in Pinte et al. (2009).

2.5. Instrument simulation

Ground-based observations are affected by atmospheric and instrumental effects such as seeing and instrumental polarization. To account for these effects, a numerical simulation including readout noise, photon noise, instrumental polarization, and different seeing conditions was performed. These simulations model the dual-beam imaging polarimeter ExPo (the Extreme Polarimeter; Rodenhuis et al. 2008 and in prep.) currently a visitor instrument at the 4.2 m William Herschel Telescope (WHT). ExPo measures linear polarization at optical wavelengths, and combines very short exposure times with the dual-beam technique to minimize systematic errors. In a standard observation, thousands of single exposure frames are combined to obtain a total exposure time of  ~15 min.

The point spread function (PSF) for ExPo is calculated as (2)where ℱ stands for Fourier transform, A and φ refer to the amplitude and phase of the incoming wavefronts, respectively. A set of 100 statistically independent wavefronts is calculated for different seeing conditions. We compute the phase of these wavefronts according to the Kolmogorov theory of turbulence (for details see Appendix A.1). Broadband PSFs are obtained from these sets of wavefronts combined with the aperture of the telescope after which they are convolved with full-resolution simulated images. Readout noise, photon noise, instrumental polarization, and guiding errors are added to these images. This set of simulated-raw images is then reduced with the ExPo Pipeline, described in detail in Canovas et al. (2011).

The characteristics of this simulation are summarized in Table A.1. The simulation of the atmosphere, telescope, and instrument is detailed in Appendix A.

Fig. 3 Comparison of the polarized intensity maps for different grain types. Upper panels are for the flared disks, lower panels for the settled disks. Left panel: homogeneous spheres; middle panels: compact, irregularly shaped grains; right panels: fluffy aggregates. For all images, the maximum size of the particles is 30 μm. Spatial scale is in arcseconds. The length of the vectors is proportional to the degree of polarization. The inclination of the disk is 55° with respect to pole on.

Fig. 4 Comparison of the polarized intensity maps for different grain sizes. All panels are for the well mixed, i.e. flaring, disk models. Upper panels are for the compact grains, lower panels for the fluffy aggregates. The maximum particle size in the size distribution increases from left to right. Spatial scale is in arcseconds. The length of the vectors is proportional to the degree of polarization. The inclination of the disk is 55° with respect to pole on. Color scale is the same as in Fig. 3.

3. Resulting full-resolution images

In this section we discuss the intrinsic polarimetric images of the disks and the effects of the different dust properties and disk parameters. All full-resolution model images are generated at a wavelength of λ = 600nm. Though detailed studies of the color behavior of the polarized images can be a useful diagnostic, the general monochromatic behavior of the images is not very sensitive to the exact wavelength. We use these single wavelength images to represent the average broadband signal, so the broadband effects are fully taken into account for the simulation of the telescope and instrument, but the full-resolution model disk images are computed at a single wavelength.

In Fig. 3 we present polarization maps for the different dust types. We show here images for our two main grain types (compact grains and fluffy aggregates) and also for compact homogeneous spherical grains. It is clear from this figure that the type of dust grain has a strong influence on the appearance of the disk. Especially in the upper part of the image, the homogeneous spherical dust grains display a behavior not seen for the more realistic dust models. The homogeneous spherical grains display here a decreasing polarized intensity towards the upper central region of the disk. This feature is seen in a much milder form in the image of the fluffy aggregates, but generally the images for the other two grain shape models show smoother behavior in this part of the image. Also, it is clear that the fluffy aggregates are fainter than the compact grain models. This is caused by the extremely forward-scattering nature of the fluffy aggregates, where only a small fraction of the light is scattered towards the observer.

In Figs. 4 and 5 we show the polarized intensity maps for the compact grains and the fluffy aggregates for different values of the maximum grain size in the size distribution. For the flared disk models (Fig. 4), it is clear that the brightness and spatial extent of the disk decreases when the particle size is increased. This is caused by the reduced opacity of the larger particles, which decreases the flaring angle of the disk. For the settled disk models we see a slightly different trend. Because the surface layer of the disk is here already in the shadow, the decrease in the height of the disk is not so important. More important in this case is the increase in the scattering albedo, causing the disk to become brighter and appear larger at first when the particle size is slightly increased. When the particle size is increased further, the decrease in the total opacity becomes dominant and the disk becomes fainter again.

Murakawa (2010) revived the discussion on the “roundabout effect” for edge on disks in near infrared polarimetric imaging. In the dark lane of an edge on disk no first-order scattering can be observed since no starlight can reach these regions in the disk. Second-order scattering causes the polarization vectors not to be aligned perpendicular to the direction towards the central star, but parallel to the dark lane of the disk. This effect was first noted by Bastien & Menard (1988), after which Whitney & Hartmann (1992) argued that detecting this feature would be very hard because of the low surface brightness in this part of the image. We also observe this effect in our optical polarization images of the most flaring edge-on disks. Figure 6 shows the polarized intensity maps for the disks seen edge-on for different grain sizes. We only show the compact grain models, but the fluffy aggregate models show similar behavior. Clearly we see the roundabout effect in the flaring disk with only small grains. As soon as the grain size is increased, or the disk is settled, the effect disappears because the dark lane becomes too narrow. The region where the effect can be observed is very weak in polarized intensity as compared to the region above and below. This implies that, in observations with a limited spatial resolution, the effect will be washed away by PSF smearing.

Fig. 5 Same as Fig. 4, but for the settled disk models.

Fig. 6 Comparison of the polarized intensity maps for different grain sizes in the nearly edge on case. All panels are for the flared disk models with compact dust grains. The settled disks display only a very thin signature. The maximum particle size in the size distribution increases from left to right. Spatial scale is in arcseconds. The length of the vectors is proportional to the degree of polarization. The inclination of the disk is 85° with respect to pole on. Color scale is the same as in Fig. 3.

Fig. 7 Same as Fig. 4, but simulated as observed with ExPo at the WHT. Color scale is logarithmic, in arbitrary flux units, and the same in all images.

Fig. 8 Same as Fig. 5, but simulated as observed with ExPo at the WHT.

4. Simulated observational images

We present images simulated for the ExPo instrument mounted at the 4.2 m WHT. The aim here is to get a feeling for the effects of atmospheric seeing, telescope diffraction, and noise on the interpretation of polarimetric images.

The final image observed through a telescope in polarization from a protoplanetary disk not only depends on the intrinsic degree of polarization caused by the grains but is also influenced by the depolarizing effects of the PSF of the central star. How strong this depolarizing effect is depends on the intensity of light scattered off the disk. Below we discuss the effects of the different grain types on the intrinsic polarization image, as well as how these different models are affected by telescope effects such as diffraction, seeing, and noise.

In Figs. 7 and 8 we show the different disk and dust types as they would be observed with ExPo at the WHT under good seeing conditions, i.e. a seeing at 500 nm of 0.8′′. There is a clear dependence of the detectibility on the particle size distribution and the vertical disk geometry. The larger and the more settled disks are much harder to detect, since they are intrinsically fainter.

The polarized intensity is the product of the total intensity and the degree of polarization. Due to the opening angle of the disk, the back side of the disk is closer to 90° scattering angle and thus has a higher degree of polarization (see also Fig. 2). This can be seen in the high-resolution images by the longer vectors in the upper parts of the images. The front side of the disk is, however, much more forward scattering and thus has a higher total intensity. Because the central star is the dominant intensity source in the observational images, the degree of polarization now is basically the polarized intensity divided by the total intensity of the central star and not by the total intensity of scattered light. Looking at the observed settled disk models (i.e. Fig. 8), it can be seen that in the images for the compact grains, the degree of polarization still dominates (also see the length of the polarization vectors), and the back side of the disk is somewhat brighter in polarized intensity. Due to the highly forward-scattering nature of the fluffy aggregates, in that case the differences in total intensity dominate, and the front side of the disk is brighter in polarized intensity and in degree of polarization.

Since we run the full ExPo data analysis pipeline as discussed in Canovas et al. (2011), we also correct for any instrumental polarization. This is done by assuming that the signal from the central resolution element, i.e. the central star, should be unpolarized. In the simulation of the instrument we do not add polarization caused by the telescope, but we include instrumental polarization due to different transmission coefficients (see Appendix). Nevertheless, the central resolution element is polarized because of scattering from the innermost regions of the disk. Since the disk is inclined, scattering by this innermost region has a non-zero, integrated polarization (see also Fig. 10). This polarization is then ‘corrected’ under the assumption that it is instrumental, and thus subtracted from the entire image. This causes the vertical depolarized region in the images where the outer disk is weak. This also causes the well oriented polarization vectors in the noise regions of the images for the cases where the disk is very weak (especially clear in Fig. 8). The amplitude of this effect strongly depends on the seeing conditions. In Fig. 9 we show two disk models under different seeing conditions. It is clear that the appearance is different and the overcorrection of the instrumental and inner disk polarization is larger in the case of worse seeing.

In Fig. 10 we show the degree of polarization of the inner resolution elements as a function of the size of the resolution element. We do this by integrating both Stokes I and Q out to a distance R from the central star. The degree of polarization shown in Fig. 10 is the ratio of these two integrated values. This shows that increasing the spatial resolution also significantly reduces the effect of polarization from the central resolution element. For images limited by an atmospheric seeing of  ~1′′, we see that the inner resolution element can be polarized by up to  ~0.8%. By adding an extreme adaptive optics system increasing the spatial resolution to 0.02′′, this degree of polarization can be reduced significantly. In the case of fluffy aggregates that have grown to significant sizes, we see that the integrated degree of polarization is very low. This is because the total scattered intensity is dominated by forward scattering, which has a very low degree of polarization. Although the polarization of the central resolution element is low, the spatially resolved degree of polarization is also extremely low. The large variations in the curves in Fig. 10 show that we do have a good diagnostic for the grain and disk properties when we have sufficient resolution and the images are properly calibrated.

4.1. Detectability

For a disk to be detectable with an imaging polarimeter at a 4 m class telescope without adaptive optics, we need the disk to be extended enough to escape the PSF of the central star and bright enough to rise above the photon and readout noise. Treating the models as observations we have determined which disk and dust parameters allow detection with this instrument. For this we assumed a total integration time of about five minutes per polarization state. The visual magnitude of the star, as derived from the input parameters in Table 2, is V = 7. All disks are assumed to be at a distance of 150 pc. The detectability of the disks is based on the signal-to-noise ratio (S/N) of the polarized intensity images. The orientation of the polarization vectors depends on this parameter, so the polarization vectors calculated for an image with poor S/N will be randomly oriented (see, for instance, Fig. 7, bottom righthand image). On the other hand, an image with high S/N will show a clear pattern on the orientation of its polarization vectors. For disks that are hardly detected, we also looked at the images for Stokes Q and U to see if a so-called “butterfly” pattern, typical of polarimetric images of circumstellar material, emerges above the noise.

The results are summarized in Table 3. We see that disks that have undergone significant grain growth are hard to detect, especially when the turbulence in the disk is not strong enough to mix up the larger grains to the disk surface layers, i.e. in the settled disk models. Also, we notice that the seeing conditions play a major role in the detectability. By eliminating this factor, an adaptive optics system will greatly increase the detection capabilities of imaging polarimeters.

Fig. 9 Comparison of the polarized intensity maps as obtained under different seeing conditions. Upper panels are for a flared disk, lower panels for a settled disk. All images use amax = 30μm and compact grains.

Fig. 10 The integrated degree of polarization as a function of the radius of integration. Upper panel is for the compact grains, lower panel for the fluffy aggregates. Black curves are for the flared disk models, grey curves for the settled disk models.

4.2. Disk gaps

Many systems have been observed, directly or indirectly, to have depleted inner disks or cleared out gaps (see e.g. Bouwman et al. 2003; Brown et al. 2007; Espaillat et al. 2008; Thalmann et al. 2010; Merín et al. 2010; Andrews et al. 2011; Verhoeff et al. 2011). One of the possible causes for this is clearing by a planet in the disk, and thus direct detection of such disk gaps is very important for constraining planet formation scenarios and timescales. To see if we can detect such disk gaps and see what it would generally look like in a polarimetric imaging observation, we ran a disk model with a cleared out gap from 120 to 180 AU and an inner disk with a depletion of a factor of ten for radii smaller than 120 AU (a parameterization very similar to the one adopted by Andrews et al. 2011). We did not aim to span parameter space for this case, but rather to provide an example of inner disk clearing. The radii of the gap were chosen to be close to what one would expect to be able to detect given the spatial resolution we had. Tests with closer in gaps show that, indeed, the radii presented here is at the edge of the detection capabilities given the simulated spatial resolution.

The images for compact grains with a_max = 30μm are shown in Fig. 11, where we show both the full-resolution images (left) and the simulated observations (right). For the instrument simulations we increased the integration time by a factor of ten to suppress the photon noise in the inner disk regions. As we can see, the depletion of the inner region is clearly detected. Comparing to the lefthand panels of Fig. 9 (which shows the same models without depleted inner disks), we see a very clear difference. We do note that the inner disk is undetectable for the flaring disk case in the simulated observations, while it is still clearly visible in the full-resolution images. This is caused by the fact that the entire inner disk is inside the central resolution element (0.8′′ seeing). As a consequence the polarization signal caused by the inner disk is considered instrumental by the data reduction pipeline and is subtracted. We also computed images where we placed the gap closer in, and we find that an inner radius of the outer disk of  ~180AU (or 1.2′′) is close to the limit where we are still able to clearly identify the disk gap in the image. For radii of the inner edge of the outer disk larger than  ~100AU, we do see significant changes to the image, but they become increasingly more difficult to identify as a clear signal from a disk gap as the radius decreases.

Fig. 11 Full resolution images (left panels) and simulated observations (right panels) of models with a depleted inner disk. The innermost regions (radii  <120 AU) are depleted by a factor of 10, while from 120 to 180AU the disk is cleared of dust. Upper panels are for a flared disk, lower panels for a settled disk. All images are using amax = 30μm and compact grains. The color scale is altered with respect to the previous images to enhance details in the bright inner disk region.

5. Discussion and conclusions

The appearance of circumstellar disks in polarized light depends heavily on the properties of both the disk and the composing dust particles. We modeled the observable signature of protoplanetary disks surrounding young, intermediate-mass stars with various dust and disk properties. We focus on the effects of grain structure and size-dependent grain settling. We also studied the effects of atmospheric blurring and observational data processing. We conclude the following.

• There is a significant difference between the appearance of a diskwith only compact grains and a disk containing fluffy aggregateddust particles.

• The commonly employed homogeneous sphere model gives significantly different results from more realistic models, because the observed signal depends strongly on the properties of these more realistic models, thus providing a potentially powerful diagnostic. The proper use of this diagnostic is not always trivial.

• Due to the extreme forward-scattering nature of fluffy aggregates, where only a small fraction of the radiation is scattered towards the observer, a disk containing these kinds of particles appears significantly dimmer in scattered light.

• The diagnostic power of imaging polarimetry is a strong function of spatial resolution. Atmospheric seeing blurs many details needed for a thorough analysis of dust properties, and thus mainly allows the global disk parameters to be deduced.

• The “roundabout effect” for edge-on disks as discussed by Murakawa (2010) for near-IR polarimetry was only found in our study for an extremely flaring disk containing only small grains. Since the effect appears in the very low-intensity regions of the image, it disappears when we consider observations with limited spatial resolution.

• The accuracy of the correction for instrumental polarization using the central star as calibration source strongly depends on the spatial resolution. In the seeing-limited case, the central resolution element can have a degree of polarization of up to  ~0.8%. This degree of polarization quickly decreases when the spatial resolution of the observations is increased.

References

Online material

Appendix A: Instrument simulation

Here we discuss the details of the instrument simulations. The instrument simulator was developed to simulate observations with the Extreme Polarimeter (ExPo), which is a visitor instrument of the 4.2 m WHT.

A.1. Wavefront generation

To generate wavefronts, we make use of a model of turbulence first proposed by Kolmogorov (1941) and later developed by Tatarski (1961) and Fried (1966). According to this model, a flat wavefront traveling through the turbulent atmosphere will modify its phase, while the change in its amplitude is negligible compared to the phase fluctuations. After it travels through the telescope, the amplitude of the wavefront is proportional to the telescope’s aperture. The power spectrum of the phase fluctuations of a wavefront traveling through the atmosphere is given by (A.1)where r₀ is Fried’s Parameter and k the wave number. The parameter r₀ represents the circular aperture over which the wavefront phase variance is equal to 1 rad² for the case of Kolmogorov turbulence. Therefore, the higher r₀, the more stable the atmosphere, or, in other words, the better the seeing is. For instance, a value of r₀ of about 10 cm produces 1 arcsec seeing at λ = 500 nm (i.e. Nightingale & Buscher 1991). The value of r₀ depends on the wavelength, following a power law: r₀ ∝ λ^6/5. Figure A.1 shows an example of an uncalibrated wavefront phase generated from Eq. (A.1).

The timescale at which one of these wavefronts remains constant, the coherence time, at optical wavelengths is only a few milliseconds (Vernin & Munoz-Tunon 1994; Law 2006; García-Lorenzo et al. 2009). The single frame PSF has two main components:

• Tilt, producing random motion of the whole image.

• Roughness, producing the observed speckle pattern.

When the exposure time is similar to the coherence time, the image motion due to the tilt component can be removed in the data processing by proper centering of each frame. Noll (1976) showed that the variance σ² of the wavefronts can be expressed in terms of the telescope diameter, D, and the Fried parameter, r₀(A.2)The standard deviation σ of the wavefront phase ϕ, expressed in wavelength units, can be obtained from the variance as σ_λ = 2πσ(A.3)In our simulation, D is the diameter of the WHT main mirror (4.2 m), and three different values of r₀, (5, 7, and 10 cm) are used to account for bad, medium, and good seeing at a reference wavelength of 500 nm, λ₅₀₀. From Eq (A.2) and the dependence of r₀ with wavelength, σ_λ can be computed at different wavelengths from our reference as (A.4)Wavefronts at different wavelengths are calibrated according to Eqs. (A.2) and (A.4). A set of wavefronts at different wavelengths is then generated for each of the three different r₀ values used in this simulation.

A.2. PSF generation

To generate a PSF, we first simulate the aperture of the WHT, taking the spiders of the telescope and the central obscuration of the main mirror into account. Once the telescope aperture and the wavefront are computed, a monochromatic PSF is calculated according to Eq. (2). A broadband PSF (PSF_bb) is calculated as the sum of monochromatic PSFs calculated over the range of 400 nm to 700 nm, in steps of 10 nm: (A.5)This PSF produces the speckle pattern obtained when observing at exposure times that are close to the coherence time (a few milliseconds). The PSF measured when observing at longer exposure times will be the sum of the short-exposure broadband PSFs. To simulate ExPo-like PSFs, we set the coherence time to 9.3 s and then generate a 28 s exposure PSF as the sum of three statistically independent short-exposure broadband PSFs: (A.6)These calculations are computed on a 2048 × 2048 pixels grid. The pixel size of this simulation is determined by the Nyquist frequency: N_ν = (λ/2D)·206265, which produces a pixel size of 0.0122 [′′/pixel], for a 4.2 m telescope and a central wavelength of 500 nm. The original 2048 × 2048 pixel grid is then binned to an smaller grid to produce simulated images at the ExPo pixel size (0.078 [′′/pixel]). Figure A.2 shows an example of a PSF for both simulated (left) and real (right) data (0.8′′ seeing, 0.028 s exposure time).

A total of 100 different PSF_ExPo are produced for each of the three different seeing conditions tested here. Figure A.3 shows an example of three different broadband PSF’s calculated for bad (left), normal (center), and good (right) seeing.

Fig. A.1 Phase of a wavefront following the Kolmogorov model of turbulence.

Fig. A.2 Left: Simulated PSF with seeing  ≈ 0.8′′, 0.028 s exposure time. Right: observed PSF under the same conditions.

A.3. Dual beam simulation

As a dual beam instrument, ExPo produces two simultaneous images with opposite polarization states that are imaged onto a CCD. We refer to these two images as left and right images. These measurements are modulated by a Ferro-electric liquid crystal (FLC), which switches between two orthogonal polarization states (A and B). Two different frames, containing four different images are produced at the end of one FLC cycle: where neither instrumental polarization nor instrumental effects are considered for the sake of simplicity, and I₀ and P₀ represent the total intensity and polarization (Stokes Q or Stokes U) of the observed target, as seeing without atmosphere and diffraction

effects, respectively. PSF_A and PSF_B are the short exposure PSF for the A and B frames, respectively. These images are given in mJy units, and they are converted to counts by considering the telescope area, exposure time, filter transmission (Johnson V filter is used here), atmospheric + instrument absorption, and CCD efficiency.

Fig. A.3 Simulated speckle pattern for different values of r0, given in centimeters.

The poidev function from the NASA IDL library is used to simulate the photon noise for each of these images. Real readout noise from ExPo measurements is added to the simulated data. Instrumental polarization is simulated by including two transmission coefficients for each of the beams (T_L,T_R) and another two coefficients to account for the FLC transmission (T_A,T_B). The value of these coefficients is listed in Table A.2.

To mimic guiding problems, a random shift with a maximum amplitude of ten pixels is applied to each of the simulated images. The final simulated images are then described by where M(x,y) represents the image-shifting function, Ph is the photon noise, and RO represents the readout noise.

We finally run the full data reduction pipeline for the simulated observations to obtain the final images presented in Figs. 7–9

All Tables

All Figures

	Fig. 1 The scattering properties of the dust particles for different values of the maximum particle size, amax. Upper panels: the normalized phase function, F11, as a function of scattering angle, θ. Lower panels: the degree of linear polarization as a function of scattering angle. The wavelength is λ = 600nm. Solid lines are for the fluffy aggregates, dashed lines for the compact grains.
In the text

	Fig. 2 Sketch of the two disk geometries. The left panel shows the flaring, well-mixed disk structure, while the right panel depicts the settled, flat disk.
In the text

Fig. 3 Comparison of the polarized intensity maps for different grain types. Upper panels are for the flared disks, lower panels for the settled disks. Left panel: homogeneous spheres; middle panels: compact, irregularly shaped grains; right panels: fluffy aggregates. For all images, the maximum size of the particles is 30 μm. Spatial scale is in arcseconds. The length of the vectors is proportional to the degree of polarization. The inclination of the disk is 55° with respect to pole on.

In the text

Fig. 4 Comparison of the polarized intensity maps for different grain sizes. All panels are for the well mixed, i.e. flaring, disk models. Upper panels are for the compact grains, lower panels for the fluffy aggregates. The maximum particle size in the size distribution increases from left to right. Spatial scale is in arcseconds. The length of the vectors is proportional to the degree of polarization. The inclination of the disk is 55° with respect to pole on. Color scale is the same as in Fig. 3.

In the text

	Fig. 5 Same as Fig. 4, but for the settled disk models.
In the text

Fig. 6 Comparison of the polarized intensity maps for different grain sizes in the nearly edge on case. All panels are for the flared disk models with compact dust grains. The settled disks display only a very thin signature. The maximum particle size in the size distribution increases from left to right. Spatial scale is in arcseconds. The length of the vectors is proportional to the degree of polarization. The inclination of the disk is 85° with respect to pole on. Color scale is the same as in Fig. 3.

In the text

	Fig. 7 Same as Fig. 4, but simulated as observed with ExPo at the WHT. Color scale is logarithmic, in arbitrary flux units, and the same in all images.
In the text

	Fig. 8 Same as Fig. 5, but simulated as observed with ExPo at the WHT.
In the text

	Fig. 9 Comparison of the polarized intensity maps as obtained under different seeing conditions. Upper panels are for a flared disk, lower panels for a settled disk. All images use amax = 30μm and compact grains.
In the text

	Fig. 10 The integrated degree of polarization as a function of the radius of integration. Upper panel is for the compact grains, lower panel for the fluffy aggregates. Black curves are for the flared disk models, grey curves for the settled disk models.
In the text

Fig. 11 Full resolution images (left panels) and simulated observations (right panels) of models with a depleted inner disk. The innermost regions (radii  <120 AU) are depleted by a factor of 10, while from 120 to 180AU the disk is cleared of dust. Upper panels are for a flared disk, lower panels for a settled disk. All images are using amax = 30μm and compact grains. The color scale is altered with respect to the previous images to enhance details in the bright inner disk region.

In the text

	Fig. A.1 Phase of a wavefront following the Kolmogorov model of turbulence.
In the text

	Fig. A.2 Left: Simulated PSF with seeing  ≈ 0.8′′, 0.028 s exposure time. Right: observed PSF under the same conditions.
In the text

	Fig. A.3 Simulated speckle pattern for different values of r0, given in centimeters.
In the text