DFT实现与python循环

Posted on

二维数组（图像）DFT变换

  • f(x,y)：图像二维数组（时域） ； F(u,v) 图像DFT变换生成的二维数组（频域）
  • 维度数N已知，时域里图像矩阵x(i,j)（对应f）已知，目标输出变换后的图像频域矩阵Z(m,n)（对应F）

code: <pre><code> dft.py（对应第一个公式）
from future import division from future import print_function

import numpy as np

def dft2(x): '''

:param x: 待分析数据,二数数组。
:return Z: 原始数据的DFT结果。
'''
N = 256
Z = np.zeros((N,N))
for m in range(N):
    for n in range(N):
        for i in range(N):
            for j in range(N):   #四层循环，先从最里层考虑，最开始先访问谁
                Z[m,n] = x[i,j]*np.exp((-2)*np.pi*1J*(m*i+n*j)/N)/N
return Z

</code></pre> main.py

<pre><code> from future import division from future import print_function

import numpy as np import matplotlib.pyplot as plt from dft import dft2

if name == 'main': N = 256 # 图片尺寸 x = np.arange(N, dtype=np.float32) y = np.arange(N, dtype=np.float32)

R = 10
theta = np.pi/6                                  # 角度可以自己修改

u = R * np.cos(theta)
v = R * np.sin(theta)
Z = np.zeros((x.size, y.size))
for i_x in range(x.size):
    for i_y in range(y.size):
        Z[i_x, i_y] = np.cos(np.pi / N * (u * x[i_x] + v * y[i_y]))

Z_dft = dft2(Z)                                  # 调用DFT

Z_dft_shift = abs(np.fft.fftshift(Z_dft))        # 将图像中的低频部分移动到图像的中心

plt.subplot(221)
plt.imshow(Z, cmap='gray')
plt.title('Original')
plt.axis('off')

plt.subplot(222)
plt.imshow(abs(Z_dft), cmap='gray')              # 展示DFT结果
plt.title('DFT2')
plt.axis('off')

plt.subplot(223)
plt.imshow(Z_dft_shift, cmap='gray')             # 展示移动低频到中心后的图像
plt.title('shift')
plt.axis('off')

plt.subplot(224)
plt.plot(Z_dft_shift[x.size//2, :])              # 取一行查看数值，位置可以自己定义
plt.title('x = {}'.format(x.size//2))

plt.show()

</code></pre>

TIPs:

  1. 循环层数 == 时间复杂度的幂数（可以相互做判据）， 如dft.py里，先访问二维矩阵x的i、j（做256✖️256次操作）,时间复杂度O(N^2)，再访问并赋值给二维矩阵Z的m、n（再做256✖️256次操作），时间复杂度变成O(N^4)，所以要写四层循环
  2. 循环顺序，先从最里层开始考虑哪个先迭代，先行后列
  3. 矩阵角度理解DFT 感觉不一定对，仅供思路提示
  4. numpy注意python和numpy里切片不一样，自带的list不能x[i,j]，只能x[i][j],numpy里的array可以
  5. 调试时，如果图像数据过大，可以resize图像维度debug（见eigenface的code）