andykee / lentil

Heart-healthy physical optics
https://andykee.github.io/lentil/
Other
14 stars 6 forks source link

Automatically interpolate phasor components to avoid wrapping #2

Closed andykee closed 2 years ago

andykee commented 4 years ago

Carl is having a problem where he's performing a propagation with tilt='angle' and npix_chip=256 but the Pupil phase has shape = (256,256). The pupil has some moderately large aberrations that are wrapping around.

It would be desirable to automatically detect when this will happen and interpolate the necessary arrays to prevent it. It should be possible to control this interpolation behavior, at least in the sense of turning it on or off so we'll need to add a new parameter to propagate(). Let's tentatively call this parameter interpolate_phasor. This parameter should default to True.

The wraparound period of the DFT is given by 1/alpha, which we should be able to compute for all Pupil to Detector propagations at the same time we cache the amplitude and phase terms in Propagate.__enter__(). I think we'll also have to (at least temporarily) increase the mask size to ensure the ptt_vector is computed correctly. Because alpha is wavelength-dependent, we'll have to use the worst case value.

Once FFT propagations are working, we'll be able to use this approach to make sure they are being performed at Nyquist or better.

In the event that interpolate_phasor = True and the wraparound period is less than npix or npix_chip (plus some buffer), we should compute the size of the amplitude, phase, and mask arrays to avoid wrapping issues and interpolate away!

We'll also have to decide on which of scipy.interp2d's interpolation methods to use. It defaults to linear but we likely want cubic at least for the phase.

andykee commented 3 years ago

Rather than forcing a user to accept no wraparound at all, it would be nicer if we could instead allow them to select how much wraparound they're willing to tolerate. In some cases, they may not care about having perfectly no wraparound since the high frequency content that is wrapping first may be drowned out by other effects like detector noise. Instead of interpolate_phasor = <bool>, maybe something like min_q = <float> would be a better parameter name. If it defaults to 0, the primary behavior would be to do nothing.

andykee commented 2 years ago

Closing in preparation for release of v0.6.0