[go: up one dir, main page]

WO2020233591A1 - 一种面向三维地表形变估计的InSAR和GNSS定权方法 - Google Patents

一种面向三维地表形变估计的InSAR和GNSS定权方法 Download PDF

Info

Publication number
WO2020233591A1
WO2020233591A1 PCT/CN2020/091273 CN2020091273W WO2020233591A1 WO 2020233591 A1 WO2020233591 A1 WO 2020233591A1 CN 2020091273 W CN2020091273 W CN 2020091273W WO 2020233591 A1 WO2020233591 A1 WO 2020233591A1
Authority
WO
WIPO (PCT)
Prior art keywords
gnss
insar
data
observations
deformation
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.)
Ceased
Application number
PCT/CN2020/091273
Other languages
English (en)
French (fr)
Inventor
胡俊
刘计洪
李志伟
朱建军
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Central South University
Original Assignee
Central South University
Priority date (The priority date 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 date listed.)
Filing date
Publication date
Application filed by Central South University filed Critical Central South University
Priority to US17/035,771 priority Critical patent/US20210011149A1/en
Publication of WO2020233591A1 publication Critical patent/WO2020233591A1/zh
Anticipated expiration legal-status Critical
Ceased legal-status Critical Current

Links

Images

Classifications

    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S13/00—Systems using the reflection or reradiation of radio waves, e.g. radar systems; Analogous systems using reflection or reradiation of waves whose nature or wavelength is irrelevant or unspecified
    • G01S13/88—Radar or analogous systems specially adapted for specific applications
    • G01S13/89—Radar or analogous systems specially adapted for specific applications for mapping or imaging
    • G01S13/90—Radar or analogous systems specially adapted for specific applications for mapping or imaging using synthetic aperture techniques, e.g. synthetic aperture radar [SAR] techniques
    • G01S13/9021—SAR image post-processing techniques
    • G01S13/9023—SAR image post-processing techniques combined with interferometric techniques
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01B—MEASURING LENGTH, THICKNESS OR SIMILAR LINEAR DIMENSIONS; MEASURING ANGLES; MEASURING AREAS; MEASURING IRREGULARITIES OF SURFACES OR CONTOURS
    • G01B15/00—Measuring arrangements characterised by the use of electromagnetic waves or particle radiation, e.g. by the use of microwaves, X-rays, gamma rays or electrons
    • G01B15/06—Measuring arrangements characterised by the use of electromagnetic waves or particle radiation, e.g. by the use of microwaves, X-rays, gamma rays or electrons for measuring the deformation in a solid
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S13/00—Systems using the reflection or reradiation of radio waves, e.g. radar systems; Analogous systems using reflection or reradiation of waves whose nature or wavelength is irrelevant or unspecified
    • G01S13/88—Radar or analogous systems specially adapted for specific applications
    • G01S13/89—Radar or analogous systems specially adapted for specific applications for mapping or imaging
    • G01S13/90—Radar or analogous systems specially adapted for specific applications for mapping or imaging using synthetic aperture techniques, e.g. synthetic aperture radar [SAR] techniques
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S19/00—Satellite radio beacon positioning systems; Determining position, velocity or attitude using signals transmitted by such systems
    • G01S19/01—Satellite radio beacon positioning systems transmitting time-stamped messages, e.g. GPS [Global Positioning System], GLONASS [Global Orbiting Navigation Satellite System] or GALILEO
    • G01S19/13—Receivers
    • G01S19/14—Receivers specially adapted for specific applications
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S19/00—Satellite radio beacon positioning systems; Determining position, velocity or attitude using signals transmitted by such systems
    • G01S19/38—Determining a navigation solution using signals transmitted by a satellite radio beacon positioning system
    • G01S19/39—Determining a navigation solution using signals transmitted by a satellite radio beacon positioning system the satellite radio beacon positioning system transmitting time-stamped messages, e.g. GPS [Global Positioning System], GLONASS [Global Orbiting Navigation Satellite System] or GALILEO
    • G01S19/42—Determining position
    • G01S19/43—Determining position using carrier phase measurements, e.g. kinematic positioning; using long or short baseline interferometry
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S19/00—Satellite radio beacon positioning systems; Determining position, velocity or attitude using signals transmitted by such systems
    • G01S19/38—Determining a navigation solution using signals transmitted by a satellite radio beacon positioning system
    • G01S19/39—Determining a navigation solution using signals transmitted by a satellite radio beacon positioning system the satellite radio beacon positioning system transmitting time-stamped messages, e.g. GPS [Global Positioning System], GLONASS [Global Orbiting Navigation Satellite System] or GALILEO
    • G01S19/42—Determining position
    • G01S19/48—Determining position by combining or switching between position solutions derived from the satellite radio beacon positioning system and position solutions derived from a further system
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S19/00—Satellite radio beacon positioning systems; Determining position, velocity or attitude using signals transmitted by such systems
    • G01S19/38—Determining a navigation solution using signals transmitted by a satellite radio beacon positioning system
    • G01S19/39—Determining a navigation solution using signals transmitted by a satellite radio beacon positioning system the satellite radio beacon positioning system transmitting time-stamped messages, e.g. GPS [Global Positioning System], GLONASS [Global Orbiting Navigation Satellite System] or GALILEO
    • G01S19/42—Determining position
    • G01S19/48—Determining position by combining or switching between position solutions derived from the satellite radio beacon positioning system and position solutions derived from a further system
    • G01S19/485—Determining position by combining or switching between position solutions derived from the satellite radio beacon positioning system and position solutions derived from a further system whereby the further system is an optical system or imaging system
    • G—PHYSICS
    • G06—COMPUTING OR CALCULATING; COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10—Complex mathematical operations
    • G06F17/16—Matrix or vector computation, e.g. matrix-matrix or matrix-vector multiplication, matrix factorization

Definitions

  • the invention relates to the field of geodetic surveying of remote sensing images, in particular to an InSAR and GNSS weighting method for three-dimensional ground deformation estimation.
  • Interferometric Synthetic Aperture Radar SAR, InSAR
  • GNSS Global Navigation Satellite System
  • InSAR technology processes two SAR images of the same area at different times (intervals ranging from several days to hundreds of days) to obtain a one-dimensional resolution unit (several meters to tens of meters) along the radar line of sight within the time interval.
  • the average deformation result, the observation accuracy is generally at the millimeter or centimeter level.
  • GNSS technology uses a ground receiver to obtain a time-continuous three-dimensional coordinate sequence, and the difference between the two moments of time can obtain the three-dimensional surface deformation at the receiver.
  • InSAR and GNSS technologies have complementary advantages in surface deformation monitoring, and provide a new perspective for obtaining high-precision, high-spatial resolution three-dimensional surface deformation.
  • InSAR and GNSS Due to the differences in the deformation observation accuracy of InSAR and GNSS and the characteristics of observation targets, accurately determining the weight ratio between the two types of observations is essential for obtaining high-precision three-dimensional surface deformation results.
  • InSAR and GNSS are very susceptible to various uncertain factors when acquiring surface deformation, such as ionosphere, atmospheric water vapor, surface vegetation cover, etc., which makes it difficult to accurately estimate the prior variance information of various observations.
  • the prior variance of GNSS is mainly obtained based on GNSS network adjustment, while the prior variance of InSAR data assumes that there is no deformation in the far-field area, and the fitting result of the semivariant variance function is used as the prior variance of the entire InSAR image.
  • InSAR observation errors are often different in space, so its weighting accuracy is limited.
  • a priori variance estimate of the observation value can also be obtained, but this method is difficult to reflect the influence of the atmospheric long-wave error in the observation value.
  • the present invention aims to solve at least one of the technical problems existing in the prior art.
  • the present invention discloses an InSAR and GNSS weighting method for three-dimensional surface deformation estimation, which includes the following steps:
  • Step 1 Using the ascending and descending InSAR data of the area to be monitored, and the GNSS data of the area to be monitored, the three-dimensional deformation d 0 of the unknown point and a certain amount of InSAR/ of the surrounding points are established based on the surface stress and strain model and the imaging geometry of the observation value. functional relationship between L i GNSS data;
  • Step 2 The value of L i K i internal observation data and lift rail and the descending rail InSAR GNSS observations and other relatively fixed weight, determine the initial weight InSAR / GNSS observations of various types of weight matrix W i;
  • Step 3 Use variance component estimation to determine the precise weight matrix between various InSAR/GNSS observations Solve the high-precision three-dimensional surface deformation d 0 based on the least squares rule;
  • Step 4 Perform the above steps 1-3 for each surface point to realize the fusion of InSAR and GNSS to estimate the high-precision three-dimensional surface deformation field.
  • step 1 further comprises, between the points of the unknown three-dimensional distortion function d 0 and a certain number of points around InSAR / GNSS data L i is:
  • P 0 means unknown point
  • Is the coefficient matrix of the surface stress-strain model
  • I is a 3 ⁇ 3 identity matrix
  • l represents the unknown parameter vector at point P 0
  • the data of up-orbit InSAR and down-orbit InSAR are all a numerical value.
  • the representative GNSS data is a 3 ⁇ 1 vector.
  • the step 2 further includes: the surface stress-strain model is a description of the physical-mechanical relationship between the three-dimensional surface deformation of adjacent points on the surface; the observation imaging geometry is the relationship between the InSAR/GNSS observation value and the three-dimensional surface deformation The geometric relationship description.
  • W i diag(W i ′) indicates that the diagonal elements are in turn the diagonal matrix of the elements in the vector W i ′.
  • the inverse distance weighted attenuation factor D 0 is determined by the following formula:
  • K′ represents the number of all GNSS stations in the entire deformation field
  • K 3 ′ represents the number of GNSS stations closest to P 0
  • K 3 ′ takes a value of 4-6
  • 'K-th GNSS sites' representing all sites and the distance K between the nearest K P 0. 3 'site of a GNSS. 3 k' of sites.
  • step 3 further includes:
  • the high-precision three-dimensional surface deformation result is obtained, that is, the first, second, and third elements of the unknown parameter vector l.
  • the present invention proposes an InSAR and GNSS weighting method for three-dimensional surface deformation estimation.
  • the method is based on the surface stress when InSAR and GNSS are fused to estimate three-dimensional surface deformation.
  • the strain model establishes the functional relationship between the InSAR/GNSS observations and the three-dimensional surface deformation of unknown points.
  • the variance component estimation algorithm is used to accurately determine the weight ratio between the two types of observations, InSAR and GNSS, and finally based on the least squares rule to achieve three-dimensional High-precision estimation of surface deformation.
  • the content of the present invention is to use the surface stress and strain model to provide redundant observations in space, so that the variance component estimation can obtain accurate InSAR/GNSS weight ratios while lacking time series data, thereby effectively improving the fusion of InSAR and GNSS to estimate the three-dimensional surface
  • the precision and universality of deformation is to use the surface stress and strain model to provide redundant observations in space, so that the variance component estimation can obtain accurate InSAR/GNSS weight ratios while lacking time series data, thereby effectively improving the fusion of InSAR and GNSS to estimate the three-dimensional surface The precision and universality of deformation.
  • Fig. 1 is a flow chart of a method of InSAR and GNSS fusion estimation of three-dimensional surface deformation based on variance component estimation of the present invention
  • FIG. 2 is a comparison diagram of the three-dimensional surface deformation field obtained by the method of the present invention and the traditional method and the original simulated three-dimensional surface deformation field;
  • Fig. 3 is a graph of InSAR simulation deformation data of ascending and descending orbits in an embodiment of the present invention.
  • Step 1 Use the ascending and descending InSAR data of the area to be monitored, and the GNSS data of the area, based on the surface stress-strain model (SM) to establish the three-dimensional surface deformation of the unknown point and a certain amount of InSAR/GNSS around the point The functional relationship between the data;
  • SM surface stress-strain model
  • H represents the unknown parameter matrix of the stress-strain model, which can be expressed as:
  • ⁇ and ⁇ represent the strain parameters and rotation parameters in the surface stress-strain model.
  • the representative up-orbit InSAR and down-orbit InSAR data are all a numerical value
  • the representative GNSS data is a 3 ⁇ 1 vector, namely Considering the geometric relationship between InSAR and GNSS observations and three-dimensional surface deformation, it can be established : Functional relationship between the three-dimensional deformation of the surface k and the point P k D
  • I is a 3 ⁇ 3 identity matrix
  • Step 2 Perform relative weighting on K i observation data within various observation values, that is, determine the initial weight matrix W i of various observation values;
  • the present invention uses the following formula to determine the initial weight of the InSAR/GNSS observation value at P k :
  • K′ represents the number of all GNSS stations in the entire deformation field.
  • K '3 represents the distance P 0 of the number of sites GNSS latest, based on experience and generally 4-6.
  • the weight ratio coefficient of the GNSS vertical observation value in equation (7) is 0.5.
  • the specific implementation process can be based on the first GNSS three-dimensional deformation value. The variance information adjusts this ratio parameter.
  • W i diag (W 'i ) denotes the diagonal elements are the vector W' i diagonal matrix elements.
  • the observation value corresponding to the minimum weight plays a negligible role in the process of solving unknown parameters. Therefore, when the method of the present invention does not consider the initial weight in the solution process GNSS stations less than the threshold W thr .
  • W thr is generally 10 -6 based on experience.
  • the number K 3 of GNSS sites used to establish a functional relationship in the first step of participation can be determined.
  • the number of various observation values should be approximately equal, that is, K 1 ⁇ K 2 ⁇ 3K 3 should be satisfied in the present invention.
  • the InSAR data of K 1 /K 2 ascending/descending orbits closest to P 0 are selected in the present invention to participate in the calculation of the unknown parameter vector l.
  • Step 3 Use variance component estimation to determine the precise weight matrix between various InSAR/GNSS observations Error in unit weight Solve high-precision three-dimensional surface deformation d 0 based on the least square rule;
  • the weight matrix of observations at this time is the optimal weight matrix. Since the initial weight matrix W i only considers the relative weights between the observation data within the same type of observations, and does not consider the weight ratios between different types of observations, the unit weights of the various observations obtained by equation (11) are The error often does not satisfy equation (12).
  • the present invention combines the idea of variance component estimation, using the following formula weights of all kinds of observed values W i is updated:
  • the high-precision three-dimensional surface deformation result can be obtained, that is, the first, second, and third elements of the unknown parameter vector l.
  • the fusion of InSAR and GNSS can be realized to estimate the high-precision three-dimensional surface deformation field.
  • Figure 2(a)-(c) are the original simulated east-west, north-south and vertical deformation data in sequence
  • Figure 2(d)- (f) In turn are the east-west, north-south and vertical deformation data obtained by the traditional method.
  • Figure 2(g)-(i) are the east-west, north-south and vertical deformation data obtained by the method of the present invention (unit : Cm)
  • Figure 3(a) is the rising orbit InSAR data
  • Figure 3(b) is the falling orbit InSAR data.
  • the triangles in the figure represent the location distribution of GNSS stations (unit: cm).
  • Simulation data description 1Simulate three-dimensional deformation fields in east-west, north-south and vertical directions in a certain area (image size 400 ⁇ 450) (as shown in Figure 2(a)-(c)); 2Combine sentinel-1A/B satellite data Calculate the ascending and descending InSAR deformation results, and the incident angle and azimuth angle of the ascending orbit data are 39.3 respectively. ,-12.2. , The incident angle and azimuth angle of the descending orbit data are 33.9 respectively. ,-167.8.
  • the inverse distance weighting method is used to amplify the prior variance of the GNSS, and the InSAR far-field data is used to fit the semivariogram to solve the prior variance of the far-field InSAR observations, and Take it as the prior variance of the entire InSAR image. Then, in the process of solving, the prior variance of InSAR and GNSS observations is used to determine the weights, and the three-dimensional surface deformation is solved under the least squares criterion.
  • this application can be provided as methods, systems, or computer program products. Therefore, this application may adopt the form of a complete hardware embodiment, a complete software embodiment, or an embodiment combining software and hardware. Moreover, this application may adopt the form of a computer program product implemented on one or more computer-usable storage media (including but not limited to disk storage, CD-ROM, optical storage, etc.) containing computer-usable program codes.
  • computer-usable storage media including but not limited to disk storage, CD-ROM, optical storage, etc.

Landscapes

  • Engineering & Computer Science (AREA)
  • Remote Sensing (AREA)
  • Radar, Positioning & Navigation (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Computer Networks & Wireless Communication (AREA)
  • Electromagnetism (AREA)
  • Mathematical Physics (AREA)
  • Mathematical Optimization (AREA)
  • Mathematical Analysis (AREA)
  • Computational Mathematics (AREA)
  • Pure & Applied Mathematics (AREA)
  • Data Mining & Analysis (AREA)
  • Theoretical Computer Science (AREA)
  • Computing Systems (AREA)
  • Algebra (AREA)
  • Databases & Information Systems (AREA)
  • Software Systems (AREA)
  • General Engineering & Computer Science (AREA)
  • Geophysics And Detection Of Objects (AREA)
  • Length Measuring Devices With Unspecified Measuring Means (AREA)

Abstract

一种面向三维地表形变估计的InSAR和GNSS定权方法，包括：步骤1：利用待监测区域升轨和降轨InSAR数据，以及待监测区域的GNSS数据，基于地表应力应变模型及观测值成像几何建立未知点三维形变d 0与周围点一定数量的InSAR/GNSS数据L i之间的函数关系；步骤2：对升轨和降轨InSAR和GNSS观测值L i内部的K i个观测数据进行相对定权，确定InSAR/GNSS各类观测值的初始权重矩阵W i；步骤3：利用方差分量估计确定InSAR/GNSS各类观测值之间的精确权重矩阵式（I），基于最小二乘准则求解高精度的三维地表形变d 0；步骤4：对每一个地表点经过上述步骤1-3实现InSAR和GNSS融合估计高精度三维地表形变场。

Description

一种面向三维地表形变估计的InSAR和GNSS定权方法 技术领域
本发明涉及遥感影像的大地测量领域，尤其涉及一种面向三维地表形变估计的InSAR和GNSS定权方法。
背景技术
合成孔径雷达干涉测量(Interferometric Synthetic Aperture Radar,SAR，InSAR)和全球导航卫星系统(Global Navigation Satellite System，GNSS)已经被广泛用于获取地震、火山、地下开采等引起的地表形变。InSAR技术对同一区域不同时间(间隔为几天至几百天)的两景SAR影像进行处理即可得到地表某分辨单元(几米至几十米)在该时间间隔内沿雷达视线向的一维平均形变结果，其观测精度一般在毫米级或厘米级。GNSS技术则是通过地面接收机获取时间连续的三维坐标序列，对两个时刻的坐标作差，即可获取接收机处的三维地表形变，其水平方向精度可达亚毫米级，垂直向精度可达毫米级。由此可见，InSAR和GNSS技术在地表形变监测方面优势互补，为获取高精度、高空间分辨率的三维地表形变提供了新视角。
由于InSAR和GNSS的形变观测精度及观测目标特点的差异性，准确确定两类观测值之间的权重比例对于获取高精度三维地表形变结果至关重要。事实上，InSAR和GNSS获取地表形变时极易受各种不确定因素影响，例如电离层，大气水汽，地表植被覆盖等，导致难以准确估计各类观测值的先验方差信息。目前，GNSS的先验方差主要根据GNSS网平差获得，而InSAR数据的先验方差，则是假定远场区域没有形变，将半变异方差函数的拟合结果作为整个InSAR影像的先验方差，进而即可实现二者之间定权。但是，InSAR观测误差在空间上往往是有差异的，因此其定权精度有限。另外，通过InSAR观测精度与相干性的经验公式，也可获取观测值的先验方差估值，但该方法难以反映观测值中大气等长波误差的影响。
发明内容
本发明旨在至少解决现有技术中存在的技术问题之一。为此，本发明公开了一种面向三维地表形变估计的InSAR和GNSS定权方法，包括以下步骤：
步骤1：利用待监测区域升轨和降轨InSAR数据，以及所述待监测区域的GNSS数据，基于地表应力应变模型及观测值成像几何建立未知点三维形变d 0与周围点一定数量的InSAR/GNSS数据L i之间的函数关系；
步骤2：对升轨和降轨InSAR和GNSS等观测值L i内部的K i个观测数据进行相对定权，确定InSAR/GNSS各类观测值的初始权重矩阵W i；
步骤3：利用方差分量估计确定InSAR/GNSS各类观测值之间的精确权重矩阵
Figure PCTCN2020091273-appb-000001
基于最小二乘准则求解高精度的所述三维地表形变d 0；
步骤4：对每一个地表点经过上述步骤1-3实现InSAR和GNSS融合估计高精度三维地表形变场。
更进一步地，所述步骤1进一步包括，所述未知点三维形变d 0与周围点一定数量的InSAR/GNSS数据L i之间函数关系为：
Figure PCTCN2020091273-appb-000002
其中，
Figure PCTCN2020091273-appb-000003
P 0表示未知点，
Figure PCTCN2020091273-appb-000004
为地表应力应变模型系数矩阵，
Figure PCTCN2020091273-appb-000005
I为3×3的单位矩阵，l代表P 0点处的未知参数向量，
Figure PCTCN2020091273-appb-000006
为InSAR/GNSS数据，且i＝1,2,3，
Figure PCTCN2020091273-appb-000007
代表的升轨InSAR、降轨InSAR数据均是一个数值，
Figure PCTCN2020091273-appb-000008
代表的GNSS数据是一个3×1的向量。
更进一步地，所述步骤2进一步包括：所述地表应力应变模型为地表临近点三维地表形变之间的物理力学关系描述；所述观测值成像几何是InSAR/GNSS观测值与三维地表形变之间的几何关系描述。
确定P k处的InSAR/GNSS观测值的初始权重：
Figure PCTCN2020091273-appb-000009
其中，
Figure PCTCN2020091273-appb-000010
表示P k处的初始权重，
Figure PCTCN2020091273-appb-000011
代表P k与P 0之间的距离，D 0代表反距离定权衰减因子；
确定各类观测值的初始权重矩阵：
W i＝diag(W i′)
其中，
Figure PCTCN2020091273-appb-000012
W i＝diag(W i′)表示对角线元素依次是向量W i′中元素的对角矩阵。
更进一步地，所述反距离定权衰减因子D 0通过下式确定：
Figure PCTCN2020091273-appb-000013
其中，K′代表整个形变场中所有GNSS站点的个数，K 3′代表距离P 0最近的GNSS站点的个数，K 3′取值4-6，
Figure PCTCN2020091273-appb-000014
代表所有K′个GNSS站点中的第k′个站点与距离P 0最近的K 3′个GNSS站点中的第k 3′个站点之间的距离。
更进一步地，所述步骤3进一步包括，
利用方差分量估计确定InSAR/GNSS各类观测值之间的精确权重矩阵
Figure PCTCN2020091273-appb-000015
及其单位权中误差
Figure PCTCN2020091273-appb-000016
基于最小二乘准则求解高精度三维地表形变d 0；
令
Figure PCTCN2020091273-appb-000017
可得：
l＝M -1N    (10)
进而根据方差分量估计算法可得：
σ 2＝Ψ -1δ   (11)
其中，
Figure PCTCN2020091273-appb-000018
为各类观测值的单位权中误差估值；Ψ为转换矩阵，δ为观测值改正数二次型向量；
通过式(13)对各类观测值权重W i进行更新：
Figure PCTCN2020091273-appb-000019
利用式(13)更新观测值权重矩阵，重新计算式(10)(11)，迭代此过程直至各类观测值单位权中误差满足
Figure PCTCN2020091273-appb-000020
之间差别小于阈值Δσ。
再根据式(10)得到高精度三维地表形变结果，即未知参数向量l的第1、2、3个元素。
更进一步地，所述转换矩阵Ψ为：
Figure PCTCN2020091273-appb-000021
更进一步地，所述观测值改正数二次型向量δ为：
Figure PCTCN2020091273-appb-000022
其中，观测值改正数v i＝B i·l-L i。
更进一步地，所述迭代此过程直至各类观测值单位权中误差满足
Figure PCTCN2020091273-appb-000023
之间差别小于阈值Δσ进一步包括：所述阈值Δσ 2＝1mm 2。
本发明与现有技术相比，取得的有益效果为：本发明提出了一种面向三维地表形变估计的InSAR和GNSS定权方法，该方法在InSAR和GNSS融合估计三维地表形变时，基于地表应力应变模型建立InSAR/GNSS观测值与未知点三维地表形变之间的函数关系,同时利用方差分量估计算法准确确定InSAR和GNSS两类观测值之间的权重比例，最后基于最小二乘准则，实现三维地表形变的高精度估计。而传统方法中需要大量时序上的InSAR/GNSS数据为方差分量估计定权提供多余观测，因此对于瞬时形变(如火山、地震等)并不适用。本发明内容则是利用地表应力应变模型在空间上提供多余观测，使得方差分量估计可在缺乏时序数据的同时也可以获取精确的InSAR/GNSS权重比例，进而有效提高了InSAR和GNSS融合估计三维地表形变的精度及普适性。
附图说明
从以下结合附图的描述可以进一步理解本发明。图中的部件不一定按比例绘制，而是将重点放在示出实施例的原理上。在图中，在不同的视图中，相同的附图标记指定对应的部分。
图1是本发明一种基于方差分量估计的InSAR和GNSS融合估计三维地表形变方法的流程图；
图2是本发明方法与传统方法得到的三维地表形变场与原始模拟三维地表形变场的对比图；
图3是本发明一实施例中的升轨和降轨InSAR模拟形变数据图。
具体实施方式
为了使本技术相关领域的人员能够更好地理解本发明，下面将对本发明的实施方案进行清楚、详细的描述。同时，在此对本发明中主要的公式符号进行说明如下：
P：点
x：点的坐标
d：三维地表形变
l：未知参数向量
B：系数矩阵
L：InSAR/GNSS观测值
W：InSAR/GNSS观测值权重
σ：InSAR/GNSS观测值单位权中误差
K：InSAR/GNSS观测值个数
D：两点间距离
V：方差
上标0/k：点的索引编号
下标i/1/2/3：InSAR/GNSS观测值类型索引编号
上标enu：与观测值相关的东西向(east-west)、南北向(north-south)及垂直向(up-down)变量
下标enu：与未知参数相关的东西向(east-west)、南北向(north-south)及垂直向(up-down)变量
实施例一
如图1所示，本实施例具体实施方式如下：
步骤1：利用待监测区域的升轨和降轨InSAR数据，以及该区域的GNSS数据，基于地表应力应变模型(Strain Model，SM)建立未知点三维地表形变与该点周围一定数量的InSAR/GNSS数据之间的函数关系；
如何确定用于建立函数关系InSAR/GNSS数据的数量将在步骤2中介绍。
假设未知点P 0的三维坐标和三维形变分别为
Figure PCTCN2020091273-appb-000024
Figure PCTCN2020091273-appb-000025
周围一点P k的三维坐标和三维形变分别为
Figure PCTCN2020091273-appb-000026
Figure PCTCN2020091273-appb-000027
那么根据地表应力应变模型有下式：
d k＝H·Δ k+d 0  (1)
其中
Figure PCTCN2020091273-appb-000028
H代表应力应变模型未知参数矩阵，可表示为：
Figure PCTCN2020091273-appb-000029
ξ和ω代表地表应力应变模型中的应变参数和旋转参数。
进而，可将式(1)写成：
Figure PCTCN2020091273-appb-000030
其中，
Figure PCTCN2020091273-appb-000031
代表地表应力应变模型系数矩阵。
Figure PCTCN2020091273-appb-000032
代表P 0点处的未知参数向量。
进一步，假设在P k处有升轨InSAR、降轨InSAR和GNSS三种数据的一种或多种，分别记为
Figure PCTCN2020091273-appb-000033
其中
Figure PCTCN2020091273-appb-000034
代表的升轨InSAR、降轨InSAR数据均是一个数值，而
Figure PCTCN2020091273-appb-000035
代表的GNSS数据是一个3×1的向量，即
Figure PCTCN2020091273-appb-000036
考虑InSAR和GNSS观测值与三维地表形变之间的几何关系，可建立
Figure PCTCN2020091273-appb-000037
与点P k处三维地表形变d k之间的函数关系：
Figure PCTCN2020091273-appb-000038
其中，
Figure PCTCN2020091273-appb-000039
I为3×3的单位矩阵，
Figure PCTCN2020091273-appb-000040
Figure PCTCN2020091273-appb-000041
分别代表获取InSAR数据时卫星的方位角和入射角。
综合式(3)和(4)，可得：
Figure PCTCN2020091273-appb-000042
其中，
Figure PCTCN2020091273-appb-000043
至此，即可建立周围点P k处的InSAR/GNSS观测值与P 0点处的未知参数向量l之间的函数关系。
假设点P 0周围有K 1个升轨InSAR、K 2个降轨InSAR和K 3个GNSS站点可用于估计未知参数向量l，则最终可得：
L＝B·l     (6)
其中，
L＝[(L 1) T,(L 2) T,(L 3) T] T,
Figure PCTCN2020091273-appb-000044
B＝[(B 1) T,(B 2) T,(B 3) T] T,
Figure PCTCN2020091273-appb-000045
步骤2：对各类观测值内部的K i个观测数据进行相对定权，即确定各类观测值的初始权重矩阵W i；
由于GNSS站点分布较为稀疏，与P 0不同距离的GNSS站点数据应赋予不同的权重。本发明利用下式来确定P k处的InSAR/GNSS观测值的初始权重：
Figure PCTCN2020091273-appb-000046
其中，
Figure PCTCN2020091273-appb-000047
代表P k与P 0之间的距离，D 0代表反距离定权衰减因子，可通过下式确定：
Figure PCTCN2020091273-appb-000048
K′代表整个形变场中所有GNSS站点的个数。K′ 3代表距离P 0最近的GNSS站点的个数，根据经验一般取4-6。
Figure PCTCN2020091273-appb-000049
代表所有K′个GNSS站点中的第k′个站点与距离P 0最近的K′ 3个GNSS站点中的第k′ 3个站点之间的距离。
值得注意的是，GNSS垂直向的观测精度往往比水平向精度低，因此式(7)中GNSS垂直向观测值的权重比例系数为0.5，具体实施过程中可根据GNSS三个维度形变值的先验方差信息调整此比例参数。
此时，即可确定各类观测值的初始权重矩阵：
W i＝diag(W′ i)   (9)
其中，
Figure PCTCN2020091273-appb-000050
W i＝diag(W′ i)表示对角线元素依次是向量W′ i中元素的对角矩阵。
同时，当一组数据中最小权重与最大权重之比小于一定阈值时，最小权重对应的观测值在未知参数解算过程中所起到的作用可以忽略不计。因此，当本发明方法在解算过程中不考虑初始权重
Figure PCTCN2020091273-appb-000051
小于阈值W thr的GNSS站点。其中，W thr根据经验一般取10 -6。
此时，即可确定参与步骤一中用于建立函数关系GNSS站点的数量K 3。为了使得方差分量估计可以更加准确的确定权重，各类观测值个数应当大致相等，即在本发明中应满足K 1≈K 2≈3K 3。基于此，本发明中选取距P 0最近的K 1/K 2个升轨/降轨的InSAR数据参与未知参数向量l的解算。
步骤3：利用方差分量估计确定InSAR/GNSS各类观测值之间的精确权重矩阵
Figure PCTCN2020091273-appb-000052
及其单位权中误差
Figure PCTCN2020091273-appb-000053
基于最小二乘准则求解高精度三维地表形变d 0；
令
Figure PCTCN2020091273-appb-000054
可得：
l＝M -1N     (10)
进而根据方差分量估计算法可得：
σ 2＝Ψ -1δ    (11)
其中，
Figure PCTCN2020091273-appb-000055
代表各类观测值的单位权中误差估值。
Figure PCTCN2020091273-appb-000056
代表转换矩阵。
Figure PCTCN2020091273-appb-000057
代表观测值改正数二次型向量。
v i＝B i·l-L i代表观测值改正数。
根据方差分量估计算法可得，当各类观测值单位权中误差近似相等时，即
Figure PCTCN2020091273-appb-000058
此时的观测值权重矩阵为最优权阵。由于初始权重矩阵W i只考虑了同一类观测值内部各个观测数据之间的相对权重，并未考虑不同类观测值之间的权重比例，因此式(11)得到的各类观测值单位权中误差往往并不满足式(12)。本发明结合方差分量估计思路，利用下式对各类观测值权重W i进行更新：
Figure PCTCN2020091273-appb-000059
利用式(13)更新观测值权重矩阵，重新计算式(10)(11)，迭代此过程直至各类观测值单位权中误差满足式(12)，即
Figure PCTCN2020091273-appb-000060
之间差别小于阈值Δσ，本发明中Δσ 2＝1mm 2。
此时，根据式(10)即可得到高精度三维地表形变结果，即未知参数向量l的第1,2,3个元素。
对每一个地表点经过上述步骤1-3即可实现InSAR和GNSS融合估计高精度三维地表形变场。
实施例二
本实施通过实验对本发明进行验证，如图2-3所示，其中，图2(a)-(c)依次为原始模拟的东西向、南北向和垂直向形变数据，图2(d)-(f)依次为传统方法得到的东西向、南北向和垂直向的形变数据，图2(g)-(i)依次为本发明方法得到的东西向、南北向和垂直向的形变数据(单位：厘米)；图3(a)为升轨InSAR数据，图3(b)为降轨InSAR数据，图中三角形代表GNSS站点的位置分布(单位：厘米)。
模拟数据描述：①在一定区域(图像大小400×450)模拟东西向、南北向及垂直向的三维形变场(如图2(a)-(c))；②结合哨兵-1A/B卫星数据的成像几何，计算升轨和降轨InSAR形变结果，其中升轨数据的入射角和方位角分别为39.3。,-12.2。，降轨数据的入射角和方位角分别为33.9。,-167.8。；③在升轨和降轨InSAR数据中分别加人方差为4mm和6mm的高斯噪声，同时在两景InSAR数据中也加入了一定量级的大气延迟误差，得到的总误差均方根分别为4.9mm和6.9mm。此时即可得到用于模拟实验的InSAR原始数据(如图3所示)。④同时，在形变场中随机选取了100个像素，作为GNSS观测站点的位置，其相应位置上原始模拟的三维形变作为GNSS观测值，同时在GNSS水平形变观测值加入方差为1mm的高斯噪声，在GNSS垂直形变观测值加入方差为2mm的高斯噪声。其GNSS站点的分布如图3中的三角形所示。
传统融合InSAR和GNSS估计三维地表形变场时，利用反距离加权的方式对GNSS先验方差进行放大，利用InSAR远场数据进行半变异函数拟合求解远场区域InSAR观测值的先验方差，并将其作为整个InSAR影像的先验方差。然后在求解过程中，利用InSAR和GNSS观测值的先验方差进行定权，在最小二乘准则下求解三维地表形变。本 次模拟实验分别利用传统方法(图2(d)-(f))和本发明方法(图2(g)-(i))对模拟数据进行三维地表形变场的解算，两种方法求解的三维地表形变场的均方根误差如表1所示。
表1 三维地表形变场残差均方根误差
Figure PCTCN2020091273-appb-000061
综合表1,图3可知，本发明提及的算法相比于传统算法可得到更为精确的三维地表形变场。
以上应用了具体个例对本发明进行阐述，只是为了帮助本领域中的普通技术人员很好的理解。在不偏离本发明的精神和范围的情况下，还可以对本发明的具体实施方式作各种推演、变形和替换。这些变更和替换都将落在本发明权利要求书所限定的范围内。
还需要说明的是，术语“包括”、“包含”或者其任何其他变体意在涵盖非排他性的包含，从而使得包括一系列要素的过程、方法、商品或者设备不仅包括那些要素，而且还包括没有明确列出的其他要素，或者是还包括为这种过程、方法、商品或者设备所固有的要素。在没有更多限制的情况下，由语句“包括一个……”限定的要素，并不排除在包括所述要素的过程、方法、商品或者设备中还存在另外的相同要素。
本领域技术人员应明白，本申请的实施例可提供为方法、系统或计算机程序产品。因此，本申请可采用完全硬件实施例、完全软件实施例或结合软件和硬件方面的实施例的形式。而且，本申请可采用在一个或多个其中包含有计算机可用程序代码的计算机可用存储介质(包括但不限于磁盘存储器、CD-ROM、光学存储器等)上实施的计算机程序产品的形式。
虽然上面已经参考各种实施例描述了本发明，但是应当理解，在不脱离本发明的范围的情况下，可以进行许多改变和修改。因此，其旨在上述详细描述被认为是例示性的而非限制性的，并且应当理解，以下权利要求(包括所有等同物)旨在限定本发明的精神和范围。以上这些实施例应理解为仅用于说明本发明而不用于限制本发明的保护范围。在阅读了本发明的记载的内容之后，技术人员可以对本发明作各种改动或修改，这些等效变化和修饰同样落入本发明权利要求所限定的范围。

Claims (8)

  1. 一种面向三维地表形变估计的InSAR和GNSS定权方法，其特征在于，包括以下步骤：
    步骤1：利用待监测区域升轨和降轨InSAR数据，以及所述待监测区域的GNSS数据，基于地表应力应变模型及观测值成像几何建立未知点三维形变d 0与周围点一定数量的InSAR/GNSS数据L i之间的函数关系；
    步骤2：对升轨和降轨InSAR和GNSS等观测值L i内部的K i个观测数据进行相对定权，确定InSAR/GNSS各类观测值的初始权重矩阵W i；
    步骤3：利用方差分量估计确定InSAR/GNSS各类观测值之间的精确权重矩阵
    Figure PCTCN2020091273-appb-100001
    基于最小二乘准则求解高精度的所述三维地表形变d 0；
    步骤4：对每一个地表点经过上述步骤1-3实现InSAR和GNSS融合估计高精度三维地表形变场。
  2. 如权利要求1所述的一种的方法，其特征在于，所述步骤1进一步包括，所述未知点三维形变d 0与周围点一定数量的InSAR/GNSS数据L i之间函数关系为：
    Figure PCTCN2020091273-appb-100002
    其中，
    Figure PCTCN2020091273-appb-100003
    P 0表示未知点，
    Figure PCTCN2020091273-appb-100004
    为地表应力应变模型系数矩阵，
    Figure PCTCN2020091273-appb-100005
    I为3×3的单位矩阵，l代表P 0点处的未知参数向量，
    Figure PCTCN2020091273-appb-100006
    为InSAR/GNSS数据，且i＝1,2,3，
    Figure PCTCN2020091273-appb-100007
    代表的升轨InSAR、降轨InSAR数据均是一个数值，
    Figure PCTCN2020091273-appb-100008
    代表的GNSS数据是一个3×1的向量。
  3. 如权利要求2所述的一种的方法，其特征在于，所述步骤2进一步包括：所述地表应力应变模型为地表临近点三维地表形变之间的物理力学关系描述；所述观测值成像几何是InSAR/GNSS观测值与三维地表形变之间的几何关系描述，
    确定P k处的InSAR/GNSS观测值的初始权重：
    Figure PCTCN2020091273-appb-100009
    其中，
    Figure PCTCN2020091273-appb-100010
    表示P k处的初始权重，
    Figure PCTCN2020091273-appb-100011
    代表P k与P 0之间的距离，D 0代表反距离定权衰减因子；
    确定各类观测值的初始权重矩阵：
    W i＝diag(W i′)
    其中，
    Figure PCTCN2020091273-appb-100012
    W i＝diag(W i′)表示对角线元素依次是向量W i′中元素的对角矩阵。
  4. 如权利要求3所述的一种的方法，其特征在于，所述反距离定权衰减因子D 0通过下式确定：
    Figure PCTCN2020091273-appb-100013
    其中，K′代表整个形变场中所有GNSS站点的个数，K 3′代表距离P 0最近的GNSS站点的个数，K 3′取值4-6，
    Figure PCTCN2020091273-appb-100014
    代表所有K′个GNSS站点中的第k′个站点与距离P 0最近的K 3′个GNSS站点中的第k 3′个站点之间的距离。
  5. 如权利要求4所述的一种的方法，其特征在于，所述3进一步包括，
    利用方差分量估计确定InSAR/GNSS各类观测值之间的精确权重矩阵
    Figure PCTCN2020091273-appb-100015
    及其单位权中误差
    Figure PCTCN2020091273-appb-100016
    基于最小二乘准则求解高精度三维地表形变d 0；
    令
    Figure PCTCN2020091273-appb-100017
    可得：
    l＝M -1N  (10)
    进而根据方差分量估计算法可得：
    σ 2＝Ψ -1δ  (11)
    其中，
    Figure PCTCN2020091273-appb-100018
    为各类观测值的单位权中误差估值；Ψ为转换矩阵，δ为观测值改正数二次型向量；
    通过式(13)对各类观测值权重W i进行更新：
    Figure PCTCN2020091273-appb-100019
    利用式(13)更新观测值权重矩阵，重新计算式(10)(11)，迭代此过程直至各类观测值单位权中误差满足
    Figure PCTCN2020091273-appb-100020
    之间差别小于阈值Δσ；
    再根据式(10)得到高精度三维地表形变结果，即未知参数向量l的第1、2、3个元素。
  6. 如权利要求5所述的一种的方法，其特征在于，所述转换矩阵Ψ为：
    Figure PCTCN2020091273-appb-100021
  7. 如权利要求6所述的一种的方法，其特征在于，所述观测值改正数二次型向量δ为：
    Figure PCTCN2020091273-appb-100022
    其中，观测值改正数v i＝B i·l-L i。
  8. 如权利要求7所述的一种的方法，其特征在于，所述迭代此过程直至各类观测值单位权中误差满足
    Figure PCTCN2020091273-appb-100023
    之间差别小于阈值Δσ进一步包括：阈值Δσ 2＝1mm 2。
PCT/CN2020/091273 2019-05-21 2020-05-20 一种面向三维地表形变估计的InSAR和GNSS定权方法 Ceased WO2020233591A1 (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
US17/035,771 US20210011149A1 (en) 2019-05-21 2020-09-29 InSAR and GNSS weighting method for three-dimensional surface deformation estimation

Applications Claiming Priority (2)

Application Number Priority Date Filing Date Title
CN201910423735.8A CN110058236B (zh) 2019-05-21 2019-05-21 一种面向三维地表形变估计的InSAR和GNSS定权方法
CN201910423735.8 2019-05-21

Related Child Applications (1)

Application Number Title Priority Date Filing Date
US17/035,771 Continuation US20210011149A1 (en) 2019-05-21 2020-09-29 InSAR and GNSS weighting method for three-dimensional surface deformation estimation

Publications (1)

Publication Number Publication Date
WO2020233591A1 true WO2020233591A1 (zh) 2020-11-26

Family

ID=67323791

Family Applications (1)

Application Number Title Priority Date Filing Date
PCT/CN2020/091273 Ceased WO2020233591A1 (zh) 2019-05-21 2020-05-20 一种面向三维地表形变估计的InSAR和GNSS定权方法

Country Status (3)

Country Link
US (1) US20210011149A1 (zh)
CN (1) CN110058236B (zh)
WO (1) WO2020233591A1 (zh)

Cited By (18)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN112835043A (zh) * 2021-01-06 2021-05-25 中南大学 一种任意方向的二维形变的监测方法
CN112986990A (zh) * 2021-02-04 2021-06-18 中国地质大学（北京） 一种大气相位改正方法及系统
CN112986993A (zh) * 2021-02-07 2021-06-18 同济大学 一种基于空间约束的InSAR形变监测方法
CN114089335A (zh) * 2021-11-16 2022-02-25 安徽理工大学 一种基于单轨道InSAR的山区开采沉陷三维变形提取方法
CN114114332A (zh) * 2021-11-03 2022-03-01 湖北理工学院 一种有效探测gnss基准站坐标时间序列不连续点的方法
CN114114256A (zh) * 2021-11-08 2022-03-01 辽宁工程技术大学 一种基于D-InSAR-GIS叠加分析技术的大面积矿区沉陷监测方法
CN114578356A (zh) * 2022-03-02 2022-06-03 中南大学 基于深度学习的分布式散射体形变监测方法、系统及设备
CN114757238A (zh) * 2022-06-15 2022-07-15 武汉地铁集团有限公司 地铁保护区变形监测的方法、系统、电子设备和存储介质
CN114966689A (zh) * 2022-05-27 2022-08-30 厦门理工学院 沿海城市时序InSAR沉降监测分析方法、装置、设备及介质
CN116049929A (zh) * 2022-10-26 2023-05-02 马培峰 一种城市建筑物风险等级InSAR评估和预测方法
CN116148852A (zh) * 2022-12-27 2023-05-23 北京理工大学 基于空时连续的北斗InSAR三维高精度形变反演方法
CN117055082A (zh) * 2023-09-01 2023-11-14 兰州交通大学 一种基于gnss时间序列的精准垂直形变提取方法
CN117109426A (zh) * 2023-08-28 2023-11-24 兰州交通大学 一种融合GNSS/InSAR观测资料的三维形变场建模方法
CN117213443A (zh) * 2023-11-07 2023-12-12 江苏省地质调查研究院 一种天地深一体化地面沉降监测网建设与更新方法
CN117607865A (zh) * 2023-10-30 2024-02-27 云南大学 一种地铁沿线形变检测方法、装置、系统以及存储介质
CN118226529A (zh) * 2024-03-28 2024-06-21 应急管理部国家自然灾害防治研究院 一种基于深度学习的InSAR同震三维形变场恢复方法
CN119395703A (zh) * 2024-11-27 2025-02-07 众智软件股份有限公司 核电站InSAR地表形变实时监测预警方法与系统
CN120405652A (zh) * 2025-03-28 2025-08-01 中国矿业大学(北京) 边坡雷达与gnss三维形变数据的联合观测融合方法及系统

Families Citing this family (50)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN110058236B (zh) * 2019-05-21 2023-04-07 中南大学 一种面向三维地表形变估计的InSAR和GNSS定权方法
CN111077525B (zh) * 2019-12-20 2022-12-27 长安大学 融合sar与光学偏移量技术的地表三维形变计算方法及系统
CN110804912B (zh) * 2020-01-06 2020-05-19 北京铁科工程检测有限公司 一种铁路线路及沿线区域形变信息的提取方法
CN111339483B (zh) * 2020-01-19 2022-03-11 武汉大学 一种基于去趋势互相关分析的gnss影像生成方法
US20230245287A1 (en) * 2020-07-20 2023-08-03 Nec Corporation Image processing device and image processing method
CN112540369A (zh) * 2020-11-27 2021-03-23 武汉大学 融合GNSS与升降轨时序InSAR的滑坡三维形变解算方法及系统
CN112711022B (zh) * 2020-12-18 2022-08-30 中国矿业大学 一种GNSS层析技术辅助的InSAR大气延迟改正方法
CN112797886B (zh) * 2021-01-27 2022-04-22 中南大学 面向缠绕相位的InSAR时序三维形变监测方法
CN113091596B (zh) * 2021-03-31 2022-01-25 中国矿业大学 一种基于多极化时序sar数据的地表形变监测方法
CN113091598B (zh) * 2021-04-06 2022-02-08 中国矿业大学 一种InSAR划定采空区建筑场地稳定性等级范围的方法
CN113096005B (zh) * 2021-04-06 2023-07-07 中国科学院生态环境研究中心 一种监测山体现今抬升速度的雷达时序差分干涉测量方法
CN113219414B (zh) * 2021-04-22 2024-04-02 桂林理工大学 一种消除卫星干涉雷达地表形变方向模糊新方法
CN113281742B (zh) * 2021-06-02 2023-07-25 西南交通大学 一种基于滑坡形变信息和气象数据的sar滑坡预警方法
CN113777606B (zh) * 2021-08-12 2023-12-26 北京理工大学 分布式geo sar三维形变反演多角度选取方法及装置
CN113723531B (zh) * 2021-09-02 2024-06-14 淮阴师范学院 面向全运行周期的矿区地表形变实时/准实时监测方法
CN113899301B (zh) * 2021-09-15 2022-07-15 武汉大学 联合gnss三维形变的区域陆地水储量变化反演方法及系统
CN114444377B (zh) * 2021-12-24 2024-06-21 北京理工大学 一种基于梯度提升机的多地面测距仪选站方法
CN115032637B (zh) * 2022-06-07 2024-06-14 安徽理工大学 一种井工开采全生命周期地表沉陷监测方法
CN115267774B (zh) * 2022-06-14 2024-12-24 深圳大学 一种城区多时相InSAR相位解缠方法、终端及存储介质
CN115201823B (zh) * 2022-07-22 2023-08-04 电子科技大学 一种利用BDS-InSAR数据融合的地表形变监测方法
CN115358311B (zh) * 2022-08-16 2023-06-16 南昌大学 地表变形监测多源数据融合处理方法
CN115455659A (zh) * 2022-08-19 2022-12-09 武汉大学 组合卫星重力场的新方法及系统
CN115574706B (zh) * 2022-10-11 2025-09-26 中国水利水电科学研究院 一种基于gnss的土石坝表面变形高精度监测方法
CN115855003B (zh) * 2022-12-01 2025-04-15 中南大学 一种InSAR辅助的GNSS监测网的优化布设方法
CN115629384B (zh) * 2022-12-08 2023-04-11 中南大学 一种时序InSAR误差的改正方法及相关设备
CN116050657B (zh) * 2023-02-28 2025-09-09 中南大学 关闭矿井的地表抬升预测方法、装置、终端设备及介质
CN116124052B (zh) * 2023-03-02 2025-09-23 中铁第四勘察设计院集团有限公司 桥梁综合变形监测系统及方法
CN116338690B (zh) * 2023-03-28 2024-01-16 中南林业科技大学 一种无模型约束的时序InSAR地形残差与地表形变估计方法
CN116485857B (zh) * 2023-05-05 2024-06-28 中山大学 一种基于多源遥感数据的高时间分辨率冰川厚度反演方法
CN117168373B (zh) * 2023-07-20 2024-07-09 中国卫通集团股份有限公司 基于卫星通导遥一体化的水库坝体形变监测系统
CN116659429A (zh) * 2023-08-01 2023-08-29 齐鲁空天信息研究院 一种多源数据高精度时序地表三维形变解算方法和系统
CN117113277A (zh) * 2023-09-05 2023-11-24 北京睿知行科技有限公司 基于传感器数据和卫星数据交互迭代的融合方法和系统
CN118053070B (zh) * 2024-01-08 2025-01-03 北京师范大学 基于双轨道雷达卫星的植被分布提取方法、装置及设备
CN118037563B (zh) * 2024-01-19 2024-11-22 中国民用航空飞行学院 一种机场跑道沉降预测方法及系统
CN117761729A (zh) * 2024-01-26 2024-03-26 中国地震局地震预测研究所 监测断层闭锁深度的gnss站点布设方法
CN118089611B (zh) * 2024-04-17 2024-08-30 东南大学 一种融合InSAR数据和物理知识的建筑三向位移监测方法及系统
CN118259280B (zh) * 2024-05-28 2024-08-06 深圳大学 联合InSAR与GNSS的填海区机场形变测评方法、系统及终端
CN118794484A (zh) * 2024-06-13 2024-10-18 国能榆林能源有限责任公司 一种基于干涉雷达的地质灾害监测预报预警系统及方法
CN118551665B (zh) * 2024-07-26 2024-10-22 兰州交通大学 一种基于InSAR和双向门控循环单元的地表形变预测方法
CN118656707B (zh) * 2024-08-16 2024-12-31 中南大学 基于InSAR和多源数据的地铁线路形变风险监测方法
CN118656684B (zh) * 2024-08-20 2024-12-10 中国科学院精密测量科学与技术创新研究院 一种基于InSAR数据和多诱发因子的形变智能预测方法
CN119064929A (zh) * 2024-08-21 2024-12-03 国能榆林能源有限责任公司 基于多源数据的煤矿采空区形变动态耦合监测方法及系统
CN119199853B (zh) * 2024-09-14 2025-05-02 中国地质大学（武汉） 一种考虑环境因素的sar时序位移同步观测值求解方法
CN119575376B (zh) * 2024-11-12 2025-06-13 鄂尔多斯市腾远煤炭有限责任公司 一种基于露天煤矿边坡形变聚集区自动提取方法
CN119936876B (zh) * 2024-12-05 2025-09-26 武汉大学 北斗辅助InSAR地表形变序列时间分辨率的优化方法
CN119902242B (zh) * 2024-12-25 2025-09-23 武汉大学 一种探测短波长形变的gnss成像方法及设备
CN119986647A (zh) * 2025-01-15 2025-05-13 湖南科技大学 InSAR与GNSS数据融合反演地表高精三维形变的建模方法
CN119986652B (zh) * 2025-01-24 2025-10-10 昆明理工大学 一种开采沉陷地表三维形变估计方法、装置、设备、介质及产品
CN120101711B (zh) * 2025-04-28 2025-08-05 中安国泰(北京)科技发展有限公司 自适应相位解缠地基sar形变监测方法及系统
CN120742316B (zh) * 2025-08-21 2025-11-11 中国矿业大学 基于总功率极化与非局域相位连接的干涉相位优化方法

Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20090237297A1 (en) * 2008-02-06 2009-09-24 Halliburton Energy Services, Inc. Geodesy Via GPS and INSAR Integration
CN104699966A (zh) * 2015-03-09 2015-06-10 中南大学 一种融合GNSS和InSAR数据获取高时空分辨率形变序列的方法
CN107102332A (zh) * 2017-05-11 2017-08-29 中南大学 基于方差分量估计与应力应变模型的InSAR三维地表形变监测方法
CN110058236A (zh) * 2019-05-21 2019-07-26 中南大学 一种面向三维地表形变估计的InSAR和GNSS定权方法

Family Cites Families (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101738620B (zh) * 2008-11-19 2012-02-15 中国农业科学院农业资源与农业区划研究所 从被动微波遥感数据amsr-e反演地表温度的方法
US8384583B2 (en) * 2010-06-07 2013-02-26 Ellegi S.R.L. Synthetic-aperture radar system and operating method for monitoring ground and structure displacements suitable for emergency conditions
US9207318B2 (en) * 2011-06-20 2015-12-08 California Institute Of Technology Damage proxy map from interferometric synthetic aperture radar coherence
CN103698750A (zh) * 2014-01-07 2014-04-02 国家卫星海洋应用中心 海洋二号卫星散射计海面风场反演方法和装置
CN104122553B (zh) * 2014-07-23 2017-01-25 中国国土资源航空物探遥感中心 一种集成多轨道、长条带CTInSAR的区域性地面沉降监测方法
CN106226767B (zh) * 2016-07-12 2018-08-21 中南大学 基于单个雷达成像几何学sar影像的矿区三维时序形变监测方法

Patent Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20090237297A1 (en) * 2008-02-06 2009-09-24 Halliburton Energy Services, Inc. Geodesy Via GPS and INSAR Integration
CN104699966A (zh) * 2015-03-09 2015-06-10 中南大学 一种融合GNSS和InSAR数据获取高时空分辨率形变序列的方法
CN107102332A (zh) * 2017-05-11 2017-08-29 中南大学 基于方差分量估计与应力应变模型的InSAR三维地表形变监测方法
CN110058236A (zh) * 2019-05-21 2019-07-26 中南大学 一种面向三维地表形变估计的InSAR和GNSS定权方法

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
HU, JUN ET AL.: "Three-Dimensional Surface Displacements From InSAR and GPS Measurements With Variance Component Estimation", 《IEEE GEOSCIENCE AND REMOTE SENSING LETTERS》, vol. 9, no. 4,, 31 July 2012 (2012-07-31), XP011448385, ISSN: 1545-598X, DOI: 20200728095337A *
LIU, JIHONG ET AL.: "A Method for Measuring 3-D Surface Deformations With InSAR Based on Strain Model and Variance Component Estimation", 《IEEE TRANSACTIONS ON GEOSCIENCE AND REMOTE SENSING》, vol. 56, no. 1,, 31 January 2018 (2018-01-31), XP55755606, ISSN: 0196-2892, DOI: 20200728095059A *

Cited By (27)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN112835043A (zh) * 2021-01-06 2021-05-25 中南大学 一种任意方向的二维形变的监测方法
CN112986990B (zh) * 2021-02-04 2023-02-17 中国地质大学（北京） 一种大气相位改正方法及系统
CN112986990A (zh) * 2021-02-04 2021-06-18 中国地质大学（北京） 一种大气相位改正方法及系统
CN112986993A (zh) * 2021-02-07 2021-06-18 同济大学 一种基于空间约束的InSAR形变监测方法
CN114114332A (zh) * 2021-11-03 2022-03-01 湖北理工学院 一种有效探测gnss基准站坐标时间序列不连续点的方法
CN114114256A (zh) * 2021-11-08 2022-03-01 辽宁工程技术大学 一种基于D-InSAR-GIS叠加分析技术的大面积矿区沉陷监测方法
CN114089335B (zh) * 2021-11-16 2022-09-06 安徽理工大学 一种基于单轨道InSAR的山区开采沉陷三维变形提取方法
CN114089335A (zh) * 2021-11-16 2022-02-25 安徽理工大学 一种基于单轨道InSAR的山区开采沉陷三维变形提取方法
CN114578356B (zh) * 2022-03-02 2024-11-08 中南大学 基于深度学习的分布式散射体形变监测方法、系统及设备
CN114578356A (zh) * 2022-03-02 2022-06-03 中南大学 基于深度学习的分布式散射体形变监测方法、系统及设备
CN114966689A (zh) * 2022-05-27 2022-08-30 厦门理工学院 沿海城市时序InSAR沉降监测分析方法、装置、设备及介质
CN114757238A (zh) * 2022-06-15 2022-07-15 武汉地铁集团有限公司 地铁保护区变形监测的方法、系统、电子设备和存储介质
CN114757238B (zh) * 2022-06-15 2022-09-20 武汉地铁集团有限公司 地铁保护区变形监测的方法、系统、电子设备和存储介质
CN116049929A (zh) * 2022-10-26 2023-05-02 马培峰 一种城市建筑物风险等级InSAR评估和预测方法
CN116049929B (zh) * 2022-10-26 2023-09-29 马培峰 一种城市建筑物风险等级InSAR评估和预测方法
US12055624B2 (en) 2022-10-26 2024-08-06 Peifeng MA Building risk monitoring and predicting based on method integrating MT-InSAR and pore water pressure model
CN116148852A (zh) * 2022-12-27 2023-05-23 北京理工大学 基于空时连续的北斗InSAR三维高精度形变反演方法
CN117109426B (zh) * 2023-08-28 2024-03-22 兰州交通大学 一种融合GNSS/InSAR观测资料的三维形变场建模方法
CN117109426A (zh) * 2023-08-28 2023-11-24 兰州交通大学 一种融合GNSS/InSAR观测资料的三维形变场建模方法
CN117055082A (zh) * 2023-09-01 2023-11-14 兰州交通大学 一种基于gnss时间序列的精准垂直形变提取方法
CN117607865A (zh) * 2023-10-30 2024-02-27 云南大学 一种地铁沿线形变检测方法、装置、系统以及存储介质
CN117213443B (zh) * 2023-11-07 2024-03-19 江苏省地质调查研究院 一种天地深一体化地面沉降监测网建设与更新方法
CN117213443A (zh) * 2023-11-07 2023-12-12 江苏省地质调查研究院 一种天地深一体化地面沉降监测网建设与更新方法
CN118226529A (zh) * 2024-03-28 2024-06-21 应急管理部国家自然灾害防治研究院 一种基于深度学习的InSAR同震三维形变场恢复方法
CN119395703A (zh) * 2024-11-27 2025-02-07 众智软件股份有限公司 核电站InSAR地表形变实时监测预警方法与系统
CN120405652A (zh) * 2025-03-28 2025-08-01 中国矿业大学(北京) 边坡雷达与gnss三维形变数据的联合观测融合方法及系统
CN120405652B (zh) * 2025-03-28 2025-09-30 中国矿业大学(北京) 边坡雷达与gnss三维形变数据的联合观测融合方法及系统

Also Published As

Publication number Publication date
CN110058236A (zh) 2019-07-26
US20210011149A1 (en) 2021-01-14
CN110058236B (zh) 2023-04-07

Similar Documents

Publication Publication Date Title
WO2020233591A1 (zh) 一种面向三维地表形变估计的InSAR和GNSS定权方法
James et al. 3‐D uncertainty‐based topographic change detection with structure‐from‐motion photogrammetry: precision maps for ground control and directly georeferenced surveys
CN112393714B (zh) 一种基于无人机航拍与卫星遥感融合的影像校正方法
Snay et al. Continuously operating reference station (CORS): history, applications, and future enhancements
CN105698764B (zh) 一种光学遥感卫星影像时变系统误差建模补偿方法及系统
Cheng et al. Making an onboard reference map from MRO/CTX imagery for Mars 2020 lander vision system
Zhang et al. Satellite SAR geocoding with refined RPC model
CN113238228B (zh) 基于水准约束的三维地表形变获取方法、系统及装置
Wu et al. Geometric integration of high-resolution satellite imagery and airborne LiDAR data for improved geopositioning accuracy in metropolitan areas
Elshambaky et al. A novel three-direction datum transformation of geodetic coordinates for Egypt using artificial neural network approach
Feng et al. A hierarchical network densification approach for reconstruction of historical ice velocity fields in East Antarctica
Tiwari et al. Multi-sensor geodetic approach for landslide detection and monitoring
Tian et al. Automatic calibration method for airborne LiDAR systems based on approximate corresponding points model
CN117572378B (zh) 基于InSAR与北斗数据的山体沉降分析方法及设备
CN118311615B (zh) 一种gnss对流层与多路径联合建模纠正方法及系统
CN117388851A (zh) 一种InSAR变形监测精度可靠估计方法
CN120559651B (zh) 一种地表三维形变反演方法、装置、设备及介质
Ashkenazi Models for controlling national and continental networks
Kwon et al. Geodetic datum transformation to the global geocentric datum for seas and islands around Korea
CN116560207A (zh) 一种用于测量世界时的方法
Santos Filho et al. Cartographic Accuracy Standard (CAS) of the digital terrain model of the digital and continuous cartographic base of the state of Amapá: case study in the city of Macapá
Grant et al. The development and implementation of New Zealand Geodetic Datum 2000
CN114924270A (zh) 基于GNSS的InSAR形变监测基准建立方法及装置
CN115407367A (zh) 一种混合星座卫星导航定位精度衰减因子估计方法
Noinak et al. Testing Horizontal Coordinate Correction Model Used for Transformation from PPP GNSS Technique to Thai GNSS CORS Network Based on ITRF2014

Legal Events

Date Code Title Description
121 Ep: the epo has been informed by wipo that ep was designated in this application

Ref document number: 20810830

Country of ref document: EP

Kind code of ref document: A1

NENP Non-entry into the national phase

Ref country code: DE

122 Ep: pct application non-entry in european phase

Ref document number: 20810830

Country of ref document: EP

Kind code of ref document: A1

122 Ep: pct application non-entry in european phase

Ref document number: 20810830

Country of ref document: EP

Kind code of ref document: A1

32PN Ep: public notification in the ep bulletin as address of the adressee cannot be established

Free format text: NOTING OF LOSS OF RIGHTS PURSUANT TO RULE 112(1) EPC (EPO FORM 1205A DATED 24/05/2022)

122 Ep: pct application non-entry in european phase

Ref document number: 20810830

Country of ref document: EP

Kind code of ref document: A1