CN116109600A - 非监督的地基微变监测雷达图像的区域变化检测方法、装置及存储介质 - Google Patents
非监督的地基微变监测雷达图像的区域变化检测方法、装置及存储介质 Download PDFInfo
- Publication number
- CN116109600A CN116109600A CN202310139606.2A CN202310139606A CN116109600A CN 116109600 A CN116109600 A CN 116109600A CN 202310139606 A CN202310139606 A CN 202310139606A CN 116109600 A CN116109600 A CN 116109600A
- Authority
- CN
- China
- Prior art keywords
- images
- correlation
- image
- change
- change detection
- 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.)
- Pending
Links
- 238000001514 detection method Methods 0.000 title claims abstract description 87
- 238000012544 monitoring process Methods 0.000 title claims abstract description 72
- 230000008859 change Effects 0.000 claims abstract description 119
- 238000000034 method Methods 0.000 claims abstract description 49
- 238000000354 decomposition reaction Methods 0.000 claims description 19
- 238000012545 processing Methods 0.000 claims description 16
- 238000009825 accumulation Methods 0.000 claims description 15
- 230000009466 transformation Effects 0.000 claims description 11
- 230000000694 effects Effects 0.000 claims description 9
- 238000005259 measurement Methods 0.000 claims description 6
- 238000004364 calculation method Methods 0.000 description 8
- 230000008569 process Effects 0.000 description 8
- 238000004458 analytical method Methods 0.000 description 7
- 238000010586 diagram Methods 0.000 description 7
- 238000003384 imaging method Methods 0.000 description 4
- 238000005065 mining Methods 0.000 description 4
- 239000000654 additive Substances 0.000 description 3
- 230000000996 additive effect Effects 0.000 description 3
- 230000008602 contraction Effects 0.000 description 3
- 238000011160 research Methods 0.000 description 3
- 238000011156 evaluation Methods 0.000 description 2
- 238000001914 filtration Methods 0.000 description 2
- 238000012986 modification Methods 0.000 description 2
- 230000004048 modification Effects 0.000 description 2
- 238000011946 reduction process Methods 0.000 description 2
- 230000000717 retained effect Effects 0.000 description 2
- 230000009286 beneficial effect Effects 0.000 description 1
- 238000010276 construction Methods 0.000 description 1
- 230000007812 deficiency Effects 0.000 description 1
- 230000002708 enhancing effect Effects 0.000 description 1
- 229910052500 inorganic mineral Inorganic materials 0.000 description 1
- 238000005305 interferometry Methods 0.000 description 1
- 230000001788 irregular Effects 0.000 description 1
- 238000003064 k means clustering Methods 0.000 description 1
- 230000014759 maintenance of location Effects 0.000 description 1
- 239000011707 mineral Substances 0.000 description 1
- 239000002689 soil Substances 0.000 description 1
Images
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T7/00—Image analysis
- G06T7/0002—Inspection of images, e.g. flaw detection
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T5/00—Image enhancement or restoration
- G06T5/70—Denoising; Smoothing
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T5/00—Image enhancement or restoration
- G06T5/90—Dynamic range modification of images or parts thereof
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06V—IMAGE OR VIDEO RECOGNITION OR UNDERSTANDING
- G06V10/00—Arrangements for image or video recognition or understanding
- G06V10/70—Arrangements for image or video recognition or understanding using pattern recognition or machine learning
- G06V10/762—Arrangements for image or video recognition or understanding using pattern recognition or machine learning using clustering, e.g. of similar faces in social networks
-
- 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/10044—Radar image
Landscapes
- Engineering & Computer Science (AREA)
- Theoretical Computer Science (AREA)
- Physics & Mathematics (AREA)
- General Physics & Mathematics (AREA)
- Computer Vision & Pattern Recognition (AREA)
- Computing Systems (AREA)
- Health & Medical Sciences (AREA)
- Artificial Intelligence (AREA)
- Quality & Reliability (AREA)
- Databases & Information Systems (AREA)
- Evolutionary Computation (AREA)
- General Health & Medical Sciences (AREA)
- Medical Informatics (AREA)
- Software Systems (AREA)
- Multimedia (AREA)
- Image Analysis (AREA)
Abstract
本公开涉及非监督的地基微变监测雷达图像的区域变化检测方法、装置及存储介质，方法包括：获取目标区域一段时间序列的图像；得到时间序列等间隔图像的相关性；基于相关性图像，获得初始差异图；对处理后的初始差异图进行聚类，获得变化检测结果。通过本公开的各实施例，克服了地基微变监测雷达噪声多、变化前后图像难以区分、差异图难获取的困难。
Description
技术领域
本公开涉及地基微波遥感领域，具体涉及一种非监督的地基微变监测雷达图像的区域变化检测方法、非监督的地基微变监测雷达图像的区域变化检装置及计算机可读存储介质。
背景技术
传统差值法、比值法、对数比值法等对变化前后两幅图像进行简单计算获取差异图，很难区分变化类与非变化类，无法获取高质量的差异图，从而影响了变化检测的结果。
发明内容
本公开意图提供一种非监督的地基微变监测雷达图像的区域变化检测方法、非监督的地基微变监测雷达图像的区域变化检装置及计算机可读存储介质，克服了地基微变监测雷达噪声多、变化前后图像难以区分、差异图难获取的困难。
根据本公开的方案之一，提供一种非监督的地基微变监测雷达图像的区域变化检测方法，包括：
获取目标区域一段时间序列的图像；
得到时间序列等间隔图像的相关性；
基于相关性图像，获得初始差异图；
对处理后的初始差异图进行聚类，获得变化检测结果。
在一些实施例中，其中，所述所述获取目标区域一段时间序列的图像，包括：
获取一组时间顺序的等间隔微变监测雷达图像。
在一些实施例中，其中，所述得到时间序列等间隔图像的相关性，包括：
获取时间顺序的等间隔的每两幅图像的相关性图像。
在一些实施例中，其中，所述基于相关性图像，获得初始差异图，包括：
基于获取的相关性图像取相邻的每两幅图像，得到相关对数比值图；
对相关对数比值图进行分组时序累加运算，以消除噪声。
在一些实施例中，其中，初始差异图的处理，包括：
小波多尺度变换对差异图进行多层分解，将初始差异图分解为低频分量与高频分量，依次对上一层得到的低频分量分解为低频分量和高频分量，直到得到预设效果的分解层数；
对分解后的高频分量进行阈值去噪和系数重构，获取处理后的初始差异图。
在一些实施例中，其中，所述进行聚类，包括：
取预设尺寸窗口遍历处理后的初始差异图的每一个像素点，并取中心像素点邻域内平均值；
对邻域内平均值进行排序，选出其中最大值与最小值，分别对邻域内所有的数据点排序，将数据值最大的数据点和最小的数据点分别作为类簇的初始聚类中心点；
输入无标签的数据，以距离作为相似度的衡量标准，根据不同数据之间的相似程度将数据聚类成不同的类簇；
将距离平方和作为准则函数，将每个数据点分配到聚类中心对应的簇中，重新计算出每一个簇中所有数据的平均值作为该簇新的聚类中心点，直至聚类中心点基本不变或者是准则函数收敛到设定的阈值。
在一些实施例中，其中，所述获得变化检测结果，包括：
对迭代计算得到的聚类中心所属的类簇进行二值化处理，得到最终的变化检测结果图。
根据本公开的方案之一，提供一种非监督的地基微变监测雷达图像的区域变化检测装置，包括：
获取模块，其配置为获取目标区域一段时间序列的图像；
变化检测模块，其配置为：
得到时间序列等间隔图像的相关性；
基于相关性图像，获得初始差异图；
对处理后的初始差异图进行聚类，获得变化检测结果。
在一些实施例中，其中，所述变化检测模块，进一步配置为：
基于获取的相关性图像取相邻的每两幅图像，得到相关对数比值图；
对相关对数比值图进行分组时序累加运算，以消除噪声。
根据本公开的方案之一，提供计算机可读存储介质，其上存储有计算机可执行指令，所述计算机可执行指令由处理器执行时，实现：
根据上述的非监督的地基微变监测雷达图像的区域变化检测方法。
本公开的各种实施例的非监督的地基微变监测雷达图像的区域变化检测方法、非监督的地基微变监测雷达图像的区域变化检装置及计算机可读存储介质，至少获取目标区域一段时间序列的图像；得到时间序列等间隔图像的相关性；基于相关性图像，获得初始差异图；对处理后的初始差异图进行聚类，获得变化检测结果，旨在基于获取的多幅图像取相同时间间隔的每两幅图像按时间顺序依次计算相关性，然后对时序相邻的每两幅相关性图像进行处理后的图像得到了非线性收缩，能够将乘性噪声转化为加性噪声，从而达到有效的去噪和增强变化和非变化区域的对比度，不同于传统的基于两幅雷达图像的变化检测，而是采用了一段时间序列的多幅图像获取差异图，克服了地基微变监测雷达噪声多、变化前后图像难以区分、差异图难获取的困难。进一步的，结合了地基微变监测雷达图像的滑坡区域在时间与空间上的特性改进了传统聚类算法，使得算法保证识别准确度的同时又减少了计算的迭代次数，简单高效。
应当理解，前面的大体描述以及后续的详细描述只是示例性的和说明性的，并非对所要求保护的本公开的限制。
附图说明
在未必按照比例绘制的附图中，不同视图中相似的附图标记可以表示相似的构件。具有字母后缀的相似附图标记或具有不同字母后缀的相似附图标记可以表示相似构件的不同实例。附图通常作为示例而非限制地图示各种实施例，并且与说明书和权利要求书一起用于解释所公开的实施例。
图1示出了本公开实施例的非监督的地基微变监测雷达图像的区域变化检测方法的一种流程图；
图2示出了传统聚类算法得到的变化检测结果图；
图3示出了本公开各实施例提供的方法获取的变化检测结果图；
图4示出了本公开实施例的非监督的地基微变监测雷达图像的区域变化检测装置的一种架构图。
具体实施方式
为了使得本公开实施例的目的、技术方案和优点更加清楚，下面将结合本公开实施例的附图，对本公开实施例的技术方案进行清楚、完整地描述。显然，所描述的实施例是本公开的一部分实施例，而不是全部的实施例。基于所描述的本公开的实施例，本领域普通技术人员在无需创造性劳动的前提下所获得的所有其他实施例，都属于本公开保护的范围。
很多地区地形复杂多样，矿产资源丰富的矿区开采过程中边坡多呈阶梯状分布，另外受现场施工活动和气候降雨等影响，易发生滑坡和崩塌等灾害，如不能及时预警和评估，存在严重安全隐患威胁了工程人员的生命安全。地基微变监测雷达通过差分干涉测量原理可以实现对易发生微小形变物体的监测，其不仅具有全天时、全天候、大范围、远距离、非接触的特点，还具有高分辨率、分钟级图像获取、设站灵活便于携带的优势，广泛应用于地质灾害的监测和预警工作。图像的变化检测是指对拍摄于不同的时间、相同的地点的多幅图像进行分析并提取变化和未变化信息的过程，利用其分析和评估灾害信息以及土表变化监测等方面发挥着重要作用。目前，研究人员利用各种算法对星载SAR图像目标区域的变化检测进行了广泛的研究。但，星载SAR图像适用于大范围的对地成像，不适用于矿区这类小范围滑坡的目标区域变化检测任务，且与星载SAR图像变化检测相比，微变监测雷达图像的变化检测研究较为缺乏，因此研究微变监测雷达图像的变化检测方法有着重要意义，不仅能为雷达监测预警矿区滑坡工作提供重要参考信息，保障矿区施工人员财产的安全，也能在微变监测雷达图像变化检测领域提供技术参考。
变化检测技术一般可以分为监督式和非监督式。监督变化检测方法往往需要准确的先验信息，非监督变化检测算法对差异图的要求较高。但，微变监测雷达由于其独特的成像原理，图像会附有大量的乘性噪声，通常导致变化前后图像先验知识难获取，变化检测任务困难，因此非监督式更适用于微变监测雷达图像的变化检测，这过程需要获取较高质量的差异图，而传统的差异图获取方法受噪声影响较大，会导致差异图的初步分类精度较低，从而影响最终的聚类结果。因此研究一种非监督的地基微变监测雷达图像的滑坡区域变化检测方法，其中要研究高质量差异图的获取以及差异图聚类的算法。
结合前文背景技术部分内容记载，本公开以实施例方式示例性记载了相应解决方案以解决现有技术中存在的缺陷，但不作为对本公开所要求专利权保护范围的限定。
作为方案之一，如图1所示，本公开的实施例提供了一种非监督的地基微变监测雷达图像的区域变化检测方法，包括：
获取目标区域一段时间序列的图像；
得到时间序列等间隔图像的相关性；
基于相关性图像，获得初始差异图；
对处理后的初始差异图进行聚类，获得变化检测结果。
针对前文所提出的问题，本公开各实施例旨在提出了一种非监督的地基微变监测雷达图像的区域变化检测方法，例如应用于滑坡区域的变化检测。该方法分为获取目标区域一段时间序列的多幅图像；计算时间序列等间隔图像的相关性，在相关性图像基础上计算相邻两幅图像的相关对数比值图并在时间序列上作累加运算获得初始差异图；利用小波多分辨率分析对差异图进行分解并对高频阈值去噪、重构小波系数获取最终差异图；用改进的聚类算法对最终差异图进行聚类获得滑坡区域变化检测结果，例如，得到一幅二值图(变化检测结果)四个部分呈现相应的变化检测结果。该方法不同于传统的基于两幅雷达图像的变化检测，而是采用了一段时间序列的多幅图像获取差异图，克服了地基微变监测雷达噪声多、变化前后图像难以区分、差异图难获取的困难，并且结合了地基微变监测雷达图像的滑坡区域在时间与空间上的特性改进了传统聚类算法，使得算法保证识别准确度的同时又减少了计算的迭代次数，简单高效。
举例来说，具体的成像过程可以实施为以如下步骤为代表的技术方案。
在一些实施方案中，本公开实施例的方法，可以为：所述获取目标区域一段时间序列的图像，包括：
获取一组时间顺序的等间隔微变监测雷达图像。
示例性的，本实施例的方法可以包括：
步骤S1：获取一组时间顺序的等间隔微变监测雷达图像I1(i,j),I2(i.j),...,In(i,j)，时序内每幅图像的采集条件不同。
这里，每幅图像是预先确定的含有目标区域的图像，具体地，目标区域为用户需要确定是否存在变化的区域，例如发生滑坡、崩塌等地质灾害之后，通过变化检测，确定A区域是否发生变化，则A区域即为目标区域。
其中，每幅图像的采集条件不同，该采集条件至少包括采集时间，即基于地基微变监测雷达获取的一组时序的等间隔雷达图像为在不同时刻采集的。
本公开的构思之一可以在于，获取时间序列的多幅地基微变监测雷达图像，要求图像大小、位置均相同，确保每幅图像的相同坐标的像素点代表的区域相同。
在一些实施方案中，本公开实施例的方法，可以为：
所述得到时间序列等间隔图像的相关性，包括：
获取时间顺序的等间隔的每两幅图像的相关性图像。
在一些实施方案中，本公开实施例的方法，可以为：
所述基于相关性图像，获得初始差异图，包括：
基于获取的相关性图像取相邻的每两幅图像，得到相关对数比值图；
对相关对数比值图进行分组时序累加运算，以消除噪声。
示例性的，本实施例的方法可以继续包括：
步骤S2：生成目标区域已发生变化的地基微变监测雷达图像的初始差异图，具体处理步骤为：
步骤S21：计算获取时间顺序的等间隔的每两幅图像的相关性图像Rt1(i,j),Rt2(i.j),...,Rtn(i,j)：
式中，In-1(i,j)和In(i,j)为经过配准滤波的时序相邻的微变监测雷达图像，Rtn(i,j)为n时刻的相关性图像，且其每个像素值在0～1之间，0表示完全不相关，1表示完全相关。
在一些实际应用场景和检测过程中，相关性可以衡量两幅图像的相似程度，但，在计算求取图像相关性的过程中会严重受到噪声的干扰，从而影响相关性图像的质量，因此仅计算两个时刻的地基微变监测雷达图像的相关性无法将噪声引起的失相干与真正滑坡引起的失相干区域区分开来。传统变化检测多基于两幅图像获取差异图，且原图变化信息明显。但，地基微变监测雷达图像由于其固有成像原理，乘性散斑噪声较多，且变化信息不明显，仅通过两幅地基微变监测雷达图像利用现有的差异图生成方法无法获取有效的差异图，噪声引起的变化信息与滑坡引起的变化信息无法区分。
有鉴于此，基于获取的多幅图像取相同时间间隔的每两幅图像按时间顺序依次计算相关性，然后对时序相邻的每两幅相关性图像进行对数比值运算获得相关对数比值图并对其进行累加运算，累加运算后的图像得到了非线性收缩，能够将乘性噪声转化为加性噪声，从而达到有效的去噪和增强变化和非变化区域的对比度。
相比于现有的变化检测方法，本公开的一些实施例中，采取时间序列上的多幅雷达图像，计算获取相关对数比值图并在时间顺序上作累加运算，生成的地基微变监测雷达图像的差异图较现有方法相比对噪声影响更具有鲁棒性，更容易聚类。
示例性的，本实施例的方法可以在前述步骤S21基础上继续包括：
步骤S22：对获取的相关性图像Rt1(i,j),Rt2(i.j),...,Rtn(i,j)取相邻的每两幅图像作对数比值运算得到相关对数比值图：
式中，RLogn(i,j)为第n时刻的相关对数比值图；
步骤S23：对得到的所有相关对数比值图进行分组时序累加运算，具体的，在时序的地基微变监测雷达图像中，形变滑坡体的相关性在时间轴上会有一段逐渐减小的平滑趋势，相关性依次降低，表示为：
Rt1＞Rt2＞…＞Rtn-1＞Rtn且0＜Rt1,Rt2,...,Rtn-1,Rtn＜1 (3)
无变化区域的相关性表示为：
R't1＞R't2＞…＞R'tn-1＞R'tn
分别对其累加运算如下：
式中，σ为噪声的相干系数值。
由于地基雷达图像固有的乘性散斑噪声的存在，在某一时刻的某一像素的相关性会发生突变的降低或升高，而滑坡区域的相关性在时序上会有一段逐渐减小的过程，通过对数累加运算能有效的消除部分噪声，使其累加后的相关性趋近于0，而形变滑坡体的相关性会被拉伸增大，增强了差异图中变化区域与未变化区域的对比度，得到较好质量的初始差异图，表示为I(i,j)。
相关对数比累加运算得到的差异图虽能消除一部分乘性散斑噪声，但地基微变监测雷达图像的噪声是随机无规律的，仍会有一些噪声点的存在，影响变化检测结果。
有鉴于此，利用小波多分辨率分析对差异图进一步的去噪。在一些实施例中，初始差异图的处理，包括：
小波多尺度变换对差异图进行多层分解，将初始差异图分解为低频分量与高频分量，依次对上一层得到的低频分量分解为低频分量和高频分量，直到得到预设效果的分解层数；
对分解后的高频分量进行阈值去噪和系数重构，获取处理后的初始差异图。
由于地基雷达图像的变化区域，例如滑坡区域，在图像中属于一个突变的现象，滑坡区域的边界和噪声集中在图像的高频部分，而整体区域在时间轴上属于图像的低频部分，小波多尺度变换既能对低频部分进行提取分析，同时也能对高频部分进行去噪分析。实际应用中，可以将高频部分视为对应为一些细节信息，可以将低频部分视为对应一些近似信息。
基于小波多尺度变换，本公开的各实施例结合地基微变监测雷达图像滑坡区域的变化特点，采用小波多尺度变换分析差异图，既能对高频分量进一步去噪，也能提取保留低频分量的滑坡区域，最后对系数进行重构获取的最终差异图，能更好的保留滑坡区域的细节信息(轮廓形状)。
示例性的，本实施例的方法可以在前述步骤S2的基础上继续包括：
步骤S3：利用小波多分辨率分析初始差异图。具体的分析步骤为：
步骤S31：小波多尺度变换对差异图进行多层分解，将差异图分解为低频分量与高频分量，再依次对上一层得到的低频分量再次分解为低频分量和高频分量，直到得到需要效果的分解层数即可。n层分解表达如下：
f(I(i,j))＝CAn+CDn+CDn-1+…+CD2+CD1 (6)
其中，I(i,j)代表差异图，CAn代表低频近似部分，CDn代表高频细节部分，n代表分解层数
具体为差异图经过小波变换后分解成四幅子图，分别是：
其中，Φ(t)为尺度函数；
为小波函数或母小波；*为垂直方向上的卷积；·表示为水平方向上的卷积。
为原图像的近似图，
为垂直和45°方向上的高频分量，
为水平和45°方向上的高频分量，
为45°方向上的高频分量。
步骤S32：对分解后的高频分量进行阈值去噪和系数重构获取最终差异图：
具体为含噪声的差异图在各尺度上进行小波分解后，图像的能量主要集中在低频分量上，而噪声信号的能量主要分布在各阶高频分量上。因而，原始差异图变化信息的小波系数绝对值较大，噪声信息小波系数的绝对值较小，在这种情况下，通过设定一个合适的阈值门限，采用阈值办法保留有用的信号系数，从而将差异图像变化信息与噪声区分开来。表示为：
其中，y为含噪声的小波变化系数；T为去噪阈值；Tsoft为软阈值滤波的收缩函数；sgn(y)为符号函数。
将最终分解层的低频系数与阈值化后的各层高频系数进行重构，抑制噪声的同时较好的保留图像的细节信息，从而得到较好质量的最终差异图。
通过上述具体实施方式，应当认为，差异图和噪声在经过小波变换后具有不同的特性，图像的能量主要集中在低频子带上，而噪声信号的能量主要分布在各阶高频子带上。因而，原始差异图信息的小波系数绝对值较大，噪声信息小波系数的绝对值较小，在这种情况下，通过设定一个合适的阈值门限，采用阈值办法保留有用的信号系数，从而将差异图像信息与噪声区分开来。通过重构阈值去噪处理后的高频系数和最终分解层的低频系数，抑制图像噪声的同时保留了较好的的细节信息，从而得到较好质量的最终差异图。
在一些实施方案中，本公开实施例的方法，可以为：
所述进行聚类，包括：
取预设尺寸窗口遍历处理后的初始差异图的每一个像素点，并取中心像素点邻域内平均值；
对邻域内平均值进行排序，选出其中最大值与最小值，分别对邻域内所有的数据点排序，将数据值最大的数据点和最小的数据点分别作为类簇的初始聚类中心点；
输入无标签的数据，以距离作为相似度的衡量标准，根据不同数据之间的相似程度将数据聚类成不同的类簇；
将距离平方和作为准则函数，将每个数据点分配到聚类中心对应的簇中，重新计算出每一个簇中所有数据的平均值作为该簇新的聚类中心点，直至聚类中心点基本不变或者是准则函数收敛到设定的阈值。
传统的K-means聚类算法是采用随机获取初始聚类中心的规则，随机选取聚类中心的规则容易增加迭代计算的次数，且选取的不同初始聚类中心可能会有不同的聚类效果，聚类结果不准确且不稳定。
本公开基于较好质量的最终差异图，进一步利用改进的聚类算法，对最终差异图进行分类，并将分类得到的每个像素进行二值化处理，确定地基微变监测雷达图像的滑坡区域。本公开涉及的各实施例提供的聚类算法对初始聚类中心选取规则进行改进，结合滑坡区域在地基微变监测雷达图像上的时间与空间特性的分析，即变化区域呈区域集中体现且数据值较未变化的区域与噪声相比数据值更高的特点，在选取初始聚类中心规则加入邻域判决项，聚类迭代次数更少，减少计算量，聚类结果更准确且稳定。
示例性的，本实施例的方法可以在前述步骤S3的基础上继续包括：
步骤S4：利用改进的聚类算法对最终差异图进行分类，并对分类得到的每个像素进行二值化处理，确定图像的变化区域，其中，具体的处理步骤为：
步骤S41：传统聚类算法随机选取聚类中心可能会增加迭代次数，增加计算量，而且不同的聚类中心往往会有不一样的聚类效果，不具有稳定性，因此对其初始聚类中心的选取进行改进。
具体的，结合本文获取的差异图的空间性与数值特性，考虑到差异图中变化区域是呈区域集中体现且数据值较未变化的区域与噪声相比数据值更高的特点，对初始聚类中心的选取规则进行改进，再选取聚类中心的过程中加入邻域均值判决项获取算法的初始聚类中心。
具体的，取r×r窗口遍历差异图的每一个像素点并取中心像素点邻域内平均值NMn，表示为：
其中，NMn表示第n个像素点的邻域均值，Nrn表示第n个像素点的r×r的邻域集合。
步骤S42：对(9)得到的邻域均值进行排序，选出其中最大值Max(NMn)与最小值Min(NMn)，然后分别对邻域内所有的数据点排序，将Max(NMn)中数据值最大的数据点和Min(NMn)中最小的数据点分别作为类簇的初始聚类中心点。
cj(1)＝Max(Max(NMn)) (10)
cj(1)＝Min(Min(NMn)) (11)
步骤S43：输入无标签的数据，以距离作为相似度的衡量标准，根据不同数据之间的相似程度将数据聚类成不同的类簇。计算像素之间的欧式距离为d(xi,cj(n))，用来计算图像中每一个像素xi到聚类中心cj的距离公式如下：
其中，xi(i＝1,2,...,n)表示属于类簇j中的像素点；聚类中心点cj(n)(j＝1,2,...,k)表示图像中的类簇，n表示迭代运算次序号，即类簇j中数据的平均值。
步骤S44：将距离平方和j(c)作为准则函数，目标是使各簇成员到其聚类中心的距离的平方和最小，将每个数据点分配到k个聚类中心对应的簇中，从而得到最好的聚类结果。目标函数表示为：
其中，c为聚类中心的集合，共有k个聚类中心，Sj(n)为类簇。
步骤S44：重新计算出每一个簇中所有数据的平均值作为该簇新的聚类中心点，计算如下：
步骤S45：重复(13)步与(14)步，直到聚类中心点基本不变或者是准则函数收敛到设定的阈值，至此算法结束。
在一些实施方案中，本公开实施例的方法，可以为：
所述获得变化检测结果，包括：
对迭代计算得到的聚类中心所属的类簇进行二值化处理，得到最终的变化检测结果图。
示例性的，本实施例的方法可以在前述步骤S45的基础上继续包括：
步骤S46：最后采用如式(15)对算法n次迭代计算的最终聚类中心cj(n)所属的类簇Sj(n)进行二值化处理，得到最终的变化检测结果图。
其中，S1(n)为变化类，S2(n)为未变化类，Img_change为变化检测结果图。
图2和图3为利用实际地基微变监测雷达图像所得到的变化检测结果图，其中，浅色区域为检测的变化区域(滑坡区域)；深色区域为未变化区域。图中坐标：Range为范围，Azimuth为方位角。
图2为传统聚类算法得到的变化检测结果图。图中圈出部分显示含有较多的噪声点被错误聚类为变化区域，且变化区域(滑坡区域)的轮廓信息模糊，细节保留度不高，变化检测的结果不准确；图3为本公开各实施例提供的方法获取的变化检测结果图。通过图3对比图2可以明显看出本公开获取的变化检测结果图保留了图像变化区域(滑坡区域)细节的同时噪声点更少，聚类结果更准确。
本实施例通过结合了地基微变监测雷达图像的滑坡区域在时间与空间上的特性改进了传统聚类算法，使得算法保证识别准确度的同时又减少了计算的迭代次数，简单高效地通过二值化处理呈现了滑坡区域变化检测结果。
作为方案之一，如图4所示，本公开的实施例提供了一种非监督的地基微变监测雷达图像的区域变化检测装置，包括：
获取模块，其配置为获取目标区域一段时间序列的图像；
变化检测模块，其配置为：
得到时间序列等间隔图像的相关性；
基于相关性图像，获得初始差异图；
对处理后的初始差异图进行聚类，获得变化检测结果。
作为一种实施方式，本公开的非监督的地基微变监测雷达图像的区域变化检测装置，可以结合前文所述各步骤，所述所述变化检测模块，进一步配置为：
基于获取的相关性图像取相邻的每两幅图像，得到相关对数比值图；
对相关对数比值图进行分组时序累加运算，以消除噪声。
作为一种实施方式，本公开的非监督的地基微变监测雷达图像的区域变化检测装置，可以结合前文所述各步骤，所述所述变化检测模块，进一步配置为：
所述进行聚类，包括：
取预设尺寸窗口遍历处理后的初始差异图的每一个像素点，并取中心像素点邻域内平均值；
对邻域内平均值进行排序，选出其中最大值与最小值，分别对邻域内所有的数据点排序，将数据值最大的数据点和最小的数据点分别作为类簇的初始聚类中心点；
输入无标签的数据，以距离作为相似度的衡量标准，根据不同数据之间的相似程度将数据聚类成不同的类簇；
将距离平方和作为准则函数，将每个数据点分配到聚类中心对应的簇中，重新计算出每一个簇中所有数据的平均值作为该簇新的聚类中心点，直至聚类中心点基本不变或者是准则函数收敛到设定的阈值。
示例性的，本实施例的非监督的地基微变监测雷达图像的区域变化检测装置的各模块以及涉及的子模块，可以进一步配置为对应实现以下步骤：
步骤S1：获取一组时间顺序的等间隔微变监测雷达图像I1(i,j),I2(i.j),...,In(i,j)，时序内每幅图像的采集条件不同。
这里，每幅图像是预先确定的含有目标区域的图像，具体地，目标区域为用户需要确定是否存在变化的区域，例如发生滑坡、崩塌等地质灾害之后，通过变化检测，确定A区域是否发生变化，则A区域即为目标区域。
其中，每幅图像的采集条件不同，该采集条件至少包括采集时间，即基于地基微变监测雷达获取的一组时序的等间隔雷达图像为在不同时刻采集的。
步骤S2：生成目标区域已发生变化的地基微变监测雷达图像的初始差异图，具体处理步骤为：
步骤S21：计算获取时间顺序的等间隔的每两幅图像的相关性图像Rt1(i,j),Rt2(i.j),...,Rtn(i,j)：
式中，In-1(i,j)和In(i,j)为经过配准滤波的时序相邻的微变监测雷达图像，Rtn(i,j)为n时刻的相关性图像，且其每个像素值在0～1之间，0表示完全不相关，1表示完全相关；
步骤S22：对获取的相关性图像Rt1(i,j),Rt2(i.j),...,Rtn(i,j)取相邻的每两幅图像作对数比值运算得到相关对数比值图：
式中，RLogn(i,j)为第n时刻的相关对数比值图；
步骤S23：对得到的所有相关对数比值图进行分组时序累加运算，具体的，在时序的地基微变监测雷达图像中，形变滑坡体的相关性在时间轴上会有一段逐渐减小的平滑趋势，相关性依次降低，表示为：
Rt1＞Rt2＞…＞Rtn-1＞Rtn且0＜Rt1,Rt2,...,Rtn-1,Rtn＜1 (18)
无变化区域的相关性表示为：
R't1＞R't2＞…＞R'tn-1＞R'tn
分别对其累加运算如下：
式中，σ为噪声的相干系数值。
由于地基雷达图像固有的乘性散斑噪声的存在，在某一时刻的某一像素的相关性会发生突变的降低或升高，而滑坡区域的相关性在时序上会有一段逐渐减小的过程，通过对数累加运算能有效的消除部分噪声，使其累加后的相关性趋近于0，而形变滑坡体的相关性会被拉伸增大，增强了差异图中变化区域与未变化区域的对比度，得到较好质量的初始差异图，表示为I(i,j)。
步骤S3：利用小波多分辨率分析初始差异图。具体的分析步骤为：
步骤S31：小波多尺度变换对差异图进行多层分解，将差异图分解为低频分量与高频分量，再依次对上一层得到的低频分量再次分解为低频分量和高频分量，直到得到需要效果的分解层数即可。n层分解表达如下：
f(I(i,j))＝CAn+CDn+CDn-1+…+CD2+CD1 (21)
其中，I(i,j)代表差异图，CAn代表低频近似部分，CDn代表高频细节部分，n代表分解层数
具体为差异图经过小波变换后分解成四幅子图，分别是：
其中，Φ(t)为尺度函数；
为小波函数或母小波；*为垂直方向上的卷积；·表示为水平方向上的卷积。
为原图像的近似图，
为垂直和45°方向上的高频分量，
为水平和45°方向上的高频分量，
为45°方向上的高频分量。
步骤S32：对分解后的高频分量进行阈值去噪和系数重构获取最终差异图：
具体为含噪声的差异图在各尺度上进行小波分解后，图像的能量主要集中在低频分量上，而噪声信号的能量主要分布在各阶高频分量上。因而，原始差异图变化信息的小波系数绝对值较大，噪声信息小波系数的绝对值较小，在这种情况下，通过设定一个合适的阈值门限，采用阈值办法保留有用的信号系数，从而将差异图像变化信息与噪声区分开来。表示为：
其中，y为含噪声的小波变化系数；T为去噪阈值；Tsoft为软阈值滤波的收缩函数；sgn(y)为符号函数。
将最终分解层的低频系数与阈值化后的各层高频系数进行重构，抑制噪声的同时较好的保留图像的细节信息，从而得到较好质量的最终差异图。
步骤S4：利用改进的聚类算法对最终差异图进行分类，并对分类得到的每个像素进行二值化处理，确定图像的变化区域，其中，具体的处理步骤为：
步骤S41：传统聚类算法随机选取聚类中心可能会增加迭代次数，增加计算量，而且不同的聚类中心往往会有不一样的聚类效果，不具有稳定性，因此对其初始聚类中心的选取进行改进。
具体的，结合本文获取的差异图的空间性与数值特性，考虑到差异图中变化区域是呈区域集中体现且数据值较未变化的区域与噪声相比数据值更高的特点，对初始聚类中心的选取规则进行改进，再选取聚类中心的过程中加入邻域均值判决项获取算法的初始聚类中心。
具体的，取r×r窗口遍历差异图的每一个像素点并取中心像素点邻域内平均值NMn，表示为：
其中，NMn表示第n个像素点的邻域均值，Nrn表示第n个像素点的r×r的邻域集合。
步骤S42：对(9)得到的邻域均值进行排序，选出其中最大值Max(NMn)与最小值Min(NMn)，然后分别对邻域内所有的数据点排序，将Max(NMn)中数据值最大的数据点和Min(NMn)中最小的数据点分别作为类簇的初始聚类中心点。
cj(1)＝Max(Max(NMn)) (25)
cj(1)＝Min(Min(NMn)) (26)
步骤S43：输入无标签的数据，以距离作为相似度的衡量标准，根据不同数据之间的相似程度将数据聚类成不同的类簇。计算像素之间的欧式距离为d(xi,cj(n))，用来计算图像中每一个像素xi到聚类中心cj的距离公式如下：
其中，xi(i＝1,2,...,n)表示属于类簇j中的像素点；聚类中心点cj(n)(j＝1,2,...,k)表示图像中的类簇，n表示迭代运算次序号，即类簇j中数据的平均值。
步骤S44：将距离平方和j(c)作为准则函数，目标是使各簇成员到其聚类中心的距离的平方和最小，将每个数据点分配到k个聚类中心对应的簇中，从而得到最好的聚类结果。目标函数表示为：
其中，c为聚类中心的集合，共有k个聚类中心，Sj(n)为类簇。
步骤S44：重新计算出每一个簇中所有数据的平均值作为该簇新的聚类中心点，计算如下：
步骤S45：重复(13)步与(14)步，直到聚类中心点基本不变或者是准则函数收敛到设定的阈值，至此算法结束。
步骤S46：最后采用如式(15)对算法n次迭代计算的最终聚类中心cj(n)所属的类簇Sj(n)进行二值化处理，得到最终的变化检测结果图。
其中，S1(n)为变化类，S2(n)为未变化类，Img_change为变化检测结果图。
具体来说，本公开的发明构思之一，旨在通过上述非监督的地基微变监测雷达图像的区域变化检测方法、非监督的地基微变监测雷达图像的区域变化检装置及计算机可读存储介质，主要根据获取目标区域一段时间序列的图像；得到时间序列等间隔图像的相关性；基于相关性图像，获得初始差异图；对处理后的初始差异图进行聚类，获得变化检测结果，旨在基于获取的多幅图像取相同时间间隔的每两幅图像按时间顺序依次计算相关性，然后对时序相邻的每两幅相关性图像进行处理后的图像得到了非线性收缩，能够将乘性噪声转化为加性噪声，从而达到有效的去噪和增强变化和非变化区域的对比度，不同于传统的基于两幅雷达图像的变化检测，而是采用了一段时间序列的多幅图像获取差异图，克服了地基微变监测雷达噪声多、变化前后图像难以区分、差异图难获取的困难。进一步的，结合了地基微变监测雷达图像的滑坡区域在时间与空间上的特性改进了传统聚类算法，使得算法保证识别准确度的同时又减少了计算的迭代次数，简单高效。
本公开涉及各实施例的有益效果至少体现在：
1)通过分析地基微变监测雷达的滑坡区域相关性在时间序列上区别于未变化区域与噪声的变化特性，设计了分组时序的相关对数比值差异图的获取算法，有效避免了地基微变监测雷达图像受乘性散斑噪声影响导致的原始图像上变化与未变化区域难区分，且仅依靠两个时相的图像无法有效的完成变化检测的难题；
2)利用小波多分辨率分析分解差异图，对差异图进一步进行去噪处理，保留图像细节部分的同时，有效避免了噪声对变化检测结果的影响；
3)通过设计改进聚类算法的初始聚类中心的选取规则，充分结合了地基微变监测雷达图像差异图的空间信息和数值特性，迭代计算次数更少，具有高效，计算量小，聚类结果稳定性更高的特点。
本公开还提供了一种计算机可读存储介质，其上存储有计算机可执行指令，所述计算机可执行指令由处理器执行时，主要实现根据上述的非监督的地基微变监测雷达图像的区域变化检测方法，至少包括：
获取目标区域一段时间序列的图像；
得到时间序列等间隔图像的相关性；
基于相关性图像，获得初始差异图；
对处理后的初始差异图进行聚类，获得变化检测结果。
以上实施例仅为本公开的示例性实施例，不用于限制本公开，本公开的保护范围由权利要求书限定。本领域技术人员可以在本公开的实质和保护范围内，对本公开做出各种修改或等同替换，这种修改或等同替换也应视为落在本公开的保护范围内。
Claims (10)
1.非监督的地基微变监测雷达图像的区域变化检测方法，包括：
获取目标区域一段时间序列的图像；
得到时间序列等间隔图像的相关性；
基于相关性图像，获得初始差异图；
对处理后的初始差异图进行聚类，获得变化检测结果。
2.根据权利要求1所述的方法，其中，所述获取目标区域一段时间序列的图像，包括：
获取一组时间顺序的等间隔微变监测雷达图像。
3.根据权利要求2所述的方法，其中，所述得到时间序列等间隔图像的相关性，包括：
获取时间顺序的等间隔的每两幅图像的相关性图像。
4.根据权利要求3所述的方法，其中，所述基于相关性图像，获得初始差异图，包括：
基于获取的相关性图像取相邻的每两幅图像，得到相关对数比值图；
对相关对数比值图进行分组时序累加运算，以消除噪声。
5.根据权利要求4所述的方法，其中，初始差异图的处理，包括：
小波多尺度变换对差异图进行多层分解，将初始差异图分解为低频分量与高频分量，依次对上一层得到的低频分量分解为低频分量和高频分量，直到得到预设效果的分解层数；
对分解后的高频分量进行阈值去噪和系数重构，获取处理后的初始差异图。
6.根据权利要求5所述的方法，其中，所述进行聚类，包括：
取预设尺寸窗口遍历处理后的初始差异图的每一个像素点，并取中心像素点邻域内平均值；
对邻域内平均值进行排序，选出其中最大值与最小值，分别对邻域内所有的数据点排序，将数据值最大的数据点和最小的数据点分别作为类簇的初始聚类中心点；
输入无标签的数据，以距离作为相似度的衡量标准，根据不同数据之间的相似程度将数据聚类成不同的类簇；
将距离平方和作为准则函数，将每个数据点分配到聚类中心对应的簇中，重新计算出每一个簇中所有数据的平均值作为该簇新的聚类中心点，直至聚类中心点基本不变或者是准则函数收敛到设定的阈值。
7.根据权利要求6所述的方法，其中，所述获得变化检测结果，包括：
对迭代计算得到的聚类中心所属的类簇进行二值化处理，得到最终的变化检测结果图。
8.非监督的地基微变监测雷达图像的区域变化检测装置，包括：
获取模块，其配置为获取目标区域一段时间序列的图像；
变化检测模块，其配置为：
得到时间序列等间隔图像的相关性；
基于相关性图像，获得初始差异图；
对处理后的初始差异图进行聚类，获得变化检测结果。
9.根据权利要求8所述的装置，其中，所述变化检测模块，进一步配置为：
基于获取的相关性图像取相邻的每两幅图像，得到相关对数比值图；
对相关对数比值图进行分组时序累加运算，以消除噪声。
10.一种计算机可读存储介质，其上存储有计算机可执行指令，所述计算机可执行指令由处理器执行时，实现：
根据权利要求1至7所述的非监督的地基微变监测雷达图像的区域变化检测方法。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202310139606.2A CN116109600A (zh) | 2023-02-17 | 2023-02-17 | 非监督的地基微变监测雷达图像的区域变化检测方法、装置及存储介质 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202310139606.2A CN116109600A (zh) | 2023-02-17 | 2023-02-17 | 非监督的地基微变监测雷达图像的区域变化检测方法、装置及存储介质 |
Publications (1)
| Publication Number | Publication Date |
|---|---|
| CN116109600A true CN116109600A (zh) | 2023-05-12 |
Family
ID=86255976
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202310139606.2A Pending CN116109600A (zh) | 2023-02-17 | 2023-02-17 | 非监督的地基微变监测雷达图像的区域变化检测方法、装置及存储介质 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN116109600A (zh) |
Cited By (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN116701848A (zh) * | 2023-08-09 | 2023-09-05 | 江苏盖亚环境科技股份有限公司 | 一种一体化装备连续探测数据处理系统 |
-
2023
- 2023-02-17 CN CN202310139606.2A patent/CN116109600A/zh active Pending
Cited By (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN116701848A (zh) * | 2023-08-09 | 2023-09-05 | 江苏盖亚环境科技股份有限公司 | 一种一体化装备连续探测数据处理系统 |
| CN116701848B (zh) * | 2023-08-09 | 2023-12-08 | 江苏盖亚环境科技股份有限公司 | 一种一体化装备连续探测数据处理系统 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| Li et al. | Complex contourlet-CNN for polarimetric SAR image classification | |
| Celik | A Bayesian approach to unsupervised multiscale change detection in synthetic aperture radar images | |
| Jamali et al. | PolSAR image classification based on deep convolutional neural networks using wavelet transformation | |
| CN102800074A (zh) | 基于轮廓波变换的sar图像变化检测差异图生成方法 | |
| CN116109600A (zh) | 非监督的地基微变监测雷达图像的区域变化检测方法、装置及存储介质 | |
| Guobin et al. | Road identification algorithm for remote sensing images based on wavelet transform and recursive operator | |
| Manjunath et al. | Brain tumor detection and classification using convolution neural network | |
| Ma et al. | Multiscale 2-D singular spectrum analysis and principal component analysis for spatial–spectral noise-robust feature extraction and classification of hyperspectral images | |
| Chen et al. | Change detection algorithm for multi-temporal remote sensing images based on adaptive parameter estimation | |
| Li et al. | A method to improve the accuracy of SAR image change detection by using an image enhancement method | |
| Zhao et al. | Change detection method of high resolution remote sensing image based on DS evidence theory feature fusion | |
| Dong et al. | Deep shearlet network for change detection in sar images | |
| Kekre et al. | SAR image segmentation using vector quantization technique on entropy images | |
| Wu et al. | Multitemporal images change detection using nonsubsampled contourlet transform and kernel fuzzy C-means clustering | |
| CN114066815B (zh) | 基于模糊c均值聚类的SAR图像变化检测方法 | |
| Deng et al. | Two-stage visual attention model guided SAR image change detection | |
| Kakadiya et al. | Shape Extraction methods for fruits: technical review | |
| Micheal et al. | Automatic mountain detection in lunar images using texture of DTM data | |
| CN113326790A (zh) | 基于异常检测思维的胶囊机器人排水管病害检测方法 | |
| Xue et al. | A new method of SAR image segmentation based on FCM and wavelet transform | |
| CN117333468B (zh) | 面向多模式时序PolSAR影像的洪涝灾害监测方法 | |
| Thakur et al. | Generalized spatial kernel based fuzzy C-means clustering algorithm for image segmentation | |
| Zhang et al. | Landslide Area Identification and Detection Method Based on Micro-Variation Monitoring Radar Images | |
| Chen et al. | Lossless Compression Method of Color Image Based on Improved Particle Swarm Optimization | |
| Zhang et al. | Image stripe noise removal based on compressed sensing |
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 |