WO2023197487A1 - 一种应用于全球导航卫星系统gnss的周跳探测与修复方法 - Google Patents

一种应用于全球导航卫星系统gnss的周跳探测与修复方法 Download PDF

Info

Publication number
WO2023197487A1
WO2023197487A1 PCT/CN2022/112718 CN2022112718W WO2023197487A1 WO 2023197487 A1 WO2023197487 A1 WO 2023197487A1 CN 2022112718 W CN2022112718 W CN 2022112718W WO 2023197487 A1 WO2023197487 A1 WO 2023197487A1
Authority
WO
WIPO (PCT)
Prior art keywords
epoch
difference
cycle slip
common
carrier phase
Prior art date
Application number
PCT/CN2022/112718
Other languages
English (en)
French (fr)
Inventor
刘欢
章浙涛
袁海军
张明凯
Original Assignee
上海海积信息科技股份有限公司
河海大学
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 上海海积信息科技股份有限公司, 河海大学 filed Critical 上海海积信息科技股份有限公司
Publication of WO2023197487A1 publication Critical patent/WO2023197487A1/zh

Links

Images

Classifications

    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO 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/00—Satellite radio beacon positioning systems; Determining position, velocity or attitude using signals transmitted by such systems
    • G01S19/01—Satellite radio beacon positioning systems transmitting time-stamped messages, e.g. GPS [Global Positioning System], GLONASS [Global Orbiting Navigation Satellite System] or GALILEO
    • G01S19/13—Receivers
    • G01S19/35—Constructional details or hardware or software details of the signal processing chain
    • G01S19/37—Hardware or software details of the signal processing chain
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO 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/00—Satellite radio beacon positioning systems; Determining position, velocity or attitude using signals transmitted by such systems
    • G01S19/38—Determining a navigation solution using signals transmitted by a satellite radio beacon positioning system
    • G01S19/39—Determining 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/42—Determining position
    • G01S19/43—Determining position using carrier phase measurements, e.g. kinematic positioning; using long or short baseline interferometry
    • G01S19/44—Carrier phase ambiguity resolution; Floating ambiguity; LAMBDA [Least-squares AMBiguity Decorrelation Adjustment] method
    • Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
    • Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02A—TECHNOLOGIES FOR ADAPTATION TO CLIMATE CHANGE
    • Y02A90/00—Technologies having an indirect contribution to adaptation to climate change
    • Y02A90/10—Information and communication technologies [ICT] supporting adaptation to climate change, e.g. for weather forecasting or climate simulation

Definitions

  • Embodiments of the present invention relate to the field of computer technology, and in particular, to a cycle slip detection and repair method, device, computing device and computer-readable storage medium applied to the Global Navigation Satellite System (GNSS).
  • GNSS Global Navigation Satellite System
  • GNSS Global Navigation Satellite System
  • the satellite signals will inevitably be affected by the observation environment, especially in complex observation environments.
  • the satellite signal is temporarily interrupted by obstacles, the lock is lost due to radio signal interference, and the GNSS reception Conditions such as loss of lock of the airborne wave phase tracking loop and low signal-to-noise ratio will lead to the cycle slip phenomenon.
  • Cycle slip refers to the jump or interruption of the entire cycle count due to the loss of lock of the satellite signal in the carrier phase measurement. Therefore, when the signal is re-tracked, the entire cycle count is incorrect. Therefore, correctly detecting and recovering cycle slips is one of the very important and must-solved problems in carrier phase measurement, and is also the key to achieving high-precision and reliable navigation and positioning. important prerequisite.
  • cycle slip detection methods include high-order difference method, polynomial fitting method, pseudorange phase combination method, etc. These methods are essentially based on geometry-free (GF) combination.
  • GF geometry-free
  • the Turbo Edit method is the most widely used cycle slip detection method. It essentially includes a combination of GF and MW. This method is also greatly affected by the ionosphere and observation noise.
  • the observation values are often noisy, resulting in the above methods being unable to accurately and effectively detect and repair cycle slips.
  • Embodiments of the present invention provide a cycle slip detection and repair method applied to a global navigation satellite system GNSS to reduce the impact of the ionosphere and observation noise on cycle slip detection and improve the accuracy and reliability of cycle slip detection.
  • embodiments of the present invention provide a cycle slip detection method applied to a global navigation satellite system GNSS, including:
  • the pseudo-range observation values of the common-view satellites determine the s2 pseudo-range three-difference observation values of the monitoring station at the t epoch; where s2 ⁇ 3;
  • the s2 carrier phase three-difference observations and the s2 pseudorange three-difference observations corresponding to the t epoch are substituted into the three-difference observation equation based on the geometric GB model to obtain 2s2 solution formulas; the three-difference observation of the GB model The equation is:
  • the least squares estimation is used to solve the 2s2 solution equations, and the solution of the s2 common-view satellites at epoch t is obtained.
  • Each integer cycle slip value if any integer cycle slip value meets the preset conditions, it is determined that the common view satellite corresponding to the integer cycle slip value has a cycle slip during the observation process.
  • the s2 carrier phase three-difference observations and the s2 pseudorange three-difference observations corresponding to epoch t are substituted into the three-difference observation equation based on the geometric GB model, and 2s2 solution equations are obtained.
  • the GB model takes into account the monitoring station and the common-view satellite. The geometric relationship between them is more tolerant to observation noise. Since the baseline length between the monitoring station and the reference station is less than the first preset threshold and the interval between epochs for obtaining observation values is less than the second preset threshold, the double-difference ionosphere between stations and satellites between epochs can be ignored. The influence of the troposphere and multipath can further simplify the three-difference observation equation of the GB model.
  • the solution before each integer cycle slip value of s2 common-view satellites at epoch t also includes:
  • the pseudo-range observation values of s1 common-view satellites received in 1 epoch determine the s1 pseudo-range three-difference observation values of the monitoring station in the t-1 epoch; the t-1 epoch is earlier than the t Epoch; the s1 common-view satellites do not have a cycle slip phenomenon or the cycle slip has been repaired in the t-1 epoch;
  • the least squares estimation is used to solve the 2s2 solution equations, and the solution of the s2 common-view satellites at epoch t is obtained Various integer cycle slip values, including:
  • the monitoring station is a static monitoring station; the b t is equal to b t-1 ; the b t-1 is the three-difference value of the satellite distance of the t-1 epoch expressed in the form of coordinate parameters. .
  • the number of parameters to be solved is further reduced, the probability that the solution formula is solvable and estimable is increased, and the number of redundant observation values in the solution formula is increased, further improving the accuracy of the cycle slip value obtained and reducing observation noise. Impact.
  • the calculated cycle slip value has higher accuracy and stability.
  • the least squares estimation is used to solve the 2s2 solution equations, and the solution is that the s2 common-view satellites are in Each integer cycle slip value of epoch t includes:
  • l is the matrix formed by the s2 carrier phase observation values and the s2 pseudorange three-difference observation values
  • E(l) is the s2 carrier phase three-difference observation values and the s2 pseudorange three-difference values.
  • B is The matrix formed by the components in the three directions of x, y and z and the wavelengths of the satellite signals of the s2 common-view satellites
  • X is the components of b t in the three directions of x, y and z and the Each of s2 common vision satellites the matrix formed;
  • the X is obtained by the following formula of least squares estimation:
  • the D3 is the variance covariance matrix of the l
  • each integer cycle slip value of the s2 common-view satellites at epoch t is obtained.
  • the above least squares method is used to solve each integer cycle slip value, which is more accurate.
  • each integer cycle slip value of the s2 common-view satellites at epoch t is obtained, including:
  • the variance and covariance matrix is determined for each cycle slip value in floating point form. According to the cycle slip value of each floating point form and the cycle slip value of each floating point form
  • the cycle slip value in the form of an integer can be obtained by using the variance covariance matrix of the value.
  • the cycle slip value in the form of an integer can reflect the occurrence of the cycle slip phenomenon. Only the cycle slip value in the form of an integer can be used to repair the cycle slip.
  • the carrier phase observation values include carrier phase observation values of multiple wavelengths; the carrier phase observation values of each wavelength determine respective corresponding carrier phase three-difference observation values;
  • the pseudorange observation values include pseudorange observation values of multiple wavelengths, and the pseudorange observation values of each wavelength determine their corresponding pseudorange three-difference observation values.
  • embodiments of the present invention also provide a cycle slip detection and repair device applied to the Global Navigation Satellite System GNSS, including:
  • the pseudo-range observation values of the common-view satellites determine the s2 pseudo-range three-difference observation values of the monitoring station at the t epoch; where s2 ⁇ 3;
  • the s2 carrier phase three-difference observations and the s2 pseudorange three-difference observations corresponding to the t epoch are substituted into the three-difference observation equation based on the geometric GB model to obtain 2s2 solution formulas; the three-difference observation of the GB model The equation is:
  • a calculation unit configured to use least squares estimation to solve the 2s2 solution equations based on the s2 carrier phase three-difference observation values and the s2 pseudorange three-difference observation values, and obtain s2 common-view satellites.
  • the determining unit is also used to:
  • the pseudo-range observation values of s1 common-view satellites received in 1 epoch determine the s1 pseudo-range three-difference observation values of the monitoring station in the t-1 epoch; the t-1 epoch is earlier than the t Epoch; the s1 common-view satellites do not have a cycle slip phenomenon or the cycle slip has been repaired in the t-1 epoch;
  • the computing unit is specifically used for:
  • the monitoring station is a static monitoring station; the b t is equal to b t-1 ; the b t-1 is the three-difference value of the satellite distance of the t-1 epoch expressed in the form of coordinate parameters. .
  • the computing unit is specifically used for:
  • l is the matrix formed by the s2 carrier phase observation values and the s2 pseudorange three-difference observation values
  • E(l) is the s2 carrier phase three-difference observation values and the s2 pseudorange three-difference values.
  • B is The matrix formed by the components in the three directions of x, y and z and the wavelengths of the satellite signals of the s2 common-view satellites
  • X is the components of b t in the three directions of x, y and z and the Each of s2 common vision satellites the matrix formed;
  • the X is obtained by the following formula of least squares estimation:
  • the D3 is the variance covariance matrix of the l
  • each integer cycle slip value of the s2 common-view satellites at epoch t is obtained.
  • the computing unit is specifically used for:
  • the computing unit is specifically used for:
  • the carrier phase observation values include carrier phase observation values of multiple wavelengths; the carrier phase observation values of each wavelength determine respective corresponding carrier phase three-difference observation values;
  • the pseudorange observation values include pseudorange observation values of multiple wavelengths, and the pseudorange observation values of each wavelength determine their corresponding pseudorange three-difference observation values.
  • an embodiment of the present invention further provides a computing device, including:
  • Memory used to store computer programs
  • a processor configured to call the computer program stored in the memory, and execute the cycle slip detection and repair method applied to GNSS listed in any of the above methods according to the obtained program.
  • embodiments of the present invention also provide a computer-readable storage medium that stores a computer-executable program, and the computer-executable program is used to cause the computer to execute any of the methods listed above. Cycle slip detection and repair methods applied to GNSS.
  • Figure 1 is a schematic flow chart of a cycle slip detection and repair method applied to a global navigation satellite system GNSS provided by an embodiment of the present invention
  • Figure 2 is a schematic diagram of an application scenario provided by an embodiment of the present invention.
  • Figure 3 is a schematic flow chart of a cycle slip detection and repair method applied to the global navigation satellite system GNSS provided by an embodiment of the present invention
  • Figure 4 is a schematic structural diagram of a cycle slip detection and repair device applied to a global navigation satellite system GNSS provided by an embodiment of the present invention
  • Figure 5 is a schematic structural diagram of a computer device provided by an embodiment of the present invention.
  • Global navigation satellite system refers to a space-based radio navigation and positioning system that can provide users with all-weather 3D coordinates, speed and time information at any location on the earth's surface or near-Earth space.
  • BDS BeiDou navigation satellite system
  • GPS United States' global positioning system
  • Galileo navigation satellite system the European Union's Galileo navigation satellite system
  • GLONASS global orbiting navigation satellite system
  • the Beidou Satellite Navigation System is a global satellite navigation system independently constructed and operated by China. It provides all-weather, all-weather, high-precision positioning, navigation and timing services to users around the world.
  • satellite navigation systems have been widely used in various fields such as aviation, car navigation, communications, surveying and mapping, and entertainment.
  • GNSS satellite signals are prone to cycle slip phenomena. How to reduce the impact of complex observation environments on cycle slip detection, thereby improving the accuracy and reliability of cycle slip detection, has become a key issue in current cycle slip detection. important topic.
  • embodiments of the present invention provide a cycle slip detection and repair method applied to the global navigation satellite system GNSS, as shown in Figure 1, including:
  • Step 101 Determine the s2 carrier phase three-difference observation values of the monitoring station at the t epoch based on the carrier phase observation values of the s2 common-view satellites received by the monitoring station at the t epoch; The received pseudo-range observation values of s2 common-view satellites determine the s2 pseudo-range three-difference observation values of the monitoring station at the t epoch; where s2 ⁇ 3;
  • Figure 2 shows an application scenario applicable to the embodiment of the present invention, including a reference satellite u, a common-view satellite s, a reference station q, and a monitoring station r.
  • a receiver is provided on the monitoring station r.
  • the monitoring station r can receive the carrier phase observation values and pseudo-range observation values of multiple common-view satellites, and send the carrier phase observation values and pseudo-range observation values of multiple common-view satellites to the computing device.
  • the computing device The carrier phase three-difference observation value of each common-view satellite is determined based on the carrier phase observation values of multiple common-view satellites, and the pseudo-range three-difference observation value of each common-view satellite is determined based on the pseudo-range observation values of multiple common-view satellites.
  • the carrier phase observation value of the common-view satellite s received by the monitoring station r in epoch t and the carrier phase observation value of the reference satellite u received by the monitoring station r in epoch t will be used.
  • the carrier phase observation value of the common-view satellite s is received, the carrier phase observation value of the reference satellite u is received by the monitoring station r in the t-1 epoch, and the carrier phase of the common-view satellite s is received by the reference station q in the t-1 epoch.
  • the specific method of determining the three-difference observation value of the carrier phase is a technical method well known to those skilled in the art and will not be described in detail here.
  • the method of determining the pseudorange three-difference observation value is the same.
  • step 103 In order to ensure that each solution formula determined later is solvable and estimable, s2 ⁇ 3 is restricted here. The specific reasons are detailed in step 103.
  • the corresponding three carrier phase triple-difference observations were determined; based on the pseudo-range observations of the three common-view satellites received by the monitoring station r, the corresponding 3 pseudorange three-difference observations.
  • Step 102 Substitute the s2 carrier phase three-difference observation values and the s2 pseudorange three-difference observation values corresponding to the t epoch into the three-difference observation equation based on the geometric GB model to obtain 2s2 solution equations.
  • the GB model takes into account the geometric relationship between the monitoring station and the common-view satellite, has a higher tolerance for observation noise, and can reduce the impact of observation noise on cycle slip detection.
  • the embodiment of the present invention optimizes the traditional GB model.
  • the three-difference observation value of the satellite distance does not take into account the geometric relationship between the monitoring station and the common vision satellite.
  • the GF linear combination method which is well known to those skilled in the art, is often used for calculation. Therefore, the obtained cycle slip value is greatly affected by the observation noise, and high-precision cycle slip detection cannot be performed.
  • a t is the design matrix of epoch t, which is used to express the geometric relationship between the monitoring station r and the common vision satellite s in epoch t
  • a t-1 is the design matrix of epoch t-1, which is used to express the geometric relationship between the monitoring station r and the common vision satellite s of epoch t-1
  • x is the satellite of epoch t expressed in the form of coordinate parameters. Distance double difference.
  • the design matrix is determined in advance based on the coordinates of the common vision satellite and the roughly known coordinates of the monitoring station (not accurate enough), and is a known value.
  • parameters expressed in “bold and italic” represent that the parameters are in matrix form.
  • the embodiment of the present invention limits the baseline length between the monitoring station and the reference station to be less than the first preset threshold, such as less than 10 kilometers, 9 kilometers, etc., so that it can be ignored Inter-epochal station star double difference ionosphere troposphere influence; and limit the sampling interval, that is, the interval between epochs for obtaining observation values is less than the second preset threshold, such as less than 5 seconds, less than 3 seconds, etc. In this way, multipathing can be ignored Impact.
  • the first preset threshold such as less than 10 kilometers, 9 kilometers, etc.
  • ⁇ (Ax) A t x t -A t-1 x t-1 .
  • the carrier phase three-difference observation value and the pseudo-range three-difference observation value can be obtained according to the satellite signal (carrier phase observation value and pseudo-range observation value) received by the monitoring station (i.e., as described in step 101), Therefore, it is a known quantity; the design matrix is a known quantity, and the wavelength of the received satellite signal is also a known quantity, so, x t , x t-1 , cycle slip value is the parameter to be solved, it is necessary to check x t , x t-1 and cycle slip value Solve (observation value noise does not need to be solved).
  • sampling interval is generally small, generally 1s-5s apart, so the magnitude of the element difference between A t and A t-1 is 10 -3 , therefore, A t -A t-1 ⁇ 0.
  • b t is the three-difference value of the satellite distance at epoch t expressed in the form of coordinate parameters.
  • the parameters to be solved are b t and cycle slip value
  • Step 103 Based on the s2 carrier phase three-difference observation values and the s2 pseudo-range three-difference observation values, use least squares estimation to solve the 2s2 solution equations, and obtain the s2 common-view satellites at t Each integer cycle slip value of the epoch; if any integer cycle slip value meets the preset conditions, it is determined that the common view satellite corresponding to the integer cycle slip value has a cycle slip during the observation process.
  • l is the matrix formed by the s2 carrier phase observation values and the s2 pseudorange three-difference observation values
  • E(l) is the s2 carrier phase three-difference observation values and the s2 pseudorange three-difference values.
  • the expectation of difference observed value B is (i.e. Gt) is a matrix formed by the components in the three directions of x, y and z and the wavelengths of the satellite signals of the s2 common-view satellites
  • X is the component of the b t in the three directions of x, y and z components and each of the s2 common-view satellites the matrix formed;
  • X is the parameter to be solved, and the X is obtained through the following formula of least squares estimation:
  • the D3 is the variance covariance matrix of the l
  • the variance covariance matrix D3 of the three-difference observation value l is obtained through the variance-covariance matrix D2 of the double-difference observation value:
  • C is the design matrix for converting double-difference observations into triple-difference observations.
  • the cycle slip value in floating point form is converted into a cycle slip value in integer form through the LAMBDA algorithm, and the cycle slip value of the s2 common-view satellites at epoch t is obtained.
  • Each integer cycle slip value is obtained.
  • Step 301 Determine the s1 carrier phase three-difference observation values of the monitoring station at the t-1 epoch based on the carrier phase observation values of s1 common-view satellites received by the monitoring station at the t-1 epoch; according to the monitoring station
  • the pseudo-range observation values of s1 common-view satellites received at the t-1 epoch determine the s1 pseudo-range three-difference observation values of the monitoring station at the t-1 epoch; the t-1 epoch is earlier than The t epoch; the s1 common-view satellites have no cycle slip phenomenon or the cycle slip has been repaired in the t-1 epoch.
  • the s1 common-view satellites corresponding to epoch t-1 and the s2 common-view satellites corresponding to epoch t may be the same or different.
  • s1 may be equal to s2, or it may be It does not mean that there is no restriction on this.
  • Step 302 Substitute s1 carrier phase three-difference observations and s1 pseudorange three-difference observations corresponding to the t-1 epoch into the three-difference observation equation of the GB model to obtain 2s1 solution equations.
  • b t-1 is the three-difference value of the satellite distance at epoch t-1 expressed in the form of coordinate parameters.
  • the number of solution formulas is 2s2, and the number of parameters to be solved is 3+s2; after the three-difference observation value of epoch t-1 is introduced, the number of solution formulas is 2 (s1+s2), the number of parameters to be solved is 3+3+s2. Therefore, in order to make the solution formula solvable, it is necessary to ensure that 2(s1+s2) ⁇ 6+s2.
  • the three-difference values of satellite distances in different epochs are equal, and the three-difference observation value of t-1 epoch is introduced, but no new parameters to be solved will be introduced.
  • the calculated cycle slip value has higher accuracy and stability.
  • the observation value can be not only a single-frequency observation value, but also a multi-frequency observation value such as dual-frequency, triple-frequency, etc.
  • the carrier phase observation value includes the carrier phase observation value of multiple wavelengths; the carrier phase observation value of each wavelength determines the corresponding carrier phase three-difference observation value; the pseudorange observation value includes the pseudorange observation value of multiple wavelengths , the pseudorange observation value of each wavelength determines the corresponding pseudorange three-difference observation value.
  • e 4 is a 4 ⁇ 1 all-1 vector
  • I s2 is the s 2- dimensional identity matrix
  • ⁇ i and ⁇ j are the wavelengths of frequency i and frequency j respectively.
  • Figure 4 exemplarily shows the structure of a device for cycle slip detection and repair of a global navigation satellite system GNSS provided by an embodiment of the present invention.
  • This structure can be applied to a global navigation satellite system.
  • GNSS cycle slip detection and repair process can be applied to a global navigation satellite system.
  • the device specifically includes:
  • Determining unit 401 used for:
  • the pseudo-range observation values of the common-view satellites determine the s2 pseudo-range three-difference observation values of the monitoring station at the t epoch; where s2 ⁇ 3;
  • the s2 carrier phase three-difference observations and the s2 pseudorange three-difference observations corresponding to the t epoch are substituted into the three-difference observation equation based on the geometric GB model to obtain 2s2 solution formulas; the three-difference observation of the GB model The equation is:
  • the calculation unit 402 is configured to use the least squares estimation to solve the 2s2 solution equations based on the s2 carrier phase three-difference observation values and the s2 pseudorange three-difference observation values, and obtain s2 common views.
  • the determining unit 401 is also used to:
  • the pseudo-range observation values of s1 common-view satellites received in 1 epoch determine the s1 pseudo-range three-difference observation values of the monitoring station in the t-1 epoch; the t-1 epoch is earlier than the t Epoch; the s1 common-view satellites do not have a cycle slip phenomenon or the cycle slip has been repaired in the t-1 epoch;
  • the computing unit 402 is specifically used for:
  • the monitoring station is a static monitoring station; the b t is equal to b t-1 ; the b t-1 is the three-difference value of the satellite distance of the t-1 epoch expressed in the form of coordinate parameters. .
  • the computing unit 402 is specifically used to:
  • l is the matrix formed by the s2 carrier phase observation values and the s2 pseudorange three-difference observation values
  • E(l) is the s2 carrier phase three-difference observation values and the s2 pseudorange three-difference values.
  • B is The matrix formed by the components in the three directions of x, y and z and the wavelengths of the satellite signals of the s2 common-view satellites
  • X is the components of b t in the three directions of x, y and z and the Each of s2 common vision satellites the matrix formed;
  • the X is obtained by the following formula of least squares estimation:
  • the D3 is the variance covariance matrix of the l
  • each integer cycle slip value of the s2 common-view satellites at epoch t is obtained.
  • the computing unit 402 is specifically used to:
  • the computing unit 402 is specifically used to:
  • the carrier phase observation values include carrier phase observation values of multiple wavelengths; the carrier phase observation values of each wavelength determine respective corresponding carrier phase three-difference observation values;
  • the pseudorange observation values include pseudorange observation values of multiple wavelengths, and the pseudorange observation values of each wavelength determine their corresponding pseudorange three-difference observation values.
  • the embodiment of the present application provides a computer device, as shown in Figure 5, including at least one processor 501, and a memory 502 connected to the at least one processor.
  • the processor is not limited in the embodiment of the present application.
  • the connection between processor 501 and memory 502 in Figure 5 is taken as an example through a bus.
  • the bus can be divided into address bus, data bus, control bus, etc.
  • the memory 502 stores instructions that can be executed by at least one processor 501.
  • at least one processor 501 can perform the above-mentioned cycle slip detection and detection applied to the global navigation satellite system GNSS. Steps of the repair method.
  • the processor 501 is the control center of the computer equipment. It can use various interfaces and lines to connect various parts of the computer equipment, and perform applications by running or executing instructions stored in the memory 502 and calling data stored in the memory 502. Cycle slip detection and repair of Global Navigation Satellite System (GNSS).
  • the processor 501 may include one or more processing units.
  • the processor 501 may integrate an application processor and a modem processor.
  • the application processor mainly processes the operating system, user interface, application programs, etc., and the modem processor
  • the debug processor mainly handles wireless communications. It can be understood that the above modem processor may not be integrated into the processor 501.
  • the processor 501 and the memory 502 can be implemented on the same chip, and in some embodiments, they can also be implemented on separate chips.
  • the processor 501 may be a general-purpose processor, such as a central processing unit (CPU), a digital signal processor, an application specific integrated circuit (Application Specific Integrated Circuit, ASIC), a field programmable gate array or other programmable logic devices, discrete gates or transistors Logic devices and discrete hardware components can implement or execute the methods, steps and logical block diagrams disclosed in the embodiments of this application.
  • a general-purpose processor may be a microprocessor or any conventional processor, etc. The steps of the methods disclosed in conjunction with the embodiments of the present application can be directly implemented by a hardware processor for execution, or can be executed by a combination of hardware and software modules in the processor.
  • the memory 502 can be used to store non-volatile software programs, non-volatile computer executable programs and modules.
  • the memory 502 may include at least one type of storage medium, for example, may include flash memory, hard disk, multimedia card, card-type memory, random access memory (Random Access Memory, RAM), static random access memory (Static Random Access Memory, SRAM), Programmable Read Only Memory (PROM), Read Only Memory (ROM), Electrically Erasable Programmable Read-Only Memory (EEPROM), Magnetic Memory, Disk , CD, etc.
  • Memory 502 is, but is not limited to, any other medium that can be used to carry or store desired program code in the form of instructions or data structures and that can be accessed by a computer.
  • the memory 502 in the embodiment of the present application can also be a circuit or any other device capable of realizing a storage function, used to store program instructions and/or data.
  • embodiments of the present invention also provide a computer-readable storage medium.
  • the computer-readable storage medium stores a computer executable program.
  • the computer executable program is used to cause the computer to execute any of the applications listed above. Methods for cycle slip detection and repair of Global Navigation Satellite System GNSS.
  • embodiments of the present application may be provided as methods, systems, or computer program products. Accordingly, the present application may take the form of an entirely hardware embodiment, an entirely software embodiment, or an embodiment that combines software and hardware aspects. Furthermore, the present application may take the form of a computer program product embodied on one or more computer-usable storage media (including, but not limited to, disk storage, CD-ROM, optical storage, etc.) having computer-usable program code embodied therein.
  • computer-usable storage media including, but not limited to, disk storage, CD-ROM, optical storage, etc.
  • These computer program instructions may also be stored in a computer-readable memory that causes a computer or other programmable data processing apparatus to operate in a particular manner, such that the instructions stored in the computer-readable memory produce an article of manufacture including the instruction means, the instructions
  • the device implements the functions specified in a process or processes of the flowchart and/or a block or blocks of the block diagram.
  • These computer program instructions may also be loaded onto a computer or other programmable data processing device, causing a series of operating steps to be performed on the computer or other programmable device to produce computer-implemented processing, thereby executing on the computer or other programmable device.
  • Instructions provide steps for implementing the functions specified in a process or processes of a flowchart diagram and/or a block or blocks of a block diagram.

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

一种应用于全球导航卫星系统GNSS的周跳探测与修复方法，包括：根据监测站在t历元接收的s2个共视卫星的载波相位观测值确定监测站在t历元的s2个载波相位三差观测值；根据监测站在t历元接收的s2个共视卫星的伪距观测值确定监测站在t历元的s2个伪距三差观测值（101）；将t历元对应的s2个载波相位三差观测值和s2个伪距三差观测值代入基于几何GB模型的三差观测方程（102）；基于s2个载波相位三差观测值和s2个伪距三差观测值，对2s2个求解式进行求解，解得s2个共视卫星在t历元的各整数周跳值；若任一整数周跳值满足预设条件，则确定整数周跳值对应的共视卫星在观测过程中发生周跳（103）。适合复杂环境下的周跳探测，得到更高精度的周跳值。

Description

一种应用于全球导航卫星系统GNSS的周跳探测与修复方法
相关申请的交叉引用
本申请要求在2022年04月13日提交中国专利局、申请号为202210387255.2、申请名称为“一种应用于全球导航卫星系统GNSS的周跳探测与修复方法”的中国专利申请的优先权，其全部内容通过引用结合在本申请中。
技术领域
本发明实施例涉及计算机技术领域，尤其涉及一种应用于全球导航卫星系统GNSS的周跳探测与修复方法、装置、计算设备及计算机可读存储介质。
背景技术
目前，全球导航卫星系统(Global Navigation Satellite System，GNSS)的卫星能够为用户提供如伪距和载波相位观测值在内的多种观测值。GNSS卫星在实时连续观测过程中，卫星信号不可避免地会受到观测环境的影响，尤其是在复杂观测环境下，例如卫星信号被障碍物挡住而暂时中断、受无线电信号干扰造成失锁、GNSS接收机载波相位跟踪环失锁、低信噪比等条件都将导致周跳现象的发生。周跳是指在载波相位测量中，由于卫星信号的失锁而导致的整周计数的跳变或中断。因此，当信号重新被跟踪后，整周计数就不正确，因此正确地探测并恢复周跳，是载波相位测量中非常重要且必须解决的问题之一，也是实现高精度且可靠的导航与定位的重要前提。
通常，针对单频观测数据，常用周跳探测方法有高次差法、多项式拟合法和伪距相位组合法等，这些方法本质上均是基于几何无关(Geometry-free，GF)组合而展开。而GF组合受电离层和观测值噪声影响较大，尤其是在复杂环境下GF组合周跳探测的性能往往欠佳且无法很好地探测出小周跳。针对双频观测数据，Turbo Edit方法是最为广泛使用的一种周跳探测方法，其本质上包含GF和MW组合，该方法同样受电离层和观测噪声的影响较大。然而，在复杂观测环境下，观测值噪声往往较大，从而导致上述这些方法均无法准确有效地探测周跳并修复周跳。
发明内容
本发明实施例提供一种应用于全球导航卫星系统GNSS的周跳探测与修复方法，用以降低电离层和观测噪声对周跳探测的影响，提高周跳探测的精度和可靠性。
第一方面，本发明实施例提供一种应用于全球导航卫星系统GNSS的周跳探测方法， 包括：
根据监测站在t历元接收的s2个共视卫星的载波相位观测值确定所述监测站在所述t历元的s2个载波相位三差观测值；根据监测站在t历元接收的s2个共视卫星的伪距观测值确定所述监测站在所述t历元的s2个伪距三差观测值；其中s2≥3；
将所述t历元对应的s2个载波相位三差观测值和s2个伪距三差观测值代入基于几何GB模型的三差观测方程，得到2s2个求解式；所述GB模型的三差观测方程为：
Figure PCTCN2022112718-appb-000001
其中，
Figure PCTCN2022112718-appb-000002
u为参考星；s为共视卫星；q为参考站；r为监测站；所述监测站与所述参考站之间的基线长度小于第一预设阈值；获取观测值的历元之间的间隔小于第二预设阈值；
Figure PCTCN2022112718-appb-000003
为监测站r在t历元针对共视卫星s的载波相位三差观测值；
Figure PCTCN2022112718-appb-000004
为监测站r在t历元针对共视卫星s的伪距三差观测值；A t为t历元的设计矩阵，用以表示t历元监测站r与共视卫星s之间的几何关系；A t-1为t-1历元的设计矩阵，用以表示t-1历元监测站r与共视卫星s之间的几何关系；x为以坐标参数的形式表示的t历元的卫地距离双差值；b t为以坐标参数的形式表示的t历元的卫地距离三差值；λ为接收的卫星信号的波长；
Figure PCTCN2022112718-appb-000005
为共视卫星s的周跳；
Figure PCTCN2022112718-appb-000006
为载波相位观测值噪声；
Figure PCTCN2022112718-appb-000007
为伪距观测值噪声；
基于所述s2个载波相位三差观测值和所述s2个伪距三差观测值，利用最小二乘估计对所述2s2个求解式进行求解，解得s2个共视卫星在t历元的各整数周跳值；若任一整数周跳值满足预设条件，则确定所述整数周跳值对应的共视卫星在观测过程中发生周跳。
将t历元对应的s2个载波相位三差观测值和s2个伪距三差观测值代入基于几何GB模型的三差观测方程，得到2s2个求解式，其中GB模型考虑了监测站与共视卫星之间的几何关系，对观测噪声的容忍度更高。由于监测站与所述参考站之间的基线长度小于第一预设阈值以及获取观测值的历元之间的间隔小于第二预设阈值，那么可忽略历元间站星双差电离层、对流层以及多路径的影响，从而可将GB模型的三差观测方程更加简化。由于卫星距离监测站较远，因此相邻历元之间的设计矩阵的差值约等于0，因此对Δ(Ax)可以进行变形，变形后减少了之后需要求解的待解参数的数量，满足最小二乘可解可估的条件，从而解决了方程病态以及不可解不可估的问题。进一步提高了对观测噪声的容忍度，更加适合于复杂环境下的周跳探测，能够得到更高精度和更高可靠性的整数周跳值。
可选地，解得s2个共视卫星在t历元的各整数周跳值之前，还包括：
根据监测站在t-1历元接收的s1个共视卫星的载波相位观测值确定所述监测站在所述t-1历元的s1个载波相位三差观测值；根据监测站在t-1历元接收的s1个共视卫星的伪距 观测值确定所述监测站在所述t-1历元的s1个伪距三差观测值；所述t-1历元早于所述t历元；所述s1个共视卫星在所述t-1历元不存在周跳现象或周跳已修复；
将所述t-1历元对应的s1个载波相位三差观测值和s1个伪距三差观测值代入所述GB模型的三差观测方程，得到2s1个求解式；
基于所述s2个载波相位三差观测值和所述s2个伪距三差观测值，利用最小二乘估计对所述2s2个求解式进行求解，解得s2个共视卫星在t历元的各整数周跳值，包括：
基于所述s2个载波相位三差观测值、所述s2个伪距三差观测值、所述s1个载波相位三差观测值和所述s1个伪距三差观测值，利用最小二乘估计对2(s1+s2)个求解式进行求解，解得s2个共视卫星在t历元的各整数周跳值。
为了进一步增加求解式中多余观测值的数量，使求解式的数量大于待解参数的数量，从而满足最小二乘可解可估的条件，降低观测噪声对周跳探测与修复的影响，进一步提高解得的周跳值的准确性和可靠性，可以进一步获取t-1历元的s1个载波相位三差观测值和s1个伪距三差观测值。最终得到2(s1+s2)个求解式，但是由于t-1历元对应的s1个共视卫星不存在周跳或周跳已修正，因此并不会引入新的关于周跳值的待解参数，解得的周跳值更加准确，对观测噪声的容忍度更高。
可选地，所述监测站为静态监测站；所述b t与b t-1相等；所述b t-1为以坐标参数的形式表示的t-1历元的卫地距离三差值。
通过添加位置域约束，进一步减少待解参数的数量，提高求解式可解可估的概率，并增加求解式中多余观测值的数量，进一步提高解得的周跳值的准确性，降低观测噪声的影响。求解出的周跳值精度和稳定性更高。
可选地，基于所述s2个载波相位三差观测值和所述s2个伪距三差观测值，利用最小二乘估计对所述2s2个求解式进行求解，解得s2个共视卫星在t历元的各整数周跳值，包括：
将所述2s2个求解式用矩阵表示成如下形式：
E(l)＝BX；
其中，l为所述s2个载波相位观测值和所述s2个伪距三差观测值形成的矩阵；E(l)为所述s2个载波相位三差观测值和所述s2个伪距三差观测值的期望；B为
Figure PCTCN2022112718-appb-000008
在x、y和z三个方向上的分量及所述s2个共视卫星的卫星信号的波长形成的矩阵；X为所述b t在x、y和z三个方向上的分量及所述s2个共视卫星的各
Figure PCTCN2022112718-appb-000009
形成的矩阵；
通过最小二乘估计的如下公式求得所述X：
Figure PCTCN2022112718-appb-000010
其中，所述D3为所述l的方差协方差矩阵；
根据所述X的矩阵的解算结果，得到所述s2个共视卫星在t历元的各整数周跳值。
采用上述最小二乘法解得各整数周跳值，更加准确。
可选地，根据所述X的矩阵的解算结果，得到所述s2个共视卫星在t历元的各整数周跳值，包括：
获取所述X的矩阵的解算结果中的s2个共视卫星在t历元的浮点形式的各周跳值；
确定各浮点形式的周跳值的方差协方差阵；
基于各浮点形式的周跳值以及各浮点形式的周跳值的方差协方差阵，将各浮点形式的周跳值转换为整数形式的周跳值，得到所述s2个共视卫星在t历元的各整数周跳值。
由于浮点解并不能真实反映各共视卫星的周跳情况，因此对浮点形式的各周跳值确定方差协方差阵，根据各浮点形式的周跳值以及各浮点形式的周跳值的方差协方差阵，得到整数形式的周跳值，整数形式的周跳值才能反映周跳现象的发生，基于整数形式的周跳值才能实现周跳的修复。
可选地，通过如下公式确定各浮点形式的周跳值的方差协方差阵，包括：
(B TD 3 -1B) -1。
可选地，所述载波相位观测值包括多种波长的载波相位观测值；每种波长的载波相位观测值确定各自对应的载波相位三差观测值；
所述伪距观测值包括多种波长的伪距观测值，每种波长的伪距观测值确定各自对应的伪距三差观测值。
第二方面，本发明实施例还提供一种应用于全球导航卫星系统GNSS的周跳探测与修复装置，包括：
确定单元，用于：
根据监测站在t历元接收的s2个共视卫星的载波相位观测值确定所述监测站在所述t历元的s2个载波相位三差观测值；根据监测站在t历元接收的s2个共视卫星的伪距观测值确定所述监测站在所述t历元的s2个伪距三差观测值；其中s2≥3；
将所述t历元对应的s2个载波相位三差观测值和s2个伪距三差观测值代入基于几何GB模型的三差观测方程，得到2s2个求解式；所述GB模型的三差观测方程为：
Figure PCTCN2022112718-appb-000011
其中，
Figure PCTCN2022112718-appb-000012
u为参考星；s为共视卫星；q为参考站；r为监测站；所述监测站与所述参考站之间的基线长度小于第一预设阈值；获取观测值的历元之间的间隔小 于第二预设阈值；
Figure PCTCN2022112718-appb-000013
为监测站r在t历元针对共视卫星s的载波相位三差观测值；
Figure PCTCN2022112718-appb-000014
为监测站r在t历元针对共视卫星s的伪距三差观测值；A t为t历元的设计矩阵，用以表示t历元监测站r与共视卫星s之间的几何关系；A t-1为t-1历元的设计矩阵，用以表示t-1历元监测站r与共视卫星s之间的几何关系；x为以坐标参数的形式表示的t历元的卫地距离双差值；b t为以坐标参数的形式表示的t历元的卫地距离三差值；λ为接收的卫星信号的波长；
Figure PCTCN2022112718-appb-000015
为共视卫星s的周跳；
Figure PCTCN2022112718-appb-000016
为载波相位观测值噪声；
Figure PCTCN2022112718-appb-000017
为伪距观测值噪声；
计算单元，用于基于所述s2个载波相位三差观测值和所述s2个伪距三差观测值，利用最小二乘估计对所述2s2个求解式进行求解，解得s2个共视卫星在t历元的各整数周跳值；若任一整数周跳值满足预设条件，则确定所述整数周跳值对应的共视卫星在观测过程中发生周跳。
可选地，所述确定单元还用于：
根据监测站在t-1历元接收的s1个共视卫星的载波相位观测值确定所述监测站在所述t-1历元的s1个载波相位三差观测值；根据监测站在t-1历元接收的s1个共视卫星的伪距观测值确定所述监测站在所述t-1历元的s1个伪距三差观测值；所述t-1历元早于所述t历元；所述s1个共视卫星在所述t-1历元不存在周跳现象或周跳已修复；
将所述t-1历元对应的s1个载波相位三差观测值和s1个伪距三差观测值代入所述GB模型的三差观测方程，得到2s1个求解式；
所述计算单元，具体用于：
基于所述s2个载波相位三差观测值、所述s2个伪距三差观测值、所述s1个载波相位三差观测值和所述s1个伪距三差观测值，利用最小二乘估计对2(s1+s2)个求解式进行求解，解得s2个共视卫星在t历元的各整数周跳值。
可选地，所述监测站为静态监测站；所述b t与b t-1相等；所述b t-1为以坐标参数的形式表示的t-1历元的卫地距离三差值。
可选地，所述计算单元，具体用于：
将所述2s2个求解式用矩阵表示成如下形式：
E(l)＝BX；
其中，l为所述s2个载波相位观测值和所述s2个伪距三差观测值形成的矩阵；E(l)为所述s2个载波相位三差观测值和所述s2个伪距三差观测值的期望；B为
Figure PCTCN2022112718-appb-000018
在x、y和z三个方向上的分量及所述s2个共视卫星的卫星信号的波长形成的矩阵；X为所述b t在x、y和z三个方向上的分量及所述s2个共视卫星的各
Figure PCTCN2022112718-appb-000019
形成的矩阵；
通过最小二乘估计的如下公式求得所述X：
Figure PCTCN2022112718-appb-000020
其中，所述D3为所述l的方差协方差矩阵；
根据所述X的矩阵的解算结果，得到所述s2个共视卫星在t历元的各整数周跳值。
可选地，所述计算单元，具体用于：
获取所述X的矩阵的解算结果中的s2个共视卫星在t历元的浮点形式的各周跳值；
确定各浮点形式的周跳值的方差协方差阵；
基于各浮点形式的周跳值以及各浮点形式的周跳值的方差协方差阵，将各浮点形式的周跳值转换为整数形式的周跳值，得到所述s2个共视卫星在t历元的各整数周跳值。
可选地，所述计算单元，具体用于：
通过如下公式确定各浮点形式的周跳值的方差协方差阵：
(B TD 3 -1B) -1。
可选地，所述载波相位观测值包括多种波长的载波相位观测值；每种波长的载波相位观测值确定各自对应的载波相位三差观测值；
所述伪距观测值包括多种波长的伪距观测值，每种波长的伪距观测值确定各自对应的伪距三差观测值。
第三方面，本发明实施例还提供一种计算设备，包括：
存储器，用于存储计算机程序；
处理器，用于调用所述存储器中存储的计算机程序，按照获得的程序执行上述任一方式所列的应用于GNSS的周跳探测与修复方法。
第四方面，本发明实施例还提供一种计算机可读存储介质，所述计算机可读存储介质存储有计算机可执行程序，所述计算机可执行程序用于使计算机执行上述任一方式所列的应用于GNSS的周跳探测与修复方法。
附图说明
为了更清楚地说明本发明实施例中的技术方案，下面将对实施例描述中所需要使用的附图作简要介绍，显而易见地，下面描述中的附图仅仅是本发明的一些实施例，对于本领域的普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。
图1为本发明实施例提供的一种应用于全球导航卫星系统GNSS的周跳探测与修复方法的流程示意图；
图2为本发明实施例提供的一种应用场景的示意图；
图3为本发明实施例提供的一种应用于全球导航卫星系统GNSS的周跳探测与修复方法的流程示意图；
图4为本发明实施例提供的一种应用于全球导航卫星系统GNSS的周跳探测与修复装置的结构示意图；
图5为本发明实施例提供的一种计算机设备的结构示意图。
具体实施方式
为使本申请的目的、实施方式和优点更加清楚，下面将结合本申请示例性实施例中的附图，对本申请示例性实施方式进行清楚、完整地描述，显然，所描述的示例性实施例仅是本申请一部分实施例，而不是全部的实施例。
基于本申请描述的示例性实施例，本领域普通技术人员在没有作出创造性劳动前提下所获得的所有其他实施例，都属于本申请所附权利要求保护的范围。此外，虽然本申请中公开内容按照示范性一个或几个实例来介绍，但应理解，可以就这些公开内容的各个方面也可以单独构成一个完整实施方式。
需要说明的是，本申请中对于术语的简要说明，仅是为了方便理解接下来描述的实施方式，而不是意图限定本申请的实施方式。除非另有说明，这些术语应当按照其普通和通常的含义理解。
本申请中说明书和权利要求书及上述附图中的术语“第一”、“第二”、“第三”等是用于区别类似或同类的对象或实体，而不必然意味着限定特定的顺序或先后次序，除非另外注明(Unless otherwise indicated)。应该理解这样使用的用语在适当情况下可以互换，例如能够根据本申请实施例图示或描述中给出那些以外的顺序实施。
此外，术语“包括”和“具有”以及他们的任何变形，意图在于覆盖但不排他的包含，例如，包含了一系列组件的产品或设备不必限于清楚地列出的那些组件，而是可包括没有清楚地列出的或对于这些产品或设备固有的其它组件。
全球导航卫星系统是指能在地球表面或近地空间的任何地点为用户提供全天候的3维坐标和速度以及时间信息的空基无线电导航定位系统。目前，全球有四大卫星导航定位系统，包括中国的北斗卫星导航系统(BeiDou navigation satellite system，BDS)、美国的全球定位系统(global positioningsystem，GPS)、欧盟的伽利略卫星导航系统(Galileo navigation satellite system，Galileo)和俄罗斯的格洛纳斯卫星导航系统(global orbitingnavigation satellite system，GLONASS)。北斗卫星导航系统是中国自主建设运行的全球卫星导航系统，为全球用户提供全天候、全天时、高精度的定位、导航和授时服务。
随着全球一体化的发展，卫星导航系统在航空、汽车导航、通信、测绘、娱乐等各个领域均具有广泛的应用。
在一些复杂的观测环境下，GNSS的卫星信号容易发生周跳现象，如何降低复杂的观测环境对周跳探测的影响，从而提高探测周跳的准确度和可靠性，成为了当前周跳探测的重要课题。
为了解决上述问题，本发明实施例提供一种应用于全球导航卫星系统GNSS的周跳探测与修复方法，如图1所示，包括：
步骤101，根据监测站在t历元接收的s2个共视卫星的载波相位观测值确定所述监测站在所述t历元的s2个载波相位三差观测值；根据监测站在t历元接收的s2个共视卫星的伪距观测值确定所述监测站在所述t历元的s2个伪距三差观测值；其中s2≥3；
图2示出了本发明实施例适用的一种应用场景，包括参考星u、共视卫星s、参考站q和监测站r，在监测站r上设置有接收机，共视卫星可以有多个，在图中仅示意出其中的一个共视卫星s。在任一历元，监测站r可以接收到多个共视卫星的载波相位观测值和伪距观测值，将多个共视卫星的载波相位观测值和伪距观测值发送至计算设备，计算设备根据多个共视卫星的载波相位观测值确定每个共视卫星的载波相位三差观测值，根据多个共视卫星的伪距观测值确定每个共视卫星的伪距三差观测值。
具体的，在计算载波相位三差观测值时会用到t历元监测站r接收到的共视卫星s的载波相位观测值、t历元监测站r接收到的参考星u的载波相位观测值、t历元参考站q接收到的共视卫星s的载波相位观测值和t历元参考站q接收到的参考星u的载波相位观测值，以及，t-1历元监测站r接收到的共视卫星s的载波相位观测值、t-1历元监测站r接收到的参考星u的载波相位观测值、t-1历元参考站q接收到的共视卫星s的载波相位观测值和t-1历元参考站q接收到的参考星u的载波相位观测值，具体的如何确定载波相位三差观测值的方式为本领域技术人员熟知的技术手段，在此不作赘述。确定伪距三差观测值的方法同理。
为了保证之后确定的各求解式可解可估，这里限制s2≥3，具体的原因在步骤103中详述。
举个例子，若s2＝3，则包括3个共视卫星。根据监测站r接收到的3个共视卫星的载波相位观测值确定了相应的3个载波相位三差观测值；根据监测站r接收的3个共视卫星的伪距观测值确定了相应的3个伪距三差观测值。
步骤102，将所述t历元对应的s2个载波相位三差观测值和s2个伪距三差观测值代入基于几何GB模型的三差观测方程，得到2s2个求解式。
GB模型相较于GF模型来说，考虑了监测站与共视卫星之间的几何关系，对观测噪声的容忍度更高，可以降低观测噪声对周跳探测的影响。本发明实施例对传统的GB模型进行了优化。
下面介绍本发明实施例得到优化后的GB模型的过程：
传统的GF模型的三差观测方程：
Figure PCTCN2022112718-appb-000021
其中，u为参考星；s为共视卫星；q为参考站；r为监测站；
Figure PCTCN2022112718-appb-000022
为监测站r在t历元针对共视卫星s的载波相位三差观测值；
Figure PCTCN2022112718-appb-000023
为监测站r在t历元针对共视卫星s的伪距三差观测值；
Figure PCTCN2022112718-appb-000024
为t历元的卫地距离三差值；λ为接收的卫星信号的波长；
Figure PCTCN2022112718-appb-000025
为共视卫星s的周跳；μ为电离层系数；
Figure PCTCN2022112718-appb-000026
为电离层延迟；
Figure PCTCN2022112718-appb-000027
为对流层延迟；
Figure PCTCN2022112718-appb-000028
为载波相位多路径误差；
Figure PCTCN2022112718-appb-000029
为载波相位观测值噪声；
Figure PCTCN2022112718-appb-000030
为伪距多路径误差；
Figure PCTCN2022112718-appb-000031
为伪距观测值噪声。
可以看出，传统的GF模型的三差观测方程中，卫地距离三差观测值并没有考虑监测站与共视卫星之间的几何关系，常用本领域技术人员熟知的GF线性组合的方式进行周跳探测，因此求得的周跳值受观测噪声的影响较大，无法进行高精度的周跳探测。
将GF模型中的
Figure PCTCN2022112718-appb-000032
以坐标参数形式表示，得到传统的GB模型的三差观测方程：
Figure PCTCN2022112718-appb-000033
其中，Δ(Ax)＝A tx t-A t-1x t-1，A t为t历元的设计矩阵，用以表示t历元监测站r与共视卫星s之间的几何关系；A t-1为t-1历元的设计矩阵，用以表示t-1历元监测站r与共视卫星s之间的几何关系；x为以坐标参数的形式表示的t历元的卫地距离双差值。其中设计矩阵为根据共视卫星的坐标和监测站的大概已知坐标(不够精确)提前确定，为已知数值。
在本发明实施例中，以“加粗斜体”表示的参数代表该参数为矩阵形式。
为了能够减少未知数的数量，提高周跳探测的准确性，本发明实施例限制监测站与参考站之间的基线长度小于第一预设阈值，例如小于10公里、9公里等，如此，可以忽略历元间站星双差电离层
Figure PCTCN2022112718-appb-000034
对流层
Figure PCTCN2022112718-appb-000035
的影响；并限制采样间隔，即获取观测值的历元之间的间隔小于第二预设阈值，例如小于5秒、小于3秒等。如此，可以忽略多路径
Figure PCTCN2022112718-appb-000036
Figure PCTCN2022112718-appb-000037
的影响。
从而，得到如下GB模型的三差观测方程：
Figure PCTCN2022112718-appb-000038
其中，Δ(Ax)＝A tx t-A t-1x t-1。
在上述GB模型中，载波相位三差观测值和伪距三差观测值可以根据监测站接收到的卫星信号(载波相位观测值和伪距观测值)得到(即，步骤101中所述)，因此为已知量；设计矩阵为已知量，接收到的卫星信号的波长也为已知量，如此，x t、x t-1、周跳值
Figure PCTCN2022112718-appb-000039
为待解参数，需要对x t、x t-1、周跳值
Figure PCTCN2022112718-appb-000040
进行求解(观测值噪声无需求解)。
为了使之后利用最小二乘法对周跳值进行求解时得到的求解式不病态，且可解可估，我们想到，可以尽可能减少待解参数的数量，因此对Δ(Ax)进行如下变形：
Figure PCTCN2022112718-appb-000041
经过研究，我们发现，采样间隔一般较小，一般相隔1s-5s，因此A t和A t-1中元素差异的量级在10 -3，因此，A t-A t-1≈0。
那么上式可进一步变为：
Figure PCTCN2022112718-appb-000042
其中，b t为以坐标参数的形式表示的t历元的卫地距离三差值。
经过上述变形后，可以看出，x t和x t-1两个待解参数合并为一个待解参数b t，在后续通过求解式求解时，只需求解b t即可，如此，保证了方程不病态，且可解可估。
最终优化后的GB模型的三差观测方程为：
Figure PCTCN2022112718-appb-000043
为了简化表示，令
Figure PCTCN2022112718-appb-000044
因此上述式子还可写作：
Figure PCTCN2022112718-appb-000045
在上述方程中，需要求解的待解参数为b t和周跳值
Figure PCTCN2022112718-appb-000046
在本发明实施例中，主要关心周跳值
Figure PCTCN2022112718-appb-000047
将在步骤101中得到的t历元对应的s2个载波相位三差观测值和s2个伪距三差观测值代入优化后的GB模型的三差观测方程，可以得到2s2个求解式，若s2等于3，则可以得到6个求解式，假设3个共视卫星分别为1、2和3，如下所示：
Figure PCTCN2022112718-appb-000048
步骤103，基于所述s2个载波相位三差观测值和所述s2个伪距三差观测值，利用最小二乘估计对所述2s2个求解式进行求解，解得s2个共视卫星在t历元的各整数周跳值；若任一整数周跳值满足预设条件，则确定所述整数周跳值对应的共视卫星在观测过程中发生周跳。
将2s2个求解式用矩阵表示成如下形式：
E(l)＝BX；
其中，l为所述s2个载波相位观测值和所述s2个伪距三差观测值形成的矩阵；E(l)为所述s2个载波相位三差观测值和所述s2个伪距三差观测值的期望；B为
Figure PCTCN2022112718-appb-000049
(即Gt)在x、y和z三个方向上的分量及所述s2个共视卫星的卫星信号的波长形成的矩阵；X为所述b t在x、y和z三个方向上的分量及所述s2个共视卫星的各
Figure PCTCN2022112718-appb-000050
形成的矩阵；
以步骤102中得到的6个求解式为例：
Figure PCTCN2022112718-appb-000051
Figure PCTCN2022112718-appb-000052
Figure PCTCN2022112718-appb-000053
综合起来即：
Figure PCTCN2022112718-appb-000054
在上式中，X为待解参数，通过最小二乘估计的如下公式求得所述X：
Figure PCTCN2022112718-appb-000055
其中，所述D3为所述l的方差协方差矩阵；
通过X的矩阵的求解结果中，可以得到3个共视卫星各自的整数周跳值
Figure PCTCN2022112718-appb-000056
和
Figure PCTCN2022112718-appb-000057
在上述式子中，待解参数有6个，求解式有6个，求解式的个数刚好等于待解参数的个数，因此求解式刚好可解可估。若设置s2＜3，例如s2＝2，则，求解式为4个，待解参数为5个，求解式不可解不可估。经过研究发现，若只获取t历元的共视卫星的载波相位值和伪距值，那么s2≥3。
在直接得到的X的矩阵的求解结果中，各周跳值可能为浮点解，此时需要采用最小二乘降相关平差(LAMBDA)算法将浮点解固定为整数。具体为：
通过双差观测值的方差协方差矩阵D2求得三差观测值l的方差协方差矩阵D3：
Figure PCTCN2022112718-appb-000058
其中，C为双差观测值转三差观测值的设计矩阵。
求得D3后，通过X的求解公式求得X中各周跳值的各浮点解，即s2个共视卫星在t历元的浮点形式的各周跳值；
通过公式(B TD 3 -1B) -1确定各浮点解的方差协方差阵；
基于各浮点解以及各浮点解的方差协方差阵，通过LAMBDA算法将各浮点形式的周跳值转换为整数形式的周跳值，得到所述s2个共视卫星在t历元的各整数周跳值。
若各整数周跳值中有非0元素，则说明该非0元素对应的共视卫星发生周跳。
可选地，为了进一步增加求解式多余观测值的数量，降低观测噪声对周跳探测的影响，进一步提高解得的周跳值的准确性，我们还可获取t-1历元的共视卫星的载波相位值和伪距值。
具体为，解得s2个共视卫星在t历元的各整数周跳值之前，还包括如下步骤，如图3所示：
步骤301，根据监测站在t-1历元接收的s1个共视卫星的载波相位观测值确定所述监 测站在所述t-1历元的s1个载波相位三差观测值；根据监测站在t-1历元接收的s1个共视卫星的伪距观测值确定所述监测站在所述t-1历元的s1个伪距三差观测值；所述t-1历元早于所述t历元；所述s1个共视卫星在所述t-1历元不存在周跳现象或周跳已修复。
之所以设置t-1历元对应的s1个共视卫星不存在周跳现象，或，存在周跳现象，但是对t-1历元对应的s1个共视卫星的周跳进行修复，是为了不引入新的待解参数，以进一步增加求解式多余观测值的数量，降低观测噪声对周跳探测的影响。
由于不同时刻能够建立联系的共视卫星不同，因此t-1历元对应的s1个共视卫星和t历元对应的s2个共视卫星可能相同，也可能不同，s1可能等于s2，也可能不等于，对此不做限制。
步骤302，将所述t-1历元对应的s1个载波相位三差观测值和s1个伪距三差观测值代入所述GB模型的三差观测方程，得到2s1个求解式。
具体得到的2s1个求解式可以参考步骤102中得到的2s2个求解式中的例子，在此不作赘述。
那么，解得s2个周跳值的方法为：
基于所述s2个载波相位三差观测值、所述s2个伪距三差观测值、所述s1个载波相位三差观测值和所述s1个伪距三差观测值，利用最小二乘估计对2(s1+s2)个求解式进行求解，解得s2个共视卫星在t历元的各整数周跳值。
通过上述方式，引入了2s1个求解式，但同时也引入了3个待解参数：b t-1x、b t-1y、b t-1z。b t-1为以坐标参数的形式表示的t-1历元的卫地距离三差值。
未引入t-1历元的三差观测值之前，求解式的数量是2s2，待解参数的数量是3+s2；引入t-1历元的三差观测值之后，求解式的数量是2(s1+s2)，待解参数的数量是3+3+s2。因此为了使求解式可解，需要保证2(s1+s2)≥6+s2。
可选地，为了进一步减少待解参数的数量，提高求解式可解可估的概率，并增加求解式的多余观测值，进一步提高解得的周跳值的准确性，降低观测噪声的影响，可以添加位置域约束，即，令监测站为静态监测站，那么b t＝b t-1，进而b t-1x＝b tx、b t-1y＝b ty、b t-1z＝b tz。
如此，不同历元的卫地距离三差值相等，引入t-1历元的三差观测值，但是不会引入新的待解参数。未引入t-1历元的三差观测值之前，求解式的数量是2s2，待解参数的数量是3+s2；引入t-1历元的三差观测值之后，求解式的数量是2(s1+s2)，待解参数的数量是3+s2。因此为了使求解式可解可估，需要保证2(s1+s2)≥3+s2。此时，即便s2＝1，s1＝1也能实现求解式可解可估，降低了对s2和s1的数量要求。求解出的周跳值精度和稳定性更高。
举个例子，通过添加位置域约束，得到的2(s1+s2)个求解式通过如下形式表示：
Figure PCTCN2022112718-appb-000059
对上述2(s1+s2)个求解式进行求解的方式同步骤103中对2s2个求解式求解的方式，在此不再赘述。
可选地，本发明实施例提供的方法中，观测值不仅可以为单频观测值，还可以为双频、三频等多频观测值。具体来说，载波相位观测值包括多种波长的载波相位观测值；每种波长的载波相位观测值确定各自对应的载波相位三差观测值；伪距观测值包括多种波长的伪距观测值，每种波长的伪距观测值确定各自对应的伪距三差观测值。
以观测值为双频观测值为例，则之前确定的求解式可表示成如下紧凑形式：
Figure PCTCN2022112718-appb-000060
其中l t-1＝[φ t-1,P t-1] T；l t＝[φ t,P t] T；x＝[b tx,b ty,b tz] T；z＝[z i,z j] T；
Figure PCTCN2022112718-appb-000061
e 4为4×1的全1向量，I s2为s 2维的单位矩阵，λ i和λ j分别为频率i和频率j的波长。
其中，
Figure PCTCN2022112718-appb-000062
Figure PCTCN2022112718-appb-000063
Figure PCTCN2022112718-appb-000064
Figure PCTCN2022112718-appb-000065
Figure PCTCN2022112718-appb-000066
Figure PCTCN2022112718-appb-000067
基于相同的技术构思，图4示例性的示出了本发明实施例提供的一种应用于全球导航 卫星系统GNSS的周跳探测与修复的装置的结构，该结构可以执行应用于全球导航卫星系统GNSS的周跳探测与修复的流程。
如图4所示，该装置具体包括：
确定单元401，用于：
根据监测站在t历元接收的s2个共视卫星的载波相位观测值确定所述监测站在所述t历元的s2个载波相位三差观测值；根据监测站在t历元接收的s2个共视卫星的伪距观测值确定所述监测站在所述t历元的s2个伪距三差观测值；其中s2≥3；
将所述t历元对应的s2个载波相位三差观测值和s2个伪距三差观测值代入基于几何GB模型的三差观测方程，得到2s2个求解式；所述GB模型的三差观测方程为：
Figure PCTCN2022112718-appb-000068
其中，
Figure PCTCN2022112718-appb-000069
u为参考星；s为共视卫星；q为参考站；r为监测站；所述监测站与所述参考站之间的基线长度小于第一预设阈值；获取观测值的历元之间的间隔小于第二预设阈值；
Figure PCTCN2022112718-appb-000070
为监测站r在t历元针对共视卫星s的载波相位三差观测值；
Figure PCTCN2022112718-appb-000071
为监测站r在t历元针对共视卫星s的伪距三差观测值；A t为t历元的设计矩阵，用以表示t历元监测站r与共视卫星s之间的几何关系；A t-1为t-1历元的设计矩阵，用以表示t-1历元监测站r与共视卫星s之间的几何关系；x为以坐标参数的形式表示的t历元的卫地距离双差值；b t为以坐标参数的形式表示的t历元的卫地距离三差值；λ为接收的卫星信号的波长；
Figure PCTCN2022112718-appb-000072
为共视卫星s的周跳；
Figure PCTCN2022112718-appb-000073
为载波相位观测值噪声；
Figure PCTCN2022112718-appb-000074
为伪距观测值噪声；
计算单元402，用于基于所述s2个载波相位三差观测值和所述s2个伪距三差观测值，利用最小二乘估计对所述2s2个求解式进行求解，解得s2个共视卫星在t历元的各整数周跳值；若任一整数周跳值满足预设条件，则确定所述整数周跳值对应的共视卫星在观测过程中发生周跳。
可选地，所述确定单元401还用于：
根据监测站在t-1历元接收的s1个共视卫星的载波相位观测值确定所述监测站在所述t-1历元的s1个载波相位三差观测值；根据监测站在t-1历元接收的s1个共视卫星的伪距观测值确定所述监测站在所述t-1历元的s1个伪距三差观测值；所述t-1历元早于所述t历元；所述s1个共视卫星在所述t-1历元不存在周跳现象或周跳已修复；
将所述t-1历元对应的s1个载波相位三差观测值和s1个伪距三差观测值代入所述GB模型的三差观测方程，得到2s1个求解式；
所述计算单元402，具体用于：
基于所述s2个载波相位三差观测值、所述s2个伪距三差观测值、所述s1个载波相位三差观测值和所述s1个伪距三差观测值，利用最小二乘估计对2(s1+s2)个求解式进行求解，解得s2个共视卫星在t历元的各周跳值。
可选地，所述监测站为静态监测站；所述b t与b t-1相等；所述b t-1为以坐标参数的形式表示的t-1历元的卫地距离三差值。
可选地，所述计算单元402，具体用于：
将所述2s2个求解式用矩阵表示成如下形式：
E(l)＝BX；
其中，l为所述s2个载波相位观测值和所述s2个伪距三差观测值形成的矩阵；E(l)为所述s2个载波相位三差观测值和所述s2个伪距三差观测值的期望；B为
Figure PCTCN2022112718-appb-000075
在x、y和z三个方向上的分量及所述s2个共视卫星的卫星信号的波长形成的矩阵；X为所述b t在x、y和z三个方向上的分量及所述s2个共视卫星的各
Figure PCTCN2022112718-appb-000076
形成的矩阵；
通过最小二乘估计的如下公式求得所述X：
Figure PCTCN2022112718-appb-000077
其中，所述D3为所述l的方差协方差矩阵；
根据所述X的矩阵的解算结果，得到所述s2个共视卫星在t历元的各整数周跳值。
可选地，所述计算单元402，具体用于：
获取所述X的矩阵的解算结果中的s2个共视卫星在t历元的浮点形式的各周跳值；
确定各浮点形式的周跳值的方差协方差阵；
基于各浮点形式的周跳值以及各浮点形式的周跳值的方差协方差阵，将各浮点形式的周跳值转换为整数形式的周跳值，得到所述s2个共视卫星在t历元的各整数周跳值。
可选地，所述计算单元402，具体用于：
通过如下公式确定各浮点形式的周跳值的方差协方差阵：
(B TD 3 -1B) -1。
可选地，所述载波相位观测值包括多种波长的载波相位观测值；每种波长的载波相位观测值确定各自对应的载波相位三差观测值；
所述伪距观测值包括多种波长的伪距观测值，每种波长的伪距观测值确定各自对应的伪距三差观测值。
基于相同的技术构思，本申请实施例提供了一种计算机设备，如图5所示，包括至少一个处理器501，以及与至少一个处理器连接的存储器502，本申请实施例中不限定处理器501与存储器502之间的具体连接介质，图5中处理器501和存储器502之间通过总线 连接为例。总线可以分为地址总线、数据总线、控制总线等。
在本申请实施例中，存储器502存储有可被至少一个处理器501执行的指令，至少一个处理器501通过执行存储器502存储的指令，可以执行上述应用于全球导航卫星系统GNSS的周跳探测与修复方法的步骤。
其中，处理器501是计算机设备的控制中心，可以利用各种接口和线路连接计算机设备的各个部分，通过运行或执行存储在存储器502内的指令以及调用存储在存储器502内的数据，从而进行应用于全球导航卫星系统GNSS的周跳探测与修复。可选的，处理器501可包括一个或多个处理单元，处理器501可集成应用处理器和调制解调处理器，其中，应用处理器主要处理操作系统、用户界面和应用程序等，调制解调处理器主要处理无线通信。可以理解的是，上述调制解调处理器也可以不集成到处理器501中。在一些实施例中，处理器501和存储器502可以在同一芯片上实现，在一些实施例中，它们也可以在独立的芯片上分别实现。
处理器501可以是通用处理器，例如中央处理器(CPU)、数字信号处理器、专用集成电路(Application Specific Integrated Circuit，ASIC)、现场可编程门阵列或者其他可编程逻辑器件、分立门或者晶体管逻辑器件、分立硬件组件，可以实现或者执行本申请实施例中公开的各方法、步骤及逻辑框图。通用处理器可以是微处理器或者任何常规的处理器等。结合本申请实施例所公开的方法的步骤可以直接体现为硬件处理器执行完成，或者用处理器中的硬件及软件模块组合执行完成。
存储器502作为一种非易失性计算机可读存储介质，可用于存储非易失性软件程序、非易失性计算机可执行程序以及模块。存储器502可以包括至少一种类型的存储介质，例如可以包括闪存、硬盘、多媒体卡、卡型存储器、随机访问存储器(Random Access Memory，RAM)、静态随机访问存储器(Static Random Access Memory，SRAM)、可编程只读存储器(Programmable Read Only Memory，PROM)、只读存储器(Read Only Memory，ROM)、带电可擦除可编程只读存储器(Electrically Erasable Programmable Read-Only Memory，EEPROM)、磁性存储器、磁盘、光盘等等。存储器502是能够用于携带或存储具有指令或数据结构形式的期望的程序代码并能够由计算机存取的任何其他介质，但不限于此。本申请实施例中的存储器502还可以是电路或者其它任意能够实现存储功能的装置，用于存储程序指令和/或数据。
基于相同的技术构思，本发明实施例还提供一种计算机可读存储介质，计算机可读存储介质存储有计算机可执行程序，计算机可执行程序用于使计算机执行上述任一方式所列的应用于全球导航卫星系统GNSS的周跳探测与修复的方法。
本领域内的技术人员应明白，本申请的实施例可提供为方法、系统、或计算机程序产品。因此，本申请可采用完全硬件实施例、完全软件实施例、或结合软件和硬件方面的实施例的形式。而且，本申请可采用在一个或多个其中包含有计算机可用程序代码的计算机可用存储介质(包括但不限于磁盘存储器、CD-ROM、光学存储器等)上实施的计算机程序产品的形式。
本申请是参照根据本申请的方法、设备(系统)、和计算机程序产品的流程图和/或方框图来描述的。应理解可由计算机程序指令实现流程图和/或方框图中的每一流程和/或方框、以及流程图和/或方框图中的流程和/或方框的结合。可提供这些计算机程序指令到通用计算机、专用计算机、嵌入式处理机或其他可编程数据处理设备的处理器以产生一个机器，使得通过计算机或其他可编程数据处理设备的处理器执行的指令产生用于实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能的装置。
这些计算机程序指令也可存储在能引导计算机或其他可编程数据处理设备以特定方式工作的计算机可读存储器中，使得存储在该计算机可读存储器中的指令产生包括指令装置的制造品，该指令装置实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能。
这些计算机程序指令也可装载到计算机或其他可编程数据处理设备上，使得在计算机或其他可编程设备上执行一系列操作步骤以产生计算机实现的处理，从而在计算机或其他可编程设备上执行的指令提供用于实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能的步骤。
显然，本领域的技术人员可以对本申请进行各种改动和变型而不脱离本申请的精神和范围。这样，倘若本申请的这些修改和变型属于本申请权利要求及其等同技术的范围之内，则本申请也意图包含这些改动和变型在内。

Claims (10)

  1. 一种应用于全球导航卫星系统GNSS的周跳探测与修复方法，其特征在于，包括：
    根据监测站在t历元接收的s2个共视卫星的载波相位观测值确定所述监测站在所述t历元的s2个载波相位三差观测值；根据监测站在t历元接收的s2个共视卫星的伪距观测值确定所述监测站在所述t历元的s2个伪距三差观测值；其中s2≥3；
    将所述t历元对应的s2个载波相位三差观测值和s2个伪距三差观测值代入基于几何GB模型的三差观测方程，得到2s2个求解式；所述GB模型的三差观测方程为：
    Figure PCTCN2022112718-appb-100001
    其中，
    Figure PCTCN2022112718-appb-100002
    u为参考星；s为共视卫星；q为参考站；r为监测站；所述监测站与所述参考站之间的基线长度小于第一预设阈值；获取观测值的历元之间的间隔小于第二预设阈值；
    Figure PCTCN2022112718-appb-100003
    为监测站r在t历元针对共视卫星s的载波相位三差观测值；
    Figure PCTCN2022112718-appb-100004
    为监测站r在t历元针对共视卫星s的伪距三差观测值；A t为t历元的设计矩阵，用以表示t历元监测站r与共视卫星s之间的几何关系；A t-1为t-1历元的设计矩阵，用以表示t-1历元监测站r与共视卫星s之间的几何关系；x为以坐标参数的形式表示的t历元的卫地距离双差值；b t为以坐标参数的形式表示的t历元的卫地距离三差值；λ为接收的卫星信号的波长；
    Figure PCTCN2022112718-appb-100005
    为共视卫星s的周跳；
    Figure PCTCN2022112718-appb-100006
    为载波相位观测值噪声；
    Figure PCTCN2022112718-appb-100007
    为伪距观测值噪声；
    基于所述s2个载波相位三差观测值和所述s2个伪距三差观测值，利用最小二乘估计对所述2s2个求解式进行求解，解得s2个共视卫星在t历元的各整数周跳值；若任一整数周跳值满足预设条件，则确定所述整数周跳值对应的共视卫星在观测过程中发生周跳。
  2. 如权利要求1所述的方法，其特征在于，解得s2个共视卫星在t历元的各整数周跳值之前，还包括：
    根据监测站在t-1历元接收的s1个共视卫星的载波相位观测值确定所述监测站在所述t-1历元的s1个载波相位三差观测值；根据监测站在t-1历元接收的s1个共视卫星的伪距观测值确定所述监测站在所述t-1历元的s1个伪距三差观测值；所述t-1历元早于所述t历元；所述s1个共视卫星在所述t-1历元不存在周跳现象或周跳已修复；
    将所述t-1历元对应的s1个载波相位三差观测值和s1个伪距三差观测值代入所述GB模型的三差观测方程，得到2s1个求解式；
    基于所述s2个载波相位三差观测值和所述s2个伪距三差观测值，利用最小二乘估计对所述2s2个求解式进行求解，解得s2个共视卫星在t历元的各整数周跳值，包括：
    基于所述s2个载波相位三差观测值、所述s2个伪距三差观测值、所述s1个载波相位 三差观测值和所述s1个伪距三差观测值，利用最小二乘估计对2(s1+s2)个求解式进行求解，解得s2个共视卫星在t历元的各整数周跳值。
  3. 如权利要求2所述的方法，其特征在于，所述监测站为静态监测站；所述b t与b t-1相等；所述b t-1为以坐标参数的形式表示的t-1历元的卫地距离三差值。
  4. 如权利要求1所述的方法，其特征在于，基于所述s2个载波相位三差观测值和所述s2个伪距三差观测值，利用最小二乘估计对所述2s2个求解式进行求解，解得s2个共视卫星在t历元的各整数周跳值，包括：
    将所述2s2个求解式用矩阵表示成如下形式：
    E(l)＝BX；
    其中，l为所述s2个载波相位观测值和所述s2个伪距三差观测值形成的矩阵；E(l)为所述s2个载波相位三差观测值和所述s2个伪距三差观测值的期望；B为
    Figure PCTCN2022112718-appb-100008
    在x、y和z三个方向上的分量及所述s2个共视卫星的卫星信号的波长形成的矩阵；X为所述b t在x、y和z三个方向上的分量及所述s2个共视卫星的各
    Figure PCTCN2022112718-appb-100009
    形成的矩阵；
    通过最小二乘估计的如下公式求得所述X：
    Figure PCTCN2022112718-appb-100010
    其中，所述D3为所述l的方差协方差矩阵；
    根据所述X的矩阵的解算结果，得到所述s2个共视卫星在t历元的各整数周跳值。
  5. 如权利要求4所述的方法，其特征在于，根据所述X的矩阵的解算结果，得到所述s2个共视卫星在t历元的各整数周跳值，包括：
    获取所述X的矩阵的解算结果中的s2个共视卫星在t历元的浮点形式的各周跳值；
    确定各浮点形式的周跳值的方差协方差阵；
    基于各浮点形式的周跳值以及各浮点形式的周跳值的方差协方差阵，将各浮点形式的周跳值转换为整数形式的周跳值，得到所述s2个共视卫星在t历元的各整数周跳值。
  6. 如权利要求5所述的方法，其特征在于，通过如下公式确定各浮点形式的周跳值的方差协方差阵，包括：
    (B TD 3 -1B) -1。
  7. 如权利要求1-6任一项所述的方法，其特征在于，所述载波相位观测值包括多种波长的载波相位观测值；每种波长的载波相位观测值确定各自对应的载波相位三差观测值；
    所述伪距观测值包括多种波长的伪距观测值，每种波长的伪距观测值确定各自对应的伪距三差观测值。
  8. 一种应用于全球导航卫星系统GNSS的周跳探测与修复装置，其特征在于，包括：
    确定单元，用于：
    根据监测站在t历元接收的s2个共视卫星的载波相位观测值确定所述监测站在所述t历元的s2个载波相位三差观测值；根据监测站在t历元接收的s2个共视卫星的伪距观测值确定所述监测站在所述t历元的s2个伪距三差观测值；其中s2≥3；
    将所述t历元对应的s2个载波相位三差观测值和s2个伪距三差观测值代入基于几何GB模型的三差观测方程，得到2s2个求解式；所述GB模型的三差观测方程为：
    Figure PCTCN2022112718-appb-100011
    其中，
    Figure PCTCN2022112718-appb-100012
    u为参考星；s为共视卫星；q为参考站；r为监测站；所述监测站与所述参考站之间的基线长度小于第一预设阈值；获取观测值的历元之间的间隔小于第二预设阈值；
    Figure PCTCN2022112718-appb-100013
    为监测站r在t历元针对共视卫星s的载波相位三差观测值；
    Figure PCTCN2022112718-appb-100014
    为监测站r在t历元针对共视卫星s的伪距三差观测值；A t为t历元的设计矩阵，用以表示t历元监测站r与共视卫星s之间的几何关系；A t-1为t-1历元的设计矩阵，用以表示t-1历元监测站r与共视卫星s之间的几何关系；x为以坐标参数的形式表示的t历元的卫地距离双差值；b t为以坐标参数的形式表示的t历元的卫地距离三差值；λ为接收的卫星信号的波长；
    Figure PCTCN2022112718-appb-100015
    为共视卫星s的周跳；
    Figure PCTCN2022112718-appb-100016
    为载波相位观测值噪声；
    Figure PCTCN2022112718-appb-100017
    为伪距观测值噪声；
    计算单元，用于基于所述s2个载波相位三差观测值和所述s2个伪距三差观测值，利用最小二乘估计对所述2s2个求解式进行求解，解得s2个共视卫星在t历元的各整数周跳值；若任一整数周跳值满足预设条件，则确定所述整数周跳值对应的共视卫星在观测过程中发生周跳。
  9. 一种计算设备，其特征在于，包括：
    存储器，用于存储计算机程序；
    处理器，用于调用所述存储器中存储的计算机程序，按照获得的程序执行权利要求1至7任一项所述的方法。
  10. 一种计算机可读存储介质，其特征在于，所述计算机可读存储介质存储有计算机可执行程序，所述计算机可执行程序用于使计算机执行权利要求1至7任一项所述的方法。
PCT/CN2022/112718 2022-04-13 2022-08-16 一种应用于全球导航卫星系统gnss的周跳探测与修复方法 WO2023197487A1 (zh)

Applications Claiming Priority (2)

Application Number Priority Date Filing Date Title
CN202210387255.2 2022-04-13
CN202210387255.2A CN116953741B (zh) 2022-04-13 2022-04-13 一种应用于全球导航卫星系统gnss的周跳探测与修复方法

Publications (1)

Publication Number Publication Date
WO2023197487A1 true WO2023197487A1 (zh) 2023-10-19

Family

ID=88328804

Family Applications (1)

Application Number Title Priority Date Filing Date
PCT/CN2022/112718 WO2023197487A1 (zh) 2022-04-13 2022-08-16 一种应用于全球导航卫星系统gnss的周跳探测与修复方法

Country Status (2)

Country Link
CN (1) CN116953741B (zh)
WO (1) WO2023197487A1 (zh)

Cited By (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN115267848A (zh) * 2022-06-27 2022-11-01 湘潭大学 一种双频周跳探测及修复方法、系统、设备及介质
CN117452463A (zh) * 2023-12-22 2024-01-26 开普勒卫星科技(武汉)有限公司 一种适用于复杂环境下单频终端的周跳探测与修复方法
CN118244307A (zh) * 2024-03-06 2024-06-25 无锡卡尔曼导航技术有限公司南京技术中心 一种复杂环境下gnss周跳探测自适应阈值确定方法

Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN108169774A (zh) * 2017-12-26 2018-06-15 北方信息控制研究院集团有限公司 支持rtppp和rtk的多模gnss单频周跳探测与修复方法
CN109799521A (zh) * 2019-03-14 2019-05-24 苏州工业园区测绘地理信息有限公司 一种bds/gps三差组合周跳探测与修复方法
CN114152961A (zh) * 2021-11-05 2022-03-08 广东汇天航空航天科技有限公司 一种导航系统的周跳处理方法和装置

Family Cites Families (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US5825326A (en) * 1996-07-09 1998-10-20 Interstate Electronics Corporation Real-time high-accuracy determination of integer ambiguities in a kinematic GPS receiver
US8077089B2 (en) * 2006-08-15 2011-12-13 Rincon Research Corporation Precision geolocation of moving or fixed transmitters using multiple observers
CN109212563B (zh) * 2017-06-29 2023-06-02 同济大学 北斗/gps三频周跳探测与修复方法
CN108196281B (zh) * 2017-11-22 2021-10-29 同济大学 一种基于位置域曲线约束的单频动态周跳探测与修复方法

Patent Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN108169774A (zh) * 2017-12-26 2018-06-15 北方信息控制研究院集团有限公司 支持rtppp和rtk的多模gnss单频周跳探测与修复方法
CN109799521A (zh) * 2019-03-14 2019-05-24 苏州工业园区测绘地理信息有限公司 一种bds/gps三差组合周跳探测与修复方法
CN114152961A (zh) * 2021-11-05 2022-03-08 广东汇天航空航天科技有限公司 一种导航系统的周跳处理方法和装置

Non-Patent Citations (3)

* Cited by examiner, † Cited by third party
Title
LI, BOFENG ET AL.: "Geometry-based cycle slip and data gap repair for multi-GNSS and multi-frequency observations", JOURNAL OF GEODESY, 15 June 2018 (2018-06-15), pages 399 - 417, XP036717253, ISSN: 0949-7714, DOI: 10.1007/s00190-018-1168-5 *
XU, TIANYANG; ZHANG, ZHETAO; HE, XIUFENG; YUAN, HAIJUN: "A Multi-Cycle Slip Detection And Repair Method for Single-Frequency GNSS Data", GEOMATICS AND INFORMATION SCIENCE OF WUHAN UNIVERSITY, CN, vol. 49, no. 3, 1 March 2024 (2024-03-01), CN, pages 465 - 472, XP009549502, DOI: 10.13203/j.whugis20210009 *
YUAN, HAIJUN ET AL.: "Real-Time Cycle Slip Detection and Repair Method for BDS-3 Five-Frequency Data", IEEE ACCESS, vol. 9, 8 April 2021 (2021-04-08), pages 51189 - 51201, XP011848412, ISSN: 2169-3536, DOI: 10.1109/ACCESS.2021.3069745 *

Cited By (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN115267848A (zh) * 2022-06-27 2022-11-01 湘潭大学 一种双频周跳探测及修复方法、系统、设备及介质
CN117452463A (zh) * 2023-12-22 2024-01-26 开普勒卫星科技(武汉)有限公司 一种适用于复杂环境下单频终端的周跳探测与修复方法
CN117452463B (zh) * 2023-12-22 2024-05-14 开普勒卫星科技(武汉)有限公司 一种适用于复杂环境下单频终端的周跳探测与修复方法
CN118244307A (zh) * 2024-03-06 2024-06-25 无锡卡尔曼导航技术有限公司南京技术中心 一种复杂环境下gnss周跳探测自适应阈值确定方法

Also Published As

Publication number Publication date
CN116953741A (zh) 2023-10-27
CN116953741B (zh) 2024-04-02

Similar Documents

Publication Publication Date Title
WO2023197487A1 (zh) 一种应用于全球导航卫星系统gnss的周跳探测与修复方法
EP2156214B1 (en) Partial search carrier-phase integer ambiguity resolution
US9035826B2 (en) Satellite differential positioning receiver using multiple base-rover antennas
CN111965673B (zh) 基于多gnss的单频精密单点定位算法的时间频率传递方法
US7961143B2 (en) Partial search carrier-phase integer ambiguity resolution
CN108802782B (zh) 一种惯导辅助的北斗三频载波相位整周模糊度求解方法
CN113466903B (zh) 一种顾及观测值系统误差的部分模糊度固定算法
WO1998029755A1 (en) Gps multipath detection method and system
CN111273327A (zh) 一种基于组合和非组合混合观测模型的精密单点定位方法
CN112285749B (zh) 全球导航卫星系统原始观测数据处理方法、装置及存储介质
CN112731496B (zh) 一种面向智能终端的gnss精密单点定位数据质量控制方法
CN113466912B (zh) 一种基于多频gnss双天线海上船舶姿态确定方法
CN112444832A (zh) 一种一机多天线接收机的周跳修复方法
US10024974B1 (en) Method and apparatus applicable to positioning in NLOS environment
CN113865592A (zh) 适用于多频gnss精密导航定位中的多路径参数化方法及存储介质
CN115327593B (zh) 基于无人机的定位方法、系统及存储介质
CN116577815A (zh) 一种多频多gnss精密单点定位方法、装置及设备
CN112987048B (zh) 一种适用于Andriod智能终端的高精度定位方法、智能终端以及存储介质
CN110941002B (zh) 一种自适应抗差的序贯最小二乘精密单点定位方法
CN109143273B (zh) 多频多模gnss周跳和数据中断的修复方法
CN116009041B (zh) 一种基于卡方检验的鲁棒自适应gnss高精度定位方法
US20230194731A1 (en) Calculating a position of one device relative to another
CN115980803B (zh) 基于双频码伪距和载波相位观测量进行伪距平滑的方法
CN115826018B (zh) 一种模糊度固定方法、装置、接收机及存储介质
CN118191890A (zh) Gnss非差非组合高精度单点定位方法及系统

Legal Events

Date Code Title Description
121 Ep: the epo has been informed by wipo that ep was designated in this application

Ref document number: 22937125

Country of ref document: EP

Kind code of ref document: A1