CN110163830A - 基于Riesz-Lap变换及PCNN的图像融合方法 - Google Patents
基于Riesz-Lap变换及PCNN的图像融合方法 Download PDFInfo
- Publication number
- CN110163830A CN110163830A CN201910318054.5A CN201910318054A CN110163830A CN 110163830 A CN110163830 A CN 110163830A CN 201910318054 A CN201910318054 A CN 201910318054A CN 110163830 A CN110163830 A CN 110163830A
- Authority
- CN
- China
- Prior art keywords
- image
- riesz
- transformation
- indicates
- indicate
- 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.)
- Granted
Links
- 230000009466 transformation Effects 0.000 title claims abstract description 113
- 238000000034 method Methods 0.000 title claims abstract description 31
- 238000013528 artificial neural network Methods 0.000 claims abstract description 34
- 230000004927 fusion Effects 0.000 claims description 24
- 238000010586 diagram Methods 0.000 claims description 19
- 238000005070 sampling Methods 0.000 claims description 18
- 238000000354 decomposition reaction Methods 0.000 claims description 13
- 238000002156 mixing Methods 0.000 claims description 13
- 238000004422 calculation algorithm Methods 0.000 claims description 6
- 239000011159 matrix material Substances 0.000 claims description 6
- 238000013507 mapping Methods 0.000 claims description 3
- 230000004044 response Effects 0.000 claims description 3
- 238000013519 translation Methods 0.000 claims description 3
- 238000005481 NMR spectroscopy Methods 0.000 abstract description 9
- 239000003814 drug Substances 0.000 abstract description 7
- 230000000007 visual effect Effects 0.000 abstract description 3
- 238000012360 testing method Methods 0.000 abstract description 2
- 238000010168 coupling process Methods 0.000 abstract 1
- 238000005859 coupling reaction Methods 0.000 abstract 1
- 230000016776 visual perception Effects 0.000 abstract 1
- 230000000694 effects Effects 0.000 description 9
- 238000011084 recovery Methods 0.000 description 8
- 230000000875 corresponding effect Effects 0.000 description 5
- 238000003384 imaging method Methods 0.000 description 4
- 241000282326 Felis catus Species 0.000 description 2
- 241000638935 Senecio crassissimus Species 0.000 description 2
- 230000008901 benefit Effects 0.000 description 2
- 238000004364 calculation method Methods 0.000 description 2
- 230000002596 correlated effect Effects 0.000 description 2
- 238000011156 evaluation Methods 0.000 description 2
- 239000000284 extract Substances 0.000 description 2
- PCHJSUWPFVWCPO-UHFFFAOYSA-N gold Chemical compound [Au] PCHJSUWPFVWCPO-UHFFFAOYSA-N 0.000 description 2
- 239000010931 gold Substances 0.000 description 2
- 229910052737 gold Inorganic materials 0.000 description 2
- 230000006872 improvement Effects 0.000 description 2
- 210000002569 neuron Anatomy 0.000 description 2
- 238000007500 overflow downdraw method Methods 0.000 description 2
- 238000013441 quality evaluation Methods 0.000 description 2
- 230000003595 spectral effect Effects 0.000 description 2
- 230000003044 adaptive effect Effects 0.000 description 1
- 230000002776 aggregation Effects 0.000 description 1
- 238000004220 aggregation Methods 0.000 description 1
- 238000004458 analytical method Methods 0.000 description 1
- 238000013459 approach Methods 0.000 description 1
- 230000015556 catabolic process Effects 0.000 description 1
- 238000006243 chemical reaction Methods 0.000 description 1
- 239000003795 chemical substances by application Substances 0.000 description 1
- 238000004883 computer application Methods 0.000 description 1
- 230000021615 conjugation Effects 0.000 description 1
- 230000007812 deficiency Effects 0.000 description 1
- 238000003745 diagnosis Methods 0.000 description 1
- 239000006185 dispersion Substances 0.000 description 1
- 235000013399 edible fruits Nutrition 0.000 description 1
- 238000002474 experimental method Methods 0.000 description 1
- 238000002844 melting Methods 0.000 description 1
- 230000008018 melting Effects 0.000 description 1
- 238000012986 modification Methods 0.000 description 1
- 230000004048 modification Effects 0.000 description 1
- 230000008569 process Effects 0.000 description 1
- 238000004445 quantitative analysis Methods 0.000 description 1
- 238000011160 research Methods 0.000 description 1
- 238000001228 spectrum Methods 0.000 description 1
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T5/00—Image enhancement or restoration
- G06T5/10—Image enhancement or restoration using non-spatial domain filtering
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T5/00—Image enhancement or restoration
- G06T5/50—Image enhancement or restoration using two or more images, e.g. averaging or subtraction
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2207/00—Indexing scheme for image analysis or image enhancement
- G06T2207/20—Special algorithmic details
- G06T2207/20016—Hierarchical, coarse-to-fine, multiscale or multiresolution image processing; Pyramid transform
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2207/00—Indexing scheme for image analysis or image enhancement
- G06T2207/20—Special algorithmic details
- G06T2207/20048—Transform domain processing
- G06T2207/20064—Wavelet transform [DWT]
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2207/00—Indexing scheme for image analysis or image enhancement
- G06T2207/20—Special algorithmic details
- G06T2207/20084—Artificial neural networks [ANN]
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2207/00—Indexing scheme for image analysis or image enhancement
- G06T2207/20—Special algorithmic details
- G06T2207/20212—Image combination
- G06T2207/20221—Image fusion; Image merging
Landscapes
- Physics & Mathematics (AREA)
- General Physics & Mathematics (AREA)
- Engineering & Computer Science (AREA)
- Theoretical Computer Science (AREA)
- Image Processing (AREA)
Abstract
本发明公开了基于Riesz‑Lap变换及PCNN的图像融合方法，Riesz变换具有分离出图像低级特征的能力，且人类视觉感知的重要信息来也来自图像的低级特征。本发明提出将待融合源图像先分解到Riesz变换域，再将Riesz变换域的图像序列进行基于分数阶拉普拉斯算子以及多重调和样条小波变换，使图像表示具有方向特性和多分辨率特性，在Riesz拉普拉斯小波基础上结合脉冲耦合神经网络所具有的全局耦合特性进行融合。经过多聚焦、医学CT和核磁共振、遥感全色及多光谱等影像多种应用场景测试，本发明得到的图像融合结果视觉效果较好，对比度和清晰度上有所提升，且运行效率高。
Description
技术领域
本发明涉及计算机应用领域，具体涉及基于Riesz-Lap变换及PCNN的图像融合方法。
背景技术
为获取某一待研究场景的图像时，由于成像时聚焦的位置不同，该场景图像的不同区域成像的清晰度也就不同，为了得到该场景更为清晰以及信息更多的图像，这就需要用到图像融合的方法去解决该类多聚焦问题。除了多聚焦融合外，还有大量其他的待解决的图像融合问题，随着技术进步成像的传感器的种类也越来越多，不同的传感器获取不同的信息能力也不同，例如遥感卫星成像设备采集到的多波段遥感图像以及全色影像与多光谱影像等，其中全色影像具有更高的空间分辨的结构信息，多光谱则包含了色彩信息，将两者融合后的图像可以更好的与实际地面信息相对应。
图像融合按图像信息描述的层级由低到高依次为像素级融合、特征级融合、决策级融合，其中像素级和特征级融合相关的研究为图像融合的重点与热点。为了更好的提取出图像结构信息，图像融合由空域进而转到变换域，变换域相对于空域直接对图像灰度值进行操作可以更容易分离出边缘细节信息，总体上在变换域基础上进行融合使图片结果更清晰对比度表现更好。变换域相关研究也大多集中在不同多分辨率变换上，有基于金字塔分解、小波变换、超小波变换等。基于小波分解的方法在局部频率分析上表现良好，相对于基于金字塔分解减少了子带系数之间相关性；超小波之一的轮廓波变换具有多方向性相对于小波变换能更好的表示自然场景中奇异的线与面等；非下采样变换(Non-subsampledContourlet Transform,NSCT)用子带分离方式避免上下采样过程造成的复原后频谱混叠，但实际上也增加了数据量和计算量。一系列的改进都有一个共同的目标就是更好的分解表示自然场景图像，相关特征进行聚集便于分离出来。Riesz变换具有保持低级特征空间的一致性，且Riesz变换后具有局部性和方向性。基于Riesz变换在计算机视觉与图像上的应用也有一些研究者在相关领域得到较理想的成果，例如对图像基于Riesz变换图像结构相似度的图像质量评价，以及基于Riesz变换分离出图像特征后图像质量评价，由于Riesz变换具有各向同性，且Riesz变换具有计算相对简便的特点，将图像转换到Riesz变换空间，再基于相位一致性检测出图像边缘，比在希尔伯特变换空间计算相位一致性更高效。
发明内容
针对现有技术的不足，本发明提供了基于Riesz-Lap变换及PCNN的图像融合方法，包括如下步骤：
步骤1，读取要处理的两张图像，分别记为图像A、图像B；
步骤2，分别对图像A、图像B进行Riesz变换，图像A、图像B变换后的图像分别记为图像C、图像D；
步骤3，对图像C、图像D分别进行拉普拉斯小波变换，得到具有多分辨率特性的图像塔形结构；
步骤4，对步骤3得到的图像采用脉冲耦合神经网络点火策略，将不同来源的待融合Riesz金字塔各层进行比较融合；
步骤5，将步骤4得到的图像用拉普拉斯金字塔逆变换回Riesz变换域(详情请见于2012年发表于计算机工程与应用期刊的《基于PCA的拉普拉斯金字塔变换融合算法研究》)；
步骤6，进行Riesz逆变换最终得到融合图像。
步骤2包括如下步骤：
步骤2-1，Riesz变换可看成是一维信号的希尔伯特变换在多维空间的扩展变换，其变换具有同希尔伯特变换类似的性质。希尔伯特变换核频域表示如公式(1)所示，
其中，为一维信号频率，为一维信号的希尔伯特变换核，j为虚数符号，为符号函数，用于取一维信号的符号，为信号频率模值；Riesz变换的频域变换核为一个二维向量(Rx,Ry)，包含两个部分，如公式(2)所示，Rx表示图像垂直方向一阶Riesz变换，Ry表示图像水平方向一阶Riesz变换，其中w为二维向量，w＝(wx,wy)，wx和wy分别为二维向量w沿垂直方向的一维向量和水平方向的一维向量；||w||为向量的模值，当wy＝0时沿着wy＝0剖分即简化为一维的希尔伯特变换；
步骤2-2，将一幅图像(即图像A或者图像B)定义为f，图像上任意一点定义为P＝(x,y)，将一幅图像左上角顶点作为二维坐标原点，每个像素作为1个单位，x、y分别是该图像上的点在垂直方向的坐标值和水平方向上的坐标值，Riesz变换核空域表示如公式(3)所示：
rx、ry分别为Riesz变换核在垂直方向的空域表示和水平方向上的空域表示。将输入的图像进行Riesz变换，如公式(4)所示，
其中，f(P)表示图像上的任意一点P，和分别为图像在垂直方向上的Riesz变换和水平方向上的Riesz变换，fr(P)表示图像在P点的Riesz变换，*表示卷积运算，图像与两个Riesz变换核在空域上分别卷积得到；
步骤2-3，在频域中图像先作傅立叶变换FT(Fourier Transform)，如公式(5)，然后分别与频域Riesz变换核相乘，如公式(6)所示：
其中，F(w)为图像傅立叶变换表示。FR(w)表示图像的频域Riesz变换，Rx表示图像垂直方向一阶Riesz变换，Ry表示图像水平方向一阶Riesz变换，和分别为图像垂直方向和水平方向的频域Riesz变换。Riesz变换的二阶变换表示如下：
其中，将一幅图像定义为f，P表示图像上的任意一点，Rx表示图像垂直方向一阶Riesz变换，Ry表示图像水平方向一阶Riesz变换，Rxx表示图像垂直方向二阶Riesz变换，Ryy表示图像水平方向二阶Riesz变换，Rxy表示图像对角线方向Riesz变换；
步骤2-4，对图像A、图像B分别执行步骤2-1～步骤2-3，图像A、图像B变换后的图像分别记为图像C、图像D。
步骤3包括如下步骤：
步骤3-1，在Riesz变换的基础上，采用拉普拉斯样条小波进行多分辨率分解：分数阶拉普拉斯算子(-Δ)a，其中α∈R+，分数阶拉普拉斯算子是2α阶的各向同性微分算子，其傅立叶域定义为：
其中，R表示实数域，是f(P)的二维傅里叶变换，f(P)表示任意一幅图像，即图像C或者图像D，其中P＝(x,y)，P表示在垂直方向坐标值为x、水平方向坐标值为y的点；
多重调和样条是与分数阶拉普拉斯算子相关的样条函数，γ阶多重样条函数φγ(x)满足：
其中Z表示整数集，φγ(x-k)表示多重样条函数φγ(x)向右移动k个单位，s(x)为样条曲线，式中s[x]＝s(x)|x＝k，s[x]为s(x)在k处的整数采样，k∈Z2，k属于整数集的平方数，对于任意满足条件的k的数都必须满足公式(9)；
φγ(x)的傅立叶响应为：
步骤3-2，根据分数阶拉普拉斯算子和多重样条函数，定义一种对称的分数阶拉普拉斯样条小波：
ψ(x)＝(-Δ)γ/2φ2γ(Dx) (11)
式中ψ(x)表示多重样条曲线，γ＞1/2，φ2γ为2γ阶多重调和样条插值算子，D为下采样矩阵，i,k分别为缩放系数和平移系数，ψ(x)定义的基函数ψi,k(x)则为：
ψi,k(x)＝|det(D)|i/2ψ(Dix-D-1k) (12)
其中，|det(D)|是求矩阵D的行列式的绝对值。对拉普拉斯小波基进行Riesz变换，得到γ阶的复Riesz拉普拉斯小波ψ′(x)，图像经过复Riesz拉普拉斯小波变换得到具有多分辨率特性的图像塔形结构：
ψ′(x)＝Rψ(x)＝R(-Δ)γ/2φ2γ(Dx) (13)。
式中，R表示Riesz变换，Rψ(x)表示对ψ(x)进行Riesz变换。
步骤4包括如下步骤：
步骤4-1，脉冲耦合神经网络(Pulse Coupled Neural Network，PCNN)是一种简化的仿猫视觉原理的神经网络，其为由若干互相连接的神经元组合成的反馈型网络，具有同步脉冲发放效应。
传统脉冲耦合神经网络PCNN的融合规则为：
其中为融合图像的小波系数，和分别为源图像A的小波系数和源图像B的小波系数，为对应小波系数的点火次数，对应小波系数的点火次数；
步骤4-2，基于Riesz变换后源图像A和源图像B的拉普拉斯小波变换各层系数分别为和i表示第i层，对得到的每层小波系数作为输入分别输入到脉冲耦合神经网络网络中，迭代Nmax次后，分别得到源图像A和源图像B每层的点火映射图MFA.i和MFB.i，取两幅图像中各层点火次数较大的像素点作为融合图像的融合系数。
步骤6中，采用乘以Riesz变换共轭系数的方法得到融合图像(详情请见于1981年发表于数学杂志的《关于二重M.Riesz变换》)。
本发明基于Riesz变换具有的保持低级特征、具有局部性、方向性等特点先将图像进行Riesz变换，将这些变换后图像再经过拉普拉斯小波变换后具有多分辨率特性，接下来采用脉冲耦合神经网络点火策略将不同来源的待融合Riesz金字塔各层进行比较融合，最后经过拉普拉斯金字塔逆变换回Riesz变换域，Riesz逆变换(此处采用乘以Riesz变换共轭系数)最终得到融合图像。
有益效果：本发明将图像用Riesz拉普拉斯样条小波金字塔分解得到一种新的表示方式，在该方式下通过脉冲耦合神经网络优化得到较好的融合图像，通过观察以及客观数量分析本发明有效且在信息量、标准差、平均梯度等参数上表现较优，且通过实验统计算法效率上亦较优。通过实验检验本发明的过程中发现本发明在脉冲耦合神经网络参数自适应选择上可以进一步优化，另外间接利用Riesz变换评价图像融合质量的方法作为融合规则中的目标函数来达到更优的融合效果亦具有研究价值。
附图说明
下面结合附图和具体实施方式对本发明做更进一步的具体说明，本发明的上述或其他方面的优点将会变得更加清楚。
图1是一阶和二阶Riesz变换图像示意图。
图2a是多聚焦Clock以脉冲耦合神经网络为融合策略的拉普拉斯金字塔分解复原示意图。
图2b是多聚焦Clock非下采样轮廓波变换分解复原示意图。
图2c是多聚焦Clock基于Riesz分数阶拉普拉斯金字塔分解复原融合示意图。
图3a是医学CT、核磁共振图像以脉冲耦合神经网络为融合策略的拉普拉斯金字塔分解复原示意图。
图3b是医学CT、核磁共振图像非下采样轮廓波变换分解复原示意图。
图3c是医学CT、核磁共振图像以脉冲耦合神经网络为融合策略的拉普拉斯金字塔分解复原示意图。
图4a是遥感全色与多光谱图像以脉冲耦合神经网络为融合策略的拉普拉斯金字塔分解复原示意图。
图4b是遥感全色与多光谱图像非下采样轮廓波变换分解复原示意图。
图4c是遥感全色与多光谱图像以脉冲耦合神经网络为融合策略的拉普拉斯金字塔分解复原示意图。
图5是基于Riesz拉普拉斯小波变换及脉冲耦合神经网络的图像融合流程图。
具体实施方式
下面结合附图及实施例对本发明做进一步说明。
如图5所示，本发明提供了基于Riesz-Lap变换及PCNN的图像融合方法，包括如下步骤：
步骤1，读取要处理的两张图像，分别记为图像A、图像B；
步骤2，分别对图像A、图像B进行Riesz变换，变换后的图像分别记为图像C、图像D；
步骤3，对图像C、图像D分别进行拉普拉斯小波变换，得到具有多分辨率特性的图像塔形结构；
步骤4，对步骤3得到的图像采用脉冲耦合神经网络点火策略，将不同来源的待融合Riesz金字塔各层进行比较融合；
步骤5，将步骤4得到的图像用拉普拉斯金字塔逆变换回Riesz变换域；
步骤6，进行Riesz逆变换最终得到融合图像。
步骤2包括如下步骤：
步骤2-1，Riesz变换可看成是一维信号的希尔伯特变换在多维空间的扩展变换，其变换具有同希尔伯特变换类似的性质。希尔伯特变换核频域表示如公式(1)所示，
其中，为一维信号频率，为一维信号的希尔伯特变换核，j为虚数符号，为符号函数，用于取一维信号的符号，为信号频率模值；Riesz变换的频域变换核为一个二维向量(Rx,Ry)，包含两个部分，如公式(2)所示，Rx表示图像垂直方向一阶Riesz变换，Ry表示图像水平方向一阶Riesz变换，其中w为二维向量w＝(wx,wy)，wx和wy分别为二维向量w沿垂直方向和水平方向的一维向量。||w||为向量的模值，当wy＝0时沿着wy＝0剖分即精简为一维的希尔伯特变换；
步骤2-2，将一幅图像(即图像A或者图像B)定义为f，图像上任意一点定义为P＝(x,y)，将一幅图像左上角顶点作为二维坐标原点，每个像素作为1个单位，x、y分别是该图像上的点在垂直方向的坐标值和水平方向上的坐标值，Riesz变换核空域表示如公式(3)所示：
rx、ry分别为Riesz变换核在垂直方向和水平方向上的空域表示。将输入的图像进行Riesz变换，如公式(4)所示，
其中，f(P)表示图像上的任意一点P，和分别为图像在垂直方向和水平方向上的Riesz变换，fr(P)表示图像在P点的Riesz变换，*表示卷积运算，图像与两个Riesz变换核在空域上分别卷积得到；
步骤2-3，在频域中图像先作傅立叶变换FT(Fourier Transform)，如公式(5)，然后分别与频域Riesz变换核相乘，如公式(6)所示：
F(w)为图像傅立叶变换表示。FR(w)表示图像的频域Riesz变换，Rx表示图像垂直方向一阶Riesz变换，Ry表示图像水平方向一阶Riesz变换，和分别为图像垂直方向和水平方向的频域Riesz变换。Riesz变换的二阶变换表示如下：
其中，将一幅图像定义为f，P表示图像上的任意一点，Rx表示图像垂直方向一阶Riesz变换，Ry表示图像水平方向一阶Riesz变换，Rxx表示图像垂直方向二阶Riesz变换，Ryy表示图像水平方向二阶Riesz变换，Rxy表示图像对角线方向Riesz变换；
一阶特征图能够很好表达图像的边缘轮廓特征，二阶特征图能够表达更复杂的特征，比如角点等，如图1所示，Rx和Rxx突出图像垂直方向轮廓，Ry和Ryy突出图像水平方向轮廓，Rxy则突出对角方向轮廓。高阶Riesz变换具有方向性。
步骤2-4，对图像A、图像B分别执行步骤2-1～步骤2-3，变换后的图像分别记为图像C、图像D。
步骤3包括如下步骤：
步骤3-1，在Riesz变换的基础上，采用拉普拉斯样条小波进行多分辨率分解：分数阶拉普拉斯算子(-Δ)a，其中α∈R+，分数阶拉普拉斯算子是2α阶的各向同性微分算子，其傅立叶域定义为：
其中，R表示实数域，是f(P)的二维傅里叶变换，f(P)表示任意一幅图像，即图像C或者图像D，其中P＝(x,y)，P表示在垂直方向坐标值为x、水平方向坐标值为y的点；
多重调和样条是与分数阶拉普拉斯算子相关的样条函数，γ阶多重样条函数φγ(x)满足：
其中Z表示整数集，φγ(x-k)表示多重样条函数φγ(x)向右移动k个单位，s(x)为样条曲线，式中s[x]＝s(x)|x＝k，s[x]为s(x)在k处的整数采样，k∈Z2，k属于整数集的平方数，对于任意满足条件的k的数都应该满足公式(9)；
φγ(x)的傅立叶响应为：
步骤3-2，根据分数阶拉普拉斯算子和多重样条函数，定义一种对称的分数阶拉普拉斯样条小波：
ψ(x)＝(-Δ)γ/2φ2γ(Dx) (11)
式中ψ(x)表示多重样条曲线，γ＞1/2，在本实施例中选取γ＝2，φ2γ为2γ阶多重调和样条插值算子，D为下采样矩阵，i,k分别为缩放系数和平移系数，在本实施例中取的i,k都为9。ψ(x)定义的基函数ψi,k(x)则为：
ψi,k(x)＝|det(D)|i/2ψ(Dix-D-1k) (12)
其中，|det(D)|是求矩阵D的行列式的绝对值。对拉普拉斯小波基进行Riesz变换，得到γ阶的复Riesz拉普拉斯小波ψ′(x)，图像经过复Riesz拉普拉斯小波变换得到具有多分辨率特性的图像塔形结构：
ψ′(x)＝Rψ(x)＝R(-Δ)γ/2φ2γ(Dx) (13)。
式中，R表示Riesz变换，Rψ(x)表示对ψ(x)进行Riesz变换。
步骤4包括如下步骤：
步骤4-1，脉冲耦合神经网络(Pulse Coupled Neural Network，PCNN)是一种简化的仿猫视觉原理的神经网络，其为由若干互相连接的神经元组合成的反馈型网络，具有同步脉冲发放效应。
传统脉冲耦合神经网络PCNN的融合规则为：
其中为融合图像的小波系数，和分别为源图像A和源图像B的小波系数，为对应小波系数的点火次数，对应小波系数的点火次数；
步骤4-2，基于Riesz变换后源图像A和源图像B的拉普拉斯小波变换各层系数分别为和i表示第i层，对得到的每层小波系数作为输入分别输入到脉冲耦合神经网络网络中，迭代Nmax次后，在本实施例中是迭代了159次，分别得到源图像A和源图像B每层的点火映射图MFA.i和MFB.i，取两幅图像中各层点火次数较大的像素点作为融合图像的融合系数。
步骤6中，采用乘以Riesz变换共轭系数的方法得到融合图像。
图像融合后的质量可以从主观上和客观上来进行评价，为了客观的检验本发明的方法，本发明选取了多组不同应用场景待融合的图片，包含多聚焦图片、医学CT和核磁共振图片、遥感全色与多光谱图片，本发明所选取的图片都是公开的标准数据，以保证对算法稳定性的检验。其中选取的多光谱与全色影像遥感图像具有多波段，故先将包含多通道的多光谱图像的红、绿、蓝(RGB)三色通道转换到色调、饱和度、亮度(HSI)模型表示，将包含主要图像结构信息的亮度分量与全色遥感影像采用本发明所提方法进行融合，再将融合新亮度分量与多光谱的带有光谱特征的色调与饱和度分量变换回红、绿、蓝模型表示。本发明待融合的图像为已经经过初步配准后的图像序列，经过以脉冲耦合神经网络为融合策略的拉普拉斯金字塔分解复原(a)、非下采样轮廓波变换分解复原(b)、Riesz分数阶拉普拉斯金字塔分解复原(c)融合方法对上述多场景应用的测试，图2a是多聚焦Clock以脉冲耦合神经网络为融合策略的拉普拉斯金字塔分解复原示意图。图2b是多聚焦Clock非下采样轮廓波变换分解复原示意图。图2c是多聚焦Clock基于Riesz分数阶拉普拉斯金字塔分解复原融合示意图。图3a是医学CT、核磁共振图像以脉冲耦合神经网络为融合策略的拉普拉斯金字塔分解复原示意图。图3b是医学CT、核磁共振图像非下采样轮廓波变换分解复原示意图。图3c是医学CT、核磁共振图像以脉冲耦合神经网络为融合策略的拉普拉斯金字塔分解复原示意图。图4a是遥感全色与多光谱图像以脉冲耦合神经网络为融合策略的拉普拉斯金字塔分解复原示意图。图4b是遥感全色与多光谱图像非下采样轮廓波变换分解复原示意图。图4c是遥感全色与多光谱图像以脉冲耦合神经网络为融合策略的拉普拉斯金字塔分解复原示意图。
主观上可见拉普拉斯金字塔分解复原、非下采样轮廓波变换以及本文融合方法都成功的融合了两幅图像。无论是针对同类多源图像还是异类多源图像，都取得了不错的效果。但是仔细观察可以发现，多聚焦的Clock图片在细节部分基于拉普拉斯金字塔分解后出现伪影较严重，非下采样轮廓波变换表现良好，而本发明的效果也有提升。而医学图像具有特殊性，基于拉普拉斯分解结果参杂伪影，对医疗诊断不利，价值不大，非下采样轮廓波变换和本发明融合效果相差不大，但本发明结果在清晰度和对比度上更好。遥感图像变换融合后拉普拉斯方法较暗，非下采样轮廓波变换和本发明融合效果都不错。为了更加全面的评价融合结果，需采用各项客观指标来进行评价。
客观上来评价图像融合效果，具有标准确定性。采用下述评价方法对不同的融合方法进行客观评价。
信息熵(Information Entropy,IE)反映了图像信息的丰富程度。
pi表示灰度值为i的像素数Ni与图像总像素数N之比，L为图像的总的灰度级。
标准差(Standard Deviation,SD)表示数据的离散程度，标准差越大表明数据越离散，在图像上反映为标准差越大对比度越大。
p(i,j)表示在图像位置(i,j)的像素灰度值，μ表示整幅图像平均像素灰度值，图像的大小为M×N。
平均梯度(Average Grades,AG)则一定程度反映了图像的清晰度。
p(i,j)表示在图像位置(i,j)的像素灰度值，图像的大小为M×N。
统计结果，其中多光谱与全色影像结果为三通道的均值结果，如表1所示。
表1
表1中LP_PCNN表示拉普拉斯金字塔分解及脉冲耦合神经网络的图像融合方法，NSCT_PCNN表示非下采样轮廓波变换分解及脉冲耦合神经网络的图像融合方法，Riesz-Lap_PCNN表示Riesz分数阶拉普拉斯金字塔分解及脉冲耦合神经网络的图像融合方法，即为本发明。信息熵、标准差、平均梯度的值与融合效果呈正相关，运行时间越短效率越高。在多聚焦和CT、核磁共振图像融合上，本发明各项指标都优于非下采样轮廓波变换方法，特别在时间上快很多，和拉普拉斯分解方法相比，在信息熵上高于该算法。在多光谱与全色图像上，本发明各项指标都高于其他算法。验证了本发明在图像细节提取与保持的优良特性。
本发明提供了基于Riesz-Lap变换及PCNN的图像融合方法，具体实现该技术方案的方法和途径很多，以上所述仅是本发明的优选实施方式，应当指出，对于本技术领域的普通技术人员来说，在不脱离本发明原理的前提下，还可以做出若干改进和润饰，这些改进和润饰也应视为本发明的保护范围。本实施例中未明确的各组成部分均可用现有技术加以实现。
Claims (5)
1.基于Riesz-Lap变换及PCNN的图像融合方法，其特征在于，包括如下步骤：
步骤1，读取要处理的两张图像，分别记为图像A、图像B；
步骤2，分别对图像A、图像B进行Riesz变换，图像A、图像B变换后的图像分别记为图像C、图像D；
步骤3，对图像C、图像D分别进行拉普拉斯小波变换，得到具有多分辨率特性的图像塔形结构；
步骤4，对步骤3得到的图像采用脉冲耦合神经网络点火策略，将不同来源的待融合Riesz金字塔各层进行比较融合；
步骤5，将步骤4得到的图像用拉普拉斯金字塔逆变换回Riesz变换域；
步骤6，进行Riesz逆变换最终得到融合图像。
2.如权利要求1所述的方法，其特征在于，步骤2包括如下步骤：
步骤2-1，希尔伯特变换核频域表示如公式(1)所示，
其中，为一维信号频率，为一维信号的希尔伯特变换核，j为虚数符号，为符号函数，用于取一维信号的符号，为信号频率模值；Riesz变换的频域变换核为一个二维向量(Rx,Ry)，包含两个部分，如公式(2)所示，Rx表示图像垂直方向一阶Riesz变换，Ry表示图像水平方向一阶Riesz变换，其中w为二维向量，w＝(wx,wy)，wx和wy分别为二维向量w沿垂直方向的一维向量和水平方向的一维向量；||w||为向量的模值，当wy＝0时沿着wy＝0剖分即简化为一维的希尔伯特变换；
步骤2-2，将一幅图像定义为f，图像上任意一点定义为P＝(x,y)，将一幅图像左上角顶点作为二维坐标原点，每个像素作为1个单位，x、y分别是该图像上的点在垂直方向的坐标值和水平方向上的坐标值，Riesz变换核空域表示如公式(3)所示：
rx、ry分别为Riesz变换核在垂直方向的空域表示和水平方向上的空域表示；将输入的图像进行Riesz变换，如公式(4)所示，
其中，f(P)表示图像上的任意一点P，和分别为图像在垂直方向上的Riesz变换和水平方向上的Riesz变换，fr(P)表示图像在P点的Riesz变换，*表示卷积运算；
步骤2-3，在频域中图像先作傅立叶变换FT，如公式(5)，然后分别与频域Riesz变换核相乘，如公式(6)所示：
其中，F(w)为图像傅立叶变换表示，FR(w)表示图像的频域Riesz变换，Rx表示图像垂直方向一阶Riesz变换，Ry表示图像水平方向一阶Riesz变换，和FRy分别为图像垂直方向的频域Riesz变换和水平方向的频域Riesz变换；Riesz变换的二阶变换表示如下：
其中，Rx表示图像垂直方向一阶Riesz变换，Ry表示图像水平方向一阶Riesz变换，Rxx表示图像垂直方向二阶Riesz变换，Ryy表示图像水平方向二阶Riesz变换，Rxy表示图像对角线方向Riesz变换；
步骤2-4，对图像A、图像B分别执行步骤2-1～步骤2-3，图像A、图像B变换后的图像分别记为图像C、图像D。
3.如权利要求2所述的方法，其特征在于，步骤3包括如下步骤：
步骤3-1，在Riesz变换的基础上，采用拉普拉斯样条小波进行多分辨率分解：分数阶拉普拉斯算子(-Δ)a，其中α∈R+，分数阶拉普拉斯算子是2α阶的各向同性微分算子，其傅立叶域定义为：
其中，R表示实数域，是f(P)的二维傅里叶变换，f(P)表示任意一幅图像，即图像C或者图像D，其中P＝(x,y)，P表示在垂直方向坐标值为x、水平方向坐标值为y的点；
多重调和样条是与分数阶拉普拉斯算子相关的样条函数，γ阶多重样条函数φγ(x)满足：
其中Z表示整数集，φγ(x-k)表示多重样条函数φγ(x)向右移动k个单位，s(x)为样条曲线，式中s[x]＝s(x)|x＝k，s[x]为s(x)在k处的整数采样，k∈Z2，k属于整数集的平方数，对于任意满足条件的k的数都必须满足公式(9)；
φγ(x)的傅立叶响应为：
步骤3-2，根据分数阶拉普拉斯算子和多重样条函数，定义一种对称的分数阶拉普拉斯样条小波：
ψ(x)＝(-Δ)γ/2φ2γ(Dx) (11)
式中ψ(x)表示多重样条曲线，γ＞1/2，φ2γ为2γ阶多重调和样条插值算子，D为下采样矩阵，i,k分别为缩放系数和平移系数，ψ(x)定义的基函数ψi,k(x)则为：
ψi,k(x)＝|det(D)|i/2ψ(Dix-D-1k) (12)
其中，|det(D)|是求矩阵D的行列式的绝对值，对拉普拉斯小波基进行Riesz变换，得到γ阶的复Riesz拉普拉斯小波ψ′(x)，图像经过复Riesz拉普拉斯小波变换得到具有多分辨率特性的图像塔形结构：
ψ′(x)＝Rψ(x)＝R(-Δ)γ/2φ2γ(Dx) (13)
式中，R表示Riesz变换，Rψ(x)表示对ψ(x)进行Riesz变换。
4.如权利要求3所述的方法，其特征在于，步骤4包括如下步骤：
步骤4-1，脉冲耦合神经网络PCNN的融合规则为：
其中为融合图像的小波系数，和分别为源图像A的小波系数和源图像B的小波系数，为对应小波系数的点火次数，对应小波系数的点火次数；
步骤4-2，基于Riesz变换后源图像A和源图像B的拉普拉斯小波变换各层系数分别为和对得到的每层小波系数作为输入分别输入到脉冲耦合神经网络网络中，迭代Nmax次后，分别得到源图像A和源图像B每层的点火映射图MFA.i和MFB.i，取两幅图像中各层点火次数较大的像素点作为融合图像的融合系数。
5.如权利要求4所述的方法，其特征在于，步骤6中，采用乘以Riesz变换共轭系数的方法得到融合图像。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201910318054.5A CN110163830B (zh) | 2019-04-19 | 2019-04-19 | 基于Riesz-Lap变换及PCNN的图像融合方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201910318054.5A CN110163830B (zh) | 2019-04-19 | 2019-04-19 | 基于Riesz-Lap变换及PCNN的图像融合方法 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN110163830A true CN110163830A (zh) | 2019-08-23 |
| CN110163830B CN110163830B (zh) | 2021-02-19 |
Family
ID=67639853
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201910318054.5A Active CN110163830B (zh) | 2019-04-19 | 2019-04-19 | 基于Riesz-Lap变换及PCNN的图像融合方法 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN110163830B (zh) |
Cited By (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN112927161A (zh) * | 2021-03-12 | 2021-06-08 | 新疆大学 | 一种多光谱遥感图像的增强方法、装置及存储介质 |
| CN113269715A (zh) * | 2021-04-08 | 2021-08-17 | 之江实验室 | 广义Bedrosian准则下的各向同性图像分解方法 |
Citations (4)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN102722877A (zh) * | 2012-06-07 | 2012-10-10 | 内蒙古科技大学 | 基于双通道pcnn的多聚焦图像融合方法 |
| CN105139371A (zh) * | 2015-09-07 | 2015-12-09 | 云南大学 | 一种基于pcnn与lp变换的多聚焦图像融合方法 |
| CN105225213A (zh) * | 2015-09-29 | 2016-01-06 | 云南大学 | 一种基于s-pcnn与拉普拉斯金字塔的彩色图像融合方法 |
| CN109102485A (zh) * | 2018-09-11 | 2018-12-28 | 中国科学院长春光学精密机械与物理研究所 | 基于nsst和自适应双通道pcnn的图像融合方法及装置 |
-
2019
- 2019-04-19 CN CN201910318054.5A patent/CN110163830B/zh active Active
Patent Citations (4)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN102722877A (zh) * | 2012-06-07 | 2012-10-10 | 内蒙古科技大学 | 基于双通道pcnn的多聚焦图像融合方法 |
| CN105139371A (zh) * | 2015-09-07 | 2015-12-09 | 云南大学 | 一种基于pcnn与lp变换的多聚焦图像融合方法 |
| CN105225213A (zh) * | 2015-09-29 | 2016-01-06 | 云南大学 | 一种基于s-pcnn与拉普拉斯金字塔的彩色图像融合方法 |
| CN109102485A (zh) * | 2018-09-11 | 2018-12-28 | 中国科学院长春光学精密机械与物理研究所 | 基于nsst和自适应双通道pcnn的图像融合方法及装置 |
Non-Patent Citations (4)
| Title |
|---|
| BO JIN 等: "Multi-modality image fusion via generalized Riesz-wavelet transformation", 《KSII TRANSACTIONS ON INTERNET AND INFORMATION SYSTEMS》 * |
| CHEN YAO 等: "Multi-sensor Image Fusion Algorithm Based on NSCT and PCNN", 《INTERNATIONAL CONFERENCE ON COMPUTATIONAL ASPECTS OF SOCIAL NETWORKS》 * |
| 李美丽 等: "基于NSCT和PCNN的红外与可见光图像融合方法", 《光电工程》 * |
| 王佺 等: "基于拉普拉斯金字塔与PCNN-SML的图像融合算法", 《计算机科学》 * |
Cited By (4)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN112927161A (zh) * | 2021-03-12 | 2021-06-08 | 新疆大学 | 一种多光谱遥感图像的增强方法、装置及存储介质 |
| CN112927161B (zh) * | 2021-03-12 | 2022-07-01 | 新疆大学 | 一种多光谱遥感图像的增强方法、装置及存储介质 |
| CN113269715A (zh) * | 2021-04-08 | 2021-08-17 | 之江实验室 | 广义Bedrosian准则下的各向同性图像分解方法 |
| CN113269715B (zh) * | 2021-04-08 | 2022-09-09 | 之江实验室 | 广义Bedrosian准则下的各向同性图像分解方法 |
Also Published As
| Publication number | Publication date |
|---|---|
| CN110163830B (zh) | 2021-02-19 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN108830796B (zh) | 基于谱空结合和梯度域损失的高光谱图像超分辨重构方法 | |
| Wang et al. | A comparative analysis of image fusion methods | |
| CN108921809B (zh) | 整体原则下基于空间频率的多光谱和全色图像融合方法 | |
| CN108780570A (zh) | 使用迭代协同滤波的图像超分辨率的***和方法 | |
| CN104867124B (zh) | 基于对偶稀疏非负矩阵分解的多光谱与全色图像融合方法 | |
| CN105719263A (zh) | 基于nsct域底层视觉特征的可见光和红外图像融合算法 | |
| CN109509160A (zh) | 一种利用逐层迭代超分辨率的分层次遥感图像融合方法 | |
| CN108765280A (zh) | 一种高光谱图像空间分辨率增强方法 | |
| CN105550989B (zh) | 基于非局部高斯过程回归的图像超分辨方法 | |
| CN109801250A (zh) | 基于adc-scm与低秩矩阵表达的红外与可见光图像融合方法 | |
| CN109993717A (zh) | 一种结合引导滤波和ihs变换的遥感图像融合方法 | |
| CN103093444A (zh) | 基于自相似性和结构信息约束的图像超分辨重建方法 | |
| CN102446351A (zh) | 一种多光谱与高分辨率全色图像融合方法研究 | |
| CN108629757A (zh) | 基于复数剪切波变换与深度卷积神经网络的图像融合方法 | |
| CN109447909A (zh) | 基于视觉显著性的红外与可见光图像融合方法及*** | |
| CN107169946B (zh) | 基于非负稀疏矩阵与超球面彩色变换的图像融合方法 | |
| CN109886908A (zh) | 红外图像与可见光图像融合方法 | |
| CN110163830A (zh) | 基于Riesz-Lap变换及PCNN的图像融合方法 | |
| CN110111276A (zh) | 基于空-谱信息深度利用的高光谱遥感影像超分辨率方法 | |
| Kotwal et al. | An optimization-based approach to fusion of hyperspectral images | |
| Zhou et al. | PAN-guided band-aware multi-spectral feature enhancement for pan-sharpening | |
| CN109166088B (zh) | 基于非降采样小波变换的双波段灰度熔池图像融合方法 | |
| Lu et al. | Pan-sharpening by multilevel interband structure modeling | |
| CN111563866B (zh) | 一种多源遥感图像融合方法 | |
| Renza et al. | Quality assessment by region in spot images fused by means dual-tree complex wavelet transform |
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 |