- Email: [email protected]
Reconstruction of two-dimensional complex amplitudes from intensity measurements
Optics Communications 225 (2003) 19–30 www.elsevier.com/locate/optcom
Reconstruction of two-dimensional complex amplitudes from intensity measurements X. Liu *, K.-H. Brenner Chair of Optoelectronics, University of Mannheim, Mannheim, Germany Received 23 January 2003; received in revised form 10 July 2003; accepted 20 July 2003
Abstract The amplitude and phase recovery of optical fields with tomographic methods employing Wigner- or ambiguity functions (WD/AF) has been demonstrated extensively for the case of one-dimensional (1D) functions. For 2D light distributions, the associated WD/AF is 4D, posing several problems. In this paper we introduce a new concept, which allows to reconstruct arbitrary 2D distributions using only 1D measurements. We reconstruct one dimension (y) of the complex light source for each position of the x-dimension. To this end, we realized a 1D propagation operator. The corresponding optical setup for measurement is shown and experimental results are presented. Ó 2003 Elsevier B.V. All rights reserved. PACS: 42.30.Rx; 42.30.wb; 42.25.Bs; 42.60.Jf Keywords: Phase retrieval; Image reconstruction; Tomography; Wave propagation; Beam characteristics
1. Introduction The utilization of intensity measurements for amplitude and phase recovery of light distributions is of considerable interest in several areas of optics. Tomographic methods employing Wigneror ambiguity functions (WD/AF) have been demonstrated extensively for the case of one-dimensional (1D) functions. For the rotation of phase space, fractional Fourier transform meth-
*
Corresponding author. Tel.: +496211812693; fax: +496211812695. E-mail address: [email protected] (X. Liu).
ods were first employed in [1]. The ambiguity function (AF) is closely related to the WDF through a 2D Fourier transform. A review about the WDF and its applications in optics and optoelectronics can be found in [2]. Previously, the WDF was used for the experimental determination of amplitude and phase of either 1D optical beams [3,4] or optical pulses [5]. The AF was used in the past as a polar display of the optical transfer function [6,7], for the description of timevarying random diffracted field [8] or as a tool for designing optical systems with high focal depth [9]. Recently, it was shown theoretically that the AF can be a more straightforward and a much faster tool for amplitude and phase retrieval than
0030-4018/$ - see front matter Ó 2003 Elsevier B.V. All rights reserved. doi:10.1016/j.optcom.2003.07.031
20
X. Liu, K.-H. Brenner / Optics Communications 225 (2003) 19–30
the WDF [10], since the retrieval algorithm is based only on Fourier transforms while in the case of the WDF a filtered back-projection Radon transform is needed, which very often is corrupted by singularities and branch points. The reconstruction of complex amplitudes with the ambiguity function was demonstrated for 1D signals in [11], suggesting that an extension to two dimensions is straightforward. Such extensions have been discussed in [3] by numerical simulation and in [11] on a proof of principle basis. Raymer et al. [3] have shown that a suitable combination of propagation distances and lens focal lengths is sufficient to recover the full 4D Wigner function. Likewise, [11] requires a set of quadratic media with different coefficients in xand y-direction. In Section 2 we point out, that a direct extension to two dimensions generally requires two independent degrees of freedom in order to measure a 4D distribution. For the reconstruction method proposed here, a practical and simple optical setup is presented. The number of measurements grows only linearly with the size of the object. The method is to our knowledge, the first practical method to reconstruct 2D wave fields using only 1D measurements. In this paper, we use the propagation properties of the 1D Ambiguity function to reconstruct a 2D complex distribution. To this end, we have experimentally realized a 1D propagation operator. The corresponding optical setup for the measurements is described and first experimental results are presented.
Aðx0 ; m0 ; uÞ ¼
u~ðmÞ ¼ uðxÞ ¼
Z Z
Z
m0 m0 2pimx0 dm; u~ m e u~ m þ 2 2 ð2Þ
uðxÞe2pixm dx
and ð3Þ
u~ðmÞe2pixm dm:
The AF is a phase-space representation of the complex signal uðxÞ. For the tomographic reconstruction of a complex signal uðxÞ from the AF, two main properties of the AF are used: the line distribution Að0; m0 ; uÞ corresponds to the FT of the intensity distribution of u. Furthermore, light propagation in the AF representation corresponds to a shear of the AF. Thus the FT of the intensity 2 juz ðxÞj at different locations along the propagation ðzÞ direction, corresponds to a line in the AF domain, Aðkzm0 ; m0 ; uÞ (Fig. 1). Therefore, we can reconstruct the AF of a complex signal uðxÞ from the samples of the intensity distribution taken at a set of different z locations (Fig. 2). Finally from Aðx0 ; m0 ; uÞ we can compute the complex signal uðxÞ. We denote u0 ðxÞ as the complex amplitude distribution at z ¼ 0. The propagation to z > 0 in the Fresnel-approximation can be treated by multiplication with a quadratic phase factor in the Fourier-domain (Eq. (4)) 2
u~z ðmÞ ¼ u~0 ðmÞeipkzm ;
ð4Þ
2. Theory of reconstruction for complex 1D signals The AF Aðx0 ; m0 ; uÞ of a complex 1D function uðxÞ is defined equivalently by Eq. (1) or Eq. (2), with the Fourier transform (FT) of uðxÞ written as u~ðmÞ. Because there are different definitions for the FT, the definition, used in this paper, is given in Eq. (3) 0
0
Aðx ; m ; uÞ ¼
Z
x0 x0 2pixm0 u xþ dx; u x e 2 2 ð1Þ
Fig. 1. (Left) juz ðxÞj2 at different locations along the propagation (z) direction. (Top right) Aðx0 ; m0 ; uÞ. (Bottom right) R 0 Aðkzm0 ; m0 ; uÞ ¼ Iz ðxÞe2pixm dx.
X. Liu, K.-H. Brenner / Optics Communications 225 (2003) 19–30
21
Fig. 2. Reconstruction of Aðx0 ; m0 ; u0 ðxÞÞ from Iz ðxÞ.
k is the wavelength of the given light source. Combining Eqs. (2) and (4), the relation between the AF of uz ðxÞ Aðx0 ; m0 ; uz Þ and the AF of u0 ðxÞ Aðx0 ; m0 ; u0 Þ is obtained Aðx0 ; m0 ; uz Þ ¼ Aðx0 kzm0 ; m0 ; u0 Þ:
ð5Þ
Fig. 3. (a) Aðx0 ; m0 ; uÞ; (b) 1. half of Aðx0 ; m0 ; uÞ; (c) 2. half of Aðx0 ; m0 ; uÞ; (d) 1. half of Aðx0 ; m0 ; u~Þ.
uz ðxÞuz ðxÞ
is related to the AF by The intensity Z 0 uz ðxÞuz ðxÞe2pixm dx ¼ Að0; m0 ; uz Þ ¼ Aðkzm0 ; m0 ; u0 Þ:
ð6Þ 2
This indicates that the FT of the intensity juz ðxÞj is a line along in the AF of the signal u0 ðxÞ. In order to obtain the complete Aðx0 ; m0 ; u0 Þ, we have to 2 implement and measure a set of intensities juz ðxÞj at different positions z. The complex amplitude u0 ðxÞ can finally be computed from Aðx0 ; m0 ; u0 Þ using the definition of the AF: Z 0 0 0 u0 ðx Þu0 ð0Þ ¼ Aðx0 ; m0 ; u0 Þepix m dm0 ; R ð7Þ 0 0 Aðx0 ; m0 ; u0 Þepix m dm0 i arg½u0 ð0Þ 0 e u0 ðx Þ ¼ : ju0 ð0Þj Since only the complex amplitude distribution is of interest, the scaling of the absolute value, ju0 ð0Þj, and the phase shift, ei arg½u0 ð0Þ , are usually irrelevant. It is, however, required that ju0 ð0Þj is nonzero in this case! By applying suitable shifts, this restriction can be avoided. In practical propagation situations, we can only 2 measure juz ðxÞj for locations z ¼ zmax , corresponding to a partial angular section of the AF. For the measurement of the remaining section of Aðx0 ; m0 ; u0 Þ we can use another property of the AF Aðm0 ; x0 ; uÞ ¼ Aðx0 ; m0 ; u~Þ:
ð8Þ
In other words, the FT results in a 90° rotation 2 for the AF, allowing to measure j~ u0 ðxÞj instead of 2 ju1 ðxÞj , as illustrated in Fig. 3.
2.1. Extension to the reconstruction of 2D signals For the reconstruction of a 2D complex signal, there are two alternative routes. One alternative is the direct extension of the 1D theory to a 2D theory, i.e., we use the definition and properties of the AF for a 2D distribution. The another alternative is to find a suitable experimental setup, by which the 1D theory of reconstruction can be applied. First we explain the extension of the 1D theory to a 2D theory and we will analyze its requirements. The definition of the AF for a 2D distribution shows that the AF is a 4D function Aðx0 ; y 0 ; m0 ; l0 ; uðx; yÞÞ Z Z x0 y0 x0 y0 ¼ u x þ ;y þ u x ;y 2 2 2 2 0
0
e2piðxm þyl Þ dx dy; 0
0
0
ð9Þ
0
Aðx ; y ; m ; l ; uðx; yÞÞ Z Z m0 l0 m0 l0 ¼ u~ m þ ; l þ u~ m ; l 2 2 2 2 0
0
e2piðx mþy lÞ dm dl; Z Z uðx; yÞ ¼ u~ðm; lÞe2piðxmþylÞ dm dl and Z Z u~ðm; lÞ ¼ uðx; yÞe2piðxmþylÞ dx dy:
ð10Þ
ð11Þ
The following properties of the 4D AF are a direct counterpart to Eqs. (5)–(7).
22
X. Liu, K.-H. Brenner / Optics Communications 225 (2003) 19–30
Propagation corresponds to a shear along the x0 and y 0 axis Aðx0 ; y 0 ; m0 ; l0 ; uz ðx; yÞÞ ¼ Aðx0 kzm0 ; y 0 kzl0 ; m0 ; l0 ; u0 ðx; yÞÞ:
ð12Þ
The 2D intensity is related to the 4D AF by Z Z 0 0 uz ðx; yÞuz ðx; yÞe2piðm xþl yÞ dx dy
3. Experimental setup for 1D light propagation
¼ Að0; 0; m0 ; l0 ; uz ðx; yÞÞ ¼ Aðkzm0 ; kzl0 ; m0 ; l0 ; uz ðx; yÞÞ:
The reconstruction steps have to be repeated for a set of different but fixed x-positions. By selecting a fixed x-position, say x0 , we are dealing with a 1D complex amplitude u0 ðyÞ. This experimental setup for this operation will be described in the next section in detail.
ð13Þ
Reconstruction of uðx; yÞ from Aðx0 ; y 0 ; m0 ; l0 ; uðx; yÞÞ RR 0 0 0 0 Aðx0 ; y 0 ; m0 ; l0 ; uðx; yÞÞepiðx m þy l Þ dm0 dl0 0 0 : uðx ; y Þ ¼ u ð0; 0Þ ð14Þ As shown in Eq. (14), the reconstruction of uðx; yÞ requires the full 4D AF. Since the shear kz in Eq. (13) is the same for the m- and l-direction, a set of different z locations is not sufficient for the reconstruction. The intensity distribution at each z location contains 2D information and the intensity measurement along the z-axis thus provides only 3D information. To reconstruct the AF of general 2D signals, an additional dimension of information is necessary. In order to acquire the full 4D information of the AF, either special symmetries are required [12], or a second dimension such as, e.g., a variable focal length is necessary. While a variation of distance is easy to implement, a setup with lenses of variable focal length is more difficult to realize. Consequently, in [3] the second degree of freedom was achieved by varying the focal length discretely in a few steps, and the distance between lenses was varied in many small steps. In this paper we introduce a new concept, which allows to reconstruct arbitrary 2D distributions using only 1D measurements. Instead of attempting to reconstruct a 4D AF, we consider the 2D wave field as a discrete series of 1D wave fields and we apply the 1D theory of reconstruction from the 2D AF. This method, however, requires to implement a 1D propagation operator, i.e., an operator which performs an exact imaging operation for one direction (x-direction) and a 1D propagation operation for the other direction (y-direction).
In the theoretical description, we require a 1D propagation operator. Such an operator would have to perform an imaging operation for one axis and an independent propagation operation for the other axis. We use the y-axis for the transformation and the x-axis as a series of parallel channels (columns), since this choice corresponds to our optical implementation. Assuming a distribution u0 ðx; yÞ as input and a control parameter z, the 1D propagation operator would have to implement this operation: Z 1 0 2 uz ðx; yÞ ¼ pffiffiffiffiffi e3pi=4 e2piz=k u0 ðx; y 0 Þeip=kzðyy Þ dy 0 : kz ð15Þ 3.1. 1D propagation operator The setup in Fig. 4 implements the desired 1D propagation operator. In the reference plane at z ¼ 0, the input u0 ðx; yÞ is located. The series of positions zP ¼ zn ¼ n dz is our control parameter. In the detector plane, the output is obtained. The goal is that, by the changing zn we can achieve propagation only for the y-direction, while for the x-axis a one-to-one mapping is performed. To analyze the optical setup, we can
Fig. 4. 1D propagation operator for 2D input distribution.
X. Liu, K.-H. Brenner / Optics Communications 225 (2003) 19–30
observe the x- and y-axis (Figs. 5 and 6) separately. This is possible, since the paraxial propagation and the optical components are x–y-separable. For the input, there is no restriction. First consider Fig. 5 for the x-axis. It performs a one-to-one mapping for every point at the input. For this system, the distance between two lenses can be varied without any effect on the intensity distribution in the output plane. Fig. 6 illustrates the operation for the y-axis. At z ¼ 4fy we have a one-to-one mapping, because this system corresponds to the classic imaging system with unit magnification. By adding a separation zP , we can implement the desired 1D propagation operator. The performance of the setup is best described with the formalism of ray-transfer matrices. As has been described in several publications [7], the action of a linear optical system in the ambiguity representation is a simple phase space coordinate transformation. Because all the operations involved (lens, propagation) are separable, we can simplify the description in terms of two 1D transformations. If x and m are the conjugate coordinate pair for the x-direction and y and l the conjugate coordinate pair for the y-direction with
Fig. 5. x view of the 1D propagation operator.
1 kx m ¼ l 2p ky
23
ð16Þ
a 1D propagation is described by the matrix 1 kz ð17Þ 0 1 and transmission through a lens is described by the matrix 1 0 ð18Þ kf1 1 acting on the phase space vector x or ly , rem spectively. Using these matrices, the optical system is described by x-axis: 1 0 Tx ¼ 4fy 4fx þzP 1 ; ð19Þ kfx2
y-axis: 1 Ty ¼ 0
kzP : 1
ð20Þ
The transformation matrix Tx is basically an identity matrix. The off-diagonal element ð4fy 4fx þ zÞ= kfx2 has only an effect on the phase. Since we measure intensity, this component has no effect. The transformation matrix Ty according to Eq. (20) is the 1D propagation matrix with inverted coordinates. Thus the setup in Fig. 4 implements a 1D propagation operator. For practical reasons, the propagation distance is limited to ð2fy fx Þ. Therefore only a finite angular segment can be measured with this system. In order to obtain a larger measuring distance, we use fx ¼ 40 mm and fy ¼ 80 mm . With this choice the measuring range is jzj 6 100 mm. In this experiment we want to reconstruct the discrete AF with dx ¼ 11 lm and dm ¼ 1=ðN dxÞ. In a diagram with k ¼ x=dx and n ¼ m=dm a longitudinal shift zmax ¼ 100 mm corresponds to a shear angle of 45° for 512 512 points. 3.2. 1D FT and propagation operator
Fig. 6. y view of the 1D propagation operator.
The propagation operator described above, can measure the intensity only for distances jzj < zmax , which in the ambiguity domain corresponds to an
24
X. Liu, K.-H. Brenner / Optics Communications 225 (2003) 19–30
angle segment of the full 2D space. In order to be able to measure the remaining angle segment, we developed a 1D FT and propagation operator, i.e., we first perform a 1D Fourier transform and then a 1D propagation operator over a distance zF . In the optical system, we require this propagation distance zF to be adjustable. Assuming a 2D complex amplitude u0 ðx; yÞ as input and the propagation distance zF as a control parameter, we have to implement an optical system which implements a one-to-one mapping for the x-axis and a sequence of Fourier transform and adjustable propagation distance for the y-axis. Fig. 7 shows a system, which realizes this requirement. To discuss its operation, again we can consider the x- and y-axis separately. Actually for the x-axis (Fig. 8) this system is identical to the 1D propagation operator. For the y-axis (Fig. 9) the system is basically a 2f -Fourier transformer, since the lens focal length and the distances fy are equal. At z ¼ 2fy we obtain u~ðx=kfy Þ. The additional separation zF allows us to combine the Fourier transform with a propagation operator. The system can also be described in terms of ray-transfer matrices. x-axis: 1 0 Tx ¼ 2fy 4fx þzF 1 ; ð21Þ kfx2
y-axis: 0 kfy 1 kzF Ty ¼ 1 0 : 0 1 kfy |fflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflffl} |fflfflfflfflfflfflffl{zfflfflfflfflfflfflffl} propagation
ð22Þ
FT
Fig. 7. 1D FT and propagation operator for 2D input distribution.
Fig. 8. x view of the 1D FT and propagation operator.
Fig. 9. y view of the 1D FT and propagation operator.
As designed, the system performs an intensity image for the x-axis and a sequence of Fourier transformation and propagation for the y-axis. For practical reasons, the propagation distance is limited to ðfy fx Þ. In order to match the Fourier-segment with the segment, obtained using the propagation operator, we use fx ¼ 10 mm and fy ¼ 40 mm, dx and dm remain the same value as for the 1D propagation operator. The intensity measuring range is between z ¼ 16 mm and z ¼ 16 mm. For N ¼ 512, the propagation distance zmax ¼ 16 mm corresponds to a shear angle of 45° in a diagram scaled with k ¼ x=dx and n ¼ m=dm. By this choice, the measurements from both operators can be combined to reconstruct the complete AF. In summary, a line uðx0 ; yÞ can be reconstructed from intensity measurements Iz ðx0 ; yÞ. Since the optical setup performs an imaging operation for the x-direction, the different x-positions can, in principle, be reconstructed in parallel. However each line has its own complex constant u ðx0 ; 0Þ (see Eq. (7)).
X. Liu, K.-H. Brenner / Optics Communications 225 (2003) 19–30
To obtain the unknown constants, one measurement of the object, along the x-direction is sufficient.
25
shown in Fig. 10, displaying the modulus in Fig. 10(a) and the phase in Fig. 10(b). 4.1. Numerical demonstration
4. Numerical reconstruction of a 2D complex wave field In order to verify the validity of the new method, we have numerically simulated the reconstruction of a 2D x–y-nonseparable complex distribution of the following form: uðx; yÞ ¼ tðx x0 ; y y0 Þmðx x0 ; y y0 Þ with the complex function pffiffiffiffiffiffiffiffiffiffiffiffiffiffi xy 2p ip 2 2 tðx; yÞ ¼ cos x þy þ 2 e Pp Pm
ð23Þ
ð24Þ
having a rotationally symmetric amplitude and a nonseparable phase distribution. The mask function is defined as
if 0 6 u < 3p=2; rect Dx rect Dy mðx; yÞ ¼ circ Dr if 3p=2 6 u < 2p; ð25Þ assuming r and u as the usual polar coordinates of x and y. Additionally, this function is shifted to the position ðx0 ; y0 Þ. The complete test function is
We have propagated one-dimensionally the test function for distances between z ¼ N dx2 =k to z ¼ N dx2 =k. N is the number of pixels, being N ¼ 256 in this example, dx is the sampling distance and k is the wavelength. For the 1D intensity in one fix column, say j, a 1D FFT was applied and the resulting transform was inserted into a 2D ambiguity function. This process was repeated for all the z-positions in the specified interval, separated by dz ¼ 2dx2 =k. Thus one half of the ambiguity function was generated. The other half was obtained by a combination of 1D Fourier transform and 1D propagation, using the same z-range. From the complete ambiguity function, one column j of the complex amplitude was reconstructed up to an unknown complex constant. By repeating this process for all columns j ¼ 0; . . . ; N 1, the complete 2D complex wave field is reconstructed up to an unknown complex constant, which is different in each column. Fig. 11 shows the modulus and the phase of the uncorrected complex amplitude ucu ðx; yÞ, which was reconstructed with the steps described
Fig. 10. (a) Modulus of u(x,y). (b) Phase of u(x,y).
26
X. Liu, K.-H. Brenner / Optics Communications 225 (2003) 19–30
previously. The suffices u and c stand for uncorrected column reconstruction. To obtain the unknown complex constants for each column, we reconstruct a row of uðx; yÞ at y ¼ 0 using the same procedure as described above, but now working on lines. From this we obtain ulu ðx; 0Þ which can be used to correct ucu ðx; yÞ. Finally we obtain u^ðx; yÞ, which should be equal to uðx; yÞ up to an overall unknown
complex constant u ð0; 0Þ. The modulus and phase of u^ðx; yÞ is shown in Fig. 12. Note that the deviations of the phase only occur, where the modulus is zero. In order to compare the reconstructed result with the original function in Fig. 10, we have multiplied our original function by u ð0; 0Þ and we display the reconstructed and original function at y ¼ 0 in Figs. 13 and 14.
Fig. 11. (a) Reconstructed modulus of uwo ðx; yÞ, without correction in x-direction. (b) Reconstructed phase of uwo ðx; yÞ, without correction in x-direction.
Fig. 12. (a) Reconstructed modulus of u^ðx; yÞ, after the correction in x-direction. (b) Reconstructed phase of u^ðx; yÞ, after correction in x-direction.
X. Liu, K.-H. Brenner / Optics Communications 225 (2003) 19–30
4.2. A simple complexity analysis In order to estimated the amount of data and operations required for the reconstruction of an arbitrary 2D complex wave field, we perform a simple complexity analysis. The reconstruction, described before using this method requires 2N measurements of a 2D intensity distribution to recover N ambiguity functions, from which in turn the uncorrected 2D wave field can be determined. For the correction we require another 2N intensity measurements. Therefore the new method requires a total of 4N measurements of 2D intensities. The
27
number of intensity frames to be recorded thus grows linear with the number of pixels N used for the reconstruction. N is also the 1D space-bandwidth product of the wave field: Dx ¼ DxDm ¼ SBP: ð26Þ N¼ dx The reconstruction of the uncorrected wave field requires 3N 2 1D FFTs of size N. For the correction, another 3N 1D FFTs of size N are required. The total number of 1D FFTs of size N is thus nFFT1 ¼ 3N ðN þ 1Þ
ð27Þ
and the number of operations thus grows with
Fig. 13. (a) Modulus of uðx; 0Þ. (b) Phase of uðx; 0Þ.
Fig. 14. (a) Modulus of the reconstructed field u^ðx; 0Þ. (b) Phase of the reconstructed field u^ðx; 0Þ.
28
nops ¼ 3N 2 ðN þ 1Þ ldN :
X. Liu, K.-H. Brenner / Optics Communications 225 (2003) 19–30
ð28Þ
A comparison of the reconstruction complexity with that of previous methods is difficult to estimate since the exact details of the procedure for reconstructing the 4D ambiguity function are not known. If a straight forward method is used, the construction of a 4D AF requires a minimum of N 2 measurements of N 2 pixel intensities. Using this method, the complete reconstruction with N ¼ 256, described above, required a calculation time of 140 s on a 1.8 GHz CPU. This includes the calculation of the propagated field intensities. For N ¼ 512 the total calculation time was 1240 s, confirming Eq. (28) quite well. 5. Experimental results In order to assess the performance of the new method, we first have verified the operation of the 1D Fresnel propagator and then we have selected as 2D input distribution first a phase-element and then an amplitude-element. 5.1. Experimental verification of the 1D propagation operator In Fig. 15, we compared the result of a numerical 1D propagation of a 2D wave field with the result of our experimental 1D propagation
operator setup. We have used a circular aperture with a diameter of 0.6 mm, which was illuminated by a plane wave. The comparison of the intensities shows clearly, that the setup performs the desired operation. 5.2. Experimental verification of the reconstruction method Since the movement of the combined mount, consisting of the cylindrical lens and the camera requires a high centering accuracy, which was difficult to obtain with standardized components, we have reconstructed only one column of the 2D input signal using our 1D operator. For the phase-element, a microlens array with a pitch of 400 lm and a focal length of 30 mm has been tested. The visible region included approximately three lenses. Fig. 16 shows the original AF using a model microlens array with the same parameters as in the experiment. Fig. 17 shows the modulus of the AF reconstructed from experimental measurements. The horizontal direction is the l-coordinate and the vertical direction is the ycoordinate. The comparison between the two figures shows, that there is a good agreement, except for the center spot. The center peak Að0; 0Þ corresponds to the total power of the input signal. We have erased this peak in the reconstructed figure, since any bias in the measurement affects only this
Fig. 15. 1D light propagation behind a circular aperture. Comparison of the intensities for the numerical (a) and optical (b) 1D propagation operator.
X. Liu, K.-H. Brenner / Optics Communications 225 (2003) 19–30
29
Fig. 18. Numerically computed phase of the lens array. Fig. 16. Numerically computed modulus of the AF.
Fig. 19. Experimentally reconstructed phase of the lens array.
Fig. 17. Experimentally reconstructed modulus of the AF.
pixel. Figs. 18 and 19 show a comparison of original and reconstructed phase curve. A practical problem in the measurement is the shift of the center during longitudinal movement. The shift occurs both in x- and y-direction. Here we only have corrected the shift in y-direction, which partly explains the deviations between reconstructed phase and ideal phase. Furthermore, there are imperfections in the real optical system, like dust and unwanted reflections, which can explain the statistic deviation.
As the amplitude-element, a rectangular aperture has been tested. Figs. 20 and 21 show a comparison of original and reconstructed AF of this rectangle. In Fig. 21 we can see another practical problem, i.e., the determination of the amplitude scale factor for the combination of data from the 1D propagation operator and the 1D Fourier operator. Because of amplitude fluctuation of the laser source, this scale factor can vary during time. Nevertheless the comparison of the reconstructed amplitude shows, that there is a good agreement with the analytic data. Fig. 22 shows the shape and diameter of the rectangular aperture and Fig. 23 shows the reconstructed
30
X. Liu, K.-H. Brenner / Optics Communications 225 (2003) 19–30
Fig. 23. Experimentally reconstructed modulus.
amplitude. The reconstructing was first performed for a plane z > 0 which was our reference plane. The propagation to z ¼ 0 was realized by shearing the experimentally reconstructed AF by 3° and then performing the reconstruction of the object amplitude according to Eq. (7). Fig. 20. Numerically computed modulus of the AF.
6. Conclusion We have described a new method to measure 2D complex wave fields. This method uses the 1D properties of the ambiguity function. The reconstruction operations are performed for one dimension while the other dimension is considered as a series of parallel optical channels. The validity of the new approach is verified by a numerical simulation and also by comparison between simulation and experiment. The experimental results are in good agreement with the numerical simulation.
References
Fig. 21. Experimentally reconstructed modulus of the AF.
Fig. 22. Numerically computed modulus.
[1] A.W. Lohmann, J. Opt. Soc. Am. A 10 (1993) 2181. [2] D. Dragoman, Prog. Opt. 37 (1997) 1. [3] M.G. Raymer, M. Beck, D.F. Alister, Phys. Rev. Lett. 72 (1994) 1137. [4] D.F. McAlister, M. Beck, L. Clarke, A. Mayer, M.G. Raymer, Opt. Lett. 20 (1995) 1181. [5] M. Beck, M.G. Raymer, I.A. Walmsley, V. Wong, Opt. Lett. 18 (1993) 2041. [6] K.-H. Brenner, A.W. Lohmann, J. Ojeda-Castaneda, Opt. Commun. 44 (1983) 323. [7] J.-P. Guigay, Opt. Commun. 26 (1978) 136. [8] R.H. Clarke, J. Opt. Soc. A 3 (1986) 1224. [9] J. Ojeda-Castaneda, L.R. Berriel-Valdos, E. Montes, Appl. Opt. 4 (1998) 790. [10] J. Tu, S. Tamura, Phys. Rev. E 55 (1997) 1946. [11] T. Alieva, F. Agullo-Lopez, Opt. Commun. 114 (1995) 161. [12] D. Dragoman, M. Dragoman, K.-H. Brenner, Appl. Opt. 41 (2002) 5521.
Reconstruction of two-dimensional complex amplitudes from intensity measurements X. Liu *, K.-H. Brenner Chair of Optoelectronics, University of Mannheim, Mannheim, Germany Received 23 January 2003; received in revised form 10 July 2003; accepted 20 July 2003
Abstract The amplitude and phase recovery of optical fields with tomographic methods employing Wigner- or ambiguity functions (WD/AF) has been demonstrated extensively for the case of one-dimensional (1D) functions. For 2D light distributions, the associated WD/AF is 4D, posing several problems. In this paper we introduce a new concept, which allows to reconstruct arbitrary 2D distributions using only 1D measurements. We reconstruct one dimension (y) of the complex light source for each position of the x-dimension. To this end, we realized a 1D propagation operator. The corresponding optical setup for measurement is shown and experimental results are presented. Ó 2003 Elsevier B.V. All rights reserved. PACS: 42.30.Rx; 42.30.wb; 42.25.Bs; 42.60.Jf Keywords: Phase retrieval; Image reconstruction; Tomography; Wave propagation; Beam characteristics
1. Introduction The utilization of intensity measurements for amplitude and phase recovery of light distributions is of considerable interest in several areas of optics. Tomographic methods employing Wigneror ambiguity functions (WD/AF) have been demonstrated extensively for the case of one-dimensional (1D) functions. For the rotation of phase space, fractional Fourier transform meth-
*
Corresponding author. Tel.: +496211812693; fax: +496211812695. E-mail address: [email protected] (X. Liu).
ods were first employed in [1]. The ambiguity function (AF) is closely related to the WDF through a 2D Fourier transform. A review about the WDF and its applications in optics and optoelectronics can be found in [2]. Previously, the WDF was used for the experimental determination of amplitude and phase of either 1D optical beams [3,4] or optical pulses [5]. The AF was used in the past as a polar display of the optical transfer function [6,7], for the description of timevarying random diffracted field [8] or as a tool for designing optical systems with high focal depth [9]. Recently, it was shown theoretically that the AF can be a more straightforward and a much faster tool for amplitude and phase retrieval than
0030-4018/$ - see front matter Ó 2003 Elsevier B.V. All rights reserved. doi:10.1016/j.optcom.2003.07.031
20
X. Liu, K.-H. Brenner / Optics Communications 225 (2003) 19–30
the WDF [10], since the retrieval algorithm is based only on Fourier transforms while in the case of the WDF a filtered back-projection Radon transform is needed, which very often is corrupted by singularities and branch points. The reconstruction of complex amplitudes with the ambiguity function was demonstrated for 1D signals in [11], suggesting that an extension to two dimensions is straightforward. Such extensions have been discussed in [3] by numerical simulation and in [11] on a proof of principle basis. Raymer et al. [3] have shown that a suitable combination of propagation distances and lens focal lengths is sufficient to recover the full 4D Wigner function. Likewise, [11] requires a set of quadratic media with different coefficients in xand y-direction. In Section 2 we point out, that a direct extension to two dimensions generally requires two independent degrees of freedom in order to measure a 4D distribution. For the reconstruction method proposed here, a practical and simple optical setup is presented. The number of measurements grows only linearly with the size of the object. The method is to our knowledge, the first practical method to reconstruct 2D wave fields using only 1D measurements. In this paper, we use the propagation properties of the 1D Ambiguity function to reconstruct a 2D complex distribution. To this end, we have experimentally realized a 1D propagation operator. The corresponding optical setup for the measurements is described and first experimental results are presented.
Aðx0 ; m0 ; uÞ ¼
u~ðmÞ ¼ uðxÞ ¼
Z Z
Z
m0 m0 2pimx0 dm; u~ m e u~ m þ 2 2 ð2Þ
uðxÞe2pixm dx
and ð3Þ
u~ðmÞe2pixm dm:
The AF is a phase-space representation of the complex signal uðxÞ. For the tomographic reconstruction of a complex signal uðxÞ from the AF, two main properties of the AF are used: the line distribution Að0; m0 ; uÞ corresponds to the FT of the intensity distribution of u. Furthermore, light propagation in the AF representation corresponds to a shear of the AF. Thus the FT of the intensity 2 juz ðxÞj at different locations along the propagation ðzÞ direction, corresponds to a line in the AF domain, Aðkzm0 ; m0 ; uÞ (Fig. 1). Therefore, we can reconstruct the AF of a complex signal uðxÞ from the samples of the intensity distribution taken at a set of different z locations (Fig. 2). Finally from Aðx0 ; m0 ; uÞ we can compute the complex signal uðxÞ. We denote u0 ðxÞ as the complex amplitude distribution at z ¼ 0. The propagation to z > 0 in the Fresnel-approximation can be treated by multiplication with a quadratic phase factor in the Fourier-domain (Eq. (4)) 2
u~z ðmÞ ¼ u~0 ðmÞeipkzm ;
ð4Þ
2. Theory of reconstruction for complex 1D signals The AF Aðx0 ; m0 ; uÞ of a complex 1D function uðxÞ is defined equivalently by Eq. (1) or Eq. (2), with the Fourier transform (FT) of uðxÞ written as u~ðmÞ. Because there are different definitions for the FT, the definition, used in this paper, is given in Eq. (3) 0
0
Aðx ; m ; uÞ ¼
Z
x0 x0 2pixm0 u xþ dx; u x e 2 2 ð1Þ
Fig. 1. (Left) juz ðxÞj2 at different locations along the propagation (z) direction. (Top right) Aðx0 ; m0 ; uÞ. (Bottom right) R 0 Aðkzm0 ; m0 ; uÞ ¼ Iz ðxÞe2pixm dx.
X. Liu, K.-H. Brenner / Optics Communications 225 (2003) 19–30
21
Fig. 2. Reconstruction of Aðx0 ; m0 ; u0 ðxÞÞ from Iz ðxÞ.
k is the wavelength of the given light source. Combining Eqs. (2) and (4), the relation between the AF of uz ðxÞ Aðx0 ; m0 ; uz Þ and the AF of u0 ðxÞ Aðx0 ; m0 ; u0 Þ is obtained Aðx0 ; m0 ; uz Þ ¼ Aðx0 kzm0 ; m0 ; u0 Þ:
ð5Þ
Fig. 3. (a) Aðx0 ; m0 ; uÞ; (b) 1. half of Aðx0 ; m0 ; uÞ; (c) 2. half of Aðx0 ; m0 ; uÞ; (d) 1. half of Aðx0 ; m0 ; u~Þ.
uz ðxÞuz ðxÞ
is related to the AF by The intensity Z 0 uz ðxÞuz ðxÞe2pixm dx ¼ Að0; m0 ; uz Þ ¼ Aðkzm0 ; m0 ; u0 Þ:
ð6Þ 2
This indicates that the FT of the intensity juz ðxÞj is a line along in the AF of the signal u0 ðxÞ. In order to obtain the complete Aðx0 ; m0 ; u0 Þ, we have to 2 implement and measure a set of intensities juz ðxÞj at different positions z. The complex amplitude u0 ðxÞ can finally be computed from Aðx0 ; m0 ; u0 Þ using the definition of the AF: Z 0 0 0 u0 ðx Þu0 ð0Þ ¼ Aðx0 ; m0 ; u0 Þepix m dm0 ; R ð7Þ 0 0 Aðx0 ; m0 ; u0 Þepix m dm0 i arg½u0 ð0Þ 0 e u0 ðx Þ ¼ : ju0 ð0Þj Since only the complex amplitude distribution is of interest, the scaling of the absolute value, ju0 ð0Þj, and the phase shift, ei arg½u0 ð0Þ , are usually irrelevant. It is, however, required that ju0 ð0Þj is nonzero in this case! By applying suitable shifts, this restriction can be avoided. In practical propagation situations, we can only 2 measure juz ðxÞj for locations z ¼ zmax , corresponding to a partial angular section of the AF. For the measurement of the remaining section of Aðx0 ; m0 ; u0 Þ we can use another property of the AF Aðm0 ; x0 ; uÞ ¼ Aðx0 ; m0 ; u~Þ:
ð8Þ
In other words, the FT results in a 90° rotation 2 for the AF, allowing to measure j~ u0 ðxÞj instead of 2 ju1 ðxÞj , as illustrated in Fig. 3.
2.1. Extension to the reconstruction of 2D signals For the reconstruction of a 2D complex signal, there are two alternative routes. One alternative is the direct extension of the 1D theory to a 2D theory, i.e., we use the definition and properties of the AF for a 2D distribution. The another alternative is to find a suitable experimental setup, by which the 1D theory of reconstruction can be applied. First we explain the extension of the 1D theory to a 2D theory and we will analyze its requirements. The definition of the AF for a 2D distribution shows that the AF is a 4D function Aðx0 ; y 0 ; m0 ; l0 ; uðx; yÞÞ Z Z x0 y0 x0 y0 ¼ u x þ ;y þ u x ;y 2 2 2 2 0
0
e2piðxm þyl Þ dx dy; 0
0
0
ð9Þ
0
Aðx ; y ; m ; l ; uðx; yÞÞ Z Z m0 l0 m0 l0 ¼ u~ m þ ; l þ u~ m ; l 2 2 2 2 0
0
e2piðx mþy lÞ dm dl; Z Z uðx; yÞ ¼ u~ðm; lÞe2piðxmþylÞ dm dl and Z Z u~ðm; lÞ ¼ uðx; yÞe2piðxmþylÞ dx dy:
ð10Þ
ð11Þ
The following properties of the 4D AF are a direct counterpart to Eqs. (5)–(7).
22
X. Liu, K.-H. Brenner / Optics Communications 225 (2003) 19–30
Propagation corresponds to a shear along the x0 and y 0 axis Aðx0 ; y 0 ; m0 ; l0 ; uz ðx; yÞÞ ¼ Aðx0 kzm0 ; y 0 kzl0 ; m0 ; l0 ; u0 ðx; yÞÞ:
ð12Þ
The 2D intensity is related to the 4D AF by Z Z 0 0 uz ðx; yÞuz ðx; yÞe2piðm xþl yÞ dx dy
3. Experimental setup for 1D light propagation
¼ Að0; 0; m0 ; l0 ; uz ðx; yÞÞ ¼ Aðkzm0 ; kzl0 ; m0 ; l0 ; uz ðx; yÞÞ:
The reconstruction steps have to be repeated for a set of different but fixed x-positions. By selecting a fixed x-position, say x0 , we are dealing with a 1D complex amplitude u0 ðyÞ. This experimental setup for this operation will be described in the next section in detail.
ð13Þ
Reconstruction of uðx; yÞ from Aðx0 ; y 0 ; m0 ; l0 ; uðx; yÞÞ RR 0 0 0 0 Aðx0 ; y 0 ; m0 ; l0 ; uðx; yÞÞepiðx m þy l Þ dm0 dl0 0 0 : uðx ; y Þ ¼ u ð0; 0Þ ð14Þ As shown in Eq. (14), the reconstruction of uðx; yÞ requires the full 4D AF. Since the shear kz in Eq. (13) is the same for the m- and l-direction, a set of different z locations is not sufficient for the reconstruction. The intensity distribution at each z location contains 2D information and the intensity measurement along the z-axis thus provides only 3D information. To reconstruct the AF of general 2D signals, an additional dimension of information is necessary. In order to acquire the full 4D information of the AF, either special symmetries are required [12], or a second dimension such as, e.g., a variable focal length is necessary. While a variation of distance is easy to implement, a setup with lenses of variable focal length is more difficult to realize. Consequently, in [3] the second degree of freedom was achieved by varying the focal length discretely in a few steps, and the distance between lenses was varied in many small steps. In this paper we introduce a new concept, which allows to reconstruct arbitrary 2D distributions using only 1D measurements. Instead of attempting to reconstruct a 4D AF, we consider the 2D wave field as a discrete series of 1D wave fields and we apply the 1D theory of reconstruction from the 2D AF. This method, however, requires to implement a 1D propagation operator, i.e., an operator which performs an exact imaging operation for one direction (x-direction) and a 1D propagation operation for the other direction (y-direction).
In the theoretical description, we require a 1D propagation operator. Such an operator would have to perform an imaging operation for one axis and an independent propagation operation for the other axis. We use the y-axis for the transformation and the x-axis as a series of parallel channels (columns), since this choice corresponds to our optical implementation. Assuming a distribution u0 ðx; yÞ as input and a control parameter z, the 1D propagation operator would have to implement this operation: Z 1 0 2 uz ðx; yÞ ¼ pffiffiffiffiffi e3pi=4 e2piz=k u0 ðx; y 0 Þeip=kzðyy Þ dy 0 : kz ð15Þ 3.1. 1D propagation operator The setup in Fig. 4 implements the desired 1D propagation operator. In the reference plane at z ¼ 0, the input u0 ðx; yÞ is located. The series of positions zP ¼ zn ¼ n dz is our control parameter. In the detector plane, the output is obtained. The goal is that, by the changing zn we can achieve propagation only for the y-direction, while for the x-axis a one-to-one mapping is performed. To analyze the optical setup, we can
Fig. 4. 1D propagation operator for 2D input distribution.
X. Liu, K.-H. Brenner / Optics Communications 225 (2003) 19–30
observe the x- and y-axis (Figs. 5 and 6) separately. This is possible, since the paraxial propagation and the optical components are x–y-separable. For the input, there is no restriction. First consider Fig. 5 for the x-axis. It performs a one-to-one mapping for every point at the input. For this system, the distance between two lenses can be varied without any effect on the intensity distribution in the output plane. Fig. 6 illustrates the operation for the y-axis. At z ¼ 4fy we have a one-to-one mapping, because this system corresponds to the classic imaging system with unit magnification. By adding a separation zP , we can implement the desired 1D propagation operator. The performance of the setup is best described with the formalism of ray-transfer matrices. As has been described in several publications [7], the action of a linear optical system in the ambiguity representation is a simple phase space coordinate transformation. Because all the operations involved (lens, propagation) are separable, we can simplify the description in terms of two 1D transformations. If x and m are the conjugate coordinate pair for the x-direction and y and l the conjugate coordinate pair for the y-direction with
Fig. 5. x view of the 1D propagation operator.
1 kx m ¼ l 2p ky
23
ð16Þ
a 1D propagation is described by the matrix 1 kz ð17Þ 0 1 and transmission through a lens is described by the matrix 1 0 ð18Þ kf1 1 acting on the phase space vector x or ly , rem spectively. Using these matrices, the optical system is described by x-axis: 1 0 Tx ¼ 4fy 4fx þzP 1 ; ð19Þ kfx2
y-axis: 1 Ty ¼ 0
kzP : 1
ð20Þ
The transformation matrix Tx is basically an identity matrix. The off-diagonal element ð4fy 4fx þ zÞ= kfx2 has only an effect on the phase. Since we measure intensity, this component has no effect. The transformation matrix Ty according to Eq. (20) is the 1D propagation matrix with inverted coordinates. Thus the setup in Fig. 4 implements a 1D propagation operator. For practical reasons, the propagation distance is limited to ð2fy fx Þ. Therefore only a finite angular segment can be measured with this system. In order to obtain a larger measuring distance, we use fx ¼ 40 mm and fy ¼ 80 mm . With this choice the measuring range is jzj 6 100 mm. In this experiment we want to reconstruct the discrete AF with dx ¼ 11 lm and dm ¼ 1=ðN dxÞ. In a diagram with k ¼ x=dx and n ¼ m=dm a longitudinal shift zmax ¼ 100 mm corresponds to a shear angle of 45° for 512 512 points. 3.2. 1D FT and propagation operator
Fig. 6. y view of the 1D propagation operator.
The propagation operator described above, can measure the intensity only for distances jzj < zmax , which in the ambiguity domain corresponds to an
24
X. Liu, K.-H. Brenner / Optics Communications 225 (2003) 19–30
angle segment of the full 2D space. In order to be able to measure the remaining angle segment, we developed a 1D FT and propagation operator, i.e., we first perform a 1D Fourier transform and then a 1D propagation operator over a distance zF . In the optical system, we require this propagation distance zF to be adjustable. Assuming a 2D complex amplitude u0 ðx; yÞ as input and the propagation distance zF as a control parameter, we have to implement an optical system which implements a one-to-one mapping for the x-axis and a sequence of Fourier transform and adjustable propagation distance for the y-axis. Fig. 7 shows a system, which realizes this requirement. To discuss its operation, again we can consider the x- and y-axis separately. Actually for the x-axis (Fig. 8) this system is identical to the 1D propagation operator. For the y-axis (Fig. 9) the system is basically a 2f -Fourier transformer, since the lens focal length and the distances fy are equal. At z ¼ 2fy we obtain u~ðx=kfy Þ. The additional separation zF allows us to combine the Fourier transform with a propagation operator. The system can also be described in terms of ray-transfer matrices. x-axis: 1 0 Tx ¼ 2fy 4fx þzF 1 ; ð21Þ kfx2
y-axis: 0 kfy 1 kzF Ty ¼ 1 0 : 0 1 kfy |fflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflffl} |fflfflfflfflfflfflffl{zfflfflfflfflfflfflffl} propagation
ð22Þ
FT
Fig. 7. 1D FT and propagation operator for 2D input distribution.
Fig. 8. x view of the 1D FT and propagation operator.
Fig. 9. y view of the 1D FT and propagation operator.
As designed, the system performs an intensity image for the x-axis and a sequence of Fourier transformation and propagation for the y-axis. For practical reasons, the propagation distance is limited to ðfy fx Þ. In order to match the Fourier-segment with the segment, obtained using the propagation operator, we use fx ¼ 10 mm and fy ¼ 40 mm, dx and dm remain the same value as for the 1D propagation operator. The intensity measuring range is between z ¼ 16 mm and z ¼ 16 mm. For N ¼ 512, the propagation distance zmax ¼ 16 mm corresponds to a shear angle of 45° in a diagram scaled with k ¼ x=dx and n ¼ m=dm. By this choice, the measurements from both operators can be combined to reconstruct the complete AF. In summary, a line uðx0 ; yÞ can be reconstructed from intensity measurements Iz ðx0 ; yÞ. Since the optical setup performs an imaging operation for the x-direction, the different x-positions can, in principle, be reconstructed in parallel. However each line has its own complex constant u ðx0 ; 0Þ (see Eq. (7)).
X. Liu, K.-H. Brenner / Optics Communications 225 (2003) 19–30
To obtain the unknown constants, one measurement of the object, along the x-direction is sufficient.
25
shown in Fig. 10, displaying the modulus in Fig. 10(a) and the phase in Fig. 10(b). 4.1. Numerical demonstration
4. Numerical reconstruction of a 2D complex wave field In order to verify the validity of the new method, we have numerically simulated the reconstruction of a 2D x–y-nonseparable complex distribution of the following form: uðx; yÞ ¼ tðx x0 ; y y0 Þmðx x0 ; y y0 Þ with the complex function pffiffiffiffiffiffiffiffiffiffiffiffiffiffi xy 2p ip 2 2 tðx; yÞ ¼ cos x þy þ 2 e Pp Pm
ð23Þ
ð24Þ
having a rotationally symmetric amplitude and a nonseparable phase distribution. The mask function is defined as
if 0 6 u < 3p=2; rect Dx rect Dy mðx; yÞ ¼ circ Dr if 3p=2 6 u < 2p; ð25Þ assuming r and u as the usual polar coordinates of x and y. Additionally, this function is shifted to the position ðx0 ; y0 Þ. The complete test function is
We have propagated one-dimensionally the test function for distances between z ¼ N dx2 =k to z ¼ N dx2 =k. N is the number of pixels, being N ¼ 256 in this example, dx is the sampling distance and k is the wavelength. For the 1D intensity in one fix column, say j, a 1D FFT was applied and the resulting transform was inserted into a 2D ambiguity function. This process was repeated for all the z-positions in the specified interval, separated by dz ¼ 2dx2 =k. Thus one half of the ambiguity function was generated. The other half was obtained by a combination of 1D Fourier transform and 1D propagation, using the same z-range. From the complete ambiguity function, one column j of the complex amplitude was reconstructed up to an unknown complex constant. By repeating this process for all columns j ¼ 0; . . . ; N 1, the complete 2D complex wave field is reconstructed up to an unknown complex constant, which is different in each column. Fig. 11 shows the modulus and the phase of the uncorrected complex amplitude ucu ðx; yÞ, which was reconstructed with the steps described
Fig. 10. (a) Modulus of u(x,y). (b) Phase of u(x,y).
26
X. Liu, K.-H. Brenner / Optics Communications 225 (2003) 19–30
previously. The suffices u and c stand for uncorrected column reconstruction. To obtain the unknown complex constants for each column, we reconstruct a row of uðx; yÞ at y ¼ 0 using the same procedure as described above, but now working on lines. From this we obtain ulu ðx; 0Þ which can be used to correct ucu ðx; yÞ. Finally we obtain u^ðx; yÞ, which should be equal to uðx; yÞ up to an overall unknown
complex constant u ð0; 0Þ. The modulus and phase of u^ðx; yÞ is shown in Fig. 12. Note that the deviations of the phase only occur, where the modulus is zero. In order to compare the reconstructed result with the original function in Fig. 10, we have multiplied our original function by u ð0; 0Þ and we display the reconstructed and original function at y ¼ 0 in Figs. 13 and 14.
Fig. 11. (a) Reconstructed modulus of uwo ðx; yÞ, without correction in x-direction. (b) Reconstructed phase of uwo ðx; yÞ, without correction in x-direction.
Fig. 12. (a) Reconstructed modulus of u^ðx; yÞ, after the correction in x-direction. (b) Reconstructed phase of u^ðx; yÞ, after correction in x-direction.
X. Liu, K.-H. Brenner / Optics Communications 225 (2003) 19–30
4.2. A simple complexity analysis In order to estimated the amount of data and operations required for the reconstruction of an arbitrary 2D complex wave field, we perform a simple complexity analysis. The reconstruction, described before using this method requires 2N measurements of a 2D intensity distribution to recover N ambiguity functions, from which in turn the uncorrected 2D wave field can be determined. For the correction we require another 2N intensity measurements. Therefore the new method requires a total of 4N measurements of 2D intensities. The
27
number of intensity frames to be recorded thus grows linear with the number of pixels N used for the reconstruction. N is also the 1D space-bandwidth product of the wave field: Dx ¼ DxDm ¼ SBP: ð26Þ N¼ dx The reconstruction of the uncorrected wave field requires 3N 2 1D FFTs of size N. For the correction, another 3N 1D FFTs of size N are required. The total number of 1D FFTs of size N is thus nFFT1 ¼ 3N ðN þ 1Þ
ð27Þ
and the number of operations thus grows with
Fig. 13. (a) Modulus of uðx; 0Þ. (b) Phase of uðx; 0Þ.
Fig. 14. (a) Modulus of the reconstructed field u^ðx; 0Þ. (b) Phase of the reconstructed field u^ðx; 0Þ.
28
nops ¼ 3N 2 ðN þ 1Þ ldN :
X. Liu, K.-H. Brenner / Optics Communications 225 (2003) 19–30
ð28Þ
A comparison of the reconstruction complexity with that of previous methods is difficult to estimate since the exact details of the procedure for reconstructing the 4D ambiguity function are not known. If a straight forward method is used, the construction of a 4D AF requires a minimum of N 2 measurements of N 2 pixel intensities. Using this method, the complete reconstruction with N ¼ 256, described above, required a calculation time of 140 s on a 1.8 GHz CPU. This includes the calculation of the propagated field intensities. For N ¼ 512 the total calculation time was 1240 s, confirming Eq. (28) quite well. 5. Experimental results In order to assess the performance of the new method, we first have verified the operation of the 1D Fresnel propagator and then we have selected as 2D input distribution first a phase-element and then an amplitude-element. 5.1. Experimental verification of the 1D propagation operator In Fig. 15, we compared the result of a numerical 1D propagation of a 2D wave field with the result of our experimental 1D propagation
operator setup. We have used a circular aperture with a diameter of 0.6 mm, which was illuminated by a plane wave. The comparison of the intensities shows clearly, that the setup performs the desired operation. 5.2. Experimental verification of the reconstruction method Since the movement of the combined mount, consisting of the cylindrical lens and the camera requires a high centering accuracy, which was difficult to obtain with standardized components, we have reconstructed only one column of the 2D input signal using our 1D operator. For the phase-element, a microlens array with a pitch of 400 lm and a focal length of 30 mm has been tested. The visible region included approximately three lenses. Fig. 16 shows the original AF using a model microlens array with the same parameters as in the experiment. Fig. 17 shows the modulus of the AF reconstructed from experimental measurements. The horizontal direction is the l-coordinate and the vertical direction is the ycoordinate. The comparison between the two figures shows, that there is a good agreement, except for the center spot. The center peak Að0; 0Þ corresponds to the total power of the input signal. We have erased this peak in the reconstructed figure, since any bias in the measurement affects only this
Fig. 15. 1D light propagation behind a circular aperture. Comparison of the intensities for the numerical (a) and optical (b) 1D propagation operator.
X. Liu, K.-H. Brenner / Optics Communications 225 (2003) 19–30
29
Fig. 18. Numerically computed phase of the lens array. Fig. 16. Numerically computed modulus of the AF.
Fig. 19. Experimentally reconstructed phase of the lens array.
Fig. 17. Experimentally reconstructed modulus of the AF.
pixel. Figs. 18 and 19 show a comparison of original and reconstructed phase curve. A practical problem in the measurement is the shift of the center during longitudinal movement. The shift occurs both in x- and y-direction. Here we only have corrected the shift in y-direction, which partly explains the deviations between reconstructed phase and ideal phase. Furthermore, there are imperfections in the real optical system, like dust and unwanted reflections, which can explain the statistic deviation.
As the amplitude-element, a rectangular aperture has been tested. Figs. 20 and 21 show a comparison of original and reconstructed AF of this rectangle. In Fig. 21 we can see another practical problem, i.e., the determination of the amplitude scale factor for the combination of data from the 1D propagation operator and the 1D Fourier operator. Because of amplitude fluctuation of the laser source, this scale factor can vary during time. Nevertheless the comparison of the reconstructed amplitude shows, that there is a good agreement with the analytic data. Fig. 22 shows the shape and diameter of the rectangular aperture and Fig. 23 shows the reconstructed
30
X. Liu, K.-H. Brenner / Optics Communications 225 (2003) 19–30
Fig. 23. Experimentally reconstructed modulus.
amplitude. The reconstructing was first performed for a plane z > 0 which was our reference plane. The propagation to z ¼ 0 was realized by shearing the experimentally reconstructed AF by 3° and then performing the reconstruction of the object amplitude according to Eq. (7). Fig. 20. Numerically computed modulus of the AF.
6. Conclusion We have described a new method to measure 2D complex wave fields. This method uses the 1D properties of the ambiguity function. The reconstruction operations are performed for one dimension while the other dimension is considered as a series of parallel optical channels. The validity of the new approach is verified by a numerical simulation and also by comparison between simulation and experiment. The experimental results are in good agreement with the numerical simulation.
References
Fig. 21. Experimentally reconstructed modulus of the AF.
Fig. 22. Numerically computed modulus.
[1] A.W. Lohmann, J. Opt. Soc. Am. A 10 (1993) 2181. [2] D. Dragoman, Prog. Opt. 37 (1997) 1. [3] M.G. Raymer, M. Beck, D.F. Alister, Phys. Rev. Lett. 72 (1994) 1137. [4] D.F. McAlister, M. Beck, L. Clarke, A. Mayer, M.G. Raymer, Opt. Lett. 20 (1995) 1181. [5] M. Beck, M.G. Raymer, I.A. Walmsley, V. Wong, Opt. Lett. 18 (1993) 2041. [6] K.-H. Brenner, A.W. Lohmann, J. Ojeda-Castaneda, Opt. Commun. 44 (1983) 323. [7] J.-P. Guigay, Opt. Commun. 26 (1978) 136. [8] R.H. Clarke, J. Opt. Soc. A 3 (1986) 1224. [9] J. Ojeda-Castaneda, L.R. Berriel-Valdos, E. Montes, Appl. Opt. 4 (1998) 790. [10] J. Tu, S. Tamura, Phys. Rev. E 55 (1997) 1946. [11] T. Alieva, F. Agullo-Lopez, Opt. Commun. 114 (1995) 161. [12] D. Dragoman, M. Dragoman, K.-H. Brenner, Appl. Opt. 41 (2002) 5521.