A Novel Error Criterion of Fundamental Matrix Based on Principal Component Analysis

Abstract

Estimating the fundamental matrix (FM) using the known corresponding points is a key step for three-dimensional (3D) scene reconstruction, and its uncertainty directly affects camera calibration and point-cloud calculation. The symmetric epipolar distance is the most popular error criterion for estimating FM error, but it depends on the accuracy, number, and distribution of known corresponding points and is biased. This study mainly focuses on the error quantitative criterion of FM itself. First, the calculated FM process is reviewed with the known corresponding points. Matrix differential theory is then used to derive the covariance equation of FMs in detail. Subsequently, the principal component analysis method is followed to construct the scalar function as a novel error criterion to measure FM error. Finally, three experiments with different types of stereo images are performed to verify the rationality of the proposed method. Experiments found that the scalar function had approximately 90% correlation degree with the Manhattan norm, and greater than 80% with the epipolar geometric distance. Consequently, the proposed method is also appropriate for estimating FM error, in which the error ellipse or normal distribution curve is the reasonable error boundary of FM. When the error criterion value of this method falls into a normal distribution curve or an error ellipse, its corresponding FM is considered to have less error and be credible. Otherwise, it may be necessary to recalculate an FM to reconstruct 3D models.

1. Introduction

With the development of stereo observation technology, the fundamental matrix (FM), which describes the geometric relationship between stereo images of the same scene, plays an increasingly important role in the field of three-dimensional (3D) remote sensing. Using FM can eliminate false corresponding points [1,2] and compute the internal and external parameters of uncalibrated images in a camera [3,4,5] and make adjacent image connections in object tracking and 3D scene reconstruction [6,7,8]. Therefore, computing FM using uncalibrated stereo images is a key step in computer vision and 3D image processing.

FM was first introduced by Longuet–Higgins as a generalization of the essential matrix described in uncalibrated images [9]. A large number of FM-estimation algorithms then emerged, which are roughly divided into linear, iterative, and robust methods [10,11]. Among them, the robust method obtains higher estimation accuracy even when corresponding points are stored in false corresponding points (outliers) [12] where the random sample consensus (RANSAC) method for eliminating outliers is widely used for FM estimation [8,13,14]. Subsequently, many scholars have proposed various improved algorithms, e.g., the LMeds method [12], L-1 method [15], MLESAC method [16], MAPSAC method [17], FSASAC method [18], PROSAC and Promeds method [2,19], DKF-RANSAC method [20], R-CNN method [21], and so on.

The effectiveness of improved algorithms depends on FM error criteria. The main error criteria [8], e.g., algebraic distance, symmetric epipolar distance, Sampson distance, and so on, can be used to estimate the FM error. Fathy [22] studied the different FM error criteria, proved that the symmetric epipolar distance is biased, and proposed the Kanatani distance as a new error criterion to remove the outliers. Li [23] proposed the 3D metric distance criterion to estimate the FM quality. However, many recent scholars have retained the epipolar geometric distance as sufficiently accurate to estimate the quality of FM and its improved algorithms [24,25,26,27]. Although the epipolar geometric distance depends on the accuracy, number, and distribution of the known corresponding points on the two uncalibrated images and is biased, it is still the most popular criterion for evaluating improved FM algorithms.

In addition, the covariance matrix can be calculated and derived by analytical or statistical methods to measure FM error [8,28,29]. Hartley pointed out that the Monte Carlo statistical method can be used to calculate the FM covariance matrix [8]. Some scholars have concretely analyzed the covariance of FMs using the eight-point algorithm [30] and the random sample consensus algorithm [31]. Combined with the covariance matrix, the Manhattan distance has been used to measure FM error [32]. There are earlier studies measuring FM error using the covariance matrix.

It is necessary to develop a novel method to estimate FM error. Principal component analysis (PCA) is the most widely used method to reduce the dimensionality of such datasets, increasing interpretability but simultaneously minimizing information loss [33,34]. This paper intends to use the PCA method to reduce the FM dimensionality and construct a scalar function to measure FM error.

Compared with the epipolar geometric distance, the proposed method may be optimized and reduced the influence the accuracy, number, and distribution of the known corresponding points on FM error. Compared with the covariance matrix, the proposed method could use a numerical value to directly measure FM error. This study aims to (1) use FM itself to construct a more direct and less biased error criterion, and (2) determine whether a certain FM calculated with the known corresponding points on stereo images is reasonable and credible.

The remainder of this paper is organized as follows. Section 2 mainly introduces the process of computing FM using two uncalibrated images. Section 3 uses matrix differential theory to derive the covariance matrix of the PCA FM method to construct a scalar function to intuitively estimate FM error. Section 4 uses experiments to verify the rationality of the proposed method. Section 5 and Section 6 are the discussion and the conclusions of this work, respectively.

2. Computing FM

In Figure 1, the camera centers C and C’, a three-space point W, and image points x and x’ lie in a common plane π. A point, x’, in the first image is transferred via the common plane π to a corresponding point x in the second image. In symbols one may write x = Fx′ where F is the fundamental matrix (FM), and the other is x′ = F⁻¹ x.

FM represents the intrinsic projective geometry between two image planes. It is independent of scene structure and only depends on the cameras’ internal and external parameters of calibrated images. However, the FM between two uncalibrated images can usually be computed from image-point correspondences alone. Specifically, FM is defined by the equation

$${x^{\prime }}^{\mathrm{T}}Fx=0$$

(1)

for any pair of corresponding points $x \leftrightarrow x’$ on two uncalibrated images. Given sufficiently many corresponding points $x_i\leftrightarrow x_i^{\prime }$, Equation (1) can be used to compute the unknown matrix F, supposing that the unknown fundamental matrix F in Equation (1), which is a 3 × 3 matrix, is

$$F=\left[\begin{array}{ccc}{F}_{11}& F_{12}& F_{13}\\ F_{21}& F_{22}& F_{23}\\ F_{31}& F_{32}& F_{33}\end{array}\right]$$

(2)

The coefficients of F in Equation (1) are easily written in terms of the known coordinates x and x′. In particular, writing $x={\left(x,y,1\right)}^{\mathrm{T}}$ and $x\prime ={\left(x\prime ,y\prime ,1\right)}^{\mathrm{T}}$, each corresponding point gives rise to one linear Equation (3) in the unknown F entries.

$$x^{\prime }x{F}_{11}+x^{\prime }y{F}_{12}+x^{\prime }{F}_{13}+y^{\prime }x{F}_{21}+y^{\prime }y{F}_{22}+y^{\prime }{F}_{23}+x{F}_{31}+y{F}_{32}+F_{33}=0$$

(3)

Supposing that the nine-vector f is made up of the entries of F in row-major order.

$$f={\left[F_{11},F_{12},F_{13},F_{21},F_{22},F_{23},F_{31},F_{32},F_{33}\right]}^{\mathrm{T}}$$

(4)

Equation (3) can be expressed as a vector inner product

$$\left[x^{\prime }x,x^{\prime }y,x^{\prime },y^{\prime }x,y^{\prime }y,y^{\prime },x,y,1\right]f=0$$

(5)

Given a set of n corresponding points, a set of linear equations can be obtained.

$$M_{n\times 9}{f}_{9\times 1}=\left[\begin{array}{ccccccccc}{x}_1{}^{\prime }{x}_1& x_1{}^{\prime }{y}_1& x_1{}^{\prime }& y_1{}^{\prime }{x}_1& y_1{}^{\prime }{y}_1& y_1{}^{\prime }& x_1& y_1& 1\\ ⋮& ⋮& ⋮& ⋮& ⋮& ⋮& ⋮& ⋮& ⋮\\ x_n{}^{\prime }{x}_n& x_n{}^{\prime }{y}_n& x_n{}^{\prime }& y_n{}^{\prime }{x}_n& y_n{}^{\prime }{y}_n& y_n{}^{\prime }& x_n& y_n& 1\end{array}\right]f=0$$

(6)

Equation (6) is a homogeneous set of equations, and matrix M must have ranked at most 8 for a non-zero f to exist. Here, the unknown vector f can only be determined up to a scale. If the data are not exact because of noise in the point coordinates, the rank of M is equal to 9, and the only solution is f = 0. To solve this question, a solution is to turn the set in Equation (6) into an inhomogeneous set of linear equations by imposing the $F_{33}=1$ condition for some entry of the vector f. In this case, Equation (6) can be converted to

$$\left\{\begin{array}{c}A\tilde{f}+L=0\\ L={\left[1, 1,\cdots ,1\right]}_{1\times n}^{\mathrm{T}}\end{array}\right.$$

(7)

where A has n rows and eight columns, and $\tilde{f}$ is an eight-element vector.

Such an equation may be solved for $\tilde{f}$ using the least-squares techniques. The specific calculation is

$$\tilde{f}=-{\left(A^{\mathrm{T}}A\right)}^{-1}{A}^{\mathrm{T}}L.$$

(8)

$\tilde{f}$ and $F_{33}$ can then be combined to form f in Equation (9). Here, f is the nine-element vector of FM concretely computed between two uncalibrated images.

$$f={\left[\begin{array}{cc}{\tilde{f}}^{\mathrm{T}}& F_{33}\end{array}\right]}^{\mathrm{T}}$$

(9)

3. Quantization of FM Error

3.1. Covariance Matrix of FM

The FMs calculated in Equations (8) and (9) involve numerical calculations and matrix-theory operations. Therefore, the following section intends to use matrix differential theory to construct the covariance of FM.

3.1.1. Matrix Differential Theory

Let $A\in F^{m\times n}$ be a matrix variable, and $A={\left(a_{ij}\right)}_{m\times n}$, $G\left(A\right)={\left(g_{kl}\right)}_{p\times q}$ be the matrix function of A. Here, g_ij (A) is a numerical function of A, and all of them are derivable at A $\left(i=1, 2,\cdots ,m; j=1, 2,\cdots , n; k=1, 2,\cdots , p; l=1, 2,\cdots , q\right)$. Thus,

$${\left[\frac{\partial G}{\partial a_{ij}}\right]}_{pm\times qn}=\nabla \otimes \Delta$$

(10)

This is the derivative of G with respect to A, denoted as $\frac{dG\left(A\right)}{dA}$. Among them, $\nabla =\left[\begin{array}{cccc}\frac{\partial G}{\partial a_{11}}& \frac{\partial G}{\partial a_{12}}& \cdots & \frac{\partial G}{\partial a_{1n}}\\ \frac{\partial G}{\partial a_{21}}& \frac{\partial G}{\partial a_{22}}& \cdots & \frac{\partial G}{\partial a_{2n}}\\ ⋮& ⋮& \ddots & ⋮\\ \frac{\partial G}{\partial a_{m1}}& \frac{\partial G}{\partial a_{m2}}& \cdots & \frac{\partial G}{\partial a_{mn}}\end{array}\right]$, $\Delta =\left[\begin{array}{cccc}\frac{\partial g_{11}}{\partial a_{ij}}& \frac{\partial g_{12}}{\partial a_{ij}}& \cdots & \frac{\partial g_{1q}}{\partial a_{ij}}\\ \frac{\partial g_{21}}{\partial a_{ij}}& \frac{\partial g_{22}}{\partial a_{ij}}& \cdots & \frac{\partial g_{2q}}{\partial a_{ij}}\\ ⋮& ⋮& \ddots & ⋮\\ \frac{\partial g_{p1}}{\partial a_{ij}}& \frac{\partial g_{p2}}{\partial a_{ij}}& \cdots & \frac{\partial g_{pq}}{\partial a_{ij}}\end{array}\right]$, and the symbol ⊗ is the Kronecker product of the matrices A and B of m × n and p × q, respectively, namely

$A\otimes B=\left[\begin{array}{cccc}{a}_{11}B& a_{12}B& \cdots & a_{1n}B\\ a_{21}B& a_{22}B& \cdots & a_{2n}B\\ ⋮& ⋮& \ddots & ⋮\\ a_{m1}B& a_{m2}B& \cdots & a_{mn}B\end{array}\right]\in F^{mp\times nq}$.

Let L be a constant matrix, $A\in F^{m\times n}$ be a matrix variable, and $D\left(\mathit{A}\right)\in F^{l\times l}, G\left(\mathit{A}\right)\in F^{p\times q},H\left(\mathit{A}\right)\in F^{q\times k}$ be matrix functions. The common rules for the differentiation of the above three types of matrices are as follows:

(1)

The derivative of the constant matrix L is a zero matrix, i.e.,

$$\frac{d\mathit{L}}{d\mathit{A}}=0$$

(11)

(2)

The derivative of the inverse matrix $D^{-1}$ is

$$\frac{d\left({\mathit{D}}^{-1}\right)}{d\mathit{A}}=-\left(I_m\otimes D^{-1}\right)\frac{d\left(\mathit{D}\right)}{d\left(\mathit{A}\right)}(I_n\otimes D^{-1})$$

(12)

(3)

The derivative of the matrix function product $G\left(\mathit{A}\right)H\left(\mathit{A}\right)$ is

$$\frac{d\left[G\left(\mathit{A}\right)H\left(\mathit{A}\right)\right]}{d\mathit{A}}=\frac{d\mathit{G}}{d\mathit{A}}\left[I_n\otimes H\left(\mathit{A}\right)\left]+\right[I_m\otimes G\left(\mathit{A}\right)\right]\frac{d\mathit{H}}{d\mathit{A}}$$

(13)

3.1.2. Covariance Matrix Derivation

In Equation (8), L is a constant vector, both A and $\tilde{f}$ are matrix functions of $X=\left[\mathrm{x},\mathrm{y},\mathrm{x}’,\mathrm{y}’\right]$. Combined with matrix differential theory, the covariance matrix of $\tilde{f}$ can then be derived.

Let $B=A^{\mathrm{T}}A$, Equation (8) can be expressed by

$$\tilde{f}=-B^{-1}{A}^{T}L.$$

(14)

The derivatives of $\tilde{f}$ and B are then described by Equation (15).

$$\left\{\begin{array}{c}\frac{d\left(\tilde{\mathit{f}}\right)}{d\mathit{X}}=-\frac{d({\mathit{B}}^{-1})}{d\mathit{X}}\left(I_n\otimes A^{\mathrm{T}}L\right)-\left(I_m\otimes B^{-1}\right)\frac{d({\mathit{A}}^{\mathrm{T}})}{d\mathit{X}}\left(I_n\otimes L\right)\\ \frac{d\mathit{B}}{d\mathit{X}}=\frac{d{\mathit{A}}^{\mathrm{T}}}{d\mathit{X}}\left(I_n\otimes A\right)+\left(I_m\otimes A^{\mathrm{T}}\right)\frac{d\mathit{A}}{d\mathit{X}}\end{array}\right.$$

(15)

Substitute Equation (12) into Equation (15) to obtain further expression of the derivatives of $\tilde{f}$.

$$\frac{d\left(\tilde{\mathit{f}}\right)}{d\mathit{X}}=\left(I_m\otimes B^{-1}\right)\frac{d\mathit{B}}{d\mathit{X}}\left(I_n\otimes B^{-1}{A}^{\mathrm{T}}L\right)-\left(I_m\otimes B^{-1}\right)\frac{d{\mathit{A}}^{\mathrm{T}}}{d\mathit{X}}\left(I_n\otimes L\right)$$

(16)

Moreover, $\frac{d\mathit{B}}{d\mathit{X}}$ is inserted into Equation (16) and the following equation is derived:

$$\frac{d\left(\tilde{\mathit{f}}\right)}{d\mathit{X}}=\left(I_m\otimes {\mathit{B}}^{-1}\right)\frac{d{\mathit{A}}^{\mathrm{T}}}{d\mathit{X}}\left(I_n\otimes A{B}^{-1}{A}^{\mathrm{T}}L\right)+\left(I_m\otimes B^{-1}{A}^{\mathrm{T}}\right)\frac{d\mathit{A}}{d\mathit{X}}\left(I_n\otimes B^{-1}{A}^{\mathrm{T}}L\right)-\left(I_m\otimes B^{-1}\right)\frac{d{\mathit{A}}^{\mathrm{T}}}{d\mathit{X}}\left(I_n\otimes L\right)$$

(17)

Herein, $\frac{d{A}^{\mathrm{T}}}{d\mathit{X}}$ and $\frac{d\mathit{A}}{d\mathit{X}}$ are calculated by Equation (10).

We know that X is a matrix of one row and four columns from the first paragraph of this section. Thus, m = 1, n = 4 in Equations (12) and (17). In addition, the subscript d is used to denote the differential matrix in this article; Equation (17) can be expressed by

$${\tilde{f}}_d=B^{-1}{\mathit{A}}_d^{\mathrm{T}}\left(I_4\otimes A{B}^{-1}{A}^{\mathrm{T}}L\right)+B^{-1}{A}^{\mathrm{T}}{\mathit{A}}_d\left(I_4\otimes B^{-1}{A}^{\mathrm{T}}L\right)-B^{-1}{A}_d^{\mathrm{T}}\left(I_4\otimes L\right)$$

(18)

Supposing that the variance $D_{\mathit{X}}$ of X is $\left[\begin{array}{c}\begin{array}{c}\begin{array}{cccc}{\sigma }_x^2& {\sigma }_{xy}& {\sigma }_{xx\prime }& {\sigma }_{xy\prime }\\ {\sigma }_{yx}& {\sigma }_y^2& {\sigma }_{yx\prime }& {\sigma }_{yy\prime }\\ {\sigma }_{x\prime x}& {\sigma }_{x\prime y}& {\sigma }_{x\prime }^2& {\sigma }_{x\prime y\prime }\\ {\sigma }_{y\prime x}& {\sigma }_{y\prime y}& {\sigma }_{y\prime x\prime }& {\sigma }_{y\prime }^2\end{array}\end{array}\end{array}\right]$, the variance of $\tilde{f}$ can then be obtained by Equation (19) using the guide to the expression of uncertainty in measurement method [35].

$$D_{\tilde{\mathit{f}}}={\tilde{\mathit{f}}}_d{D}_{\mathit{X}}{\tilde{\mathit{f}}}_d^{\mathrm{T}}$$

(19)

The variance of $f$ can then be calculated following Equation (20).

$$D_{\mathit{f}}=\left[\begin{array}{cc}{D}_{\tilde{\mathit{f}}}& 0\\ 0& 0\end{array}\right]$$

(20)

Some studies have pointed out [32,36] that the variance of a matrix is the same as the variance of its row-first vectors. Therefore, the covariance matrix $D_{\mathit{F}}$ of FM can also be expressed by the variance $D_{\mathit{f}}$, i.e.,$D_{\mathit{F}}=D_{\mathit{f}}$.

3.2. Scalar Function of FM Error

The 3 × 3 FM is a multivariate of nine elements, which is difficult to comprehend the FM error represented by 9 × 9 covariance matrix. Scalar functions in matrix analysis are often used to describe the one-to-one mapping relationship between matrices and real-valued functions [36]. Therefore, this paper uses the PCA method to reduce FM dimensionality and construct a scalar function to measure FM error.

According to the eigen decomposition principle, the covariance matrix $D_{\mathit{F}}$ can be divided into eigenvalues ${\mathsf{\lambda }}_i$ and eigenvectors $u_i$, where $i=1, 2, \dots ,9$. Next, the eigenvalues ${\mathsf{\lambda }}_i$ are arranged in descending order, and the first m eigenvalues are then used to calculate the cumulative contribution rate following Equation (21). Herein, the cumulative contribution rate is the information preservation degree of the original data by the newly generated components, which is usually required to be >85%. When $\eta$ is greater than 85%, the j value at this time is considered as the m value.

$$\eta \left(m\right)=\frac{{{\displaystyle \sum }}_{j=1}^m{\mathsf{\lambda }}_j}{{{\displaystyle \sum }}_{i=1}^{i=9}{\mathsf{\lambda }}_i}, m\le 9$$

(21)

The eigenvectors corresponding to the first m eigenvalues are extracted to form a transformation matrix $T^{\mathrm{T}}$:

$$T^{\mathrm{T}}=\left(u_1, u_2,\dots , u_m\right)$$

(22)

The vector $Y$, calculated by Equation (23), is composed of the most important m principal components of FM, and replaces the original nine-element FM.

$$Y={\left(Y_1,Y_2,\cdots ,Y_m\right)}^{\mathrm{T}}=T\times f$$

(23)

Herein, the $Y_1$ and $Y_2$ coordinate axis are the original direction with the largest and the second largest variance of FM, and so on. It is well known for Equation (23) that the correlation coefficients between Y elements are all 0, the scalar function Y_F can be constructed by the weight ${\alpha }_i$, where ${\alpha }_i$ is the ratio of ${\mathsf{\lambda }}_i$ to the sum ${{\displaystyle \sum }}_{i=1}^{i=9}{\mathsf{\lambda }}_i$ of the eigenvalues, as is shown in Equation (24).

$$Y_F={\alpha }_1{Y}_1+{\alpha }_2{Y}_2+\dots +{\alpha }_i{Y}_i, {\alpha }_i={\mathsf{\lambda }}_i/{{\displaystyle \sum }}_{i=1}^{i=9}{\mathsf{\lambda }}_i$$

(24)

The FM describes the mapping relationship between two uncalibrated images, and its rank is 2 [8]. Therefore, the scale of the FM value is uncertain. In order to unify the scale of FMs in a pair of stereo images, the min-max normalization method is utilized in the proposed method to convert the scalar function Y_F to the ratio R_F. Supposing that the boundary (the lower bound $Y_{\mathit{F}}^L$ and the upper bound $Y_{\mathit{F}}^U$) of Y_F can be calculated or counted, the ratio R_F can be calculated by Equation (25) to determine the quality of FM. When the ratio R_F is in the range $\left[0, 1\right]$, its corresponding FM is considered credible. Otherwise, when R_F > 1, it is considered as a low-quality FM.

$$R_F=\frac{Y_{\mathit{F}}-Y_{\mathit{F}}^L}{Y_{\mathit{F}}^U-Y_{\mathit{F}}^L}$$

(25)

Alternatively, when the lower bound of Y is known, the FM error direction can be represented by an vector arrow. In Figure 2b, the difference vector $\Delta Y$ between a certain Y and the lower bound of Y is the FM error direction in 2D space.

3.3. Review of the Proposed Method

As we know, FM represents a transformation from the coordinate system centered at C’ to the other centered at C (Figure 1). Therefore, the FM F in Equation (1) is the linear transformation function between known corresponding points, and its eigenvectors are FM directions in 2D image space (Figure 3a). The covariance matrix $D_F$ can then be derived using matrix differentiation theory, which describes the FM error direction in 9D FM space in Figure 3b; herein, the eigenvector corresponding to the maximum eigenvalue of $D_F$ is the maximum variance direction of FM error. The FM can then be converted into Y or Y_F by the PCA method. Figure 3 is the specific conversion process from nine-element FM to two-element Y in a new 2D space and to one-element Y_F in a new 1D space.

When the boundary of Y or Y_F is counted or calculated, the FM error can be represented by the ratio R_F and the difference vector $\Delta Y$. The calculation process of FM error in the proposed method is shown in Algorithm 1.

Algorithm 1. The calculated steps of FM error.

The FM can be calculated with known corresponding points on two uncalibrated images. The covariance matrix of FM can be calculated by Equation (20), and its eigenvectors can then be obtained by eigendecomposition. The m-element vector Y can be calculated by Equation (23). Here, m < 9. The scalar function YF can be calculated by Equation (24), which represents the value of the maximum variance direction of FM. When the boundary of Y or YF is known, the error value and the error direction of FM can then be measured by calculating the ratio RF and the difference vector ΔY, respectively.

When the specific ratio $R_{\mathit{F}}$ is within the range of [0, 1], its corresponding FM is considered reasonable and credible. Conversely, it may be necessary to recalculate an FM to evaluate the intrinsic projective geometry relationship between two uncalibrated images.

4. Experiment

A large number of experiments found that the cumulative contribution rate of the first two eigenvalues for FM covariance matrix, which is calculated by Equation (23), accounted for > 85%. Therefore, the first two principal components Y₁ and Y₂ are used to form the vector Y = (Y₁, Y₂) and to construct the scalar function Y_F. Therefore, the Y and Y_F could be calculated mainly by Equation (26) in the following experiments.

$$\left\{\begin{array}{c}{Y}_{\mathit{F}}={\lambda }_1{Y}_1+{\lambda }_2{Y}_2\\ Y_1=u_{1}f\\ Y_2=u_{2}f\end{array}\right.$$

(26)

Herein,$f$ is the nine-vector of FM, ${\lambda }_1$ and ${\lambda }_2$ are the largest and sublargest eigenvalues decomposed by the covariance matrix D_F. Additionally, $u_1$ and $u_2$ are eigenvectors corresponding to the eigenvalues ${\lambda }_1$ and ${\lambda }_2$, respectively.

4.1. Data Source

Three groups of stereo images in this work were selected with different image acquisition methods for testing, and their corresponding points were extracted using SURF and the nearest-neighbor search algorithms. The test1 (Figure 4(a1,b1)) was with normal close-digital stereo imagery of the Tsinghua School, Beijing, China, which was published by the Institute of Automation of the Chinese Academy of Sciences, Beijing, China. In addition, 835 pairs of corresponding points are shown in Figure 4(c1) with + symbols, of which about 250 pairs are outliers. The test2 (Figure 4(a2,b2)) was actual digital drone images, which were taken by measuring using a SONY ILCE-6000 digital camera at an altitude of 300 m. Consequently, 1962 pairs are shown in Figure 4(c2), of which about 550 pairs are outliers. The test3 (Figure 4(a3,b3)) are actual digital aerial images, which were taken using the digital aerial UltraCam camera at an altitude of 3600 m. Figure 4(c3) shows 3656 pairs and 150 pairs are outliers.

4.2. Correctness of Scalar Function

It is often desirable to use an element or scalar function to refer to multiple variables, in which the norm is an important method to study least-squares solutions and matrix perturbations [36]. There are many types of norms, e.g., the Euclidean norm, row addition norm, column addition norm, spectral norm, and the Mahalanobis norm, in which Euclidean and Mahalanobis norms are metric functions that measure matrix similarity. Previous research has pointed out that the Manhattan distance in matrix theory could be used to evaluate FM error [32] and was expressed as

$${\parallel \mathit{F}\parallel }_{{\mathit{D}}_{\mathit{f}}}=\sqrt{\mathrm{tr}\left(F^{\mathrm{T}}\times D_{\mathit{f}}\times F\right)}$$

(27)

Therefore, this paper uses the Mahalanobis norm ${\parallel \mathit{F}\parallel }_{{\mathit{D}}_{\mathit{f}}}$ as the truth value to verify the correctness of the scalar function Y_F.

4.2.1. Experiment Process

The Mahalanobis norm computed by Equation (27) can be used as the truth value to verify the derived value in Equation (26). Figure 5 is a verification flowchart and its verification process is as follows:

1.

The RANSAC algorithm is run n times to obtain n R_i sets of corresponding points after eliminating outliers, with which the covariance matrix D_Fi can then be calculated by Equation (20), where $i=1, 2, \dots ,n$.

2.

The Mahalanobis norm ${\parallel \mathit{F}\parallel }_{{\mathit{D}}_{\mathit{f}}}$ can be calculated following Equation (27), which is considered the truth value of FM error.

3.

The scalar function Y_F can also be calculated using Equation (26), which is a measured value of FM error.

4.

The correlation of the truth and measured errors can be calculated. When the correlation coefficient is large enough, the scalar function derived in this article is considered correct.

4.2.2. Experiment Result

The RANSAC algorithm was run 100 times to obtain 100 R_i sets using known corresponding points in Figure 4(c1,c2,c3), and the covariance matrix D_Fi then was calculated by Equation (20), here,$i=1, 2, \dots ,100$. The scalar function Y_F and Mahalanobis norm ${\parallel \mathit{F}\parallel }_{{\mathit{D}}_{\mathit{f}}}$ of FM were then calculated following Equations (26) and (27), and their distribution (the first and second) and correlation (the third) scattergrams are shown in Figure 6.

A linear relationship exists between the scalar function Y_F and the Mahalanobis norm ${\parallel \mathit{F}\parallel }_{{\mathit{D}}_{\mathit{f}}}$. By linear fitting, the correlation coefficients are 0.951, 0.935, and 0.933 in tests 1, 2, and 3, respectively. Therefore, Y_F and ${\parallel \mathit{F}\parallel }_{{\mathit{D}}_{\mathit{f}}}$ have great correlation, >90% in all three tests. The Mahalanobis norm as the truth error can represent the difference of FMs, thus, the scalar function Y_F derived in Equation (26) can also be utilized in this study.

4.3. Comparison with the Existing Method

In the past few years, although the 3D metric distance was proposed as a new error criterion to estimate the FM quality, it also has a similar performance in accuracy compared to the symmetric epipolar distance [23]. Therefore, the symmetric epipolar distance, as the most popular error criterion, is chosen to verify the proposed method. Specifically, the symmetric epipolar distance E_d is the squared distance between a point’s epipolar line and the corresponding point in the other image (computed for both points of the match), averaged over all N corresponding points.

$$E_d=\frac{1}{N}{\displaystyle \sum }_{i=1}^N\left(d{\left(x_i^{\prime },F{x}_i\right)}^2+d{\left(x_i,F^{\mathrm{T}}{x}_i^{\prime }\right)}^2\right)$$

(28)

where $d\left(x_i,F^{\mathrm{T}}{x}_i^{\prime }\right)$ is the distance (in pixels) between a point $x_i$ and a line $F^{\mathrm{T}}{x}_i^{\prime }$.

4.3.1. Experiment Process

The RANSAC algorithm can be run multiple times to obtain many FMs. Y_F and $E_d$ could then be calculated with Equations (26) and (28), the correlation of which could also be fitted by the linear regression method. The comparison flowchart (Figure 7) of Y_F and $E_d$ is designed as follows:

1.

The RANSAC algorithm is run n times to obtain n R_i sets of corresponding points after eliminating errors, with which the covariance matrix D_Fi can then be calculated by Equation (20), where $i=1, 2, \dots ,n$.

2.

The symmetric epipolar distance $E_d$ can be calculated following Equation (28), which is considered as the compared value of FM error.

3.

The scalar function Y_F can be calculated using Equation (26), which is a measured value of FM error.

4.

After the correlation of Steps 2 and 3 is calculated, the difference between the proposed method and the existing method can be analyzed.

4.3.2. Experiment Result

The RANSAC algorithm was run 100 times to obtain 100 FMs using known corresponding points in Figure 4(c1,c2,c3). Following this, 100 Y_F and 100 $E_d$ values were calculated by Equations (26) and (28), respectively. Figure 8 shows the distribution and correlation scattergrams of Y_F and $E_d$.

In this experiment, three groups of scattergrams in Figure 8(a1,a2,a3) were calculated with known corresponding points in Figure 4(c1,c2,c3), and in Figure 8(b1,b2,b3) with the corresponding points after removing the false points. Comparing Figure 8(a1,a2,a3) with Figure 8(b1,b2,b3), it is found that the linear distribution in Figure 8 (b1,b2,b3) is a part in Figure 8(a1,a2,a3), as the rectangle in Figure 8(a1,a2,a3) points from Figure 8(b1,b2,b3). By linear fitting, the correlation coefficients between Y_F and E_d are 0.274, 0.126, and 0.861 in Figure 8(a1,a2,a3), and 0.842, 0.839, and 0.946 in Figure 8(b1,b2,b3). Relative to Figure 8(b1,b2,b3), the correlation coefficients of Figure 8(a1,a2,a3) are all lower, especially those of Figure 8(a1) and Figure 8(a2). The specific reason is as follows: There are about 250 pairs in test1, about 550 pairs in test2, and about 150 pairs of false matches in test3, and they, respectively, account for about 30%, 28%, and 4% of all corresponding points in Figure 4(c1,c2,c3). Tests 1 and 2 in Figure 8(a1,a2) contain many false corresponding points, which makes it difficult to eliminate all the false points using the RANSAC algorithm, so that the R set in Figure 7 contains several false corresponding points, calculating an outlier of E_d, and resulting in less correlation between Y_F and E_d.

Therefore, false corresponding points affect the linear relationship of Y_F and E_d. When there are many false points in an R set, the E_d calculated by Equation (28) may be wrong or a large error, resulting in the incorrect estimation of FM. However, the covariance $D_F$ of FM in this article is calculated with the least-squares principle, which ignores the influence of the false corresponding points on the whole. Therefore, compared with the symmetric epipolar distance, the proposed method relies less on the accuracy, number, and distribution of known corresponding points, and is more robust and more reliable for evaluating FM error.

4.4. Application of This Method

With this proposed method, the FM and its covariance matrix can be calculated with known corresponding points on uncalibrated stereo images, and the vector Y and scalar function Y_F can be calculated by Equation (26). When the boundary of Y_F or Y can be counted or calculated for a group of stereo images, the ratio $R_F$ and the vector arrow $\Delta Y$ can be calculated to estimate the error value and the error direction of a certain FM in this group of stereo images.

Figure 9 shows the statistical boundaries (broken lines) of Y_F and Y based on three groups of stereo images in Figure 4. In this experiment, the RANSAC algorithm was run 100 times to obtain 100 Y and 100 Y_F values, most of which are clustered together and form the normal distribution and the error ellipse. For a certain FM estimated with some corresponding points in Figure 4(c1,c2,c3), the Y and Y_F can be calculated.

With the statistical boundaries (broken lines) of Y_F and Y in Figure 9, the ratio $R_F$ and the vector arrow $\Delta Y$ can then be calculated. When the Y_F or Y falls into the normal distribution curve or inside the error ellipse, its ratio $R_F$ is less than 1 and its vector arrow is inside the error ellipse, as shown in the red dot A. This means that the FM corresponding to the red dot A falls into the normal distribution and the error ellipse, and is considered to have less error and be more credible. On the contrary, the red dot B falls outside the normal distribution and its vector arrow goes beyond the error ellipse. We can say that the FM corresponding to the red dot B may be low quality, and it is necessary to recalculate an FM to evaluate the intrinsic projective geometry relationship between two uncalibrated images.

5. Discussion

FM represents the intrinsic projective geometry between two image planes, and is a multivariate matrix of nine elements, in which is difficult to measure the FM error represented by 9 × 9 covariance matrix. The PCA is the most widely used method to reduce the dimensionality of such datasets. Therefore, this study used the PCA method to reduce FM dimensionality and constructed the scalar function Y_F to measure FM error.

In matrix theory, the Manhattan distance as a standard can be used to measure the difference of FMs. This study used the Mahalanobis norm ${\parallel \mathit{F}\parallel }_{{\mathit{D}}_{\mathit{f}}}$ as the truth value to verify the correctness of the scalar function Y_F. The experiment in Section 4.2 found that there was a linear relationship between Y_F and ${\parallel \mathit{F}\parallel }_{{\mathit{D}}_{\mathit{f}}}$, and there was a high linear correlation in three groups of tests. The Mahalanobis norm, as the truth value, can represent the difference of FMs, the scalar function Y_F derived in Equation (26) can also be used to measure the FM error.

The symmetric epipolar distance is the most popular error criterion. This study analyzed the relationship between the scalar function Y_F and the symmetric epipolar distance E_d in Section 4.3. Compared with tests 1 and 2, test 3 had a higher linear relationship in Figure 8. The reason may be that there are many false corresponding points in tests 1 and 2, the E_d calculated with these corresponding points may be wrong or involve a large error, resulting in the incorrect estimation of FM. Therefore, false corresponding points affect E_d and the linear correlation coefficient of Y_F and E_d. Compared with E_d, the proposed method relies less on the accuracy, location, and number of known corresponding points, and is more robust and more reliable to evaluate FM error.

It is worth noting that both the Mahalanobis norm ${\parallel \mathit{F}\parallel }_{{\mathit{D}}_{\mathit{f}}}$ and the symmetric epipolar distance E_d are error criteria greater than 0. When ${\parallel \mathit{F}\parallel }_{{\mathit{D}}_{\mathit{f}}}$ or E_d is larger than 0, it means the FM error is larger. However, the experiments in this article found that the ranges of Y_F calculated in the three tests were different in Figure 6 and Figure 8, in which some Y_F values were less than 0. Although Y_F is linear with respect to ${\parallel \mathit{F}\parallel }_{{\mathit{D}}_{\mathit{f}}}$ or E_d, the minimum value of Y_F is uncertain. Therefore, it is difficult to measure the FM error for a specific Y_F. In order to solve this problem, the min-max normalization method is utilized to transform the scalar function Y_F to the ratio R_F. In Section 4.4, the boundary of Y_F or Y can be counted using the different corresponding points of the same stereo images multiple times. As is shown in Figure 9, the boundary of Y_F is a normally distributed curve in 1D space, and the boundary of Y is an error ellipse in 2D space. When the Y_F or Y corresponding to a specific FM falls into a normal distribution curve or inside an error ellipse, the ratio R_F is less than 1, and the FM corresponds to the Y_F, Y and R_F are then considered to have less error and be credible. Otherwise, the FM corresponding to the ratio R_F > 1 may be low quality.

In summary, three groups of experiments verified the rationality of the proposed method and obtained results as follows: (1) The proposed method is highly fitted to the Mahalanobis norm of FM and can be used to replace the nine-element FM. (2) Compared with the symmetric epipolar distance, the proposed method relies less on the accuracy, location, and number of known corresponding points, and is more robust and more reliable in evaluating FM error. (3) Through statistical analysis, it was found that the FM error is normally distributed and forms an error ellipse. When counted or calculated, the boundary of Y_F or Y, the ratio R_F, and the vector $\Delta Y$ can be calculated to estimate the error value and error direction of FM. When the ratio ${\mathrm{R}}_F$ is within [0, 1], the FM corresponding to $R_F$ is considered reasonable. Conversely, it may be a low-quality FM.

6. Conclusions

This study accomplished two main objectives: one was to derive the covariance matrix by matrix differential theory and covariance-propagation theory, and the second was to use the PCA method to construct the scalar function to estimate FM error. Specifically, the FM error in this paper was estimated as follows:

• Compute the covariance matrix using the known corresponding points extracted from two uncalibrated images.
• Following the PCA method, decompose the covariance matrix and then construct the vector Y and the scalar function Y_F.
• Count or calculate the boundary of Y and Y_F, and then calculate the ratio R_F to estimate FM error_. Here, the FM corresponding to $R_F\in \left[0, 1\right]$ is considered reasonable. When R_F > 1, it is considered a low-quality FM.

The main significance of the proposed method is to determine whether a specific FM falls into a normal distribution or inside an error ellipse by calculating the ratio R_F. When the R_F is in the range [0, 1], FM is then considered to be inside a normal distribution or an error ellipse, and it is credible to describe the intrinsic projective geometry between the two uncalibrated images. Conversely, it may be necessary to recalculate an FM for computing camera parameters, making adjacent image connections, and reconstructing 3D models.

The boundary of the error criterion in this paper was obtained by statistical methods in experimental application. They will be studied and calculated by theoretical methods in the future.