An Improved Synthesis Phase Unwrapping Method Based on Three-Frequency Heterodyne

An improved three-frequency heterodyne synthesis phase unwrapping method is proposed to improve the measurement accuracy through phase difference and phase sum operations. This method can reduce the effect of noise and increase the equivalent phase frequency. According to the distribution found in the phase difference calculation process, the Otsu segmentation is introduced to judge the phase threshold. The equivalent frequency obtained from the phase sum is more than those of all projected fringe patterns. In addition, the appropriate period combinations are also studied. The simulations and related experiments demonstrate the feasibility of the proposed method and the ability to improve the accuracy of the measurement results further.


Introduction
Fringe projection profilometry (FPP) is an optical three-dimensional (3D) profile measurement technique [1][2][3][4], which is a non-contact active measurement method to obtain the 3D profile of an object by projecting and collecting fringe patterns. With the development of digital projectors and sensor devices, this method achieves high accuracy at a low cost. Therefore, this method has many applications in reverse engineering, mechanical assembly, biomedicine, heritage conservation, and other fields [5][6][7].
The core content in structured light 3D measurement is obtaining the continuous phase of the object [8,9]. The phase acquisition is divided into Fourier transform profilometry (FTP) [10] and phase-shifting profilometry (PSP) [11]. For FTP, the fundamental frequency component is filtered out in the frequency domain after the collected image is transformed by the Fourier transform. Then, the Fourier inverse transform is used to obtain the phase. This method only needs one fringe image [12,13] and has good performance in high-speed measurement, but the measurement accuracy of this method is low. However, according to the principle of PSP, the phase is obtained by projecting multiple groups of phase-shifting fringe patterns to the object and performing the point-to-point operations on the collected fringes. This method effectively avoids the influence between adjacent points and has high measurement accuracy [14]. Since these methods extract the phase by tangent calculation, the phase values range from −π to π. Therefore, a phase unwrapping process is required to obtain a continuous phase.
Among the phase unwrapping methods, they can be broadly classified into two categories, spatial phase unwrapping (SPU) [15] and temporal phase unwrapping (TPU) [16]. The SPU [17][18][19] algorithms unwrap the phase by the phase value of adjacent pixels. However, one problem with this SPU algorithm is that the phase errors can spread to other locations. In contrast, the TPU algorithm does not have this problem [20][21][22][23]. The basic idea of this algorithm is to make the frequency of fringes change with time, and the fringe

Principle
The FPP measurement system comprises a projector, a camera, a processing unit (computer), and a working plane. In a 3D measurement process, a computer generates pre-designed sinusoidal fringe patterns, which are then projected onto the object's surface by a projector. The camera captures the deformed fringes modulated by the object, then uses the phase extraction formula to obtain the wrapped phase. The schematic diagram is shown in Figure 1.

Wrapped Phase Extraction
In N-step PSP, the projected fringe patterns

Wrapped Phase Extraction
In N-step PSP, the projected fringe patterns I p n (n = 1, . . . , N; N ≥ 3) can be denoted as I p n (x p , y p ) = a p + b p cos[2π f x p − 2πn/N], (1) where (x p , y p ) is the pixel coordinate of the projector, a p is the average intensity, b p is the amplitude, f is the frequency of the projected fringe ( f = 1 T , T is the period of the projected fringe), n represents the phase-shifting index, and N is the phase-shifting step number [14].
The intensity of deformed fringe patterns captured by a camera can be expressed as where (x, y) is the pixel coordinate of the camera, A(x, y) is the average intensity relating to the pattern brightness and background illumination, B(x, y) is the intensity modulation relating to the pattern contrast and surface reflectivity, and ϕ(x, y) is the corresponding wrapped phase which can be extracted by [14] ϕ(x, y) = tan −1 ∑ N n=1 I n (x, y) sin(2πn/N) ∑ N n=1 I n (x, y) cos(2πn/N) . ( It can be noted that there are three unknown quantities, A(x, y), B(x, y), and ϕ(x, y), in Equation (2), so at least three fringe patterns are needed to calculate the wrapped phase ϕ(x, y). Due to the different projected fringe frequencies used in the TFH method, the wrapped phase of each frequency needs to be obtained for further unwrapping.

Principle of Phase Unwrapping
The phase values calculated by Equation (3) range from −π to π and they will have a phase jump of 2π. To obtain a continuous phase, the wrapped phase plus integer multiple of 2π is needed, which is given by Φ(x, y) = ϕ(x, y) + 2πk(x, y), where ϕ(x, y) is the wrapped phase, Φ(x, y) is the unwrapped phase, and k(x, y) is the fringe order. The main task of phase unwrapping is to calculate the correct fringe orders k(x, y). In the traditional TFH phase unwrapping algorithm, three groups of fringes with frequencies of f 1 , f 2 , and f 3 ( are projected onto the target object. Then, the wrapped phase ϕ 1 , ϕ 2 , and ϕ 3 are obtained by Equation (3). The synthetic phases ϕ 12d , ϕ 23d , and ϕ 123d can be expressed, respectively, as where subscript d refers to phase difference operation. The equivalent frequencies of synthetic phases meet The obtained synthetic phase is equivalent to the extract phase from equivalent frequency fringes. The equivalent phase with a larger period is obtained by phase difference operation. If the fringe period is chosen appropriately, the synthetic phase ϕ 123 with Sensors 2022, 22, 9388 4 of 16 a period covering the whole measurement field can be obtained after two phase difference operations. Φ(x, y) = ϕ 123d (x, y) is the continuous phase that helps calculate the fringe orders.
The fringe order k of ϕ h and continuous phase Φ h can be obtained by Equation (6) as where round[·] is the symbol indicating the nearest integer to the value. f l and f h are the frequencies of low-frequency and high-frequency fringes, ϕ h is the high-frequency wrapped phase, and Φ l and Φ h are the low-frequency and high-frequency continuous phases [29]. By repeatedly applying the process of recovering the high-frequency phase with the low-frequency phase, the final continuous phase with the highest frequency is obtained. For the traditional TFH phase unwrapping algorithm, it is easy to misjudge the phase when the phase values are closer due to the noise effect in the phase difference calculation process, which will lead to the phase jump with the value of 2π. In addition, the reduction of the equivalent frequency will also affect the measurement accuracy. Figure 2 shows the schematic diagram of the proposed algorithm. Both phase difference and phase sum operations are introduced to obtain synthetic phases, which could conduct the final continuous phase. Firstly, three groups of fringes with different periods are projected, and the periods meet T 1 < T 2 < T 3 . Similar to the traditional method, by calculating the phase difference with Equation (5), the phase differences ϕ 12d and ϕ 23d can be obtained. process, which will lead to the phase jump with the value of 2π . In addition, the reduction of the equivalent frequency will also affect the measurement accuracy. Figure 2 shows the schematic diagram of the proposed algorithm. Both phase difference and phase sum operations are introduced to obtain synthetic phases, which could conduct the final continuous phase. Firstly, three groups of fringes with different periods are projected, and the periods meet 1 2 3 T T T   . Similar to the traditional method, by calculating the phase difference with Equation (5), the phase differences 12d  and 23d  can be obtained. After this first phase difference calculation, the phase with the larger equivalent period is achieved. However, during the second phase difference calculation for 123d  , the positions with very small phase differences are sensitive to noise. The possible noise will make the phase fluctuate, which is easy to cause misjudgment, and eventually lead to phase jumps. To investigate how to obtain robust phase difference results in this process, the lateral change of the two synthetic phases was observed after the first phase difference calculation, and it was found that the phase difference absolute value 12

The Proposed Algorithm
 is mainly divided into two parts, (1) and (2), respectively, as shown in Figure 3a. After this first phase difference calculation, the phase with the larger equivalent period is achieved. However, during the second phase difference calculation for ϕ 123d , the positions with very small phase differences are sensitive to noise. The possible noise will make the phase fluctuate, which is easy to cause misjudgment, and eventually lead to phase jumps. To investigate how to obtain robust phase difference results in this process, the lateral change of the two synthetic phases was observed after the first phase difference calculation, and it was found that the phase difference absolute value |ϕ 12d − ϕ 23d | is mainly divided into two parts, (1) and (2), respectively, as shown in Figure 3a.
positions with very small phase differences are sensitive to noise. The possible noise will make the phase fluctuate, which is easy to cause misjudgment, and eventually lead to phase jumps. To investigate how to obtain robust phase difference results in this process, the lateral change of the two synthetic phases was observed after the first phase difference calculation, and it was found that the phase difference absolute value 12  is mainly divided into two parts, (1) and (2), respectively, as shown in Figure 3a. It can be observed from Figure 3a that the distribution of phase difference in (1) and (2) is discontinuous, and the value in the (1) region is smaller than that in the (2) region. There is a steep jump in the change of phase difference value. In addition, the statistics of It can be observed from Figure 3a that the distribution of phase difference in (1) and (2) is discontinuous, and the value in the (1) region is smaller than that in the (2) region. There is a steep jump in the change of phase difference value. In addition, the statistics of phase differences are shown in Figure 3b. It also can be found that the difference is mainly divided into two parts, (3) and (4). Part (3) corresponds to the region (1), which has smaller values, and part (4) corresponds to the region (2), which has larger values. It can be clearly noticed that there is a gap between part (3) and part (4) from Figure 3b.
Based on the above findings, the threshold can be set to carry out the correct phase difference operation in the second phase difference calculation process. This algorithm uses the Otsu thresholding segmentation algorithm to calculate the intermediate thresholds of part (3) and part (4). The Otsu method is an adaptive threshold segmentation algorithm, and it can quickly calculate the interclass data threshold by the maximum variance between two data classes after segmentation. Therefore, the threshold T h can be obtained by where OTSU(·) is the Otsu operation. By comparing the phase difference absolute value and the threshold value, the improved second phase difference calculation can be expressed as where ϕ id is the wrapped phase difference, f id is the equivalent frequency, and i = 12, 23, 123. The phase difference operation can obtain the phases of larger equivalent periods. When selecting an appropriate combination of periods, the period of the synthetic phase ϕ 123d will cover the whole measurement field, while Φ 123d = ϕ 123d is a continuous phase, which could help calculate the fringe orders. The phase sum operation in the proposed method is also introduced to obtain the synthetic phases with higher equivalent frequency. The phases ϕ 12s and ϕ 23s after the first phase sum calculation are obtained by where subscript s refers to the phase sum operation. After the first phase sum calculation, all the phase values ranging from 0 to 2π are obtained. If these results are used directly for the second phase sum calculation, the phase sum in higher frequencies cannot be obtained correctly. Therefore, the first phase sum results need to be shifted between −π and π by where ϕ is is the phase sum after the shift. Then, the ϕ 123s is the result of the second phase sum calculation, which is expressed as where f is is the equivalent frequency of phase sum and i = 12, 23, 123. The phase sum operation can obtain the phases of smaller equivalent periods. Afterward, the phase sum ϕ 123s is restored to a continuous phase using phase difference ϕ 123d . The phase calculation process is shown in Figure 4.
where is f is the equivalent frequency of phase sum and 12, 23,123 i  . The phase sum operation can obtain the phases of smaller equivalent periods.
Afterward, the phase sum 123s  is restored to a continuous phase using phase difference 123d  . The phase calculation process is shown in Figure 4.
where i  and j  are the lower frequency and higher frequency continuous phase, j  is the higher frequency wrapped phase ( 123 ,12 ,1,12 is the fringe order. The sequence of the whole phase recovery is as follows: The phase is unwrapped from low to high equivalent frequencies step by step.

Mathematical Derivation and Analysis
The sensitivity gain G (between 123d  and 123s  ) of the three-frequency method can be calculated according to the two-frequency method [34][35][36].   Figure 4c, ϕ 123d is a continuous phase (Φ 123d = ϕ 123d ), while ϕ 123s is a wrapped phase. The fringe order k s of ϕ 123s and the final continuous phase Φ 123s can be obtained by where Φ i and Φ j are the lower frequency and higher frequency continuous phase, ϕ j is the higher frequency wrapped phase (i = 123d, 12d, 1, 12s, j = 12d, 1, 12s, 123s), k is the fringe order. The sequence of the whole phase recovery is as follows: The phase is unwrapped from low to high equivalent frequencies step by step.

Mathematical Derivation and Analysis
The sensitivity gain G (between Φ 123d and Φ 123s ) of the three-frequency method can be calculated according to the two-frequency method [34][35][36]. Obtaining Sensors 2022, 22, 9388 7 of 16 As the T 1 , T 2 , and T 3 are greater than zero, the G is always greater than one. The phase sum Φ 123s is G = T 123d /T 123s times more sensitive than Φ 123d . Next, the SNRs of phase difference and phase sum will be discussed. In practice, the phases Φ 1n (x, y), Φ 2n (x, y), and Φ 3n (x, y) are corrupted by additive white Gaussian noise n 1 (x, y), n 2 (x, y), and n 3 (x, y) as where n 1 , n 2 , and n 3 (omit the coordinate notation (x, y) henceforth) are uncorrelated samples with a variance of σ 2 , and h(x, y) (h(x, y) = Φ(x, y)T/2π) is the height contained in the phase. The phase difference Φ 123dn and phase sum Φ 123sn with noise are Then, the SNRs for Φ 123dn and Φ 123sn are where the (x, y) ∈ Ω is the two-dimensional region where the fringe data are well-defined. Since n 1 , n 2 , and n 3 are generated by the same Gaussian zero-mean stationary stochastic process, then the average energies of n 1 − n 2 , n 1 + n 2 ; n 3 − n 2 , and n 3 + n 2 are equal [34][35][36][37]. The SNR gain between Φ 123sn and Φ 123dn is Thus, Φ 123sn (x, y) has G 2 times higher SNR than Φ 123dn (x, y). A phase sum operation can obtain higher G and SNR through mathematical analysis.
In addition, the noise-induced phase error is also Gaussian distributed and is assumed to have a variance of σ 2 Φ according to the study of Ref [29]. After the phase unwrapping is performed, the phase is extended from 2π to 2π f t ( f t is the total number of fringes in the fringe pattern). If the already expanded continuous phase is scaled down to the scale [−π, π), the variance of the phase error is equivalent to being reduced by a factor of f t 2 .
Equation (19) shows that three factors can be used to reduce the impact of noise, but the most convenient way to reduce phase error for a system that already exists is increasing f t . Compared with the DFH method, the TFH method can select the fringes with higher frequency to calculate the phases because it performs the phase difference twice, and it can measure a larger field. This paper introduces the phase sum operation with higher sensitivity gain and SNR into the TFH method. Combining the TFH method with the phase sum operation can obtain better measurement results through mathematical derivation and theoretical analysis.

Phase Calculation from Fringe Patterns with Noise Added
To verify the whole process and the anti-noise performance, random noise with an SNR of 29.8 dB is added to the simulated fringes. The simulated image size is 512 × 512 pixels, and the computer generates three groups of sinusoidal fringe patterns (T 1 = 21, T 2 = 23, T 3 = 25). The phase is calculated by the four-step PSP algorithm. The phase calculation process with noise is shown in Figure 5.
ing t f . Compared with the DFH method, the TFH method can select the fringes with higher frequency to calculate the phases because it performs the phase difference twice, and it can measure a larger field. This paper introduces the phase sum operation with higher sensitivity gain and SNR into the TFH method. Combining the TFH method with the phase sum operation can obtain better measurement results through mathematical derivation and theoretical analysis.

Phase Calculation from Fringe Patterns with Noise Added
To verify the whole process and the anti-noise performance, random noise with an SNR of 29.8 dB is added to the simulated fringes. The simulated image size is 512 512 pixels  , and the computer generates three groups of sinusoidal fringe patterns . The phase is calculated by the four-step PSP algorithm. The phase calculation process with noise is shown in Figure 5.  In the first phase difference calculation process, the statistics of the phase difference absolute value |ϕ 12d − ϕ 23d | with noise are calculated and shown in Figure 6.  It can be observed that the phase difference values after adding noise still remain in the found distribution, which could be divided into two parts. The threshold is    It can be observed that the phase difference values after adding noise still remain in the found distribution, which could be divided into two parts. The threshold is T h = 3.1 using the Otsu method, and it will be used for the second phase difference calculation. The comparison of improved phase difference calculation results and the traditional method are shown in Figure 7.  It can be observed that the phase difference values after adding noise still remain in the found distribution, which could be divided into two parts. The threshold is 3.1 h T  using the Otsu method, and it will be used for the second phase difference calculation. The comparison of improved phase difference calculation results and the traditional method are shown in Figure 7.    Figure 7a shows the results of traditional phase difference calculation. At the position where the phase value is close, due to the effect of noise, the calculated continuous phase has multiple jumps of 2π. Figure 7b demonstrates the result of the improved method proposed in this work, and it can effectively eliminate the phase jumps. To further study the improvement compared with other methods, one of the improved DFH methods [39], the traditional TFH method, and the proposed method are compared. To ensure that each method can finally obtain a continuous phase, the fringe periods of the improved DFH method are set as T 1 = 31 and T 2 = 32, and the fringe periods of the traditional TFH method and the proposed method are set as T 1 = 21, T 2 = 23, and T 3 = 25. All methods are employed to measure a virtual plane, and the same level of random noise is added to fringe patterns. The plane phase of each method is calculated and shown in Figure 8. Comparing these three methods shows that the improved DFH method and the traditional TFH method have some steep phase jumps, and the improved DFH method has more small phase fluctuations in the linear region. In general, using the proposed method can effectively improve the accuracy of the measurement and reduce the phase jumps caused by noise.

Effect of Fringe Period Selection on Phase Calculation
Then, the selection of the fringe period is discussed and the effect of various combinations on the measurement accuracy is studied. It can be predicted that not all combinations are feasible. The combination must meet as , , where m is the image resolution in the fringe coding direction. Since there are three periods that need to be determined step by step, we first fix 3 T in the range of 20 to 80, then fix 2 T in the range of 1 to 3 1 T  , and finally fix 1 T in the range of 1 to 2 1 T  . The random noise is also added to the simulation, and the root mean square error (RMSE) of different combinations is calculated after determining each period. When the maximum pe- 1920 pixels pixels  . The error bar graph for various period combinations is shown in Figure 9. Comparing these three methods shows that the improved DFH method and the traditional TFH method have some steep phase jumps, and the improved DFH method has more small phase fluctuations in the linear region. In general, using the proposed method can effectively improve the accuracy of the measurement and reduce the phase jumps caused by noise.

Effect of Fringe Period Selection on Phase Calculation
Then, the selection of the fringe period is discussed and the effect of various combinations on the measurement accuracy is studied. It can be predicted that not all combinations are feasible. The combination must meet as where m is the image resolution in the fringe coding direction. Since there are three periods that need to be determined step by step, we first fix T 3 in the range of 20 to 80, then fix T 2 in the range of 1 to T 3 − 1, and finally fix T 1 in the range of 1 to T 2 − 1. The random noise is also added to the simulation, and the root mean square error (RMSE) of different combinations is calculated after determining each period. When the maximum period T 3 = 40, 50, 60, 70, we select the T 1 with the smallest error after T 2 is fixed to draw the error bar graph of Φ 123s , which is labeled as (T 2 , T 1 ). To avoid the random fluctuation caused by a single calculation as much as possible, the mean value and peak-valley value are obtained by repeating each combination 20 times. The image size is 1920 pixels × 1920 pixels. The error bar graph for various period combinations is shown in Figure 9.  According to Figure 10, it can be known that those combinations with larger RMSE in Figure 9 have a wider range of fluctuations and errors. This leads to larger phase errors and reduces the measurement accuracy during the phase calculation. Therefore, selecting the appropriate fringe period combination is necessary before conducting the experiments.
To obtain the optimal period combination, those points with lower mean and peakvalley values are selected, as shown in Figure 11. It can be seen that the optimal combination of periods shows a linear trend. Therefore, it can be fitted by a linear function, Equation (21). According to Figure 9, it is known that when the value of T 3 is larger, there are more stable period combinations that meet Equation (20). When the values of T 1 , T 2 , and T 3 are relatively close, the error value and fluctuation are smaller, while the error in the place where the period difference is is larger, and the fluctuation of the calculation results is also larger. Those period combinations are sensitive to noise in the measurement process. Therefore, that should be avoided. To investigate the reasons for the large errors and fluctuations generated by some combinations, the three combinations with the periods of T 1 = 15, T 2 = 23, T 3 = 50, T 1 = 44, T 2 = 47, T 3 = 50, and T 1 = 46, T 2 = 48, T 3 = 50 are conducted without adding noise to reflect the impact of different combinations simply. The errors of the Φ 123d and Φ 123s corresponding to the three combinations are shown in Figure 10.
According to Figure 9, it is known that when the value of 3 T is larger, there are more stable period combinations that meet Equation (20). When the values of 1 T , 2 T , and 3 T are relatively close, the error value and fluctuation are smaller, while the error in the place where the period difference is is larger, and the fluctuation of the calculation results is also larger. Those period combinations are sensitive to noise in the measurement process. Therefore, that should be avoided. To investigate the reasons for the large errors and fluctuations generated by some combinations, the three combinations with the periods of  According to Figure 10, it can be known that those combinations with larger RMSE in Figure 9 have a wider range of fluctuations and errors. This leads to larger phase errors and reduces the measurement accuracy during the phase calculation. Therefore, selecting the appropriate fringe period combination is necessary before conducting the experiments.
To obtain the optimal period combination, those points with lower mean and peakvalley values are selected, as shown in Figure 11. It can be seen that the optimal combination of periods shows a linear trend. Therefore, it can be fitted by a linear function, Equation (21). According to Figure 10, it can be known that those combinations with larger RMSE in Figure 9 have a wider range of fluctuations and errors. This leads to larger phase errors and reduces the measurement accuracy during the phase calculation. Therefore, selecting the appropriate fringe period combination is necessary before conducting the experiments.
To obtain the optimal period combination, those points with lower mean and peakvalley values are selected, as shown in Figure 11. It can be seen that the optimal combination of periods shows a linear trend. Therefore, it can be fitted by a linear function, Equation (21).
Using the fitting function, the values of 1 k , 2 k , 3 k , and 4 k can be obtained, respectively, as Before the fringe projection, the maximum period needs to be first determined for the measurement, and then the periods of the other two groups are obtained by function calculation. The accuracy of the measurement results can be further improved through this process.

Experiment
A fringe projection measurement system is set up for the experiment, as shown in Figure 12. The system includes a CCD camera, a digital projector, a high-precision motorized linear translation stage, a checkerboard, and a computer. The CCD camera is the digital camera IMAVISION MER-231-41GM-P of the Mercury series from Daheng Imaging, Using the fitting function, the values of k 1 , k 2 , k 3 , and k 4 can be obtained, respectively, as Before the fringe projection, the maximum period needs to be first determined for the measurement, and then the periods of the other two groups are obtained by function calculation. The accuracy of the measurement results can be further improved through this process.

Experiment
A fringe projection measurement system is set up for the experiment, as shown in Figure 12. The system includes a CCD camera, a digital projector, a high-precision motorized linear translation stage, a checkerboard, and a computer. The CCD camera is the digital camera IMAVISION MER-231-41GM-P of the Mercury series from Daheng Imaging, with a resolution of 1920 × 1200. The motorized linear translation stage is a GCD-203300M from Daheng Optics, with an accuracy of 0.001 mm. The projector is an Epson CH-TW5600, with a resolution of 1920 × 1080. The size of the checkerboard square is 15 mm. Under the premise of satisfying the period relationship of Equation (20), the maximum period is selected to be 3 60 T  , so that the continuous phases can be obtained for fair and reliable comparison when other methods are calculated, which will be discussed Under the premise of satisfying the period relationship of Equation (20), the maximum period is selected to be T 3 = 60, so that the continuous phases can be obtained for fair and reliable comparison when other methods are calculated, which will be discussed later. According to the period optimization method, the other two groups of fringe periods are T 2 = 58 and T 1 = 56. Firstly, the smooth surface continuous objects are measured, as shown in Figure 13. Under the premise of satisfying the period relationship of Equation (20), the maximum period is selected to be 3 60 T  , so that the continuous phases can be obtained for fair and reliable comparison when other methods are calculated, which will be discussed later. According to the period optimization method, the other two groups of fringe periods are 2 58 T  and 1 56 T  . Firstly, the smooth surface continuous objects are measured, as shown in Figure 13.  Figure 13a shows the deformed fringe patterns, and Figure 13b illustrates the wrapped phases calculated by PSP. Figure 13c describes the phase difference maps. After two phase difference operations, the continuous phase covering the whole measurement field is obtained. Figure 13d shows the phase sum maps, and the equivalent frequency after two phase sum operations is much more than that of all projected fringes. Figure 13e demonstrates the histogram graph of 12 , and it still conforms to the found distribution. The threshold obtained by the Otsu segmentation algorithm is 2.51. Figure 13f illustrates the final continuous phase 123s  , obtained by the introduced phase recovery sequence. The calibration is performed using a translation stage and a checkerboard to present the results in the world coordinate system. The translation stage moves with a known height multiple times to calculate the phase-to-height parameter [40]. The calibrated height volume is 100 mm. Twenty-five checkerboard images are captured and calculated for camera internal, external, and distortion parameters by Zhang's camera calibration technique [41]. After calibration, the subsequent experimental results are converted to the world coordinate system.  Figure 13a shows the deformed fringe patterns, and Figure 13b illustrates the wrapped phases calculated by PSP. Figure 13c describes the phase difference maps. After two phase difference operations, the continuous phase covering the whole measurement field is obtained. Figure 13d shows the phase sum maps, and the equivalent frequency after two phase sum operations is much more than that of all projected fringes. Figure 13e demonstrates the histogram graph of |ϕ 12d − ϕ 23d |, and it still conforms to the found distribution. The threshold obtained by the Otsu segmentation algorithm is 2.51. Figure 13f illustrates the final continuous phase Φ 123s , obtained by the introduced phase recovery sequence.
The calibration is performed using a translation stage and a checkerboard to present the results in the world coordinate system. The translation stage moves with a known height multiple times to calculate the phase-to-height parameter [40]. The calibrated height volume is 100 mm. Twenty-five checkerboard images are captured and calculated for camera internal, external, and distortion parameters by Zhang's camera calibration technique [41]. After calibration, the subsequent experimental results are converted to the world coordinate system.
To make the error comparison fairer and more reliable, a 24-step phase-shifting plus multi-frequency method [28] is used to obtain the ground truth in the experiment. The periods of the projected fringe are determined by T i = m 2 i−1 (i = 1, 2, 3, . . .), where the projector's resolution is m × n, and m and n are the horizontal and vertical resolutions of the projector, respectively.
Then, the comparison of the proposed method, the traditional TFH method, and the improved DFH method is performed. To ensure that these three methods can obtain continuous phases, and make their period values close to proceed with a fairer comparison, the periods of DFH are selected as T 2 = 68 and T 1 = 66, while the combination of the proposed method and the traditional TFH method are selected as T 3 = 60, T 2 = 58, and T 1 = 56. The target of this combination with a smaller period value is to avoid f id > m (i = 12, 23) after the first phase difference calculation, where m is the horizontal resolution of the projector (m = 1920 in this research). Additionally, the periods of the multi-frequency method are determined by T i = 1920 2 i−1 (i = 1, 2, 3, 4, 5, 6), and the continuous phase is solved with T 6 = 60 as the ground truth. In this case, each method can satisfy its own period relationship.
A ceramic standard gauge block with a height of 20 mm was measured. A linear area was selected on its upper surface to compare these three methods, as shown in Figure 14. T   as the ground truth. In this case, each method can satisfy its own period relationship.
A ceramic standard gauge block with a height of 20 mm was measured. A linear area was selected on its upper surface to compare these three methods, as shown in Figure 14. As shown in Figure 14, the result of the proposed method is closer to the actual height, and the RMSE is smaller. A part of a gourd and a statue were measured to compare the three methods further. The comparison between these methods was performed, as shown in Figure 15.  As shown in Figure 14, the result of the proposed method is closer to the actual height, and the RMSE is smaller. A part of a gourd and a statue were measured to compare the three methods further. The comparison between these methods was performed, as shown in Figure 15. T   as the ground truth. In this case, each method can satisfy its own period relationship.
A ceramic standard gauge block with a height of 20 mm was measured. A linear area was selected on its upper surface to compare these three methods, as shown in Figure 14. As shown in Figure 14, the result of the proposed method is closer to the actual height, and the RMSE is smaller. A part of a gourd and a statue were measured to compare the three methods further. The comparison between these methods was performed, as shown in Figure 15.  Analyzing the results in Figure 15, the proposed method has fewer jump points than the other two methods. To show more clearly the differences between several methods, their Euclidean distance (ED = (x − x ) 2 + (y − y ) 2 + (z − z ) 2 , where (x, y, z) and (x , y , z ) are the coordinates of two space points, respectively) maps were calculated from the ground truth, as shown in Figure 16. Figure 16 shows that the overall RMSE of the proposed method is minimal, and there are fewer jump points. The above comparison proves that the proposed method can further improve measurement accuracy and reduce jump errors.
The optimization method of the fringe period is verified. According to the above conclusion, the optimal fringe periods are T 1 = 56, T 2 = 58, T 3 = 60, and the results are compared with several other periods' results, as shown in Figure 17.
The optimal combination has much fewer jumps and is almost the same as the ground truth. The error of other combinations is more prominent, especially combination 2. To more clearly observe the difference between the results of several combinations and the ground truth, the ED maps were calculated as shown in Figure 18.
Analyzing the results in Figure 15, the proposed method has fewer jump points than the other two methods. To show more clearly the differences between several methods, their Euclidean distance ( , , x y z and   , , x y z    are the coordinates of two space points, respectively) maps were calculated from the ground truth, as shown in Figure 16.  Figure 16 shows that the overall RMSE of the proposed method is minimal, and there are fewer jump points. The above comparison proves that the proposed method can further improve measurement accuracy and reduce jump errors.
The optimization method of the fringe period is verified. According to the above conclusion, the optimal fringe periods are 1 2 3 56, 58, 60 T T T    , and the results are compared with several other periods' results, as shown in Figure 17.    Figure 16 shows that the overall RMSE of the proposed method is minimal, and there are fewer jump points. The above comparison proves that the proposed method can further improve measurement accuracy and reduce jump errors.
The optimization method of the fringe period is verified. According to the above conclusion, the optimal fringe periods are 1 2 3 56, 58, 60 T T T    , and the results are compared with several other periods' results, as shown in Figure 17.  The optimal combination has much fewer jumps and is almost the same as the ground truth. The error of other combinations is more prominent, especially combination 2. To more clearly observe the difference between the results of several combinations and the ground truth, the ED maps were calculated as shown in Figure 18. As shown in Figure 18, the ED map and RMSE of the optimal combination are the smallest, while the error of combination 2 is the largest, and the error level of combination 3, 4, and 5 is between them. Figure 18 shows that the optimal combination has the least phase jumps and a higher measurement accuracy. This experiment proves that the opti- As shown in Figure 18, the ED map and RMSE of the optimal combination are the smallest, while the error of combination 2 is the largest, and the error level of combination 3, 4, and 5 is between them. Figure 18 shows that the optimal combination has the least phase jumps and a higher measurement accuracy. This experiment proves that the optimal period selection criterion can obtain better measurement results.
The above experiments can show that the proposed method is feasible and has a high anti-noise ability. By comparing the results obtained by several methods, it can be seen that the proposed method has better measurement accuracy, also with the optimal fringe period combination.

Conclusions
An improved TFH synthesis phase unwrapping algorithm is proposed in this research. The phase sum operation is introduced into the TFH algorithm, and the phase sum and phase difference operations are improved. In the process of phase difference calculation, according to the found distribution of phase difference value, the Otsu method is used to calculate the threshold to help phase judgment and reduce the effect of noise. In the process of phase sum operation, the phase with a much higher frequency than those of projected fringes can be realized. The continuous phase difference is used to assist in phase sum unwrapping. The improved phase sum and phase difference operation jointly achieve the effect of accuracy improvement and error reduction. The fringe period selection method is deduced based on the distribution of error values for different period combinations, and better results can be obtained using the optimal combination. Simulation and experimental results show the feasibility and anti-noise ability of the proposed method, which can improve the accuracy of measurement results to a certain extent.