CN103679711B - A kind of remote sensing satellite linear array push sweeps optics camera outer orientation parameter calibration method in-orbit - Google Patents

A kind of remote sensing satellite linear array push sweeps optics camera outer orientation parameter calibration method in-orbit Download PDF

Info

Publication number
CN103679711B
CN103679711B CN201310632138.9A CN201310632138A CN103679711B CN 103679711 B CN103679711 B CN 103679711B CN 201310632138 A CN201310632138 A CN 201310632138A CN 103679711 B CN103679711 B CN 103679711B
Authority
CN
China
Prior art keywords
reference mark
pixi
point
error
pointing vector
Prior art date
Application number
CN201310632138.9A
Other languages
Chinese (zh)
Other versions
CN103679711A (en
Inventor
马灵霞
郝胜勇
王剑
邹同元
郭翠翠
马狄
张荞
丁火平
纪强
李宝明
Original Assignee
航天恒星科技有限公司
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by 航天恒星科技有限公司 filed Critical 航天恒星科技有限公司
Priority to CN201310632138.9A priority Critical patent/CN103679711B/en
Publication of CN103679711A publication Critical patent/CN103679711A/en
Application granted granted Critical
Publication of CN103679711B publication Critical patent/CN103679711B/en

Links

Abstract

A kind of remote sensing satellite linear array push sweeps optics camera outer orientation parameter calibration method in-orbit, and step is: (1) gathers high-precision control point information; (2) obtain the auxiliary datas such as reference mark place camera internal position element, and track, set up the strict imaging geometry model in each reference mark place; (3) calculate the positioning error at each reference mark place, excluding gross error reference mark, obtain preliminary calibration reference mark information; (4) set up camera outer orientation element geometry calibration mathematical model, calculate each reference mark process opinion pointing vector and actual pointing vector; (5) set up three axle rotation error equations between theoretical pointing vector and actual pointing vector, equation is carried out linearization process; (6) method of least squares resolves error equation parameter, obtains demarcating outer parameter. After the present invention utilizes the outer parameter of demarcation to carry out systematic error compensation, camera is stablized without control positioning precision.

Description

A kind of remote sensing satellite linear array push sweeps optics camera outer orientation parameter calibration method in-orbit

Technical field

The present invention relates to a kind of remote sensing satellite linear array push-broom type optics camera outer orientation parameter calibration method in-orbit, belong to remote sensing satellite image processing technology field, for the optimization of the positioning precision in-orbit lifting of the linear push-broom type camera that remote sensing satellite carries.

Background technology

Along with the continuous transmitting of domestic and international remote sensing satellite, the resolving power of remote sensing satellite improves constantly, and the features such as the big area of satellite remote sensing, short period, low cost become the dominant direction of remote sensing development. Ended for the end of the year 2012, China has transmitted multi-series many Optical remote satellites such as environment, resource, remote sensing, but a large amount of domestic remotely-sensed datas is not fully used, and traces it to its cause, geometric positioning accuracy is low is the important factor of subsequent applications effect of a restriction remotely-sensed data. External commercial optical remote sensing satellite geometric positioning accuracy reaches the precision of rice level, and such as Geoeye, WorldView etc., and in the remote sensing satellite that China has launched, resource 02 star positioning precision is 7Km, resource No. two 03 star positioning precision 200m. Domestic remote sensing satellite is launch the Accurate Calibration that camera is not carried out outer orientation element by rear ground without a major reason of control positioning precision difference.

Linear array push-broom type optics camera is the main flow of current Optical remote satellite sensor, and linear array push-broom type optics camera pushes ahead dynamic imaging with satellite along the track pre-defined, and forms a width bidimensional image one by one after scanning. In certain photography moment, row corresponding on image with taken the photograph ground and be there is strict central projection relation. By ground survey error, the impact launching the factors such as rear environment impulse force, the outer orientation parameter after remote sensing satellite transmitting changes so that it seriously declines without control geometric positioning accuracy, lacks the direct demarcation means to satellite in orbit outer orientation element in floor treatment.

Summary of the invention

The technical problem that the present invention solves is: overcome the deficiencies in the prior art, a kind of remote sensing satellite linear array push is provided to sweep optics camera outer orientation parameter calibration method in-orbit, utilize a reference mark, scape image collection accurate ground, the rigorous geometry model of camera can be revised, resolve accurate outer orientation parameter, can be used for promoting the nothing control geometric positioning accuracy that camera obtains image, support the geometry calibration of remote sensing satellite linear array push-broom type optics camera.

The technical scheme of the present invention is: a kind of remote sensing satellite linear array push sweeps optics camera outer orientation parameter calibration method in-orbit, comprises the following steps:

1) calibration image I is being treated1, control image I2On utilize the automatic matching algorithm in reference mark, gather reference mark information, described reference mark information is that T registering control points is to { PixI1i, PixI2i, wherein i=1,2,3 ..., T); Record each dominating pair of vertices { PixI1i, PixI2iTreat calibration image I1Picture point coordinate (i, j) and the corresponding control image unique point position coordinate (lat, lon, height) under WGS84 system of coordinates;

2) camera internal position element and the auxiliary data at each reference mark place is obtained, when described auxiliary data comprises the track of satellite, attitude, row; Set up the rigorous geometry model at each reference mark place X G Y G Z G = X s Y s Z s + uM 3 M 2 M 1 M 0 tan p s i X - tan p s i Y 1 ; Wherein [XS,YS,ZS] inscribe the position of satellite in agreement geocentric coordinates system constantly for row that each reference mark is corresponding, i.e. track position, reference mark place (PX, PY, PZ); [XG,YG,ZG] coordinate in agreement geocentric coordinates system of terrain object point corresponding to the corresponding picture unit in each reference mark; PsiX, psiY be respectively each reference mark correspondence image picture unit primary optical axis unit vector and satellite body system of coordinates X-axis, Y-axis angle; U is scale factor; M0For satellite body system of coordinates is relative to the installation matrix of camera, obtain by ground survey before satellite is launched; M1For the row that each reference mark is corresponding inscribes satellite constantly to track system of coordinates rotation matrix, it is made up of the attitude angle measured on star; M2For the row that each reference mark is corresponding inscribes track constantly to J2000.0 system of coordinates rotation matrix; M3For inscribing J2000.0 to WGS84 system of coordinates rotation matrix time this;

3) according to step 2) rigorous geometry model that obtains, calculate the positioning error size and direction that obtain each reference mark place, preliminary excluding gross error reference mark, obtains initial control point pair;

4) camera outer orientation element geometry calibration mathematical model is set up, according to step 3) initial control point pair that obtains, calculate and obtain each reference mark process opinion pointing vector and actual pointing vector;

41) in rigorous geometry model, introduce error matrix M, set up camera outer orientation element geometry calibration mathematical model

X G Y G Z G = X s Y s Z s + uM 3 M 2 M 1 M 0 tan p s i X - tan p s i Y 1 ;

42) to arbitrary group of reference mark, element (psiX, psiY) calculation control point place, the interior orientation optical axis at reference mark place is utilized to point to vector V in camera coordinates system, as the actual pointing vector V at this reference mark place after being normalizedActual;

43) to arbitrary group of reference mark, utilize reference mark ground point location coordinate (lat, lon, height) calculation control point place optical axis points to vector V in camera coordinates system, as the theoretical pointing vector V at this reference mark place after being normalizedTheoretical;

Wherein X G Y G Z G = ( a + h e i g h t ) × c o s ( l o n ) × c o s ( l a t ) ( a + h e i g h t ) × s i n ( l o n ) × c o s ( l a t ) ( b + h e i g h t ) × sin ( l a t ) , A is earth semi-major axis, and b is earth minor semi-axis;

5) set up three axle rotation error equations between theoretical pointing vector and actual pointing vector, this three axles rotation error equation is carried out linearization process, obtains system of linear equations;

51) set error matrix M as three axle rotating orthogonal matrixes, set up the matrix conversion model between theoretical pointing vector and actual pointing vector:

Wherein X Y Z For theoretical pointing vector, x y z For actual pointing vector;

52) model is changed, obtain the nonlinear equation of the theoretical each parameter of pointing vector with each parameter of actual pointing vector;

53) nonlinear equation is carried out Taylor expansion at initial point (0,0,0) place, obtain the linear-apporximation equation of nonlinear equation;

54) for each group of theoretical pointing vector and actual pointing vector, calculate linear equation parameter, obtain system of linear equations;

6) method of least squares iteration process of solution 5 is adopted) system of linear equations that obtains, obtain and demarcate outer parameter.

Step 1) in gather reference mark information concrete grammar be:

11) calibration original image I is treated based on SIFT algorithm1Carry out feature point extraction, obtain M unique point PixI1i(i=1,2,3 ..., M), M is positive integer; Record the SIFT feature vector of each unique point;

12) based on SIFT algorithm, high precision is controlled image I2Carry out feature point extraction, obtain N number of unique point PixI2i(i=1,2,3 ..., N), N is positive integer; Record the SIFT feature vector at each unique point place;

13) adopt Europe formula distance the unique point of two width images to be mated as similarity measurement criterion, obtain T registering control points to { PixI1i, PixI2i(i=1,2,3 ..., T).

Step 2) in set up each reference mark place the concrete grammar of rigorous geometry model be:

21) according to dominating pair of vertices { PixI1i, PixI2iImage coordinate (i, j), when obtaining the row of its correspondence;

22) Lagrange's interpolation algorithm interpolation calculation dominating pair of vertices { PixI is utilized1i, PixI2iCorresponding row time track position (PX, PY, PZ, VX, VY, VZ);

23) attitude quaternion passed down by satellite carries out ordinate transform process, utilizes Lagrange algorithm interpolation calculation dominating pair of vertices { PixI1i, PixI2iCorresponding row time moment place's camera relative to the three-axis attitude angle (Roll, Pitch, Yaw) of track system of coordinates;

24) according to camera laboratory measurements or interior orientation element geometry calibration result, calculation control point is to { PixI1i, PixI2iCorresponding CCD visits unit's optical axis and points to angle (psiX, psiY);

25) according to step 21)-step 24) result that obtains, set up the rigorous geometry model at reference mark place, for certain on the image that the linear array push-scanning image camera a certain moment obtains a bit, build the rigorous geometry model of remote sensing image;

Step 3) in the concrete acquisition methods of initial control point pair be:

31) according to step 2) rigorous geometry model that obtains, calculate each dominating pair of vertices { PixI1i, PixI2iIn treat calibration imaging point PixI1iInitial ground point location coordinate { lon ', lat ' };

32) according to each dominating pair of vertices { PixI1i, PixI2iMiddle control imaging point PixI2iLocation coordinate information, calculation control point PixI1iInitial fix error { �� xi, �� yi};

33) calculation control point is to { PixI1i, PixI2iPositioning error mean value { �� x, �� y};

34) to each dominating pair of vertices, judge its initial fix error �� xi, �� yi} whether �� x �� 20% �� x, within the scope of �� y �� 20% �� y}, will not { dominating pair of vertices within the scope of �� x �� 20% �� x, �� y �� 20% �� y} be rejected;

35) reject error control point to rear, obtain K group initial control point to { PixI1i, PixI2i(i=1,2,3 ..., K).

The present invention's useful effect compared with prior art is:

(1) remote sensing satellite linear array push-broom type optics camera rigorous geometry model and auxiliary data characteristic are conducted in-depth analysis by the present invention, the auxiliary data optimization process such as devise when comprising the information acquisition of reference mark, image high-precision ground, track camera, attitude, row, parameter calibration flow process outside strict geometry calibration mathematical model is set up, desirable pointing vector calculates with actual pointing vector, triaxial error establishing equation and equation linearization process, least square resolve the steps such as outer orientation parameter high precision, meet satellite geometry calibration demand in-orbit.

(2) in the information extraction of reference mark, high precision ground, for features such as the rotation of basic image and high-precision control point image, yardsticks, have employed the Control point extraction based on SIFT algorithm and coupling, it is to increase the accuracy of Control point extraction.

(3) in whole calibration process, it is contemplated that fault-tolerant based on the initial control point of rigorous geometry model, registering control points to screening, is ensured to gather the exactness of dominating pair of vertices by the principle that utilizes image each point positioning error basically identical.

(4) when setting up geometry calibration mathematical model, introduce three axle rotating orthogonal matrixes as error matrix, after error equation linearizing, the basis reducing the difficulty that error equation resolves well ensures that linearized stability equation is to the matching similarity of true error.

(5) have employed the iteration error based on least square when resolving outer orientation parameter and reject design, it is possible to further excluding gross error reference mark, it is to increase the precision of parameter calibration.

Accompanying drawing explanation

Fig. 1 is the inventive method schema.

Embodiment

Below in conjunction with accompanying drawing 1, the specific embodiment of the present invention is further described in detail:

1. treating calibration image I1, high precision control image I2On utilize the automatic matching algorithm in reference mark to gather high-precision control point to information,

(1) calibration original image I is treated based on SIFT algorithm1Carry out feature point extraction, obtain m unique point PixI1i(i=1,2 ..., m), record the SIFT feature vector of each unique point;

SIFT feature matching algorithm characterizes as follows:

(1.1) characteristic point position coordinate and place yardstick is determined. Set up image gaussian pyramid, 26 neighborhoods in pyramid scale space detect extreme value, if certain point (x, y) is maximum in 26 neighborhoods of this layer of pyramid scale space and upper and lower two layers or during minimum value, defining this point is the unique point of image under this yardstick L.

(1.2) unique point direction parameter calculates. Sample in the neighborhood window centered by unique point (x, y), with the gradient direction of statistics with histogram neighborhood territory pixel. The scope of histogram of gradients is divided into 36 directions, and every 10 degree represent a direction. The direction at definition histogram peak place is the direction parameter of this unique point.

(1.3) SIFT feature vector description is calculated. Around unique point center, get the window of 8*8, window is cut into the sub-window of 2*2. Every sub-window calculates the gradient orientation histogram in 8 directions, adds up the cumulative value in each direction, as the directional information of this sub-window. 4 sub-windows are the final 32 dimensional characteristics vectors generating this unique point place after having added up.

(2) based on SIFT algorithm, high precision is controlled image I2Carry out feature point extraction, obtain n unique point PixI2i(i=1,2 ..., n), record the SIFT feature vector at each unique point place;

(3) adopt Europe formula distance the unique point of two width images to be mated as similarity measurement criterion, obtain M registering control points to { PixI1i, PixI2i(i=1,2 ..., M);

(3.1) for treating calibration original image I1Arbitrary reference mark PixI1i, calculate its SIFT feature vector and high precision control image I2Europe formula distance between n the proper vector that upper extraction obtains. For two vectorial l1(x1,x2,��,xn)��l2(y1,y2,��,yn), its Europe formula distance characterizes as follows:

Δ l = Σ i = 0 n ( x i - y i ) 2

(3.2) minimum value in m Europe formula distance is calculated, the unique point PixI at minimum value place2jIt is and treats calibration original image I1Unique point PixI1iMatching characteristic point.

(4) to each matching characteristic point to { PixI1i, PixI2j, calibration original image I treated in record1Unique point PixI1iPlace's picture point planimetric coordinates (sample, line), and high precision control image I2Unique point PixI2jPlace's three-dimensional location coordinates (lat, lon, height), as the high-precision control point information collected.

2. the auxiliary data such as when obtaining interior orientation, place, each reference mark element, track, attitude, row, sets up the rigorous geometry model at each reference mark place.

(1) data such as when resolving track in calibration video imaging time range, attitude, row in the raw data passed from satellite;

(2) to arbitrary dominating pair of vertices Pointsi{ sample, line, lat, lon, height}, according to the image coordinate (sample, line) of dominating pair of vertices, obtain the photography time scanTime of its correspondence;

Directly resolving in the auxiliary data that photography time corresponding to reference mark can pass from satellite, imaging time corresponding to the capable auxiliary data of line is photography time corresponding to this reference mark.

(3) satellite orbital position (PX, PY, PZ, VX, VY, VZ) of Lagrange's interpolation algorithm interpolation calculation photography time scanTime is utilized; Satellite passes orbital data according under certain frequency, and therefore, before and after the satellite orbital position needs utilization photography moment that photography time scanTime is corresponding, the orbital data of several groups carries out interpolation calculation. The present invention adopts Lagrange's interpolation algorithm, utilizes the satellite orbital position of three groups of orbital data calculating photography times before and after photography time.

(3.1) from first group of orbital data, judge the relation between the generation time of this group orbital data, the generation time of next group orbital data and scanTime, if scanTime is greater than the generation time of i-th group of orbital data, the generation time being simultaneously less than the i-th+1 group orbital data, then recording i is the distance nearest orbital data sequence number of photography time.

(3.2) utilize the i-th-1, i, i+1 group orbital data, calculate track position and the speed in photography moment based on Lagrange algorithm. Lagrange's interpolation algorithm is expressed as follows:

For the function table (x of known y=f (x)i,f(xi)) (i=0,1 ..., n), at [xo,xn] arbitrary x in scope, have:

P n ( x ) = Σ i = 0 n ( y i * Π j = 0 , j ≠ i n x - x j x i - x j )

(4) utilize Lagrange algorithm interpolation calculation photography time scanTime camera relative to the three-axis attitude angle (Roll, Pitch, Yaw) of track system of coordinates;

(5) according to camera laboratory measurements, the sample the CCD reading image coordinate (sample, the line) correspondence of reference mark dominating pair of vertices visits optical axis sensing angle (psiX, psiY) of unit.

(6) rigorous geometry model at reference mark place is set up, for arbitrary reference mark, the strict imaging geometry model that the strict geometry imaging model Model linear array push of data construct remote sensing image sweeps camera when utilizing its interior orientation, track, attitude, row etc. as follows described in.

X G Y G Z G = X s Y s Z s + uM 3 M 2 M 1 M 0 tan p s i X - tan p s i Y 1

Wherein:

[XS,YS,ZS] for inscribing the position of satellite in agreement geocentric coordinates system time this, i.e. track position, reference mark place (PX, PY, PZ);

[XG,YG,ZG] it is the coordinate of terrain object point corresponding to this picture unit in agreement geocentric coordinates system.

PsiX, psiY be respectively image corresponding camera picture unit primary optical axis unit vector and satellite body system of coordinates X-axis, Y-axis angle

U scale factor.

M0For satellite body system of coordinates is relative to the installation matrix of camera, obtain by ground survey before satellite is launched;

M1For inscribing satellite time this to track system of coordinates rotation matrix, it is made up of the attitude angle measured on star.

M2For this moment lower railway is to J2000.0 system of coordinates rotation matrix, it is made up of dragon's head right ascension, orbital inclination, argument etc.

M3For inscribing J2000.0 to WGS84 system of coordinates rotation matrix time this, precession of the equinoxes correction, nutating correction, Greenwich sidereal time correction and pole need to be carried out and move correction.

3. calculating positioning error size and the direction at each reference mark place, preliminary excluding gross error reference mark, obtains initial control point pair

(1) according to the rigorous geometry model set up, based on the auxiliary data such as when interior orientation, place, each reference mark element, track, attitude, row, calculate the initial ground point location coordinate (lon ', lat ') of original picture point (sample, line) in each reference mark

(2) image spot location (lon, lat, height) coordinate information is controlled according to each reference mark high precision, the initial fix error (�� xi, �� yi) of scaling system picture point. For each reference mark,

(2.1) by initial ground point location coordinate (lon ', lat ') after utm projection, its planimetric coordinates (x ', y ') is obtained

(2.2) by reference mark imaging point coordinate (lon, lat, height) after utm projection, its planimetric coordinates (x, y) is obtained

(2.3) definition initial fix error (�� xi, �� yi)=(x-x ', y-y ')

(3) calculation control point is to positioning error mean value (�� x, �� y)

(4) to each dominating pair of vertices, judge that its initial fix error { whether in [�� x �� 20% �� x, �� y �� 20% �� y] scope, if not within the scope of this, rejected by �� xi, �� yi} by this reference mark. After rejecting error control point, obtain initial control point pair, for setting up outer orientation element geometry calibration mathematical model

4. set up camera outer orientation element geometry calibration mathematical model, calculate each reference mark process opinion pointing vector and actual pointing vector.

(1) set up camera outer orientation element geometry calibration mathematical model, in strict imaging geometry model, introduce error matrix M, as follows

X G Y G Z G = X s Y s Z s + uM 3 M 2 M 1 M 0 tan p s i X - tan p s i Y 1

(2) to arbitrary group of reference mark, the actual pointing vector V at original image picture point (sample, line) place is calculatedActual:

Element (psiX, psiY) calculation control point place, the interior orientation optical axis at reference mark place is utilized to point to vector V in camera coordinates system, as the actual pointing vector V at this reference mark place after being normalizedActual;

(3) to arbitrary group of reference mark, the theoretical pointing vector V at calculation control point placeTheoretical:

(3.1) utilizing reference mark ground point location coordinate (lat, lon, height) to calculate the three axle positions (XG, YG, ZG) of ground millet cake in agreement geocentric coordinates system, calculation formula is as follows:

X G Y G Z G = ( a + h e i g h t ) × c o s ( l o n ) × c o s ( l a t ) ( a + h e i g h t ) × s i n ( l o n ) × c o s ( l a t ) ( b + h e i g h t ) × sin ( l a t )

(3.2) the three-axis attitude angle (Roll, Pitch, Yaw) at reference mark place is utilized to calculate by camera coordinates system to the rotation matrix M1 of track system of coordinates

(3.3) utilize the satellite orbital position (PX, PY, PZ, VX, VY, VZ) at reference mark place to calculate the track six roots of sensation number of satellite under J2000.0 system of coordinates, thus calculate the rotation matrix M2 of track system of coordinates to J2000.0 system of coordinates

(3.4) utilize the imaging time at reference mark place, calculate the rotation matrix M3 of J2000.0 system of coordinates to WGS84 system of coordinates

(3.5) the theoretical pointing vector V at calculation control point placeTheoretical, calculation formula is:

5. set up three axle rotation error equations between theoretical pointing vector and actual pointing vector, equation is carried out linearization process

(1) set error matrix M as three axle rotating orthogonal matrixes, set up the matrix conversion model between theoretical pointing vector and actual pointing vector; Medial error matrix M of the present invention is characterized by around X, Y, the orthogonal matrix that the angle that Z coordinate axle rotates is formed, as follows:

Wherein:

(2) model is changed, obtain the nonlinear equation of the theoretical each parameter of pointing vector with each parameter of actual pointing vector, as follows:

(3) nonlinear equation is carried out Taylor expansion at initial point (0,0,0) place, obtain the linear-apporximation equation of nonlinear equation.

Wherein:

a 1 = - X Y Z 2

a 2 = 1 + X 2 Z 2

a 3 = - Y Z

c 1 = x - X Z

b 1 = - ( 1 + Y 2 Z 2 )

b 2 = X Y Z 2

b 3 = X Z

c 2 = y - Y Z

(4) for each group of theoretical pointing vector and actual pointing vector, calculate linear equation parameter, obtain system of linear equations

Utilize the system of linear equations between three axle rotation error establishing equations its theoretical pointing vector and actual pointing vector, for N number of reference mark, system of linear equations can be obtained as follows:

6. method of least squares iteration resolves error equation parameter, obtains demarcating outer parameter

(6.1) the error equation M1 set up for N number of reference mark, least-squares iteration resolves equation parameter

For the large scale sparse linear equations Ax=b that error equation calculates out, adopt iterative method solving equation numerical solution. For equation Ax=b, set up its recursion formula:

xi+1=(ATA+E)-1ATb+(AtA+E)-1xi

Set up least square solution error equation:

��=xi+1-xi=(ATA+E)-1ATb+(AtA+E)-1xi-xi

The orthogonal each parameter of error matrix M is all similar to and equals 0, therefore for x composes initial value (0,0,0).

Sequence vector x is resolved according to recursion formula0,x1,x2,...,xn, utilize minimum mean-square error equation solver ��1,��2,...,��n, until ��k�� min, nowIt is the least square solution of this equation.

(6.2) using the least square solution of equation M1 as preliminary outer orientation parameter calibration result, substitute into the strict imaging geometry model at each reference mark place, calculate the calibration residual error at each reference mark place, the reference mark that deleted residual is big.

(6.3) utilizing M reference mark after rejecting excessive residual error, set up error equation M2, least-squares iteration resolves equation parameter. The least square solution calculated is as outer orientation parameter calibration result. Calibration result is applied in Image correction in remote sensing algorithm, strict imaging geometry model parameter is revised, it is possible to significantly improve the nothing control geometric positioning accuracy of this camera remote sensing image product.

As shown in table 1-table 3, the outer orientation parameter calibration result that table 1 calculates for the present invention, table 2 is not for using VRSS-1 satellite geometry positioning precision measuring result before the present invention, and table 3 carries out the geometric positioning accuracy measuring result of VRSS-1 satellite after calibration in-orbit for utilizing the present invention. VRSS-1 satellite outer orientation parameter is being rolled and all there occurs bigger change in pitching two as can be seen from Table 1, reach 0.06 measurement level, average orbit altitude (official's nominal value is 640 kilometers) measuring and calculating according to VRSS-1 satellite, the outer parameter that rotating direction-0.06239 is spent can cause ground to be about 640000 �� tan (-0.06239)=696.1 meter without control positioning precision; Pitch orientation 0.06804 is outside one's consideration, and without controlling, positioning precision is about 640000 �� tan (0.06804)=760.07 meter on parameter inducible ground; The inducible overall error of rolling pitch orientation isBefore adopting the present invention, the actual measurement of satellite is 1073.18 meters of (mean value) magnitudes without control positioning precision as can be seen from Table 2, it is seen that what outer orientation parameter caused is consistent with measured value substantially without control positioning error. As can be seen from Table 3, after adopting the present invention that satellite is carried out outer orientation parameter calibration, the actual measurement of satellite is 36.31 meters of (mean value) magnitudes without control positioning precision, and the nothing control positioning precision of this satellite is improve 29 times by the present invention.

Table 1

Error matrix parameter Calibration result (unit: degree) Roll -0.062391107351712 Pitch 0.068045533894513 Yaw 0.000843235759485

Table 2

Table 3

The part that the present invention does not elaborate belongs to techniques well known.

Claims (4)

1. a remote sensing satellite linear array push sweeps optics camera outer orientation parameter calibration method in-orbit, it is characterised in that comprise the following steps:
1) calibration image I is being treated1, control image I2On utilize the automatic matching algorithm in reference mark, gather reference mark information, described reference mark information is that T registering control points is to { PixI1i, PixI2i, wherein i=1,2,3 ..., T; Record each dominating pair of vertices { PixI1i, PixI2iTreat calibration image I1Picture point coordinate (i, j) and the corresponding control image unique point position coordinate (lat, lon, height) under WGS84 system of coordinates;
2) camera internal position element and the auxiliary data at each reference mark place is obtained, when described auxiliary data comprises the track of satellite, attitude, row; Set up the rigorous geometry model at each reference mark placeWherein [XS,YS,ZS] inscribe the position of satellite in agreement geocentric coordinates system constantly for row that each reference mark is corresponding, i.e. track position, reference mark place (PX, PY, PZ); [XG,YG,ZG] coordinate in agreement geocentric coordinates system of terrain object point corresponding to the corresponding picture unit in each reference mark; PsiX, psiY be respectively each reference mark correspondence image picture unit primary optical axis unit vector and satellite body system of coordinates X-axis, Y-axis angle; U is scale factor; M0For satellite body system of coordinates is relative to the installation matrix of camera, obtain by ground survey before satellite is launched; M1For the row that each reference mark is corresponding inscribes satellite constantly to track system of coordinates rotation matrix, it is made up of the attitude angle measured on star; M2For the row that each reference mark is corresponding inscribes track constantly to J2000.0 system of coordinates rotation matrix; M3For inscribing J2000.0 to WGS84 system of coordinates rotation matrix time this;
3) according to step 2) rigorous geometry model that obtains, calculate the positioning error size and direction that obtain each reference mark place, preliminary excluding gross error reference mark, obtains initial control point pair;
4) camera outer orientation element geometry calibration mathematical model is set up, according to step 3) initial control point pair that obtains, calculate and obtain each reference mark process opinion pointing vector and actual pointing vector;
41) in rigorous geometry model, introduce error matrix M, set up camera outer orientation element geometry calibration mathematical model
42) to arbitrary group of reference mark, element (psiX, psiY) calculation control point place, the interior orientation optical axis at reference mark place is utilized to point to vector V in camera coordinates system, as the actual pointing vector V at this reference mark place after being normalizedActual;
43) to arbitrary group of reference mark, utilize reference mark ground point location coordinate (lat, lon, height) calculation control point place optical axis points to vector V in camera coordinates system, as the theoretical pointing vector V at this reference mark place after being normalizedTheoretical;
WhereinA is earth semi-major axis, and b is earth minor semi-axis;
5) set up three axle rotation error equations between theoretical pointing vector and actual pointing vector, this three axles rotation error equation is carried out linearization process, obtains system of linear equations;
51) set error matrix M as three axle rotating orthogonal matrixes, set up the matrix conversion model between theoretical pointing vector and actual pointing vector:
WhereinFor theoretical pointing vector,For actual pointing vector;
52) model is changed, obtain the nonlinear equation of the theoretical each parameter of pointing vector with each parameter of actual pointing vector;
53) nonlinear equation is carried out Taylor expansion at initial point (0,0,0) place, obtain the linear-apporximation equation of nonlinear equation;
54) for each group of theoretical pointing vector and actual pointing vector, calculate linear equation parameter, obtain system of linear equations;
6) method of least squares iteration process of solution 5 is adopted) system of linear equations that obtains, obtain and demarcate outer parameter.
2. a kind of remote sensing satellite linear array push according to claim 1 sweeps optics camera outer orientation parameter calibration method in-orbit, it is characterised in that: step 1) in gather reference mark information concrete grammar be:
11) calibration original image I is treated based on SIFT algorithm1Carry out feature point extraction, obtain M unique point PixI1i, i=1,2,3 ..., M, M are positive integer; Record the SIFT feature vector of each unique point;
12) based on SIFT algorithm, high precision is controlled image I2Carry out feature point extraction, obtain N number of unique point PixI2i, i=1,2,3 ..., N, N are positive integer; Record the SIFT feature vector at each unique point place;
13) adopt Europe formula distance the unique point of two width images to be mated as similarity measurement criterion, obtain T registering control points to { PixI1i, PixI2i, i=1,2,3 ..., T.
3. a kind of remote sensing satellite linear array push according to claim 1 sweeps optics camera outer orientation parameter calibration method in-orbit, it is characterised in that: step 2) in set up each reference mark place the concrete grammar of rigorous geometry model be:
21) according to dominating pair of vertices { PixI1i, PixI2iImage coordinate (i, j), when obtaining the row of its correspondence;
22) Lagrange's interpolation algorithm interpolation calculation dominating pair of vertices { PixI is utilized1i, PixI2iCorresponding row time track position (PX, PY, PZ, VX, VY, VZ);
23) attitude quaternion passed down by satellite carries out ordinate transform process, utilizes Lagrange algorithm interpolation calculation dominating pair of vertices { PixI1i, PixI2iCorresponding row time moment place's camera relative to the three-axis attitude angle (Roll, Pitch, Yaw) of track system of coordinates;
24) according to camera laboratory measurements or interior orientation element geometry calibration result, calculation control point is to { PixI1i, PixI2iCorresponding CCD visits unit's optical axis and points to angle (psiX, psiY);
25) according to step 21)-step 24) result that obtains, set up the rigorous geometry model at reference mark place, for certain on the image that the linear array push-scanning image camera a certain moment obtains a bit, build the rigorous geometry model of remote sensing image.
4. a kind of remote sensing satellite linear array push according to claim 1 sweeps optics camera outer orientation parameter calibration method in-orbit, it is characterised in that: step 3) in the concrete acquisition methods of initial control point pair be:
31) according to step 2) rigorous geometry model that obtains, calculate each dominating pair of vertices { PixI1i, PixI2iIn treat calibration imaging point PixI1iInitial ground point location coordinate { lon ', lat ' };
32) according to each dominating pair of vertices { PixI1i, PixI2iMiddle control imaging point PixI2iLocation coordinate information, calculation control point PixI1iInitial fix error { �� xi, �� yi};
33) calculation control point is to { PixI1i, PixI2iPositioning error mean value { �� x, �� y};
34) to each dominating pair of vertices, judge its initial fix error �� xi, �� yi} whether �� x �� 20% �� x, within the scope of �� y �� 20% �� y}, will not { dominating pair of vertices within the scope of �� x �� 20% �� x, �� y �� 20% �� y} be rejected;
35) reject error control point to rear, obtain K group initial control point to { PixI1i, PixI2i, i=1,2,3 ..., K.
CN201310632138.9A 2013-11-29 2013-11-29 A kind of remote sensing satellite linear array push sweeps optics camera outer orientation parameter calibration method in-orbit CN103679711B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201310632138.9A CN103679711B (en) 2013-11-29 2013-11-29 A kind of remote sensing satellite linear array push sweeps optics camera outer orientation parameter calibration method in-orbit

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201310632138.9A CN103679711B (en) 2013-11-29 2013-11-29 A kind of remote sensing satellite linear array push sweeps optics camera outer orientation parameter calibration method in-orbit

Publications (2)

Publication Number Publication Date
CN103679711A CN103679711A (en) 2014-03-26
CN103679711B true CN103679711B (en) 2016-06-01

Family

ID=50317162

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201310632138.9A CN103679711B (en) 2013-11-29 2013-11-29 A kind of remote sensing satellite linear array push sweeps optics camera outer orientation parameter calibration method in-orbit

Country Status (1)

Country Link
CN (1) CN103679711B (en)

Families Citing this family (10)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN105374009A (en) * 2014-10-22 2016-03-02 航天恒星科技有限公司 Remote sensing image splicing method and apparatus
CN104613942B (en) * 2015-01-30 2017-01-18 北京林业大学 Analytic method of image pair relative elements of exterior orientation
CN105180963B (en) * 2015-07-22 2018-02-16 北京航空航天大学 Unmanned plane telemetry parameter modification method based on online calibration
CN105115477B (en) * 2015-07-27 2017-12-15 上海卫星工程研究所 To the spaceborne method for solving of ground point target push-scanning image task parameters
CN105136164B (en) * 2015-08-13 2019-04-05 航天恒星科技有限公司 Consider the comprehensive staring imaging emulation moved of satellite and method for evaluating quality and device
CN106595600B (en) * 2016-12-23 2019-04-30 航天恒星科技有限公司 The stereo mapping attitude of satellite measures the compensation method of system low frequency aberration and system
CN107085856A (en) * 2017-04-06 2017-08-22 上海航天测控通信研究所 A kind of in-orbit high-precision real-time location method based on optical image
CN108663535A (en) * 2018-05-30 2018-10-16 北京市遥感信息研究所 A method of air speed is estimated based on single scape high-resolution optical remote sensing image
CN109029931A (en) * 2018-08-02 2018-12-18 北京空间机电研究所 A kind of remote sensor pointing accuracy on-orbit calibration device and method
CN109931925B (en) * 2019-03-12 2019-12-27 中国人民解放军军事科学院国防科技创新研究院 Method for optimizing and estimating spinning attitude of space rolling satellite based on sequence image axis

Citations (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN103034994A (en) * 2012-11-20 2013-04-10 江南大学 Feature point calibration camera parameter in general case of camera shooting by using scale invariant feature transform (sift)

Family Cites Families (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
KR100912715B1 (en) * 2007-12-17 2009-08-19 한국전자통신연구원 Method and apparatus of digital photogrammetry by integrated modeling for different types of sensors

Patent Citations (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN103034994A (en) * 2012-11-20 2013-04-10 江南大学 Feature point calibration camera parameter in general case of camera shooting by using scale invariant feature transform (sift)

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
基于线阵光学图像的运动参数测量技术及其应用研究;赵竹新;《CNKI博士学位论文全文库》;20120501;10-17 24-34页 *

Also Published As

Publication number Publication date
CN103679711A (en) 2014-03-26

Similar Documents

Publication Publication Date Title
JP2009145314A (en) Digital photogrammetry by integrated modeling of different types of sensors, and its device
Li Potential of high-resolution satellite imagery for national mapping products
Zhang et al. An Unmanned Aerial Vehicle‐Based Imaging System for 3D Measurement of Unpaved Road Surface Distresses 1
Hu et al. Understanding the rational function model: methods and applications
GREJNER‐BRZEZINSKA Direct exterior orientation of airborne imagery with GPS/INS system: Performance analysis
Li et al. The global image of the Moon obtained by the Chang’E-1: Data processing and lunar cartography
Fraser et al. Sub-metre geopositioning with Ikonos GEO imagery
CN103983254B (en) The motor-driven middle formation method of a kind of novel quick satellite
CN101893440B (en) Celestial autonomous navigation method based on star sensors
Li et al. Rigorous photogrammetric processing of HiRISE stereo imagery for Mars topographic mapping
CN102506824B (en) Method for generating digital orthophoto map (DOM) by urban low altitude unmanned aerial vehicle
CN101344391B (en) Lunar vehicle posture self-confirming method based on full-function sun-compass
FR2953940A1 (en) Method for geo-referencing an image area
CN103674063B (en) A kind of optical remote sensing camera geometric calibration method in-orbit
CN102607526B (en) Target posture measuring method based on binocular vision under double mediums
CN103954283B (en) Inertia integrated navigation method based on scene matching aided navigation/vision mileage
CN104501779A (en) High-accuracy target positioning method of unmanned plane on basis of multi-station measurement
CN101762273A (en) Autonomous optical navigation method for soft landing for deep space probe
CN101852623B (en) On-track calibration method for internal element of satellite optical remote sensing camera
CN102168972A (en) RPC-based method for improving and calibrating block adjustment of three-linear array three-dimensional satellite
CN101520511A (en) Method for formation configuration of distributed satellites with synthetic aperture radars
Li et al. Impact of imaging geometry on 3D geopositioning accuracy of stereo IKONOS imagery
CN104931022A (en) Satellite image three-dimensional area network adjustment method based on satellite-borne laser height measurement data
CN104215258B (en) Method and system for measuring precision of angle measurement of vehicle theodolite
CN103616028B (en) A kind of starlight refraction autonomous navigation of satellite method based on single star sensor

Legal Events

Date Code Title Description
PB01 Publication
SE01 Entry into force of request for substantive examination
C10 Entry into substantive examination
GR01 Patent grant
C14 Grant of patent or utility model