Self-calibration using two same circles

Self-calibration using two same circles

Optics & Laser Technology 44 (2012) 1924–1933 Contents lists available at SciVerse ScienceDirect Optics & Laser Technology journal homepage: www.els...

807KB Sizes 2 Downloads 57 Views

Optics & Laser Technology 44 (2012) 1924–1933

Contents lists available at SciVerse ScienceDirect

Optics & Laser Technology journal homepage: www.elsevier.com/locate/optlastec

Self-calibration using two same circles Feipeng Da n, Qin Li 1, Hu Zhang 2, Xu Fang 3 Southeast University, School of Automation, No. 2, Si Pai Lou, Nanjing 210096, PR China

a r t i c l e i n f o

a b s t r a c t

Article history: Received 2 November 2011 Received in revised form 14 February 2012 Accepted 14 February 2012 Available online 3 March 2012

Camera calibration is a fundamental step in 3D reconstruction and computer vision. Considering the easy accessibility of two same circles, a method for self-calibration based on two same circles is proposed. The proposed method does not need any prior knowledge and known camera parameters. By taking three photos of the object containing two same circles from different views, the close solution for intrinsic and extrinsic parameters are first obtained by the invariance of circular points and common tangent of two circles. Then the solution is refined by the nonlinear optimization. This method could get the intrinsic parameters as well as the extrinsic parameters without complicated matching process. Extensive results show the accuracy, robustness and wide applications of the proposed method. & 2012 Elsevier Ltd. All rights reserved.

Keywords: Camera calibration Self-calibration Circular points

1. Introduction Camera calibration is a significant task in computer vision since the intrinsic parameters and extrinsic parameters of the camera are indispensable for 3D reconstruction. Some traditional calibration methods based on circular target are proposed due to its accuracy and robustness [1–3]. Compared with the traditional calibration methods, self-calibration is a flexible method. Until now, many results have been obtained on the self-calibration based on the circular target. Some researchers introduce the selfcalibration methods based on concentric circles [4–7]. Meng et al. utilize a special temple which contains a circle and a bundle of lines crossing the center of the circle to obtain the intrinsic parameters [8]. Colombo et al. propose a new method based on coaxial circles [9]. Furthermore, some methods first obtain the image of infinity line by parallel lines or orthogonal lines and then get the image of circular points (ICPs) by solving the intersecting points of the infinity line and the ellipses transformed by the concentric circles. Finally, the camera parameters could be obtained by the ICPs [10–12]. Wu et al. propose two calibration methods based on two parallel circles and two coplanar circles [13,14]. The intrinsic parameters could be solved, but the solution for extrinsic parameters was not given. Chen et al. propose a calibration method from analytic geometry to solve both intrinsic parameters and extrinsic parameters [15]. But unfortunately, the

n

Corresponding author. Tel.: þ86 13 801585093; fax: þ86 25 83793000. E-mail addresses: [email protected] (F. Da), [email protected] (Q. Li), [email protected] (H. Zhang), [email protected] (X. Fang). 1 Tel.: þ86 15950505939; fax: þ86 25 83793000. 2 Tel.: þ86 15365182747; fax: þ 86 25 83793000. 3 Tel.: þ86 15950516283; fax: þ86 25 83793000. 0030-3992/$ - see front matter & 2012 Elsevier Ltd. All rights reserved. doi:10.1016/j.optlastec.2012.02.013

equations derived by analytic geometry are too complicated and some intrinsic parameters are needed to be known beforehand. Gurdjos et al. analyze the Euclidean structure of parallel circles by ICPs, which makes a solid foundation for calibration based on parallel circles and coplanar circles [16]. Then Zheng et al. utilize ICPs and quadric enveloping lines to calibrate the camera [17]. This method could avoid the complicated equations derived by analytic geometry in Chen’s method. However, for the projection of circle center the rotated matrix is needed in addition. Due to the perspective projection deviations [18], the true projection of circle center is not the fitting center of the ellipse projected by the circle. The solution for the true projection of circle center is not given. Furthermore, this method also needs some known intrinsic parameters. The basic idea of all the above methods is the invariabilities of ICPs and absolute quadric curve. Due to the highly symmetric trait of the circle, these calibration methods could only obtain the intrinsic parameters. However, the extrinsic parameters are indispensable for 3D reconstruction. All the present methods could get simultaneously the intrinsic parameters and extrinsic parameters only under the condition that one part of intrinsic parameters are known or some other geometric objects (such as lines) are needed. Due to tangent invariance in perspective projection, the common tangents of two circles have been used to obtain camera parameters. In this paper, the property of ICPs and the common tangents of two same circles are utilized to obtain intrinsic parameters and extrinsic parameters of the camera. The only requirement is that there are two same circles in the calibration template or geometric object. In the calibration the relative positions of the two same circles do not need to be known beforehand. Moreover, any prior information about the camera

F. Da et al. / Optics & Laser Technology 44 (2012) 1924–1933

is not needed too. Because two same circles are quite easy to get and it could be extracted easily in practice (such as CD, pop can and bowl cover), therefore the proposed method could make the whole calibration process more flexible and simple.

1925

From simple algebraic knowledge, we have tan y ¼

1 e2iy 1 : i e2iy þ 1

ð5Þ

Substituting Eqs. (3) and (4) into Eq. (5), we obtain



2. Background 2.1. Camera parameters, absolute quadric curve and circular point In projective transformation, a point in 3D world coordinate could be denoted as M¼[X,Y,Z]T and the homogeneous coordinate ~ ¼ ½X,Y,Z,tT . When the point is on the infinite of this point is M plane, t is equal to zero. Otherwise, t is equal to one. After projective projection, the corresponding point in 2D image coordinate is denoted as m¼[u,v]T and the homogeneous coordi~ ¼ ½u,v,1T . The relationship between the point in 3D nate is m world coordinate and its corresponding point in 2D image coordinate is ~ ~ ¼ K½R,TM: lm

ð1Þ

where 2 fx 6 K ¼4 0

fy

7 v0 5

0

0

1

s

u0

3

is the intrinsic parameter matrix, R is the rotation matrix in the extrinsic parameters, T is the translation vector in the extrinsic parameters. In matrix K, s denotes the obliquity factor, fx and fy are the effective focal lengths in x axis and y axis, respectively, (u0,v0) is the principal point. Absolute conic consists of points which are on the infinite plane. The point [X,Y,Z]T on the absolute quadric curve satisfies ~ TM ~ ¼ 0. According to Eq. (1), the image of X2 þY2 þZ2 ¼0, i.e., M absolute quadric curve should satisfy ~ ¼ 0: ~ T K T K 1 m m

ð2Þ

Circular point is the intersection of line at infinity and the circle. There is one pair of circular points in a plane: I(1,i,0,0), J(1,  i,0,0) which satisfies ( x2 þ y2 þz2 ¼ 0 : t¼0 This illustrates that circular point is on the absolute quadric curve. If the ICPs mI, mJ are known, the intrinsic matrix K could be obtained by Eq. (2).

ð6Þ

This is the Laguerre Theorem [19]. According to this theorem and cross-ratio invariance in projective projection, if the ICPs and the images of two common intersecting lines are known, the actual angle between the two lines could be obtained.

3. Calibration procedures The whole calibration procedures are shown in Fig. 1. 3.1. Compute the intrinsic parameters From projective geometry, circular point is the intersecting point of line at infinity and the circle and the ICPs are the intersecting points of the image of line at infinity and image of the circle. For two coplanar circles, every circle could determine one pair of ICPs. Since there is only one pair of circular points on a plane, in the ideal case, ICPs determined by the two circles are the same and the ICPs are on the image of line at infinity. From the above analysis, it can be concluded that ICPs are included in the intersecting points of the images of two circles. Wu et al. have used this idea to obtain the ICPs [17]. Due to the fact that a standard circle is changed to an ellipse under the projective transformation, the ICPs are a pair of conjugate points. Since two separate coplanar circles are used in this paper, the line which passes the circular points should not lie among the two circles. The line is shown as L2 in Fig. 2(a). Therefore, we choose the method for determining circular points as follows. As shown in Fig. 2(a), in the world coordinate system, the intersecting points of two circles are two pairs of conjugate points. Line L1 and line L2 are determined by every pair of the conjugate points. O1 and O2 are the

Extract the two ellipses and obtain the ellipse coefficients Obtain the coordinates of the ICPs Compute the Intrinsic Parameters Determine Inner Tangent Commons and External Tangent Commons

2.2. Isotropic line and Laguerre theorem Isotropic line is a virtual line which passes the circular point. Obviously, the intersection of isotropic line and line at infinity is the circular point. Let two common tangents be l1, l2, the slope of l1, l1 be l1, l2 and the angle between them be y. Let isotropic line be m1,m2 and the slope of m1,m2 be i, þi, respectively. m is the temporary variable. Then we have sin ðl1 ,m1 Þ=sin ðl2 ,m1 Þ m ¼ CRðl1 ,l2 ,m1 ,m2 Þ ¼ sin ðl1 ,m2 Þ=sin ðl1 ,m2 Þ ¼

1 ln m: 2i

Recover the Angle between the Two Inner Tangent Commons Obtain Key Points Determine the Angle Solve the Homography Matrix

1 þi 1lþ2 l1ll12

ðl1 þ iÞðl2 iÞ ð1 þ l1 l2 Þ þ iðl2 l1 Þ ¼ ¼ : ðl2 þ iÞðl1 iÞ ð1 þ l1 l2 Þiðl2 l1 Þ 1i 1lþ2 l1ll12

ð3Þ

Solve the External Camera Parameters and Optimize the Rotated Matrix

According to the straight-line angle formula, we obtain tan y ¼

l2 l1 : 1 þ l1 l2

Non-linear optimization ð4Þ Fig. 1. Flow chart of the algorithm.

1926

F. Da et al. / Optics & Laser Technology 44 (2012) 1924–1933

Zw

Zw

v

Yw

Xw C1

O1

L1

P1

L

u

C

L2 C2

O2

L 1' C 1'

' 2

O

' 2

line_out 1 line_in 1 line_in 2

C1

' 2

C2'

v

Yw

Xw

u C2

O

line_out1' line_in 2' line_in 1' O'

O2

O1

' 1

P

O1'

C1'

line_out 2

O1'

O2'

line_out 2'

Fig. 3. Common tangent of the two circles and their projections. Fig. 2. Schematic diagram for circular points.

center of circle C1 and center of circle C2, respectively. P1 is the intersecting point of line O1 O2 and line L1. From the Fig. 2(a), point P1 lies between point O1 and point O2, while the intersecting point of line O1 O2 and line L2 is beyond that. Fig. 2(b) shows the result of Fig. 2(a) after projective transformation. O01 and O02 are the image point of O1 and O2 after projective transformation, respectively. Here, the fitting centers of ellipses which are projected by the two 0 0 same circles are denoted as O~ 1 and O~ 2 , which are not shown in Fig. 2(b). According to the invariance in projective projection, the line L02 which passes ICPs should not lie between the two ellipses. Therefore, it could be thought that the intersecting point of line L02 and line O01 O02 should lie besides the area which is between point O01 0 0 and point O02 . Since both the fitting centers (O~ 1 , O~ 2 ) and the 0 0 projective centers (O1 , O2 ) should lie in the circles, it can be 0 0 concluded that the intersecting point of line L02 and line O~ 1 , O~ 2 0 0 should lie besides the area which is between point O~ 1 and point O~ 2 . In this paper, two ellipses are extracted in every image. Then the edge points of ellipses are fitted by the least square method and the ellipse coefficients {Aei,Bei,Cei,Dei,Eei}, i¼1,2 (i denotes the serial number of the two ellipses) are obtained. The intersecting point (u, v) could be solved by following equation set: ( u2 þ Ae1 uv þ Be1 v2 þC e1 u þ De1 v þ Ee1 ¼ 0, ð7Þ u2 þ Ae2 uv þ Be2 v2 þC e2 u þ De2 v þ Ee2 ¼ 0: Then the ICPs mI and mJ could be obtained according to the above discussions. If we take three photos, three pairs of ICPs are got. Therefore, the intrinsic parameter matrix K can be obtained from Eq. (2). 3.2. Compute the extrinsic parameters Common tangents are exclusively determined by two circles or ellipses. According to the tangency invariance in the projective projection, common tangents of two ellipses are utilized to obtain the extrinsic parameters. Step 1: Determine inner tangent commons and external tangent commons Let the equation of common line be I1v¼ I2uþI3. I1, I2 and I3 are the coefficients of the line. From the nature of tangency, combine the common line equation and two ellipse equations, then the common line can be solved by following equations: 8 2 2 > < u þ Ae1 uv þ Be1 v þC e1 u þ De1 v þ Ee1 ¼ 0, 2 u þ Ae2 uv þ Be2 v2 þC e2 u þ De2 v þ Ee2 ¼ 0, ð8Þ > : I1 v ¼ I2 u þ I3 : Four groups of solutions are obtained and each group represents each common line. I1, I2 and I3 are the coefficients we want to obtain. They are the slope and the intercept of lines that are obtained from Eq. (8). But the question is how to determine which common line is inner tangent common or external tangent 0 common. We choose the following method. In Fig. 3(b), {line_out1, 0 0 0 line_out2, line_in1,line_in2} are the common tangents of the two ellipses. In the world coordinate system, as shown in Fig. 3(a), the

P4

v

'

C 2'

u

P3 C1'

O2'

'

P2

P1

O1'

P5

'

'

O'

P6

'

'

Pvsp Fig. 4. Schematic diagram for key points.

intersecting point of line O1O2 and external tangent common does not lie between the two circles, while the intersecting point of line O1O2 and inner tangent common lies between the two circles. Therefore, it can be concluded from the invariance in projective projection that if the intersecting point of common line of ellipse and line O01 O02 lies between O01 and O02 , this common line is the inner tangent common; otherwise this common line is the external tangent common. (Here, similar to the discussions mentioned at the end of the second paragraph in this section, 0 0 O01 and O02 can be replaced by O~ 1 and O~ 2 , respectively.) Step 2: Recover the angle between the two inner tangent commons Let the origin of world coordinate system, which is denoted as O shown in Fig. 3(a), be the middle of two circle centers. Due to the two same circles, the intersecting point of two inner tangent commons (line_in1, line_in2) is at O. As shown in Fig. 3(b), from the knowledge of projective geometry, line_in01 and line_in02 are the 0 0 images of two non-isotropic lines and O mI and O mJ are the images of two isotropic lines, where mI, mJ are ICPs. According to Laguerre Theorem and cross-ratio invariance in projective projection, the angle between two inner tangent commons (line_in1, line_in2) is



1 ln ðcrossðline_in01 ,line_in02 ,O0 mI ,O0 mJ ÞÞ, 2i

ð9Þ

where cross(  ) denotes the cross-ratio in projective projection. Step 3: Obtain key points The following steps need some information about O0 , P 01 , P 02 , O01 0 and O02 , which are shown in Fig. 4. Point O is the intersecting point of the two inner tangent commons. P 01 and P 02 are the intersecting points of the line O01 O02 and the ellipse C 01 and C 02 , respectively. These points could be obtained from the definition of intersecting point. Due to the center deviations in projective projection [18], the true projective centers O01 and O02 are not the fitting centers. The method for locating the true projective centers is described as follows. As shown in Fig. 4, P 03 and P05 are the intersecting point of the two external tangent commons and the ellipse C 01 , respectively. P04 and P 06 are the intersecting points of the two external tangent commons and the ellipse C 02 , respectively. From the geometric

F. Da et al. / Optics & Laser Technology 44 (2012) 1924–1933

prosperity of two same circles and the invariance in projective projection, O01 and O02 should be the intersecting points of line O01 O02 and line P 03 P 05 and P 04 P 06 . In the world coordinate system, the original image of line O01 O02 is parallel with the two external tangent commons of two circles. After projective transformation, line O01 O02 intersects the two external tangent commons of two ellipses at an infinite point Pvsp. At the same time, line O01 O02 should pass though O0 due to the two same circles. Therefore, line O01 O02 is equal to line O0 Pvsp. From above-mentioned statement, the point Pvsp which is the intersecting point of the two external tangent commons in two ellipses can be found out and the point 0 O which is the intersecting point of the two inner tangent commons in two ellipses can be found out. Then, the true projective center of ellipse C 01 and C 02 are the interesting point of the line P03 P05 and P 04 P 06 with line O0 Pvsp, respectively. Step 4: Determine the angle The recovered angle in Step 2 by straight-line angle formula is either an acute angle or a right angle. When the world coordinate is recovered, two possible cases will occur. As shown in Fig. 5, for two sets of same circles ({O1, O2} and {O3, O4}), the inner tangency commons are the same as well as the angle between the inner tangency commons. However, the geometric configurations for the two sets of same circles are different. For distinguishing these two cases, the included angle of the two inner tangency commons which does not lie between the two circles is denoted as y, which is shown in Fig. 6. P1 is the intersecting point of line segment O1 O2 and circle C1. P2 is the intersecting point of line segment O1 O2 and circle C2. Let the length of radius in the circle be 1 and the distance between O and P2 be x. Obviously, with the increasing x, the angle pffiffiffi pffiffiffi y increases. When x ¼ 21, y ¼ p=2, therefore, when x 4 21, y 4 p=2. According to the cross-ratio of O1 P2 P1 O2, we have crossðO1 P2 P1 O2 Þ ¼

1 : 4xðx þ 1Þ

1927

From Eq. (10), cross(O1P2P1O2) is an increasing function when x40. Therefore, it can be concluded that

y 4 p=2 if crossðO1 P 2 P 1 O2 Þ 4 0:4268, y r p=2 else:

ð11Þ

From the cross-ratio invariance in the projective projection, we have

y 4 p=2 if crossðO01 P 02 P 01 O02 Þ 4 0:4268, y r p=2 else:

ð12Þ

From the above discussions, the true angle is the difference of

p and the obtained angle in Step 2 when y is larger than p/2. Otherwise, the true angle is the obtained angle in Step 2. Step 5: Solve the homography matrix The homography matrix could be solved by point-to-point or line-to-line. Considering that line is more robust than point, the line-to-line method is selected for the homography matrix to improve the robustness of the algorithm. Let the origin in the world coordinate system, which is denoted as point O in Fig. 3(a), be the intersecting point of the two external tangent commons. According to the angle y obtained above, two external tangency commons in the world coordinate system are line_out 1 ¼ ½ 0 line_out 2 ¼ ½ 0

r

1 T

r

1 T ,

ð13Þ

where r denotes the radius of the circle. Generally, r is taken as the physical measurement unit in the world coordinate system. So two inner tangency commons in the world coordinate system

ð10Þ Table 1 Averages of intrinsic parameters under different noise levels.

O3

NL

fx

fy

s

u0

v0

0 0.4 0.8 1.2 1.6 2.0

1260.000 1256.741 1252.185 1256.232 1249.999 1234.772

1320.0 1316.731 1311.680 1316.717 1309.629 1294.266

0.0 0.926 1.923 2.833 4.673 6.635

600.0 597.854 594.100 598.474 588.669 580.910

500.0 496.909 493.772 497.167 487.064 476.949

C2

C1 O1

O

O2 Table 2.1 Averages of extrinsic parameters under different noise levels.

O4

Fig. 5. Two forms of the same angle.

Fig. 6. Judging for the angle.

NL

a

b

g

0

25.7143

15.0000

 13.8462

0.4

25.6983

14.9479

 13.8134

0.8

25.7064

14.9427

 13.7953

1.2

25.7087

14.9520

 13.7882

1.6

25.6301

14.8285

 13.7291

2.0

25.6164

14.6207

 13.5509

T1 [3.0000 3.0000 10.0000] [3.0140 3.0223 9.9716] [3.0257 3.0219 9. 9595] [3.0118 3.0204 9.9715] [3.0726 3.0919 9.9016] [3.1262 3.1646 9.7744]

1928

F. Da et al. / Optics & Laser Technology 44 (2012) 1924–1933

8

25

7 6 Standard deviation of s

Standard deviation of u0,v0

20

15

10

5 4 3 2

5 1 0

0 0

0.2

0.4

0.6

0.8

1

1.2

1.4

1.6

1.8

2

0

0.2

0.4

0.6

Noise level (pixel)

0.8

1

1.2

1.4

1.6

1.8

2

1.4

1.6

1.8

2

Noise level (pixel)

0.014

25

0.012

Standard deviation of θ

15

10

0.01

0.008

0.006

0.004 5 0.002

0

0 0.2

0.4

0.6

0.8

1

1.2

1.4

1.6

1.8

2

0

0.2

0.4

0.6

Noise level (pixel)

0.8

1

1.2

Noise level (pixel)

Fig. 7. Standard deviations for intrinsic parameters and y under different noise levels.

0.4

0.1

α β γ

0.35

Standard deviation of translation vector

0

Standard deviation of α, β,γ

Standard deviation of fx,fy

20

0.3 0.25 0.2 0.15 0.1 0.05 0

The first element The second element The third element

0.09 0.08 0.07 0.06 0.05 0.04 0.03 0.02 0.01 0

0

0.2

0.4

0.6

0.8

1

1.2

Noise level (pixel)

1.4

1.6

1.8

2

0

0.2

0.4

0.6

0.8

1

1.2

Noise level (pixel)

Fig. 8. Standard deviations for extrinsic parameters and y under different noise levels.

1.4

1.6

1.8

2

F. Da et al. / Optics & Laser Technology 44 (2012) 1924–1933

are line_in1 ¼ ½ tanððpyÞ=2Þ line_in2 ¼ ½ tanððpyÞ=2Þ

1 1

0

T

ð14Þ

0 T ,

In the image coordinate system, let the common tangents of the two ellipses be line_in01 , line_in02 , line_out01 and line_out02 ,which could be obtained by ellipse coefficients. The line transformation in projective projection is

Assume that the model plane is on Z ¼0 of the world coordinate system. From the pin-hole model, we have   H ¼ sK r 1 r 2 t : ð17Þ Let bi, i¼1,2,3 be the column vector of K  1H. According to Eq. (17) and property of rotated matrix, the external camera parameters are 8   detððb1 =:b1 :Þ, 7 ðb2 =:b2 :Þ,ððb1 b2 Þ=:b1 :U:b2 :ÞÞðb1 b2 Þ > , < R ¼ 7 :bb1 : , 7 :bb2 : , 7 :b :U:b : 1

0

l ¼ HT l,

ð15Þ

where H is the homography, l denotes the homogenous coordinate of the practical line in the world coordinate system, l0 denotes the homogenous coordinate of corresponding line in the image coordinate system. According to the four common tangents of the two ellipses and Eq. (15), we can obtain 8 > line_in01 ¼ HT  line_in1 > > > > < line_in0 ¼ HT  line_in2 2 : T 0 > line_out  line_out1 > 1 ¼H > > > : line_out0 ¼ HT  line_out 2 2

2

1

ð18Þ where the sign 7 could be judged by the positions of the camera and the model. It is related with Eq. (15). In Eq. (15) the l0 is assigned as (k0 , 1,d0 ) from k0 x yþd0 ¼0, while l is assigned as (  k,1, d) from y kx  d ¼0. Then the sign¼  1. If l is assigned as (k,  1,d), then the sign¼ 1. Due to the noises in the data, the obtained matrix R does not meet the property of rotated matrix. The optimal rotated matrix R^ of matrix R is the closet matrix to given matrix R in the Frobenius normal number meanings, that is, 2

ð19Þ

R^

T At the same time, R^ satisfies R^ R^ ¼ I. Then, we have

H ¼ ðline_out 01 ,line_out02 ,line_in01 ÞT diagðe1 ,e2 ,e3 Þ T

ðline_out 1 ,line_out2 ,line_in1 Þ ,

ð16Þ

T 2 T ^ ^ ^ ¼ traceððRRÞ ðRRÞÞ ¼ 3 þ traceðRT RÞ2traceðR^ RÞ: :RR: F

T

½ðline_out1 ,line_out2 ,line_in1 Þ

1

line_in2 j

½ðline_out01 ,line_out02 ,line_in01 Þ1 line_in02 j

:

Here, [vector]j is the jth component of vector. Step 6: Solve the external camera parameters and optimize the rotated matrix

Eq. (19) is equivalent to maxR^ ðtraceðR^ RÞÞ. Matrix R could be degenerated as UDVT by the SVD method, where D ¼ diag (s1,s2,s3), s1, s2 and s3 denote the singular values of R, U and ^ then V are orthogonal matrices. Assume Z ¼ V T RU, T T T traceðR^ RÞ ¼ traceðR^ UDV T Þ ¼ traceðV T R^ UDÞ ¼ traceðZDÞ

¼ Table 2.2 Averages of y under different noise levels.

3 X 1

zii si r

3 X

si :

1

Therefore, y

0 96.379

0.4 96.381

0.8 96.375

ð20Þ

According to Eq. (19) and Eq. (20), it could be found that

where

NL

2

b > : t ¼ 7 :b31 : ,

^ min :RR: : F

Then from matrix computation, we have

ej ¼

1929

1.2 96.383

1.6 96.374

2.0 96.373

T

maxðtraceðR^ RÞÞ ¼ traceðDÞ: R^

Zhang’s template

The proposed template Fig. 9. Actual images.

1930

F. Da et al. / Optics & Laser Technology 44 (2012) 1924–1933 T

Assuming V T R^ U ¼ Z ¼ I, we can obtain that R^ ¼ UV T :

ð21Þ

Step 7: Optimum objective function In practice, since the camera parameters are bound together to the projective process of the actual object, the optimization for camera parameters could improve the accuracy of camera calibration. According to the world coordinate system set up above, N points are selected uniformly on the two circles, which are denoted as ðX wi ,Y wi Þ, i ¼ 1, 2    N. Then the image coordinates of these points could be obtained by the world coordinate of the chosen points and the initial camera parameters. The image coordinates are denoted as ðuij ,vij Þ, i ¼ 1, 2    N, j ¼ 1, 2, where j is the serial number of two circles. Since the point (uij,vij) should be on the according ellipse, the following optimum objective function is chosen. N_image X

N X 2 X

k

i¼1j¼1

ð:uij 2 þAej *uij vij þ Bej vij 2 þ C ej uij þ Dej vij þ Eej :Þ2 ,

ð22Þ

where N_image is the number of photos. It should be noted that Levenberg–Marquardt algorithm is used here for non-linear optimization [20]. Since the radius of the circle is the physical measurement unit, the obtained translation vector in extrinsic parameters of the camera has a scale factor difference with the true value.

4. Simulation experiments In the following experiments, the world coordinates of two same circles are set as (xþ3)2 þy2 ¼4 and (x  3)2 þy2 ¼4. The

angle y shown in Fig. 6 is y ¼ p2  a sin ð2=3Þ ¼ 96:37941. The inner camera parameter matrix is 2 3 2 3 1260 0 600 f x s u0 6 7 6 7 1320 500 5: K ¼ 4 0 f y v0 5 ¼ 4 0 0 0 1 0 0 1 The first extrinsic camera parameters are a1 ¼25.71431,

b1 ¼151, g1 ¼ 13.84621, T1 ¼[3,3,10]T. The second extrinsic camera parameters are a2 ¼361, b2 ¼ 16.36371, g2 ¼ 22.51, T2 ¼[2,2,10]T. The third extrinsic camera parameters are a3 ¼ 91,b3 ¼181, g3 ¼ 301, T3 ¼[1,2,10]T. According to the above camera parameters, three photos are generated and every image resolution ratio is 1600 pixeln1600 pixel. Then, Gauss noises with mean 0 and standard deviation s are added to the image points. We varied the noise level from 0 to 2 pixels. For each noise level, we have performed 200 independent tests. The averages of obtained camera parameters are shown in Table 1 and Table 2.1. The standard deviations of the obtained camera parameters are shown in Figs. 7 and 8 (Due to space limitations, the first extrinsic camera parameters are only given in the extrinsic camera parameters.). However the accuracy of y influences the extrinsic camera parameters, and averages and deviations of y are shown in Table 2.2 and Fig. 7(d), respectively, where NL denotes the noise level. The unit of a, b, g and y is the degree and the unit of other variables is pixel. From the experimental results, it can be seen that the proposed method has good accuracy and robustness. At the same time, high accuracy of y ensures the good performances of extrinsic parameters.

5. Real experiment

Table 3 Comparison for intrinsic parameters.

5.1. Data Comparisons with Zhang’s method [21]

Method

fx

fy

s

u0

v0

Zhang’s method Proposed method

2578.0350 2553.2098

2571.4527 2551.3671

 3.0192  3.8318

686.9230 674.8289

488.4222 469.9054

Since Zhang’s method belongs to the traditional calibration method and has high precision, we validate the performance of the proposed method by comparing with Zhang’s method.

Table 4.1 Comparison for extrinsic parameters of the first image. Method

a1

b1

g1

T1

Zhang’s method After standardization Proposed method After standardization

85.1780 4.8220 4.4917 4.4917

 10.9165 10.9165 10.8138 10.8138

191.1374 11.1374 11.3823 11.3823

[  0.7621  0.8017 7.0600]T [  0.1079  0.1136 1.0000]T [0.2776 0.0920 14.0330]T [0.0198 0.0066 1.0000]T

Table 4.2 Comparison for extrinsic parameters of the second image. Method

a2

b2

g2

T2

Zhang’s method After standardization Proposed method After standardization

80.7090 9.2910 8.9793 8.9793

 28.1374 28.1374 28.0467 28.0467

 166.4264 13.5736 13.9816 13.9816

[  1.1401 [  0.1615 [  0.5532 [  0.0394

 0.7788 6.8031]T  0.1103 0.9636]T 0.1293 14.0015]T 0.0092 0.9978]T

Table 4.3 Comparison for extrinsic parameters of the third image. Method

a3

b3

g3

T3

Zhang’s method After standardization Proposed method After standardization

91.3620  1.3620  1.7422  1.7422

19.3543  19.3543  19.8531  19.8531

189.8685 9.8685 10.4399 10.4399

[  0.6677  0.8200 7.7147]T [  0.0946  0.1161 1.0927]T [0.2332 0.0900 14.4992]T [0.0166 0.0064 1.0332]T

F. Da et al. / Optics & Laser Technology 44 (2012) 1924–1933

We first display Zhang’s template and proposed template on LCD, and then rotate the LCD. For each rotation, we use the camera to capture one photo for the two templates. The photos are shown in Fig. 9 and the used camera is UNIQ UP-1800 with image resolution ratio 1030  1380 pixel. Finally, we compute the camera parameters by Zhang’s method and the proposed method, respectively. The intrinsic parameters are shown in Table 3. From the results in Table 3, it can be concluded that the results by Zhang’s method and by the proposed methods are very close. The inclination factor by Zhang’s method and the proposed method are 3.0192 and  3.8318, respectively. The angle between the abscissa and ordinate of the camera by Zhang’s method is f ¼ a tanðf x =absðsÞÞ ¼ 89:93291, while the angle between the abscissa and ordinate of the camera by the proposed method is f ¼ a tanðf x =absðsÞÞ ¼ 89:91401. Both results are very close to 901

Table 5 Comparison for translation vector. Difference of the vectors

Value

T1_zh-T1_pro T2_zh-T2_pro T3_zh-T3_pro

[  0.1277  0.1201 0.0000]T [  0.1221  0.1195  0.0341]T [  0.1112  0.1225 0.0595]T

1931

which is the angle between the abscissa and ordinate of the ideal camera. Since the origin and y-axis in the world coordinate system and the scale factor of the translation vector defined in Zhang’s method are different from that defined in the proposed method, to compare the results with these two methods, the extrinsic parameters in the world coordinate system by Zhang’s method need to change into the world coordinate system by the proposed method, and at the same time the translation vector by the two methods needs to normalize. The special operation methods are following: The y-axis defined in Zhang’s method is opposite to that defined in the proposed method, therefore, the three rotated angle should be reversed accordingly. Since the standard range of a is [0,p/2], the difference between the value of a solved by Zhang’s method and p/2 is taken as the normalized value for a. Similarly, since the standard range of b and g are [0,p], the difference between the value of b and g solved by Zhang’s method and p is taken as the normalized value for b and g, respectively. Furthermore, there is a scale factor between the translation vectors by the two methods. The third weight of each translation vector should be quantized to 1. The results of the extrinsic camera parameters are shown in Table 4.1–4.3. The results show that the rotated angles by the two methods are very close. Furthermore, considering difference of translation vectors by the two methods is the difference of the

100

100

50

50

0

0

-50

-50

-100

-100 40

80 60

20

80 40

40

20 0

40 0

60

20 -20

20

-40

0

8.8

8.85

8.75

8.8

8.7

8.75

8.65

8.7

8.6

8.65

8.55

8.6

8.5

8.55

8.45

-20

0

-40

8.5 0

50

100

150

200

250

300

0

50

100

150

200

250

Fig. 10. Results of 3D reconstruction. (a) Zhang’s method, (b) the proposed method, (c) Zhang’s method and (d) the proposed method.

300

1932

F. Da et al. / Optics & Laser Technology 44 (2012) 1924–1933

origins defined by the two methods, the differences of translation vectors by the two methods for the three photos should be equivalent. At the same time, since the origins defined by the two methods are coplanar, the third components of the three vector differences should be close to zero. So we compute the differences of translation vectors by the two methods for the three photos which are shown in Table 5. From Table 5, the three vector differences are nearly equal. Especially, the third components of the three vector differences are close to zero. All the above experimental data show that the camera parameters obtained by the proposed method are almost the same to the camera parameters by Zhang’s method. This illustrates that the proposed method could get good camera parameters.

0.0746 Unit. It could be concluded that the proposed method has a little advantage over Zhang’s method in planar fitting error and it is close to Zhang’s method in consistency of the neighbor points’ distance. 5.3. Practicability To illustrate the practicability of the proposed method, three common kinds of circular targets (CD, pop can and bowl cover) are taken by the camera from three different views as three photos which are shown in Fig. 11. These photos are then processed by the proposed method and the final calibration results are shown in Table 6.1 and 6.2.

5.2. Reconstruction comparisons with Zhang’s method 6. Conclusion According to the angular points of the Zhang’s template extracted in the two former images shown in Fig. 9, 3D reconstruction results are set up by the camera parameters which are obtained by Zhang’s method and the proposed method. 3D reconstruction results are shown in Fig. 10(a) and (b). Since the physical measurement units are different, the used camera parameters are normalized which is mentioned above. Then, we compute the planar fitting error from the 3D reconstructed information. The error by Zhang’s method is 0.0864 Unit and the error by the proposed method is 0.0789 Unit, where Unit denotes the physical measurement unit. In Fig. 10(a) and (b), there are translations of the coordinates of points. This is because the origins of the world coordinate systems by the two methods are different. Furthermore, we compute the distances between every two neighbor points. The results are shown in Fig. 10(c) and (d). The mean value and standard deviation by Zhang’s method are 8.6250 Unit and 0.0713 Unit. The mean value and standard deviation by the proposed method are 8.6859 Unit and

CD 1

Pop can 1

Bowl cover 1

A novel self-calibration method using two same circles is proposed in this paper. This method is simple and flexible.

Table 6.1 Calibration results of intrinsic parameters.

CD 1 CD 2 CD 3 PC 1 PC 2 PC 3 BC 1 BC 2 BC 3

fx

fy

u0

v0

s

2458.3171

2457.6924

708.0534

514.1853

6.2982

2605.2464

2614.3241

716.2536

527.3701

8.0115

2454.6873

2459.7186

676.2873

536.5353

3.3769

CD 2

CD 3

Pop can 2

Pop can 3

Bowl cover 2

Bowl cover 3

Fig. 11. Three common circular targets.

F. Da et al. / Optics & Laser Technology 44 (2012) 1924–1933

References

Table 6.2 Calibration results of extrinsic parameters.

a

CD 1 CD 2 CD 3 PC 1 PC 2 PC 3 BC 1 BC 2 BC 3

29.7187 14.7846  4.1940  50.8794 60.8270  9.4625 1.2566 17.2931 0.9530

b

 21.1967 0.8661 20.0295 11.4533 27.4083  41.3984 11.4985  16.8883  35.0617

g

16.0020 18.2225 11.2160 1.1908  1.8448 9.4483 15.0244 14.8126 21.5268

1933

T

[  0.2410 0.1297 10.9038]T [0.0166 0.0580 10.1890]T [  0.0661 0.1536 10.2496]T [0.1318  0.2569 16.1973]T [0.0676 0.7827  18.1204]T [  1.1617  1.3035 16.1339]T [0.0456 0.0225 12.4987]T [  0.1916 0.0928 13.1516]T [  0.6697  0.2114 14.6015]T

Compared to the traditional calibration methods such as Zhang’s method, this method does not need any special and fine-made calibration-targets as well as the complicated process for image coordinates matching. Comparing with the other self-calibration methods based on circular targets, this method does not need any other geometric objects and prior information about the camera parameters. Extensive experiments validate the precision and robustness of this method. At the same time, since the target with two same circles is the common geometric object in the nature and especially quite easy for drawing, this method has great values in practice. It should be noted that the edge location of the two ellipses in this paper are pixel-level. Sub-pixel location for the edge points should benefit to obtain the ellipse coefficients and finally improve the calibration results.

Acknowledgments This work was partially funded by National Natural Science Foundation of P.R. China (51175081, 61107001), Natural Science Foundation of JiangSu Province (BK2010058).

[1] Heikkila J. Geometric camera calibration using circular control points. IEEE Transactions on Pattern Analysis and Machine Intelligence 2000;22:1066–77. [2] Asari KV, Kumar S, Radhakrishnan D. A new approach for nonlinear distortion correction in endoscopic images based on least squares estimation. IEEE Transactions on Medical Imaging 1999;18:345–54. [3] Kannala J, Brandt. S. A generic camera model and calibration method for conventional, wide-angle and fish-eye lenses. IEEE Transactions on Pattern Analysis and Machine Intelligence 2006;28:1335–40. [4] Kim JS, Grudjos P, Kweon IS. Geometric and algebraic constraints of projected concentric circles and their application to camera calibration. IEEE Transactions on Pattern Analysis and Machine Intelligence 2005;27:637–42. [5] Kim JS, Kim HW, Kweon IS. A camera calibration method using concentric circles for vision applications. 2002; 23–5. [6] Abad F, Camahort E, Vivo R. Camera calibration using two concentric circles. 2004; 688–96. [7] Fremont V, Chellali R. Direct camera calibration using two concentric circles from a single view. 2002; 93–8. [8] Meng XQ, Hu. ZY. A new easy camera calibration technique based on circular points. Pattern Recognition 2003;36:1155–64. [9] Colombo C, Comanducci D, Bimbo AD. Camera calibration with two arbitrary coaxial circles. 2006; 265–76. [10] Chen YS, Horace HS. Planar metric rectification by algebraically estimating the image of the absolute conic. 2004; 88–91. [11] Wang GH, Wua J, Ji ZQ. Single view based pose estimation from circle or parallel lines. Pattern Recognition Letters 2008;29:977–85. [12] Zhong H, Mai F, Hung YS. Camera calibration using circle and right angles. 2006; 646–9. [13] Wu YH, Zhu H. Camera calibration from the quasi-affine invariance of two parallel circles. 2004; 190–202. [14] Wu YH, Li X, Wu FC. Coplanar circles, quasi-affine invariance and calibration. Image and Vision Computing 2006;24:319–26. [15] Chen Q, Wu H, Wada T. Camera Calibration with two arbitrary coplanar circles. 2004; 521–32. [16] Gurdjos P, Sturm P, Wu YH. Euclidean structure from N 4 ¼ 2 parallel circles: theory and algorithms. 2006; 238–52. [17] Zheng Y, Liu Y. Camera calibration using one perspective view of two arbitrary coplanar circles. Optical Engineering 2008;47:067203. [18] Heikkila J, Silven O. A four-step camera calibration procedure with implicit image correction. 1997; 1106–12. [19] Hartley R, Zisserman A. Multiple view geometry in computer vision. Cambridge University; 2000. [20] More. J. The Levenberg-Marquardt algorithm implementation and theory. 1977; 770636. [21] Zhang ZY. A flexible new technique for camera calibration. IEEE Transactions on Pattern Analysis and Machine Intelligence 2000;22:1330–4.