CC BY
64
COHERENT FIELD PHASE RETRIEVAL USING A PHASE ZERNIKE FILTER
V.V. Kotlyar, S.N. Khonina, V.A. Soifer, Y Wang *, D. Zhao*
Image Processing Systems Institute, Russian Academy of Sciences *) Beijing Institute of Technology, Beijing 100081, China, e-mail:ytwang@sun.ihep.ac.cn
Abstract
Aberrations of the coherent wavefront are analyzed using a phase Zemike filter. Developed iterative methods allow us to design a filter that decomposes the analyzed light field into a set of diffraction orders with amplitudes proportional to the circular Zernike polynomials. We also apply the algorithm to the calculation of the light field phase from measurements of the modules of decomposition coefficients. Operation of a 25-channel filter is simulated.
1. Introduction
Retrieval of the coherent light field phase is a topical problem in digital data processing. We cannot directly measure the light field phase but have to determine it indirectly, via light intensity measurements. For example, the wavefront of the light field can be reconstructed from an interferogram [1], from measurements of the intensity distribution of the spatial Fourier-spectrum [2]. A Shack-Hartmann wavefront sensor consisting of an array of equivalent pinholes or a microlens matrix [3] has also found use for phase retrieval. One can also reconstruct the phase using amplitude and phase filters capable of decomposing the light field in an orthogonal basis [4,5]. Paper [6] deals with the theory of spatial filters, called modans, intended for analyzing (selecting) the spatial modes of laser light. The number of basis terms that can provide the effective analysis of light fields is minimized using a basis matched to the field under analysis. For example, Gauss-Laguerre and Gauss-Hermite light modes are analyzed in [7] using phase spatial filters decomposing the light field in terms of orthogonal Laguerre and Hermite polynomials.
In the present paper, we analyze wavefront aberrations using a phase spatial filter. Based on an approach related to the expansion of field in an orthogonal basis [4-7] we propose that orthogonal circular Zernike polynomial [8] should be used as a basis. Note, however, that it turns out to be impossible to design a spatial filter that would be able to produce an intensity distribution proportional to the coefficients expansion of the sought light field phase in terms of Zernike polynomials. In other words, we cannot solve the formulated problem using a purely optical approach. We will have to do additional calculations.
In this paper we consider the decomposition of the entire complex amplitude, rather than the phase field, into a series in terms of Zernike polynomials. In this case, an intensity proportional to the field decomposition coefficients will be generated in the spatial Fourier-spectrum plane. Next, the measured coefficient modules are used for calculating the argument of the light field complex amplitude.
We also consider iterative algorithms for calculating the phase of a Zernike filter and for calculating the phase of the light field complex amplitude. Results of numerical simulation are discussed.
2. Circular Zernike polynomials
There is complete set of orthogonal functions with angular harmonics in the circle of radius r0. These are the circular Zernike polynomials [8]:
^nm (r,P) = AnmK ( r)exp( imP), (1)
where
A =
nm
In +1
(2)
(n-m) / 2
r:(r) = i (-1)p(n - p)b
p=0
n + m n - m
p!(—2— p)!(^2— *
n-2 p
(3)
x (-)
r0
r: (r) are the radial Zernike polynomials:
Rnm (r) = R:(r), R± (ro) =1,
Rll+1 (r) = 0,R22k+1 (r) = 0,|m < n,Ro0 (r) = 1,
(r, p) are polar coordinates.
The decomposition of the light complex amplitude E(r, p) into a series in terms of (1) is given by
E(r,p) =1 I Cnm ^nm (r,p)
0 r
j r: (r )r: (r )rdr=—°+- ,5
2(n +1)
np
• 0
Cnm = j j E(r, PY^lm (r, P)rdrdP .
(4)
(5)
(6)
0 0
In the plane of a spatial Fourier-spectrum which can be generated using a spherical lens with focal length f the light field complex amplitude F(p,d) will take the form
k
F (p^) —7X
2nf
r0 2n
. (7)
j j E(r,p)exp[ -i—rpcos(p-d) }-drdp 0 0 f where k=2n/A is the wavenumber of light, 2 is the wavelength of light , and (p, d) are polar coordinates. Based on Eq. (4), we can represent the decomposition of the light field in (7) into a series in terms of Zernike polynomials in Eq. (1) as:
k
F(p,d) = - I I (-i)mCnm X
J n=0m=-n 2n —
X Anme'm0 j R: (r ) Jm (- rp)rdr
o J
(8)
nr0
In deriving Eq. (8) we employed an integral representation of the Bessel function of first kind and of m-th order:
i m 2n
Jm (x) = — j exp[- ix cos t + imt\}t.
2n 0
We can find the integral in Eq. (8) in an explicit form [8]
X
X
n=0 m=-n
so n
* —
Wnm (p) = j R: (r) Jm (-frp)rdr ■■
o J
= ( i)(n-m)/2 „2 Jn+1(—f r0p)
=<_1) r0 (—f "V)
(9)
From Eq.(9) it is seen that at n>0 the complex amplitude is equal to zero at the central point p=0:
0,n > 0
Wnm (p = 0) = j 22 r02,n = 0.
(10)
Therefore, at n>0 the intensity distribution of basis diffraction orders in the Fourier plane will have a like-ring structure.
3. An algorithm for designing a phase Zernike filter
In Fig.1 is shown an optical setup (Zernike’s analyzer) demonstrating the use of a phase Zernike filter for analysis of the wavefront with complex amplitude E(r, p) . Analogous to the wavefront Shack-Hartmann sensor [3], the Zernike filter (ZF) is placed immediately in the plane of the wavefront under study. A spherical lens L of focal length f is place beside it. In the rear focal plane of the lens L we place an array of photoreceiver PA coupled to the computer PC. The transmission function of the ZF should be the phase one
r(r, p) = exp[iS(r, p)] (11)
and is sought-for in the form
N n
T(r,p) =I I^nm (r,p) X
n=0m=-n
X exp|i—f [nm cos(p- 6nm ) + Vn
(12)
where (prm, pm are the vectors of carrier spatial frequencies in polar coordinates and vnm are the task’s free parameters that should be fitted in such a manner that Eq. (12) be the correct equality.
ZF L
Fig.1 An optical configuration of a Zernike analyzer: ZF is a Zernike filter, L is a spherical lens, PA is an array of photosensors, and PC is a computer.
Provided the effective spatial separation in the Fourier plane of separate basis diffraction orders, we can consider the functions
^m (r,p)exp[ikf Vn
X cos(p- 6nm ) + Vnm ]
(13)
as being nearly orthogonal and compute the parameters in Eq. (12) by
|r0 2n
j j exp[(r,p)]nm (r,p)X
0 0 .
X exp[- ikf ^rpnm cos(p- dnm )rdrdp]}
The proposed iterative algorithm for the computation of the filter phase function S(r, 6) is based on the consecutive computation of the sums (12) and the integrals (13) using the FFT algorithm, with some constraints imposed. The estimate of the ZF phase function in the (k+1)-th step has the form:
Sk+1 (^ p) = argi I I ^m (r, p) :
rpnm cos(p- 6nm ) + V
(14)
(k) nm
where Vnm is the estimate of free parameters in the k ■
(15)
th iteration, with the parameters being deduced using Eq. (13) from the preceding phase estimate Sk (r,p). For the convergence of similar iterative algorithms the reader is referred to [7].
The Zernike filter of Eq. (12) produces the spatial separation of certain coefficients Cnm from the expansion of the field E(r, p) in Eq. (4). If such a filter r(x,y) is placed near a spherical lens and is illuminated by a light wave of amplitude E(x,y), the light intensity at the focal plane points u = anm and v = Pnm will be approximately proportional to the squared modulus of the expansion coefficients Cnm:
—f-1 j j E(x, yT(x, y) x
-o -o
x exp[- ikf(ux + vy)]xdy «
N n
~I I Cnm exp(i Vnm ) X
n=0m=-n
X5(u -anm , V - P nm )
where (anm , Pnm) are the vectors of spatial carrier frequencies in Cartesian coordinates. The more effectively is performed the spatial separation of basis beams propagating at different angles, the more accurate will be the approximate equality (15).
4. An algorithm of the light field phase retrieval
After the light intensity proportional to the squared modules of coefficients of the expansion (4)
i = |c I2 (16)
nm | nm | v '
has been measured at discrete points of the Fourier-plane (Fig.1), we will have to do additional calculations to find the light field phase
Q(r ,p) = arg E (r,p). (17)
For this purpose, we may use an iterative algorithm similar to the algorithm of Eqs.(13) and (14). In
r,
n=0m=-n
this case, the estimate of the light field phase in the (k+1)-th iteration will be given by
q—+1(r ,p)=arg^I
Nn
(18)
xT,
,(r, p)exp[ rim]}
,(k)
where Y nm are the free parameters in the k -th iteration obtained from the equation
Y(m = arg| j jexp[iQk (r,p)]
0 0
(19)
xT,
t(r ,p)rdrdp}
where Qk (r,p) is the desired phase estimation the k -th iteration.
5. Peculiarities of analysing aberrations using a Zernike filter
Since wavefront aberrations met with in optical systems are described by even functions in an azimuth angle p [8] , we can represent the wave field E(r, p) in the form
f N n ]
E(r,p) = expji I I BnmRm (r)cos(mp)J. (20)
L n=0m=-n J
Because of this, instead of the general expansion of Eq.(4) we may use an expansion in even functions
E(r,p) =I I Cnm T nm (r,p).
Tnm (r,p) = £
n + 1
m^—TR:(r)cos(mp).
(21)
(22)
where sm = {2,m ^ 0;1,m = 0} .
With small aberrations, the relation between the decomposition coefficients Bnm and Cnm is linear:
1 + iB00 =
C0,
iB = r
Cn
(23)
With arbitrary aberrations, the relation between Bnm and Cnm is nonlinear and after measuring the modules |C nm | we need to use the algorithm of Eqs. (18)
and (19) in order that we may the phase Q(r,p) in Eq. (17). After this, we can find the coefficients of wave aberrations Bmn using Eq. (20).
Note that because r0° (r) = 1 the Zernike polynomial basis contains a unit as a term, which means that if we illuminate a Zernike filter by a plane wave of amplitude E(r,p)=1, we find that in the frequency plane only one coefficient of the expansion (4) will be non-zero:
|C001 ^ 0 .
From Eq. (9) it also follows that diffraction orders corresponding to the basis functions with different numbers m , but with the same numbers n , will have similar diffraction patterns (of like-ring structure at n > 0 ) in the Fourier-plane:
\Wnm (p)| = r0
|Jn+1(kf-Vp p)| (kf ~'r0p)
(24)
6. Numerical examples
The simulation parameters are: 256 pixels on radius r and 256 pixels on angle 6, ro = 1 mm, k = 104 mm-1, f = 100 mm. We designed a 25-channel Zernike filter forming the basis diffraction orders (n,m) at m<8 and n < 8 that propagate at some angles to the optical axis. Shown in Fig. 2 are: (a) a half-tonic phase of the Zernike filter (black colour - phase 0 and white colour -phase 2n), (b) 25 diffraction orders generated in the lens frequency plane (negative), (c) the configuration of the numbers (n,m) distribution between orders. The filter was assumed to be illuminated by a plane wave. In this case, the analyzer splits the incident plane beam into 25 beams having nearly equivalent energies. All diffraction orders account for over 80% of the entire illuminating beam energy. From Fig.2b it is seen that the zero intensity is found in the central points of all orders in the Fourier-plane, expect for the zero order (0,0). This means that the illuminated wavefront has not aberrations.
Figure 3 shows the result of operation of the same 25-channel Zernike filter illuminated by a beam consisting of three terms with the same weight and the numbers (n,m): (2,0)+(5,3)+(7,7). In Fig.3 are shown: (a) the illuminating beam intensity, (b) the filter phase, (c) the diffraction pattern in the Fourier-plane. Checking the table of order numbers (Fig.2c) from Fig.3c it is seen that the intensity is non-zero (black points in Fig. 3 c) at the central points of the diffraction orders with numbers (2,0), (5,3), and (7,7). The Table gives relative intensity values in the vicinity of central points for all orders (the size of the vicinity is 3x3 pixels). As can be seen from the Table, if weight coefficients of the terms entering into the illuminating beam amplitude are equal to each other at the Zernike analyzer (Fig.1) input, they have different values at the output:
|C20|2 = 0.82,| C53|2 = 1,| C77|2 = 0.68.
In addition, the coefficients are also non-zero in the other measurement channels. Note that the greatest values among the coefficients, which were, absent at the input but are present at the output of the Zernike analyzer (Fig.1) are as follows:
|C64|2 = 0.22,| C66|2 = 0.21,| C00I2 =
= 0.22,| C42|2 = 0.16,| C73|2 = 0.15
X
n=0m=-n
n=0 m=-n
b)
7 , 1 6 , 6 6 4 6 , 2 6 0
7 , 3 3 , 1 2 2 2 , o 5 5
7 , 5 3 , 3 0 0 1 , 1 5 3
7 , 7 4 , 0 4 2 4 , 4 5 , 1
8 ,0 8 ,2 8 4 8 ■ & 8 8
c)
Fig.2 (a) Half-tonic phase of a 25-channel Zernike filter, (b) diffraction pattern in the frequency Fourier
plane under illumination by a plane wave, and (c) a
table of correspondence between diffraction orders and the numbers of Zernike polynomials.
Table .
Relative intensity in the vicinities of central points for basis diffraction orders in the Fourier plane when the Zernike filter (Fig.3b) is illuminated by a superposition of three beams: (2,0)+(5,3)+(7,7)
(7,1) 0.02 (6,6) 0.21 (6,4) 0.22 (6,2) 0.01 (6,0) 0.00
(7,3) 0.15 (3,1) 0.02 (2,2) 0.03 (2,0) 0.82 (5,5) 0.04
(7,5) 0.03 (3,3) 0.09 (0,0) 0.22 (1,1) 0.05 (5,3) 1.00
(7,7) 0.68 (4,0) 0.10 (4,2) 0.16 (4,4) 0.00 (5,1) 0.02
(8,0) 0.02 (8,2) 0.00 (8,4) 0.04 (8,6) 0.05 (8,8) 0.12
b)
c)
Fig.3 (a) The input light field amplitude distribution for a Zernike analyzer consisting of three terms numbered (2,0)+(5,3)+(7,7), (b) the phase of a 25-channel Zernike filter corresponds to that shown in Fig.2a, (c) the diffraction pattern in the frequency Fourier plane at the analyzer output.
Thus, we can see that at the Zernike analyzer output the squared modules of coefficients are measured with a relative error of about 20%. We can reduce this error in two ways: via the improvement of calculating the Zernike filter’s phase ( note that numerically simulated squared modules of the functions (24) in the Fourier plane, see Fig.2b, are different from the reference (24) by 20% on the average) and via the reduction of the vicinity magnitude (3x3 pixels) in the diffraction order centres in which the coefficient modules are measured.
7. Conclusions
In this paper we have obtained the following results:
• an iterative algorithm for the calculating phase functions of Zernike filters has been developed. The ZF are intended for the analysis of wavefront aberrations and
are diffractive optical elements that produce in the frequency plane the spatial separation of different circular Zernike polynomials contained in an incident light field;
• an iterative algorithm has been produced and is able to reconstruct the light field phase from measurements of a small number of intensity pixels at the Zernike analyzer output, with the pixels being proportional to the squared module of coefficients in the expansion of a given field amplitude in terms of circular Zernike polynomials;
• operation of a 25-channel Zernike filter has been numerical simulated.
Acknowledgments
This work was supported by Russian Foundation of Fundamental Research (N96-01-10021).
References
1. T. Yatagai, S. Nakadata, M. Idesawa, H.Saito. Automatic fringe analysis using digital image processing techniques. Opt. Eng., 1982, v.21, no.2, pp.432-435.
2. J.R. Fienup. Phase retrieval algorithms. A comparison. Appl. Opt., 1982, v.21, no. 15, pp.2758-2769.
3. G. Artzner. Microlens arrays for Shack-Hartmann wavefront sensors. Opt. Eng., 1992, v.31, no.6, pp. 1311-1322.
4. M.A. Golub, A.M. Prokhorov, I.N. Sisakian, V.A. Soifer. Synthesis of spatial filter for investigation of the transverse mode composition of coherent radiation. Sov. J. Quant. Electr,, 1982, v.12, no.9, pp.1208-1209.
5. V.V Kotlyar. Decomposition of coherent field by orthogonal basis. Computer Optics, MCNTI, Moscow, 1989, issue 5, pp.31-33.
6. V.A. Soifer, M.A. Golub. Laser beam mode selection by computer generated holograms, CRC Press, Boca Raton, 1994.
7. V.A. Soifer, V.V. Kotlyar, L.L. Doskolovich. Iterative methods for diffractive optical elements computation, Taylor and Francis, London, 1997.
8. Born M., Wolf E., Principles of Optics, Pergamon Press, Oxford, 1968.
