METHOD OF FULLY AUTONOMOUS GEOMETRIC CALIBRATION FOR LINEAR-ARRAY REMOTE SENSING SATELLITES

Abstract

A method of fully autonomous geometric calibration for linear-array remote sensing satellite (LARSS) based on the joint observation for stars and earth by satellite, with the support of satellite's high maneuverability is proposed. This invention realizes the full-link processing from data acquisition to internal and external calibration. Based on the ultra-high attitude stability and agile maneuverability, this invention designs a joint observation mode for the star and the earth, which is suitable for autonomous geometric calibration. With the joint observations, this invention achieves the external calibration through the star observations acquired in the solar shadow area, and achieves the internal calibration through the ground overlapping images acquired in the solar illumination area. Therefore, the high-precision geometric imaging model of the LARSS would be restored by the method, under the condition without using the ground calibration sites.

Claims

1. A method of fully autonomous geometric calibration for LARSS relying on a joint observation for stars and earth by satellites, and the method comprising: obtaining star images in a solar shadow region using satellite's agile maneuvers and extracting star positions from the star images; constructing an on-orbit geometric external calibration model by introducing generalized installation angles into a rigorous star imaging model; estimating the generalized installation angles using the extracted sequential star observations, in which a collinearity relationship between imagery point and star point is used to construct an adjustment model, and a least squares adjustment algorithm is adopted for parameter estimation; obtaining ground overlapping images in a solar illumination region based on satellite's agile mobility and an overlap requirement, and identity dense corresponding imagery points from the overlapping images; constructing an internal calibration model by introducing a fitted viewing angle model into a rigorous earth imaging model; and estimating coefficients of the fitted viewing angle model using coplanar constraints between the overlapping images.

2. The method of claim 1, wherein a systematic geometric error in full imaging link of a linear-array remote sensing satellite (LARSS) is calibrated using the joint observation for stars and the ground overlapping images, instead of reference data of calibration sites, which overcomes shortcomings of traditional methods due to their strong dependence on ground calibration sites, reducing calibration cost and improving calibration timeliness.

3. The method of claim 1, wherein constructing an on-orbit geometric external calibration model by introducing a generalized installation angles into a rigorous star imaging model comprises expressing geometric errors of external orientation parameters using the generalized installation angles, determining an optical aberration correction angle, establishing a transform quaternion to correct optical aberration in the rigorous star imaging model, establishing the geometric external calibration model using the rigorous star imaging model.

4. The method of claim 3, wherein an optical aberration is corrected by a matrix in the rigorous star imaging model; by separating the optical aberration error in the star imaging model, it is ensured that the installation angles determined by the star-based external calibration is consistent with that in the earth imaging model of internal calibration, so as to ensure the complementarity of internal and external calibration without coupling.

5. The method of claim 1, wherein generalized installation angles are estimated using sequential star observations, and a least squares adjustment algorithm.

6. The method of claim 1, wherein overlap degree between the overlapping images collected by the same Charge-coupled Device (CCD) needs to be between 45% and 75%, and the optimal overlap degree is 65%.

7. The method of claim 1, wherein constructing an internal calibration model comprises expressing internal orientation errors of a camera using the fitted viewing angles of CCD detectors, constructing the rigorous earth imaging model, and constructing the internal calibration model by introducing the fitted viewing angle model into the rigorous earth imaging model.

8. The method of claim 1, wherein estimating coefficients of the fitted viewing angle model using coplanar constraints between overlapping images comprises constructing an adjustment model for internal calibration, performing absolute internal calibrations CCD by CCD, and performing an integrated relative internal calibration of all CCDs.

9. The method of claim 8, wherein absolute geometric distortions of each CCD and the relative geometric distortions among CCDs are both calibrated; therefore, all CCDs are registered under the same external installation angles.

10. The method of claim 8, wherein the adjustment model for internal calibration is established through using the common ground coordinates as the connection to express the coplanar condition of corresponding imagery points.

11. The method of claim 8, wherein absolute internal calibrations are performed CCD by CCD using the corresponding imagery points, and constant terms in the viewing angle model of each CCD are not estimated, because of their independence from coplanar conditions; due to the correlation between estimated parameters in internal and external calibration, not calculating constant term does not affect calibration accuracy, it can be considered that the errors of these constant terms have been compensated in external calibration.

12. The method of claim 8, wherein the corresponding imagery points between overlapping images of adjacent CCDs are used to perform relative internal calibration of all CCDs, and relative internal calibration parameters for all CCDs are estimated together through an overall adjustment, to suppress error accumulation among CCDs.

13. The method of claim 12, wherein a CCD is selected as a reference CCD, and the constant terms in the viewing angle models of all non-reference CCDs are estimated based on the reference CCD, to ensure the accuracy and stability of the overall adjustment.

14. A computing device for executing the method of claim 1, the device comprising a processor, and a memory that comprises instructions that, when executed by the processor, cause the processor to perform acts comprising: reading the collected star images, the ground overlapping imagery, the attitude auxiliary data, the orbit auxiliary data, and the imaging time auxiliary data into the memory, conducting the star points extraction and corresponding imagery points matching, conducting the external calibration and the internal calibration according to method of the disclosure, and then outputting the estimated accurate installation angles and coefficient of the viewing angle model into the memory.

Description

BRIEF DESCRIPTION OF THE DRAWINGS

[0014] FIG. 1 is a flowchart of the method of external calibration using star observations and the internal calibration using ground overlapping images;

[0015] FIG. 2 is a schematic diagram of acquiring star observations based on the agile mobility of the LARSS;

[0016] FIG. 3 is a schematic diagram of collecting ground overlapping images in the same orbit circle based on the agile mobility of the LARSS;

[0017] FIG. 4 is a diagram of the distributions of ground overlapping images and their corresponding imagery points; and

[0018] FIG. 5 is a schematic diagram of CCD distortion in different stages after absolute and relative calibration in internal geometric calibration.

DETAILED DESCRIPTION

[0019] FIG. 1 is a flowchart that illustrates the method of fully autonomous geometric calibration for the LARSS based on the joint observation of stars and earth by satellites. A further detailed description of the method is given below for each step in the embodiments.

[0020] Step 1: Obtain the star images in the solar shadow region based on satellite agile maneuvers, and extracting star positions from star images.

[0021] FIG. 2 is the schematic diagram showing the satellite from earth imaging to star imaging. Since the number, distribution and magnitude of stars vary from one celestial zone to another, these conditions of stars in the different celestial zones needs to be considered when collecting star images. The whole celestial region is divided into zones according to satellite camera's field of view, and the zone with suitable star number, distribution and magnitude is selected as the observation zone according to the imaging sensitivity of the camera. By adjusting the satellite's attitude, the sequential star-viewing images are obtained.

[0022] Then, the precise star points are extracted by image processing methods such as image denoising, binarization, edge extraction and center-of-mass fitting with the aid of the navigation catalog to obtain the right ascension and declination of each star point.

[0023] Step 2: Construct the on-orbit geometric external calibration model based on the generalized installation angles.

[0024] This step further includes the following sub-steps:

[0025] Step 2.1: According to the right ascension α and declination δ of a star point, its star observation vector in the celestial coordinate system can be determined by v.sub.star=(cos α cos δ, sin α cos δ, sin δ).sup.T, and then based on the co-linearity between the star observation vector and the satellite viewing vector v.sub.sat, a rigorous geometric imaging model of star imaging for the LARSS image can be established as follows:

[00001]$\begin{array}{cc}\left(\begin{array}{c}\tan {\phi }_s\\ \tan {\phi }_y\\ 1\end{array}\right)=\mu R_{\mathrm{Body}}^{\mathrm{Cam}}\left(\phi ,\omega ,\kappa \right)R_{J2000}^{\mathrm{Body}}{R}_{\mathrm{Aber}}\left[\begin{array}{c}\cos \alpha \cos \delta \\ \sin \mathrm{\alpha cos}\delta \\ \sin \delta \end{array}\right]& \left(1\right)\end{array}$

[0026] In which, φ.sub.x and φ.sub.y are the viewing-angle of the image point in the across and along the CCD directions within the camera. μ is the scaling factor; R.sub.J 2000.sup.Body is the rotation matrix from the J2000 coordinate system to the satellite body coordinate system, which is obtained by the attitude of the satellite; R.sub.Aber is the correction matrix for optical aberration, which can be determined by the optical aberration angle θ and the rotation vector n, the determination of the correction matrix R.sub.Aber as follows:

[00002]$\begin{array}{cc}{R}_{\mathrm{Aber}}={\left[\begin{array}{ccc}{p}_1^2-p_2^2-p_3^2+p_0^2& 2\left(p_1{p}_2+p_3{p}_0\right)& 2\left(p_1{p}_3-p_2{p}_0\right)\\ 2\left(p_1{p}_2-p_3{p}_0\right)& -p_1^2+p_2^2-p_3^2+p_0^2& 2\left(p_2{p}_3+p_1{p}_0\right)\\ 2\left(p_1{p}_3+p_2{p}_0\right)& 2\left(p_2{p}_3-p_1{p}_0\right)& -p_1^2-p_2^2+p_3^2+p_0^2\end{array}\right]}^T& \left(2\right)\end{array}$$\begin{array}{cc}{p}_0=n_x\cos \left(\frac{\theta }{2}\right),p_1=n_x\sin \left(\frac{\theta }{2}\right),p_2=n_y\sin \left(\frac{\theta }{2}\right),p_3=n_z\sin \left(\frac{\theta }{2}\right)& \left(3\right)\end{array}$$\begin{array}{cc}\theta =\nu \sin \beta /c& \left(4\right)\end{array}$$\begin{array}{cc}n=\left(n_x,n_y,n_z\right)=\frac{v_{\mathrm{star}}\times v_{\mathrm{sat}}}{.Math.v_{\mathrm{star}}\times v_{\mathrm{sat}}.Math.}& \left(5\right)\end{array}$

[0027] In which, P=(p.sub.0, p.sub.1, p.sub.2, p.sub.3) is the transformation quaternion; v is the velocity of the satellite; c is the speed of light; β is the included angle between the star observation vector v.sub.star and the satellite viewing vector v.sub.sat.

[0028] By separating the optical aberration error in the star imaging model, it is ensured that the installation angle determined by the star-based external calibration is consistent with that in the earth imaging model of internal calibration, so as to ensure the complementarity of the internal and external calibration without coupling.

[0029] Step 2.2: The external calibration model constructed in this invention is shown in Eq. (6). R.sub.Body.sup.Cam is a generalized installation matrix from the satellite body coordinate system to the camera coordinate system, and all the errors of the external orientation parameters are unified into this matrix for compensation, which is determined by the three camera installation angles (φ, ω, κ), as follows:

[00003]$\begin{array}{cc}{R}_{\mathrm{Body}}^{\mathrm{Cam}}=\left[\begin{array}{ccc}\cos \phi & 0& \sin \phi \\ 0& 1& 0\\ -\sin \phi & 0& \cos \phi \end{array}\right]\left[\begin{array}{ccc}1& 0& 0\\ 0& \cos \omega & -\sin \omega \\ 0& \sin \omega & \cos \omega \end{array}\right]\left[\begin{array}{ccc}\cos \kappa & -\sin \kappa & 0\\ \sin \kappa & \cos \kappa & 0\\ 0& 0& 1\end{array}\right]& \left(6\right)\end{array}$

[0030] Step 3: Estimate the external calibration parameters based on sequential star observations.

[0031] This step further includes the following sub-steps:

[0032] Step 3.1: The least squares adjustment is used to estimate the external calibration parameters, and the adjustment model for parameters solution is first established based on the geometric calibration model. The constructed adjustment model (G.sub.x, G.sub.y) is as follows:

[00004]$\begin{array}{cc}\{\begin{array}{c}{G}_x=\stackrel{\_}{X}-\stackrel{\_}{Z}.Math.\tan {\phi }_x\\ G_y=\stackrel{\_}{Y}-\stackrel{\_}{Z}.Math.\tan {\phi }_y\end{array}& \left(7\right)\end{array}$

[0033] In which,

[00005]$\begin{array}{cc}\left(\begin{array}{c}\stackrel{\_}{X}\\ \stackrel{\_}{Y}\\ \stackrel{\_}{Z}\end{array}\right)=\mu R_{\mathrm{Body}}^{\mathrm{Cam}}\left(\phi ,\omega ,\kappa \right)R_{J2000}^{\mathrm{Body}}{R}_{\mathrm{Aber}}\left[\begin{array}{c}\cos \alpha \cos \delta \\ \sin \mathrm{\alpha cos}\delta \\ \sin \delta \end{array}\right]& \left(8\right)\end{array}$

[0034] For each star point, the adjustment equation used for the parameters solution is established according to Eq. (7) based on its imaging time, attitude and current calibration parameters.

[0035] Step 3.2: Based on the adjustment equation, the error equation for i-th star point can be constructed by model linearization:

v.sub.E.sup.i=A.sub.ix−L  (9)

[0036] In which, x=[dφ, dω, dκ].sup.T is the correction vector for the installation angles; A.sub.i is the coefficient matrix of the error equation for a star point; L.sub.i is the corresponding constant vector of the error equation, as follows:

[00006]$\begin{array}{cc}{A}_i=\left[\begin{array}{ccc}\frac{\partial G_x}{\partial \phi }& \frac{\partial G_x}{\partial \omega }& \frac{\partial G_x}{\partial \kappa }\\ \frac{\partial G_y}{\partial \phi }& \frac{\partial G_y}{\partial \omega }& \frac{\partial G_y}{\partial \kappa }\end{array}\right],L_i={\left[\begin{array}{c}-G_x\\ -G_y\end{array}\right]}_i& \left(10\right)\end{array}$

[0037] Based on the least squares theory, the estimation of x can be deduced:

[00007]$\begin{array}{cc}x={\left(\underset{i=1}{\stackrel{m}{.Math.}}{A}_i^T{A}_i\right)}^T\left(\underset{i=1}{\stackrel{m}{.Math.}}{A}_i^T{L}_i\right)& \left(11\right)\end{array}$

[0038] In which, m is the number of star points.

[0039] The estimation of the external parameters is an iterative process, the current values of the parameters are updated according to the estimated corrections, and used as inputs in next solution. The iteration ends when the difference between the results of the two successive solutions is less than the limit threshold.

[0040] Step 4: Obtain the ground overlapping images in the solar illumination region based on satellite's agile mobility, and identity dense corresponding imagery points from overlapping images.

[0041] FIG. 3 is a schematic diagram of a LARSS acquiring the overlapping images. When the satellite orbits to the solar illumination region, the satellite begins to collect the overlapping images using its maneuvering ability in three phases. In “T-FWD” phase, an image covering the earth's surface is obtained. Then the satellite adjusts its attitude according to the overlap requirements in the calibration (Motorization phase). Finally, in the third “B-FWD” phase, another overlapping image satisfying overlapping degree is obtained.

[0042] FIG. 4 is the schematic diagram showing the overlap condition of the acquired two scenes overlapping images through three stages. Since the linear-array satellite camera is usually composed of multiple CCDs (it is assumed that there are 6 CCDs in FIG. 4), and the distortions of each CCD are different, each CCD needs to have different parameters to be calibrated in the internal calibration, and each CCD needs to be processed separately. Therefore, the two images acquired by each CCD needs to have a certain overlap degree T.sub.CCD, which is required in this invention to be between 45% and 75%, and the best T.sub.CCD is 65%. The relationship between the overall overlap degree T and that of each CCD is as follows:

[00008]$\begin{array}{cc}T=1-\frac{1-T_{CCD}}{n}& \left(12\right)\end{array}$

[0043] In which, n is the number of CCDs.

[0044] FIG. 4 also shows the distribution of identified corresponding imagery points in the overlapping area of the two scenes of images. In order to limit the effect of the random jitter of the attitude fitting error that can bring nonlinear error to the internal calibration, the corresponding imagery points are identified in a shorter section of the line (along the orbit) direction of the overlapping images.

[0045] Step 5: Construct the internal calibration model based on the coplanar constraints between overlapping images.

[0046] This step further includes the following sub-steps:

[0047] Step 5.1: Obtain the imaging time of the corresponding imagery points according to their imagery line number, and then interpolate the attitude and the orbit parameters according to the imaging time. Therefore, the rigorous imaging model for ground observation is established based on the relationship among the imagery point, ground point and projection center, the model is as follows:

[00009]$\begin{array}{cc}\left(\begin{array}{c}\tan {\phi }_x\\ \tan {\phi }_y\\ 1\end{array}\right)=\mu R_{\mathrm{Body}}^{\mathrm{Cam}}\left(\phi ,\omega ,\kappa \right)R_{J2000}^{\mathrm{Body}}{R_{\mathrm{WGS}84}^{J2000}\left[\begin{array}{c}{X}_g-X_{\mathrm{gps}}\\ Y_g-Y_{\mathrm{gps}}\\ Z_g-Z_{\mathrm{gps}}\end{array}\right]}_{\mathrm{WGS}84}& \left(13\right)\end{array}$

[0048] In which, R.sub.Body.sup.Cam is still the rotation matrix from the satellite body coordinate system to the camera coordinate system, determined by the above external calibration; R.sub.WGS84.sup.J2000 is the rotation matrix from the WGS84 coordinate system to the J2000 coordinate system, determined by the ephemeris parameters at the moment of imaging; R.sub.J2000.sup.Body is the same as that in the imaging model for the stars; (X.sub.gps, Y.sub.gps, Z.sub.gps) represents the coordinates of the GPS antenna phase center in the WGS84 coordinate system, obtained by the GPS receiver on the satellite; (X.sub.g, Y.sub.g, Z.sub.g) represents the ground 3D coordinates of the corresponding imagery points in the WGS84 coordinate system. The transformation between the geographic coordinates (Lat, Lon, Hei) (latitude, longitude, elevation) and 3D coordinates is as follows:

[00010]$\begin{array}{cc}\{\begin{array}{c}{X}_g=\left(N+\mathrm{Hei}\right)\cos \mathrm{Lat}\cos \mathrm{Lon}\\ Y_g=\left(N+\mathrm{Hei}\right)\cos \mathrm{Lat}\sin \mathrm{Lon}\\ Z_g=\left(N\left(1-e^2\right)+\mathrm{Hei}\right)\sin \mathrm{Lat}\end{array}& \left(14\right)\end{array}$

[0049] In which, N is the radius of the Earth's curvature in prime vertical, and e is the first eccentricity of the Earth's ellipsoid.

[0050] Step 5.2: For the viewing of each CCD detector of a linear-array camera, its viewing angle in the camera coordinate system can be expressed by (φ.sub.x, φ.sub.y) accurately. According to the distortion characteristics of the LARS S cameras, the viewing angle (φ.sub.x, φ.sub.y) of each detector can be fitted using two cubic polynomials, and then the internal calibration model based on the fitted viewing angle is constructed:

[00011]$\begin{array}{cc}\{\begin{array}{c}\tan {\phi }_x=a_0+a_{1}s+a_2{s}^2+a_3{s}^3\\ \tan {\phi }_y=b_0+b_{1}s+b_2{s}^2+b_3{s}^3\end{array}& \left(15\right)\end{array}$

[0051] In which, s is the number of CCD detector, and (a.sub.0, a.sub.1, a.sub.2, a.sub.3, b.sub.0, b.sub.1, b.sub.2, b.sub.3) are the coefficients of cubic polynomials, also the internal calibration parameters need to be determined.

[0052] By introducing the constructed internal calibration model of each CCD into the rigorous geometric imaging model, the on-orbit geometric calibration model is obtained, in which each CCD has its own set of internal calibration parameters to be solved.

[0053] Step 6: Stepwise estimation of the internal calibration parameters based on coplanar constraints between overlapping images.

[0054] This step further includes the following sub-steps:

[0055] FIG. 5 is a schematic diagram of CCD distortion in different stages after absolute and relative calibration in internal geometric calibration. Before internal geometric calibration, the CCDs of the camera have various distortions with inconsistent properties and sizes (rotational distortion, scaling distortion, optical distortion, etc., or comprehensive distortion coupled with these distortions), resulting in low image geometric quality and accuracy. After the absolute calibration, the absolute distortion of each CCD is well corrected, but there is still relative geometric distortion among CCDs, resulting in the inability of accurate splicing among CCDs. Further, the relative distortion is corrected in the relative calibration, to improve the geometric splicing accuracy among the CCDs.

[0056] Step 6.1: The least squares adjustment is used to solve the parameters in absolute internal calibration, so the adjustment model (16) for parameters solution is first established based on the calibration model.

[00012]$\begin{array}{cc}\{\begin{array}{c}{F}_x=\overline{U}-\overline{W}.Math.\tan {\phi }_x\\ F_y=\overline{V}-\overline{W}.Math.\tan {\phi }_y\end{array}& \left(16\right)\end{array}$

[0057] In which,

[00013]$\begin{array}{cc}\left(\begin{array}{c}\overline{U}\\ \overline{V}\\ \overline{W}\end{array}\right)=\mu R_{\mathrm{Body}}^{\mathrm{Cam}}\left(\phi ,\omega ,\kappa \right)R_{J2000}^{\mathrm{Body}}{R_{\mathrm{WGS}84}^{J2000}\left[\begin{array}{c}{X}_g-X_{\mathrm{gps}}\\ Y_g-Y_{\mathrm{gps}}\\ Z_g-Z_{\mathrm{gps}}\end{array}\right]}_{\mathrm{WGS}84}& \left(17\right)\end{array}$

[0058] It should be noted that, the models of different CCDs are determined by the respective viewing angle models in the adjustment model.

[0059] Step 6.2: On the basis of the external calibration, the above estimated external calibration parameters are considered as true values in internal calibration. For each CCD, its internal calibration parameters are estimated separately using the dense corresponding imagery points on its overlapping images. Since the constant term (a.sub.0, b.sub.0) is independent of the coplanar constraints between the images, the translation error of the constant term parameters is allocated to the external calibration for compensation, rather than in the internal calibration. Therefore, only the other higher order parameters are calculated in the CCD-by-CCD absolute internal calibration. According to the constructed adjustment model, the error equation can be constructed for each pair of corresponding imagery point by linearizing the model:

[00014]$\begin{array}{cc}\{\begin{array}{c}{v}_{\mathrm{IH}}^{\mathrm{Fi}}=B_i^{F}y+C_i^F{t}_i-R_i^F\\ v_{\mathrm{IH}}^{\mathrm{Bi}}=B_i^{B}y+C_i^B{t}_i-R_i^B\end{array}& \left(18\right)\end{array}$

[0060] In which, v.sub.IH.sup.Fi and v.sub.IG.sup.Bi are the correction vectors corresponding to the imagery points on the T-FWD and T-BWD images, respectively, y=[da.sub.1, da.sub.2, da.sub.3, db.sub.1, db.sub.2, db.sub.3].sup.T is the correction vectors of the calibration parameters for a CCD; t.sub.i=[dLat,dLon].sub.i.sup.T is the correction vectors of the ground plane coordinates of each corresponding imagery point, and its ground elevation is interpolated from the DEM of the image coverage area, the ground flat area should be selected for internal calibration as far as possible, so the open source DEM could usually meet the calibration accuracy requirements; the matrices and B.sub.i.sup.B are the partial derivative coefficient matrices on the calibration parameters in the error equations of the T-FWD and T-BWD image points, respectively; the matrices C.sub.i.sup.F and C.sub.i.sup.B are the partial derivative coefficient matrices on the ground plane coordinates in the error equation for the T-FWD and T-BWD imagery points, respectively; R.sub.i.sup.F and R.sub.i.sup.B are the constant vectors in the error equation for the T-FWD and T-BWD imagery points, respectively. Taking the image point on the T-FWD image as an example, the specific form of each matrix in the error equation is as follows:

[00015]$\begin{array}{cc}{B}_i^F={\left[\begin{array}{cccccc}\frac{\partial F_x^F}{\partial a_1}& \frac{\partial F_x^F}{\partial a_2}& \frac{\partial F_x^F}{\partial a_3}& 0& 0& 0\\ 0& 0& 0& \frac{\partial F_y^F}{\partial b_1}& \frac{\partial F_y^F}{\partial b_2}& \frac{\partial F_y^F}{\partial b_3}\end{array}\right]}_i{C}_i^F={\left[\begin{array}{cc}\frac{\partial F_x^F}{\partial \mathrm{Lat}}& \frac{\partial F_x^F}{\partial \mathrm{Lon}}\\ \frac{\partial F_y^F}{\partial \mathrm{Lat}}& \frac{\partial F_y^F}{\partial \mathrm{Lon}}\end{array}\right]}_i{R}_i^F={\left[\begin{array}{c}-F_x^F\\ -F_y^F\end{array}\right]}_i& \left(19\right)\end{array}$

[0061] Therefore, the estimation of y can be deduced according to the least squares adjustment:

[00016]$\begin{array}{cc}y={\left(\underset{i=1}{\overset{k}{.Math.}}\left(B_i^T{B}_i-B_i^T{C_i\left(C_i^T{C}_i\right)}^{-1}{C}_i^T{B}_i\right)\right)}^{-1}\left(\underset{i=1}{\overset{k}{.Math.}}\left(B_i^T{R}_i-B_i^T{C_i\left(C_i^T{C}_i\right)}^{-1}{C}_i^T{R}_i\right)\right)& \left(20\right)\end{array}$

[0062] In which, k is the number of corresponding imagery points on the overlapping image of this CCD, and

[00017]$B_i=\left[\begin{array}{c}{B}_i^F\\ B_i^B\end{array}\right],C_i=\left[\begin{array}{c}{C}_i^F\\ C_i^B\end{array}\right],R_i=\left[\begin{array}{c}{R}_i^F\\ R_i^B\end{array}\right],$

[0063] Similarly, the internal calibration solution is also an iteration process, and the current internal calibration parameters are updated according to the estimated corrections in each iteration, and be used as the inputs in the next iteration. The iterative solution ends when the estimated corrections in two consecutive iterations is less than a threshold. Given that the process of solving the high order internal calibration parameters is the same for each CCD, so the solutions for other CCDs are not repeated here.

[0064] Step 6.3: To ensure the geometric splicing accuracy among the CCDs, it is necessary to register the constant terms of all CCD under the same external calibration parameters. The corresponding imagery points between the overlapping images of adjacent CCDs are used as the observations to perform relative internal calibration. A CCD is chosen as the reference CCD, and its constant term is regarded as true values and not calculated. The constant terms of all non-reference CCDs are estimated with the reference CCD as the benchmark.

[0065] Using the above constructed adjustment model to construct the error equation, the specific form of this equation is similar to Eq. (18), the difference is that the parameters to be solved here are the constant terms of the viewing angle model of all non-reference CCDs, as follows:

[00018]$\begin{array}{cc}\{\begin{array}{c}{v}_{\mathrm{IL}}^{\mathrm{Fi}}=D_i^{F}z+E_i^F{t}_i-H_i^F\\ v_{\mathrm{IL}}^B=D_i^{B}z+E_i^B{t}_i-H_i^B\end{array}& \left(21\right)\end{array}$

[0066] In which, z=[dz.sub.1, dz.sub.2 . . . dz.sub.n].sup.T is the corrections of constant term for all non-reference CCDs, where dz=[da.sub.0, db.sub.0].sub.i and n is the number of non-reference CCDs; t.sub.i is still the correction vector of the ground plane coordinates of each corresponding imagery point; and Ds are the partial derivative matrices on the constant term of calibration parameters in the corresponding error equations, respectively; E.sub.i.sup.F and E.sub.i.sup.B are the partial derivative matrices on ground plane coordinates in the corresponding error equations, respectively; H.sub.i.sup.F and H.sub.i.sup.B are the constant vectors in the corresponding error equations, respectively.

[0067] Finally, the estimation of z can be deduced according to the least squares adjustment:

[00019]$\begin{array}{cc}z={\left(\underset{i=1}{\overset{k}{.Math.}}\left(D_i^T{D}_i-D_i^T{E_i\left(E_i^T{E}_i\right)}^{-1}{E}_i^T{D}_i\right)\right)}^{-1}\left(\underset{i=1}{\overset{k}{.Math.}}\left(D_i^T{H}_i-D_i^T{E_i\left(E_i^T{E}_i\right)}^{-1}{E}_i^T{H}_i\right)\right)& \left(22\right)\end{array}$

[0068] In which, λ represents the number of corresponding imagery points used for the calculation, and

[00020]$D_i=\left[\begin{array}{c}{D}_i^F\\ D_i^B\end{array}\right],E=\left[\begin{array}{c}{E}_i^F\\ E_i^B\end{array}\right],H_i=\left[\begin{array}{c}{H}_i^F\\ H_i^B\end{array}\right],$

[0069] The solution of relative internal calibration is also an iteration process, and the current calibration parameters are updated according to the estimated corrections in each iteration, and be used as the inputs in next iteration. The iterative solution ends when the estimated corrections in two consecutive iterations is less than a threshold.

[0070] In summary, the method achieves fully autonomous geometric calibration of the LARSS aiming at correcting the systematic geometric errors in the imaging model, based on the joint observation for stars and earth by satellites, under the condition without ground calibration sites. The external calibration is achieved using the star observations to determine the exact viewing axis direction of the camera. The internal calibration is achieved through the ground overlapping images to realize the compensation of absolute and relative geometric distortion within the camera. Therefore, this method is suitable for satellites with various designs such as multi-CCDs splicing and multi-bands registration and so on.