Open Access
Add to BibSonomy
Abstract
In phase-shifting fringe projection profilometry, fringe harmonics caused by device nonlinearities as well as other factors may badly ruin the measurement results. Generally, the used phase-shifting algorithm enables restraint of effects of harmonics below a certain order depending on the number of phase shifts. When reducing the number of phase shifts for efficiency, high order harmonics will affect the phase-measuring results because of aliasing caused by insufficient sampling rate. To overcome this issue, this paper suggests a non-filtering technique operating in frequency domain, that enables improvement of measurement accuracy by eliminating effects of high order harmonics. With this technique, the phase-shifting algorithm is restated as a process of retrieving the fundamental complex fringes from the phase-shifted fringe patterns. Implementing a Fourier transform to this calculated complex fringe pattern, the actual fundamental signals and the aliased harmonics have their own lobes with separated peaks in the frequency domain. We reconstruct each order of the aliased harmonics by exploiting their relations with the fundamental signals and then estimate their magnitudes by using the spectral peaks. Instead of directly filtering the fringe spectrum, we subtract spectra of the harmonics from Fourier transform of the just calculated complex fringes, so that the Fourier spectrum of the fundamental fringes without harmonics is recovered through an iterative operation. Further, the phase map is measured accurately. Simulation and experimental results confirm that this proposed method can significantly suppress effects of fringe harmonics. Meanwhile, by taking advantage of non-filtering, it effectively preserves the edges and details of the measured surfaces from being blurred.
© 2023 Optica Publishing Group under the terms of the Optica Open Access Publishing Agreement
1. Introduction
In phase-shifting fringe projection profilometry [1–3], the nonlinearities of devices, especially of the projector, induce harmonics in fringe profile, thus ruining the measurement results by inducing ripple-like artifacts in them [4]. The most straightforward technique for solving this problem is to reshape the fringe profile in advance. For this purpose, a photometric calibration [5,6] or a trial experiment [7–9] is usually implemented before measurement to determine the input-output response curves (or functions) of the devices. Using these response curves allows us to actively compensate for the nonlinearity errors in the early stage of measurement by generating nearly perfect sinusoidal fringes. When the calibration data are not available or the device nonlinearities vary over time [10], however, we have to design or enhance the phase measuring algorithms for suppressing effects of the device nonlinearities as far as possible.
In the recent decades, many phase-measuring algorithms have been developed, that enable restraining effects of the fringe harmonics on the measurement results. For example, the synchronous detection algorithm [11], as one of the earliest phase-shifting algorithms and the most popularly used algorithm today, uses uniform phase shifts and is immune to up to the $({N - 2} )$th order of fringe harmonics when the number of phase shifts is $N$[12,13]. Similarly, most phase-shifting algorithms behave like a digital band-pass filter exactly retrieving fundamental fringe signals. With them, the errors are mainly sourced from harmonic aliasing caused by the insufficient sampling rates. Generally, increasing the number of phase shifts is helpful for restraining effects of more harmonics, at the expense of increased time duration for image capturing thus taking more risks of object motions [14] and illumination fluctuations [15]. Because fringe patterns of different frequencies have harmonics staggered over the patterns, it is possible to eliminate these fringe harmonics by capturing fringe patterns of multiple frequencies [16,17]. This type of methods is efficient in juggling temporal phase unwrapping and harmonics elimination, but they still require capturing a number of fringe patterns so that sufficiently many independent equations are available for determine coefficients of certain harmonics. Selecting non-uniform phase shifts enables us to avoid aliasing when sampling few fringe patterns, in which case the harmonics can be decoupled from one another by solving a complicated system of nonlinear equations [18–20]. The most difficult case occurs when only few, e.g., three or four, fringe patterns having uniform phase shifts are available. In this situation, the standard method of removing the high order harmonics generally involves a procedure of processing each fringe pattern using a low-pass [21] or band-pass [22] filter, having a disadvantage of blurring edges and details of the measured object.
To overcome the limitations of the existing techniques, this paper suggests, to the best of our knowledge, a novel non-filtering method that works with few fringe patterns and enables improvement of measurement accuracy by eliminating effects of high order harmonics. In doing it, we retrieve a complex fringe pattern first by using the standard phase-shifting algorithm, and then calculate Fourier spectrum of this complex fringe pattern, in which the fundamental signals and the aliased harmonics have separated lobe peaks. Next, we reconstruct each order of the aliased harmonics by exploiting their relations with the fundamental signals and then estimate their magnitudes according to the ratios of peak heights of spectra. Instead of directly filtering the fringe spectrum, we subtract spectra of the harmonics from Fourier transform of the just calculated complex fringes, so that the spectrum of the fundamental fringes without harmonics is recovered iteratively. Further, the phase map is measured accurately. Simulation and experimental results have confirmed that this proposed method can significantly suppress effects of fringe harmonics. Meanwhile, by taking advantage of non-filtering, it effectively preserves the edges and details of the measured surfaces from being blurred.
2. Background
2.1 Synchronous detection phase-shifting algorithm and phase-to-depth conversion
Phase-shifting fringe projection profilometry works with a measurement system shown in Fig. 1. It is mainly composed of a projector and a camera. The projector casts a sequence of sinusoidal fringe patterns onto the measured object and correspondingly the camera grabs the deformed fringe patterns. Analyzing the deformed patterns results in their phase maps from which the object depths are further reconstructed.
In this measurement, capturing the sequence of phase-shifted fringe patterns corresponds to a procedure of sampling temporal fringe signals. If the used devices are perfectly linear, the temporal fringe signal for each pixel $({x,y} )$ has the general form
where $a({x,y} )$ is the background intensity and $b({x,y} )$ is the modulation at $({x,y} )$; $\phi ({x,y} )$ is the phase that has a mapping relationship with depth variation; and t represents time and ${\omega _0}{\; }$represents temporal frequency of this fringe signal.Assuming sampling is uniform and the sampling period is
the $n$th sample is represented with where N is the number of samples in a $2\pi $ period and $n = 0,1, \cdots N - 1$. Here, we use $\hat{I}[{x,y,\; n} ]$ with its variables enclosed in square brackets denoting a temporal discrete signal, which is distinguished from the continuous one $\hat{I}({x,y,t} )$. In this case, the phase is calculated by using the synchronous detection algorithm [11]. The underlying principle is to retrieve the fundamental complex fringesIn measurement practices, however, the phase measuring results using this synchronous detection algorithm may have errors. First, the nonlinearities of the used devices (especially of the projector) and some other factors make the fringe signal non-sinusoidal, or in other words induce harmonics in the fringes. Second, the limited sampling rate leads to aliasing. We shall address these issues in the next subsection.
When the fringe phases are calculated, we convert them into the object depths. For this purpose, we should calibrate the measurement system first. There have been many techniques developed for performing this task [23–25]. Here, we employ a simple implicit mapping function between the unwrapped phase $\mathrm{\Phi }({x,y} )$ and the depth $h({x,y} )$
whose coefficients ${D_1}({x,y} )$, ${D_2}({x,y} )$, ${D_3}({x,y} )$ and ${D_4}({x,y} )$ are determined through a simple calibration method based on cross-ratio invariance [23]. Through Eq. (6), the phase errors induce ripple-like artifacts in the profile of the depth map.2.2 Errors caused by aliasing of harmonics
Because of the device nonlinearities and other factors, the observed fringe signals are usually not perfectly sinusoidal but contains infinite harmonics, in which case, Eq. (1) becomes [26]
We considered the principle of synchronous detection phase-shifting algorithm afresh. It is to calculate the integration in Eq. (9) for the fundamental term ${c_1}({x,y} )$ and then calculate its phase angle, i.e., $\textrm{arg}[{{c_1}({x,y} )} ]$. In practical measurement, we must implement a numerical integration in Eq. (4) for calculating ${c_1}({x,y} )$ from discrete samples. In this procedure, aliasing caused by inadequate sampling rate induces errors.
Using the sampling period in Eq. (2), the $n$th sample of the fringe signal in Eq. (8) is
Considering its $m$th term, if m and k satisfy $m - k = Nl$ with l being an integer, we have ${c_m}({x,y} )\exp [{j({2\pi m/N} )n} ]= {c_m}({x,y} )\exp [{j({2\pi k/N} )n} ]$. This phenomenon means that, in this case, aliasing occurs between the $m$th and $k$th terms. As a result, in the N phase-shifted fringe patterns captured by the camera, i.e., in $\hat{I}[{x,y,n} ]$ for $n = 0,1 \cdots N - 1$, only N frequency components are observable, namely,
Specifically, what we calculated by using Eq. (4) is
Equation (13) reveals that, with the $N$-step synchronous detection algorithm, the harmonics of orders $1 + lN$ ($l ={\pm} 1, \pm 2, \cdots $) affect the phase measuring results. To depict the errors caused by aliasing of harmonics, we conducted a simulation here by predefining a phase map and then measuring it using the three-step synchronous detection algorithm. Figure 2(a) shows the first one of the simulated three fringe patterns. In these fringe patterns, harmonics were introduced by setting a gamma of 1.55 for the projector. The background intensities and modulations were simulated to be Gaussian-shaped with a 50% decrease at corners. Figures 2(b) and 2(c), respectively, show the wrapped and unwrapped phase maps reconstructed from the simulated fringe patterns. Figure 2(d) shows the same phase map but the carrier has been removed, from which we observed ripple-like artifacts covering the phase map. Figure 2(e) exhibits the phase errors obtained by subtracting the predefined phases from the calculated ones. According to Eq. (13), these ripple-like errors are caused by aliasing of some special harmonics, e.g., the $- 2$nd and $4$th orders.
Fig. 2. A simple simulation about the effects of fringe harmonics on phase measuring when using three-step phase-shifting. (a) The first one of phase-shifted fringe patterns. (b) Wrapped phases calculated from phase-shifted fringe patterns including (a). (c) Unwrapped phases. (d) Calculated phase map with the carrier removed. (e) Phase errors. In (b)-(e), the unit is radian.
Because aliasing occurs in the early stage of measurement when capturing the phase-shifted fringe patterns, its effects cannot be removed in the temporal and temporal frequency domains. Specially, when using $3$-, $4$-, and $5$-step algorithms, the harmonics in aliasing are listed in Table 1. Form these results, we know that an effective solution to the aliasing issue is to increase the sampling rate, i.e., increasing the number of phase shifts. If only few fringe patterns are available, however, removal of the errors caused by aliasing is a challenging task. Note that, regarding Table 1, the similar results have been reported in early literature [12,13]. We repeat listing them here just for the convenience of deducing the non-filtering algorithm in the next section.
Table 1. Harmonics in aliasing with the synchronous detection phase-shifting algorithm
3. Non-filtering algorithm for harmonics elimination
3.1 Separability of harmonics in spatial frequency domain
As the aforementioned, the device nonlinearities may produce harmonics in the fringes with the phase-shifting fringe projection technique, but these harmonics do not necessarily affect the phase measuring results. The influences of harmonics are mainly dependent on the sampling, and in other words, on the selected phase shifts. If some harmonics are not aliased with the fundamental fringe signal, it is possible to remove their effect by designing a temporal phase-shifting algorithm. Using uniform phase shifts, however, harmonic aliasing is inevitable. Since aliasing is caused by time sampling, its induced errors cannot be removed in the temporal domain or in the temporal frequency domain. In order to solve this problem, we analyzed the separability of harmonics in spatial frequency domain.
Note that, in fringe projection technique, the fringe phase $\phi ({x,y} )$ has a linear carrier, so that each term in Eq. (13) is very similar to an amplitude modulated signal and different terms have different carrier frequencies. By calculating spatial Fourier transform of Eq. (13), we have its frequency domain representation
Because the terms in Eq. (15) have different linear carriers, they appear in frequency domain like a row of peaks separate from one another.
Following from Fig. 2, we continued the simulation. Using Eq. (4), the calculated fundamental term ${\hat{c}_1}({x,y} )$ does not equal to the true fundamental term ${c_1}({x,y} )$. By performing Fourier transform of ${\hat{c}_1}({x,y} )$, we got its frequency domain representation ${\hat{C}_1}({u,v} )$. Figure 3(a) shows the magnitudes of ${\hat{C}_1}({u,v} )$ as a function of frequencies. In Fig. 3(a), the highest peak corresponds to the spectrum of the true fundamental term ${c_1}({x,y} )$. On its both sides, we observe a sequence of separate peaks along $u$-axis (i.e., along the direction perpendicular to the fringes). These peaks correspond to different harmonics. Typically, the peaks corresponding harmonics of orders −2 and 4, have the second and third largest heights, respectively, and are conspicuous in Fig. 3(a).
Fig. 3. The simulation of implementing the proposed non-filtering method. (a) The spectral magnitudes of the calculated fundamental term having harmonics. (b) The spectral magnitudes with the harmonics subtracted using the proposed method. (c) Wrapped phases calculated from (b). (d) Unwrapped phases of (c). (e) Phase map without carrier. (f) Phase errors. In (c)-(f), the unit is radian.
3.2 Principle of fringe harmonics subtraction
The previous subsection reveals that the aliased harmonics sourced from time sampling appear like a row of separate peaks in the spatial frequency domain. To remove these harmonics, the most popularly used technique is to process ${\hat{C}_1}({u,v} )$ using a band-pass filter, so that the true fundamental term ${C_1}({u,v} )$ is extracted. Instead of filtering ${\hat{C}_1}({u,v} )$, some techniques directly process the captured fringe patterns $\hat{I}[{x,y,n} ]$ using a low- or band-pass filter. Doing so will deliver the same result because the used filters are linear. In measurement, however, these filtering-based techniques have a common drawback of inducing blurring in the measurement results. The reason is that the terms in Eq. (15) may have relatively wide bandwidths, especially when measuring objects having complex shapes. In this case, the fundamental term and the harmonics have separate spectral peaks, but their spectra still overlap at the bottom between the peaks. To solve this problem, we suggest a non-filtering technique here by estimating and subtracting the harmonics.
Using Eqs. (13) and (14), we get the estimates of ${c_1}({x,y} )$ and $\phi ({x,y} )$, i.e., ${\hat{c}_1}({x,y} )$ and $\hat{\phi }({x,y} )$, respectively. Substituting them into Eq. (16), we have the estimate of ${\varPsi _{1 + lN}}({u,v} )$ as
Note that ${C_{1 + lN}}({u,v} )$ and ${\varPsi _{1 + lN}}({u,v} )$ have the very similar profile. They are different just in magnitude by a scale ${d_{1 + lN}}$. They have peaks at the same position in the frequency domain. For determine the scale ${d_{1 + lN}}$, we seek the peak position of $|{{{\hat{\varPsi }}_{1 + lN}}({u,v} )} |$, namely $({{u_p},{v_p}} )$ as
and further the scale is determined byFrom ${\hat{C}_1}({u,v} )$, subtracting the estimated harmonics yields the fundamental term more accurate than ${\hat{C}_1}({u,v} )$, that is
Here, the harmonics above the $({1 + LN} )$th order are truncated because they are small and omittable. When ${C_1}({u,v} )$ is obtained from Eq. (20), ${c_1}({x,y} )$ is recalculated by using an inverse Fourier transform
and then the more accurate phase $\phi ({x,y} )$ is calculated by using Eq. (5).This proposed algorithm requires performing a Fourier transform. Generally, methods involving Fourier transform, e.g., those using a low- or band-pass filter to process fringes, are considered to have a deduced measurement resolution. The reason for this is the same as that with Fourier transform profilometry [27,28]. When measuring objects having complex shapes, spectra of different components are not separate but partially overlapping. Spectral leakage causes the spectrum of the fringe component to be spread across other frequencies. Therefore, it is difficult to design a perfect filter to separate them. Our proposed algorithm does not require a filtering processing, so that spectral overlap and leakage have a much less impact on its processing results. As a result, the fundamental fringe component can be separated with a much higher accuracy.
3.3 Iterative steps
In the previous subsection, because ${\hat{\varPsi }_{1 + lN}}({u,v} )$ and ${\hat{d}_{1 + lN}}{\; }$ are estimates rather than accurate values, the recalculated ${c_1}({x,y} )$ and $\phi ({x,y} )$ still contains residual errors. To improve the accuracy, we suggested an iterative method whose procedure is summarized as follows with the numbers of iterations being indicated by using the superscripts in parentheses.
Step 1: Determine the initial values for the fundamental terms ${c_1}({x,y} )$ and phases $\phi ({x,y} )$. Using Eqs. (13) and (14), we calculate ${\hat{c}_1}({x,y} )$ and $\hat{\phi }({x,y} )$, and take them as the initial values $\hat{c}_1^{(0 )}({x,y} )$ and ${\hat{\phi }^{(0 )}}({x,y} )$, respectively. By implementing Fourier transform of $\hat{c}_1^{(0 )}({x,y} )$, we have $\hat{C}_1^{(0 )}({u,v} )$.
Step 2: Update ${c_1}({x,y} )$ and $\phi ({x,y} )$. If the $i$th iterative results $\hat{c}_1^{(i )}({x,y} )$ and ${\hat{\phi }^{(i )}}({x,y} )$ are obtained, we substitute them into Eq. (17) for calculating $\hat{\varPsi }_{1 + lN}^{(i )}({u,v} )$ (for $l ={-} L, \cdots , - 1,1, \cdots L$). We seek the peak position for each $|{\hat{\varPsi }_{1 + lN}^{(i )}({u,v} )} |$ and determine $\hat{d}_{1 + lN}^{(i )}$ through Eq. (19). By substituting $\hat{\varPsi }_{1 + lN}^{(i )}({u,v} )$ and $\hat{d}_{1 + lN}^{(i )}$ (for $l ={-} L, \cdots , - 1,1, \cdots L$) into the right-hand side of Eq. (20), we calculate the new spectrum of fundamental term, i.e., $\hat{C}_1^{({i + 1} )}({u,v} )$. Implementing an inverse Fourier transform of $\hat{C}_1^{({i + 1} )}({u,v} )$, we calculate $\hat{c}_1^{({i + 1} )}({x,y} )$ and finally its phases ${\hat{\phi }^{({i + 1} )}}({x,y} )= \arg [{\hat{c}_1^{({i + 1} )}({x,y} )} ]$. This step always uses $\hat{C}_1^{(0)}$ on the right sides of Eqs. (19) and (20).
Step 3: Repeat running Step 2 until the algorithm converges.
Using this method, we continued the simulation in Fig. 3. Figure 3(a) shows the magnitudes of ${\hat{C}_1}({u,v} )$ as a row of peaks corresponding to the fundamental fringe signal and harmonics. By using the iterative procedure just presented, the harmonics were estimated and subtracted from the spectrum. The result is given in Fig. 3(b) where the peaks of harmonics have been eliminated. Using this resulting spectrum of fundamental fringes, the wrapped fringe phases were recalculated and shown in Fig. 3(c). Figure 3(d) shows the phase-unwrapping results. Figure 3(e) gives the phase with their carrier removed. Figure 3(f) shows the phase errors. By comparing Fig. 3(f) and Fig. 2(e), we find that the errors caused by the aliased harmonics have been suppressed significantly.
4. Further numerical simulations
In the above sections, we have performed a numerical simulation in Figs. 2 and 3. Its purpose is just to explain the principle and operation of the proposed non-filtering technique. In this section, we made further simulations in order to examine performances of this technique. In these simulations, the fringe parameters were set to be the same as those in Figs. 2 and 3, but the measured phase map was assumed to have edges in its profile.
When the number of phase shifts is three, the simulation results are shown in Fig. 4, with its three columns corresponding to different means to process the fringe patterns. The leftmost column shows the results of simply using three-step phase-shifting algorithm without correcting the phase errors. In it, Fig. 4(a) is the first one of phase-shifted fringe patterns. Figure 4(d) shows its calculated phase map, where the ripple artifacts caused by harmonics are observably large. Figure 4(g) gives the errors of the calculated phase map. These errors have a structure parallel to the fringes and have spatial frequencies three times higher than the fringe frequencies, which are typically caused by the $- 2$nd and $4$th order harmonics. The middle column was obtained by using a filtering method. Below the fringe pattern in Fig. 4(b), Fig. 4(e) shows the calculated phase map. For getting it, we retrieved the fundamental complex fringes using Eq. (13) and then filtered them using a band-pass filter. From this phase map, we can observe that the error artifacts at most points have been removed. At the same time, however, the edges of phase profile get blurred. For this reason, large errors appear near edges as we have seen from the error map in Fig. 4(h). The rightmost column gives the results of using our proposed non-filtering method. Its phase map in Fig. 4(f) and error map in Fig. 4(i) illustrate that the harmonics-caused error artifacts have been eliminated and the edges of the phases were well preserved.
Fig. 4. Further numerical simulation results when using three-step phase-shifting. The columns, from left to right, show the results of simply using three-step phase-shifting algorithm without correcting the phase errors, using filtering technique, and using the proposed non-filtering method, respectively. The panels in each column, from top to bottom, show the first one of the phase-shifted fringe patterns, the calculated phase map in radians without carrier, and the phase errors in radians, respectively.
When using four-step phase-shifting, we repeated doing the simulation. The simulation results are shown in Fig. 5, with its layout being the same as that of Fig. 4. These results are very similar to those in Fig. 4. The biggest differences appear in the leftmost column, when simply using four-step phase-shifting algorithm without correcting the errors. In this case, the phase errors are mainly caused by the $- 3$rd and $5$th order harmonics, and hence the artifacts, as shown in Figs. 5(d) and 5(g), have decreased amplitudes and frequencies four times higher than the fringe frequencies. In the middle column, it is clearly visible in Figs. 5(e) and 5(h) that the edges of phases have been blurred by the filtering method. In the rightmost column, when using the proposed technique, the errors were removed and the edges were perfectly protected as exhibited in Figs. 5(f) and 5(i).
Fig. 5. Further numerical simulation results when using four-step phase-shifting. The layout is the same as that of Fig. 4.
With the synchronous detection algorithm, the calculated fundamental term ${\hat{c}_1}({x,y} )$ in Eq. (13) may have infinite harmonics, but their amplitudes rapidly decline as their orders increase in both positive and negative directions. We truncated high order terms in Eq. (20) when estimating true fundamental term ${c_1}({x,y} )$. Therefore, it is necessary to investigate effects of truncation on accuracies of the proposed method. We implemented the simulation many times by subtracting different number of harmonics and showed the corresponding RMS values of residual errors in Fig. 6. In it, Figs. 6(a) and 6(b) were obtained with three-step and four-step phase- shifting, respectively. These simulation results show that, with the proposed technique, subtracting more harmonics is helpful for reducing residual errors but these subtracted harmonics, as their orders become higher, have weaker (even negligible) effects on the results.
Fig. 6. Simulation results of RMS errors versus the number of subtracted harmonics. (a) RMS errors when using three-step phase-shifting. (b) RMS errors when using four-step phase-shifting.
Another concern is regarding convergence of the proposed algorithm. In Figs. 7(a) and 7(b), we plotted the curves of the RMS errors against the number of iterations when using three-step and four-step phase-shifting, respectively. It is shown from them that the proposed algorithm rapidly converges in the first several iterations. As the number of iterations increases, the RMS values of residual errors approach certain amounts depending on noise levels. In fact, after the algorithm converges, the final residual errors are mainly caused by random noise rather than by the fringe harmonics. Especially, under the noise-free condition, the RMS phase errors should approach 0 as the number of iterations increases.
Fig. 7. Simulation results of the RMS errors converging with the number of iterations under different noise conditions. (a) Results of using three-step phase-shifting. (b) Results of using four-step phase-shifting.
5. Experiments and discussions
In this section, we experimentally verified the feasibility of the proposed technique by measuring a practical object. The measurement system has presented in Section 2 by Fig. 1. It consists of a camera (AVT G192B, 1600 × 1200) and a projector (Liying LDLP-500).
When the number of phase shifts is three and the relative phase step is $2\pi /3$ radians, the first one of the phase-shifted fringe patterns, each of which has a size 791 × 791 pixels, is shown in Fig. 8(a). Simply using the three-step phase-shifting algorithm, we calculated the fringe phases as shown in Fig. 8(b) whose values are wrapped in the range from $- \pi $ to $\pi $ radians. Figure 8(c) is the unwrapped phase map. In this work, for reliably unwrapping phases, we excluded invalid pixels by thresholding the calculated fringe modulations [29] and then unwrapped the phase maps using a spatial method based on geometric constraints [30]. It is worth noting that phase measuring results in this work are independent of the used phase-unwrapping techniques. Temporal phase-unwrapping methods [31] can also be used for the same purpose at the expense of capturing more fringe patterns. Figure 8(d) gives the phases with their carrier removed, from which ripple-like artifacts caused by the fringe harmonics are apparently visible.
Fig. 8. Experimental results of measuring a practical object. (a) The first one of captured phase-shifted fringe patterns. (b) Wrapped phases calculated from (a) by using three-step phase-shifting algorithm. (c) Unwrapped phases. (d) Phase map with the carrier removed. (e) Phase errors. In (b)-(e), the unit is radian.
For quantitively evaluating errors of the calculated phases, we additionally measured the same object using 12-step phase-shifting, by which effects of the harmonics below the tenth order were eliminated and influence of noise was also suppressed. Therefore, its measurement accuracy is much higher than those with few phase shifts, and its measured phase map was taken as a benchmark to evaluate measurement accuracies of different methods in our experiments. In Fig. 8, the phase errors of simply using three-step algorithm were calculated by subtracting from its result the benchmark phase map. These phase errors are shown in Fig. 8(e), from which we observed the errors have frequencies three times higher than the fringes, which are typically induced by the $- 2$nd and $4$th order harmonics.
We used the proposed non-filtering technique to remove the phase errors. Figure 9 summarizes the procedure. We calculated the fundamental complex term from the fringe patterns by using Eq. (13). Its Fourier transform has magnitudes in Fig. 9(a). In it, the highest peak corresponds to the true fundamental term. On its both sides along $v$-axis (perpendicular to the fringe direction), some small peaks corresponding harmonics of orders $- 2$, $4$ and so forth are visible. [In Fig. 9(a), there are four additional peaks appearing in the four quadrants far away from $u$- and $v$-axes. They are caused by the pixel structure of the projector.] Following the steps in Section 3.3, by iteratively estimating and subtracting spectra of harmonics, Fig. 9(b) shows the finally obtained spectrum corresponding to the fundamental term. In comparison with Fig. 9(a), the spectral lobes of the harmonics have been removed and at the same time the useful high-frequency components are preserved in this spectrum. After performing an inverse Fourier transform, the phases of the fundamental fringes were calculated and show in Fig. 9(c) with their values within the range from $- \pi $ to $\pi $ radians. Figure 9(d) shows the unwrapped phase map. Figure 9(e) shows the phase map without carrier. Figure 9(f) gives the phase error map obtained by subtracting the benchmark phases mentioned above. From these phase-measuring results, it is evident that the harmonics-caused errors have been significantly suppressed by using our proposed non-filtering technique.
Fig. 9. Experimental results of implementing the proposed non-filtering method. (a) The spectral magnitudes of the calculated fundamental term having harmonics. (b) The spectral magnitudes with the harmonics subtracted using the proposed method. (c) Wrapped phases calculated from (b). (d) Unwrapped phases of (c). (e) Phase map without carrier. (f) Phase errors. In (c)-(f), the unit is radian.
Figure 9 shows the phase measuring result with the newly proposed algorithm. As mentioned in Section 2.1, the phase errors degrade measurement accuracies by inducing ripple-like artifacts in the profile of the depth map. For more accurately evaluating effectiveness of the proposed method in improving measurement accuracies, we converted the phase map into depth map and compared it with the results of other techniques. Figure 10(a) is the depth map reconstructed from the phases calculated by simply using three-step algorithm without correcting the harmonics-caused errors. Figure 10(b) is the depth map with the fringe harmonics were eliminated by using a filtering method. Figure 10(c) shows the result of that the projector nonlinearities were compensated for through a photometric calibration. Figure 10(d) displays the depth map with the errors were removed using the proposed non-filtering method. Below them, Figs. 10(e) through 10(h), respectively, are their error maps obtained by subtracting the benchmark depths measured using a 12-step phase-shifting technique.
Fig. 10. Experimental results of measuring object depths using three-step phase-shifting. (a) The depth map reconstructed simply using three-step phase-shifting algorithm without correcting the harmonics-induced errors. (b) The depth map with the fringe harmonics removed using a filtering method. (c) The depth map with the projector nonlinearities compensated for through a photometric calibration. (d) The phase map obtained with the effects of harmonics suppressed using the proposed non-filtering method. (e)-(h) The errors of (a)-(d), respectively. In (a)-(h), the unit is millimeter.
For further confirming the performance of this proposed method in restraining effects of harmonics when capturing very few patterns, we changed the number of phase shifts as four and measured the same object once again. The procedure was the same as that with three-step phase shifting. We just give the results in Fig. 11. In it, the used methods and the layout are the same as those in Fig. 10.
Fig. 11. Experimental results of measuring object depths using four-step phase-shifting. The layout is the same as that of Fig. 10.
From Figs. 10(e) and 11(e), we observe that the measurement errors with a $N$-step phase-shifting algorithm appear as ripple-like artifacts having frequencies N times higher than the frequencies. In Figs. 10(f) and 11(f), the filtering method enables removing harmonics-caused artifacts but it induces blurring in the calculated depth maps, especially at edges. In Figs. 10(g) and 11(g), the errors have been suppressed because a photometric calibration was implemented in advance before doing the measurements. The error maps in Figs. 10(h) and 11(h) demonstrate that the proposed non-filtering technique enabled us to eliminate effects of fringe harmonics depending on very few fringe patterns in the absence of priori knowledge or calibration data of the projector.
To display the residual artifacts more clearly, we selected a cross-section perpendicular to the fringes from the reference plane board behind the measured object. Figure 12 plots the error distributions of different methods along this cross-section. When simply using phase-shifting algorithm, the errors, depending on the number of phase shifts, have relatively large fluctuations. When using the filtering method, the reconstructed depths became smooth but large errors appeared near the boundaries of the plane. With the photometric calibration method, the errors were suppressed. When using the proposed method, the error artifacts were suppressed but a photometric calibration was avoided.
Fig. 12. Measurement errors along a cross-section of the reference plane. (a) The errors when the number of phase shifts is three. (b) The errors when the number of phase shifts is four. In the figure legend, Methods 1 through 4 are the methods simply using phase-shifting algorithm, with the harmonics removed using a filtering method, with the projector nonlinearities compensated for through a photometric calibration, and proposed in this paper, respectively.
For quantitively evaluating the accuracies the used methods achieve, Table 2 lists the RMS values of the errors in the above experiments. Note that, when using the filtering method, extremely large errors occurred at edges and the image boundaries. In the third column of Table 2, we calculated the RMS errors of the filtering method and gave the results of two cases, i.e., including and excluding edge points, on the left and right sides of the slash sign, respectively. The data in Table 2 shows that the proposed method achieved higher accuracies over others. Note that the data in Table 2 cannot be used to measure the accuracy limit of this proposed method. In fact, this method is an accuracy-enhancing method. In other words, if more precision devices are used, this method enables further improving their accuracy.
Table 2. The RMS depth errors (mm) with different phase measuring methods
Figure 13 investigates the convergence of the proposed method with its horizontal and vertical axes denoting the number of iterations and the RMS depth errors, respectively. It is shown that this algorithm has converged in the first several iterations to certain values depending on noise levels. In other words, after the algorithm converges, the residual errors in a depth map are mainly caused by random noise. Naturally, using four-step phase-shifting have a final RMS error a little bit smaller than that of using three-step phase-shifting. A related issue is regarding computational efficiency of the proposed method. In these experiments, we used a laptop with Intel Core i5-7300HQ CPU and 16GB memory. When simply using phase-shifting algorithm without correcting the errors, it took 0.51 and 0.57 seconds for processing three- and four-step phase-shifted fringe patterns, respectively. When the proposed method was implemented, these computational time durations became 1.52 and 1.58 seconds, respectively. It has a satisfying computational efficiency meeting demands from most ordinary applications.
Fig. 13. With the proposed non-filtering method, the RMS depth errors decrease as the number of iterations increases.
Besides device nonlinearities, random noise is a main error-inducing factor. In the above experiments, after implementing the proposed technique, the residual errors are predominantly caused by noise. Noise level in a captured fringe pattern partially depends on the used camera and projector. Therefore, using more advanced devices is helpful for enhancing the signal-to-noise ratio (SNR) of images. By use of temporal randomness of noise, capturing and averaging a number of frames is the most effective approach to reduce the noise variance, at an expense of increased time duration for image capturing. Because noise induced by complicated physical factors is generally treated as a Gaussian variable according to the central limit theorem, spatial low-pass filter is usually used to smooth it, by which details in the image may be blurred depending on the frequency response of the filter. In addition, illumination may slightly fluctuate during measurement thus inducing errors. From Figs. 10(h) and 11(h), we observe some residual errors which are mainly caused by illumination fluctuations. This error-inducing factor cannot be modeled as fringe harmonics, and therefore is difficult to remove using a phase-shifting algorithm, but this issue can be overcome by using the method based on the fringe histograms [15].
The proposed algorithm enables restraining effects of the aliased harmonics on phase measuring, and the aliasing phenomenon is related to the selection of phase shifts. It has been known that, if the phase shifts, together with $2\pi $, have no common divisors, aliasing is avoided. In this case, the harmonic-caused phase errors can be removed through alternate iterative least-squares algorithms [18,19]. These algorithms fail in dealing with the situation of that the phase shifts are uniform over a $2\pi $ period, because harmonic aliasing occurs. Our newly proposed method works with uniform phase shifts and enables suppressing effects of harmonic aliasing. In fact, by use of separability of harmonics in spatial frequency domain, this method in principle is still valid when the phase shifts are nonuniform. In this case, more peaks having smaller space between them will appear in the frequency domain and a more sophisticated method is required to deal with them.
6. Conclusion
In this paper, we have proposed a non-filtering method, operating in frequency domain, that enables us to eliminate fringe harmonics caused by projector nonlinearities in the phase-shifting fringe projection profilometry. Firstly, it calculates the fundamental complex term from the phase-shifted fringe patterns and implements Fourier transform. By use of separability of fringe harmonics in the spatial frequency domain, it then estimates magnitudes of spectral peaks of the harmonics. Iteratively subtracting the estimated spectra of harmonics enables us to get the fundament complex fringes without influences of harmonics, and finally the phase map is accurately recovered. Simulation and experimental results have demonstrated that this proposed method provides some advantages over others. Typically, it does not require calibrating the projector nonlinearities in advance thus having a higher flexibility in practical measurements. It significantly suppresses effects of fringe harmonics by capturing very few fringe patterns. Meanwhile, by taking advantage of non-filtering, it effectively preserves the edges and details of the measured surfaces from being blurred.
Funding
National Natural Science Foundation of China (51975345).
Disclosures
The authors declare no conflicts of interest.
Data availability
Data underlying the results presented in this paper are not publicly available at this time but may be obtained from the authors upon reasonable request.
References
1. V. Srinivasan, H. C. Liu, and M. Halioua, “Automated phase-measuring profilometry of 3-D diffuse objects,” Appl. Opt. 23(18), 3105–3108 (1984). [CrossRef]
2. S. Zhang, “Absolute phase retrieval methods for digital fringe projection profilometry: A review,” Opt. Lasers Eng. 107, 28–37 (2018). [CrossRef]
3. C. Zuo, S. Feng, L. Huang, T. Tao, W. Yin, and Q. Chen, “Phase shifting algorithms for fringe projection profilometry: a review,” Opt. Lasers Eng. 109, 23–59 (2018). [CrossRef]
4. H. Guo, H. He, and M. Chen, “Gamma correction for digital fringe projection profilometry,” Appl. Opt. 43(14), 2906–2914 (2004). [CrossRef]
5. C. R. Coggrave and J. M. Huntley, “High-speed surface profilometer based on a spatial light modulator and pipeline image processor,” Opt. Eng. 38(9), 1573–1581 (1999). [CrossRef]
6. S. Kakunai, T. Sakamoto, and K. Iwata, “Profile measurement taken with liquid-crystal grating,” Appl. Opt. 38(13), 2824–2828 (1999). [CrossRef]
7. T. Hoang, B. Pan, D. Nguyen, and Z. Wang, “Generic gamma correction for accuracy enhancement in fringe-projection profilometry,” Opt. Lett. 35(12), 1992–1994 (2010). [CrossRef]
8. Z. Li and Y. Li, “Gamma-distorted fringe image modeling and accurate gamma correction for fast phase measuring profilometry,” Opt. Lett. 36(2), 154–156 (2011). [CrossRef]
9. S. Ma, C. Quan, R. Zhu, L. Chen, B. Li, and C. J. Tay, “A fast and accurate gamma correction based on Fourier spectrum analysis for digital fringe projection profilometry,” Opt. Commun. 285(5), 533–538 (2012). [CrossRef]
10. B. Li, Y. Wang, J. Dai, W. Lohry, and S. Zhang, “Some recent advances on superfast 3D shape measurement with digital binary defocusing techniques,” Opt. Lasers Eng. 54, 236–246 (2014). [CrossRef]
11. J. H. Bruning, D. R. Herriott, J. E. Gallagher, D. P. Rosenfeld, A. D. White, and D. J. Brangccio, “Digital wavefront measuring interferometer for testing optical surfaces and lenses,” Appl. Opt. 13(11), 2693–2703 (1974). [CrossRef]
12. K. A. Stetson and W. R. Brohinsky, “Electro-optic holography and its application to hologram interferometry,” Appl. Opt. 24(21), 3631–3637 (1985). [CrossRef]
13. H. Guo and M. Chen, “Fourier analysis of the sampling characteristics of the phase-shifting algorithm,” Proc. SPIE 5180, 437–444 (2004). [CrossRef]
14. Z. Liu, P. C. Zibley, and S. Zhang, “Motion-induced error compensation for phase shifting profilometry,” Opt. Express 26(10), 12632–12637 (2018). [CrossRef]
15. Y. Lu, R. Zhang, and H. Guo, “Correction of illumination fluctuations in phase-shifting technique by use of fringe histograms,” Appl. Opt. 55(1), 184–197 (2016). [CrossRef]
16. C. Jiang, S. Xing, and H. Guo, “Fringe harmonics elimination in multi-frequency phase-shifting fringe projection profilometry,” Opt. Express 28(3), 2838–2856 (2020). [CrossRef]
17. Y. Lee, Y. Mao, and Z. Chen, “Fast combined-frequency phase extraction for phase shifting profilometry,” Opt. Express 30(25), 45288–45300 (2022). [CrossRef]
18. H. Zhu and H. Guo, “Anti-aliasing phase reconstruction via a non-uniform phase-shifting technique,” Opt. Express 30(3), 3835–3853 (2022). [CrossRef]
19. H. Zhu and H. Guo, “Alternate iterative least-squares algorithm based on nonuniform phase shifting for suppressing nonlinearity errors in fringe projection profilometry,” IEEE Trans. Instrum. Meas. 71, 1–13 (2022). [CrossRef]
20. H. Zhu and H. Guo, “Joint least-squares algorithm correcting phase-shift errors and detector nonlinearity simultaneously in phase-shifting interferometry,” Appl. Opt. 62(7), 1829–1839 (2023). [CrossRef]
21. Y. Hu, J. Xi, J. Chicharo, and Z. Yang, “Improved Three-step Phase Shifting Profilometry Using Digital Fringe Pattern Projection,” International Conference on Computer Graphics, Imaging and Visualisation (CGIV'06), Sydney, NSW, Australia, 2006, pp. 161–167.
22. S. Ma, R. Zhu, C. Quan, B. Li, C. J. Tay, and L. Chen, “Blind phase error suppression for color-encoded digital fringe projection profilometry,” Opt. Commun. 285(7), 1662–1668 (2012). [CrossRef]
23. R. Zhang and H. Guo, “Depth recovering method immune to projector errors in fringe projection profilometry by use of cross-ratio invariance,” Opt. Express 25(23), 29272–29286 (2017). [CrossRef]
24. H. Guo, “Least-squares calibration method for fringe projection profilometry,” Opt. Eng. 27(9), 1171–1187 (2004). [CrossRef]
25. S. Xing and H. Guo, “Iterative calibration method for measurement system having lens distortions in fringe projection profilometry,” Opt. Express 28(2), 1177–1196 (2020). [CrossRef]
26. S. Xing and H. Guo, “Directly recognizing and removing the projector nonlinearity errors from a phase map in phase-shifting fringe projection profilometry,” Opt. Commun. 435, 212–220 (2019). [CrossRef]
27. M. Takeda and K. Mutoh, “Fourier transform profilometry for the automatic measurement of 3-D object shapes,” Appl. Opt. 22(24), 3977–3982 (1983). [CrossRef]
28. X. Su and Q. Zhang, “Dynamic 3-D shape measurement method: A review,” Opt. Lasers Eng. 48(2), 191–204 (2010). [CrossRef]
29. M. Zhao, L. Huang, Q. Zhang, X. Su, A. Asundi, and Q. Kemao, “Quality-guided phase unwrapping technique: comparison of quality maps and guiding strategies,” Appl. Opt. 50(33), 6214–6224 (2011). [CrossRef]
30. Y. An, J. Hyun, and S. Zhang, “Pixel-wise absolute phase unwrapping using geometric constraints of structured light system,” Opt. Express 24(16), 18445 (2016). [CrossRef]
31. C. Zuo, L. Huang, M. Zhang, Q. Chen, and A. Asundi, “Temporal phase unwrapping algorithms for fringe projection profilometry: A comparative review,” Opt. Lasers Eng. 85, 84–103 (2016). [CrossRef]
