CN117452454A - GNSS receiver timing solution calculation method, device, equipment and storage medium - Google Patents

GNSS receiver timing solution calculation method, device, equipment and storage medium Download PDF

Info

Publication number
CN117452454A
CN117452454A CN202311544584.4A CN202311544584A CN117452454A CN 117452454 A CN117452454 A CN 117452454A CN 202311544584 A CN202311544584 A CN 202311544584A CN 117452454 A CN117452454 A CN 117452454A
Authority
CN
China
Prior art keywords
matrix
abnormal
values
state
timing solution
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.)
Pending
Application number
CN202311544584.4A
Other languages
Chinese (zh)
Inventor
彭竞
付栋
龚航
马明
臧文驰
陈飞强
于美婷
陈曦
刘婷
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
National University of Defense Technology
Original Assignee
National University of Defense Technology
Priority date (The priority date 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 date listed.)
Filing date
Publication date
Application filed by National University of Defense Technology filed Critical National University of Defense Technology
Priority to CN202311544584.4A priority Critical patent/CN117452454A/en
Publication of CN117452454A publication Critical patent/CN117452454A/en
Pending legal-status Critical Current

Links

Classifications

    • GPHYSICS
    • G01MEASURING; TESTING
    • G01SRADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S19/00Satellite radio beacon positioning systems; Determining position, velocity or attitude using signals transmitted by such systems
    • G01S19/01Satellite radio beacon positioning systems transmitting time-stamped messages, e.g. GPS [Global Positioning System], GLONASS [Global Orbiting Navigation Satellite System] or GALILEO
    • G01S19/13Receivers
    • G01S19/35Constructional details or hardware or software details of the signal processing chain
    • G01S19/37Hardware or software details of the signal processing chain
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01SRADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S19/00Satellite radio beacon positioning systems; Determining position, velocity or attitude using signals transmitted by such systems
    • G01S19/38Determining a navigation solution using signals transmitted by a satellite radio beacon positioning system
    • G01S19/39Determining a navigation solution using signals transmitted by a satellite radio beacon positioning system the satellite radio beacon positioning system transmitting time-stamped messages, e.g. GPS [Global Positioning System], GLONASS [Global Orbiting Navigation Satellite System] or GALILEO
    • G01S19/40Correcting position, velocity or attitude

Landscapes

  • Engineering & Computer Science (AREA)
  • Radar, Positioning & Navigation (AREA)
  • Remote Sensing (AREA)
  • Signal Processing (AREA)
  • Computer Networks & Wireless Communication (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Position Fixing By Use Of Radio Waves (AREA)

Abstract

The application relates to a GNSS receiver timing solution calculation method, device, equipment and storage medium. The method comprises the following steps: according to abnormal conditions of observed quantity, a timing solution and abnormal quantity thereof are calculated based on a closed-loop observer, and then the abnormal quantity is deducted from the timing solution, so that a correct result is recovered, the time service reliability is effectively improved, and the time service precision is ensured. The method can be realized by only improving the algorithm of the navigation resolving layer of the receiver, and has good feasibility and popularization value. The problem of calculating reliable timing solutions by using abnormal observables in a complex environment can be solved.

Description

GNSS receiver timing solution calculation method, device, equipment and storage medium
Technical Field
The present disclosure relates to the field of satellite navigation technologies, and in particular, to a method, an apparatus, a device, and a storage medium for calculating a timing solution of a GNSS receiver.
Background
A global satellite navigation system (Global Navigation Satellite System, GNSS) receiver synchronizes the local time of a user to satellite navigation system time by receiving satellite signals, calculating the clock difference and Zhong Piao (collectively referred to as timing solutions) of the local clock relative to the system clock, and further tracing to coordinated universal time, so that the GNSS receiver can provide wide-area high-precision time service for the user. Because the environment where the receiver is located is complex, shielding of satellite signals by tall buildings and trees in cities, multipath caused by refraction, artificially generated in-band radio frequency signals and the like can possibly cause abnormal changes of navigation signal parameters, so that abnormal observables are obtained, and in most cases, users can distinguish whether the observables are credible or not, but cannot recover correct timing solutions under the condition that the observables are not credible. Therefore, how to calculate reliable timing solutions using anomaly observables is a key issue for trusted time service.
However, the current state of the art of GNSS trusted time service is under development. In recent years, students have been publishing related research results, but there are many limitations. The existing method can only judge whether the timing solution is credible or not and can not recover the correct result, and can only autonomously watch time or use other time service means after judging that the timing solution is not credible.
Disclosure of Invention
In view of the foregoing, it is desirable to provide a GNSS receiver timing solution calculation method, apparatus, device, and storage medium capable of recovering a correct timing solution calculation result.
A GNSS receiver timing solution calculation method, the method comprising:
step 1: acquiring state quantity, observed quantity and abnormal vector of an observation satellite; establishing a state equation and an observation equation of the GNSS receiver according to the state quantity, the observed quantity, the abnormal vector, the state transition matrix, the coefficient matrix and the observation matrix; constructing a closed-loop observer model by using a state equation and an observation equation of the GNSS receiver and a plurality of estimated vectors;
step 2: initializing a closed-loop observer model, and calculating the values of all matrixes according to priori knowledge; all the matrixes comprise a state transition matrix, a coefficient matrix, an observation matrix and a constant matrix;
step 3: establishing matrix inequality of a plurality of positive definite matrixes and real matrixes to be solved according to the system stability condition, and substituting values of a state transition matrix and values of an observation matrix into the matrix inequality to obtain values of the positive definite matrixes and the real matrixes; calculating the value of the feedback gain matrix by using the values of the positive definite matrix and the real matrix; substituting the values of the feedback gain matrix into a closed-loop observer model to calculate to obtain the observed quantity of each satellite;
step 4: detecting whether the observed quantity of all satellites is normal according to a preset judging standard, if so, removing abnormal satellites and then transferring to a step 9; otherwise, executing the step 5;
step 5: judging whether all satellites are abnormal for the first time, if so, executing the step 6, and if not, turning to the step 8;
step 6: calculating the deviation between the timing solution calculated by each channel of the previous epoch and the timing solution estimated by the observer;
step 7: performing low-pass filtering on the historical solution value of the abnormal vector to obtain a data set with noise filtered; fitting the frequency drift of the local clock of the receiver according to the data set after noise filtering;
step 8: correcting observed quantities of all the channels before all the abnormal channels are used for calculating the timing solution;
step 9: step 4 is executed to carry out observational quantity judgment according to the channel estimation state quantity after the anomaly correction and the anomaly vector, and if no normal satellite is available, the step 11 is executed; otherwise, executing the step 10;
step 10: taking the timing solution state quantity obtained by performing difference calculation on the state quantity and the abnormal vector in the step 9 as a timing solution result of the current epoch, and turning to the step 4 to perform calculation of the next epoch;
step 11: correcting the state quantity and the observed quantity according to the frequency drift to obtain the corrected state quantity and the corrected observed quantity;
step 12: and (4) taking the corrected state quantity in the step (11) as a timing solution calculation result, and turning to the step (4) to calculate the next epoch.
In one embodiment, constructing a closed-loop observer model using state equations and observation equations of the GNSS receiver and a plurality of estimated vectors includes:
constructing a closed-loop observer model as using state equations and observation equations and multiple estimated vectors of a GNSS receiver
Wherein x is m =x-d represents the corrected timing solution state quantity; x= [ t ] u ,f u ] T Is the state quantity, t u And f u Respectively representing clock differences and Zhong Piao;for the observed quantity ρ n And->Pseudo-range and pseudo-range rate respectively representing nth satelliteN is the number of satellites observed; d= [ d ] 1 ,d 2 ] T Is an abnormal vector +.>Is x m Is used for estimating the vector of the vector; />Is an estimated vector of the observed quantity y; />Is an estimated vector of the anomaly quantity d; />And->Represents a feedback gain matrix, A represents a state transition matrix, B 2 Representing a coefficient matrix, C representing an observation matrix, C l Representing a constant matrix.
In one embodiment, calculating the values of all matrices based on a priori knowledge includes:
calculating the values of all matrices as based on a priori knowledge
Wherein A represents a state transition matrix, B 2 Representing a coefficient matrix, C representing an observation matrix, C l Represents a constant matrix, k represents the current epoch, c represents the speed of light, T s Representing the step size of the discretization,and f s n Respectively representing the calculation result from ephemerisSatellite position, velocity, satellite clock bias, and Zhong Piao; p is p u And v u Representing the position and velocity of the receiver, respectively, as a priori information,to initialize exception vector, ">For the corrected state quantity +.>For observational quantity +.>N is the number of satellites observed.
In one embodiment, a matrix inequality of a plurality of positive definite matrixes and real matrixes to be solved is established according to the system stability condition, and values of a state transition matrix and values of an observation matrix are substituted into the matrix inequality to obtain the values of the positive definite matrixes and the real matrixes; constructing a feedback gain matrix using values of the positive definite matrix and the real matrix, comprising:
establishing matrix inequality of a plurality of positive definite matrixes and real matrixes to be solved as follows according to system stability conditions
{P,Q,M}={P T ,Q T ,M T }>0
Wherein, A represents a state transition matrix, C represents an observation matrix, G, P, Q represents different positive definite matrices, M represents a real matrix, and superscript T represents transposition operation;
constructing feedback gain matrices L and Γ as using the values of the positive and real matrices
In one embodiment, the process of correcting the observed quantity of each channel comprises the following steps:
correcting observed quantity of each channel according to deviation between timing solution calculated by each channel of previous epoch and timing solution estimated by observer to obtain corrected observed quantity as
Wherein,representing the deviation between the timing solution calculated for each channel of the previous epoch and the timing solution estimated by the observer, k representing the current epoch, k 0 A previous epoch, y, representing all channel anomalies n (k) Representing the observed quantity of the current channel.
In one embodiment, the process of correcting the state quantity and the observed quantity according to the frequency drift to obtain the corrected state quantity and the corrected observed quantity includes:
the effect of frequency drift is subtracted by:
wherein,and->The clock difference estimated in the step 9 and the Zhong Piao abnormal quantity, f drift Representing frequency drift, T s Representing a discretization step size;
then calculates the accumulated anomaly vector according to the following formulaThe form is as follows:
in one embodiment, the slaveAnd->Less->Obtaining corrected state quantity->And observed quantity->
Wherein,representing the estimated state quantity ∈ ->Representing the estimated observables.
A GNSS receiver timing solution computing device, the device comprising:
a closed-loop observer model module is constructed and used for acquiring the state quantity, the observed quantity and the abnormal vector of the observation satellite; establishing a state equation and an observation equation of the GNSS receiver according to the state quantity, the observed quantity, the abnormal vector, the state transition matrix, the coefficient matrix and the observation matrix; constructing a closed-loop observer model by using a state equation and an observation equation of the GNSS receiver and a plurality of estimated vectors;
the observed quantity calculation module is used for initializing a closed-loop observer model and calculating the values of all matrixes according to priori knowledge; all the matrixes comprise a state transition matrix, a coefficient matrix, an observation matrix and a constant matrix; establishing matrix inequality of a plurality of positive definite matrixes and real matrixes to be solved according to the system stability condition, and substituting values of a state transition matrix and values of an observation matrix into the matrix inequality to obtain values of the positive definite matrixes and the real matrixes; calculating the value of the feedback gain matrix by using the values of the positive definite matrix and the real matrix; substituting the values of the feedback gain matrix into a closed-loop observer model to calculate to obtain the observed quantity of each satellite;
a timing solution calculation module, configured to step 4: detecting whether the observed quantity of all satellites is normal according to a preset judging standard, if so, removing abnormal satellites and then transferring to a step 9; otherwise, executing the step 5; step 5: judging whether all satellites are abnormal for the first time, if so, executing the step 6, and if not, turning to the step 8; step 6: calculating the deviation between the timing solution calculated by each channel of the previous epoch and the timing solution estimated by the observer; step 7: performing low-pass filtering on the historical solution value of the abnormal vector to obtain a data set with noise filtered; fitting the frequency drift of the local clock of the receiver according to the data set after noise filtering; step 8: correcting observed quantities of all the channels before all the abnormal channels are used for calculating the timing solution; step 9: step 4 is executed to carry out observational quantity judgment according to the channel estimation state quantity after the anomaly correction and the anomaly vector, and if no normal satellite is available, the step 11 is executed; otherwise, executing the step 10; step 10: taking the timing solution state quantity obtained by performing difference calculation on the state quantity and the abnormal vector in the step 9 as a timing solution result of the current epoch, and turning to the step 4 to perform calculation of the next epoch; step 11: correcting the state quantity and the observed quantity according to the frequency drift to obtain the corrected state quantity and the corrected observed quantity; step 12: and (4) taking the corrected state quantity in the step (11) as a timing solution calculation result, and turning to the step (4) to calculate the next epoch.
A computer device comprising a memory storing a computer program and a processor which when executing the computer program performs the steps of:
step 1: acquiring state quantity, observed quantity and abnormal vector of an observation satellite; establishing a state equation and an observation equation of the GNSS receiver according to the state quantity, the observed quantity, the abnormal vector, the state transition matrix, the coefficient matrix and the observation matrix; constructing a closed-loop observer model by using a state equation and an observation equation of the GNSS receiver and a plurality of estimated vectors;
step 2: initializing a closed-loop observer model, and calculating the values of all matrixes according to priori knowledge; all the matrixes comprise a state transition matrix, a coefficient matrix, an observation matrix and a constant matrix;
step 3: establishing matrix inequality of a plurality of positive definite matrixes and real matrixes to be solved according to the system stability condition, and substituting values of a state transition matrix and values of an observation matrix into the matrix inequality to obtain values of the positive definite matrixes and the real matrixes; calculating the value of the feedback gain matrix by using the values of the positive definite matrix and the real matrix; substituting the values of the feedback gain matrix into a closed-loop observer model to calculate to obtain the observed quantity of each satellite;
step 4: detecting whether the observed quantity of all satellites is normal according to a preset judging standard, if so, removing abnormal satellites and then transferring to a step 9; otherwise, executing the step 5;
step 5: judging whether all satellites are abnormal for the first time, if so, executing the step 6, and if not, turning to the step 8;
step 6: calculating the deviation between the timing solution calculated by each channel of the previous epoch and the timing solution estimated by the observer;
step 7: performing low-pass filtering on the historical solution value of the abnormal vector to obtain a data set with noise filtered; fitting the frequency drift of the local clock of the receiver according to the data set after noise filtering;
step 8: correcting observed quantities of all the channels before all the abnormal channels are used for calculating the timing solution;
step 9: step 4 is executed to carry out observational quantity judgment according to the channel estimation state quantity after the anomaly correction and the anomaly vector, and if no normal satellite is available, the step 11 is executed; otherwise, executing the step 10;
step 10: taking the timing solution state quantity obtained by performing difference calculation on the state quantity and the abnormal vector in the step 9 as a timing solution result of the current epoch, and turning to the step 4 to perform calculation of the next epoch;
step 11: correcting the state quantity and the observed quantity according to the frequency drift to obtain the corrected state quantity and the corrected observed quantity;
step 12: and (4) taking the corrected state quantity in the step (11) as a timing solution calculation result, and turning to the step (4) to calculate the next epoch.
A computer readable storage medium having stored thereon a computer program which when executed by a processor performs the steps of:
step 1: acquiring state quantity, observed quantity and abnormal vector of an observation satellite; establishing a state equation and an observation equation of the GNSS receiver according to the state quantity, the observed quantity, the abnormal vector, the state transition matrix, the coefficient matrix and the observation matrix; constructing a closed-loop observer model by using a state equation and an observation equation of the GNSS receiver and a plurality of estimated vectors;
step 2: initializing a closed-loop observer model, and calculating the values of all matrixes according to priori knowledge; all the matrixes comprise a state transition matrix, a coefficient matrix, an observation matrix and a constant matrix;
step 3: establishing matrix inequality of a plurality of positive definite matrixes and real matrixes to be solved according to the system stability condition, and substituting values of a state transition matrix and values of an observation matrix into the matrix inequality to obtain values of the positive definite matrixes and the real matrixes; calculating the value of the feedback gain matrix by using the values of the positive definite matrix and the real matrix; substituting the values of the feedback gain matrix into a closed-loop observer model to calculate to obtain the observed quantity of each satellite;
step 4: detecting whether the observed quantity of all satellites is normal according to a preset judging standard, if so, removing abnormal satellites and then transferring to a step 9; otherwise, executing the step 5;
step 5: judging whether all satellites are abnormal for the first time, if so, executing the step 6, and if not, turning to the step 8;
step 6: calculating the deviation between the timing solution calculated by each channel of the previous epoch and the timing solution estimated by the observer;
step 7: performing low-pass filtering on the historical solution value of the abnormal vector to obtain a data set with noise filtered; fitting the frequency drift of the local clock of the receiver according to the data set after noise filtering;
step 8: correcting observed quantities of all the channels before all the abnormal channels are used for calculating the timing solution;
step 9: step 4 is executed to carry out observational quantity judgment according to the channel estimation state quantity after the anomaly correction and the anomaly vector, and if no normal satellite is available, the step 11 is executed; otherwise, executing the step 10;
step 10: taking the timing solution state quantity obtained by performing difference calculation on the state quantity and the abnormal vector in the step 9 as a timing solution result of the current epoch, and turning to the step 4 to perform calculation of the next epoch;
step 11: correcting the state quantity and the observed quantity according to the frequency drift to obtain the corrected state quantity and the corrected observed quantity;
step 12: and (4) taking the corrected state quantity in the step (11) as a timing solution calculation result, and turning to the step (4) to calculate the next epoch.
According to the GNSS receiver timing solution calculation method, device, equipment and storage medium, firstly, the timing solution and the abnormal quantity thereof are calculated based on the closed loop observer according to the abnormal condition of the observed quantity, then the abnormal quantity is deducted from the timing solution, so that a correct result is recovered, the timing solution abnormality caused by the observed quantity disturbance can be corrected in real time, the frequency drift influence caused by the clock with poor precision can be restrained, no extra hardware is needed, the integration is only developed on the navigation solution layer, the precision and the credibility of GNSS timing are ensured, the timing solution with reliable timing calculation by using the abnormal observed quantity is ensured, the problem that the timing solution is not credible in the complex environment is solved, and the method can be realized by only improving the algorithm of the receiver navigation solution layer, has good feasibility and popularization value and is a feasible economic and efficient method.
Drawings
FIG. 1 is a flow chart illustrating a method for calculating a timing solution of a GNSS receiver according to an embodiment;
FIG. 2 is a block diagram illustrating an apparatus for calculating a timing solution of a GNSS receiver according to an embodiment;
FIG. 3 is an internal block diagram of a computer device in one embodiment.
Detailed Description
In order to make the objects, technical solutions and advantages of the present application more apparent, the present application will be further described in detail with reference to the accompanying drawings and examples. It should be understood that the specific embodiments described herein are for purposes of illustration only and are not intended to limit the present application.
In one embodiment, as shown in fig. 1, a GNSS receiver timing solution calculation method is provided, including the following steps:
step 1: acquiring state quantity, observed quantity and abnormal vector of an observation satellite; establishing a state equation and an observation equation of the GNSS receiver according to the state quantity, the observed quantity, the abnormal vector, the state transition matrix, the coefficient matrix and the observation matrix; constructing a closed-loop observer model by using a state equation and an observation equation of the GNSS receiver and a plurality of estimated vectors;
because the environment where the receiver is located is complex, the satellite signals are blocked, multipath caused by refraction, artificially generated in-band radio frequency signals and the like are likely to cause disturbance to observables, the invention aims to simultaneously estimate the timing solution and abnormal quantity thereof, so that the correct timing solution is recovered when the observables are not reliable, and the calculation process of the timing solution of the GNSS receiver can be modeled as the following system state equation and observation equation according to the basic theory of a control system:
x(k+1)=Ax(k)+B 2 d(k)
y(k)=Cx(k)+c l (k)+ξ(k)
wherein k represents an epoch; c represents the speed of light; x= [ t ] u ,f u ] T Is a state vector, t u And f u Respectively representing clock differences and Zhong Piao;for observing the vector ρ n And->The pseudo range and the pseudo range rate of the nth satellite are respectively represented, and N is the number of the observed satellites; d= [ d ] 1 ,d 2 ] T Is an anomaly vector, d 1 And d 2 Abnormal amounts respectively representing clock skew and Zhong Piao; a represents a state transition matrix; b (B) 2 Representing a coefficient matrix; c represents an observation matrix; c l Representing a constant matrix; ζ represents a noise vector.
Step 2: initializing a closed-loop observer model, and calculating the values of all matrixes according to priori knowledge; all the matrixes comprise a state transition matrix, a coefficient matrix, an observation matrix and a constant matrix;
by correcting the closed-loop observer corresponding to the system, the observer capable of simultaneously estimating the timing solution and the abnormal state quantity thereof can be obtained, and the values of each matrix in the system model are calculated according to the set and known parameters:
wherein the discretization step length T s Setting to 1 second; position p may be determined by a stationary receiver or other means u And velocity v u Then c l (k) It can be calculated entirely from ephemeris.
At the same time, anomaly vectorNormally initialized to zero vector, modified state vector +.>And an observation vectorThe values can be assigned by the timing solution and the observed quantity in the last period, and the observer can finally converge to the correct result as long as the initial value is in a reasonable range.
Step 3: establishing matrix inequality of a plurality of positive definite matrixes and real matrixes to be solved according to the system stability condition, and substituting values of a state transition matrix and values of an observation matrix into the matrix inequality to obtain values of the positive definite matrixes and the real matrixes; calculating the value of the feedback gain matrix by using the values of the positive definite matrix and the real matrix; substituting the values of the feedback gain matrix into a closed-loop observer model to calculate to obtain the observed quantity of each satellite;
substituting the matrix value obtained by initializing the system model into a matrix inequality to solve the values of a positive definite matrix and a real matrix, and calculating the values of the positive definite matrix and the real matrix to obtain the value of a feedback gain matrix; substituting the values of the feedback gain matrix into a closed-loop observer model to calculate the observed quantity of each satellite.
Step 4: detecting whether the observed quantity of all satellites is normal according to a preset judging standard, if so, removing abnormal satellites and then transferring to a step 9; otherwise, executing the step 5;
and judging the observed quantity state. Detecting whether the observed quantity of all satellites is normal, if so, excluding the observed quantity of abnormal satellites, and switching to S9 only by using a normal satellite list; otherwise, executing the next step.
Step 5: judging whether all satellites are abnormal for the first time, if so, executing the step 6, and if not, turning to the step 8;
and judging an abnormal state. After all satellites are abnormal, only abnormal observables can be used, but abnormal satellites cannot be directly eliminated, so that whether all satellites are abnormal for the first time or not needs to be judged, and the next step is executed; otherwise, go to S8.
Step 6: calculating the deviation between the timing solution calculated by each channel of the previous epoch and the timing solution estimated by the observer;
when all satellite channels are identified as anomalous, the deviation between the timing solution calculated for each channel and the correct timing solution estimated by the closed-loop observer needs to be calculated:
wherein k is 0 A previous epoch that represents all channel anomalies;a timing solution calculated from the pseudorange/pseudorange rate observation equation for channel n under known receiver position conditions.
The deviation is used to calculate the accumulated error caused in the later epoch.
Step 7: performing low-pass filtering on the historical solution value of the abnormal vector to obtain a data set with noise filtered; fitting the frequency drift of the local clock of the receiver according to the data set after noise filtering;
because of the different precision of the local clocks used by GNSS receivers, the worse clocks have significant frequency drift, which affects the anomaly vector estimated by the closed loop observer, and therefore the value of the frequency drift needs to be calculated by fitting. Due to abnormal vectorsThe history of (2) will have a lot of high frequency noise, which needs to be doneLow-pass filtering to obtain the filtered abnormal quantity d f To fit the frequency drift:
f drift =(1 T 1) -1 1 T d f
wherein 1 is d f A co-dimensional unit vector.
Step 8: correcting observed quantities of all the channels before all the abnormal channels are used for calculating the timing solution;
before all outlier channels are used to calculate the timing solution, the observed deviations accumulated from the previous epoch need to be corrected:
step 9: step 4 is executed to carry out observational quantity judgment according to the channel estimation state quantity after the anomaly correction and the anomaly vector, and if no normal satellite is available, the step 11 is executed; otherwise, executing the step 10;
executing the estimation process of the closed-loop observer, and selecting the next operation according to the abnormal state of the S4 satellite: when no normal satellite is available, correcting abnormal quantity based on estimated timing solution, namely turning to S11; otherwise, the calculation result of the closed-loop observer can be directly output at S10. The accuracy of the timing solution is improved by correcting timing solution anomalies caused by observed quantity disturbances in real time.
Step 10: taking the timing solution state quantity obtained by performing difference calculation on the state quantity and the abnormal vector in the step 9 as a timing solution result of the current epoch, and turning to the step 4 to perform calculation of the next epoch;
the abnormal-state-free quantity obtained by calculation of S9As a result of the timing solution of the present epoch, and go to S4 to perform calculation of the next epoch.
Step 11: correcting the state quantity and the observed quantity according to the frequency drift to obtain the corrected state quantity and the corrected observed quantity;
by correcting the state quantity and the observed quantity by using the frequency drift, the influence of the frequency drift caused by the clock with poor precision can be restrained, and the calculation precision of the state quantity and the observed quantity is improved.
Step 12: and (4) taking the corrected state quantity in the step (11) as a timing solution calculation result, and turning to the step (4) to calculate the next epoch.
According to the GNSS receiver timing solution calculation method, firstly, the timing solution and the abnormal quantity thereof are calculated based on the closed loop observer according to the abnormal condition of the observed quantity, then the abnormal quantity is deducted from the timing solution, so that a correct result is recovered, the timing solution abnormality caused by the observed quantity disturbance can be corrected in real time, the frequency drift influence caused by a clock with poor precision can be restrained, no extra hardware is needed, the integration is only developed on a navigation solution layer, the precision and the credibility of GNSS timing are guaranteed, the timing safety of a great infrastructure is guaranteed, the reliable timing solution is calculated by using the abnormal observed quantity, the problem that the timing solution is not credible in a complex environment is solved, the method can be realized by only improving the algorithm of the receiver navigation solution layer, and the method has good feasibility and popularization value and is a feasible economical and efficient method.
In one embodiment, constructing a closed-loop observer model using state equations and observation equations of the GNSS receiver and a plurality of estimated vectors includes:
constructing a closed-loop observer model as using state equations and observation equations and multiple estimated vectors of a GNSS receiver
Wherein x is m =x-d represents the corrected timing solution state quantity; x= [ t ] u ,f u ] T Is the state quantity, t u And f u Respectively representing clock differences and Zhong Piao;for the observed quantity ρ n And->The pseudo range and the pseudo range rate of the nth satellite are respectively represented, and N is the number of the observed satellites; d= [ d ] 1 ,d 2 ] T Is an abnormal vector +.>Is x m Is used for estimating the vector of the vector; />Is an estimated vector of the observed quantity y; />Is an estimated vector of the anomaly quantity d; />And->Represents a feedback gain matrix, A represents a state transition matrix, B 2 Representing a coefficient matrix, C representing an observation matrix, C l Representing a constant matrix.
In one embodiment, calculating the values of all matrices based on a priori knowledge includes:
calculating the values of all matrices as based on a priori knowledge
Wherein A represents a state transition matrix, B 2 Representing a coefficient matrix, C representing an observation matrix, C l Represents a constant matrix, k represents the current epoch, c represents the speed of light, T s Representing the step size of the discretization,and->Respectively representing satellite position, velocity, satellite clock difference and Zhong Piao calculated according to ephemeris; p is p u And v u Representing the position and velocity of the receiver, respectively, as a priori information,to initialize exception vector, ">For the corrected state quantity +.>For observational quantity +.>N is the number of satellites observed.
In one embodiment, a matrix inequality of a plurality of positive definite matrixes and real matrixes to be solved is established according to the system stability condition, and values of a state transition matrix and values of an observation matrix are substituted into the matrix inequality to obtain the values of the positive definite matrixes and the real matrixes; constructing a feedback gain matrix using values of the positive definite matrix and the real matrix, comprising:
establishing matrix inequality of a plurality of positive definite matrixes and real matrixes to be solved as follows according to system stability conditions
{P,Q,M}={P T ,Q T ,M T }>0
Wherein, A represents a state transition matrix, C represents an observation matrix, G, P, Q represents different positive definite matrices, M represents a real matrix, and superscript T represents transposition operation;
constructing feedback gain matrices L and Γ as using the values of the positive and real matrices
/>
In one embodiment, the process of correcting the observed quantity of each channel comprises the following steps:
correcting observed quantity of each channel according to deviation between timing solution calculated by each channel of previous epoch and timing solution estimated by observer to obtain corrected observed quantity as
Wherein,representing the deviation between the timing solution calculated for each channel of the previous epoch and the timing solution estimated by the observer, k representing the current epoch, k 0 A previous epoch, y, representing all channel anomalies n (k) Representing the observed quantity of the current channel.
In one embodiment, the process of correcting the state quantity and the observed quantity according to the frequency drift to obtain the corrected state quantity and the corrected observed quantity includes:
the effect of frequency drift is subtracted by:
wherein,and->The clock difference estimated in the step 9 and the Zhong Piao abnormal quantity, f drift Representing frequency drift, T s Representing a discretization step size;
then calculates the accumulated anomaly vector according to the following formulaThe form is as follows:
in one embodiment, the slaveAnd->Less->Obtaining corrected state quantity->And observed quantity->
Wherein,representing the estimated state quantity ∈ ->Representing the estimated observables.
It should be understood that, although the steps in the flowchart of fig. 1 are shown in sequence as indicated by the arrows, the steps are not necessarily performed in sequence as indicated by the arrows. The steps are not strictly limited to the order of execution unless explicitly recited herein, and the steps may be executed in other orders. Moreover, at least some of the steps in fig. 1 may include multiple sub-steps or stages that are not necessarily performed at the same time, but may be performed at different times, nor do the order in which the sub-steps or stages are performed necessarily performed in sequence, but may be performed alternately or alternately with at least a portion of other steps or sub-steps of other steps.
In one embodiment, as shown in fig. 2, there is provided a GNSS receiver timing solution calculating apparatus, including: a closed-loop observer model module 202, an observed quantity calculation module 204, and a timing solution calculation module are constructed, wherein:
a closed-loop observer model module is constructed and used for acquiring the state quantity, the observed quantity and the abnormal vector of the observation satellite; establishing a state equation and an observation equation of the GNSS receiver according to the state quantity, the observed quantity, the abnormal vector, the state transition matrix, the coefficient matrix and the observation matrix; constructing a closed-loop observer model by using a state equation and an observation equation of the GNSS receiver and a plurality of estimated vectors;
the observed quantity calculation module is used for initializing a closed-loop observer model and calculating the values of all matrixes according to priori knowledge; all the matrixes comprise a state transition matrix, a coefficient matrix, an observation matrix and a constant matrix; establishing matrix inequality of a plurality of positive definite matrixes and real matrixes to be solved according to the system stability condition, and substituting values of a state transition matrix and values of an observation matrix into the matrix inequality to obtain values of the positive definite matrixes and the real matrixes; calculating the value of the feedback gain matrix by using the values of the positive definite matrix and the real matrix; substituting the values of the feedback gain matrix into a closed-loop observer model to calculate to obtain the observed quantity of each satellite;
a timing solution calculation module, configured to step 4: detecting whether the observed quantity of all satellites is normal according to a preset judging standard, if so, removing abnormal satellites and then transferring to a step 9; otherwise, executing the step 5; step 5: judging whether all satellites are abnormal for the first time, if so, executing the step 6, and if not, turning to the step 8; step 6: calculating the deviation between the timing solution calculated by each channel of the previous epoch and the timing solution estimated by the observer; step 7: performing low-pass filtering on the historical solution value of the abnormal vector to obtain a data set with noise filtered; fitting the frequency drift of the local clock of the receiver according to the data set after noise filtering; step 8: correcting observed quantities of all the channels before all the abnormal channels are used for calculating the timing solution; step 9: step 4 is executed to carry out observational quantity judgment according to the channel estimation state quantity after the anomaly correction and the anomaly vector, and if no normal satellite is available, the step 11 is executed; otherwise, executing the step 10; step 10: taking the timing solution state quantity obtained by performing difference calculation on the state quantity and the abnormal vector in the step 9 as a timing solution result of the current epoch, and turning to the step 4 to perform calculation of the next epoch; step 11: correcting the state quantity and the observed quantity according to the frequency drift to obtain the corrected state quantity and the corrected observed quantity; step 12: and (4) taking the corrected state quantity in the step (11) as a timing solution calculation result, and turning to the step (4) to calculate the next epoch.
For specific limitations on the GNSS receiver timing solution calculation means, reference may be made to the above limitations on the GNSS receiver timing solution calculation method, and no further description is given here. The above-described modules in the GNSS receiver timing solution computing device may be implemented in whole or in part by software, hardware, or a combination thereof. The above modules may be embedded in hardware or may be independent of a processor in the computer device, or may be stored in software in a memory in the computer device, so that the processor may call and execute operations corresponding to the above modules.
In one embodiment, a computer device is provided, which may be a terminal, and the internal structure of which may be as shown in fig. 3. The computer device includes a processor, a memory, a network interface, a display screen, and an input device connected by a system bus. Wherein the processor of the computer device is configured to provide computing and control capabilities. The memory of the computer device includes a non-volatile storage medium and an internal memory. The non-volatile storage medium stores an operating system and a computer program. The internal memory provides an environment for the operation of the operating system and computer programs in the non-volatile storage media. The network interface of the computer device is used for communicating with an external terminal through a network connection. The computer program is executed by a processor to implement a GNSS receiver timing solution calculation method. The display screen of the computer equipment can be a liquid crystal display screen or an electronic ink display screen, and the input device of the computer equipment can be a touch layer covered on the display screen, can also be keys, a track ball or a touch pad arranged on the shell of the computer equipment, and can also be an external keyboard, a touch pad or a mouse and the like.
It will be appreciated by those skilled in the art that the structure shown in fig. 3 is merely a block diagram of some of the structures associated with the present application and is not limiting of the computer device to which the present application may be applied, and that a particular computer device may include more or fewer components than shown, or may combine certain components, or have a different arrangement of components.
Those skilled in the art will appreciate that implementing all or part of the above described methods may be accomplished by way of a computer program stored on a non-transitory computer readable storage medium, which when executed, may comprise the steps of the embodiments of the methods described above. Any reference to memory, storage, database, or other medium used in the various embodiments provided herein may include non-volatile and/or volatile memory. The nonvolatile memory can include Read Only Memory (ROM), programmable ROM (PROM), electrically Programmable ROM (EPROM), electrically Erasable Programmable ROM (EEPROM), or flash memory. Volatile memory can include Random Access Memory (RAM) or external cache memory. By way of illustration and not limitation, RAM is available in a variety of forms such as Static RAM (SRAM), dynamic RAM (DRAM), synchronous DRAM (SDRAM), double Data Rate SDRAM (DDRSDRAM), enhanced SDRAM (ESDRAM), synchronous Link DRAM (SLDRAM), memory bus direct RAM (RDRAM), direct memory bus dynamic RAM (DRDRAM), and memory bus dynamic RAM (RDRAM), among others.
The technical features of the above embodiments may be arbitrarily combined, and all possible combinations of the technical features in the above embodiments are not described for brevity of description, however, as long as there is no contradiction between the combinations of the technical features, they should be considered as the scope of the description.
The above examples merely represent a few embodiments of the present application, which are described in more detail and are not to be construed as limiting the scope of the invention. It should be noted that it would be apparent to those skilled in the art that various modifications and improvements could be made without departing from the spirit of the present application, which would be within the scope of the present application. Accordingly, the scope of protection of the present application is to be determined by the claims appended hereto.

Claims (10)

1. A method for calculating a timing solution for a GNSS receiver, the method comprising:
step 1: acquiring state quantity, observed quantity and abnormal vector of an observation satellite; establishing a state equation and an observation equation of the GNSS receiver according to the state quantity, the observed quantity, the abnormal vector, the state transition matrix, the coefficient matrix and the observation matrix; constructing a closed-loop observer model by using a state equation and an observation equation of the GNSS receiver and a plurality of estimated vectors;
step 2: initializing the closed-loop observer model, and calculating the values of all matrixes according to priori knowledge; all the matrixes comprise a state transition matrix, a coefficient matrix, an observation matrix and a constant matrix;
step 3: establishing matrix inequality of a plurality of positive definite matrixes and real matrixes to be solved according to the system stability condition, and substituting values of a state transition matrix and values of an observation matrix into the matrix inequality to obtain values of the positive definite matrixes and the real matrixes; calculating the value of the feedback gain matrix by using the values of the positive definite matrix and the real matrix; substituting the values of the feedback gain matrix into the closed-loop observer model to calculate the observed quantity of each satellite;
step 4: detecting whether the observed quantity of all satellites is normal according to a preset judging standard, if so, removing abnormal satellites and then transferring to a step 9; otherwise, executing the step 5;
step 5: judging whether all satellites are abnormal for the first time, if so, executing the step 6, and if not, turning to the step 8;
step 6: calculating the deviation between the timing solution calculated by each channel of the previous epoch and the timing solution estimated by the observer;
step 7: performing low-pass filtering on the historical solution value of the abnormal vector to obtain a data set with noise filtered; fitting the frequency drift of the local clock of the receiver according to the data set after noise filtering;
step 8: correcting observed quantities of all the channels before all the abnormal channels are used for calculating the timing solution;
step 9: step 4 is executed to carry out observational quantity judgment according to the channel estimation state quantity after the anomaly correction and the anomaly vector, and if no normal satellite is available, the step 11 is executed; otherwise, executing the step 10;
step 10: taking the timing solution state quantity obtained by performing difference calculation on the state quantity and the abnormal vector in the step 9 as a timing solution result of the current epoch, and turning to the step 4 to perform calculation of the next epoch;
step 11: correcting the state quantity and the observed quantity according to the frequency drift to obtain the corrected state quantity and the corrected observed quantity;
step 12: and (4) taking the corrected state quantity in the step (11) as a timing solution calculation result, and turning to the step (4) to calculate the next epoch.
2. The method of claim 1, wherein constructing a closed-loop observer model using state equations and observation equations and a plurality of estimated vectors of the GNSS receiver comprises:
constructing a closed-loop observer model by using a state equation and an observation equation of the GNSS receiver and a plurality of estimated vectors as follows
Wherein x is m =x-d represents the corrected timing solution state quantity; x= [ t ] u ,f u ] T Is the state quantity, t u And f u Respectively representing clock differences and Zhong Piao;for the observed quantity ρ n And->The pseudo range and the pseudo range rate of the nth satellite are respectively represented, and N is the number of the observed satellites; d= [ d ] 1 ,d 2 ] T Is an abnormal vector +.>Is x m Is used for estimating the vector of the vector; />Is an estimated vector of the observed quantity y; />Is an estimated vector of the anomaly quantity d; />And->Represents a feedback gain matrix, A represents a state transition matrix, B 2 Representing a coefficient matrix, C representing an observation matrix, C l Representing a constant matrix.
3. The method of claim 1, wherein calculating values of all matrices based on a priori knowledge comprises:
calculating the values of all matrices as based on a priori knowledge
Wherein A represents a state transition matrix, B 2 Representing a coefficient matrix, C representing an observation matrix, C l Represents a constant matrix, k represents the current epoch, c represents the speed of light, T s Representing the step size of the discretization,and f s n Respectively representing satellite position, velocity, satellite clock difference and Zhong Piao calculated according to ephemeris; p is p u And v u Representing the position and velocity of the receiver, respectively, as a priori information,to initialize exception vector, ">For the corrected state quantity +.>For observational quantity +.>N is the number of satellites observed.
4. The method according to claim 1, wherein matrix inequalities of a plurality of positive definite matrices and real matrices to be solved are established according to system stability conditions, and values of a state transition matrix and values of an observation matrix are substituted into the matrix inequalities to obtain values of the positive definite matrices and the real matrices; constructing a feedback gain matrix using the values of the positive definite matrix and the real matrix, including:
establishing matrix inequality of a plurality of positive definite matrixes and real matrixes to be solved as follows according to system stability conditions
{P,Q,M}={P T ,Q T ,M T }>0
Wherein, A represents a state transition matrix, C represents an observation matrix, G, P, Q represents different positive definite matrices, M represents a real matrix, and superscript T represents transposition operation;
constructing feedback gain matrices L and Γ as using the values of the positive and real matrices
5. The method of claim 1, wherein correcting the observed quantity of each channel comprises:
correcting observed quantity of each channel according to deviation between timing solution calculated by each channel of previous epoch and timing solution estimated by observer to obtain corrected observed quantity as
Wherein,representing the deviation between the timing solution calculated for each channel of the previous epoch and the timing solution estimated by the observer, k representing the current epoch, k 0 A previous epoch, y, representing all channel anomalies n (k) Representing the observed quantity of the current channel.
6. The method of claim 1, wherein correcting the state quantity and the observed quantity based on the frequency drift, the process of obtaining the corrected state quantity and the observed quantity comprises:
the effect of frequency drift is subtracted by:
wherein,and->The clock difference estimated in the step 9 and the Zhong Piao abnormal quantity, f drift Representing frequency drift, T s Representing a discretization step size;
then calculates the accumulated anomaly vector according to the following formulaThe form is as follows:
7. the method of claim 6, wherein the method further comprises:
from the slaveAnd->Less->Obtaining corrected state quantity->And observed quantity->
Wherein,representing the estimated state quantity ∈ ->Representing the estimated observables.
8. A GNSS receiver timing solution computing device, the device comprising:
a closed-loop observer model module is constructed and used for acquiring the state quantity, the observed quantity and the abnormal vector of the observation satellite; establishing a state equation and an observation equation of the GNSS receiver according to the state quantity, the observed quantity, the abnormal vector, the state transition matrix, the coefficient matrix and the observation matrix; constructing a closed-loop observer model by using a state equation and an observation equation of the GNSS receiver and a plurality of estimated vectors;
the observed quantity calculation module is used for initializing the closed-loop observer model and calculating the values of all matrixes according to priori knowledge; all the matrixes comprise a state transition matrix, a coefficient matrix, an observation matrix and a constant matrix; establishing matrix inequality of a plurality of positive definite matrixes and real matrixes to be solved according to the system stability condition, and substituting values of a state transition matrix and values of an observation matrix into the matrix inequality to obtain values of the positive definite matrixes and the real matrixes; calculating the value of the feedback gain matrix by using the values of the positive definite matrix and the real matrix; substituting the values of the feedback gain matrix into the closed-loop observer model to calculate the observed quantity of each satellite;
a timing solution calculation module, configured to step 4: detecting whether the observed quantity of all satellites is normal according to a preset judging standard, if so, removing abnormal satellites and then transferring to a step 9; otherwise, executing the step 5; step 5: judging whether all satellites are abnormal for the first time, if so, executing the step 6, and if not, turning to the step 8; step 6: calculating the deviation between the timing solution calculated by each channel of the previous epoch and the timing solution estimated by the observer; step 7: performing low-pass filtering on the historical solution value of the abnormal vector to obtain a data set with noise filtered; fitting the frequency drift of the local clock of the receiver according to the data set after noise filtering; step 8: correcting observed quantities of all the channels before all the abnormal channels are used for calculating the timing solution; step 9: step 4 is executed to carry out observational quantity judgment according to the channel estimation state quantity after the anomaly correction and the anomaly vector, and if no normal satellite is available, the step 11 is executed; otherwise, executing the step 10; step 10: taking the timing solution state quantity obtained by performing difference calculation on the state quantity and the abnormal vector in the step 9 as a timing solution result of the current epoch, and turning to the step 4 to perform calculation of the next epoch; step 11: correcting the state quantity and the observed quantity according to the frequency drift to obtain the corrected state quantity and the corrected observed quantity; step 12: and (4) taking the corrected state quantity in the step (11) as a timing solution calculation result, and turning to the step (4) to calculate the next epoch.
9. A computer device comprising a memory and a processor, the memory storing a computer program, characterized in that the processor implements the steps of the method of any of claims 1 to 7 when the computer program is executed.
10. A computer readable storage medium, on which a computer program is stored, characterized in that the computer program, when being executed by a processor, implements the steps of the method of any of claims 1 to 7.
CN202311544584.4A 2023-11-17 2023-11-17 GNSS receiver timing solution calculation method, device, equipment and storage medium Pending CN117452454A (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202311544584.4A CN117452454A (en) 2023-11-17 2023-11-17 GNSS receiver timing solution calculation method, device, equipment and storage medium

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202311544584.4A CN117452454A (en) 2023-11-17 2023-11-17 GNSS receiver timing solution calculation method, device, equipment and storage medium

Publications (1)

Publication Number Publication Date
CN117452454A true CN117452454A (en) 2024-01-26

Family

ID=89587420

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202311544584.4A Pending CN117452454A (en) 2023-11-17 2023-11-17 GNSS receiver timing solution calculation method, device, equipment and storage medium

Country Status (1)

Country Link
CN (1) CN117452454A (en)

Similar Documents

Publication Publication Date Title
US9513129B2 (en) Low authority GPS aiding of navigation system for anti-spoofing
Li et al. PPP/INS tightly coupled navigation using adaptive federated filter
US9939532B2 (en) Heading for a hybrid navigation solution based on magnetically calibrated measurements
US20100141515A1 (en) Position tracking device and method
US8543266B2 (en) Modified Kalman filter for generation of attitude error corrections
JP2013511038A (en) Detection and correction of anomalous measurements and determination of ambiguity in a global navigation satellite system receiver.
CN112711039B (en) Time synchronization attack detection and correction method and device based on optimal estimation
CN106093991A (en) A kind of fuzziness quick recovery method for GNSS location and system
CN110376620A (en) Real-time clock deviation forecasting procedure, device and computer equipment
US8949027B2 (en) Multiple truth reference system and method
JP2018531371A6 (en) Integer resolution of multi-time carrier phase GNSS
JP2018531371A (en) Integer resolution of multi-time carrier phase GNSS
CN114647178B (en) Automatic atomic clock calibration method and system based on Beidou and ground reference transmission
Gao et al. Double-channel sequential probability ratio test for failure detection in multisensor integrated systems
Dai et al. Real-time precise orbit determination for BDS satellites using the square root information filter
CN109212560B (en) Method and apparatus for compensating frequency inaccuracy of frequency generator such as oscillator
CN116399351A (en) Vehicle position estimation method
KR20180092791A (en) Method and apparatus for improving position-velocity solution in gnss receivers
CN117452454A (en) GNSS receiver timing solution calculation method, device, equipment and storage medium
CN117055323A (en) Star-based precise time service method and system based on Beidou/Galileo system fusion
CN103727947A (en) BDS (Beidou Navigation System) and GIS (Geographic Information System) deep coupling location method and system based on UKF (Unscented Kalman Filter)
CN115060275A (en) Navigation information optimization method for multiple inertial navigation devices
Chen et al. A factor set-based GNSS fault detection and exclusion for vehicle navigation in urban environments
CN114001730A (en) Fusion positioning method and device, computer equipment and storage medium
Yun et al. Reducing the computation time in the state chi-square test for IMU fault detection

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