A large-height-difference GNSS deformation monitoring method and system

Abstract

The application provides a large-elevation-difference GNSS deformation monitoring method and system, and relates to the technical field of engineering surveying, which comprises obtaining first information and second information, the first information comprising meteorological parameter information of GNSS reference stations and monitoring stations, and the second information comprising observation information of the GNSS reference stations and the monitoring stations, the observation information comprising pseudo-range and phase observation values; the first information is sent to a preset meteorological parameter model for calculation to obtain third information, the third information comprising tropospheric zenith total delay of the GNSS reference stations and tropospheric zenith total delay of the monitoring stations; a double-difference observation equation is constructed based on the second information and the third information to obtain a GNSS relative positioning observation equation; the double-difference observation equation is solved, and deformation estimation data of GNSS monitoring points are determined based on a Kalman filtering algorithm and a calculation result of the double-difference observation equation. The application realizes high-precision deformation monitoring calculation under a large-elevation-difference condition.

Description

TECHNICAL FIELD

[0001] The present application relates to the field of engineering surveying, in particular to a GNSS deformation monitoring method and system for large height difference. BACKGROUND

[0002] The GNSS deformation monitoring adopts a relative positioning mode, that is, a GNSS monitoring station is arranged at a monitoring point, a GNSS reference station is arranged at a nearby stable position, and a GNSS relative positioning method is used to solve the relative position of the GNSS monitoring station relative to the GNSS reference station. However, when the height difference between stations is large, the tropospheric delay error between stations will be significantly large. Therefore, a method and system capable of calculating the prior tropospheric delay error by using the measured meteorological parameters at all times, and taking the residual tropospheric delay error as a parameter to estimate in the solving process, are needed to realize high-precision deformation monitoring solution under large height difference. SUMMARY

[0003] The present application aims to provide a GNSS deformation monitoring method and system for large height difference to improve the above problems. In order to achieve the above purpose, the technical scheme adopted by the present application is as follows:

[0004] In one aspect, the present application provides a GNSS deformation monitoring method for large height difference, comprising:

[0005] obtaining first information and second information, the first information comprising meteorological parameter information of the GNSS reference station and the monitoring station, the meteorological parameter information comprising temperature parameter information, pressure parameter information and water vapor pressure parameter information, the second information comprising observation information of the GNSS reference station and the monitoring station, the observation information comprising pseudorange and phase observation value;

[0006] sending the first information to a preset meteorological parameter model for calculation respectively to obtain third information, the third information comprising tropospheric zenith total delay of the GNSS reference station and tropospheric zenith total delay of the monitoring station;

[0007] constructing a double-difference observation equation based on the second information and the third information to obtain a GNSS relative positioning observation equation;

[0008] solving the double-difference observation equation, and determining deformation estimation data of the GNSS monitoring point based on the calculation result of the Kalman filtering algorithm and the double-difference observation equation.

[0009] In another aspect, the present application also provides a GNSS deformation monitoring system for large height difference, comprising:

[0010] The acquisition unit is used for acquiring first information and second information, the first information comprises meteorological parameter information of GNSS reference stations and monitoring stations, the meteorological parameter information comprises temperature parameter information, pressure parameter information and water vapor pressure parameter information, and the second information comprises observation information of the GNSS reference stations and the monitoring stations, the observation information comprises pseudo-range and phase observation value;

[0011] The first calculation unit is used for transmitting the first information to a preset meteorological parameter model respectively for calculation, obtaining third information, and the third information comprises tropospheric zenith total delay of the GNSS reference stations and tropospheric zenith total delay of the monitoring stations;

[0012] The second calculation unit is used for constructing a double-difference observation equation based on the second information and the third information, obtaining a GNSS relative positioning observation equation;

[0013] The third calculation unit is used for solving the double-difference observation equation, and determining deformation estimation data of a GNSS monitoring point based on a Kalman filtering algorithm and a calculation result of the double-difference observation equation.

[0014] The beneficial effects of the present application are as follows:

[0015] The present application is used for arranging actual measurement meteorological stations at GNSS reference stations and monitoring stations, calculating prior tropospheric delay error by using actual measurement meteorological parameters when solving, further taking residual tropospheric delay error as a parameter, estimating in a solving process, and realizing high-precision deformation monitoring solution under a large height difference. The present application improves the accuracy of displacement estimation by modeling and calculating tropospheric zenith total delay, and further improves the stability and real-time performance of displacement estimation by applying a Kalman filtering algorithm.

[0016] Other features and advantages of the present application will be described in the following description, and some will become apparent from the description, or will be learned from practice of the present application. The purpose and other advantages of the present application can be achieved and obtained by the structure specially pointed out in the written description, claims, and drawings. BRIEF DESCRIPTION OF DRAWINGS

[0017] In order to more clearly illustrate the technical solutions of the embodiments of the present application, the following will briefly introduce the drawings needed to be used in the embodiments. It should be understood that the following drawings only show some embodiments of the present application, and therefore should not be regarded as limiting the scope. For those skilled in the art, other related drawings can also be obtained without creative labor on the basis of these drawings.

[0018] Figure 1 The flow chart of the large height difference GNSS deformation monitoring method described in the embodiments of the present application;

[0019] Figure 2 A large height difference GNSS deformation monitoring system structure diagram described in the embodiment of the present application.

[0020] In the figure: 701, acquisition unit; 702, first calculation unit; 7021, first calculation subunit; 7022, second calculation subunit; 7023, third calculation subunit; 703, second calculation unit; 704, third calculation unit; 7041, first processing subunit; 7042, fourth calculation subunit; 7043, fifth calculation subunit; 7044, sixth calculation subunit; 7045, second processing subunit; 7046, third processing subunit; 7047, fourth processing subunit; 7048, fifth processing subunit; 7049, sixth processing subunit; 70410, seventh processing subunit; 70411, eighth processing subunit; 70471, ninth processing subunit; 70472, seventh calculation subunit; 70473, tenth processing subunit. DETAILED DESCRIPTION

[0021] In order to make the objects, technical solutions and advantages of the embodiments of the present application clearer, the technical solutions in the embodiments of the present application will be described clearly and completely below with reference to the drawings in the embodiments of the present application. Obviously, the described embodiments are some but not all of the embodiments of the present application. The components of the embodiments of the present application described and shown in the drawings herein can be arranged and designed in various different configurations. Therefore, the following detailed description of the embodiments of the present application provided in the drawings is not intended to limit the scope of the claimed present application, but only represents selected embodiments of the present application. Based on the embodiments in the present application, all other embodiments obtained by those of ordinary skill in the art without creative work are within the scope of protection of the present application.

[0022] It should be noted that: similar reference numerals and letters represent similar items in the following drawings, therefore, once an item is defined in one drawing, it does not need to be further defined and explained in subsequent drawings. Meanwhile, in the description of the present application, the terms “first”, “second” and the like are only used to distinguish description, and cannot be understood as indicating or implying relative importance.

[0023] Embodiment 1:

[0024] The embodiment provides a large height difference GNSS deformation monitoring method.

[0025] Referring to Figure 1 , the method includes steps S1, S2, S3 and S4.

[0026] Step S1, obtaining first information and second information, the first information includes meteorological parameter information of GNSS reference station and monitoring station, the meteorological parameter information includes temperature parameter information, pressure parameter information and water vapor pressure parameter information, the second information includes observation information of GNSS reference station and monitoring station, the observation information includes pseudo-range and phase observation value;

[0027] It can be understood that the temperature parameter information in this step can be measured by a temperature sensor, the pressure parameter information can be measured by a barometer, and the water vapor pressure parameter information can be measured by a humidity sensor; the pseudo-range observation value is an estimate of the time delay of satellite signal propagation, while the phase observation value provides more accurate phase information, which can be used to solve the relative phase of the satellite signal. GNSS receiver receives signals from multiple satellites, measures the propagation time and phase information of these signals. These data are recorded and used for subsequent calculations to estimate the displacement of the monitoring station. Through these data, the displacement estimation can more accurately consider the influence of atmospheric factors, and the reliability of the displacement estimation of the monitoring station is improved.

[0028] Step S2, sending the first information to a preset meteorological parameter model for calculation to obtain third information, the third information includes the zenith tropospheric delay of the GNSS reference station and the zenith tropospheric delay of the monitoring station;

[0029] It can be understood that this step considers the influence of meteorological parameters and provides accurate calculation of the zenith tropospheric delay, and accurate zenith tropospheric delay calculation helps to improve the accuracy and reliability of displacement estimation. In this step, step S2 includes step S21, step S22 and step S23.

[0030] Step S21, calculating the tropospheric static delay based on the first information to obtain fourth information, the fourth information includes the tropospheric static delay of the GNSS reference station and the monitoring station;

[0031] The calculation formula of the tropospheric static delay calculation in this step is as follows:

[0032]

[0033] Wherein, ZHD is the tropospheric static delay, P is the atmospheric pressure at the reference station or the monitoring station, φ is the latitude of the reference station or the monitoring station; H is the site elevation of the reference station or the monitoring station.

[0034] Step S22, calculating the tropospheric wet delay based on the first information to obtain fifth information, the fifth information includes the tropospheric wet delay of the GNSS reference station and the monitoring station;

[0035] The calculation formula of the tropospheric wet delay calculation in this step is as follows:

[0036]

[0037] wherein ZWD is the tropospheric wet delay, T is the air temperature at the reference station or monitoring station, and e is the water vapor pressure at the reference station or monitoring station.

[0038] Step S23, tropospheric zenith total delay calculation is performed based on the first sub-information and the second sub-information, to obtain third information, wherein the third information comprises tropospheric zenith total delays of the GNSS reference station and the monitoring station.

[0039] The calculation formula of the tropospheric zenith total delay calculation in this step is as follows:

[0040] ZTD = ZHD + ZWD

[0041] wherein ZTD is the tropospheric zenith total delay, ZWD is the tropospheric wet delay, and ZHD is the tropospheric hydrostatic delay.

[0042] Step S3, a double-difference observation equation is constructed based on the second information and the third information, to obtain a GNSS relative positioning observation equation.

[0043] It can be understood that the double-difference observation equation is constructed by using the second information and the third information in this step. This equation describes the distance difference between different satellites and the influence of the tropospheric delay on signal propagation, so the double-difference observation equation is named as the GNSS relative positioning observation equation in this step, and the formula of the double-difference observation equation is as follows:

[0044]

[0045]

[0046] wherein, is a double-difference phase observation value, is a double-difference satellite-ground distance, is a prior double-difference tropospheric delay, λ f is a wavelength, is a double-difference ambiguity parameter, is a double-difference phase observation error; is a double-difference pseudo-range observation value, is a double-difference pseudo-range observation error, is a remaining tropospheric delay error to be estimated.

[0047] Step S4, the double-difference observation equation is solved, and deformation estimation data of a GNSS monitoring point is determined based on a Kalman filtering algorithm and a calculation result of the double-difference observation equation.

[0048] It is understood that the deformation estimation data in this step refers to the displacement estimation data of the GNSS monitoring points. Specifically, by combining the Kalman filtering algorithm with the calculation results of the double-difference observation equation, a more accurate displacement estimation of the GNSS monitoring points is achieved. The Kalman filter possesses temporal and adaptive properties, enabling real-time updates to the displacement estimation and improving its stability and accuracy. Step S4 includes steps S41, S42, S43, S44, S45, and S46.

[0049] Step S41: Define the position, velocity and clock difference of the GNSS monitoring point based on the second information to obtain the fourth information, which includes state vector data, state transition matrix and system noise covariance matrix;

[0050] It is understandable that the state vector in this step includes position, velocity, and clock error information, which will be updated in the iteration of Kalman filtering to provide accurate displacement estimation.

[0051] Step S42: Initialize the state vector data and the state estimation error covariance matrix, and use the state estimation error covariance matrix and the fourth information to send a preset prediction equation for calculation to predict the fifth information, which includes the state vector data and the state estimation error covariance matrix of the next time step.

[0052] Understandably, this step initializes the state vector and the state estimation error covariance matrix, and then uses them for state prediction. This helps establish the initial state for the Kalman filter algorithm, facilitating subsequent observation updates and displacement estimation. The pre-defined prediction equation is as follows:

[0053] X k =F k,k-1 *X k-1

[0054]

[0055] Among them, P k,k-1 Let F represent the state estimation error covariance matrix. k.k-1 Let P represent the state transition matrix. k-1 This represents the covariance of the state estimation error at time step k-1. Q represents the transpose of the state transition matrix. k-1 Let X represent the system noise covariance matrix at time step k-1. k Let X represent the state vector at time step k. k-1 This represents the state vector at time step (k-1).

[0056] Step S43, solving the double-difference observation equation to obtain the sixth information, the sixth information including an observation vector, an observation matrix and an observation noise covariance matrix;

[0057] It can be understood that after solving the double-difference observation equation, the observation vector is obtained, which includes pseudo-range and phase observation values from GNSS satellites. Among them, the observation matrix and the observation noise covariance matrix help to associate these observation values with the state vector for subsequent Kalman filter observation update. This helps to improve the accuracy and reliability of displacement estimation.

[0058] Step S44, performing Kalman gain calculation based on the fifth information and the sixth information to obtain a Kalman gain value;

[0059] It can be understood that this step determines how to balance the prior prediction and the actual observation value by calculating the Kalman gain to minimize the estimation error.

[0060] Step S45, if the Kalman gain value is within a preset threshold range, iteratively updating the state vector data and the state estimation error covariance matrix to obtain the state vector data and the state estimation error covariance matrix at each time step;

[0061] It can be understood that this step ensures the stability and accuracy of state estimation by iteratively updating the state vector and the state estimation error covariance matrix, while avoiding unnecessary updates to improve computational efficiency. The formula for iteratively updating the state vector data and the state estimation error covariance matrix is as follows:

[0062]

[0063] y=Z k -H k *X k,k-1

[0064] X k =X k,k-1 +K k *y

[0065] P k =(I-K k *H k )*P k,k-1

[0066] Where K k represents the Kalman gain value at the kth time step, P k,k-1 state estimation error covariance matrix, H k represents the observation matrix at the kth time step, represents the transpose of the observation matrix at the kth time step, R kRkdenotes the observation noise covariance matrix at the kth time step, y denotes the observation residual, Z k Rkdenotes the observation noise covariance matrix at the kth time step, y denotes the observation residual, Z k,k-1 Rkdenotes the observation noise covariance matrix at the kth time step, y denotes the observation residual, Z k Rkdenotes the observation noise covariance matrix at the kth time step, y denotes the observation residual, Z k Rkdenotes the observation noise covariance matrix at the kth time step, y denotes the observation residual, Z

[0067] Step S46, determining the deformation estimation data of the GNSS monitoring point based on the state vector data and the state estimation error covariance matrix of each time step.

[0068] It can be understood that this step determines the estimated position of the GNSS monitoring point based on the state vector data and the state estimation error covariance matrix of each time step, and after step S46, steps S47, S48, S49, S410 and S411 are further included.

[0069] Step S47, performing error judgment on the preset historical estimated position of the GNSS monitoring point and the historical actual position of the GNSS monitoring point, and obtaining an abnormal threshold range based on the error judgment result;

[0070] It can be understood that this step establishes an abnormal threshold range by comparing the historical estimated position and the actual position, which is used to judge whether the estimated position is abnormal, wherein step S47 includes steps S471, S472 and S473.

[0071] Step S471, performing clustering processing on all the error judgment results based on a K-means algorithm to obtain at least one cluster set, each cluster set containing at least one error judgment result;

[0072] Step S472, calculating a threshold range corresponding to each cluster set based on all the cluster sets and a Lya-punov criterion;

[0073] Step S474, analyzing all the threshold ranges, and taking the minimum threshold range formed by all the threshold ranges as the threshold range for judging error as abnormal.

[0074] It can be understood that this step determines the range of abnormal data by using a clustering algorithm and a Lya-punov criterion, so that the abnormal detection is more accurate and reliable, which helps to improve the quality and reliability of the displacement estimation of the GNSS monitoring point.

[0075] Step S48, performing abnormal type labeling on the preset historical estimated position of the GNSS monitoring point based on the abnormal threshold range, to obtain labeled abnormal estimation data;

[0076] It can be understood that this step calibrates the abnormal estimation data to further detect and analyze the abnormality.

[0077] Step S49, constructing a CART decision tree based on the CART algorithm and the calibrated abnormal estimation data, wherein the CART decision tree is randomly pruned and the constant of the CART decision tree is determined to obtain at least one untrained sub-decision tree;

[0078] It can be understood that this step helps to establish the basis of the abnormal class recognition model by constructing the CART decision tree and performing random pruning, so as to more accurately determine the abnormality of the estimated position.

[0079] Step S410, obtaining an optimal sub-decision tree based on the untrained sub-decision tree and the Gini index calculation method, and obtaining the abnormal class recognition model based on the optimal sub-decision tree, wherein the abnormal class recognition model comprises the optimal sub-decision tree and its corresponding target constant;

[0080] It can be understood that this step helps to improve the accuracy and reliability of abnormal detection by calculating the optimal sub-decision tree and constructing the abnormal class recognition model, so as to identify different types of abnormal situations.

[0081] Step S411, sending the estimated position of the GNSS monitoring point to the trained abnormal class recognition model for judgment, and determining whether the estimated position of the GNSS monitoring point is abnormal based on the judgment result of the abnormal class recognition model.

[0082] It can be understood that this step judges whether the estimated position of the GNSS monitoring point is abnormal by using the trained abnormal class recognition model.

[0083] Embodiment 2:

[0084] As shown in Figure 2 The present embodiment provides a GNSS deformation monitoring system for large height difference, as shown in Figure 2 The system comprises an acquisition unit 701, a first calculation unit 702, a second calculation unit 703 and a third calculation unit 704.

[0085] The acquisition unit 701 is configured to acquire first information and second information, wherein the first information comprises meteorological parameter information of GNSS reference stations and monitoring stations, and the meteorological parameter information comprises temperature parameter information, pressure parameter information and water vapor pressure parameter information; the second information comprises observation information of the GNSS reference stations and the monitoring stations, and the observation information comprises pseudorange and phase observation value.

[0086] The first computing unit 702 is configured to send the first information to a preset meteorological parameter model respectively to obtain third information, wherein the third information comprises tropospheric zenith total delay of the GNSS reference station and tropospheric zenith total delay of the monitoring station.

[0087] The first computing unit 702 comprises a first computing subunit 7021, a second computing subunit 7022 and a third computing subunit 7023.

[0088] The first computing subunit 7021 is configured to perform tropospheric static delay calculation based on the first information to obtain fourth information, wherein the fourth information comprises tropospheric static delay of the GNSS reference station and the monitoring station.

[0089] The second computing subunit 7022 is configured to perform tropospheric wet delay calculation based on the first information to obtain fifth information, wherein the fifth information comprises tropospheric wet delay of the GNSS reference station and the monitoring station.

[0090] The third computing subunit 7023 is configured to perform tropospheric zenith total delay calculation based on the first sub-information and the second sub-information to obtain the third information, wherein the third information comprises tropospheric zenith total delay of the GNSS reference station and the monitoring station.

[0091] The second computing unit 703 is configured to construct a double-difference observation equation based on the second information and the third information to obtain a GNSS relative positioning observation equation.

[0092] The third computing unit 704 is configured to solve the double-difference observation equation and determine deformation estimation data of the GNSS monitoring point based on a Kalman filtering algorithm and a calculation result of the double-difference observation equation.

[0093] The third computing unit 704 comprises a first processing subunit 7041, a fourth computing subunit 7042, a fifth computing subunit 7043, a sixth computing subunit 7044, a second processing subunit 7045 and a third processing subunit 7046.

[0094] The first processing subunit 7041 is configured to define position, velocity and clock error of the GNSS monitoring point based on the second information to obtain fourth information, wherein the fourth information comprises state vector data, state transition matrix and system noise covariance matrix.

[0095] The fourth computing subunit 7042 is configured to initialize the state vector data and the state estimation error covariance matrix, and use the state estimation error covariance matrix and the fourth information to send a preset prediction equation to perform calculation to predict fifth information, wherein the fifth information comprises state vector data and state estimation error covariance matrix of the next time step.

[0096] The fifth calculation sub-unit 7043 is configured to solve the double-difference observation equation to obtain sixth information, wherein the sixth information comprises an observation vector, an observation matrix, and an observation noise covariance matrix.

[0097] The sixth calculation sub-unit 7044 is configured to perform Kalman gain calculation based on the fifth information and the sixth information to obtain a Kalman gain value.

[0098] The second processing sub-unit 7045 is configured to, if the Kalman gain value is within a preset threshold range, iteratively update the state vector data and the state estimation error covariance matrix to obtain state vector data and state estimation error covariance matrix of each time step.

[0099] The third processing sub-unit 7046 is configured to determine the deformation estimation data of the GNSS monitoring point based on the state vector data and the state estimation error covariance matrix of each time step.

[0100] The third processing sub-unit 7046 further comprises a fourth processing sub-unit 7047, a fifth processing sub-unit 7048, a sixth processing sub-unit 7049, a seventh processing sub-unit 70410, and an eighth processing sub-unit 70411.

[0101] The fourth processing sub-unit 7047 is configured to perform error judgment on a preset historical estimated position of a GNSS monitoring point and a historical actual position of the GNSS monitoring point, and obtain an abnormal threshold range based on the error judgment result.

[0102] The fourth processing sub-unit 7047 comprises a ninth processing sub-unit 70471, a seventh calculation sub-unit 70472, and a tenth processing sub-unit 70473.

[0103] The ninth processing sub-unit 70471 is configured to perform clustering processing on all the error judgment results based on a K-means algorithm to obtain at least one cluster set, wherein each cluster set comprises at least one error judgment result.

[0104] The seventh calculation sub-unit 70472 is configured to calculate a threshold range corresponding to each cluster set based on all the cluster sets and a Lyaupunov criterion.

[0105] The tenth processing sub-unit 70473 is configured to analyze all the threshold ranges, and take a minimum threshold range formed by all the threshold ranges as a threshold range for judging error as abnormal.

[0106] The fifth processing sub-unit 7048 is configured to perform abnormal type labeling on a preset historical estimated position of a GNSS monitoring point based on the abnormal threshold range to obtain labeled abnormal estimation data.

[0107] The sixth processing subunit 7049 is used to construct a CART decision tree based on the CART algorithm and the calibrated anomaly estimation data, wherein the CART decision tree is randomly pruned and the constants of the CART decision tree are determined to obtain at least one untrained sub-decision tree.

[0108] The seventh processing subunit 70410 is used to obtain the optimal sub-decision tree based on the untrained sub-decision tree and the Gini index calculation method, and to obtain the anomaly category recognition model based on the optimal sub-decision tree. The anomaly category recognition model includes the optimal sub-decision tree and its corresponding target constant.

[0109] The eighth processing subunit 70411 is used to send the estimated position of the GNSS monitoring point to the trained anomaly category identification model for judgment, and determine whether the estimated position of the GNSS monitoring point is abnormal based on the judgment result of the anomaly category identification model.

[0110] It should be noted that the specific methods by which each module performs operations in the system described in the above embodiments have been described in detail in the embodiments related to the method, and will not be elaborated here.

[0111] The above description is merely a preferred embodiment of the present invention and is not intended to limit the invention. Various modifications and variations can be made to the present invention by those skilled in the art. Any modifications, equivalent substitutions, improvements, etc., made within the spirit and principles of the present invention should be included within the scope of protection of the present invention.

[0112] The above description is merely a specific embodiment of the present invention, but the scope of protection of the present invention is not limited thereto. Any variations or substitutions that can be easily conceived by those skilled in the art within the technical scope disclosed in the present invention should be included within the scope of protection of the present invention. Therefore, the scope of protection of the present invention should be determined by the scope of the claims.

Claims

1. A method for GNSS deformation monitoring with large height difference, characterized in that, The method comprises the following steps: acquiring first information and second information, the first information comprising meteorological parameter information of GNSS reference stations and monitoring stations, the meteorological parameter information comprising temperature parameter information, pressure parameter information and water vapor pressure parameter information, the second information comprising observation information of GNSS reference stations and monitoring stations, the observation information comprising pseudo-range and phase observation values; sending the first information to a preset meteorological parameter model for calculation to obtain third information, the third information comprising zenith tropospheric total delay of GNSS reference stations and monitoring stations; constructing a double-difference observation equation based on the second information and the third information to obtain a GNSS relative positioning observation equation; solving the double-difference observation equation and determining GNSS monitoring point deformation estimation data based on a Kalman filtering algorithm and a calculation result of the double-difference observation equation; wherein, sending the first information to a preset meteorological parameter model for calculation to obtain third information comprises: performing tropospheric static delay calculation based on the first information to obtain first sub-information, the first sub-information comprising tropospheric static delay of GNSS reference stations and monitoring stations; performing tropospheric wet delay calculation based on the first information to obtain second sub-information, the second sub-information comprising tropospheric wet delay of GNSS reference stations and monitoring stations; performing tropospheric zenith total delay calculation based on the first sub-information and the second sub-information to obtain third information, the third information comprising zenith tropospheric total delay of GNSS reference stations and monitoring stations; wherein, solving the double-difference observation equation and determining GNSS monitoring point deformation estimation data based on a Kalman filtering algorithm and a calculation result of the double-difference observation equation comprises: defining the position, velocity and clock error of the GNSS monitoring point based on the second information to obtain fourth information, the fourth information comprising state vector data, state transition matrix and system noise covariance matrix; initializing the state vector data and state estimation error covariance matrix, and sending the state estimation error covariance matrix and the fourth information to a preset prediction equation for calculation to predict fifth information, the fifth information comprising state vector data and state estimation error covariance matrix of the next time step; solving the double-difference observation equation to obtain sixth information, the sixth information comprising observation vector, observation matrix and observation noise covariance matrix; performing Kalman gain calculation based on the fifth information and the sixth information to obtain Kalman gain value; if the Kalman gain value is within a preset threshold range, iteratively updating the state vector data and state estimation error covariance matrix in the fifth information to obtain state vector data and state estimation error covariance matrix of each time step; determining the GNSS monitoring point deformation estimation data based on the state vector data and state estimation error covariance matrix of each time step; wherein, after determining the GNSS monitoring point deformation estimation data based on the state vector data and state estimation error covariance matrix of each time step, the method further comprises: The preset GNSS monitoring point historical estimated position is compared with the GNSS monitoring point historical actual position to determine the error, and an abnormal threshold range is obtained based on the error determination result; An abnormal type of the preset GNSS monitoring point historical estimated position is labeled based on the abnormal threshold range, and labeled abnormal estimated data is obtained; A CART decision tree is constructed based on the CART algorithm and the labeled abnormal estimated data, wherein the CART decision tree is randomly pruned and the constant of the CART decision tree is determined to obtain at least one untrained sub-decision tree; An optimal sub-decision tree is obtained based on the untrained sub-decision tree and the Gini index calculation method, and an abnormal category recognition model is obtained based on the optimal sub-decision tree, wherein the abnormal category recognition model comprises the optimal sub-decision tree and a corresponding target constant; The GNSS monitoring point estimated position is sent to the trained abnormal category recognition model for judgment, and whether the GNSS monitoring point estimated position is abnormal is determined based on the judgment result of the abnormal category recognition model.

2. The large-height-difference GNSS deformation monitoring method according to claim 1, characterized in that , and an abnormal threshold range is obtained based on the error determination result, comprising: All the error determination results are clustered based on the K-means algorithm to obtain at least one cluster set, and each cluster set comprises at least one error determination result; The threshold range corresponding to each cluster set is calculated based on all the cluster sets and the Lyaupin criterion; All the threshold ranges are analyzed, and the minimum threshold range formed by all the threshold ranges is taken as the threshold range for determining the error as abnormal.

3. A GNSS deformation monitoring system for large height differences, characterized in that comprising: An acquisition unit is configured to acquire first information and second information, wherein the first information comprises meteorological parameter information of GNSS reference stations and monitoring stations, the meteorological parameter information comprises temperature parameter information, pressure parameter information, and water vapor pressure parameter information, and the second information comprises observation information of the GNSS reference stations and the monitoring stations, the observation information comprises pseudo-range and phase observation values; A first calculation unit is configured to send the first information to a preset meteorological parameter model for calculation to obtain third information, wherein the third information comprises zenith total delay of the GNSS reference stations and zenith total delay of the monitoring stations; A second calculation unit is configured to construct a double-difference observation equation based on the second information and the third information to obtain a GNSS relative positioning observation equation; A third calculation unit is configured to solve the double-difference observation equation and determine deformation estimation data of GNSS monitoring points based on a Kalman filtering algorithm and a calculation result of the double-difference observation equation; The first calculation unit comprises: A first calculation sub-unit is configured to perform tropospheric static delay calculation based on the first information to obtain first sub-information, wherein the first sub-information comprises tropospheric static delay of the GNSS reference stations and the monitoring stations; A second calculation sub-unit is configured to perform tropospheric wet delay calculation based on the first information to obtain second sub-information, wherein the second sub-information comprises tropospheric wet delay of the GNSS reference stations and the monitoring stations; The third calculation subunit is configured to calculate tropospheric zenith total delay based on the first sub-information and the second sub-information to obtain third information, wherein the third information comprises tropospheric zenith total delay of a GNSS reference station and tropospheric zenith total delay of a monitoring station; The third calculation unit comprises: The first processing subunit is configured to define position, velocity and clock bias of the GNSS monitoring point based on the second information to obtain fourth information, wherein the fourth information comprises state vector data, state transition matrix and system noise covariance matrix; The fourth calculation subunit is configured to initialize the state vector data and state estimation error covariance matrix, and send the state estimation error covariance matrix and the fourth information to a preset prediction equation for calculation to predict fifth information, wherein the fifth information comprises state vector data and state estimation error covariance matrix of a next time step; The fifth calculation subunit is configured to solve the double-difference observation equation to obtain sixth information, wherein the sixth information comprises observation vector, observation matrix and observation noise covariance matrix; The sixth calculation subunit is configured to calculate Kalman gain based on the fifth information and the sixth information to obtain a Kalman gain value; The second processing subunit is configured to iteratively update the state vector data and state estimation error covariance matrix in the fifth information if the Kalman gain value is within a preset threshold range to obtain state vector data and state estimation error covariance matrix of each time step; The third processing subunit is configured to determine deformation estimation data of the GNSS monitoring point based on the state vector data and state estimation error covariance matrix of each time step; The third calculation unit further comprises: The fourth processing subunit is configured to judge error between a preset historical estimated position of the GNSS monitoring point and a historical actual position of the GNSS monitoring point, and obtain an abnormal threshold range based on the error judgment result; The fifth processing subunit is configured to label abnormal types of the preset historical estimated position of the GNSS monitoring point based on the abnormal threshold range to obtain labeled abnormal estimation data; The sixth processing subunit is configured to construct a CART decision tree based on a CART algorithm and the labeled abnormal estimation data, wherein the CART decision tree is randomly pruned and a constant of the CART decision tree is determined to obtain at least one untrained sub-decision tree; The seventh processing subunit is configured to obtain an optimal sub-decision tree based on the untrained sub-decision tree and a Gini index calculation method, and obtain an abnormal category recognition model based on the optimal sub-decision tree, wherein the abnormal category recognition model comprises the optimal sub-decision tree and a corresponding target constant thereof; The eighth processing subunit is configured to send the estimated position of the GNSS monitoring point to the trained abnormal category recognition model for judgment, and determine whether the estimated position of the GNSS monitoring point is abnormal based on a judgment result of the abnormal category recognition model.

4. The large-height-difference GNSS deformation monitoring system according to claim 3, characterized in that, The fourth processing subunit comprises: a ninth processing subunit configured to perform clustering processing on all of the error judgment results based on a K-means algorithm to obtain at least one cluster set, each cluster set containing at least one error judgment result; a seventh calculation subunit configured to calculate a threshold range corresponding to each cluster set based on all of the cluster sets and a Laplace criterion; a tenth processing subunit configured to analyze all of the threshold ranges and take a minimum threshold range formed by all of the threshold ranges as a threshold range for judging errors as abnormal.

Images