mat_fft.m

附录 光盘说明 本书附赠的光盘包括各章节实例的设计工程与源码

%-------------------------------------------------------------------------%
%--     时间抽取32点基－2fft的matlab验证程序。                          --%
%--         作者：Xudong                                                --%
%--         时间：2004年5月4日                                          --%
%-------------------------------------------------------------------------%
clear;
clc;

PI=3.14159265357;
DWIDTH=8;
POINTS=32;

%i=sqrt(-1);
%indata=[39+58*i, 37+58*i, 36+5*i,  33+5*i,  36+69*i, 22+99*i, 22+65*i, 22+67*i,...%--仿真输入数据（波形文件中的数据）。
%        22+58*i, 33+53*i, 32+51*i, 32+54*i, 32+23*i, 35+21*i, 52+14*i, 96+25*i,...
%        55+22*i, 25+34*i, -20+21*i,22+20*i, 22+21*i, 22+58*i, 58+98*i, 55+56*i,...
%        69+68*i, 5+57*i,  57+59*i, 58+63*i, 5+33*i,  6+32*i,  9+31*i,  7+38*i];

%indata=[indata_inv(1) indata_inv(17) indata_inv(9)  indata_inv(25) indata_inv(5) indata_inv(21) indata_inv(13) indata_inv(29)...
%        indata_inv(3) indata_inv(19) indata_inv(11) indata_inv(27) indata_inv(7) indata_inv(23) indata_inv(15) indata_inv(31)...
%        indata_inv(2) indata_inv(18) indata_inv(10) indata_inv(26) indata_inv(6) indata_inv(22) indata_inv(14) indata_inv(30)...
%        indata_inv(4) indata_inv(20) indata_inv(12) indata_inv(28) indata_inv(8) indata_inv(24) indata_inv(16) indata_inv(32)];
%        indata_inv是前面的indata。

%--下面的程序用于产生仿真所需的输入信号（正弦信号）。
%                  ^
%                  |       *
%                  |    *     *
%                  |  *         *
%                  | *            *               
%                --|*--------------*--------------*----->
%                  |0 1 ...         *            *31
%                  |                  *         *
%                  |                    *     *
%                  |                       *

MULT_F=3;
OMIGA=2*PI/(POINTS-1)*MULT_F;
%indata=round(sin(0:OMIGA*1:2*PI*MULT_F)*2^(DWIDTH-1));
N=[0:1:31];
%indata=round(sin(OMIGA*N)*2^(DWIDTH-1));
indata=round(sin(OMIGA*N)*2^(6-1));%--考虑A/D动态范围。
for i=1:POINTS %--产生仿真波形中的输入数组（要用补码来表示）。
    if indata(i)<0
        indata_comp(i)=2^DWIDTH+indata(i); %--取负数的补码。
    elseif indata(i)==2^(DWIDTH-1)  %--避免对“+1”归一化时的错误。
        indata_comp(i)=indata(i)-1;
    else
        indata_comp(i)=indata(i);
    end;
end;

fftout_mat=fft(indata);%--Malab计算出的FFT理论输出。
fftout_xil=fftout_mat/(2^log2(POINTS));%--因为在设计中为了防止每级蝶形运算的数据逸出，
                                           %--所以每级输出数据缩小了2倍。
%fftout_xil_sample=abs([fftout_xil(1) fftout_xil(17) fftout_xil(9)  fftout_xil(25) fftout_xil(5) fftout_xil(21) fftout_xil(13) fftout_xil(29) ...
%                   fftout_xil(3) fftout_xil(19) fftout_xil(11) fftout_xil(27) fftout_xil(7) fftout_xil(23) fftout_xil(15) fftout_xil(31)]);   
           
absout=abs(fftout_xil);%--理论FFT结果的幅度值。
thetaout=angle(fftout_xil);%--理论FFT结果相位值（弧度表示）。
%--与thetaout=angle(fftout_mat)效果是一样的。

thetaout_xil=thetaout*2^(DWIDTH-3);%--按Xilin的CORDIC核中定义的数据格式缩放后的数值（波形中的理论数值，便于比较用）。

%thetaout_xil_sample=[thetaout_xil(1) thetaout_xil(17) thetaout_xil(9)  thetaout_xil(25) thetaout_xil(5) thetaout_xil(21) thetaout_xil(13) thetaout_xil(29) ...
%                   thetaout_xil(3) thetaout_xil(19) thetaout_xil(11) thetaout_xil(27) thetaout_xil(7) thetaout_xil(23) thetaout_xil(15) thetaout_xil(31)];

thetaout_du=thetaout/PI*180;%--度数表示的FFT相位理论输出结果。

%--下面的程序是为了将Modelsim中的仿真结果变换成角度的表示，以便更好的分析误差。
%--xilsimtheta_r=[  4850 -15678 -16190 -18165 -20950 -23310 -24708 -23941 -21972 ...
%--               -19181 -19927 -21461 -22926 -22415 -20438 -19181 -20182 -21972 ...
%--               -21972 -22415 -21972 -20438 -20183 -20950 -21972 -21460 -20438 ...
%--               -20182 -20182 -21972 -21972 -21972 -19181 -20182 -20182 -21460 ...
%--               -22415 -20438 -20438 -19181 -20182 -20950 -21460 -20438 -20438 ...
%--               -18421 -20182 -20950 -21460 -21972 -19181 -18421    0      0   ...
%--                  0     0      0       0      0      0      0      0      0   ... 
%--                  0];%--输出结果原始值（从波形中得来）。
            
% xilsimtheta_du=xilsimtheta_r/8192/PI*180; %--输出结果的度数表示。     
% abs_xilsimtheta_du=abs(xilsimtheta_du);
% abs_thetaout_du=abs(thetaout_du);%--理论值的绝对值。
% det_du=abs(abs_xilsimtheta_du-abs_thetaout_xil_du);%--仿真结果（FFT输出相位）的绝对值误差。

%--图形化输出。
subplot(2,1,1);
stem(indata,'*-r');
title('输入数据时域值','color','b','FontSize',8);

subplot(2,1,2);
stem(abs(fftout_xil),'*-r');
title('FFT输出理论值','color','b','FontSize',8);