Clustering-driven residue filter for profile measurement system

Jun Jiang; Jun Cheng; Ying Zhou; Guang Chen

doi:10.1364/JOSAA.28.000214

Journal of the Optical Society of America A
Vol. 28,
Issue 2,
pp. 214-221
(2011)
•https://doi.org/10.1364/JOSAA.28.000214

Clustering-driven residue filter for profile measurement system

Jun Jiang, Jun Cheng, Ying Zhou, and Guang Chen

Open Access

Get PDF
Email
Share
Get Citation
Copy Citation Text
Jun Jiang, Jun Cheng, Ying Zhou, and Guang Chen, "Clustering-driven residue filter for profile measurement system," J. Opt. Soc. Am. A 28, 214-221 (2011)

Export Citation
- BibTex
- Endnote (RIS)
- HTML
- Plain Text
Citation alert
Save article
Spotlight Summary

More Like This

Unsupervised-clustering-driven noise-residue filter for phase images
Jun Jiang, et al.
Appl. Opt. 49(11) 2143-2150 (2010)

Noise-residue filtering of interferometric phase images
Giovanni Nico
J. Opt. Soc. Am. A 17(11) 1962-1974 (2000)

Dynamic phase measurement in shearography by clustering method and Fourier filtering
Yuanhao Huang, et al.
Opt. Express 19(2) 606-615 (2011)

Related Topics
Table of Contents Category
- Image Processing
Optics & Photonics Topics
?

The topics in this list come from the Optics and Photonics Topics applied to this article.

About this Article
History
- Original Manuscript: September 28, 2010
- Manuscript Accepted: December 2, 2010
- Published: January 26, 2011
Virtual Issues
February 18, 2011 Spotlight on Optics

Abstract

The profile measurement system is widely used in industrial quality control, and phase unwrapping (PU) is a key technique. An algorithm-driven PU is often used to reduce the impact of noise-induced residues to retrieve the most reliable solution. However, measuring speed is lowered due to the searching of optimal integration paths or correcting of phase gradients. From the viewpoint of the rapidity of the system, this paper characterizes the noise-induced residues, and it proposes a clustering-driven residue filter based on a set of directional windows. The proposed procedure makes the wrapped phases included in the filtering window have more similar values, and it groups the correct and noisy phases into individual clusters along the local fringe direction adaptively. It is effective for the tightly packed fringes, and it converts the algorithm-driven PU to the residue-filtering-driven one. This improves the operating speed of the 3D reconstruction significantly. The tests performed on simulated and real projected fringes confirm the validity of our approach.

1. INTRODUCTION

Surface profile measurement by noncontact optical methods has been extensively applied to automated optical inspection, defect detection of surface shape, and solid modeling. In such a context, active methods based on the projection of structured light are very popular and commercially successful. A typical fringe projection profilometry (FPP) system consists of a projection unit, an image acquisition module, and a processing module [1, 2]. The implementation, as illustrated in Fig. 1, involves (i) projecting a set of structured patterns (e.g., sinusoidal fringes) onto the observed object surface, (ii) recording the fringe images that are phase modulated by the object’s height distribution, (iii) resolving the phase modulation using a fringe analysis technique and a proper phase unwrapping (PU) algorithm to acquire continuous phase distribution, and (iv) calibrating the system for mapping the continuous phase distribution to real-world 3D coordinates.

Common to most applications of FPP is the fact that the interesting physical information, e.g., shape, deformation, and displacement, is related to the absolute phase ϕ. However, only phase ψ of modulo- $2 π$ , the so-called wrapped phase, can be measured limited to the interval $[- π, π)$ . Thus, in order to recover the continuous phase distribution, a 2D PU algorithm has to be used. Formally, PU is defined as given the wrapped phase $ψ \in [- π, π)$ , find the absolute phase $ϕ \in [- \infty, \infty)$ , and the relationship between these two kinds of phases is as follows:

ψ = W (ϕ) = ϕ - 2 π 〈 \frac{ϕ}{2 π} 〉,

where W is the wrapping operator and

〈 \cdot 〉

rounds its argument to the closest integer. In fact, PU is mathematically an ill-posed problem if no further constraint is taken into account, because W is surjection rather than bijection. Thus, a constraint is given that the absolute value of the gradient of two arbitrary adjacent absolute phases is less than π [3]. However, PU is still a problem due to the existence of residues, which are the potential source of phase-error propagation [4]. Thus, the problem for PU is to minimize the impact of those residue-induced false jumps on the estimated unwrapped phase image. There are two ways to realize this purpose: algorithm driven and residue filtering driven. The existing algorithm-driven PUs for 3D reconstruction can be classified into four major categories: (i) path-following [5, 6], (ii) minimum-norm [7, 8, 9], (iii) model-based [10], and (iv) Bayesian-based [11, 12]. The advantage of these algorithm-driven approaches is that they minimize the number of phase jumps larger than π between neighboring points of the estimated unwrapped phases through searching optimal integration paths or correcting the phase gradients to retrieve the most reliable solution. The disadvantage is the computational burden of the optimization procedure. For example, the computational complexity of the minimum spanning tree method is

O (N \log_{2} N)

[13], the computational complexity of the branch-cut method is

O (N^{2})

, and the computational complexity of the network-flow method reaches

O (N^{3})

In comparison, the motivation of the residue-filtering- driven approach is to identify residues and eliminate them in order to provide an irrotational field that can be integrated along any path. This makes 3D reconstruction faster with similar accuracy. Capanni et al. [14] proposed an adaptive median filter based on some histogram properties of the filtering window to identify noisy pixels in a wrapped phase image. However, this method introduces a smoothing impact on noise-free phases, and noisy phases cannot be located in areas with relatively dense noise. Lee et al. [15] presented a sigma filter as an adaptive extension of the averaging filter. It averages only selected phase values within the $k_{1} \times k_{2}$ neighborhood along the directions of phase fringes. This filter, however, results in an undesired phase smoothing in residue-free areas. Nico [16] was confined to the description of configurations of residues, and designed a special filter. This method is based on the observation of the noisy phases inducing adjacent phase inconsistencies. A combinatorial analysis, however, shows that the filter designed is very complicated, owing to various configurations of the adjacent residues. Fornaro et al. [17] addressed a maximum-likelihood PU based on the phase difference, and the noise is postprocessed in the unwrapped phase image using a median filter. Recently, Qian et al. [18] proposed a windowed filter in the Fourier domain, and the filtered amplitude is used as a quality map. Hitherto, no new method can remove the residues thoroughly within the minimum operating time.

In this paper, we analyze the configurations of noise- induced residues, introduce clustering into the problem of residue filtering, design a set of directional windows to involve more similar phases, and propose a feasible filter to remove residues without the need of searching optimal integration paths or correcting the phase gradients. The remainder of the paper is organized as follows. Section 2 briefly introduces the principle of the FPP system. Section 3 analyzes the characteristics of residues. Section 4 addresses our motivations and the proposed filter. Section 5 demonstrates the validity of our approach in the profile measurement system with arbitrarily arranged devices. Section 6 presents the conclusion.

2. PRINCIPLE OF FRINGE PROJECTION PROFILOMETRY

In the FPP system, a set of structured fringe patterns is projected on the object surface. The captured images are phase modulated by the object height distribution. Thus, the task is to analyze the deformed fringes to determine the profile of the measured object. The main steps consist of fringe analysis, PU, and system calibration.

2A. Fringe Analysis

In the fringe projection technique, an object shape is evaluated through a phase distribution extracted from the captured fringes. In a general FPP system, the illumination is divergent and a nonlinear carrier will be introduced. Figure 2 shows the optical geometry of the FPP system. When a sinusoidal fringe is projected onto a 3D diffusing object, its deformed image can be expressed as

I (x, y) = A (x, y) + B (x, y) \cos [2 n (x, y) π + ψ (x, y) + φ],

where

(x, y)

is the coordinate on the reference plane,

A (x, y)

is the reflected intensity of the object surface,

B (x, y)

is the fringe contrast, φ denotes the phase-shifting increment,

ψ (x, y) \in [- π, π)

, and

n (x, y) \in Z

. The task is to obtain

ψ (x, y)

from captured fringe images. According to the principle of phase-shifting technique, when N measurements

I_{1}

I_{2}, \dots, I_{N}

of I are made with a phase increment

\frac{2 π}{N}

, we can acquire

ψ (i, j) = \arctan \frac{\sum_{k = 1}^{N} I_{k} (i, j) \sin (\frac{2 k π}{N})}{\sum_{k = 1}^{N} I_{k} (i, j) \cos (\frac{2 k π}{N})},

where

(i, j)

is the pixel coordinate on the imaging plane. It is notable that although the fringe spacing is varied in the x direction, the phase-shifting increment is constant for any point

(x, y)

on the object.

2B. Phase Unwrapping

In practice, the presence of shadows, low modulations, fringe discontinuities, and noise makes PU difficult. In the next sections, we will analyze how to eliminate the impact of noise-induced residues on PU to improve the running speed of the system. Given that all the residues are removed, PU can be operated in the wrapped phase image using Eqs. (4, 5) to acquire a continuous phase image:

ϕ (i, j) = ϕ (1, 1) + \sum_{i^{'} = 1}^{i - 1} Δ_{1} ϕ (i^{'}, 1) + \sum_{j^{'} = 1}^{j - 1} Δ_{2} ϕ (i, j^{'}),

Δ_{1} ϕ (i, j) = ψ_{1} (i, j) = Δ_{1} ψ (i, j) + 2 π n_{1} (i, j), Δ_{2} ϕ (i, j) = ψ_{2} (i, j) = Δ_{2} ψ (i, j) + 2 π n_{2} (i, j),

where

ϕ (i, j)

is the unwrapped phase,

Δ_{1} ϕ (i, j)

Δ_{2} ϕ (i, j)

are the discrete partial derivatives of neighboring unwrapped phases,

Δ_{1} ψ (i, j)

Δ_{2} ψ (i, j)

are the discrete partial derivatives of neighboring wrapped phases, and

n_{1} (i, j)

n_{2} (i, j) \in Z

are selected to make

ψ_{1} (i, j)

ψ_{2} (i, j) \in [- π, π)

, respectively.

2C. System Calibration

In Fig. 2, the unwrapped phase $ϕ_{B}$ at point B must equate to the phases at points Q and A, respectively, which is the essence of calibration for an arbitrarily arranged FPP system. Du and Wang [19] derived a mathematical description of the out-of-plane height distribution. The analytic equation is as follows:

h = z_{Q} = \frac{1 + c_{1} ϕ_{B} + (c_{2} + c_{3} ϕ_{B}) i_{B} + (c_{4} + c_{5} ϕ_{B}) j_{B}}{d_{o} + d_{1} ϕ_{B} + (d_{2} + d_{3} ϕ_{B}) i_{B} + (d_{4} + d_{5} ϕ_{B}) j_{B}},

where h or

z_{Q}

is the out-of-plane height of point Q on the object;

ϕ_{B}

is the unwrapped phase of point B on the imaging plane;

(i_{B}, j_{B})

is the pixel coordinate of B; and

c_{1} - c_{5}

and

d_{0} - d_{5}

are constants produced by the system geometric param eters, including the lens focuses of the camera and projector, the relative positions among the CCD image plane, projection plane, and reference plane.

To calculate the out-of-plane height h, the unknown coefficients must be determined first. A proper way is to use the nonlinear least-squares algorithm, e.g., the Levenberg–Marquardt method. The least-squares error can be expressed as

E = \sum_{k = 1}^{m} [\frac{1 + c_{1} ϕ_{k} + (c_{2} + c_{3} ϕ_{k}) i_{k} + (c_{4} + c_{5} ϕ_{k}) j_{k}}{d_{o} + d_{1} ϕ_{k} + (d_{2} + d_{3} ϕ_{k}) i_{k} + (d_{4} + d_{5} ϕ_{k}) j_{k}} - z_{k}^{g}]^{2},

where

z^{g}

denotes the heights of the gauge blocks and m is the number of points on the gauge blocks used in the calculation. A larger m generally yields a higher accuracy. Two gauge objects with uniform but different heights are required for calibrating the system, because using a single gauge block with a uniform height would produce indeterminate solutions for nonlinear equations.

3. CHARACTERISTICS OF NOISE-INDUCED RESIDUES

3A. Localization of Residues

We assume that the projected fringes are sampled sufficiently, in that, Nyquist sampling theorem is satisfied, and no abrupt shape change occurs on the object. A general method utilizes the irrotational property to check phase inconsistency [8], and a residue is located by evaluating the sum of the phase gradients counterclockwise (or clockwise) around each set of four adjacent points:

Δ (i, j) = (- 1 / 2 π) [ψ_{1} (i, j + 1) - ψ_{1} (i, j) - ψ_{2} (i + 1, j) + ψ_{2} (i, j)],

where

Δ (i, j) \in {- 1, 0, + 1}

ψ_{1} (i, j)

and

ψ_{2} (i, j)

are explained in Eq. (5). For a given

2 \times 2

pixel neighborhood, if

Δ (i, j)

is different from zero, the

2 \times 2

area is called a residue in this paper. If the summation takes the value

+ 1

, it is referred to as positive residue (positive polarity); if the summation takes the value

- 1

, it is referred to as negative residue (negative polarity). Figure 3 shows a phase image and its corresponding residue image.

3B. Configurations of Residues

The noise-induced residue is caused by the random fluctuation of phases due to noise in the fringe images. Of the four phase gradients in Eq. (8), only three are independent. For a given $2 \times 2$ pixel neighborhood, it can be shown that the summation defined in Eq. (8) can be different from zero only if an odd number, one or three, of its four phase gradients round to a nonzero value. Some triads cannot occur, and only triads satisfying the following rules are observed: (i) two consecutive phase gradients greater than π in modulus must be of an opposite sign and (ii) two nonconsecutive phase gradients greater than π in modulus must share the same sign, and the third one must have the opposite sign. As Jiang and Cheng [20] depicted, most of the residues have only one phase gradient greater than π in modulus. From a quantitative viewpoint, the adjacent residues are the majority (more than 70% [16]). In this paper, we classify the noise-induced residues into four configurations: (i) a couple of adjacent opposite-sign residues in the form of left–right or up–down, (ii) a couple of adjacent opposite-sign residues in the form of diagonal, (iii) a couple of disjoint opposite-sign residues, and (iv) more than two residues joined together. Figure 4 depicts the four categories of residues. Besides, Karout et al. [21] indicate that a single residue only occurs close to the borders of the phase image—the simple fact of their border location causing their opposite-sign residue to lie outside the field of measurement.

4. CLUSTERING-DRIVEN RESIDUE FILTER

4A. Motivations

The operating speed of the PU algorithm determines the running speed of a FPP system. Although the algorithm-driven PUs can minimize the adverse impact of residues, it is very time consuming to search integration paths or correct the phase gradients. Thus, the performance of the system might be reduced drastically. From the viewpoint of rapidity, the path-independent PU based on residue filtering is superior to the algorithm-driven one. Consequently, the main task reduces to filtering the residues effectively.
Because residues are the potential source of phase-error propagation in the process of PU, we only have to filter the phases within a residue area rather than all the phases in the wrapped phase image. In an actual FPP system, residues are the minority relative to correct phases, and it is reasonable to regard these noisy phases inducing residues as outliers or anomalies locally. From the viewpoint of anomaly detection, the correct phases belong to large and dense clusters, while the noisy phases (or anomalies) belong to small or sparse clusters. Thus, we formulate the problem of residue filtering as a simple problem of 1D data clustering, as illustrated in Fig. 5.
For steep slope areas, the fringe is tightly packed. Squared windows, such as the one for the boxcar filter, will contain more dissimilar phases. This will destroy the continuity of fringes and make the PU incorrect. Thus, we want to design directional windows to make the wrapped phases within the filtering window have more similar values. This makes filtering more effective and preserves the details of the fringes.

4B. Identification of the Filtering Window

To filter the noise-induced residues, the local wrapped phases in the operating window have to be unwrapped, because these discontinuous phases are improper to measure their similarity, and we are not able to identify the local fringe orientation as Lee et al. [15] depicted. We should choose an effective PU to unwrap the pixels within the window to continuous phases reliably.

We use 16 directional windows, as shown in Fig. 6. The window size is set to $9 \times 9$ empirically, and the windows enclose the four categories of residues as shown in Fig. 4. Only the white pixels in windows are included for computation. The window to be selected should be approximately parallel to the local fringe. In the implementation, these windows are convolved with the phase image. The center of the adjacent residue area serves as the center of the windows. We take advantage of local statistics to analyze the similarities of local unwrapped phases, and corresponding variances are computed in all 16 windows individually using Eqs. (9, 10):

v = \frac{1}{N} \sum_{i = 1}^{N} (x_{i} - \bar{x})^{2},

\bar{x} = \frac{1}{N} \sum_{i = 1}^{N} x_{i},

where N is the number of the local unwrapped phases,

x_{i}

is the ith phase, v and

\bar{x}

are the variance and mean of the N phases. Note that the phases, coincided with the white pixels in windows, do not include the pixels in the adjacent residue area. The window, with the minimum variance, is selected as the filtering window.

4C. Unsupervised Clustering

Because we are discussing a clustering problem, the k-means and the single-linkage (SL) algorithms are selected as candidates for their excellent performance. For the problem of 1D data clustering, the computational complexity of the SL algorithm is $O (n^{2})$ with the n samples, and the complexity of the k-means algorithm is $O (n c T)$ with c clusters and T iterations [22]. Because the number of clusters of an arbitrary residue area cannot be definite, SL, with unsupervised attribute, is selected as the basis of the proposed noise-residue-filtering algorithm rather than k-means clustering.

4D. Proposed Filtering Algorithm

Locate the residues in the wrapped phase image. Residues are tested using Eq. (8). An accompanying flag matrix, the same size as the corresponding wrapped phase image, is introduced to mark the locations of the residues. The $2 \times 2$ pixel neighborhood forming a residue is marked “1” (white), and the other pixels are marked “0” (black).
Identify the filtering window. The window size is set to $9 \times 9$ empirically. We choose Costantini’s approach [8] to unwrap the pixels locally within the windows to continuous phases. The sixteen directional windows in Fig. 6 convolve with a residue area, and the corresponding variances are computed individually. The window, with the minimum variance, is selected as the filtering window, as described in Subsection 4B.
Locate the noisy phases and remove them. The wrapped correspondences of the local unwrapped phases in the filtering window serve as the input of the unsupervised classifier. The SL algorithm is used to classify the selected wrapped phases ${ψ_{1}, ψ_{2}, \dots}$ in the filtering window into different groups ${c_{1}, c_{2}, \dots}$ . The phases within the minimal cluster $c_{\min}$ [23] are regarded as the noisy phases based on the assumption of an unsupervised anomaly detection. The assumption is that normal data instances belong to large and dense clusters, while anomalies either belong to small or sparse clusters [24]. Thus, the noisy phases are replaced by the median phase within the maximal cluster $c_{\max}$ [25].
Judge whether there is still any residue. If any residue occurs, go back to step 1 until all the residues are removed.

After all the residues have been removed, $Δ_{1} ϕ (i, j)$ and $Δ_{2} ϕ (i, j)$ in Eq. (5) can be estimated by $ψ_{1} (i, j)$ and $ψ_{2} (i, j)$ reliably. Then we can unwrap the filtered wrapped phase image using Eq. (4) to acquire consistent unwrapped phase image.

5. EXPERIMENTS

In this section, the filter is used as a preprocessing step of PU to remove noise-induced residues. The filtering performance is evaluated by applying the filter to synthetic and real phase signals. We demonstrate the capability of the proposed filter in three aspects:

processing tightly packed fringes with dense noise- induced residues effectively,
removing noisy phases inducing residues without reducing the detailed information, and
substituting the residue-filtering-driven PU for the algorithm-driven one reliably.

We test our algorithm on simulated and real phase images. The simulated data allow a quantitative validation of the method. We compare the proposed filter with the complex average filter [26] and the local histogram-based filter [14] from the viewpoints of residue reduction and fringe preservation, while we use real data to compare the path-independent PU based on the proposed filter with the algorithm-driven methods [6, 7, 8] from the viewpoint of rapidity. The accuracy is also demonstrated after calibrating the FPP system in a less controlled situation. Besides, the reason why we choose these comparative algorithms is that they are representative and widely used in practice.

We simulate four $256 \times 256$ phase-shifted fringes as the object to be constructed. The fringe images are contaminated by dense random noise, and the corresponding phase image is shown in Fig. 7a. In the proposed algorithm, the residues are filtered iteratively, because new residues may occur when previous residues are removed. This means that a wrapped phase, not inducing a residue, is not always a noise-free one. Figure 7b is the rewrapped phase reconstruction using the proposed filter. The fringes are clear and unbroken, indicating that the proposed filter removes the noise-induced residues thoroughly and preserves the fringe boundaries effectively.

Figure 7c illustrates the rewrapped phase fringe based on the complex average filter. Obviously, area A is blurred, indicating that there is residue-induced erroneous unwrapping. Using a local histogram-based filter, there are similar phenomena at B and C, as shown in Fig. 7d. Each of the three discussed algorithms applies a filtering window; however, the results are quite different. The fundamental reason is that the selected window used in the proposed algorithm varies adaptively according to the local orientation of the fringe. The phases within this window are classified into correct phase groups and noisy phase groups accurately. The noisy phases are removed according to the assumption of unsupervised anomaly detection. On the contrary, although the complex average filter is feasible to filter random noise in fringe images, the selection of window size may impact the preservation of the fringe edges and the accuracy of the PU, particularly for the case of narrow spacing fringes. On the other hand, if the noise is dense, the shape of the local histogram may be distorted and a clear peak may not be distinguished, so that it becomes difficult to decide whether there is a phase jump crossing the window, and the rate of false detection increases consequently. Besides, the selection of threshold σ (see [14]) may impact the filtering effect significantly. If σ is too small, the window will smooth the noise-free phases, and if σ is too big, the noisy phases inducing residues cannot be identified. Thus, the criterion to determine optimal σ must be constructed beforehand.

We also verify the rapidity of path-independent PU based on the proposed filter using real phase signals. The surface of a container is selected as the object, because the corresponding measurement is relatively complex in natural environment. This can verify our filter convictively. Our FPP system is shown in Fig. 8. Four phase-shifted sinusoidal fringe images are projected on the surface of a container, and the corresponding deformed fringes are captured by the camera successively. All the captured images share the same size $1024 \times 768$ , and the corresponding measured area is $450 mm \times 320 mm$ . The actual size of the observed surface is shown in Fig. 9. Figure 10 illustrates the process of 3D reconstruction of a partial surface of the container. The number of residues reduces from 589 to 0 within three iterations. The PU result using our filter is shown in Fig. 10d. Obviously, no error propagation occurs, and the time consumed is only $6.5 s$ . Although the algorithm-driven methods can acquire the same result, the operation for PU is time consuming. From the viewpoint of rapidity, the sequence in ascending order is PUs based on minimum-cost-flow, branch-cut, least-squares, and SL, as shown in Table 1. Table 2 lists the computational complexity of the major procedures for separate approaches. The problem of minimum cost flow is equivalent to the problem of linear programming, and this procedure is the most time consuming. For the branch-cut algorithm, flood filling is used to prevent the integration path from crossing the branch cuts. The least-squares-driven PU is a computationally efficient algorithm. However, the least-squares procedure tends to spread the errors rather than containing them within a limited set of points. In comparison with the above three methods, the most time-consuming procedure is the locating residues for an SL-driven PU. In Table 2, N is the pixel number of an image; T is the iteration number of residue filtering, $T ≪ N$ ; $n_{r}$ is the number of residues in each iteration, $n_{r} ≪ N$ ; and n is the point number within the filtering window, $n ≪ N$ .

To confirm further that the estimated unwrapped phase image is almost correct, we will calibrate the profile measurement system for mapping the unwrapped phases to real world 3D coordinates. In order to avoid an undetermined solution for the problem of nonlinear least squares, two gauge blocks, with uniform top heights of 9.115 and $13.085 mm$ , are employed. Figure 11 shows one projected fringe image on the two gauge blocks. When the Levenberg–Marquard algorithm is used, we should set a suitable initial value to prevent the solution from converging to the local minimum. In Eq. (7), let the content within the brackets equal zero, for $k = 1, 2, \dots, m$ , $m > 11$ , and then we get

{(\begin{matrix} - ϕ_{1} & \dots & - ϕ_{m} \\ - i_{1} & \dots & - i_{m} \\ - ϕ_{1} i_{1} & \dots & - ϕ_{m} i_{m} \\ - j_{1} & \dots & - j_{m} \\ - ϕ_{1} j_{1} & \dots & - ϕ_{m} j_{m} \\ z_{1}^{g} & \dots & z_{m}^{g} \\ ϕ_{1} z_{1}^{g} & \dots & ϕ_{m} z_{m}^{g} \\ i_{1} z_{1}^{g} & \dots & i_{m} z_{m}^{g} \\ ϕ_{1} i_{1} z_{1}^{g} & \dots & ϕ_{m} i_{m} z_{m}^{g} \\ j_{1} z_{1}^{g} & \dots & j_{m} z_{m}^{g} \\ ϕ_{1} j_{1} z_{1}^{g} & \dots & ϕ_{m} j_{m} z_{m}^{g} \end{matrix})}^{T} (\begin{matrix} c_{1} \\ c_{2} \\ c_{3} \\ c_{4} \\ c_{5} \\ d_{0} \\ d_{1} \\ d_{2} \\ d_{3} \\ d_{4} \\ d_{5} \end{matrix}) = (\begin{matrix} 1 \\ 1 \\ 1 \\ 1 \\ 1 \\ 1 \\ 1 \\ 1 \\ 1 \\ 1 \\ 1 \end{matrix}) .

Equation (11) can be described as

Ax = b .

From Eq. (12), we have

x = (A^{T} A)^{- 1} A^{T} b .

Thus, we can estimate an initial

x_{0}

using Eq. (13). Substituting

x_{0}

into the Levenberg–Marquardt algorithm, we obtain

x = (- 2.59 \times 10^{- 3}, 3.03 \times 10^{- 5}, 2.80 \times 10^{- 8}, 1.42 \times 10^{- 3}, 4.39 \times 10^{- 8}, 2.78 \times 10^{- 3}, 4.58 \times 10^{- 6}, 2.23 \times 10^{- 7}, - 5.33 \times 10^{- 10}, 2.88 \times 10^{- 6}, 1.89 \times 10^{- 9})^{T}

. Finally, we retrieve the height map using Eq. (6). The height map is shown in Fig. 12. After repeated measurements, we list the accuracy analysis of the height measurement of the top surface in Table 3:

h_{s}

is the actual value,

h_{\max}

is the maximum measured value,

h_{\min}

is the minimum measured value,

\bar{h}

is the average measured value,

σ_{\max}

is the maximum deviation, and

\bar{σ}

is the average deviation. Consequently, the error of repeated measurements is less than 2%.

6. CONCLUSION

In terms of the characteristics of noise-induced residues, a fast clustering-driven residue filter is proposed, and it is embedded in the arbitrarily arranged FPP system. It uses a set of directional windows to include more similar phases, and it classifies the noisy phases inducing residues and the other correct phases within the filtering window into separate clusters accurately. The noisy phases are removed based on the assumption of anomaly detection. The filter is designed to make the measurement system substitute the residue-filtering-driven PU for the algorithm-driven one. Experimental results demonstrated the filter’s effectiveness in residue reduction, fringe preservation, and time savings. After all the residues are removed, the error of repeated measurements is less than 2%.

ACKNOWLEDGMENTS

The authors would like to thank Dacheng Tao and Wei Bian from Nanyang Technological University for their constructive comments, and Bruno Luong from FOGALE Nanotech for his help with compiling the program of “2D phase unwrapping algorithm based on network flow.”

The work described in this paper is supported by the National Natural Science Foundation of China (NSFC) (grant 60806050), the Knowledge Innovation Program of the Chinese Academy of Sciences (CAS) (KGCX2-YW-156, KGCX2-YW-154), the Key Laboratory of Robotics and Intelligent System of Guangdong Province (2009A060800016), the Shenzhen Technology Project (JC200903160416A), and the Shenzhen Nanshan Research Project (2009016).

Table 1. Comparison of Rapidity between the Proposed Algorithm and Several Famous Phase Unwrappings

View Table | View all tables in this article

Table 2. Comparison of Computational Complexity between the Proposed Algorithm and Several Famous Phase Unwrappings

View Table | View all tables in this article

Table 3. Analysis of the Measured Height Data (Unit: mm)

View Table | View all tables in this article

Fig. 1 Work flow in FPP: (a) projection of sinusoidal fringes, (b) image acquisition, (c) fringe analysis and PU, and (d) phase-to-height conversion.

Download Full Size | PDF

Fig. 2 Optical geometry for fringe analysis: the reference plane $o x y$ , the projection plane $o^{p} x^{p} y^{p}$ , and the imaging plane $o^{c} x^{c} y^{c}$ are arbitrarily arranged; Q represents an arbitrary point on the object; A and B indicate the original fringe point projected at Q and the imaging point of Q respectively; $O_{p}$ and $O_{c}$ denote the lens centers of the projector and the camera, respectively.

Download Full Size | PDF

Fig. 3 (a) $256 \times 256$ phase image and (b) corresponding residue image of (a); the white and black flags denote positive and negative residues, respectively.

Download Full Size | PDF

Fig. 4 Configurations of residues: (a) a couple of adjacent opposite-sign residues in the form of up–down, (b) a couple of adjacent opposite-sign residues in the form of diagonal, (c) a couple of disjoint opposite-sign residues, and (d) two couples of adjacent opposite- sign residues joined together. $g (\cdot)$ denotes one pixel in the adjacent residue area.

Download Full Size | PDF

Fig. 5 This is a clustering example. A local window contains a pair of adjacent residues in the form of the diagonal in a wrapped phase image. The seven phases in this pair of residues are shown by a bigger and bold style. The seven phases can be classified into three groups easily based on the knowledge of clustering. In fact, $- 0.0008$ is the noisy phase inducing this pair of residues.

Download Full Size | PDF

Fig. 6 Sixteen directional windows for the filtering of noisy phases. Only the white pixels are included in the computation.

Download Full Size | PDF

Fig. 7 (a) Corresponding phase image of a set of four $256 \times 256$ fringe images contaminated by random noise, (b) the rewrapped phase fringe based on the proposed filter, (c) the rewrapped phase reconstruction using complex average filter with area A blurred, and (d) the rewrapped phase reconstruction using local histogram based filter with areas B and C blurred.

Download Full Size | PDF

Fig. 8 Schematic diagram of our FPP system. There are four modules: fringe projection, image acquisition, fringe analysis, and high-precision motion.

Download Full Size | PDF

Fig. 9 Actual size of the observed surface of a container.

Download Full Size | PDF

Fig. 10 Process of 3D reconstruction: (a) one fringe image of a partial surface of the container, (b) the wrapped phase image of (a), (c) the residue image corresponding to (b), and (d) the reconstruction result using PU based on the proposed filter.

Download Full Size | PDF

Fig. 11 Projected fringe image on the two gauge blocks.

Download Full Size | PDF

Fig. 12 Height map of a partial surface of the measured container.

Download Full Size | PDF

1. F. Chen, G. M. Brown, and M. Song, “Overview of three- dimensional shape measurement using optical methods,” Opt. Eng. 39, 10–22 (2000). [CrossRef]

2. S. S. Gorthi and P. Rastogi, “Fringe projection techniques: whither we are?” Opt. Lasers Eng. 48, 133–140 (2010). [CrossRef]

3. K. Itoh, “Analysis of the phase unwrapping problem,” Appl. Opt. 21, 2470–2470 (1982). [CrossRef] [PubMed]

4. J. Jiang, J. Cheng, and B. Luong, “Unsupervised-clustering- driven noise-residue filter for phase images,” Appl. Opt. 49, 2143–2150 (2010). [CrossRef] [PubMed]

5. H. A. Zebker and Y. P. Lu, “Phase unwrapping algorithms for radar interferometry: residue-cut, least-squares, and synthesis algorithms,” J. Opt. Soc. Am. A 15, 586–598 (1998). [CrossRef]

6. A. Baldi, “Phase unwrapping by region growing,” Appl. Opt. 42, 2498–2505 (2003). [CrossRef] [PubMed]

7. D. C. Ghiglia and L. A. Romero, “Minimum $L^{p}$ -norm two- dimensional phase unwrapping,” J. Opt. Soc. Am. A 13, 1999–2013 (1996). [CrossRef]

8. M. Costantini, “A novel phase unwrapping method based on network programming,” IEEE Trans. Geosci. Remote Sensing 36, 813–821 (1998). [CrossRef]

9. M. Hubig, S. Suchandt, and N. Adam, “A class of solution- invariant transformations of cost functions for minimum cost flow phase unwrapping,” J. Opt. Soc. Am. A 21, 1975–1987 (2004). [CrossRef]

10. Z. P. Liang, “A model-based method for phase unwrapping,” IEEE Trans. Med. Imaging 15, 893–897 (1996). [CrossRef] [PubMed]

11. J. M. N. Leitão and M. A. T. Figueiredo, “Absolute phase image reconstruction: a stochastic nonlinear filtering approach,” IEEE Trans. Image Process. 7, 868–882 (1998). [CrossRef]

12. G. Nico, G. Palubinskas, and M. Datcu, “Bayesian approaches to phase unwrapping: theoretical study,” IEEE Trans. Signal Process. 48, 2545–2556 (2000). [CrossRef]

13. L. An, Q. S. Xiang, and S. Chavez, “A fast implementation of the minimum spanning tree method for phase unwrapping,” IEEE Trans. Med. Imaging 19, 805–808 (2000). [CrossRef] [PubMed]

14. A. Capanni, L. Pezzati, D. Bertani, M. Cetica, and F. Francini, “Phase-shifting speckle interferometry: a noise reduction filter for phase unwrapping,” Opt. Eng. 36, 2466–2472 (1997). [CrossRef]

15. J. S. Lee, K. P. Papathanassiou, T. L. Ainsworth, M. R. Grunes, and A. Reigber, “A new technique for noise filtering of SAR interferometric phase images,” IEEE Trans. Geosci. Remote Sensing 36, 1173 (1998). [CrossRef]

16. G. Nico, “Noise-residue filtering of interferometric phase images,” J. Opt. Soc. Am. A. 17, 1962–1974 (2000). [CrossRef]

17. G. Fornaro, A. Pauciullo, and E. Sansosti, “Phase difference-based multichannel phase unwrapping,” IEEE Trans. Image Process. 14, 960–972 (2005). [CrossRef] [PubMed]

18. K. M. Qian, W. Gao, and H. Wang, “Windowed Fourier-filtered and quality-guided phase-unwrapping algorithm,” Appl. Opt. 47, 5408 (2008). [CrossRef]

19. H. Du and Z. Y. Wang, “Three-dimensional shape measurement with an arbitrarily arranged fringe projection profilometry system,” Opt. Lett. 32, 2438–2440 (2007). [CrossRef] [PubMed]

20. J. Jiang and J. Cheng, “Noise-residue filtering based on unsupervised clustering for phase unwrapping,” in Proceedings of 5th International Symposium on Visual Computing, M. Z. Nashed, ed. (Springer, 2009), pp. 719–727.

21. S. A. Karout, M. A. Gdeisat, D. R. Burton, and M. J. Lalor, “Residue vector, an approach to branch-cut placement in phase unwrapping: theoretical study,” Appl. Opt. 46, 4712–4727 (2007). [CrossRef] [PubMed]

22. P. Berkhin, Survey of Clustering Data Mining Techniques (Springer, 2002).

23. $c_{\min}$ is the cluster in which the number of members is smallest.

24. V. Chandola, A. Banerjee, and V. Kumar, “Anomaly detection: a survey,” ACM Comput. Surv. 41, 1–58 (2009). [CrossRef]

25. $c_{\max}$ is the cluster in which the number of members is biggest.

26. D. C. Ghiglia and M. D. Pritt, Two-Dimensional Phase Unwrapping: Theory, Algorithms, and Software (Wiley, 1998).

Previous Article Next Article

Cited By

Optica participates in Crossref's Cited-By Linking service. Citing articles from Optica Publishing Group journals and other participating publishers are listed here.

Alert me when this article is cited.

Fig. 1 Work flow in FPP: (a) projection of sinusoidal fringes, (b) image acquisition, (c) fringe analysis and PU, and (d) phase-to-height conversion.

View in Article | Download Full Size | PDF

Fig. 2 Optical geometry for fringe analysis: the reference plane

o x y

, the projection plane

o^{p} x^{p} y^{p}

, and the imaging plane

o^{c} x^{c} y^{c}

are arbitrarily arranged; Q represents an arbitrary point on the object; A and B indicate the original fringe point projected at Q and the imaging point of Q respectively;

O_{p}

and

O_{c}

denote the lens centers of the projector and the camera, respectively.

View in Article | Download Full Size | PDF

Fig. 3 (a)

256 \times 256

phase image and (b) corresponding residue image of (a); the white and black flags denote positive and negative residues, respectively.

View in Article | Download Full Size | PDF

g (\cdot)

denotes one pixel in the adjacent residue area.

View in Article | Download Full Size | PDF

- 0.0008

is the noisy phase inducing this pair of residues.

View in Article | Download Full Size | PDF

Fig. 6 Sixteen directional windows for the filtering of noisy phases. Only the white pixels are included in the computation.

View in Article | Download Full Size | PDF

Fig. 7 (a) Corresponding phase image of a set of four

256 \times 256

fringe images contaminated by random noise, (b) the rewrapped phase fringe based on the proposed filter, (c) the rewrapped phase reconstruction using complex average filter with area A blurred, and (d) the rewrapped phase reconstruction using local histogram based filter with areas B and C blurred.

View in Article | Download Full Size | PDF

Fig. 8 Schematic diagram of our FPP system. There are four modules: fringe projection, image acquisition, fringe analysis, and high-precision motion.

View in Article | Download Full Size | PDF

Fig. 9 Actual size of the observed surface of a container.

View in Article | Download Full Size | PDF

View in Article | Download Full Size | PDF

Fig. 11 Projected fringe image on the two gauge blocks.

View in Article | Download Full Size | PDF

Fig. 12 Height map of a partial surface of the measured container.

View in Article | Download Full Size | PDF

Tables (3)

Table 1 Comparison of Rapidity between the Proposed Algorithm and Several Famous Phase Unwrappings

View Table | View all tables in this article

Table 2 Comparison of Computational Complexity between the Proposed Algorithm and Several Famous Phase Unwrappings

View Table | View all tables in this article

Table 3 Analysis of the Measured Height Data (Unit: mm)

View Table | View all tables in this article

Equations (13)

Equations on this page are rendered with MathJax. Learn more.

ψ = W (ϕ) = ϕ - 2 π 〈 \frac{ϕ}{2 π} 〉,

I (x, y) = A (x, y) + B (x, y) \cos [2 n (x, y) π + ψ (x, y) + φ],

ψ (i, j) = \arctan \frac{\sum_{k = 1}^{N} I_{k} (i, j) \sin (\frac{2 k π}{N})}{\sum_{k = 1}^{N} I_{k} (i, j) \cos (\frac{2 k π}{N})},

ϕ (i, j) = ϕ (1, 1) + \sum_{i^{'} = 1}^{i - 1} Δ_{1} ϕ (i^{'}, 1) + \sum_{j^{'} = 1}^{j - 1} Δ_{2} ϕ (i, j^{'}),

Δ_{1} ϕ (i, j) = ψ_{1} (i, j) = Δ_{1} ψ (i, j) + 2 π n_{1} (i, j), Δ_{2} ϕ (i, j) = ψ_{2} (i, j) = Δ_{2} ψ (i, j) + 2 π n_{2} (i, j),

h = z_{Q} = \frac{1 + c_{1} ϕ_{B} + (c_{2} + c_{3} ϕ_{B}) i_{B} + (c_{4} + c_{5} ϕ_{B}) j_{B}}{d_{o} + d_{1} ϕ_{B} + (d_{2} + d_{3} ϕ_{B}) i_{B} + (d_{4} + d_{5} ϕ_{B}) j_{B}},

E = \sum_{k = 1}^{m} [\frac{1 + c_{1} ϕ_{k} + (c_{2} + c_{3} ϕ_{k}) i_{k} + (c_{4} + c_{5} ϕ_{k}) j_{k}}{d_{o} + d_{1} ϕ_{k} + (d_{2} + d_{3} ϕ_{k}) i_{k} + (d_{4} + d_{5} ϕ_{k}) j_{k}} - z_{k}^{g}]^{2},

Δ (i, j) = (- 1 / 2 π) [ψ_{1} (i, j + 1) - ψ_{1} (i, j) - ψ_{2} (i + 1, j) + ψ_{2} (i, j)],

v = \frac{1}{N} \sum_{i = 1}^{N} (x_{i} - \bar{x})^{2},

\bar{x} = \frac{1}{N} \sum_{i = 1}^{N} x_{i},

{(\begin{matrix} - ϕ_{1} & \dots & - ϕ_{m} \\ - i_{1} & \dots & - i_{m} \\ - ϕ_{1} i_{1} & \dots & - ϕ_{m} i_{m} \\ - j_{1} & \dots & - j_{m} \\ - ϕ_{1} j_{1} & \dots & - ϕ_{m} j_{m} \\ z_{1}^{g} & \dots & z_{m}^{g} \\ ϕ_{1} z_{1}^{g} & \dots & ϕ_{m} z_{m}^{g} \\ i_{1} z_{1}^{g} & \dots & i_{m} z_{m}^{g} \\ ϕ_{1} i_{1} z_{1}^{g} & \dots & ϕ_{m} i_{m} z_{m}^{g} \\ j_{1} z_{1}^{g} & \dots & j_{m} z_{m}^{g} \\ ϕ_{1} j_{1} z_{1}^{g} & \dots & ϕ_{m} j_{m} z_{m}^{g} \end{matrix})}^{T} (\begin{matrix} c_{1} \\ c_{2} \\ c_{3} \\ c_{4} \\ c_{5} \\ d_{0} \\ d_{1} \\ d_{2} \\ d_{3} \\ d_{4} \\ d_{5} \end{matrix}) = (\begin{matrix} 1 \\ 1 \\ 1 \\ 1 \\ 1 \\ 1 \\ 1 \\ 1 \\ 1 \\ 1 \\ 1 \end{matrix}) .

Ax = b .

x = (A^{T} A)^{- 1} A^{T} b .

Algorithm	Computing Time (s)
Branch cut	187
Least squares	13
Min. cost flow	494
SL driven	6.5

Algorithm	Computational Complexity
Branch cut	residue	branch cut	integration
Branch cut	$O (N)$	$O (N)$	$O (N^{2})$
Least squares	discrete cosine transform	inverse discrete cosine transform	integration
Least squares	$O (N \log N)$	$O (N \log N)$	$O (N)$
Min. cost flow	linear programming		integration
Min. cost flow	$O (N^{3})$		$O (N)$
SL driven	residue	window	clustering	integration
SL driven	$O (T N)$	$O (T n_{r} n^{3})$	$O (T n_{r} n^{2})$	$O (N)$

Algorithm	Computing Time (s)
Branch cut	187
Least squares	13
Min. cost flow	494
SL driven	6.5

Algorithm	Computational Complexity
Branch cut	residue	branch cut	integration
Branch cut	$O (N)$	$O (N)$	$O (N^{2})$
Least squares	discrete cosine transform	inverse discrete cosine transform	integration
Least squares	$O (N \log N)$	$O (N \log N)$	$O (N)$
Min. cost flow	linear programming		integration
Min. cost flow	$O (N^{3})$		$O (N)$
SL driven	residue	window	clustering	integration
SL driven	$O (T N)$	$O (T n_{r} n^{3})$	$O (T n_{r} n^{2})$	$O (N)$

Abstract

1. INTRODUCTION

2. PRINCIPLE OF FRINGE PROJECTION PROFILOMETRY

2A. Fringe Analysis

2B. Phase Unwrapping

2C. System Calibration

3. CHARACTERISTICS OF NOISE-INDUCED RESIDUES

3A. Localization of Residues

3B. Configurations of Residues

4. CLUSTERING-DRIVEN RESIDUE FILTER

4A. Motivations

4B. Identification of the Filtering Window

4C. Unsupervised Clustering

4D. Proposed Filtering Algorithm

5. EXPERIMENTS

6. CONCLUSION

ACKNOWLEDGMENTS

Cited By

Figures (12)

Tables (3)

Equations (13)

Journal of the Optical Society of America A