ss/deconvolving ringmap

[ssiegelx]

refactor(ringmapmaker): refactor DeconvolveHybridM

• Move most of the contents of DeconvolveHybridM into
  DeconvolveHybridMBase.

• Create subclasses that enable the user to deconvolve either the
  driftscan analytical beam model or an external beam model.

• Create subclasses that enable a Tikhonov scheme for deconvolution
  that uses a single regularisation parameter, or Weiner scheme that
  assumes some model for the m-mode power spectrum of the signal.

• DeconvolveHybridMBase now computes the el-dependent EW dirty beam
  and saves it to the output RingMap container if requested.

• DeconvolveHybridMBase now normalizes by the el-dependent
  peak of the dirty beam so that point source peak flux is
  preserved through the map making process.

• Create beam.py module and CreateBeamStream task that takes a
  CommonModeGridBeam container and converts it to a HybridVisStream
  container that can be m-mode transformed and used in the
  deconvolving ringmapmakers.

• Add an el axis to the weight dataset in the RingMap container,
  since the noise after deconvolution will be el dependent.

• Propagate noise variance through the deconvolution process.
  Create RADependentWeight task that re-introduces the RA dependence
  to the weights that was lost in the round-trip m-mode transform.
  Applying this procedure to Gaussian noise realizations recovers
  the standard deviation of the noise at the <10% level.

Branches off deconvolving-ringmap (PR #122)

[lgtm-com]

This pull request introduces 1 alert and fixes 1 when merging 537a0e4 into 0c4b397 - view on LGTM.com

new alerts:

• 1 for Unused import

fixed alerts:

• 1 for Unnecessary pass

[wulfda02]

@ssiegelx Here is what I have been using for the sky covariance:

def S(f,m):
    dfse = 2 * np.power(f/600., -3.5) * np.power(m, -1.5)
    psrc = 2e-3 * np.power(f/600., -2)
    return np.sqrt(dfse**2 + psrc**2)

Attached, I've plotted this function over the m-transform of a deconvolved ringmap (absolute square, so units of Jy^2). The different colors are different frequencies, which have been artificially shifted so they don't lay on top of each other. Note that this is not a fit---just tweaked by eye.

sky_covar.pdf

The dashed lines correspond to the equivalent prior inv_SN = 1e-6 that's currently in the code, assuming that the variance is ~1e-8 Jy^2. Based on these plots, I think we are over-estimating the signal to noise at high m with the default prior. Hopefully that's consistent with what @jrs65 has found---that a larger inv_SN works better.

[jrs65]

@ssiegelx I didn't find that I needed to tweak the SN parameter very much. It looks like I was using 1e-5 in my runs.

[ssiegelx]

Thanks @wulfda02 and @jrs65.

I generated maps using a range of values for inv_s2n. This pdf file summarizes the results:

inv_s2n_noise_optimization_driftscan_compressed.pdf

Page 2 shows the expected thermal noise as a function of inv_s2n (relative to the noise at the default inv_s2n = 1e-6). Increasinginv_s2n from 1e-6 to 1e-5 results in a 30% reduction in the standard deviation for XX pol and 12.5% for YY pol. I think this explains a decent chunk of the discrepancy between the noise in the ringmap stack versus the beam former stack.

Increasing inv_s2n from 1e-5 to 1e-4 does not reduce the noise further by much (roughly 7.5% for XX pol and 2.5% YY pol). So 1e-5 seems like a good number.

The maps at a single frequency are shown on pg 3. They don't look very different to me. However pg 4 and pg 5 show a slice through two different point sources. You can see that the ability to deconvolve the beam degrades as you increase inv_s2n, as evident by the ratio of the peak response to the grating lobe response.

[lgtm-com]

This pull request introduces 1 alert and fixes 1 when merging 78c08f8 into 6135579 - view on LGTM.com

new alerts:

• 1 for Unused import

fixed alerts:

• 1 for Unnecessary pass

[wulfda02]

I'm not really qualified to asses changes to core.io and core.task, but everything else looks good to me. My two comments are minor accounting details, not something that will break anything.

[jrs65]

@ssiegelx I've forced pushed a branch with the conflicts resolved so it's easier to see what's changed. In case that's gone wrong at all the original commit hash is 78c08f8914f4831ed3f14f4da550d1b0f26d8f37

[jrs65]

@ssiegelx this looks great. Just a few comments, the biggest being that Wiener is incorrectly spelled everywhere (currently Weiner is used)

[ssiegelx]

I think this PR is ready. I've addressed the comments and also add several features that I had in my local version, specifically: use the prior Dallas suggested for the Wiener ringmap maker, add capability to apply a window to shape synthesized beam, and fix a bug in how the number of RA samples was being calculated.

It was failing the ci tests because cora now requires python >= 3.7. I changed the GitHub workflow in the same way as was done in this commit. However, it still says that it need to run the tests for python 3.6 and 3.8, but never does. Any ideas how to fix this? Should draco's minimum python requirement in setup.py also be bumped up to 3.7?

[jrs65]

@ssiegelx this looks good to me. I think I beat you with merging a commit to get the CI working again, so you may need to rebase to get rid of any conflict.

I wondered about the minimum version in setup.py, in theory there's nothing about draco itself that can't support Python 3.6, it's just the dependency. I think maybe we should bump it though to try and keep things a bit more uniform across the repos.