© ESO, 2016

1. Introduction

The sky at high Galactic latitude and frequencies above about 100 GHz is dominated by thermal emission from the Galactic interstellar medium, specifically arising from dust grains of about 0.1 μm in size. Asymmetrical dust grains align with the Galactic magnetic field to produce polarized emission. This polarized submillimetre emission has been measured from ground-based and balloon-borne telescopes (e.g., Hildebrand et al. 1999; Benoît et al. 2004; Ponthieu et al. 2005; Vaillancourt 2007; Matthews et al. 2014). The observed polarization relates to the nature, size, and shape of dust grains and the mechanisms of alignment, discussed for example by Draine (2004) and Martin (2007). It also probes the structure of the Galactic magnetic field, which is an essential component of models of Galactic dust polarization (Baccigalupi 2003; Fauvet et al. 2011, 2012; O’Dea et al. 2012; Jaffe et al. 2013; Delabrouille et al. 2013).

The polarized emission from dust is also of interest in the context of foregrounds (Tucci et al. 2005; Dunkley et al. 2009a; Gold et al. 2011) to the cosmic microwave background (CMB). On angular scales between 10′ and a few tens of degrees, cosmological B-mode polarization signals may be present that were imprinted during the epoch of inflation. The discovery of a primordial B-mode polarization signature is a major scientific goal of many CMB experiments. These include ground-based experiments (ACTPol, Niemack et al. 2010; BICEP2, BICEP2 Collaboration 2014a; Keck-array, Staniszewski et al. 2012; POLARBEAR, Arnold et al. 2010; QUBIC, Ghribi et al. 2014; QUIJOTE, Rubiño-Martín et al. 2010; and SPTpol, Austermann et al. 2012), stratospheric balloon missions (EBEX, Grainger et al. 2008; and SPIDER, Fraisse et al. 2013), and the ESA Planck¹ satellite (Tauber et al. 2010). Accurate assessment and, if necessary, subtraction of foreground contamination is critical to the measurement of CMB E- and B-mode polarization because the expected signals from inflation and late-time reionization are expected to be small.

Planck has measured the all-sky dust polarization at 353 GHz, where the dust emission dominates over other polarized signals. These data have been presented in a first set of publications in which the focus was on the structure of the Galactic magnetic field and the characterization of dust polarization properties (Planck Collaboration Int. XIX 2015; Planck Collaboration Int. XX 2015; Planck Collaboration Int. XXI 2015; Planck Collaboration Int. XXII 2015). Here, we use the Planck polarized data to compute the and power spectra of dust polarization over the multipole range 40 <ℓ< 600, on large fractions of the sky away from the Galactic plane. We also investigate dust polarization in sky patches at high Galactic latitude with sizes comparable to those surveyed by ground-based CMB experiments. We derive statistical properties of dust polarization from these spectra, characterizing the shape of the spectra and their amplitude with respect to both the observing frequency and the mean dust intensity of the sky region over which they are computed. We verify that these properties hold in low-column-density patches at high Galactic latitude and we explore statistically the potential existence of “clean” patches on the sky that might be suitable for cosmology.

Our analysis of dust polarization is relevant to the present generation of CMB polarization observations, and to the design of future experiments. It gives a statistical description of Galactic dust polarization, providing input for the modelling of Galactic dust as part of component separation methods and for CMB polarization likelihood analysis parameterization.

The BICEP2 collaboration has recently reported a significant detection of the B-mode power spectrum around the expected angular scale of the recombination bump, namely a few degrees (BICEP2 Collaboration 2014a,b). Their analysis was based on dust polarization models that predicted subdominant contamination of their B-mode signal by dust polarization. We use information from our detailed analysis of Planck polarization data at 353 GHz to assess the potential dust contamination.

The paper is organized as follows. In Sect. 2, we present the Planck HFI polarization data used in this work and describe the general properties of the polarization maps in terms of emission components and systematic effects. In Sect. 3, we describe our method for computing the dust and angular power spectra, including the selected science regions of interest on the sky. We assess and compare the two methods we use to compute the power spectra in Appendix A. In Sect. 4, we present power spectra of dust polarization for multipoles ℓ> 40, computed with high signal-to-noise ratio (S/N) on large fractions of the sky, and characterize their shape and amplitude. Complementary angular power spectra involving temperature and polarization, and , and cross-polarization,   are given in Appendix B. We extend this analysis to smaller sky patches at high Galactic latitude in Sect. 5. In Appendix C we discuss some complementary aspects of this analysis of patches. These results are used specifically in Sect. 5.3 to build a map of the expected dust contamination of the power spectrum at 150 GHz and ℓ = 80. In Sect. 6 we present a study of the polarized dust emission in the vicinity of the BICEP2 region. Systematic effects relating to the Planck angular power spectrum estimates are assessed in Appendix D and we discuss the decorrelation of the dust signal between frequencies in Appendix E. Section 7 summarizes the main conclusions and discusses the implications of this work for future CMB experiments².

2. Planck polarization maps

2.1. Planck data

The Planck Collaboration recently released the Planck satellite nominal mission temperature data and published a set of papers describing these data and their cosmological interpretation (e.g., Planck Collaboration I 2014; Planck Collaboration XVI 2014). These results are based on the data from the two instruments on board the satellite (Low Frequency Instrument, LFI, Mennella et al. 2011, and High Frequency Instrument, HFI, Planck HFI Core Team 2011). The data processing of the nominal mission data (Surveys 1 and 2, 14 months) was summarized in Planck Collaboration II (2014) and Planck Collaboration VI (2014).

Planck HFI measures the linear polarization at 100, 143, 217, and 353 GHz (Rosset et al. 2010). The properties of the detectors (sensitivity, spectral response, noise properties, beams, etc.) are described in detail in Lamarre et al. (2010) and their in-flight performance is reported in Planck HFI Core Team (2011), Planck Collaboration VII (2014), Planck Collaboration VIII (2014), Planck Collaboration IX (2014), and Planck Collaboration X (2014), while Planck Collaboration VI (2014) describes the general processing applied to the data to measure polarization. In this paper, we make use of full-mission (Surveys 1 to 5, 30 months, Planck Collaboration I 2014) polarization maps of Planck HFI (internal data release “DX11d”), projected into the HEALPix pixelization scheme (Górski et al. 2005). This is one of the first publications to use these maps, which will be described in the Planck cosmology 2015 release.

To compute polarization angular power spectra, we use Q and U maps at 100, 143, 217, and 353 GHz. Specifically, we calculate power spectra using the so-called “Detector-Set” maps (hereafter “DetSets”), constructed using two subsets of polarization sensitive bolometers (PSBs) at a given frequency (see Table 3 in Planck Collaboration VI 2014). Each DetSet polarization map is constructed using data from two pairs of PSBs, with the angle between the two PSBs in a pair being 90°, and the angle between pairs being 45°. In this paper we concentrate on the Q and U maps at 353 GHz. The Stokes Q and U maps at lower frequencies (100, 143, and 217 GHz) are only used to determine the spectral energy distribution (SED) of the dust emission in polarization.

To quantify systematic effects, we also use maps made from other data subsets (Planck Collaboration VI 2014). We use the ring halves (hereafter “HalfRing”), where the approximately 60 circles performed for each Planck telescope ring (also called a stable pointing period) are divided into two independent subsets of 30 circles. We use observational years (hereafter “Years”), consisting of Surveys 1 and 2 on the one hand and Surveys 3 and 4 on the other, to build two additional maps with independent noise.

The Planck maps we use are in thermodynamic units (K_CMB). To characterize the SED of the dust emission in polarization we express the data as the specific intensity (such as that for Stokes I dust emission, I_d(ν)) at the Planck reference frequencies, using the conversion factors and colour corrections from Planck Collaboration IX (2014)³. For the average dust SED at intermediate Galactic latitudes, the colour correction factor is 1.12 at 353 GHz (see Table 3 in Planck Collaboration Int. XXII 2015).

In addition to these basic products, a Planck CO map from Planck Collaboration XIII (2014), the so-called “type 3” map, and the Planck 857 GHz map, are also used in the selection of the large intermediate latitude analysis regions (see Sect. 3.3.1).

2.2. Emission contributions to the Planck HFI polarization maps

2.2.1. Polarized thermal dust emission

Thermal dust emission is partially linearly polarized (e.g., Hildebrand et al. 1999; Benoît et al. 2004; Ponthieu et al. 2005; Vaillancourt 2007). It is the dominant polarized foreground signal in the high frequency Planck bands (Tucci et al. 2005; Dunkley et al. 2009b; Fraisse et al. 2009; Fauvet et al. 2011; Planck Collaboration Int. XXII 2015).

Dust polarization arises from alignment of non-spherical grains with the interstellar magnetic field (e.g., Hildebrand 1988; Draine 2004; Martin 2007). The structure of the dust polarization sky has already been described using maps of the polarization fraction (p) and angle (ψ) derived from the Planck HFI 353 GHz data (Planck Collaboration Int. XIX 2015; Planck Collaboration Int. XX 2015). The map of p shows structure on all scales, with polarization fractions ranging from low (less than 1%) to high values (greater than 18%). Planck Collaboration Int. XIX (2015) and Planck Collaboration Int. XX (2015) report an anti-correlation between p and the local dispersion of ψ, which indicates that variations in p arise mainly from depolarization associated with changes in the magnetic field orientation within the beam, rather than from changes in the efficiency of grain alignment.

Planck Collaboration Int. XXII (2015) showed that the SED of polarized dust emission over the four Planck HFI frequencies from 100 to 353 GHz is consistent with a modified blackbody emission law of the type I_d(ν) ∝ ν^β_dB_ν(T_d), with spectral index β_d = 1.59 for T_d = 19.6 K ⁴, and where B_ν is the Planck function. About 39% of the sky at intermediate Galactic latitudes was analysed⁵. Among 400 circular patches with 10° radius (equivalent to a sky fraction ) the 1σ dispersion of β_d was 0.17 for constant T_d = 19.6 K. We scale this uncertainty on β_d to larger sky areas by using the factor . This is a conservative choice because this uncertainty includes the effects of noise in the data and so is an upper limit to the true regional variations of β_d on this scale. This polarization spectral index can be compared to variations in the spectral index for the intensity SED. For that quantity the S/N of the data is higher than for polarization and Planck Collaboration Int. XXII (2015) report a dispersion of 0.07 (1σ) over the same sized circular patches. Planck Collaboration Int. XVII (2014) extend this analysis for intensity to high Galactic latitudes in the southern Galactic cap, using the dust-H i correlation to separate the faint emission of dust from the anisotropies of the cosmic infrared background, and find a dispersion of about 0.10 in . We expect spectral variations to be correlated in polarization and intensity, unless the dust emission has a significant component that is unpolarized.

2.2.2. CMB

The CMB temperature anisotropies have been measured with unprecedented accuracy by the Planck Collaboration (Planck Collaboration I 2014; Planck Collaboration XV 2014), and preliminary Planck polarization results have been demonstrated to be in very good agreement with the cosmology inferred from temperature measurements (Planck Collaboration I 2014; Planck Collaboration XVI 2014).

For , the ΛCDM concordance model has been shown to be a very good fit to all the available data (including preliminary Planck results at ℓ ≳ 50; see Barkats et al. 2014, for a recent compendium). For 353 GHz data at small angular scales (ℓ ≳ 400), the E-mode CMB polarization is comparable to the power of dust polarization at high Galactic latitudes.

The CMB B-mode power, even for the highest primordial tensor perturbation models, is negligible with respect to the dust polarization at 353 GHz. Since no reliable published CMB polarization maps are available, we have chosen not to remove the CMB polarization from the Planck HFI Q and U maps. Nevertheless, because the CMB E-mode polarization is significant with respect to the dust at 353 GHz at high multipoles (and even at lower multipoles for the lower frequencies), when studying the Planck HFI bands we subtract the Planck best-fit ΛCDM model (Col. 2 of Table 2 in Planck Collaboration XVII 2014) from the dust power spectra, paying the price of an increased error due to sample variance. No CMB is removed in this work when computing the dust spectra.

2.2.3. Synchrotron emission

Synchrotron emission is known to be significantly polarized (up to 75% for typical relativistic electron spectra, Rybicki & Lightman 1979). Since its specific intensity scaling with frequency follows a power law with a spectral index close to −3 (Gold et al. 2011; Macellari et al. 2011; Fuskeland et al. 2014), synchrotron polarized emission is expected to be subdominant in the Planck HFI channels in general and negligible at 353 GHz (Tucci et al. 2005; Dunkley et al. 2009a; Gold et al. 2011; Fauvet et al. 2011; Fuskeland et al. 2014; Planck Collaboration Int. XXII 2015). Hence, we neither subtract nor mask any synchrotron contribution before estimating the angular power spectra of dust polarization. The justification of this assumption will be demonstrated below by studying the frequency dependence of the polarized dust power between 100 and 353 GHz (but see also Appendix D.4).

2.2.4. Polarized point sources

Radio sources have been shown to have a fractional polarization of a few percent (e.g., Battye et al. 2011; Massardi et al. 2013).Their contribution to the polarization angular power spectra in the Planck HFI bands is expected to be negligible at low and intermediate multipoles (Battye et al. 2011; Tucci & Toffolatti 2012). Upper limits have been set on the polarization of infrared galaxies, and their contribution to the polarization power spectra is also expected to be negligible (e.g., Seiffert et al. 2007). However, the brightest of the polarized point sources can be responsible for ringing in the angular power spectra estimation, and therefore need to be masked (see Sect. 3.3.1).

2.2.5. CO emission

The first three carbon monoxide (CO) Galactic emission lines at 115 GHz (J = 1 → 0), 230 GHz (J = 2 → 1), and 345 GHz (J = 3 → 2), contribute significantly to the power in the Planck HFI bands at 100, 217, and 353 GHz, respectively (Planck Collaboration IX 2014). The Planck data were used to produce the first all-sky maps of Galactic CO emission (Planck Collaboration XIII 2014). It is known that CO emission can be intrinsically polarized (Goldreich & Kylafis 1982; Li & Henning 2011). Furthermore, CO emission can induce spurious polarization because of the differences in spectral transmission at the CO frequencies between Planck HFI detectors (Planck Collaboration IX 2014). For these reasons, we mask CO-emitting regions (Sect. 3.3.1). Outside this mask, as has been shown in Planck Collaboration XIII (2014) by comparing the Planck CO maps to high Galactic latitude ground-based CO observations (Hartmann et al. 1998; Magnani et al. 2000), the CO emission is negligible in the Planck channels (lower than one fourth of each channel noise rms at 95% CL). We will check that our polarization analysis is not contaminated by CO by examining the frequency dependence of the polarized angular power spectra of the dust (see Sect. 4.5).

2.3. Systematics of the Planck HFI polarization maps

The first CMB polarization results from Planck were presented in Planck Collaboration I (2014) and Planck Collaboration XVI (2014). The EE power spectrum at ℓ> 200 was found to be consistent with the cosmological model derived from temperature anisotropies. Systematic effects in the data have so far limited the use of Planck HFI polarization data on large angular scales. The polarization systematics in the 2013 data were discussed and estimated in Planck Collaboration VI (2014; see the power spectra shown in their Fig. 27). The same data were used in the Planck Galactic polarization papers (Planck Collaboration Int. XIX 2015; Planck Collaboration Int. XX 2015; Planck Collaboration Int. XXI 2015; Planck Collaboration Int. XXII 2015) and there the low brightness regions of the 353 GHz sky were masked.

In this paper, we use a new set of Planck polarization maps for which the systematic effects have been significantly reduced. Corrections will be fully described and applied in the Planck 2015 cosmology release, as well as the remaining systematic effects that we describe briefly here.

Two main effects have been corrected in the time-ordered data prior to mapmaking. The correction for the nonlinearity of the analogue-to-digital converters was improved, and we have also corrected the data for very long time constants that were not previously identified (see Planck Collaboration VI 2014; Planck Collaboration X 2014). After these corrections, the main systematic effects over the multipole range 40 <ℓ< 600 relevant to this analysis result from leakage of intensity into polarization maps.Effects arising from polarization angle and polarization efficiency uncertainties have been shown to be second order (Planck Collaboration VI 2014). The leakage effect can be expressed as (1)where ν is the frequency, and , the leakage pattern for bolometer b, is multiplied by the different s leakage source maps and their associated scaling coefficients . The leakage patterns are fully determined by the scanning strategy. They represent the cumulative result of all the systematic effects that lead to a leakage of intensity to polarization. In the Planck HFI bands, there are three main sources of intensity to polarization leakage: (i) monopole differences between detectors not corrected by data destriping (the intensity source term is constant over the sky); (ii) bolometer inter-calibration mismatch (the source term is the full intensity map, including the CMB dipole and the Galactic emission); and (iii) a dust spectral mismatch term. The spectral bandpass varies from one bolometer to another for a given band. As the bolometer gain is calibrated on the CMB dipole, the differential gains on the dust emission produce the bandpass mismatch term (for which the source term is the dust intensity map).

Results from a global fit of the Q and U maps with these three leakage terms are used to quantify the leakage. This fit yields estimates of the scaling coefficients at each frequency ν in Eq. (1), which allows us to compute angular power spectra of the leakage terms in Sect. 4. We point out that the fit captures any emission in the maps that has a pattern on the sky similar to one of the leakage patterns. These global fit maps are used to assess the level of systematics, but are not removed from the data.

An independent and complementary estimate of systematic effects in the data is provided by the null tests that we can build from the DetSets, HalfRings, and Years data subsets (see Sect. 2.1). These null tests are a good way of determining the level of any systematics other than intensity to polarization leakage. These are pursued in Sect. 4.1 and Appendices C.1 and D.3.

3. Computation of angular power spectra of polarized dust emission

3.1. Methods

We can use the Planck data to compute the polarization angular power spectra ( and ) of the polarized dust emission within selected sky regions. Even if the statistical properties of the dust emission on the sky might not be entirely captured by a two-point function estimator, because the scope of this paper is to assess the level of dust in the framework of CMB data analysis, we follow the approximation that is generally made when processing such data, i.e., that a large fraction of the information is contained in power spectra.

On an incomplete sky, a polarization field can be divided into three classes of modes: pure E-modes; pure B-modes; and “ambiguous” modes, which are a mixture of the true E- and B-modes (Bunn et al. 2003).

The ambiguous modes represent a cross-talk between E and B, which is often referred to as “E-to-B leakage” for the CMB (because for the CMB ). Methods used to estimate the CMB angular power spectrum for polarization account for and correct analytically for the incomplete sky coverage. However, the presence of the ambiguous modes yields a biased estimate of the variance of the spectra, unless so-called “pure” power spectrum estimators are used (Smith 2006). For dust polarization, the power in E- and B-modes are comparable, and we do not expect a significant variance bias.

The two specific approaches we use are Xpol, our main method, which we describe as a “classical” pseudo-C_ℓ estimator⁶ and, for comparison, Xpure, a pure pseudo-C_ℓ estimator. Since they give similar results for the present study, as we demonstrate below, the former is chosen as our main method because it is computationally less expensive.

3.1.1. Xpol

Xpol is an extension of the Xspect method (Tristram et al. 2005) to polarization. Xspect computes the pseudo power spectra and corrects them for incomplete sky coverage, filtering effects, and pixel and beam window functions. Correction for incomplete sky coverage is performed using a Master-like algorithm (Hivon et al. 2002), consisting of the inversion of the mode-mode coupling matrix M_ℓℓ′ that describes the effect of the partial sky coverage; M_ℓℓ′ is computed directly from the power spectrum of the mask that selects the data in the analysis region of interest (Sect. 3.3). The HEALPix pixel window functions (Górski et al. 2005) are used to correct for pixelization effects. For the beams we use the Planck HFI individual detector beam transfer functions described in Planck Collaboration VII (2014).

Xpol estimates the error bars analytically, without requiring Monte Carlo simulations. Using simulated data comprising inhomogeneous white noise and a Gaussian map with a dust power spectrum, we have checked that this analytical estimate is not biased for multipoles ℓ> 40. The analytical error bars combine the contributions from instrumental noise and sample variance. The Gaussian approximation of the sample variance is (2)where X = { E,B }, f_sky is the retained sky fraction (which can be if the sky field is apodized), and Δℓ_bin is the size of the multipole bin ℓ_bin. To estimate the error relevant to the polarized dust signal measured within a given region, we subtract quadratically this estimate of the contribution from sample variance. The performance of Xpol on Gaussian simulations that have EE and BB dust-like angular power spectra is presented in Appendix A.

3.1.2. Xpure

Xpure is a numerical implementation of the pure pseudo-spectral approach described and validated in Grain et al. (2009). The method is optimized for computing CMB B-mode power spectra over small sky patches. It uses a suitably chosen sky apodization that vanishes (along with its first derivative) at the edges of the patch in order to minimize the effects of E-to-B leakage. For the estimation of the angular power spectra of the Planck data we compute the cross-correlation of two different DetSets. Uncertainties are obtained by performing Monte Carlo inhomogeneous white noise simulations, considering the diagonal terms of the pixel-pixel covariance for the two data sets.

We have applied this algorithm to Planck 353 GHz maps and simulations to estimate the B-mode power spectrum of the Galactic thermal dust emission in regions ranging from 1% to 30% of the sky. We have used Xpure here as a cross-check for the robustness of the Xpol method presented in the previous section. Validation and comparison of the performance of the two algorithms is presented in Appendix A on simulations and in Appendix D.1 on application to the data.

3.2. Computing cross-spectra

To avoid a bias arising from the noise, we compute all the Planck power spectra from cross-correlations of DetSets maps (see Sect. 2). The noise independence of the two DetSet maps at a given frequency was quantified in Planck Collaboration XV (2014) and the resulting level of noise bias in the cross-power spectra between them has been shown to be negligible. The cross-power spectrum at a given frequency ν is (3)where and are the two independent DetSet 1 and DetSet 2 maps at the frequency ν. The Planck cross-band spectrum between the frequencies ν and ν′ is (4)or equivalently the cross-spectrum between ν and ν′ of the averaged frequency maps .

As we noted in Sect. 2.2.2, from each computed spectrum we subtract the Planck best-fit ΛCDM model (Planck Collaboration XVI 2014), i.e., the theoretical model obtained from the temperature data fit, while spectra are kept unaltered.

3.3. Selection of regions

To measure the dust polarization power spectra with high S/N, we select six large regions, the analysis regions of interest at intermediate Galactic latitude, which have effective coverage of the sky from 24 to 72% (see Sect. 3.3.1)⁷. For statistical studies at high Galactic latitude, we compute spectra on a complete set of smaller regions or patches (Sect. 3.3.2), similar in size to the patches observed in typical CMB experiments.

3.3.1. Large regions

For selection of all of the large regions, we used the Planck CO map from Planck Collaboration XIII (2014), smoothed to a 5° resolution, to mask the sky wherever the CO line brightness I_CO ≥ 0.4 K km s^-1 8. This mask is shown in Fig. 1. The complement to this mask by itself defines a preliminary region that retains a sky fraction f_sky = 0.8.

We then mask the sky above successively lower thresholds of I₈₅₇ in the Planck 857 GHz intensity map, smoothed to a 5° resolution, chosen such that together with the CO mask we select five more preliminary regions that retain f_sky from 0.7 to 0.3 in steps of 0.1. These six regions are shown in Fig. 1.

To avoid power leakage, these six masks are then apodized by convolving with a 5° FWHM Gaussian which alters the window function by gradually reducing the signal towards the edges of the retained regions and thus lowers the effective retained sky coverage. The value is simply defined as the mean sky coverage of the window function map.

Finally, we mask data within a radius 2σ_beam of point sources selected from the Planck Catalogue of Compact Sources (PCCS, Planck Collaboration XXVIII 2014) at 353 GHz.The selected sources have S/N> 7 and a flux density above 400 mJy⁹. Selection of spurious infrared sources in bright dust-emitting regions is avoided by using contamination indicators of infrared cirrus listed in the PCCS description (Planck Collaboration XXVIII 2014). This point-source masking is done to prevent the brightest polarized sources from producing ringing in the power spectrum estimation, while avoiding the removal of dust emitting regions and their statistical contribution to the angular power spectra. The details of this source selection will be presented in the Planck 2014 release papers. The edges of the masks around point sources were apodized with a 30′ FWHM Gaussian, further reducing the retained net effective sky coverage.

In combination these masking and apodization procedures result in six large retained (LR) regions, which we label using the percentage of the sky retained (the net effective fractional sky coverages, , are listed in Table 1), e.g., LR72 for the largest region and LR24 for the smallest.

Table 1 also lists other properties of the regions, including ⟨I₃₅₃⟩, the mean specific intensity at 353 GHz within the region in MJy sr^-1, and N_{H i}, the mean H i column density in units of 10²⁰ cm^-2 computed on the LAB H i survey data cube (Kalberla et al. 2005).

3.3.2. Small patches at high Galactic latitude

To examine the statistics of the angular power spectra at high Galactic latitude, we also analyse the Planck polarization maps within patches with a size similar to those of typical ground-based and balloon-borne CMB experiments (e.g., QUIET Collaboration et al. 2012; Hanson et al. 2013; BICEP2 Collaboration 2014b; Ade et al. 2014; Naess et al. 2014). Specifically, we consider 400 deg² circular areas (radius 11°̣3) centred on the central pixel positions of the HEALPix N_side = 8 grid that have Galactic latitude | b | > 35°. This results in 352 such patches. These are apodized with a 2° FWHM Gaussian, which reduces the retained sky fraction to for each patch. For these small patches, we do not mask point sources as was done in selecting the large regions, since we want to preserve the same for each mask, nor is the CO mask needed, since we use these masks only on the 353 GHz high Galactic latitude data. We note that for N_side = 8, pixels have an area of about 54 deg² and characteristic centre-to-centre spacing of about 7°̣4. Therefore, this HEALPix grid oversamples the sky relative to this patch size, so that the patches overlap and are not independent.

4. Dust polarized angular power spectra at intermediate Galactic latitude

In this section, we quantify the 353 GHz dust polarization in the power spectrum domain, achieving a high S/N through the use of LR regions. For convenience we present results for and , where .

4.1. Description of the spectra

Using Xpol we have computed and from the two DetSets at 353 GHz, as a function of the multipole ℓ in the range 40–600¹⁰. These represent the first measurements of the thermal dust and power spectra on large fractions of the sky for ℓ> 40 (see Ponthieu et al. 2005; Gold et al. 2011, for earlier related studies).

In Fig. 2 we present the results for f_sky = { 0.3,0.5,0.7 }. The amplitudes of the spectra increase with increasing f_sky because the polarized emission is brighter on average when more sky is retained. We point out that the LR regions presented in Sect. 3.3.1 overlap because they form a nested set. However, the power spectra derived from them are almost independent of each other because each spectrum is dominated by the brightest areas in the corresponding region, namely the parts closest to the Galactic plane and hence the areas in which each LR region differs from those nested inside it.

The and spectra are characterized by a power-law dependence on multipole ℓ over the range ℓ = 40 to 600; furthermore, the slope is similar for different regions that retain from 24 to 63% of the sky.

Figure 2 also shows the power spectrum computed from the Planck 2013 best-fit ΛCDM model of the CMB temperature data (Planck Collaboration XVI 2014), and the much lower expectation for from the CMB model with primordial gravitational waves with amplitude r = 0.2. At 353 GHz, the angular power spectra of the dust are about 3−4 orders of magnitude larger than the CMB model at ℓ = 30, 1−2 orders of magnitude larger at ℓ = 100, and about the same order of magnitude as the CMB at ℓ> 300. At 353 GHz, the angular power spectra for dust are much greater than the CMB model power spectrum for all ℓ values in Fig. 2. The dust power spectra are larger than the r = 0.2 CMB spectrum by 4−5 orders of magnitude at ℓ = 30, and by 3−4 orders of magnitude at ℓ = 100. At ℓ = 500, where the lensing of CMB anisotropies is the dominant contribution to the CMB model spectrum, the dust is still 2−3 orders of magnitude higher.

As discussed in Appendix A, we do not expect any significant bias, or E-to-B leakage, from the computation of the dust angular power spectra using Xpol. Figure 2 also includes the and spectra at 353 GHz, computed from our estimate of the leakage terms from intensity to polarization (discussed in Sect. 2.3), for the same three LR regions. The dust spectra are much higher than the corresponding spectra for the leakage, which represent the main systematic effects over the ℓ range of interest for this work. The largest contamination of the dust signal by leakage is about 3.5% at ℓ = 50, for both the and spectra. Because we consider that our estimate of the leakage maps is conservative, we conclude that contamination of the dust and power spectra by systematic effects amounts to a maximum of 4% at ℓ = 50, and is less at higher multipoles. Therefore, we have not corrected the Planck data for intensity-to-polarization leakage in this work.

Finally, in Fig. 2 we present for the three LR regions the absolute value of the null-test spectra anticipated in Sect. 2.3, here computed from the cross-spectra of the HalfRing/DetSet differences, i.e., (353_DS1,HR1−353_DS1,HR2) /2 × (353_DS2,HR1 − 353_DS2,HR2)/2. These and spectra show a behaviour that is close to what is expected from a white-noise dominated (thus ℓ²) spectrum. The amplitudes of these error estimates are consistent with the noise expectations; there is no evidence for any effects of systematics.

For completeness, in Appendix B we present a further quantification of the power spectrum of thermal dust emission at 353 GHz via the spectra involving temperature and polarization, and , and the polarization cross-spectrum, .

4.2. Power-law fit

To assess the apparent power-law dependence quantitatively, we made a χ² fit to the spectra at 353 GHz using the form , where X ∈ { E,B }. For and we fit the 22 band-powers in the range 60 <ℓ< 500. For both power spectra we restricted the fit to ℓ> 60 to avoid a possible bias from systematic effects in the data, and also because the angular power spectra of the dust polarization exhibit more spatial variation than is expected for a Gaussian random field, particularly on large angular scales.

The exponents of the power-law fits, α_XX, are plotted in Fig. 3 for each of the six LR regions identified with f_sky. All exponents are consistent with constant values of α_EE = −2.41 ± 0.02 and of α_BB = −2.45 ± 0.03. While there is a slight indication of a steeper slope for than for , hereafter we adopt the mean exponent −2.42 ± 0.02. This exponent is consistent with the value α_TT fitted to the Planck 353 GHz dust intensity power spectra in this range of ℓ (Planck Collaboration XV 2014; Planck Collaboration Int. XXII 2015) on similar-sized regions at intermediate Galactic latitude, but slightly flatter than the α_TT fitted at higher ℓ and higher frequency (Miville-Deschênes et al. 2007, 2010; Planck Collaboration XV 2014).

For the fits with fixed exponent α_EE,BB = −2.42, the values of the χ² (with number of degrees of freedom, N_d.o.f. = 21) are listed in Table 1 for the six LR regions. For the spectra, the χ² values range from 26.3 (probability to exceed, PTE = 0.2) to 44.8 (PTE = 0.002), with a trend for a quality of fit that degrades with increasing f_sky. For the spectra, the χ² values range from 14.0 (PTE = 0.87) to 22.1 (PTE = 0.39), with a trend related to f_sky.

Possible explanations for this difference between and spectrum shape descriptions, are the chance correlation between dust and CMB polarization, which is not taken into account in the subtraction of the CMB spectrum, and the increasing S/N degrading the overall quality of the fit when going from to and for from f_sky = 0.3 to f_sky = 0.8 as the amplitude of the dust polarized signal increases. We stress that while the power law in ℓ is a good general description of the shape of the dust polarized power spectra, a full characterization would have to consider the detailed features that can be seen in Fig. 2.

4.3. Amplitude dependence on ⟨I₃₅₃⟩

Similarly to what has been done for intensity power spectra in, e.g., Gautier et al. (1992) and Miville-Deschênes et al. (2007), we investigate how the amplitude of the dust polarization power spectrum scales with the dust intensity. To quantify the dust emission at 353 GHz, we use the model map derived from a modified blackbody fit to the Planck data at ν ≥ 353 GHz and IRAS at λ = 100 μm, presented in Planck Collaboration XI (2014). This map is corrected for zodiacal emission and the brightest extragalactic point sources are subtracted. The Galactic reference offsets of the underlying IRAS and Planck data were obtained through a method based on correlation with 21 cm data from the LAB H i survey (Kalberla et al. 2005) integrated in velocity, effectively removing the CIB monopole in the model map. The mean dust intensity, ⟨I₃₅₃⟩, listed for each LR region in Table 1, ranges from 0.068 to 0.227 MJy sr^-1 for increasing . The mean column density calculated from the LAB survey data is also listed in Table 1, with values ranging from 1.65 × 10²⁰ cm^-2 to 6.05 × 10²⁰ cm^-2.

For all of the LR regions, we fit the and spectra with a power law in ℓ, using the fixed exponent α_EE,BB = − 2.42, over the ℓ-ranges defined in Sect. 4.2. The amplitudes A^EE derived from these fits are listed in Table 1; the A^BB amplitudes can be retrieved from the A^BB/A^EE ratio. These amplitudes are plotted as a function of ⟨I₃₅₃⟩ in Fig. 4 after normalization by the maximum value found for the largest region (LR72). We fit the empirical dependence of these amplitudes on ⟨I₃₅₃⟩ as a power law of the form A^XX(⟨I₃₅₃⟩)= K_XX⟨I₃₅₃⟩ ^ϵ_XX where X ∈ {E,B}. The two fitted exponents are quite similar, ϵ_EE = 1.88 ± 0.02 and ϵ_BB = 1.90 ± 0.02. The exponent that we find for polarization is close to the one observed in the diffuse interstellar medium for the dust intensity, consistent with , where ⟨I_ν⟩ is the mean value of the dust specific intensity (Miville-Deschênes et al. 2007). Values close to 2 are expected, because we compute angular power spectra, which deal with squared quantities.

Although the data points roughly follow this ⟨I₃₅₃⟩^1.9 dependence, the empirical law fails to fully describe individual dust amplitudes (e.g., the estimate is off by about 20% for on LR33). The scaling can help to asses the order of magnitude of the dust polarization level on a specific region, but is not a substitute for actually characterizing the polarized angular power spectra.

4.4. Amplitude of relative to

We examine the ratio of the amplitudes of the fitted power laws found in Sect. 4.3. The A^BB/A^EE ratios are listed in Table 1, and are plotted for different values of f_sky in Fig. 5. For all of the LR regions, we observe more power in the dust spectrum than in . All ratios are consistent with a value of A^BB/A^EE = 0.52 ± 0.03, significantly different from unity, over various large fractions of the intermediate latitude sky.

This result is not taken into account in existing models of polarized microwave dust emission that have been developed to test component separation methods. For example, we have computed the and spectra over the LR regions for the Planck Sky Model (Delabrouille et al. 2013) and the model of O’Dea et al. (2012); for both models and all LR regions we find a ratio A^BB/A^EE close to 1. However, these two models are based on a very simplified picture of the Galactic magnetic field geometry and assumptions on how the polarized emission depends on it. Further insight into the structure of the dust polarization sky is required to account for the observed ratio.

4.5. Amplitude dependence on frequency

Finally, we explore the frequency dependence of the amplitude of the angular power spectra. We compute the and angular power spectra from the Q and U DetSet maps at 100, 143, 217, and 353 GHz (see Sect. 3.2). From these four sets of polarization maps, we compute ten power spectra: 100 × 100; 100 × 143; 100 × 217; 100 × 353; 143 × 143; 143 × 217; 143 × 353; 217 × 217; 217 × 353; and 353 × 353.

The ten angular cross-power spectra are consistent with a power law in ℓ, with the exponent α_EE,BB = −2.42 measured at 353 GHz (Sect. 4.2). Therefore, to each of these spectra we fit the amplitudes of a power-law function that has a fixed exponent α_EE,BB = − 2.42, in the range 40 <ℓ< 500, for and . As an illustration of the quality of the fit, for the smallest region (LR24, f_sky = 0.3) the averages and dispersions of the χ² (21 degrees of freedom) of the fits are and for the ten cross-frequency spectra.

To compare the frequency dependence of the results of the fits to that expected from the SED for dust polarization from Planck Collaboration Int. XXII (2015), we converted the fitted amplitudes A^EE,BB from μK to units of (MJy sr^-1)², taking into account the Planck colour corrections¹¹. For all regions, we examined the frequency dependence by plotting the amplitudes normalized to unity at 353 GHz versus the effective frequency¹². A representative example is shown in Fig. 6 for the smallest region, LR24 (f_sky = 0.3).

For all of the LR regions the frequency dependence found is in good agreement with the square of the adopted dust SED, which is a modified blackbody spectrum having β_d = 1.59 and T_d = 19.6 K (Planck Collaboration Int. XXII 2015). We note that no fit was performed at this stage and that the square of the adopted dust SED goes through the 353 GHz data point. However, if we fit the amplitude of the dust frequency dependence, with fixed β_d and T_d, the χ² (N_d.o.f. = 9) is 13.1 for EE (PTE = 0.16) and 3.2 for BB (PTE = 0.96). This good agreement supports the assumption made in the present work about the faintness of the synchrotron and CO emission at high Galactic latitude (see Sect. 2.2). The residuals to the fit do not show any evidence for excess power at 100 GHz, which might arise from polarized synchrotron emission. This result is consistent with the study of synchrotron polarization at high Galactic latitudes by Fuskeland et al. (2014), and confirmed in Appendix D.4. Furthermore, we do not see any excess at either 100 or 217 GHz, such as could arise from leakage and/or polarization associated with CO line emission in these bands (Sect. 2.2.5).

5. Statistical study of the dust power spectra at high Galactic latitudes

In Sect. 4, statistical properties of dust polarization were derived from angular power spectra computed with the Planck353 GHz data on large fractions of the sky ( from 24% to 72%). Most of the CMB experiments target fields at high Galactic latitude with a smaller size than this. So now we evaluate the dust B-mode power in such patches.

We perform a statistical analysis computing the and power spectra of the 353 GHz Planck polarization data for high Galactic latitude circular patches of size 400 deg², similar to the size of the fields observed in ground-based and balloon-borne CMB polarization experiments (see Sect. 3.3.2). We verify that the empirical scaling law of the power spectra amplitudes as a function of the mean intensity, derived in Sect. 4.3, holds on small independent fractions of the sky that are not nested.

5.1. Data processing

We compute the dust and angular power spectra at 353 GHz from the cross-spectra of the two DetSets with independent noise in this channel, using Xpol on all of the 352 high Galactic latitude (| b | > 35°) patches defined in Sect. 3.3.2. As shown in Appendix A, we do not expect any significant cut-sky leakage from E polarization to B polarization for the dust in such patches. The and spectra are computed in the range 40 <ℓ< 370, using top-hat binning in the intervals defined between multipoles 40, 70, 110, 160, 220, 290, and 370.

To each of the 352 and spectra we fit the power law in ℓ presented in Sect. 4.2 over the range 40 <ℓ< 370 with a slope fixed to α_EE,BB = −2.42 in order to compute the amplitudes A^EE,BB and their associated errors. Since we are interested in quantifying the possible contamination by dust polarization for ground-based CMB experiments, here we express the amplitude of the fitted power laws at the position of the recombination bump at ℓ = 80 (). In Appendix C.1, using a null test, we show that the potential systematic effects in the data do not affect the computation of the and spectrum amplitudes.

5.2. Results

For patches covering 1% or more of the sky, we do not observe any significant departure in the shapes of the and spectra from the results presented in Sect. 4.2, even if a large dispersion of the power-law exponent (attributable to the low S/N for some patches) is observed. The / ratios are on average consistent with what has been derived on the larger LR regions in Sects. 4.4. These ratios are presented in Appendix C.2.

The fitted dust and amplitudes (A^EE and A^BB) at ℓ = 80 are presented in Fig. 7 in units of μK at 353 GHz and as a function of the mean dust intensity ⟨I₃₅₃⟩ of each patch.

We see a clear correlation of the dust and amplitudes with ⟨I₃₅₃⟩. The data points are consistent with the same ⟨I₃₅₃⟩ dependence of the scaling law found in Sect. 4.3, A^EE,BB ∝⟨I₃₅₃⟩^1.9, indicating that this empirical law also applies in the faintest regions of the sky and holds reasonably well for the description of the amplitudes for any approximately 1% patch of the sky at latitudes above 35°.

The amplitudes of the two empirical scaling laws are adjusted to the data points by computing the median over all the patches of A^EE,BB/⟨I₃₅₃⟩^1.9. The amplitudes that we derive from Fig. 7 are in agreement with the amplitudes reported for the LR regions in Sect. 4.3 and Table 1.

The dust brightness ⟨I₃₅₃⟩ in the cleanest 400 deg² patches (⟨I₃₅₃⟩ down to 0.038 MJy sr^-1) is about a factor of 2 lower than the value listed in Table 1 for the LR24 region (0.068 MJy sr^-1). Applying the empirical scaling derived in Sect. 4.3, the expected level of the dust polarization angular power spectrum is a factor of (0.068/0.038)^1.9 = 3 higher on 24% of the sky than on the faintest patches.

A dispersion around this empirical power law for the patches is observed. The shaded regions in Fig. 7 show the expected scatter in this relation coming from the noise, together with Gaussian sample variance. The former is computed from Planck inhomogeneous noise Monte Carlo simulations presented in Appendix C.1 and the latter is computed given the and the binning used to determine the spectra (see Sect. 5.1). Thus, the general trend of the dust and amplitudes among these patches follows the empirical scaling in ⟨I₃₅₃⟩, but the inhomogeneous nature of the polarized dust emission, including variations in magnetic field orientation and grain alignment, is responsible for the large dispersion around it, larger than expected for a stationary Gaussian process plus instrumental noise.

Even at low ⟨I₃₅₃⟩, the statistical error on the amplitude at ℓ = 80 for such small patches is at least 7.5 μK (3σ).

5.3. Optimizing the search for a primordial B-mode signal at the recombination peak

In Sect. 5.1 the 353 GHz A^BB amplitudes were obtained for the 352 patches of 400 deg² at Galactic latitudes | b | > 35°. These amplitudes are expressed at ℓ = 80, the approximate position of the maximum of the CMB tensor B-mode recombination peak. In order to extrapolate these amplitudes to 150 GHz, an observing frequency typical of CMB experiments searching for B-modes, we use the dust SED of Planck Collaboration Int. XXII (2015) presented in Sect. 2.2.1, which is a modified blackbody spectrum with β_d = 1.59 ± 0.17 and T_d = 19.6 K. The conversion factor from the Planck 353 GHz band to the 150 GHz central frequency for this extrapolation, taking into account the Planck colour corrections and the statistical uncertainty of β_d, is in the SED, and this squared in .

These extrapolated estimates are divided by the value of the r = 1 primordial tensor CMB spectrum at ℓ = 80, 6.71 × 10^-2 μK, to express the estimated power in units that we denote r_d. Because the CMB primordial tensor B-mode power scales linearly with r ¹³, a value r_d = 0.1 would mean that the expected contamination from dust at ℓ = 80 is equal to the amplitude of the primordial tensor CMB for r = 0.1. For each of these estimates we also compute σ(r_d), the quadratic sum of the fit errors on A^BB and the above uncertainty from the extrapolation to 150 GHz. We note that the fitted amplitudes A^BB for five of these patches are negative¹⁴, but are consistent with r_d = 0 at 1σ. The r_d values vary from −0.17 to more than 10 and σ(r_d) ranges from as small as 5.6 × 10^-2 to larger than 10. Taking the smallest value of σ(r_d) we see that Planck measurements of the dust spectra for such small patches have, at best, a statistical uncertainty of 0.17 (3 σ) in r_d units.

To reveal the spatial dependence over the high latitude sky, a map of the absolute value of r_d is shown in Fig. 8. Each computed value of r_d is given in a pixel at the position of the centre of each of the 352 patches defined in Sect. 3.3.2 on an N_side = 8HEALPix map. The accompanying map of σ(r_d) is also presented.

We can see in Fig. 8 that there is a high latitude region in the southern Galactic hemisphere for which r_d is quite low. This region is also associated with a small estimated uncertainty. For example, the six HEALPix N_side = 8 pixels numbered 741, 742, 754, 755, 762, and 763 (in the HEALPix “ring” ordering scheme) have r_d = 0.053 ± 0.096, 0.027 ± 0.098, −0.062 ± 0.052, −0.020 ± 0.127, 0.057 ± 0.122, and −0.031 ± 0.121, respectively. These pixels are located around Galactic coordinates l = −30°, b = −70°.

We stress that the expectation for this low level of dust contamination is valid only for these particular patches, including their positions, sizes, shapes, and apodizations. In addition, because we found the amplitudes of the dust spectra associated with these patches based on a power-law fit, our estimate does not take into account possible features in the power spectra that might alter the precise value of dust contamination. Nevertheless, there are clearly some patches that appear to be “optimal," i.e., cleaner than the others. It needs to be emphasized, however, that finding the cleanest areas of the polarized sky for primordial B-mode searches cannot be accomplished accurately using the Planck total intensity maps alone.

6. BB angular power spectrum of dust in the BICEP2 field

In this section, we use the Planck data and the results presented above to assess the dust polarization in the field observed by the BICEP2 experiment (BICEP2 Collaboration 2014a,b). As above, finding at 150 GHz involves two steps, measuring the power spectrum at 353 GHz and then extrapolating the amplitude to 150 GHz.

6.1. An approximation to the observed BICEP2 field

To define a field similar to and representative of the actual BICEP2 field, we carried out the following steps (at the HEALPix N_side = 2048 pixelization): (i) we constructed a mask M by filling the inside of the BICEP2 deep-field outline¹⁵ with 1 and the outside with 0; (ii) we took the complement of this mask, M′ = 1−M; (iii) inside M′ we computed the distance to the border using the HEALPix “iprocess_mask” procedure; (iv) we smoothed this “distance map” with a 7°̣5 FWHM Gaussian; (v) we apodized, computing the 7°̣5 FWHM Gaussian weight from the distance map to obtain M′′; and (vi) we took the complement of M′′ to be our definition of the BICEP2 field (i.e., M_B2 = 1−M′′). The resulting field, which in the rest of this paper we will refer to as M_B2, has (689 deg²) and ⟨I₃₅₃⟩= 0.060 MJy sr^-1. Its is larger than the BICEP2 373 deg² deep-field region, but M_B2 is similar to the BICEP2 inverse noise variance map presented in BICEP2 Collaboration (2014a), even if it extends further in declination but less in right ascension. In Appendix D.2 we show that the main results in this section do not depend significantly on the definition of this field. However, we note that the various applications of filtering in the BICEP2 pipeline mean that their observed spatial modes are a subset of those present within their nominal observed region.

6.2. Statistical estimate of the dust B-mode level in the BICEP2 field

We have seen in Sect. 5.2 that the dust spectra amplitudes statistically follow an empirical scaling law of the mean dust intensity ⟨I₃₅₃⟩ for the patch on which they are computed. We can use this scaling law to assess the most probable dust value for M_B2. For this field, which has a mean dust intensity of ⟨I₃₅₃⟩ = 0.060 MJy sr^-1, the expected value at 353 GHz is (13.4 ± 0.26) μK at ℓ = 80, taking into account the uncertainty on the fitted ⟨I₃₅₃⟩^1.9 amplitude from Sect. 4.3.

In order to extrapolate this value to the BICEP2 observing frequency, 150 GHz, we proceed as in Sect. 5.3, using the same typical dust SED from Planck Collaboration Int. XXII (2015). For patches the size of M_B2 the expected dispersion of β_d (Sect. 2.2.1) is 0.11, which introduces an uncertainty in the extrapolation. The extrapolation, unlike in that section, is made to 150 GHz taking into account the BICEP2 bandpass¹⁶. The conversion factor from the Planck 353 GHz band to the 150 GHz BICEP2 band, taking into account both Planck and BICEP2 colour corrections, as well as the statistical error on β_d, is thus in the SED, and therefore this value squared in .

The resulting expected dust amplitude at 150 GHz is  μK at ℓ = 80. The uncertainty quoted here applies to the most probable value for patches with a similar ⟨I₃₅₃⟩, but it ignores the intrinsic dispersion around this value (see Sect. 5.2). Therefore, a direct measurement of the dust spectrum in M_B2 is required to complete the assessment of its polarized dust level.

In Sect. 5.3, we have localized the amplitudes from the analysis of patches on a map. The pixels falling inside the approximate BICEP2 deep-field region shown in Fig. 8 give a mean value of r_d = 0.207, i.e., an expected dust power of = 1.39 × 10^-2 μK at 150 GHz and ℓ = 80; we note that in this region the Planck estimate of the dust contamination is significantly higher, by a factor of about 2, than it is in some more optimal patches. Nevertheless, a more accurate estimate of the dust polarization amplitude in the BICEP2 observed field requires the direct computation of the angular power spectrum on a similar region.

6.3. BB angular power spectrum of dust in the BICEP2 field

We now compute the spectrum of the dust at 150 GHz on M_B2 presented in Sect. 6.1, using the Planck 353 GHz data. For this purpose, we compute the 353 GHz spectrum by cross-correlating the two DetSets that have independent noise at this frequency in the three ℓ-bins defined by the intervals between multipoles 40, 120, 250, and 400. The ℓ-bin sizes were increased with respect to the previous sections in order to increase the S/N, especially around the range of the CMB recombination bump. Appendix D.1 confirms that the result does not depend on the method of computing the power spectrum.

This power spectrum is extrapolated to 150 GHz as in Sect. 6.2, with an extrapolation uncertainty estimated from the inferred dispersion of β_d. Our final estimate of the spectrum is presented in Fig. 9, together with its 1σ error budget. For the first bin, ℓ = 40−120, the expected level of dust polarized , as extrapolated to 150 GHz, is 1.32 × 10^-2 μK (Fig. 9). The statistical error, estimated from Monte Carlo simulations of inhomogeneous Planck noise (presented in Appendix A for this particular binning), is ± 0.29 × 10^-2 μK, so that the dust spectrum is statistically detected at 4.5σ in this broad ℓ bin.

In order to assess the potential contribution from systematics, we have computed the dust spectrum on M_B2 on different subsets of the data and performed null tests, which are presented in Appendix D.3. In this lowest bin of ℓ, we do not observe any departure from what is allowed by noise. Nevertheless, we stress that below the noise level our cross-spectra could be subject to a positive or negative bias due to systematic effects. For example, if instead of taking the DetSets cross-spectra (as we have done throughout this paper) we take the mean value computed from the DetSets, HalfRings, and Years cross-spectra (presented in Appendix D.3), the statistical significance of our measurement is decreased from 4.5σ to 3.6σ.

The uncertainty coming from the M_B2 definition (presented in Appendix D.2) is 0.04 × 10^-2 μK for this bin, thus much less than the statistical error. For this reason, it is not added to the error budget. However, the spectral extrapolation to 150 GHz adds an additional uncertainty (+ 0.28,−0.24) × 10^-2 μK to the estimated power in M_B2, added in quadrature in Fig. 9.

The expected value in this lowest-ℓ bin from direct computation of the power spectrum on M_B2, as shown in Fig. 9, is lower than (but consistent with) the statistical expectation from the analysis of the 352 high Galactic latitude patches presented in Sects. 5.2 and 6.2. This indicates that M_B2 is not one of the outliers of Fig. 7 and therefore its dust B-mode power is well represented by its mean dust intensity through the empirical scaling law ⟨I₃₅₃⟩^1.9.

These values of the amplitude in the ℓ range of the primordial recombination bump are of the same magnitude as those reported by BICEP2 Collaboration (2014b). Our results emphasize the need for a dedicated joint analysis of the B-mode polarization in this region incorporating all pertinent observational details of the Planck and BICEP2 data sets, which is in progress.

6.4. Frequency dependence

We complement the power spectrum analysis of the 353 GHz map with Planck data at lower frequencies. As in the analysis in Sect. 4.5, we compute the frequency dependence of the BB power measured by Planck at HFI frequencies in the BICEP2 field, using the patch M_B2 as defined in Sect. 6.1.

We compute on M_B2 the Planck auto- and cross-power spectra from the three Planck HFI bands at 100, 143, 217, and 353 GHz, using the two DetSets with independent noise at each frequency, resulting in ten angular power spectra (100 × 100, 100 × 143, 100 × 217, 100 × 353, 143 × 143, 143 × 217, 143 × 353, 217 × 217, 217 × 353, and 353 × 353), constructed by combining the cross-spectra as presented in Sect. 3.2. We use the same multipole binning as in Sect. 6.3. To each of these spectra, we fit the amplitude of a power law in ℓ with a fixed exponent α_BB = −0.42 (see Sect. 4.2). In Fig. 10 we plot these amplitudes as a function of the effective frequency from 143 to 353 GHz, in units of sky brightness squared, as in Sect. 4.5. Data points at effective frequencies below 143 GHz are not presented, because the dust polarization is not detected at these frequencies. An upper limit on the synchrotron contribution at 150 GHz from the Planck LFI data is given in Appendix D.4.

We can see that the frequency dependence of the amplitudes of the Planck HFI spectra is in very good agreement with a squared dust modified blackbody spectrum having β_d = 1.59 and T_d = 19.6 K (Planck Collaboration Int. XXII 2015). We note that this emission model was normalized only to the 353 GHz point and that no global fit has been performed. Nevertheless, the χ² value from the amplitudes relative to this model is 4.56 (N_d.o.f. = 7). This shows that dust dominates in the specific M_B2 region defined where these cross-spectra have been computed. This result emphasizes the need for a dedicated joint Planck-BICEP2 analysis.

7. Conclusions

We have presented the first nearly all-sky statistical analysis of the polarized emission from interstellar dust, focussing mostly on the characterization of this emission as a foreground contaminant at frequencies above 100 GHz. Our quantitative analysis of the angular dependence of the dust polarization relies on measurements at 353 GHz of the and (alternatively and ) angular power spectra for multipoles 40 <ℓ< 500. At this frequency only two polarized components are present: dust emission; and the CMB, which is subdominant in this multipole range. We have found that the statistical, spatial, and spectral distribution properties can be represented accurately by a simple model over most of the sky, and for all frequencies at which Planck HFI measures polarization.

• The angular power spectra and at 353 GHz are well fit by power laws in ℓ with exponents consistent with α_EE,BB = −2.42 ± 0.02, for sky fractions ranging from 24% to 72% for the LR regions used.

• The amplitudes of and in the LR regions vary with mean dust intensity at 353 GHz, ⟨I₃₅₃⟩, roughly as ⟨I₃₅₃⟩^1.9.

• The frequency dependence of the dust and from 353 GHz down to 100 GHz, obtained after removal of the prediction from the Planck best-fit CMB model (Planck Collaboration XVI 2014), is accurately described by the modified blackbody dust emission law derived in Planck Collaboration Int. XXII (2015), with β_d = 1.59 and T_d=19.6 K.

• The ratio between the amplitudes of the two polarization power spectra is , which is not consistent with the simplest theoretical models.

• Dust and spectra computed for 352 high Galactic latitude 400 deg² patches satisfy the above general properties at 353 GHz and have the same frequency dependence.

We have shown that the Planck determination of the 353 GHz dust polarization properties is unaffected by systematic errors for ℓ> 40. This enables us to draw the following conclusions relevant for CMB polarization experiments aimed at the detection of primordial CMB tensor B-modes.

• Extrapolating the Planck 353 GHz spectra computed on the 400 deg² circular patches at high Galactic latitude to 150 GHz shows that we expect significant contamination by dust over most of the high Galactic latitude sky in the ℓ range of interest for detecting a primordial spectrum.

• Even for the cleanest of these regions, the Planck statistical error on the estimate of amplitude at ℓ = 80 for such small regions is at best 0.17 (3σ) in units of r_d.

• Our results show that subtraction of polarized dust emission will be essential for detecting primordial B-modes at a level of r = 0.1 or below.

• There is a significant dispersion of the polarization amplitude for a given dust total intensity. Choices of the cleanest areas of the polarized sky cannot be made accurately using the Planck total intensity maps alone.

• At present, component separation, or template cleaning, can best be done with the Planck HFI 353 GHz data, but the accuracy of such cleaning is limited by Planck noise in small fields. Ground-based or balloon-borne experiments should include dust channels at high frequency. Alternatively, if they intend to rely on the Planck data to remove the dust emission, they should optimize the integration time and area so as to have a similar S/N level for the CMB and dust power spectra.

Turning specifically to the part of the sky mapped by the BICEP2 experiment, our analysis of the M_B2 region indicates the following results.

• Over the multipole range 40 <ℓ< 120, the Planck 353 GHz power spectrum extrapolated to 150 GHz yields a value 1.32 × 10^-2 μK, with statistical error ± 0.29 × 10^-2 μK and a further uncertainty (+0.28, −0.24) × 10^-2 μK from the extrapolation. This value is comparable in magnitude to the BICEP2 measurements at these multipoles that correspond to the recombination bump.

• The frequency dependence of across the Planck bands is consistent with the typical SED of dust polarization (Planck Collaboration Int. XXII 2015).

• Assessing the dust contribution to the B-mode power measured by the BICEP2 experiment requires a dedicated joint analysis with Planck, incorporating all pertinent observational details of the two data sets, such as masking, filtering, and colour corrections.

• We have identified regions in which the dust polarization amplitude may be significantly lower, by about a factor of 2, than in the BICEP2 observing region.

---

Acknowledgments

The development of Planck has been supported by: ESA; CNES and CNRS/INSU-IN2P3-INP (France); ASI, CNR, and INAF (Italy); NASA and DoE (USA); STFC and UKSA (UK); CSIC, MICINN, J.A., and RES (Spain); Tekes, AoF, and CSC (Finland); DLR and MPG (Germany); CSA (Canada); DTU Space (Denmark); SER/SSO (Switzerland); RCN (Norway); SFI (Ireland); FCT/MCTES (Portugal); and PRACE (EU). A description of the Planck Collaboration and a list of its members, including the technical or scientific activities in which they have been involved, can be found at http://www.rssd.esa.int/index.php?project=PLANCK&page=Planck_Collaboration. Some of the results in this paper have been derived using the HEALPix package. The research leading to these results has received funding from the European Research Council under the European Union’s Seventh Framework Programme (FP7/2007-2013)/ERC grant agreement No. 267934.

References

Appendix A: Power spectrum estimator performance

Using Xpol (see Sect. 3.1.1) we have estimated the dust polarization and angular power spectrum, not on the full sky, but on particular cuts of the sky, such as the LR24 region defined in Sect. 3.3.1 and the much smaller M_B2 patch defined in Sect. 6.1. To validate the performance of Xpol on such cuts we compare the results of this algorithm to those from Xpure (see Sect. 3.1.2), using simulated data in both cases.

Simulating the dust polarization for assessment of the performance of polarization angular power spectra estimators is an important issue. One either relies on statistically isotropic, Gaussian simulations with polarization power spectra similar to those of the dust, or alternatively one uses non-Gaussian, anisotropic simulations. The problem in the first case is that the hypotheses of Gaussianity and statistical isotropy might apply only partially to a dust polarization map. The problem in the second case is that one has to produce a non-Gaussian and anisotropic simulation, which has a defined input angular power spectrum, in order to be able to characterize what is retrieved in the output. In this section, we will use statistically isotropic, Gaussian simulations to explore the performance of our algorithms for recovering the dust polarization angular power spectra on a given masked sky. This choice is motivated by the fact that cosmologists have tended to make this assumption when analysing the CMB anisotropies from a sky contaminated by non-Gaussian and anisotropic processes.

We generated 1000 full sky Q and U Gaussian map simulations from and power spectra having , the typical dust power spectrum shape that we measured in Sect. 4.2. Additionally, we gave to these spectra the B/E amplitude hierarchy that we have measured in Sect. 4.4, i.e., . For each of these simulated dust polarized signal Q and U maps, we computed the and angular power spectra on the cut sky LR24 region and the much smaller M_B2 patch. We note that the structure in the simulations is driven by the assumed power spectra and the random phases drawn for each mode and is not correlated with the geometry of the mask.

The results are shown in Fig. A.1. For each region, we give the angular power spectra recovered by Xpol and Xpure from the simulations, together with the input power spectrum and the sample variance associated with each region.

We can see that for both methods there is no overall bias in the recovered angular power spectra, demonstrating that there is no substantial leakage from E to B. At multipoles above 40, the recovered values agree with the input values and the maximum excursions (5%) seen in some bins are consistent with the expected accuracy for this number of simulations. Moreover, Xpol and Xpure give very similar results in terms of the mean value and uncertainty.In the multipole bin ℓ ≃ 40 on M_B2 the small offsets for both Xpol and Xpure (−4% and 8%, respectively) might be hinting at a limitation of the estimation procedure in this multipole bin because of the size of the patch.

A related comparative assessment of the two methods on actual Planck data is presented in Appendix D.1.

Appendix B: Complementary dust spectra at 353 GHz: , , and

Using the same procedure as described in Sects. 3 and 4 for the dust polarized angular power spectra for the LR regions at intermediate Galactic latitude, we have computed the 353 GHz dust spectra involving temperature and polarization, and , and the cross-spectrum for polarization, . This completes the entire set of polarization-related power spectra. Spectra such as and necessarily vanish for cosmic fields that satisfy statistical isotropy and parity invariance. However, the emission from our Galaxy satisfies neither of these and so the and spectra are generally expected to be non-zero.

In Fig. B.1 we present the results for f_sky = { 0.3,0.5,0.7 }; we note that the Planck best-fit CMB power spectrum was removed, as in Sect. 3.2. Again, the amplitudes of the spectra increase with increasing f_sky because the polarized emission is brighter on average when more sky is retained. The amplitudes of the spectra are about a factor of 3 higher than those reported for in Fig. 2. For the largest region there is evidence of a power-law dependence on multipole ℓ over the range ℓ = 40 to 600, with a slope compatible with that found for and in Sect. 4.2 (α_TE = −2.37, −2.43, and −2.47 for f_sky = 0.3, 0.5, and 0.7, respectively). The fact that the slopes of the spectra are similar to those found in Sect. 4.2 for the dust polarization and in Planck Collaboration Int. XXII (2015) for dust intensity, indicates that the TE cross-correlation is dominated by the dust correlations in this multipole range, and not by the other components in the intensity maps (e.g., CIB or point sources). We also show the power spectrum computed from the Planck 2013 best-fit ΛCDM model of the CMB temperature data (Planck Collaboration XVI 2014). As in Fig. 2 for , the CMB model is below the 353 GHz angular power spectra at low ℓ and about the same order of magnitude at ℓ> 300.

Similarly, the results for are shown in Fig. B.2, but now on a linear scale. For the two largest regions (f_sky = 0.5 and 0.7), we detect a significant signal in the range ℓ = 60−130 that increases with the sky area. In the same multipole range, the dust is only marginally detected for the f_sky = 0.3 region. At lower and higher ℓ, the individual binned measurements of on all the regions are compatible with zero.

Finally, we show the results for in Fig. B.3. There is a marginal detection of the spectrum in the range ℓ = 60−130 for the f_sky = 0.7 mask. For the other masks and at other multipoles, the spectra are compatible with zero.

Appendix C: Notes on the analysis of the high Galactic latitude patches

Appendix C.1: Assessment of Planck polarization systematic uncertainties

To assess the level of systematic effects in the results presented in Sect. 5, we use a null-test analysis (see also Appendix D.3). We repeated the statistical analysis for the 353 GHz and cross-spectra on the high-Galactic latitude patches. The steps of the analysis are exactly the same as in Sect. 5.1, except that instead of using the 353 GHz 353_DS1 × 353_DS2 cross-power spectra, we use the (353_DS1,HR1−353_DS1,HR2)/2 × (353_DS2,HR1−353_DS2,HR2)/2 half-difference cross-spectra. These cross-spectra are expected to be consistent with zero for uncorrelated noise and subdominant systematics.

The computed dust and amplitudes A^EE,BB (at ℓ = 80 computed on each of these patches) for the null-test spectra, in μK at 353 GHz, are shown in Fig. C.1 as a function of ⟨I₃₅₃⟩. We also report in this figure the expected noise level computed from Monte Carlo simulations of Planck inhomogeneous white noise.

We can see that the fitted amplitudes for the null-test spectra on the 352 high Galactic latitude patches are scattered around zero, without any obvious bias, for the entire range of ⟨I₃₅₃⟩; this is the expectation for half-difference cross-spectra. Taking into account the individual error bars, there is no patch that has a null-test spectrum amplitude inconsistent with the ± 3σ range of the instrumental noise dispersion.

We thus conclude from this analysis that the results presented in Sect. 5 are not significantly affected by instrumental systematics.

Appendix C.2: Dust B-mode amplitudes compared to E-mode

In Sect. 5, we computed the 353 GHz dust and amplitudes on the 352 high Galactic latitude patches defined in Sect. 3.3.2. We present in Fig. C.2 the ratio of the amplitudes of these spectra, A^BB/A^EE as a function of the mean dust intensity ⟨I₃₅₃⟩ of the patch on which it was computed.

These ratios are scattered around the inverse noise variance weighted ⟨A^BB/A^EE⟩ = 0.52 ± 0.19 value. This value is consistent with the one computed in Sect. 4.4 on the LR regions. Nevertheless, it can be seen from Fig. C.2 that there are patches with a significantly different A^BB/A^EE ratio, e.g., some smaller values for ⟨I₃₅₃⟩ in the range 0.08−0.1 MJy sr^-1.

We additionally report in Fig. C.2 the A^BB/A^EE value determined for the M_B2 spectra computed in Sect. 6.3. The A^EE and A^BB amplitudes are computed as they are in Sect. 5, fitting a power law in ℓ to the and spectra with a fixed exponent α_EE,BB = −2.42. Thus, for M_B2, we find that under these assumptions A^BB/A^EE = 0.70 ± 0.19. This value is higher than the average for such fields, but consistent with the mean expected value.

Appendix D: Systematic effects relating to the power spectrum estimate in the BICEP2 field

Appendix D.1: Dust spectrum measurements in the BICEP2 field using Xpol and Xpure

We present in Fig. D.1 the difference of the power spectra computed by Xpol and Xpure in M_B2 (). The processing is the same as described in Sect. 6.3.

The Xpol and Xpure angular power spectra are consistent with each other in all three ℓ bins. The maximum difference observed is in the ℓ = 40–120 bin, where the Xpure spectrum is 10% higher than the Xpol spectrum presented in Sect. 6.3. Moreover, the Xpol recovered error bars and the Xpure error bars from Monte Carlo simulations of Planck inhomogeneous noise, also presented in Fig. D.1, are consistent.

We conclude from this comparison that the results presented in Sect. 6.3 do not depend significantly on the method used for estimating the dust spectrum.

Appendix D.2: Impact of uncertainties in the specification of the BICEP2 field

Our estimation of the dust angular power spectrum in M_B2 can depend on the detailed definition of the field. To explore the dependence of our results on the extent and specification of our approximation of the BICEP2 field, we repeat our computation of the dust power spectrum on the 353 GHz data extrapolated to 150 GHz, as in Sect. 6.3, for slight variations in the definition of this field.

For this purpose, we start from the mask M defined in Sect. 6.1, which consists of zeros and ones. We increase its size (with a smoothing and a cut) in order to have an outline that extends exactly 1° further in every direction. This defines (“b” for binary, since it contains zeros and ones like M). We then decrease the size of (with smoothing and cut) by exactly 1° in every direction to obtain . After that is shrunk in the same way to obtain and we proceed iteratively down to . Then to each mask we apply exactly the same procedure for apodizing that we applied to M to produce M_B2 in Sect. 6.1 to yield the final M_i field. The effective area () is 805 deg² for M₅, 677 deg² for M₄, 564 deg² for M₃, 460 deg² for M₂, and 365 deg² for M₁. We note that because of the forward and backward smoothing, M_B2 defined in Sect. 6.1 is similar, but not exactly equal to M₄.

We compute the 353 GHz power spectrum on all five fields using Xpol and extrapolate the results to 150 GHz with the frequency dependence from Planck Collaboration Int. XXII (2015), as in Sect. 6.3. We show in Fig. D.2 the difference between these spectra and the spectrum computed on M_B2 presented in Sect. 6.3. We also show the error budget for M_B2, centred on zero because we are looking at differences between power spectra.

The recovered scatter is very low in the bin with ℓ = 40–120 (the standard deviation of the mean value among the fields is 0.04 × 10^-2 μK), much lower than the statistical uncertainty in this particular bin. In the two other bins, the scatter is more important, but still consistent with the statistical uncertainty.

We conclude from this analysis that the main result presented in Sect. 6.3 does not depend on our precise definition of M_B2.

Appendix D.3: Assessment of the Planck HFI systematics in the BICEP2 field

To assess the potential contribution of systematic effects to the dust 353 GHz spectrum computed on M_B2, we repeat the analysis presented in Sect. 6.3 on null tests coming from subtracting different Planck 353 GHz data subsets.

The cross-DetSets, cross-HalfRings, and cross-Years spectra are computed for the 353 GHz data extrapolated to 150 GHz. The differences of these cross-spectra (DetSets minus Years, Years minus HalfRings, and HalfRings minus DetSets) are presented in Fig. D.3. In the first ℓ bin (ℓ = 40−120), the three differences of the dust power spectra are consistent with zero, given the statistical uncertainty. They exhibit a dispersion consistent with what is expected from statistical noise. In the other bins the dispersion is larger, but still consistent with the noise. We conclude that the systematic errors are subdominant.

In addition to these data cuts, we create null-test spectra on M_B2 from the four data subsets at 353 GHz involving detector sets and ring halves. We compute the cross-DetSet spectra between their two HalfRing half-differences (353_{DS1, HR1}−353_{DS1, HR2})/2 × (353_{DS2, HR1}−353_{DS2, HR2})/2 and the cross-HalfRings between the two DetSet half-differences (353_{DS1, HR1}−353_{DS2, HR1})/2 × (353_{DS1, HR2}−353_{DS2, HR2})/2. These null-test spectra are also shown in Fig. D.3. We can see that these spectra are consistent with zero in all the ℓ bins, given the statistical and systematic uncertainties.

We list in Table D.1 the χ² values for the null hypothesis, as well as the PTE for these null tests (plus additional values not presented in Fig. D.3). All these tests are consistent with the null hypothesis.

Appendix D.4: Synchrotron contribution

An estimate can be made of the level of a synchrotron contribution in M_B2 by a simple extrapolation from the lowest frequency Planck LFI channels using a reasonable range of spectral indices. Starting with the Planck LFI 28.4 GHz Q and U maps masked to retain M_B2, a pure pseudo-C_ℓ power spectrum estimator was used to calculate the B-mode band-powers . The noise bias and significance was estimated from 100 FFP8 noise simulations¹⁷, which show the 28.4 GHz power spectrum to be completely within the expected ± 2σ noise envelope. We therefore use the 3σ upper limit on the scatter of the simulations as the maximum allowable synchrotron signal at 28.4 GHz.

A spectral index of −3.12 ± 0.04 for synchrotron emission at high Galactic latitudes was found by Fuskeland et al. (2014) using TT plots of WMAP 9 yr data. Adopting a spectral index of −3.1 for a direct extrapolation to 150 GHz, we find 3σ upper limits on the possible contribution of synchrotron B-mode power at ℓ = 45 and ℓ = 74 to be 2.1%, expressed relative to the primordial spectrum with r = 0.2. This estimate is consistent with the values from Fuskeland et al. (2014). A more detailed and complete analysis of diffuse polarized synchrotron emission will appear in a forthcoming Planck paper.

Appendix E: Dust polarization decorrelation

Extrapolation of the polarized dust signal observed at 353 GHz to lower frequencies (relevant to the CMB) hinges on the assumption that the pattern of dust polarization on the sky is independent of frequency, with only the amplitude changing, scaling as the SED. However, because the dust emission is the sum of many contributions from grains of different composition, size, shape, and alignment, which will change along the line of sight, there can in principle be some decorrelation of the signal from frequency to frequency (e.g. Tassis & Pavlidou 2014). We have seen in Sect. 4.5 that the dust polarization power spectra amplitudes are compatible with the typical dust SED, indicating that if there is such a decorrelation it has to be relatively weak. Nevertheless, in this Appendix, we will assess this possible effect quantitatively.

If there is a decorrelation of the dust polarization with frequency, one expects the cross-spectrum between two frequencies to be smaller than the geometric mean of the two auto-spectra. Hence, we compute the ratio between the cross- and auto-spectra for 217 and 353 GHz, where we have the sensitivity to quantify the effect. We compute at these two frequencies and for the cross-spectrum and fit a power-law amplitude A^BB in each case, with a slope fixed to α_BB = −2.42, as in Sect. 4.2.

The ratios are shown in Fig. E.1 as a function of ⟨I₃₅₃⟩  for the 352 patches used in Sect. 5. Although there is a large dispersion at low ⟨I₃₅₃⟩  one can see that the ratios are centred on unity, indicating that on average there is no decorrelation; among the 353 patches, we find a mean decorrelation ratio of .

In Fig. E.1 we also show this ratio for the six LR regions of Sect. 4 and the M_B2 region of Sect. 6, finding no evidence for a significant decorrelation; for the LR regions, we specifically find .

All Tables

All Figures