数字信号处理实验二.doc

实验二信号的分析与处理综合实验 一、实验目的目的：综合运用数字信号处理的理论知识进行信号的采样，重构，频谱分析和滤波器的设计，通过理论推导得出相应结论，再利用Matlab作为编程工具进行计算机实现，从而加深对所学知识的理解，建立概念二、基本要求1．掌握数字信号处理的基本概念、基本理论和基本方法；2．学会MATLAB的使用，掌握MATLAB的程序设计方法；3．掌握用MATLAB设计简单实验验证采样定理的方法；4．掌握在Windows环境下语音信号采集的方法；5．学会用MATLAB对信号进行频谱分析；6．掌握MATLAB设计FIR和IIR数字滤波器的方法；三、实验内容1． 利用简单正弦信号设计实验验证采样定理： （1） Matlab产生离散信号的方法，作图的方法，以及基本运算操作（2） 对连续正弦信号以不同的采样频率作采样（3） 对采样前后信号进行傅立叶变换，并画频谱图（4） 分析采样前后频谱的有变化，验证采样定理掌握画频谱图的方法，深刻理解采样频率，信号频率，采样点数，频率分辨率等概念2． 真实语音信号的采样重构：录制一段自己的语音信号，并对录制的信号进行采样；画出采样前后语音信号的时域波形和频谱图；对降采样后的信号进行插值重构，滤波，恢复原信号。

1） 语音信号的采集（2） 降采样的实现（改变了信号的采样率）（3） 以不同采样率采样后，语音信号的频谱分析（4） 采样前后声音的变化（5） 对降采样后的信号进行插值重构，滤波，恢复原信号3． 带噪声语音信号的频谱分析（1） 设计一频率已知的噪声信号，与实验2中原始语音信号相加，构造带噪声信号（2） 画出原始语音信号和加噪声后信号，以及它们的频谱图（3） 利用频谱图分析噪声信号和原语音信号的不同特性4． 对带噪声语音信号滤波去噪：给定滤波器性能指标，采样窗函数法或双线性变换设计滤波器，并画出滤波器的频率响应；然后用自己设计的滤波器对采样的语音信号进行滤波，画出滤波后信号的时域波形和频谱，并对滤波前后的信号进行对比，分析信号的变化；回放语音信号；（1） 分析带噪声信号频谱，找出噪声所在的频率段（2） 利用matlab中已有的滤波器滤波（3） 根据语音信号特点，自己设计滤波器滤波（4） 比较各种滤波器性能（至少四种），选择一种合适的滤波器将噪声信号滤除（5） 回放语音信号，比较滤波前后声音的变化四、实验原理参考《数字信号处理》教材《数字信号处理的MATLAB实现》万永革 编著五、主要实验仪器及材料微型计算机、Matlab。

六、实验步骤1．设计一简单正弦信号，通过改变采样率观察采样前后的信号变化例如：假设有一振幅为1，频率为10Hz，相位为0.3的模拟信号，即，用0.01s的采样间隔（采样频率为100Hz）来表示原始信号（注意：实际上模拟信号不能用离散值表示，此处为了在计算机上表示，用采样率非常高的离散信号表示模拟信号）分别以5Hz，10Hz（每秒采样10次，即采样间隔为0.1s），20Hz，40Hz，80Hz，200Hz对原始信号进行采样，画出采样前后的信号，并画出其频谱图，对比前后的变化，验证采样定理1）可以用t=0:1/fs:9/f; 取9个周期，通过改变采样率，自动改变采样点数2）也可以通过设置dt1(采样间隔)，已知采样点数n1，t1=n1*dt1，如图所示，采样率为40Hz时的原始信号，采样过程和采样后的信号时域图和频谱图，可见，当采样率大于原始信号频率的两倍时，采样前后信号频率基本不发生变化，信号不失真图2-1 采样的过程图2-2 原始信号和采样后信号的频谱2．语音信号的采集利用windows下的录音机或其他软件，录制一段自己的话音，时间控制在1秒左右然后在MATLAB软件平台下，利用函数wavread对语音信号进行采样，记住采样频率和采样点数。

通过wavread函数的使用，要求理解采样频率、采样位数等概念wavread函数调用格式：y=wavread(file)，读取file所规定的wav文件，返回采样值放在向量y中[y,fs,nbits]=wavread(file)，采样值放在向量y中，fs表示采样频率（Hz），nbits表示采样位数y=wavread(file,N)，读取前N点的采样值放在向量y中y=wavread(file,[N1,N2])，读取从N1点到N2点的采样值放在向量y中3．对真实语音信号的采样、重构1）降采样：利用windows下的录音机录制的音频采样率是固定的fs（=22050），可以选择以下函数实现对语音信号的降采样y=x((1:N:length(x))); %对原始信号每隔N个点取一位，即采样率变为原来的1/Ny=resample(yn,L,M); %采样率变为原来的L/M倍y=downsample(yn,N); %%采样率变为原来的1/N倍改变采样率为原来的1/2倍，1/4倍，1/20倍，1/50倍，1/100倍等，分别画出降采样前后的信号波形和频谱图，分析采样前后信号的变化图 2-3 采样率为原来的1/2时的信号波形频谱图图 2-4 采样率为原来的1/10时信号波形频谱图在MATLAB中，函数sound可以对声音进行回放。

其调用格式：sound(x,fs,bits)；通过调用此函数，感觉采样前后声音的变化2）重构原信号：降采样后，信号的采样率和采样点数同时变化如采样率变为原来的1/2，即对原始信号每隔一个点采样如果要恢复原始信号，即信号长度和采样率须变为原来同样大小所以，必须对降采样后信号进行插值重构具体过程参见附件1中例子对采样后的真实语音信号进行插值重构，滤波，恢复原始信号画出插值前后信号的波形以及频谱图，并将重构后信号与原始信号进行比较如，对采样率降为原来1/5的降采样后信号插值重构，结果如下图所示图2-5 采样率为原来的1/5时的波形频谱图图2-6 插值后的信号波形频谱图图2-7 低通滤波器的频率响应图图2-8 滤波后的波形频谱图调用sound函数感受插值后的声音，发现会有高频的噪声经过低通滤波器之后，高频噪声被滤出但是，因为之前的原始信号经过了降采样，所以插值后的效果一定不如原始声音也可用用wavwrite函数将经过处理的语音信号保存下来，调用格式为如wavwrite(y,fs,bits,'sound.wav')，其中，y为所要保存的语音信号，fs为其采样率，bits为采样位数，'sound.wav'是保存的语音信号的文件名。

4．语音信号的频谱分析要求：首先画出语音信号的时域波形；然后对语音信号进行频谱分析，在MATLAB中，可以利用函数fft对信号进行快速付立叶变换，得到信号的频谱特性；从而加深对频谱特性的理解[x,fs,bits] =wavread('e:\sound1.wav',[1024 5120]); %读取[1024 5120]段%语音信号数据；如果不指定所读取读取语音信号数据的长度，则读取整段语音信号X=fft(x，4096); % fftshift(fft(x))N=length(x);k=0:N-1;subplot(221);plot(k,x);title('原始信号波形');subplot(222);plot((-N/2:N/2-1)/N*fs,abs(X));title('原始信号频谱');5．对原始语音信号加噪声：设计频率已知的噪声信号或者用自然噪声信号加在原始语音信号上，构建带噪声信号如对原始信号加上频率为5000的正弦波噪声信号，程序如下：[x,fs,bits] =wavread('sound1.wav'); %读取原始语音信号X=abs(fftshift(fft(x))); %原始信号频谱c1=0.01*sin (2*pi*5000*k/fs); %构建频率为5000Hz的正弦波噪声信号yn=x+c1'; %构建带噪声信号Yn=abs(fftshift(fft (yn))); %求带噪声信号频谱subplot(321);plot(x);title('原始信号波形');subplot(322);plot((-N/2:N/2-1)/length(k)*fs,X);title('原始信号频谱');subplot(323);plot(yn);title('带噪声信号波形');subplot(324);plot((-N/2:N/2-1)/length(k)*fs,Yn);title('带噪声信号频谱');运行程序，结果如图2-9所示。

从频谱图可以看出，原始信号频谱段集中在0~5000Hz之间，主要频率集中在0Hz附近的低频部分加上噪声信号后，在5000Hz处有个幅值非常大的的高频成份，即以上所加的高频正弦波噪声信号图2-9 加噪声前后信号的波形频谱图6．设计数字滤波器和画出频率响应根据分析所得的原始信号的频谱和噪声信号频谱特点，给出有关滤波器的性能指标首先用窗函数法或者最优化法设计高通，低通，带通，带阻滤波器，在MATLAB中，可以利用函数fir1，firls设计FIR滤波器；然后在用双线性变换法或脉冲响应法设计上面几种滤波器，在MATLAB中，可以利用函数butte、cheby1和ellip设计IIR滤波器；最后，利用MATLAB中的函数freqz画出各滤波器的频率响应具体方法参加附件3种滤波器设计的步骤和实例如，根据以上语音信号的特点给出有关IIR滤波器的性能指标：1）低通滤波器性能指标，fp=4500，fc=6500，Rs=100，Rp=1fp：通带截至频率；fc：阻带截至频率；Rs：通带波纹；Rp：阻带波纹）2）带阻滤波器性能指标，fp1=4800 Hz，fp2=5200 Hz，fc1=4600 Hz，fc2=5400 Hz，Rs=30dB，Rp=1dB。

[fp1 fp2]：阻带截至频率；[fc1 fc2]：通带截至频率）程序如下：%椭圆带阻滤波器figure,wp=[4800 5200]*2/fs;ws=[4600 5400]*2/fs; %通带和阻带边界频率Rp=1;Rs=30;Nn=128; %通带波纹和阻带衰减以及绘制频率特性的数据点数[NN,Wn]=ellipord(wp,ws,Rp,Rs); %求取数字滤波器的最小阶数和归一化截至频率[b,a]=ellip(NN,Rp,Rs,Wn,'stop'); %设计滤波器freqz(b,a,512,fs); title('椭圆带阻滤波器');图 2-10 椭圆带阻滤波器的频率响应%双线性变换法设计的椭圆低通滤波器：figure,fp=4500;fc=6500;Rs=100;Rp=1;fs=22050;wc=2*fc/fs;wp=2*fp/fs;[n,wn]=ellipord(wp,wc,Ap,As);[b,a]=ellip(n,Rp,Rs,wn);freqz(b,a,512,fs);若设计FIR滤波器，要阻止5000Hz的高频信号通过，即设计如下滤波器特性，% 最优化设计法 设计带阻滤波器figure,n=100; %滤波器阶数f=[ 0 4600*2/fs 4800*2/fs 5400*2/fs 5600*2/fs 1 ]; %频率向量 a=[ 1 1 0。