06频率域图像增强

背景知识

  • 线性系统

X 1 ( t ) − > Y 1 ( t ) X 2 ( t ) − > Y 2 ( t ) X_1(t) -> Y_1(t)\\ X_2(t)->Y_2(t) X1(t)>Y1(t)X2(t)>Y2(t)

则:
X 1 ( t ) + X 2 ( t ) − > Y 1 ( t ) + Y 2 ( t ) a X 1 ( t ) − > a Y 1 ( t ) X_1(t)+X_2(t)->Y_1(t)+Y_2(t)\\ aX_1(t)->aY_1(t) X1(t)+X2(t)>Y1(t)+Y2(t)aX1(t)>aY1(t)

  • 平移不变性

对于某系统,有 x ( t ) − > y ( t ) x(t)->y(t) x(t)>y(t) 当输入信号平移 T T T 时, x ( t − T ) − > y ( t − T ) x(t-T)->y(t-T) x(tT)>y(tT) 则称该系统具有平移不变性。

  • 卷积

对于一个线性系统的输入f(x)和输出g(x),他们之间必然存在这样的关系:

g ( x ) = ∫ − ∞ + ∞ f ( τ ) h ( x − τ ) g(x)=\int_{-\infin}^{+\infin}f(\tau)h(x-\tau) g(x)=+f(τ)h(xτ)

其中 h ( x ) h(x) h(x)是这个系统的单位冲激响响应。该积分也称为卷积积分,记为: g ( x ) = f ( x ) ∗ h ( x ) g(x)=f(x)*h(x) g(x)=f(x)h(x)

离散一维卷积:
f ( x ) ∗ h ( x ) = 1 M ∑ m = 0 M − 1 f ( m ) h ( x − m ) f(x)*h(x)=\frac{1}{M}\sum_{m=0}^{M-1}f(m)h(x-m) f(x)h(x)=M1m=0M1f(m)h(xm)

二维卷积:
f ( x , y ) ∗ g ( x , y ) = ∫ − ∞ + ∞ f ( u , v ) h ( x − u , y − v ) d u d v f(x,y)*g(x,y) = \int_{-\infin}^{+\infin}f(u,v)h(x-u,y-v)dudv f(x,y)g(x,y)=+f(u,v)h(xu,yv)dudv

离散二维卷积:
f ( x , y ) ∗ g ( x , y ) = 1 M N ∑ m = 0 M − 1 ∑ n = 0 N − 1 f ( m , n ) h ( x − m , y − n ) f(x,y)*g(x,y) = \frac{1}{MN}\sum_{m=0}^{M-1}\sum_{n=0}^{N-1}f(m,n)h(x-m,y-n) f(x,y)g(x,y)=MN1m=0M1n=0N1f(m,n)h(xm,yn)

傅里叶变换

一维傅里叶变换:
F ( u ) = ∫ − ∞ + ∞ f ( x ) e − j 2 π u x d x F(u)=\int_{-\infin}^{+\infin}f(x)e^{-j2\pi ux}dx F(u)=+f(x)ej2πuxdx
逆变换:
f ( x ) = ∫ − ∞ + ∞ F ( u ) e j 2 π u x d u f(x)=\int_{-\infin}^{+\infin}F(u)e^{j2\pi ux}du f(x)=+F(u)ej2πuxdu

f(x)分解为不同频率的成分,F(u)代表频率为u的成分的强度

函数f(x)的傅立叶变换后一般是一个复量,它可以用下式表示:
F ( u ) = R ( u ) + j I ( u ) F ( u ) = ∣ F ( u ) ∣ e j ϕ ( u ) ∣ F ( u ) ∣ = R 2 ( u ) + I 2 ( u ) , f ( x ) 的 幅 度 谱 ϕ ( u ) = a r t a n [ I ( u ) R ( u ) ] , f ( x ) 的 相 位 谱 P ( u ) = ∣ F ( u ) ∣ 2 = P 2 ( u ) + I 2 ( u ) , f ( x ) 的 功 率 谱 F(u)=R(u)+jI(u)\\ F(u) = |F(u)|e^{j\phi(u)}\\ |F(u)|=\sqrt{R^2(u)+I^2(u)} , f(x)的幅度谱\\ \phi(u)=artan[\frac{I(u)}{R(u)}], f(x)的相位谱\\ P(u)=|F(u)|^2=P^2(u)+I^2(u),f(x)的功率谱 F(u)=R(u)+jI(u)F(u)=F(u)ejϕ(u)F(u)=R2(u)+I2(u) ,f(x)ϕ(u)=artan[R(u)I(u)],f(x)P(u)=F(u)2=P2(u)+I2(u),f(x)

二维连续傅里叶变换及其逆变换:
F ( u , v ) = ∫ − ∞ + ∞ ∫ − ∞ + ∞ f ( x , y ) e − j 2 π ( u x + v y ) d x d y F(u,v)=\int_{-\infin}^{+\infin}\int_{-\infin}^{+\infin}f(x,y)e^{-j2\pi (ux+vy)}dxdy F(u,v)=++f(x,y)ej2π(ux+vy)dxdy

f ( x , y ) = ∫ − ∞ + ∞ ∫ − ∞ + ∞ F ( u , v ) e j 2 π ( u x + v y ) d u d v f(x,y)=\int_{-\infin}^{+\infin}\int_{-\infin}^{+\infin}F(u,v)e^{j2\pi (ux+vy)}dudv f(x,y)=++F(u,v)ej2π(ux+vy)dudv

二维傅里叶变换的幅度谱和相位谱为:
∣ F ( u , v ) ∣ = R 2 ( u , v ) + I 2 ( u , v ) , 幅 度 谱 ϕ ( u , v ) = a r t a n [ I ( u , v ) R ( u , v ) ] , 相 位 谱 |F(u,v)|=\sqrt{R^2(u,v)+I^2(u,v)} , 幅度谱\\ \phi(u,v)=artan[\frac{I(u,v)}{R(u,v)}], 相位谱\\ F(u,v)=R2(u,v)+I2(u,v) ,ϕ(u,v)=artan[R(u,v)I(u,v)],

离散傅里叶变换 D F T DFT DFT (Discrete Fourier Transform)

一维离散傅立叶变换:
F ( u ) = ∑ k = 0 N − 1 f ( k ) e − j 2 π u k N , u = 0 , 1 , 2 , . . . , N − 1 f ( k ) = 1 N ∑ u = 0 N − 1 F ( u ) e j 2 π u k N , k = 0 , 1 , 2 , . . . , N − 1 F(u)=\sum_{k=0}^{N-1}f(k)e^{-\frac{j2\pi uk}{N}},u=0,1,2,...,N-1\\ f(k)=\frac{1}{N}\sum_{u=0}^{N-1}F(u)e^{\frac{j2\pi uk}{N}},k=0,1,2,...,N-1\\ F(u)=k=0N1f(k)eNj2πuk,u=0,1,2,...,N1f(k)=N1u=0N1F(u)eNj2πuk,k=0,1,2,...,N1
二维离散傅立叶变换:
F ( u , v ) = ∑ x = 0 M − 1 ∑ y = 0 N − 1 f ( x , y ) e − j 2 π ( u x M + v y N ) u = 0 , 1 , 2 , . . . , M − 1 v = 0 , 1 , 2 , . . . , N − 1 f ( x , y ) = 1 M N ∑ u = 0 M − 1 ∑ v = 0 N − 1 F ( u , v ) e j 2 π ( u x M + v y N ) x = 0 , 1 , 2 , . . . , M − 1 y = 0 , 1 , 2 , . . . , N − 1 F(u,v)=\sum_{x=0}^{M-1}\sum_{y=0}^{N-1}f(x,y)e^{-j2\pi(\frac{ ux}{M}+\frac{vy}{N})}\\ u=0,1,2,...,M-1\\ v=0,1,2,...,N-1\\ f(x,y)=\frac{1}{MN}\sum_{u=0}^{M-1}\sum_{v=0}^{N-1}F(u,v)e^{j2\pi(\frac{ ux}{M}+\frac{vy}{N})}\\ x=0,1,2,...,M-1\\ y=0,1,2,...,N-1\\ F(u,v)=x=0M1y=0N1f(x,y)ej2π(Mux+Nvy)u=0,1,2,...,M1v=0,1,2,...,N1f(x,y)=MN1u=0M1v=0N1F(u,v)ej2π(Mux+Nvy)x=0,1,2,...,M1y=0,1,2,...,N1
二维离散傅立叶变换特性:

  • 可分离性
  • 平移性
  • 线性
  • 旋转特性
  • 均值性
  • 微分特性

图像傅立叶变换的物理意义是得到图像的不同频率成分的分布情况。


频率域平滑滤波器

思想:边缘和其它尖锐变化(如噪声)在傅立叶变换中对于高频成分,因此平滑可通过衰减图像傅立叶变换的高频成分,保留低频成分来实现。

理想低通滤波器

一个二维的理想低通过滤器( I L P F ILPF ILPF)的系统函数满足(是一个分段函数):
H ( u , v ) = { 1 D ( u , v ) ≤ D 0 0 D ( u , v ) > D 0 D 0 为 截 止 频 率 , D ( u , v ) 为 距 离 函 数 D ( u , v ) = ( u 2 + v 2 ) 1 2 H(u,v) = \begin{cases} 1&D(u,v)\leq D_0\\ 0&D(u,v)>D_0\\ \end{cases}\\ D_0 为截止频率,D(u,v)为距离函数 D(u,v)=(u^2+v^2)^\frac{1}{2} H(u,v)={10D(u,v)D0D(u,v)>D0D0D(u,v)D(u,v)=(u2+v2)21
D 0 D_0 D0 半径内的频率分量无损通过,圆外的频率分量会被滤除。若滤除的高频分量中含有大量的边缘信息,会发生图像边缘模糊现象。钝化的图像被一种非常严重的振铃效果(理想低通滤波器的一种特性)所影响。

Code

import matplotlib.pyplot as plt
import numpy as np
import cv2 as cv

# 二维离散傅里叶变换
"""
def fft2(img):
    cp = np.copy(img).astype(np.complex)
    temp = np.copy(img)
    M,N = cp.shape
    def fun(u,v,f):
        value = 0
        for x in range(M):
            for y in range(N):
                exp = np.exp(np.complex(0,-2*np.pi*(u*x/M+v*y/N)))
                value += f[x][y]*exp
        return value
    for u in range(M):
        for v in range(N):
            cp[u][v] = fun(u,v,temp)
            print(cp[u][v])
    return cp
"""
# 理想低通滤波
def ilpf(f,D0):
    cp = np.copy(f)
    r,c = cp.shape
    H = np.zeros((r,c))
    for i in range(r):
        for j in range(c):
            if (i-r/2)**2+(j-c/2)**2<= D0**2:
                H[i][j] = 1

    return H*cp

def show(img,title,sub):
    plt.subplot(sub[0],sub[1],sub[2])
    plt.imshow(img,cmap="gray")
    plt.axis('off') # 去掉坐标轴
    plt.title(title)

if __name__=='__main__':
    img = cv.imread('10.png',0)

    # 傅里叶变换
    f = np.fft.fftn(img)

    # 坐标原点平移到中心
    fshift = np.fft.fftshift(f)

    # 避免频谱范围过大影响视觉效果
    fimg = np.log(1+np.abs(fshift))


    D = [5,15,30,80,230]
    i = 1
    for D0 in D:
        # 傅里叶逆变换
        ishift = np.fft.ifftshift(ilpf(fshift,D0))
        re = np.fft.ifftn(ishift)
        reimg = np.abs(re)
        show(reimg,'D0=%d'%D0,(2,3,i))
        i+=1
    show(img,'Original img',(2,3,6))
    plt.show()

在这里插入图片描述

B u t t e r w o r t h Butterworth Butterworth 低通滤波器

一个截止频率在与原点距离为 D 0 D0 D0 n n n B u t t e r w o r t h Butterworth Butterworth 低通过滤器( B L P F BLPF BLPF)的变换函数如下:

H ( u , v ) = 1 1 + [ D ( u , v ) D 0 ] 2 n H(u,v)=\frac{1}{1+[\frac{D(u,v)}{D_0}]^{2n}} H(u,v)=1+[D0D(u,v)]2n1
Code

if __name__=='__main__':    
    x = np.linspace(0,10,200)
    H = []
    for i in range(1,5):
        H.append(1/(1+np.power(x/5.,2*i)))
    plt.plot(x, H[0])
    plt.plot(x, H[1])
    plt.plot(x, H[2])
    plt.plot(x, H[3])
    plt.legend(['n=1','n=2','n=3','n=4'])
    plt.grid()
    plt.ylabel('H(u,v)')
    plt.xlabel('D(u,v)')
    plt.show()

在这里插入图片描述

# 巴特沃斯低通滤波器
def blpf(f,D0,n):
    """
    :param f: 傅里叶变换后的幅度谱
    :param D0: 截止距离
    :param n: 阶数
    :return: 新的幅度谱
    """
    cp = np.copy(f)
    r,c = cp.shape
    H = np.zeros((r, c))
    for i in range(r):
        for j in range(c):
            D = np.sqrt((i-r/2)**2+(j-c/2)**2)
            H[i][j] = 1/(1+np.power(D/D0,2*n))
    return H*cp

在这里插入图片描述

  • 在任何经二阶 B L P F BLPF BLPF 处理过的图像中都没有明显的振铃效果,这是过滤器在低频和高频之间的平滑过渡的结果。
  • 同样的截止频率下,比理想低通保留相对较多高频成分,所以模糊程度减少,但是去噪效果不如理想低通。
高斯低通滤波器

H ( u , v ) = e x p ( − D 2 ( u , v ) 2 D 0 2 ) H(u,v)=exp(-\frac{D^2(u,v)}{2D_0^2}) H(u,v)=exp(2D02D2(u,v))

[外链图片转存失败,源站可能有防盗链机制,建议将图片保存下来直接上传(img-3WI5B2RL-1584523520945)(image-20200318154429995.png)]

高斯低通滤波器 ( G L P F GLPF GLPF) 有更加平滑的过渡带,比 B L P F BLPF BLPF 相比衰减更快。无振铃,比 B L F P BLFP BLFP 去噪效果好。

Code

def glpf(f,D0):
    cp = np.copy(f)
    r,c = cp.shape
    H = np.zeros((r, c))
    for i in range(r):
        for j in range(c):
            D = (i-r/2)**2+(j-c/2)**2
            H[i][j] = np.exp(-D/(2*D0**2))

    return H*cp

在这里插入图片描述

频率域锐化滤波器

基本思想:低通滤波器的反操作,即 H h p ( u , v ) = 1 − H l p ( u , v ) H_{hp}(u,v)=1-H_{lp}(u,v) Hhp(u,v)=1Hlp(u,v)。与低通滤波器相对应,频率域内常用的高通滤波器有3种:

  1. 理想高通滤波器

在这里插入图片描述

D 0 = 15 D0=15 D0=15 D 0 = 30 D0=30 D0=30 时,振铃问题十分明显。

  1. 巴特沃斯高通滤波器

在这里插入图片描述

  1. 高斯高通滤波器

在这里插入图片描述

高通滤波器的问题:

  • 低频成分被严重地消弱了,图像失去主体成分,只剩下细节。
  • 改进措施
    • 加一个常数到变换函数 H ( u , v ) + A H(u,v) + A H(u,v)+A 这种方法被称为 高频提升
    • 为了解决变暗的趋势,在变换结果图像上再进行一次直方图均衡化。这种方法被称为 后过滤处理

高频提升滤波:
f h b = A f ( x , y ) − f l p ( x , y ) = ( A − 1 ) f ( x , y ) + f h p ( x , y ) H h p = ( A − 1 ) + H h p ( u , v ) f_{hb}=Af(x,y)-f_{lp}(x,y)\\ =(A-1)f(x,y)+f_{hp}(x,y)\\ H_{hp}=(A-1)+H_{hp}(u,v) fhb=Af(x,y)flp(x,y)=(A1)f(x,y)+fhp(x,y)Hhp=(A1)+Hhp(u,v)
高频加强滤波:在高通滤波器函数前乘以系数 b b b,再加上 a a a 倍的原始图像。
H h p e ( u , v ) = a + b H h p ( u , v ) , a ≥ 0 , b > a 当 b = 1 时 , 高 频 加 强 转 换 为 高 频 提 升 滤 波 H_{hpe}(u,v)=a+bH_{hp}(u,v), a\geq0,b>a\\ 当b=1时,高频加强转换为高频提升滤波 Hhpe(u,v)=a+bHhp(u,v),a0,b>ab=1


  • 改进措施
    • 加一个常数到变换函数 H ( u , v ) + A H(u,v) + A H(u,v)+A 这种方法被称为 高频提升
    • 为了解决变暗的趋势,在变换结果图像上再进行一次直方图均衡化。这种方法被称为 后过滤处理

高频提升滤波:
f h b = A f ( x , y ) − f l p ( x , y ) = ( A − 1 ) f ( x , y ) + f h p ( x , y ) H h p = ( A − 1 ) + H h p ( u , v ) f_{hb}=Af(x,y)-f_{lp}(x,y)\\ =(A-1)f(x,y)+f_{hp}(x,y)\\ H_{hp}=(A-1)+H_{hp}(u,v) fhb=Af(x,y)flp(x,y)=(A1)f(x,y)+fhp(x,y)Hhp=(A1)+Hhp(u,v)
高频加强滤波:在高通滤波器函数前乘以系数 b b b,再加上 a a a 倍的原始图像。
H h p e ( u , v ) = a + b H h p ( u , v ) , a ≥ 0 , b > a 当 b = 1 时 , 高 频 加 强 转 换 为 高 频 提 升 滤 波 H_{hpe}(u,v)=a+bH_{hp}(u,v), a\geq0,b>a\\ 当b=1时,高频加强转换为高频提升滤波 Hhpe(u,v)=a+bHhp(u,v),a0,b>ab=1


  • 0
    点赞
  • 0
    收藏
    觉得还不错? 一键收藏
  • 0
    评论

“相关推荐”对你有帮助么?

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值