dsp 语音信号处理代码+报告.zip

云炬网络

于 2024-06-10 21:05:27 发布

阅读量869

点赞数 7

本文链接：https://blog.csdn.net/m0_66999594/article/details/139582188

版权

《数字信号处理》

课程大作业

班级：

学号：

姓名：

设计目的
学会MATLAB的使用，掌握MATLAB程序设计方法；
掌握使用MATLAB对语言信号采集的方法；
掌握数字信号处理的基本概念、基本理论和基本方法；
掌握MATLAB设计IIR数字滤波器的设计方法；
学会用MATLAB对信号进行分析和处理。
设计要求
使用手机录制一段自我介绍音频，命名为“.wav”。
利用Matlab软件对“.wav”音频信号进行数字信号采样，分别对采样后的信号进行时/频域分析，并提供仿真图（包括语音信号的时域波形和频谱图）和分析说明；
设计数字滤波器，对“学号.wav”音频信号分别进行低通滤波（截止频率4000Hz，阻带衰减20dB，过渡带宽0.1）与高通滤波（截止频率3000Hz，阻带衰减40dB，过渡带宽0.1），提供时/频域仿真图（包括所设计滤波器的损耗函数曲线以及滤波后的频谱图）和分析说明；
给“.wav”添加一段幅度大小为0.1，频率f为9000Hz的正弦噪声，合成一段有噪声的音频，统一命名为“noisy_signal.wav”。对添加噪声后的信号进行时/频域分析，并提供仿真图（包括语音信号的时域波形和频谱图）和分析说明；
设计合理的数字滤波器，滤去音频信号中的噪声，提供时/频域仿真图（包括所设计滤波器的损耗函数曲线以及滤波后的频谱图）和分析说明，并将数字滤波后的数字信号转换成wav格式音频文件，统一命名为denoised.wav。
设计步骤、实现方法、结果分析
时域采样
设计流程
利用‘audioread’函数对录制的‘xx.wav’音
频文件进行采样，调用格式为
[voice,fs]=audioread('xx.wav');%声音读取
采样后的采样值保存为‘voice’，‘fs’为采样频率，得到‘fs=48000’;voice中包含两个声道，将两个声道的音频分别提取出来，保存为‘v1’‘v2’。对‘voice’‘v1’‘v2’进行快速傅里叶变换，并绘制出时域波形和频谱。
实现方式
%% 声音读取
[voice,fs]=audioread('xx.wav');%声音读取
v1=voice(:,1); %抽取第一声道
v2=voice(:,2); %抽取第二声道
%sound(voice,fs); %立体声音播放
%sound(v1,fs); %播放第一声道
%sound(v2,fs); %播放第二声道
%% 傅里叶变换
fft1=fft(voice,N)/N; %voice%傅里叶变换
mag1=abs(fftshift(fft1)); %取模
f1=(-N/2:1:N/2-1)'*fs/N; %双边变换
fft2=fft(v1,N)/N; %v1%傅里叶变换
mag2=abs(fftshift(fft2)); %取模
f2=(-N/2:1:N/2-1)'*fs/N; %双边变换
fft3=fft(v2,N)/N; %v2%傅里叶变换
mag3=abs(fftshift(fft3)); %取模
f3=(-N/2:1:N/2-1)'*fs/N; %双边变换
%%
figure(1); %原始语音信号时域频域波形
subplot(311);plot(t,voice);
title('立体音频信号时域波形');
xlabel('t/s');
ylabel('幅度');
subplot(312);plot(t,v1);
title('第一声道音频信号时域波形');
xlabel('t/s');
ylabel('幅度');
subplot(313);plot(t,v2);
title('第二声道音频信号时域波形');
xlabel('t/s');
ylabel('幅度');
figure(123)
subplot(311);plot(f1,mag1);
title('立体音频信号频域波形');
xlabel('f/Hz');
ylabel('幅度');
subplot(312);plot(f2,mag2);
title('第一声道音频信号频域波形');
xlabel('f/Hz');
ylabel('幅度');
subplot(313);plot(f3,mag3);
title('第二声道音频信号频域波形');
xlabel('f/Hz');
ylabel('幅度');
结果分析

图3.1.1

图3.1.2

从时域波形和频谱中我们可以发现，‘voice’是‘v1’从

时域波形和频谱中我们可以发现，‘voice’是‘v1’‘v2’叠加形成的，‘v1’‘v2’的时域波形和频谱相差并不是很大，因此在之后只对第一声道‘v1’处理以及分析，得到的结论都是一致的。对‘v1’的时域波形和频谱进行分析可以发现，采样时间为12s,声音的频率主要集中在0—5000 Hz频段附近，对应的归一化数字频率范围约为 0-0.208pi rad。

数字低通滤波器设计
设计流程
确定数字低通滤波器指标
截止频率 f = 4000Hz ；
阻带最小衰减 20dB ；
过渡带宽 0.1π 。
确定窗函数类型和长度
阻带最小衰减为20dB,查表后选择使用矩形窗，矩形窗的
过渡带宽度精确值为1.8/N，过渡带指标为0.1，故窗函数长度为18，即N=18；
生成滤波器的单位脉冲响应
利用hn=fir1(M,wc,'ftype',window)得到单位脉冲响应；
绘制频率响应曲线

利用‘freqz’函数，直接绘制出损耗函数曲线和相位特

性曲线。

实现方式
%% 矩形窗低通滤波器
wc1=4000*2*pi/fs; %截至频率
rs_low=20; %矩形窗阻带衰减
DB=0.1*pi; %过渡带
M=19; %矩形窗长度
hn=fir1(M-1,wc1/pi,'low',boxcar(M)); %矩形窗单位脉冲响应
figure(10) %矩形窗
freqz(hn,1); %直接使用freqz
v_low=filter(hn,1,v1);
%sound(v_low,fs); %播放经过低通滤波器声音
结果分析

图3.2.1

图3.2.2

从图3.2.1中可以看到阻带衰减约为20dB，过渡带约为

0.1π，数字截至频率约为0.166π，即为4000Hz;

从图3.2.2中可以看到4000Hz以上的频带已经被滤除；

播放经过该滤波器的音频，听到声音低沉，但更加清楚了，说明频率高的声音已经被滤除掉了；

综上所述，数字低通滤波器设计基本符合指标。

数字高通滤波器设计
设计流程
确定数字高通滤波器指标
截止频率 f = 4000Hz；
阻带最小衰减 40dB；
过渡带宽 0.1π；
确定窗函数类型和长度
阻带最小衰减为40dB,查表后选择使用汉宁窗，汉宁窗的
过渡带宽度精确值为6.2/N，过渡带指标为0.1，故窗函数长度为62，根据第一类线性相位FIRDF的约束条件，设计高通滤波器N应当为奇数，所以取N=63；
生成滤波器的单位脉冲响应
利用hn=fir1(M,wc,'ftype',window)得到单位脉冲响应；
绘制频率响应曲线

利用‘freqz’函数，直接绘制出损耗函数曲线和相位特

性曲线

实现方式
%%汉宁窗高通滤波器
wc2=3000*2*pi/fs; %截至频率
rs_high=40; %汉宁窗阻带衰减
DB=0.1*pi; %过渡带
Q0=ceil(6.2*pi/DB); %汉宁窗所需长度
Q=Q0+mod(Q0+1,2); %保证汉宁窗长度为奇数
hm=fir1(Q-1,wc2/pi,'high',hanning(Q));%汉宁窗冲激序列
figure(11)
freqz(hm,1); %汉宁窗
v_high=filter(hm,1,v1);
%sound(v_high,fs); %播放经过低通滤波器声音
结果分析

图3.3.1

图3.3.2

从图3.3.1中可以看到阻带衰减为44dB>40dB，过渡带为

0.1π，数字截至频率约为0.1797π，约为4000Hz;

从图3.2.2中可以看到0—4000Hz的频带已经被滤除；

播放经过该滤波器的音频，发现声音很轻，很尖锐，这是因为大部分音频都是低频率，已经被滤除了，只留下少部分高频声音；

综上所述，数字高通滤波器设计基本符合指标。

合成噪声
设计流程
生成一个频率为9000Hz的正弦信号作为单频噪声，转置后加到原始音频信号上；
实现方式
%% 加入正弦噪声
f=9000; %正弦噪声频率
N=length(voice);
t=(0:N-1)/fs;
%noise=0.01*randn(N,1); %加随机噪声
noisy=(0.1*sin(2*pi*f*t))'; %正弦噪声信号
v_n=v1+noisy; %加噪后语音信号
结果分析

图3.4.1

由图3.4.1可以看到，加入幅度为0.1，频率为9000Hz的

正弦噪声后，在时域波形中音频信号幅度在原有基础上发生了-0.1-0.1的规律变化，在频谱上则可以更加清楚的发现，在9000Hz处，有幅度为0.05的冲激。

播放加噪后的声音，发现声音中有一段持续的很尖锐的噪声，说明加噪成功；

综上所述，加噪后音频信号实际的变化符合理论值。

数字滤波器设计
设计方法
1. 设计目标
  滤除9000Hz处正弦噪声；经之前的分析可知，音频主要集
  中在0—5000Hz,因此选择使用低通滤波器，并不会过于损害原
  始音频的效果。由于本次设计中已经使用窗函数法设计过一次低
  通滤波器，因此这里采用等波纹最佳逼近法。
2. 确定低通滤波器设计指标
  通带衰减为1dB
  阻带衰减60dB
  逼近通带[0，7999]
  逼近阻带[8999，fs]
3. 生成滤波器的单位脉冲响应
  调用‘remezord’‘remez’分别估算等波纹最佳逼近指标
  以及生成滤波器的单位脉冲响应
4. 绘制频率响应曲线
  利用‘freqz’函数，直接绘制出损耗函数曲线和相位特
  性曲线

实现方式
%% 滤除9000hz正弦噪声数字滤波器设计
f=[7999,8999]; %边界（模拟）频率
m=[1,0];
rp=1; %通带衰减
rs=60; %阻带衰减
dat1=(10^(rp/20)-1)/(10^(rp/20)+1);%通带振荡波纹幅度
dat2=10^(-rs/20); %阻带振荡波纹幅度
rip=[dat1,dat2];
[M,fo,mo,w]=remezord(f,m,rip,fs);
M=M+1;
Hn=remez(M,fo,mo,w);
figure(100)
freqz(Hn,1);
v_n_low=filter(Hn,1,v_n);
%sound(v_n_low,fs); %播放经过低通滤波器的加噪语音信号
结果分析

图3.5.1

图3.5.2

由图3.5.1可知阻带衰减为60dB，过渡带约为0.1Π，

0.375pi位于阻带；

由图3.5.2可知9000Hz以上频带的音频全部被滤除；

播放经过该滤波器的音频，可以听到基本上已经没有噪声了，说明滤波器设计成功；

综上所述，该数字低通滤波器的设计达到了设计目的。

心得体会

在本次DSP大作业设计过程中我学习到了很多东西，主要包括两个方面；一方面是自己利用完成大作业的过程中，对理论知识进行梳理，在Matlab中加以实践，加深了理解，也更加有利于我对课程的理解，尤其是在将模拟信号转换成离散信号和各种滤波器的设计，包括数字低通滤波器、数字高通滤波器，以及不同的设计方法，包括窗函数法、频率采样法以及等波纹最佳逼近法等等；另一方面是自己在制作大作业的过程中，整个过程完全独立完成，增强了自己发现问题、分析问题、解决问题的能力，还对同学进行了一定的帮助。以下是我在本次大作业设计过程中遇到的问题、解决方案以及一些我的本次作业的独特之处：

在本次大作业设计过程中，由于对于matlab的使用有一定的基
础，因此基础编程包括绘制图形、快速傅里叶变换、正弦噪声信号的产生都没有问题，但是在对音频信和噪声进行叠加时遇到了问题，问题是叠加后得到了一个n*n的矩阵，直接把程序搞崩溃掉；后来在工作区内发现两个信号一个是n*1矩阵，一个是1*n矩阵，也就是一个是行向量，一个是列向量，于是对正弦噪声信号进行了转置，这个问题就得到了解决。
第二个问题是在滤波器设计过程中，对于滤波器的设计方法、滤

波器的指标、滤波器的实现都不是很熟悉，这是由于我在学习的过程中

只知其然不知其所以然，于是就下苦功夫仔细的学习课本的第五、六、七章，在理解了原理后对课本的例程进行研究，初步设计了滤波器，在网上查找了一些资料，对滤波器进行了改进，并学习了一些与滤波器有关的函数，最后，在又复习了一段时间后对滤波器有了更加深刻的理解于是又重新设计了滤波器，最终完成了设计。在这个过程中持续了很长时间，我由此发觉，一个设计过程应但是不断的推翻、改进、完善、创新的过程，需要不断的学习，不断的进步。

3、第三个是我帮助别人解决的一个问题，问题是别人在进行傅里叶变换时，得到的结果并不是理论上的幅度，我由于在之前的学习中学习了进行傅里叶变换，发现他在进行傅里叶变换时没有除以采样频率Fs；但是，解决了问题后，我并不能做出很好的解释，就只好让应用走到了理论的前头。在复习时，我对课本进行了仔细的研究，在课本的第四章第五节“模拟信号的频谱分析”中发现这么一句话，“模拟信号时域采样得到N点采样序列，经过DFT（FFT），在乘以T，就是模拟信号的傅里叶在频域的采样”。T就是采样间隔，也就是采样频率的倒数，这也就解释了为什么要除以采样频率Fs，我们在设计的过程中，对频域采样，然后DFT，却把DFT就当成了实际模拟信号的频谱。发现这句话的时候有一种豁然开朗的感觉，真的是“书读百遍，其异自现”。

4、在编写matlab代码时，我还注意了交互问题，因为老师再检查时不知道或者短时间很难找到在哪里播放声音，因此我简单的设计了一个选择交互，可以让老师在检查时选择要验听的音频。

最后，由于本次大作业设计赶上了期末复习，所以设计的过程中略显匆忙，滤波器的指标不是很完美，但是都基本符合标准，达到了锻炼自己的目的。最后再次感谢老师在这一学期的辛勤付出，2019新年快乐。