CN106226767B - 基于单个雷达成像几何学sar影像的矿区三维时序形变监测方法 - Google Patents
基于单个雷达成像几何学sar影像的矿区三维时序形变监测方法 Download PDFInfo
- Publication number
- CN106226767B CN106226767B CN201610546270.1A CN201610546270A CN106226767B CN 106226767 B CN106226767 B CN 106226767B CN 201610546270 A CN201610546270 A CN 201610546270A CN 106226767 B CN106226767 B CN 106226767B
- Authority
- CN
- China
- Prior art keywords
- time
- mining area
- sinking
- sar
- insar
- 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.)
- Active
Links
- 238000005065 mining Methods 0.000 title claims abstract description 64
- 238000000034 method Methods 0.000 title claims abstract description 58
- 238000012544 monitoring process Methods 0.000 title claims abstract description 22
- 238000003384 imaging method Methods 0.000 claims description 19
- 239000011159 matrix material Substances 0.000 claims description 11
- 230000000007 visual effect Effects 0.000 claims description 3
- 230000001360 synchronised effect Effects 0.000 claims description 2
- 241001269238 Data Species 0.000 abstract 1
- 238000007796 conventional method Methods 0.000 abstract 1
- 238000005516 engineering process Methods 0.000 description 7
- 230000008569 process Effects 0.000 description 7
- 238000010586 diagram Methods 0.000 description 4
- 238000004088 simulation Methods 0.000 description 4
- 230000000694 effects Effects 0.000 description 2
- 238000005259 measurement Methods 0.000 description 2
- 238000012986 modification Methods 0.000 description 2
- 230000004048 modification Effects 0.000 description 2
- 238000012545 processing Methods 0.000 description 2
- 208000004350 Strabismus Diseases 0.000 description 1
- 230000001186 cumulative effect Effects 0.000 description 1
- 238000011156 evaluation Methods 0.000 description 1
- 230000006872 improvement Effects 0.000 description 1
- 238000005305 interferometry Methods 0.000 description 1
- 238000004519 manufacturing process Methods 0.000 description 1
- 238000013508 migration Methods 0.000 description 1
- 230000005012 migration Effects 0.000 description 1
- 230000000644 propagated effect Effects 0.000 description 1
- 230000001902 propagating effect Effects 0.000 description 1
- 238000012502 risk assessment Methods 0.000 description 1
- 230000001629 suppression Effects 0.000 description 1
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
- 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/904—SAR modes
Landscapes
- Engineering & Computer Science (AREA)
- Remote Sensing (AREA)
- Radar, Positioning & Navigation (AREA)
- Physics & Mathematics (AREA)
- Electromagnetism (AREA)
- Computer Networks & Wireless Communication (AREA)
- General Physics & Mathematics (AREA)
- Radar Systems Or Details Thereof (AREA)
Abstract
本发明公开了一种基于单个雷达成像几何学SAR影像的矿区三维时序形变监测方法,利用单个雷达成像几何学SAR数据生成可用InSAR干涉对;并生成矿区地表多时域差分下沉观测值;建立时间相邻SAR影像期间的下沉速度与多时域下沉观测值之间的观测方程,并使用加权最小二乘法求解下沉速度;使用解出的下沉速度计算矿区地表在SAR影像获取时间的时序下沉;基于矿区地表东西、南北方向水平移动与下沉在对应方向的梯度成比例的关系,使用计算的矿区时序下沉估计出东西、南北方向的时序形变。本发明实现了仅利用单个雷达成像几何学SAR数据监测矿区地表三维时序形变,大大降低了传统方法中需要三个或以上的不同雷达成像几何学SAR数据的苛刻要求。
Description
技术领域
本发明涉及一种基于单个雷达成像几何学SAR影像的矿区三维时序形变监测方法。
背景技术
矿区地表三维时序形变监测对于提前评估和控制地下开采导致的潜在地质灾害和建构筑物损坏有着重要作用。利用传统大地测量技术(如GPS、水准测量等)监测矿区地表三维时序形变不仅成本高、效率低,且监测范围较小。合成孔径雷达干涉测量(Interferometric Synthetic Aperture Radar,简称InSAR)是一种能够获得地表厘米甚至毫米级形变的遥感技术。其基本原理就是通过对两幅或以上的合成孔径雷达(SyntheticAperture Radar,简称SAR)卫星影像进行差分干涉处理,从相位差中提取厘米甚至毫米级的雷达视线方向(line-of-sight,LOS)形变值。相对于传统的大地测量技术而言,InSAR具有全天候、大范围覆盖、高精度、低成本等优势。
时序InSAR技术为高级InSAR技术。该技术基于同一雷达成像几何学的多时域SAR影像实现了地表时序形变监测。然而,由于当前的雷达传感器均为斜视,所以时序InSAR技术获取的地表时序形变是沿着雷达LOS方向的,,其为地表真实三维形变按照雷达成像几何学的投影。因此,若想从一维LOS向时序形变中分解出地表三维时序形变,则至少需要三个不同雷达成像几何学的多时域SAR影像。然而,由于当前可用SAR卫星较少,且矿区大形变梯度容易导致InSAR干涉对失相关而不可用。所以,找到覆盖同一个矿区同一时间段内的三个以上的不同雷达成像几何学的多时域SAR数据几乎不可能。这种局限也大大制约了InSAR技术在矿区三维时序形变获取中的应用。
鉴于矿区地表水平移动与下沉梯度之间为比例关系,之前提出了“一种利用单个InSAR干涉对获取矿区地表三维形变场的方法”(专利号:CN201210440875)(简称“SIP”)。其仅利用单个InSAR干涉对即可获取矿区地表三维形变场,大大降低了传统方法对SAR数据的苛刻要求,拓展了InSAR技术在矿区的应用前景。然而,由于该方法仅能获取矿区地表在组成单个InSAR干涉对的两景SAR影像期间的地表三维形变,而不能获取三维时序形变。因此,该方法的应用前景受到了很大的制约,比如无法实现矿区地表动态地质灾害风险评估。其次,该方法中估计的下沉值的误差在传播到水平移动的过程中容易被放大,从而降低其获取的三维形变的精度。
发明内容
本发明提供了一种基于单个雷达成像几何学SAR影像的矿区三维时序形变监测方法,其目的在于克服传统InSAR时序三维形变获取方法至少需要三个不同雷达成像几何学SAR数据的苛刻要求,降低了监测成本、提高了监测精度。
一种基于单个雷达成像几何学SAR影像的矿区三维时序形变监测方法,包括以下几个步骤:
步骤1:利用单个雷达成像几何学SAR影像生成覆盖待监测矿区的可用InSAR干涉对;
步骤2:利用现有的SIP方法生成待监测矿区多时域差分下沉观测值ΔW:
ΔW=[ΔW1,ΔW2,…,ΔWm]T;
其中,m表示可用InSAR干涉对的个数;
所述多时域差分下沉观测值指一系列与可用InSAR干涉对时间同步的下沉值;
步骤3:建立矿区地表各点在两相邻SAR影像时刻间的平均下沉速度V=[V1,V2,…,Vn]T与多时域差分下沉观测值ΔW=[ΔW1,ΔW2,…,ΔWm]T之间的函数关系:Bm×n·Vn×1=ΔWm×1;
其中,SAR影像个数为n+1;B为m×n的可用InSAR干涉对辅影像和主影像获取时间差系数矩阵;
对于B的任意第k行,第IMk个元素前的所有元素均为0,从第IMk到第ISk-1个元素,依次为:第ISk-1个元素以后的所有元素均为0;
其中,IMk和ISk分别为生成第k个可用InSAR干涉对的主、辅影像获取时间索引,根据InSAR干涉对的组成情况获得;
本质上,即为在和之间沉降的累计相加;
对于任意的第k个InSAR干涉对,其差分下沉观测值为ΔWk,主影像的时间为辅影像时间为
用速度表示该观测值为(即为在和之间沉降的累计相加)。在系数矩阵Bm×n中,只有从第IMk个元素开始才有值,一直到ISk-1结束,在这一行里,其他元素全部为0:
即其余元素均为0;
在系数矩阵Bm×n中的第k行的具体表达形式如下:
步骤4:求解地表各点在时间相邻SAR影像时刻的平均下沉速度V:
V=[BT·B]-1·[BT·ΔW];
步骤5:估计矿区地表时序下沉W,并利用时序下沉计算获得东西、南北方向的时序水平移动,完成待监测矿区三维时序形变监测;
其中,Wk表示在tk时刻矿区地表的时序下沉值,tl表示时间,Vl表示tl-1至tl期间的平均下沉速度;
其中,E(i,j)和N(i,j)分别表示像素坐标为(i,j)的点在tk时刻的东西、南北方向的水平移动,k=1,2,…,n;b为待监测矿区水平移动系数,H为开采采深,tanβ为主要影响角正切,RE和RN分别为LOS向形变图在东西和南北方向的空间分辨率。
所述求解地表各点在时间相邻SAR影像时刻的平均速度V时采用加权最小二乘法求解:
V=[BT·Q·B]-1·[BT·Q·ΔW]
其中,Q地表各点多时域差分下沉值的加权值,
Q=diag[1000·γ1 3,1000·γ2 3,…,1000·γm 3],diag表示主对角矩阵;γ为各干涉对在该点的相干性。
所述SIP方法是一种利用单个InSAR干涉对获取矿区地表三维形变场的方法,基于矿区地表水平移动与下沉梯度的比例关系,利用单个InSAR干涉对获取的一维雷达视线向形变获取了矿区地表在该时间段内沿着垂直、东西和南北方向的三维差分形变。
有益效果
本发明公开了一种基于单个雷达成像几何学SAR影像的矿区三维时序形变监测方法,利用单个雷达成像几何学SAR数据生成可用InSAR干涉对;利用现有的基于单个InSAR干涉对的矿区地表三维形变获取方法处理每一个可用InSAR干涉对,从而生成矿区地表多时域下沉观测值;建立时间相邻SAR影像期间的下沉速度与多时域下沉观测值之间的观测方程,并使用加权最小二乘法求解下沉速度;使用解出的下沉速度计算矿区地表在SAR影像获取时间的时序下沉;基于矿区地表东西、南北方向水平移动与下沉在对应方向的梯度成比例的关系,使用计算的矿区时序下沉估计出东西、南北方向的时序形变。本发明实现了仅利用单个雷达成像几何学SAR数据监测矿区地表三维时序形变,大大降低了传统方法中需要三个或以上的不同雷达成像几何学SAR数据的苛刻要求,该方法构思巧妙,数据制约少,监测结果准确有效,大大拓宽了InSAR技术的应用前景。通过加权最小二乘和大量的观测值,有效地抑制了InSAR视线向形变误差传播到时序下沉中。此外,本发明利用误差抑制后的时序下沉计算东西、南北方向的时序水平移动,其比直接使用SPI方法估计获得的结果的精度高。因而能有效地降低了水平移动中的误差传播,提高了三维时序形变的精度和可靠性。
因此,本方法其对于拓宽InSAR应用空间,降低矿区三维时序形变监测成本,提高InSAR矿区三维形变监测精度有着重要意义。此外,其对于指导矿区安全生产、预警矿区地表地质灾害以及生态环境保护也起着重要作用。
附图说明
图1本发明的所述方法的流程示意图;
图2表示走向长壁工作面模拟示意图;
图3表示垂直Ws、东西Es和南北Ns方向的三维时序形变的模拟示意图,其中Ws0=0,Es0=0,Ns0=0;
图4表示InSAR监测LOS向形变模拟示意图,其中,dLOS(ti-tj)表示该干涉对由获取时间为ti和tj的两景SAR影像生成;
图5为采用本发明所述方法估计的垂直W、东西E和南北N方向三维时序形变示意图,其中,W0=0,E0=0,N0=0。
具体实施方式
下面将结合附图1对本发明做进一步的说明。
步骤1:利用单个雷达成像几何学SAR影像生成可用InSAR干涉对;
假设有n+1幅覆盖待监测矿区的单个雷达成像几何学SAR影像,其获取时间分别为(t0,t1,…,tn)(其中t0<t1<…<tn)。
设定InSAR干涉对的时间基线和空间基线阈值,并利用SAR影像生成时间基线和空间基线均小于对应阈值的可用InSAR干涉对。假定可用InSAR干涉对的数量m,且覆盖整个时域过程(即从时间t0到tn期间均有InSAR干涉对覆盖)。令IM和IS分别表示组成m个可用InSAR干涉对的时间索引,即
其中,且比如第k个可用InSAR干涉对由时间为分别为t4和t5的两幅主、辅SAR影像生成,则IMk=4且ISk=5。
步骤2:生成多时域差分下沉观测值;
随后,利用现有的方法“一种利用单个InSAR干涉对获取矿区地表三维形变场的方法”(以下简称SIP方法)处理每一个可用InSAR干涉对,得到m个沿着垂直方向的多时域差分下沉观测值ΔW=[ΔW1,ΔW2,…,ΔWm]T。
步骤3:建立矿区地表各点在两相邻SAR影像时刻间的平均下沉速度V=[V1,V2,…,Vn]T与多时域差分下沉观测值ΔW=[ΔW1,ΔW2,…,ΔWm]T之间的观测方程;
令像素坐标为(i,j)的地表点在时间相邻SAR影像期间的平均下沉速度为V(i,j)=[V1(i,j),V2(i,j),…,Vn(i,j)]T,即
其中,Wk(i,j)和Wk-1(i,j)表示点(i,j)在tk和tk-1时刻相对于t0时刻(即W0(i,j)=0)的时序下沉值。因此,对于在点(i,j)在SAR影像获取时刻的时序沉降值可表示为:
从式(3)可以看出,若想获取矿区地表点(i,j)处的时序沉降值,其在时间相邻SAR影像期间的平均速度V(i,j)=[V1(i,j),V2(i,j),…,Vn(i,j)]T首先需要被估计。也就是说,对于矿区地表任意一点,其有n个未知数需要估计。理论上,若能够组建至少n个关于未知数的线性无关方程组,即可解出这n个未知数,从而估计出矿区该点在SAR影像获取时刻的下沉值。在矿区地表点(i,j)上,第k个SIP方法获取的下沉观测值ΔWk(i,j)可以表示为:
从式(4)中可以看出,对于每个SIP方法获取的差分下沉值ΔWk(i,j)(k=1,2,…,m)均可建立如式(4)的线性方程。若可用InSAR干涉对覆盖整个时域过程(即m≥n),即可解出点(i,j)的n个未知数。为了便于表示,将式(4)表示为矩阵形式:
Bm×n·V(i,j)n×1=ΔW(i,j)m×1 (5)
其中B为m×n的系数矩阵,其形式取决于可用InSAR干涉对的组成,对于任意第k行的形式为:
其余元素均为0。比如ΔW1(i,j)=W1(i,j)-W0(i,j)、ΔW2(i,j)=W3(i,j)-W0(i,j),则式(5)可表示为:
假设有n+1景SAR影像,其获取时间按照先后顺序排序为t=[t1,t2,…,tn],生成了m个可用InSAR干涉对。令IM=[IM1,IM2,…,IMm]为InSAR干涉对的主影像时间索引,IS=[IS1,IS2,…,ISm]为辅影像时间索引。所述时间索引为影像获取时间的下标。
1)假设干涉对1是时间为t1和t2的SAR影像生成,则ΔW1为t1和t2之间的差分下沉观测值,其等于ΔW1=v1·(t2-t1);
2)假设干涉对2是时间为t1和t3的SAR影像生成,则ΔW2为t1和t3之间的差分下沉观测值,其等于ΔW1=v1·(t2-t1)+v2·(t3-t2)(两个时间段内的沉降相加);
3)假设干涉对3是时间为t2和t4的SAR影像生成,则ΔW3为t2和t4之间的差分下沉观测值,其等于ΔW3=v2·(t3-t2)+v3·(t4-t3)(两个时间段内的沉降相加);
由于第1-3个InSAR干涉对分别有时间为t1-t2,t1-t3和t2-t4的SAR影像生成,所以就有IM1=1,IM2=1,IM3=2,而IS1=2,IS2=3,IS3=4;
把上面三个假设写成矩阵的形式即为:B·V=ΔW
步骤4:求解地表各点在时间相邻SAR影像时刻的平均速度;
式(6)有m个观测方程和n个未知数,若可用InSAR干涉对覆盖了整个时域过程(即m≥n),则可利用最小二乘法解出n个未知数。然而,由于不同InSAR干涉对的噪声水平不同,所以SIP方法获取的下沉观测值ΔWk(k=1,2,…,M)的精度也是不同的。
为了尽可能抑制SIP方法获取的差分下沉值ΔWk中误差传播到估计的未知数中,本发明选用加权最小二乘法估计未知数V(i,j)。
考虑到相干性γ是评价InSAR干涉对获取的LOS向形变精度的重要参考指标,本发明设计了基于相干性的加权函数Q,即:
Qk(i,j)=[10·γk(i,j)]3,(k=1,2,…,m) (7)
其中γk(i,j)为第k个InSAR干涉对在点(i,j)处的相干性,Qk(i,j)为ΔWk的加权值。基于该加权函数,即解出未知数的加权最小二乘解:
V(i,j)=[BTQ(i,j)B]-1·[BTQ(i,j)ΔW(i,j)] (8)
式中,Q(i,j)是对角线元素为Q1(i,j),Q2(i,j),…,QM(i,j)的主对角矩阵。
步骤5:估计矿区地表时序下沉;
在解出点(i,j)处的下沉速度V(i,j)后,即可利用公式(4)计算出该点的时序沉降值。对于矿区其它地表点重复以上步骤,即可估计出整个待监测矿区在SAR影像获取时刻的地表时序下沉值W=[W1,W2,…,Wn]。
步骤6:计算东西、南北方向的时序水平移动;
鉴于矿区地表水平在东西E、南北方向的水平移动N与对应方向的下沉梯度之间是存在比例关系;即点(i,j)在tk(k=1,2,…,n)时刻的东西E(i,j)、南北方向的水平移动N(i,j)按以下公式表示:
式中,b为待监测矿区水平移动系数,H为开采采深,tanβ为主要影响角正切,RE和RN分别为LOS向形变图在东西和南北方向的空间分辨率。
至此,矿区地表在SAR影像获取时刻沿着东西E=[E1,E2,…,En]、南北N=[N1,N2,…,Nn]和垂直三个方向W=[W1,W2,…,Wn]的三维时序形变。
本发明使用模拟数据进行可行性和精度论证,假设有一走向长壁开采工作面,其开采深度为700m,倾向宽度为150m,倾角0°,至西向东以平均速度3.3m/天推进(如图2所示,其中虚线表示各SAR影像获取时刻工作面的推进位置)。假设覆盖该矿区的可用SAR数据为6景升轨的ALOS PALSAR数据(单一雷达成像几何学),其时间间隔均为46天,即(t0,t1,…,t5)=(0,46,92,138,184,230天)。采用快速拉格朗日分析软件FLAC3D(Fast LagrangianAnalysis of Continua)模拟该工作面在SAR影像期间的地表在垂直、东西和南北方向的三维时序形变(其结果如图3所示)Ws=[Ws0,Ws1,Ws2,…,Ws5],Es=[Es0,Es1,Es2,…,Es5],Ns=[Ns0,Ns1,Ns2,…,Ns5]。由于在t0=0时刻,地下开采刚开始,所以Ws0=0,Es0=0,Ns0=0。通过分析三维时序形变后发现在该地质采矿条件下的走向和倾向水平移动系数为b=0.3,主要影响角正切为tanβ=1.6。
将6景SAR影像中任意两景组成InSAR干涉对,并设置时间基线阈值为180天。剔除时间基线大于时间基线阈值的InSAR干涉对,得到12对可用InSAR干涉对(干涉对构成如图4所示)。随后,按照可用InSAR干涉对的主辅影像时间,并将对应的FLAC3D模拟的三维形变投影到LOS方向(如图4所示),从而生成12个模拟的InSAR监测的LOS向形变。比如,第一个可用InSAR干涉对有第一景t0和第三景t2SAR影像组成,基于FLAC3D模拟的三维时序形变可计算出0天到92天之间地表的差分三维形变,即ΔW1=Ws1-Ws3,ΔE1=Es1-Es3和ΔN1=Ns1-Ns3。将三维差分形变投影至LOS方向即可获得模拟的InSAR干涉对监测的LOS向形变dLOS(t0-t2),即:
式中θ为雷达入射角,αh为雷达飞行方位角,在本实施例中θ=38.7°,αh=349.7°。
假设LOS向干涉图中各点的相干性相等,且均为0.6,即γ=0.6。利用本发明的方法处理12个模拟的LOS向形变图,从而得到在6景SAR影像获取时刻矿区地表的三维时序形变W=[W0,W1,…,W5],Es=[E0,E1,…,E5],Ns=[N0,N1,…,N5],其结果如图5所示。通过对比本发明所述方法获取的三维时序形变和模拟的三维时序形变后得出:两者吻合较好,且在垂直、东西和南北方向的中误差分别为0.006、0.007和0.009m。该结果表明本发明所述方法是可行的,且结果较为可靠。
以上所述仅为本发明的优选实施例而已,并不用于限制本发明,对于本领域的技术人员来说,本发明可以有各种更改和变化。凡在本发明的精神和原则之内,所作的任何修改、等同替换、改进等,均应包含在本发明的保护范围之内。
Claims (3)
1.一种基于单个雷达成像几何学SAR影像的矿区三维时序形变监测方法,其特征在于,包括以下几个步骤:
步骤1:利用单个雷达成像几何学SAR影像生成覆盖待监测矿区的可用InSAR干涉对;
步骤2:利用现有的SIP方法生成待监测矿区多时域差分下沉观测值ΔW:
ΔW=[ΔW1,ΔW2,…,ΔWm]T;
其中,m表示可用InSAR干涉对的个数;ΔW1,ΔW2,…,ΔWm分别表示第1到第m个可用InSAR干涉对生成的差分下沉观测值;
所述多时域差分下沉观测值指一系列与可用InSAR干涉对时间同步的下沉值;
步骤3:建立矿区地表各点在两相邻SAR影像时刻间的平均下沉速度V=[V1,V2,…,Vn]T与多时域差分下沉观测值ΔW=[ΔW1,ΔW2,…,ΔWm]T之间的函数关系:Bm×n·Vn×1=ΔWm×1;
其中,SAR影像个数为n+1;V1,V2,…,Vn表示矿区地表各点依次在n+1个影像中,两相邻SAR影像时刻间的平均下沉速度;B是一个维度为m×n,由可用InSAR干涉对辅影像和主影像获取时间差组成的系数矩阵;
对于B的任意第k行,第IMk个元素前的所有元素均为0,从第IMk到第ISk-1个元素,依次为:第ISk-1个元素以后的所有元素均为0;
其中,IMk和ISk分别为生成第k个可用InSAR干涉对的主、辅影像获取时间索引,根据InSAR干涉对的组成情况获得;
步骤4:求解地表各点在相邻SAR影像时刻间的平均下沉速度V:V=[BT·B]-1·[BT·ΔW];
步骤5:估计矿区地表时序下沉W,并利用时序下沉计算获得东西、南北方向的时序水平移动,完成待监测矿区三维时序形变监测;
其中,Wk表示在tk时刻矿区地表的时序下沉值,tl表示时间,Vl表示在tl-1至tl期间的平均下沉速度;
其中,E(i,j)和N(i,j)分别表示像素坐标为(i,j)的点在tk时刻的东西、南北方向的水平移动,k=1,2,…,n;b为待监测矿区水平移动系数,H为开采采深,tanβ为主要影响角正切,RE和RN分别为LOS向形变图在东西和南北方向的空间分辨率。
2.根据权利要求1所述的方法,其特征在于,所述求解地表各点在时间相邻SAR影像时刻的平均速度V时采用加权最小二乘法求解:
V=[BT·Q·B]-1·[BT·Q·ΔW]
其中,Q表示地表各点多时域差分下沉值的加权矩阵,Q=diag[1000·γ1 3,1000·γ2 3,…,1000·γm 3],diag表示主对角矩阵;γ为各干涉对在该点的相干性。
3.根据权利要求1或2所述的方法,其特征在于,所述SIP方法是一种利用单个InSAR干涉对获取矿区地表三维形变场的方法,基于矿区地表水平移动与下沉梯度的比例关系,利用单个InSAR干涉对获取的一维雷达视线向形变获取了矿区地表在该时间段内沿着垂直、东西和南北方向的三维差分形变。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201610546270.1A CN106226767B (zh) | 2016-07-12 | 2016-07-12 | 基于单个雷达成像几何学sar影像的矿区三维时序形变监测方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201610546270.1A CN106226767B (zh) | 2016-07-12 | 2016-07-12 | 基于单个雷达成像几何学sar影像的矿区三维时序形变监测方法 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN106226767A CN106226767A (zh) | 2016-12-14 |
CN106226767B true CN106226767B (zh) | 2018-08-21 |
Family
ID=57520564
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN201610546270.1A Active CN106226767B (zh) | 2016-07-12 | 2016-07-12 | 基于单个雷达成像几何学sar影像的矿区三维时序形变监测方法 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN106226767B (zh) |
Families Citing this family (7)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN106767380B (zh) * | 2017-01-19 | 2019-04-19 | 中南大学 | 一种基于两景sar强度影像的矿区地表大量级三维形变估计方法 |
CN107144213A (zh) * | 2017-06-29 | 2017-09-08 | 中南大学 | 基于sar强度影像的矿区大量级三维时序形变估计方法及装置 |
CN108763822B (zh) * | 2018-06-15 | 2022-03-08 | 安徽理工大学 | 一种基于沉陷监测的煤矿采空区空间几何特征精准识别方法 |
CN109738892B (zh) * | 2019-01-24 | 2020-06-30 | 中南大学 | 一种矿区地表高时空分辨率三维形变估计方法 |
CN110058236B (zh) * | 2019-05-21 | 2023-04-07 | 中南大学 | 一种面向三维地表形变估计的InSAR和GNSS定权方法 |
CN111076704B (zh) * | 2019-12-23 | 2022-05-20 | 煤炭科学技术研究院有限公司 | 一种利用insar精确解算采煤沉陷区地表下沉量的方法 |
CN114777633B (zh) * | 2022-04-02 | 2024-05-14 | 山西省煤炭地质勘查研究院有限公司 | 一种关闭矿区阶段性形变监测及分析方法 |
Family Cites Families (3)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN102927934B (zh) * | 2012-11-07 | 2015-01-28 | 中南大学 | 一种利用单个InSAR干涉对获取矿区地表三维形变场的方法 |
CN103091675B (zh) * | 2013-01-11 | 2014-07-30 | 中南大学 | 一种基于InSAR技术的矿区开采监测方法 |
CN103822598B (zh) * | 2014-02-26 | 2016-05-25 | 北京理工大学 | 地基sar在时间去相关严重区域的形变监测方法 |
-
2016
- 2016-07-12 CN CN201610546270.1A patent/CN106226767B/zh active Active
Also Published As
Publication number | Publication date |
---|---|
CN106226767A (zh) | 2016-12-14 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN106226767B (zh) | 基于单个雷达成像几何学sar影像的矿区三维时序形变监测方法 | |
CN106526590B (zh) | 一种融合多源sar影像工矿区三维地表形变监测及解算方法 | |
CN102927934B (zh) | 一种利用单个InSAR干涉对获取矿区地表三维形变场的方法 | |
Perski et al. | InSAR analyses of terrain deformation near the Wieliczka Salt Mine, Poland | |
Samieie-Esfahany et al. | On the effect of horizontal deformation on InSAR subsidence estimates | |
Tamburini et al. | Retrieving surface deformation by PSInSAR™ technology: A powerful tool in reservoir monitoring | |
Casson et al. | Contribution of multi-temporal remote sensing images to characterize landslide slip surface‒Application to the La Clapière landslide (France) | |
AU2021439678A9 (en) | Method for extracting three-dimensional surface deformation by combining unmanned aerial vehicle doms and satellite-borne sar images | |
CN106767380B (zh) | 一种基于两景sar强度影像的矿区地表大量级三维形变估计方法 | |
Hu et al. | Kalman-filter-based approach for multisensor, multitrack, and multitemporal InSAR | |
Montazeri et al. | Three-dimensional deformation monitoring of urban infrastructure by tomographic SAR using multitrack TerraSAR-X data stacks | |
CN108983232B (zh) | 一种基于邻轨数据的InSAR二维地表形变监测方法 | |
Tobita et al. | 3‐D surface deformation of the 2000 Usu eruption measured by matching of SAR images | |
Ren et al. | Calculating vertical deformation using a single InSAR pair based on singular value decomposition in mining areas | |
Wang et al. | 3D coseismic deformations and source parameters of the 2010 Yushu earthquake (China) inferred from DInSAR and multiple-aperture InSAR measurements | |
Ng et al. | Estimating horizontal and vertical movements due to underground mining using ALOS PALSAR | |
Wang et al. | A method of monitoring three-dimensional ground displacement in mining areas by integrating multiple InSAR methods | |
Zhao et al. | Pre-, co-, and post-rockslide analysis with ALOS/PALSAR imagery: A case study of the Jiweishan rockslide, China | |
CN112233232B (zh) | 一种基于单轨InSAR观测的三维地壳形变转换方法 | |
CN105824022A (zh) | 一种电网不良地质体三维形变监测方法 | |
Sarychikhina et al. | Application of satellite SAR interferometry for the detection and monitoring of landslides along the Tijuana-Ensenada Scenic Highway, Baja California, Mexico | |
Huang et al. | Time-series SBAS pixel offset tracking method for monitoring three-dimensional deformation in a mining area | |
Ittycheria et al. | Time series analysis of surface deformation of Bengaluru city using Sentinel-1 images | |
Liang et al. | Mapping surface deformation over Tatun volcano group, northern Taiwan using multitemporal insar | |
Suncar et al. | Deformations of a rapidly moving landslide from high-resolution optical satellite imagery |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
C06 | Publication | ||
PB01 | Publication | ||
SE01 | Entry into force of request for substantive examination | ||
SE01 | Entry into force of request for substantive examination | ||
GR01 | Patent grant | ||
GR01 | Patent grant |