zhangfft.m

给出了用MATLAB实现FFT算法

function y=fuyufft(x)

m=nextpow2(x);  % NEXTPOW2(N) returns the first P such that 2^P >= abs(N). 
N=2^m;
if length(x)<N
    x=[x,zeros(1,N-length(x))];% if datanumber of x is larger than 2^m ,then make up the blank space by zeros
end
a=[0:1:N-1];
a=dec2bin(a);  % DEC2BIN Convert decimal integer to a binary string.
               % DEC2BIN(D) returns the binary representation of D as a string.
a=fliplr(a);   %Flip matrix in left/right direction.   eg: X = 1 2 3     becomes  3 2 1
a=bin2dec(a)+1;
y=x(a);        %make the right sequence
tep=y;

for a=0:m-1    %   m可以认为是m级堞形运算   a表示现在在做第几级堞形运算
        c=2^a; %   c可以认为是堞形运算两元素间隔的元素数，比如说N＝8，则c=1,2,4
        d=a+1;
        
        for b=1:2*c:N   
            for e=0:1:c-1
                tep(b+e)=y(b+e)+y(b+e+c)*exp(-i*2*pi*e/2^d);
                tep(b+e+c)=y(b+e)-y(b+e+c)*exp(-i*2*pi*e/2^d);
            end
        end
        y=tep;
end        

%算法可以这么理解：有三个循环变量  a, b, e    观察比如8点的DFT的流图 共有三级堞形运算(a表示现在在做第几级堞形运算; ) ，第一级有4个(N/(2*c))独立的，不交叉的堞形块，
%每个是2点DFT; 第二级有两个独立的，每个是4点的DFT。每个独立的堞形块中运算的变量的递增都是1，这里是变量e表示。b是用来表示从一个堞形块到另一个堞形块元素变量的递增