理想陷波滤波器
陷波滤波器可以阻止或允许以某一频率为中心的邻域内的频率通过,所以本质上还是带阻滤波器或带通滤波器,可以分别称为
带阻陷波滤波器和带通陷波滤波器。
理想的带阻陷波滤波器传递函数为:
带阻陷波滤波器通常成对工作,具有两个频率中心:
反演后,得到一个陷波带通滤波器。
Numpy陷波滤波:
import cv2 as cv
import numpy as np
import seaborn as sns
import matplotlib.pyplot as plt
img = cv.imread('./images/R-C.png')
img = cv.cvtColor(img, cv.COLOR_BGR2RGB)
img_gray = cv.cvtColor(img, cv.COLOR_RGB2GRAY)
测试用例:
理想的陷波带通、带阻滤波:
def notch_filter(img_gray, u0=0, v0=0, d0=50, ftype='pass'):
# 以频谱左上角为坐标原点
dft = cv.dft(img_gray.astype('float32'),flags = cv.DFT_COMPLEX_OUTPUT)
dft_shift = np.fft.fftshift(dft)
m, n, _ = dft_shift.shape
mask = np.zeros_like(dft_shift)
x_arr = np.concatenate([np.arange(m).reshape(m, 1)], axis=1)
y_arr = np.concatenate([np.arange(n).reshape(1, n)], axis=0)
dist1 = np.sqrt((x_arr - u0)**2 + (y_arr - v0)**2)
u1 = m - u0
v1 = n - v0
dist2 = np.sqrt((x_arr - u1)**2 + (y_arr - v1)**2)
mask[dist1 <= d0] = 1
mask[dist2 <= d0] = 1
if ftype != 'pass':
mask = 1 - mask
bpf_dft_shift = dft_shift * mask
magnitude_spectrum = cv.magnitude(bpf_dft_shift[:,:,0], bpf_dft_shift[:,:,1])
log_magnitude_spectrum = 20*np.log(magnitude_spectrum+1)
bpf_dft = np.fft.ifftshift(bpf_dft_shift)
img_ = cv.idft(bpf_dft)
img_bpf = cv.magnitude (img_[:,:,0],img_[:,:,1])
return img_bpf, log_magnitude_spectrum
img_1, log_1 = notch_filter(img_gray, u0=0, v0=0, d0=150, ftype='pass')
img_2, log_2 = notch_filter(img_gray, u0=0, v0=0, d0=150, ftype='stop')
img_3, log_3 = notch_filter(img_gray, u0=150, v0=250, d0=150, ftype='pass')
img_4, log_4 = notch_filter(img_gray, u0=150, v0=250, d0=150, ftype='stop')
巴特沃斯陷波带阻滤波器
def bw_notch_filter(img_gray, u0=0, v0=0, d0=50, N=1, ftype='pass'):
# 以频谱左上角为坐标原点
dft = cv.dft(img_gray.astype('float32'),flags = cv.DFT_COMPLEX_OUTPUT)
dft_shift = np.fft.fftshift(dft)
m, n, _ = dft_shift.shape
x_arr = np.concatenate([np.arange(m).reshape(m, 1)], axis=1)
y_arr = np.concatenate([np.arange(n).reshape(1, n)], axis=0)
dist1 = np.sqrt((x_arr - u0)**2 + (y_arr - v0)**2)
dist2 = np.sqrt((x_arr + u0)**2 + (y_arr + v0)**2)
mask = 1 / (1. + ((d0**2)/(np.multiply(dist1, dist2)+0.00001))**N)
if ftype == 'pass':
mask = 1 - mask
bpf_dft_shift = dft_shift * mask.reshape(m, n, 1)
magnitude_spectrum = cv.magnitude(bpf_dft_shift[:,:,0], bpf_dft_shift[:,:,1])
log_magnitude_spectrum = 20*np.log(magnitude_spectrum+1)
bpf_dft = np.fft.ifftshift(bpf_dft_shift)
img_ = cv.idft(bpf_dft)
img_bpf = cv.magnitude (img_[:,:,0],img_[:,:,1])
return img_bpf, log_magnitude_spectrum
测试效果:
高斯带阻陷波滤波器
def gaussian_notch_filter(img_gray, u0=0, v0=0, d0=50, ftype='pass'):
# 以频谱左上角为坐标原点
dft = cv.dft(img_gray.astype('float32'),flags = cv.DFT_COMPLEX_OUTPUT)
dft_shift = np.fft.fftshift(dft)
m, n, _ = dft_shift.shape
x_arr = np.concatenate([np.arange(m).reshape(m, 1)], axis=1)
y_arr = np.concatenate([np.arange(n).reshape(1, n)], axis=0)
dist1 = np.sqrt((x_arr - u0)**2 + (y_arr - v0)**2)
dist2 = np.sqrt((x_arr + u0)**2 + (y_arr + v0)**2)
mask = 1 - np.exp(np.multiply(dist1, dist2)/(d0**2)*(-0.5))
if ftype == 'pass':
mask = 1 - mask
bpf_dft_shift = dft_shift * mask.reshape(m, n, 1)
magnitude_spectrum = cv.magnitude(bpf_dft_shift[:,:,0], bpf_dft_shift[:,:,1])
log_magnitude_spectrum = 20*np.log(magnitude_spectrum+1)
bpf_dft = np.fft.ifftshift(bpf_dft_shift)
img_ = cv.idft(bpf_dft)
img_bpf = cv.magnitude (img_[:,:,0],img_[:,:,1])
return img_bpf, log_magnitude_spectrum
img_1, log_1 = gaussian_notch_filter(img_gray, u0=400, v0=400, d0=100, ftype='pass')
img_2, log_2 = gaussian_notch_filter(img_gray, u0=400, v0=400, d0=100, ftype='stop')
版权声明:本文为博主今晚打佬虎原创文章,版权归属原作者,如果侵权,请联系我们删除!
原文链接:https://blog.csdn.net/u014281392/article/details/123029131