Issue 
A&A
Volume 594, October 2016
Planck 2015 results



Article Number  A22  
Number of page(s)  24  
Section  Cosmology (including clusters of galaxies)  
DOI  https://doi.org/10.1051/00046361/201525826  
Published online  20 September 2016 
Planck 2015 results
XXII. A map of the thermal SunyaevZeldovich effect
^{1} APC, AstroParticule et Cosmologie, Université Paris Diderot, CNRS/IN2P3, CEA/lrfu, Observatoire de Paris, Sorbonne Paris Cité, 10 rue Alice Domon et Léonie Duquet, 75205 Paris Cedex 13, France
^{2} Aalto University Metsähovi Radio Observatory and Dept of Radio Science and Engineering, PO Box 13000, 00076 Aalto, Finland
^{3} African Institute for Mathematical Sciences, 6–8 Melrose Road, 7945 Muizenberg, Cape Town, South Africa
^{4} Agenzia Spaziale Italiana Science Data Center, via del Politecnico snc, 00133 Roma, Italy
^{5} AixMarseille Université, CNRS, LAM (Laboratoire d’Astrophysique de Marseille) UMR 7326, 13388 Marseille, France
^{6} Astrophysics Group, Cavendish Laboratory, University of Cambridge, J J Thomson Avenue, Cambridge CB3 0HE, UK
^{7} Astrophysics & Cosmology Research Unit, School of Mathematics, Statistics & Computer Science, University of KwaZuluNatal, Westville Campus, Private Bag X54001, 4000 Durban, South Africa
^{8} Atacama Large Millimeter/submillimeter Array, ALMA Santiago Central Offices, Alonso de Cordova 3107, Vitacura, 763 0355 Casilla, Santiago, Chile
^{9} CITA, University of Toronto, 60 St. George St., Toronto, ONM5S 3H8, Canada
^{10} CNRS, IRAP, 9 Av. colonel Roche, BP 44346, 31028 Toulouse Cedex 4, France
^{11} CRANN, Trinity College, Dublin 2, Ireland
^{12} California Institute of Technology, Pasadena, CA 91125, USA
^{13} Centre for Theoretical Cosmology, DAMTP, University of Cambridge, Wilberforce Road, Cambridge CB3 0WA, UK
^{14} Centro de Estudios de Física del Cosmos de Aragón (CEFCA), Plaza San Juan, 1, Planta 2, 44001 Teruel, Spain
^{15} Computational Cosmology Center, Lawrence Berkeley National Laboratory, Berkeley, CA 94720, USA
^{16} DSM/Irfu/SPP, CEASaclay, 91191 GifsurYvette Cedex, France
^{17} DTU Space, National Space Institute, Technical University of Denmark, Elektrovej 327, 2800 Kgs. Lyngby, Denmark
^{18} Département de Physique Théorique, Université de Genève, 24 quai E. Ansermet, 1211 Genève 4, Switzerland
^{19} Departamento de Física, Universidad de Oviedo, Avda. Calvo Sotelo s/n, 33003 Oviedo, Spain
^{20} Department of Astrophysics/IMAPP, Radboud University Nijmegen, PO Box 9010, 6500 GL Nijmegen, The Netherlands
^{21} Department of Physics & Astronomy, University of British Columbia, 6224 Agricultural Road, Vancouver, British Columbia, Canada
^{22} Departmentof Physics and Astronomy, Dana and David Dornsife College of Letter, Arts and Sciences, University of Southern California, Los Angeles, CA 90089, USA
^{23} Department of Physics and Astronomy, University College London, London WC1E 6BT, UK
^{24} Department of Physics, Florida State University, Keen Physics Building, 77 Chieftan Way, Tallahassee, Florida, USA
^{25} Department of Physics, Gustaf Hällströmin katu 2a, University of Helsinki, 00014 Helsinki, Finland
^{26} Department of Physics, Princeton University, Princeton, NJ 08544, USA
^{27} Department of Physics, University of California, Santa Barbara, CA 93106, USA
^{28} Department of Physics, University of Illinois at UrbanaChampaign, 1110 West Green Street, Urbana, Illinois, USA
^{29} Dipartimento di Fisica e Astronomia G. Galilei, Università degli Studi di Padova, via Marzolo 8, 35131 Padova, Italy
^{30} Dipartimento di Fisica e Scienze della Terra, Università di Ferrara, via Saragat 1, 44122 Ferrara, Italy
^{31} Dipartimento di Fisica, Università La Sapienza, P.le A. Moro 2, 00185 Roma, Italy
^{32} Dipartimento di Fisica, Università degli Studi di Milano, via Celoria 16, 20133 Milano, Italy
^{33} Dipartimento di Fisica, Università degli Studi di Trieste, via A. Valerio 2, 34127 Trieste, Italy
^{34} Dipartimento di Matematica, Università di Roma Tor Vergata, via della Ricerca Scientifica 1, 00133 Roma, Italy
^{35} Discovery Center, Niels Bohr Institute, Blegdamsvej 17, 2100 Copenhagen, Denmark
^{36} Dpto. Astrofísica, Universidad de La Laguna (ULL), 38206 La Laguna, Tenerife, Spain
^{37} European Southern Observatory, ESO Vitacura, Alonso de Cordova 3107, Vitacura, Casilla 19001, Santiago, Chile
^{38} European Space Agency, ESAC, Planck Science Office, Camino bajo del Castillo, s/n, Urbanización Villafranca del Castillo, 28692 Villanueva de la Cañada, Madrid, Spain
^{39} European Space Agency, ESTEC, Keplerlaan 1, 2201 AZ Noordwijk, The Netherlands
^{40} Facoltà di Ingegneria, Università degli Studi eCampus, via Isimbardi 10, 22060 Novedrate (CO), Italy
^{41} Gran Sasso Science Institute, INFN, viale F. Crispi 7, 67100 L’ Aquila, Italy
^{42} HGSFP and University of Heidelberg, Theoretical Physics Department, Philosophenweg 16, 69120 Heidelberg, Germany
^{43} Helsinki Institute of Physics, Gustaf Hällströmin katu 2, University of Helsinki, 00014 Helsinki, Finland
^{44} ICTP South American Institute for Fundamental Research, Instituto de Física Teórica, Universidade Estadual Paulista, 01140070 São Paulo, Brazil
^{45} INAF–Osservatorio Astronomico di Padova, Vicolo dell’Osservatorio 5, 35122 Padova, Italy
^{46} INAF–Osservatorio Astronomico di Roma, via di Frascati 33, 00040 Monte Porzio Catone, Italy
^{47} INAF–Osservatorio Astronomico di Trieste, via G.B. Tiepolo 11, 40127 Trieste, Italy
^{48} INAF/IASF Bologna, via Gobetti 101, 40129 Bologna, Italy
^{49} INAF/IASF Milano, via E. Bassini 15, 20133 Milano, Italy
^{50} INFN, Sezione di Bologna, via Irnerio 46, 40126 Bologna, Italy
^{51} INFN, Sezione di Roma 1, Università di Roma Sapienza, P.le Aldo Moro 2, 00185 Roma, Italy
^{52} INFN, Sezione di Roma 2, Università di Roma Tor Vergata, via della Ricerca Scientifica 1, 00185 Roma, Italy
^{53} INFN/National Institute for Nuclear Physics, via Valerio 2, 34127 Trieste, Italy
^{54} IPAG: Institut de Planétologie et d’Astrophysique de Grenoble, Université Grenoble Alpes, IPAG; CNRS, IPAG, 38000 Grenoble, France
^{55} Imperial College London, Astrophysics group, Blackett Laboratory, Prince Consort Road, London, SW7 2AZ, UK
^{56} Infrared Processing and Analysis Center, California Institute of Technology, Pasadena, CA 91125, USA
^{57} Institut Néel, CNRS, Université Joseph Fourier Grenoble I, 25 rue des Martyrs, Grenoble, France
^{58} Institut Universitaire de France, 103 bd SaintMichel, 75005 Paris, France
^{59} Institut d’Astrophysique Spatiale, CNRS (UMR 8617) Université ParisSud 11, Bâtiment 121, 91405 Orsay, France
^{60} Institut d’Astrophysique de Paris, CNRS (UMR 7095), 98bis boulevard Arago, 75014 Paris, France
^{61} Institute of Astronomy, University of Cambridge, Madingley Road, Cambridge CB3 0HA, UK
^{62} Institute of Theoretical Astrophysics, University of Oslo, Blindern, 0371 Oslo, Norway
^{63} Instituto de Astrofísica de Canarias, C/Vía Láctea s/n, La Laguna, 38205 Tenerife, Spain
^{64} Instituto de Física de Cantabria (CSICUniversidad de Cantabria), Avda. de los Castros s/n, 39005 Santander, Spain
^{65} Istituto Nazionale di Fisica Nucleare, Sezione di Padova, via Marzolo 8, 35131 Padova, Italy
^{66} Jet Propulsion Laboratory, California Institute of Technology, 4800 Oak Grove Drive, Pasadena, CA 31109, USA
^{67} Jodrell Bank Centre for Astrophysics, Alan Turing Building, School of Physics and Astronomy, The University of Manchester, Oxford Road, Manchester, M13 9PL, UK
^{68} Kavli Institute for Cosmology Cambridge, Madingley Road, Cambridge, CB3 0HA, UK
^{69} LAL, Université ParisSud, CNRS/IN2P3, 91898 Orsay, France
^{70} LAPTh, Univ. de Savoie, CNRS, BP 110, 74941 AnnecyleVieux, France
^{71} LERMA, CNRS, Observatoire de Paris, 61 avenue de l’Observatoire, 75014 Paris, France
^{72} Laboratoire AIM, IRFU/Service d’Astrophysique − CEA/DSM − CNRS − Université Paris Diderot, Bât. 709, CEASaclay, 91191 GifsurYvette Cedex, France
^{73} Laboratoire Traitement et Communication de l’Information, CNRS (UMR 5141) and Télécom ParisTech, 46 rue Barrault, 75634 Paris Cedex 13, France
^{74} Laboratoire de Physique Subatomique et Cosmologie, Université GrenobleAlpes, CNRS/IN2P3, 53 rue des Martyrs, 38026 Grenoble Cedex, France
^{75} Laboratoire de Physique Théorique, Université ParisSud 11 & CNRS, Bâtiment 210, 91405 Orsay, France
^{76} Lebedev Physical Institute of the Russian Academy of Sciences, Astro Space Centre, 84/32 Profsoyuznaya st., GSP7, 117997 Moscow, Russia
^{77} MaxPlanckInstitut für Astrophysik, KarlSchwarzschildStr. 1, 85741 Garching, Germany
^{78} National University of Ireland, Department of Experimental Physics, Maynooth, Co. Kildare, Ireland
^{79} Niels Bohr Institute, Blegdamsvej 17, 2100 Copenhagen, Denmark
^{80} Optical Science Laboratory, University College London, Gower Street, WC1E 6 BT London, UK
^{81} SBITPLPPC, EPFL, 1015 Lausanne, Switzerland
^{82} SISSA, Astrophysics Sector, via Bonomea 265, 34136 Trieste, Italy
^{83} School of Physics and Astronomy, Cardiff University, Queens Buildings, The Parade, Cardiff, CF24 3AA, UK
^{84} Sorbonne UniversitéUPMC, UMR7095, Institut d’Astrophysique de Paris, 98bis boulevard Arago, 75014 Paris, France
^{85} Space Research Institute (IKI), Russian Academy of Sciences, Profsoyuznaya Str, 84/32, 117997 Moscow, Russia
^{86} Space Sciences Laboratory, University of California, Berkeley, CA 94720, USA
^{87} Special Astrophysical Observatory, Russian Academy of Sciences, Nizhnij Arkhyz, Zelenchukskiy region, 369167 KarachaiCherkessian Republic, Russia
^{88} SubDepartment of Astrophysics, University of Oxford, Keble Road, Oxford OX1 3RH, UK
^{89} Theory Division, PHTH, CERN, 1211 Geneva 23, Switzerland
^{90} UPMC Univ Paris 06, UMR 7095, 98bis boulevard Arago, 75014 Paris, France
^{91} Université de Toulouse, UPSOMP, IRAP, 31028 Toulouse Cedex 4, France
^{92} University Observatory, Ludwig Maximilian University of Munich, Scheinerstrasse 1, 81679 Munich, Germany
^{93} University of Granada, Departamento de Física Teórica y del Cosmos, Facultad de Ciencias, 18071 Granada, Spain
^{94} University of Granada, Instituto Carlos I de Física Teórica y Computacional, 18071 Granada, Spain
^{95} Warsaw University Observatory, Aleje Ujazdowskie 4, 00478 Warszawa, Poland
^{⋆}
Corresponding author: B. Comis email: comis@lpsc.in2p3.fr, M. Remazeilles, email: mathieu.remazeilles@manchester.ac.uk
Received: 5 February 2015
Accepted: 8 March 2016
We have constructed allsky Compton parameters maps, ymaps, of the thermal SunyaevZeldovich (tSZ) effect by applying specifically tailored component separation algorithms to the 30 to 857 GHz frequency channel maps from the Planck satellite. These reconstructed ymaps are delivered as part of the Planck 2015 release. The ymaps are characterized in terms of noise properties and residual foreground contamination, mainly thermal dust emission at large angular scales, and cosmic infrared background and extragalactic point sources at small angular scales. Specific masks are defined to minimize foreground residuals and systematics. Using these masks, we compute the ymap angular power spectrum and higher order statistics. From these we conclude that the ymap is dominated by tSZ signal in the multipole range, 20 <ℓ< 600. We compare the measured tSZ power spectrum and higher order statistics to various physically motivated models and discuss the implications of our results in terms of cluster physics and cosmology.
Key words: largescale structure of Universe / cosmology: observations
© ESO, 2016
1. Introduction
This paper, one of a set associated with the 2015 release of data from the Planck^{1} mission, describes the Planck Compton parameter (y) map, which is part of the Planck 2015 data release.
The thermal SunyaevZeldovich (tSZ) effect (Sunyaev & Zeldovich 1972) is produced by the inverse Compton scattering of cosmic microwave background (CMB) photons by hot electrons along the line of sight and in particular in clusters of galaxies. The tSZ effect has proved to be a major tool for studying the physics of clusters of galaxies, and also structure formation in the Universe. Catalogues of clusters of galaxies selected via the tSZ effect have become available in the last few years, including, for example, those from the Planck satellite (Planck Collaboration VIII 2011; Planck Collaboration XXIX 2014), the Atacama Cosmology Telescope (ACT, Hasselfield et al. 2013), and the South Pole Telescope (SPT, Reichardt et al. 2013; Bleem et al. 2015). These catalogues and their associated sky surveys have been used to study the physics of clusters of galaxies (Planck Collaboration XII 2011; Planck Collaboration XI 2011; Planck Collaboration X 2011) and their cosmological implications (Planck Collaboration XX 2014; Benson et al. 2013; Das et al. 2013; Wilson et al. 2012; Mak & Pierpaoli 2012).
The study of cluster number counts and their evolution with redshift using tSZ selected catalogues is an important cosmological test (Carlstrom et al. 2002; Dunkley et al. 2013; Benson et al. 2013; Planck Collaboration XX 2014; Hou et al. 2014). This study is complemented by the measurement of the tSZ effect power spectrum (Komatsu & Seljak 2002), for which no explicit measurement of cluster masses is required. Lowmass clusters, thus fainter in tSZ, which can not be detected individually, also contribute statistically to the measured signal (Battaglia et al. 2010; Shaw et al. 2010). Another complementary approach (as pointed out in RubiñoMartín & Sunyaev 2003) is to study the higherorder statistics of the tSZ signal, and in particular the skewness or, equivalently, the bispectrum. The bispectrum of the tSZ effect signal is dominated by massive clusters at intermediate redshifts (Bhattacharya et al. 2012), for which highprecision Xray observations exist. This contrasts with the power spectrum, for which most of the signal comes from the lower mass and higher redshift groups and clusters (e.g. Trac et al. 2011). Thus, theoretical uncertainties in the tSZ bispectrum are expected to be significantly smaller than those in the estimation of the tSZ power spectrum. Therefore, combined measurements of the power spectrum and the bispectrum can be used to distinguish the contribution to the power spectrum from different cluster masses and redshift ranges. However, contamination from point sources (RubiñoMartín & Sunyaev 2003; Taburet et al. 2010) and other foregrounds (Planck Collaboration XX 2014) needs to be dealt with carefully.
As shown in the Planck 2013 results (Planck Collaboration XX 2014), the allsky coverage and unprecedented wide frequency range of Planck allowed us to produce allsky tSZ Compton parameter maps, also referred to as ymaps. From these maps we can obtain an accurate measurement of the tSZ power spectrum at intermediate and large angular scales, for which the tSZ fluctuations are almost insensitive to the cluster core physics. Furthermore, the expected nonGaussianity properties of the tSZ signal can be studied using higher order statistical estimators, such as the skewness, and the bispectrum. The ymap can be also used to extract the tSZ signal on regions centred at cluster positions and in particular to perform stacking analysis. Finally, this map may also be cross correlated with other cosmological probes (see for example Ma et al. 2015; Van Waerbeke et al. 2014; Hill & Spergel 2014) as a consistency test.
In this paper we construct revised tSZ allsky maps from the individual Planck frequency maps. With respect to the Planck Compton parameter map in the Planck 2013 results (Planck Collaboration XX 2014), we have extended the analysis to the full mission data set and performed an indepth characterization of its statistical properties. This extended analysis allows us to deliver this map as part of the Planck 2015 data release. The paper is structured as follows. Section 2 describes the Planck data used to compute the tSZ allsky maps and the simulations used to characterize them. Section 3 presents the reconstruction of the Planck allsky Compton parameter map. In Sect. 4 we discuss the validation of the reconstructed ymap in the pixel domain including signal and noise characterization. Section 5 describes the power spectrum analysis. Crosschecks using higher order statistics are presented in Sect. 6. The cosmological interpretation of the results is discussed in Sect. 7, and we finally present our conclusions in Sect. 8.
2. Data and simulations
2.1. The Planck data
This paper is based on the Planck’s full mission (Planck Collaboration I 2014; Planck Collaboration I 2016), corresponding to five and eight fullsky surveys for the High Frequency Instrument (HFI) and Low Frequency Instrument (LFI) data (Planck Collaboration II 2016; Planck Collaboration VII 2016; Planck Collaboration VIII 2016), respectively. The Planck channel maps are provided in HEALPix (Górski et al. 2005) pixelization scheme at N_{side} = 2048. A noise map is associated with each channel map. This map is obtained from the half difference of maps made from the first and second halves of each stable pointing period (also called “ring”). In the halfdifference maps the astrophysical emission cancels out, which makes them a good representation of the statistical instrumental noise. These halfdifference maps are used to estimate the noise in the final Compton parameter map. In addition, survey maps, which are also available for each channel, will be used to estimate possible residual systematic effects in the ymap.
For the purpose of this paper we approximate the Planck effective beams by circular Gaussians, the FWHM estimates of which are given in Table 1 for each frequency channel.
2.2. Simulations
We also use simulated Planck frequency maps obtained from the full focal plane (FFP) simulations (Planck Collaboration XII 2016), which are described in the Planck Explanatory supplement (Planck Collaboration 2013).
These simulated maps include the most relevant sky components at microwave and millimetre frequencies, based on foregrounds from the Planck sky model (PSM, Delabrouille et al. 2013), namely CMB anisotropies, thermal SZ effect, diffuse Galactic emissions (synchrotron, freefree, thermal and spinning dust and CO), radio and infrared point sources, and the clustered cosmic infrared background (CIB). The tSZ signal was constructed using hydrodynamical simulations of clusters of galaxies up to redshift 0.3. For higher redshifts pressure profilebased simulations of individual clusters of galaxies randomly drawn on the sky were added. The noise in the maps was obtained from realizations of Gaussian random noise in the time domain and therefore accounts for noise inhomogeneities in the maps.
Conversion factors for tSZ Compton parameter y to CMB temperature units and the FWHM of the beam of the Planck channel maps.
The simulation set also includes Monte Carlo noiseonly realizations for each Planck channel map. These will also be used to estimate the noise properties in the final ymap.
3. Reconstruction of the allsky tSZ maps
The thermal SZ Compton parameter (Sunyaev & Zeldovich 1972) in a given direction, , can be expressed as (1)where k_{B} is the Boltzmann contsant, m_{e} the electron mass, σ_{T} the Thomson crosssection, ds the distance along the line of sight, , and n_{e} and T_{e} are the electron number density and temperature.
In CMB temperature units the tSZ effect contribution to the Planck maps for a given observational frequency ν is given by (2)Neglecting relativistic corrections, g(ν) = x coth(x/ 2)−4, with x = hν/ (k_{B}T_{CMB}). Table 1 presents the conversion factors for Compton parameter to CMB temperature, K_{CMB}, for each frequency channel, after integrating over the instrumental bandpass (Planck Collaboration IX 2014).
3.1. Reconstruction methods
The tSZ effect signal in the Planck frequency maps is subdominant with respect to the CMB and other foreground emissions. By contrast to the CMB emission, the tSZ effect from galaxy clusters is spatially localized and leads to a highly nonGaussian signal. CMBoriented componentseparation methods (Planck Collaboration XII 2014) are not adequate to recover the tSZ signal. We therefore use specifically tailored component separation algorithms to reconstruct the tSZ signal from the Planck frequency channel maps, as in Planck Collaboration XX (2014). These algorithms rely on the spatial localization of the different astrophysical components and on their spectral diversity to separate them. As in Planck Collaboration XX (2014) we consider here the MILCA (Modified Internal Linear Combination Algorithm, Hurier et al. 2013) and NILC (Needlet Independent Linear Combination, Remazeilles et al. 2011) methods. Both are based on the Internal Linear Combination (ILC) approach that searches for the linear combination of the input maps that minimizes the variance of the final reconstructed map, under the constraint of offering unit gain to the component of interest (here the tSZ effect, whose frequency dependence is known). Both algorithms have been extensively tested on simulated Planck data.
For both methods, the Planck HFI maps from 100 to 857 GHz are used convolved to a common resolution of 10′. In the case of NILC we also use the LFI data at large angular scales (ℓ< 300). Similarly, for both methods the 857 GHz map, which traces the thermal dust emission on large angular scales, is only used for multipoles ℓ< 300 to minimize residuals from IR point sources and clustered CIB emission. The final ymaps have a resolution of 10′.
Fig. 1 Window functions corresponding to the spectral localization of the NILC and MILCA algorithms. For NILC (solid line) there are 10 Gaussian window functions defining 10 needlet scales. MILCA (dashed line) uses 11 overlapping Gaussian windows. 
3.1.1. NILC
In the multicomponent extensions of NILC (Delabrouille et al. 2009; Remazeilles et al. 2011), initially developed to extract the CMB, the weights for component separation (i.e., covariances) are computed independently in domains of a needlet decomposition (in the spherical wavelet frame). The needlet decomposition provides localization of the ILC filters both in pixel and in multipole space, allowing us to deal with local contamination conditions varying both in position and in scale. We imposed constraints to remove the CMB contamination and preserve the tSZ effect. To avoid strong foreground effects, part of the Galactic plane (corresponding to about 2% of the sky) was masked before applying NILC to the Planck frequency maps.
The localization in multipole space is achieved by using ten Gaussian window functions { h_{j}(ℓ) } _{1 ≤ j ≤ 10} as bandpass filters^{2}, socalled “needlet bands”, allowing for smooth localization in ℓ (see Fig. 1, solid lines). NILC performs a weighted linear combination of the bandpass filtered Planck maps for each needlet scale independently. The localization in the spatial domain is achieved by defining scaledependent zones over the sky on which the covariance matrices and ILC weights are computed. More precisely, the pixel domain on which the covariance is computed is defined from the smoothing of the product of the relevant needlet maps with a symmetric Gaussian window in pixel space whose FWHM depends on the needlet scale considered. This avoids artificial discontinuities at pixel edges. The localization reduces the number of modes on which the statistics is computed; this may be responsible for an “ILC bias” due to chance correlations between SZ and foregrounds (Delabrouille et al. 2009). At the coarsest scale in particular (first needlet band), the area of the spatial localization must be large enough to counterbalance the lack of modes in multipole space due to spectral localization. In practice, the zones for spatial localization are not predefined but their area is automatically adjusted to the needlet scale considered. The ILC bias b_{ILC} is related to the number of channels N_{ch} and to the number of modes N_{m}Delabrouille et al. (2009) as (3)This offers the possibility to adapt N_{m} and N_{ch} in order to control the ILC bias to a set threshold, as discussed in Remazeilles et al. (2013): given both the number of channels and the number of spectral modes in ℓspace at the needlet scale considered, our NILC algorithm consistently computes the number of spatial modes (similarly, the zone area for spatial localization) required to control the ILC bias.
Fig. 2 Reconstructed Planck allsky Compton parameter maps for NILC (top) and MILCA (bottom) in orthographic projections. The apparent difference in contrast observed between the NILC and MILCA maps comes from differences in the residual foreground contamination and from the differences in the filtering applied for display purposes to the original Compton parameter maps. For the MILCA method, filtering out low multipoles significantly reduces the level of foreground emission in the final ymap. The wavelet basis used in the NILC method was tailored for tSZ extraction. For details see Planck Collaboration XXII (2016). 
Fig. 3 A small region of the reconstructed Planck allsky Compton parameter maps for NILC (left) and MILCA (right) at intermediate Galactic latitudes in the southern sky centred at (0°,−45°) in Galactic coordinates. The colour scale is in units of y × 10^{6}. 
Fig. 4 Inscan and crossscan contributions in the NILC (top line) and MILCA (bottom line) ymaps in Compton parameter units times 10^{6}. From left to right we present the original ymaps, and their in and cross scan contributions for a small region at intermediate Galactic latitudes in the southern sky centred at (0°,−45°) in Galactic coordinates. 
3.1.2. MILCA
MILCA (Hurier et al. 2013), when applied to the extraction of the tSZ signal, also uses two spectral constraints: preservation of the tSZ signal (the tSZ spectral signature discussed above is assumed); and removal of the CMB contamination in the final SZ map, making use of the well known spectrum of the CMB. In addition, to compute the weights of the linear combination, we have used the extra degrees of freedom in the linear system to minimize residuals from other components (two degrees of freedom) and from the noise (two additional degrees). The noise covariance matrix is estimated from the halfdifference maps described in Sect. 2.1. To improve the efficiency of the MILCA algorithm, weights are allowed to vary as a function of multipole ℓ, and are computed independently on different sky regions. We have used 11 filters in ℓ space as shown in Fig. 1 (dashed lines). These filters have an overall transmission of unity, except for ℓ< 8. For these large angular scales we have used a Gaussian filter to reduce foreground contamination. To ensure sufficient spatial localization for each required resolution, the size of the independent sky regions was adapted to the multipole range. We used a minimum of 12 regions at low resolution and a maximum of 3072 regions at high resolution.
3.2. Reconstructed Compton parameter map
Figure 2 shows the reconstructed Planck allsky Compton parameter map for NILC (top panel) and MILCA (bottom panel). For display purposes, the maps are filtered using the procedure described in the first paragraph of Sect. 6.1. Clusters appear as positive sources, the Coma cluster and Virgo supercluster clearly visible near the north Galactic pole. The Galactic plane is masked in both maps, leaving 67% of the sky. Residual Galactic contamination is also visible as diffuse positive structures in the MILCA ymap. We can also observe a granular structure in the NILC ymap that corresponds to an excess of noise at large angular scales, as discussed in Sect. 4.2.
Weaker and more compact clusters are visible in the zoomed region of the Southern cap, shown in Fig. 3. Strong Galactic and extragalactic radio sources show up as negative bright spots on the maps and were masked prior to any scientific analysis, as discussed below in Sect. 4.4.1. We can again observe residual Galactic contamination around the edges of the masked area, which is more important for MILCA. Finally, we note in the NILC and MILCA ymaps the presence of systematic residuals along the scanning direction that show up as stripes. These are the consequence of stripes in the original Planck maps (see Planck Collaboration VIII 2016, for a detailed discussion). We discuss the level of residual stripes in Sect. 4.1.
In addition to the full Compton parameter maps, we also produce the socalled first and last (“F” and “L” hereafter) Compton parameter maps from the first and second halves of the survey rings (i.e., pointing periods). These maps are used for power spectrum estimation in Sect. 5 as well as to estimate the noise properties in the ymaps (see Sect. 4.2).
3.3. Comparison to other Compton parameter maps in the literature
Van Waerbeke et al. (2014) and Hill & Spergel (2014) have used the Planck nominal data to reconstruct a Compton parameter map over a large fraction of the sky. As for the present paper, they use an ILC approach, imposing spectral constraints to have unitary gain for the SZ component and null contribution from CMB, but they do not use the spectral and spatial separation discussed above. Van Waerbeke et al. (2014) consider only the four HFI channels between 100 and 353 GHz and force a null contribution from dust emission. Hill & Spergel (2014) also include the 545 GHz map, while using the 857 GHz channel only as a template for dust emission: a flux cut is imposed in order to build a mask that keeps only the 30% of the sky. They also apply a point source mask (radio and IR) before computing the ILC coefficients, reducing the final sky fraction to about the 25%. Furthermore, Van Waerbeke et al. (2014) and Hill & Spergel (2014) aim mainly at studying the crosscorrelation with the gravitational lensing mass map from the CanadaFranceHawaii Telescope Lensing Survey (CFHTLenS) and the publiclyreleased Planck CMB lensing potential map, respectively. Their ymap reconstruction methods are tailored to fulfil these objectives and as a consequence they give larger overall foreground contamination, and in the case of Hill & Spergel (2014) also a significantly smaller sky fraction.
4. Pixel space analysis
4.1. Stripes in the ymap
As discussed in Sect. 3.2, residual stripes are visible in the NILC and MILCA ymaps. These stripes are mainly due to inscan direction residual systematics after subtraction of an offset for each Planck stable pointing period (Planck Collaboration VIII 2016). To study how these stripes contaminate the final ymaps we have decomposed them in their inscan and crossscan direction contributions. We first convert the ymaps from Galactic to Ecliptic coordinates for which the scan direction corresponds mainly to a fixed longitude value. Secondly, we apply a Galactic mask to the Ecliptic ymaps and decompose them into spherical harmonics. Third, we select the inscan (crossscan) direction components by nullifying the spherical harmonic coefficients, a_{ℓ,m}, for  m  >ℓ/ 2 ( m  < = ℓ/ 2). Finally, we construct maps from those transformed coefficients and convert back to Galactic coordinates. We have chosen to mask the brightest 40% of the sky in the 857 GHz Planck map to keep Galactic ringing residuals negligible and a sufficiently large fraction of the sky for the analysis.
Figure 4 shows (from left to the right) the original, inscan, and crossscan ymaps for NILC (top) and MILCA (bottom) in the southern sky region presented in Fig. 3. We note that the ymaps look noisier since they are not filtered using the procedure described in the first paragraph of Sect. 6.1. The stripes are apparent in both the original and inscan ymaps, as expected. We find that the ratio of the rms of the in and crossscan maps is 1.16 (1.17) for NILC (MILCA), consistent with residual stripe contamination (for Gaussian noise only maps we find a ratio of 1). Here we have measured an overall increase of the rms of the maps due to stripe contamination; a more detailed estimate of the effect in terms of tSZ power spectrum is presented later in Sect. 5.2.1.
The inscan and crossscan decomposition method modifies the cluster signal significantly due to ringing effects – see positive and negative patterns around clusters in Fig. 4. To explore this effect we have also applied the inscan and crossscan decomposition method to the simulated Compton parameter map for the detected and confirmed clusters of galaxies in the Planck catalogue (Planck Collaboration XXIX 2014), which is presented in Sect. 5.3. We find that those negative and positive patterns are also present in the decomposed maps and the ratio between the rms of the inscan and crossscan maps is 1, as expected. We finally stress that the inscan and crossscan decomposition of the ymaps does not preserve the positiveness of the tSZ signal. As a consequence we can not use them to estimate the contamination of the stripe systematic effect in the analysis presented in Sect. 6.
Fig. 5 Top: standard deviation maps for the NILC (top) and MILCA (middle) ymaps corresponding to the inhomogenous noise contribution computed from the halfdifference of the halfring maps. Bottom: angular power spectrum of the homogenous noise contribution for the NILC (solid grey line) and MILCA (dotted black line) ymaps (see main text for details). 
4.2. Noise distribution on the ymap
The Planck maps present a highly nonhomogeneous structure for the noise that is mainly due to the inhomogeneous scanning strategy. Such complicated structure is propagated into the ymap and needs to be considered for further analysis. Here, we have chosen to describe the noise structure in the ymap by: (A) a variance map, which will capture the inhomogeneity of the noise; and (B) the angular power spectrum of a homogenous noise map, which is obtained after correcting for inhomogeneity. The variance map and homogeneous noise power spectrum have been obtained using two different methods, which we now describe.

1.
Halfdifference. The halfdifference of the halfring ymaps is used to obtain an estimate of the noise in the NILC and MILCA ymaps. From this halfdifference map we compute an estimate of the noise variance at lower resolution by squaring and downgrading it to Healpix N_{side} = 128. The homogenous noise contribution is obtained by dividing the halfdifference map by the variance map after upgrading it to Healpix N_{side} = 2048.

2.
Simulationsbased. 100 noise only realizations from the simulations described in Sect. 2.2 are used to compute an estimate of the variance per pixel at full resolution, which is then averaged to Healpix N_{side} = 128 resolution. The homogenous noise contribution map is then computed as above.
The two methods give consistent results at high Galactic latitudes, but differ at low Galactic latitudes, where we expect larger foreground residuals that show up in the halfdifference method. This is visible in the top panel of Fig. 5 where we display variance maps obtained from the halfdifference method for NILC (left) and MILCA (right). For a wide range of analyses with the ymap (e.g., computing the radial profile, surface density or total flux of a cluster, stacking of faint sources, detection of shocks, and study of the ICM) an estimate of the overall uncertainties per pixel (including foreground contribution) is required. Thus, the halfdifference maps are released along with the ymap. We also show in the bottom panel of Fig. 5 the power spectrum of the homogeneous noise contribution. We observe a significant large angular scale (low multipole) component in the homogeneous noise contribution. This low multipole component is significantly larger for the NILC ymap than for the MILCA one.
4.3. Foreground contamination and masking
As discussed above, a residual Galactic foreground contribution is observed at low Galactic latitudes in both the NILC and MILCA ymaps, while being more important for the MILCA ymap (see Sect. 5.2.1 for a more quantitative comparison). This contribution is mainly induced by thermal dust emission, as discussed in detail in Planck Collaboration XXIX (2014), and is therefore mainly associated with the Galactic plane. When dealing with individual objects, or with the stacking of faint objects, we recommend that it is account for by using the variance maps discussed above. For statistical analyses, such as those presented in Sects. 5 and 6, specific masks need to be defined and these are fully discussed in those sections.
4.4. tSZ signal from resolved sources
4.4.1. Point sources
Point source contamination is an important issue for the cosmological interpretation of the Planck Compton parameter map.
In the reconstructed tSZ maps radio sources will appear as negative peaks, while infrared sources will show up as positive peaks, mimicking the signal from clusters. To avoid contamination from these sources we introduce a point source mask (PSMASK, hereafter). This mask is the union of the individual frequency pointsource masks discussed in Planck Collaboration XXVIII (2014). The reliability of this mask was verified by looking at the 1D Probability Density Function (PDF) of the pixel signal (as will be discussed in Sect. 6.1). We found that a more robust mask can be obtained by enlarging the mask size around the strongest radio sources in order to minimize their contribution.
Alternatively, we have also performed a blind search for negative sources in the ymaps using the MHW2 algorithm (LópezCaniego et al. 2006). We have detected 997 and 992 negative sources for NILC and MILCA respectively. We find that 67 and 42 (for NILC and MILCA, respectively) of those detected sources are not masked by the PSMASK. However, most of them (54 for for NILC and 36 for MILCA) are false detections made by the algorithm in regions surrounding very strong positive sources (i.e., galaxy clusters). For those detected as additional negative sources, the PSMASK has been updated to account for them, even though the strongest are already excluded when applying the 50% Galactic mask used for the cosmological analysis below.
For infrared sources, estimating the efficiency of the masking is hampered by the tSZ signal itself. The residual contamination from point sources is discussed in Sects. 5.2 and 6. It is also important to note that the PSMASK may exclude some clusters of galaxies. This is particularly true in the case of clusters with strong central radio sources, such as the Perseus cluster (see Planck Collaboration XXIX 2014).
4.4.2. Blind search for clusters
Following Planck Collaboration XXI (2014) a blind search for tSZ (positive) sources has alsobeen performed on the allsky NILC and MILCA ymaps. We use the IFCAMEX (MHW2, GonzálezNuevo et al. 2006; LópezCaniego et al. 2006) and the single frequency matched filter (MF, Melin et al. 2006) methods. The detected sources with signaltonoise ratio S/N> 4 have been matched to the Planck cluster catalogue (PSZ2, Planck Collaboration XXVII 2016). We have associated the detected sources to PSZ2 objects if the distance between their positions is smaller than 10′ (the resolution of the SZ allsky maps). The MHW2 algorithm finds 1018 and 1522 candidates for NILC and MILCA, respectively, out of which 500 and 457 correspond to objects present within the PSZ2 catalogue. For NILC (MILCA) we have 41 (25) positions that are not masked by the PSMASK and that do not correspond to PSZ2 clusters. However 14 (5) are located in regions excluded by the mask used to build the PSZ2 catalogue. With the MF method we have 1472 and 1502 candidates for NILC and MILCA, out of which 1107 and 1096 correspond to objects present within the PSZ2 catalogue (867 and 835 corresponding to already validated objects, respectively). This implies good agreement between the Planck cluster sample and the detected sources in the ymaps.
This agreement is improved by taking the union of the MHW2 and MF catalogues. In this case, and considering only validated sources in the PSZ2 catalogue (1070 sources), we find 907 and 870 matches with blind detected cluster candidates for NILC and MILCA, respectively. We have also performed a visual inspection of the NILC and MILCA ymaps at the positions of validated PSZ2 sources for which we find no counterpart, and here we find evidence of an excess of signal. For most of these sources low S/N and/or foreground contamination can explain why they are not detected blindly in the ymap. This is consistent with the fact that we expect blind detection methods in the ymap to be less sensitive than multifrequency searches for compact sources (Melin et al. 2012). In some cases these nondetected sources are extended but with a relatively high S/N (between 6 and 12). We note that the MHW2 and MF methods are tuned to detect pointlike and compact sources primarily, and thus they may fail for extended ones.
We present in the top panel of Fig. 6 a comparison of the flux measured for these common sources in the NILC ymap with that derived from the PSZ2. The quantity Y_{5r500} corresponds to the integrated signal within a radius equal to 5 ×r_{500}^{3}. We observe good agreement between the two. Most of the observed outliers (low ymap flux with respect to the PSZ2 one) correspond to cluster candidates detected close to masked point sources, or in regions with strong Galactic contamination. We also show in the bottom panel a comparison of the flux measured with the NILC and MILCA ymaps. We find good agreement between the two, with three outliers corresponding to sources lying either in a region badly affected by radio and/or infrared point sources or for which we observe large uncertainties in the estimated characteristic radius of the cluster candidate. For both maps (NILC and MILCA), the average S/N of the recovered clusters with the MHW2 algorithm is around 10, while S/N of PSZ2 clusters only found with the matched filter technique is, as expected, lower (about 6).
Fig. 6 Top: comparison between the measured tSZ flux reported in the Planck cluster sample and that estimated directly on the NILC ymap for the blindlydetected common sources. Bottom: as above but for the NILC and MILCA measured fluxes. Outliers in the figures are discussed in the text. 
4.4.3. Maps of selected clusters
We have performed a more detailed analysis for some of the cluster candidates in the PSZ2 sample. In Fig. 7 we present Compton parameter maps centred at the positions of three of the newly discovered Planck clusters with S/N of 9.3 (left), 6.2 (middle) and 4.6 (right), both for MILCA (top row) and NILC (middle row). In the bottom row of the figure we present the radial profiles of these clusters as obtained from the MILCA (black lines) and NILC (blue) maps. We find that the profiles are consistent between MILCA and NILC. For comparison we also show the radial profile of a 10′ Gaussian beam. We observe that even for low S/N and compact clusters we are able to detect extended emission. Thus, these maps could be used to expand the pressure profile analysis presented in Planck Collaboration Int. V (2013) to fainter and higher redshift clusters.
One of the major outcomes of the NILC and MILCA ymaps is the possibility to study diffuse faint emission between clusters, as well as the emission from the cluster outskirts. We show in Fig. 8 some wellknown merging systems including the Shapley supercluster, and the A3395A3391 and the A339A401 interacting systems (see Planck Collaboration Int. VIII 2013, for a detailed description). We observe that the clusters themselves, their outskirts, and the intercluster emission are well reconstructed in the NILC and MILCA ymaps (which show consistent results). The quality of these maps will permit analyses similar to what is presented in Planck Collaboration Int. VIII (2013), with a significant increase of the S/N.
4.5. SZ signal below the noise level
We use the catalog of 132 684 clusters of galaxies identified from the Sloan Digital Sky Survey III (Wen et al. 2012) in order to quantify the tSZ signal below the noise level in the ymap. This catalogue provides estimates of cluster redshift, richness (N_{200}), and characteristic radius (r_{200}).
We focus on unresolved groups and clusters for which the richness is 8 ≤ N_{200} ≤ 100. These objects are expected to be at S/N well below 1 in the ymap and so their direct detection is not possible. However, their average tSZ signal can be detected using a stacking approach on the ymaps. Figure 9 shows the stacked maps (obtained with the IAS stacking library, Bavouzet 2008; Béthermin et al. 2010) for six richness intervals with N_{200} varying from 8 to 100, on patches of 4° × 4°, for the full sky NILC ymap. For the stacking we exclude all the positions for which point sources are found within a radius of 10′ from the centre of the object. The noise of the stacked maps scales, as expected, with , and the average signal increases for increasing richness. Then, when considering a sufficiently large number of objects, we are able to significantly detect the average SZ signal, even for small groups (Fig. 9). We obtain consistent results for the MILCA and NILC ymaps.
Fig. 7 Compton parameter MILCA (top row) and NILC (middle row) maps for a selection of PSZ2 cluster candidates with signal to noise ratio of 9.3, 6.2 and 4.6 from left to right. The maps are centred at the positions of the clusters in Galactic coordinates; (40°,75°), (25.1°,−78.6°) and (2.4°,69.6°), respectively. The colour scale is in units of y × 10^{6}. The bottom row presents the cluster radial profiles for MILCA (black) and NILC (blue). The beam profile is shown as a blue dashed line. 
Fig. 8 Compton parameter maps in Galactic coordinates of well known merging systems for MILCA (left column) and for NILC (right column). The colour scale is in units of y × 10^{6}. 
For the N_{200} intervals reported in Fig. 9 we have estimated the total stacked fluxes (Y_{5r500}) as the integrated signal within a radius 5 × r_{500}, with r_{500} = 0.7 ⟨ r_{200} ⟩. Here ⟨ r_{200} ⟩ is the mean of the r_{200} reported by Wen et al. (2012) for clusters belonging to each considered subsample. The average values and associated errors have been obtained using a bootstrap approach. For this, we have constructed and stacked cluster samples obtained by randomly replacing sources from the original sample, so that each of them contains a number of clusters equal to the initial one. By repeating the process several times, we have determined the statistical properties of the population being stacked. In Fig. 10 (left panel) we report Y_{5r500} as a function of richness (N_{200}) for all the objects with z ≤ 0.42. In fact higher redshift clusters present in this catalogue may be biased to lower richness because of incompleteness of member galaxies, as detailed in Wen et al. (2012). Following Planck Collaboration XII (2011), we have fitted a power law of the form (4)with E(z)^{2} = (H(z) /H_{0})^{2}, D_{A} the angular diameter distance of the cluster, and Y_{20} a normalization parameter.
By considering only N_{200} ≤ 60 (for which the number of clusters in each richness bin is >1000), we find a slope α = 1.85 ± 0.20, which is consistent with the scaling obtained in Planck Collaboration XII (2011). We have verified that the scaling obtained is insensitive to whether we limit our cluster sample to high Galactic latitude objects or not, and that a change in the integration radius only affects the normalization Y_{20}. Using Y_{500} = Y_{5r500}/ 1.79 (Planck Collaboration XXIX 2014), we obtain Y_{20} = (16.5 ± 1.2) × 10^{5} arcmin^{2}. However, for independent cluster data sets, different choices are made for the fainter magnitude of the member galaxies contributing to N_{200}. This may affect the richness associated with an object of given mass and complicate the comparison of the constant term Y_{20} for analyses that are based on different cluster data sets.
The massrichness scaling relation, (5)can then be used to also derive the Y_{500}–M_{tot} scaling (see Fig. 10, right panel), by fitting a power law of the form (6)Ford et al. (2014) compares different results for the scaling of the mass with richness, in the form of Eq. (5). In particular they discuss the results obtained by Wen et al. (2012), for a subsample of clusters with already known masses, and the results of Covone et al. (2014), from weak lensing mass measurements of 1176 clusters of the catalogue used here. Assuming β = 1.2 ± 0.1 and M_{0} = (1.1 ± 0.11) × 10^{14}M_{⊙} (Wen et al. 2012; Covone et al. 2014,and a 10% error) in Eq. (5), we compute via a Monte Carlo approach the average M_{200} and associated uncertainties for each Y_{5r500} bin. From this we find B = 1.92 ± 0.42 for the Y_{5r500}−M_{200} scaling (Eq. (6)). This is consistent with the selfsimilar value (5/3) and with what is expected for this mass range according to simulations (Sembolini et al. 2014). Through a weak lensing analysis of more than 18 000 clusters in the CanadaFranceHawaii Telescope Lensing Survey (CFHTLenS), Ford et al. (2014) obtain a steeper calibration for the massrichness relation, with β = 1.58 ± 0.34. For higher β, the Y_{5r500}−M_{200} slope B becomes lower, but still consistent with the value obtained before. Again comparison of M_{0} and Y_{0} is complicated by the different choices of the fainter limit on galaxies contributing to N_{200}. By exploring the range 0.2 ≲ z ≲ 0.9, Ford et al. (2014) finds no evolution with redshift. Thus, we do not expect our result to be biased by the z ≤ 0.42 cut.
5. Angular power spectrum
We now consider the validation of the Compton parameter maps at the power spectrum level.
5.1. Power spectrum estimation
To estimate the power spectrum of the tSZ signal we use the XSPECT method (Tristram et al. 2005), initially developed for the crosscorrelation of independent detector maps. XSPECT uses standard MASTERlike techniques (Hivon et al. 2002) to correct for the beam convolution and the pixelization, as well as the modecoupling induced by masking foreground contaminated sky regions.
Fig. 9 4° × 4° average maps for different ranges in N_{200} from 8 to 100. The colour scale is in units of y × 10^{6}. 
Fig. 10 Integrated tSZ signal as a function of cluster richness (left) and total mass (right). The black points correspond to the average signal obtained for the richness bins considered in Fig. 9. The red line represents the corresponding bestfit power law (Eqs. (4) and (6), respectively). Considering z ≤ 0.42, there are 13 814 objects for 8 <N_{200} ≤ 10, 37 250 for 10 <N_{200} ≤ 20, 7458 for 20 <N_{200} ≤ 30, 2069 for 30 <N_{200} ≤ 40, and 1133 for 40 <N_{200} ≤ 60. 
In the following, all the spectra will use a common multipole binning scheme, which was defined in order to minimize the correlation between adjacent bins at low multipoles and to increase the S/N at high multipole values. We will also only consider angular crosspower spectra between the ymaps reconstructed from the first (F) and second (L) halves of the data. This allows us to avoid the bias induced by the noise in the angular auto power spectrum, the correction of which would require a large number of noise simulations. The angular crossspectra are named in the following “NILC F/L” and “MILCA F/L”. Error bars in the spectrum are computed analytically from the autopower and crosspower spectra of the pairs of maps, as described in Tristram et al. (2005). We do not consider here higher order terms in the power spectrum variance. All of our Compton parameter maps assume a circular Gaussian beam of 10′ FWHM. The additional filtering at large angular scales in the MILCA Compton parameter maps is also accounted for and deconvolved when producing power spectra.
Fig. 11 Angular crosspower spectra of the Planck NILC F/L (left) and MILCA F/L (right) reconstructed Compton parameter maps for different Galactic masks corresponding to 30% (cyan), 40% (green), 50% (blue), and 60% (pink) of the sky. For comparison we also show MILCANILC F/L (red) and NILC70 F/L (black) on 50% of the sky. See text for details. 
5.2. Foreground contamination
In Planck Collaboration XXI (2014) we identified the dominant foreground contributions to the angular power spectrum of the reconstructed ymaps using the FFP6 simulated data set. We have repeated here the same analysis on updated simulated data sets (see Sect. 2.2), for which the description of foreground components has been improved according to our current knowledge (Planck Collaboration XII 2016). We find no major difference between the two analyses and therefore we do not repeat the discussion here. At large angular scales the dominant foreground contribution is the Galactic thermal dust emission, as we have already discussed in Sect. 4.3. At intermediate and small angular scales the major contribution comes from the clustered CIB emission. Equally important at small angular scales are the residual contributions from radio and IR sources.
5.2.1. Lowmultipole contribution
Assuming that at large angular scales the Compton parameter maps are mainly affected by diffuse Galactic dust emission, we have tested several Galactic masks by imposing flux cuts on the Planck 857 GHz channel intensity map. In particular we investigated masking out 30%, 40%, 50%, and 60% of the sky. The edges of these masks have been apodized to limit ringing effects on the reconstruction of the angular power spectrum. We stress that to compute power spectra we used the XSPECT algorithm (Tristram et al. 2005), which has been specifically developed to account for mode coupling. We have proved using simulations that we have no bias in the final power spectrum, and so we do not expect that any differences between MILCA and NILC could be induced by the different sizes of the masked area. Figure 11 shows the power spectra for NILC (left) and MILCA (right). We observe that the MILCA F/L large angular scale power decreases when imposing a more severe masking (larger fraction of the sky is masked out). A similar effect, but significantly smaller, is also observed for the MILCA F/L cross–power spectrum. The MILCA ymap is significantly more contaminated by thermal dust emission at multipoles below 30. However, we also observe for NILC F/L an extra noise contribution at large angular scales, as discussed in Sect. 4.2. These two problems limit the reliability of the results at multipoles below 30.
To extend the measurement of the angular power spectrum of the tSZ emission to multipoles below 30 we have considered two options: (1) as in Planck Collaboration XXI (2014) we apply a more severe Galactic mask (30% of the sky is masked) before the computation of the NILC weights to produce the ymap (NILC70); and (2) compute the crosscorrelation of NILC and MILCA ymaps. Considering 50% of the sky, we show in Fig. 11 the angular power spectrum of the NILC70 ymap (black) and the crossspectrum of the NILC first half and MILCA second half (NILCMILCA F/L) ymaps. We observe that the two are compatible within error bars in the multipole range 10 <ℓ< 1500. As expected the cross spectrum NILCMILCA F/L shows larger error bars.
Fig. 12 Comparison of the tSZ angular power spectrum estimated from the crosspowerspectrum of the NILC (left) and MILCA (right) F/L maps (black), with the expected angular power spectrum of the confirmed clusters in the Planck Cluster Sample (green solid line). The angular crosspower spectrum between the NILCA and MILCA Compton parameter maps and the simulated detected cluster map is shown in red. The correlation between the reconstructed ymap and the simulated detected cluster map, to which an arbitrary rotation has been applied, is plotted in grey (negative values dashed–lines for the error bars). 
Although the NILC70 ymap seems to be the best choice in terms of power spectrum estimation, it results in a significant reduction of the available sky area for other kind of studies as those presented in Sect. 4. Furthermore, it is difficult to accurately estimate the ultimate residual foreground contribution at very large angular scales. Because of this, and to preserve the coherence of the delivered products and the analysis presented in this paper, we have chosen the angular cross–power spectrum of NILCMILCA F/L as a baseline. This is obviously a more conservative choice in terms of noise–induced uncertainties. The NILCMILCA F/L angular cross–power spectrum bandpowers and uncertainties are further discussed in Sect. 7. Using the inscan and crossscan ymaps presented in Sect. 4.1 we find that stripe contamination accounts for less than 10% of the total signal in the NILCMILCA F/L cross angular power spectrum. We have enlarged the error bars to account for this systematic effect.
5.2.2. Highmultipole contribution
At small angular scales the measured tSZ power spectrum is affected by residual foreground contamination, coming from clustered CIB emission as well as radio and IR point sources. They show up in the MILCANILC F/L cross–power spectrum (see Fig. 11) as an excess of power at large multipoles.
To deal with those residuals we adopt the same strategy as in Planck Collaboration XXI (2014). We define physically motivated models of the angular power spectrum of the foreground components for each observation channel, including cross–correlations between channels. In contrast to Planck Collaboration XXI (2014) we also account for the cross correlation between the clustered CIB and the tSZ emission. A detailed description of this cross correlation, as well as of the clustered CIB model is presented in Planck Collaboration XXIII (2016). For the radio and IR point source models we refer to Planck Collaboration XXI (2014).
We use the outputs of Planck Collaboration XVIII (2011) and Planck Collaboration XXX (2014) for the clustered CIB modelling. For the six Planck HFI frequencies considered in this paper, the clustered CIB model consists of six autopower spectra and 24 crosspower spectra. For frequencies above 217 GHz, these spectra are fitted in Planck Collaboration XXX (2014) to the measured CIB, consistently with Planck Collaboration XVIII (2011). The model is extrapolated at 100 and 143 GHz following Béthermin et al. (2012) and Planck Collaboration XVIII (2011). The uncertainties in the clusteredCIB model are mainly due to the crosscorrelation coefficients that relate the crosspower spectra to the autopower spectra. Following Planck Collaboration XXX (2014) we consider 5% global uncertainties on those coefficients.
We use the Béthermin et al. (2012) model to compute the starforming dusty galaxy contribution. Finally, we use the Tucci et al. (2011) model, fitted to the Planck ERCSC (Planck Collaboration Int. VII 2013), for extragalactic radio sources. Notice that these models are also used for the study of the clustered CIB with Planck (Planck Collaboration XXX 2014).
Using the models described above we generate Gaussian realizations of the foreground contribution for each HFI frequency channel between 100 and 857 GHz. Note that the LFI channels are only used at large angular scales. We apply the MILCA or NILC weights to these simulated maps. From these simulations we find that the cross–correlation between the CIB and tSZ contribution can be neglected to first order with respect to the others and will therefore not be considered hereafter. Uncertainties on the parameters describing the foreground models were also propagated using simulations. We find that the clustered CIB model uncertainties might be as large as 50% in amplitude. In addition, we notice that the amplitude of the point source models can vary significantly depending on the point source mask applied. These uncertainties are taken into account hereafter. The amplitude of the residual foreground models are jointly fitted with the cosmologydependent tSZ model, as detailed in Sect. 7.1.
5.3. Contribution of resolved clusters to the tSZ power spectrum
We simulate the expected Compton parameter map for the detected and confirmed clusters of galaxies in the Planck catalogue (Planck Collaboration XXIX 2014) from their measured integrated Compton parameter Y_{5r500} and redshift z. We assume hydrostatic equilibrium and an Arnaud et al. (2010) pressure profile. The green solid line in Fig. 12 shows the power spectrum of this simulated map. Figure 12 also shows the crosspower spectra of the NILC and MILCA F/L maps (in black). We also compute the crosspower spectrum of the simulated cluster map and the Planck reconstructed NILC and MILCA Compton parameter maps. This is shown in red in the figure. Here again, the signal is consistent with the expected power spectrum of the confirmed Planck clusters of galaxies. These results show that the tSZ signal from clusters is preserved in the ymaps.
6. Higher–order statistics
The power spectrum analysis presented above only provides information on the 2point statistics of the Compton parameter distribution over the sky. An extended characterization of the field can be performed by studying the higherorder moments in the 1D PDF of the map, or by measuring 3point statistics, i.e., the bispectrum.
Fig. 13 1D PDF of the Planckymap before (red) and after (orange) masking the PSZ2 clusters, and of the halfdifference map (black) for the NILC (left) and MILCA (right) methods. We also show for comparison the 1D PDF of the simulated PSZ2 cluster map (dark blue). 
6.1. 1D PDF analysis
Following Planck Collaboration XXI (2014), we perform an analysis of the 1D PDF of the NILC and MILCA reconstructed Compton parameter maps. For the tSZ effect we expect an asymmetric distribution with a significantly positive tail (RubiñoMartín & Sunyaev 2003). We thus focus on the asymmetry of the distribution and its unnormalized skewness. First, we filter the maps in order to enhance the tSZ signal with respect to foreground contamination and noise. We then apodize the combined PSMASK and Galactic mask to avoid residual point source ringing. We follow the approach of Wilson et al. (2012) and use a filter in harmonic space. For each multipole ℓ, the filter is constructed as the ratio between the angular power spectrum of the expected tSZ signal (, obtained from the simulations in Sect. 2.2) and the power spectrum of the halfdifference y maps (, see Sect. 4.2), such that . We smooth this filter function in multipole space using a 21point square kernel and normalize it to one. Notice that this filter only selects the multipole range for which the tSZ signal is large with respect to the noise, and thus, it does not modify significantly the nonGaussianity properties of the ymaps. Furthermore, we have found that the filter used here behaves better than the more traditionally used Wiener filter, because it is less affected by pointsource ringing. Following this procedure, the 1D PDF of the filtered Compton parameter maps, P(y), is computed from the histogram of the pixels. As for the power spectrum, several Galactic masks have been considered in order to tests the robustness of the results.
Figure 13 shows the 1D PDF for the NILC (left) and MILCA (right) Compton parameter map in red when masking 50% of the sky. They correspond to the convolution of the 1D PDF of the different components in the map: the tSZ effect; foregrounds; and the noise. The 1D PDFs of the NILC and MILCA ymaps clearly show three distinct contributions: a Gaussian central part that is slightly wider than the noise contribution, as expected from the halfdifference map 1D PDF (black curve); a small negative tail, corresponding most likely to residual radio sources; and a positive tail corresponding mainly to the tSZ signal as observed for the PSZ2 cluster simulated map (dark blue). In comparison to Planck Collaboration XXI (2014) we now find a better agreement between the NILC and MILCA results. This is mainly due to the improved masking of radio sources presented in Sect. 4.4.1. Finally we also show the 1D PDF of the reconstructed ymaps after masking the PSZ2 clusters (orange). We observe that most of the tSZ is removed, indicating that the 1D PDF of the ymap is dominated by detected clusters.
A simple analysis of the measured 1D PDF can be performed by considering the asymmetry of the distribution: (7)where y_{p} is the peak value of the normalized PDF P(y). In addition, the nonGaussianity of the positive tail can be quantified by (8)with G(y) the expected PDF if fluctuations were only due to noise. The latter can be obtained from the halfdifference ymaps. See Hill et al. (2014) for a similar analysis conducted on ACT data.
Masking 60% of the sky, we find A = 0.218 and Δ = 0.10 for the NILC Compton parameter map. Equivalently, for the MILCA Compton parameter map we find A = 0.223 and Δ = 0.11. Uncertainties on A and Δ are mainly given by the accuracy to which measure the Gaussian noise distribution and are below 1 %. These results are in agreement and consistent with a positive tail in the 1D PDF, confirming the tSZ nature of the signal. Actually, we find similar values for the FFP6 simulations, A = 0.30 and Δ = 0.13 as presented in Planck Collaboration XX (2014). The agreement between the NILC and MILCA results worsens when reducing the masked area, as a consequence of a larger foreground contribution in the MILCA ymap.
6.2. Bispectrum
Since the SZ signal is nonGaussian, significant statistical information is contained in the bispectrum, complementary to the power spectrum (RubiñoMartín & Sunyaev 2003; Bhattacharya et al. 2012). Results on SPT data have been obtained by Crawford et al. (2014) and George et al. (2015). The bispectrum also provides an additional statistic assessing the compatibility of the NILC and MILCA reconstructed Compton parameter maps, as well as their reliability in terms of foreground contamination.
We therefore analyse the bispectra of the NILC and MILCA maps. The estimation method is essentially the same as for the 2013 results (Planck Collaboration XX 2014), briefly recapped here. We use the binned bispectrum estimator described in Bucher et al. (2010) and Lacasa et al. (2012), which is also used for the Planck primordial nonGaussianity analysis (Planck Collaboration XXIV 2014; Planck Collaboration XVII 2016). We use a multipole bin size Δℓ = 64 and a maximum multipole ℓ_{max} = 2048 for the analysis, working at a resolution N_{side} = 1024 to reduce computing time. The NILC and MILCA maps are masked with the combination of PSMASK, described in Sect. 4.4.1, and a Galactic mask at 50, 60, or 70%, described in Sect. 5.2.1 (in the rest of this section we will simply denote as X% mask the combination of PSMASK and the Galactic mask at X%). The bestfit monopole and dipole outside the mask are finally removed before estimation.
An important part of the pipeline is then to correct for the bias introduced by masking. To this end, we compute the ratio of the fullsky and masked sky bispectra, on highly nonGaussian simulations with a tSZlike bispectrum and 10′ resolution. This ratio is used to correct the measured bispectra and flag unreliable (ℓ_{1},ℓ_{2},ℓ_{3}) configurations. Specifically we flag configurations where the ratio is different by more than 25% from the naive expectation f_{sky}B(ℓ_{1}) B(ℓ_{2}) B(ℓ_{3}), where B(ℓ) is the Gaussian 10′ beam.
For both NILC and MILCA, we find that the bispectra computed on the 50, 60 and 70% masks are consistent. This indicates that there is no detectable residual galactic contamination in these bispectra. However we did find Galactic contamination on less aggressive Galactic masks (when we mask less than 50% of the sky), specifically positive Galactic dust. As Galactic dust is highly nonGaussian, we do not recommend the measurement of higher order statistics using Galactic masks smaller than 50%. In the following we adopt the 50% mask as a baseline, since it leaves the most sky available for estimation and minimizes masking effects in the measurement.
Figure 14 shows the obtained bispectra as a function of multipole for the NILC and MILCA Compton parameter maps. We observe good agreement between the bispectra of the two maps, and the bispectral behaviour is consistent with that expected from a tSZ signal (see, e.g., Lacasa 2014, chapter 5). We further compare these measurements with the bispectrum of the simulated map for the PSZ2 clusters, which is presented in blue in Fig. 14. We observe good agreement between the bispectra of the NILC and MILCA and that of the PSZ2 clusters. We therefore conclude that the observed bispectrum in the ymap is dominated by detected clusters.
Fig. 14 Bispectra of the NILC (green) and MILCA (red) ymaps for four different configurations (equilateral, orthogonal, flat and squeezed), compared with the bispectrum of the projected map of the PSZ2 clusters (blue). ± 1σ uncertainties are indicated as black dotted lines. 
Finally, in Fig. 14 are shown the ± 1σ uncertainties of the measurements, in black dotted lines. The error bars were computed in a similar manner to that of the 2013 results (Planck Collaboration XX 2014), and Appendix A.3 gives more detailed discussion.
With a detection per configuration at an average significance of 3.5σ, and a total significance of around 60σ, the Planck data thus provide a high quality measurement of the nonGaussianity of the thermal SunyaevZel’dovich signal, with contamination from foregrounds below the uncertainty levels when masking more than 50 % of the sky.
7. Cluster physics and cosmology
7.1. Power spectrum analysis
7.1.1. tSZ power spectrum modelling
As a measure of structure growth, the tSZ power spectrum can provide independent constraints on cosmological parameters. As shown by Komatsu & Seljak (2002), the power spectrum of the tSZ effect is highly sensitive to the normalization of the matter power spectrum, commonly parameterized by the rms of the z = 0 mass distribution on 8 h^{1} Mpc scales, σ_{8}, and to the total amount of matter Ω_{m}. We expect the tSZ power spectrum to also be sensitive to other cosmological parameters, e.g., the baryon density parameter Ω_{b}, the Hubble contant H_{0}, and the primordial special index n_{s}. For reasonable external priors on those parameters, however, the variations are expected to be negligible with respect to those introduced by changes in Ω_{m} and σ_{8} and so they are not considered here.
Following Planck Collaboration XXI (2014) we consider here a twohalo model for the tSZ power spectrum, which is fully described in Appendix A.1. This model accounts for both intrahalo (1halo term) and interhalo (twohalos) correlations. Following Eq. (A.6), the tSZ spectrum is computed using the 2halo model, the Tinker et al. (2008) mass function, and the Arnaud et al. (2010) pressure profile. In particular, we use the numerical implementation presented in Taburet et al. (2009, 2010, 2011), and integrate in redshift from 0 to 3 and in mass from 10^{13}M_{⊙} to 5 × 10^{15}M_{⊙}. Our model allows us to compute the tSZ power spectrum at the largest angular scales and is consistent with the tSZ spectrum presented in Efstathiou & Migliaccio 2012 (EM12), which was used as a template in the CMB cosmological analysis in Planck Collaboration XVI (2014) and Planck Collaboration XI (2016). We also include the mass bias parameter, b, which accounts for bias between the observationally deduced (M^{obs}) and true (M^{true}) mass of the cluster (see Planck Collaboration XX 2014, for details), such that M^{obs} = (1−b)M^{true}.
Cosmological parameter results are very sensitive to the mass bias and in particular we expect σ_{8} and Ω_{m} to be strongly degenerate with b (Planck Collaboration XXI 2014). By contrast to Planck Collaboration XXIV (2016), for which external priors in the mass bias have been used, we consider here only two distinct values, b = 0.2 and b = 0.4. The former corresponds to the average value that numerical simulations seem to indicate (see Fig. A.2 in Planck Collaboration XX 2014). On the other hand, the latter value for the bias alleviates the inconsistency with the constraints derived from the analysis based on primary CMB anisotropies (see Planck Collaboration XX 2014).
This value is, however, larger than that obtained by mass comparison of clusters present in both the Planck cosmology sample (Planck Collaboration XX 2014) and the Weighing the Giants (WtG von der Linden et al. 2014) project. Even if studies based on lensing mass measurements still provide different and inconsistent results for the cluster mass calibration, their number and their accuracy has improved dramatically in the recent past (e.g., Mahdavi et al. 2013, for the CCCP project, Umetsu et al. 2014, for CLASH; Israel et al. 2014, for 400d WL; Ford et al. 2014 for CFHTLenS; Gruen et al. 2014 for WISCy) and they are expected to provide more useful constraints in the near future.
7.1.2. Maximum likelihood analysis
As in Planck Collaboration XXI (2014), cosmological constraints are obtained from a fit to the tSZ power spectrum. Following the discussion of Sect. 5.2 we take the NILCMILCA F/L crosspower spectrum (black dots in Fig. 15) as a reference and limit the analysis to 50% of the sky to minimize foreground residuals. In terms of astrophysical signal we consider a fourcomponent model: tSZ; clustered CIB; radio point sources; and infrared point sources. We restrict the analysis to multipoles ℓ> 10 so that we can neglect the residual thermal dust contamination (see Sect. 5.2.1). For ℓ> 2000 the total signal in the tSZ map is dominated by correlated noise, which is also accounted for in the fit. Because of this correlated noise and the expected high value of foreground contamination at the smallest scales we limit the fit to multipoles ℓ< 1411. Finally, the observed ymap power spectrum, , is modelled as (9)Here is the tSZ power spectrum (in red in Fig. 15), is the clustered CIB power spectrum (in green), and and are the infrared and radio source power spectra, respectively (in cyan and in blue). is an empirical model for the high multipole correlated noise.
Fig. 15 NILC  MILCA F/L crosspower spectrum (black) compared to the power spectra of the physically motivated foreground models. The foregrounds considered are: clustered CIB (green line); infrared sources (cyan line); and radio sources (blue line). Additionally, the bestfit tSZ power spectrum model presented in Sect. 7.1 is also plotted as a solid red line. 
Foreground contamination is modelled following Sect. 5.2.2. As discussed there, the main uncertainties in the residual power spectrum translate into up to 50% uncertainty in the amplitude of the clustered CIB. We have not considered in this analysis the CIBtSZ crosscorrelation that was shown to be negligible in terms of the power spectrum. The amplitude of the IR and radio pointsource contribution will depend on the exact Galactic mask used for the analysis. However, we expect the shape of the their power spectra to remain the same. We thus allow for a variation of the normalization amplitudes for the clustered CIB, A_{CIB}, and for the point sources, A_{IR} and A_{rad}.
We assume a Gaussian approximation for the likelihood function. Bestfit values and uncertainties are obtained using an adapted version of the CosmoMC algorithm (Lewis & Bridle 2002). Only σ_{8} and Ω_{m} are allowed to vary. All other cosmological parameters are fixed to their bestfit values, as obtained in Table 2 of Planck Collaboration XVI (2014). The normalization amplitudes, A_{CIB}, A_{rad}, and A_{IR}, considered as nuisance parameters, are allowed to vary between 0 and 3. For the range of multipoles considered here, the tSZ angular power spectrum varies approximately as . The results are thus presented in terms of this parameter combination.
7.2. Bestfit parameters and tSZ power spectrum
Figure 16 presents the 2D and 1D likelihood distributions for the cosmological parameter combination σ_{8}(Ω_{m}/ 0.28)^{3 / 8} and for the foreground nuisance parameters. We present the results obtained assuming a mass bias of 0.2 (black) and 0.4 (red). We obtain very similar values for the nuisance parameters in both cases. In particular the bestfit values for a mass bias of 0.2 are , , and . However, there is a significant shift in the value of σ_{8}(Ω_{m}/ 0.28)^{3 / 8} as one would expect (Planck Collaboration XX 2014). In the case of a mass bias of 0.2 we have , while for a mass bias of 0.4 we have . Notice that these values are obtained in a specific framework, with all other cosmological parameters being fixed and a fiducial fixed model used for the signals. Relaxing this framework would likely weaken the constraints presented here, as discussed below.
Fig. 16 2D and 1D likelihood distributions for the combination of cosmological parameters σ_{8}(Ω_{m}/ 0.28)^{3 / 8}, and for the foreground parameters A_{Rad.PS}, A_{CIB} and A_{IR.PS}. We show the 68.3% and 95.4% contours. The red and black contours correspond to a fixed mass bias of 0.2 and 0.4, respectively. 
Marginalized bandpowers of the angular power spectrum of the Planck tSZ Compton parameter map (in dimensionless (ΔT/T)^{2} units), statistical and foreground errors, and bestfit tSZ power spectrum and number count models (also dimensionless).
Figure 15 shows the NILCMILCA F/L angular crosspower spectrum before correcting (black dots) for foreground contribution. We also show the bestfit foreground models: clustered CIB (green line); and radio (blue line) and IR (cyan line) point sources. The statistical (thick line) and total (statistical plus foreground, thin line) are also shown. The bestfit tSZ power spectrum is presented as a solid red line. We conclude that the NILCMILCA F/L angular crosspower spectrum is dominated by tSZ for multipoles ℓ< 700, and by foreground contribution for multipoles ℓ> 1200. We also note that for the bestfit model the radio pointsources contribution seems to be negligible with respect to the IR one. This is not a realistic result and it is most probably explained by the strong degeneracy observed between the radio and IR pointsource amplitude (see Fig. 16).
Finally we present in Fig. 17 the NILCMILCA F/L angular crosspower spectrum after correcting for foreground contributions. Uncertainties account for statistical and systematic errors, as well as for uncertainties in the foreground subtraction. The marginalized bandpowers and uncertainties are also presented in Table 2. We note that foregroundinduced uncertainties dominate at multipoles ℓ> 100. Bandpowers for the bestfit model for the angular tSZ power spectrum are also given for comparison. We also show in Fig. 17 th tSZ power spectrum estimates at high multipoles obtained in CMB oriented analyses by the Atacama Cosmology Telescope (ACT) and the South Pole Telescope (SPT, George et al. 2015). The black line shows the tSZ power spectrum template (EM12, Efstathiou & Migliaccio 2012) used in the Planck CMB cosmological analysis (Planck Collaboration XVI 2014; Planck Collaboration XI 2016) assuming the bestfit amplitude A_{tSZ} and the grey region ±2σ uncertainties from Planck Collaboration XI (2016). A direct comparison of our results to the estimates of the SZ power spectrum amplitude by ACT, SPT and Planck Collaboration XI (2016) is difficult as we need on the one hand to account for foreground contamination in the Planck SZ data and on the other hand to marginalise over all uncertainties on the SZ power spectrum amplitude. Accounting for foreground uncertainties as in Table 2 we find that χ^{2} of the bestfit value is not modified significantly by including the ACT and SPT data.
7.2.1. Dependence on cluster physics
As discussed in Planck Collaboration XXI (2014), we also expect the tSZ power spectrum amplitude to be sensitive to the physics of clusters of galaxies. To explore this dependence we have considered a set of predicted tSZ spectra for various physical models. In Fig. 18 we compare these models to the foregroundcleaned Planck tSZ power spectrum derived above, as well as to the Atacama Cosmology Telescope (ACT) and the South Pole Telescope (SPT, George et al. 2015) power spectrum estimates. We consider the predictions derived from hydrodynamical simulations (Battaglia et al. 2010, 2012), from Nbody simulations plus semianalytical models (Trac et al. 2011, TBO2) and from analytical calculations (Shaw et al. 2010). These models were originally computed for the set of cosmological parameters in Hinshaw et al. (2012) with σ_{8} = 0.8 and have been rescaled in amplitude to our bestfit value for . We note that there is some dispersion in the predicted amplitudes and shapes of the tSZ power spectrum. These differences reflect the range of methodologies and assumptions used both in the physical properties of clusters and in the technical details of the computation. The latter includes differences in the redshift ranges and also in the mass intervals probed by the limited sizes of the simulation boxes of the hydrodynamical simulations. Analytical predictions are also sensitive to the model ingredients, such as the mass function, mass bias and scaling relations adopted.
Fig. 17 NILC–MILCA F/L crosspower spectrum after foreground subtraction (red points), compared to the Atacama Cosmology Telescope (ACT; cyan dot) and the South Pole Telescope (SPT; orange, George et al. 2015) power spectrum estimates. The black line shows the tSZ power spectrum template (EM12, Efstathiou & Migliaccio 2012) used in the Planck CMB cosmological analysis (Planck Collaboration XVI 2014; Planck Collaboration XI 2016) with its best fit amplitude A_{tSZ} (Planck Collaboration XI 2016). The grey region allows comparison with the ±2σ interval. 
We see from Fig. 18 that the models presented above (the tSZ template for CMB analyses, plus the Battaglia et al. 2012; Shaw et al. 2010; and TBO2 models) provide reasonable fits to the data for multipoles above 200. For lower multipoles the Shaw et al. 2010 and TBO2 models are not consistent with the data.
We have also performed a simplified likelihood analysis to evaluate the uncertainties on cosmological parameters induced by the uncertainties in the modelling of the cluster physics. We replace our own model of the tSZ power spectrum by the models discussed above and recompute σ_{8}(Ω_{m}/ 0.28)^{3 / 8}, A_{CIB}, A_{rad} and A_{IR} from a simple linear fit to the NILCMILCA F/L crosspower spectrum. In the case of a mass bias of 0.2, we obtain values for σ_{8}(Ω_{m}/ 0.28)^{3 / 8} between 0.77 and 0.80, which lie within the 1σ uncertainties (0.03) presented above.
In the case of our fiducial model (see Appendix A.1), we can also consider uncertainties in the parameters describing the scaling relations allowing us to relate the observed tSZ flux to the mass of the cluster for a given redshift. Following Eq. (7) in Planck Collaboration XXVIII (2014) the main parameters to be considered are the mass bias b, the overall amplitude Y_{∗} and the scaling slope β. As discussed above, the mass bias is fully degenerate with σ_{8}. Similar conclusions can be drawn for Y_{∗}, which is expected to be known at the percent level (see Table 1 in Planck Collaboration XXVIII 2014) and therefore it is subdominant with respect to the uncertainties in the mass bias. Although the uncertainties in the slope of the scaling relation are relatively large, we have checked that they lead to negligible uncertainties on cosmological parameters.
Fig. 18 tSZ power spectrum for existing models in the literature. The NILCMILCA F/L crosspower spectrum after foreground correction (black dots) is compared to the Atacama Cosmology Telescope (ACT; cyan dot) and the South Pole Telescope (SPT; orange, George et al. 2015) power spectrum estimates. We also show the tSZ power spectrum models from hydrodynamic simulations (Battaglia et al. 2012, blue), from Nbody simulations plus semianalytical dust and gas models (Trac et al. 2011, cyan; TBO2), and from analytical calculations (Shaw et al. 2010, green). 
7.3. Higherorder statistics
7.3.1. Skewness measurements
The skewness of the 1D PDF distribution, / can also be used to derive constraints on cosmological parameters. Following Wilson et al. (2012) and Planck Collaboration XXI (2014) we have chosen a hybrid approach, by computing the skewness of the filtered Compton parameter maps outside the 50% sky mask. In particular, we have computed the skewness of the Planck data Compton parameter maps ⟨ y^{3} ⟩, and of the halfdifference maps .
Using the models presented in Appendix A we can show that the unnormalized skewness of the tSZ fluctuation, scales approximately as , whereas the amplitude of the bispectrum scales as with α = 11–12, as shown by Bhattacharya et al. (2012). In the following we do not consider the dependency of the bispectrum and the unnormalized skewness on other cosmological parameters, since such dependencies are expected to be significantly lower than for σ_{8} (Bhattacharya et al. 2012).
We derive constraints on σ_{8} by comparing the measured unnormalized skewness and bispectrum amplitudes with those obtained from simulations of the tSZ effect. The tSZ contribution was obtained from a hybrid simulation including a hydrodynamic component for z< 0.3, plus extra individual clusters at z> 0.3, and with σ_{8} = 0.789. This approach is strongly limited by systematic uncertainties and the details of the theoretical modelling (see Hill & Sherwin 2013). Uncertainties due to foreground contamination are computed using the simulations and are accounted for in the final error bars.
We obtain σ_{8} = 0.77 for NILC and σ_{8} = 0.78 for MILCA. Combining the two results and considering model and foreground uncertainties, we obtain σ_{8} = 0.78 ± 0.02 (68% C.L.). Notice that the reported uncertainties are mainly dominated by foreground contamination. However the model uncertainties only account for the expected dependence of the unnormalized skewness upon σ_{8}, as shown in Appendix A. As was also the case in Wilson et al. (2012), we have neglected the dependence on other cosmological parameters. We have also not considered any uncertainties coming from the combination of the hydrodynamical and individual cluster simulations. Because of these limitations, our error bars might be underestimated.
7.3.2. Fit of the 1D PDF distribution
We also derive constraints on σ_{8} by fitting the 1D PDF obtained in Sect. 6.1. Here, we follow the formalism described in Hill et al. (2014) that evaluates the tSZ 1D PDF theoretically integrating across individual cluster contributions. We specifically use the Tinker et al. (2008) mass function and the Arnaud et al. (2010) pressure profile. The latter is normalized following the Y–M scaling relations described in Planck Collaboration XX (2014), and considering a mass bias parameter of b = 0.2. All cosmological parameters are fixed to the Planck 2015 CMB analysis bestfit values, and we fit only for σ_{8}.
We note that the Hill et al. (2014) formalism explicitly neglects any effects due to overlapping clusters along the line of sight. For this reason, and given the uncertainties in the modelling of the foreground residuals, the bestfit solution to the observed 1D PDF is fitted in the region that is dominated by nonoverlapping clusters and where the noise and foregrounds contributions are minimal. In our case, we use y> 4.5 × 10^{6}.
The confidence limits on σ_{8} are obtained from a maximum likelihood approach, in which the likelihood has a multivariate Gaussian shape with a covariance matrix which only depends on σ_{8}. This covariance matrix is evaluated numerically, accounting only for Poisson terms (both for pixeltopixel covariance and due to the correlations introduced by the cluster’s yprofile). We obtain σ_{8} = 0.77 ± 0.02 (68% C.L.) both for NILC and MILCA maps, including statistical and systematic uncertainties.
Fig. 19 Marginalized likelihood distribution for σ_{8}(Ω_{m}/ 0.28)^{3 / 8} for tSZ– and CMB–based analyses. We represent the tSZ power spectrum analysis results assuming a mass bias, b, of 0.2 (red) and 0.4 (orange), the cluster number count analysis results (green; Planck Collaboration XXIV 2016), and the combined Planck CMB and BAO analysis (Planck Collaboration XIII 2016) with (cyan) and without (blue) extra lensing constraints. 
7.4. Comparison to other Planck cosmological probes
We have shown in Fig. 17 that the amplitude of the tSZ power spectrum measured in this paper and from the Planck CMB analysis are in good agreement. However, the Planck 2013 results (Planck Collaboration XX 2014) have shown some tension between CMB and tSZ derived constraints on σ_{8} for wide range of experiments including Planck. Figure 19 shows the marginalized likelihood distribution for σ_{8}(Ω_{m}/ 0.28)^{3 / 8}, as obtained from the combined Planck CMB and BAO analysis (Planck Collaboration XIII 2016) with and without lensing constraints. We also presents the results obtained from the Planck 2015 cluster number count analysis assuming a mass bias of 0.2 (Planck Collaboration XXIV 2016), and for the tSZ power spectrum analysis in this paper assuming a mass bias of 0.2 and 0.4. We observe that assuming a mass bias of 0.2 the two tSZ analyses are in good agreement, but the tension with the CMB measurements remains. Agreement would be better if we adopted a larger mass bias that increases the value of σ_{8}(Ω_{m}/ 0.28)^{3 / 8}. Furthermore, we find that including lensing constraints leads to smaller values of σ_{8}(Ω_{m}/ 0.28)^{3 / 8}, which as a consequence, are in better agreement with the tSZ results. From this we can conclude that the Large Scale Structure (clusters and lensing) data would have a marginal preference for small values of the mass bias.
8. Summary and conclusions
Because of its wide frequency coverage, from 30 to 857 GHz, the Planck satellite mission is particularly well suited for the measurement of the thermal SunyaevZeldovich effect. Working with the Planck frequency channel maps from 30 to 857 GHz, we have reconstructed the tSZ signal over the full sky using tailored component separation methods.
We tested and validated the Planckymaps extensively and characterized them in terms of noise and foreground contamination. As expected the noise in the ymap is inhomogeneous and can be characterized by pixeldependent variance and homogeneous correlated Gaussian noise. Foreground contamination by thermal dust emission is found to be important at large angular scales. Additional foreground contamination is due to radio and IR point sources, for which a mask is provided. In terms of tSZ signal we find good agreement between the flux of blindly detected clusters in the ymap and that measured for clusters in the Planck cluster sample. Furthermore, we find that the sensitivity of the ymap is sufficient to detect faint and diffuse structures such as bridges between merging clusters. Moreover, we have proved, via a stacking analysis, that the very low signaltonoise regions in the ymap preserve the tSZ signal even for small galaxy groups (tens of galaxies).
After accounting for foreground contribution, mainly thermal dust emission at large angular scales, and clustered CIB and point sources at small angular scales, we have derived from the ymap the tSZ angular power spectrum in the multipole range from 9 <ℓ< 1411. The results are compatible with previous Planck measurements (Planck Collaboration XXI 2014) and extend significantly the multiple range giving for the first time access to the 2halo term contribution. The cosmological analysis of the tSZ power spectrum allows us to set constraints on cosmological parameters representing the matter content in the Universe, mainly σ_{8} and Ω_{m}. These constraints are consistent with those obtained from cluster number counts (Planck Collaboration XXIV 2016) and in soft tension with those derived from CMB analysis (Planck Collaboration XIII 2016).
The analysis of the nonGaussian properties of the ymap using the 1D PDF, the unnormalized skewness and the bispectrum of the map have confirmed the tSZ nature of the signal.
The Planckymaps and additional ancillary data (noise variance maps, foreground masks and ILC weights) are made available to the public for the Planck 2015 release (see Appendix C for details). These ymaps are expected to be useful in a wide range of astrophysical and cosmological analyses with clusters. For any of these analyses, and depending on the scientific goals, the inhomogeneous properties of the noise, and the systematics and foreground contamination should be taken into account in different ways, as described in this paper. Regions masked by the point source mask should never be used. In the case of pixelbased analyses, quality flags can be defined by combining the information from the variance map and the various foreground masks. For power spectrum, cross correlation, and higher order statistic analyses, we note that the ymaps present significant foreground contamination that needs to be taken into account both by masking highly contaminated regions (mainly the Galactic plane region) and by using adequate foreground models to which the ILC weights are applied. Taking these necessary precautions, the Planckymaps will prove a very useful tool for further studies.
Planck (http://www.esa.int/Planck) is a project of the European Space Agency (ESA) with instruments provided by two scientific consortia funded by ESA member states and led by Principal Investigators from France and Italy, telescope reflectors provided through a collaboration between ESA and a scientific consortium led and funded by Denmark, and additional contributions from NASA (USA).
Note that in Planck Collaboration XXI (2014) cosine window functions were used.
Acknowledgments
The Planck Collaboration acknowledges the support of: ESA; CNES and CNRS/INSUIN2P3INP (France); ASI, CNR, and INAF (Italy); NASA and DoE (USA); STFC and UKSA (UK); CSIC, MINECO, J.A., and RES (Spain); Tekes, AoF, and CSC (Finland); DLR and MPG (Germany); CSA (Canada); DTU Space (Denmark); SER/SSO (Switzerland); RCN (Norway); SFI (Ireland); FCT/MCTES (Portugal); 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/planckcollaboration.
References
 Arnaud, M., Pratt, G. W., Piffaretti, R., et al. 2010, A&A, 517, A92 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Battaglia, N., Bond, J. R., Pfrommer, C., Sievers, J. L., & Sijacki, D. 2010, ApJ, 725, 91 [NASA ADS] [CrossRef] [Google Scholar]
 Battaglia, N., Bond, J. R., Pfrommer, C., & Sievers, J. L. 2012, ApJ, 758, 75 [NASA ADS] [CrossRef] [Google Scholar]
 Bavouzet, N. 2008, Theses, Université Paris Sud, Paris XI [Google Scholar]
 Benson, B. A., de Haan, T., Dudley, J. P., et al. 2013, ApJ, 763, 147 [NASA ADS] [CrossRef] [Google Scholar]
 Béthermin, M., Dole, H., Beelen, A., & Aussel, H. 2010, A&A, 512, A78 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Béthermin, M., Daddi, E., Magdis, G., et al. 2012, ApJ, 757, L23 [NASA ADS] [CrossRef] [Google Scholar]
 Bhattacharya, S., Nagai, D., Shaw, L., Crawford, T., & Holder, G. P. 2012, ApJ, 760, 5 [NASA ADS] [CrossRef] [Google Scholar]
 Bleem, L. E., Stalder, B., de Haan, T., et al. 2015, ApJS, 216, 27 [NASA ADS] [CrossRef] [Google Scholar]
 Bucher, M., Tent, B. V., & Carvalho, C. S. 2010, MNRAS, 407, 2193 [NASA ADS] [CrossRef] [Google Scholar]
 Carlstrom, J. E., Holder, G. P., & Reese, E. D. 2002, ARA&A, 40, 643 [NASA ADS] [CrossRef] [Google Scholar]
 Covone, G., Sereno, M., Kilbinger, M., & Cardone, V. F. 2014, ApJ, 784, L25 [NASA ADS] [CrossRef] [Google Scholar]
 Crawford, T. M., Schaffer, K. K., Bhattacharya, S., et al. 2014, ApJ, 784, 143 [NASA ADS] [CrossRef] [Google Scholar]
 Das, S., Louis, T., Nolta, M. R., et al. 2013,1301, 1037 [Google Scholar]
 Delabrouille, J., Cardoso, J.F., Jeune, M. L., et al. 2009, A&A, 493, 835 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Delabrouille, J., Betoule, M., Melin, J.B., et al. 2013, A&A, 553, A96 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Diego, J. M., & Majumdar, S. 2004, MNRAS, 352, 993 [NASA ADS] [CrossRef] [Google Scholar]
 Dunkley, J., Calabrese, E., Sievers, J., et al. 2013, 1301, 776 [Google Scholar]
 Durret, F., Laganá, T. E., Haider, M. 2011, A&A, 529, A38 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Efstathiou, G., & Migliaccio, M. 2012, MNRAS, 423, 2492 [NASA ADS] [CrossRef] [Google Scholar]
 Ford, J., Hildebrandt, H., Van Waerbeke, L., et al. 2014, MNRAS, 439, 3755 [NASA ADS] [CrossRef] [Google Scholar]
 George, E. M., Reichardt, C. L., Aird, K. A., et al. 2015, ApJ, 799, 177 [NASA ADS] [CrossRef] [Google Scholar]
 GonzálezNuevo, J., Argüeso, F., LópezCaniego, M., et al. 2006, MNRAS, 369, 1603 [NASA ADS] [CrossRef] [Google Scholar]
 Górski, K. M., Hivon, E., Banday, A. J., et al. 2005, ApJ, 622, 759 [NASA ADS] [CrossRef] [Google Scholar]
 Gruen, D., Seitz, S., Brimioulle, F., et al. 2014, MNRAS, 442, 1507 [NASA ADS] [CrossRef] [Google Scholar]
 Hasselfield, M., Hilton, M., Marriage, T. A., et al. 2013, 1301, 816 [Google Scholar]
 Hill, J. C., & Sherwin, B. D. 2013, PRD, 87, 23527 [Google Scholar]
 Hill, J. C., & Spergel, D. N. 2014, J. Cosmol. Astropart. Phys., 2, 30 [Google Scholar]
 Hill, J. C., Sherwin, B. D., Smith, K. M., et al. 2014, ArXiv eprints [arXiv:1411.8004] [Google Scholar]
 Hinshaw, G., Larson, D., Komatsu, E., et al. 2012, 1212, 5226 [Google Scholar]
 Hivon, E., Górski, K. M., Netterfield, C. B., et al. 2002, ApJ, 567, 2 [NASA ADS] [CrossRef] [Google Scholar]
 Hou, Z., Reichardt, C. L., Story, K. T., et al. 2014, ApJ, 782, 74 [NASA ADS] [CrossRef] [Google Scholar]
 Hurier, G., MaciasPerez, J. F., & Hildebrandt, S. 2013, A&A, 558, A118 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Israel, H., Reiprich, T. H., Erben, T., et al. 2014, A&A, 564, A129 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Komatsu, E., & Kitayama, T. 1999, ApJ, 526, L1 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
 Komatsu, E., & Seljak, U. 2002, MNRAS, 336, 1256 [NASA ADS] [CrossRef] [Google Scholar]
 Lacasa, F. 2014, ArXiv eprints [arXiv:1406.0441] [Google Scholar]
 Lacasa, F., Aghanim, N., Kunz, M., & Frommert, M. 2012, MNRAS, 421, 1982 [NASA ADS] [CrossRef] [Google Scholar]
 Lewis, A., & Bridle, S. 2002, Phys. Rev. D, 66, 103511 [NASA ADS] [CrossRef] [Google Scholar]
 LópezCaniego, M., Herranz, D., GonzálezNuevo, J., et al. 2006, MNRAS, 370, 2047 [NASA ADS] [CrossRef] [Google Scholar]
 Ma, Y.Z., Van Waerbeke, L., Hinshaw, G., Hojjati, A., & Scott, D. 2015, JCAP, 09, 046 [NASA ADS] [CrossRef] [Google Scholar]
 Mahdavi, A., Hoekstra, H., Babul, A., et al. 2013, ApJ, 767, 116 [NASA ADS] [CrossRef] [Google Scholar]
 Mak, D. S. Y., & Pierpaoli, E. 2012, Phys. Rev. D, 86, 123520 [NASA ADS] [CrossRef] [Google Scholar]
 Melin, J., Bartlett, J. G., & Delabrouille, J. 2006, A&A, 459, 341 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Melin, J.B., Aghanim, N., Bartelmann, M., et al. 2012, A&A, 548, A51 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Mo, H. J., & White, S. D. M. 1996, MNRAS, 282, 347 [NASA ADS] [CrossRef] [Google Scholar]
 Planck Collaboration. 2013, The Explanatory Supplement to the Planck 2013 results, http://pla.esac.esa.int/pla/index.html (ESA) [Google Scholar]
 Planck Collaboration VIII. 2011, A&A, 536, A8 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Planck Collaboration X. 2011, A&A, 536, A10 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Planck Collaboration XI. 2011, A&A, 536, A11 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Planck Collaboration XII. 2011, A&A, 536, A12 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Planck Collaboration XVIII. 2011, A&A, 536, A18 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Planck Collaboration I. 2014, A&A, 571, A1 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Planck Collaboration IX. 2014, A&A, 571, A9 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Planck Collaboration XII. 2014, A&A, 571, A12 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Planck Collaboration XVI. 2014, A&A, 571, A16 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Planck Collaboration XX. 2014, A&A, 571, A20 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Planck Collaboration XXI. 2014, A&A, 571, A21 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Planck Collaboration XXIV. 2014, A&A, 571, A24 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Planck Collaboration XXVIII. 2014, A&A, 571, A28 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Planck Collaboration XXIX. 2014, A&A, 571, A29 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Planck Collaboration XXX. 2014, A&A, 571, A30 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Planck Collaboration I. 2016, A&A, 594, A1 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Planck Collaboration II. 2016, A&A, 594, A2 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Planck Collaboration III. 2016, A&A, 594, A3 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Planck Collaboration IV. 2016, A&A, 594, A4 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Planck Collaboration V. 2016, A&A, 594, A5 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Planck Collaboration VI. 2016, A&A, 594, A6 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Planck Collaboration VII. 2016, A&A, 594, A7 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Planck Collaboration VIII. 2016, A&A, 594, A8 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Planck Collaboration IX. 2016, A&A, 594, A9 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Planck Collaboration X. 2016, A&A, 594, A10 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Planck Collaboration XI. 2016, A&A, 594, A11 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Planck Collaboration XII. 2016, A&A, 594, A12 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Planck Collaboration XIII. 2016, A&A, 594, A13 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Planck Collaboration XIV. 2016, A&A, 594, A14 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Planck Collaboration XV. 2016, A&A, 594, A15 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Planck Collaboration XVI. 2016, A&A, 594, A16 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Planck Collaboration XVII. 2016, A&A, 594, A17 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Planck Collaboration XVIII. 2016, A&A, 594, A18 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Planck Collaboration XIX. 2016, A&A, 594, A19 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Planck Collaboration XX. 2016, A&A, 594, A20 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Planck Collaboration XXI. 2016, A&A, 594, A21 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Planck Collaboration XXII. 2016, A&A, 594, A22 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Planck Collaboration XXIII. 2016, A&A, 594, A23 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Planck Collaboration XXIV. 2016, A&A, 594, A24 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Planck Collaboration XXV. 2016, A&A, 594, A25 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Planck Collaboration XXVI. 2016, A&A, 594, A26 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Planck Collaboration XXVII. 2016, A&A, 594, A27 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Planck Collaboration XXVIII. 2016, A&A, 594, A28 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Planck Collaboration Int. V. 2013, A&A, 550, A131 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Planck Collaboration Int. VII. 2013, A&A, 550, A133 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Planck Collaboration Int. VIII. 2013, A&A, 550, A134 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Reichardt, C. L., Stalder, B., Bleem, L. E., et al. 2013, ApJ, 763, 127 [NASA ADS] [CrossRef] [Google Scholar]
 Remazeilles, M., Delabrouille, J., & Cardoso, J.F. 2011, MNRAS, 410, 2481 [NASA ADS] [CrossRef] [Google Scholar]
 Remazeilles, M., Aghanim, N., & Douspis, M. 2013, MNRAS, 430, 370 [NASA ADS] [CrossRef] [Google Scholar]
 RubiñoMartín, J. A., & Sunyaev, R. A. 2003, MNRAS, 344, 1155 [NASA ADS] [CrossRef] [Google Scholar]
 Sembolini, F., De Petris, M., Yepes, G., et al. 2014, MNRAS, 440, 3520 [NASA ADS] [CrossRef] [Google Scholar]
 Shaw, L. D., Nagai, D., Bhattacharya, S., & Lau, E. T. 2010, ApJ, 725, 1452 [NASA ADS] [CrossRef] [Google Scholar]
 Sunyaev, R. A., & Zeldovich, Y. B. 1972, Comments on Astrophysics and Space Physics, 4, 173 [Google Scholar]
 Taburet, N., Aghanim, N., Douspis, M., & Langer, M. 2009, MNRAS, 392, 1153 [NASA ADS] [CrossRef] [Google Scholar]
 Taburet, N., Douspis, M., & Aghanim, N. 2010, MNRAS, 404, 1197 [NASA ADS] [Google Scholar]
 Taburet, N., HernándezMonteagudo, C., Aghanim, N., Douspis, M., & Sunyaev, R. A. 2011, MNRAS, 418, 2207 [NASA ADS] [CrossRef] [Google Scholar]
 Tinker, J., Kravtsov, A. V., Klypin, A., et al. 2008, ApJ, 688, 709 [NASA ADS] [CrossRef] [Google Scholar]
 Trac, H., Bode, P., & Ostriker, J. P. 2011, ApJ, 727, 94 [NASA ADS] [CrossRef] [Google Scholar]
 Tristram, M., MacíasPérez, J. F., Renault, C., & Santos, D. 2005, MNRAS, 358, 833 [NASA ADS] [CrossRef] [Google Scholar]
 Tucci, M., Toffolatti, L., de Zotti, G., & MartínezGonzález, E. 2011, A&A, 533, A57 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Umetsu, K., Medezinski, E., Nonino, M., et al. 2014, ApJ, 795, 163 [NASA ADS] [CrossRef] [Google Scholar]
 Van Waerbeke, L., Hinshaw, G., & Murray, N. 2014, Phys. Rev. D, 89, 023508 [NASA ADS] [CrossRef] [Google Scholar]
 von der Linden, A., Allen, M. T., Applegate, D. E., et al. 2014, MNRAS, 439, 2 [NASA ADS] [CrossRef] [Google Scholar]
 Wen, Z. L., Han, J. L., & Liu, F. S. 2012, ApJS, 199, 34 [NASA ADS] [CrossRef] [Google Scholar]
 Wilson, M. J., Sherwin, B. D., Hill, J. C., et al. 2012, Phys. Rev. D, 86, 122005 [NASA ADS] [CrossRef] [Google Scholar]
Appendix A: Modelling the expected tSZ signal
Appendix A.1: tSZ power spectrum
The representation of the ymap in spherical harmonics, Y_{ℓm}, reads (A.1)Thus, its angular power spectrum is given by (A.2)Note that is a dimensionless quantity here, like y.
As in Planck Collaboration XX (2014), the tSZ power spectrum is modelled using a 2halo model, to account both for intrahalo and interhalo correlations: (A.3)Following Komatsu & Seljak (2002) the 1halo term reads (A.4)where dV_{c}/ (dzdΩ) is the comoving volume per unit redshift and solid angle and n(M,z)dM dV_{c}/ (dzdΩ) is the probability of having a galaxy cluster of mass M at a redshift z in the direction dΩ. The quantity is the 2D Fourier transform on the sphere of the 3D radial profile of the Compton yparameter of individual clusters, (A.5)where x = r/r_{s}, ℓ_{s} = D_{A}(z) /r_{s}, r_{s} is the scale radius of the 3D pressure profile, D_{A}(z) is the angular diameter distance to redshift z, and P_{e} is the electron pressure profile.
The 2halo term (Komatsu & Kitayama 1999; Diego & Majumdar 2004; Taburet et al. 2011) is given by (A.6)where P(k,z) is the 3D matter power spectrum at redshift z. The quantity B(M,z) is the timedependent linear bias factor that relates the matter power spectrum, P(k,z), to the power spectrum of the cluster correlation function. Following Komatsu & Kitayama (1999, see also Mo & White 1996 we adopt B(M,z) = 1 + (ν^{2}(M,z)−1) /δ_{c}(z), where ν(M,z) = δ_{c}(M) /D(z)σ(M), σ(M) is the presentday rms mass fluctuation, D(z) is the linear growth factor, and δ_{c}(z) is the threshold overdensity of spherical collapse.
Finally, we use the Tinker et al. (2008) mass function, dn(M,z) / dM, including an observedtotrue mass bias b, as discussed in detail in Planck Collaboration XX (2014), and we model the SZ Compton parameter using the pressure profile of Arnaud et al. (2010).
Appendix A.2: Nth moment of the tSZ field
Assuming a Poisson distribution (1halo term) of the clusters on the sky, and neglecting clustering between clusters the Nth moment of the tSZ signal (Komatsu & Kitayama 1999; Wilson et al. 2012; Planck Collaboration XX 2014) reads (A.7)where y(θ,M,z) is the integrated Compton parameter along the line of sight for a cluster of mass M at redshift z.
Appendix A.3: Bispectrum
The angular bispectrum is given by (A.8)where the angleaveraged quantity in the fullsky limit can be written as (A.9)and satisfies the conditions m_{1} + m_{2} + m_{3} = 0, ℓ_{1} + ℓ_{2} + ℓ_{3} = even, and ℓ_{i}−ℓ_{j} ≤ ℓ_{k} ≤ ℓ_{i} + ℓ_{j}, for the Wigner 3j function in brackets. Assuming a Poissonian spatial distribution of the clusters, as above, the bispectrum reads (Bhattacharya et al. 2012)
Appendix B: Bispectrum cosmic variance
Following Lacasa (2014, chapter 2), the bispectrum cosmic variance is composed of a Gaussian term, a bispectrum × bispectrum term, a spectrum × trispectrum term and a connected 6point term. Due to the lack of model or measurement of the trispectrum and 6point function, we neglected the last two terms. Note that they are, however, expected to yield a subdominant contribution. Thus we have two contributions (for full sky)

Gaussian cosmic variance: (B.1)where (B.2)and where C_{ℓ} is the auto power spectrum of the Compton parameter map, thus containing the noise contribution.

Bispectrum × bispectrum cosmic variance (B.3)This is the only term which gives offdiagonal contributions to the covariance matrix.
For our purpose, this cosmic variance is multiplied by f_{sky} and binned to the appropriate binning scheme.
We also consider systematic errors induced by foreground residuals or masking effects. We estimate systematic errors due to component separation uncertainties from the half difference of the NILC and MILCA bispectra. Masking effects are normally corrected for using simulations, which may under or overestimate leakage from large to small scales. We thus take a conservative ± 25% error on the debiasing ratio, consistent with the fact that the selected configurations have a ratio within ± 25 % of f_{sky}B(ℓ_{1}) B(ℓ_{2}) B(ℓ_{3}). This error is most likely a conservative overestimate.
Appendix C: Products
In the following we list the ymap related products delivered in the Planck 2015 data release^{4}:

fullsky MILCA and NILC ymaps for the full mission and for the first (F) and second (L) halves of Planck stable pointing period in Compton parameter units (see Sect. 3.2);

ILC weights per filter and per frequency used for the reconstruction of the MILCA and NILC full mission ymaps, as described in Sect. 3.2;

variance map accounting for the nonhomogeneous coverage and power spectrum of the correlated homogeneous counterpart, , for the MILCA and NILC full mission ymaps in Compton parameter units (see Sect. 4.2);

point source masks including known radio and IR sources, as described in Sect. 4.4.1;

Galactic masks used in the analyses presented in Sects. 5 and 6;
All Tables
Conversion factors for tSZ Compton parameter y to CMB temperature units and the FWHM of the beam of the Planck channel maps.
Marginalized bandpowers of the angular power spectrum of the Planck tSZ Compton parameter map (in dimensionless (ΔT/T)^{2} units), statistical and foreground errors, and bestfit tSZ power spectrum and number count models (also dimensionless).
All Figures
Fig. 1 Window functions corresponding to the spectral localization of the NILC and MILCA algorithms. For NILC (solid line) there are 10 Gaussian window functions defining 10 needlet scales. MILCA (dashed line) uses 11 overlapping Gaussian windows. 

In the text 
Fig. 2 Reconstructed Planck allsky Compton parameter maps for NILC (top) and MILCA (bottom) in orthographic projections. The apparent difference in contrast observed between the NILC and MILCA maps comes from differences in the residual foreground contamination and from the differences in the filtering applied for display purposes to the original Compton parameter maps. For the MILCA method, filtering out low multipoles significantly reduces the level of foreground emission in the final ymap. The wavelet basis used in the NILC method was tailored for tSZ extraction. For details see Planck Collaboration XXII (2016). 

In the text 
Fig. 3 A small region of the reconstructed Planck allsky Compton parameter maps for NILC (left) and MILCA (right) at intermediate Galactic latitudes in the southern sky centred at (0°,−45°) in Galactic coordinates. The colour scale is in units of y × 10^{6}. 

In the text 
Fig. 4 Inscan and crossscan contributions in the NILC (top line) and MILCA (bottom line) ymaps in Compton parameter units times 10^{6}. From left to right we present the original ymaps, and their in and cross scan contributions for a small region at intermediate Galactic latitudes in the southern sky centred at (0°,−45°) in Galactic coordinates. 

In the text 
Fig. 5 Top: standard deviation maps for the NILC (top) and MILCA (middle) ymaps corresponding to the inhomogenous noise contribution computed from the halfdifference of the halfring maps. Bottom: angular power spectrum of the homogenous noise contribution for the NILC (solid grey line) and MILCA (dotted black line) ymaps (see main text for details). 

In the text 
Fig. 6 Top: comparison between the measured tSZ flux reported in the Planck cluster sample and that estimated directly on the NILC ymap for the blindlydetected common sources. Bottom: as above but for the NILC and MILCA measured fluxes. Outliers in the figures are discussed in the text. 

In the text 
Fig. 7 Compton parameter MILCA (top row) and NILC (middle row) maps for a selection of PSZ2 cluster candidates with signal to noise ratio of 9.3, 6.2 and 4.6 from left to right. The maps are centred at the positions of the clusters in Galactic coordinates; (40°,75°), (25.1°,−78.6°) and (2.4°,69.6°), respectively. The colour scale is in units of y × 10^{6}. The bottom row presents the cluster radial profiles for MILCA (black) and NILC (blue). The beam profile is shown as a blue dashed line. 

In the text 
Fig. 8 Compton parameter maps in Galactic coordinates of well known merging systems for MILCA (left column) and for NILC (right column). The colour scale is in units of y × 10^{6}. 

In the text 
Fig. 9 4° × 4° average maps for different ranges in N_{200} from 8 to 100. The colour scale is in units of y × 10^{6}. 

In the text 
Fig. 10 Integrated tSZ signal as a function of cluster richness (left) and total mass (right). The black points correspond to the average signal obtained for the richness bins considered in Fig. 9. The red line represents the corresponding bestfit power law (Eqs. (4) and (6), respectively). Considering z ≤ 0.42, there are 13 814 objects for 8 <N_{200} ≤ 10, 37 250 for 10 <N_{200} ≤ 20, 7458 for 20 <N_{200} ≤ 30, 2069 for 30 <N_{200} ≤ 40, and 1133 for 40 <N_{200} ≤ 60. 

In the text 
Fig. 11 Angular crosspower spectra of the Planck NILC F/L (left) and MILCA F/L (right) reconstructed Compton parameter maps for different Galactic masks corresponding to 30% (cyan), 40% (green), 50% (blue), and 60% (pink) of the sky. For comparison we also show MILCANILC F/L (red) and NILC70 F/L (black) on 50% of the sky. See text for details. 

In the text 
Fig. 12 Comparison of the tSZ angular power spectrum estimated from the crosspowerspectrum of the NILC (left) and MILCA (right) F/L maps (black), with the expected angular power spectrum of the confirmed clusters in the Planck Cluster Sample (green solid line). The angular crosspower spectrum between the NILCA and MILCA Compton parameter maps and the simulated detected cluster map is shown in red. The correlation between the reconstructed ymap and the simulated detected cluster map, to which an arbitrary rotation has been applied, is plotted in grey (negative values dashed–lines for the error bars). 

In the text 
Fig. 13 1D PDF of the Planckymap before (red) and after (orange) masking the PSZ2 clusters, and of the halfdifference map (black) for the NILC (left) and MILCA (right) methods. We also show for comparison the 1D PDF of the simulated PSZ2 cluster map (dark blue). 

In the text 
Fig. 14 Bispectra of the NILC (green) and MILCA (red) ymaps for four different configurations (equilateral, orthogonal, flat and squeezed), compared with the bispectrum of the projected map of the PSZ2 clusters (blue). ± 1σ uncertainties are indicated as black dotted lines. 

In the text 
Fig. 15 NILC  MILCA F/L crosspower spectrum (black) compared to the power spectra of the physically motivated foreground models. The foregrounds considered are: clustered CIB (green line); infrared sources (cyan line); and radio sources (blue line). Additionally, the bestfit tSZ power spectrum model presented in Sect. 7.1 is also plotted as a solid red line. 

In the text 
Fig. 16 2D and 1D likelihood distributions for the combination of cosmological parameters σ_{8}(Ω_{m}/ 0.28)^{3 / 8}, and for the foreground parameters A_{Rad.PS}, A_{CIB} and A_{IR.PS}. We show the 68.3% and 95.4% contours. The red and black contours correspond to a fixed mass bias of 0.2 and 0.4, respectively. 

In the text 
Fig. 17 NILC–MILCA F/L crosspower spectrum after foreground subtraction (red points), compared to the Atacama Cosmology Telescope (ACT; cyan dot) and the South Pole Telescope (SPT; orange, George et al. 2015) power spectrum estimates. The black line shows the tSZ power spectrum template (EM12, Efstathiou & Migliaccio 2012) used in the Planck CMB cosmological analysis (Planck Collaboration XVI 2014; Planck Collaboration XI 2016) with its best fit amplitude A_{tSZ} (Planck Collaboration XI 2016). The grey region allows comparison with the ±2σ interval. 

In the text 
Fig. 18 tSZ power spectrum for existing models in the literature. The NILCMILCA F/L crosspower spectrum after foreground correction (black dots) is compared to the Atacama Cosmology Telescope (ACT; cyan dot) and the South Pole Telescope (SPT; orange, George et al. 2015) power spectrum estimates. We also show the tSZ power spectrum models from hydrodynamic simulations (Battaglia et al. 2012, blue), from Nbody simulations plus semianalytical dust and gas models (Trac et al. 2011, cyan; TBO2), and from analytical calculations (Shaw et al. 2010, green). 

In the text 
Fig. 19 Marginalized likelihood distribution for σ_{8}(Ω_{m}/ 0.28)^{3 / 8} for tSZ– and CMB–based analyses. We represent the tSZ power spectrum analysis results assuming a mass bias, b, of 0.2 (red) and 0.4 (orange), the cluster number count analysis results (green; Planck Collaboration XXIV 2016), and the combined Planck CMB and BAO analysis (Planck Collaboration XIII 2016) with (cyan) and without (blue) extra lensing constraints. 

In the text 
Current usage metrics show cumulative count of Article Views (fulltext article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.
Data correspond to usage on the plateform after 2015. The current usage metrics is available 4896 hours after online publication and is updated daily on week days.
Initial download of the metrics may take a while.