理解傅里叶变换

理解傅里叶变换

0 预备知识

0.1 正交基

如何表示二维坐标系下的一个向量?

使用二维坐标(x,y)。再具体一点,向量应该表示为xi+yjx\cdot\overrightarrow{i}+y\cdot\overrightarrow{j}

其中i\overrightarrow{i}j\overrightarrow{j}表示该二维坐标系下的基向量,其长度为单位长度,且ij1\overrightarrow{i}\cdot\overrightarrow{j}\neq 1.

这是基向量的简单定义。正交基就是一组相互垂直正交的基向量,也就是ii=1\overrightarrow{i}\cdot\overrightarrow{i}= 1ij=0\overrightarrow{i}\cdot\overrightarrow{j}=0

0.2 傅里叶级数

一句话概括:周期性函数f(t)f(t)可以写作一系列正余弦函数的和。

f(t)=a02+ansin(nωt+ϕn)=a02+ansin(nωt)+bncos(nωt)f(t)=\frac{a_0}{2}+\sum a_n\sin(n\omega t+\phi_n)=\frac{a_0}{2}+\sum a_n\sin(n\omega t)+\sum b_n\cos(n\omega t)

其中n是整数表示n倍的频率,ana_nbnb_n是幅值。可以看出f(t)f(t)由三个基向量函数表示分别是常函数1、sin(nωt)\sin(n\omega t)cos(nωt)\cos(n\omega t)

快速理解:正余弦函数是周期性函数,周期性函数相加仍然是周期性函数。那么反过来,一个周期性函数是不是可以由多个周期性函数表示呢?后人得证,可以。因此将正余弦函数作为正交基来表示周期性函数。

那么为什么选择正余弦函数作为正交基?

  • 相同频率、相位情况下,正余弦函数在一个周期内是正交的。这意味着这组正交基可以任意改变频率,而不改变它们正交的性质。从而表示不同频域的信号。
  • 正余弦函数是实数函数且相对简单

为什么说正余弦函数是正交基,正交基不应该是向量吗?

通常情况下,我们将正交基(orthogonal basis)用于描述向量空间中的一组基向量,使得这些基向量两两之间的内积为零。但是在泛函分析中,我们也可以将正交基推广到函数空间中,特别是在描述傅里叶级数和傅里叶变换中。

在泛函分析中,我们考虑的不再是有限维向量空间,而是无穷维函数空间。在函数空间中,我们可以将函数视为向量,因此我们可以将一组函数视为该空间的基。这些函数之间的内积定义为函数之间的积分,而不再是向量的点积。

两个函数f(x)f(x)g(x)g(x)在区间[a,b]上内积(点积)可以表示为:

<f,g>=abf(x)g(x)dx<f,g>=\int_a^b f(x)\cdot g(x)\text{d}x

如果一组函数在某个特定区间内的内积等于零,则它们被称为正交的。

因此,虽然正交基的概念最初是从向量空间中衍生出来的,但在泛函分析中,我们也可以将其应用于函数空间,并将满足正交性质的函数称为正交基。在傅里叶变换中,正弦函数和余弦函数就是函数空间中的一组正交基,它们用于将函数分解为频域上的正弦和余弦分量。

1 一维傅里叶变换

上面提到傅里叶级数只说了周期性函数f(t)f(t)可以写作一系列正余弦函数的和。

那么非周期的函数怎么表示?可以看作周期无穷大的周期函数

1.1 欧拉公式

eiθ=cosθ+isinθe^{i\theta}=\cos\theta+i\sin\theta

将动图中的t理解为角度θ\theta更为合适,然后令θ=ωt\theta=\omega t,其中tt时间,ω\omega表示角速度也可以理解为频率。

为什么ω\omega表示角速度也可以理解为频率?从动图中可以看出eiθe^{i\theta}在实数轴和虚数轴的投影分别是sin和cos函数,角速度越快,sin和cos的周期越短,频率越快。所以角速度和频率等价。换言之T=2π/ωT=2\pi/\omega,这里带入加速度和频率是同一个意思。更严谨的说,这里的ω\omega应该称为角频率。而频率f=ω/2πf=\omega/2\pi

引入欧拉公式的目的:eiωte^{i\omega t}表示一组正余弦函数,不同的ω\omega表示不同组的正余弦函数。从而通过eiωte^{i\omega t}对非周期函数进行解耦操作,也就是分解出非周期函数中频率与eiωte^{i\omega t}相同的部分。

具体怎么操作?

内积

1.2 傅里叶变换

B站首发!草履虫都能看懂的【傅里叶变换】讲解,清华大学李永乐老师教你如何理解傅里叶变换,辨清美颜和变声原理,!!_哔哩哔哩_bilibili

将目标函数f(t)f(t)ejωte^{-j\omega t}做内积。已知f(t)f(t)是由一系列不同ω\omegaejωte^{-j\omega t}组成的函数。将f(t)f(t)与某个特定ω\omegaejωte^{-j\omega t}做内积,f(t)f(t)如果不含该ω\omega分量,则f(t)f(t)与该特定ω\omegaejωte^{-j\omega t}正交,内积为零。

F^T(ω)=+f(t)ejωtdt{=0,ω0,ω\hat{F}_T(\omega)=\int_{-\infty}^{+\infty}f(t)\cdot e^{-j\omega t}\text{d}t \begin{cases} =0&, 不含\omega分量\\ \neq 0&, 含\omega分量 \end{cases}

该公式就是傅里叶变换的过程,对f(t)f(t)在不同频率上解耦。F^T(ω)\hat{F}_T(\omega)是傅里叶变换后的复数,在实部表示幅值、虚部表示相位。

为什么说“f(t)f(t)如果不含该ω\omega分量,则f(t)f(t)与该特定ω\omegaejωte^{-j\omega t}正交”?

f(t)f(t)可以表达为不同频率的复指数函数ejωte^{-j\omega t}的线性组合。

那么”f(t)f(t)与该特定ω\omegaejωte^{-j\omega t}正交“就变为”不同频率的复指数函数ejωte^{-j\omega t}与某一特定频率的不同频率的复指数函数ejωte^{-j\omega t}做正交“。

而不同ω\omega的复指数函数ejωte^{-j\omega t}之间相互正交,为什么?

+ejω1tejω2tdt=+ej(ω1ω2)tdt\int_{-\infty}^{+\infty}e^{-j\omega_1 t}e^{j\omega_2 t}\text{d}t=\int_{-\infty}^{+\infty}e^{-j(\omega_1-\omega_2)t}\text{d}t

如果ω\omega不同,则积分的对象是一个周期函数,周期函数在无穷域上积分为0。

如果ω\omega相同,则积分对象是一个常数,在无穷域上积分是无穷大。

为什么上述公式的ω2\omega_2的符号是正的?

当我们谈论两个复指数函数的内积时,实际上是要考虑一个函数与另一个函数的复共轭的乘积。而ejω2te^{-j\omega_2 t}的复共轭就是ejω2te^{j\omega_2 t}。所以傅里叶变换实际是f(t)f(t)ejωte^{j\omega t}做内积,但是ejωte^{j\omega t}是复指数函数,所以带入公式中就变成了ejωte^{j\omega t}的复共轭ejωte^{-j\omega t}

f(t)f(t)可以表达为不同频率的复指数函数ejωte^{j\omega t}的线性组合。如果f(t)f(t)含有ω\omega分量,那么+f(t)ejωtdt\int_{-\infty}^{+\infty}f(t)\cdot e^{-j\omega t}\text{d}t是不是可以看作+aejωtejωtdt\int_{-\infty}^{+\infty}a\cdot e^{-j\omega t}\cdot e^{-j\omega t}\text{d}t?但是这样的计算结果不就是无穷大吗?

理论上确实是这样的。

但是实际中我们通常处理的是有限时间内的信号或是周期性信号。因此,实际计算中,积分是在有限区间或周期上进行的。

因此,当f(t)f(t)包含ω\omega频率分量时,傅里叶变换在频率ω\omega处的计算结果是该频率分量的幅度和相位的复数表示,而不是无限大。这个结果提供了信号在特定频率ω\omega上的特征,例如幅度a|a|和相位arg(a)\arg(a)

傅里叶逆变换写为:

f(t)=+F^T(ω)ejωtdωf(t)=\int_{-\infty}^{+\infty}\hat{F}_T(\omega)\cdot e^{j\omega t}\text{d}\omega

以上,为最简单的傅里叶变换理解。公式存在些许错误,为了方便理解,不进行纠正。

1.3 希尔伯特空间坐标变换

卷积神经网络的底层是傅里叶变换,傅里叶变换的底层是希尔伯特空间坐标变换_哔哩哔哩_bilibili

傅里叶变换的底层理论是希尔伯特空间坐标变换。因此,从希尔伯特空间这个角度重新理解傅里叶变换。

先从二维曲线与高维向量的联系说起。

对于二维坐标系下的曲线,可以写作一系列点的集合C={(xi,yi)xi,yiR,iZ}C=\{(x_i,y_i)|x_i,y_i\in\R,i\in\Z\}。如果将xix_i作为数组的索引,yiy_i作为索引对应的值,那么CC可以看作是一个一维的长数组。现在有一个高维空间的点,该点的维度与该数组的长度一致,该点的坐标对应该数组的值。那么该高维点与二维曲线等价。将点变为向量,也就是说,一个高维向量与一个二维曲线是等价的。

介绍高维向量与二维曲线的联系有什么用?

傅里叶变换就是高维向量在希尔伯特空间的一种坐标变换。

通俗一点就是将高维向量用一个新的坐标系表示,再进一步说,就是使用一组基向量dn\overrightarrow{d}_n来表示这个高维向量。

那么怎么表示呢?就是将向量投影到基向量dn\overrightarrow{d}_n上获得un\overrightarrow{u}_n

公式化表示:

|\overrightarrow{u}_n|\cdot\overrightarrow{d}_n &=&|\overrightarrow{f}|\cdot|\overrightarrow{d}_n|\cdot\cos(\theta)\cdot\frac{\overrightarrow{d}_n}{|\overrightarrow{d}_n|}\\ &=&\langle\overrightarrow{f},\overrightarrow{d}_n\rangle\cdot\hat{d}_n\\

这个公式是在一个无穷维的空间里讨论一个无穷维的向量。如果在这个无穷维的空间中,向量的内积是完备的,那么这个空间就是希尔伯特空间。相当于将欧几里得空间变为了无穷维。

将向量内积<><\cdot>这个符号展开,也就是向量间逐元素相乘再相加,公式就变为了:

undn=+f(t)dn(t)dtd^n|\overrightarrow{u}_n|\cdot\overrightarrow{d}_n=\int_{-\infty}^{+\infty}f(t)d_n(t)\text{d}t\cdot\hat{d}_n

上图是希尔伯特空间坐标变换的通用写法。

前文提到,傅里叶变换就是高维向量在希尔伯特空间的一种坐标变换。

当变量nn变为角频率ω\omega,基向量dn\overrightarrow{d}_n变为ejωte^{j\omega t}时,这个变换就是傅里叶变换,见下图。

进一步的说,基向量d^n\hat{d}_n同样是一个无穷维的向量。因此,它也可以表示一个二维曲线。对于标准正交坐标系,它对应的是一个宽度极小的峰值为1的方波。对于ejωte^{j\omega t},则是一组正余弦曲线,见1.1节动图。

使用ejωte^{j\omega t}作为基向量,获得的特征,也就是F(n)F(n)表示的是全局特征。

但是有的时候,全局的特征并不好。因为图像的一些细微改变,都会使特征获得非常大的变化。

因此,怎么获得一些局部的特征呢?

这需要先考虑为什么使用ejωte^{j\omega t}作为基向量获得的是全局特征。因为ejωte^{j\omega t}本身的作用域就是全局的,从负无穷到正无穷。因此,一个简单的想法就是将ejωte^{j\omega t}分段,也就是加一个窗口。

这种有窗口范围的函数应该怎么表示呢?

其中ss表示窗口位置,[1,1][-1,1]表示窗口的范围。定义窗口范围相当于一个基向量分为了许多个分段的基向量,从而增加基向量个数。

但是这种分段函数不是处处可微的,因此使用高斯函数进一步修改。

其中aa是方差,表示窗口的范围,ss是期望,表示窗口的位置。

使用gag_a窗口函数表示的基向量,如下图所示

而这个使用高斯窗口分段基向量的变换称为Gabor变换。

此外,窗口函数还有其他的类型。

如果窗口函数是指数函数,那么这个变换称为拉普拉斯变换。

如果窗口函数随频率动态变换,那么这个变换称为小波变换。

以上,是从希尔伯特空间角度出发的理解。

1.4 傅里叶变换的幅度和相位

在一维傅里叶变换中,对于频谱F(ω)F(\omega)需要明确的三个参数是频率ω\omega、振幅AA和相位ϕ\phi

频率、相位、时间、基向量和目标信号的关系如下图所示:

对信号f(t)f(t)进行傅里叶变换可以获得复数形式的频谱F(ω)F(\omega)

F(ω)=+f(t)ejωtdtF(\omega)=\int_{-\infty}^{+\infty}f(t)\cdot e^{-j\omega t}\text{d}t

这里,F(ω)F(\omega)是一个复数形式,它可以表示为:

F(ω)=A(ω)ejϕ(ω)F(\omega)=A(\omega)e^{j\phi(\omega)}

其中,A(ω)A(\omega)是振幅谱,ϕ(ω)\phi(\omega)是相位谱。

更方便理解的求解方法是,将F(ω)F(\omega)分解为实部和虚部:

F(ω)=Re(F(ω))+jIm(F(ω))F(\omega)=\text{Re}(F(\omega))+j\text{Im}(F(\omega))

而振幅就是Re(F(ω))\text{Re}(F(\omega))

相位ϕ(ω)\phi(\omega) 可以通过实部和虚部计算得出:

ϕ(ω)=tan1(Im(F(ω))Re(F(ω)))\phi(\omega)=\tan^{-1}(\frac{\text{Im}(F(\omega))}{\text{Re}(F(\omega))})

那么F(ω)F(\omega)怎么分解为实部和虚部呢?

已知F(ω)=+f(t)ejωtdtF(\omega)=\int_{-\infty}^{+\infty}f(t)\cdot e^{-j\omega t}\text{d}t,根据欧拉公式,ejωte^{-j\omega t}可以表示为cos(ωt)jsin(ωt)\cos(\omega t)-j\sin(\omega t)

因此,将积分拆解为两部分:

Re(F(ω))=+f(t)cos(ωt)dt\text{Re}(F(\omega))=\int_{-\infty}^{+\infty}f(t)\cdot \cos(\omega t)\text{d}t

Im(F(ω))=+f(t)sin(ωt)dt\text{Im}(F(\omega))=-\int_{-\infty}^{+\infty}f(t)\cdot \sin(\omega t)\text{d}t

2 二维傅里叶变换

形象理解二维傅里叶变换 - 知乎 (zhihu.com)

一维傅里叶变换是以多组正余弦函数为基向量进行变换,这些正余弦函数实际是不同频率、相位、振幅的波。各种各样的波叠在一起形成了实际的信号。

一维信号通过波来表示,那么二维图像通过什么来表示?答案是复平面波ej2π(ux+vy)e^{j2\pi(ux+vy)}

一维傅里叶变换公式:

F(ω)=+f(t)ejωtdtF(\omega)=\int_{-\infty}^{+\infty}f(t)\cdot e^{-j\omega t}\text{d}t

类比一维傅里叶变换公式,将时序信号自变量tt变为空域图像的像素坐标(x,y)(x,y),将一维波复平面波ej2πte^{j2\pi t}变为二维复平面波ej2π(ux+vy)e^{j2\pi(ux+vy)}。那么,二维傅里叶变换公式写为:

F(u,v)=++f(x,y)ej2π(ux+vy)dxdyF(u,v)=\int_{-\infty}^{+\infty}\int_{-\infty}^{+\infty}f(x,y)\cdot e^{-j2\pi(ux+vy)}\text{d}x\text{d}y

2.1 二维频率域K-Space

对于正弦平面波,可以这样理解,在一个方向上存在一个正弦函数,在法线方向上将其拉伸。前面说过三个参数可以确定一个一维的正弦波。哪几个参数可以确定一个二维的正弦平面波呢?答案是四个,其中三个和一维的情况一样(频率 𝑤𝑤 ,幅度 𝐴𝐴 ,相位 𝜑),但是具有相同这些参数的平面波却可以有不同的方向 n\overrightarrow{n} 。如下图所示:

类比一维中,幅度和相位可以用一个复数表示,它可以作为我们存储的内容。但是还有两个:一个频率一个方向。这时想到向量是有方向的,也是有长度的。所以我们用一个二维的矩阵的来保存分解之后得到的信息。这个矩阵就是K空间。(一般用k来表示空间频率,单位是1/m)

什么意思呢?就是说一个二维矩阵点 (𝑢,𝑣)(𝑢,𝑣)代表这个平面波的法向量 n\overrightarrow{n},这个向量的模u2+v2\sqrt{u^2+v^2}代表这个平面波的频率ω\omega,这个点里面保存的内容复数就是此平面波的幅度和相位。下面这个图很好的体现了这一点:

再如下面这个图片,中心低频贡献了图像的主体,周围高频提供图像的细节和边缘

因此,k空间的每一个位置存储的数代表了所在位置复平面波在图像中占多少成分,我们就可以用每个系数X所代表的平面波相加得到原来的图像,也就是下图。所以k空间和对应图像储存的信息含量是一样的,只不过表现形式不同,或者说基不同。

2.2 K-Space的一些特性

2.2.1 FOV、分辨率与K空间的关系

在数字图像中,数据都是离散的。也就涉及到采样的问题,和一维一样,如果采样率过低,k空间就会混叠。同时在k空间中采样过低,图像也会混叠。

FOV和分辨率在k空间和图像中是相反的关系。也就是:

Δx=12kxmax,Δk=12xxmax=1FOV\Delta x=\frac{1}{2*k_{xmax}},\Delta k=\frac{1}{2*x_{xmax}}=\frac{1}{FOV}

2.2.2 旋转不变性

从平面波的角度很容易理解,旋转没有改变平面波的幅度相位,只是将所有的平面波都旋转了一个角度。二维傅里叶变换中,实空间旋转多少,频率空间也会相应旋转多少。这其实是高维傅里叶变换缩放定理的一种特殊情况。(连续的是可以证明的,离散的涉及插值 ,不一定完全准确)

3 Gabor变换

【Gabor滤波器】提取图像纹理(Python、C++两种实现)_gabor滤波提取图像的纹理特征向量-CSDN博客

Gabor滤波器_gabor函数 wavelength-CSDN博客

标准的Gabor滤波器及Log_Gabor滤波器的实现、解析、速度优化及其和Halcon中gen_gabor的比较。_log-gabor滤波器-CSDN博客

从傅里叶(Fourier)变换到伽柏(Gabor)变换再到小波(Wavelet)变换_gabor变换-CSDN博客

Gabor是一个用于边缘提取的线性滤波器,其频率和方向表达与人类视觉系统类似,能够提供良好的方向选择和尺度选择特性,而且对于光照变化不敏感,因此十分适合纹理分析。

图中第一行是脊椎动物的视觉响应,第二行是Gabor滤波器的响应,第三行是两者之差,可以看到,二者相差极小。

基于以上特性,Gabor滤波器被广泛应用于人脸识别的预处理。

我们知道,傅里叶变换可以将信号从时域转换到频域,但无法获得频谱中不同频率之间的先后关系。
然而实际应用中我们更多的关心信号局部范围内的的特性,Gabor和小波变换突破了傅里叶变换的局限性

Gabor变换是D.Gabor于1946年提出的,为了提取傅里叶变换的局部信息,引入了时间局部化的窗函数(把信号划分成许多小的时间间隔,用傅里叶变换分析每一个间隔)。因此Gabor变换又称为窗口傅里叶变换(短时傅里叶变换)。

3.1 Gabor滤波器公式推导

首先看一下Gabor函数长什么样,再进行推导:

g(x,y;λ,θ,ϕ,σ,γ)=e(x2+γ2y22σ2)ei(2πfx+ϕ)g(x,y;\lambda,\theta,\phi,\sigma,\gamma)=e^{\Big( -\frac{x'^2+\gamma^2 y'^2}{2\sigma^2}\Big)}\cdot e^{i(2\pi f x'+\phi)}

Gabor滤波器是使用Gabor函数进行滤波。Gabor函数是二维正弦函数和高斯函数的叠加。

基本上,在空域他的形状就是一些有间隔的白色过度条,在频域,则基本为两处白色亮点,如下图所示:

二维正弦函数写作:

f(x,y)=Acos(kxx+kyy+ϕ)=Acos(kr+ϕ)f(x,y)=A\cos(k_x x+k_y y+\phi)=A\cos(\mathbf{k}\cdot\mathbf{r}+\phi)

其中AA是振幅、kxk_xkyk_y分别是在x和y方向上的波数,ϕ\phi是相位。k\mathbf{k}是波矢向量、r=(x,y)\mathbf{r}=(x,y)是方向向量

波矢k=(kx,ky)=k(cos(θ),sin(θ))\mathbf{k}=(k_x,k_y)=k(\cos(\theta),\sin(\theta)),其中θ\theta是波的方向,kk是波矢的模长。

波矢顾名思义,波矢的模长是什么东西?使用波长来理解,就是2π2\pi内有多少个波,写作

k=2πλk=\frac{2\pi}{\lambda}

实际的Gabor函数中没有显式表达λ\lambda,类比傅里叶变换,Gabor使用频率ff来定性的表达波长λ\lambda

λ=vf\lambda=\frac{v}{f}

vv是波速,但是在空域傅里叶变换没有波速这个概念,所以可以将波速忽略,定性的理解为频率ff与波长λ\lambda成反比,带入二维正弦波函数中,整理为:

f(x,y)=Acos(2πf(xcos(θ)+ysin(θ))+ϕ)f(x,y)=A\cos(2\pi f(x \cos(\theta)+y \sin(\theta))+ \phi)

二维高斯函数写为:

g(x,y)=e(x2+y22σ2)g(x,y)=e^{\Big( -\frac{x^2+y^2}{2\sigma^2}\Big)}

再次说明,Gabor函数是二维正余弦波函数与高斯函数的叠加,所以Gabor函数可以写为:

g(x,y)=e(x2+y22σ2)Acos(2πf(xcos(θ)+ysin(θ))+ϕ)g(x,y)=e^{\Big( -\frac{x^2+y^2}{2\sigma^2}\Big)}\cdot A\cos(2\pi f(x \cos(\theta)+y \sin(\theta))+ \phi)

这距离实际的Gabor函数还差了一个空间纵横比γ\gamma,该参数是用来调整高斯函数的形状,也就是相对x,y的变化比例。空间纵横比γ\gamma等于1,则是圆形,否则是椭圆。将空间纵横比γ\gamma加入公式中,写为:

g(x,y)=e(x2+γ2y22σ2)Acos(2πf(xcos(θ)+ysin(θ))+ϕ)g(x,y)=e^{\Big( -\frac{x^2+\gamma^2 y^2}{2\sigma^2}\Big)}\cdot A\cos(2\pi f(x \cos(\theta)+y \sin(\theta))+ \phi)

其中高斯函数因为加入了空间纵横比γ\gamma而产生了方向性,而高斯函数的方向应该与波函数的方向保持一致,也就是同样应该旋转θ\theta角,因此高斯函数变为了:

g(x,y)=exp((xcos(θ)+ysin(θ))2+γ2(xsin(θ)+ycos(θ))22σ2)=e(x2+γ2y22σ2)g(x,y)=\exp\Big( -\frac{(x \cos(\theta)+y \sin(\theta))^2+\gamma^2 (x \sin(\theta)+y \cos(\theta))^2}{2\sigma^2}\Big)=e^{\Big( -\frac{x'^2+\gamma^2 y'^2}{2\sigma^2}\Big)}

其中x=xcos(θ)+ysin(θ)x'=x \cos(\theta)+y \sin(\theta)y=xsin(θ)+ycos(θ)y'=-x \sin(\theta)+y \cos(\theta),表示旋转后的坐标(这里是极坐标表示,与笛卡尔坐标系的方向相反)。重新带入Gabor函数,得到:

g(x,y)=e(x2+γ2y22σ2)Acos(2πfx+ϕ)g(x,y)=e^{\Big( -\frac{x'^2+\gamma^2 y'^2}{2\sigma^2}\Big)}\cdot A\cos(2\pi fx'+ \phi)

使用Gabor函数进行傅里叶变换操作的话,Gabor函数应该包含一组基,根据欧拉公式:

cosω+isinω=eiω\cos\omega+i\sin\omega=e^{i\omega}

则,上面的公式只表示了实数部分,将实部和虚部整合到一起,Gabor函数再次修改为:

g(x,y)=Ae(x2+γ2y22σ2)ei(2πfx+ϕ)g(x,y)=A e^{\Big( -\frac{x'^2+\gamma^2 y'^2}{2\sigma^2}\Big)}\cdot e^{i(2\pi fx'+ \phi)}

也就是:

g(x,y)=Aexp(x2+γ2y22σ2)exp(i(2πfx+ϕ))g(x,y)=A \exp\Big( -\frac{x'^2+\gamma^2 y'^2}{2\sigma^2}\Big)\cdot \exp(i(2\pi fx'+ \phi))

其中:

  • x=xcos(θ)+ysin(θ)x'=x \cos(\theta)+y \sin(\theta)
  • y=xsin(θ)+ycos(θ)y'=-x \sin(\theta)+y \cos(\theta)
  • θ\theta为波函数的方向
  • γ\gamma为空间纵横比,用来调整高斯函数的形状,等于1,则是圆形,否则是椭圆
  • σ\sigma为高斯函数标准差,用来调整高斯函数的范围,也就是滤波器的大小
  • ff为波函数的频率
  • ϕ\phi为相位偏移

此时的Gabor函数是中心处于原点的样子,假设x0,y0x_0,y_0是Gabor函数的中心点,则应该做如下修改:

  • x=(xx0)cos(θ)+(yy0)sin(θ)x'=(x-x_0) \cos(\theta)+(y-y_0) \sin(\theta)
  • y=(xx0)sin(θ)+(yy0)cos(θ)y'=-(x-x_0) \sin(\theta)+(y-y_0) \cos(\theta)

3.2 Gabor滤波器的图像理解

图2对比了不同参数的Gabor滤波器,从第1到第5行λ\lambda取3、6、9、12、15,从第1到第8列θ\theta0,π8,2π8,3π8,4π8,5π8,6π8,7π80,\frac{\pi}{8},\frac{2\pi}{8},\frac{3\pi}{8},\frac{4\pi}{8},\frac{5\pi}{8},\frac{6\pi}{8},\frac{7\pi}{8},其他参数保持不变:ϕ=0,σ=2π,γ=0.5\phi=0,\sigma=2\pi,\gamma=0.5

使用上述滤波器对以下两张人脸进行滤波

得到以下两个结果:

3.3 Gabor变换

Gabor 特征总结 (mengqi92.github.io)

Gabor transform - Wikipedia

回顾一维傅里叶变换的形式:

F(ω)=+f(t)ejωtdtF(\omega)=\int_{-\infty}^{+\infty}f(t)\cdot e^{-j\omega t}\text{d}t

Gabor变换是D.Gabor于1946年提出的,为了提取傅里叶变换的局部信息,引入了时间局部化的窗函数(把信号划分成许多小的时间间隔,用傅里叶变换分析每一个间隔)。因此Gabor变换又称为窗口傅里叶变换(短时傅里叶变换)。

也就是说Gabor变换就是在傅里叶变换的基础上加入了”窗口“这个概念。

怎么给基函数ejωte^{-j\omega t}加上窗口?使用一个函数g(t)g(t)与之相乘,这个函数在窗口范围内具有值,范围外为0,从而将基函数进行窗口划分。可参考1.3节。

于是窗口傅里叶变换写为:

F(ω)=+f(t)g(t)ejωtdtF(\omega)=\int_{-\infty}^{+\infty}f(t)\cdot g(t)\cdot e^{-j\omega t}\text{d}t

而Gabor变换就是使用高斯函数来作为窗口函数,一维高斯函数写作:

g(t)=exp(t22σ2)g(t)=\exp\Big(-\frac{t^2}{2\sigma^2} \Big)

窗口函数需要确定窗口的函数,所以使用高斯函数的期望来确定高斯函数中心位置:

g(ts)=exp((ts)22σ2)g(t-s)=\exp\Big(-\frac{(t-s)^2}{2\sigma^2} \Big)

其中σ\sigma是方差,表示窗口的范围,ss是期望,表示窗口的位置。

那么二维Gabor变换类比1.3和3.1节,自行推导。

3.4 Log-Gabor滤波器

Python OpenCV实现Log Gabor滤波器(由LGHD描述符扩展)_log-gabor-CSDN博客

Log Gabor filter - 维基百科,自由的百科全书 (wikipedia.org)

LGHD(Log-Gabor Histogram Descriptor)
描述符的思想是基于高频分量分布的描述符对于不同的非线性强度变化具有鲁棒性。可以理解为虽然非线性强度差异可以影响图像,但是场景中所包含对象的形状的整体外观仍趋于保持恒定。Log-Gabor 滤波器是根据Gabor滤波器衍生而来,Gabor变换能达到时频局部化的目的,它能够在整体上提供信号的全部信息又能提供在任意局部时间内信号变化剧烈程度的信息。

Log-Gabor滤波器
Log-Gabor 滤波器的表达是从Gabor衍生而来的,相比于Gabor滤波器,它的传递函数是对数频率尺度下的高斯函数,始终没有直流分量,在图像处理时可以不受亮度的影响。
并经研究Log-Gabor函数更符合哺乳动物的视觉观察系统,有更优的纹理提取、目标检测效果。

3.4.1 一维Log-Gabor滤波器

Log-Gabor滤波器具体是怎么操作的?就是在高斯函数的基础上使用了对数log,具体表示为:

g(f)=exp((log(f)log(f0))22(log(σ)log(f0))2)=exp((log(f/f0))22(log(σ/f0))2)g(f)=\exp\Big(-\frac{(\log(f)-\log(f_0))^2}{2(\log(\sigma) - \log(f_0))^2} \Big)=\exp\Big(-\frac{(\log(f/f_0))^2}{2(\log(\sigma/f_0))^2} \Big)

但是根据高斯函数的定义,分母应该是2(log(σ))22(\log(\sigma))^2而不是2(log(σ)log(f0))22(\log(\sigma) - \log(f_0))^2

是的,这样看确实应该是2(log(σ))22(\log(\sigma))^2。但是窗口函数或者滤波函数应该在函数中心f0f_0进行移动的情况下保证对于不同的尺度(带宽σ\sigma)具有相似的响应。见下图意会

左上图表示Gabor函数在频率为0时具有非零响应,这意味着什么?频率为0也就是直流分量,通常值图像的平均亮度或者亮度偏移。

频率为0处具有非零响应,意味着,无论Gabor函数的范围、位置的值是多少,Gabor函数始终对图像的平均亮度具有响应值。

而Log-Gabor最直接的目的就是使Gabor去除这个性质,也就是将频率放入对数域。这样做的好处就是将函数的响应值变为了左下图的样子。在将频率进行对数拉伸(右下图),发现Log-Gabor在对数频率域具有相似的响应值。

好的,回到刚才那个问题,但是根据高斯函数的定义,分母应该是2(log(σ))22(\log(\sigma))^2而不是2(log(σ)log(f0))22(\log(\sigma) - \log(f_0))^2,为什么?

观测左下图,随着中心频率的增加,Log-Gabor函数的带宽是增加的。所以Log-Gabor函数的带宽需要使用2(log(σ)log(f0))22(\log(\sigma) - \log(f_0))^2定义,使得在中心频率f0f_0移动的同时使得相同的σ\sigma具有相似的响应值。

3.4.2 二维Log-Gabor滤波器

二维Log-Gabor滤波器在一维Log-Gabor滤波器的基础上加入了方向滤波。具体为:

G(f,θ)=exp((log(f/f0))22(log(σf/f0))2)exp((θθ0)22σθ2)G(f,\theta)=\exp\Big(-\frac{(\log(f/f_0))^2}{2(\log(\sigma_f/f_0))^2} \Big)\exp\Big( -\frac{(\theta-\theta_0)^2}{2\sigma_\theta^2}\Big)

方向滤波与高斯函数无异,θ0\theta_0表示方向滤波的主方向(方向中心),σθ\sigma_\theta表示方向范围。见下图意会:

图a是Log-Gabor函数的频域图(Gabor函数中心不会有空洞),图b是方向滤波的频域图,图c是两者的合并,也就是二维Log-Gabor滤波器。

这是另一张参考图像。

3.5 补充问题

Gabor 函数与普通高斯函数在滤波中的区别?

Gabor 函数与普通高斯函数在滤波中的应用有一些区别,这主要体现在以下几个方面:

  1. 空间-频率分解

    • Gabor 滤波器:Gabor 滤波器既具有空间域特性,又具有频率域特性。它在空间域中类似于高斯函数,但同时还有频率选择性,能够对图像进行局部特征提取,并且在频域中有明确定义的中心频率和带宽。
    • 普通高斯函数:普通高斯函数主要用于平滑图像或信号,它只在空间域中定义,没有频率域特性。普通高斯函数的作用是模糊图像,去除噪声或平滑图像。
  2. 频域特性

    • Gabor 滤波器:Gabor 滤波器在频域中具有明确的中心频率和带宽,可以对不同频率的图像特征进行选择性地响应。
    • 普通高斯函数:普通高斯函数在频域中的响应相对均匀,不具有明确的频率选择性,对所有频率成分的影响相似。
  3. 方向选择性

    • Gabor 滤波器:Gabor 滤波器可以具有方向选择性,即可以设计为在特定方向上更敏感,这使得它在处理具有方向性结构的图像时表现更好。
    • 普通高斯函数:普通高斯函数在各个方向上的响应相同,没有方向选择性。
  4. 边缘响应

    • Gabor 滤波器:由于其正弦波的特性,Gabor 滤波器对图像中的边缘和纹理具有更好的响应,可以用于边缘检测和纹理分析。
    • 普通高斯函数:普通高斯函数主要用于平滑图像,对边缘和纹理的响应相对较弱。

综上所述,Gabor 函数与普通高斯函数相比,在滤波中具有更丰富的特性和应用场景,尤其适用于对图像的局部特征进行分析和提取。

4 相位一致性

相位一致性的基本原理及应用问题-CSDN博客

从相位一致性理解傅里叶变换中的相位 Phase Congruency - 知乎 (zhihu.com)

相位一致性的理解及两种新的相位一致性模型 - 百度文库 (baidu.com)

4.1 理解相位和幅度对于空域图像的关系

上面三节介绍了空域图像的频域变换。进行频域变换后有几个参数是需要明确的。对于二维傅里叶变换来说,分别是幅度、相位、频率、方向。对于Log Gabor变换来说,分别是幅度、相位、频率中心、方向中心、频率尺度、方向尺度。

那么这里面的幅度和相位和空域图像的关系是什么?

上图是对两幅毫不相干的图像进行傅里叶变换,分别获得其幅度图和相位图,将两者的相位和幅度进行交换,然后逆傅里叶变换获得两个新图。可以看出,相位携带着图像的信息,也就是说图像的轮廓和结构特征主要是由相位提供的。而幅度携带着图像的能量,也就是图像的能量分布(颜色值分布)。

如何理解”图像的轮廓和结构特征主要是由相位提供的“这句话?

Morrone和Owens在1987年的论文中认为,在图像中傅里叶分量最大相位的点处可以感知特征,傅里叶分量在方波的阶跃点以及三角波的波峰和波谷处都是同相的,并且该性质在尺度上是稳定的。