Understanding PolSpice results

Questions:

I am getting weird or negative C(ℓ) values for TT (or EE or BB) auto-spectrum. Is there something wrong?
Why do I get spurious EE and BB spectra on cut sky polarized maps?
Can I make Spice run faster, especially in the Monte-Carlo analysis of many maps?
Can I read arbitrary extensions and/or columns from a FITS file?
What is the 'double precision' compilation?

Back to PolSpice Main Page Known Bugs FAQ

Answers:

Q: I am getting weird or negative C(ℓ) values for TT (or EE or BB) auto-spectrum. Is there something wrong?

A: This is not un-expected, especially if you are analysing a small fraction of the sky (ie f_sky ≪ 1).
PolSpice returns an un-biased estimator of the C(ℓ), whose average (on many different random realisations of the same model) converges toward the true spectrum. However, as for any estimator, the C(ℓ) returned for a particular sky will oscillate around this expected value, with a variance inversely proportional to f_sky. Since PolSpice does not impose the positivity of the auto-spectra, this may lead to weird or negative auto-spectrum values (or large positive values at places where a vanishing result was expected).
To mitigate those fluctuations one can:

make sure that apodizesigma and thetamax are set to values close to the diameter of the sky of patch being analyzed.
For instance: apodizesigma ~ thetamax ~ 180 deg when f_sky > 0.5,
apodizesigma ~ thetamax ~ 20 deg when f_sky =0.01.
Tweaking their values will change the amplitude and period of the spurious oscillations.
The relation between the C_Spice(ℓ) produced by PolSpice and the "true" underlying C_true(ℓ) is in average given by
C_Spice(ℓ₁) = Σ_ℓ₂ K(ℓ₁,ℓ₂)C_true(ℓ₂)
where the kernels K, made available for temperature and polarization with kernelsfileout, depend on
apodizetype, apodizesigma, nlmax, thetamax (but not on the signal map nor on the pixel mask or weightfile).
See Eqs. (56-60), (62) and (91) of Chon et al (2004) and the FITS header of kernelsfileout for details;
regressing out the best fit monopole and dipole (subav and subdipole options) is always a good idea,
especially when some pixels are masked, or down weighted;
use an apodized weightfile instead of a boolean maskfile, to reduce the cross-talk between multipoles;
The F90 facility process_mask or its IDL wrapper iprocess_mask.pro can be used to apodize a mask.
and finally, since it is not possible to measure properly individual multipoles in the presence of small f_sky
one should bin the output C(ℓ), with roughly Δℓ ~ 1/f_sky.
The Healpix/IDL routine bin_llcl.pro or the python routine bin_llcl.py can be used to do so.

Q: Why do I get spurious EE and BB spectra on cut sky polarized maps?

A: When analyzing polarized cut-sky maps (ie with a maskfile and/or weightfile and with polarization=YES),

it is important to set decouple=YES to reduce the E/B coupling induced by the pixel masking;
decreasing the tolerance keyword (which quantifies the relative error accepted on numerical integrals)
from its default value of 1.e-5 can also improve the EE and BB spectra accuracy,
at the cost of longer compute times.

In the case of Planck-like apodized masks (f_sky ~ 0.5),
it was found on ideal simulations with a Planck-like input model and vanishing BB,
that setting tolerance=1.e-6 improved the accuracy of the recovered EE at larger multipoles,
and that of BB at all scales, with a negligible cost in CPU or wall-clock time.

Q: Can I make Spice run faster, especially in the Monte-Carlo analysis of many maps?

A: Version v03-04-01 (of Aug 2017) of PolSpice has been optimized for speed in high multipole (nlmax) analysis.
However, when analyzing a polarized map, one of the most costly steps of PolSpice remains the computation of a normalization function for the produced TE (and TB) spectra, which depends on the parameters apodizetype, apodizesigma, nlmax, and thetamax.
A key ingredient of this calculation is the TE kernel (or the sum over its rows) which can be exported and re-imported with a combination of windowfileout and windowfilein in all PolSpice versions, or the less disk-space hungry tenormfileout and tenormfilein in version 03-04-01 and above.
Therefore, when many (pairs of) polarized maps are analyzed with the parameters apodizetype, apodizesigma, nlmax, and thetamax set to a fixed value, as in the case for Monte-Carlo simulations of a given observation, it may save a lot of CPU time to set,

In all PolSpice versions:
on the first run: windowfilein=NO, windowfileout= a FITS_file_name1,
on all subsequent runs: windowfileout=NO, windowfilein= same FITS_file_name1.
The FITS file thus generated will have a size in MB of 67x(nlmax/2048)²;
In version v03-04-01 and above:
on the first run: tenormfilein=NO, tenormfileout= a FITS_file_name2,
on all subsequent runs: tenormfileout=NO, tenormfilein= same FITS_file_name2.
The FITS file thus generated will have a size in kB of 60x(nlmax/2048).

Q: Can I read arbitrary extensions and/or columns from a FITS file?

A: The extension and columns to be read can be specified in mapfile*, maskfile*, weightfile*, listmapfiles*, listmapweights* and beam_file*.

To read columns 2,3 and 6 (1 based), named for instance TEMP, SQ, X from file map.fits, you can use
mapfile=map.fits[col 2,3,6] or
mapfile=map.fits[col teMp, Sq, X] or
mapfile=map.fits[col te*, s*, X*].
To read extension 2 (1 based), named for instance SKY from file map.fits, you can use
mapfile=map.fits[2] or
mapfile=map.fits[SKY]
These filters can be combined, eg
mapfile=map.fits[2][col 2,3,6]

These CFITSIO features require PolSpice to be compiled with a recent version of Healpix (3.10 or more).

Back to Top

Q: What is the double precision compilation?

A: Because of the multithreading (= parallel execution) used in PolSpice and the aggressive optimization performed by current compilers, two consecutive runs with the same input data and parameters can see their numerical operations done in different order, which, because of round-off errors, may lead to slightly different final results.
Starting with version v03-05-01, PolSpice is by default compiled in 'Double Precision', allowing for smaller round-off errors and better reproducibility of the results, at the cost of larger memory foot-prints, and a bit longer computation times.
Compilation in the original 'Single Precision' setting remains possible (see the INSTALL document in PolSpice top directory).
This change in precision concerns only large arrays used as intermediate variables, so that whatever the precision chosen, the input maps, masks and weights can be read indifferently from single or double precision FITS files, and the ouput C(ℓ) will always be in double precision FITS files.

Back to Top

Updated 2018-11-05

For unlisted problems, contact Eric Hivon (hivon at iap dot fr)
Back to PolSpice Main Page Known Bugs FAQ