智能图像处理：Python和OpenCV实现-课件 第六章 图像频域滤波

第6章图像频域滤波频域滤波

频率域是指从函数的频率角度出发分析函数，和频率域相对的是时间域。简单说就是如果从时间域分析信号时，时间是横坐标，振幅是纵坐标。而在频率域分析的时候则是频率是横坐标，振幅是纵坐标。在频率域对图像进行滤波处理的目的一方面是因为滤波在频率域更为直观，它可以解释空间域滤波的某些性质；另一方面也是因为在频率域中可以指定滤波器，进行滤波，然后在空间域使用在频域中滤波的结果在空间域进行逆滤波，得到在空间域中的图像，这时的图像就是已经滤掉噪声杂波的干净图像。

相比于图像的空间域处理，频域图像处理有许多优点，它可以通过频域成分的特殊性质完成一些空间域中无法完成的任务，而且频域图像处理更利于信号处理的解释，对于滤波过程产生的某些效果有着更直观的解释。目录CONTENTS傅里叶变换01低通滤波02高通滤波03带通和带阻滤波04同态滤波056.1傅里叶变换傅里叶变换傅里叶是18世纪法国的一位伟大的数学家，他最大的贡献在于指出任何周期函数都可以表示为不同频率的正弦或余弦之和的形式。无论函数有多复杂，只要它是周期性的，并且满足一定的数学条件，就一定可以用这样的正弦或余弦和的形式来表示。甚至在有些情况下，非周期函数也可以用正弦或余弦和的形式来表示。用傅里叶变换表示的函数特征可以完全通过傅里叶反变换来重建，而不会丢失任何信息。傅里叶变换的核心贡献在于如何求出每种正弦波和余弦波的频率，以及在给定每种正弦波和余弦波的比例系数时可以恢复出原始信号。6.1.1Numpy中的傅里叶变换图像二维傅里叶变换公式如6-1所示：

图像二维傅里叶逆变换公式如6-2所示：其中图像长M、宽N。f(x,y)表示时域图像，F(u,v)表示频域图像。x的范围为[0,M-1]，y的范围为[0,N-1]。

在Numpy库中有FFT（快速傅里叶变换）函数来实现傅立叶变换。二维的傅里叶变换函数语法格式为：

np.fft.fft2(src,n=None,axis=-1,norm=None)其中第一个参数是输入图像，即灰度图像。第二个参数是可选的，它决定输出数组的大小。如果它大于输入图像的尺寸，则在计算FFT之前用零填充输入图像。如果小于输入图像，将裁切输入图像。如果未传递任何参数，则输出数组的大小将与输入的大小相同。其余参数选择默认值即可。例6.1利用python的numpy库实现图像傅里叶变换及反变换。importcv2importnumpyasnpfrommatplotlibimportpyplotasplt#读取图像img=cv2.imread('d:\pics\lena.jpg',0)#傅里叶变换fft_img=np.fft.fft2(img)fft_shift=np.fft.fftshift(fft_img)fft_res=np.log(np.abs(fft_shift))#傅里叶逆变换ifft_shift=np.fft.ifftshift(fft_shift)ifft_img=np.fft.ifft2(ifft_shift)ifft_img=np.abs(ifft_img)#显示图像plt.subplot(131),plt.imshow(img,'gray')plt.title('OriginalImage'),plt.axis('off')plt.subplot(132),plt.imshow(fft_res,'gray')plt.title('FourierImage'),plt.axis('off')plt.subplot(133),plt.imshow(ifft_img,'gray')plt.title('InverseFourierImage'),plt.axis('off')plt.show()Numpy的FFT变换与逆变换图像6.1.2

OpenCV中的傅里叶变换OpenCV提供了离散傅里叶变换函数cv2.dft()和离散傅里叶逆变换函数cv2.idft()。傅里叶变换函数cv2.dft()语法格式为：dft=cv2.dft(np.float32(img),flags=cv2.DFT_COMPLEX_OUTPUT)其中第一个参数输入图像应转换为np.float32格式。第二个参数flags=cv2.DFT_COMPLEX_OUTPUT表示进行傅里叶变化的方法。输出参数dft有两个通道，一个通道是实部，另一个通道是虚部。傅里叶逆变换函数cv2.idft()的语法格式为：

iimg=cv2.idft(dft)求傅里叶逆变换后图像的幅值函数cv2.magnitude()语法格式为：res2=cv2.magnitude(iimg[:,:,0],iimg[:,:,1])表6-1flags标识符列表标识符名称意义DFT_INVERSE用一维或二维逆变换代替默认的正向变换。DFT_SCALE缩放比例标识符，输出的结果都会以I/N进行缩放，通常会结合DFT_INVERSE一起使用。DFT_ROWS对输入矩阵的每行进行正向或反向的变换，此标识符可以在处理多种矢量的时候用于减小资源开销，这些处理常常是三维或高维变换等复杂操作。DFT_COMPLEX_OUTPUT进行一维或二维实数数组正变换。这样的结果虽然是复数阵列，但拥有复数的共轭对称性，所以可以被写成一个拥有同样尺寸的实数阵列。DFT_REAL_OUTPUT进行一维或二维复数数组反变换。这样的结果通常是一个大小相同的|复矩阵。如果输入的矩阵有复数的共轭对称性(比如是一个带有DFT_COMPLEX_OUTPUT标识符的正变换结果)，便会输出实矩阵。例6.2利用OpenCV提供了函数cv2.dft()和函数cv2.idft实现傅里叶变换和逆变换。importcv2importnumpyasnpfrommatplotlibimportpyplotasplt#读取图像img=cv2.imread('d:\pics\lena.jpg',0)#傅里叶变换dft=cv2.dft(np.float32(img),flags=cv2.DFT_COMPLEX_OUTPUT)dftshift=np.fft.fftshift(dft)res1=20*np.log(cv2.magnitude(dftshift[:,:,0],dftshift[:,:,1]))#傅里叶逆变换ishift=np.fft.ifftshift(dftshift)iimg=cv2.idft(ishift)res2=cv2.magnitude(iimg[:,:,0],iimg[:,:,1])#显示图像plt.subplot(131),plt.imshow(img,'gray')plt.title('OriginalImage'),plt.axis('off')plt.subplot(132),plt.imshow(res1,'gray')plt.title('FourierImage'),plt.axis('off')plt.subplot(133),plt.imshow(res2,'gray')plt.title('InverseFourierImage'),plt.axis('off')plt.show()OpenCV中的FFT变换与逆变换图像6.2低通滤波低通滤波低通滤波是将频域图像中的高频部分滤除而通过低频部分。图像的边缘和噪声对应于频域图像中的高频部分，而低通滤波的作用即是减弱这部分的能量，从而达到图像平滑去噪的目的。对一幅图像使用低通滤波器，可以使图像比原始图像少尖锐的细节部分而突出平滑过渡部分，典型效果就是有效控制图像的模糊程度。常用的低通滤波器有理想低通滤波器、巴特沃斯低通滤波器和高斯低通滤波器。6.2.1理想低通滤波器最简单的低通滤波器是理想低通滤波器，基本思想是给定一个频率阈值，将高于该阈值的所有部分设置为0，而低于该频率的部分保持不变。理想低通滤波器的传递函数如下所示：

其中D0表示通带的半径，D(u,v)表示频域点(u,v)到频域图像原点的距离，称为截止频率。

例6.3对一幅图像进行不同范围频域的理想低通滤波，观察效果。importcv2importnumpyasnpfrommatplotlibimportpyplotasplt#第一步读入图像img=cv2.imread('d:\pics\lena.jpg',0)#第二步：进行数据类型转换img_float=np.float32(img)#第三步：使用cv2.dft进行傅里叶变换dft=cv2.dft(img_float,flags=cv2.DFT_COMPLEX_OUTPUT)#第四步：使用np.fft.fftshift将低频转移到图像中心dft_center=np.fft.fftshift(dft)#第五步：定义掩模：生成的掩模中间为1，周围为0#求得图像的中心点位置crow,ccol=int(img.shape[0]/2),int(img.shape[1]/2)#设置掩模区域为40*40的正方形mask=np.zeros((img.shape[0],img.shape[1],2),np.uint8)mask[crow-20:crow+20,ccol-20:ccol+20]=1#设置掩模区域为100*100的正方形mask1=np.zeros((img.shape[0],img.shape[1],2),np.uint8)mask1[crow-50:crow+50,ccol-50:ccol+50]=1#第六步：将掩模与傅里叶变化后图像相乘，保留中间部分mask_img=dft_center*maskmask1_img=dft_center*mask1#第七步：使用np.fft.ifftshift将低频移动到原来的位置img_idf=np.fft.ifftshift(mask_img)img1_idf=np.fft.ifftshift(mask1_img)#第八步：使用cv2.idft进行傅里叶的反变换img_idf=cv2.idft(img_idf)img1_idf=cv2.idft(img1_idf)#第九步：使用cv2.magnitude转换为空间域img_idf=cv2.magnitude(img_idf[:,:,0],img_idf[:,:,1])img1_idf=cv2.magnitude(img1_idf[:,:,0],img1_idf[:,:,1])#第十步：输出图像plt.subplot(131),plt.title('OriginImage')plt.imshow(img,cmap='gray'),plt.axis('off')plt.subplot(132),plt.title('Lowpassmask=20')plt.imshow(img_idf,cmap='gray'),plt.axis('off')plt.subplot(133),plt.title('Lowpassmask=50')plt.imshow(img1_idf,cmap='gray'),plt.axis('off')plt.show()程序运行结果如图所示。从结果看，图中左边为原图像，中间图像为在傅里叶变换频域上截取40ⅹ40正方形的低通区域，逆变换后的得到的图像，可以看出图像变得模糊、丢失了许多细节、有振铃现象，这也是理想低通滤波的特点。右边图像为在傅里叶变换频域上截取100ⅹ100正方形的低通区域，逆变换后的得到的图像，这时图像看起来与原图像基本差不多了，这也是对图像进行有损压缩的原理。不同掩模下的低通滤波图像6.2.2巴特沃斯低通滤波（Butterworthfilter）

例6.4对一幅图像进行巴特沃斯低通滤波，观察结果。importcv2importnumpyasnpdefcombine_images(images):#滤波后图像与频域图组合在一起

shapes=np.array([mat.shapeformatinimages])rows=np.max(shapes[:,0])copy_imgs=[cv2.copyMakeBorder(img,0,rows-img.shape[0],0,0,cv2.BORDER_CONSTANT,(0,0,0))forimginimages]returnnp.hstack(copy_imgs)deffft(img):#傅里叶变换

rows,cols=img.shape[:2]nrows=cv2.getOptimalDFTSize(rows)#得到傅里叶最优尺寸大小

ncols=cv2.getOptimalDFTSize(cols)nimg=np.zeros((nrows,ncols))nimg[:rows,:cols]=imgfft_mat=cv2.dft(np.float32(nimg),flags=cv2.DFT_COMPLEX_OUTPUT)returnnp.fft.fftshift(fft_mat)deffft_image(fft_mat):log_mat=cv2.log(1+cv2.magnitude(fft_mat[:,:,0],fft_mat[:,:,1]))cv2.normalize(log_mat,log_mat,0,255,cv2.NORM_MINMAX)returnnp.uint8(np.around(log_mat))defifft(fft_mat):#逆傅里叶变换

f_ishift_mat=np.fft.ifftshift(fft_mat)img_back=cv2.idft(f_ishift_mat)img_back=cv2.magnitude(*cv2.split(img_back))cv2.normalize(img_back,img_back,0,255,cv2.NORM_MINMAX)returnnp.uint8(np.around(img_back))deffft_distances(m,n):u=np.array([iifi<=m/2elsem-iforiinrange(m)],dtype=np.float32)v=np.array([iifi<=m/2elsem-iforiinrange(m)],dtype=np.float32)v.shape=n,1ret=np.sqrt(u*u+v*v)returnnp.fft.fftshift(ret)defBWfilter(rows,cols,d0,n):#巴特沃斯低通滤波

duv=fft_distances(*fft_mat.shape[:2])filter_mat=1/(1+np.power(duv/d0,2*n))filter_mat=cv2.merge((filter_mat,filter_mat))returnfilter_matdefdo_filter(_=None):d0=cv2.getTrackbarPos('D0',filter_win)n=cv2.getTrackbarPos('n',filter_win)filter_mat=BWfilter(fft_mat.shape[0],fft_mat.shape[1],d0,n)filtered_mat=filter_mat*fft_matimg_back=ifft(filtered_mat)cv2.imshow(image_win,combine_images([img_back,fft_image(filter_mat)]))if__name__=='__main__':img=cv2.imread('d:\pics\lena.jpg',0)rows,cols=img.shape[:2]filter_win='FilterParameters'image_win='ButterworthLowPassFilteredImage'dWindow(filter_win)dWindow(image_win)cv2.createTrackbar('D0',filter_win,20,min(rows,cols)//4,do_filter)cv2.createTrackbar('n',filter_win,1,5,do_filter)fft_mat=fft(img)do_filter()cv2.resizeWindow(filter_win,512,20)cv2.waitKey(0)cv2.destroyAllWindows()当D0=30、n=2时程序运行结果如图6-4所示。可以通过调整滤波器参数D0和n，看到滤波器参数D0和n的变换对图像滤波的影响效果。当滤波半径D0越大，滤波效果越好，图像越清晰；当滤波器阶数越大，图像变得越模糊。因为巴特沃斯低通滤波器在高、低频间的过渡平滑，因此没有出现明显的振铃效应。巴特沃斯低通滤波器参数D0、n调整巴特沃斯低通滤波后的结果巴特沃斯低通滤波器6.3高斯低通滤波（GaussianLowPassFilter）

例6.5对一幅图像进行高斯低通滤波，观察结果。importcv2importnumpyasnpdefcombine_images(images):shapes=np.array([mat.shapeformatinimages])rows=np.max(shapes[:,0])copy_imgs=[cv2.copyMakeBorder(img,0,rows-img.shape[0],0,0,cv2.BORDER_CONSTANT,(0,0,0))forimginimages]returnnp.hstack(copy_imgs)deffft(img):rows,cols=img.shape[:2]nrows=cv2.getOptimalDFTSize(rows)ncols=cv2.getOptimalDFTSize(cols)nimg=np.zeros((nrows,ncols))nimg[:rows,:cols]=imgfft_mat=cv2.dft(np.float32(nimg),flags=cv2.DFT_COMPLEX_OUTPUT)returnnp.fft.fftshift(fft_mat)deffft_image(fft_mat):log_mat=cv2.log(1+cv2.magnitude(fft_mat[:,:,0],fft_mat[:,:,1]))cv2.normalize(log_mat,log_mat,0,255,cv2.NORM_MINMAX)returnnp.uint8(np.around(log_mat))defifft(fft_mat):f_ishift_mat=np.fft.ifftshift(fft_mat)img_back=cv2.idft(f_ishift_mat)img_back=cv2.magnitude(*cv2.split(img_back))cv2.normalize(img_back,img_back,0,255,cv2.NORM_MINMAX)returnnp.uint8(np.around(img_back))deffft_distances(m,n):u=np.array([iifi<=m/2elsem-iforiinrange(m)],dtype=np.float32)v=np.array([iifi<=m/2elsem-iforiinrange(m)],dtype=np.float32)v.shape=n,1ret=np.sqrt(u*u+v*v)returnnp.fft.fftshift(ret)defglpfilter(rows,cols,d0,n):#高斯低通滤波

duv=fft_distances(*fft_mat.shape[:2])filter_mat=np.exp(-(duv*duv)/(2*d0*d0))filter_mat=cv2.merge((filter_mat,filter_mat))returnfilter_matdefdo_filter(_=None):d0=cv2.getTrackbarPos('D0',filter_win)n=cv2.getTrackbarPos('n',filter_win)filter_mat=glpfilter(fft_mat.shape[0],fft_mat.shape[1],d0,n)filtered_mat=filter_mat*fft_matimg_back=ifft(filtered_mat)cv2.imshow(image_win,combine_images([img_back,fft_image(filter_mat)]))if__name__=='__main__':img=cv2.imread('d:\pics\lena.jpg',0)rows,cols=img.shape[:2]filter_win='FilterParameters'image_win='GaussianLowPassFilteredImage'dWindow(filter_win)dWindow(image_win)cv2.createTrackbar('D0',filter_win,20,min(rows,cols)//4,do_filter)cv2.createTrackbar('n',filter_win,1,5,do_filter)fft_mat=fft(img)do_filter()cv2.resizeWindow(filter_win,512,20)cv2.waitKey(0)cv2.destroyAllWindows()当D0=30、n=3时程序运行结果如图6-5所示。可以通过调整滤波器参数D0、n，看到D0、n的变换对图像的影响效果。高斯低通滤波器参数D0、n调整高斯低通滤波后的结果高斯低通滤波器6.3高通滤波高通滤波高通滤波是使高频通过而使低频衰减的滤波，高通滤波器能削弱一幅图像的傅里叶变换里的低频分量，而增强高频分量（或保持它们不变）。被高通滤波的图像比原始图像少灰度级的平滑过渡而突出边缘等细节部分。6.3.1理想高通滤波只保留图像的高频部分，滤波过程与理想低通滤波一样，构造的掩模中间为0，边缘为1，然后与傅里叶变化后的图像结合，保留高频部分，去除低频部分如下所示：其中D0表示通带的半径，D(u,v)表示频域点(u,v)到频域图像原点的距离，称为截止频率。

例6.对一幅图像进行范围大小不同的频域高通滤波，使用掩模只保留高通部分，观察效果。importcv2importnumpyasnpfrommatplotlibimportpyplotasplt#第一步读入图片img=cv2.imread('d:\pics\lena.jpg',0)#第二步：进行数据类型转换img_float=np.float32(img)#第三步：使用cv2.dft进行傅里叶变换dft=cv2.dft(img_float,flags=cv2.DFT_COMPLEX_OUTPUT)#第四步：使用np.fft.fftshift将低频转移到图像中心dft_center=np.fft.fftshift(dft)#第五步：定义掩模：生成的掩模中间为0周围为1#求得图像的中心点位置crow,ccol=int(img.shape[0]/2),int(img.shape[1]/2)#设置掩模区域为10*10的正方形mask=np.ones((img.shape[0],img.shape[1],2),np.uint8)mask[crow-5:crow+5,ccol-5:ccol+5]=0

#设置掩模区域为60*60的正方形mask1=np.ones((img.shape[0],img.shape[1],2),np.uint8)mask1[crow-30:crow+30,ccol-30:ccol+30]=0#第六步：将掩模与傅里叶变化后图像相乘，保留中间部分mask_img=dft_center*maskmask1_img=dft_center*mask1#第七步：使用np.fft.ifftshift将低频移动到原来的位置img_idf=np.fft.ifftshift(mask_img)img1_idf=np.fft.ifftshift(mask1_img)#第八步：使用cv2.idft进行傅里叶的反变换img_idf=cv2.idft(img_idf)img1_idf=cv2.idft(img1_idf)#第九步：使用cv2.magnitude转换为空间域img_idf=cv2.magnitude(img_idf[:,:,0],img_idf[:,:,1])img1_idf=cv2.magnitude(img1_idf[:,:,0],img1_idf[:,:,1])#第十步：显示结果图像plt.subplot(131),plt.title('OriginImage')plt.imshow(img,cmap='gray'),plt.axis('off')plt.subplot(132),plt.title('Highpassmask=5')plt.imshow(img_idf,cmap='gray'),plt.axis('off')plt.subplot(133),plt.title('Highpassmask=30')plt.imshow(img1_idf,cmap='gray'),plt.axis('off')plt.show()程序运行结果如图所示。左图为原始图像，中图为10ⅹ10像素的频域高通滤波后的图像，右图为30ⅹ30像素的频域高通滤波后的图像。从图中可以看出高频部分对图像的重要性，减少了高频部分，图像的边缘变得模糊，图像不清晰。理想高通滤波后的图像6.3.2巴特沃斯高通滤波巴特沃斯高通滤波的函数表达式如下所示：其中n称为Butterworth高通滤波器的阶数。从函数表达式可知，巴特沃斯高通滤波没有理想高通滤波器那么剧烈。在巴特沃斯高通滤波中，阶数越高，滤波器过度越剧烈，振铃现象越明显。

例6.7对一幅图像进行巴特沃斯高通滤波，调节滤波半径和阶数，观察效果。importcv2importnumpyasnpdefcombine_images(images):shapes=np.array([mat.shapeformatinimages])rows=np.max(shapes[:,0])copy_imgs=[cv2.copyMakeBorder(img,0,rows-img.shape[0],0,0,cv2.BORDER_CONSTANT,(0,0,0))forimginimages]returnnp.hstack(copy_imgs)deffft(img):rows,cols=img.shape[:2]nrows=cv2.getOptimalDFTSize(rows)ncols=cv2.getOptimalDFTSize(cols)nimg=np.zeros((nrows,ncols))nimg[:rows,:cols]=imgfft_mat=cv2.dft(np.float32(nimg),flags=cv2.DFT_COMPLEX_OUTPUT)returnnp.fft.fftshift(fft_mat)deffft_image(fft_mat):log_mat=cv2.log(1+cv2.magnitude(fft_mat[:,:,0],fft_mat[:,:,1]))cv2.normalize(log_mat,log_mat,0,255,cv2.NORM_MINMAX)returnnp.uint8(np.around(log_mat))defifft(fft_mat):f_ishift_mat=np.fft.ifftshift(fft_mat)img_back=cv2.idft(f_ishift_mat)img_back=cv2.magnitude(*cv2.split(img_back))cv2.normalize(img_back,img_back,0,255,cv2.NORM_MINMAX)returnnp.uint8(np.around(img_back))deffft_distances(m,n):u=np.array([iifi<=m/2elsem-iforiinrange(m)],dtype=np.float32)v=np.array([iifi<=m/2elsem-iforiinrange(m)],dtype=np.float32)v.shape=n,1ret=np.sqrt(u*u+v*v)returnnp.fft.fftshift(ret)defBhpfilter(rows,cols,d0,n):duv=fft_distances(rows,cols)duv[rows//2,cols//2]=0.000001filter_mat=1/(1+np.power(d0/duv,2*n))filter_mat=cv2.merge((filter_mat,filter_mat))returnfilter_matdefdo_filter(_=None):d0=cv2.getTrackbarPos('D0',filter_win)n=cv2.getTrackbarPos('n',filter_win)filter_mat=Bhpfilter(fft_mat.shape[0],fft_mat.shape[1],d0,n)filtered_mat=filter_mat*fft_matimg_back=ifft(filtered_mat)cv2.imshow(image_win,combine_images([img_back,fft_image(filter_mat)]))if__name__=='__main__':img=cv2.imread('d:\pics\lena.jpg',0)rows,cols=img.shape[:2]filter_win='FilterParameters'image_win='ButterworthHighPassFilteredImage'dWindow(filter_win)dWindow(image_win)cv2.createTrackbar('D0',filter_win,20,min(rows,cols)//4,do_filter)cv2.createTrackbar('n',filter_win,1,5,do_filter)fft_mat=fft(img)do_filter()cv2.resizeWindow(filter_win,512,20)cv2.waitKey(0)cv2.destroyAllWindows()当D0=20、n=2时程序运行结果如图所示。巴特沃斯高通滤波器参数D0、n调整巴特沃斯高通滤波后的图像巴特沃斯高通滤波器6.3.3高斯高通滤波高斯高通滤波的函数表达式如下所示：其中D0表示通带的半径，D(u,v)表示频域点(u,v)到频域图像原点的距离。由于高斯高通滤波器的过渡是非常平坦的，不会产生振铃现象。

例6.8对一幅图像进行高斯高通滤波，调节D0和n，观察效果。importcv2importnumpyasnpdefcombine_images(images):shapes=np.array([mat.shapeformatinimages])rows=np.max(shapes[:,0])copy_imgs=[cv2.copyMakeBorder(img,0,rows-img.shape[0],0,0,cv2.BORDER_CONSTANT,(0,0,0))forimginimages]returnnp.hstack(copy_imgs)deffft(img):rows,cols=img.shape[:2]nrows=cv2.getOptimalDFTSize(rows)ncols=cv2.getOptimalDFTSize(cols)nimg=np.zeros((nrows,ncols))nimg[:rows,:cols]=imgfft_mat=cv2.dft(np.float32(nimg),flags=cv2.DFT_COMPLEX_OUTPUT)returnnp.fft.fftshift(fft_mat)deffft_image(fft_mat):log_mat=cv2.log(1+cv2.magnitude(fft_mat[:,:,0],fft_mat[:,:,1]))cv2.normalize(log_mat,log_mat,0,255,cv2.NORM_MINMAX)returnnp.uint8(np.around(log_mat))defifft(fft_mat):f_ishift_mat=np.fft.ifftshift(fft_mat)img_back=cv2.idft(f_ishift_mat)img_back=cv2.magnitude(*cv2.split(img_back))cv2.normalize(img_back,img_back,0,255,cv2.NORM_MINMAX)returnnp.uint8(np.around(img_back))deffft_distances(m,n):u=np.array([iifi<=m/2elsem-iforiinrange(m)],dtype=np.float32)v=np.array([iifi<=m/2elsem-iforiinrange(m)],dtype=np.float32)v.shape=n,1ret=np.sqrt(u*u+v*v)returnnp.fft.fftshift(ret)defghpfilter(rows,cols,d0,n):duv=fft_distances(*fft_mat.shape[:2])filter_mat=1-np.exp(-(duv*duv)/(2*d0*d0))filter_mat=cv2.merge((filter_mat,filter_mat))returnfilter_matdefdo_filter(_=None):d0=cv2.getTrackbarPos('D0',filter_win)n=cv2.getTrackbarPos('n',filter_win)filter_mat=ghpfilter(fft_mat.shape[0],fft_mat.shape[1],d0,n)filtered_mat=filter_mat*fft_matimg_back=ifft(filtered_mat)cv2.imshow(image_win,combine_images([img_back,fft_image(filter_mat)]))if__name__=='__main__':img=cv2.imread('d:\pics\lena.jpg',0)rows,cols=img.shape[:2]filter_win='FilterParameters'image_win='GaussianHighPassFilteredImage'dWindow(filter_win)dWindow(image_win)cv2.createTrackbar('D0',filter_win,20,min(rows,cols)//4,do_filter)cv2.createTrackbar('n',filter_win,1,5,do_filter)fft_mat=fft(img)do_filter()cv2.resizeWindow(filter_win,512,20)cv2.waitKey(0)cv2.destroyAllWindows()当D0=10、n=2时程序运行结果如图所示。高斯低通滤波器参数D0、n调整高斯高通滤波后图像高斯高通滤波器6.4带通和带阻滤波带通和带阻滤波在频域内进行滤波时，可能遇在某一范围的频域受到了杂波的干扰，需要对其滤除掉，这时就用到了带通或带阻滤波。6.4.1带通滤波带通滤波器是只保留某一范围区域的频率带，频率范围外的信息过滤掉，选择图像的部分信息。常用的带通滤波器包括理想带通滤波器、巴特沃斯带通滤波器、高斯带通滤波器。假设BW代表带宽，D0代表带宽的径向中心，则三种带通滤波器可以表示为：理想带通滤波器：

例6.9对一幅图像分别进行理想带通滤波器、巴特沃斯带通滤波器、高斯带通滤波器频域的带通滤波，调整参数，观察效果。importcv2importnumpyasnpdefcreateBPFilter(shape,center,bandCenter,bandWidth,lpType,n=2):rows,cols=shape[:2]r,c=np.mgrid[0:rows:1,0:cols:1]c-=center[0]r-=center[1]d=np.sqrt(np.power(c,2.0)+np.power(r,2.0))lpFilter_matrix=np.zeros(shape,np.float32)iflpType==0:#理想带通滤波器

lpFilter=np.copy(d)lpFilter[:,:]=1lpFilter[d>(bandCenter+bandWidth/2)]=0lpFilter[d<(bandCenter-bandWidth/2)]=0eliflpType==1:#巴特沃斯带通滤波器

lpFilter=1.0-1.0/(1+np.power(d*bandWidth/(d-pow(bandCenter,2)),2*n))eliflpType==2:#高斯带通滤波器lpFilter=np.exp(-pow((d-pow(bandCenter,2))/(d*bandWidth),2))lpFilter_matrix[:,:,0]=lpFilterlpFilter_matrix[:,:,1]=lpFilterreturnlpFilter_matrixdefstdFftImage(img_gray,rows,cols):fimg=np.copy(img_gray)fimg=fimg.astype(np.float32)#1.图像矩阵乘以（-1）^(r+c),中心化

forrinrange(rows):forcinrange(cols):if(r+c)%2:fimg[r][c]=-1*img_gray[r][c]img_fft=fftImage(fimg,rows,cols)returnimg_fftdeffftImage(img_gray,rows,cols):#离散傅里叶变换

rPadded=cv2.getOptimalDFTSize(rows)cPadded=cv2.getOptimalDFTSize(cols)imgPadded=np.zeros((rPadded,cPadded),dtype=np.float32)imgPadded[:rows,:cols]=img_grayimg_fft=cv2.dft(imgPadded,flags=cv2.DFT_COMPLEX_OUTPUT)returnimg_fftdefgraySpectrum(fft_img):#幅度、频谱

real=np.power(fft_img[:,:,0],2.0)imaginary=np.power(fft_img[:,:,1],2.0)amplitude=np.sqrt(real+imaginary)spectrum=np.log(amplitude+1.0)spectrum=cv2.normalize(spectrum,0,1,norm_type=cv2.NORM_MINMAX,\dtype=cv2.CV_32F)spectrum*=255returnamplitude,spectrumdefnothing(args):passif__name__=="__main__":img_file=r"d:\pics\lena.jpg"img_gray=cv2.imread(img_file,0)#1.快速傅里叶变换

rows,cols=img_gray.shape[:2]img_fft=stdFftImage(img_gray,rows,cols)amplitude,_=graySpectrum(img_fft)minValue,maxValue,minLoc,maxLoc=cv2.minMaxLoc(amplitude)#中心化后频谱的最大值在图像中心位置处

dWindow("tracks")max_radius=np.sqrt(pow(rows,2)+pow(cols,2))cv2.createTrackbar("BandCenter","tracks",0,int(max_radius),nothing)cv2.createTrackbar("BandWidth","tracks",0,int(max_radius),nothing)cv2.createTrackbar("Filtertype","tracks",0,2,nothing)whileTrue:#2.构建带通滤波器

bandCenter=cv2.getTrackbarPos("BandCenter","tracks")bandWidth=cv2.getTrackbarPos("BandWidth","tracks")lpType=cv2.getTrackbarPos("Filtertype","tracks")nrows,ncols=img_fft.shape[:2]ilpFilter=createBPFilter(img_fft.shape,maxLoc,bandCenter,bandWidth,lpType)#3.带通滤波器滤波

img_filter=ilpFilter*img_fft_,gray_spectrum=graySpectrum(img_filter)#观察滤波器的变化

#4.傅里叶反变换，取实部进行裁剪,并去除中心化

img_ift=cv2.dft(img_filter,flags=cv2.DFT_INVERSE+cv2.DFT_REAL_OUTPUT+cv2.DFT_SCALE)ori_img=np.copy(img_ift[:rows,:cols])forrinrange(rows):forcinrange(cols):if(r+c)%2:ori_img[r][c]=-1*ori_img[r][c]#截断高低值

ifori_img[r][c]<0:ori_img[r][c]=0ifori_img[r][c]>255:ori_img[r][c]=255ori_img=ori_img.astype(np.uint8)#5.带通滤波器的输出

cv2.imshow("Origin_img",img_gray)cv2.imshow("Filter_spectrum",gray_spectrum)cv2.imshow("BWF_image",ori_img)key=cv2.waitKey(1)ifkey==27:breakcv2.destroyAllWindows()理想滤波参数例6.9程序运行结果如图所示。理想带通滤波器理想带通滤波后图像巴特沃斯滤波参数巴特沃斯带通滤波器巴特沃斯带通滤后图像高斯滤波参数高斯带通滤波器高斯带通滤波后图像6.4.2带阻滤波

例6.10对一幅图像分别进行理想带阻滤波器、巴特沃斯带阻滤波器、高斯带阻滤波器频域的带阻滤波，调整参数，观察效果。importcv2importnumpyasnpdefcreateBRFilter(shape,center,bandCenter,bandWidth,lpType,n=2):rows,cols=shape[:2]r,c=np.mgrid[0:rows:1,0:cols:1]c-=center[0]r-=center[1]d=np.sqrt(np.power(c,2.0)+np.power(r,2.0))lpFilter_matrix=np.zeros(shape,np.float32)iflpType==0:#理想带阻滤波器

lpFilter=np.copy(d)lpFilter[:,:]=0lpFilter[d>(bandCenter+bandWidth/2)]=1lpFilter[d<(bandCenter-bandWidth/2)]=1eliflpType==1:#巴特沃斯带阻滤波器

lpFilter=1.0/(1+np.power(d*bandWidth/(d-pow(bandCenter,2)),2*n))eliflpType==2:#高斯带阻滤波器lpFilter=1-np.exp(-pow((d-pow(bandCenter,2))/(d*bandWidth),2))lpFilter_matrix[:,:,0]=lpFilterlpFilter_matrix[:,:,1]=lpFilterreturnlpFilter_matrixdefstdFftImage(img_gray,rows,cols):fimg=np.copy(img_gray)fimg=fimg.astype(np.float32)#图像矩阵乘以（-1）^(r+c),中心化

forrinrange(rows):forcinrange(cols):if(r+c)%2:fimg[r][c]=-1*img_gray[r][c]img_fft=fftImage(fimg,rows,cols)returnimg_fftdeffftImage(img_gray,rows,cols):rPadded=cv2.getOptimalDFTSize(rows)cPadded=cv2.getOptimalDFTSize(cols)imgPadded=np.zeros((rPadded,cPadded),dtype=np.float32)imgPadded[:rows,:cols]=img_grayimg_fft=cv2.dft(imgPadded,flags=cv2.DFT_COMPLEX_OUTPUT)returnimg_fftdefgraySpectrum(fft_img):real=np.power(fft_img[:,:,0],2.0)imaginary=np.power(fft_img[:,:,1],2.0)amplitude=np.sqrt(real+imaginary)spectrum=np.log(amplitude+1.0)spectrum=cv2.normalize(spectrum,0,1,norm_type=cv2.NORM_MINMAX,dtype=cv2.CV_32F)spectrum*=255returnamplitude,spectrumdefnothing(args):passif__name__=="__main__":img_file=r"d:\pics\lena.jpg"img_gray=cv2.imread(img_file,0)#1.快速傅里叶变换

rows,cols=img_gray.shape[:2]img_fft=stdFftImage(img_gray,rows,cols)amplitude,_=graySpectrum(img_fft)minValue,maxValue,minLoc,maxLoc=cv2.minMaxLoc(amplitude)#中心化后频谱的最大值在图片中心位置处

dWindow("tracks")max_radius=np.sqrt(pow(rows,2)+pow(cols,2))cv2.createTrackbar("BandCenter","tracks",0,int(max_radius),nothing)cv2.createTrackbar("BandWidth","tracks",0,int(max_radius),nothing)cv2.createTrackbar("Filtertype","tracks",0,2,nothing)whileTrue:#2.构建带阻滤波器

bandCenter=cv2.getTrackbarPos("BandCenter","tracks")bandWidth=cv2.getTrackbarPos("BandWidth","tracks")lpType=cv2.getTrackbarPos("Filtertype","tracks")nrows,ncols=img_fft.shape[:2]ilpFilter=createBRFilter(img_fft.shape,maxLoc,bandCenter,bandWidth,lpType)#3.带阻滤波器滤波

img_filter=ilpFilter*img_fft_,gray_spectrum=graySpectrum(img_filter)#观察滤波器的变化

#4.傅里叶反变换，并取实部进行裁剪,并去中心化

img_ift=cv2.dft(img_filter,flags=cv2.DFT_INVERSE+\cv2.DFT_REAL_OUTPUT+cv2.DFT_SCALE)ori_img=np.copy(img_ift[:rows,:cols])ori_img[ori_img<0]=0ori_img[