本发明涉及一种算法，具体是一种基于Curvelet变换的彩色图像融合算法。
背景技术：
图像是人们从客观世界获取信息的重要来源，图像的清晰度的好坏影响着后期的图像处理与分析。在采用彩色CCD或者CMOS数码相机采集某一场景中的多个目标时,由于CCD或者CMOS图像传感器的景深问题和聚焦点的不同,无法保证场景中的所有的目标都是清晰的。在不改变图像采集系统的情况下，为了使整个场景中的目标是清晰的，可以采用多传感器图像融合的方法来解决。
多传感器信息融合亦称多源信息融合，是指对来自多个传感器的数据进行多级别、多方面、多层次的处理和综合，从而获得更可靠、更丰富、更精确、更有意义的信息。多种文献对信息融合的定义大致可概括为：对按时序获得的多传感器观测数据在一定准则下加以自动分析、优化综合，以完成所需的准确状态估计和决策任务而进行的处理过程。
可见光图像与红外图像融合是多传感器图像融合研究中的一个重要部分，它的一个主要应用领域就是基于图像的战场侦察。由于可见光成像传感器只敏感于目标场景的反射，而与目标场景的热对比度无关，因此符合人的视觉，具有直观、清晰的特点，能通过图像一目了然地观察到敌方前沿地形、阵地布设、武器装备、兵力部署和调动。而红外成像传感器是靠探测目标与背景间的热辐射差异来识别目标的，因而具有特殊的识别伪装能力，有利于目标隐蔽情况下或夜间和恶劣气候条件下探测目标信息，能避免外界因素的干扰。监控侦察的任务不仅是要发现、识别目标，而且应确定目标所在的准确位置，要求获取的图像具有较高的清晰度。红外与可见光图像的融合将有利于综合红外图像较好的目标指示特性和可见光图像的清晰场景信息，因此研究红外与可见光图像的融合很有意义。另外，可见光成像传感器对目标场景的亮度变化敏感，能较好的显示对比度、形状和纹理细节等，但它却受光照条件的约束。红外成像传感器反映目标和场景的红外辐射特性，可以昼夜两用，但缺乏对目标场景纹理细节的表现。基于可见光和红外图像的互补性，可以提取并通过一定算法组合两幅同一场景可见光和红外图像的各自特征，得到一幅融合图像，可以用在比如军事领域中的隐藏目标检测当中。
以往的可见光图像和红外图像都是单色图像，其最大缺陷是信噪比低且图像缺乏体视感（深度感），不利于目标的探测。利用人眼彩色视觉的高分辨力和高灵敏度的特性，实现图像彩色融合可提高探测和识别目标能力，具有明显的战场效果，能增强场景理解、突出目标，有利于在隐藏、伪装和迷惑的军用背景下更快更精确地探测目标。另外，将融合图像显示成简单、适合人眼观察的自然形式，可明显改善人眼的识别性能，减小操作者的疲劳感。但是，由于彩色图像融合的色彩往往与人类本身的视觉不同，过于鲜明的彩色图像对比易造成观察者的疲劳和烦恼，不适宜长时间的密集观察,甚至难以获得理想的效果。因此，能够更好的描述图像的细节和边缘等几何特征，寻找接近人眼视觉感的彩色图像融合算法成为研究重点。
近几十年来，国际上在图像信息融合的不同层次上开展了大量的模型和算法研究。但是由于应用环境（目标与场景）不同，图像的特性多种多样，对各种任务的要求也各不相同，因此图像融合算法也各种各样。目前尚未发现在各种任务要求下融合效果全优的算法。所使用的彩色图像融合算法大体可分为三类：灰度融合后采用伪彩色编码融合，基于RGB色空间的彩色融合，基于HIS色空间的彩色融合。
从当前国内外图像融合技术的应用来看，国外发达国家由于起步较早，投入较大，在融合系统的融合算法和实际应用的研究上处于领先地位，迄今为止已研制出上百种极具应用价值的图像融合系统。最早投入使用的图像融合系统是美国陆地资源卫星(LANDSTAR)。英国以通用II类组件为基础，研制出了具有图像融合处理功能的双波段热像仪。美国TI公司于1995年底从美国夜视与电子传感器管理局(NVESD)获得了将红外与三代微光图像融合设计集成到先进直升机驾驶(AHP)传感器系统的合同，信号处理由TMS320C30 DSP完成。美欧多套大型战区级传感器信息融合演示验证系统，其中相当重要的部分就是武器平台上分布式的图像融合装置。英美联合研制的“追踪者”战术侦察车是将热成像仪、电视摄像仪以及激光测距仪等多个传感器进行融合利用。美国波音公司航空电子飞行实验室2000年成功地演示验证了联合攻击机航空电子系统的多传感器信息融合技术和功能。法国研 制的Vp，MPIRMB红外全景监视系统，美国研制的凝视型N/RRAF-43，扫描型N/RRAF-44红外警戒接收机及AN/ALQF-153，154，156等多普勒雷达探测器，美国和加拿大联合研制的AN/SARF-8红外系统及荷兰研制的单、双波段IR SCAN等系统都成功地应用了图像融合技术。
众所周知，在小波变换出现的近20年间，它在信号处理中得到了广泛的应用，基于小波变换的图像融合的方法被广泛采用。然而，对于二维图像处理，常用的二维小波变换的基函数是各向同性的，变换系数的局部模极大值只能反映出这个小波系数出现的位置是“过”边缘的，而无法表达“沿”边缘的信息，这也正是在基于小波的一系列处理方法，如图像去噪、图像融合以及图像压缩等应用中，均不可避免地在图像边缘和细节位置引入一定程度模糊的原因所在。针对小波变换的上述缺点，Donoho等人于1999年提出了第一代Curvelet变换理论，该理论是Ridgelet理论衍生而来，是由一种特殊的滤波过程和多尺度Ridgelet变换组合而成的。与小波变换不同，除了尺度和位移参量，Curvelet还增加了一个方向参量，对图像的边缘，如曲线、直线等几何特征的表达更加优于小波，因而具有更好的应用前景。
国内对图像融合技术研究起步较晚。目前，上海交通大学，中科院上海技术物理研究所，中科院遥感所，武汉测绘大学等单位正在此领域进行一些相关研究。南京理工大学和北京理工大学分别进行了一些探索性的研究，在基于伪彩色的图像融合领域都取得了有益的成果。
技术实现要素：
本发明的目的在于提供一种基于Curvelet变换的彩色图像融合算法，以解决上述背景技术中提出的问题。
为实现上述目的，本发明提供如下技术方案：
一种基于Curvelet变换的彩色图像融合算法，包括以下步骤：(1) 选择Curvelet 的紧框架作为原子模型，设定该框架中的离散尺度参数j、离散方向参数θ和位移参数 k 的取值范围，即-10 ≤ j ≤ 6，θ ∈ [0，π)，k ＝ (k1，k2)，其中 k1 为原子中心沿 x 轴的位移，k2 为原子中心沿 y 轴的位移，k1，k2 ∈ [0，n)，n 表示图像块的边长；
(2)在上述设定的各参数范围内，离散尺度参数 j，其间隔为 1，离散方向参数 θ，其间隔为离散位移参数k1，k2，其间隔为 1，从而形成了 Curvelet 冗余字典；
将输入图像分块，通过正交匹配追踪 OMP 算法求每个图像块在 Curvelet 冗余字典下的稀疏分解系数向量，组合所有稀疏分解系数向量得到稀疏矩阵，将该稀疏矩阵与 Curvelet 冗余字典相乘，得到输入图像的稀疏表示结果。
作为本发明进一步的方案：其中步骤(1) 所述的Curvelet 紧框架由三元组 (j，θ，k) 表征的 Curvelet 函数 来表示，其中 j 为离散尺度参数，θ 为离散方向参数和 k 为位移参数，即：
其中，Da 为尺度算子， a ＝ 2-2j ；Rθ 为 θ 角度的旋转算子，
离散位移参数 k1，k2 均为整数；(x，y) 为图像中像素点的坐标值；
函数是小波母函数：其中t为自由变量；将上述Curvelet 的紧框架通过点积的形式简化，即：
其中，X ＝ 22j xcosθ-22jysinθ-k1，Y ＝ 2j xsinθ+2jycosθ-k2
作为本发明再进一步的方案：其中步骤(3) 所述的通过正交匹配追踪 OMP 算法求每个图像块在Curvelet冗余字典下的稀疏分解系数向量，按如下步骤进行： (3a) 令初始残差信号为原始信号，初始迭代次数t ＝0，根据所需稀疏分解的精确度确定信号稀疏分解中需要选取最大原子个数 m 和信号能量阈值 ξstop ；
(3b) 计算原始信号或残差信号与 Curvelet 冗余字典中所有原子的内积，并从冗余字典中选择使内积最大的原子；
(3c) 利用 Gram-Schmidt 正交化方法对内积最大的原子进行正交化处理；
(3d) 用残差信号减去最大内积值与内积最大原子的乘积来更新残差信号；
(3e) 确定是否停止迭代：如果迭代次数t＜m 或者步骤 (3d) 中更新后的残差信号能量小于ξstop，则停止迭代，得到相应的稀疏分解系数向量，否则t增加 1，进行下一次迭代，转至步骤 (3b)。
与现有技术相比，本发明的有益效果是：本发明计算复杂度低，运行时间短；能更好地捕捉到图像中的曲线奇异性，从而提高了稀疏表示图的质量和视觉效果。
具体实施方式
下面对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。
本发明实施例中，一种基于Curvelet变换的彩色图像融合算法，包括以下步骤：(1) 选择Curvelet 的紧框架作为原子模型，设定该框架中的离散尺度参数j、离散方向参数θ和位移参数 k 的取值范围，即-10 ≤ j ≤ 6，θ ∈ [0，π)，k ＝ (k1，k2)，其中 k1 为原子中心沿 x 轴的位移，k2 为原子中心沿 y 轴的位移，k1，k2 ∈ [0，n)，n 表示图像块的边长；(3)在上述设定的各参数范围内，离散尺度参数 j，其间隔为 1，离散方向参数 θ，其间隔为 离散位移参数 k1，k2，其间隔为 1，从而形成了 Curvelet 冗余字典；将输入图像分块，通过正交匹配追踪 OMP 算法求每个图像块在 Curvelet 冗余字典下的稀疏分解系数向量，组合所有稀疏分解系数向量得到稀疏矩阵，将该稀疏矩阵与 Curvelet 冗余字典相乘，得到输入图像的稀疏表示结果。
其中步骤(1) 所述的Curvelet 紧框架由三元组 (j，θ，k) 表征的 Curvelet 函数 来表示，其中 j 为离散尺度参数，θ 为离散方向参数和 k 为位移参数，即：
其中，Da 为尺度算子， a ＝ 2-2j ；Rθ 为 θ 角度的旋转算子，
离散位移参数 k1，k2 均为整数；(x，y) 为图像中像素点的坐标值；
函数是小波母函数：其中 t 为自 由变量；将上述 Curvelet 的紧框架通过点积的形式简化，即：
其中，X ＝ 22j xcosθ-22jysinθ-k1，Y ＝ 2j xsinθ+2jycosθ-k2。
其中步骤(3) 所述的通过正交匹配追踪 OMP 算法求每个图像块在Curvelet冗余字典下的稀疏分解系数向量，按如下步骤进行： (3a) 令初始残差信号为原始信号，初始迭代次数t ＝0，根据所需稀疏分解的精确度确定信号稀疏分解中需要选取最大原子个数 m 和信号能量阈值 ξstop ；
(3b) 计算原始信号或残差信号与 Curvelet 冗余字典中所有原子的内积，并从冗余字典中选择使内积最大的原子；
(3c) 利用 Gram-Schmidt 正交化方法对内积最大的原子进行正交化处理；
(3d) 用残差信号减去最大内积值与内积最大原子的乘积来更新残差信号；
(3e) 确定是否停止迭代：如果迭代次数t＜m 或者步骤 (3d) 中更新后的残差信号能量小于ξstop，则停止迭代，得到相应的稀疏分解系数向量，否则t增加 1，进行下一次迭代，转至步骤 (3b)。