zonca / pico-simulations

Simulations for PICO using TOAST
0 stars 2 forks source link

Run PICO 4pi simulations with conviqt #24

Closed zonca closed 6 years ago

zonca commented 6 years ago

Using hpc4cmb/toast#241 I can run some simulations with the real beams provided by Karl.

Please crosscheck the configuration:

zonca commented 6 years ago

Should I also run the same simulations with no beam convolution, or with a gaussian beam?

ksyoung commented 6 years ago

That seems reasonable to me, but Brendan or Shaul should have the final word.

Running it with a simple beam (either no beam or a gaussian) sounds like a good idea to me. This would be an important cross check. I'm mildly for no convolution, unless you think that would add some confusing artifacts.

zonca commented 6 years ago

actually as we are just rescanning the input sky, I'll provide the input PySM map for comparison, so we get the "no convolution" map directly. As another cross-check, I'll run the simulations with conviqt but with a symmetric gaussian beam.

bcrill commented 6 years ago

I like Andrea's suggestion - it would be good to compare the output with those two cases. Thanks guys.

zonca commented 6 years ago

ok, I created the input alm and nside 512 map, see https://gist.github.com/zonca/c9ad8250857bc1d69a529b919fc3c911, they are available at:

/global/cscratch1/sd/zonca/pico-simulations/conviqt_sim/input/
zonca commented 6 years ago

@ksyoung simulation with real beams should be ready soon, I have realized I don't have a symmetric gaussian beam with the same fwhm of the real beams (6.2 arcmin) in alm format, I think they would be useful for comparison, can you please produce them?

zonca commented 6 years ago

all the scripts I am using are in #26, if @ksyoung or @bcrill would like to review them, it would be useful to spot any mistake.

zonca commented 6 years ago

@ksyoung sorry the real beams are 6.2 arcmin (I wrote 9.2 arcmin before)

zonca commented 6 years ago

I have the first results for the real beams for the 10h simulation, see:

/global/cscratch1/sd/zonca/pico/conviqt_sim/201806_boresight_1pix_conviqt_realbeams_10h

I checked the maps and they look reasonable, for this channel the beam is very small, consider that the average size of a 512 pixel is about 7 arcmin. It looks like the larger effect is ringing, see the comparison below.

pysm i diff

zonca commented 6 years ago

messed up the reference frame of the input, I'll fix an rerun

zonca commented 6 years ago

fix #27

zonca commented 6 years ago

@ksyoung @bcrill I got the results of the 10 hours simulation, the maps are at NERSC in:

/global/cscratch1/sd/zonca/pico/conviqt_sim/201806_boresight_1pix_conviqt_realbeams_10h

See the plots here: https://gist.github.com/zonca/7ec107d45c0a1265a43c9f4486750e20

I'll proceed with 1 year simulations.

zonca commented 6 years ago

I wanted to do an intermediate test before running 1 year, so I ran 30 days, see here the maps:

/global/cscratch1/sd/zonca/pico/conviqt_sim/201806_boresight_1pix_conviqt_realbeams_30d

pico_realbeams_30d

it shows some ringing. I'll go ahead with the 1 year to test it, but please let me know if that ringing is unexpected.

zonca commented 6 years ago

@ksyoung @bcrill got 1 year!

/global/cscratch1/sd/zonca/pico/conviqt_sim/201806_boresight_1pix_conviqt_realbeams_1y

See the executed notebook with plots here: https://gist.github.com/zonca/398b99eed03b1a9f52df264ab5c38222

1year_gaussian

1year

1year_diff

Please checkout the maps and let me know if you find anything wrong. I'll switch and work on the calibration simulations.

bcrill commented 6 years ago

The ringing from the Galaxy is too big. I doubt that it's a pipeline issue but more likely to be related to the beam map we're using for convolution. @ksyoung are you available again to take a look at this? Are there any sharp edges in the "real" beam map? We might be able to add a slight amount of apodization if that's the case and it should knock down the ringing.

ksyoung commented 6 years ago

I'll be back full time starting next week. I'll think about this in the mean time. I'm not sure if these beam maps have the center pixel replaced with an average or with zero. If the center pixel is zero maybe that's causing ringing. I'll look into it.

ksyoung commented 6 years ago

Hi @bcrill, I'm back in the office and can work on this again. I'm not sure what the issue is but my best idea is the following. The beam maps have ~ 1 deg pixels so maybe there is a sharp edge in the sense that the center pixel is 1 deg x 1 deg and the next pixel out is ~ 40 dB down. So this may effectively look like a 1 deg tophat beam, which would give ringing. Basically the GRASP output is too coarse of a representation of the few arcmin beam.

@zonca are you using the beam files "4pibeam[x,y]feed_center_150GHz_alm.fit" or the ones titled "4pi_centermask"? The centermask files have the central pixel replaced with 0's which may solve the tophat problem. It should also make the galaxy disappear since there is no main beam. This was done at Jacques suggestion to look at far sidelobe effects. Since the center pixel = 0 there still may be ringing effects, although they should be smaller.

@zonca, do you have an easy way to plot the beam alm files? Comparing the masked and unmasked versions may be helpful. The right thing to do may be to infill the center pixel with the average value of it's neighbours to make a beam with no discontinuities. This is a bit more work and invites me to make mistakes, so let's see if it is necessary first.

zonca commented 6 years ago

@ksyoung I ran 4pi_beam_[x,y]feed_center_150GHz_alm.fits, should I run the same simulations with the centermask as well?

the best way to visualize the alm is to convert them back to a map and then visualize the map.

ksyoung commented 6 years ago

@zonca I'm looking at the converted alm beams and they don't seem to quite make sense, so let's wait to run new simulations until I figure that out. Should be only a couple days.

ksyoung commented 6 years ago

@zonca @bcrill, I figured out the issue. The previous alms calculated from the beam maps were at too coarse of a resolution to reconstruct the far sidelobes. The maps from alms up to ell = 90 are completely dominated by ringing so the far sidelobes aren't even visible. I've figured out the ell and m range I need and am reproducing some GRASP beams at the proper resolution and recalculating the alms. I should have new beams to run by Monday, and hopefully they'll be ready Friday.

zonca commented 6 years ago

I got the new beams from @ksyoung and ran 30d and 1y simulations, they are available at NERSC at:

@bcrill see plots below:

30d 1y

please check them out and send me some feedback

ksyoung commented 6 years ago

Thanks @zonca! @bcrill, these look like what I would expect. The galaxy is broadly smoothed (by the 5 deg mask) and at low power, because the beam has been cut from about 65 dB to -15 dB. Below are some replots of the same map with different colorbar limits. You can see some structure at high galactic latitude that is probably due to the sidelobe pickup. We should compare to what real structure is there though.

Input galaxy, copied from earlier in this conversation: input_galaxy_from_zonca

1 yr sim at uK scale: beam_scan_1yr_pico

1 yr sim galaxy masked out at +- 30 deg galactic lat: beam_scan_1yr_pico_masked_30d_gal_lat

@zonca can you also send me the input map? Or make "/global/cscratch1/sd/zonca/pico/input_conviqt" where it lives have group access for mp107? Thanks1!

zonca commented 6 years ago

Fixed permissions

ksyoung commented 6 years ago

Thanks, @zonca! How is polarization handled? The output Q, U maps have max and min values of +-10^-19 mK. I expected Q,U values to scale similar to I values, so the Q,U maps would have amplitudes around 10^-6 mK. I get that estimate by noting the galactic I map input is ~ 10^5 mK and the I output is ~ 1 mK. The Q input is ~ .1 mK, so if it drops similar to the I map (just from the effect of masking the main beam) I get 10^-6. This may be completely the wrong scaling approach though. I haven't thought about this sort of thing before.

zonca commented 6 years ago

does the beam have also polarization components?

ksyoung commented 6 years ago

Yes, the beam is defined in GRASP as co-pol and cross-pol. Then the pipeline converts that to Stokes I,Q,U. See the input and output maps below. The reconstruction isn't as good for QU as for I, which I didn't check before, but it likely is still fine. Pre-alm fits. I, Q, U are upper left, upper right, lower left absolute value of maps, so sign of QU is lost combined_mainb_fullsky_iqu

reconstructed from alms: I, Q, U are upper left, upper right, lower left absolute value of maps, so sign of QU is lost from_alms_fullsky_400pxmask5d_l720_iqu

I don't know how the conviqt code deals with QU though. It seems like given the QU beam maps above convolving with the galaxy should give QU signals at the fraction of a uK level. But there may be something I'm missing.
Thanks!

ksyoung commented 6 years ago

For completeness here are the input galaxy IQU maps and the final IQU maps. @bcrill does this make sense to you? The final QU maps seem to not have enough signal. But maybe there are some effects I'm not thinking about. Input galaxy: conviqt_in_201810_smooth1deg_iqu conviqt out: conviqt_out_201810_smooth1d_iqu

keskitalo commented 6 years ago

Hi.

The polarization looks anomalously low. What are you using for the polarization efficiency? Are you using high enough expansion order to sample the main beam? It would also be helpful to smooth the input map with a representative Gaussian beam to have a rough idea about how the beam-convolved map should look like.

Conviqt assumes that the input beam expansion reflects the actual polarization efficiency. For a perfect co-polar beam in the Pxx frame, the intensity and co-polar beams should be identical and both normalized to 0.5 over the sky. The cross-polar beam should vanish.

ksyoung commented 6 years ago

There is some issue in the conversion of a stokes parameter map to alms by the Planck LevelS code. Q and U are being reconstructed incorrectly. Here's a dummy example using a 10 deg FWHM gaussian which is 100% crosspol (linear y polarized in GRASP language, or 100% +Q).

Beam in stokes parameters from 'grasp2stokes': 10d_gauss_stokes_20181011 Note the I and Q maps are equal, U = 0.

Beam reconstructed from alms after 'beam2alm' 10d_ypol_gauss_from_alms_20181011 The Q, U maps have not been reconstructed properly and have a quadrature sign flip.

Thoughts?

keskitalo commented 6 years ago

That looks more like a rotation between two beam reference frames. Are you by any chance rotating the beams to Dxx at some point?

ksyoung commented 6 years ago

I got in touch with Mark Ashdown and it sounds like there is a rotation happening in the LevelS code that I didn't understand. @keskitalo, Can you point me to a document that defines the Dxx, Pxx and similar coordinate systems from Planck?

zonca commented 6 years ago

@ksyoung, emailed you a doc about Planck reference frames

ksyoung commented 6 years ago

I talked to Mark Ashdown and solved the issue with the dummy beams. The 'beam2alm' code converts to spherical polarization basis which causes the quadrupole shape. When converting back to the GRASP defined basis I get, 10d_ypol_gauss_from_alms_20181015 as expected.

This is also true for the PICO beams. I get the correct reconstruction of Q and U when I convert to the correct polarization basis.

@zonca A few simple possibilities but quite unlikely . . .

More possible:

zonca commented 6 years ago

@ksyoung

ksyoung commented 6 years ago

@zonca I went back and looked at the full sky maps from the simulation with the unmasked beam and those also show no polarization response. This makes me wonder if there's something wrong in the orientations of the two beams.

Did the simple gaussian case you did, pico-sim issue 17 include a polarized input? If so the output maps should be just a smoothed version of the QU inputs, since the 10 deg Gaussian beams I sent you were perfectly polarized. If not can you do a similar test with that same I map and a Q map that is I / 10? Or some other simple case. Or scan the PySM input galaxy with the 10deg Gaussian.

zonca commented 6 years ago

ok @ksyoung as you suggested I ran with 4 detectors, and it looks like this fixed the issue, the polarization is now in the tenths of microK. Is it expected that Q is significantly larger?

See the output maps in: /global/cscratch1/sd/zonca/pico/conviqt_sim/201810_boresight_1pix4det_conviqt_realbeams_mask5deg_1y

pico_4det_i pico_4det_q pico_4det_u

ksyoung commented 6 years ago

Thanks Andrea! This does indeed look more like what we'd expect. And the Galaxy being mostly Q is also what we expect.

zonca commented 6 years ago

Simulation completed, details about the run at: https://github.com/zonca/pico-simulations/tree/master/conviqt_sim/runs/201810_boresight_1pix4det_conviqt_realbeams_mask5deg_1y