simfft.m

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

%-------------------------------------------------------------------------%
%--     固定缩小比例的FFT算法Matlab验证程序。                             --%
%--         作者：Xudong                                                --%
%--         时间：2004年5月6日                                          --%
%------------------------------------------------------------------------%
function [g1out,g2out,g3out,g4out,g5out]=simfft(indata);
%--输出数据：g1out...g5out都是8位有符号数（-128～+127）；
%--输入数据：indata也是8位有符号数（-128～+127）；

PI=3.14159265357;
DataWidth=8; %--旋转因子的精度（数据位宽）。
PointDot=32;   %--时间抽取基－2FFT的点数。

WnDot=PointDot/2; %--旋转因子数组元素个数。
theta=[0:1/WnDot*PI:(WnDot-1)/WnDot*PI];
cos_ary=cos(theta);
sin_ary=sin(theta);

wn_r_t=round(cos_ary*2^(DataWidth-1)); %--四舍五入取整。
wn_i_t=-round(sin_ary*2^(DataWidth-1));%--此处要加“－”号，为便于后面处理。

for i=1:WnDot %--产生wn_r和wn_i数组（不用补码来表示）。
    if wn_r_t(i)==2^(DataWidth-1)  %--避免对“+1”归一化时的错误。
        wn_r(i)=wn_r_t(i)-1;
    else
        wn_r(i)=wn_r_t(i);
    end;
    if wn_i_t(i)==2^(DataWidth-1)
        wn_i(i)=wn_i_t(i)-1;
    else
        wn_i(i)=wn_i_t(i);
    end;    
end;
wn=wn_r+wn_i*sqrt(-1);%--构建复数形式的旋转因子。

for i=1:PointDot %--产生仿真波形中的输入数组（不要用补码来表示）。
    if indata(i)==2^(DataWidth-1)  %--避免对“+1”归一化时的错误。
        indata_comp(i)=indata(i)-1;
    else
        indata_comp(i)=indata(i);
    end;
end;
indata_bitrev=bitrevorder(indata_comp);%--输入数据位反转。
%--第一级输出。
for i=1:2:PointDot
        [g1out(i),g1out(i+1)]=butfly(indata_bitrev(i),indata_bitrev(i+1),wn(1));
end;
%--第二级输出。
for i=1:2%--外循环：每组蝶形结的个数。
    for j=1:4:PointDot%--内循环：每组蝶形结需要循环的次数。
        [g2out(i-1+j),g2out(i-1+j+2)]=butfly(g1out(i-1+j),g1out(i-1+j+2),wn(8*(i-1)+1));
    end;
end;
%--第三级输出。
for i=1:4
    for j=1:8:PointDot
        [g3out(i-1+j),g3out(i-1+j+4)]=butfly(g2out(i-1+j),g2out(i-1+j+4),wn(4*(i-1)+1));
    end;
end;
%--第四级输出。
for i=1:8
    for j=1:16:PointDot
        [g4out(i-1+j),g4out(i-1+j+8)]=butfly(g3out(i-1+j),g3out(i-1+j+8),wn(2*(i-1)+1));
    end;
end;
%--第五级（最后一级）输出。
for i=1:16
        [g5out(i),g5out(i+16)]=butfly(g4out(i),g4out(i+16),wn(i));
end;
%--经过验证第一级和第五级的输出完全符合ISE与Modelsim中的结果；
%--幅度误差与Matlab中算出的理论值也相差很小，可是相位误差却很大，原因有待进一步分析。
%--模值小的点相位误差大，而模值大的点，例如仿真正弦信号频率点处的相位误差居然也很大。

%-------------------------------------------------------------------------%
%--幅度理论值：
%Columns 1 through 8 
%
%         0    1.3212    4.3405   61.9818    7.7471    4.0671    2.8353    2.3614
%
%  Columns 9 through 16 
%
%    1.9445    1.8409    1.5481    1.4261    1.3912    1.2839    1.2585    1.3036
%
%  Columns 17 through 24 
%
%    1.2500    1.3036    1.2585    1.2839    1.3912    1.4261    1.5481    1.8409
%
%  Columns 25 through 32 
%
%    1.9445    2.3614    2.8353    4.0671    7.7471   61.9818    4.3405    1.3212

%--“abs(g5out)”幅度值（与Modelsim中的输出值不完全吻合，主要是CORDIC模块的误差，
    %但是g5out与Modelsim中算的FFT值是完全吻合的）：
%Columns 1 through 8 
%
%    4.0000    2.8284    2.0000   58.5235    8.0623    6.4031    4.4721    4.4721
%
%  Columns 9 through 16 
%
%    2.8284    3.6056    2.8284    2.8284    1.4142    2.2361    2.2361    2.2361
%
%  Columns 17 through 24 
%
%         0    2.0000    2.0000    2.2361    1.0000    2.2361    2.8284    2.8284
%
%  Columns 25 through 32 
%
%    2.8284    2.2361    2.8284    4.4721    7.6158   61.6847    3.6056    1.0000

%--“g5out”与“fftout_xil”比较：
%--g5out
%Columns 1 through 8 
%
%  -4.0000            -2.0000 - 2.0000i        0 + 2.0000i  17.0000 +56.0000i  -4.0000 - 7.0000i  -4.0000 - 5.0000i  -2.0000 - 4.0000i  -2.0000 - 4.0000i
%
%  Columns 9 through 16 
%
%  -2.0000 - 2.0000i  -3.0000 - 2.0000i  -2.0000 - 2.0000i  -2.0000 - 2.0000i  -1.0000 - 1.0000i  -2.0000 - 1.0000i  -2.0000 - 1.0000i  -2.0000 - 1.0000i
%
%  Columns 17 through 24 
%
%        0            -2.0000            -2.0000            -1.0000 + 2.0000i        0 + 1.0000i  -2.0000 + 1.0000i  -2.0000 + 2.0000i  -2.0000 + 2.0000i
%
%  Columns 25 through 32 
%
%  -2.0000 + 2.0000i  -1.0000 + 2.0000i  -2.0000 + 2.0000i  -2.0000 + 4.0000i  -3.0000 + 7.0000i  18.0000 -59.0000i   2.0000 - 3.0000i        0 - 1.0000i
%--fftout_xil
% Columns 1 through 8 
%
%        0             0.1295 - 1.3148i   0.8468 - 4.2571i  17.9924 -59.3129i  -2.9647 + 7.1574i  -1.9172 + 3.5868i  -1.5752 + 2.3575i  -1.4980 + 1.8254i
%
%  Columns 9 through 16 
%
%  -1.3750 + 1.3750i  -1.4230 + 1.1678i  -1.2872 + 0.8601i  -1.2577 + 0.6723i  -1.2853 + 0.5324i  -1.2286 + 0.3727i  -1.2343 + 0.2455i  -1.2973 + 0.1278i
%
%  Columns 17 through 24 
%
%  -1.2500            -1.2973 - 0.1278i  -1.2343 - 0.2455i  -1.2286 - 0.3727i  -1.2853 - 0.5324i  -1.2577 - 0.6723i  -1.2872 - 0.8601i  -1.4230 - 1.1678i
%
%  Columns 25 through 32 
%
%  -1.3750 - 1.3750i  -1.4980 - 1.8254i  -1.5752 - 2.3575i  -1.9172 - 3.5868i  -2.9647 - 7.1574i  17.9924 +59.3129i   0.8468 + 4.2571i   0.1295 + 1.3148i

%--由此可见这是算法本身的缺陷导致了后面计算相位时的误差是如此的大。下面要分析产生如此大的相位误差的原因所在，以及能够采用的改进方法。
   %ffttheorytest.m文件用于分析相位误差的原因所在。2004年5月7日上午。
%--原因已经查明：是Wi中的符号差了一个负号的原因。修正后的结果在信号频点处的幅度和相位误差都已经达到了很好的精度，在没有信号的频点处
  %的相位误差仍然很大，这是因为这些点的模值太小，相位受截断误差的影响是不可避免的。2004年5月7日14:16。