Two-step random phase retrieval approach based on Gram-Schmidt orthonormalization and Lissajous ellipse fitting method

: To overcome the phase shift error in phase shifting interferometry, a two-step random phase retrieval approach based on Gram-Schmidt (GS) orthonormalization and Lissajous ellipse fitting (LEF) method (GS&LEF) is proposed. It doesn’t need pre-filtering, and can obtain relatively accurate phase distribution with only two phase shifted interferograms and less computational time. It is suitable for different background intensity, modulation amplitude distributions and noises. Last but not least, this method is effective for circular, straight or complex fringes. The simulations and experiments verify the correctness and feasibility of GS&LEF.


Introduction
Phase-shifting interferometry (PSI) has been widely used in the high precision phase measurement, such as optical surface testing, deformation measurement [1][2][3].In order to obtain accurate measured phase, the outstanding phase-shifting algorithm (PSA) is essential.For the standard PSA, its accuracy depends on the accuracy of the phase shift [4][5][6], which should be a special constant (e.g.π/2).However, the practical phase shift is difficult to equal to the theoretical phase shift due to the phase shift error caused by the miscalibration of piezotransducer, vibrational error, air turbulence in the working environment, instability of the laser frequency, and so on [7][8][9].
In recent years, to suppress the phase shift error, kinds of random PSAs have been developed [10][11][12][13][14][15][16][17], they can be divided to two types.The first type is multi-step random PSA, where the number of the phase shifted interferograms is greater than or equal to 3. While it can obtain the tested phase more easily, the phase error may be introduced due to the multiinterferograms, and it spends more time on the measurement and calculation [18].introduced an advanced iterative algorithm (AIA) to extract phase distribution from randomly phase shifted interferograms, it is based on a least-squares iterative procedure.The algorithm provides stable convergence and accurate phase extraction with as few as three interferograms, even when the phase shifts are completely random.However, this approach consumes a lot of computational time [19][20][21][22][23][24][25].proposed a series of PSAs based on principal component analysis (PCA), which can fast and easily extract the phase distribution from randomly phase shifted interferograms.The PCA is an efficient technique for phase extraction by converting a set of possibly correlated variables into a set of values of uncorrelated variables.The main drawback of this approach is that it cannot determine the global sign of the measured phase.
The second type is the two-step random PSA [26].presented a two-step demodulation based on the Gram-Schmidt orthonormalization method (GS), where phase-shift value is random and can be any value.The main drawback is that it requires subtracting the DC term by filtering, and filtering introduces error.In [11], Farrell and Player utilized Lissajous figures and ellipse fitting to calculate the phase difference between two interferograms, but it also needs pre-filtering and the correction result is not accurate if the intensity distribution is nonuniform.In [27], Liu et al. proposed a PSA which can simultaneously extract the tested phase and phase shift from only two interferograms using Lissajous figure and ellipse fitting technology, however, the two interferograms used in this algorithm need to be filtered by the Hilbert-Huang pre-filtering.Most two-step PSAs need pre-filtering to subtract the background intensity, the pre-filtering costs more time and may introduce extra errors, but two-step PSAs will save the measuring time since they only need two interferograms.
Recently, we proposed a random two-step PSA based on Lissajous ellipse fitting and least squares technologies [16], although this algorithm uses only two interferograms to extract the relatively accurate tested phase distribution and unknown phase shift without pre-filtering, the least squares technologies are time consuming.
For two-step random PSAs, it is difficult to obtain the high-precision phase distribution with less time because of the pre-filtering or the DC term subtraction.Hence, to save time and increase the accuracy, the research of two-step PSA without pre-filtering is essential.
In this paper, we will discuss the fast and accurate two-step phase retrieval approach with unknown phase shift.Section 2 presents the principle and process of the proposed PSA based on Gram-Schmidt orthonormalization and Lissajous ellipse fitting method (GS&LEF).In Section 3 the simulation of GS&LEF is discussed, and the comparison with GS is performed.Section 4 evaluates the novel algorithm with the experimental data.The conclusion is finally drawn in Section 5.

Principles
In PSI, the intensity distribution of the phase shifted interferograms can be expressed as: For orthonormalizing two vectors { 1 u , 2 u }, there are three simple steps according to the GS method [26].First, we take one of the vectors and normalize it: Then, we orthogonalize u 2 with respect to the 1 u  vector, subtracting its projection as At last, we obtain 2 u  by normalizing 2 û 2 2 2 ˆˆ.
Note that, in the above equations, ⋅ and , ⋅ ⋅ respectively represent the 2-norm and the inner product.
According to the GS method outlined below, we orthonormalizing 1 I and 2 I as: ( ) where 1   1 . Then, following Eq.( 3), we can obtain 2 Î as: N and y N corresponds to the image columns and rows, respectively.Finally, we obtain 2 I  by Eqs. ( 4) and ( 6) ( ) ( ) where sin cos 1 11) and ( 12) can be rewritten as ( ) ( ) According to Eq. ( 13), a general ellipse function can be obtained ( ) A general conic function can be also expressed by the following second order polynomial: F ax bxy cy dx fy g For an ellipse, Eq. ( 15) needs to meet the conditions of 0 F = and 2 4 0 b ac − < .In the following, we will use the Lissajous ellipse fitting (LEF) method to extract the real phase.At first, we plot an approximate ellipse with N * as the x coordinate and D * as the y coordinate, secondly, the coefficients of the ellipse function can be obtained by the least squares fitting, thirdly, the semi-major amplitude x a , semi-minor amplitude y a , the center offset 0 x , 0 y and the ellipse orientation angle θ with respect to the x-axis can be calculated by Fourthly, we transform the ellipse to a perfect circle centered at the origin using Eq. ( 17).
where c N and c D are the numerator and denominator of ( ) tan ϕ after using the LEF method, the transformation matrix is where θ θ ′ = − , x y r a a = .
After the transformation, the real phase can be finally calculated by Eq. ( 19)

Simulation
To express the performance of the proposed method GS&LEF, we perform several simulations, moreover, we compare the proposed method with the well-evaluated two-step method GS.In the following, we assume that the phase distribution ( ) ( ) where 1 1 , 1 1 x y − ≤ ≤ − ≤ ≤ , the relative phase shift between the first and second interferograms is 1 rad, the image size is 401*401.Figure 1(a) shows the theoretical phase distribution (PV = 31.416rad, RMS = 6.656 rad).
In the first simulation, the background intensity and modulation amplitude distributions of the two interferograms are uniform, 1 1, 0.9, 0.9, 0.8 In the third simulation, the background intensity and modulation amplitude are non-uniform, , , e x p 0 .0 2 a x y a x y x y In the fourth simulation, both the fluctuation and non-uniformity of the background intensity and modulation amplitude exist, such as The fifth simulation is similar to the first simulation with 1 The results of six different simulations are shown in Fig. 2, N represents the index of the simulation, the first and third rows show the phase distributions calculated by GS and GS&LEF, and the phase error distributions are displayed in the second and fourth 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 six simulations, that is to say, the two methods are both robust for different situations.However, the phase error distributions are different for the different simulations and methods.
From the first simulation, we can see that the phase error is approximately equal to zero for GS&LEF when the background intensity and modulation amplitude are perfect, but the phase error also exists for GS (RMS = 0.1610 rad) because of the filtering error and approximation.In the second simulation, the RMS phase error for GS is also very large (RMS = 0.1667 rad), and the phase error for GS&LEF is a very small value with a RMS value of 2.4492e-6 rad, we can get the conclusion that GS&LEF can suppress the effect of background intensity and modulation amplitude fluctuation.The third simulation is a little complex because the background intensity and modulation amplitude distributions are non-uniform for each pixel.The non-uniformity causes some errors for both methods, but the phase error of GS&LEF (RMS = 0.0090 rad) is far smaller than that of GS (RMS = 0.1596 rad).The fourth simulation is more complex than the third simulation, both the fluctuation and non-uniformity of the background intensity and modulation amplitude exist, but the RMS phase errors for third and fourth simulations are similar, the RMS phase errors are 0.1604 rad and 0.0081 rad for GS and GS&LEF respectively.From the first to fourth simulations, we can know that, for GS, the RMS phase error is similar because the filtering error is larger than other errors, such as the non-uniformity between different pixels and fluctuation error between different interferograms, and for GS&LEF, the effect of the non-uniformity is larger than the fluctuation since most of the fluctuation error between different interferograms can be suppressed by the LEF process.
From the fifth simulation, we begin to discuss the effect of noise, we only add the noise to the interferograms, and the background intensity and modulation amplitude are perfect, but we found that the phase error is increasing for both methods, especially for GS&LEF, the RMS phase error which is 0.1281 rad is more than 10 times of the third simulation, that is to say, LEF process only can correct the piston, it cannot correct the noise, and the RMS phase error of GS (0.1701 rad) is a little bigger than that of the third simulation, but the effect of the filtering error is also largest.The sixth simulation 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 for both methods, but the RMS phase error of GS (0.1883 rad) is also larger than that of GS&LEF (0.1449 rad).
From Table 1, we can see that the accuracy of GS&LEF is higher than that of GS for each simulation since there is no pre-filtering in GS&LEF, and the LEF process can extract the accurate phase distribution.In addition, the computational time of GS&LEF is less than that of GS since the pre-filtering spends more time.
Figure 3 shows the ellipses before and after using LEF method for GS&LEF in different simulations, we can see that, after using the LEF method, the ellipse not centered at the origin was transformed to an approximate circle centered at the origin as the above theoretical description.For the front four simulations, whether the ellipse before using the LEF method or the circle after using the LEF method, the curve is smooth, but the curve is not smooth for the last two simulations since the noise exists, and the LEF method cannot remove this effect, hence, it would be best to suppress the noise before using the proposed method to further increase the accuracy.We know that the phase shift value and the level of noise are important to the PSAs, hence, it is necessary to discuss the phase error due to different phase shifts and levels of noises for the proposed method and compare it with GS.We uniformly set the distribution of the phase shift in the range of 0.4 rad to 2.6 rad, and we study the effect of different levels of noises for every phase shift value, other simulated conditions are same as the fourth simulation discussed above, both the fluctuation and non-uniformity of the background intensity and modulation amplitude exist, the RMS phase error curves of GS and GS&LEF for the different phase shifts with the different levels of noises are plotted in Fig. 4, from Fig. 4, we can conclude that: 1) for GS&LEF, the RMS phase error is decreasing with the decrease of the noise, but for GS, there is no relationship between the RMS phase error and the level of noise because the pre-filtering error is the main error; 2) for GS-LEF, the RMS phase errors are relatively stable and small when the SNR of noise is more than 50dB, they are less than 0.03 rad, hence we can ignore it, but for GE&LEF, the phase error can't be ignored when the SNR of noise is more than 50dB; 3) when the phase shift is equal to 0.4 rad and the SNR of noise is 20dB, for GS&LEF, the phase error is relatively large, that is to say, the GS&LEF is invalid in this situation since the phase shift value is too small, the phase error of GS is also large in this situation, but it is valid; 4) when the phase shift is between 0.6 rad and 2.6 rad, both GS and GS&LEF are valid, the RMS phase errors for GS&LEF and GS are respectively less than 0.25 rad and 0.3 rad; 5) when the phase shift is between 1.0 rad and 2.0 rad, the RMS phase error of GS&LEF is less than that of GS, for other range, sometimes the accuracy of GS is higher, hence the best range of phase shift for GS&LEF is between 1 rad and 2 rad; 6) of course, like other two-step PSAs, if the phase shift is close to 0 or π rad, a large phase error also exists, so the phase shifts of both methods should be deviated from 0 or π rad to achieve high accuracy.Finally, we study the effects induced by different phase distributions to further verify the effectiveness of GS&LEF and compare it with GS.We respectively simulate the plane wavefront with 4 x ϕ π = and the complex wavefront ( ) , where we use the function peaks in Matlab, other conditions are same as the sixth simulation.For the plane wavefront, the fringes are straight, as shown in Fig. 5(a), and the theoretical phase distribution is shown in Fig. 5(b).Figure 6 represents the simulated results of plane wavefront using GS and GS&LEF.For the complex wavefront, the fringes are the asymmetrical complex fringes, as shown in Fig. 7(a), the complex phase distribution is drawn in Fig. 7(b), and the simulated results are shown in Fig. 8.According to the above simulations, the calculated results are as follows:1) for the plane wavefront, the RMS phase errors of GS and GS&LEF are respectively 0.1798 rad and 0.1561 rad, and the computational time of GS and GS&LEF are 4.85 s and 3.49 s; 2) for the complex wavefront, the RMS phase error of GS and GS&LEF are respectively 0.2074 rad and 0.1545 rad, and the computational time of GS and GS&LEF are 4.93 s and 3.52 s.The results indicate that the proposed method is suitable for circular, straight and complex fringes, and the accuracy of GS&LEF is higher than that of GS for all the phase distributions, furthermore, for all the phase distributions, GS&LEF costs less time than GS.Based on the above simulations, the conclusions of the proposed GS&LEF can be summarized as: 1) It can obtain the accurate tested phase distribution by only two phase shifted interferograms without the pre-filtering; 2) it can save the computational time; 3) the phase shift can be random, when the large noise exist, the best range of phase shift is between 1.0 rad and 2.0 rad; 4) whether the circular, straight or complex fringes, the proposed method is valid, and the accuracy of GS&LEF is higher than that of GS.

Demonstration with experimental data
To verify the feasibility and outstanding performance of the proposed approach, we perform the experiments with the different fringes.Firstly, four phase shifted interferograms with circular fringes are collected to perform the phase retrieval by GS and GS&LEF, the size of the interferograms is 301*301, and the phase shifts are 0, π/2, π and 3π/2 respectively.The phase extracted by standard 4-step PSA is set as the reference phase due to its high accuracy.One of the interoferograms is shown in Fig. 9(a), Fig. 9(b) shows the reference phase distribution, and the phase distributions extracted by GS and GS&LEF are drawn in Figs.9(c) and 9(d).The differences between the reference phase and the phase extracted by GS and GS&LEF are shown in Figs.9(e) and 9(f), the RMS values of the differences are respectively 0.1453 rad and 0.0817 rad, further verifying the accuracy of GS&LEF is higher than that of GS.Moreover, the computational time of GS&LEF (3.52 s) is less than that of GS (4.37 s).For GS&LEF, the ellipses before and after using LEF method are plotted in Figs.9(g) and 9(h).Then, to verify the effectiveness of GS&LEF, the phase shifted interferograms with the straight and relatively complex fringes which are asymmetrical are also collected, and the comparison for GS and GS&LEF are performed as the circular fringes.The size of the interferograms with the straight fringes is 401*401, and the size of the interferograms with the complex fringes is 301*301, other conditions are same as the circular fringes.Figures 10  and 11 show the results of the straight and complex fringes, we can see that, both GS and GS&LEF 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 GS and GS&LEF are 0.3317 rad and 0.3110 rad, and the computational time for GS and GS&LEF are 5.67 s and 4.48 s respectively.And, for the complex fringes, the RMS values of the differences for GS and GS&LEF are 0.1069 rad and 0.1054 rad, and the computational time are 4.27 s and 3.64 s.For these two kinds of fringes, the accuracy of GS&LEF is a little higher than that of GS, the differences between the reference phase and the phase extracted by GS and GS&LEF are similar since there are large noises in the interferograms.For the computational time, we get the conclusion the same as the circular fringes.Through the above experiments, we verify that, for both the circular, straight and complex fringes, the proposed GS&LEF without pre-filtering can obtain relatively accurate result with less computational time by only two interferograms.(e) and (f) the differences between the reference and extracted phase distributions by GS and GS&LEF, (g) and (h) the ellipses before and after using LEF method.

Conclusion
In this paper, we have proposed a two-step phase retrieval approach based on Gram-Schmidt orthonormalization and Lissajous ellipse fitting method with random phase shift, the GS method is performed firstly, then the tested phase distribution are obtained by the LEF process.We have compared this proposed method with GS by the simulated data and experimental data.The proposed method can directly obtain the tested phase with no prefiltering and less computational time, and it can achieve high accuracy with different background intensity, modulation amplitude distributions and noises.In addition, the best range of phase shift which is between 1.0 rad and 2.0 rad is given.Finally, the proposed algorithm is effective for the circular, straight or complex fringes.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 arbitrary phase shift.
m = 1,2,...,M represents the image index with M the total number of phase shifted interferograms, M is set to 2. ( ) background intensity, the modulation amplitude and the measured phase.m δ represents the phase shift between interferograms.For simplicity, we define 1 the following discussion.
the two phase shifted interferograms with circular fringes are shown in Figs. 1 (b) and 1 (c).In the second simulation, there are fluctuations in the background intensity and modulation amplitude distributions between different interferograms, we set 1 2 1 2 only add 20dB Gauss noise generated by the awgn function of Matlab to the two interferograms.In the last simulation, we add noise with SNR of 20dB to the fourth simulation.

Fig. 1 .
Fig. 1.Simulated phase distribution and two phase shifted interferograms.(a) The theoretical phase distribution, (b) and (c) the first interferogram and the second interferogram.

Fig. 2 .
Fig. 2. The calculated phase distributions and phase error distributions by GS and GS&LEF in different simulations.

Fig. 3 .Fig. 4 .
Fig.3.The ellipses before and after using LEF method for GS&LEF in different simulations.

Fig. 6 .
Fig. 6.Simulated results of straight fringes using GS and GS&LEF.(a) and (b) the extracted phase distributions using GS and GS&LEF, (c) and (d) the phase error distributions after using GS and GS&LEF.

Fig. 8 .
Fig. 8. Simulated results of complex fringes using GS and GS&LEF.(a) and (b) the extracted phase distributions using GS and GS&LEF, (c) and (d) the phase error distributions after using GS and GS&LEF.

Fig. 9 .
Fig. 9. Experimental results of the circular fringes.(a) One of the phase shifted interferograms, (b) the reference phase distribution obtained by 4-step PSA (PV = 38.4319rad, RMS = 7.9988 rad), (c) and (d) the extracted phase distributions obtained by GS (PV = 37.7590 rad, RMS = 7.9190 rad) and GS&LEF (PV = 38.2341rad, RMS = 7.9983 rad), (e) and (f) the differences between the reference and extracted phase distributions by GS and GS&LEF, (g) and (h) the ellipses before and after using LEF method.

Fig. 10 .
Fig. 10.Experimental results of the straight fringes.(a) One of the phase shifted interferograms, (b) the reference phase distribution obtained by 4-step PSA (PV = 21.6492rad, RMS = 6.0171 rad), (c) and (d) the extracted phase distributions obtained by GS (PV = 24.0863rad, RMS = 6.0483 rad) and GS&LEF (PV = 22.1756 rad, RMS = 6.0483 rad), (e) and(f) the differences between the reference and extracted phase distributions by GS and GS&LEF, (g) and (h) the ellipses before and after using LEF method.

Fig. 11 .
Fig. 11.Experimental results of the complex fringes.(a) One of the phase shifted interferograms, (b) the reference phase distribution obtained by 4-step PSA (PV = 135.5400rad, RMS = 30.5750rad), (c) and (d) the extracted phase distributions obtained by GS (PV = 135.1772rad, RMS = 30.5159rad) and GS&LEF (PV = 135.9246rad, RMS = 30.5455rad),(e) and (f) the differences between the reference and extracted phase distributions by GS and GS&LEF, (g) and (h) the ellipses before and after using LEF method.