CN104835159A - Digital image correction method for continuous variable-focal-length optical imaging system - Google Patents

Abstract

The invention discloses a digital image correction method for a continuous variable-focal-length optical imaging system. The digital image correction method comprises the steps of acquiring planar template images which are photographed by an image camera in a plurality of discrete focal lengths, performing off-line calibration on a distortion parameter which corresponds with each focal length of a variable-focal-length lens; performing curve fitting on the distortion parameter which corresponds with each focal length for obtaining a fitting formula or establishing a distortion parameter lookup table; according to the actual operation focal length of the image camera, obtaining a lens distortion parameter which corresponds with the actual operation focal length through the distortion parameter lookup table or calculation by the fitting formula; and constructing a projection conversion relationship from an image camera coordinate system to a map coordinate system according to position gesture data and the lens distortion parameter in imaging of the image camera, and performing resampling on the pixel brightness value after coordinate conversion for obtaining an orthographic projection image after correction on the squint deformation and lens distortion. The digital image correction method realizes a purpose of simultaneous correction for squint trapezoidal distortion and variable-focal-length nonlinear distortion in the variable-focal-length imaging system.

Description

For the method for correcting digital images of zoom optical imaging system

Technical field

The present invention relates to digital image processing field, especially, relate to a kind of digital picture distortion correction method for zoom optical imaging system.

Background technology

Zoom optical system has the advantage of not lose objects in the transfer process of visual field, can search for fast-moving target and catch, realize the function of Large visual angle search target and small field of view resolution target preferably, be widely used in aeroplane photography, Tracking and Measurment.Due to system manufacture and rigging error, the nonlinear distortion characteristic of zoom lens changes with the change of focal length, simultaneously due to the change in aspect and video camera optical axis orientation, aerial camera is with the photographic imagery over the ground of certain angle of inclination, create complicated geometry deformation, have a strong impact on the precision utilizing its image to carry out in real time splicing, also cause that the target photographed rotates, convergent-divergent, to be even out of shape, greatly increase and continual and steady difficulty of following the tracks of is carried out to target.

For the keystone distortion that stravismus photography causes, bearing calibration has the geometric correction method (multinomial model, rational function model, support vector machine etc.) based on empirical model, the method that corrects and the geometric correction method (collinearity equation method, projective transform method) based on conformation model is carried out with reference picture registration, wherein the latter is not owing to needing manually gather ground control point data or provide reference picture, be convenient to adopt embedded system automatically to realize, obtain in the geometry correction of remote sensing images and pay close attention to widely and apply.Current most Method of Remote Sensing Image Geometric Correction generally only corrects the geometry deformation that stravismus photography causes, less for itself distortion of long focal length lens, the correction accuracy needed for engineer applied can be reached, and zoom lens is larger in the lens distortion of short focus state, image border place aberration rate is up to 2% ~ 5%, and the nonlinear distortion characteristic of zoom lens changes with the change of focal length, therefore all must correct stravismus distortion and camera lens nonlinear distortion.

For the nonlinear distortion of imaging system, bearing calibration is mainly divided into two classes: a class is camera marking method, comprise Tsai based on the two-step approach of radial constraint and the Zhang Zhengyou scaling method based on plane target, these class methods are based on camera imaging model, consider in video camera simultaneously, outer parameter, solving precision is high, but calculated amount is large, in Optimizing Search process in distortion model parameter and video camera, outer parameter couples easily causes solution procedure do not restrain or converge on local minimum, and usually need to gather the coordinate that several target images extract a large amount of calibration point, because the nonlinear distortion of zoom lens changes with the change of focal length, for the particular focal length state of system, general needs re-start demarcation, complex operation is consuming time, another kind of is based on projective geometry unchangeability or invariant (geometrical perspective Projective invariance, cross ratio invariability principle, line segment slope, vanishing point and plane restriction etc.) non-measurement bearing calibration, point Main Basis photographic subjects with common trait still has the character of this feature in ideal image, true picture is found and demarcates the point meeting special characteristic, foundation take distortion parameter as the linear functional of optimum solution, the method of Optimizing Search is adopted to solve distortion parameter, this kind of bearing calibration will obtain the structural information of scene in advance, to under circumstances not known condition (as scout, the application such as information) image that obtains not necessarily is suitable for, such as, in image, edge feature is not significantly (as grassland, desert, ocean etc.) or when comprising curve in curve and some real worlds that lens distortion causes simultaneously, utilize the distortion parameter estimated result that the optimizing search method of geometrical perspective Projective invariance possibly cannot be implemented or lead to errors.

There is no the method that simultaneously can correct the nonlinear distortion of imaging system and stravismus keystone distortion at present, therefore existing bearing calibration is difficult to the imaging demand meeting zoom optical system.

Summary of the invention

The invention provides a kind of method for correcting digital images for zoom optical imaging system, to solve the technical matters that existing method for correcting image cannot correct the nonlinear distortion of imaging system and stravismus keystone distortion simultaneously.

The technical solution used in the present invention is as follows:

For a method for correcting digital images for zoom optical imaging system, comprising:

Acquisition camera is at the plane template image of some discrete focal length value shootings, and off-line calibration goes out distortion parameter corresponding to each focal length of zoom lens;

The distortion parameter corresponding to each focal length carries out curve fitting and obtains fitting formula or set up distortion parameter look-up table;

According to the real work focal length value of video camera, searched by distortion parameter look-up table and to obtain or fitting formula calculates lens distortion parameter corresponding to real work focal length value;

The Transformation Relation of Projection of map coordinates system is tied to according to position and attitude data during video camera imaging and lens distortion parametric configuration camera coordinates, by photographed images re-projection in map coordinates system, resampling is carried out to the pixel brightness value after coordinate transform and obtains the orthographic projection images after correcting stravismus distortion and lens distortion.

Further, distortion parameter comprises: distortion factor k ₁with center of distortion coordinate (u ₀, v ₀).

Further, some discrete focal length values according to the corresponding selection of the focal-distance tuning range of zoom lens, the corresponding at least one secondary plane template image of each discrete focal length value,

The distortion parameter that each focal length is corresponding adopts distortion model parameter estimation algorithm to solve and obtains.

Further, distortion model parameter estimation algorithm solves distortion parameter and comprises the following steps:

Determine the span of corresponding distortion parameter according to the distortion performance of general zoom lens, comprise distortion factor k ₁span and center of distortion coordinate (u ₀, v ₀) span;

The edge pixel point of detection plane template image, obtains corresponding edge image;

According to the span of distortion parameter, corresponding step-length is selected to obtain distortion parameter combination;

To the correcting image often organizing distortion parameter edge calculation image, and calculate correction back edge pixel gradient, obtain coordinate and the gradient direction of each edge pixel point;

Calculate ballot to each edge pixel point in correcting image, the distortion parameter of trying to achieve ballot sum maximal value corresponding is optimum value.

Further, gone out the distortion parameter under some discrete focal length values by described distortion model parameter estimation algorithm off-line calibration, the distortion parameter corresponding to each focal length carries out curve fitting and obtains in fitting formula, distortion factor k ₁comprise with the funtcional relationship of focal distance f:

Work as f ₁≤ f≤f ₂time,

k ₁(f)＝ρ·f ²+σ·f+τ，

Work as f ₂≤ f≤f ₃time,

$k_1\left(f\right)=\frac{\mathrm{\δ}}{{(f+\mathrm{\η})}^2}+\mathrm{\ψ},$

Parameter ρ, σ, τ, δ, η, ψ, f in formula ₁, f ₂, f ₃value different with different camera lens, the method combined with curve by off-line calibration is obtained.

Further, gone out the distortion parameter under some discrete focal length values by described distortion model parameter estimation algorithm off-line calibration, the distortion parameter corresponding to each focal length carries out curve fitting and obtains in fitting formula, center of distortion coordinate (u ₀, v ₀) with the pass of focal distance f be:

u ₀(f)＝λ ₁·f+μ ₁，v ₀(f)＝λ ₂·f+μ ₂,

Parameter lambda in formula ₁, μ ₁, λ ₂, μ ₂value different with different camera lens, the method combined with curve by off-line calibration is obtained.

Further, the real work focal length value of video camera is detected by the position transducer on the focus adjusting mechanism on video camera and obtains.

Further, the information that position and attitude data during video camera imaging are detected by the angular transducer of GPS locator data and video camera inside calculates.

Further, video camera is aviation zoom camera.

The present invention has following beneficial effect:

The present invention is used for the method for correcting digital images of zoom optical imaging system, adopt the method that on ground off-line calibration and machine, on-line correction combines to correct stravismus keystone distortion and zoom lens's nonlinear distortion simultaneously, avoid in the correcting distorted process of general camera marking method the shortcoming needing to repeat each focal length state loaded down with trivial details camera calibration consuming time, break away from again non-measurement bearing calibration to image cathetus, the dependence of the geometric properties such as vanishing point, decrease the number of processes to pel data, define the Processing Algorithm settled at one go, at correction accuracy, correction rate and automaticity are significantly improved.The present invention adopts the mode that on ground off-line calibration and machine, on-line correction combines, do not need carry out registration with reference picture and gather ground control point data, greatly reduce correction time, the orthogonal projection position calculating each pixel has concurrency, field programmable gate array (field programming gate arrays can be mapped to easily, FPGA) in, meet the processing demands of Real-time System by the advantage of hardware parallel computation in speed, thus realize automatic on-line real time correction on machine.

Except object described above, feature and advantage, the present invention also has other object, feature and advantage.Below with reference to figure, the present invention is further detailed explanation.

Accompanying drawing explanation

The accompanying drawing forming a application's part is used to provide a further understanding of the present invention, and schematic description and description of the present invention, for explaining the present invention, does not form inappropriate limitation of the present invention.In the accompanying drawings:

Fig. 1 is preferred embodiment of the present invention aviation zoom camera geometric correction of imagery process flow diagram;

Fig. 2 is that preferred embodiment of the present invention aviation zoom camera is as geometry correction schematic diagram;

Fig. 3 is distortion coefficients of camera lens k ₁with the matched curve of focal distance f;

Fig. 4 is center of distortion coordinate (u ₀, v ₀) with the matched curve of focal distance f.

Description of reference numerals:

G: video camera projection centre;

F _c: camera coordinate system;

F _b: carrier aircraft coordinate system;

F _v: carrier aircraft geographic coordinate system;

M: map coordinates system;

A, b, c, d: 4 angle points of original image;

O: principal point or center of distortion;

U, v: two coordinate axis of original image pixels coordinate system, the row of u axle mark original image, the row of v axle mark original image;

X, y: two coordinate axis of original image physical coordinates system, x-axis is parallel with u axle, and y-axis is parallel with v axle;

T: object point t ＇ actual image point on the image plane;

T ₀: object point t ＇ ideal image point on the image plane;

A ', b ', c ', d ': the subpoint of original image 4 angle points in map coordinates system;

T ＇: a certain object point;

O ＇: the principal point o subpoint in map coordinates system;

U ＇, v ＇: two coordinate axis of correcting image pixel coordinate system, the row of u ＇ axle mark correcting image, the row of v ＇ axle mark correcting image;

O _m: map coordinates system initial point;

X _m, y _m, z _m: three coordinate axis of map coordinates system.

Embodiment

Below in conjunction with accompanying drawing, embodiments of the invention are described in detail, but the multitude of different ways that the present invention can be defined by the claims and cover is implemented.

The preferred embodiments of the present invention provide a kind of method for correcting digital images for zoom optical imaging system, nonlinear distortion that zoom lens causes and stravismus can be corrected to photograph the stravismus keystone distortion caused simultaneously, the inventive method adopts the mode that on ground off-line calibration and machine, on-line correction combines, according to geometrical perspective Projective invariance principle, utilize one-parameter to remove formula model calibrates some discrete focal length value lower plane template images distortion factor and center of distortion coordinate by variable step-size search method, set up distortion coefficients of camera lens, Mathematical Fitting formula between center of distortion coordinate and focal length or construct distortion parameter look-up table, according to this fitting formula or the lens distortion parameter obtained under arbitrary continuation focal length of tabling look-up, introduce aircraft position, attitude and the video camera optical axis point to the factors such as orientation, structure camera coordinates is tied to the Transformation Relation of Projection of map coordinates system, by aerial image re-projection in map coordinates system, resampling is carried out to the pixel brightness value after coordinate transform and obtains the orthographic projection images after correcting stravismus keystone distortion and nonlinear distortion.This method avoid in the correcting distorted process of general camera marking method the shortcoming needing to repeat each focal length state loaded down with trivial details camera calibration consuming time, break away from again the dependence of non-measured bearing calibration to geometric properties such as image cathetus, vanishing points, decrease the number of processes to pel data, define the Processing Algorithm settled at one go, correction accuracy, correction rate and automaticity are significantly improved.

The preferred embodiment of the present invention is imaged as example with aviation zoom camera, and adopt operating system to be WINDOWS XP, software platform is MATLAB 2012b, and processor is second generation Duo i5-4210M.

As shown in Figure 1, the aviation zoom camera image deformation bearing calibration of the specific embodiment of the invention, comprises the following steps:

Step 1: the plane template image that acquisition camera is taken under some discrete focal length values, goes out distortion factor and the center of distortion coordinate of zoom lens at ground off-line calibration.

First the plane template image taken under different focal of acquisition camera, image size is 1024 × 768pixels, adopt the plane checkerboard pattern of 35 × 35 grids as plane template, pattern magnitude is 350mm × 350mm, focal length of camera scope is from 5mm to 100mm, select 13 discrete focal length values: 5.8, 10, 15, 20, 25, 30, 40, 50, 60, 70, 80, 90, 100mm, respectively a width template image is taken to each focal length value, distortion model parameter estimation algorithm is below adopted to solve distortion coefficients of camera lens corresponding to this focal length and center of distortion coordinate, double counting is averaged as parameter estimation result 20 times, wherein distortion model parameter estimation algorithm comprises the following steps:

Step 1-1: determine distortion factor k ₁with center of distortion coordinate (u ₀, v ₀) scope, if image size is w × h pixels, D is image diagonal length, distortion factor k ₁generally be positioned at [-1/D ², 1/D ²] in scope, center of distortion is generally positioned at the rectangular area of 0.1w × 0.1h size near picture centre, obtains u ₀∈ [0.45w, 0.55w], v ₀∈ [0.45h, 0.55h].

Step 1-2: adopt the edge pixel point in Canny edge detector detection original distortion image, obtain corresponding edge image, wherein the threshold value of Canny operator is determined according to the certain percentage of image gradient value.

Step 1-3: according to above-mentioned distortion parameter k ₁, u ₀, v ₀span, select corresponding step-length δ k respectively ₁, δ u ₀, δ v ₀, obtain distortion parameter combination (k ₁ ⁱ, u ₀ ^j, v ₀ ^t) as follows: $k_1^i=-1/D^2+i\·{\mathrm{\δk}}_1,u_0^j=0.45w+j\·{\mathrm{\δu}}_0,v_0^t=0.45h+t\·{\mathrm{\δv}}_0,$ Wherein i=1,2 ..., N ₁, j=1,2 ..., N ₂, t=1,2 ..., N ₃, δ k ₁=2/ (N ₁d ²), δ u ₀=0.1w/N ₂, δ v ₀=0.1h/N ₃.

Step 1-4: adopt above-mentioned various distortion parameter combination (k according to formula (1) ~ (2) ₁ ⁱ, u ₀ ^j, v ₀ ^t) correcting image of edge calculation image, and according to formula (3) ~ (4) calculate correct after the gradient of each edge pixel point, obtain the coordinate (u of each edge pixel point after correcting _n, v _n) and gradient direction α (u _n, v _n) as follows.

x _d＝(u _d-u ₀)d _x，y _d＝(v _d-v ₀)d _y， (1)

$\left.\begin{array}{c}{u}_n=u_0+x_d/(1+k_1{x_d}^2+k_1{y_d}^2)/d_x\\ v_n=v_0+y_d/(1+k_1{x_d}^2+k_1{y_d}^2)/d_y\end{array}\right\},---\left(2\right)$

$\left.\begin{array}{c}{G}_u=(I_{v_n,u_n+1}-I_{v_n,u_n}+I_{v_n+1,u_n+1}-I_{v_n+1,u_n})/2\\ G_v=(I_{v_n+1,u_n}-I_{v_n,u_n}+I_{v_n+1,u_n+1}-I_{v_n,u_n+1})/2\end{array}\right\},---\left(3\right)$

α(u _n,v _n)＝arctan(G _v/G _u). (4)

D in formula _x, d _yrepresent the physical size of single pixel, unit is μm, (u _d, v _d), (x _d, y _d) be respectively the distortion pixel coordinate of picture point and physical coordinates, (u _n, v _n) for correcting the pixel coordinate of after image point, I represents image brightness values, G _u, G _vfor the edge image brightness value after correction is at (u _n, v _n) first-order partial derivative at place.

Step 1-5: the HOUGH conversion calculating correcting image, tries to achieve the Nl bar straight-line segment edge that top n HOUGH converter unit peak value is corresponding, and with the distance dist (q) of initial point and direction β (q), q=1,2 ..., N.

Step 1-6: ballot is calculated to each edge pixel point in correcting image: if the gradient direction α (u of this pixel _n, v _n) to differ with direction β (q) of q article of straight line and be less than a certain threshold value δ _α, such as δ _α=2 °, calculate the distance of this point and q article of straight line:

d _q＝|u _ncos(β(q))+v _nsin(β(q))-dist(q)|， (5)

If d _qbe less than a certain threshold value δ _d, such as δ _d=2pixels, calculates the ballot value votes=1/ (1+d of this point _q), calculate the ballot sum of all edge pixels point, try to achieve the distortion parameter k that ballot sum maximal value is corresponding ₁ ⁽⁰⁾, u ₀ ⁽⁰⁾, v ₀ ⁽⁰⁾for optimum value:

$\max \left\{\underset{q=1}{\overset{N}{\mathrm{\Σ}}}\mathrm{votes}(\mathrm{dist}\left(q\right),\mathrm{\β}\left(q\right),k_1^i,u_0^j,v_0^t)\right\}.---\left(6\right)$

In formula for adopting distortion parameter (k ₁ ⁱ, u ₀ ^j, v ₀ ^t) in correcting image each edge pixel point to the ballot value of Nl bar straight line.

In order to distortion factor k ₁with center of distortion u ₀, v ₀estimate more accurately, preferably, also comprise step 1-7: respectively at [k ₁ ⁽⁰⁾-δ k ₁, k ₁ ⁽⁰⁾+ δ k ₁], [u ₀ ⁽⁰⁾-δ u ₀, u ₀ ⁽⁰⁾+ δ u ₀], [v ₀ ⁽⁰⁾-δ v ₀, v ₀ ⁽⁰⁾+ δ v ₀] in scope, select step-size in search to be 1/N in step 1-3 ₁, 1/N ₂, 1/N ₃, repeat step 1-2 to step 1-7, until k ₁hunting zone be less than 10 ^-10, algorithm terminates, now corresponding coefficient k ₁, u ₀, v ₀for the optimum value after optimization.

Step 2: the distortion parameter corresponding to each focal length carries out curve fitting or set up distortion parameter look-up table.

The distortion parameter corresponding to each focal length carries out curve fitting, and as Fig. 3, shown in 4, can find distortion factor k from Fig. 3 ₁following funtcional relationship is there is with f:

As 5.8≤f≤20mm

k ₁(f)＝ρ·f ²+σ·f+τ， (7)

ρ=1.369 × 10 in formula ^-10, σ=-1 × 10 ^-9, τ=-1.108 × 10 ^-7

As 20≤f≤100mm

$k_1\left(f\right)=\frac{\mathrm{\δ}}{{(f+\mathrm{\η})}^2}+\mathrm{\ψ},---\left(8\right)$

δ=-5.245 × 10 in formula ^-5, η=2.768, ψ=2.564 × 10 ^-8.

As can be seen from Figure 4, when focal length variations, center of distortion is located substantially on straight line, the pixel coordinate u of center of distortion ₀, v ₀can be similar to the relation of focal distance f and represent with following formula:

u ₀(f)＝λ ₁·f+μ ₁，v ₀(f)＝λ ₂·f+μ ₂. (9)

λ in formula ₁=0.5722, μ ₁=500.4904, λ ₂=0.2897, μ ₂=363.5341.

Step 3: during actual on-line correction, video camera real work focal length value is obtained according to respective sensor measurement on video camera focus adjusting mechanism, lens distortion parameter corresponding to this focal length value is calculated by look-up method or according to the fitting formula between distortion parameter with focal length, elements of exterior orientation is calculated again according to position and attitude data during video camera imaging, structure camera coordinates is tied to the Transformation Relation of Projection of map coordinates system, by aerial image re-projection in map coordinates system, resampling is carried out to the pixel brightness value after coordinate transform and obtains the orthographic projection images after correcting stravismus distortion and lens distortion.

The center of distortion coordinate figure structure look-up table utilizing different focal corresponding, actual timing obtains the real work focal length of video camera by sensor measurement, obtaining the center of distortion coordinate that adjacent two focal length values of this focal length value are corresponding by tabling look-up, adopting linear interpolation method to calculate center of distortion coordinate (u corresponding to this focal length value ₀, v ₀), utilize fitting formula (7) ~ (8) to calculate distortion factor value k corresponding to this focal length value ₁.

Elements of exterior orientation when calculating video camera imaging according to the camera position of the angular transducer record of boat appearance measuring system, GPS and video camera inside and attitude information, structure camera coordinates is tied to the Transformation Relation of Projection of map coordinates system: as shown in Figure 2, when not considering lens distortion, object point t ＇, ideal image point t ₀with video camera projection centre G point-blank, meet pin-hole imaging model; Due to the impact of lens distortion, object point in the corresponding picture point of the plane of delineation from ideal position t ₀move to distortional point position t, the map reference (x of object point t ＇ _m, y _m) to the pixel coordinate (u of corresponding picture point t _d, v _d) transformational relation derive as follows.

The coordinate system situation that the present embodiment relates to is described as follows.

Camera coordinate system (F _c): initial point is video camera projection centre G, x _caxle, y _caxle respectively with image pixel coordinates system u axle (row of marking image, unit is pixel), v axle (row of marking image, unit is pixel) is parallel and direction is consistent; The initial point of image physical coordinates system o-xy is positioned at intersection point and the principle point location of camera optical axis and the plane of delineation, x-axis, y-axis respectively with u axle, v axle is parallel and direction is consistent, this coordinate system is in units of m or mm.

Carrier aircraft coordinate system (F _b): initial point is boat appearance measuring system barycenter, and general boat appearance measuring system is installed on the level reference of photoelectric platform, and boat appearance measuring system barycenter and video camera projection centre, apart from very little, can be similar to and think that both overlap, x _baxle is 0 ° of direction of boat appearance measuring system, y _baxle is 90 ° of directions of boat appearance measuring system, z _baxle is determined by right-hand screw rule, distance lambda between the video camera that the position angle Θ that photoelectric platform interior angle scrambler exports and angular altitude Ψ and range finder using laser export and field of view center target ₁it is this coordinate system relative.

Carrier aircraft geographic coordinate system (F _v): initial point is positioned at boat appearance measuring system barycenter, is NED (North East Down) coordinate system, the carrier aircraft course angle β that boat appearance measuring system exports, and angle of pitch ε and roll angle γ are this coordinate systems relative.

, only there is a translational movement with carrier aircraft geographic coordinate system in map coordinates system (m): be also NED coordinate system, if video camera projection centre G is at the subpoint o of ground level during shooting the 1st width image _mfor map coordinates system initial point, translational movement when taking the i-th width image between Two coordinate system is the coordinate [x of carrier aircraft position in map coordinates system that this moment GPS exports _miy _miz _mi] ^t.

Principle without loss of generality, supposes that the orthographic projection images of original image is positioned at map coordinates system z _m=0 plane, calculates object point t ＇ (x _m, y _m) corresponding ideal picture point t ₀physical coordinates (x _n, y _n) as follows:

$s\left[\begin{array}{c}{x}_n\\ y_n\\ 1\end{array}\right]=C_{\mathrm{cm}}\left[\begin{array}{c}{x}_m\\ y_m\\ 0\\ 1\end{array}\right]=\left[\begin{array}{cccc}{m}_{11}& m_{12}& m_{13}& m_{14}\\ m_{21}& m_{22}& m_{23}& m_{24}\\ m_{31}& m_{32}& m_{33}& m_{34}\end{array}\right]\left[\begin{array}{c}{x}_m\\ y_m\\ 0\\ 1\end{array}\right]=\left[\begin{array}{ccc}{m}_{11}& m_{12}& m_{14}\\ m_{21}& m_{22}& m_{24}\\ m_{31}& m_{32}& m_{34}\end{array}\right]\left[\begin{array}{c}{x}_m\\ y_m\\ 1\end{array}\right],---\left(10\right)$

C in formula _cmfor video camera elements of exterior orientation, represent the transformational relation being tied to camera coordinate system from map reference, through the conversion of map coordinates system, carrier aircraft geographic coordinate system, carrier aircraft coordinate system, camera coordinate system four coordinate systems, need represent with following formula:

$C_{\mathrm{cm}}=\left[\begin{array}{cccc}{c}_{\mathrm{\Θ}}{s}_{\mathrm{\Ψ}}& s_{\mathrm{\Ψ}}{s}_{\mathrm{\Θ}}& -c_{\mathrm{\Ψ}}& 0\\ -s_{\mathrm{\Θ}}& c_{\mathrm{\Θ}}& 0& 0\\ c_{\mathrm{\Θ}}{c}_{\mathrm{\Ψ}}& s_{\mathrm{\Θ}}{c}_{\mathrm{\Ψ}}& s_{\mathrm{\Ψ}}& 0\end{array}\right]\left[\begin{array}{cccc}{c}_{\mathrm{\ϵ}}{c}_{\mathrm{\β}}& c_{\mathrm{\ϵ}}{s}_{\mathrm{\β}}& -s_{\mathrm{\ϵ}}& 0\\ c_{\mathrm{\β}}{s}_{\mathrm{\γ}}{s}_{\mathrm{\ϵ}}-c_{\mathrm{\γ}}{s}_{\mathrm{\β}}& c_{\mathrm{\γ}}{c}_{\mathrm{\β}}+s_{\mathrm{\γ}}{s}_{\mathrm{\ϵ}}{s}_{\mathrm{\β}}& c_{\mathrm{\ϵ}}{s}_{\mathrm{\γ}}& 0\\ s_{\mathrm{\γ}}{s}_{\mathrm{\β}}+c_{\mathrm{\γ}}{c}_{\mathrm{\β}}{s}_{\mathrm{\ϵ}}& c_{\mathrm{\γ}}{s}_{\mathrm{\ϵ}}{s}_{\mathrm{\β}}-c_{\mathrm{\β}}{s}_{\mathrm{\γ}}& c_{\mathrm{\γ}}{c}_{\mathrm{\ϵ}}& 0\\ 0& 0& 0& 1\end{array}\right]\left[\begin{array}{cccc}1& 0& 0& -x_{\mathrm{mi}}\\ 0& 1& 0& -y_{\mathrm{mi}}\\ 0& 0& 1& -z_{\mathrm{mi}}\\ 0& 0& 0& 1\end{array}\right],---\left(11\right)$

C=cos (*) in formula, s=sin (*), Θ, Ψ is respectively camera coordinate system relative to the position angle of carrier aircraft coordinate system and angular altitude, ε, β, γ is respectively the angle of pitch of the relative carrier aircraft geographic coordinate system of carrier aircraft coordinate system (NED coordinate system), course angle, roll angle, [x _miy _miz _mi] ^tfor the coordinate of carrier aircraft position in map coordinates system, the GPS location by carrier aircraft calculates.

First obtained longitude L, latitude M and the geodetic height H of carrier aircraft position by GPS measurement, calculate the coordinate of carrier aircraft at WGS-84 (WorldGeodetic System 1984) the earth's core rectangular coordinate system in space:

$P_e=\left[\begin{array}{c}{x}_e\\ y_e\\ z_e\end{array}\right]=\left[\begin{array}{c}(N+H)\cos M\cos L\\ (N+H)\cos M\sin L\\ [N(1-e^2)+H]\sin M\end{array}\right],---\left(12\right)$

In formula: semimajor axis of ellipsoid a=6378137.0m, semiminor axis of ellipsoid b=6356752.0m, ellipsoid first excentricity $e=\sqrt{(a^2-b^2)/a^2},$ Ellipsoid radius of curvature in prime vertical $N=a/\sqrt{1-e^2{\sin }^{2}M}.$

The coordinate P of carrier aircraft position in map coordinates system _m=[x _miy _miz _mi] ^tcan be calculated as follows:

P _m＝[x _miy _miz _mi] ^T＝R _m/e(P _e-P _e,ref)， (13)

P in formula _{e, ref}for map coordinates system initial point o _mcoordinate in the rectangular coordinate system in space of the earth's core, is calculated as follows:

$P_{e,\mathrm{ref}}=\left[\begin{array}{c}{x}_{e1}\\ y_{e1}\\ z_{e1}\end{array}\right]=\left[\begin{array}{c}(N_1+H_1)\cos M_1\cos L_1\\ (N_1+H_1)\cos M_1\sin L_1\\ [N_1(1-e^2)+H_1]\sin M_1\end{array}\right],---\left(14\right)$

L in formula ₁, M ₁, H ₁for map coordinates system initial point o _mlongitude, latitude and geodetic height, ellipsoid radius of curvature in prime vertical $N_1=a/\sqrt{1-e^2{\sin }^2{M}_1}.$

R in formula (13) _m/efor from the earth's core rectangular coordinate system in space to the rotation matrix of map coordinates system, be calculated as follows:

$R_{m/e}=\left[\begin{array}{ccc}-\cos L_1\sin M_1& -\sin L_1\sin M_1& \cos M_1\\ -\sin L_1& \cos L_1& 0\\ -\cos L_1\cos M_1& -\cos M_1\sin L_1& -\sin M_1\end{array}\right],---\left(15\right)$

General lens distortion is based on radial distortion form, and tangential distortion and the thin prism very I that distorts is ignored, and adopts one-parameter to describe zoom lens's distortion, ideal image point t except formula model here ₀following relation is there is with the physical coordinates of distortion picture point t:

x _n＝x _d/(1+k ₁r _d ²)，y _n＝y _d/(1+k ₁r _d ²)， (16)

K in formula ₁for coefficient of radial distortion, obtain by above-mentioned lens distortion parameter estimation algorithm, (x _d, y _d) for distortion picture point t physical coordinates, for it is to the Euclidean distance of center of distortion, order formula (16) can be converted into r _dquadratic equation with one unknown:

$r_d^2-r_d/\left(k_1{r}_n\right)+1/k_1=0,---\left(17\right)$

Separate above-mentioned system of equations and obtain r _dshown in (18), wherein k ₁corresponding pincushion distortion during >0, k ₁corresponding barrel distortion during <0

Solve r _dafter, calculate object point t ＇ (x according to following formula _m, y _m) pixel coordinate (u of corresponding picture point t in fault image _d, v _d) as follows

[u _dv _d1] ^T＝A[(r _d/r _n)x _n(r _d/r _n)y _n1] ^T. (19)

In like manner, the pixel coordinate (u of distortion picture point t _d, v _d) to the map reference (x of corresponding object point t ＇ _m, y _m) conversion as follows: first by the pixel coordinate (u of picture point t _d, v _d) calculate its physical coordinates (x _d, y _d), be shown below

[x _dy _d1] ^T＝A ^-1[u _dv _d1] ^T， (20)

Then object point t ＇ corresponding picture point t in desirable orthoscopic image is calculated according to formula (16) ₀physical coordinates (x _n, y _n), the map reference (x of object point t ＇ is finally calculated according to following formula _m, y _m)

$\left[\begin{array}{c}{x}_m\\ y_m\\ 1\end{array}\right]=s{\left[\begin{array}{ccc}{m}_{11}& m_{12}& m_{14}\\ m_{21}& m_{22}& m_{24}\\ m_{31}& m_{32}& m_{34}\end{array}\right]}^{-1}\left[\begin{array}{c}{x}_n\\ y_n\\ 1\end{array}\right].---\left(21\right)$

The Transformation Relation of Projection of map coordinates system is tied to, by original image 4 angular coordinate a (u according to camera coordinates ₂₁, v ₂₁), b (u ₂₂, v ₂₂), c (u ₂₃, v ₂₃), d (u ₂₄, v ₂₄) by formula (20), (16), (21) project in earth-fixed co-ordinate system and go, and obtain 8 coordinate figure a ' (x _m1, y _m1), b ' (x _m2, y _m2), c ' (x _m3, y _m3), d ' (x _m4, y _m4), as shown in Figure 2, and x is pressed to these 8 coordinate figures _mand y _mtwo set of coordinates ask its minimum value (x respectively _min, y _min) and maximal value (x _max, y _max).

Divide ground grid, the floor measurements d of each pixel of definition output image _xmand d _ym, the total line number row and the total columns col that obtain output calibration image are respectively

$\left.\begin{array}{c}\mathrm{row}=\mathrm{fix}[(y_{\max }-y_{\min })/d_{\mathrm{ym}}]+1\\ \mathrm{col}=\mathrm{fix}[(x_{\max }-x_{\min })/d_{\mathrm{xm}}]+1\end{array}\right\},---\left(22\right)$

In formula, fix represents bracket function.

For each pixel coordinate (u ', v ') in correcting image, calculate its map reference (x _m, y _m)

x _m＝x _min+(u'-1)d _xm，y _m＝y _min+(v'-1)d _ym. (23)

Wherein u '=1,2 ..., col, v '=1,2 ..., row.

According to formula (10), (18), (19), by map reference (x _m, y _m) calculate corresponding original image pixels point coordinate (u _d, v _d).

Bilinear interpolation algorithm is adopted to calculate original image pixels point coordinate (u _d, v _d) gray-scale value, by this gray-scale value assignment to correcting image (u ', v ') pixel, namely

In formula, q represents the gray matrix of correcting image, and p represents the gray matrix of original image, and h, w are respectively the ranks coordinate range of original image.

According to above-mentioned steps, calculate the gray-scale value of each picture point in correcting image successively, until all picture points calculate complete, obtain the correcting image that a width eliminates zoom lens's distortion and stravismus keystone distortion.

The foregoing is only the preferred embodiments of the present invention, be not limited to the present invention, for a person skilled in the art, the present invention can have various modifications and variations.Within the spirit and principles in the present invention all, any amendment done, equivalent replacement, improvement etc., all should be included within protection scope of the present invention.

Claims (9)

1. for a method for correcting digital images for zoom optical imaging system, it is characterized in that, comprising:

Acquisition camera is at the plane template image of some discrete focal length value shootings, and off-line calibration goes out distortion parameter corresponding to each focal length of zoom lens;

The distortion parameter corresponding to each focal length carries out curve fitting and obtains fitting formula or set up distortion parameter look-up table;

According to the real work focal length value of described video camera, searched by described distortion parameter look-up table and to obtain or described fitting formula calculates lens distortion parameter corresponding to described real work focal length value;

The Transformation Relation of Projection of map coordinates system is tied to according to position and attitude data during video camera imaging and described lens distortion parametric configuration camera coordinates, by photographed images re-projection in map coordinates system, resampling is carried out to the pixel brightness value after coordinate transform and obtains the orthographic projection images after correcting stravismus distortion and lens distortion.

2. the method for correcting digital images for zoom optical imaging system according to claim 1, is characterized in that,

Described distortion parameter comprises: distortion factor k ₁with center of distortion coordinate (u ₀, v ₀).

3. the method for correcting digital images for zoom optical imaging system according to claim 2, is characterized in that,

Described some discrete focal length values according to the corresponding selection of the focal-distance tuning range of described zoom lens, the described plane template image of the corresponding at least one pair of each described discrete focal length value,

The distortion parameter that described each focal length is corresponding adopts distortion model parameter estimation algorithm to solve and obtains.

4. the method for correcting digital images for zoom optical imaging system according to claim 3, is characterized in that,

Described distortion model parameter estimation algorithm solves described distortion parameter and comprises the following steps:

Determine the span of corresponding described distortion parameter according to the distortion performance of general zoom lens, comprise distortion factor k ₁span and center of distortion coordinate (u ₀, v ₀) span;

Detect the edge pixel point of described plane template image, obtain corresponding edge image;

According to the span of described distortion parameter, corresponding step-length is selected to obtain distortion parameter combination;

To the correcting image often organizing distortion parameter edge calculation image, and calculate correction back edge pixel gradient, obtain coordinate and the gradient direction of each edge pixel point;

Calculate ballot to each edge pixel point in described correcting image, the distortion parameter of trying to achieve ballot sum maximal value corresponding is optimum value.

5. the method for correcting digital images for zoom optical imaging system according to claim 4, is characterized in that,

Gone out the distortion parameter under some discrete focal length values by described distortion model parameter estimation algorithm off-line calibration, the distortion factor corresponding to each focal length carries out curve fitting and obtains in fitting formula, described distortion factor k ₁comprise with the funtcional relationship of focal distance f:

Work as f ₁≤ f≤f ₂time,

k ₁(f)＝ρ·f ²+σ·f+τ，

Work as f ₂≤ f≤f ₃time,

$k_1\left(f\right)=\frac{\mathrm{\&delta;}}{{(f+\mathrm{\&eta;})}^2}+\mathrm{\&psi;},$

Parameter ρ, σ, τ, δ, η, ψ, f in formula ₁, f ₂, f ₃value different with different camera lens, the method combined with curve by off-line calibration is obtained.

6. the method for correcting digital images for zoom optical imaging system according to claim 4, is characterized in that,

Gone out the distortion parameter under some discrete focal length values by described distortion model parameter estimation algorithm off-line calibration, the center of distortion coordinate corresponding to each focal length carries out curve fitting and obtains in fitting formula, center of distortion coordinate (u ₀, v ₀) with the pass of focal distance f be:

u ₀(f)＝λ ₁·f+μ ₁，v ₀(f)＝λ ₂·f+μ ₂,

Parameter lambda in formula ₁, μ ₁, λ ₂, μ ₂value different with different camera lens, the method combined with curve by off-line calibration is obtained.

7. the method for correcting digital images for zoom optical imaging system according to claim 1, is characterized in that,

The real work focal length value of described video camera is detected by the position transducer on the focus adjusting mechanism on described video camera and obtains.

8. the method for correcting digital images for zoom optical imaging system according to claim 1, is characterized in that,

The information that position and attitude data during described video camera imaging are detected by the angular transducer of GPS locator data and described video camera inside calculates.

9. the method for correcting digital images for zoom optical imaging system according to claim 1, is characterized in that,

Described video camera is aviation zoom camera.