Can we develop phase retrieval code that can work with less than Nyquist data? We'll clearly be unable to estimate high frequency abberations (high frequency zernike modes) but we should be able to estimate the lower order ones which could still be useful for characterizing commercial / core scopes with hard to change magnifications.
In [1]:
%pylab inline
import pyotf.otf
In [ ]:
class HanserPSFSubsample(object):
"""A class defining the pupil function and its closely related methods.
Based on the following work
[(1) Hanser, B. M.; Gustafsson, M. G. L.; Agard, D. A.; Sedat, J. W.
Phase-Retrieved Pupil Functions in Wide-Field Fluorescence Microscopy.
Journal of Microscopy 2004, 216 (1), 32–48.](dx.doi.org/10.1111/j.0022-2720.2004.01393.x)
[(2) Hanser, B. M.; Gustafsson, M. G. L.; Agard, D. A.; Sedat, J. W.
Phase Retrieval for High-Numerical-Aperture Optical Systems.
Optics Letters 2003, 28 (10), 801.](dx.doi.org/10.1364/OL.28.000801)
"""
def __init__(self, wl, na, ni, res, size, **kwargs):
"""zrange : array-like
An alternate way to specify the z range for the calculation
must be expressed in the same units as wavelength
"""
nyquist_res = 1 / (2 * self.na / self.wl) / 2
self.ratio = np.ceil(res / nyquist_res).astype(int)
self.res = res
self.size = size
self.nyquist = pyotf.otf.HanserPSF(wl, na, ni, res / ratio, size * ratio, **kwargs)
# do a general getter/setter that redirects to wrapped class for everything
# aside from res
def _gen_kr(self):
"""Internal utiltiy to generate coordinate system and other internal
parameters"""
k = fftfreq(self.size, self.res)
kxx, kyy = np.meshgrid(k, k)
self._kr, self._phi = cart2pol(kyy, kxx)
# kmag is the radius of the spherical shell of the OTF
self._kmag = self.ni / self.wl
# because the OTF only exists on a spherical shell we can calculate
# a kz value for any pair of kx and ky values
self._kz = psqrt(self._kmag**2 - self._kr**2)
def _gen_pupil(self):
"""Generate an ideal pupil"""
kr = self._kr
# define the diffraction limit
# remember we"re working with _coherent_ data _not_ intensity,
# so drop the factor of 2
diff_limit = self._na / self._wl
# return a circle of intensity 1 over the ideal passband of the
# objective make sure data is complex
return (kr < diff_limit).astype(complex)
def _calc_defocus(self):
"""Calculate the defocus to apply to the base pupil"""
kz = self._kz
return np.exp(2 * np.pi * 1j * kz *
self.zrange[:, np.newaxis, np.newaxis])
def _gen_psf(self, pupil_base=None):
"""An internal utility that generates the PSF
Kwargs
------
pupil_base : ndarray
provided so that phase retrieval algorithms can hook into this
method.
NOTE: that the internal state is created with fftfreq, which creates
_unshifted_ frequences"""
self._PSFa = self.nyquist.PSFa
# Because the _attribute_changed() method sets all the internal OTFs and
# PSFs None we can recalculate them only when needed
@property
def OTFa(self):
if self._OTFa is None:
self._OTFa = easy_fft(self.PSFa, axes=(1, 2, 3))
return self._OTFa
@property
def PSFa(self):
if self._PSFa is None:
self._gen_psf()
return self._PSFa
In [2]:
import scipy.ndimage as ndi
In [3]:
ndi.zoom?
In [ ]: