生物医学信号处理大作业解析.docx
《生物医学信号处理大作业解析.docx》由会员分享，可在线阅读，更多相关《生物医学信号处理大作业解析.docx（17页珍藏版）》请在冰豆网上搜索。
生物医学信号处理大作业解析
生物医学信号处理大作业
题目：
基于matlab的语音信号处理
学生姓名：
学号：
专业：
学院：
精密仪器与光电子工程学院
作业要求
录制自己的一段语音：
“天津大学精密仪器与光电子工程学院,Collegeofprecisioninstrumentandopto-electronicsengineering,biomedicalengineering”，时间控制在15秒到30秒左右；利用wavread函数对自己的语音进行采样，记住采样频率。
（1）求原始语音信号的特征频带（比如谱峰位置）：
可以分别对一定时间间隔内，求功率谱（傅里叶变换结果取模的平方）并画出功率谱。
（2）根据语音信号频谱特点，设计FIR或IIR滤波器，分别画出滤波器幅频和相频特性曲线。
说明滤波器特性参数。
用设计的滤波器对信号滤波，画出滤波前后信号的频谱图。
用sound函数回放语音信号,说明利用高通/低通/带通滤波后的效果，不同特征频带被滤除后分别有什么效果。
（3）求出特征频段语音信号随时间变化的曲线（每隔一定时间求一次功率谱，连接成曲线，即短时FFT）。
（4）选做：
语谱图：
横轴为时间，纵轴为频率，灰度值大小表示功率谱值的大小。
（提示，可以采用spectrogram函数）
（5）选做：
分析自己的语音频谱特点，比如中英文发音的区别。
基于matlab的语音信号处理
摘要:
对录制的语音信号进行采样,分析其时域波形和频谱图。
给定数字滤波器的性能指标,采用窗函数法和双线性变换法设计数字滤波器,并对语音信号进行滤波,得到滤波前后的信号幅频响应。
通过对比高通、低通两种滤波处理结果,简单而有效地论证了两种数字滤波器在语音信号处理上的异同。
并进一步求出特征频段语音信号随时间变化的曲线，分析自身语音信号的特点。
关键词:
MATLAB数字滤波器语音信号
SpeechSignalProcessingbyDigitalFilterBasedonMATLAB
Abstract：
Time-domainwaveformandfrequencyspectrumoftherecordedspeechsignalsareanalyzedbysampling.Theperformanceindexesofdigitalfiltersaregiven.Twomethodsofwindowfunctionandbilineartransformationareusedtodesignthedigitalfilters.Thespeechsignalisfilteredbythefilters,andthenmagnitude-frequencyresponsesofthesignalbeforeandafterfilteringarereceived.Theadvantagesoftwodigitalfilters（filterlowpassandfilterhighpass）inspeechsignalprocessingaredemonstratedbycomparingdifferentmethodsforfilteringsimplyandeffectively.Formore,weareabletofigureoutthetimecurvesofcharacteristicbandsofspeechsignalandthenanalysesthecharacterofourownspeechsignals.
Keywords:
MATLAB,digitalfilter,speechsignal
为了进一步观察和确定语音信号的频谱特征，首先分别画出每秒的频带特征。
具体代码如下所示：
clear
[yfull,fs,nbits]=wavread（'D:
\meng'）;
y=yfull（:
1）;
win=blackman（fs）;
yfft=zeros（15,fs）;
ypower=zeros（15,fs）;
forii=0:
14
yfft（ii+1,1:
fs）=abs（fft（y（1+ii*fs:
（ii+1）*fs）.*win,fs））;
ypower（ii+1,1:
fs）=yfft（ii+1,1:
fs）.^2;
figure（ii+1）;
f=1:
3000;5
plot（f,ypower（ii+1,1:
3000））;
end
由此得到十五秒中每一秒的强度分布（注：
横坐标为频率，纵坐标为强度）
第1秒第2秒
第3秒第4秒
第5秒第6秒
第7秒第8秒
第9秒第10秒
第11秒第12秒
第13秒第14秒
第15秒
FIR滤波器的设计：
FIR滤波器的设计方法有窗函数法、频率取样法、最佳等波纹法等。
我选用了窗函数法设计FIR滤波器，窗函数设计法是一种通过截短,用有限长脉冲响应序列来逼近无限长非因果序列。
根据所要设计的数字滤波器的性能指标要求，采用窗函数设计FIR滤波器的步骤如下:
（1）根据最小阻带衰减As要求,选择窗函数类型。
（2）根据截止频率、滤波器的类型、最大通带衰减Ap、最小阻带衰减As以及采样频率,利用MATLAB阶次估计函数获得窗函数的参数。
（3）利用MATLAB设计加窗有限冲激响应滤波器函数fir1设计FIR数字滤波器。
（4）通过计算数字滤波器的频率响应函数freqz绘出FIR数字滤波器的频率响应曲线。
常用窗函数的性能指标如表1所示。
其中,凯泽窗[2]是一种适应性较强的窗,窗函数中的B是一个可自由选择的参数,它可以同时调整主瓣宽度与旁瓣电平。
因而,改变B值就可对主瓣宽度与旁瓣衰减进行选择。
由于滤波器指标中没有对主瓣宽度做要求,所以B值的确定由最小阻带衰减As决定。
根据以上十五个图像可以看出所有图像在200Hz至1200Hz处均有谱峰，因此选择这两个参数进行归一处理。
并利用凯瑟窗（学号29对应于凯瑟窗的设计）进行加窗处理设计滤波器，并进行滤波。
具体程序如下所示:
%----------归一化频率------------
w1=200/fs*pi;
w2=1200/fs*pi;
%----------凯瑟窗设计滤波器----------
N=40;
beta=6;
w=kaiser（N,beta）;
h=fir1（N-1,[w1w2],w）;
[h1,w3]=freqz（h,1）;
figure（16）;
plot（w3/pi,20*log10（abs（h1）））;
gridon;
axis（[0,1,-80,10]）;
xlabel（'归一化频率/pi'）;
ylabel（'幅度/dB'）;
%------------滤波---------------
yfilt=filter（h,1,y）;
figure（17）;
f=0:
fs/778239:
fs;
plot（f,yfilt）;
得到如下所示的波形：
下面分别用低通，高通（带通已在凯瑟窗的设计中实现）进行滤波，画出滤波前后信号的频谱图并进行相应的比较和分析。
具体程序如下所示：
%---------------低通滤波器的设计及滤波---------------
Fl=400;Fh=800;
wp=（Fl/fs）*2*pi;
ws=（Fh/fs）*2*pi;
OmegaP=2*fs*tan（wp/2）;
OmegaS=2*fs*tan（ws/2）;
[n,Wn]=buttord（OmegaP,OmegaS,1,25,'s'）;
[b,a]=butter（n,Wn,'s'）;
%freqs（b,a）%设计模拟滤波器
[bz,az]=bilinear（b,a,fs）;%映射为数字滤波器参数
[hl,wl]=freqz（bz,az）;
fl=wl/pi*fs./2;
figure（18）;
plot（fl,20*log10（abs（hl）））;
title（'低通滤波器'）
grid;
yfiltlowpass=filter（bz,az,y）;
figure（19）;
plot（yfiltlowpass）;
title（'通过低通滤波器后的波形'）;
%---------------高通滤波器的设计及滤波----------------
wc=2*fs*tan（2*pi*600/（2*fs））;
wr=2*fs*tan（2*pi*200/（2*fs））;
[Nh,wnh]=cheb1ord（wc,wr,0.8,20,'s'）;%估计滤波器阶数
[Bh,Ah]=cheby1（Nh,0.8,wnh,'high','s'）;
[numh,denh]=bilinear（Bh,Ah,1000）;
[hh,wh]=freqz（numh,denh）;
fh=wh/pi*500;
figure（20）;
plot（fh,20*log10（abs（hh）））;
grid;
xlabel（'频率/Hz'）;
ylabel（'幅度/dB'）;
title（'高通滤波器'）;
yfilthighpass=filter（numh,denh,y）;
figure（21）;
plot（yfilthighpass）;
title（'通过高通滤波器后的波形'）;
通过试听可以发现，利用凯瑟窗设计的滤波器的效果最为理想，既达到了声音听起来更加清晰，又有效地保持了声音的特性，所有文字都可以清晰的听到。
通过低通滤波器后，声音变得比较浑厚，听起来有些发闷，但是还是比较接近原声的，我认为原因是滤掉了高频成分，所以音调听起来比较低。
而通过使用高通滤波器，声音明显变小变细，而且变化是极为明显的，无论怎样改变参数，声音都比较小，高通滤波后只能听到少许的杂音。
接下来，求出特征频段语音信号随时间变化的曲线（每隔一定时间求一次功率谱，连接成曲线，即短时FFT）。
程序集图像如下图所示：
%----------短时傅里叶变换-------------
yfiltfft=zeros（166,44100）;
yfiltfftwidth300_350=zeros（1,166）;
yfiltfftwidth580_700=zeros（1,166）;
yfiltfftwidth1000_1050=zeros（1,166）;
forjj=0:
165
yfiltfft（jj+1,1:
44100）=abs（fft（yfilt（jj*4410+1:
44100+jj*4410）,44100））;
yfiltfftwidth300_350（1,jj+1）=mean（yfiltfft（jj+1,300:
350））^2;
yfiltfftwidth580_700（1,jj+1）=mean（yfiltfft（jj+1,580:
700））^2;
yfiltfftwidth1000_1050（1,jj+1）=mean（yfiltfft（jj+1,1000:
1050））^2;
end
figure（22）;
plot（yfiltfftwidth300_350）;
figure（23）;
plot（yfiltfftwidth580_700）;
figure（24）;
plot（yfiltfftwidth1000_1050）;
对于选做题目，做出了语谱图
%----------语谱图--------------
spectrogram（yfilt,hamming（fs）,50,512,fs）
对于语谱图的相关知识：
语谱图就是语音频谱图，一般是通过处理接收的时域信号得到频谱图，因此只要有足够时间长度的时域信号就可。
（时间长度为保证频率分辨率）专业点讲，那是频谱分析视图，如果针对语音数据的话，叫语谱图。
语谱图的横坐标是时间，纵坐标是频率，坐标点值为语音数据能量。
由于是采用二维平面表达三维信息，所以能量值的大小是通过颜色来表示的，颜色深，表示该点的语音能量越强。
之所以采用二维平面来表示三维信息是由于一下原因：
在数字技术发展以前，人们可视化研究语音数据的方法是把数据通过频率滤波器，然后各个频率的数据驱动相应的类似针式打印的设备按频率高低顺序记录在一卷纸上，信号的强弱由记录在纸上的灰度来表示。
记录纸按照一定的速度旋转，即相当于在不同的时间里记录下语音数据。
语谱图还分为窄带语谱图和宽带语谱图。
从语谱图上，我们也可以看出很多的信息：
我们可以观察语音不同频段的信号强度随时间的变化情况。
从以上的图中可以看到明显的、一条条横方向的条纹，我们称为“声纹”，有很多应用。
条纹的地方实际是颜色深的点聚集的地方，随时间延续，就延长成条纹，也就是表示语音中频率值为该点横坐标值的能量较强，在整个语音中所占比重大，那么相应影响人感知的效果要强烈得多。
而一般语音中数据是周期性的，所以，能量强点的频率分布是频率周期的，即存在300Hz强点，则一般在n*300Hz点也会出现强点，所以我们看到的语谱图都是条纹状的。
尽管客观人发声器官的音域是有限度的，即一般人发声最高频率4000Hz，但是，由于我们数字分析频率时，采用的是算法实现的，一般是FFT，所以其结果是由采样率决定的，即尽管是上限为4000Hz的语音数据，如果采用16Khz的采样率来分析，则仍然可以在4000Hz以上的频段发现有数据分布，则可以认为是算法误差，非客观事实。
由我的语谱图可以很清楚的看到我说话的特点：
开始的第一个音强度在整句话中大致是相同的，后面的语句随着不同的语言而有所差异，中间部分的语气是比较弱的。
进一步分析自身的语音频谱特点：
总的来讲，我的中文发音和英文发音差距不是很大，只是英文的语调起伏相对较大（表现为纵坐标强度）。
同时，我说话的频率基本集中在200Hz到1200Hz的区域。
接下来，可以再找一个声音信号（男生）分析男女生说话的差异、分析男声和女声的差别，实验代码如下：
y1=wavread（'nansheng.wav'）;
Y1=y1（:
1）;
Fs=44100;%采样频率%
yt1=fft（y1）;%傅里叶变换%
df=Fs/length（yt1）;
Fx=df*（0:
length（yt1）-1）;%将横轴变为频率轴%
subplot（1,1,1）;%subplot将图像画在一张图上%
plot（y1）;title（'语音时域波形'）;%画语音信号的时域波形%
subplot（1,1,1）;
plot（Fx,abs（yt1））;axis（[0500002000]）;title（'频谱图'）;xlabel（'frequency/Hz'）;
%频谱图%
y1=wavread（‘nvsheng.wav'）;
Y1=y1（:
1）;
Fs=44100;%采样频率%
yt1=fft（y1）;%傅里叶变换%
df=Fs/length（yt1）;
Fx=df*（0:
length（yt1）-1）;%将横轴变为频率轴%
subplot（1,1,1）;%subplot将图像画在一张图上%
plot（y1）;title（'语音时域波形'）;%画语音信号的时域波形%
subplot（1,1,1）;
plot（Fx,abs（yt1））;axis（[0500002000]）;title（'频谱图'）;xlabel（'frequency/Hz'）;
%频谱图%
以下为所得图形：
由此可见，男生的声音信号强度比女生的大，而且大的十分明显。
在频率方面，男生的音域要比女生的广。
从男生、女生时域波形看出：
男生说话时的杂音比较大，女生说话的声音更加清晰。
同时，通过基于MATLAB和短时频域分析，并结合实际操作的分析，能够得出:
一．长窗具有较高的频率分辨率，但具有较低的时间分辨率。
从一个周期到另一个周期，共振峰是要发生变化的，这一点即使从语音波形上也能够看出来。
然而，如果采用较长的窗，这种变化就模糊了，因为长窗起到了时间上的平均作用。
二．短窗的频率分辨率低，但具有较高的时间分辨率。
采用短窗时，能够从短时频谱中提取出共振峰从一个周期到另一个周期所发生的变化。
当然，激励源的谐波结构也从短时频谱上消失了。
三．在对语音信号进行短时傅里叶分析时，窗长需要折衷考虑。
一方面，短窗具有较好的时间分辨率因而能够提取出语音信号中的短时变化;但另一方面，损失了频率分辨率。
还应当注意到，语音信号的基音周期是有一个覆盖范围的（例如，女性和儿童的基音周期就比成年男性短得多），因此，窗宽的选择应当考虑到这个因素。
四．我所采用的凯瑟窗更接近于带通特性，对于处理本人的特性是比较合适的。
升级会员 