2. 北京遥感信息研究所，北京 100192
2. Beijing Institute of Remote Sensing, Beijing 100192, China
1 引 言
对原始灰度图像进行彩色化处理，能够大大提高遥感图像的可判读性，更容易提取所需的信息。
对于遥感影像彩色化，传统方法有假彩色增强、真彩色增强、伪彩色增强和IHS变换方法等，这些方法都能够有效提高遥感图像的可分辨性。但是假彩色和伪彩色技术所得到的增强图像反映的并不是场景的真实色彩，所以对于非专业图像判读人员来说往往难以理解和接受[1]。真彩色技术要求彩色合成时与遥感成像时的太阳光谱成份相一致，这种条件难以满足。IHS方法能提高目标探测能力，对同一景物进行多光谱观察，可以揭示出目标在单一传感器观察时无法探测或不易探测的特征。但这种方法如果应用于不同数据源、不同传感器、不同时段的数据,在图像配准时经常会出现矛盾。
近年来，遥感彩色化增强领域的研究，主要是针对某种特定类型图像在上述传统方法的基础上加以改进。例如，文献[2]提出的地物波谱数据辅助的SPOT影像模拟真彩色方法，该方法从遥感波谱物理特性出发，进行地物波谱库辅助下的遥感影像真彩色模拟方法研究，旨在利用地物波谱库中丰富的光谱信息实现SPOT影像蓝光波段模拟。文献[3]提出的基于多方向多源合成孔径雷达数据融合的假彩色正射影像制作，利用不同侧视方向、不同分辨率、多极化的SAR影像数据，融合制作信息量丰富，能够满足应用要求的SAR假彩色正射影像。
文献[4]基于人类的视觉感知研究，提出了一种新的颜色空间，称为lαβ颜色空间。其中各通道之间相关性很小，因此在不同的颜色通道之间进行不同的操作时，不会受其他通道的影响或者影响其他通道。文献[5]提出一种利用lαβ颜色空间在两幅彩色图像间进行颜色传递的算法。文献[6]在此基础上提出一种从彩色图像到灰度图像间进行颜色传递的算法：对于目标图像上的每一个像素，都从源图像上找一个与它匹配的像素，它使用一个像素的邻域的统计值来指导匹配过程。找到匹配像素后，就将该像素的颜色信息传输到目标图像但保留目标图的亮度信息。
国内外学者提出了一些针对该方法的改进与应用。文献[7]提出将颜色传递技术应用于夜视图像，取得了较好的效果。文献[8]提出基于颜色传递的单波段热成像自然感彩色化处理算法，该算法应用于几种场景类型的单波段红外图像，得到了色彩自然感较好的彩色化热图像。文献[9]为克服当前偏振图像融合方法存在的不足，提出了一种基于颜色迁移和聚类分割的偏振图像近自然彩色融合方法。
由于遥感图像的复杂性和颜色传递的局限性，目前该技术一般局限于夜视图像、单波段热图像等相对较简单的图像类型。本文根据遥感图像的特点，将颜色传递改进后尝试用于全色遥感图像彩色化增强，所采用的参考彩色图像不需要与目标图像配准，只要图像特征、风格相近似，就可以对纹理结构相对较简单的全色遥感图像实现全自动的彩色化增强。本文所作的主要改进如下。
1.1 复数小波包降噪预处理
遥感图像往往会受到噪声的污染。全色图像如果不预先进行降噪的预处理，直接转化彩色后噪声会更明显，影响到最终的增强效果。灰度图像的噪声也会影响和干扰后续的图像分割等操作。笔者采用分频带高斯模型，利用复数小波包，采用软阈值方法进行小波降噪预处理。
1.2 多尺度图像分割
利用多尺度图像分割技术[10]，对遥感图像进行区域分割。将图像分割为若干子区域后，在子区域之间利用直方图进行匹配，在匹配好的子区域之间搜索像素与颜色传递。由于原算法是对图像中所有像素进行全局搜索与匹配，而分割后只需要对子区域进行搜索与匹配，因而大大降低了匹配像素中的搜索代价。同时还提高了颜色传递的准确度，彩色化增强效果更好。
本文算法的主要框架为：首先采用复数小波包工具，对目标灰度图像使用软阈值技术进行滤波降噪；利用多尺度归一化切分的分割方法对图像进行分割；将目标图像和源图像转换到lαβ颜色空间；通过亮度直方图，匹配子区域；在对应的子区域内，对各个子区域基于邻域的统计值进行亮度匹配；取对应的颜色值赋予对应的区域内的像素；最后将结果图像转换回RGB空间。算法流程如下图所示。
 |
| 图 1 流程图 Fig. 1 Flow chart |
大多数情况下，遥感图像噪声现为孤立的离散点，没有空间关联性，均可以视为高斯分布，因此本文考虑的噪声模型以高斯为主[11,12]。
2.1 降噪方法的选择
小波变换是目前使用最多的变换域降噪方法。该方法不仅能够在不同分辨率下根据信号和噪声分布的特点降噪，还可以很好地刻画信号的非平稳特征。小波变换的关键是阈值和小波基的选择。
阈值可分为硬阈值和软阈值两种类型。硬阈值能够更好地保存边缘信息，更接近实际情况，但噪声抑制效果较差；软阈值具有连续性，降噪效果好，获得的结果更加平滑，因此视觉上更加自然，容易接受，但图像边缘信息模糊。软阈值往往会应用到细节少的分解级别(小波子带)中，而硬阈值用在其他情形[13,14,15,16]。由于本文采用的是分频带高斯模型，所以对于小波子带来说，更适合采用软阈值方法。
小波基可选的种类很多，如单小波、多小波、小波包等。小波基的选择对降噪的效果影响很大，对于不同类型的图像，最优的小波基是不同的[17]。由于复数小波包变换兼备复数小波变换和小波包变换的优点，具有更好的方向选择特性和更为细致的高频信息，并且具有较强的去除遥感图像中彩色噪声的能力，所以本文选择复数小波包作为小波基[18]。
2.2 降噪的算法及过程
对于输入的遥感图像，本文使用统一的分频带高斯模型。在这个模型中，每个频带取决于所选择的小波基。对于复数小波包分解，即是每个复数小波包频带。
不妨假设X是输入图像，Z是由X添加上彩色噪声P获得的带噪的观测图像，如公式(1)所示
在每个频带中假设噪声是一个高斯分布，均值为0，方差为σk2。其中方差σk2可以使用仿真的方法由P的统计参数估计出来。
在估计获得各频带的模型参数之后，首先通过一个单尺度的快速估计方法(如Wiener或者RHEA)获得估计值ξ′，再通过软阈值滤波修正获得预估值e2
然后利用预估值e2，分频带计算所有ξ的值
式中，z是观测值Z在对应分解域下的值。
最后再利用复数小波包逆变换将ξ还原到空间域上，得到结果X，实现遥感图像的降噪[19,20,21]。
3 多尺度归一化切分
遥感图像分割，属于自然环境图像分割的一种，相对于室内场景或工业环境等结构性环境而言，要复杂很多。就像植被、水系和交通等要素，仅依靠灰度信息不足以对图像进行分割，还需要各种物体的物理及概念层次的精细描述[22,23]。
对于本文来说，图像分割的目的不是用于特征提取及分析，要求相对较低，主要是能够区分不同灰度、颜色和纹理；但由于遥感图像尺寸巨大而纹理复杂，所以寻找一种能够利用颜色和纹理信息得到较好结果、并且速度较快的分割方法，是非常重要的。
1997年，Shi和Malik提出归一化切分的分割方法，将分割问题形式化为图分割问题，利用颜色、纹理或者结合使用颜色和纹理及其他特性进行分割，将问题化解为求特征向量和特征值的问题。
文献[10]对归一化切分进行改进，对图像的多重比例尺采用并行计算，没有迭代，计算高效，能够分割大尺寸的图像。用这个方法能够准确分割以前因为尺寸太大而无法有效分割的图像，算法框架如下。
给定一个p×q的图像I，对于s=1，2，…，S(S为缩放比例)：
按照
采样，像素点i∈C，图像Is以规则的格网映射自Is-1，这里ρ是采样因子。
计算约束
计算Is的仿射Wsc，小半径为r，利用在比例s的图像边。
通过(Wsc,Cs－1,s)s，计算W、C，如下式所示
计算投影Q
这里D是一个对角线矩阵
计算V，是矩阵
的前k个特征向量；计算
，使之离散化。
图 2是采用该方法对卫星图像分割的示例。
 |
| 图 2 图像分割示例 Fig. 2 Example of image segmentation |
首先将彩色源图像和目标灰度图都转换到lαβ颜色空间；然后利用亮度直方图来匹配最接近的子区域；然后在建立好对应关系的子区域之间搜索最匹配的像素，将颜色值α、β值一一赋给目标灰度图，再转换为RGB颜色空间，就得到了目标图像对应的彩色图。对于较复杂的图像，可以通过人工交互的方式来指定对应的传输区域。具体步骤如下：
(1) 根据文献[4]提出的转换方法，将彩色源图像和目标灰度图像从RGB颜色空间转换到lαβ颜色空间。
(2) 对分割后的子区域利用亮度直方图进行匹配。先计算目标图像与源图像子区域的亮度直方图；然后利用直方图交叉法，使目标图像中的每一个子区域在源图像中找到最接近的子区域，建立起对应关系。
(3) 在对应好的子区域中搜索最匹配的像素，然后赋予颜色值。像素的匹配原则是根据像素的邻域亮度均值和标准差，邻域的尺寸可选为5×5或7×7；当查找到最匹配的像素时，将匹配像素的颜色值α、β赋给目标图像中的像素，而像素亮度值l保持不变。
(4) 将结果图像从lαβ颜色空间转换到RGB颜色空间[4]。
5 试验与分析
本文选取不同风格的彩色源图像和灰度目标图像，按照图 3所示的算法框架与流程进行试验。图 4采用的是SPOT5的全色灰度图像，空间分辨率为2.5m，尺寸为256像素×256像素。图 5像采用的是CNES(法国国家空间研究中心)提供的Nimes地区全色灰度图[20]，裁切的尺寸为256像素×256像素，添加了0.5的噪声方差。用于参考的真彩色图像其空间分辨率为0.6m。试验的中间过程与最终结果的对比如图 4—图 5所示。
 |
| 图 3 算法框架与流程 Fig. 3 Algorithm overview |
 |
| 图 4 试验的中间过程与最终结果对比1 Fig. 4 Procedure and comparison of rsults 1 |
 |
| 图 5 试验的中间过程与最终结果对比2 Fig. 5 Procedure and comparison of results 2 |
由图 4—图 5中的试验结果可以看出：基于图像分割和色彩传递相结合的方法只需要一个特征、风格相近的参考彩色图像，无须配准，就可以对纹理结构相对较简单的遥感图像实现全自动的彩色化增强。通过利用多尺度图像分割的方法，在分割后的子区域之间进行亮度匹配和颜色传递，能够减少像素间的颜色传递误差；与软阈值的复数小波包降噪预处理结合在一起，能明显增强遥感图像的细节分辨能力。彩色化的结果图像无论是在美观程度，还是目标的可识别度都有了明显增强。
本文还对传统的数据融合算法与本文方法进行了对比试验。图 6(a)为全色灰度影像(空间分辨率为2m，尺寸为256像素×256像素)，图 6(b)为采用数据融合方法得到的彩色图像(与分辨率为10m的多光谱图像采用pansharp算法得到的融合结果)，图 6(c)为本文算法的结果。
 |
| 图 6 传统数据融合与本文结果对比 Fig. 6 Comparison of results |
通过与传统数据融合算法的对比，可以发现本文提出的算法有两个优点：①其参考彩色源图像不需要与目标图像配准；②较好地将复数小波包图像预处理与颜色传递相结合，对图像的细节部分有着明显的复原与增强。
由于原算法是对图像中所有像素(假设像素数为N)进行全局搜索与匹配，则搜索代价为O(Nlog(N))；而分割后只需要对子区域(假设子区域内的像素数为n)进行搜索与匹配，搜索代价为O(nlog(n))，因为n<<N，所以大大降低了匹配像素中的搜索代价。本算法额外增加的时间代价包含小波降噪和图像分割，这两种数据预处理本文均采用了高效的算法，基于复数小波包的降噪与去模糊算法的复杂度为O(N)[19]，多尺度的图像分割算法的复杂度也保持为线性[10]。因此整体的时间代价主要取决于图像的尺寸大小。
试验中，发现这种方法存在以下局限性：①目前该方法只能对纹理相对简单、亮度反差较大的全色遥感图像实现全自动的彩色化，对于亮度值接近但属性不同或者纹理复杂的图像还需要人工指定颜色分配；②参考图像尽可能挑选特征、风格相近的图像，最终彩色化的效果很大程度上依赖于参考图像的选择。
6 结 论
本文根据遥感图像的特点，尝试利用复数小波包降噪与去模糊、多尺度图像分割和颜色传递技术相结合的方法，来实现无须配准的全自动的遥感图像彩色化增强。试验结果改善了全色遥感图像的细节现力和视觉感受。本文的方法是对遥感光学图像传统数据融合方法之外的一种探索和补充。
针对所提方法的局限性，未来可对遥感图像彩色化作如下改进：①尝试基于样例的彩色化技术[24]，以期实现更准确的颜色传递；②对图像的分割方法和像素的统计值匹配及映射方式作进一步的研究与改进；③进一步研究其他不同类型颜色空间的转换方法。
| [1] | ZHANG Yan, WANG Yu. Remote Sensing Image Colorization Using Optimization[J]. Science Paper Online, 2008,3(2): 114-119. (张砚,王昱. 基于优化的遥感影像彩色化[J]. 中国科技论文在线,2008,3(2): 114-119.) |
| [2] | ZHU Changming, LUO Jiancheng, SHEN Zhanfeng, et al. Approach of SPOT True Color Image Simulation Assisted by Spectral Library Data[J]. Acta Geodaetica et Cartographica Sinica, 2010, 39(2):169-174. (朱长明,骆剑承,沈占锋,等. 地物波谱数据辅助的SPOT影像模拟真彩色方法研究[J].测绘学报,2010, 39(2):169-174.) |
| [3] | ZHANG Jixian, WEI Jujie, ZHAO Zheng, et al. Color Orthophoto Map Generation Based on Multi-direction and Multi-polarization SAR Data Fusion[J]. Acta Geodaetica et Cartographica Sinica, 2011, 40(3):276-282. (张继贤,魏钜杰,赵争,等. 基于多方向多源合成孔径雷达数据融合的假彩色正射影像制作[J]. 测绘学报,2011, 40(3):276-282.) |
| [4] | RUDERMAN D L, CRONIN T W, CHIAO C C. Statistics of Cone Responses to Natural Images: Implications for Visual Coding[J]. Optical Society of America, 1998, 15(8): 2036-2045. |
| [5] | REINHARD E, ASHIKMIN M, GOOCH B, et al. Color Transfer between Images[J]. IEEE Computer Graphics and Applications, 2001, 21(5): 34-40. |
| [6] | TOMIHISA W, MICHAEL A, KLAUS M. Transferring Color to Greyscale Images[C]//Proceedings of the SIGGRAPH.San Antonio:[s.n.], 2002: 277-280. |
| [7] | TOET A. Natural Color Mapping for Multiband Night Vision Imagery[J]. Information Fusion, 2003, 4(3): 155-166. |
| [8] | ZHAO Yuanmeng, WANG Lingxue, JIN Weiqi, et al. Single-band Infrared Image Colorization Algorithm based on Color Transfer[J]. Acta Optica Sinica, 2009, 29(3): 654-658.(赵源萌, 王岭雪, 金伟其,等. 基于色彩传递的单波段热图像彩色化处理算法[J]. 光学学报, 2009, 29(3): 654-658.) |
| [9] | ZHOU Pucheng, ZHANG Hongkun, XUE Mogen. Polarization Image Fusion Method Using Color Transfer and Clustering-Based Segmentation[J]. Acta Photonica Sinica, 2011, 40(1): 149-153. (周浦城, 张洪坤, 薛模根. 基于颜色迁移和聚类分割的偏振图像融合方法[J]. 光子学报, 2011, 40(1): 149-153.) |
| [10] | COUR T, BENEZIT F, SHI J. Spectral Segmentation with Multiscale Graph Decomposition[J]. Computer Vision and Pattern Recogniton, 2005, 2:1124-1131. |
| [11] | SUN Yankui. Wavelets Analysis and Application[M]. Beijing: China Machine Press,2005. (孙延奎. 小波分析及其应用[M]. 北京:机械工业出版社,2005.) |
| [12] | NI Lin. Wavelet Translation and Image Processing[M]. Hefei: Chinese Science and Technology University Press, 2010. (倪林. 小波变换与图像处理[M]. 合肥:中国科学技术大学出版社,2010.) |
| [13] | CHANG S G, VETTERLI M. Spatial Adaptive Wavelet Thresholding for Image Denoising[C]// Proceedings of ICIP.[S.l.]:ICIP, 1997. |
| [14] | DONOHO D.Denoising by Soft Thresholding[J]. IEEE Transaction, 1995, 41:613-627. |
| [15] | YANG Yuanxi. Properties of the Adaptive Filtering for Kinematic Positioning[J]. Acta Geodaetica et Cartographica Sinica, 2003, 32(3): 189-192.(杨元喜. 动态定位自适应滤波解的性质[J]. 测绘学报,2003,32(3):189-192.) |
| [16] | YANG Yuanxi, CUI Xianqiang. Influence Functions of Colored Noises on Kinematic Positioning:Taking the AR Model of First Class as an Example[J]. Acta Geodaetica et Cartographica Sinica, 2003, 32(1):6-10. (杨元喜,崔先强. 动态定位有色噪声影响函数:以一阶AR模型为例[J].测绘学报,2003,32(1):6-10.) |
| [17] | JIANG Sanping. Image Denoising Based on Wavelets Translation[M]. Beijing: National Defense Industry Press, 2009. (姜三平. 基于小波变换的图像降噪[M]. 北京: 国防工业出版社, 2009.) |
| [18] | JALOBEANU A, BLANC-FERAUD L, ZERUBIA J. Satellite Image Deblurring Using Complex Wavelet Packets[J]. International Journal of Computer Vision, 2003,51(3): 205-217. |
| [19] | JALOBEANU A, ZERUBIA J, BLANC-FERAUD L. Bayesian Estimation of Blur and Noise in Remote Sensing Imaging[C]//Blind Image Deconvolution: Theory and Applications.London:CRC Press, 2007:239- 271. |
| [20] | ZHANG Yan, MAN Yiyun. Satellite Image Adaptive Restoration Using Periodic Plus Smooth Image Decomposition and Complex Wavelet Packet Transforms[J]. Tsinghua Science and Technology, 2012, 17(3): 337-343. |
| [21] | ZHANG Yan, LI Xianying, MAN Yiyun. Remote Sensing Images Upsampling Based on Projection onto Convex Sets and Complex Wavelet Packet Transfer[J].Chinese Journal of Computers, 2011, 34(3): 482-488. (张砚, 李先颖, 满益云. 基于凸集投影法和复数小波包域的遥感图像上采样研究[J]. 计算机学报, 2011, 34(3): 482-488.) |
| [22] | LINDA G S, GEORGE C S. Computer Vision[M]. ZHAO Qingjie, QIAN Fang, CAI Lidong,Translated. Beijing: China Machine Press, 2005. (LINDA G S, GEORGE C S. 计算机视觉[M]. 赵清杰,钱芳,蔡利栋,译. 北京:机械工业出版社,2005.) |
| [23] | XIAO Pengfeng, FENG Xuezhi, ZHAO Shuhe, et al. Segmentation of High-resolution Remotely Sensed Imagery Based on Phase Congruency[J]. Acta Geodaetica et Cartographica Sinica, 2007,36(2):146-151,186. (肖鹏峰,冯学智,赵书河,等. 基于相位一致的高分辨率遥感图像分割方法[J]. 测绘学报,2007,36(2):146-151,186.) |
| [24] | IRONY R, COHEN D, LISCHINSKI D. Colorization by Example[C]//Proceedings of Eurographics Symposium on Rendering. Konstanz:[s.n.],2005. |