Python中的DFT矩阵

Question

[最简单的方法是在python中获取2-D DFT的DFT matrix？我在numpy.fft中找不到此类功能。谢谢！

Answer 1

import numpy as np
def DFT_matrix(N):
    i, j = np.meshgrid(np.arange(N), np.arange(N))
    omega = np.exp( - 2 * pi * 1J / N )
    W = np.power( omega, i * j ) / sqrt(N)
    return W

EDIT对于2D FFT矩阵，可以使用以下代码：
x = np.zeros(N, N) # x is any input data with those dimensions
W = DFT_matrix(N)
dft_of_x = W.dot(x).dot(W)

Answer 2

import scipy as sp

def dftmtx(N):
    return sp.fft(sp.eye(N))

如果您知道更快的方式（可能会更复杂），我们将不胜感激。

只是使其与主要问题更相关-您也可以使用numpy来做到这一点：
import numpy as np

dftmtx = np.fft.fft(np.eye(N))

[当我对它们两个进行基准测试时，我的印象都比较快，但是我尚未完全完成，而且是在某个时间之前，所以请不要相信我。

这里是python中FFT实现的很好的参考资料：http://nbviewer.ipython.org/url/jakevdp.github.io/downloads/notebooks/UnderstandingTheFFT.ipynb而是从速度的角度来看，但是在这种情况下，我们实际上可以看到有时它也具有简单性。 

Answer 3

def DFT_matrix_2d(N):
    i, j = np.meshgrid(np.arange(N), np.arange(N))
    A=np.multiply.outer(i.flatten(), i.flatten())
    B=np.multiply.outer(j.flatten(), j.flatten())
    omega = np.exp(-2*np.pi*1J/N)
    W = np.power(omega, A+B)/N
    return W

Answer 4

dftmtx = lambda N: np.fft.fft(np.eye(N))

您可以使用dftmtx（N）进行调用。示例：

In [62]: dftmtx(2)
Out[62]: 
array([[ 1.+0.j,  1.+0.j],
       [ 1.+0.j, -1.+0.j]])

Answer 5

具有16点DFT矩阵的示例：

>>> import scipy.linalg
>>> import numpy as np
>>> m = scipy.linalg.dft(16)

验证单一属性，音符矩阵未缩放，因此16*np.eye(16)：

>>> np.allclose(np.abs(np.dot( m.conj().T, m )), 16*np.eye(16))
True

对于2D DFT矩阵，在我们处理矩阵代数的时候，这只是张量积或特别是Kronecker积的问题。

>>> m2 = np.kron(m, m) # 256x256 matrix, flattened from (16,16,16,16) tensor

现在我们可以为它提供平铺的可视化效果，它是通过将每行重新排列成一个正方形块来完成的

>>> import matplotlib.pyplot as plt
>>> m2tiled = m2.reshape((16,)*4).transpose(0,2,1,3).reshape((256,256))
>>> plt.subplot(121)
>>> plt.imshow(np.real(m2tiled), cmap='gray', interpolation='nearest')
>>> plt.subplot(122)
>>> plt.imshow(np.imag(m2tiled), cmap='gray', interpolation='nearest')
>>> plt.show()

结果（真实部分和影像部分分开：]

如您所见，它们是2D DFT基本函数
Link到文档

Answer 6

M = 16
N = 16
X = np.random.random((M,N)) + 1j*np.random.random((M,N))
Y = np.fft.fft2(X)
W = np.zeros((M*N,M*N),dtype=np.complex)
​
hold = []
for m in range(M):
    for n in range(N):
        hold.append((m,n))
​
for j in range(M*N):
    for i in range(M*N):
        k,l = hold[j]
        m,n = hold[i]
        W[j,i] = np.exp(-2*np.pi*1j*(m*k/M + n*l/N))

np.allclose(np.dot(W,X.ravel()),Y.ravel())
True

如果要将归一化更改为正交，则可以除以1 / sqrt（MN），或者如果要进行逆变换，只需更改指数中的符号即可。