Issue 
A&A
Volume 534, October 2011



Article Number  A88  
Number of page(s)  19  
Section  Astronomical instrumentation  
DOI  https://doi.org/10.1051/00046361/201116871  
Published online  10 October 2011 
Iterative destriping and photometric calibration for PlanckHFI, polarized, multidetector mapmaking
^{1}
Laboratoire de l’Accélérateur Linéaire, Université ParisSud, CNRS/IN2P3, 91898 Orsay, France
email: tristram@lal.in2p3.fr
^{2}
Laboratoire Astroparticule et Cosmologie, Université Paris7 Denis Diderot, CNRS/IN2P3, Paris, France
Received: 11 March 2011
Accepted: 20 July 2011
We present an iterative scheme designed to recover calibrated I, Q, and U maps from PlanckHFI data using the orbital dipole due to the satellite motion with respect to the Solar System frame. It combines a map reconstruction, based on a destriping technique, juxtaposed with an absolute calibration algorithm. We evaluate systematic and statistical uncertainties incurred during both these steps with the help of realistic, Plancklike simulations containing CMB, foreground components and instrumental noise, and assess the accuracy of the sky map reconstruction by considering the maps of the residuals and their spectra. In particular, we discuss destriping residuals for polarization sensitive detectors similar to those of PlanckHFI under different noise hypotheses and show that these residuals are negligible (for intensity maps) or smaller than the white noise level (for Q and U Stokes maps), for ℓ > 50. We also demonstrate that the combined level of residuals of this scheme remains comparable to those of the destripingonly case except at very low ℓ where residuals from the calibration appear. For all the considered noise hypotheses, the relative calibration precision is on the order of a few 10^{4}, with a systematic bias of the same order of magnitude.
Key words: cosmology: observations / cosmic background radiation / techniques: photometric / methods: data analysis / techniques: polarimetric
© ESO, 2011
1. Introduction
The Planck mission^{1}, launched on May 14 2009, is a thirdgeneration satellite dedicated to observations of cosmic microwave background (CMB) anisotropies after COBE^{2} and WMAP^{3}. Its primary objectives are to measure the fullsky CMB anisotropies down to the cosmic variance limit reaching beyond ℓ ~ 2000 in temperature and ℓ ~ 1500 in Emode polarization. Other scientific goals include indepth studies of the Galactic emissions, extraction of catalogs of galaxy clusters and extragalactic sources and searches for nonGaussianity.
Planck observes the sky in nine frequency bands using two instruments: LFI with three channels centered at 30, 44, and 70 GHz and HFI with six channels at 100, 143, 217, 353, 545, and 857 GHz with an angular resolution in the HFI CMBdominated bands, 217 and 143 GHz, of from 5 to 10 arcmin. A mapmaking procedure designed to recover the sky maps from the noisy timeordered measurements is a necessary step down the data analysis pipeline of the Planck data needed to ensure a delivery of high quality sky images in the observed bands.
Mapmaking algorithms have been extensively tested and compared with each other within the Planck framework (see e.g. Ashdown et al. 2007a,b) assuming some standardized, and unavoidably, idealized circumstances. This study complements these earlier works describing an iterative procedure to combine mapmaking and CMB orbitaldipolebased calibration procedures and investigating in detail its precision and its dependence on lowfrequency correlated noise due to drifts and artifacts present in HFIlike timelines. The mapmaking technique used in this paper is destriping (e.g., Burigana et al. 1999; Delabrouille et al. 1998; Maino et al. 1999, 2002; Revenu et al. 2000; Keihanen et al. 2004), which is designed to mitigate the longterm correlation of the HFI bolometers and approach the white noise limit in the map domain. The algorithm first projects the timeordered data (TOD) for each stable pointing period of the Planck satellite separately onto the sky resulting in a set of overlapping ringlike maps. The relative offsets of these (one per ring), i.e. ones that lead to a coherent fullsky map, are then derived using a maximumlikelihood method that marginalizes over the sky signal. The sky maps are produced by projecting the offsetcorrected TOD onto the sky. Calibration is done iteratively using offsets determined from a previous step to constrain the next approximation of the calibration coefficients using the maximum likelihood algorithm. We also address the validity of the destriping hypothesis for different models of Fourier noise spectrum and estimate the effect of iteration on the calibration for the highprecision polarized mapmaking of Planck.
This paper is organized as follows. We first describe the algorithm for destriping and calibration in Sects. 2 and 3. Section 4 is devoted to Polkapix, the implementation of those algorithms (ringmaking, offset determination, calibration and projection) within the PlanckHFI Data Processing Center infrastructure. Simulations are detailed in Sect. 5. The results of destriping and calibration are described in Sects. 6 and 7, respectively. Finally, the results of the combined method solving for offsets and gain by iteration are discussed in Sect. 8.
2. Mapmaking
2.1. The mapmaking problem
PlanckHFI bolometers measure the brightness of the sky in a given direction convolved with an instrumental beam. The timeordered data vector, d, may therefore be modeled as the sum of the signal from the convolved sky T and the noise n(1)The pointing matrix A, of size N_{s} × N_{p}, relates each time sample s to the corresponding pixel p in the sky. It can typically be represented as very sparse. For an axisymmetric beam response, the “smearing” and “pointing” operations commute, and one can solve directly for the beamconvolved map, assuming that the matrix A contains only one (three for polarizationsensitive detectors) nonnull values in each row, as each sample is sensitive to only one pixel of the convolved sky. The data model, Eq. (1) can then be written as (2)where (I_{p},Q_{p},U_{p}) are beamsmoothed Stokes parameters in pixel p, ε is the polarization efficiency (equal to zero for the total power experiments), and ψ_{t} is the angle of the detector’s polarization direction, with respect to the polarization basis in the pixel p, at the time t.
For arbitrary beam patterns, one may either need to deal with nontrivial, positiondependent smoothing in the final maps, while retaining the simplicity of the matrix A, or to try to mitigate these effects with the help of one of the proposed beamdeconvolution techniques (e.g. Arnau & Saez 2000; Burigana & Sáez 2003; Armitage & Wandelt 2004), which typically require more complex pointing matrices.
After preprocessing of the TOD, noise can typically be considered as Gaussian and piecewise stationary so that all the statistical information of the noise is contained in its covariance matrix N = ⟨ nn^{T} ⟩ for which each stationary block is a symmetric Toeplitz matrix (Golub & van Loan 1996). Nevertheless, the instrumental noise is usually not white and contains some low frequency components that, if not accounted for properly, are projected to the sky as large stripes following the scan pattern.
2.2. Maximum likelihood (ML) methods
The most general solution to the mapmaking problem is obtained by maximizing the likelihood of the data given a noise model (Wright 1996; Tegmark 1997; Stompor et al. 2002). The solution, resulting in a minimum variance map, is given by the wellknown generalized least squares (GLS) equations, which provide an estimate of the sky signal and its covariance In practice, the inversion or even the calculation of the N_{p} × N_{p} matrix A^{T}N^{1}A is impossible for large datasets such as Planck, where the number of sky pixels, N_{p}, can be as large as many millions. Thus, Eq. (3) is usually solved using iterative methods such as the preconditioned conjugate gradient (PCG) (Golub & van Loan 1996) and capitalizes on fast Fourier transforms (FFT) to perform the Toeplitz matrixvector products (Doré et al. 2001; de Gasperis et al. 2005; Cantalupo et al. 2010).
2.3. Destriping methods
Socalled “destripers” attempt to simplify the general mapmaking problem described above for Plancklike scanning strategy when the solution would require too many resources. The destriping technique for the CMB mapmaking has been investigated in minute detail for the Planck satellite (Burigana et al. 1999; Delabrouille et al. 1998; Maino et al. 1999, 2002; Revenu et al. 2000; Keihanen et al. 2004; Stompor & White 2004). This technique takes advantage of the detectors observing large circles on the sky when the Planck telescope is spinning. Each circle is observed several times (of order 50) consecutively. Averaging over these circles enables us to compile datasets of higher signaltonoise ratio. It has been shown by Janssen et al. (1996) that the low frequency noise can be represented by a uniform offset on a given baseline, corresponding in the case of Planck to a stable pointing period called a ring. However, it can also be adapted to groundbased or balloonborne experiments (Patanchon et al. 2008; Sutton et al. 2010).
In the destriping approach, the noise is approximated as a lowfrequency component represented by the offsets x unfolded to timeordered data by the offset “pointing” matrix Γ and remaining Gaussian noise n characterized by a correlation matrix, N(5)The matrix Γ assigns to every consecutive piece of the time stream data the respective offset amplitude and its elements are therefore either 1 or 0. The maximumlikelihood estimates of the offset amplitudes, x, can be found from the timeordered data, d, by solving (6)where (7)and the covariance matrix for offsets reads (8)We note that, unlike Eq. (6), which provides an unbiased, though potentially nonoptimal, estimate of the offset amplitudes even if the assumed and actual properties of the instrumental noise differ, this last equation holds only if our data model given by Eq. (5) and the assumed noise correlation, N, are both correct. This is worth emphasizing as in the destriper applications one usually assumes that the remaining noise, n, is uncorrelated and typically piecewise white, and thus the matrix N is diagonal. This is indeed a key assumption from which the performance advantage of the destriper technique stems (Ashdown et al. 2009). Once the offsets, x, are estimated, we can calculate the sky map, by coadding samples into pixels after first subtracting the estimated offsets from the TODs, i.e., (9)Again such an estimate of the sky signal will be unbiased even if the simplifying assumptions about the noise, n, are adopted. If these are correct, or sufficiently accurate, then the error covariance matrix of the destriper map can be calculated as (10)Equations (8) and (10) can clearly provide at least some insights into the error correlation structure of the recovered maps and even be suitable for at least some of their statistical analyses, even if the reminder noise, n, is not strictly speaking uncorrelated and for this reason their calculation is included in the software described here and discussed later on.
As shown in Poutanen et al. (2006) and Ashdown et al. (2007b), the baseline length used to define the offsets need not be tied to the length of a Planck ring. Including priors on the low frequency noise, the destriping algorithm is equivalent to the GLS algorithm in the short baseline limit. In the Planck case, the duration of the stable pointing period is not constant over the time so that the baseline lengths vary slightly from 40 to 65 min. However, as we show later, this generalization of the Janssen et al. (1996) prescription does not affect the destriper performances.
3. Photometric calibration
The bolometer signal measured through Ibiasing is proportional to the small variation in the incoming power from the sky. To express the measurement in sky temperature units, one has to determine a gain per detector based on a known source in the sky. For low frequency channels (roughly between 20 and 300 GHz) for which the CMB signal is sufficiently high, the CMB dipole is usually used as a primary calibrator for experiments with large sky coverage. This is because it is only marginally affected by pointing errors and beam uncertainties, is a stronger signal than CMB anisotropies (by a factor 100) but not bright enough to cause nonlinearities in the detectors, and has the same electromagnetic spectrum. At higher frequency (typically above 300 GHz), calibration of data on Galactic signal, based on ancillary data (e.g. FIRAS maps Finkbeiner et al. 1999), is preferred. For smaller sky coverage, calibration is usually derived from point source objects of known flux.
The CMB dipole is induced by the Doppler effect of the relative motion of the satellite with respect to the last scattering surface. The solar system motion with respect to the last scattering surface (referred to as the solar dipole) is the dominant component of the satellite velocity. A residual contribution (called hereafter the orbital dipole) is induced by the motion of the satellite with respect to the solar system. The solar dipole can be considered as sky stationary during the observations and is thus projected onto the sky as an ℓ = 1 component with an amplitude measured by WMAP of 3.355 ± 0.008 mK (Hinshaw et al. 2005). Relativistic corrections to the solar dipole produce second order anisotropies at higher multipoles with amplitudes proportional to β^{ℓ} and more importantly couple both the dipole components as we discuss later on. Though the orbital dipole flux is typically one order of magnitude lower than the solar dipole, it is time dependent, and its timevariability is precisely determined by the satellite velocity.
The orbital dipole signal is modulated on a one year period. As it is not projected onto the sky, the ringordered data, d, for each sample s is modeled as (11)where A is the pointing matrix relating ring pixels to those of the sky, and t_{orbital} is the timedependent orbital dipole signal. The unknown are the sky temperature T (including the solar dipole) and the gain g.
As t_{orbital} is known, the calibration problem is linear and can be solved directly in the same way as the destriping problem. For instance, the maximumlikelihood estimate of the coefficient g can be obtained by marginalizing for T over the sky solving the equation (12)where (13)and subsequently used to estimate the sky map.
For Planck, using polarizationsensitive detectors to solve such a calibration problem for a single detector at a time turns out to be insufficient, as the mapmaking problem is either degenerate if polarization is explicitly modeled, or the answer biased, if polarization is neglected (see Sect. 7.1). However including multiple detectors in the calibration problem in Eq. (12) leads to a set of nonlinear equations. To avoid this problem, our algorithm performs only a single detector gain estimation at any given time, although it does this for a few detectors in parallel, the calibrated data of which are then used to determine offset estimates. The latter are subsequently used for the calibration estimates resulting in an iterative procedure.
4. Implementation
In the context of Planck data analysis, several destriper codes have been developed, such as MADAM (Keihanen et al. 2005; Keihänen et al. 2010) and Springtide (Ashdown et al. 2007a). In this paper, we present a modular code called Polkapix based on four separate steps: the ringmaking, the offset determination, the photometric calibration, and the projection of the rings onto a map. The offset determination and the calibration are solved iteratively. We note that explicit ringmaking is not necessary in destriper codes, although we use it not only as a convenient means of compressing the data but also as intermediate stage products that are useful for monitoring the systematics. In addition, and for the same reason, our implementation allows for internal, thus potentially different, sky resolutions for the offset and calibration determination independently of the one used for the final map projection.
4.1. Ringmaking
Rings are partial sky maps produced via a projection onto the sky of each single pointing period separately. These rings therefore provide a compressed and higher signaltonoise ratio rendition of the original TODs. The length of each stable pointing period varies between 40 and 65 min (with a mean of 45 min) in the Planck scanning strategy. During each period, noise is assumed to be white, with the low frequency part folded in as an overall ring offset, and consequently the rings are computed using a simple noiseweighted projection procedure. The resulting noise in the ring domain is thus (nearly) white, modulo the offset, and mostly uncorrelated from one ring to another.
For calibration purposes, as in Eq. (12), the orbital dipole has to be averaged in the same way for each pixel of the ring. This is done using the pointing direction and the satellite velocity calculated at each sample. As the dispersion of the observed direction falling into each ring pixel is very low (below 10^{3}), the averaging effects of pixelization are found to be negligible and no loss of accuracy is incurred as a result.
To avoid introducing any additional binning of the data, we choose a sky pixelization as a basis for this ring making (HEALPix, Górski et al. 2005). Our ring structure is therefore called HPR for HEALPix Rings in the following. As the ring are produced via simple projection assuming white noise only, samples that are flagged as bad, are simply not included in the ring making procedure, without any need for additional processing, such as gapfilling (e.g., Stompor et al. 2002).
4.2. Offsets determination
The destriping algorithm solves Eq. (6) for the offsets x (one offset per ring) taking as an input a set of HPR, d, computed at the ringmaking stage. We also impose the external constraint ∑ x = const. to break the degeneracy in the offset determination and set the arbitrary constant to 0 so that the mean value of the TODs (unphysical in any case) is conserved. Though we typically use the same underlying pixelization, HEALPix (Górski et al. 2005), for the ring and sky signals, these may have different resolutions for the intensity and polarization parts. The resolution should be high enough for the sky signal to be considered constant across the pixel, thus typically lower than the characteristic instrumental beam scale. However, a tradeoff should be made between the offsets uncertainty related to the level of noise in each pixel (which decreases for larger pixels) and the correlation between offsets (demanding smaller pixels). Moreover, to satisfy the first criterion above it is typically necessary to mask the inner part of the Galaxy, while estimating the offsets, as the galactic gradients are often too strong, leading to a poor estimate of the sky signal. These questions are discussed in more detail in Sect. 6.
Two solvers of Eq. (6) have been implemented here: the full inversion approach, allowing for, and requiring, the estimation of the offset covariance matrix and a (much faster) iterative method through conjugate gradient. The code has been implemented within the HFI Data Processing Center infrastructure, and is fully parallelized as far as both the workload and memory consumption are concerned. It is able to deal simultaneously with data from several detectors at the same frequency, as needed, for instance, to reconstruct polarization sky. As an example, the destriping is performed in a few minutes for 4 detectors on 32 processors using the full simulation described in Sect. 5. The reconstruction of the polarization for each pixel requires the inversion of each threebythree elements of the (A^{T}N^{1}A) matrix, some of which may be illconditioned for pixels with a insufficient number of bolometer orientations (see Ashdown et al. 2007a). We considered only pixels for which the condition number for this threebythree matrix is lower than 10^{3}.
The code described in this work was used to produce the Planck Early Release (Planck Collaboration 2011) and is one of the workhorses of the PlanckHFI Data Processing Centers. The implemented algorithm is overall very close to that of Springtide as discussed in detail in Ashdown et al. (2007a). We checked that indeed both these implementations produce the same offsets on a simple simulation of 4 HFI detectors at 143 GHz and including only the CMB signal (see Sect. 5 for details on simulations). The residual maps, i.e., a difference between noiseless input and output maps, give very similar rms with relative differences being: ~10^{6} for I and ~5 × 10^{5} for Q and U. The rms offsets differences agree within a relative error of ~5 × 10^{8}. For a detailed comparison of mapmaking methods, we refer the reader to the papers of the PlanckCTP Working Group (Ashdown et al. 2007a,b, 2009).
4.3. Dipole calibration
The calibration module solves Eq. (12) for gain g for each of the included bolometers individually. As an input, it takes precomputed rings corrected for the offset estimates, x, as derived earlier. The code is parallelized so that the gain estimation is achieved in fewer then two minutes on 16 processors. We use an underlying sky pixelization for intensity based on HEALPix and set its resolution as discussed in Sect. 7. The calibration error in the fitted gain is then estimated as (14)We note that this equation is an approximation neglecting the errors due to an offset determination. However, as we use multiple detectors for the latter, the approximation is expected to be very good.
The calibration algorithm takes advantage of the orbital dipole not being fixed on the sky, unlike the solar dipole. In practice, relativistic corrections couple solar and orbital dipoles creating a additional nonstationary signal. Though the latter is three orders of magnitude (β = 10^{3}) below the orbital dipole signal, we have to include the relativistic correction in the calibrator. We measure a relative bias of 6 × 10^{6} to the recovered gain when using the orbital dipole, instead of the nonstationary signal (with coupling corrections). Given the statistical error (~5 × 10^{5}, see Sect. 7.2), we conclude that the orbital dipole can be safely used as a calibrator for PlanckHFI data even if, for the rest of the study, we use the exact nonstationary calibrator constructed as the difference between the CMB dipole (including solar, orbital, and relativistic corrections) and the solar dipole.
The procedure may be susceptible to pixelization effects caused by signals with large intrapixel variations such as those found in the Galactic plane or for point sources. We evaluate the impact of this effect using a Galactic mask in Sect. 7.1. The calibration code can also provide a constraint on the gain due to the data corresponding to each sky pixel separately, as well as an estimate of its uncertainty. The overall calibration constraint is then a weighted average of the these pixelspecific values, with the weight being its estimated uncertainty. Thanks to this facility, we can identify the sky areas that contribute most significantly to the final result of the calibration procedure.
4.4. Rings projection
Once the offsets have been estimated, they are subtracted from each ring, which removes the low frequency noise component from the data. A “corrected” map is then obtained by simply coadding ringpixel amplitudes corresponding to the same sky pixels and weighted by the estimated ringpixel white noise estimate. We note that though the resolution of the final map does not have to be the same as that of the rings it should be typically no higher than that. In this paper, we always project the final fullsky maps on high resolution (nside = 1024), whereas we study the impact of degrading the resolution on the underlying maps in both destriping and calibration.
5. Simulations
Timeordered data (TOD) were simulated for 380 days of observations corresponding to slightly more than two sky surveys. The sampling frequency for HFI was set to 180 Hz, so that the total number of samples is about 5.6 × 10^{9} per bolometer. For 12 detectors and the addition of the pointing information, the total volume of the data amounts to 700 Gb per simulation. We restrict our discussion to four bolometers only per channel, which allows for a polarization reconstruction while limiting the overall data size to one more readily manageable. Nevertheless, we rescale the noise level so that noise in the final map is comparable to that expected for a complete Plancklike channel map. We focus on two CMB channels centered at 143 and 217 GHz. Simulations were done using the LevelS simulation codes (Reinecke et al. 2006) ported into the DPCHFI infrastructure, and in particular the Planck I/O library. Beam effects were simulated assuming Gaussian circular beams.
5.1. Scanning strategy
We used a Plancklike cycloidal scanning strategy. The satellite was assumed to be spinstabilized and rotates on its axis once per minute. The spin axis follows a cycloidal path on the sky by stepwise displacements of approximately 2 arcmin every 48 min. The stable pointing period between two repointings is not constant in order to take into account the spin axis not remaining in the Ecliptic plane but following a cycloid path (Fig. 1). Between the repointings, the spin axis nutates with a mean amplitude of 2 arcsec. Detectors scan the sky following almost great circles as they point at 83.8° away from the spin axis (Fig. 2). We also include small variations in the spin rate (rms 2.16 arcsec/s).
For the polarization studies, we simulated four bolometers at 143 GHz coming from two horns of the PlanckHFI focal plane (1431 and 1433). Each horn groups two detectors (labeled a and b) sharing the same pointing with a polarized grid orientation rotated by 90°. The two horns follow the same path on the sky but are 2° apart from each other. Their relative orientation of the polarized grids in each horn differ by 45°. During the 12.5 months of the mission, the four detectors observe each pixel on the sky with multiple polarization orientation in a 3.44 arcmin resolution map (HEALPix nside = 1024). Therefore we are able to determine the three Stokes components in each pixel. In addition, we simulated one total power, 217 GHz bolometer.
Fig. 1 Dwelling time of the simulated stable pointing periods for the 12.5 monthlong observation. The variation in the period length reflects that the spin axis does not remain in the Ecliptic plane but follows a cycloid path. 
Fig. 2 A map of hit counts per 3.4 arcmin pixels (Nside = 1024) in Galactic coordinates corresponding to 12.5 months of observation for four detectors. The sky is covered at least twice. Redundancy of the revisits is maximal around the Ecliptic poles. The conspicuous sshaped band with a higher hit count corresponds to the beginning of the third survey. 
5.2. Sky signals
The simulated signal combines the CMB and diffuse Galactic foregrounds from the Planck Sky Model^{4} at a resolution of 1.7 arcmin (HEALPix nside = 2048). The CMB anisotropies are generated using CMBFAST (Seljak & Zaldarriaga 1996) with the WMAP fiveyear bestfit model without lensing (Ω_{m} = 0.257,Ω_{b} = 0.044,Ω_{Λ} = 0.743,h = 0.719,τ = 0.1,n_{s} = 0.963,r = 0.02). Both solar and orbital dipoles are included. The Galactic foreground signals include thermal and spinning dust, synchrotron radiation, and freefree scattering in intensity. Polarized signals for thermal dust and synchrotron are also added. More details about the sky model used for this work can be found in Leach et al. (2008) and Betoule et al. (2009).
The maps were smoothed using a symmetric Gaussian beam with FWHM of 7.05 and 4.72 arcmin at 143 and 217 GHz, respectively. TODs were then generated using LevelS codes (Reinecke et al. 2006) extracting the signal from the input sky maps according to the scan strategy. No bolometer time constant was included.
5.3. Noise properties
Streams of 1/f^{α}noise were generated at the TOD level using the algorithm described in Plaszczynski (2007). We selected one of two values of the slope (α = 1 and 2) and two knee frequencies (f_{knee} = 0.1 and 0.01 Hz), with a cutoff frequency f_{min} = 10^{5} Hz below which noise becomes white. In this approach, a white Gaussian noise realization is generated and the long range correlations are obtained through a linear digital filter. For a given seed of the generator, one can directly study the effect of the pure 1/f^{α}noise component using the difference of reconstructed maps with and without the correlation turned on. For Monte Carlo simulations, we need a significant number of noise realizations, which renders impractical the writing of the TOD to disk. We therefore directly project the noise at HPR level using a precomputed HPR pointing matrix.
The thermal fluctuation expected for PlanckHFI detectors are usually modeled by a 1/f^{2} Fourier spectrum with a low knee frequency. The noise power spectra of HFI bolometers are more accurately described by a 1/f spectrum (α = 1) with f_{knee} ≃ 0.1 Hz plus white noise (Planck Collaboration 2011). We set the TOD white noise level in order to match the expected noise level of the foreseen PlanckHFI, 143 GHz frequency map. Thus, the noise rms per sample was set to 589 μK in thermodynamic CMB scale (Planck Collaboration 2005). The noise covariance matrix in the map domain is given by Eq. (4). If correlated errors are small, the covariance matrix for (I, Q, U) at pixel p reads (15)where the sums extend over all samples of all considered detectors falling into the pixel p, and we have assumed that the polarization efficiency, Eq. (2), is the same for all detectors. The second step results from the cancellation of the offdiagonal contributions coming from two detectors of each horn, which thus have their polarizers rotated by 90 degrees with respect to each other. As the polarization efficiencies are not all the same and the polarizers are not rotated by exactly 90 degrees, the second step is only an approximation with the offdiagonal elements being on the order of and , where Δϵ denotes a difference between the polarization efficiencies of two detectors of the same horn and Δα the deviation from orthogonality. The offdiagonal terms are further suppressed for sky pixels, which are observed by the same horn multiple times with different attack angles (the correlation is lower than 5% for our scanning strategy). Consequently, the matrices M_{p} are usually diagonally dominated for sufficiently wellobserved pixels, with the diagonal elements being approximately given by Eq. (15). Given the hit counts per pixel n_{p} (Fig. 2) and the characteristics of the HFI detectors that we used (ε = 0.83, 0.85, 0.84, 0.90), the averaged white noise per pixel we expect in our simulations is therefore 15.4 μK_{CMB} for I and 25.6 μK_{CMB} for Q and U maps at 143 GHz. In the following, we compare the level of residuals with those values.
6. Destriping parameters and performances
6.1. Systematic studies
We now investigate the various systematic effects possibly affecting the offset determination. We use a simulated sky with CMB and Galactic signal at 143 GHz for four bolometers, allowing for intensity and polarization reconstruction. For each bolometer, we generate 100 pure white noise realizations. HPR are constructed at high resolution (nside = 2048) and the final fullsky map is projected at nside = 1024. We first built the signal map using coaddition in each pixel. We also compiled the white noise map for each realization without destriping. For each noise realization, we constructed the residual map subtracting the signal+white noise map from the destriped map. We then computed fullsky angular power spectra for each residual map using the Xpol (Tristram 2006) code, based on the S^{2}HAT library^{5} (Hupca et al. 2010; Szydlarski et al. 2011) and averaged them over the simulated realizations.
First we investigate the effect of large signal variations within the pixels on the offset determination. These pixels can be found mostly in the Galactic plane so we remove the pixels with the highest Galactic signal using a mask (based on FIRAS intensity maps) in the offset determination. We gradually increase the maskedout sky fraction from 0% to 25%. We use the destriping module to estimate the offsets using an internal sky resolution of 6.9 arcmin (HEALPix nside = 512). The results are shown in Fig. 3. We can see that only temperature residuals are sensitive to the masking and that a 5% mask is sufficient to suppress pixel effects due to strong subpixel signal gradients caused by the Galactic signals. Consequently for the following, we use an internal 5% Galactic mask at 143 GHz for offset determination.
Subsequently, we assess the impact of the assumed sky resolution on the offset determination procedure. We do this in two steps. First, we neglect the polarization signal while determining the offsets, and consider several internal resolutions from nside = 64 to 2048 in the offset determination. Figure 4 shows the power spectra of the destriping residuals for a pure white noise. In terms of temperature, at low internal resolutions (nside lower than 256), residuals due to the pixelization effect dominate at low multipoles over the instrumental noise level. Residual power spectra in polarization are unaffected by changes in the internal resolution for the temperature and are consistently higher at lowℓs than the noise level. We attribute this last observation to no polarization being included during the offset computation step.
Fig. 3 Noise residual fullsky power spectra with various internal Galactic mask applied to the offset determination. Residual maps are at nside = 1024. Spectra are smoothed over five multipoles. Signal spectra are plotted in red. As a reference, we also plot the CMB WMAP fiveyear fiducial model (black) and a total power of the simulated sky signal, CMB+Galaxy+noise, (red curve). From left to right: temperature, Emode, and Bmode. 
Fig. 4 Noise residual fullsky power spectra with various internal resolutions (nside = 64, 128, 256, 512, 1024) and neglecting polarization in the offset determination. Residual maps are at nside = 1024. Spectra are smoothed over five multipoles. Signal spectra are plotted in red. As a reference, we also plot the CMB WMAP5yr fiducial model (black). The horizontal, dashed lines show the predicted noise level. From left to right: temperature, Emode, and Bmode. 
Finally, we compare the destriping performances when the polarization is included (for various internal sky resolutions) in the offset determination. For this comparison, the internal resolution for intensity is fixed at 3.4 arcmin (HEALPix nside = 1024). As shown in Fig. 5, including the polarization reduces the level of the residuals of the E and B mode polarized spectra. For polarization, at high resolution (nside greater than 512), we found residuals no larger than the white noise level over the entire ℓ range, including the low multipole range (ℓ < 20). At higher multipoles, a higher resolution seems to produce a slight increase in the residuals but we understand this is negligible given the white noise level. As far as the intensity is concerned, a low resolution for polarization induces strong residuals in the intensity map for multipoles lower than ℓ = 100.
Fig. 5 Noise residual fullsky power spectra when using different internal sky resolution for polarization in the offset determination. Internal resolution for intensity is fixed to nside = 1024. Residual maps are at nside = 1024. Spectra are smoothed over five multipoles. Signal spectra are plotted in red. As a reference, we also plot the CMB WMAP fiveyear fiducial model (black). From left to right: temperature, Emode, and Bmode. 
From these studies, we conclude that polarization must be taken into account during the destriping, even if this requires more data (for Planck at least four bolometers) leading to an increase in the CPU time as well as the introduction of possible systematics due to detector mismatches (not discussed here). Moreover, we have demonstrated that a resolution of 6.9 arcmin (HEALPix nside = 512) for both intensity and polarization ensures that the residuals remain negligible relative to white noise level.
6.2. Statistical uncertainties and offset covariance matrix
The statistical uncertainties in the offset determination are quantified by the offsetcovariance matrix, , which in general depends not only on the scanning strategy but also instrumental characteristics, both of which enter via the pointing matrices A and Γ, Eq. (8), as well as noise properties of the detectors encoded in the matrix N. However, this is typically the scan pattern, which plays a dominant role in determining the overall structure of the noise covariance (e.g. Stompor & White 2004; Efstathiou 2007), and conversely most of the major features of the covariance can be traced back to the details of the scanning. Below we discuss the results of the numerical computations of , examples of which are displayed in Figs. 6–9, emphasizing this connection. To this end, we first rewrite Eq. (8) more explicitly as (16)where, is the sky error covariance, given in general by Eq. (10), and which for our purpose in this section will be approximated as either diagonal or blockdiagonal for the unpolarized and polarized detectors respectively as given by the first term on the righthand side (rhs) of Eq. (10).
The first term on the rhs of Eq. (16) describes the offset variance derived assuming that the sky signal is known. The effect of the marginalization over the unknown sky performed in the presence of the crossings between the rings is then quantified by the second term. We emphasize that no correction is explicitly applied in the above equation to account for the singularity of the problem (see Sect. 4.2). Though this is sufficient for the qualitative discussion presented here, it does have to be treated in the actual numerical work. In the following, we discuss separately two specific cases of interest: one assuming unpolarized total power detectors and the other polarizationsensitive ones. For simplicity, we assume the noise variance per sample, , to be timeindependent and the same for all considered detectors.
Fig. 6 Left: diagonal of the offset covariance matrix for various internal sky resolutions used for the offset recovery. Black line shows the first term of Eq. (17). Right: offset covariance matrix profiles centered on the diagonal showing the ringtoring covariance for various internal sky resolution. 
Fig. 7 Offset covariance matrix without (left), with a Galactic mask (middle) and the difference (right). 
Fig. 8 Binned offset covariance matrices. left: Polkapix estimate, right: covariance for 10 000 MC. The bin width is 100 rings. 
Fig. 9 Offset correlation matrix for four polarizationsensitive bolometers at 143 GHz without (left) and with (right) polarization reconstruction in the destriping. 
6.2.1. Total power detectors
For total power detectors, the entries of the sky pointing matrix, A, are either 0 or 1, with the indices of the nonzero elements, (t, p), defining which pixel p has been observed at the time, t. The matrix, C_{m} is approximately diagonal with the elements equal to the estimated noise level for each pixel. Elements of the matrix, Γ^{t} N^{1} A are then given by a number of times each pixel was observed within each offset baseline, i.e., a ring, weighted by the sample noise. In this special case, we can therefore rewrite Eq. (16) as (17)Hereafter we use n_{t}, n_{r} to denote the total number of samples and rings per detector; n_{p} the total number of pixels observed by all detectors; n_{obs}(p,r,d) the number of observations of a pixel p within a pointing period (ring), r, for a detector, d, where n_{obs}(p) ≡ ∑ _{r,d} n_{obs}(p,r,d) is the total number of its observations; and n_{obs}(r,d) the number of samples in a pointing period r. Moreover, r and r′ denote two rings corresponding to the pointing periods of two detectors, d and d′, respectively, and a symbol p ∈ r ∩ r′ refers to the pixels observed during both pointing periods, r and r′. We note that the second term on the rhs is always nonnegative and vanishes if r ∩ r′ = ∅. It therefore increases the offset determination uncertainty and leads to a correlation, rather than an anticorrelation, of the offsets. This reflects that the offset determination precision is always slightly lower whenever the sky signal has to be estimated simultaneously. This has an easy intuitive explanation: if one of two offsets corresponding to two rings, which cross on the sky, turns out to be overestimated, the sky signal in pixels common to both rings will be on average underestimated, and the offset of the second ring will tend to compensate for this and therefore be overestimated as well. Nevertheless, the anticorrelations are unavoidably present in the offset correlation matrices. This can be seen in Figs. 7 and 8 and arises because of the condition we impose on the sum of all offsets in order to break the problem singularity. This effectively requires that the sum of all the elements in each matrix column or row vanishes, hence that at least some, or more typically many, elements of the correlation matrix are indeed negative.
We first consider the case with r = r′ and d = d′. Equation (17) then reads (18)Given that typically n_{obs}(r,d) ~ n_{t}/n_{r} ≫ 1, the second term in the parenthesis is strongly suppressed. However, the magnitude of the sum over the pixels observed in the period, r, can be large enough to compensate for this, rendering this term important. To illustrate this, we first observe that (19)and consequently that the sum in Eq. (18) is larger, whenever varies more rapidly from a pixel to a pixel along the ring, or, in other words, whenever the distribution of the observation within the pointing period, r, is more inhomogeneous. In contrast, whenever both n_{obs}(p,r,d) and n_{obs}(p) are homogeneous, thus pixelindependent, the second term is given as (20)where n_{p}(r,d) is the number of pixels observed in the period, r. The offset covariance is then close to its minimum, where it is dominated by the first term, with the second term contribution nearly negligible thanks to an implicitly implied perfectly linked network of rings on the sky.
At the other extreme, all samples of the ring r are located in a single pixel, p. The second term is then roughly given by (21)as n_{obs}(p) ≥ n_{obs}(p,r,d) = n_{t}/n_{r} and the equality is realized whenever p is observed only within the considered pointing period. We thus conclude that the second term, though typically smaller, may be comparable in magnitude to the first one. The important consequence of these considerations is that whenever the inhomogeneity of pixel observations for some pointing period increases the diagonal elements of the offsetoffset covariance matrix tend to increase. This explains the results displayed in Fig. 6, which display a noticeable increase in the offset variances as a result of a decrease in pixel size, which in turn enhances the observation inhomogeneity. This has important practical implications when selecting an appropriate pixel size during the offset determination step (see Sect. 6.1). We note that for large pixels the second term is unimportant to the observation distribution among the pixels reaching the limit as discussed earlier in the context of Eq. (20).
For two different rings, i.e., r ≠ r′ or d ≠ d′, the first term in the expression for the covariance is absent and therefore the latter is determined by the second term. Its magnitude is now driven by two effects, the first one – as in the previous case – related to the observations inhomogeneity, which in general tends to increase the strength of the correlations, and the second one related to the number of pixels in common between the two rings. The latter effect is however the one that drives the largescale appearance of the offset correlation matrix, Fig. 7. Consequently, the matrix is diagonal dominated with the correlations decaying away from the diagonal as the number of crossing pixels between the two rings decreases. The decay rate as well as the autocorrelation change with the ring number as does the Planck scan pattern, during which the satellite spin axis follows a cycloid in the Ecliptic coordinates. Consequently, the width of the diagonal band is found to vary along the diagonal. As Planck observes a full sky in six months of observations, with the end of the six month corresponding roughly to ring 5000, the overlaps between the rings separated by ~5000 are again enhanced leading to a significant level of the corresponding correlations. This again happens for the ring separation of on the order of 10 000, when the third fullsky survey starts. We note that the pattern of sky scanning changes between the first and second survey, owing to the different satellite tilt at which it is preformed, and consequently the offset correlations at this lag are somewhat more smeared out than those around the main diagonal. As the third scan follows the first nearly exactly, the correlations at the lag of ~10 000 are expected to resemble those at the diagonal very closely as indeed already seen to some extent in the figure.
Owing to the two competing effects mentioned above, the impact of the pixel size on the offdiagonal correlations is more complex than on the offset autocorrelations as they tend to change the correlation amplitude in an opposite way. Hence, though the correlation length will decrease when any pixel is shared on average by fewer rings, the changes in the correlation amplitude will depend on the fine details of the scanning strategy and in general will be neither simple nor easy to predict. This is indeed confirmed by the tests presented in the right panel of Fig. 6.
The contribution of the second term in Eq. (17) is sensitive to the pixels that are retained for the offset determination. This is true in particular when only a few pixels are shared between two rings, or whenever wellobserved pixels are removed. In this context, the pixel selection is vitally important to practical applications as discussed in more detail in Sect. 6.1. Here we only illustrate it in Fig. 7 by contrasting two covariance matrices, obtained with all sky pixels included (left panel) or with those corresponding to the discarded Galactic plane (middle).
For two rings that are different on the sky, whether they correspond to pointing periods of a single or two different detectors is completely irrelevant. This is because of our assumptions about lack of timedomain noise correlations and the same value of the noise variance. This, together with all the detectors being assumed to scan the sky in the same way, explains the periodic structure of the four detector offset covariance shown in Fig. 9. We indeed see that all the offdiagonal blocks are essentially identical, while the only difference between any of the diagonal and any of the offdiagonal blocks is caused by their diagonal values and is a result of the nonvanishing contribution to the first term in the case of the diagonal blocks. Our assumptions about the scanning imply that n_{obs}(p,r,d) is independent of the detector number, rendering the offdiagonal blocks symmetric with respect to their diagonal.
We also checked our calculation of the offset covariance matrix from Polkapix using 10 000 MonteCarlo simulations of pure white noise. We simulated a simplified sky with CMB and Galactic emissions for a single PlanckHFI detector at 217 GHz. We then generated 10 000 white noise realizations. For each simulation, we estimated the offsets using a 5% Galactic mask. We found a very good agreement between the MC standard deviation and the covariance matrix as shown in Fig. 8.
6.2.2. Polarizationsensitive detectors
For polarizationsensitive detectors, the pointing matrix, , is more complex, Eq. (2), with elements, that can be both positive and negative. The sky covariance can be approximated, as before, by the first term of Eq. (10). This would in general render a blockdiagonal structure with each block of a size 3 × 3 characterizing singlepixel correlations of three Stokes parameters. However, in our case (see Sect. 5 and in particular Eq. (15)) where the data of four detectors are analyzed simultaneously, it can again be treated as diagonal with the diagonal elements corresponding to Q and U parameters always being a factor 2/ε^{2} different than those for the total intensity, I. Hereafter, for brevity we omit the polarization efficiency ε. Assuming as before that the white noise is stationary over the course of the entire set of observations for which an rms value, , is identical for all detectors, we can now specialize Eq. (16) to the case at hand obtaining (22)Here, (23)We note that for each ring, r, and each detector, d, the averages above are done separately and therefore the arguments used in simplifying Eq. (15) do not apply here and the terms with cosines and sines do not vanish trivially in the Planck case. In contrast, they have an important and easily discernibly impact on the offset covariance matrix as shown in Fig. 9.
For Planck scans, the angle between the polarizer orientation and the celestial coordinate system for a given pixel p and during a ring r does not change significantly from one to another pixel crossing as they are fixed with respect to the scanning direction. The averages of the trigonometric functions are therefore well approximated by the cosines and sines of this specific angle. Equation (22) can then be simplified yielding (24)Consequently, for any single ring, i.e., r = r′ and d = d′, the extra terms due to polarization are always positive, and, twice as large as the third, angleindependent term on the rhs. For the two corresponding rings belonging to two detectors in the same horn, e.g., with the relative polarizer orientation at 90 degrees, the third term is equal to − 2 forcing the corresponding element of the offset covariance to be negative. If the two detectors above are chosen from two different horns, i.e., their polarizers are either at 45 or 135 degrees with respect to each other, the corresponding rings of the detectors are again positively correlated. Moreover, in the latter case and for two different rings, r and r′ the cosine function in Eq. (24) turns into a sine and the contribution is asymmetric with respect to an exchange of r ↔ r′ giving rise to an asymmetry of the corresponding offdiagonal blocks. This is not so for the same horn detector blocks as the correlations are given by a cosine function, symmetric with respect to ring exchanges. All these observations and in particular a dichotomy between identical and different horn detectors are indeed confirmed by our numerical calculations, Fig. 9.
As in the unpolarized case discussed earlier, we can track the pattern of each block to the specific, singledetector scan features, such as the beginning of the second (r ~ 5000) and the third (r ~ 10 000) full sky survey, leading to an overlap between their respective rings with those of the previous surveys and thus to an enhancement in the strength of the correlations.
6.2.3. Offset covariance – a recap
We briefly summarize all salient points of the discussion presented above. We have found that in the Planck case the offset covariance matrix is in general nondiagonal and has a rather complex structure, major features of which are due to the pattern of the Planck scan. The computation of this matrix requires rather lengthy and involved numerical calculations, although quick intuitive insights can be obtained in a semianalytic way (see also, e.g., Stompor & White 2004; Efstathiou 2007). Those can be not only useful in devising new, efficient scanning strategies but also provide a useful crosscheck of the numerical results. We have also tested the latter by means of extensive Monte Carlo simulations and found good agreement.
We have demonstrated that the offset variance and their correlations depend on the observation (in)homogeneity and thus on the assumed pixel size and the number of crossing points between the rings, i.e., sky areas covered by the experiment within a single pointing period. We have also found that the overall structure of the covariance depends on whether the detectors are polarizationsensitive, an effect that, for a Plancklike scanning, has to be accounted for if a highprecision offset determination is desired. This last conclusion seconds one of the conclusions obtained in Sect. 6.1.
Fig. 10 Residuals spectra after destriping the signal plus noise in the three cases of (1) pure white noise (blue), low frequency noise plus white noise with (α = 2, f_{knee} = 0.1 Hz), (2) (α = 2, f_{knee} = 0.01 Hz) and (α = 1, f_{knee} = 0.1 Hz) compared to the simulation power spectra (red), and (3) CMB WMAP fiveyear fiducial model (black) on 100 simulations. Upper line, from left to right: temperature, Emode, and Bmode. Bottom line, from left to right: crosscorrelation TE, TB, and EB. 
6.3. Residuals from low frequency noise
In the previous sections, we have described the systematic effects related to the destriping algorithm presented in Sect. 4. In particular, we have shown that the pixelization effects linked to the strong variation in the sky signal within a pixel are negligible for the offset estimation. However, we have found that a significant level of residuals arise whenever the polarization signal is neglected or the resolution adopted for the estimation of the underlying sky signal is reduced. In this section, we study the residuals for various levels of correlated noise. For each bolometer, we add to the simulated signal (containing CMB and Galactic emissions) 100 noise realizations for three cases: pure white noise, correlated noise with α = 2 for two knee frequencies (f_{knee} = 0.1 and 0.01 Hz), and α = 1 for a knee frequency of 0.1 Hz. We reduce the systematic bias using a 5% internal Galactic mask and an internal sky resolution of 6.9 arcmin (HEALPix nside = 512), as shown in Sect. 6.1.
Figure 10 shows the power spectra of the residuals computed for the three selected cases in both temperature and polarization. The level of the residuals is determined by the level of lowfrequency noise introduced into the TOD. Thus, residuals increase with both the knee frequency f_{knee} and the slope α. We show that lowfrequency noise residuals are negligible for f_{knee} ≤ 0.01 Hz when α = 2. For a higher f_{knee}, a residual at low frequency biased the spectrum from a pure whitenoise power spectrum at low multipoles (ℓ < 50). In the case of α = 1, the residuals are lower by one order of magnitude.
Apart from the small deviations coming from the realistic scanning strategy used in this study, we found results that are fully compatible with the analytical predictions of Efstathiou (2007). In particular, we found that the power spectra of the residuals scale like σ^{2} (where σ is the white noise level in map domain) in a similar way to the noise power spectra itself. This implies that the results discussed in this section do not depend on the value of the white noise level but only on the characteristics of the correlated noise (through the f_{knee} and α parameters of the noise Fourier spectrum).
7. Calibration results
7.1. Systematics
The calibration method is based on the assumption that when observing the same pixel of the sky at different times, the only variation in signal is due to the orbital dipole (assuming that the lowfrequency noise is perfectly handled previously by the destriping). Any other effect leading to a difference between two different measurements of the sky power in one pixel, e.g. due to polarization because of a change in the bolometer orientation or a large signal variation inside the pixel) will result in a bias in the gain reconstruction. We performed a simulation including only CMB and Galactic emission intensities and estimate the gain for several internal resolutions and various internal mask sizes. The masks are based on the gradient of the Galactic emission, which is the component inducing the strongest intrapixel variations. We have found a bias smaller than 10^{5} for all considered internal resolutions as soon as the 1% most variable Galactic plane areas were masked. As shown in Fig. 11, the gain is biased only for the lowest resolution when the Galaxy is not masked at all (g ≤ 6 × 10^{5}). We conclude that in the calibration case, intrapixel variations induce very small systematics. To match the destriping settings, we fix hereafter the internal resolution to nside = 512.
Fig. 11 Evolution of the calibration bias due to the sky resolution with the masked fraction. Error bars are statistical estimates from Polkapix. 
We also evaluate the bias in the gain reconstruction caused by the polarized sky signals from the Galaxy and CMB. We study three cases: CMB including intensity and polarization, Galactic emission with only polarized signal, and CMB + Galaxy (intensity and polarization). We reconstruct the gain applying several internal masks based on the intensity of polarization of the Galactic component. Figure 12 shows the bias of the reconstructed gain for the three simulations using each mask. We first note that the CMB polarization induces a bias (typically ~5 × 10^{5}) that is constant with respect to the percentage of sky masked but depends on the CMB realization and the bolometer characteristics. For a polarized Galactic signal, we find a stronger bias of up to 4 × 10^{4}, which decreases with the mask size, reaching zero for a 20% masked fraction. The statistical error bar increases with the size of the mask because of the smaller sky coverage.
Fig. 12 Evolution of the calibration bias due to the polarized signal and the masked fraction (blue: simulation of CMB; orange: simulation of Galactic emissions; red: simulation of CMB + Galaxy). Error bars are statistical estimates from Polkapix. 
Bias in the gain reconstruction leads to a misestimation of the amplitude of the recovered intensity map. Its effects are however more dramatic for the polarization. If the bias affects all the detectors at the selected frequency in the same way, in the case of PlanckHFI, this would not result in a leakage of the total intensity into the polarization. However, in our case, the bias depends on the detector and it will cause a leakage of the temperature signal (dominated by the solar dipole) into Q and U Stokes parameters. Figure 13 illustrates the effect of this leakage on residual Q and U maps for a gain error of 5 × 10^{5}.
Fig. 13 Residual maps of the Q and U Stokes parameters for a calibration uncertainty of 5 × 10^{5}. The main features on large scales are related to the CMB dipole leaking into the polarization. 
7.2. Gain statistical error
The uncertainty in the gain reconstruction is estimated from Eq. (14). We have verified this estimation using 1000 Monte Carlo simulations including only the CMB dipoles together with pure white noise for one bolometer at 143 GHz. We apply the 20% Galactic mask used in the previous section. Figure 14 shows the distribution of the gain bias (the relative difference between the reconstructed and simulated gains) for the 1000 simulations compared to a Gaussian with FWHM derived from Eq. (14). The standard deviation of the MC is fully compatible with the analytical error bar given by the number of simulations. The statistical error does not depend on the input signal or noise characteristics but only on the white noise level and the orbital dipole signal. For the four bolometers at 143 GHz, statistical errors are ~4.5 × 10^{5}.
Fig. 14 Gain bias distribution for 1000 simulations including CMB dipoles and pure white noise relative to Polkapix posterior PDF. A mask covering 20% of the sky is applied prior to calibration. 
7.3. Calibration efficiency on the sky
From the calibration algorithm, we derive the contribution from each sky pixel to the global gain bias for a simulation including either CMB signal only or CMB and Galactic emissions (Fig. 15). Around the Ecliptic poles, this contribution displays large variations because the dipole has a lower signaltonoise ratio. In these regions, the line of sight is orthogonal to the satellite speed direction, which remains close to the Ecliptic plane, leading to a small amplitude orbital dipole signal. We note also that whenever Galactic emissions are included the contributions from the Galactic plane pixels become very significant, which given that the sky signal in this area is notoriously difficult to estimate with sufficient fidelity, could give rise to a bias, as indeed already pointed out in Fig. 12. The impact of the Galactic plane is magnified by the weights, which need to be applied to the contribution of each pixel and are shown in Fig. 16. These define the sensitivity to the orbital dipole for each pixel. In particular, we see that the sensitivity is much higher in the Galactic plane than at the Ecliptic poles, hence the former not the latter is more likely to be important to the gain estimation. This observation justifies the use of a 20% Galactic mask and explains the observations made in Sect. 7.1.
Fig. 15 Maps of the contribution of each sky pixel to the global gain bias obtained for CMBonly (left) and CMB+galaxy (right) simulations. 
Fig. 16 Map of the weight applied to the contribution per pixel toward the global gain illustrating the sensitivity to orbital dipole. The highlevel sensitivity band matchs the overlap with the beginning of the third survey (see also Fig. 2). 
8. Solving for offsets and gain by iteration
Precise knowledge of the gains is required for destriping in order to precisely estimate the sky (I,Q,U) and the offsets, and remove the nonstationary component (the orbital dipole). In contrast, calibration cannot be performed before the removal of the lowfrequency noise, i.e. the subtraction of the offsets from the data. This leads to the following iterative pipeline scheme. At each iteration step, we calibrate the detectors independently and then estimate offsets on a multidetector basis. We repeat these operations until the relative difference between two consecutive gain estimations for each detector is smaller than 10^{7}.
8.1. Characterization of the convergence
To check the convergence of the iterations for both gain and offset reconstruction, we present in Fig. 17 an example of the decrease in the gain bias with respect to the number of iterations for one of our simulations of signal (CMB and Galactic emissions) plus pure white noise. The initial condition corresponds to random gain factors generated with 5% rms of the simulated value. We reach the convergence level in about five iterations. Figure 18 shows, for the same simulation, the rms of the offset residuals with respect to the iteration number.
Fig. 17 Bias of the reconstructed gains versus an iteration number in percent for a simulation with CMB, Galactic emission and pure white noise. The iterations are stopped when the relative differences between consecutive gains reaches 10^{7}. 
Fig. 18 Rms of the offset residuals with respect to iteration number in K_{CMB} for a simulation with CMB, Galactic emission, and pure white noise. 
Fig. 19 Gain biases for the four bolometers averaged over the Monte Carlo simulations, considering each dataset (from left to right: pure white noise without destriping, pure white noise after iterations, and low frequency noise with (α = 2,f_{knee} = 0.01), (α = 2,f_{knee} = 0.1), and (α = 1,f_{knee} = 0.1)). Error bars are derived from Monte Carlo simulations. dash line: systematics estimated in Sect. 7.1. 
Fig. 20 Average fullsky spectra of noise residual maps after the destriping and calibration pipeline for the datasets: pure white noise, and low frequency noise plus white noise with f_{knee} = 0.1 Hz and 0.01 Hz. The effect of the calibration uncertainty is represented by the spectra of the residual map constructed using estimated gains (blue dashed line). The simulation power spectra and CMB WMAP5yr fiducial model are also plotted for comparison. Upper line, from left to right: temperature, Emode, Bmode. Bottom line, from left to right: crosscorrelation TE, TB, and EB. 
Fig. 21 Residual maps of the Q and U Stokes parameters after iteration on a simulation of signal (CMB+Galactic emissions) plus low frequency noise (α = 1,f_{knee} = 0.1). Maps are degraded at 27.5 arcmin (HEALPix nside = 128). Largescale residuals are compatible with the orbital dipole leakage as shown in Fig. 13. Stripes are destriping residuals. 
8.2. Gain and sky map accuracy
We now present results from the analysis of simulations including four bolometers at 143 GHz (see Sect. 5 for details). The sky signal is the combination of CMB (including the dipoles) and Galactic emissions to which we have added different noise TOD of pure white noise, or correlated noise with either (α = 2,f_{knee} = 0.01 or 0.1 Hz) or (α = 1,f_{knee} = 0.1 Hz). One hundred noise realizations were generated in each case. Table 1 describes the results of the gain reconstruction for one bolometer (1431a). As for the destriping (see Sect. 6.3), the calibration precision is significantly worse in this case (α = 2,f_{knee} = 0.1 Hz). However, even in this case, the relative precision is on the order of a few 10^{4}.
Gain biases and errors for bolometer 1431a, for the four simulated noises: pure white noise, correlated noise with (α = 2,f_{knee} = 0.01 or 0.1 Hz), and (α = 1,f_{knee} = 0.1 Hz). Statistical error is constant (see Sect. 7.2). Polarized signals introduce a systematic bias which we reduce to 5.3 × 10^{5} by removing 20% of the brightest Galactic regions (see Sect. 7.1).
Figure 19 summarizes the calibration reconstruction accuracy for the four bolometers at 143 GHz. The gain biases are all consistent with zero given the statistical and systematic errors. The first set corresponds to gains reconstructed without destriping, for the pure white noise simulations. The following points represent the results of the pipeline for each dataset. The averaged gain for each bolometer is stable when changing the noise properties. The systematic bias, on the order of 5 × 10^{5}, results from the anisotropies of the CMB polarization, as explained in Sect. 7.1.
To evaluate the quality of the sky map reconstruction from our pipeline, we used a similar procedure to that adopted for the destriping analysis presented in Sect. 6. Figure 20 presents a compilation of averaged power spectra from residual maps for each datasets. The residual spectra have an additional excess power at low ℓ (below ~20) compared to those shown in Fig. 10. We attribute this excess to the calibration uncertainty that induces a leakage, of roughly 1 μK, of the orbital dipole signal into polarization. We computed the power spectra of the difference in the pure signal I, Q, U maps built with the recovered gains and the input maps (Fig. 13). As shown in Fig. 20, it matches the increase in power in the residual power spectra.
Figure 21 shows an example of the residual maps after iterations for one of the MC simulation. We clearly see the largescale structures correlated with the orbital dipole residuals together with the stripes produced by the destriping errors.
9. Conclusion
We have presented an iterative scheme for mapmaking, by means of the destriping and photometric calibration of PlanckHFI data. This method has been fully implemented and tested within the HFI Data Processing Center. We have used simulated datasets for 143 and 217 GHz HFI detectors, under various noise hypotheses to evaluate the performances and derive the impact of systematics.
We have first set the parameters of the destriping to minimize the systematic effects that could bias the offset determination (foregrounds signal, sky pixelization). We have shown that for correlated noise with a knee frequency up to 0.1 Hz, destriping residuals are below the white noise level except for ℓ < 50 for temperature and polarization. We found a significant bias for polarization at a knee frequency above 0.01 Hz. For a Plancklike Fourier spectrum presented by Planck Collaboration (2011), we found residuals below the white noise level above ℓ = 5.
We have then studied the performances of a gain reconstruction scheme based on the orbital dipole signal. We have demonstrated that in the same noise hypothesis as previously, we are able to reconstruct the photometric calibration with a statistical error of 5 × 10^{5}. For pure white noise, we have also evaluated the systematic uncertainty to be 5 × 10^{5}, which is dominated by CMB polarization anisotropies that are not modeled in our singlebolometer calibration scheme.
Finally, we presented the results when solving for offsets and gain by iteration. With low frequency noise, the gain bias stays constant but the uncertainties increase up to a few 10^{4} for correlated noise with (α = 2,f_{knee} = 0.1 Hz). The calibration uncertainty induces a further worsening of the destriping residuals concentrated at very low multipoles (ℓ < 20). Altogether, this scheme could be fruitfully applied to the PlanckHFI flight data, provided they satisfy our very basic noise hypotheses.
Acknowledgments
We acknowledge the use of the Planck Sky Model, developed by the Component Separation Working Group (WG2) of the Planck Collaboration. We thank also the CTP working group for fruitful discussions on mapmaking algorithms. We acknowledge use of CCIN2P3 facilities. The authors would like to thank F. Couchot for a detailed review of our manuscript.
References
 Armitage, C., & Wandelt, B. D. 2004, Phys. Rev. D, 70, 123007 [NASA ADS] [CrossRef] [Google Scholar]
 Arnau, J. V., & Saez, D. 2000, New Astron., 5, 121 [NASA ADS] [CrossRef] [Google Scholar]
 Ashdown, M. A. J., Baccigalupi, C., Balbi, A., et al. 2007a, A&A, 467, 761 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Ashdown, M. A. J., Baccigalupi, C., Balbi, A., et al. 2007b, A&A, 471, 361 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Ashdown, M. A. J., Baccigalupi, C., Bartlett, J. G., et al. 2009, A&A, 493, 753 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Betoule, M., Pierpaoli, E., Delabrouille, J., Le Jeune, M., & Cardoso, J.F. 2009, A&A, 503, 691 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Burigana, C., Malaspina, M., Mandolesi, N., et al. 1999, Int. Rep. TeSRE/CNR, 198 [arXiv:astroph/9906360] [Google Scholar]
 Burigana, C., & Sáez, D. 2003, A&A, 409, 423 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Cantalupo, C. M., Borrill, J. D., Kisner, T. S., Jaffe, A. H., & Stompor, R. 2010, ApJS, 187, 212 [NASA ADS] [CrossRef] [Google Scholar]
 de Gasperis, Balbi, A., Cabella, P., Natoli, P., & Vittorio, N. 2005, A&A, 436, 1159 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Delabrouille, J. 1998, A&AS, 127, 555 [Google Scholar]
 Doré, O., Teyssier, R., Bouchet, F. R., Vibert, D., & Prunet, S. 2001, A&A, 374, 358 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Efstathiou, G. 2007, MNRAS, 380, 1621 [NASA ADS] [CrossRef] [Google Scholar]
 Finkbeiner, D. P., Davis, M., & Schlegel, D. J. 1999, ApJ, 524, 867 [NASA ADS] [CrossRef] [Google Scholar]
 Golub, G. H., & van Loan, C. F. 1996, Matrix Computations 3rd edn. (John Hopkins University Press) [Google Scholar]
 Górski, K. M., Hivon, E., Banday, A. J., et al. 2005, ApJ, 622, 759 [NASA ADS] [CrossRef] [Google Scholar]
 Hinshaw, G., Weiland, J. L., Hill, R. S., et al. 2009, ApJS, 180, 225 [NASA ADS] [CrossRef] [Google Scholar]
 Hupca, I. O., Falcou, J., Grigori, L., & Stompor, R. 2010, INRIA technical report 7409 [arXiv:1010.1260] [Google Scholar]
 Janssen, M. A., S. Douglas, W. Martin, et al. 1996, unpublished [arXiv:astroph/9602009] [Google Scholar]
 Keihänen, E., KurkiSuonio, H., Poutanen, T., et al. 2004, A&A, 428, 287 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Keihänen, E., KurkiSuonio, H., & Poutanen, T. 2005, MNRAS, 360, 390 [NASA ADS] [CrossRef] [Google Scholar]
 Keihänen, E., Keskitalo, R., KurkiSuonio, H., Poutanen, T., & Sirviö, A.S. 2010, A&A, 510, A57 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Leach, S., Cardoso, J.F., Baccigalupi, C., et al. 2008, A&A, 491, 597 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Maino, D., Burigana, C., Maltoni, M., et al. 1999, A&AS, 140, 383 [Google Scholar]
 Maino, D., Burigana, C., Gorski, K. M., Mandolesi, N., & Bersanelli, M. 2002, A&A, 387, 356 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Patanchon, G., Ade, P. A. R., Bock, J. J., et al. 2008, ApJ, 681, 708 [NASA ADS] [CrossRef] [Google Scholar]
 Planck Collaboration 2005, ESA pub. ESASCI(2005)/01 [Google Scholar]
 Planck Collaboration 2011, A&A, 536, A4 [arXiv:astroph/1101.2039] [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Plaszczynski, S. 2007, Fluctuation and Noise Letters, 7 [Google Scholar]
 Poutanen, T., de Gasperis, G., Hivon, E., et al. 2006, A&A, 449, 1311 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Revenu, B., Kim, A., Ansari, R., et al. 2000, A&AS, 142, 499 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Reinecke, M., Dolag, K., Hell, R., Bartelmann, M., & Ensslin, T. 2006, A&A, 445, 373 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Seljak, U., & Zaldarriaga, M. 1996, ApJ, 469, 437 [NASA ADS] [CrossRef] [Google Scholar]
 Sutton, D., Zuntz, J. A., Ferreira, P. G., et al. 2010, MNRAS, 407, 1387 [NASA ADS] [CrossRef] [Google Scholar]
 Stompor, R., et al. 2002, Phys. Rev. D, 65, 022003 [Google Scholar]
 Stompor, R., & White, M. 2004, A&A, 419, 783 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Szydlarski, M., Esterie, P., Falcou, J., Grigori, L., & Stompor, R. 2011, INRIA technical report 7635 [arXiv:1106.0159] [Google Scholar]
 Tegmark, M. 1997, Phys. Rev. D, 55, 5895 [Google Scholar]
 Tristram, M. 2006, PoS (CMB2006), 063 [Google Scholar]
 Wright, E. L. 1996 [arXiv:astroph/9612006] [Google Scholar]
All Tables
Gain biases and errors for bolometer 1431a, for the four simulated noises: pure white noise, correlated noise with (α = 2,f_{knee} = 0.01 or 0.1 Hz), and (α = 1,f_{knee} = 0.1 Hz). Statistical error is constant (see Sect. 7.2). Polarized signals introduce a systematic bias which we reduce to 5.3 × 10^{5} by removing 20% of the brightest Galactic regions (see Sect. 7.1).
All Figures
Fig. 1 Dwelling time of the simulated stable pointing periods for the 12.5 monthlong observation. The variation in the period length reflects that the spin axis does not remain in the Ecliptic plane but follows a cycloid path. 

In the text 
Fig. 2 A map of hit counts per 3.4 arcmin pixels (Nside = 1024) in Galactic coordinates corresponding to 12.5 months of observation for four detectors. The sky is covered at least twice. Redundancy of the revisits is maximal around the Ecliptic poles. The conspicuous sshaped band with a higher hit count corresponds to the beginning of the third survey. 

In the text 
Fig. 3 Noise residual fullsky power spectra with various internal Galactic mask applied to the offset determination. Residual maps are at nside = 1024. Spectra are smoothed over five multipoles. Signal spectra are plotted in red. As a reference, we also plot the CMB WMAP fiveyear fiducial model (black) and a total power of the simulated sky signal, CMB+Galaxy+noise, (red curve). From left to right: temperature, Emode, and Bmode. 

In the text 
Fig. 4 Noise residual fullsky power spectra with various internal resolutions (nside = 64, 128, 256, 512, 1024) and neglecting polarization in the offset determination. Residual maps are at nside = 1024. Spectra are smoothed over five multipoles. Signal spectra are plotted in red. As a reference, we also plot the CMB WMAP5yr fiducial model (black). The horizontal, dashed lines show the predicted noise level. From left to right: temperature, Emode, and Bmode. 

In the text 
Fig. 5 Noise residual fullsky power spectra when using different internal sky resolution for polarization in the offset determination. Internal resolution for intensity is fixed to nside = 1024. Residual maps are at nside = 1024. Spectra are smoothed over five multipoles. Signal spectra are plotted in red. As a reference, we also plot the CMB WMAP fiveyear fiducial model (black). From left to right: temperature, Emode, and Bmode. 

In the text 
Fig. 6 Left: diagonal of the offset covariance matrix for various internal sky resolutions used for the offset recovery. Black line shows the first term of Eq. (17). Right: offset covariance matrix profiles centered on the diagonal showing the ringtoring covariance for various internal sky resolution. 

In the text 
Fig. 7 Offset covariance matrix without (left), with a Galactic mask (middle) and the difference (right). 

In the text 
Fig. 8 Binned offset covariance matrices. left: Polkapix estimate, right: covariance for 10 000 MC. The bin width is 100 rings. 

In the text 
Fig. 9 Offset correlation matrix for four polarizationsensitive bolometers at 143 GHz without (left) and with (right) polarization reconstruction in the destriping. 

In the text 
Fig. 10 Residuals spectra after destriping the signal plus noise in the three cases of (1) pure white noise (blue), low frequency noise plus white noise with (α = 2, f_{knee} = 0.1 Hz), (2) (α = 2, f_{knee} = 0.01 Hz) and (α = 1, f_{knee} = 0.1 Hz) compared to the simulation power spectra (red), and (3) CMB WMAP fiveyear fiducial model (black) on 100 simulations. Upper line, from left to right: temperature, Emode, and Bmode. Bottom line, from left to right: crosscorrelation TE, TB, and EB. 

In the text 
Fig. 11 Evolution of the calibration bias due to the sky resolution with the masked fraction. Error bars are statistical estimates from Polkapix. 

In the text 
Fig. 12 Evolution of the calibration bias due to the polarized signal and the masked fraction (blue: simulation of CMB; orange: simulation of Galactic emissions; red: simulation of CMB + Galaxy). Error bars are statistical estimates from Polkapix. 

In the text 
Fig. 13 Residual maps of the Q and U Stokes parameters for a calibration uncertainty of 5 × 10^{5}. The main features on large scales are related to the CMB dipole leaking into the polarization. 

In the text 
Fig. 14 Gain bias distribution for 1000 simulations including CMB dipoles and pure white noise relative to Polkapix posterior PDF. A mask covering 20% of the sky is applied prior to calibration. 

In the text 
Fig. 15 Maps of the contribution of each sky pixel to the global gain bias obtained for CMBonly (left) and CMB+galaxy (right) simulations. 

In the text 
Fig. 16 Map of the weight applied to the contribution per pixel toward the global gain illustrating the sensitivity to orbital dipole. The highlevel sensitivity band matchs the overlap with the beginning of the third survey (see also Fig. 2). 

In the text 
Fig. 17 Bias of the reconstructed gains versus an iteration number in percent for a simulation with CMB, Galactic emission and pure white noise. The iterations are stopped when the relative differences between consecutive gains reaches 10^{7}. 

In the text 
Fig. 18 Rms of the offset residuals with respect to iteration number in K_{CMB} for a simulation with CMB, Galactic emission, and pure white noise. 

In the text 
Fig. 19 Gain biases for the four bolometers averaged over the Monte Carlo simulations, considering each dataset (from left to right: pure white noise without destriping, pure white noise after iterations, and low frequency noise with (α = 2,f_{knee} = 0.01), (α = 2,f_{knee} = 0.1), and (α = 1,f_{knee} = 0.1)). Error bars are derived from Monte Carlo simulations. dash line: systematics estimated in Sect. 7.1. 

In the text 
Fig. 20 Average fullsky spectra of noise residual maps after the destriping and calibration pipeline for the datasets: pure white noise, and low frequency noise plus white noise with f_{knee} = 0.1 Hz and 0.01 Hz. The effect of the calibration uncertainty is represented by the spectra of the residual map constructed using estimated gains (blue dashed line). The simulation power spectra and CMB WMAP5yr fiducial model are also plotted for comparison. Upper line, from left to right: temperature, Emode, Bmode. Bottom line, from left to right: crosscorrelation TE, TB, and EB. 

In the text 
Fig. 21 Residual maps of the Q and U Stokes parameters after iteration on a simulation of signal (CMB+Galactic emissions) plus low frequency noise (α = 1,f_{knee} = 0.1). Maps are degraded at 27.5 arcmin (HEALPix nside = 128). Largescale residuals are compatible with the orbital dipole leakage as shown in Fig. 13. Stripes are destriping residuals. 

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