CN102998681A - High-frequency clock error estimation method of satellite navigation system - Google Patents

Abstract

The invention discloses a high-frequency clock error estimation method of a satellite navigation system. The method includes acquiring observed data of a tracking station; testing the reliability of the observed data; rejecting outliers from the observed data after reliability test by using linear combination of the observation quantity, determining an initial ambiguity, detecting and repairing the cycle slip, and obtaining codes and phase observed values; acquiring a satellite obit parameter, an earth polar motion parameter, a troposphere parameter and a station coordinate parameter; obtaining clock error estimated values of a satellite and a receiver according to original observed data, the satellite obit parameter, the earth polar motion parameter, the troposphere parameter and the station coordinate parameter; and interpolating the obtained clock error estimation by using the code and the phase observed value to obtain the high-frequency clock error estimation of the satellite navigation system. The method has the advantages of high reliability and accuracy and capability of improving the high-frequency clock error estimation accuracy of the satellite navigation system.

Description

A kind of high frequency clock bias estimation method of satellite navigation system

Technical field

The present invention relates to satellite navigation positioning technical field, particularly relate to a kind of high frequency clock bias estimation method of satellite navigation system.

Background technology

Satellite clock correction refers to satellite clock time and the difference of navigational system between the standard time that the satellite clock frequency drift causes.The navigational system observed quantity be take the clock of satellite and receiver frequently signal obtain as benchmark, user's bearing accuracy and its precision of obtainable satellite clock correction and track closely related, and the precision closely related (also can be called a kind of coupled relation) that the estimated accuracy of clock correction and track are estimated.

With GPS(Global Positioning System, GPS) Navsat is example, at present the gps satellite clock correction precision as the broadcast ephemeris of navigation message broadcast only is 7ns, can't satisfy the requirement of decimeter grade consumer positioning, and IGS(The International GNSS Service, international GPS service tissue) the quick and supper-fast clock correction of gps satellite clock correction of centrales prerequisite confession only has the 15min sampling interval, and the gps satellite clock correction of the 30s that the CODE center provides also must postpone 14-17 hour, sometimes is difficult to satisfy the requirement of engineering.In addition on the one hand, if want to obtain high frequency, high-precision satellite clock correction product, observation data is more so in theory, the satellite clock correction precision of resolving is higher, if but the ground observation data are too many, the so corresponding parameter of resolving also will be more, and then the travelling speed of computing machine is had higher requirement.Most important one side is if the IGS center because some natures or politics cause can't provide some gps satellite product, will certainly exert an influence to the precise orbit determination of picture based on spaceborne GPS low orbit satellite so.

Summary of the invention

The technical problem to be solved in the present invention is a kind of high frequency clock bias estimation method of satellite navigation system, and the satellite clock correction that can effectively solve high-precision 30s or higher frequency is obtained hard problem.

For solving the problems of the technologies described above, the invention provides a kind of high frequency clock bias estimation method of satellite navigation system, comprising:

S101 obtains the observation data of tracking station;

S102, the check observation the reliability of the data;

S103 to utilize the linear combination of observed quantity through the observation data of certificate authenticity, rejects wild value, determines initial blur level, surveys and repair cycle slip, through obtaining code and phase observations value after the smoothing processing;

S104 obtains satellite orbit parameter, polar motion of globe parameter, troposphere parameter, survey station coordinate parameters;

S105, orbit parameter, polar motion of globe parameter, troposphere parameter, the survey station coordinate parameters of the satellite that the original observed data that step S103 is obtained and step S104 obtain; Be brought in the non-poor pseudorange and phase observations value error equation of having eliminated ionosphere effect, the selection reference clock is treated and is estimated that clock imposes restriction, and carries out unified adjustment and processes, and obtains the clock bias estimation value of satellite and receiver;

S106 utilizes code and phase observations value among the step S103, and the clock bias estimation that step S105 is obtained carries out interpolation, obtains the high frequency clock bias estimation of satellite navigation system.

Further, eliminated non-poor phase place and the Pseudo-range Observations of ionosphere effect, the observed reading error equation is as follows:

$v_{k,\mathrm{\Φ}}^j\left(i\right)={\mathrm{\Δt}}_k\left(i\right)-{\mathrm{\Δt}}^j\left(i\right)+{\mathrm{\ρ}}_k^j\left(i\right)/C+\mathrm{\δ}{\mathrm{\ρ}}_{k,\mathrm{trop}}^j\left(i\right)/C+\mathrm{\λ}\·N_k^j/C+{\mathrm{\ϵ}}_{k,\mathrm{\Φ}}^j\left(i\right)-\mathrm{\λ}\·{\mathrm{\Φ}}_k^j\left(i\right)/C---\left(4\right)$

$v_{k,p}^j\left(i\right)={\mathrm{\Δt}}_k\left(i\right)-{\mathrm{\Δt}}^j\left(i\right)+{\mathrm{\ρ}}_k^j\left(i\right)/C+\mathrm{\δ}{\mathrm{\ρ}}_{k,\mathrm{trop}}^j\left(i\right)/C+{\mathrm{\ϵ}}_{k,p}^j\left(i\right)-P_k^i\left(i\right)/C---\left(5\right)$

In formula (4), (5): k is for surveying station number, and j is for defending asterisk, and i is corresponding epoch of observation, and C is the light velocity in the vacuum; Δ t _k(i) be receiver clock correction, Δ t ^j(i) be satellite clock correction,

Be the tropospheric delay impact;

Be not modeled error effect;

Be the blur level parameter;

Be pseudorange and the phase combination observed reading without ionosphere effect; Be its observational error, λ is the wavelength without ionosphere combination observation value; For the satellite position of signal x time to the geometric distance between the signal receiver location time of reception.

Further, the pseudorange of clock correction and phase error equation can be expressed as:

$V\left(i\right)=\left(\begin{array}{c}{V}_p\left(i\right)\\ V_{\mathrm{\Φ}}\left(i\right)\end{array}\right)=\left(\begin{array}{c}{A}_p\\ A_{\mathrm{\Φ}}\end{array}\right)X\left(i\right)+\left(\begin{array}{c}{L}_p\left(i\right)\\ L_{\mathrm{\Φ}}\left(i\right)\end{array}\right)---\left(6\right)$

In the formula (6): $A_p=\left[\begin{array}{cc}\underset{(m\·n)\×(2m+n)}{A}& \underset{(m\·n)\×(m\·n)}{0}\end{array}\right];$ $A_{\mathrm{\Φ}}=\left[\begin{array}{cc}\underset{(m\·n)\×(2m+n)}{A}& \underset{(m\·n)\×(m\·n)}{\mathrm{\λ}/C\·I}\end{array}\right];$ If I is n rank unit matrixs, A is the design matrix when determining receiver and satellite clock correction, be (mn) (m+n) matrix, V (i) is the observational error of i observation position, L (i) is the free term of i observation position, Φ and p represent respectively phase place and pseudorange observation, and X (i) is i parameter to be asked, V (i), L _p(i), L _Φ(i) be known quantity; M base station quantity, n is the quantity of synchronous tracking satellite.

Further, utilize a receiver clock in the tracking network as reference clock, the clock correction precision of reference clock is better than 10 ^-6S; The relative clock correction of all satellites is zero-mean in the same navigational system.

Further, reference clock is hydrogen atomic clock.

Further, step S106 specifically comprises:

S1061 estimates based on the clock that the code observed reading is carried out;

S1062 estimates based on the clock that the phase observations value is carried out;

S1063 carries out the clock combination, obtains the high frequency clock bias estimation of Navsat and receiver.

Beneficial effect of the present invention is as follows:

The present invention can improve the estimated accuracy of satellite navigation system high frequency clock correction, has reliability height, precision advantages of higher.

Description of drawings

Fig. 1 is the process flow diagram of the high frequency clock bias estimation method of a kind of satellite navigation system in the embodiment of the invention.

Embodiment

Below in conjunction with accompanying drawing and embodiment, the present invention is further elaborated.Should be appreciated that specific embodiment described herein only in order to explain the present invention, does not limit the present invention.

To the estimation of the satellite clock correction difference by the parameter processing mode, two kinds of thinkings can be arranged: the one, timing parameter and other parameter, be put into one such as satellite orbit parameter, survey station coordinate, troposphere parameter etc., simultaneously estimation; Another kind is to resolve first non-epoch of the parameter such as satellite orbit parameter, survey station coordinate parameters, they is calculated and fix, and then resolves parameter epoch, by the number that reduces parameter epoch, resolves at last the clock correction parameter.For these two kinds of methods, if the less the first parameter that adopts of parameter epoch, the second method substep, carry out stage by stage resolving of parameter, can be used for resolving the situation of many epoch of parameter, and the method can be more effective in efficient, and therefore the embodiment of the invention adopts second method when the clock correction parameter calculation.

Because being the relative time between survey station and the satellite, the GPS observed reading postpones.Therefore, can not determine simultaneously all satellites and receiver clock correction, must first by the clock correction (receiver clock or satellite clock) of fixing a certain reference clock, then determine the relative clock correction of other receiver and satellite, but must guarantee that the clock correction precision of reference clock is better than 10 ^-6S.Clock correction and absolute clock correction are of equal value for user's positioning result relatively, and namely the systematic bias of clock correction can be absorbed by receiver user clock correction in user's location model fully relatively, and does not affect user's bearing accuracy.

As shown in Figure 1, the embodiment of the invention relates to a kind of high precision, high reliability high frequency Navsat clock bias estimation method of satellite navigation system, comprising:

Step S101 obtains tracking station's observation data;

The original observed data that ground tracking station provides is obtained in download, comprises the broadcast ephemeris of RINEX observation data and satellite.

Step S102, the check observation the reliability of the data;

Reliability index is applicable to yard observed reading and the phase observations value is tested, namely

$\mathrm{\γ}=\underset{i=1}{\overset{v}{\mathrm{\Σ}}}\underset{j=1}{\overset{v}{\mathrm{\Σ}}}{\left(\mathrm{PR}\right)}_{\mathrm{ij}}/\underset{i=1}{\overset{v}{\mathrm{\Σ}}}\underset{j=1}{\overset{v}{\mathrm{\Σ}}}{P}_{\mathrm{ij}}---\left(1\right)$

In the formula (1), γ is that system is in the degree of reiability of i observation position; For the code observed reading, this moment v=1, γ is larger, reliability is better, γ is less, reliability is poorer; P is the weight of observation battle array, and R is orthogonal intersection cast shadow matrix.If the component of redundant observation is approximately equal to zero, the reliability of i observation position is very poor so, and the redundant observation component is less, and system is poorer in the reliability of i observation position; J is the satellite numbering.Internal Reliability index MDB and extrinsic reliabilities MDE can be expressed as

${\left|\mathrm{MDB}\right|}_{c\left(p\right)}={\mathrm{\λ}}_0^{\frac{1}{2}}{\mathrm{\σ}}_0{\left[\mathrm{\γ}\underset{i=1}{\overset{v}{\mathrm{\Σ}}}\underset{j=1}{\overset{v}{\mathrm{\Σ}}}{P}_{\mathrm{ij}}\right]}^{-\frac{1}{2}}---\left(2\right)$

MDE _c(p)=λ ₀(γ ^-1-1) （3）

In formula (2), (3), c (p) expression code or phase observations value; λ ₀Be the non-centrality parameter under certain Probability Condition, it is the function of abandoning true probability, the power of test and model parameter (rough error) number to be detected.These three parameters are once selecting λ ₀It is exactly a constant; Usually getting these three parameters is 0.001,0.80 and 1, then λ ₀Be 17.07.When the reliability of the different observation positions that compare same system, only be concerned about their relative ratio, constant λ ₀Inoperative.On the other hand, the selection of abandoning true probability and power of test parameter generally is information rule of thumb, and certain subjectivity and arbitrariness are arranged.σ ₀Be the middle error of observed reading, v is the Continuous Observation number to a certain amount.P, R implication are the same.

Step S103, the data pre-service obtains code and phase observations value;

Step S102 has been carried out the data of certificate authenticity, if the non-constant of reliability is done to reject and processed, threshold value is selected according to the actual observation experience.The data that obtain for the S102 step are further utilized the linear combination of observed quantity, reject wild value, determine initial blur level, survey and also repair cycle slip, through obtaining code and phase observations value after the smoothing processing.By the satellite clock correction in Pseudo-range Observations and the broadcast ephemeris, adopt the least-squares estimation receiver clock correction of each epoch.

Step S104, the track of satellite, the isoparametric calculating of earth rotation;

Utilize the pretreated code of step S103 and phase observations value to form difference information, find the solution satellite orbit parameter, polar motion of globe parameter, troposphere parameter, survey station coordinate parameters etc.Can certainly obtain from the IGS center forecast part of supper-fast track product, can find the solution orbit parameter at that time, directly arrive next step.

Step S105, the satellite clock bias estimation of low sampling rate;

Step S103 is processed orbit parameter, polar motion of globe parameter, troposphere parameter, the survey station coordinate parameters of the satellite that the original observed data obtain and step S104 obtain and be used as given value, be brought in the non-poor pseudorange and phase observations value error equation of having eliminated ionosphere effect, the selection reference clock estimates that clock imposes restriction in order to treat, carry out unified adjustment and process, obtain the clock bias estimation value of high-precision satellite and receiver.

In Clock Bias is estimated, general non-poor phase place and the Pseudo-range Observations of having eliminated ionosphere effect that adopt, the observed reading error equation is as follows:

$v_{k,\mathrm{\Φ}}^j\left(i\right)={\mathrm{\Δt}}_k\left(i\right)-{\mathrm{\Δt}}^j\left(i\right)+{\mathrm{\ρ}}_k^j\left(i\right)/C+\mathrm{\δ}{\mathrm{\ρ}}_{k,\mathrm{trop}}^j\left(i\right)/C+\mathrm{\λ}\·N_k^j/C+{\mathrm{\ϵ}}_{k,\mathrm{\Φ}}^j\left(i\right)-\mathrm{\λ}\·{\mathrm{\Φ}}_k^j\left(i\right)/C---\left(4\right)$

$v_{k,p}^j\left(i\right)={\mathrm{\Δt}}_k\left(i\right)-{\mathrm{\Δt}}^j\left(i\right)+{\mathrm{\ρ}}_k^j\left(i\right)/C+\mathrm{\δ}{\mathrm{\ρ}}_{k,\mathrm{trop}}^j\left(i\right)/C+{\mathrm{\ϵ}}_{k,p}^j\left(i\right)-P_k^i\left(i\right)/C---\left(5\right)$

In formula (4), (5): k is for surveying station number, and j is for defending asterisk, and i is corresponding epoch of observation, and C is the light velocity in the vacuum;

Be receiver clock correction, Δ t ^j(i) be satellite clock correction,

Be the tropospheric delay impact;

Be not modeled error effect;

Be the blur level parameter;

Be pseudorange and the phase combination observed reading without ionosphere effect;

Be its observational error, λ is the wavelength without ionosphere combination observation value;

For the satellite position of signal x time to the geometric distance between the signal receiver location time of reception.

If m base station followed the tracks of n satellite synchronously, the pseudorange of clock correction and phase error equation can be expressed as:

$V\left(i\right)=\left(\begin{array}{c}{V}_p\left(i\right)\\ V_{\mathrm{\Φ}}\left(i\right)\end{array}\right)=\left(\begin{array}{c}{A}_p\\ A_{\mathrm{\Φ}}\end{array}\right)X\left(i\right)+\left(\begin{array}{c}{L}_p\left(i\right)\\ L_{\mathrm{\Φ}}\left(i\right)\end{array}\right)---\left(6\right)$

In the formula (6):

$A_p=\left[\begin{array}{cc}\underset{(m\·n)\×(2m+n)}{A}& \underset{(m\·n)\×(m\·n)}{0}\end{array}\right];$ $A_{\mathrm{\Φ}}=\left[\begin{array}{cc}\underset{(m\·n)\×(2m+n)}{A}& \underset{(m\·n)\×(m\·n)}{\mathrm{\λ}/C\·I}\end{array}\right];$ If I is n rank unit matrixs, A is the design matrix when determining receiver and satellite clock correction, be (mn) (m+n) matrix, V (i) is the observational error of i observation position, L (i) is the free term of i observation position, Φ and p represent respectively phase place and pseudorange observation, and X (i) is i parameter to be asked, V (i), L _p(i), L _Φ(i) be known quantity.If find the solution satellite clock correction parameter take following formula as observation equation, normal equation is unusual, in order to find the solution the clock correction parameter, must introduce a reference clock, ask other receiver clock and satellite clock with respect to the clock correction of reference clock.When introducing reference clock, the present invention is directed to ground receiver clock and satellite clock and adopted two kinds of measures: (1) utilizes a receiver clock in the tracking network as reference clock, but prerequisite is to guarantee reference clock reliable and stable (for example hydrogen atomic clock).(2) the relative clock correction of all satellites is the zero-mean condition in the same navigational system of supposition.This spline equation just can be separated the relative clock correction of satellite.

Step S106, high precision, high frequency satellite clock correction are obtained;

The code and the phase observations value that re-use among the step S103 are carried out interpolation to low frequency satellite clock correction.Common interpolating method is owing to be subject to the impact of data time sequence itself, and the interpolation precision is relatively relatively poor.In order to obtain higher interpolation precision, this step relates to a kind of method of special satellite clock correction interpolation, and this algorithm is a kind of high precision and efficient clock correction interpolating method; Have following characteristics: 1, use the method to carry out satellite clock correction interpolation and can make the clock correction result of each satellite more stable, the saltus step of very large or even several ns can not occur, can improve the interpolation precision; 2, for the decline that has than the precision after the clock correction interpolation of small part satellite a little, but the magnitude of impact is smaller, and is generally high than using general interpolating method precision; 3, along with the increase of website number, clock correction interpolation as a result trend can be smooth-out, and corresponding precision also can improve, and the steadily saltus step of relatively large level can not occur.

Take gps satellite as example, in order to guarantee this algorithm complexity, must use the satellite clock correction of gps satellite track, polar motion of globe parameter, survey station coordinate and troposphere parameter and low sampling rate, and these parameters two poor processing of gps satellite precise orbit determination in step S104 obtained in the non-difference data processing procedure of neutralization procedure S105 medium and low frequency clock correction, according to these information, utilize that the interpolating method of this step comes that interpolation gets to the high frequency clock error correction, can carry out the high frequency clock error correction interpolation of 30s and 5s sampling interval with this algorithm.

In order to introduce this interpolating method, the present invention describes this algorithm as an example of 300s clock correction interpolation 30s clock correction example.Parameter estimation generally is based on least square method and carries out.This interpolating method is divided into following a few step: 1, carry out clock based on the code observed reading and estimate; 2, carry out clock bias estimation based on the phase difference observed reading; 3, select the Reference clock row clock combination of going forward side by side.It should be stressed that at this if want to obtain the clock correction of high sampling rate, original observed data must have identical sample frequency with the high frequency clock correction of trying to achieve at last.In step S105, for fear of the too much parameter of finding the solution, we have carried out the vacuate processing to original observation data, the raw data frequency of so resulting clock correction frequency and vacuate identical (the original observed samples rate of for example participating in calculating is 300s, and the clock correction sampling rate that calculates so also is 300s).So when finding the solution high frequency clock correction, the raw data sample frequency must improve, just can try to achieve the satellite clock correction of 30s or higher sampling interval.

Step S1061 estimates based on the clock that the code observed reading is carried out;

Receiver clock correction must with the gps system time synchronized, this is to be reduced in below the 1mm for the satellite that guarantees to be caused by receiver clock correction and the geometric distance error between the receiver, utilize a code observed reading to carry out when synchronous, for each epoch, estimate respectively receiver clock error correction and satellite clock error correction.

Utilize the observation data (data sampling rate 30s) after step S103 processes to make up GPS code observation equation, be expressed as follows (unit: m, formula of reduction, and do not have error term):

$p_{\mathrm{fk}}^j\left(t_i\right)={\mathrm{\ρ}}_k^j\left(t_i\right)-C\·{\mathrm{\δ}}^j\left(t_i\right)+C\·{\mathrm{\δ}}_k\left(t_i\right)+{\mathrm{\ρ}}_{\mathrm{fk},\mathrm{iono}}^j\left(t_i\right)+{\mathrm{\ρ}}_{\mathrm{fk},\mathrm{trop}}^j\left(t_i\right)+C\·b^j---\left(7\right)$

Wherein,

It is the code pseudorange of observation;

It is oblique distance; δ ^jIt is the satellite clock error correction; δ _kIt is the receiver clock error correction; It is ionospheric refraction;

It is troposphere reflection; b ^jIt is the code deviation (DCB) of difference; C is the light velocity; Subscript k represents receiver, and f represents frequency, and subscript j refers to satellite, t _iRefer to epoch of observation.

The satellite orbit parameter, polar motion of globe parameter, survey station coordinate parameters, troposphere parameter and the DCB that obtain with step S104 are as known number, and observation equation is so:

$p_{\mathrm{fk}}^j\left(t_i\right)=-C\·{\mathrm{\δ}}^j\left(t_i\right)+C\·{\mathrm{\δ}}_k\left(t_i\right)+{\mathrm{\ρ}}_{\mathrm{fk},\mathrm{iono}}^j\left(t_i\right)---\left(8\right)$

If make up

With

Deion layer LC combination

So only remaining satellite and receiver clock error correction:

${\mathrm{Cp}}_k^j\left(t_i\right)=-C\·{\mathrm{\δ}}^j\left(t_i\right)+C\·{\mathrm{\δ}}_k\left(t_i\right)---\left(9\right)$

Based on this class observation equation, to all relevant satellite and receivers, set up respectively the normal equation (Normal Equation, NEQ) of each epoch.Simultaneously, the step that in clock bias estimation, has comprised screening, owing to estimate a parameter with mass data, so receiver and satellite clock error correction algorithm for estimating are very sane, and can come that by iteration underproof data are carried out mark according to the size of residual error and process, and in iteration subsequently, reject.

Step S1062 estimates based on the clock that the phase observations value is carried out;

Generation for the high frequency clock error correction, the phase observations value is main observed reading, contrast with the code inspection process, utilize observation data (code and phase observations value after step S103 processes, data sampling rate 30s) makes up GPS phase observations equation and be expressed as follows (unit: m, formula of reduction, and do not have error term):

${\mathrm{\φ}}_{\mathrm{fk}}^j\left(t_i\right)={\mathrm{\ρ}}_k^j\left(t_i\right)-C\·{\mathrm{\δ}}^j\left(t_i\right)+C\·{\mathrm{\δ}}_k\left(t_i\right)-{\mathrm{\ρ}}_{\mathrm{fk},\mathrm{iono}}^j\left(t_i\right)+{\mathrm{\ρ}}_{\mathrm{fk},\mathrm{trop}}^j\left(t_i\right)+{\mathrm{\λ}}_f\·N_k^j---\left(10\right)$

Wherein,

It is the phase observations value; λ _fIt is the wavelength of frequency f;

It is blur level; The phase observations equation is compared with the code observation equation, and the difference of existence is: the blur level item

With code deviation information (DCB), and the opposite in sign of ionospheric refraction.In order to eliminate the blur level item in the phase observations equation, epoch subsequently carried out the phase observations computing:

${\mathrm{\Δ\φ}}_{\mathrm{fk}}^j(t_{i+1},i)={\mathrm{\φ}}_{\mathrm{fk}}^j\left(t_{i+1}\right)-{\mathrm{\φ}}_{\mathrm{fk}}^j\left(t_i\right)---\left(11\right)$

In addition, use the step same with the code inspection process, the deion layer LC that obtains combination

The phase observations equation as follows:

$\mathrm{\Δ}{\mathrm{C\φ}}_k^j(t_{i+1},i)=-C\·\mathrm{\Δ}{\mathrm{\δ}}^j(t_{i+1},i)+C{\·\mathrm{\Δ\δ}}_k(t_{i+1},i)---\left(12\right)$

The satellite of estimating with single Epoch Algorithm and the difference result of receiver clock error correction value are as the priori value of phase difference process.Therefore, only estimated the corrected value of priori clock correction difference.If occur inconsistent situation between code and the phase observations value, so just reject corresponding phase observations value, be equivalent to utilize step S1061 that the result that step S1062 calculates is retrained; Adjust normal equation with the zero-mean condition of the satellite clock correction corrected value that participates in calculating simultaneously.From the variance-covariance matrix of each of difference, obtain every satellite j's epoch With every receiver k's Take the high frequency clock error correction of satellite and receiver as condition, in computing subsequently, can use these values.

Step S1063 carries out the clock combination, obtains the high frequency clock bias estimation of Navsat and receiver;

In order to obtain the corrected value of clock correction difference, among the code observed reading that obtains according to step S1061 or the step S105 from high precision clock correction time series (such as 5min clock correction solution) the satellite that obtains and the clock bias estimation value of receiver, be respectively every satellite and receiver and carry out the clock correction combination; In order to obtain the clock error correction value of " absolute " phase place clock correction difference, be fixed the interpolation of clock error correction with these clock correction difference informations.

With among the step S1062 as the phase place clock correction difference of observed reading

I=1, L, n-1(n are numbers epoch) set up a new observation equation.With known clock error correction value

As pseudo-observed reading, it is brought on the clock correction collection of 5min solution, form new observation equation group, as follows:

Wherein, δ (t _i) be the phase place clock correction of i observation,

With

Be clock correction Difference Terms that subtracts last after the phase place clock correction,

Be fixing clock error correction, 1≤i≤n.Ignore the correlativity between Difference Solution epoch, and dimension is d=(n-1)+n _Fix* (n-1)+n _FixThe power battle array be diagonal matrix (n _FixBe the clock correction solution interval quantity of lower epoch).From upper one epoch Difference Solution posteriority RMS value obtain pseudo-observed reading

Power.

$\mathrm{P\Δ\δ}\left(t_{i+1,i}\right)=P_{i,i+1}=\frac{{\mathrm{\σ}}_0^2}{{\mathrm{\σ}}_{\mathrm{\Δ\δ}\left(t_{i+1,i}\right)}^2},i=1,...,n---\left(14\right)$

Wherein, σ ₀It is the pseudo range difference of i observation

Error in the weight unit of correspondence, P is the power battle array, for fixing known clock error correction value, according to weighing battle array information to pseudo-observed reading Retrain.

At last, based on formula (13) and constraint condition (14), try to achieve the efficient solution of system of equations, the satellite of the high reliability that namely will try to achieve, high-precision 30s or higher sampling interval and receiver clock correction.

Table 1 has provided the gps satellite clock bias estimation result who utilized three days 30s sampling interval of 40 IGS reference stations calculating in 2012 the 197th day to the 199th day, and with result of calculation and the comparison of CODE center 30s clock correction, as can be seen from Table 1, the poor mean value of satellite clock mistake of 31 gps satellites estimation in three days is about 0.55ns, RMS is within 0.3ns, test findings shows that satellite clock correction and the CODE center 30s clock correction of using this method processing to obtain compare thus, and precision can reach about 0.2ns.

Table 1GPS satellite clock correction estimated result statistical form

As can be seen from the above-described embodiment, the present invention can improve the estimated accuracy of satellite navigation system high frequency clock correction, has reliability height, precision advantages of higher.

Although be the example purpose, the preferred embodiments of the present invention are disclosed, it also is possible those skilled in the art will recognize various improvement, increase and replacement, therefore, scope of the present invention should be not limited to above-described embodiment.

Claims (6)

1. the high frequency clock bias estimation method of a satellite navigation system is characterized in that, comprising:

S101 obtains the observation data of tracking station;

S102, the check observation the reliability of the data;

S103 to utilize the linear combination of observed quantity through the observation data of certificate authenticity, rejects wild value, determines initial blur level, surveys and repair cycle slip, through obtaining code and phase observations value after the smoothing processing;

S104 obtains satellite orbit parameter, polar motion of globe parameter, troposphere parameter, survey station coordinate parameters;

S105, orbit parameter, polar motion of globe parameter, troposphere parameter, the survey station coordinate parameters of the satellite that the original observed data that step S103 is obtained and step S104 obtain; Be brought in the non-poor pseudorange and phase observations value error equation of having eliminated ionosphere effect, the selection reference clock is treated and is estimated that clock imposes restriction, and carries out unified adjustment and processes, and obtains the clock bias estimation value of satellite and receiver;

S106 utilizes code and phase observations value among the step S103, and the clock bias estimation that step S105 is obtained carries out interpolation, obtains the high frequency clock bias estimation of satellite navigation system.

2. the high frequency clock bias estimation method of satellite navigation system as claimed in claim 1 is characterized in that, has eliminated non-poor phase place and the Pseudo-range Observations of ionosphere effect, and the observed reading error equation is as follows:

$v_{k,\mathrm{\Φ}}^j\left(i\right)={\mathrm{\Δt}}_k\left(i\right)-{\mathrm{\Δt}}^j\left(i\right)+{\mathrm{\ρ}}_k^j\left(i\right)/C+\mathrm{\δ}{\mathrm{\ρ}}_{k,\mathrm{trop}}^j\left(i\right)/C+\mathrm{\λ}\·N_k^j/C+{\mathrm{\ϵ}}_{k,\mathrm{\Φ}}^j\left(i\right)-\mathrm{\λ}\·{\mathrm{\Φ}}_k^j\left(i\right)/C---\left(4\right)$

$v_{k,p}^j\left(i\right)={\mathrm{\Δt}}_k\left(i\right)-{\mathrm{\Δt}}^j\left(i\right)+{\mathrm{\ρ}}_k^j\left(i\right)/C+\mathrm{\δ}{\mathrm{\ρ}}_{k,\mathrm{trop}}^j\left(i\right)/C+{\mathrm{\ϵ}}_{k,p}^j\left(i\right)-P_k^i\left(i\right)/C---\left(5\right)$

In formula (4), (5): k is for surveying station number, and j is for defending asterisk, and i is corresponding epoch of observation, and C is the light velocity in the vacuum; Δ t _k(i) be receiver clock correction, Δ t ^j(i) be satellite clock correction,

Be the tropospheric delay impact;

Be not modeled error effect;

Be the blur level parameter;

Be pseudorange and the phase combination observed reading without ionosphere effect;

Be its observational error, λ is the wavelength without ionosphere combination observation value;

For the satellite position of signal x time to the geometric distance between the signal receiver location time of reception.

3. the high frequency clock bias estimation method of satellite navigation system as claimed in claim 2 is characterized in that, the pseudorange of clock correction and phase error equation can be expressed as:

$V\left(i\right)=\left(\begin{array}{c}{V}_p\left(i\right)\\ V_{\mathrm{\Φ}}\left(i\right)\end{array}\right)=\left(\begin{array}{c}{A}_p\\ A_{\mathrm{\Φ}}\end{array}\right)X\left(i\right)+\left(\begin{array}{c}{L}_p\left(i\right)\\ L_{\mathrm{\Φ}}\left(i\right)\end{array}\right)---\left(6\right)$

In the formula (6): $A_p=\left[\begin{array}{cc}\underset{(m\·n)\×(2m+n)}{A}& \underset{(m\·n)\×(m\·n)}{0}\end{array}\right];$ $A_{\mathrm{\Φ}}=\left[\begin{array}{cc}\underset{(m\·n)\×(2m+n)}{A}& \underset{(m\·n)\×(m\·n)}{\mathrm{\λ}/C\·I}\end{array}\right];$ If I is n rank unit matrixs, A is the design matrix when determining receiver and satellite clock correction, be (mn) (m+n) matrix, V (i) is the observational error of i observation position, L (i) is the free term of i observation position, Φ and p represent respectively phase place and pseudorange observation, and X (i) is i parameter to be asked, V (i), L _p(i), L _Φ(i) be known quantity; M base station quantity, n is the quantity of synchronous tracking satellite.

4. the high frequency clock bias estimation method of satellite navigation system as claimed in claim 3 is characterized in that, utilizes a receiver clock in the tracking network as reference clock, and the clock correction precision of reference clock is better than 10-6s; The relative clock correction of all satellites is zero-mean in the same navigational system.

5. the high frequency clock bias estimation method of satellite navigation system as claimed in claim 4 is characterized in that, reference clock is hydrogen atomic clock.

6. the high frequency clock bias estimation method of satellite navigation system as claimed in claim 1 is characterized in that step S106 specifically comprises:

S1061 estimates based on the clock that the code observed reading is carried out;

S1062 estimates based on the clock that the phase observations value is carried out;

S1063 carries out the clock combination, obtains the high frequency clock bias estimation of Navsat and receiver.