Digital holographic moiré for the direct and simultaneous estimation of strain and slope fields

Rishikesh Kulkarni; Pramod Rastogi

doi:10.1364/OE.22.023192

1. Introduction

In many applications especially in structural analysis and in non-destructive testing, the estimation of phase derivatives is of crucial importance since these are directly related to the strain and slope fields developed on the object surface due to deformation. While various fringe analysis techniques have been developed for the estimation of single phase derivative components corresponding to slopes [1], the problem of the simultaneous estimation of multiple phase derivatives has remained almost unaddressed so far. To perform such measurements, a new class of technique based on digital holographic moiré has recently emerged [2]. Although this technique allows the simultaneous estimation of multiple phase derivatives from a single moiré fringe pattern, it restricts the estimation of phase derivatives to the second order polynomial approximation in a given segment of the interference field.

In this paper, we propose a new method for the direct and simultaneous estimation of multiple phase derivatives corresponding to strain and slope fields from a single moiré fringe pattern in digital holographic moiré. The interference field in a given row/column contains multiple signal components and is modeled as a spatially-varying autoregressive (SVAR) process. The spatially-varying coefficients of this model are estimated by approximating them as the linear combination of the predefined linearly independent basis functions. The poles of the transfer function corresponding to the SVAR model are derived using the estimated coefficients which provide the accurate estimation of multiple phase derivatives. It is important to note that the unwrapping and phase differentiation operations are not required for the proposed technique as it directly provides the unwrapped estimates of the phase derivatives.

2. Theory

The schematic of the experimental set-up of digital holographic moiré for the multiple phase derivative estimation is shown in Fig. 1. A single reference beam directly illuminating the CCD and two object beams illuminating the object are derived from a single color laser source. The object illumination is symmetrical to the normal to the object surface. According to the analysis performed in [3], with the proposed experimental set-up, the sum and the difference of the phase derivatives provide the measurements of the slope and strain fields, respectively. Digital holographic moiré follows the same scheme as the traditional digital holographic interferometry, that is to say, a hologram of the object in its initial state is recorded followed by the recording of another hologram of the object in its deformed state.

Fig. 1: Schematic of the experimental set-up.

Download Full Size | PDF

The object deformation introduces a change in the phase of the wavefields scattered by the object. This change in phase is proportional to the object deformation. An interference field is generated by the conjugate multiplication of the complex wavefields computed by the numerical reconstruction of the two holograms. Let us represent this interference field as I[n, m],

I [n, m] = e^{j ψ_{1} [n, m]} + e^{j ψ_{2} [n, m]} + ς [n, m],

where,

j = \sqrt{- 1}

, n and m represent the pixels along the rows and columns, respectively; ψ₁[n, m] and ψ₂[n, m] are the interference phases associated with the two object beams; ς[n, m] is the complex additive white Gaussian noise. However, the approach proposed in this paper requires the addition of carrier frequency in one arm of the interferometer in order to achieve the spectral separation of the signal components. We propose to introduce the carrier frequency into the system prior to the recording of the hologram of the object in its deformed state. In order to concretize this, we block one of the two beams falling on the object and record a hologram. A second hologram is recorded after providing a tilt to the beam illuminating the object surface. The complex wavefields scattered by the object are computed by numerical reconstruction of these two holograms. An interference field corresponding to the carrier frequency is generated by the conjugate multiplication of the two wavefields. Note that the hologram recording for the carrier frequency estimation is required to be performed only once for a given experimental set-up. Considering the introduction of the carrier frequency in the second arm of the interferometer, Eq. (1) becomes,

I [n, m] = e^{j ψ_{1} [n, m]} + e^{j ψ_{2}^{'} [n, m]} + ς [n, m],

where, ψ′₂[n, m] = ψ₂[n, m] + ω_nn + ω_mm with (ω_n, ω_m) as the carrier frequencies.

The interference phase derivatives can be computed either with respect to m or n. The analysis is carried out to compute the interference phase derivatives with respect to n by considering the interference field in a given column m at a time. Note that a similar analysis holds with respect to m as well, where the interference field in a given row n is considered at a time. The interference field in a given column m is represented as,

I [n] = e^{j ψ_{1} [n]} + e^{j ψ_{2}^{'} [n]} + ς [n] .

Since the interference field in column m is essentially a 1D signal, the interference phase derivatives with respect to n can be written as,

\begin{array}{l} {\dot{ψ}}_{1} [n] = \frac{\partial ψ_{1} [n]}{\partial n} \\ {\dot{ψ}}_{2}^{'} [n] = \frac{\partial ψ_{2}^{'} [n]}{\partial n} . \end{array}

The formulation of the interference field given in Eq. (3) can be looked upon as a signal consisting of two complex exponential frequency modulated components embedded in noise. Consequently, according to [4], I[n] can be expressed as a spatially-varying autoregressive (SVAR) process of order P as,

I [n] = \sum_{p = 1}^{P} a_{p} [n] I [n - p] + ς [n],

where, a_p[n] are the spatially-varying coefficients of the difference equation. The spatially-varying transfer function corresponding to SVAR model can be expressed as,

H (z; n) = \frac{1}{1 - \sum_{p = 1}^{P} a_{p} [n] z^{- p}} .

The roots of the polynomial in the denominator formed by the spatially-varying coefficients a_p[n] represent the poles of the transfer function at each value of n. Under zero noise condition, the SVAR model with model order P = 2 accurately represents the spatial variation of the interference field. In this case, the arguments of the two spatially-varying poles provide the estimates of the interference phase derivatives ψ̇₁[n] and ψ̇′₂[n], respectively. However, in the presence of noise, the interference field cannot be accurately modeled as SVAR with model order of P = 2. Consequently, the SVAR model of I[n] with overdetermined model order P > 2 has been suggested [4]. Accordingly, there are two pole trajectories associated with the two signal components which are different than the remaining P − 2 pole trajectories associated with the noise. The pole trajectories corresponding to the signal lie along the unit circle, whereas the pole trajectories corresponding to the noise are randomly distributed inside or outside the unit circle. This property is utilized to identify the poles associated with the signal components at each value of n.

In order to obtain the accurate estimation of the poles associated with the signal components, the spatially-varying coefficients a_p[n] have to be estimated accurately. To do so, we utilize the method proposed in [5], wherein the spatially-varying coefficients are approximated as a weighted linear combination of a predefined set of linearly independent basis functions. That is,

a_{p} [n] = \sum_{k = 0}^{K} a_{p k} β_{k} [n],

where, β[n] is a set of linearly independent basis functions with the basis dimension K and β₀[n] = 1; a_pk is a set of constant coefficients of size P(K + 1) which is required to be estimated. Essentially, this approach reduces the difficult problem of estimation of spatially-varying coefficients to the relatively simple task of the estimation of a certain number of constant coefficients. Substituting the expression of a_p[n] from Eq. (6) in Eq. (4) we obtain,

I [n] = \sum_{p = 1}^{P} \sum_{k = 0}^{K} a_{p k} β_{k} [n] I [n - p] + ς [n] .

In fact, the noise contributes to error in the estimation of I[n] which is computed using the spatially-varying coefficients and preceding samples of I. We write the estimate of I[n] as,

\hat{I} [n] = \sum_{p = 1}^{P} \sum_{k = 0}^{K} a_{p k} β_{k} [n] I [n - p] .

Consequently, the error in the estimation of I[n] can be written as,

ς [n] = I [n] - \hat{I} [n] .

The sum of the squares of the error in the estimation is computed as,

{‖ ς ‖}^{2} = {| ς [1] |}^{2} + {| ς [2] |}^{2} + \dots + {| ς [N] |}^{2} .

It can be understood by observing Eqs. (6)–(10) that the minimization of ‖ς‖² is basically dependent on the optimum values of constant coefficients a_pk. To compute these coefficients, let us first define the dot product of two complex valued 1D vectors x and y each of length N as,

(x, y) = \sum_{n = 1}^{N} x [n] y^{*} [n],

where, ‘*’ represents the complex conjugate. The sum of the squares of the error in the estimation can now be expressed as,

{‖ ς ‖}^{2} = (ς, ς) .

Substituting the expression of ς from Eq. (9) in Eq. (11),

\begin{array}{l} {‖ ς ‖}^{2} & = (ς, ς) \\ = {‖ I ‖}^{2} - 2 (I, \hat{I}) + {‖ \hat{I} ‖}^{2} . \end{array}

To compute the expression on the right hand side of Eq. (12), we note the following definitions given in [5],

\begin{matrix} w_{k p} & = & {β_{k} [n] I [n - p]} \\ s_{k l} (p, q) & = & (w_{k p} - w_{l q}) \\ s_{0 l} & = & {(s_{0 l} (0, 1), \dots, s_{0 l} (0, P))}^{⊤}, \end{matrix}

where, the superscript “⊤” represents the transpose operation; k, l = 0,...,K and p, q = 1,..,P; Let S_kl represent the P × P matrix with s_kl (p, q) as its (pq)-th element. A Hermitian matrix S of size P(K + 1) × P(K + 1) is built with S_kl as its (kl)-th sub-block. The constant coefficient matrix A of size P(K + 1) × 1 is defined as A = (a₁₀, a₂₀,..., a_P0,..., a_1K,..., a_PK)^⊤. Further, a column matrix s₀ of size P(K + 1) × 1 is defined with the ordered vectors s_0l as its elements. Using these definitions and Eq. (12) we obtain,

{‖ ς ‖}^{2} = {‖ I ‖}^{2} - 2 A^{⊤} s_{0} + A^{⊤} S A .

From Eq. (13), we observe that the absolute minimum of ‖ς‖² is achieved with the optimum value of A when,

S A = s_{0} .

Subsequently, the vector A is computed as, A = S⁻¹s₀. Substituting the estimated coefficients from A in Eq. (6), we compute the spatially-varying coefficients a_p[n].

The selection of the basis functions and the basis dimension plays an important role in the accuracy of coefficient estimation. We have used polynomial basis and Fourier basis functions in our study. The use of other basis functions such as Chebyshev, Legendre, wavelet, etc. have also been mentioned in the literature. The polynomial and Fourier basis functions are defined as,

Polynomial basis: $β_{k} [n] = {(\frac{n}{N})}^{k}$
Fourier basis: $β_{k} [n] = {\begin{array}{l} cos \frac{k π n}{2 N} & k even \\ sin \frac{(k + 1) π n}{2 N} & k odd \end{array}$

The selection of basis dimension depends upon the extent of spatial variation of the coefficients of the difference equation. To capture the rapid variations of coefficients, large basis dimension, i.e., large number of basis functions are required. On the other hand, smaller basis dimension suffices if the coefficients exhibit slow spatial variation. Once the spatially-varying coefficients are estimated, the poles corresponding to the transfer function of the SVAR model are computed and the arguments of the poles associated with the signal components provide the estimates of the phase derivatives. This procedure is repeated in all the columns to obtain the complete 2D phase derivative estimates.

3. Simulation and experimental results

To evaluate the performance of the proposed method, we performed a numerical analysis using the 2D simulated moiré fringe pattern of size 513 × 513 shown in Fig. 2(a) generated with two simulated interference phases. A carrier frequency was added in one of the signal components. Initially, in the noiseless condition, the SVAR model order P was set to be equal to 2. The proposed method was implemented to estimate the interference phase derivatives. The root mean square errors (RMSEs) in the estimation are computed using the above mentioned basis functions as shown in Table 1. The error values indicate that the SVAR model with P = 2 accurately models the two component complex exponential signal. Next, the same fringe pattern was simulated with a signal to noise ratio of 15 dB. In this case, the RMSEs in the estimation of phase derivatives shown in Table 2 suggest that the SVAR model with P = 2 fails to accurately model the noisy signal and results in erroneous phase derivative estimates. Consequently, the SVAR model with P = 4 and P = 6 was considered. The RMSEs in the estimation of phase derivatives were computed and are given in Table 2. A simple median filtering was applied to the estimated phase derivatives. The trajectories of the estimated poles computed in the case of noiseless and noisy condition in the column m = 255 are shown in Figs. 2(b) and 2(c), respectively. It can be observed that the pole trajectories associated with the signal components lie along the unit circle, whereas the pole trajectories associated with the noise components lie inside and outside of the unit circle.

Fig. 2: (a) Moiré fringe pattern. Pole trajectories in the case of (b) noiseless signal (c) noisy signal for m = 255.

Download Full Size | PDF

Table 1:. Errors in the estimation of phase derivatives (in radians/pixel) from a noiseless signal using different basis functions and P = 2

View Table | View all tables in this article

Table 2:. Errors in the estimation of phase derivatives (in radians/pixel) from a noisy signal using different basis functions and model order

View Table | View all tables in this article

The phase derivative estimates computed using the Fourier basis and basis dimension K = 15 with model order P = 6 are shown in Figs. 3(a) and 3(b). The carrier frequency was estimated with K = 0. The required phase derivative corresponding to second arm of the interferometer is computed as,

{\dot{ψ}}_{2} [n] = \frac{\partial ψ_{2}^{'} [n]}{\partial n} - ω_{n} .

Fig. 3: Estimated (a) phase derivative 1 (b) phase derivative 2. Error in the estimation of (c) phase derivative 1 (d) phase derivative 2.

Download Full Size | PDF

The errors in the phase derivative estimation are shown in Figs. 3(c) and 3(d). The RMSE values in the estimation of phase derivatives given in Table 2 indicate that the proposed method accurately estimates the phase derivatives with the over-determined model order in the presence of noise.

The experimental validation of the proposed method is performed using a digital holographic moiré set-up. A reference beam and two object beams were derived from a laser source of wavelength 532 nm. A circular object of diameter 5 cm clamped on its edges was used. A carrier frequency was added in one of the object beam arms by tilting a mirror along its path. The object was subjected to out-of-plane deformation by a point load applied at its center. The in-plane deformation was generated by rotating the object along an axis perpendicular to its plane. The resulting moiré fringe pattern after the median filtering of the interference field is shown in Fig. 4(a). The linear fringe pattern corresponding to the carrier frequency added to the second object beam is shown in Fig. 4(b). The proposed method was implemented with the SVAR model order P = 6 and polynomial basis with the basis dimension K = 7. The carrier frequency was estimated with K = 0. The estimated phase derivatives are shown in Figs. 5(a) and 5(b). The sum and the difference of the estimated phase derivatives are shown in Figs. 5(c) and 5(d), respectively.

Fig. 4: (a) Moiré fringe pattern. (b) linear fringe pattern corresponding to the carrier frequency.

Download Full Size | PDF

Fig. 5: Estimated (a) phase derivative 1 (b) phase derivative 2 (c) Sum of phase derivatives (d) Difference of phase derivatives.

Download Full Size | PDF

The accuracy of the phase derivative estimation can further be improved by adaptive selection of the model order and the basis dimension. The computation of ‖ς‖² is performed with the different combinations of the selected values of the model order and basis dimension. The pair for which the absolute minimum of ‖ς‖² is achieved, can be considered as the optimum values of model order and basis dimension. The fringe analysis technique based on the filtering of signal components in the Fourier domain as proposed in [6] assumes that the spectrum of the signal components lie either around the zero frequency or the carrier frequency. This assumption does not hold in the practical situations where the in-plane deformation normally introduce a shift in the spectrum of the signal components. This assumption is not made in the proposed method because it does not require the filtering of the signal components in Fourier domain.

4. Conclusion

A new technique based on the spatially-varying autoregressive (SVAR) modeling of the interference field is proposed for the accurate estimation of the in-plane and out-of-plane phase derivatives from a single moiré fringe pattern in a digital holographic moiré set-up. A method for the spatially-varying coefficient estimation is developed, and it is shown that the over-determination of the SVAR model order allows for the accurate modeling of the interference field in the presence of high noise. A criteria for the accurate selection of spatially-varying poles corresponding to the signal components is developed to obtain the accurate phase derivative estimates. The proposed method does not involve the spectral filtering of the signal components and possesses high noise robustness. Moreover, the carrier frequency present in one of the signal components is also computed by the proposed method. The simulation and experimental results show that the proposed method has a strong potential to become an important tool in fringe analysis for the direct estimation of multiple phase derivatives.

References and links

1. G. Rajshekhar and P. Rastogi, “Fringe analysis: premise and perspectives,” Opt. Laser Eng. 50, iii–x (2012). [CrossRef]

2. R. Kulkarni, S. S. Gorthi, and P. Rastogi, “Measurement of in-plane and out-of-plane displacements and strains using digital holographic moiré,” J. Mod. Opt. 61, 755–762 (2014). [CrossRef]

3. G. Rajshekhar, S. S. Gorthi, and P. Rastogi, “Simultaneous measurement of in-plane and out-of-plane displacement derivatives using dual-wavelength digital holographic interferometry,” Appl. Opt. 50, H16–H21 (2011). [CrossRef] [PubMed]

4. A. A. Beex and P. Shan, “Time-varying prony method for instantaneous frequency estimation at low SNR,” in Proceedings of IEEE International Symposium on Circuits and Systems (IEEE, 1999), pp. 5–8.

5. L. A. Liporace, “Linear estimation of nonstationary signals,” J. Acoust. Soc. Am. 58, 1288 (1975). [CrossRef] [PubMed]

6. G. Rajshekhar, S. SivaGorthi, and P. Rastogi, “Simultaneous multidimensional deformation measurements using digital holographic moiré,” Appl. Opt. 50, 4189–4197 (2011). [CrossRef] [PubMed]

Estimation Error	Basis functions
Estimation Error	Polynomial	Fourier
Phase derivative 1	0.0004	0.0004
Phase derivative 2	0.0006	0.0006

Estimation Error	Basis functions
	Polynomial			Fourier

	P = 2	P = 4	P = 6	P = 2	P = 4	P = 6
Phase derivative 1	0.0535	0.0056	0.0026	0.0532	0.0056	0.0026
Phase derivative 2	0.0535	0.0054	0.0031	0.0534	0.0051	0.0031

Estimation Error	Basis functions
Estimation Error	Polynomial	Fourier
Phase derivative 1	0.0004	0.0004
Phase derivative 2	0.0006	0.0006

Estimation Error	Basis functions
	Polynomial			Fourier

	P = 2	P = 4	P = 6	P = 2	P = 4	P = 6
Phase derivative 1	0.0535	0.0056	0.0026	0.0532	0.0056	0.0026
Phase derivative 2	0.0535	0.0054	0.0031	0.0534	0.0051	0.0031

Digital holographic moiré for the direct and simultaneous estimation of strain and slope fields

Abstract

1. Introduction

2. Theory

3. Simulation and experimental results

4. Conclusion

References and links

Cited By

Figures (5)

Tables (2)

Equations (20)

Optics Express