© ESO, 2014

1. Introduction

This paper, one of a set associated with the 2013 release of data from the Planck¹ mission (Planck Collaboration I 2014), describes the component separation techniques applied to the Planck data to produce maps of the cosmic microwave background (CMB) temperature anisotropies (see Fig. 1) and of diffuse foregrounds.

The sky at millimetre and sub-millimetre wavelengths contains a wealth of cosmological and astrophysical information. Accessing it is an inversion process, known as component separation, to extract the sources of emission contributing to a set of maps observed at different frequencies. Planck gives us a powerful data set to unlock new information in this manner by observing the entire sky from 30 to 857 GHz in nine frequency bands at higher angular resolution and sensitivity than its predecessors. Accurate and detailed component separation is a central objective of the mission.

We divide the foregrounds into two distinct categories: diffuse emission from the Galaxy and compact sources. The Galactic foregrounds are the principal source of contamination of the CMB on large angular scales, with fluctuation power decreasing roughly as a power law towards higher multipoles (Bennett et al. 2003). They are dominated by synchrotron, free-free and anomalous microwave emission (AME, ascribed to spinning dust grains) at frequencies below 70 GHz, and by rotational line emission from carbon monoxide (CO) molecules and thermal dust emission at frequencies above 100 GHz. Extragalactic foregrounds, on the other hand, dominate the small-scale contamination of the CMB. They arise from discrete, individually detectable compact sources and the collective emission from unresolved radio and infrared (IR) sources, and also from the Sunyaev-Zeldovich (SZ) effect in galaxy clusters (Planck Collaboration XVIII 2011; Planck Collaboration XXVIII 2014; Planck Collaboration XXIX 2014).

In the Planck analyses, these foregrounds are dealt with in a variety of ways. At the power spectrum and likelihood level, the extragalactic foregrounds are modelled with parameterized power spectra, appropriate to their statistical isotropy, over regions restricted to low Galactic emission (Planck Collaboration XV 2014). Component separation as described in the present paper aims at removing Galactic emission to produce CMB maps covering the largest possible sky area for studies of the large-scale properties and higher-order statistics of the CMB. In addition, this component separation provides a reconstruction of the diffuse emission from our Galaxy. Detailed studies of specific extragalactic foregrounds, such as the cosmic infrared background (CIB; Planck Collaboration XVIII 2014) and the diffuse SZ signal (Planck Collaboration XXI 2014), employ methods tailored to their particular needs.

Building on previous work (Leach et al. 2008), we approach CMB extraction with a philosophy designed to ensure robustness by applying four distinct algorithms based on two different methodologies. The first avoids any assumptions concerning the foregrounds and relies solely on a minimum variance criterion for the data component possessing a blackbody spectrum (i.e., the CMB), while the second methodology relies on parametric modelling of the foregrounds in either real or harmonic space. We evaluate the performance of these component separation algorithms through detailed simulations, and we examine the robustness of the recovered CMB maps by comparing them, their power spectra, and their resulting cosmological constraints. As a diagnostic, we also briefly examine their higher-order statistics.

The CMB results presented in this work serve a number of applications. We use the real-space modelling to produce a clean CMB map and power spectra on large angular scales, where diffuse Galactic emission is the main contaminant, to construct the likelihood function at low multipoles; this is then combined with the high multipole likelihood function that models extragalactic foregrounds with power spectra (Planck Collaboration XV 2014). The high resolution CMB maps are used as a check on primary cosmological constraints (see below), for lensing studies (Planck Collaboration XVII 2014), studies of the integrated Sachs-Wolfe ISW effect (Planck Collaboration XIX 2014), of the isotropy of the CMB (Planck Collaboration XXIII 2014), of non-Gaussian statistics (Planck Collaboration XXIV 2014), in searches for topological defects (Planck Collaboration XXV 2014), and for examination of the geometry and topology of the Universe (Planck Collaboration XXVI 2014).

In addition, we present maps of diffuse Galactic emission divided into low- and high-frequency components, as well as a molecular CO component. We judge the adequacy of this reconstruction through simulations and by comparison with known properties of the diffuse Galactic foregrounds.

The paper is organized as follows. In Sect. 2 we discuss the expected sources of sky emission over the Planck frequency interval and how they are modelled. Then in Sect. 3 we detail the overall approach and introduce the four component separation methods. In Sect. 4 we present the Planck data set and pre-processing procedure, and we describe our simulations. This is followed by a presentation of the derived CMB maps and their characterization in Sect. 5. Section 6 is dedicated to power spectra and cosmological parameter constraints obtained from these maps, and Sect. 7 to studies of higher-order statistics. Section 8 presents a reconstruction of the diffuse Galactic foregrounds, and Sect. 9 concludes. We relegate details of the algorithms to appendices.

2. The sky at Planck frequencies

The properties of Galactic emission vary significantly across the Planck frequency range from 30 to 857 GHz. At frequencies below 70 GHz, the dominant radiation processes are: synchrotron emission from cosmic ray electrons interacting with the Galactic magnetic field (e.g., Haslam et al. 1982; Reich & Reich 1988; Broadbent et al. 1989; Davies et al. 1996; Platania et al. 2003; Bennett et al. 2003; Gold et al. 2011); thermal Bremsstrahlung (or free-free emission) from electron-electron and electron-ion scattering (e.g., Banday et al. 2003; Dickinson et al. 2003; Davies et al. 2006; Ghosh et al. 2012; Planck Collaboration Int. XII 2013; Planck Collaboration XX 2011); and AME from dust grains (Kogut 1996; Leitch et al. 1997; Banday et al. 2003; Lagache 2003; de Oliveira-Costa et al. 2004; Finkbeiner et al. 2004; Davies et al. 2006; Bonaldi et al. 2007; Dobler & Finkbeiner 2008; Miville-Deschênes et al. 2008; Ysard et al. 2010; Gold et al. 2011; Planck Collaboration XX 2011), possibly due to their rotational emission (Draine & Lazarian 1998; Ali-Haïmoud et al. 2009; Ysard & Verstraete 2010; Hoang & Lazarian 2012). Over the frequency range covered by Planck, both synchrotron and free-free spectra are well approximated by power laws in antenna temperature, T_B ∝ ν^β, with the synchrotron index, β_synch, ranging from − 3.2 to − 2.8 (Davies et al. 1996) and the free-free index, β_ff, lying between − 2.2 and − 2.1. Less is known about the AME spectrum, but spinning dust models with a spectrum peaking at frequencies below 20 GHz (in antenna temperature units) adequately describe current observations². Above the peak, the spectrum appears consistent with a power-law (e.g., Banday et al. 2003; Davies et al. 2006; Dobler & Finkbeiner 2008; Ghosh et al. 2012). In addition to these three, the existence of a fourth low-frequency foreground component, known as the “Galactic haze”, has been claimed, possibly due to a hard-spectrum synchrotron population near the Galactic centre (e.g., Finkbeiner 2004; Dobler & Finkbeiner 2008; Pietrobon et al. 2012; Planck Collaboration Int. IX 2013).

At frequencies higher than 100 GHz, thermal dust emission dominates over most of the sky and is commonly described by a modified blackbody spectrum with power-law emissivity, ϵ_ν ∝ ν^β_d, and temperature, T_d. Both the temperature and spectral index, β_d, vary spatially. Prior to Planck, the best-fitting single component dust model had a temperature T_d ≈ 18 K and spectral index β_d ≈ 1.7 (Finkbeiner et al. 1999; Bennett et al. 2003; Gold et al. 2011), although there is evidence of flattening of the spectral index from around 1.8 in the far-infrared to 1.55 in the microwave region (Finkbeiner et al. 1999; Planck Collaboration Int. XIV 2014), the interpretation of which is still under study.

In addition to these diffuse Galactic components, extragalactic emission contributes at Planck frequencies. In particular, a large number of radio and far-infrared (FIR; Planck Collaboration XIII 2011) galaxies, clusters of galaxies and the CIB (Planck Collaboration XVIII 2011) produce a statistically isotropic foreground, with frequency spectra well approximated by models similar to those applicable to the Galactic foregrounds (modified blackbody spectra, power laws, etc.). Except for a frequency-dependent absolute offset, which may be removed as part of the overall offset removal procedure, these extragalactic components are therefore typically absorbed by either the low-frequency or thermal dust components during component separation. No special treatment is given here to extragalactic foregrounds, beyond the masking of bright objects. Dedicated scientific analyses of these sources are described in detail in Planck Collaboration XVIII (2011), Planck Collaboration XXVIII (2014), and Planck Collaboration XXIX (2014). In the Planck likelihood, extragalactic sources are modelled in terms of power spectrum templates at high ℓ (Planck Collaboration XV 2014).

Other relevant sources include emission from molecular clouds, supernova remnants, and compact H ii regions inside our own Galaxy, as well as the thermal and kinetic SZ effects, due to inverse Compton scattering of CMB photons off free electrons in ionized media. Planck provides new and important information on all these processes, as described both in the following and in the companion papers Planck Collaboration XIII (2014), Planck Collaboration XXI (2014), and Planck Collaboration XI (2014). In particular, Planck’s frequency range, angular resolution and sensitivity make it a powerful probe of thermal dust, resulting in new and tight constraints on dust temperature and emissivity. The same frequencies also allow extraction of the first ever full-sky maps of the emission resulting from the CO J = 1 → 0, J = 2 → 1 and J = 3 → 2 rotational transitions at 115, 230 and 345 GHz, respectively (Planck Collaboration XIII 2014).

The focus of this paper is to reconstruct the CMB anisotropies over a large sky fraction, exploiting only the Planck frequency bands. We also present a detailed reconstruction of the thermal dust emission at high frequencies, as well as CO emission lines. At low frequencies and over the region used for CMB analysis, the total foreground contribution is well approximated by a single power law (see Sect. 8). We therefore model the sum of all low-frequency foregrounds by a power law with spatially varying spectral index whose numerical value in any pixel results from the influence of the dominant foreground component at that location. The full analysis of diffuse foregrounds, using ancillary data to resolve the individual components at low frequencies, will be presented in a forthcoming publication.

3. Approach to component separation

The rich content of the Planck data encourages application of several component separation techniques. We consider four, as summarized in Table 1, which we classify according to one of two different general methodologies. The first makes minimal assumptions concerning the foregrounds and seeks only to minimize the variance of the CMB, i.e., the sky component possessing a blackbody spectrum. We implement this approach with a needlet (wavelet on the sphere) version of the internal linear combination (ILC) algorithm (NILC; Delabrouille et al. 2009), and also with a template-based method to remove foreground contamination from the CMB-dominant bands. These foreground templates are constructed from the lowest and highest frequency channels (Fernández-Cobos et al. 2012, Spectral Estimation Via Expectation Maximization, SEVEM).

The second methodology uses parametric modelling of the foregrounds. In our real space implementation, we explore model parameters through Bayesian parameter estimation techniques, fitting a parametric signal model per pixel (Commander; Eriksen et al. 2006, 2008); a similar implementation is presented by Stompor et al. (2009). To estimate spectral indices robustly in pixel space, this procedure requires identical angular resolution across all frequencies included in the analysis, and is therefore limited in resolution by the 30 GHz LFI channel. However, this is sufficient to generate the low-resolution CMB map and power spectrum samples required for the low multipole part of the Planck likelihood function for cosmological parameters (Planck Collaboration XV 2014). To produce full resolution maps, we use the resulting low-resolution spectral parameter samples to solve for the component amplitudes, in an extension to the method known as Ruler (we refer to the combined method as Commander-Ruler, or C-R). In our fourth technique, we implement a CMB-oriented parametric approach that fits the amplitude and spectral parameters of CMB and foregrounds in the harmonic domain (Spectral Matching Independent Component Analysis, SMICA; Cardoso et al. 2008).

Details of each algorithm are given in the appendices. We now turn to their application to the data and evaluate their performance using simulations.

4. Data, simulations and masks

We use the data set from the first 15.5 months of Planck observations, corresponding to 2.6 sky surveys, from both the Low Frequency Instrument (LFI) and High Frequency Instrument (HFI). The primary inputs for component separation are the frequency channel maps, including half-ring maps, bandpasses, and beam characteristics; a full description of these products is given in Planck Collaboration II (2014) and Planck Collaboration VI (2014). No special corrections are made for zodiacal light emission (ZLE; Planck Collaboration VI 2014) in the analyses presented here. The ZLE is not stationary on the sky, since it depends on Planck’s position and scanning strategy. Therefore the frequency maps contain a projected version of the emission averaged over the nominal mission. Despite this, a series of exploratory analyses showed that our algorithms naturally correct for this component within their existing model space. It was also found that larger CMB residuals were induced when applying a correction based on a ZLE model than when applying no correction, most likely due to uncertainties in the model itself.

To evaluate and validate our algorithms, we analyse a large suite of realistic simulations, the so-called full focal plane (FFP) simulations, based on detailed models of the instrument and sky. The version used for this data release is denoted FFP6, and is described in Planck Collaboration (2013). The simulation procedure generates time streams for each detector, incorporating the satellite pointing, the individual detector beams, bandpasses, noise properties, and data flags, and then produces simulated frequency channel maps through the mapmaking process. For the input sky, we use the Planck Sky Model (PSM), which includes the CMB, diffuse Galactic emission (synchrotron, free-free, thermal dust, AME, and molecular CO lines), and compact sources (thermal and kinetic SZ effects, radio sources, infrared sources, the CIB, and ultra-compact H ii regions). The pre-launch version of the PSM is described by Delabrouille et al. (2013), and has been modified for the present work as described in Planck Collaboration (2013). Each FFP data set consists of three parts: the simulated observations, Monte Carlo realizations of the CMB, and Monte Carlo realizations of the instrumental noise.

For both the data and the simulations, we reconstruct the CMB and foregrounds from the full frequency channel maps and the corresponding half-ring maps, which are made from the data in the first half or second half of each stable pointing period. The half-ring maps can be used to obtain an estimate of the noise in each channel by taking half of the difference between the two maps, thereby normalizing the noise level to that of the full map. This is referred to as the half-ring half-difference (HRHD) map. The signals fixed to the sky will be cancelled leaving only the noise contribution. The HRHD map can be treated as a realization of the same underlying noise processes and it can be used to estimate the power spectrum, and other properties, of the noise. If there are noise correlations between the half-ring maps, then the estimates of the noise properties thus obtained can be biased. This is the case for HFI channels; the cosmic ray glitch removal (Planck Collaboration VI 2014; Planck Collaboration X 2014) induces correlations that lead to the noise power spectrum being underestimated by a few percent at high ℓ when using the HRHD maps.

Prior to processing the data through each component separation pipeline, we define masks for the point sources and bright Galactic regions. Point source masking is based on the source catalogues obtained by filtering the input sky maps with the Mexican Hat wavelet 2 (MHW2) filter and applying a 4σ threshold for the LFI bands and a 5σ threshold for the HFI bands (Planck Collaboration XVIII 2011; Planck Collaboration XXVIII 2014). The mask radius of each source is different for the LFI and HFI. Due to the large beam size of LFI channels, we define a variable masking radius for each source according to its signal-to-noise ratio (S/N) as $\mathit{r}\mathrm{=}\mathrm{(}\mathrm{2}\log \hspace{0.17em}\mathrm{(}\mathit{A}\mathit{/}\mathit{m}\mathrm{)}{\mathrm{)}}^{\mathrm{1}\mathit{/}\mathrm{2}}\mathit{/}\mathrm{(}\mathrm{2}\sqrt{\mathrm{2}\log \mathrm{2}}\mathrm{)}\hspace{0.17em}\mathrm{\times }\hspace{0.17em}\mathit{FWHM}$, where r is the radius, A is the S/N, and m is the maximum amplitude (given in units of the background noise level) allowed for the tail of unmasked point sources; we set m = 0.1, which is a compromise between masking the source tails and minimizing the number of masked pixels. For HFI, the mask radius around each source is 1.27 × FWHM, using the average FWHM obtained from the effective beams.

A basic set of Galactic masks is defined as follows. We subtract a CMB estimate from the 30 and 353 GHz maps, mask point sources, and smooth the resulting maps by a Gaussian with FWHM of 5°. We then threshold and combine them, generating a series of masks with different amounts of available sky. The resulting combined Galactic (CG) masks, shown in Fig. 2, correspond to sky fractions of 20, 40, 60, 70, 75, 80, 90, 97, and 99%, and are named CG20, etc.

5. CMB Maps

We begin the discussion of our results by presenting the foreground-cleaned CMB maps. These maps are shown in Fig. 1 for each of the four component separation algorithms. Already from this figure it is clear that the wide frequency coverage and high angular resolution of Planck allow a faithful reconstruction of the CMB field over most of the sky. The fluctuations appear visually consistent with the theoretical expectation of a Gaussian and isotropic signal everywhere except inside a small band very close to the Galactic plane³

Each CMB map is accompanied by its own confidence mask outside which the corresponding solution is considered statistically robust, shown in Fig. 3; for a definition of each mask, see Appendices A–D. Accepted sky fractions are 75, 93, 76, and 89%, respectively, for Commander-Ruler, NILC, SEVEM, and SMICA. These masks are denoted CS-CR75, CS-NILC93, CS-SEVEM76, CS-SMICA89, respectively. The union of the confidence masks accepts 73% of the sky and is denoted U73. It is adopted as the default mask for evaluation purposes in this paper.

In addition to the CMB maps from the full data set, the half-ring frequency maps have been processed by each algorithm to provide half-ring CMB maps. They are used to provide estimates of the instrumental noise contribution to the foreground-cleaned maps in the power spectrum analysis (see Sect. 6). The algorithms were also used to process Monte Carlo simulations: 1000 realizations of the CMB and 1000 realizations of noise. They are not used in the analyses presented in this paper, but are used by Planck Collaboration XXIII (2014) and Planck Collaboration XXIV (2014).

The beam transfer functions of the foreground-cleaned CMB maps have been estimated for each algorithm, as shown in Fig. 4. The angular resolution of the NILC, SEVEM, and SMICA maps corresponds to a Gaussian beam with FWHM of 5′. The difference between SEVEM and NILC/SMICA is due to their different treatment of the HEALPix⁴ pixel window function (Górski et al. 2005). The deviation of NILC beam from a Gaussian shape at ℓ > 2800 is caused by the last needlet window (see Appendix B). Commander-Ruler has a larger beam, because it is defined explicitly as a weighted average of frequency maps in pixel space. Its resolution is equivalent to a Gaussian beam with FWHM of approximately . The beam transfer functions have been computed assuming the best-fit beam transfer function for each frequency channel, and the uncertainties in the latter have not been propagated to these estimates.

In Fig. 5 we show the standard deviation per pixel among the four foreground-cleaned CMB maps downgraded to N_side = 128, and in Fig. 6 we show all pairwise difference maps. Typical differences at high Galactic latitudes are smaller than 5 μK. Considering the difference maps in more detail, it is clear that the Commander-Ruler map is the most different from the other three, whereas NILC and SMICA are the most similar. This is not completely unexpected, because while Commander-Ruler uses only frequencies between 30 and 353 GHz in its solution, the other three codes additionally include the dust-dominated 545 and 857 GHz maps.

This difference in data selection may explain some of the coherent structures seen in Fig. 6. In particular, the most striking large-scale feature in the difference maps involving Commander-Ruler is a large negative band roughly following the ecliptic plane. This is where the ZLE (Planck Collaboration VI 2014) is brightest. Since the ZLE is also stronger at high frequencies, having a spectrum close to that of thermal dust, it is possible that this pattern may be an imprint of residual ZLE either in the Commander-Ruler map, or in all of the other three maps. Both cases are plausible. The Commander-Ruler solution may not have enough high-frequency information to distinguish between ZLE and normal thermal dust emission, and, by assuming a thermal dust spectrum for the entire high-frequency signal at 353 GHz, over-subtracts the ZLE at lower frequencies. It is also possible that the other three CMB solutions have positive ZLE residuals from extrapolating the high-frequency signal model from 857 GHz to the CMB frequencies. Without an accurate and detailed ZLE model, it is difficult to distinguish between these two possibilities. It is of course also possible that the true explanation is in fact unrelated to ZLE, and the correlation with the ecliptic plane is accidental. In either case, it is clear that the residuals are small in amplitude, with peak-to-peak values typically smaller than 10 μK, of which by far the most is contained in a quadrupole aligned with the ecliptic. This provides additional evidence that residual ZLE is not important for the CMB power spectrum and cosmological parameter estimation, although some care is warranted when using these maps to study the statistics of the very largest angular scales (e.g., Planck Collaboration XXIII 2014); checking consistency among all four maps for a given application alleviates much of this concern.

We end this section by showing in Fig. 7 a set of residual maps derived by analysing the FFP6 simulation with exactly the same analysis approaches as applied to the data. It is evident that SMICA produces the map with lowest level of residuals. Considering the morphology in each case, we see that the main contaminant for Commander-Ruler is under-subtracted free-free emission, while for both NILC and SEVEM it is over-subtracted thermal dust emission, and for SMICA it is under-subtracted thermal dust emission. However, at high latitudes and outside the confidence masks, the residuals are generally below a few μK in amplitude. It is also worth noting that each algorithm has been optimized (in terms of model definition, localization parameters, etc.) for the data, and the same configuration was subsequently used for the FFP6 simulations without further tuning. The simulations presented here therefore provide a conservative estimate of the residuals in the data. This is also reflected in the fact that the differences between CMB reconstructions for the FFP6 simulations are larger than those found in the data. See Appendix E for further details.

6. Power spectrum and cosmological parameters

In this section we evaluate the foreground-cleaned maps in terms of CMB power spectra and cosmological parameters. Our purpose in doing this is to show that the maps are consistent with the high-ℓ likelihood obtained from the cross-spectrum analysis of detector set and frequency maps in Planck Collaboration XV (2014), and with the cosmological parameters derived from them in Planck Collaboration XVI (2014). This also establishes the consistency between Planck’s cosmological constraints and studies of the large-scale structure and higher order statistics of the CMB.

6.1. Power spectra

Figure 8 shows the power spectra of the foreground-cleaned CMB maps and the corresponding HRHD maps, evaluated using the U73 mask with a 30′ cosine apodization. The spectra have been corrected for the effect of the mask and the beam transfer function of each algorithm has been deconvolved. The spectra of the HRHD maps give an estimate of the instrumental noise contribution to the power spectrum of the cleaned map. The correlations between the HFI half-ring frequency maps are inherited by the half-ring CMB maps that use them as input. At small angular scales, the CMB solution comes almost entirely from data in the HFI channels, and therefore the spectrum of the CMB HRHD maps is also biased low.

At small angular scales, the effective noise levels of NILC, SEVEM, and SMICA are very similar, and lower than that of Commander-Ruler. The last has larger noise because it operates entirely in pixel space and therefore applies the same weights to all multipoles. It cannot take advantage of the changing signal-to-noise ratio of the frequency channels with angular scale.

We can estimate the contribution of residual foregrounds to the foreground-cleaned CMB maps by making use of the FFP6 simulations. In addition to processing the simulated frequency maps, the maps of the individual input sky components were processed by the algorithms after fixing their parameters or weights to the values obtained from the “observed” maps. Figure 9 shows the power spectra of the simulated FFP6 components, in this case CMB, noise and the sum of the foreground components. The top panel shows the spectra computed using the union mask derived from the simulation with a 30′ cosine apodization. The total foreground contribution becomes comparable to the CMB signal at ℓ ≈ 2000. The bottom panel shows the same computed with an apodized point source mask applied to the maps (i.e., no diffuse masking, although this mask does removes a large part of the Galactic plane). The residual foreground contribution is larger at all angular scales, but still it only becomes comparable to the CMB signal at ℓ ≈ 1800 in the worst case. For both masks, SMICA has the smallest residual foreground contamination at large angular scales, which is also demonstrated in Fig. 7. A more detailed examination of the contribution of the individual foreground components to the power spectrum is in Appendix E.

6.2. Likelihood and cosmological parameters

We estimate the binned power spectra with XFaster (Rocha et al. 2010, 2011) and determine cosmological parameter constraints using a correlated Gaussian likelihood. Parameter constraints are derived using a Metropolis-Hastings Markov Chain Monte Carlo sampler. To speed up this process, we additionally use PICO (Parameters for the Impatient COsmologist; Fendt & Wandelt 2008), a tool which interpolates the CMB power spectra and matter power spectra as a function of cosmological parameters.

6.2.1. Model and methods

We compute the power spectrum for each foreground-cleaned map over the multipole range 2 ≤ ℓ ≤ 2500, while parameter constraints are derived using only 70 ≤ ℓ ≤ 2000; as shown in Appendix E through simulations, modelling errors become non-negligible between ℓ = 2000 and 2500. For parameter estimation, we adopt a standard six-parameter ΛCDM model, and impose an informative Gaussian prior of τ = 0.0851 ± 0.014, since polarization data are not included in this analysis.

While the foreground-cleaned maps should have minimal contamination from diffuse Galactic emission, they do contain significant contamination from unresolved extragalactic sources. These contributions are most easily modelled in terms of residual power spectra, therefore we marginalize over the corresponding parameters at the power spectrum level. To the six ΛCDM parameters, describing the standard cosmology, we add two foreground parameters, A_ps, the amplitude of a Poisson component (and hence constant, C_ℓ = A_ps), and A_cl, the amplitude of a clustered component with shape D_ℓ = ℓ(ℓ + 1)C_ℓ/ 2π ∝ ℓ ^0.8. Both are expressed in terms of D_ℓ at ℓ = 3000 in units of  μK².

The power spectrum calculation is based on the half-ring half-sum (HRHS) and HRHD CMB maps (see Sect. 5); the latter is used to estimate the noise bias in the power spectra extracted from the HRHS maps. From these, we calculate the pseudo-spectra, $\begin{array}{c}{\mathrm{˜}}_{}\\ {\mathit{C}}_{\mathit{\ell }}\end{array}$ and $\begin{array}{c}{\mathrm{˜}}_{}\\ {\mathit{N}}_{\mathit{\ell }}\end{array}$ (Hivon et al. 2002), respectively, after applying the U73 mask. These are used as inputs to XFaster together with the beam transfer functions provided by each method (see Fig. 4).

To avoid aliasing of power from large to small scales, which would add an offset between the signal-plus-noise and noise pseudo-spectra at high ℓ, we use the apodized version of the U73 mask. The known mismatch in the noise level between the spectra due to the correlation between the half-ring maps is not explicitly corrected. It is left to be absorbed into the two foreground parameters.

Using the pseudo-spectra and XFaster, we then reconstruct an estimate of the power spectrum of each foreground-cleaned HRHS map, removing the noise bias as estimated from the corresponding HRHD map. To this end we apply an iterative scheme starting from a flat spectrum model. The result is a binned power spectrum and the associated Fisher matrix, which are then used to construct the likelihood, approximated here by a correlated Gaussian distribution.

To study consistency in the low-ℓ range, we fit a two-parameter q–n (amplitude-tilt) model relative to the Planck best-fit ΛCDM model on the form, ${\mathit{C}}_{\mathit{\ell }}\mathrm{=}\mathit{q}\mathrm{(}\mathit{\ell }\mathit{/}{\mathit{\ell }}_{\mathrm{pivot}}{\mathrm{)}}^{\mathit{n}}{\mathit{C}}_{\mathit{\ell }}^{\mathrm{bf}}$, using a pixel-space likelihood for maps smoothed to 6° FWHM; see Planck Collaboration XV (2014) for further algorithmic details.

6.2.2. Results

We perform the power spectrum and parameter estimation analysis for both the data and the FFP6 simulations described in Sect. 4. The results for the latter are given in Appendix E.

Figure 10 shows estimates of the angular power spectrum for each foreground-cleaned map, with the uncertainties given by the Fisher matrix. The parameter summary given in Fig. 11 shows the parameter constraints derived using multipoles between ℓ = 70 and 2000, and compares these to results obtained with the CamSpec and Plik likelihoods (Planck Collaboration XV 2014).

Differences in the power spectra at high ℓ are mostly absorbed by the two-parameter foreground model, rendering consistent cosmological parameters. For example, the high-ℓ power excess seen in the Commander-Ruler map is well-fitted in terms of residual point sources, which makes intuitive sense, considering the lower angular resolution of this map (see Sect. 5). However, the ΛCDM parameter uncertainties derived from the four codes are very consistent. This indicates that most of the cosmological information content above ℓ ≥ 1500 is degenerate with the extragalactic foreground model, and a more sophisticated foreground treatment is required in order to recover significant cosmological parameter constraints from these scales. Beyond this, deviations among cosmological parameters are small and within 1σ for all methods and most of the parameters. Further, the parameters derived from the four foreground-cleaned CMB maps are in good agreement with those obtained by CamSpec and Plik using cross-spectra; departures are well within 1σ for most parameters.

Inspecting the differences between the best-fit models derived from the four foreground-cleaned maps and from CamSpec plotted in Fig. 12, we find that the relative residuals are within 40 μK² for all multipole ranges, and smaller than 20 μK² at high ℓ. This can be compared to the corresponding residuals for the FFP6 simulation shown in Appendix E.

The likelihood used for this analysis does not take into account some systematic effects that will affect our foreground-cleaned CMB maps, such as relative calibration uncertainties between the frequency channel maps used to construct them, or their beam uncertainties. These effects are accounted for in the likelihoods in Planck Collaboration XV (2014). We have also adopted a very simple two-parameter model for the residual extragalactic foregrounds. Despite these limitations, the four CMB maps yield cosmological parameters in agreement with the cross-spectrum based likelihoods for a basic six-parameter ΛCDM model. Thus we can be confident that the CMB maps are consistent with the power spectrum analysis.

Before concluding this section, we show in Fig. 13 the results from a two-parameter fit of an amplitude-tilt model to each of the four foreground-cleaned maps, downgraded to 6° and repixelized at an N_side = 32 grid. Clearly, the maps are virtually identical on large angular scales measured relative to cosmic variance, with any differences being smaller than 0.1σ in terms of cosmological parameters. However, it is worth noting that the best-fit model, (q,n) = (1,0), is in some tension with the low-ℓ spectrum, at about 1.7σ in this plot. The same tension between large and small angular scales is observed in Planck Collaboration XV (2014) and Planck Collaboration XVI (2014) with higher statistical significance using the full Planck likelihood. Irrespective of physical interpretation, the calculations presented here demonstrate that these low-ℓ features are robust with respect to component separation techniques.

7. Higher-order statistics

The foreground-cleaned CMB maps presented in this paper are used as inputs for most Planck analyses of higher-order statistics, including non-Gaussianity studies (Planck Collaboration XXIV 2014), studies of statistical isotropy (Planck Collaboration XXIII 2014), gravitational lensing by large-scale structure (Planck Collaboration XVII 2014), and of the ISW effect (Planck Collaboration XXIV 2014). In this section we provide a summary of the non-Gaussianity and gravitational lensing results.

7.1. Non-Gaussianity

Primordial non-Gaussianity is typically constrained in terms of the amplitude, ${\mathit{f}}_{\mathrm{NL}}^{\mathrm{local}}$, of the quadratic corrections to the gravitational potential, as well as by means of the three-point correlation function based on different triangle configurations. The results from these calculations for the foreground-cleaned CMB maps are presented in Planck Collaboration XXIV (2014). After subtraction of the lensing-ISW correlation contribution, the final result is ${\mathit{f}}_{\mathrm{NL}}^{\mathrm{local}}\mathrm{=}\mathrm{2.7}\hspace{0.17em}\mathrm{\pm }\hspace{0.17em}\mathrm{5.8}$, as estimated from the SMICA map using the KSW bispectrum estimator (Komatsu et al. 2005), consistent within 1σ with results from other methods and foreground-cleaned maps.

Uncertainties are evaluated by means of the FFP6 simulations, and potential biases are studied using both Gaussian and non-Gaussian CMB realizations. In particular, when a detectable level of primordial non-Gaussianity (${\mathit{f}}_{\mathrm{NL}}^{\mathrm{local}}\mathrm{=}\mathrm{20.4075}$) is injected into the FFP6 simulations, each foreground-cleaned map yielded a positive detection within 2σ of the expected value, recovering values of ${\mathit{f}}_{\mathrm{NL}}^{\mathrm{local}}\mathrm{=}\mathrm{8.8}\mathrm{\pm }\mathrm{8.6}$, 19.0 ± 7.5, 11.1 ± 7.6 and 19.7 ± 7.4 for Commander-Ruler, NILC, SEVEM, SMICA, respectively. We see that NILC and SMICA demonstrate the best recovery of the injected non-Gaussianity, and we favoured the latter for non-Gaussian studies for its faster performance over NILC. The foreground-cleaned CMB maps presented in this paper do not provide significant evidence of a non-zero value of ${\mathit{f}}_{\mathrm{NL}}^{\mathrm{local}}$, and realistic simulations show that the component separation methods do not suppress real non-Gaussian signatures within expected uncertainties. The implications of these results in terms of early Universe physics are discussed in the relevant papers (Planck Collaboration XXIV 2014; Planck Collaboration XXII 2014).

7.2. Gravitational lensing by large-scale structure

Gravitational lensing by the intervening matter imprints a non-Gaussian signature in the CMB, which allows the reconstruction of the gravitational potential integrated along the line of sight to the last scattering surface. In Planck Collaboration XVII (2014), this effect has been detected at a high significance level (greater than 25σ) using the Planck temperature maps. Specifically, the lensing induced correlations between the total intensity and its gradients have been used to reconstruct a nearly full sky map of the lensing potential φ, which has been used for further studies on Planck data, including the detection of a non-zero correlation with the ISW (Planck Collaboration XXIV 2014; Planck Collaboration XIX 2014) and other tracers of large-scale structure (notably, significant correlation with the CIB is reported in Planck Collaboration XVIII 2014), as well as the estimate of the power spectrum of the lensing potential and the associated likelihood. The latter was constructed using a simple minimum variance combination of the 143 and 217 GHz maps on about 70% of the sky, as well as subtracting dust contamination using the 857 GHzPlanck channel as a template (Planck Collaboration XVII 2014). These lensing results have improved the cosmological constraints from Planck (Planck Collaboration XVI 2014).

The foreground-cleaned CMB maps described in Sect. 5 were used to perform a lensing extraction on a larger sky fraction, reaching about 87% of the sky. We found the lensing power spectrum to be in good agreement with the one obtained using the minimum variance combination, i.e., the signal agrees within 1σ in the majority of the angular domain bins, and is characterized by an equivalent uncertainty. The foreground-cleaned maps were further exploited on the baseline 70% sky fraction for assessing the robustness of the main reconstruction against the foreground contamination (Planck Collaboration XVII 2014).

We show that the component separation algorithms presented in this paper do not bias the lensing reconstruction in the case of the large sky fraction considered here. We consider FFP6 simulations including noise and lensed CMB signal, propagated through each of the component separation algorithms described in Sect. 3. We perform a lensing potential reconstruction in the pixel domain based on the CMB maps processed by the four component separation methods using the metis algorithm described in Planck Collaboration XVII (2014). This method uses the quadratic estimator presented in Okamoto & Hu (2003), which corrects for the mean-field bias caused by extra sources of statistical anisotropy in addition to the CMB.

For each method, we combine the masks of CO regions, nearby galaxies and compact objects as defined in Planck Collaboration XVII (2014), with the CG90 mask described in Sect. 4. This procedure results in masks with sky fractions f_sky = 0.836,0.851,0.850,0.846 for Commander-Ruler, NILC  SEVEM, and SMICA, respectively.

We estimate the lensing potential power spectrum, ${\mathit{C}}_{\mathit{L}}^{\mathit{\phi \phi }}$, following the methodology described in Planck Collaboration XVII (2014). It consists of a pseudo-C_ℓ estimate based on a highly-apodized version of the lensing potential reconstruction, which has an effective available sky fraction f_sky,2 = 0.648,0.690,0.686,0.683 for Commander-Ruler, NILC, SEVEM and SMICA, respectively. The band-power reconstructions in 17 bins in the range 2 ≤ ℓ ≤ 1025 are plotted in Fig. 14, as well as the residuals relative to the theoretical lens power spectrum. All algorithms achieved an unbiased estimation of the underlying lensing power spectrum, with χ² = 10.58,17.34,18.54,15.30, for Commander-Ruler, NILC, SEVEM, and SMICA respectively, with 17 degrees of freedom. The associated probability-to-exceed (PTE) values are 83%, 36%, 29%, 50%.

The power spectrum estimates are in remarkable agreement with each other. However, the Commander-Ruler solution has significantly larger uncertainties, as expected from its lower signal-to-noise ratio to lensing due to its larger beam. These results on simulated foreground-cleaned CMB maps demonstrate that the component separation algorithms do not alter the lensing signal, and this provides a strategy for achieving a robust lensing reconstruction on the largest possible sky coverage. The foreground-cleaned maps have been used in Planck Collaboration XVII (2014) to obtain lensing potential estimates on 87% of the sky.

8. Foreground components

In this section we consider the diffuse Galactic components, and present full-sky maps of thermal dust and CO emission, as well as a single low-frequency component map representing the sum of synchrotron, AME, and free-free emission. Our all-sky CO map is a “type 3” product as presented in Planck Collaboration XIII (2014). To assess the accuracy of these maps, we once again take advantage of the FFP6 simulation. The Commander-Ruler method used in the following is described in Appendix A and consists of a standard parametric Bayesian MCMC analysis at low angular resolution, followed by a generalized least-squares solution for component amplitudes at high resolution.

8.1. Data selection and processing

We only use the seven lowest Planck frequencies, from 30 to 353 GHz. The two highest channels have significantly different systematic properties than the lower frequency bands, for instance concerning calibration, ZLE, and noise correlations, and they are more relevant to thermal dust and CIB studies than to the present CMB analysis. Studies of specific foregrounds (CO, thermal dust, CIB etc.) using all Planck frequencies as well as ancillary data are discussed in companion papers Planck Collaboration VI (2014), Planck Collaboration XXI (2014), Planck Collaboration XI (2014), and additional future publications will consider extensions to AME, synchrotron and free-free emission.

In order to obtain unbiased estimates of the spectral parameters across all frequency bands, each map is downgraded from its native resolution to a common angular resolution of 40′ and repixelized at N_side = 256, a limit imposed by the LFI 30 GHz channel. Once the spectral indices have been determined, we re-estimate the component amplitudes at native Planck resolution (see Appendix A).

Although the smoothing operation introduces noise correlations between pixels, we model the noise of the smoothed maps as uncorrelated white noise with an effective standard deviation, σ(p), for each pixel p. This approximation does not bias the final solution, because the analysis is performed independently for each pixel. However, it is important to note that correlations between pixels are not taken into account in this analysis. The effective noise uncertainty, σ(p), is estimated using realistic noise simulations downgraded in the same way as the data. The measured instrumental bandpasses are taken into account by integrating the emission laws over the bandpass for each component at each Monte Carlo step in the analysis.

The monopole (zero-point) of each frequency map is not constrained by Planck, but is rather determined by post-processing, and associated with a non-negligible uncertainty (see Table 5 of Planck Collaboration I 2014). In addition, each frequency map includes a significant monopole contribution from isotropic extragalactic sources and CIB fluctuations not traced by local Galactic structure, ranging from less than about 10–20 μK at 70 GHz to several hundreds of  μK at 353 GHz. Finally, the effective dipole in each map is associated with significant uncertainty due to the large kinematic CMB dipole. In order to prevent these effects from introducing modelling errors during component separation, they must be fit either prior to or jointly with the Galactic parameters. Unfortunately, when allowing free spectral parameters per pixel, there is a near-perfect degeneracy among the offsets, the foreground amplitudes and the spectral indices, and in order to break this degeneracy, it is necessary to reduce the number of spectral degrees of freedom.

We adopt the method described by Wehus et al. (in prep.) for this purpose, which has the additional advantage of making minimal assumptions about the foreground spectra. In short, this method uses linear regression between data from CMB-subtracted maps evaluated on pixels falling within each large N_side = 8 pixel to estimate the relative offsets, m₁ and m₂, between any two maps at each position on the sky. Each regression provides a constraint of the form m₁ = am₂ + b, where a and b are the slope and offset, respectively, and where each value of m_i consists of the sum of both a monopole and a dipole term evaluated at that position. The individual monopoles and dipoles can then be reconstructed by measuring a and b in different regions of the sky, exploiting spatial variations in spectral indices, and solving jointly for two monopoles and dipoles, including constraints from all positions. To minimize degeneracies, a positivity prior is imposed on the fit, such that statistically significant negative pixels are heavily penalized. For 44 and 70 GHz, we retain the dipole values determined during the mapmaking process, and do not attempt to fit them.

The resulting complete set of monopole and dipole values is listed in Table 2. As a cross-check, we performed a dedicated Commander run in which we fitted for the dipole at 353 GHz, together with the foreground amplitudes and spectral indices, and only found sub- μK differences. This channel is by far the most problematic in our data set in terms of offset determination, because of the very bright dust emission at this frequency. As a result, there is a large relative uncertainty between the zero-level of the dust amplitude map and the 353 GHz channel offset not accounted for in the following analyses. However, the sum of the two terms is well determined, and a potential error in either therefore does not compromise the quality of the other signal components (e.g., CMB and low-frequency components). A similar comment applies between the offset at 30 GHz and the zero-level of the low-frequency component, although at a significantly lower level. These degeneracies can only be broken by including additional observations at lower and higher frequencies, respectively, and this will be addressed in a future publication.

8.2. Component models and priors

Our model for the low-resolution CMB analysis includes four independent physical components: CMB; “low-frequency” emission; CO emission; and thermal dust emission. It can be written schematically in the form $\begin{array}{ccc}{\mathit{s}}_{\mathit{\nu }}\mathrm{\left(}\mathit{p}\mathrm{\right)}& \mathrm{=}& {\mathit{A}}_{\mathrm{CMB}}\mathrm{\left(}\mathit{p}\mathrm{\right)}\mathrm{+}{\mathit{A}}_{\mathrm{lf}}\mathrm{\left(}\mathit{p}\mathrm{\right)}{\left(\frac{\mathit{\nu }}{{\mathit{\nu }}_{\mathrm{0}\mathit{,}\mathrm{lf}}}\right)}^{{\mathit{\beta }}_{\mathrm{lf}}\mathrm{\left(}\mathit{p}\mathrm{\right)}}\\ & & \end{array}$(1)where A_i(p) denotes the signal amplitude for component i at pixel p, ν_0,i is the reference frequency for each component, and ν refers to frequency. (Note that for readability, integration over bandpass, as well as unit conversions between antenna, flux density and thermodynamic units, is suppressed in this expression.) Thus, each component is modelled with a simple frequency spectrum parameterized in terms of an amplitude and a small set of free spectral parameters (a power-law index for the low-frequency component, and an emissivity index and temperature for the thermal dust component); no spatial priors are imposed. One goal of the present analysis is to understand how well this simple model captures the sky signal in terms of effective components over the considered frequency range, and we exploit the FFP6 simulation (see Sect. 4) for this purpose.

In order to take into account the effect of bandpass integration, each term in the above model is evaluated as an integral over the bandpass as described in Sect. 3 of Planck Collaboration IX (2014), and converted internally to thermodynamic units. Accordingly, the reference frequencies in Eq. (1) are computed as effective integrals over the bandpass, such that the amplitude map, A_i(p), corresponds to the foreground map observed by the reference detector, i.e., after taking into account the bandpass. In order to minimize degeneracies between the different signal components, the reference band for a given component is set to the frequency at which its relative signal-to-noise ratio is maximized.

The foreground model defined in Eq. (1) is motivated by prior knowledge about the foreground composition over the CMB frequencies as outlined in Sect. 2, as is our choice of priors. In addition to the Jeffreys prior⁵ (Eriksen et al. 2008), we adopt Gaussian priors on all spectral parameters with centre values and widths attempting to strike a balance between prior knowledge and allowing the data to find the optimal solution. Where needed, we have also run dedicated analyses, either including particular high signal-to-noise ratio subsets of the data or using a lower resolution parameterization to increase the effective signal-to-noise in order to inform our prior choices. We now consider each foreground component in turn, and note in passing that the CMB component, by virtue of being a blackbody signal, is given by a constant in thermodynamic temperature units.

We approximate the low-frequency component by a straight power law in antenna temperature with a free spectral index per pixel, and adopt a prior of β = −3 ± 0.3 (this is the index in terms of antenna temperature). This choice is determined by noting that the prior is in practice only relevant at high Galactic latitudes where the signal-to-noise ratio is low and the dominant foreground component is expected to be synchrotron emission; in the signal-dominated and low-latitude AME and free-free regions, the data are sufficiently strong to render the prior irrelevant. For validation purposes, we have also considered minor variations around this prior, such as β = −2.9 ± 0.3 and β = −3.05 ± 0.2, finding only small differences in the final solutions. The reference band for the low-frequency component is set to 30 GHz, where the low-frequency foreground signal peaks. The final low-frequency amplitude map is provided in units of thermodynamic microkelvin.

The CO emission is modelled in terms of a single line ratio for each frequency. Specifically, the CO amplitude is normalized to the 100 GHz band, and defined in units of K km s^-1 (Planck Collaboration XIII 2014). The amplitude at other frequencies is determined by a single multiplicative factor relative to this, with a numerical value of 0.595 at 217 GHz and 0.297 at 353 GHz; all other frequencies are set to zero. These values are obtained from a dedicated CO analysis that includes only high signal-to-noise ratio CO regions covering a total of 0.5% of the sky. The derived values are in good agreement with those presented by Planck Collaboration XIII (2014).

Thermal dust emission is modelled by a one-component modified blackbody emission law with a free emissivity spectral index, β_d, and dust temperature, T_d, per pixel. However, since we only include frequencies below 353 GHz, the dust temperature is largely unconstrained in our fits, and we therefore impose a tight prior around the commonly accepted mean value of T_d = 18 ± 0.05 K. The only reason we do not fix it completely to 18 K is to allow for modelling errors near the Galactic centre. The dust emissivity prior is set to β_d = 1.6 ± 0.3, where the mean is determined by a dedicated run fitting for a single best-fit value for the high-latitude sky, where the prior is relevant. The reference band for the thermal dust component is 353 GHz, and the final map is provided in units of megajansky per steradian.

8.3. Results and validation

The output of the Bayesian component separation algorithm is a set of samples drawn from the joint posterior distribution of the model parameters, as opposed to a single well-defined value for each. For convenience, we summarize this distribution in terms of posterior mean and standard deviation maps, computed over the sample set, after rejecting a short burn-in phase. The goodness-of-fit is monitored in terms of the χ² per pixel. Although convenient, it is, however, important to note that this description does not provide a comprehensive statistical representation of the full posterior distribution, which is intrinsically non-Gaussian. One should be careful about making inferences in the low signal-to-noise regime based on this simplified description.

The low-resolution Commander posterior mean amplitude maps are shown in Fig. 15 for the low-frequency, CO, and thermal dust components, and the spectral index maps in Fig. 16. The associated χ² map is plotted in Fig. 17. Note that because we are sampling from the posterior instead of searching for the maximum-likelihood point, the expected number of degrees of freedom is equal to N_band = 7 in this plot, not N_band − N_par.

Several features can be seen here, foremost of which is that the Galactic plane is strikingly obvious, with χ² values exceeding 10⁴ for seven degrees-of-freedom in a few pixels. This is not surprising, given the very simplified model at low frequencies (i.e., a single power law accounting for AME, synchrotron, and free-free emission), as well as the assumption of a nearly constant dust temperature of 18 K. Second, there is an extended region of moderately high χ² roughly aligned with great circles going through the ecliptic poles, indicating the presence of correlated noise in the scanning rings not accounted for in our white noise model.

Based on these, and other considerations, it is clear that parts of the sky must be masked before proceeding to CMB power spectrum and likelihood analyses. This masking process is discussed at greater length in Planck Collaboration XV (2014), and results in different masks for specific applications. The goal of our present discussion is to evaluate the adequacy of the mask adopted for low-ℓ likelihood analysis (L87), which is based on the fits presented here. This mask removes 13% of the sky, and is derived from a combination of χ² and component amplitude thresholding.

Figure 18 compares the high-resolution Ruler solution to the low-resolution Commander solution for CMB, CO and thermal dust on a particularly strong CO complex near the Fan region, centred on Galactic coordinates (l,b) = (110°,15°).

For validation purposes, we analyse the simulations described in Sect. 4 in the same way as the real data, including monopole and dipole determination, CO line ratio estimation and spectral index estimation. Individual component maps at each observed frequency are available from the simulation process, and used for direct comparison with the reconstructed products.

In Fig. 19 we show the differences between the recovered and input component maps at their respective reference frequencies. The boundary of the 13% Commander mask is traced by the white contours, and a best-fit monopole and dipole have been subtracted from each difference map. All difference maps are shown in units of thermodynamic  μK. The top panel of Fig. 20 gives the error histograms outside the masked region for each component, normalized to the respective estimated standard deviation; if the recovered solution has both correct mean and standard deviation, these histograms should match a Gaussian distribution with zero mean and unit variance, indicated by the dashed black line. Conversely, a significant bias would be visible as a horizontal shift in this plot, while under-estimation of the errors would result in too wide a distribution and vice versa. The bottom panel shows the fractional error (i.e., the error divided by the true input value) for all pixels with signal above 5σ; the fractional error is not a useful quantity for noisier signals.

The difference maps in Fig. 19 display significant errors in the Galactic plane. For the low-frequency component, the residuals are dominated by free-free emission, while for thermal dust the dominant contaminant is CO emission. However, outside the mask the residuals are small, and, at least for the low-frequency and CO components, the spatial characteristics appear similar to instrumental noise. This is more clear in the histograms shown in the top panel of Fig. 20; the mean and standard deviations are δ_lf = 0.01 ± 1.12, δ_CO = 0.00 ± 0.87, and δ_td = 0.00 ± 2.01, respectively, for the low-frequency, CO and thermal dust components. There is no evidence of bias outside the mask in any component, and the error estimates are accurate to 12 and 13% for the low-frequency and CO components. Note, though, that the estimated error for the CO component is actually larger than the true uncertainty, suggesting that the white noise approximation for the 100 GHz channel overestimates the true noise. This can occur if the correlated instrumental noise is important in regions where there is no significant CO emission. Locally re-scaling the white noise to account for spatially varying correlated noise would correct this effect.

For the thermal dust component, on the other hand, the error is underestimated by a factor of 2. The explanation for this is most easily seen from the lower panel of Fig. 19. This map is dominated by isotropic CIB fluctuations, rather than instrumental noise. Because these fluctuations have a slightly different spectrum than the dominant Galactic dust emission, and the model does not account for a separate CIB component, the error on the Galactic component is underestimated. When using the Galactic map presented here for detailed analysis near the noise limit, taking into account these residual fluctuations is essential, and the effective noise per pixel should be increased by a factor of 2.

As clearly seen in Fig. 19, the residuals inside the mask are highly significant in a strict statistical sense. However, as seen in the bottom panel of Fig. 20, they are relatively small in terms of fractional errors. Specifically, the three histograms have means and standard deviations of f_lf = 0.00 ± 0.10, f_CO = −0.03 ± 0.10, and f_td = 0.00 ± 0.06, respectively, for the low-frequency, CO and thermal dust components. The largest bias is observed for the CO component, for which the absolute amplitude is biased by 3%. The bias in the low-frequency and thermal dust components is negligible, and the fractional uncertainties are 10 and 6%, respectively. This confirms that approximating the sum of the three low-frequency components by a single power-law over the Planck frequency bands is reasonable; if modelling errors dominated, one would expect to see a significant bias in the resulting amplitude.

In order to validate the spectral parameters, we show in Fig. 21 histograms of the normalized residuals for each foreground component evaluated at its two leading sub-dominant frequencies (i.e., at 44 and 70 GHz for the low-frequency component; at 217 and 353 GHz for the CO component; and at 143 and 217 GHz for the thermal dust component). The means and standard deviations of these distributions are: δ_lf(44 GHz) = −0.41 ± 1.98 and δ_lf(70 GHz) = −0.34 ± 2.04 for the low-frequency component; δ_CO(217 GHz) = 0.10 ± 0.84 and δ_CO(353 GHz) = 0.51 ± 1.00 for the CO component; and δ_td(143 GHz) = −0.02 ± 1.53 and δ_td(217 GHz) = −0.13 ± 1.87 for the thermal dust component. As expected, the effect of modelling errors is more significant at the sub-dominant frequencies than at the pivot frequencies, when measured in terms of statistical uncertainties, since the foreground signal is weaker and the confusion with the other components relatively larger. Nevertheless, we see that the absolute bias is at most 0.5σ for the CO component at 353 GHz, while the thermal dust bias is negligible even at 143 GHz. The estimated uncertainties are generally underestimated by up to a factor of two due to these modelling errors.

Finally, the efficiency of the adopted foreground model for CMB analysis is quantified in Appendix E in terms of power spectrum residuals and cosmological parameter estimation.

To summarize, we find that the simplified model, defined by Eq. (1), provides a good fit to the realistic FFP6 simulation for most of the sky. Absolute residuals are small, and the amplitude uncertainty estimates are accurate to around 12%, except for the thermal dust component for which unmodelled CIB fluctuations are important. Further, we find that the real Planck data behave both qualitatively and quantitatively very similarly to the FFP6 simulation, suggesting that this approach also performs well on the real sky.

8.4. Interpretation and comparison with other results

The maps shown in Figs. 15 and 16 provide a succinct summary of the average foreground properties over the Planck frequency range. We now consider their physical interpretation and compare them to products from alternative methods.

First, the top panel of Fig. 22 shows a difference map between the thermal dust map⁶ derived in the present paper using Planck frequencies between 30 and 353 GHz, and the one determined from only the three highest Planck frequencies (353 to 857 GHz) and the 100  μm IRIS map by Planck Collaboration XI (2014). The large-scale features in this map are dominated by the different CO, zodiacal light and offset modelling approaches adopted by the two pipelines. Specifically, while the Commander solution jointly fits for CO and thermal dust, the high-frequency solution assumes CO to be negligible at these frequencies; while Commander analyzes raw Planck maps with ZLE still present, the high-frequency analysis considers maps that have been explicitly corrected for zodiacal light using the Planck Collaboration VI 2014 model (which itself has a poorly defined zero-level); and while the Commander solution fits internally for spurious monopoles and dipoles by regression as described in Sect. 8.1, the high-frequency solution determines these through correlations with external HI observations. The middle panel of Fig. 22 shows the difference between the two solutions after subtracting residual monopoles, dipoles and zodiacal light components. The differences are now below 0.01 MJy sr^-1 everywhere at high latitudes. This is to be compared with the corresponding thermal dust amplitude map in the bottom panel of Fig. 15, which ranges between 0 and 2.5 MJy sr^-1. Inside the Galactic plane, the differences are dominated by residuals due to different CO modelling, seen as solid blue colours in Fig. 22; however, even in this region the differences are smaller than 5% of the amplitude. Thus, the agreement between the two solutions is excellent, despite their very different choice of data and processing methodology. The bottom panel shows the T–T plot between the two maps, excluding any pixel for which the Commander CO amplitude at 100 GHz is larger than 1 K km s^-1.

From Fig. 16, we see that the dust emissivity ranges between 1.3 and 1.7 for most of the sky, in good agreement with the mean value of β_d = 1.56 ± 0.02 derived by Planck Collaboration Int. XXII (2014) by averaging over 10° disks and including additional ancillary data. Considering only the pixels with a posterior distribution width that is a third of the prior width (i.e., σ(β_d) < 0.1), we find a mean value of 1.49. The two exceptions are a large region of shallow indices northeast of the Galactic centre, and steep indices near the Galactic plane. The former region corresponds to a part of the sky with low dust emission, where we expect the spectral index to be sensitive to both monopole and dipole residuals, as well as instrumental systematics, such as correlated 1 /f noise. The latter appears to be particularly pertinent here because the shallow index region at least partially traces the Planck scanning strategy; as a result, the systematic error on the spectral index in this region is considerable. The main systematic uncertainty connected to the region of steep indices around the Galactic plane is confusion with CO emission.

The CO map shown in Figs. 15 and 18 is discussed in greater detail in Planck Collaboration XIII (2014). A distinct advantage of this particular solution over available alternatives is its high signal-to-noise ratio per pixel, which is achieved by reducing all information into a single value per pixel. Consequently, this map serves as a unique tool for follow-up CO observations. However, the assumption of a constant line ratio over the full sky may lead to a significant systematic uncertainty on CO amplitude per pixel.

Finally, the spectral index map for the low-frequency component shown in Fig. 16 can be used to determine the dominant low-frequency component (synchrotron, free-free or AME) as a function of position on the sky. To illustrate this connection, we once again take advantage of the FFP6 simulation for which we know the amplitude of each low-frequency component per pixel. In the top panel of Fig. 23, we use this information to make a “dominant component map”; dark blue indicates that synchrotron emission is strongest at a given pixel, light blue that free-free is strongest, and orange that spinning dust (AME) dominates. In the bottom panel, we show our derived power-law index map from the same simulation. As expected, the correspondence between the power-law index and the dominant component is very strong, implying that the spectral index map can be used to trace the individual components. In particular, we see that an index below about − 3.3 reflects the presence of a component consistent with a spinning dust model peaking below 20 GHz over the Planck frequency range, while an index higher than around − 2.3 signals the importance of free-free emission. Intermediate values typically indicate synchrotron emission, although it should be noted that the signal-to-noise ratio at very high latitudes is low and the results are therefore prior-driven in these regions.

Returning to the spectral index map shown in Fig. 16, we see a good correspondence between the real data and the simulation. Features present in the simulation also appear in the data. For instance, we see that the spectral index in the so-called Fan Region (i.e., near Galactic coordinates (l,b) = (90°,20°)) is low in both cases, and this alone provides strong evidence for the presence of AME. Further, the AME spectral index is consistent with the spinning dust interpretation. This power-law index map may be used to identify particular AME regions for follow-up observations. Finally, we note that regions known for strong free-free emission, such as the Gum Nebula or Zeta Ophiuchi, have spectral indices close to − 2.1 or − 2.2, as expected.

9. Conclusions

We have tested four component separation algorithms on the Planck frequency maps to produce clean maps of the CMB anisotropies over a large area of sky. These CMB maps are used for studies of statistics and isotropy (Planck Collaboration XXIII 2014), primordial non-Gaussianity (Planck Collaboration XXIV 2014), gravitational lensing (Planck Collaboration XVII 2014), the ISW effect (Planck Collaboration XIX 2014), cosmic geometry and topology (Planck Collaboration XXVI 2014), searches for cosmic defects from primordial phase transitions (Planck Collaboration XXV 2014), as well as an integral part of the low-ℓPlanck likelihood (Planck Collaboration XV 2014). Two of the methods, one using internal foreground templates (SEVEM) and the other an ILC in needlet space (NILC), are non-parametric, extracting the CMB map by minimizing the variance of the total contamination. The other two methods fit models of the foregrounds to clean the CMB of their emission. One fits a parametric model in real space (C-R) and one fits a non-parametric in the harmonic domain (SMICA).

All four methods have been demonstrated to work well both on real and simulated data, and to yield consistent results. Nevertheless, there are differences between the methods, making them more or less suitable for specific applications. For instance, Commander-Ruler allows a joint parametric foreground estimation and CMB power spectrum estimation, with full propagation of foreground uncertainties to cosmological parameters, but is limited to a lower angular resolution than the other codes. This method has therefore been selected for the low-ℓPlanck likelihood (Planck Collaboration XV 2014) and to produce astrophysical component maps (Sect. 8), while it is sub-optimal for applications requiring full angular resolution, e.g., gravitational lensing reconstruction or estimation of primordial non-Gaussianity. For these purposes, we use the three higher-resolution maps. We take SMICA to be the leading method, based on its superior performance on the FFP6 simulation, where it has be shown to have the lowest residual foreground contamination at large scales and to preserve primordial non-Gaussianity. When subjecting foreground-cleaned Planck maps to scientific analysis, we use the other two or three maps, as appropriate, to assess the uncertainties inherent in the choice of methods and the assumptions they make. Indeed, this is the main purpose for presenting four different CMB solutions to the general community.

The CMB anisotropies are robustly recovered over a large fraction (73%) of the sky and down to small angular scales, reaching to multipoles ℓ ≈ 2000. We characterize the CMB maps with angular power spectra and cosmological parameter constraints. Parameter constraints from these maps are consistent with those from the Planck likelihood function based on cross-spectra and large sky cuts (Planck Collaboration XV 2014). This agreement supports the robustness of both our component separation methodology and cosmological parameter constraints.

The real-space parametric fits of Commander-Ruler enable us to characterize the diffuse Galactic foregrounds. We parameterize them with a low-frequency power-law component, representing the sum of synchrotron, free-free, and AME emission, a high-frequency modified blackbody spectrum describing thermal dust emission, and a molecular CO component. Using only the Planck data from 30 to 353 GHz, we fit for the amplitude and spectral parameters of the three foregrounds and the CMB simultaneously at each pixel of a 40-arcmin resolution map. The spectral parameters are the low-frequency component power-law exponent and the modified blackbody emissivity power-law exponent; the CO line ratios are spatially fixed. These parameters give us the source mixing matrix, which we then use in a direct inversion to deduce the component amplitudes at higher resolution. Through Gibbs sampling, we obtain realizations drawn from the full posterior distribution of possible foreground and CMB solutions, giving us a powerful ability to statistically characterize our results.

Our in-depth analysis of the recovered CMB anisotropies is unprecedented for component separation studies, concerning both the accuracy of cosmological parameter constraints, and studies of early Universe physics and structure formation through gravitational lensing. On the other hand, the complex nature of the foreground emission over such a large frequency range limits us to the use of relatively simple methods when analysing Planck data alone. An extensive study in combination with other probes of Galactic foregrounds will be presented in forthcoming papers. In particular, the separation of individual components at low frequencies requires the use of ancillary data, for example, from the WMAP and radio surveys.

---

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, JA 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); 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.sciops.esa.int/index.php?project=planck&page=Planck_Collaboration. The authors acknowledge the support provided by the Advanced Computing and e-Science team at IFCA. This work made use of the COSMOS supercomputer, part of the STFC DiRAC HPC Facility. Some of the results in this paper have been derived using the HEALPix package.

References

Appendix A: Physical parametrization

The Commander-Ruler (C-R) approach implements Bayesian component separation in pixel space, fitting a parametric model to the data by sampling the posterior distribution for the model parameters. For computational reasons, the fit is performed in a two-step procedure: first, both foreground amplitudes and spectral parameters are found at low-resolution using Markov Chain Monte Carlo (MCMC)/Gibbs sampling algorithms (Jewell et al. 2004; Wandelt et al. 2004; Eriksen et al. 2004, 2007, 2008); second, the amplitudes are recalculated at high resolution by solving the generalized least squares system (GLSS) per pixel with the spectral parameters fixed to their values from the low-resolution run.

For the CMB-oriented analysis presented in this paper, we only use the seven lowest Planck frequencies, i.e., from 30 to 353 GHz. We first downgrade each frequency map from its native angular resolution to a common resolution of 40 arcminutes and repixelize at HEALPixN_side = 256. Second, we set the monopoles and dipoles for each frequency band as described in Sect. 8.1, using a method that locally records spectral indices (Wehus et al., in prep.). We approximate the effective instrumental noise as white with a root mean square (rms) per pixel given by the Planck scanning pattern and an amplitude calibrated by smoothing simulations of the instrumental noise, including correlations, to the same resolution. For the high-resolution analysis, the important pre-processing step is the upgrading of the effective low-resolution mixing matrices to full Planck resolution: this is done by repixelizing from N_side = 256 to 2048 in harmonic space, ensuring that potential pixelization effects from the low-resolution map do not introduce sharp boundaries in the high-resolution map.

Our model for the data, a map d_ν at frequency ν, consists of a linear combination of N_c astrophysical components, plus instrumental noise, $\begin{array}{ccc}{{d}}_{\mathit{\nu }}\mathrm{=}\sum _{\mathit{i}\mathrm{=}\mathrm{1}}^{{\mathit{N}}_{\mathrm{c}}}{\mathrm{F}}_{\mathit{\nu }}^{\mathit{i}}\mathrm{\left(}\mathit{\theta }\mathrm{\right)}\mathrm{·}{{A}}^{\mathit{i}}\mathrm{+}{{n}}_{\mathit{\nu }}\mathit{,}& & \end{array}$(A.1)where Aⁱ denotes a sky map vector containing the foreground amplitude map for component i normalized at a reference frequency, and ${\mathrm{F}}_{\mathit{\nu }}^{\mathit{i}}\mathrm{\left(}{\mathit{\theta }}^{\mathit{i}}\mathrm{\right)}$ is a diagonal matrix describing the spectral emission law for component i as a function of frequency and which depends on a (small) set of spectral parameters, θ. The CMB signal is included in the sum and, as a special case, it may be represented either in harmonic or pixel space, depending on whether the main goal of the analysis is CMB power spectrum analysis or component separation. The former representation is used for the Commander-based low-ℓ likelihood presented in Planck Collaboration XV (2014), while the latter is used for the foreground fits presented in Sect. 8 of this paper.

Bayes theorem specifies the posterior distribution for the model parameters, $\begin{array}{ccc}\mathit{P}\mathrm{\left(}{{A}}^{\mathit{i}}\mathit{,\theta }\mathrm{\right|}{d}\mathrm{)}\mathrm{\propto }\mathrm{ℒ}\mathrm{(}{{A}}^{\mathit{i}}\mathit{,\theta }\mathrm{\left)}\mathit{P}\mathrm{\right(}{{A}}^{\mathit{i}}\mathit{,\theta }\mathrm{)}\mathit{,}& & \end{array}$(A.2)where ℒ(Aⁱ,θ) = P(Aⁱ,θ | d) is a Gaussian likelihood of observing data, d, given model parameters (Aⁱ,θ), and the prior P(Aⁱ,θⁱ) depends on the application. In this paper, the prior on spectral indices is a product of a Jeffreys prior and physical priors, as detailed in Sect. 8.2; no priors are imposed on the foreground amplitudes.

In the low-resolution Commander analysis, we exploit a Gibbs sampler to map out the posterior distribution (Eriksen et al. 2008), adopting the following minimal two-step scheme: $\begin{array}{ccc}{A}& \mathrm{\leftarrow }\mathit{P}\mathrm{\left(}{A}\mathrm{\right|}\mathit{\theta ,}{d}\mathrm{)}\mathrm{;}& \\ \mathit{\theta }& \mathrm{\leftarrow }\mathit{P}\mathrm{\left(}\mathit{\theta }\mathrm{\right|}{A}\mathit{,}{d}\mathrm{)}\mathit{.}& \end{array}$The first conditional distribution is a multivariate Gaussian, while the second distribution does not have an analytic form and must be mapped out numerically.

The high-resolution Ruler analysis maximizes the foreground amplitude conditional in Eq. (A.3) numerically by solving the generalized least squares system $\begin{array}{ccc}{{A}}^{\mathit{i}}\mathrm{\left(}\mathit{\theta }\mathrm{\right)}\mathrm{=}(\sum _{\mathit{\nu }}{\mathrm{F}}_{\mathit{\nu }}^{\mathit{i,}\mathrm{T}}{\mathrm{N}}_{\mathit{\nu }}^{-1}{\mathrm{F}}_{\mathit{\nu }}^{\mathit{i}}{)}^{-1}\sum _{\mathit{\mu }}{\mathrm{F}}_{\mathit{\mu }}^{\mathit{i,}\mathrm{T}}{\mathrm{N}}_{\mathit{\mu }}^{-1}{{d}}_{\mathit{\mu }}\mathrm{\equiv }\sum _{\mathit{\mu }}\mathrm{W}\mathrm{(}\mathit{\theta }{\mathrm{)}}_{\mathit{\mu }}^{\mathit{i}}{{d}}_{\mathit{\mu }}\mathit{,}& & \end{array}$(A.5)were N_ν is the noise covariance matrix (assumed to be diagonal) of the νth channel, Fⁱ ≡ Fⁱ(θ) and we have neglected the different angular resolutions of the channels. The posterior marginal average for the high-resolution amplitude maps is then given by $\mathrm{⟨}{{A}}^{\mathit{i}}\mathrm{⟩}\mathrm{=}{}^{\mathrm{\int }}{A}^{\mathit{i}}\mathrm{\left(}\mathit{\theta }\mathrm{\right)}\mathit{P}\mathrm{\left(}\mathit{\theta }\mathrm{\right)}\mathrm{d}\mathit{\theta }\mathrm{\simeq }\mathrm{\left(}{\sum }_{\mathit{\theta ,\nu }}\mathrm{W}\mathrm{\right(}\mathit{\theta }{\mathrm{)}}_{\mathit{\nu }}^{\mathit{i}}{{d}}_{\mathit{\nu }}\mathrm{)}\mathit{/}{\mathit{N}}_{\mathrm{sample}}\mathrm{\equiv }{\sum }_{\mathit{\nu }}{\mathrm{W}}_{\mathit{\nu }}^{\mathit{i}}{{d}}_{\mathit{\nu }}$, a sum over the N_sample samples of the spectral parameters θ.

Once the channel weights, W_ν, have been computed, processing a large number of simulations requires negligible computational resources. This feature has been extensively used for computation of the effective beam of Ruler maps: FFP6 CMB simulations for the 30 to 353 GHz channels are combined according to W_ν and the effective beam transfer function is found as ${\mathit{b}}_{\mathit{\ell }}^{\mathrm{2}}\mathrm{\equiv }\mathrm{⟨}{\mathit{C}}_{\mathit{\ell }}^{\mathrm{out}}\mathit{/}\mathrm{\left(}{\mathit{C}}_{\mathit{\ell }}^{\mathrm{inp}}{\mathit{w}}_{\mathit{\ell }}^{\mathrm{2}}\mathrm{\right)}\mathrm{⟩}$. Here, ${\mathit{C}}_{\mathit{\ell }}^{\mathrm{inp}}$ is the power spectrum of the input simulation before convolution with the instrumental beam, w_ℓ is the HEALPix pixel window function, and the average is taken over the set of simulations. Missing pixels are set to 0 when computing C_ℓ in the above expression. The low number (~500) of missing pixels in the data renders the impact of such a choice negligible at ℓ< 2000. A similar procedure is used for defining the effective beam of the non-CMB components.

The above algorithm produces a set of samples drawn from the posterior distribution, as opposed to a direct estimate of individual component amplitudes or spectral parameters. While this sample set provides a statistically complete representation of the posterior, it is non-trivial to visualize or to compare the distribution with external data. For convenience, we therefore summarize the distribution in terms of mean and standard deviation maps for each component. We emphasize, however, that the distribution is significantly non-Gaussian, and when searching for features in the maps at low signal-to-noise levels, one must take into account the exact distribution.

Finally, the Commander-Ruler confidence mask (see Sect. 5) is primarily defined by the product of the CG80 mask and the point source mask described in Sect. 4. We additionally remove any pixels excluded by the 13% Commander likelihood mask described by Planck Collaboration XV (2014); however, this is almost entirely included within the CG80 mask, and this step therefore has very little impact on the final result. To complete the mask, we remove any pixels for which the high-resolution Ruler CMB map, smoothed to 40 arcmin, differs by more than 3σ from the low-resolution Commander CMB map, which can happen due to spatial spectral variations on pixel scales.

Appendix B: Internal linear combination

NILC is a method to extract the CMB (or any component with known spectral behaviour) by applying the ILC technique to multi-channel observations in needlet space, that is, with weights that are allowed to vary over the sky and over the full multipole range.

The ability to linearly combine input maps varying over the sky and over multipoles is called localization. In the needlet framework, harmonic localization is achieved using a set of bandpass filters defining a series of scales, and spatial localization is achieved at each scale by defining zones over the sky. The harmonic localization adopted here uses nine spectral bands covering multipoles up to ℓ = 3200 (see Fig. B.1). The spatial localization depends on the scale. At the coarsest scale, which includes the multipoles of lowest degree, we use a single zone (no localization), while at the finest scales (which include the highest multipoles) the sky is partitioned into 20 zones (again, see Fig. B.1).

The NILC method amounts to computing an ILC in each zone of each scale, allowing the ILC weights to adapt naturally to the varying strength of the other components as a function of position and multipole. A complete description of the basic NILC method can be found in Delabrouille et al. (2009).

In the present work, however, there is an important difference in the processing of the coarsest scale. Since the coarsest scale of the NILC filter is not localized, a plain NILC map would be equivalent to a pixel-based ILC for all the multipoles of that scale. This procedure, however, is known to be quite susceptible to the “ILC bias”, due to chance correlations between the CMB and foregrounds. In order to mitigate this effect, the (single) covariance matrix which determines the ILC coefficients at the coarsest scale is not computed as a pixel average, but is rather estimated in the harmonic domain as an average over spherical harmonic coefficients using a spectral weight which equalizes the power of the CMB modes (based on a fiducial spectrum). This can be shown significantly to decrease the large scale errors.

In practice, our NILC processing depends on several implementation choices, as follows:

• Input channels: in this work, the NILC algorithm isapplied to all Planck channelsfrom 44 to 857 GHz, omitting onlythe 30 GHz channel.

• Pre-processing of point sources: identical to the SMICA pre-processing (see Appendix D).

• Masking and inpainting: the NILC CMB map is actually produced in a three-step process. In a first step, the NILC weights are computed from covariance matrices evaluated using a Galactic mask removing about 2% of the sky (and is apodized at 1°). In a second step, those NILC weights are applied to needlet coefficients computed over the complete sky (except for point source masking/subtraction), yielding a NILC CMB estimate over the full sky (except for the point source mask). In short, the weights are computed over a masked sky but are applied to a full sky (up to point sources). In a final step, the pixels masked due to point source processing are replaced by the values of a constrained Gaussian realization (inpainting).

• Spatial localization: the boundaries of the zones used for spatial localisation (shown at Fig. B.1) are obtained as iso-level curves of a low resolution map of Galactic emission.

• Beam and transfer function: as in the SMICA processing, the resolution of the input maps is adjusted so that they have a beam of 5′. Therefore the resulting CMB map is automatically synthesized with an effective Gaussian beam of 5′, according to the unbiased nature of the ILC.

• Using SMICA recalibration: in our current implementation, the NILC solution uses the values determined by SMICA for the CMB spectrum, given in Eq. (D.6).

Appendix C: Template fitting

The original SEVEM algorithm produced clean CMB maps at several frequencies through template fitting, followed by an estimation of the CMB power spectrum from these clean maps using a method based on the Expectation Maximization algorithm (Martınez-Gonzalez et al. 2003; Leach et al. 2008; Fernández-Cobos et al. 2012). From this power spectrum, a multi-frequency Wiener-filtered CMB map was produced. For the present work, only the first step of the method, producing clean CMB maps at different frequencies, is considered. In addition, two of these clean maps are optimally combined to produce a final CMB map.

The templates used for cleaning are internal, i.e., they are constructed from Planck data, avoiding the need for external data sets, which usually complicate the analyses and may introduce inconsistencies. In the cleaning process, no assumptions about the foregrounds or noise levels are needed, rendering the technique very robust. The fitting can be done in real or wavelet space (using a fast wavelet adapted to the HEALPix pixelization; Casaponsa et al. 2011) to properly deal with incomplete sky coverage. For expediency, however, we fill in the small number of unobserved pixels at each channel with the mean value of their neighbouring pixels before applying SEVEM.

We construct our templates by subtracting two close Planck frequency channel maps, after first smoothing them to a common resolution, and converting to CMB temperature units if necessary, to ensure that the CMB signal is properly removed. A linear combination of the templates, t_j, is then subtracted from (hitherto unused) map, d, to produce a clean CMB map at that frequency. This is done either in real or wavelet space (i.e., scale by scale) at each position on the sky, $\begin{array}{ccc}{\mathit{T}}_{\mathrm{c}}\mathrm{\left(}{x}\mathit{,\nu }\mathrm{\right)}\mathrm{=}\mathit{d}\mathrm{\left(}{x}\mathit{,\nu }\mathrm{\right)}\mathrm{-}\sum _{\mathit{j}\mathrm{=}\mathrm{1}}^{{\mathit{n}}_{\mathit{t}}}{\mathit{\alpha }}_{\mathit{j}}{\mathit{t}}_{\mathit{j}}\mathrm{\left(}{x}\mathrm{\right)}\mathit{,}& & \end{array}$(C.1)where n_t is the number of templates. If the cleaning is performed in real space, the α_j coefficients are obtained by minimizing the variance of the clean map, T_c, outside a given mask. When working in wavelet space, the cleaning is performed in the same way at each wavelet scale independently (i.e., the linear coefficients depend on the scale). Although we exclude very contaminated regions during the minimization, the subtraction is performed for all pixels and, therefore, the cleaned maps cover the full-sky (although foreground residuals are expected to be present in the excluded areas).

An additional level of flexibility may also be considered: the linear coefficients can be fixed over the full sky, or in several regions. The regions are then combined in a smooth way, by weighting the pixels at the boundaries to reduce possible discontinuities in the clean maps.

Since the method is linear, we may easily propagate the noise properties to the final CMB map. Moreover, it is very fast and permits the generation of thousands of simulations to characterize the statistical properties of the outputs, a critical need for many cosmological applications. The final CMB map retains the angular resolution of the original frequency map.

There are several configurations possible for SEVEM, depending the number of frequency maps to be cleaned or the number of templates used in the fitting. Note that the production of clean maps at different frequencies is of great interest in order to test the robustness of the results. Therefore, to define the best strategy, one needs to find a compromise between the number of maps that can be cleaned independently and the number of templates that can be constructed.

In particular, we have cleaned the 143 GHz and 217 GHz maps using four templates constructed as the difference of the following Planck channels (smoothed to a common resolution): (30−44), (44−70), (545−353), and (857−545). For simplicity, the two maps have been cleaned in real space, since there was no significant improvement found when using wavelets, especially at high latitude. In order to take into account the different spectral behaviour of the foregrounds at low and high Galactic latitudes, we considered two independent sky regions, using different sets of coefficients (see Table C.1 for the values of the linear coefficients for the main considered region). The first region corresponds to the brightest 3% of Galactic emission, while the second region is defined by the remaining 97% of the sky (see Sect. 4 for a detailed description of our Galactic mask construction). For the first region, the coefficients are actually estimated over the complete sky (we find that this is better than performing the minimization on only the brightest 3% of the sky, where the CMB is very sub-dominant), while for the second region, we exclude the bright 3% sky fraction, point sources detected at any frequency and those pixels which have not been observed in all channels.

Note that, for consistency, we have used the same configuration (four templates, cleaning in real space, two regions) for the analysis of the FFP6 simulations. However, we find that this simple configuration produces more contaminated CMB maps than for the data (although the region outside the confidence mask still has low contamination), indicating some differences between the foreground level in the data and in simulations. Therefore, conclusions derived from the FFP6 results for SEVEM should be taken with caution, since they are expected to provide overestimated residuals.

Our final CMB map was constructed by combining the 143 and 217 GHz maps by weighting the maps in harmonic space, taking into account the noise level, the resolution, and a rough estimation of the foreground residuals of each map (obtained from realistic simulations). This final map has a resolution corresponding to a Gaussian beam of 5′ FWHM.

Moreover, additional clean CMB maps (at frequencies 44, 70, 100, and 353 GHz) were also produced using different combinations of templates. In particular, to clean the 100 GHz map, we used the same templates and regions as for 143 and 217 GHz. This allowed us to produce three (almost) independent clean maps, in the sense that none of the three maps to be cleaned was used to construct the templates. For 44, 70, and 353 GHz, different combinations of templates were used, and the linear coefficients were chosen to be the same over the full sky. They were obtained by minimizing the variance of the map outside the same mask as that used to clean the central frequency maps on the largest region. The templates and the corresponding linear coefficients used for each of the considered frequencies are given in Table C.1.

The SEVEM clean frequency maps have been used in analyses of the isotropy and statistics of the CMB (Planck Collaboration XXIII 2014) and to obtain cosmological constraints from the ISW effect (Planck Collaboration XIX 2014). In particular, clean maps from 44 to 353 GHz were used for the stacking analysis presented in Planck Collaboration XIX (2014), while frequencies from 70 to 217 GHz were used for consistency tests in Planck Collaboration XXIII (2014).

The confidence mask provided for SEVEM is constructed by combining four types of selected regions. In particular, it excludes zones with high residuals identified through the subtraction of SEVEM clean maps at different frequencies, as well as the sources detected at all Planck channels using the Mexican hat wavelet algorithm (Planck Collaboration XXVIII 2014). These point sources are masked with holes of varying radius, according to the flux of each source. The pixels that are not observed by all channels are also masked. Finally, to ensure that the area left outside the mask is statistically robust, we also exclude from the analysis the brightest 20% of Galactic emission, leaving a useful area of around 76%. This provides a conservative mask for CMB analysis; however, we point out that smaller masks could also be used in specific applications, such as the lensing potential reconstruction described in Sect. 7).

Appendix D: Spectral matching

SMICA (Spectral Matching Independent Component Analysis) reconstructs a CMB map as a linear combination in harmonic space of N_chan input frequency maps with weights that depend on multipole ℓ. Given the N_chan × 1 vector x_ℓm of spherical harmonic coefficients for the input maps, it computes coefficients ŝ_ℓm for the CMB map as $\begin{array}{ccc}\mathit{ŝ}\mathit{\ell m}\mathrm{=}{{w}}_{\mathit{\ell }}^{†}{{x}}_{\mathit{\ell m}}\mathit{,}& & \end{array}$(D.1)where the N_chan × 1 vector w_ℓ containing the multipole-dependent weights is chosen to give unit gain to the CMB with minimum variance. This is achieved with $\begin{array}{ccc}{{w}}_{\mathit{\ell }}\mathrm{=}\frac{{\mathrm{C}}_{\mathit{\ell }}^{-1}{a}}{{{a}}^{†}{\mathrm{C}}_{\mathit{\ell }}^{-1}{a}}\mathit{,}& & \end{array}$(D.2)where vector a is the spectrum of the CMB evaluated at each channel (allowing for possible inter-channel re-calibration factors) and C_ℓ is the N_chan × N_chan spectral covariance matrix of x_ℓm. Taking C_ℓ in Eq. (D.2) to be the sample spectral covariance matrix $\begin{array}{c}{\mathrm{􏽢}}^{}\\ {\mathrm{C}}_{\mathit{\ell }}\end{array}$ of the observations, $\begin{array}{ccc}\begin{array}{c}\mathrm{􏽢}\\ {\mathrm{C}}_{\mathit{\ell }}\end{array}\mathrm{=}\frac{\mathrm{1}}{\mathrm{2}\mathit{\ell }\mathrm{+}\mathrm{1}}\sum _{\mathit{m}}{{x}}_{\mathit{\ell m}}{{x}}_{\mathit{\ell m}}^{†}\mathit{,}& & \end{array}$(D.3)would implement a simple harmonic-domain ILC, similar to Tegmark et al. (2003). At the largest scales, we instead use a model, C_ℓ(θ), and determine the covariance matrix to be used in Eq. (D.2) by fitting C_ℓ(θ) to $\begin{array}{c}{\mathrm{􏽢}}^{}\\ {\mathrm{C}}_{\mathit{\ell }}\end{array}$. This is done in the maximum likelihood sense for stationary Gaussian fields, that is, the best fit matrices, , are obtained for $\begin{array}{ccc}\mathit{\theta ̂}\mathrm{=}\arg \underset{\mathit{\theta }}{\min }\sum _{\mathit{\ell }}\mathrm{(}\mathrm{2}\mathit{\ell }\mathrm{+}\mathrm{1}\mathrm{)}\left(\begin{array}{c}\mathrm{􏽢}\\ {\mathrm{C}}_{\mathit{\ell }}\end{array}{\mathrm{C}}_{\mathit{\ell }}\mathrm{(}\mathit{\theta }{\mathrm{)}}^{-1}\mathrm{+}\log \begin{array}{c}\det \end{array}{\mathrm{C}}_{\mathit{\ell }}\mathrm{(}\mathit{\theta }\mathrm{)}\right)\mathit{.}& & \end{array}$(D.4)Equations (D.1)–(D.4) summarize the basic principles of SMICA; its actual operation depends on a choice for the spectral model C_ℓ(θ), and on several implementation-specific details, which we briefly describe below.

We model the data as a superposition of CMB, noise and foregrounds. The latter are not parametrically modelled; instead, we represent the total foreground emission by d templates with arbitrary frequency spectra, angular spectra and correlations. In the spectral domain, this is equivalent to modelling the covariance matrices as $\begin{array}{ccc}{\mathrm{C}}_{\mathit{\ell }}\mathrm{\left(}\mathit{\theta }\mathrm{\right)}\mathrm{=}{a}{{a}}^{†}\hspace{0.17em}{\mathit{C}}_{\mathit{\ell }}\mathrm{+}\mathrm{A}{\mathrm{P}}_{\mathit{\ell }}{\mathrm{A}}^{†}\mathrm{+}{\mathrm{N}}_{\mathit{\ell }}\mathit{,}& & \end{array}$(D.5)where C_ℓ is the angular power spectrum of the CMB, A is a N_chan × d matrix, P_ℓ is a positive definite d × d matrix, and N_ℓ is a diagonal matrix representing the noise power spectra of the data. The parameter vector θ contains all or part of the quantities in (D.5).

The decomposition D.5 reflects the fact that CMB, foregrounds and noise are independent components of the signal. Thus, SMICA is an ICA (independent component analysis) method. It operates by matching the observations $\begin{array}{c}{\mathrm{􏽢}}^{}\\ \mathrm{C}\end{array}$ to the spectral model (D.5) using the criterion (D.4).

The maximal flexibility in a SMICA fit of model (D.5) is obtained with all the parameters free, that is without any constraint on the spectrum C_ℓ, on the diagonal entries of N_ℓ, on a, or on A and P_ℓ. One would ideally fit all those parameters (except for obvious degeneracies, like that between a scale factor in a and the overall normalization of the CMB spectrum C_ℓ) over the whole multipole range. In practice, this turns out to be too difficult, given the large dynamic range both over the sky and over multipoles. We resort to a pragmatic three-step approach in which the criterion (D.4) is minimized by first fitting a, then A, and finally the linear parameters C_ℓ and N_ℓ. Each fit is conducted over the multipole ranges and the sky fraction most appropriate for the parameter of interest, as follows.

We first estimate the CMB spectral law a by fitting all model parameters (that is, without constraint) over a clean fraction of sky (f_sky = 40 %) in the range 100 ≤ ℓ ≤ 680 where the signal is CMB-dominated in most of the channels and the beam window functions are accurately known. In this fit, which is done over a clean part of the sky, we use a foreground emission matrix, A, with only four columns. From this step, we only retain the best fit value for vector a. In the second step, we estimate the foreground emissivity by fixing a to its value from the previous step and fitting all the other parameters over a large fraction of sky (f_sky = 97%) in the range 4 ≤ ℓ ≤ 150 where the signal is dominated by the Galactic emission in all channels. From this second step, we retain the best fit value for the matrix A which, again, is adjusted without constraint other than having d = 6 columns. In the last step, we fit all power spectrum parameters: we fix a and A to their previously found values and fit for C_ℓ and P_ℓ at each ℓ.

Note that the first step (fitting a) amounts to re-calibrating the input maps on the basis of CMB anisotropies. For the maps in thermodynamics units, we find $\begin{array}{ccc}\begin{array}{c}\\ & & \\ & & \end{array}& & \end{array}$(D.6)The value at 857 GHz is not accurately recovered by SMICA, so we have set a₈₅₇ = 1. Since the norm of a is degenerate with a global scale factor for the CMB angular spectrum, it can only be recovered by SMICA up to a scale factor. This degeneracy is fixed here by taking a₁₄₃ = 1. The re-calibration step could have been omitted, since $\hat{{a}}$ is very close the unit vector. However, we found that using $\hat{{a}}$ improved the behavior of SMICA over using a = [ 1,...,1 ].

Before describing implementation details, we explain how SMICA deals with the varying resolution of the input channels, since the discussion thus far assumed that all input maps had the same resolution. Since SMICA works in the harmonic domain, it is a simple matter to account for the beam transfer function, b_i(ℓ), of the ith input map. The CMB sky multipoles s_ℓm contribute s_ℓma_ib_i(ℓ)p_i(ℓ) to the harmonic coefficient ${\mathit{x}}_{\mathit{\ell m}}^{\mathit{i}}$ of the ith map (where p_i(ℓ) is the pixel window function for the HEALPix map at ${\mathit{N}}_{\mathrm{side}}^{\mathit{i}}$). Therefore, in order to produce a final CMB map at 5′ resolution, close to the highest resolution of Planck, we only need to work with input spherical harmonics re-beamed to 5′; that is, to apply SMICA on vectors ${\begin{array}{c}{\mathrm{˜}}_{}\\ {x}\end{array}}_{\mathit{\ell m}}$ with entries $\stackrel{}{{{\mathit{x}}_{\mathit{\ell m}}^{\mathit{i}}}_{\mathrm{˜}}}\mathrm{=}{\mathit{x}}_{\mathit{\ell m}}^{\mathit{i}}{\mathit{b}}_{\mathrm{5}}\mathrm{\left(}\mathit{\ell }\mathrm{\right)}\mathit{/}{\mathit{b}}_{\mathit{i}}\mathrm{\left(}\mathit{\ell }\mathrm{\right)}\mathit{/}{\mathit{p}}_{\mathit{i}}\mathrm{\left(}\mathit{\ell }\mathrm{\right)}$, where b₅(ℓ) is a 5′ Gaussian beam function. By construction, SMICA then produces an CMB map with an effective Gaussian beam of 5′ (without the pixel window function).

We now give further details on the actual implementation of SMICA:

• Inputs: SMICA uses all nine Planck frequency channels from 30 to 857 GHz, harmonically transformed up to ℓ = 4000.

• Pre-processing of point sources: SMICA is applied on input maps in which point sources are subtracted or masked. We start by fitting the PCCS point sources with S/N> 5 to a Gaussian shape where the source amplitude is estimated, together with its position and a constant factor representing the background variance. If the fit is successful (χ² ≤ 2), the fitted point source is removed from the map; otherwise it is masked in all channels and the hole is inpainted by a simple diffusive filling process. This is done at all frequencies except 545 and 857 GHz, where all point sources with S/N> 7.5 are masked and inpainted.

• Beams: when the harmonic coefficients of the input maps are re-beamed at 5′, we do not apply exactly the expression $\stackrel{}{{{\mathit{x}}_{\mathit{\ell m}}^{\mathit{i}}}_{\mathrm{˜}}}\mathrm{=}{\mathit{x}}_{\mathit{\ell m}}^{\mathit{i}}{\mathit{b}}_{\mathrm{5}}\mathrm{\left(}\mathit{\ell }\mathrm{\right)}\mathit{/}{\mathit{b}}_{\mathit{i}}\mathrm{\left(}\mathit{\ell }\mathrm{\right)}\mathit{/}{\mathit{p}}_{\mathit{i}}\mathrm{\left(}\mathit{\ell }\mathrm{\right)}$ mentioned above, because the factor 1 /b_i(ℓ) would diverge at high ℓ for the lowest resolution input channels. That may not be a problem in infinite precision arithmetic, but would lead to matrices $\begin{array}{c}{}^{\mathrm{􏽢}}\\ \mathrm{C}\end{array}\mathrm{\left(}\mathit{\ell }\mathrm{\right)}$ with extremely large condition numbers. Instead, we re-beam with the factor 1 /b_i(ℓ) replaced by min(1 /b_i(ℓ),1000). The re-beaming of the CMB modes then is no longer perfect, but this is of course irrelevant because the thresholding occurs in a regime where the signal is completely dominated by the noise, so that the contribution of the corresponding channel is already highly attenuated by the SMICA weights (as shown in Fig. D.1).

  Fig. D.1 Weights, wℓ, given by SMICA to the input maps, after they are re-beamed to 5′and expressed in KRJ, as a function of multipole. Top panel: linear scale; bottom panel: the absolute value of the weights on a logarithmic scale.

  Fig. D.2 Contribution of each input channel to the noise in the SMICA map.

• Masking: in practice, SMICA operates on a masked sky, the mask being applied after the point source processing. The mask is obtained by thresholding a heavily smoothed version of the point source mask. The threshold is chosen to leave about 97% of the sky. Because of the heavy smoothing, the mask has smooth contours and is only sensitive to large aggregates of point sources: the masked areas mostly lie in the Galactic plane, but include also a few bright regions like the Large Magellanic Cloud.

• Inpainting: the SMICA map used in this paper has no real power in the masked region described above. However, for convenience, an inpainted SMICA map has also been produced by replacing the masked pixels with a constrained Gaussian realization obtained by the method of Benoit-Lévy et al. (2013). That map appears in Planck Collaboration I (2014).

• Binning: in our implementation, we use binned spectra.

• Processing at fine scales: since there is little point trying to model the spectral covariance at high multipoles, because the sample estimate is sufficient, SMICA implements a simple harmonic ILC at ℓ > 1500; that is, it applies the filter (D.2) with ${\mathrm{C}}_{\mathit{\ell }}\mathrm{=}\begin{array}{c}{}^{\mathrm{􏽢}}\\ {\mathrm{C}}_{\mathit{\ell }}\end{array}$.

• Confidence mask: a confidence mask (Fig. 3) is provided with SMICA, constructed in the following way. The SMICA CMB map is bandpass filtered through a spectral window v(ℓ) = exp [−((ℓ − 1700)/200)²/2 ]. The result is squared and smoothed at 2° resolution, yielding a map of the (bandpassed) variance of the CMB map. That variance is corrected for the noise contribution by subtracting the variance map for the noise obtained by the same procedure applied to the SMICA HRHD map. If the SMICA map contained only CMB and noise, the variance map would have a uniform value ∑ _ℓv(ℓ)²b₅(ℓ)²C(ℓ)(2ℓ + 1)/4π = 31.3 μK² over the sky. The confidence map is obtained by thresholding the noise-corrected variance map at 70 μK².

Viewed as a filter, SMICA can be summarized by the weights w_ℓ applied to each input map as a function of multipole. In this sense, SMICA is strictly equivalent to co-adding the input maps after convolution by specific axisymmetric kernels directly related to the corresponding entry of w_ℓ.

The SMICA weights used here are shown in Fig. D.1 (for input maps in units of K_RJ). We see, in particular, the (expected) progressive attenuation of the lowest resolution channels with increasing multipole. Figure D.2 shows the contribution of each input channel to the noise in the SMICA CMB map as a function of multipole. The spectral noise contribution from channel i is simply obtained as w(ℓ)²N_i(ℓ), where w_i(ℓ) is the ith entry of the weight vector w(ℓ) and N_i(ℓ) is the angular spectrum of the ith noise map.

More details about SMICA are given in Cardoso et al. (2008), as well as in applications to the analysis of WMAP (Patanchon et al. 2005) and Archeops data (Tristram et al. 2005). An application to the measurement of the tensor-to-scalar ratio using CMB B-modes is discussed in Betoule et al. (2009). Within the Planck collaboration, SMICA is used to define the Plik high-ℓ likelihood (Planck Collaboration XV 2014), but physical models of foreground emission are used there instead of the non-parametric foreground model used here. SMICA is also used to cross-check the HFI calibration (see Planck Collaboration VIII 2014).

Appendix E: FFP6 simulation results

Here we study the performance of the component separation algorithms on the FFP6 simulation described in Sect. 4, providing additional information beyond that in the body of the paper. Much of the analysis presented here mirrors that shown for the data in Sects. 5 and 6, allowing a direct comparison between the two analyses.

Appendix E.1: CMB maps

First, we show in Fig. E.1 the standard deviation between the four foreground-cleaned FFP6 maps, similar to that shown in Fig. 5 for the data. Figure E.2 shows all pairwise differences between the same maps, mirroring Fig. 6 for the data. These two plots highlight an important point concerning the FFP6 analysis already mentioned in Sect. 5, namely that in near-Galactic regions, where the foregrounds are important, the internal differences between the four algorithms are larger in the FFP6 simulation than in the real data. This is due to the fact that each component separation algorithm has been optimized using the real data in terms of model definition, localization, etc. Then, the same models have been used for the FFP6 simulation without change. Only the parameters within those models are refitted to the new data set. This implies, in fact, that we expect each method to perform better on the data than the simulations in terms of absolute residuals, to the extent that the simulation matches the real sky. In other words, the FFP6 simulation provides a conservative estimate of the residual errors in the real data.

Appendix E.2: Power spectrum residuals from individual components

In Fig. E.3 we show the residual effect of some of the individual components on the foreground-cleaned CMB map. The thermal dust emission, CIB fluctuations, point sources, and noise have been processed individually with each algorithm. All other components (free-free, synchrotron, spinning dust, CO, thermal SZ, and kinetic SZ) are shown as a single, composite residual component.

Appendix E.3: CMB power spectra and cosmological parameters

We assess the performance of our component separation techniques by evaluating cosmological constraints from the foreground-cleaned CMB maps derived from the FFP6 simulation.

Figure E.4 shows the estimates of the angular power spectra of the CMB maps. Figure E.5 compares the cosmological parameters derived from the four foreground-cleaned CMB maps, together with CamSpec⁷ and Plik, to the input (theoretical) parameters for different ℓ-ranges. The parameter space is defined by the same model applied to the real data in Sect. 6, including six ΛCDM and two foreground parameters. All deviations from input parameters are small and within 1σ up to ℓ = 2000, verifying that all methods work well in this multipole range. However, for ℓ_max = 2500 we start to see significant shifts, e.g., for Ω_bh² and n_s. Further, the point source foreground parameter, A_PS, reaches large values, implying that assumptions concerning the high-ℓ foreground model become important. For these reasons, we consider ℓ_max = 2000 as the maximum recommended ℓ range for these maps in the current data release.

Still, the overall agreement is excellent between all codes and all ℓ ranges. In particular, we see that differences in the bandpower spectra at high ℓ between the different codes are mostly absorbed by the two-parameter foreground model. For instance, the Commander-Ruler band power spectrum has more power at high ℓ due to noise or residual point sources, but this excess is well fitted by the two-parameter foreground model, and mostly interpreted in terms of a residual point source component; this is expected, given the lower angular resolution of this map.

As mentioned above, n_s and A_s are to some extent sensitive to ℓ_max. These parameters are degenerate with the foreground parameters. This may suggest that our C_ℓ foreground templates deviate more from the shape of the Poissonian and clustered component in the CMB map. This is a limitation of the simple foreground templates used here. To properly describe the foreground residuals in the reconstructed maps, we should use a foreground power spectrum template tailored to each method. For instance, such templates may be constructed by processing simulated foreground maps though each of the four pipelines. The templates are then given by the pseudo-C_ℓ of each of the processed foreground map. However, our analysis shows that the current simple model provides accurate results when restricting the analysis to ℓ_max = 2000.

Figure E.6 shows the best-fit power spectrum residuals for the CMB map, CamSpec and Plik relative to the input CMB ΛCDM model estimated up to ℓ = 2000. These plots show that the residuals of the CMB map-based best-fit models are comparable to the CamSpec and Plik residuals, and smaller than 40 μK² for most of the ℓ range, with larger deviations observed for CamSpec at ℓ ~ 200. At higher ℓs the residuals are smaller than 10 μK² for both approaches, all showing similar trends. Thus, both the map- and spectrum-based likelihoods recover input parameters reasonably well, with the latter yielding slightly larger deviations from the best-fit model of the input CMB realization.

All Tables

All Figures