CN111708095B - Satellite gravitational field inversion method and system based on bidirectional integration - Google Patents
Satellite gravitational field inversion method and system based on bidirectional integration Download PDFInfo
- Publication number
- CN111708095B CN111708095B CN202010446920.1A CN202010446920A CN111708095B CN 111708095 B CN111708095 B CN 111708095B CN 202010446920 A CN202010446920 A CN 202010446920A CN 111708095 B CN111708095 B CN 111708095B
- Authority
- CN
- China
- Prior art keywords
- satellite
- gravitational field
- orbit
- reverse
- field inversion
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Active
Links
Classifications
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
- G01V7/00—Measuring gravitational fields or waves; Gravimetric prospecting or detecting
- G01V7/02—Details
- G01V7/06—Analysis or interpretation of gravimetric records
Landscapes
- Physics & Mathematics (AREA)
- Life Sciences & Earth Sciences (AREA)
- General Life Sciences & Earth Sciences (AREA)
- General Physics & Mathematics (AREA)
- Geophysics (AREA)
- Radio Relay Systems (AREA)
- Navigation (AREA)
Abstract
The invention relates to a satellite gravitational field inversion method and system based on 'bidirectional integration'. The satellite gravitational field inversion method based on the 'bidirectional integration' comprises the following steps: acquiring an initial epoch and an end epoch of an orbit in a gravitational field of a satellite to be inverted; acquiring a satellite gravitational field inversion model; and inputting the initial epoch and the final epoch into the satellite gravitational field inversion model to finish inversion of the satellite gravitational field. According to the satellite gravitational field inversion method and system based on the bidirectional integration, the satellite gravitational field is inverted by adopting the satellite gravitational field inversion model, so that the influence of system errors on the satellite gravitational field inversion result can be greatly reduced, and the satellite gravitational field inversion precision is improved.
Description
Technical Field
The invention relates to the technical field of satellite gravitational field inversion, in particular to a satellite gravitational field inversion method and system based on 'bidirectional integration'.
Background
The earth's gravitational field is a fundamental physical field that reflects the distribution, movement, and variation of earth's mass. The precise structure and the space-time variation of the earth gravity field are precisely determined, and the precise structure and the space-time variation have extremely important roles and significance in the research of military national defense construction, national economy, space science, earth science and related disciplines. The research of the earth's gravitational field is a fundamental task and is also a research core and hotspot in the field of geodetics. Satellite gravity measurement opens up an epoch of the gravitational field of human detection earth and is considered to be a high-efficiency gravity detection technology with the highest value and application prospect in the beginning of the 21 st century.
The gravity satellite data is utilized to invert the earth gravity field, namely, a satellite gravity field inversion model is constructed, which is a core task of the earth gravity satellite. Common methods for inverting the earth's gravitational field using gravity satellite data are divided into time domain methods and space domain methods. The time domain method mainly comprises the following steps: kaula linear perturbation method, dynamic method, baseline method, energy method, acceleration method, short arc value method, etc., which is suitable for the resolution of CHAMP satellites (Challenging Minisatellite Payload, challenging small satellite payload satellites) and GRACE (Gravity Recovery and Climate Experiment ) satellites. The airspace method is suitable for the data calculation of GOCE satellites (GravityField and Stead-state Ocean Circulation Explore, gravitational field and static ocean current detection satellites) to calculate the gravitational field of the earth. Among the numerous satellite gravitational field resolution methods, the dynamic method is a classical method that is widely adopted.
When solving the earth satellite gravitational field inversion model by adopting a dynamic method, solving Newton's motion equation and variation equation is involved, orbit integration is needed first, then a normal equation is constructed based on residual information, and gravitational field solution is further carried out. The numerical integration is to directly perform numerical integration on a satellite motion equation and a variation equation, and gradually calculate satellite positions, speeds and state transition matrixes at any moment by taking satellite positions and speeds of reference epochs as initial values.
The numerical integration method gives an approximation at discrete steps, and there are unavoidable errors, in addition to the initial value errors, truncation errors and rounding errors. Errors (including initial errors and rounding errors) generated by a certain step are propagated down, and errors are accumulated, so that the corresponding numerical method is stable only when the accumulation of errors is controlled.
Orbit integration is in fact a process of extrapolating future satellite orbits by integrating satellite equations of motion based on a set of initial values. In the extrapolation process, the longer the calculation interval is, the larger the error accumulated by integration is due to the influence of error propagation. The unidirectional error accumulation effect of the orbit integration makes the normal equation information obtained based on the error propagation law worse, resulting in larger estimated bit coefficient deviation.
Therefore, aiming at the situation of unidirectional transmission containing the systematic error, the satellite gravitational field inversion method capable of greatly weakening the influence of the systematic error on the satellite gravitational field inversion result and improving the gravitational field inversion precision is a technical problem to be solved in the field.
Disclosure of Invention
The invention aims to provide a satellite gravitational field inversion method and system based on 'bidirectional integration', which are used for greatly weakening the influence of system errors on satellite gravitational field inversion results and improving satellite gravitational field inversion accuracy.
In order to achieve the above object, the present invention provides the following solutions:
a satellite gravitational field inversion method based on 'two-way integration', comprising:
acquiring an initial epoch and an end epoch of an orbit in a gravitational field of a satellite to be inverted;
acquiring a satellite gravitational field inversion model;
and inputting the initial epoch and the final epoch into the satellite gravitational field inversion model to finish inversion of the satellite gravitational field.
Optionally, the construction process of the satellite gravitational field inversion model specifically includes:
acquiring gravity satellite actual measurement data; the gravity satellite actual measurement data comprises: orbit data, accelerometer data, star sensor data, and inter-satellite ranging data; the track data includes: an initial epoch and an end epoch of the track;
correcting the initial forward track value of the track by using the initial epoch as an integral initial value and adopting a forward numerical integration method to obtain a corrected forward track initial value;
correcting the initial value of the reverse track of the track by using the last epoch as an integral initial value and adopting a reverse numerical integration method to obtain a corrected initial value of the reverse track;
the corrected forward orbit initial value is taken as an integral initial value, and a forward numerical integration method is adopted to obtain a forward reference orbit;
determining a forward reference inter-satellite variability according to the forward reference orbit;
taking the corrected initial value of the reverse track as an integral initial value, and adopting a reverse numerical integration method to obtain a reverse reference track;
determining a reverse reference inter-satellite distance variability according to the reverse reference orbit;
and constructing a satellite gravitational field inversion model according to the forward reference orbit, the forward reference inter-satellite space variability, the reverse reference orbit and the reverse reference inter-satellite space variability.
Optionally, the constructing a satellite gravitational field inversion model according to the forward reference orbit, the forward reference inter-satellite variability, the backward reference orbit and the backward reference inter-satellite variability specifically includes:
constructing a forward observation equation according to the forward reference orbit and the forward reference inter-satellite distance variability;
determining a forward satellite gravitational field inversion model according to the forward observation equation;
constructing a reverse observation equation according to the reverse reference orbit and the reverse reference inter-satellite distance variability;
determining an inverse model of the inverse satellite gravitational field according to the inverse observation equation;
and determining a satellite gravitational field inversion model according to the forward satellite gravitational field inversion model and the reverse satellite gravitational field inversion model.
Optionally, before obtaining the initial epoch and the final epoch of the orbit in the gravitational field of the satellite to be inverted, the method further includes:
collecting gravity satellite actual measurement data in the satellite gravity field to be inverted; the gravity satellite actual measurement data comprises: orbit data, accelerometer data, star sensor data, and inter-satellite ranging data; the track data includes: an initial epoch and an end epoch of the track;
preprocessing the collected gravity satellite actual measurement data.
A satellite gravitational field inversion system based on "two-way integration", comprising:
the data acquisition module is used for acquiring an initial epoch and an end epoch of an orbit in a gravitational field of a satellite to be inverted;
the model acquisition module is used for acquiring a satellite gravitational field inversion model;
and the gravity field inversion module is used for inputting the initial epoch and the final epoch into the satellite gravity field inversion model to finish inversion of the satellite gravity field.
Optionally, the system further comprises:
the measured data acquisition module is used for acquiring the measured data of the gravity satellite; the gravity satellite actual measurement data comprises: orbit data, accelerometer data, star sensor data, and inter-satellite ranging data; the track data includes: an initial epoch and an end epoch of the track;
the forward track initial value correction module is used for correcting the forward track initial value of the track by using the initial epoch as an integral initial value and adopting a forward numerical integration method to obtain a corrected forward track initial value;
the reverse orbit initial value correction module is used for correcting the reverse orbit initial value of the orbit by adopting a reverse numerical integration method by taking the last epoch as an integral initial value to obtain a corrected reverse orbit initial value;
the forward reference orbit determination module is used for obtaining a forward reference orbit by taking the corrected forward orbit initial value as an integral initial value and adopting a forward numerical integration method;
the forward reference inter-satellite variability determining module is used for determining a forward reference inter-satellite variability according to the forward reference orbit;
the reverse reference orbit determination module is used for obtaining a reverse reference orbit by taking the corrected initial value of the reverse orbit as an integral initial value and adopting a reverse numerical integration method;
the reverse reference inter-satellite distance variable rate determining module is used for determining a reverse reference inter-satellite distance variable rate according to the reverse reference orbit;
the satellite gravitational field inversion model construction module is used for constructing a satellite gravitational field inversion model according to the forward reference orbit, the forward reference inter-satellite variability, the reverse reference orbit and the reverse reference inter-satellite variability.
Optionally, the satellite gravitational field inversion model building module specifically includes:
a forward observation equation construction unit for constructing a forward observation equation according to the forward reference orbit and the forward reference inter-satellite variability;
the forward satellite gravitational field inversion model determining unit is used for determining a forward satellite gravitational field inversion model according to the forward observation equation;
the reverse observation equation construction unit is used for constructing a reverse observation equation according to the reverse reference orbit and the reverse reference inter-satellite space variability;
the inverse satellite gravitational field inversion model determining unit is used for determining an inverse satellite gravitational field inversion model according to the inverse observation equation;
the satellite gravitational field inversion model determining unit is used for determining a satellite gravitational field inversion model according to the forward satellite gravitational field inversion model and the reverse satellite gravitational field inversion model.
Optionally, the system further comprises:
the data acquisition module is used for acquiring the gravity satellite actual measurement data in the satellite gravity field to be inverted; the gravity satellite actual measurement data comprises: orbit data, accelerometer data, star sensor data, and inter-satellite ranging data; the track data includes: an initial epoch and an end epoch of the track;
the preprocessing module is used for preprocessing the acquired gravity satellite actual measurement data.
According to the specific embodiment provided by the invention, the invention discloses the following technical effects:
according to the satellite gravitational field inversion method and system based on the bidirectional integration, after the initial epoch and the final epoch of the orbit in the satellite gravitational field to be inverted are acquired, inversion of the satellite gravitational field can be completed by adopting the satellite gravitational field inversion model, so that influence of system errors on satellite gravitational field inversion results can be greatly reduced, and satellite gravitational field inversion accuracy is improved.
Drawings
In order to more clearly illustrate the embodiments of the present invention or the technical solutions of the prior art, the drawings that are needed in the embodiments will be briefly described below, it being obvious that the drawings in the following description are only some embodiments of the present invention, and that other drawings may be obtained according to these drawings without inventive effort for a person skilled in the art.
FIG. 1 is a flow chart of a satellite gravitational field inversion method based on "two-way integration" provided by an embodiment of the invention;
FIG. 2 is a schematic diagram of the principle of "bi-directional integration" in an embodiment of the present invention;
FIG. 3 is a flow chart of a gravity field inversion method of a low-rail gravity satellite GRACE in an embodiment of the present invention;
FIGS. 4 a-4 c are graphs of track residuals before correcting initial data values for GRACE tracks for 24H long arc segments in an embodiment of the present invention;
FIGS. 5 a-5 c are graphs of corrected track residuals after "bi-directional integration" of GRACE track data for 24H long arc segments in an embodiment of the present invention;
FIGS. 6 a-6 c are graphs of track residuals before correcting initial data values for GRACE tracks for 6H short arc segments in an embodiment of the present invention;
FIGS. 7 a-7 c are graphs of corrected track residuals after "bi-directional integration" of GRACE track data for a 6H short arc segment in an embodiment of the present invention;
fig. 8 is a schematic structural diagram of a satellite gravitational field inversion system based on "two-way integration" according to an embodiment of the present invention.
Detailed Description
The following description of the embodiments of the present invention will be made clearly and completely with reference to the accompanying drawings, in which it is apparent that the embodiments described are only some embodiments of the present invention, but not all embodiments. All other embodiments, which can be made by those skilled in the art based on the embodiments of the invention without making any inventive effort, are intended to be within the scope of the invention.
The invention aims to provide a satellite gravitational field inversion method and system based on 'bidirectional integration', which are used for greatly weakening the influence of system errors on satellite gravitational field inversion results and improving satellite gravitational field inversion accuracy.
The basic idea of the satellite gravitational field inversion method and system based on 'bidirectional integration' provided by the invention is as follows: when the orbit numerical integration is carried out, the opposite-direction integration is added, and the inversion solution is calculated to obtain higher precision by utilizing the two integrations with opposite directions to more accurately integrate the orbit.
In order that the above-recited objects, features and advantages of the present invention will become more readily apparent, a more particular description of the invention will be rendered by reference to the appended drawings and appended detailed description.
Fig. 1 is a flowchart of a satellite gravitational field inversion method based on "bidirectional integration", which is provided in an embodiment of the present invention, as shown in fig. 1, and the satellite gravitational field inversion method based on "bidirectional integration" includes:
s1, acquiring an initial epoch and an end epoch of an orbit in a gravitational field of a satellite to be inverted.
S2, acquiring a satellite gravitational field inversion model.
S3, inputting the initial epoch and the final epoch into the satellite gravitational field inversion model to finish inversion of the satellite gravitational field.
The construction process of the satellite gravitational field inversion model specifically comprises the following steps:
and acquiring the gravity satellite actual measurement data. The measured data includes: the initial epoch and the last epoch of the track.
And correcting the forward track initial value of the track by using the initial epoch as an integral initial value and adopting a forward numerical integration method to obtain a corrected forward track initial value.
And correcting the initial value of the reverse track of the track by using the last epoch as an integral initial value and adopting a reverse numerical integration method to obtain a corrected initial value of the reverse track.
And taking the corrected forward orbit initial value as an integral initial value, and obtaining a forward reference orbit by adopting a forward numerical integration method.
And determining a forward reference inter-satellite variability according to the forward reference orbit.
And taking the corrected initial value of the reverse track as an integral initial value, and obtaining a reverse reference track by adopting a reverse numerical integration method.
And determining a reverse reference inter-satellite distance variability according to the reverse reference orbit.
And constructing a satellite gravitational field inversion model according to the forward reference orbit, the forward reference inter-satellite space variability, the reverse reference orbit and the reverse reference inter-satellite space variability.
The construction process of the satellite gravitational field inversion model mainly adopts a construction idea based on 'bidirectional integration'. The principle of "two-way integration" is schematically shown in fig. 2, wherein forward integration takes the initial epoch of the track as an initial value, spans one step forward in time, integrates once, and corrects it according to the initial condition. The backward integration takes the last epoch of the track as an initial value, and the backward time is one step length, and the backward integration is corrected according to the initial condition. In the example of 'two-way integration' in fig. 2, in order to demonstrate the effect, the invention intercepts a 100-minute track arc segment as a real track, respectively takes two ends of the arc segment as initial epochs of forward and reverse integration, the integration step length is 5 minutes, the integration duration is 50 minutes, and forward and reverse integration is carried out. The integration step length and the effectiveness of the integration interval are main factors influencing the integration precision due to the accumulation characteristic of the integration errors, and practice shows that the integration step length of 5 seconds has better integration precision. Therefore, in the subsequent gravity field inversion process, the invention preferentially selects 5 seconds as an integral step length.
In order to further improve accuracy of inversion of the gravitational field, the process of constructing a satellite gravitational field inversion model according to the forward reference orbit, the forward reference inter-satellite variability, the reverse reference orbit and the reverse reference inter-satellite variability specifically includes:
and constructing a forward observation equation according to the forward reference orbit and the forward reference inter-satellite space variability.
And determining a forward satellite gravitational field inversion model according to the forward observation equation.
And constructing a reverse observation equation according to the reverse reference orbit and the reverse reference inter-satellite distance variability.
And determining an inversion model of the reverse satellite gravitational field according to the reverse observation equation.
And carrying out weighted fusion processing on the forward satellite gravitational field inversion model and the reverse satellite gravitational field inversion model to obtain a satellite gravitational field inversion model.
The specific process of constructing the satellite gravitational field inversion model is described below by taking a mode of bidirectional numerical integration by using a KSG integrator as an example.
When a KSG integrator is used for bidirectional numerical integration, the forward and reverse integration formulas have sign differences.
The whole process of constructing the satellite gravitational field inversion model is as follows:
in satellite gravitational field inversion solution, parameters to be estimated include: double star state vectors (position, velocity), double star accelerometer bias and scale factors, and gravitational field model bit coefficients. The gravity field inversion solution order is N, and the satellite position vector is assumed to beThe velocity vector is +.> Is a vector to be estimated, and
wherein, the liquid crystal display device comprises a liquid crystal display device,
x, y, z represent three axes of the coordinate system, respectively, as three components of the position. v x ,v y ,v z Is a three-component speed. b x ,b y ,b z Is a three component deviation. k (k) x ,k y ,k z Is a scale three component. C (C) nm 、S nm All are gravity field model bit coefficients, n is model order, m is model number, then m=0, coefficient C n0 Called harmonic term coefficients, coefficient C when n=m nm 、S nm Called fan harmonic term coefficients, coefficient C when n.noteq.m nm 、S nm Called the field harmonic term coefficient.
The corresponding state vectorVector initial value to be estimated +.>Is a partial derivative of phi>Vector initial value to be estimated +.>Partial derivative of>The method comprises the following steps:
assuming an integration step size of h, the integrator order is i (a 14 th order KSG integrator is actually selected, i.e., i=14). Because the KSG integrator adopts a central iteration initialization process (the process is not the key point of the patent and is not specifically described), i node values are needed to be calculated in integration startingWherein->The integral vector isAnd->Respectively carrying out forward integration and reverse integration, wherein the calculation formula of each node in the integration starting stage is as follows:
1) During forward integration, the initial epoch of the known observation arc segment is selected as the integration initial epoch, and the initial time is recorded as t 0 The initial value of the integral vector is recorded asAnd->Then:
wherein, the liquid crystal display device comprises a liquid crystal display device,for the start time integration right function +.>At t k Integrating vector of time,/->At t k The integral right function of time (i.e. t k The satellite experiences a perturbation acceleration).
2) When in backward integration, the last epoch of the observation arc segment is selected as the integration starting time, and the initial time is recorded as t 0 The initial value of the integral vector is recorded asAnd->Then:
from the above deductions, the coefficients in the KSG formula of the forward integration and the backward integration are identical, and when the method is implemented, h in the forward integration formula can be replaced by-h, namely the backward integration formula.
The partial derivative of the initial value of the parameter to be estimated can be obtained through the numerical integration, the orbit residual error and the inter-satellite variable rate residual error can be calculated and obtained by utilizing satellite orbit data and inter-satellite variable rate residual error, and after the orbit residual error and the inter-satellite variable rate residual error are obtained, the observation equation of the parameter to be estimated (comprising the initial state parameter, accelerometer deviation, scale parameter and gravitational field bit coefficient) of a single epoch can be established by the following formula:
Δr i an orbit residual of satellite i (i=a or B),For inter-satellite variability residual, +.>R is the inter-satellite variability reference value i0 ,B i0 ,K i0 Representing initial state parameters, accelerometer bias and scale parameters of satellite i (i=a or B), respectively, beta being the gravitational field coefficient, +.>Solving, < > in the course of integration of the equations of orbital motion>For the partial derivative of the inter-satellite variability to the initial state parameters and the gravitational field coefficients, the partial derivative is composed of +.>Linear combination.
The two observation equations can be arranged to obtain respective error equations, and the single epoch error equations are uniformly expressed as follows:
wherein V is i To observe errors, L i Is track residual error or inter-satellite variability residual error, A i The coefficient matrix of the error equation is formed by combining the partial derivatives, deltaX is the variation of the parameter to be estimated, and the value of i is from 1 to the number of single arc segment epochs. In practical calculation, a plurality of epochs are generally used as a group to form a matrix for calculation, so that the error equation matrix of the epochs is as follows:
V＝A·ΔX-L (12)
wherein V is an error vector, L is a residual vector, and A is a coefficient matrix. The equation that can be derived from the error equation is:
wherein A is T The transpose of the coefficient matrix A, N is the normal left matrix, B is the normal right matrix, and W is the parameter weight matrix to be estimated (usually taken as the identity matrix).
In actual calculation, according to the characteristics of parameters to be estimated, the normal equation is expressed in blocks:
x L 、x G representing local parameters (initial state parameters of A star and B star, accelerometer bias parameters) and global parameters (accelerometer scale parameters of A star and B star, gravitational field coefficients), respectively, N 11 ，N 22 Is a sub-matrix corresponding to local and global parameters,is normal right vector.
After the orbit observation values GNV1B of the A star and the B star are arranged by an error equation (formula 12), an orbit observation equation (formula 9) can form a normal equation
Similarly, each arc segment inter-satellite ranging observation KBRR (K-band range-rate, K-band inter-satellite ranging range rate) may form a normal equation
In combination with the track data and the KBRR data, the following relationship exists:
wherein, the liquid crystal display device comprises a liquid crystal display device,
σ kbrr for KBRR observation accuracy, sigma orb The track observation precision is obtained. In the calculation process, the KBRR observation precision is 0.2nm/s, and the orbit observation precision is 2cm.
The normal equation is solved by adopting a conventional solution, and the normal equation about the global parameter can be obtained by eliminating the local parameter of the normal equation
The above is described as
Solving the reduced post-equation to obtain the global parameter variation as
Will Δx G The equation below is replaced by the local parameter variation delta x which can be solved L
The above is a forward solution for a single arc segment. In practice, to obtain a high-precision resolving result according to satellite observation characteristics, a plurality of arc segment data are generally accumulated and resolved, a time-varying gravitational field is generally resolved in units of months, and a static gravitational field needs to use observation data of longer arc segments.
In order to ensure the accuracy and the integrity of the data, the method provided by the invention further comprises the following steps before the initial epoch and the final epoch of the orbit in the gravitational field of the satellite to be inverted are acquired:
collecting gravity satellite actual measurement data in the satellite gravity field to be inverted; the gravity satellite actual measurement data comprises: orbit data, accelerometer data, star sensor data, and inter-satellite ranging data.
Preprocessing the collected gravity satellite actual measurement data, specifically preprocessing actual measurement data (the actual measurement data comprises the initial epoch, the final epoch, accelerometer observation data, star sensor observation data, inter-satellite ranging observation data and the like) by coarse difference rejection, downsampling and the like, and completing data preparation.
In another embodiment of the present invention, a low-rail gravity satellite GRACE is taken as an example, and the satellite gravity field inversion method based on "bidirectional integration" provided by the present invention is specifically described.
The inversion process of the gravity field of the low-orbit gravity satellite GRACE by adopting the satellite gravity field inversion method based on the 'bidirectional integration' provided by the invention specifically comprises the following steps:
step 1: four main load 1B-level data (mainly orbit data, inter-satellite ranging data, accelerometer data and star sensor data) were collected for the low-rail gravity satellite GRACE (Gravity Recovery and Climate Experiment ). Preprocessing the collected four main load 1B-level data (including coarse difference rejection, error correction, space-time reference unification, interpolation deficiency, downsampling and the like), and forming two groups of data sets of 24H and 6H (namely, the arc lengths are 24 hours and 6 hours respectively) arc segments. Each set of data contains level 1B data for four main loads, orbit data, inter-satellite ranging data, accelerometer data, and satellite sensor data, GNV1B, KBR1B, ACC B and SCA1B, respectively. The two sets of data represent different arc lengths for comparison of the bi-directional integration effect at different arc lengths. Integration of different arc lengths, track error accumulation is different. Theoretically, the longer the integral arc length, the greater the error accumulation and the rapid accumulation.
Step 2: and (3) taking the initial epoch of the GRACE satellite orbit arc segment GNV1B data obtained in the step (1) as an initial value, adopting a KSG integrator (taking 5s as a step length, and carrying out forward integration on arc segments by arc segment to obtain a forward reference orbit and a partial derivative, and carrying out initial value correction on the orbit arc segment initial epoch by arc segment, wherein the step length is not specifically described in the follow-up, and the default integral step length is 5 s).
The orbit initial value correction only uses an orbit observation equation, and the parameters to be estimated are only satellite state vectors, and the specific correction process is as follows: 1) Track residuals are obtained using the forward reference track and the track. 2) Only the partial derivatives with respect to the state vector are extracted from the partial derivatives, forming a coefficient matrix. 3) And constructing an orbit normal equation by using the orbit residual error and the coefficient matrix. 4) And solving an algorithm equation to obtain the initial value correction quantity of the orbit. 5) And performing initial track value correction by using the initial track value correction amount.
Step 3: and (3) taking the last epoch of the GRACE satellite orbit arc segment GNV1B data obtained in the step (1) as an initial value, adopting a KSG integrator to carry out reverse integration on arc segments by arc segments, carrying out initial value correction on the last epoch of the orbit arc segment by arc segments, wherein the correction process is consistent with the forward orbit initial value correction, and only observing data are inconsistent (the observing data are GNV1B data, but the initial epoch of GNV1B is taken as the initial value for integration during the forward initial value correction, and the last epoch of GNV1B is taken as the initial value for integration during the reverse initial value correction).
Step 4: and taking the initial epoch of the corrected arc segment as an initial value, performing forward integration to obtain a forward reference orbit and a partial derivative, and constructing a forward reference inter-satellite distance variable rate based on the forward reference orbit and a corresponding orbit, wherein the construction process is to utilize the difference of the double-satellite forward reference orbit to form the reference inter-satellite distance variable rate.
Step 5: and taking the end epoch of the corrected arc segment as an initial value, carrying out backward integration to obtain a backward reference orbit and a partial derivative, and constructing a backward reference inter-satellite variability based on the backward reference orbit and a corresponding orbit, wherein the construction process is consistent with the forward reference inter-satellite variability construction process.
Step 6: in the fusion calculation stage, two fusion calculation strategies can be adopted, and specifically include:
a) Bit coefficient fusion solution
And constructing a normal equation by using the forward reference orbit and the forward inter-satellite variability, and acquiring a forward satellite gravitational field inversion model according to the normal equation. And constructing a normal equation by using the reverse reference orbit and the reverse reference inter-satellite space variability and solving to obtain an inversion model of the gravity field of the reverse satellite. And obtaining and fusing the two groups of model weighted fusion to obtain the satellite gravitational field inversion model. The whole bit coefficient fusion calculation process is shown in fig. 3.
The weighted fusion strategy adopts a spectrum combination method based on the gravity field model bit coefficient. Two ways are adopted to estimate gravity field bit coefficients, and a linear unbiased estimation model is expressed as follows:
wherein, the liquid crystal display device comprises a liquid crystal display device,for the fused model bit coefficients, beta 1 For positively estimated bit coefficients, beta 2 For bit coefficients estimated in reverse, W 1 For forward estimating the weight of the bit coefficient, W 2 Is the weight of the bit coefficient estimated in the reverse direction. According to the least square principle, W is taken because the same set of observation data is used in two modes, which are neglected 1 ＝W 2 ＝0.5。
b) French equation fusion solution
Construction using forward integration orbit and forward inter-star variabilityAnd (3) constructing a normal equation by using the reverse integral orbit and the reverse inter-satellite variability, and carrying out weighted fusion and solution on the two sets of normal equations to obtain a fused satellite gravitational field inversion model, as shown by a dot-dash line in fig. 3. The model here refers to a set of data results that are actually a data file containing gravitational field spherical harmonic coefficients C nm 、S nm In the actual application process, the model data file is used as input, and different business processes can be performed.
In order to verify the specific effect of the satellite gravitational field inversion method based on the bidirectional integration, the invention adopts GRACE satellite data of 12 months in 2005 to carry out the bidirectional integration inversion resolving experiment. Wherein the numerical integrator employs a KSG integrator of 14 th order (KSG is a fixed order, fixed step, linear, multi-step integrator suitable for second order differential equation, the first name of which is taken from the last name of three scholars researching it: krogh, shampine and golden), the integration step length is 5s, the integration arc length comprises two groups of 24H and 6H, and the mechanical model employed in the orbit integration process is shown in the following Table 1:
TABLE 1
Test results
Forward integration and backward integration are carried out on GRACE track data of 24H long arc segments, track residuals before initial value correction are shown in fig. 4 a-4 c, and track residuals after iterative correction of the track initial values are shown in fig. 5 a-5 c. Forward integration and backward integration are carried out on GRACE track data of the 6H short arc section, track residuals before initial value correction are shown in fig. 6 a-6 c, and track residuals after the track initial value is subjected to iterative correction are shown in fig. 7 a-7 c. The above results indicate that: initial value correction can effectively reduce initial track errors, forward integration and backward integration can obtain reference tracks with higher precision, and unidirectional accumulation of systematic errors in the integration process is effectively reduced, so that transmission of errors in the gravity field inversion process is reduced, and gravity field inversion precision is improved.
Compared with the prior art, the invention has the following advantages:
1. the bidirectional integration is provided, so that the accumulation of unidirectional errors of track integration can be effectively avoided, and the initial track calculation precision is improved.
2. The bidirectional integration method is applied to the gravitational field calculation, so that the accumulation and the diffusion of integration errors can be effectively inhibited, the system errors are eliminated, and the gravitational field inversion accuracy is improved.
In addition, with respect to the satellite gravitational field inversion method based on the "two-way integration" provided above, the present invention correspondingly provides a satellite gravitational field inversion system based on the "two-way integration", as shown in fig. 8, the system includes: a data acquisition module 1, a model acquisition module 2 and a gravity field inversion module 3.
The data acquisition module 1 is used for acquiring an initial epoch and an end epoch of an orbit in a gravitational field of a satellite to be inverted. The model acquisition module 2 is used for acquiring a satellite gravitational field inversion model. The gravity field inversion module 3 is used for inputting the initial epoch and the final epoch into the satellite gravity field inversion model to complete inversion of the satellite gravity field.
To improve inversion and accuracy, the system further comprises: the system comprises an actual measurement data acquisition module, a forward orbit initial value correction module, a reverse orbit initial value correction module, a forward reference orbit determination module, a forward reference inter-satellite space variability determination module, a reverse reference orbit determination module, a reverse reference inter-satellite space variability determination module and a satellite gravitational field inversion model construction module.
The measured data acquisition module is used for acquiring the measured data of the gravity satellite. The measured data includes: the initial epoch and the last epoch of the track. And the forward track initial value correction module is used for correcting the forward track initial value of the track by adopting a forward numerical integration method by taking the initial epoch as an integral initial value to obtain a corrected forward track initial value. And the reverse track initial value correction module is used for correcting the reverse track initial value of the track by adopting a reverse numerical integration method by taking the last epoch as an integral initial value to obtain a corrected reverse track initial value. The forward reference orbit determination module is used for obtaining a forward reference orbit by taking the corrected forward orbit initial value as an integral initial value and adopting a forward numerical integration method. The forward reference inter-satellite variability determination module is used for determining a forward reference inter-satellite variability according to the forward reference orbit. The back reference orbit determination module is used for obtaining a back reference orbit by taking the corrected back orbit initial value as an integral initial value and adopting a back numerical integration method. The reverse reference inter-satellite variability determination module is used for determining a reverse reference inter-satellite variability according to the reverse reference orbit. The satellite gravitational field inversion model construction module is used for constructing a satellite gravitational field inversion model according to the forward reference orbit, the forward reference inter-satellite variability, the reverse reference orbit and the reverse reference inter-satellite variability.
The satellite gravitational field inversion model construction module specifically comprises: the system comprises a forward observation equation construction unit, a forward satellite gravitational field inversion model determination unit, a reverse observation equation construction unit, a reverse satellite gravitational field inversion model determination unit and a satellite gravitational field inversion model determination unit.
The forward observation equation construction unit is used for constructing a forward observation equation according to the forward reference orbit and the forward reference inter-satellite variability. And the forward satellite gravitational field inversion model determining unit is used for determining a forward satellite gravitational field inversion model according to the forward observation equation. The reverse observation equation construction unit is used for constructing a reverse observation equation according to the reverse reference orbit and the reverse reference inter-satellite distance variability. The inverse satellite gravitational field inversion model determining unit is used for determining an inverse satellite gravitational field inversion model according to the inverse observation equation. The satellite gravitational field inversion model determining unit is used for determining a satellite gravitational field inversion model according to the forward satellite gravitational field inversion model and the reverse satellite gravitational field inversion model.
For the accuracy and integrity of the acquired data, the system further comprises: and the data acquisition module and the preprocessing module.
The data acquisition module is used for acquiring the gravity satellite actual measurement data in the satellite gravity field to be inverted; the gravity satellite actual measurement data comprises: orbit data, accelerometer data, star sensor data, and inter-satellite ranging data. The preprocessing module is used for preprocessing the acquired gravity satellite actual measurement data.
In the present specification, each embodiment is described in a progressive manner, and each embodiment is mainly described in a different point from other embodiments, and identical and similar parts between the embodiments are all enough to refer to each other. For the system disclosed in the embodiment, since it corresponds to the method disclosed in the embodiment, the description is relatively simple, and the relevant points refer to the description of the method section.
The principles and embodiments of the present invention have been described herein with reference to specific examples, the description of which is intended only to assist in understanding the methods of the present invention and the core ideas thereof; also, it is within the scope of the present invention to be modified by those of ordinary skill in the art in light of the present teachings. In view of the foregoing, this description should not be construed as limiting the invention.
Claims (6)
1. The satellite gravitational field inversion method based on 'two-way integration' is characterized by comprising the following steps:
acquiring an initial epoch and an end epoch of an orbit in a gravitational field of a satellite to be inverted;
acquiring a satellite gravitational field inversion model;
inputting the initial epoch and the final epoch into the satellite gravitational field inversion model to finish inversion of the satellite gravitational field;
the construction process of the satellite gravitational field inversion model specifically comprises the following steps:
acquiring gravity satellite actual measurement data; the gravity satellite actual measurement data comprises: orbit data, accelerometer data, star sensor data, and inter-satellite ranging data; the track data includes: an initial epoch and an end epoch of the track;
correcting the initial forward track value of the track by using the initial epoch as an integral initial value and adopting a forward numerical integration method to obtain a corrected forward track initial value;
correcting the initial value of the reverse track of the track by using the last epoch as an integral initial value and adopting a reverse numerical integration method to obtain a corrected initial value of the reverse track;
the corrected forward orbit initial value is taken as an integral initial value, and a forward numerical integration method is adopted to obtain a forward reference orbit;
determining a forward reference inter-satellite variability according to the forward reference orbit;
taking the corrected initial value of the reverse track as an integral initial value, and adopting a reverse numerical integration method to obtain a reverse reference track;
determining a reverse reference inter-satellite distance variability according to the reverse reference orbit;
and constructing a satellite gravitational field inversion model according to the forward reference orbit, the forward reference inter-satellite space variability, the reverse reference orbit and the reverse reference inter-satellite space variability.
2. The method for inverting a satellite gravitational field based on 'two-way integration' according to claim 1, wherein the constructing a satellite gravitational field inversion model according to the forward reference orbit, the forward reference inter-satellite variability, the backward reference orbit and the backward reference inter-satellite variability specifically comprises:
constructing a forward observation equation according to the forward reference orbit and the forward reference inter-satellite distance variability;
determining a forward satellite gravitational field inversion model according to the forward observation equation;
constructing a reverse observation equation according to the reverse reference orbit and the reverse reference inter-satellite distance variability;
determining an inverse model of the inverse satellite gravitational field according to the inverse observation equation;
and determining a satellite gravitational field inversion model according to the forward satellite gravitational field inversion model and the reverse satellite gravitational field inversion model.
3. The method for inverting a satellite gravitational field based on "two-way integration" according to claim 1, wherein before obtaining the initial epoch and the final epoch of the orbit in the satellite gravitational field to be inverted, further comprises:
collecting measured data of a gravity satellite in the gravitational field of the satellite to be inverted; the gravity satellite actual measurement data comprises: orbit data, accelerometer data, star sensor data, and inter-satellite ranging data; the track data includes: an initial epoch and an end epoch of the track;
preprocessing the collected measured data of the gravity satellite.
4. A satellite gravitational field inversion system based on "two-way integration", comprising:
the data acquisition module is used for acquiring an initial epoch and an end epoch of an orbit in a gravitational field of a satellite to be inverted;
the model acquisition module is used for acquiring a satellite gravitational field inversion model;
the gravity field inversion module is used for inputting the initial epoch and the final epoch into the satellite gravity field inversion model to finish inversion of the satellite gravity field;
the system further comprises:
the measured data acquisition module is used for acquiring the measured data of the gravity satellite; the gravity satellite actual measurement data comprises: orbit data, accelerometer data, star sensor data, and inter-satellite ranging data; the track data includes: an initial epoch and an end epoch of the track;
the forward track initial value correction module is used for correcting the forward track initial value of the track by using the initial epoch as an integral initial value and adopting a forward numerical integration method to obtain a corrected forward track initial value;
the reverse orbit initial value correction module is used for correcting the reverse orbit initial value of the orbit by adopting a reverse numerical integration method by taking the last epoch as an integral initial value to obtain a corrected reverse orbit initial value;
the forward reference orbit determination module is used for obtaining a forward reference orbit by taking the corrected forward orbit initial value as an integral initial value and adopting a forward numerical integration method;
the forward reference inter-satellite variability determining module is used for determining a forward reference inter-satellite variability according to the forward reference orbit;
the reverse reference orbit determination module is used for obtaining a reverse reference orbit by taking the corrected initial value of the reverse orbit as an integral initial value and adopting a reverse numerical integration method;
the reverse reference inter-satellite distance variable rate determining module is used for determining a reverse reference inter-satellite distance variable rate according to the reverse reference orbit;
the satellite gravitational field inversion model construction module is used for constructing a satellite gravitational field inversion model according to the forward reference orbit, the forward reference inter-satellite variability, the reverse reference orbit and the reverse reference inter-satellite variability.
5. The satellite gravitational field inversion system based on "two-way integration" according to claim 4, wherein said satellite gravitational field inversion model construction module specifically comprises:
a forward observation equation construction unit for constructing a forward observation equation according to the forward reference orbit and the forward reference inter-satellite variability;
the forward satellite gravitational field inversion model determining unit is used for determining a forward satellite gravitational field inversion model according to the forward observation equation;
the reverse observation equation construction unit is used for constructing a reverse observation equation according to the reverse reference orbit and the reverse reference inter-satellite space variability;
the inverse satellite gravitational field inversion model determining unit is used for determining an inverse satellite gravitational field inversion model according to the inverse observation equation;
the satellite gravitational field inversion model determining unit is used for determining a satellite gravitational field inversion model according to the forward satellite gravitational field inversion model and the reverse satellite gravitational field inversion model.
6. The satellite gravitational field inversion system based on "two-way integration" of claim 4, said system further comprising:
the data acquisition module is used for acquiring the gravity satellite actual measurement data in the satellite gravity field to be inverted; the gravity satellite actual measurement data comprises: orbit data, accelerometer data, star sensor data, and inter-satellite ranging data; the track data includes: an initial epoch and an end epoch of the track;
the preprocessing module is used for preprocessing the acquired gravity satellite actual measurement data.
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202010446920.1A CN111708095B (en) | 2020-05-25 | 2020-05-25 | Satellite gravitational field inversion method and system based on bidirectional integration |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202010446920.1A CN111708095B (en) | 2020-05-25 | 2020-05-25 | Satellite gravitational field inversion method and system based on bidirectional integration |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN111708095A CN111708095A (en) | 2020-09-25 |
| CN111708095B true CN111708095B (en) | 2023-07-28 |
Family
ID=72537340
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202010446920.1A Active CN111708095B (en) | 2020-05-25 | 2020-05-25 | Satellite gravitational field inversion method and system based on bidirectional integration |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN111708095B (en) |
Family Cites Families (5)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN102393535B (en) * | 2011-07-20 | 2013-06-05 | 中国科学院测量与地球物理研究所 | Satellite gravity inversion method based on two-star energy interpolation principle |
| CN103076640B (en) * | 2013-01-17 | 2015-03-18 | 中国科学院测量与地球物理研究所 | Method for inverting earth gravitational field by using variance-covariance diagonal tensor principle |
| CN103513294B (en) * | 2013-09-29 | 2016-05-18 | 清华大学 | A kind of low star tracking satellite gravitational field measurement performance Analytic Calculation Method |
| CN108267792B (en) * | 2018-04-13 | 2019-07-12 | 武汉大学 | Building global gravitational field model inversion method |
| CN108845366A (en) * | 2018-06-04 | 2018-11-20 | 中国人民解放军61540部队 | The system of harmonic analysis model and its modeling method of Satellite gravity field tensor diagonal line three-component inverting earth gravitational field |
-
2020
- 2020-05-25 CN CN202010446920.1A patent/CN111708095B/en active Active
Also Published As
| Publication number | Publication date |
|---|---|
| CN111708095A (en) | 2020-09-25 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN108225337B (en) | Star sensor and gyroscope combined attitude determination method based on SR-UKF filtering | |
| CN111156987B (en) | Inertia/astronomy combined navigation method based on residual compensation multi-rate CKF | |
| Searcy et al. | Magnetometer-only attitude determination using novel two-step Kalman filter approach | |
| CN106767837A (en) | Based on the spacecraft attitude method of estimation that volume quaternary number is estimated | |
| CN107729592B (en) | Time-varying structure modal parameter identification method based on generalized subspace tracking | |
| CN113375677B (en) | Spacecraft speed fixing method based on pulsar observation | |
| CN107246883A (en) | A kind of Rotating Platform for High Precision Star Sensor installs the in-orbit real-time calibration method of matrix | |
| CN104121927A (en) | Inertial measurement unit calibration method applicable to low-accuracy no-azimuth-reference single-axis transposition equipment | |
| CN108896456B (en) | Aerosol extinction coefficient inversion method based on feedback type RBF neural network | |
| CN103123487B (en) | A kind of spacecraft attitude determination method | |
| CN111044082B (en) | Gyro error parameter on-orbit rapid calibration method based on star sensor assistance | |
| CN111708095B (en) | Satellite gravitational field inversion method and system based on bidirectional integration | |
| CN109557594B (en) | Gravity reference graph time-varying correction method and system based on gravity abnormal time variation | |
| CN110554443A (en) | Method for determining earth gravity field based on carrier phase observed value and point acceleration method | |
| CN111198402B (en) | Earth gravity field model modeling method based on orbit mask differential operator | |
| CN115235513A (en) | Inertial navigation correction method based on satellite navigation pseudo range and pseudo range rate | |
| Zhao et al. | Unscented Kalman filtering for SINS attitude estimation | |
| CN108981753B (en) | Ground alignment method, system and storage medium based on multi-inertial-unit information constraint | |
| Shou et al. | Micro-satellite attitude determination: using Kalman filtering of magnetometer data approach | |
| CN112013849A (en) | Autonomous positioning method and system for surface ship | |
| CN114608587B (en) | Template-independent single-pulse-satellite spacecraft orbit determination method | |
| Gao et al. | Performance comparison of two filtering approaches for INS rapid transfer alignment | |
| CN116659483A (en) | On-orbit autonomous astronomical navigation method using multiple sensors | |
| Wang et al. | Application of a Sigma-point Kalman filter for alignment of MEMS-IMU | |
| Lijun et al. | Moving horizon estimation and attitude optimization algorithm for MEMS gyroscope array |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant |