GNILC templates processing

[giuspugl]

Summary

This PR addresses several open issues related to the new models of thermal dust.

https://github.com/galsci/pysm/issues/85, https://github.com/galsci/pysm/issues/84 ,https://github.com/galsci/pysm/issues/83, https://github.com/galsci/pysm/issues/69

Small scales in the IQU amplitude templates

In https://github.com/galsci/pysm/pull/74 we fit separately for TT, EE and BB fit a power law spectrum in different multipole ranges. However, as we get flatter spectral indices for polarization spectra, this will yield to injecting smaller angular scales in TT whose power at given multipole is smaller than EE and BB for all To avoid this, in this PR we force EE and BB small scales to follow the fitted TT power law. We also remark that in this way we get EB ratio closer to ~2.

Validate outputs estimating NaMaster power spectra on both input and outputs

we have expanded our validation by means of 3 figures of merit:

1. TT,EE, BB power spectra with Namaster
2. B-to-E ratio
3. Polarization fraction

Moreover, we have considered 4 masks to evaluate the quality of the small scales injected for different .

Another small difference with previous analysis is that we estimated the spectra on binned equally-spaced multipoles. We choose respectively for the 4 masks.

In notebook you can find the other validation plots.

Spectral parameters with small scales

Finally, in the last few cells of the notebook I've finalized a procedure to inject smaller angular scales to the spectral parameters.

• we use the Planck GNILC maps of Bd and Td.

• We inject smaller angular scales to the maps by extrapolating the power law fitted from the GNILC spectral parameter maps

• Smaller angular scales are modulated similarly as the intensity amplitude map

• the multipoles where the fit is performed are different given the observed spectra . In any case we don't fit beyond ell=400, which is consistent with the TT analysis above

• given the fact that we inject smaller angular scales with a steeper spectral index than TT

  

we don't expect to injecti small scale noise when rescaling at frequencies orders of magnitude lower or larger than the reference one ( 353 GHz).

[review-notebook-app[bot]]

Check out this pull request on 

See visual diffs & provide feedback on Jupyter Notebooks.

---

Powered by ReviewNB

[zonca]

thanks @giuspugl there are 2 notebooks in the pull request:

• expansion_logpoltens.ipynb is similar to what I had before with your improvements, so I can use it to generate the templates and create the docs
• what about gnilc_gaussian_scales.ipynb?

[giuspugl]

• what about gnilc_gaussian_scales.ipynb?

yeah sorry you can neglect this one .

[zonca]

@giuspugl is the BK15_region_Gal_apo.fits mask available publicly?

[giuspugl]

it doesn't seem so: https://lambda.gsfc.nasa.gov/product/suborbit/bicep2_prod_table.cfm

[zonca]

thanks @giuspugl it is only a mask I think I can share it. Can you please make sure it is readable and send me the location at NERSC?

[zonca]

thanks @giuspugl I have been able to reexecute your notebook, but results are a bit different.

Here is the original notebook you added in this pull request:

https://github.com/giuspugl/pysm/blob/206429fc58512b06c8781e9cfeac8041317f92c9/docs/preprocess-templates/expansion_logpoltens.ipynb

here is my clean execution, I have deleted all spectra in npz format so I'm sure I'm not loading an old result:

https://nbviewer.org/gist/zonca/feab8cf34805fddd95c2b0fb340f0c04

Some possible differences:

• I see you load masks HFI_Mask_GalPlane-apo2_2048_R2.00_GAL040.fits, have you done any processing, or are they the same as extension GAL040 of the mask file HFI_Mask_GalPlane-apo2_2048_R2.00.fits?

Plots

[giuspugl]

thanks @zonca,

have you done any processing, or are they the same as extension GAL040 of the mask file HFI_Mask_GalPlane-apo2_2048_R2.00.fits?

no those are the same masks as the one released in the file HFI_Mask_GalPlane-apo2_2048_R2.00.fits of Planck legacy archive

I think you are estimating the spectra within the same mask. I would change codecell #58 and #59 of your notebook with :

output_planck_mask= hp.ud_grade(planck_masks[k], nside_out=nside_out)  
planck_mask= planck_masks[k] 

[zonca]

@giuspugl changed that, plots are still different, see https://nbviewer.org/gist/zonca/a288a1f6964f8987e286723cba63e794

[zonca]

@giuspugl I executed the notebook again, see

https://nbviewer.org/gist/zonca/ff5176b4d6245ebe3ddd657706724bc3

unfortunately spectra are still different:

can you please put a cleanly executed version of your notebook with no pre-existing npz spectra on Gist so I can compare cell-to-cell and find differences?

[giuspugl]

Seems there is an inconsistency in the spectra estimated between @zonca and me . As a test case, i suggest to run this notebook https://gist.github.com/giuspugl/c5eb030329bc9aa12629ed0e8cd4a8d9 and compare the results.

[zonca]

thanks @giuspugl it seems results are consistent:

https://gist.github.com/b8ab166bf0f314793223849376c3e426

[zonca]

thanks @giuspugl so it looks like the same script running in a notebook gives a good result. The same running in a SLURM job on a Haswell node gives all NaN. I will file a bug report to NaMaster and see if we can track it down.

What should we do with the dust model? You were actually getting good results with a SLURM job and we both were getting bad results in a Notebook.

[zonca]

I am worried about dtypes of the arrays, I'll make some tests.

[zonca]

@giuspugl it was the dtypes. if I force the read maps to be float64, NaMaster gives consistent results. Filed an issue about this https://github.com/LSSTDESC/NaMaster/issues/149.

https://gist.github.com/zonca/958629d26af6794492bb5b85f0373095/revisions#diff-b6147ca046f048ae39d8725113fd7907fa84d7b472341b7bbc3f043f1354e95c

Now, not sure if this was affecting your script.

[zonca]

@zonca will test @giuspugl's scripts on Cori in @giuspugl 's scratch in extending_gnilc_dust/:

make_spectra.py
namaster_spectr.slrm

Test with healpy 1.14.x

[zonca]

@zonca will test @giuspugl's scripts on Cori in @giuspugl 's scratch in extending_gnilc_dust/:

make_spectra.py
namaster_spectr.slrm

Test with healpy 1.14.x

I don't think this test is very useful, the objective is to understand what is wrong in the notebook, so I'd rather directly focus on the notebook.

So I am executing the notebook carefully cell by cell cross-checking all steps. I am using the notebook is this branch, which is derived from: https://gist.github.com/giuspugl/d36a0e529a3b4a89c19bce5b58c9494a

@giuspugl first thing I noticed is in cell 53:

planck_masks ={ k:    np.ma.masked_equal(m ,  1 ).mask  for k,m in planck_masks.items() } 

you are turning the apodized masks into binary masks, what is the reason? doesn't it mess up the spectra with NaMaster?

[zonca]

ok, I commented out that line and spectra now agree with @giuspugl and seem fine, light color is input, darker color is output:

[zonca]

fully executed notebook: https://nbviewer.org/gist/zonca/515355d88ca0105e36e935ea28e9930c

Let's clarify #100 first, then I can go ahead and finalize the model.

[zonca]

Updated executed notebook:

https://nbviewer.org/gist/zonca/b0ec6cc14fa4c22ded497cc445a7ea46

This has been run at 2048, the production runs to create the templates is the same notebook just executed at 8192.

See also in this pull request the notebook gnilc_generate_map.ipynb which generates the templates from large scales and spectrum of small scales.

[giuspugl]

you are turning the apodized masks into binary masks, what is the reason? doesn't it mess up the spectra with NaMaster?

thanks a lot @zonca that was the bug!!

[zonca]

squashed to 1 commit to avoid having large notebooks in the history of the repo. Only the final execution of the notebook is merged into the docs.