CN110120046A - 一种融合dem、光学遥感和形变信息的潜在滑坡识别方法 - Google Patents

一种融合dem、光学遥感和形变信息的潜在滑坡识别方法 Download PDF

Info

Publication number
CN110120046A
CN110120046A CN201910237520.7A CN201910237520A CN110120046A CN 110120046 A CN110120046 A CN 110120046A CN 201910237520 A CN201910237520 A CN 201910237520A CN 110120046 A CN110120046 A CN 110120046A
Authority
CN
China
Prior art keywords
landslide
dem data
dem
area
data
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Granted
Application number
CN201910237520.7A
Other languages
English (en)
Other versions
CN110120046B (zh
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.)
Changan University
Original Assignee
Changan 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 Changan University filed Critical Changan University
Priority to CN201910237520.7A priority Critical patent/CN110120046B/zh
Publication of CN110120046A publication Critical patent/CN110120046A/zh
Application granted granted Critical
Publication of CN110120046B publication Critical patent/CN110120046B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

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/9004—SAR image acquisition techniques
    • G01S13/9005—SAR image acquisition techniques with optical processing of the SAR signals
    • 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
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F18/00—Pattern recognition
    • G06F18/20—Analysing
    • G06F18/24—Classification techniques
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00—Image analysis
    • G06T7/10—Segmentation; Edge detection
    • G06T7/11—Region-based segmentation
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06V—IMAGE OR VIDEO RECOGNITION OR UNDERSTANDING
    • G06V20/00—Scenes; Scene-specific elements
    • G06V20/10—Terrestrial scenes
    • G06V20/13—Satellite images
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00—Indexing scheme for image analysis or image enhancement
    • G06T2207/10—Image acquisition modality
    • G06T2207/10032—Satellite or aerial image; Remote sensing
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00—Indexing scheme for image analysis or image enhancement
    • G06T2207/10—Image acquisition modality
    • G06T2207/10032—Satellite or aerial image; Remote sensing
    • G06T2207/10036—Multispectral image; Hyperspectral image
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00—Indexing scheme for image analysis or image enhancement
    • G06T2207/20—Special algorithmic details
    • G06T2207/20048—Transform domain processing
    • G06T2207/20064—Wavelet transform [DWT]
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00—Indexing scheme for image analysis or image enhancement
    • G06T2207/30—Subject of image; Context of image processing
    • G06T2207/30181—Earth observation

Landscapes

  • Engineering & Computer Science (AREA)
  • Remote Sensing (AREA)
  • Radar, Positioning & Navigation (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Computer Networks & Wireless Communication (AREA)
  • Electromagnetism (AREA)
  • Data Mining & Analysis (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Evolutionary Biology (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Artificial Intelligence (AREA)
  • Bioinformatics & Cheminformatics (AREA)
  • Bioinformatics & Computational Biology (AREA)
  • Signal Processing (AREA)
  • Evolutionary Computation (AREA)
  • General Engineering & Computer Science (AREA)
  • Astronomy & Astrophysics (AREA)
  • Multimedia (AREA)
  • Pit Excavations, Shoring, Fill Or Stabilisation Of Slopes (AREA)
  • Testing Or Calibration Of Command Recording Devices (AREA)

Abstract

本发明属于潜在滑坡识别领域，提供了一种融合DEM、光学遥感和形变信息的潜在滑坡识别方法，具体包括：对研究区域内的SAR影像处理得到重采样后的形变速率图；获取研究区域的光学遥感影像数据和DEM数据，根据得到的研究区域的光学遥感影像数据和DEM数据，得到多个对象；根据重采样后的形变速率图和多个对象，处理得到利用面向对象分类技术获取的疑似滑坡区；选取研究区域内已知滑坡区和已知非滑坡区，得到优选阈值对应的疑似滑坡区域，与利用面向对象分类技术获取的疑似滑坡区进行结合，得到潜在滑坡区域；本发明操作简单，得到的疑似滑坡区域在空间上更加连续完整，有效减少了滑坡区域的错判漏判。

Description

一种融合DEM、光学遥感和形变信息的潜在滑坡识别方法
技术领域
本发明属于潜在滑坡识别领域，具体涉及一种融合DEM、光学遥感和形变信息的潜在滑坡识别方法。
背景技术
滑坡是一种常见的工程地质问题导致的自然灾害，现在进行滑坡的早期识别主要是基于常规野外图幅调查或者常规工程勘察，这些方法主要依赖于专业人员的经验，不仅耗时费力，覆盖面有限，还需要投入大量的人力资源。而光学遥感解译技术识别滑坡，尽管其覆盖范围大，且可以通过人机交互进行识别滑坡；但是，常规的光学遥感调查往往只能提供半定量滑坡识别成果，存在漏判、错判现象，自动化程度低，且多数是针对滑坡之后的定位，难以对潜在的滑坡进行早期识别，导致很多灾害难以得到有效的防治。
SAR影像不仅含有强度信息，还包含相位信息，InSAR技术可以获取厘米级甚至毫米级的地表形变，从而提高滑坡识别与监测的可靠性与正确率；但是它只能获取一维形变而且存在阴影、顶底倒置和透视缩短等现象，因此也会存在对滑坡的漏判、错判现象，特别是对于空间尺度较小的滑坡，基于InSAR形变直接进行识别往往不成功。高精度DEM数据可以客观定量的分析精细尺度下的地表特征进行滑坡识别，但也存在漏判、错判现象。InSAR技术和光学遥感融合的方法中，需要构建复杂的分类规则集，且其是在大范围形变的基础上探测疑似滑坡，然后以单个疑似滑坡为单位进行提取，不利于对大范围的潜在滑坡进行快速识别，仍旧错判、漏判的现象。
发明内容
针对现有技术中存在错判、漏判的缺点，本申请的目的在于，提供了一种融合DEM、光学遥感和形变信息的潜在滑坡识别方法。
为了实现上述目的，本发明采取以下技术方案予以实现：
一种融合DEM、光学遥感和形变信息的潜在滑坡识别方法，具体包括以下步骤：
步骤1、对研究区域内的SAR影像进行处理，得到研究区域的形变速率图；获取研究区域的光学遥感影像数据和DEM数据，所述光学遥感影像数据包括全色波段数据和4波段多光谱数据；通过DEM数据对光学遥感影像数据进行预处理，得到融合的光学遥感影像图；所述DEM指数字高程模型；
步骤2、对研究区域的形变速率图定义投影和重采样，得到重采样后的形变速率图；对所述DEM数据定义投影，得到具有投影信息的DEM数据；采用多尺度分割方法对具有投影信息的DEM数据和融合的光学遥感影像图中的4波段多光谱数据整体进行分割，得到多个对象；
还包括以下步骤：
步骤3、将重采样后的形变速率图与得到的多个对象进行叠加，得到融合的光学遥感影像图中各像元的形变速率；通过融合的光学遥感影像图中各像元的形变速率计算得到多个对象的形变速率；采用阈值分类的方法对多个对象的形变速率进行处理，得到疑似滑坡区；
步骤4、在研究区域内分别选取已知滑坡区和已知非滑坡区，分别沿已知滑坡区和已知非滑坡区的边界，对具有投影信息的DEM数据进行剪裁，得到已知滑坡区的DEM数据和已知非滑坡区的DEM数据；给定多个尺度因子，计算多个尺度因子下已知滑坡区的DEM数据的小波系数的方差和多个尺度因子下已知非滑坡区的DEM数据的小波系数的方差；通过多个尺度因子下已知非滑坡区的DEM数据的小波系数的方差，对多个尺度因子下已知滑坡区DEM数据的小波系数的方差进行归一化处理，得到多个尺度因子下已知滑坡区的DEM数据归一化后的小波系数方差；判断得到的多个尺度因子下已知滑坡区的DEM数据归一化后的小波系数方差的大小，得到特征尺度因子；采用特征尺度因子对DEM数据进行小波变换，得到DEM数据中每个节点在多个整数尺度因子下的小波系数的平方和；
步骤5、设定多个给定阈值，判断多个整数尺度因子下DEM数据中每个节点的小波系数的平方和与多个给定阈值的大小，得到多个给定阈值对应的疑似滑坡区域；采用误差分析方法对得到的多个给定阈值对应的疑似滑坡区域进行处理，得到优选阈值；从多个给定阈值对应的疑似滑坡区域中选取优选阈值对应的疑似滑坡区域，将优选阈值对应的疑似滑坡区域与疑似滑坡区进行交运算，得到潜在滑坡区域。
进一步的，所述步骤3的具体操作如下：
将重采样后的形变速率图与得到的多个对象进行叠加，得到融合的光学遥感影像图中各像元的形变速率；分别计算多个对象中每个对象的所有像元形变速率的平均值，将得到的平均值分别作为多个对象的形变速率；采用阈值分类的方法对多个对象的形变速率进行分类，得到阈值范围内的形变速率，从多个对象中提取阈值范围内的形变速率所对应的对象，将提取的阈值范围内的形变速率所对应的对象作为疑似滑坡区。
进一步的，所述步骤4的具体操作如下：
步骤4.1、在研究区域内分别选取已知滑坡区和已知非滑坡区，分别沿已知滑坡区和已知非滑坡区的边界，对具有投影信息的DEM数据进行剪裁，分别得到已知滑坡区的DEM数据和已知非滑坡区的DEM数据；给定多个尺度因子，计算得到多个尺度因子下已知滑坡区的DEM数据的小波系数和多个尺度因子下已知非滑坡区的DEM数据的小波系数；
步骤4.2、计算多个尺度因子下已知滑坡区的DEM数据的小波系数的方差，和多个尺度因子下已知非滑坡区的DEM数据的小波系数的方差；
步骤4.3、通过多个尺度因子下已知非滑坡区的DEM数据的小波系数的方差，对多个尺度因子下已知滑坡区的DEM数据的小波系数的方差进行归一化处理，得到多个尺度因子下已知滑坡区的DEM数据归一化后的小波系数方差；将多个尺度因子下已知滑坡区域DEM数据归一化后小波系数方差中的最大值的1/2作为截断值；提取多个尺度因子下已知滑坡区域DEM数据归一化小波系数方差中大于截断值的已知滑坡区域DEM数据归一化小波系数方差，将大于截断值的已知滑坡区域DEM数据归一化小波系数方差所对应的尺度因子作为特征尺度因子，所述特征尺度因子中包括多个整数尺度因子；
步骤4.4、采用特征尺度因子中的多个整数尺度因子对DEM数据进行小波变换，得到多个整数尺度因子下的DEM数据的小波系数；通过多个整数尺度因子下的DEM数据的小波系数，计算DEM数据中每个节点在多个整数尺度因子下的小波系数的平方和。
进一步的，所述步骤5的具体操作如下：
步骤5.1、设定多个给定阈值，对DEM数据中的每个节点和多个给定阈值分别进行编号，得到DEM数据节点序列和给定阈值序列；将DEM数据节点序列中的第一个节点作为当前节点，将给定阈值序列中的第一个给定阈值作为当前阈值；
步骤5.2、从DEM数据中每个节点在多个整数尺度因子下的小波系数的平方和中选取当前节点在多个整数尺度因子下的小波系数的平方和；
当当前节点在多个整数尺度因子下得到的小波系数的平方和大于等于当前阈值时，将当前节点作为当前阈值对应的疑似滑坡区域节点；将当前节点的下一个节点作为当前节点，重复执行与步骤5.2相同的操作，直至判断完最后一个节点结束，得到当前阈值对应的疑似滑坡区域的多个节点；根据得到的当前阈值的疑似滑坡区域的多个节点构成得到当前阈值对应的疑似滑坡区域；
当当前节点在多个整数尺度因子下得到的小波系数的平方和小于当前阈值时，表示当前节点属于非滑坡区域；将当前节点的下一个节点作为当前节点，重复执行与步骤5.2相同的操作，直至判断完DEM数据中所有节点操作完成时操作结束；
步骤5.3、将当前阈值的下一个给定阈值作为当前阈值，重复执行与步骤5.2相同的操作，直至判断完最后一个给定阈值结束，得到多个给定阈值对应的疑似滑坡区域；
步骤5.4、采用误差分析方法对得到的多个给定阈值对应的疑似滑坡区域进行处理，得到优选阈值；从多个给定阈值对应的疑似滑坡区域中选取优选阈值对应的疑似滑坡区域，将优选阈值对应的疑似滑坡区域与疑似滑坡区进行交运算，得到潜在滑坡区域。
进一步的，所述步骤4.1中，采用公式(1)计算得到多个尺度因子下已知滑坡区的DEM数据的小波系数，采用公式(2)计算得到多个尺度因子下已知非滑坡区的DEM数据的小波系数；
其中，(x,y)表示已知滑坡区域的DEM数据的节点，x和y均大于等于0；zfld(x,y)为已知滑坡区域的DEM数据中节点为(x,y)的高程，zfld(x,y)为实数；Cfld(si,a,b)为zfld(x,y)经过小波变换得到的小波系数，Cfld(si,a,b)∈R；si表示第i个尺度因子，si>0；为已知滑坡区的DEM数据中节点为(x,y)在第i个尺度因子下的小波函数值，a、b均表示小波变换的位移，a和b均为实数；
其中，(m,n)表示已知非滑坡区的DEM数据中的节点，m和n均大于等于0；zunfld(m,n)为已知非滑坡区的DEM数据中节点为(m,n)的高程，zunfld(m,n)为实数；Cunfld(si,a,b)为zunfld(m,n)经过小波变换得到的小波系数，Cunfld(si,a,b)∈R；si表示第i个尺度因子，si>0；为已知非滑坡区的DEM数据中节点为(m,n)在第i个尺度因子下的小波函数值，a、b均表示小波变换的位移，a和b均为实数。
进一步的，所述步骤4.2的具体操作如下：
通过多个尺度因子下已知滑坡区的DEM数据的小波系数，采用公式(3)计算得到多个尺度因子下已知滑坡区的DEM数据的小波系数的方差；
其中，VCWT_fld(si)表示第i个尺度因子时已知滑坡区域的DEM数据的小波系数的方差，VCWT_fld(si)≥0；Na×Nb表示DEM数据中的节点个数，Na,Nb分别表示DEM数据中每一行和每一列的节点个数，Na,Nb均大于0；
通过多个尺度因子下已知非滑坡区的DEM数据的小波系数，采用公式(4)计算得到多个尺度因子下已知非滑坡区的DEM数据的小波系数的方差；
其中，VCWT_unfld(si)表示第i个尺度因子时已知非滑坡区域的DEM数据的小波系数的方差，VCWT_unfld(si)≥0；Na,Nb表示DEM数据中的节点个数，Na,Nb分别表示DEM数据中每一行和每一列的节点个数，Na,Nb均大于0。
进一步的，所述步骤4.3中，通过多个尺度因子下已知非滑坡区的DEM数据的小波系数方差，采用公式(5)对多个尺度因子下已知滑坡区DEM数据的小波系数方差进行归一化处理，得到多个尺度因子下已知滑坡区的DEM数据归一化后的小波系数方差；
VCWT_norm(si)＝VCWT_fld(si)/VCWT_unfld(si) (5)
其中，VCWT_norm(si)表示第i个尺度因子的已知滑坡区的DEM数据的归一化小波系数方差，VCWT_fld(si)表示第i个尺度因子的已知滑坡区的DEM数据的小波系数的方差，VCWT_unfld(si)表示第i个尺度因子的已知非滑坡区的DEM数据的小波系数的方差，VCWT_norm(si),VCWT_fld(si)，VCWT_unfld(si)均大于等于0。
进一步的，所述步骤4.4中，采用特征尺度因子中的多个整数尺度因子对DEM数据进行小波变换，具体采用公式(6)计算得到多个整数尺度因子下的DEM数据的小波系数；
其中，(xj,yj)表示DEM数据中第j个节点，且j>0，xj与yj均为实数；z(xj,yj)为DEM数据中第j个节点的高程，z(xj,yj)为实数；C(s'i,a,b)为z(xj,yj)经过小波变换得到的小波系数，C(s'i,a,b)∈R；s'i表示特征尺度因子中第i个整数尺度因子，s'i>0，i>0；为DEM数据中第j个节点在第i个整数尺度因子下的小波函数值，且a、b均表示小波变换的位移，且a,b均为实数。
进一步的，所述步骤1的具体操作如下：
步骤1.1、对研究区域内的SAR卫星采集得到的SAR影像，采用小基线集技术对SAR影像进行处理，得到研究区域的的形变速率图；
步骤1.2、通过遥感卫星获取覆盖研究区域的光学遥感影像数据，所述光学遥感影像数据包括分辨率为0.5m的全色波段数据和分辨率为2m的4波段多光谱数据，所述4波段包括红、绿、蓝和近红外波段；
步骤1.3、利用无人机获取研究区域的分辨率为3m的DEM数据，利用光学遥感影像数据自带的RPC文件和RPC模型得到完整的RPC模型；利用DEM数据和完整的RPC模型分别对光学遥感影像数据进行正射校正，得到正射光学遥感影像，所述RPC指有理多项式系数；
步骤1.4、采用NNDiffuse Pan Sharpening算法对正射光学遥感影像进行融合，得到融合的光学遥感影像图。
进一步的，所述步骤2的具体操作如下：
步骤2.1、采用地图投影方法对研究区域的形变速率图定义投影，得到具有投影信息的形变速率图，使具有投影信息的形变速率图与融合的光学遥感影像图的投影坐标系统一致；采用三次卷积的方法对具有投影信息的形变速率图进行重采样，得到重采样后的形变速率图，使重采样后的形变速率图与融合的光学遥感影像图的分辨率保持一致；
步骤2.2、采用地图投影方法对所述DEM数据定义投影，得到具有投影信息的DEM数据，使具有投影信息的DEM数据与融合的光学遥感影像图的投影坐标系统保持一致；
步骤2.3、将具有投影信息的DEM数据作为一个波段数据，采用多尺度分割方法，对具有投影信息的DEM数据和融合的光学遥感影像图中的4波段多光谱数据整体进行分割，得到多个对象，所述对象是指具有同质性的像元集合。
与现有技术相比，本发明的有益效果如下：
1、本发明融合了高精度DEM、光学遥感和InSAR形变信息自动化识别潜在滑坡区域，以形变信息作为一特征因子参与光学遥感对潜在滑坡的识别提取，得到利用面向对象分类技术获取的疑似滑坡区，并利用二维连续小波变化基于高精度DEM自动绘制研究区域的疑似滑坡区域，将利用面向对象分类技术获取的疑似滑坡区与优选阈值对应的疑似滑坡区域进行融合，得到潜在滑坡区域，能够有效避免沉降区域，减少错判现象；
2、本发明无需野外调查，对于人员的地质专业知识要求较低，对于人员难以到达的区域也能进行研究；相较于传统的方法，本发明不需要建立规则集，更加简便、客观，自动化程度更高，得到的结果空间上更加连续、完整，有利于进一步的滑坡分析，为及时预防滑坡灾害的发生提供技术支撑。
附图说明
图1为本发明流程图；
图2为基于光学遥感影像及InSAR形变信息并利用OBIA技术获取的疑似滑坡的结果图；
图3为利用投影后的DEM数据通过小波变换自动绘制的疑似滑坡区域的结果图；
图4为本发明提取的潜在滑坡区域的结果图。
具体实施方式
如图1所示，本发明提供了一种融合DEM、光学遥感和形变信息的潜在滑坡识别方法，包括以下步骤：
步骤1、对研究区域内的SAR影像进行处理，得到研究区域的形变速率图；获取研究区域的光学遥感影像数据和DEM数据，通过DEM数据对光学遥感影像数据进行预处理，得到融合的光学遥感影像；所述DEM指数字高程模型；所述光学遥感影像数据包括全色波段数据和4波段多光谱数据；
步骤2、对研究区域的形变速率图定义投影和重采样，得到重采样后的形变速率图；对所述DEM数据定义投影，得到具有投影信息的DEM数据；采用多尺度分割方法对具有投影信息的DEM数据和融合的光学遥感影像中的4波段多光谱数据整体进行分割，得到多个对象；
步骤3、将重采样后的形变速率图与得到的多个对象进行叠加，得到融合的光学遥感影像图中各像元的形变速率；通过融合的光学遥感影像图中各像元的形变速率计算得到多个对象的形变速率根据；采用阈值分类的方法对多个对象的形变速率进行处理，得到疑似滑坡区(即利用面向对象分类技术获取的疑似滑坡区)；
步骤4、在研究区域内分别选取已知滑坡区和已知非滑坡区，分别沿已知滑坡区和已知非滑坡区的边界，对具有投影信息的DEM数据进行剪裁，得到已知滑坡区的DEM数据和已知非滑坡区的DEM数据；给定多个尺度因子，计算多个尺度因子下已知滑坡区的DEM数据的小波系数的方差和多个尺度因子下已知非滑坡区的DEM数据的小波系数的方差；通过多个尺度因子下已知非滑坡区的DEM数据的小波系数方差，对多个尺度因子下已知滑坡区DEM数据的小波系数方差进行归一化处理，得到多个尺度因子下已知滑坡区的DEM数据归一化后的小波系数方差；判断得到的多个尺度因子下已知滑坡区的DEM数据归一化后的小波系数方差的大小，得到特征尺度因子；采用特征尺度因子对DEM数据进行小波变换，得到DEM数据中每个节点在多个整数尺度因子下的小波系数的平方和；
步骤5、设定多个给定阈值，判断多个整数尺度因子下DEM数据中每个节点的小波系数的平方和与多个给定阈值的大小，得到多个给定阈值对应的疑似滑坡区域；采用误差分析方法对得到的多个给定阈值对应的疑似滑坡区域进行处理，得到优选阈值；从多个给定阈值对应的疑似滑坡区域中选取优选阈值对应的疑似滑坡区域，将优选阈值对应的疑似滑坡区域与疑似滑坡区进行交运算，得到潜在滑坡区域。
本发明融合了高精度DEM、光学遥感和InSAR形变信息自动化识别潜在滑坡区域，得到利用面向对象分类技术获取的疑似滑坡区；并利用二维连续小波变化基于高精度DEM自动绘制研究区域的疑似滑坡区域，得到优选阈值对应的疑似滑坡区域，将利用面向对象分类技术获取的疑似滑坡区与优选阈值对应的疑似滑坡区域进行融合，得到潜在滑坡区域，能够有效避免沉降区域，减少错判现象；本发明无需野外调查，对于人员的地质专业知识要求较低，对于人员难以到达的区域也能进行研究；相较于传统的方法，不需要建立规则集，更加简便、客观，自动化程度更高，得到的结果空间上更加连续、完整，有利于进一步的滑坡分析，为及时预防滑坡灾害的发生提供技术支撑。
具体的，所述步骤1的具体操作如下：
步骤1.1、对研究区域内的SAR卫星采集得到的SAR影像，采用SBAS(小基线集，Small Baseline Subsets，简称SBAS)技术对SAR影像进行处理，得到研究区域的的形变速率图。
步骤1.2、通过遥感卫星获取覆盖研究区域的光学遥感影像数据，所述光学遥感影像数据包括分辨率为0.5m的全色波段数据和分辨率为2m的4波段多光谱数据，所述4波段包括红、绿、蓝和近红外波段；
步骤1.3、利用无人机获取研究区域的分辨率为3m的DEM数据，利用光学遥感影像数据自带的RPC文件和RPC模型得到完整的RPC模型；利用DEM数据和完整的RPC模型分别对光学遥感影像数据中的全色波段数据与4波段多光谱数据进行正射校正，得到包括全色波段数据和4波段多光谱数据的正射光学遥感影像，所述RPC指有理多项式系数；
步骤1.4、采用NNDiffuse Pan Sharpening算法对包括全色波段数据和4波段多光谱数据的正射光学遥感影像进行融合，得到融合的光学遥感影像图。
该方式通过采集SAR影像、研究区域的光学遥感影像数据和DEM数据，处理得到研究区域的形变速率图与融合的光学遥感影像图。
具体的，所述步骤2的具体操作为：
步骤2.1、采用地图投影方法对研究区域的形变速率图定义投影，得到具有投影信息的形变速率图，使具有投影信息的形变速率图与融合的光学遥感影像图的投影坐标系统一致；采用三次卷积的方法对具有投影信息的形变速率图进行重采样，得到重采样后的形变速率图，使重采样后的形变速率图与融合的光学遥感影像图的分辨率保持一致；
步骤2.2、采用地图投影方法对所述DEM数据定义投影，得到具有投影信息的DEM数据，使具有投影信息的DEM数据与融合的光学遥感影像图的投影坐标系统保持一致；
步骤2.3、将具有投影信息的DEM数据作为一个波段数据，采用多尺度分割方法，对具有投影信息的DEM数据和融合的光学遥感影像图中的4波段多光谱数据整体进行分割，得到多个对象，所述对象是指具有同质性的像元集合。
该方式中通过对DEM数据及研究区域的形变速率图进行处理，得到多个对象及重采样后的形变速率图；能够采用形变信息这一要素进行分类，提取疑似滑坡区域。
具体的，所述步骤3的具体操作为：
将重采样后的形变速率图与得到的多个对象进行叠加，得到融合的光学遥感影像图中各像元的形变速率；分别计算多个对象中每个对象的所有像元的形变速率的平均值，将得到的平均值作为多个对象的形变速率；设置形变速率的两个阈值m、n，其中m<0,n>0，通过阈值分类的方法提取得到多个对象中对象的形变速率k在k≤m或k≥n的范围内的对象，将提取的对象作为利用OBIA技术获取的疑似滑坡区；
所述OBIA表示Object-based Image Analysis，面向对象分类。
该方式中，只利用形变信息这一要素进行分类，设置形变阈值直接提取疑似滑坡区域，不需要计算一系列的特征要素(地形特征、纹理特征、光谱特征等)，既不需要专家经验来决定分类特征要素，也不需要建立十分复杂的规则集对非滑坡地类进行逐一剔除，对工作人员的解译经验及专业知识没有苛刻的要求，操作步骤简洁明了，自动化程度高，十分适合于在大范围内进行疑似滑坡的识别和提取。
具体的，所述步骤4的具体操作为：
步骤4.1、在研究区域内分别选取已知滑坡区和已知非滑坡区，分别沿已知滑坡区和已知非滑坡区的边界，对具有投影信息的DEM数据进行剪裁，分别得到已知滑坡区的DEM数据和已知非滑坡区的DEM数据；给定多个尺度因子，采用公式(1)计算得到多个尺度因子下已知滑坡区的DEM数据的小波系数，采用公式(2)计算得到多个尺度因子下已知非滑坡区的DEM数据的小波系数；
其中，(x,y)表示已知滑坡区域的DEM数据的节点，x和y均大于等于0；zfld(x,y)为已知滑坡区域的DEM数据中节点为(x,y)的高程，zfld(x,y)为实数；Cfld(si,a,b)为zfld(x,y)经过小波变换得到的小波系数，Cfld(si,a,b)∈R；si表示第i个尺度因子，si>0；为已知滑坡区的DEM数据中节点为(x，y)在第i个尺度因子下的小波函数值，a、b均表示小波变换的位移，a和b均为实数；
其中，(m,n)表示已知非滑坡区的DEM数据中的节点，m和n均大于等于0；zunfld(m,n)为已知非滑坡区的DEM数据中节点为(m,n)的高程，zunfld(m,n)为实数；Cunfld(si,a,b)为zunfld(m,n)经过小波变换得到的小波系数，Cunfld(si,a,b)∈R；si表示第i个尺度因子，si>0；为已知非滑坡区的DEM数据中节点为(m,n)在第i个尺度因子下的小波函数值，a、b均表示小波变换的位移，a和b均为实数；
步骤4.2、通过多个尺度因子下已知滑坡区的DEM数据的小波系数，采用公式(3)计算得到多个尺度因子下已知滑坡区的DEM数据的小波系数的方差；通过多个尺度因子下已知非滑坡区的DEM数据的小波系数，采用公式(4)计算得到多个尺度因子下已知非滑坡区的DEM数据的小波系数的方差；
其中，VCWT_fld(si)表示第i个尺度因子时已知滑坡区域的DEM数据的小波系数的方差，VCWT_fld(si)≥0；Na,Nb表示DEM数据中的节点个数，Na,Nb分别表示DEM数据中每一行和每一列的节点个数，Na,Nb均大于0；
其中，VCWT_unfld(si)表示第i个尺度因子时已知非滑坡区域的DEM数据的小波系数的方差，VCWT_unfld(si)≥0；Na,Nb表示DEM数据中的节点个数，Na,Nb分别表示DEM数据中每一行和每一列的节点个数，Na,Nb均大于0；
步骤4.3、通过多个尺度因子下已知非滑坡区的DEM数据的小波系数方差，采用公式(5)将多个尺度因子下已知滑坡区DEM数据的小波系数方差进行归一化处理，得到多个尺度因子下已知滑坡区的DEM数据归一化后的小波系数方差，及多个尺度因子下已知滑坡区的DEM数据归一化后的小波系数方差与尺度因子之间的对应关系；将多个尺度因子下已知滑坡区域DEM数据归一化后小波系数方差中的最大值的1/2作为截断值，提取多个尺度因子下已知滑坡区域DEM数据归一化小波系数方差中大于截断值的已知滑坡区域DEM数据归一化小波系数方差，将大于截断值的已知滑坡区域DEM数据归一化小波系数方差所对应的尺度因子作为特征尺度因子，所述特征尺度因子中包括多个整数尺度因子；
VCWT_norm(si)＝VCWT_fld(si)/VCWT_unfld(si) (5)
其中，VCWT_norm(si)表示第i个尺度因子的已知滑坡区的DEM数据的归一化小波系数方差，VCWT_fld(si)表示第i个尺度因子的已知滑坡区的DEM数据的小波系数的方差，VCWT_unfld(si)表示第i个尺度因子的已知非滑坡区的DEM数据的小波系数的方差，VCWT_norm(si),VCWT_fld(si)，VCWT_unfld(si)均大于等于0。
步骤4.4、根据特征尺度因子中的多个整数尺度因子，采用公式(6)对DEM数据进行小波变换，得到多个整数尺度因子下的DEM数据的小波系数；通过多个整数尺度因子下的DEM数据的小波系数，计算得到DEM数据中每个节点在多个整数尺度因子下的小波系数的平方和；
其中，(xj,yj)表示DEM数据中第j个节点，且j>0，xj与yj均为实数；z(xj,yj)为DEM数据中第j个节点的高程，z(xj,yj)为实数；C(s'i,a,b)为z(xj,yj)经过小波变换得到的小波系数，C(s'i,a,b)∈R；s'i表示特征尺度因子中第i个整数尺度因子，s'i>0，i>0；为DEM数据中第j个节点在第i个整数尺度因子下的小波函数值，且a、b均表示小波变换的位移，且a,b均为实数。
具体的，所述步骤5具体操作如下：
步骤5.1、设定多个给定阈值，对DEM数据中的每个节点和多个给定阈值分别进行编号，得到DEM数据节点序列和给定阈值序列；将DEM数据节点序列中的第一个节点作为当前节点，将给定阈值序列中的第一个给定阈值作为当前阈值；
步骤5.2、从DEM数据中每个节点在多个整数尺度因子下的小波系数的平方和中选取当前节点在多个整数尺度因子下的小波系数的平方和；
当当前节点在多个整数尺度因子下得到的小波系数的平方和大于等于当前阈值时，将当前节点作为当前阈值对应的疑似滑坡区域节点；将当前节点的下一个节点作为当前节点，重复执行与步骤5.2相同的操作，直至判断完最后一个节点结束，得到当前阈值对应的疑似滑坡区域的多个节点；根据得到的当前阈值的疑似滑坡区域的多个节点构成得到当前阈值对应的疑似滑坡区域；
当当前节点在多个整数尺度因子下得到的小波系数的平方和小于当前阈值时，表示当前节点属于非滑坡区域；将当前节点的下一个节点作为当前节点，重复执行与步骤5.2相同的操作，直至判断完DEM数据中所有节点操作完成时操作结束；
步骤5.3、将当前阈值的下一个给定阈值作为当前阈值，重复执行与步骤5.2相同的操作，直至判断完最后一个给定阈值结束，得到多个给定阈值对应的疑似滑坡区域；执行步骤5.4；
步骤5.4、根据研究区域已知的滑坡编目图，采用误差分析方法对得到的多个给定阈值对应的疑似滑坡区域进行处理，得到优选阈值；从多个给定阈值对应的疑似滑坡区域中选取优选阈值所对应的疑似滑坡区域；将优选阈值所对应的疑似滑坡区域与疑似滑坡区进行交运算，得到潜在滑坡区域。
步骤4与步骤5中，基于高精度的DEM数据，利用二维连续小波变换进行疑似滑坡区域的绘制，得到优选阈值对应的疑似滑坡区域，即研究区域的疑似滑坡区域；高精度的DEM数据能够对滑坡的精细地表特征进行客观定量的分析；将利用面向对象分类技术获取的疑似滑坡区与优选阈值对应的疑似滑坡区域进行融合，得到潜在滑坡区域，能够有效避免沉降区域，减少错判现象。
实施例
本发明实验数据采用了真实的TerraSAR-X降轨数据和WorldView-02光学遥感影像数据。其中，TerraSAR-X降轨数据为19景，WorldView-02光学遥感影像数据为一景，包括分辨率为0.5m的全色数据与分辨率为2m的多光谱数据，采集日期为2015年11月29日；
实验过程：
步骤1、对研究区域内的SAR影像(TerraSAR-X降轨数据)进行处理，得到研究区域的形变速率图；获取研究区域的光学遥感影像数据和DEM数据，所述光学遥感影像数据包括全色波段数据和4波段多光谱数据；根据DEM数据对光学遥感影像数据进行预处理，得到融合的光学遥感影像图；所述DEM指数字高程模型；
步骤2、对研究区域的形变速率图定义投影和重采样，得到重采样后的形变速率图；对所述DEM数据定义投影，得到具有投影信息的DEM数据；采用多尺度分割方法对具有投影信息的DEM数据和融合的光学遥感影像图中的4波段多光谱数据整体进行分割，得到多个对象；
步骤3、将重采样后的形变速率图与得到的多个对象进行叠加，得到融合的光学遥感影像图中各像元的形变速率；通过融合的光学遥感影像图中各像元的形变速率计算多个对象的形变速率；采用阈值分类的方法对多个对象的形变速率进行处理，得到疑似滑坡区；
步骤4、在研究区域内分别选取已知滑坡区和已知非滑坡区，根据已知滑坡区和已知非滑坡区的边界，对具有投影信息的DEM数据进行剪裁，得到已知滑坡区的DEM数据和已知非滑坡区的DEM数据；给定多个尺度因子，计算多个尺度因子下已知滑坡区的DEM数据的小波系数的方差和多个尺度因子下已知非滑坡区的DEM数据的小波系数的方差；通过多个尺度因子下已知非滑坡区的DEM数据的小波系数的方差，对多个尺度因子下已知滑坡区DEM数据的小波系数的方差进行归一化处理，得到多个尺度因子下已知滑坡区的DEM数据归一化后的小波系数方差；判断得到的多个尺度因子下已知滑坡区的DEM数据归一化后的小波系数方差的大小，得到特征尺度因子；采用特征尺度因子对DEM数据进行小波变换，得到DEM数据中每个节点在多个整数尺度因子下的小波系数的平方和；
步骤5、设定多个给定阈值，判断多个整数尺度因子下DEM数据中每个节点的小波系数的平方和与多个给定阈值的大小，得到多个给定阈值对应的疑似滑坡区域；采用误差分析方法对得到的多个给定阈值对应的疑似滑坡区域进行处理，得到优选阈值；从多个给定阈值对应的疑似滑坡区域中选取优选阈值对应的疑似滑坡区域，将优选阈值对应的疑似滑坡区域与疑似滑坡区进行交运算，最终的得到如图4中用黑线绘制的潜在滑坡区域。
利用光学遥感影像及InSAR形变信息设置形变阈值的方法得到如图2中用黑线绘制的疑似滑坡区域；利用投影后的DEM数据通过小波变换自动绘制得到如图3中用黑线绘制的疑似滑坡区域。
通过图2、图4的对比可以看出，InSAR技术监测得到的形变不仅包含有潜在滑坡的形变，也包含有地面沉降的形变。仅利用光学遥感影像进行分割得到对象后利用形变信息提取得到的疑似滑坡区域会将沉降区包含在内，这是明显的错判现象；通过图3、图4的对比可以看出，基于高精度DEM数据利用小波变换进行滑坡区域的自动绘制所获取的疑似滑坡区域，可以有效避免沉降区域，但其覆盖范围广泛，存在大量的错判现象。
从图4中可以看出，本发明融合了高精度DEM数据、光学遥感影像和InSAR形变信息自动化识别潜在滑坡区域，以形变信息作为一特征因子参与光学遥感影像对潜在滑坡的识别提取，并利用二维连续小波变化基于高精度DEM数据自动绘制研究区域的疑似滑坡区域，以二者的交集作为最终的潜在滑坡区域，与滑坡编目图中识别出的滑坡具有很高的一致性，并可以在大范围内自动化识别潜在滑坡区域，不需要建立规则集，更加简便、客观，自动化程度更高，得到的结果空间上更加连续、完整，有利于进一步的滑坡分析，为及时预防滑坡灾害的发生提供技术支撑。
以上公开的仅为本发明的具体实施例，但是，本发明实施例并非局限于此，任何本领域的技术人员能思之的变化都应落入本发明的保护范围。

Claims (10)

1.一种融合DEM、光学遥感和形变信息的潜在滑坡识别方法，具体包括以下步骤：
步骤1、对研究区域内的SAR影像进行处理，得到研究区域的形变速率图；获取研究区域的光学遥感影像数据和DEM数据，所述光学遥感影像数据包括全色波段数据和4波段多光谱数据；通过DEM数据对光学遥感影像数据进行预处理，得到融合的光学遥感影像图；所述DEM指数字高程模型；
步骤2、对研究区域的形变速率图定义投影和重采样，得到重采样后的形变速率图；对所述DEM数据定义投影，得到具有投影信息的DEM数据；采用多尺度分割方法对具有投影信息的DEM数据和融合的光学遥感影像图中的4波段多光谱数据整体进行分割，得到多个对象；
其特征在于，还包括以下步骤：
步骤3、将重采样后的形变速率图与得到的多个对象进行叠加，得到融合的光学遥感影像图中各像元的形变速率；通过融合的光学遥感影像图中各像元的形变速率计算得到多个对象的形变速率；采用阈值分类的方法对多个对象的形变速率进行处理，得到疑似滑坡区；
步骤4、在研究区域内分别选取已知滑坡区和已知非滑坡区，分别沿已知滑坡区和已知非滑坡区的边界，对具有投影信息的DEM数据进行剪裁，得到已知滑坡区的DEM数据和已知非滑坡区的DEM数据；给定多个尺度因子，计算多个尺度因子下已知滑坡区的DEM数据的小波系数的方差和多个尺度因子下已知非滑坡区的DEM数据的小波系数的方差；通过多个尺度因子下已知非滑坡区的DEM数据的小波系数的方差，对多个尺度因子下已知滑坡区DEM数据的小波系数的方差进行归一化处理，得到多个尺度因子下已知滑坡区的DEM数据归一化后的小波系数方差；判断得到的多个尺度因子下已知滑坡区的DEM数据归一化后的小波系数方差的大小，得到特征尺度因子；采用特征尺度因子对DEM数据进行小波变换，得到DEM数据中每个节点在多个整数尺度因子下的小波系数的平方和；
步骤5、设定多个给定阈值，判断多个整数尺度因子下DEM数据中每个节点的小波系数的平方和与多个给定阈值的大小，得到多个给定阈值对应的疑似滑坡区域；采用误差分析方法对得到的多个给定阈值对应的疑似滑坡区域进行处理，得到优选阈值；从多个给定阈值对应的疑似滑坡区域中选取优选阈值对应的疑似滑坡区域，将优选阈值对应的疑似滑坡区域与疑似滑坡区进行交运算，得到潜在滑坡区域。
2.如权利要求1所述的融合DEM、光学遥感和形变信息的潜在滑坡识别方法，其特征在于，所述步骤3的具体操作如下：
将重采样后的形变速率图与得到的多个对象进行叠加，得到融合的光学遥感影像图中各像元的形变速率；分别计算多个对象中每个对象的所有像元形变速率的平均值，将得到的平均值分别作为多个对象的形变速率；采用阈值分类的方法对多个对象的形变速率进行分类，得到阈值范围内的形变速率，从多个对象中提取阈值范围内的形变速率所对应的对象，将提取的阈值范围内的形变速率所对应的对象作为疑似滑坡区。
3.如权利要求1所述的融合DEM、光学遥感和形变信息的潜在滑坡识别方法，其特征在于，所述步骤4的具体操作如下：
步骤4.1、在研究区域内分别选取已知滑坡区和已知非滑坡区，分别沿已知滑坡区和已知非滑坡区的边界，对具有投影信息的DEM数据进行剪裁，分别得到已知滑坡区的DEM数据和已知非滑坡区的DEM数据；给定多个尺度因子，计算得到多个尺度因子下已知滑坡区的DEM数据的小波系数和多个尺度因子下已知非滑坡区的DEM数据的小波系数；
步骤4.2、计算多个尺度因子下已知滑坡区的DEM数据的小波系数的方差，和多个尺度因子下已知非滑坡区的DEM数据的小波系数的方差；
步骤4.3、通过多个尺度因子下已知非滑坡区的DEM数据的小波系数的方差，对多个尺度因子下已知滑坡区的DEM数据的小波系数的方差进行归一化处理，得到多个尺度因子下已知滑坡区的DEM数据归一化后的小波系数方差；将多个尺度因子下已知滑坡区域DEM数据归一化后小波系数方差中的最大值的1/2作为截断值；提取多个尺度因子下已知滑坡区域DEM数据归一化小波系数方差中大于截断值的已知滑坡区域DEM数据归一化小波系数方差，将大于截断值的已知滑坡区域DEM数据归一化小波系数方差所对应的尺度因子作为特征尺度因子，所述特征尺度因子中包括多个整数尺度因子；
步骤4.4、采用特征尺度因子中的多个整数尺度因子对DEM数据进行小波变换，得到多个整数尺度因子下的DEM数据的小波系数；通过多个整数尺度因子下的DEM数据的小波系数，计算DEM数据中每个节点在多个整数尺度因子下的小波系数的平方和。
4.如权利要求1所述的融合DEM、光学遥感和形变信息的潜在滑坡识别方法，其特征在于，所述步骤5的具体操作如下：
步骤5.1、设定多个给定阈值，对DEM数据中的每个节点和多个给定阈值分别进行编号，得到DEM数据节点序列和给定阈值序列；将DEM数据节点序列中的第一个节点作为当前节点，将给定阈值序列中的第一个给定阈值作为当前阈值；
步骤5.2、从DEM数据中每个节点在多个整数尺度因子下的小波系数的平方和中选取当前节点在多个整数尺度因子下的小波系数的平方和；
当当前节点在多个整数尺度因子下得到的小波系数的平方和大于等于当前阈值时，将当前节点作为当前阈值对应的疑似滑坡区域节点；将当前节点的下一个节点作为当前节点，重复执行与步骤5.2相同的操作，直至判断完最后一个节点结束，得到当前阈值对应的疑似滑坡区域的多个节点；根据得到的当前阈值的疑似滑坡区域的多个节点构成得到当前阈值对应的疑似滑坡区域；
当当前节点在多个整数尺度因子下得到的小波系数的平方和小于当前阈值时，表示当前节点属于非滑坡区域；将当前节点的下一个节点作为当前节点，重复执行与步骤5.2相同的操作，直至判断完DEM数据中所有节点操作完成时操作结束；
步骤5.3、将当前阈值的下一个给定阈值作为当前阈值，重复执行与步骤5.2相同的操作，直至判断完最后一个给定阈值结束，得到多个给定阈值对应的疑似滑坡区域；
步骤5.4、采用误差分析方法对得到的多个给定阈值对应的疑似滑坡区域进行处理，得到优选阈值；从多个给定阈值对应的疑似滑坡区域中选取优选阈值对应的疑似滑坡区域，将优选阈值对应的疑似滑坡区域与疑似滑坡区进行交运算，得到潜在滑坡区域。
5.如权利要求3所述的融合DEM、光学遥感和形变信息的潜在滑坡识别方法，其特征在于，所述步骤4.1中，采用公式(1)计算得到多个尺度因子下已知滑坡区的DEM数据的小波系数，采用公式(2)计算得到多个尺度因子下已知非滑坡区的DEM数据的小波系数；
其中，(x,y)表示已知滑坡区域的DEM数据的节点，x和y均大于等于0；zfld(x,y)为已知滑坡区域的DEM数据中节点为(x,y)的高程，zfld(x,y)为实数；Cfld(si,a,b)为zfld(x,y)经过小波变换得到的小波系数，Cfld(si,a,b)∈R；si表示第i个尺度因子，si>0；为已知滑坡区的DEM数据中节点为(x,y)在第i个尺度因子下的小波函数值，a、b均表示小波变换的位移，a和b均为实数；
其中，(m,n)表示已知非滑坡区的DEM数据中的节点，m和n均大于等于0；zunfld(m,n)为已知非滑坡区的DEM数据中节点为(m,n)的高程，zunfld(m,n)为实数；Cunfld(si,a,b)为zunfld(m,n)经过小波变换得到的小波系数，Cunfld(si,a,b)∈R；si表示第i个尺度因子，si>0；为已知非滑坡区的DEM数据中节点为(m,n)在第i个尺度因子下的小波函数值，a、b均表示小波变换的位移，a和b均为实数。
6.如权利要求3所述的融合DEM、光学遥感和形变信息的潜在滑坡识别方法，其特征在于，所述步骤4.2的具体操作如下：
通过多个尺度因子下已知滑坡区的DEM数据的小波系数，采用公式(3)计算得到多个尺度因子下已知滑坡区的DEM数据的小波系数的方差；
其中，VCWT_fld(si)表示第i个尺度因子时已知滑坡区域的DEM数据的小波系数的方差，VCWT_fld(si)≥0；Na,Nb表示DEM数据中的节点个数，Na,Nb分别表示DEM数据中每一行和每一列的节点个数，Na,Nb均大于0；
通过多个尺度因子下已知非滑坡区的DEM数据的小波系数，采用公式(4)计算得到多个尺度因子下已知非滑坡区的DEM数据的小波系数的方差；
其中，VCWT_unfld(si)表示第i个尺度因子时已知非滑坡区域的DEM数据的小波系数的方差，VCWT_unfld(si)≥0；Na,Nb表示DEM数据中的节点个数，Na,Nb分别表示DEM数据中每一行和每一列的节点个数，Na,Nb均大于0。
7.如权利要求3所述的融合DEM、光学遥感和形变信息的潜在滑坡识别方法，其特征在于，所述步骤4.3中，通过多个尺度因子下已知非滑坡区的DEM数据的小波系数方差，采用公式(5)对多个尺度因子下已知滑坡区DEM数据的小波系数方差进行归一化处理，得到多个尺度因子下已知滑坡区的DEM数据归一化后的小波系数方差；
VCWT_norm(si)＝VCWT_fld(si)/VCWT_unfld(si) (5)
其中，VCWT_norm(si)表示第i个尺度因子的已知滑坡区的DEM数据的归一化小波系数方差，VCWT_fld(si)表示第i个尺度因子的已知滑坡区的DEM数据的小波系数的方差，VCWT_unfld(si)表示第i个尺度因子的已知非滑坡区的DEM数据的小波系数的方差，VCWT_norm(si),VCWT_fld(si)，VCWT_unfld(si)均大于等于0。
8.如权利要求3所述的融合DEM、光学遥感和形变信息的潜在滑坡识别方法，其特征在于，所述步骤4.4中，采用特征尺度因子中的多个整数尺度因子对DEM数据进行小波变换，具体采用公式(6)计算得到多个整数尺度因子下的DEM数据的小波系数；
其中，(xj,yj)表示DEM数据中第j个节点，且j>0，xj与yj均为实数；z(xj,yj)为DEM数据中第j个节点的高程，z(xj,yj)为实数；C(s'i,a,b)为z(xj,yj)经过小波变换得到的小波系数，C(s'i,a,b)∈R；s'i表示特征尺度因子中第i个整数尺度因子，s'i>0，i>0；为DEM数据中第j个节点在第i个整数尺度因子下的小波函数值，且a,b均表示小波变换的位移，且a,b均为实数。
9.如权利要求1所述的融合DEM、光学遥感和形变信息的潜在滑坡识别方法，其特征在于，所述步骤1的具体操作如下：
步骤1.1、对研究区域内的SAR卫星采集得到的SAR影像，采用小基线集技术对SAR影像进行处理，得到研究区域的的形变速率图；
步骤1.2、通过遥感卫星获取覆盖研究区域的光学遥感影像数据，所述光学遥感影像数据包括分辨率为0.5m的全色波段数据和分辨率为2m的4波段多光谱数据，所述4波段包括红、绿、蓝和近红外波段；
步骤1.3、利用无人机获取研究区域的分辨率为3m的DEM数据，利用光学遥感影像数据自带的RPC文件和RPC模型得到完整的RPC模型；利用DEM数据和完整的RPC模型分别对光学遥感影像数据进行正射校正，得到正射光学遥感影像，所述RPC指有理多项式系数；
步骤1.4、采用NNDiffuse Pan Sharpening算法对正射光学遥感影像进行融合，得到融合的光学遥感影像图。
10.如权利要求1所述的融合DEM、光学遥感和形变信息的潜在滑坡识别方法，其特征在于，所述步骤2的具体操作如下：
步骤2.1、采用地图投影方法对研究区域的形变速率图定义投影，得到具有投影信息的形变速率图，使具有投影信息的形变速率图与融合的光学遥感影像图的投影坐标系统一致；采用三次卷积的方法对具有投影信息的形变速率图进行重采样，得到重采样后的形变速率图，使重采样后的形变速率图与融合的光学遥感影像图的分辨率保持一致；
步骤2.2、采用地图投影方法对所述DEM数据定义投影，得到具有投影信息的DEM数据，使具有投影信息的DEM数据与融合的光学遥感影像图的投影坐标系统保持一致；
步骤2.3、将具有投影信息的DEM数据作为一个波段数据，采用多尺度分割方法，对具有投影信息的DEM数据和融合的光学遥感影像图中的4波段多光谱数据整体进行分割，得到多个对象，所述对象是指具有同质性的像元集合。
CN201910237520.7A 2019-03-27 2019-03-27 一种融合dem、光学遥感和形变信息的潜在滑坡识别方法 Active CN110120046B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201910237520.7A CN110120046B (zh) 2019-03-27 2019-03-27 一种融合dem、光学遥感和形变信息的潜在滑坡识别方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201910237520.7A CN110120046B (zh) 2019-03-27 2019-03-27 一种融合dem、光学遥感和形变信息的潜在滑坡识别方法

Publications (2)

Publication Number Publication Date
CN110120046A true CN110120046A (zh) 2019-08-13
CN110120046B CN110120046B (zh) 2022-09-16

Family

ID=67520628

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201910237520.7A Active CN110120046B (zh) 2019-03-27 2019-03-27 一种融合dem、光学遥感和形变信息的潜在滑坡识别方法

Country Status (1)

Country Link
CN (1) CN110120046B (zh)

Cited By (10)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN110927723A (zh) * 2019-11-11 2020-03-27 中国地质环境监测院 毫米波雷达泥石流智能监测预警系统与方法
CN111047616A (zh) * 2019-12-10 2020-04-21 中国人民解放军陆军勤务学院 一种遥感图像滑坡目标约束性主动轮廓特征提取方法
CN111077525A (zh) * 2019-12-20 2020-04-28 长安大学 融合sar与光学偏移量技术的地表维形变计算方法及系统
CN111308468A (zh) * 2019-11-27 2020-06-19 北京东方至远科技股份有限公司 一种基于In SAR技术的形变风险区域自动识别的方法
CN111398958A (zh) * 2020-04-03 2020-07-10 兰州大学 一种确定地面沉降与黄土填挖方区建筑高度相关性的方法
CN112036460A (zh) * 2020-08-24 2020-12-04 河海大学 一种识别量化控制泉流量潜在因素的方法
CN113192086A (zh) * 2021-05-11 2021-07-30 中国自然资源航空物探遥感中心 地质灾害隐患变形强度分布图的生成方法和存储介质
CN113219912A (zh) * 2021-03-31 2021-08-06 成都飞机工业（集团）有限责任公司 一种基于Multi-Agent的数控加工柔性制造系统加工过程预警方法
CN113887515A (zh) * 2021-10-28 2022-01-04 中国自然资源航空物探遥感中心 一种基于卷积神经网络的遥感滑坡识别方法及系统
CN114820552A (zh) * 2022-05-11 2022-07-29 中国地质环境监测院(自然资源部地质灾害技术指导中心) 一种利用光学卫星立体影像检测滑坡位移场的方法

Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20090161944A1 (en) * 2007-12-21 2009-06-25 Industrial Technology Research Institute Target detecting, editing and rebuilding method and system by 3d image
EP2392943A1 (en) * 2010-06-03 2011-12-07 Ellegi S.r.l. Synthetic-aperture radar system and operating method for monitoring ground and structure displacements suitable for emergency conditions
CN105989322A (zh) * 2015-01-27 2016-10-05 同济大学 一种基于高分辨率遥感影像的多指标融合滑坡检测方法
CN107132539A (zh) * 2017-05-03 2017-09-05 中国地质科学院探矿工艺研究所 一种基于小基线集的时间序列InSAR的滑坡早期识别方法

Patent Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20090161944A1 (en) * 2007-12-21 2009-06-25 Industrial Technology Research Institute Target detecting, editing and rebuilding method and system by 3d image
EP2392943A1 (en) * 2010-06-03 2011-12-07 Ellegi S.r.l. Synthetic-aperture radar system and operating method for monitoring ground and structure displacements suitable for emergency conditions
CN105989322A (zh) * 2015-01-27 2016-10-05 同济大学 一种基于高分辨率遥感影像的多指标融合滑坡检测方法
CN107132539A (zh) * 2017-05-03 2017-09-05 中国地质科学院探矿工艺研究所 一种基于小基线集的时间序列InSAR的滑坡早期识别方法

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
丁辉等: "基于多特征面向对象区域滑坡现象识别", 《遥感技术与应用》 *
陈天博等: "无人机遥感数据处理与滑坡信息提取", 《地球信息科学学报》 *

Cited By (15)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN110927723A (zh) * 2019-11-11 2020-03-27 中国地质环境监测院 毫米波雷达泥石流智能监测预警系统与方法
CN110927723B (zh) * 2019-11-11 2021-01-01 中国地质环境监测院 毫米波雷达泥石流智能监测预警系统与方法
CN111308468A (zh) * 2019-11-27 2020-06-19 北京东方至远科技股份有限公司 一种基于In SAR技术的形变风险区域自动识别的方法
CN111047616A (zh) * 2019-12-10 2020-04-21 中国人民解放军陆军勤务学院 一种遥感图像滑坡目标约束性主动轮廓特征提取方法
CN111077525A (zh) * 2019-12-20 2020-04-28 长安大学 融合sar与光学偏移量技术的地表维形变计算方法及系统
CN111077525B (zh) * 2019-12-20 2022-12-27 长安大学 融合sar与光学偏移量技术的地表三维形变计算方法及系统
CN111398958A (zh) * 2020-04-03 2020-07-10 兰州大学 一种确定地面沉降与黄土填挖方区建筑高度相关性的方法
CN112036460A (zh) * 2020-08-24 2020-12-04 河海大学 一种识别量化控制泉流量潜在因素的方法
CN112036460B (zh) * 2020-08-24 2022-08-30 河海大学 一种识别量化控制泉流量潜在因素的方法
CN113219912A (zh) * 2021-03-31 2021-08-06 成都飞机工业（集团）有限责任公司 一种基于Multi-Agent的数控加工柔性制造系统加工过程预警方法
CN113219912B (zh) * 2021-03-31 2022-03-15 成都飞机工业（集团）有限责任公司 一种基于Multi-Agent的数控加工柔性制造系统加工过程预警方法
CN113192086A (zh) * 2021-05-11 2021-07-30 中国自然资源航空物探遥感中心 地质灾害隐患变形强度分布图的生成方法和存储介质
CN113192086B (zh) * 2021-05-11 2022-01-28 中国自然资源航空物探遥感中心 地质灾害隐患变形强度分布图的生成方法和存储介质
CN113887515A (zh) * 2021-10-28 2022-01-04 中国自然资源航空物探遥感中心 一种基于卷积神经网络的遥感滑坡识别方法及系统
CN114820552A (zh) * 2022-05-11 2022-07-29 中国地质环境监测院(自然资源部地质灾害技术指导中心) 一种利用光学卫星立体影像检测滑坡位移场的方法

Also Published As

Publication number Publication date
CN110120046B (zh) 2022-09-16

Similar Documents

Publication Publication Date Title
CN110120046B (zh) 一种融合dem、光学遥感和形变信息的潜在滑坡识别方法
CN110929607B (zh) 一种城市建筑物施工进度的遥感识别方法和系统
Ke et al. A review of methods for automatic individual tree-crown detection and delineation from passive remote sensing
CN107451982B (zh) 一种基于无人机影像的高郁闭度林分树冠面积获取方法
CN108492260B (zh) 基于张量投票耦合霍夫变换的地质线性体提取方法
CN103729848B (zh) 基于光谱显著性的高光谱遥感图像小目标检测方法
CN111027547A (zh) 一种针对二维图像中的多尺度多形态目标的自动检测方法
CN103218787B (zh) 多源异构遥感影像控制点自动采集方法
CN107977661B (zh) 基于fcn与低秩稀疏分解的感兴趣区域检测方法
Gueguen et al. Toward a generalizable image representation for large-scale change detection: Application to generic damage analysis
CN113610070A (zh) 一种基于多源数据融合的滑坡灾害识别方法
CN105069818A (zh) 一种基于图像分析的皮肤毛孔识别方法
Tarsha Kurdi et al. Automatic filtering and 2D modeling of airborne laser scanning building point cloud
Touati et al. A reliable mixed-norm-based multiresolution change detector in heterogeneous remote sensing images
CN104951765B (zh) 基于形状先验信息和视觉对比度的遥感图像目标分割方法
CN117423010A (zh) 一种基于遥感数据的河湖划界识别监测方法
CN109544608B (zh) 一种无人机图像采集特征配准方法
CN105678720A (zh) 一种全景拼接判断图像匹配方法及装置
CN115063698A (zh) 斜坡地表变形裂缝自动识别与信息提取方法及系统
CN115019163A (zh) 基于多源大数据的城市要素识别方法
CN109461176A (zh) 高光谱图像的光谱配准方法
WO2024222610A1 (zh) 基于深度卷积网络的遥感影像变化检测方法
Vats et al. Terrain-Informed Self-Supervised Learning: Enhancing Building Footprint Extraction from LiDAR Data with Limited Annotations
Huang et al. Morphological building index (MBI) and its applications to urban areas
Mahphood et al. Virtual first and last pulse method for building detection from dense LiDAR point clouds

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant