Motion compensation and noise tolerance in phase-shifting digital in-line holography

doi:10.1364/OE.14.004286

« Show journal navigation

Motion compensation and noise tolerance in phase-shifting digital in-line holography

Michael D. Stenner and Mark A. Neifeld »View Author Affiliations

Optics Express, Vol. 14, Issue 10, pp. 4286-4299 (2006)
http://dx.doi.org/10.1364/OE.14.004286

View Full Text Article

Acrobat PDF (1398 KB)

Journal Search
Article Lookup

Article Tools

Citations

Alert me when this article is cited
Export Citation/Save

Article Navigation

▸ Article Outline
– Abstract
1. Introduction
2. Technique
3. System design and performance limitations
4. Conclusions
– References and links

Article Tools
Full Text
Article Info
References
Cited By
Figures
Metrics
Download PDF

Viewing Equations in the
HTML Article

Optimal Viewing Recommendations

Full Text
Article Info
References (6)
Cited By
Figures (11)
Metrics

Abstract

We present a technique for phase-shifting digital in-line holography which compensates for lateral object motion. By collecting two frames of interference between object and reference fields with identical reference phase, one can estimate the lateral motion that occurred between frames using the cross-correlation. We also describe a very general linear framework for phase-shifting holographic reconstruction which minimizes additive white Gaussian noise (AWGN) for an arbitrary set of reference field amplitudes and phases. We analyze the technique’s sensitivity to noise (AWGN, quantization, and shot), errors in the reference fields, errors in motion estimation, resolution, and depth of field. We also present experimental motion-compensated images achieving the expected resolution.

1. Introduction

Holographic imaging offers many advantages over conventional imaging, including the ability to measure the full complex electric field and to reconstruct it in three dimensions. Recently, the development and proliferation of CCD and CMOS devices have driven holographic imaging, like traditional imaging, more and more toward digital sensing. In holographic imaging, where the detected intensity is not a simple isomorphism of the object intensity, some reconstruction operation is required to extract the image from the detected signal. In digital holography, this operation is usually performed computationally [1

T. Kreis, Handbook of Holographic Interferometry (Wiley-VCH, Weinheim, 2005).

One problem that plagues conventional holographic imaging is the need to separate the three scattered fields that occur during reconstruction; the reference field, the reconstructed object field, and the conjugate field. Traditionally, this is accomplished by angularly separating the reference and object fields during the recording stage so that the three fields will be physically separate upon reconstruction. The result of this technique is that it records not only the interference from structure on the object field, but also a spatial carrier, or grating, from the angle between the two fields. If the fields are to be separated upon reconstruction, then the spatial carrier must be of higher spatial frequency than the object field itself contains, and so recording the spatial carrier requires greater sensor resolution than would be required to simply record the object field. The additional resolution is used to store the phase information of the field [2

J. W. Goodman, Introduction to Fourier Optics , 3rd ed. (Roberts & Company, Greenwood Village, 2004).

Alternately, phase-shifting holographic imaging can be used to record the phase information in an in-line geometry by acquiring multiple frames of intensity data [1

T. Kreis, Handbook of Holographic Interferometry (Wiley-VCH, Weinheim, 2005).

, 3

I. Yamaguchi and T. Zhang, “Phase-shifting digital holography,” Opt. Lett. 22, 1268–1270 (1997). [CrossRef] [PubMed]

, 4

S. Lai, B. King, and M. A. Neifeld, “Wave front reconstruction by means of phase-shifting digital in-line holography,” Opt. Comm. 173, 155–160 (2000). [CrossRef]

]. Each frame is acquired using a different phase shift (and possibly amplitude) in the plane wave reference beam. Because there is little or no angle between the object and reference beams, no resolution is sacrificed detecting a spatial carrier.

One obvious problem with phase-shifting holographic imaging is that object motion during the multi-frame acquisition will lead to errors in the reconstruction. If the object motion is limited to fixed-velocity lateral translation only, then the motion can be estimated by calculating the cross-correlation of two frames acquired with identical reference beams. The two fields will be identical except for lateral translation, which can be extracted from the peak of the cross-correlation. This translation can then be used to estimate the translation of each frame, and a motion-compensated reconstruction can be computed.

In this paper, we present in detail the process described above. We begin by describing the physical setup, motion compensation, and reconstruction process. We include several experimental images and also discuss the limitations of this technique, including noise performance and sensitivity to various forms of error.

2. Technique

In this section, we describe the process of acquiring the raw data and reconstructing the field in the object plane from it. We begin by describing the physical apparatus and the process of data acquisition. We then describe a technique for reconstructing the field in the object plane from the raw images, which occurs in three stages; motion compensation, sensor field reconstruction, and Fresnel-Kirchoff back-propagation. Finally, we present several experimental images acquired using this technique.

2.1. Physical apparatus and data acquisition

The holographic imaging technique described herein can be used with a very simple physical apparatus, as shown in Fig. 1. In our setup, the object is illuminated in transmission with coherent light. The diffracted light from the object is then combined with a plane reference field using a beam splitter, and the combined field is observed using a CCD camera. The phase of the reference beam can be controlled with a liquid crystal phase modulator, which is in turn controlled by a computer. The object is placed on a mechanical stage, which allows for controlled object motion.

The holographic data is acquired by recording M+1 raw sensor frames, where M≥3. For the j-th frame, the reference field has complex amplitude R _j . The first and last frames use the same reference field, R ₁=R_M ₊₁.

2.2. Motion compensation

The first stage of image reconstruction from the raw data is motion compensation. The total lateral motion between the first (j=1) and last (j=M+1) frames can be estimated using the cross correlation of those frames. Because those two frames use the same reference field, they differ only due to the motion of the object. That is, I ₁(x,y)=I_M ₊₁(x+Δx,y+Δy) where Δx and Δy are the distances by which the object has moved in the x and y directions, respectively. In this ideal limit, the cross-correlation

X_{1, M + 1} (x', y') = I_{1} ⋆ I_{M + 1} = ℱ^{- 1} [ℱ {(I_{1})}^{*} ℱ (I_{M + 1})]

(1)

will have a peak at x′=Δx,y′=Δy, where ℱ and ℱ^-1 indicate Fourier transform and inverse Fourier transform, respectively.

Fig. 1. Experimental setup for motion-compensated digital in-line holography. The object is transmissive, three-dimensional, and movable with a translation stage. The setup is effectively a Mach-Zehnder interferometer with an attenuator and phase shifter in one arm, and the object in the other arm.

Once Δx and Δy have been estimated, the frames can be un-shifted according to

{\bar{I}}_{j} (x, y) = I_{j} (x + \frac{j - 1}{M} Δ x, y + \frac{j - 1}{M} Δ y),

(2)

which assumes the object is moving laterally with constant velocity, and that the original frames are acquired at a constant rate. Here, the intensities Ī_j (x,y) represent the intensity patterns that would have been recorded if the object had remained at the same position (it’s position for frame j=1) for all frames.

For the transmissive geometry described in Fig. 1, the reconstruction is largely insensitive to axial motion, as described in Subsection 3.5. However, for a reflective geometry, axial motion will introduce additional phase shifts which will result in phase errors in the reconstruction. For translations large compared to λ, such phase shifts can become nearly random, making reconstruction difficult. There are, however, techniques for estimating the phase shifts from the image data itself [1

T. Kreis, Handbook of Holographic Interferometry (Wiley-VCH, Weinheim, 2005).

, 5

K. Larkin, “A self-calibrating phase-shifting algorithm based on the natural demodulation of two-dimensional fringe patterns,” Opt. Express 9, 236–253 (2001), http://www.opticsinfobase.org/abstract.cfm?URI=oe-9-5-236. [CrossRef] [PubMed]

], which may be useful in such cases.

2.3. Reconstruction

Having compensated for the object motion, the next step is to reconstruct the object field in the sensor plane f(x,y) from the intensity measurements Ī_j . Many mathematical formalisms for such reconstruction have been developed [1

T. Kreis, Handbook of Holographic Interferometry (Wiley-VCH, Weinheim, 2005).

], but we present here a very general framework allowing for linear reconstruction of the electric field using M≥3 frames, each with a different reference field. Having served its purpose in the motion compensation, frameM+1 is no longer needed. For each frame j, the recorded intensity can be written as

{\bar{I}}_{j} = {∣ R_{j} ∣}^{2} + {∣ f ∣}^{2} + 2 (R_{j}^{'} f' + R_{j}^{″} f ″),

(3)

where R_j (f) is the complex amplitude of the reference (object) field at the sensor. A single prime indicates the real part of a quantity, and a double prime indicates the imaginary part. Here, the spatial dependence of Ī_j and f_j have been suppressed for clarity. The reference field R_j is a plane wave and so has no spatial dependence. This equation holds for all pixels independently. If the values of R_j are known, then we can define the quantity

D_{j} = {\bar{I}}_{j} - {∣ R_{j} ∣}^{2} = {∣ f ∣}^{2} + 2 (R_{j}^{'} f' + R_{j}^{″} f ″) .

(4)

For M exposures, this can be written in matrix form as

D = (2 R + 1 f^{⊤}) f,

(5)

where ^T indicates transpose and

D = [\begin{matrix} D_{1} \\ D_{2} \\ ⋮ \end{matrix}], R = [\begin{matrix} R_{1}^{'} & R_{1}^{″} \\ R_{2}^{'} & R_{2}^{″} \\ ⋮ & ⋮ \end{matrix}], 1 = [\begin{matrix} 1 \\ 1 \\ ⋮ \end{matrix}], f = [\begin{matrix} f' \\ f ″ \end{matrix}] .

(6)

In general, this system of equations is nonlinear, but it can be linearized if a matrix P can be found such that both

P R_{p} = T,

(7)

where R_{p} = [\begin{matrix} R_{1}^{'} & R_{1}^{″} & 1 \\ R_{2}^{'} & R_{2}^{″} & 1 \\ ⋮ & ⋮ & ⋮ \end{matrix}] and T = [\begin{matrix} 1 & 0 & 0 \\ 0 & 1 & 0 \end{matrix}] .

(8)

Solutions for P exist as long an M≥3 and R is rank-2. As discussed in Subsection 3.2, a good choice is P=

{TR}_{p}^{+}

, where ⁺ indicates the Moore-Penrose pseudo-inverse. Once P is found, the field at the sensor plane can be calculated as f=PD/2.

The equations presented above describe reconstruction of a single pixel. Extending this technique to multiple pixels requires only that D becomes an M×N matrix and f becomes 2×N, where each column represents one of N pixels.

2.4. Back-propagation

With the object field at the sensor plane f(x,y) calculated, estimating the field at the object plane g(x,y) is simply a matter of Fresnel-Kirchoff back-propagation according to

\hat{g} (x, y) = f (x, y) * K (x, y) = ℱ^{- 1} [ℱ (f (x, y)) ℱ (K (x, y))],

(9)

where ∗ denotes convolution and K(x,y) is the Fresnel-Kirchoff diffraction kernel

K (x, y) = \frac{e^{ikz}}{i λ z} e^{\frac{ik}{2 z} (x^{2} + y^{2})} .

(10)

Even in the absence of noise, the field ĝ(x,y) is an estimate of g(x,y) because the finite-aperture diffraction process acts as a filter, as discussed below.

Figure 2 shows several experimental images of a sample object. The images are acquired with M=4 and the reference fields are R,iR,-R,-iR,R. The illumination is provided by a He:Ne laser operating at 632.8 nm and the raw images are captured on a CCD with 765×511 pixels, each 9µm×9µm. Each image is 500×500 pixels and represents a 4.5 mm×4.5 mm area. Figure 2(a) shows a traditional (non-holographic) image of the object acquired using a 75 mm focal-length lens in a one-to-one imaging configuration. For the holographic imaging, the object was moved horizontally by approximately 100µm between each frame. The intensity at the object plane (z=107 mm), reconstructed without motion compensation is shown in Fig. 2(b). Reconstruction from the same data with motion compensation yields Fig. 2(c). For comparison, a single raw captured frame of the holographic data is shown in Fig. 2(d).

Fig. 2. Several images of an USAF-1951 resolution target. (a) shows a traditional image of the resolution target. (b) shows a holographic reconstruction of the intensity at the plane of the horizon ally-moving target without motion compensation. (c) shows the same reconstruction with motion compensation. (d) shows a single raw frame of the holographic data.

By varying z in the back-propagation calculation, it is possible to “focus” on different object depths computationally from a single acquisition, as shown in Fig. 3. Both Figs. 3(a) and 3(b) are reconstructions of the same holographic data. The top and bottom targets were placed at different distances from the sensor.

3. System design and performance limitations

3.1. Optical design

In the simplest configuration, as in Fig. 1, there is no lens between the object and sensor. In this case, the resolution limits of the system are dominated by the sensor pixel size and traditional diffraction limitations. The pixel size presents an object-space resolution limit because the finest spatial frequency which can be sampled on the sensor has period 2p, where p is the pixel separation. In the lensless case, such a fringe pattern on the sensor is generated by a grating of the same period at the object.

The other major limit on optical resolution is the traditional diffraction limit. It may not be obvious at first that such limitations apply for the holographic case, but the limitations can be understood simply by considering a localized grating of period Δ on the object, as illustrated in Fig. 4. Such a grating will scatter incident light at angles θ _m =sin^-1(mλ/Δ). If we insist that at least two scattering orders (including m=0) fall on the sensor regardless of the position of the localized grating, then we require N_p ≥2λd/Δ, where N is the linear number of pixels and d is the distance between the object and sensor. The diffraction-limited object resolution is then

r_{d} = \frac{Δ}{2} \approx λ \frac{d}{N_{p}} .

(11)

This result is extremely similar to Rayleigh’s criterion, which provides r_R ≈1.22λd/(N_p ).

Fig. 3. Two reconstructions at different depths from a single holographic data set. (a) shows a reconstruction at z=107 mm. The top target is “in focus” and the bottom target is not. (b) Using the same data but reconstructing at z=244 mm, the reverse is true.

Fig. 4. An illustration of how diffraction limits resolution for holographic reconstruction. So that the localized grating can be reconstructed, at least two diffraction orders (including m=0) must be recorded by the sensor. This limits the minimum grating size that can be resolved.

A reasonable system design approach is to equate these two resolution limits, suggesting d _opt≈

N_{p}^{2}

/λ. For example, given the camera used to acquire the data presented in Figs. 2 and 3, the ideal working distance is d _opt=511×(9×10^-6)²/(632.8×10^-9)=0.065 m, yielding an achievable diffraction-limited resolution of r_d =9µm, the pixel size. A larger distance (0.107 m) was used above due to mechanical constraints, resulting in a predicted resolution of r_d =15µm. The observed resolution is approximately 17µm in the vertical direction and 25µm in the horizontal direction, suggesting that the horizontal motion led to a reduction in resolution by approximately one pixel, as described in Subsection 3.5.

Fig. 5. Portion of experimental setup when using a relay lens. The holographic sensor observes the image, allowing control of the object resolution by adjusting the magnification.

Using the design approach described above, the object resolution is completely controlled by the sensor pixel size. In order to change the object resolution, one can use relay optics to magnify the object as shown in Fig. 5. The holographic sensor then observes the optical image rather than the object itself. Because Fresnel-Kirchoff propagation can go either backward or forward, the image may fall beyond the sensor. Also, the image can be either real or virtual. In these terms, the image magnification is M_i =s_i /s_o and the resolution becomes

r_{d} \approx \frac{1}{M_{i}} λ \frac{d}{N_{p}} .

(12)

Of course, the lens must support this resolution. Specifically, its diameter D_l must satisfy D_l ≥λs_i /p. Also, the plane reference field must either be mixed with the scattered object light after the lens (as shown in Fig. 5) or the reference field must be prepared such that it has planar wavefronts as it strikes the sensor.

3.2. Noise and dynamic range

In any imaging system, noise and dynamic range are concerns. Due to the post-processing of the raw image data, the noise properties of this holographic system are not obvious. Here, we describe the effects of additive white Gaussian noise (AWGN) from the sensor, quantization noise, and shot noise.

For AWGN with the same distribution in each raw frame, it can be shown that the noise variance on the reconstructed field f̂ at the sensor is given by

σ_{fg}^{2} = σ_{fg}^{' 2} + σ_{fg}^{″ 2}

(13)

= \frac{1}{4} \sum_{j = 1}^{M} [{(P_{j}^{'} σ_{Ig})}^{2} + {(P_{j}^{″} σ_{Ig})}^{2}]

(14)

= \frac{1}{4} σ_{Ig}^{2} {∥ P ∥}_{F}^{2},

(15)

where

σ_{Ig}^{2}

is the noise variance of measurement I_j and ‖P‖ _F is the Frobenius norm of P. The real and imaginary noise will in general be correlated, but Eq. (13) holds because the noise adds as n_f =n′_f +in″_f and so |n_f |²=

{n'}_{f}^{2}

{n ″}_{f}^{2}

, where n_f is the actual noise on the field, and n′_f and n″_f are the real and imaginary parts, respectively. From Eq. (15), it is clear that to minimize the noise in f̂, one must minimize ‖P‖ _F . Note that the third column in R _p only serves to force the row-sums of P to be zero, and the actual value in that column is unimportant. In the limit that it goes to infinity, the values in the third row of

R_{p}^{+}

approach zero and ‖

R_{p}^{+}

‖ _F =‖P‖ _F . Therefore, because the Moore-Penrose pseudo-inverse is minimum-norm, choosing P=

{TR}_{p}^{+}

provides the minimum-norm solution of PR _p =T and the minimum-AWGNlinear solution for the field f.

The effects of quantization error are very similar to those of AWGN. For typical signals, the quantization noise takes the form of additive uniformly distributed noise. The quantization noise variance at the sensor can be written as

σ_{Iq}^{2} = {(\frac{I_{max}}{(2 \sqrt{3}) 2^{b}})}^{2},

(16)

where I _max is the maximum measurable intensity, b is the number of bits, and the quantity 2√3 comes from the assumption of a uniform intensity distribution. This assumption need only be valid on the scale of a quantization level I _max/2 ^b . Calculating the quantization noise variance on f̂ is very similar to the AWGN case, leading to

σ_{fq}^{2} = \frac{1}{4} σ_{Iq}^{2} {∥ P ∥}_{F}^{2} .

(17)

The shot noise calculation is less simple than either of the cases presented above, largely because shot noise is not additive. The shot noise variance observed on a raw pixel is proportional to the intensity of the light falling on that pixel, with a proportionality constant α which depends on the wavelength and sensor properties. The resulting noise variance is given by

σ_{fs}^{2} = \frac{α}{4} \sum_{j = 1}^{M} ({P'}_{j}^{2} + {P ″}_{j}^{2}) [{∣ R_{j} ∣}^{2} + {∣ f ∣}^{2} + 2 (R_{j}^{'} f' + R_{j}^{″} f ″)],

(18)

which simplifies considerably in many common cases. Consider the specific case of

R_{j} = R e^{i (j - 1) \frac{2 π}{M}},

(19)

in which the combined AWGN, quantization, and shot noise contributions take the form

σ_{f}^{2} = σ_{fg}^{2} + σ_{fq}^{2} + σ_{fs}^{2} = \frac{σ_{Ig}^{2}}{{MR}^{2}} + \frac{σ_{Iq}^{2}}{{MR}^{2}} + \frac{α}{M} (1 + \frac{{∣ f ∣}^{2}}{R^{2}}) .

(20)

The corresponding signal-to-noise ratio (SNR) is given by

{SNR}_{f} = \frac{∣ f ∣}{σ_{f}} = {[\frac{{MR}^{2} {∣ f ∣}^{2}}{α (R^{2} + {∣ f ∣}^{2}) + σ_{fg}^{2} + σ_{fq}^{2}}]}^{1 ⁄ 2} .

(21)

Equations (20) and (21) may be slightly misleading because they assume that the energy devoted to each frame remains constant as the number of frames increases. If the energy from the object E_f =M|f|² and from the reference E_R =MR ² is held fixed, then the noise variance becomes

σ_{E}^{2} = \frac{σ_{Ig}^{2}}{E_{R}} + \frac{σ_{Iq}^{2}}{E_{R}} + \frac{α}{M} (1 + \frac{E_{f}}{E_{R}})

(22)

and the SNR becomes

{SNR}_{E} = \frac{∣ f ∣}{σ_{f}} = {[\frac{E_{f} E_{R}}{α (E_{f} + E_{R}) + M σ_{fg}^{2} + M σ_{fq}^{2}}]}^{1 ⁄ 2} .

(23)

Fig. 6. Simulated SNR _E as a function of sensor SNR for linear (solid) and nonlinear (dotted) reconstruction. The upper (high SNR _E ) lines are for M=3 and the lower (low SNR _E ) lines are for M=6. The region to the left is dominated by AWGN and the region to the right by shot noise.

As expected, the SNR _E contribution from shot noise does not depend on M, whereas the contribution from detector AWGN and quantization noise decreases with increasing M. Therefore, the best noise performance is achieved with small M.

Figure 6 shows a plot of simulated SNR _E vs. sensor SNR=I _max/σ _Ig for two values of M using both the linear inversion technique described above and also a nonlinear technique. The nonlinear approach chooses f to minimize the magnitude of (2R+1f ^T)f-D at each pixel. The simulation includes both AWGN and shot noise. In the high-sensor-noise region on the left, the performance is dominated by sensor noise. In that limit, the linear dependence of SNR _E on the sensor SNR dominates and the 1/√M dependence is clearly visible. Also, we see that there is a slight improvement (less than 1 dB) in SNR _E achieved by using nonlinear reconstruction. At high sensor SNR, performance is dominated by shot noise. In that regime, nonlinear reconstruction continues to provide improved SNR _E , but the 1/√M dependence is lost.

Equations (21) and (23) clearly show that it is always advantageous to use large field strengths for both the reference and objects fields. As illustrated in Fig. 7, best performance is achieved by maximizing the dynamic range, given by 4f _max R _max, where f _max and R _max are the maximum amplitudes of the object and reference fields, respectively. However, the sensor places a practical upper limit on the field strengths because the total intensity may not exceed I _max. Therefore, one must impose the constraint

f_{max} + R_{max} \leq \sqrt{I_{max}}

. As an example, given an object field with fixed strength

f_{max} = \sqrt{I_{max}} ⁄ 4

, one should choose

R_{max} = \sqrt{I_{max}} - f_{max} = 3 \sqrt{I_{max}} ⁄ 4

yielding a dynamic range of 3I _max/4. If both object and reference field strengths can be adjusted, then one should choose

f_{max} = R_{max} = \sqrt{I_{max}} ⁄ 2

All of the preceding noise analysis is for the field on a single pixel of the sensor. The final field calculation also includes the Fresnel-Kirchoff back-propagation. The primary effect of this back-propagation is to “color” the AWGN according to the magnitude of the amplitude transfer function (ATF) |H(k_x ,k_y )|=|ℱ[K(x,y)]|. This filtering process can also be thought of as introducing spatial correlations. Fig. 8 shows plots of the ATF as a function of k_x for a 500×500 sensor with 9µm pixels for several reconstruction distances z. In all cases, the high frequencies are attenuated, but for larger z the coloring effect is more extreme. The vertical lines mark the spatial frequencies k_x =1/(2r_d ), where r_d is the resolution derived in Eq. (11). Figure 9 shows several simulated images corresponding to the same conditions. For larger z, the noise is less grainy, suggesting weaker high-frequency components. Obviously, the object field also experiences the same low-pass filtering, leading to worse resolution and a blurry image.

Fig. 7. An illustration of the intensities that may fall on the sensor given an object field of maximum amplitude f _max and maximum reference field amplitude R _max. The central region between the dashed lines of height 4f _max R _max represents the part of the sensor’s full dynamic range that is actually used.

Fig. 8. Amplitude transfer function (ATF) as a function of spatial frequency k_x for several values of z. In each case, the image is 500×500 9µm pixels. The vertical lines show the spatial frequency k_x =1/(2r_d ) corresponding to the predicted resolution at that value of z.

The shot noise will likely already have spatial structure because it depends on the object field intensity |f(x,y)|²). Therefore, calculating its final spatial correlations is difficult, although like the AWGN, it will be (additionally) colored by the ATF.

Fig. 9. Several simulated images corresponding to the ATF plots in Fig. 8. The top row of images shows the entire 500×500 for each value of z, whereas the lower images are 120×120 pixels, showing the central region. Both the noise and the object field are smoothed at increasing z.

3.3. Errors in Rj

The previous section discussed the result of errors in the measurements D due to noise. It is also likely that the actual values of the reference fields will differ from the desired values. The impact of such error is easily seen by considering an estimated reference matrix R̂, the corresponding inversion matrix P̂, and the resulting estimated field f̂. Equation (5) can be solved using the same measurements D with either the ideal or estimated versions of R, P, and f. By equating both versions, one finds

\hat{f} = [I + (\hat{P} - P) R] f

(24)

= (I + E) f,

(25)

where E is a 2×2 error matrix. This matrix is constant for all pixels. In general, it is difficult to predict the value of this matrix given arbitrary R and arbitrary error. However, we present a few common cases here.

If the measured reference fields relate to the true fields as R̂=(1+ε)R, then E=-ε I/(1+ε) and the resulting field is also scaled in amplitude. If all reference fields have a constant phase error R̂_j =e^iϕR_j , then the resulting field measurement is f̂=e- ^iϕf.

If the reference fields take the form of Eq. (19) and R̂_j =(1+ε_j )R_j with small Gaussian-distributed ε_j , then the RMS of the error matrix is given approximately by

rms [E] = [\begin{matrix} 1 & 0 \\ 0 & 1 \end{matrix}] \sqrt{\frac{2}{M}} σ_{ε} (real ε_{j}),

(26)

= [\begin{matrix} 0 & 1 \\ 1 & 0 \end{matrix}] \sqrt{\frac{2}{M}} σ_{ε} (imaginary ε_{j}),

(27)

= [\begin{matrix} 1 & 1 \\ 1 & 1 \end{matrix}] \sqrt{\frac{1}{M}} σ_{ε} (complex ε_{j}),

(28)

where

σ_{ε}^{2}

is the variance of ε_j . We see that the variance of the error in f depends only on the variance of the error in R_j . Also, amplitude error in R_j tends to create amplitude error in f, and phase error creates phase error. However, these hold only in the limit that σ_ε≪1; for larger σ_ε, phase and amplitude errors begin to mix.

Fig. 10. Illustration of how the resolving power with a defocus of Δz is the same for holographic reconstruction and traditional imaging with a lens of focal length F=z/2 and diameter D=N_p .

3.4. Depth of field and reconstruction

Because this technique measures the electric field f (x,y) at the sensor and then reconstructs the field g(x,y) at any position z, it is useful to define two quantities: the depth of field (DOF) and depth of reconstruction (DOR). The DOF has a similar meaning to its use in traditional imaging; it is the range over which an object remains (approximately) in focus for a single reconstruction. In contrast, the DOR is the range of z over which a given f(x,y) can be used to reconstruct g(x,y) with good resolution.

The DOR has largely been calculated in Subsection 3.1, where we found that the diffraction-limited resolution varies approximately as λz/(N_p ) for the lensless case. Therefore, a factor of two loss in resolution occurs at z=2

N_{p}^{2}

/λ. There is no loss of resolution at z<

N_{p}^{2}

/λ; the resolution is limited by the pixel size in that domain. However, one must use an appropriate diffraction calculation [2

J. W. Goodman, Introduction to Fourier Optics , 3rd ed. (Roberts & Company, Greenwood Village, 2004).

] and it may be necessary to upsample the image to avoid aliasing in the diffraction kernel.

The depth of field is less obvious. The resolving power of this system for an object placed at z+Δz, where z is the reconstruction depth, is identical to a traditional lens with focal length F=z/2 and diameter D=N_p , as illustrated in Fig. 10. This is because the lens and propagation through distance z ₀ introduce the same quadratic phase term as back-propagation by a distance -z0 [2

J. W. Goodman, Introduction to Fourier Optics , 3rd ed. (Roberts & Company, Greenwood Village, 2004).

]. Therefore, the DOF is roughly given by twice the Rayleigh range of the equivalent lens system,

DOF = 2 z_{R} = \frac{π}{2} {(\frac{r_{d}}{λ})}^{2} λ = \frac{π}{2} λ {(\frac{d}{N_{P}})}^{2} .

(29)

Over this range, the resolution remains less than approximately √2rd. For the experimental results presented in Fig. 2, we find DOF=0.56 mm.

3.5. Motion-related error

While the technique presented in Subsection 2.2 compensates for constant-velocity linear lateral motion, deviation from that assumption has several effects on the resulting image. The most obvious is the traditional blur that results from an object moving a distance t _ev_o long compared to the pixel size p, where t_e is the exposure time and v_o is the object velocity. However, this form of blur affects each frame identically and results in final-image blur that is very similar to the blur in a conventional image. However, in a traditional image, the intensity I(x,y) is blurred, whereas in this case, the field f̂(x,y) (and ĝ(x,y)) is blurred. That is, a traditional images is blurred according to Ī(x,y)=I(x,y)∗B(x,y), where B(x,y) is a blur function. Because the blur is applied to each I_j (x,y) identically, the result is

\hat{f} (x, y) = f (x, y) * (x, y) and

(30)

\hat{g} (x, y) = [f (x, y) * B (x, y)] * K (x, y)

(31)

= B (x, y) * [f (x, y) * K (x, y)] .

(32)

Fig. 11. MSE of a sample image as a function of defocus (solid) and axial shifts (dashed). The shifts were performed for M=3,6,9, with increasing error for larger M because the overall shift is larger. For defocus, the object was placed at z=z ₀+Δz. For the shifted images, each frame was acquired with the object at z_j =z ₀+[j-(M+1)/2]Δ _z . In all cases, z ₀=65 mm and the image size 500×500 9µm pixels.

More pervasive are errors in the velocity estimation or deviation from constant-velocity linear motion. Either of these will result in motion-compensation errors such that the individual frames are not properly re-aligned, leading to blur-like errors on the length scale of the frame misalignment. Frame misalignment on the order of the pixel size are inevitable in the simplest approach, where one simply uses the location of the largest value of the cross-correlation to calculate the shift, and then shifts each frame by the nearest integer number of pixels. One can, estimate sub-pixel shifts to reduce this error further [6

R. C. Hardie, K. J. Barnard, J. G. Bognar, E. E. Armstrong, and E. A. Watson, “High-resolution image reconstruction from a sequence of rotated and translated frames and its application to an infrared imaging system,” Opt. Eng. 37, 247–260 (1998). [CrossRef]

Another important consideration for velocity estimation is that the first and last frames should be, as much as possible, shifted versions of each other. Specifically, a fixed background in addition to the moving object, or large parts of the object falling out of the field of view (or entering it) between frames will confuse the cross-correlation and lead to errors in the velocity estimation. The algorithm is quite robust, but when it fails, it often fails catastrophically, leaping to a peak far from the correct one. For example, a fixed background and moving object will result in two strong peaks which may be far apart, and the algorithm will tend to “hop” between them.

As mentioned above, the process described here does not estimate or compensate for axial motion. However, axial translation introduces errors similar to that produced by simple defocus. Figure 11 shows the simulated mean square error (MSE) of a sample image as a function of defocus (solid) and axial shifts (dashed). Several values of M are shown and because total shift increases with M, the MSE is also larger for larger M. This figure demonstrates that in order to avoid losing resolution as a result of axial motion, one must keep the axial translation small compared to the DOF.

4. Conclusions

We have presented a technique for motion-compensated phase-shifting digital in-line holography based on velocity estimation with an additional frame and the cross-correlation. We have also described a new linear reconstruction approach which minimizes AWGN for an arbitrary set of reference fields.

Acknowledgments

The authors would like to thank George Barbastathis for many helpful discussions.

References and links

1.	T. Kreis, Handbook of Holographic Interferometry (Wiley-VCH, Weinheim, 2005).
2.	J. W. Goodman, Introduction to Fourier Optics , 3rd ed. (Roberts & Company, Greenwood Village, 2004).
3.	I. Yamaguchi and T. Zhang, “Phase-shifting digital holography,” Opt. Lett. 22, 1268–1270 (1997). [CrossRef] [PubMed]
4.	S. Lai, B. King, and M. A. Neifeld, “Wave front reconstruction by means of phase-shifting digital in-line holography,” Opt. Comm. 173, 155–160 (2000). [CrossRef]
5.	K. Larkin, “A self-calibrating phase-shifting algorithm based on the natural demodulation of two-dimensional fringe patterns,” Opt. Express 9, 236–253 (2001), http://www.opticsinfobase.org/abstract.cfm?URI=oe-9-5-236. [CrossRef] [PubMed]
6.	R. C. Hardie, K. J. Barnard, J. G. Bognar, E. E. Armstrong, and E. A. Watson, “High-resolution image reconstruction from a sequence of rotated and translated frames and its application to an infrared imaging system,” Opt. Eng. 37, 247–260 (1998). [CrossRef]

OCIS Codes
(090.2880) Holography : Holographic interferometry
(110.6880) Imaging systems : Three-dimensional image acquisition

ToC Category:
Holography

History
Original Manuscript: February 27, 2006
Revised Manuscript: April 27, 2006
Manuscript Accepted: April 27, 2006
Published: May 15, 2006

Citation
Michael D. Stenner and Mark A. Neifeld, "Motion compensation and noise tolerance in phase-shifting digital in-line holography," Opt. Express 14, 4286-4299 (2006)
http://www.opticsinfobase.org/oe/abstract.cfm?URI=oe-14-10-4286

Sort: Author | Year | Journal | Reset

References

T. Kreis, Handbook of Holographic Interferometry (Wiley-VCH, Weinheim, 2005).
J. W. Goodman, Introduction to Fourier Optics, 3rd ed. (Roberts & Company, Greenwood Village, 2004).
I. Yamaguchi and T. Zhang, "Phase-shifting digital holography," Opt. Lett. 22, 1268-1270 (1997). [CrossRef] [PubMed]
S. Lai, B. King, and M. A. Neifeld, "Wave front reconstruction by means of phase-shifting digital in-line holography," Opt. Commun. 173, 155-160 (2000). [CrossRef]
K. Larkin, "A self-calibrating phase-shifting algorithm based on the natural demodulation of two-dimensional fringe patterns," Opt. Express 9, 236-253 (2001). [CrossRef] [PubMed]
R. C. Hardie, K. J. Barnard, J. G. Bognar, E. E. Armstrong, and E. A. Watson, "High-resolution image reconstruction from a sequence of rotated and translated frames and its application to an infrared imaging system," Opt. Eng. 37, 247-260 (1998). [CrossRef]

Cited By

Alert me when this paper is cited

OSA is able to provide readers links to articles that cite this paper by participating in CrossRef's Cited-By Linking service. CrossRef includes content from more than 3000 publishers and societies. In addition to listing OSA journal articles that cite this paper, citing articles from other participating publishers will also be listed.

Click here to see a list of articles that cite this paper