A geometric positioning processing method based on a linear array push-broom imaging satellite

Abstract

The application belongs to the technical field of satellite geometric positioning, and specifically discloses a geometric positioning processing method based on a linear array push-broom imaging satellite, which comprises the following steps: constructing a bias matrix and a pointing angle polynomial model, the bias matrix being used for correcting the pointing of the imaging light of a detector element, and the pointing angle polynomial model being used for simulating the process that the pointing angle changes with the column of the detector element image; the pointing angle is the angle between the imaging light of the detector element and the direction perpendicular to the orbit direction and the direction along the orbit; based on the bias matrix and the pointing angle polynomial model, a rigorous geometric positioning model is adjusted to obtain a joint calibration model, and the rigorous geometric positioning model is used for representing the geometric relationship between an image point and a ground point. Through the application, the interior and exterior orientation elements are checked, so that the rigorous geometric positioning model after adjustment can more accurately represent the geometric relationship between the image point and the ground point, and the geometric positioning precision of the linear array push-broom imaging satellite is improved.

Description

TECHNICAL FIELD

[0001] The present application belongs to the technical field of satellite geometric positioning, and more particularly relates to a geometric positioning processing method based on a linear array push-broom imaging satellite. BACKGROUND

[0002] Satellite in-orbit geometric calibration is a key link to fully play the geometric performance of the satellite. The geometric positioning model can provide geographic position information for satellite images, and is a key for subsequent measurement, positioning and other applications of remote sensing images. A general rigorous geometric positioning model can be constructed for a linear array push-broom imaging satellite, and the geometric relationship between an image point and a real ground point can be established through the model. However, the construction process of the positioning model involves the combination of many internal and external orientation element information, and the errors of the external orientation elements and the internal orientation elements will affect the positioning accuracy of the geometric positioning model. How to improve the geometric positioning accuracy based on the linear array push-broom imaging satellite is a technical problem to be solved in the field. SUMMARY

[0003] In view of the defects of the prior art, the purpose of the present application is to improve the geometric positioning accuracy based on the linear array push-broom imaging satellite.

[0004] To achieve the above purpose, in a first aspect, the present application provides a geometric positioning processing method based on a linear array push-broom imaging satellite, which comprises:

[0005] constructing a bias matrix and a pointing angle polynomial model, the bias matrix being used for correcting the pointing of the imaging light of the detector element, and the pointing angle polynomial model being used for simulating the process that the pointing angle changes with the image column of the detector element, the pointing angle being the angle of the imaging light of the detector element relative to the direction perpendicular to the orbit direction and the direction along the orbit direction;

[0006] adjusting the rigorous geometric positioning model based on the bias matrix and the pointing angle polynomial model, and obtaining a joint calibration model, the rigorous geometric positioning model being used for representing the geometric relationship between the image point and the ground point.

[0007] In a second aspect, the present application provides a geometric positioning processing device based on a linear array push-broom imaging satellite, which comprises:

[0008] a construction module, which is used for constructing a bias matrix and a pointing angle polynomial model, the bias matrix being used for correcting the pointing of the imaging light of the detector element, and the pointing angle polynomial model being used for simulating the process that the pointing angle changes with the image column of the detector element, the pointing angle being the angle of the imaging light of the detector element relative to the direction perpendicular to the orbit direction and the direction along the orbit direction;

[0009] an obtaining module, which is used for adjusting the rigorous geometric positioning model based on the bias matrix and the pointing angle polynomial model, and obtaining a joint calibration model, the rigorous geometric positioning model being used for representing the geometric relationship between the image point and the ground point.

[0010] Thirdly, this application provides an electronic device, comprising: at least one memory for storing a program; and at least one processor for executing the program stored in the memory, wherein when the program stored in the memory is executed, the processor is configured to execute the method described in the first aspect or any possible implementation thereof.

[0011] Fourthly, this application provides a computer-readable storage medium storing a computer program that, when run on a processor, causes the processor to perform the method described in the first aspect or any possible implementation thereof.

[0012] Fifthly, this application provides a computer program product that, when run on a processor, causes the processor to perform the method described in the first aspect or any possible implementation thereof.

[0013] Overall, the technical solutions conceived in this application have the following beneficial effects compared with the prior art:

[0014] (1) By adjusting the rigorous geometric positioning model based on the bias matrix and the pointing angle polynomial model, the adjusted rigorous geometric positioning model can be used as a joint calibration model. The bias matrix in the joint calibration model can be used to correct the pointing of the imaging rays of the probe to realize the calibration of the external orientation elements. The pointing angle polynomial model in the joint calibration model can be used to simulate the process of the pointing angle changing with the image column of the probe to realize the calibration of the internal orientation elements. By calibrating the internal and external orientation elements, the adjusted rigorous geometric positioning model can more accurately represent the geometric relationship between the image points and the ground points, thereby improving the geometric positioning accuracy of linear array pushbroom imaging satellites.

[0015] (2) Construct a corresponding geometric calibration model based on the error sources and error characteristics of positioning accuracy, and obtain control point data. The correction amount, by refining the geometric positioning model, can effectively improve the geometric positioning accuracy of the image. While improving positioning accuracy, the consistency accuracy between multispectral images is also considered. Joint calibration is performed using high-precision corresponding points between panchromatic and multispectral data to obtain... The correction amount is used to modify the geometric positioning model and ensure the positioning consistency of the geometric positioning model for images in each band.

[0016] (3) For satellites with obvious attitude jitter, commonly used attitude models such as polynomial models are difficult to simulate accurately. However, the on-board attitude sampling frequency is high, and the penalized spline function model has a good fitting effect. It can accurately simulate the on-board attitude, eliminate the intersection error of corresponding points with obvious trends between spectral bands, and effectively ensure the positioning accuracy of the rigorous geometric positioning model. Attached Figure Description

[0017] Figure 1 This is a flowchart illustrating the geometric positioning processing method based on a linear array pushbroom imaging satellite provided in this application embodiment;

[0018] Figure 2 This is a schematic diagram of linear array push-broom imaging provided in an embodiment of this application;

[0019] Figure 3 This is a schematic diagram of CCD rotation error provided in an embodiment of this application;

[0020] Figure 4 This is a schematic diagram of the camera probe position provided in an embodiment of this application;

[0021] Figure 5 This is a schematic diagram of satellite attitude data information provided in an embodiment of this application;

[0022] Figure 6 This is a scatter plot of the dynamic changes in satellite attitude angles provided in the embodiments of this application;

[0023] Figure 7 This is a polynomial model attitude fitting diagram provided in the embodiments of this application;

[0024] Figure 8 This is a schematic diagram of the positioning error of the corresponding point provided in the embodiment of this application. Detailed Implementation

[0025] To make the objectives, technical solutions, and advantages of this application clearer, the following detailed description is provided in conjunction with the accompanying drawings and embodiments. It should be understood that the specific embodiments described herein are merely illustrative and not intended to limit the scope of this application.

[0026] In the embodiments of this application, the terms "exemplary" or "for example" are used to indicate that something is an example, illustration, or description. Any embodiment or design that is described as "exemplary" or "for example" in the embodiments of this application should not be construed as being more preferred or advantageous than other embodiments or design. Specifically, the use of the terms "exemplary" or "for example" is intended to present the relevant concepts in a specific manner.

[0027] In the description of the embodiments of this application, unless otherwise stated, "multiple" means two or more, for example, multiple processing units means two or more processing units, multiple elements means two or more elements, etc.

[0028] The embodiments of this application are described below with reference to the accompanying drawings.

[0029] Figure 1This is a flowchart illustrating the geometric positioning processing method based on linear array pushbroom imaging satellites provided in this application embodiment, as shown below. Figure 1 As shown, the method includes the following steps S101 and S102.

[0030] Step S101: Construct the bias matrix and the pointing angle polynomial model. The bias matrix is ​​used to correct the pointing of the imaging ray of the probe. The pointing angle polynomial model is used to simulate the process of the pointing angle changing with the probe image column. The pointing angle is the angle of the probe's imaging ray relative to the direction perpendicular to the track and along the track.

[0031] Step S102: Based on the bias matrix and pointing angle polynomial model, adjust the rigorous geometric positioning model to obtain the joint calibration model. The rigorous geometric positioning model is used to characterize the geometric relationship between image points and ground points.

[0032] Understandably, by adjusting the rigorous geometric positioning model based on the bias matrix and the pointing angle polynomial model, the adjusted rigorous geometric positioning model can be used as a joint calibration model. The bias matrix in the joint calibration model can be used to correct the pointing of the imaging rays of the probe to calibrate the exterior orientation elements. The pointing angle polynomial model in the joint calibration model can be used to simulate the process of the pointing angle changing with the probe image column to calibrate the interior orientation elements. By calibrating the interior and exterior orientation elements, the adjusted rigorous geometric positioning model can more accurately represent the geometric relationship between image points and ground points, thereby improving the geometric positioning accuracy of linear pushbroom imaging satellites.

[0033] The geometric positioning processing method based on linear array pushbroom imaging satellites provided in this application is illustrated below with several examples.

[0034] (1) Construction of positioning model and analysis of positioning error of linear array pushbroom imaging camera.

[0035] Linear pushbroom imaging satellites employ linear pushbroom imaging cameras. A linear pushbroom imaging camera is a device that uses a linear image sensor for imaging, capturing image information through line-by-line scanning. During imaging, the photosensitive elements on the sensor sequentially expose the scene, converting light signals into electrical signals to form one-dimensional image data. As the camera or the scanned object moves, the accumulated data line by line is eventually combined to form a complete two-dimensional image.

[0036] Geometric positioning models provide geographic location information for satellite imagery, which is crucial for subsequent measurement and positioning applications of remote sensing imagery. The function of a geometric positioning model is to establish the mapping relationship between the coordinates of points on the image and their actual geographic coordinates. During the imaging process, each pixel on the sensor can be simplified into a series of points; therefore, geometric positioning models corresponding to different satellite imagery can be established based on the structural characteristics of different sensors.

[0037] As the most accurate positioning model, the rigorous geometric model is constructed strictly according to the geometric process of satellite imaging. This process requires obtaining the real parameter information of the satellite's sensors when it is in orbit, and it also requires the satellite's real-time orbit, attitude and other information to establish the relationship between coordinate systems.

[0038] (1-1) Construction of geometric positioning model of linear array pushbroom imaging camera.

[0039] (1-1-1) Principle of rigorous (imaging) geometric model;

[0040] For satellites employing a linear pushbroom imaging mode, the received electromagnetic wave signals from ground objects are converted into electrical signals by a photosensor, recorded, and stored. These electrical signals are then converted into digital signals and transmitted to the ground receiving station via a data transmission system. The linear pushbroom imaging camera captures images line by line as the satellite moves, as illustrated in the diagram below. Figure 2 As shown.

[0041] During satellite operation, remote sensing images are stitched together from line images as the satellite moves. At the moment of each line image imaging, the imaging payload CCD (Charge-Coupled Device) probe, the optical perspective center, and the target location are collinear. Therefore, a rigorous geometric model can be established based on the collinearity equation using the internal orientation elements of the satellite sensor and the actual external orientation elements measured by the satellite platform during operation to achieve precise positioning of remote sensing images.

[0042] (1-1-2) Construction of a rigorous geometric model;

[0043] The general rigorous geometric positioning model is as follows:

[0044] (1);

[0045] In the formula, The transformation matrix from the camera coordinate system to the body coordinate system was obtained through measurement before satellite launch. Represents the camera's principal distance. The field angle of the CCD array along the track. The position corresponding to the main viewing axis (the point perpendicular to the CCD linear array through the principal point). Indicates the position of each pixel in the CCD array, subscript Represents the index of a pixel. For the size of the probe, It is the offset of the camera's origin in the body coordinate system. It is the transformation matrix for converting the body coordinate system to the geocentric inertial coordinate system. The J2000 (Julian 2000) coordinate system is the internationally used geocentric inertial coordinate system. In this application, the attitude quaternion data measured by the satellite using the star sensor is the attitude of the satellite in the body coordinate system relative to the attitude in the J2000 coordinate system. The matrix can be determined by using the quaternions parsed from the satellite auxiliary data. Let t be the transformation matrix for converting the J2000 coordinate system to the WGS84 (World Geodetic System 1984) coordinate system at time t; This refers to the coordinates of the origin of the GPS (Global Positioning System) module in the body coordinate system. Let t be the position vector of the phase center of the GPS module in the WGS84 coordinate system at time t; This indicates the coordinates of a pixel in the satellite image in the WGS84 coordinate system, where m is the scaling factor.

[0046] The origin of the camera coordinate system is located at the center of the camera projection. The Z-axis of the camera coordinate system is the principal optical axis of the camera and points positively toward the focal plane. The Y-axis of the camera coordinate system is parallel to the direction of the CCD array, and the X-axis of the camera coordinate system roughly points toward the direction of satellite flight. The three axes are aligned according to the right-hand coordinate system rule.

[0047] The satellite's body coordinate system is a coordinate system fixed to the satellite. Typically, the satellite's center of mass is taken as the origin O, and the three principal axes of inertia are taken as the X, Y, and Z axes. The OZ axis is positive, pointing from the center of mass to the ground; the OX axis is positive, pointing in the direction of the satellite's flight; and the OY axis is determined by the right-hand coordinate system rules.

[0048] This formula is a general and rigorous geometric positioning model for linear array pushbroom imaging satellites. Through this model, the geometric relationship between image points and real ground points can be established, which is the basis for geometric processing of remote sensing data.

[0049] (1-2) Positioning model error analysis;

[0050] The construction of the localization model involves combining numerous interior and exterior orientation elements. Analyzing the sources of error throughout the imaging process is crucial for improving the localization accuracy of the image and the consistency accuracy across multiple spectral bands. This also forms the basis for subsequent localization model optimization. By analyzing the process of constructing the geometric model, the errors can be subdivided into interior orientation element errors, exterior orientation element errors, etc.

[0051] (1-2-1) External orientation element error analysis;

[0052] In light of the characteristics of linear array pushbroom imaging satellites, the exterior orientation elements refer to the three linear elements used to describe the spatial coordinates of the imaging center and the three angular elements used to express the spatial attitude at the instant a line of images is formed. The position of the imaging center can be indirectly obtained through the GPS module carried on the satellite, while the attitude of the satellite image can be obtained by measuring data from the star sensor.

[0053] The orbital position error, i.e. the position error of the imaging center, can be decomposed into three directions: along the orbit, perpendicular to the orbit, and radial. According to the characteristics of the linear array pushbroom imaging mode, the errors along the orbit and perpendicular to the orbit are translation errors, while the image point offset caused by the radial error is a translation error and a scaling error. Based on the on-board high-precision orbit determination system, the scaling error caused by the radial error is relatively small. Therefore, all three types of errors can be equivalent to translation errors.

[0054] The errors of the onboard attitude system can also be classified into roll angle error, pitch angle error and yaw angle error according to the direction of rotation. The roll angle error is similar to the pitch angle error and can be equivalent to the translation error. The error caused by the yaw angle can be equivalent to the CCD rotation error.

[0055] From satellite launch to in-orbit operation, errors can occur in camera installation, star sensors, and other components due to external environmental influences. This is the main source of equipment installation errors. Translational errors caused by displacement during this process can be equivalent to orbital deviation errors, and installation angle errors can be equivalent to attitude measurement errors.

[0056] (1-2-2) Internal Orientation Element Error Analysis;

[0057] The internal orientation element error mainly includes principal point offset error, satellite camera principal distance error, detector size error, CCD rotation error, and error caused by lens distortion. The internal orientation element error reflects the degree of internal distortion of the satellite.

[0058] The principal point offset error can be equivalent to the translation error along the track and perpendicular to the track. The principal distance error will cause the image point coordinates to be scaled, which can be equivalent to the proportional error along the track and perpendicular to the track. The image point offset caused by the probe size error can be equivalent to the proportional error.

[0059] like Figure 3 As shown, during the satellite's operation in orbit, external forces cause the CCD's placement to be non-perpendicular to the satellite's flight direction. The resulting offsets in the x and y directions are:

[0060] (2);

[0061] In the formula: For rotation angle, Center of rotation This represents the y-coordinate of a pixel.

[0062] Errors caused by deformation of the camera lens surface or displacement of its position are collectively referred to as errors caused by lens distortion. This distortion is mainly divided into radial distortion and eccentric distortion. Radial distortion is caused by the surface curvature error of the lens element, while eccentric distortion refers to the distortion caused by the lenses not being kept on a level plane. The errors in the two directions caused by radial distortion are:

[0063] (3);

[0064] in: This indicates radial distortion. , and For polynomials with fixed coefficients, For variables in a polynomial, Like a dot Coordinates in the camera coordinate system , Represents the principal distance of the image point (camera principal distance).

[0065] The error caused by eccentric distortion is:

[0066] (4);

[0067] in, , and represents the coefficients of the equation to be solved.

[0068] (2) Verification of internal and external orientation elements.

[0069] (2-1) Exterior orientation element check;

[0070] The purpose of exterior orientation element calibration is to eliminate systematic errors in data such as orbit and attitude. Orbit measurement errors can be treated as equivalent to attitude measurement errors, and load installation errors are equivalent to attitude measurement errors. Therefore, in the exterior orientation element calibration process, these three can be unified by introducing an offset matrix to correct the imaging rays during the imaging process, ensuring that the rays point to their true positions, thus eliminating systematic errors. The offset matrix is ​​an orthogonal rotation matrix formed by the angle between the actual observation direction and the theoretical observation direction of a pixel in the body coordinate system. Represents the bias matrix. The angles by which the imaging rays rotate around the x-axis, y-axis, and z-axis of the body coordinate system represent these rotations. , , , As shown below:

[0071] (5);

[0072] To obtain the true direction of light rays in the body coordinate system, in practical applications, remote sensing satellite imagery and control imagery data are matched to obtain control points. The geographic coordinates of the control points are then used to calculate the true direction of light rays in the body coordinate system. The deviation between the direction of light rays from the image point in the body coordinate system and the ray line calculated from the control point is calculated, and the direction of light rays is corrected using an offset matrix. According to formula (1), the corrected form after adding the offset matrix is:

[0073] (6);

[0074] Formula (6) is the external orientation element compensation model.

[0075] (2-2) Internal orientation element check;

[0076] Satellites can use three CCDs to stitch together images. The internal orientation errors introduced by multiple CCDs mainly include translation and rotation errors caused by the positional deviation of the CCD devices. They also include internal orientation element errors such as principal point offset error, principal distance error, and probe size error. These errors are highly correlated, and some errors cannot be modeled and eliminated. Therefore, it is difficult to establish a distortion model to achieve internal orientation element calibration by using the distortion caused by various errors.

[0077] In fact, the combined error caused by the interior orientation elements can be regarded as the influence on the true ray direction of each imaging element, which causes the offset of the true position of the image point. In practical applications, high-precision corresponding points of the control image and the satellite can be obtained as ground control points, and the true orientation of the pixel can be restored by using the control points.

[0078] The relationship between the camera sensor and the imaging center in the camera coordinate system is as follows: Figure 4 As shown, the imaging beam of the probe... Decomposed into directions perpendicular to the orbit and along the orbit, the angles of the imaging ray relative to these two directions are respectively... Then it's like a point. Coordinates in the camera coordinate system These two angles can be represented as: , This represents the principal distance of the image points (camera principal distance). However, due to the presence of coarse errors in the matching process, using the pointing angle model to calibrate each pixel cannot account for local smoothness features. This leads to significant errors in the pointing angle of the detector obtained using the coarse error solution. Furthermore, the number of control points required to independently solve the pointing angle for each pixel is too large. Therefore, a polynomial is used to approximate the process of the pointing angle changing with the detector image series (which can be represented by the x-coordinate of the image points). Typically, a polynomial with a maximum degree of 5 can adequately simulate the errors of various interior orientation elements. The pointing angle simulated using the polynomial method can be expressed as:

[0079] (7);

[0080] in, Represents the x-coordinate of an image point;

[0081] The geometric model can now be rewritten as follows:

[0082] (8);

[0083] Combining the aforementioned internal and external orientation element models, the satellite joint calibration model is as follows:

[0084] (9).

[0085] (3) Analysis of satellite imaging geometric parameter characteristics and model construction;

[0086] The construction of geometric positioning models for different satellites requires knowledge of the satellite's attitude and orbit information during imaging. The imaging mode of linear array pushbroom imaging satellites necessitates acquiring the time, attitude, and orbit corresponding to each row of images. Due to the limited measurement frequency of onboard attitude and orbit measurement equipment, the orbit and attitude data acquired onboard are discrete, and generally, it is not possible to directly obtain the data information corresponding to each row. Therefore, it is necessary to perform error analysis on the attitude, orbit, and timing characteristics of different satellites, and use the most suitable interpolation model based on its data characteristics. This allows the model to effectively recover the exterior orientation element information of the satellite during imaging, thereby improving the positioning accuracy of the rigorous geometric model.

[0087] This application analyzes the characteristics of the corresponding imaging geometric parameters of satellite orbit, travel time, and attitude auxiliary data, and constructs corresponding orbit, travel time, and attitude models.

[0088] (3-1) Track characteristic analysis and model construction

[0089] To provide a more intuitive analysis of the actual data, this application uses real data transmitted from the satellite for analysis. The satellite uses the Global Positioning System to record its position and velocity information in the WGS-84 system at a specific moment.

[0090] Based on the on-board auxiliary data format, the corresponding orbital data information can be parsed. By analyzing the variation patterns of the original orbital data of the calibration scene image with time as the variable, a scatter plot of the X, Y, and Z positions in the WGS-84 system can be obtained from the original on-board data, along with the X-velocity in the WGS-84 system. Y speed Z-speed The scatter plot shows that the satellite's orbital position remains relatively stable over time without significant fluctuations, and its velocity also tends to stabilize. This indicates that the satellite's orbital position is relatively stable during its operation. Based on the smooth characteristics of the scatter plot, a polynomial model can be used to fit the orbital data, thus providing the orbital information at any given time.

[0091] The orbital polynomial fitting model is shown below:

[0092] (10);

[0093] In the formula: Represents the coefficients of each polynomial. Represents the moment of imaging.

[0094] (3-2) Streaming time characteristic analysis and model construction;

[0095] The satellite records the travel time information and interpolates the travel time on the satellite before transmitting it to the ground.

[0096] The on-board time data obtained from the parsing of the original on-board data is line-by-line time, so the image behavior can be used as an independent variable to analyze the change of the imaging time of each line with the line number.

[0097] The satellite's in-flight imaging time follows the in-flight number linearly, and the in-flight integration time is within the range of 0.000148 seconds to 0.000151 seconds. The in-flight integration time is stable, with the largest jump in integration time being less than 0.000003 seconds. Linear interpolation methods can be used to model the in-flight time.

[0098] (3-3) Attitude characteristics analysis and model construction.

[0099] (3-3-1) Attitude characteristics analysis;

[0100] The satellite's attitude defines the rotation matrix between the Earth-fixed coordinate system and the orbital coordinate system. After extended Kalman filtering, it is represented by a quaternion q in the J2000 coordinate system. During its on-orbit operation, the satellite calculates the quaternion values ​​using onboard star sensors and transmits them to the auxiliary data. Therefore, the initial auxiliary data, after parsing, yields quaternions, as detailed below. Figure 5 As shown. However, quaternions do not actually have geometric meaning, so it is impossible to analyze the satellite's attitude by simply observing the changes in quaternions. In order to better analyze the changes in the satellite's attitude, quaternions are converted to Euler angles with geometric meaning, and the changes in the satellite's attitude are analyzed by observing the Euler angles.

[0101] For satellites exhibiting high-frequency attitude jitter (high-frequency attitude jitter type satellites), the changes in onboard attitude can also be analyzed by observing Euler angles. The converted pitch, roll, and yaw angles are analyzed with sampling time as the independent variable, yielding results such as... Figure 6 The results are shown.

[0102] Figure 6 Subgraphs a, b, and c Figure Three The three plots represent scatter plots showing the changes in pitch, roll, and yaw angles over time. All three plots show a clear overall upward or downward trend, and the changes in the three angles do not exhibit a smooth characteristic but rather show significant fluctuations, especially in the roll and yaw angles, where attitude jitter is very pronounced. Among the three angles, the yaw angle has the highest oscillation frequency, approximately 0.25 Hz.

[0103] Commonly used pose models include polynomial models based on Euler angles. For example... Figure 7 The image shown is a graph illustrating the simulation results of the polynomial model on the original data. From... Figure 7 As can be seen, the third-order and fourth-order polynomials do not fit the original data well and do not reflect the obvious oscillation effect in the original data. Therefore, it is not advisable to use polynomials as the basis of the attitude model.

[0104] analyze Figure 6 It can be observed that the yaw angle has the highest oscillation frequency, approximately 0.25 Hz, while the onboard attitude sampling frequency is 4 Hz. If a correct attitude model is selected, the attitude oscillation process can be simulated relatively accurately. Therefore, this application proposes a penalized spline function model as the attitude model to simulate the onboard attitude for satellites with high-frequency attitude jitter.

[0105] (3-3-2) Penalized spline function model;

[0106] For attitude data with high oscillation frequency and large oscillation amplitude, common attitude models cannot simulate it accurately. In order to better solve this problem, spline function models can be used to solve related issues.

[0107] Spline functions are composed of polynomials over subspaces with continuity conditions. Assume there are n+1 points. These points satisfy These points are called nodes. For a node... Specify an integer On these nodes Subspline function The following two conditions must be met:

[0108] a) In each interval superior, It is a number of times polynomials;

[0109] b) In superior have The continuous derivative of order 1.

[0110] Common spline functions include cubic splines, which assume the existence of points The cubic spline function is actually... Spline functions that pass through these points and satisfy the above conditions.

[0111] Based on the conditions of cubic spline functions, the following conditions can be derived:

[0112] (11);

[0113] From the definition of spline functions, we know that in each small interval... To determine the spline function, we need to identify 4 coefficients across n intervals. Based on the continuity condition of the second derivative, we can determine 3n-3 conditions. Since the cubic spline function passes through all points, we can determine n+1 conditions, for a total of 4n-2 conditions. Therefore, we need 2 more conditions to determine the spline function. Typically, it will be in the range. Two conditions are added to the boundary to calculate all undetermined coefficient values. However, due to noise in the measurement data, the attitude jitter frequency is high, and overfitting may occur when using cubic spline functions to simulate the attitude data.

[0114] Considering the presence of noise, we hope to use a loss function to better balance the relationship between noise and pose fitting accuracy. The form of the loss function is as follows:

[0115] (12);

[0116] In the formula The basis matrix represents the spline function. Here is the parameter matrix of the spline function model. represents the undetermined coefficients for determining the spline function, and 'a' represents the attitude recording value. This is a penalty factor. Let be a symmetric penalized positive definite matrix. In the formula... The penalty method of Eilers and Marx is used to determine the matrix by calculating the second derivative or difference of the basis functions, and constructing a symmetric penalty positive definite matrix. For example, if spline functions are used as basis functions, they can be constructed using the second derivative of the splines. Simultaneously, B-spline functions are chosen as the basis functions of the spline functions to determine the above equation. And the penalty factor It plays an important role in balancing noise control and attitude fitting accuracy.

[0117] Wahba proposed a generalized cross-validation method to estimate... , The estimated value is The minimum value, The definition of is:

[0118] (13);

[0119] in, It is an (n+1)×(n+1) matrix. The definition is as follows:

[0120] (14);

[0121] This method allows us to determine the penalty factor, thereby completing the construction of the penalty spline function.

[0122] Understandably, the above example analyzes the attitude changes of satellites experiencing high-frequency attitude jitter during on-orbit operation. Based on onboard data, significant attitude jitter was observed. Considering the relatively high stability of onboard attitude measurement components, random errors caused by measurement issues were ruled out. For attitudes with significant jitter, commonly used attitude models such as polynomial models are difficult to simulate accurately. However, the high sampling frequency of onboard attitude allows for the simulation of onboard attitude using precise models. This application uses a penalized spline function as the basis for establishing the attitude model. Real data transmitted from the satellite was used to compare different attitude models with the penalized spline model, and the accuracy of the attitude model was reflected by the multi-spectral band registration accuracy. Experiments demonstrate that the penalized spline function model has a good fitting effect, eliminating intersection errors of corresponding points with significant trends between spectral bands, thus ensuring the positioning accuracy of the satellite's product data.

[0123] (4) Solving the joint calibration model;

[0124] During joint calibration, the obtained matching data includes control points matched between panchromatic and multispectral images and control images (reference images used in Geographic Information Systems (GIS) or image registration); corresponding points matched between multispectral and panchromatic images; corresponding points between different CCDs in panchromatic images; and corresponding points between different CCDs in multispectral images. Panchromatic images (images captured over a wide band, typically covering the entire visible light range) have higher resolution and can match high-precision control points, thus achieving better calibration results compared to multispectral images. However, the imaging time interval between multispectral and panchromatic images is short, allowing for high-precision matching of corresponding points. Therefore, adjusting the interior orientation elements of the multispectral image using corresponding points between panchromatic and multispectral images yields a more significant effect. For example, a satellite can generate five images of the same area within a short time interval, including a panchromatic image, red band image B1, green band image B2, blue band image B3, and infrared band image B4.

[0125] (4-1) Control point equations (standard equations fitted by correct, error-free control points);

[0126] The joint calibration model utilizes bias matrices and a polynomial model to eliminate errors in the interior and exterior orientation elements. The equations are constructed using control data from matching panchromatic and multispectral imagery with control imagery:

[0127] (15);

[0128] In the formula: These represent the three-dimensional coordinates of the control points matched with the control images in the panchromatic band, multispectral red band, and control image, respectively. , This represents the coefficients of the pointing angle polynomial model. The variable in this formula is... and polynomial coefficients. Let: .

[0129] Formula (15) can be simplified to:

[0130] (16);

[0131] in:

[0132] , (17);

[0133] According to formula (17), after eliminating m, let:

[0134] (18);

[0135] Linearizing the above equation, we get:

[0136] (19);

[0137] In the formula Represent The correction amount. The above process is the control point solution process in the joint calibration model.

[0138] (4-2) Equation of small intersection angle connection point (In the field of remote sensing, especially in satellite image processing, "small intersection angle" usually refers to the small angle between the lines of sight (LOS) of two sensors).

[0139] like Figure 8 As shown, if satellites A and B image the same ground feature, the image points will be located at... place, This represents the imaging location corresponding to satellite A. Let B be the imaging position corresponding to satellite B. At this point, the imaging equation between the two points can be expressed as:

[0140] (20);

[0141] Subscript and Used to distinguish Corresponding image points and The corresponding image point, and Indicates the pointing angle in the x and y directions.

[0142] If the imaging geometry of these two satellites perfectly matches the actual satellite configuration and the elevation of the ground point S is known, then the geographic coordinates of the two points calculated based on the imaging model should be the same. That is...

[0143] (twenty one);

[0144] Elevation is typically obtained using a Digital Elevation Model (DEM). However, due to the limited accuracy of DEMs, errors caused by elevation are difficult to avoid. Even if the imaging geometry of two satellites is completely correct, the geographic coordinates of the two points will still not be equal. The error caused by elevation can be approximated as:

[0145] (twenty two);

[0146] In the formula: and The attitude angles of the two images taken consecutively; This represents the elevation error. According to formula (22). Depend on as well as and Control, if it can be made and When the satellite's attitude angles are very close to those of the two images taken, the intersection error of corresponding points caused by elevation can be effectively eliminated. In this case, the positioning error of the corresponding points can directly reflect the interior and exterior orientation element errors of the two bands during imaging.

[0147] The panchromatic and multispectral imaging CCD elements of the satellite can be mounted on the same camera and the distance between CCDs of different bands is very small. When the satellite pushbrooms the image, the time interval between imaging the same ground object by different bands is short. It can be assumed that the attitude angles of the two images are very close. At this time, the intersection error of the corresponding points caused by the elevation error can be effectively controlled.

[0148] Because the imaging attitude angles of two points in different bands are similar, the Earth ellipsoid equation needs to be introduced when solving for the 3D coordinates of corresponding points. There are two common methods for solving for the 3D coordinates of ground points where corresponding points intersect using the ellipsoid equation: iterative transformation based on DEM; and stereo intersection method with small intersection angle.

[0149] (4-2-1) Iterative transformation based on DEM;

[0150] Based on the Earth ellipsoid model, the three-dimensional coordinates of any point on the ellipsoid surface are... The following relationship should be satisfied:

[0151] (twenty three);

[0152] In the formula: and These represent the major and minor axes of the ellipsoid, respectively. This represents the altitude of a point on the ground. For the WGS84 ellipsoid, we have... m, m. In the formula This can be obtained through DEM, thus establishing it. and , The relationship between them.

[0153] According to formula (20), it can be determined that , , and The relationship at this time

[0154] (twenty four);

[0155] In the formula Since all values ​​are known, the first expression of formula (20) can be simplified to:

[0156] (25);

[0157] Substituting this equation into the ellipsoid equation, we get

[0158] (26);

[0159] The only unknown in this formula is... Therefore, by solving the quadratic equation in one variable, we can obtain... At this time, Substituting back into formula (25), a three-dimensional coordinate can be obtained. .make = As the initial value for elevation, Substitute the values ​​into the DEM and interpolate to obtain the corresponding elevation values. In formula (32) Updated to Calculation yields new Calculate the new elevation values ​​again in the DEM. Repeat the preceding process iteratively until the two consecutive iterations are completed. If the difference is less than a certain threshold, then the calculated value is... That is, the three-dimensional coordinates corresponding to that point.

[0160] Calculate using the above methods respectively , Two points ( and The geographic coordinates corresponding to two image points are used to distinguish them. The average value is the geographic coordinate of the ground point corresponding to the two points.

[0161] (4-2-2) Small intersection angle solid intersection method;

[0162] Based on the principle of stereo intersection, when multiple satellite images converge on the ground with imaging rays from different angles, the three-dimensional coordinates of ground points can be calculated without relying on a DEM.

[0163] However, because the time interval between panchromatic and multispectral imaging is extremely short and they are located on the same camera, the attitude angles of the two images are very close when imaging corresponding points. This makes the elevation solution in the three-dimensional coordinates unstable, so it is also necessary to introduce the Earth ellipsoid model. The three-dimensional coordinates of corresponding points are solved by adjustment method using the intersection model.

[0164] Transform the first expression in formula (26) (the expression to the left of the first plus sign in formula (26) into:

[0165] (27);

[0166] Will , All are known values, let Then we have:

[0167] (28);

[0168] Substituting the above equation into formula (33), we get:

[0169] (29);

[0170] Eliminate according to the above formula Then establish the equation and let . and Only in China , , Since it is an unknown, and according to the equation of the ellipsoid:

[0171] (30);

[0172] At this point, the equation contains only two independent variables. , ,right and After expanding using a Taylor series and retaining the first-order terms, we obtain:

[0173] (31);

[0174] in:

[0175] ;

[0176] According to formula (30), the following can be calculated:

[0177] (32);

[0178] At this point, the error equation can be constructed according to formula (31) to solve for the three-dimensional coordinates of the corresponding points.

[0179] (4-2-3) Solving for the small intersection angle connection point of the joint calibration model;

[0180] Combining formula (26), the corresponding points of panchromatic and multispectral matching, the corresponding points between different panchromatic CCDs, and the corresponding points between different multispectral CCDs are combined. Similarly, the bias matrix and pointing angle polynomial model are used to eliminate the errors of interior and exterior orientation elements. Then the equation can be rewritten as:

[0181] (33);

[0182] In the formula This represents the corresponding points matched between the panchromatic band and a specific band in the multispectral spectrum. Since the imaging attitude angles of corresponding points in different bands are approximately the same, the intersection error caused by elevation is controlled within a certain range. The intersection error between two points reflects the exterior and interior orientation errors of the two images. Therefore, the geographic coordinates of the corresponding points in each band of the panchromatic and multispectral spectrum can be calculated as initial values. The three-dimensional coordinates of the corresponding points are not control data and inherently contain errors; therefore, this model requires solving... In addition to the polynomial coefficients, the geographic coordinates of the corresponding points also need to be determined. The amount of correction.

[0183] The corresponding points matched from the CCD overlapping areas of images in each band have the same characteristics as the corresponding points between spectral bands. The imaging attitude angles are approximately the same, and the error equations are constructed in the same way as those for corresponding points between spectral bands.

[0184] Equation (33) is transformed similarly to the control point equation, let:

[0185] (34);

[0186] In the formula:

[0187] , , … This represents the values ​​of the bias matrix.

[0188] Linearizing the above equation yields:

[0189] (35);

[0190] In the formula Represent The correction amount. This process is the solution process for the small intersection angle connection points in the joint calibration model.

[0191] (4-3) Solving for joint calibration parameters

[0192] The error equation is constructed as follows:

[0193] (36);

[0194] Error vector (or observation vector): Represents the difference between observed values ​​and model predictions. In remote sensing, this typically refers to the difference between the actual coordinates of a ground control point and the coordinates predicted by a geometric model. Each element is an error term.

[0195] The coefficient matrix (design matrix) contains coefficients related to the model parameters, used to translate changes in parameters into changes in observations. Each row corresponds to an observation, and each column corresponds to a model parameter. In the design matrix, different combinations of rows and columns can represent different geometric relationships and physical models.

[0196] : Parameter vector, which contains the model parameters that need to be estimated, such as the sensor position, attitude parameters, and coordinates of ground points. Each element is a model parameter.

[0197] The observation vector (or constant term vector) contains the actual observed values, usually referring to the coordinates of ground control points. Each element is an observation.

[0198] The weight matrix is ​​a diagonal matrix where the elements on the diagonal represent the precision or weight of each observation. In remote sensing data processing, different observations may have different reliability or precision, and therefore their weights in the error equation will also differ. Used to adjust the effect of different observations on the solution in the least squares method.

[0199] In the formula The coefficient matrix is ​​represented as follows:

[0200] (37);

[0201] in The attitude angle coefficient matrix specifically includes:

[0202] (38);

[0203] in These represent the two equations established using each pair of observation data. According to formula (34), let:

[0204] ;

[0205] but The result of each term is obtained from the following formula (39):

[0206]

[0207] After substituting the initial values ​​into formula (39), the values ​​of each coefficient can be obtained. Since each CCD in each band of the panchromatic and multispectral systems is located on the same camera, the same bias matrix is ​​used to correct the exterior orientation elements for each band.

[0208] In formula (37), These represent the pointing angle coefficient matrices for each band of the panchromatic and multispectral spectra, respectively. For example, The matrix is ​​specifically as follows:

[0209] (40);

[0210] in These represent the two equations established for each pair of observation data. According to the formula...

[0211] ;

[0212] Substituting into formula (34), we get:

[0213] (41);

[0214] but Each item, after calculation, is:

[0215] , (42);

[0216] The solution method is the same. same.

[0217] In formula (37), The matrix representing the connection point coefficients is expressed as:

[0218] (43);

[0219] in Let represent the three-dimensional coordinates of the ground point corresponding to the nth pair of points with the same name. According to formula (43), let: Later it was found that:

[0220] (44);

[0221] Based on the ellipsoid model, It can be done Then the only independent variables in the formula are Coefficient matrix The value was calculated to be:

[0222] (45);

[0223] In the formula:

[0224] ;

[0225] in and These represent the major and minor axes of the ellipsoid, respectively. This represents the altitude of a point on the ground. For the WGS84 ellipsoid, we have... m, m.

[0226] In formula (45), the coefficient matrix Corresponding adjustment parameters for:

[0227] (46);

[0228] in Represents the pointing angle model coefficients for the nth band. This represents the geographic coordinates corresponding to the nth pair of points with the same name.

[0229] In formula (45) This represents the initial value calculated based on each observation. ), This represents the observation weight matrix.

[0230] According to the least squares principle, the calculation is as follows: In the formula, all values ​​are known quantities. The initial values ​​are continuously updated iteratively using the calculated initial values. The result obtained when the initial values ​​are less than a certain threshold is the joint calibration parameter.

[0231] Understandably, this application analyzes the error sources and characteristics of interior and exterior orientation elements in the construction process of a linear array pushbroom imaging model. Based on the satellite's own orbit, attitude, and timing characteristics, it analyzes and establishes corresponding orbit, attitude, and timing models. It addresses the attitude jitter problem caused by satellite platform instability and proposes a method based on a penalized spline function model to reduce model errors caused by insufficient model accuracy, ensuring the positioning accuracy of a rigorous model. Based on the error sources and characteristics of positioning accuracy, a corresponding geometric calibration model is constructed. Control data is used to improve the geometric positioning accuracy of the images. While improving positioning accuracy, the consistency accuracy between multispectral images is considered. High-precision corresponding points between panchromatic and multispectral data are used for joint calibration to correct the geometric positioning model, ensuring the positioning consistency of the geometric positioning models for each band of images, and providing accuracy support for the final application of the satellite.

[0232] The geometric positioning processing device based on linear array pushbroom imaging satellite provided in this application is described below. The geometric positioning processing device based on linear array pushbroom imaging satellite described below can be referred to in correspondence with the geometric positioning processing method based on linear array pushbroom imaging satellite described above.

[0233] This application provides a geometric positioning processing device based on a linear array pushbroom imaging satellite, comprising: a construction module and an acquisition module. Wherein:

[0234] The module is used to build the bias matrix and the pointing angle polynomial model. The bias matrix is ​​used to correct the pointing of the imaging ray of the probe. The pointing angle polynomial model is used to simulate the process of the pointing angle changing with the probe image column. The pointing angle is the angle of the probe's imaging ray relative to the direction perpendicular to the track and along the track.

[0235] The acquisition module is used to adjust the rigorous geometric positioning model based on the bias matrix and the pointing angle polynomial model, and to acquire the joint calibration model. The rigorous geometric positioning model is used to characterize the geometric relationship between image points and ground points.

[0236] It is understood that the detailed functional implementation of each of the above units / modules can be found in the description in the aforementioned method embodiments, and will not be repeated here.

[0237] It should be understood that the above-described device is used to execute the methods in the above embodiments. The implementation principle and technical effect of the corresponding program modules in the device are similar to those described in the above methods. The working process of the device can be referred to the corresponding process in the above methods, and will not be repeated here.

[0238] Based on the methods in the above embodiments, this application provides an electronic device that may include a processor, a communications interface, a memory, and a communication bus, wherein the processor, communications interface, and memory communicate with each other via the communication bus. The processor may invoke logical instructions stored in the memory to execute the methods in the above embodiments.

[0239] Furthermore, the logical instructions in the aforementioned memory can be implemented as software functional units and, when sold or used as independent products, can be stored in a computer-readable storage medium. Based on this understanding, the technical solution of this application, in essence, or the part that contributes to the prior art, or a portion of the technical solution, can be embodied in the form of a software product. This computer software product is stored in a storage medium and includes several instructions to cause a computer device (which may be a personal computer, server, or network device, etc.) to execute all or part of the steps of the methods described in the various embodiments of this application.

[0240] Based on the methods in the above embodiments, this application provides a computer-readable storage medium storing a computer program that, when run on a processor, causes the processor to execute the methods in the above embodiments.

[0241] Based on the methods in the above embodiments, this application provides a computer program product that, when run on a processor, causes the processor to execute the methods in the above embodiments.

[0242] It is understood that the various numerical designations used in the embodiments of this application are merely for the convenience of description and are not intended to limit the scope of the embodiments of this application.

[0243] Those skilled in the art will readily understand that the above description is merely a preferred embodiment of this application and is not intended to limit this application. Any modifications, equivalent substitutions, and improvements made within the spirit and principles of this application should be included within the scope of protection of this application.

Claims

1. A geometric positioning processing method based on linear array pushbroom imaging satellites, characterized in that, include: An offset matrix and a pointing angle polynomial model are constructed. The offset matrix is ​​used to correct the pointing of the imaging ray of the probe. The pointing angle polynomial model is used to simulate the process of the pointing angle changing with the probe image column. The pointing angle is the angle of the probe's imaging ray relative to the direction perpendicular to the track and along the track. Based on the bias matrix and pointing angle polynomial model, the rigorous geometric positioning model is adjusted to obtain the joint calibration model. The rigorous geometric positioning model is used to characterize the geometric relationship between image points and ground points. Also includes: Based on spline functions and the following loss function, construct a penalized spline function model: ； The on-board attitude data was simulated using a penalized spline function model; Where W represents the parameter matrix of the spline function model, B represents the basis matrix of the spline function, c represents the undetermined coefficients of the spline function, and a represents the attitude recording value. Let P represent the penalty factor, and let P represent the positive definite penalty matrix.

2. The geometric positioning processing method based on linear array pushbroom imaging satellite according to claim 1, characterized in that, The method of adjusting the rigorous geometric positioning model based on the bias matrix and pointing angle polynomial model to obtain the joint calibration model includes obtaining the joint calibration model through the following formula: ； in, This indicates the coordinates of the pixels in the satellite image in the WGS84 coordinate system. express The position vector of the phase center of the GPS module in the WGS84 coordinate system. This represents the proportionality coefficient. for The transformation matrix for converting the J2000 coordinate system to the WGS84 coordinate system at time t. for The transformation matrix from the body coordinate system to the geocentric inertial coordinate system at time t. Represents the bias matrix. This represents the transformation matrix from the camera coordinate system to the body coordinate system; Represents a pointing angle polynomial model. Representing pixels coordinate, and The subscripts represent the coefficients of the pointing-angle polynomial model. This indicates the degree of the polynomial.

3. The geometric positioning processing method based on linear array pushbroom imaging satellites according to claim 2, characterized in that, The bias matrix is ​​specifically: ； in, This represents the rotation angle of the imaging ray of the probe around the x-axis of the body coordinate system. This represents the rotation angle of the imaging ray of the probe around the y-axis of the body coordinate system. This represents the rotation angle of the imaging ray of the probe around the z-axis of the body coordinate system.

4. The geometric positioning processing method based on linear array pushbroom imaging satellite according to claim 3, characterized in that, Also includes: Substitute the control points into the joint calibration model to construct the control point equations: ； By solving the control point equations, we can obtain... The amount of correction; in, The three-dimensional coordinates of the control points are represented by either matching the panchromatic image with the control image or matching the multispectral image with the control image.

5. The geometric positioning processing method based on linear array pushbroom imaging satellite according to claim 4, characterized in that, Also includes: Substitute the corresponding points into the joint calibration model to construct the equations for the corresponding points: ； By solving the equations of the corresponding points, we can obtain... The amount of correction; in, This represents the three-dimensional coordinates of a corresponding point, which is obtained by matching a panchromatic image with a band in a multispectral image.

6. A geometric positioning processing device based on a linear array pushbroom imaging satellite, characterized in that, include: The module is used to build the bias matrix and the pointing angle polynomial model. The bias matrix is ​​used to correct the pointing of the imaging ray of the probe. The pointing angle polynomial model is used to simulate the process of the pointing angle changing with the probe image column. The pointing angle is the angle of the probe's imaging ray relative to the direction perpendicular to the track and along the track. The acquisition module is used to adjust the rigorous geometric positioning model based on the bias matrix and the pointing angle polynomial model, and to acquire the joint calibration model. The rigorous geometric positioning model is used to characterize the geometric relationship between image points and ground points. Also includes: Based on spline functions and the following loss function, construct a penalized spline function model: ； The on-board attitude data was simulated using a penalized spline function model; Where W represents the parameter matrix of the spline function model, B represents the basis matrix of the spline function, c represents the undetermined coefficients of the spline function, and a represents the attitude recording value. Let P represent the penalty factor, and let P represent the positive definite penalty matrix.

7. An electronic device, characterized in that, include: At least one memory for storing computer programs; At least one processor is configured to execute a program stored in the memory, wherein when the program stored in the memory is executed, the processor is configured to perform the method as described in any one of claims 1-5.

8. A computer-readable storage medium storing a computer program, characterized in that, When the computer program is run on the processor, it causes the processor to perform the method as described in any one of claims 1-5.

9. A computer program product, characterized in that, When the computer program product is run on a processor, the processor causes the processor to perform the method as described in any one of claims 1-5.

Images