CN103745474A - Image registration method based on inertial sensor and camera - Google Patents

Abstract

The invention puts forward an image registration method based on an inertial sensor and a camera. The camera and the inertial sensor are synchronously controlled to operate to simultaneously acquire the color image and depth image of the same scene and the three-axis motion angular velocity and three-axis motion acceleration of the camera; Kalman filtering is utilized to carry out filtering on the three-axis motion angular velocity and the three-axis motion acceleration, the Eulerian dynamics equation is used for converting the three-axis acceleration into a Eulerian angle, the uniformly accelerated motion formula is used for converting the three-axis angular velocity into three-axis displacement, and thereby an external reference matrix of the camera is obtained; by utilizing the affine transformation relation between the color image and the depth image and the proportional relation between depth image grey value and pixel depth information, the depth information of angular points in the depth image is solved; an OFCM (Orthogonal Frequency-division Multiplexing) model is used for calculation to obtain image background optical flow; the RANSAC (Random Sample Consensus) algorithm is used for completing image registration. The calculation load of the method is low, both the real-timeness and the precision are high, and the method is applicable to hand-held devices and mobile devices.

Description

Based on the method for registering images of inertial sensor and video camera

Technical field

The invention belongs to digital image processing field, be specifically related to a kind of method for registering images based on inertial sensor and video camera.

Background technology

The key issue of image registration techniques is the estimation to camera motion, generally by camera calibration, completes.In prior art, the method for camera self-calibration has:

The people such as Li Lichun have proposed based on 8 method basis matrix three-dimensional reconstruction algorithm [Li Lichun, Qiu Zhiqiang, king roc etc. based on match measure, weighting solves the three-dimensional reconstruction algorithm [J] of basis matrix. computer utility, 2007,27 (10): 2530-2533], the method complexity and calculated amount are very large, are difficult to realize the real-time of image processing.

The people such as Zhang Yang the have proposed FOE modelling verification translation effect [Zhang Yang of camera motion, Kai Sele Vios, Bai Sen. the moving object detection algorithm based on mobile camera in visible-range. computer vision and the .IEEE of pattern-recognition operating room international conference, 2006,130-131], the method is not considered the internal reference matrix of rotatablely moving of video camera and video camera, and precision is not high.

The scaling method that the people such as Zheng Hong have proposed based on active vision is mainly to utilize the hardware platform of the known motions such as The Cloud Terrace that camera motion information [Zheng Hong is provided, Liu Zhenqiang. the monopod video camera self-calibration [J] based on accurate model. robot, 2013.5,35 (3): 1994-2013], the method is because the equipment such as The Cloud Terrace are larger, although can realize well real-time, be difficult to spread in handheld device and mobile device, use face is limited.

Summary of the invention

The present invention proposes a kind of method for registering images based on inertial sensor and video camera, and the method calculated amount is little, and real-time and precision are all higher, can be applicable to handheld device and mobile device.

In order to solve the problems of the technologies described above, the invention provides a kind of method for registering images based on inertial sensor and video camera, comprise the following steps:

The coloured image I of Same Scene and the three-axis moving angular velocity [w of J, depth image I ' and J ' and video camera are obtained in step 1, synchro control video camera and inertial sensor work simultaneously _x, w _y, w _z] and three-axis moving acceleration [a _x, a _y, a _z];

The angle point m of step 2, extraction coloured image I ₁(x ₁, y ₁);

Step 3, utilize the three axis angular rate [ws of Kalman filtering to video camera _x, w _y, w _z] and 3-axis acceleration [a _x, a _y, a _z] carry out filtering, obtain Kalman Filter Estimation three axis angular rates with Kalman Filter Estimation 3-axis acceleration use euler dynamical equations formula will

convert Eulerian angle [α to _k, β _k, γ _k], use equation of uniformly accelerated motion will

convert three-shaft displacement [t to _xk, t _yk, t _zk], thereby obtain video camera, join rotation matrix R and the translation matrix T in matrix [R, T] outward;

Step 4, utilize the affine transformation relationship of coloured image and depth image, solve the angle point m of coloured image I ₁(x ₁, y ₁) coordinate position m in depth image I ' and J ' _1ir(x ₁, y ₁) and m _2ir(x ₁, y ₁); Utilize the proportionate relationship of depth image gray-scale value and pixel depth information to solve angle point m ₁(x ₁, y ₁) depth information u in depth image I ' and J ' ₁and u ₂;

The depth information u of the coloured image that step 5, the rotation matrix R that step 3 is obtained and translation matrix T and step 4 obtain ₁, u ₂bring OFCM model into, calculate the background light stream m that obtains image J ₂(x ₂, y ₂), OFCM model as shown in Equation (1),

$m_2=\frac{u_1}{u_2}{\mathrm{KRK}}^{-1}{m}_1+\frac{1}{u_2}\mathrm{KT}---\left(1\right)$

In formula (1), the internal reference matrix that K is video camera, K ^-1for the inverse matrix of K;

The background light stream m of step 6, use RANSAC color image I ₂carry out Robust estimation, by background light stream m ₂the light stream being produced by moving object of middle existence and mismatch quasi-optical stream and weed out, thereby complete the image registration of coloured image I and coloured image J, then utilize the correct background light stream of the coloured image I retaining to carry out initialization to coloured image J, carry out the registration of successive image.

Main inventive principle:

1, the present invention has proposed the relational model (OFCM) of image background light stream and camera motion in conjunction with the rigid body variation relation of theorem in Euclid space, and this model mode of mainly deriving is as follows:

Supposing that colored CCD is pin-hole model, only there is the relative motion with colored CCD in the background in scene, and the foundation of camera coordinate system is O centered by the photocentre of colored CCD, and the straight line at optical axis place is Z axis, then sets up X, Y-axis according to left hand principle.

Before supposing camera motion, camera coordinate system overlaps with world coordinate system, and the coordinate of any point P in space background meets so:

P _W＝P _C (2)

In formula (2), P _wfor the world coordinate system coordinate of space point, P _wfor the camera coordinate system coordinate of space point, the homogeneous coordinates m of the imaging point on space background point P and its imaging plane ₁(x ₁, y ₁, 1) and meet formula (3):

u ₁m ₁＝KP _C (3)

In formula (3), m ₁(x ₁, y ₁) be the imaging point on imaging plane, while utilizing model to calculate, adopt Harris angle point to replace; u ₁can be approximately m ₁(x ₁, y ₁) depth information of corresponding spatial point; K is the internal reference matrix of video camera, considers the relation that in image, pixel coordinate is and physical coordinates is, can draw the five-parameter model of video camera internal reference matrix, as shown in Equation (4),

$K=\left[\begin{array}{ccc}{f}_x& s& u_0\\ 0& f_y& v_0\\ 0& 0& 1\end{array}\right]---\left(4\right)$

In formula (4),

wherein, dx is the physical size of pixel in x direction, and dy is the physical size of pixel in y direction, [u ₀, v ₀] be the coordinate of principal point in pixel coordinate system, s is camera distortion coefficient.

Suppose that video camera exists the translation motion of any direction and rotatablely moves, translation motion can be used translation vector T=[t _x, t _y, t _z] ^trepresent, rotatablely move and can represent with the rotation matrix R of 3*3, now the coordinate of space background point P in camera coordinate system meets:

P' _C＝R×P _W+T＝R×P _C+T (5)

In formula (5), rotation matrix R can be with the Eulerian angle α of theorem in Euclid space, and beta, gamma represents, as the formula (6),

$R=\left[\begin{array}{ccc}\cos \mathrm{\γ}\cos \mathrm{\β}& \cos \mathrm{\γ}\sin \mathrm{\β}\sin \mathrm{\α}-\sin \mathrm{\γ}& \cos \mathrm{\γ}\sin \mathrm{\β}\cos \mathrm{\α}-\sin \mathrm{\α}\sin \mathrm{\γ}\\ \sin \mathrm{\γ}\cos \mathrm{\β}& \sin \mathrm{\α}\sin \mathrm{\β}\sin \mathrm{\γ}+\cos \mathrm{\γ}\cos \mathrm{\α}& \sin \mathrm{\γ}\sin \mathrm{\β}\cos \mathrm{\α}+\cos \mathrm{\γ}\sin \mathrm{\α}\\ -\sin \mathrm{\β}& \cos \mathrm{\β}\sin \mathrm{\α}& \cos \mathrm{\β}\cos \mathrm{\α}\end{array}\right]---\left(6\right)$

After camera motion, the homogeneous coordinates m of the imaging point on spatial point P and its imaging plane ₂(x ₂, y ₂, 1) and meet relation:

u ₂m ₂＝KP' _C＝K(R×P _C+T) (7)

According to formula (3) and formula (5), can obtain the motion light stream of camera motion front and back spatial point P on imaging plane, i.e. the final form of OFCM:

$m_2=\frac{u_1}{u_2}{\mathrm{KRK}}^{-1}{m}_1+\frac{1}{u_2}\mathrm{KT}---\left(8\right)$

Wherein, u ₁and u ₂can be approximated to be the depth information of spatial point P in camera coordinate system, the depth image obtaining by video camera calculates.Outer ginseng matrix [R, T] is obtained and is processed by hardware platform and draws, m ₁can obtain by Harris Corner Detection.Degree of depth CCD and the colored CCD of supposing video camera are all operated under the condition of fixed focal length f, and this time image plane and object plane are one to one.Therefore, the internal reference matrix K of colored CCD is an invariable amount in whole image acquisition procedures, can be through once demarcating after obtaining and use as experience matrix.

Advantage: the camera model existing now has FOE model, MOF model, FOE model and MOF model are not all considered the internal reference matrix of video camera, and FOE model is only supported the translation motion of video camera, therefore under camera complex motion, cannot keep higher precision.The internal reference matrix that the OFCM model that the present invention proposes comprises video camera, supports three axle translation motions and three axles of video camera to rotatablely move, and can overcome scene brightness sudden change.Therefore, OFCM model, in camera complex motion and scene brightness sudden change situation, can keep high-precision image registration.

2, the hardware platform based on inertial sensor and video camera:

The present invention utilizes OFCM model to carry out the background light stream of computed image, and from OFCM model, the calculating of image background light stream need to be known outer ginseng matrix [R, T] and the depth information u corresponding to image pixel of video camera ₁, u ₂.The present invention adopts hardware platform to complete video camera to join obtaining of matrix and image pixel depth information outward.

Hardware platform of the present invention comprises three parts: control device FPGA, inertial sensor and video camera, video camera can adopt Kinect video camera, inertial sensor comprises gyroscope and accelerometer, Kinect video camera comprises that depth image obtains CCD and coloured image obtains CCD, the minimum system of FPGA and inertial sensor is all fixed on depth image and obtains on CCD, inertial sensor be positioned at Kinect video camera directly over, gyroscope is for obtaining the three-axis moving angular velocity [w of video camera _x, w _y, w _z], accelerometer is for obtaining the three-axis moving acceleration [a of video camera _x, a _y, a _z], Kinect video camera is for color image shot I, J and depth image I ', J '.After hardware platform powers on, FPGA is by the shutter work of gating pulse control Kinect video camera, thereby complete obtaining of coloured image and depth image, FPGA is by SPI protocol synchronization control Kinect video camera and inertial sensor work, corresponding one by one to guarantee the camera motion of image sequence that this step is obtained and collection.

Advantage: the image registration algorithm calculated amount based on this model is little, and system real time is good; Because inertial sensor has, volume is little, precision advantages of higher, obtains the equipment of camera motion compare with The Cloud Terrace etc., can popularize better in handheld device and mobile device.

3, the foundation of the Kalman filter model based on inertial sensor:

There is zero migration, scale factor error and random noise equal error in the camera motion collecting due to inertial sensor, therefore needs the camera motion to collecting to carry out filtering and noise reduction processing.

When using the classical formulas of Kalman filtering, need to set up according to the performance of inertial sensor and parameter error model and the state transitions relational expression of inertial sensor:

z′ _k=z _k(M+1)+B _f+U (9)

z _k=Ψ _k,k-1z _k-1+V (10)

Wherein, z ' _kthe actual observed value of etching system during for k, z _kthe actual value of etching system during for k, M is scale factor error, B _ffor zero point drift, V and U are random Gaussian white noise, Ψ _{k, k-1}for state-transition matrix.

According to formula (9) binding time sequence analysis, solve the scale factor error M of inertial sensor, zero point drift B _f, and white Gaussian noise V and U; According to formula (10) binding time sequence analysis, solve the state-transition matrix Ψ between parameter _{k, k-1}.

The result of formula formula (9) and formula (10) calculating is brought into the classical formulas (11) of Kalman filtering, the parameter in k moment is carried out to filtering, obtain estimated value

$\left\{\begin{array}{c}{z}_{k,k-1}={\mathrm{\Ψ}}_{k,k-1}{z}_{k-1}\\ \widehat{z_k}=z_{k,k-1}+K_k\left[{z^{\′}}_k(M+1)z_{k,k-1}\right]\\ P_{k,k-1}={\mathrm{\Ψ}}_{k,k-1}{P}_{k-1}{{\mathrm{\Ψ}}_{k,k-1}}^T+Q_{k-1}\\ K_k=P_{k,k-1}{(M+1)}^T{P}_k^{-1}\\ P_k=[I-K_k(M+1)]P_{k,k-1}\end{array}\right.---\left(11\right)$

Wherein, the noise variance that Q is system, K _kfor kalman gain, P _kfor state covariance, w _{k, k-1}for state discreet value.

Compared with prior art, its remarkable advantage is in the present invention:

(1) because the estimation of camera parameters is by inertial sensor and video camera Real-time Obtaining, wherein the gyroscope in inertial sensor and accelerometer obtain respectively rotational motion parameter and the translation motion parameter of video camera; Video camera comprises that depth image obtains CCD (degree of depth CCD) and coloured image obtains CCD (colored CCD), and the two can be in conjunction with the depth information that obtains coloured image.Although therefore this method is based on hardware platform, because inertial sensor has, volume is little, precision advantages of higher, after video camera binding, can popularize well in handheld device and mobile device.

(2) from photography geometric angle, in conjunction with the rigid body variation relation of theorem in Euclid space, derive the relational model (optical flow and camera ' s motion) of image background light stream and camera motion, be denoted as OFCM model.The internal reference matrix that this model comprises video camera, and support three axle translation motions and three axles of video camera to rotatablely move; Meanwhile, the in the situation that of jump in brightness, still there is higher precision, compared to FOE model and light stream matching method, had larger breakthrough.

Accompanying drawing explanation

Fig. 1 is the method for registering images process flow diagram that the present invention is based on inertial sensor and video camera.

Fig. 2 is the image registration design sketch that uses the inventive method to obtain in scene jump in brightness situation.

Fig. 3 is the image registration design sketch that uses the inventive method to obtain when video camera generation translation motion.

Fig. 4 is the image registration design sketch that uses the inventive method to obtain when video camera rotates motion.

Fig. 5 is in video camera generation translation motion and the image registration design sketch that uses while rotatablely moving the inventive method to obtain.

Embodiment

As shown in Figure 1, the present invention is based on the method for registering images of inertial sensor and video camera, comprise the following steps

Step 1, in hardware platform, synchro control video camera and inertial sensor work are obtained the coloured image I of Same Scene and the kinematic parameter of J, depth image I ' and J ' and video camera, i.e. three-axis moving angular velocity [w simultaneously _x, w _y, w _z] and three-axis moving acceleration [a _x, a _y, a _z], be implemented as follows:

Hardware platform of the present invention comprises three parts: control device FPGA, inertial sensor and video camera, inertial sensor comprises gyroscope and accelerometer, video camera comprises that depth image obtains CCD and coloured image obtains CCD, for color image shot I, J and depth image I ', J ', the minimum system of FPGA and inertial sensor is all fixed on depth image and obtains CCD above, and gyroscope is for obtaining the three-axis moving angular velocity [w of video camera _x, w _y, w _z], accelerometer is for obtaining the three-axis moving acceleration [a of video camera _x, a _y, a _z].After hardware platform powers on, FPGA is by the shutter work of gating pulse control video camera, thereby complete obtaining of coloured image and depth image, FPGA is by SPI protocol synchronization control video camera and inertial sensor work, corresponding one by one to guarantee the camera motion of image sequence that this step is obtained and collection.

Step 2, utilize Harris Corner Detection to extract the angle point m of coloured image I ₁(x ₁, y ₁).Be implemented as follows:

2.1 calculate pixel p[x in coloured image I, y] correlation matrix M, computing method as shown in formula (1),

$M=\underset{x,y}{\mathrm{\Σ}}w(x,y)\left[\begin{array}{cc}{I}_x^2& I_x{I}_y\\ I_x{I}_y& I_y^2\end{array}\right]=w(x,y)\⊗\left[\begin{array}{cc}{I}_x^2& I_x{I}_y\\ I_x{I}_y& I_y^2\end{array}\right]---\left(2\right)$

In formula (1), (x, y) is the image coordinate of pixel p, and w (x, y) is that window size is Gauss's weighting function of w*w, I _xrepresent the derivative of line direction, I _yrepresent the derivative of column direction;

2.2 calculating pixel point m ₁(x ₁, y ₁) angle point response R (x, y), computing method as shown in formula (2),

R(x,y)＝detM-k*(traceM) ² (2)

In formula (2), detM represents the value of correlation matrix M determinant, and traceM represents the mark of correlation matrix M, and k is constant, is generally 0.04-0.06;

2.3 judge pixel m ₁(x ₁, y ₁) whether be angle point, determination methods is shown in formula (3), if pixel m ₁(x ₁, y ₁) angle point response R (x, y) meet formula (3), pixel m ₁(x ₁, y ₁) be angle point; If pixel m ₁(x ₁, y ₁) angle point response R (x, y) do not meet formula (11), by m ₁(x ₁, y ₁) give up;

[R(x,y)≥R(x′,y′)]∩[R(x,y)≥0.01R _max] (3)

In formula (3), and R (x ', y ') represent in window w*w except pixel m ₁(x ₁, y ₁) outside the angle point response of other pixel, R _maxfor edge feature frame I _fin the response of maximum angle point.

Step 3, utilize the three axis angular rate [ws of Kalman filtering to video camera _x, w _y, w _z] and 3-axis acceleration [a _x, a _y, a _z] carry out filtering, obtain Kalman Filter Estimation three axis angular rates

with Kalman Filter Estimation 3-axis acceleration use euler dynamical equations formula will

convert Eulerian angle [α to _k, β _k, γ _k], use equation of uniformly accelerated motion will

convert three-shaft displacement [t to _xk, t _yk, t _zkthereby] obtain video camera join rotation matrix R and the translation matrix T in matrix [R, T] outward;

3.1 utilize the kinematic parameter of Kalman filtering to video camera, i.e. three axis angular rate [w _x, w _y, w _z] and 3-axis acceleration [a _x, a _y, a _z] carry out filtering, obtain Kalman Filter Estimation three axis angular rates with Kalman Filter Estimation 3-axis acceleration

concrete account form as follows:

3.1.1 error model and the state transitions relational expression of according to the performance of inertial sensor and kinematic parameter, setting up inertial sensor, binding time sequence analysis solves scale factor error M, the zero point drift B of inertial sensor _f, the state-transition matrix Ψ between white Gaussian noise U, kinematic parameter _{k, k-1}and white Gaussian noise V;

Error model as shown in Equation (4), shown in state transitions relational expression formula (5),

z′ _k=z _k(M+1)+B _f+U (4)

z _k=Ψ _k,k-1z _k-1+V (5)

In formula (4) and formula (5), z represents three axis angular rate [w _x, w _y, w _z] and 3-axis acceleration [a _x, a _y, a _z] in a certain kinematic parameter, z ' _kfor k moment inertial sensor is measured the measured value of a certain kinematic parameter obtaining, z _kfor the actual value of k a certain kinematic parameter of moment, M is scale factor error, B _ffor zero point drift, V and U are random Gaussian white noise, Ψ _{k, k-1}for state-transition matrix;

According to formula (4) binding time sequence analysis, solve the scale factor error M of inertial sensor, zero point drift B _f, and white Gaussian noise V and U; According to formula (5) binding time sequence analysis, solve the state-transition matrix Ψ between kinematic parameter _{k, k-1};

3.1.2 the kinematic parameter z of the classical formulas (6) that the result of step 3.1.1 being asked for is brought Kalman filtering into the k moment _kcarry out filtering, obtain the estimated value of this kinematic parameter

$\left\{\begin{array}{c}{z}_{k,k-1}={\mathrm{\Ψ}}_{k,k-1}{z}_{k-1}\\ \widehat{z_k}=z_{k,k-1}+K_k\left[{z^{\′}}_k(M+1)z_{k,k-1}\right]\\ P_{k,k-1}={\mathrm{\Ψ}}_{k,k-1}{P}_{k-1}{{\mathrm{\Ψ}}_{k,k-1}}^T+Q_{k-1}\\ K_k=P_{k,k-1}{(M+1)}^T{P}_k^{-1}\\ P_k=[I-K_k(M+1)]P_{k,k-1}\end{array}\right.---\left(6\right)$

In formula (6), the noise variance that Q is system, K _kfor kalman gain, P _kfor state covariance, z _{k, k-1}for state discreet value.

According to the method described above respectively to three-axis moving angular velocity [w _x, w _y, w _z] and three-axis moving acceleration [a _x, a _y, a _z] carry out Kalman filtering, obtain with

3.2 use euler dynamical equations by filtered k moment three-axis moving angular velocity w _xk, w _yk, w _zkconvert the Eulerian angle α in European coordinate system to _k, β _k, γ _k, by the α of gained _k, β _k, γ _kbring in rotation matrix R, complete the calculating of colored CCD rotation matrix R;

Described euler dynamical equations as shown in formula (7),

$\left[\begin{array}{c}{\stackrel{\·}{\mathrm{\α}}}_k\\ {\stackrel{\·}{\mathrm{\β}}}_k\\ {\stackrel{\·}{\mathrm{\γ}}}_k\end{array}\right]=\left[\begin{array}{ccc}-\sin {\mathrm{\α}}_k/\tan {\mathrm{\β}}_k& {\cos \mathrm{\α}}_k/\tan {\mathrm{\β}}_k& 1\\ {\cos \mathrm{\α}}_k& \sin {\mathrm{\α}}_k& 0\\ \sin {\mathrm{\α}}_k/\sin {\mathrm{\β}}_k& -\cos {\mathrm{\α}}_k/\sin {\mathrm{\β}}_k& 0\end{array}\right]*\left[\begin{array}{c}{w}_{\mathrm{xk}}\\ w_{\mathrm{yk}}\\ w_{\mathrm{zk}}\end{array}\right]---\left(7\right)$

In formula (7), for Eulerian angle α _k, β _k, γ _kfirst order derivative, the parameter α that the k-1 moment is tried to achieve _k-1, β _k-1, γ _k-1as the initial value of formula (7);

Described rotation matrix R is as shown in formula (8):

$R=\left[\begin{array}{ccc}\cos {\mathrm{\γ}}_k\cos {\mathrm{\β}}_k& \cos {\mathrm{\γ}}_k\sin {\mathrm{\β}}_k\sin {\mathrm{\α}}_k-\sin {\mathrm{\γ}}_k& \cos {\mathrm{\γ}}_k\sin {\mathrm{\β}}_k\cos {\mathrm{\α}}_k-\sin {\mathrm{\α}}_k\sin {\mathrm{\γ}}_k\\ \sin {\mathrm{\γ}}_k\cos {\mathrm{\β}}_k& \sin {\mathrm{\α}}_k\sin {\mathrm{\β}}_k\sin {\mathrm{\γ}}_k+\cos {\mathrm{\γ}}_k\cos {\mathrm{\α}}_k& \sin {\mathrm{\γ}}_k\sin {\mathrm{\β}}_k\cos {\mathrm{\α}}_k+\cos {\mathrm{\γ}}_k\sin {\mathrm{\α}}_k\\ -\sin {\mathrm{\β}}_k& \cos {\mathrm{\β}}_k\sin {\mathrm{\α}}_k& \cos {\mathrm{\β}}_k\cos {\mathrm{\α}}_k\end{array}\right]---\left(8\right)$

The speed displacement formula of 3.3 uses as shown in formula (9) is by 3-axis acceleration a _xk, a _yk, a _zkconvert three-shaft displacement t to _xk, t _yk, t _zk.Because the frequency acquisition of camera motion is higher, therefore the translation motion of video camera can be approximately to uniformly accelerated motion, computing formula is as shown in formula (9):

$\left\{\begin{array}{c}{t}_k=v_{k}T+\frac{1}{2}{a}_k{T}^2\\ v_k=v_{k-1}+a_{k-1}\×T\end{array}\right.----\left(9\right)$

In formula (9), t _krepresent three-shaft displacement t _xk, t _yk, t _zkin the displacement of any change in coordinate axis direction, the shooting cycle that T is video camera, a _k-1for the axial acceleration of k-1 moment respective coordinates, v _k-1and v _kbe respectively the initial velocity in k-1 moment and k moment.The present invention is to the initial velocity 0 of fixed system initial time, i.e. v ₀=0.By 3-axis acceleration a _xk, a _yk, a _zkbring respectively formula (9) into, can calculate the translation vector T=[t that obtains three axle translation motions _xk, t _yk, t _zk].

Step 4, utilize the affine transformation relationship of coloured image and depth image, solve respectively the angle point m of coloured image I ₁(x ₁, y ₁) coordinate position m in depth image I ' and J ' _1ir(x ₁, y ₁) and m _2ir(x ₁, y ₁); Utilize the proportionate relationship of gray-scale value and pixel depth information in depth image to solve angle point m ₁(x ₁, y ₁) depth information u in depth image I ' and J ' ₁and u ₂, be implemented as follows:

4.1 utilize the affine transformation relationship of coloured image and depth image, solve respectively the angle point m of coloured image I ₁(x ₁, y ₁) coordinate position m in depth image I ' and J ' _1ir(x ₁, y ₁) and m _2ir(x ₁, y ₁), be implemented as follows:

According to direct linear transformation's algorithm, calculate the affine transformation relationship of coloured image and depth image, concrete as shown in formula (10):

m _ir＝Hm (10)

In formula (10), m _irfor the pixel coordinate in depth image, m is the pixel coordinate in coloured image.By the angle point m of coloured image I ₁(x ₁, y ₁) bring formula (9) into and solve angle point m ₁(x ₁, y ₁) coordinate position m in depth image I ' and J ' _1ir(x ₁, y ₁) and m _2ir(x ₁, y ₁).

4.2 utilize the proportionate relationship of gray-scale value and pixel depth information in depth image to solve coloured image angle point m ₁(x ₁, y ₁) depth information u in depth image I ' and J ' ₁and u ₂, concrete computation process is as follows

The rule-of-thumb relation that the proportionate relationship of gray-scale value and pixel depth information can propose by reverse engineering, as shown in formula (11):

u＝l×m (11)

In formula (11), m (x, y) is the gray-scale value that pixel coordinate (x, y) is located, and u is that pixel coordinate (x, y) is located corresponding depth value, and l is scale-up factor; The m that step 4.2 is solved _1ir(x ₁, y ₁), m _2ir(x ₁, y ₁) bring formula (11) into and solve corresponding depth information u ₁and u ₂.

The depth information u of the coloured image of matrix R and T and step 4 acquisition is joined in step 5, the motion that step 3 is obtained outward ₁, u ₂bring OFCM model into, calculate the background light stream m that obtains image J ₂(x ₂, y ₂), be implemented as follows:

The depth information u of the coloured image of matrix R and T and step 4 acquisition is joined in the motion that step 3 is obtained outward ₁, u ₂bring OFCM model as shown in Equation (12) into, calculate and obtain k frame angle point m ₁background light stream m in k+1 two field picture ₂,

$m_2=\frac{u_1}{u_2}{\mathrm{KRK}}^{-1}{m}_1+\frac{1}{u_2}\mathrm{KT}---\left(12\right)$

In formula (12), the internal reference matrix that K is video camera, K ^-1for the inverse matrix of K.

The background light stream m of step 6, use RANSAC color image I ₂carry out Robust estimation, by background light stream m ₂the light stream being produced by moving object of middle existence and mismatch quasi-optical stream and weed out, thereby complete the image registration of coloured image I and coloured image J, then utilize the correct background light stream of the coloured image I retaining to carry out initialization to coloured image J, carry out the registration of successive image.Be implemented as follows:

The RANSAC computation process of described coloured image I is:

The current all angle point m of 6.1 initialization ₁(x ₁, y ₁), m ₂(x ₂, y ₂) be intra-office point;

6.2 utilize DLT algorithm to calculate affine transformation matrix H ₁, according to formula (13) in conjunction with angle point m ₁(x ₁, _y1) estimate angle point estimated value m ₂' (x ₂', y ₂'):

m ₂′=H ₁m ₁ (13)

6.3 calculate angle point m according to formula (14) ₂(x ₂, y ₂) and angle point estimated value m ₂' (x ₂', y ₂') Euclidean distance d, if d is less than threshold value d ', by angle point m ₁(x ₁, y ₁), m ₂(x ₂, y ₂) be set to intra-office point (being background light stream), if d is more than or equal to threshold value d ', by angle point m ₁(x ₁, y ₁), m ₂(x ₂, y ₂) being set to point not in the know, point not in the know is moving object light stream or mismatches quasi-optical stream.When completing after a RANSAC calculating, the current intra-office point of initialization is point not in the know, and repeating step 6.2 and 6.3, when no longer producing new intra-office point or reaching regulation multiplicity, finishes current RANSAC algorithm,

d＝sqrt((x ₂-x ₂′) ²+(y ₂-y ₂′) ²) (14)

Describedly coloured image J is carried out to initialized process be: first travel through all angle points of this frame, if this angle point overlaps with the correct background light stream that coloured image I retains, the initial intra-office point of calculating using this angle point as RANSAC.Repeating step 6.1,6.2,6.3, completes the estimation of RANSAC frame by frame.

According to above-mentioned steps, constantly by the moving object light stream in image or mismatch quasi-optical stream and weed out; Utilize the angle point light stream of RANSAC algorithm extraction, complete the registration of image.

Beneficial effect of the present invention can further illustrate by following experiment:

The main control device that this experiment adopts is the FPGA that the model of ALTERA company is EP2C5Q208C6, sensor comprises two parts: the inertial sensor that model is ADIS16405 and Kinect video camera, and wherein the gyroscope in inertial sensor and accelerometer obtain respectively rotational motion parameter and the translation motion parameter of video camera; Kinect video camera comprises that depth image obtains CCD (degree of depth CCD) and coloured image obtains CCD (colored CCD), and the two can be in conjunction with the depth information that obtains coloured image.Because inertial sensor is of a size of 35mm*30mm, so hardware platform volume is little, can popularize in mobile device.

Utilize Kinect collected by camera to obtain size to be coloured image and the depth image of 640*480.According to the real time analysis to PC end algorithm, the image acquisition frequency of hardware platform is set the present embodiment and the parameter acquisition frequency of inertial sensor is 40 frames/s.The window of Harris Corner Detection is chosen for 5*5, and the calculation times of RANSAC algorithm is set to 30 times.

For the effect of OFCM model of verifying that better the inventive method proposes, this experiment has verified that respectively Kinect video camera is under jump in brightness, random translation motion, Random-Rotation motion, random translation rotation compound movement condition, use side of the present invention to carry out image registration to Same Scene, obtain respectively the design sketch as shown in Fig. 2, Fig. 3, Fig. 4, Fig. 5.In Fig. 2, Fig. 3, Fig. 4, Fig. 5, the left-half of every width figure is that left figure, right half part are right figure, and figure interval, left and right is 10 frames, and wherein left figure is k two field picture, and right figure is k+10 two field picture.From Fig. 2, Fig. 3, Fig. 4, Fig. 5, find out, the inventive method can complete the registration of image accurately.

Table 1 is under above-mentioned different cameras moving condition, PC end image registration statistical analysis of performance.According to table 1, the indices of mode of motion of the same race is carried out to longitudinal comparison, when can find the inventive method in the generation translation motion of Kinect video camera and rotatablely moving, image registration has all obtained good result, and average treatment precision can reach 69.5%; The indices of the situation to Different Exercise Mode hypograph registration carries out lateral comparison, can find the inventive method to the registration accuracy of Kinect video camera single movement the registration accuracy higher than hybrid motion, this error is mainly to be caused by the error stack of gyroscope and accelerometer; Just process real-time, the about 22.92ms of average handling time on PC, temporal difference is mainly that the probability of Kinect video camera hybrid motion generation Mismatching point becomes large, therefore RANSAC algorithm is estimated relatively consuming time; According to the time Estimate of PC end, consider the actual difficulty that hardware gathers, the synchronous acquisition frequency that the Image Acquisition frame frequency of Kinect video camera and the parameter of inertial sensor are set is 30Hz.

Table 1PC end image registration performance statistics table

Claims (8)

1. the method for registering images based on inertial sensor and video camera, is characterized in that, comprises the following steps:

The coloured image I of Same Scene and the three-axis moving angular velocity [w of J, depth image I ' and J ' and video camera are obtained in step 1, synchro control video camera and inertial sensor work simultaneously _x, w _y, w _z] and three-axis moving acceleration [a _x, a _y, a _z];

The angle point m of step 2, extraction coloured image I ₁(x ₁, y ₁);

Step 3, utilize the three axis angular rate [ws of Kalman filtering to video camera _x, w _y, w _z] and 3-axis acceleration [a _x, a _y, a _z] carry out filtering, obtain Kalman Filter Estimation three axis angular rates

with Kalman Filter Estimation 3-axis acceleration ; Use euler dynamical equations formula will

convert Eulerian angle [α to _k, β _k, γ _k], use equation of uniformly accelerated motion will

convert three-shaft displacement [t to _xk, t _yk, t _zk], thereby obtain video camera, join rotation matrix R and the translation matrix T in matrix [R, T] outward;

Step 4, utilize the affine transformation relationship of coloured image and depth image, solve the angle point m of coloured image I ₁(x ₁, y ₁) coordinate position m in depth image I ' and J ' _1ir(x ₁, y ₁) and m _2ir(x ₁, y ₁); Utilize the proportionate relationship of depth image gray-scale value and pixel depth information to solve angle point m ₁(x ₁, y ₁) depth information u in depth image I ' and J ' ₁and u ₂;

The depth information u of the coloured image that step 5, the rotation matrix R that step 3 is obtained and translation matrix T and step 4 obtain ₁, u ₂bring OFCM model into, calculate the background light stream m that obtains image J ₂(x ₂, y ₂), OFCM model as shown in Equation (1),

$m_2=\frac{u_1}{u_2}{\mathrm{KRK}}^{-1}{m}_1+\frac{1}{u_2}\mathrm{KT}---\left(1\right)$

In formula (1), the internal reference matrix that K is video camera, K ^-1for the inverse matrix of K;

The background light stream m of step 6, use RANSAC color image I ₂carry out Robust estimation, by background light stream m ₂the light stream being produced by moving object of middle existence and mismatch quasi-optical stream and weed out, thereby complete the image registration of coloured image I and coloured image J, then utilize the correct background light stream of the coloured image I retaining to carry out initialization to coloured image J, carry out the registration of successive image.

2. the method for registering images based on inertial sensor and video camera as claimed in claim 1, it is characterized in that, in step 1, FPGA is by SPI protocol synchronization control video camera and inertial sensor, inertial sensor comprises gyroscope and accelerometer, video camera comprises that depth image obtains CCD and coloured image obtains CCD, and the minimum system of FPGA and inertial sensor is all fixed on depth image and obtains on CCD.

3. the method for registering images based on inertial sensor and video camera as claimed in claim 1, is characterized in that, uses Harris Corner Detection to extract the angle point m of coloured image I in step 2 ₁(x ₁, y ₁), computation process is as follows:

3.1 calculate pixel p[x in coloured image I, y] correlation matrix M, computing method as shown in formula (2),

$M=\underset{x,y}{\mathrm{\Σ}}w(x,y)\left[\begin{array}{cc}{I}_x^2& I_x{I}_y\\ I_x{I}_y& I_y^2\end{array}\right]=w(x,y)\⊗\left[\begin{array}{cc}{I}_x^2& I_x{I}_y\\ I_x{I}_y& I_y^2\end{array}\right]---\left(2\right)$

In formula (2), (x, y) is the image coordinate of pixel p, and w (x, y) is that window size is Gauss's weighting function of w*w, I _xrepresent the derivative of line direction, I _yrepresent the derivative of column direction;

3.2 calculating pixel point m ₁(x ₁, y ₁) angle point response R (x, y), computing method as shown in formula (3),

R(x,y)＝detM-k*(traceM) ² (3)

In formula (3), detM represents the value of correlation matrix M determinant, and traceM represents the mark of correlation matrix M, and k is constant, is generally 0.04-0.06;

3.3 judge pixel m ₁(x ₁, y ₁) whether be angle point, determination methods is shown in formula (4), if pixel m ₁(x ₁, y ₁) angle point response R (x, y) meet formula (4), pixel m ₁(x ₁, y ₁) be angle point; If pixel m ₁(x ₁, y ₁) angle point response R (x, y) do not meet formula (4), by m ₁(x ₁, y ₁) give up;

[R(x,y)≥R(x',y')]∩[R(x,y)≥0.01R _max] (3)

In formula (4), R (x', y') represents in window w*w except pixel m ₁(x ₁, y ₁) outside the angle point response of other pixel, R _maxfor edge feature frame I _fin the response of maximum angle point.

4. the method for registering images based on inertial sensor and video camera as claimed in claim 1, is characterized in that, in step 3, and described acquisition Kalman Filter Estimation three axis angular rates

with Kalman Filter Estimation 3-axis acceleration computation process be:

4.1 set up error model and the state transitions relational expression of inertial sensor according to the performance of inertial sensor and kinematic parameter, and binding time sequence analysis solves scale factor error M, the zero point drift B of inertial sensor _f, the state-transition matrix Ψ between white Gaussian noise U, kinematic parameter _{k, k-1}and white Gaussian noise V;

Error model as shown in Equation (4), shown in state transitions relational expression formula (5),

z' _k＝z _k(M+1)+B _f+U (4)

z _k＝Ψ _k,k-1z _k-1+V (5)

In formula (4) and formula (5), z represents three axis angular rate [w _x, w _y, w _z] and 3-axis acceleration [a _x, a _y, a _z] in a certain kinematic parameter, z ' _kfor k moment inertial sensor is measured the measured value of a certain kinematic parameter obtaining, z _kfor the actual value of k a certain kinematic parameter of moment;

4.2 results that step 4.1 is asked for are brought the kinematic parameter z of Kalman filtering formula to the k moment into _kcarry out filtering, obtain the estimated value of this kinematic parameter

kalman filtering formula as shown in formula (6),

$\left\{\begin{array}{c}{z}_{k,k-1}={\mathrm{\Ψ}}_{k,k-1}{z}_{k-1}\\ \widehat{z_k}=z_{k,k-1}+K_k\left[{z^{\′}}_k(M+1)z_{k,k-1}\right]\\ P_{k,k-1}={\mathrm{\Ψ}}_{k,k-1}{P}_{k-1}{{\mathrm{\Ψ}}_{k,k-1}}^T+Q_{k-1}\\ K_k=P_{k,k-1}{(M+1)}^T{P}_k^{-1}\\ P_k=[I-K_k(M+1)]P_{k,k-1}\end{array}\right.---\left(6\right)$

In formula (6), the noise variance that Q is system, K _kfor kalman gain, P _kfor state covariance, z _{k, k-1}for state discreet value.

5. the method for registering images based on inertial sensor and video camera as claimed in claim 1, is characterized in that, in step 3, the computation process of described acquisition video camera rotation matrix R and translation matrix T is:

5.1 use euler dynamical equations by filtered k moment three-axis moving angular velocity w _xk, w _yk, w _zkconvert the Eulerian angle α in European coordinate system to _k, β _k, γ _k, by the α of gained _k, β _k, γ _kbring in rotation matrix R, complete the calculating of colored CCD rotation matrix R;

Described euler dynamical equations as shown in formula (7),

$\left[\begin{array}{c}{\stackrel{\·}{\mathrm{\α}}}_k\\ {\stackrel{\·}{\mathrm{\β}}}_k\\ {\stackrel{\·}{\mathrm{\γ}}}_k\end{array}\right]=\left[\begin{array}{ccc}-\sin {\mathrm{\α}}_k/\tan {\mathrm{\β}}_k& {\cos \mathrm{\α}}_k/\tan {\mathrm{\β}}_k& 1\\ {\cos \mathrm{\α}}_k& \sin {\mathrm{\α}}_k& 0\\ \sin {\mathrm{\α}}_k/\sin {\mathrm{\β}}_k& -\cos {\mathrm{\α}}_k/\sin {\mathrm{\β}}_k& 0\end{array}\right]*\left[\begin{array}{c}{w}_{\mathrm{xk}}\\ w_{\mathrm{yk}}\\ w_{\mathrm{zk}}\end{array}\right]---\left(7\right)$

In formula (7),

for Eulerian angle first order derivative, the parameter α that the k-1 moment is tried to achieve _k-1, β _k-1, γ _k-1as the initial value of formula (7);

Described rotation matrix R is as shown in formula (8):

$R=\left[\begin{array}{ccc}\cos {\mathrm{\γ}}_k\cos {\mathrm{\β}}_k& \cos {\mathrm{\γ}}_k\sin {\mathrm{\β}}_k\sin {\mathrm{\α}}_k-\sin {\mathrm{\γ}}_k& \cos {\mathrm{\γ}}_k\sin {\mathrm{\β}}_k\cos {\mathrm{\α}}_k-\sin {\mathrm{\α}}_k\sin {\mathrm{\γ}}_k\\ \sin {\mathrm{\γ}}_k\cos {\mathrm{\β}}_k& \sin {\mathrm{\α}}_k\sin {\mathrm{\β}}_k\sin {\mathrm{\γ}}_k+\cos {\mathrm{\γ}}_k\cos {\mathrm{\α}}_k& \sin {\mathrm{\γ}}_k\sin {\mathrm{\β}}_k\cos {\mathrm{\α}}_k+\cos {\mathrm{\γ}}_k\sin {\mathrm{\α}}_k\\ -\sin {\mathrm{\β}}_k& \cos {\mathrm{\β}}_k\sin {\mathrm{\α}}_k& \cos {\mathrm{\β}}_k\cos {\mathrm{\α}}_k\end{array}\right]---\left(8\right)$

The speed displacement formula of 5.2 uses as shown in formula (9) is by 3-axis acceleration a _xk, a _yk, a _zkconvert three-shaft displacement t to _xk, t _yk, t _zk,

$\left\{\begin{array}{c}{t}_k=v_{k}T+\frac{1}{2}{a}_k{T}^2\\ v_k=v_{k-1}+a_{k-1}\×T\end{array}\right.----\left(9\right)$

In formula (9), t _krepresent three-shaft displacement t _xk, t _yk, t _zkin the displacement of any change in coordinate axis direction, the shooting cycle that T is video camera, a _k-1for the axial acceleration of k-1 moment respective coordinates, v _k-1and v _kbe respectively the initial velocity in k-1 moment and k moment, giving the initial velocity of fixed system initial time is 0.

6. the method for registering images based on inertial sensor and video camera as claimed in claim 1, it is characterized in that, in step 4, according to direct linear transformation's algorithm, calculate the affine transformation relationship of coloured image and depth image, affine transformation relationship is as shown in formula (10):

m _ir＝Hm (10)

In formula (10), m _irfor the pixel coordinate in depth image, m is the pixel coordinate in coloured image.

7. the method for registering images based on inertial sensor and video camera as claimed in claim 1, it is characterized in that, the rule-of-thumb relation proposing by reverse engineering in step 4 obtains the proportionate relationship of gray-scale value and pixel depth information, and the proportionate relationship of gray-scale value and pixel depth information is as shown in formula (11):

u＝l×m (11)

In formula (11), m (x, y) is the gray-scale value that pixel coordinate (x, y) is located, and u is that pixel coordinate (x, y) is located corresponding depth value, and l is scale-up factor.

8. the method for registering images based on inertial sensor and video camera as claimed in claim 1, is characterized in that, in step 6,

The background light stream m of described use RANSAC color image I ₂the computation process of carrying out Robust estimation is:

The current all angle point m of 8.1 initialization ₁(x ₁, y ₁), m ₂(x ₂, y ₂) be intra-office point;

8.2 utilize DLT algorithm to calculate affine transformation matrix H ₁, according to formula (12) in conjunction with angle point m ₁(x ₁, y ₁) estimate angle point estimated value m ₂' (x ₂', y ₂'):

m ₂′=H ₁m ₁ (12)

8.3 calculate angle point m according to formula (13) ₂(x ₂, y ₂) and angle point estimated value m ₂' (x ₂', y ₂') Euclidean distance d, if d is less than threshold value d ', by angle point m ₁(x ₁, y ₁), m ₂(x ₂, y ₂) be set to intra-office point (being background light stream), if d is more than or equal to threshold value d ', by angle point m ₁(x ₁, y ₁), m ₂(x ₂, y ₂) be set to point not in the know; When completing after a RANSAC calculating, the current intra-office point of initialization is point not in the know, and repeating step 8.2 and 8.3, when no longer producing new intra-office point or reaching regulation multiplicity, finishes current RANSAC algorithm,

d＝sqrt((x ₂-x ₂') ²+(y ₂-y ₂') ²) (13)

Describedly coloured image J is carried out to initialized process be: first travel through all angle points of this frame, if this angle point overlaps with the correct background light stream that coloured image I retains, the initial intra-office point calculating this angle point as RANSAC, method described in repeating step 8.1,8.2,8.3, completes the estimation of RANSAC frame by frame.

Images