© J.-M. Delouis et al. 2019

1. Introduction

The development of the SRoll global solution has been initiated within the Planck High Frequency Instrument (HFI) Consortium, attempting to reduce the systematic effects that impact the large-scale polarization of the HFI frequency maps at 100, 143, 217, and 353 GHz. The SRoll algorithm (hereafter SRoll1) and its performance when applied to data and simulations, were first described in Planck Collaboration Int. XLVI (2016, hereafter LowEll paper). The efficiency of SRoll led to significant improvements in the measurement of the cosmological reionization optical depth. The basically unchanged algorithm with marginal improvements, described in Planck Collaboration III (2019, hereafter L03 paper), was used for building the 2018 Legacy HFI frequency maps. Those maps are referred to in this paper as HFI2018 maps.

This paper presents the second generation of this algorithm, called SRoll2, developed beyond the Planck Collaboration one, to further reduce the main residual systematic effects left in the HFI2018 maps, firstly the second-order correction of the analog-to-digital convertor nonlinearity (ADCNL); secondly the temperature-to-polarization leakage from the polarization efficiency and imperfect orientation knowledge of the bolometer;thirdly the temperature-to-polarization leakage from bandpass mismatch of foregrounds and associated foreground modeling; fourthly the transfer function associated with very long time constants.

This paper is organized as follows: Sect. 2 presents the SRoll2 algorithm; Sect. 3 presents the resulting SRoll2 maps; Sect. 4 presents the characterization of those maps and estimation of the systematic effect residuals.

2. The SRoll2 solution

As in SRoll1, the improved SRoll2 is a global solution algorithm, estimating instrumental systematic effect parameters from the sky data, as described in detail in the LowEll, and L03 papers. It works on the data in a single frequency band and produces Stokes I, Q, and U maps, using all polarization-sensitive bolometers (PSBs) at this frequency.

SRoll2 is better at handling systematic effects than SRoll1, in particular the ADCNL effect. The nonlinearity of the analog-to-digital convertors (ADC) was measured in the ground tests, as well as in flight (during the warm phase of the Planck HFI when the dilution cooler 3He tank was empty), and corrected in the first module of the Time Ordered Information (TOI) data processing. The accuracy of this correction was not good enough for the CMB channels, and left residuals which were then corrected empirically in the mapmaking process as a time varying linear detector chain gain. The 353 GHz modulated signals were more spread in analog-to-digital units (ADU) making the correction much smaller. Starting from the same TOIs, but not applying this first-order linear correction, SRoll2 builds a full nonlinear correction per detector and per stable pointing observation (also called “ring”) in which, for each ADU of the signal, a nonlinear correction of the ADC transfer function model is fitted iteratively to the difference between the data and the data model at each iteration. As the HFI bolometers have never shown in tests gain variations other than the one associated with the temperature drift of the 0.1 K-plate (which is corrected for), we attribute to the ADC any residual apparent gain variation, and we adopt a detector response stable in time during the mission. This is a major change which should coherently correct several systematic effects associated with the nonlinearity of the ADC residuals.

The SRoll2 ADCNL correction is computed as a function of the ADU value of the demodulated signal based on a set of spline functions (Sect. 2.3.2) which models semi-empirically the ADCNL residuals in each pointing period. This is done at the ring level (and not at the timeline level), reducing the computational requirements in the mapmaking process to an acceptable value. The number of spline functions needed is much smaller than the number of ADUs to correct because the electronic noise and modulation smooth the signal over a large number of ADUs. The ADCNL spline functions thus replace the correction as an empirical time variation of the detector response, and they also account for the nonlinear distortion of strong signals. This decomposition in spline functions varies in time¹ to empirically adapt to the slowly changing shape of the electronic modulation of the signal, which can be considered stable during one ring.

This ADCNL determination method converges if the initial estimated frequency map used is close enough to the exact solution. This is why the ADCNL computation is always done after a computation of corrections of all other systematic effects, including the calibrated response, to produce an estimated frequency map based on the last minimization of the difference between the data rings (called HPRs for HEALPix (Górski et al. 2005) rings) and the model, providing the best correction available. The ADCNL residual semi-empirical correction is fitted as a function of the demodulated signal (not going through a full model of the ADC and the modulation electronic signal). The ADCNL correction is degenerate with the calibrated response variation. Therefore, it is not possible to fit the ADCNL correction simultaneously with the calibration which is the result of the frequency map solution from the previous step. If our assumption that the apparent gain variations are only coming from the ADCNL and not from the bolometer itself is right, the calibration response should converge to a stable one when the ADCNL correction converges. We show in Sect. 2.3.2 that this is happening.

Figure 1 shows the SRoll2 algorithm scheme composed of two steps: the first step (called Step#1) is similar to the SRoll1 algorithm, but with improvements that increase the convergence speed. It builds polarized frequency maps by optimizing mapmaking parameters for absolute bolometer calibration, polarization efficiencies, and orientations, minimizing the discrepancies between bolometer signals within a frequency band, and minimizing the difference with the data.

Fig. 1. Schematic view of SRoll2 iterations for building the Stokes I, Q, and U maps for one frequency band. The top box contains the initial improvements of the maps. Step#1 algorithm uses data HPRs to consistently compute the maps, the detector response, and the systematic effects correction parameters other than ADCNL. Step#2 computes the ADCNL correction that is propagated to the subsequent iteration. In the middle box, Step#1 produces maps used as input for Step#2 that computes the ADCNL correction. This whole process is then iterated 19 times. The bottom box is again a Step#1 to build the final maps.

In order to obtain a model for the ADCNL computation, Step#1 must be first run to produce frequency maps as input for the loop. Simulations have verified that one Step#1 iteration before computing the ADCNL correction is enough to get convergence.

The second step (called Step#2) is new: it starts from the polarized map generated in Step#1 and refines, for each bolometer, the ADCNL parameters which are the spline function weights used to model the ADC correction.

This is contained in the central box of Fig. 1. In about 20 iterations, the ADCNL converges, and the final maps are obtained by returning to Step#1.

This section details the SRoll2 algorithm: the first subsection explains the data model hypotheses, the second subsection details the Step#1 that computes the calibrated detector chain response, but also the systematic correction parameters at the exception of the ADCNL ones, the third subsection describes the Step#2 that fits the ADCNL model parameters. The last subsections address map projection, future improvements, and simulations.

2.1. Improved scheme

Similarly to the SRoll1 data model (see Eq. (1) of L03), the SRoll2 data model for bolometer signal is given by Eq. (1) where indices are:

• b for the bolometer, up to n_bolo;

• r for the stable pointing period (ring), up to n_ring;

• p for the sky pixel;

• h for bins of spin frequency harmonics (up to n_harm), labeled as bin_h = 1 for the first harmonic, bin_h = 2 for harmonics 2 and 3, bin_h = 3 for harmonics 4 to 7, and bin_h = 4 for harmonics 8 to 15;

• f for the polarized foreground, up to n_comp (dust and free-free).

In this paper, “response” refers to the detector chain response of a given bolometer to a CMB signal (in ADU/μK), calibrated on the orbital dipole signal.

(1)

where:

• g_b is the calibrated response of a bolometer chain;

• M_b, r, p is the measured bolometer total signal;

• A_b, r, p = [1,ρ_bcos(2ϕ_b, r, p),ρ_bsin(2ϕ_b, r, p)] is the pointing vector giving the observed pixel for a given bolometer in a given ring;

• S_p = [I_p,Q_p,U_p] is the sky signal in pixel p after subtraction of the orbital dipole assumed to be known with an amplitude invariant in time. I_p, Q_p, and U_p represent the common sky maps seen by all bolometers (excluding the Solar dipole)²;

• ρ_b is the ground-measured polarization efficiency;

• ϕ_b, r, p is the ground-measured detector polarization angle for bolometer b with respect to the north-south axis;

• V_{b, r, p, h} is the spatial template of the empirical time-transfer function computed from timeline convolution with long time constant (e.g., > 1 s);

• γ_b, h is the empirical time-transfer function complex amplitude;

• C_p, f are the input foreground component spatial templates which are stable for a given SRoll2 run;

• δL_b, f is the bolometer bandpass foreground color-correction coefficient relative to the bolometers average bandpass, for foreground f, that is, for each foreground component, we set . The component brightness map from the frequency bandpass averaged over the bolometers is part of the frequency sky signal S_p;

• is the total CMB dipole signal (sum of solar and orbital dipoles), with being the template for the solar dipole with a constant direction and amplitude and being the template of the orbital dipole;

• and are respectively the total dipole and the Galactic signal integrated over the far sidelobes (FSL);

• O_b, r is the offset per pointing period r used to model the 1/f noise, and we set , since the Planck data provide no information on the monopole;

• N_b, r, p is the pixel white noise (electronic and photon noises), with variance ;

• ℵ_b, r is the signal distortion induced by the ADCNL residual effect following the TOI correction for the bolometer b, at ring r.

The ℵ_b, r term, new with respect to SRoll1 data model, is added to describe the ADCNL residual systematic effect as a function of the signal level. This model and the method to fit its parameters from the data redundancy is the key of the correction of the ADCNL systematic effect, including both the time variability and the nonlinearity of strong signals.

2.2. Step#1: computing the bolometer response and associated leakage coefficients

The Step#1 process is similar to the SRoll1 algorithm but dropping the empirical time varying gain which was a first-order correction of the ADCNL systematic effect. Thus, SRoll2 models the detector response as a constant value. In SRoll2, the ADCNL systematic effect is cleaned in the Step#2 process (see Sect. 2.3).

The Step#1 of SRoll2 is an upgraded version of the SRoll1 algorithm that improves the response convergence, and decreases the number of required iterations. Solving for a single response per bolometer for the whole mission necessarily involves solving a nonlinear least-squares equation. Like SRoll1, SRoll2 uses an iterative scheme to solve for the response g_b. We therefore introduce g_b, n where n stands for the iteration number. At iteration n, we set:

(2)

where δg_b, n is the difference between the responses g_b, n and the response g_b associated with SRoll2 inputs (e.g., the foreground templates) and the instrumental model.

The goal is to iteratively fit δg_b, n. In order to reduce the complexity of Eq. (1), we first remove the low-amplitude signal contribution of the FSLs; this does not induce any visible remaining signature signal in specific null tests. We also remove the best estimate of ℵ_b, r, n, initially set to ℵ_{b, r, n = 0} = 0, and then iteratively computed in Step#2 (see Sect. 2.3). The reduced bolometer total signal becomes:

(3)

We define the systematic effects ϵ_b, r, p, including detector time constant amplitudes as:

(4)

where (δL_b, f C_p, f) is the correction to the coefficient intensity to polarization leakage due to the bandpass mismatch between the polarized detectors with respect to (L_f C_p, f) which is the averaged foreground (f) signal over all bolometers within a frequency. (L_f C_p, f) is part of the signal S_p. Those systematic effects are most often much smaller than the signal. To avoid extra leakage from frequency bandpass mismatch and/or beam anisotropy related to strong signal brightness gradients, a mask is used to ensure that ϵ_b, r, p is small enough on the sky used to compute the systematic effect parameters. This is particularly important at 353 GHz where the Galactic plane has strong gradients. The masks are the same as those used in L03, and they retain 86.2%, 85.6%, 84.6%, 86.1% of the sky for the 100, 143, 217, and 353 GHz maps, respectively. Then, for all the used pixels, we have:

(5)

Using Eqs. (1), (2), and (5), Eq. (3) now becomes:

(6)

The term should converge toward zero after several iterations. This is only possible if is known. Equation (5) gives the value of ϵ_b, r, p, and it can be verified to be negligible. Then, in Eq. (6), (δg_b, n/g_b) can be solved only knowing S_p and . The quantity () is modeled and refined iteratively.

The Step#1 algorithm is based on the redundancy of different observations of the signal in pixel p, S_p, made by the same detector at different times, and different detectors in the same pointing period, to determine the response of each bolometer. The key quantity for extracting the response is the observation ring vector residual model R_b, r, p for detector b and pointing period r:

(7)

where:

• T_b, r, p defined as the projection matrix from frequency map to rings (we define A_p as the pointing matrix built from all A_b, r, p for the pixel p);

• is the matrix built from all weighted by the corresponding within one pixel p;

• , γ_h, nV_p, h, δL_f, nC_p, f, O_p, n, and N_p are the matrices respectively built from observations of pixel p by all g_b, n , γ_b, h, nV_{b, r, p, h}, L_b, f, nC_p, f, O_b, r, n, and N_b, r, p, each being weighted by the corresponding (e.g., ).

It is important to note that retrieval of the SRoll2 bandpass mismatch parameter L_f is based on measurements from different detectors in one frequency band giving coherent results. In practice, SRoll2 only cleans the difference in coefficients between each bolometer and their average. The average of the corrections is therefore set to zero in order to make the system invertible.

The model used to fit the parameters we want to measure (e.g., the gain correction δg_b, n), is to minimize the difference for each pixel of the signal for a given bolometer and ring with the signal from the average of all bolometers and rings. The quantity to minimize is defined as³:

(8)

with being a model for . The response value converges if

(9)

In SRoll1, was taken as the total dipole only. Although η was indeed small for 100 GHz, this was not true at 353 GHz where the Galactic signal was not taken into account, leading to a slow convergence. In SRoll2, the is replaced by the ring signal extracted from 1° smoothed HFI2018 frequency maps. The consistency of this signal model with the real sky brightness is much better than 1% at all frequencies, thus η ≪ 0.01, and the convergence speed is spectacularly improved as illustrated in Fig. 2.

Fig. 2. Average on all bolometers of the response correction at each iteration of the 353 GHz map when using only the total dipole (in red) and using the sky model to fit the response (in blue).

After each Step#1, a solution for the ADCNL model is obtained by a Step#2 algorithm. The calibrated response is thus slightly modified by the change of the ADCNL, and, at the subsequent iteration, the response converges toward another value inside the Step#1.

At each Step#1, we minimize the following χ²:

(10)

In Eq. (10), we introduce the g_b, n terms in the weight at the denominator to compensate the g_b, nN_b, r, p term present in the numerator. Without this response term in the weights, the χ² minimization would be biased towards a smaller value of g_b, n.

The χ² minimization through a positive and symmetric matrix inversion is performed with a parallel conjugate gradient. This provides an estimate of , and also , , and .

Subsequently, the cleaned and calibrated signal ℳ_b, r, p is built using

(11)

and a new map can be computed from the clean and calibrated timeline:

(12)

2.3. Step#2: correcting for the ADCNL systematic effect

The ADCNL systematic effect is due to the imperfect knowledge of the analog-to-digital signal conversion. The essential SRoll2 improvement is to correct for this effect by fitting a model of the ADCNL in the data as done for the other instrumental systematic effects. This subsection describes the model, and how we extract its parameters from the data.

2.3.1. ADCNL model

The bolometer analog signal is sampled at 7.2 kHz, and modulated at 90 Hz. The analog-to-digital conversion is performed at 7.2 kHz, and a 180 Hz digital signal is computed as the sum over 40 samples (called “fast samples”) for each half modulation (Lamarre et al. 2010). One TOI sample is described as:

(13)

where:

• m is the fast sample index within one half-modulation;

• S_b, r, p is the signal and photon noise;

• E_m is the positive/negative modulation electronic signal at 90 Hz, initiated as a square but distorted by the stray entrance capacitance after demodulation;

• N_m is the electronic noise per sample, assumed to be a white noise with σ ≈ 4 digital units;

• ℶ(S_b, r, p+E_m+N_m) is the bias introduced by the ADCNL for a given level of ADC input electronic signal (S_b, r, p + E_m + N_m) at the fast sampling rate;

• is the ADCNL bias at the TOI sampling frequency.

In order to remove the ADCNL bias from the measured value , we have to compute Δ_b, r, p. To avoid biasing the ADCNL correction by the noise, we do not compute the Δ_b, r, p function for each TOI sample, but average over a large enough number of samples (typically 50 samples) to make the noise negligible.

For a period t, while the transients of the modulation E_m are not changing significantly, the average bias is:

(14)

where

• is the distribution of the E_m function over the 40 samples around the mean of the E_m;

• ⊗ is a convolution;

• N_m is the Gaussian distribution of the electronic noise;

• (S_b, r, p + ⟨E_m⟩) at each sample p for the ring r is computed using the mean of the known digital 180-Hz value sample at the pixel p.

In Eq. (14), the convolution of the signal by the noise and the electronic modulation leads to ⟨Δ_{b, r, p, t}⟩ being smoothed over approximately several hundred of ADUs.

As demonstrated by simulations, the ADCNL correction for the demodulated signal can be fitted on a limited number of spline bases in ADU level (called 1D spline). When all other systematic effects are removed, the residuals between the sky signal and the data model are attributed to the ADCNL residuals averaged over all observations with this bolometer for which the demodulated signal is in the same range of ADUs. The modulation signal varies slowly and is stable on one ring.

When needed, an additional time-dependant description (called 2D spline) with a limited number of stable correction ranges of rings is introduced. Figure 3 illustrates that, for the bolometer 100-4 b, the 1D spline catches all patterns of the ADCNL residuals. On the contrary, the 2D spline catches small time-dependant effects for a few bolometers, like the bolometer 143-1a, especially at the beginning of the mission. This is due to the variability of the modulation signal.

Fig. 3. Simulated ADCNL residuals (color coded) for the ADU level on the ADC (ordinate axis) along the time, expressed as ring number. Two bolometers are shown: left column for the 100-4b and right column for the 143-1a. First row: input of the simulation, second and third rows: their fitted 1D or 2D spline corrections without noise. Typical calibration is ≃6 μK/ADU.

We note that the function () has a linear component that is degenerate with the linear part of the response. This is the reason for not computing this function in Step#1.

The electronic modulation changes during the mission. In particular, the 4-K lines vary along the mission. Electrical coupling between the ⁴He-JT cooler drive and the bolometer readout appears as very narrow lines (called “4-K lines”) in the power spectrum of the TOI at harmonics of 10 Hz. These are removed from the TOI (see, e.g., Planck Collaboration VII 2016). Those lines modify the electronic modulation, and thus slightly modify the effective ADCNL residuals for a given ring. A given pixel is observed by one bolometer several times within a ring, but with different 4-K line phases when the spacecraft spin frequency is not synchronized with the 4-K line frequency. In order to avoid the effect that would be induced when they are synchronized, the rings for which a 4-K line is a harmonic of the spacecraft spin frequency were not used to build the HFI maps to avoid corruption (Planck Collaboration VI 2014). Thus, the mean ADCNL residual for a given ring is only affected by the mean of the 4-K line variation which, not being synchronous with the sky, tends to cancel out. Therefore, the 4-K line variation is a marginal effect at the ring level and is described as a long time variation between pointing periods.

Finally, the 4-K lines only induce a small variance increase: for the most intense lines and the most sensitive to 4-K lines, bolometer 143-3a displays twice the variance of the electronic noise, otherwise around 25% of the electronic noise. For the 143-3a, the amplitude is multiplied by four during the mission, leading to a smoother ADC residual response at the end of the mission. Therefore, the only issue is to sample well enough the ADCNL function of ADU and time to catch all the nonlinearity variability. In the following section, we show a test on simulation of one of the worst bolometers to determine how many spline parameters are needed to describe the ℶ function.

2.3.2. Fitting the ADCNL model parameters

An estimated residual signal model δℳ_b, r, p is computed using the parameters and the cleaned and calibrated signal from the previous Step#1:

(15)

We identify this residual as being due to the ADCNL correction which is the systematic effect not treated in the first step due to its degeneracy with the gain, and we define:

(16)

where:

• 𝒜_b, r, p is the value in ADU of the signal M_b, r, p.

• a_b is the bolometer linear response of the ADU that is an adjustable addition to the response. The computed detector response in the first Step#1 iteration is affected by the full ADCNL residual which is not corrected yet. Even in the further iterations, the gain is still converging to its final value. To allow the Step#2 to retrieve the ADCNL function not affected by small detector gain residual, it is necessary to include a linear part of the ADCNL which should converge to zero. To avoid any divergence of the algorithm, this a_b is computed again at each iteration.

• ℵ_b, r(𝒜_b, r, p) is the model of the ADCNL for the bolometer b at the sample ring r for the corresponding 𝒜_b, r, p.

• δS_p is the residual signal that is due to the imperfect estimation of all parameters because the ADCNL correction, which is taken to be null in the first iteration, has not yet reached the best estimation at this stage of the algorithm.

Because of the noise and electronic modulation smoothing, the effective ADCNL mismatch can be expressed as a smooth function ℶ of ADU for each bolometer, and is not a pure ADCNL correction(see Sect. 2.3.1). Thus the iteration, in the Step#2, on the weights of the spline description of ℵ_b, r defined in Eq. (17), converges when δS_p becomes small enough. The overall gain iteration converges when the signal variations between two iterations become sufficiently small. Equation (16) shows that ℵ_b, r is fully determined when δS_p = 0.

We describe ℵ_b, r as a spline function:

(17)

where

• t is the index of the spline weight varying with time;

• s is the spline index of the bolometer b;

• w_{t, b, s, n} is the spline weight at the itteration n;

• α_b, s(𝒜_b, r, p) are the spline functions of ADU level for the bolometer b to correct the whole signal around level M_b, r, p, therefore around the corresponding ADU level;

• β_b, t are the spline functions of time.

The iteration, in Step#2, on the weights of the spline description of ℵ_b, r converges when δS_p becomes small enough. The overall gain iteration converges when the signal variations between two iterations become sufficiently small.

For a bolometer, the model for ADCNL residual can be represented by different numbers of spline functions for different frequencies. This is due to the difference in digital range for the signals: for 100, 217, and 353 GHz, the ADCNL is described as 32 third-degree splines in ADU, where at 143 GHz time variation is needed, and ADCNL are described as 8 third degree splines in ADU and 128 splines in time.

From Eqs. (15), (16), and (17), and using again the redundancies in the data, we minimize the following χ² including a implicit binning in ADU due to the smoothness of the splines:

(18)

and we obtain a new estimation of , and an estimation of the totality of w_{t, b, s, n}. changes because the ADCNL correction based on the ADU level that follows the very slow focal plane thermal evolution is degenerate with the very low frequency noise. After this Step#2, a new iteration, with a new gain and a new ADCNL model, can start.

Figure 4 shows the residual between the simulated HPR data and the model given the number of iterations for a simulation without noise. It demonstrates that ADCNL effect can be cleaned very well, as the convergence, measured by the residual spectrum at ℓ = 3, falls below the noise for 10 gain iteration. We therefore limit the number of iterations to 20. Increasing the number of iterations does not improve the results. These ADCNL residuals describe the level of uncertainty reached if one ignores all other systematics at low multipoles which are shown to be smaller in Sect. 4.5.

Fig. 4. EE power spectra for difference between two successive iterations of residual maps from a noiseless simulation including only the ADCNL effect and SRoll2 mapmaking for increasing number of iterations at 143 GHz. 20 iterations are enough to reach a good enough convergence.

2.4. Map projection

At the last iteration, only the Step#1 is performed, and the final map is given by Eq. (12). The resulting SRoll2 frequency maps are presented Fig. 6 and discussed in the following section.

2.5. Future improvements

Using a nonlinear correction of the ADCNL as a substitute for the gain variation in time, the improvements of the polarization leakage and the transfer function are exhibited in Sect. 4 where we show that these systematic effects are under control.

It is also important to note that CO template accuracy is critical as CO has a strong impact on 100 GHz map computation, and a lesser impact at 217 and 353 GHz. At all these frequencies, CO maps can be built from the different single detector response to CO lines, either from ground measurements or using large-enough CO maps from radioastronomy observations as was done in L03. Such CO templates for the two main isotopologues used in this analysis are described in L03 characterization, but are not used in the corresponding HFI2018 maps. The changes that could be brought to CO templates are discussed in Sect. 4.1. At this level of sensitivity, no other molecular line is significant.

Further improvements will be made to reduce the large-scale detector noise due to the undetected cosmic ray glitch tails that now dominate at low multipoles for some frequencies.

2.6. Simulations

SRoll2 simulations, hereafter SRoll2Sims, built for the characterization of SRoll2 data, follow the same schematic as the one described in L03 for the FFP10 simulations. The differences are improved inputs: firstly updated thermal dust template as the “fiducial” FFP10 one (unlike the 300 FFP10 iterations ones); secondly solar dipole as ; thirdly inclusion of a kinetic Doppler quadrupole whose amplitudes, in μK, are 0.84 at 100 GHz, 1.00 at 143 GHz, 1.35 at 217 GHz, and 2.1 at 353 GHz.

All the other simulation pipeline parameters (scanning strategy, optical beams, spectral transmissions, polarization efficiencies, conversion factors, thermal baseline, noise realizations and auto-correlation, ADC nonlinearities, etc.) are unchanged. Besides the inputs, the most significant improvement with respect to the FFP10 simulations is the use of the SRoll2 algorithm at the mapmaking step.

As the sky simulations do not perfectly mimic data, a final post processing is needed to add a small correlated white noise between detectors. Having no correlated noise between bolometers in the timeline simulation allows for a much greater level of computational parallelization.

Figure 5 illustrates the quality of those simulations compared to the data for the half-mission null test. Overall, simulations show good consistency with the data. We see a low frequency noise rise, much reduced with respect to Fig. 17 of L03 below ℓ = 50, as expected. At 143 GHz, at low multipoles, there is a small discrepancy related to the mask used to compute the spectra. This lack of low multipole noise in the simulations is mainly in the temperature spectra while the polarization noise and systematic effect residuals are well modeled. At 353 GHz the simulation spectra show a higher level, demonstrating an overestimation of systematic effects at ℓ <  10. The largest 353 GHz systematic effect is the dust signal distortion which cannot be modeled properly since the sky model is not representative of the real data at very large scale.

Fig. 5. Noise and systematic residuals in EE, BB, and EB spectra, for difference maps of the half-mission null test binned by Δℓ=10. Spectra are rescaled to the full mission, and also rescaled to full sky from the unmasked sky fraction fsky = 0.43 used. Data spectra are represented by thick blue lines, and the average of simulations by thin black lines. The blue band shows one ±σ equivalent of the statistical distribution of the simulations. The linear y-axis scales are adapted to show the full data range on each panel.

3. Results

3.1. SRoll2 frequency maps

Figure 6 presents the SRoll2 “full mission” maps at 100, 143, 217, and 353 GHz maps, the total dipole being removed. Along with those full mission maps, we build half data maps by splitting the data in two sets for which both give a complete set of frequency maps. We can then build null tests with differences of the two sets of maps in a pair to obtain residuals in which the uncorrelated noise and the signal distortion due to systematic effects not identical in the pair will appear, providing information about the noise and systematic effect residuals.

Fig. 6. Full sky, full mission SRoll2 maps at 100 to 353 GHz (in rows), for Stokes I, Q, and U (in columns).

The half-mission (hm) null test is sensitive to the slow temporal variation of instrumental effects (e.g., the ADCNL effect). The detector-set (detset) null test separates the four PSB detectors, at 100, 143, 217, and 353 GHz, into two detector sets. The difference of detset maps can be used as a very useful test of systematics arising from the detector chain specificities not fully accounted for (e.g., calibration or bandpass differences). The odd-even ring (oe) null test, using the difference of two pointing periods that just follow each other, one odd and the other even in the numbering of pointing periods. The difference should be the closest to maps built with the detector noise only. The odd and even six-month surveys are scanning nearly the same ring in opposite directions six months apart. This is very useful to test how well time constants (including very long ones) or asymmetrical FSL systematic effects are corrected. The caveat with this null test is that a map built with only odd (or even) surveys will not be able to properly deal with the polarized systematic effects due to the reduced redundancy in polarization angles. We thus use odd-even survey differences with the systematic effects estimated on the full mission and then build the half maps taking only the HPRs from odd (or even) surveys.

3.2. Comparison of HFI2018 and SRoll2 frequency maps

Here we discuss the temperature maps. The noise in the SRoll2 143 and 217 GHz temperature maps is higher than the one in the HFI2018 maps. This is expected as SRoll2 only uses 8 PSBs where SRoll1 was using both PSBs and SWBs, namely 11 bolometers at 143 GHz and 12 at 217 GHz.

The solar dipole directions and amplitude were measured with high accuracy in L03. The bias in amplitude for 100 simulations of SRoll1 was 1.5 × 10⁻⁴ for the CMB channel, thus about 1σ measured on the dispersion on the 100 simulations, and corrected. This marginally significant correction was applied to maintain coherence with the Planck 2015 determination. In the present analysis, the solar dipole was not extracted again. The solar dipole amplitude residuals show differences, for each frequency, between HFI2018 and SRoll2 maps estimated using a sky fraction of f_sky = 0.43: −0.87 μK (−2.4 × 10⁻⁴) at 100 GHz; −0.62 μK (−1.7 × 10⁻⁴) at 143 GHz; +0.79 μK (2.1 × 10⁻⁴) at 217 GHz.

The negative average bias on the three CMB frequencies is −2 × 10⁻⁴, and reduces to −0.5 × 10⁻⁴ when SRoll1 bias is removed. In summary, SRoll2 does not show any sign of bias with respect to the L03 best value without debiasing SRoll1.

Dipole biases due to SRoll1 miscalibration residuals were estimated to be ±1 μK at CMB frequencies (see L03). For the HFI2018 and SRoll2 simulations, the input absolute calibration was taken to be the same, and thus no calibration solar dipole residual should be present. The solar dipole amplitude differences are at a comparable level as the uncertainties may be driven by other systematic effects than the calibration errors found for the L03. Systematic effects like ADCNL modify the dipole amplitude. The foreground removal, assuming a constant spectral energy distribution (SED) on the whole sky, should be performed using the gradient of SED found in L03 which could induce a dipole error ±1 μK. Thus, the amplitude of the residual solar dipoles, at CMB frequencies, between the HFI2018 and SRoll2 maps should not be interpreted as a calibration mismatch but as being due to a sum of comparable errors in foreground templates and other weak instrumental systematic effects. This is illustrated in Fig. 19. The higher difference at 353 GHz by an order of magnitude (−13 μK or −3.6 × 10⁻³) than at CMB frequencies, is at least partly explained by the uncertainties on the polar efficiencies and the improvement of the time transfer function model.

The 1 μK accuracy on the solar dipole amplitude is thus a limit of SRoll1 and SRoll2 algorithms. The only way to improve previous dipole estimation is to introduce a full component separation in the mapmaking to avoid bias introduced by degeneracy between systematics and foreground dipoles and taking into account the foreground SED variations at large scales.

Figure 7 shows the difference between the HFI2018 and the SRoll2 frequency maps in polarization. These difference maps are expected to show mostly the patterns associated with the systematic effects better handled in the SRoll2 processing.

Fig. 7. HFI2018 minus SRoll2 difference maps at 100 to 353 GHz in Stokes Q (left column) and U (right column), smoothed by a 5° beam.

At 100 GHz, in comparison with 143 GHz where there is no CO line contribution, the improvement in the CO component can be seen in the Taurus and Orion molecular clouds below the Galactic disk, and in the Ophiuchus cloud above the Galactic disk in the central region. This is very likely due to the improved CO full sky map templates, including the ¹³CO, not included before. These are built using the bandpass mismatch coefficients extracted from a specific run using, as input templates, the CO maps of the Taurus molecular cloud from Goldsmith et al. (2008). The use of this small but very well measured CO region to determine better bandpass mismatch is described in Sect. 3.1.3 of L03.

For the dust component, residuals should appear in the same regions as the CO, if the dust was not properly corrected. These are not apparent at 217 GHz, and appear marginally at 353 GHz. This is consistent with the SRoll2 input dust template being the same as the SRoll1 one.

At 100, 143, and 217 GHz, the striping scanning patterns (banana-like patterns) are associated with ADCNL systematic effect residuals. At 353 GHz, those stripes are not visible as the ADCNL effect residuals are much smaller than other instrumental systematic effects due to a wider spread of the signal range on the ADC.

At all frequencies, improvements show up in polarization, thanks to SRoll2 specific polarized parameters adjustments.

Figure 8 shows the 5° smoothed polarized map. For foreground cleaning at 100 GHz, no synchrotron cleaning is performed. After cleaning the dust signal using a single all-sky SED emission law approximation to propagate the 353 GHz polarized information to other frequencies obtained by the correlation factors between CMB frequencies and the 353 GHz after removing the 100 GHz to make it fully dominated by dust. The maps shown are expected to be consistent with noise, signal and foreground or instrumental residuals. We note that the polarized CMB brightness is very small at large scale (less than 0.5 μK). HFI2018 maps had striping patterns along the scanning direction that become much smaller in the SRoll2 maps demonstrating the effectiveness of the algorithm in reducing large-scale instrumental systematic effect residuals.

Fig. 8. Comparison between 100 to 217 GHz HFI2018 and SRoll2Q and U maps, cleaned by the 353 GHz maps smoothed by a 5° beam.

Even at 217 GHz, the differences between the HFI2018 and SRoll2 maps show less large structures for the U signal. Nevertheless, the 217 GHz map dust cleaning using a single coefficient for the whole sky to be applied to the 353 GHz map cannot remove the dust well enough at large scales. In the Planck Legacy release, to get consistency on the solar dipole measured at all frequencies from 100 to 545 GHz, it was necessary to introduce, for the dust removal using the 857 GHz map, an SED variation, on very large scales (dipoles and quadrupoles). It was noted that a similar correction is also indicated by the frequency map intensity ratios (Figs. 20 and 21 of L03). The proper cleaning of the dust foreground at 217 GHz needs to use more parameters at large scales to project the 353 GHz signal to take into account the dust temperature distribution variation over the sky.

Figure 9 shows the spatial pattern consistency between the HFI2018 and SRoll2 maps difference at 100 and 143 GHz and the same for a simulation of the instrumental systematic effects better corrected in SRoll2. The single realization chosen shows similar pattern and amplitude when compared to the data demonstrating the effectiveness of the SRoll2 simulations in capturing the main systematic effects.

Fig. 9. Polarization intensity map difference between SRoll1 and SRoll2 algorithms at 100 and 143 GHz (in columns), applied to the data (first row) and to one realization of the end-to-end simulations (second row).

3.3. Power spectra of null test maps

Figure 10 shows two sets of 100 × 143 power spectra at large scales, for half-mission, detsets, and odd-even rings null tests, EE and BB, for HFI2018 and for SRoll2. Data points are overplotted, and simulation average and error bars are shown. Three hundred simulations are available for the HFI2018, and 500 are built for SRoll2 processing. The very large scales ℓ <  5 show significant deviations with respect to zero both in the data and in the simulations for HFI2018; these are highly improved by SRoll2. HFI2018 simulations show higher variance than SRoll2 at very large scales due to ADCNL correction only to first order. SRoll2 data do not show the obvious biases which were seen and predicted by the simulations for the 2018 release version of the data.

Fig. 10. EE (upper 3 plots) and BB (lower 3 plots) power spectra of difference of null test maps for: half mission, detsets and odd-even rings, for HFI2018 (orange and blue) and SRoll2 (red and green). Data are overplotted on top of the simulations (average and dispersion).

The two data distributions are consistent with the simulations for both data sets as shown by the PTEs for two multipole ranges ℓ = 2 − 10, and ℓ = 2 − 29 in Table 1. Measurements of cosmological parameters (e.g., τ) on the reionization peak using SRoll2 maps will have a negligible bias with respect to the error budget.

Furthermore, the half-mission map difference power spectrum, which is related to the level of ADCNL residuals and noise, is now consistent with the detector noise dominating the odd-even rings map difference. This test demonstrates that the level of ADCNL residuals is now lower than the 1/f detector noise.

The very low multipole variances are substantially reduced by more than a factor of two for ℓ <  5 between HFI2018 and SRoll2. B-modes error bars, assuming r = 0, are computed for the 100 × 143 GHz power spectra, and shown in Fig. 11.

Fig. 11. Comparison between HFI2018 and SRoll2 100 × 143 B-modes variances, computed analyzing the respective E2E simulations.

Figure 12 uses a suite of maps, for both HFI2018 and SRoll2 data, built from half split-data sets, namely detsets, half missions, and odd-even rings. The figure shows EE and BB cross-spectra of the two halves of the split maps showing signal, and the power spectra of the difference of those maps showing noise plus residuals of systematic effects. This is another sensitive and quantitative estimate of the level of improvements in the SRoll2 processing over the SRoll1 one.

Fig. 12. EE (left column) and BB (right column) power spectra of HFI2018 and SRoll2, based on half split data sets (detests, half-mission, and odd-even rings) maps at 100, 143, 217, and 353 GHz. Auto-spectra of the difference maps show the noise plus systematic effect residuals, rescaled to full mission equivalent and corrected for sky fraction used here (43%). The cross-spectra between the same maps are displayed and show the signal. The binning is: δℓ=1 for 2 ≤ ℓ <  30; δℓ=5 for 30 ≤ ℓ <  50; δℓ=10 for 50 ≤ ℓ <  160; δℓ=20 for 160 ≤ ℓ <  1000; and δℓ=100 for ℓ >  1000.

White noises measured at ℓ = 100 − 300 are compatible with HFI2018 data. The low-frequency noise plus systematic residuals in the range ℓ = 6–15 is now compatible with the detector 1/f noise. This is an indication that systematic effects have all been brought to a significantly lower level. This is discussed in Sect. 4.

The odd-even ring difference uses the same systematic effect-removal parameters based on the whole mission and should give the measurement closest to the one limited by the detector noise. The half mission difference gives a very good estimate of the ADCNL distortion of the signal varying on long time scales for which the induced residuals were clearly seen at 143 GHz in the HFI2018 maps. At 143 GHz, such residuals have almost disappeared in the SRoll2 data. This again advocates for future improvements concentrating on reducing the large-scale detector noise (see Sect. 2.5).

At other frequencies, where the half mission differences were more similar to the other null test in the HFI2018 data, the reduction of the residuals is of course smaller.

3.4. Caveats when using the SRoll2 maps

Although the SRoll2 algorithm concentrates on improving the large-scale mapmaking, several systematic effects at small scales have also been reduced by the indirect effect of the improvements discussed above through the degeneracies between the two sets of systematic effects.

Several instrumental effects, already identified in L03, which could improve the maps have not been implemented in SRoll2. Firstly the Beam anisotropy: one should be aware that the beam leakage may correlate temperature and polarized signals at small scales. Secondly the 1/f detector noise at ℓ ≤ 100 shown to be dominated by undetected glitches: no attempt has yet been made to clean this Gaussian noise, and this dominates over the white noise. Thirdly the 4-K line cleaning has not been improved and known features are still present in the maps. Finally the sub-pixel effect of the CO template that affects a few pixels in regions of very strong gradients in the Galactic center and molecular cloud cores is still present. As the 143 GHz map is of course not affected, the 100 and 217 GHz polarized signal in such regions should not be used.

We note that, for the SRoll2 frequency maps, the foregrounds bandpass mismatch coefficients L_f have been cleaned. For this, all detectors have been brought to a common response for each foreground which is the bandpass average of all detectors in the band. Foreground maps cannot be computed from the difference between any single bolometer map⁴ built from the cleaned HPRs.

We also note that the SRoll2 data processing includes nonlinear steps that make the data much more complex than a signal with a simple additive residual corrections on top of a Gaussian noise. As for all previous Planck releases, it is mandatory to validate any use of these data by testing against the results of simulations.

Furthermore, the simulation sky model has a limited representativeness due to the adopted simple foreground emission laws and is also incomplete (e.g., the ¹³CO emission is not simulated).

4. Characterization

4.1. Foreground residuals

To model foreground induced systematic effects, the SRoll1 and SRoll2 algorithms both rely on sky templates to fit their correction parameters from the data redundancies, and then for the implementation of the corrections. The quality of the frequency maps produced depends on the representativity of these input foreground templates with respect to the sky ones. The main foreground templates for the HFI data (¹²CO, ¹³CO, and dust) are used to correct for bandpass mismatch leakage in polarization. In SRoll2, the all-sky templates are the dust thermal and the free-free emissions, based on the Commander Planck 2015 component separation products. The ¹²CO and ¹³CO line emission all-sky templates are built with the leakage coefficients obtained from the Taurus molecular cloud maps from Goldsmith et al. (2008) as described in L03, and used in the characterization of the HFI2018 maps.

Here, we develop a test of the impact of the input foreground templates which are not modified within one SRoll2 run on the final maps. This is done by implementing an iterative approach to improve the accuracy of the foreground templates outside SRoll2. This test is based on a limited characterization which itself is based on an addition to the foreground templates, and a simple model with two foregrounds: the 100 GHz map, assumed to be dominated by the CO bandpass leakage, and the 353 GHz map, dominated by the dust bandpass leakage.

In the first iteration of this toy model, we degrade the input CO and dust foreground templates by adding a colored noise map with the same power spectrum as the foreground but with random phases. This makes the initial CO foreground leakage template significantly different from the CO input maps shown in L03 to be very good by testing them within the Taurus region mapped by Goldsmith et al. (2008). A similar degradation is done for the input I, Q, U dust maps obtained by component separation methods and used as the input foreground component maps in this simulation. Subsequently, we perform a full SRoll2 run to build a set of frequency maps. We then build new foreground templates. The CO templates are built using the bolometer efficiency coefficients improved in the first SRoll2 run, while the dust map is simply the difference between the 353 GHz and the 100 GHz maps, after being cleaned from the new CO templates. We then use these new foreground templates in the following iteration of SRoll2 run. Subsequently, we test that the residuals of the foreground leakage converge to the solution obtained with the nondegraded input templates.

The implementation is done by integrating the two CO template maps in an extended definition of the pointing vector A:

(19)

The first run starts with the degraded templates. Then, in the following runs, the values of the leakage coefficients and are given by the previous iteration of the SRoll2 fit. Using Eqs. (11), (12), and (19), we build the frequency maps and the two new CO leakage template maps. At each step, the 100 GHz single bolometer HPRs are used to compute both the ¹²CO and the ¹³CO templates. The 100 and 353 GHz frequency maps are then used to build the new dust template. This procedure is then iterated. We show in Fig. 13 that four runs, one starting with the best estimates of foreground template and three iterations starting with the degraded templates, lead to a set of foreground template maps showing decreasing differences between iterations. This demonstrates that, despite the strongly degraded quality of the initial template, the procedure provides a method able to quickly build much better component templates for both CO and dust.

Fig. 13. Differences between two runs of SRoll2, one using the best available CO and dust templates, and another using a degraded version of these initial ones (as described in the text) designated as iteration 1. We display the foreground templates evolution: 12CO, 13CO (two top rows) and I, U, Q dust template (three lower rows). The first column is the difference between the best templates (iteration 0) and the run with the degraded templates (iteration1). SRoll2 is then run with foreground templates built from the sky maps and leakage coefficients from the iteration 1 displayed in the middle column. The results from a third iteration are displayed in the right column.

To be more quantitative on the convergence, Fig. 14 shows the EE power spectra of the difference of the 100 and 353 GHz polarized frequency maps built using the nominal input templates and the maps obtained with the degraded ones. The two panels show in purple the starting “bad” input templates. At 100 GHz, the first iteration (in blue), shows, for ℓ <  100, a similar level as the foreground power spectrum purposely added to distort the associated template. The second and third iterations show much smaller differences showing quantitatively the fast convergence of the templates down to differences much smaller than the noise. The convergence toward zero of the dust foreground error at 353 GHz is almost reached as soon as the first iteration.

Fig. 14. Power spectra of the initial templates and then differences between two successive iterations: first iteration in cyan, second in green, and third in red. The differences becomes quickly much smaller than the noise (black line) that is estimated from half mission difference computed with fsky = 0.43 and realigned at the full mission level.

We thus demonstrate here on a simple model with two frequencies and two foregrounds, that, after a few iterations, the foreground templates produced by SRoll2 after adding very different initial templates (random phases) converge to the initial one. This does not demonstrate that the component templates found by this method are the optimal component maps that could be built, but that they are sufficiently good to clean the polarized data from foreground bandpass mismatch leakage to a level much lower than the noise. This provides a path to integrate component separation to the mapmaking, by combining all detectors from different frequencies in the same mapmaking run with the caveat that each foreground component should be negligible in at least one frequency band.

4.2. Transfer function, polarization angles, and efficiencies

4.2.1. The SRoll2 transfer function

The transfer function evaluates how much of the sky signal is either absorbed by the mapmaking into a systematic effect correction, or distorted by the processing. This cannot be done by pure simulation as the models are not statistically reliable models for systematic effects and foregrounds. We thus run SRoll2 with the sky data to which we add a white noise map in I, Q, and U of comparable rms. We subtract from these maps the ones obtained with the sky data alone, and build the power spectra from the resulting map. The ratio of these power spectra to the power spectra of the additional white noise simulated signal gives the transfer function in conditions close to those of the sky. Figure 15 shows these transfer functions measured on the auto-power-spectra at 100 GHz and 143 GHz. This result demonstrates that the signal distortion is of the order of or smaller than 10⁻³ at ℓ = 3–10, and a few 10⁻⁴ above. We also note a very small suppression of power at the level of few 10⁻⁴ in the 143 GHz power spectra. We conclude that there is no evidence for any transfer function effect to be corrected for.

Fig. 15. Transfer function of the SRoll2 100 and 143 GHz maps. Thick line shows the average of ten simulations while the envelop shows the peak-to-peak dispersion.

4.2.2. Thermal time constants

The bolometer thermal chain includes the copper bolometer housing and feed horns, and its interface with the bolometer plate. These have very long time constants (the main time constant of the HFI bolometers is a few milliseconds long). When testing the instrument, the heaters of the temperature control system were used to inject power into the bolometer plate and measure the time delay of the bolometers reaction. Those show a distribution of time constants not associated in propagation in the plate of order 10–30 s. The previous mapmaking algorithms, including SRoll1, could not correct them well. This needed to be improved to reach the detector noise at the largest scales. The shorter time constants identified in the ground tests before flight (up to 2 s) were already corrected in the timeline processing in the previous versions of the data processing. In the HFI2018 maps, the very long time constants were corrected through a temporal transfer function only modeled in three ranges of harmonics of the spin frequency. Within each range, the time transfer function was assumed to be constant. The best null test for evaluating time constant correction is the comparison of odd and even surveys which scan the sky in opposite directions. At 353 GHz, this null test showed (see L03) striping not aligned with the scanning pattern (the so called “zebra effect”), and systematic calibration discrepancies between odd and even surveys. This description of the transfer function does not properly correct the very long time constant systematic effect, although Fig. 16 shows that replacing the description of the transfer function in steps with a continuous complex transfer function description associated with a single time constant (of the order of 25 s) removes the zebra effect.

Fig. 16. Intensity map of the odd-even survey difference at 353 GHz from HFI2018 (top) and SRoll2 (bottom) after changing the model of the real part of the transfer function.

In SRoll2, a new description of the complex time transfer function has been added. The imaginary part for dipoles is the critial one, and we take a similar approach as the Planck release to correct the phase shift only for dipoles.

The real part affects the transfer function for spectra, and is adjusted to a function with a single time constant⁵ which accounts for the more complex behavior. Due to the sky observation period of 60 s, it is not possible to properly identify any time constant longer than 10 s. Figure 16 shows that using this empirical model made of two independent sets of transfer functions corrects the striping. It also corrects the odd-even calibration discrepancy between odd and even surveys. This confirms that the SRoll1 solution was not the right one, and that just adding a single very long time constant was not accounting for the complexity of the heat transfer from the bolometers to the bolometer plate.

Figure 17 shows that the odd-even survey pattern seen in the response variation of the HFI2018 353 GHz data (shown Fig. 39 of L03) has been cleaned in SRoll2 using the new time constant model. Furthermore, the calibration dispersion between surveys is 2.6 × 10⁻⁴ for HFI2018, and 5 × 10⁻⁴ for SRoll2 data. This leads to an estimate of the uncertainty when averaging over all detectors of 2.2 × 10⁻⁴ for SRoll2 for which the noise is dominating. This value is in agreement with the estimate of the detector noise given in Table 7 of L03. The reduced dispersion between bolometers in each survey is clearly smaller than the noise, showing that there was an overfitting, discussed below.

Fig. 17. Residual solar dipole amplitudes with respect to the average dipole, for each 353 GHz PSB, shown per survey.

Figure 18 shows the difference between intensity maps from the two single bolometers within a single PSB pair at 353 GHz, which demonstrates the need for this procedure. The two lower panels show the cases with and without time constant correction. They show the reduction of the residuals brought by SRoll2. In the top panel, the map difference for the L03 mapmaking shows a lower large-scale noise at high latitude than in the SRoll2 solution (bottom panel), which is also again the sign of overfitting.

Fig. 18. Difference of two single-bolometer intensity maps (353-3a and 353-5b) smoothed at 5°. Top: HFI2018; middle: SRoll2 without cleaning the long time constant; bottom: SRoll2 with cleaning the long time constant.

The overfitting is due to the degeneracy of the low-frequency 1/f detector noise with the very long time constant. To avoid the overfitting, the very long time constant is fitted, within a single PSB pair, for bolometer a on a template built from bolometer b. This removes the degeneracy between the very long time constant and the 1/f noise. This method is justified because the very long time constants, due to the bolometer housing interface with the plate, are common to the a and b bolometers.

4.2.3. Detector polarization angles and efficiencies

Although the polarization angle and efficiency errors are partially degenerate in their effect on the polarized maps, a global error on the overall angle of the focal plane would induce EB and TB correlation which are observed to be very low. This is discussed and confirmed for the SRoll2 maps in the following section. The absolute value of the polarization efficiency is degenerate with the polarized map calibration (L03). The polarization angles are well measured and induce a smaller error than the polarization efficiency uncertainties. Furthermore, L03 demonstrated that the ground measurements of the bolometer polarization efficiency were not accurate, especially for SWB. Without an external calibrator, it is not possible to clean the induced absolute polarization mismatch. It is nevertheless possible to compute for each detector a relative polarization efficiency with respect to the average for the frequency band detector set. The following equation introduces the polarization parameter error in the Eq. (1) data model:

(20)

where

• T_ρ = [Q_pcos(2ϕ_b, r, p)+U_psin(2ϕ_b, r, p)],

• T_ϕ = 2ρ_b[U_pcos(2ϕ_b, r, p)−Q_psin(2ϕ_b, r, p)].

Therefore, by fitting both T_ρ and T_ϕ templates, SRoll2 measures the relative errors on polarization efficiency and angle. As mentioned above, T_ρ is degenerate with the polarized signal calibration, and, in order to avoid removing polarization signal from the final map, the mean of all δρ_b within a frequency is set to zero. With T_ϕ being degenerate with the leakage from E to B modes, a similar hypothesis is to set the mean δϕ_b, r, p to zero.

Using end-to-end simulations, we obtain the following uncertainties for the polarization efficiency per bolometer: 2% at 100 GHz, 0.8% at 143 GHz, 0.3% at 217 GHz, and 0.3% at 353 GHz. This shows that no improvement of the polarization efficiency knowledge can be obtained. Thanks to the increase of the dust emission with the frequency, the polarized signal becomes much higher than the noise, and we can thus improve its knowledge on the 2% level quoted in L03 at 217 and 353 GHz.

Figure 19 shows the efficiency of the SRoll2 method in cleaning residual relative polarization efficiency with respect to the ground measurements average on all 353 GHz PSB detectors seen on single-bolometer HFI maps. The residuals after correction (in blue) are the HFI2018 ones, as shown in Fig. 35 of L03. The red relative residuals from SRoll2 are consistent with the 0.3% quoted above, thus confirming the nearly complete correction of the polarization efficiencies with a reduced dispersion of 0.6%, compatible with the error bars.

Fig. 19. Relative polarization efficiency with respect to ground-based measurements, extracted from single-bolometer maps for the 353 GHz PSBs. This figure compares the relative polarization efficiency from HFI2018 (in blue) with polarization efficiency cleaning (in red). The error bars are obtained by testing the same procedure on simulations.

4.3. Global focal plane rotation

We calibrate the effects of a possible global focal plane rotation by adding, to each simulation, the same constant polarization angle ψ to all bolometers at all frequencies. This generates E-to-B leakage, and we can then infer an absolute polarization angle correction as the one minimizing this leakage. We use the mean of the EB cross-spectra for ℓ ∈ [60, 167] as an indicator of the leakage, and we verify that it behaves linearly with ψ. This multipole range is a best compromise driven by the fact that, at lower multipoles, dust has unknown E-to-B leakage, and at higher multipoles, there is unknown beam leakage. We then produce two runs (one shifted by ψ₁ = −0.6°, and one not shifted with respect to the ground-based measurements ψ₀ = 0°) of 100 simulations, and we use a linear interpolation to find the value of ψ for which the leakage is zero.

Figure 20 shows the data and the dispersion of the simulations. A similar analysis has previously been performed (see LowEll) but without estimating the error from end-to-end simulations. It is nevertheless very important to go through the mapmaking process to take into account the other instrumental systematic effects that could generate leakage from E to B (e.g., beam errors).

Fig. 20. Distribution of the polarization absolute angle shift on 100 simulations without any shift in the input. The corresponding value for the data is overplotted (red line). There is no detection of any bias of the absolute angular shift in frequency average. The data value at 100 GHz shows 2σ detection, too small to be used as a correction.

The distribution of the simulations shows that there is no detectable focal-plane-rotation-angle bias introduced by SRoll2. Finally, this analysis demonstrates that the detected absolute polarization angle shift is consistent with the distribution of 100 simulations without shift, and therefore the SRoll2 maps do not have to be corrected for any absolute angle rotation.

4.4. ADCNL removal

The level of the ADCNL residuals cannot be checked only on the null test (see Sect. 3.3). To demonstrate the efficiency of the SRoll2 algorithm in removing the ADCNL systematic effect, we use the end-to-end simulations. We compare the SRoll1 approach when the ADCNL was corrected by a linear time varying gain, with the SRoll2 approach in which the full ADCNL of each detector is modeled by spline functions for each ring. The spline weights are fitted using the redundancy in the input maps.

Figure 21 shows the residual maps obtained taking the difference between input and output of an end-to-end simulation containing only dipole, sky signal, systematic effects, and ADC electronic noise. In the left column, SRoll1 is used to process the simulated data, while in the right one, we adopt the SRoll2 method. At CMB frequencies, the large-scale residuals present in the SRoll1 computation are much lower in the SRoll2 approach, leading to residuals below 0.5 μK outside the Galactic plane. At 353 GHz, the ADCNL effect is not dominant, and the 1/f noise dominates.

Fig. 21. SRoll1 (left column) and SRoll2 (right column) polarization intensity maps obtained by subtracting, from the output sky map, the input sky map from a simulation including a nonlinear ADC. The simulation input map contains dipole, sky signal, and systematic effects but electronic noise only to reveal effects lower than the full noise.

The central part of the Galactic plane is not improved as the signal gradient is too large and induces strong sub-pixel effects which dominate. The level of the residual is however very small (below 1%) compared to the Galactic plane signal.

The top panel of Fig. 22 shows the EE power spectrum of the SRoll2 residual maps presented in the right column of Fig. 21. At 217 GHz, the residual power left in the maps by SRoll2 is below the noise given by simplified FFP9 simulations (Planck Collaboration XII 2016) for ℓ >  10 but dominates at ℓ <  5. At 100 and 143 GHz, the ADCNL residuals are below the FFP9 noise down to ℓ = 3.

Fig. 22. EE pseudo spectra (auto-spectra in the top panel, cross-spectra in the bottom panel) evaluated for 100 GHz, 143 GHz, and 217 GHz on the simulations shown in Fig. 21 with a symmetric Galactic cut of 20°, 66% of the sky. The black solid line corresponds to a EE power spectrum with τ = 0.055 and 1010As = 21.14.

The bottom panel of Fig. 22 shows the ADCNL residuals for cross-power spectra which fall significantly below the auto spectra, showing that those residuals are weakly correlated between frequencies, and well below the fiducial cosmological signal. The much lower level of these cross spectra demonstrates that cross frequency spectra should be used for cosmology studies based on power spectra.

Figure 23 shows half-mission null test power spectra of five SRoll2 simulations with and without ADCNL. This figure shows that SRoll2 cleans the ADCNL effect in such a way that the levels of the power spectra of the residual maps are comparable for simulations including or not ADCNL. This demonstrates the ability of SRoll2 to clean the ADCNL residual at the detector noise level. The half-mission null test is no longer affected by the ADCNL residuals and any detected residuals in this null test are caused by systematic effects other than the ADCNL.

Fig. 23. Power spectrum of simulated half-mission null test. The plots in the left column are computed using SRoll1 while those in the right column are computed with SRoll2. The two lower frequencies are shown in the two cases when the ADCNL is included in the input and when it is not. These plots demonstrate the efficiency of the SRoll2 method to remove the ADCNL residuals while the previous method left some artefacts at very large scale.

In conclusion, this section shows that the ADCNL correction is very efficient and produces maps cleaned from the associated polarization leakage with residuals much lower than the noise amplitude, especially if we use cross-frequency spectra. This test confirms that large-scale residuals are now dominated by the 1/f detector noise.

4.5. All systematics

Finally, it is important to identify the main systematic effects that impact the data for each frequency. Figure 24, which can be compared to Fig. 52 of L03, shows the level of each systematic effect for each frequency. For all frequencies, at large scales, the 1/f detector noise becomes the main limitation, while this was not the case with previous versions of the Planck HFI data. A classical 1/f cleaning that uses offsets on timescales shorter than one minute has not been included in SRoll2 because we verified that it removes some signal due to insufficient redundancy in the Planck scanning strategy. Implementing this would require end-to-end simulations to assess the trade-off of reducing the 1/f noise versus the knowledge of the induced amplitude transfer function. For the three lower frequencies, the main systematic effect below the 1/f noise is the ADCNL residual. At 353 GHz the main systematic effect is the distortion of the dust emission by the polar efficiency and angle errors. Another very important improvement on the previous mapmaking is the very small level of leakage residuals, thanks to a better 353 GHz map due to the identification of the degeneracy between the gain determination and the very long time transfer functions.

Fig. 24. Polarization power spectra in , showing the noise (red) and the main systematic residuals: ADCNL remaining after cleaning (purple), polarization efficiency (dark blue), bandpass-mismatch leakage (light blue), leakage from calibration mismatch (orange), and the sum of all these (green). SRoll2 residual empirical transfer function (gray), and CO template subpixel effect (turquoise) have not been included in the sum. The fiducial CMB power spectrum is shown as a black lines. At 353 GHz, the systematic residuals are dominated by the polarization efficiency knowledge that has been fitted at the noise level for the lowest multipoles. At other frequencies, ADCNL is still the instrumental systematic effect that has the highest residuals.

This demonstrates that SRoll2 maps are the best that can be achieved before removing the 1/f detector noise. This cannot be carried out before the systematic effects at large scales can be brought down to a level where the degeneracies between 1/f noise and systematic effects do not remove signal. To improve the map at large angular scales the next generation of mapmaking should undertake an improved model reducing the dimensionality of the 1/f noise parameters to be computed and make it compatible with the information brought by the redundancy of the data.

Another limitation, that is especially acute at 353 GHz, lies in our ability to build statistically realistic simulations of the interstellar dust. The next generation of maps should integrate the component-separation tools inside the mapmaking and include the higher frequencies of 545 and 857 GHz (although they are not polarized) in order to provide the best dust foreground template. In that case, the frequency maps will be limited by the residual after removal of the foreground components which is expected to be dominated by the CMB signal after enough galactic masking (the central part of the galactic disk is probably not sufficiently sampled in frequency, or spatially, to be decomposed).

5. Conclusions

This paper presents the SRoll2 mapmaking method and the resulting frequency maps using the Planck HFI data. This method corrects or improves the correction of the instrumental systematic effects identified earlier and described in L03 but not fully cleaned in the Planck collaboration legacy release. The main goal of SRoll2 is to provide polarized maps for studies requiring the largest scale. This paper demonstrates and quantifies the efficiency of this method in correcting for systematic effects, and provides frequency maps at the detector noise level for all multipoles from 2 to 100 when masking the Galactic plane. The SRoll2 maps are free of systematic errors at large scale but, for higher multipoles, the maps need to be cleaned from other instrumental effects, such as the beam anisotropies or the 4-K line residuals, to reach a level lower than the noise.

We demonstrate that this new mapmaking approach provides polarized CMB maps with a variance at ℓ <  6 reduced by a factor larger than two with respect to L03, and without instrumental bias outside the Galactic disk. Furthermore, significant improvements are made for multipoles up to 100. This new data set can therefore be used to improve cosmological parameter measurements based on the re-ionization peak at very low multipoles (τ for E modes, and r and lensing for B-modes). SRoll2 maps (data and simulations) are available online⁶.

The future Planck HFI data mapmaking process should integrate component separation and systematic effect correction inside the mapmaking to avoid the degeneracies induced by the foreground bandpass leakage, for example.

Acknowledgments

This work is part of the Bware project, partially supported by CNES. The program was granted access to the HPC resources of CINES (http://www.cines.fr) under the allocation 2017-A0030410267 made by GENCI (http://www.genci.fr). The authors acknowledge the heritage of the Planck HFI consortium regarding data, software, knowledge, and the IAP hosted HFI DPC computing facility. LP acknowledges the support of the CNES postdoctoral program.

References

All Tables

All Figures

Fig. 1. Schematic view of SRoll2 iterations for building the Stokes I, Q, and U maps for one frequency band. The top box contains the initial improvements of the maps. Step#1 algorithm uses data HPRs to consistently compute the maps, the detector response, and the systematic effects correction parameters other than ADCNL. Step#2 computes the ADCNL correction that is propagated to the subsequent iteration. In the middle box, Step#1 produces maps used as input for Step#2 that computes the ADCNL correction. This whole process is then iterated 19 times. The bottom box is again a Step#1 to build the final maps.

In the text

	Fig. 2. Average on all bolometers of the response correction at each iteration of the 353 GHz map when using only the total dipole (in red) and using the sky model to fit the response (in blue).
In the text

	Fig. 3. Simulated ADCNL residuals (color coded) for the ADU level on the ADC (ordinate axis) along the time, expressed as ring number. Two bolometers are shown: left column for the 100-4b and right column for the 143-1a. First row: input of the simulation, second and third rows: their fitted 1D or 2D spline corrections without noise. Typical calibration is ≃6 μK/ADU.
In the text

	Fig. 4. EE power spectra for difference between two successive iterations of residual maps from a noiseless simulation including only the ADCNL effect and SRoll2 mapmaking for increasing number of iterations at 143 GHz. 20 iterations are enough to reach a good enough convergence.
In the text

Fig. 5. Noise and systematic residuals in EE, BB, and EB spectra, for difference maps of the half-mission null test binned by Δℓ=10. Spectra are rescaled to the full mission, and also rescaled to full sky from the unmasked sky fraction fsky = 0.43 used. Data spectra are represented by thick blue lines, and the average of simulations by thin black lines. The blue band shows one ±σ equivalent of the statistical distribution of the simulations. The linear y-axis scales are adapted to show the full data range on each panel.

In the text

	Fig. 6. Full sky, full mission SRoll2 maps at 100 to 353 GHz (in rows), for Stokes I, Q, and U (in columns).
In the text

	Fig. 7. HFI2018 minus SRoll2 difference maps at 100 to 353 GHz in Stokes Q (left column) and U (right column), smoothed by a 5° beam.
In the text

	Fig. 8. Comparison between 100 to 217 GHz HFI2018 and SRoll2Q and U maps, cleaned by the 353 GHz maps smoothed by a 5° beam.
In the text

	Fig. 9. Polarization intensity map difference between SRoll1 and SRoll2 algorithms at 100 and 143 GHz (in columns), applied to the data (first row) and to one realization of the end-to-end simulations (second row).
In the text

	Fig. 10. EE (upper 3 plots) and BB (lower 3 plots) power spectra of difference of null test maps for: half mission, detsets and odd-even rings, for HFI2018 (orange and blue) and SRoll2 (red and green). Data are overplotted on top of the simulations (average and dispersion).
In the text

	Fig. 11. Comparison between HFI2018 and SRoll2 100 × 143 B-modes variances, computed analyzing the respective E2E simulations.
In the text

Fig. 12. EE (left column) and BB (right column) power spectra of HFI2018 and SRoll2, based on half split data sets (detests, half-mission, and odd-even rings) maps at 100, 143, 217, and 353 GHz. Auto-spectra of the difference maps show the noise plus systematic effect residuals, rescaled to full mission equivalent and corrected for sky fraction used here (43%). The cross-spectra between the same maps are displayed and show the signal. The binning is: δℓ=1 for 2 ≤ ℓ <  30; δℓ=5 for 30 ≤ ℓ <  50; δℓ=10 for 50 ≤ ℓ <  160; δℓ=20 for 160 ≤ ℓ <  1000; and δℓ=100 for ℓ >  1000.

In the text

Fig. 13. Differences between two runs of SRoll2, one using the best available CO and dust templates, and another using a degraded version of these initial ones (as described in the text) designated as iteration 1. We display the foreground templates evolution: 12CO, 13CO (two top rows) and I, U, Q dust template (three lower rows). The first column is the difference between the best templates (iteration 0) and the run with the degraded templates (iteration1). SRoll2 is then run with foreground templates built from the sky maps and leakage coefficients from the iteration 1 displayed in the middle column. The results from a third iteration are displayed in the right column.

In the text

	Fig. 14. Power spectra of the initial templates and then differences between two successive iterations: first iteration in cyan, second in green, and third in red. The differences becomes quickly much smaller than the noise (black line) that is estimated from half mission difference computed with fsky = 0.43 and realigned at the full mission level.
In the text

	Fig. 15. Transfer function of the SRoll2 100 and 143 GHz maps. Thick line shows the average of ten simulations while the envelop shows the peak-to-peak dispersion.
In the text

	Fig. 16. Intensity map of the odd-even survey difference at 353 GHz from HFI2018 (top) and SRoll2 (bottom) after changing the model of the real part of the transfer function.
In the text

	Fig. 17. Residual solar dipole amplitudes with respect to the average dipole, for each 353 GHz PSB, shown per survey.
In the text

	Fig. 18. Difference of two single-bolometer intensity maps (353-3a and 353-5b) smoothed at 5°. Top: HFI2018; middle: SRoll2 without cleaning the long time constant; bottom: SRoll2 with cleaning the long time constant.
In the text

	Fig. 19. Relative polarization efficiency with respect to ground-based measurements, extracted from single-bolometer maps for the 353 GHz PSBs. This figure compares the relative polarization efficiency from HFI2018 (in blue) with polarization efficiency cleaning (in red). The error bars are obtained by testing the same procedure on simulations.
In the text

	Fig. 20. Distribution of the polarization absolute angle shift on 100 simulations without any shift in the input. The corresponding value for the data is overplotted (red line). There is no detection of any bias of the absolute angular shift in frequency average. The data value at 100 GHz shows 2σ detection, too small to be used as a correction.
In the text

	Fig. 21. SRoll1 (left column) and SRoll2 (right column) polarization intensity maps obtained by subtracting, from the output sky map, the input sky map from a simulation including a nonlinear ADC. The simulation input map contains dipole, sky signal, and systematic effects but electronic noise only to reveal effects lower than the full noise.
In the text

	Fig. 22. EE pseudo spectra (auto-spectra in the top panel, cross-spectra in the bottom panel) evaluated for 100 GHz, 143 GHz, and 217 GHz on the simulations shown in Fig. 21 with a symmetric Galactic cut of 20°, 66% of the sky. The black solid line corresponds to a EE power spectrum with τ = 0.055 and 1010As = 21.14.
In the text

Fig. 23. Power spectrum of simulated half-mission null test. The plots in the left column are computed using SRoll1 while those in the right column are computed with SRoll2. The two lower frequencies are shown in the two cases when the ADCNL is included in the input and when it is not. These plots demonstrate the efficiency of the SRoll2 method to remove the ADCNL residuals while the previous method left some artefacts at very large scale.

In the text

Fig. 24. Polarization power spectra in , showing the noise (red) and the main systematic residuals: ADCNL remaining after cleaning (purple), polarization efficiency (dark blue), bandpass-mismatch leakage (light blue), leakage from calibration mismatch (orange), and the sum of all these (green). SRoll2 residual empirical transfer function (gray), and CO template subpixel effect (turquoise) have not been included in the sum. The fiducial CMB power spectrum is shown as a black lines. At 353 GHz, the systematic residuals are dominated by the polarization efficiency knowledge that has been fitted at the noise level for the lowest multipoles. At other frequencies, ADCNL is still the instrumental systematic effect that has the highest residuals.

In the text