Issue 
A&A
Volume 677, September 2023



Article Number  A87  
Number of page(s)  16  
Section  Cosmology (including clusters of galaxies)  
DOI  https://doi.org/10.1051/00046361/202245834  
Published online  11 September 2023 
Thermal Sunyaev–Zeldovich measurements and cosmic infrared background leakage mitigation combining upcoming groundbased telescopes
^{1}
ArgelanderInstitut für Astronomie, Universität Bonn, Auf dem Hügel 71, 53121 Bonn, Germany
email: mcharmetant@astro.unibonn.de
^{2}
Deutsches Zentrum für Luft und Raumfahrt e.V. (DLR) Projektträger, JosephBeuysAllee 4, 53113 Bonn, Germany
Received:
31
December
2022
Accepted:
5
June
2023
Context. The Fred Young Submillimeter Telescope (FYST) and the Simons Observatory Large Aperture Telescope (SO LAT) will deliver unprecedented highresolution measurements of microwave sky emissions. Notably, one of those microwave sky emissions, the thermal Sunyaev–Zeldovich (tSZ) signal, is an essential probe for cluster astrophysics and cosmology. However, an obstacle to its measurement is contamination by the cosmic infrared background (CIB), especially at high frequencies.
Aims. Our goal is to assess the detection and purity of tSZ power spectrum measurements from these two telescopes. We demonstrate that FYST’s highfrequency coverage helps lower CIB contamination and improves signal detection.
Methods. We simulated the various components of the microwave sky at the frequencies, sensitivities, and beam sizes of the upcoming SO LAT and FYST telescopes using fullsky Hierarchical Equal Area isoLatitude Pixelisation (HEALPix) map templates from the Websky simulations and the Python Sky Model (PySM). We used a mapbased internal linear combination (ILC) and a constrained ILC (CILC) to extract the tSZ signal and compute residual noises to assess CIB contamination and signal recovery.
Results. We find that the CIB’s residual noise power spectrum in the ILCrecovered tSZ is lowered by ∼35% on average over the scales ℓ ∈ [500, 5000] when SO LAT and FYST are combined compared to when SO LAT is used alone. We find that when using CILC to deproject CIB, the combined abilities of SO LAT and FYST offer a large ℓ ∈ [1800, 3500] window in which the recovered tSZ power spectrum is not noise dominated.
Key words: cosmic background radiation / cosmology: observations / galaxies: clusters: general / methods: data analysis / submillimeter: diffuse background
© The Authors 2023
Open Access article, published by EDP Sciences, under the terms of the Creative Commons Attribution License (https://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.
This article is published in open access under the Subscribe to Open model. Subscribe to A&A to support open access publication.
1. Introduction
The cosmic microwave background (CMB) photons we observe have been influenced by their interaction with the baryonic matter as they travel through the Universe. These interactions modify the CMB photons and are called secondary anisotropies of the CMB. The study of these secondary anisotropies allows for the extraction of information about the largescale structure (LSS) of the Universe as well as about the primary anisotropies, the matter content of the Universe. One of the secondary anisotropies, the thermal Sunyaev–Zeldovich (tSZ) effect, is due to the inverse Compton scatting of the CMB photons by the hot electrons in the intracluster medium (ICM) of galaxy clusters. This effect is a crucial tool for cluster detection and was first used by the South Pole Telescope (SPT; Staniszewski et al. 2009) and the Atacama Cosmology Telescope (ACT) Collaborations (Menanteau et al. 2010) and later by the Planck Collaboration (Planck Collaboration XXVII 2016) to construct cluster catalogues. The tSZ effect can be used to measure cluster properties, such as temperature (Erler et al. 2018). It is also crucial for cosmology, as it has a strong dependence on matter clustering σ_{8} and matter density (Komatsu & Seljak 2002; Planck Collaboration XXII 2016; Bolliet 2018). By fitting the tSZ power spectrum extracted from data with models of the signal and its contaminants, assuming a Gaussian likelihood function, and applying a Markov chain Monte Carlo (MCMC) analysis with a fixed mass bias b, one can constrain these cosmological parameters (Planck Collaboration XXI 2014). This implies that a small improvement in the tSZ power spectra statistical significance or a reduction of the bias due to residual contamination from other signals will lead to a significant improvement in measurements of σ_{8} and Ω_{m}. Therefore, a clean, uncontaminated tSZ map is essential.
The cosmic infrared background (CIB) is the emission from stars and galaxies absorbed and reemitted by dust in the infrared domain. It is the main contaminant in tSZ reconstruction and is also a tracer of the matter distribution and the LSS (Hurier 2015). The CIB contamination, or leakage, originates from the spatial correlation between CIB and tSZ and is also due to the component separation techniques, such as internal linear combination (ILC; Bennett et al. 2003; Eriksen et al. 2004; Park et al. 2007; Kim et al. 2008; Delabrouille et al. 2009), not being able to perfectly separate the signal from the many other astrophysical emissions present in the sky. Previous studies on CIB leakage in the Planck tSZ map made by Planck Collaboration XXII (2016) showed that the amplitude of the leakage increases when the CIB originates from higher redshift populations. The dependence of the leakage on the redshift of the CIB population comes from the fact that the separation methods of components primarily focus on reducing the main contaminant, which at high frequencies is Galactic thermal dust, and is thus less efficient at reducing the highz CIB contaminants. This could lead to a bias in tSZ estimation for highz clusters, as mentioned in Vikram et al. (2017). Another study by Alonso et al. (2018) showed that the error bars on the tSZ profile are not much improved by the deprojection (i.e. nulling) of the CIB with more constrained component separation methods such as constrained internal linear combination (CILC; Remazeilles et al. 2011). Dust contamination can reduce the number of detected galaxy clusters by up to ∼9% within z ∈ [0.3, 0.8] (Melin et al. 2018). Similarly, Zubeldia et al. (2023) demonstrated that CIB contamination in the tSZ signal can lead to ∼20% fewer clusters being detected and proposed a publicly available code to produce CIBfree cluster catalogues.
In the past three decades, several space missions, such as the COsmic Background Explorer (COBE; Mather et al. 1994), the Wilkinson Microwave Anisotropy Probe (WMAP; Bennett et al. 2003), and Planck (Planck Collaboration I 2014), have observed the sky at microwave frequencies and extracted fullsky maps of the CMB, tSZ, CIB, and other astrophysical components. Upcoming groundbased experiments, such as those using the Simons Observatory Large Aperture Telescope (SO LAT; Ade et al. 2019) and the Fred Young Submillimeter Telescope (FYST; CCATprime Collaboration 2023), will allow for unprecedented measurements of the tSZ and CIB thanks to their high spatial resolution and high sensitivities.
In this work, our focus is on quantifying the benefits of combining future experiments, such as SO LAT and FYST to reduce the CIB leakage in ILCreconstructed tSZ maps. A cleaner tSZ map will allow better constraints on the cosmological parameters σ_{8} and Ω_{m} to be drawn, more accurate cluster counts and detection, and more accurate measurements of the properties of clusters, such as their temperature. Without mitigation, CIB leakage along the line of sight (LoS) could systematically bias measurements of the tSZ effect and lead to higher tSZ measurements, or it could lower such measurements by filling the characteristic tSZ decrement at frequencies below 217 GHz (Schaan et al. 2021). Obtaining less contaminated signals by combining microwave experiments is common and will be part of the observation strategies of future experiments. For example, previous studies, such as Ade et al. (2019), have shown that SO LAT will achieve precise measurement of the tSZ power spectrum when combined with Planck data. As Planck observations are not contaminated by the atmosphere, this combination allows SO LAT to mitigate such contamination on large scales. Similarly, Melin et al. (2021) constructed a bigger cluster catalogue by combining Planck and the South Pole Telescope (SPTSZ) data. Our work also encourages the combination of all available data from previous and upcoming CMB experiments.
This paper is structured as follows: In Sect. 2, we formalise the tSZ and CIB description in our study. In Sect. 3, we describe the modelling of highresolution maps of the different microwave sky emissions. In particular, we first explain how we generated Galactic foregrounds, the extragalactic components, and their processing. Then, we detail the generation of noises and the handling of beams and masks. This is followed by a conclusion of the output of the sky modelling. In Sect. 4, we summarise the methods used to extract the tSZ effect, starting with the ILC and followed by CILC. Next, we explain beam convolution handling and the treatment of fullsky maps, and we provide an overview of the methods used. We present a detailed description of our results and an overview in Sect. 5. Section 6 contains a summary and discussion of the abovementioned results.
2. Formalism
2.1. The thermal Sunyaev–Zeldovich effect
The thermal tSZ effect, predicted by Sunyaev & Zeldovich (1970, 1972), is the inverse Compton scattering of CMB photons on the hot electrons present in the ICM. A consequence of this effect is the distortion of the CMB blackbody spectrum, due to the shift of the low frequency (corresponding to low energies) CMB photons to higher frequencies (corresponding to higher energies). The temperature change of the CMB due to the tSZ effect can be expressed as:
where T_{CMB} = 2.72548 ± 0.00057 K (Fixsen 2009) is the temperature of the CMB and y is the dimensionless Comptony parameter measuring the electron pressure integrated over the LoS proper distance l,
where σ_{T} is the Thomson crosssection, m_{e} the electron mass, and P_{e}(l) = n_{e}(l)k_{B}T_{e}(l) the electron pressure composed of the electron density n_{e}, the Boltzmann constant k_{B}, and the electron temperature T_{e}. In Eq. (1), f(x_{ν}) is the spectral shape of the tSZ (Sunyaev & Zeldovich 1972, 1980):
where x_{ν} = hν/k_{B}T_{CMB}, h being the Planck constant and ν the frequency. The frequency spectrum of the tSZ is characterised by a decrement in the intensity change for frequencies lower than ν ∼ 217 GHz and an increment at higher frequencies. To express the frequency variations induced by the tSZ effect in terms of an intensity variation, we use the conversion factor to pass from units of thermodynamic temperature (K_{CMB}; Planck Collaboration IX 2014) in Eq. (1), to intensity units (Jy sr^{−1}),
where B_{ν}(T) is the Planck blackbody spectrum at a fixed frequency ν, I_{0} = 2(k_{B}T_{CMB})^{3}/(hc)^{2} ≈ 270 MJy sr^{−1}, c is the speed of light, and
Therefore, the tSZ intensity change is given by
More details on the Sunyaev–Zeldovich (SZ) effects can be found in the review by Mroczkowski et al. (2019). In this study, we place ourselves in a simplistic approach and do not take into account the relativistic correction of the tSZ effect (Wright 1979; Itoh et al. 1998; Pointecouteau et al. 1998; Chluba & Sunyaev 2012), even though Planck data have shown that this leads to an underestimation of the amplitude of its power spectrum (Remazeilles et al. 2018).
2.2. The cosmic infrared background
The CIB is the extragalactic integrated emission produced by all stars and galaxies during their formation and throughout their life. The CIB frequency spectrum peaks in the infrared regime because the majority of sources are far and their light is redshifted by the expansion of the Universe but also because some of this light is absorbed by dust and reemitted into the infrared regime. This existence of the CIB was predicted by Partridge & Peebles (1967) and first detected by Puget et al. (1996), Fixsen et al. (1998), Hauser et al. (1998). Contrary to the tSZ, one frequency spectrum is not sufficient to fully model it. The CIB spectral energy distribution (SED) depends on the redshift of sources. It can be approximated as a modified blackbody spectrum:
where A_{dust} is the amplitude of the dust emission, β is the dust emissivity spectral index, and T_{dust} is the temperature of the dust grains. The CIB is a tracer of the entire history of star formation in the Universe. One characteristic of the CIB of particular importance for our study is its nonzero correlation with the tSZ signal. Indeed, some of the unresolved dusty starforming galaxies constituting the CIB are found in galaxy clusters. Therefore, the spatial correlation between tSZ and CIB is nonzero, and it was measured in Reichardt et al. (2012, 2021), Dunkley et al. (2013), George et al. (2015), Choi et al. (2020b).
3. Modelling the microwave sky
We set up a pipeline called Skymodel to create mock highresolution maps of the microwave sky that would be observed in future experiments. Our pipeline is based on Hierarchical Equal Area isoLatitude Pixelisation (HEALPix) template maps of the microwave sky from existing simulations as well as the most recent observations of the microwave sky. The microwave sky components are usually separated into two categories: the emissions coming from the Galaxy, called Galactic foregrounds, and the extragalactic components. Our Skymodel deals differently with these two categories, using an existing Python package for the Galactic foregrounds and HEALPix maps from mock simulations for the extragalactic ones.
3.1. Galactic foregrounds
Galactic foregrounds are the various emissions originating from the Milky Way and are thus in the foreground of extragalactic emissions. We considered the following Galactic foregrounds: Galactic synchrotron; freefree emission; anomalous microwave emission (AME), which is supposedly produced by spinning dust (Dickinson et al. 2018); and thermal emission from Galactic dust grains.
All of the abovementioned components were generated using the Python Sky Model (PySM; Thorne et al. 2017). The PySM^{1} package is based on Galactic foreground template maps obtained from Planck Collaboration XXV (2016) that were extracted using the COMMANDER algorithm (Eriksen et al. 2004, 2008; Planck Collaboration IX 2015; Planck Collaboration IV 2020) on Planck 2015 raw data, data from WMAP, and 408 MHz allsky survey by Haslam et al. (1981, 1982), reprocessed by Remazeilles et al. (2015). The original data maps were at a HEALPix (Górski et al. 2005) pixel resolution of pix_{size} = 13.74′ (N_{side} = 256). The resolution of a HEALPix pixel is defined as the square root of the pixel area. The PySM upgrades them to 6.9′ (N_{side} = 512) and adds smallscale fluctuations generated randomly by extrapolating the power spectra of those maps (MivilleDeschênes et al. 2007). We generated each foreground component at the desired frequency using the basic spectral model of PySM (i.e. “s1”, “f1”, “d1”, and “a1”). These maps were then upgraded to a higher N_{side} = 4096 and smoothed with a 10′ full width at half maximum (FWHM) Gaussian to mitigate pixelation artefacts (see Fig. 1).
Fig. 1. Projections of the PySM sky at 353 GHz composed of freefree, synchrotron, AME, and Galactic dust. The upper row shows cutouts around the Galactic south pole, and the lower row shows cutouts around the Galactic centre. Each cutout is composed of 600 × 600 pixels of the size 0.86′, making 8.6° ×8.6° maps. The left column shows cutouts from the original PySM sky, which feature pixelation artefacts due to their pixel size of 6.9′ (N_{side} = 512). The right column shows modified cutouts obtained by first updating the original map to a higher resolution, 0.86′ (N_{side} = 4096), by oversampling using the HEALPy function ‘ud_grade()’ and then smoothing with a Gaussian kernel of a width that is close the original pixel size of the PySM maps (FWHM = 10′). This reprocessing of the PySM maps removes the pixelation artefacts while only marginally lowering the spatial resolution. 
3.2. The extragalactic components
The extragalactic components include the CMB, tSZ, CIB, and the kinematic Sunyaev–Zeldovich (kSZ) effect. These components were generated at various frequencies using template maps provided by the Websky extragalactic CMB mocks (Stein et al. 2020), which are HEALPix fullsky maps of the extragalactic sky at N_{side} = 4096 (pix_{size} = 0.86′). The largescale structure of the Websky simulations was obtained through a peakpatch and 2ndorder Lagrangian perturbation theory (2LPT) simulation with 12288^{3} particles and a box size of 15.4 h^{−1} Gpc. The authors assumed a flat Λ cold dark matter (ΛCDM) cosmology with h = 0.68, Ω_{cdm} = 0.261, Ω_{b} = 0.049, Ω_{Λ} = 0.69, n_{s} = 0.965, and σ_{8} = 0.81. The simulations can be downloaded from their website^{2}. An overview of the provided template maps, their frequencies, and units can be seen in Table 1.
Summary of the base extragalactic templates of the Websky simulations.
To generate the tSZ map, Stein et al. (2020) projected the pressure profiles from hydrodynamical simulations (Battaglia et al. 2012) onto the massPeak Patch halo catalogue. They did not consider field contributions, and the minimum halo mass considered to compute their pressure profile was ∼1 × 10^{13} M_{⊙}. Their tSZ power spectrum is in good agreement with the Planck data (Stein et al. 2020). Websky provides allsky, dimensionless Comptony maps (linked to the tSZ effect through Eqs. (1) or (6)), which can be multiplied by the tSZ spectrum in Eqs. (3) and (5) to obtain frequency maps of the tSZ. The tSZ signal was simulated by Websky in its nonrelativistic limit.
The CIB was generated using the halo occupation distribution model proposed by Shang et al. (2012) and with the parameters given in Viero et al. (2013) to fit the CIB power spectrum with the Herschel observations. Within this model, the restframe SED of a source is a modified blackbody at low frequencies and a power law at high frequencies. The CIB’s SED depends on the source redshift, halo mass, and frequency at which it is observed. The Websky halo catalogue is populated with central and satellite galaxies. The field contribution was not considered. The Websky simulations generated CIB maps at the Planck HighFrequency Instrument (HFI) frequencies (see Table 1), and in particular, the 545 GHz CIB power spectrum was normalised to fit the Planck 545 GHz CIB power spectrum (Planck Collaboration I 2014) at ℓ = 500.
Because of the redshift dependence of the CIB distribution, a decorrelation was expected between the CIB at different frequency channels. The decorrelation was determined in the Websky simulations (Stein et al. 2020) by computing the Pearson productmoment correlation coefficients on the power spectrum of the CIB maps, evaluated at the Planck HFI frequency channels. The correlation matrix between the different CIB maps provided by the Websky simulations is illustrated in Fig. 2. The Websky CIB decorrelation values are in agreement with Planck Collaboration I (2014) and Lenz et al. (2019).
Fig. 2. Correlation matrix between the Websky CIB maps provided at the Planck HFI frequencies. The axis shows the frequencies of the maps in GHz, and the dimensionless colour bar shows the correlation coefficient. 
In this work, we extrapolated the CIB model to arbitrary frequencies by fitting each pixel of the existing Websky CIB maps with the following modified blackbody:
where ν_{0} = 353 GHz. This was done using the Scipy “curve_fit()” function to derive three allsky parameter value maps for the variables A_{CIB}, T_{dust}, and β. These fullsky maps were injected in Eq. (8) to produce Webskylike CIB intensity maps at new frequencies.
In the first approximation, we assumed this exact modified black body across frequencies, which is not fully realistic since the CIB is composed of different sources at different redshifts, and their SED is not the same. The Websky simulation includes crosscorrelation between the tSZ and CIB, which is in agreement with the Planck data results (see Fig. 11 of Stein et al. 2020).
For the CMB component, the Websky simulations provided us with two files containing random realisations of the primary T/Q/U a_{ℓm} coefficients, one being unlensed and one being lensed. The lensed file was obtained using the unlensed a_{ℓm} file and the CMB lensing convergence from 0 < z < 1100 in the Born approximation. From those a_{ℓm}, we could reconstruct a lensed and an unlensed CMB map in μK_{CMB}. Using those maps and the conversion factor (Eq. (4)), we could obtain a CMB map at any desired frequency. In the first approximation, we used the unlensed CMB. For the kSZ component, the situation is analogous to the tSZ one. The Websky templates provided frequencyindependent maps of the kSZ, which can be seen as a ‘y_{kSZ}map’. We multiplied the maps by the conversion factor in Eq. (4) to get a kSZ map at any frequency. The Websky simulations do currently not provide radio point source (RPS) templates, but the Websky team will add them in a future release. This extragalactic component is therefore not accounted for in our study. The different extragalactic fullsky map templates were all converted to the same unit, that is, either temperature units (K_{CMB}) or intensity units (Jy sr^{−1}) units, using Eq. (4).
3.3. Future experiments
We simulated the microwave sky as will be seen by two future groundbased telescopes. One of these telescopes is the FYST of the Cerro Chajnantor Atacama Telescope prime (CCATprime) Collaboration (CCATprime Collaboration 2023). It is an upcoming 6 m diameter telescope that will be located at 5600 m on the mountain of Cerro Chajnantor in Chile. Due to its very low precipitable water vapour, this exceptional site will provide one of the best atmospheric transmissions for groundbased microwave observations to date (Bustos et al. 2014). The crossedDragone optical design proposed by Niemack (2016) will deliver a highthroughput, wideFoV telescope capable of rapidly scanning large areas of the sky. The first light is expected in 2024.
The second telescope of interest will be part of the Simons Observatory (SO; Ade et al. 2019). The SO will be composed of three refracting 0.4 m small aperture telescopes (SATs) and one crossDragone 6 m large aperture telescope (LAT). The four instruments will be located at 5300 m on Cerro Toco in the Atacama Desert in Chile.
We studied what benefits may result from the combination of FYST with the SO LAT. To do so, we used PySM and the Websky mock sky to simulate the microwave sky at frequencies and sensitivities, meaning the smallest signal a telescope can detect above the random background noise, as well as the beam sizes of future observations conducted with these telescopes, see Table 2.
Frequencies, beam sizes in FWHM, and baseline instrumental sensitivities for SO LAT and FYST.
3.4. Noises, masks, and beams
Our Skymodel pipeline includes two types of noise: white noise, to reproduce the instrumental noise, and atmospheric red noise, which is present for groundbased experiments. The noise implementation is based on the SO Collaboration (Ade et al. 2019) noise curves forecasted and adapted for CCATprime by Choi et al. (2020a):
where ℓ_{knee} = 1000 and α_{knee} = −3.5. The values for N_{red} and N_{white} for the SO LAT and FYST are included in the Skymodel. For this reason, the atmospheric noise model can only be computed at valid SO LAT frequencies or FYST frequencies (see Table 2) and takes the form of a series of power spectra, as defined by Eq. (9) and shown in Fig. 3. Using those power spectra, the Skymodel generates fullsky map realisations of the atmospheric noise or atmospheric and instrumental noise. This is done through the HEALPy function ‘synfast()’, which generates map realisations from the power spectrum while assuming the field to be Gaussian. This was done for each frequency band. In reality, neighbouring frequency bands will observe similar regions of the sky at the same time. Thus, the atmospheric noise is partially correlated. We adopted a 90% correlation between each pair of neighbouring frequency channels as suggested in Ade et al. (2019). Fixing the ‘synfast()’ seed allowed for the production of a 100% correlated map for two neighbouring channels. Another seed and a correction factor can be applied to the power spectrum to produce a downweighted correlated map and an uncorrelated map of one of the two channels. This new map and the fully correlated map can be linearly combined to produce a 90% correlated map using
Fig. 3. SO LAT (blue lines) and FYST (red lines) sum of the atmospheric red and instrumental white noise spectra based on the noise model presented by Ade et al. (2019) and adopted by Choi et al. (2020a). The noise curves have been beam corrected. 
where M is the final map at the frequency ν, which will be 90% correlated with its neighbouring channel. The HEALPy synfast() function is represented by h, and is the atmospheric noise curve, which is 100% correlated when produced from the same seed as its neighbouring channel and 0% when produced from another seed.
The frequency bands are broad and extended over some frequency range, and their transmission is not perfect but varies over their extent. In the first approximation, we considered the frequency bands to be delta functions with perfect transmission.
In addition to the maps of astrophysical emission and noise, the Skymodel uses survey masks. The masks are stored at N_{side} = 256 and were upgraded to N_{side} = 4096 in order to correspond to the Websky map resolution. In this study, we used a mask that excludes the Galactic plane with a sky fraction of f_{sky} = 0.6 (see Fig. 4). This mask was derived in Erler et al. (2018) by masking the 40% brightest pixels of a Planck Galactic emission map. This is an initial approximation, as in reality SO LAT and FYST will observe a smaller fraction of the sky.
Fig. 4. Galactic foreground mask derived from Planck data by Erler et al. (2018; f_{sky} ≈ 0.6). The white part is the masked region, and the blue part is what was used in our analysis. Overplotted is the HEALPix tessellation scheme for N_{side} = 4. 
The Skymodel pipeline also simulates the effect of telescope beams using the Gaussian FWHMs given in Table 2. The product of the Skymodel is a fullsky HEALPix map containing all the Galactic foreground components, the extragalactic components, instrumental white noise, and atmospheric red noise at a given frequency and assuming a given beam size (see Appendix D).
4. Method to extract thermal Sunyaev–Zeldovich map
To study the power spectrum measurement of the tSZ effect for the SO LAT and FYST telescopes, a tSZ map needs to be extracted from the simulated microwave sky. The sky is composed of various Galactic and extragalactic emissions that are considered to be contaminants of our signal of interest. Component separation methods permit the recovery of a signal from the multicomponent microwave sky. Internal linear combination (ILC) is a multifrequency component separation technique that was first used by Bennett et al. (2003) to extract CMB maps from WMAP data and later refined by Eriksen et al. (2004). Multiple flavours of ILC exist. In their paper, Tegmark et al. (2003) suggested applying it in Fourier space rather than in map space. Needlet ILC (NILC; Delabrouille et al. 2009) processes the maps in both spatial and harmonic space by proposing to filter each multifrequency map by a set of ℓdependent window functions called needlets in order to apply multiple ILCs on different multipole intervals. This makes the ILC weights ℓdependent, allowing them to take advantage of the scaledependent distribution of astrophysical components to minimise the noise and better extract the signal of interest. The Modified ILC Algorithm (MILCA; Hurier et al. 2013) differentiates itself by allowing multiple constraints on the Astrophysical emissions in order to reduce their contamination. Scalediscretised directional wavelet ILC (SILC; Rogers et al. 2016) filters the map with ℓ windows that are direction dependent.
In this work, we used ILC and CILC methods in a mapbased space. It is also possible to use these methods directly in Fourier space, which is computationally faster but does not allow for spatial separation of signals. Fourier space ILC/CILC only has information on the power spectra of the component in ℓspace and no spatial information on the localisation of the signal over the fullsky map.
4.1. Internal linear combination
Internal linear combination is widely used in tSZ studies (Bennett et al. 2003; Fornazier et al. 2022). Its main advantage is that it is a socalled “blind” method because only the SED of the component of interest must be known, and no information on the contaminants is needed. The ILC supposes that (1) the maps are composed of a linear combination of astrophysical components and noise, and (2) the astrophysical components are uncorrelated.
A set of N_{obs} observed map m_{i}(p) at a frequency i and a pixel p can be expressed as
where a_{i} is the amplitude of the frequency spectrum of the signal of interest s at each frequency i. The noise at a given frequency i and pixel p is represented as n_{i}(p).
To retrieve the signal of interest s, an estimator is defined as
where ω_{i} is the ILC weight at the frequency i. The ILC weights are defined such that they follow two conditions. One of them is that in order to obtain an unbiased estimate, the ILC weights must have a unit response to the frequency spectrum of the signal of interest s characterised by a_{i} at each frequency i (see Eq. (13)):
The other condition is that the minimisation of the presence of the unwanted astrophysical contaminants and the instrumental noise in the reconstructed map of the signal s is achieved through the condition of minimum variance of the reconstructed map, VAR(). As shown in Eriksen et al. (2004), because the components are uncorrelated, we have VAR(, where is the value of the covariance of the observed maps between the frequencies j and k. Thus, the minimum variance condition is given by
The minimum variance under the constraint of Eq. (14) can be achieved using the Lagrange multiplier method to solve the above equation. The solution gives the following ILC weights:
where W is the vector containing the ILC weights, and A is the vector containing the frequency spectrum information of the signal of interest and is called ‘the mixing vector’. The covariance matrix of the set of maps is represented by . The map containing our signal of interest is given by
From that, we can define the residual noise (r_{N}), that is, the noise left over by the ILC in the reconstructed map. This noise can be defined as the difference between the ILCrecovered map and the true signal:
The contaminants are not absent from the recovered signal map. Moreover, the low variance contaminants might not be reduced in the resulting map. This becomes a problem when assumption two of the ILC is not fully verified and our target signal is spatially correlated with another emission, such as for the tSZ and CIB. This is the case between CIB and tSZ. The remaining CIB signal in the map can bias the ILC estimate of the tSZ. One way to reduce this effect is through CILC.
4.2. Constrained internal linear combination
The CILC method was first introduced by Remazeilles et al. (2011). It adds one extra constraint to the ILC approach. If the spectral shape of one of the contaminants is known, it can be taken into account in the modelling such that instead of Eq. (11), it is supposed that
where s(p) is the signal of interest and q(p) is the contaminant we want to suppress. Respectively, a_{i} and b_{i} are the intensity of the frequency spectrum evaluated at a frequency i of the signal and contaminant. All the remaining astrophysical contaminants and the instrumental noise at a given frequency i and pixel p are contained by n_{i}(p). The CILC weights are defined such that
The CILC weights have a unit response to the SED of the signal of interest s and a null response to the SED of the astrophysical contaminant q we want to suppress so that this contaminant is absent from the reconstructed signal.
The condition of minimum variance (see Eq. (14)) stays the same, and only the expression of the weights is changed to take into account the two conditions imposed by the CILC,
where W is the vector containing the CILC weights, B is the mixing vector of the contaminant, A is the mixing vector of the signal, and is the covariance matrix of the set of maps. Those weights are then applied to the set of multifrequency maps using Eq. (16) to recover a map containing the signal of interest free from the nulled contaminant and with minimum contamination from other astrophysical signals and instrumental noise. The CILC method can be generalised to null more than one contaminant if their SED is known. However, these additional constraints reduce the degree of freedom available for the condition of minimum variance to be verified (Hurier et al. 2013), thus generating increased noise for the CILC compared to the ILC. Similar to the ILC, the residual noise can be defined as the difference between the CILCrecovered map and the signal of interest (see Eq. (17)), except that the weights are given by Eq. (20).
4.3. Beam convolution
The set of multifrequency maps was smoothed with a Gaussian kernel to the lowest frequency channel resolution so that they could be linearly combined by the ILC. No meaningful physical information is available on scales smaller than the resolution of each map. When linearly combining the maps, the upper limit of meaningful physical information was therefore given by the lowest resolution map.
More sophisticated component separation methods, such as NILC, allows the original resolution of each frequency map to be kept. Indeed, the NILC weights are a function of the multipole ℓ for each frequency map and can therefore be set to zero when the scale exceeds the resolution of the map. Setting the weights to 0 above a given ℓ, allows the lowerresolution maps to contribute to the signal reconstruction only up to their beamlimited scale, leaving the higherresolution maps to be the contributors of signal reconstruction on smaller scales.
We note that when simulating the microwave sky, only the astrophysical components and atmospheric noise are limited by the instrument beam. Instrumental noise is not beam limited. The inclusion of the beam was done using HEALPy to smooth the maps a_{ℓ, m} in the multipole space and multiply them with a ℓdependent Gaussian kernel. Our set of multifrequency maps can be expressed as
where are the a_{ℓ, m} coefficient of the smoothed map at the frequency i. The Gaussian beam in the multipole space at the frequency i is indicated as b_{i}(l). Next, a_{ℓ, m}(s) is the a_{ℓ, m} decomposition of component of interest map, and a_{ℓ, m}(c_{i}) represents all astrophysical contaminants and the atmospheric noise at a frequency i. Finally, a_{ℓ, m}(n_{i}) is the instrumental noise a_{ℓ, m} at the frequency i. In general, to apply an ILC, all maps have to be smoothed down to the lowest resolution. This is necessary to avoid adding information from other higher resolution, higher frequencies maps that are smaller than the beam size of some of the maps into the recovered map. The smoothing of the maps to the lowest map resolution was done using HEALPy, and this process can be formalised as
where are the a_{ℓ, m} coefficient of the map at a frequency i smoothed down to the lowest resolution. The Gaussian kernel at this desired ‘lowest’ resolution is indicated at b^{f}. The consequence of the necessary smoothing to the lowest resolution is an increase of the instrumental white noise contribution in each frequency map by a factor b_{f}(ℓ)/b_{i}(ℓ) because for all the frequency channels, except the lowest one, we have b_{f}(ℓ) > b_{i}(ℓ).
4.4. Fullsky treatment
The microwave fullsky maps are composed of many different astrophysical signals. Some of those signals are spatially localised. The Galactic foregrounds are mainly localised around the Galactic plane. When applying the mapbased ILC to a fullsky map, instead of applying it straightforwardly to the multifrequency set of maps, thus getting only one set of weights for the full sky and equally treating regions that are strongly contaminated by Galactic foregrounds and regions that are less contaminated, we tessellated the sky into different zones and applied an ILC on each of these zones separately. This resulted in a set of different weights for each region so that they are always optimised for the level of contamination present in their zone. This has been done previously in various studies (Tegmark et al. 2003). The most common tessellation approach is to separate the sky into zones of Galactic contamination.
To tessellate the sky, we used the native nested HEALPix tessellation (see Fig. 4). This nested scheme tessellates the sky originally into 12 equalarea fields. The sides of each of those fields can be further divided to create more equalarea cells of smaller sizes. Using HEALPix nested scheme to create subfields (Górski et al. 2005), we chose to tessellate the sky using N_{side} = 4, which gives fields, each covering an area of ≈215 deg^{2}. In general, tessellation can create border effects along the limits of the fields. In Eriksen et al. (2004), the border effects were attenuated by smoothing the ILC weights. In this work, we do not correct the border effects.
The Galactic plane was masked (see Fig. 4) with a f_{sky} ≈ 0.6 Galactic mask map to reduce foreground contamination and better extract the tSZ effect. To compute power spectra on a masked map, we used PyMaster (Alonso et al. 2019), which corrects the power spectra for the masking of part of the sky. In this work, we do not apodise the masks, as it is computationally very intensive at our high resolution and does not make a significant difference in the resulting power spectrum (see Fig. A.2). The full treatment which includes, creation, tessellation, masking and ILC is performed at the map level. Only when a tSZ fullsky map was extracted was the power spectrum of the map computed using spherical harmonics.
4.5. Debiasing
The instrumental white noise biasing of the tSZ power spectrum extraction was estimated for a 4000 h observation time (Choi et al. 2020b). Therefore, the mean of this noise could be computed and removed to debias the tSZ reconstruction. This was done by generating random realisation of the instrumental white noise for each frequency at the sensitivities given in Table 2 and computing the projected noise by multiplying the ILC weights from the tSZ reconstruction with the noise maps and summing them. The power spectrum of this noise was then subtracted from the ILCextracted Comptony (y_{ILC}) power spectrum. All the following power spectrum results were debiased for instrumental white noise.
Similarly, as the CMB power spectrum is wellknown, we could generate random realisation from the Websky CMB power spectrum to debias for the CMB residual noise. In the first approximation, this study presents an idealised case where the CMB power spectrum is perfectly known and its residual is removed. With real data, this CMB debiasing will be limited by the uncertainties on the recovered CMB power spectrum.
4.6. Overview of steps
In this section, we describe our tSZ extraction method and illustrate it in Fig. 5.
Fig. 5. Scheme of the method we use in this paper to extract the tSZ signal. The input of the sky model is template maps coming from the Websky simulations and PySM. The sky model processes them to produce maps of the different sky components at a given frequency, sensitivity, and beam. This set of multifrequency maps {M_{i1}, …, M_{in}} was smoothed to a common resolution and tessellated, and the Comptony signal was extracted by the ILC. The PyMaster algorithm was used on the resulting masked map to extract the power spectrum. 
First, using the Skymodel based on the Websky simulation and PySM, we generated multifrequency maps of the microwave sky at the frequencies, sensitivities, and beam sizes of SO LAT and FYST (see Table 2). The set of maps is in the HEALPix format and uses intensity MJy sr^{−1} units.
Second, the simulated maps were all smoothed down to the resolution of the lowest frequency channel used in our analysis.
Third, a Galactic mask that excludes 40% of the sky (see Fig. 4) was applied to the set of multifrequency maps.
Fourth, each map was tessellated into smaller equalarea regions using the HEALPix nested format.
Fifth, an ILC was applied to each multifrequency set of the tessellated area to obtain the Comptony patches. The patches were reassembled to obtain the fullsky Comptony map. The covariance matrix was computed only on the unmasked data by our ILC.
Sixth, the power spectrum of the recovered Comptony map was computed using PyMaster and debiasing, and the residual noises were computed using Eq. (17).
5. Results
To quantify how much CIB could leak into the reconstructed tSZ maps for future experiments and access the detectability of the tSZ power spectrum, we used the microwave sky model detailed in Sect. 3. We simulated a set of multifrequency maps at the frequencies, sensitivities, and beam sizes for SO LAT and FYST (see Table 2) with atmospheric noise. We applied an ILC or CILC component separation method to extract a tSZ map and looked at the ratios between the different noises and signal power spectra. We restricted ourselves to the scales ℓ ∈ [500, 5000]. We chose this lower limit to avoid the cosmic variance and the upper limit because of the 2.2′ resolution of our HEALPix maps, which allowed us to resolve scales up to ℓ ∼ 4909.
5.1. Internal linear combination results
By construction, the weights of an ILC must verify the first condition, which is to have a unit response to the SED of the tSZ effect (see Eq. (13)). This first condition means that the elementwise product of weights (ω_{i}) times the mixing vector (a_{i}) is an indicator of the fractional contribution of each frequency band to the recovery of the tSZ signal. It can be negative since the frequency spectrum of the tSZ effect is negative below 217 GHz. The weights attributed by the ILC can themselves be negative for optimisation purposes, with the only constraint being that the sum over all frequencies must be one. The value of the weights is given by the mean of the 192 weights (one per tessellated field), and the error bars are given by the standard deviation of those weights. The top panel of Fig. 6 shows that for SO LAT alone (orange dots), the 93 GHz channel is the main contributor to the tSZ reconstruction by a factor ω_{93}a_{93} ≈ 1.4. This channel overestimates the tSZ signal, and the SO LAT 145 GHz channel was used by the ILC to correct this. As a result, the 145 GHz channel contributes roughly by a factor ω_{145}a_{145} ≈ −0.4 to the tSZ reconstruction. This result is not surprising, as those frequency bands are located at the tSZ decrement that is characteristic of the signal in the microwave sky and where the values of the mixing vector are negative. When SO LAT and FYST are combined (SO LAT+FYST; blue dots), the SO LAT 93 GHz and 145 GHz channels remain the main contributors to the tSZ reconstruction and have similar values, ω_{93}a_{93} ≈ 1.45 and ω_{145}a_{145} ≈ −0.45. Other channels from SO LAT and FYST do not contribute to the reconstruction. When FYST is alone (green dots), the main contributor to the tSZ reconstruction is the 280 GHz channel, by a factor ω_{280}a_{280} ≈ 1.1, and the 350 GHz channel, by ω_{350}a_{350} ≈ 0.05. The 220 GHz and 405 GHz channels were used to mitigate the over evaluation of the tSZ signal with the respective factors ω_{220}a_{220} ≈ −0.05 and ω_{405}a_{405} ≈ −0.1.
Fig. 6. Internal linear combination weights for SO LATlike (orange) and FYSTlike (green) experiments and SO LAT+FYST (blue) when the sky contains all extragalactic components (see Table 1), Galactic foregrounds, and atmospheric noise. Top figure: fractional contribution of each frequency channel to the recovered ymap. Bottom figure: frequency channels that contribute to maximising the tSZ signal and minimising the noise. The error bars show the standard deviation around the mean weight of all of the 192 tessellated fields for each frequency map. 
The second condition that has to be verified by the ILC weights is to ensure that the variance of the recovered map has been minimised (see Eq. (14)). The bottom panel of Fig. 6 shows the values of the weights. In combination with the top panel, we concluded that the weights associated with some frequency bands that do not contribute to the tSZ signal reconstruction but have a high absolute value contribute to the noise minimisation. In particular, the 220 GHz channel located close to the tSZ null is often attributed a high weight by the ILC even though it does not contribute to the tSZ reconstruction. This is because it is the only map containing all the noise contaminants and almost no signal. It is thus the perfect map for noise minimisation. But this is not the case in Fig. 6, except for when FYST is alone (green dots), because of the presence of atmospheric noise in our mock sky. The strong atmospheric noise contaminates the 220 GHz band, making it less exploitable by the ILC (see the comparison with Fig. A.1 for a case without atmospheric noise). Figure 6 shows that in the presence of all the microwave emissions and atmospheric noise, the SO LAT frequency channels that probe the tSZ decrement, 93 GHz and 145 GHz, drive both the tSZ reconstruction and the noise minimisation of the tSZ signal.
However, the weights do not tell us anything about the CIB contamination or the detectability of the tSZ signal. For that, we turned to the power spectrum and the residual ILC noise power spectra. The ILC output is a map containing the Comptony signal and residual noises (Eq. (17)). Using PyMaster, the power spectra of this map and the residual noises can be computed and compared to the inputpure Comptony map power spectra from the Websky simulations. We defined as the power spectra of the quantity X or of the residual noise r_{X} of the quantity X. For SO LAT alone (see Fig. 7 topleft panel) the dominant noise residual on small scales is the CIB residual (r_{CIB}), which starts to dominate around ℓ ∼ 3700. The cumulative foreground residual (r_{f}) and the kSZ residual (r_{kSZ}) are both subdominants to the input signal by around one order of magnitude over the scales ℓ ∈ [500, 5000]. The atmospheric noise residual (r_{Atmo}) dominates the input signal for ℓ < 1400 and ℓ > 4500.
Fig. 7. Power spectra of the ILC residual noises compared to the input Comptony power spectrum from the Websky simulations (red). The simulated sky contains all extragalactic components (see Table 1), Galactic foregrounds, and atmospheric noise. Top panels: power spectra of various quantities. The grey curves show the ILC noises residual power spectra, denoted by r_{X} (X being the kSZ); the cumulative Galactic foregrounds (f); and the atmospheric noise (Atmo). The black curve shows the power spectrum of the CIB residual noise. The blue curve shows the sum of all grey and black curves, that is, the total ILC noise residual debiased for instrumental white noise and the CMB. The topleft panel shows the results for the simulated SO LAT data, while the topright panel shows the results for the simulated data for SO LAT and FYST combined. Bottom panel: ratio of the input Comptony power spectrum (red) and the total noise residual (blue) for SO LAT and SO LAT+FYST. All power spectra have been beam corrected and bin averaged over a window Δℓ = 50. 
In the combined case of SO LAT+FYST (see right panel of Fig. 7), the situation is similar to the case where SO LAT is alone except that the atmospheric noise residual (r_{Atmo}) is higher when FYST is present. This could be due to higher atmospheric noise being present in the highfrequency channels of FYST (see Fig. 3). The CIB noise residual power spectrum is significantly lower for SO LAT+FYST. It starts to dominate above the input signal around ℓ ∼ 4200, while for SO LAT alone it was around ℓ ∼ 3700. We concluded that adding FYST to SO LAT reduces the CIB noise residual.
The lower panel of Fig. 7 shows the extent the input Comptony power spectrum is dominated by the total ILC noise residual as a function of the scale. The total noise residual (r_{tot}) power spectrum dominates over the input signal on nearly all scales except around ℓ ∈ [2100, 2700] for SO LAT alone. This is because of our debiasing of the CMB on those critical scales. Overall, SO LAT alone has a better input/noise ratio on scales ℓ < 2500 because of its slightly lower atmospheric noise residual, and SO LAT+FYST does better than SO LAT alone on scales ℓ > 3000 because of its lower CIB noise residual.
5.2. Constrained internal linear combination results
Constrained internal linear combination can be used to ‘null’ the contribution of a component of the microwave sky using its frequency spectrum. Nullifying several components was tested.
Before any debiasing, the CMB is the main contaminant of the recovered Comptony power spectra on large scales, so it seemed natural to deproject (i.e. ‘null’) it. However, the additional constraints on the ILC weights (see Eq. (18)) led to a smaller solution space for the weights to extract the signal, which resulted in a ‘noise penalty’, that is, higher noise in the recovered signal map. In this particular case of a mock sky containing all Galactic foregrounds and extragalactic emission, excluding RPSs, and applying a CILC to deproject CMB, the ILC weights are so constrained to null CMB that they become less efficient at reducing the second main contaminant, the atmospheric noise contamination, which blows up. As a consequence, the total residual noise dominates the input signal even more, by up to one order of magnitude on small scales (see Appendix B). Deprojecting the CMB is therefore not a recommended option for our study. This effect was already reported in Fig. 36 of Ade et al. (2019), but in their study, the blowing up of the noise was mitigated on large scales by the addition of Planck data.
Because of its highfrequency coverage, we expected FYST to better probe the CIB, which is the dominant contaminant above ν > 300 GHz, and therefore better target CIB contamination and remove it. Deprojecting the CIB using CILC seemed the natural option to take full advantage of the potential constraining power of the FYST highfrequency channels. We note that a perfect deprojection of the CIB is not possible because its signal is not composed of only one SED. We adopted the modified blackbody spectrum given in Eq. (8), with A_{CIB} = 1, T_{dust} = 10 K, and β = 1.6. Parameter values were calculated by computing the average for each of the parameter maps derived in Sect. 3.2.
The top panel of Fig. 8 shows that the CIB noise residual (r_{CIB}) power spectrum is reduced by the CILC deprojection for both SO LAT and SO LAT+FYST, but it is not zero because the CIB is composed of multiples spectra and only one is used for the deprojection. It also shows that the atmospheric noise residual (r_{Atmo}) power spectrum increases for SO LAT and decreases for the combined SO LAT+FYST when the CIB is deprojected (see Fig. 7 for comparison). This can be explained by the higher degrees of freedom available with the higher number of frequency channels for SO LAT+FYST, making the CILC noise penalty have a less penalising effect in comparison to when SO LAT is alone. The lower panel of the figure shows that SO LAT+FYST offers a large ℓ ∈ [1700, 3500] window over which the Comptony signal is not noise dominated, while for SO LAT alone, the signal remains noise dominated (exclusively by the atmospheric noise). For SO LAT alone or for SO LAT+FYST combined, compared to the ILC case, the ratio between the input target signal and the resulting noise is not much higher in general (Fig. 7). This was also seen when applying our pipeline to Planck data, as the recovered Comptony map amplitude is roughly the same as the ILC and CILC. This is not surprising because we also saw that the ILC and CILC weights are very similar. This observation of the CILC not providing a much lower noise bias Comptony map compared to the ILC was already reported by Alonso et al. (2018).
Fig. 8. Power spectra of the CILC residual noises when deprojecting CIB compared to the input Comptony power spectrum from the Websky simulations (red). The simulated sky contains all extragalactic components (see Table 1), Galactic foregrounds, and atmospheric noise. Top panels: power spectra of various quantities. The grey curves show the ILC noise residual power spectra, denoted by r_{X} (X being the kSZ); the Galactic foregrounds (f); and the atmospheric noise (Atmo). The black curve shows the power spectrum of the CIB residual noise. The blue curve shows the sum of all grey and black curves, that is, the total ILC noise residual debiased for instrumental white noise and the CMB. The topleft panel shows the results for the simulated SO LAT data, while the topright panel shows the results for the simulated data for SO LAT and FYST combined. Bottom panel: ratio between the input Comptony power spectrum (red) and the total noise residual (blue) for SO LAT and SO LAT+FYST. All power spectra have been beam corrected and bin averaged over a window Δℓ = 50. 
5.3. Overview
Figure 9 gives an overview of the relative amplitude changes for the power spectrum of residual noises of each astrophysical component when combining SO LAT and FYST and shows a comparison to SO LAT alone in the two previous cases applying an ILC (in blue) and a CILC deprojecting the CIB (in orange). To see how much SO LAT+FYST reduces residual noise power spectra in comparison to SO LAT alone, we took the ratio of the residual power spectra of each component in both cases:
Fig. 9. Value of the ratios of the power spectra of X for SO LAT+FYST compared to SO LAT(). The ILCextracted Comptony signal (y_{ILC}) and r_{X} represent the residual noise of the quantity X. The shaded regions are the density distribution of the ratios over the ℓ range. The bar is the mean averaged over the ℓ range of the ratio distribution. In blue is the ILC case (see Sect. 5.1) and in orange is the CILC case when CIB is deprojected (see Sect. 5.2). 
where is the power spectrum of the component X or of the residual noise r_{X} for SO LAT+FYST and is the same quantity but for SO LAT alone. For the extracted Comptony map (y_{ILC}), the blue violin plot shows a mean improvement of around 4% for SO LAT+FYST compared to SO LAT alone but with a large scatter. Especially at ℓ < 2800, SO LAT performs better (see bottom panel of Fig. 7). When deprojecting the CIB (orange violin plot), we moved to a mean improvement of around 8% in the ILCextracted Comptony map for SO LAT+FYST compared to SO LAT alone. In this case, the scatter is also much smaller than before because of the high number of degrees of freedom offered by the combined SO LAT+FYST frequency bands, and the CILC noise penalty is much lower for SO LAT+FYST compared to SO LAT alone, making the combination strictly superior to extract a less biased tSZ power spectrum in the CILC case.
For the residual CIB noise (r_{CIB}), the power spectrum is on average around 35% lower when combining SO LAT and FYST in the ILC case. In the CILC case, it drops to ∼7%, which might not be surprising because the CILC focusses on nulling as much CIB as possible no matter the frequencies. The addition of FYST does reduce the CIB residual even more but not as much as without the null constraint. Another reason for the drop in improvement when combining SO LAT and FYST in the CILC case might be the higher atmospheric noise at the high frequencies where the CIB signal dominates, undermining the advantage of FYST to constrain and remove CIB.
In the ILC case, the CMB residual noise (r_{CMB}) power spectrum is reduced by ∼22% when combining SO LAT with FYST compared to SO LAT alone. Because the ILC cannot separate CMB and kSZ due to their identical spectral frequency distribution, the residual kSZ noise (r_{kSZ}) is also suppressed when combining SO LAT and FYST, on average by around 18%. When deprojecting CIB with a CILC, adding FYST to SO LAT also reduces the CMB and kSZ residual noises but only by ∼5% over all scales on average.
The sum of all galactic foreground noise residuals is defined as the cumulative foreground residual (r_{f}). Its noise power spectrum is improved on average by ∼10% for SO LAT+FYST compared to SO LAT in the ILC case. In the CILC case, the cumulative foreground residual power spectrum is reduced by ∼30% when combining FYST with SO LAT. This could be due to the additional constraint on the CIB dust spectrum, nulling part of the Galactic dust contribution, which is more easily nulled with FYST highfrequency coverage.
The FYST highfrequency channels are more contaminated by the atmosphere. For this reason, the atmospheric red noise residual (r_{Atmo}) power spectrum is ∼30% higher when SO LAT and FYST are combined rather than when SO LAT is alone in the ILC case. When deprojecting the CIB, the atmospheric noise residual becomes slightly lower ∼1% for SO LAT+FYST compared to SO LAT alone. This is due to the noise penalty generated by the additional constraint in the CILC worsening the atmospheric noise for LAT alone compared to SO LAT+FYST. Because SO LAT+FYST benefits from additional frequency bands and therefore more degrees of freedom, it is more resilient to the inevitable noise increase caused by additional spectral constraints.
6. Summary and discussion
Using fullsky HEALPix template maps from the Websky simulations (Stein et al. 2020) and PySM, we created a pipeline that allowed us to simulate the highresolution maps of the different microwave sky emissions at any beam, frequency, and sensitivity. Based on the Choi et al. (2020a) noise modelling, our pipeline also simulates atmospheric red noise at the locations of the Simons Observatory and FYST.
We simulated the sky that would be seen during experiments at these locations and used a mapbased ILC and CILC to extract back the thermal Sunyaev–Zeldovich (tSZ) effect from this multifrequency, multicomponent sky. To minimise Galactic foreground contamination, we used a Galactic mask derived from Planck (see Fig. 4). To optimise our ILC weights in order to take advantage of the spatial distribution of some of the microwave emissions, we tessellated the sky using the HEALPix equalarea nested tessellation scheme. We made use of CILC to deproject the CIB by using a simple modified blackbody spectrum to represent the CIB frequency variation and constrain the ILC weights to nullify it. We then used PyMaster to compute the power spectrum and its associated residual noise power spectra, which are the noise coming from all the other microwave sky emissions left by ILC and CILC in the recovered map.
We note that our study does have certain shortcomings that should be taken into account. For example, the Websky simulations do not yet offer a map of RPSs, and these extragalactic components are therefore absent from our study. The RPSs are an important contaminant of the tSZ signal, as they are located in galaxy clusters. If RPSs are not considered, they can fill up the tSZ decrement and bias an estimation towards a lower estimated tSZ signal. A solution to this problem is to mask the RPSs, but this also lowers the estimated tSZ power spectrum, as it removes parts of the tSZ signal. This effect especially impacts measurements of the tSZ signal below 90 GHz, which we do not study here. A study from Holder (2002) showed that RPS subtraction can lead to an underestimation of the tSZ power spectrum of between ∼30% and ∼50%, depending on the scale. We did not consider the relativistic correction of the tSZ effect, which can lead to a significant underestimation of the signal (Remazeilles et al. 2018). This correction is left to future studies. Regarding the method, our tessellation scheme can lead to border effects that were reported by Eriksen et al. (2004) but not accounted for. Our pipeline was tested on Planck data with and without tessellation and only a marginal (i.e. < 5%) difference could be seen in the power spectrum amplitude. Applying ILC in mapbased space also has the disadvantage of the necessary resolution downgrade of all maps to the lower beam channel, which is 2.2′ for SO LAT, thus reducing the benefits of FYST’s higher resolution channels. This also prevented us from using Planck data that are essential to lowering atmospheric and CMB contamination on large scales. This was addressed by performing CMB debiasing of the residual CMB noise in the extracted Comptony map. A Fourier space or needlet ILC implementation is left for future work.
We find that when using an ILC to extract tSZ, the combination of SO LAT with FYST reduces the CIB residual noise power spectrum by ∼35%, on average, on scales ℓ ∈ [500, 5000] when compared to SO LAT alone. For the CMB and kSZ components, the addition of FYST to SO LAT reduces the residual noise by ∼22% on average, with a strong scale dependence for the CMB, as the addition of FYST mainly helped on large scales. For the kSZ component, the reduction observed when adding FYST to SO LAT is not scale dependent. Adding FYST also helps reduce the Galactic foreground contamination by, on average, ∼10% on scales ℓ ∈ [500, 5000]. However, FYST’s highfrequency channels are more affected by atmospheric noise due to poorer transmission. Therefore, the residual noise power spectrum coming from the atmosphere is on average 30% higher when SO LAT and FYST are combined compared to when SO LAT is alone. This greatly affects the total residual noise budget, lowering the benefits of adding FYST to SO LAT for the tSZ recovery to only an average ∼4%, and on larger scales, SO LAT alone performs better. We find that over ℓ ∈ [2100, 2700], the ILCextracted tSZ power spectrum is not noise dominated for SO LAT alone. But on smaller scales, it is extremely biased, and its noise budget is dominated by the CIB. When SO LAT is combined with FYST, the CIB noise residual dominates much later ℓ ∼ 4200, but because the atmospheric noise is higher, there is no window over which the ILCextracted tSZ power spectrum is not noise dominated.
We note that the CIB residual noise reduction that FYST brings to SO LAT is not something that SO LAT could achieve on its own by simply having a longer survey time. This was tested by reducing the instrumental white noise in SO LAT bands by two, and we observed a reduction of the ILC residual instrumental noise in the reconstructed Comptony map, but the other noise components remained the same. This ‘boosted’ version of SO LAT alone did not have a lower CIB residual power spectrum compared to the planned SO LAT or SO LAT plus FYST combination.
Because tSZ and CIB have a nonzero correlation in the Websky simulations and because FYST’s high frequencies probe the CIB better, we used CILC to deproject the CIB component, modelling it with a simple onecomponent modified blackbody spectrum. Compared to the ILC case, adding FYST to SO LAT lowers the CIB noise residual by only ∼7%. The CMB and kSZ noise residuals are also reduced by ∼5% when combining SO LAT and FYST when using CILC. However, the cumulative foregrounds are suppressed by ∼30%, and atmospheric noise is a bit lower for the SO LAT and FYST combination. Despite a lower reduction of the CIB residual noise when combining SO LAT and FYST, the CIB deprojection case offers a large ℓ ∈ [1800, 3500] window where the CILCextracted Comptony map is not noise dominated after instrumental noise and CMB debiasing.
Acknowledgments
The authors would like to thank Kaustuv Basu for his guidance and advice throughout the project. Srinivasan Raghunathan, Jacques Delabrouille, and Frank Bertoldi for their time, valuable insights, and numerous helpful discussions. The referee for useful comments and numerous improvements to the draft. The editor Laura Pentericci for her advice and help on the submission process. Michael Coronado for his language editing work which significantly helped improve the manuscript. The CCATprime Collaboration for their feedback on this work and Colin Hill for his help. Laila Linke for her useful comments on this draft. We furthermore acknowledge support from the International MaxPlanck Research School (IMPRS) and the BonnCologne Graduate School of Physics and Astronomy (BCGS). The simulations of the CIB used in this paper were developed by the Websky Extragalactic CMB Mocks team, with the continuous support of the Canadian Institute for Theoretical Astrophysics (CITA), the Canadian Institute for Advanced Research (CIFAR), and the Natural Sciences and Engineering Council of Canada (NSERC), and were generated on the GPC supercomputer at the SciNet HPC Consortium. SciNet is funded by the Canada Foundation for Innovation under the auspices of Compute Canada, the Government of Ontario, Ontario Research Fund – Research Excellence, and the University of Toronto.
References
 Ade, P., Aguirre, J., Ahmed, Z., et al. 2019, JCAP, 2019, 056 [Google Scholar]
 Alonso, D., Colin Hill, J., Hložek, R., et al. 2018, Phys. Rev. D, 97, 093514 [Google Scholar]
 Alonso, D., Sanchez, J., Slosar, A., et al. 2019, MNRAS, 484, 4127 [NASA ADS] [CrossRef] [Google Scholar]
 Battaglia, N., Bond, J. R., Pfrommer, C., et al. 2012, ApJ, 758, 75 [NASA ADS] [CrossRef] [Google Scholar]
 Bennett, C. L., Bay, M., Halpern, M., et al. 2003, ApJ, 583, 1 [NASA ADS] [CrossRef] [Google Scholar]
 Bolliet, B. 2018, arXiv eprints [arXiv:1806.04786] [Google Scholar]
 Bustos, R., Rubio, M., Otárola, A., et al. 2014, PASP, 126, 1126 [NASA ADS] [CrossRef] [Google Scholar]
 CCATprime Collaboration (Aravena, M., et al.) 2023, ApJS, 264, 7 [NASA ADS] [CrossRef] [Google Scholar]
 Chluba, J., & Sunyaev, R. A. 2012, MNRAS, 419, 1294 [Google Scholar]
 Choi, S. K., Austermann, J., Basu, K., et al. 2020a, J. Low Temp. Phys., 199, 1089 [NASA ADS] [CrossRef] [Google Scholar]
 Choi, S. K., Hasselfield, M., Ho, S.P. P., et al. 2020b, JCAP, 2020, 045 [CrossRef] [Google Scholar]
 Delabrouille, J., Cardoso, J.F., Le Jeune, M., et al. 2009, A&A, 493, 835 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Dickinson, C., AliHaïmoud, Y., Barr, A., et al. 2018, New Astron. Rev., 80, 1 [CrossRef] [Google Scholar]
 Dunkley, J., Calabrese, E., Sievers, J., et al. 2013, JCAP, 2013, 025 [CrossRef] [Google Scholar]
 Eriksen, H. K., Banday, A. J., Górski, K. M., et al. 2004, ApJ, 612, 633 [NASA ADS] [CrossRef] [Google Scholar]
 Eriksen, H. K., Jewell, J. B., Dickinson, C., et al. 2008, ApJ, 676, 10 [Google Scholar]
 Erler, J., Basu, K., Chluba, J., et al. 2018, MNRAS, 476, 3360 [NASA ADS] [CrossRef] [Google Scholar]
 Fixsen, D. J. 2009, ApJ, 707, 916 [Google Scholar]
 Fixsen, D. J., Dwek, E., Mather, J. C., et al. 1998, ApJ, 508, 123 [NASA ADS] [CrossRef] [Google Scholar]
 Fornazier, K. S. F., Abdalla, F. B., Remazeilles, M., et al. 2022, A&A, 664, A18 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 George, E. M., Reichardt, C. L., Aird, K. A., et al. 2015, ApJ, 799, 177 [Google Scholar]
 Górski, K. M., Hivon, E., Banday, A. J., et al. 2005, ApJ, 622, 759 [Google Scholar]
 Haslam, C. G. T., Klein, U., Salter, C. J., et al. 1981, A&A, 100, 209 [NASA ADS] [Google Scholar]
 Haslam, C. G. T., Salter, C. J., Stoffel, H., et al. 1982, A&AS, 47, 1 [Google Scholar]
 Hauser, M. G., Arendt, R. G., Kelsall, T., et al. 1998, ApJ, 508, 25 [Google Scholar]
 Holder, G. P. 2002, ApJ, 580, 36 [NASA ADS] [CrossRef] [Google Scholar]
 Hurier, G. 2015, A&A, 575, A11 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Hurier, G., MacíasPérez, J. F., & Hildebrandt, S. 2013, A&A, 558, A118 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Itoh, N., Kohyama, Y., Nozawa, S., et al. 1998, J. Phys. Condens. Matter, 10, 11273 [NASA ADS] [CrossRef] [Google Scholar]
 Kim, J., Naselsky, P., & Christensen, P. R. 2008, Phys. Rev. D, 77, 103002 [CrossRef] [Google Scholar]
 Komatsu, E., & Seljak, U. 2002, MNRAS, 336, 1256 [NASA ADS] [CrossRef] [Google Scholar]
 Lenz, D., Doré, O., & Lagache, G. 2019, ApJ, 883, 75 [Google Scholar]
 Mather, J. C., Cheng, E. S., Cottingham, D. A., et al. 1994, ApJ, 420, 439 [Google Scholar]
 Melin, J.B., Bartlett, J. G., Cai, Z.Y., et al. 2018, A&A, 617, A75 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Melin, J.B., Bartlett, J. G., Tarrío, P., et al. 2021, A&A, 647, A106 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Menanteau, F., González, J., Juin, J.B., et al. 2010, ApJ, 723, 1523 [Google Scholar]
 MivilleDeschênes, M. A., Lagache, G., Boulanger, F., et al. 2007, A&A, 469, 595 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Mroczkowski, T., Nagai, D., Basu, K., et al. 2019, Space Sci. Rev., 215, 17 [Google Scholar]
 Niemack, M. D. 2016, Appl. Opt., 55, 1686 [Google Scholar]
 Park, C.G., Park, C., & Gott, J. R., III 2007, ApJ, 660, 959 [NASA ADS] [CrossRef] [Google Scholar]
 Partridge, R. B., & Peebles, P. J. E. 1967, ApJ, 148, 377 [NASA ADS] [CrossRef] [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 XXI. 2014, A&A, 571, A21 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Planck Collaboration IX. 2015, A&A, 594, A9 [Google Scholar]
 Planck Collaboration XXII. 2016, A&A, 594, A22 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Planck Collaboration XXV. 2016, A&A, 594, A25 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Planck Collaboration XXVII. 2016, A&A, 594, A27 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Planck Collaboration IV. 2020, A&A, 641, A4 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Pointecouteau, E., Giard, M., & Barret, D. 1998, A&A, 336, 44 [NASA ADS] [Google Scholar]
 Puget, J.L., Abergel, A., Bernard, J.P., et al. 1996, A&A, 308, L5 [NASA ADS] [Google Scholar]
 Reichardt, C. L., Shaw, L., Zahn, O., et al. 2012, ApJ, 755, 70 [Google Scholar]
 Reichardt, C. L., Patil, S., Ade, P. A. R., et al. 2021, ApJ, 908, 199 [NASA ADS] [CrossRef] [Google Scholar]
 Remazeilles, M., Delabrouille, J., & Cardoso, J.F. 2011, MNRAS, 410, 2481 [Google Scholar]
 Remazeilles, M., Dickinson, C., Banday, A. J., et al. 2015, MNRAS, 451, 4311 [NASA ADS] [CrossRef] [Google Scholar]
 Remazeilles, M., Bolliet, B., Rotti, A., & Chluba, J. 2018, MNRAS, 483, 3459 [Google Scholar]
 Rogers, K. K., Peiris, H. V., Leistedt, B., et al. 2016, MNRAS, 460, 3014 [NASA ADS] [CrossRef] [Google Scholar]
 Schaan, E., Ferraro, S., Amodeo, S., et al. 2021, Phys. Rev. D, 103, 063513 [NASA ADS] [CrossRef] [Google Scholar]
 Shang, C., Haiman, Z., Knox, L., & Oh, S. P. 2012, MNRAS, 421, 2832 [NASA ADS] [CrossRef] [Google Scholar]
 Staniszewski, Z., Ade, P. A. R., Aird, K. A., et al. 2009, ApJ, 701, 32 [CrossRef] [Google Scholar]
 Stein, G., Alvarez, M. A., Bond, J. R., et al. 2020, JCAP, 2020, 012 [CrossRef] [Google Scholar]
 Sunyaev, R. A., & Zeldovich, Y. B. 1970, Comments Astrophys. Space Phys., 2, 66 [NASA ADS] [Google Scholar]
 Sunyaev, R. A., & Zeldovich, Y. B. 1972, Comments Astrophys. Space Phys., 4, 173 [NASA ADS] [EDP Sciences] [Google Scholar]
 Sunyaev, R. A., & Zeldovich, Y. B. 1980, ARA&A, 18, 537 [NASA ADS] [CrossRef] [Google Scholar]
 Tegmark, M., de OliveiraCosta, A., & Hamilton, A. J. 2003, Phys. Rev. D, 68, 123523 [Google Scholar]
 Thorne, B., Dunkley, J., Alonso, D., & Næss, S. 2017, MNRAS, 469, 2821 [Google Scholar]
 Viero, M. P., Wang, L., Zemcov, M., et al. 2013, ApJ, 772, 77 [Google Scholar]
 Vikram, V., Lidz, A., & Jain, B. 2017, MNRAS, 467, 2315 [NASA ADS] [Google Scholar]
 Wright, E. L. 1979, ApJ, 232, 348 [Google Scholar]
 Zubeldia, I., Chluba, J., Battye, R., et al. 2023, MNRAS, 522, 5123 [NASA ADS] [CrossRef] [Google Scholar]
Appendix A: Planck data
A.1. Internal linear combination and constrained internal linear combination weights
We applied our pipeline to Planck HighFrequency instrument (HFI) Data Release 2^{3} (DR2; Planck Collaboration XXII 2016) in order to see how well we could recover the true y using an ILC. We also considered a CILC to deproject CIB by constraining the weights so that they have a null response to the frequency spectrum of the CIB (see Eq. (20)). However, the CIB deprojection is not perfect, because the CIB frequency spectrum depends on the population of the sources and their redshift.
Looking at the ILC weights gave us information about which channel contributes the most to extracting the tSZ signal and minimising the noise. These were the two conditions that the weights had to satisfy (see Eq. (13) and (14)). The ILC weights ponderated by their respective frequency mixing vector value showed the fractional contribution of each frequency channel to the recovered y. The fractional contribution of the ∼220 GHz band was null because the value of the mixing vector at this frequency is null.
For a mapbased ILC applied to Planck data (see blue points in Fig. A.1), we concluded from the weight that the 143 GHz channel is the one that contributes the most, with a factor of ω_{143}a_{143} ≈ 0.7, to the reconstructed final map as well as to minimising the noise. The 217 GHz channel does not contribute to the final reconstructed y map because the tSZ signal is null at this frequency. However, its high ω_{217} value makes it the main contributor to noise minimisation. This was expected, as the tSZ signal is absent at this frequency and the map contains all the noise contaminants. Next are the contributions of the 353 GHz channel and the 545 GHz channel. The contribution values of those channels are consistent with previous studies (Erler et al. 2018). In the case where a CILC is applied to Planck data to deproject the CIB (orange points, Fig. A.1), we saw that the weights are almost identical to the ILC ones. This is expected, as Alonso et al. (2018) also found that deprojecting CIB with a CILC does not result in much less CIB contamination in the reconstructed tSZ profile, arguing that this shows the ILC already does a sufficient job of cleaning the CIB.
Fig. A.1. Internal linear combination (blue) and CILC (orange) weights with CIB deprojected on Planck data. Top figure: Fractional contribution of each frequency channel to the recovered ymap. Bottom figure: Frequency channels that contribute to maximising the tSZ signal and minimising the noise. The error bars show the standard deviation around the mean weight of all of the 192 tessellated fields for each frequency map. 
A.2. Processing of the maps
We tested the consistency and effect of each of the postprocessing steps applied on the map (see Fig. A.2). The blue curve in the figure shows the ILCextracted power spectrum computed with a Galactic mask and without tessellation (see Fig. 4), while the lighter blue dotted curve shows the same but with apodisation. The black curve shows the case with the tessellation of the sky and the ILC being applied to each rhombus separately. The grey dotted line shows when apodisation was also applied. For our particular case of tSZ extraction, the mask apodisation only lowers the power spectra of the maps by less than 5%. Comparing the power spectrum with tessellation (black) and without (blue) showed that tessellation helps reduce the noise on large scales, up to ℓ ≈ 2000, by a factor of approximately two. At ℓ = 2500, it increases the noise by around 15%. This can be understood as the size of the rhombus being small enough to break down regions highly contaminated by unwanted astrophysical components on large scales and to weigh down such regions in order to better recover the tSZ signal. But the rhombus size might be too big to optimally deal with smaller region noise. Moreover, with 192 total rhombuses over the sky, we might be creating border effects on smaller scales and thus increasing the noise residual on those scales. The scales (ℓ > 2200) are also affected by the ∼9.66’ beam.
Fig. A.2. Power spectra of the Comptony map using our ILC pipeline on Planck DR2 with different configurations. In blue is the power spectrum of the extracted map with a simple f_{sky} = 0.6 Galactic mask. In black is the same thing but with a tessellation of the sky. The dotted lines (blue and black) indicate the inclusion of apodisation in each case. 
Appendix B: Constrained ILC: Deprojecting the cosmic microwave background
When applying an ILC to extract back a Comptony map, the CMB is the dominant noise residual (r_{CMB}) on large scales (see Fig. B.1). Using CILC to deproject the CMB, the residual CMB noise is suppressed by three orders of magnitude. However, the additional constraints imposed on the weights to deproject the CMB led to a smaller solution space and therefore increased noise. This led to a one order of magnitude increase of the atmospheric residual noise (r_{Atmo}), which worsens the estimate of the extracted Comptony power spectrum, especially on small scales, by more than one order of magnitude. This is commented on extensively in Sec. 5.2 and was tested using our Skymodel pipeline (see Sec. 3).
Fig. B.1. Power spectra of ILC/CILC (full lines) and their noise residuals (dotted lines). Specifically, the orange curve is the power spectrum of the ILCextracted Comptony map. The lighter orange lines are the noise residuals. Specifically, they are the ILC residual CMB noise (r_{CMB}) and r_{Atmo} for the atmospheric noise. The blue lines represent the case with CILC where the CMB is deprojected. The solid, darker blue curve is the power spectrum of the CILCextracted Comptony map. The lighter blue curves are the noise residuals. 
Appendix C: ILC weights
In our study, we tesselated the sky into 192 rhombi of equal area using the HEALPix nested scheme. This sky separation allowed us to better adapt to the spatially dependant contaminant of some Astrophysical emissions, such as the Galactic foregrounds. The ILC we applied on each zone resulted in one weight per frequency and per rhombus. An example of a full sky weights map produced and used by our ILC can be seen in Fig. C.1.
Fig. C.1. Maps of the weights for each patch of the sky where the ILC was applied for the 860 GHz frequency band. The weights were computed from simulations that included all Galactic foregrounds and extragalactic components (see Table 1) at the sensitivities, beam sizes, and frequencies of SO LAT and FYST (see Table 2) and with 90% correlated atmospheric red noise generated using the noise curves given by Choi et al. (2020a), (see Fig. 3). 
Appendix D: Simulated maps
In this work, we used the Skymodel pipeline based on PySM and the Websky simulations to generate mock maps of the microwave sky at the frequencies, sensitivities, and instrument beam size of the upcoming SO LAT and FYST (see Table 2). The atmospheric red noise was generated using the noise curves of Choi et al. (2020a) and the red noise of the neighbouring frequency channels of a given telescope were correlated at 90% using Eq. 10. The resulting simulated microwave sky maps used can be seen in Fig. D.1.
Fig. D.1. Mock HEALPy maps of the microwave sky simulated using the Skymodel pipeline (see Sec. 3). All Galactic foregrounds and extragalactic components (see Table 1) at the sensitivities, beam sizes, and frequencies of SO LAT and FYST (see Table 2) and the 90% correlated atmospheric red noise between neighbouring channels generated using the noise curves given by Choi et al. (2020a) were included (see Fig. 3). 
All Tables
Frequencies, beam sizes in FWHM, and baseline instrumental sensitivities for SO LAT and FYST.
All Figures
Fig. 1. Projections of the PySM sky at 353 GHz composed of freefree, synchrotron, AME, and Galactic dust. The upper row shows cutouts around the Galactic south pole, and the lower row shows cutouts around the Galactic centre. Each cutout is composed of 600 × 600 pixels of the size 0.86′, making 8.6° ×8.6° maps. The left column shows cutouts from the original PySM sky, which feature pixelation artefacts due to their pixel size of 6.9′ (N_{side} = 512). The right column shows modified cutouts obtained by first updating the original map to a higher resolution, 0.86′ (N_{side} = 4096), by oversampling using the HEALPy function ‘ud_grade()’ and then smoothing with a Gaussian kernel of a width that is close the original pixel size of the PySM maps (FWHM = 10′). This reprocessing of the PySM maps removes the pixelation artefacts while only marginally lowering the spatial resolution. 

In the text 
Fig. 2. Correlation matrix between the Websky CIB maps provided at the Planck HFI frequencies. The axis shows the frequencies of the maps in GHz, and the dimensionless colour bar shows the correlation coefficient. 

In the text 
Fig. 3. SO LAT (blue lines) and FYST (red lines) sum of the atmospheric red and instrumental white noise spectra based on the noise model presented by Ade et al. (2019) and adopted by Choi et al. (2020a). The noise curves have been beam corrected. 

In the text 
Fig. 4. Galactic foreground mask derived from Planck data by Erler et al. (2018; f_{sky} ≈ 0.6). The white part is the masked region, and the blue part is what was used in our analysis. Overplotted is the HEALPix tessellation scheme for N_{side} = 4. 

In the text 
Fig. 5. Scheme of the method we use in this paper to extract the tSZ signal. The input of the sky model is template maps coming from the Websky simulations and PySM. The sky model processes them to produce maps of the different sky components at a given frequency, sensitivity, and beam. This set of multifrequency maps {M_{i1}, …, M_{in}} was smoothed to a common resolution and tessellated, and the Comptony signal was extracted by the ILC. The PyMaster algorithm was used on the resulting masked map to extract the power spectrum. 

In the text 
Fig. 6. Internal linear combination weights for SO LATlike (orange) and FYSTlike (green) experiments and SO LAT+FYST (blue) when the sky contains all extragalactic components (see Table 1), Galactic foregrounds, and atmospheric noise. Top figure: fractional contribution of each frequency channel to the recovered ymap. Bottom figure: frequency channels that contribute to maximising the tSZ signal and minimising the noise. The error bars show the standard deviation around the mean weight of all of the 192 tessellated fields for each frequency map. 

In the text 
Fig. 7. Power spectra of the ILC residual noises compared to the input Comptony power spectrum from the Websky simulations (red). The simulated sky contains all extragalactic components (see Table 1), Galactic foregrounds, and atmospheric noise. Top panels: power spectra of various quantities. The grey curves show the ILC noises residual power spectra, denoted by r_{X} (X being the kSZ); the cumulative Galactic foregrounds (f); and the atmospheric noise (Atmo). The black curve shows the power spectrum of the CIB residual noise. The blue curve shows the sum of all grey and black curves, that is, the total ILC noise residual debiased for instrumental white noise and the CMB. The topleft panel shows the results for the simulated SO LAT data, while the topright panel shows the results for the simulated data for SO LAT and FYST combined. Bottom panel: ratio of the input Comptony power spectrum (red) and the total noise residual (blue) for SO LAT and SO LAT+FYST. All power spectra have been beam corrected and bin averaged over a window Δℓ = 50. 

In the text 
Fig. 8. Power spectra of the CILC residual noises when deprojecting CIB compared to the input Comptony power spectrum from the Websky simulations (red). The simulated sky contains all extragalactic components (see Table 1), Galactic foregrounds, and atmospheric noise. Top panels: power spectra of various quantities. The grey curves show the ILC noise residual power spectra, denoted by r_{X} (X being the kSZ); the Galactic foregrounds (f); and the atmospheric noise (Atmo). The black curve shows the power spectrum of the CIB residual noise. The blue curve shows the sum of all grey and black curves, that is, the total ILC noise residual debiased for instrumental white noise and the CMB. The topleft panel shows the results for the simulated SO LAT data, while the topright panel shows the results for the simulated data for SO LAT and FYST combined. Bottom panel: ratio between the input Comptony power spectrum (red) and the total noise residual (blue) for SO LAT and SO LAT+FYST. All power spectra have been beam corrected and bin averaged over a window Δℓ = 50. 

In the text 
Fig. 9. Value of the ratios of the power spectra of X for SO LAT+FYST compared to SO LAT(). The ILCextracted Comptony signal (y_{ILC}) and r_{X} represent the residual noise of the quantity X. The shaded regions are the density distribution of the ratios over the ℓ range. The bar is the mean averaged over the ℓ range of the ratio distribution. In blue is the ILC case (see Sect. 5.1) and in orange is the CILC case when CIB is deprojected (see Sect. 5.2). 

In the text 
Fig. A.1. Internal linear combination (blue) and CILC (orange) weights with CIB deprojected on Planck data. Top figure: Fractional contribution of each frequency channel to the recovered ymap. Bottom figure: Frequency channels that contribute to maximising the tSZ signal and minimising the noise. The error bars show the standard deviation around the mean weight of all of the 192 tessellated fields for each frequency map. 

In the text 
Fig. A.2. Power spectra of the Comptony map using our ILC pipeline on Planck DR2 with different configurations. In blue is the power spectrum of the extracted map with a simple f_{sky} = 0.6 Galactic mask. In black is the same thing but with a tessellation of the sky. The dotted lines (blue and black) indicate the inclusion of apodisation in each case. 

In the text 
Fig. B.1. Power spectra of ILC/CILC (full lines) and their noise residuals (dotted lines). Specifically, the orange curve is the power spectrum of the ILCextracted Comptony map. The lighter orange lines are the noise residuals. Specifically, they are the ILC residual CMB noise (r_{CMB}) and r_{Atmo} for the atmospheric noise. The blue lines represent the case with CILC where the CMB is deprojected. The solid, darker blue curve is the power spectrum of the CILCextracted Comptony map. The lighter blue curves are the noise residuals. 

In the text 
Fig. C.1. Maps of the weights for each patch of the sky where the ILC was applied for the 860 GHz frequency band. The weights were computed from simulations that included all Galactic foregrounds and extragalactic components (see Table 1) at the sensitivities, beam sizes, and frequencies of SO LAT and FYST (see Table 2) and with 90% correlated atmospheric red noise generated using the noise curves given by Choi et al. (2020a), (see Fig. 3). 

In the text 
Fig. D.1. Mock HEALPy maps of the microwave sky simulated using the Skymodel pipeline (see Sec. 3). All Galactic foregrounds and extragalactic components (see Table 1) at the sensitivities, beam sizes, and frequencies of SO LAT and FYST (see Table 2) and the 90% correlated atmospheric red noise between neighbouring channels generated using the noise curves given by Choi et al. (2020a) were included (see Fig. 3). 

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.