What is the effect of reprojecting data circles as those of Fig. 3 (click here) on sky maps? The answer depends on the statistical properties of the noise along the "data circles", on the scanning strategy, and on the way that the relative "offsets" of the data circles are readjusted with respect to one another.
It is clear from looking at Fig. 3 (click here) and from the above discussion that we should first find a way to adjust the relative values of the offsets Ai of the different data circles before reprojecting on the sky. Even for very low values of the knee frequency, there is no hope to keep the signal from diverging on time-scales of a few months. Just setting the average level of each circle to zero is not good enough, since there is no reason the average level of the useful signal on each circle should be the same, and the only measurement we have access to is signal+noise.
For PLANCK, the observing strategy ensures that, in addition to the short term redundancy provided by "spin-chopping", there is also redundancy at long periods. In particular, as each data circle (obtained from 120 superimposed scans) crosses 3360 other such circles in two spots for a 1 year mission, it is quite natural, as a first order method for suppressing low-frequency noise effects, to try to readjust the average relative offset of the data circles by imposing that all the differences between signal over the same sky pixels but observed with different rotation axis of the satellite (at different times) are minimised. One key requirement on the scanning strategy for this method to work properly is that any given circle do not cut all the other circles in just a small number of pixels, in which case the adjusted value for the offset of that circle would depend drastically on the realization of the (white part of the) noise on the pixels used for the adjustment. For PLANCK, the scanning strategy is such that all the pixels of a circle are seen by at least one other circle, and thus all of them are used for the readjustment.
This method allows one to estimate the drifts due to low-frequency noise independently of the useful signal from the sky, since only differences of the total (signal+noise) at times where the beam is pointing on the same direction of the sky are used.
What should the accuracy of this method be? For a simple scanning strategy where the spin axis of the satellite is kept in the ecliptic, all circles play the same role (total symmetry), and thus , the rms of the error on readjusting the offset for circle i, should be the same for all circles. If we denote by the rms of the noise per pixel on a single circle obtained by averaging 120 scans, the accuracy of the determination of the offsets will be of the order of , where Nj is the number of independent pixels per circle.
It would be possible to optimise the scanning so that there be no preferred direction for this residual striping (which is anyway small, as for the PLANCK HFI Nj is about 2000 to 4500, depending on the channel), but it does not seem worth compromising the monitoring of other systematics and noise in the process.
Note that if we could build an ideal experiment with no low-frequency noise, the value for that we would estimate just from sample variance on a circle of Nj pixels would be , and thus no significant additional striping should be added by this inversion.
In order to check the above assertions about the accuracy of offset readjusting, and to evaluate the effect of low frequency noise on PLANCK maps, a low-resolution version of the expected PLANCK data has been simulated. For the moment, the simulation of the complete data stream for one detector at the actual PLANCK resolution is out of reach of our computers. However, the conclusions obtained with a degraded resolution can be scaled to the actual resolution.
In this simulation, a vector of spin-axis positions on the sky corresponding to a scanning strategy of PLANCK is generated. Each position of the spin-axis is distant from the previous one by a step equal to the resolution at which the simulation is performed, and for each such position, a vector of beam positions on a circle of radius 70 on the sky is generated. For each beam position, a single one degree by one degree pixel on a simulated map of the sky is selected. A two-dimensional set of data corresponding to these beam positions is generated. Here we assume that the knee frequency of the low frequency noise is small enough that there be no significant low-frequency contribution to the noise along individual data circles. We simulate the effect of low frequency drifts by adding to each data circle some offset Ai.
The actual data taking process is such that the zero level of the measurements is frequently readjusted, in order to avoid the saturation of the detectors due to slow drifts in the signal. In the end, the average value of the signal on each circle is not known precisely, but certainly does not drift by more than a couple orders of magnitude (extremely conservative assumption!).
The actual inversion is performed by the following least mean square method. The linear system of equations on the constants Ai to be determined is obtained by requiring a function S(A1, ...,An) to be minimal. Here we take:
where is the estimated difference of offsets obtained by comparing the measurements of circle i2 and circle i1 on pixel p, n(p) is the number of times pixel p is seen by the experiment, and wi1,i2 is a weight attributed to measurement . Since pixel p contributes terms in the sum, we set wi1,i2 = 2/(n(p)-1) (the total weight of the contribution of the n(p) measurements at pixel p is thus proportional to n(p)).
We get the linear system to be inverted by writing that all partial derivatives of S with respect to all constants Ai should vanish. This system is degenerated, since if a set of constants A1, ...,An is solution, the same set with a constant K added to all the Ai is also solution. This degeneracy is lifted by adding the requirement that the mean of the constants Ai vanish.
We first use the nominal scanning strategy of PLANCK, with a spin axis kept rigorously anti-solar, and check that the accuracy of the inversion does not depend on the input offsets for the circles, but does improve with the resolution. For each simulation, the data of the mission are stored as a table of Ni scan circles of Nj samples each. For a simulation with a 1 degree resolution, Ni = 421 and Nj = 339 for a 14 month mission and a 140 scan diameter. We evaluate the amount of striping by calculating the mean value of the signal on each scan circle, and computing the square root of the variance of the collection of Ni numbers so obtained. This is denoted as .
For a 1 resolution simulation, we generate a element table of Gaussian random numbers with rms 1. This simulates the white noise of the mission. To this noise we add offsets which simulate the low frequency drifts. Various types of offsets were considered: a) no offset, in order to test that the method does not add striping where there is none; b) a linear drift of 100i/Ni, where i is the index of the circle and Ni the number of circles; c) a sinusoidal drift of ; d) a sequence of random offsets with an rms of 100 (hundred times the rms of the noise); e) a random walk of rms 10 per step; f) in order to check the intrinsic precision of the inversion, we try to recover the offsets of the random walk without adding any white noise.
In order to be able to compare the performance of the method in all cases above, let us consider one specific realization of the white part of the noise here, and add to it the offsets described above to get five different "noise signals". In all cases a) through e), the value of after the inversion is 0.0412, which proves that the accuracy does not depend on the offsets to be corrected for to any significant level. It is interesting to realize also that if there is no offset, the expected value of before any inversion just from sample variance is 0.0543 (and on the particular noise realization we used here it happened to be 0.0522), so that the method does not only remove the striping due to drifts in the offsets due to low-frequency noise if any, but does even readjust to some extent the variations in the average level of the circles due to the sample variance of the average of the Nj points on a circle. This is anecdotic, of course, but does show that no additional noise or bias is generated in the readjusting process - in fact, for a 14 month mission the method partly suppresses low frequency components of the white noise itself. It is worth stressing that the success of the method is due for a large part to the fact that all the points of any circle contribute to the evaluation of the value of its offset, and not only two points at the North and South ecliptic poles. The precision of the inversion (performed in double-precision) is evaluated from the results of case f) above, for which we get a totally negligible residual striping of .
In order to check that the accuracy does indeed depend on the resolution, we repeat the simulation of a) through e) for resolutions of 2 and . In these cases, the value of after the inversion is 0.0549 and 0.0800 respectively, independent of the original stripes again, and below the value we get from sample variance on the original data. Note that the weight given to each term in the least-square sum above is critical to reach this accuracy. For instance, if we put the same weight to each term, too much importance is given to points on the circles that are close to the north and south ecliptic poles (because they are "seen" by many more circles than the ones in the ecliptic), and some residual striping of the order of magnitude of the white noise per pixel remains, because the value of the offset recovered depends most on the realization of the white noise at points near the poles.
We tried this inversion scheme with all kinds of offsets, with various kinds of scanning strategy and different mission durations (from 3 to 14 months). In all cases the inversion works extremely well. To be a little more illustrative, we show images of a noise generated according to method e) above, for a 1 year mission, at the resolution of , reprojected on the sky, with no destriping treatment (Fig. 6 (click here)) and after destriping (Fig. 7 (click here)). For this simulation we used a scanning strategy for which the spin axis has been made to oscillate sinusoidally out of the ecliptic with an amplitude of 15 and a frequency of 8 oscillations per year, in order to maximise the sky coverage. Note the different scales for the amplitude of the structures that can be seen. No striping whatsoever remains on the map after the treatment.
Figure 6: Map of reprojected noise and offsets if no processing of the signal is done. The map is heavily striped. This should not be interpreted as a noise map for PLANCK, as at least some offset readjustment can be made before reprojecting. For this data set, the average level of circles has been left to drift some hundred sigma away!
Figure 7: Map of reprojected noise and offsets after signal-preserving destriping (see text). Compare with Fig. 6 (click here) (note the different colour scale!)
Whereas it has been made clear that our method is satisfactory to remove the striping due to low-frequency noise off maps obtained by PLANCK in the case where the knee frequency of the noise is low enough that all significant low-frequency noise appears in the form of different offsets for different circles on the sky, it is worth mentioning that if there is significant low-frequency noise power at frequencies higher than the spinning frequency then the method above may not be sufficient anymore: After relative offsets are readjusted, the next effect of low-frequency noise is to add very low level drifts along individual circles, as shown in Fig. 4 (click here) and Fig. 5 (click here), and one should investigate whether or not these fluctuations may cause a problem.
It is not easy to generate degraded resolution maps that preserve exactly both the visual aspect (i.e. do we see striping or not) and the statistical properties of the noise for real 1/f noise on all scales. For instance, the rms amplitude of low-frequency noise relative to white noise depends on the sampling frequency (the variance of white noise is proportional to the sampling frequency, and the variance of 1/f noise diverges proportionnaly to the logarithm of the sampling frequency). For these simulations, we decide to scale both the white noise and low frequency noise rms to their rms values for 10 arcminute pixels, independently of the size of the pixel of the simulation (here 1 square degree pixels). The next problem is that for 1 degree resolution simulations, instead of one circle every two hours, we generate one circle every day. However, we wish to preserve the possible correlations between consecutive circles of the simulation. Thus, we have to face the problem that either the correlations between consecutive circles is right, or the relative value of the offset between circles taken at a 1 day time interval is right. As the paragraph above demonstrates that relative levels of the circles do not play a role in the accuracy of the readjustment of the offsets, and as offsets will not drift too far away because they will be readjusted by the instrument anyway, I decided, in this next part of the simulations, to preserve the correlation between consecutive circles rather than their relative offsets.
Thus, a time sequence for the noise is generated by the following method: for a 1 year mission at a 1 degree resolution, we need 361 circles (so as to re-observe the first circle at the end of the mission) of 339 points each. Each of these circles is obtained from averaging 120 scans, and thus we need to generate a noise of almost 15 million data points. We wish that these data points preserve the relative rms amplitudes of low-frequency noise to white noise of the full 10 arcminute resolution maps, i.e. we do not smooth the noise down to the lower resolution, but rather affect to each 1 degree by 1 degree pixel the noise of its "central" 10 arcminute by 10 arcminute pixel, which corresponds to reprojecting on a 1 degree resolution map a subset of the 10 arcminute resolution pixels. This method preserves best the relative rms of the low-frequency residual noise and the white noise part of true full resolution maps.
Generating a 15 million sample dataset with a "true" 1/f spectrum is out of the reach of the computer I used. In order to generate the 1/f noise, an under-sampled (by a factor 32) white noise with a unit variance was generated, converted to Fourier space (by FFT) where it was multiplied by (with the appropriate rescaling to take into account the effect of under-sampling). Then we compute the inverse FFT and interpolate between the samples, add consecutive circles by packs of 120 to get a array of data points (corresponding to 361 circles), and add to it a randomly generated array of data points with a Gaussian statistic and a rms of . We then re-scale everything conveniently by multiplying by . We check that the relative rms values and the visual aspect of the circles obtained in this way are correct by comparing with circles obtained from fully fast-sampled simulated low frequency noise of a few circles only.
First, let us generate such a noise for nominal low-frequency noise parameters for the PLANCK SURVEYOR High Frequency Instrument, i.e. and , and readjust the offsets by the method of the previous sub-section. Again we use the scanning strategy for which the spin axis has been made to oscillate sinusoidally out of the ecliptic with an amplitude of 15 and a frequency of 8 oscillations per year. The resulting output noise map is shown in Fig. 8 (click here). No striping whatsoever is apparent on the map, and the increase of the rms of the noise is about 0.47%. In order to check the effect of the out-of-ecliptic motion, we do the simulation with the same realization of the noise but with a nominal anti-solar spin axis. In this case the increase in the rms of the noise is 0.49% (no significant difference). Repeated simulations with different noise realizations show that the small difference is always in favour of sinusoidal out-of-ecliptic motion.
Figure 8: Map of reprojected simulated 1/f noise, , , after correction by simple offset readjustment. The very-low level residual striping can not be distinguished by eye. This is a 1-degree resolution projection of what the processed PLANCK all-sky noise could look like. The relative standard deviations of reprojected white noise and reprojected residual striping are those of a 10 arcminute resolution mission (see text)
In their 1996 paper, Janssen et al. computed the maximum noise increase along a scan circle. Adapted to our notations, their formula can be written as
where m is the number of independent data points along one circle. For PLANCK, m = 2034 for a 10 arcminute resolution, sec., Hz, and we get . The average additional noise that we find by simulations is about half the maximum noise increase along one circle they predict (for diametrically opposed pixels), so our result is fully consistent with theirs. This is consistent with most of the additional noise being due to small drifts along individual circles.
In order to give an estimate of the order of magnitude of striping and the statistical properties of the noise, we show, in Fig. 9 (click here), plots of cuts through the noise map of Fig. 8 (click here). On both panels of Fig. 9 (click here), total noise is represented with diamonds, and the component of noise due to residual low-frequency drifts (obtained by computing the difference between the output map of Fig. 8 (click here) and the map obtained by simple reprojection on the sky of the white part of the noise used for the simulation only) as a plain line. In these plots, although one point only is plotted per degree, the spread in represented points is typical of what we would get in a 10 arcminute resolution map, not a smoothed map with a 1 degree resolution. The top panel corresponds to a vertical cut in the middle of the map, and the bottom panel to an horizontal cut in the middle of the map (thus perpendicular to the expected striping, if any). In both cases the total contribution of striping to the total noise is very small. The structure of residual low-frequency noise is not similar in both directions, which can be understood from the direction of scans. The small structure on the plot of residual low-frequency noise in the middle of the top panel is characteristic of very-low level striping. Finally, it is obvious from the distribution of the total noise that some regions of the sky are integrated more than others (as the region around ecliptic longitude and ecliptic latitude ).
Figure 9: Cuts through the noise map of Fig. 8 (click here). The top panel is a vertical cut at 0 ecliptic longitude, and the bottom panel an horizontal cut along the ecliptic plane. The plain line represents residual striping, obtained by computing the difference between the map of reprojected total noise after inversion, and the map of reprojected white part of the noise only. Diamonds correspond to the total reprojected noise. The structure due to striping is very small compared to the total noise, and the noise r.m.s. increase due to this residual striping is 0.47% on a 10 arcminute resolution map of the sky
Now we want to investigate what would happen if the low frequency noise were much worse than expected. The purpose of the following simulations, using parameters that are unrealistically pessimistic for the PLANCK HFI, is to show that the conclusions of the previous paragraphs do not depend drastically on noise assumptions, and that even unforeseen instabilities can be monitored quite well with PLANCK. Here we assume that some unwanted temperature fluctuations (for instance) generate low frequency noise with the parameters of (ten times the nominal!) and (same parameters as were used to generate the plots of Figs. 2 (click here)-5 (click here)). The resulting increase in the noise standard deviation of the map (scaled to 10 arcminute pixel size as explained above) is 7.5%. This is not much, but it is significantly larger than the value of 0.5% obtained with nominal HFI noise. Because the excess power is not white, it could be distinguished on the maps, especially at degraded resolution. For instance, within the framework of standard assumptions of physical cosmology this could cause trouble for the estimation of the properties of primordial supra-horizon fluctuations, which are important for constraining the models of inflation. It could also be annoying in pattern recognition methods looking for discontinuities generated by cosmic strings. Finally, some optimal foreground separation methods rely on statistical methods which use prior knowledge of the spectrum (as a function of scale) of the various astrophysical components (i.e. cirrus clouds, primary CMB, free-free, synchrotron, etc...).
In order to get rid of this residual striping, we can try a more sophisticated treatment. The idea again is to find a method which does not depend on the real signal, and thus estimates low-frequency components by using signal differences on common pixels. To do so, we can adapt the above method: instead of fitting just one constant for each data circle, we fit a function with more parameters, so that along a scan circle i we may write
In the equation above, j indexes the samples along the data circle, and the functions fk(j) are vectors of a basis of functions on which to decompose the noise ni(j). Typically, the set of functions fk can be sines and cosines (Fourier modes), or polynomials, or other well-chosen functions. The sum we want then to minimise is:
where j1(p) and j2(p) index the samples on circles i1 and i2 respectively for which pixel p on the sky is observed.
In this framework, for instance, we can take advantage of the fact that the most interesting property of low-frequency noise is that it does not have significant high frequency power. Because of that, low frequency noise itself can be estimated by sampling it at a much lower sampling frequency than the true signal. Thus, it seems to be a good idea to use as a basis of functions for noise estimation along one circle a set of a few "top-hat" functions, corresponding each to a constant on a fraction of a circle only.
Using the same realization of low-frequency and white noise, we inverted the data set by adjusting more than one constant for each circle. We do it for two constants (i.e. one for each half-circle), three constants (one for each third of a circle), and four constants (one for each fourth of a circle). We performed the inversion also using through degree polynomials.
Table 1 (click here) gives the noise rms increase on maps in all the cases discussed above. Maps, too numerous to be shown here, can be provided by the author upon request. It is clear that the residual striping can be reduced substantially by this method (and the more so at high resolution, as more points are available to estimate individual parameters of the fitting functions fk, whereas no more functions are needed at high resolution than at low resolution to estimate low-frequency structures in the noise). 4 constants per circle instead of one is the best that could be done because of computational limitations. For full-resolution data sets, with a good computer, going to 10 constants or so should be possible and should improve the fit significantly. This method, used on very pessimistic noise here for illustrative purposes, could also be applied to destripe further full-resolution maps obtained with nominal PLANCK HFI noise, if one wished to do so.
|1 constant||+ 7.5%||+ 7.5%|
|2 constants||+ 3.3%||+ 2.9%|
|3 constants||+ 3.2%||+ 2.3%|
|4 constants||+ 2.9%||+ 6.3%|
For each noise spectrum there must be an optimal set of functions fk to use. For instance, a degree polynomial is better than two constants per circle for very steep noise spectra, as most of the low-frequency noise contribution comes from very low frequencies. For each noise spectrum, there is also an optimal number of functions fk to use, as the more functions one uses, the less constraints one gets on each of the functions. These optimal solutions are yet to be found.