Simulando curvas de taxa de erro - M-PAM.

Notar que a curva 2-PAM "descola" da curva teórica quando a relação sinal-ruído é maior que 12 dB.

Obs: essa postagem é exclusiva para "meus" ex-alunos da disciplina de Sistema de Comunicação II.

A taxa de erro de símbolo para uma constelação M-PAM que sofre a ação do ruído guassino pode ser expressa por:

Podemos usar um código Scilab para gerar as curvas simuladas e comparar com o desempenho teórico. O resultado é o gráfico de abertura desta postagem.

O código Scilab (não otimizado e com poucos comentários) segue abaixo.

function fq=Fq(z)
    s2 = sqrt(2);
    fq = 0.5*erfc(z/s2);
endfunction

function mp=mpam(Ns, Im)
   xx = rand(Ns,1,'uniform');
   if Im == 4 then
   mp = zeros(Ns,1)-3;
   for k=1:Ns
       if (xx(k)>0.25)&(xx(k)<0.50) then mp(k)=-1; end;
       if (xx(k)>0.50)&(xx(k)<0.75) then mp(k)=1; end;
       if (xx(k)>0.75) then mp(k)=3; end;
   end
   end;
   if Im == 8 then
   mp = zeros(Ns,1)-7;
   for k=1:Ns
       if (xx(k)>0.125)&(xx(k)<0.250) then mp(k)=-5; end;
       if (xx(k)>0.250)&(xx(k)<0.375) then mp(k)=-3; end;
       if (xx(k)>0.375)&(xx(k)<0.500) then mp(k)=-1; end;
       if (xx(k)>0.500)&(xx(k)<0.625) then mp(k)=1; end;
       if (xx(k)>0.625)&(xx(k)<0.750) then mp(k)=3; end;
       if (xx(k)>0.750)&(xx(k)<0.825) then mp(k)=5; end;
       if (xx(k)>0.825) then mp(k)=7; end;
   end
   end;
endfunction

function mx=dmod(ax, Im)
   if Im == 4 then
   mx = 0*ax-3;
   for k=1:max(size(mx))
       if (ax(k)>-2)&(ax(k)<0) then mx(k)=-1; end;
       if (ax(k)>0)&(ax(k)<2) then mx(k)=1; end;
       if (ax(k)>2) then mx(k)=3; end;
   end
   end;
   if Im == 8 then
   mx = 0*ax-7;
   for k=1:max(size(mx))
       if (ax(k)>-6)&(ax(k)<-4) then mx(k)=-5; end;
       if (ax(k)>-4)&(ax(k)<-2) then mx(k)=-3; end;
       if (ax(k)>-2)&(ax(k)<0) then mx(k)=-1; end;
       if (ax(k)>0)&(ax(k)<2) then mx(k)=1; end;
       if (ax(k)>2)&(ax(k)<4) then mx(k)=3; end;
       if (ax(k)>4)&(ax(k)<6) then mx(k)=5; end;
       if (ax(k)>6) then mx(k)=7; end;
   end
   end;
endfunction

close;
clc;

N = 42000;  // quanto maior o N, melhor a simulação

a2 = sign(rand(N,1,'normal'));
a4 = mpam(N,4);
a8 = mpam(N,8);

P2 = round(variance(a2));
P4 = round(variance(a4));
P8 = round(variance(a8));

disp([P2, P4, P8]);
n = rand(N,1,'normal');  // Pn = 1;

p=1;
for dB=0:18
    A = sqrt(P2/(10^(dB/10)));
    an = a2 + A*n;
    an_d = sign(an);
    erros = abs(a2 - an_d)/2;
    n_erros2(p) = sum(erros)/N;
    p = p + 1;
end

p=1;
for dB=0:18
    A = sqrt(P4/(10^(dB/10)));
    an = a4 + A*n;
    an_d = dmod(an,4);
    erros = abs(a4 - an_d)/2;
    for k=1:N
        if erros(k)>1 then erros(k)=1;
        end
    end
    n_erros4(p) = sum(erros)/N;
    p = p + 1;
end

p=1;
for dB=0:18
    A = sqrt(P8/(10^(dB/10)));
    an = a8 + A*n;
    an_d = dmod(an,8);
    erros = abs(a8 - an_d)/2;
    for k=1:N
        if erros(k)>1 then erros(k)=1;
        end
    end
    n_erros8(p) = sum(erros)/N;
    p = p + 1;
end

Am = 2;
dB = 0:18;
A2 = 2*sqrt((1.0)./(10.^(dB/10)));
ne2 = Fq(Am./A2);

A2 = 2*sqrt((5.0)./(10.^(dB/10)));
ne4 = 2*(1-1/4)*Fq(Am./A2);

A2 = 2*sqrt((21.0)./(10.^(dB/10)));
ne8 = 2*(1-1/8)*Fq(Am./A2);

n_erros2 = n_erros2 + 1e-21;  // para poder colocar em escala log
n_erros4 = n_erros4 + 1e-21; // para poder colocar em escala log

plot(dB,ne2,dB,ne4,dB,ne8,dB,n_erros2,'-o',...
dB,n_erros4,'->',dB,n_erros8,'-x'); xgrid;
a = gca();
a.log_flags = "nln"; // eixo y em escala log

xlabel('SRN - dB'); ylabel('Taxa de erro');
legend('2-PAM teórico','4-PAM teórico','8-PAM teórico',...
'2-PAM simulado','4-PAM simulado','8-PAM simulado');