© ESO, 2016

1. Introduction

This paper, one of a set associated with the 2015 release of data from the Planck¹ satellite, presents maps of the cosmic microwave background (CMB) anisotropies derived from the full Planck data set comprising 50 months of observations from the Low Frequency Instrument (LFI) and 29 months from the High Frequency Instrument (HFI; Planck Collaboration I 2016). This analysis updates the temperature-only analysis of the first 15.5 months of Planck observations discussed in Planck Collaboration XII (2014), and presents the first CMB polarization maps derived from Planck observations.

Much of the Planck analysis effort since the 2013 data release has revolved around understanding and reducing instrumental systematic uncertainties. As summarized in Planck Collaboration I (2016), this work has been highly successful, reducing the net power from systematic errors in the HFI CMB channels by almost two orders of magnitude on large angular scales (Planck Collaboration XI 2016). The main contributions to these improvements have come from better temperature-to-polarization leakage modelling, reduced analogue-to-digital conversion (ADC) errors, and improved modelling of very long time constants (VLTCs; Planck Collaboration VI 2016; Planck Collaboration VIII 2016). With these improvements, the Planck observations are now sufficiently free of instrumental artefacts to allow a robust determination of the CMB polarization anisotropies on intermediate and small angular scales, covering multipoles ℓ ≳ 20. However, as described both in this paper, and in Planck Collaboration VIII (2016) and Planck Collaboration XI (2016), residual systematics are still not negligible compared to the CMB signal on the very largest scales (ℓ ≲ 20), and these modes are therefore removed from the current maps using a high-pass filter. The rate of progress is still excellent, however, and updated all-scale maps with low large-scale systematics are expected to be released in the near future.

For temperature, the most significant improvement in the Planck 2015 analysis pipeline is calibration based on the orbital CMB dipole rather than the solar dipole (Planck Collaboration I2016). This change, combined with a better understanding of both the Planck beams and transfer functions (Planck Collaboration VI 2016; Planck Collaboration VIII 2016), has reduced the uncertainties in gain calibration, and the agreement between LFI, HFI, and the Wilkinson Microwave Anisotropy Probe (WMAP) has improved to the level of the uncertainties. The calibration of the 70, 100, and 143 GHz channels determined from the CMB temperature power spectrum around the first acoustic peak agrees to within a few tenths of a percent. A similar comparison between 70 GHz and WMAP V-band, and between 100 GHz and WMAP W band shows agreement at the 0.5% level. For further details see Sect. 5.4, Table 4, and Fig. 6 of Planck Collaboration I (2016).

The component separation efforts of the Planck 2015 release are summarized in three papers. The current paper is dedicated to CMB extraction, and presents the main Planck 2015 CMB maps in both temperature and polarization. Planck Collaboration X (2016) addresses astrophysical component separation as implemented by Commander (Eriksen et al. 2004, 2008), a Bayesian component separation algorithm, and presents a global model of the microwave temperature sky ranging from 408 MHz to 857 GHz, including detailed maps of synchrotron, free-free, spinning dust, thermal dust, and CO emission. In addition, a few minor components (line emission around 90 GHz and the thermal Sunyaev-Zeldovich effect near the Coma and Virgo clusters) are included in the model, as are instrumental parameters in the form of calibration, bandpasses, monopoles, and dipoles. The corresponding polarization model includes only synchrotron and thermal dust emission. Planck Collaboration XXV (2016) presents a detailed analysis of the foregrounds below about 100 GHz in both temperature and polarization. For detailed descriptions of the various foreground components relevant for microwave component separation, see either Sect. 2 of Planck Collaboration XII (2014) or Sect. 4 of Planck Collaboration X (2016).

The foreground amplitude relative to CMB polarization is such that effective foreground suppression is required for almost any cosmological analysis, but the optimal approach depends sensitively on the topic in question. For instance, because the fluctuation power of diffuse polarized foregrounds decays as a power-law in multipole moment ℓ (Page et al. 2007; Gold et al. 2011; Planck Collaboration X 2016; Planck Collaboration Int. XXX 2016), it is of greater importance for high-ℓ CMB power spectrum and likelihood estimation to minimize noise sensitivity and to marginalize over unresolved point sources, than it is to model diffuse foregrounds with high accuracy. In this case, it is more convenient to parameterize the residual foregrounds in terms of power spectrum models, and to marginalize over these in terms of a few global parameters, than to marginalize over a large number of per-pixel foreground parameters. The Planck 2015 CMB likelihood therefore employs cross-spectra coupled to simple harmonic space foreground modelling (Planck Collaboration XI 2016), rather than the detailed foreground modelling described in this paper.

Similarly, because of the low – but non-negligible – level of residual instrumental systematics on large angular scales in the Planck 2015 data, the low-ℓ polarization likelihood also implements a special purpose cleaning algorithm, in terms of a simple template fit including only the cleanest 30, 70, and 353 GHz channels (Planck Collaboration VI 2016; Planck Collaboration VIII 2016; Planck Collaboration XI 2016). This approach is similar to that adopted for the 9-yr WMAP likelihood (Bennett et al. 2013), and allows for easy propagation of uncertainties from correlated noise in terms of full pixel-pixel noise covariance matrices.

Other applications, however, such as gravitational lensing and integrated Sachs-Wolfe (ISW) reconstructions (Planck Collaboration XV 2016; Planck Collaboration XXI 2016), constraints on isotropy and statistics (Planck Collaboration XVI 2016), searches for primordial non-Gaussianity (Planck Collaboration XVII 2016), and constraints on global geometry and topological defects (Planck Collaboration XVIII 2016), require actual CMB maps, and these are all based on the products described in this paper.

As in 2013, we apply four complementary CMB component separation algorithms to the Planck 2015 sky maps. In alphabetical order, these are: (1) Commander, a parametric pixel-based Bayesian CMB Gibbs sampler (Eriksen et al. 2004, 2008); (2) NILC, a needlet-based internal linear combination method (Basak & Delabrouille 2012, 2013); (3) SEVEM, which implements linear template fitting based on internal templates in pixel space (Leach et al. 2008; Fernández-Cobos et al. 2012); and (4) SMICA, a semi-blind spectral-matching algorithm fully defined in harmonic space (Cardoso et al. 2008; Planck Collaboration XII 2014). These codes were all applied to CMB temperature reconstruction in Planck Collaboration XII (2014), and have now been extended to polarization, as described in Appendices A–D. In addition, each algorithm is applied to several subsets of the full data set, including half-ring, half-mission, and yearly data splits (Planck Collaboration I 2016). Comparing the resulting maps, both between algorithms and data splits, provides a good understanding of instrumental and algorithmic uncertainties.

The paper is organized as follows. In Sect. 2 we describe the Planck 2015 data selection and pre-processing. In Sect. 3 we briefly review the component separation methods, deferring mathematical details to Appendices A–D. In Sect. 4 we present the derived CMB maps. In Sect. 5 we quantify the residual foreground emission present in the maps by cross-correlation with foreground templates. In Sect. 6 we present angular power spectra and corresponding cosmological parameters. In Sects. 7 and 8 we consider higher-order statistics and gravitational lensing. In Sect. 9 we summarize the main features and limitations of these maps, and provide recommendations on their application. We conclude in Sect. 10.

2. Data selection and pre-processing

In this paper we use the full-mission Planck data (Planck Collaboration VI 2016; Planck Collaboration VIII 2016) and accompanying simulations (Planck Collaboration XII 2016). The CMB temperature maps are derived using all nine frequency channels, from 30 to 857 GHz. The CMB polarization maps are derived from the seven frequency channels sensitive to polarization, from 30 to 353 GHz. In most cases, the component separation methods use frequency channel maps as input, with the exception of the temperature analysis performed by Commander, which uses maps from subsets of detectors in each frequency channel, as specified in Table 1 of Planck Collaboration X (2016). Commander additionally uses the WMAP and the 408 MHz sky maps (Bennett et al. 2013; Haslam et al. 1982) for low-resolution temperature reconstruction, feeding into the low-ℓ temperature likelihood; see Appendix A for further details.

The primary foreground-reduced CMB maps are derived from full-mission maps, maximizing signal-to-noise ratio and minimizing destriping errors. In addition, a number of data splits are analysed to enable internal consistency checks and to make estimates of the properties of the CMB maps. These data splits can also be used in analyses where more than one map is required as input. For the purposes of component separation, each subset of the data must contain the same combination of frequency channels as the full-mission data set. The following data splits have been analysed:

• maps from the first and second half of each pointing period (“half rings”; HR1 and HR2);

• maps from odd and even years, consisting of year 1+3 maps for the LFI channels plus year 1 maps for the HFI channels, and year 2+4 maps for the LFI channels plus year 2 maps for the HFI channels (“years”; YR1 and YR2);

• maps from the first and second half of the mission, consisting of year 1+2 maps for the LFI channels plus half-mission 1 maps for the HFI channels, and year 3+4 maps for the LFI channels plus half-mission 2 for the HFI channels (“half mission”; HM1 and HM2).

The HFI maps in the half-ring data split have spurious correlations between them (Planck Collaboration VIII 2016). If these maps are used to estimate the noise level in the power spectrum, the correlations cause the estimate to be biased low. This was already seen in the analysis of the nominal-mission CMB maps in the 2013 release (Planck Collaboration XII 2014), and the same effect is seen in the full-mission maps. The odd- and even-year data split does not include all of the data, since it omits HFI survey 5; as a consequence, the HFI maps have more missing pixels than in the half-mission data split. For these reasons, the half-mission maps have been used as the primary data split for assessing the properties of the CMB maps and for further analysis.

2.1. Data

The frequency maps are described in detail in Planck Collaboration II (2016) for LFI and Planck Collaboration VII (2016) for HFI. The maps have several features that are relevant to component separation. We summarize them here, and the reader is recommended to consult the references for further details.

• Monopole and dipole contributions from the CMB, CIB, and other astrophysical components are estimated and removed during mapmaking (Planck Collaboration VI 2016; Planck Collaboration VIII 2016). This has an effect on component separation, and each method treats monopoles and dipoles in a different way, as described in Appendices A–D.

• The HFI maps are corrected for zodiacal light emission (ZLE) by subtracting a model of the emission at the ring level during mapmaking (Planck Collaboration VIII 2016). This differs from the treatment in the 2013 component separation (Planck Collaboration XII 2014), where the ZLE model was not subtracted.

• Leakage from intensity to polarization due to bandpass mismatches between detectors is estimated from the measurements of the bandpasses before launch (Planck Collaboration IV 2016; Planck Collaboration VI 2016; Planck Collaboration VIII 2016). These estimates are subtracted from the polarized maps.

• The HFI maps at 100, 143, and 217 GHz are renormalized in order to correct for far-sidelobe effects in the calibration (Planck Collaboration VII 2016).

• Missing pixels are filled with the average values of pixels in the surrounding area. This area is defined as being within a radius of 1° for LFI maps and HFI frequency maps, and within a radius of 1.̊5 for HFI detector subset maps.

Point source catalogues and masks have been provided as input to the component separation process. Construction of the PCCS2 catalogue is described in Planck Collaboration XXVI (2016). Masks have been constructed based on these catalogues. For each frequency channel, the intensity source mask removes all sources detected down to a S/N threshold of 4 for the LFI channels and 5 for the HFI channels. For each polarized frequency channel, the polarization source mask removes all sources detected at 99.99% significance or greater. The details of the masking procedures used by each component separation method are described in Appendices A–D.

2.2. Simulations

To validate our results, we analyse realistic simulations of the Planck data set called the full focal plane 8 (FFP8) simulations. They are based on detailed models of the instrument and sky, and are described in full in Planck Collaboration XII (2016). We summarize their contents here.

2.2.1. CMB

The CMB was simulated using an input Λ cold dark matter (ΛCDM) model based on the 2013 cosmological parameter results (Planck Collaboration XVI 2014). The fiducial simulation contains no primordial tensor modes or primordial non-Gaussianity. However, four variants of the same CMB realization have been produced that include non-zero values of the tensor-to-scalar ratio and non-Gaussianity of a local type.

2.2.2. Foregrounds

The Planck sky model (PSM) has been used to simulate the foreground components. The intensity part of the simulation includes all astrophysical components that were identified in the 2013 release. The diffuse components that are relevant at low frequencies consist of synchrotron, free-free, and spinning dust. At high frequencies, CO, thermal dust, and CIB are included. The foreground modelling has been improved in this version of the simulations. In particular, Planck 353 GHz data were used to improve the frequency scaling of the dust emission. In polarization, the diffuse foregrounds are synchrotron and thermal dust. The extragalactic emission from radio and infrared sources has been simulated in intensity and polarization, and the SZ effect from clusters of galaxies has been included in intensity.

2.2.3. Simulated observations

Time-ordered data (TOD) for each detector were simulated using the satellite pointing, and the individual detector beams, bandpasses, noise properties, and data flags. The same mapmaking used for the data is used to generate maps from the simulated TOD. All of the maps from subsets of the data have also been generated from the simulations.

Two versions of the maps are available, with and without bandpass mismatch leakage. The latter is simulated using the average bandpass for all detectors in a frequency channel, eliminating the leakage effect. The version of the maps without bandpass leakage is considered in this paper.

In addition to the fiducial maps, a set of 10 000 Monte Carlo realizations of CMB and noise has been generated. These realizations are intended to be used to assess the uncertainties on the results.

2.2.4. Mismatch between simulations and data

In analysing the simulations, a number of deficiencies have become evident. First, the amplitude of the CMB component does not match that of the data. This is the (expected) consequence of the fact that the CMB model for the simulations was specified before the recalibration of the Planck data between the 2013 release and the present one was completed. This mismatch can be mitigated by increasing the amplitude of the CMB simulations by 1.3% when comparing them to the data. Second, the noise properties of the simulated maps do not precisely match those of the data. This does not significantly affect the analysis of the CMB temperature maps, since they are signal-dominated. However, it does affect the analysis of CMB polarization maps because they are more noisy. The noise mismatch appears to be scale-dependent, since the adjustment of the amplitude of the noise simulations to match the data depends on the resolution of the maps. This is explored in Sect. 7.

Any analysis that relies on simulations to estimate the uncertainties of a result from the CMB polarization will be limited by these mismatches. Despite this, many analyses are possible, including those using cross-spectrum, cross-correlation, or stacking techniques.

3. Component separation methods

The four methods used by Planck to separate the CMB from diffuse foreground emission were described in detail in Planck Collaboration XII (2014). They are representative of the main approaches to component separation developed in recent years. The methods can be divided into two types. The first type assumes only knowledge of the blackbody spectrum of the CMB, and the foregrounds are removed by combining the multi-frequency data to minimize the variance of the CMB component. The second type constructs an explicit parameterized model of the CMB and foregrounds with an associated likelihood, and the CMB component is obtained by maximizing or sampling from the posterior distribution of the parameters. Either type of method may be implemented in the map domain or in the harmonic domain. We recall briefly their main features and comment on their application to polarization data. Descriptions of the changes in each algorithm since 2013 are given in the appendices.

• Commander (Eriksen et al. 2006, 2008) is a Bayesian parametric method that works in the map domain. Both the CMB and foregrounds are modelled using a physical parameterization in terms of amplitudes and frequency spectra, so the method is well suited to perform astrophysical component separation in addition to CMB extraction (Planck Collaboration X 2016). The joint solution for all components is obtained by sampling from the posterior distribution of the parameters given the likelihood and a set of priors. To produce a high-resolution CMB map, the separation is performed at multiple resolutions with different combinations of input channels. The final CMB map is obtained by combining these solutions in the spherical harmonic domain. This obviates the need for the Ruler step that was used in 2013 to extend the Commander solution to high resolution. A low-resolution version of the separation is used to construct the temperature power spectrum likelihood for large angular scales, as described in Planck Collaboration XI (2016). We note that Commander employs detector and detector set maps rather than full frequency maps, and excludes some specific detector maps judged to have significant systematic errors. Thus, the selection of data is not identical between Commander and the other three methods. For further details regarding data selection and processing, see Appendix A.

• NILC (Delabrouille et al. 2009) is an implementation of internal linear combination (ILC) that works in the needlet (wavelet) domain. The input maps are decomposed into needlets at a number of different angular scales. The ILC solution for the CMB is produced by minimizing the variance at each scale. This has the advantage that the weights used to combine the data can vary with position on the sky and also with angular scale. The solutions are then combined to produce the final CMB map.

• SEVEM (Fernández-Cobos et al. 2012) is an implementation of the template-cleaning approach to component separation that works in the map domain. Foreground templates are typically constructed by differencing pairs of maps from the low- and high-frequency channels. The differencing is done in order to null the CMB contribution to the templates. These templates are then used to clean each CMB-dominated frequency channel by finding a set of coefficients to minimize the variance of the map outside of a mask. Thus SEVEM produces multiple foreground-cleaned frequency channel maps. The final CMB map is produced by combining a number of the cleaned maps in harmonic space.

• SMICA (Cardoso et al. 2008) is a non-parametric method that works in the spherical harmonic domain. Foregrounds are modelled as a small number of templates with arbitrary frequency spectra, arbitrary power spectra, and arbitrary correlation between the components. The solution is obtained by minimizing the mismatch of the model to the auto- and cross-power spectra of the frequency channel maps. From the solution, a set of weights is derived to combine the frequency maps in the spherical harmonic domain to produce the final CMB map. Maps of the total foreground emission in each frequency channel can also be produced. In the analysis performed for the 2013 release (Planck Collaboration XII 2014), SMICA was the method that performed best on the simulated temperature data.

3.1. Extension to polarization

The methods described above were applied to Planck temperature data for the 2013 release (Planck Collaboration XII 2014), and they have been extended to operate on polarization data in the present work. A key distinction between the methods is the choice of operating domain. Two of the methods, Commander and SEVEM, operate in the map domain, so it is most natural for them to perform the polarized component separation on the Q and U maps. The other two methods, NILC and SMICA, operate in the harmonic or needlet domain. An intrinsic part of the transform of polarized maps to these domains is the decomposition of Q and U into E- and B-modes, which is accomplished by using spherical harmonic transforms on the full sky. Thus these two methods perform their separation directly on E and B.

3.2. Outputs

In addition to producing CMB maps, each method provides “confidence” masks to define the region of the sky in which the CMB solution is trusted in temperature and polarization. The procedure each method uses to create the mask is described in Appendices A–D. The confidence masks are used to define masks for further analysis of the data. Two of the pipelines, Commander and SMICA  also produce foreground products, described in Planck Collaboration X (2016).

The first 1000 Monte Carlo realizations of CMB and noise have been propagated through the four pipelines. This has been done twice, once using the parameters derived from the data, and once using the parameters derived from the fiducial FFP8 maps, to provide a set of simulations to accompany each data set.

The methods produce maps at different resolutions, as described in Appendices A−D. The products have been brought to a standard resolution to compare them and for distribution. Standard resolution temperature maps have a Gaussian beam of 5′ FWHM and HEALPix² (Górski et al. 2005) resolution N_side = 2048. Standard resolution polarization maps have a Gaussian beam of 10′ FWHM and HEALPix resolution N_side = 1024. If maps are produced at higher resolution then they are downgraded to these standard resolutions. The downgrading procedure for maps is to decompose them into spherical harmonics (T, or E and B, as appropriate) on the full sky at the input HEALPix resolution. The spherical harmonic coefficients, a_ℓm, are convolved to the new resolution using (1)where b_ℓ is the beam transfer function, p_ℓ is the HEALPix pixel window function, and the “in” and “out” superscripts denote the input and output resolutions. They are then synthesized into a map directly at the output HEALPix resolution. Masks are downgraded in a similar way. The binary mask at the starting resolution is first downgraded like a temperature map. The smooth downgraded mask is then thresholded by setting pixels where the value is less than 0.9 to zero and all others to unity, to make a binary mask. This has the effect of enlarging the mask to account for the smoothing of the signal. In addition to the standard resolution products, the maps, masks, and Monte Carlo realizations have been downgraded to lower resolutions for analyses that need them, using the above procedure.

The polarization maps and Monte Carlo realizations have been decomposed into E- and B-mode maps, and downgraded to lower resolutions too, for analyses that work on E and B directly. The decomposition is done on the full-sky maps using spherical harmonic transforms. The CMB maps from the real-space methods, Commander and SEVEM, are inpainted before doing the decomposition. Both the standard Commander and SEVEM CMB maps are inpainted inside their corresponding confidence masks using a constrained realization (see Appendix A for details). For both methods, only the CMB maps are inpainted, not the Monte Carlo realizations, since this would be too computationally expensive. The other two methods, NILC and SMICA, work on E- and B-modes, so it is possible to make E and B maps directly from their outputs in addition the the standard Q and U maps.

4. CMB maps

In this section we present and discuss the component-separated CMB maps in temperature and polarization. The maps are available to download from the Planck Legacy Archive (PLA)³. For temperature, we compare the maps to those extracted from the nominal mission data in 2013. We also compare the maps from different methods to assess their consistency, and we use the FFP8 simulations to assess the accuracy of the methods and the robustness of the solutions. Throughout the discussion, we make use of appropriate masks in order to highlight differences at intermediate and high Galactic latitudes, ignoring the plane of the Galaxy, where differences are much higher due to the complexity of the foreground signal and its dominance over the CMB.

4.1. Temperature maps

Temperature confidence masks produced by the methods have been used to make combined masks for further analysis of the data. The first mask is constructed as the union of the Commander, SEVEM, and SMICA confidence masks. The NILC mask is not included in the union because it removes a significantly smaller fraction of the sky. This union mask has f_sky = 77.6%. We refer to it as UT78 and adopt it as the preferred mask for analysing the temperature maps. It is shown on the left in Fig. 1. An extended version of the mask has been constructed by adding to the UT78 mask those pixels where the the standard deviation between the four temperature maps is greater than 10 μK. This mask has f_sky = 76.1%, and we refer to it as UTA76. A union mask has been created for the FFP8 simulations in the same way as for the data. It has f_sky = 73.5%, and we refer to it as FFP8-UT74.

The CMB temperature maps produced by the four methods are shown in Fig. 2. No obvious differences are seen in these maps, and they appear visually consistent outside the mask. An important assessment of the robustness and consistency of the CMB T component separation solutions is provided in Fig. 3, which shows the differences between the Planck 2013 and 2015 maps for each method. Several interesting features may be seen in these differences, most of which correspond directly to a better understanding of the systematic uncertainties in the new maps.

Starting with Commander, the most striking features are large-scale swaths tracing the Planck scanning strategy with a peak-to-peak amplitude of approximately 10 μK. This pattern is very similar to that originally pointed out by Larson et al. (2015), who found this by subtracting the 9-yr WMAP ILC map (Bennett et al. 2013) from a template-cleaned version of the Planck 2013 100 GHz map. Similar patterns are also seen in the Commander residual maps shown in Fig. 2 of Planck Collaboration X (2016), corresponding to detector maps that are rejected from the new 2015 analysis. These structures are primarily due to two effects, namely destriping errors from bandpass mismatch between detectors and far sidelobe contamination (Planck Collaboration VIII 2016; Planck Collaboration X 2016). By rejecting particularly susceptible channels in the updated Commander analysis, these errors are greatly reduced in the new map.

Turning to the other three difference maps in Fig. 3, we see that the residuals are internally very similar, but quite different from the Commander residuals. Clear traces of the ZLE are apparent, which is explained by the fact that it was not subtracted in the mapmaking in 2013, but is subtracted in the updated processing (Planck Collaboration VI 2014; Planck Collaboration VIII 2016). Commander, on the other hand, is less sensitive to residual ZLE, for the following two reasons. First, in 2013 it used channels only up to 353 GHz, which are less affected by the ZLE than the higher frequencies. Second, by virtue of fitting independent thermal dust spectral parameters (index and dust temperature) per pixel, it can efficiently absorb the ZLE in the thermal dust component. However, some of the ZLE may still be observed in the Commander differences; remnants of the “red arcs” typically seen in the second and fourth quadrant of the sky are just visible in the difference of the Commander solutions, while being very evident in all other cases. In addition, a quadrupole component aligned with the CMB dipole (Galactic coordinates (l,b) = (264°,48°); Planck Collaboration I 2016) can be seen. This is in part associated with the relativistic Doppler quadrupole. In 2013, the HFI data processing did not subtract this contribution, whereas the frequency-independent component of the quadrupole is removed from the 2015 maps. The difference map thus contains the frequency-independent part of this quadrupole, with small additional contributions due to frequency-dependent residuals arising due to the change of frequency weights between 2013 and 2015. In a reference frame aligned with the CMB dipole, this corresponds to a 4–5 μK amplitude for the a₂₀ component.

The residuals seen in Fig. 3 are small compared to the typical CMB anisotropies, with features mostly constrained to below 5–10 μK, with a distinct large-scale pattern. In particular, these small differences are completely negligible for power spectrum and cosmological parameter estimation. The only cosmological application for which some care is needed is the study of large scale isotropy (Planck Collaboration XXIII 2014). An updated isotropy analysis of the new sky maps is presented in Planck Collaboration XVI (2016); no significant differences are found compared to the 2013 results. However, for studies of the CMB quadrupole−octopole alignment or 2-point correlation function, correcting for the residual Doppler quadrupole may be desirable.

The residual from the frequency-dependent component of the Doppler quadrupole has been estimated using the model described in Planck Collaboration XII (2016). The amplitude of the residual is found by propagating simulated maps of the frequency-dependent part of the kinematic quadrupole through each of the component separation pipelines. The resulting maps give an estimate of the residual effect in the CMB maps after the removal of the frequency-independent component, and may be applied to them as a further correction for the kinematic quadrupole effect. They have been made publicly available through the PLA as part of the ancillary data accompanying the CMB maps.

Table 1 summarizes the main properties of the CMB maps derived both from the data and from the fiducial set of FFP8 simulations. We evaluate standard deviations in two cases, corresponding to high (FWHM 10′, N_side = 1024) and intermediate resolution (FWHM 160′, N_side = 64). The values in parentheses are standard deviations of half-mission half-difference (HMHD) maps, and they give an estimate of the level of uncertainties due to instrumental noise and systematic effects. The same quantities are given for the FFP8 simulations. At this level, results show good consistency for both data and simulations. The SEVEM maps have the lowest standard deviation, as measured by the HMHD maps, at small and intermediate angular scales.

Pairwise differences between all four maps are shown in Fig. 4, after smoothing to 80′ FWHM and downgrading to N_side = 128. As expected, differences are largest close the edge of the mask, where the absolute foreground level is the highest. Comparing these with the corresponding maps from 2013 shown in figure 6 of Planck Collaboration XII (2014), and noting that the new colour bar spans a range that is 4 times narrower than the previous one, we see that the internal agreement between the four methods is substantially better in the new maps, typically by about a factor of 2.

Figure 5 shows the differences between the FFP8 outputs and the input CMB map. The residuals are smallest for NILC and largest for SEVEM. However, we note that the foreground model adopted for the simulations was chosen to be more complex than the real sky, in order to explicitly probe modelling errors. In particular, the simulated thermal dust frequency spectrum exhibits a strong positive (and spatially dependent) curvature at low frequency, which is neither captured in the parametric models adopted by Commander, nor easily modelled by the spatial templates adopted by SEVEM. This additional complexity makes it hard to draw any strong conclusions about the performance on the real data, for which the thermal dust spectrum may be very well approximated by a simple one-component modified blackbody with a nearly spatially constant spectral index (Planck Collaboration X 2016).

4.2. Polarization maps

We now turn our attention to the foreground-reduced CMB polarization maps. As discussed extensively in Planck Collaboration I (2016), Planck Collaboration VI (2016), and Planck Collaboration VIII (2016), the residual systematics in the Planck polarization maps have been dramatically reduced compared to 2013, by as much as two orders of magnitude on large angular scales. Nevertheless, on angular scales greater than 10°, corresponding to ℓ ≲ 20, systematics are still non-negligible compared to the expected cosmological signal. Different combinations of input frequency channels for the component separation have been explored in order to mitigate the polarization residuals. However, it was not possible, for this data release, to fully characterize the large-scale residuals from the data or from simulations. Therefore the CMB polarization maps provided in the current release have been high-pass filtered to remove the large angular scales. This has been implemented by applying a cosine filter to the E and B spherical harmonic coefficients of the maps. This filter is defined as (2)and we have used ℓ₁ = 20 and ℓ₂ = 40. The same filtering has been applied to the FFP8 fiducial maps and to the Monte Carlo simulations.

As we did for temperature, individual polarization confidence masks derived for each method have been used to make combined masks for further analysis. Our first polarization mask is simply the union of the Commander, SEVEM, and SMICA confidence masks, which has f_sky = 77.6% and we refer to it as UP78. However, the 1-point statistics analysis summarized in Sect. 7.1, revealed significant point source contamination in the Commander, NILC, and SMICA maps using this mask. SEVEM was not affected by this problem because it applies an inpainting technique to remove the brightest point sources (see Appendix C for further details). For this reason, two extended versions of the mask were created, the first by excluding in addition the pixels where the standard deviation between the CMB maps, averaged in Q and U, exceeds 4 μK. This mask has f_sky = 76.7%, and we refer to it as UPA77. The second was made by additionally excluding from the union mask the polarized point sources detected in each frequency channel. It has f_sky = 77.4%, and we refer to it as UPB77. This mask is shown on the right side of Fig. 1, and we adopt this as the preferred polarization mask, since it is physically better motivated than UPA77. Also, although it keeps a larger fraction of the sky than UPA77, it is sufficient to alleviate the point source contamination.

Masks have been made in a similar way for the FFP8 simulations. The union mask has f_sky = 76.3%, and we refer to it as FFP8-UP76. An extended mask that also excludes polarized point sources has f_sky = 75.7%, and we refer to it as FFP8-UPA76. FFP8-UPA76 is the preferred mask for FFP8 polarization analysis.

Figures 6 and 7 show the high-pass filtered Q and U Stokes parameter CMB maps after applying the UPB77 mask. The maps are shown at full resolution, and are thus dominated by instrumental noise, except in the regions at the ecliptic poles where integration time is greatest. Visually, the methods operating in the harmonic (NILC and SMICA) and spatial (Commander and SEVEM) domains are more similar to each other than the other methods.

In order to have better visual insight at the map level, in Fig. 8 we show a 20° × 20° patch of the high-pass-filtered Commander polarization maps centred on the north ecliptic pole. The top row shows the full-mission Q and U maps. Note the characteristic “+” pattern in Q and “×” pattern in U; this is the expected signal for a pure E mode signal. To make this point more explicit, the middle row shows the same map decomposed into E and B components. There is a clear asymmetry between them, with E having visibly more coherent power than B, again as expected for an E-dominated signal. Finally, the third row shows the half-ring half-difference (HRHD) E and B maps, illustrating the noise level in the full-mission maps. Comparing the middle and bottom rows, there is clearly an E-mode excess in the full-mission map, whereas the corresponding full-mission B-mode map is consistent with the HRHD B-mode map. In addition, the HRHD E-mode map is also consistent with both the full-mission and HRHD B-mode maps, suggesting that all are consistent with instrumental noise.

Pairwise differences between the four polarization maps are shown in Figs. 9 and 10. The NILC and SMICA solutions appear closest to each other. The regions of the sky that are most affected by differences appear to be those with a higher noise level, as may be seen by comparing to Figs. 6 and 7. For completeness, Figs. 11 and 12 show the differences between the FFP8 outputs and the input CMB map. The differences show a pattern similar to that of the noise, though with higher amplitude with respect to pairwise differences of solutions from data, possibly reflecting again the enhanced complexity of the simulated sky with respect to real data.

The combination of high-pass filtering and noise makes visual comparison of these maps difficult. The rms summary provided in Table 1 is more informative in this respect. Comparing the HMHD rms values listed in parentheses, we see that NILC and SMICA have the lowest effective polarization noise levels at high angular resolution, with rms values that are roughly 20% lower than those observed for Commander and SEVEM. One plausible explanation for this difference is the angular resolution adopted for the fitting process by the four methods; whereas Commander and SEVEM perform the polarization analysis at 10′ FWHM resolution, NILC and SMICA adopt a 5′ FWHM resolution. When comparing the rms values at an angular resolution of 10′, as presented in Table 1, the maps from the latter two methods are smoothed by post-processing to a lower resolution, whereas the maps from the two former codes are not.

This effect is not relevant at intermediate angular scales, for instance at 160′ FWHM, as shown in the bottom section of Table 1. On these angular scales, we see that the situation among the codes is reversed, and Commander provides a 20% lower effective noise than the other three codes.

5. Correlation with external templates

Correlation of the CMB maps with templates of foreground emission provides a first diagnostic of residual contamination in the maps. We compute the correlation coefficient r between two maps x and y as (3)where the index i runs over the N_pix pixels observed with the common mask, ⟨ x ⟩ = ∑ _ix_i/N_pix, σ_x = [ ∑ _i(x_i− ⟨ x ⟩ )²/ (N_pix−1)] ^{1 / 2}, and similarly for y. Both maps and templates are smoothed to FWHM 1° and downgraded to N_side = 256 before computing the correlation, and the analysis is performed separately on temperature and polarization maps.

The foreground templates considered for temperature are: the 408 MHz radio survey of Haslam et al. (1982); the velocity-integrated CO map of Dame et al. (2001); the full-sky Hα template of Finkbeiner (2003); and the Planck 857 GHz channel map. The uncertainty in the value of r due to chance correlations between foregrounds and the cleaned maps is estimated by computing the correlations between the templates and the 1000 simulations of CMB and noise provided by each method. For polarization, the only template we consider is one constructed from the WMAP 9-yr maps. The template is made by smoothing the K and Ka band maps to FWHM 1° and differencing them to remove the CMB contribution. This produces a template containing the low-frequency polarized foreground emission and noise. The correlation analysis is done twice, once with the original maps and templates, and once with a high-pass filtered version. The resulting coefficient factors and 1σ uncertainties are shown in Table 2.

For temperature, the results for all methods are compatible with no correlation within 1σ. From this, we conclude that there are no significant temperature foreground residuals with the same morphology as the templates in the map, to a precision set by cosmic variance. For polarization, the analysis of unfiltered maps shows that SEVEM and SMICA are compatible with no correlation at the 1σ level, Commander has a moderate level of correlation at the 2σ level, and for NILC we find a correlation around 4.5σ. For high-pass filtered maps (labelled by HPF), the level of correlation is reduced for Commander and NILC to about 1.5σ, while it increases for SEVEM and SMICA to about 4 and 2.5σ, respectively. These levels may be due to accidental correlations of the map with the template that are not taken into account when computing the uncertainty. Neither the noise in the template nor the systematics in the maps at large angular scales are modelled, the latter being important for the unfiltered case. In the filtered case, the signal is reduced relative to the noise in the template, which could exacerbate spurious correlations. From this analysis, it is difficult to draw strong conclusions about the residual contamination in the polarized CMB maps.

6. Power spectra and cosmological parameters

In this section we evaluate the foreground-cleaned maps in terms of CMB power spectra and cosmological parameters. We emphasize that the Planck 2015 parameter results are not based on the high-resolution foreground-cleaned CMB maps presented in this paper, but are instead derived from the likelihood described in Planck Collaboration XIII (2016). That likelihood combines the low-resolution Commander temperature map derived from Planck, WMAP, and 408 MHz with a template-cleaned LFI 70 GHz polarization map in a pixel-based low-ℓ likelihood, and adopts a cross-spectrum-based estimator for the high-ℓ temperature and polarization likelihood. The high-ℓ likelihood, called Plik, relies on a careful choice of masks along with templates and modelling, all in the power spectrum domain, to reduce the contribution from diffuse Galactic and discrete Galactic and extragalactic foreground emission. The templates and source models are marginalized over when estimating cosmological parameters.

The parameter estimation methodology that we use here is primarily a tool to evaluate the quality of the high-resolution CMB maps and to assess overall consistency with the Planck 2015 likelihood. We start with the foreground-cleaned CMB maps and masks described in Sect. 4. The maps have been cleaned from diffuse foregrounds and, to a varying extent, from extragalactic foregrounds. We construct simplified templates for the residual extragalactic foregrounds that are marginalized over when estimating cosmological parameters.

While the Planck 2015 likelihood takes into account calibration and beam uncertainties, we have not done so in this analysis. Two of the methods, Commander and SMICA, fit for the relative calibration between frequency channels, but the uncertainties from this process are not propagated into the maps. The other two methods, NILC and SEVEM, assume that the frequency channels are perfectly calibrated. None of the four methods propagate the uncertainties in the beam transfer functions into the CMB maps.

We are interested in assessing the relative quality of the CMB maps, for which it is more important to assess the spread of parameters between methods and as a function of angular scale, rather than to provide absolute numerical values. However, since the CMB maps we describe here are the basis for the analysis of the statistical isotropy of the CMB, primordial non-Gaussianity, and gravitational lensing by large-scale structure, it is both important and reassuring that the parameter values that we find are reasonably close to those of Planck Collaboration XIII (2016).

Figure 13 shows the TT and EE power spectra of the foreground-cleaned maps, masked with an apodized versions of the UT78 mask in temperature and the UP78 mask in polarization. The power spectra of the half-mission half-sum (HMHS) data are shown as thick lines; those of the half-mission half-difference (HMHD) data are shown as thin lines. The HMHS spectra contain signal and noise, whereas the HMHD spectra contain only noise and some residual systematic effects. The variations in the temperature noise spectra at low to intermediate multipoles are caused by the component separation methods optimizing the trade-off between foreground signal and noise as a function of scale. At high multipoles, the same spectra are smooth because the relative weighting of the frequency channels is set by the noise levels. The breaks in the Commander noise spectra are caused by the hybridization of maps at different resolutions to make the final map.

Figure 14 shows the TT and EE power spectra of the CMB signal in the maps estimated using XFaster (Rocha et al. 2010, 2011) from the HMHS and HMHD maps. The first panel shows the four TT power spectra derived from the CMB maps, the Plik spectrum derived from the frequency channel maps, and the best-fit ΛCDM power spectrum derived from the Planck 2015 likelihood using multipoles up to ℓ = 2500. The second panel shows the TT residuals after subtraction of the best-fit theoretical spectrum. The third and fourth panels show the same information for the EE spectra.

The component separation pipelines are principally designed to remove diffuse foreground emission from the CMB, so residuals from unresolved compact sources remain in the maps at small angular scales. For the TT spectrum, these residuals are modelled by using the FFP8 simulations to construct extragalactic source templates as a function of multipole for each method individually, and marginalizing over the corresponding amplitudes during parameter estimation. These templates, however, are only as good as the inputs to the simulations. If the model of the extragalactic components is not correct, this will translate into an error in the shape of the templates. For the EE spectrum, we do not have a detailed model for the residuals, so we assume that it has a shape D_ℓ ∝ ℓ², and marginalize over its amplitude during cosmological parameter estimation. Figure 14 shows the power spectra after subtracting the best-fit model for high−ℓ extragalactic source residuals for TT and EE.

There are differences between the power spectra determined from the CMB maps and the spectrum determined by Plik from the frequency maps. The differences are most noticeable around the first peak in the TT power spectrum, where the CMB map spectra have deficits of power that range from 50–100μK² compared to the Plik spectrum. The spectrum of the SEVEM map is the closest to the Plik spectrum.

Some of the differences can be ascribed to sample variance from the different sky coverage and the different combination of input data. As noted above, the temperature spectra from the CMB maps are determined using the UT78 mask, which has f_sky = 69% after apodization, whereas the Plik spectrum is determined from a combination of frequency-channel spectra using more conservative masks with f_sky = 66% at 100 GHz, f_sky = 57% at 143 GHz, and f_sky = 47% at 217 GHz. A similar situation obtains in polarization, where the CMB map spectra are determined using the UP78 mask with f_sky = 76% after apodization, but the Plik spectrum is obtained from the frequency maps using masks with f_sky = 70% at 100 GHz, f_sky = 50% at 143 GHz, and f_sky = 41% at 217 GHz.

Not all of the observed differences in the spectra can be explained by sample variance. There are differences between the spectra of the maps from the four component separation methods. Sample variance from sky coverage does not contribute to them because the spectra have all been determined using the same mask. One possible source of the differences is the relative calibration of the frequency channels in the component separation. If there is a mismatch in the assumed relative calibration between the channels, or an error in calibration is introduced in the internal recalibration by a method that does it, then the CMB map will not have the expected properties. The relative calibration of each frequency map propagates through the pipelines in an ℓ-dependent way, so, in general, an error in relative calibration produces an error in the effective beam transfer function of the CMB map. There is some evidence to support this hypothesis. During the revision of this paper it was found that the recalibration applied by SMICA to the 44 GHz channel was erroneously large. If the 44 GHz channel is not recalibrated relative to the other channels, then the power spectrum of the resulting SMICA map is much closer to that of the SEVEM map.

Another potential source of differences between the methods is the modelling of the diffuse foregrounds. The models used by Commander and NILC allow the properties of the foregrounds to vary across the sky, whereas SEVEM and SMICA assume they are fixed. It has been assumed that the component separation methods remove the diffuse foregrounds perfectly, leaving only the (compact) extragalactic foregrounds to be modelled in the parameter analysis. If this is not the case, then the level of the residuals will vary between methods, and it must be taken into account in the subsequent analysis.

Cosmological parameters from the component-separated maps in both temperature and polarization are determined using XFaster power spectra, coupled to CosmoMC (Lewis & Bridle 2002) using a correlated Gaussian likelihood. Specifically, we include multipoles between ℓ_min = 50 and ℓ_max, where ℓ_max = 1000, 1500, or 2000 for temperature, and 1000 or 1500 for polarization. We adopt a standard 6-parameter ΛCDM model, and, since low-ℓ data are not used in the likelihood, impose an informative Gaussian prior of τ = 0.070 ± 0.006. As mentioned above, we construct foreground templates in TT for each CMB map by propagating the simulated full-sky FFP8 foreground maps through the respective pipeline and estimating the resulting power spectra, normalized to some pivotal multipole.

The resulting cosmological parameters are summarized in Fig. 15 for both TT and EE. Starting with the temperature cases, we first observe good overall internal agreement between the four component separation methods, with almost all differences smaller than 1σ within each multipole band. Second, we also observe good agreement with the best-fit Planck 2015 ΛCDM model derived from the likelihood, since most of the differences are within 1σ. The notable exception is the power spectrum amplitude, A_se^{− 2τ}, which is systematically low at about the 2σ level for ℓ_max = 1000 for all methods.

In more detail, however, there is some evidence of internal tensions at the 1–2σ level, most clearly seen in Ω_ch² and A_se^{− 2τ}. For these two parameters, there are almost 1σ shifts for NILC, SEVEM, and SMICA going from ℓ_max = 1000 to 1500, and again from ℓ_max = 1500 to 2000. Commander appears to be somewhat more robust with respect to multipole range than the other three methods. Although some variation is indeed expected by statistical variation alone, the combination of the shapes of the power spectrum differences seen in Fig. 14 and the systematic parameter trends suggest systematic uncertainties at the 1−2σ level due to various unmodelled residuals, as discussed above. For this reason, we do not recommend using the component-separated maps for cosmological parameter estimation at this time; for that purpose the Planck likelihood method is preferred, since it shows much better stability with respect to multipole range (Planck Collaboration XI 2016).

For polarization, the results are more ambiguous, with fluctuations relative to the temperature prediction beyond 2σ. Note, however, that the cosmic variance contributions to the temperature and polarization parameters are essentially independent, and the two estimates are therefore only expected to agree statistically, not point-by-point. Nevertheless, SEVEM in particular appears to show evidence of larger deviations than expected in polarization for several parameters, and the θ parameter for Commander in polarization shows some hints of tension with respect to the corresponding temperature estimate.

The differences in the parameters are driven by the differences in the power spectra. In the absence of a model of errors in the relative calibration and the diffuse foreground residuals, they are absorbed in part by the extragalactic foreground model or by the cosmological parameters. However, the remaining differences are at a low level, so overall we find good consistency between the four different component separation methods in both temperature and polarization. In the next data release, we will investigate the effect of the relative calibration of the frequency channels and the modelling of the residual diffuse emission, and we will improve the extragalactic foreground model.

7. Higher-order statistics

We now consider the higher-order statistics of the CMB maps in the form of 1-point statistics (variance, skewness, and kurtosis), N-point correlation functions, and primordial non-Gaussianity (f_NL). We focus in particular on the polarization maps and their degree of consistency with the FFP8 simulations. The temperature results are described in Planck Collaboration XVI (2016), Planck Collaboration XVII (2016), and Planck Collaboration XVIII (2016).

7.1. 1-point statistics

The variance, skewness, and kurtosis of the maps, and the pre-processing steps needed to compute them, are described in detail in Monteserín et al. (2008), Cruz et al. (2011), and Planck Collaboration XVI (2016). The procedure is to normalize the data, d_p, in pixel p by its expected rms, . The rms is calculated from 1000 FFP8 realizations (Planck Collaboration XI 2016), for both CMB anisotropies and instrumental noise. To the extent that represents an accurate description of the data variance, and both the sky signal and the instrumental noise are Gaussian-distributed quantities, will be Gaussian distributed with zero mean and unit variance. In temperature, the variance of the instrumental noise is subtracted in order to determine the variance of the CMB signal in the data. Once the variance of the CMB signal is estimated, it is used to extract the skewness and kurtosis from the normalized map. This procedure is well established, and provides a direct test for the presence of residual foregrounds and of CMB Gaussianity.

The temperature analysis reveals that the Planck 2015 CMB maps and the FFP8 fiducial CMB maps are fully compatible with the Monte Carlo simulations, therefore they can be used for further statistical analyses. For more details about the temperature results and the FFP8 validation analysis please refer to Appendix E.2 and Planck Collaboration XVI (2016).

For polarization, since the Q and U maps are not rotationally invariant, we consider the polarized intensity . P is not Gaussian-distributed with zero mean, and its skewness and kurtosis are non-vanishing; however, by comparing the data with the Monte Carlo ensemble, we can test for residual foregrounds and quantify the performance of the component separation methods. Noise in the P maps is complicated; rather than trying to remove the noise contribution (as is done in the case of temperature), we compare the P maps with the Monte Carlo simulations of CMB and noise together. Table 3 gives the resulting lower-tail probabilities, that is, the percentage of simulations that show a lower variance, skewness, or kurtosis than the P maps. This is done at three different resolutions (N_side = 1024, 256, and 64) for both the data maps (using mask UPB77) and the fiducial FFP8 maps (using mask FFP8-UPA76).

At lower resolutions, Table 3 shows good agreement between the fiducial CMB maps and the Monte Carlo simulations for all methods. At high resolution, skewness, and kurtosis are mostly in agreement with the simulations, with NILC and SMICA showing slightly high values for the kurtosis (lower tail probabilities of 98.5% and 99.4%, respectively). However, there is excess variance in the maps at high resolution, 2−3.5σ away from the mean of the simulations, with Commander deviating the least.

Using the individual components of the FFP8 fiducial maps (CMB, noise, and foregrounds), we are able to quantify the contributions to the statistics of the high-resolution maps from each component separately. Figure 16 compares the variance (left column), skewness (middle column), and kurtosis (right column) extracted from the FFP8 Monte Carlo simulation P maps (histograms) with those extracted from the FFP8 fiducial realization of: the sum of CMB and noise (brown); the sum of CMB, noise, and thermal dust (light teal); the sum of CMB, noise, and radio point sources (violet); and the sum of CMB, noise, and all foregrounds (pink). The last case defines the lower tail probabilities in Table 3. We have also investigated other components, but found that their contributions are negligible.

We see that the only case that is fully compatible with the Monte Carlo ensemble for all methods is that of CMB and noise. In particular, adding the thermal dust component increases the variance slightly outside the acceptable range for both NILC and SMICA  and the same effect becomes even stronger when adding the rest of the foreground components for these methods. Commander is only slightly affected by foregrounds by this measure, and remains within the acceptable range even for the full foreground model, while SEVEM is an intermediate case. For the kurtosis, we find relatively high sensitivity to radio source residuals (the violet and pink lines coincide), but very low sensitivity to diffuse foregrounds. We conclude that the anomalous statistics seen in analysis of the high-resolution component-separated FFP8 CMB maps are due to the foreground components, which could plausibly be caused by the additional complexity of the FFP8 foreground model with respect to the real sky. The Monte Carlo simulations are compatible with the CMB and noise components of the fiducial map, as they should be by design.

Figure 17 shows the variance of the Monte Carlo simulations (from which the lower tail probabilities are summarized in the top half of Table 3) compared with that of the Planck 2015 CMB P maps. The maps analysed are the sum of CMB and noise signal. Since the CMB is practically the same for all methods, the different average values of the variance of the Monte Carlo simulations are a reflection of somewhat higher noise on small angular scales in the SEVEM and Commander maps than in the SMICA and NILC maps. Looking at the data (pink bars), we immediately see that they do not match the simulations. No simulation has a variance as high as the data for any method at any resolution. As already mentioned, this discrepancy is due to an underestimation of the noise in the FFP8 simulations (Planck Collaboration XII 2016). In the variance column of Table 3, in parentheses, we report the excess variance of the data in percent relative to the mean variance of the Monte Carlo simulations. The excess is 3−4% at N_side = 1024, increasing to 10−20% at N_side = 64. Moreover, we see that there are some extreme values for both the skewness and kurtosis, in particular at low resolution for Commander and at high resolution for NILC. However, as long as the variance in the simulations is inaccurate, it is difficult to decide whether residual foregrounds, a true non-Gaussian feature in the map, or noise underestimation is the cause. For now, we simply conclude that the simulations are inconsistent with the data, and this effectively prevents studies of higher-order statistics of this kind for the polarization maps. For temperature, the same is not true because of the much higher signal-to-noise ratio, although care is warranted even then when probing into the noise-dominated regime above ℓ ≈ 1500−2000.

We also performed the 1-point analysis on the Planck 2015 polarization maps using the UP78 mask. In this case, the values of the skewness and kurtosis for SMICA, NILC, and Commander were significantly affected by the presence of point sources. Conversely, the results for SEVEM were consistent with the expected distribution because this method applies an inpainting technique to remove the signal from the brightest point sources (see Appendix C for more details). This motivated the construction of the UPB77 mask, which excluded the brightest point sources detected in polarization, and strongly alleviated this problem for the other three methods, as can be seen in the results in Table 3. The same behaviour was also found for the FFP8 fiducial maps.

7.2. N-point correlation functions

Real-space N-point correlation functions are a useful diagnostic of the statistics of CMB maps, complementary to harmonic analyses. In this section we describe their application to the Planck 2015 CMB polarization maps. Results for the FFP8 CMB maps are given in Appendix E.3. Details of their application to temperature maps may be found in Planck Collaboration XVI (2016).

For observed fields X measured in a fixed relative orientation on the sky, the N-point correlation function is defined as (4)where the unit vectors span an N-point polygon on the sky. Assuming statistical isotropy, N-point functions are functions only of the geometrical configuration of the N-point polygon. In the case of the CMB, the fields X correspond to ΔT and two Stokes parameters Q and U describing the linear polarization of the radiation in direction . In standard CMB conventions, Q and U are defined with respect to the local meridian of the spherical coordinate system of choice. However, Q and U form a spin-2 field and depend on a rotational coordinate system transformation (Leahy et al. 2010). To obtain coordinate-system-independent N-point correlation functions the Stokes parameters are rotated with respect to a local coordinate system defined by the centre of mass of the polygon (see Gjerløw et al. 2010 for details). The Stokes parameters in this new “radial” system are denoted by Q_r and U_r.

Given rotationally invariant quantities X ∈ {ΔT,Q_r,U_r}, the correlation functions are estimated by simple product averages over all sets of N pixels fulfilling the geometric requirements set by θ₁,...,θ_2N−3, which characterize the shape and size of the polygon, (5)The pixel weights are introduced to reduce noise or mask boundary effects. Masks set weights to 1 for included pixels and 0 for excluded pixels.

The shapes of the polygons selected for this analysis are pseudo-collapsed and equilateral configurations for the 3-point function, and a rhombic configuration for the 4-point function, comprising two equilateral triangles sharing a common side. The 4-point function is only computed in the analysis of temperature maps. We use the same definition of pseudo-collapsed as in Eriksen et al. (2005), that is, an isosceles triangle where the length of the baseline falls within the second bin of the separation angles. The length of the longer edges of the triangle, θ, parameterizes its size. Analogously, in the case of the equilateral triangle and rhombus, the size of the polygon is parameterized by the length of the edge, θ. We note that these functions are chosen because of ease of implementation, not because they are better suited for testing Gaussianity than other configurations. In the following, all results refer to the connected 4-point function.

We analyse the Planck 2015 CMB temperature and high-pass filtered polarization maps at resolution FWHM 160′, N_side = 64. We used the downgraded version of the UP78 mask in the analysis.

The N-point functions are used to test the quality of the CMB estimates derived from the Planck data, and are shown in Fig. 18. We show the differences between the N-point functions for the high-pass filtered CMB maps and the corresponding mean values estimated from the 1000 Monte Carlo simulations. The probabilities of obtaining values of the χ² statistic for the Planck fiducial ΛCDM model at least as large as for the CMB map are given in Table 4. The results of the analysis of the temperature maps can be found in Planck Collaboration XVI (2016).

The results for the data deviate significantly from the Monte Carlo simulations for almost all N-point functions involving at least two polarization fields. The smallest deviation is seen for the Commander map. The N-point functions for Monte Carlo simulations have smaller variance than for the data. By comparing HMHD maps with Monte Carlo simulations of noise corresponding to the CMB estimates, we find that the amplitude of the noise is underestimated by 18%. After adjusting the Monte Carlo simulations to compensate, the N-point functions are more consistent with data. However, it is difficult to draw conclusions about the Gaussianity of the polarization maps from these results, since the N-point functions themselves are used to estimate the mismatch between the noise Monte Carlo simulations and the data.

7.3. Primordial non-Gaussianity

Primordial non-Gaussianity is often measured in terms of the amplitude, , of the quadratic corrections to the gravitational potential, as well as by means of the 3-point correlation function based on different triangle configurations. The results from these calculations for the foreground-cleaned CMB maps are presented in Planck Collaboration XVII (2016). Compared to the previous data release, we can now include both temperature and polarization bispectra in the analysis. We thus consider bispectrum f_NL estimates obtained from temperature and polarization data only, as well as the full constraint from all eight possible TTT, TTE, EET, and EEE combinations.

Results obtained from application of the Komatsu-Spergel-Wandelt (KSW) estimator (Komatsu et al. 2005) to the CMB maps after subtraction of the lensing-ISW correlation (see below) are listed in Table 5 for various geometrical configurations that have been considered (Planck Collaboration XVII 2016). It is interesting to evaluate the impact of polarization data on the estimation of f_NL measurements. By considering only the temperature bispectrum for the SMICA map, we obtain , while the polarization alone yields . Uncertainties are evaluated by means of Gaussian FFP8 simulations. We find consistency between pipelines at the 1σ level. The results confirm and extend to polarization the absence of evidence of non-Gaussianity of primordial origin estimated through the bispectrum. The performance of the methods is also tested using Gaussian and non-Gaussian FFP8 simulations, showing that SMICA and SEVEM give the results closest to the inputs in both cases (see Sect. 7.3 of Planck Collaboration XVII 2016).

An interesting case to consider is that of the ISW-lensing bispectrum (Planck Collaboration XXI 2016). Unlike with primordial shapes, there is a specific prediction for the expected amplitude of the ISW-lensing three-point signal in a given cosmological model. We can therefore use this shape to check whether the expected level of NG is recovered, verifying that different component separation methods neither spuriously add nor remove any signal, at least in the squeezed limit where this shape is peaked. This is indeed the case. By normalizing the ISW-lensing shape in such a way as to have an f_NL amplitude of 1 in the best-fit model, we recover f_NL(SMICA) = 0.85 ± 0.2 from the full analysis including all bispectra, f_NL(SMICA) = 0.6 ± 0.3 from temperature alone, and f_NL(SMICA) = 4.7 ± 6.0 from polarization alone.

The constraint from polarization alone is, as expected, much looser than the one from temperature alone. However, polarization does not have a negligible impact on the final combined measurement, mostly due to contributions coming from TTE configurations. Adding polarization reduces the final uncertainty by about 30%, as well as moving the recovered amplitude parameter closer to its expected value of 1. The implications of these results in terms of the physics of the early Universe, as well as the study of many additional shapes, are discussed in Planck Collaboration XVII (2016), which explains the algorithmic details and all results from the procedure summarized here, and in Planck Collaboration XX (2016), which discusses the implications for inflationary physics.

8. Gravitational lensing by large-scale structure

Gravitational lensing from intervening matter imprints a non-Gaussian signature in the CMB temperature and polarization maps, which in turn can be exploited to extract the gravitational potential integrated along the line of sight back to the surface of last scattering. Planck accurately measures the lensing potential over most of the sky (Planck Collaboration XV 2016). As it remaps the CMB polarization, lensing partially transforms primordial E-modes into B-modes, resulting in a secondary B-mode spectrum peaking at around ℓ = 1000. By forming a weighted product of the E-modes from component-separated CMB maps and the reconstructed lensing potential, it is possible to generate a map of the expected lensing-induced B-modes of CMB polarization. Cross-correlating this lensing B-mode template with the total observed B-modes provides an indirect measurement of the lensing B-mode power spectrum (Planck Collaboration PIP116, in prep.; Planck Collaboration XV 2016). Lensing B-mode template maps were synthesized for the four component separation methods considered in this paper, using the Stokes parameter maps and the lensing potential reconstruction from the intensity map, as discussed in Planck Collaboration PIP116 (in prep.). A common mask was generated by combining the union polarization mask and the lensing potential 80% mask of Planck Collaboration PIP116 (in prep.). After apodization using a cosine function over 3°, the resulting mask preserves an effective sky fraction of about 60%.

The lensing B-mode power spectra for Commander, NILC, SEVEM, and SMICA, which have been measured by cross-correlating with the corresponding foreground cleaned polarization maps, are shown in Fig. 19. The four CMB polarization solutions lead to consistent lensing B-mode power spectrum measurements within 1σ over the entire probed multipole range. In addition, the χ² relative to the 2015 PlanckΛCDM base model (Planck Collaboration XIII 2016) of the lensing B-mode band powers in eight multipole bins are 12.2, 9.0, 8.0, and 7.0, with corresponding probability-to-exceed values of 14, 34, 44, and 54% for Commander, NILC, SEVEM, and SMICA, respectively, which indicate the agreement of the Planck lensing B-mode signal with theoretical expectations. We measure amplitude fits with respect to the Planck base model of

corresponding to 10, 11, 10, and 12σ detections of lensing B-modes for Commander, NILC, SEVEM, and SMICA, respectively.

9. Summary and recommendations

We now summarize the results, and provide a critical analysis of the applicability of the derived maps for cosmological purposes.

Starting with the temperature case, we have shown that the four CMB maps are in excellent agreement overall. The amplitudes of the pairwise difference maps are smaller than 5 μK over most of the sky on large angular scales, and the high-ℓ power spectra agree to around 1σ. Correspondingly, the differences with respect to the Planck 2013 maps are typically smaller than 10 μK over most of the sky, and the morphology of the differences is well understood in terms of improved treatment of systematic errors in the 2015 analysis. We conclude that the Planck 2015 temperature maps provide a more accurate picture of the CMB sky than the 2013 temperature maps, having both lower noise and lower levels of systematics, and we expect these maps to find the same cosmological applications as the previous generation of maps. However, we emphasize that these maps are not cleaned of high-ℓ foregrounds, such as extragalactic point source or SZ emission. Residuals lead to biases in cosmological parameters at the 2σ level, beyond ℓ ≃ 2000. Cosmological analyses using small angular scales must therefore take care to marginalize over such foregrounds, as appropriate.

For polarization, the situation is significantly more complicated, due to two different problems. First, a low level of residual systematics on large angular scales prevents a faithful CMB polarization reconstruction for multipoles ℓ ≲ 20. These modes are therefore removed by high-pass filtering in the current maps. Any cosmological analysis of these maps must take into account the corresponding transfer function in order to avoid biases. Second, due to the current noise mismatch between the FFP8 simulations and the data, we strongly caution against using the polarization maps provided here for any cosmological analysis that depends sensitively on the assumed noise level. Nevertheless, the maps should prove useful for a number of other applications that do not require detailed noise simulations, for instance estimation of cross-spectra, cross-correlations, or stacking analyses, and we therefore release the maps to the public despite these limitations.

Considering the four component separation methods in greater detail, the results may be distinguished according to two criteria, namely data selection and basis functions. While Commander performs data selection at the detector (or detector-set) map level, rejecting potentially problematic maps, the other three methods employ frequency channel maps, and thereby maximize the signal-to-noise ratio. Likewise, while Commander and SEVEM perform their analyses in pixel space, NILC and SMICA perform all operations in harmonic space. These distinctions can explain many of the qualitative differences discussed in the previous sections.

We make the following recommendations regarding the use of the four maps. First and foremost, we strongly recommend that any cosmological analysis based on these maps consider all four maps in parallel, in order to assess the impact of specific choices of implementation and modelling. To be considered robust, no results should depend strongly on the specifics of a given component-separation algorithm. Considering specific details, we generally consider Commander to be the preferred solution on large and intermediate angular scales, due to its somewhat lower large-scale effective polarization noise (Table 1), weaker cross-correlation with the high-pass filtered WMAP K−Ka band synchrotron template (Table 2), lower N-point correlation function fluctuations (Fig. 18 and Table 4), and weaker cosmological parameter dependence on ℓ_max, suggesting less internal tension between low, intermediate, and high multipoles (Fig. 15). In addition, the method is able to propagate uncertainties from the input maps to final products by means of Monte Carlo samples drawn from the full posterior. For these reasons, we adopt the Commander solution for the low-ℓPlanck 2015 temperature likelihood (Planck Collaboration XI 2016).

However, at high multipoles the Commander solution exhibits a significantly higher effective point source amplitude than the other three maps, due to the exclusion of frequencies below 217 GHz. For temperature, the lowest residual high-ℓ foregrounds are instead obtained by SMICA, as shown in Fig. 14. As a result, as in 2013, we confirm our preference for the SMICA map for analyses that require full-resolution observations in temperature, such as f_NL (first three rows in Table 5 and Planck Collaboration XVII 2016) or lensing reconstruction (Fig. 19 and Planck Collaboration XV 2016). SEVEM is also a very good choice for temperature, since it provides the map with the lowest level of noise at a wide range of scales, as well as a smooth noise power spectrum as measured by the HMHD maps (see Table 1 and Fig. 14). It also performs equally well as SMICA with regard to the estimation of f_NL.

In polarization, NILC and SMICA perform equivalently at high multipoles (Fig. 14). The NILC polarization maps yield measurements of f_NL that are most consistent with zero (Table 5). The NILC and SMICA analyses also provide an effective mapping of the weights of the Planck frequencies in the needlet/harmonic domains: a given weight tends to be high if a given frequency channel, in a given band in the harmonic domain, is relevant for foreground cleaning; on the other hand, the higher the statistical noise at a given frequency, the lower the associated weight. For NILC and SMICA, the weights for T-, E-, and B-modes are shown in Figs. B.2 and D.1, respectively. SEVEM provides the most stability with respect to the effect of bright point sources in polarization (see Sect. 7.1). This is due to the inpainting procedure applied to these sources, which significantly reduces the effect of this contaminant. Therefore, SEVEM could be a suitable choice for those analyses in polarization that cannot easily deal with the presence of point source holes in a mask.

Finally, we note that the SEVEM approach is unique in its ability to provide independent CMB estimates in a number of frequency channels. For analyses that benefit significantly from, or even require, such information, SEVEM is the only meaningful choice. Specific examples include various isotropy estimators (Planck Collaboration XVI 2016), the integrated Sachs-Wolfe stacking analysis (Planck Collaboration XXI 2016), Doppler boosting (Planck Collaboration XXVII 2014), and Rayleigh scattering analyses (Lewis 2013).

10. Conclusions

We have presented four different foreground-reduced CMB maps in both temperature and polarization, derived from the Planck 2015 observations. These maps are based on the full Planck data, including a total of 50 months of LFI observations and 29 months of HFI observations. The temperature component of these maps represents the most accurate description of the CMB intensity sky published to date. In the polarization component, the characteristic E-mode signal expected in a standard ΛCDM model is easily discernible over the full sky. Corresponding astrophysical foreground products are described in Planck Collaboration X (2016).

The CMB maps presented here are the direct result of the detailed analyses of systematic errors described in Planck Collaboration VI (2016) and Planck Collaboration VIII (2016), which led to an effective reduction of systematic errors by almost two orders of magnitude in power on large angular scales in polarization compared to 2013. However, despite these improvements, the polarization systematic errors in the Planck 2015 data set are not yet negligible in several frequencies on the very largest scales. Multipoles below ℓ = 20 are therefore suppressed by low-pass filtering in the current component-separated CMB maps.

Additionally, as already noted in Planck Collaboration XII (2016), we observe a mismatch in the effective noise amplitude of 10–20% when comparing the latest generation of Planck simulations (FFP8) with the data. Considering that the current temperature sky maps are signal-dominated up to ℓ ≈ 2000, this noise mismatch is of little practical importance for cosmological analyses based on temperature observations, except on the very smallest scales. As in 2013, the temperature maps presented here are therefore used for a wide range of important applications, including large-scale temperature likelihood estimation (Planck Collaboration XI 2016), gravitational lensing (Planck Collaboration XV 2016), studies of isotropy and statistics (Planck Collaboration XVI 2016), primordial non-Gaussianity (Planck Collaboration XVII 2016), and non-trivial cosmological topologies (Planck Collaboration XVIII 2016). For high-ℓ power spectrum and likelihood estimation, we recommend the cross-spectrum based methods described in Planck Collaboration XI (2016), primarily due to difficulties in establishing sufficiently accurate models of unresolved extragalactic high-ℓ foregrounds for the maps presented here. Cosmological parameters derived by temperature power spectra from these maps have been compared with the results of the Planck 2015 likelihood analysis (Planck Collaboration XI 2016) and found to agree at the 1σ level.

For polarization, the noise mismatch is not negligible, and we therefore do not yet recommend using the provided maps for cosmological studies that require a highly accurate noise model. The polarization maps presented here may still be very useful for many important cosmological applications, including cross-correlation and cross-spectrum based analyses, and we therefore release the maps despite the current noise mismatch. Analysis of the higher-order statistics of these maps has been performed within the current framework of precision assessment and presented in this paper. The bispectrum analysis, including E-mode polarization, applied to the results of the four component separation methods, gives evidence of a vanishing non-Gaussian signal for three geometrical configurations, namely local, orthogonal, and equilateral. The B-modes derived from the current polarization maps have been cross-correlated with the predicted lensing B-modes from the measured E signal and the lensing potential measured independently from temperature in Planck Collaboration XV (2016); the result is found to be in excellent agreement among the four component separation methods, as well as with the prediction of the Planck best-fit cosmology.

On the basis of these encouraging results, intense work is continuing to reduce the large-scale polarization systematics to negligible levels, as well as to resolve the noise simulation mismatch, and good progress is being made. Updated products will be published as soon as this work has reached a successful completion.

---

Acknowledgments

The Planck Collaboration acknowledges the support of: ESA; CNES, and CNRS/INSU-IN2P3-INP (France); ASI, CNR, and INAF (Italy); NASA and DoE (USA); STFC and UKSA (UK); CSIC, MINECO, 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); ERC and PRACE (EU). A description of the Planck Collaboration and a list of its members, indicating which technical or scientific activities they have been involved in, can be found at http://www.cosmos.esa.int/web/planck/planck-collaboration. Some of the results in this paper have been derived using the HEALPix package.

References

Appendix A: Bayesian parametric fitting

Commander (Eriksen et al. 2004, 2008) fits a physical model to a set of observations within a standard Bayesian parametric framework, defined by a set of explicit physical parameters and priors. The code can be operated in two modes, either employing Gibbs sampling to map out the full parameter posterior, or using iterative non-linear searches to derive the maximum-likelihood solution; the only implementational difference between the two is whether to sample from or maximize the conditional posterior in each Gibbs step. All maps presented in this paper are derived in the maximum-likelihood mode, and uncertainties are evaluated through simulations.

Commander forms the core of the Planck 2015 foreground-targeted diffuse component separation efforts, and the corresponding results are described in full detail in Planck Collaboration X (2016). In this section we only summarize the most relevant steps for CMB-oriented analysis.

Appendix A.1: Intensity

In 2013, the Commander CMB temperature solution was derived using only the seven lowest Planck frequency maps between 30 and 353 GHz, adopting a simple four-component signal model, including CMB, CO, modified blackbody thermal dust, and a single power-law low-frequency component (Planck Collaboration XII 2014). The only instrumental parameters included in the analysis were monopole and dipole corrections. In the current release, significantly more data are included in the analysis, and the astrophysical and instrument models have been expanded to account for more effects. Specifically, a total of 32 maps are considered in the analysis, including 21 Planck detector and detector-set maps, 10 WMAP differencing assembly maps, and a 408 MHz low-frequency survey map (Haslam et al. 1982). This wide frequency range allows us to fit separately for synchrotron, free-free, spinning dust, thermal dust, and CMB, as well as individual CO transitions at 115, 230, and 345 GHz, a common line emission component in the 94 and 100 GHz channels (primarily HCN), and thermal Sunyaev-Zeldovich emission for the Coma and Virgo clusters. On the instrumental side, we now fit for both calibration and bandpass uncertainties, in addition to monopoles and dipoles (Planck Collaboration X 2016).

However, the current Commander implementation requires uniform angular resolution across frequencies in order to estimate spectral parameters correctly. This implies that all frequency maps must be smoothed to the resolution of the lowest resolution channel before analysis. In the 2013 Planck release, this problem was partially solved by first determining spectral parameters at low resolution using Commander, as described above, and then solving for the component amplitudes from full resolution data using a so-called Ruler step, leading to the Commander-Ruler hybrid. In the current release, we adopt a Commander-only multi-stage approach, in which the system is solved at four different angular resolutions, using different subsets of the data, but each with internally coherent angular resolutions. At each higher-resolution step, all frequency channels with lower angular resolution are dropped from the analysis, and the model is simplified as required by the new frequency subset. We consider the following four cases in our analysis:

• 1.

  Planck+WMAP+408 MHz at1°FWHM resolution, with a full foreground and instrumental model as described above. A corresponding mask is defined through a χ² evaluation, as described in Planck Collaboration X (2016), and allows 93% of the sky for analysis. This solution represents the most complete astrophysical model of the microwave frequencies on large angular scales available as of today, and forms the basis of the Planck 2015 temperature likelihood at large angular scales (Planck Collaboration XI 2016).

• 2.

  Planck 30–857 GHz at 40′FWHM resolution, with a single power-law foreground model at low frequencies to account simultaneously for synchrotron, free-free and spinning dust emission. All global instrumental parameters are fixed to the values obtained in the first step. As in the first step, the mask is defined in terms of the per-pixel χ², but with a higher (i.e., less restrictive) threshold, resulting in a total sky fraction of 98.4%; see Fig. A.1. The low-frequency foreground sector of this model is identical to that adopted for the Planck 2013 component separation analysis.

• 3.

  Planck 143–857 GHz at FWHM resolution, including only CMB, CO, and thermal dust emission. The dust emissivity is fitted in each pixel, while the dust temperature is fixed to the low-resolution solution. The mask is defined as the product of the per-pixel χ² at the native resolution, the 40′ FWHM mask, and an external point source mask, and it admits a sky fraction of 95%.

• 4.

  Planck 217–857 GHz at ~ 5′FWHM resolution, fitting only the CMB and thermal dust emission amplitude per pixel. All other parameters are fixed to their lower-resolution values. The highest-resolution mask is defined in the same way as the FWHM one, but with an additional cut on the CO amplitude of 0.5 K km s^-1 (see Planck Collaboration X 2016), since CO is no longer part of the foreground model. This mask admits a sky fraction of 82.4%.

A comparison of the various Planck-only masks are shown in Fig. A.1, while the likelihood mask is shown in Planck Collaboration X (2016).

Given these single-resolution maps, we hybridize the three Planck-only solutions (ranging between 40 and 5′ FWHM) into a single map as follows. First, we replace the masked regions at each resolution with a constrained Gaussian realization (see Eriksen et al. 2004 for details), drawn from P(a_CMB | C_ℓ,d), to suppress ringing and Galactic leakage effects. Then, we combine the three Planck-only maps in harmonic space using the following cosine apodization scheme, (A.1)where is the beam transfer function and is the pixel window function for a particular resolution, and the cosine apodization weights are given by as illustrated in Fig. A.2. Figure A.3 shows the contributions from the three different resolutions and their sum.

The result is a single full-resolution CMB map with variable effective sky fraction as a function of multipole. For multipoles below ℓ = 200, a total of 98.4% of the sky are formally derived from real measurements, while the remaining 1.6% are filled with a Gaussian constrained realization. For multipoles between ℓ = 300 and 1000, 95% of the sky are derived from real data, while at multipoles above ℓ = 1200, only 82% of the sky are derived from real data. However, the fact that the map is data-driven in a particular pixel for a given angular scale does not necessarily guarantee proper goodness-of-fit, and we therefore recommend application of the most conservative 82% mask for high-precision science. For analyses requiring maximum sky coverage, we instead recommend the likelihood map, which supports 93% of the sky, albeit at lower angular resolution. Figure A.4 shows the difference between the Planck-only and the Planck+WMAP+408 MHz Commander CMB solutions at 1° FWHM resolution. The rms difference between these two maps is 2.7  μK outside the common UT78 mask.

Appendix A.2: Polarization

The Commander CMB polarization map is derived in an analogous manner to the temperature solution, but relying on Planck observations alone. At low resolution, we derive a 40′ map from all frequencies between 30 and 353 GHz, including CMB, synchrotron and thermal dust in the signal model. At high resolution, we derive a 10′ map from frequencies between 100 and 353 GHz, including only CMB and thermal dust. The two maps are hybridized into one 10′ map using the same formalism as in Eq. (A.1), but with only one hybridization operator covering ℓ = 200−300.

As discussed in Planck Collaboration X (2016), several spectral index models have been explored for polarized synchrotron and thermal dust emission, including: (1) spatially constant; (2) smoothed over large angular scales; (3) based on the temperature model. The main conclusion, however, is that the current polarization data have too low signal-to-noise ratio for both synchrotron and thermal dust indices to discriminate between the three models at a statistically significant level. Allowing such additional degrees of freedom only increases the degeneracies between the various components without improving the overall fit, and also leaves the solution sensitive to large-scale residual systematics. For now, we therefore adopt the temperature-derived spectral parameters for both synchrotron and thermal dust for the primary Commander CMB polarization map.

For polarization masking purposes, we construct a mask from the product of the thresholded low-resolution χ² map and the CO amplitude map. The same mask is applied at both the 40 and 10′ FWHM resolutions.

Appendix B: Internal linear combination in needlet space

The goal of NILC is to estimate the CMB from multi-frequency observations while minimizing the contamination from Galactic and extragalactic foregrounds, as well as instrumental noise. The method makes a linear combination of the data from the input maps with minimum variance using a basis of spherical wavelets called needlets (Narcowich et al. 2006). Due to their unique properties, needlets allow localized filtering in both pixel space and harmonic space. Localization in pixel space allows the weights of the linear combination to adapt to local conditions of foreground contamination and noise, whereas localization in harmonic space allows the method to favour foreground rejection on large scales and noise rejection on small scales. Needlets permit the weights to vary smoothly on large scales and rapidly on small scales, which is not possible by cutting the sky in zones prior to processing.

The NILC pipeline (Basak & Delabrouille 2012, 2013) is applicable to scalar fields on the sphere, hence we work separately on maps of temperature and the E- and B-modes of polarization. The decomposition of input polarization maps into E and B is done on the full sky. At the end, the CMB Q and U maps are reconstructed from the E and B maps.

Prior to applying NILC, all of the input maps are convolved or deconvolved in harmonic space to a common resolution corresponding to a Gaussian beam of 5′ FWHM. Each map is then decomposed into a set of needlet coefficients. For each scale j, needlet coefficients of a given map are stored in the form of a single HEALPix map. The filters used to compute filtered maps are shaped as follows For each scale j, the filter has compact support between the multipoles and with a peak at (see Table B.1 and Fig. B.1). The needlet coefficients are computed from these filtered maps on HEALPix pixels with N_side equal to the smallest power of 2 larger than .

In order to show the contribution of the various frequency channels to the final CMB map at different needlet bands, we compute the full sky average of needlet weights for each frequency channel and needlet band. Figure B.2 shows that most of the contribution to the reconstructed CMB maps comes from the 143 GHz and 217 GHz channels. In the low-ℓ needlet bands, the contribution from 143 GHz is large compared to that from 217 GHz. However, due to better angular resolution, the 217 GHz channel contributes more than the 143 GHz channel in the highest ℓ needlet bands. In the intermediate needlet bands, the contributions from these two channels are comparable. It is interesting to note the contribution from the LFI in the lowest ℓ bands; in intensity, the 70 GHz channel serves mostly for foreground removal, while for polarization it contributes to the CMB solution with a positive weight, and very similarly between E and B. We stress again that the lowest ℓ modes have been filtered out in the results presented here, and therefore the low ℓ results will require further investigation.

Appendix B.1: Masking

The confidence masks for NILC for intensity and polarization have been generated following a procedure similar to that used by SMICA, but adopting a different parameterization.

For intensity, the NILC CMB map is filtered through a spectral window (B.1)The result is then squared and smoothed with a Gaussian circular beam with FWHM 120′. The variance map obtained in this way is then corrected for the noise contribution by subtracting the variance map for the noise obtained in the same way from the NILC HRHD map. The confidence map is obtained by thresholding the noise-corrected variance map at 73.5 μK². For polarization, the polarized intensity is obtained from the NILC outputs, and is filtered through a circularly symmetric Gaussian window function of FWHM 30′. The result is then squared and smoothed with a Gaussian of FWHM 210′. The resulting variance map is corrected for the noise contribution by following the same procedure used for intensity. The confidence map is obtained by thresholding at 6.75 μK². The resulting masks are shown in Fig. B.3.

Appendix C: Template fitting

The SEVEM method (Leach et al. 2008; Fernández-Cobos et al. 2012) aims to produce clean CMB maps for several frequency channels by using internal template fitting. The templates are constructed from the Planck data, typically as the subtraction of two close Planck frequency channels to remove the CMB signal. The beam sizes of the maps are equalized before subtraction. The cleaning is achieved simply by subtracting a linear combination of the templates t_j(x) from the data, with coefficients α_j obtained by minimizing the variance outside a given mask, (C.1)Here n_t is the number of templates used, and T_c(x,ν) and d(x,ν) correspond to the cleaned and raw maps at frequency ν, respectively. The same expression applies for T, Q, or U.

The cleaned frequency maps are then combined in harmonic space, taking into account the noise level, resolution, and, for temperature, a rough estimate of the foreground residuals of each cleaned channel, to produce a final CMB map at the required resolution.

Appendix C.1: Implementation for temperature

For temperature, we followed a similar procedure to that used for the Planck 2013 release, i.e., the 100, 143, and 217 GHz maps are cleaned using four templates constructed from the six remaining frequency channels. A few differences have been implemented in the current pipeline with respect to the previous work: the use of a single coefficient over the whole sky for each template (instead of defining two regions), the use of inpainting to reduce contamination from sources, and the use of the 857 GHz channel as a template (instead of 857−545). We note that the other three templates (30−44, 44−70, and 353−217) are the same as for the previous release.

The six frequency channels used to construct templates (30 to 70 GHz and 353 to 857 GHz) are inpainted at the position of sources detected by the Mexican hat wavelet (MHW) algorithm (Planck Collaboration XXVI 2016). The size of the holes to be inpainted is determined by taking into account the beam size of the channel as well as the flux density of each source. We do a simple diffusive inpainting, which fills one pixel with the mean value of the neighbouring pixels in an iterative way. To avoid inconsistencies when subtracting two channels, each map is inpainted on the sources detected in both that map and on the other map used to construct the template. For example, to construct the (30−44) template, both maps are inpainted in the positions of the sources detected at 30 and 44 GHz. This reduces significantly the contamination from compact sources in the templates.

Once they have been inpainted, the maps are smoothed to a common resolution and then subtracted. To construct the first three templates, the first channel in the subtraction is smoothed with the beam of the second map and vice versa. For the 857 GHz template, we simply filter the map with the 545 GHz beam (this is done for comparison with the 857−545 template from the 2013 pipeline, where the 857 GHz map was smoothed by the 545 GHz beam). We note that the coefficients used to multiply this template are typically of order 10^-4, so the level of CMB signal introduced by this template in the final cleaned map is negligible. We take advantage of this to drop subtraction of the 545 GHz map, as was done in 2013, nominally to remove the CMB signal. This simplifies the method and also reduces the noise in this template.

The 100, 143, and 217 GHz maps are then cleaned by subtracting a linear combination of the four templates. The coefficients of the linear combination (Table C.1) are obtained by minimizing the variance outside an analysis mask. The main difference with respect to the 2013 release is that we have used the same coefficients for the whole sky (instead of dividing it into two regions), since this simplifies the procedure without affecting the quality of the reconstruction (other than on those pixels very close to the Galactic centre that need to be masked in any case). This analysis mask covers the 1% of the sky with the brightest emission, as well as sources detected at all frequency channels. Once the maps are cleaned, each is inpainted on the source positions detected at that (raw) channel. Then, the MHW algorithm is run again, now on the cleaned maps. A relatively small number of new sources are found and are also inpainted at each channel. The resolution of the cleaned map is the same as that of the raw map. We note that no assumptions about the noise or foregrounds are made in order to construct the single-frequency cleaned maps.

Finally, the SEVEM CMB map is constructed by combining the cleaned and inpainted 143 and 217 GHz maps. In the combination, the maps are weighted taking into account the noise, resolution, and a rough estimate (obtained from realistic simulations) of the foreground residuals in each map. The weights are shown in Fig. C.1. The resolution of this map corresponds to a Gaussian beam of FWHM 5′ and HEALPix resolution N_side = 2048, with a maximum multipole ℓ_max = 4000.

The same procedure, including the full inpainting process for point sources, is applied when running the pipeline on FFP8 simulations.

Appendix C.2: Implementation for polarization

To clean the polarization maps, a procedure similar to the one used for the temperature maps is applied to the Q and U maps independently. However, given that narrower frequency coverage is available for polarization, the templates and maps to be cleaned are different. In particular, we clean the 70, 100, and 143 GHz maps using three templates. Since the signal-to-noise ratio is lower, at 100 and 143 GHz, the clean maps are produced at N_side = 1024, with resolution corresponding to a Gaussian beam of FWHM 10′ and a maximum multipole ℓ_max = 3071. At 70 GHz, the map is produced at its native resolution.

The first step of the pipeline is to inpaint the positions of the sources detected using the MHW algorithm in those channels that will be used to construct templates. Similarly to the temperature case, for a given template (constructed as the difference of two frequency channels), the inpainting is performed for all of the sources detected in any of the channels. We note that the inpainting is performed in the frequency maps at their native resolution.

These inpainted maps are then used to construct a total of four templates. To trace the synchrotron emission, we construct a 30−44 template. For the dust emission, the following templates are used: 353−217 (smoothed to 10′ resolution); 217−143 (used to clean 70 and 100 GHz); and 217−100 (to clean 143 GHz). These last two templates are constructed at 1° resolution, since an additional smoothing becomes necessary in order to increase the signal-to-noise ratio of the template. Otherwise, the estimated coefficients are driven by the noise and the cleaning is less efficient. Since fewer frequency channels are available in polarization, it becomes necessary to use the maps that are to be cleaned as part of one of the templates. Therefore the 100 GHz map is used to clean the 143 GHz frequency channel and vice versa, making the clean maps less independent than in the temperature case.

These templates are then used to clean the non-inpainted 70 (at its native resolution), 100 (at 10′ resolution), and 143 GHz maps (also at 10′). The corresponding linear coefficients (listed in Table C.2) are estimated independently for Q and U by minimizing the variance of the clean maps outside a mask that covers compact sources and the 3% of sky with the brightest Galactic emission. Once the maps have been cleaned, inpainting of the sources detected at each map is carried out. The size of the holes to be inpainted takes into account the additional smoothing of the 100 and 143 GHz maps. Similarly to temperature, the same inpainting processing is applied to the point source positions when running the pipeline on the FFP8 simulations. The 100 and 143 GHz clean maps are then combined in harmonic space, using a full-sky E and B decomposition, to produce the final CMB maps for the Q and U components with a Gaussian beam of FWHM 10′ and HEALPix resolution N_side = 1024. Each map is weighted taking into account its noise level at each multipole.

Before applying the post-processing high-pass filter to the cleaned Q and U polarization maps, we inpaint the region with the brightest Galactic residuals (5% of the sky) with the same simple algorithm used for point source holes. This is done to avoid introducing ringing around the Galactic centre when the maps are filtered.

In addition, E- and B-mode maps are constructed from the clean Q and U maps. In this case, prior to performing the decomposition, the region of the Q and U maps defined by the SEVEM confidence mask is filled with a Gaussian-constrained realization (Eriksen et al. 2004).

Appendix C.3: Masks

In temperature, the SEVEM confidence mask is produced by thresholding differences obtained between different CMB reconstructions. In particular, we construct the difference map between the clean 217 and 143 GHz maps at a resolution of FWHM 30′ and N_side = 256. The brightest pixels of this map (and their direct neighbours) are successively removed from the CMB combined map (at the same resolution) and the dispersion of the CMB combined map calculated. If a sufficient number of pixels is removed, the dispersion of the CMB map goes down and reaches a plateau, indicating that convergence has been achieved. The removed pixels constitute the mask. This mask is then smoothed with a Gaussian beam of 1° to avoid sharp edges, and upgraded to full resolution. The same procedure is repeated for the other two difference maps: the cleaned 143−100 map; and the difference of two clean CMB combined maps, whose linear coefficients have been obtained by minimizing the variance outside two different masks. Finally, the three masks are multiplied in order to produce our final confidence mask, which leaves a suitable sky fraction of approximately 85%. Residual monopoles and dipoles outside this mask are subtracted for the single frequency and combined cleaned maps.

For polarization, the clean combined map is downgraded to a resolution of FWHM 90′ and N_side = 128. The dispersion at each pixel is estimated from a circle centred in the pixel being considered. Those pixels with a dispersion above a given threshold are included in the mask, which is then smoothed with a Gaussian beam of FWHM 90′ to avoid sharp edges and upgraded to the required N_side. Finally, this mask is multiplied by a mask customized to cover the CO emission, in order to discard those pixels contaminated by this foreground component due to the bandpass leakage. An additional 1% of the sky is added to the final mask to remove those pixels most affected by the high-pass filtering that is subsequently applied to the cleaned maps. The final mask allows for a useful fraction of the sky of approximately 80%. The 2015 SEVEM masks are shown in Fig. C.2.

Appendix D: Spectral matching

Spectral Matching Independent Component Analysis (SMICA) is a semi-blind component separation algorithm that operates in harmonic space. CMB maps are synthesized from spherical harmonic coefficients, s_ℓm, obtained as weighted linear combinations of the coefficients of N_chan input maps, (D.1)where x_ℓm is the N_chan × 1 vector of the spherical harmonic coefficients of the input maps and w_ℓ is the N_chan × 1 vector of weights. Here and in the following, for the sake of exposition, a real-valued spherical harmonic basis is used (it is always possible to do so) so that all the quantities are real. In practice, we use the standard complex basis, but this is immaterial. The spectral weights used to produce the SMICA maps are designed to minimize the total foreground and noise contamination at each multipole, under the constraint that the resulting map has a well defined effective beam window function, that of a Gaussian beam with 5′ FWHM.

In theory, such weights are given by (D.2)where the N_chan × 1 vector a is the frequency spectrum of the CMB and the N_chan × N_chan spectral covariance matrix R_ℓ contains in its (i,j) entry the cross-power spectrum at multipole ℓ between the input frequency maps i and j, at 5′ resolution. In practice, the spectral covariance matrices must be estimated from the data. The natural sample estimate (D.3)can be used, possibly after some binning, to replace R_ℓ in the weight formula, Eq. (D.2). This works well at high ℓ because a large number of modes are averaged in Eq. (D.3), so this estimate has low variance. At larger scales, however, it is necessary to constrain the spectral covariance matrices in order to obtain reliable estimates and mitigate the effects of chance correlations. For that purpose, SMICA uses a semi-parametric model of these matrices.

The CMB, the foregrounds, and the noise are independent processes, so the spectral covariance matrices, after beam correction, can be decomposed into (D.4)We assume that the noise is uncorrelated between frequency channels, therefore the noise contribution is diagonal. The signal parts are modelled as (D.5)where C_ℓ is the CMB power spectrum, matrix F is N_chan × d, and P_ℓ is a d × d positive definite matrix. A SMICA fit consists of fitting the model of Eqs. (D.4) and (D.5) to empirical covariance matrices by minimizing the spectral matching criterion The fitted values of R_ℓ can be seen as regularized versions of their empirical counterparts , to be used in Eq. (D.2). If no constraints are imposed on matrices F and P_ℓ, except that the latter is positive definite, the foreground model is equivalent to d sky templates with arbitrary frequency spectra (represented by the columns of F) and arbitrary power spectra and correlations between them (represented by P_ℓ). Ultimately, the foreground contribution is controlled by a single parameter, the rank d of the foreground model.

More details regarding the principles of SMICA can be found in Cardoso et al. (2008). We now describe its extension to polarization data.

There are several options to extend SMICA to polarization. We choose to obtain the CMB Q and U maps through a joint processing of the E- and B-modes. The CMB U and Q maps are synthesized from their E- and B-modes, and , obtained as linear combinations (D.6)where the N_chan × 1 vectors and are the spherical harmonic coefficients of the input maps and W_ℓ is a 2N_chan × 2 matrix of weights. The optimal weights are obtained as a simple generalization of the single-field case, (D.7)where R_ℓ now is a 2N_chan × 2N_chan covariance matrix. The weights defined by Eq. (D.7) can be safely obtained at high ℓ by using the sample covariance matrices, possibly after binning them. At low multipoles, some regularization via modelling is in order, as for the temperature analysis. We use a natural extension of temperature Eq. (D.5), in the form where F_E and F_B are N_chan × d matrices and P_ℓ is a 2d × 2d matrix. For this release, no constraints are imposed on those matrices, except that each P_ℓ must be positive definite. Hence, as in temperature, SMICA fits a foreground model representing d polarized templates with arbitrary frequency spectra, power spectra, and correlations.

Appendix D.1: Implementation for temperature

The production of the CMB temperature map is mostly unchanged from 2013 (see Planck Collaboration XII 2014 for details). Here, we recall the 3-step fitting strategy adopted in 2013 and mention the changes made for this release.

• First fit: recalibration. A preliminary and independent SMICA fit is done with power spectra estimated over a clean fraction of the sky, including a as a free parameter. This can be understood as a recalibration procedure; the estimated value of a is kept fixed in the following steps. The number of foreground templates, d, is now 5, whereas it was 4 in 2013. The 30 GHz channel is not recalibrated, unlike in 2013.

• Second fit: foreground emission. The foreground emission, matrix F, is estimated in a second SMICA fit with a kept fixed at the value found in the first step. Spectra are estimated over a large fraction of the sky (f_sky = 97%) and the fit is made over the range 4 ≤ ℓ ≤ 150. For the data, this step is the same as in 2013. For FFP8, foreground emission is captured using seven templates, whereas six templates were used for the FFP6 simulations in 2013.

• Third fit: power spectra. The frequency spectra captured by vector a and matrix F are kept fixed at the values found in steps 1 and 2; we fit only the signal power spectra C_ℓ and P_ℓ and the frequency channel noise power spectra. Spectral covariance matrices are computed over 97% of the sky and the fit includes all multipoles up to ℓ = 1500. This step is the same as in 2013.

The weights determined for the temperature maps are shown at the top of Fig. D.1. It should be noted that the CMB map is synthesized from spherical harmonic coefficients which have been set explicitly to zero for ℓ = 0 and ℓ = 1. Therefore, the SMICA CMB map has no monopole or dipole components over the full sky.

Appendix D.2: Implementation for polarization

All Planck polarized channels are used to produce the polarized CMB maps. The process is easier in polarization than it is in temperature because less precision is required due to the lower signal-to-noise ratio and also because the foregrounds appear to have a simpler structure. In particular, we do not preprocess the frequency maps by subtracting or masking point sources, as is done to the temperature maps.

The SMICA fit in polarization is conducted with the same parameters as in temperature, but with two differences. First, the recalibration step is omitted, and instead we use the CMB frequency spectrum (vector a) determined from temperature maps. Second, the foreground model comprises six polarized templates (as for temperature) but the matrices F_E and F_B are fitted in the second step over the multipole range 4 ≤ ℓ ≤ 50. The weights determined for the polarization data are shown at the bottom of Fig. D.1.

In order to mitigate spectral leakage from E to B and from the Galactic plane onto other regions of the sky, we do not compute spherical transforms directly from masked Q and U maps. Instead, for each pair of input Q and U maps, we first produce full-sky E and B maps to which an apodized Galactic mask is applied. It is from those masked E and B maps that spherical harmonic coefficients are computed for the estimation of the spectra and cross-spectra going into and for the synthesis of the final CMB map.

Appendix D.3: Masks

Confidence masks, shown in Fig. D.2, are built using the procedure described in Planck Collaboration XII (2014), with the following changes. For temperature, we apply a bandpass filter to the CMB map that is then squared and smoothed at 3.̊5 (it was 2° in 2013). The confidence mask is obtained by thresholding the resulting map of local power. The threshold is determined by visual inspection to be 50 μK² (it was 70 μK² in 2013). The resulting mask is then enlarged by multiplication by a Galactic mask covering 10% of the sky. For polarization, the confidence mask is obtained by a procedure similar to temperature, but the CMB Q and U maps are low-pass filtered (rather than bandpassed) using a Gaussian beam with 30′ FWHM. They are then squared and smoothed to 3.̊5 resolution. Any area where the resulting map is above 5 μK² is excluded from the confidence mask. In addition to that, we exclude pixels that are less than 7° away from the Galactic equator.

Appendix E: FFP8 simulations

In this Appendix we provide a compendium of analyses evaluated from the FFP8 data set, corresponding directly to those performed on the Planck 2015 data in the main text.

Appendix E.1: CMB map differences

Figures 4, 9, and 10 show pairwise difference maps between any of the Commander, NILC, SEVEM, and SMICA maps derived from the Planck 2015 data for each of the three Stokes parameters, I, Q, and U. In Figs. E.1–E.3 we show the same quantities evaluated from the FFP8 simulation set.

Overall, the relative differences are of similar magnitude in the simulated data set as in the data, although they have somewhat different dominant morphology. The high-latitude differences are in general somewhat weaker in FFP8 than in the data, while the low-latitude differences are somewhat stronger. This is primarily due to the rather complicated foreground model adopted for thermal dust emission in the FFP8 simulations (see Planck Collaboration XII 2016 for details). In short, thermal dust emission is modelled in the FFP8 simulations by a sum of two modified blackbody components, one of which has a temperature-dependent spectral index, β(T_d). This was motivated by the results presented by Planck Collaboration Int. XIV (2014); however, as shown by the updated analysis in Planck Collaboration X (2016), there is no evidence for steepening from the Planck data alone. Only when the IRAS 100μm data are included in the fit is any such effect seen. At such high frequencies the thermal dust physics is far more complicated, and the overall calibration problem more difficult. For the present discussion, it is sufficient to note that the thermal dust model adopted for the FFP8 simulations is significantly more complicated in terms of frequency dependence than the observed sky, and includes significant spectral curvature below 353 GHz that is not seen in the actual data. None of the component separation codes account for this additional curvature, and this results in the strong features near the Galactic plane seen in some of the pairwise difference maps shown in Fig. 4.

Appendix E.2: 1-point statistics of the total intensity FFP8 maps

Next, we consider the 1-point statistics of the FFP8 fiducial map in temperature, using the FFP8-UT74 mask (Table E.1). These statistics were defined in Sect. 7.1, and corresponding polarization results were shown in Table 3. We see that all the methods are in good agreement with the Monte Carlo simulations at all resolutions considered, and the differences among codes are small.

Appendix E.3: The real-space N-point correlation functions for FFP8 maps

Finally, we present the real-space N-point correlation functions derived from the FFP8 simulations, complementing the analysis in Sect. 7.2. As we did for the data, we analyse the FFP8 CMB temperature and high-pass filtered polarization estimates downgraded to N_side = 64 and smoothed with a Gaussian kernel of 160′ FWHM. For the polarization analysis we employ a low resolution version of the common polarization mask, and for the temperature analysis we employ a low resolution version of the common temperature mask.

The N-point functions for the FFP8 simulations are shown in Figs. E.4, plotted in terms of differences between the N-point functions for the high-pass filtered fiducial FFP8 map and the mean value estimated from 1000 Monte Carlo simulations. The probabilities of obtaining larger values for the χ² statistic, compared to the Planck fiducial ΛCDM model, are tabulated in Table E.2.

Overall, there is good agreement between the recovered FFP8 CMB estimates and the corresponding Monte Carlo ensembles. The main outlier is NILC, for which the fluctuations are larger than expected for all 3-point polarization functions.

Appendix E.4: CMB Power spectra and residuals from individual components

Here we discuss the recovery of the power spectrum in the FFP8 simulations, and assess the impact of foregrounds in the CMB solutions presented in this paper. CMB angular power spectra from the foreground cleaned CMB maps from FFP8 are shown in Fig. E.5, compared to the input for TT and EE. The performance on simulations is comparable to that on the data, as may be seen by comparison with Fig. 14.

Figures E.6 and E.7 show the residual effect of combined individual components on the foreground-cleaned CMB maps. These plots give an estimate of how the residuals of the foreground components in the CMB solutions presented in this paper compare in power and are distributed in ℓ. The FFP8 simulated skies of the individual components are propagated through the various component-separation pipelines, and their power spectra are calculated using the FFP8-UT74 sky mask. The individual components are those simulated in FFP8 (Sect. 2.2 and Planck Collaboration XII 2016) and include in TT dust, CIB fluctuations, unresolved sources (PS), and the sum of the remaining components (“other”). In EE, they include synchrotron, dust, noise, and the sum of the remaining components. With the exception of the lowest multipoles, the dominant contaminant to the TT signal is noise. For EE, the residuals from dust dominate over the other foregrounds, but remain sub-dominant with respect to instrumental noise at all multipoles.

All Tables

All Figures