数字图像处理-傅里叶变换系列编程

在这里插入图片描述

问题

问题1

通过计算一维傅里叶变换实现图像二维快速傅里叶变换
在这里插入图片描述
代码实现

def dft2D(f):h,w = f.shapeF = np.zeros(f.shape, dtype=complex)for i in range(h):F[i,:] = np.fft.fft(f[i,:])for i in range(w):F[:,i] = np.fft.fft(F[:,i])return F

问题2

图像二维快速傅里叶逆变换
在这里插入图片描述
代码实现

def idft2D(F):h,w = F.shapeF1 = np.conj(F)f = np.zeros(F1.shape,dtype=complex)for i in range(h):f[i,:] = np.fft.fft(F1[i,:])for i in range(w):f[:,i] = np.fft.fft(f[:,i])f = f/(h*w)f = np.conj(f)f = np.abs(f)return f

问题3

测试图像二维快速傅里叶变换与逆变换
代码

import cv2
import numpy as np
#归一化
def normalization(data):_range = np.max(data) - np.min(data)return (data - np.min(data)) / _rangedef dft2D(f):h,w = f.shapeF = np.zeros(f.shape, dtype=complex)for i in range(h):F[i,:] = np.fft.fft(f[i,:])for i in range(w):F[:,i] = np.fft.fft(F[:,i])return Fdef idft2D(F):h,w = F.shapeF1 = np.conj(F)f = np.zeros(F1.shape,dtype=complex)for i in range(h):f[i,:] = np.fft.fft(F1[i,:])for i in range(w):f[:,i] = np.fft.fft(f[:,i])f = f/(h*w)f = np.conj(f)f = np.abs(f)return fif __name__ == '__main__':f1 = cv2.imread("rose512.tif",cv2.IMREAD_GRAYSCALE)f = normalization(f1)F = dft2D(f)g = idft2D(F)d = np.abs(f-g)cv2.imshow('original', f1)cv2.imshow('g', g)cv2.imshow('d', d)cv2.waitKey(0)cv2.destroyAllWindows()

结果展示
在这里插入图片描述

问题4

计算图像的中心化二维快速傅里叶变换与谱图像

问题分析

首先合成矩形图像，测试结果如下图所示，其长为 60，宽为 10。
在这里插入图片描述

代码

import numpy as np
import cv2
import Q3if __name__ == '__main__':img = np.zeros((512,512))img[512//2-30:512//2+30,512//2-5:512//2+5]=1cv2.imshow('img', img)F1 = Q3.dft2D(img)F1 = Q3.normalization(np.abs(F1))img_c = np.zeros(img.shape)cv2.imshow('F1',  F1)for i in range (img.shape[0]):for j in range (img.shape[1]):img_c[i,j] = img[i,j]*((-1)**(i+j))F2 = Q3.dft2D(img_c)cv2.imshow('F2', Q3.normalization(np.abs(F2)))F3 = np.log(1+np.abs(F2))cv2.imshow('F3',  Q3.normalization(F3))cv2.waitKey(0)cv2.destroyAllWindows()

结果展示
在这里插入图片描述

选做题

计算其他 5 幅图像的二维快速傅里叶变换： house.tif, house02.tif, lena_gray_512.tif, lunar_surface.tif, characters_test_pattern.tif。注意，有些图像的尺寸不是 2 的整数次幂，需要进行相应的像素填补处理。如果图像有多个通道可以选择其中的一个通道进行计算。

问题分析
在这里插入图片描述
代码实现

import numpy as np
import cv2
import Q3def isPower(n):if n < 1:return Falsei = 1while i <= n:if i == n:return Truei <<= 1return Falsedef pad(img):if isPower(img.shape[0]) and isPower(img.shape[1]):return imgh = w = 0i = 1while True:if 2**i>=img.shape[0] and h==0:h = 2**iif 2**i>=img.shape[1] and w==0:w = 2**iif h!=0 and w!=0:breaki=i+1img2 = np.zeros((h,w))img2[(h-img.shape[0])//2:h-(h-img.shape[0])//2,(w-img.shape[1])//2:w-(w-img.shape[1])//2] = imgreturn img2def center(f):img_c = np.zeros(f.shape)for i in range(f.shape[0]):for j in range(f.shape[1]):img_c[i, j] = f[i, j] * ((-1) ** (i + j))return img_cif __name__ == '__main__':f1 = cv2.imread("house.tif", cv2.IMREAD_GRAYSCALE)f2 = cv2.imread("house02.tif", cv2.IMREAD_GRAYSCALE)f3 = cv2.imread("lena_gray_512.tif",cv2.IMREAD_GRAYSCALE)f4 = cv2.imread("lunar_surface.tif", cv2.IMREAD_GRAYSCALE)f5 = cv2.imread("Characters_test_pattern.tif", cv2.IMREAD_GRAYSCALE)f11 = pad(Q3.normalization(f1))F11 = Q3.dft2D(f11)F11 = Q3.normalization(np.abs(F11))f12 = center(f11)F12 = Q3.dft2D(f12)F13 = np.log(1 + np.abs(F12))cv2.imshow('orign_house', cv2.resize(Q3.normalization(f11), (f11.shape[0]//2,f11.shape[1]//2))	)cv2.imshow('F1_house',cv2.resize(Q3.normalization(F11) , (F11.shape[0]//2,F11.shape[1]//2)))cv2.imshow('F2_house', cv2.resize(Q3.normalization(np.abs(F12)), (F12.shape[0]//2,F12.shape[1]//2)))cv2.imshow('F3_house',cv2.resize( Q3.normalization(F13), (F13.shape[0]//2,F13.shape[1]//2)))f21 = pad(Q3.normalization(f2))F21 = Q3.dft2D(f21)F21 = Q3.normalization(np.abs(F21))f22 = center(f21)F22 = Q3.dft2D(f22)F23 = np.log(1 + np.abs(F22))cv2.imshow('orign_house02', cv2.resize(Q3.normalization(f21), (f21.shape[0]//2,f21.shape[1]//2))	)cv2.imshow('F1_house02',cv2.resize(Q3.normalization(F21) , (F21.shape[0]//2,F21.shape[1]//2)))cv2.imshow('F2_house02',cv2.resize(Q3.normalization(np.abs(F22)), (F22.shape[0]//2,F22.shape[1]//2)))cv2.imshow('F3_house02',cv2.resize( Q3.normalization(F23), (F23.shape[0]//2,F23.shape[1]//2)))f31 = pad(Q3.normalization(f3))F31 = Q3.dft2D(f31)F31 = Q3.normalization(np.abs(F31))f32 = center(f31)F32 = Q3.dft2D(f32)F33 = np.log(1 + np.abs(F32))cv2.imshow('orign_lena_gray_512',cv2.resize(Q3.normalization(f31), (f31.shape[0]//2,f31.shape[1]//2))	)cv2.imshow('F1_lena_gray_512',  cv2.resize(Q3.normalization(F31) , (F31.shape[0]//2,F31.shape[1]//2)))cv2.imshow('F2_lena_gray_512',  cv2.resize(Q3.normalization(np.abs(F32)), (F32.shape[0]//2,F32.shape[1]//2)))cv2.imshow('F3_lena_gray_512', cv2.resize( Q3.normalization(F33), (F33.shape[0]//2,F33.shape[1]//2)))f41 = pad(Q3.normalization(f4))F41 = Q3.dft2D(f41)F41 = Q3.normalization(np.abs(F41))f42 = center(f41)F42 = Q3.dft2D(f42)F43 = np.log(1 + np.abs(F42))cv2.imshow('orign_lunar_surface', cv2.resize(Q3.normalization(f41), (f41.shape[1]//2,f41.shape[0]//2))	)cv2.imshow('F1_lunar_surface', cv2.resize(Q3.normalization(F41) , (F41.shape[1]//2,F41.shape[0]//2)))cv2.imshow('F2_lunar_surface', cv2.resize(Q3.normalization(np.abs(F42)), (F42.shape[1]//2,F42.shape[0]//2)))cv2.imshow('F3_lunar_surface',cv2.resize( Q3.normalization(F43), (F43.shape[1]//2,F43.shape[0]//2)))f51 = pad(Q3.normalization(f5))F51 = Q3.dft2D(f51)F51 = Q3.normalization(np.abs(F51))f52 = center(f51)F52 = Q3.dft2D(f52)F53 = np.log(1 + np.abs(F52))cv2.imshow('orign_Characters_test_pattern', cv2.resize(Q3.normalization(f51), (f51.shape[1]//2,f51.shape[0]//2))	)cv2.imshow('F1_Characters_test_pattern', cv2.resize(Q3.normalization(F51) , (F51.shape[1]//2,F51.shape[0]//2)))cv2.imshow('F2_Characters_test_pattern', cv2.resize(Q3.normalization(np.abs(F52)), (F52.shape[1]//2,F52.shape[0]//2)))cv2.imshow('F3_Characters_test_pattern',cv2.resize( Q3.normalization(F53), (F53.shape[1]//2,F53.shape[0]//2)))cv2.waitKey(0)cv2.destroyAllWindows()