CN100588269C - Camera array calibration method based on matrix decomposition - Google Patents

Abstract

The invention discloses a camera array calibration method based on matrix decomposition, which comprises the following steps of 1) numbering a camera and taking photos of a calibration plate, and alsonumbering the obtained images according to the taking-photo sequence and the corresponding camera; 2) calculating homographic matrix of the camera by using the actual physic coordinate data of characteristics points on the calibration plate and the image data obtained from taking the photos; 3) constructing and measuring the matrix by using the calculated homographic matrix, the corresponding camera and the number of the calibration plate, and adjusting a proportional factor of the homographic matrix; 4) choosing a corresponding method to conduct matrix decomposition to measurement matrix W with data gap according to the arrangement condition of the camera array minus concentration or sparsity; 5) adding measurement restriction and calculating the internal and the external parameters of the camera; 6) conducting nonlinear optimization to the obtained data. The invention is suitable for camera array arranged in any way and can exactly obtain the internal and the external parameters ofthe camera.

Description

Camera array calibration method based on matrix decomposition

Technical field

The invention belongs to technical field of computer multimedia, particularly stereoscopic vision, 3DTV (three-dimensional television), the camera array calibration technology among the IBR (based on the rendering technique of image).

Background technology

Camera calibration is asked for the process of the imaging geometry model of video camera exactly, and the parameter of this geometrical model is exactly a camera parameters.Inner geometry optical parametric (inner parameter) A comprising video camera _iAnd video camera comprises spin matrix R with respect to the position relation (external parameter) of world coordinate system _iWith translation vector t _iA wherein _iForm is as follows:

$A_i=\left(\begin{array}{ccc}{\mathrm{\α}}_i& c_i& u_{0i}\\ 0& {\mathrm{\β}}_i& v_{0i}\\ 0& 0& 1\end{array}\right)$

For the demarcation of single camera, existing at present a lot of very ripe algorithms comprise: linear approach, nonlinear optimization method, two-step method, biplane method etc.

But along with Development of Multimedia Technology, single camera little by little can not satisfy people's demand.No matter be that the 3DTV (three-dimensional television) that rising at present or the drafting by multiple image obtain the more IBR technology of high quality graphic, all need a plurality of video cameras, even video camera array worked simultaneously.Therefore, how accurate, simply obtain the inside and outside parameter of each video camera in the video camera array, and and then obtain the spatial information of object thus, depth information becomes video camera array key issue in application content.

At present, video camera array there is the comparison extensive studies abroad.MIT once set up a 3-Droom experimental situation in early days, by from different angles the house interior scene being taken simultaneously, realized the three-dimensional reconstruction of inner scene.The video camera array that Stanford university is set up has about more than 100 video cameras, and array is realized different application demands and image processing according to the difference of the dense degree of arranging.In addition, the video camera array of CMU university is a video camera array that can dispose automatically and move, thereby makes whole system can obtain better effect when using the IBR technology to carry out scene rendering.To the application of similar above several so large-scale video camera arrays, if the method that we also use single camera to demarcate one by one will cause: on the one hand, the complexity of algorithm is very high; On the other hand, the precision of system also is difficult to reach requirement.

The demarcation of video camera array, existing method mainly contains: the expansion of (1) Zhang Zhengyou algorithm.At first video camera is demarcated one by one, asked for the confidential reference items of video camera.Utilize video camera and the relation of the mutual alignment between the scaling board that back draws then, obtain correlation between video camera and the video camera by nonlinear optimization; (2) Peter Sturm, Andrei Zaharescu, Henrik

Also successively proposed to demarcate separately under the good situation at video camera Deng the people, the anglec of rotation rank of matrix 3 that utilizes video camera and surface plate to demarcate between the thing retrains, and asks for the anglec of rotation and displacement between video camera and the video camera by matrix decomposition; (3) the homography matrix structure between video camera and the plane scaling board that utilizes that ToshioUeshiba proposes is measured matrix (Measurement Matrix), and this measurement matrix satisfies order 4 constraints simultaneously, thereby will measure matrix decomposition is camera parameters matrix and scaling board parameter matrix, comes video camera array is demarcated.Compare with first method, have better consistency between the camera parameters that method (3) is tried to achieve demarcation.This consistency will make the virtual view reconstruction of back in carrying out depth extraction and 3DTV have higher accuracy.Compare with second method, method (3) no longer requires when asking for the locus attitude of video camera video camera to demarcate restriction accurately in advance, makes calibration process more flexible.

But in the method for Toshio Ueshiba, suppose that the visual field between the video camera is basic the coincidence, when in the calibration process scaling board of diverse location being taken pictures in other words, these scaling boards must be seen simultaneously by all video cameras.Yet, huger when our video camera array, perhaps the video camera array arrangement is more sparse, more special (as: video monitoring system, the 3-D ROOM of MIT) time, we just can not guarantee that all position calibration plates can both be photographed by all video cameras.Even if video camera array is than the situation of comparatively dense, we are difficult to also accomplish that each position calibration plate can both be photographed by all video cameras, and picture quality can reach the requirement of demarcating.

When video camera is taken less than certain scaling board, or the picture quality of this scaling board is can not reach requirement the time, we just can't obtain the homography matrix between this video camera and the scaling board, and this moment, the Toshio method adopted SVD to decompose just no longer feasible to measuring matrix.In order to decompose, need do corresponding adjustment in addition, order 4 constraints when measuring matrix decomposition with assurance to the scale factor of homography matrix to measuring matrix.But when a certain location position plate was not in the visual field at certain video camera, the homography matrix of this disappearance might have influence on the adjustment of other homography matrix scale factor, at this moment measured matrix and just no longer satisfied order 4 constraints, and the decomposition of matrix is also no longer feasible.

Summary of the invention

The objective of the invention is in order to overcome the deficiency of existing calibration technique, it is huger fully to solve video camera array, or the distribution of video camera relatively disperses and situation that the camera coverage that causes differs greatly, and a kind of camera array calibration method based on matrix decomposition is provided

Camera array calibration method based on matrix decomposition of the present invention may further comprise the steps:

1) take pictures with camera number and to scaling board, simultaneously with the image of gained also according to order and the corresponding camera number thereof of taking pictures;

2) view data of utilizing the actual physics coordinate data of characteristic point on the scaling board and taking pictures to obtain is asked for the homography matrix of video camera

3) homography matrix of trying to achieve is measured matrix according to the numbering structure of its corresponding video camera and scaling board

4) according to the intensive or sparse arranging situation of video camera array, adjust the scale factor of homography matrix, and to the measurement matrix of obliterated data

Carry out matrix decomposition;

5) add the tolerance constraint, find the solution the inside and outside parameter of video camera;

6) the gained inside and outside parameter is carried out nonlinear optimization.

Among the present invention, the step of homography matrix that the view data that the said actual physics coordinate data of utilizing characteristic point on the scaling board and taking pictures obtains is asked for video camera is as follows:

A) imaging model of video camera is suc as formula (1)

Wherein, the coordinate of X representation space point in world coordinate system, (u v) projects to coordinate in the camera review coordinate system, A for this spatial point _iThe confidential reference items matrix of representing i video camera, R _i, t _iBe rotation and the displacement of camera coordinate system in world coordinate system,

Refer to that both members differs a scale factor, the transformational relation on the scaling board between characteristic point coordinate and the world coordinate system is suc as formula (2)

More than in two formulas, P _iBe camera parameters matrix, Q ^jLocation parameter matrix for scaling board;

B) the actual coordinate data of characteristic point on view data and the known scaling board are done standardization processing, and ask for the homography matrix of video camera by formula (3) ${\tilde{H}}_i^j={\mathrm{\λ}}_i^j{H}_i^j,$

H wherein _i ^jRepresent the true homography matrix between i video camera and j the scaling board image.

Among the present invention, the said method that the homography matrix of trying to achieve is measured matrix according to the numbering structure of its corresponding video camera and scaling board is, with step 2) in the video camera that obtains and the homography matrix between the scaling board

Arrange by mode shown in the formula (4) according to subscript on it, obtain measuring matrix

Be the label of " rower " expression video camera of matrix, the label of " row mark " expression scaling board image, each All be the matrix of (3 * 3),

$\tilde{W}=\left[\begin{array}{ccc}\widetilde{H_1^1}& ...& \widetilde{H_1^J}\\ .& & .\\ .& & .\\ .& & .\\ \widetilde{H_I^1}& ...& \widetilde{H_I^J}\end{array}\right]---\left(4\right)$

Wherein

Homography matrix that the expression Practical Calculation is come out and real homography matrix H _i ^jDiffer a proportionality factors lambda _i ^jThat is: ${\tilde{H}}_i^j={\mathrm{\λ}}_i^j{H}_i^j.$

Among the present invention, said intensive or sparse arranging situation according to video camera array is adjusted the scale factor of homography matrix, and the measurement matrix of obliterated data is carried out matrix decomposition, and its step is as follows respectively by intensive or sparse arranging situation:

A) when the video camera array dense arrangement

(1) scale factor of adjustment homography matrix

Adjust ${\tilde{H}}_i^j={\mathrm{\λ}}_i^j{H}_i^j$ In proportionality factors lambda _i ^j, method is as follows:

Concerning l and f video camera, they are respectively about the relevant homography matrix N of r and p scaling board image _Lf ^RpForm is as follows:

$N_{\mathrm{lf}}^{\mathrm{rp}}=H_f^r{\left(H_l^r\right)}^{-1}{H}_l^p{\left(H_f^p\right)}^{-1}$

The relevant homography matrix that has scale factor that then calculates For:

${\tilde{N}}_{\mathrm{lf}}^{\mathrm{rp}}={\tilde{H}}_f^r{\left({\tilde{H}}_l^r\right)}^{-1}{\tilde{H}}_l^p{\left({\tilde{H}}_f^p\right)}^{-1}=\left({\mathrm{\λ}}_f^r{\mathrm{\λ}}_l^p\right)/\left({\mathrm{\λ}}_l^r{\mathrm{\λ}}_f^p\right)N_{\mathrm{lf}}^{\mathrm{rp}}---\left(5\right)$

Having a value is u _Ll ^RpDouble characteristic root suc as formula (6)

$u_{\mathrm{ll}}^{\mathrm{rp}}=\left({\mathrm{\λ}}_f^r{\mathrm{\λ}}_l^p\right)/\left({\mathrm{\λ}}_l^r{\mathrm{\λ}}_f^p\right)---\left(6\right)$

Measurement matrix with the step 3) structure

In element press row and column ordering, make the shooting function of its first three rows correspondence photograph all scaling boards, the scaling board of first three columns correspondence can be taken simultaneously by all video cameras, and first video camera in the measurement matrix of fixing this moment and first scaling board be as the reference video camera of the relevant homography matrix of calculating with reference to scaling board, and then first video camera with first scaling board with the relevant homography matrix of j scaling board image of i video camera is: ${\tilde{N}}_{l1}^{\mathrm{<figure-callout\; id="1"\; label="r"\; filenames="C2008101212380029H1.png,C2008101212380029H2.png"\; state="\{\{state\}\}">r</figure-callout>}1}=u_{l1}^{r1}{N}_{l1}^{\mathrm{<figure-callout\; id="1"\; label="r"\; filenames="C2008101212380029H1.png,C2008101212380029H2.png"\; state="\{\{state\}\}">r</figure-callout>}1},$ And u _L1 ^R1By

Double characteristic root calculate, the note u _I1 ^J1Be u _i ^j,

Right

In each element do following adjustment,

Order ${\tilde{H}}_i^j=u_i^j{\tilde{H}}_i^j=u_i^j{\mathrm{\&lambda;}}_i^j{H}_i^j=({\mathrm{\&lambda;}}_i^j{\mathrm{\&lambda;}}_1^j/{\mathrm{\&lambda;}}_1^1)P_i{Q}^j---\left(7\right)$

Then adjusted

For:

$\tilde{W}=\left[\begin{array}{ccc}\widetilde{H_1^1}& ...& \widetilde{H_1^J}\\ .& & .\\ .& & .\\ .& & .\\ \widetilde{H_I^1}& ...& \widetilde{H_I^J}\end{array}\right]$

(2) the measurement matrix of obliterated data

Decomposition

Undertaken by following three steps:

1) be 4 to an order, size is the measurement matrix of 3I * 3J The 3I of given any initial value * 4 matrix P ₀

2) at P _K-1Under the known situation, by formula $B{(j,:)}^T=\left(\tilde{W}{(:,j)}^{T}A\right){\left(A^{T}A\right)}^{-1}(j=1\&CenterDot;\&CenterDot;\&CenterDot;3J)$ Ask for 3J * 4 matrix Q _k, and to Q _kColumn vector do orthogonalization process: Q _k=Q _k* N _k, make

Minimum, the ij representative is right

The factor summation that exists;

3) at Q _kUnder the known situation, by formula $A{(i,:)}^T=\left(\tilde{W}(i,:)B\right){\left(B^{T}B\right)}^{-1}(i=1\&CenterDot;\&CenterDot;\&CenterDot;3I)$ Ask for 3I * 4 matrix P _k, make

Minimum, the ij representative is right

The factor summation that exists;

4) repeat above 2), 3) two steps are up to P _kQ _k ^TConvergence, and, after the contrary ordering of Q, obtain measuring matrix to P

Decomposed form suc as formula (8):

$\tilde{W}=(1/{\mathrm{\&lambda;}}_1^1)\left[\begin{array}{c}{\mathrm{\&lambda;}}_1^1{\mathrm{<figure-callout\; id="1"\; label="P"\; filenames="C2008101212380029H1.png,C2008101212380029H2.png"\; state="\{\{state\}\}">P</figure-callout>}}_1\\ .\\ .\\ .\\ {\mathrm{\&lambda;}}_I^1{P}_I\end{array}\right]\left[\begin{array}{ccc}{\mathrm{\&lambda;}}_1^1{Q}^{1^T}& ...& {\mathrm{\&lambda;}}_1^J{Q}^{J^T}\end{array}\right]---\left(8\right)$

When iteration does not restrain or error

When surpassing threshold value, change iterative initial value P ₀Again iteration gets final product.

B) when the arrangement of video camera array is sparse

(1) the measurement matrix to constructing in the step 3)

Reorder, make

The shooting function of first three rows correspondence photograph the scaling board image of maximum numbers, the scaling board of first three columns correspondence can be photographed simultaneously by the video camera of maximum numbers, and first video camera in the measurement matrix of fixing this moment and first scaling board be as the reference video camera of the relevant homography matrix of calculating with reference to scaling board, and then first video camera with first scaling board with the relevant homography matrix of j scaling board image of i video camera is: ${\tilde{N}}_{l1}^{\mathrm{<figure-callout\; id="1"\; label="r"\; filenames="C2008101212380029H1.png,C2008101212380029H2.png"\; state="\{\{state\}\}">r</figure-callout>}1}=u_{l1}^{r1}{N}_{l1}^{\mathrm{<figure-callout\; id="1"\; label="r"\; filenames="C2008101212380029H1.png,C2008101212380029H2.png"\; state="\{\{state\}\}">r</figure-callout>}1},$ And u _L1 ^R1By

Double characteristic root calculate, the note u _I1 ^J1Be u _i ^j, right

Middle first three rows and first three columns are not empty

With

Those corresponding row and columns

Do following adjustment,

Order ${\tilde{H}}_i^j=u_i^j{\tilde{H}}_i^j=u_i^j{\mathrm{\&lambda;}}_i^j{H}_i^j=({\mathrm{\&lambda;}}_i^1{\mathrm{\&lambda;}}_1^j/{\mathrm{\&lambda;}}_1^1)P_i{Q}^j---\left(7\right)$

Then adjusted

For:

$\tilde{W}=\left[\begin{array}{ccc}\widetilde{H_1^1}& ...& \widetilde{H_1^J}\\ .& & .\\ .& & .\\ .& & .\\ \widetilde{H_I^1}& ...& \widetilde{H_I^J}\end{array}\right]$

(2) the measurement matrix of obliterated data

Decomposition

Undertaken by following step:

1) be 4 to an order, size is the matrix of 3I * 3J

The 3I of given any initial value * 4 matrix P ₀

2) establish each

Scale factor

3) utilize current λ _i ^j, make up and measure matrix

4) at P _K-1Under the known situation, by formula $Q_k{(j,:)}^T=\left(\tilde{W}{(:,j)}^T{P}_{k-1}\right){\left({P_{k-1}}^T{P}_{k-1}\right)}^{-1}(j=1\&CenterDot;\&CenterDot;\&CenterDot;3J)$ Ask for 3J * 4 matrix Q _k, and to Q _kColumn vector do orthogonalization process: Q _k=Q _k* N _k, make

Minimum, the ij representative is right

The factor summation that exists;

5) at Q _kUnder the known situation, by formula $P_k{(i,:)}^T=\left(\tilde{W}(i,:)Q_k\right){\left({Q_k}^T{Q}_k\right)}^{-1}(i=1\&CenterDot;\&CenterDot;\&CenterDot;3I)$ Ask for 3I * 4 matrix P _k, make

Minimum, the ij representative is right

The factor summation that exists;

6) adjust

Middle first three rows and first three are classified as empty

With

Those corresponding row and columns

Scale factor as follows: ${\mathrm{\&lambda;}}_i^j=\mathrm{\&Sigma;}{P}_i(3,:)\&times;Q_j{(3,:)}^T,$ P wherein _i(3 :) represent P _iThe third line, Q _j(3 :) represent Q ^jThe 3rd row;

7) repeating step 3) to 6) up to

Convergence, and to P, after the contrary ordering of Q, the decomposed form that obtains is suc as formula (8):

$\tilde{W}=(1/{\mathrm{\&lambda;}}_1^1)\left[\begin{array}{c}{\mathrm{\&lambda;}}_1^1{\mathrm{<figure-callout\; id="1"\; label="P"\; filenames="C2008101212380029H1.png,C2008101212380029H2.png"\; state="\{\{state\}\}">P</figure-callout>}}_1\\ .\\ .\\ .\\ {\mathrm{\&lambda;}}_I^1{P}_I\end{array}\right]\left[\begin{array}{ccc}{\mathrm{\&lambda;}}_1^1{Q}^{1^T}& ...& {\mathrm{\&lambda;}}_1^J{Q}^{J^T}\end{array}\right]---\left(8\right)$

When iteration does not restrain or error

When surpassing threshold value, change iterative initial value P ₀Again iteration gets final product.

Among the present invention, the constraint of said interpolation tolerance, the step of inside and outside parameter of finding the solution video camera is as follows:

A) to P, Q goes standardization, obtains measuring matrix Decomposed form is as follows:

$\tilde{W}\&cong;\left[\begin{array}{c}{\mathrm{<figure-callout\; id="1"\; label="P"\; filenames="C2008101212380029H1.png,C2008101212380029H2.png"\; state="\{\{state\}\}">P</figure-callout>}}_1^{\&prime;}\\ .\\ .\\ .\\ P_I^{\&prime;}\end{array}\right]\left[\begin{array}{ccc}{\mathrm{<figure-callout\; id="1"\; label="Q"\; filenames="C2008101212380029H1.png,C2008101212380029H2.png"\; state="\{\{state\}\}">Q</figure-callout>}}^{1^{\&prime;}}& ...& Q^{J^{\&prime;}}\end{array}\right]$

B) ask for the video camera confidential reference items

Real camera parameters matrix P and actual P ' relation as (9) formula that obtains of decomposing:

$P_i^{\&prime;}\&cong;P_{i}T;$ $Q^{j^{\&prime;}}\&cong;T^{-1}{Q}^j---\left(9\right)$

Wherein T is a nonsingular matrix, can determine as follows:

Concerning first video camera, R ₁Be unit matrix, t ₁Each component be 0, then ${\mathrm{<figure-callout\; id="1"\; label="P"\; filenames="C2008101212380029H1.png,C2008101212380029H2.png"\; state="\{\{state\}\}">P</figure-callout>}}_1=A_1{R_1}^T[I-t_1]=\left[\begin{array}{cc}{\mathrm{<figure-callout\; id="1"\; label="K"\; filenames="C2008101212380029H1.png,C2008101212380029H2.png"\; state="\{\{state\}\}">K</figure-callout>}}_1& 0\end{array}\right],$ With the matrix P ' after decomposing season ₁Has form P ' ₁=[I 0] is so can get $T=\left[\begin{array}{cc}{K}_1^{-1}& 0\\ h^T& h\end{array}\right],$ H wherein ^TBe one 3 * 1 vector, h is a scale factor,

If the scale factor in the formula (9) is β ^j, then formula (9) can be written as: ${\mathrm{TQ}}^{j^{\&prime;}}={\mathrm{\&beta;}}^j{Q}^j,$

Q to each scaling board correspondence ^j, by the p in the formula (2), the orthogonality of q draws formula (10):

$p^{j^{\&prime;}T}{\mathrm{\&omega;}}_1{p}^{j^{\&prime;}}=q^{j^{\&prime;}T}{\mathrm{\&omega;}}_1{q}^{j^{\&prime;}}={\left({\mathrm{\&beta;}}^j\right)}^2,p^{j^{\&prime;}T}{\mathrm{\&omega;}}_1{q}^{j^{\&prime;}}=0---\left(10\right)$

Wherein ${\mathrm{\&omega;}}_1=K_1^{-T}{K}_1^{-1},$ Obtain the quadrature equation of first all scaling boards of photographing of shooting function by formula (10), and all equations that superpose, obtain K ₁Thereby:

${\mathrm{\&beta;}}^j=\sqrt{(p^{j^{\&prime;}T}{\mathrm{\&omega;}}_1{p}^{j^{\&prime;}}+q^{j^{\&prime;}T}{\mathrm{\&omega;}}_1{q}^{j^{\&prime;}})/2}---\left(11\right)$

Again, each scaling board all satisfies equation (12):

$\left[\begin{array}{cc}{h}^T& h\end{array}\right]Q^{j^{\&prime;}}=\left[\begin{array}{ccc}0& 0& {\mathrm{\&beta;}}^j\end{array}\right]---\left(12\right)$

The above-mentioned equation (12) that superposes, [h ^TH], and then obtain transformation matrix T, so first camera parameters matrix P is arranged ₁=P ' ₁T ^-1, the parameter of other video camera also can be by P _i=P ' _iT ^-1Obtain.

C) ask for the outer ginseng of video camera

External parameters of cameras is tried to achieve as follows:

$P_i=A_i{R}_i^T[I-t_i]=\left[M\right|b]---\left(13\right)$

So: $R_i^{\&prime;}=M^T{A}_i^{-T},$ And to R ' _iBe and try to achieve spin matrix R after SVD decomposes _iSuc as formula (14),

$R_i^{\&prime;}={\mathrm{USV}}^T\&DoubleRightArrow;R_i={\mathrm{UV}}^T---\left(14\right)$

Displacement between video camera: t _i=-M ^-1 _iB (15).

Among the present invention, said step of the gained inside and outside parameter being carried out nonlinear optimization:

Carry out nonlinear optimization by formula (16):

$\min \underset{i=1}{\overset{I}{\mathrm{\&Sigma;}}}\underset{j=1,i}{\overset{J}{\mathrm{\&Sigma;}}}\underset{n}{\overset{\mathrm{Nj}}{\mathrm{\&Sigma;}}}{\left|\right|m_{\mathrm{in}}^j-\hat{m}(K_i,R_i,t_i,Q^j,M_n^j)\left|\right|}^2---\left(16\right)$

In the formula,

Be with a certain characteristic point M on the scaling board _n ^jAgain project to coordinate in the image coordinate system, m by the camera model of trying to achieve _In ^jBe M _n ^jPosition coordinates in actual photographed image.

The present invention is applicable to the video camera array of arranging in any way.When the quality that has well solved or some image very big when the camera coverage difference can not reach requiring of camera calibration, correctly be decomposed into the implementation method of camera parameters matrix and scaling board parameter matrix with measuring matrix, and then obtain the inside and outside parameter of video camera accurately.In addition, by the observation to the matrix decomposition form, we have also proposed effectively to adjust by iteration the scale factor of video camera homography matrix, guarantee to measure rank of matrix 4 constraints, and optimize calibration result.Compare with present existing other camera array calibration method, the consistency that the present invention has between the camera parameters of trying to achieve is better, the precision height, and noise robustness is strong, realizes simply, easily, plurality of advantages such as is widely used.

Description of drawings

Fig. 1 is the schematic diagram with the video camera array of compact mode arrangement;

Fig. 2 is the schematic diagram of the video camera array arranged in sparse mode;

Fig. 3 is camera parameters matrix and scaling board parameter matrix schematic diagram;

Fig. 4 is the focal length of camera error of calculation and noise factor graph of a relation;

Fig. 5 is the video camera central point u0 error of calculation and noise factor graph of a relation;

Fig. 6 is the video camera central point v0 error of calculation and noise factor graph of a relation;

Fig. 7 is the video camera central point u0 error of calculation and noise factor graph of a relation;

Fig. 8 is the video camera central point v0 error of calculation and noise factor graph of a relation;

Embodiment

Further specify the present invention below in conjunction with embodiment.

Camera array calibration method based on matrix decomposition may further comprise the steps:

1. hypothesis has six video cameras that scaling board is taken pictures (as shown in Figure 1) at 24 diverse locations, simultaneously with the image of gained also according to order and the corresponding camera number thereof of taking pictures;

We will represent camera number with subscript i in the expression of back, and subscript j represents the numbering of scaling board image.Each video camera photographs the image of 4 width of cloth scaling boards at least simultaneously.

2. utilize the actual physics coordinate data of characteristic point on the scaling board and the view data that obtains of taking pictures is asked for the homography matrix of video camera

A) imaging model of video camera is suc as formula (1)

Wherein, the coordinate of X representation space point in world coordinate system, (u v) projects to coordinate in the camera review coordinate system, A for this spatial point _iThe confidential reference items matrix of representing i video camera, R _i, t _iBe rotation and the displacement of camera coordinate system in world coordinate system, Refer to that both members differs a scale factor.Transformational relation on the scaling board between characteristic point and the world coordinate system is suc as formula (2)

P in the formula (2) ^j, q ^jBe respectively the unit vector of scaling board coordinate system, d ^jBe the distance between scaling board and the world coordinate system, more than in two formulas, P _iBe camera parameters matrix, Q ^jLocation parameter matrix for scaling board; Its concrete physical significance can be with reference to Fig. 2, and among the figure, x, y are the scaling board coordinate system, x _i, y _i, z _iBe the camera coordinate system of video camera i, x ₀, y ₀, z ₀Overlap for the camera coordinate system of first video camera and with world coordinate system.

B) scaling board of our use has 12 row, 12 row, totally 144 characteristic points, to this 144 characteristic point in the imager coordinate data on actual coordinate on the scaling board (x y) and the video camera thereof (after u v) does standardization processing, (normalized step please refer to document: " In defense ofthe 8-point algorithm " .R.Hartley.In Proc.5th International Conference on Computer Vision, pages 1064-1070, Boston, MA, June 995; ) ask for the homography matrix of video camera by formula (3) ${\tilde{H}}_i^j={\mathrm{\&lambda;}}_i^j{H}_i^j,$

H wherein _i ^jRepresent the true homography matrix between i video camera and j the scaling board image.

3. the homography matrix of trying to achieve is measured matrix according to the numbering structure of its corresponding video camera and scaling board

With step 2) in the video camera that obtains and the homography matrix between the scaling board

Press mode shown in the formula according to subscript on it and arrange, obtain measuring matrix

Be the label of " rower " expression video camera of matrix, the label of " row mark " expression scaling board.Attention: each

It all is the matrix of (3 * 3).

$\tilde{W}=\left[\begin{array}{ccc}\stackrel{~}{{\mathrm{<figure-callout\; id="1"\; label="H"\; filenames="C2008101212380029H1.png,C2008101212380029H2.png"\; state="\{\{state\}\}">H</figure-callout>}}_1^1}& ...& \widetilde{H_1^J}\\ .& & .\\ .& & .\\ .& & .\\ \widetilde{H_I^1}& ...& \widetilde{H_I^J}\end{array}\right]---\left(4\right)$

Wherein

Homography matrix that the expression Practical Calculation is come out and real homography matrix H _i ^jDiffer a proportionality factors lambda _i ^jThat is: ${\tilde{H}}_i^j={\mathrm{\&lambda;}}_i^j{H}_i^j.$

Will be because of scaling board not in the camera coverage scope and imponderable homography matrix

Treat by empty matrix, in Fig. 1, video camera 1 be can't see the scaling board of position at scaling board 2 places, then

Be empty matrix, do not fill out any data.

4. according to the intensive or sparse arranging situation of video camera array, adjust the scale factor of homography matrix, and the measurement matrix of obliterated data is carried out matrix decomposition;

According to H in the formula (3) and P, the relation of Q, measure matrix W really and can decompose as follows:

Promptly to the measurement matrix matrix W of a 3I * 3J (I represents the number of video camera, and J represents the number of the scaling board image taken), can be decomposed into size and be 3J * 4 scaling board parameter matrix Q for 3I * 4 camera parameters matrix P and size.

A) when the video camera array dense arrangement

Here intensive is meant to have at least a shooting function to photograph all position calibration plates, and at least on a position this scaling board can be photographed simultaneously by all video cameras simultaneously.This kind situation relatively is applicable to the video camera array that plane parallel is arranged, shown in figure one.The video camera array of the video camera array of Stanford university and CMU university all belongs to this type.

(1) scale factor of adjustment homography matrix

Video camera number I is 6 in this example, and the image number J of the scaling board of shooting is 24, as shown in Figure 1.

Adjust ${\tilde{H}}_i^j={\mathrm{\&lambda;}}_i^j{H}_i^j$ In proportionality factors lambda _i ^j, method is as follows:

The relevant homography matrix of definition: concerning l and f video camera, they are about the relevant homography matrix N of r and p scaling board image _Lf ^RpForm is as follows:

$N_{\mathrm{lf}}^{\mathrm{rp}}=H_f^r{\left(H_l^r\right)}^{-1}{H}_l^p{\left(H_f^p\right)}^{-1}$

The relevant homography matrix that has scale factor that then calculates

For:

${\tilde{N}}_{\mathrm{lf}}^{\mathrm{rp}}={\tilde{H}}_f^r{\left({\tilde{H}}_l^r\right)}^{-1}{\tilde{H}}_l^p{\left({\tilde{H}}_f^p\right)}^{-1}=\left({\mathrm{\&lambda;}}_f^r{\mathrm{\&lambda;}}_l^p\right)/\left({\mathrm{\&lambda;}}_l^r{\mathrm{\&lambda;}}_f^p\right)N_{\mathrm{lf}}^{\mathrm{rp}}---\left(5\right)$

Because N _Lf ^RpHas a unit double root, therefore

Having a value is u _Lf ^RpDouble characteristic root suc as formula (6)

$u_{\mathrm{ll}}^{\mathrm{rp}}=\left({\mathrm{\&lambda;}}_f^r{\mathrm{\&lambda;}}_l^p\right)/\left({\mathrm{\&lambda;}}_l^r{\mathrm{\&lambda;}}_f^p\right)---\left(6\right)$

Measurement matrix with the step 3) structure

In element press row and column ordering, make the shooting function of its first three rows correspondence photograph all scaling boards, the scaling board of first three columns correspondence can be taken simultaneously by all video cameras, and first video camera in the measurement matrix of fixing this moment and first scaling board be as the reference video camera of the relevant homography matrix of calculating with reference to scaling board, and then first video camera with first scaling board with the relevant homography matrix of j scaling board image of i video camera is: ${\tilde{N}}_{i1}^{\mathrm{<figure-callout\; id="1"\; label="j"\; filenames="C2008101212380029H1.png,C2008101212380029H2.png"\; state="\{\{state\}\}">j</figure-callout>}1}=u_{i1}^{j1}{N}_{i1}^{\mathrm{<figure-callout\; id="1"\; label="j"\; filenames="C2008101212380029H1.png,C2008101212380029H2.png"\; state="\{\{state\}\}">j</figure-callout>}1},$ And u _I1 ^J1Can pass through

Double characteristic root calculate, the note u _I1 ^J1Be u _i ^j,

Right

In each element do following adjustment,

Order ${\tilde{H}}_i^j=u_i^j{\tilde{H}}_i^j=u_i^j{\mathrm{\&lambda;}}_i^j{H}_i^j=({\mathrm{\&lambda;}}_i^j{\mathrm{\&lambda;}}_1^j/{\mathrm{\&lambda;}}_1^1)P_i{Q}^j---\left(7\right)$

Then adjusted

For:

$\tilde{W}=\left[\begin{array}{ccc}\stackrel{~}{{\mathrm{<figure-callout\; id="1"\; label="H"\; filenames="C2008101212380029H1.png,C2008101212380029H2.png"\; state="\{\{state\}\}">H</figure-callout>}}_1^1}& ...& \widetilde{H_1^J}\\ .& & .\\ .& & .\\ .& & .\\ \widetilde{H_I^1}& ...& \widetilde{H_I^J}\end{array}\right]$

Under the dense arrangement mode, for

In disappearance element (as: H _m ⁿ), do not adjust.

(2) to the measurement matrix of obliterated data

Carry out matrix decomposition.

Undertaken by following three steps:

1) be 4 to an order, size is the measurement matrix of 3I * 3J

The 3I of given any initial value * 4 matrix P ₀

2) at P _K-1Under the known situation, by formula $B{(j,:)}^T=\left(\tilde{W}{(:,j)}^{T}A\right){\left(A^{T}A\right)}^{-1}(j=1\&CenterDot;\&CenterDot;\&CenterDot;3J)$ Ask for 3J * 4 matrix Q _k, and to Q _kColumn vector do orthogonalization process: Q _k=Q _k* N _k, make

Minimum, the ij representative is right

The factor summation that exists;

3) at Q _kUnder the known situation, by formula $A{(i,:)}^T=\left(\tilde{W}(i,:)B\right){\left(B^{T}B\right)}^{-1}(i=1\&CenterDot;\&CenterDot;\&CenterDot;3I)$ Ask for 3I * 4 matrix P _k, make

Minimum, the ij representative is right

The factor summation that exists;

4) repeat above 2), 3) two steps are up to P _kQ _k ^TConvergence, and, after the contrary ordering of Q, obtain measuring matrix to P

Decomposed form suc as formula (8):

$\tilde{W}=(1/{\mathrm{\&lambda;}}_1^1)\left[\begin{array}{c}{\mathrm{\&lambda;}}_1^1{\mathrm{<figure-callout\; id="1"\; label="P"\; filenames="C2008101212380029H1.png,C2008101212380029H2.png"\; state="\{\{state\}\}">P</figure-callout>}}_1\\ .\\ .\\ .\\ {\mathrm{\&lambda;}}_I^1{P}_I\end{array}\right]\left[\begin{array}{ccc}{\mathrm{\&lambda;}}_1^1{Q}^{1^T}& ...& {\mathrm{\&lambda;}}_1^J{Q}^{J^T}\end{array}\right]---\left(8\right)$

When iteration does not restrain or error

When surpassing threshold value, change iterative initial value P ₀Again iteration gets final product.

B) when the arrangement of video camera array is sparse

Here sparse refers to when the arrangement of video camera more scattered, and the perhaps great disparity relatively of the anglec of rotation between the video camera causes the very big situation of visual field difference of video camera.Mainly use as video monitoring system the 3D ROOM of MIT etc.As shown in Figure 3,8 video cameras form a circle, and at this moment, no matter where scaling board is placed in, all will exist certain several video camera can not photograph the front of this scaling board, then corresponding homography matrix between these video cameras and the scaling board Just can not obtain.

(1) scale factor of adjustment homography matrix

Video camera number I is 8 in this example, and the image number J of the scaling board of shooting is 60, as shown in Figure 3, in the process of taking pictures, should make as far as possible that the number of the scaling board front plate that each shooting function photographs is more even.

To the measurement matrix of constructing in the step 3) Reorder, make

The shooting function of first three rows correspondence photograph the scaling board image of maximum numbers, the scaling board of first three columns correspondence can be photographed simultaneously by the video camera of maximum numbers, and first video camera in the measurement matrix of fixing this moment and first scaling board be as the reference video camera of the relevant homography matrix of calculating with reference to scaling board, and then first video camera with first scaling board with the relevant homography matrix of j scaling board image of i video camera is: ${\tilde{N}}_{l1}^{\mathrm{<figure-callout\; id="1"\; label="r"\; filenames="C2008101212380029H1.png,C2008101212380029H2.png"\; state="\{\{state\}\}">r</figure-callout>}1}=u_{l1}^{r1}{N}_{l1}^{\mathrm{<figure-callout\; id="1"\; label="r"\; filenames="C2008101212380029H1.png,C2008101212380029H2.png"\; state="\{\{state\}\}">r</figure-callout>}1},$ And u _L1 ^R1By

Double characteristic root calculate, the note u _I1 ^J1Be u _i ^j,

Right

Middle first three rows and first three columns are not empty With Those corresponding row and columns Do following adjustment,

Order ${\tilde{H}}_i^j=u_i^j{\tilde{H}}_i^j=u_i^j{\mathrm{\&lambda;}}_i^j{H}_i^j=({\mathrm{\&lambda;}}_i^1{\mathrm{\&lambda;}}_1^j/{\mathrm{\&lambda;}}_1^1)P_i{Q}^j---\left(7\right)$

Because of the restriction of sparse arrangement, to what can't adjust according to formula (7)

Scale factor, we utilize the measurement matrix

Pro forma characteristics obtain in the matrix decomposition iterative process.

(2) the measurement matrix of obliterated data

Decomposition

Undertaken by following step:

1) be 4 to an order, size is the matrix of 3I * 3J

The 3I of given any initial value * 4 matrix P ₀

2) establish each

Scale factor ${\mathrm{\&lambda;}}_i^j=1;$

3) utilize current λ _i ^j, make up and measure matrix

4) at P _K-1Under the known situation, by formula $Q_k{(j,:)}^T=\left(\tilde{W}{(:,j)}^T{P}_{k-1}\right){\left({P_{k-1}}^T{P}_{k-1}\right)}^{-1}(j=1\&CenterDot;\&CenterDot;\&CenterDot;3J)$ Ask for 3J * 4 matrix Q _k, and to Q _kColumn vector do orthogonalization process: Q _k=Q _k* N _k, make Minimum, the ij representative is right

The factor summation that exists;

5) at Q _kUnder the known situation, by formula $P_k{(i,:)}^T=\left(\tilde{W}(i,:)Q_k\right){\left({Q_k}^T{Q}_k\right)}^{-1}(i=1\&CenterDot;\&CenterDot;\&CenterDot;3I)$ Ask for 3I * 4 matrix P _k, make Minimum, the ij representative is right The factor summation that exists;

6) adjust

Middle first three rows and first three are classified as empty

With

Those corresponding row and columns

Proportionality factors lambda _i ^jAs follows: ${\mathrm{\&lambda;}}_i^j=\mathrm{\&Sigma;}{P}_i(3,:)\&times;Q_j{(3,:)}^T,$ P wherein _i(3 :) represent P _iThe third line, Q _j(3 :) represent Q ^jThe 3rd row;

7) repeating step 3) to 6) up to

Convergence, and to P, after the contrary ordering of Q, the decomposed form that obtains is suc as formula (8):

$\tilde{W}=(1/{\mathrm{\&lambda;}}_1^1)\left[\begin{array}{c}{\mathrm{\&lambda;}}_1^1{\mathrm{<figure-callout\; id="1"\; label="P"\; filenames="C2008101212380029H1.png,C2008101212380029H2.png"\; state="\{\{state\}\}">P</figure-callout>}}_1\\ .\\ .\\ .\\ {\mathrm{\&lambda;}}_I^1{P}_I\end{array}\right]\left[\begin{array}{ccc}{\mathrm{\&lambda;}}_1^1{Q}^{1^T}& ...& {\mathrm{\&lambda;}}_1^J{Q}^{J^T}\end{array}\right]---\left(8\right)$

When iteration does not restrain or error

When surpassing threshold value, change iterative initial value P ₀Again iteration gets final product.

5. add the tolerance constraint, find the solution the inside and outside parameter of video camera

A) to P, Q goes standardization, obtains measuring matrix

Decomposed form is as follows:

$\tilde{W}\&cong;\left[\begin{array}{c}{\mathrm{<figure-callout\; id="1"\; label="P"\; filenames="C2008101212380029H1.png,C2008101212380029H2.png"\; state="\{\{state\}\}">P</figure-callout>}}_1^{\&prime;}\\ .\\ .\\ .\\ P_I^{\&prime;}\end{array}\right]\left[\begin{array}{ccc}{\mathrm{<figure-callout\; id="1"\; label="Q"\; filenames="C2008101212380029H1.png,C2008101212380029H2.png"\; state="\{\{state\}\}">Q</figure-callout>}}^{1^{\&prime;}}& ...& Q^{J^{\&prime;}}\end{array}\right]$

B) ask for the video camera confidential reference items

Real camera parameters matrix P and actual P ' relation as (9) formula that obtains of decomposing:

$P_i^{\&prime;}\&cong;P_{i}T;$ $Q^{j^{\&prime;}}\&cong;T^{-1}{Q}^j---\left(9\right)$

Wherein T is a nonsingular matrix, can determine as follows:

Concerning first video camera, R ₁Be unit matrix, t ₁Each component be 0, then ${\mathrm{<figure-callout\; id="1"\; label="P"\; filenames="C2008101212380029H1.png,C2008101212380029H2.png"\; state="\{\{state\}\}">P</figure-callout>}}_1=A_1{R_1}^T[I-t_1]=\left[\begin{array}{cc}{\mathrm{<figure-callout\; id="1"\; label="K"\; filenames="C2008101212380029H1.png,C2008101212380029H2.png"\; state="\{\{state\}\}">K</figure-callout>}}_1& 0\end{array}\right],$ With the matrix P ' after decomposing season ₁Has form P ' ₁=[I 0] is so can get $T=\left[\begin{array}{cc}{K}_1^{-1}& 0\\ h^T& h\end{array}\right],$ H wherein ^TBe one 3 * 1 vector, h is a scale factor.

If the scale factor in the formula (9) is β ^j, then formula (9) can be written as: ${\mathrm{TQ}}^{j^{\&prime;}}={\mathrm{\&beta;}}^j{Q}^j.$

Q to each scaling board correspondence ^j, by the p in the formula (2), the orthogonality of q draws formula (10):

$p^{j^{\&prime;}T}{\mathrm{\&omega;}}_1{p}^{j^{\&prime;}}=q^{j^{\&prime;}T}{\mathrm{\&omega;}}_1{q}^{j^{\&prime;}}={\left({\mathrm{\&beta;}}^j\right)}^2,p^{j^{\&prime;}T}{\mathrm{\&omega;}}_1{q}^{j^{\&prime;}}=0---\left(10\right)$

Wherein ${\mathrm{\&omega;}}_1=K_1^{-T}{K}_1^{-1},$ Obtain the quadrature equation of first all scaling boards of photographing of shooting function by formula (10), and all equations that superpose, obtain K ₁Thereby:

${\mathrm{\&beta;}}^j=\sqrt{(p^{j^{\&prime;}T}{\mathrm{\&omega;}}_1{p}^{j^{\&prime;}}+q^{j^{\&prime;}T}{\mathrm{\&omega;}}_1{q}^{j^{\&prime;}})/2}---\left(11\right)$

Again, each scaling board all satisfies equation (12):

$\left[\begin{array}{cc}{h}^T& h\end{array}\right]Q^{j^{\&prime;}}=\left[\begin{array}{ccc}0& 0& {\mathrm{\&beta;}}^j\end{array}\right]---\left(12\right)$

The above-mentioned equation (12) that superposes, [h ^TH], and then receive transformation matrix T, so first camera parameters matrix P is arranged ₁=P ' ₁T ^-1, the parameter of other video camera also can be by P _i=P ' _iT ^-1Obtain.

C) ask for the outer ginseng of video camera

The outer ginseng of video camera is asked for as follows:

$P_i=A_i{R}_i^T[I-t_i]=\left[M\right|b]---\left(13\right)$

So have: $R_i^{\&prime;}=M^T{A}_i^{-T},$ And to R ' _iBe and try to achieve spin matrix R after SVD decomposes _iSuc as formula (14),

$R_i^{\&prime;}={\mathrm{USV}}^T\&DoubleRightArrow;R_i={\mathrm{UV}}^T---\left(14\right)$

Displacement between video camera: t _i=-M ^-1 _iB.(15)

6. nonlinear optimization

For further improving the accuracy of the camera interior and exterior parameter that calculates, these parameters are carried out nonlinear optimization by formula (16):

$\min \underset{i=1}{\overset{I}{\mathrm{\&Sigma;}}}\underset{j=1,i}{\overset{J}{\mathrm{\&Sigma;}}}\underset{n}{\overset{\mathrm{Nj}}{\mathrm{\&Sigma;}}}{\left|\right|m_{\mathrm{in}}^j-\hat{m}(K_i,R_i,t_i,Q^j,M_n^j)\left|\right|}^2---\left(16\right)$

In the formula,

Be with a certain characteristic point M on the scaling board _n ^jAgain project to coordinate in the image coordinate system, m by the camera model of trying to achieve _In ^jBe M _n ^jPosition coordinates in actual photographed image.

Fig. 4 is that method of the present invention is carried out error performance figure after the emulation to Fig. 8.Employed video camera array comprises 6 video cameras in the emulation, and wherein, the inside and outside parameter of video camera is provided with as follows in each experiment:

1) external parameters of cameras setting: video camera and video camera at the interval of x direction at (300mm, evenly distribute 600mm), y, the interval of z direction is at (5mm, evenly distribute 5mm), each video camera is about world coordinate system (x, y, z) anglec of rotation of reference axis also evenly distributes between (0 °, 18 °).

2) camera intrinsic parameter setting: it is 790 that each focus of camera meets average, and variance is 20 normal distribution, and the center point coordinate of video camera is (300,300).

In the emulation, scaling board moves freely before video camera array 500mm, and it is about (z) anglec of rotation of each reference axis evenly distributes between (0 °, 60 °) for x, y.Simultaneously, in order to test the signal to noise ratio of the inventive method, in the experiment, giving each picture point stack average is 0, be the Gaussian noise that step-length is gone forward one by one with 0.2 the variance from 0 to 2, and each noise factor done to test for 100 times is averaged the result.

Fig. 4 is the relation curve of focus of camera sum of errors noise factor, only listed first respectively among the figure, three, five video cameras are mean error and the average focus error of accumulative total of all video cameras and the relation curve between the noise factor of focal length separately, wherein the x axle is the noise superimposed coefficient value, and the y axle is the focal length of camera error amount.

Fig. 5 is the relation curve of the sum of errors noise factor of video camera central point abscissa u0, only listed first respectively among the figure, three, five video cameras separately the average u0 error of accumulative total of the mean error of abscissa u0 and all video cameras with the relation curve of the noise factor in the emulation, wherein the x axle is the noise superimposed coefficient value, and the y axle is a video camera abscissa u0 error amount.

Fig. 6 is the relation curve of the sum of errors noise factor of video camera central point ordinate v0, only listed first respectively among the figure, three, five video cameras separately the average v0 error of accumulative total of the mean error of ordinate v0 and all video cameras with the relation curve of the noise factor in the emulation, wherein the x axle is the noise superimposed coefficient value, and the y axle is a video camera ordinate v0 error amount.

Fig. 7 is the relation curve that the mutual distance and position of video camera concerns the sum of errors noise factor, only listed second respectively among the figure, three, five video cameras separately with the mutual distance and position error of first video camera and all video cameras with the accumulative total mean error of the mutual distance and position error of first video camera and the relation curve of the noise factor in the emulation, wherein the x axle is the noise superimposed coefficient value, and the y axle is the mutual range error value of video camera.

Fig. 8 is the relation curve of the mutual anglec of rotation sum of errors of video camera noise factor, only listed second respectively among the figure, three, five video cameras separately with the anglec of rotation error of first video camera and all video cameras with the accumulative total mean error of the anglec of rotation error of first video camera and the relation curve of the noise factor in the emulation, wherein the x axle is the noise superimposed coefficient value, and the y axle is the mutual anglec of rotation error amount of video camera.

By Fig. 4 to error performance curve shown in Figure 8 as seen, utilize each camera interior and exterior parameter accuracy of the inventive method gained higher, noise resisting ability is strong, can satisfy the requirement that video camera array is demarcated.

Claims (6)

1. based on the camera array calibration method of matrix decomposition, it is characterized in that may further comprise the steps:

1) take pictures with camera number and to scaling board, simultaneously with the image of gained also according to order and the corresponding camera number thereof of taking pictures;

2) view data of utilizing the actual physics coordinate data of characteristic point on the scaling board and taking pictures to obtain is asked for the homography matrix of video camera

3) homography matrix of trying to achieve is measured matrix according to the numbering structure of its corresponding video camera and scaling board

4) according to the intensive or sparse arranging situation of video camera array, adjust the scale factor of homography matrix, and to the measurement matrix of obliterated data Carry out matrix decomposition;

5) add the tolerance constraint, find the solution the inside and outside parameter of video camera;

6) the gained inside and outside parameter is carried out nonlinear optimization.

2. the camera array calibration method based on matrix decomposition according to claim 1, it is characterized in that the said actual physics coordinate data of utilizing characteristic point on the scaling board and the view data that obtains of taking pictures to ask for the step of homography matrix of video camera as follows:

A) imaging model of video camera is suc as formula (1)

Wherein, the coordinate of X representation space point in world coordinate system, (u v) projects to coordinate in the camera review coordinate system, A for this spatial point _iThe confidential reference items matrix of representing i video camera, R _i, t _iBe rotation and the displacement of camera coordinate system in world coordinate system, Refer to that both members differs a scale factor,

Transformational relation on the scaling board between characteristic point coordinate and the world coordinate system is suc as formula (2)

More than in two formulas, P _iBe camera parameters matrix, Q _jBe the location parameter matrix of scaling board, (x y) represents the actual coordinate of characteristic point on the scaling board;

B) the actual coordinate data of characteristic point on view data and the known scaling board are done standardization processing, and ask for the homography matrix of video camera by formula (3) ${\tilde{H}}_i^j={\mathrm{\&lambda;}}_i^j{H}_i^j,$

H wherein _i ^jRepresent the true homography matrix between i video camera and j the scaling board image.

3. the camera array calibration method based on matrix decomposition according to claim 1, it is characterized in that the said method that the homography matrix of trying to achieve is measured matrix according to the numbering structure of its corresponding video camera and scaling board is, with step 2) in the video camera that obtains and the homography matrix between the scaling board

Arrange by mode shown in the formula (4) according to subscript on it, obtain measuring matrix

Be the label of " rower " expression video camera of matrix, the label of " row mark " expression scaling board image, each

All be the matrix of (3 * 3),

$\tilde{W}=\left[\begin{array}{ccc}\widetilde{H_1^1}& \&CenterDot;\&CenterDot;\&CenterDot;& \widetilde{H_1^J}\\ \&CenterDot;& & \&CenterDot;\\ \&CenterDot;& & \&CenterDot;\\ \&CenterDot;& & \&CenterDot;\\ \widetilde{H_I^1}& \&CenterDot;\&CenterDot;\&CenterDot;& \widetilde{H_I^J}\end{array}\right]---\left(4\right)$

Wherein

Homography matrix that the expression Practical Calculation is come out and real homography matrix H _i ^jDiffer a proportionality factors lambda _i ^j, that is: ${\tilde{H}}_i^j={\mathrm{\&lambda;}}_i^j{H}_i^j,$ I represents the number of video camera, and J represents the number of the scaling board image taken.

4. the camera array calibration method based on matrix decomposition according to claim 2, it is characterized in that said intensive or sparse arranging situation according to video camera array, adjust the scale factor of homography matrix, and the measurement matrix of obliterated data carried out matrix decomposition, its step is as follows respectively by intensive or sparse arranging situation:

A) when the video camera array dense arrangement

(1) scale factor of adjustment homography matrix

Adjust ${\tilde{H}}_i^j={\mathrm{\&lambda;}}_i^j{H}_i^j$ In proportionality factors lambda _i ^j, method is as follows:

Concerning l and f video camera, they are respectively about the relevant homography matrix N of r and p scaling board image _Lf ^RpForm is as follows:

$N_{\mathrm{lf}}^{\mathrm{rp}}=H_f^r{\left(H_l^r\right)}^{-1}{H}_l^p{\left(H_f^p\right)}^{-1}$

The relevant homography matrix that has scale factor that then calculates

For:

${\tilde{N}}_{\mathrm{lf}}^{\mathrm{rp}}={\tilde{H}}_f^r{\left({\tilde{H}}_l^r\right)}^{-1}{\tilde{H}}_l^p{\left({\tilde{H}}_f^p\right)}^{-1}=\left({\mathrm{\&lambda;}}_f^r{\mathrm{\&lambda;}}_l^p\right)/\left({\mathrm{\&lambda;}}_l^r{\mathrm{\&lambda;}}_f^p\right)N_{\mathrm{lf}}^{\mathrm{rp}}---\left(5\right)$

Having a value is u _Ll ^RpDouble characteristic root suc as formula (6)

$u_{\mathrm{ll}}^{\mathrm{rp}}=\left({\mathrm{\&lambda;}}_f^r{\mathrm{\&lambda;}}_l^p\right)/\left({\mathrm{\&lambda;}}_l^r{\mathrm{\&lambda;}}_f^p\right)---\left(6\right)$

Measurement matrix with the step 3) structure

In element press row and column ordering, make the shooting function of its first three rows correspondence photograph all scaling boards, the scaling board of first three columns correspondence can be taken simultaneously by all video cameras, and first video camera in the measurement matrix of fixing this moment and first scaling board be as the reference video camera of the relevant homography matrix of calculating with reference to scaling board, and then first video camera with first scaling board with the relevant homography matrix of j scaling board image of i video camera is: ${\tilde{N}}_{l1}^{r1}=u_{l1}^{r1}{N}_{l1}^{r1},$ And u _L1 ^RlBy

Double characteristic root calculate, the note u _I1 ^J1Be u _i ^j,

Right

In each element do following adjustment,

Order ${\tilde{H}}_i^j=u_i^j{\tilde{H}}_i^j=u_i^j{\mathrm{\&lambda;}}_i^j{H}_i^j=({\mathrm{\&lambda;}}_i^1{\mathrm{\&lambda;}}_1^j/{\mathrm{\&lambda;}}_1^1)P_i{Q}^j---\left(7\right)$

Then adjusted

For:

$\tilde{W}=\left[\begin{array}{ccc}\widetilde{H_1^1}& \&CenterDot;\&CenterDot;\&CenterDot;& \widetilde{H_1^J}\\ \&CenterDot;& & \&CenterDot;\\ \&CenterDot;& & \&CenterDot;\\ \&CenterDot;& & \&CenterDot;\\ \widetilde{H_I^1}& \&CenterDot;\&CenterDot;\&CenterDot;& \widetilde{H_I^J}\end{array}\right]$

(2) the measurement matrix of obliterated data

Decomposition

Undertaken by following three steps:

1) be 4 to an order, size is the measurement matrix of 3I * 3J

The 3I of given any initial value * 4 matrix P ₀

2) at P _K-1Under the known situation, by formula $B{(j,:)}^T=\left(\tilde{W}{(:,j)}^{T}A\right){\left(A^{T}A\right)}^{-1},(j=1\&CenterDot;\&CenterDot;\&CenterDot;3J)$ Ask for 3J * 4 matrix Q _k, and to Q _kColumn vector do orthogonalization process: Q _k=Q _k* N _k, make

Minimum, the ij representative is right

The factor summation that exists;

3) at Q _kUnder the known situation, by formula $A{(i,:)}^T=\left(\tilde{W}(i,:)B\right){\left(B^{T}B\right)}^{-1},(i=1\&CenterDot;\&CenterDot;\&CenterDot;3I)$ Ask for 3I * 4 matrix P _k, make

Minimum, the ij representative is right

The factor summation that exists;

4) repeat above 2), 3) two steps are up to P _kQ _k ^TConvergence, and, after the contrary ordering of Q, obtain measuring matrix to P

Decomposed form suc as formula (8):

$\tilde{W}=(1/{\mathrm{\&lambda;}}_1^1)\left[\begin{array}{c}{\mathrm{\&lambda;}}_1^1{P}_1\\ \&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\\ {\mathrm{\&lambda;}}_I^1{P}_I\end{array}\right]\left[\begin{array}{ccc}{\mathrm{\&lambda;}}_1^1{Q}^{1^T}& \&CenterDot;\&CenterDot;\&CenterDot;& {\mathrm{\&lambda;}}_1^J{Q}^{J^T}\end{array}\right]---\left(8\right)$

When iteration does not restrain or error When surpassing threshold value, change iterative initial value P ₀Again iteration gets final product;

B) when the arrangement of video camera array is sparse

(1) the measurement matrix to constructing in the step 3)

Reorder, make

The shooting function of first three rows correspondence photograph the scaling board image of maximum numbers, the scaling board of first three columns correspondence can be photographed simultaneously by the video camera of maximum numbers, and first video camera in the measurement matrix of fixing this moment and first scaling board be as the reference video camera of the relevant homography matrix of calculating with reference to scaling board, and then first video camera with first scaling board with the relevant homography matrix of j scaling board image of i video camera is: ${\tilde{N}}_{l1}^{r1}=u_{l1}^{r1}{N}_{l1}^{r1},$ And u _L1 ^RlBy

Double characteristic root calculate, the note u _I1 ^J1Be u _i ^j, right Middle first three rows and first three columns are not empty

With

Those corresponding row and columns Do following adjustment,

Order ${\tilde{H}}_i^j=u_i^j{\tilde{H}}_i^j=u_i^j{\mathrm{\&lambda;}}_i^j{H}_i^j=({\mathrm{\&lambda;}}_i^1{\mathrm{\&lambda;}}_1^j/{\mathrm{\&lambda;}}_1^1)P_i{Q}^j---\left(7\right)$

Then adjusted

For:

$\tilde{W}=\left[\begin{array}{ccc}\widetilde{H_1^1}& \&CenterDot;\&CenterDot;\&CenterDot;& \widetilde{H_1^J}\\ \&CenterDot;& & \&CenterDot;\\ \&CenterDot;& & \&CenterDot;\\ \&CenterDot;& & \&CenterDot;\\ \widetilde{H_I^1}& \&CenterDot;\&CenterDot;\&CenterDot;& \widetilde{H_I^J}\end{array}\right]$

(2) the measurement matrix of obliterated data

Decomposition

Undertaken by following step:

1) be 4 to an order, size is the matrix of 3I * 3J

The 3I of given any initial value * 4 matrix P ₀

2) establish each

Scale factor

3) utilize current λ _i ^j, make up and measure matrix

4) at P _K-1Under the known situation, by formula $Q_k{(j,:)}^T=\left(\tilde{W}{(:,j)}^T{P}_{k-1}\right){\left({P_{k-1}}^T{P}_{k-1}\right)}^{-1},(j=1\&CenterDot;\&CenterDot;\&CenterDot;3J)$ Ask for 3J * 4 matrix Q _k, and to Q _kColumn vector do orthogonalization process: Q _k=Q _k* N _k, make

Minimum, the ij representative is right

The factor summation that exists;

5) at Q _kUnder the known situation, by formula $P_k{(i,:)}^T=\left(\tilde{W}(i,:)Q_k\right){\left({Q_k}^T{Q}_{,}\right)}^{-1},(i=1\&CenterDot;\&CenterDot;\&CenterDot;3I)$ Ask for 3I * 4 matrix P _k, make

Minimum, the ij representative is right

The factor summation that exists;

6) adjust

Middle first three rows and first three are classified as empty

With

Those corresponding row and columns

Scale factor as follows: ${\mathrm{\&lambda;}}_i^j=\mathrm{\&Sigma;}{P}_i(3,:)\&times;Q_j{(3,:)}^T,$ P wherein _i(3 :) represent P _iThe third line, Q _j(3 :) represent Q ^jThe 3rd row;

7) repeating step 3) to 6) up to

Convergence, and to P, after the contrary ordering of Q, the decomposed form that obtains is suc as formula (8):

$\tilde{W}=(1/{\mathrm{\&lambda;}}_1^1)\left[\begin{array}{c}{\mathrm{\&lambda;}}_1^1{P}_1\\ \&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\\ {\mathrm{\&lambda;}}_I^1{P}_I\end{array}\right]\left[\begin{array}{ccc}{\mathrm{\&lambda;}}_1^1{Q}^{1^T}& \&CenterDot;\&CenterDot;\&CenterDot;& {\mathrm{\&lambda;}}_1^J{Q}^{J^T}\end{array}\right]---\left(8\right)$

When iteration does not restrain or error

When surpassing threshold value, change iterative initial value P ₀Again iteration gets final product.

5. the camera array calibration method based on matrix decomposition according to claim 2 is characterized in that the constraint of said interpolation tolerance, and the step of inside and outside parameter of finding the solution video camera is as follows:

A) to P, Q goes standardization, obtains measuring matrix Decomposed form is as follows:

$\tilde{W}\&cong;\left[\begin{array}{c}{P}_1^{\&prime;}\\ \&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\\ P_I^{\&prime;}\end{array}\right]\left[\begin{array}{ccc}{Q}^{1^{\&prime;}}& \&CenterDot;\&CenterDot;\&CenterDot;& Q^{J^{\&prime;}}\end{array}\right]$

B) ask for the video camera confidential reference items

Real camera parameters matrix P and actual P ' relation as (9) formula that obtains of decomposing:

$P_i^{\&prime;}\&cong;P_{i}T;$ $Q^{j^{\&prime;}}\&cong;T^{-1}{Q}^j---\left(9\right)$

Wherein T is a nonsingular matrix, can determine as follows:

Concerning first video camera, R ₁Be unit matrix, t ₁Each component be 0, P then ₁=A ₁R ₁ ^T[I-t ₁]=[K ₁0], with the matrix P after decomposing season ₁' have a form P ₁'=[, I 0], so can get $T=\left[\begin{array}{cc}{K}_1^{-1}& 0\\ h^T& h\end{array}\right],$ H wherein ^TBe one 3 * 1 vector, h is a scale factor,

If the scale factor in the formula (9) is β ^j, then formula (9) can be written as: TQ ^j'=β ^jQ ^j,

Q to each scaling board correspondence ^j, by the p in the formula (2), the orthogonality of q draws formula (10):

p ^j′Tω ₁p ^j′＝q ^j′Tω ₁q ^j′＝(β ^j) ²，p ^j′Tω ₁q ^j′＝0 (10)

Wherein ${\mathrm{\&omega;}}_1=K_1^{-T}{K}_1^{-1},$ Obtain the quadrature equation of first all scaling boards of photographing of shooting function by formula (10), and all equations that superpose, obtain K ₁Thereby:

${\mathrm{\&beta;}}^j=\sqrt{(p^{j^{\&prime;}T}{\mathrm{\&omega;}}_1{p}^{j^{\&prime;}}+q^{j^{\&prime;}T}{\mathrm{\&omega;}}_1{q}^{j^{\&prime;}})/2}---\left(11\right)$

Again, each scaling board all satisfies equation (12):

[h ^T h]Q ^j′＝[0 0 β ^j] (12)

The above-mentioned equation (12) that superposes, [h ^TH], and then obtain transformation matrix T, so first camera parameters matrix P is arranged ₁=P ₁' T ^-1, the parameter of other video camera also can be by P _i=P _i' T ^-1Obtain.

C) ask for the outer ginseng of video camera

External parameters of cameras is tried to achieve as follows:

$P_i=A_i{R}_i^T[I-t_i]=\left[M\right|b]---\left(13\right)$

So: $R_i^{\&prime;}=M^T{A}_i^{-T},$ And to R _i' be and try to achieve spin matrix R after SVD decomposes _iSuc as formula (14),

$R_i^{\&prime;}={\mathrm{USV}}^T\&DoubleRightArrow;R_i={\mathrm{UV}}^T---\left(14\right)$

Displacement between video camera: t _i=-M ^-1 _iB (15).

6. a kind of camera array calibration method based on matrix decomposition according to claim 2 is characterized in that said step of the gained inside and outside parameter being carried out nonlinear optimization:

Carry out nonlinear optimization by formula (16):

$\min \underset{i=1}{\overset{I}{\mathrm{\&Sigma;}}}\underset{j=1,i}{\overset{J}{\mathrm{\&Sigma;}}}\underset{n}{\overset{\mathrm{Nj}}{\mathrm{\&Sigma;}}}{\left|\right|m_{\mathrm{in}}^j-\hat{m}(K_i,R_i,t_i,Q^j,M_n^j)\left|\right|}^2---\left(16\right)$

In the formula,

Be with a certain characteristic point M on the scaling board _n ^jAgain project to coordinate in the image coordinate system, m by the camera model of trying to achieve _In ^jBe M _n ^jPosition coordinates in actual photographed image.

Images