彩色图像的超复数傅氏变换和超复数互相关的快速算法

摘要

本发明属于图像信号处理技术领域，具体为一种彩色图像的超复数傅氏变换和超复数互相关的快速算法。超复数把彩色图像作为一个矢量整体进行处理，与传统算法相比，能更好地描述图像的色彩关联，超复数互相关已广泛应用到彩色图像处理的多个领域。目前在彩色图像处理中主要应用“分解变换方法”和“分解垂直平行向量法”，来进行超复数傅氏变换和超复数互相关的快速计算。本发明通过把超复数的信号处理按实部和各个虚部展开，利用传统的快速傅氏变换FFT工具分别计算，再把对应的单位虚向量还原，从而为彩色图像的超复数傅氏变换和超复数互相关提出了一种新的快速算法。分析表明，本算法比目前现有的快速方法更简单易行，且计算量更小。

说明书

技术领域

本发明属于图像信号处理技术领域，具体涉及一种彩色图像的超复数傅氏变换和超复数互相关的快速算法。

背景技术

过去的一个世纪以来，互相关作为一种信号匹配技术，几乎应用到所有的数字信号处理领域，但是其中大多数都局限于标量信号处理。而很多信号处理的应用都涉及到多维矢量，过去的算法常常是对各维信号进行标量处理，然后再简单合成，这样往往达不到很好的信号处理效果。近年来，超复数的应用就是针对多维矢量信号处理问题的。在彩色图像处理中，超复数在多维空间上把彩色图像作为一个矢量整体进行描述，提供了更准确的色彩信息。如今，超复数互相关技术已经应用到彩色图像的图像配准^[1]、光谱分析^[2]、图像数据压缩和边缘检测^[3]等多个领域。因为超复数的一些特殊的数学性质(如：不满足乘法交换律)，所以，过去在使用超复数傅氏变换来快速计算超复数互相关时，需要应用超复数的运算规则重新进行超复数傅氏变换算法，而不能应用传统的复数快速傅氏变换(FFT)工具，很不方便。

T.A.Ell等人提出了对超复数傅氏变换进行分解变换的方法^[4]，把它变成了两个传统复数傅氏变换对，这样就可以快速计算超复数傅氏变换。对于超复数互相关，由于不满足相关定律，所以T.A.Ell等人又提出一种方法^[5]，通过把超复数傅氏变换分解成平行和垂直的分量，来计算超复数互相关。目前在彩色图像处理中主要应用T.A.Ell等人的快速算法。本发明提出了一种彩色图像超复数处理的新的快速算法，把超复数按实部和各个虚部展开，分别进行传统的快速傅氏变换，然后把对应的虚向量还原，从而实现彩色图像的超复数傅氏变换和超复数互相关的快速计算。本发明的算法比目前现有的快速算法更简单易行，且计算量更小。

发明内容

本发明的目的在于提供一种简单易行的计算彩色图像的二维超复数傅氏变换和超复数互相关的快速方法。

超复数把彩色图像作为一个矢量整体进行描述，因而能更好地描述图像的色彩关联。设(m，n)为彩色图像中像素的坐标，则彩色图像RGB(R、G、B分别表示红、绿、蓝分量)模型可以表示为如下的无实部的纯超复数：f(m，n)＝R(m，n)i+G(m，n)j+B(m，n)k。

其中，i、j、k为超复数的虚数单位。

设四元超复数为：q(m，n)＝a(m，n)+b(m，n)i+c(m，n)j+d(m，n)k，超复数乘法满足结合律和分配律，不满足交换律。超复数的二维傅氏右变换^[6]为：

$>F^R(v,u)=\frac{1}{\sqrt{\mathrm{MN}}}\underset{m=0}{\overset{M-1}{\Sigma }}\underset{n=0}{\overset{N-1}{\Sigma }}(a+\mathrm{bi}+\mathrm{cj}+\mathrm{dk})e^{-\mu 2\pi (\frac{\mathrm{mv}}{M}+\frac{\mathrm{nu}}{N})}---\left(1\right)$>

M、N为图像尺寸，即M为像素行坐标m的总数，N为像素列坐标n的总数。

其中μ为单位虚向量，在彩色图像处理中可以取强度图像矢量$>\mu =(i+j+k)/\sqrt{3.}$>在超复数图像空间，任何色彩向量与μ轴的平行分量代表该向量的亮度，与该轴的垂直分量代表该向量的色度^[4]。把(1)式中的指数项移到左边时得到结果是不同的，对应的左变换形式为F^L。

根据分配律，本发明算法把(1)式展开成四个傅氏变换项，例如，其中第二个傅氏变换项为：$>i\left(\frac{1}{\sqrt{\mathrm{MN}}}\underset{m=0}{\overset{M-1}{\Sigma }}\underset{n=0}{\overset{N-1}{\Sigma }}\left({\mathrm{be}}^{-\mu 2\pi (\frac{\mathrm{mv}}{M}+\frac{\mathrm{nu}}{N}}\right)\right),$>b为实数，而μ²＝i²＝-1。设real()表示取复数的实部，image()表示取复数的虚部，我们用传统实数傅氏变换(RFT)来计算上式：上式＝i(real(b_RFT)+μ·imag(b_RFT))。这样，就利用RFT得到了超复数傅氏右变换：

F^R(v，u)＝(real(a_RFT)+μ·imag(a_RFT))+i(real(b_RFT)+μ·imag(b_RFT))    (2)

+j(real(c_RFT)+μ·imag(c_RFT))+k(real(d_RFT)+μ·imag(d_RFT))

类似的，可利用实数的快速傅氏逆变换(IRFT)得到超复数傅氏右逆变换。

F^-R(v，u)＝(real(A_IRFT)+μ·imag(A_IRFT))+i(real(B_IRFT)+μ·imag(B_IRFT)    (3)

+j(real(C_IRFT)+μ·imag(C_IRFT))+k(real(D_IRFT)+μ·imag(D_IRFT))

而要得到超复数的傅氏左变换F^L和左逆变换F^-L，只要把(2)式和(3)式中左乘i，j，k变成右乘i，j，k，即可。

两幅彩色图像的超复数互相关的形式为：

$>\mathrm{cr}(m,n)=f*g=\underset{x=0}{\overset{M-1}{\Sigma }}\underset{y=0}{\overset{N-1}{\Sigma }}f(x,y)\overline{g(x-m,y-n)}---\left(4\right)$>

f、g分为两幅彩色图像的纯超复数模型，为复数g的共轭复数，下同。

它比传统互相关更能体现图像的色彩关联，它表达了两帧图像色彩之间的映射和旋转^[1]。超复数乘法不满足交换律，所以彩色图像超复数互相关不满足相关定律，不能由两幅图像各自的傅氏变换相乘的逆变换得到。

对于RGB彩色图像，我们把式(4)的超复数互相关展开得：

$>\mathrm{cr}(m,n)=f*g=\underset{x=0}{\overset{M-1}{\Sigma }}\underset{y=0}{\overset{N-1}{\Sigma }}[(f_R(x,y)i+f_G(x,y)j+f_B(x,y)k)---\left(5\right)$>

$>(-g_R(x-m,y-n)i-g_G(x-m,y-n)j-g_B(x-m,y-n)k\left)\right]$>

$>\mathrm{cr}(m,n)=\underset{x=0}{\overset{M-1}{\Sigma }}\underset{y=0}{\overset{N-1}{\Sigma }}{[f}_R(x,y)g_R(x-m,y-n)+f_G(x,y)g_G(x-m,y-n)+f_B(x,y)g_B(x-m,y-n)]$>

$>+i\underset{x=0}{\overset{M-1}{\Sigma }}\underset{y=0}{\overset{N-1}{\Sigma }}[f_B(x,y)g_G(x-m,y-n)-f_G(x,y)g_B(x-m,y-n)]---\left(6\right)$>

$>+j\underset{x=0}{\overset{M-1}{\Sigma }}\underset{y=0}{\overset{N-1}{\Sigma }}[f_R(x,y)g_B(x-m,y-n)-f_B(x,y)g_R(x-m,y-n)]$>

$>+k\underset{x=0}{\overset{M-1}{\Sigma }}\underset{y=0}{\overset{N-1}{\Sigma }}[f_G(x,y)g_R(x-m,y-n)-f_R(x,y)g_G(x-m,y-n)]$>

对上述的每一个求和项，可以先用传统实数傅氏变换RFT，再进行复数傅氏逆变换IFFT快速求解。以第2个求和项为例：

这样，我们就利用传统的快速傅氏变换FFT工具，实现了一种很简单易行的彩色图像超复数互相关的快速算法。

由可分析可见，本发明提出的彩色图像的超复数傅氏变换和超复数互相关的快速算法，是把彩色图像的超复数信号处理按实部及各个虚部展开分别进行传统FFT变换，并把对应的单位虚向量还原。其中：

彩色图像的超复数傅氏变换的快速算法，是用纯超复数形式表示彩色图像的RGB模型，把彩色图像的超复数傅氏变换展开成四个傅氏变换项，分别进行传统实数傅氏变换RFT，然后把对应的单位虚向量还原，从而实现了二维超复数傅氏变换的快速运算。

彩色图像的超复数互相关的快速算法，是用纯超复数形式表示彩色图像的RGB模型，把彩色图像的超复数互相关展开成九个求和项，对每一求和项先用传统实数傅氏变换RFT，再进行复数傅氏逆变换IFFT快速求解，从而实现了二维超复数互相关的快速运算。

本发明计算彩色图像的超复数傅氏变换的快速计算的步骤归纳如下：

(1)把彩色图像的RGB模型表示为纯超复数形式：

f(m，n)＝R(m，n)i+G(m，n)j+B(m，n)k，其中，(m，n)为彩色图像中像素的坐标；

(2)分别求出彩色图像的每一色的实数傅氏变换：R_RFT，G_RFT和B_RFT；

(3)计算彩色图像的超复数傅氏变换：

左变换：F^L(v，u)＝(real(R_RFT)+μ·imag(R_RFT))i+(real(G_RFT)+μ·imag(G_RFT))j

+(real(B_RFT)+μ·imag(B_RFT))k

右变换：F^R(v，u)＝i(real(R_RFT)+μ·imag(R_RFT))+j(real(G_RFT)+μ·imag(G_RFT))

+k(real(B_RFT)+μ·imag(B_RFT))

其中，μ为强度图像矢量，$>\mu =(i+j+k)/\sqrt{3.}$>

本发明计算彩色图像的超复数互相关的快速计算的步骤归纳如下：

(1)设两幅彩色图像为f(x，y)和g(x，y)，把它们分别表示为纯超复数形式：f(x，y)＝f_R(x，y)i+f_G(x，y)j+f_B(x，y)k，g(x，y)＝g_R(x，y)i+g_G(x，y)j+g_B(x，y)k；

(2)分别求出f和g的每一色的实数傅氏变换：f_RRFT，f_GRFT，f_BRFT，g_RRFT，g_GRFT和g_BRFT；

(3)彩色图像f和g在(m，n)处的超复数互相关为$>\mathrm{cr}(m,n)=f*g=g\underset{x=0}{\overset{M-1}{\Sigma }}\underset{y=0}{\overset{N-1}{\Sigma }}f(x,y)\overline{g(x-m,y-n)},$>把它按照(6)式展开，从而可以通过快速实数傅氏变换RFT和快速复数傅氏逆变换IFFT来计算超复数互相关：

$>\mathrm{cr}(m,n)=\mathrm{IFFT}[f_{\mathrm{RRFT}}\overline{g_{\mathrm{RRFT}}}+f_{\mathrm{GRFT}}\overline{g_{\mathrm{GRFT}}}+f_{\mathrm{BRFT}}\overline{g_{\mathrm{BRFT}}}]+i·\mathrm{IFFT}[f_{\mathrm{BRFT}}\overline{g_{\mathrm{GRFT}}}-f_{\mathrm{GRFT}}\overline{g_{\mathrm{BRFT}}}]$>

$>+j·\mathrm{IFFT}[f_{\mathrm{RRFT}}\overline{g_{\mathrm{BRFT}}}-f_{\mathrm{BRFT}}\overline{g_{\mathrm{RRFT}}}]+k·\mathrm{IFFT}[f_{\mathrm{GRFT}}\overline{g_{\mathrm{RRFT}}}-f_{\mathrm{RRFT}}\overline{g_{\mathrm{GRFT}}}]$>

具体实施方式

我们以对两幅彩色图像进行图像配准来描述本发明的具体实施方式。

1.设χ＝(x，y)为彩色图像中像素的坐标，把两幅彩色图像f和g中每一像素的RGB模型分别表示为纯超复数形式：f(x，y)＝f_R(x，y)i+f_G(x，y)j+f_B(x，y)k和g(x，y)＝g_R(x，y)i+g_G(x，y)j+g_B(x，y)k。

2.分别求出彩色图像f和g的每一色的实数傅氏变换：f_RRFT，f_GRFT，f_BRFT，g_RRFT，g_GRFT和g_BRFT。

3.计算f和g超复数互相关。我们用g来匹配f，设待匹配图像为g(χ+ξ)，其中ξ＝(m，n)表示匹配点在图像中的相对移动，则两幅图像在(m，n)的超复数互相关cr(m，n)，如(4)式所示，按R，G，B(红、绿、蓝)三色展开如(6)式所示，这样可以通过RFT和IFFT来快速计算超复数互相关：

$>\mathrm{cr}(m,n)=\mathrm{IFFT}[f_{\mathrm{RRFT}}\overline{g_{\mathrm{RRFT}}}+f_{\mathrm{GRFT}}\overline{g_{\mathrm{GRFT}}}+f_{\mathrm{BRFT}}\overline{g_{\mathrm{BRFT}}}]+i·\mathrm{IFFT}[f_{\mathrm{BRFT}}\overline{g_{\mathrm{GRFT}}}-f_{\mathrm{GRFT}}\overline{g_{\mathrm{BRFT}}}]$>

$>+j·\mathrm{IFFT}[f_{\mathrm{RRFT}}\overline{g_{\mathrm{BRFT}}}-f_{\mathrm{BRFT}}\overline{g_{\mathrm{RRFT}}}]+k·\mathrm{IFFT}[f_{\mathrm{GRFT}}\overline{g_{\mathrm{RRFT}}}-f_{\mathrm{RRFT}}\overline{g_{\mathrm{GRFT}}}]$>

4.通过寻找超复数互相关的最大值来进行图像配准。找出使超复数互相关cr(m，n)最大的(m，n)，即为彩色图像f和g的配准点。

目前普遍采用“分解变换方法”^[4]和“分解垂直平行向量法”^[5]来快速计算超复数傅氏变换和超复数互相关，我们通过与它们的比较来阐述本发明的技术效果。

“分解变换方法”^[4]通过另取三个正交的单位虚向量μ₁，μ₂，μ₃，变换四元超复数的表达式为：q＝(a′+μ₁b′)+μ₂(c′+μ₁d′)＝q₁+μ₂q₂，其中a′＝a，$>b^{\prime }=(b+c+d)/\sqrt{3,}$>c′＝(c-d)/，$>d^{\prime }=(2b-c-d)/\sqrt{6,}$>从而把超复数傅氏变换，转变成2次复数FFT运算：

F^R(v，u)＝Q₁(v，u)+μ₂Q₂(v，u)    (7)

＝real(q_1FFT)+μ₁·imag(q_1FFT)+μ₂[real(q_2FFT)+μ₁·imag(q_2FFT)]

对于彩色图像RGB模型，因为无实部，所以计算Q₁(v，u)只要进行b′的实数傅氏变换(RFT)：Q₁(v，u)＝real(b′_RFT)+μ₁·imag(b′_RFT)。

“分解变换方法”计算(7)式共需要2次复数FFT和10MN次实数乘法。由于计算二维FFT需要MN·log₂(MN)次实数乘法[7]，而计算“分解变换”a′，b′，c′，d′，需要3MN次实数乘法，所以“分解变换法”计算超复数傅氏变换的计算量为：

MN·(2log₂(MN)+13)    (8)

类似分析可得，该方法进行超复数傅氏逆变换的计算量相同。

对于彩色图像RGB模型，由于二维RFT算法复杂度为(0.5MN)·log₂(MN)^[7]，所以采用“分解变换法”计算彩色图像的超复数傅氏变换的计算量为：

MN·(1.5log₂(MN)+13)    (9)

彩色图像进行超复数傅氏变换后的结果是四元数，所以用“分解变换法”对其频域的结果进行超复数傅氏逆变换的计算量仍然为(8)式。

若采用本发明方法，计算(μ·实数)只要进行1次实数乘法，所以，计算一个超复数傅氏变换(或逆变换)，需要4次RFT运算和4MN次乘法，计算量为：

MN·(2log₂(MN)+4)    (10)

对于彩色图像RGB模型，因为无实部，则对其时域进行傅氏变换(或逆变换)只要3次RFT运算和3MN次乘法。本发明方法的计算量为：

MN·(1.5log₂(MN)+3)    (11)

用本发明方法对彩色图像经过变换后在频域内的结果进行超复数傅氏变换(或逆变换)的计算量仍然为(10)式。

本发明的算法不需要单位虚向量替换和分解变换，要比T.A.Ell等人的“分解变换法”更加简单易行，而且计算量更小。在进行图像处理时，若图像尺寸为M·N＝320×240，log₂(320×240)＝16.2288，则本发明方法与“分解变换法”的计算量的比值为：

或者：

对于快速计算彩色图像超复数互相关的“分解垂直平行向量法”^[5]，它把纯四元数q分解为与另一纯四元数p平行和垂直的向量：

q_⊥＝(q+pqp)/²，q_⊥⊥p；q_‖＝(q-pqp)/²，q_‖‖p   (14)

通过把超复数傅氏变换分解成和强度图像矢量μ平行和垂直的分量，来计算超复数互相关^[1，5，8]，以下给出两种形式：

$>\mathrm{cr}(m,n)=F^{-R}\left(F^R(v,u)\overline{{G_{\left|\right|}}^R(v,u)}\right)+F^R\left(F^R(v,u)\overline{{G_{\perp }}^R(v,u)}\right)---\left(15\right)$>

$>\mathrm{cr}(m,n)F^{-R}(F^R(v,u)\overline{{G_{\left|\right|}}^L(v,u)}+F^{-R}(v,u)\overline{{G_{\perp }}^L(v,u)})---\left(16\right)$>

由于计算pqp需要3×3+4×3＝21次实数乘法，所以一次分解垂直平行向量需要21次实数乘法。(15)式一共要进行2次彩色图像的超复数傅氏变换，对频域结果进行傅氏变换和傅氏逆变换各1次，2MN次四元数乘法和MN次“分解向量”运算，而两个四元数相乘要进行4×4＝16次实数乘法，综上得，(15)式计算量为：

MN(2×16+21)+MN(7log₂(MN)+52)

                                  (17)

＝MN·(7log₂(MN)+105)

(16)式要进行2次彩色图像的超复数傅氏变换，1次彩色图像的超复数逆傅氏变换，1次对频域结果的傅氏逆变换，2MN次四元数乘法和MN次“分解向量”运算，所以(16)式计算量为：

MN(2×16+21)+MN(6·5log₂(MN)+52)

                                     (18)

＝MN·(6.5log₂(MN)+105)

若采用本发明方法，两个复数相乘要进行2×2＝4次实数乘法。所以求解两幅彩色图像的超复数互相关，需要6次RFT，4次IFFT和9MN次复数乘法，计算量为：

MN(7log₂(MN)+36)    (19)

本发明方法比“分解垂直平行向量法”简单得多，计算量也更少，若图像的尺寸为M·N＝320×240，则两种方法算法计算量的比值为：

形式一：

形式二：

通过以上分析可得，本发明提出的快速算法比目前的计算超复数傅氏变换和超复数互相关快速算法的快速算法更加简单易行而且计算量更小，对于320×240的彩色图像，分析的结果表明：本发明提出的算法分别比目前快速算法的计算量减少了20％和30％。

参考文献

[1]C E Moxey，S J Sangwine，T A Ell.Hypercomplex correlation techniques for vector images[J].IEEETrans on Signal Processing，2003，51(7)：1941-1953.

[2]Ja-Han Chang，Soo-Chang Pei，Jian-Jiun Ding.2D quaternion Fourier spectral analysis and its applications[C].Proceedings of IEEE International Symposium on Circuits and Systems(ISCAS′04)，23-26 May2004，3：III-241-4.

[3]S J Sangwine.Color image edge detector based on QCV[J].IEE Electron.Letter，1998，34(10)：969-971.

[4]T A Ell，S J Sangwine.Decomposition of 2D hypercomplex Fourier transforms into pairs of complexFourier transforms[C].European Signal Processing Conference(EUSIPCO)，2000，151-154.

[5]T A Ell，and S J Sangwine.Hypercomplex Wiener-Khintchine theorem with application to color imagecorrelation[C].IEEE International Conference on Image Processing(ICIP2000)，Sept 10-13 2000，2：792-795.

[6]S J Sangwine，and T.A.Ell，The discrete Fourier transform of a colour image[A].J M Blackledge，M JTurner.Image processing II mathematical methods，algorithms and applications[C]，Chichester，U.K.：Ellis Horwood，2000：430-441.

[7]P Duhamel.Implementation of split-radix FFT algorithms for complex，real and real-symmetric data[J].IEEE Trans on Acoust，Speech，Signal Processing，1986，ASSP-34：285-295.

[8]C E Moxey，T A Ell，S J Sangwine.Hypercomplex operators and vector correlation[C].European SignalProcessing Conference(EUSIPCO)，Toulouse，France，Sept 03-06，2002.