Time comparison method and device, electronic equipment and storage medium

Abstract

The invention provides a time comparison method and device, electronic equipment and a storage medium, and relates to the technical field of satellite navigation time service, and the method comprises the steps: extracting an observation value of the inclined path total delay of an RDSS signal through a fixed elevation penetration path of the RDSS signal; decomposing the theoretical value of the oblique path total delay of the RDSS signal into a zenith delay variable and a horizontal gradient variable; minimizing the difference between the theoretical value and each observation value to obtain an optimal solution of troposphere vertical path delay and an optimal solution of a horizontal gradient vector in the zenith direction, and constructing a target constraint function of troposphere delay in a double-difference carrier phase observation equation of the RNSS signal; and calculating an estimated value of the double-difference phase ambiguity based on the double-difference carrier phase observation equation and the target constraint function, and calculating the time deviation based on the estimated value of the double-difference phase ambiguity. According to the method, high precision close to that of a TWSTFT method can be kept, and the advantages of low cost and high availability of a GNSS PPP method can be achieved.

Description

TECHNICAL FIELD

[0001] The present application relates to the field of satellite navigation and timing technology, and in particular to a time comparison method and device, an electronic device and a storage medium. BACKGROUND

[0002] High-precision time comparison technology has penetrated into many aspects of national science and technology, economy, military and social life. With the development of science and technology, the requirement for time comparison accuracy is also getting higher and higher.

[0003] Currently, the commonly used high-precision long-distance time comparison methods mainly include the following.

[0004] 1) Satellite two-way time frequency transfer (TWSTFT) method, which sends and receives time signals through the same communication satellite between two ground stations, and eliminates the common path error by calculation, thereby realizing high-precision comparison of clocks at two places. TWSTFT method is one of the internationally recognized highest precision long-distance time transfer technologies, and its time comparison uncertainty is usually better than 1 nanosecond, which is often used as the "gold standard" for verifying other methods. However, the TWSTFT method relies on expensive ground station equipment and dedicated satellite communication resources, and the system construction and operation cost is extremely high, which greatly limits its large-scale popularization and application.

[0005] 2) Global satellite navigation system precise point positioning (GNSS PPP) method, which uses the precise satellite ephemeris and clock error products provided by the global continuous operation reference station network, and the user end single receiver can correct the observation value through complex error model (such as ionospheric delay, tropospheric delay, etc.), to realize high-precision absolute positioning and time recovery. The GNSS PPP method breaks away from the dependence on common-view satellites in traditional common-view method, and has high flexibility, but its final accuracy depends heavily on the accuracy of tropospheric delay modeling, and there is a convergence process of tens of minutes to several hours, which makes it difficult to realize high-precision real-time or near-real-time time transfer.

[0006] Specifically, the GNSS PPP method usually relies on the carrier phase observation value of the Radio Navigation Satellite Service (RNSS) signal. However, this method is extremely affected by the time-varying tropospheric delay, and in the long baseline scenario, this error can reach centimeter level (equivalent to tens of nanoseconds), and the ambiguity resolution efficiency is low and difficult to fix, which becomes the main bottleneck restricting its accuracy and reliability.

[0007] 3) China's global satellite navigation system-Beidou satellite navigation system provides not only RNSS service, but also characteristic radio determination satellite service (RDSS), which can provide communication service and time transfer by the method of two-way pseudorange comparison. However, the accuracy of this method is limited by the pseudorange observation noise (typical value 10 ns), which is difficult to meet the demand of precise time and frequency transfer of nanosecond level.

[0008] In addition, the existing PPP technology widely uses the troposphere model based on experience (such as Saastamoinen, UNB3m, etc.) for delay correction. However, this kind of model does not fully consider the temporal and spatial variation of water vapor content, especially the change of horizontal gradient, resulting in a measured residual greater than 5 cm, and the error will further increase under bad weather conditions, introducing significant systematic bias.

[0009] Therefore, it is urgent to develop a new long-distance time comparison method which can not only maintain high accuracy close to the TWSTFT method, but also has the advantages of low cost and high availability of GNSS PPP method. SUMMARY

[0010] The application provides a time comparison method, device, electronic equipment and storage medium, which can not only maintain high accuracy close to the TWSTFT method, but also has the advantages of low cost and high availability of GNSS PPP method.

[0011] The application provides a time comparison method, comprising: extracting an observation value of the total slant path delay of the satellite radio determination service signal by using the fixed elevation penetration path of the satellite radio determination service signal; decomposing a theoretical value of the total slant path delay of the satellite radio determination service signal into a zenith delay variable and a horizontal gradient variable; the zenith delay variable includes the troposphere vertical path delay in the zenith direction, and the horizontal gradient variable includes a horizontal gradient vector; minimizing the difference between the theoretical value and each observation value to obtain an optimal solution of the troposphere vertical path delay in the zenith direction and an optimal solution of the horizontal gradient vector; using the optimal solution of the troposphere vertical path delay in the zenith direction and the optimal solution of the horizontal gradient vector as a prior value to construct a target constraint function of the troposphere delay in a double-difference carrier phase observation equation of the satellite radio navigation service signal; calculating an estimated value of double-difference phase ambiguity based on the double-difference carrier phase observation equation and the target constraint function, and calculating a time deviation based on the estimated value of the double-difference phase ambiguity.

[0012] According to the time comparison method provided by the application, the total slant path delay of the satellite radio navigation service signal is extracted by using the fixed elevation angle penetration path of the satellite radio determination service signal, and the method comprises the following steps: Obtaining a plurality of pseudo-range observation values of the satellite radio determination service signal; For each pseudo-range observation value, the geometric distance between the corresponding satellite and receiver and the pseudo-range coarse synchronization clock difference of the satellite radio navigation service signal are stripped from the pseudo-range observation value to obtain the observation value of the total slant path delay of the satellite radio determination service signal.

[0013] According to the time comparison method provided by the application, the difference between the theoretical value and each observation value is minimized to obtain the optimal solution of the tropospheric vertical path delay in the zenith direction and the optimal solution of the horizontal gradient vector, and the method comprises the following steps: Constructing a target function for minimizing the square of the Euclidean norm of the difference between each observation value and the theoretical value; Solving the target function by using the least square method to obtain the optimal solution of the tropospheric vertical path delay in the zenith direction and the optimal solution of the horizontal gradient vector.

[0014] According to the time comparison method provided by the application, the optimal solution of the tropospheric vertical path delay in the zenith direction and the optimal solution of the horizontal gradient vector are used as the prior value to construct the target constraint function of the tropospheric delay in the double-difference carrier phase observation equation of the satellite radio navigation service signal, and the method comprises the following steps: Obtaining an observation vector composed of double-difference carrier phase observation values of the satellite radio navigation service signal at a plurality of epochs; Based on the double-difference hardware delay, the double-difference phase ambiguity, the optimal solution of the tropospheric vertical path delay in the zenith direction and the optimal solution of the horizontal gradient vector, a parameter vector to be estimated is constructed; Based on the observation equation, a Jacobian matrix of the observation vector relative to the parameter vector is calculated; Obtaining a tropospheric delay vector composed of tropospheric delays at a plurality of epochs; Obtaining a tropospheric prior value obtained by inversion of the satellite radio determination service signal; Based on the observation vector, the parameter vector to be estimated and the Jacobian matrix, a first regularization term is constructed; Based on the tropospheric delay vector, the tropospheric prior value and the constraint weight, a second regularization term is constructed; Based on the first regularization term and the second regularization term, the target constraint function is constructed.

[0015] According to the time comparison method provided by the application, the first regular term is constructed based on the observation vector, the parameter vector to be estimated and the Jacobian matrix, and the method comprises the following steps: The first regular term is constructed to minimize the square of the Euclidean norm of the difference between the product of the observation vector and the parameter vector to be estimated and the Jacobian matrix.

[0016] According to the time comparison method provided by the application, the second regular term is constructed based on the tropospheric delay vector, the tropospheric prior value and the constraint weight, and the method comprises the following steps: The second regular term is constructed to be the product of the square of the Euclidean norm of the difference between the tropospheric delay vector and the tropospheric prior value and the constraint weight.

[0017] The application further provides a time comparison device, which comprises: The extraction module is used for extracting observation values of slant path total delays of satellite radio determination service signals by using fixed elevation penetration paths of the satellite radio determination service signals. The decomposition module is used for decomposing theoretical values of slant path total delays of the satellite radio determination service signals into zenith delay variables and horizontal gradient variables; the zenith delay variables comprise tropospheric vertical path delays in zenith directions, and the horizontal gradient variables comprise horizontal gradient vectors. The solving module is used for minimizing the difference between the theoretical values and the observation values to obtain optimal solutions of the tropospheric vertical path delays in the zenith directions and optimal solutions of the horizontal gradient vectors. The construction module is used for constructing a target constraint function of tropospheric delays in a double-difference carrier phase observation equation of the satellite radio navigation service signals by using the optimal solutions of the tropospheric vertical path delays in the zenith directions and the optimal solutions of the horizontal gradient vectors as prior values. The calculation module is used for calculating an estimated value of double-difference phase ambiguities based on the double-difference carrier phase observation equation and the target constraint function, and calculating a time bias based on the estimated value of the double-difference phase ambiguities.

[0018] The application further provides an electronic device, which comprises a memory, a processor and a computer program stored in the memory and capable of being executed on the processor, and the processor realizes the time comparison method according to any one of the above-mentioned methods when executing the computer program.

[0019] The application further provides a non-transient computer readable storage medium, which stores a computer program, and the computer program is executed by a processor to realize the time comparison method according to any one of the above-mentioned methods.

[0020] The application further provides a computer program product comprising a computer program which, when executed by a processor, implements the time comparison method according to any one of the above.

[0021] The time comparison method, device, electronic equipment and storage medium provided by the application firstly extracts the observation value of the total slant path delay of the satellite radio determination service signal by using the fixed elevation angle penetration path of the satellite radio determination service signal, can realize independent measurement of the troposphere delay, break through the limitation of the traditional model relying on meteorological parameters, avoid introducing significant systematic bias, and thus improve the accuracy of time comparison; the theoretical value of the total slant path delay of the satellite radio determination service signal is decomposed into a zenith delay variable and a horizontal gradient variable; the zenith delay variable includes the troposphere vertical path delay in the zenith direction, and the horizontal gradient variable includes a horizontal gradient vector; then, the difference between the theoretical value and each observation value is minimized to obtain the optimal solution of the troposphere vertical path delay in the zenith direction and the optimal solution of the horizontal gradient vector; the optimal solution of the troposphere vertical path delay in the zenith direction and the optimal solution of the horizontal gradient vector are used as prior values to construct a troposphere delay target constraint function in a double-difference carrier phase observation equation of the satellite radio navigation service signal, and the estimation value of the double-difference phase ambiguity is calculated based on the double-difference carrier phase observation equation and the target constraint function; first, the modification of the troposphere delay is constrained by the satellite radio determination service signal measurement, which can improve the modification accuracy of the troposphere delay; second, the target constraint function takes into account the horizontal gradient variation, which can reduce the error and thus improve the accuracy of time comparison; third, the use of the target constraint function can reduce the ambiguity search range, thereby improving the success rate of ambiguity fixing, and thus improving the accuracy and reliability of time comparison, and can reduce the dimension of the to-be-estimated parameters, thereby improving the ambiguity solution efficiency; fourth, the time comparison is not directly performed using the satellite radio determination service signal measurement, but a troposphere correction model (i.e., the double-difference carrier phase observation equation and the target constraint function) is established using the measurement value, and parameter estimation is performed through multiple epoch observations, which is equivalent to smoothing processing, thereby suppressing the influence of pseudorange noise, and can meet the demand of precise time and frequency transfer in the order of nanoseconds; fifth, finally, the time bias is calculated based on the estimation value of the double-difference phase ambiguity, which can realize the solution of the sub-nanosecond time bias between receivers. Moreover, the RDSS and RNSS dual-mode receiver has a lower equipment price than the TWSTFT method, which can reduce the system construction and operation cost and is conducive to large-scale popularization and application. Therefore, the time comparison method of the application can not only maintain high precision close to the TWSTFT method, but also has the advantages of low cost and high usability of the GNSS PPP method. BRIEF DESCRIPTION OF DRAWINGS

[0022] In order to make the technical solutions in the present application or prior art clearer, the accompanying drawings needed in the embodiments or prior art description will be briefly introduced. Obviously, the accompanying drawings in the following description are only some embodiments of the present application, and other drawings can be obtained by those of ordinary skill in the art without any creative effort.

[0023] Figure 1 is a flowchart of a time comparison method provided by an embodiment of the present application.

[0024] Figure 2 is a block diagram of a system provided by an embodiment of the present application.

[0025] Figure 3 is a flowchart of a joint solution of a time comparison method provided by an embodiment of the present application.

[0026] Figure 4 is a timing diagram of two sets of clock differences provided by an embodiment of the present application.

[0027] Figure 5 is an Allen variance diagram of two sets of clock differences provided by an embodiment of the present application.

[0028] Figure 6 is a structural diagram of a time comparison device provided by an embodiment of the present application.

[0029] Figure 7 is a structural diagram of an electronic device provided by an embodiment of the present application. DETAILED DESCRIPTION

[0030] In order to make the technical solutions in the present application or prior art clearer, the accompanying drawings needed in the embodiments or prior art description will be briefly introduced. Obviously, the accompanying drawings in the following description are only some embodiments of the present application, and other drawings can be obtained by those of ordinary skill in the art without any creative effort.

[0031] The technical solutions of the present application will be described below with reference to the accompanying drawings. Figures 1-4 The time comparison method of the present application will be described below.

[0032] Please refer to Figure 1 , Figure 1 is a flowchart of a time comparison method provided by an embodiment of the present application. As shown in Figure 1 , the method comprises the following steps 101-105.

[0033] Step 101, an observation value of a slant path total delay of a satellite radio determination service signal is extracted by using a fixed elevation angle penetration path of the satellite radio determination service signal.

[0034] Specifically, by using the physical characteristics of the fixed elevation penetration path of the satellite radio determination service (RDSS) signal, the independent measurement of the slant total delay (STD) of the RDSS signal is realized for the first time, and the observation value of the slant total delay of the RDSS signal is obtained, which can break through the limitation of the traditional model relying on meteorological parameters.

[0035] In an embodiment, step 101 comprises: obtaining a plurality of pseudo-range observation values of the satellite radio determination service signal; for each pseudo-range observation value, stripping the geometric distance between the corresponding satellite and the receiver and the pseudo-range coarse synchronization clock difference of the satellite radio navigation service signal from the pseudo-range observation value to obtain the observation value of the slant total delay of the satellite radio determination service signal.

[0036] Specifically, the observation value of the slant total delay of the RDSS signal can be extracted by the following expression (1): (1) Wherein, represents the observation value of the slant total delay of the RDSS signal, represents the pseudo-range observation value of the RDSS signal, represents the geometric distance between the corresponding satellite and the receiver (calculated by the RDSS positioning coordinates), represents the speed of light, represents the pseudo-range coarse synchronization clock difference of the satellite radio navigation service (RNSS) signal (accuracy 5ns).

[0037] Step 102, decompose the theoretical value of the slant total delay of the satellite radio determination service signal into zenith delay variable and horizontal gradient variable; the zenith delay variable includes the zenith direction troposphere vertical path delay, and the horizontal gradient variable includes the horizontal gradient vector.

[0038] Specifically, the decomposition can be performed by the following expression (1): (2) Wherein, the left side of the plus sign is the zenith component variable, and the right side of the plus sign is the horizontal gradient component variable, represents the theoretical value of the slant total delay of the RDSS signal, represents the zenith direction troposphere vertical path delay, and represents the components of the horizontal gradient in the north and east directions, represents the satellite azimuth, Satellite fixed elevation angle (value range: 30°-40°, typical value in China is about 35°).

[0039] Step 103, minimizing the difference between the theoretical value and each observation value to obtain the optimal solution of the tropospheric vertical path delay in zenith direction and the optimal solution of the horizontal gradient vector.

[0040] In an embodiment, step 103 comprises: constructing an objective function that minimizes the square of the Euclidean norm of the difference between each observation value and the theoretical value; and solving the objective function by using the least square method to obtain the optimal solution of the tropospheric vertical path delay in zenith direction and the optimal solution of the horizontal gradient vector.

[0041] Specifically, the objective function is constructed by the following expression (3): (3) wherein, represents the horizontal gradient vector, .

[0042] The optimal solution of and the optimal solution of are obtained by solving expression (3) by using the least square method.

[0043] Step 104, using the optimal solution of the tropospheric vertical path delay in zenith direction and the optimal solution of the horizontal gradient vector as the prior value, constructing the tropospheric delay objective constraint function in the double-difference carrier phase observation equation of the satellite radio navigation service signal.

[0044] Specifically, the double-difference carrier phase observation equation of the RNSS signal is shown in expression (4): (4) wherein, represents the physical parameter after the double-difference (inter-satellite and inter-station difference) processing, represents the wavelength of the RNSS signal, represents the carrier phase observation value of the RNSS signal, represents the geometric distance between the satellite and the antenna phase center of the station, represents the hardware delay, represents the phase ambiguity, represents the tropospheric delay, which is the part that can be constrained by the prior value, represents the measurement noise.

[0045] In an embodiment, step 104 comprises sub-steps 1041-1048.

[0046] Step 1041, obtaining the double-difference carrier phase observation value of the satellite radio navigation service signal at multiple epochs an observation vector consisting of the ionosphere-free carrier phase observations (in meters) of the RNSS signals , a double-difference carrier phase observation (in cycles, i.e. number of periods of the electromagnetic wave) of the RNSS signals.

[0047] Step 1042, based on the double-difference hardware delays , double-difference phase ambiguities , optimal solution of the tropospheric vertical path delay in zenith direction and the horizontal gradient vector constructing the parameter vector to be estimated , .

[0048] Step 1043, based on the observation equation calculating the observation vector with respect to the parameter vector the Jacobian matrix (i.e. the matrix of partial derivatives).

[0049] Step 1044, obtaining the tropospheric delay vector consisting of the tropospheric delays at a plurality of epochs .

[0050] Step 1045, obtaining the tropospheric a priori values from the satellite radio determination service signals .

[0051] Step 1046, constructing a first regularizer based on the observation vector , the parameter vector to be estimated and the Jacobian matrix .

[0052] Optionally, step 1046 comprises constructing a first regularizer that minimizes the square of the Euclidean norm of the difference between the observation vector and the product of the parameter vector to be estimated and the Jacobian matrix . .

[0053] Step 1047, constructing a second regularizer based on the tropospheric delay vector , the tropospheric a priori values and a constraint weight .

[0054] Optionally, step 1047 comprises constructing a second regularizer that is the product of the square of the Euclidean norm of the difference between the tropospheric delay vector and the tropospheric a priori values and the constraint weight .​ .

[0055] Step 1048, constructing the target constraint function based on the first regularization term and the second regularization term.

[0056] Optionally, the target constraint function is constructed by the following expression (5): (5) wherein, is the Jacobian matrix relative to calculated according to the expression (4), and optionally, = 0.7 as an empirical value of the constraint weight.

[0057] Step 105, calculating the estimated value of the double-difference phase ambiguity based on the double-difference carrier phase observation equation and the target constraint function, and calculating the time bias based on the estimated value of the double-difference phase ambiguity.

[0058] Specifically, the estimated value of the double-difference phase ambiguity can be calculated based on the double-difference carrier phase observation equation and the target constraint function by using the classic LAMBDA method .

[0059] The tropospheric prior value obtained by RDSS inversion is embedded into the double-difference carrier phase observation equation as a prior constraint, which can reduce the dimension of the estimated parameters by 50% and improve the ambiguity resolution efficiency.

[0060] The time bias is calculated by the following expression (6): (6) wherein, represents the time bias.

[0061] In a specific implementation, the time comparison precision calculated by the time comparison method of the embodiment reaches 0.98 ns (3000 km baseline), which meets the sub-nanosecond level requirement of quantum communication, gravitational wave detection, etc.

[0062] The time comparison method provided by the embodiment of the application firstly extracts the observation value of the total slant path delay of the satellite radio determination service signal by using the fixed elevation angle penetration path of the satellite radio determination service signal, can realize independent measurement of the troposphere delay, breaks through the limitation of the traditional model relying on meteorological parameters, avoids introducing significant systematic bias, and thus improves the accuracy of time comparison; decomposes the theoretical value of the total slant path delay of the satellite radio determination service signal into a zenith delay variable and a horizontal gradient variable; the zenith delay variable includes the troposphere vertical path delay in the zenith direction, and the horizontal gradient variable includes a horizontal gradient vector; then, the difference between the theoretical value and each observation value is minimized to obtain the optimal solution of the troposphere vertical path delay in the zenith direction and the optimal solution of the horizontal gradient vector; the optimal solution of the troposphere vertical path delay in the zenith direction and the optimal solution of the horizontal gradient vector are used as prior values to construct a troposphere delay target constraint function in a double-difference carrier phase observation equation of the satellite radio navigation service signal, and the estimated value of the double-difference phase ambiguity is calculated based on the double-difference carrier phase observation equation and the target constraint function. First, the modification of the troposphere delay is constrained by the satellite radio determination service signal measurement, which can improve the modification accuracy of the troposphere delay. Second, the target constraint function takes into account the horizontal gradient variation, which can reduce errors and thus improve the accuracy of time comparison. Third, the use of the target constraint function can reduce the ambiguity search range, thereby improving the success rate of ambiguity fixing and further improving the accuracy and reliability of time comparison, and can reduce the dimension of the estimated parameters, thereby improving the efficiency of ambiguity solution. Fourth, the time comparison is not directly performed using the satellite radio determination service signal measurement, but a troposphere correction model (i.e., the double-difference carrier phase observation equation and the target constraint function) is established using the measurement value, and parameter estimation is performed through multiple epoch observations, which is equivalent to smoothing processing, thereby suppressing the influence of pseudorange noise and meeting the demand for precise time and frequency transfer of the order of nanoseconds. Fifth, finally, the time bias is calculated based on the estimated value of the double-difference phase ambiguity, which can realize the solution of the sub-nanosecond time bias between receivers. Moreover, the RDSS and RNSS dual-mode receiver has a lower equipment price than the TWSTFT method, which can reduce the system construction and operation cost and is conducive to large-scale popularization and application. Therefore, the time comparison method of the embodiment of the application can not only maintain high precision close to the TWSTFT method, but also has the advantages of low cost and high availability of the GNSS PPP method.

[0063] The following will be combined Figures 2-4 with the data of a certain region in a certain time range in 2021 to verify the accuracy (including accuracy and stability) of the time comparison method provided by the embodiment of the application.

[0064] As shown in Figure 2 , the time comparison system comprises an STD extraction module, a joint solver and a time bias output interface.

[0065] For the STD extraction module, as shown in Figure 2 and 3 , first, input the multiple pseudorange observations of the RDSS signal . Specifically, the RDSS signal of the Geostationary Earth Orbit (GEO) satellite (PRN C01, Beidou C01 satellite) is received by the RDSS receiver, and the signal is transmitted to multiple satellites and then retransmitted by the satellites to the master station for measurement. Therefore, a set of measurements can have multiple satellite measurements at the same time, i.e., multiple pseudorange observations of the RDSS signal . Then, based on the RDSS signal, the elevation data sent by the barometer / digital elevation model (DEM), and the satellite ephemeris, the user position (i.e., the RDSS positioning coordinates) is calculated, and the geometric distance between the satellite and the receiver is calculated based on the RDSS positioning coordinates . Next, input the pseudorange coarse synchronization clock difference of the RNSS signal . Finally, the observation value of the slant path total delay of the RDSS signal is extracted by expression (1) .

[0066] For the joint solver, as shown in Figure 3 , after the observation value of the slant path total delay of the RDSS signal is extracted , the search space is set: , , the least squares method is used to solve expression (3) to solve the tropospheric parameters: 、 、 . The tropospheric prior value is injected, i.e., the prior constraint. The double-frequency RNSS carrier phase is collected, and the double-difference carrier phase observation value of the RNSS signal is input , and the double-difference carrier phase observation equation of the RNSS signal is constructed as shown in expression (4). The target constraint function as shown in expression (5) is solved to obtain the estimated value of the double-difference phase ambiguity , and the time bias is calculated by expression (6) .

[0067] For the time bias output interface, the time bias is output.

[0068] In the case of using the RDSS signal for assistance, the success rate of fixing the double-difference phase ambiguity is increased from 89% to 93%, and the convergence time is shortened from 34 minutes to 26 minutes.

[0069] In a specific implementation, the clock difference result of the two-site integrated inter-satellite link and the C-band inter-station time synchronization is used as a reference benchmark for verification. Figure 4 is a timing diagram of two sets of clock differences provided by the embodiment of the application. As shown in Figure 4 , the blue color is the clock difference result calculated by the time comparison method of the embodiment of the application, the red color is the clock difference result of the C-band inter-station time synchronization, the horizontal coordinate is the day of year in 2021, and the vertical coordinate is the clock difference result (diff Clk). From Figure 4 , it can be seen that the clock difference result (diff Clk) calculated by the time comparison method of the embodiment of the application is consistent with the clock difference of the reference benchmark, and the root mean square error of the clock difference is 0.96 ns, which can verify the accuracy of the time comparison method of the embodiment of the application. Figure 5 is a diagram of Allan deviation of two sets of clock differences provided by the embodiment of the application. As shown in Figure 5 , the blue color is the clock difference result calculated by the time comparison method of the embodiment of the application, the red color is the clock difference result of the C-band inter-station time synchronization, the horizontal coordinate is the averaging time, and the vertical coordinate is the Allan deviation. From Figure 5 , it can be seen that the clock difference result calculated by the time comparison method of the embodiment of the application has better frequency stability (Frequency Stability) than the clock difference result of the C-band inter-station time synchronization, especially in the short-term stability aspect.

[0070] The time comparison device provided by the application is described below. The time comparison device described below can be correspondingly referred to the time comparison method described above.

[0071] Please refer to Figure 6 , Figure 6 is a structural diagram of the time comparison device provided by the embodiment of the application. As shown in Figure 6 , the device can include: an extraction module 10 configured to extract an observation value of a slant path total delay of a satellite radio determination service signal by using a fixed elevation angle penetration path of the satellite radio determination service signal; a decomposition module 20 configured to decompose a theoretical value of the slant path total delay of the satellite radio determination service signal into a zenith delay variable and a horizontal gradient variable; the zenith delay variable includes a zenith direction troposphere vertical path delay, and the horizontal gradient variable includes a horizontal gradient vector; a solving module 30 configured to minimize the difference between the theoretical value and each of the observation values to obtain an optimal solution of the zenith direction troposphere vertical path delay and an optimal solution of the horizontal gradient vector; The constructing module 40 is configured to construct a target constraint function of tropospheric delay in a double-difference carrier phase observation equation of the satellite radio navigation service signal by using the optimal solution of the tropospheric vertical path delay in the zenith direction and the optimal solution of the horizontal gradient vector as prior values. The calculating module 50 is configured to calculate an estimated value of a double-difference phase ambiguity based on the double-difference carrier phase observation equation and the target constraint function, and calculate a time bias based on the estimated value of the double-difference phase ambiguity.

[0072] In an embodiment, the extracting module 10 is specifically configured to: Obtain a plurality of pseudo-range observation values of the satellite radio determination service signal. For each pseudo-range observation value, strip the geometric distance between a corresponding satellite and a receiver and a pseudo-range coarse synchronization clock bias of the satellite radio navigation service signal from the pseudo-range observation value to obtain an observation value of the slant path total delay of the satellite radio determination service signal.

[0073] In an embodiment, the solving module 30 is specifically configured to: Construct a target function of minimizing the square of the Euclidean norm of the difference between each observation value and the theoretical value; Solve the target function by using a least square method to obtain the optimal solution of the tropospheric vertical path delay in the zenith direction and the optimal solution of the horizontal gradient vector.

[0074] In an embodiment, the constructing module 40 comprises: A first obtaining unit configured to obtain an observation vector composed of double-difference carrier phase observation values of the satellite radio navigation service signal at a plurality of epochs; A first constructing unit configured to construct a parameter vector to be estimated based on a double-difference hardware delay, a double-difference phase ambiguity, the optimal solution of the tropospheric vertical path delay in the zenith direction and the optimal solution of the horizontal gradient vector; A calculating unit configured to calculate a Jacobian matrix of the observation vector with respect to the parameter vector based on an observation equation; A second obtaining unit configured to obtain a tropospheric delay vector composed of tropospheric delays at a plurality of epochs; A third obtaining unit configured to obtain a tropospheric prior value obtained by inversion of the satellite radio determination service signal; A second constructing unit configured to construct a first regularization term based on the observation vector, the parameter vector to be estimated and the Jacobian matrix; A third constructing unit configured to construct a second regularization term based on the tropospheric delay vector, the tropospheric prior value and a constraint weight; A fourth constructing unit is configured to construct the target constraint function based on the first regular term and the second regular term.

[0075] In an embodiment, the second constructing unit is specifically configured to: construct the first regular term that minimizes a square of a Euclidean norm of a difference between the observation vector and a product of the parameter vector to be estimated and the Jacobian matrix.

[0076] In an embodiment, the third constructing unit is specifically configured to: construct the second regular term that is a product of a square of a Euclidean norm of a difference between the tropospheric delay vector and the tropospheric prior value and the constraint weight.

[0077] Figure 7 An example of a schematic diagram of a physical structure of an electronic device is shown in Figure 7 As shown, the electronic device can include a processor 710, a communications interface 720, a memory 730, and a communications bus 740, wherein the processor 710, the communications interface 720, and the memory 730 can communicate with each other through the communications bus 740. The processor 710 can invoke a logical instruction in the memory 730 to execute a time comparison method, which includes: extracting an observation value of a slant path total delay of a satellite radio navigation service signal using a fixed elevation penetration path of the satellite radio navigation service signal; decomposing a theoretical value of the slant path total delay of the satellite radio navigation service signal into a zenith delay variable and a horizontal gradient variable; the zenith delay variable includes a tropospheric vertical path delay in a zenith direction, and the horizontal gradient variable includes a horizontal gradient vector; minimizing a difference between the theoretical value and each observation value to obtain an optimal solution of the tropospheric vertical path delay in the zenith direction and an optimal solution of the horizontal gradient vector; using the optimal solution of the tropospheric vertical path delay in the zenith direction and the optimal solution of the horizontal gradient vector as a prior value to construct a target constraint function of a tropospheric delay in a double-difference carrier phase observation equation of a satellite radio navigation service signal; calculating an estimated value of a double-difference phase ambiguity based on the double-difference carrier phase observation equation and the target constraint function, and calculating a time bias based on the estimated value of the double-difference phase ambiguity.

[0078] Moreover, the logic instructions in the memory 730 described above can be implemented in the form of software functional units and sold or used as independent products, and can be stored in a computer readable storage medium. Based on such understanding, the technical solutions of the present application essentially or the part that contributes to the prior art or part of the technical solutions can be embodied in the form of a software product. The computer software product is stored in a storage medium, and includes several instructions for making a computer device (which can be a personal computer, a server, or a network device, etc.) execute all or part of the steps of the methods described in the various embodiments of the present application. The aforementioned storage medium includes: a U disk, a mobile hard disk, a read-only memory (ROM, Read-Only Memory), a random access memory (RAM, Random Access Memory), a magnetic disk or an optical disk, and various media that can store program codes.

[0079] In another aspect, the present application also provides a computer program product, which comprises a computer program, the computer program can be stored on a non-transitory computer readable storage medium, and the computer program can be executed by a processor to enable a computer to execute the time comparison method provided by the above-mentioned methods. The method comprises: extracting an observation value of a slant path total delay of a satellite radio determination service signal by using a fixed elevation angle penetration path of the satellite radio determination service signal; decomposing a theoretical value of the slant path total delay of the satellite radio determination service signal into a zenith delay variable and a horizontal gradient variable; the zenith delay variable comprises a zenith direction troposphere vertical path delay, and the horizontal gradient variable comprises a horizontal gradient vector; minimizing the difference between the theoretical value and each observation value to obtain an optimal solution of the zenith direction troposphere vertical path delay and an optimal solution of the horizontal gradient vector; using the optimal solution of the zenith direction troposphere vertical path delay and the optimal solution of the horizontal gradient vector as a prior value to construct a target constraint function of a troposphere delay in a double-difference carrier phase observation equation of a satellite radio navigation service signal; calculating an estimated value of a double-difference phase ambiguity based on the double-difference carrier phase observation equation and the target constraint function, and calculating a time deviation based on the estimated value of the double-difference phase ambiguity.

[0080] In yet another aspect, the present application also provides a non-transitory computer readable storage medium having stored thereon a computer program, which, when executed by a processor, implements the time comparison method provided by any of the above methods, and the method comprises: extracting an observation value of a slant path total delay of a satellite radio determination service signal by using a fixed elevation angle penetration path of the satellite radio determination service signal; decomposing a theoretical value of the slant path total delay of the satellite radio determination service signal into a zenith delay variable and a horizontal gradient variable; the zenith delay variable comprises a zenith direction troposphere vertical path delay, and the horizontal gradient variable comprises a horizontal gradient vector; minimizing a difference between the theoretical value and each observation value to obtain an optimal solution of the zenith direction troposphere vertical path delay and an optimal solution of the horizontal gradient vector; using the optimal solution of the zenith direction troposphere vertical path delay and the optimal solution of the horizontal gradient vector as a priori value to construct a target constraint function of a troposphere delay in a double-difference carrier phase observation equation of a satellite radio navigation service signal; calculating an estimated value of a double-difference phase ambiguity based on the double-difference carrier phase observation equation and the target constraint function, and calculating a time bias based on the estimated value of the double-difference phase ambiguity.

[0081] The device embodiments described above are merely illustrative, wherein the units described as separate components can or can not be physically separate, and the components shown as units can or can not be physical units, i.e., can be located in one place or distributed on multiple network units. Part or all of the modules can be selected to achieve the purpose of the embodiment according to actual needs. Those skilled in the art can understand and implement without creative labor.

[0082] From the above description of the embodiments, those skilled in the art can clearly understand that each embodiment can be realized by means of software and a necessary general hardware platform, and of course can also be realized by hardware. Based on such understanding, the above technical solutions, essentially or in terms of contribution to the prior art, can be embodied in the form of a software product, which can be stored in a computer readable storage medium, such as a ROM / RAM, a magnetic disk, an optical disk, etc., and includes a number of instructions to make a computer device (which can be a personal computer, a server, or a network device, etc.) execute the methods described in each embodiment or some parts of the embodiments.

[0083] Finally, it should be noted that: the above embodiments are only used to illustrate the technical solutions of the present application, and not to limit them; although the present application has been described in detail with reference to the foregoing embodiments, those skilled in the art should understand that they can still modify the technical solutions recorded in the foregoing embodiments, or make equivalent replacements to some technical features; and these modifications or replacements do not make the corresponding technical solutions deviate from the spirit and scope of the technical solutions of the embodiments of the present application.

Claims

1. A time comparison method, characterized by, The method comprises the following steps: extracting an observation value of total slant path delay of the satellite radio navigation service signal by using a fixed elevation penetration path of the satellite radio determination service signal; decomposing a theoretical value of total slant path delay of the satellite radio navigation service signal into a zenith delay variable and a horizontal gradient variable; the zenith delay variable comprises tropospheric vertical path delay in zenith direction, and the horizontal gradient variable comprises a horizontal gradient vector; minimizing a difference between the theoretical value and each observation value to obtain an optimal solution of the tropospheric vertical path delay in the zenith direction and an optimal solution of the horizontal gradient vector; using the optimal solution of the tropospheric vertical path delay in the zenith direction and the optimal solution of the horizontal gradient vector as a prior value to construct a target constraint function of tropospheric delay in a double-difference carrier phase observation equation of the satellite radio navigation service signal; calculating an estimated value of double-difference phase ambiguity based on the double-difference carrier phase observation equation and the target constraint function, and calculating a time bias based on the estimated value of the double-difference phase ambiguity.

2. The time comparison method of claim 1, wherein, The method for extracting the observation value of total slant path delay of the satellite radio determination service signal by using the fixed elevation penetration path of the satellite radio determination service signal comprises the following steps: obtaining a plurality of pseudo-range observation values of the satellite radio determination service signal; for each pseudo-range observation value, stripping the corresponding geometric distance between the satellite and the receiver and the pseudo-range coarse synchronization clock difference of the satellite radio navigation service signal from the pseudo-range observation value to obtain the observation value of total slant path delay of the satellite radio determination service signal.

3. The time comparison method according to claim 1 or 2, characterized in that, The method for minimizing the difference between the theoretical value and each observation value to obtain the optimal solution of the tropospheric vertical path delay in the zenith direction and the optimal solution of the horizontal gradient vector comprises the following steps: constructing a target function for minimizing the square of the Euclidean norm of the difference between each observation value and the theoretical value; solving the target function by using a least square method to obtain the optimal solution of the tropospheric vertical path delay in the zenith direction and the optimal solution of the horizontal gradient vector.

4. The time comparison method according to claim 1 or 2, characterized in that, The method for using the optimal solution of the tropospheric vertical path delay in the zenith direction and the optimal solution of the horizontal gradient vector as a prior value to construct a target constraint function of tropospheric delay in a double-difference carrier phase observation equation of the satellite radio navigation service signal comprises the following steps: obtaining an observation vector composed of double-difference carrier phase observation values of the satellite radio navigation service signal at a plurality of epochs; constructing a parameter vector to be estimated based on double-difference hardware delay, double-difference phase ambiguity, the optimal solution of the tropospheric vertical path delay in the zenith direction and the optimal solution of the horizontal gradient vector; calculating a Jacobian matrix of the observation vector relative to the parameter vector based on an observation equation; obtaining a tropospheric delay vector composed of tropospheric delays at a plurality of epochs; obtaining a tropospheric prior value obtained by inversion of the satellite radio determination service signal; constructing a first regularization term based on the observation vector, the parameter vector to be estimated and the Jacobian matrix; constructing a second regularization term based on the tropospheric delay vector, the tropospheric prior value and a constraint weight; constructing the target constraint function based on the first regular term and the second regular term.

5. The time comparison method of claim 4, wherein, The constructing the first regular term based on the observation vector, the parameter vector to be estimated and the Jacobian matrix comprises: constructing the first regular term that minimizes a square of a Euclidean norm of a difference between the observation vector and a product of the parameter vector to be estimated and the Jacobian matrix.

6. The time comparison method of claim 4, wherein, The constructing the second regular term based on the troposphere delay vector, the troposphere prior value and a constraint weight comprises: constructing the second regular term that is a product of a square of a Euclidean norm of a difference between the troposphere delay vector and the troposphere prior value and the constraint weight.

7. A time comparison device, characterized by comprise: an extracting module configured to extract observation values of slant path total delays of satellite radio determination service signals by using fixed elevation angle penetration paths of the satellite radio determination service signals; a decomposing module configured to decompose theoretical values of the slant path total delays of the satellite radio determination service signals into zenith delay variables and horizontal gradient variables; the zenith delay variables comprise troposphere vertical path delays in zenith directions, and the horizontal gradient variables comprise horizontal gradient vectors; a solving module configured to minimize a difference between the theoretical values and the observation values to obtain optimal solutions of the troposphere vertical path delays in the zenith directions and optimal solutions of the horizontal gradient vectors; a constructing module configured to construct a target constraint function of troposphere delays in double difference carrier phase observation equations of the satellite radio navigation service signals by using the optimal solutions of the troposphere vertical path delays in the zenith directions and the optimal solutions of the horizontal gradient vectors as prior values; a calculating module configured to calculate an estimated value of double difference phase ambiguities based on the double difference carrier phase observation equations and the target constraint function, and calculate a time bias based on the estimated value of the double difference phase ambiguities.

8. An electronic device comprising a memory, a processor, and a computer program stored on the memory and running on the processor, characterized in that, The processor implements the time comparison method according to any one of claims 1 to 6 when executing the computer program. 9.A non-transitory computer-readable storage medium having stored thereon a computer program, characterized in that, The computer program, when executed by a processor, implements the time comparison method according to any one of claims 1 to 6.

10. A computer program product comprising a computer program, characterized in that, The computer program, when executed by a processor, implements the time comparison method according to any one of claims 1 to 6.

Images