GNSS landslide geological disaster monitoring and analysis method based on virtual reference station

Abstract

A GNSS landslide geological disaster monitoring and analysis method based on a virtual reference station (VRS) includes: Step 1, continuously operating reference stations determination, namely continu

Description

TECHNICAL FIELD

[0001] The present invention relates to the technical field of remote sensing disaster identification, and in particular, to a GNSS landslide geological disaster monitoring and analysis method based on a virtual reference station (VRS). BACKGROUND

[0002] Landslide geological disaster monitoring mainly relies on mass monitoring and prevention as well as professional monitoring. In terms of professional monitoring, it primarily utilizes Global Navigation Satellite System (GNSS) technology and sensor technologies, such as crack meters, tension meters, and rain gauges. GNSS technology offers advantages such as all-weather capability, high precision, automation, and no need for inter-visibility. The relative baseline solution precision can reach 10-9, making it an important technical means for landslide monitoring and early warning.

[0003] At present, GNSS landslide monitoring mainly adopts the "1+N" mode, that is, 1 reference station and N monitoring stations. By solving the baselines between the reference station located in a geologically stable area and the monitoring stations within a landslide body, the displacement of the landslide body can be obtained. However, the “1+N” mode requires the construction of at least one reference station near each monitoring area, leading to problems such as high construction costs, difficulty in site selection, unstable benchmarks, and resource waste. SUMMARY

[0004] The technical problem to be solved by the present invention is to provide a GNSS landslide geological disaster monitoring and analysis method based on a VRS, so as to solve the problems in existing GNSS landslide monitoring such as high construction costs, difficulty in site selection, unstable benchmarks, and resource waste.

[0005] To solve the above technical problem, the present invention adopts the following technical solution:

[0006] A GNSS landslide geological disaster monitoring and analysis method based on a VRS comprises the following steps:

[0007] Step 1, continuously operating reference station (CORS) network determination and site selection of monitoring stations: comprising determining CORS network coverage in a project area, and selecting landslide feature monitoring stations to be monitored within a CORS coverage area for the installation of monitoring equipment;

[0008] Step2, solution of a multi-system VRS network: based on the VRS, the key to landslide monitoring is to generate effective VRS observation values near the monitoring stations and perform short baseline solutions between the VRS observation values and monitoring station observation values to obtain virtual reference station observation values of the monitoring stations;

[0009] Step 3, quality inspection of the virtual observation values: inspecting and evaluating the quality of the virtual observation values obtained by the VRS by means of physical reference station data; and

[0010] Step 4, data processing and result analysis: solving the observation values with different durations, statistically calculating accuracy, and statistically evaluating the accuracy of solution results in north (N), east (E), and up (U) directions.

[0011] Step 1, the CORS network determination and site selection of monitoring stations comprises determining the coverage area of the CORS network in the project area, and determining the site selection of the landslide feature monitoring stations.

[0012] Step 2, a solution method of VRS observation values comprises three parts: firstly, solving baselines in a CORS network and fixing double-difference integer ambiguities; secondly, using BeiDou observation values of each baseline, the double-difference integer ambiguities, and accurate position coordinates of reference stations to calculate ionospheric delay, tropospheric delay, and a combined error of each baseline; and thirdly, determining differential correction information at the VRS through interpolation based on the coordinates of the VRS, so as to provide real-time virtual reference services for a geological disaster monitoring terminal.

[0013] A solution process for the integer ambiguities in the solution of the VRS observation values comprises:

[0014] 1) forming a Melboume-Wiibbena combination, namely an M-W combination, by means of a carrier-phase observation value and a P-code pseudo-range, and calculating a wide-lane integer ambiguity, wherein a specific formula as follows:

[0015] where superscripts p and q represent an observation satellite p and a reference satellite q, the subscript WL represents a GNSS wide-lane combined observation value, subscripts i and j represent reference stations i and j, X represents a carrier wavelength, and respectively represent frequencies of a first frequency point LI and a second frequency point L2 of GNSS observation values, AV represents a double-difference operator, N represents the integer AW pq = ambiguity, P represents a pseudo-range observation value, and (p represents a carrier phase observation value;

[0016] 2) performing an ionosphere-free combination of the carrier phase observation value and the pseudo-range observation value, and solving to obtain a real number solution of an ionosphere-free ambiguity, wherein an observation equation is as follows: 2 • ^(ppq = ^ppq - 2 • ^Npq + ^Troppq + £pq j y y y y y = + ^Troppq +

[0017] where p represents a distance from a monitoring station to a satellite, Trop represents a tropospheric delay error, e represents carrier-phase observation value noise, and y represents pseudo-range observation value noise;

[0018] Since the GLONASS uses a frequency-division multiple-access signal structure, a GLONASS double-difference ambiguity is not an integer in a double-difference observation equation in which the unit is distance, the GLONASS double-difference ambiguity needs to be converted into a double-difference integer ambiguity and a single-difference integer ambiguity, and the observation equation is as follows: U, ,pq =Wppq-2, -AV7V ,^+(2, -2, )-A2V / glop T gioIJ r IJ glop gloij v glop gloqg gloij + ^Troppq+spq ^Pglpq = ^ppq + ^Troppq + y™

[0019] where Aglo, <pglo, Pglo, and lVglo respectively represent a wavelength, a phase observation value, a pseudo-range observation value, and an integer ambiguity of GLONASS observation values respectively, / ? represents an observation satellite, and q represents a reference satellite;

[0020] 3) using the wide-lane integer ambiguity solved in Step 1), and the real number solution of the ionosphere-free ambiguity and a corresponding variance-covariance matrix solved in Step 2) to obtain a real number solution of a narrow-lane ambiguity and a corresponding variance-covariance matrix, and using a least-square ambiguity decorrelation adjustment (LAMBDA) algorithm to fix a narrow-lane integer phase ambiguity, and then calculating integer phase ambiguities of GNSS satellites at LI and L2 frequencies.

[0021] The calculation process for the ionospheric delay, the tropospheric delay and the combined error comprises:

[0022] during the calculation of errors of a multi-system VRS, dividing the errors into ionospheric delay, tropospheric delay and a combined error including second-order ionospheric delay, an orbit error, a multipath error, etc., and performing calculations separately; after fixing the integer ambiguity of each baseline, calculating the above three errors based on the GNSS observation value of each baseline, the double-difference integer ambiguities and the accurate position coordinates of the reference stations, wherein specific formulas are as follows: ^lono = (—^_)[(^AV^ - ^AW2) + (4AWj - ^AW2)] / 1 / 2 ^Trop = (MFw(ekm) -MFw «)) • ZWDm - (MFv -MFv• ZWDn + (MFh(4) ~MFh «)) • ZHDm - (MF„ (ek) —MFH «)) • ZHDn AV Other = AVp -(^ AV A VA^) - Milano - Ml Trop

[0023] whereAV / ono represents double-difference ionospheric delay, represents a carrier observation value of the frequency point LI, represents a carrier observation value of the frequency point L2, and MF represents a tropospheric mapping function; ZWD represents tropospheric wet delay; ZHD represents tropospheric hydrostatic delay; Other represents a combined error; and ^2 represent wavelengths of LI and L2 observation values, and ^2 represent integer ambiguities of the LI and L2 observation values, e represents an elevation angle of an observation value, AV represents a double-difference operator, subscripts H and W represent tropospheric dry delay and tropospheric wet delay respectively, subscripts n and m represent a reference station and a monitoring station respectively, and superscripts s and k represent a reference satellite and an observation satellite respectively.

[0024] A specific process for determining differential correction information at the VRS through interpolation based on the coordinates of the VRS is as follows:

[0025] For the monitoring terminal of the multi-system VRS network, the key to achieving real-time accurate positioning lies in whether distance-related errors can be accurately obtained. Due to the strong spatio-temporal correlation of an ionospheric delay error, a tropospheric delay error, and a combined error, and considering that a vertical distance between a rover station and the reference station is much smaller than a horizontal distance, an inverse distance weighting method is used to interpolate the three types of errors: the ionospheric delay error, the tropospheric delay error, and the combined error in the east and north direction vectors of the baselines separately. Observation variances corresponding to interpolation coefficients aA and a2 according to the east and north direction sectors, respectively, are as follows: L=B■X+^

[0026] where L represents an observation value matrix, B represents a design matrix, X represents a to-be-estimated parameter matrix, and △ represents a residual matrix, wherein a weight matrix P of the observation equation is as follows: Err Err, L = . 2 _Errn_ ; ’AXj Aff B= A72 _^n EYn_, x= 1 a2

[0027] AX, and A^ represent eastward and northward distances from each point on baseline 1 in a topocentric coordinate system to a reference point, AX2 and Ay2 represent eastward and northward distances from each point on baseline 2 in the topocentric coordinate system to the reference point, AXn and A / n represent eastward and northward distances from each point on baseline n in the topocentric coordinate system to the reference point;

[0028] the interpolation coefficients^ and a2 corresponding to the east and north direction vectors can be solved based on the least-squares principle, and the errors of the rover station are then interpolated based on the following formula: Errv

[0029] where Err represents a baseline error, AX and AK represent the east and north direction vectors of the baselines, D represents a baseline length, ar and u2 represent the interpolation coefficients corresponding to the east and north direction vectors, Errv. AXV and A / v represent an interpolation error of the rover station, and the east and north direction vectors of baselines from the rover station to a master station.

[0030] Step 3, quality inspection of observation values: inspecting and evaluating the quality of the VRS observation values by means of physical reference station data comprises three steps: firstly, selecting same-named monitoring stations and determining inspection indicators; secondly, using GNSS preprocessing software Anubis to perform preprocessing and quality inspections on VRS data and physical reference station data of the same-named monitoring stations; and finally, statistically analyzing data quality of each indicator under both monitoring modes to evaluate the VRS data quality.

[0031] The selecting same-named monitoring stations and determining inspection indicators comprises determining quality inspection indicators for the VRS data and the physical reference station data of the same-named monitoring stations, where the inspection indicators comprise a signal-to-noise ratio (SNR), a multipath effect, and a cycle slip ratio.

[0032] The method further comprises performing preprocessing and quality inspections on the VRS data and the physical reference station data of the same-named monitoring stations respectively, specifically using the GNSS preprocessing software Anubis to process the VRS data and the physical reference station data of the same-named monitoring stations and to determine values of the quality inspection indicators.

[0033] The method further comprises statistically analyzing and evaluating the VRS data quality, specifically performing comparative analysis on values of the quality inspection indicators such as the SNR, the multipath effect and the cycle slip ratio of the VRS data and the physical reference station data of the same-named monitoring stations, and evaluating feasibility of VRS solution data and determining the data quality.

[0034] According to the GNSS landslide geological disaster monitoring and analysis method based on the VRS provided by the present invention, by means of the landslide geological disaster monitoring method based on the VRS, the CORSs for satellite navigation and positioning are utilized to generate VRS data in a monitoring area, replacing a physical reference station for landslide monitoring to reduce the cost of GNSS landslide monitoring. Meanwhile, the VRS is not constrained by geological conditions and distances, thus meeting the requirements for wide-range coverage, universality, and adaptability to complex scenarios in landslide monitoring. DESCRIPTION OF THE DRAWINGS

[0035] The present invention is further described with reference to the accompanying drawings and embodiments:

[0036] FIG. 1 is a flow chart of a method according to the present invention;

[0037] FIG. 2 is a schematic diagram of coordinates of a VRS and a corresponding monitoring station according to an embodiment; and

[0038] FIG. 3 shows data solution results in an embodiment. DETAILED DESCRIPTION OF EMBODIMENTS

[0039] The technical solution of the present invention will be described in detail in combination with the accompanying drawings and the specific embodiments.

[0040] A GNSS landslide geological disaster monitoring and analysis method based on a VRS includes the following steps:

[0041] Step 1, CORS network determination and site selection of monitoring stations: including determining CORS network coverage in a project area, and selecting landslide feature monitoring stations to be monitored within a CORS coverage area for the installation of monitoring equipment, wherein CORS refers to continuously operating reference stations;

[0042] Step2, solution of a multi-system VRS network: based on the VRS, the key to landslide monitoring is to generate an effective VRS near the monitoring stations and perform short baseline solutions between the VRS and the monitoring stations to obtain displacement of the monitoring stations;

[0043] Step 3, quality inspection of virtual observation values: inspecting and evaluating the quality of VRS observation values by means of physical reference station data: firstly, selecting same-named monitoring stations and determining inspection indicators; secondly, using GNSS preprocessing software Anubis to perform preprocessing and quality inspections on VRS data and physical reference station data of the same-named monitoring stations; and finally, statistically analyzing data quality of each indicator under both monitoring modes to evaluate the VRS data quality.

[0044] Step 4, data processing and result analysis: solving the observation values with different durations, statistically calculating accuracy, and statistically evaluating the accuracy of solution results in north (N), east (E), and up (U) directions.

[0045] Step 1, the CORS network determination and site selection of monitoring stations includes determining the coverage area of the CORS network in the project area, and determining the site selection of the landslide feature monitoring stations.

[0046] Step 2 (which is completed by a computer), a solution method of VRS observation values includes three parts: firstly, solving baselines in a BeiDou reference station network and fixing double-difference integer ambiguities; secondly, using BeiDou observation values of each baseline, the double-difference integer ambiguities, and accurate position coordinates of reference stations to calculate ionospheric delay, tropospheric delay, and a combined error of each baseline; and thirdly, determining differential correction information at the VRS through interpolation based on the coordinates of the VRS, so as to provide high-accuracy real-time virtual reference services for a geological disaster monitoring terminal.

[0047] A solution process for the integer ambiguities in the solution of the VRS observation values includes: 1) 1) forming a Melbourne-Wtibbena combination, namely the M-W combination, by means of AV? / pq = a carrier-phase observation value and a P-code pseudo-range, and calculating a wide-lane integer ambiguity, where a specific formula as follows: y+ / 2 A-n

[0048] where superscripts p and q represent an observation satellite p and a reference satellite q, the subscript WL represents a GNSS wide-lane combined observation value, subscripts i and j represent reference stations i and j, X represents a carrier wavelength, and ^2 respectively represent frequencies of a first frequency point LI and a second frequency point L2 of GNSS observation values, AV represents a double-difference operator, N represents the integer ambiguity, P represents a pseudo-range observation value, and (p represents a carrier phase observation value;

[0049] 2) performing an ionosphere-free combination of the carrier phase observation value and the pseudo-range observation value, and solving to obtain a real number solution of an ionosphere-free ambiguity, where an observation equation is as follows: A. = AVpf - A • AW” + WTrop™ + s™ WPf = AVpf + AV Prop™ +

[0050] where p represents a distance from a monitoring station to a satellite, Trop represents a tropospheric delay error, e represents carrier-phase observation value noise, and y represents pseudo-range observation value noise;

[0051] Since the GLONASS uses a frequency-division multiple-access signal structure, a GLONASS double-difference ambiguity is not an integer in a double-difference observation equation in which the unit is distance, the GLONASS double-difference ambiguity needs to be converted into a double-difference integer ambiguity and a single-difference integer ambiguity, and therefore, the observation equation is as follows: -AV^; M =AVpM-2; -A, )-AN,q glop ^gloy glop glop v glop gloq >gloij + AVTroppq+spq = Wp + Wlnp +

[0052] where Aglo, <pglo, Pglo, and Aglo respectively represent a wavelength, a phase observation value, a pseudo-range observation value, and an integer ambiguity of GLONASS observation values respectively, / ? represents an observation satellite, and q represents a reference satellite;

[0053] 3) using the wide-lane integer ambiguity solved in Step 1), and the real number solution of the ionosphere-free ambiguity and a corresponding variance-covariance matrix solved in Step 2) to obtain a real number solution of a narrow-lane ambiguity and a corresponding variance-covariance matrix, and using a least-square ambiguity decorrelation adjustment (LAMBDA) algorithm to fix a narrow-lane integer phase ambiguity, and then calculating integer phase ambiguities at LI and L2 frequencies.

[0054] The calculation process (which is completed by a computer) of the ionospheric delay, the tropospheric delay and the combined error includes:

[0055] during the calculation of errors of the multi-system VRS, dividing the errors into ionospheric delay, tropospheric delay and a combined error including second-order ionospheric delay, an orbit error, a multipath error, etc., and performing calculations separately; after fixing the integer ambiguity of each baseline, calculating the above three errors based on the GNSS observation value of each baseline, the double-difference integer ambiguities and the accurate position coordinates of the reference stations, where specific formulas are as follows: f2 AV / owo = ( / 2 2)[(^AV^ - ^AV^) + (4AVAj -^AWJ] f ~ fi ^Trop = (MFW(ekm) -MFW «)) • ZWDm - (MFW (ef-MFW«)) • ZWD„ + WFH(4) ~MFh «)) • ZHDm - (MFH (ek) —MFH «)) • ZHDn AV Other = AVp -(^ AV AW,) - £NIono - Trop

[0056] whereAV / ono represents double-difference ionospheric delay, represents a carrier observation value of the frequency point LI, represents a carrier observation value of the frequency point L2, and MF represents a tropospheric mapping function; ZWD represents tropospheric wet delay; ZHD represents tropospheric hydrostatic delay; Other represents a combined error; and ^2 represent wavelengths of LI and L2 observation values, and ^2 represent integer ambiguities of the LI and L2 observation values, e represents an elevation angle of an observation value, AV represents a double-difference operator, subscripts H and W represent tropospheric dry delay and tropospheric wet delay respectively, subscripts n and m represent a reference station and a monitoring station respectively, and superscripts s and k represent a reference satellite and an observation satellite respectively.

[0057] The specific process (which is completed by a computer) of determining differential correction information at the VRS through interpolation based on the coordinates of the VRS is as follows:

[0058] For the monitoring terminal of the multi-system VRS network, the key to achieving real-time accurate positioning lies in whether distance-related errors can be accurately obtained. Due to the strong spatio-temporal correlation of an ionospheric delay error, a tropospheric delay error, and a combined error, and considering that a vertical distance between a rover station and the reference station is much smaller than a horizontal distance, an inverse distance weighting method is used to interpolate the three types of errors: the ionospheric delay error, the tropospheric delay error, and the combined error in the east and north direction vectors of the baselines separately. Observation variances corresponding to interpolation coefficients aA and a2 according to the east and north direction sectors, respectively, are as follows: L=B-X+E

[0059] where L represents an observation value matrix, B represents a design matrix, X represents a to-be-estimated parameter matrix, and △ represents a residual matrix, wherein a weight matrix P of the observation equation is as follows: _1 D 0 0 • Er^ Err, L = . 2 _Errn EX, ex. AX EY n n a, X = 1 a2

[0060] AXt and A^ represent eastward and northward distances from each point on baseline 1 in a topocentric coordinate system to a reference point, AX2 and A / 2 represent eastward and northward distances from each point on baseline 2 in the topocentric coordinate system to the reference point, EXn and EYn represent eastward and northward distances from each point on baseline n in the topocentric coordinate system to the reference point;

[0061] the interpolation coefficientsu! and a2 corresponding to the east and north direction vectors can be solved based on the least-squares principle, and the errors of the rover station are then interpolated based on the following formula: Errv=ai-!EXv+a2-^v

[0062] where Err represents a baseline error, AX and AT represent the east and north direction vectors of the baselines, D represents a baseline length, aA and a2 represent the interpolation coefficients corresponding to the east and north direction vectors, Errv, AXV and ATV represent an interpolation error of the rover station, and the east and north direction vectors of baselines from the rover station to a master station.

[0063] The method is based on the network real-time kinematic (NRTK) technology of the CORS, and is capable of generating a VRS at any location within the coverage area of a reference station network. The VRS is not constrained by geological conditions or distance, thereby reducing the cost of GNSS landslide monitoring, enabling wide-range landslide monitoring, and meeting the requirements for wide coverage, universality, and adaptability to complex scenarios in landslide monitoring.

[0064] Step 3, quality inspection of observation values: inspecting and evaluating the quality of VRS observation values by means of physical reference station data includes three steps: firstly, selecting same-named monitoring stations and determining inspection indicators; secondly, using GNSS preprocessing software Anubis to perform preprocessing and quality inspections on VRS data and physical reference station data of the same-named monitoring stations; and finally, statistically analyzing data quality of each indicator under both monitoring modes to evaluate the VRS data quality.

[0065] The selecting same-named monitoring stations and determining inspection indicators includes determining quality inspection indicators of the VRS data and the physical reference station data of the same-named monitoring stations, where the inspection indicators include a signal-to-noise ratio (SNR), a multipath effect, and a cycle slip ratio.

[0066] The method further includes performing preprocessing and quality inspections on the VRS data and the physical reference station data of the same-named monitoring stations respectively, specifically using the GNSS preprocessing software Anubis to process the VRS data and the physical reference station data of the same-named monitoring stations and to determine values of the quality inspection indicators.

[0067] The method further comprises statistically analyzing and evaluating the VRS data quality, specifically performing comparative analysis on values of the quality inspection indicators such as the SNR, the multipath effect and the cycle slip ratio of the VRS data and the physical reference station data of the same-named monitoring stations, and with reference to the technical specification “Technical Specification for Reference Station Construction and Acceptance of BeiDou Ground-Based Augmentation System (GBT 39772.2-2021) Part 2: Acceptance Specification”, evaluating feasibility of VRS solution data and determining the data quality.

[0068] Step 4, data processing and result analysis: solving the observation values with different durations, statistically calculating accuracy, and statistically evaluating the accuracy of solution results in north (N), east (E), and up (U) directions.

[0069] Embodiment:

[0070] To verify the present invention, data from an XX project is selected as an illustrative example. In this project, station JD38N is selected as a geological disaster monitoring station. Differential data from a reference station is obtained by solving both a VRS and a traditional physical reference station, respectively, and the monitoring data are compared. A schematic diagram of the baselines is shown in FIG. 2.

[0071] The data solved from 16:00:00 on Sep. 16, 2022, to 16:00:00 on Sep. 18, 2022, is selected for accuracy comparison and analysis. A resolution strategy is set as dual-system (GPS + BDS) and multi-frequency differential data processed by using a 4-hour Kalman filter, with a solution interval of 5 min. Coordinate curves of the sets of solved results in the east, north, and up directions are shown in FIG. 3.

[0072] By calculating standard deviations (STD) of coordinate variations in three-dimensional directions for both sets of results, the internal consistency accuracy of positioning results obtained from a VRS mode and a physical reference station mode can be obtained. The statistical results are as follows: Physical reference station STD( / m) VRS STD( / m) East 0.0216 0.0179 North 0.0124 0.0093 Up 0.0597 0.0455

[0073] To verify the external consistency accuracy of the two sets of results, this project uses a single-day static baseline solution as a true value of a baseline vector, and the external consistency accuracy of the two sets of VRS positioning results is calculated based on this. The statistical results are shown as follows: Physical reference station RMS( / m) VRS RMS( / m) East 0.0236 0.0182 North 0.0124 0.0096 Up 0.0621 0.0465

[0074] The present invention aims to remotely monitor landslide displacement and trends at landslide sites through sensors and satellite navigation wireless signal transmission. The smaller the difference between a remotely calculated accuracy value and a monitoring value from an on-site physical reference station, the higher the accuracy. The comparison results shown in FIG. 3 present the fitting between the actual coordinate deviations in the three-dimensional directions (east, north, and up) and the variances calculated by this application. It can be seen from the data in the above chart that the results of deformation monitoring of the differential data in the two modes both have centimeter-level accuracy.

[0075] The actual on-site environmental conditions are as follows:

[0076] (1) The baselines are relatively long. In this experiment, the selected monitoring station is approximately 5 km away from the VRS and has a significant height difference.

[0077] (2) Both modes only support the GPS and the Beidou-2 system when generating differential correction data, and do not accommodate satellites of the Beidou-3 system for the time being, resulting in a relatively small number of available satellites during the solution.

[0078] From the above, it can be seen that the currently used method can still achieve centimeter-level accuracy under the relatively poor actual field conditions (1. long baselines and large height difference; 2. few satellites), which indicates that the method of the present invention has good practicability and easy popularization, and the method has good performance.

[0079] In addition, the above results also show that the monitoring accuracy of the VRS mode is close to that of the physical reference station mode. Under the VRS mode, meteorological constraint conditions are incorporated. Atmospheric errors including ionospheric delay and tropospheric delay have significant spatial correlations. Due to the access of AK25 station data, atmospheric delay information with stronger spatial correlation is fused in the construction process for an atmospheric error model, which greatly enhances the accuracy of atmospheric error correction values of the network rtk, and further, the accuracy of the VRS is already very close to that of the physical reference station mode.

Claims

1. A GNSS landslide geological disaster monitoring and analysis method based on a VRS, comprising the following steps:Step 1, CORS network determination and site selection of monitoring stations: comprising determining CORS network coverage in a project area, and selecting landslide feature monitoring stations to be monitored within a CORS coverage area for the installation of monitoring equipment;Step 2, solution of a multi-system VRS network: generating an effective VRS near the monitoring stations, and performing short baseline solutions between the VRS and the monitoring stations to obtain virtual observation values of the monitoring stations;Step 3, quality inspection of the virtual observation values: inspecting and evaluating the quality of the virtual observation values obtained by the VRS by means of physical reference station data; andStep 4, data processing and result analysis: solving the virtual observation values obtained by the VRS with different durations, statistically calculating accuracy, and statistically evaluating the accuracy of solution results in north (N), east (E), and up (U) directions.

2. The GNSS landslide geological disaster monitoring and analysis method based on a VRS according to claim 1, wherein in Step 1, the CORS network determination and site selection of monitoring stations comprises determining a CORS network coverage area in the project area, and determining the site selection of the landslide feature monitoring stations.

3. The GNSS landslide geological disaster monitoring and analysis method based on a VRS according to claim 1, wherein in Step 2, a solution method of VRS observation values comprises three parts: firstly, solving baselines in a CORS network and fixing double-difference integer ambiguities; secondly, using BeiDou observation values of each baseline, the double-difference integer ambiguities, and accurate position coordinates of reference stations to calculate ionospheric delay, tropospheric delay, and a combined error of each baseline; and thirdly, determining differential correction information at the VRS through interpolation based on the coordinates of the VRS, so as to provide real-time virtual reference services for a geological disaster monitoring terminal.

4. The GNSS landslide geological disaster monitoring and analysis method based on a VRS according to claim 3, wherein a solution process for the integer ambiguities in the solution of the VRS observation values comprises:1) forming a Melboume-Wiibbena combination, namely an M-W combination, by means of a carrier-phase observation value and a P-code pseudo-range, and calculating a wide-lane integer ambiguity, wherein a specific formula is as follows:\ / ;-av / ;"+ / 2-avp,” / -Av^-y-Av^^ y+ / 2 A-n5where superscripts p and q represent an observation satellite p and a reference satellite q, subscript WL represents a GNSS wide-lane combined observation value, subscripts i and j represent reference stations i and j, X represents a carrier wavelength, &and ^2 respectively represent frequencies of a first frequency point L1 and a second frequency point L2 of GNSS observation values, AV represents a double-difference operator, N represents the integer ambiguity, P represents a pseudo-range observation value, and <p represents a carrier phase observation value;2) performing an ionosphere-free combination of the carrier phase observation value and the pseudo-range observation value, and solving to obtain a real number solution of an ionosphere-free ambiguity, wherein an observation equation is as follows:2 • AVr / ?w = ^ppq - 2 • ^Npq + ^Troppq + spq j y y y y y^Ppq = + ^Troppq +5where / ^represents a distance from a monitoring station to a satellite, Trop represents a tropospheric delay error, e represents carrier-phase observation value noise, and y represents pseudo-range observation value noise;wherein since GLONASS uses a frequency-division multiple-access signal structure, a GLONASS double-difference ambiguity is not an integer in a double-difference observation equation in which the unit is distance, the GLONASS double-difference ambiguity needs to be converted into a double-difference integer ambiguity and a single-difference integer ambiguity, and thus an observation equation is as follows:f / l, -AW, M = AV■^N,pg+a, -2, YkN,gglop Tgloy ry glop gloij V glop gloqg gloij+ ^Troppg+£pq^Pglopg = ^ppg +^Troppg+ypg 5where Aglo, <pglo, Pglo, and Aglo represent a wavelength, a phase observation value, a pseudo-range observation value, and an integer ambiguity of GLONASS observation valuesrespectively, / ? represents an observation satellite, and q represents a reference satellite; and3) using the wide-lane integer ambiguity solved in Step 1), and the real number solution of the ionosphere-free ambiguity and a corresponding variance-covariance matrix solved in Step 2) to obtain a real number solution of a narrow-lane ambiguity and a corresponding variance-covariance matrix, and using a least-square ambiguity decorrelation adjustment (LAMBDA) algorithm to fix a narrow-lane integer phase ambiguity, and then calculating integer phase ambiguities of GNSS satellites at LI and L2 frequencies.

5. The GNSS landslide geological disaster monitoring and analysis method based on a VRS according to claim 4, wherein a calculation process for the ionospheric delay, the tropospheric delay and the combined error comprises:during the calculation of errors of a multi-system VRS, dividing the errors into ionospheric delay, tropospheric delay and a combined error comprising second-order ionospheric delay, an orbit error and a multipath error, and performing calculations separately; after fixing the integer ambiguity of each baseline, calculating the above three errors based on the GNSS observation value of each baseline, the double-difference integer ambiguities and the accurate position coordinates of the reference stations, wherein specific formulas are as follows:^lono = (■)[(^AV^ -^AV^j + ^AV^ -^AWJ]^Trop = (MFw(ekm)’?WDm -(MFw(ek)-MF^’?WDn + <MFh(e^)-MFH(esm))• ZHDm -(MFH(ek)-MFH«))• ZHDnAV Other = A Vp - (^ AV A Wj) - ^lono - AV TropwhereAV / ono represents double-difference ionospheric delay, represents a carrierobservation value of the frequency point LI, represents a carrier observation value of the frequency point L2, and MF represents a tropospheric mapping function; ZWD represents tropospheric wet delay; ZHD represents tropospheric hydrostatic delay; Other represents a combined error; and ^2 represent wavelengths of LI and L2 observation values, and ^2 represent integer ambiguities of the L1 and L2 observation values, e represents an elevation angle of an observation value, AV represents a double-difference operator, subscripts H and W represent tropospheric dry delay and tropospheric wet delay respectively, subscripts n and m represent a reference station and a monitoring station respectively, and superscripts s and k represent a reference satellite and an observation satellite respectively.

6. The GNSS landslide geological disaster monitoring and analysis method based on a VRS according to claim 5, wherein a specific process for determining differential correction information at the VRS through interpolation based on the coordinates of the virtual reference station comprises:using an inverse distance weighting method to interpolate three types of errors, namely an ionospheric delay error, a tropospheric delay error and a combined error, of a rover station according to eastward and northward direction vectors of the baselines, wherein observation variances a, and u2 °f interpolation coefficients corresponding to the east and north direction vectors are given as follows:L = B-X+k5where L represents an observation value matrix, B represents a design matrix, X represents a to-be-estimated parameter matrix, and △ represents a residual matrix, whereinII 3 l-j > >> ••• xc । ।aX = 1a2 5AX1 and A^ represent eastward and northward distances from each point on baseline 1 in a topocentric coordinate system to a reference point, AX2 and Ay2 represent eastward and northward distances from each point on baseline 2 in the topocentric coordinate system to the reference point, AXn and ATn represent eastward and northward distances from each point on baseline n in the topocentric coordinate system to the reference point; andsolving to obtain the interpolation coefficients^ and a2 corresponding to the east and north direction vectors based on a least-squares principle; and then interpolating the errors of the rover station based on the following formula:Errv5where Err represents a baseline error, AX and AT represent the east and north direction vectors of the baselines, D represents a baseline length, a, and a2 represent the interpolation coefficients corresponding to the east and north direction vectors, Errv, AXV and XYV represent an interpolation error of the rover station, and east and north direction vectors of baselines from the rover station to a master station.

7. The GNSS landslide geological disaster monitoring and analysis method based on a VRS according to claim 1, wherein the quality inspection of the virtual observation values: inspecting and evaluating the quality of the VRS observation values by means of physical reference station data in Step 3 comprises three steps: firstly, selecting same-named monitoring stations and determining inspection indicators; secondly, using GNSS preprocessing software Anubis to perform preprocessing and quality inspections on VRS data and physical reference station data of the same-named monitoring stations; and finally, statistically analyzing data quality of each indicator under both monitoring modes to evaluate the VRS data quality.

8. The GNSS landslide geological disaster monitoring and analysis method based on a VRS according to claim 7, wherein the selecting same-named monitoring stations and determining inspection indicators comprises determining quality inspection indicators for the VRS data and the physical reference station data of the same-named monitoring stations, wherein the inspection indicators comprise a signal-to-noise ratio (SNR), a multipath effect, and a cycle slip ratio.

9. The GNSS landslide geological disaster monitoring and analysis method based on a VRS according to claim 7, further comprising performing preprocessing and quality inspections on the VRS data and the physical reference station data of the same-named monitoring stations respectively, specifically using the GNSS preprocessing software Anubis to process the VRS data and the physical reference station data of the same-named monitoring stations and to determine values of the quality inspection indicators.

10. The GNSS landslide geological disaster monitoring and analysis method based on a VRS according to claim 7, further comprising: statistically analyzing and evaluating the VRS data quality: performing comparative analysis on values of the quality inspection indicators, i.e., the SNR, the multipath effect and the cycle slip ratio of the VRS data and the physical reference station data of the same-named monitoring stations, and evaluating feasibility of VRS solution data and determining the data quality.

Images