基于2ask和74汉明码的通信系统实现

Posted 2021-07-03 小瓜瓜的小呱呱

是用matlab进行的仿真，将所有function和主函数放在了一起。

先对无汉明编码的数据进行2ask传输，载波是s2，加噪声，然后相干解调，计算误码率，并且画出每个步骤的波形；

然后对原始数据线74汉明编码然后再2ask传输，载波s1，同样的步骤，只不过多了一个解汉明码并且纠错，然后计算误码率。

最终对两条路线分别画出信噪比和误码率的曲线对比图。代码中低通带通滤波器还有一些零零散散的是参考的其他博主的文章，具体记不清了，就不放链接了，代码和结果图如下

clc
clear all
i=40;%40个码元
j=4000;
M=linspace(0,40,j);
I=70;
J=7000;
T=linspace(0,70,I);
K=linspace(0,70,J);
KK=linspace(0,40,j);

fc=1;%载波频率
fs=2000;
fm=i/10;%码元速率
b_exten=K;
dm=K;
c=T;
d=T;
a_exten=M;
x=(rand(1,i));
a=round(x);  %生成随机序列a
for n=1:40
    if a(n)==0
        for m=j/i*(n-1)+1:j/i*n  
            a_exten(m)=0;
        end
    else
        for m=j/i*(n-1)+1:j/i*n
            a_exten(m)=1;
        end
    end
end

%% 将a编码为b，b是经过74汉明码编码后的矩阵
b = reshape(a,i/4,[]); 
Q=[1 1 1;1 1 0;1 0 1;0 1 1];
G=[eye(4) Q];
U=b*G;
b=mod(U,2);


%% 将b展示，不重要
c = reshape(b',1,70);
for v=1:70
    if c(v)==0
        for m=J/I*(v-1)+1:J/I*v 
            b_exten(m)=0;
        end
    else
        for m=J/I*(v-1)+1:J/I*v
            b_exten(m)=1;
        end
    end
end

%% 载波
s1=cos(2*pi*fc*K);
s2=cos(2*pi*fc*KK);
%% 调制
e_2ask=b_exten.*s1;%载波
e_2askk=a_exten.*s2;

P=1;
countP=0;
for p=-50:20
    P=p+51;
    res(P)=0;
    ress(P)=0;
    Pe(P)=erfc(sqrt((power(10,p/10))/4))/2;
    for tt=1:200
i=40;
j=4000;
M=linspace(0,40,j);
I=70;
J=7000;
T=linspace(0,70,I);
K=linspace(0,70,J);
KK=linspace(0,40,j);

%% 给已调信号加噪声

rp  = 0.1;                   % Passband ripple
rs  = 50;                    % Stopband ripple
fs  = 200;                  % Sampling frequency
f   = [10 15 25 30];     % Cutoff frequencies
aq   = [0 1 0];               % Desired amplitudes
dev             = [10^(-rs/20) (10^(rp/20)-1)/(10^(rp/20)+1)  10^(-rs/20)]; 
[n,fo,ao,w]     = firpmord(f,aq,dev,fs);
PassFreqcoeff   = firpm(n,fo,ao,w);     % 计算得到FIR带通滤波器的系数
Snr             = p;
[SignalAddNoise,Noisenor2] = SignalAddPassBandNoise(e_2ask,Snr,PassFreqcoeff);
[SignalAddNoisee,Noisenor22] = SignalAddPassBandNoise(e_2askk,Snr,PassFreqcoeff);
e_2asknoise=SignalAddNoise;%加入噪声
e_2askknoise=SignalAddNoisee;
%% 相干解调
cm=e_2asknoise.*s1;
cm=cm-mean(cm);%因为是单极性波形，还有直流分量，应去掉
ccm=e_2askknoise.*s2;
ccm=ccm-mean(ccm);
[f,af] = T2F(K,cm);%通过低通滤波器
[K,cm] = Ipf(f,af,2*fm);
[f,af] = T2F(KK,ccm);%通过低通滤波器
[KK,ccm] = Ipf(f,af,2*fm);

%% 74编码的抽样判决
for m=0:69 %i=10   i是码元个数
    if (cm(1,m*100+50)+0.5)<0.5%100是1个码元的长度，+50就是正好每次都定位到每个码元中部
        d(m+1)=0;
        for j=m*100+1:(m+1)*100
            dm(1,j)=0; 
        end
    else
        d(m+1)=1;
        for j=m*100+1:(m+1)*100
            dm(1,j)=1;%否则判为1
        end
    end
end

%% 未经编码的抽样判决
for m=0:39 
    if ccm(1,m*100+50)<0
        dd(m+1)=0;
        for j=m*100+1:(m+1)*100
            ddm(1,j)=0; 
        end
    else
        dd(m+1)=1;
        for j=m*100+1:(m+1)*100
            ddm(1,j)=1;
        end
    end
end
%% 纠错 
eo=reshape(d,7,10); %e是解调后的序列转换成的10行7列矩阵
e=eo';
ee=e;
for i=1:10
    ss1(i)=xor(xor(e(i,3),xor(e(i,1),e(i,2))),e(i,5));
    ss2(i)=xor(xor(e(i,4),xor(e(i,1),e(i,2))),e(i,6));
    ss3(i)=xor(xor(e(i,4),xor(e(i,1),e(i,3))),e(i,7));
    ssum=4*ss1(i)+2*ss2(i)+ss3(i);
    if(ssum==3)
        ee(i,4)=~e(i,4);
        else if(ssum==5)
            ee(i,3)=~e(i,3);                
        else if(ssum==6)
            ee(i,2)=~e(i,2);    
        else if(ssum==7)
            ee(i,1)=~e(i,1);
        end
        end
        end 
    end
    eee(i,1)=ee(i,1);
    eee(i,2)=ee(i,2);
    eee(i,3)=ee(i,3);
    eee(i,4)=ee(i,4);
    ssum=0;
end
eeee=reshape(eee,1,40);

%% 计算经过编码后的误码率
count=0;

for m=1:40 
    if eeee(m)~=a(m)
        count=count+1;
    end
end
y(P,tt)=count/40;
res(P)=res(P)+y(P,tt);

%% 计算未经编码的误码率
countt=0;
for m=1:40 
    if dd(m)~=a(m)
        countt=countt+1;
    end
end
yy(P,tt)=countt/40;
ress(P)=ress(P)+yy(P,tt);
end
res(P)=res(P)/200;
ress(P)=ress(P)/200;
countP=countP+1;
end


%% 画图
figure(1)
subplot(421);
plot(M,a_exten)
title('初始序列');
axis([0,40,-1,2]);
subplot(422);
plot(K,b_exten)
title('74汉明编码后');
axis([0,70,-1,2]);
subplot(423);
plot(s1);
title('载波信号s1');
axis([0,7000,-1,1]);
subplot(424);
plot(K,e_2ask);
title('已调信号1');
 axis([0,70,-1.5,1.5]);
 subplot(425);
plot(K,e_2asknoise);
title('加入噪声s1的信号');
axis([0,70,-2,2]);
subplot(426);
plot(K,cm);
title('相干解调后波形');
 axis([0,70,-0.5,1]);
subplot(427);
plot(K,dm);
axis([0,70,-1,2]);
title('抽样判决后波形')

figure(2)
subplot(421);
plot(M,a_exten)
title('初始序列');
axis([0,40,-1,2]);
subplot(422);
plot(s1);
title('载波信号s2');
axis([0,4000,-1,1]);
subplot(423);
plot(KK,e_2askk);
title('已调信号2');
 axis([0,40,-1.5,1.5]);
 subplot(424);
plot(KK,e_2askknoise);
title('加入噪声s2的信号');
axis([0,40,-2,2]);
subplot(425);
plot(KK,ccm);
title('相干解调后波形');
 axis([0,40,-0.5,1]);
subplot(426);
plot(KK,ddm);
axis([0,40,-1,2]);
title('抽样判决后波形')

figure(3)
plot(res)
hold on
plot(ress)
hold on
plot(Pe)



function [f,sf]= T2F(t,st)
%利用FFT计算信号的频谱并与信号的真实频谱的抽样比较。
%脚本文件T2F.m定义了函数T2F，计算信号的傅立叶变换。
%Input is the time and the signal vectors,the length of time must greater
%than 2
%Output is the frequency and the signal spectrum
dt = t(2)-t(1);
T=t(end);
df = 1/T;
N = length(st);
f=-N/2*df : df : N/2*df-df;
sf = fft(st);
sf = T/N*fftshift(sf);
end

function [t,st]=Ipf(f,sf,B)
%此功能使用低通滤波器过滤输入数据
%Inputsf频率样本
% sf输入数据部分样本
% B矩形低通滤波器的低通带宽
%Outputst时间样本
% st输出数据时间样本
df = f(2)-f(1);
T = 1/df;
hf = zeros(1,length(f));%全零矩阵
bf = [-floor( B/df ): floor( B/df )] + floor( length(f)/2 );
hf(bf)=1;
yf=hf.*sf;
[t,st]=F2T(f,yf);
st = real(st);
end

function [SignalAddNoise,Noisenor2] = SignalAddPassBandNoise(SignalInput,Snr,PassFreqcoeff)
%{
    函数功能说明：给一段信号添加带限高斯白噪声；
    输入参数说明：
    SignalInput：   需要添加噪声的信号
    Snr:            添加噪声后的信噪比
    PassFreqcoeff:  带通滤波器系数
%}
    %% 输入参数判断
    % 输入的序列必须是N×1或者1×N;
    [m,n] = size(SignalInput);
    if (m~=1) && (n~=1)
        error('输入信号SignalInput必须为向量形式！')
    elseif m == 1
        SendSignalNum   = n;
        TempSignalInput = SignalInput;
    else
        SendSignalNum   = m;
        TempSignalInput = SignalInput.';
    end
    
    %% 给信号添加噪声
    Noise           = randn(1, 2 * SendSignalNum);               % 产生噪声点数是信号噪声的两倍
    NoiseFilterd    = filter(PassFreqcoeff, 1,  Noise);
    NoiseExtra      = NoiseFilterd(SendSignalNum - SendSignalNum/2:SendSignalNum + SendSignalNum/2 - 1);
    ENoise          = NoiseExtra * NoiseExtra';                 % 计算噪声的能量
    ESignal         = TempSignalInput * TempSignalInput';       % 计算信号的能量
    Noisenorl       = NoiseExtra / sqrt(ENoise);                % 噪声归一化
    Snr_ratio       = 10^(-Snr/20);
    Noisenor2       = Noisenorl * Snr_ratio * sqrt(ESignal);    % 得到带限的高斯白噪声
    SignalAddNoise  = TempSignalInput + Noisenor2;
 
end