2.1. Projection Models

A projection model is a mapping from the 3D object point P(X,Y,Z) to its corresponding 2D image point p(x,y). A pinhole camera is always described by the perspective projection model. If we denote the focal length of the camera as f, the perspective projection could be formulated as

In addition to the projection model, imaging system design usually involves many other considerations, e.g., size, cost, and uniform illumination. Assembly defects are also introduced during camera manufacture. These will make the actual output of cameras deviates from the ideal projection model to varying degrees. This deviation is often treated as distortion. In classical camera calibration methods, the most commonly used distortion model is the Brown–Conrady polynomial model [21]:

Here, [xd,yd]T is the actual position of the image point, [xu,yu]T is the position of the corresponding point without distortion, and r denotes the distance from [xu,yu]T to the distortion center c=[xc,yc]T, which is also the principle point here. Further, kn is the nth radial distortion coefficient. When the distortion is small, usually two coefficients k1 and k2 are sufficient. We can say that the distortion of a camera is a mapping from [xu,yu]T to [xd,yd]T and can be characterized by a set of parameters λ={xc,yc,k1,k2}.

The FOV of a pinhole camera is restricted by the size of the sensor. To take the entire scene of a 180∘ FOV into a limited sensor, fish-eye lenses are designed to severely bend the incident ray coming from an ultra-wide angle and therefore cannot be described properly by the linear model given by (1). A two-stage procedure using a sphere (say, a viewing sphere), whose center coincides with the optical center of the camera and whose radius is equal to the camera's focal length, is introduced to describe a nonlinear projection model [22]. The procedure is demonstrated in Figure 1.

Figure 1

Two-stage procedure for describing a projection model. A 3D point is first projected onto the viewing sphere, followed by another perspective projection from some point on the optical axis onto the image plane.

First, a 3D point P is projected perspectively through the center O onto point Ps on the sphere. Then Ps is projected to a 2D point p on the image plane π from another point O′ on the optical axis ON. The position of O′ along ON defines a series of projection models. Here O′ coincides with N in Figure 1. A projection model can be clearly represented as a mapping r=ρ(α), where α is the incident angle of ray emitted from P, and r is the radial distance as mentioned. The representative models below have been widely discussed:

Note that the perspective projection can also be described by this procedure. In Ref. [23], the author prefers to stereographic projection rather than other alternative ones because it can represent a wider FOV while preserving a suitable range of shape properties of the projected objects, including the circularity, the local symmetry, the shape of small regions, and the angle at which two curves intersect. These properties are very useful in applications such as intelligent surveillance and automatic navigation. Thus, in this paper, the stereographic projection model is chosen for our estimation method. It is exactly the model shown in Figure 1.

2.2. Geometric Invariants Under Stereographic Projection

Simple geometric targets, such as straight lines or spheres, are often utilized in the distortion estimation methods because the shape of the images of these targets is known under a certain projection model. For stereographic projection, this is shown in Figure 2.

Figure 2

Images of a line and a sphere under the stereographic projection. The projection of a space line on the viewing sphere, as well as the projection of a sphere (dotted line), is a circle. Circles on the viewing sphere maintain their shape when projected onto the image plane from the northern pole of the viewing sphere (solid line on the image plane).

When projected onto the viewing sphere from its center O, the image of a space line L is a great circle determined by the sphere surface and the plane OL. On the other hand, the image of a sphere S on the surface of the viewing sphere is a small circle. Let (nx,ny,nz,d0) be the plane containing the small circle. Here, (nx,ny,nz) is the normal vector of the plane, and d0 is the distance from O to the plane. Obviously, the great circle is actually a special case of the small circle in which d0 is zero. When projected from N onto the image plane π, the function of a point on the image of S can be derived as

We set d0=0 to obtain the formula of the image of L:

Notice that both (4) and (5) are the general equations of a circle. So the conclusion is that the image of a space line section under the stereographic projection is a circular arc (a portion of a circle), and the image of a sphere is a circle. These are the main geometric invariants used in our estimation method.

3.1. Generate the Ideal Fish-Eye Image

The common method to generate images of a virtual scene is backward ray-tracing. For fish-eye camera simulation, it is the inverse procedure of the projection model described in last section. For a given 2D point p on π, a direction vector starting from O is found. A ray in this direction is traced until it hits a surface in the scene at a point P. The value of p can be computed according to the properties of P (material, color, orientation, lighting, etc.). Let the polar coordinate of point p is denoted by (r, ϕ), where r is the radial distance, and ϕ is determined by

From the stereographic projection function given in (3), the incident angle α of the light ray corresponding to p is

Actually, in order to achieve a more realistic rendering effect, the tracing may not stop at P. It will go between objects in the scene until it reaches a light source. This is a tedious work, and difficult to implement. So we let the 3D engine to do this job. This is a so called postprocessing method, because it is based on the prerendered images.

We define 6 identical pinhole cameras in the virtual scene, all positioned at the center of the viewing sphere. The cameras are arranged orthogonally. If the orientation of the fish-eye camera is (0, 0, 1), then the orientations of the 6 pinhole cameras are (0, 0, 1), (0, 0, −1), (0, 1, 0), (0, −1, 0), (1, 0, 0), (−1, 0, 0), respectively. The cameras all have an FOV of 90 degrees in both horizontal and vertical directions. The deployment of the cameras is shown in Figure 3.

Figure 3

The viewing sphere and the deployment of the six pinhole cameras. The cameras are moved from their original position (i.e., the center of the viewing sphere) along the axis to make the figure clearer.

In this way, the images rendered through these six cameras form a cube-map. In our simulation, only five cameras are sufficient. The camera looking backward is omitted because we just need an FOV of 180 degrees in front of our fish-eye camera. OP→ is transformed into the coordinate frames of the five cameras to see if it hits one of the image planes at a point m(xm,ym), as shown in Figure 4.

Figure 4

Mapping from the fish-eye image plane to the five orthogonal cameras. The incident light ray is back-traced and hits one of the faces of the cube-map at a point.

This creates a mapping between the simulated fish-eye camera and the five cameras. This mapping doesn't change from frame to frame in real-time rendering, thus can be stored in a mapping table. Each camera is assigned an ID, and for a pixel in the output image, its corresponding entry in the lookup table is like (ID,xm,ym). The pixel value can then be obtained from the position in the pointed image. The rendered outputs of the five cameras and the simulated ideal fish-eye image is given in Figure 5.

Figure 5

The rendered outputs of the five cameras and the simulated ideal fish-eye image.

3.2. Image Distortion

Figure 6(a) shows the ideal image of several spheres generated by the method described in last section. According to the property of the stereographic projection, the image of a sphere is a circle. However, it is not the case in a real fish-eye image as shown in Figure 6(b).

Figure 6

(a) The simulated ideal fish-eye image of several spheres. (b) Image of spheres captured by a real fish-eye camera.

In this section, we detail our distortion estimation method. The estimation of the distortion center and the distortion coefficients are performed in a hierarchically iterative manner because it has been shown that including the distortion center in the nonlinear optimization may make the system unstable with the presence of noise [21]. In the outer loop, feature points extracted from images of parallel lines are undistorted with the latest estimated distortion parameters and then used to update the distortion center. If the new center is close enough to the current one, the algorithm terminates. Otherwise, in the inner loop, with the latest distortion center, a genetic algorithm (GA) module is adopted to optimize the distortion coefficients based on the geometric invariants discussed in Section 2.2.

3.2.1. Estimation of the distortion center

We estimate the distortion center with the help of the vanishing points. A vanishing point is the intersection of a group of parallel straight lines under the stereographic projection model. A straight line in space can be given as a point direction form equation as

L(t)=qt+b,(9)

where q=[qX,qY,qZ]T is the direction vector, and b=[bX,bY,bZ]T is a fixed point on the line. Further, t is a parameter indicating any given point on the line. By combining (1) and (9), we have the projection of a straight line on the perspective image plane:

x(t)y(t)=−fqZt+bZqXt+bXqYt+bY.(10)

Let r(t) denotes the radial distance from [x(t),y(t)]T to the distortion center. If we assume that the distortion center coincides with the origin [0,0]T, then we have r(t)=x2(t)+y2(t).

From the perspective and stereographic projection functions in (3), the radial distance of the corresponding point [x′(t),y′(t)]T under the stereographic projection is

r′(t)=2ftanarctanr(t)f2.(11)

Rearranging gives

r′(t)=2ff2+r2(t)−fr(t).(12)

The projection of a given 3D point under different projection models always moves along the radial direction. Therefore, the projected 2D points are related by their radial distance:

x′(t)y′(t)=r′(t)r(t)x(t)y(t).(13)

Inserting (10) and (12) into (13) gives

x′(t)y′(t)=1−D(qXt+bX)2+(qYt+bY)2⋅2f⋅qXt+bXqYt+bY.(14)

where D=(qXt+bX)2+(qYt+bY)2+(qZt+bZ)2.

The limit of (14) as t→∞ then gives two vanishing points, v1 and v2, of a group of parallel lines projected onto the stereographic image plane, i.e.,

xvyv=limt→±∞x′(t)y′(t)=1∓qX2+qY2+qZ2qX2+qY2⋅2f⋅qXqY,(15)

If we consider the two points as two vectors from the distortion center, examining the angle ϕ between these two vectors gives

ϕ=arccosv1⋅v2v1v2=π,(16)

which means that the two vanishing points lie in opposite directions from the distortion center. Thus, the line joining v1 and v2 must pass through the distortion center. The position of the distortion center is determined as follows:

1. Extract points. An edge detection algorithm such as Sobel or Canny is used to extract the points on the image of lines. As shown in Figure 7(a), the extracted points are grouped into several sets:

   PL=lji=(xji,yji)j=1mii=1n,
   where for each i=1,…,n the set {lji}j=1mi is a collection of extracted points belonging to the image of the ith space line.

2. Undistort the points and do circle fitting. PL is undistorted according to (2) with the currently available distortion parameters, if they exist. Circles are fitted to each set of the resulting points P¯L in a least square manner. Let the circle fitted to the point set {l¯ji}j=1mi be Ci. Ideally, all {Ci}i=1n will converge exactly at a pair of vanishing points. However, this is not the case because of the error introduced by the distortion and measurement noise, as shown in Figure 7(b). Thus, the estimation of the vanishing points needs a nonlinear optimization.

3. Initialize the vanishing points. An initial guess is required for the optimization. The intersections of every two circles in {Ci}i=1n are calculated to be the candidates of initial vanishing points. All these candidates are examined to pick the best one. To do this, we do circle fitting again to P¯L. For a pair of candidate vanishing points (v̂1,v̂2) that is being examined, the circle Ĉi fitted to the point set {l¯ji}j=1mi is forced to pass through both of these two points. Ĉi is then adjusted to minimize the error function:

   ϵi=∑j=1mix¯ji−x0i2+y¯ji−y0i2−Ri2,(17)
   where [xi0,yi0]T and Ri are the center and radius of Ĉi respectively. Summing the minimized fitting errors calculated from every point set in P¯L gives
   E=1∑i=1nϵi.(18)

   E is taken as the score of the pair (v̂1,v̂2). The pair with the highest score is chosen as the initial position of the vanishing points.

4. Optimize the vanishing points. The initial positions of v1 and v2 are then refined with the Levenberg–Marquardt algorithm. Eq. (18) is used again as the objective function. To ensure the stable convergence of the optimization, the two vanishing points are refined alternately. When v1 is being refined, the position of v2 is fixed. Once v1 converges, v1 is fixed and v2 is refined. This repeats until both v1 and v2 converge. Figure 7(c) shows the final position of v1 and v2.

5. Estimate the distortion center. Several other pairs of vanishing points are determined from other pictures of the parallel lines in the same way. The lines joining each pair of vanishing points intersect with each other, and the centroid of the intersections is located as the distortion center, as shown in Figure 7(d).

Figure 7

Estimation of the distortion center: (a) Points on the images of lines are extracted and grouped into several sets. (b) The fitting circles have many pairs of intersections, and the one with the lowest score is chosen as the initial position of the vanishing points. (c) The initial points are optimized and then linked. (d) At least three pairs of vanishing points need to be estimated to determine the distortion center.

3.2.2. Estimation of the distortion coefficients

In this subsection, the geometric invariants discussed in Section 2.2 are used to estimate the distortion coefficients of the fish-eye cameras. We need to define an object function by which, when minimized, the distortion parameters are optimized to put the points extracted from the image of the geometric targets back onto given shapes. Our object function is defined via the images of spheres as well as straight lines. As we will show later in Section 4, spheres are more appropriate to serve as the geometric target for the distortion estimation of fish-eye cameras, especially with the presence of measurement noise.

The first two steps are similar to the ones described in last subsection. Pictures of the target objects (lines in space or spheres) are captured, and edge points Ptar are extracted. Ptar is then undistorted with a set of distortion parameters λ̂={x̂c,ŷc,k̂1,k̂2}, and circles are fitted to the resulting points P¯tar.

Measurement noise in the extracted points will be amplified nonlinearly due to the distortion if the objective function is defined in the space of the undistorted points [24]. Thus, in this paper, we define the object function in the space of distorted image points to ensure the robustness of the estimation method with the presence of noise. As shown in Figure 8, pji is one of the points in Ptar, p¯ji is the corresponding point in P¯tar, and Ci is the circle fitted to {p¯ji}j=1mi. We search for a point p̂ji near pji which, after being undistorted with λ̂, will lie exactly on Ci at p̂ji. Then the object function is defined as follows:

e=∑i=1n∑j=1mipji−p̂ji2.(19)

Figure 8

Setup for determining objective function.

The GA is used here to minimize (19). A simple GA using single-point crossover and single-point mutation is sufficient for our problem. Any other nonlinear optimization algorithm could also be used here. Once the distortion parameters are determined, the distorted position of the pixels in the ideal image can be found, and the mapping relationships between the output image and the cube-map recorded in the mapping table need to be adjusted accordingly.

4.1. Experiment on Fully Synthetic Points

In this experiment, we show how the number of the target objects affects the estimation result with the presence of noise of different levels. The effectiveness of lines and spheres as the target is also compared. A fish-eye camera using an 800*800 CCD sensor is considered. With a focal length of 160 pixels, the whole FOV of 180∘ can be fitted within the sensor under the stereographic projection. The ideal distrotion-free image of parallel space lines and spheres are shown in Figure 9(a). We assume that all the lines in space are parallel to the image plane, which will not affect the objectiveness of the result, so the vanishing points are at the angular position of 90∘, i.e. on the edge of a FOV of 180∘. The images of spheres is relatively simple, only circles with random positions and radii. Points evenly sampled on these arcs, denoted as PLG and PSG, are grouped like PL.

Figure 9

Illustration of the experiment using fully synthetic data. (a) Points on the ideal image. (b) Points in (a) are distorted. (c) Distorted points are contaminated by Gaussian random noise.

Next, PLG and PSG are distorted with known parameters to obtain PLD and PSD, as shown in Figure 9(b). Here, for simplicity, we let the distortion center coincide with the center of the image, i.e. [xc,yc]T=[0,0]T. Figure 9(c) shows that an independent, zero mean, zero-correlated, two-dimensional Gaussian noise is added to PLD and PSD. The points with noise are denoted as PLD′ and PSD′.

In a single run of the estimation, a group of three PLD′s is generated to serve as PL which is used to estimate the distortion center. A PSD′ or another PLD′ is generated to serve as Ptar. In order to make the experiment more rigorous, we generate another PSD(notice, not PSD′ here) other than that used in the estimation to evaluate the estimation result. This PSD is undistorted by the estimated distortion parameters to get P¯SD. The average distance d between P¯SD and its corresponding ground truth PSG can intuitively reflect the accuracy of our method. So a test set used here is PLD′∗3,PSD′(orPLD′),PSD.

We vary the noise level w in increment of 1 pixel starting from 1 pixel, and vary the number of the targets n under different noise levels. For each (w,n) pair, the proposed method is performed on 100 random test sets. Figure 10 reports the inter-quartile range (IQR) of d when PSD′ is used as Ptar.

Figure 10

The inter-quartile range (IQR) of is fixed at 0.

Although becoming larger with the increasing noise level, d is controlled no more than 3.5 pixels, even when w reaches up to 5 pixels. When w is lower than 2 pixels, our method achieves sub-pixel accuracy. We can see that when more spheres are involved, the performance of the method is enhanced. Both the accuracy and the stability of the algorithm are increased when the number of the spheres reaches 6. The distortion is a polynomial of the radial distance r. The estimation of the distortion parameters is essentially a process of polynomial fitting. More targets will provide more image points which can cover a larger range of r. This improves the stability of the polynomial fitting. These results demonstrate that our method is robust to noise when the spheres are employed.

On the other hand, when a group of lines is used as the target, the result is not so good, as show in Figure 11. Even if the number of the lines increases, the improvement of the result is limited. The main reason for this is that the image of a space line under stereographic projection is only a portion of a circle, and circle fitting with points lying on a portion of a circle is an ill-conditioned problem. When a circle is fitted to these points, owing to the stronger fitting capacity of quadratics and the lack of geometric constraints, the objective function (19) has a great chance of reaching the minimum even if the undistorted points are very far from the ground truth. This will cause the iteration to converge to the wrong parameters. However, the image of a sphere is itself a complete circle; thus, the problem described above will not exist when spheres are used.

Figure 11

The inter-quartile range (IQR) of is set to 1.

4.2. Experiment on Simulated Images

The complete process of the method is carried out in this experiment, including the extraction of the image points. The experimental setup is the same as that described in Section 4.1 except that the ground truth of the points doesn't exist. The parameters of the distortion are still known. Therefore, in this section, we directly examine the estimated parameters. To obtain the simulated images, two virtual scenes are created by the OGRE engine, as shown in Figure 12(a) and 12(b).

Figure 12

(a) and (b) Scenes containing a plane textured by black-and-white stripes and some spheres, respectively. (c) and (d) Distorted fish-eye images of (a) and (b), rendered by the modified OGRE engine.

Scene 1 contains a plane target textured by some stripes to provide parallel lines. Scene 2 contains several identical spheres placed at random positions. We don't care about the interval between the lines and the size of the spheres. We have modified the engine to render blurry images. A Gaussian Blur with a radius of 2 pixels is applied here, and the coefficients k1 and k2 are set to 3E−6 and 6E−13, respectively. Three images of scene 1 and one image of scene 2 comprise a test set. The proposed method is performed 100 times on a single test set to investigate its stability. The estimated distortion center is c=−0.047±0.89,−0.036±0.88T, which is relatively stable. Figure 13 shows the results of the distortion coefficients.

Figure 13

(a) and (b) Results of 100 runs of our method with the same test set. Both is fixed at 0.

The result shows that when both of the coefficients are estimated, the value of k2 fluctuates, and so does k1. On the other hand, when k2 is fixed to 0, the stability of the estimation is enhanced. This explains why the stability of the algorithm is improved when only k1 is estimate, especially when the noise level increases, as shown in Figure 10. Thus, it is recommended that k2 should be omitted when an actual estimation is performed.

4.3. Experiment on Real Images

In this subsection, we give a quantitative evaluation of our simulation method. The error between the relative distortion of the real and the simulated image is measured.

The device used in the experiment is an inexpensive fish-eye camera. The size of its sensor is 704*576 in pixels. The distortion parameters are estimated from the pictures of the targets, as shown in Figure 14.

Figure 14

Images captured by a real fish-eye camera used to estimate the distortion.

The estimation result is c=25.82,−6.83T and (k1,k2)=(−1.61E−6,2.5E−13).

A chessboard pattern is captured with the target device, and the corresponding simulated image is given, as shown in Figure 15.

Figure 15

An image of a chessboard pattern captured by the real fish-eye camera and the corresponding simulated image.

Corners are extracted from both of the images. The ratio between the radial distance of the corners in the observed image and their distortion-free positions is a relative measure of distortion. The similarity between the real and the simulated image can then be determined by comparing this ratio, as shown in Figure 16.

Figure 16

The ratio between the ideal and distorted position of the corners in the real and simulated image and the relative error calculated between the two images.

For points near the image center (radial distance smaller than 60%), the relative error is very small. The error increases with the radial distance, which could have been expected due to the limited number of terms used in the distortion polynomial. At the edge of the images, the error between the real and the simulated images is about 1.2%. This is a good enough result.

4.4. Validation of Practicality

In this subsection, the practicality of our method is validated through an actual application of vehicle omniview system. In a common vehicle omniview system, there are often four fish-eye cameras mounted: one in the front of the vehicle, one in the back of the vehicle, one in the left outside mirror, and one in the right outside mirror. These four simulated fish-eye images are shown in Figure 17(a). The imaging system simulated above is used here. The ground is textured with chessboard pattern to clarify the result. The inverse perspective transformation is adopted to transform the image parts containing the ground into top views according to the parameters of the cameras. The top views are then fused into one single omniview image, as shown in Figure 17(b).

Figure 17

(a) The simulated images from the four fish-eye cameras mounted on the vehicle. (b) The omniview image generated from the images in (a).

As we can see, compared with the sides of the vehicle, the top views of the ground in front and rear of the vehicle are more blurred and mismatched. The main reason of this is the mounting angle of the cameras. Both of the side cameras look straight downward, while the cameras on both ends of the vehicle raise their line of sight in order to avoid being occlude by the bumpers. The image parts of ground in these images are closer to the edge, where the distortion is more severe. The quality of the top view transformed from these images is poor. To solve this, we can modify the algorithm, change the camera mounting position, or even consider other cameras. In this way, problems can be found before the assembling of the actual system and the cost is saved.