Exponential fringe pattern projection approach to gamma-independent phase computation without calibration for gamma nonlinearity in 3D optical metrology

This paper presents a method that expresses the fringe pattern as an exponential function and a mathematical model for gamma-independent phase computation. The method was compared to: (i) conventional phase measurement without nonlinearity correction, and (ii) conventional gamma correction by pattern pre-distortion based on an input-to-projector camera-output look-up table. The pre-distorted and exponential methods achieved large reduction in error compared to conventional computation with no gamma correction. The advantage of the exponential method is that no system gamma nonlinearity calibration procedure or information is required. This reduces optical system setup before measurement and permits easier use of off-the-shelf projectors.


Introduction
Digital Fringe Projection (DFP) is a commonly used method for three-dimensional (3D) surface shape measurement. DFP systems essentially consist of a projector that projects fringe patterns onto an object surface and a camera that captures images of the deformed patterns on the object surface. In phase-shifting methods, multiple phase-shifted patterns and captured images are used, and a phase-shifting algorithm is employed to compute a phase map from the captured images. This phase map contains the depth information of the object, from which the surface shape is calculated, after first applying a phase unwrapping algorithm [1] [2] to produce an unwrapped phase map.
Among the various sources of error in DFP measurement systems, projector nonlinearity, commonly referred to as gamma nonlinearity, has been known to be a major source of measurement error [3]. Projector gamma is the nonlinear response to the grey values input to the projector, purposely incorporated in off-the-shelf digital data projectors as a visual enhancement for the nonlinear sensitivity of human vision to intensity. The nonlinear mapping of the projector input to captured-image intensity causes errors in the computed phase map which leads to errors in the 3D surface height map. Correction methods have been developed to either: correct the captured image intensities and then compute the phase with the corrected images, correct the phase computed from the captured images with distorted intensities, or recover the depth from distorted phase without correcting phase or intensity.
Correction of image intensities, either after image capture (passive approach) or before pattern projection (active approach), are common gamma correction methods [4]. The projector-input -camera-output relationship has been determined by inputting uniform gray level patterns to the projector, projecting the patterns onto a white flat plate, capturing the images as output, and determining the input-output relationship [5][6][7] either as a gamma function [8], polynomial [9], constrained cubic spline [10], or look-up-table (LUT) [11,12].
Once the gamma curve is determined prior to measurement, an inverse gamma function can be applied to the gray levels before input to the projector [13,14] by first pre-distorting (correcting) the projected pattern using a projector inverse intensity transfer function [11,15], using a combined projector-camera nonlinearity function [8], or superimposing specific harmonics to adjust pattern amplitude [16]. While the above methods are effective in reducing errors in the computed phase and height, they need prior input-output mapping that requires the projection and capture of numerous patterns and images, respectively.
Rather than determining the projector-camera intensity mapping and correcting the intensities, other methods use the gamma-distorted intensities, but correct the phase computed from the captured images. Phase error compensation has been applied using a precomputed LUT of stored phase errors [17,18]. This compensation has been effective in the three-step phase shifting algorithm [18] and is generic for any phase-shifting method [3]. While these methods are useful, these approaches typically require prior knowledge of the gamma nonlinearity function or ideal phase map based on projection and capture of numerous patterns and images, in order to apply the correction.
Other methods of handling system gamma nonlinearity do not use direct intensity or phase mapping. Gamma error has been reduced by eliminating the higher order harmonics of the computed wrapped phase map [19][20][21]. Based on the period of phase residuals being equal to the number of phase shifts of the captured fringe patterns, a mathematical model that iteratively eliminates the phase error has been used [22]. In a similar approach, a cubic spline was utilized to estimate the ideal phase from the computed phase using a smoothing parameter based on a planarity constraint on the reference phase [23]. This method is suitable for continuous surfaces. In another approach, phase is computed without an explicit gamma mapping procedure. The gamma and phase were computed in a minimization of the energy of the higher order harmonics of the captured fringe patterns for every pixel [24]. Based on [24], the phase error amplitude was estimated to compensate for the phase map error using both active and passive methods [25]. While the results of these methods were excellent, a high number of phase-shifted patterns are required for calibration. A comparison of the cumulative intensity distribution functions of captured fringe images and of an ideal sinusoidal pattern [26], and matching the histogram of captured fringe patterns with those of the ideal fringe pattern [27] have been used to estimate the gamma value. Blind gamma correction without requiring calibration [28] uses bicoherence as a measure to determine the correlation between different higher order harmonics, and bicoherence minimization to solve for gamma. Gamma nonlinearity compensation without a specific nonlinearity model or calibration data involved optimization by minimizing a ratio of high to low frequencies [29]. While the above methods have worked well, they involve iterative processes or high computational cost to correct for phase error.
Still other methods of handling system gamma nonlinearity recover depth from gammadistorted phase without correcting phase or intensity. As an example, the Generalized Analysis Model (GAM) treats the deformed fringe pattern on the object surface as a shifted version of the pattern on the reference plane, where an Inverse Function Shift Estimation (IFSE) retrieves the object depth directly from the spatial shift, independent of projector gamma [30]. This method requires a prior calibration process to determine the shift-to-depth conversion.
The above techniques of handling system gamma nonlinearity either require a prior mapping such as projector-input to camera-output, additional information (such as an ideal phase map) that often leads to greater acquisition time, iterative processes, or high computational cost. This paper presents a novel passive method of handling system gamma nonlinearity that expresses the fringe pattern as an exponential function, and develops a mathematical model for gamma-independent phase computation. The method greatly reduces the effect of projector nonlinearity without requiring any prior or post knowledge of the projector's nonlinear behavior, and without involving iterative processes or high computational cost. This greatly reduces optical system setup and preparation before measurement, permitting easier use of off-the-shelf consumer projectors.

Conventional fringe pattern projection
A fringe pattern is conventionally generated by a cosine function with the desired pitch and amplitude. A vertical fringe pattern in phase shifting DFP is expressed as follows: − the fringe pattern index, N the number of phase shifts, p indicates projector, f is the spatial frequency (inverse of pitch) of the fringe pattern, p B x y is the intensity modulation (amplitude), and (x,y) are the pattern coordinates (column, row).
The nonlinear input-to-projector to camera-output mapping can be modelled by a gamma function. After pattern projection and image acquisition, the captured image can thus be modeled as follows: where γ is the gamma function, the background and amplitude parameters are similar to those in Eq. (1), c indicates camera, φ is the wrapped phase between -π and π, and (x,y) are the image coordinates. In the case of four phase shifts, N = 4, and the wrapped phase map ( , ) x y φ is calculated by:

Exponential fringe pattern projection
The novel approach of this paper is to take advantage of the nonlinear behaviour of the projector expressed as a gamma function, and restate the formulation of fringe patterns using an exponential function: To include the effect of projector gamma nonlinearity, Eq. (4) is modified as follows: A plot of conventional and exponential fringe patterns based on Eqs. (1) and (4), respectively, and simulated camera-captured fringe images: conventional with gamma distortion and exponential with gamma distortion, based on Eqs. (2) and (5), respectively, are shown for one pattern column (Fig. 1) to illustrate the difference in intensity between the conventional and exponential patterns, and the effect of gamma distortion, where gamma γ is set to 2.2.
Equation (6) can be expressed as: where ( ) For phase computation based on the four-step phaseshifting algorithm, the four phase-shifted patterns are: These equations can be simplified as follows:

Simulation
To investigate the effectiveness of the exponential fringe-pattern projection method, a simulation of fringe projection on a flat plate with image capture was performed. First, a set of conventional and a set of exponential phase-shifted fringe patterns (both four step with equal phase shifts, and 16-pixel pitch) were generated. Then, the fringe patterns of both sets (conventional and exponential) were distorted by a gamma function acquired by projector gamma calibration [9]. Afterward, white noise with zero mean and 1% variance in gray values were added to both fringe pattern sets (conventional and exponential). This was done using eleven different amounts of added noise, from 0% to 100% of image pixels, in increments of 10%. Finally, the phase maps using conventional and exponential patterns were calculated using Eqs. (3) and (19), respectively, and both were compared to the ideal (undistorted noise-free) phase map. The root mean square error (RMSE) of the phase for both image sets (conventional and exponential) uncorrected for gamma nonlinearity for various amounts of added noise is shown in Table 1. The uncorrected conventional fringe patterns (images) had errors of approximately 0.04 rad for all amounts of noise added, while errors for the exponential fringe patterns (images) were negligible at 0%, and 0.003 to 0.015 rad from 10% to 100% noise added. The uncorrected exponential fringe pattern (image) errors were thus much lower than for the uncorrected conventional fringe patterns (images), from nearly three times lower at 100% added noise, to several orders of magnitude lower at 0% added noise. The phase computed with exponential fringe images was more sensitive to image noise, compared to the phase computed with conventional fringe images; however, sensitivity is a less dominant factor, since the highest error for exponential fringe images, was lower than the smallest error for conventional fringe images.

Implementation in real phase measurement
The method of phase computation based on exponential fringe patterns was implemented to demonstrate its effectiveness in a real optical measurement environment. The optical setup consisted of a LCD projector (Panasonic PT-AE7000U) and monochrome CCD camera (Basler avA1000-100gm, 1024 x 1024 resolution), with approximately 10 deg cameraprojector angle and 2.5 m camera-projector-plane to object distance. Three sets of phase-shifted fringe patterns (all four-step and equal phase shifts) were used for experiments: (i) conventional, (ii) conventional pre-distorted (corrected) by an input-toprojector camera-output look-up table, determined from prior gamma calibration, and (iii) exponential. In three separate tests, the three sets of patterns were generated with pitches ranging from 16 to 2048 pixels, using 2 b , b = [4,11], and then projected onto a flat plate mounted on a translation stage that was translated to 21 known positions (depths) from 0 to 100 mm in 5 mm intervals. For each projected pattern (3 sets, 4 phase shifts, 8 pitches, 21 plate positions), 10 images were captured and averaged for noise suppression [31]. For each plate position, a wrapped phase map was calculated for each fringe pitch using: (i) images from the conventional patterns and Eq. (3), (ii) images from pre-distorted (corrected) conventional patterns and Eq. (3), and (iii) images from the exponential patterns and Eq. (19), and multiple wavelength phase unwrapping [32] was performed using the eight wrapped phase maps (8 pitches) to generate one unwrapped phase map for each set (i, ii, iii). A total of 63 experiments (3 sets x 21 plate positions) were thus performed.
The gamma function for the optical measurement system and the inverse gamma function (Fig. 3) show one-to-one input to output mapping in the intensity range 45 to 225. This range was thus used as the input intensity range for generating the pre-distorted (corrected) patterns. A comparison of the camera-captured images of the conventional, pre-distorted (corrected) and exponential fringe patterns are shown in Fig. 4.
The phase maps computed from conventional fringe patterns, and exponential fringe patterns before and after outlier removal are shown in Fig. 5. The salt-and-pepper noise in the exponential phase map (Fig. 5(b)) (mainly in the bottom left and top right corners) is due to ambiguity in Eqs. (16) and (17), as discussed earlier (spikes in phase in Fig. 2). While the error occurs in all columns based on simulated input data (not shown in figures), the error does not occur in each column of the phase map based on the real images, possibly due to camera intensity quantization and some projector defocus. Outlier removal based on value deviation from the mean of all data was performed to obtain the result in Fig. 5(c).
To evaluate the exponential fringe method, phase errors were calculated for the exponential fringe method and compared to the errors for the conventional (uncorrected) and pre-distorted (corrected) fringe methods. The phase error at each pixel was defined as the difference between the calculated phase and expected ideal (linear) phase. The proposed exponential fringe method performed nearly as well as the pre-distorted (corrected) fringe method, both of which had greatly reduced phase errors compared to conventional phase calculation without gamma correction. This is shown in the plot of phase errors of one column of the phase maps for conventional, pre-distorted (corrected), and exponential methods (Fig. 6).     The phase RMSE for the region-of-interest (ROI) 400x400 pixels, centred in the image, was computed for each of the 21 positions of the flat plate over the range 0 to 100 mm with 5 mm intervals as follows. For each of the three unwrapped phase maps (conventional, predistorted (corrected), and exponential), a plane was fit to the unwrapped phase map and considered as the ideal phase map. The phase RMSE at each position (Fig. 7) was based on the differences between the phase at all pixels in the ROI for each phase map and the ideal phase map. The ROI was the same as that used for the gamma calibration and pre-distorting patterns. Over the 100 mm depth range, the phase RMSE varied from 0.0418 to 0.0602 rad for conventional, 0.0097 to 0.0093 rad for pre-distorted (corrected), and 0.0117 to 0.0115 rad for exponential methods. The mean phase RMSE and standard deviation (SD) of phase RMSE for conventional, pre-distorted (corrected), and exponential methods over all plate positions are shown in Table 2. Both the pre-distorted (corrected) and exponential methods achieved large reduction in error (81% and 77%, respectively) compared to the conventional phase computation with no gamma correction. The difference between the new exponential method and the pre-distorted (corrected) method was small (0.0019 rad). The small error difference may be partly due to the assumption of a gamma function in Eq. (5) for the exponential method. The exponential method might be more sensitive to noise than the predistorted (corrected) method because of the lower fringe-pattern intensity of the exponential method (Fig. 4). The standard deviation based on multiple measurements was low for the predistorted and exponential methods compared to the conventional method. Phase maps for a measured mask show no observable difference in phase between exponential and predistorted (corrected) methods (Fig. 8). Both phase maps are smooth while the conventional method has ripples. The advantage of the proposed exponential method is that no gamma calibration procedure is required and no information of the projector and camera nonlinearity is required. This greatly reduces optical system setup and preparation before measurement, permitting easier use of off-the-shelf consumer projectors. Furthermore, phase is computed without involving iterative processes or high computational cost. The exponential method has been demonstrated for four-step phase shifting. Application to other fringe steps requires further investigation.

Conclusion
The method of exponential fringe-pattern projection greatly reduced the gamma nonlinearity compared to the conventional method with no nonlinearity correction, and the new method performed nearly as well as the pre-distorted (corrected) fringe method using a LUT. While the exponential fringe-pattern projection method had slightly higher phase error than the predistorted (corrected) method, the exponential method has the advantage of not requiring any prior knowledge of the camera-projector nonlinear gamma behaviour, and thus does not require any optical-system gamma calibration procedure.

Funding
This research has been funded by the Ministry of Science, Research and Technology of I.R. Iran, and the Natural Sciences and Engineering Research Council of Canada (NSERC).