galsci / pysm

PySM 3: Sky emission simulations for Cosmic Microwave Background experiments
https://pysm3.readthedocs.io/
BSD 3-Clause "New" or "Revised" License
36 stars 23 forks source link

GNILC dust templates with polarization fraction tensor #72

Closed zonca closed 3 years ago

zonca commented 3 years ago

View and comment on the notebook at: https://app.reviewnb.com/healpy/pysm/pull/72/files/

The plan as suggested by Jacques Delabrouille is to add small scales fluctuations directly in needlet space, and Mathieu Remazeilles is working on it.

However I think it is useful to have a simpler model to compare to, so the plan here is to fit for a high-ell slope to the GNILC maps and then add small scale fluctuations. The difference with PySM 2 is that now I am using the log of the polarization fraction tensor as suggested by Andrei Frolov, it is based on polarization fraction, so we get the polarization fluctuations automatically scale with amplitude of the large scale emission:

image

review-notebook-app[bot] commented 3 years ago

Check out this pull request on  ReviewNB

See visual diffs & provide feedback on Jupyter Notebooks.


Powered by ReviewNB

zonca commented 3 years ago

Masking

Following what was done for PySM 2, we are masking the brighter parts of the galaxy and fit on the rest of the sky. So I am using a union of the Planck R3 intensity and polarization masks.

image

Is this the right approach, should we instead try a smaller mask? or use a small dust-dominated patch?

zonca commented 3 years ago

Resolution of the GNILC maps

image

Considering masking, here the sky fractions at each resolution in arcmin:

0,25.90%
10,0.00%
15,2.68%
20,6.65%
30,14.71%
60,44.62%
80,5.43%

For temperature Mathieu recommends to use the R2 release:

in order to be consistent with Q and U, the I map in this 2018 release was constructed without using the unpolarized channels 545 and 857 GHz, hence losing efficiency in cleaning CIB in the GNILC dust I map.

Which in fact has way better resolution: image

0,25.90%
5,39.50%
9,14.13%
12,7.19%
15,11.10%
22,2.17%

The issue is, if we are using the polarization fraction tensor, how bad is to use maps at different resolutions?

zonca commented 3 years ago

Fitting using GNILC R3

In the notebook attached to this PR I am loading the R3 maps (also for intensity), building the log of the polarization fraction tensor and fitting a slope to that.

I select manually a range of ell [100-300] which is clean from large scale fluctuations and before emission starts to flatten due to noise (I believe).

image image image

And that is exactly the problem, what section of the spectrum do we want to use for fitting?

brandonshensley commented 3 years ago

A few quick thoughts, but hopefully others can chime in:

  1. The masking looks good. Things get complicated (and possibly qualitatively different) too near the Galactic plane, so best to mask it even though there is a lot of signal there.
  2. Hopefully Andrei has thoughts re mixing resolutions. Definitely sounds non-ideal, but it would also be a shame to smooth the I maps.
  3. Is the reason that the TT spectrum flattens at the same ell as EE and BB just that you are using the polarization tensor formalism? It looks like the PySM2 model used ells well below 100 for their fits (see Figure 5). Are the differences due to the maps used or the polarization tensor approach?

Thank you for putting this together!

zonca commented 3 years ago

A few quick thoughts, but hopefully others can chime in:

  1. The masking looks good. Things get complicated (and possibly qualitatively different) too near the Galactic plane, so best to mask it even though there is a lot of signal there.

very good thanks

  1. Hopefully Andrei has thoughts re mixing resolutions. Definitely sounds non-ideal, but it would also be a shame to smooth the I maps.

ok wrote him email about this

  1. Is the reason that the TT spectrum flattens at the same ell as EE and BB just that you are using the polarization tensor formalism? It looks like the PySM2 model used ells well below 100 for their fits (see Figure 5). Are the differences due to the maps used or the polarization tensor approach?

I think it is due to the polarization tensor approach, down below see the plot of temperature only. Also notice that in the polarization tensor approach, the y axis is already logaritmic so I am plotting using semilogx instead of loglog as in the PySM 2 plots and the previous GNILC II only analysis. Is it possible there is something wrong in how I implemented that?

Figure 5 from the PySM 2 paper

image

Fit on GNILC 80 armin map II

This is from my previous analysis at: https://nbviewer.jupyter.org/github/zonca/pysm/blob/gnilc_dust/docs/preprocess-templates/planck_gnilc_dust_I.ipynb

image

zonca commented 3 years ago

Power spectrum of the GNILC R2 II map

image

brandonshensley commented 3 years ago

Do we understand the features in the spectrum above ell of 1000? Or do we not use that part anyway so it doesn't matter?

zonca commented 3 years ago

I have a new notebook with lots of improvements, will post on Monday

zonca commented 3 years ago

ok, I have fixed the issue using the polarization fraction tensor, now the fit looks better (now properly in log-log):

image image image

I changed the mask, I noticed that GNILC had already sources removed, so chose a large galactic mask (GAL060 from HFI PR2):

image

This was due to the fact that when I want to remove the small scales from the GNILC temperature map, then when I apply the mask and run NaMaster, I see residual power at high-ell, this is improved if I use a mask which has no point sources and I don't apply too steep a cut (I am using a sigmoid so I can tweak the steepness of the cut). Even now I can see the effect (please send suggestions on how to improve here, I would like to apply a steeper cut):

image

However the residuals are finally small enough that the output spectra look fine (before they showed a step in the II output spectrum around ell = 1400 because the residual small scale signal from the large scale signal map was added to the actual simulated gaussian small scales:

image image image

The full notebook with the analysis is at:

https://nbviewer.jupyter.org/gist/zonca/d00d3d2e17d902e5aac6f80e4c9f20f4

If anyone has suggestions on how to reduce the leaking of power to high ell, I would use a sharper transition between the input map and the gaussian realization, otherwise I think this is acceptable and I can produce the final templates for inspection from other people.

giuspugl commented 3 years ago

This is great! Are the spectra shown here in real units or in the pol. fract. tensor ones?

zonca commented 3 years ago

all of them in pol frac tensor, I'll also add a summary plot in standard IQU at the end of the notebook

zonca commented 3 years ago

TE

During the CMB-S4 meeting there was a mention of PySM 3 not having TE in the galactic models.

I made a first attempt at the fit:

image

the problem is that when I try to add the small scales generated using this fit, the polarization maps have way too much power and are completely off:

image

I'll keep investigating but it would be useful to have a general opinion about this and an answer to the 2 questions above.

b-thorne commented 3 years ago

Hi Andrea, I've been taking another look at this model, and was wondering about something. In this model we care about the polarization fraction, and so it is important to get the zero level of the intensity maps right. The GNILC maps do not seem to remove the CIB monopole, and get the overall zero-level contribution from dust slightly wrong (according to section 2.2 of 1807.06212).

I've taken a look at the effect of removing at least the CIB monopole from the GNILC intensity map, and combining this with gaussian realizations of Frolov (q, u) to produce (Q, U) with power tuned to the level in the BK maps (4.7 uK at ell=80). Whether or not you include the monopole in the intensity map seems to have a large qualitative effect on the simulations for a given power level, at these higher latitudes.

The small scale realizations you're doing here are slightly different, and I'm not sure yet what effect this has, if any. Have you already considered this?

zonca commented 3 years ago

no, haven't considered this at all, can you please share the code you used to fix the monopole level? or modify my notebook to incorporate that?

zonca commented 3 years ago

current issues slides: https://us02st1.zoom.us/web_client/ehjzr5/html/externalLinkPage.html?ref=https://docs.google.com/presentation/d/19Tcg_jccTgSg8KtZgf5PT0R2MVNFdl5_PM5mJjnN0EM/edit?usp=sharing

zonca commented 3 years ago

@bthorne93 just a reminder to send me details about fixing the monopole level.

b-thorne commented 3 years ago

Based on Section 2.2 of Planck 2018: XII, the contributions to the monopole in the 353 GHz GNILC intensity map are CIB and dust. They remove the CIB contribution, and apply some process to improve on the estimate of the dust monopole the maps already contained.

The numbers for each contribution can be found in Table 12 of Planck 2018: III, which I've pasted below. This seems to indicate there is a relevant (but very small) contribution from zodiacal emission as well, which is not discussed in XII.

The simplest first thing to do is just remove the CIB contribution, for which I've been using the number in Table 12: 0.13 MJy/sr at 353 GHz. Perhaps we also need to remove zodiacal emission, and perhaps we also need to refine the dust monopole model as was done in XII.

Screen Shot 2021-04-26 at 10 27 15 AM

brandonshensley commented 3 years ago

I'd be inclined not to adjust for a Galactic offset a la Planck XII 2018. They estimate it to be 63+/-40 uK_CMB at 353 GHz, and I've noticed that it is also sensitive to the HI dataset used which I don't think was explored in their treatment of systematics. But CIB (and maybe even Zodi) monopole corrections, yes.

zonca commented 3 years ago

Resolution with GAL080 mask

image

0,27.70%
10,0.04%
15,0.44%
20,4.26%
30,16.09%
60,45.93%
80,5.54%

image

0,27.70%
5,36.91%
9,14.51%
12,7.34%
15,11.29%
22,2.25%
zonca commented 3 years ago

ok, I have made a few changes:

Need to do:

the full notebook executed is at:

https://nbviewer.jupyter.org/gist/zonca/b398f85f701088ba60e3e66a20e0f48f

Spectra fitting plots

image image image

Fitted slopes are very close to zero:

'TT': -1.001
'EE': -0.053
'BB': -0.012
zonca commented 3 years ago

Thanks @brandonshensley and @seclark for the chat today, some notes:

zonca commented 3 years ago

see #74