CN111652826B - 基于Wallis滤波+直方图匹配的有重多/高光谱遥感图像匀色方法 - Google Patents

Abstract

基于Wallis滤波+直方图匹配的有重多/高光谱遥感图像匀色方法，本发明涉及遥感图像匀色方法。本发明的目的是为了解决现有图像的匀色方法只注重了整体效果而忽略了图像局部特性，处理后可能使待匀色图像数据出现大的偏差，尤其是在图像交叠区域出现大的偏差的问题。过程为：一、获得M组同波段灰度图像；二、选取参考图像；三、获取新的待匀色图像；四、对参考和新的待匀色图像做直方图匹配全局匀色；五、将匀色后的图像作为参考图像，重复执行三至四，直至第m组同波段灰度图像全部匀色处理完，进行拼接；六、重复二至五获得所有M组拼接好的灰度图像，将所有M组拼接好的灰度图像合成新多/高光谱遥感图像。本发明用于遥感图像匀色领域。

Description

基于Wallis滤波+直方图匹配的有重多/高光谱遥感图像匀色方法

技术领域

本发明涉及多/高光谱遥感影像图像匀色方法。

背景技术

无论是通过何种手段或者方式获取得到的遥感数据，均会存在一些外在的因素比如传感器因素、人为因素以及天气条件等原因，导致不能将获取的遥感影像被直接广泛应用，这些影像需要科研人员经过一定的处理校正，比如图像的几何或者辐射校正、图像的匀色处理等，影像只有经过这样的后期处理才能被应用到各大领域中。而图像的匀色这一处理需要研究人员经过长期不懈地坚持努力，才能让前期获取的遥感影像达到视觉效果上良好、影像的色彩以及亮度尽可能跟地物本身相似、影像中的地物纹理清晰等目的。

近年来，图像的匀色算法逐渐地引起了国内学者的高度重视。在建立大范围的空间影像数据库时，要求同一场景的影像经过拼接后是一幅影像之间基本没有色彩和亮度偏差的无缝的大影像。现有的软件INPHO是在过渡区域对拼接缝的进行处理来消除影像之间的色彩差异，但是如果在影像的整体色彩不一样的情况下，拼接的结果也就难以达到预期的效果，从而无法保证一幅大影像整体色彩的一致性。对此，我们可以采用先对影像进行匀色处理，比较常见的、简单快速的影像色彩一致性处理方法主要有：基于直方图的图像匀色处理及改进方法、基于信息熵的图像调整、基于均值方差的图像调整等方法。

但上述方法只注重了整体效果而忽略了图像局部特性，处理后可能使待匀色图像数据出现大的偏差，尤其是在图像交叠区域。为此提出一种Wallis滤波作为局部匀色的补充算法，针对有重叠区域的多/高光谱遥感图像匀色。

发明内容

本发明的目的是为了解决现有图像的匀色方法只注重了整体效果而忽略了图像局部特性，处理后可能使待匀色图像数据出现大的偏差，尤其是在图像交叠区域出现大的偏差的问题，而提出基于Wallis滤波+直方图匹配的有重多/高光谱遥感图像匀色方法。

基于Wallis滤波+直方图匹配的有重多/高光谱遥感图像匀色方法具体过程为：

步骤一、获得多/高光谱遥感图像，对多/高光谱遥感图像进行分波段处理，获得M组同波段灰度图像，第m组同波段灰度图像为X1、X2、X3、…、Xi、…、Xn；m＝1,2,...,M；i＝1,2,...,n；

步骤二、通过统计第m组同波段灰度图像X1、X2、X3、…、Xi、…、Xn的信息量及地物丰富度，选取出参考图像作为参考；

步骤三、取两幅图像，一幅是步骤二得到的第m组同波段灰度图像X1、X2、X3、…、Xi、…、Xn的参考图像Xi，另一幅是与参考图像有交叠的同波段灰度图像Xj；

对交叠区域进行分块处理，获得子块，设置子块取整数的平方；

将与参考图像有交叠的同波段灰度图像和参考图像的交叠区域做局部匀色；

用局部匀色图像替代原待匀色图像与参考图像的交叠区域，替换交叠区域后的原待匀色图像作为新的待匀色图像；

步骤四、对步骤三的参考图像和新的待匀色图像做直方图匹配全局匀色；

步骤五、将匀色处理后的图像作为参考图像，重复执行步骤三至步骤四，直至第m组同波段灰度图像全部匀色处理完，将匀色处理完的第m组同波段灰度图像拼接；

步骤六、重复步骤二至步骤五获得所有M组拼接好的灰度图像，将所有M组拼接好的灰度图像合成新多/高光谱遥感图像。

本发明的有益效果为：

本发明实现了利用选取好的参考图像对待匀色图像进行匀色处理的改进算法，该算法结合了Wallis滤波的局部匀色效果与直方图匹配的全局匀色能力。首先利用Wallis滤波的局部处理优势应对有重叠区域的多/高光谱遥感图像匀色。提出了一种由局部到整体的一整套匀色方案，有效的处理了有重多/高光谱遥感图像的色彩和亮度不均问题。该方法受到Wallis滤波算法在增强影像对比度的启发，通过将待匀色图像的灰度均值和标准差映射到参考图像的均值和标准差，再不断调整图像标准差扩展系数与图像的亮度扩展系数使达到局部匀色最佳效果，这样不仅通过Wallis滤波增强了图像对比度、对噪声进行了一定程度的抑制，还为直方图匹配提供了更有利的数据基础，提高全局匹配性能和匀色效果。解决现有图像的匀色方法只注重了整体效果而忽略了图像局部特性，处理后可能使待匀色图像数据出现大的偏差，尤其是在图像交叠区域出现大的偏差的问题。

为了验证本发明所提出的方法的性能，针对一组有重多/高光谱遥感图像数据进行验证，GF2_160628与GF2_150906数据为高分2号采集的山东某村庄多光谱遥感图像，每张图像大小为1536×1536像素，有四个波段。实验结果验证了本发明提出的基于Wallis滤波+直方图匹配的有重区域多/高光谱遥感图像匀色算法的有效性。在使用评判标准平均梯度法的情况下，本发明的平均梯度指数达到1.4154明显高于传统方法1.0912。

附图说明

图1是传统匀色法的实现流程示意图；

图2是本发明的实现流程示意图；

图3是一组原始拼接图像图，其中左上图为待匀色图像，右下图为参考图像；

图4是一组传统匀色法结果图，其中左上图为直方图匹配匀色图像，右下图为参考图像；

图5是一组改进算法结果图，其中左上图为Wallis滤波+直方图匹配匀色图像，右下图为参考图像；

图6a表示待匀色图像重叠区域的直方图对比评价图；

图6b表示直方图匹配匀色图像重叠区域的直方图对比评价图；

图6c表示Wallis滤波+直方图匹配匀色图像重叠区域的直方图对比评价图；

图6d表示参考图像重叠区域的直方图对比评价图。

具体实施方式

具体实施方式一：本实施方式基于Wallis滤波+直方图匹配的有重多/高光谱遥感图像匀色方法具体过程为：

本发明以Wallis滤波与传统直方图匹配算法为基础，从局部优化到整体匀色的角度出发，充分利用参考图像的信息丰富度，对待匀色图像进行匀色处理。本发明的目的是为了使得处理后的遥感影像的在整体色调上保持一致，整幅图像亮度均匀，对比度适合，达到增强影像视觉感官上好的效果，满足研究人员和后期处理人员对地面物题的解译和分析，提出了一种基于Wallis滤波+直方图匹配的有重区域多/高光谱遥感图像匀色方法。

步骤一、获得多/高光谱遥感图像，对多/高光谱遥感图像进行分波段处理(鉴于多/高光谱遥感图像的多波段特性，对多/高光谱遥感图像进行分波段处理，M个波段，就分M组，每组为一个波段)，获得M组同波段灰度图像(一个波段是一个灰度图)，第m组同波段灰度图像为X1、X2、X3、…、Xi、…、Xn；m＝1,2,...,M；i＝1,2,...,n；

步骤二、通过统计第m组同波段灰度图像X1、X2、X3、…、Xi、…、Xn的信息量及地物丰富度(地物种类)，选取出参考图像作为参考(参考价值Q为最大值的图像)；

步骤三、取两幅图像，一幅是步骤二得到的第m组同波段灰度图像X1、X2、X3、…、Xi、…、Xn的参考图像Xi，另一幅是与参考图像有交叠的同波段灰度图像Xj(X1、X2、X3、…、Xi、…、Xn中的一个)；

对交叠区域进行分块处理，获得子块，设置子块取整数的平方为优，例如2*2、3*3……n*n，使得每个子块的面积保持在50*50～100*100像素范围内，子块过多会导致分区明显，过少会形容模糊；

将与参考图像有交叠的同波段灰度图像和参考图像的交叠区域做局部匀色(参考图像Xi和与参考图像有交叠的同波段灰度图像Xj的交叠区域做局部匀色，其中Xi为参考图像，Xj是局部待匀色图像，处理后的局部匀色结果图像也可作为新的参考图像)；

用局部匀色图像替代原待匀色图像与参考图像的交叠区域，替换交叠区域后的原待匀色图像作为新的待匀色图像；

步骤四、对步骤三的参考图像和新的待匀色图像做直方图匹配全局匀色；

步骤五、将匀色处理后的图像作为参考图像，重复执行步骤三至步骤四，直至第m组同波段灰度图像全部匀色处理完，将匀色处理完的第m组同波段灰度图像拼接；

步骤六、重复步骤二至步骤五获得所有M组拼接好的灰度图像，将所有M组拼接好的灰度图像合成新多/高光谱遥感图像。

具体实施方式二：本实施方式与具体实施方式一不同的是：所述步骤二中通过统计第m组同波段灰度图像X1、X2、X3、…、Xi、…、Xn的信息量及地物丰富度(地物种类)，选取出参考图像作为参考(参考价值Q为最大值的图像)；具体过程为：

步骤二一、信息熵表示第m组同波段灰度图像X1、X2、X3、…、Xi、…、Xn包含信息量的多少；

对于第m组同波段灰度图像X1、X2、X3、…、Xi、…、Xn，认为第m组同波段灰度图像X1、X2、X3、…、Xi、…、Xn每幅图像中各灰度值是相互独立的样本，单幅图像(第m组同波段灰度图像X1、X2、X3、…、Xi、…、Xn中一幅)中各灰度值所占比例为p＝{p₁、p₂、p₃…p_n}；

利用以下公式分别计算第m组同波段灰度图像X1、X2、X3、…、Xi、…、Xn中每幅图像的信息熵：

式中，i表示第m组同波段灰度图像Xi的某灰度级；p_i表示第m组同波段灰度图像Xi的对应灰度级出现的概率；n表示第m组同波段灰度图像Xi的灰度级个数；

步骤二二、按各占百分之五十的权重统计第m组同波段灰度图像X1、X2、X3、…、Xi、…、Xn中每幅图像的参考价值，计算公式如下

式中，Q为参考价值，H为第m组同波段灰度图像Xi的信息熵，H_max为第m组同波段灰度图像X1、X2、X3、…、Xi、…、Xn中信息熵中最大值，N为第m组同波段灰度图像Xi的包含地物种类(农田、房屋、河流等)，N_max为第m组同波段灰度图像X1、X2、X3、…、Xi、…、Xn中包含地物种类最大值；

参考价值Q取最大值时所指的灰度图像Xi为第m组同波段灰度图像X1、X2、X3、…、Xi、…、Xn的参考图像；

步骤二都是在同波段图像里操作；

其它步骤及参数与具体实施方式一相同。

具体实施方式三：本实施方式与具体实施方式一或二不同的是：所述步骤三中取两幅图像，一幅是步骤二得到的第m组同波段灰度图像X1、X2、X3、…、Xi、…、Xn的参考图像Xi，另一幅是与参考图像有交叠的同波段灰度图像Xj(X1、X2、X3、…、Xi、…、Xn中的一个)；

对交叠区域进行分块处理，获得子块，设置子块取整数的平方为优，例如2*2、3*3……n*n，使得每个子块的面积保持在50*50～100*100像素范围内，子块过多会导致分区明显，过少会形容模糊；

将与参考图像有交叠的同波段灰度图像和参考图像的交叠区域做局部匀色(参考图像Xi和与参考图像有交叠的同波段灰度图像Xj的交叠区域做局部匀色，其中Xi为参考图像，Xj是局部待匀色图像，处理后的局部匀色结果图像也可作为新的参考图像)；

用局部匀色图像替代原待匀色图像与参考图像的交叠区域，替换交叠区域后的原待匀色图像作为新的待匀色图像；具体过程为：

对参考图像和与参考图像有交叠的同波段灰度图像的交叠区域分别进行分块处理(分块数一致)，获得子块，设置子块为整数的平方为优，(例如2*2、3*3……n*n，)使得每个子块的面积保持在50*50～100*100像素范围内，子块过多会导致分区明显，过少会形容模糊；(比如参考图像和与参考图像有交叠的同波段灰度图像的交叠区域500×500，所以各划成5×5个子块)；

交叠区域内的参考图像每个子块对同波段灰度图像的每个子块进行Wallis滤波，完成交叠区域内的参考图像每个子块(参考图像)对同波段灰度图像每个子块(待匀色图像)(同波段灰度图像为与参考图像有交叠的交叠区域内的同波段灰度图像)的局部匀色(子块一对一匀色)；过程为：

式中，c为图像标准差扩展系数，其取值范围为0≤c≤1；b为图像的亮度扩展系数，其取值范围为0≤b≤1；当b取值靠近1时，待局部匀色图像的均值将向m_f靠近，当b取值靠近0时，待局部匀色图像的均值将向m_g靠近；m_f为参考图像第a个子块的灰度均值，m_g为待局部匀色图像第a个子块的灰度均值，s_f为参考图像第a个子块的灰度标准差，s_g为待局部匀色图像第a个子块的灰度标准差，g(x,y)为待局部匀色图像的第a个子块，g′(x,y)为局部匀色图像的第a个子块；不同子块的m_f、m_g、s_f、s_g值可能不同；

根据局部匀色图像的每个子块g′(x,y)得到交叠区域内与参考图像有交叠的同波段灰度图像的局部匀色图像；

用局部匀色图像替代原待匀色图像与参考图像的交叠区域，替换交叠区域后的原待匀色图像作为新的待匀色图像。

其它步骤及参数与具体实施方式一或二相同。

具体实施方式四：本实施方式与具体实施方式一至三之一不同的是，所述c取值为0.7～0.8，b取值为0.6～0.7在实际操作过程中效果更理想。

其它步骤及参数与具体实施方式一至三之一相同。

具体实施方式五：本实施方式与具体实施方式一至四之一不同的是，所述步骤四中对步骤三的参考图像和新的待匀色图像做直方图匹配全局匀色；具体过程为：

步骤四一、分别对新的待匀色图像和参考图像的交叠部分进行直方图均衡化；过程为：

设新的待匀色图像的(x,y)处的灰度值为h(x,y)，对新的待匀色图像进行归一化处理，得到归一化处理后新的待匀色图像r：

式中，h_min为h(x,y)中的最小值，h_max为h(x,y)中的最大值；

设参考图像的(x′,y′)处的灰度值为h(x′,y′)，对参考图像进行归一化处理，得到归一化处理后参考图像z：

式中，h′_min为h(x′,y′)中的最小值，h′_max为h(x′,y′)中的最大值；

设归一化处理后新的待匀色图像r的直方图分布为P_r(r)，均衡化结果为P_s(s)，其中：

s＝T[r]

式中，T[r]为P_r(r)的累积分布函数；

同理设归一化处理后参考图像z的直方图分布为P_z(z)，均衡化结果为P_v(v)，其中：

v＝G[z]

式中，G[z]为P_z(z)的累积分布函数；

步骤四二、根据均衡化结果的一致性，认为

s＝v

由公式s＝T[r]、v＝G[z]和s＝v，得到一种建立在r与z之间的映射关系：

z＝G^-1[T[r]]

根据这种映射算法将新的待匀色图像进行匀色处理。

其它步骤及参数与具体实施方式一至四之一相同。

采用以下实施例验证本发明的有益效果：

实施例一：

本实施例一种基于目标增强的高空间分辨率遥感图像场景分类算法具体是按照以下步骤制备的：

实验所用数据是一组有重多/高光谱遥感图像数据：GF2_160628与GF2_150906数据为高分2号采集的山东某村庄多光谱遥感图像，每张图像大小为1536×1536像素，有四个波段。附图3为原始拼接图像，可见色彩与亮度差异大，拼接交界处有明显边界；附图4为采用传统方法即直方图匹配算法实现的匀色拼接效果图，与原始图像相较有一定改善，但拼接重合处依然边界明显，例如农田色彩与参考图像依然相去甚远，过度极不自然；附图5为本发明算法，该算法融合了Wallis滤波算法与直方图匹配算法的优势，统筹局部与全局匀色，实现了较好的匀色效果，同时也消除了边界效应。附图6a、6b、6c、6d是一组直方图对比，对比对象是各图像重叠区域的灰度值，可见在数据上本发明方法同样具有更贴近于参考图像的性能优势。

表1为一组评判标准为平均梯度法的情况下各图像对比，本发明的平均梯度指数达到1.4154明显高于传统方法1.0912，体现了信息量的丰富度，在后期的图像处理等实际应用中给工程人员带来更大的操作空间，便于后期利用。综上，实验结果验证了本发明提出的基于Wallis滤波+直方图匹配的有重区域多/高光谱遥感图像匀色算法的有效性。

表1利用平均梯度法总体评估信息量对比

图像类别	参考图像	待匀色图像	直方图匹配	Wallis滤波+直方图匹配
					平均梯度指数	1.2662	0.9480	1.0912	1.4154

本发明还可有其它多种实施例，在不背离本发明精神及其实质的情况下，本领域技术人员当可根据本发明作出各种相应的改变和变形，但这些相应的改变和变形都应属于本发明所附的权利要求的保护范围。

Claims (4)

1.基于Wallis滤波+直方图匹配的有重多/高光谱遥感图像匀色方法，其特征在于：

步骤一、获得多/高光谱遥感图像，对多/高光谱遥感图像进行分波段处理，获得M组同波段灰度图像，第m组同波段灰度图像为X1、X2、X3、…、Xi、…、Xn；m＝1,2,...,M；i＝1,2,...,n；

步骤二、通过统计第m组同波段灰度图像X1、X2、X3、…、Xi、…、Xn的信息量及地物丰富度，选取出参考图像作为参考；

步骤三、取两幅图像，一幅是步骤二得到的第m组同波段灰度图像X1、X2、X3、…、Xi、…、Xn的参考图像Xi，另一幅是与参考图像有交叠的同波段灰度图像Xj；

对交叠区域进行分块处理，获得子块，设置子块取整数的平方；

将与参考图像有交叠的同波段灰度图像和参考图像的交叠区域做局部匀色；

用局部匀色图像替代原待匀色图像与参考图像的交叠区域，替换交叠区域后的原待匀色图像作为新的待匀色图像；

步骤四、对步骤三的参考图像和新的待匀色图像做直方图匹配全局匀色；

步骤五、将匀色处理后的图像作为参考图像，重复执行步骤三至步骤四，直至第m组同波段灰度图像全部匀色处理完，将匀色处理完的第m组同波段灰度图像拼接；

步骤六、重复步骤二至步骤五获得所有M组拼接好的灰度图像，将所有M组拼接好的灰度图像合成新多/高光谱遥感图像；

所述步骤二中通过统计第m组同波段灰度图像X1、X2、X3、…、Xi、…、Xn的信息量及地物丰富度，选取出参考图像作为参考；具体过程为：

步骤二一、信息熵表示第m组同波段灰度图像X1、X2、X3、…、Xi、…、Xn包含信息量的多少；

对于第m组同波段灰度图像X1、X2、X3、…、Xi、…、Xn，认为第m组同波段灰度图像X1、X2、X3、…、Xi、…、Xn每幅图像中各灰度值是相互独立的样本，单幅图像中各灰度值所占比例为p＝{p₁、p₂、p₃…p_n}；

利用以下公式分别计算第m组同波段灰度图像X1、X2、X3、…、Xi、…、Xn中每幅图像的信息熵：

式中，i表示第m组同波段灰度图像Xi的某灰度级；p_i表示第m组同波段灰度图像Xi的对应灰度级出现的概率；n表示第m组同波段灰度图像Xi的灰度级个数；

步骤二二、按各占百分之五十的权重统计第m组同波段灰度图像X1、X2、X3、…、Xi、…、Xn中每幅图像的参考价值，计算公式如下

式中，Q为参考价值，H为第m组同波段灰度图像Xi的信息熵，H_max为第m组同波段灰度图像X1、X2、X3、…、Xi、…、Xn中信息熵中最大值，N为第m组同波段灰度图像Xi的包含地物种类，N_max为第m组同波段灰度图像X1、X2、X3、…、Xi、…、Xn中包含地物种类最大值；

参考价值Q取最大值时所指的灰度图像Xi为第m组同波段灰度图像X1、X2、X3、…、Xi、…、Xn的参考图像。

2.根据权利要求1所述基于Wallis滤波+直方图匹配的有重多/高光谱遥感图像匀色方法，其特征在于：所述步骤三中取两幅图像，一幅是步骤二得到的第m组同波段灰度图像X1、X2、X3、…、Xi、…、Xn的参考图像Xi，另一幅是与参考图像有交叠的同波段灰度图像Xj；

对交叠区域进行分块处理，获得子块，设置子块取整数的平方；

将与参考图像有交叠的同波段灰度图像和参考图像的交叠区域做局部匀色；

用局部匀色图像替代原待匀色图像与参考图像的交叠区域，替换交叠区域后的原待匀色图像作为新的待匀色图像；具体过程为：

对参考图像和与参考图像有交叠的同波段灰度图像的交叠区域分别进行分块处理，获得子块，设置子块为整数的平方；

交叠区域内的参考图像每个子块对同波段灰度图像的每个子块进行Wallis滤波，完成交叠区域内的参考图像每个子块对同波段灰度图像每个子块的局部匀色；过程为：

式中，c为图像标准差扩展系数，其取值范围为0≤c≤1；b为图像的亮度扩展系数，其取值范围为0≤b≤1；当b取值靠近1时，待局部匀色图像的均值将向m_f靠近，当b取值靠近0时，待局部匀色图像的均值将向m_g靠近；m_f为参考图像第a个子块的灰度均值，m_g为待局部匀色图像第a个子块的灰度均值，s_f为参考图像第a个子块的灰度标准差，s_g为待局部匀色图像第a个子块的灰度标准差，g(x,y)为待局部匀色图像的第a个子块，g′(x,y)为局部匀色图像的第a个子块；

根据局部匀色图像的每个子块g′(x,y)得到交叠区域内与参考图像有交叠的同波段灰度图像的局部匀色图像；

用局部匀色图像替代原待匀色图像与参考图像的交叠区域，替换交叠区域后的原待匀色图像作为新的待匀色图像。

3.根据权利要求2所述基于Wallis滤波+直方图匹配的有重多/高光谱遥感图像匀色方法，其特征在于：所述c取值为0.7～0.8，b取值为0.6～0.7。

4.根据权利要求3所述基于Wallis滤波+直方图匹配的有重多/高光谱遥感图像匀色方法，其特征在于：所述步骤四中对步骤三的参考图像和新的待匀色图像做直方图匹配全局匀色；具体过程为：

步骤四一、分别对新的待匀色图像和参考图像的交叠部分进行直方图均衡化；过程为：

设新的待匀色图像的(x,y)处的灰度值为h(x,y)，对新的待匀色图像进行归一化处理，得到归一化处理后新的待匀色图像r：

式中，h_min为h(x,y)中的最小值，h_max为h(x,y)中的最大值；

设参考图像的(x′,y′)处的灰度值为h(x′,y′)，对参考图像进行归一化处理，得到归一化处理后参考图像z：

式中，h′_min为h(x′,y′)中的最小值，h′_max为h(x′,y′)中的最大值；

设归一化处理后新的待匀色图像r的直方图分布为P_r(r)，均衡化结果为P_s(s)，其中：

s＝T[r]

式中，T[r]为P_r(r)的累积分布函数；

同理设归一化处理后参考图像z的直方图分布为P_z(z)，均衡化结果为P_v(v)，其中：

v＝G[z]

式中，G[z]为P_z(z)的累积分布函数；

步骤四二、根据均衡化结果的一致性，认为

s＝v

由公式s＝T[r]、v＝G[z]和s＝v，得到一种建立在r与z之间的映射关系：

z＝G^-1[T[r]]

根据这种映射算法将新的待匀色图像进行匀色处理。

Images