Accurate phase retrieval algorithm based on linear correlation in self-calibration phase-shifting interferometry with blind phase shifts

Journal Pre-proof Accurate phase retrieval algorithm based on linear correlation in self-calibration phase-shifting interferometry with blind phase shifts Yao Li, Yihui Zhang, Yongying Yang, Chen Wang, Yuankai Chen, Jian Bai

PII: DOI: Reference:

S0030-4018(20)30192-9 https://doi.org/10.1016/j.optcom.2020.125612 OPTICS 125612

To appear in:

Optics Communications

Received date : 12 December 2019 Revised date : 2 February 2020 Accepted date : 23 February 2020 Please cite this article as: Y. Li, Y. Zhang, Y. Yang et al., Accurate phase retrieval algorithm based on linear correlation in self-calibration phase-shifting interferometry with blind phase shifts, Optics Communications (2020), doi: https://doi.org/10.1016/j.optcom.2020.125612. This is a PDF file of an article that has undergone enhancements after acceptance, such as the addition of a cover page and metadata, and formatting for readability, but it is not yet the definitive version of record. This version will undergo additional copyediting, typesetting and review before it is published in its final form, but we are providing this version to give early visibility of the article. Please note that, during the production process, errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain. © 2020 Published by Elsevier B.V.

Highlights

Journal Pre-proof Highlights 1. Any prior knowledge doesn’t need to be measured in advance, and the restrictions on background intensity are removed. 2. The proposed method is applicable to the different kinds of interferograms. 3. It has an excellent phase extraction accuracy especially when the phase shift distributes unevenly.

lP repro of

4. It has no complicated mathematical approximation and transformation thus to be implemented

Jou

rna

simply and easily.

Journal Pre-proof *Manuscript Click here to view linked References

lP repro of

Accurate phase retrieval algorithm based on linear correlation in self-calibration phase-shifting interferometry with blind phase shifts

Yao Li a, Yihui Zhang b, Yongying Yang a,*, Chen Wang a, Yuankai Chen a, Jian Bai a a

State Key Laboratory of Modern Optical Instrumentation, College of Optical Science and Engineering, Zhejiang University, Hangzhou 310027, China b School of Mechatronics Engineering, Henan University of Science and Technology, Luoyang 471003, China *Corresponding author: [email protected]

Abstract: The self-calibration algorithm is a powerful phase map reconstruction technique for phaseshifting interferometry with random phase shifts. The existing phase reconstruction algorithms has some restrictions on preconditions including the uniformity of background and the well distributed phase shifts. In this paper, a simple, applicable and accurate phase retrieval algorithm based on linear correlation (LCA) with non- iterative characteristic is proposed to improve this issue. Without pre-filtering, the algorithm searches for optimal linear combination coefficients of the difference intensity map by maximizing the linear correlation, and then phase shifts can be solved by these searched coefficients. Finally, according to the obtained phase shifts, the phase map is readily retrieved by using the matrix form of the linear equations. The whole phase retrieval process has no complicated mathematical transformation, hence it is simply and easily implemented. Good performance of the proposed LCA are verified by reconstructed phase maps from multiple kinds of simulation and experiment interferograms with blind phase shifts, along with the comparisons to that of the popular algorithms. Especially, the proposed algorithm has an excellent accuracy of phase reconstruction when the phase shift distributes unevenly. Keywords: Self-calibration algorithm; Phase map reconstruction; Phase-shifting interferometry; Linear correlation. 1. Introduction

rna

Phase-shifting interferometry (PSI), as a powerful phase extraction technique with high accuracy and good repeatability, has been widely applied in optical testing, microscopic measurement, in situ wavefront measurement for adaptive optics and other fields [1]. In PSI, several algorithms, such as the least-squares algorithm [2, 3], the Carré algorithm [4] and the hariharan algorithm [5], have been developed to reconstruct the measured wavefront from multiple interferograms. These classical phase-shifting algorithms (PSA) require a priori knowledge of phase shifts or demand linear distribution of phase shifts. However, accurate phase shift is difficult to control due to the miscalibration and the nonlinearity of the phase shifter, vibrational error, air turbulence, and the instability of the laser frequency, etc. Comparatively, the selfcalibrating PSAs with blind and arbitrary phase shifts [6-11] are less dependent on the stability and accuracy of the phase shifter and insensitive to environmental conditions. Up to date, several self-calibrating PSAs using different numbers of interferograms have been reported. For multi-interferograms, the advanced iterative algorithm (AIA) [9] and principle component analysis (PCA) algorithm [10, 11] are two of the most popular techniques for accurate phase extraction. The AIA based on the least-squares method (LSM) is judged as the criterion in this field, while it is great timeconsuming because of multi-interferograms acquisition and algorithm convergence. The latter algorithm distinguishes uncorrelated sine and cosine terms from the captured intensity maps, and its non-iterative nature makes it low computational cost. But if the high-accuracy phase reconstruction is required in the PCA method, it is necessary to satisfy the condition that the phase shift is uniformly distributed in the range of [0, 2π] and the fringes in interferograms exceed one. Besides, Hao et al [12] proposed a statistical method. It determines the background and modulation amplitudes pixel by pixel, which is much timeconsuming. Cai et al. [13, 14] exploited a series of methods for phase reconstruction, in which either the amplitude of reference wave was assumed to be constant or the reference amplitude and the object amplitude were recorded prior to the measurement. Zhang et al. [15] proposed a three-step PSI based on the difference map normalization and diagonal vector normalization method on the assumption of more than

Jou

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65

Journal Pre-proof

rna

lP repro of

one fringe in interferograms. With the same assumption, the GS3 method without filtering was demonstrated well for three-frame PSI [16]. Recently much more attentions have been paid for two-frame PSA. Based on Lissajous figure and ellipse fitting technology [17], Liu et al. [18] realized the simultaneous extraction for measured wavefront and phase shifts from two interferograms processed by Hilbert-Huang pre-filtering, while this algorithm needs uniform intensity distribution. Cheng et al. [19] researched a fast analytic method by solving a quartic polynomial equation with an approximation mentioned above, but this technique also uses a highpass Guassian filter to suppress the direct-current (DC) background. In addition, two-step PSA based on Gram-Schmidt orthonormalization method was proposed in Ref. [20] for retrieving random phase distributing in [0, 2π] except π. Tian et al. [21] adopted global optimization by the difference evolution algorithm to search the phase shifts, which is insensitive to the noises and defects and is able to deal with small fringes numbers. Nevertheless, it is time consuming owing to the complicated calculation. Also an algorithm, proposed by the same authors, evaluated the phase shifts in a local mask [22] and suppressed the background using a high-pass Gaussian filter. Similar to the method in Refs. [18, 19], these algorithms still need to subtract DC background term by pre-filtering, hence the residual background would influence the accuracy of the phase extraction. Overall, almost all two-frame PSAs have their own restrictions on preconditions, such as the requirement for the constant of background and the singular phase shift π rad. These rigid requirements narrow the application in optical testing. Compared with two-frame PSAs, multi-frame PSAs, such as AIA and PCA, are capable of removing above mentioned restrictions, meanwhile, for all we known, they can also achieve high-accuracy phase extraction. But in self-calibration PSAs, less computational time is also required simultaneously when maintaining these properties. For this reason, a non-iterative phase retrieval algorithm based on linear correlation (LCA) with blind phase shifts is developed in this paper, which achieves the accurate, efficient and robust phase reconstruction, especially for uneven phase shifts. In the proposed LCA, the phase shifts are calculated based on searching procedure of making a merit function maximum, as derived merit function in Refs. [23-27]. Guo et al. [23] utilized an objective function based on cross power spectrum, and a gradient-guided search strategy is used for minimizing this objective function so as to estimate the phase shifts. Cheng et al. [26] searched the solution of phase shifts to obtain the minimum coefficient of variation of the modulation amplitude. For our merit function of LCA, we search for the linear combination coefficients of difference intensity map, defined as b and c, to maximum the linear correlation, so that the phase shifts can be solved by these searched coefficients. And the calculation procedure of the proposed LCA is simple, fast and without complicated mathematical approximation and transformation, so the accurate and time-saving phase reconstruction can be attained. Additionally, this algorithm is suitable for the phase extraction from interferograms with blind phase shift and without any prior information. And its accuracy for phase reconstruction is insensitive to the distribution of background and modulation amplitudes. In this paper, the proposed LCA will be discussed in principle and analyzed in simulation and experiment. Section 2 presents the basic principle and calculation process of proposed LCA; the numerical simulations for the phase reconstruction are conducted by LCA, along with that of AIA and PCA method for comparisons in several cases in Section3; the experiment results are presented to demonstrate the validity of the proposed algorithm in Section 4; Section 5 presents a conclusion.

Jou

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65

2. Principles

2.1 Basic principles of the LCA The proposed LCA are described in detail below. Assuming that the background amplitude Aj and the modulation amplitude B j are the invariant in each frame, the intensity distribution of the obtained interferograms can be expressed by [28] I n , j  Aj  B j cos( j   n ).

(1)

where the subscript n and j represent the nth phase-shifted interferogram and the pixel position in the obtained interferogram, respectively;  j is the phase distribution that requires to be extracted, usually only a function of the pixel;  n denotes the unknown phase shift between the nth interferogram and the initial interferogram. In the initial interferogram, we assume  0  0 generally.

Journal Pre-proof

In this method, four interferograms are required, assuming they are initial interferogram I 0, j , first phase-shifted interferogram I1, j , followed by phase-shifted interferograms I 2, j and I 3, j in sequence. By

lP repro of

subtracting I n , j (n  1, 2,3) from I 0, j , Aj can be eliminated. Then, we can get the difference intensity map Dn, j  I 0, j  I n, j  2B j sin(

Then we set  j   j 

n 2

)sin( j 

n 2

),(n  1,2,3).

(2)

  1 1  (n  1,2,3) , the Eq. (2) can be expressed by  and  n  n 2 2 2 Dn, j  2B j sin(

n 2

)cos   j   n .

(3)

Note that modulation terms B j sin( n / 2) have amplitude inconsistency because of different phase shifts

 n . Hence to ensure jth pixel has the same modulation amplitude in each Dn , j , the vector Dn , j is simplified with the normalization method as Dn , j

Dˆ n , j =



B j cos( j   n )

m

B

Dn , j , Dn , j

j 1

2

j

.

(4)

cos ( j   n ) 2

where ,  represents the inner product, and m is the total number of pixels in each Dn , j . And if we have m

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

B j 1

Then we have m

B j 1

2 j

2 j

 cos   j   sin   j   0 [29].

 cos 2  j  cos 2   j   n  

  B 2j  cos  j  cos   j   n    cos  j  cos   j   n   m

j 1

rna

m            4   B 2j  cos   j  n   sin   j  n   cos  n   sin  n  2 2    2  2 j 1

(5)

m        2  sin  n    B 2j  cos   j  n   sin   j  n  2 2   j 1 0

We can rewrite Eq. (5) as

m

Jou

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65

Then let Cn 



m

j 1

B j 1

2

j

cos 2 ( j ) 

m

B j 1

2 j

cos 2 ( j   n )

(6)

B j 2 cos2 ( j   n ) , where the phase shift  1 is zero. So simplify Eq. (6) with Cn , we can

get C1 C2  1 and C1 C3  1 . According to Eq. (4), we can see that the normalized difference intensity map is the cosine of ( j   n ) . However what we desired to get are two terms: the sine and cosine of  j , so in view of this case we naturally expand the Eq. (4) given by

Dˆ1, j =Tc , j .

(7)

C Dˆ 2, j = 1 [Tc , j cos( 2 )  Ts, j sin( 2 )]  Tc , j cos( 2 )  Ts, j sin( 2 ). C2

(8)

Journal Pre-proof

C Dˆ 3, j = 1 [Tc , j cos( 3 )  Ts, j sin( 3 )]  Tc , j cos( 3 )  Ts, j sin( 3 ). C3 where Ts , j  B j sin( j ) / C1 and Tc , j  B j cos( j ) / C1 .

lP repro of

(9)

Up to now, although these three equations are derived, both Ts , j term and Tc , j term still cannot be solved in case of unknown phase shifts. However interesting, through the linear combination of Eqs. (7) ~ (9), such as y  b  Dˆ  Dˆ and y  Dˆ , there is a unique constant b, so that y and y can satisfy 1, j

1, j

3, j

2, j

2, j

1, j

2, j

the linear relationship y1, j =  y2, j , where the coefficients of Ts , j and Tc , j are proportional. So we can get cos( 2 )sin( 3 )  cos( 3 ). sin( 2 ) above, we infer two matrices b

According

to

the

analysis

Y1  [ y1,1 , y1,2 , y1,3  y1,m ]

(10) and

Y2  [ y2,1 , y2,2 , y2,3  y2,m ] are linearly related. The correlation coefficient ρ between Y1 and Y2 is given by



cov(Y1 , Y2 ) . std(Y1 )  std(Y2 )

where cov(, ) and std() represent the covariance and standard deviation, respectively, defined as

(11)

m

cov(Y1,Y2 )   ( y1, j  y1 )( y2, j  y2 ).

(12)

j 1

m

( y

std(Y ) 

j 1

j

 y )2 .

(13)

where y is the sample mean of Y . Our task is to search for a constant b that maximizes the correlation coefficient ρ between Y1 and Y2 and is immune to the phase data distribution. Similarly, for another linear combination of y  c  Dˆ  Dˆ and y  Dˆ , we continue to search 1, j

2, j

3, j

2, j

1, j

rna

for the constant c, making maximum correlations between variables y1, j and y 2, j . Then we have c

sin( 3 ) . sin( 2 )

(14)

where  2 , 3 [0, ] . Therefore, according to Eqs. (10) and (14), the unknown phase shift  2 and  3 can be solved, they are  1  c2  b2  .  2bc 

(15)

 c2 1 b    .  2b 2b 2 

(16)

 2  arc cos 

Jou

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65

 3  arc cos 

Hence, after getting the values of  2 and  3 , the terms Ts , j and Tc , j can be solved by using the matrix form of Eqs. (7) ~ (9) as

 Dˆ1, j  0    T   c, j sin  2     Dˆ 2, j   Ts , j   ˆ  sin  3    D3, j  Finally, the measured phase distribution can be reconstructed by  1  cos  2   cos  3

 Ts , j  j  tan 1  T  c, j

  . 

(17)

(18)

Journal Pre-proof

lP repro of

The reconstructed phase  j has a constant deviation (  j   j  (1   ) 2 ) from the true phase  j . It would be fitted to the piston term of the Zernike polynomial. Generally, the Zernike polynomial is applied to fit the phase distribution in the optical testing [1]. As far as we know, the piston aberration would not affect the reconstructed phase distribution [1]. 2.2 Calculation process of the LCA

Different from other methods, the proposed LCA utilizes the maximum linear correlation of the combined intensity maps to obtain a function of the phase shifts. For clarity, as shown in Fig. 1(a), we summarize the entire calculation process of the proposed LCA in detail: Step 1. Obtain four frames blind phase-shifted interferograms, and then reshape each frame into one column with size of m. Step 2. Subtract I n, j  n  1,2,3 from I 0, j to get the difference intensity vector Dn , j with no background intensity. And then normalize the vector Dn , j .

Step 3.1. Firstly, determine two sets of linear combinations of the difference intensity vector Dn , j . One of them is y1, j  b  Dˆ1, j  Dˆ 3, j and y2, j  Dˆ 2, j , then search for a constant b to maximize |ρ|. The search algorithm is described in detail in Fig. 1(b).

Step 3.2. And the same process is conducted again for the other set y1, j  c  Dˆ 2, j  Dˆ 3, j and y2, j  Dˆ1, j . Step 4. Utilize b and c determined in the previous step to solve the phase shifts  2 and  3 in Eqs. (15) and (16). Step 5. Substitute the solved phase shifts  2 and  3 into Eqs. (7) ~ (9), then calculate the sine term Ts , j and cosine term Tc , j by using the matrix form of Eqs. (7) ~ (9).

rna

Step 6. Finally extract wrapped phase map  j with the Eq. (18).

Jou

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65

Fig. 1. (a) Flow chart of the entire calculation process of the proposed LCA. (b) A detailed flow chart of step 3.1 with respect to searching for a constant b to satisfy the maximum correlation condition.

As to the search process in Fig. 1(a), we take the step 3.1 for an example to describe in Fig. 1(b) in great detail. Set a set of initial input variable b [bmin , bmax ] , and evenly sample n points at interval  . Then we calculate the correlation coefficient ρ of variable b by using Eq. (10). As shown in Fig. 2(a), the distribution k of the correlation coefficients has an extreme value  max , which corresponds an input variable value bpk , and under this condition, the correlation matrices Y1 and Y2 are optimally linearly related. Finally judge k 1 k   max |  (  is the predefined minimum) is satisfied. If so, the best whether the merit function |  max

optimal value b is determined to solve the phase shift. If not, set bmin  bpk   and bmax  bpk   , and then repeat the search process until the merit function is fulfilled.

Fig. 2. The relationship between the correlation coefficients ρ and (a) the input variable b and (b) the input variable c.

According to the procedures in Fig. 1(b), four blind phase-shifted interferograms were simulated with different signal-to-noise ratios (SNRs). In the simulation, we searched for the optimal coefficients (that is b and c) and the maximum correlation coefficient ρ. Fig. 2 shows that the correlation coefficient ρ varies with the input variable under different Gauss noise conditions. From Fig. 2, it can be inferred that there is always an extreme point of the correlation coefficient no matter how much noise is added into the interferograms. The searched extreme points are recorded in Table 1. As can be seen from Table 1, in the case of in Fig. 2(a) and (b), the extreme value ρ of the linear correlation coefficients are all above 0.99, and the changes in the optimum values of the corresponding b and c are negligible. In Fig. 2(a), although the extreme value of ρ is 0.9988 when noise with an SNR of 35dB is added, the corresponding optimum b is 9.6845, which differs by only 0.8% from that of an SNR of 60dB. In summary, as a key procedure of the proposed algorithm, the searched coefficients (that is b and c) which satisfy the linear correlation are insensitive to the noise. Table 1. The comparisons of the extreme value of coefficients ρ and corresponding optimal values of b and c in Fig. 2 when Gauss noises with SNR of 35dB, 45dB and 60dB are added into the interferograms. Figure 2(a)

Figure 2(b)

SNR 35dB

Optimal value b 9.6845

Extreme value ρ 0.9988

Optimal value c -9.4992

Extreme value ρ 0.9988

45dB

9.6043

0.9999

-9.4591

0.9999

60dB

9.6043

1.0000

-9.4591

1.0000

rna

3. Simulation validation

In order to verify the effectiveness, universality and robustness of the proposed LCA method, a series of computational numerical simulations were performed on four interferograms with a pixel size of 256×256. In these interferograms, it is assumed that the non-uniform background amplitude and modulation amplitude are A  1.3  exp[0.02  ( x 2  y 2 )], B  1.2  exp[0.02  ( x 2  y 2 )].

(19) respectively, where 1  x  1 and 1  y  1 . In addition, White Gaussian noise with a SNR of 35 dB is added to these interferograms. In the following simulation, all numerical processing was performed by MATLAB R2017a on a computer (Intel Core i3−8100, 3.6 GHz). In the first circular interferogram simulation, the phase distribution is defined as   5  ( x 2  y 2 ) . To demonstrate the validity of the approximation in Eq. (6), we calculated Cn using the simulation phase data with phase shifts in the case 2 of Table 2. After being calculated, the values of C1, C2 and C3 are respectively 215.290, 216.136 and 215.340, which are approximately equal to each other. The values of the ratio C1/C2 and C1/C3 are approximately equal to 0.996 and 1.000, respectively. Hence, the approximation operation used in proposed method is verified to be effective. Table 2 summarizes the reconstructed phase error root-mean-square (RMS) and processing time of the proposed LCA, along with that of PCA and AIA for comparison in three representative cases. The true phase shift is [0, 0.57, 0.78, 4.15] rad in the case 1, [0, 0.57, 2.5, 2.9] rad in the case 2 and [0, 1.45, 3.14, 4.77] rad in the case 3. The results clearly show that in all cases, the RMS values of AIA are the most accurate, but due to multiple iterations, its processing time is much longer than the other two algorithms. For example, in the case 2, the phase errors RMS extracted by proposed LCA, PCA and AIA are 0.0289 rad, 0.0397 rad and 0.0246 rad, respectively, and their processing time are 0.035s, 0.029s and 36.536s, respectively. In addition, we can see that the PCA method has the

Jou

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65

lP repro of

Journal Pre-proof

Journal Pre-proof

lP repro of

shortest processing time, but its accuracy is worse than that of LCA, especially in the case 1, the proposed LCA has a residual phase error RMS of 0.0240 rad, which is much smaller than that of PCA of 0.0339 rad. It is also worth noting that the PCA algorithm requires extra time to correct the global sign of the measured phase. Therefore, it can be concluded that the proposed LCA approach has an ability of accurate phase reconstruction like AIA, as well as fast processing time. Additionally, when the phase shifts are well distributed between 0 and 2π in the case 3, the PCA method with 0.0136 rad of phase error RMS has a better performance than the LCA with 0.0147 rad of that. However the LCA can have better results in the case of uneven phase shifts distribution, such as the case 1 and the case 2. Table 2. Comparisons of the simulated phase error RMS and processing time between the proposed LCA, PCA and AIA method for four-frame phase-shifted interferograms. Case 1a RMS (rad)

a

LCA

0.0240

PCA

0.0339

AIA

0.0215

Case 2b

Case 3c

Time (s)

RMS (rad)

Time (s)

RMS (rad)

Time (s)

0.049

0.0289

0.035

0.0147

0.039

0.033

0.0397

0.029

0.0136

0.032

9.30

0.0246

36.536

0.0106

19.285

The true phase shift in the case 1 is [0 0.57 0.78 4.15] rad. bThe true phase shift in the case 2 is [0 0.57 2.5 2.9] rad. cThe true phase

shift in the case 3 is [0 1.45 3.14 4.77] rad.

rna

Fig. 3(a)-(c) show three unwrapped phase maps retrieved by the proposed LCA, PCA and AIA method, respectively, when phase shift is the case 1 in Table 2. It can be seen that the three reconstructed phase maps are very similar. To better compare the performance of the these three algorithms, we subtracted the true phase map from the reconstructed phase maps, and the corresponding residual phase error maps are shown in Fig. 3(d)-(f). We can see that there is a significant distinction between the phase error maps. Meanwhile, it can be seen from the error histogram and quantitative comparison in Fig. 3(g)-(i) that the number of pixels in the histograms is mainly concentrated in the near zero phase error region from the LCA and the AIA. In addition, the phase reconstruction error RMS of the proposed LCA is 0.0240 rad, which is only 0.0025 rad higher than AIA, while, the RMS difference between PCA and AIA is 0.0124 rad. Obviously, the proposed LCA is superior to the PCA method in the reconstruction accuracy.

Jou

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65

Fig. 3. Comparisons of the phase retrieval simulation with circular interferograms when the phase shift is the case 1 in Table 2. (a), (b) and (c) are unwrapped phase maps (unit: rad) reconstructed by the proposed LCA, PCA and AIA, respectively, whose phase error maps are (d), (e) and (f). (g), (h) and (i) are the error histograms from the proposed LCA, PCA and AIA method respectively.

Journal Pre-proof

lP repro of

In order to examine the applicability and universality of the proposed LCA, simple straight and complex ghost interferograms with a size of 256×256 pixels were simulated. Similar to the circular interferogram, we calculate the value of Cn to examine whether the approximation operations in Eq. (6) are satisfied with different application. In the case of straight interferograms, the values of C1/C2 and C1/C3 are 1.000 and 0.996 respectively. And for complex ghost interferograms, the ratio values are 1.001 and 1.011, which are approximately equal to one. And the reconstruction results in two cases are shown in Fig. 4. and Fig. 5, respectively. Similar to the circular interferogram, the proposed LCA method still performs well for both straight interferograms and ghost interferograms. For these three different forms of interferograms, the residual phase error maps reconstructed by the LCA, as shown in Fig. 3(d), Fig. 4(a) and Fig. 5(d), have the RMS of 0.024 rad, 0.0257 rad and 0.0293 rad, respectively, and they are very close to each other. Also the retrieved phase maps of the proposed LCA are more accurate than that of PCA for all kinds of interferograms, and the LCA takes less computation time than the AIA. It is proved that the proposed method is accurate, applicable and universal.

rna

Fig. 4. (a), (b) and (c) are the reconstructed residual phase error maps by the proposed LCA, PCA and AIA method, respectively, along with the error histograms (d) (e) and (f).

Jou

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65

Fig. 5. (a), (b) and (c) are phase maps reconstructed by the proposed LCA, PCA and AIA using simulated complex ghost interferograms, (d), (e) and (f) are the corresponding residual phase errors distribution relative to the true phase maps. All units of the phase are radians. (g), (h) and (i) are the error histograms respectively.

Journal Pre-proof

lP repro of

In principle, the proposed LCA has a limit on the fringe number in the interferograms. To quantify and clarify this, we processed a series of simulated interferograms with different fringe numbers of 2kπ rad, as shown in Fig. 6. As can be seen from Fig. 6(a1) and (a2), when the fringe number is less than 0.8, the LCA is inferior in accuracy to the other two methods due to the residual error from the approximation in Eqs. (5) and (6). However, as the fringe number is further increased, the performance of the proposed LCA becomes stable and better, while the phase error of the PCA method is unstable and larger than the proposed method. Furthermore, when the fringe numbers are more than 3, the proposed LCA has very similar accuracy to the phase reconstruction of the AIA. Therefore, when the fringe number exceeds 0.8, the proposed LCA can achieve good accuracy, but if the measurement demands higher accuracy, it is necessary to tilt the reference mirror in the experiment to adjust fringe numbers to 3 or more. In addition, as shown in Fig. 6(b1) and (b2), circle, straight and complex ghost interferograms were also simulated to verify the universality of the proposed method. When dealing with fringes less than 0.8, it can be seen that LCA has significant better performance on the straight interferograms than the other two forms. And when the fringe numbers are more than 1.2, the proposed LCA has the similar accuracy for these three forms of interferograms.

rna

Fig. 6. Performance comparisons of processing interferograms with different fringe numbers. (a1) and (a2) are the simulation results for the circular fringes by the proposed LCA, PCA and AIA. (a2) is the enlarged view of (a1) in which the fringe number is more than 0.7. (b1) and (b2) are the simulation results for different forms of interferograms by the proposed LCA. (b2) is the enlarged view of (b1) in which the fringe number is more than 0.8.

Finally, we further analyze the robustness of the proposed LCA as the SNR of the interferogram increases from 25dB to 70dB, comparing with the robustness of the PCA and the AIA, as depicted in Fig. 7. And Table 3 gives the detailed results for two cases, one for high-contrast circular interferograms and the other for low-contrast circular interferograms. It is concluded that the three methods have the same sensitivity to noise. Besides the proposed LCA has a better accuracy than the PCA method in the case of larger noise. It is notable that when the contrast of the interferogram is poor, the residual error RMS of three methods will increase over the whole range of noise. And when the noise is almost zero, the proposed LCA has a residual phase error of about 0.008 rad due to the assumptions in Eqs. (5) and (6). However, in actual experiment, there is always some inevitable noise generated by the detector, making the residual error negligible.

Jou

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65

Fig. 7. The phase error RMS simulated by the three algorithms in the case of (a) high-contrast interferograms and (b) low-contrast interferograms with different levels of noise.

Journal Pre-proof

Table 3. The phase error RMS simulated by the three algorithms in the case of high-contrast interferograms and low-contrast interferograms with different levels of noise. SNR in high-contrast case 25dB

30 dB

35 dB

LCA

0.0872

0.0452

0.0239

PCA

0.0910

0.0565

0.0339

AIA

0.0679

0.0383

0.0215

RMS

40 dB

45 dB

50 dB

55 dB

60 dB

65 dB

70 dB

No noise

0.0158

0.0111

0.0092

0.0085

0.0083

0.0082

0.0082

0.0081

0.0188

0.0121

0.0067

0.0042

0.0027

0.0021

0.0018

0.0017

0.0121

0.0068

0.0038

0.0022

0.0014

0.0009

0.0007

0.0005

lP repro of

RMS

SNR in low-contrast case

25dB

30 dB

35 dB

LCA

0.4886

0.2577

0.1539

PCA

0.4673

0.2809

0.1756

AIA

0.4730

0.2408

0.1298

4. Experiment validation

40 dB

45 dB

50 dB

55 dB

60 dB

65 dB

70 dB

No noise

0.0941

0.0485

0.0272

0.0165

0.0115

0.0093

0.0085

0.0082

0.1041

0.0580

0.0343

0.0175

0.0113

0.0062

0.0043

0.0017

0.0721

0.0411

0.0232

0.0134

0.0082

0.0056

0.0047

0.0039

rna

Verification experiments are carried out by the proposed LCA, along with AIA and PCA for comparisons. The experiment setup is the same as in the point diffraction interferometry system of Ref. [30] in our previous research, which uses a nearly ideal spherical wave generated by a pinhole as reference wave. A PZT is applied to generate phase shifts of the spherical mirror under test. Besides, by adjusting a five-axis lens holder of the test mirror, we can capture straight, circular and complex interferograms on the CCD. In our first experiment, multiple straight interferograms were used to extract the reference phase map. Fig. 8(a) and (b) show the peak-to-valley (PV) and RMS convergence of the reconstructed phase maps by the AIA method. The results indicate that phase error has a satisfactory convergence when processed interferograms are no less than 19 frames. Therefore, as shown in Fig. 8(c), the reference phase map is reconstructed by the AIA method with 19 frames interferograms, and its PV value and RMS value are 28.518 rad and 7.102 rad, respectively. After determining the reference phase, we utilize these three algorithms to reconstruct the measured phase from four interferograms. From Fig. 9(b) (e) and (h), we can see that all the retrieved phase maps are similar in an RMS value of 7.071 rad. Therefore, for better comparison, the residual phase error maps relative to the reference phase map in Fig. 8(c) are calculated, as shown in Fig. 9(c) (f) and (i). As can be seen from Fig. 9(c) (f) and (i) , the reconstruction phase error RMS value by the proposed LCA is 0.0296 rad, which is close to that of 0.0303 rad by the AIA method, and better than that of 0.0433 rad by the PCA. The corresponding processing time of phase reconstruction is 28.49s for AIA, 0.219s for PCA and 0.274s for the proposed LCA. Besides, the phase shift  n calculated by the proposed LCA method is [0, 0.110, 0.949] rad, which is very close to the phase shift [0, 0.126, 0.971] rad from AIA method. These statistical data verify the consistency of the accuracy between the proposed LCA and the AIA, meanwhile the proposed LCA takes much less time than the AIA. Therefore, it can be concluded that our proposed method has good phase map retrieval performance.

Jou

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65

Fig. 8. (a) PV and (b) RMS convergence of the extracted phase map by the AIA method. (c) The reference unwrapped phase map (unit: rad) reconstructed from 19 interferograms by the AIA method.

Fig. 9. The result comparisons between these three methods based on four experimental interferograms in the left figure. (a), (d) and (g) are the wrapped phase maps reconstructed by AIA, PCA and proposed LCA method, respectively. (b), (e) and (h) are the corresponding unwrapped phase map. For better comparison, (c), (f) and (i) are the reconstructed phase error maps relative to the reference phase map, respectively.

rna

In addition, the circular interferograms were also obtained in the experiment. In the case of circular interferograms, the results calculated by these three methods are shown in Fig. 10. The residual error RMS from AIA, PCA and proposed LCA are respectively 0.0585 rad, 0.0714 rad and 0.0617 rad. In addition, the phase shift  n calculated by the proposed LCA and AIA method are [0, 0.083, 1.431] rad and [0, 0.085, 1.460] rad, respectively, which are very close to each other. It is concluded that the proposed LCA is also applicable to the circular interferograms, and compared with PCA, it has phase reconstruction accuracy closer to AIA. In the case of complex interferograms, the reconstructed wrapped phase map by AIA and proposed method are shown in Fig. 11(a1), respectively. To examine their consistency, we unwrapped the reconstructed phase maps shown in Fig. 11 (b1). It can been seen from Fig. 11 (b2) that even if the edge of the interferogram has a low contrast, the phase map still can be restored well by the proposed LCA method. The PV value and RMS value of the unwrapped phase map reconstructed by proposed LCA method are 95.2862 rad and 17.5375 rad, respectively, which are extremely close to that of phase map (with PV value of 95.2324 rad and RMS value of 17.5378 rad) from AIA method. Fig. 11 (b3) shows the difference phase map between the reconstructed phase by LCA and AIA method, whose RMS value is only 0.0389 rad. Besides, the phase shift  n calculated by the proposed LCA is [0, 0.404, 1.433] rad, which is similar to the phase shift [0, 0.442, 1.464] rad from AIA method. These statistics verify the consistency between the phase maps from the proposed LCA and AIA methods for both circular and complex interferograms.

Jou

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65

lP repro of

Journal Pre-proof

rna

Fig. 10. Circular interferograms experimental result comparisons between these three methods. (a) The reference phase map reconstructed by the AIA method from 19 interferograms. (b) Experimental interferograms for phase retrieval. (c), (d) and (e) are the unwrapped phase maps reconstructed by AIA, PCA and proposed LCA method, respectively. For better comparison, (f), (g) and (h) are the corresponding residual error maps relative to the reference phase map.

Fig. 11. The experimental result with the complex interferograms in the top left corner figure (size: 623×623). (a1) and (b1) are the reconstructed wrapped phase maps from AIA and proposed LCA method, respectively. For better comparison, the corresponding reconstructed unwrapped phase maps are shown in (a2) and (b2), respectively. (b3) is the difference map between the phase map reconstructed by AIA and proposed LCA method.

To demonstrate the performance of the proposed LCA method more intuitively, the residual RMS errors and processing times are summarized in Table 4, which are similar to the simulation results. Consequently, it is verified that the proposed LCA method have acceptable accuracy and computation time for processing different shapes of interferograms.

Jou

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65

lP repro of

Journal Pre-proof

Table 4. Comparisons of the residual phase error RMS and processing time between AIA, PCA and LCA method. Straight interferograms

Circular interferograms

Complex interferograms

RMS (rad)

Time (s)

RMS (rad)

Time (s)

RMS (rad)

Time (s) 30.56

AIA

0.0303

28.49

0.0585

29.184

-

PCA

0.0433

0.219

0.0714

0.267

-

-

LCA

0.0296

0.274

0.0617

0.314

0.0389

0.382

5. Conclusions

In conclusion, we have proposed a non-iterative phase retrieval algorithm based on linear correlation (LCA) in self-calibration interferometry with blind phase shifts, which is simple, accurate and efficient. This method can calculate the unknown phase shifts and reconstruct the phase map simultaneously by searching for the optimal linear combination coefficients b and c to satisfy the maximum linear correlation condition. During the whole phase reconstruction process, any prior knowledge doesn’t need to be

Journal Pre-proof

lP repro of

measured in advance and the restrictions on the background intensity are removed. Also it doesn’t demand that the phase shift is well distributed in the range of [0, 2π] like the PCA method. We have tested the proposed method through simulated and experiment interferograms, and compared its performance with the AIA and PCA methods. It is demonstrated that the proposed method can achieve high accuracy and take less time simultaneously, especially in the case of uneven phase shift distribution. Although it still requires more than one fringes, it can be easily implemented in the experiment. Acknowledgment

This work is supported by the National Natural Science Foundation of China (61627825, 61875173). References

3. 4. 5. 6. 7. 8. 9. 10. 11. 12. 13. 14. 15. 16. 17. 18.

19. 20. 21. 22. 23. 24. 25. 26. 27. 28.

D. Malacara, Optical Shop Testing, 3rd ed. (John Wiley & Sons, Inc., 2007). 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) (1974) 2693-2703. J. E. Greivenkamp, Generalized Data Reduction for Heterodyne Interferometry, Opt. Eng. 23(4) (1984) 234-350. P. Carré. Installation et utilisation du comparateur photoélectrique et interférentiel du Bureau International des Poids et Mesures, Metrologia 2(1) (1966) 13-23. P. Hariharan, B. F. Oreb, and T. Eiju, Digital phase-shifting interferometry: a simple error-compensating phase calculation algorithm, Appl. Opt. 26(13) (1987) 2504-2506. 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) (2011) 20483-20492. K. Okada, A. Sato, and J. Tsujiuchi. Simultaneous calculation of phase distribution and scanning phase shift in phase shifting interferometry, Opt. Com. 84(3-4) (1991) 118-124. G. Han and S. Kim, Numerical correction of reference phases in phase-shifting interferometry by iterative least-squares fitting, Appl. Opt. 33(31) (1994) 7321-7325. Z. Wang and B. Han, Advanced iterative algorithm for phase extraction of randomly phase-shifted interferograms, Opt. Lett. 29(14) (2004) 1671-1673. J. Vargas, J. A. Quiroga, and T. Belenguer, Analysis of the principal component algorithm in phase-shifting interferometry, Opt. Lett. 36(12) (2011) 2215-2217. J. Vargas, J. A. Quiroga, and T. Belenguer, Phase-shifting interferometry based on principal component analysis, Opt. Lett. 36(8) (2011) 1326-1328. Q. Hao, Q. Zhu, and Y. Hu, Random phase-shifting interferometry without accurately controlling or calibrating the phase shifts, Opt. Lett. 34(8) (2009) 1288-1290. 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) (2003) 1808-1810. 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) (2004) 183-185. Y. Zhang, X. Tian, and R. Liang, Three-step random phase retrieval approach based on difference map normalization and diamond diagonal vector normalization, Opt. Express 26(22) (2018) 29170-29182. 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) (2014) 5147–5154. C. T. Farrell and M. A. Player, Phase step measurement and variable step algorithms in phase-shifting interferometry, Meas. Sci. Technol. 3(10) (1992) 953–958. 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. Opt. 18(10) (2016) 105604. Z. Cheng and D. Liu, Fast and accurate wavefront reconstruction in two-frame phase-shifting interferometry with unknown phase step, Opt. Lett. 43(13) (2018) 3033-3036. 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) (2012) 443–445. C. Tian and S. Liu, Demodulation of two-shot fringe patterns with random phase shifts by use of orthogonal polynomials and global optimization, Opt. Express 24(4) (2016) 3202-3215. C. Tian and S. Liu, Phase retrieval in two-shot phase-shifting interferometry based on phase shift estimation in a local mask, Opt. Express 25(18) (2017) 21673-21683. H. Guo, Y. Yu, and M. Chen, Blind phase shift estimation in phase-shifting interferometry, J. Opt. Soc. Am. A 24(1) (2007) 2533. P. Gao, B. Yao, N. Lindlein, K. Mantel, I. Harder, and E. Geist, Phase-shift extraction for generalized phase-shifting interferometry, Opt. Lett. 34(22) (2009) 3553-3555. 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) (2006) 1966-1968. Z. Cheng and D. Liu, Universal phase reconstruction approach of self-calibrating phase-shifting interferometry, Opt. Lett. 44(15) (2019) 3857-3860. H. Ding and H. Guo, Estimating phase shifts from three fringe patterns by use of cross spectrum, Appl. Opt. 56(4) (2017) 916927. M. Born and E. Wolf, Principles of Optics, 7th ed. (Cambridge University, 1999).

rna

1. 2.

Jou

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65

Journal Pre-proof

rna

lP repro of

29. J. Vargas, J Antonio Quiroga, C. O. S. Sorzano, J. C. Estrada, J. M. Carazo. Two-step demodulation based on the Gram-Schmidt orthonormalization method. Opt. Lett. 37(3) (2012) 443-447. 30. Y. Li, Y. Yang, C. Wang, Y. Chen, Y. Zhang, and J. Bai, Comprehensive design and calibration of an even aspheric quarterwave plate for polarization point diffraction interferometry, Appl. Opt. 57(8) (2018) 1789-1799.

Jou

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65

*Credit Author Statement

Journal Pre-proof CRediT authorship contribution statement Yao Li: Writing — original draft, Conceptualization, Methodology, Validation. Yihui Zhang: Validation, Supervision. Yongying Yang: Conceptualization, Methodology, Validation, Supervision.

Jou

rna

lP repro of

Yuankai Chen: Validation. Chen Wang: Validation. Jian Bai: Supervision.

Declaration of Interest Statement

Journal Pre-proof Declaration of interest statement The authors declare that they have no known competing financial interests or personal relationships

Jou

rna

lP repro of

that could have appeared to influence the work reported in this paper.