一种数学形态学滤波结合卡方变换的多时相遥感影像变化检测方法技术

本发明专利技术公开了一种数学形态学滤波结合卡方变换的多时相遥感影像变化检测方法，首先计算多时相多光谱影像的差异影像，在此基础上，计算影像上每一个点的CST值，根据致信水平获取阈值，得到初步的变化检测结果，然后再对该初步结果进行众数滤波(嵌入空间信息)，并根据滤波结果重新计算非变化区域的均值和方差矩阵。重复上述过程直到检测结果没有变化为止。其中检测过程中置信水平的选择是通过伪训练样本集来选择的。在最优的置信水平基础上，获取最终的变化检测结果。本发明专利技术在迭代估计卡方变换的统计参数过程中，加入数学形态学来对检测结果进行空间约束，通用性好，且提高了检测精度。

A multi-phase remote sensing image change detection method based on mathematical morphological filtering and chi square transform

The invention discloses a mathematical morphology filtering with chi square transform of multitemporal remote sensing image change detection method, first calculate the difference of multitemporal multispectral images of the image, based on the calculated image of every point of the CST value, according to the level of access to the letter to the threshold, the preliminary change detection results, and then the preliminary results of mode filtering (embedded spatial information), and calculate the non regional variation in mean and variance matrix according to the filtering results. Repeat the above process until the test results are not changed. The selection of confidence level in the detection process is selected by the pseudo training sample set. On the basis of the optimal confidence level, the final change detection results are obtained. In the process of iteratively estimating the statistical parameters of the chi square transformation, the present invention adds mathematical morphology to the spatial constraint of the test results, which has good universality and improves the detection accuracy.

【技术实现步骤摘要】
一种数学形态学滤波结合卡方变换的多时相遥感影像变化检测方法
本专利技术涉及一种数学形态学滤波结合卡方变换的多时相遥感影像变化检测方法，属于多时相光学遥感影像变化检测

技术介绍
随着多时相遥感数据的不断积累以及空间数据库的相继建立，如何从这些遥感数据中提取和检测变化信息已成为遥感科学和地理信息科学的重要研究课题。根据同一区域不同时相的遥感影像，可以提取城市、环境等动态变化的信息，为资源管理与规划、环境保护等部门提供科学决策的依据。遥感影像的变化检测就是从不同时期的遥感数据中，定量地分析和确定地表变化的特征与过程。各国学者从不同的角度和应用研究提出了许多有效的检测算法，如变化矢量分析法(ChangeVectorAnalysis，CVA)、基于FuzzyC-means(FCM)的聚类方法等。其中，传统的基于卡方变换(Chi-SquaredTransform,CST)的多时相光学遥感变化检测，先计算差异影像的均值和方差矩阵，然后再基于置信水平，确定变化检测的阈值，进而得到变化检测结果。该类技术中，常规的CST的不足是仅使用多时相高分辨率差异影像的光谱信息，没有利用空间信息。Shi等提出的基于空间约束的CST结合MRF模型(SCCSTMRF)的非监督变换检测方法通过众数(mode)滤波将空间信息嵌入到检测过程中，提高了检测的精度[AiyeShi,ChaoWang,ShaohongShen,FengchenHuang,andZhenliMa.Unsupervisedchangedetectionofmultispectralimagesbasedonspatialconstraintchi-squaredtransformandMarkovrandomfieldmodel.JournalofAppliedRemoteSensing,2016,10(4),046028:1-18.]。该空间约束的CST方法的不足是(1)容易引起较大的虚检率。(2)计算较复杂，对CST的结果须要采用MRF模型进一步精炼。
技术实现思路
针对现有技术存在的不足，本专利技术目的是提供一种数学形态学滤波结合卡方变换的多时相遥感影像变化检测方法，在迭代估计卡方变换的统计参数过程中，加入数学形态学来对检测结果进行空间约束，通用性好，且提高了检测精度。为了实现上述目的，本专利技术是通过如下的技术方案来实现：本专利技术的一种数学形态学滤波结合卡方变换的多时相遥感影像变化检测方法，包括以下几个步骤：步骤1：输入两时相的高分辨率光学遥感影像，分别记为：X1和X2；步骤2：利用遥感软件对X1和X2进行影像配准，分为几何粗校正和几何精校正；步骤3：利用多元变化检测方法对X1和X2进行辐射归一化校正；步骤4：对输入的多时相高分辨率光学遥感影像分别计算(是现有的计算方法，即简单的波段相减，此处不再赘述)每一个波段的差异影像，最后将每一个波段的差异影像进行组合(就是简单的组合，为现有技术。比如将独立的3个波段影像，组合为一个3波段的单独影像，此处不再赘述)，获取多时相差异影像Dx，DX＝X1-X2；步骤5：计算Dx的模值，记为XM；并利用Bayes原理，基于最大期望算法获取最优分割阈值T，设定Dx的模值动态范围，并依据动态范围的百分比确定伪训练样本集；步骤6：设定置信水平1-α的搜索范围和搜索步长；步骤7：按照迭代步长增加置信水平；步骤8：针对设定的置信水平，计算非变化区域的均值矢量m和方差矩阵∑，并计算卡方值；步骤9：对每一个置信水平计算相应的阈值，并根据阈值获取初步的检测结果M0；步骤10：设定数学形态学滤波结构单元B的形状和大小，对步骤9的检测结果M0进行数学形态学开操作滤波，记滤波的结果为M1：式中，表示图像的开操作；步骤11：利用步骤10的结果M1，确定影像的变化区域和非变化区域；步骤12：判读本次迭代的检测结果与前一次检测结果是否有改变，如果没有改变，回到步骤7；如果有改变，计算DX中相对应于非变化区域的均值矢量和方差矩阵，回到步骤8；步骤13：判读置信水平是否达到最大值；如果是，终止迭代；如果否，回到步骤7；步骤14：针对每一个置信水平，判断(判断方法为现有技术，精度的指标可以选为分类的正确率)伪训练样本集的精度，选择精度最高的置信水平，并输出对应的变化检测结果。步骤2中，对于所述几何粗校正，利用ENVI4.8软件实现，具体操作步骤为：步骤2.1：显示基准影像和待校正影像；步骤2.2：采集地面控制点GCPs，所述地面控制点GCPs应均匀分布在整幅图像内，地面控制点GCPs的数目至少大于等于9；步骤2.3：按照计算同名控制点的均方根误差，其中P1i，分别是前一时相和后一时相影像的同名控制点，S表示同名控制点的数目；步骤2.4：选择两时相影像的同名控制点间的匹配模型为二次多项式模型；步骤2.5：采用双线性插值法进行重采样输出。上述双线性插值法具体的方法如下：若求未知函数f在点P＝(x,y)的f(x，y)值，假设已知函数f在Q11＝(x1,y1),Q12＝(x1,y2),Q21＝(x2,y1)及Q22＝(x2,y2)四个点的值；如果选择一个坐标系统使得这四个点的坐标分别为(0,0)、(0,1)、(1,0)和(1,1)，那么双线性插值公式就可以表示为：f(x,y)≈f(0,0)(1-x)(1-y)+f(1,0)x(1-y)+f(0,1)(1-x)y+f(1,1)xy(1)。步骤2中，所述几何精校正具体的方法如下：将经过所述几何粗校正的多光谱遥感影像数据，利用三角剖分法进行几何精校正；所述三角剖分法为，采用逐点插入法构建Delaunay三角网，对每一个三角形，利用其三个顶点的行列号与其对应的基准影像同名点的地理坐标，来确定(具体的确定方法为现有技术，此处不再赘述)该三角形内部的仿射变换模型参数，对待校正影像进行纠正，得到校正后的遥感影。步骤3中，具体的方法如下：首先，找到两期影像各波段亮度值的一个线性组合，得到变化信息增强的差异影像，通过分割阈值确定变化区域和未变化区域；然后，通过未变化区域对应的两时相像元对的映射方程，完成相对辐射校正。步骤5中，将|XM-T|≤δ的区域作为伪训练样本集，其中δ的选择为XM动态范围的15％，最优分割阈值T的计算过程及伪训练样本集的构建如下：1)假设XM影像上未变化类ωn和变化类ωc服从如下的高斯分布，即：未变化类的均值和方差为mn和σn，变化类的均值和方差为mc和σc，X表示差异影像的模值，ωi表示类别，可以为非变化类ωn，也可以为变化类ωc；2)采用最大期望算法估计mn、σn、mc和σc这四个参数，其中，I和J分别表示影像的行数和列数，t表示迭代次数，pt+1(ωn)表示第t+1次迭代后，类别ωn的先验概率，表示非变化类在t次迭代后的均值，表示非变化类在t+1次迭代后的方差，公式中的p为均值概率；3)依据Bayes最小误差准则，求解变化矢量幅值图像XM的分割阈值T，在高斯分布的情况下，等同于求解下式：4)根据最大期望算法所估计的阈值T，伪训练样本集的构建分为如下两个部分：其中，未变化类伪训练集样本为：变化类伪训练样本集为：步骤6中，所述置信水平1-α的搜索范围为0.95-0.999，搜索步长为0.001。步骤8中，第一本文档来自技高网...

【技术保护点】
一种数学形态学滤波结合卡方变换的多时相遥感影像变化检测方法，其特征在于，包括以下几个步骤：步骤1：输入两时相的高分辨率光学遥感影像，分别记为：X1和X2；步骤2：利用遥感软件对X1和X2进行影像配准，分为几何粗校正和几何精校正；步骤3：利用多元变化检测方法对X1和X2进行辐射归一化校正；步骤4：对输入的多时相高分辨率光学遥感影像分别计算每一个波段的差异影像，最后将每一个波段的差异影像进行组合，获取多时相差异影像Dx，DX＝X1‑X2；步骤5：计算Dx的模值，记为XM；并利用Bayes原理，基于最大期望算法获取最优分割阈值T，设定Dx的模值动态范围，并依据动态范围的百分比确定伪训练样本集；步骤6：设定置信水平1‑α的搜索范围和搜索步长；步骤7：按照迭代步长增加置信水平；步骤8：针对设定的置信水平，计算非变化区域的均值矢量m和方差矩阵∑，并计算卡方值；步骤9：对每一个置信水平计算相应的阈值，并根据阈值获取初步的检测结果M0；步骤10：设定数学形态学滤波结构单元B的形状和大小，对步骤9的检测结果M0进行数学形态学开操作滤波，记滤波的结果为M1：

【技术特征摘要】
1.一种数学形态学滤波结合卡方变换的多时相遥感影像变化检测方法，其特征在于，包括以下几个步骤：步骤1：输入两时相的高分辨率光学遥感影像，分别记为：X1和X2；步骤2：利用遥感软件对X1和X2进行影像配准，分为几何粗校正和几何精校正；步骤3：利用多元变化检测方法对X1和X2进行辐射归一化校正；步骤4：对输入的多时相高分辨率光学遥感影像分别计算每一个波段的差异影像，最后将每一个波段的差异影像进行组合，获取多时相差异影像Dx，DX＝X1-X2；步骤5：计算Dx的模值，记为XM；并利用Bayes原理，基于最大期望算法获取最优分割阈值T，设定Dx的模值动态范围，并依据动态范围的百分比确定伪训练样本集；步骤6：设定置信水平1-α的搜索范围和搜索步长；步骤7：按照迭代步长增加置信水平；步骤8：针对设定的置信水平，计算非变化区域的均值矢量m和方差矩阵∑，并计算卡方值；步骤9：对每一个置信水平计算相应的阈值，并根据阈值获取初步的检测结果M0；步骤10：设定数学形态学滤波结构单元B的形状和大小，对步骤9的检测结果M0进行数学形态学开操作滤波，记滤波的结果为M1：式中，表示图像的开操作；步骤11：利用步骤10的结果M1，确定影像的变化区域和非变化区域；步骤12：判读本次迭代的检测结果与前一次检测结果是否有改变，如果没有改变，回到步骤7；如果有改变，计算DX中相对应于非变化区域的均值矢量和方差矩阵，回到步骤8；步骤13：判读置信水平是否达到最大值；如果是，终止迭代；如果否，回到步骤7；步骤14：针对每一个置信水平，判断伪训练样本集的精度，选择精度最高的置信水平，并输出对应的变化检测结果。2.根据权利要求1所述的数学形态学滤波结合卡方变换的多时相遥感影像变化检测方法，其特征在于，步骤2中，对于所述几何粗校正，利用ENVI4.8软件实现，具体操作步骤为：步骤2.1：显示基准影像和待校正影像；步骤2.2：采集地面控制点GCPs，所述地面控制点GCPs应均匀分布在整幅图像内，地面控制点GCPs的数目至少大于等于9；步骤2.3：按照计算同名控制点的均方根误差，其中P1i，分别是前一时相和后一时相影像的同名控制点，S表示同名控制点的数目；步骤2.4：选择两时相影像的同名控制点间的匹配模型为二次多项式模型；步骤2.5：采用双线性插值法进行重采样输出。3.根据权利要求2所述的数学形态学滤波结合卡方变换的多时相遥感影像变化检测方法，其特征在于，所述双线性插值法具体的方法如下：若求未知函数f在点P＝(x,y)的f(x，y)值，假设已知函数f在Q11＝(x1,y1),Q12＝(x1,y2),Q21＝(x2,y1)及Q22＝(x2,y2)四个点的值；如果选择一个坐标系统使得这四个点的坐标分别为(0,0)、(0,1)、(1,0)和(1,1)，那么双线性插值公式就可以表示为：f(x,y)≈f(0,0)(1-x)(1-y)+f(1,0)x(1-y)+f(0,1)(1-x)y+f(1,1)xy(1)。4.根据权利要求1所述的数学形态学滤波结合卡方变换的多时相遥感影像变化检测方法，其特征在于，步骤2中，所述几何精校正具体的方法如下：将经过所述几何粗校正的多光谱遥感影像数据，利用三角剖分法进行几何精校正；所述三角剖分法为，采用逐点插入法构建Delaunay三角网，对每一个三角形，利用其三个顶点的行列号与其对应的基准影像同名点的地理坐标，来确定该三角形内部的仿射变换模型参数，对待校正影像进行纠正，得到校正后的遥感影。5.根据权利要求1所述的数学形态学滤波结合卡方变换的多时相遥感影像变化检测方法，其特征在于，步骤3中，具体的方法如下：首先，找到两期影像各波段亮度值的一个线性组合，得到变化信息增强的差异影像，通过分割阈值确定变化区域和未变化区域；然后，通过未变化区域对应的两时相像元对的映射方程，完成相对辐射校正。6.根据权利要求1所述的数学形态学滤波结合卡方变换的多时相遥感影像变化检测方法，其特征在于，步骤5中，将|XM-T|≤δ的区域作为伪训练样本集，其中δ的选择为XM动态范围的15％，最优分割阈值T的计算过程及伪训练样本集的构建如下：1)假设XM影像上未变化类ωn和变化...

【专利技术属性】
技术研发人员：石爱业，许磊，
申请(专利权)人：河海大学，
类型：发明
国别省市：江苏,32