信号分析实验报告

1、 信号分析与处理实验报告班级: 姓名：学号：一、实验目的（1）掌握计算序列的离散傅立叶变换（DFT）的方法；（2）掌握实现时间抽取快速傅立叶变换（FFT）的编程方法；（3）复习复数序列的运算方法；（4）掌握设计IIR数字滤波器的冲激响应不变法和双线性变换法；（5）掌握IIR数字滤波器的实现方法。2、 实验原理1、快速傅利叶变换的软件实现（1）程序输入元素的数目必须为2的整数次幂，即N为2的M次幂，整个运算需要M级蝶形运算。（2）输入序列按二进制码位倒置排列，输出序列按自然顺序排列。（3）输出数据占用输入数据的存储单元。（4）每一级含N/2个基本蝶形运算。（5）第L级中有N/2的L次幂个群，群与

2、群间隔为2的L次幂。（6）同一级中各个群的系数W分布相同，第L级的群中有2的L-1次幂个系数。（7）处于第L级的群的系数是(p=1,2,3，,)（8）对于第L级的蝶形运算，输入数据的间隔为2的L-1次幂。2、 码位倒置程序流程图3、 蝶形运算程序流程图 4、 IIR数字滤波器的设计以教材第九章无限冲激响应数字滤波器的例子说明软件实现方法。设滤波器为级联型结构，其系统函数如下： 式中 显然，该数字滤波器为无限冲激响应数字滤波器，其数字网络和软件流程图分别如图所示。依据软件流程图可以编制出基于C语言的实现程序。应用该程序可以计算上述数字滤波器输出序列的前N个数据。若要求该数字滤波器为具有输入输出端

3、口的实时数字滤波器，则可以将有限循环改为无限循环，并取消输出数据的长度限制。这可通过修改for循环语句来实现，即for( ; ; )然后，将输入输出数据的操作对象由数组修改为端口。y(n)y2(n)-0.375340.166760.08338z -1z -1y1(n)0.08338x(n)1.0541-0.234220.945920.166760.08338z -1z -11.314320.16676z -1z -1-0.714890.083380.083380.08338IIR数字滤波器级联型结构图5、IIR数字滤波器级联型结构实现程序流程图x(n-2)=x(n-1), x(n-1)=x(n

4、)y2(n-2)= y2(n-1), y2(n-1)=y2(n)n=n+1n<N ？No赋初值：x(-1)=0, x(-2)=0, y1(-1)=0, y1(-2)=0, y2(-1)=0, y2(-2)=0, y(-1)=0, y(-2)=0,y1(n-2)= y1(n-1), y1(n-1)=y1(n)n =0输入x(n)y1(n)=1.31432y1(n-1)-0.71489y1(n-2)+0.08338x(n)+0.16676x(n-1)+0.08338x(n-2)y2(n)=1.0541y2(n-1)-0.37534y2(n-2)+0.08338y1(n)+0.16676y1(

5、n-1)+0.08338y1(n-2)y(n)=0.94592y(n-1)-0.23422y(n-2)+0.08338y2(n)+0.16676y2(n-1)+0.08338y2(n-2)Yes结束3、 实验内容与步骤1、 为了验证所设计数字滤波器的实际功能，首先将下列分别位于通带、过渡带和阻带的三个模拟信号数字化 V ， V ， V ，式中：为采样间隔；T为记录时间长度。，数字化后的序列为： 2、 将x(n)作为激励信号输入到所设计数字滤波器的输入端，计算输出端的响应信号y(n)；3、 计算x(n)的快速傅立叶变换X(k)，和y(n)的快速傅立叶变换Y(k)；4、 用MATLAB画出x(n)

6、、y(n)、X(k)和Y(k)的图像并做比较。四、设计程序#include<math.h>#include<stdio.h>#define swap(a,b) T=(a); (a)=(b); (b)=T; /*宏定义的交换函数*/void fft(float A,float B,unsigned M) /*快速傅利叶变换程序，数组A、B分别是待变换序列 的实部和虚部*/ unsigned long N,I,J,K,L,LE,LE1,P,Q,R; float Wr,Wi,W1r,W1i,WTr,WTi,theta,Tr,Ti,T; N=1<<M; /*N=2M

7、 是序列的总长度*/ J=0; for(I=0;I<N-1;I+) /*码位倒置*/ if(J>I) swap(AI,AJ);swap(BI,BJ); K=N>>1; while(K>=2&&J>=K) J-=K; K>>=1; J+=K; for(L=1;L<=M;L+) /*外层循环由级数L控制，执行M次*/ LE=1<<L; /* LE=2L是群间隔*/ LE1=LE/2; /* LE1=2L-1是每个群的系数W数目*/ Wr=1.0; Wi=0.0; theta=(-1)*3.1415926536/LE1

8、; W1r=cos(theta); W1i=sin(theta); for(R=0;R<LE1;R+) /*中层循环由群系数控制，执行次*/ for(P=R;P<N-1;P+=LE) /*R是群系数的编号，P、Q是基本蝶形运 算两个输入数据在数组中的编号，循环每 次完成同一个系数的蝶形运算*/ Q=P+LE1; Tr=Wr*AQ-Wi*BQ; Ti=Wr*BQ+Wi*AQ; /*Tr、Ti是的实部和虚部*/ AQ=AP-Tr; /*即*/ BQ=BP-Ti; AP+=Tr; /*即*/ BP+=Ti; WTr=Wr; /*Wr、Wi是的实部和虚部*/ WTi=Wi; Wr=WTr*

9、W1r-WTi*W1i; /*用和差化积公式求的实虚部Wr、Wi*/ Wi=WTr*W1i+WTi*W1r; return; void IIRDF(float A,unsigned long N) /*IIR滤波器程序，A是采样序列，N是序列长度*/ unsigned long n; float x3=0,0,0,y13=0,0,0,y23=0,0,0,y3=0,0,0; for(n=0;n<N;n+) x0=An; y10=1.31432*y11-0.71489*y12+0.08338*x0+0.16676*x1+0.08338*x2;/*第一级*/ x2=x1; x1=x0; y20

10、=1.0541*y21-0.37534*y22+0.08338*y10+0.16676*y11+0.08338*y12; /*第二级*/ y12=y11; y11=y10; y0=0.94592*y1-0.23422*y2+0.08338*y20+0.16676*y21+0.08338*y22; /*第三级*/ y22=y21; y21=y20; y2=y1; y1=y0; An=y0;main() int i, N, M=6; float A1024,B1024,C1024,D1024,E1024,F1024,pi=3.1415927; N=1<<M; for (i=0;i<

11、;N;i+) Ai=4*sin(pi/8*i)+4*sin(pi/4*i)+4*sin(pi/2*i); /*Ai、Bi分别表示实部和虚部*/ Bi=0; Ci=Ai; printf("滤波前序列x(n)的值: n"); for(i=0;i<N;i+) printf("%f ",Ai); IIRDF(A,N); /*将x(n)用设计的低通滤波器滤波*/ printf("nn滤波后序列y(n)的值: n"); for(i=0;i<N;i+) printf("%f ",Ai); fft(C,B,M); /*

12、将滤波前的序列x(n)做FFT*/ printf("nn对滤波前序列x(n)进行FFT得X(k)的幅值: n"); for(i=0;i<N;i+) Di=sqrt(Ci*Ci+Bi*Bi); printf("%f ",Di); fft(A,B,M); /*将滤波后的序列y(n)做FFT*/ printf("nn对滤波后序列y(n)进行FFT得Y(k)的幅值: n"); for(i=0;i<N;i+) Ei=sqrt(Ai*Ai+Bi*Bi); printf("%f ",Ei); 五、实验结果与数据处理根

13、据C语言输出的结果，用MATLAB的stem函数画出x(n)、y(n)、X(k)和Y(k)的图像。MATLAB程序：x=0:1:63;y1=0.000000 8.359161 6.828427 2.523945 4.000000 4.867091 -1.171574 -5.297693 0.000002 5.297694 1.171571 -4.867092 -3.999999 -2.523946 -6.828429 -8.359159 0.000005 8.359162 6.828424 2.523944 4.000002 4.867089 -1.171577 -5.297693 0.000

14、005 5.297694 1.171568 -4.867093 -3.999998 -2.523947 -6.828432 -8.359158 0.000010 8.359164 6.828422 2.523944 4.000003 4.867088 -1.171581 -5.297693 0.000009 5.297694 1.171564 -4.867095 -3.999996 -2.523947 -6.828434 -8.359156 0.000015 8.359166 6.828420 2.523943 4.000005 4.867087 -1.171584 -5.297692 0.0

15、00012 5.297695 1.171561 -4.867095 -3.999995 -2.523948 -6.828436 -8.359155;figure, subplot(221),stem(x,y1),title('滤波前序列x(n)的值');x=0:1:63;y2=0.000000 0.004846 0.049092 0.236617 0.728676 1.627885 2.840505 4.054023 4.844191 4.839596 3.907749 2.294106 0.557038 -0.743326 -1.377211 -1.531653 -1.630

16、350 -2.046460 -2.840957 -3.640089 -3.789852 -2.771374 -0.592484 2.143054 4.504156 5.645319 5.186960 3.413738 1.136026 -0.765726 -1.811080 -2.082345 -2.043807 -2.197339 -2.744945 -3.406926 -3.552500 -2.626251 -0.571409 2.067082 4.389311 5.548735 5.142144 3.423888 1.181404 -0.713346 -1.774680 -2.07195

17、3 -2.056172 -2.221020 -2.767230 -3.419286 -3.552813 -2.617826 -0.560111 2.075908 4.392835 5.547058 5.137420 3.418879 1.178197 -0.713980 -1.773221 -2.069581;subplot(222),stem(x,y2),title('滤波后序列y(n)的值');x=0:1:63;y3=0.000001 0.000009 0.000012 0.000021 127.999985 0.000012 0.000012 0.000029 127.9

18、99954 0.000028 0.000012 0.000010 0.000008 0.000015 0.000022 0.000052 127.999985 0.000061 0.000034 0.000025 0.000019 0.000016 0.000013 0.000013 0.000018 0.000011 0.000014 0.000011 0.000017 0.000009 0.000013 0.000010 0.000014 0.000010 0.000013 0.000009 0.000016 0.000011 0.000014 0.000011 0.000025 0.00

19、0013 0.000013 0.000016 0.000022 0.000025 0.000034 0.000061 127.999985 0.000052 0.000022 0.000015 0.000014 0.000010 0.000012 0.000028 127.999954 0.000029 0.000012 0.000012 127.999969 0.000021 0.000012 0.000009;subplot(223),stem(x,y3),title('对滤波前序列x(n)进行FFT得X(k)的幅值');x=0:1:63;y4=22.439137 22.2

20、89307 21.832138 21.043167 113.883820 18.160545 15.252707 9.832901 49.362049 1.986152 1.817431 1.849492 1.803223 1.718803 1.627275 1.541272 1.221894 1.397443 1.339100 1.288442 1.244439 1.206240 1.173103 1.144456 1.119829 1.098827 1.081161 1.066581 1.054885 1.045951 1.039638 1.035882 1.034635 1.035880 1.039637 1.045949 1.054895 1.066580 1.081162 1.098828 1.119828 1.144456 1.173103 1.206239 1.244431 1.288443 1.