150 lines
4.5 KiB
Matlab
150 lines
4.5 KiB
Matlab
clear;
|
|
close all;
|
|
|
|
Fe = 10000; % fréquence d'échantillonage (Hz)
|
|
Te = 1/Fe; % période d'échantillonage (s)
|
|
N = 50000; % nombre de bits envoyés
|
|
Rb = 2000; % débit binaire
|
|
Fp = 2000; % fréquence porteuse
|
|
|
|
load('bits.mat')
|
|
%bits = randi([0, 1], 1, N); % bits envoyés
|
|
|
|
M = 2^2; % signal QPSK
|
|
Ts = log2(M)/Rb; % période symbole
|
|
Ns = floor(Ts/Te); % facteur suréchantillonage
|
|
T = (0:N*Ns/log2(M)-1) * Te; % échelle temporelle
|
|
|
|
fc = 1.25*Fe; % fréquence de coupure du lowpass
|
|
ordre = 100; % ordre du filtre
|
|
lowpass = (2*fc/Fe)*sinc(2*fc*[-(ordre-1)*Te/2 : Te : (ordre-1)*Te/2]); % coefficients du filtre passe-bas
|
|
|
|
alpha = 0.35; % roll off du racine de cosinus surélevé
|
|
span = 10; % largeur du racine de cosinus surélevé
|
|
h = rcosdesign(alpha, span, Ns); % mise en forme: racine de cosinus surélevé
|
|
h_c = [ 1 zeros(1, Ns-1) ]; % propagation/canal: dirac
|
|
h_r = h; % réception: racine de cosinus surélevé
|
|
|
|
g = conv(conv(h, h_c), h_r); % réponse impulsionnelle globale
|
|
n0 = 1; % déterminé en traçant le diagramme de l'oeil, ou g
|
|
|
|
|
|
|
|
X_grp = reshape(bits, log2(M), N/log2(M)).'; % on regroupe les bits
|
|
X_ma = X_grp * 2 - 1; % mapping
|
|
X_co = X_ma(:,1) + 1j * X_ma(:,2); % symboles complexes
|
|
|
|
X_kr = kron(X_co.', [1 zeros(1, Ns-1)]); % Suréchantillonnage
|
|
|
|
X_env = [ X_kr zeros(1, span*Ns/2) ]; % on rajoute des zéros pour le retard
|
|
X_env = filter(h, 1, X_env); % enveloppe complexe
|
|
X_env = X_env(span*Ns/2 + 1 : end); % on dégage les zéros
|
|
|
|
figure;
|
|
preview = 250;
|
|
plot(T(1:preview), real(X_env(1:preview)));
|
|
title("Partie réelle de l'enveloppe complexe du signal");
|
|
xlabel("Temps (s)");
|
|
ylabel("Amplitude");
|
|
|
|
figure;
|
|
pwelch(real(X_env)); % on trace la DSP de l'enveloppe
|
|
title("DSP de l'enveloppe complexe associé");
|
|
|
|
figure;
|
|
oeil = reshape(real(X_env), Ns, length(X_env)/Ns);
|
|
plot(oeil);
|
|
title("Diagramme de l'oeil, signal non bruité");
|
|
xlabel("Temps (s)");
|
|
ylabel("Amplitude");
|
|
|
|
|
|
|
|
X_port = X_env .* exp( 2*pi*j*Fp*T ); % transposition en fréquence
|
|
|
|
figure;
|
|
pwelch(real(X_port)); % on trace la DSP du signal modulé sur porteuse
|
|
title("DSP du signal modulé sur porteuse");
|
|
|
|
X_reel = real(X_port); % passage en réel
|
|
|
|
figure;
|
|
plot(T(1:preview), X_reel(1:preview));
|
|
title("Signal réel transposé en fréquence");
|
|
xlabel("Temps (s)");
|
|
ylabel("Amplitude");
|
|
|
|
X_tran = filter(h_c, 1, X_reel); % passage dans le canal
|
|
|
|
|
|
|
|
P_x = mean(abs(X_reel).^2);
|
|
EbN0_db = linspace(0, 8, 256);
|
|
EbN0 = 10.^(EbN0_db./10);
|
|
TEBs = [];
|
|
TESs = [];
|
|
%figure;
|
|
for e=EbN0
|
|
|
|
sigma2_x = P_x * Ns / (2 * log2(M) * e); % calcul de sigma^2
|
|
X_bruit = X_tran + sqrt(sigma2_x) * randn(1, length(X_tran)); % signal bruité
|
|
|
|
X_cos = X_bruit .* cos( 2*pi*Fp*T ); % on prend la partie en phase du signal
|
|
X_cos = [ X_cos zeros(1, ordre/2) ]; % on rajoute des zéros
|
|
X_cos = filter(lowpass, 1, X_cos); % passe bas
|
|
X_cos = X_cos(ordre/2 + 1:end); % on dégage les zéros
|
|
|
|
X_sin = X_bruit .* sin( 2*pi*Fp*T ); % on prend la partie en quadrature du signal
|
|
X_sin = [ X_sin zeros(1, ordre/2) ]; % on rajoute des zéros
|
|
X_sin = filter(lowpass, 1, X_sin); % passe bas
|
|
X_sin = X_sin(ordre/2 + 1:end); % on dégage les zéros
|
|
|
|
X_base = X_cos - 1j * X_sin; % on reconsitue les symboles
|
|
X_base = -2 * X_base; % sinon ça marche pas ?
|
|
|
|
%plot(X_env)
|
|
%hold
|
|
%plot(X_base)
|
|
|
|
X_recu = [ X_base zeros(1, span*Ns/2) ]; % on rajoute des zéros pour le retard
|
|
X_recu = filter(h_r, 1, X_recu); % signal reçu
|
|
X_recu = X_recu(span*Ns/2 + 1:end); % on dégage les zéros
|
|
|
|
X_sample = X_recu( n0 : Ns : N*Ns/log2(M) ); % signal échantillonné
|
|
|
|
X_demod_re = real(X_sample) > 0; % detection par seuil des réels
|
|
X_demod_im = imag(X_sample) > 0; % detection par seuil des imaginaires
|
|
|
|
recu = double( [ X_demod_re' X_demod_im' ] ); % on reconstitue les bits à partir des symboles
|
|
TES = mean((bi2de(recu) - bi2de(X_grp)).^2); % on calcule l'erreur symbole
|
|
|
|
recu = reshape(recu', size(bits)); % on transforme les symboles en bits
|
|
TEB = mean((recu - bits).^2); % on calcule l'erreur binaire
|
|
|
|
TESs = [ TESs TES ];
|
|
TEBs = [ TEBs TEB ];
|
|
|
|
%plot(X_sample, '+');
|
|
%hold;
|
|
%plot(X_co, '*');
|
|
%hold;
|
|
%axis square;
|
|
%grid;
|
|
%axis([-2, 2, -2, 2]);
|
|
%title(sprintf("Eb/N0 = %.2fdB, erreur = %.2f%s", 10*log10(e), TEB*100, "%"));
|
|
%legend("symboles reçus", "symboles théoriques");
|
|
%drawnow;
|
|
|
|
end
|
|
|
|
TEB_theorique = qfunc(sqrt(2*EbN0));
|
|
|
|
figure;
|
|
semilogy(EbN0_db, TESs, '+');
|
|
hold;
|
|
semilogy(EbN0_db, TEBs, '+');
|
|
plot(EbN0_db, TEB_theorique);
|
|
title("TEB = f(E_b/N_0)");
|
|
xlabel("E_b/N_0 (dB)");
|
|
ylabel("TEB");
|
|
legend("TES numérique", "TEB numérique", "TEB théorique"); |