brandondube / prysm

physical optics: integrated modeling, phase retrieval, segmented systems, polynomials and fitting, sequential raytracing...
https://prysm.readthedocs.io/en/stable/
MIT License
257 stars 44 forks source link

Zernikes? #22

Closed ljpguimen closed 4 years ago

ljpguimen commented 4 years ago

Hi Brandon,

Thanks for sharing this nice code. I have an issue when transforming some random phase wavefront data into Zernike polynomials:

x, y, phase = np.linspace(-1,1,128), np.linspace(-1,1,128), 10*np.random.rand(128,128) fz = FringeZernike(x=x, y=y, phase=phase, z_unit='um') fz.magnitudes

The outputs are all 0. I'm expecting some positive numbers for this random phase - where am I doing wrong?

brandondube commented 4 years ago

Are you trying to fit zernikes to random phase, or make a surface out of random zernike coefficients?

To fit random phase:

import numpy as np

from prysm.zernike import zernikefit
from prysm.geometry import circle

# you need to punch out the corners or the zernike polynomials will not be orthogonal
# prysm does LSQ fitting so it doesn't really care, and with a random phase screen
# the output itself is probably not meaningful, but to be complete...
z = np.random.rand(128,128)
mask = circle(128)
z[mask == 0] = np.nan

# just coefs = if residual=False (the default)
# these are 'orthonormalized' zernikes with meaning of coefficients are RMS
# value of that term, if you want the "unnormalized" ones, use norm=False
coefs, residual = zernikefit(z, terms=47, norm=True, residual=True, map_='fringe')

To generate phase from random zernike terms:

import numpy as np

from prysm import NollZernike, FringeZernike, ANSI1TermZernike

nz = NollZernike(np.random.rand(36))  # 36 terms

fz = NollZernike(np.random.rand(48))  # 48 terms

ansi = ANSI1TermZernike(np.random.rand(1000), samples=512) # 1000 terms -- the higher radial order will require more samples to satisfy Nyquist requirements

The .magnitudes property calculates the magnitude of each Zernike term, which is the L2 norm or RSS of terms for which radial order n is the same and azimuthal order m has the same absolute value. If you just want the coefficient vector, it will be stored in .coefs, which will be zeros in that case. Specifying the z kwarg will set the phase (.plot2d() to see, for example) but you are actually just bypassing the Zernike classes altogether and using arguments to the prysm.Pupil constructor.

Did you find the User's guide helpful? Actually looking at it now, I can see it is outdated and actually has some things wrong, and does not even include the ANSI Zernike implementation. Perhaps you could open a PR to improve it?

ljpguimen commented 4 years ago

Thanks! I'm fitting Zernikes to a specific phase. Your first section solves the problem.

I'll open a PR once I get a chance.

brandondube commented 4 years ago

There are more parameters than I showed, please see help(zernikefit) or find it in the API reference. In particular, you can fairly elegantly change the normalization radius by modifying the radial coordinates you give it. Or handle shifted coordinate frames, and so on.

because it is a LSQ fit and not something like Gram-Schmidt, the orthogonality remains over a unit disk. I'm going to close this issue for now, please re-open if you see the need.