MATLAB 中的二维离散傅立叶变换实现
2D Discrete Fourier Transform implementation in MATLAB
我正在使用矩阵乘法在 Matlab 中实现二维离散傅立叶变换。
我意识到这可以是一个可分离的操作,所以我为一维 DFT 创建了一个矩阵,并将其与输入图像的列相乘,然后与图像的行相乘。
但是,得到的2D DFT的值与使用MATLAB中内置函数计算的DFT(即fft2
)有很大差异。因此,在执行逆 DFT 重建图像时,结果图像没有正确重建(即它与原始图像不同,但如果我使用 fft2
函数,它是相同的)。
这是我写的代码。
%% signal is a matrix of MxN size
function res=myDFT(signal)
signal=double(signal);
l=size(signal,1);
x=[1:l];
%% u and x are matrices of MxM size
[x u]=meshgrid(x,x);
M1=l-1;
pre_dft=exp(1i*(-2*pi)./M1)/sqrt(M1);
pre_dft=(pre_dft.^(u.*x));
%the below matrix will be multiplied with the rows of the signal
post_dft=pre_dft;
% res is the resultant DFT of the signal matrix
% 1D DFT matrix is first multiplied with columns and then the
% rows of signal matrix
res=pre_dft*signal*post_dft;
end
如果有人能指出任何对编辑我的代码有用的东西或指出我的理论理解中的缺陷,我将不胜感激。
好的,您有一些我们需要修复的错误...如果您希望它的工作方式与 fft2
在 MATLAB 中的工作方式完全相同。
错误 #1 - 错误的归一化因子
fft
和最终 fft2
在计算变换时有 没有 归一化因子。在计算 DFT 矩阵的第一部分时,您可以去掉该除法语句。另外,DFT矩阵的第一部分,在指数中,需要除以信号的总长度,而不是长度减1。
错误 #2 - 完成 DFT 矩阵最后一部分的幂运算错误
你计算了一个从 1 到 l
的 meshgrid
,但是幂运算要求你 从 0
开始。因此,在进行幂计算之前,您需要将u
和x
减1。
错误 #3 - 将 FFT 应用于行和列
正如您已经注意到的,2D FFT 是可分离的,可以通过先对行执行 1D FFT,然后对列执行来完成。不幸的是,您正在执行的操作不正确。 post_dft
需要 转置 ,因为您要对列应用生成的中间操作。
因此,通过我上面提到的所有修复,您更正后的代码是:
function res=myDFT(signal)
signal=double(signal);
l=size(signal,1);
x=[1:l];
[x, u]=meshgrid(x,x);
%// Error #1
pre_dft=exp(1i*(-2*pi)./l); %// Change
%// Error #2
pre_dft=(pre_dft.^((u-1).*(x-1))); %// Change
%// Error #3
post_dft = pre_dft.'; %// Change
res = pre_dft*signal*post_dft;
end
用一些随机数据测试以上内容,并与 fft2
:
进行比较
rng(123123);
in = rand(7);
out1 = myDFT(in);
out2 = fft2(in);
out1
包含更正的自定义实现,而 out2
包含 MATLAB 的 fft2
算法的结果。我们得到:
>> out1
out1 =
Columns 1 through 4
26.2182 + 0.0000i -1.3805 + 1.0956i 2.2881 - 0.4435i -0.8005 + 1.5133i
-1.3067 + 0.3236i -0.5703 - 1.3884i -0.7127 + 0.1303i 3.2689 - 0.5995i
-0.6136 - 2.0731i -0.2776 + 1.2926i 0.1587 - 1.4504i -0.5305 + 1.4018i
-3.0703 + 0.3715i -1.5094 - 0.0216i 0.2573 - 3.2934i -1.1100 + 1.4180i
-3.0703 - 0.3715i -0.9940 + 0.1328i 1.7269 + 0.2915i 0.2814 + 2.2212i
-0.6136 + 2.0731i 1.8351 + 1.0629i 1.2891 + 1.7418i 0.4402 + 1.8756i
-1.3067 - 0.3236i -0.4974 - 0.9678i -2.2419 + 2.0839i -1.6844 + 0.9781i
Columns 5 through 7
-0.8005 - 1.5133i 2.2881 + 0.4435i -1.3805 - 1.0956i
-1.6844 - 0.9781i -2.2419 - 2.0839i -0.4974 + 0.9678i
0.4402 - 1.8756i 1.2891 - 1.7418i 1.8351 - 1.0629i
0.2814 - 2.2212i 1.7269 - 0.2915i -0.9940 - 0.1328i
-1.1100 - 1.4180i 0.2573 + 3.2934i -1.5094 + 0.0216i
-0.5305 - 1.4018i 0.1587 + 1.4504i -0.2776 - 1.2926i
3.2689 + 0.5995i -0.7127 - 0.1303i -0.5703 + 1.3884i
>> out2
out2 =
Columns 1 through 4
26.2182 + 0.0000i -1.3805 + 1.0956i 2.2881 - 0.4435i -0.8005 + 1.5133i
-1.3067 + 0.3236i -0.5703 - 1.3884i -0.7127 + 0.1303i 3.2689 - 0.5995i
-0.6136 - 2.0731i -0.2776 + 1.2926i 0.1587 - 1.4504i -0.5305 + 1.4018i
-3.0703 + 0.3715i -1.5094 - 0.0216i 0.2573 - 3.2934i -1.1100 + 1.4180i
-3.0703 - 0.3715i -0.9940 + 0.1328i 1.7269 + 0.2915i 0.2814 + 2.2212i
-0.6136 + 2.0731i 1.8351 + 1.0629i 1.2891 + 1.7418i 0.4402 + 1.8756i
-1.3067 - 0.3236i -0.4974 - 0.9678i -2.2419 + 2.0839i -1.6844 + 0.9781i
Columns 5 through 7
-0.8005 - 1.5133i 2.2881 + 0.4435i -1.3805 - 1.0956i
-1.6844 - 0.9781i -2.2419 - 2.0839i -0.4974 + 0.9678i
0.4402 - 1.8756i 1.2891 - 1.7418i 1.8351 - 1.0629i
0.2814 - 2.2212i 1.7269 - 0.2915i -0.9940 - 0.1328i
-1.1100 - 1.4180i 0.2573 + 3.2934i -1.5094 + 0.0216i
-0.5305 - 1.4018i 0.1587 + 1.4504i -0.2776 - 1.2926i
3.2689 + 0.5995i -0.7127 - 0.1303i -0.5703 + 1.3884i
我觉得不错!
我正在使用矩阵乘法在 Matlab 中实现二维离散傅立叶变换。
我意识到这可以是一个可分离的操作,所以我为一维 DFT 创建了一个矩阵,并将其与输入图像的列相乘,然后与图像的行相乘。
但是,得到的2D DFT的值与使用MATLAB中内置函数计算的DFT(即fft2
)有很大差异。因此,在执行逆 DFT 重建图像时,结果图像没有正确重建(即它与原始图像不同,但如果我使用 fft2
函数,它是相同的)。
这是我写的代码。
%% signal is a matrix of MxN size
function res=myDFT(signal)
signal=double(signal);
l=size(signal,1);
x=[1:l];
%% u and x are matrices of MxM size
[x u]=meshgrid(x,x);
M1=l-1;
pre_dft=exp(1i*(-2*pi)./M1)/sqrt(M1);
pre_dft=(pre_dft.^(u.*x));
%the below matrix will be multiplied with the rows of the signal
post_dft=pre_dft;
% res is the resultant DFT of the signal matrix
% 1D DFT matrix is first multiplied with columns and then the
% rows of signal matrix
res=pre_dft*signal*post_dft;
end
如果有人能指出任何对编辑我的代码有用的东西或指出我的理论理解中的缺陷,我将不胜感激。
好的,您有一些我们需要修复的错误...如果您希望它的工作方式与 fft2
在 MATLAB 中的工作方式完全相同。
错误 #1 - 错误的归一化因子
fft
和最终 fft2
在计算变换时有 没有 归一化因子。在计算 DFT 矩阵的第一部分时,您可以去掉该除法语句。另外,DFT矩阵的第一部分,在指数中,需要除以信号的总长度,而不是长度减1。
错误 #2 - 完成 DFT 矩阵最后一部分的幂运算错误
你计算了一个从 1 到 l
的 meshgrid
,但是幂运算要求你 从 0
开始。因此,在进行幂计算之前,您需要将u
和x
减1。
错误 #3 - 将 FFT 应用于行和列
正如您已经注意到的,2D FFT 是可分离的,可以通过先对行执行 1D FFT,然后对列执行来完成。不幸的是,您正在执行的操作不正确。 post_dft
需要 转置 ,因为您要对列应用生成的中间操作。
因此,通过我上面提到的所有修复,您更正后的代码是:
function res=myDFT(signal)
signal=double(signal);
l=size(signal,1);
x=[1:l];
[x, u]=meshgrid(x,x);
%// Error #1
pre_dft=exp(1i*(-2*pi)./l); %// Change
%// Error #2
pre_dft=(pre_dft.^((u-1).*(x-1))); %// Change
%// Error #3
post_dft = pre_dft.'; %// Change
res = pre_dft*signal*post_dft;
end
用一些随机数据测试以上内容,并与 fft2
:
rng(123123);
in = rand(7);
out1 = myDFT(in);
out2 = fft2(in);
out1
包含更正的自定义实现,而 out2
包含 MATLAB 的 fft2
算法的结果。我们得到:
>> out1
out1 =
Columns 1 through 4
26.2182 + 0.0000i -1.3805 + 1.0956i 2.2881 - 0.4435i -0.8005 + 1.5133i
-1.3067 + 0.3236i -0.5703 - 1.3884i -0.7127 + 0.1303i 3.2689 - 0.5995i
-0.6136 - 2.0731i -0.2776 + 1.2926i 0.1587 - 1.4504i -0.5305 + 1.4018i
-3.0703 + 0.3715i -1.5094 - 0.0216i 0.2573 - 3.2934i -1.1100 + 1.4180i
-3.0703 - 0.3715i -0.9940 + 0.1328i 1.7269 + 0.2915i 0.2814 + 2.2212i
-0.6136 + 2.0731i 1.8351 + 1.0629i 1.2891 + 1.7418i 0.4402 + 1.8756i
-1.3067 - 0.3236i -0.4974 - 0.9678i -2.2419 + 2.0839i -1.6844 + 0.9781i
Columns 5 through 7
-0.8005 - 1.5133i 2.2881 + 0.4435i -1.3805 - 1.0956i
-1.6844 - 0.9781i -2.2419 - 2.0839i -0.4974 + 0.9678i
0.4402 - 1.8756i 1.2891 - 1.7418i 1.8351 - 1.0629i
0.2814 - 2.2212i 1.7269 - 0.2915i -0.9940 - 0.1328i
-1.1100 - 1.4180i 0.2573 + 3.2934i -1.5094 + 0.0216i
-0.5305 - 1.4018i 0.1587 + 1.4504i -0.2776 - 1.2926i
3.2689 + 0.5995i -0.7127 - 0.1303i -0.5703 + 1.3884i
>> out2
out2 =
Columns 1 through 4
26.2182 + 0.0000i -1.3805 + 1.0956i 2.2881 - 0.4435i -0.8005 + 1.5133i
-1.3067 + 0.3236i -0.5703 - 1.3884i -0.7127 + 0.1303i 3.2689 - 0.5995i
-0.6136 - 2.0731i -0.2776 + 1.2926i 0.1587 - 1.4504i -0.5305 + 1.4018i
-3.0703 + 0.3715i -1.5094 - 0.0216i 0.2573 - 3.2934i -1.1100 + 1.4180i
-3.0703 - 0.3715i -0.9940 + 0.1328i 1.7269 + 0.2915i 0.2814 + 2.2212i
-0.6136 + 2.0731i 1.8351 + 1.0629i 1.2891 + 1.7418i 0.4402 + 1.8756i
-1.3067 - 0.3236i -0.4974 - 0.9678i -2.2419 + 2.0839i -1.6844 + 0.9781i
Columns 5 through 7
-0.8005 - 1.5133i 2.2881 + 0.4435i -1.3805 - 1.0956i
-1.6844 - 0.9781i -2.2419 - 2.0839i -0.4974 + 0.9678i
0.4402 - 1.8756i 1.2891 - 1.7418i 1.8351 - 1.0629i
0.2814 - 2.2212i 1.7269 - 0.2915i -0.9940 - 0.1328i
-1.1100 - 1.4180i 0.2573 + 3.2934i -1.5094 + 0.0216i
-0.5305 - 1.4018i 0.1587 + 1.4504i -0.2776 - 1.2926i
3.2689 + 0.5995i -0.7127 - 0.1303i -0.5703 + 1.3884i
我觉得不错!