SRoll2: an improved mapmaking approach to reduce large-scale systematic effects in the Planck High Frequency Instrument legacy maps

This paper describes an improved mapmaking approach with respect to the one used for the Planck High Frequency Instrument 2018 Legacy release. The algorithm SRoll2 better corrects the known instrumental effects that still affected mostly the polarized large-angular-scale data by distorting the signal, and/or leaving residuals observable in null tests. The main systematic effect is the nonlinear response of the onboard analog-to-digital convertors that was cleaned in the Planck HFI Legacy release as an empirical time-varying linear detector chain response which is the first-order effect. The SRoll2 method fits the model parameters for higher-order effects and corrects the full distortion of the signal. The model parameters are fitted using the redundancies in the data by iteratively comparing the data and a model. The polarization efficiency uncertainties and associated errors have also been corrected based on the redundancies in the data and their residual levels characterized with simulations. This paper demonstrates the effectiveness of the method using end-to-end simulations, and provides a measure of the systematic effect residuals that now fall well below the detector noise level. Finally, this paper describes and characterizes the resulting SRoll2 frequency maps using the associated simulations that will eventually be released to the community.


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.

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 analogto-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 1 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. 1 In most of the cases, the amplitude of the difference in ADCNL residuals between successive rings is less than 10 −2 µK; in cases of solar flares, this difference can reach 1 µK for a few rings.  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.

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 freefree). 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.
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, ρ b cos(2φ b,r,p ), ρ b sin(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) 2 ; -ρ 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 ; -D tot r,p is the total CMB dipole signal (sum of solar and orbital dipoles), with D sol p being the template for the solar dipole with a constant direction and amplitude and D orb r,p being the template of the orbital dipole; -F dip b,r,p and F gal b,r,p 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 n bolo b=1 n ring r=1 O b,i = 0, 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,p 2 ; 2 (·) refers to the scalar product.
-ℵ 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: 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 M b,r,p becomes: We define the systematic effects b,r,p , including detector time constant amplitudes as: 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: Using Eqs. (1), (2), and (5), Eq. (3) now becomes: The term δg b,n g b should converge toward zero after several iterations. This is only possible if (A b,r,p · S p + D tot r,p + b,r,p ) 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 D tot r,p . The quantity (A b,r,p · S p + D tot r,p ) 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: where: -T b,r,p defined as the projection matrix A b,r,p A p A p −1 A p from frequency map to rings (we define A p as the pointing matrix built from all A b,r,p for the pixel p); -δg n g p A p · S p + D orb p is the matrix built from all weighted by the corresponding 1/σ b 2 within one pixel p; -M p , γ h,n V p,h , δL f,n C p, f , O p,n , and N p are the matrices respectively built from observations of pixel p by all g b,n M b,r,p , γ b,h,n V b,r,p,h , L b, f,n C p, f , O b,r,n , and N b,r,p , each being weighted by the corresponding 1/σ 2 b (e.g., O p ≡ b r O b,r,n,p∈r /σ 2 b ). 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 3 : withS tot r,p being a model for (A p · S p + D tot p ). The response value converges if In SRoll1,S r,p 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, theS tot r,p 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.
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 χ 2 : In Eq. (10), we introduce the g b,n terms in the weight at the denominator to compensate the g b,n N b,r,p term present in the numerator. Without this response term in the weights, the χ 2 minimization would be biased towards a smaller value of g b,n . The χ 2 minimization through a positive and symmetric matrix inversion is performed with a parallel conjugate gradient. This provides an estimate ofδg b,n , and alsoÕ b,r ,γ b,h , and δL b, f . Subsequently, the cleaned and calibrated signal M b,r,p is built using and a new mapS p can be computed from the clean and calibrated timeline: 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.

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: 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; -∆ b,r,p = 40 m=1 S b,r,p + E m + N m is the ADCNL bias at the TOI sampling frequency. In order to remove the ADCNL bias from the measured value M b,r,p , 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: 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.
We note that the function (Ẽ⊗N m ⊗ ) 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 4 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.

Fitting the ADCNL model parameters
An estimated residual signal model δM b,r,p is computed using the parameters and the cleaned and calibrated signal from the previous Step#1: 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: where: 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 A b,r,p is the model of the ADCNL for the bolometer b at the sample ring r for the corresponding A 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: 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 (A 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 χ 2 including a implicit binning in ADU due to the smoothness of the splines: and we obtain a new estimation ofÕ b,i,n , and an estimation of the totality of w t,b,s,n .Õ b,i 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.

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.

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.
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. 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.

SRoll2 frequency maps
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 A38, page 7 of 19 A&A 629, A38 (2019) estimated on the full mission and then build the half maps taking only the HPRs from odd (or even) surveys.

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 −4 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 −4 ) at 100 GHz; −0.62 µK (−1.7 × 10 −4 ) at 143 GHz; +0.79 µK (2.1 × 10 −4 ) at 217 GHz.
The negative average bias on the three CMB frequencies is −2 × 10 −4 , and reduces to −0.5 × 10 −4 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 −3 ) 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.
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 13 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. 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. 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. 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.

Power spectra of null test maps
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. 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 crossspectra 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.
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.
A38, page 12 of 19 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 4 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 13 CO emission is not simulated).

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 ( 12 CO, 13 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 12 CO and 13 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: The first run starts with the degraded templates. Then, in the following runs, the values of the leakage coefficients L b, 12 CO and L b, 13 CO 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 12 CO and the 13 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.
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.
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 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: 12 CO, 13 CO (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. 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 f sky = 0.43 and realigned at the full mission level. from different frequencies in the same mapmaking run with the caveat that each foreground component should be negligible in at least one frequency band.

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 −3 at = 3-10, and a few 10 −4 above. We also note a very small suppression of power at the level of few 10 −4 in the 143 GHz power spectra. We conclude that there is no evidence for any transfer function effect to be corrected for.

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.
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 5 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 −4 for HFI2018, and 5 × 10 −4 for SRoll2 data. This leads to an estimate of the uncertainty when averaging over all detectors of 2.2 × 10 −4 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. 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.
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.

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 T B 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: where -T ρ = Q p cos(2φ b,r,p ) + U p sin(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 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. 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.

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 ψ 1 = −0.6 • , and one not shifted with respect to the ground-based measurements ψ 0 = 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).
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.

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.
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. 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 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.
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.

All systematics
Finally, it is important to identify the main systematic effects that impact the data for each frequency. Figure 24, 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 10 10 A s = 21.14.
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 endto-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. 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).

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. . Polarization power spectra in D EE , 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.
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 6 .
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.