CN109655849A - A kind of PPP quickly positions convergent method - Google Patents

Abstract

The present invention provides a kind of PPP quickly to position convergent method, using EKF the disposal of gentle filter, reduces the randomized jitter in convergence process；The location data of preceding 30 epoch is positioned using PPP, classified with mean cluster, acquire the cluster centre of each cluster, cluster centre is averaging, the convergent Fixed Initial Point of PPP is obtained, then convergence calculating is carried out to PPP positioning result with Hatch algorithm, the PPP of proposition of the invention quickly positions convergent method, compared to the previous way by pseudorange One-Point Location or Differential positioning, meets separate unit receiver i.e. and can reach the requirement of cm grades of positioning accuracy；PPP proposed by the present invention quickly positions convergent method, compared to the way that conventional fine One-Point Location resolves, saves for 90% time, has greatly reinforced PPP method real-time and practicability in the calibration of outfield, has been quickly obtained cm grades of positioning accuracies.

Description

A kind of PPP quickly positions convergent method

Technical field

The invention belongs to field of satellite navigation, are related to a kind of localization method.

Background technique

Satellite navigation has many advantages, such as round-the-clock, wide covering and low cost, has become most important in navigation system lead Boat means.In satellite navigation and positioning, mainly there are two kinds of positioning methods of pseudorange One-Point Location and Differential positioning.Pseudorange One-Point Location As long as positioning can be realized in a receiver, there is good real-time, single positioning accuracy is generally tens meters, can be well Meet rough position location requirement, but required precision can not meet the actual demands such as mapping, detection, calibration.GPS accurate one-point is fixed Position technology (PPP) is the completely new station-keeping mode of one kind developed in recent years, using the IGS precise ephemeris provided and accurate clock Difference, user receive ephemeris and observation data using separate unit GPS dual-frequency receiver, and subsequent resolve may be implemented at thousands of square kilometres Or even any position high accuracy positioning in global range.But slowly due to PPP positioning convergence, generally require 15 minutes~30 points Clock can just converge to cm class precision, greatly reduce the efficiency in mapping, detection, calibration, seriously constrain PPP technology in work Application in journey practice.

Summary of the invention

For overcome the deficiencies in the prior art, the present invention provides a kind of method of PPP fast convergence, is smoothly filtered using EKF Wave processing, reduces the randomized jitter in convergence process；The location data that preceding 30 epoch is positioned using PPP, is classified with mean cluster, The cluster centre for acquiring each cluster is averaging cluster centre, obtains the convergent Fixed Initial Point of PPP.Hatch algorithm is used again, Convergence calculating is carried out to PPP positioning result, PPP static immobilization convergence time is substantially reduced, cm can be converged in 2 minutes Grade positioning accuracy, saved for 90% time, has greatly reinforced PPP method real-time and practicability in the calibration of outfield, has obtained cm Grade positioning accuracy.

The technical solution adopted by the present invention to solve the technical problems the following steps are included:

(1) precise ephemeris and precise clock correction for obtaining IGS website orientation in real time obtain broadcast ephemeris and observation data

(2) it using GPS satellite ephemeris and satellite clock correction and dual frequency carrier phase observation, is carried out using non-differential mode type High-precision point positioning,

P(L_i)=ρ+c (dt-dT)+d_orb+d_trop+d_ion/Li+d_mpath/P(Li)+ε(P(L_i))

Wherein, Φ (L_i) it is carrier phase observation data on wave band Li, i=1,2；P(L_i) it is L_iOn Pseudo-range Observations；ρ It is geometric distance of the satellite to receiver,(X_s,Y_s,Z_s) it is satellite launch The coordinate of moment t, (x, y, z) are the receiver coordinates of signal time of reception T=t+ ρ/c；C is the light velocity；Dt is satellite clock correction；dT It is receiver clock-offsets；d_orbIt is satellite orbital error；d_tropIt is tropospheric delay；d_ion/LiFor the ionosphere delay on Li；λ_iIt is Li Wavelength；N_iIt is the integer ambiguity of Li；It is the initial phase of receiver oscillator；It is satellite oscillation The initial phase of device；d_mpath/P(Li)It is the multipath effect of the pseudo-range measurements on Li；d_mpath/φ(Li)It is the carrier phase on Li The multipath effect (m) of measured value；ε () is measurement noise；

(3) precise ephemeris and precise clock correction for using IGS, by tropospheric delay d_tropBe expressed as tropospheric delay zpd with The product of mapping function M；It is combined using double frequency without ionosphere observation to eliminate the single order in ionosphere and influence, the mathematics of One-Point Location Model

(4) fuzziness and coordinate accumulation are carried out using Extended Kalman filter (EKF)；

Construct system equation x_k=(X_k,Y_k,Z_k,zpd_k,N_k), X_k,Y_k,Z_kIt represents k reception seat in the plane to set, zpd_kWhen for k The tropospheric delay at quarter, N_kFor the fuzziness at k moment；

Construction measurement equation y_k=(P_IF,k,Φ_IF,k)=H (x_k)x_k, P_IF,kThe k moment is represented without ionosphere Pseudo-range Observations, Φ_IF,kThe k moment is represented without ionosphere carrier phase observable, H (x_k) it is the coefficient matrix for measuring equation；

Construct system estimation equation Representing the k-1 moment is System estimated value,Represent k moment system estimation value, K_kFor gain matrix, yk is measurement observation；

Construct system estimation equation covarianceI is unit matrix, P_k(-) is k- 1 momentCovariance, P_k(+) is the k momentCovariance；

Construct gain equation

(5) PPP location data is set as X={ x₁,x₂,…,x_n, for each sample point x_i, dot density z_i=1/min ({d_i,j), wherein 1≤i≤n, d_i,jIndicate sample x_iWith sample x_jEuclidean distance；

Correct the regulatory factor of the distance of FCMWherein α_jIt is sample x_jClassification home identity, work as sample x_j∈S_j, i.e. x_jBelong to subset S_j, then α_j=1；Otherwise α_j=0；Whole point where regulatory factor reflects sample in classification Density information.

Modifying factor is obtained, then the distance metric of sample and cluster centre is modified toThen objective function

Construction Lagrange function obtains cluster centre to objective function minimizingAnd degree of membership

Validity Index based on distance correctionWherein,

In PPP positioning convergence process, data are continuously increased, if new points are x_new, x_newTo each class center Distance D_i=d (x_new,v_i), threshold valueWork as D_iWhen > maxdist, x_newIndividually it is classified as one Class, using the data point as new a kind of cluster centre；Work as D_iWhen < maxdist, x_newIt is included into the i-th class；Work as D_i< Maxdist and D_jWhen < maxdist, then the i-th class and jth class are merged, by x_newIt is included into and merges in class；It repeats above operation, directly Until meeting the variance condition of convergence；

(6) PPP is positioned using Hatch algorithm and carries out fast convergence

X (i) is certain location data of PPP positioning,For the average value of cluster centre；

By constantly restraining, when | δ (i) |-| δ (i-1) |≤10^-6, obtain PPP positioning result.

The beneficial effects of the present invention are:

(1) PPP proposed by the present invention quickly positions convergent method, compared in the past by pseudorange One-Point Location or difference The way of positioning meets separate unit receiver i.e. and can reach the requirement of cm grades of positioning accuracy；

(2) PPP proposed by the present invention quickly positions convergent method, does compared to what conventional fine One-Point Location resolved Method saved for 90% time, has greatly reinforced PPP method real-time and practicability in the calibration of outfield, has been quickly obtained cm grades and determines Position precision.

Detailed description of the invention

Fig. 1 is that Static Precise Point Positioning resolved data obtains schematic diagram；

Fig. 2 is the convergence graph of conventional PPP positioning X-direction；

Fig. 3 is the convergence graph of conventional PPP positioning Y-direction；

Fig. 4 is the convergence graph of conventional PPP positioning Z-direction；

Fig. 5 is the convergence graph of quick PPP positioning X-direction；

Fig. 6 is the convergence graph of quick PPP positioning Y-direction；

Fig. 7 is the convergence graph of quick PPP positioning Z-direction.

Specific embodiment

Present invention will be further explained below with reference to the attached drawings and examples, and the present invention includes but are not limited to following implementations Example.

The present invention obtains the Precise Orbit and precise clock correction of IGS website orientation in real time online, calculates PPP using mean cluster Fixed Initial Point is restrained, then convergence process is carried out with Hatch algorithm again, accelerates the convergence time of PPP, be quickly obtained cm grades of positioning Precision.Particular content is as follows:

(1) the online Precise Orbit and precise clock correction for obtaining IGS website orientation in real time

With reference to Fig. 1, the real-time One-Point Location of GPS resolves the real-time accurate ephemeris and precise clock correction for needing to obtain the website IGS, Increase communication module before receiver positioning calculation module, transmitted by mobile data, obtains the website IGS precise ephemeris in real time And clock deviation, broadcast ephemeris and observation data are obtained by receiver antenna, then by real-time accurate ephemeris, precise clock correction, broadcast Ephemeris and observation data are sent into receiver positioning calculation module and are resolved.

(2) double frequency Pseudo-range Observations and carrier phase observation data are constructed

Using high-precision GPS satellite ephemeris and satellite clock correction and dual frequency carrier phase observation, using non-differential mode type High-precision point positioning is carried out, equation is as follows:

P(L_i)=ρ+c (dt-dT)+d_orb+d_trop+d_ion/Li+d_mpath/P(Li)+ε(P(L_i)) (2)

Φ(L_i) it is carrier phase observation data (m) on Li, i=1,2；P(L_i) it is L_iOn Pseudo-range Observations (m)；ρ is Geometric distance of the satellite to receiver, (X_s,Y_s,Z_s) be satellite launch moment t coordinate, (x, y, z) is signal time of reception T The receiver coordinate of=t+ ρ/c,C is the light velocity (m/s)；Dt is satellite clock Poor (s)；DT is receiver clock-offsets (s)；d_orbIt is satellite orbital error (m)；d_tropIt is tropospheric delay (m)；d_io_n/LiFor on Li Ionosphere delay (m)；λ_iIt is the wavelength (m) of Li；N_iIt is the integer ambiguity (week) of Li；It is receiver oscillator Initial phase；It is the initial phase of satellite oscillator；d_mpath/P(Li)It is the multipath of the pseudo-range measurements on Li Effect (m)；d_mpath/φ(Li)It is the multipath effect (m) of the carrier-phase measurement on Li；ε () is measurement noise.

(3) double frequency is calculated without the non-subtractive combination observation in ionosphere

Using the precise satellite ephemeris and Clock Bias of IGS, by tropospheric delay d_tropIt is expressed as tropospheric delay The product of zpd and mapping function M.It is combined using double frequency without ionosphere observation to eliminate the single order in ionosphere and influence, One-Point Location Mathematical model will become:

(4) it is based on the positioning of Extended Kalman filter (EKF)

Fuzziness and coordinate accumulation, fast convergence are carried out using Extended Kalman filter (EKF).Equationof structure is as follows:

System equation:

x_k=(X_k,Y_k,Z_k,zpd_k,N_k) (5)

X_k,Y_k,Z_kIt represents k reception seat in the plane to set, zpd_kFor the tropospheric delay at k moment, N_kFor the fuzziness at k moment.

Measure equation:

y_k=(P_IF,k,Φ_IF,k)=H (x_k)x_k (6)

P_IF,kThe k moment is represented without ionosphere Pseudo-range Observations, Φ_IF,kThe k moment is represented without ionosphere carrier phase observable, H (x_k) For the coefficient matrix for measuring equation.

System estimation equation:

K-1 moment system estimation value is represented,Represent k moment system estimation value, K_kFor gain matrix, y_kFor Measure observation.

System estimation equation covariance:

I is unit matrix, P_k(-) is the k-1 momentCovariance, P_k(+) is the k momentCovariance.

Gain equation:

(5) initial convergence point is calculated using means clustering algorithm (FCM)

If PPP location data is X={ x₁,x₂,…,x_n, for each sample point x_i, dot density expression formula are as follows:

z_i=1/min ({ d_i,j})1≤i≤n (10)

Wherein d_i,jIndicate sample x_iWith sample x_jEuclidean distance.

Correct the regulatory factor of the distance of FCM are as follows:

Wherein α_jIt is sample x_jClassification home identity, as sample x_j∈S_j, i.e. x_jBelong to subset S_j, then α_j=1；Otherwise α_j=0.Whole dot density information where regulatory factor reflects sample in classification.

Modifying factor is obtained, then the distance metric of sample and cluster centre is corrected are as follows:

Then objective function are as follows:

Cluster centre v then can be obtained to objective function minimizing in construction Lagrange function_iWith degree of membership u_i,j:

Validity Index V based on distance correction_KwonAre as follows:

Wherein,

According to the experience of Cluster Classification, if m=2,Iterated conditional ε=10^-6.It is subordinate to random number initialization Spend matrix (u_i.j), c=2:c is calculated according to formula (10)~(16)_maxCorresponding V_hw(U, V, c) value, works as V_hw(U, V, c) takes minimum When value, preferable clustering number c is obtained.Then subordinated-degree matrix is calculated by formula (16)~(21) again and update cluster centre, until two Secondary iteration meets | (u_i,j)_t+1-(u_i,j)_t|≤ε, then completion of classifying.

In PPP positioning convergence process, data are continuously increased, if new points are x_new, x_newTo each class center Distance are as follows: D_i=d (x_new,v_i), threshold value1) work as D_iWhen > maxdist, then x_newIt is single It solely is classified as one kind, using the data point as new a kind of cluster centre；2)D_iWhen < maxdist, then x_newIt is included into the i-th class In；3) work as D_i< maxdist and D_jWhen < maxdist, then the i-th class and jth class are merged, by x_newIt is included into and merges in class.It repeats The above operation, until meeting the variance condition of convergence.

(5) PPP is positioned using Hatch algorithm and carries out fast convergence

X (i) is certain location data of PPP positioning,For the average value of cluster centre.

By constantly restraining, when | δ (i) |-| δ (i-1) |≤10^-6, available PPP positioning result.

The technical proposal for solving the technical problem of the invention is by the online precision for obtaining IGS website orientation in real time Track and precise clock correction calculate PPP using mean cluster and restrain Fixed Initial Point, then carry out convergence process with Hatch algorithm again, add The speed convergence time of PPP, is quickly obtained cm grades of positioning accuracies.Specific step is as follows:

(1) the online Precise Orbit and precise clock correction for obtaining IGS website orientation in real time

With reference to Fig. 1, real-time One-Point Location resolving needs to obtain the website IGS real-time accurate ephemeris and clock deviation, fixed in receiver Position increases communication module before resolving module, is transmitted by mobile data, obtains the ultrafast ephemeris in the website IGS and clock deviation in real time, Broadcast ephemeris and observation data by receiver antenna obtain, then by these observed quantities feeding receiver positioning calculation module into Row resolves.

(2) dual-frequency code and carrier phase observable are constructed

When using double-frequency GPS receiver, code and carrier phase observable between receiver and satellite on L1 and L2 can be indicated Are as follows:

P(L_i)=ρ+c (dt-dT)+d_orb+d_trop+d_ion/Li+d_mpath/P(Li)+ε(P(L_i)) (10)

(3) double frequency is calculated without the non-subtractive combination observation in ionosphere

Ionosphere delay error is one of the main error of un-differenced phase precise point positioning precision, generation apart from shadow It rings, when satellite elevation angle is lower, in some instances it may even be possible to reach rice up to a hundred.It is combined using double frequency without ionosphere un-differenced observation to eliminate electricity The single order of absciss layer influences, and observational equation becomes:

Using the precise satellite ephemeris and Clock Bias of IGS, by tropospheric delay d_tropIt is expressed as tropospheric delay The product of zpd and its mapping function M.The mathematical model of One-Point Location will become:

(4) fuzziness and coordinate are calculated

Fuzziness and coordinate accumulation, fast convergence are carried out using Extended Kalman filter (EKF).Equationof structure is as follows:

x_k=(x_k,y_k,z_k,zpd_k,N_k) (16)

y_k=(P_IF,k,Φ_IF,k)=H (x_k)x_k (17)

(5) initial convergence point is calculated using means clustering algorithm (FCM)

If PPP location data is X={ x₁,x₂,…,x_n, for each sample point x_i, dot density expression formula are as follows:

z_i=1/min ({ d_i,j})1≤i≤n (21)

Wherein d_i,jIndicate sample x_iWith sample x_jEuclidean distance.

Correct the regulatory factor of the distance of FCM are as follows:

Wherein α_jIt is sample x_jClassification home identity, as sample x_j∈S_j, i.e. x_jBelong to subset S_j, then α_j=1；Otherwise α_j=0.Whole dot density information where regulatory factor reflects sample in classification.

Modifying factor is obtained, then the distance metric of sample and cluster centre is corrected are as follows:

Then objective function are as follows:

Cluster centre v then can be obtained to objective function minimizing in construction Lagrange function_iWith degree of membership u_i,j:

Validity Index V based on distance correction_KwonAre as follows:

Wherein,

According to the experience of Cluster Classification, if m=2,Iterated conditional ε=10^-6.It is subordinate to random number initialization Spend matrix (u_i.j), c=2:c is calculated according to formula (21)~(27)_maxCorresponding V_hw(U, V, c) value, works as V_hw(U, V, c) takes minimum When value, preferable clustering number c is obtained.Then subordinated-degree matrix is calculated by formula (16)~(21) again and update cluster centre, until two Secondary iteration meets | (u_i,j)_t+1-(u_i,j)_t|≤ε, then completion of classifying.

In PPP positioning convergence process, data are continuously increased, if new points are x_new, x_newTo each class center Distance are as follows: D_i=d (x_new,v_i), threshold value1) work as D_iWhen > maxdist, then x_newIt is single It solely is classified as one kind, using the data point as new a kind of cluster centre；2)D_iWhen < maxdist, then x_newIt is included into the i-th class In；3) work as D_i< maxdist and D_jWhen < maxdist, then the i-th class and jth class are merged, by x_newIt is included into and merges in class.It repeats The above operation, until meeting the variance condition of convergence.

(5) PPP is positioned using Hatch algorithm and carries out fast convergence

δ²(i)=((n-2) δ²(i-1)+(X(i)-X(i))²)/(i-1) (19)

X (i) is certain location data of PPP positioning,For the average value of cluster centre.

By constantly restraining, when | δ (i) |-| δ (i-1) |≤10^-6, available PPP positioning result.

Below by taking the calibrating and positioning result of somewhere as an example, the specific embodiment of calibrating and positioning method of the invention is carried out Illustrate, and with conventional fine One-Point Location Comparative result.Using method proposed by the present invention, positioned in certain practical calibration Specific step is as follows:

(1) receiver is set up on known point

Receiver antenna is set up in known point R (- 3050576.283,4551789.819,3192861.287), receives GPS Satellite-signal, 2 receivers are observed using common antenna mode, hourly observation data and positioning result.

(2) online real-time accurate track and precise clock correction and positioning calculation

It is transmitted by mobile data, receiver obtains the ultrafast ephemeris in the website IGS and clock deviation in real time, passes through receiver day Line obtains broadcast ephemeris and observation data, these observed quantities feeding receiver positioning calculation module is resolved, is obtained in real time The position of the calibration point.

(3) routine PPP positioning and quick PPP positioning

After continuous observation, point coordinate that conventional PPP finally restrains be P (- 3050576.179, 4551789.876,3192861.307), the coordinate of quick PPP anchor point be Q (- 3050576.218,4551789.857, 3192861.335)。

(4) position error compares

With known point R (- 3050576.283,4551789.819,3192861.287) for true value, more conventional PPP positioning With quick PPP positioning error, subtracting each other can obtain: the error that conventional PPP is located on tri- directions XYZ be respectively as follows: (0.034, 0.057,0.020)；The error that quick PPP is located on tri- directions XYZ is respectively as follows: (0.065,0.038,0.048).By difference Known to value: conventional PPP positioning and quick PPP position error are restrained by the positioning of a period of time at cm grades, have all reached cm's Positioning accuracy.

(5) convergence time compares

With reference to Fig. 2,3,4 as can be seen that routine PPP can converge to cm grades by about 15min, and reference Fig. 5,6, 7 can be seen that quick PPP can converge to cm in 2min, substantially reduce convergence time.

Claims (1)

1. a kind of PPP quickly positions convergent method, it is characterised in that include the following steps:

(1) precise ephemeris and precise clock correction for obtaining IGS website orientation in real time obtain broadcast ephemeris and observation data

(2) it using GPS satellite ephemeris and satellite clock correction and dual frequency carrier phase observation, is carried out using non-differential mode type high-precision One-Point Location is spent,

P(L_i)=ρ+c (dt-dT)+d_orb+d_trop+d_ion/Li+d_mpath/P(Li)+ε(P(L_i))

Wherein, Φ (L_i) it is carrier phase observation data on wave band Li, i=1,2；P(L_i) it is L_iOn Pseudo-range Observations；ρ is to defend Star to receiver geometric distance,(X_s,Y_s,Z_s) it is the satellite launch moment The coordinate of t, (x, y, z) are the receiver coordinates of signal time of reception T=t+ ρ/c；C is the light velocity；Dt is satellite clock correction；DT is to connect Receipts machine clock deviation；d_orbIt is satellite orbital error；d_tropIt is tropospheric delay；d_ion/LiFor the ionosphere delay on Li；λ_iIt is the wave of Li It is long；N_iIt is the integer ambiguity of Li；It is the initial phase of receiver oscillator；It is satellite oscillator Initial phase；d_mpath/P(Li)It is the multipath effect of the pseudo-range measurements on Li；d_mpath/φ(Li)It is the carrier phase measurement on Li The multipath effect (m) of value；ε () is measurement noise；

(3) precise ephemeris and precise clock correction for using IGS, by tropospheric delay d_tropIt is expressed as tropospheric delay zpd and mapping letter The product of number M；It is combined using double frequency without ionosphere observation to eliminate the single order in ionosphere and influence, the mathematical model of One-Point Location

(4) fuzziness and coordinate accumulation are carried out using Extended Kalman filter (EKF)；

Construct system equation x_k=(X_k,Y_k,Z_k,zpd_k,N_k), X_k,Y_k,Z_kIt represents k reception seat in the plane to set, zpd_kFor the k moment Tropospheric delay, N_kFor the fuzziness at k moment；

Construction measurement equation y_k=(P_IF,k,Φ_IF,k)=H (x_k)x_k, P_IF,kThe k moment is represented without ionosphere Pseudo-range Observations, Φ_IF,k The k moment is represented without ionosphere carrier phase observable, H (x_k) it is the coefficient matrix for measuring equation；

Construct system estimation equation Etching system is estimated when representing k-1 Evaluation,Represent k moment system estimation value, K_kFor gain matrix, y_kTo measure observation；

Construct system estimation equation covarianceI is unit matrix, P_k(-) is the k-1 momentCovariance, P_k(+) is the k momentCovariance；

Construct gain equation

(5) PPP location data is set as X={ x₁,x₂,…,x_n, for each sample point x_i, dot density z_i=1/min ({ d_i,j), Wherein, 1≤i≤n, d_i,jIndicate sample x_iWith sample x_jEuclidean distance；

Correct the regulatory factor of the distance of FCMWherein α_jIt is sample x_jClassification home identity, as sample x_j∈ S_j, i.e. x_jBelong to subset S_j, then α_j=1；Otherwise α_j=0；Whole dot density where regulatory factor reflects sample in classification Information；

Modifying factor is obtained, then the distance metric of sample and cluster centre is modified toThen objective function

Construction Lagrange function obtains cluster centre to objective function minimizingAnd degree of membership

Validity Index based on distance correctionWherein,

In PPP positioning convergence process, data are continuously increased, if new points are x_new, x_newTo each class center away from From D_i=d (x_new,v_i), threshold valueWork as D_iWhen > maxdist, x_newIt individually is classified as one kind, Using the data point as new a kind of cluster centre；Work as D_iWhen < maxdist, x_newIt is included into the i-th class；Work as D_i< Maxdist and D_jWhen < maxdist, then the i-th class and jth class are merged, by x_newIt is included into and merges in class；It repeats above operation, directly Until meeting the variance condition of convergence；

(6) PPP is positioned using Hatch algorithm and carries out fast convergence

X (i) is certain location data of PPP positioning,For the average value of cluster centre；

By constantly restraining, when | δ (i) |-| δ (i-1) |≤10^-6, obtain PPP positioning result.