CN117523404B - 一种堆土场动态变化监测方法、装置、终端及存储介质 - Google Patents

一种堆土场动态变化监测方法、装置、终端及存储介质 Download PDF

Info

Publication number
CN117523404B
CN117523404B CN202410022651.4A CN202410022651A CN117523404B CN 117523404 B CN117523404 B CN 117523404B CN 202410022651 A CN202410022651 A CN 202410022651A CN 117523404 B CN117523404 B CN 117523404B
Authority
CN
China
Prior art keywords
image
time sequence
storage yard
time
phase
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
Application number
CN202410022651.4A
Other languages
English (en)
Other versions
CN117523404A (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.)
Shenzhen University
Original Assignee
Shenzhen 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 Shenzhen University filed Critical Shenzhen University
Priority to CN202410022651.4A priority Critical patent/CN117523404B/zh
Publication of CN117523404A publication Critical patent/CN117523404A/zh
Application granted granted Critical
Publication of CN117523404B publication Critical patent/CN117523404B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • 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/188—Vegetation
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01C—MEASURING DISTANCES, LEVELS OR BEARINGS; SURVEYING; NAVIGATION; GYROSCOPIC INSTRUMENTS; PHOTOGRAMMETRY OR VIDEOGRAMMETRY
    • G01C11/00—Photogrammetry or videogrammetry, e.g. stereogrammetry; Photographic surveying
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01C—MEASURING DISTANCES, LEVELS OR BEARINGS; SURVEYING; NAVIGATION; GYROSCOPIC INSTRUMENTS; PHOTOGRAMMETRY OR VIDEOGRAMMETRY
    • G01C5/00—Measuring height; Measuring distances transverse to line of sight; Levelling between separated points; Surveyors' levels
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F16/00—Information retrieval; Database structures therefor; File system structures therefor
    • G06F16/20—Information retrieval; Database structures therefor; File system structures therefor of structured data, e.g. relational data
    • G06F16/29—Geographical information databases
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T17/00—Three dimensional [3D] modelling, e.g. data description of 3D objects
    • G06T17/05—Geographic models
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2200/00—Indexing scheme for image data processing or generation, in general
    • G06T2200/04—Indexing scheme for image data processing or generation, in general involving 3D image data
    • Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
    • Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02A—TECHNOLOGIES FOR ADAPTATION TO CLIMATE CHANGE
    • Y02A90/00—Technologies having an indirect contribution to adaptation to climate change
    • Y02A90/10—Information and communication technologies [ICT] supporting adaptation to climate change, e.g. for weather forecasting or climate simulation

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Remote Sensing (AREA)
  • Theoretical Computer Science (AREA)
  • Multimedia (AREA)
  • Databases & Information Systems (AREA)
  • Geometry (AREA)
  • Software Systems (AREA)
  • Radar, Positioning & Navigation (AREA)
  • General Health & Medical Sciences (AREA)
  • Health & Medical Sciences (AREA)
  • Data Mining & Analysis (AREA)
  • General Engineering & Computer Science (AREA)
  • Computer Graphics (AREA)
  • Image Processing (AREA)

Abstract

本发明提供一种堆土场动态变化监测方法、装置、终端及存储介质，所述方法包括：获取第一光学遥感卫星采集的时序多光谱遥感影像，得到时序植被状况数据和时序土地类型并确定潜在堆土场；获取无人机采集的现状高清影像集，得到现状数字高程模型和现状数字正射影像，获取历史多时相数字正射影像，获取时序合成孔径雷达影像，得到历史多时相数字高程模型、地表形变速率图和时序沉降曲线；生成变化趋势信息；确定堆土场，生成堆土场的动态变化监测结果。本发明通过结合多光谱遥感影像、无人机拍摄的现状高清影像集、历史多时相数字正射影像和时序合成孔径雷达影像对堆土场进行动态变化监测，可有效提高堆土场监测结果的精度。

Description

一种堆土场动态变化监测方法、装置、终端及存储介质
技术领域
本发明涉及遥感技术领域，尤其涉及的是一种堆土场动态变化监测方法、装置、终端及存储介质。
背景技术
随着城市化和工业化进程的加速，基础设施建设迅猛发展，同时也不断产生大量废弃土壤和建筑材料。有效处理和安全存放这些废弃物已经成为环境管理中的一个重要问题。
传统的监测方法，如现场调查和样本分析，均较为耗时耗力。随着卫星遥感技术的不断进步，利用卫星影像对堆土场进行监测已成为一种新的技术手段。卫星遥感的优势在于其能够提供持续的、大范围的监测能力，对于难以接近或面积广阔的地区尤为有用。然而，尽管现有的卫星遥感技术可以覆盖广阔的地域，提供地表的宏观视图，但它们的空间分辨率较低，这导致一些细节信息难以被有效捕捉，使得堆土场监测结果的精度较低。
因此，现有技术存在缺陷，有待改进与发展。
发明内容
本发明要解决的技术问题在于，针对现有技术的上述缺陷，提供一种堆土场动态变化监测方法、装置、终端及存储介质，旨在解决现有技术中堆土场监测结果的精度较低问题。
本发明解决技术问题所采用的技术方案如下：
第一方面，本发明实施例提供一种堆土场动态变化监测方法，所述方法包括：
获取第一光学遥感卫星采集的监测区域的时序多光谱遥感影像，以确定所述监测区域的时序植被状况数据和时序土地类型数据，并基于所述时序植被状况数据或所述时序土地类型数据确定所述监测区域内的潜在堆土场；
获取无人机采集的所述潜在堆土场的现状高清影像集以得到现状数字高程模型和现状数字正射影像、获取第二光学遥感卫星提供的所述潜在堆土场的历史多时相数字正射影像、以及获取星载平台提供的所述潜在堆土场的时序合成孔径雷达影像以得到历史多时相数字高程模型、地表形变速率图和时序沉降曲线；
基于所述时序植被状况数据、所述时序土地类型数据、所述现状数字高程模型、所述现状数字正射影像、所述历史多时相数字正射影像、所述历史多时相数字高程模型、所述地表形变速率图和所述时序沉降曲线得到所述潜在堆土场的变化趋势信息，将所述变化趋势信息和预设条件进行比对；
若所述变化趋势信息符合预设条件，将所述潜在堆土场确定为堆土场，利用GIS平台对所述时序植被状况数据、所述时序土地类型数据、所述现状数字高程模型、所述现状数字正射影像、所述历史多时相数字正射影像、所述历史多时相数字高程模型、所述地表形变速率图和所述时序沉降曲线进行融合分析，得到所述堆土场的动态变化监测结果。
在一种实施方式中，所述获取第一光学遥感卫星采集的监测区域的时序多光谱遥感影像，以确定所述监测区域的时序植被状况数据和时序土地类型数据，包括：
获取第一光学遥感卫星采集的监测区域的时序多光谱遥感影像，利用预设的坐标投影转换公式，将所述时序多光谱遥感影像从原始投影坐标系转换为第一目标坐标系，生成多个基于时间序列的第一多光谱遥感影像；
提取并组合每个所述第一多光谱遥感影像中预设波段的多个单一波段图像，得到每个所述第一多光谱遥感影像对应的多波段图像；
提取并组合每个所述第一多光谱遥感影像中的红光谱波段图像、绿光谱波段图像和蓝光谱波段图像，得到每个所述第一多光谱遥感影像对应的真彩色图像；
对每个所述多波段图像和每个所述真彩色图像均进行辐射定标、大气校正和边界区域裁剪，得到多个第一波段图像和多个第一真彩色图像；
对全部所述第一波段图像进行处理，得到所述监测区域的时序植被状况数据；
利用最大似然估计监督学习算法对全部所述第一真彩色图像进行处理，得到所述监测区域的时序土地类型数据。
在一种实施方式中，所述基于所述时序植被状况数据或所述时序土地类型数据确定所述监测区域内的潜在堆土场，包括：
基于所述时序植被状况数据或所述时序土地类型数据，确定所述监测区域内的突变区域；
将所述突变区域作为潜在堆土场。
在一种实施方式中，所述获取无人机采集的所述潜在堆土场的现状高清影像集以得到现状数字高程模型和现状数字正射影像，包括：
获取无人机采集的所述潜在堆土场的现状高清影像集，所述现状高清影像集中包括多个现状高清影像；
对每个所述现状高清影像均进行匀光匀色处理和影像畸变纠正，得到多个第一现状高清影像；
利用空中三角测量技术对每个所述第一现状高清影像进行特征点提取和多视影像联合平差，得到多个第二现状高清影像；
基于每个所述第二现状高清影像，利用密集匹配技术，得到所述潜在堆土场对应的密集点云；
基于所述密集点云，得到现状数字高程模型和现状数字正射影像。
在一种实施方式中，所述获取星载平台提供的所述潜在堆土场的时序合成孔径雷达影像以得到历史多时相数字高程模型、地表形变速率图和时序沉降曲线，包括：
获取星载平台提供的所述潜在堆土场的时序合成孔径雷达影像，所述时序合成孔径雷达影像由每个月所述星载平台拍摄的合成孔径雷达影像对组成，所述合成孔径雷达影像对是两个拍摄时间间隔为12天的合成孔径雷达影像；
按照预设规则对多个所述合成孔径雷达影像对进行选择，得到多个目标合成孔径雷达影像对，确定每个所述目标合成孔径雷达影像对的主影像和从影像，将每个所述从影像与对应的所述主影像进行配准，得到每个所述从影像对应的第一配准影像；
对每对所述主影像和所述第一配准影像进行干涉相位计算和相干性估计，生成多个干涉图，以及每个干涉图对应的相干图，对每个所述干涉图去除平地相位误差，生成多个第一干涉图，利用Goldstein滤波平滑每个所述第一干涉图，得到多个第二干涉图，基于每个所述相干图，利用最小费用流算法对应解缠每个所述第二干涉图，得到多个解缠后的相位图；
将全部所述相位图中的相位信息向高程转化，得到全部所述相位图中每个像素的高程信息，将全部所述相位图中每个像素的所述高程信息从雷达影像坐标系转换到第二目标坐标系，得到全部所述相位图中每个像素的第一高程信息；
采用双线性内插法对每个所述相位图中每个像素的所述第一高程信息进行重采样，得到历史多时相数字高程模型。
在一种实施方式中，所述获取星载平台提供的所述潜在堆土场的时序合成孔径雷达影像之后，还包括：
对所述潜在堆土场的所述时序合成孔径雷达影像进行前置滤波处理，得到第一时序合成孔径雷达影像，选取所述第一时序合成孔径雷达影像中成像效果符合预设标准的影像作为第一主影像，在所述时序合成孔径雷达影像中去除所述第一主影像，得到第一从影像集，所述第一从影像集中包括多个第一从影像；
将每个所述第一从影像与所述第一主影像进行配准，得到每个所述第一从影像对应的第二配准影像；
采用小基线原则将所述第一主影像和每个所述第二配准影像进行组合，生成多个第三干涉图，以及每个所述第三干涉图对应的第一相干图，所述第一相干图中包含对应的所述第三干涉图中不同点的相干性信息，对每个所述第三干涉图依次进行去平地相位误差、Goldstein滤波以及去地形相位后，得到多个差分干涉图，对每个所述差分干涉图进行多视处理，得到多个第一差分干涉图；
确定每个所述第一相干图中的高相干点，利用每个所述第一相干图中的全部所述高相干点，对每个所述第一相干图对应的每个所述第一差分干涉图进行相位解缠，得到多个第二差分干涉图，所述第二差分干涉图中包含每个所述高相干点对应的高相干区域内全部像素的初始解缠相位；
对每个所述第二差分干涉图估算对流层延迟相位，得到每个所述第二差分干涉图对应的区域对流层改正图，所述区域对流层改正图用于表征所述第二差分干涉图中每个像素的对流层延迟相位；
利用每个所述区域对流层改正图对每个所述区域对流层改正图对应的每个所述第二差分干涉图进行对流层延迟相位消除，得到多个第三差分干涉图，所述第三差分干涉图中包含每个所述高相干点对应的高相干区域内全部像素的第二解缠相位，所述第二解缠相位为对所述高相干区域内的每个像素的所述初始解缠相位消除对应的所述对流层延迟相位所得；
获取每个所述第三差分干涉图对应的所述第一主影像的成像时间和所述第一从影像的成像时间，基于每个所述第三差分干涉图对应的所述第一主影像的成像时间和所述第一从影像的成像时间，得到多个时间差矩阵；
基于每个所述第三差分干涉图中每个所述高相干点对应的高相干区域内全部像素的所述第二解缠相位和每个所述第三差分干涉图对应的所述时间差矩阵，得到时序合成孔径雷达影像在视线方向的形变速率，基于所述形变速率得到地表形变速率图和时序沉降曲线。
在一种实施方式中，所述对每个所述第二差分干涉图进行对流层延迟校正，得到每个所述第二差分干涉图对应的区域对流层改正图，包括：
获取所述潜在堆土场及所述潜在堆土场预设范围内的周边的气象信息，基于所述气象信息得到多个气象数据；
获取所述潜在堆土场中预设地点的地表高程和映射系数，所述映射系数用于将天顶对流层延迟投影至卫星视线方向；
将各个所述气象数据、所述地表高程和所述映射系数输入对流层延迟校正模型，得到所述潜在堆土场中每个所述预设地点的对流层延迟相位；
利用插值方法对所述潜在堆土场中全部所述预设地点的所述对流层延迟相位进行插值处理，得到每个所述第二差分干涉图对应的区域对流层改正图。
第二方面，本发明实施例还提供一种堆土场动态变化监测装置，包括：
潜在堆土场确定模块，用于获取第一光学遥感卫星采集的监测区域的时序多光谱遥感影像，以确定所述监测区域的时序植被状况数据和时序土地类型数据，并基于所述时序植被状况数据或所述时序土地类型数据确定所述监测区域内的潜在堆土场；
遥感成果生成模块，用于获取无人机采集的所述潜在堆土场的现状高清影像集以得到现状数字高程模型和现状数字正射影像、获取第二光学遥感卫星提供的所述潜在堆土场的历史多时相数字正射影像、以及获取星载平台提供的所述潜在堆土场的时序合成孔径雷达影像以得到历史多时相数字高程模型、地表形变速率图和时序沉降曲线；
变化趋势信息生成模块，用于基于所述时序植被状况数据、所述时序土地类型数据、所述现状数字高程模型、所述现状数字正射影像、所述历史多时相数字正射影像、所述历史多时相数字高程模型、所述地表形变速率图和所述时序沉降曲线得到所述潜在堆土场的变化趋势信息，将所述变化趋势信息和预设条件进行比对；
动态变化监测结果生成模块，用于若所述变化趋势信息符合预设条件，将所述潜在堆土场确定为堆土场，利用GIS平台对所述时序植被状况数据、所述时序土地类型数据、所述现状数字高程模型、所述现状数字正射影像、所述历史多时相数字正射影像、所述历史多时相数字高程模型、所述地表形变速率图和所述时序沉降曲线进行融合分析，得到所述堆土场的动态变化监测结果。
第三方面，本发明实施例还提供一种终端，所述终端包括：存储器、处理器及存储在所述存储器上并可在所述处理器上运行的堆土场动态变化监测程序，所述堆土场动态变化监测程序被所述处理器执行时实现如上所述的堆土场动态变化监测方法的步骤。
第四方面，本发明实施例还提供一种计算机可读存储介质，所述计算机可读存储介质存储有堆土场动态变化监测程序，所述堆土场动态变化监测程序能够被执行以用于实现如上所述的堆土场动态变化监测方法的步骤。
本发明的有益效果：本发明实施例通过获取第一光学遥感卫星采集的时序多光谱遥感影像，得到时序植被状况数据和时序土地类型并确定潜在堆土场；获取无人机采集的现状高清影像集，得到现状数字高程模型和现状数字正射影像，获取历史多时相数字正射影像，获取时序合成孔径雷达影像，得到历史多时相数字高程模型、地表形变速率图和时序沉降曲线；生成变化趋势信息；确定堆土场，生成堆土场的动态变化监测结果。本发明通过结合多光谱遥感影像、无人机拍摄的现状高清影像集、历史多时相数字正射影像和时序合成孔径雷达影像对堆土场进行动态变化监测，可有效提高堆土场监测结果的精度。
附图说明
图1是本发明中堆土场动态变化监测方法较佳实施例的流程图。
图2是本发明中GIS平台进行融合分析的示意图。
图3是本发明中堆土场动态变化监测装置较佳实施例的功能原理框图。
图4是本发明的终端原理框图。
具体实施方式
为使本发明的目的、技术方案及优点更加清楚、明确，以下参照附图并举实施例对本发明进一步详细说明。应当理解，此处所描述的具体实施例仅用以解释本发明，并不用于限定本发明。
随着城市化和工业化进程的加速，基础设施建设迅猛发展，同时也不断产生大量废弃土壤和建筑材料。有效处理和安全存放这些废弃物已经成为环境管理中的一个重要问题。随着填埋场滑坡事故的多次发生，迫切需要动态掌握堆土场的状态，及早发现并制止违法堆土行为，以更好地保护自然环境和人员安全。
传统的监测方法，如现场调查和样本分析，均较为耗时耗力。随着卫星遥感技术的不断进步，利用卫星影像对堆土场进行监测已成为一种新的技术手段。卫星遥感的优势在于其能够提供持续的、大范围的监测能力，对于难以接近或面积广阔的地区尤为有用。然而，尽管现有的卫星遥感技术可以覆盖广阔的地域，提供地表的宏观视图，但它们的空间分辨率较低，这导致一些细节信息难以被有效捕捉，使得堆土场监测结果的精度较低。
针对现有技术的上述缺陷，本发明提供一种堆土场动态变化监测方法、装置、终端及存储介质，所述方法包括：获取第一光学遥感卫星采集的时序多光谱遥感影像，得到时序植被状况数据和时序土地类型并确定潜在堆土场；获取无人机采集的现状高清影像集，得到现状数字高程模型和现状数字正射影像，获取历史多时相数字正射影像，获取时序合成孔径雷达影像，得到历史多时相数字高程模型、地表形变速率图和时序沉降曲线；生成变化趋势信息；确定堆土场，生成堆土场的动态变化监测结果。本发明通过结合多光谱遥感影像、无人机拍摄的现状高清影像集、历史多时相数字正射影像和时序合成孔径雷达影像对堆土场进行动态变化监测，可有效提高堆土场监测结果的精度。
请参见图1，本发明实施例所述的堆土场动态变化监测方法包括如下步骤：
步骤S100、获取第一光学遥感卫星采集的监测区域的时序多光谱遥感影像，以确定所述监测区域的时序植被状况数据和时序土地类型数据，并基于所述时序植被状况数据或所述时序土地类型数据确定所述监测区域内的潜在堆土场。
具体地，多光谱影像是一种通过将成像技术和光谱测量技术相结合而获取的图像，它不仅包括二维空间信息，还包含随波长分布的光谱辐射信息。本发明中，所述第一光学遥感卫星为Landsat-8卫星。所述时序多光谱遥感影像为一系列基于时间序列的多光谱遥感影像。所述时序植被状况数据为一系列基于时间序列的植被状况数据，所述时序土地类型数据为一系列基于时间序列的土地类型数据。所述土地类型数据包括耕地、林地、草地、沙地、裸土地、建筑体、水域等。本发明利用时序多光谱遥感影像，能够确定监测区域内潜在的堆土场。除了能够识别正规的堆土场外，还能有效识别出隐藏在河道、山谷等隐蔽处的违规堆土场，为堆土场的监管提供了有效的指导。
在一种实现方式中，所述获取第一光学遥感卫星采集的监测区域的时序多光谱遥感影像，以确定所述监测区域的时序植被状况数据和时序土地类型数据，包括：
获取第一光学遥感卫星采集的监测区域的时序多光谱遥感影像，利用预设的坐标投影转换公式，将所述时序多光谱遥感影像从原始投影坐标系转换为第一目标坐标系，生成多个基于时间序列的第一多光谱遥感影像；
提取并组合每个所述第一多光谱遥感影像中预设波段的多个单一波段图像，得到每个所述第一多光谱遥感影像对应的多波段图像；
提取并组合每个所述第一多光谱遥感影像中的红光谱波段图像、绿光谱波段图像和蓝光谱波段图像，得到每个所述第一多光谱遥感影像对应的真彩色图像；
对每个所述多波段图像和每个所述真彩色图像均进行辐射定标、大气校正和边界区域裁剪，得到多个第一波段图像和多个第一真彩色图像；
对全部所述第一波段图像进行处理，得到所述监测区域的时序植被状况数据；
利用最大似然估计监督学习算法对全部所述第一真彩色图像进行处理，得到所述监测区域的时序土地类型数据。
具体地，当获取到时序多光谱遥感影像后，先进行坐标转换，生成多个第一多光谱遥感影像。然后提取每个第一多光谱遥感影像中B2-B7波段的波段图像，将B2-B7波段的波段图像进行组合，生成多个多波段图像。后续将基于多个多波段图像得到时序植被状况数据。提取每个第一多光谱遥感影像中的红光谱波段图像、绿光谱波段图像和蓝光谱波段图像，将红光谱波段图像、绿光谱波段图像和蓝光谱波段图像进行组合，生成多个真彩色图像，后续将基于多个多彩色图像得到时序土地类型数据。
当生成多波段图像和真彩色图像后，对全部多波段图像和全部真彩色图像进行辐射定标、大气校正和边界区域裁剪。
辐射定标过程为把多波段图像和真彩色图像中每个像素的地物灰度值转化为辐射亮度，其辐射定标公式表示为：
（1），
其中，是波段λ的辐射亮度、/>为像素灰度值、/>为增益参数、/>为偏移参数。增益参数和偏移参数均为已知参数，像素灰度值可从影像中得到。
然后，采用快速大气校正模型(QUick Atmospheric Correction, QUAC)进行大气校正，其公式表示为：
（2），
其中，为从影像中得到的端元光谱的反射率，/>为端元光谱的平均反射率，/>表示图像中的端元个数。利用快速大气校正模型对全部多波段图像和全部真彩色图像处理，进一步提升图像质量。
可以理解的，通过辐射定标、大气校正这两个步骤可逐步提高多波段图像和真彩色图像的质量。通过边界区域裁剪这个步骤可裁剪掉边界或无效区域，提高后续图像处理的效率。
当进行完上述的预处理过程后，需要对第一波段图像进行处理，得到监测区域的时序植被状况数据。时序植被状况数据包括植被指数（NDVI）和植被覆盖度(FVC)。利用植被指数计算公式计算植被指数，利用植被覆盖度公式计算植被覆盖度。
植被指数计算公式表示为：
（3），
其中，为植被指数，/>为B5波段，/>为B4波段，即红波段。植被指数（NDVI）数值越接近1表示植被越繁茂，植被指数（NDVI）数值越接近0则表示越缺乏植被。
对所有第一波段图像的全部像素进行计算得到区域植被指数分布图。
植被覆盖度公式表示为：
（4），
其中，为植被覆盖度，/>为植被指数，/>为区域植被指数分布图中累积频率为5%的区域对应的植被指数（NDVI），视为裸土区的植被指数（NDVI），/>为区域植被指数分布图中累积频率为95%的区域对应的植被指数（NDVI），代表植被区像素的植被指数（NDVI）。植被覆盖度（FVC）越接近1表示地表被植被覆盖的比例越高。植被指数（NDVI）和植被覆盖度（FVC）可以相互补充，提供更准确和全面的植被信息。
在计算监测区域的时序土地类型数据时，先获取包含类别标签的训练样本，训练样本为经过感兴趣区域编辑（ROI）处理的真彩色图像中的像素光谱值。类别标签包括耕地、林地、草地、沙地、裸土地、建筑体、水域等。然后，使用最大化样本的似然函数来确定最优的分类决策边界。
最大化样本的似然函数的公式可表示为：
（5），
其中，L为似然函数，是训练样本，/>是概率密度函数，表示在给定参数/>的情况下，观测到数据点/>的概率。/>是模型参数。通过大量的训练样本来不断调整参数/>，使得似然函数最大化。当得到最大化似然函数后，获取全部第一真彩色图像中的全部像素光谱值，利用公式（5）对全部第一真彩色图像中的全部像素光谱值进行处理，得到监测区域的时序土地类型数据。
在一种实现方式中，所述基于所述时序植被状况数据或所述时序土地类型数据确定所述监测区域内的潜在堆土场，包括：
基于所述时序植被状况数据或所述时序土地类型数据，确定所述监测区域内的突变区域；
将所述突变区域作为潜在堆土场。
具体的，堆填土方会导致裸土大面积出现，同时伴随着短期内植被指数下降。本发明中，突变区域基于时序植被状况数据确定或基于时序土地类型数据确定。
在一种实现方式中，将植被指数（NDVI）下降区间超过第一预设区间的区域作为突变区域。
在一种实现方式中，将植被覆盖度（FVC）下降区间超过第二预设区间的区域作为突变区域。
在一种实现方式中，将土地类型从耕地、林地、水域向沙地或裸土地变化的区域作为突变区域。
请参见图1，本发明实施例所述的堆土场动态变化监测方法还包括如下步骤：
步骤S200、获取无人机采集的所述潜在堆土场的现状高清影像集以得到现状数字高程模型和现状数字正射影像、获取第二光学遥感卫星提供的所述潜在堆土场的历史多时相数字正射影像、以及获取星载平台提供的所述潜在堆土场的时序合成孔径雷达影像以得到历史多时相数字高程模型、地表形变速率图和时序沉降曲线。
具体地，无人机被广泛应用于实景三维建模，具有高机动、低成本和灵活部署的优势，在小区域的实景影像重建以及数字正射影像（DOM）、数字高程模型（DEM）制作上，相比于卫星遥感、点云等测量建模方式具有更高的精度、分辨率和更真实的纹理映射等优势。然而，由于无人机的续航距离较短，故暂不支持区域级的航拍建模，因此，基于无人机建模的观测仅适用于对已知堆土场进行的精细观察。对于尚未开展定期航拍的堆土场（例如尚未发现的违填堆土场、已封场停用的堆土场等），由于缺乏历史时期的无人机航拍影像记录，无法对堆土场过往的堆填变化过程进行分析。本发明通过无人机来采集潜在堆土场的现状高清影像集，并基于现状高清影像集生成现状数字高程模型（现状DEM）和现状数字正射影像（现状DOM），可有效提高后续动态变化监测结果的精度。
在一种实现方式中，所述获取无人机采集的所述潜在堆土场的现状高清影像集以得到现状数字高程模型和现状数字正射影像，包括：
获取无人机采集的所述潜在堆土场的现状高清影像集，所述现状高清影像集中包括多个现状高清影像；
对每个所述现状高清影像均进行匀光匀色处理和影像畸变纠正，得到多个第一现状高清影像；
利用空中三角测量技术对每个所述第一现状高清影像进行特征点提取和多视影像联合平差，得到多个第二现状高清影像；
基于每个所述第二现状高清影像，利用密集匹配技术，得到所述潜在堆土场对应的密集点云；
基于所述密集点云，得到现状数字高程模型和现状数字正射影像。
具体地，本发明利用具有实时动态测量（RTK）定位功能的多旋翼无人机倾斜摄影得到潜在堆土场的现状影像集。可使用仿地飞行路径或五航线法拍摄潜在堆土场，在此不做限定。使用仿地飞行路径进行拍摄可减少采集时间，提高处理效率和三维重建质量，使用五航线法进行拍摄需保证航向重叠率不低于80%，旁向重叠率不低于70%。拍摄过程中记录拍摄路径，以保证后续拍摄的拍摄路径和之前的拍摄路径保持一致。
当获取到现状高清影像集后，需要对每个现状高清影像进行匀光匀色处理和影像畸变纠正，以提高影像质量，得到第一现状高清影像。
影像畸变纠正过程的公式可表示为：
（6），
其中，为像点坐标校正值，是未知量，/>为影像中像点的坐标，为像主点坐标，r是像点和像主点间距离，k1、k2、k3为径向畸变参数，p1、p2为切向畸变参数。其中，/>和/>可从影像中获取，k1、k2、k3和p1、p2可预先进行校准得到。
经过上述的影像预处理过程后，还需考虑无人机倾斜摄影带来的仿射变形误差，需要继续利用空中三角测量技术对第一现状高清影像进行处理，为后续的三维重建提供数据支撑。利用空中三角测量技术对每个所述第一现状高清影像进行特征点提取和多视影像联合平差，得到多个第二现状高清影像。此部分技术为现有技术，在此不再赘述。
当得到多个第二现状高清影像后，采用SGM半全局匹配算法快速提取并自动匹配全部第二现状高清影像中的同名点，通过多视前方交会恢复物方三维点坐标，生成潜在堆土场对应的密集点云。
基于密集点云，构建多稀疏度三角网格，得到数字表面模型（DSM），将数字表面模型作为现状数字高程模型（现状DEM），将多稀疏度三角网格封装为三维白膜模型，分析三维白膜模型中的任意三角形与相应范围影像的单一映射关系，获取所有三角形的纹理范围，完成纹理与三维白膜模型的配准和映射，得到潜在堆土场的三维模型，对三维模型进行正射纠正，得到现状数字正射影像（现状DOM）。
本发明中基于无人机的遥感技术对潜在堆土场进行三维重建，对比于控制点可实现厘米级别的平面和高程误差，在不使用外部控制点平差的情况下也可将高程中误差控制在1米以内，满足堆土场堆填过程大尺度变化监测需求。
在一种实现方式中，所述第二光学遥感卫星为SPOT-7卫星或WV-3卫星。SPOT-7卫星和WV-3卫星均可提供历史时期的数字正射影像（DOM）。
在一种实现方式中，所述获取星载平台提供的所述潜在堆土场的时序合成孔径雷达影像以得到历史多时相数字高程模型、地表形变速率图和时序沉降曲线，包括：
获取星载平台提供的所述潜在堆土场的时序合成孔径雷达影像，所述时序合成孔径雷达影像由每个月所述星载平台拍摄的合成孔径雷达影像对组成，所述合成孔径雷达影像对是两个拍摄时间间隔为12天的合成孔径雷达影像；
按照预设规则对多个所述合成孔径雷达影像对进行选择，得到多个目标合成孔径雷达影像对，确定每个所述目标合成孔径雷达影像对的主影像和从影像，将每个所述从影像与对应的所述主影像进行配准，得到每个所述从影像对应的第一配准影像；
对每对所述主影像和所述第一配准影像进行干涉相位计算和相干性估计，生成多个干涉图，以及每个干涉图对应的相干图，对每个所述干涉图去除平地相位误差，生成多个第一干涉图，利用Goldstein滤波平滑每个所述第一干涉图，得到多个第二干涉图，基于每个所述相干图，利用最小费用流算法对应解缠每个所述第二干涉图，得到多个解缠后的相位图；
将全部所述相位图中的相位信息向高程转化，得到全部所述相位图中每个像素的高程信息，将全部所述相位图中每个像素的所述高程信息从雷达影像坐标系转换到第二目标坐标系，得到全部所述相位图中每个像素的第一高程信息；
采用双线性内插法对每个所述相位图中每个像素的所述第一高程信息进行重采样，得到历史多时相数字高程模型。
具体地，所述星载平台为哨兵1A号（Sentinel-1A）。所述历史多时相数字高程模型反映了地表在不同时间段的高程变化情况。
本发明中，获取的时序合成孔径雷达影像由每个月所述星载平台拍摄的合成孔径雷达影像对组成，合成孔径雷达影像对是两个拍摄时间间隔为12天的合成孔径雷达影像。当获取了时序合成孔径雷达影像后，需要按照预设规则对多个所述合成孔径雷达影像对进行选择，得到多个目标合成孔径雷达影像对。
在一种实现方式中，预设规则为每一个季度选取一个合成孔径雷达影像对。可以理解的，堆土场的监测是个长期的过程，通过选取出目标合成孔径雷达影像对，可以减少数据的处理量，同时也可实现堆土场的监测。
当得到目标合成孔径雷达影像对后，先确定目标合成孔径雷达影像对中的主影像和从影像。
在一种实现方式中，将目标合成孔径雷达影像对中拍摄时间靠前的一景合成孔径雷达影像作为主影像，将目标合成孔径雷达影像对中拍摄时间靠后的一景合成孔径雷达影像作为从影像。
确定好主从影像后，对主影像和从影像进行前置滤波后，确定配准的配准控制点，然后采用二次多项式拟合主影像和从影像的同名点之间的几何变换模型来实现粗配准，在此基础上，对配准控制点的值进行内插，再利用最大频谱法实现亚像素精配准，通过两次配准，将主影像和从影像对齐到一个统一的参考坐标系中，配准的精度达到0.1个像素级。本发明通过两次配准，可以有效保证干涉图的生成质量和高程反演精度。经过两次配准后，将配准后生成的配准图像进行灰度重采样，得到第一配准影像。然后对每对主影像和第一配准影像进行干涉相位计算和相干性估计。
干涉相位计算的公式表示为：
（7），
其中，表示进行干涉计算的两个影像的干涉相位集合，里面包含了两个影像中每个像素的干涉相位，/>表示待计算的第一影像中每个像素的相位，/>表示待计算的第二影像中每个像素的相位，u表示第一影像和第二影像共轭相乘后的复数幅角值，/>和分别表示u的虚部和实部。在本发明中，可根据每对主影像和第一配准影像计算出u，进而得到每对主影像和第一配准影像的干涉相位集合。
相干性计算公式表示为：
（8），
其中，表示相干性，N为估计相干性时所取的像素空间，u1、u2表示两个零均值高斯复数随机信号，n为求和的索引。利用公式（8）可求出每对主影像和对应的第一配准影像的相干性。
对每对主影像和第一配准影像进行干涉相位计算和相干性估计，生成多个干涉图，以及每个干涉图对应的相干图，对每个所述干涉图去除平地相位误差，生成多个第一干涉图，利用Goldstein滤波平滑每个所述第一干涉图，得到多个第二干涉图，基于每个所述相干图，利用最小费用流算法对应解缠每个所述第二干涉图，得到多个解缠后的相位图。
当得到多个解缠后的相位图后，将全部相位图的相位信息向高程转化，此过程的公式可表示为：
（9），
其中，z表示高程，是解缠后的相位，/>、/>分别是水平基线和铅垂基线，/>是入射角，/>是波长，r1表示哨兵1A号（Sentinel-1A）卫星上的传感器到地面像素的斜距。
经过公式（9）的处理，可以把每个解缠后的相位图中的每个像素的相位信息向高程转化，得到每个解缠后的相位图中每个像素的高程信息，然后对每个像素的高程信息进行坐标转换，得到全部解缠后的相位图中每个像素的第一高程信息。
采用双线性内插法对每个解缠后的相位图中每个像素的所述第一高程信息进行重采样，得到历史多时相数字高程模型。
以往通过雷达干涉测量技术（InSAR）生成数字高程模型（DEM）时，多采用时间基线短（<1天）、空间基线较长（>2km）的两景影像进行干涉从而保证较高的数字高程模型（DEM）精度。但这对于合成孔径雷达（SAR）卫星平台的成像模式要求较高，使得使用成本较高，难以长期开展。无人机生成的数字高程模型(DEM)精度较高，但其航程和续航能力有限，操作成本较高。
本发明提出可以将时间上相隔12天的两景Sentinel-1A（属于哨兵1号卫星系列的一部分，A表示第一颗卫星）单视复数（SLC）影像（空间基线<500m）用于数字高程模型（DEM）制作，从而观察堆土场高程变化。
根据实验，本发明使用时间间隔12天的两景Sentinel影像制作的数字高程模型（DEM）与同月基于无人机生成的数字高程模型(DEM)间高程平均误差为3米，最大误差约为6米。而堆土场封场前后高度差可达近百米，因此本发明的基于雷达干涉测量技术（InSAR）生成的米级精度数字高程模型（DEM）足以反映不同时期堆土场高程变化趋势。此部分的处理可监测潜在堆土场的大尺度形变。
在一种实现方式中，除了潜在堆土场的大尺度形变的监测之外，由于堆填的弃土大多比较松散，在长期施工扰动下，容易发生沉降变形进而引起塌方，故还需要对地表的小尺度形变进行监测。故在获取星载平台提供的所述潜在堆土场的时序合成孔径雷达影像之后需要进一步处理，以生成地表形变速率图和时序沉降曲线。
首先，对所述潜在堆土场的时序合成孔径雷达影像进行前置滤波处理，得到第一时序合成孔径雷达影像，选取第一时序合成孔径雷达影像中成像效果符合预设标准的影像作为第一主影像，举例说明，可选择成像效果最好的合成孔径雷达影像作为第一主影像。在第一时序合成孔径雷达影像中去除第一主影像，得到第一从影像集，第一从影像集中包括多个第一从影像。
然后，将每个第一从影像与第一主影像进行配准，得到每个第一从影像对应的第二配准影像。
采用小基线原则将第一主影像和每个第二配准影像进行组合，生成多个第三干涉图，以及每个所述第三干涉图对应的第一相干图，所述第一相干图中包含对应的所述第三干涉图中不同点的相干性信息。
其中，使用tA时刻和tB时刻（tB>tA）基于时序合成孔径雷达影像生成的第j幅（j=1，2，……M）第三干涉图上的干涉相位可以表示为：
（10），
具体地，表示tB时刻的观测相位，可从tB时刻的时序合成孔径雷达影像中得到，表示tA时刻的观测相位，可从tA时刻的时序合成孔径雷达影像中得到。基于公式（10）可求出每个第三干涉图中每个像素的干涉相位，对每个第三干涉图依次进行去平地相位误差、Goldstein滤波以及去地形相位后，得到多个差分干涉图，对每个所述差分干涉图进行多视处理，得到多个第一差分干涉图。
此时相位仍是纠缠状态，需要选取高相干点进行解缠。确定每个第一相干图中的高相干点，基于每个所述第一相干图中的全部所述高相干点，利用最小费用流解缠算法（Minimum Cost Flow Algorithm）对每个第一相干图对应的每个第一差分干涉图进行相位解缠后，得到多个第二差分干涉图，第二差分干涉图中包含每个所述高相干点对应的高相干区域内全部像素的初始解缠相位。
然后，考虑对流层延迟对监测结果的影响，先对每个第二差分干涉图估算对流层延迟相位，得到每个第二差分干涉图对应的区域对流层改正图，区域对流层改正图用于表征第二差分干涉图中每个像素的对流层延迟相位，利用每个区域对流层改正图对与每个区域对流层改正图对应的每个第二差分干涉图进行对流层延迟相位进行消除，得到多个第三差分干涉图，所述第三差分干涉图中包含每个所述高相干点对应的高相干区域内的全部像素的第二解缠相位，第二解缠相位为对高相干区域内的每个像素的所述初始解缠相位消除对应的所述对流层延迟相位所得。
对于解缠后的第j幅（j=1，2，……M）第三差分干涉图上的高相干点对应的高相干区域的像素的第二解缠相位/>可表示为：
（11），
其中，表示像素的方位向坐标，/>表示像素的距离向坐标，/>为信号的中心波长，表示tB时刻相较于t0时刻卫星视线方向的累计形变量，/>表示tA时刻相较于t0时刻卫星视线方向的累计形变量。/>表示第二差分干涉图中的残余地形相位，为大气延迟相位，/>表示去相干噪声。可以理解的，/>、和/>均为误差相位，需要去除。
具体地，由于在生成第一差分干涉图之前已经进行了多视和Goldstein滤波处理，已经大幅削弱了噪声，可视为已经去除了相干噪声。
在不考虑大气延迟相位和残余地形相位对后续计算的影响的前提下，可将公式（11）表示成如下矩阵形式：
（12），
其中，B表示基于第一主影像的成像时间和第一从影像的成像时间之间的时间间隔得到的时间差矩阵，为M×N阶矩阵。时间差矩阵中的每个元素Bij表示以i为第一主影像（i<M）、j为第一从影像(j<N)，二者的成像时间间隔。表示第三差分干涉图中全部高相干点对应的高相干区域内全部像素的第二解缠相位，第二解缠相位是对高相干区域内的每个像素的所述初始解缠相位消除对应的所述对流层延迟相位所得。v表示时序合成孔径雷达影像在视线方向的形变速率，为未知量。
利用公式（12）进行奇异值分解（SVD）反演之后，通过时间上的高通滤波（High-pass filter，HP）和空间上的低通滤波（Low-pass filter，LP）估计出大气延迟相位和残余地形相位，并进行去除。此时的第三差分干涉图中全部高相干点对应的高相干区域内全部像素的第二解缠相位（）中不包含大气延迟相位和残余地形相位，然后利用公式（12）再次进行奇异值分解（SVD）反演，得到时序合成孔径雷达影像在视线方向的形变速率v。在各时段内对时序合成孔径雷达影像在视线方向的形变速率v进行积分变换，最终得到地表形变速率图和时序沉降曲线(时间序列形变量)。
在一种实现方式中，所述对每个所述第二差分干涉图进行对流层延迟校正，得到每个所述第二差分干涉图对应的区域对流层改正图，包括：
获取所述潜在堆土场及所述潜在堆土场预设范围内的周边的气象信息，基于所述气象信息得到多个气象数据；
获取所述潜在堆土场中预设地点的地表高程和映射系数，所述映射系数用于将天顶对流层延迟投影至卫星视线方向；
将各个所述气象数据、所述地表高程和所述映射系数输入对流层延迟校正模型，得到所述潜在堆土场中每个所述预设地点的对流层延迟相位；
利用插值方法对所述潜在堆土场中全部所述预设地点的所述对流层延迟相位进行插值处理，得到每个所述第二差分干涉图对应的区域对流层改正图。
具体地，气象信息可从当地气象站或其他公开平台获取，在此不做限制。常规的大气延迟相位消除操作通常假设大气具有在空间上低频变化、在时间上高频变化的特点，并通过时空滤波去除大气延迟相位。但是堆土场通常位于山坳等地形起伏区域，大气压力、温度和季节性水汽含量等变化产生的对流层延迟并未在时空滤波过程中有效去除。故本发明引入了考虑大气压力、温度以及水汽三大气象要素的对流层延迟校正模型来进行对流层延迟的消除。
对流层延迟校正模型可表示为：
（13），
具体的，表示某时刻某个预设地点的对流层延迟；z0表示预设地点的地表高程，是已知数据，/>表示对流层顶部高程，/>代表以为单位的地表的大气压力值，/>是对流层平均重力加速度，e是水汽压力，T是温度，/>是将天顶对流层延迟投影至卫星视觉方向（LOS）的映射系数，是通过标准映射函数基于气象站点的地理位置、观测时间以及卫星的方位和高度角进行计算所得，为已知数据，dz表示积分函数，N表示折射率。剩余参数均为常数，Rd=287.05J/（kg·K），k1=0.776K/Pa，k2=0.233K/Pa，k3=3.75×103K2/Pa。
通过公式（13）可以得到潜在堆土场中某个预设地点在某时刻的对流层延迟，在tA时刻，对流层延迟表示为，在tB时刻，对流层延迟表示为/>。所以在第二差分干涉图中，某个预设地点的对流层延迟改正值可表示为/>。当求出全部预设地点的对流层延迟改正值后，利用插值方法对潜在堆土场中全部预设地点的对流层延迟改正值进行插值处理，得到每个第二差分干涉图对应的区域对流层改正图。本发明通过将对流层延迟相位进行消除，可进一步提高堆土场监测结果的精度。
请参见图1，本发明实施例所述的堆土场动态变化监测方法还包括如下步骤：
步骤S300、基于所述时序植被状况数据、所述时序土地类型数据、所述现状数字高程模型、所述现状数字正射影像、所述历史多时相数字正射影像、所述历史多时相数字高程模型、所述地表形变速率图和所述时序沉降曲线得到所述潜在堆土场的变化趋势信息，将所述变化趋势信息和预设条件进行比对。
具体地，通过叠加时序植被状况数据、时序土地类型数据、现状数字高程模型、现状数字正射影像、历史多时相数字正射影像、历史多时相数字高程模型、地表形变速率图和时序沉降曲线这些数据进行综合研判，得到潜在堆土场的变化趋势信息，综合多类型数据可以得到较为准确的变化趋势信息，从而更合理地安排有限的巡查资源。
请参见图1，本发明实施例所述的堆土场动态变化监测方法还包括如下步骤：
步骤S400、若所述变化趋势信息符合预设条件，将所述潜在堆土场确定为堆土场，利用GIS平台对所述时序植被状况数据、所述时序土地类型数据、所述现状数字高程模型、所述现状数字正射影像、所述历史多时相数字正射影像、所述历史多时相数字高程模型、所述地表形变速率图和所述时序沉降曲线进行融合分析，得到所述堆土场的动态变化监测结果。
具体地，请参见图2，GIS平台融合了多个数据源的数据进行分析，包括时序植被状况数据、时序土地类型数据、现状数字高程模型（现状DEM）、现状数字正射影像（现状DOM）、历史多时相数字正射影像(历史多时相DOM)、历史多时相数字高程模型(历史多时相DEM)、地表形变速率图和时序沉降曲线，可提取出不同时段下堆土场的面积、体积、高程、厚度、坡度、形变等几何信息以及土地类型、植被覆盖等属性信息，进而反演出堆土场在各个时段的状态信息，作为堆土场监控结果。后续可进一步结合土地利用规划可确定该堆土场是否存在压占林地、耕地等环境破坏行为；结合气象等外部因素分析影响堆土场沉降状态的因素；结合历史多时序数字高程模型（历史多时序DEM）和现状数字高程模型（现状DEM）计算当前堆放弃土的厚度、体积和堆填速率等。
综上所述，本发明首先利用多源遥感影像确定出监测区域的潜在堆土场，然后利用无人机遥感技术获取潜在堆土场的现状高清影像，无人机的高清拍摄可以提供准确的现状数字高程模型（现状DEM）和现状数字正射影像（现状DOM），同时，在历史数据的获取上，获取第一光学遥感卫星提供的时序植被状况数据和时序土地类型数据，获取第二光学遥感卫星提供的历史多时相数字正射影像(历史多时相DOM)，以及获取星载平台提供的时序合成孔径雷达影像，得到历史多时相数字高程模型(历史多时相DEM)、地表形变速率图和时序沉降曲线。无人机提供了潜在堆土场的现状情况，而第一光学遥感卫星、第二光学遥感卫星和星载平台提供的数据可得到潜在堆土场的历史情况，不但在时间维度上可以互为补充，在感知内容上可也互为补充，完成堆土场的监测。同时，本发明结合大尺度变化监测和地表小尺度变化监测进行互补，可以分析堆土场在不同状态、不同尺度下的变化。这些互补的数据和技术共同实现堆土场的综合监测，最终获得全面的堆土场的动态变化监测结果。
在一种实施例中，如图3所示，基于上述堆土场动态变化监测方法，本发明还相应提供了一种堆土场动态变化监测装置，包括：
潜在堆土场确定模块100，用于获取第一光学遥感卫星采集的监测区域的时序多光谱遥感影像，以确定所述监测区域的时序植被状况数据和时序土地类型数据，并基于所述时序植被状况数据或所述时序土地类型数据确定所述监测区域内的潜在堆土场；
遥感成果生成模块200，用于获取无人机采集的所述潜在堆土场的现状高清影像集以得到现状数字高程模型和现状数字正射影像、获取第二光学遥感卫星提供的所述潜在堆土场的历史多时相数字正射影像、以及获取星载平台提供的所述潜在堆土场的时序合成孔径雷达影像以得到历史多时相数字高程模型、地表形变速率图和时序沉降曲线；
变化趋势信息生成模块300，用于基于所述时序植被状况数据、所述时序土地类型数据、所述现状数字高程模型、所述现状数字正射影像、所述历史多时相数字正射影像、所述历史多时相数字高程模型、所述地表形变速率图和所述时序沉降曲线得到所述潜在堆土场的变化趋势信息，将所述变化趋势信息和预设条件进行比对；
动态变化监测结果生成模块400，用于若所述变化趋势信息符合预设条件，将所述潜在堆土场确定为堆土场，利用GIS平台对所述时序植被状况数据、所述时序土地类型数据、所述现状数字高程模型、所述现状数字正射影像、所述历史多时相数字正射影像、所述历史多时相数字高程模型、所述地表形变速率图和所述时序沉降曲线进行融合分析，得到所述堆土场的动态变化监测结果。
需要说明的是，前述对堆土场动态变化监测方法实施例的解释说明也适用于本实施例的堆土场动态变化监测装置，此处不再赘述。
基于上述实施例，本发明还提供了一种终端，其原理框图可以如图4所示。上述终端包括通过系统总线连接的处理器、存储器、网络接口以及显示屏。其中，该终端的处理器用于提供计算和控制能力。该终端的存储器包括非易失性存储介质、内部存储器。该非易失性存储介质存储有操作系统和堆土场动态变化监测程序。该内部存储器为非易失性存储介质中的操作系统和堆土场动态变化监测程序的运行提供环境。该终端的网络接口用于与外部的终端通过网络连接通信。该堆土场动态变化监测程序被处理器执行时实现上述任意一种堆土场动态变化监测方法的步骤。该终端的显示屏可以是液晶显示屏或者电子墨水显示屏。
本领域技术人员可以理解，图4中示出的原理框图，仅是与本发明方案相关的部分结构的框图，并不构成对本发明方案所应用于其上的终端的限定，具体的终端可以包括比图中所示更多或更少的部件，或者组合某些部件，或者具有不同的部件布置。
在一个实施例中，提供了一种终端，上述终端包括存储器、处理器以及存储在上述存储器上并可在上述处理器上运行的堆土场动态变化监测程序，上述堆土场动态变化监测程序被上述处理器执行时实现本发明实施例提供的任意一种堆土场动态变化监测方法的步骤。
本发明实施例还提供一种计算机可读存储介质，上述计算机可读存储介质上存储有堆土场动态变化监测程序，上述堆土场动态变化监测程序被处理器执行时实现本发明实施例提供的任意一种堆土场动态变化监测方法的步骤。
应理解，上述实施例中各步骤的序号大小并不意味着执行顺序的先后，各过程的执行顺序应以其功能和内在逻辑确定，而不应对本发明实施例的实施过程构成任何限定。
所属领域的技术人员可以清楚地了解到，为了描述的方便和简洁，仅以上述各功能单元、模块的划分进行举例说明，实际应用中，可以根据需要而将上述功能分配由不同的功能单元、模块完成，即将上述装置的内部结构划分成不同的功能单元或模块，以完成以上描述的全部或者部分功能。实施例中的各功能单元、模块可以集成在一个处理单元中，也可以是各个单元单独物理存在，也可以两个或两个以上单元集成在一个单元中，上述集成的单元既可以采用硬件的形式实现，也可以采用软件功能单元的形式实现。另外，各功能单元、模块的具体名称也只是为了便于相互区分，并不用于限制本发明的保护范围。上述系统中单元、模块的具体工作过程，可以参考前述方法实施例中的对应过程，在此不再赘述。
在上述实施例中，对各个实施例的描述都各有侧重，某个实施例中没有详述或记载的部分，可以参见其它实施例的相关描述。
本领域普通技术人员可以意识到，结合本文中所公开的实施例描述的各实例的单元及算法步骤，能够以电子硬件、或者计算机软件和电子硬件的结合来实现。这些功能究竟是以硬件还是软件方式来执行，取决于技术方案的特定应用和设计约束条件。专业技术人员可以对每个特定的应用来使用不同法来实现所描述的功能，但是这种实现不应认为超出本发明的范围。
在本发明所提供的实施例中，应该理解到，所揭露的装置/终端设备和方法，可以通过其它的方式实现。例如，以上所描述的装置/终端设备实施例仅是示意性的，例如，上述模块或单元的划分，仅为一种逻辑功能划分，实际实现时可以由另外的划分方式，例如多个单元或组件可以结合或者可以集成到另一个系统，或一些特征可以忽略，或不执行。
以上所述实施例仅用以说明本发明的技术方案，而非对其限制；尽管参照前述实施例对本发明进行了详细的说明，本领域的普通技术人员应当理解；其依然可以对前述各实施例所记载的技术方案进行修改，或者对其中部分技术特征进行等同替换；而这些修改或者替换，并不是相应技术方案的本质脱离本发明各实施例技术方案的精神和范围，均应包含在本发明的保护范围之内。

Claims (10)

1.一种堆土场动态变化监测方法，其特征在于，所述方法包括：
获取第一光学遥感卫星采集的监测区域的时序多光谱遥感影像，以确定所述监测区域的时序植被状况数据和时序土地类型数据，并基于所述时序植被状况数据或所述时序土地类型数据确定所述监测区域内的潜在堆土场；
获取无人机采集的所述潜在堆土场的现状高清影像集以得到现状数字高程模型和现状数字正射影像、获取第二光学遥感卫星提供的所述潜在堆土场的历史多时相数字正射影像、以及获取星载平台提供的所述潜在堆土场的时序合成孔径雷达影像以得到历史多时相数字高程模型、地表形变速率图和时序沉降曲线；
基于所述时序植被状况数据、所述时序土地类型数据、所述现状数字高程模型、所述现状数字正射影像、所述历史多时相数字正射影像、所述历史多时相数字高程模型、所述地表形变速率图和所述时序沉降曲线得到所述潜在堆土场的变化趋势信息，将所述变化趋势信息和预设条件进行比对；
若所述变化趋势信息符合预设条件，将所述潜在堆土场确定为堆土场，利用GIS平台对所述时序植被状况数据、所述时序土地类型数据、所述现状数字高程模型、所述现状数字正射影像、所述历史多时相数字正射影像、所述历史多时相数字高程模型、所述地表形变速率图和所述时序沉降曲线进行融合分析，得到所述堆土场的动态变化监测结果。
2.根据权利要求1所述的堆土场动态变化监测方法，其特征在于，所述获取第一光学遥感卫星采集的监测区域的时序多光谱遥感影像，以确定所述监测区域的时序植被状况数据和时序土地类型数据，包括：
获取第一光学遥感卫星采集的监测区域的时序多光谱遥感影像，利用预设的坐标投影转换公式，将所述时序多光谱遥感影像从原始投影坐标系转换为第一目标坐标系，生成多个基于时间序列的第一多光谱遥感影像；
提取并组合每个所述第一多光谱遥感影像中预设波段的多个单一波段图像，得到每个所述第一多光谱遥感影像对应的多波段图像；
提取并组合每个所述第一多光谱遥感影像中的红光谱波段图像、绿光谱波段图像和蓝光谱波段图像，得到每个所述第一多光谱遥感影像对应的真彩色图像；
对每个所述多波段图像和每个所述真彩色图像均进行辐射定标、大气校正和边界区域裁剪，得到多个第一波段图像和多个第一真彩色图像；
对全部所述第一波段图像进行处理，得到所述监测区域的时序植被状况数据；
利用最大似然估计监督学习算法对全部所述第一真彩色图像进行处理，得到所述监测区域的时序土地类型数据。
3.根据权利要求2所述的堆土场动态变化监测方法，其特征在于，所述基于所述时序植被状况数据或所述时序土地类型数据确定所述监测区域内的潜在堆土场，包括：
基于所述时序植被状况数据或所述时序土地类型数据，确定所述监测区域内的突变区域；
将所述突变区域作为潜在堆土场。
4.根据权利要求1所述的堆土场动态变化监测方法，其特征在于，所述获取无人机采集的所述潜在堆土场的现状高清影像集以得到现状数字高程模型和现状数字正射影像，包括：
获取无人机采集的所述潜在堆土场的现状高清影像集，所述现状高清影像集中包括多个现状高清影像；
对每个所述现状高清影像均进行匀光匀色处理和影像畸变纠正，得到多个第一现状高清影像；
利用空中三角测量技术对每个所述第一现状高清影像进行特征点提取和多视影像联合平差，得到多个第二现状高清影像；
基于每个所述第二现状高清影像，利用密集匹配技术，得到所述潜在堆土场对应的密集点云；
基于所述密集点云，得到现状数字高程模型和现状数字正射影像。
5.根据权利要求1所述的堆土场动态变化监测方法，其特征在于，所述获取星载平台提供的所述潜在堆土场的时序合成孔径雷达影像以得到历史多时相数字高程模型、地表形变速率图和时序沉降曲线，包括：
获取星载平台提供的所述潜在堆土场的时序合成孔径雷达影像，所述时序合成孔径雷达影像由每个月所述星载平台拍摄的合成孔径雷达影像对组成，所述合成孔径雷达影像对是两个拍摄时间间隔为12天的合成孔径雷达影像；
按照预设规则对多个所述合成孔径雷达影像对进行选择，得到多个目标合成孔径雷达影像对，确定每个所述目标合成孔径雷达影像对的主影像和从影像，将每个所述从影像与对应的所述主影像进行配准，得到每个所述从影像对应的第一配准影像；
对每对所述主影像和所述第一配准影像进行干涉相位计算和相干性估计，生成多个干涉图，以及每个干涉图对应的相干图，对每个所述干涉图去除平地相位误差，生成多个第一干涉图，利用Goldstein滤波平滑每个所述第一干涉图，得到多个第二干涉图，基于每个所述相干图，利用最小费用流算法对应解缠每个所述第二干涉图，得到多个解缠后的相位图；
将全部所述相位图中的相位信息向高程转化，得到全部所述相位图中每个像素的高程信息，将全部所述相位图中每个像素的所述高程信息从雷达影像坐标系转换到第二目标坐标系，得到全部所述相位图中每个像素的第一高程信息；
采用双线性内插法对每个所述相位图中每个像素的所述第一高程信息进行重采样，得到历史多时相数字高程模型。
6.根据权利要求5所述的堆土场动态变化监测方法，其特征在于，所述获取星载平台提供的所述潜在堆土场的时序合成孔径雷达影像之后，还包括：
对所述潜在堆土场的所述时序合成孔径雷达影像进行前置滤波处理，得到第一时序合成孔径雷达影像，选取所述第一时序合成孔径雷达影像中成像效果符合预设标准的影像作为第一主影像，在所述时序合成孔径雷达影像中去除所述第一主影像，得到第一从影像集，所述第一从影像集中包括多个第一从影像；
将每个所述第一从影像与所述第一主影像进行配准，得到每个所述第一从影像对应的第二配准影像；
采用小基线原则将所述第一主影像和每个所述第二配准影像进行组合，生成多个第三干涉图，以及每个所述第三干涉图对应的第一相干图，所述第一相干图中包含对应的所述第三干涉图中不同点的相干性信息，对每个所述第三干涉图依次进行去平地相位误差、Goldstein滤波以及去地形相位后，得到多个差分干涉图，对每个所述差分干涉图进行多视处理，得到多个第一差分干涉图；
确定每个所述第一相干图中的高相干点，利用每个所述第一相干图中的全部所述高相干点，对每个所述第一相干图对应的每个所述第一差分干涉图进行相位解缠，得到多个第二差分干涉图，所述第二差分干涉图中包含每个所述高相干点对应的高相干区域内全部像素的初始解缠相位；
对每个所述第二差分干涉图估算对流层延迟相位，得到每个所述第二差分干涉图对应的区域对流层改正图，所述区域对流层改正图用于表征所述第二差分干涉图中每个像素的对流层延迟相位；
利用每个所述区域对流层改正图对每个所述区域对流层改正图对应的每个所述第二差分干涉图进行对流层延迟相位消除，得到多个第三差分干涉图，所述第三差分干涉图中包含每个所述高相干点对应的高相干区域内全部像素的第二解缠相位，所述第二解缠相位为对所述高相干区域内的每个像素的所述初始解缠相位消除对应的所述对流层延迟相位所得；
获取每个所述第三差分干涉图对应的所述第一主影像的成像时间和所述第一从影像的成像时间，基于每个所述第三差分干涉图对应的所述第一主影像的成像时间和所述第一从影像的成像时间，得到多个时间差矩阵；
基于每个所述第三差分干涉图中每个所述高相干点对应的高相干区域内全部像素的所述第二解缠相位和每个所述第三差分干涉图对应的所述时间差矩阵，得到时序合成孔径雷达影像在视线方向的形变速率，基于所述形变速率得到地表形变速率图和时序沉降曲线。
7.根据权利要求6所述的堆土场动态变化监测方法，其特征在于，所述对每个所述第二差分干涉图进行对流层延迟校正，得到每个所述第二差分干涉图对应的区域对流层改正图，包括：
获取所述潜在堆土场及所述潜在堆土场预设范围内的周边的气象信息，基于所述气象信息得到多个气象数据；
获取所述潜在堆土场中预设地点的地表高程和映射系数，所述映射系数用于将天顶对流层延迟投影至卫星视线方向；
将各个所述气象数据、所述地表高程和所述映射系数输入对流层延迟校正模型，得到所述潜在堆土场中每个所述预设地点的对流层延迟相位；
利用插值方法对所述潜在堆土场中全部所述预设地点的所述对流层延迟相位进行插值处理，得到每个所述第二差分干涉图对应的区域对流层改正图。
8.一种堆土场动态变化监测装置，其特征在于，包括：
潜在堆土场确定模块，用于获取第一光学遥感卫星采集的监测区域的时序多光谱遥感影像，以确定所述监测区域的时序植被状况数据和时序土地类型数据，并基于所述时序植被状况数据或所述时序土地类型数据确定所述监测区域内的潜在堆土场；
遥感成果生成模块，用于获取无人机采集的所述潜在堆土场的现状高清影像集以得到现状数字高程模型和现状数字正射影像、获取第二光学遥感卫星提供的所述潜在堆土场的历史多时相数字正射影像、以及获取星载平台提供的所述潜在堆土场的时序合成孔径雷达影像以得到历史多时相数字高程模型、地表形变速率图和时序沉降曲线；
变化趋势信息生成模块，用于基于所述时序植被状况数据、所述时序土地类型数据、所述现状数字高程模型、所述现状数字正射影像、所述历史多时相数字正射影像、所述历史多时相数字高程模型、所述地表形变速率图和所述时序沉降曲线得到所述潜在堆土场的变化趋势信息，将所述变化趋势信息和预设条件进行比对；
动态变化监测结果生成模块，用于若所述变化趋势信息符合预设条件，将所述潜在堆土场确定为堆土场，利用GIS平台对所述时序植被状况数据、所述时序土地类型数据、所述现状数字高程模型、所述现状数字正射影像、所述历史多时相数字正射影像、所述历史多时相数字高程模型、所述地表形变速率图和所述时序沉降曲线进行融合分析，得到所述堆土场的动态变化监测结果。
9.一种终端，其特征在于，所述终端包括：存储器、处理器及存储在所述存储器上并可在所述处理器上运行的堆土场动态变化监测程序，所述堆土场动态变化监测程序被所述处理器执行时实现如权利要求1-7任意一项所述的堆土场动态变化监测方法的步骤。
10.一种计算机可读存储介质，其特征在于，所述计算机可读存储介质存储有堆土场动态变化监测程序，所述堆土场动态变化监测程序能够被执行以用于实现如权利要求1-7任意一项所述的堆土场动态变化监测方法的步骤。
CN202410022651.4A 2024-01-08 2024-01-08 一种堆土场动态变化监测方法、装置、终端及存储介质 Active CN117523404B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202410022651.4A CN117523404B (zh) 2024-01-08 2024-01-08 一种堆土场动态变化监测方法、装置、终端及存储介质

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202410022651.4A CN117523404B (zh) 2024-01-08 2024-01-08 一种堆土场动态变化监测方法、装置、终端及存储介质

Publications (2)

Publication Number Publication Date
CN117523404A CN117523404A (zh) 2024-02-06
CN117523404B true CN117523404B (zh) 2024-03-22

Family

ID=89746156

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202410022651.4A Active CN117523404B (zh) 2024-01-08 2024-01-08 一种堆土场动态变化监测方法、装置、终端及存储介质

Country Status (1)

Country Link
CN (1) CN117523404B (zh)

Families Citing this family (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN117784075B (zh) * 2024-02-27 2024-05-10 中国铁路设计集团有限公司 一种融合雷达遥感和空间分析的岩溶发育强度分级方法

Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN107421514A (zh) * 2017-06-21 2017-12-01 广东中冶地理信息股份有限公司 一种垃圾填埋场范围变化监测方法
CN111723757A (zh) * 2020-06-28 2020-09-29 北京隆普智能科技有限公司 垃圾填埋场的监控方法及系统
CN113033381A (zh) * 2021-03-23 2021-06-25 生态环境部卫星环境应用中心 基于遥感数据的矿山和固废填埋场自动提取方法和装置
WO2021227423A1 (zh) * 2020-05-13 2021-11-18 深圳大学 基于动态基线的InSAR数字高程模型构建方法及系统

Patent Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN107421514A (zh) * 2017-06-21 2017-12-01 广东中冶地理信息股份有限公司 一种垃圾填埋场范围变化监测方法
WO2021227423A1 (zh) * 2020-05-13 2021-11-18 深圳大学 基于动态基线的InSAR数字高程模型构建方法及系统
CN111723757A (zh) * 2020-06-28 2020-09-29 北京隆普智能科技有限公司 垃圾填埋场的监控方法及系统
CN113033381A (zh) * 2021-03-23 2021-06-25 生态环境部卫星环境应用中心 基于遥感数据的矿山和固废填埋场自动提取方法和装置

Non-Patent Citations (3)

* Cited by examiner, † Cited by third party
Title
Xiaoqiong Qin et al..Stability assessment of high-speed railway using advanced InSAR technique.《2016 IEEE International Geoscience and Remote Sensing Symposium》.2016,第6010-6013页. *
徐十锋 等.无人机遥感技术在麻地梁煤矿中的精度分析及应用.《采矿技术》.2022,第22卷(第1期),第155-158页. *
李清泉 等. 海岸带地理环境遥感监测综述.《遥感学报》.2016,第1216-1229页. *

Also Published As

Publication number Publication date
CN117523404A (zh) 2024-02-06

Similar Documents

Publication Publication Date Title
CN113205475B (zh) 基于多源卫星遥感数据的森林高度反演方法
Kolecka et al. Assessment of the accuracy of SRTM C-and X-Band high mountain elevation data: A case study of the Polish Tatra Mountains
Sharma et al. Potential of airborne LiDAR data for terrain parameters extraction
CN113340191B (zh) 时间序列干涉sar的形变量测量方法及sar系统
CN117523404B (zh) 一种堆土场动态变化监测方法、装置、终端及存储介质
CN109100719B (zh) 基于星载sar影像与光学影像的地形图联合测图方法
CN110174673B (zh) 一种利用时序接力干涉图叠加高效减弱大气相位影响的方法
CN114820552A (zh) 一种利用光学卫星立体影像检测滑坡位移场的方法
CN111666896A (zh) 一种基于线性融合模型的遥感影像时空融合方法
Sefercik et al. Validation of Sentinel-1A and AW3D30 DSMs for the metropolitan area of Istanbul, Turkey
Pendyala et al. Dense DSM and DTM point cloud generation using CARTOSAT-2E satellite images for high-resolution applications
Sefercik et al. DEM extraction in urban areas using high-resolution TerraSAR-X imagery
Recla et al. From Relative to Absolute Heights in SAR-based Single-Image Height Prediction
Barrile et al. Analysis of hydraulic risk territories: comparison between LIDAR and other different techniques for 3D modeling
Ajayi et al. Modelling 3D Topography by comparing airborne LiDAR data with Unmanned Aerial System (UAS) photogrammetry under multiple imaging conditions
Deo et al. Evaluation of interferometric SAR DEMs generated using TanDEM-X data
Van Alphen A Technique-Based Approach to Structure-from-Motion: Applications to Human-Coastal Environments
Agrawal et al. Geospatial data for the himalayan region: requirements, availability, and challenges
Sohn et al. Shadow-effect correction in aerial color imagery
Wei et al. 3D digital elevation model generation
Torma et al. Sentinel-2 Images and Finnish Corine Land Cover Classification
Lombardi et al. Accuracy of high resolution CSK interferometric Digital Elevation Models
Piermattei et al. Analysis of glacial and periglacial processes using structure from motion.
Oshio et al. Generating DTM from DSM Using a Conditional GAN in Built-Up Areas
CN118230258B (zh) 一种顾及干涉基线优化的滑坡灾害早期识别方法

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