用DFT矩阵法实现傅里叶变换不采用MATLAB的fft与ifft函数

一、序言

本文主要实现以下两个功能。

1.探究傅里叶变换FFT得到的心得，并且分别采用矩阵法以及DFT公式法实现MATLAB函数的FFT功能

2.对FFT变换功率倍数变化，以及幅值倍数变化问题，以及归一化做出完整的解释，并给出代码。

二、DFT与FFT变换的详细解释

其实这个网站上已经很多了，他的原理和公式我就不再叙述了，大家要搞清楚四个变换就好了，FT  连续时间傅里叶变换  

DTFT离散时间傅里叶变换

DFS傅里叶序列变换

DFT （FFT）

这是关于DSP信号处理的，大家可以自己去百度搜下他们的区别。总结就是下面两句话

1.信号在时域如果是连续的，频域就是非周期

2.时域离散信号 对应频域周期信号  

反过来也一样。

因此给大家一些结论：

1.傅里叶变换：是一种频域分析方法，可以让信号从时域转化为频域。

2.信号如果满足狄里克雷 (Dirichlet)条件，信号可以由无穷个基频以及高频频率正弦波或者余弦波的叠加。

3.为什么要变成离散信号呢，因为连续信号计算机没法处理，只能采样为离散信号，因此这就是数字信号处理这本书再讲的一个事情吧。

4.此外DFT与FFT是一样的，只不过FFT是一个算法方便计算出来DFT而已。

三、正文开始

针对开篇的两个问题：1.矩阵法如何实现2.FFT后功率是多少倍呢

1.矩阵法以及DFT公式法实现FFT

1.1矩阵法

首先matlab的函数fft(x,Nfft),Nfft为进行fft变换的点数，一般为2的n次方。

%% matlab DFT矩阵法实现

Nfft=64;
x=rand(1,Nfft);
%% FFT 特点 FFT后会对信号能量增加N倍
%IIFT特点 能量会降低N倍
%方法一、数据乘以dft矩阵矩阵为行向量
% y=dftmtx(Nfft)*x.'
%方法二、dft矩阵矩阵乘以矩阵x为列向量
y=x*dftmtx(Nfft);
y1=fft(x,Nfft);
norm(y-y1) %验证
%关于系数的问题 与能量有关。 % 
y2=y1/sqrt(Nfft);
sum(abs(x).^2)
sum(abs(y1).^2)

1.2.结论

函数dftmtx()是用来获得DFT矩阵的

可以通过两个方法：

1.y=x*dftmtx(Nfft) 获得傅里叶变换，此时信号是列向量

2.第二个是y=dftmtx(Nfft)*x.' 此时信号是行向量，

3.能量的问题：FFT 后能量会增加N倍，N是fft的点数，因此fft后整个序列应该除以sqrt(N)

有一点大家一定要注意：matlab对复数矩阵的转置是有点问题的，如果需要转置的向量x是个1x64的复数向量， 直接用x'语法 转置得到的是他的共轭向量的转置。例如x=[1+i 1-i 2+i] 直接x'

ans=[1-i; 1+i;2-i ] 发现虽然变成列向量了，但是虚数都取反了，正确的做法应该是用

conj(x')语法

2.IFFT矩阵法

%% IIFT特点 能量会降低N倍


x1=ifft(y1,Nfft); 
% 方法二 数据乘以矩阵转置
x2=real(y1*dftmtx(Nfft)'/Nfft);
% norm(x1-x2) %验证 MATLAB 与 矩阵法
% norm(x-x1) %验证 矩阵法 与原x数据
sum(abs(x).^2)
sum(abs(x2).^2)

1.2.结论

函数dftmtx()是用来获得DFT矩阵的

1.能量的问题：IFFT 后能量会降低N倍，N是fft的点数，因此fft后整个序列应该乘以sqrt(N)

3.公式法

这个是直接用DFT和IDFT的公式得到傅里叶变换序列的 for循环就好，

%% matlab DFT公式法实现

Nfft=64;
x=rand(1,Nfft);
%% FFT 特点 FFT后会对信号能量增加N倍
%IIFT特点 能量会降低N倍
% fft法
y1=fft(x,Nfft);
%公式法
DFT_x=0;
n=0:Nfft-1;
for k=0:Nfft-1
    temp=x(k+1)*exp(-1j*2*pi*n*k/Nfft);
DFT_x=DFT_x+temp;
end
y=DFT_x;
norm(y-y1) %验证
%关于系数的问题 与能量有关。 % 
y2=y1/sqrt(Nfft);
sum(abs(x).^2)
sum(abs(y1).^2)