Matlab系统建模实践

语音数字化量化噪声的改善（脉冲编码调制PCM）

一、主要任务

（1） 采集语音信号（模拟信号），并画出其时域图与频域图
（2） 对信号进行滤波处理
（3） 利用非均匀量化方法对信号进行编码
（4） 对编码过的信号进行解码，还原语言信号，并画出其时域图与频域图，并和原始信号进行比较
（5） 设计GUI界面，包括原始信号采集按钮（显示时域图与频域图）、滤波按钮（显示时域图与频域图）、编码按钮、解码按钮（显示时域图与频域图）

二、相关内容（略）

（1） 语音滤波降噪

（2） 均匀量化与非均匀量化

（3） 脉冲编码调制PCM

三、设计思路

（1） 采集音频并分析

调用matlab中读取音频文件的函数audioread采集音频信号，对这个信号加入噪声，噪声的频率为12.5KHz再通过快速傅里叶变换fft函数做出频域图。

[husheng0,fs]=audioread('H:\husheng_matlab0.mp3');%单声道
%sound(hushung0,fs);%播放
T=1/fs;%采样周期
l=length(husheng0);%信号长度
t=(0:l-1)*T;%时域范围
tt=(1:l);
noise=0.02*cos(12500*2*pi/fs*tt');%噪声
husheng=husheng0+noise;%加噪声
husheng1=fft(husheng,l)*T;
H=fftshift(husheng1);
S=abs(H);
F=2*pi*fs;%频率
FF=(-F/2+(0:l-1)*F/l)/(2*pi);

（2） 设计IIR低通滤波器

这次采用IIR巴特沃斯低通滤波器。因为人声的频率在300Hz到3.4KHz左右，所以设计了两种不同的技术指标的滤波器（10KHz和5KHz），以10KHz的为例，通带边界频率Fp为2000Hz，阻带边界频率Fs为1000Hz，通带最大衰减为2dB，阻带最小衰减为80dB。

Rp=2;Rs=80;%10khz
%Rp=1;Rs=60;%5khz
%Ft=48000;%我的test参数
%Ft=44100;%周杰伦的参数
Ft=8000;Fp=1000;Fs=2000;%备选参数1  10khz
%Fp=2000;Fs=4000; %参数2   5khz
Wp=2*pi*Fp/Ft;
Ws=2*pi*Fs/Ft; %求出待设计的模拟滤波器的边界频率
[N,Wn]=buttord(Wp,Ws,Rp,Rs, 's'); %低通滤波器的阶数和截止频率
[B,A]=butter(N,Wn, 's'); %S域频率响应的参数即：滤波器的传输函数
[Bz,Az]=bilinear(B,A,0.5); %利用双线性变换实现频率响应S域到Z域的变换

（3）音频滤波并分析

采用一维数字滤波器函数，调用上面产生的低通滤波器对信号进行滤波。

husheng2=filter(Bz,Az,husheng);%滤波处理

再对滤波后的信号进行快速傅里叶变换，做出时域图和频域图。

husheng3=fft(husheng2);
%sound(husheng2,fs);%播放一下
%滤波后时域频域图
H1=fftshift(husheng3);
S1=abs(H1);

（4） 编码与解码

编码与解码采用A律13折线法，均做成函数，这样子方便在后来的过程中进行调用。并且为了保证信号幅度的完整性，解码时加入一个归一化操作，即mas的值（信号幅度的最大值），以便解码恢复到原来的幅度。

%非均匀量化对信号编码
function code = PCMcoding(S)
z=sign(S);
    MaxS=max(abs(S));
    S=abs(S/MaxS);
    Q=2048*S;
    code=zeros(length(S),8);
    
    for i=1:length(S)
        if (Q(i)>=128)&&(Q(i)<=2048)
            code(i,2)=1;
        end
        if (Q(i)>32)&&(Q(i)<128)||(Q(i)>=512)&&(Q(i)<=2048)
            code(i,3)=1
        end
        if (Q(i)>=16)&&(Q(i)<32)||(Q(i)>=64)&&(Q(i)<128)||(Q(i)>=256)&&(Q(i)<512)||(Q(i)>=1024)&&(Q(i)<=2048)
            code(i,4)=1;
        end
    end
    
    N=zeros(length(S));
    
    for i=1:length(S)
        N(i)=bin2dec(num2str(code(i,2:4)))+1;
    end
    
    。。。。。
 
end

%对信号解码
function s = PCMdecoding(encode,max)
encode=(reshape(encode',8,length(encode)/8))';
    l=size(encode,1);
    a=[0,16,32,64,128,256,512,1024];
    b=[1 1 2 4 8 16 32 64];
    c=[0 1.5:15.5];
    for i=1:l
        x=encode(i,1);
。。。。。。。
        if Y==0
            k(i)=a(T)/2048;
        else
            k(i)=(a(T)+b(T)*c(Y))/2048;
        end
        if x==0
            s(i)=-k(i);
        else
            s(i)=k(i);
        end
    end
    s = s.*max;
 
end

（5） GUI界面设计

GUI界面设计共有四个按钮，分别控制：采集信号（时域图与频域图）、滤波器（滤波器特性曲线）、滤波后的信号（时域图与频域图）、编码、解码（时域图与频域图）。
点击按钮“原音频”进行信号的采集与分析，之后点击“滤波器”，做出所采用的低通滤波器的特性曲线，再点击“滤波”，做出滤波后的信号的时域图与频域图，接着点击“编码”对滤波后的信号进行编码，等待编码完成后将会显示“编码完成”点击“解码”，即做出解码后的信号的时域图与频域图。

四、结果与不足之处

（1）实验的结果符合预期的要求，能够完成对音频信号进行较低失真度的编码与解码。

（2）不足之处
在对信号进行采集的时候，时域坐标不符合实际，即时域坐标是采样点的个数，这个在请教老师之后及时进行了改正的。另外音频信号大多有左右声道，本次实践默认为单声道，缺少实际问题的多方面考虑，可以再添加一个对单声道和立体声的判断，在之后的编码与解码过程中也对单声道和双声道的问题进行考虑，做出两者不同的编码与解码。另外还可以加入失真度分析，即加入一个失真度测量算法，计算编码与解码过程中的失真度大小。

五、参考资料

（1） 《基于MATLAB的PCM编码解码实现》CSDN博主「virlaser」，原文链接：https://blog.csdn.net/weixin_44041884/article/details/86354999
（2） 《基于MATLAB的音频信号低通滤波实验》CSDN博主「和嘻泥」，原文链接：https://blog.csdn.net/qq_43160569/article/details/86544542
（3） 《如何用MATLAB设计GUI》CSDN博主「qq_26895463」，原文链接：https://blog.csdn.net/qq_26895463/article/details/81002338