CN117739797A - Beidou/GNSS-based multi-time-scale deformation monitoring method - Google Patents
Beidou/GNSS-based multi-time-scale deformation monitoring method Download PDFInfo
- Publication number
- CN117739797A CN117739797A CN202311517400.5A CN202311517400A CN117739797A CN 117739797 A CN117739797 A CN 117739797A CN 202311517400 A CN202311517400 A CN 202311517400A CN 117739797 A CN117739797 A CN 117739797A
- Authority
- CN
- China
- Prior art keywords
- representing
- station
- difference
- stations
- monitoring
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Granted
Links
- 238000012544 monitoring process Methods 0.000 title claims abstract description 154
- 238000000034 method Methods 0.000 title claims abstract description 35
- 238000007781 pre-processing Methods 0.000 claims abstract description 39
- 238000007667 floating Methods 0.000 claims abstract description 21
- 238000001914 filtration Methods 0.000 claims abstract description 13
- 238000001514 detection method Methods 0.000 claims abstract description 10
- 238000004364 calculation method Methods 0.000 claims abstract description 9
- 239000005433 ionosphere Substances 0.000 claims abstract description 7
- 239000011159 matrix material Substances 0.000 claims description 24
- 238000003908 quality control method Methods 0.000 claims description 17
- 238000012545 processing Methods 0.000 claims description 16
- 238000012937 correction Methods 0.000 claims description 11
- 238000013461 design Methods 0.000 claims description 8
- 230000008030 elimination Effects 0.000 claims description 6
- 238000003379 elimination reaction Methods 0.000 claims description 6
- 239000013585 weight reducing agent Substances 0.000 claims description 6
- 239000005436 troposphere Substances 0.000 claims description 4
- 238000003491 array Methods 0.000 claims description 3
- 230000035772 mutation Effects 0.000 abstract description 6
- 230000009286 beneficial effect Effects 0.000 abstract description 3
- 238000012360 testing method Methods 0.000 description 12
- 230000007547 defect Effects 0.000 description 5
- 238000006073 displacement reaction Methods 0.000 description 5
- 238000005516 engineering process Methods 0.000 description 2
- 238000012986 modification Methods 0.000 description 2
- 230000004048 modification Effects 0.000 description 2
- 238000012216 screening Methods 0.000 description 2
- 238000012300 Sequence Analysis Methods 0.000 description 1
- 230000037237 body shape Effects 0.000 description 1
- 238000004891 communication Methods 0.000 description 1
- 238000010276 construction Methods 0.000 description 1
- 230000007812 deficiency Effects 0.000 description 1
- 238000010586 diagram Methods 0.000 description 1
- 230000000694 effects Effects 0.000 description 1
- 230000000873 masking effect Effects 0.000 description 1
- 230000003287 optical effect Effects 0.000 description 1
- 238000002407 reforming Methods 0.000 description 1
- 230000011218 segmentation Effects 0.000 description 1
- 230000004304 visual acuity Effects 0.000 description 1
Landscapes
- Position Fixing By Use Of Radio Waves (AREA)
Abstract
The invention discloses a Beidou/GNSS-based multi-time scale deformation monitoring method, which is used for performing cycle slip and gross error detection, satellite position and clock error calculation and monitoring station rough coordinate calculation based on Beidou/GNSS original observation data of a reference station and a monitoring station and broadcast ephemeris information; constructing an ionosphere weighting single difference equation between stations, estimating a floating solution of the coordinates of the monitoring station by adopting a Kalman filtering algorithm, and obtaining an instantaneous fixed solution of the coordinates of the monitoring station after partial ambiguity is fixed; storing the pre-processing data; when K hours are accumulated, constructing a multi-epoch station star double-difference observation equation according to the first preprocessing data to obtain a monitoring station coordinate hour scale fixed solution; and taking the floating solution and the estimated value variance thereof as the observed value and the weight inverse thereof when accumulating for one day, and estimating the coordinate space solution of the monitoring station by using least square. The invention can provide high-frequency, low-frequency and long-trend multi-time scale deformation monitoring information, and is beneficial to quick response of large mutation and forecast of potential micro deformation.
Description
Technical Field
The invention relates to the field of deformation monitoring, in particular to a Beidou/GNSS-based multi-time-scale deformation monitoring method.
Background
The GNSS technology is one of important means for monitoring the body shape of important natural or artificial structures such as slopes, dams, bridges and the like due to the characteristics of low cost, high precision, all weather and the like. The GNSS deformation monitoring technology usually adopts a "1+n" working mode, i.e. a reference station and N monitoring stations, wherein the reference station is generally built on a bedrock with stable geological environment, and the monitoring stations are built on a monitoring deformation body. On the basis, common errors such as satellite clock errors, orbit errors and the like are eliminated through short base line difference, and most of spatial correlation errors such as ionosphere delay, troposphere delay and the like are weakened, so that the relative position of the monitoring station relative to the reference station is obtained, and the three-dimensional displacement of the monitoring station, namely the surface deformation of the monitoring point, is further obtained according to the time sequence of the relative position. In the prior art, although the method has a certain effect on millimeter-level slow deformation monitoring which adopts a high-performance resolving platform and has better GNSS observation data quality, the method has the following defects for large-scale multi-scale deformation monitoring:
1) The solution amount is large. In order to obtain the millimeter-level precision relative position, the observation values of enough epochs need to be collected, and the baseline is solved through integral least square, so that the method comprises more original observation values and parameters to be estimated, and the data preprocessing and the parameter estimation are time-consuming, so that high-frequency deformation monitoring is difficult to realize even if a sliding window is used;
2) Masking large mutations. The premise of the integral least squares solution of the base line is that the coordinates of the arc segment monitoring station are approximately unchanged, and in the actual environment, the monitoring body is greatly deformed due to natural or artificial factors such as rainfall, earthquake or construction, and the integral least squares solution of the base line can average the mutation to each epoch, and the resolving interval is larger, so that the larger mutation of the monitoring body cannot be reflected rapidly;
3) A single time scale. The method has the advantages that one set of program only outputs the displacement information of the surface of the monitoring body with one time scale, and the coordinate value with millimeter-level precision is generally 2-4 hours, so that the method can be used for monitoring slow micro deformation and is difficult to reflect high-frequency abrupt change or longer (one day) trend deformation;
4) Poor resistance. The lack of an efficient quality control strategy, i.e. the lack of detection and rejection of cycle slip, gross errors, etc. in GNSS observations, will affect the reliability of deformation monitoring.
Disclosure of Invention
In order to solve the technical problems, the embodiment of the invention provides a multi-time scale deformation monitoring method based on Beidou/GNSS.
The embodiment of the invention provides a multi-time scale deformation monitoring method based on Beidou/GNSS, which comprises the following steps:
based on Beidou/GNSS original observation data and broadcast ephemeris information of a reference station and a monitoring station obtained in real time, cycle slip and rough detection, satellite position and clock error calculation and monitoring station sketch coordinate calculation are carried out;
based on the rough coordinates of the monitoring stations and the satellite positions, constructing a design matrix, a residual error matrix and a weight matrix of a single difference equation among ionosphere weighted stations, estimating floating solutions of the coordinates of the monitoring stations in real time by adopting a Kalman filtering algorithm with additional quality control according to the design matrix, the residual error matrix and the weight matrix, and carrying out partial ambiguity fixation on the floating solutions to obtain instantaneous fixed solutions of the coordinates of the monitoring stations as high-frequency deformation monitoring information;
storing satellite altitude angle, observation value residual error, ambiguity cycle slip information and satellite-to-monitoring station direction cosine as first preprocessing data in real time, and taking the floating solution of the monitoring station coordinates and estimation variances thereof as second preprocessing data; wherein the satellite altitude is determined by the satellite position and the observation residual is determined by the residual array;
when the time for storing the first preprocessing data is accumulated to K hours, constructing a multi-epoch station star double-difference observation equation according to the first preprocessing data, and carrying out least square estimation and partial ambiguity fixation according to the multi-epoch station star double-difference observation equation to obtain a monitoring station coordinate hour scale fixed solution as low-frequency deformation monitoring information; k is greater than or equal to 1;
and when the time for storing the second preprocessing data is accumulated to one day, taking the floating solution of the monitoring station coordinates and the estimated value variance thereof as an observed value and a weight inverse array thereof respectively, and estimating the monitoring station coordinates and the day solution by utilizing least square according to the observed value and the weight inverse array thereof to be taken as long-trend deformation monitoring information.
Further, the partial ambiguity fixing is performed by an LAMBDA algorithm; the instantaneous fixed solution of the coordinates of the monitoring station is obtained through epoch-by-epoch real-time solution, and the precision is in the centimeter level and is used for representing burst large deformation.
Further, the multi-epoch station star double difference observation equation is constructed by the following steps:
for a plurality of continuous epochs, caching first preprocessing data of each epoch in real time;
under the condition that first preprocessing data of N epochs are acquired and accumulated for K hours, a reference star is selected, the number of ambiguity to be estimated is determined, and a multi-epoch station star double difference observation equation is constructed based on the reference star, the number of ambiguity to be estimated and the first preprocessing data of the N epochs; wherein N is a preset epoch number threshold.
Further, when the time of storing the second preprocessing data is accumulated to one day, the floating solution and the estimated variance of the monitoring station coordinate are respectively used as an observed value and a weight inverse array thereof, and the monitoring station coordinate solution is estimated by using least square according to the observed value and the weight inverse array thereof, and the method comprises the following steps:
for a plurality of continuous epochs, caching second preprocessing data of each epoch in real time;
under the condition that second preprocessing data of M epochs are acquired, constructing a baseline solution quadratic adjustment based on the observed values of the M epochs and the weight inverse arrays thereof, and solving the baseline solution quadratic adjustment by using least square estimation to obtain a monitoring station coordinate day solution; wherein M is the minimum epoch number of the represented time which is more than or equal to one day.
Further, the cycle slip and gross error detection specifically includes:
removing coarse difference data contained in the Beidou/GNSS original observation data;
and marking and initializing cycle slip data contained in the Beidou/GNSS original observation data.
Further, the ionospheric weighted inter-station single difference equation is determined by the following formula (1) -formula (2):
wherein a represents a reference station and b represents a monitoring station; i represents an epoch number, j represents a frequency number;a pseudo-range observation representing a single difference between stations corresponding to satellites with a satellite number s, +.>Representing the phase observations corresponding to satellites with the satellite number s; />Representing an inter-site Shan Chawei ground clearance approximation; />Receiver clock for indicating single difference between stationsDifference;representing zenithal troposphere wet delay projection function, τ ab (i) Represents zenith tropospheric wet delay; />(i) Represents the single difference ionospheric bias delay, mu, between stations j ＝(λ j /λ 1 ) 2 A scale factor representing ionospheric bias delay, where lambda j The phase wavelength of the j-th frequency point is represented; />Representing the code deviation of a single difference receiver between stations; />Representing the phase deviation of the single-difference receiver between stations; />A station star double difference phase ambiguity relative to the reference station a and the first reference star r; /> Pseudo-range observation noise representing single difference between stations, +.>Representing inter-station single-difference phase observation noise; />The direction cosine from the monitoring station to the satellite is represented; dx (dx) b (i)＝[dx b (i) dy b (i) dz b (i)] T Representing the coordinate corrections of the monitoring station X, Y, Z in three directions; dt (dt) ab (i) Representing the single difference receiver clock difference between the original stations; d, d ab,1 Representing the reference frequency point between stationsCode bias of single difference receiver, d ab,j Representing code deviation of a single-difference receiver between j frequency point stations; delta ab,j Representing phase deviation of a single-difference receiver between j frequency point stations; />Representing the single differential phase ambiguity between the stations of the reference star r.
Further, the multi-epoch station star double difference observation equation is determined by the following formula (3):
wherein a represents a reference station and b represents a monitoring station; i represents an epoch number, j represents a frequency number; q represents the second reference star, and,the coordinate corrections representing the three directions of the monitoring station X, Y, Z are an hour arc segment solution,representing the single difference phase ambiguity between stations relative to reference station a and second reference star q, < >>A pseudo-range observation representing a single difference between stations corresponding to satellites with a satellite number s, +.>Pseudo-range observation value representing single difference between stations corresponding to second reference star, < >>Representing inter-station single difference phase observations corresponding to satellites with satellite numbers s, +.>Representing the single difference phase observation value between stations corresponding to the second reference star, < >>The directional cosine of the monitoring station b to the satellite s; lambda (lambda) j Representing the phase wavelength.
Further, the baseline solution quadratic deviation is determined by the following formula (4):
wherein 1= [1 1 1 1]Representing an all 1 vector;the coordinate corrections representing the three directions of the monitoring station X, Y, Z are solved for days; i denotes an epoch number.
Further, the quality control includes the steps of:
performing residual size check on the data obtained after the processing of the Kalman filtering algorithm;
if the residual error exceeds a preset residual error threshold value, carrying out weight reduction or elimination processing on an observed value with the largest residual error in the data, carrying out Kalman filtering processing on the data subjected to the weight reduction or elimination processing, and repeating the steps until the residual error is smaller than the preset residual error threshold value, and outputting the data as the floating solution.
In summary, the invention has the following beneficial effects:
1) In terms of resolving pressure: the existing hour scale multi-epoch integral least square method has large resolving capacity and long time consumption, but the embodiment of the invention adopts a lightweight Kalman filtering algorithm, has small resolving pressure from epoch to epoch, and can carry out 1-10 Hz high-frequency resolving; meanwhile, the first preprocessing data can be further processed, so that the resolving efficiency is higher compared with the prior art;
2) In terms of deformation multiscale monitoring: the prior art is mainly focused on millimeter-scale deformation monitoring with a single hour scale, and is easy to submerge for abrupt changes with a centimeter level or more, and difficult to respond in time; for the defect of insufficient potential trend deformation monitoring of the sky scale, early deformation warning is difficult to carry out. According to the embodiment of the invention, the positions of the monitoring stations in the epoch-by-epoch scale, the hour scale and the day scale can be estimated in parallel, so that the surface displacement of the monitoring points can be reflected more comprehensively;
3) In terms of satellite observation data quality control: in the prior art, the data is less subjected to pre-test screening and rejecting, the comprehensive post-test quality control is lacked, and the small gross error of pre-test missed judgment is difficult to find. When the embodiment of the invention is used for processing the main thread, the data is subjected to strict pre-test and post-test quality control, so that the accuracy of system deformation monitoring can be ensured.
Drawings
FIG. 1 is a flow chart of an embodiment of a Beidou/GNSS-based multi-time scale deformation monitoring method provided by the invention;
fig. 2 is a schematic diagram of a multi-time scale deformation monitoring method based on beidou/GNSS according to an embodiment of the present invention.
Detailed Description
The following description of the embodiments of the present invention will be made clearly and completely with reference to the accompanying drawings, in which it is apparent that the embodiments described are only some embodiments of the present invention, but not all embodiments. All other embodiments, which can be made by those skilled in the art based on the embodiments of the invention without making any inventive effort, are intended to be within the scope of the invention.
In the description of this application, the terms "first," "second," "third," and the like are used for descriptive purposes only and are not to be construed as indicating or implying relative importance or implicitly indicating the number of technical features indicated. Thus, a feature defining "a first", "a second", "a third", etc. may explicitly or implicitly include one or more such feature. In the description of the present application, unless otherwise indicated, the meaning of "a plurality" is two or more.
In the description of the present application, it should be noted that, unless explicitly specified and limited otherwise, the terms "mounted," "connected," and "connected" are to be construed broadly, and may be either fixedly connected, detachably connected, or integrally connected, for example; can be mechanically or electrically connected; can be directly connected or indirectly connected through an intermediate medium, and can be communication between two elements. The specific meaning of the terms in this application will be understood by those of ordinary skill in the art in a specific context.
In the description of the present application, it should be noted that all technical and scientific terms used herein have the same meaning as commonly understood by one of ordinary skill in the art unless defined otherwise. The terminology used in the description of the present invention is for the purpose of describing particular embodiments only and is not intended to be limiting of the invention, as the particular meaning of the terms described above in this application will be understood to those of ordinary skill in the art in the specific context.
Referring to fig. 1, a flow chart of an embodiment of a multi-time scale deformation monitoring method based on beidou/GNSS provided by the present invention includes steps S1 to S5, specifically as follows:
s1, performing cycle slip and gross error detection, satellite position and clock error calculation and monitoring station sketch coordinate calculation based on Beidou/GNSS original observation data and broadcast ephemeris information of a reference station and a monitoring station, which are acquired in real time;
s2, based on the rough coordinates of the monitoring stations and the satellite positions, constructing a design matrix, a residual error matrix and a weight matrix of a single difference equation among ionosphere weighted stations, estimating a floating solution of the coordinates of the monitoring stations in real time by adopting a Kalman filtering algorithm with additional quality control according to the design matrix, the residual error matrix and the weight matrix, and fixing partial ambiguity of the floating solution to obtain an instantaneous fixed solution of the coordinates of the monitoring stations as high-frequency deformation monitoring information;
s3, storing satellite altitude angle, observation value residual error, ambiguity cycle slip information and satellite-to-monitoring station direction cosine as first preprocessing data in real time, and taking the floating solution of the monitoring station coordinates and estimation variances thereof as second preprocessing data; wherein the satellite altitude is determined by the satellite position and the observation residual is determined by the residual array;
s4, when the time for storing the first preprocessing data is accumulated to K hours, constructing a multi-epoch station star double difference observation equation according to the first preprocessing data, and carrying out least square estimation and partial ambiguity fixation according to the multi-epoch station star double difference observation equation to obtain a monitoring station coordinate hour scale fixation solution as low-frequency deformation monitoring information; k is greater than or equal to 1;
and S5, when the time for storing the second preprocessing data is accumulated to one day, taking the floating solution of the monitoring station coordinates and the estimated value variance thereof as an observed value and a weight inverse array thereof respectively, and estimating the monitoring station coordinates and the day solution by utilizing least square according to the observed value and the weight inverse array thereof as long-trend deformation monitoring information.
In specific implementation, referring to fig. 2, S1-S2 in this embodiment are used as the main thread therein, that is, it is to be understood that this embodiment may be applied alone or in combination with the auxiliary thread S4 and/or S5, where the combination with the auxiliary thread is exemplified by the following other embodiments.
Illustratively, the floating solution is fixed in partial ambiguity by using the LAMBDA algorithm.
The method further comprises the steps of: and carrying out deformation early warning based on the instantaneous fixed solution of the coordinates of the monitoring station, the fixed solution of the hour scale of the coordinates of the monitoring station and/or the weather solution of the coordinates of the monitoring station.
For example, the Beidou/GNSS raw observations include Beidou (BDS) observations, GPS observations, GLONASS observations, galileo observations, QZSS observations, and/or IRNSS observations.
In an alternative embodiment, the partial ambiguity fixing is performed by the LAMBDA algorithm; the instantaneous fixed solution of the coordinates of the monitoring station is obtained through epoch-by-epoch real-time solution, and the precision is in the centimeter level and is used for indicating burst large deformation. Therefore, the instantaneous fixed solution of the coordinates of the monitoring station can be used as the first type of deformation monitoring information, namely high-frequency deformation monitoring information.
In an alternative embodiment, the multi-epoch station star double difference observation equation is constructed by:
for a plurality of continuous epochs, caching first preprocessing data of each epoch in real time;
under the condition that first preprocessing data of N epochs are acquired and accumulated for K hours, a reference star is selected, the number of ambiguity to be estimated is determined, and a multi-epoch station star double difference observation equation is constructed based on the reference star, the number of ambiguity to be estimated and the first preprocessing data of the N epochs; wherein N is a preset epoch number threshold.
The fixed solution of the coordinate hour scale of the monitoring station is a sliding arc section solution of a plurality of hours, is smoother and has millimeter precision, and can be used for representing smaller deformation trend, so that the fixed solution is used as second type deformation monitoring information, namely low-frequency deformation monitoring information.
In specific implementation, referring to fig. 2, the embodiment may be a first auxiliary thread, in which, the first auxiliary thread is started to calculate the coordinates of the millimeter-scale monitoring station in an hour scale, when the number of the calendar elements and the calculation interval meet the set threshold, the first auxiliary thread is used to collect the preprocessing data such as the single-difference pseudo-range/phase observation value residual array between stations, satellite position, satellite clock difference, satellite altitude angle, ambiguity cycle skip segmentation information, and the cosine of the monitoring station in the sanitation direction, etc.
In an optional embodiment, when the time of storing the second preprocessed data is accumulated to one day, the floating solution and the estimated variance of the monitoring station coordinate are respectively used as an observed value and a weight inverse array thereof, and the least square is used to estimate the monitoring station coordinate space solution according to the observed value and the weight inverse array thereof, which comprises:
for a plurality of continuous epochs, caching second preprocessing data of each epoch in real time;
under the condition that second preprocessing data of M epochs are acquired, constructing a baseline solution quadratic adjustment based on the observed values of the M epochs and the weight inverse arrays thereof, and solving the baseline solution quadratic adjustment by using least square estimation to obtain a monitoring station coordinate day solution; wherein M is the minimum epoch number of the represented time which is more than or equal to one day.
In specific implementation, referring to fig. 2, the embodiment may be a second auxiliary thread, in which, the correction estimation value and the variance of the dynamic coordinate of the monitoring station are collected from epoch to epoch, when the day of preprocessing data is collected, the second auxiliary thread is started, the day scale sub-millimeter level monitoring station coordinate is resolved, and time sequence analysis is performed.
In an alternative embodiment, the cycle slip and gross error detection specifically includes:
removing coarse difference data contained in the Beidou/GNSS original observation data;
and marking and initializing cycle slip data contained in the Beidou/GNSS original observation data.
In an alternative embodiment, the ionospheric weighted inter-station single difference equation is determined by the following equation (1) -equation (2):
wherein a represents a reference station and b represents a monitoring station; i represents an epoch number, j represents a frequency number;a pseudo-range observation representing a single difference between stations corresponding to satellites with a satellite number s, +.>Representing the phase observations corresponding to satellites with the satellite number s; />Representing an inter-site Shan Chawei ground clearance approximation; />Representing the single difference receiver clock difference between stations;representing zenithal troposphere wet delay projection function, τ ab (i) Represents zenith tropospheric wet delay; />(i) Represents the single difference ionospheric bias delay, mu, between stations j ＝(λ j /λ 1 ) 2 A scale factor representing ionospheric bias delay, where lambda j The phase wavelength of the j-th frequency point is represented; />Representing the code deviation of a single difference receiver between stations; />Representing the phase deviation of the single-difference receiver between stations; />A station star double difference phase ambiguity relative to the reference station a and the first reference star r; /> Pseudo-range observation noise representing single difference between stations, +.>Representing inter-station single-difference phase observation noise; />The direction cosine from the monitoring station to the satellite is represented; dx (dx) b (i)＝[dx b (i) dy b (i) dz b (i)] T Representing the coordinate corrections of the monitoring station X, Y, Z in three directions; dt (dt) ab (i) Representing the single difference receiver clock difference between the original stations; d, d ab,1 Representing code deviation of single-difference receiver between base frequency point stations, d ab,j Representation ofCode deviation of a single-difference receiver between j-th frequency point stations; delta ab,j Representing phase deviation of a single-difference receiver between j frequency point stations; />Representing the single differential phase ambiguity between the stations of the reference star r.
In specific implementation, a current epoch design matrix and a residual matrix can be constructed according to the formula (1), and an observation value weight matrix can be constructed according to a height angle weighting method.
It can be understood that, because of the linear correlation between some parameters in the original observation equation, the design matrix has rank deficiency, and cannot estimate all parameters, so that the parameters need to be reformed, the interested parameters are reserved, and a full rank function model is obtained, and meanwhile, the method adoptsThe method comprises the steps of representing parameters to be estimated after parameter reforming, wherein an estimated form of the reformed parameters relative to original parameters can be represented as a formula (2), and by combining with a Kalman filter and adopting a cycle slip, coarse difference detection and removal and quality control strategy, the real-time dynamic estimation of the relative coordinates of a monitoring station can be realized.
In an alternative embodiment, the station star double difference multiple epoch baseline solution function model is determined by the following equation (3):
wherein a represents a reference station and b represents a monitoring station; i represents an epoch number, j represents a frequency number; q represents the second reference star, and,the coordinate corrections representing the three directions of the monitoring station X, Y, Z are an hour arc segment solution,representing the single difference phase ambiguity between stations relative to reference station a and second reference star q, < >>A pseudo-range observation representing a single difference between stations corresponding to satellites with a satellite number s, +.>Pseudo-range observation value representing single difference between stations corresponding to second reference star, < >>Representing inter-station single difference phase observations corresponding to satellites with satellite numbers s, +.>Representing the single difference phase observation value between stations corresponding to the second reference star, < >>The directional cosine of the monitoring station b to the satellite s; lambda (lambda) j Representing the phase wavelength.
It should be noted that, different from the first reference star being a preset reference star, the second reference star may be selected successively based on the sliding arc segment; from dx b (i)，And regarding the positions of the monitoring stations at all epoch moments in one arc segment as unchanged, namely, the coordinate correction values of all epoch monitoring stations are a group of parameters.
The observation value in the embodiment is obtained after the quality control processing, wherein the satellite arc section with cycle slip is segmented, and the integral least square can be directly carried out, so that the method has stronger robust difference resistance.
In an alternative embodiment, the baseline solution quadratic adjustment is determined by the following equation (4):
wherein 1= [1 1 ] 1]Representing an all 1 vector;the coordinate corrections representing the three directions of the monitoring station X, Y, Z are solved for days; i denotes an epoch number.
In an alternative embodiment, the quality control comprises the steps of:
performing residual size check on the data obtained after the processing of the Kalman filtering algorithm;
if the residual error exceeds a preset residual error threshold value, carrying out weight reduction or elimination processing on an observed value with the largest residual error in the data, carrying out Kalman filtering processing on the data subjected to the weight reduction or elimination processing, and repeating the steps until the residual error is smaller than the preset residual error threshold value, and outputting the data as the floating solution.
Embodiment one:
referring to fig. 2, multi-time scale, multi-precision level monitoring station displacement information is extracted in parallel using multithreading, wherein: 1) The main thread adopts a lightweight Kalman filtering algorithm, so that the problem that the high-frequency monitoring is difficult to carry out due to large resolving power in the existing method is solved, the centimeter-level position solution of the monitoring station can be obtained in real time, and the large mutation high-frequency monitoring of the monitoring point is realized; 2) The main thread adopts an ionosphere weighting inter-station single difference differential model, so that not only is inverse distance weighting constraint carried out on a base line to realize reliable estimation of the position of a monitoring station, but also the pretreatment information such as observation value single difference residual information and the like can be directly provided for the auxiliary thread, and the data pretreatment pressure of the auxiliary thread is reduced; 3) The main thread adopts double guarantees of pre-test gross error, cycle slip detection removal and post-test quality control, and after pretreatment, intermediate information is transmitted to the auxiliary thread for processing, so that the defect of the existing method in data quality control is overcome, and the anti-poor performance of position solution of a monitoring station is improved; 4) The main thread realizes instantaneous scale centimeter-level large mutation monitoring alarm, the auxiliary thread realizes hour scale millimeter-level micro deformation monitoring and early warning, the auxiliary thread realizes day scale sub-millimeter-level trend deformation monitoring and early warning, and the defect of single hour scale monitoring in the existing method is overcome.
In summary, the invention has the following beneficial effects:
1) In terms of resolving pressure: the existing hour scale multi-epoch integral least square method has large resolving capacity and long time consumption, but the embodiment of the invention adopts a lightweight Kalman filtering algorithm, has small resolving pressure from epoch to epoch, and can carry out 1-10 Hz high-frequency resolving; meanwhile, the first preprocessing data can be further processed, so that the resolving efficiency is higher compared with the prior art;
2) In terms of deformation multiscale monitoring: the prior art is mainly focused on millimeter-scale deformation monitoring with a single hour scale, and is easy to submerge for abrupt changes with a centimeter level or more, and difficult to respond in time; for the defect of insufficient potential trend deformation monitoring of the sky scale, early deformation warning is difficult to carry out. According to the embodiment of the invention, the positions of the monitoring stations in the epoch-by-epoch scale, the hour scale and the day scale can be estimated in parallel, so that the surface displacement of the monitoring points can be reflected more comprehensively;
3) In terms of satellite observation data quality control: in the prior art, the data is less subjected to pre-test screening and rejecting, the comprehensive post-test quality control is lacked, and the small gross error of pre-test missed judgment is difficult to find. When the embodiment of the invention is used for processing the main thread, the data is subjected to strict pre-test and post-test quality control, so that the accuracy of system deformation monitoring can be ensured.
From the above description of the embodiments, it will be clear to those skilled in the art that the present invention may be implemented by means of software plus necessary hardware platforms, but may of course also be implemented entirely in hardware. With such understanding, all or part of the technical solution of the present invention contributing to the background art may be embodied in the form of a software product, which may be stored in a storage medium, such as ROM/RAM, a magnetic disk, an optical disk, etc., including several instructions for causing a computer device (which may be a personal computer, a server, or a network device, etc.) to perform the method described in the embodiments or some parts of the embodiments of the present invention.
While the foregoing is directed to the preferred embodiments of the present invention, it will be appreciated by those skilled in the art that changes and modifications may be made without departing from the principles of the invention, such changes and modifications are also intended to be within the scope of the invention.
Claims (9)
1. The multi-time-scale deformation monitoring method based on Beidou/GNSS is characterized by comprising the following steps of:
based on Beidou/GNSS original observation data and broadcast ephemeris information of a reference station and a monitoring station obtained in real time, cycle slip and rough detection, satellite position and clock error calculation and monitoring station sketch coordinate calculation are carried out;
based on the rough coordinates of the monitoring stations and the satellite positions, constructing a design matrix, a residual error matrix and a weight matrix of a single difference equation among ionosphere weighted stations, estimating floating solutions of the coordinates of the monitoring stations in real time by adopting a Kalman filtering algorithm with additional quality control according to the design matrix, the residual error matrix and the weight matrix, and carrying out partial ambiguity fixation on the floating solutions to obtain instantaneous fixed solutions of the coordinates of the monitoring stations as high-frequency deformation monitoring information;
storing satellite altitude angle, observation value residual error, ambiguity cycle slip information and satellite-to-monitoring station direction cosine as first preprocessing data in real time, and taking the floating solution of the monitoring station coordinates and estimation variances thereof as second preprocessing data; wherein the satellite altitude is determined by the satellite position and the observation residual is determined by the residual array;
when the time for storing the first preprocessing data is accumulated to K hours, constructing a multi-epoch station star double-difference observation equation according to the first preprocessing data, and carrying out least square estimation and partial ambiguity fixation according to the multi-epoch station star double-difference observation equation to obtain a monitoring station coordinate hour scale fixed solution as low-frequency deformation monitoring information; k is greater than or equal to 1;
and when the time for storing the second preprocessing data is accumulated to one day, taking the floating solution of the monitoring station coordinates and the estimated value variance thereof as an observed value and a weight inverse array thereof respectively, and estimating the monitoring station coordinates and the day solution by utilizing least square according to the observed value and the weight inverse array thereof to be taken as long-trend deformation monitoring information.
2. The beidou/GNSS based multi-time scale deformation monitoring method according to claim 1, wherein the partial ambiguity fixing is performed by an LAMBDA algorithm; the instantaneous fixed solution of the coordinates of the monitoring station is obtained through epoch-by-epoch real-time solution, and the precision is in the centimeter level and is used for representing burst large deformation.
3. The Beidou/GNSS multi-time scale deformation monitoring method of claim 1, wherein the multi-epoch station star double difference observation equation is constructed by the following steps:
for a plurality of continuous epochs, caching first preprocessing data of each epoch in real time;
under the condition that first preprocessing data of N epochs are acquired and accumulated for K hours, a reference star is selected, the number of ambiguity to be estimated is determined, and a multi-epoch station star double difference observation equation is constructed based on the reference star, the number of ambiguity to be estimated and the first preprocessing data of the N epochs; wherein N is a preset epoch number threshold.
4. The beidou/GNSS based multi-time scale deformation monitoring method of claim 1, wherein when the time of storing the second preprocessing data is accumulated to one day, the floating solution of the monitoring station coordinates and the estimated variance thereof are respectively used as an observed value and an inverse matrix thereof, and the least square is used to estimate the monitoring station coordinates day solution according to the observed value and the inverse matrix thereof, which comprises:
for a plurality of continuous epochs, caching second preprocessing data of each epoch in real time;
under the condition that second preprocessing data of M epochs are acquired, constructing a baseline solution quadratic adjustment based on the observed values of the M epochs and the weight inverse arrays thereof, and solving the baseline solution quadratic adjustment by using least square estimation to obtain a monitoring station coordinate day solution; wherein M is the minimum epoch number of the represented time which is more than or equal to one day.
5. The multi-time scale deformation monitoring method based on Beidou/GNSS of any one of claims 1-4, wherein the cycle slip and gross error detection specifically comprises:
removing coarse difference data contained in the Beidou/GNSS original observation data;
and marking and initializing cycle slip data contained in the Beidou/GNSS original observation data.
6. The Beidou/GNSS based multi-time scale deformation monitoring method of any one of claims 1-4, wherein the ionosphere weighted inter-station single difference equation is determined by the following formulas (1) - (2):
wherein a represents a reference station and b represents a monitoring station; i represents an epoch number, j represents a frequency number;a pseudo-range observation representing a single difference between stations corresponding to satellites with a satellite number s, +.>Representing the phase observations corresponding to satellites with the satellite number s; />Representing an inter-site Shan Chawei ground clearance approximation; />Representing the single difference receiver clock difference between stations; />Representing zenithal troposphere wet delay projection function, τ ab (i) Represents zenith tropospheric wet delay; />(i) Represents the single difference ionospheric bias delay, mu, between stations j ＝(λ j /λ 1 ) 2 A scale factor representing ionospheric bias delay, where lambda j The phase wavelength of the j-th frequency point is represented;representing the code deviation of a single difference receiver between stations; />Representing the phase deviation of the single-difference receiver between stations; />A station star double difference phase ambiguity relative to the reference station a and the first reference star r; /> Representing single-difference pseudorange observation noise between stations,representing inter-station single-difference phase observation noise; />The direction cosine from the monitoring station to the satellite is represented; dx (dx) b (i)＝[dx b (i) dy b (i) dz b (i)] T Representing the coordinate corrections of the monitoring station X, Y, Z in three directions; dt (dt) ab (i) Representing the single difference receiver clock difference between the original stations; d, d ab,1 Representing code deviation of single-difference receiver between base frequency point stations, d ab,j Representing code deviation of a single-difference receiver between j frequency point stations; delta ab,j Representing phase deviation of a single-difference receiver between j frequency point stations; />Representing the single differential phase ambiguity between the stations of the reference star r.
7. The multi-time scale deformation monitoring method according to any one of claims 1-4, wherein the multi-epoch station star double difference observation equation is determined by the following formula (3):
wherein a represents a reference station and b represents a monitoring station; i represents an epoch number, j represents a frequency number; q represents the second reference star, and,coordinate correction for three directions representing the monitoring station X, Y, Z hour arc solution, +.>Representing the single difference phase ambiguity between stations relative to reference station a and second reference star q, < >>A pseudo-range observation representing a single difference between stations corresponding to satellites with a satellite number s, +.>Representing single difference pseudorange observations between stations corresponding to a second reference star,representing satellitesSingle-difference phase observation value between stations corresponding to satellites with number s, < >>Representing the single difference phase observation value between stations corresponding to the second reference star, < >>The directional cosine of the monitoring station b to the satellite s; lambda (lambda) j Representing the phase wavelength.
8. The beidou/GNSS based multi-time scale deformation monitoring method of claim 4 wherein the baseline solution quadratic adjustment is determined by the following formula (4):
wherein 1= [1 1 1 1]Representing an all 1 vector;the coordinate corrections representing the three directions of the monitoring station X, Y, Z are solved for days; i denotes an epoch number.
9. The multi-time scale deformation monitoring method based on beidou/GNSS as claimed in any of the claims 1-4, wherein said quality control comprises the steps of:
performing residual size check on the data obtained after the processing of the Kalman filtering algorithm;
if the residual error exceeds a preset residual error threshold value, carrying out weight reduction or elimination processing on an observed value with the largest residual error in the data, carrying out Kalman filtering processing on the data subjected to the weight reduction or elimination processing, and repeating the steps until the residual error is smaller than the preset residual error threshold value, and outputting the data as the floating solution.
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202311517400.5A CN117739797B (en) | 2023-11-15 | 2023-11-15 | Beidou/GNSS-based multi-time-scale deformation monitoring method |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202311517400.5A CN117739797B (en) | 2023-11-15 | 2023-11-15 | Beidou/GNSS-based multi-time-scale deformation monitoring method |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN117739797A true CN117739797A (en) | 2024-03-22 |
| CN117739797B CN117739797B (en) | 2024-07-12 |
Family
ID=90282103
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202311517400.5A Active CN117739797B (en) | 2023-11-15 | 2023-11-15 | Beidou/GNSS-based multi-time-scale deformation monitoring method |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN117739797B (en) |
Citations (8)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN108981559A (en) * | 2018-08-28 | 2018-12-11 | 郑州信大先进技术研究院 | Real-time deformation monitoring method and system based on Beidou ground strengthening system |
| CN109059751A (en) * | 2018-09-10 | 2018-12-21 | 中国科学院国家授时中心 | A kind of deformation data monitoring method and system |
| WO2021237804A1 (en) * | 2020-05-29 | 2021-12-02 | 湖南联智科技股份有限公司 | Infrastructure structure deformation monitoring method based on beidou high-precision positioning |
| CN114912551A (en) * | 2022-07-18 | 2022-08-16 | 中国铁路设计集团有限公司 | GNSS and accelerometer real-time fusion algorithm for bridge deformation monitoring |
| CN115235330A (en) * | 2022-07-27 | 2022-10-25 | 长安大学 | Deformation monitoring data processing method and deformation monitoring device |
| CN116256774A (en) * | 2023-03-24 | 2023-06-13 | 淮北矿业股份有限公司 | Mining area railway line deformation monitoring method based on Beidou satellite positioning |
| CN116336931A (en) * | 2022-10-17 | 2023-06-27 | 南京凌远时空科技有限公司 | GNSS-based real-time dam deformation monitoring method and system |
| CN116718153A (en) * | 2023-08-07 | 2023-09-08 | 成都云智北斗科技有限公司 | Deformation monitoring method and system based on GNSS and INS |
-
2023
- 2023-11-15 CN CN202311517400.5A patent/CN117739797B/en active Active
Patent Citations (8)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN108981559A (en) * | 2018-08-28 | 2018-12-11 | 郑州信大先进技术研究院 | Real-time deformation monitoring method and system based on Beidou ground strengthening system |
| CN109059751A (en) * | 2018-09-10 | 2018-12-21 | 中国科学院国家授时中心 | A kind of deformation data monitoring method and system |
| WO2021237804A1 (en) * | 2020-05-29 | 2021-12-02 | 湖南联智科技股份有限公司 | Infrastructure structure deformation monitoring method based on beidou high-precision positioning |
| CN114912551A (en) * | 2022-07-18 | 2022-08-16 | 中国铁路设计集团有限公司 | GNSS and accelerometer real-time fusion algorithm for bridge deformation monitoring |
| CN115235330A (en) * | 2022-07-27 | 2022-10-25 | 长安大学 | Deformation monitoring data processing method and deformation monitoring device |
| CN116336931A (en) * | 2022-10-17 | 2023-06-27 | 南京凌远时空科技有限公司 | GNSS-based real-time dam deformation monitoring method and system |
| CN116256774A (en) * | 2023-03-24 | 2023-06-13 | 淮北矿业股份有限公司 | Mining area railway line deformation monitoring method based on Beidou satellite positioning |
| CN116718153A (en) * | 2023-08-07 | 2023-09-08 | 成都云智北斗科技有限公司 | Deformation monitoring method and system based on GNSS and INS |
Non-Patent Citations (1)
| Title |
|---|
| 冯威;黄丁发;李萌;张熙;严丽;: "高频GPS双差残差模型监测强震地表运动", 地球物理学报, no. 09, 15 September 2013 (2013-09-15) * |
Also Published As
| Publication number | Publication date |
|---|---|
| CN117739797B (en) | 2024-07-12 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN108317949B (en) | RTK high-precision differential positioning deformation monitoring system and method | |
| CN108981559B (en) | Real-time deformation monitoring method and system based on Beidou foundation enhancement system | |
| Bilitza et al. | The international reference ionosphere today and in the future | |
| DE102011076604A1 (en) | GNSS Atmospheric Estimation with Ambiguity Mooring | |
| CN108196284B (en) | GNSS network data processing method for fixing single-difference ambiguity between satellites | |
| CN113358017B (en) | Multi-station cooperative processing GNSS high-precision deformation monitoring method | |
| Baldysz et al. | Comparison of GPS tropospheric delays derived from two consecutive EPN reprocessing campaigns from the point of view of climate monitoring | |
| CN113253314B (en) | Time synchronization method and system between low-orbit satellites | |
| WO2023197714A1 (en) | Gnss multi-path error reducing method suitable for dynamic carrier platform | |
| CN103592657A (en) | Method for realizing single-mode RAIM (Receiver Autonomous Integrity Monitoring) under small number of visible satellites based on assistance of clock correction | |
| CN111735380A (en) | Method for extracting dynamic deflection of high-speed rail bridge in real time by using accelerometer to assist GNSS | |
| CN114002727A (en) | Differential speed measurement method, system, device and medium between non-combined carrier phase epochs | |
| Urschl et al. | Validating ocean tide loading models using GPS | |
| CN117194928B (en) | GNSS-based geographic deformation monitoring system | |
| CN117739797B (en) | Beidou/GNSS-based multi-time-scale deformation monitoring method | |
| CN115980317B (en) | Foundation GNSS-R data soil moisture estimation method based on corrected phase | |
| Tang et al. | GNSS localization propagation error estimation considering environmental conditions | |
| Wang et al. | Analysis of GNSS-R Code-Level Altimetry using QZSS C/A, L1C, and BDS B1C signals and their Combinations in a Coastal Experiment | |
| Li et al. | Study on Precise Point Positioning based on combined GPS and GLONASS | |
| Tserolas et al. | The Western Crete geodetic infrastructure: Long-range power-law correlations in GPS time series using Detrended Fluctuation Analysis | |
| Kashani et al. | The impact of the ionospheric correction latency on long-baseline instantaneous kinematic GPS positioning | |
| CN113324627A (en) | Buoy-based sea level observation and prediction method | |
| Marques et al. | Shoreline monitoring by gnss-ppp aiming to attendance the law 14.258/2010 from Pernambuco state, Brazil | |
| Enjolras et al. | Using altimetry waveform data and ancillary information from SRTM, Landsat, and MODIS to retrieve river characteristics | |
| Malinowski | Accuracy of Precise Point Positioning with GPS and GLONASS Satellites Constellation Using Online Web Services |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| GR01 | Patent grant |