CN104965208A - Ephemeris for simulation of satellite navigation signals, and almanac obtaining method - Google Patents

Abstract

The invention provides an ephemeris for simulation of satellite navigation signals, and an almanac obtaining method, and belongs to the technical field of satellite navigation. The ephemeris for the simulation of the satellite navigation signals, and the almanac obtaining method provided by the invention are applied to one from a known almanac or an ephemeris, do not directly obtain the other one through acquisition of satellite orbit information when the other one is needed but obtain the other one through mutual conversion between the two, and can realize accurate conversion between ephemeris/almanac parameters.

Description

A kind of for the ephemeris of satellite navigation signals simulation and the acquisition methods of almanac

Technical field

The invention belongs to technical field of satellite navigation, particularly relate to a kind of ephemeris for satellite navigation signals simulation or almanac conversion method.

Background technology

The position of Real-time Obtaining satellite, speed and clock correction correction information are the preconditions that satellite navigation and location system positions, and the position of satellite, speed and clock correction correction can be calculated by satellite broadcasting ephemeris and corresponding clock correction corrected parameter, therefore, satellite ephemeris and corresponding clock correction parameters are indispensable for satellite navigation system.Broadcast ephemeris is as the one of ephemeris, calculated by the ground handling control center station of satellite navigation system, be transmitted to user by satellite with the form of radio news program, be mainly used in real-time navigation location, wherein comprise 6 orbital trackings, ephemeris corrects parameter with reference to moment and 9 perturbations.Compared with broadcast ephemeris, almanac only provides the primary track parameter of each satellite, and almanac is with reference to moment and two clock error correction numbers, and therefore, the precision of almanac is lower, but effective time increases greatly.In satellite navigation simulation, the acquisition of ephemeris and almanac parameters is very important.At present, the algorithm carrying out ephemeris/almanac matching according to satellite orbit parameter is comparatively perfect, but between ephemeris/almanac, the method for conversion mutually is not yet suggested, and in satellite navigation simulation, needs to obtain ephemeris and almanac parameters, for artificial actual satellite-signal.When adopting the outside mode imported to obtain ephemeris or ephemeris time, if both import simultaneously, needing the unitarity considering effective time, being sometimes not easy to realize, in this case, this method can being utilized to carry out the mutual conversion of ephemeris and almanac.What the present invention relates to is a kind of ephemeris/almanac accurate transformation method for satellite navigation simulation.

Summary of the invention

For solving the problem, the invention provides a kind of ephemeris for satellite navigation signals simulation or almanac conversion method, realizing the accurate transformation of ephemeris/almanac parameters.

Ephemeris for satellite navigation signals simulation of the present invention and the acquisition methods of almanac, be applied in known almanac or ephemeris, needs another time directly do not obtain by gathering satellite-orbit information, and by mutual conversion acquisition therebetween; If A is known terms, B is unknown term, when ephemeris is A, almanac is B, and when almanac is A, ephemeris is B, and described conversion comprises the following steps:

Step 1, the ground utilizing the known terms in m moment to calculate each moment is Satellite position admittedly;

Step 2, sets up state equation X=X (X₀, t₀, t) with observation equation Y=Y (X, t)=Y (X₀, t₀, t);

Wherein, X is the parameter vector of the B of t, X₀for t₀the parameter vector of the B in moment, t₀for the ephemeris reference time, Y is an observation column vector containing m observed quantity, and m is greater than or equals the number of parameters of A, and m observed quantity is Satellite position corresponding to the ground in m moment admittedly;

Step 3, with y and H for intermediate quantity sets up X₀optimal estimationexpression formulautilize the mode of iteration according toobtain the optimal estimation of iteration each timethe optimal estimation of i-th iteration is designated asand willas x_i/0, according to x_i/0=X₀-X_i/0relation obtain X₀;

Wherein, y=Y-Y (X_i/0, t₀, t),x_i/0for intermediate quantity, X_i/0be the X of i-th iteration₀iterative initial value;

The content of iteration comprises:

During the 1st iteration, initial value is X_1/0, according toobtain current iteration x_1/0optimal estimationwherein, X_1/0for setting value;

During the 2nd iteration, initial value isand according toobtain current iteration x_2/0optimal estimation

During i-th iteration, initial valueand according toobtain current iteration x_i/0optimal estimationif current iteration meetsthen stop interative computation;

Wherein, σ_ibe the statistical error of i-th iteration, ε is the minimum arranged according to accuracy requirement;

If the I time stops iteration, then willand X_i/0substitute into x_i/0=X₀-X_i/0obtain X₀;

Step 4, the X that step 3 is obtained₀as the parameter vector value of B, be the ephemeris/almanac parameters after conversion.

Embodiment

Need to obtain ephemeris parameter or almanac parameters in satellite navigation simulation, ephemeris parameter is relatively more comprehensive, thus accuracy is high, but its effective time is short; Almanac parameters then number ratio is less, thus accuracy is not high, but its effective time is longer.Ephemeris parameter or almanac parameters respectively have its merits and demerits, need to go to obtain target component according to the actual demand of satellite navigation simulation.When needing another kind of parameter in prior art in known a kind of parameter, need to adopt the outside mode imported to obtain, need the unitarity of the effective time of consideration two kinds of parameters during importing, inconvenience realizes; Also likely exterior storage is complete, causes the parameter that can not obtain needs.And also need ephemeris and almanac when dummy source generating literary composition simultaneously, and for some reason may ephemeris or almanac obtain less than, at this time can only remove to obtain by multiple acquisition mode the ephemeris or almanac that need.Just based on this technical matters, the application is when known parameters, and the technical scheme proposing the application converts the parameter of needs to according to known parameters, and method is simple, be convenient to realize, and accuracy is higher.

Matching ephemeris or ephemeris time, need not many groups satellite orbit data in the same time, therefore, need to calculate respective satellite orbital data according to certain time interval in almanac or ephemeris effective time.

When known ephemeris information, with the ephemeris reference moment for benchmark, reasonable time interval should be chosen, calculates ephemeris with reference to the satellite orbit data information in a period of time before and after the moment.It should be noted that during the access time of interval that the time that the satellite orbital position calculated is corresponding needs to comprise the value of almanac with reference to the moment in order to make fitting result more accurate.For meeting this requirement, for the GPS of the U.S. (GPS), consider that its ephemeris every two hours upgrades once, ephemeris is with reference to moment t_oe(in seconds) be 7200 integral multiple, and the reference moment t of almanac_oa(in seconds) be 4096 integral multiple, comprise t to make selected satellite orbit data_oethe satellite orbit data information in moment, the selected time interval (in seconds) should be the minimum common factor of 4096 and 7200 or can be divided exactly by its minimum common factor, and specific implementation is shown in step one.

Step one, according to known ephemeris parameter utilize ephemeris user algorithm to calculate be Satellite position admittedly, or according to known almanac parameters utilize almanac user algorithm to calculate be Satellite position admittedly;

If step 1 A is ephemeris, B is almanac, then utilize ephemeris user algorithm to calculate be Satellite position admittedly:

Step 101: the satellite mean angular velocity n inscribed when calculating ephemeris reference₀:

$n_0=\sqrt{\frac{\mathrm{\μ}}{A^3}}$(formula 1)

Wherein, μ is Gravitational coefficient of the Earth, and A is semi-major axis of orbit; For different satellite navigation systems, μ needs to correspond to different coordinate systems and carries out value, and in the Beidou satellite navigation system of China, μ is the Gravitational coefficient of the Earth of CGCS2000 coordinate system, and value as shown in Equation 2;

μ=3.986004418 × 10¹⁴m³/ s²(formula 2)

$A={\left(\sqrt{A}\right)}^2$(formula 3)

Step 102, the satellite mean angular velocity n utilizing step 101 to obtain₀ephemeris is revised with reference to moment satellite mean angular velocity, obtains revised mean angular velocity n:

N=n₀+ Δ n (formula 5)

Wherein, Δ n is the parameter value announced in ephemeris parameter;

Step 103, calculates current time t and ephemeris reference time t_oemistiming t_k:

T_k=t-t_oe(formula 4)

Step 104, utilizes the mean anomaly M of reference time in ephemeris parameter₀, step 102 obtain revised mean angular velocity n and step 103 obtain mistiming t_kcalculate the mean anomaly M of current time_k:

M_k=M₀+ nt_k(formula 6)

Step 105, the mean anomaly M that integrating step 104 obtains_kthe method of iteration is utilized to calculate the eccentric anomaly E of current time_k:

E_k=M_k-e sin E_k(formula 7)

This alternative manner is prior art, and iterative formula is E_{k (i+1)}=M_k-e sin E_{k (i)}, iteration initial value E_{k (1)}given, if current E_{k (i)}do not meet setting demand, then ask E according to iterative formula_{k (i+1)}if, E_{k (i+1)}acquisition meets setting demand, then stop iteration, current E_{k (i+1)}be the eccentric anomaly E of current time_k.

E in above formula is the excentricity in ephemeris parameter.

Step 106, utilizes the eccentric anomaly E that step 105 obtains_kcalculate the true anomaly f of current time_k:

$f_k=\mathrm{ta}{n}^{-1}\left\{\frac{\sqrt{1-e^2}\sin E_k/(1-e\cos E_k)}{(\cos E_k-e)/(1-e\cos E_k)}\right\}$(formula 8)

Step 107, the true anomaly f utilizing step 106 to obtain_kcalculate the latitude argument of current time

(formula 9)

W in above formula is the argument of perigee in ephemeris parameter.

Step 108, the latitude argument utilizing step 107 to obtainconsider that second order zonal harmonics perturbation corrects latitude argument, satellite radius vector and orbit inclination, obtain latitude argument correction member δ μ_k, radial correction member δ r_kand orbit inclination correction member δ i_kas shown in the formula:

(formula 10)

C in above formula_us, C_ucbe respectively the amplitude of the sine of latitude argument, cosine mediation correction member, C_rs, C_rcbe respectively the amplitude of the sine of orbit radius, cosine mediation correction member, C_is, C_icbe respectively the amplitude of the sine of orbit inclination, cosine mediation correction member.

Step 109, the latitude argument correction member δ μ utilizing step 108 to obtain_k, radial correction member δ r_kand orbit inclination correction member δ i_kcalculate the latitude argument μ after correcting_k, satellite radius vector r_kwith orbit inclination i_k:

(formula 11)

IDOT in above formula is the orbit inclination rate of change in ephemeris parameter, i₀for the orbit inclination announced in ephemeris.

Step 110, the latitude argument μ after the correction utilizing step 109 to obtain_k, satellite radius vector r_kcalculate the coordinate (x' of satellite in orbital coordinate system_k, y'_k):

x′_k＝r_kcosμ_k

Y '_k=r_ksin μ_k(formula 12)

Step 111, the mistiming t utilizing step 103 to obtain_kthe geodetic longitude Ω of calculating observation moment ascending node_k:

${\mathrm{\Ω}}_k={\mathrm{\Ω}}_0+(\stackrel{\·}{\mathrm{\Ω}}-{\stackrel{\·}{\mathrm{\Ω}}}_e)t_k-{\stackrel{\·}{\mathrm{\Ω}}}_e{t}_{\mathrm{oe}}$(formula 13)

In above formula, Ω₀for ephemeris in ephemeris parameter is with reference to the right ascension of ascending node in moment,for the right ascension of ascending node rate of change in ephemeris parameter,for earth rotation speed, for different satellite navigation systems,need to correspond to different coordinate systems and carry out value, in the Beidou satellite navigation system of China,for the earth rotation speed of CGCS2000 coordinate system, value as shown in Equation 14.

${\stackrel{\·}{\mathrm{\Ω}}}_e=7.2921150\×10^{-5}\mathrm{rad}/s$(formula 14)

Step 112, utilizes the orbit inclination i that step 109 obtains_k, the satellite that obtains of step 110 is at the coordinate (x' of orbital coordinate system_k, y'_k) and step 111 obtain geodetic longitude Ω_kcalculating ground is Satellite position (x admittedly_k, y_k, z_k):

X_k=x '_kcos Ω_k-y '_ksin Ω_kcos i_ky_k=x '_ksin Ω_k+ y '_kcos Ω coi_k(formula 15) z_k=y '_ksin i_k

When known almanac information, should with the almanac reference moment for benchmark, choose reasonable time interval, calculate almanac with reference to the satellite orbit data information in a period of time before and after the moment, require similar to known ephemeris time for choosing of the time interval, the time that namely calculated satellite orbital position is corresponding needs to comprise the value in ephemeris reference moment.Specific implementation step is shown in step 2.

Step 2, if A is almanac, B is ephemeris, then utilize almanac user algorithm to calculate be Satellite position admittedly:

Step 201: the satellite mean angular velocity n inscribed when calculating almanac reference₀:

$n_0=\sqrt{\frac{\mathrm{\μ}}{A^3}}$(formula 16)

Wherein, μ is Gravitational coefficient of the Earth, and A is semi-major axis of orbit; For different satellite navigation systems, μ needs to correspond to different coordinate systems and carries out value, and in the Beidou satellite navigation system of China, μ is the Gravitational coefficient of the Earth of CGCS2000 coordinate system, and value as shown in Equation 2.

Step 202, calculates current time t and almanac reference time t_oamistiming t_k^*:

T_k^*=t-t_oa(formula 17)

Step 203, utilizes the mean anomaly M of reference time in almanac parameters₀^*, step 201 obtain mean angular velocity n₀and the mistiming t that step 202 obtains_k^*calculate the mean anomaly M of current time_k^*:

${M_k}^{*}={M_0}^{*}+n_0^{*}{t_k}^{*}$(formula 18)

Step 204, the mean anomaly M that integrating step 203 obtains_k^*the method of iteration is utilized to calculate the eccentric anomaly E of current time_k^*:

E_k^*=M_k^*-e^*sinE_k^*(formula 19)

E in above formula^*for the excentricity in almanac parameters.

Step 205, utilizes the eccentric anomaly E that step 204 obtains_k^*calculate the true anomaly of current time

$f_k^{*}=\mathrm{ta}{n}^{-1}\left\{\frac{\sqrt{1-{\left(e^{*}\right)}^2}\sin {E_k}^{*}/(1-e^{*}\cos {E_k}^{*})}{(\cos {E_k}^{*}-e^{*})/(1-e^{*}\cos {E_k}^{*})}\right\}$(formula 20)

Step 206, utilizes the true anomaly that step 205 obtainscalculate the latitude argument of current time

(formula 21)

W in above formula^*for the argument of perigee in almanac parameters.

Step 207, calculates satellite radius vector

$r_k^{*}=A(1-e^{*}\cos E_k^{*})$(formula 22)

Step 208, calculates the orbit inclination i of almanac with reference to the moment:

$i=i_0^{*}+\mathrm{\δ}$(formula 23)

δ i in above formula is the reduction at the track reference inclination angle of almanac reference time in almanac parameters,for the orbit inclination announced in almanac.

Step 209, the latitude argument utilizing step 206 to obtainwith the satellite radius vector that step 207 obtainscalculate the coordinate (x' of satellite in orbital coordinate system_k, y'_k):

(formula 24)

Step 210, the geodetic longitude of the calculating observation moment ascending node utilizing step 209 to obtain

${\mathrm{\Ω}}_k^{*}={\mathrm{\Ω}}_0^{*}+({\stackrel{\·}{\mathrm{\Ω}}}^{*}-{\stackrel{\·}{\mathrm{\Ω}}}_e^{*})t_k^{*}-{\stackrel{\·}{\mathrm{\Ω}}}_e^{*}{t}_{\mathrm{oa}}$(formula 25)

In above formula,for almanac in almanac parameters is with reference to the right ascension of ascending node in moment,for the right ascension of ascending node rate of change in almanac parameters,for earth rotation speed, for different satellite navigation systems,need to correspond to different coordinate systems and carry out value, in the Beidou satellite navigation system of China,for the earth rotation speed of CGCS2000 coordinate system, value as shown in Equation 26.

${\stackrel{\·}{\mathrm{\Ω}}}_e^{*}=7.2921150\×10^5\mathrm{rad}/s$(formula 26)

Step 211, the geodetic longitude that the orbit inclination i step 210 utilizing step 208 to obtain obtainswith the coordinate (x' that step 209 obtains_k, y'_k) be Satellite position (x admittedly with calculating_k, y_k, z_k):

$\begin{array}{c}{x}_k=x_k^{\′}\cos {\mathrm{\Ω}}_k^{*}-y_k^{\′}\sin {\mathrm{\Ω}}_k^{*}\cos i\\ y_k=x_k^{\′}\sin {\mathrm{\Ω}}_k^{*}+y_k^{\′}\cos {\mathrm{\Ω}}_k^{*}\cos i\\ z_k=y_k^{\′}\sin i\end{array}$(formula 27)

Utilize the not many groups satellite orbit data in the same time obtained by ephemeris computation, can almanac parameters matching be carried out; Utilize the not many groups satellite orbit data in the same time calculated by almanac, can ephemeris parameter matching be carried out.Detailed process is shown in step 2.

Step 2, ephemeris/almanac parameters matching

Step 21, the ground according to obtaining in step one is Satellite position (x' admittedly_k, y'_k) and ephemeris/almanac parameters set up state equation (formula 28) and observation equation (formula 29):

X=X (X₀, t₀, t) (formula 28)

Y=Y (X, t)=Y (X₀, t₀, t) (formula 29)

Wherein, X is quantity of state, be again t corresponding for solve for parameter vector X=X (X₀, t₀, t), X=X (X₀, t₀, t) represent that X and t is relevant, X₀for t₀corresponding solve for parameter vector, t₀for the ephemeris reference time, and t₀∈ t.Y is an observation column vector (m is greater than and equals 15) containing m observed quantity, and an observed quantity corresponds to the location components (x of Navsat t in ground is admittedly_k, y_k, z_k), represent Y and X Y=Y (X, t) relevant to t again, and due to t₀∈ t, so Y=Y (X, t)=Y (X₀, t₀, t).Because state equation and observation equation are nonlinear equation, the valuation problem of ephemeris/almanac parameters is the Least Squares Estimating problem of nonlinear system, and row iteration of nonlinear equation linearization being gone forward side by side solves.

When fitting parameter is ephemeris parameter, X₀for ephemeris parameter, as shown in Equation 30; .

$X_0={(\sqrt{A},e,i,{\mathrm{\Ω}}_0,w,M_0,\mathrm{\Δn},\stackrel{\·}{\mathrm{\Ω},\mathrm{IDOT},C_{\mathrm{us}},C_{\mathrm{uc}},C_{\mathrm{is}},C_{\mathrm{ic}},C_{\mathrm{rs}},C_{\mathrm{rc}}})}^T$(formula 30)

When fitting parameter is almanac parameters, X₀for almanac parameters, as shown in Equation 31; t₀for the almanac reference time.

$X_0={(\sqrt{A},e^{*}\mathrm{\δi},{\mathrm{\Ω}}_0^{*},w^{*},{M_0}^{*},{\stackrel{\·}{\mathrm{\Ω}}}^{*},t_{\mathrm{oa}})}^T$(formula 31)

Step 22, by formula 29 at X_i/0place launches:

$Y=Y(X_{i/0},t_0,t)+{\left(\frac{\∂Y}{\∂X_0}\right)}_{X_0=X_{i/0}}(X_0-X_{i/0})+O\left[{(X_0-X_{i/0})}^2\right]$(formula 32)

X_i/0for X₀for iterative initial value, it is setting value.Ο () is higher order term more than second order and second order.

Step 23, order:

X₀=X₀-X_i/0(formula 33)

Y=Y-Y (X_i/0, t₀, t) (formula 34)

$H={\left(\frac{\∂Y}{\∂X_0}\right)}_{X_0=X_{i/0}}{\left(\frac{\∂Y}{\∂X}\frac{\∂X}{\∂X_0}\right)}_{X_0=X_{i/0}}$(formula 35)

Step 24, the higher order term that the formula 32 in step 22 that (formula 33) (formula 34) in step 23 and (formula 35) substituted into is omitted in formula 32 more than second order can obtain:

Y=Hx₀(formula 36)

Step 25, solves formula 36 according to Least Squares Estimating principle, can obtain x₀optimal estimation:

${\hat{x}}_0={\left(H^{T}H\right)}^{-1}{H}^{T}y$(formula 37)

Step 26, during the 1st iteration, initial value is X_1/0, obtain x according to formula 37₀optimal estimationduring the 2nd iteration, initial valueand obtain x according to formula 37₀optimal estimationduring i-th iteration, initial valueand obtain x according to formula 37₀optimal estimationif now meet formula 39 to stop calculating, and willas solve for parameter vector value X_ovalue

$\frac{|{\mathrm{\&sigma;}}_{i+1}-{\mathrm{\&sigma;}}_i|}{{\mathrm{\&sigma;}}_i}<\mathrm{\&epsiv;}$(formula 39)

In above formula, σ is statistical error in iterative process, and ε is the minimum (usually getting 0.01) arranged according to accuracy requirement:

$\mathrm{\&sigma;}=\sqrt{\frac{{(y-H{x}_0)}^T(y-H{x}_0)}{m-n}}$(formula 40)

In above formula, n is the dimension of quantity of state.The fitting result that step 27 obtains-solve for parameter vector value X_obe final required ephemeris/almanac parameters.

Certainly; the present invention also can have other various embodiments; when not deviating from the present invention's spirit and essence thereof; those of ordinary skill in the art are when making various corresponding change and distortion according to the present invention, but these change accordingly and are out of shape the protection domain that all should belong to the claim appended by the present invention.

Claims (4)

1. the ephemeris for satellite navigation signals simulation and the acquisition methods of almanac, it is characterized in that, be applied in known almanac or ephemeris, needs another time directly do not obtain by gathering satellite-orbit information, and by mutual conversion acquisition therebetween; If A is known terms, B is unknown term, when ephemeris is A, almanac is B, and when almanac is A, ephemeris is B, and described conversion comprises the following steps:

Step 1, the ground utilizing the known terms in m moment to calculate each moment is Satellite position admittedly;

Step 2, sets up state equation X=X (X₀, t₀, t) with observation equation Y=Y (X, t)=Y (X₀, t₀, t);

Wherein, X is the parameter vector of the B of t, X₀for t₀the parameter vector of the B in moment, t₀for the ephemeris reference time, Y is an observation column vector containing m observed quantity, and m is greater than or equals the number of parameters of A, and m observed quantity is Satellite position corresponding to the ground in m moment admittedly;

Step 3, with y and H for intermediate quantity sets up X₀optimal estimationexpression formulautilize the mode of iteration according toobtain the optimal estimation of iteration each timethe optimal estimation of i-th iteration is designated asand willas x_i/0, according to x_i/0=X₀-X_i/0relation obtain X₀;

Wherein, y=Y-Y (X_i/0, t₀, t),x_i/0for intermediate quantity, X_i/0be the X of i-th iteration₀iterative initial value;

The content of iteration comprises:

During the 1st iteration, initial value is X_1/0, according toobtain current iteration x_1/0optimal estimationwherein, X_1/0for setting value;

During the 2nd iteration, initial valueand according toobtain current iteration x_2/0optimal estimation

During i-th iteration, initial valueand according toobtain current iteration x_i/0optimal estimationif current iteration meetsthen stop interative computation;

Wherein, σ_ibe the statistical error of i-th iteration, ε is the minimum arranged according to accuracy requirement;

If the I time stops iteration, then willand X_i/0substitute into x_i/0=X₀-X_i/0obtain X₀;

Step 4, the X that step 3 is obtained₀as the parameter vector value of B, be the ephemeris/almanac parameters after conversion.

2. as claimed in claim 1 for ephemeris or the almanac conversion method of satellite navigation signals simulation, it is characterized in that, if A is ephemeris, then by ephemeris reference moment t in the parameter of ephemeris_oebe the minimum common factor of 4096 and 7200 or can be divided exactly by its minimum common factor.

If A is almanac, then almanac reference time t in the parameter of almanac_oabe the minimum common factor of 4096 and 7200 or can be divided exactly by its minimum common factor.

3., as claimed in claim 1 for ephemeris or the almanac conversion method of satellite navigation signals simulation, it is characterized in that, ε is 0.01.

4., as claimed in claim 1 for ephemeris or the almanac conversion method of satellite navigation signals simulation, it is characterized in that,