disp("請輸入一個128點序列");for ii=1:128 %使用者可以自由輸入序列x(ii) = input(["x(",num2str(ii),")="]);end%整體運用原位計算m=nextpow2(x);N=2^m; % 求x的長度對應的2的最低冪次mif length(x)<N x=[x,zeros(1,N-length(x))]; % 若x的長度不是2的冪，補零到2的整數冪end nxd=bin2dec(fliplr(dec2bin([1:N]-1,m)))+1; % 求1:2^m數列序號的倒序y=x(nxd); % 將x倒序排列作為y的初始值for mm=1:m % 將DFT作m次基2分解,從左到右，對每次分解作DFT運算,共做m級蝶形運算，每一級都有2^(mm-1)個蝶形結Nz=2^mm;u=1; % 旋轉因子u初始化為WN^0=1WN=exp(-i*2*pi/Nz); % 本次分解的基本DFT因子WN=exp(-i*2*pi/Nz)for j=1:Nz/2 % 本次跨越間隔內的各次蝶形運算,在進行第mm級運算時需要2^(mm-1)個 蝶形for k=j:Nz:N % 本次蝶形運算的跨越間隔為Nz=2^mmkp=k+Nz/2; % 蝶形運算的兩個因子對應單元下標的關係t=y(kp)*u; % 蝶形運算的乘積項y(kp)=y(k)-t; % 蝶形運算y(k)=y(k)+t; % 蝶形運算endu=u*WN; % 修改旋轉因子,多乘一個基本DFT因子WNendendyy1=fft(x) %自己編的FFT跟直接呼叫的函式運算以後的結果進行對比
disp("請輸入一個128點序列");for ii=1:128 %使用者可以自由輸入序列x(ii) = input(["x(",num2str(ii),")="]);end%整體運用原位計算m=nextpow2(x);N=2^m; % 求x的長度對應的2的最低冪次mif length(x)<N x=[x,zeros(1,N-length(x))]; % 若x的長度不是2的冪，補零到2的整數冪end nxd=bin2dec(fliplr(dec2bin([1:N]-1,m)))+1; % 求1:2^m數列序號的倒序y=x(nxd); % 將x倒序排列作為y的初始值for mm=1:m % 將DFT作m次基2分解,從左到右，對每次分解作DFT運算,共做m級蝶形運算，每一級都有2^(mm-1)個蝶形結Nz=2^mm;u=1; % 旋轉因子u初始化為WN^0=1WN=exp(-i*2*pi/Nz); % 本次分解的基本DFT因子WN=exp(-i*2*pi/Nz)for j=1:Nz/2 % 本次跨越間隔內的各次蝶形運算,在進行第mm級運算時需要2^(mm-1)個 蝶形for k=j:Nz:N % 本次蝶形運算的跨越間隔為Nz=2^mmkp=k+Nz/2; % 蝶形運算的兩個因子對應單元下標的關係t=y(kp)*u; % 蝶形運算的乘積項y(kp)=y(k)-t; % 蝶形運算y(k)=y(k)+t; % 蝶形運算endu=u*WN; % 修改旋轉因子,多乘一個基本DFT因子WNendendyy1=fft(x) %自己編的FFT跟直接呼叫的函式運算以後的結果進行對比