1、数字信号处理实习报告中国地质大学（武汉）数字信号处理上机实习学生姓名： 班 级： 学 号： 指导老师： 题目一 离散卷积计算1、实验题目设线性时不变（LTI）系统的冲激响应为h(n)，输入序列为x(n)1、h(n)=(0.8)n，0n4； x(n)=u(n)-u(n-4)2、h(n)=(0.8)nu(n), x(n)=u(n)-u(n-4)3、h(n)=(0.8)nu(n), x(n)=u(n)求以上三种情况下系统的输出y(n)，显示输入和输出波形。2、实验目的1.理解和掌握离散卷积计算；2.学习如何用Mtalab实现离散卷积计算。三、算法设计离散卷积定义为：1、， (a) 当时，；(b) 当
2、时, (0.8)n ；(c) 当时， (0.8)n；(d) 当时，；2、， (a) 当时，；(b) 当时, (0.8)n； (c) 当时， (0.8)n； (d) 当时， (0.8)n；(e) 当时，；3、， (a) 当时，； (b) 当时, (0.8)n； (c) 当时， (0.8)n；(d) 当时，；4、程序分析所用到的函数：（1）y=conv（x.,h）：卷积运算函数，计算;（2）n1=0:4：n1取04；（3）subplot(m,n,p):subplot()函数是将多个图画到一个平面上的工具。其中，m表示是图排成m行，n表示图排成n列，也就是整个figure中有n个图是排成一行的，一共
3、m行，如果m=2就是表示2行图。p表示图所在的位置，p=1表示从左到右从上到下的第一个位置。（4）title(content ):title()函数的功能是为当前坐标系添加标题“content ”。五、程序设计n=0:4;h=0.8.n;x=1 1 1 1subplot(331);stem(x); title(x(n);y=conv(x,h);subplot(332);stem(h); title(h(n);subplot(333);stem(y); title(y(n); n=0:40;h=0.8.n;x=1 1 1 1subplot(334);stem(x); title(x(n);y=c
4、onv(x,h);subplot(335);stem(h); title(h(n);subplot(336);stem(y); title(y(n); n=0:40;h=0.8.n;x=zeros(1,0),ones(1,40);subplot(337);stem(x); title(x(n);y=conv(x,h);subplot(338);stem(h); title(h(n);subplot(339);stem(y); title(y(n);六、运行结果图中从左至右三列依次对应、及卷积结果题目二 离散傅立叶变换1、实验题目设有离散序列 分析下列三种情况下的幅频特性。（1） 采集数据长度N
5、=16，分析16点的频谱，并画出幅频特性。采集数据长度N=16，并补零到64点，分析其频谱，并画出幅频特性。（2） 采集数据长度N=64，分析64点的频谱，并画出幅频特性。观察三幅不同的幅频特性图，分析和比较它们的特点及形成原因。2、实验目的1、理解掌握DFT及FFT算法；2、利用FFT算法计算信号的频谱。三、算法设计当抽样数N=2M时，以下为蝶形算法图。 1、 当N=2M时，则要进行M次分解，即进行M级蝶形单元的计算；2、按自然顺序输入，输出是码位倒置；3、每一级包含N/2个基本蝶形运算；4、第L级有2L-1个蝶群，蝶群间隔为N/2L-1；四、程序分析所用到的函数：（1）fft（）：求x的一
6、维傅里叶变换，计算X（k）；使用方法Xk=fft(xn)，例如：N=8; n=0:N-1; xn=4 3 2 6 7 8 9 0; Xk=fft(xn)（2）stem（）:显示函数图像。（3）xlabel():横坐标轴的名称；（4）ylabel():纵坐标轴的名称；五、程序设计%x(n)=cos(0.48n)+cos(0.52n)n=0:15; %产生序列x(n) 取16点x=cos(0.48*pi*n)+cos(0.52*pi*n);subplot(231);stem(n,x);title(采集数据长度N=16); %显示x(n)x1k=fft(x,16);%进行16点傅里叶变换subplo
7、t(234);stem(0:15,abs(x1k),.); %显示X（k）xlabel(n);ylabel(X1(k); title(16点傅立叶变换频谱); x=x zeros(1,48); %取16点，补零到64点x2k=fft(x,64);%进行64点傅里叶变换subplot(232);stem(x);title(采集N=16,补零到64);subplot(235);stem(0:63,abs(x2k),.);xlabel(n);ylabel(X2(k);title(64点傅立叶变换频谱); n=0:63; %产生序列x(n)，取64点x=cos(0.48*pi*n)+cos(0.52*
8、pi*n);subplot(233);stem(n,x);title(采集数据长度N=64);x3k=fft(x);%做64点傅里叶变换subplot(236);stem(0:63,abs(x3k),.);xlabel(n);ylabel(X3(k);title(64点傅立叶变换频谱);六、运行结果七结果分析 N点DFT的频谱分辨率是2 /N。一节指出可以通过补零观察到更多的频点，但是这并不意味着补零能够提高真正的频谱分辨率。这是因为xn 实际上是x(t) 采样的主值序列，而将xn补零得到的xn 周期延拓之后与原来的序列并不相同，也不是x(t) 的采样。因此是不同离散信号的频谱。对于补零至M点
9、的x的DFT，只能说它的分辨率2 /M仅具有计算上的意义，并不是真正的、物理意义上的频谱。频谱分辨率的提高只能通过提高采样频率实现。第三幅图形取样点数较多，分辨率最高，而第一幅图形则显得较为稀疏，主要是因为取样点数太少。八.专题实习心得离散傅里叶变换是一种快速算法，由于有限长序列在其频域也可离散化为有限长序列，因此离散傅里叶变换在数字信号处理中是非常有用的。DFT是重要的变换，在分析有限长序列的有用工具、信号处理的理论上有重要意义、运算方法上起核心作用，谱分析、卷积、相关都可以通DFT在计算机上实现。通过编程实践体会到了时域、频域信号的对应关系，也对采样频率的含义有了深刻的认识，同时也加深了对
10、采样信号频谱周期性的理解。 题目三 IIR滤波器的设计1、实验题目1、设计一个Butterworth数字低通滤波器，设计指标如下： 通带截止频率：0.2，幅度衰减不大于1分贝 阻带截止频率：0.3，幅度衰减大于15分贝2、分析不同滤波器的特点和结果。3、编程设计实现IIR滤波器。二、实验目的1.理解和掌握不同IIR滤波器的性质、特点。2.掌握用Matlab设计butterworth滤波器和椭圆滤波器的方法。三、算法设计1、滤波器类型Butterworth滤波器： Butterworth滤波器的特点是在通带内的频率特性是平坦的,并且随着频率的增加而衰减。Butterworth滤波器又是最简单的滤
11、波器。N阶低通Butterworth滤波器的幅度平方函数为:椭圆滤波器：椭圆滤波器在通带和阻带内都是等纹波振荡。椭圆滤波器的特性函数为：2、变换方法（a）冲激响应不变法冲激响应不变法的基本原理是从滤波器的冲激响应出发，对模拟滤波器冲激响应h(t)进行取样，所得到的离散序列h(nT)作为数字滤波器的单位取样响应。H(z)是由H(s)通过下式的对应关系得到。（b）双线性变换是在所得到满足性能指标要求的模拟滤波器的基础上,通过变换,从而得到相应的数字滤波器。4、程序分析所用到的函数：（1）N,wn=buttord(Wp,Ws,ap,as,s)：阶数，截止频率=buttord（通带截止频率，阻带截止频
12、率，通带衰减，阻带衰减，模拟）；（2）num,den=bilinear(b,a,fs)：双线性变换可选频率预畸变，在IIR数字滤波器设计中有把模拟的原型滤波器转变成相应的数字滤波器，即通过已知S变换中模拟滤波器的系数求数字滤波器的系数。其中有两种方法来进行这种变换，一种是激冲不变法，另一种是双线性Z变换法。bilinear 函数便是双线性Z变换法。函数中的b和a是原模拟滤波器S变换中的分子和分母系数，num和den是数字滤波器的系数。（3）grid on：增加主要网格线为当前轴；（4）h,w=freqz(num,den)：求频率响应。函数的输出：a.滤波器的频率响应H（N点）b.频率向量W(N
13、点，且单位为弧度)，频率向量W是均匀分布在滤波器的上半区，即：0：pi，这些点上的频率响应都将通过此函数计算出来。；滤波器的系数：分子为num，分母为den。 五、程序设计clc;Rs=15; Rp=1;Wp1=0.2*pi;Ws1=0.3*pi;Op1=tan(Wp1/2); %预畸变Op2=tan(Ws1/2);N,Wn=ellipord(Op1,Op2,Rp,Rs,s); %阶数及截止频率num1,den1=ellip(N,Rp,Rs,Wn,s); %设计椭圆低通数字滤波器，返回传输函数num,den=bilinear(num1,den1,0.5); %双线性变换subplot(221)
14、;zplane(num,den); %绘制零极点图title(零极点图)w=0:pi/256:pi;h=freqz(num,den,w);g=abs(h);g1=angle(h); subplot(223);plot(w/pi,g);gridaxis(0 1 0 2);xlabel(频率); ylabel(幅度); title(幅频响应);subplot(224);plot(w/pi,g1);grid onaxis(0 1 -5 5);xlabel(频率); ylabel(相位); title(相频响应); %测试1t=1:300;I=sin(0.1*pi*t)+sin(0.55*pi*t)+
15、sin(0.9*pi*t);%设计正弦波figure;subplot(311);plot(I);title(原信号x(n)=sin(0.1*pi*t)+sin(0.55*pi*t)+sin(0.9*pi*t);A=filter(num,den,I);%正弦波通过滤波器subplot(312);plot(A);title(滤波后y(n);w=(0:255)/256*(0.5/2);subplot(313)plot(w,abs(S(1:256) SF(1:256);title(傅里叶变换后频谱图对比); %测试2t=1:300;I=sin(0.1*pi*t)+sin(0.8*pi*t)+sin(0
16、.95*pi*t);%设计正弦波figure;subplot(311);plot(I);title(原信号x(n)=sin(0.1*pi*t)+sin(0.8*pi*t)+sin(0.95*pi*t);A=filter(num,den,I);%正弦波通过滤波器subplot(312);plot(A);title(滤波后y(n);w=(0:255)/256*(0.5/2);subplot(313)plot(w,abs(S(1:256) SF(1:256);title(傅里叶变换后频谱图对比);6、运行结果（1）零极点图、幅频响应、相频响应(2)输入为I=sin(0.1*pi*t)+sin(0.5
17、5*pi*t)+sin(0.9*pi*t)时的滤波结果(3)输入为I=sin(0.1*pi*t)+sin(0.8*pi*t)+sin(0.95*pi*t)时的滤波结果7、结果分析由幅频特性可知，椭圆滤波器在通带内的频率特性是平坦的,并且随着频率的增加而衰减。通带截止频率：0.2，幅度衰减不大于1分贝，阻带截止频率：0.3，幅度衰减大于15分贝。正弦信号在经过IIR滤波器滤波后，由傅里叶变换后的频谱图可看出高频信号被滤除，低频信号被保留了下来。8、专题实习心得通过本专题的学习，熟悉和巩固了Butterworth数字低通滤波器和椭圆滤波器的设计方法和原理，实现滤波器设计的有关经典算法，更重要的是熟
18、练掌握使用MATLAB语言设计各种要求的数字滤波器。这一部分的内容相对来说是有些难度了，做起来花费的精力也多了一些，不过对数字信号处理内容的掌握上又加深了一层。题目四 FIR滤波器的设计一、实验内容选取合适窗函数设计一个线性相位FIR低通滤波器，使它满足如下性能指标：通带截止频率：p=0.5,通带截止频率处的衰减不大于3分贝；阻带截止频率：s=0.66，阻带衰减不小于40分贝二、实验目的1、掌握用窗函数法设计FIR滤波器的原理和方法。2、熟悉线性相位滤波器特性。3、了解各种窗函数对滤波器特性的影响。三、原理及算法概要1、算法通过其通带截止频率p与阻带截止频率s算出其过渡带的宽度与滤波器的长度，
19、从而得到理想滤波器的截止频率，根据所要求的理想滤波器，得到hd（n）。由于其通带截止频率处的衰减不大于3分贝与阻带衰减不小于40分贝，我选择最接近的汉宁窗，最后调用函数h=hd.*win 及freqz(h,1,512)得到实际汉宁窗的响应和实际滤波器的幅度响应。2、原理概要利用窗函数法设计FIR 滤波器，FIR滤波器的最大特点是其相位特性可以设计为严格的线性，而其幅值可以任意设置，这样输出波形就不会相位失真。理想低通滤波器的单位取样响应hd（n）是无限长的，所以要用一个有限长的因果序列h（n）进行逼近，最有效的方法是截断hd（n），即用有限长的窗函数w（n）来截取hd（n），表示为h（n）=h
20、d（n）w（n）。为获得线性相位的FIR滤波器，h（n）必须满足中心对称条件，序列h（n）应有一定的延迟，且=（N-1）/2频率响应逼近hd（ejw）的FIR滤波器，最简单的窗函数为矩形窗： 1 n(N-1)/2加窗后的频谱加窗后使实际频响偏离理想频响，影响主要有两个方面：（1） 通带和阻带之间存在过渡带，过渡带宽度取决于窗函数频响的主瓣宽度。（2） 通带和阻带区间有纹波，这是由窗函数的旁瓣引起的，旁瓣越多，纹波越多。增加窗函数的宽度N，其主瓣宽度减小，但不改变旁瓣的相对值。为了改善滤波器的性能，要求窗函数的主瓣宽度尽可能窄，以获得较窄的过渡带；旁瓣衰减尽可能大，数量尽可能大，从而改善纹波状况
21、，使实际频响H（ej）更好地逼近理想频响Hd（ej）。除了矩形窗外，一般还可以采用以下几种窗函数 汉宁窗： 海明窗布来克曼窗4、程序设计 clc; wp=0.5*pi;ws=0.66*pi;wd=ws-wp;%hanning窗N_hann=ceil(8*pi/wd); wc=(wp+ws)/2;if mod(N_hann,2)=0 N_hann=N_hann+1;endhann_win=hanning(N_hann);b_hann=fir1(N_hann-1,wc/pi,hann_win);freq_axis=0:pi/512:pi-pi/512;freqn=0:511/512;H_hann=
22、freqz(b_hann,1,512);figure(1);subplot(311);plot(freqn,20*log10(abs(H_hann);xlabel(归一化频率w/pi);ylabel(幅度谱（dB）);title(采用hanning窗设计fir滤波器幅度谱);subplot(312);plot(freqn,angle(H_hann);xlabel(归一化频率w/pi);ylabel(相位谱);title(采用hanning窗设计fir滤波器相位谱);subplot(313);stem(b_hann);xlabel(n);ylabel(h(n);title(单位脉冲响应序列);%
23、hamming窗N_hamm=ceil(8*pi/wd);wc=(wp+ws)/2;if mod(N_hamm,2)=0 N_hamm=N_hamm+1;endhamm_win=hamming(N_hamm);b_hamm=fir1(N_hamm-1,wc/pi,hamm_win);H_hamm=freqz(b_hamm,1,512);figure(2); subplot(311);plot(freqn,20*log10(abs(H_hamm);xlabel(归一化频率w/pi);ylabel(幅度谱（dB）);title(采用hamming窗设计fir滤波器幅度谱);subplot(312)
24、;plot(freqn,angle(H_hamm);xlabel(归一化频率w/pi);ylabel(相位谱);title(采用hamming窗设计fir滤波器相位谱);subplot(313);stem(b_hann);xlabel(n);ylabel(h(n);title(单位脉冲响应序列);%N*2N2_hamm=2*N_hamm;wc=(wp+ws)/2;if mod(N2_hamm,2)=0 N2_hamm=N2_hamm+1;end%hamm_win=hamming(N2_hamm);b_hamm=fir1(N2_hamm-1,wc/pi);H_hamm=freqz(b_hamm,
25、1,512);figure(3); subplot(311);plot(freqn,20*log10(abs(H_hamm);xlabel(归一化频率w/pi);ylabel(幅度谱（dB）);title(窗长增大一倍，采用hamming窗设计fir滤波器幅度谱);subplot(312);plot(freqn,angle(H_hamm);xlabel(归一化频率w/pi);ylabel(相位谱);title(窗长增大一倍，采用hamming窗设计fir滤波器相位谱);subplot(313);stem(b_hann);xlabel(n);ylabel(h(n);title(单位脉冲响应序列)
26、;5程序运行结果（1）采用hanning窗设计FIR滤波器（2）采用hamming窗设计FIR滤波器(3)窗长增大一倍，采用hamming窗设计FIR滤波器6结果分析上图为选用汉宁窗和汉明窗的FIR 滤波器的相频特性、幅频特性及单位取样响应hd（n）。由运行得出的幅频特性可知，满足设计指标通带截止频率：p=0.5；阻带截止频率：s=0.66。得出的窗主瓣宽度较宽，旁瓣宽度较窄，也合乎要求。七.专题实习心得在本实验中，选用汉宁窗即可达到要求。窗函数法设计FIR数字滤波器是傅立叶变换的典型运用。通过上机实习，加深对数字滤波器的理解，并提高使用傅立叶变换分析解决实际问题的能力。这一阶段的实习让我体会
27、到了能读懂程序和自己编写程序的差距是很大的，理论与实际相结合是很重要的，只有理论知识是远远不够的，只有把所学的理论知识与实践相结合起来，经过大量的亲身实践后才能够得到更深的掌握，最终达到融会贯通的效果。题目五一、实验内容录制一段自己的语音信号，并对录制的信号进行采样；画出采样后语音信号的时域波形和频谱图；给定滤波器的性能指标，采用窗函数法或双线性变换设计滤波器，并画出滤波器的频率响应；然后用自己设计的滤波器对采集的语音信号进行滤波，画出滤波后信号的时域波形和频谱，并对滤波前后的信号进行对比，分析信号的变化；回放语音信号；最后，用MATLAB设计一信号处理系统界面。二、实验目的1、掌握用滤波器处
28、理语音信号的原理和方法。2、掌握用matlab设计界面。三、原理及算法概要调用函数function pushbutton1_Callback(hObject, eventdata, handles)实现一个信号处理系统界面。选择button1时，用双线性变换法设计滤波器来对信号进行处理。读取语音信号，对语音信号进行f=8000的频率进行采样，调用函数y1=fft(x1,2048)对所采集的点做2048点FFT变换。先设计butterworth模拟滤波器，再用双线性变换法实现模拟滤波器到数字滤波器的转换。最后调用函数f1=filter(bz,az,x2)对加了噪声的语音信号进行滤波，得到滤波后的
29、频谱图。4程序% - Executes on button press in pushbutton1.function pushbutton1_Callback(hObject, eventdata, handles)% hObject handle to pushbutton1 (see GCBO)% eventdata reserved - to be defined in a future version of MATLAB% handles structure with handles and user data (see GUIDATA)clc;fs=8000; %?8000x1=wavread(D:matlab2011workshua.wav);t=(0:length(x1)-1)/8000;y1=fft(x1,2048); %?2048?FFT?f=fs*(0:1023)/2048;plot(handles.axes1,t,x1); %?grid on;axis tight;title(handles.axes1,?);xlabel(handles.axes1,time(s)
