《数字信号处理》课程设计报告

数字信号处理实验报告附件实验报告实验名称数字信号课程设计课程名称数字信号处理姓名汪潇成绩班级电信1501学号0903150129日期2017/9/19地点备注：一，课程设计题目和题目设计要求；课程设计题目：设计一个数字信号处理演示系统，该系统包含信号生成，频谱分析，滤波器设计、数字滤波和输出信号分析5个主要模块。题目设计要求：各模块的具体功能要求如下：⑴信号生成多频正弦：输入由多个不同频率正弦信号叠加组合而成的模拟信号公式，指定采样频率（Hz）以及采样点数，动态生成该信号的采样序列作为测试信号，用于理论验证滤波器性能。⑵频谱分析使用FFT对产生的测试信号进行频域变换，展示其幅频、相频特性，指定需要滤除或保留的频带，通过选择滤波器类型（IIR/FIR）确定对应的滤波器（低通、高通、带通、带阻）技术指标。⑶滤波器设计根据IIR/FIR数字滤波器技术指标设计滤波器，生成相应的滤波器系数，并展示对应的滤波器幅频（衰减）、相频特性。①IIRDF设计：使用双线性变换法，可选择滤波器基型（巴特沃斯/切比雪夫型）；②FIRDF设计：使用窗口法，可选择窗口类型。数字滤波根据设计的滤波器系数，对测试信号进行滤波，得到滤波后信号。①IIRDF：要求通过差分方程迭代实现滤波，未知初值置零处理；②FIRDF：要求通过快速卷积实现滤波。可以选择使用重叠相加或重叠保留法进行卷积运算，并动态展示卷积运算的详细过程。⑸输出信号分析展示滤波后信号的时域波形、幅频与相频特性，分析是否满足滤波要求。对同一滤波要求，根据输出信号频谱，对比分析各类滤波器的差异。二，设计思想和系统功能结构及功能说明；设计思想：（1）首先利用matlabGUI合理设计出程序界面，界面应包括原信号频率的文本输入部分，滤波器参数的文本输入部分，原信号及其频谱的显示部分，滤波器特性的显示部分，滤波后信号及其频谱的显示部分，以及一些相关的按钮控件；（2输入原信号频率f1,f2,f3,f4,输入采样频率，采样点数和FFT点数，对于输入原信号yt=100*sin(2*pi*f1*t)+100*sin(2*pi*f2*t)+100*sin(2*pi*f3*t)+100*sin(2*pi*f4*t)进行fft变换，将幅频相频显示出来，再根据幅频特性选择合适的滤波器进行滤波，显示滤波器特性和滤波后信号的特性，同时将理性滤波信号与滤波后信号进行对比，观察滤波效果并验证滤波器特性系统功能结构及功能说明：（1）原信号产生部分:对于任意输入的4个频率，输出yt=100*sin(2*pi*f1*t)+100*sin(2*pi*f2*t)+100*sin(2*pi*f3*t)+100*sin(2*pi*f4*t)的原始信号并绘制在坐标轴上；（2）抽样及fft变换部分输入采样频率Fs,采样点数对原信号进行t=n*T(T=1/Fs)的替换进行抽样,得到抽样序列yn,输入fft点数n0，对yn进行n0点FFT变换，利用abs()和angel()函数绘制采样信号的幅频相频特性曲线；（3）IIR滤波器设计部分IIR数字滤波器冲激响应无限宽，与模拟滤波器相匹配。因此可采取在模拟滤波器设计的基础上通过变换的方法设计IIR滤波器，即通过已知的模拟滤波器系统的系统函数H(s)来设计数字滤波器的系统函数H(z)，主要通过脉冲响应不变法或双线性变换法完成s平面到Z平面的变换。通过典型的模拟滤波器可实现一定参数要求的数字滤波器。IIR数字滤波器的设计具体步骤：1，按照既定规则将数字滤波器技术指标转为模拟低通滤波器技术指标。2，以转换后的技术指标为依据设计模拟低通滤波器H(s)（传递函数)。3，再按照一定规则将H(s)转换成H(z)（数字滤波器的传递函数）。若设计的是高通、带通或带阻滤波器，还需进行以下步骤：将高通、带通或带阻数字滤波器的技术指标转换为低通模拟滤波器的技术指标,然后设计出低通H(s),再将H(s)转换为H(z)。根据滤波器指标要求和滤波器类型输入通带（下）截止频率，阻带（下）截止频率，通带上截止频率，阻带上截止频率，以及通带衰减，阻带衰减，结合采样频率Fs得到归一化指标wp1,wp2,ws1,ws2,直接调用buttord和butter函数设计出相应的巴特沃斯数字滤波器；调用cheb1ord和cheby1函数设计相应的切比雪夫数字滤波器。调用格式如下：[N2,wc2]=buttord(wp2,ws2,ap,as);[B2,A2]=butter(N2,wc2);[N2,wc2]=cheb1ord(wp2,ws2,ap,as);[B2,A2]=cheby1(N2,ap,wc2);其中wp,ws为相应数字频率指标对pi的归一化频率，ap,as为通阻带衰减，N为滤波器阶数，wc为3dB截止频率，B为滤波器系数的分子多项式，A为滤波器系数的分母多项式。（4）IIR滤波器输出部分对于设计好的IIR滤波器，为显示其滤波器特性，调用freq函数。调用格式如下：[H,w]=freqz(B,A,N)；其中B和A分别为离散系统的系统函数分子、分母多项式的系数向量，返回量H则包含了离散系统频响在0~pi范围内N个频率等分点的值（其中N为正整数），w则包含了范围内N个频率等分点。调用默认的N时，其值是512。并利用abs()和angel()函数滤波器的幅频相频特性曲线。对于滤波信号要求使用差分迭代方法进行滤波，所以调用函数filter调用格式如下：Y=filter(B,A,X)；其中输入X为滤波前序列，Y为滤波结果序列，B/A提供滤波器系数，B为分子，A为分母然后对滤波后信号再次进行fft并输出其幅频相频特性曲线（5）FIR滤波器设计部分设计原理理同IIR滤波器部分，并且对于六种FIR滤波器分别调用boxcar，矩形窗bartlett，三角窗hanning，汉宁窗hamming，哈明窗blackman，布莱克曼窗kaiser，凯赛贝塞尔窗进行对应数字滤波器的设计，调用格式为：bt=ws-wp;long=ceil(B/bt);hn=fir1(long-1,wc,boxcar(long));其中ws,wp为数字频率对Fs的归一化，B为滤波器对应的过渡带宽，long为窗口长度，此外在设计高通带阻滤波器时要保证long为奇数，这里用n0=ceil(B/bt);long=n0+mod(n0+1,2);来保证long的合理性。（6）FIR滤波器输出部分对于设计好的FIR滤波器，为显示其滤波器特性，调用FFT函数，调用格式同上。对于滤波信号要求使用快速卷积方法进行滤波，所以调用函数fftfilt调用格式如下：ywt=fftfilt(hn,yn,long);其中hn是设计出来的窗函数，yn是要求滤波的序列，long是窗函数长度，ywt是滤波输出；然后对滤波后信号再次进行fft并输出其幅频相频特性曲线。三，设计中关键部分的详细描述和介绍，采用流程图描述关键模块的设计思路；1，原信号的输出及采样输出开始开始输入原信号频率输入原信号频率调用函数进行采样输入采样信息调用函数进行采样输入采样信息调用函数进行频谱分析调用函数进行频谱分析输出特性曲线输出特性曲线结束结束2，滤波器特性及滤波输出开始开始输入滤波器技术指标输入滤波器技术指标指标对Fs归一化指标对Fs归一化调用函数设计滤波器参数调用函数设计滤波器参数调用函数输出滤波器特性调用函数输出滤波器特性调用函数对原信号滤波调用函数对原信号滤波输出滤波后信号输出滤波后信号调用函数输出滤波后信号特性调用函数输出滤波后信号特性结束结束四，总结，包括设计过程中遇到的问题和解决方法，心得体会等；遇到的问题有：1，用户GUI界面设计不合理；2，GUI相关控件回调函数的输入参数的处理不当；3，对相关系统函数如FFT，FILTER的错误使用；4，对不同滤波器设计指标的错误理解；5，对输出波形纵横坐标范围的处理不当；解决方法：合理的安排相关控件及面板，保证GUI界面设计的合理性与实用性；理解分析相关控件回调函数的作用，合理的设置输入参数；查阅相关资料了解部分系统函数的正确使用格式查阅相关资料了解各种滤波器设计指标的差别，如IIR滤波器设计函数与FIR滤波器设计函数的输入参数是不一样的；使用坐标控制语句实现对坐标轴的控制达到良好的输出效果；本次课程设计首先考察了对数字信号处理相关理论的掌握以及实际应用的设计能力，其次考察了一部分matlab程序设计的能力；结合课本知识来理解分析课程设计的要求，明白为了满足要求需求要掌握哪些方面的理论；通过本次课程设计，我基本掌握了简单滤波器的设计思想及验证方法，了解了数字信号处理技术在实际中的部分应用，加深了对于数字信号处理技术这套理论工具的理解，最后还锻炼了个人的问题分析能力和动手操作能力；五，参考文献；[1]高西全，丁玉美.数字信号处理（第三版）[M],西安：西安电子科技大学出版社，2013.[2]高西全，丁玉美.数字信号处理（第三版）学习指导[M],西安：西安电子科技大学出版社，2013.[3]刘卫国，MATLAB程序设计与应用（第二版）[M],北京：高等教育出版社，2012.六，程序源代码清单。附件1：实验程序部分源代码function%原信号输出及采样输出函数output=origin_output(f1,f2,f3,f4,n0,Fs,fft_n,ax1,ax2,ax3,ax4,ax5)F=[f1,f2,f3,f4];max_f=max(F);T=1/Fs;long=100/max_f;%x轴长度t=linspace(0,long,800);yt=100*sin(2*pi*f1*t)+100*sin(2*pi*f2*t)+100*sin(2*pi*f3*t)+100*sin(2*pi*f4*t);plot(ax1,t,yt);title(ax1,'源信号')xlabel(ax1,'t');ylabel(ax1,'yt');axis(ax1,[0,30/max_f,-inf,inf]);%显示范围%{cla(ax4);plot(ax4,t,yt);hold(ax4,'on');axis(ax4,[0,30/max_f,-inf,inf]);cla(ax5);plot(ax5,t,yt);hold(ax5,'on');axis(ax5,[0,30/max_f,-inf,inf]);%}n=0:n0-1;yn=100*sin(2*pi*f1*n*T)+100*sin(2*pi*f2*n*T)+100*sin(2*pi*f3*n*T)+100*sin(2*pi*f4*n*T);%n点采样，默认500N=fft_n;%FFT，默认512点fft_f=[0:N-1]*Fs/N;%频率分辨率Fs/NY=fft(yn,N);fft_w=fft_f*2/Fs;plot(ax2,fft_f,2*abs(Y));title(ax2,'信号幅频特性');xlabel(ax2,'f/hz');axis(ax2,[0,Fs/2,-inf,inf]);plot(ax3,fft_w*pi,angle(Y)/pi);title(ax3,'信号相频特性');axis(ax3,[0,2,-inf,inf]);ylabel(ax3,'pi');xlabel(ax3,'w/pi');function%iir滤波器特性输出函数（低通高通）output=iirfilter_character_2(h1,h2,h3,h4,Fs,ax1,ax2,a,b)s1=get(h1,'String');%获取字符if(isempty(s1))s1='0';end;s2=get(h2,'String');if(isempty(s2))s2='0.1';end;s3=get(h3,'String');if(isempty(s3))s3='0';end;s4=get(h4,'String');if(isempty(s4))s4='60';end;fp=str2num(s1);ap=str2num(s2);fs=str2num(s3);as=str2num(s4);wp2=2*fp/Fs;ws2=2*fs/Fs;fk=linspace(0,Fs,2000);wk=2*fk/Fs;if(a==1&&b==1)[N2,wc2]=buttord(wp2,ws2,ap,as);[B2,A2]=butter(N2,wc2);elseifa==1&&b==2[N2,wc2]=buttord(wp2,ws2,ap,as);[B2,A2]=butter(N2,wc2,'high');elseifa==2&&b==1[N2,wc2]=cheb1ord(wp2,ws2,ap,as);[B2,A2]=cheby1(N2,ap,wc2);elseifa==2&&b==2[N2,wc2]=cheb1ord(wp2,ws2,ap,as);[B2,A2]=cheby1(N2,ap,wc2,'high');end;[hk,w]=freqz(B2,A2,wk);%[H，W]=freqz（b，a，n）返回n点复频响应矢量H和n点的频率向量wplot(ax1,w/pi/2*Fs,20*log10(abs(hk)));%axis(ax1,[0,1.5*fs,-inf,inf]);title(ax1,'IIR滤波器幅频特性');ylabel(ax1,'幅度/dB');xlabel(ax1,'f/hz');axis(ax1,[0,Fs/2,-100,inf]);grid(ax1,'on');plot(ax2,w,angle(hk)/pi);title(ax2,'IIR滤波器相频特性');ylabel(ax2,'pi');xlabel(ax2,'w/pi');grid(ax2,'on');axis(ax2,[0,2,-inf,inf]);function%iir滤波器滤波输出函数（低通高通）output=iirfilter_output_2(max_f,Fs,f1,f2,f3,f4,h5,h6,h7,h8,ax1,a,b,ax2,ax3,n0,fft_n)s5=get(h5,'String');if(isempty(s5))s5='0';end;s6=get(h6,'String');if(isempty(s6))s6='0.1';end;s7=get(h7,'String');if(isempty(s7))s7='0';end;s8=get(h8,'String');if(isempty(s8))s8='60';end;fp=str2num(s5);ap=str2num(s6);fs=str2num(s7);as=str2num(s8);wp2=2*fp/Fs;ws2=2*fs/Fs;T=1/Fs;n=0:n0-1;yn=100*sin(2*pi*f1*n*T)+100*sin(2*pi*f2*n*T)+100*sin(2*pi*f3*n*T)+100*sin(2*pi*f4*n*T);if(a==1&&b==1)[N2,wc2]=buttord(wp2,ws2,ap,as);[B,A]=butter(N2,wc2);elseifa==1&&b==2[N2,wc2]=buttord(wp2,ws2,ap,as);[B,A