信号处理IIRFIRFFT算法滤波器.doc
信号处理课程设计题1：自相关与FFT得分一、（25分）已知某正弦波信号混有白噪声信号（正弦波信号频率=学号后两位/100，其他参数任选）。对上述信号进行自相关和FFT变换。1、设定正选信号的频率为79/100=0.79HZ，抽样频率为100HZ；2、设定N(0,0.25)高斯白噪声，及噪声功率为0.25W；要求：1.给出Matlab 源程序及说明；2.给出原始波形、自相关函数图形和FFT变换后的图形；3.根据结果分析信号当中混有的主要频率成分；4.对比说明上述两种信号分析方法各有什么特点。1.答：源程序及说明fs=100;fc=0.79;%频率为学号/100的带有白噪声的正弦信号x=(0:1/fs:2);n=201;y1=sin(2*pi*fc*x); %原正弦信号，频率为10a=0;b=0.5; %均值为a，方差为b2subplot(3,2,1);plot(x,y1,'r');title('y=sin(20pi*x)');ylabel('y');xlabel('x/20pi');grid;y2=a+b*rand(1,n); %均匀白噪声subplot(3,2,2);plot(x,y2,'r');title('N(0，0.25)的均匀白噪声');ylabel('y');xlabel('x/20pi');grid;y=y1+y2; %加入噪声之后的信号subplot(3,2,3);plot(x,y,'r');title('叠加了均匀白噪声的sinx');ylabel('y');xlabel('x/20pi');grid;FY=fft(y); %傅里叶变换得出频谱函数FY1=fftshift(FY); %频谱校正f=(0:200)*fs/n-fs/2;subplot(3,2,4);plot(f,abs(FY1),'r');title('函数频谱图');ylabel('F(jw)');xlabel('w');grid;2.答：Butterworth数字滤波器的模拟滤波器的传递函数：数字滤波器的转移函数为：根据以下运行结果：bs2 =9.9565*10-4,0,0,0，as2 = 1.0000,0.0946,.0233,.0010，bz2 = 0.0009,0.0027,0.0027,0.0009，az2 = 1.0000, -2.7423,2.5786, -0.8291Chebyshev字滤波器的模拟滤波器的传递函数：数字滤波器的转移函数为：3.答：根据第三个图,-1和0.8Hz为主要成分4.答：自相关函数：1）在0点的值最大；之后变小，2）若信号中有周期成分，则自相关函数也有周期性，且不衰减！如：正弦信号的自相关函数为余弦函数；3）若信号中无周期成分，自相关函数一般衰减到均方值（未去直流）或0(在信号中去掉直流成分) FFT变换：大大的减少了运算量题2: IIR 滤波器设计得分二、（25分）请设计一个Butterworth及Chebyshev低通IIR数字滤波器，其具体技术指标为：通带截止频率为50Hz，对应衰减3dB, 阻带截止频率为80Hz，对应的衰减为20dB, 采样频率为1000Hz。要求：1.给出设计滤波器的Matlab源程序及说明；2.根据运行结果分别给出模拟滤波器的传递函数及数字滤波器的转移函数并画出相应的数字滤波器频率特性图；3.用所设计的滤波器对题1中的信号进行滤波处理，给出源程序及滤波前后信号的图形；4.利用仿真分析改变指标（衰减或者截止频率），上述结果有什么变化。程序：function ButterworthFs=1000;%根据采样频率将滤波器边界频率进行转换Wp=50;Ws=80;Rp=3;Rs=20;Fs=1000;Nn=256;%调用freqz所用频率的点数N,Wn=buttord(Wp,Ws,Rp,Rs,'s');%模拟滤波器的最小阶数z,p,k=buttap(N);Bap,Aap=zp2tf(z,p,k);%将零极点增益形式转换成传递函数形式b,a=lp2lp(Bap,Aap,Wn);%经频率转换bz,az=impinvar(b,a,Fs);%脉冲响应不变法得到的传递函数figure(1)H,f=freqz(bz,az,Nn,Fs) %绘制数字滤波器的幅频特性和相频特性subplot(2,1,1),plot(f,20*log10(abs(H);xlabel('f');ylabel('振幅');subplot(2,1,2),plot(f,180/pi*unwrap(angle(H);xlabel('f');ylabel('相位');figure(2) N=1000,Fs=20;%数据长度和采样频率n=0:N-1;t=n/Fs;%时间序列lag=50;randn('state',0);x=sin(2*0.79*pi*t)+0.3*randn(1,length(t);%频率是学号/100的带有白噪声的正弦信号subplot(2,1,1),plot(t,x);xlabel('t');ylabel('x(t)');title('u');grid on;y=filtfilt(bz,az,x);%用函数filtfilt对输入信号进行滤波subplot(2,1,2),plot(t,y);,title('y（t）');,xlabel('t')程序2：function ChebyshevWp=50;Ws=80;Rp=3;Rs=20;Fs=1000;Nn=256;%调用freqz所用频率的点数wp=50*2/Fs;ws=80*2/Fs;N,Wn=cheb1ord(wp,ws,Rp,Rs);%求最小阶数和归一化截止频率b,a=cheby1(N,Rp,Wn);%按最小阶数，通带波纹和截止频率figure(1)H,f=freqz(b,a,Nn,Fs) %绘制滤波器的频率特征subplot(2,1,1),plot(f,20*log10(abs(H);xlabel('f');ylabel('振幅');subplot(2,1,2),plot(f,180/pi*unwrap(angle(H);xlabel('f');ylabel('相位');figure(2)N=1000,Fs=20;n=0:N-1;t=n/Fs;%时间序列lag=50;%延迟样点数randn('state',0);x=sin(2*0.79*pi*t)+0.3*randn(1,length(t);%频率为学号/100的带有白噪声的正弦信号subplot(2,1,1),plot(t,x);xlabel('t');ylabel('x(t)');title('u');grid on;y=filtfilt(b,a,x);%用函数filtfilt对输入信号进行滤波subplot(2,1,2),plot(t,y);,title('y（t）');,xlabel('t')Butterworth Chebyshev 改变衰减阻带从10增加到25由此可得切比雪夫型的滤波器性能较好。如果阻带的截止频率变小，会使滤波的性能变好题3：FIR滤波器设计得分三、（25分）请用窗函数法（矩形窗及汉明窗）及Parks-McClellan方法分别设计FIR低通数字滤波器，其具体技术指标为：c=0.3，窗长度N取学号后两位。(若学号小于10，N可以取学号乘10)要求：1.给出设计滤波器的Matlab源程序及说明；2.给出所设计窗函数的幅频特性及滤波器的频率特性图；3.利用仿真分析改变参数（窗长度N或者通带及阻带纹波）对滤波器频率特性有何影响；4.根据结果对比总结上述各种FIR滤波器的特点。程序如下：function hamming11Wc=0.3*pi;N=79;%窗长度N取学号后两位Nf=512;%窗函数复数频率特性的数据点数b=fir1(N,Wc/pi,hamming(N+1);%哈明窗H,f=freqz(b,1,Nf,50) % 求出频率响应，采样频率为50Hz，subplot(2,1,1),plot(f,20*log10(abs(H);xlabel('频率');ylabel('振幅');subplot(2,1,2),plot(f,180/pi*unwrap(angle(H);xlabel('频率');ylabel('相位');figure(2)b=fir1(N,Wc/pi,boxcar(N+1);%矩形窗H,f=freqz(b,1,512,50) % 求出频率响应，采样频率为50Hz，subplot(2,1,1),plot(f,20*log10(abs(H);xlabel('频率');ylabel('振幅');subplot(2,1,2),plot(f,180/pi*unwrap(angle(H);xlabel('频率');ylabel('相位');汉明窗矩形窗： 程序2：function ParksMcClellanclf;n=20;%滤波器阶数f=0 0.4 0.5 1;a=1 1 0 0;b=remez(n,f,a);%采用remez设计滤波器h,w=freqz(b);%计算滤波器的频率响应plot(w/pi,abs(h);Parks-McClellan图像改变窗体长度从10到70；图像汉明窗矩形窗可以看出床窗的长度越窄，滤波器的效果越好，上述图形分别是根据窗函数法（矩形窗及汉明窗）及Parks-McClellan方法分别设计FIR低通数字滤波器，根据结果可以看出如果是设计同一阶数的滤波器使用Parks-McClelland的优化方法较窗函数设计更加优越，能取得更好的衰耗特性，使用Parks-McClelland的优化方法可以获得最佳的频率特性，而且一定程度上还具有通带和阻带平坦过渡带不宽等优点，即所设计出的滤波器滤波效果更好。 利用窗函数方法设计低通滤波器时，不同的窗函数所设计的滤波器效果不一样，而且阶数的大小对不同的窗函数所达到的设计效果也有影响，阶数较低时，阻带特性不满足设计要求，只有当阶数较高时，使用汉明窗窗口法设计基本可以达到所需的阻带衰耗要求。题4：信号处理的应用报告得分补偿零点的量化我们所提出的补偿零点量化法，在对第二级子模块进行量化时，会考虑第一子模块的量化误差影响该方法的关键点是让量化后级联结构的频率响应与初始未量化的直接结构的频率响应理想匹配。采用补偿零点进行量化时，先基于上一节描述的简单量化方法，把第一级级联模块C1（n）量化为C1(n)，把增益因子k量化为k 。与简单量化方法不同的是，接下来我们并不对C2（n）进行量化，而是计算一个最佳的Ccomp(n)，使得Ccomp(n)的目标是要补偿由C1(n)量化引起的性能损失，因此Ccomp(n)的计算过程。量化前级联滤波器的传递函数H(z)=kC1(z)C2(z)式中k为增益因子。采用补偿零点方法得到的半量化滤波器传递函数为：HCOMP(z)=kC1(z)Ccomp(z)使以上两式相等kC1(z)Ccomp(z)=kC1(z)C2(z)，这样就可以求解得到Ccomp(n)1、根据给定的直接型非量化系数h（n），计算量化前级联结构的各系数c1（n）、c2（n）和k。2、利用简单量化方法，将c1（n）量化为c1(n),k量化为k。3、选择M个唯一的正品率点。当滤波器具有对称结构时，M=【（N-2）/2】，否则M=【（N-1）/2】，其中N为第二级联子模块c2（n）的长度。4、基于式（3）求解M个未知的系数Ccomp（n）。当滤波器具有对称结构时，通过将步骤3中选定的M个正频率分量代入式（3）来进行求解。反之，当滤波器为非对称结构时，需要将M个正频率分量和对应的M个负频率分量一同代入式（3）进行求解。5、利用简单量化方法对浮点型系数Ccomp（n）进行量化得到c2(n)，期间将使用步骤2执行后剩余的T。6、由公式H(z)=kC1(z)C2(z),计算补偿零点量化滤波器的直接系数
