Expand this Topic clickable element to expand a topic
Skip to content
Optica Publishing Group

Efficient vision ray calibration of multi-camera systems

Open Access Open Access

Abstract

Vision ray calibration provides imaging properties of cameras for application in optical metrology by identifying an independent vision ray for each sensor pixel. Due to this generic description of imaging properties, setups of multiple cameras can be considered as one imaging device. This enables holistic calibration of such setups with the same algorithm that is used for the calibration of a single camera. Obtaining reference points for the calculation of independent vision rays requires knowledge of the parameters of the calibration setup. This is achieved by numerical optimization which comes with high computational effort due to the large amount of calibration data. Using the collinearity of reference points corresponding to individual sensor pixels as the measure of accuracy of system parameters, we derived a cost function that does not require explicit calculation of vision rays. We analytically derived formulae for gradient and Hessian matrix of this cost function to improve computational efficiency of vision ray calibration. Fringe projection measurements using a holistically vision ray calibrated system of two cameras demonstrate the effectiveness of our approach. To the best of our knowledge, neither any explicit description of vision ray calibration calculations nor the application of vision ray calibration in holistic camera system calibration can be found in literature.

© 2021 Optical Society of America under the terms of the OSA Open Access Publishing Agreement

1. Introduction

Optical metrology techniques based on the utilization of incoherent light such as the fringe projection technique (FPT) [1] and phase measuring deflectometry (PMD) [2,3] utilize cameras as sensors for three-dimensional topography measurement. Absolute measurement requires camera calibration to establish a metric relation between the image plane and real-world coordinates [4,5]. Hence, accuracy of measured quantities depends directly on the quality of the camera calibration. Due to its inherent importance not only in optical metrology but in computer vision in general, the topic of camera calibration has been widely discussed in literature [49]. Many applications obtain depth utilizing information from different perspectives provided by stereo- or multi-camera systems. FPT with stereo cameras does not require a calibration of the projector unit since height information is determined by triangulation of camera vision rays [10]. In PMD, stereo vision is often applied to solve the inherent ambiguity between height and slope information [11]. In addition to the calibration of imaging properties of the individual cameras, such applications also require precise knowledge of relative positions and orientations of the imaging devices.

Most calibration techniques utilize the photogrammetric camera model, assuming a pinhole camera and accounting for distortions by means of polynomial functions [47]. Knowledge of the parameters of the photogrammetric model provides a vision ray, pointing towards the source of received light, for each point on the image plane. Model parameters are typically determined from images of a known calibration target captured from multiple perspectives. Since a limited number of reference points suffices for calibration, classical marker plates are commonly used for this purpose. However, this technique can only account for image distortion with low spatial frequencies. Furthermore, the assumption of a common point of all vision rays, the camera pinhole, states an avoidable approximation and also renders the photogrammetric calibration inapplicable to certain imaging systems such as fisheye cameras [12].

An alternative approach for camera calibration is to avoid any model assumptions by determination of an independent vision ray for each sensor pixel [1316]. Hence, this generic camera calibration also determines strongly localized image distortions [16,17]. Furthermore, since there are no inherent specifications regarding sensor pixel positions and vision ray orientations, setups of multiple cameras can be considered as a single imaging system, enabling holistic calibration without the introduction of additional requirements. Determination of independent vision rays necessitates a calibration procedure providing reference points for each sensor pixel. This requires either interpolation of reference points or application of a calibration target providing continuous spatial coding of its surface.

Vision ray calibration (VRC) utilizes liquid crystal displays (LCDs) as targets for this model-free camera calibration [16]. Application of phase shifting technique (PST) [18], evaluating images of sinusoidal fringe patterns shifted by constant phase angles, provides continuous spatial coding of the surface of LCDs. Using FPT for qualification, it was demonstrated that VRC can achieve superior accuracy compared to photogrammetric camera calibration due to providing a truer description of imaging properties [17,19]. Further improvements were accomplished by the introduction of an extended model of the display utilized as calibration target [20]. However, the inherent advantage of VRC being able to holistically calibrate multiple cameras at once, has not yet been properly capitalized on. Furthermore, we found a lack of explicit documentation regarding the implementation of VRC. This seems adverse since the calculation of the required numerical optimization can be very time-consuming if the mathematical problem at hand is posed unfavorably.

With this work, we contribute a detailed mathematical description of VRC and present a cost function as well as analytically derived formulae for its gradient and Hessian matrix. This provides the means to perform the required numerical optimization for the determination of system parameters with high computational efficiency.

This paper is structured as follows: We will first introduce the quantities obtained by calibration measurements as well as those describing the setup. Based on these quantities, we derive a cost function for numerical optimization of system parameters that does not require explicit calculation of vision rays. After demonstrating the calculation of derivatives of reference points with respect to individual system parameters, we present the analytical derivations of gradient and Hessian matrix of our cost function. The subsequent experimental section demonstrates the application of our holistic VRC and provides FPT measurements for qualification of accuracy.

2. Theory

2.1 Local reference points

Utilizing PST [18], an LCD serves as target for our calibration. Evaluation of multiple images of phase shifted sinusoidal fringe patterns yields unique spatial coding of the LCD surface and constitutes a calibration measurement. We perform $N$ calibration measurements using different relative positions and orientations between our imaging system with $M$ sensor pixels and the calibration target. For conciseness, we assume that each sensor pixel $m$ observes a point on the surface of the calibration target in each measurement $n$.

Since the surface of our calibration target provides continuous spatial coding, we receive $N$ two-dimensional reference points $(x_{m,n}^l,y_{m,n}^l )^T$ for each sensor pixel $m$. The superscript $l$ refers to these coordinates being defined within the local coordinate system of the calibration target. Knowledge of the shape of the calibration target provides corresponding $z_{m,n}^l$-values given as a function of $x_{m,n}^l$ and $y_{m,n}^l$. Let us assume an approximately plane calibration target with small flatness deviations. We can utilize a bivariate polynomial function for surface parameterization by defining

$$z_{m,n}^l = \sum_{p,q} c_{pq} \cdot (x_{m,n}^l)^p (y_{m,n}^l)^q.\\$$

The polynomial coefficients $c_{pq}$ correspond to the powers $p,q$ of $x,y$-coordinates. Such a parameterization is advisable when using LCDs in optical metrology since these typically exhibit flatness deviation in the range of about 1 mm [21] which will affect the accuracy of the calibration if ignored. It has been demonstrated that such a surface parameterization of LCD significantly reduces systematic errors in FTP [20] and PMD [22]. The selection of polynomial terms for surface parameterization is a choice of the user. For conciseness of later notations, we introduce the surface parameter vector $\vec {C}$ containing all polynomial coefficients $c_{pq}$ corresponding to the polynomial terms of choice.

2.2 Global reference points

For each measurement $n$ exists a transformation of local coordinates $(x_{m,n}^l,y_{m,n}^l,z_{m,n}^l)^T$ into global coordinates $(x_{m,n},y_{m,n},z_{m,n})^T$ given by

$$\begin{pmatrix} x_{m,n} \\ y_{m,n} \\ z_{m,n} \end{pmatrix} = R^{(z)}(\gamma_n) R^{(y)}(\beta_n) R^{(x)}(\alpha_n) \begin{pmatrix} x_{m,n}^l \\ y_{m,n}^l \\ z_{m,n}^l \end{pmatrix} + \begin{pmatrix} t_n^{(x)} \\ t_n^{(y)} \\ t_n^{(z)} \end{pmatrix}.$$

Here, $R^{(x)}$, $R^{(y)}$ and $R^{(z)}$ are common 3-by-3 rotation matrices with respect to the coordinate axis marked in the superscript and $\alpha _n$, $\beta _n$ and $\gamma _n$ are the rotation angles corresponding to a certain measurement $n$. Note that other conventions regarding the order of application of rotation matrices are as well possible. Translations are denoted as $t_n^{(x)}$, $t_n^{(y)}$ and $t_n^{(z)}$ with respect to the coordinate axis marked in the superscript. Therefore, a set of six parameters $\vec {T}_n = (\alpha _n,\beta _n,\gamma _n,t_n^{(x)},t_n^{(y)},t_n^{(z)})$ describes the pose of the calibration target for each measurement $n$.

2.3 Vision rays

Assuming knowledge of pose parameters $\vec {T}_n$ for all measurements $n$ and of surface parameters $\vec {C}$, there exists a known point-cloud $(x_{m,n},y_{m,n},z_{m,n})^T$ for each pixel $m$ of the imaging system. This enables the determination of independent visions rays by linear regression. Introducing the convention that start points lie in the plane $z=0$ of the global coordinate system, each vision ray

$$\vec{r}_m(z) = \begin{pmatrix} x_{m,0} \\ y_{m,0} \\ 0 \end{pmatrix} + \begin{pmatrix} u_m \\ v_m \\ 1 \end{pmatrix} \cdot z; \quad z \geq 0$$
can be described by four parameters: $u_m$ and $v_m$ are the slopes in $x$- and $y$-direction of the vision ray of pixel $m$ of the imaging system while $(x_{m,0},y_{m,0},0)^T$ describes the corresponding startpoint. Determination of these parameters for all pixels yields the characterization of the imaging system and constitutes the calibration result. Figure 1 illustrates the relation of local reference points $(x_{m,n}^l,y_{m,n}^l )^T$, global reference points $(x_{m,n},y_{m,n},z_{m,n} )^T$ and vision rays given by pose parameters $\vec {T}_n$ and surface parameters $\vec {C}$ of the calibration setup.

 figure: Fig. 1.

Fig. 1. Schematic representation of VRC quantities using three calibration measurements. An LCD acts as calibration target and is observed from varying perspectives. (a) Three reference points (denoted green, red and blue) corresponding to an arbitrary pixel $m$ of the imaging system are shown. These reference points are the positions on the LCD surface that were observed by pixel $m$ in these measurements. Using surface parameters $\vec {C}$ and pose parameters $\vec {T}_n$ describing the shape of the LCD surface and its relative pose towards the imaging system for every measurement enable transformation into a global coordinate system. (b) A line fit through the global reference points yields the vision ray $\vec {r}_m$ corresponding to the pixel $m$.

Download Full Size | PDF

2.4 Cost function

For each sensor pixel $m$, we define the error value

$$\Delta_m = \sum_n \left|\left| \begin{pmatrix} x_{m,n} \\ y_{m,n} \\ z_{m,n} \end{pmatrix} - \vec{r}_m(z_{m,n}) \right|\right|^2 = \sum_n \Delta_{m,n}^{(x)^2} + \sum_n \Delta_{m,n}^{(y)^2}$$
as the sum of squared absolute values of deviations between reference points in global coordinates $(x_{m,n},y_{m,n},z_{m,n})^T$ and the intersection points of our vision rays $\vec {r}_m$ with the $x$-$y$-plane at $z_{m,n}$ over all phase measurements $n$. Due to this choice of definition, only $x$- and $y$ components contribute to the error value $\Delta _m$ since the $z$-component of $\vec {r}_m(z_{m,n})$ is equal to $z_{m,n}$ as can be seen in Eq. (3). This allows the introduction of $x$-deviations $\Delta _{m,n}^{(x)}$ and $y$-deviations $\Delta _{m,n}^{(y)}$ of vision rays of pixels $m$ towards individual reference points given by
$$\Delta_{m,n}^{(x)} = x_{m,n} - u_m z_{m,n} - x_{m,0}, \quad \Delta_{m,n}^{(y)} = y_{m,n} - v_m z_{m,n} - y_{m,0}.$$

Figure 2 illustrates the relations between reference points (circles), the fitted vision ray (blue line) and the vision ray deviations (black lines) given in Eq. (5). As shown by Draper and Smith [23], one can eliminate start points in Eq. (5) by using centered coordinates

$$\hat{x}_{m,n} = {x}_{m,n} - \langle x_m \rangle_n, \quad \hat{y}_{m,n} = {y}_{m,n} - \langle y_m \rangle_n, \quad \hat{z}_{m,n} = {z}_{m,n} - \langle z_m \rangle_n.$$

 figure: Fig. 2.

Fig. 2. Illustration of $x$-deviations $\Delta _{m,n}^{(x)}$ and $y$-deviations $\Delta _{m,n}^{(y)}$ of a vision ray (blue line) for a single pixel $m$ of the imaging system and three calibration measurements. Deviations between reference points $(x_{m,n},y_{m,n},z_{m,n})^T$ and the fitted vision ray $\vec {r}_m$ are measured along the corresponding axis of the coordinate system within $z$-planes given by $z_{m,n}$ reference points coordinate values. The three colors correspond to those of the points in Fig. 1.

Download Full Size | PDF

Here, $\langle \cdots \rangle _n$ indicates average values with respect phase measurements $n$

$$\langle\xi_m\rangle_n = \frac{1}{N} \sum_n \xi_{m,n},$$
where $N$ is the number of phase measurements contributing reference point information. With these, $x$- and $y$-deviations are given by [23]
$$\Delta_{m,n}^{(x)} = \hat{x}_{m,n} - u_m\cdot\hat{z}_{m,n}, \quad \Delta_{m,n}^{(y)} = \hat{y}_{m,n} - v_m\cdot\hat{z}_{m,n}.$$

Furthermore, we can compute the slope of individual vision rays using [24]

$$u_m = \frac{ \langle \hat{x}_m \hat{z}_m \rangle_n }{ \langle \hat{z}_m^2 \rangle_n }, \quad v_m = \frac{ \langle \hat{y}_m \hat{z}_m \rangle_n }{ \langle \hat{z}_m^2 \rangle_n }.$$

Using Eqs. (8) and (9) for the calculation of vision ray deviations does not require explicit line regression for determination of vision ray parameters. Error values $\Delta _m$ constitute a measure for the collinearity of reference points corresponding to individual sensor pixels $m$. The sum

$$\Delta = \sum_m \Delta_m$$
over all these error contributions provides a measure for the overall system accuracy. This system error value $\Delta$ depends on pose parameter values in $\vec {T}$ and surface coefficients in $\vec {C}$ due to the relation of local and global reference points given by Eqs. (1) and (2). Hence, Eq. (10) can be used as a cost function for the determination of unknown system parameters. Identifying the minimum of $\Delta$ with respect to system parameter values yields the best description of the calibration setup and provides the vision rays of the imaging system constituting the calibration result.

2.5 Determination of system parameters

We now define our system parameter vector $\vec {S}$ by writing all system parameters, except for those of one reference pose $\vec {T}_1$, into the vector

$$\vec{S} = (\vec{T}_2,\ldots, \vec{T}_n, \vec{C}).$$

We achieve calibration of our imaging system by determination of those system parameter values yielding the minimum of our error function given by Eq. (10) using numerical optimization. Parameter values of a reference pose $\vec {T}_1$ must be considered constant for the cost function to have a unique minimum. This is due to the cost function quantifying collinearity of reference points $(x_{m,n},y_{m,n},z_{m,n})^T$ for all pixel $m$ of the imaging system. Collinearity can be achieved anywhere in space, depending only on relative poses of the calibration target. As a consequence, vision rays resulting from the calibration are defined within an arbitrary coordinate system depending on the initial guess of this reference pose. However, by post-processing of vision ray data, it is easily possible to transform vision ray data into natural coordinate systems of individual imaging sensors [16,25]. It should be noted that the application of a VRC imaging systems in optical metrology typically only requires knowledge regarding relative orientations, rendering this step optional. Software development tools such as Matlab provide several functions for numerical optimization. We use a general-purpose trust-region algorithm as provided by the Matlab function fminunc [26]. This algorithm requires an analytical formula for the gradient of the cost function. Providing an analytical formula for the Hessian matrix of the cost function is not required since the algorithm can numerically determine its eigenvalues. However, calculating the elements of the Hessian matrix analytically significantly reduces computation time. Due to the continuous spatial coding of the calibration target in VRC, the amount of reference point data is very large compared to photogrammetric camera calibration with marker plates. Therefore, computational efficient numerical optimization is beneficial to make best use of the available data. The entries of gradient $G_k$ and of the Hessian matrix $H_{j,k}$ are determined by calculating first and second order derivatives of the cost function given by Eq. (10) with respect to all system parameters $S_k$, where $S_k$ is the $k$th entry of the system parameter vector $\vec {S}$:

$$G_k = \frac{ \partial\Delta }{ \partial S_k }, \quad H_{j,k} = \frac{ \partial^2\Delta }{ \partial S_j \partial S_k }.$$

2.6 Derivatives of reference points

The dependence of the cost function given by Eq. (10) originates from the transformation of local reference points $(x_{m,n}^l,y_{m,n}^l)^T$ into global reference points $(x_{m,n},y_{m,n},z_{m,n})^T$. Determination of gradient and Hessian matrix of our cost function therefore requires calculation of derivatives of global reference points with respect to the entries $S_k$ of the system parameter vector $\vec {S}$. These derivatives can be derived from Eqs. (1) and (2). Regarding derivatives with respect to translation parameters, one obtains

$$\frac{\partial}{\partial t_i^{(x)}} \begin{pmatrix} x_{m,n} \\ y_{m,n} \\ z_{m,n} \end{pmatrix} = \delta_{i,n} \begin{pmatrix} 1 \\ 0 \\ 0 \end{pmatrix}, \quad \frac{\partial}{\partial t_i^{(y)}} \begin{pmatrix} x_{m,n} \\ y_{m,n} \\ z_{m,n} \end{pmatrix} = \delta_{i,n} \begin{pmatrix} 0 \\ 1 \\ 0 \end{pmatrix}, \quad \frac{\partial}{\partial t_i^{(z)}} \begin{pmatrix} x_{m,n} \\ y_{m,n} \\ z_{m,n} \end{pmatrix} = \delta_{i,n} \begin{pmatrix} 0 \\ 0 \\ 1 \end{pmatrix}.$$

The Kronecker delta $\delta _{i,n}$ results from the fact that derivatives with respect to a pose parameter corresponding to the $i$th measurement only affect reference points of that measurement. The same applies to derivatives with respect to rotation angles:

$$\begin{aligned} \frac{\partial}{\partial \alpha_i} \begin{pmatrix} x_{m,n} \\ y_{m,n} \\ z_{m,n} \end{pmatrix} &= \delta_{i,n} R^{(z)}(\gamma_i) R^{(y)}(\beta_i) R^{'(x)}(\alpha_i) \begin{pmatrix} x_{m,n}^l \\ y_{m,n}^l \\ z_{m,n}^l \end{pmatrix}, \\ \frac{\partial}{\partial \beta_i} \begin{pmatrix} x_{m,n} \\ y_{m,n} \\ z_{m,n} \end{pmatrix} &= \delta_{i,n} R^{(z)}(\gamma_i) R^{'(y)}(\beta_i) R^{(x)}(\alpha_i) \begin{pmatrix} x_{m,n}^l \\ y_{m,n}^l \\ z_{m,n}^l \end{pmatrix},\\ \frac{\partial}{\partial \gamma_i} \begin{pmatrix} x_{m,n} \\ y_{m,n} \\ z_{m,n} \end{pmatrix} &= \delta_{i,n} R^{'(z)}(\gamma_i) R^{(y)}(\beta_i) R^{(x)}(\alpha_i) \begin{pmatrix} x_{m,n}^l \\ y_{m,n}^l \\ z_{m,n}^l \end{pmatrix}. \end{aligned}$$

Here, $R^{'(x)}$, $R^{'(y)}$ and $R^{'(z)}$ are the derivatives of the corresponding rotation matrices. Derivatives with respect to polynomial coefficients $c_{pq}$ of the surface parameterization of the calibration target only affect local $z$-coordinates $z_{m,n}^l$. Using Eq. (1) one obtains

$$\begin{aligned} \begin{pmatrix} x_{m,n} \\ y_{m,n} \\ z_{m,n} \end{pmatrix} &= \frac{\partial}{\partial c_{ij}} \left( R^{(z)}(\gamma_n) R^{(y)}(\beta_n) R^{(x)}(\alpha_n) \begin{pmatrix} x_{m,n}^l \\ y_{m,n}^l \\ \sum_{p,q} c_{pq} \cdot (x_{m,n}^l)^p (y_{m,n}^l)^q \end{pmatrix} \right) \\ &= (x_{m,n}^l)^i (y_{m,n}^l)^j R^{(z)}(\gamma_n) R^{(y)}(\beta_n) R^{(x)}(\alpha_n) \begin{pmatrix} 0 \\ 0 \\ 1 \end{pmatrix}. \end{aligned}$$

For conciseness, we will omit the second order derivatives required for the calculation of the Hessian matrix here. Since first order derivatives with respect to translations are constant, as shown in Eq. (13), second order derivatives that include one derivation with respect to translation will always be zero. Equation (15) shows that first order derivatives with respect to surface coefficients only depend on rotational pose parameters. Hence, second order derivatives including one derivation with respect to a surface coefficient are zero unless including the derivation with respect to a rotational pose parameter. In that case, the derivative can be determined from Eq. (14) using the scheme applied in Eq. (15). One can easily determine second order derivatives with respect to two rotational pose parameters from Eq. (14) using the same scheme as with the calculation of the first order derivatives.

2.7 Gradient of the cost function

From Eqs. (4), (10), and (12), we find that the entries $G_k$ of the gradient of our cost function can be separated into $x$- and $y$-contributions:

$$\begin{aligned} G_k &= \frac{\partial \Delta}{\partial S_k} = \sum_m \frac{\partial \Delta_m}{\partial S_k} = \sum_{m,n}\frac{\partial \Delta_{m,n}^{(x)^2}}{\partial S_k} + \sum_{m,n}\frac{\partial \Delta_{m,n}^{(y)^2}}{\partial S_k} \\ &= 2 \sum_{m,n} \Delta_{m,n}^{(x)} \frac{\partial \Delta_{m,n}^{(x)}}{\partial S_k} + 2 \sum_{m,n} \Delta_{m,n}^{(y)} \frac{\partial \Delta_{m,n}^{(y)}}{\partial S_k}\\ &= G_k^{(x)} + G_k^{(y)}. \end{aligned}$$

As shown in Eqs. (13), (14), and (15), we can calculate the derivatives of reference points with respect to system parameters. Hence, the task at hand is to find expressions for the terms Eq. (16) which only contain these derivatives.

Due to the definition of $\Delta _{m,n}^{(x)}$ and $\Delta _{m,n}^{(y)}$ in Eq. (8), the only difference between the expressions for $x$- and $y$-contributions of the gradient in Eq. (16) are the centered reference point coordinate components $\hat {x}_{m,n}$ and $\hat {y}_{m,n}$ and the vision ray slope components $u_m$ and $v_m$. By introducing $\xi$ as a placeholder for either $x$ or $y$ respectively, calculating

$$\frac{1}{2} G_k^{(\xi)} = \frac{1}{2} \sum_{m,n} \frac{\partial \Delta_{m,n}^{(\xi)^2}}{\partial S_k} = \sum_{m,n} \Delta_{m,n}^{(\xi)} \frac{\partial \Delta_{m,n}^{(\xi)}}{\partial S_k}$$
provides the solution for Eq. (16).

The calculation of error contributions $\Delta _{m,n}^{(\xi )}$ in Eq. (17) uses centered reference point coordinates. Derivatives of centered coordinates, required to solve Eq. (17), can be derived from Eq. (6), yielding

$$\frac{\partial \hat{\xi}_{m,n}}{\partial S_k} = \frac{\partial \xi_{m,n}}{\partial S_k} - \frac{\partial \langle \xi_m \rangle_n}{\partial S_k}, \quad \frac{\partial \hat{z}_{m,n}}{\partial S_k} = \frac{\partial z_{m,n}}{\partial S_k} - \frac{\partial \langle z_m \rangle_n}{\partial S_k}.$$

Solving Eq. (17), we will make use of the fact that the sum of $x$- and $y$-deviation of a single vision ray over all reference points

$$\sum_n \Delta_{m,n}^{(\xi)} = 0$$
is always zero, which can be derived from Eq. (8). Furthermore, the sum of products of centered $z$-coordinates with corresponding vision ray deviations over all reference points
$$\begin{aligned} \sum_n \hat{z}_{m,n} \Delta_{m,n}^{(\xi)} &= \sum_n \hat{z}_{m,n} (\hat{\xi}_{m,n} - \mu_m \hat{z}_{m,n}) \\ &= N \langle \hat{\xi}_m \hat{z}_m \rangle_n - N \mu_m \langle \hat{z}_m^2\rangle_n\\ &= N \langle \hat{\xi}_m \hat{z}_m \rangle_n - N \frac{ \langle \hat{\xi}_{m,n} \hat{z}_{m,n} \rangle_n } { \langle \hat{z}_m^2 \rangle_n } \langle \hat{z}_m^2 \rangle_n\\ &= 0 \end{aligned}$$
also always yields zero. We introduced $\mu _m$ in Eq. (20) as a placeholder for vision ray slopes $u_m$ and $v_m$ in $x$- and $y$-direction corresponding to our $x$-$y$-placeholder $\xi$. Utilizing these conditions, a solution for the gradient components in Eq. (17), containing only derivatives of reference point coordinates to be calculated using Eqs. (13), (14), and (15), can be found:
$$\begin{aligned} \frac{1}{2} G_k^{(\xi)} &= \sum_{m,n} \Delta_{m,n}^{(\xi)} \frac{\partial}{\partial S_k} (\hat{\xi}_{m,n} - \mu_m \hat{z}_{m,n}) \\ &= \sum_m \left( \sum_n \Delta_{m,n}^{(\xi)} \frac{\partial \hat{\xi}_{m,n}}{\partial S_k} - \mu_m \sum_n \Delta_{m,n}^{(\xi)} \frac{\partial \hat{z}_{m,n}}{\partial S_k} - \frac{\partial \mu_m}{\partial S_k} \sum_n \Delta_{m,n}^{(\xi)} \hat{z}_{m,n} \right)\\ &= \sum_m \Bigg( \sum_n \Delta_{m,n}^{(\xi)} \frac{\partial \xi_{m,n}}{\partial S_k} - \frac{\partial \langle \xi_m \rangle_n}{\partial S_k} \sum_n \Delta_{m,n}^{(\xi)}\\ &\qquad- \mu_m \sum_n \Delta_{m,n}^{(\xi)} \frac{\partial z_{m,n}}{\partial S_k} + \mu_m \frac{\partial \langle z_m \rangle_n}{\partial S_k} \sum_n \Delta_{m,n}^{(\xi)} \Bigg)\\ &= \sum_{m,n} \Delta_{m,n}^{(\xi)} \frac{\partial \xi_{m,n}}{\partial S_k} - \sum_{m,n} \mu_m \Delta_{m,n}^{(\xi)} \frac{\partial z_{m,n}}{\partial S_k}\\ &= \sum_{m,n} \Delta_{m,n}^{(\xi)} \left( \frac{\partial \xi_{m,n}}{\partial S_k} - \mu_m \frac{\partial z_{m,n}}{\partial S_k} \right). \end{aligned}$$

We obtain the expression in line one of Eq. (21) by inserting the definition of $\Delta _{m,n}^{(\xi )}$ from Eq. (8) into Eq. (17). Using the definition for derivatives of centered coordinates from Eq. (6) yields the term $\frac {\partial \mu _m}{\partial S_k} \sum _n \Delta _{m,n}^{(\xi )} \hat {z}_{m,n}$ in line two, which is zero as shown in Eq. (20). Next, we replace the remaining derivatives of centered coordinates with the definition in Eq. (6) and obtain two terms containing $\sum _n \Delta _{m,n}^{(\xi )}$ which is zero due to Eq. (19). Inserting Eq. (21) into Eq. (16) provides the formula for the gradient vector of our cost function

$$G_k = \frac{\partial \Delta}{\partial S_k} = 2 \sum_{m,n} \Delta_{m,n}^{(x)} \left( \frac{\partial x_{m,n}}{\partial S_k} -u_m \frac{\partial z_{m,n}}{\partial S_k} \right) +2 \sum_{m,n} \Delta_{m,n}^{(y)} \left( \frac{\partial y_{m,n}}{\partial S_k} -v_m \frac{\partial z_{m,n}}{\partial S_k} \right)$$
which contains only derivatives of reference point coordinates.

2.8 Hessian matrix of the cost function

Let $S_j$ and $S_k$ be the $j$th and the $k$th component of our system parameter vector $\vec {S}$. To calculate the entries $H_{j,k}$ of the Hessian matrix of our cost function, we have to determine the second order derivatives of our cost function with respect to all combinations of $j$ and $k$ of system parameter indices. As with the derivation of gradient entries, the goal is to find a formula only containing derivatives of reference point coordinates. Here, we can use the formula for the entries of the gradient vector $\vec {G}$ given by Eq. (22) to calculate the entries $H_{j,k}$ of the Hessian matrix as derivatives of the gradient:

$$H_{j,k} = \frac{\partial^2 \Delta}{\partial S_j \partial S_k} = \frac{\partial G_k}{\partial S_j} = \frac{\partial G_k^{(x)}}{\partial S_j} + \frac{\partial G_k^{(y)}}{\partial S_j} = H_{j,k}^{(x)} + H_{j,k}^{(y)}.$$

As introduced in the derivation of the gradient of the cost function, we will again substitute coordinate components $x_{m,n}$ and $y_{m,n}$ with $\xi _{m,n}$ and vision ray slope components $u_m$ and $v_m$ with $\mu _m$. Furthermore, we substitute the bracket terms in Eq. (22) with

$$A_{m,n}^{(\xi)}(S_k) = \left( \frac{\partial \xi_{m,n}}{\partial S_k} - \mu_m \frac{\partial z_{m,n}}{\partial S_k} \right).$$

Using Eqs. (22), (23), and (24), this leaves us with the task of solving

$$\frac{1}{2} H_{j,k}^{(\xi)} = \frac{1}{2} \frac{\partial G_k^{(\xi)}}{\partial S_j} = \frac{\partial}{\partial S_j} \sum_{m,n} \Delta_{m,n}^{(\xi)} A_{m,n}^{(\xi)}(S_k)$$
until only derivatives of reference point coordinates remain.

On behalf of conciseness, we will omit an explicit derivation for solving Eq. (25) as we did for the gradient of the cost function in Eq. (21). The procedure is substantial but straight forward. Since both, $\Delta _{m,n}^{(\xi )}$ and $A_{m,n}^{(\xi )}(S_k)$ depend on system parameters $\vec {S}$, the first step is to apply product rule. Afterwards, one has to substitute the identities of bracket terms $A_{m,n}^{(\xi )}(S_k)$, vision ray residuals $\Delta _{m,n}^{(\xi )}$, vision ray slopes $\mu _m$ and derivatives of centered coordinates $\frac {\partial \hat {\xi }_{m,n}}{\partial S_k}$ given by Eqs. (8), (9), (18), and (24) respectively. The derivation of vision ray slopes $\mu _m$ yields derivatives of average values of squared and mixed centered coordinates

$$\frac{\partial \langle \hat{z}^2 \rangle_n}{\partial S_j} = 2 \left\langle \hat{z} \frac{\partial z}{\partial S_j} \right\rangle_n, \quad \frac{\partial \langle \hat{\xi} \hat{z} \rangle_n}{\partial S_j} = \left\langle \frac{\partial \xi}{\partial S_j} \hat{z} + \frac{\partial z}{\partial S_j} \hat{x} \right\rangle_n.$$

Unfortunately, in contrast to the derivation of our gradient formula Eq. (21), we do not encounter terms yielding zero. As our final result for solving Eq. (25), we obtain

$$\begin{aligned} \frac{1}{2} H_{j,k}^{(\xi)} &= \sum_{m,n} \Delta_{m,n}^{(\xi)} \left( \frac{\partial^2 \xi_{m,n}}{\partial S_j \partial S_k} -\mu_m \frac{\partial^2 z_{m,n}}{\partial S_j \partial S_k} \right) \\ &- \sum_m \frac{1}{\sum_n \hat{z}_{m,n}^2} \left( \sum_n \left( A_{m,n}^{(\xi)}(S_j) \hat{z}_{m,n} + \Delta_{m,n}^{(\xi)} \frac{\partial z_{m,n}}{\partial S_j} \right) \right)\\ & \qquad \cdot \left( \sum_m \left( A_{m,n}^{(\xi)}(S_k) \hat{z}_{m,n} + \Delta_{m,n}^{(\xi)} \frac{\partial z_{m,n}}{\partial S_k} \right) \right) \\ &- \frac{1}{N} \sum_m \left( \sum_n A_{m,n}^{(\xi)}(S_j) \right) \cdot \left( \sum_n A_{m,n}^{(\xi)}(S_k) \right) + \sum_{m,n} A_{m,n}^{(\xi)}(S_j) \cdot A_{m,n}^{(\xi)}(S_k) \end{aligned}$$
which contains only derivatives of reference points. Hence, by inserting Eq. (27) into Eq. (23), we have found an analytical formula for the calculation of the Hessian matrix of our cost function defined in Eq. (10).

3. Experimental validation

3.1 Vision ray calibration

Having described the relation of measured local reference points to global reference points, we can determine vision rays of imaging systems based on calibration measurements and assumed system parameter values. Using our cost function, one can perform numerical optimization of assumed system parameters to obtain best calibration results. Gradient and Hessian matrix of our cost function can be explicitly calculated, enabling a computational efficient optimization procedure. The calibration can be applied to determine vision rays of a single camera as well as for imaging systems of multiple cameras since no distinctions between individual sensor pixels and corresponding optical elements are made. The latter case states a holistic calibration of all cameras of the imaging system. For demonstration and qualification of our approach, we perform a Vision Ray Calibration of an imaging system consisting of two cameras and a Liquid Crystal Display (LCD). The geometry of the setup and the technical details of the components are given in Fig. 3.

 figure: Fig. 3.

Fig. 3. Representation of the calibration setup as determined by vision ray calibration. 20 planes illustrate the relative pose of the LCD surface used as calibration target. Blue lines indicate the fields of view of both cameras. We use two AVT Manta G-223 B cameras with a resolution of 2048x1088 pixels using Ricoh FL-CC2514A-2M objectives with a focal length of 25 mm. Both cameras were mounted on a CFRP pipe. The calibration target is a common consumer Samsung LU28E590DS/EN Liquid-Crystal Display with a resolution of 3840x2160 pixels.

Download Full Size | PDF

We perform 20 calibration measurements using different relative poses between imaging system and calibration target. To acquire an initial guess for these poses we use Zhang’s calibration algorithm [4] on measured local reference point data from PST of one camera, assuming an ideal pinhole camera without distortions and a perfectly plane calibration target. We reliably obtain a good estimation with low computational effort by applying downsampling by a factor of about 10,000, using only every 100th pixel of each line and column of the 2048×1088 pixel matrix of the sensor. For the polynomial parameterization of the LCD surface in the main step of the Vision Ray Calibration, we obtain good results using polynomial terms up to the 5th order. Systematic investigations on how the choice of polynomial terms affects the calibration result are pending. In any case, one must omit polynomial terms corresponding to $c_{00}$, $c_{01}$ and $c_{10}$ [see Eq. (1)]. These represent offset and tilt of the calibration target, which are already described by the pose parameters in $\vec {T}_n$. Overall, we obtain 33 polynomial coefficients which we initialize with zeros. In total, 147 system parameters, 6 pose parameters for each measurement, except for one used as reference, plus the 33 polynomial coefficients, must be determined by numerical optimization. For this demonstration, we apply downsampling by a factor of about 400 to the reference point data of each camera, utilizing every 20th pixel of each line and column of the sensor matrix. This yields 226,600 reference points, 5665 for each camera and each measurement. The optimization takes 7 iterations and about 40 seconds. Using Matlabs fminunc-function [26], there are several exit criteria. The one that was met is the step tolerance criterium which was set to $10^{-11}$, stating an absolute limit due to the accuracy of numerical precision. Practically this means that further optimization was prevented due to numerical accuracy while evaluating the cost function. For the subsequent calculation of vision rays for each pixel of both cameras, utilizing the system parameters determined by the optimization, we use all measurement data without any downsampling or interpolation. This step takes another about 14 seconds and completes the calibration process. All calculations were performed with Matlab Version R2019b Update 7 on a desktop computer with Windows 10 Pro Version 1909, an Intel Core i5-7400 CPU and 16 GB of RAM.

The calibration result provides poses and shape of the calibration target and a vision ray for each sensor pixel in a common coordinate system. The vision ray data implicitly contains the information of relative positions and orientation of both cameras. This is sufficient for application in optical metrology but inconvenient for graphical representation of the calibration setup. Therefore, we determine natural coordinate systems for the vision rays of the individual cameras as described in [25] and corresponding transformations. Using this additional information, we can finally display the graphical representation of the calibration setup shown in Fig. 3. Here, 20 planes mark the poses of the calibration target in relation to the cameras. Blue lines indicate the fields of view of the cameras.

3.2 Fringe projection measurements

For qualification of the calibration result, we perform Fringe Projection measurements using the vision ray calibrated camera setup and a digital projector of type Optoma ML750ST for structured illumination. A spherical sample with curvature radius 37.5947 mm serves as reference surface. The sample was qualified by a Mitutoyo Crysta Apex C7106 coordinate measurement machine. Surface points measured with the coordinate measurement machine result in a standard deviation of 6.5 µm with respect to the sphere radius. From Fringe Projection measurements, absolute shape data is calculated based on calibration results using (a) the initial guess of the system parameters without any numerical optimization, (b) system parameters after only one iteration of numerical optimization, and (c) system parameters after the complete optimization procedure.

Figures 4 and 5 show our results by means of color-coded residuals in normal direction to the corresponding spherical fit. Visible high-frequency microstructures of the illustrated residuals are real surface features of the sample. Figure 4 shows the residuals with respect the best fit of a sphere with reference radius 37.5947 mm. Residuals in Fig. 4(a) set the scale of color-code in all subimages, visualizing the decline of low-frequency systematic measurement errors due to the numerical optimization. Standard deviation (STD) is reduced from 25.0 µm in (a) to 14.6 µm in (b) due to only one iteration of numerical optimization. Using the final calibration result achieved after seven iterations, the STD is reduced to 8.8 µm.

 figure: Fig. 4.

Fig. 4. Measured shape deviations of a spherical sample with radius 37.5947 mm using FPT with different system calibrations. Deviations in normal direction with respect to the fit sphere surface. The dark dot in the center of the sphere is not an artefact, but a surface flaw of sample. The accuracy of the measured shape increases due to the optimization of system parameters during the holistic camera calibration. (a) Deviations using initial guess of system parameters. STD: 25.0 µm. (b) Deviations using system parameters received after first iteration of VRC. STD: 14.6 µm. (c) Deviations using system parameters from full VRC-process (7 iterations). STD: 8.8 µm.

Download Full Size | PDF

 figure: Fig. 5.

Fig. 5. Measured shape deviations of the spherical sample with respect to a best-fit sphere using FPT with different system calibrations. Measurement data, used in Fig. 4 before, yields different residuals since no constraints to the radius of the best-fit sphere were made. Deviations in normal direction with respect to the fit sphere surface. Deviations decrease due to the optimization of system parameters during the holistic camera calibration. (a) Deviations using initial guess of system parameters. STD: 12.9 µm, radius: 37.7086 mm. (b) Deviations using system parameters received after first iteration of VRC. STD: 10.4 µm, radius: 37.6481 mm. (c) Deviations using system parameters from full VRC-process (7 iterations). STD: 8.5 µm, radius: 37.5846 mm.

Download Full Size | PDF

Residuals displayed in Fig. 5 are derived from the same measurement data, using best-fit spheres without constraints on the sphere radius. Colormaps are individually adjusted based on minimum and maximum residual values. Shape data in (a), obtained by using the initial guess of system parameters to determine vision rays, yields a radius of 37.7086 mm with STD of 12.9 µm of residuals. Performing one iteration of numerical optimization reduced the best-fit sphere radius to 37.6481 mm and STD of residuals to 10.4 µm as shown in (b). Completing numerical optimization in (c) yields a best-fit sphere radius of 37.5846 mm and an STD of 8.5 µm. Comparison of visible low-frequency structures in the residual images and of STD-values shown in Fig. 5 demonstrates that the numerical optimization of system parameters in Vision Ray Calibration does not only increase accuracy of measured radii but also reduces systematical deviations with respect to the spherical shape for our sample.

Our result for the radius of the spherical sample of 37.5846 mm does not quite match reference radius 37.5947 mm. However, we find that the sum of standard deviations of residuals from measurements with fringe projection and coordinate measurement machine of 15.0 µm is smaller than the difference in obtained sphere radii of 10.1 µm. Furthermore, comparison of residuals displayed in Figs. 4(c) and 5(c) shows that the difference of 10.1 µm of fitted sphere radius affects the STD of residuals by only 0.3 µm.

Examining the effect of the number of iterations during optimization on the calibration result shows that we have chosen a far too severe exit criterium. This is demonstrated in Fig. 6 using the radius of the best-fit sphere and the standard deviation of residuals received from spherical fits using the reference radius. After only two iterations of the numerical optimization, the accuracy of the measurements of our sample has converged close to that of the coordinate measurement machine. Based on our calibration measurements, performing more than three iterations yields no significant improvements.

 figure: Fig. 6.

Fig. 6. Convergence of (a) best-fit sphere radius and (b) standard deviation (STD) of residuals with respect to the reference radius against the number iterations for the determination of system parameters during calibration of the Fringe Projection measurement setup. Corresponding reference values (blue lines) were obtained using a coordinate measurement machine. The results show that we achieve measurement accuracy close to that of the coordinate measurement machine after only 2 iterations of optimization during system calibration.

Download Full Size | PDF

4. Conclusion

Vision ray calibration utilizes a generic model for the description of imaging devices. An independent vision ray is determined for each sensor pixel without any model assumptions regarding the imaging system. Hence, a setup of multiple cameras can be considered as one sensor unit. This comes with the advantage that an algorithm for vision ray calibration of a single camera will provide full holistic calibration when applied to such a setup. Subsequent determination of relative positions and orientations of cameras can be completely omitted since this information is implicitly contained within the obtained vision ray data. This renders vision ray calibration exceptionally suited for application in optical metrology techniques based on triangulation such as Fringe Projection Technique.

Compared to model-based photogrammetric camera calibration using marker plates, VRC requires processing of large amounts of data due to the continuous spatial coding of the calibration target, yielding a reference point for each sensor pixel in each measurement. It is therefore desirable to find numerically efficient calibration algorithms, minimizing computational costs while making best use of available data. As demonstrated in this work, a cost functions for the necessary numerical optimization does not require the numerically expensive explicit calculation of vision rays. Additionally, avoiding explicit calculation of vision rays furthermore allows analytical derivation of gradient and Hessian matrix of the cost function for efficient optimization of system parameters.

By performing holistic VRC to a stereo fringe projection metrology system, we demonstrated that our algorithm handles large amounts of data in good time. In our experiments, we achieved full calibration in less than one minute computation time while utilizing 20 measurements and every 400th pixel of each 2-megapixel camera in the optimization procedure. Using a common desktop computer, calculation of system parameters by means of numerical optimization required a computation time of about 40 s. Additional 14 s were required for the subsequent calculation of vision ray parameters for all camera pixels. Evaluation of the measurement of a spherical reference surface shows that the numerical optimization during calibration in fact improves accuracy of the metrology systems. Performing full calibration, we obtained an accuracy of the shape measurement close to that of a coordinate measurement machine with our fringe projection setup.

We would like to point out that this work focusses on the demonstration of vision ray calibration and the algorithms involved, providing the means for efficient implementation. The formulae for gradient and Hessian matrix we derived do not contain any constraints regarding system parameters. Hence, further system parameters can be introduced in order to improve or adapt the calibration procedure.

Currently, we cannot provide quantitative comparison of accuracy of vision ray calibration with photogrammetric calibration techniques. Such an investigation requires the development of a holistic photogrammetric calibration algorithm. Also, the influence of the calibration target must be considered since vision ray calibration relies on a display for this purpose. These investigations will be addressed in future work.

Funding

Deutsche Forschungsgemeinschaft (289307220, 411170139).

Acknowledgments

The authors gratefully acknowledge funding from the Deutsche Forschungsgemeinschaft (DFG) for the projects "Sichtstrahl Plus" (Vision ray plus) and "MultiDeflect" (Holistic multi camera deflectometry). J. Bartsch acknowledges valuable discussions with Mr. M. Kalms and technical assistance provided by Mr. R. Klattenhoff and Mr. C. Kapitza, enabling our investigations.

Disclosures

The authors declare no conflict 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. J. Burke, T. Bothe, W. Osten, and C. F. Hess, “Reverse engineering by fringe projection,” Proc. SPIE 4778, 312–324 (2014). [CrossRef]  

2. M. C. Knauer, J. Kaminski, and G. Häusler, “Phase Measuring Deflectometry: a new approach to measure specular free-form surface,” Proc. SPIE 5457, 366–377 (2014). [CrossRef]  

3. T. Bothe, W. Li, C. Kopylow, and W. P. Jüptner, “High-resolution 3D shape measurement on specular surfaces by fringe reflection,” Proc. SPIE 5457, 411–422 (2014). [CrossRef]  

4. Z. Zhang, “A flexible new technique for camera calibration,” in Proceedings of IEEE Transactions on Pattern Analysis and Machine Intelligence22 (IEEE, 2000), pp. 1330–1334.

5. T. Luhmann, “Erweiterte Verfahren zur geometrischen Kamerakalibrierung in der Nahbereichsphotogrammetrie,” in: Deutsche Geodätische Kommission, Reihe C, Nr. 645 (Verlag der Bayerischen Akademie der Wissenschaften in Kommission beim Verlag CH Beck, 2010).

6. A. S. Poulin-Girard, S. Thibault, and D. Laurendeau, “Influence of camera calibration conditions on the accuracy of 3D reconstruction,” Opt. Express 24(3), 2678–2686 (2016). [CrossRef]  

7. W. Li, S. Shan, and H. Liu, “High-precision method of binocular camera calibration with a distortion model,” Appl. Opt. 56(8), 2368–2377 (2017). [CrossRef]  

8. Z. Liu, Q. Wu, S. Wu, and X. Pan, “Flexible and accurate camera calibration using grid spherical images,” Opt. Express 25(13), 15269–15285 (2017). [CrossRef]  

9. L. R. Ramírez-Hernández, J. C. Rodríguez-Quiñonez, M. J. Castro-Toscano, D. Hernández-Balbuena, W. Flores-Fuentes, R. Rascón-Carmona, L. Lindner, and O. Sergiyenko, “Improve three-dimensional point localization accuracy in stereo vision systems using a novel camera calibration method,” International Journal of Advanced Robotic Systems 17(1), 172988141989671 (2020). [CrossRef]  

10. D. Li, “Fast Phase-based Stereo Matching Method for 3D Shape Measurement,” in Proceedings of 2010 International Symposium on Optomechatronic Technologies (IEEE, 2010), pp. 1–5.

11. S. B. Werling, “Deflektometrie zur automatischen Sichtprüfung und Rekonstruktion spiegelnder Oberflächen,” in: Schriftenreihe Automatische Sichtprüfung und Bildverarbeitung Band 3, edited by J. Beyerer, (KIT Scientific Publishing, 2011).

12. S. Ramalingam and P. Sturm, “A Unifying Model for Camera Calibration,” in Proceedings of IEEE Transactions on Pattern Analysis and Machine Intelligence39 (IEEE, 2017), pp. 1309–1319.

13. M. D. Grossberg and S. K. Nayar, “A General Imaging Model and a Method for Finding its Parameters,” in Proceedings of 8th IEEE International Conference on Computer Vision (IEEE, 2001), pp. 108–115.

14. P. Sturm and S. Ramalingam, “A Generic Calibration Concept – Theory and Algorithms,” INRIA - Institut national de recherche en informatique et en automatique, RR-5058 (2003).

15. S. Ramalingam, P. Sturm, and S. K. Lodha, “Theory and Experiments towards Complete Generic Calibration,” INRIA - Institut national de recherche en informatique et en automatique, RR-5562 (2005).

16. T. Bothe, W. Li, M. Schulte, C. Kopylow, R. B. Bergmann, and W. P. Jüptner, “Vision ray calibration for the quantitative geometric description of general imaging and projection optics in metrology,” Appl. Opt. 49(30), 5851 (2010). [CrossRef]  

17. W. Li, M. Schulte, T. Bothe, C. Kopylow, N. Köpp, and W. Jüptner, “Beam based calibration for optical imaging device,” in Proceedings of 3DTV Conference (IEEE, 2007), pp. 1–4.

18. C. Zuo, S. Feng, L. Huang, T. Tao, W. Yin, and Q. Chen, “Phase shifting algorithms for fringe projection profilometry: A review,” Optics and Lasers in Engineering 109, 23–59 (2018). [CrossRef]  

19. W. Li, T. Bothe, M. Schulte, C. Kopylow, and W. Jüptner, “Sichtrahlkalibrierung für optisch abbildende Systeme,” in: DGaO-Proceedings 2008 – ISSN: 1614-8436.

20. T. Reh, W. Li, J. Burke, and R. B. Bergmann, “Improving the Generic Camera Calibration technique by an extended model of calibration display,” J. Europ. Opt. Soc. Rap. Public. 9, 14044 (2014). [CrossRef]  

21. M. Fischer, M. Petz, and R. Tutsch, “Evaluation of LCD monitors for deflectometric measurement systems,” Proc. SPIE 7726, 77260V (2010). [CrossRef]  

22. J. Bartsch, M. Kalms, and R. B. Bergmann, “Improving the Calibration of phase measuring deflectometry by a polynomial representation of the display shape,” J. Europ. Opt. Soc. Rap. Public. 15(1), 20 (2019). [CrossRef]  

23. N. R. Draper and H. Smith, Applied Regression Analysis, Vol. 326 (John Wiley & Sons, 1998), pp. 15–45.

24. N. R. Draper and H. Smith, Applied Regression Analysis, Vol. 326 (John Wiley & Sons, 1998), p. 24.

25. M. H. U. Prinzler, J. Bartsch, M. Kalms, and R. B. Bergmann, “Metric for comparison of generic camera calibration,” Proc. SPIE 10750, 107500I (2018). [CrossRef]  

26. MathWorks, “fminunc,” from the web: https://www.mathworks.com/help/optim/ug/fminunc.html, last visited on February 05, 2021.

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.

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.


Figures (6)

Fig. 1.
Fig. 1. Schematic representation of VRC quantities using three calibration measurements. An LCD acts as calibration target and is observed from varying perspectives. (a) Three reference points (denoted green, red and blue) corresponding to an arbitrary pixel $m$ of the imaging system are shown. These reference points are the positions on the LCD surface that were observed by pixel $m$ in these measurements. Using surface parameters $\vec {C}$ and pose parameters $\vec {T}_n$ describing the shape of the LCD surface and its relative pose towards the imaging system for every measurement enable transformation into a global coordinate system. (b) A line fit through the global reference points yields the vision ray $\vec {r}_m$ corresponding to the pixel $m$ .
Fig. 2.
Fig. 2. Illustration of $x$ -deviations $\Delta _{m,n}^{(x)}$ and $y$ -deviations $\Delta _{m,n}^{(y)}$ of a vision ray (blue line) for a single pixel $m$ of the imaging system and three calibration measurements. Deviations between reference points $(x_{m,n},y_{m,n},z_{m,n})^T$ and the fitted vision ray $\vec {r}_m$ are measured along the corresponding axis of the coordinate system within $z$ -planes given by $z_{m,n}$ reference points coordinate values. The three colors correspond to those of the points in Fig. 1.
Fig. 3.
Fig. 3. Representation of the calibration setup as determined by vision ray calibration. 20 planes illustrate the relative pose of the LCD surface used as calibration target. Blue lines indicate the fields of view of both cameras. We use two AVT Manta G-223 B cameras with a resolution of 2048x1088 pixels using Ricoh FL-CC2514A-2M objectives with a focal length of 25 mm. Both cameras were mounted on a CFRP pipe. The calibration target is a common consumer Samsung LU28E590DS/EN Liquid-Crystal Display with a resolution of 3840x2160 pixels.
Fig. 4.
Fig. 4. Measured shape deviations of a spherical sample with radius 37.5947 mm using FPT with different system calibrations. Deviations in normal direction with respect to the fit sphere surface. The dark dot in the center of the sphere is not an artefact, but a surface flaw of sample. The accuracy of the measured shape increases due to the optimization of system parameters during the holistic camera calibration. (a) Deviations using initial guess of system parameters. STD: 25.0 µm. (b) Deviations using system parameters received after first iteration of VRC. STD: 14.6 µm. (c) Deviations using system parameters from full VRC-process (7 iterations). STD: 8.8 µm.
Fig. 5.
Fig. 5. Measured shape deviations of the spherical sample with respect to a best-fit sphere using FPT with different system calibrations. Measurement data, used in Fig. 4 before, yields different residuals since no constraints to the radius of the best-fit sphere were made. Deviations in normal direction with respect to the fit sphere surface. Deviations decrease due to the optimization of system parameters during the holistic camera calibration. (a) Deviations using initial guess of system parameters. STD: 12.9 µm, radius: 37.7086 mm. (b) Deviations using system parameters received after first iteration of VRC. STD: 10.4 µm, radius: 37.6481 mm. (c) Deviations using system parameters from full VRC-process (7 iterations). STD: 8.5 µm, radius: 37.5846 mm.
Fig. 6.
Fig. 6. Convergence of (a) best-fit sphere radius and (b) standard deviation (STD) of residuals with respect to the reference radius against the number iterations for the determination of system parameters during calibration of the Fringe Projection measurement setup. Corresponding reference values (blue lines) were obtained using a coordinate measurement machine. The results show that we achieve measurement accuracy close to that of the coordinate measurement machine after only 2 iterations of optimization during system calibration.

Equations (27)

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

z m , n l = p , q c p q ( x m , n l ) p ( y m , n l ) q .
( x m , n y m , n z m , n ) = R ( z ) ( γ n ) R ( y ) ( β n ) R ( x ) ( α n ) ( x m , n l y m , n l z m , n l ) + ( t n ( x ) t n ( y ) t n ( z ) ) .
r m ( z ) = ( x m , 0 y m , 0 0 ) + ( u m v m 1 ) z ; z 0
Δ m = n | | ( x m , n y m , n z m , n ) r m ( z m , n ) | | 2 = n Δ m , n ( x ) 2 + n Δ m , n ( y ) 2
Δ m , n ( x ) = x m , n u m z m , n x m , 0 , Δ m , n ( y ) = y m , n v m z m , n y m , 0 .
x ^ m , n = x m , n x m n , y ^ m , n = y m , n y m n , z ^ m , n = z m , n z m n .
ξ m n = 1 N n ξ m , n ,
Δ m , n ( x ) = x ^ m , n u m z ^ m , n , Δ m , n ( y ) = y ^ m , n v m z ^ m , n .
u m = x ^ m z ^ m n z ^ m 2 n , v m = y ^ m z ^ m n z ^ m 2 n .
Δ = m Δ m
S = ( T 2 , , T n , C ) .
G k = Δ S k , H j , k = 2 Δ S j S k .
t i ( x ) ( x m , n y m , n z m , n ) = δ i , n ( 1 0 0 ) , t i ( y ) ( x m , n y m , n z m , n ) = δ i , n ( 0 1 0 ) , t i ( z ) ( x m , n y m , n z m , n ) = δ i , n ( 0 0 1 ) .
α i ( x m , n y m , n z m , n ) = δ i , n R ( z ) ( γ i ) R ( y ) ( β i ) R ( x ) ( α i ) ( x m , n l y m , n l z m , n l ) , β i ( x m , n y m , n z m , n ) = δ i , n R ( z ) ( γ i ) R ( y ) ( β i ) R ( x ) ( α i ) ( x m , n l y m , n l z m , n l ) , γ i ( x m , n y m , n z m , n ) = δ i , n R ( z ) ( γ i ) R ( y ) ( β i ) R ( x ) ( α i ) ( x m , n l y m , n l z m , n l ) .
( x m , n y m , n z m , n ) = c i j ( R ( z ) ( γ n ) R ( y ) ( β n ) R ( x ) ( α n ) ( x m , n l y m , n l p , q c p q ( x m , n l ) p ( y m , n l ) q ) ) = ( x m , n l ) i ( y m , n l ) j R ( z ) ( γ n ) R ( y ) ( β n ) R ( x ) ( α n ) ( 0 0 1 ) .
G k = Δ S k = m Δ m S k = m , n Δ m , n ( x ) 2 S k + m , n Δ m , n ( y ) 2 S k = 2 m , n Δ m , n ( x ) Δ m , n ( x ) S k + 2 m , n Δ m , n ( y ) Δ m , n ( y ) S k = G k ( x ) + G k ( y ) .
1 2 G k ( ξ ) = 1 2 m , n Δ m , n ( ξ ) 2 S k = m , n Δ m , n ( ξ ) Δ m , n ( ξ ) S k
ξ ^ m , n S k = ξ m , n S k ξ m n S k , z ^ m , n S k = z m , n S k z m n S k .
n Δ m , n ( ξ ) = 0
n z ^ m , n Δ m , n ( ξ ) = n z ^ m , n ( ξ ^ m , n μ m z ^ m , n ) = N ξ ^ m z ^ m n N μ m z ^ m 2 n = N ξ ^ m z ^ m n N ξ ^ m , n z ^ m , n n z ^ m 2 n z ^ m 2 n = 0
1 2 G k ( ξ ) = m , n Δ m , n ( ξ ) S k ( ξ ^ m , n μ m z ^ m , n ) = m ( n Δ m , n ( ξ ) ξ ^ m , n S k μ m n Δ m , n ( ξ ) z ^ m , n S k μ m S k n Δ m , n ( ξ ) z ^ m , n ) = m ( n Δ m , n ( ξ ) ξ m , n S k ξ m n S k n Δ m , n ( ξ ) μ m n Δ m , n ( ξ ) z m , n S k + μ m z m n S k n Δ m , n ( ξ ) ) = m , n Δ m , n ( ξ ) ξ m , n S k m , n μ m Δ m , n ( ξ ) z m , n S k = m , n Δ m , n ( ξ ) ( ξ m , n S k μ m z m , n S k ) .
G k = Δ S k = 2 m , n Δ m , n ( x ) ( x m , n S k u m z m , n S k ) + 2 m , n Δ m , n ( y ) ( y m , n S k v m z m , n S k )
H j , k = 2 Δ S j S k = G k S j = G k ( x ) S j + G k ( y ) S j = H j , k ( x ) + H j , k ( y ) .
A m , n ( ξ ) ( S k ) = ( ξ m , n S k μ m z m , n S k ) .
1 2 H j , k ( ξ ) = 1 2 G k ( ξ ) S j = S j m , n Δ m , n ( ξ ) A m , n ( ξ ) ( S k )
z ^ 2 n S j = 2 z ^ z S j n , ξ ^ z ^ n S j = ξ S j z ^ + z S j x ^ n .
1 2 H j , k ( ξ ) = m , n Δ m , n ( ξ ) ( 2 ξ m , n S j S k μ m 2 z m , n S j S k ) m 1 n z ^ m , n 2 ( n ( A m , n ( ξ ) ( S j ) z ^ m , n + Δ m , n ( ξ ) z m , n S j ) ) ( m ( A m , n ( ξ ) ( S k ) z ^ m , n + Δ m , n ( ξ ) z m , n S k ) ) 1 N m ( n A m , n ( ξ ) ( S j ) ) ( n A m , n ( ξ ) ( S k ) ) + m , n A m , n ( ξ ) ( S j ) A m , n ( ξ ) ( S k )
Select as filters


Select Topics Cancel
© Copyright 2024 | Optica Publishing Group. All rights reserved, including rights for text and data mining and training of artificial technologies or similar technologies.