quinta-feira, 28 de novembro de 2019

Probabilidade de erro nas modulações 16 QAM e 16 PSK: simulação usando o Scilab

Probabilidade de erro de símbolo para as constelações 16-QAM e 16-PSK. Note que as constelações PSK e QAM não possuem a mesma energia na figura acima.
Em geral, um sistema de comunicação apresenta algum ruído. Esse ruído pode afetar a integridade dos dados recebidos e induzir a erros de recepção. Uma forma de comparar diferentes técnicas de modulação é justamente verificar qual é menos sensível ao ruído.

Nesta postagem comparamos a taxa de erro de símbolo das constelações 16-PSK e 16-QAM. Verificamos as curvas de probabilidade de erro experimentais (simulação) e teóricas usando um código Scilab. As curvas de erros teóricas são dadas por:
\begin{equation} P_{M-PSK} \leq 2Q\left(\sin(\pi/M)\sqrt{\frac{2E_s}{N_0}}\right) \end{equation} \begin{equation} P_{M-QAM} = 4\left(1 - \frac{1}{\sqrt{M}}\right)Q\left(\sqrt{\frac{3E_s}{(M-1)N_0}}\right) - 4\left(1 - \frac{1}{\sqrt{M}}\right)^2 Q^2\left(\sqrt{\frac{3E_s}{(M-1)N_0}}\right) \end{equation}

O código Scilab (com poucos comentários) que compara as modulações 16-PSK e 16-QAM é:

clc; xdel(winsid());  

function y=Qf(x)
    y = 0.5*erfc(x/sqrt(2));
endfunction

function y=dpsk(x)   /// demodulação 16 PSK
    y = 0*x;
    agg = 0:15;
    c = exp(2*%i*%pi*agg/M);
    tam = max(size(x));
    for k=1:tam
        dist = abs(x(k) - c);
        [aa,bb] = min(dist);
        y(k) = c(bb);
    end
endfunction

function y=dqam(x)   /// demodulação 16 QAM
    y = 0*x;
    tam = max(size(x))
    for k=1:tam
        xr = real(x(k));
        xi = imag(x(k));
        
        if xr<(-2) then xr = -3; end;
        if ((-2)<xr)&(xr<0) then xr = -1; end;
        if ((0)<xr)&(xr<2) then xr = 1; end;
        if xr>(2) then xr = 3; end;
        
        if xi<(-2) then xi = -3; end;
        if ((-2)<xi)&(xi<0) then xi = -1; end;
        if ((0)<xi)&(xi<2) then xi = 1; end;
        if xi>(2) then xi = 3; end;
        
        y(k) = xr + %i*xi;
    end
endfunction

function y=caltaxa(x)
    tam = max(size(x));
    y = 0;
    for k=1:tam
        if x(k)>0 then y = y + 1; end;
    end
    y = y/tam;
endfunction

//// curva teórica:
M = 16;
EsNo = 1:0.1:60;
Pepsk = 2*Qf(sin(%pi/M)*sqrt(2*EsNo));

px = (1-1/sqrt(M))*Qf(sqrt(3*EsNo/(M-1)));
Peqam = 4*px - 4*px.*px;
snrdb = 10*log10(EsNo);

subplot(2,1,1); plot(snrdb,Pepsk, snrdb, Peqam); 
title('Curva teórica - 16QAM e 16PSK. Taxa de erro de símbolo','fontsize',4);
legend('16-PSK','16-QAM',3);
aa = gca();
aa.log_flags = "nln";

////////  curva experimental:
N = 15000;
qx = (sign(rand(1,4*N,'n'))+1)/2;
p=1;
spsk = zeros(1,N);
for k=1:N
    b = qx(p:p+3);
    ag = b(1) + 2*b(2) + 4*b(3) + 8*b(4);
    p = p + 4;
    spsk(k) = exp(2*%i*%pi*ag/M);
end
qx = sign(rand(1,N,'n')) + 2*sign(rand(1,N,'n'));
qy = %i*(sign(rand(1,N,'n')) + 2*sign(rand(1,N,'n')));
sqam = qx + qy;

/// Ruido:
sqamn = sqam + 0.01*(rand(1,N,'n') + %i*rand(1,N,'n'));
subplot(2,2,3); plot(real(sqamn),imag(sqamn),'.');
title('16-QAM');
/// Ruido:
spskn = spsk + 0.01*(rand(1,N,'n') + %i*rand(1,N,'n'));
subplot(2,2,4); plot(real(spskn),imag(spskn),'.');
title('16-PSK');

Esmpsk = 1;
Esmqam = 10;
Ap = 0.5*sqrt(2)/0.9;
Aq = 0.5*sqrt(2)*sqrt(Esmqam)/0.9;
vtp = [];  // vetor de taxa de erro PSK
vsp = [];  // vetor de SNRdb
vtq = [];  // vetor de taxa de erro QAM
vsq = [];  // vetor de SNRdb
janH=waitbar('Calculando, aguarde um pouco!');
for k=1:20
    //// PSK:
    Ap = Ap*0.9;
    r = Ap*(rand(1,N,'n') + %i*rand(1,N,'n'));
    vr = variance(r);
    snrp = 10*log10(Esmpsk/vr);
    snpsk = spsk + r;
    sd = dpsk(snpsk);
    erros = abs(sd-spsk);
    taxap = caltaxa(erros);
    disp([taxap,snrp]);
    vtp = [vtp, taxap];
    vsp = [vsp, snrp];

    /// QAM:
    Aq = Aq*0.9;
    r = Aq*(rand(1,N,'n') + %i*rand(1,N,'n'));
    vr = variance(r);
    snrq = 10*log10(Esmqam/vr);
    sn = sqam + r;
    sd = dqam(sn);
    erros = abs(sd-sqam);
    taxaq = caltaxa(erros);
    disp([taxaq,snrq]);
    vtq = [vtq, taxaq];
    vsq = [vsq, snrq];
        
    vez = k/20;
    waitbar(vez,'Ainda calculando ...',janH);
end;
close(janH);

figure;
ff = gcf();
plot(snrdb,Pepsk,snrdb,Peqam,vsp,vtp,'r-*',vsq,vtq,'m->');
title('Curvas teórica e experimental - 16PSK. Taxa de erro de símbolo',...
'fontsize',4);
ff.background = 8;
aa = gca();
aa.log_flags = "nln";
aa.background = 8;
legend('Curva teórica 16PSK','Curva teórica 16QAM','Pontos experimentais PSK','Pontos experimentais QAM',3);
bb = aa.children;
bb(1).font_size = 3;
xgrid();


Que gera o seguintes gráfico:

Nenhum comentário:

Postar um comentário