基-4FFT和分裂基FFT的MATLAB仿真实现,与Walsh变换比较

HarbinInstituteofTechnology数字信号处理作业院系：电信学院通信工程系姓名：学号：哈尔滨工业大学作业一：在MATLAB上完成4096点的基-4FFT和分裂基FFT的运算程序一、4096点的基-4FFT1、假设被分析信号为一个四频信号，234sinsinsinsin2345xtttt对该信号在时间t上每隔1Ts采样一次，共取4096个采样点，得到长度为4096N的序列xn，对离散序列做基-4FFT。2、基-4FFT的MATLAB仿真结果在对基-4FFT的结果进行仿真时，为了验证其正确性，使用了MATLAB中自带的fft函数对分析信号进行基-2FFT变换，与基-4FFT变换的结果进行对比。图1-14096点的基-4FFT的频谱图图1-1是MATLAB中基-4FFT的仿真结果，从图中明显可以看出，四频信号x的基-4FFT计算出的频谱图与基-2FFT计算出的频谱完全相同，证明了基-4FFT程序正确。3、基-4FFT的MATLAB程序05001000150020002500300035004000-4-2024分析信号的时域波形050010001500200025003000350040000100020003000matlab中基-2FFT下的频谱050010001500200025003000350040000100020003000基-4FFT下的频谱clc;clear;N=4096;%DFT的点数%———生成四频的分析信号————————t=0:4095;x=sin(1*pi/2*t)+sin(2*pi/3*t)+sin(3*pi/4*t)+sin(4*pi/5*t);%四频的信号subplot(3,1,1);plot(x);axis([04096-44]);title('分析信号的时域波形');subplot(3,1,2);%绘制基-2FFT的频谱图,用于比较plot(abs(fft(x)));axis([0409603000]);title('matlab中基-2FFT下的频谱');M=log(N)/log(4);%基4FFT的分解级数，M=log4(N)Wn=exp(-2j*pi/N);%旋转因子temp=zeros(1,N);%中间临时数组，用于存储各级蝶形运算计算结果%————————四进制逆序————————n=0:N-1;screen=ones(1,N);n=bitor(bitand(n,screen*hex2dec('cccc'))/4,bitand(n,screen*hex2dec('3333'))*4);n=bitor(bitand(n,screen*hex2dec('f0f0'))/16,bitand(n,screen*hex2dec('0f0f'))*16);n=bitor(bitand(n,screen*hex2dec('ff00'))/256,bitand(n,screen*hex2dec('00ff'))*256);n=n/4^(8-M)+1;forn1=1:Ntemp(n(n1))=x(n1);endx=temp;%forL=1:M%运算级之间的循环group_cont_2=4^(M-L);%第L级数据分组数group_cont_1=4^(M-L+1);%第L-1级数据分组数group_interval_2=4^L;%第L级组间数据间隔个数，也是组内数据个数group_interval_1=4^(L-1);%第L-1级组间数据间隔个数，也是组内数据个数G=group_cont_2-1;%分组上限K=group_interval_1-1;%组内数据上限forg=0:G%每一级中包含的组循环，遍历每一组，计算各组中的数据fork0=0:K%遍历每一组中的所有数据，计算次级数据k=k0+g*group_interval_2+1;%每组数据中第一个数据序号m=group_cont_2*k0;%每一级所乘旋转因子的指数因子k1=k;%每组数据中四个数据的序号k2=k1+group_interval_1;k3=k2+group_interval_1;k4=k3+group_interval_1;X1=x(k1);X2=Wn^m*x(k2);X3=Wn^(2*m)*x(k3);X4=Wn^(3*m)*x(k4);temp(k1)=X1+X2+X3+X4;%递推公式计算结果temp(k2)=X1-1j*X2-X3+1j*X4;temp(k3)=X1-X2+X3-X4;temp(k4)=X1+1j*X2-X3-1j*X4;endendx=temp;%将temp中临时存储的第L级计算结果赋值给x%作为次级运算的输入end%—————绘制4FFT下的频谱—————————————subplot(3,1,3);plot(abs(x));axis([0409603000]);title('基-4FFT下的频谱');二、4096点的分裂基FFT1、假设被分析信号与前面相同，也是四频信号，234sinsinsinsin2345xtttt对该信号每隔1秒取4096个采样点，得到长度为4096N的序列xn，对离散序列做分裂基FFT。2、分裂基FFT的MATLAB仿真结果在对分裂基FFT的结果进行仿真时，为了验证其正确性，使用MATLAB中自带的fft函数对分析信号进行基-2FFT变换，与分裂FFT变换的结果进行对比。图1-24096点的分裂基FFT的频谱图图1-2是MATLAB中分裂基FFT的仿真结果，从图中明显可以看出，四频信号x的分裂基FFT计算出的频谱图与基-2FFT计算出的频谱完全相同，证明了分裂基FFT程序正确。3、分裂基FFT的MATLAB程序%————时域抽取法分裂基FFT算法——————clear;clc;t=0:(N-1);x=sin(1*pi/2*t+sin(2*pi/3*t+sin(3*pi/4*t+sin(4*pi/5*t;%四频输入信号subplot(3,1,1);plot(x);axis([04096-33]),title('分析信号的时域波形');y=fft(x,4096);%绘制基-2FFT的频谱图,用于验证基-4FFT结果的正确性subplot(3,1,2);plot(abs(y));axis([0409603000]);title('matlab中基-2FFT下的频谱');05001000150020002500300035004000-4-2024分析信号的时域波形050010001500200025003000350040000100020003000matlab中基-2FFT下的频谱050010001500200025003000350040000100020003000分裂基FFT下的频谱ticN=4096;%FFT点数M=log(N)/log(4);%FFT变换级数%——————倒序运算——————————————J=1;N1=N-1;forI=1:N1ifIJT=x(I);x(I)=x(J);x(J)=T;endK=N/2;whileKJJ=J-K;K=K/2;endJ=J+K;end%—————进行第一级的2点FFT运算——————IS=1;step=4;whileISNforn=IS:step:N%2点DFTa=x(n)+x(n+1);b=x(n)-x(n+1);x(n)=a;%原位计算x(n+1)=b;endIS=2*step-1;step=4*step;end%—————进行后面M-1级倒L形分裂基蝶形运算——————N2=2;forK=2:MN2=N2*2;%每一级DFT点数N4=N2/4;%每一级对应的蝶形数目，一个基本碟形运算为倒L形forJ=1:N4W_J_N2=exp(-j*2*pi*(J-1)/N2);%旋转因子；W_3J_N2=exp(-j*2*pi*3*(J-1)/N2);%旋转因子；IS=J;step=2*N2;whileISNforn=IS:step:N-1a=x(n)+x(n+2*N4)*W_J_N2+x(n+3*N4)*W_3J_N2;b=x(n+N4)-j*(x(n+2*N4)*W_J_N2-x(n+3*N4)*W_3J_N2);c=x(n)-x(n+2*N4)*W_J_N2-x(n+3*N4)*W_3J_N2;d=x(n+N4)+j*(x(n+2*N4)*W_J_N2-x(n+3*N4)*W_3J_N2);x(n)=a;x(n+N4)=b;x(n+2*N4)=c;x(n+3*N4)=d;endIS=2*step-N2+J;step=4*step;endendend%——————绘制分裂基FFT频谱图——————————subplot(3,1,3);plot(abs(x));axis([0409603000]);title('分裂基FFT下的频谱');toc作业二：比较基-2FFT，基-4FFT，分裂基FFT和Walsh变换的运算时间FFT的复杂度主要来源于计算过程中与旋转因子的乘法，基-2FFT的每个基本蝶形中需乘1次旋转因子，运算量为2log2NN；基-4FFT的每个基本蝶形中只需3次乘旋转因子，运算量为43log4NN，分裂基FFT充分结合了基-2FFT和基-4FFT的优点，对序列的偶序号项和奇数项项分别使用基-2FFT和基-4FFT，按照理论分析，分裂基FFT的运算量最小，基-4FFT次之，基-2FFT运算量最大。变换基底后有得到Walsh变换，可以由时域转化到列率域来分析。为了比较几种FFT和Walsh变换的运算量，统计几种变换的运算时间并对比。使用MATLAB中的tic/toc命令，tic用在程序的开始，作用是启动一个计时器，然后在程序尾部放一个toc，表示终止计时器，并返回tic启动以来的总时间，单位为秒（s）。1、FFT的运算时间由于软件函数库中自带的基-2fft函数和Walsh变换的fwht函数是预先编译好的内部代码，执行效率非常高，所以就会出现基-2FFT的运算时间远少于基-4FFT。例如，在做4096点基-4FFT的仿真时，4096点基-2FFT和基-4FFT的运算时间分别是：Elapsedtimeis0.013441seconds.（基-2FFT）Elapsedtimeis0.032516seconds.（基-4FFT）为了能够公平地比较几种变换的运算时间，在这里给出了一个自编的基-2FFT程序，如下所示：%基2FFT算法的MATLAB程序closeall;clearall;clc;%————生成分析信号—————t=0:4095;x=sin(1*pi/2*t+sin(2*pi/3*t+sin(3*pi/4*t+sin(4*pi/5*t;%四频输入信号subplot(2,1,1),plot(x);axis([04095-22]);title('分析信号的时域波形');tic%进行参数设置以及初始化N=4096;m=log2(N);n=0:N-1;nT=bin2dec(fliplr(dec2bin(n,m)))+1;y=x(nT);%进行基2FFT变换运算formm=1:m;nz=2^mm;u=1;wn=exp(-i*2*pi/nz);forj=1:nz/2;fork=j:nz:N;kp=k+nz/2;t=y(kp)*u;y(kp)=y(k)-t;y(k)=y(k)+t;endu=u*wn;endend%绘制基2FFT运算结果subplot(2,1,2);plot(abs(y));axis([0409603000]);title('基2FFT计算出的频谱');toc由以上程序得到的基-2FFT运算时间为：Elapsedtimeis0.129479seconds.分裂基FFT的运算时间可以直接由前面的分裂基FFT程序计算出来，运行一次程序得到