.. zernpol documentation master file, created by sphinx-quickstart on Tue Oct 19 09:43:11 2021. You can adapt this file completely to your liking, but it should at least contain the root `toctree` directive. zernpol documentation ===================== .. toctree:: :maxdepth: 3 :caption: Contents This simple package offers tools to handle Zernike Polynoms and the different indexing systems (as the Noll system for instance). The package offers object and non-object oriented functions. Install ------- From pip: :: > pip install zernpol From sources: :: > git clone https://gricad-gitlab.univ-grenoble-alpes.fr/guieus/zernpol.git > cd zernpol > python setup.py install Usage ----- zernpol +++++++ The :func:`zernpol.zernpol` function is almost all whats one needs. It creates one or a list of :class:`zernpol.Zernpol` object which is basically a 2 tuple with some useful properties and methods from various input parameters. In the following all examples are equivalents: :: >>> from zernpol import zernpol, Zernpol >>> zernpol( (2,-2) ) Zernpol(2, -2) :: >>> Zernpol(2,-2) Zernpol(2, -2) :: >>> z = zernpol(5) # 5 is the index number system here in the default Noll system Zernpol(2, -2) :: >>> z = zernpol(5, 'Noll') Zernpol(2, -2) :: >>> z = zernpol(3, 'Ansi') Zernpol(2, -2) :: >>> z = Zernpol.from_noll(5) Zernpol(2, -2) :: >>> z = zernpol('astig_0') Zernpol(2, -2) The :func:`zernpol.zernpol` function accept also iterable, as list for instance, and will return a list of :class:`zernpol.Zernpol` : :: >>> from zernpol import zernpol >>> [z.name for z in zernpol( range(1,10) )] # default is Noll index system ['piston', 'tip', 'tilt', 'defocus', 'astig_0', 'astig_45', 'comax', 'comay', 'trifoll_0'] Indexing system +++++++++++++++ The used system can be change from a string or thanks to the :class:`zernpol.ZIS` (for Zernike Indexing system) containing functions to convert back and forward zernike polynom coefficients to indexes or indexes to zernike polynom coefficients :: >>> from zernpol import zernpol, ZIS >>> zernpol( range(1,5), ZIS.Noll) [Zernpol(0, 0), Zernpol(1, 1), Zernpol(1, -1), Zernpol(2, 0)] >>> zernpol( range(1,5), ZIS.Ansi) # zernpol( range(1,5), 'Ansi') also works [Zernpol(1, -1), Zernpol(1, 1), Zernpol(2, -2), Zernpol(2, 0)] Four Indexing systems exists (defined on `wikipedia `_): - :class:`zernpol.ZIS.Noll` - :class:`zernpol.ZIS.Ansi` - :class:`zernpol.ZIS.Fringe` - :class:`zernpol.ZIS.Wyant` These classes can also be simply used to convert index to zernike polynoms or the contrary :: >>> from zernpol import ZIS >>> ZIS.Noll.i2z(5) Zernpol(2, -2) >>> ZIS.Noll.z2i((2,-2)) 5 The default indexing system, used when the second optional argument of :func:`zernpol.zernpol` is None, is contained inside the :class:`zernpol.ZIS` as :class:`zernpol.ZIS.Default` parameters. One can change it with cautions at the highest level of an application: :: >>> from zernpol import ZIS, zernpol >>> ZIS.Default = ZIS.Ansi >>> zernpol(3) Zernpol(2, -2) :: >>> ZIS.Default = ZIS.Noll >>> zernpol(3) Zernpol(1, -1) Build the polynom +++++++++++++++++ The :meth:`zernpol.Zernpol.func` and :meth:`zernpol.Zernpol.func_cart` are used to build the zernike polynomials from polar and cartesian coordinate respectively. :: >>> from zernpol import zernrange, zernpol_func_cart >>> from matplotlib.pylab import plt >>> import numpy as np >>> X,Y = np.meshgrid( np.linspace(-1,1,50), np.linspace(-1,1,50)) >>> fig, axs = plt.subplots(3,4) >>> for zp,ax in zip(zernrange(1,13),axs.flat): >>> ax.imshow(zp.func_cart(X,Y)) >>> ax.axis('off') .. image:: img/zernikes_grid_example.png The :func:`zernpol.zernpol_func` and :func:`zernpol.zernpol_func_cart` function do the same but in a vectorialised fashion optimised for array of zernike polynome coefficients :: >>> from zernpol import zernpol, zernpol_func_cart >>> import numpy as np >>> X,Y = np.meshgrid( np.linspace(-1,1,50), np.linspace(-1,1,50)) >>> Z = zernpol_func_cart( zernpol(range(1,13)), X, Y) >>> Z.shape (12, 50, 50) :: >>> from zernpol import zernpol_func, zernrange >>> import numpy as np >>> from matplotlib.pylab import plt >>> r,theta = np.linspace(0,1,50), np.pi/4. >>> plt.plot(r, zernpol_func( zernrange(1,21), r,theta).T ) >>> plt.gca().set( title="20 first zernike modes at pi/4", xlabel="r") .. image:: img/zernikes_cutpi4_example.png .. note:: On the example above ``zernrange(1,21)`` is similar to ``zernpol(range(1,21))`` By default the zernike polynoms are normalised to rms = 1 :: >>> from zernpol import zernpol, zernpol_func_cart >>> import numpy as np >>> X,Y = np.meshgrid( np.linspace(-1,1,2000), np.linspace(-1,1,2000)) >>> Z = zernpol_func_cart( zernpol(range(1,13)), X, Y) >>> rms = np.sqrt( np.nanmean(Z**2, axis=(1,2)) ) >>> rms array([1. , 0.99998983, 0.99998983, 0.99997967, 1.00002123, 0.99993809, 0.99996952, 0.99996952, 0.9999695 , 0.9999695 , 0.99995938, 0.99989017]) Pupil +++++ For wavefront analysis one have to project the zernike polynomes over a specific portion of an 2d image. For this :mod:`zernpol` provides two things the :func:`zernpol.zernpol_pupil` (and its :meth:`zernpol.Zernpol.func_pupil` homolog) function and the :class:`zernpol.PupilMask`. :func:`zernpol.zernpol_pupil` is using :class:`zernpol.PupilMask` to project the zernike polynomes in an image. :class:`zernpol.PupilMask` has several parameters: - :attr:`zernpol.PupilMask.mask` the 2d boolean mask where True determine illuminated phase screen pixels - :attr:`zernpol.PupilMask.radius` is the true radius (in fraction of pixel) of the the pupil where zernike polynoms will be projected - :attr:`zernpol.PupilMask.center` true (x0,y0) coordinates of the pupil center - :attr:`zernpol.PupilMask.angle` (default 0.0) the angle of projection of zernike polynomials - :attr:`zernpol.PupilMask.flip` 2 tuple of 1 or -1 to determine any flip in the zernike polynoms for instance (-1,1) if a flip in x direction only, (1,1) is no flip at all. They are several ways to create a :class:`zernpol.PupilMask` and probably one will want to custom function to application (for instance from an Influence Function matrix). However :mod:`zernpol` offer basic ways : :: >>> from zernpol import PupilDisk, zernpol >>> from matplotlib.pylab import plt >>> disk = PupilDisk(36.0) # A disk with a 36mm radius (the unit does not mater, see below) >>> mask = disk.make_mask( [200,200], scale=0.34, center=[98,87] ) # scale is in User Unit/pixel >>> plt.imshow(zernpol(4).func_pupil( mask )) >>> plt.show() .. image:: img/zernikes_pupil_example.png The :func:`zernpol.func_pupil` is the vectorialized alternative : :: >>> from zernpol import PupilDisk, zernpol, zernpol_pupil >>> mask = PupilDisk(36.0).make_mask( [200,200], scale=0.34, center=[98,87] ) >>> zernpol_pupil( list(zernrange(1,15)), mask ).shape (14, 200, 200) The array returned can have a lot of NaN values when the image is large and the illuminated spot is small. The `inpupil_only` keyword allows to return only phases of illuminated pixel in a flat vector (per zernikes). This allow to reduce the memory footprint of a matrix of zernike polynoms and reduce computation time when multiplying or inverting matrix. Following the last example: :: >>> zernikes = zernpol_pupil( list(zernrange(1,15)), mask, inpupil_only=True ) >>> zernikes.shape (14, 8789) One can however reconstruct the images easily with the :meth:`zernpol.PupilMask.reconstruct` method :: >>> mask.reconstruct(zernikes).shape (14, 200, 200) The simple reverse operation can be done easily as well (to remove NaN of a mesured phase screen for instance). :: >>> mask.deconstruct( mask.reconstruct(zernikes) ).shape (14, 8789) A :class:`zernpol.PupilMask` can be build with :func:`zernpol.mask_from_image` and from a 2d image where guess of the projected zernike mode radius is made from the illuminated mask of the pupil and the its center from the baricenter of illuminated pixels. :: >>> import numpy as np >>> from zernpol import mask_from_image >>> x,y = np.meshgrid( np.arange(200), np.arange(200)) >>> z = np.random.random( (200,200) ) >>> z[ np.sqrt( (x-98.2)**2+(y-87.6)**2 )>52.3 ] = np.nan >>> mask = mask_from_image(z) >>> mask.radius, mask.center (51.75, (98.18354283054003, 87.58088919925513)) As seen above the pixel digitalisation will make the radius and center accurate at ~1/2 pixel The example bellow demonstrate the effect of the flip and the angle parameters: :: >>> from zernpol import PupilDisk, zernpol >>> from matplotlib.pylab import plt >>> import numpy as np >>> disk = PupilDisk(80) >>> ax = plt.subplot(221) >>> plt.imshow( zernpol(2).func_pupil( disk.make_mask([100,100], angle=0.0, flip=(1,1) ))) >>> ax.set_title("a=0.0 flip=(1,1)") >>> ax = plt.subplot(222) >>> plt.imshow( zernpol(2).func_pupil( disk.make_mask([100,100], angle=np.pi/4, flip=(1,1) ))) >>> ax.set_title("a=$pi/4$ flip=(1,1)") >>> ax = plt.subplot(223) >>> plt.imshow( zernpol(2).func_pupil( disk.make_mask([100,100], angle=0.0, flip=(-1,1) ))) >>> ax.set_title("a=0.0 flip=(-1,1)") >>> ax = plt.subplot(224) >>> plt.imshow( zernpol(2).func_pupil( disk.make_mask([100,100], angle=np.pi/4, flip=(1,-1) ))) >>> ax.set_title("a=$pi/4$ flip=(1,-1)") >>> plt.show() .. image:: img/zernikes_flip_example.png Misc features +++++++++++++ A latex equation representing the zernike polynomial can be reached : :: >>> from zernpol import zernpol >>> z = zernpol(5) >>> z.latex '\\sqrt{6} r^{2} \\ \\sin{\\,2\\theta}' Is equivalent to: :: >>> from zernpol import latex_formula >>> latex_formula((2, -2)) '\\sqrt{6} r^{2} \\ \\sin{\\,2\\theta}' :: >>> from zernpol import zernpol >>> from matplotlib.pylab import plt >>> for i,z in enumerate(zernpol(range(1,24))): >>> y = 24-i >>> plt.text(0.7 ,y, "$"+z.latex+"$") >>> plt.text(0.15,y, z.long_name) >>> plt.text(0.01,y, z) >>> plt.gca().set( ylim=(-1, i+1)) >>> plt.gca().axis('off') .. image:: img/zernikes_formula_example.png for the fun of it, a pyramide of zernike polynoms :: >>> from zernpol import zernpol_func_cart >>> from matplotlib.pylab import plt >>> x,y = np.meshgrid(np.linspace(-1,1,100), np.linspace(-1,1,100)) >>> ns = range(0,8) >>> f,axs = plt.subplots(8,16) >>> for a in axs.flat: a.axis('off') >>> for n in ns: >>> for m in range(-n,n+1,2): >>> axs[n, 7+m].imshow( zernpol_func_cart((n,m),x,y) ) >>> axs[n, 7+m].set_title((n,m)) >>> plt.gcf().set_size_inches( 14, 7 ) .. image:: img/zernikes_pyramid_example.png Module Indexes ============== See also Alphabetical :ref:`genindex` .. automodule:: zernpol :members: zernpol, Zernpol, ZIS, zernpol_pyramid, latex_formula, zernpol_norm, zernpol_func, zernpol_func_cart, zernrange, zernpol_to_noll, noll_to_zernpol, zernpol_to_ansi, ansi_to_zernpol, zernpol_to_fringe, fringe_to_zernpol, zernpol_to_wyant, wyant_to_zernpol, zernpol_info, cart2pol, pol2cart, PupilMask, PupilDisk, mask_from_image