第三个数字信号设计.docx
《第三个数字信号设计.docx》由会员分享，可在线阅读，更多相关《第三个数字信号设计.docx（20页珍藏版）》请在冰豆网上搜索。
第三个数字信号设计
广西工学院
数字信号处理课程设计
说明书
设计题目语音信号的处理与滤波
系别计算机工程系
专业班级通信061班
学生姓名李凤梅
学号200600402033
指导教师周坚和
日期2009年1月2日
目录
一．设计目的……………………………………………3
二．设计内容……………………………………………3
三．设计原理……………………………………………3
四．具体实现……………………………………………5
1.录制声音……………………………………………5
2.设计FIR滤波器……………………………………6
3.将声音信号送入滤波器滤波………………………10
4.回放语音信号………………………………………16
5.男女语音信号的频谱分析…………………………16
6.噪声叠加和滤除……………………………………19
五．设计一个系统界面…………………………………23
六．结果分析……………………………………………25
七．课程设计总结………………………………………26
八．参考文献……………………………………………27
一．设计目的
综合运用本课程的理论知识进行频谱分析以及滤波器设计，通过理论推导得出相应结论，并利用MATLAB作为工具进行实现，从而复习巩固课堂所学的理论知识，提高对所学知识的综合应用能力，并从实践上初步实现对数字信号的处理。
二．设计内容
录制一段个人自己的语音信号，并对录制的信号进行采样；画出采样后语音信号的时域波形和频谱图；给定滤波器的性能指标，采用窗函数法和双线性变换法设计滤波器，并画出滤波器的频率响应；然后用自己设计的滤波器对采集的信号进行滤波，画出滤波后信号的时域波形和频谱，并对滤波前后的信号进行对比，分析信号的变化；回放语音信号；换一个与你性别相异的人录制同样一段语音内容，分析两段内容相同的语音信号频谱之间有什么特点；再录制一段同样长时间的背景噪声叠加到你的语音信号中，分析叠加前后信号频谱的变化，设计一个合适的滤波器，能够把该噪声滤除；最后，设计一个信号处理系统界面。
三．设计原理
1.在Matlab软件平台下，利用函数wavrecord（），wavwrite（）,wavread（）,wavplay（）对语音信号进行录制，存储，读取，回放。
2.用y0=fft（x）对采集的信号做快速傅立叶变换，并用[h1,w]=freqz（h）进行DTFT变换。
3.掌握FIRDF线性相位的概念，即线性相位对
及零点的约束，了解四种FIRDF的频响特点。
4.在Matlab中，FIR滤波器利用函数fftfilt对信号进行滤波。
5.抽样定理
连续信号经理想抽样后时域、频域发生的变化（理想抽样信号与连续信号频谱之间的关系）
理想抽样信号能否代表原始信号、如何不失真地还原信号即由离散信号恢复连续信号的条件（抽样定理）
理想采样过程描述：
时域描述：
频域描述：
利用傅氏变换的性质，时域相乘频域卷积，若
则有
与
的关系：
理想抽样信号的频谱是连续信号频谱的周期延拓,重复周期为s（采样角频率）。
如果：
即连续信号是带限的，且信号最高频率不超过抽样频率的二分之一，则可不失真恢复。
奈奎斯特采样定理：
要使实信号采样后能够不失真还原，采样频率必须大于信号最高频率的两倍：
四．具体实现
1．录制一段声音
在MATLAB中用wavrecord、wavread、wavplay对声音进行录制、读取、回放。
程序如下：
fs=8000;%抽样频率
x=wavrecord（3*fs,fs）;%按回车键，录入声音
wavwrite（x,fs,'lfm.wav'）;%保存录音
[x,fs,bits]=wavread（'lfm.wav'）;%播放原始信号
wavplay（x,fs）;%播放原始信号
N=length（x）;%返回采样点数
df=fs/N;%采样间隔
n1=1:
N/2;
f=[（n1-1）*（2*pi/N）]/pi;%频带宽度
figure
（2）;
subplot（2,1,1）;
plot（x）;%录制信号的时域波形
title（'原始信号的时域波形'）;%加标题
ylabel（'幅值/A'）;%显示纵坐标的表示意义
grid;%加网格
y0=fft（x）;%快速傅立叶变换
figure
（2）;
subplot（2,1,2）;
plot（f,abs（y0（n1）））;%原始信号的频谱图
title（'原始信号的频谱图'）;%加标题
xlabel（'频率w/pi'）;%显示横坐标表示的意义
ylabel（'幅值'）;%显示纵坐标表示的意义
title（'原始信号的频谱图'）;%加标题
grid;%加网格
2．设计FIR滤波器
（1）．低通滤波器
Fs=18000;
fb=1000;fc=1200;As=100;Ap=1;
wb=2*pi*fb/Fs;wc=2*pi*fc/Fs;
deltaw=wc-wb;
ws=（wc+wb）/2;
N1=ceil（2*pi*6.4/deltaw）;
N0=N1+mod（N1+1,2）;
hdwindow=ideallp（ws,N0）;
wdwindow=kaiser（N0,10.056）;
hr=hdwindow.*wdwindow';
[h1,w]=freqz（hr,1）;
subplot（2,2,1）;
plot（w/pi,20*log10（abs（h1）））;
title（'凯泽窗低通滤波器频率响应'）
xlabel（'\omega/\pi'）;ylabel（'绝对振幅'）;
grid;
（2）高通滤波器
Fs=12000;%自定义的抽样频率
fc=4800;fb=5000;As=100;Ap=1;%给定的滤波器性能指标
wb=2*fb*pi/Fs;wc=2*pi*fc/Fs;%转化为数字频率
deltaw=wb-wc;%通带截止频率与阻带截止频率的距离
ws=（wb+wc）/2;%过渡带宽
N2=ceil（2*pi*6.4/deltaw）;%求抽样长度
N0=N2+mod（N2+1,2）;%确保抽样长度为奇数
hdwindow=ideallp（pi,N0）-ideallp（ws,N0）;%理想脉冲响应
wdwindow=kaiser（N0,10.056）;%凯泽窗
hr=wdwindow.*hdwindow';%凯泽窗单位抽样响应
[h2,w]=freqz（hr,1）;%对hr进行DTFT变换
subplot（2,2,2）;
plot（w/pi,20*log10（abs（h2）））;%画凯泽窗高通滤波器频率响应图
title（'凯泽窗高通滤波器频率响应'）%标题
xlabel（'\omega/\pi'）;%X轴名称
ylabel（'绝对振幅'）;%Y轴名称
grid;
（3）．带通滤波器
fs=24000;%抽样频率
fp1=1200;fp2=3000;fsc1=1000;fsc2=3200;
Ap=1;As=100;%给定的滤波器性能指标
Wp1=2*pi*fp1/fs;Wp2=2*pi*fp2/fs;%转化为数字频率
Wsc1=2*pi*fsc1/fs;Wsc2=2*pi*fsc2/fs;%转化为数字频率
Wc1=（Wp1+Wsc1）/2;Wc2=（Wp2+Wsc2）/2;%过渡带宽
deltaw=min（（Wp1-Wsc1）,（Wsc2-Wp2））;%最小通带截止频率与阻带截止频率的距离
N3=ceil（2*6.4*pi/deltaw）;%求抽样长度
N0=N3+mod（N3+1,2）;%确保抽样长度为奇数
hdwindow=ideallp（Wc2,N0）-ideallp（Wc1,N0）;%理想脉冲响应
wdwindow=kaiser（N0,10.056）;%凯泽窗
hr=wdwindow.*hdwindow';%凯泽窗单位抽样响应
[h3,w]=freqz（hr,1）;%对hr进行DTFT变换
subplot（2,2,3）;
plot（w/pi,20*log10（abs（h3）））;%画凯泽窗带通滤波器频率响应图
title（'凯泽窗带通滤波器频率响应'）%标题
xlabel（'\omega/\pi'）;%X轴名称
ylabel（'绝对振幅'）;%Y轴名称
grid;%显示网格
图形如下：
分析：
根据老师给出各滤波器的性能指标：
（1）低通滤波器性能指标：
fb＝1000Hz，fc＝1200Hz，As＝100dB，Ap＝1dB。
（2）高通滤波器性能指标：
fc＝4800Hz，fb＝5000HzAs＝100dB，Ap＝1dB。
（3）带通滤波器性能指标：
fp1＝1200Hz，fp2＝3000Hz，fsc1＝1000Hz，fsc2＝3200Hz，As＝100dB，Ap＝1dB。
得出的FIR滤波器图形
3．将声音信号送入滤波器滤波
[x,fs,bits]=wavread（'lfm.wav'）;%播放原始信号
wavplay（x,fs）;%播放原始信号
N=length（x）;%返回采样点数
df=fs/N;%采样间隔
n1=1:
N/2;
f=[（n1-1）*（2*pi/N）]/pi;%频带宽度
figure（4）;
subplot（4,2,1）;
plot（x）;%录制信号的时域波形
title（'原始信号的时域波形'）;%加标题
ylabel（'幅值/A'）;%显示纵坐标的表示意义
grid;%加网格
y0=fft（x）;%快速傅立叶变换
figure（4）;
subplot（4,2,3）;
plot（f,abs（y0（n1）））;%原始信号的频谱图
title（'原始信号的频谱图'）;%加标题
xlabel（'频率w/pi'）;%显示横坐标表示的意义
ylabel（'幅值'）;%显示纵坐标表示的意义
title（'原始信号的频谱图'）;%加标题
grid;%加网格
（1）．低通滤波器滤波
fs=8000;
beta=10.056;
wc=2*pi*1000/fs;
ws=2*pi*1200/fs;
width=ws-wc;
wn=（ws+wc）/2;
n=ceil（12.8*pi/width）;
h=fir1（n,wn/pi,'band',kaiser（n+1,beta））;
[h1,w]=freqz（h）;
ys=fftfilt（h,x）;%信号送入滤波器滤波，ys为输出
fftwave=fft（ys）;%将滤波后的语音信号进行快速傅立叶变换
figure（4）;
subplot（4,2,2）;%在四行两列的第二个窗口显示图形
plot（ys）;%信号的时域波形
title（'低通滤波后信号的时域波形'）;%加标题
xlabel（'频率w/pi'）;
ylabel（'幅值/A'）;%显示标表示的意义
grid;%网格
figure（4）;
subplot（4,2,4）;%在四行两列的第四个窗口显示图形
plot（f,abs（fftwave（n1）））;%绘制模值
xlabel（'频率w/pi'）;
ylabel（'幅值/A'）;%显示标表示的意义
title（'低通滤波器滤波后信号的频谱图'）;%标题
grid;%加网格
wavplay（ys,8000）;%播放滤波后信号
（2）．高通滤波器滤波
fs=8000;
beta=10.056;
ws=2*5000/fs;
wc=2*4800/fs;
width=ws-wc;
wn=（ws+wc）/2;
n=ceil（12.8*pi/width）;
h=fir1（n,wn/pi,'high',kaiser（n+2,beta））;
[h1,w]=freqz（h）;
ys=fftfilt（h,x）;%将信号送入高通滤波器滤波
figure（4）;
subplot（4,2,5）;%在四行两列的第五个窗口显示图形
plot（ys）;%信号的时域波形
xlabel（'频率w/pi'）;
ylabel（'幅值/A'）;%显示标表示的意义
title（'高通滤波后信号的时域波形'）;%标题
ylabel（'幅值/A'）;%显示纵坐标的表示意义
grid;%网格
fftwave=fft（ys）;%将滤波后的语音信号进行快速傅立叶变换
figure（4）;
subplot（4,2,7）;%在四行两列的第七个窗口显示图形
plot（f,abs（fftwave（n1）））;%绘制模值
axis（[01050]）;
xlabel（'频率w/pi'）;
ylabel（'幅值/A'）;%显示标表示的意义
title（'高通滤波器滤波后信号的频谱图'）;%标题
grid;%加网格
wavplay（ys,8000）;%播放滤波后信号
（3）．带通滤波器
fs=8000;
beta=10.056;
wc1=2*pi*1000/fs;
wc2=2*pi*3200/fs;
ws1=2*pi*1200/fs;
ws2=2*pi*3000/fs;
width=ws1-wc1;
wn1=（ws1+wc1）/2;
wn2=（ws2+wc2）/2;
wn=[wn1wn2];
n=ceil（12.8/width*pi）;
h=fir1（n,wn/pi,'band',kaiser（n+1,beta））;
[h1,w]=freqz（h）;
ys1=fftfilt（h,x）;%将信号送入高通滤波器滤波
figure（4）;
subplot（4,2,6）;%在四行两列的第六个窗口显示图形
plot（ys1）;%绘制后信号的时域的图形
title（'带通滤波后信号的时域波形'）;%加标题
xlabel（'频率w/pi'）;
ylabel（'幅值/A'）;%显示纵坐标表示的意义
grid;%网格
fftwave=fft（ys1）;%对滤波后的信号进行快速傅立叶变换
figure（4）;
subplot（4,2,8）;%在四行两列的第八个窗口显示图形
plot（f,abs（fftwave（n1）））;%绘制模值
axis（[01050]）;
xlabel（'频率w/pi'）;
ylabel（'幅值/A'）;%显示标表示的意义
title（'带通滤波器滤波后信号的频谱图'）;%加标题
grid;%网格
wavplay（ys1,8000）;%播放滤波后信号
图形如下：
分析：
三个滤波器滤波后的声音与原来的声音都发生了变化。
其中低通的滤波后与原来声音没有很大的变化，其它两个都又明显的变化
4.回放语音信号
在滤波后加wavpla（）函数对语音进行了回放
5.男女语音信号的频谱分析
（1）．录制一段异性的声音进行频谱分析
fs=8000;%抽样频率
y=wavrecord（3*fs,fs）;%按回车键，录入声音
wavwrite（y,fs,'nansheng.wav'）;%保存到指定的位置
（2）．分析男女声音的频谱
[x,fs,bits]=wavread（'lfm.wav'）;%播放原始信号,解决落后半拍
wavplay（x,fs）;%播放原始信号
N=length（x）;%返回采样点数
df=fs/N;%采样间隔
n1=1:
N/2;
f=[（n1-1）*（2*pi/N）]/pi;%频带宽度
figure（4）;
subplot（2,2,1）;
plot（x）;%录制信号的时域波形
title（'原始女生信号的时域波形'）;%加标题
ylabel（'幅值/A'）;%显示纵坐标的表示意义
grid;%加网格
y0=fft（x）;%快速傅立叶变换
figure（4）;
subplot（2,2,2）;
plot（f,abs（y0（n1）））;%原始信号的频谱图
title（'原始女生信号的频谱图'）;%加标题
xlabel（'频率w/pi'）;%显示横坐标表示的意义
ylabel（'幅值'）;%显示纵坐标表示的意义
grid;%加网格
[y,fs,bits]=wavread（'nansheng.wav'）;%对语音信号进行采样
wavplay（y,fs）;%播放原始信号
N=length（y）;%返回采样点数
df=fs/N;%采样间隔
n1=1:
N/2;
f=[（n1-1）*（2*pi/N）]/pi;%频带宽度
subplot（2,2,3）;
plot（y）;%录制信号的时域波形
title（'原始男生信号的时域波形'）;%加标题
ylabel（'幅值/A'）;%显示纵坐标的表示意义
grid;%加网格
y0=fft（y）;%快速傅立叶变换
subplot（2,2,4）;%在四行两列的第三个窗口显示图形
plot（f,abs（y0（n1）））;%原始信号的频谱图
title（'原始男生信号的频谱图'）;%加标题
xlabel（'频率w/pi'）;%显示横坐标表示的意义
ylabel（'幅值'）;%显示纵坐标表示的意义
grid;%加网格
（3）男女声音的频谱图
分析：
就时域图看，男生的时域图中振幅比女生的高，对于频谱图女生的高频成分比较多
6.噪声的叠加和滤除
（1）录制一段背景噪声
fs=8000;%抽样频率
x1=wavrecord（3*fs,fs）;%按回车键，录入声音
wavwrite（x1,fs,’zaosheng.wav'）;%保存到指定的位置
（2）．对噪声进行频谱的分析
[x1,fs,bits]=wavread（'zaosheng.wav'）;%对语音信号进行采样
wavplay（x1,fs）;%播放噪声信号
N=length（x1）;%返回采样点数
df=fs/N;%采样间隔
n1=1:
N/2;
f=[（n1-1）*（2*pi/N）]/pi;%频带宽度
figure（5）;
subplot（3,2,1）;
plot（x1）;%信号的时域波形
title（'噪声信号的时域波形'）;
grid;
ylabel（'幅值/A'）;
y0=fft（x1）;%快速傅立叶变换
figure（5）;
subplot（3,2,2）;
plot（f,abs（y0（n1）））;%噪声信号的频谱图
ylabel（'幅值'）;
title（'噪声信号的频谱图'）;
（3）．将噪声和自己录制的声音叠加并分析叠加后的频谱
噪声的录制与叠加
fs=8000;
[x,fs,bits]=wavread（'lfm.wav'）;%对录入信号进行采样
[x1,fs,bits]=wavread（'zaosheng.wav'）;%对噪声信号进行采样
yy=x+x1;%将两个声音叠加
叠加信号的频谱分析：
wavplay（yy,fs）;%播放叠加后信号
N=length（yy）;%返回采样点数
df=fs/N;%采样间隔
n1=1:
N/2;
f=[（n1-1）*（2*pi/N）]/pi;%频带宽度
figure（5）;
subplot（3,2,3）;
plot（yy,'LineWidth',2）;%信号的时域波形
title（'叠加信号的时域波形'）;
xlabel（'时间/t'）;
ylabel（'幅值/A'）;
grid;
y0=fft（yy）;%快速傅立叶变换
figure（5）;
subplot（3,2,4）;
plot（f,abs（y0（n1）））;%叠加信号的频谱图
title（'叠加信号的频谱图'）;
xlabel（'频率w/pi'）;
ylabel（'幅值/db'）;
grid;
（4）．设计一个合适的滤波器将噪声滤除
fs=18000;%采样频率
Wp=2*1000/fs;%通带截至频率
Ws=2*2000/fs;%阻带截至频率
Rp=1;%最大衰减
Rs=100;%最小衰减
[N,Wn]=buttord（Wp,Ws,Rp,Rs）;%buttord函数（n为阶数，Wn为截至频率）
[num,den]=butter（N,Wn）;%butter函数（num为分子系数den为分母系数）
[h,w]=freqz（num,den）;%DTFT变换
ys=filter（num,den,yy）;%信号送入滤波器滤波，ys为输出
fftwave=fft（ys）;%将滤波后的语音信号进行快速傅立叶变换
figure（5）;
subplot（3,2,5）;
plot（ys）;%信号的时域波形
title（'低通滤波后信号的时域波形'）;%加标题
ylabel（'幅值/A'）;%显示标表示的意义
grid;%网格
figure（5）;
subplot（3,2,6）;
plot（f,abs（fftwave（n1）））;%绘制模值
title（'低通滤波器滤波后信号的频谱图'）;%标题
xlabel（'频率w/pi'）;
ylabel（'幅值/A'）;%显示标表示的意义
grid;%加网格
wavplay（ys,8000）;%播放滤波后信号
grid;
图形如下：
五、设计一个系统界面
1.GUIDE的控制程序：
%主体界面
h=figure（'Visible','off'）;
set（h,'color',[0.60.21],'position',[100100350400],...
'name','语音信号滤波分析',...
'NumberTitle','off','pointer','cross',...
'Units','pixels','Resize','on'）;
set（h,'Visible','on'）;
%控件
A=uicontrol（'parent',h,'style','pushbutton','Units','points','position',[7025010025],...
'String','原始信号','callback','OriginalSignal'）;
B=uicontrol（'parent',h,'style','pushbutton','Units','points','position',[7021010025],...
'String','滤波器的幅频特性','callback','Filter1'）;
C=uicontrol（'parent',h,'style','pushbutton','Units','points','position',[7017010025],...
'String','信号通过滤波器滤波','callback','lvbo'）;
E=uicontrol（'parent',h,'style','pushbutton','Units','points','position',[7013010025],...
'String','男女声音信号频谱比较','callback','n
升级会员 