基于matlab的语音信号PCM编码-2PSK传输仿真

主要内容：采用PCM编码，2PSK结合simulink实现语音信号的传输仿真。（注：代码为本人参考部分csdn文章改写）封面为部分运行结果。

PCM，2PSK相关知识就不多说了，直接贴代码吧：

%% 产生信号
load chirp % matlab自带语音信号
y1=[y,y];
x=y1(1:20000); %取前20000个采样点
sound(x,Fs);
%% PCM编码
x1=x/0.8.*2048;
yy=pcm_encode(x1);
time=0.0001:0.0001:16;
time_workspace=time';
pcm_code1=yy';
workspace_var=[time_workspace,pcm_code1];
figure;
subplot(1,1,1);
stem(yy(1:80),'.');
title('PCM编码后的波形');
figure;
subplot(1,1,1);
plot(x);
title('原始语音信号');
%% 基带信号求反
%由于PSK中的是双极性信号，因此对上面所求单极性信号取反来与之一起构成双极性码
st1=yy(1:160000);
t=1:1:160000;
st2=t;
j=160000;
for k=1:j
    if st1(k)>=1
        st2(k)=0;
    else
        st2(k)=1;
    end
end
figure;
subplot(1,1,1)
stem(st2(1:80),'.');
title('基带信号反码');
axis([0,80,0,1]);
pcm_code2=st2';
%% 译码
demodata=simout(2:160001);
zz=pcm_decode(demodata,0.8);
figure;
subplot(1,1,1);
plot(zz);
title('译码的语音信号');
sound(zz,Fs);
figure;
plot(x,'b');
hold on
plot(zz,'r');
legend('原语音信号','恢复的语音信号');
title('语音信号对比');
% sound(zz,Fs);
%% 误码计算
code1=code(2:160001);
[number1,ratio1]=symerr(code1,demodata);
%% PCM编码函数
function y = pcm_encode( x )
y=zeros(length(x),8); %存储矩阵（全零）
z=sign(x); %判断x的正负
x=abs(x);%取绝对值
%%段落码判断段区间的取值范围为前开后闭区间
for k=1:length(x)
    %符号位的判断
    if z(k)>0
        y(k,1)=1;
    elseif z(k)<0
        y(k,1)=0;
    end
    if x(k)>128 && x(k)<=2048 %在第五段与第八段之间，段位码第一位都为“1”
        y(k,2)=1;
    end
    if (x(k)>32 && x(k)<=128) || (x(k)>512 && x(k)<=2048)
        y(k,3)=1; %在第三四七八段内，段位码第二位为“1”
    end
    if (x(k)>16&&x(k)<=32)||(x(k)>64&&x(k)<=128)||(x(k)>256&&x(k)<=512)||(x(k)>1024&&x(k)<=2048)
        y(k,4)=1; %在二四六八段内，段位码第三位为“1”
    end
end
%段内码判断程序
N=zeros(1,length(x));
for k=1:length(x)
    N(k)=y(k,2)*4+y(k,3)*2+y(k,4)+1; %找到x位于第几段
end
a=[0,16,32,64,128,256,512,1024]; %量化间隔
b=[1,1,2,4,8,16,32,64]; %除以16，得到每段的最小量化间隔
for m=1:length(x)
    q=ceil((x(m)-a(N(m)))/b(N(m))); %求出在段内的位置
    if q==0
        y(m,(5:8))=[0,0,0,0]; %如果输入为零则输出“0”
    else
        k=num2str(dec2bin(q-1,4)); %编码段内码为二进制
        y(m,5)=str2double(k(1));
        y(m,6)=str2double(k(2));
        y(m,7)=str2double(k(3));
        y(m,8)=str2double(k(4));
    end
end
%将N行8列矩阵转换为1行8*N列的矩阵
y=y';
y=reshape(y,1,length(x)*8);
end
%% PCM译码函数
function x=pcm_decode(y,max)
%将1行8*N列的矩阵转换为N行8列矩阵
y=reshape(y,8,length(y)/8);
y=y';%PCM译码
n=size(y,1); %求出输入码组的个数
a=[0,16,32,64,128,256,512,1024]; %段落起点值
b=[1,1,2,4,8,16,32,64]; %每段的最小量化间隔
for k=1:n
    t1=y(k,1); %取符号
    t2=y(k,2)*4+y(k,3)*2+y(k,4)+1; %判断段落位置
    t3=y(k,5)*8+y(k,6)*4+y(k,7)*2+y(k,8); %判断段内位置
    if t3==0 %段内码为零时
        m(k)=(a(t2)+1+0.5*b(t2))/2048*max;
    else
        m(k)=(a(t2)+b(t2)*t3+0.5*b(t2))/2048*max; %还原出量化后的电平值
    end        %判断符号位
    if t1==0
        x(k)=-m(k);
    else
        x(k)=m(k);
    end
end
end

部分运行结果：

 

 

代码需要结合simlink模型使用，可根据相关知识写一个，不复杂。代码中从工作区导入的变量名：workspace_var，使用两个from workspace模块分别设置：[time_workspace,pcm_code1]，[time_workspace,pcm_code2];从simulink导出的数据变量名：simout，使用to workspace模块。

注意事项：电脑配置低的simulink运行可能比较慢，我设置的仿真时间16，如果需要改这个时间，则相对的需要改代码中的simout(2:160001)，否则可能出现索引数组超出范围。

可以去我上传的资源中直接下载simulink文件，和上面代码配套。