Method for generating troposphere delay products combining gnss and polar orbiting remote sensing satellites

Abstract

The application relates to the field of high-precision positioning technology and discloses a troposphere delay product generation method combining GNSS and polar orbit remote sensing satellites, which comprises the following steps: calculating the difference between real-time GNSS true value points and a remote sensing satellite troposphere delay background field to generate a group of differential anchor points; after reliability arbitration is performed on the differential anchor points, a continuous troposphere differential correction field is generated through a terrain-constrained interpolation algorithm, the weight of the interpolation algorithm is based on the effective path distance considering the terrain obstacle penalty, and finally the correction field is superimposed on the background field; through generation and application of a dynamic correction field based on physical deviation, the method replaces the direct mathematical fitting of different source data in the prior art, the interpolation method of terrain constraint is used, space propagation of correction information can follow the restriction of real physical barriers, and unrealistic correction results in complex terrain areas are avoided.

Description

TECHNICAL FIELD

[0001] The application relates to a troposphere delay product generation method combining GNSS and polar orbit remote sensing satellites, and belongs to the technical field of high-precision positioning. BACKGROUND

[0002] To achieve centimeter-level or even more accurate positioning, it is necessary to effectively correct the delay caused by atmospheric refraction in the signal propagation path, and the troposphere zenith total delay is one of the main error sources. In order to obtain continuous troposphere delay products in space and time, a common technical method in the industry is to combine the high-precision point values of GNSS reference stations distributed sparsely on the ground, which are inversed in real time, with the background field obtained by polar orbit remote sensing satellites, which has spatial continuity but a time delay of several hours, and then perform data fusion through mathematical algorithms such as Kriging interpolation or weighted average.

[0003] However, when this technical method is applied to scenarios such as unmanned aerial vehicle (UAV) automated logistics in hilly areas, which have higher requirements for positioning continuity and reliability, its inherent limitations are exposed. The effectiveness of this method depends on a core assumption that the troposphere changes are a spatially continuous and slow stationary field. However, in areas with sharp changes in topography, the troposphere water vapor distribution often presents spatial jumps due to local circulation, and under dynamic weather conditions, the temporal evolution of the troposphere delay is also nonlinear. This theoretical assumption deviates from the physical reality, resulting in insufficient accuracy of the troposphere products generated in the above-mentioned scenarios, which in turn causes the positioning results of the UAV to jump when crossing different micro-terrain regions, frequently triggering its obstacle avoidance or re-planning system.

[0004] Further analysis shows that the above limitations are due to two aspects of the processing logic of the existing technical method. Firstly, the pure mathematical interpolation logic cannot consider real geographical barriers such as mountains or coastlines as physical constraints in the interpolation process, resulting in the incorrect propagation of correction information to physically disconnected areas. This process is inconsistent with the real physical process of troposphere water vapor distribution. Secondly, the existing method treats data sources with different time delays and physical characteristics equally in mathematics, so it cannot effectively remove the systematic bias introduced by time delay in the background field, limiting the reliability of the final product to the update frequency of the background field data. Therefore, how to establish a new technical method that can break away from the assumption of a stationary troposphere field, and in the process of generating high spatial and temporal resolution troposphere products, can fully utilize the spatial coverage advantage of remote sensing satellites, and truly reflect the blocking effect of complex terrain and the real-time dynamic changes of troposphere delay, has become a technical problem to be solved by the present application. SUMMARY

[0005] The application provides a troposphere delay product generation method combining GNSS and polar orbit remote sensing satellites, which mainly aims to solve the problem of the existing static fusion mode with the core of post-weighting, which cannot generate physically reliable and high-precision troposphere delay products in space and time continuity when facing complex terrain and dynamic weather because it depends on the internal assumption of the smooth change of troposphere delay.

[0006] To achieve the above-mentioned purpose, the application provides a troposphere delay product generation method combining GNSS and polar orbit remote sensing satellites, which comprises the following steps:

[0007] Step a: obtaining a remote sensing satellite troposphere delay background field covering a target area and real-time obtained troposphere delay true value points at at least two GNSS sites in the target area;

[0008] Step b: calculating the difference between the troposphere delay true value points and the values of the remote sensing satellite troposphere delay background field at the corresponding positions at each GNSS site to generate a set of differential anchor points;

[0009] Step c: performing reliability arbitration on the set of differential anchor points to identify and correct outlier anchor points;

[0010] Step d: based on the differential anchor points subjected to reliability arbitration, generating a continuous troposphere differential correction field covering the target area by a terrain-constrained interpolation algorithm, wherein the distance for the weight of the terrain-constrained interpolation algorithm in calculating the correction value at any position is an effective path distance, and the effective path distance takes into account the terrain obstacle penalty;

[0011] Step e: superimposing the troposphere differential correction field on the remote sensing satellite troposphere delay background field to generate the final troposphere delay product.

[0012] Preferably, the terrain-constrained interpolation algorithm is an inverse distance weighting method, and the calculation of the effective path distance is based on digital elevation model data, and if a path must pass through an elevation obstacle identified according to the digital elevation model data between two geographical points, the calculation distance of the path is increased by a penalty coefficient defined in the system.

[0013] Preferably, the reliability arbitration in step c comprises: temporarily removing each target differential anchor point to be arbitrated from the set of differential anchor points to form a neighbor anchor point subset, calculating a virtual differential value at the position of the target differential anchor point based on the neighbor anchor point subset by using the terrain-constrained interpolation algorithm, and comparing the real value of the target differential anchor point with the virtual differential value, if the absolute value of the difference between the two values exceeds a statistical threshold determined according to the historical troposphere delay statistical characteristics of the target area, the target differential anchor point is determined as an outlier anchor point, and the virtual differential value is used to replace it.

[0014] Preferably, after step b, there is further a step of physical genesis state discrimination of the differential anchor, the discrimination step comprising: at the GNSS site, synchronously obtaining the carrier phase high frequency variance of its received GNSS signal and the signal signal-to-noise ratio fluctuation variance, constructing a artifact discrimination index wherein, is the signal signal-to-noise ratio fluctuation variance, is the carrier phase high frequency variance, and when performing step d, if the artifact discrimination index exceeds a discrimination threshold, then assigning a suppressive contribution weight to the differential anchor, otherwise, assigning a normal contribution weight to the differential anchor.

[0015] Preferably, there is further an adaptive processing mode switching mechanism based on the spatial density of the differential anchors, the mechanism comprising: before performing step c, performing a spatial density evaluation of the differential anchors within the target region, if the spatial density is lower than a data abundance threshold, then discontinuing the performance of steps c to e, and instead performing the following time gradient feedforward step: based on the historical archived remote sensing satellite troposphere delay background field data, calculating the time evolution gradient of each grid point within the target region, using at least one differential anchor within the target region or its neighboring region to globally calibrate the time evolution gradient field, and based on the calibrated time evolution gradient, extrapolating from the nearest historical remote sensing period to the current time to generate the final troposphere delay product.

[0016] Preferably, there is further a step of generating an accompanying positioning quality indicator field, the step comprising: defining the absolute value of the residual between the true value of each target differential anchor and its corresponding virtual differential value calculated in the reliability arbitration step as the local tropospheric instability index of the position, based on the local tropospheric instability indices of all GNSS site positions, and using a terrain-constrained interpolation algorithm to generate a continuous positioning quality indicator field covering the target region, and outputting the positioning quality indicator field together with the final troposphere delay product.

[0017] Preferably, the elevation obstacle is defined as a continuous region where the terrain slope or elevation change rate exceeds an obstacle threshold.

[0018] Preferably, the statistical threshold is dynamically determined according to the historical troposphere delay statistical characteristics of the target region and the average terrain path distance between GNSS sites within the region.

[0019] Preferably, the differential anchor is assigned a suppressive contribution weight, specifically multiplying its weight in the inverse distance weighted method by a suppressive factor less than 1, the value of the suppressive factor being positively correlated with the size of the artifact discrimination index .

[0020] Preferably, the final tropospheric delay product is a double-layer data package, the data layer one is the final tropospheric delay product, and the data layer two is a positioning quality index field, which is used for risk assessment and decision-making for downstream high-precision positioning applications.

[0021] Compared with the prior art, the present application has the following advantages:

[0022] 1. A set of differential anchor points is generated by calculating the difference between real-time high-precision GNSS true value points and obsolete remote sensing background fields, and a tropospheric differential correction field is generated based on the set of differential anchor points. The method changes the idea of directly mathematically fusing two kinds of data with different temporalities and physical properties in the prior art. Instead of focusing on the absolute values of the two kinds of data, the method focuses on the deviation between them. This deviation captures the sum of the real atmospheric changes since the remote sensing satellite passed by and the inherent deviation of the remote sensing model, thereby converting the problem from numerical fitting to dynamic correction of physical deviation, making the finally generated tropospheric delay product more self-consistent and reliable in physical sense.

[0023] 2. In the present application, a terrain-constrained interpolation algorithm is used to generate the tropospheric differential correction field. The interpolation weight is based on the effective path distance considering terrain obstacle penalty, rather than the simple straight-line distance. This way makes the propagation of correction information respect the physical reality of the discontinuous distribution of tropospheric water vapor caused by real geographical barriers such as mountains or coastlines, avoiding the incorrect crossing of physical barriers in complex terrain areas by traditional interpolation methods, which produces interpolation results that do not conform to the actual situation. Especially in hilly areas, the method can provide more detailed and accurate tropospheric delay correction for high-precision positioning, reducing the problems of flight efficiency reduction and energy consumption increase caused by positioning jump.

[0024] 3. The present application establishes an internal data quality self-consistent verification mechanism without external reference source by performing reliability arbitration on the differential anchor points before generating the correction field. The arbitration uses the idea of leave-one-out cross-validation to construct a virtual true value for each anchor point using its neighbor anchor points. By comparing the difference between the real value and the virtual value, the method identifies and corrects the wild points caused by instantaneous faults of GNSS stations. This self-purification ability avoids local pollution of the correction field caused by false differential anchor points, and improves the fault tolerance of the entire technical process to input data quality fluctuations and the stability of the final product from a methodological point of view.

[0025] 4. The application also provides a physical cause state discrimination mechanism, which can distinguish whether the differential anchor point is caused by real water vapor activity or non-water vapor factors such as ice crystal artifacts that cannot be distinguished by remote sensing detectors, by synchronously analyzing the carrier phase high frequency variance of GNSS signals and the signal signal-to-noise ratio fluctuation variance of two different physical meaning accompanying parameters. This multi-dimensional information collaborative judgment method improves the ability of the system to understand the physical causes behind the remote sensing data, avoids the error correction of remote sensing artifacts as real atmospheric changes, and prevents the introduction of systematic deviations that are inconsistent with physical reality in the final product.

[0026] 5. The application generates an accompanying positioning quality index field in parallel by reusing the correction residual generated in the reliability arbitration process. This design converts the residual information originally used as intermediate calculation waste into a second product that can represent the local tropospheric temporal and spatial stability and has a clear physical meaning. The final output is no longer a single tropospheric delay value, but a decision-level information package containing the value and its confidence evaluation, providing a decision basis for dynamic adjustment of behavior strategies for downstream risk-sensitive applications such as autonomous driving or unmanned aerial vehicles. BRIEF DESCRIPTION OF DRAWINGS

[0027] Fig. 1 It is the overall flowchart of the tropospheric delay product generation method of the application.

[0028] Fig. 2 It is the logical data flowchart of the multi-source data fusion and double-layer product generation of the application.

[0029] Fig. 3 It is the full-chain service architecture diagram of the application from data collection to terminal application. DETAILED DESCRIPTION

[0030] To make the purpose, technical scheme and advantages of the application clearer, the specific embodiments of the application will be described in detail below, but those skilled in the art should understand that the technical scheme of the application can be modified or replaced equivalently without departing from the spirit and scope of the application, and these modifications and equivalent replacements should be covered within the protection scope of the application.

[0031] The application discloses a troposphere delay product generation method combining GNSS and polar orbit remote sensing satellites, and mainly comprises the following stages: data acquisition, differential anchor point generation and arbitration, terrain constraint correction field construction, and final product synthesis and output, wherein data flow starts from acquisition of multi-source data, a set of differential anchor points representing physical deviation is generated and arbitrated, a troposphere differential correction field is constructed by using an interpolation algorithm considering geographical barriers, and the correction field is superimposed on a remote sensing background field to generate a corrected troposphere delay product, in the application of providing high-precision positioning services for unmanned aerial vehicles and other users in complex terrain areas such as hills or coastal areas, a technical challenge is caused by the troposphere delay background field provided by the polar orbit remote sensing satellite which is systematically deviated from the current real atmospheric state due to time delay; in order to solve the challenge, the application is configured to: in a target area, acquire a remote sensing satellite troposphere delay background field covering the area, and simultaneously acquire real troposphere delay points on at least two GNSS sites arranged in the area in real time, at the position of each GNSS site , differential calculation is performed to generate differential anchor points , and the calculation formula is ; the generated differential anchor point set , the value of which represents the sum of real atmospheric changes and inherent deviation of the remote sensing model from the remote sensing satellite passing to the current time at the positions of the sites, considering that a single GNSS site may generate instantaneous outliers due to local interference, thereby affecting the accuracy of the correction field, the system performs a reliability arbitration procedure: for each target differential anchor point to be arbitrated in the set, the target differential anchor point is logically removed from the anchor point set to form a neighbor anchor point subset composed of the remaining anchor points; then, a virtual differential value at the position of the target differential anchor point is calculated based on the neighbor anchor point subset by calling a terrain constraint interpolation algorithm described below, then, the residual absolute value of the real value and the virtual differential value is calculated, and the residual is compared with a statistical threshold value, the statistical threshold value is dynamically determined according to the historical troposphere delay statistical characteristics of the target area and the average terrain path distance between the GNSS sites in the area, so as to adapt to the atmospheric change characteristics and site density in different regions, if the residual absolute value exceeds the dynamic threshold value, the target differential anchor point is determined to be an outlier anchor point, and the corresponding virtual differential value is used to replace it, and then used for subsequent correction field generation, in order to further distinguish whether the differential anchor point is caused by real water vapor change or by ice crystals in high-altitude thin cirrus, the system also performs physical cause state discrimination, and simultaneously acquires carrier phase high-frequency variance and signal-to-noise ratio fluctuation variance of signals at each GNSS site , and construct an artifact discrimination index , and construct an artifact discrimination index The calculated value is compared with a discrimination threshold value, if exceeds the threshold value, it indicates that the differential anchor point at this location is likely to be affected by non-water vapor factors, accordingly, in the subsequent interpolation calculation, an inhibitory contribution weight is assigned to this anchor point, this weight is realized by multiplying its standard weight in the inverse distance weighted method by an inhibition factor less than 1 and positively correlated with the

[0032] After obtaining a set of arbitrated and discriminated differential anchor points, the system generates a continuous tropospheric differential correction field covering the target area through terrain-constrained interpolation algorithm; the interpolation weight used in this algorithm is based on the effective path distance, rather than the simple straight-line distance; the calculation of the effective path distance is based on the digital elevation model DEM data, when calculating the path between any two points, if the path must pass through a recognized elevation obstacle, i.e. a continuous area with a terrain slope or elevation rate exceeding a pre-set obstacle threshold, a system-defined penalty coefficient will be applied to the calculated distance of this path; in this way, the spatial propagation of correction information follows the constraints of physical barriers.

[0033] The specific value of the system adjustment coefficient in the physical cause state discrimination mechanism follows an optimization calibration procedure based on historical data; this procedure first constructs an independent verification data set containing known classification labels, i.e. pure weather samples and artifact samples, from historical GNSS signal data and contemporaneous meteorological data, then introduces a comprehensive index for evaluating classification performance . Fractions, which are defined as Among them, accuracy The proportion of correctly identified artifact samples out of all samples identified as artifacts is denoted as Recall, which is the proportion of correctly identified artifact samples out of all real artifact samples. Then, within a pre-defined parameter scanning interval, different... The value corresponding to The score, and will eventually make The score reaches the maximum value The numerical values ​​are determined as the final engineering parameters for deployment in this area; the data abundance threshold required for the adaptive processing mode switching mechanism and the positive margin in the terrain-constrained interpolation algorithm are also considered. Its calibration process is directly related to the service level requirements and basic data attributes of specific application scenarios.

[0034] Example 1: In a continuously operating automated logistics delivery operation using unmanned aerial vehicles (UAVs), the flight path needs to cross a ridge with a significant elevation difference from a distribution center located in a valley to reach a designated location on the other side. The positioning correction service used for this flight path relies on a conventional tropospheric delay product that combines ground GNSS station data and remote sensing satellite data. When the UAV crosses the ridge peak along the predetermined path, the three-dimensional coordinates calculated by its onboard positioning system experience a momentary jump. This jump triggers the obstacle avoidance and flight path replanning procedures of the flight control system, causing the aircraft to hover and perform energy-intensive attitude adjustments. Ultimately, this results in the delivery mission failing to meet timeliness standards and the return battery power falling below the warning value. To address this situation, the tropospheric delay product generation method of this invention is deployed in the operating area. The system acquires the remote sensing satellite tropospheric delay background field at the same time as the aforementioned scenario. and real-time tropospheric delay truth points of several GNSS stations deployed on the plains on the other side of the valley and ridge. Then, a set of differential anchor points was first generated. In terms of technical processing, this step no longer performs direct mathematical fitting of the absolute numerical fields of the two heterogeneous sources, but instead turns to solving a physical difference field that characterizes the systematic deviation of the background field. Subsequently, the system performs reliability arbitration on the group of difference anchor points, identifies and corrects the erroneous anchor point corresponding to one of the GNSS stations that has generated outliers due to local signal interference.

[0035] In the link of generating the troposphere differential correction field, the differential anchor point generation and the terrain constraint interpolation algorithm are interdependent. The former provides the interpolation process with the correction information to be propagated, which has a clear physical meaning. The latter provides the propagation path for the information, which conforms to the physical law. In the calculation of the correction value at the ridge vertex, the interpolation algorithm considers the effective path distance that is punished by the terrain obstacle as the weight basis, so that the weight contribution of the differential anchor point located in the valley to the ridge vertex is greatly attenuated due to the effective path distance between the two points. The correction value at this position is mainly determined by the differential anchor points on the plain side. This mechanism not only utilizes the spatial coverage of the remote sensing background field, but also realizes the physical correction of local details through the combination of GNSS anchors and terrain constraint interpolation. When the unmanned aerial vehicle performs the same route task again using the troposphere delay product generated by the application, the positioning solution result remains continuous in the entire flight profile. During the flight over the ridge, the positioning coordinates do not jump, and the flight attitude is stable. The final task flight time and energy consumption meet the expected standard operation.

[0036] Example 2: To verify the effectiveness of the technical solution of the application in practical application, especially the improvement effect of the troposphere delay product precision in the complex terrain area, a comparative test is specially designed and performed. The test purpose is to introduce a control group using the conventional interpolation method in the industry, to quantitatively evaluate the precision performance of the application method relative to the conventional method at different terrain feature points under the same input data conditions. The test platform is selected in a coastal mountain area, which has both mountainous areas with sharp changes in elevation and relatively flat coastal plains in the area, and a high-precision GNSS reference station network is deployed to provide continuous observation data. The test data sources include the troposphere delay background field acquired by a polar orbit remote sensing satellite, covering the entire test area The real-time troposphere delay true value points acquired by the GNSS reference station network in the test area , and a 30-meter resolution digital elevation model data covering the area, for constructing the test group and the control group, a part of the GNSS network is selected as the basis data source for differential anchor point generation, and the remaining stations of the network are used as independent precision checking points, the value provided by the checking points is regarded as the reference true value of the position and does not participate in the product generation process of the two methods. The test group adopts the technical solution method A of the application, and the control group adopts a weighted average interpolation scheme method B based on straight-line distance. The two schemes use the same remote sensing background field and differential anchor point basis data source.

[0037] In the test parameter setting, for the terrain constraint interpolation algorithm in the scheme A of the present application, the key parameter terrain obstacle penalty coefficient is set to balance the effective isolation of physical barriers and the smooth adaptation to general terrain undulations. If the parameter value is too small, the correction information cannot cross the main mountain range, and if the parameter value is too large, unnecessary gradient changes may occur on the secondary terrain. The setting rule of the coefficient is to make the path penalty value applied by it greater than the maximum tropospheric delay difference statistically calculated from regional historical data and possibly occurring when crossing the main terrain obstacles. In the coastal mountain scenario of the present experiment, the coefficient is set to a non-limiting example value. After the test is started, the method A and the method B generate the tropospheric delay products covering the entire test area in parallel within a test period of 48 hours, and at each independent accuracy check point, the tropospheric delay values generated by the two methods are compared with the reference true value. At the coastal plain check point CK-01 with relatively flat terrain, the root mean square errors of method A and method B are 4.1 mm and 5.3 mm respectively, which are close to each other. However, at the check point with sharp terrain features, the accuracy difference between the two methods is significant. At the check point CK-04 located on the ridge line, the root mean square error of method B increases to 25.4 mm, while the error of method A is 7.1 mm. Similarly, at the check point CK-03 inside the valley, the error of method B is 18.5 mm, and the error of method A is 6.2 mm. The error difference is due to the fact that the interpolation logic based on straight-line distance used by method B performs mathematical averaging of the differential anchor point information on both sides of the valley and the plain across the terrain, while the terrain constraint interpolation algorithm in method A prevents the propagation of such information across physical barriers due to the design of the effective path distance, thereby obtaining a correction result more consistent with the physical process. The test data show that the technical scheme provided by the present application generates a tropospheric delay product with smaller error in complex terrain areas under the condition of using the same input data source as a conventional scheme. The reduction of the error is due to the compliance of the terrain constraint interpolation algorithm with the physical law of water vapor distribution, thereby avoiding the model distortion problem of the conventional scheme due to the inability to consider the real geographical barriers.

[0038] Embodiment 3: The present embodiment is combined with Figs. 1 to 3 The method for generating the tropospheric delay product of the joint GNSS and polar orbit remote sensing satellite is described as follows. Fig. 1As shown, the flow starts with data acquisition, i.e. obtaining troposphere delay true value points from GNSS site real-time data and obtaining troposphere delay background field from remote sensing satellite data, then computing GNSS true value and remote sensing background field bias through generating differential anchor points step, and performing reliability arbitration on differential anchor points to identify and correct outlier anchor points, the correction residual output by the arbitration step is used to generate accompanying positioning quality index, and the arbitrated anchor points are used to generate correction information through generating troposphere differential correction field step and applying terrain-constrained interpolation algorithm under the data support of digital elevation model DEM, the generated correction field is superimposed to the original background field through product synthesis step, together with the accompanying positioning quality index field, it constitutes the final output double-layer data package, which contains data layer one: final troposphere delay product and data layer two: positioning quality index field.

[0039] As shown, Fig. 2 The external entity: GNSS network provides troposphere delay true value points and the external entity: remote sensing satellite system provides troposphere delay background field, both of which are sent to module 1.0 to generate differential anchor points to output the original differential anchor point set, which is processed by module 2.0 reliability arbitration to generate arbitrated differential anchor points and correction residual, the former is input into module 3.0 to generate differential correction field together with elevation data from data storage D1: digital elevation model DEM, and the latter is input into module 4.0 to generate positioning quality index field together with elevation data, the output of troposphere differential correction field and positioning quality index field from modules 3.0 and 4.0, together with the troposphere delay background field initially provided by the remote sensing satellite system, are finally synthesized into a double-layer data package in module 5.0 to synthesize the final product, and delivered to the external entity: downstream high-precision positioning application.

[0040] As shown, Fig. 3 The node polar orbit remote sensing satellite provides remote sensing background field data through satellite downlink, and the GNSS delay true value points are provided by the real-time data acquisition module in the node GNSS ground station network through the ground network, both of which are collected to the node data and algorithm processing center, which processes the troposphere delay product generation system inside to generate the troposphere delay product, and provides it to the node user terminal such as: unmanned aerial vehicle through data broadcasting service, and the positioning correction service client in the user terminal receives the product for positioning correction.

[0041] Example 4: Before deploying the method of the present application in a new operational area, a set of standardized offline calibration procedures are performed to match the key parameters of the internal algorithms to the tropospheric and terrain environmental characteristics of the area, providing a set of initial values with engineering setting basis for the automated operation of the system; for the statistical threshold value in the reliability arbitration mechanism, the calibration process is as follows: obtain the historical tropospheric delay true value point time series data of all available GNSS stations in the target operational area for at least one year, and the digital elevation model data of the area, for any station in the network, traverse all historical time points, take other stations at that time as the neighbor anchor point subset, and use the terrain constraint interpolation algorithm to calculate a virtual differential value sequence at the target station position, then, subtract the real value sequence of the station from the virtual value sequence to obtain a residual sequence, and calculate the root mean square value of the residual sequence , traverse all stations in the network and repeat the above steps, finally average all the calculated by the stations, and the three times value is set as the statistical threshold value of the area, which is used to determine the outlier anchor point in online operation.

[0042] For the artifact discrimination index in the physical cause state discrimination mechanism , the discrimination threshold value and the inhibitory contribution weight required are determined as follows: collect historical GNSS signal data of the target area, and combine with the same period meteorological satellite cloud image data, screen out two types of sample data sets, one is the sample under pure weather conditions, and the tropospheric delay change is dominated by water vapor activity, and the other is the sample covered by the thin cirrus cloud composed of ice crystals at high altitude, calculate the exponential distribution of the two types of sample sets respectively, the discrimination threshold value is set at the numerical point that can best distinguish the two distributions, for the calculation of the inhibitory contribution weight, when the exponent of a certain differential anchor point exceeds the discrimination threshold value, its weight in interpolation calculation will be multiplied by an inhibitory factor , the factor follows the formula , where is a system adjustment coefficient greater than 0, used to control the rate of weight drop.

[0043] For the terrain obstacle threshold value and the penalty coefficient required by the terrain constraint interpolation algorithm, the calibration procedure is as follows: based on the digital elevation model data of the target area, calculate the slope value distribution histogram of the whole area, and set the slope value at the 95th percentile point of the distribution as the terrain obstacle threshold value, to identify the main physical barriers such as ridge lines in the area, then, for a typical terrain obstacle identified by the threshold value, calculate the average straight line path length of the machine crossing the obstacle , and the average shortest effective path distance around the obstacle to make the algorithm choose detour when routing, the value of the penalty coefficient is set as wherein, is a small positive value margin; by performing the above procedure, all the key parameters in the method of the present application are given quantified values with clear setting basis matching the environmental characteristics of the specific work area before deployment.

[0044] In the embodiment 5, when the method of the present application is applied to a large area with uneven distribution of GNSS sites, in order to calibrate the data abundance threshold relied on by the internal adaptive processing mode switching mechanism, a simulation analysis procedure needs to be performed in advance. The procedure first selects a sub-area with the densest distribution of GNSS sites in the area, and generates a high-precision troposphere delay reference field using all the site data in the sub-area. Then, the system enters a simulation iteration cycle. In each iteration, one anchor point is removed from the set of differential anchor points used to generate the correction field in a predetermined sequence to simulate the gradual reduction of site density.

[0045] After each iteration, the system regenerates a test field using the terrain-constrained interpolation method based on the remaining anchor point set, and calculates the root mean square error of the test field relative to the aforementioned high-precision reference field. Through continuous iteration, the system records the complete data of the root mean square error changing with the spatial density of anchor points, identifies the point where the error growth rate first exceeds the preset change rate threshold as the performance inflection point, and sets the spatial density value corresponding to the inflection point as the data abundance threshold. This procedure directly links the trigger condition of mode switching to a quantifiable output product accuracy change rate; wherein, the determination of the data abundance threshold starts from a clear application scenario tolerable error increment which is a performance indicator preset according to the risk control strategy of the downstream such as autonomous navigation of unmanned aerial vehicles, based on which and the root mean square error of the reference field , a performance warning root mean square error threshold is calculated , and the spatial density of anchor points corresponding to the first time when the root mean square error of the test field reaches is directly determined as the data abundance threshold. At the same time, the positive margin in the terrain obstacle penalty is bound with the grid resolution of the digital elevation model (DEM), and its value is determinedly set as to provide a non-arbitrary decision bias based on the physical properties of the data source for path optimization.

[0046] Embodiment 6: To ensure that the method of the present application can be clearly defined and guaranteed in service level in uninterrupted and high reliability requirements applications, a set of fault tolerance and service degradation operation procedures covering the main data source interruption scenarios need to be configured before the system is formally put into operation. The procedures define the logical path that should be executed when the system monitors the upstream data link interruption, resulting in the complete absence of critical input data flow within the predetermined update period.

[0047] The procedures include the following: when the system monitors the new remote sensing satellite troposphere delay background field When the system monitors the new remote sensing satellite troposphere delay background field When the system monitors the new remote sensing satellite troposphere delay background field When the system monitors the new remote sensing satellite troposphere delay background field

[0048] It is obvious to those skilled in the art that the present application is not limited to the details of the above exemplary embodiments, and can be implemented in other specific forms without departing from the spirit or essential characteristics of the present application.

[0049] Finally, it should be noted that the above embodiments are only used to illustrate the technical solutions of the present application and are not limiting. Although the present application has been described in detail with reference to the preferred embodiments, it should be understood by those skilled in the art that the technical solutions of the present application can be modified or replaced equivalently without departing from the spirit and scope of the technical solutions of the present application.

Claims

1. A method for troposphere delay product generation combining GNSS and polar orbiting remote sensing satellites, characterized in that, The method comprises the following steps: Step a, obtaining a troposphere delay background field of remote sensing satellites covering a target area, and troposphere delay true value points obtained in real time at at least two GNSS sites within the target area; Step b, at the location of each GNSS site, calculating the difference between the troposphere delay true value point and the value of the troposphere delay background field of remote sensing satellites at the corresponding location to generate a set of differential anchor points; Step c, performing reliability arbitration on the set of differential anchor points to identify and correct outlier anchor points; Step d, based on the differential anchor points that have passed reliability arbitration, generating a continuous troposphere differential correction field covering the target area by a terrain-constrained interpolation algorithm, the distance on which the weight of the terrain-constrained interpolation algorithm depends is the effective path distance, which takes into account the terrain obstacle penalty; Step e, superimposing the troposphere differential correction field on the troposphere delay background field of remote sensing satellites to generate the final troposphere delay product; Wherein the terrain-constrained interpolation algorithm is the inverse distance weighted method, and the calculation of the effective path distance is based on digital elevation model data, and between two geographical points, if the path must pass through an elevation obstacle identified according to the digital elevation model data, the calculated distance of the path is increased by a penalty coefficient defined within the system; The elevation obstacle is defined as a continuous area with a terrain slope or an altitude change rate exceeding an obstacle threshold.

2. The method according to claim 1, wherein, The reliability arbitration in step c comprises: temporarily removing each target differential anchor point to be arbitrated from the set of differential anchor points to form a neighbor anchor point subset, calculating a virtual differential value at the location of the target differential anchor point based on the neighbor anchor point subset using the terrain-constrained interpolation algorithm, and comparing the true value of the target differential anchor point with the virtual differential value, if the absolute value of the difference between the two values exceeds a statistical threshold determined according to the statistical characteristics of the historical troposphere delay of the target area, the target differential anchor point is determined to be an outlier anchor point, and the virtual differential value is used to replace it.

3. The method according to claim 1, wherein, comprising a step of physical cause state discrimination of the differential anchor, the discrimination step comprising: at the GNSS site, synchronously obtaining the carrier phase high frequency variance of its received GNSS signal and the signal SNR fluctuation variance, constructing a pseudorange discrimination index wherein, is the signal SNR fluctuation variance, is the carrier phase high frequency variance, and when performing step d, if the pseudorange discrimination index exceeds a discrimination threshold, assigning a suppressive contribution weight to the differential anchor, otherwise, assigning a normal contribution weight to the differential anchor.

4. The method of claim 1, wherein the method further comprises: It also comprises an adaptive processing mode switching mechanism based on the spatial density of differential anchor points, which comprises: before step c is executed, evaluating the spatial density of differential anchor points in the target area, if the spatial density is lower than a data abundance threshold, steps c to e are abandoned, and the following time gradient feedforward step is executed: calculating the time evolution gradient of each grid point in the target area based on the historical archived troposphere delay background field data of remote sensing satellites, using at least one differential anchor point in the target area or its adjacent area to overall calibrate the time evolution gradient field, and based on the calibrated time evolution gradient, deriving from the nearest historical remote sensing period to the current time to generate the final troposphere delay product.

5. The method of claim 2, wherein the method further comprises: The method further comprises a step of generating a concomitant positioning quality index field, which comprises: defining the absolute value of the residual between the true value of each target differential anchor point and its corresponding virtual differential value calculated in the reliability arbitration step as the local tropospheric instability index of the position, generating a continuous positioning quality index field covering the target area based on the local tropospheric instability indexes of all GNSS site positions, and generating a continuous positioning quality index field covering the target area using a terrain-constrained interpolation algorithm, and outputting the positioning quality index field together with the final tropospheric delay product.

6. The method of claim 2, wherein the method further comprises: The statistical threshold is dynamically determined according to the historical tropospheric delay statistical characteristics of the target area and the average terrain path distance between the GNSS sites in the area.

7. The method of claim 3, wherein the method further comprises: An inhibitory contribution weight is assigned to the differential anchor, specifically multiplying its weight in the inverse distance weighting method by an inhibitory factor smaller than 1, the value of the inhibitory factor being positively correlated with the size of the artifact discrimination index .

Images