通过matlab实现语音信号输入抽样量化转换.docx
《通过matlab实现语音信号输入抽样量化转换.docx》由会员分享,可在线阅读,更多相关《通过matlab实现语音信号输入抽样量化转换.docx(14页珍藏版)》请在冰豆网上搜索。
![通过matlab实现语音信号输入抽样量化转换.docx](https://file1.bdocx.com/fileroot1/2023-2/6/e03b0bf1-4524-4892-ba5f-f4b12238409c/e03b0bf1-4524-4892-ba5f-f4b12238409c1.gif)
通过matlab实现语音信号输入抽样量化转换
公司内部档案编码:
[OPPTR-OPPT28-OPPTL98-OPPNN08]
通过matlab实现语音信号输入抽样量化转换
通过matlab实现
1、语言信号的抽样、量化(A/D);
2、数字语言信号的存储;
3、一维随机噪声的产生;
4、IIR、FIR低通滤波器设计
IIR数字滤波器的设计方法a冲激响应不变法b双线性变换法
FIR数字滤波器设计Kaiser(凯塞窗)海明窗布莱克曼窗
5附上fr,和ideal_lp,
Fs=8000;
time=input('输入录音时间(整数秒):
');
fprintf('按回车开始录音%ds',time);
pause;
disp('录音中......');
x=audiorecorder(Fs,8,1);
recordblocking(x,time);
disp('录音结束\n');
fprintf('按回车回放录音');
pause;
play(x);
myrecording=getaudiodata(x);
audiowrite('lll.wav',myrecording,Fs);
y=audioread('lll.wav');
figure;
plot(y);
title('录音语音信号');
time=10;
N=Fs*time;
u=fft(y,N);
figure;
subplot(2,2,1);
plot(abs(u));
title('N点FFT幅度谱');
grid;
subplot(2,2,2);
plot(y);
title('语音信号');
grid;
L=length(y);
noise=0.05*randn(L,1);
yn=y+noise;
audioplayer(yn,Fs);
N=Fs*time;
m=fft(yn,N);
subplot(2,2,3);
plot(abs(m));
title('加噪N点FFT幅度谱');
grid;
subplot(2,2,4);
plot(yn);
title('语音信号加噪声');
grid;
以上是基础部分
紧接双线性变换iir
Ft=8000;
Fp=1000;
Fs=1200;
wp=2*pi*Fp/Ft;
ws=2*pi*Fs/Ft;
fp=2*Ft*tan(wp/2);
fs=2*Fs*tan(ws/2);
[n,Wn]=buttord(wp,ws,1,50,'s');
[Z,P,K]=buttap(n);
[Bap,Aap]=zp2tf(Z,P,K);
[b,a]=lp2lp(Bap,Aap,Wn);
[bz,ba]=bilinear(b,a,1);
[H,W]=freqz(bz,ba,64);
W1=W/pi;
H1=abs(H);
H2=20*log10(H1);
figure;
subplot(2,1,1);
stem(W1,H1);
grid;
xlabel('频率');
ylabel('幅度');
subplot(2,1,2);
stem(W1,H2);
grid;
xlabel('频率');ylabel('幅度(db)');
n=length(yn);
Y=fft(yn,n);
z22=filter(bz,ba,yn);
sound(z22);
m22=fft(z22);
figure;
subplot(2,2,1);
plot(abs(Y));
title('双线性IIR滤波前信号频谱');
grid;
subplot(2,2,2);
plot(abs(m22));
title('双线性IIR滤波后信号频谱');
grid;
subplot(2,2,3);
plot(z22);
title('双线性IIR滤波前信号波形');
grid;
subplot(2,2,4);
plot(y);
title('双线性IIR滤波后号波形');
grid;
凯赛窗fir
Ft=5000;
Fp=1000;
Fs=1200;
wp=2*Fp/Ft;
ws=2*Fs/Ft;
DB=ws-wp;%计算过渡带宽度
Rs=49;
beta=0.5842*(Rs-21)^0.4+0.07886*(Rs-21);%计算凯塞窗的控制参数
M=ceil((Rs-8)/2.285/DB);%计算凯塞窗所需阶数M
wc=(wp+ws)/2/pi;%计算理想低通滤波器通带截止频率
hn=fir1(M,wc,kaiser(M+1,beta));%调用firl函数计算低通FIRDF的h(n)
figure;
plot(hn);
xlabel('频率/Hz');ylabel('幅值');
title('数字滤波器幅频响应|H(ejOmega)|')
figure;
freqz(hn,1,512);
yy=fftfilt(hn,yn);
yy1=fft(yy,L);
yyy=abs(yy1);
figure;
subplot(2,2,1);
plot(abs(m));
title('加噪N点FFT幅度谱');
grid;
subplot(2,2,2);
plot(yn);
title('语音信号加噪声');
grid;subplot(2,2,3);
plot(yy);
title('凯赛处理后的语音信号波形图');
grid;subplot(2,2,4);
plot(yyy);
title('凯赛处理后的FFt频谱图');
grid;
到这里是个完整的程序,如果需要进行部分改变请从下面选取响应部分进行替换
冲击响应不变法iir
Ft=8000;
Fp=1000;
Fs=1200;
wp=2*pi*Fp/Ft;
ws=2*pi*Fs/Ft;
T=1;
[n,Wn]=buttord(wp,ws,1,50,'s');%计算模拟滤波器阶数N和截止频率wc
[B,A]=butter(n,Wn,'s');%计算相应模拟滤波器系统函数
[bz,ba]=impinvar(B,A,T);%用脉冲响应不变法将模拟转化成数字滤波器
figure;
freqz(bz,ba,512);
n=length(yn);
Y=fft(yn,n);
z22=filter(bz,ba,yn);
sound(z22);
m22=fft(z22);
figure;
subplot(2,2,1);
plot(abs(Y));
title('冲击响应不变法IIR滤波前信号频谱');
grid;
subplot(2,2,2);
plot(abs(m22));
title('冲击响应不变法IIR滤波后信号频谱');
grid;
subplot(2,2,3);
plot(z22);
title('冲击响应不变法IIR滤波前信号波形');
grid;
subplot(2,2,4);
plot(y);
title('冲击响应不变法IIR滤波后号波形');
grid;
海明窗fir
Ft=5000;
Fp=1000;
Fs=1200;
wp=2*Fp/Ft;
ws=2*Fs/Ft;
tr_width=ws-wp;
disp('海明窗设计低通滤波器参数:
');
M=ceil(6.6*pi/tr_width)+1;
disp(['滤波器的长度为',num2str(M)]);
n=0:
M-1;
wc=(ws+wp)/2;%理想LPF的截止频率
hd=ideal_lp(wc,M);
w_ham=(hamming(M))';
hn=hd.*w_ham;
[db,mag,pha,gfd,w]=fr(hn,1);
delta_w=2*pi/1000;
Rp=-(min(db(1:
1:
wp/delta_w+1)));%求出实际通带波动
disp(['实际带通波动为',num2str(Rp)]);
As=-round(max(db(ws/delta_w+1:
1:
501)));%求出最小阻带衰减
disp(['最小阻带衰减为',num2str(As)]);
figure;
subplot(2,2,1);
plot(n,hd);
title('理想冲击响应');
axis([0,M-1,-0.1,0.3]);
ylabel('hd(n)');
subplot(2,2,2);
plot(n,w_ham);
title('海明窗');
axis([0,M-1,0,1.1]);
ylabel('w(n)');
subplot(2,2,3);
plot(n,hn);
title('实际冲激响应');
axis([0,M-1,-0.1,0.3]);
xlabel('n');
ylabel('h(n)');
subplot(2,2,4);
plot(w/pi,db);
title('幅度响应(db)');
axis([0,1,-100,10]);
grid;
xlabel('以pi为单位的频率');
ylabel('分贝数');
yy=fftfilt(hn,yn);
yy1=fft(yy,L);
yyy=abs(yy1);
figure;
subplot(2,2,1);
plot(abs(m));
title('加噪N点FFT幅度谱');
grid;
subplot(2,2,2);
plot(yn);
title('语音信号加噪声');
grid;subplot(2,2,3);
plot(yy);
title('海明窗处理后的语音信号波形图');
grid;subplot(2,2,4);
plot(yyy);
title('海明窗处理后的FFt频谱图');
grid;
布莱克曼窗fir
Ft=5000;
Fp=1000;
Fs=1200;
wp=2*Fp/Ft;
ws=2*Fs/Ft;
r_width=ws-wp;
disp('布莱克曼窗设计低通滤波器的参数:
');
M=ceil(11.0*pi/tr_width)+1;
disp(['滤波器的长度为',num2str(M)]);
n=0:
M-1;%理想LPF的截止频率?
wc=(ws+wp)/2;
hd=ideal_lp(wc,M);
w_bla=(blackman(M))';
hn=hd.*w_bla;
[db,mag,pha,gfd,w]=fr(hn,1);
delta_w=2*pi/1000;
Rp=-(min(db(1:
1:
wp/delta_w+1)));%求出实际通带波动
disp(['实际带通波动为',num2str(Rp)]);
As=-round(max(db(ws/delta_w+1:
1:
501)));%求出最小阻带衰减
disp(['最小阻带衰减db',num2str(As)]);
figure;
subplot(2,2,1);
plot(n,hd);
title('理想冲击响应');
axis([0,M-1,-0.1,0.3]);
ylabel('hd(n)');
subplot(2,2,2);
plot(n,w_bla);
title('布莱克曼窗');
axis([0,M-1,0,1.1]);
ylabel('w(n)');
subplot(2,2,3)
plot(n,hn);
title('实际冲激响应');
axis([0,M-1,-0.1,0.3]);
xlabel('n');
ylabel('hn(n)');
subplot(2,2,4);
plot(w/pi,db);
title('幅度响应(db)');
axis([0,1,-100,10]);
grid;
xlabel('以pi为单位的频率');
ylabel('分贝数');
yy=fftfilt(hn,yn);
yy1=fft(yy,L);
yyy=abs(yy1);
figure;
subplot(2,2,1);
plot(abs(m));
title('加噪N点FFT幅度谱');
grid;
subplot(2,2,2);
plot(yn);
title('语音信号加噪声');
grid;subplot(2,2,3);
plot(yy);
title('布莱克曼处理后的语音信号波形图');
grid;subplot(2,2,4);
plot(yyy);
title('布莱克曼处理后的FFt频谱图');
grid;
以下2个需要单独保存
function[db,mag,pha,gfd,w]=fr(b,a)%N阶差分方程所描述的系统频响函数的m函数文件fr.m
%求解系统响应
%db为相位振幅(db)
%mag为绝对振幅
%pha为相位响应
%grd为群延时
%w为频率样本点矢量
%b为Ha(z)分析多项式系数(对FIR而言,b=h)
%a为Hz(z)分母多项式系数(对FIR而言,a=1)
[H,w]=freqz(b,a,1000,'whole');
H=(H(1:
501))';
w=(w(1:
501))';
mag=abs(H);
db=20*log10((mag+eps)/max(mag));
pha=angle(H);
gfd=grpdelay(b,a,w);
functionhd=ideal_lp(wc,M)%理想低通滤波器单位抽样响应hd(n)的m函数文件ideal.m
%IdealLowPassfiltercomputation
%--------------------------------
%[hd]=ideal_lp(wc,M);
%hd=idealimpulseresponsebetween0toM-1
%wc=cutofffrequencyinradians
%M=lengthoftheidealfilter
%
alpha=(M-1)/2;
n=(0:
1:
(M-1));
m=n-alpha+eps;%addsmallestnumbertoavoidividedbyzero
hd=sin(wc*m)./(pi*m);