CN109726649B - 遥感影像云检测方法、系统及电子设备 - Google Patents

遥感影像云检测方法、系统及电子设备 Download PDF

Info

Publication number
CN109726649B
CN109726649B CN201811537640.0A CN201811537640A CN109726649B CN 109726649 B CN109726649 B CN 109726649B CN 201811537640 A CN201811537640 A CN 201811537640A CN 109726649 B CN109726649 B CN 109726649B
Authority
CN
China
Prior art keywords
threshold
image
component image
remote sensing
principal component
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN201811537640.0A
Other languages
English (en)
Other versions
CN109726649A (zh
Inventor
刘军
付华联
陈劲松
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Shenzhen Institute of Advanced Technology of CAS
Original Assignee
Shenzhen Institute of Advanced Technology of CAS
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Shenzhen Institute of Advanced Technology of CAS filed Critical Shenzhen Institute of Advanced Technology of CAS
Priority to CN201811537640.0A priority Critical patent/CN109726649B/zh
Publication of CN109726649A publication Critical patent/CN109726649A/zh
Application granted granted Critical
Publication of CN109726649B publication Critical patent/CN109726649B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Landscapes

  • Image Analysis (AREA)
  • Image Processing (AREA)

Abstract

本申请涉及一种遥感影像云检测方法、系统及电子设备。所述方法包括:步骤a:对多波段遥感影像进行主成分分析,提取所述多波段遥感影像的第一主成分图像和第二主成分图像,并进行自适应S曲线增强;步骤b:利用至少两种阈值方法分别得到所述自适应S曲线增强后的第一主成分图像和第二主成分图像的至少两种二值化阈值结果;步骤c:利用集成阈值法分别对所述第一主成分图像和第二主成分图像的至少两种二值化阈值结果进行集成,并将所述第一主成分图像和第二主成分图像的阈值集成结果进行合并,得到所述遥感影像的云检测结果。本申请有效改善云检测中的错检和漏检情况,提高云检测效率,提高遥感影像的利用率。

Description

遥感影像云检测方法、系统及电子设备
技术领域
本申请属于遥感影像云检测技术领域,特别涉及一种遥感影像云检测方法、系统及电子设备。
背景技术
随着成像技术的不断发展,遥感图像获取的渠道越来越多。在风云遥感影像获取过程中,受云雾干扰等因素的影响,导致原地物光谱失真,影响遥感产品与影像判读,对信息提取造成很大的影响。正确地分离遥感图像中的有云像元和无云像元,对于天气预报、气象灾害的预防、温度反演、救援及生态环境的监控等都有重大的影响。因此,在遥感影像的使用过程中,对受云层遮挡的遥感图像进行有效的云检测是遥感数据修复过程中需解决的首要问题。
国内外利用遥感图像的光谱信息来进行云检测已有多年的历史,目前常见的云检测方法有:
1)阈值法;阈值法包括:光谱结合阈值法、频率结合阈值法。光谱结合阈值法主要利用云在第二主成分图像具有强反射的特性,该类算法对阈值的敏感程度较高,同一卫星数据因时间、天气等原因,检测阈值将发生巨大变化,加大了此类方法的局限性。频率结合阈值法主要利用云的低频特性,通过小波分析、傅里叶变换等方法获取影像低频数据进行云检测,但由于受地面低频信息干扰,通常采用多层小波变换排除,这大大降低了云检测效率。同时,普通的单通道阈值法面对大范围以及复杂下垫面情况下的云检测效果较差;多通道多阈值法的阈值在不同环境下阈值波动较大导致阈值选取困难,并且对先验知识要求高,一些动态阈值法对于额外附加信息要求高,实时性低,算法运行速度较慢。
2)纹理分析法;纹理分析法利用云与地面纹理特征差异,常以分块子图为单位,结合二阶矩、分形维数、灰度共生矩阵和多次双边滤波进行纹理特征计算,该类方法需要提前获得可靠云特征区间才能保证分类的精度,效率较低。
3)统计学方法;统计学方法主要分为统计方程和聚类分析法。统计方程法利用样本数据建立模拟公式计算云的反射率或者亮温来进行云检测,聚类分析法是根据不同地物类型的像元观测值存在明显的差别的原理实现云检测,在样本量较大时,要获得聚类结论有一定困难,需要人为干预,极大影响检测效率。
4)综合智能法;综合智能法主要包括人工神经网络、支持向量机和模糊逻辑算法等。综合智能法在实现的过程中需要获取大量的训练样本,对分类特征的选取要求较高,针对不同数据需要重新选取样本,导致效率低下。
另外,现有的云检测方法在卫星影像云检测方面仍然会有将部分晴空、陆地、海洋检测为云以及部分漏检情况,云检测效率较低。
发明内容
本申请提供了一种遥感影像云检测方法、系统及电子设备,旨在至少在一定程度上解决现有技术中的上述技术问题之一。
为了解决上述问题,本申请提供了如下技术方案:
一种遥感影像云检测方法,包括以下步骤:
步骤a:对多波段遥感影像进行主成分分析,提取所述多波段遥感影像的第一主成分图像和第二主成分图像,并对所述第一主成分图像和第二主成分图像分别进行自适应S曲线增强;
步骤b:利用至少两种阈值方法分别得到所述自适应S曲线增强后的第一主成分图像和第二主成分图像的至少两种二值化阈值结果;
步骤c:利用集成阈值法分别对所述第一主成分图像和第二主成分图像的至少两种二值化阈值结果进行集成,并将所述第一主成分图像和第二主成分图像的阈值集成结果进行合并,得到所述遥感影像的云检测结果。
本申请实施例采取的技术方案还包括:在所述步骤a中,所述对第一主成分图像和第二主成分图像分别进行自适应S曲线增强具体为:设遥感影像的最大灰度值为r,最小灰度值为s,当r与s的和与2取余计算的结果等于0时,令初始阈值为T0=(r+s)/2,否则令初始阈值为T0=(r+s-1)/2,将图G的像素值G(x,y)小于T0的集合记为f1(x,y),图G的像素值G(x,y)大于等于T0的集合记为f2(x,y),S曲线拉伸的像素值F(x,y)计算公式为:
Figure GDA0002723055130000031
上述公式中,k为拉伸系数。
本申请实施例采取的技术方案还包括:在所述步骤b中,所述至少两种阈值方法为十种阈值方法,所述十种阈值方法分别为:大津法、分块大津法、局部动态阈值法、全局阈值与局部阈值结合、Wellner自适应阈值、最小误差法、双峰法、迭代阈值法、最大熵阈值法和固定阈值分割法。
本申请实施例采取的技术方案还包括:在所述步骤c中,所述利用集成阈值法分别对所述第一主成分图像和第二主成分图像的至少两种二值化阈值结果进行集成,并将所述第一主成分图像和第二主成分图像的阈值集成结果进行合并,得到所述遥感影像的云检测结果具体为:使用投票的方式作为结合策略,投票系数δ决定结合策略的结合程度;设所述十种阈值方法得到的云检测结果为(F1,F2,......Fp),Fp为第p种阈值方法对应的云检测结果,使用δ序列(δ12,......,δi)作为投票系数,得到云检测结果
Figure GDA0002723055130000041
如果每个像素标记得到的投票数超过投票系数δi,则检测结果为云,否则为非云;最后,选择最佳的
Figure GDA0002723055130000042
值为集成阈值的云检测结果。
本申请实施例采取的另一技术方案为:一种遥感影像云检测系统,包括:
波段提取模块:用于对多波段遥感影像进行主成分分析,提取所述多波段遥感影像的第一主成分图像和第二主成分图像;
云特征增强模块:用于对所述第一主成分图像和第二主成分图像分别进行自适应S曲线增强;
二值化模块:用于利用至少两种阈值方法分别得到所述自适应S曲线增强后的第一主成分图像和第二主成分图像的至少两种二值化阈值结果;
阈值集成模块:用于利用集成阈值法分别对所述第一主成分图像和第二主成分图像的至少两种二值化阈值结果进行集成,并将所述第一主成分图像和第二主成分图像的阈值集成结果进行合并,得到所述遥感影像的云检测结果。
本申请实施例采取的技术方案还包括:所述云特征增强模块对第一主成分图像和第二主成分图像分别进行自适应S曲线增强具体为:设遥感影像的最大灰度值为r,最小灰度值为s,当r与s的和与2取余计算的结果等于0时,令初始阈值为T0=(r+s)/2,否则令初始阈值为T0=(r+s-1)/2,将图G的像素值G(x,y)小于T0的集合记为f1(x,y),图G的像素值G(x,y)大于等于T0的集合记为f2(x,y),S曲线拉伸的像素值F(x,y)计算公式为:
Figure GDA0002723055130000051
上述公式中,k为拉伸系数。
本申请实施例采取的技术方案还包括:所述至少两种阈值方法为十种阈值方法,所述十种阈值方法分别为:大津法、分块大津法、局部动态阈值法、全局阈值与局部阈值结合、Wellner自适应阈值、最小误差法、双峰法、迭代阈值法、最大熵阈值法和固定阈值分割法。
本申请实施例采取的技术方案还包括:所述阈值集成模块利用集成阈值法分别对所述第一主成分图像和第二主成分图像的至少两种二值化阈值结果进行集成,并将所述第一主成分图像和第二主成分图像的阈值集成结果进行合并,得到所述遥感影像的云检测结果具体为:使用投票的方式作为结合策略,投票系数δ决定结合策略的结合程度;设所述十种阈值方法得到的云检测结果为(F1,F2,......Fp),Fp为第p种阈值方法对应的云检测结果,使用δ序列(δ12,......,δi)作为投票系数,得到云检测结果
Figure GDA0002723055130000052
如果标记得到的投票数超过投票系数δi,则检测结果为云,否则为非云;最后,选择最佳的
Figure GDA0002723055130000053
值为集成阈值的云检测结果。
本申请实施例采取的又一技术方案为:一种电子设备,包括:
至少一个处理器;以及
与所述至少一个处理器通信连接的存储器;其中,
所述存储器存储有可被所述一个处理器执行的指令,所述指令被所述至少一个处理器执行,以使所述至少一个处理器能够执行上述的遥感影像云检测方法的以下操作:
步骤a:对多波段遥感影像进行主成分分析,提取所述多波段遥感影像的第一主成分图像和第二主成分图像,并对所述第一主成分图像和第二主成分图像分别进行自适应S曲线增强;
步骤b:利用至少两种阈值方法分别得到所述自适应S曲线增强后的第一主成分图像和第二主成分图像的至少两种二值化阈值结果;
步骤c:利用集成阈值法分别对所述第一主成分图像和第二主成分图像的至少两种二值化阈值结果进行集成,并将所述第一主成分图像和第二主成分图像的阈值集成结果进行合并,得到所述遥感影像的云检测结果。
相对于现有技术,本申请实施例产生的有益效果在于:本申请实施例的遥感影像云检测方法、系统及电子设备通过PCA主成分提取、S曲线云特征增强及集成阈值方法进行基于自适应S曲线增强和集成阈值的遥感影像云检测,能够有效的改善云检测中的错检和漏检情况,提高云检测效率,并提高了遥感影像的利用率。
附图说明
图1是本申请实施例的遥感影像云检测方法的流程图;
图2为基于灰度直方图的S曲线云特征值增强示意图;
图3是本申请实施例的遥感影像云检测系统的结构示意图;
图4是本申请实施例提供的遥感影像云检测方法的硬件设备结构示意图;
图5为第一实施例的风云影像云检测示意图,其中,图5(a)原图,图5(b)为手工标记,图5(c)为OTSU,图5(d)本申请的方法;
图6为第二实施例的Landsat8影像云检测示意图,其中,图6(a)原图,图6(b)为手工标记,图6(c)为OTSU,图6(d)本申请的方法;
图7为第三实施例的Sentinal-2影像云检测示意图,其中,图7(a)原图,图7(b)为手工标记,图7(c)为OTSU,图7(d)本申请的方法。
具体实施方式
为了使本申请的目的、技术方案及优点更加清楚明白,以下结合附图及实施例,对本申请进行进一步详细说明。应当理解,此处所描述的具体实施例仅用以解释本申请,并不用于限定本申请。
请参阅图1,是本申请实施例的遥感影像云检测方法的流程图。本申请实施例的遥感影像云检测方法包括以下步骤:
步骤100:通过PCA(PrincipalComponents Analysis,主成分分析)技术对多波段遥感影像进行主成分分析,提取多波段遥感影像的第一主成分图像和第二主成分图像;
步骤100中,PCA的基本思想是将原来具有一定相关性的指标X1,X2,...,Xp,重新组合成一组较少个数的互不相关的综合指标Fm来代替原来的指标。Fm既能最大程度的反映原变量Xp所代表的信息,又能保证新指标之间保持相互无关。
具体的波段提取方式为:设F1表示原变量的第一个线性组合所形成的主成分指标,即F1=a11X1+a21X2+...+ap1Xp,由数学相关知识可知,每一个主成分所提取的信息量可以用方差来度量,其方差var(F1)越大,表示F1包含的信息越多。通常第一主成分图像F1所含的信息量最大,因此在所有的线性组合中选取的F1应该是X1,X2,...,Xp的所有线性组合中方差最大的,即第一主成分图像为F1,如第一主成分图像不足以代表原来p个指标的信息,再考虑选取第二个主成分指标F2来反映原数据信息,且Cov(F1,F2)=0,即F2与F1独立、不相关。所以F2是与F1不相关的X1,X2,...,Xp的所有线性组合中方差最大的,以此类推构造出F1,F2,......,Fm为原变量指标X1,X2,...,Xp第一、第二、......、第m个主成分。
Figure GDA0002723055130000081
公式(1)中,Fi与Fj互不相关,即Cov(Fi,Fj)=0,并有Var(Fi)=ai’Σai,其中Σ为X的协方差矩阵;Fm是F1,F2,......,Fm-1都不相关的X1,X2,...,Xp所有线性组合中方差最大者;F1,F2,......,Fm(m≤p)为构造的新指标,即原变量指标的第一、第二、......、第m个主成分。
步骤200:利用基于灰度直方图的自适应S曲线对第一主成分图像和第二主成分图像中的云特征进行增强处理;
步骤200中,基于灰度直方图的S曲线云特征值增强的示意图如图2所示。设遥感影像的最大灰度值为r,最小灰度值为s,当r与s的和与2取余计算的结果等于0时,令初始阈值为T0=(r+s)/2,否则令初始阈值为T0=(r+s-1)/2。例如,图G的像素值G(x,y)小于T0的集合记为f1(x,y),图G的像素值G(x,y)大于等于T0的集合记为f2(x,y),S曲线拉伸的像素值F(x,y)计算公式如下:
Figure GDA0002723055130000082
公式(2)中,k为拉伸系数,k值越大,拉伸越明显,本申请实施例中,取k=2,具体可根据实际应用进行设定。S曲线拉伸能使得含云像素在遥感影像上更亮,非云像素在遥感影像上表现更暗,有助于遥感影像的后续处理。
步骤300:利用至少两种阈值方法分别得到第一主成分图像和第二主成分图像的二值化阈值结果;
步骤300中,由于单个阈值方法得到的云检测二值结果不够准确,采用至少两种阈值方法进行二值化,可以提高云检测二值结果的准确度。本申请实施例中,至少两种阈值方法为十种阈值方法,可以理解,阈值方法的数量可以根据实际应用进行设定。具体地,本申请实施例的十种阈值方法分别为:
1、大津法
大津法算法是按图像的灰度特性,将图像分成背景和前景两部分,分别得到前景灰度平均值和前景像素点占总像素点的概率以及背景灰度均值和背景占总像素点的概率、整幅图像的均值,计算方差。因方差是灰度分布均匀性的一种度量,背景和前景之间的类间方差越大,说明构成图像的两部分的差别越大,当部分前景错分为背景或部分背景错分为前景都会导致图像的两部分差别变小。因此,类间方差最大的分割意味着错分概率最小。
2、分块大津法
分块大津法即是将图像分成相同大小的不同块,针对每一块使用大津法阈值处理,分块大津法能够更好的保留局部特征,使得细节部分更加明显。
3、局部动态阈值法
局部动态阈值法是根据像素的邻域块的像素值分布来确定该像素位置上的二值化阈值。每个像素位置处的二值化阈值不是固定不变的,而是由其周围邻域像素的分布来决定的。亮度较高的图像区域的二值化阈值通常会较高,而亮度较低的图像区域的二值化阈值则会相适应地变小。不同亮度、对比度、纹理的局部图像区域将会拥有相对应的局部二值化阈值。
4、全局阈值与局部阈值结合
首先为全局阈值选择一个初始估计值Th(图像的平均灰度),用Th分割图像,分割完成产生两种像素:G1像素由灰度值大于Th的像素组成,G2像素由小于等于T的像素组成。分别计算G1像素的平均灰度值m1和G2像素的平均灰度值m2,并计算m1和m2的均值作为新的阈值,最后重复前述步骤,直到连续迭代中的Th值间的差为零。局部阈值采用同样的原理,将局部阈值和全局阈值结合起来能比较好的得到局部二值化结果。
5、Wellner自适应阈值
Wellner自适应阈值方法首先遍历图像,将图像的所有行假设成一个行向量,对于每个像素计算一个移动的平均值,如果某个像素明显的低于这个平均值,则设置为黑色,否则设置为白色。
假设pn为图像中位于点n处的像素,fs(n)是点n处后s个像素的总和,最后的图像是T(n)为1(白色)或0(黑色),则依赖于其是否比其前s个像素的平均值的百分之t的暗,公式如下:
Figure GDA0002723055130000101
6、最小误差法
最小误差法的思想是假设灰度图像由目标和背景组成,且目标和背景满足一混合高斯分布,计算目标和背景的均值、方差,根据最小分类误差思想得到的最小误差目标函数,取目标函数最小时的阈值即为最佳阈值。最后按此阈值将图像二值化。
7、双峰法
双峰法图像分割是一种简单的分割算法,双峰法图像二值化就是根据双峰法得到的阈值对图像进行二值化。在直方图中有两个山峰状的图像分布,山峰的顶记为Hmax1和Hmax2,他们对应的灰度值分别为T1和T2,双峰法图像分割的思想就是找到图像两个山峰之间的谷地最低值,即在[T1,T2]的灰度范围内寻找阈值Th,使其满足对应的像素数目最少,表现在图像上就是高度最低,用Th对图像进行二值化。
8、迭代阈值法
迭代法图像二值化的算法思想是首先初始化一个阈值Th,然后按照某种策略通过迭代不断更新这一阈值,直到满足给定的约束条件为止。其基本步骤如下:首先对于一幅图像,假设当前像素为f(x,y),设定一阈值Th,跟据当前阈值,循环f(x,y),将图像分为两类像素的集合A,B;其次分别计算A,B集合的像素均值μA和μB;更新阈值Th为μA和μB的均值;最后判断当前计算阈值与上次计算阈值的差值是否满足约束条件,即两次阈值差值小于一约束值Th,若小于,则当前阈值Th即为所求最佳阈值,否则,继续求取A,B的像素均值μA和μB
9、最大熵阈值法
一维最大熵法图像分割就是利用图像的灰度分布密度函数定义图像的信息熵,通过优化一定的熵准则得到熵最大时对应的阈值,从而进行图像分割的方法。算法的基本过程首先对于一幅灰度图像,灰度范围为[0,L-1],分别求取图像的最小灰度级min和最大灰度级max;其次按照熵的公式求取灰度t对应的熵值E(t);最后计算t从最小灰度min到最大灰度max之间不同灰度级所对应的熵值E(t),求取E(t)最大时所对应的灰度级t,该灰度级即为所求的阈值Th
10、固定阈值分割
固定阈值分割即人为设定一个阈值Th,设定的阈值需要根据经验,参考前文的阈值,计算前述9种方法的平均阈值并设置为固定阈值,当图像当前像素小于该固定阈值,将该像素设为0,否者设为1,人为的设定阈值需要图像的灰度分布范围。
步骤400:利用集成阈值法分别对第一主成分图像和第二主成分图像的至少两种二值化阈值结果进行集成,并将第一主成分图像和第二主成分图像的二值化阈值集成结果进行合并,得到遥感影像的云检测结果;
步骤400中,集成阈值法是基于多种阈值集成的一种云检测方法,相对于单个阈值方法,本申请通过一种结合策略将多个阈值继续集成,所得到的云检测结果可以比较好的反映出云检测结果。
在基于自适应S曲线增强和集成阈值的遥感影像云检测过程中,使用投票的方式作为结合策略。影响阈值检测的参数主要为投票系数δ,δ决定了结合策略的结合程度。前述的十种阈值方法得到的云检测结果为(F1,F2,......Fp),Fp为第p种阈值方法对应的云检测结果。为了找到更好的δ值,使用δ序列(δ12,......,δi)作为投票系数,得到云检测结果
Figure GDA0002723055130000121
如果标记得到的投票数超过投票系数δi,则检测结果为云,否则为非云。最后,选择最佳的
Figure GDA0002723055130000122
值为集成阈值的云检测结果。
对于风云影像,由于第二主成分图像无法检测出无光照区域的云,而第一主成分图像能够比较全面的检测出有光照和无光照区域的云,如果可见光和第一主成分图像的阈值集成结果对应位置上至少有一个检测为云则将该位置标记为云,否则为非云。对于非风云影像,将第一主成分图像和第二主成分图像的阈值结果进行合成,取两个波段的云检测结果的并集作为云检测结果。
请参阅图3,是本申请实施例的遥感影像云检测系统的结构示意图。本申请实施例的遥感影像云检测系统包括波段提取模块、云特征增强模块、二值化模块和阈值集成模块。
波段提取模块:用于通过PCA(PrincipalComponents Analysis,主成分分析)技术对多波段遥感影像进行主成分分析,提取多波段遥感影像的第一主成分图像和第二主成分图像;其中,PCA的基本思想是将原来具有一定相关性的指标X1,X2,...,Xp,重新组合成一组较少个数的互不相关的综合指标Fm来代替原来的指标。Fm既能最大程度的反映原变量Xp所代表的信息,又能保证新指标之间保持相互无关。
具体的波段提取方式为:设F1表示原变量的第一个线性组合所形成的主成分指标,即F1=a11X1+a21X2+...+ap1Xp,由数学相关知识可知,每一个主成分所提取的信息量可以用方差来度量,其方差var(F1)越大,表示F1包含的信息越多。通常第一主成分图像F1所含的信息量最大,因此在所有的线性组合中选取的F1应该是X1,X2,...,Xp的所有线性组合中方差最大的,即第一主成分图像为F1,如第一主成分图像不足以代表原来p个指标的信息,再考虑选取第二个主成分指标F2来反映原数据信息,且Cov(F1,F2)=0,即F2与F1独立、不相关。所以F2是与F1不相关的X1,X2,...,Xp的所有线性组合中方差最大的,以此类推构造出F1,F2,......,Fm为原变量指标X1,X2,...,Xp第一、第二、......、第m个主成分。
Figure GDA0002723055130000131
公式(1)中,Fi与Fj互不相关,即Cov(Fi,Fj)=0,并有Var(Fi)=ai’Σai,其中Σ为X的协方差矩阵;Fm是F1,F2,......,Fm-1都不相关的X1,X2,...,Xp所有线性组合中方差最大者;F1,F2,......,Fm(m≤p)为构造的新指标,即原变量指标的第一、第二、......、第m个主成分。
云特征增强模块:用于利用基于灰度直方图的自适应S曲线对第一主成分图像和第二主成分图像中的云特征进行增强处理;其中,基于灰度直方图的S曲线云特征值增强的示意图如图2所示。设遥感影像的最大灰度值为r,最小灰度值为s,当r与s的和与2取余计算的结果等于0时,令初始阈值为T0=(r+s)/2,否则令初始阈值为T0=(r+s-1)/2。例如,图G的像素值G(x,y)小于T0的集合记为f1(x,y),图G的像素值G(x,y)大于等于T0的集合记为f2(x,y),S曲线拉伸的像素值F(x,y)计算公式如下:
Figure GDA0002723055130000141
公式(2)中,k为拉伸系数,k值越大,拉伸越明显,本申请实施例中,取k=2,具体可根据实际应用进行设定。S曲线拉伸能使得含云像素在遥感影像上更亮,非云像素在遥感影像上表现更暗,有助于遥感影像的后续处理。
二值化模块:用于利用至少两种阈值方法分别得到第一主成分图像和第二主成分图像的至少两种二值化阈值结果;由于单个阈值方法得到的云检测二值结果不够准确,采用至少两种阈值方法进行二值化,可以提高云检测二值结果的准确度。本申请实施例中,至少两种阈值方法为十种阈值方法,可以理解,阈值方法的数量可以根据实际应用进行设定。具体地,本申请实施例的十种阈值方法分别为:
1、大津法
大津法算法是按图像的灰度特性,将图像分成背景和前景两部分,分别得到前景灰度平均值和前景像素点占总像素点的概率以及背景灰度均值和背景占总像素点的概率、整幅图像的均值,计算方差。因方差是灰度分布均匀性的一种度量,背景和前景之间的类间方差越大,说明构成图像的两部分的差别越大,当部分前景错分为背景或部分背景错分为前景都会导致图像的两部分差别变小。因此,类间方差最大的分割意味着错分概率最小。
2、分块大津法
分块大津法即是将图像分成相同大小的不同块,针对每一块使用大津法阈值处理,分块大津法能够更好的保留局部特征,使得细节部分更加明显。
3、局部动态阈值法
局部动态阈值法是根据像素的邻域块的像素值分布来确定该像素位置上的二值化阈值。每个像素位置处的二值化阈值不是固定不变的,而是由其周围邻域像素的分布来决定的。亮度较高的图像区域的二值化阈值通常会较高,而亮度较低的图像区域的二值化阈值则会相适应地变小。不同亮度、对比度、纹理的局部图像区域将会拥有相对应的局部二值化阈值。
4、全局阈值与局部阈值结合
首先为全局阈值选择一个初始估计值Th(图像的平均灰度),用Th分割图像,分割完成产生两种像素:G1像素由灰度值大于Th的像素组成,G2像素由小于等于T的像素组成。分别计算G1像素的平均灰度值m1和G2像素的平均灰度值m2,并计算m1和m2的均值作为新的阈值,最后重复前述步骤,直到连续迭代中的Th值间的差为零。局部阈值采用同样的原理,将局部阈值和全局阈值结合起来能比较好的得到局部二值化结果。
5、Wellner自适应阈值
Wellner自适应阈值方法首先遍历图像,将图像的所有行假设成一个行向量,对于每个像素计算一个移动的平均值,如果某个像素明显的低于这个平均值,则设置为黑色,否则设置为白色。
假设pn为图像中位于点n处的像素,fs(n)是点n处后s个像素的总和,最后的图像是T(n)为1(白色)或0(黑色),则依赖于其是否比其前s个像素的平均值的百分之t的暗,公式如下:
Figure GDA0002723055130000161
6、最小误差法
最小误差法的思想是假设灰度图像由目标和背景组成,且目标和背景满足一混合高斯分布,计算目标和背景的均值、方差,根据最小分类误差思想得到的最小误差目标函数,取目标函数最小时的阈值即为最佳阈值。最后按此阈值将图像二值化。
7、双峰法
双峰法图像分割是一种简单的分割算法,双峰法图像二值化就是根据双峰法得到的阈值对图像进行二值化。在直方图中有两个山峰状的图像分布,山峰的顶记为Hmax1和Hmax2,他们对应的灰度值分别为T1和T2,双峰法图像分割的思想就是找到图像两个山峰之间的谷地最低值,即在[T1,T2]的灰度范围内寻找阈值Th,使其满足对应的像素数目最少,表现在图像上就是高度最低,用Th对图像进行二值化。
8、迭代阈值法
迭代法图像二值化的算法思想是首先初始化一个阈值Th,然后按照某种策略通过迭代不断更新这一阈值,直到满足给定的约束条件为止。其基本步骤如下:首先对于一幅图像,假设当前像素为f(x,y),设定一阈值Th,跟据当前阈值,循环f(x,y),将图像分为两类像素的集合A,B;其次分别计算A,B集合的像素均值μA和μB;更新阈值Th为μA和μB的均值;最后判断当前计算阈值与上次计算阈值的差值是否满足约束条件,即两次阈值差值小于一约束值Th,若小于,则当前阈值Th即为所求最佳阈值,否则,继续求取A,B的像素均值μA和μB
9、最大熵阈值法
一维最大熵法图像分割就是利用图像的灰度分布密度函数定义图像的信息熵,通过优化一定的熵准则得到熵最大时对应的阈值,从而进行图像分割的方法。算法的基本过程首先对于一幅灰度图像,灰度范围为[0,L-1],分别求取图像的最小灰度级min和最大灰度级max;其次按照熵的公式求取灰度t对应的熵值E(t);最后计算t从最小灰度min到最大灰度max之间不同灰度级所对应的熵值E(t),求取E(t)最大时所对应的灰度级t,该灰度级即为所求的阈值Th
10、固定阈值分割
固定阈值分割即人为设定一个阈值Th,设定的阈值需要根据经验,参考前文的阈值,计算前述9种方法的平均阈值并设置为固定阈值,当图像当前像素小于该固定阈值,将该像素设为0,否者设为1,人为的设定阈值需要图像的灰度分布范围。
阈值集成模块:用于利用集成阈值法分别对第一主成分图像和第二主成分图像的至少两种二值化阈值结果进行集成,并将第一主成分图像和第二主成分图像的二值化阈值集成结果进行合并,得到遥感影像的云检测结果;其中,集成阈值法是基于多种阈值集成的一种云检测方法,主要思想是单个阈值方法得到的云检测二值结果不够准确,由多个阈值通过一种结合策略得到的云检测结果可以比较好的反映出云检测结果。
在遥感影像云检测过程中,使用投票的方式作为结合策略。影响阈值检测的参数主要为投票系数δ,δ决定了结合策略的结合程度。前述的十种阈值方法得到的云检测结果为(F1,F2,......Fp),Fp为第p种阈值方法对应的云检测结果。为了找到更好的δ值,使用δ序列(δ12,......,δi)作为投票系数,得到云检测结果
Figure GDA0002723055130000181
如果标记得到的投票数超过投票系数δi,则检测结果为云,否则为非云。最后,选择最佳的
Figure GDA0002723055130000182
值为集成阈值的云检测结果。
对于风云影像,由于第二主成分图像无法检测出无光照区域的云,而第一主成分图像能够比较全面的检测出有光照和无光照区域的云,如果可见光和第一主成分图像的阈值集成结果对应位置上至少有一个检测为云则将该位置标记为云,否则为非云。对于非风云影像,将第一主成分图像和第二主成分图像的阈值结果进行合成,取两个波段的云检测结果的并集作为云检测结果。
图4是本申请实施例提供的遥感影像云检测方法的硬件设备结构示意图。如图4所示,该设备包括一个或多个处理器以及存储器。以一个处理器为例,该设备还可以包括:输入系统和输出系统。
处理器、存储器、输入系统和输出系统可以通过总线或者其他方式连接,图4中以通过总线连接为例。
存储器作为一种非暂态计算机可读存储介质,可用于存储非暂态软件程序、非暂态计算机可执行程序以及模块。处理器通过运行存储在存储器中的非暂态软件程序、指令以及模块,从而执行电子设备的各种功能应用以及数据处理,即实现上述方法实施例的处理方法。
存储器可以包括存储程序区和存储数据区,其中,存储程序区可存储操作系统、至少一个功能所需要的应用程序;存储数据区可存储数据等。此外,存储器可以包括高速随机存取存储器,还可以包括非暂态存储器,例如至少一个磁盘存储器件、闪存器件、或其他非暂态固态存储器件。在一些实施例中,存储器可选包括相对于处理器远程设置的存储器,这些远程存储器可以通过网络连接至处理系统。上述网络的实例包括但不限于互联网、企业内部网、局域网、移动通信网及其组合。
输入系统可接收输入的数字或字符信息,以及产生信号输入。输出系统可包括显示屏等显示设备。
所述一个或者多个模块存储在所述存储器中,当被所述一个或者多个处理器执行时,执行上述任一方法实施例的以下操作:
步骤a:对多波段遥感影像进行主成分分析,提取所述多波段遥感影像的第一主成分图像和第二主成分图像,并对所述第一主成分图像和第二主成分图像分别进行自适应S曲线增强;
步骤b:利用至少两种阈值方法分别得到所述自适应S曲线增强后的第一主成分图像和第二主成分图像的至少两种二值化阈值结果;
步骤c:利用集成阈值法分别对所述第一主成分图像和第二主成分图像的至少两种二值化阈值结果进行集成,并将所述第一主成分图像和第二主成分图像的阈值集成结果进行合并,得到所述遥感影像的云检测结果。
上述产品可执行本申请实施例所提供的方法,具备执行方法相应的功能模块和有益效果。未在本实施例中详尽描述的技术细节,可参见本申请实施例提供的方法。
本申请实施例提供了一种非暂态(非易失性)计算机存储介质,所述计算机存储介质存储有计算机可执行指令,该计算机可执行指令可执行以下操作:
步骤a:对多波段遥感影像进行主成分分析,提取所述多波段遥感影像的第一主成分图像和第二主成分图像,并对所述第一主成分图像和第二主成分图像分别进行自适应S曲线增强;
步骤b:利用至少两种阈值方法分别得到所述自适应S曲线增强后的第一主成分图像和第二主成分图像的至少两种二值化阈值结果;
步骤c:利用集成阈值法分别对所述第一主成分图像和第二主成分图像的至少两种二值化阈值结果进行集成,并将所述第一主成分图像和第二主成分图像的阈值集成结果进行合并,得到所述遥感影像的云检测结果。
本申请实施例提供了一种计算机程序产品,所述计算机程序产品包括存储在非暂态计算机可读存储介质上的计算机程序,所述计算机程序包括程序指令,当所述程序指令被计算机执行时,使所述计算机执行以下操作:
步骤a:对多波段遥感影像进行主成分分析,提取所述多波段遥感影像的第一主成分图像和第二主成分图像,并对所述第一主成分图像和第二主成分图像分别进行自适应S曲线增强;
步骤b:利用至少两种阈值方法分别得到所述自适应S曲线增强后的第一主成分图像和第二主成分图像的至少两种二值化阈值结果;
步骤c:利用集成阈值法分别对所述第一主成分图像和第二主成分图像的至少两种二值化阈值结果进行集成,并将所述第一主成分图像和第二主成分图像的阈值集成结果进行合并,得到所述遥感影像的云检测结果。
为了验证本申请方法的可行性和有效性,以下实施例利用多幅包含云的遥感影像进行了测试,同时与手工标记云检测以及OTSU(最大类间方差法)进行了比较,具体如图5、图6和图7所示,图5为第一实施例的风云影像云检测示意图,其中,图5(a)原图,图5(b)为手工标记,图5(c)为OTSU,图5(d)本申请的方法。图6为第二实施例的Landsat8影像云检测示意图,其中,图6(a)原图,图6(b)为手工标记,图6(c)为OTSU,图6(d)本申请的方法。图7为第三实施例的Sentinal-2影像云检测示意图,其中,图7(a)原图,图7(b)为手工标记,图7(c)为OTSU,图7(d)本申请的方法。为客观评价各算法的结果图像的对比度,采用命中率(Probability of Detection,POD)、误报率(False Alarm Ratio,FAR)和临界成功指数(Critical Success Index,CSI)来进行评价,即:
Figure GDA0002723055130000211
公式(4)中,NH表示FY-2G的云检测结果和手工标记云检测都为云的像素点数,NM表示FY-2G云检测结果中无云而手工标记云检测有云的像素点数,NF表示FY-2G云检测结果中有云而手工标记云检测中无云的像素点数。POD越高,则检测的正确率越高;FAR越低则,检测的错误率越低;临界成功指数能够反映检测结果接近真实值的综合度量。
评价结果如下表1所示:
表1
Figure GDA0002723055130000212
由表1可知单一的OTSU方法对于云检测不一定都能得到较好的结果,本申请的POD值总是最高的,效果最好。即表明本申请能有效的检测出云,且各个指标的平均值表明,本申请具有最高的POD、CSI值,最低的FAR值。各个指标的标准差旨在评估在不同情况下方法的稳定性,结果表明本申请具有更稳定的云检测结果。
本申请实施例的遥感影像云检测方法、系统及电子设备通过PCA主成分提取、S曲线云特征增强及集成阈值方法进行基于自适应S曲线增强和集成阈值的遥感影像云检测,适用于所有类型的遥感影像,能够有效的改善云检测中的错检和漏检情况,提高云检测效率,并提高了遥感影像的利用率。
对所公开的实施例的上述说明,使本领域专业技术人员能够实现或使用本申请。对这些实施例的多种修改对本领域的专业技术人员来说将是显而易见的,本申请中所定义的一般原理可以在不脱离本申请的精神或范围的情况下,在其它实施例中实现。因此,本申请将不会被限制于本申请所示的这些实施例,而是要符合与本申请所公开的原理和新颖特点相一致的最宽的范围。

Claims (7)

1.一种遥感影像云检测方法,其特征在于,包括以下步骤:
步骤a:对多波段遥感影像进行主成分分析,提取所述多波段遥感影像的第一主成分图像和第二主成分图像,并对所述第一主成分图像和第二主成分图像分别进行自适应S曲线增强;
步骤b:利用至少两种阈值方法分别得到所述自适应S曲线增强后的第一主成分图像和第二主成分图像的至少两种二值化阈值结果;
步骤c:利用集成阈值法分别对所述第一主成分图像和第二主成分图像的至少两种二值化阈值结果进行集成,并将所述第一主成分图像和第二主成分图像的阈值集成结果进行合并,得到所述遥感影像的云检测结果;
在所述步骤a中,所述对第一主成分图像和第二主成分图像分别进行自适应S曲线增强具体为:设遥感影像的最大灰度值为r,最小灰度值为s,当r与s的和与2取余计算的结果等于0时,令初始阈值为T0=(r+s)/2,否则令初始阈值为T0=(r+s-1)/2,将图G的像素值G(x,y)小于T0的集合记为f1(x,y),图G的像素值G(x,y)大于等于T0的集合记为f2(x,y),S曲线拉伸的像素值F(x,y)计算公式为:
Figure FDA0002991946400000011
上述公式中,k为拉伸系数。
2.根据权利要求1所述的遥感影像云检测方法,其特征在于,在所述步骤b中,所述至少两种阈值方法为十种阈值方法,所述十种阈值方法分别为:大津法、分块大津法、局部动态阈值法、全局阈值与局部阈值结合、Wellner自适应阈值、最小误差法、双峰法、迭代阈值法、最大熵阈值法和固定阈值分割法。
3.根据权利要求2所述的遥感影像云检测方法,其特征在于,在所述步骤c中,所述利用集成阈值法分别对所述第一主成分图像和第二主成分图像的至少两种二值化阈值结果进行集成,并将所述第一主成分图像和第二主成分图像的阈值集成结果进行合并,得到所述遥感影像的云检测结果具体为:使用投票的方式作为结合策略,投票系数δ决定结合策略的结合程度;设所述十种阈值方法得到的云检测结果为(F1,F2,......Fp),Fp为第p种阈值方法对应的云检测结果,使用δ序列(δ12,......,δi)作为投票系数,得到云检测结果
Figure FDA0002991946400000021
如果每个像素标记得到的投票数超过投票系数δi,则检测结果为云,否则为非云;最后,选择最佳的
Figure FDA0002991946400000022
值为集成阈值的云检测结果。
4.一种遥感影像云检测系统,其特征在于,包括:
波段提取模块:用于对多波段遥感影像进行主成分分析,提取所述多波段遥感影像的第一主成分图像和第二主成分图像;
云特征增强模块:用于对所述第一主成分图像和第二主成分图像分别进行自适应S曲线增强;
二值化模块:用于利用至少两种阈值方法分别得到所述自适应S曲线增强后的第一主成分图像和第二主成分图像的至少两种二值化阈值结果;
阈值集成模块:用于利用集成阈值法分别对所述第一主成分图像和第二主成分图像的至少两种二值化阈值结果进行集成,并将所述第一主成分图像和第二主成分图像的阈值集成结果进行合并,得到所述遥感影像的云检测结果;
所述云特征增强模块对第一主成分图像和第二主成分图像分别进行自适应S曲线增强具体为:设遥感影像的最大灰度值为r,最小灰度值为s,当r与s的和与2取余计算的结果等于0时,令初始阈值为T0=(r+s)/2,否则令初始阈值为T0=(r+s-1)/2,将图G的像素值G(x,y)小于T0的集合记为f1(x,y),图G的像素值G(x,y)大于等于T0的集合记为f2(x,y),S曲线拉伸的像素值F(x,y)计算公式为:
Figure FDA0002991946400000031
上述公式中,k为拉伸系数。
5.根据权利要求4所述的遥感影像云检测系统,其特征在于,所述至少两种阈值方法为十种阈值方法,所述十种阈值方法分别为:大津法、分块大津法、局部动态阈值法、全局阈值与局部阈值结合、Wellner自适应阈值、最小误差法、双峰法、迭代阈值法、最大熵阈值法和固定阈值分割法。
6.根据权利要求5所述的遥感影像云检测系统,其特征在于,所述阈值集成模块利用集成阈值法分别对所述第一主成分图像和第二主成分图像的至少两种二值化阈值结果进行集成,并将所述第一主成分图像和第二主成分图像的阈值集成结果进行合并,得到所述遥感影像的云检测结果具体为:使用投票的方式作为结合策略,投票系数δ决定结合策略的结合程度;设所述十种阈值方法得到的云检测结果为(F1,F2,......Fp),Fp为第p种阈值方法对应的云检测结果,使用δ序列(δ12,......,δi)作为投票系数,得到云检测结果
Figure FDA0002991946400000032
如果标记得到的投票数超过投票系数δi,则检测结果为云,否则为非云;最后,选择最佳的
Figure FDA0002991946400000033
值为集成阈值的云检测结果。
7.一种电子设备,包括:
至少一个处理器;以及
与所述至少一个处理器通信连接的存储器;其中,
所述存储器存储有可被所述一个处理器执行的指令,所述指令被所述至少一个处理器执行,以使所述至少一个处理器能够执行上述1至3任一项所述的遥感影像云检测方法的以下操作:
步骤a:对多波段遥感影像进行主成分分析,提取所述多波段遥感影像的第一主成分图像和第二主成分图像,并对所述第一主成分图像和第二主成分图像分别进行自适应S曲线增强;
步骤b:利用至少两种阈值方法分别得到所述自适应S曲线增强后的第一主成分图像和第二主成分图像的至少两种二值化阈值结果;
步骤c:利用集成阈值法分别对所述第一主成分图像和第二主成分图像的至少两种二值化阈值结果进行集成,并将所述第一主成分图像和第二主成分图像的阈值集成结果进行合并,得到所述遥感影像的云检测结果。
CN201811537640.0A 2018-12-15 2018-12-15 遥感影像云检测方法、系统及电子设备 Active CN109726649B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201811537640.0A CN109726649B (zh) 2018-12-15 2018-12-15 遥感影像云检测方法、系统及电子设备

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201811537640.0A CN109726649B (zh) 2018-12-15 2018-12-15 遥感影像云检测方法、系统及电子设备

Publications (2)

Publication Number Publication Date
CN109726649A CN109726649A (zh) 2019-05-07
CN109726649B true CN109726649B (zh) 2021-08-24

Family

ID=66297598

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201811537640.0A Active CN109726649B (zh) 2018-12-15 2018-12-15 遥感影像云检测方法、系统及电子设备

Country Status (1)

Country Link
CN (1) CN109726649B (zh)

Families Citing this family (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN110096454B (zh) * 2019-05-15 2021-06-01 武昌理工学院 一种基于非易失存储器的遥感数据快存方法
CN110796677B (zh) * 2019-10-29 2022-10-21 北京环境特性研究所 一种基于多波段特性的卷云虚警源检测方法
CN110942452B (zh) * 2019-11-21 2020-11-06 自然资源部国土卫星遥感应用中心 一种基于多时相热红外遥感影像的云检测方法
CN112634382B (zh) * 2020-11-27 2024-03-19 国家电网有限公司大数据中心 一种非自然对象的图像识别、替换方法和装置
CN113298839B (zh) * 2021-05-18 2023-08-15 中国矿业大学 一种sar影像半自动阈值提取方法

Family Cites Families (9)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US7505604B2 (en) * 2002-05-20 2009-03-17 Simmonds Precision Prodcuts, Inc. Method for detection and recognition of fog presence within an aircraft compartment using video images
US7414706B2 (en) * 2004-12-22 2008-08-19 Northrop Grumman Corporation Method and apparatus for imaging a target using cloud obscuration prediction and detection
GB2537681B (en) * 2015-04-24 2018-04-25 Univ Oxford Innovation Ltd A method of detecting objects within a 3D environment
CN104933415B (zh) * 2015-06-25 2017-12-29 华中科技大学 一种实时可见光遥感影像云区检测方法
CN105354865B (zh) * 2015-10-27 2018-01-26 武汉大学 多光谱遥感卫星影像自动云检测方法及系统
CN107563397A (zh) * 2016-06-30 2018-01-09 中国电力科学研究院 一种卫星云图中云团运动矢量的确定方法
CN106295498B (zh) * 2016-07-20 2019-04-16 湖南大学 光学遥感图像目标区域检测装置与方法
CN106909902B (zh) * 2017-03-01 2020-06-05 北京航空航天大学 一种基于改进的层次化显著模型的遥感目标检测方法
CN108648184A (zh) * 2018-05-10 2018-10-12 电子科技大学 一种遥感图像高空卷云的检测方法

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
一种基于集成学习技术的图像分割算法的研究;罗会兰 等;《江西理工大学学报》;20120630;第33卷(第3期);第56-61页 *
用双通道动态阈值对GMS-5图像进行自动云检测;刘希 等;《应用气象学报》;20050831;第16卷(第4期);434-435 *

Also Published As

Publication number Publication date
CN109726649A (zh) 2019-05-07

Similar Documents

Publication Publication Date Title
CN109740639B (zh) 一种风云卫星遥感影像云检测方法、系统及电子设备
CN109726649B (zh) 遥感影像云检测方法、系统及电子设备
CN107229918B (zh) 一种基于全卷积神经网络的sar图像目标检测方法
CN109902715B (zh) 一种基于上下文聚合网络的红外弱小目标检测方法
CN106650812B (zh) 一种卫星遥感影像的城市水体提取方法
CN107808138B (zh) 一种基于FasterR-CNN的通信信号识别方法
CN107038416B (zh) 一种基于二值图像改进型hog特征的行人检测方法
CN110569782A (zh) 一种基于深度学习目标检测方法
CN105117736B (zh) 基于稀疏深度堆栈网络的极化sar图像分类方法
CN108399430B (zh) 一种基于超像素和随机森林的sar图像舰船目标检测方法
CN110889843B (zh) 基于最大稳定极值区域的sar图像舰船目标检测方法
CN108829711B (zh) 一种基于多特征融合的图像检索方法
WO2021118463A1 (en) Defect detection in image space
CN103366373B (zh) 基于模糊相容图的多时相遥感影像变化检测方法
CN107392887A (zh) 一种基于同质像素点转化的异质遥感图像变化检测方法
CN116168240A (zh) 基于注意力增强的任意方向密集舰船目标检测方法
CN116664565A (zh) 一种光伏太阳能电池片的隐裂检测方法及系统
CN110689060B (zh) 一种基于聚合特征差异学习网络的异源图像匹配方法
CN111091580B (zh) 一种基于改进ResNet-UNet网络的立木图像分割方法
CN110310263B (zh) 一种基于显著性分析和背景先验的sar图像居民区检测方法
CN109241865B (zh) 一种弱对比度交通场景下的车辆检测分割算法
CN113284066B (zh) 遥感影像自动云检测方法和装置
CN115033721A (zh) 基于大数据的图像检索方法
CN112233042A (zh) 一种含非合作目标的大场景sar图像快速生成方法
CN112215104A (zh) 一种基于超像素分割的海冰提取方法、装置和设备

Legal Events

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