高分三号SAR影像双阈值变化检测

高分三号SAR影像双阈值变化检测 Dual-thresholds change detection in GF-3 SAR images 1 2 3 first-author 崔斌, 1989年生，男, 博士研究生, 研究方向为SAR图像处理，变化检测。E-mail: cuibinwing@126.com 崔斌, 1989年生，男, 博士研究生, 研究方向为SAR图像处理，变化检测。E-mail: cuibinwing@126.com cuibinwing@126.com 2 corresp 张永红，1973年生，男，研究员，研究方向为SAR图像处理，土地利用与覆盖变化。E-mail: yhzhang@casm.ac.cn 张永红，1973年生，男，研究员，研究方向为SAR图像处理，土地利用与覆盖变化。E-mail: yhzhang@casm.ac.cn yhzhang@casm.ac.cn 1 2 武汉大学 测绘学院, 武汉 430079 武汉大学 测绘学院, 武汉 430079 School of Geodesy and Geomatics, Wuhan University, Wuhan 430079, China School of Geodesy and Geomatics, Wuhan University, Wuhan 430079, China 中国测绘科学研究院, 北京 100830 中国测绘科学研究院, 北京 100830 Chinese Academy of Surveying and Mapping, Beijing 100830, China Chinese Academy of Surveying and Mapping, Beijing 100830, China 城市空间信息工程北京市重点实验室, 北京 100038 城市空间信息工程北京市重点实验室, 北京 100038 Beijing Key Laboratory of Urban Spatial Information Engineering, Beijing 100038, China Beijing Key Laboratory of Urban Spatial Information Engineering, Beijing 100038, China 双阈值合成孔径雷达SAR(Synthetic Aperture Radar)变化检测算法具有在发现变化区域的同时还能确定地表发生后向散射变化类型的优点。针对广义高斯双阈值最小误差法D-GKIT(Dual Generalized Kittler and Illingworth Thresholding)在进行阈值选取时直方图中不同类别像素灰度级重叠严重时，分割结果容易在尖峰单侧选取出双阈值而导致无法正确分割差异图的问题，本文提出一种结合归一化最大类间方差和广义高斯最小误差法GKIT(Generalized Kittler and Illingworth Thresholding)的双阈值SAR变化检测方法。首先，提出以归一化最大类间方差值作为灰度级重叠程度的判别参数，确定阈值的选取顺序及两个候选区间；然后，利用GKIT在候选区间内进行分割，获取单侧阈值及非变化类拟合函数；最后，提出利用非变化类拟合函数更新后的直方图作为另一侧阈值选取基础进行分割，得到对应分割阈值。以宁波地区高分三号(GF-3)SAR卫星影像作为试验研究数据，结果表明：本文方法能较好地解决灰度级重叠时D-GKIT无法进行正确分割的问题，具有良好的变化检测效果和更强的鲁棒性且达到了利用研究区数据验证利用GF-3号SAR卫星影像进行变化检测研究可行性的目的。 Compared with the single threshold segment method in SAR change detection, the dual thresholds segment method can simultaneously identify the change areas and confirm the change types. Although D-GKIT shows a superior performance, a strongly overlapping gray level is observed in the histogram of the difference image, thereby inaccurately identifying double thresholds on the same side of the peak. In this paper, we apply a dual-thresholds method combined with normalized maximal between-class variance and GKIT test on GF-3 images to verify its feasibility of our method and the ability of change detection ability.First, the normalized maximal between-class variance values of two sides surrounding the peak in the histogram are taken as the degrees of the overlapping gray level, and then the thresholds selection sequence and the candidate intervals are confirmed. Second, the side at which the gray level lightly overlaps is segmented by GKIT, and the threshold and the fitting function of the unchanged class are obtained. Third, the fitting function of the unchanged class is used to replace the corresponding part in the origin histogram to form a new histogram that is subsequently segmented to obtain the threshold in the second candidate interval. Finally, the two thresholds are applied on the difference image to obtain the final change result. The experiment on GF-3 SAR images reveals that the performance on our proposed method outperforms D-GKIT and can deal well with the overlapping gray level overlapped in the histogram of the difference image. The confusion matrix of the results for various local areas in the change image also shows that the proposed method has been slightly influenced by the overlapping gray level overlapped and obtains generally good results. Therefore, the feasibility of our method and its change detection ability by using GF-3 images are verified. We propose a method based on the normalized maximal between-class variance and GKIT to segment a difference image by applying dual thresholds in SAR change detection. The effectiveness of the proposed method and its change detection ability by using GF-3 images are validated by the experiment results. 遥感 SAR变化检测 双阈值分割 高分三号(GF-3) GKIT 归一化最大类间方差 remote sensing SAR change detection dual thresholds segment GF-3 GKIT normalized maximal between-class variance 1　引　言 1 遥感影像变化检测与传统方法进行地表覆盖变化检测相比，具有方便高效，节省外业工作等优势。随着SAR影像质量和数量不断提升，以及自身不受天气和光照影响，主动成像的特性( 郝洪美 等，2012 郝洪美 等，2012 )，使SAR变化检测在城市扩张研究( Hu和Ban，2014 Hu和Ban，2014 ； Yousif和Ban等，2014 Yousif和Ban等，2014 )，舰船探测( 张小强 等，2015 张小强 等，2015 )，森林监测( Hame等，1998 Hame等，1998 ； Pantze等，2014 Pantze等，2014 )和灾害检测( Whittle等，2012 Whittle等，2012 )等方面有着重要的利用价值。 SAR变化检测重点和难点( Xiong等，2012 Xiong等，2012 )主要集中在差异图生成和差异图分割两方面，围绕这两方面国内外学者做了大量研究。对于差异图生成，由于对数比方法LR (Log-Ratio)可以将乘性噪声转为加性噪声，且对数运算可以增强非变化类与变化类的对比度( 公茂果 等，2016 公茂果 等，2016 )，使其成为使用最广泛的差异图生成方法。此外还有结合纹理特征的差异图算法( Gong 等，2014 Gong 等，2014 )， Xiong等(2012) Xiong等(2012) 提出的邻域似然比算法，利用小波融合构造差异图等算法( Ma等，2012 Ma等，2012 )。差异图分割方面，直观高效的阈值法使用较为广泛，最初的研究主要集中在单阈值分割，包括传统的OTSU法( Otsu，1979 Otsu，1979 )，最小误差法KI(Kittler and Illingworth)(Kittler 和Illingworth，1986)，最大熵法等，近些年一些结合统计分布信息的阈值分割算法也取得了较好的效果，比如针对广义Gamma模型的自适应KI阈值分割方法( 高丛珊 等，2010 高丛珊 等，2010 )，针对广义高斯GG (Generalized Gaussian)模型的GKIT算法( Bazi 等，2005 Bazi 等，2005 )以及在此基础上增加了模型种类的结合KI准则的算法( Moser和Serpico，2006 Moser和Serpico，2006 )等。 为了在SAR变化检测时能同时确定出变化区域及变化类型，为确定地表覆盖改变情况打下基础，不少学者开始研究利用双阈值算法对差异图进行分割。主要有在假定不同类别像素服从高斯分布模型的基础上，采用期望最大化算法EM(Expectation Maximization)进行双阈值分割( 黄世奇等，2010 黄世奇等，2010 )。 Bazi等(2006) Bazi等(2006) 进一步完善的结合Hessian矩阵检验阈值个数的GKIT算法先根据Hessian矩阵性质对阈值个数进行判断，单阈值利用GKIT算法，双阈值则利用广义高斯双阈值最小误差法D-GKIT(Dual Generalized Kittler and Illingworth Thresholding)进行双阈值分割。 胡召玲(2013) 胡召玲(2013) 将D-GKIT算法进行了独立使用，用来确定在明显存在两种变化像素情况下的SAR变化检测双阈值分割，在假设3种类别像素服从GG分布的基础上，结合KI准则，确定准则函数最小时对应的双阈值，将GKIT推广为使用范围更广泛的双阈值分割算法。 Hu和Ban (2014) Hu和Ban (2014) 针对LR差异图变化类和非变化类像素数量相差过大导致直方图中只存在一个明显尖峰致使直方图分割困难的问题，提出在阈值分割前先判断直方图形态的方法：若检验结果直方图中存在两个波峰，则利用GKIT分割；若检测结果呈现单峰特性，则人工选取差异图中一些区域作为样本生成新直方图来确定分割阈值，方法虽然得到了可能正确的阈值，但其过程复杂且需人工干预，结果受人工选取区域的影响。 D-GKIT是使用较为广泛的SAR变化检测双阈值分割算法，但由于其变化类和非变化类的统计参数都是有偏估计，当不同类别像素灰度级重叠严重时，即发生差异图直方图中一侧变化类C1的灰度级拟合为GG分布的准则函数值大于将另一侧变化类C2灰度级拟合为两个GG分布的准则函数值的情况时，其盲目的从全局选取最小误差将导致双阈值均选取在尖峰同一侧的错误情况。特别是当差异图中两种不同类型的变化像素数量相差较大，差异图整体差异度较小时更容易出现这一错误情况。因此，本文提出一种结合归一化最大类间方差和GKIT的双阈值分割算法。首先对一维高斯核平滑后的差异图直方图进行两侧归一化最大类间方差计算，判断两侧灰度级重叠程度，确定阈值选取顺序和候选区间；之后利用GKIT对重叠程度较小一侧在候选区间内进行阈值选取，得到对应单侧阈值和非变化类GG分布拟合函数；最后利用非变化类GG分布拟合函数更新原始直方图以消除一侧变化像素灰度级对其他灰度级拟合选取阈值时的影响，得到另一侧对应阈值。通过利用国产GF-3号浙江宁波地区两时相SAR影像数据进行试验研究，在验证了本文算法有效性的同时也证明了利用GF-3号SAR卫星影像进行变化检测研究的可行性。 2　算法介绍 1 本文算法的整体流程如 图1 图1 所示，大致分为影像预处理，差异图生成和阈值分割3部分。 2.1　差异图生成 2 采用窗口大小为3的增强LEE滤波对多视后的影像进行滤波。因LR差异图近似满足GG分布，能直观反应两时相间增强与减弱变化，且具有将图像固有的乘性噪声转化为加性噪声的优点。故在本文中采用LR算法构造差异图 <math display="block" id="M1"><mi>L</mi><mo>=</mo><mi mathvariant="normal">l</mi><mi mathvariant="normal">o</mi><mi mathvariant="normal">g</mi><mfrac><mrow><msub><mrow><mi>E</mi></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msub></mrow><mrow><msub><mrow><mi>E</mi></mrow><mrow><mn mathvariant="normal">1</mn></mrow></msub></mrow></mfrac></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165217&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165214&type= 15.57866669 9.05933285 式中， <math id="M2"><msub><mrow><mi>E</mi></mrow><mrow><mn mathvariant="normal">1</mn></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165223&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165220&type= 3.47133350 3.72533321 和 <math id="M3"><msub><mrow><mi>E</mi></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165230&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165226&type= 3.47133350 3.72533321 为滤波后两时相影像， <math id="M4"><mi>L</mi></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165238&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165234&type= 1.86266661 2.62466669 为LR差异图。 2.2　阈值选取 2 GKIT对服从GG分布且只含单一变化类型的差异图有良好分割效果。D-GKIT虽然可以进行双阈值分割，但分割结果受灰度级重叠影响。通过判断直方图两侧归一化最大类间方差值来确定变化像素灰度级与非变化像素灰度级间的重叠程度，归一化最大类间方差值大代表这一侧受灰度级重叠影响小，优先选取重叠度较小一侧进行分割可将拟合时一侧变化类灰度级对非变化类灰度级产生的影响降低。 2.2.1　直方图平滑 3 原始直方图局部频率震荡会引起出现局部波峰和波谷，影响两侧归一化类间方差计算。为在不影响直方图基本形态前提下，减小直方图局部震荡对类间方差的影响，本文先对直方图进行平滑处理。 一维高斯核函数平滑能有效减少灰度级频率不稳定带来的局部峰值和谷值( Tsai，1995 Tsai，1995 )，减小对类间方差计算的影响。1维高斯核函数为 <math display="block" id="M5"><mi>F</mi><mo stretchy="false">(</mo><mi>k</mi><mo>,</mo><mi>σ</mi><mo stretchy="false">)</mo><mo>=</mo><mfrac><mrow><mn mathvariant="normal">1</mn></mrow><mrow><mi>σ</mi><mroot><mrow><mn mathvariant="normal">2</mn><mi mathvariant="normal">π</mi></mrow><mrow/></mroot></mrow></mfrac><msup><mrow><mi mathvariant="normal">e</mi></mrow><mrow><mo>-</mo><mfrac><mrow><msup><mrow><mi>k</mi></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msup></mrow><mrow><mn mathvariant="normal">2</mn><msup><mrow><mi>σ</mi></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msup></mrow></mfrac></mrow></msup></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165245&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165242&type= 33.35866928 9.90599918 式中， <math id="M6"><mi>k</mi></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165253&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165248&type= 1.69333339 2.62466669 为灰度级， <math id="M7"><mi>σ</mi></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165261&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165258&type= 2.45533323 2.62466669 为窗口宽度。窗口大小为3时， F (-1)=0.2661, F (0)=0.5478, F (1)=0.2661( Tsai，1995 Tsai，1995 )，1维高斯核卷积如 式(3) 式(3) ， <math id="M8"><mo>⊗</mo></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165275&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165266&type= 3.72533321 2.87866688 表示卷积 <math display="block" id="M9"><mi>H</mi><mo stretchy="false">(</mo><mi>k</mi><mo>,</mo><mi>σ</mi><mo stretchy="false">)</mo><mo>=</mo><msub><mrow><mi>H</mi></mrow><mrow><mn mathvariant="normal">0</mn></mrow></msub><mo stretchy="false">(</mo><mi>k</mi><mo stretchy="false">)</mo><mo>⊗</mo><mi>F</mi><mo stretchy="false">(</mo><mi>k</mi><mo>,</mo><mi>σ</mi><mo stretchy="false">)</mo></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165284&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165279&type= 40.13199997 4.65666676 <math id="M10"><msub><mrow><mi>H</mi></mrow><mrow><mn mathvariant="normal">0</mn></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165292&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165288&type= 3.64066648 3.72533321 为原始直方图，离散情况有平滑式(4) <math display="block" id="M11"><mi>H</mi><mo stretchy="false">(</mo><mi>k</mi><mo>,</mo><mi>w</mi><mo stretchy="false">)</mo><mo>=</mo><munderover><mo>∑</mo><mrow><mi>s</mi><mo>=</mo><mo>-</mo><mfenced open="[" close="]" separators="|"><mrow><mfrac><mrow><mi>w</mi></mrow><mrow><mn mathvariant="normal">2</mn></mrow></mfrac></mrow></mfenced></mrow><mrow><mi>s</mi><mo>=</mo><mfenced open="[" close="]" separators="|"><mrow><mfrac><mrow><mi>w</mi></mrow><mrow><mn mathvariant="normal">2</mn></mrow></mfrac></mrow></mfenced></mrow></munderover><mrow><msub><mrow><mi>H</mi></mrow><mrow><mn mathvariant="normal">0</mn></mrow></msub></mrow><mo stretchy="false">(</mo><mi>k</mi><mo>+</mo><mi>s</mi><mo stretchy="false">)</mo><mo>⊗</mo><mi>F</mi><mo stretchy="false">(</mo><mi>s</mi><mo stretchy="false">)</mo></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165303&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165296&type= 49.95333481 14.73199940 式中，[·]表示取整, <math id="M12"><mi>w</mi></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165313&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165307&type= 2.28600001 2.62466669 为平滑窗口大小。 <math id="M13"><mi>H</mi></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165322&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165318&type= 2.53999996 2.62466669 为平滑后直方图，作为归一化最大类间方差的计算基础。 2.2.2　归一化最大类间方差 3 最大类间方差最早由 Otsu(1979) Otsu(1979) 提出，最初用以影像分割。随后 Demirkaya和Asyali(2004) Demirkaya和Asyali(2004) 将类间方差概念用于直方图形态检验，在分别给定类内方差和类间方差定义后，又进一步给定归一化类间方差的计算形式，以反应不同类别像素灰度级间的差异程度。基于此，本文用来判断直方图两侧不同类别像素灰度级间的重叠程度，归一化最大类间方差越大则灰度级重叠度越小。假设平滑后直方图中最大值出现在非变化类的中间位置(真实LR差异值0附近)，此时直方图两侧均包含一部分非变化类灰度级和一种变化类的灰度级。以直方图尖峰灰度级将直方图分为两部分，分别计算两侧归一化最大类间方差，值较大一侧变化类与非变化类灰度级重叠程度较小，即为强侧，优先计算阈值，另一侧为弱侧。 类内方差 <math id="M14"><msubsup><mrow><mi>σ</mi></mrow><mrow><mi>w</mi></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msubsup><mo stretchy="false">(</mo><mi>k</mi><mo stretchy="false">)</mo></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165328&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165325&type= 8.80533314 4.82599974 计算公式为 <math display="block" id="M15"><msubsup><mrow><mi>σ</mi></mrow><mrow><mi>w</mi></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msubsup><mo stretchy="false">(</mo><mi>k</mi><mo stretchy="false">)</mo><mo>=</mo><msub><mrow><mi>p</mi></mrow><mrow><mn mathvariant="normal">1</mn></mrow></msub><mo stretchy="false">(</mo><mi>k</mi><mo stretchy="false">)</mo><msubsup><mrow><mi>σ</mi></mrow><mrow><mn mathvariant="normal">1</mn></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msubsup><mo stretchy="false">(</mo><mi>k</mi><mo stretchy="false">)</mo><mo>+</mo><msub><mrow><mi>p</mi></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msub><mo stretchy="false">(</mo><mi>k</mi><mo stretchy="false">)</mo><msubsup><mrow><mi>σ</mi></mrow><mrow><mn mathvariant="normal">2</mn></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msubsup><mo stretchy="false">(</mo><mi>k</mi><mo stretchy="false">)</mo></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165335&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165332&type= 50.96933365 4.82599974 类间方差 <math id="M16"><msubsup><mrow><mi>σ</mi></mrow><mrow><mi>b</mi></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msubsup><mo stretchy="false">(</mo><mi>k</mi><mo stretchy="false">)</mo></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165344&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165340&type= 8.63599968 4.82599974 计算公式为 <math display="block" id="M17"><msubsup><mrow><mi>σ</mi></mrow><mrow><mi>b</mi></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msubsup><mo stretchy="false">(</mo><mi>k</mi><mo stretchy="false">)</mo><mo>=</mo><msub><mrow><mi>p</mi></mrow><mrow><mn mathvariant="normal">1</mn></mrow></msub><mo stretchy="false">(</mo><mi>k</mi><mo stretchy="false">)</mo><msub><mrow><mi>p</mi></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msub><mo stretchy="false">(</mo><mi>k</mi><mo stretchy="false">)</mo><mo stretchy="false">(</mo><msub><mrow><mi>m</mi></mrow><mrow><mn mathvariant="normal">1</mn></mrow></msub><mo stretchy="false">(</mo><mi>k</mi><mo stretchy="false">)</mo><mo>-</mo><msub><mrow><mi>m</mi></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msub><msup><mrow><mo stretchy="false">(</mo><mi>k</mi><mo stretchy="false">)</mo><mo stretchy="false">)</mo></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msup></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165350&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165347&type= 55.71066666 4.82599974 归一化类间方差 <math id="M18"><mi>G</mi><mo stretchy="false">(</mo><mi>k</mi><mo stretchy="false">)</mo></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165358&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165354&type= 7.19666672 3.47133350 计算方式为 <math display="block" id="M19"><mi>G</mi><mo stretchy="false">(</mo><mi>k</mi><mo stretchy="false">)</mo><mo>=</mo><mfrac><mrow><msubsup><mrow><mi>σ</mi></mrow><mrow><mi>b</mi></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msubsup><mo stretchy="false">(</mo><mi>k</mi><mo stretchy="false">)</mo></mrow><mrow><msubsup><mrow><mi>σ</mi></mrow><mrow><mi>w</mi></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msubsup><mo stretchy="false">(</mo><mi>k</mi><mo stretchy="false">)</mo><mo>+</mo><msubsup><mrow><mi>σ</mi></mrow><mrow><mi>b</mi></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msubsup><mo stretchy="false">(</mo><mi>k</mi><mo stretchy="false">)</mo></mrow></mfrac></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165367&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165361&type= 33.95133591 8.97466660 针对灰度级为 <math id="M20"><mi>N</mi></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165375&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165370&type= 2.53999996 2.62466669 ，分割灰度级为 <math id="M21"><mi>t</mi></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165382&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165379&type= 1.10066664 2.62466669 的两类别灰度级的频率和 <math id="M22"><msub><mrow><mi>p</mi></mrow><mrow><mi>i</mi></mrow></msub><mo stretchy="false">(</mo><mi>t</mi><mo stretchy="false">)</mo><mo>,</mo><mi>i</mi><mo>=</mo><mo stretchy="false">(</mo><mn mathvariant="normal">1,2</mn><mo stretchy="false">)</mo></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165388&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165385&type= 19.98133469 4.57200003 ，均值 <math id="M23"><msub><mrow><mi>m</mi></mrow><mrow><mi>i</mi></mrow></msub><mo stretchy="false">(</mo><mi>t</mi><mo stretchy="false">)</mo><mo>,</mo><mi>i</mi><mo>=</mo><mo stretchy="false">(</mo><mn mathvariant="normal">1,2</mn><mo stretchy="false">)</mo></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165394&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165391&type= 20.99733353 4.57200003 ,和方差 <math id="M24"><msubsup><mrow><mi>σ</mi></mrow><mrow><mi>i</mi></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msubsup><mo stretchy="false">(</mo><mi>t</mi><mo stretchy="false">)</mo><mo>,</mo><mi>i</mi><mo>=</mo><mo stretchy="false">(</mo><mn mathvariant="normal">1,2</mn><mo stretchy="false">)</mo></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165400&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165397&type= 21.16666603 4.82599974 ，如式 (8) (8) , (10) (10) 所示 <math display="block" id="M25"><msub><mrow><mi>p</mi></mrow><mrow><mn mathvariant="normal">1</mn></mrow></msub><mo stretchy="false">(</mo><mi>t</mi><mo stretchy="false">)</mo><mo>=</mo><munderover><mo>∑</mo><mrow><mi>k</mi><mo>=</mo><mn mathvariant="normal">1</mn></mrow><mrow><mi>t</mi></mrow></munderover><mrow><mi>H</mi><mo stretchy="false">(</mo><mi>k</mi><mo stretchy="false">)</mo></mrow><mo>,</mo><msub><mrow><mi>p</mi></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msub><mo stretchy="false">(</mo><mi>t</mi><mo stretchy="false">)</mo><mo>=</mo><munderover><mo>∑</mo><mrow><mi>k</mi><mo>=</mo><mi>t</mi><mo>+</mo><mn mathvariant="normal">1</mn></mrow><mrow><mi>N</mi></mrow></munderover><mrow><mi>H</mi><mo stretchy="false">(</mo><mi>k</mi><mo stretchy="false">)</mo></mrow></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165407&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165403&type= 51.47733307 8.29733276 <math display="block" id="M26"><mtable><mtr><mtd><msub><mrow><mi>m</mi></mrow><mrow><mn mathvariant="normal">1</mn></mrow></msub><mo stretchy="false">(</mo><mi>t</mi><mo stretchy="false">)</mo><mo>=</mo><munderover><mo>∑</mo><mrow><mi>k</mi><mo>=</mo><mn mathvariant="normal">1</mn></mrow><mrow><mi>t</mi></mrow></munderover><mrow><mrow><mrow><mi>k</mi><mo>×</mo><mi>H</mi><mo stretchy="false">(</mo><mi>k</mi><mo stretchy="false">)</mo></mrow><mo>/</mo><mrow><msub><mrow><mi>p</mi></mrow><mrow><mn mathvariant="normal">1</mn></mrow></msub></mrow></mrow></mrow><mo stretchy="false">(</mo><mi>t</mi><mo stretchy="false">)</mo></mtd></mtr><mtr><mtd><msub><mrow><mi>m</mi></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msub><mo stretchy="false">(</mo><mi>t</mi><mo stretchy="false">)</mo><mo>=</mo><munderover><mo>∑</mo><mrow><mi>k</mi><mo>=</mo><mi>t</mi><mo>+</mo><mn mathvariant="normal">1</mn></mrow><mrow><mi>N</mi></mrow></munderover><mrow><mrow><mrow><mi>k</mi><mo>×</mo><mi>H</mi><mo stretchy="false">(</mo><mi>k</mi><mo stretchy="false">)</mo></mrow><mo>/</mo><mrow><msub><mrow><mi>p</mi></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msub></mrow></mrow></mrow><mo stretchy="false">(</mo><mi>t</mi><mo stretchy="false">)</mo></mtd></mtr></mtable></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165413&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165409&type= 44.11133194 18.37266731 <math display="block" id="M27"><mtable><mtr><mtd><msubsup><mrow><mi>σ</mi></mrow><mrow><mn mathvariant="normal">1</mn></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msubsup><mo stretchy="false">(</mo><mi>t</mi><mo stretchy="false">)</mo><mo>=</mo><munderover><mo>∑</mo><mrow><mi>k</mi><mo>=</mo><mn mathvariant="normal">1</mn></mrow><mrow><mi>t</mi></mrow></munderover><mrow><mrow><mrow><mo stretchy="false">(</mo><mi>k</mi><mo>-</mo><msub><mrow><mi>m</mi></mrow><mrow><mn mathvariant="normal">1</mn></mrow></msub><msup><mrow><mo stretchy="false">(</mo><mi>t</mi><mo stretchy="false">)</mo><mo stretchy="false">)</mo></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msup><mo>×</mo><mi>H</mi><mo stretchy="false">(</mo><mi>k</mi><mo stretchy="false">)</mo></mrow><mo>/</mo><mrow><msub><mrow><mi>p</mi></mrow><mrow><mn mathvariant="normal">1</mn></mrow></msub></mrow></mrow></mrow><mo stretchy="false">(</mo><mi>t</mi><mo stretchy="false">)</mo></mtd></mtr><mtr><mtd><msubsup><mrow><mi>σ</mi></mrow><mrow><mn mathvariant="normal">2</mn></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msubsup><mo stretchy="false">(</mo><mi>t</mi><mo stretchy="false">)</mo><mo>=</mo><munderover><mo>∑</mo><mrow><mi>k</mi><mo>=</mo><mi>t</mi><mo>+</mo><mn mathvariant="normal">1</mn></mrow><mrow><mi>N</mi></mrow></munderover><mrow><mrow><mrow><mo stretchy="false">(</mo><mi>k</mi><mo>-</mo><msub><mrow><mi>m</mi></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msub><msup><mrow><mo stretchy="false">(</mo><mi>t</mi><mo stretchy="false">)</mo><mo stretchy="false">)</mo></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msup><mo>×</mo><mi>H</mi><mo stretchy="false">(</mo><mi>k</mi><mo stretchy="false">)</mo></mrow><mo>/</mo><mrow><msub><mrow><mi>p</mi></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msub></mrow></mrow></mrow><mo stretchy="false">(</mo><mi>t</mi><mo stretchy="false">)</mo></mtd></mtr></mtable></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165425&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165417&type= 60.79066849 18.37266731 计算最大值两侧的归一化类间方差为 <math display="block" id="M28"><mtable><mtr><mtd><msub><mrow><mi>G</mi></mrow><mrow><mi mathvariant="normal">l</mi><mi mathvariant="normal">e</mi><mi mathvariant="normal">f</mi><mi mathvariant="normal">t</mi></mrow></msub><mo stretchy="false">(</mo><msub><mrow><mi>k</mi></mrow><mrow><mn mathvariant="normal">1</mn></mrow></msub><mo stretchy="false">)</mo><mo>=</mo><mi mathvariant="normal">m</mi><mi mathvariant="normal">a</mi><mi mathvariant="normal">x</mi><mo stretchy="false">(</mo><mi>G</mi><mo stretchy="false">(</mo><msub><mrow><mi>k</mi></mrow><mrow><mn mathvariant="normal">1</mn></mrow></msub><mo stretchy="false">)</mo><mo stretchy="false">)</mo><mo>,</mo><mo stretchy="false">(</mo><mn mathvariant="normal">1</mn><mo>≤</mo><msub><mrow><mi>k</mi></mrow><mrow><mn mathvariant="normal">1</mn></mrow></msub><mo>≤</mo><msub><mrow><mi>h</mi></mrow><mrow><mi mathvariant="normal">m</mi><mi mathvariant="normal">a</mi><mi mathvariant="normal">x</mi></mrow></msub><mo stretchy="false">)</mo></mtd></mtr><mtr><mtd><msub><mrow><mi>G</mi></mrow><mrow><mi mathvariant="normal">r</mi><mi mathvariant="normal">i</mi><mi mathvariant="normal">g</mi><mi mathvariant="normal">h</mi><mi mathvariant="normal">t</mi></mrow></msub><mo stretchy="false">(</mo><msub><mrow><mi>k</mi></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msub><mo stretchy="false">)</mo><mo>=</mo><mi mathvariant="normal">m</mi><mi mathvariant="normal">a</mi><mi mathvariant="normal">x</mi><mo stretchy="false">(</mo><mi>G</mi><mo stretchy="false">(</mo><msub><mrow><mi>k</mi></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msub><mo stretchy="false">)</mo><mo stretchy="false">)</mo><mo>,</mo><mo stretchy="false">(</mo><msub><mrow><mi>h</mi></mrow><mrow><mi mathvariant="normal">m</mi><mi mathvariant="normal">a</mi><mi mathvariant="normal">x</mi></mrow></msub><mo>≤</mo><msub><mrow><mi>k</mi></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msub><mo>≤</mo><mi>N</mi><mo stretchy="false">)</mo></mtd></mtr></mtable></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165430&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165427&type= 57.91199875 9.82133293 式中， <math id="M29"><msub><mrow><mi>G</mi></mrow><mrow><mi mathvariant="normal">l</mi><mi mathvariant="normal">e</mi><mi mathvariant="normal">f</mi><mi mathvariant="normal">t</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165435&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165434&type= 4.99533367 3.72533321 和 <math id="M30"><msub><mrow><mi>G</mi></mrow><mrow><mi mathvariant="normal">r</mi><mi mathvariant="normal">i</mi><mi mathvariant="normal">g</mi><mi mathvariant="normal">h</mi><mi mathvariant="normal">t</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165442&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165438&type= 6.09599972 3.72533321 为左侧和右侧的归一化最大类间方差， <math id="M31"><msub><mrow><mi>k</mi></mrow><mrow><mn mathvariant="normal">1</mn></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165447&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165444&type= 2.87866688 3.72533321 和 <math id="M32"><msub><mrow><mi>k</mi></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165451&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165450&type= 2.87866688 3.72533321 分别为两侧归一化最大类间方差对应的灰度级， <math id="M33"><msub><mrow><mi>h</mi></mrow><mrow><mi mathvariant="normal">m</mi><mi mathvariant="normal">a</mi><mi mathvariant="normal">x</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165456&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165454&type= 5.24933338 3.72533321 为频率最大值对应的灰度级。 2.2.3　阈值选取区间 3 比较 <math id="M34"><msub><mrow><mi>G</mi></mrow><mrow><mi mathvariant="normal">l</mi><mi mathvariant="normal">e</mi><mi mathvariant="normal">f</mi><mi mathvariant="normal">t</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165460&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165457&type= 4.99533367 3.72533321 和 <math id="M35"><msub><mrow><mi>G</mi></mrow><mrow><mi mathvariant="normal">r</mi><mi mathvariant="normal">i</mi><mi mathvariant="normal">g</mi><mi mathvariant="normal">h</mi><mi mathvariant="normal">t</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165464&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165461&type= 6.09599972 3.72533321 ，较大一侧优先进行GKIT阈值选取。将两侧归一化最大类间方差对应灰度级 <math id="M36"><msub><mrow><mi>k</mi></mrow><mrow><mn mathvariant="normal">1</mn></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165468&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165466&type= 2.87866688 3.72533321 和 <math id="M37"><msub><mrow><mi>k</mi></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165472&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165470&type= 2.87866688 3.72533321 作为阈值候选区间的起始和终止灰度级。左侧阈值候选区间为 [1, <math id="M38"><msub><mrow><mi>k</mi></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165477&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165475&type= 2.87866688 3.72533321 ]，只包含后向散射减弱变化和非变化像素，右侧为[ <math id="M39"><msub><mrow><mi>k</mi></mrow><mrow><mn mathvariant="normal">1</mn></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165481&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165479&type= 2.87866688 3.72533321 , <math id="M40"><mi>N</mi></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165487&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165484&type= 2.53999996 2.62466669 ]，只包含增强和非变化像素。与保持直方图基本形态不变的基础上计算两侧归一化最大类间方差不同，为尽可能保证阈值选取精确，阈值选取都是在原始直方图 <math id="M41"><msub><mrow><mi>H</mi></mrow><mrow><mn mathvariant="normal">0</mn></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165490&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165488&type= 3.64066648 3.72533321 上进行。 2.2.4　确定强侧阈值 3 不同于D-GKIT，本文将同时确定双侧阈值优化为分别求取强弱两侧分割阈值，避免了变化与非变化类灰度级重叠引起的单侧双阈值错误情况发生。 KI准则的基本思想是将像元划分为隶属度较高的类别，求取准则函数最小时的对应阈值。 GG分布函数为 <math display="block" id="M42"><mi>p</mi><mo stretchy="false">(</mo><mi>k</mi><mo stretchy="false">|</mo><msub><mrow><mi>C</mi></mrow><mrow><mi>i</mi></mrow></msub><mo stretchy="false">)</mo><mo>=</mo><msub><mrow><mi>a</mi></mrow><mrow><mi>i</mi></mrow></msub><msup><mrow><mi mathvariant="normal">e</mi></mrow><mrow><mo>-</mo><mo stretchy="false">(</mo><msub><mrow><mi>b</mi></mrow><mrow><mi>i</mi></mrow></msub><msup><mrow><mfenced open="|" close="|" separators="|"><mrow><mi>k</mi><mo>-</mo><msub><mrow><mi>m</mi></mrow><mrow><mi>i</mi></mrow></msub></mrow></mfenced><mo stretchy="false">)</mo></mrow><mrow><msub><mrow><mi>β</mi></mrow><mrow><mi>i</mi></mrow></msub></mrow></msup></mrow></msup><mo>,</mo><mi>i</mi><mo>∈</mo><mo stretchy="false">(</mo><mi>u</mi><mo>,</mo><mi>c</mi><mo stretchy="false">)</mo></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165494&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165492&type= 43.85733414 6.68866634 式中， <math id="M43"><mi>p</mi><mo stretchy="false">(</mo><mi>k</mi><mo stretchy="false">|</mo><msub><mrow><mi>C</mi></mrow><mrow><mi>i</mi></mrow></msub><mo stretchy="false">)</mo></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165497&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165495&type= 10.07533360 4.57200003 表示差异图中某一类 <math id="M44"><msub><mrow><mi>C</mi></mrow><mrow><mi>i</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165500&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165498&type= 3.04799986 3.72533321 的概率密度函数， <math id="M45"><msub><mrow><mi>a</mi></mrow><mrow><mi>i</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165503&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165502&type= 2.62466669 3.72533321 和 <math id="M46"><msub><mrow><mi>b</mi></mrow><mrow><mi>i</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165506&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165505&type= 2.37066674 3.72533321 是常数，可以通过均值 <math id="M47"><msub><mrow><mi>m</mi></mrow><mrow><mi>i</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165509&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165508&type= 3.47133350 3.72533321 和方差 <math id="M48"><msubsup><mrow><mi>σ</mi></mrow><mrow><mi>i</mi></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msubsup></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165512&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165511&type= 3.64066648 4.06400013 和形状参数 <math id="M49"><msub><mrow><mi>β</mi></mrow><mrow><mi>i</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165516&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165514&type= 2.62466669 4.48733330 求出，如 式(13) 式(13) 所示 <math display="block" id="M50"><msub><mrow><mi>a</mi></mrow><mrow><mi>i</mi></mrow></msub><mo>=</mo><mfrac><mrow><msub><mrow><mi>b</mi></mrow><mrow><mi>i</mi></mrow></msub><msub><mrow><mi>β</mi></mrow><mrow><mi>i</mi></mrow></msub></mrow><mrow><mn mathvariant="normal">2</mn><mi mathvariant="normal">Γ</mi><mfenced separators="|"><mrow><mfrac bevelled="true"><mrow><mn mathvariant="normal">1</mn></mrow><mrow><msub><mrow><mi>β</mi></mrow><mrow><mi>i</mi></mrow></msub></mrow></mfrac></mrow></mfenced></mrow></mfrac><mo>,</mo><msub><mrow><mi>b</mi></mrow><mrow><mi>i</mi></mrow></msub><mo>=</mo><mfrac><mrow><mn mathvariant="normal">1</mn></mrow><mrow><msub><mrow><mi>σ</mi></mrow><mrow><mi>i</mi></mrow></msub></mrow></mfrac><mroot><mrow><mfrac><mrow><mi mathvariant="normal">Γ</mi><mfenced separators="|"><mrow><mfrac bevelled="true"><mrow><mn mathvariant="normal">3</mn></mrow><mrow><msub><mrow><mi>β</mi></mrow><mrow><mi>i</mi></mrow></msub></mrow></mfrac></mrow></mfenced></mrow><mrow><mi mathvariant="normal">Γ</mi><mfenced separators="|"><mrow><mo stretchy="false">(</mo><mfrac bevelled="true"><mrow><mn mathvariant="normal">1</mn></mrow><mrow><msub><mrow><mi>β</mi></mrow><mrow><mi>i</mi></mrow></msub></mrow></mfrac><mo stretchy="false">)</mo></mrow></mfenced></mrow></mfrac></mrow><mrow/></mroot></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165519&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165518&type= 55.37200165 19.13466644 准则函数 J 为 <math display="block" id="M51"><mi>J</mi><mo stretchy="false">(</mo><mi>T</mi><mo stretchy="false">)</mo><mo>=</mo><munderover><mo>∑</mo><mrow><mi>k</mi><mo>=</mo><mn mathvariant="normal">0</mn></mrow><mrow><mi>N</mi></mrow></munderover><mrow><msub><mrow><mi>H</mi></mrow><mrow><mn mathvariant="normal">0</mn></mrow></msub><mo stretchy="false">(</mo><mi>k</mi><mo stretchy="false">)</mo><mi>c</mi><mo stretchy="false">(</mo><mi>k</mi><mo>,</mo><mi>T</mi><mo stretchy="false">)</mo></mrow></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165522&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165521&type= 34.54399872 8.29733276 T 为分割阈值，代价函数 <math id="M52"><mi>c</mi><mo stretchy="false">(</mo><mi>k</mi><mo>,</mo><mi>T</mi><mo stretchy="false">)</mo></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165525&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165524&type= 9.48266697 3.47133350 为 <math display="block" id="M53"><mi>c</mi><mo stretchy="false">(</mo><mi>k</mi><mo>,</mo><mi>T</mi><mo stretchy="false">)</mo><mo>=</mo><mo>-</mo><mn mathvariant="normal">2</mn><mi mathvariant="normal">l</mi><mi mathvariant="normal">n</mi><mo stretchy="false">(</mo><mi>P</mi><mo stretchy="false">(</mo><msub><mrow><mi>C</mi></mrow><mrow><mi>i</mi></mrow></msub><mo stretchy="false">|</mo><mi>k</mi><mo>,</mo><mi>T</mi><mo stretchy="false">)</mo><mo stretchy="false">)</mo><mo>,</mo><mi>i</mi><mo>∈</mo><mo stretchy="false">(</mo><mi>u</mi><mo>,</mo><mi>c</mi><mo stretchy="false">)</mo></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165529&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165527&type= 51.13866806 4.57200003 式中， <math id="M54"><mi>P</mi><mo stretchy="false">(</mo><msub><mrow><mi>C</mi></mrow><mrow><mi>i</mi></mrow></msub><mo stretchy="false">|</mo><mi>k</mi><mo>,</mo><mi>T</mi><mo stretchy="false">)</mo></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165532&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165531&type= 14.13933372 4.57200003 是变化类与非变化类的后验概率。通过GKIT可在候选区间内选取合理强侧阈值 <math id="M55"><msub><mrow><mi>T</mi></mrow><mrow><mi mathvariant="normal">s</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165537&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165535&type= 3.04799986 3.72533321 。假设强侧为右侧，则阈值选取区间为[ <math id="M56"><msub><mrow><mi>k</mi></mrow><mrow><mn mathvariant="normal">1</mn></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165541&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165538&type= 2.87866688 3.72533321 , <math id="M57"><mi>N</mi></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165545&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165543&type= 2.53999996 2.62466669 ]，非变化类的拟合函数为 <math id="M58"><mi>p</mi><mi>d</mi><msub><mrow><mi>f</mi></mrow><mrow><mi mathvariant="normal">r</mi><mi mathvariant="normal">s</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165550&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165548&type= 6.18066692 4.57200003 。若强侧为左侧，则阈值选取区间为[1, <math id="M59"><msub><mrow><mi>k</mi></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165553&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165551&type= 2.87866688 3.72533321 ]，非变化类的拟合函数为 <math id="M60"><mi>p</mi><mi>d</mi><msub><mrow><mi>f</mi></mrow><mrow><mi mathvariant="normal">l</mi><mi mathvariant="normal">s</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165557&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165555&type= 6.01133347 4.57200003 。 2.2.5　更新直方图 3 对弱侧进行阈值分割时，强侧变化像素灰度级依然会影响其他灰度级GG函数的拟合，利用上步中获取的非变化类拟合函数 <math id="M61"><mi>p</mi><mi>d</mi><msub><mrow><mi>f</mi></mrow><mrow><mi mathvariant="normal">r</mi><mi mathvariant="normal">s</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165562&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165559&type= 6.18066692 4.57200003 中的[ <math id="M62"><msub><mrow><mi>h</mi></mrow><mrow><mi mathvariant="normal">m</mi><mi mathvariant="normal">a</mi><mi mathvariant="normal">x</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165566&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165564&type= 5.24933338 3.72533321 , <math id="M63"><mi>N</mi></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165570&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165568&type= 2.53999996 2.62466669 ]或 <math id="M64"><mi>p</mi><mi>d</mi><msub><mrow><mi>f</mi></mrow><mrow><mi mathvariant="normal">l</mi><mi mathvariant="normal">s</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165573&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165571&type= 6.01133347 4.57200003 中的[1, <math id="M65"><msub><mrow><mi>h</mi></mrow><mrow><mi mathvariant="normal">m</mi><mi mathvariant="normal">a</mi><mi mathvariant="normal">x</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165576&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165575&type= 5.24933338 3.72533321 ]替换原始直方图中对应区间值形成新直方图 <math id="M66"><msub><mrow><mi>H</mi></mrow><mrow><mi mathvariant="normal">n</mi><mi mathvariant="normal">e</mi><mi mathvariant="normal">w</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165580&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165579&type= 5.75733376 3.72533321 ，将直方图一侧变为标准GG分布，以减小对其他灰度级分割选取阈值的影响。 直方图左侧为强侧时有式(16) <math display="block" id="M67"><mtable><mtr><mtd><msub><mrow><mi>H</mi></mrow><mrow><mi mathvariant="normal">n</mi><mi mathvariant="normal">e</mi><mi mathvariant="normal">w</mi></mrow></msub><mfenced open="[" close="]" separators="|"><mrow><mn mathvariant="normal">1</mn><mo>,</mo><msub><mrow><mi>h</mi></mrow><mrow><mi mathvariant="normal">m</mi><mi mathvariant="normal">a</mi><mi mathvariant="normal">x</mi></mrow></msub></mrow></mfenced><mo>=</mo><mi>p</mi><mi>d</mi><msub><mrow><mi>f</mi></mrow><mrow><mi mathvariant="normal">l</mi><mi mathvariant="normal">s</mi></mrow></msub><mfenced open="[" close="]" separators="|"><mrow><mn mathvariant="normal">1</mn><mo>,</mo><msub><mrow><mi>h</mi></mrow><mrow><mi mathvariant="normal">m</mi><mi mathvariant="normal">a</mi><mi mathvariant="normal">x</mi></mrow></msub></mrow></mfenced></mtd></mtr><mtr><mtd><msub><mrow><mi>H</mi></mrow><mrow><mi mathvariant="normal">n</mi><mi mathvariant="normal">e</mi><mi mathvariant="normal">w</mi></mrow></msub><mo stretchy="false">(</mo><msub><mrow><mi>h</mi></mrow><mrow><mi mathvariant="normal">m</mi><mi mathvariant="normal">a</mi><mi mathvariant="normal">x</mi></mrow></msub><mo>,</mo><mi>N</mi><mo stretchy="false">]</mo><mo>=</mo><msub><mrow><mi>H</mi></mrow><mrow><mn mathvariant="normal">0</mn></mrow></msub><mo stretchy="false">(</mo><msub><mrow><mi>h</mi></mrow><mrow><mi mathvariant="normal">m</mi><mi mathvariant="normal">a</mi><mi mathvariant="normal">x</mi></mrow></msub><mo>,</mo><mi>N</mi><mo stretchy="false">]</mo></mtd></mtr></mtable></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165583&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165582&type= 39.36999893 11.51466751 直方图右侧为强侧的情况对应式(17) <math display="block" id="M68"><mtable><mtr><mtd><msub><mrow><mi>H</mi></mrow><mrow><mi mathvariant="normal">n</mi><mi mathvariant="normal">e</mi><mi mathvariant="normal">w</mi></mrow></msub><mo stretchy="false">[</mo><mn mathvariant="normal">1</mn><mo>,</mo><msub><mrow><mi>h</mi></mrow><mrow><mi mathvariant="normal">m</mi><mi mathvariant="normal">a</mi><mi mathvariant="normal">x</mi></mrow></msub><mo stretchy="false">)</mo><mo>=</mo><msub><mrow><mi>H</mi></mrow><mrow><mn mathvariant="normal">0</mn></mrow></msub><mo stretchy="false">[</mo><mn mathvariant="normal">1</mn><mo>,</mo><msub><mrow><mi>h</mi></mrow><mrow><mi mathvariant="normal">m</mi><mi mathvariant="normal">a</mi><mi mathvariant="normal">x</mi></mrow></msub><mo stretchy="false">)</mo></mtd></mtr><mtr><mtd><msub><mrow><mi>H</mi></mrow><mrow><mi mathvariant="normal">n</mi><mi mathvariant="normal">e</mi><mi mathvariant="normal">w</mi></mrow></msub><mfenced open="[" close="]" separators="|"><mrow><msub><mrow><mi>h</mi></mrow><mrow><mi mathvariant="normal">m</mi><mi mathvariant="normal">a</mi><mi mathvariant="normal">x</mi></mrow></msub><mo>,</mo><mi>N</mi></mrow></mfenced><mo>=</mo><mi>p</mi><mi>d</mi><msub><mrow><mi>f</mi></mrow><mrow><mi mathvariant="normal">r</mi><mi mathvariant="normal">s</mi></mrow></msub><mfenced open="[" close="]" separators="|"><mrow><msub><mrow><mi>h</mi></mrow><mrow><mi mathvariant="normal">m</mi><mi mathvariant="normal">a</mi><mi mathvariant="normal">x</mi></mrow></msub><mo>,</mo><mi>N</mi></mrow></mfenced></mtd></mtr></mtable></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165586&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165585&type= 41.06333542 10.41399956 2.2.6　确定弱侧阈值 3 对弱侧进行阈值分割时，新直方图只包含非变化类和弱侧变化像素灰度级，避免了强侧变化像素灰度级的影响。根据阈值候选区间利用GKIT对新直方图 <math id="M69"><msub><mrow><mi>H</mi></mrow><mrow><mi mathvariant="normal">n</mi><mi mathvariant="normal">e</mi><mi mathvariant="normal">w</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165589&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165588&type= 5.75733376 3.72533321 在候选区间上进行分割得到弱侧分割阈值 <math id="M70"><msub><mrow><mi>T</mi></mrow><mrow><mi mathvariant="normal">w</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165592&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165590&type= 3.47133350 3.72533321 。由此得到了两侧分割阈值 <math id="M71"><msub><mrow><mi>T</mi></mrow><mrow><mi mathvariant="normal">w</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165594&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165593&type= 3.47133350 3.72533321 和 <math id="M72"><msub><mrow><mi>T</mi></mrow><mrow><mi mathvariant="normal">s</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165597&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165596&type= 3.04799986 3.72533321 。 3　试验与分析 1 3.1　试验数据 2 GF-3作为中国首颗高分辨率SAR卫星，是高分对地观测系统中唯一的民用微波遥感卫星，具有高分辨率，多成像模式和大成像幅宽等特点，在各邻域都有着重要利用价值( 张庆君，2017 张庆君，2017 )。其成像基本信息如 表1 表1 。本文以GF-3宁波地区的强度影像作为双阈值变化检测的研究对象，试验区数据基本信息如 表2 表2 。利用SAR影像处理软件GDEMSI完成强度提取，配准，多视，裁剪等预处理流程，配准精度为亚像素级，最终试验中采用定标后未做地理编码的强度影像，像素深度分别为[0,3159]和[0,3146] 。为使影像符合人类视觉系统，使方位向和距离向分辨率相近，对原始影像采用5×4多视处理，处理后的分辨率均在13.5 m左右。预处理完成后试验影像大小为4071×3753。 图2 图2 (a)和 图2 图2 (b)为前后时相影像。 2 2 成像模式 3 分辨率/m 2 成像幅宽/km 2 极化方式 标称 方位向 距离向 标称 范围 6 条带成像模式 超精细条带 3 3 2.5—5.0 30 30 可选单极化 精细条带1 5 5 4—6 50 50 可选双极化 精细条带2 10 10 8—12 100 95—110 可选双极化 标准条带 25 25 15—30 130 95—150 可选双极化 全极化条带1 8 8 6—9 30 20—35 全极化 全极化条带2 25 25 15—30 40 35—50 全极化 2 聚束模式 1 1.0—1.5 0.9—2.5 10×10 10×10 可选单极化 3 扫描成像模式 窄幅扫描 50 50—60 30—60 300 300 可选双极化 宽幅扫描 100 100 50—110 500 500 可选双极化 全球观测成像模式 500 500 350—700 650 650 可选双极化 2 波成像模式 10 10 8—12 5×5 5×5 全极化 2 扩展入射角模式 低入射角 25 25 15—30 130 120—150 可选双极化 高入射角 25 25 20—30 80 70—90 可选双极化 基本信息 前时相 后时相 获取时间 2016-11-14 2017-03-10 成像模式 数据级别 精细条带(FSI) L1A(SLC) 精细条带(FSI) L1A(SLC) 轨道方式 降轨 降轨 极化方式 HH HH 图像尺寸 21525×16285 21560×16884 入射角 42.918° 42.883° 分辨率 2.868×2.248/m 2.863×2.248/m 侧视 右视 右视 1,2; 试验区内包括海洋，湖泊，农田，裸地，建筑区，港口等。沿海地区水产养殖，舰船位移和地表覆盖及土地含水量变化等是造成区域内地表覆盖变化的主要原因，且变化类型及差异程度也不甚相同。 3.2　差异图生成 2 在生成差异图之前，采用3×3的增强Lee滤波对原始图像进行滤波处理。由于时间间隔较短(4个月)，发生剧烈变化的区域较少，造成了直方图不同类别间灰度级发生重叠。而且研究采用的影像像幅较大，理论上随着样本(像素)个数的增加，差异图直方图应与GG分布有更好的拟合效果。差异图直方图及平滑处理后结果如 图3 图3 ，像素深度为[-9.49，9.50]，灰度范围较小，所以本文参考8位图的方式将原始影像范围拉伸为256个区间进行阈值选取。 3.3　阈值选取及结果分析 2 D-GKIT阈值为145和154，根据 Bazi等(2006) Bazi等(2006) ，利用2维曲面分析其准则函数形态，如 图4 图4 。利用最大值灰度级( <math id="M73"><msub><mrow><mi>h</mi></mrow><mrow><mi mathvariant="normal">m</mi><mi mathvariant="normal">a</mi><mi mathvariant="normal">x</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165624&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165623&type= 5.24933338 3.72533321 =125)( 图4 图4 中红线和紫线)将阈值平面分为4个象限：第1象限： T 1， T 2小于 <math id="M74"><msub><mrow><mi>h</mi></mrow><mrow><mi mathvariant="normal">m</mi><mi mathvariant="normal">a</mi><mi mathvariant="normal">x</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165627&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165626&type= 5.24933338 3.72533321 ；第2象限： T 1大于 <math id="M75"><msub><mrow><mi>h</mi></mrow><mrow><mi mathvariant="normal">m</mi><mi mathvariant="normal">a</mi><mi mathvariant="normal">x</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165630&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165629&type= 5.24933338 3.72533321 ， T 2小于 <math id="M76"><msub><mrow><mi>h</mi></mrow><mrow><mi mathvariant="normal">m</mi><mi mathvariant="normal">a</mi><mi mathvariant="normal">x</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165633&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165632&type= 5.24933338 3.72533321 ；第3象限： T 1， T 2大于 <math id="M77"><msub><mrow><mi>h</mi></mrow><mrow><mi mathvariant="normal">m</mi><mi mathvariant="normal">a</mi><mi mathvariant="normal">x</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165636&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165635&type= 5.24933338 3.72533321 ；第4象限： T 1小于 <math id="M78"><msub><mrow><mi>h</mi></mrow><mrow><mi mathvariant="normal">m</mi><mi mathvariant="normal">a</mi><mi mathvariant="normal">x</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165639&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165638&type= 5.24933338 3.72533321 ， T 2大于 <math id="M79"><msub><mrow><mi>h</mi></mrow><mrow><mi mathvariant="normal">m</mi><mi mathvariant="normal">a</mi><mi mathvariant="normal">x</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165642&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165640&type= 5.24933338 3.72533321 。试验中准则函数 J 最小值( 图4 图4 中紫圈位置)投影位于第3象限中，证明其能够得到2个分割阈值。但考虑两阈值均出现在直方图尖峰的一侧，明显与正确分割时最小值应出现在第4象限不符，显示出D-GKIT在进行双阈值分割时盲目的从全局选取双阈值而忽略阈值应出现在正确区间的不足。 图5 图5 是D-GKIT函数拟合结果，由于直方图左侧变化类灰度级与非变化类灰度级重叠严重，使 图5 图5 (a)中明显是非变化类的拟合结果被错误判定为左侧变化，而 图5 图5 (c)将误认为的“非变化类”的区间进行了强行拟合，且右侧拟合结果如 图5 图5 (b)的拟合效果也较差，证明D-GKIT在这种情况下无法得到正确的分割阈值。 1,2,3; 利用归一化最大类间方差进行形态判断，具体形态参数如 表3 表3 ，左侧归一化最大类间方差为0.4099小于右侧值0.5273，证明直方图右侧受不同类别灰度级重叠影响较小，应优先进行阈值选取。 参数 左侧 右侧 G 0.4099 0.5273 k 119 134 右侧的阈值候选区间为[119,255]，GKIT分割得到阈值为142，利用拟合的非变化类GG函数代替原始直方图中的右侧区间形成新直方图，弱侧阈值候选区间为[1,134] ，分割得到阈值为99。 图6 图6 (b)为优先选取阈值的右侧的拟合结果，右上方空白处为非变化类拟合结果； 图6 图6 (a)为左侧拟合结果，左上方空白处为针对新生成的直方图中非变化类在区间内的拟合结果。通过给定阈值候选区间，在灰度级重叠情况下利用更新后直方图在候选区间进行阈值选取，避免了类似于D-GKIT因全局选取阈值而出现的错误分割结果，能够稳定得到可靠的分割阈值。 1,2; 表4 表4 为D-GKIT和本文算法的拟合参数。D-GKIT左侧拟合结果和本文算法中非变化类相差不大，原因是左侧变化类与非变化类灰度级重叠严重，D-GKIT无法对这一侧进行较好的分割。 2 方法 a B β m σ 3 D-GKIT 左侧 0.11 0.2106 1.1 124.93 5.58 不变 0.1313 0.2353 3.2 148.01 2.56 右侧 0.1021 0.2459 0.7 162.58 10.92 4 本文算法 右侧 0.2057 0.5207 0.7 147.89 6.01 不变 0.1601 0.1991 1.2 125.29 5.39 新不变 0.1051 0.1948 1.28 125.37 5.16 左侧 0.0846 0.19 0.8 90.09 11.29 由于LR对弱反射区域变化敏感度高(灰度值从1变为2和从100变化为200对数比值相同)，使水面容易产生对数比值虚高的虚警情况，尤其未经拉伸的原始影像更为明显，这也是LR的主要缺点之一。而研究区内存在大量水域，因此本文通过设定两时相间差值阈值( T =10)的方法对检测结果进行处理，消除水域间比值虚高带来的影响。即当两时相影像中对应像素绝对值之差小于10时， <math id="M80"><mi>L</mi></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165687&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165686&type= 1.86266661 2.62466669 中的对应像素在变化结果图中均标记为未变化，如式(18) <math display="block" id="M81"><mfenced open="|" close="|" separators="|"><mrow><msub><mrow><mi>E</mi></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msub><mo>-</mo><msub><mrow><mi>E</mi></mrow><mrow><mn mathvariant="normal">1</mn></mrow></msub></mrow></mfenced><mo>≤</mo><mn mathvariant="normal">10</mn></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165689&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21165688&type= 21.92866707 5.16466665 对差异图进行3×3的中值滤波后，最终检测结果中未发生变化的像素灰度值标记为灰色，后向散射减弱和增强像素分别标记为绿色和红色。去除对数比值虚高后的最终检测结果如 图7 图7 。时间间隔4个月，区域内发生的变化主要有(1)沿海滩涂变化；(2)海面舰船位置变化；(3)地表土地利用类型变化。目视结果本算法能避免灰度级重叠引起的单侧双阈值错分情况，而D-GKIT缺陷明显。选取3个具有代表性区域进行对比，如 图8 图8 第1列是沿海滩涂变化区域，第2列是海面舰船位移变化区域，第3列为土地覆盖变化区域。因为单侧双阈值的错误情况，D-GKIT无法正确得到减弱变化像素检测结果，增强和非变化区域也出现检测结果不稳定的错分情况，而本文算法基本正确的检测出了大致变化情况。 1,2; 1;2;3; 选取3个样本区域利用混淆矩阵对本文算法进行精度评价，因D-GKIT错分严重，不再对其进行精度评价。精度评价区域的地表真实变化结果均参考邻近时间光学影像解译获得。区域一( 图8 图8 )是位于 图7 图7 中左下方的宁波市北仑区沿海滩涂区域，大小为6468像素。沿海滩涂的利用使两时相影像间既有增强变化也有减弱变化。利用本文算法进行检测，混淆矩阵精度评价结果如 表5 表5 ，减弱、无变化和增强变化的检测正确率分别为56.1%，98.43%和92.1%，而D-GKIT算法由于只在灰度级重叠度较小的右侧选取了两个阈值，而将没有发生变化的像素都误判为后向散射减弱像素且将增强变化错分为非变化和增强变化两种类别。区域二( 图8 图8 )是位于 图7 图7 中左上方北仑区的海面区域，大小为2500像素。因靠近港口，船舶往来密集，由于船舶位移造成两种类型变化。精度评价结果如 表6 表6 ，减弱、无变化和增强的检测精度分别为68.75%，97.43%和61.76，D-GKIT没有检测到驶出本区域的船只。区域三( 图8 图8 )为 图7 图7 中右上方镇海区农田区域，大小为2500像素。相隔4个月，农作物生长及含水量变化是造成本区域内发生增强和减弱变化的主要原因。混淆矩阵精度评价结果如 表7 表7 ，减弱，非变化和增强的检测正确率分别为80%，99.96%和76.19%，相比于本文算法，D-GKIT完全没有检测出后向散射减弱的变化像素。 减弱变化 无变化 增强变化 减弱变化 56.1 43.9 0 无变化 0.04 98.43 1.53 增强变化 0 7.90 92.1 减弱变化 无变化 增强变化 减弱变化 68.75 31.25 0 无变化 0.5 97.43 2.07 增强变化 0 38.24 61.76 减弱变化 无变化 增强变化 减弱变化 80 20 0 无变化 0 99.96 0.04 增强变化 0 23.81 76.19 通过对不同类型数据进行结果对比，并进行精度评定。证明相比于D-GKIT算法，本文算法能较好的应对灰度级重叠情况且具有较好检测精度。通过利用GF-3数据进行试验，达到了验证利用国产GF-3号卫星影像进行SAR变化检测研究的目的。 4　结　论 1 针对D-GKIT在进行SAR变化检测阈值分割时因直方图灰度级重叠而无法正确分割的问题，本文提出了一种结合归一化最大类间方差和GKIT的双阈值分割算法。通过对1维高斯核函数平滑后的差异图直方图尖峰两侧进行归一化最大类间方差计算，判断两侧不同类别像素灰度级间的重叠程度，确定阈值选取顺序及候选区间，再利用GKIT针对原始直方图和更新后的直方图在不同区间上分别进行分割，得到两侧的分割阈值。利用宁波地区GF-3影像将本文算法与D-GKIT进行试验并对比分析，结果表明，本文方法在不同类别灰度级发生重叠时依然具有良好的分割可靠性和较高的检测精度，一方面证明了本文方法可靠有效，另一方面也验证了利用GF-3 SAR影像进行变化检测研究的可行性。 本文算法的创新点在于：(1)提出以归一化最大类间方差值作为不同类别间灰度级重叠程度的判决因素，并以此在1维高斯核平滑后的直方图上确定阈值选取顺序和候选区间，将阈值选取从全局优化为两合理区间，既减少了运算量也基本杜绝了灰度级重叠引发的单侧双阈值错误情况的发生；(2)提出以非变化类拟合函数更新后的直方图作为重叠较大一侧的阈值选取基础，减小了一侧变化像素对其他灰度级分割选取阈值时的影响，提高了检测精度，避免出现了全局选取阈值的盲目性。并且可以预见的是，随着图像像幅的增加，差异图直方图将更加趋近于一个标准的“尖峰”型GG分布，灰度级重叠程度会更明显，本文方法相对于D-GKIT算法的优势会更突出。 分别获取直方图两侧的阈值能较好地解决灰度级重叠引起的单侧双阈值错误现象，但算法依然在存在一些不足。例如没有利用到图像邻域信息，易受孤点噪声干扰；另外，对于变化类灰度级近似于隐含在非变化类 (变化差异小，变化像素少) 的情况，无法进行较好分割，这是由LR差异图不同类别间像素灰度级分离度较小的原因造成的，可以通过改进差异图生成方式，增大不同类别间的差异程度来应对这一不足，这也是未来的研究重点。 参考文献( References ) Bazi Y , Bruzzone L and Melgani F . 2005 . An unsupervised approach based on the generalized gaussian model to automatic change detection in multitemporal SAR images . IEEE Transactions on Geoscience and Remote Sensing , 43 ( 4 ): 874 - 887 [ DOI: 10.1109/TGRS.2004.842441 http://dx.doi.org/10.1109/TGRS.2004.842441 ] An unsupervised approach based on the generalized gaussian model to automatic change detection in multitemporal SAR images Bazi Y F , Bruzzone L and Melgani F . 2006 . Automatic identification of the number and values of decision thresholds in the log-ratio image for change detection in SAR images . IEEE Geoscience and Remote Sensing Letters , 3 ( 3 ): 349 - 353 [ DOI: 10.1109/LGRS.2006.869973 http://dx.doi.org/10.1109/LGRS.2006.869973 ] Automatic identification of the number and values of decision thresholds in the log-ratio image for change detection in SAR images Demirkaya O and Asyali M H . 2004 . Determination of image bimodality thresholds for different intensity distributions . Signal Processing : Image Communication , 19 ( 6 ): 507 - 516 [ DOI: 10.1016/j.image.2004.04.002 http://dx.doi.org/10.1016/j.image.2004.04.002 ] Determination of image bimodality thresholds for different intensity distributions Gao C S , Zhang H , Wang C and Wu F . 2010 . SAR change detection based on generalized Gamma distribution divergence and auto-threshold segmentation . Journal of Remote Sensing , 14 ( 4 ): 710 - 724 SAR change detection based on generalized Gamma distribution divergence and auto-threshold segmentation 高丛珊 , 张红 , 王超 , 吴樊 . 2010 . 广义Gamma模型及自适应KI阈值分割的SAR图像变化检测 . 遥感学报 , 14 ( 4 ): 710-724 [ DOI: 10.11834/jrs.20100407 http://dx.doi.org/10.11834/jrs.20100407 ] 广义Gamma模型及自适应KI阈值分割的SAR图像变化检测 Gong M G , Li Y , Jiao L C , Jia M and Su L Z . 2014 . SAR change detection based on intensity and texture changes . ISPRS Journal of Photogrammetry and Remote Sensing , 93 : 123 - 135 [ DOI: 10.1016/j.isprsjprs.2014.04.010 http://dx.doi.org/10.1016/j.isprsjprs.2014.04.010 ] SAR change detection based on intensity and texture changes Gong M G , Su L Z , Li H and Liu J . 2016 . A survey on change detection in synthetic aperture radar imagery . Journal of Computer Research and Development , 53 ( 1 ): 123 - 137 A survey on change detection in synthetic aperture radar imagery 公茂果 , 苏临之 , 李豪 , 刘嘉 . 2016 . 合成孔径雷达影像变化检测研究进展 . 计算机研究与发展 , 53 ( 1 ): 123-137 [ DOI: 10.7544/issn1000-1239.2016.20150662 http://dx.doi.org/10.7544/issn1000-1239.2016.20150662 ] 合成孔径雷达影像变化检测研究进展 Hame T , Heiler I and Miguel-Ayanz J S . 1998 . An unsupervised change detection and recognition system for forestry . International Journal of Remote Sensing , 19 ( 6 ): 1079 - 1099 [ DOI: 10.1080/014311698215612 http://dx.doi.org/10.1080/014311698215612 ] An unsupervised change detection and recognition system for forestry Hao H M , Zhang Y H , Shi H Y and Huang J B . 2012 . Application of test statistic method in fully polarimtric SAR change detection . Journal of Remote Sensing , 16 ( 3 ): 520 - 532 Application of test statistic method in fully polarimtric SAR change detection 郝洪美 , 张永红 , 石海燕 , 黄金波 . 2012 . 统计假设检验方法在全极化SAR变化检测中的应用 . 遥感学报 , 16 ( 3 ): 520-532 [ DOI: 10.11834/jrs.20120384 http://dx.doi.org/10.11834/jrs.20120384 ] 统计假设检验方法在全极化SAR变化检测中的应用 Hu H T and Ban Y F . 2014 . Unsupervised change detection in multitemporal SAR images over large urban areas . IEEE Journal of Selected Topics in Applied Earth Observations and Remote Sensing , 7 ( 8 ): 3248 - 3261 [ DOI: 10.1109/JSTARS.2014.2344017 http://dx.doi.org/10.1109/JSTARS.2014.2344017 ] Unsupervised change detection in multitemporal SAR images over large urban areas Hu Z L . 2013 . An unsupervised change detection approach based on KI dual thresholds under the generalized gauss model assumption in SAR images . Acta Geodaetica et Cartographica Sinica , 42 ( 1 ): 116 - 122 An unsupervised change detection approach based on KI dual thresholds under the generalized gauss model assumption in SAR images 胡召玲 . 2013 . 广义高斯模型及KI双阈值法的SAR图像非监督变化检测 . 测绘学报 , 42 ( 1 ): 116 - 122 广义高斯模型及KI双阈值法的SAR图像非监督变化检测 Huang S Q , Liu D Z , Hu M X and Wang S C . 2010 . Multi-temporal SAR image change detection technique based on wavelet transform . Acta Geodaetica et Cartographica Sinica , 39 ( 2 ): 180 - 186 Multi-temporal SAR image change detection technique based on wavelet transform 黄世奇 , 刘代志 , 胡明星 , 王仕成 . 2010 . 基于小波变换的多时相SAR图像变化检测技术 . 测绘学报 , 39 ( 2 ): 180 - 186 基于小波变换的多时相SAR图像变化检测技术 Kittler J and Illingworth J . 1986 . Minimum error thresholding . Pattern Recognition , 19 ( 1 ): 41 - 47 [ DOI: 10.1016/0031-3203(86)90030-0 http://dx.doi.org/10.1016/0031-3203(86)90030-0 ] Minimum error thresholding Ma J J , Gong M G and Zhou Z Q . 2012 . Wavelet fusion on ratio images for change detection in SAR images . IEEE Geoscience and Remote Sensing Letters , 9 ( 6 ): 1122 - 1126 [ DOI: 10.1109/LGRS.2012.2191387 http://dx.doi.org/10.1109/LGRS.2012.2191387 ] Wavelet fusion on ratio images for change detection in SAR images Moser G and Serpico S B . 2006 . Generalized minimum-error thresholding for unsupervised change detection from SAR amplitude imagery . IEEE Transactions on Geoscience and Remote Sensing , 44 ( 10 ): 2972 - 2982 [ DOI: 10.1109/TGRS.2006.876288 http://dx.doi.org/10.1109/TGRS.2006.876288 ] Generalized minimum-error thresholding for unsupervised change detection from SAR amplitude imagery Otsu N . 1979 . A threshold selection method from gray-level histograms . IEEE Transactions on Systems, Man, and Cybernetics , 9 ( 1 ): 62 - 66 [ DOI: 10.1109/TSMC.1979.4310076 http://dx.doi.org/10.1109/TSMC.1979.4310076 ] A threshold selection method from gray-level histograms Pantze A , Santoro M and Fransson J E S . 2014 . Change detection of boreal forest using bi-temporal ALOS PALSAR backscatter data . Remote Sensing of Environment , 155 : 120 - 128 [ DOI: 10.1016/j.rse.2013.08.050 http://dx.doi.org/10.1016/j.rse.2013.08.050 ] Change detection of boreal forest using bi-temporal ALOS PALSAR backscatter data Tsai D M . 1995 . A fast thresholding selection procedure for multimodal and unimodal histograms . Pattern Recognition Letters , 16 ( 6 ): 653 - 666 [ DOI: 10.1016/0167-8655(95)80011-H http://dx.doi.org/10.1016/0167-8655(95)80011-H ] A fast thresholding selection procedure for multimodal and unimodal histograms Whittle M , Quegan S , Uryu Y , Stüewe M and Yulianto K . 2012 . Detection of tropical deforestation using ALOS-PALSAR: a Sumatran case study . Remote Sensing of Environment , 124 : 83 - 98 [ DOI: 10.1016/j.rse.2012.04.027 http://dx.doi.org/10.1016/j.rse.2012.04.027 ] Detection of tropical deforestation using ALOS-PALSAR: a Sumatran case study Xiong B L , Chen J M and Kuang G Y . 2012 . A change detection measure based on a likelihood ratio and statistical properties of SAR intensity images . Remote Sensing Letters , 3 ( 3 ): 267 - 275 [ DOI: 10.1080/01431161.2011.572093 http://dx.doi.org/10.1080/01431161.2011.572093 ] A change detection measure based on a likelihood ratio and statistical properties of SAR intensity images Yousif O and Ban Y F . 2014 . Improving SAR-based urban change detection by combining MAP-MRF classifier and nonlocal means similarity weights . IEEE Journal of Selected Topics in Applied Earth Observations and Remote Sensing , 7 ( 10 ): 4288 - 4300 [ DOI: 10.1109/JSTARS.2014.2347171 http://dx.doi.org/10.1109/JSTARS.2014.2347171 ] Improving SAR-based urban change detection by combining MAP-MRF classifier and nonlocal means similarity weights Zhang Q J . 2017 . System design and key technologies of the GF-3 satellite . Acta Geodaetica et Cartographica Sinica , 46 ( 3 ): 269 - 277 System design and key technologies of the GF-3 satellite 张庆君 . 2017 . 高分三号卫星总体设计与关键技术 . 测绘学报 , 46 ( 3 ): 269-277 [ DOI: 10.11947/j.AGCS.2017.20170049 http://dx.doi.org/10.11947/j.AGCS.2017.20170049 ] 高分三号卫星总体设计与关键技术 Zhang X Q , Xiong B L and Kuang G Y . 2015 . A ship target discrimination method based on change detection in SAR imagery . Journal of Electronics and Information Technology , 37 ( 1 ): 63 - 70 A ship target discrimination method based on change detection in SAR imagery 张小强 , 熊博莅 , 匡纲要 . 2015 . 一种基于变化检测技术的SAR图像舰船目标鉴别方法 . 电子与信息学报 , 37 ( 1 ): 63-70 [ DOI: 10.11999/JEIT140143 http://dx.doi.org/10.11999/JEIT140143 ] 一种基于变化检测技术的SAR图像舰船目标鉴别方法 10.11834/jrs.20208179 1007-4619（2020）01-0001-10 崔斌,张永红,闫利,魏钜杰.2020.高分三号SAR影像双阈值变化检测.遥感学报,24(1): 1-10Cui B,Zhang Y H,Yan L and Wei J J. 2020. Dual-thresholds change detection in GF-3 SAR images. Journal of Remote Sensing(Chinese), 24(1): 1-10[DOI： 10.11834/jrs.20208179 http://dx.doi.org/10.11834/jrs.20208179 ] 崔斌,张永红,闫利,魏钜杰.2020.高分三号SAR影像双阈值变化检测.遥感学报,24(1): 1-10Cui B,Zhang Y H,Yan L and Wei J J. 2020. Dual-thresholds change detection in GF-3 SAR images. Journal of Remote Sensing(Chinese), 24(1): 1-10[DOI：10.11834/jrs.20208179] CUI Bin,ZHANG Yonghong,YAN Li,Wei Jujie. 2020. Dual-thresholds change detection in GF-3 SAR images. Journal of Remote Sensing(Chinese). 24(1): 1-10 CUI Bin,ZHANG Yonghong,YAN Li,Wei Jujie. 2020. Dual-thresholds change detection in GF-3 SAR images. Journal of Remote Sensing(Chinese). 24(1): 1-10 国家自然科学基金(编号：41271430，41801284)；中国测绘科学研究院基本科研业务费项目(编号：7771624,7771813)；城市空间信息工程北京市重点实验室经费资助项目(编号: 2017201) Supported by National Natural Science Foundation of China (No.41271430, No.41801284) 2018-04-10 2020-01-07 2021-10-18 2018-09-07 遥感学报 1 24