The deflectometry technique is used for a surface measurement where the local slopes are optically measured and the surface is reconstructed using an integration procedure. This approach has several advantages over a direct height measurement.^{1, 2} The resultant measurement of this technique is the directional derivative of a wavefront *W*, which is defined by

## 1

[TeX:] \documentclass[12pt]{minimal}\begin{document} \begin{eqnarray} {\phi}^{k}({\bf r}) = \triangledown W \cdot {\bf v}_{k} {+} \eta^{k} = \cos \varphi_{k} \frac{\partial W( {\bf r})}{\partial x} + \sin \varphi_{k} \frac{\partial W({\bf r})}{\partial y} + {\eta}^{k},\!\!\!\!\!\!\nonumber\\ \end{eqnarray} \end{document} $$\begin{array}{c}\hfill {\phi}^{k}\left(\mathbf{r}\right)=\u25bfW\xb7{\mathbf{v}}_{k}+{\eta}^{k}=\mathrm{cos}{\varphi}_{k}\frac{\partial W\left(\mathbf{r}\right)}{\partial x}+\mathrm{sin}{\varphi}_{k}\frac{\partial W\left(\mathbf{r}\right)}{\partial y}+{\eta}^{k},\end{array}$$**r**= (

*x*,

*y*) is a position in an

*M*×

*N*rectangular grid of pixels,

**v**

_{k}= (cos φ, sin φ) is the measurement direction of a sensor,

*k*= 1, …,

*K*with

*K*⩾ 2, and η

^{k}represents the measurement and processing errors.

A common way to integrate the derivatives is to use a global integration method,^{3, 4} where most of the integration techniques only use two orthogonal directions in the measurement, that is ϕ^{1} and ϕ^{2} with **v**_{1} = (1, 0) and **v**_{2} = (0, 1). Recently, in Ref. 5 a procedure was proposed that uses multiple directional derivatives, which reduces the noise contribution in the reconstruction of the wavefront. This consists of acquiring *k*-number of directional derivatives with different measurement direction, and the estimated wavefront *W*_{e} is computed as the minimizer of the following cost functional:

## 2

[TeX:] \documentclass[12pt]{minimal}\begin{document} \begin{eqnarray} U(W) &=& \sum _{k}\sum _{(x,y)\in L} \lbrace \left[W(x+1,y) -W(x,y)\right]\cos \varphi _{k}\nonumber \\[-5pt] && + [W(x,y+1) -W(x,y)]\sin \varphi _{k} - \phi ^{k}(x,y) \rbrace ^{2},\nonumber\\ \end{eqnarray} \end{document} $$\begin{array}{ccc}\hfill U\left(W\right)& =& \sum _{k}\sum _{(x,y)\in L}\{\left[W(x+1,y)-W(x,y)\right]\mathrm{cos}{\varphi}_{k}\hfill \\ & & +[W(x,y+1)-W(x,y)]\mathrm{sin}{\varphi}_{k}-{\phi}^{k}{(x,y)\}}^{2},\hfill \end{array}$$*L*defines the region with valid data.

One drawback of the above cost functional is the computational time employed to found the minimizer of the functional. An alternative to reduce the processing time is to perform the integration through faster methods, that is to perform the integration in the frequency domain. To our knowledge, in the literature there are only two very similar works related to the integration of directional derivatives in frequency domain;^{6, 7} both methods search an approximation between the measurement data and the directional derivative model in the least-squares sense. One problem with this least-squares approach is its sensitivity to a high level of noise, even when these methods use a large number of derivatives. This situation is illustrated in Fig. 1: (a) shows a synthetic wavefront used as reference, (b) shows the wavefront reconstructed with the method reported in Ref. 5, and (c) shows the wavefront reconstructed with the method reported in Ref. 7. In all cases, a Gaussian noise with a signal-to-noise ratio equal to −20 db was added to 10 directional derivatives employed in the reconstructions. As one can see, the reconstructed wavefront shows artifacts due to noise despite the large number of derivatives utilized on the methods.

To overcome this drawback and compute a fast reconstruction, we add to the cost functional given in Eq. 2, a term which assumes that the wavefront is continuous.^{8} Thus, the wavefront *W*_{e} is estimated as the minimizer of our proposed regularized cost functional

## 3

[TeX:] \documentclass[12pt]{minimal}\begin{document} \begin{eqnarray} J(W) &=& \sum _{k}\int\!\!\int |\triangledown W\cdot {\bf v}^{k} - \phi ^{k} |^{2} {\rm d}x\,{\rm d}y \nonumber \\ && + \lambda \int\!\!\int (|W_{xx}|^{2} + |W_{xy}|^{2} + |W_{yy}|^{2} )\, {\rm d}x\,{\rm d}y, \end{eqnarray} \end{document} $$\begin{array}{ccc}\hfill J\left(W\right)& =& \sum _{k}\int \int {|\u25bfW\xb7{\mathbf{v}}^{k}-{\phi}^{k}|}^{2}\mathrm{d}x\phantom{\rule{0.16em}{0ex}}\mathrm{d}y\hfill \\ & & +\lambda \int \int \left(\right|{W}_{xx}{|}^{2}+|{W}_{xy}{|}^{2}+|{W}_{yy}{|}^{2})\phantom{\rule{0.16em}{0ex}}\mathrm{d}x\phantom{\rule{0.16em}{0ex}}\mathrm{d}y,\hfill \end{array}$$^{9}

The minimization of the above cost function is performed by means of Fourier transform theory:^{10} Taking the Fourier transform and using its properties, we obtain the minimizer of Eq. 3 expressed as

## 4

[TeX:] \documentclass[12pt]{minimal}\begin{document} \begin{equation} {\mathfrak F}\lbrace W_{e}\rbrace = {\displaystyle {i\;\sum \limits _{k} [ {\mathfrak F}\lbrace \phi ^{k}\rbrace ({\bf q}\cdot {\bf v}_{k})] \over \sum \limits _{k}[({\bf q}\cdot {\bf v}_{k})^{2} ] + \lambda \;||{\bf q}||^{2}}}, \end{equation} \end{document} $$\mathfrak{F}\left\{{W}_{e}\right\}={\displaystyle \frac{i\phantom{\rule{0.28em}{0ex}}\sum _{k}\left[\mathfrak{F}\left\{{\phi}^{k}\right\}(\mathbf{q}\xb7{\mathbf{v}}_{k})\right]}{\sum _{k}\left[{(\mathbf{q}\xb7{\mathbf{v}}_{k})}^{2}\right]+{\lambda \phantom{\rule{0.28em}{0ex}}\left|\right|\mathbf{q}\left|\right|}^{2}}},$$**q**= (

*u*,

*v*) is the position vector on frequency domain, and [TeX:] $ {\mathfrak F}\lbrace \cdot \rbrace$ $\mathfrak{F}\{\xb7\}$ denotes the Fourier transform. This expression becomes the same reported in Ref. 10 for the case of

*K*= 2 with

**v**

_{1}= (1, 0) and

**v**

_{2}= (0, 1), and for the case of λ = 0, Eq. 4 becomes the same expression reported in Ref. 7. For practical purposes, one has to consider that the discrete Fourier transform always needs full fields of valid data. The missing data, for instance due to occlusions in the measurement, has to be extrapolated by use of a Gershberg-type algorithm or an alternative technique.

^{11, 12}

The performance of the proposed technique is illustrated by two experiments. All numerical experiments were made in a 2.53 GHz Pentium Dual Core PC with 8 GB of main memory using Ubuntu 9.04 as the operative system. The algorithm used to compute the discrete Fourier transform was the FFTW library.^{13} The first experiment was a numerical estimation using the synthetic wavefront shown in Fig. 1. As mentioned above, a Gaussian noise with a signal-to-noise ratio equal to −20 db was added to 10 directional derivatives, which were generated by use of Eq. 1; the direction **v**_{k} was equally spaced from 0° to 180°, where every derivative is a 512 × 512 matrix. The resultant wavefront is shown in Fig. 1 using λ = 20. As one can observe, the artifacts in the reconstructed wavefront are notably reduced. The mean square error (MSE) obtained in each reconstruction is shown on Table 1.

## Table 1

Mean square error for the reconstructed wavefronts.

Method | MSE | Time employed (ms) |
---|---|---|

Ref. 5 | 0.243157 | 9011 |

Ref. 7 | 0.320633 | 472 |

Eq. 4, λ = 20 | 0.147356 | 472 |

The second experiment was the reconstruction of a progressive lens surface. An example of the fringe pattern used in the experiment is shown in Fig. 2. Five directional derivatives were acquired using Eq. 1 with **v**_{k} equally spaced from 0° to 180° with steps of 30°, where every image has a resolution of 2048 × 1552 pixels. The fringe patterns were processed using the procedure reported in Ref. 14. The computational time employed for the reconstruction was nearly 879 ms using the proposed technique with λ = 5. In the case of the procedure of Ref. 5, the time employed was 272,275 ms. The resultant reconstructions are shown in Figs. 2 and 2. As one can observe, both reconstructions are similar where the uncertainty of both measurements are about 0.04 mm. However, there is a huge difference in the computational process employed for the reconstruction.

The performance of the cost functional proposed here, Eq. 3 and their minimizer given in Eq. 4, was proved with two fast and accurate reconstructions using a large number of directional derivatives, as it was shown in the results described above. One extra advantage of the proposed method is its feasibility to be implemented on a dedicated hardware for processing in real-time, which is the aim of future research.

## Acknowledgments

The authors were supported in part by grants of the Consejo Nacional de Ciencia y Tecnologia (CONACyT).