Accurate and fast two-step phase shifting algorithm based on principle component analysis and Lissajous ellipse fitting with random phase shift and no pre-filtering

Yu Zhang; Xiaobo Tian; Rongguang Liang

doi:10.1364/OE.27.020047

1. Introduction

With the development of advanced optical manufacturing technology, the advanced optical testing technology has also been developed, and the interferometry is an easy and effective optical testing tool which is always seen as a testing standard. Since the optical phase distribution can be easily extracted by several interferograms, the phase shifting interferometry (PSI) has been widely used in optical measurement [1–3]. The performance of PSI mainly depends on the interferometer, environment and phase shifting algorithm (PSA), for the fixed interferometer and environment, outstanding PSA can improve the performance of PSI. An outstanding PSA should have two advantages which are high accuracy and timesaving feature, outstanding PSA could suppress the different kinds of errors to improve the accuracy, such as the miscalibration of piezo-transducer (PZT), detector error, vibrational error, air turbulence in the working environment, instability of the laser, etc., moreover, if it is easy to implement and only needs less phase shifted interferograms, it will save time, lastly, the PSA with random phase shift will be insensitive to the phase shift error due to the miscalibration of PZT, vibrational error, air turbulence, instability of the laser frequency, hence, a fast and accurate random PSA is essential for the high-quality optical testing with interferometry.

In recent years, many random PSAs have been developed, it can be divided into the iterative and non-iterative PSAs. Generally, the accuracy of the iterative PSA is relatively high, but it costs more time because of the iterative operation, and the non-iterative PSA spends less time than the iterative PSA, but the accuracy may be not as high as the iterative PSA. For the optical metrology, especially for the in situ metrology, the instantaneity of PSA is very important, hence, only a small number of iterative PSAs have been developed [4–6], most scientists were committed to research the non-iterative random PSA with high accuracy.

In 1992, Farrell and Player [7] utilized Lissajous figures and ellipse fitting to calculate the phase difference between two interferograms, and in 2016, Liu et al. [8] proposed a PSA which can simultaneously extract the tested phase and phase shift from only two interferograms using Lissajous figure and ellipse fitting technology, but these two algorithms both need pre-filtering and the non-uniform intensity distribution will affect the accuracy. From 2003 to 2014, Cai et al. [9–17] proposed a series of statistical algorithms which can extract the phase shifts and tested phase, however, most of these algorithms need to know the intensities of object and reference, and many approximations in the PSAs will affect the accuracy. From 2011 to 2017 [18–24], proposed a series of PSAs based on principal component analysis (PCA) which is an efficient technique for phase extraction by converting a set of possibly correlated variables into a set of values of uncorrelated variables, but it needs more than three interferograms with the phase shift well distributed between 0 and 2π because it needs to subtract relatively accurate mean-background intensity, hence, the more the interferograms, the higher the accuracy is, but more interferograms will cost more acquisition time and computational time, it is difficult to obtain the high accuracy and high speed simultaneously. In 2012 [25], presented a two-step demodulation based on the Gram-Schmidt orthonormalization method (GS2), it requires subtracting the DC term by filtering before performing GS2. In 2014 [26], proposed an advanced GS method called GS3, it needs three phase shifted interferograms, and the major advantage of this method is that it performs well when the phase shift is close to π as most two-step algorithms become invalid in this situation.

From the above literatures, we found a phenomenon that many PSAs need pre-filtering before the phase extraction, the pre-filtering process will cost more time, and it may affect the accuracy, in addition, some PSAs which don’t need pre-filtering may need more than three interferograms, they will also spend more time. Generally, the PSA with more interferograms has high accuracy and low speed, and the PSA with less interferograms has high speed and low accuracy, it is difficult to obtain the high accuracy and speed simultaneously. To balance the computational time and accuracy, the research of non-iterative random PSA with less interferograms and no pre-filtering is essential.

In this paper, we will discuss a fast and accurate two-step PSA with random phase shift. Section 2 presents the principle and process of the proposed PSA based on principal component analysis and Lissajous ellipse fitting (PCA&LEF). In Section 3 the simulation of PCA&LEF is discussed, and the comparison of PCA&LEF with GS is performed. Section 4 evaluates the novel PCA&LEF with the experimental data. The conclusion is finally drawn in Section 5.

2. Principles

The intensity expressions of the phase shifted interferograms are

I_{n} (x, y) = a_{n} (x, y) + b {}_{n}{(x, y)} \cos (φ (x, y) + δ_{n}) + ξ_{n} (x, y) .

where

I_{n} (x, y)

is the n^th phase shifted interferogram with size of

N_{x} \times N_{y}

, n = 1,2,...,N represents the image index with N the total number of phase shifted interferograms, N is set to 2,

a_{n} (x, y)

and

b_{n} (x, y)

respectively represent the background intensity and modulation amplitude of the phase shifted interferograms,

φ (x, y)

is the tested phase,

ξ_{n} (x, y)

is the noise, and

δ_{n}

is the phase shift. For convenience, the spatial coordinate has been omitted in the following.

Equation (1) can be rewritten as

I_{n} = α_{n} I_{c} + β_{n} I_{s} + κ_{n}

where

α_{n} = \cos (δ_{n})

,

β_{n} = - \sin (δ_{n})

,

I_{c} = b_{n} \cos (φ)

,

I_{s} = b_{n} \sin (φ)

and

κ_{n} = a_{n} + ξ_{n}

.

For PCA, the first step is filtering the background intensity $a_{n}$ by subtracting the average of all the phase shifted interferograms, however, the background intensity can be well eliminated only when the phase shift is well distributed between 0 and 2π and the number of the phase shifted interferograms is large enough, when the phase shift is “randomly” distributed and the number of phase shifted interferograms is small, the phase extracted by PCA will be not accurate, and especially when there are only two phase shifted interferograms, the real phase can't be extracted by PCA since the background intensity can’ be eliminated in this situation. Hence, we design a new method which can extract the accurate phase without the background intensity filtering by PCA, and it only needs two randomly phase shifted interferograms.

From Eq. (2), we can see that the intensity of the phase shifted interferogram can be expressed as a linear combination of two signals, the background intensity and noise term. Then, we can rewrite Eq. (2) as

I = Q Γ + ℜ .

where

Q = [p, q]

with size of

N_{x} N_{y} \times 2

, p and q are column vectors with size of

N_{x} N_{y} \times 1

whose elements are taken columnwise from

I_{c}

and

I_{s}

,

Γ = {[α, β]}^{T}

with size of

2 \times N

,

α

and

β

are the column vectors with size of

N \times 1

whose elements are taken columnwise from

α_{n}

and

β_{n}

,

ℜ

is the background intensity and noise matrix with size of

N_{x} N_{y} \times N

,where the n^th column is taken columnwise from

κ_{n}

, lastly,

I

is a matrix with size of

N_{x} N_{y} \times N

,where the n^th column is taken columnwise from

I_{n}

.

The covariance matrix C can be expressed as

C = I^{T} I = {(Q Γ + ℜ)}^{T} (Q Γ + ℜ) = Γ^{T} Q^{T} Q Γ + Γ^{T} Q^{T} ℜ + ℜ^{T} Q Γ + ℜ^{T} ℜ \approx Γ^{T} Q^{T} Q Γ + ℜ^{T} ℜ

The product of two uncorrelated matrixes-

Γ^{T} Q^{T} ℜ

and

ℜ^{T} Q Γ

can be ignored because they are significantly smaller than

Γ^{T} Q^{T} Q Γ

and

ℜ^{T} ℜ

.

Γ Γ^{T} = [\begin{matrix} α \\ β \end{matrix}] [\begin{matrix} α & β \end{matrix}] = (\begin{matrix} {‖ α ‖}^{2} & 〈 α, β 〉 \\ 〈 α, β 〉 & {‖ β ‖}^{2} \end{matrix}) .

where

‖ \cdot ‖

and

〈 \cdot, \cdot 〉

represent the norm and inner product,

{‖ α ‖}^{2} = \sum_{n = 1}^{N} \cos^{2} (δ_{n})

,

{‖ β ‖}^{2} = \sum_{n = 1}^{N} \sin^{2} (δ_{n})

, and

〈 α, β 〉 = - \sum_{n = 1}^{N} \cos (δ_{n}) \sin (δ_{n})

.

Note that, $Γ Γ^{T}$ is real and symmetric matrix, it can be diagonalized as $Γ Γ^{T} = P_{Γ}^{T} D_{Γ} P_{Γ}$ , where $D_{Γ}$ and $P_{Γ}$ are diagonal and orthogonal matrices.

Γ Γ^{T} = P_{Γ}^{T} D_{Γ} P_{Γ} = (P_{Γ}^{T} D_{Γ}^{1 / 2} \hat{Γ}) {({\hat{Γ}}^{T} {(D_{Γ}^{1 / 2})}^{T} P_{Γ})}^{T} .

_$\hat{Γ}$is a new matrix that

\hat{Γ} {\hat{Γ}}^{T} = E

, where

E

is unit matrix. Hence, according to Eq. (6), we can get the expression of

Γ

as

Γ = P_{Γ}^{T} D_{Γ}^{1 / 2} \hat{Γ} .

The new matrix

\hat{Γ}

can be expressed as

\hat{Γ} = D_{Γ}^{- 1 / 2} P_{Γ} Γ .

The

D_{Γ}

is given by

D_{Γ} = (\begin{matrix} λ_{1} & 0 \\ 0 & λ_{2} \end{matrix}) .

with

λ_{1, 2} = \frac{({‖ α ‖}^{2} + {‖ β ‖}^{2}) \pm \sqrt{{({‖ α ‖}^{2} - {‖ β ‖}^{2})}^{2} + 4 {〈 α \cdot β 〉}^{2}}}{2} .

and

P_{Γ}

is given by

P_{Γ} = (\begin{matrix} \frac{(λ_{1} - {‖ α ‖}^{2})}{\sqrt{{〈 α \cdot β 〉}^{2} + {(λ_{1} - {‖ α ‖}^{2})}^{2}}} & \frac{〈 α \cdot β 〉}{\sqrt{{〈 α \cdot β 〉}^{2} + {(λ_{2} - {‖ β ‖}^{2})}^{2}}} \\ \frac{〈 α \cdot β 〉}{\sqrt{{〈 α \cdot β 〉}^{2} + {(λ_{1} - {‖ α ‖}^{2})}^{2}}} & \frac{(λ_{2} - {‖ β ‖}^{2})}{\sqrt{{〈 α \cdot β 〉}^{2} + {(λ_{2} - {‖ β ‖}^{2})}^{2}}} \end{matrix}) .

Additionally,

Q^{T} Q

can be expressed as

Q^{T} Q = (\begin{matrix} p \\ q \end{matrix}) (\begin{matrix} p & q \end{matrix}) = (\begin{matrix} {‖ p ‖}^{2} & 〈 p \cdot q 〉 \\ 〈 p \cdot q 〉 & {‖ q ‖}^{2} \end{matrix}) = \sum_{N_{x} \times N_{y}} (\begin{matrix} b^{2} \cos^{2} (φ) & b^{2} \cos (φ) \sin (φ) \\ b^{2} \cos (φ) \sin (φ) & b^{2} \sin^{2} (φ) \end{matrix}) .

If we have more than one fringe in the interferograms, we can use the approximation

\sum_{N_{x} \times N_{y}} b^{2} \cos (φ) \sin (φ) \approx 0.

and

\sum_{N_{x} \times N_{y}} b^{2} \cos^{2} (φ) \approx \sum_{N_{x} \times N_{y}} b^{2} \sin^{2} (φ) \approx σ .

Hence, Eq. (12) can be rewritten as

Q^{T} Q = D_{Q} \approx σ (\begin{matrix} 1 & 0 \\ 0 & 1 \end{matrix}) .

According to Eqs. (4), (7) and (15), we have

C \approx {\hat{Γ}}^{T} D_{Γ}^{1 / 2} P_{Γ} D_{Q} P_{Γ}^{T} D_{Γ}^{1 / 2} \hat{Γ} + ℜ^{T} ℜ .

Since

D_{Q}

is an approximate diagonal matrix, and

P_{Γ}^{T} P_{Γ} = P_{Γ} P_{Γ}^{T} = E

, Eq. (16) can be rewritten as

C \approx {\hat{Γ}}^{T} (D_{Γ} D_{Q} + D_{ℜ}) \hat{Γ} .

where

ℜ^{T} ℜ = {\hat{Γ}}^{T} {(ℜ {\hat{Γ}}^{T})}^{T} (ℜ {\hat{Γ}}^{T}) \hat{Γ} = {\hat{Γ}}^{T} D_{ℜ} \hat{Γ} .

According to Eq. (7), Eq. (3) can be rewritten as

I \approx Q (P_{Γ}^{T} D_{Γ}^{1 / 2} \hat{Γ}) + ℜ = (Q P_{Γ}^{T}) D_{Γ}^{1 / 2} \hat{Γ} + ℜ .

PCA is a technique from statistics for reducing an image or data set that transforms a number of uncorrelated images into the smallest number uncorrelated images called the principle components. Since the covariance matrix C is a real and symmetric matrix, it can be diagonalized as

C = U^{T} D U .

where U and D are orthogonal and diagonal matrices.

The principle components of the interferograms are given by

Z = I U^{T} .

where Z is matrix with size of

N_{x} N_{y} \times N

, and its column vectors z_n are the principle components.

The phase will be extracted by PCA using the first and second components (z₁ and z₂), that corresponds to the highest eigenvalues, as Eq. (22).

φ = arc \tan (\frac{I_{s}}{I_{c}}) = arc \tan (\frac{z_{2}}{z_{1}}) .

According to Eqs. (17) and (20), we can state that U and D correspond to

\hat{Γ}

and

D_{Γ} D_{Q} + D_{ℜ}

respectively. Then, we can rewrite Eq. (21) as

z_{i} = I {\hat{Γ}}^{T} = {((Q P_{Γ}^{T}) D_{Γ}^{1 / 2} + ℜ \hat{Γ})}_{i} = {(\hat{Q} D_{Γ}^{1 / 2} + ℜ \hat{Γ})}_{i}, i = 1, 2.

where

\hat{Q} = Q P_{Γ}^{T}

.

We know that $P_{Γ}^{T}$ is a 2 × 2 orthogonal matrix, so it can be expressed as

P_{Γ}^{T} = (\begin{matrix} \cos (θ) & \sin (θ) \\ - \sin (θ) & \cos (θ) \end{matrix}) .

Then

\begin{array}{l} \hat{Q} = Q P_{Γ}^{T} = (\begin{matrix} b \cos (φ) & b \sin (φ) \end{matrix}) (\begin{matrix} \cos (θ) & \sin (θ) \\ - \sin (θ) & \cos (θ) \end{matrix}) . \\ = (\begin{matrix} b \cos (φ) \cos (θ) - b \sin (φ) \sin (θ) & b \cos (φ) \sin (θ) + b \sin (φ) \cos (θ) \end{matrix}) \\ = (\begin{matrix} b \cos (φ + θ) & b \sin (φ + θ) \end{matrix}) = (\begin{matrix} \hat{p} & \hat{q} \end{matrix}) \end{array}

Then,

\begin{array}{l} z_{1} = \hat{p} λ_{1}^{1 / 2} + {(ℜ \hat{Γ})}_{1} = b \cos (φ + θ) λ_{1}^{1 / 2} + {(ℜ \hat{Γ})}_{1} . \\ z_{2} = \hat{q} λ_{2}^{1 / 2} + {(ℜ \hat{Γ})}_{2} = b \sin (φ + θ) λ_{2}^{1 / 2} + {(ℜ \hat{Γ})}_{2} \end{array}

Then we can obtain

\begin{array}{l} \cos (φ + θ) = \frac{z_{1} - {(ℜ \hat{Γ})}_{1}}{b λ_{1}^{1 / 2}} . \\ \sin (φ + θ) = \frac{z_{2} - {(ℜ \hat{Γ})}_{2}}{b λ_{2}^{1 / 2}} \end{array}

Because

\sin^{2} (φ + θ) + \cos^{2} (φ + θ) = 1

, and we use X and Y instead of z₂ and z₁, then Eq. (27) can be rewritten as

{(\frac{X - x_{0}}{a_{x}})}^{2} + {(\frac{Y - y_{0}}{a_{y}})}^{2} = 1.

Note that Eq. (28) is just an ellipse equation,

a_{x} = b λ_{2}^{1 / 2}, a_{y} = b λ_{2}^{1 / 2}, x_{0} = {(ℜ \hat{Γ})}_{2}, y_{0} = {(ℜ \hat{Γ})}_{1} .

Equation (28) can be expanded as a general conic function:

\frac{1}{a_{x}^{2}} X^{2} + \frac{1}{a_{y}^{2}} Y^{2} - 2 \frac{x_{0}}{a_{x}^{2}} X - 2 \frac{y_{0}}{a_{y}^{2}} Y + \frac{x_{0}^{2}}{a_{x}^{2}} + \frac{y_{0}^{2}}{a_{y}^{2}} - 1 = 0.

A general conic function can be also expressed by the following second order polynomial:

F = c x^{2} + d x y + e y^{2} + f x + g y + h .

For an ellipse, Eq. (31) needs to meet the conditions of

F = 0

and

d^{2} - 4 c e < 0

. In the following, the real phase distribution will be obtained by the Lissajous ellipse fitting (LEF) method. Firstly, by plotting z₂ against z₁ in a Cartesian coordinate, a Lissajous ellipse which is just the representation of Eq. (28) will be created, and it will be easy to calculate the conic coefficients of Eq. (31) by the least squares algorithm, then, the semi-major amplitude

a_{x}

, semi-minor amplitude

a_{y}

, the center offset

x_{0}

and

y_{0}

can be calculated as

\begin{array}{l} a_{x} = \sqrt{2 \frac{c g^{2} + e f^{2} + h d^{2} - d f g - 4 c e h}{(d^{2} - 4 c e) (\sqrt{{(c - e)}^{2} + d^{2}} - (c + e))}}, a_{y} = \sqrt{2 \frac{c g^{2} + e f^{2} + h d^{2} - d f g - 4 c e h}{(d^{2} - 4 c e) (- \sqrt{{(c - e)}^{2} + d^{2}} - (c + e))}} . \\ x_{0} = \frac{2 e f - d g}{d^{2} - 4 c e}, y_{0} = \frac{2 c g - d f}{d^{2} - 4 c e} \end{array}

Lastly, according to Eqs. (27) and (28), the phase

Φ

can be easily calculated as

Φ = φ + θ = \tan^{- 1} (\frac{X - x_{0}}{Y - y_{0}} \cdot \frac{a_{y}}{a_{x}}) .

We know that there is only a constant

θ

between between

φ

and

Φ

, which doesn’t affect the whole phase distribution, hence we can use

Φ

to express the tested phase distribution.

Based on the principles of PCA and LEF, we propose a novel method called principle components analysis and Lissajous ellipse fitting algorithm (PCA&LEF), it can extract the phase distribution from two randomly phase shifted interferograms without background intensity filtering, it removes the restriction that PCA needs more than three interferograms with well distributed phase shifts to subtract relatively accurate mean-background intensity. In the following, we will introduce the process of the proposed method in detail:

1) Generate a matrix I with size of $N_{x} N_{y} \times N$ , where the n^th column is taken columnwise from the intensity of n^th phase shifted interferograms $I_{n}$ ;
2) calculate the covariance matrix C by equation $C = I^{T} I$ ;
3) calculate the orthogonal matrix $\hat{Γ}$ including the eigenvectors of the covariance matrix C;
4) obtain the first and second principle components (z₁ and z₂) which corresponds to the highest eigenvalues by $z_{i} = I {\hat{Γ}}^{T}$ ;
5) plot an approximate ellipse with z₂ as the x coordinate and z₁ as the y coordinate;
6) calculate the semi-major amplitude $a_{x}$ , semi-minor amplitude $a_{y}$ , the center offset $x_{0}$ and $y_{0}$ of the Lissajous ellipse using LEF by Eq. (32);
7) calculate the phase distribution using Eq. (33).

The whole procedure of the PCA&LEF is illustrated in Fig. 1.

Fig. 1 Flow chart of PCA&LEF.

Download Full Size | PDF

3. Simulation

To verify the effectiveness of the method proposed above, we perform a series of numerical simulations, and compare it with GS which is a well-evaluated two-step random PSA, note that, the Hilbert-Huang pre-filtering will be performed before using GS. In the following, all computations are performed with the CPU of Intel(R) Core(TM) i7-6700 and the 8 GB memory, and we use the Matlab software for coding; Code 1 of PCA&LEF is provided in Ref [31].

We perform PCA&LEF and GS to process two phase shifted interferograms with different situations. In the following, the tested phase is set as $φ = N_{f} π (x^{2} + y^{2})$ , in which $N_{f} = 5$ is the fringe number in the interferogram. Figure 2(a) shows the theoretical phase distribution. In situation 1, the background intensity and modulation amplitude distributions of the two interferograms are uniform, $a_{1} = a_{2} = 1$ , $b_{1} = b_{2} = 1$ . In situation 2, the background intensity and modulation amplitude are set as $a_{i} (x, y) = N_{a} \exp [- 0.02 (x^{2} + y^{2})]$ and $b_{i} (x, y) = N_{b} \exp [- 0.02 (x^{2} + y^{2})]$ respectively, both the fluctuation and non-uniformity of the background intensity and modulation amplitude exist, hence, $N_{a}$ of the 1st and 2nd interferograms are set as 1 and 0.95, $N_{b}$ of the 1st and 2nd interferograms are set as 0.9 and 0.85. In situation 3, we add noise with SNR of 20dB which is generated by the function awgn in Matlab to situation 2. With the above parameters setting, two simulated phase shifted interferograms with the size of 401 × 401 and the phase shift of 2 rad in situation 3 are generated, as shown in Figs. 2(b) and 2(c).

Fig. 2 Simulated phase distribution and two phase shifted interferograms. (a) The theoretical phase distribution (PV = 31.4159 rad, RMS = 6.6561 rad), (b) and (c) the first interferogram and the second interferogram.

Download Full Size | PDF

Figure 3 shows the results of three different situations, N represents the index of the situation, the 1st and 3rd rows show the phase distributions calculated by PCA&LEF and GS, and the phase error distributions are displayed in the 2nd and 4th rows. The detailed RMS phase errors and computational time are represented in Table 1. We can see that the phase distributions are similar for the two different methods in the three situations, that is to say, the two methods are both robust for different situations. However, the phase error distributions are different for the different situations and methods.

Fig. 3 The phase distributions and phase error distributions calculated by PCA&LEF and GS in different situations.

Download Full Size | PDF

Table 1. The RMS phase errors and computational time of PCA&LEF and GS in different situations

View Table | View all tables in this article

From Table1, we can see that the accuracy of PCA&LEF is higher than that of GS in different situations since there is no filtering in PCA&LEF and LEF can suppress some errors to extract the accurate phase distribution. Moreover, the computational time of PCA& LEF is less than that of GS because the filtering before performing GS costs more time. Through the three situations, we found that the RMS phase errors are similar for GS because the main error of GS is filtering error, however, for the RMS phase errors of PCA&LEF, there are big differences in different situations, the RMS phase error is close to zero in situation 1, it means that PCA&LEF can obtain absolutely accurate result when the background intensity and modulation amplitude are perfect, and the RMS phase error of PCA&LEF in situation 2 is larger than that in situation 1, but it is only 0.0127 rad, that is to say, PCA&LEF can partly suppress the effect of the background intensity and modulation amplitude fluctuation and non-uniformity, lastly, situation 3 is most complex, the non-uniformity between different pixels, fluctuation error between different interferograms and noise are all added to the interferograms, the mixed errors cause the largest phase error, the RMS phase error which is 0.1307 rad is more than 10 times of situation 2, that is to say, LEF can’t suppress the effect of noise.

From Fig. 1, we know that there are 7 steps for PCA&LEF, to further study the computational time of PCA&LEF, we respectively calculate the computational time of every step in three different situations, as shown in Table 2. We found that the computational time of the front five steps are relatively small, especially the front four steps, the computational time can be ignored. Step 6 which is the LEF process costs most time, although LEF process spends many time, it increases the accuracy and avoids the pre-filtering, if the pre-filtering is used, it will cost more time, and the accuracy will be decreased, GS is an example. Moreover, step 7 costs more time because the unwrapped process is the main part of the phase calculation, and the unwrapped process is unavoidable for every PSA. From the above study, we can conclude that PCA hardly costs time, LEF process and phase calculation cost more time in PCA&LEF.

Table 2. The computational time of every step for PCA&LEF in different situations

View Table | View all tables in this article

Figure 4 shows the ellipses before and after using LEF for PCA&LEF in different situations, we can see that, before using LEF, the approximate ellipse with X as the x coordinate and Y as the y coordinate is not centered at the origin, after using LEF, the ellipse was transformed an approximate circle with $(X - x_{0}) / a_{x}$ as the x coordinate and $(Y - y_{0}) / a_{y}$ as the y coordinate centered at the origin. For situations 1 and 2, whether the ellipse before using LEF or the circle after using LEF, the curve is smooth, but the curve is not smooth for the last situation since the noise exists, and LEF cannot remove this effect.

Fig. 4 The ellipses before and after using LEF for PCA&LEF in different situations.

Download Full Size | PDF

In order to analyze the effect of different phase shifts and different levels of noises for two methods, we calculate the RMS phase errors of PCA&LEF and GS with different phase shifts in 6 different situations, the range of phase shift is between 0.2 rad and 6.2 rad, and the situations 1, 2 and 3 are same as the above simulation, we additionally discuss the different levels of noises in situation 3, such as 30 dB, 40 dB and 50 dB. The results are presented in Fig. 5. Figure 5(a) shows the results of PCA&LEF, we can see that, the more complex the situation, the larger the RMS phase error is, for situation 1, the RMS phase errors approximate to zero for most different phase shifts except for 3.2 rad which is close to π rad, for situation 2, the RMS phase errors are larger than that in situation 1 due to the fluctuation and non-uniformity of the background intensity and modulation amplitude, for situation 3, the RMS phase error is increasing with the increase of the SNR of noise since LEF can’t suppress the noise, moreover, for situation1, 2 and 3 (50dB), PCA&LEF is workable for the different phase shifts, however, for situation 3 with larger noise, PCA&LEF sometimes doesn’t work in some special phase shifts which is close to 0 rad, π rad and 2π rad since PCA&LEF can’t suppress the noise, and the noise affect the fitting of Lissajous ellipse, the larger the noise, the smaller the working range of phase shift is, for situation 3 (20dB), 0 rad to 0.6 rad, 2.8rad to 3.4 rad, and 5.8 rad to 6.2 rad are invalid range, we found a common feature that the nearer the phase shift to 0 rad, π rad and 2π rad, the larger the RMS phase error is. Figure 5(b) shows the results of GS, the RMS phase errors of GS have no concern with the different situations because the filtering error is the largest error, and for different situations, they have the common invalid range of phase shift, they are 2.8 rad to 3.8 rad, 6.0 rad to 6.2 rad, hence, GS is not always valid with random phase shift even through the experimental environment and interferometer are perfect.

Fig. 5 RMS phase errors of PCA&LEF and GS with different phase shifts in different situations.

Download Full Size | PDF

Moreover, we also compare PCA&LEF and GS in the above 6 situations, the compared results are shown in Fig. 6, we can see that, in situation 1, the RMS phase errors with random phase shift are far smaller than that of GS, and in situations 2 and 3, the results are similar as situation 1 except 3.2 rad, in situations 4 and 5, most RMS phase errors of PCA&LEF are also smaller than that of GS except the phase shifts which is close to 0 rad, π rad and 2π rad, lastly, in situation 6, the orders of magnitudes for two methods are similar, that is to say, the effect of the noise for PCA&LEF and the effect of the filtering error for GS are similar. If the relatively high accuracy is demanded, it would be best to suppress the noise before using the proposed method, and it's best to choose a phase shift which is far away from 0 rad, π rad and 2π rad to further increase the accuracy.

Fig. 6 The compared results of PCA&LEF and GS in different situations.

Download Full Size | PDF

4. Demonstration with experimental data

In order to verify the performance of the proposed method, three groups of experiments are performed to do the phase retrieval by the proposed method and GS. The experimental setup is Twyman-Green interferometer with 4D camera which is a kind of synchronous phase-shifting interferometer (SPSI) [27], four phase-shifted interferograms with the phase shifts 0, π/2, π and 3π/2 can be extracted from a single image snapshotted by the 12 bit polarization camera PolarCam with the pixel number of 1208 × 1348 and the pixel size of 7.4 µm from 4D Technology, Inc [28–30]. Moreover, we test different objects to capture the circular, straight and complex fringes. For the first experiment, the phase shifted interferograms with the circular fringes are collected, the size of the interferograms is 301 × 301, and the phase extracted by standard 4-step PSA is set as the reference phase due to its high accuracy. One of the interferograms is shown in Fig. 7(a), Fig. 7(b) shows the reference phase distribution, and the phase distributions extracted by PCA&LEF and GS are drawn in Figs. 7(c) and 7(d). The differences between the reference phase and the phase extracted by PCA&LEF and GS are shown in Figs. 7(e) and 7(f), the RMS values of the differences are respectively 0.0581 rad and 0.1805 rad, further verifying the accuracy of PCA&LEF is higher than that of GS. Moreover, the computational time of PCA&LEF (1.2458 s) is less than that of GS (1.8681 s), actually, the computational time of GS is only 1s, but the Hilbert-Huang pre-filtering before performing GS costs long time which is 1.40s, and the filtering process introduces the extra phase error to the result. For PCA&LEF, the ellipses before and after using LEF are plotted in Figs. 7(g) and 9(h).

Fig. 7 Experimental results of the circular fringes . (a) One of the phase shifted interferograms, (b) the reference phase distribution extracted by 4-step PSA (PV = 24.9105rad, RMS = 5.0283 rad), (c) and (d) the phase distributions extracted by PCA&LEF (PV = 25.1416 rad, RMS = 5.0314 rad) and GS (PV = 25.1539 rad, RMS = 4.9495 rad), (e) and (f) the differences between the reference and phase distributions extracted by PCA&LEF and GS, (g) and (h) the ellipses before and after using LEF.

Download Full Size | PDF

Then, the second and third experiments with the straight and relatively complex fringes are performed, the size of the interferograms with the straight fringes is 401 × 401, and the size of the interferograms with the complex fringes is 201 × 201, other conditions are same as the above circular fringes. Figures 8 and 9 show the results of the straight and complex fringes, we can see that, both PCA&LEF and GS are effective for the different fringes. Moreover, for the straight fringes, the RMS values of the differences between the reference phase and the phase extracted by PCA&LEF and GS are 0.0453 rad and 0.2392 rad, and the computational time of PCA&LEF and GS are 1.6769 s and 3.2142 s respectively. And, for the complex fringes, the RMS values of the differences for PCA&LEF and GS are 0.0832 rad and 0.1116 rad, and the computational time are 0.9674 s and 1.2804 s. For these two kinds of fringes, we get the same conclusion as the circular fringes in regard to the accuracy and computational time. Through the above experiments, we verify that, both PCA&LEF and GS are suitable for the circular, straight and complex fringes, and the proposed PCA&LEF without pre-filtering can obtain relatively accurate result with less computational time by only two interferograms.

Fig. 8 Experimental results of the straight fringes. (a) One of the phase shifted interferograms, (b) the reference phase distribution extracted by 4-step PSA (PV = 22.2217rad, RMS = 6.0449 rad), (c) and (d) the phase distributions extracted by PCA&LEF (PV = 22.1837 rad, RMS = 6.0488 rad) and GS (PV = 24.0806 rad, RMS = 6.0492 rad), (e) and (f) the differences between the reference and phase distributions extracted by PCA&LEF and GS, (g) and (h) the ellipses before and after using LEF.

Download Full Size | PDF

Fig. 9 Experimental results of the complex fringes . (a) One of the phase shifted interferograms, (b) the reference phase distribution extracted by 4-step PSA (PV = 40.4129rad, RMS = 8.4209 rad), (c) and (d) the phase distributions extracted by PCA&LEF (PV = 40.7467 rad, RMS = 8.4380 rad) and GS (PV = 24.0806 rad, RMS = 6.0492 rad), (e) and (f) the differences between the reference and phase distributions extracted by PCA&LEF and GS, (g) and (h) the ellipses before and after using LEF.

Download Full Size | PDF

Finally, we study the computational time of every step for PCA&LEF with the different fringes, as shown in Table 3. The conclusions are same as the simulation, for the different fringes, steps 6 and 7 also cost more time than other steps, and the computational time of every step are different for the different fringes because the sizes of the interferograms are different.

Table 3. The computational time of every step for PCA&LEF with the different fringes

View Table | View all tables in this article

5. Conclusion

In this paper, we present a PSA based on principal component analysis and Lissajous ellipse fitting, PCA is firstly used without subtracting the mean-background intensity, then the LEF process is performed to extract the real phase distribution. We have compared PCA&LEF and well-evaluated GS by the simulated and experimental data. The proposed method can achieve high accuracy with only two randomly phase shifted interferograms and no pre-filtering, and it can directly obtain the tested phase with less computational time. Then, it removes the restriction that PCA needs more than three interferograms with well distributed phase shifts to subtract relatively accurate mean-background intensity. In addition, it is suitable for different background intensity, modulation amplitude distributions and noises, and it can obtain absolutely accurate phase distribution when the background intensity and modulation amplitude are perfect, also, it can partly suppress the effect of the imperfect background intensity and modulation amplitude. Lastly, if the higher accuracy is requested, the noise is best to be suppressed, and it's best to choose a phase shift which is far away from 0 rad, π rad and 2π rad. The simulations and experiments demonstrate the validity of the proposed method. In summary, this proposed method is a power tool for the phase retrieval with random phase shift.

Funding

National Natural Science Foundation of China (NSFC) (11304034); Department of Science and Technology of Jilin Province (20190701018GH); Education Department of Jilin Province (JJKH20190691KJ); State Key Laboratory of Applied Optics.

References

1. D. Malacara, Optical Shop Testing, 3rd ed. (John Wiley & Sons, Inc., 2007).

2. J. H. Bruning, D. R. Herriott, J. E. Gallagher, D. P. Rosenfeld, A. D. White, and D. J. Brangaccio, “Digital wavefront measuring interferometer for testing optical surfaces and lenses,” Appl. Opt. 13(11), 2693–2703 (1974). [CrossRef] [PubMed]

3. C. Tian and S. Liu, “Two-frame phase-shifting interferometry for testing optical surfaces,” Opt. Express 24(16), 18695–18708 (2016). [CrossRef] [PubMed]

4. Z. Wang and B. Han, “Advanced iterative algorithm for phase extraction of randomly phase-shifted interferograms,” Opt. Lett. 29(14), 1671–1673 (2004). [CrossRef] [PubMed]

5. J. Xu, Q. Xu, and L. Chai, “Iterative algorithm for phase extraction from interferograms with random and spatially nonuniform phase shifts,” Appl. Opt. 47(3), 480–485 (2008). [CrossRef] [PubMed]

6. Y. C. Chen, P. C. Lin, C. M. Lee, and C. W. Liang, “Iterative phase-shifting algorithm immune to random phase shifts and tilts,” Appl. Opt. 52(14), 3381–3386 (2013). [CrossRef] [PubMed]

7. C. T. Farrell and M. A. Player, “Phase step measurement and variable step algorithms in phase-shifting interferometry,” Meas. Sci. Technol. 3(10), 953–958 (1992). [CrossRef]

8. F. Liu, J. Wang, Y. Wu, F. Wu, M. Trusiak, K. Patorski, Y. Wan, Q. Chen, and X. Hou, “Simultaneous extraction of phase and phase shift from two interferograms using Lissajous figure and ellipse fitting technology with Hilbert-Huang prefiltering,” J. Optics-UK 18(10), 105604 (2016). [CrossRef]

9. L. Z. Cai, Q. Liu, X. L. Yang, and Y. R. Wang, “Sensitivity adjustable contouring by digital holography and a virtual reference wavefront,” Opt. Commun. 221(1-3), 49–54 (2003). [CrossRef]

10. L. Z. Cai, Q. Liu, and X. L. Yang, “Phase-shift extraction and wave-front reconstruction in phase-shifting interferometry with arbitrary phase steps,” Opt. Lett. 28(19), 1808–1810 (2003). [CrossRef] [PubMed]

11. L. Z. Cai, Q. Liu, and X. L. Yang, “Generalized phase-shifting interferometry with arbitrary unknown phase steps for diffraction objects,” Opt. Lett. 29(2), 183–185 (2004). [CrossRef] [PubMed]

12. X. F. Meng, L. Z. Cai, X. F. Xu, X. L. Yang, X. X. Shen, G. Y. Dong, and Y. R. Wang, “Two-step phase-shifting interferometry and its application in image encryption,” Opt. Lett. 31(10), 1414–1416 (2006). [CrossRef] [PubMed]

13. X. F. Xu, L. Z. Cai, X. F. Meng, G. Y. Dong, and X. X. Shen, “Fast blind extraction of arbitrary unknown phase shifts by an iterative tangent approach in generalized phase-shifting interferometry,” Opt. Lett. 31(13), 1966–1968 (2006). [CrossRef] [PubMed]

14. X. F. Xu, L. Z. Cai, Y. R. Wang, X. L. Yang, X. F. Meng, G. Y. Dong, X. X. Shen, and H. Zhang, “Generalized phase-shifting interferometry with arbitrary unknown phase shifts: Direct wave-front reconstruction by blind phase shift extraction and its experimental verification,” Appl. Phys. Lett. 90(12), 121124 (2007). [CrossRef]

15. X. F. Xu, L. Z. Cai, Y. R. Wang, X. F. Meng, W. J. Sun, H. Zhang, X. C. Cheng, G. Y. Dong, and X. X. Shen, “Simple direct extraction of unknown phase shift and wavefront reconstruction in generalized phase-shifting interferometry: algorithm and experiments,” Opt. Lett. 33(8), 776–778 (2008). [CrossRef] [PubMed]

16. X. F. Xu, L. Z. Cai, Y. R. Wang, and R. S. Yan, “Direct phase shift ectraction and wavefront reconstruction in two-step generalized phase-shifting interferometry,” J. Opt. 12(1), 015301 (2010). [CrossRef]

17. R. S. Yan, L. Z. Cai, and X. F. Meng, “Correction of wave-front retrieval errors caused by the imperfect collimation of reference beam in phase-shifting interferometry,” Optik (Stuttg.) 125(2), 601–605 (2014). [CrossRef]

18. J. Vargas, J. A. Quiroga, and T. Belenguer, “Phase-shifting interferometry based on principal component analysis,” Opt. Lett. 36(8), 1326–1328 (2011). [CrossRef] [PubMed]

19. J. Vargas, J. A. Quiroga, and T. Belenguer, “Analysis of the principal component algorithm in phase-shifting interferometry,” Opt. Lett. 36(12), 2215–2217 (2011). [CrossRef] [PubMed]

20. J. Deng, K. Wang, D. Wu, X. Lv, C. Li, J. Hao, J. Qin, and W. Chen, “Advanced principal component analysis method for phase reconstruction,” Opt. Express 23(9), 12222–12231 (2015). [CrossRef] [PubMed]

21. J. Xu, W. Jin, L. Chai, and Q. Xu, “Phase extraction from randomly phase-shifted interferograms by combining principal component analysis and least squares method,” Opt. Express 19(21), 20483–20492 (2011). [CrossRef] [PubMed]

22. K. Yatabe, K. Ishikawa, and Y. Oikawa, “Improving principal component analysis based phase extraction method for phase-shifting interferometry by integrating spatial information,” Opt. Express 24(20), 22881–22891 (2016). [CrossRef] [PubMed]

23. K. Yatabe, K. Ishikawa, and Y. Oikawa, “Simple, flexible, and accurate phase retrieval method for generalized phase-shifting interferometry,” J. Opt. Soc. Am. A 34(1), 87–96 (2017). [CrossRef] [PubMed]

24. K. Yatabe, K. Ishikawa, and Y. Oikawa, “Hyper ellipse fitting in subspace method for phase-shifting interferometry: Practical implementation with automatic pixel selection,” Opt. Express 25(23), 29401–29416 (2017). [CrossRef]

25. J. Vargas, J. A. Quiroga, C. O. Sorzano, J. C. Estrada, and J. M. Carazo, “Two-step demodulation based on the Gram-Schmidt orthonormalization method,” Opt. Lett. 37(3), 443–445 (2012). [CrossRef] [PubMed]

26. H. Wang, C. Luo, L. Zhong, S. Ma, and X. Lu, “Phase retrieval approach based on the normalized difference maps induced by three interferograms with unknown phase shifts,” Opt. Express 22(5), 5147–5154 (2014). [CrossRef] [PubMed]

27. Y. Zhang, X. Tian, and R. Liang, “Fringe-print-through error analysis and correction in snapshot phase-shifting interference microscope,” Opt. Express 25(22), 26554–26566 (2017). [CrossRef] [PubMed]

28. M. Novak, J. Millerd, N. Brock, M. North-Morris, J. Hayes, and J. Wyant, “Analysis of a micropolarizer array-based simultaneous phase-shifting interferometer,” Appl. Opt. 44(32), 6861–6868 (2005). [CrossRef] [PubMed]

29. D. Wang and R. Liang, “Simultaneous polarization Mirau interferometer based on pixelated polarization camera,” Opt. Lett. 41(1), 41–44 (2016). [CrossRef] [PubMed]

30. https://www.4dtechnology.com/products/polarimeters/polarcam/.

31. https://doi.org/10.6084/m9.figshare.8195426.

	N	1	2	3
RMS Phase error (rad)	PCA&LEF	3.4949 × 10⁻⁶	0.0127	0.1370
RMS Phase error (rad)	GS	0.1471	0.1491	0.1449
Computational time (s)	PCA&LEF	1.6307	1.6116	1.6058
Computational time (s)	GS	3.1117	3.0439	3.1555

Time (s)	Situation 1	Situation 2	Situation 3
Step 1	0.003110	0.002576	0.003138
Step 2	0.000586	0.000558	0.000563
Step 3	0.000633	0.000558	0.000581
Step 4	0.007928	0.005904	0.005288
Step 5	0.074457	0.074653	0.073545
Step 6	1.017182	1.002472	1.001089
Step 7	0.526767	0.524898	0.521588
Total	1.6307	1.6116	1.6058

Time (s)	Circular fringes	Straight fringes	Complex fringes
Step 1	0.001714	0.004360	0.002139
Step 2	0.000477	0.000651	0.000354
Step 3	0.000663	0.000648	0.000749
Step 4	0.004215	0.008880	0.003593
Step 5	0.051740	0.050191	0.050170
Step 6	0.768652	0.915817	0.504136
Step 7	0.418376	0.696384	0.406307
Total	1.2458	1.6769	0.9674

	N	1	2	3
RMS Phase error (rad)	PCA&LEF	3.4949 × 10⁻⁶	0.0127	0.1370
RMS Phase error (rad)	GS	0.1471	0.1491	0.1449
Computational time (s)	PCA&LEF	1.6307	1.6116	1.6058
Computational time (s)	GS	3.1117	3.0439	3.1555

Time (s)	Situation 1	Situation 2	Situation 3
Step 1	0.003110	0.002576	0.003138
Step 2	0.000586	0.000558	0.000563
Step 3	0.000633	0.000558	0.000581
Step 4	0.007928	0.005904	0.005288
Step 5	0.074457	0.074653	0.073545
Step 6	1.017182	1.002472	1.001089
Step 7	0.526767	0.524898	0.521588
Total	1.6307	1.6116	1.6058

Accurate and fast two-step phase shifting algorithm based on principle component analysis and Lissajous ellipse fitting with random phase shift and no pre-filtering

Abstract

1. Introduction

2. Principles

3. Simulation

4. Demonstration with experimental data

5. Conclusion

Funding

References

Supplementary Material (1)

Cited By

Figures (9)

Tables (3)

Equations (33)

Optics Express