CD_AWGN.m

   1 numSymbs = 10000;
   2 M = 4;
   3
   4 Rsym = 2.5e10; % symbol rate (sym/sec)
   5
   6 rolloff = 0.25;
   7 span = 6; % filter span
   8 sps = 4; % samples per symbol
   9
  10 fs = Rsym * sps; % sampling freq (Hz)
  11 Tsamp = 1 / fs;
  12
  13 t = (0 : 1 / fs : numSymbs / Rsym + (1.5 * span * sps - 1) / fs)';
  14
  15 EbN0_db = 0:0.2:14;
  16 EbN0 = 10 .^ (EbN0_db ./ 10);
  17
  18 Es = 1;
  19 Eb = Es / log2(M);
  20 N0 = Eb ./ EbN0;
  21
  22 EsN0 = EbN0 .* log2(M);
  23 EsN0_db = 10 .* log10(EsN0);
  24
  25 plotlen = length(EbN0);
  26
  27 ber = zeros(1, plotlen);
  28
  29 data = randi([0 M - 1], numSymbs, 1);
  30 modData = pskmod(data, M, 0, 'gray');
  31 x = txFilter(modData, rolloff, span, sps);
  32
  33 %% Simulate chromatic dispersion
  34 D = 20; % ps / (nm km)
  35 lambda = 1550; % nm
  36 z = 10; % km
  37
  38 xCD = chromaticDispersion(x, D, lambda, z, Tsamp);
  39 xCD = normalizeEnergy(xCD, numSymbs, 1);
  40
  41
  42 for i = 1:plotlen
  43   snr = EbN0_db(i) + 10 * log10(log2(M)) - 10 * log10(sps);
  44   noiseEnergy = 10 ^ (-snr / 10);
  45
  46   y = awgn(xCD, snr, 'measured');
  47
  48   yCDComp = CDCompensation(y, D, lambda, z, Tsamp);
  49
  50   r = rxFilter(yCDComp, rolloff, span, sps);
  51   %% normalize energy
  52   %r = normalizeEnergy(r, numSymbs, 1 + noiseEnergy);
  53
  54   rSampled = r(sps*span/2+1:sps:(numSymbs + span/2) * sps);
  55   demodData = pskdemod(rSampled, M, 0, 'gray');
  56
  57   [bitErrors, ber(i)] = biterr(data, demodData);
  58 end
  59
  60 figure(1);
  61 clf;
  62
  63 %% Plot simulated results
  64 semilogy(EbN0_db, ber, 'r', 'LineWidth', 2);
  65 hold on;
  66
  67 theoreticalPSK(EbN0_db, M, 'b', 'LineWidth', 1);
  68 legend({'CD + AWGN + CD compensation', 'AWGN only'}, 'Location', 'southwest');
  69
  70 title(strcat(num2str(M), '-PSK with chromatic dispersion and compensation'));
  71 grid on;
  72 xlabel('$E_b/N_0$ (dB)');
  73 ylabel('BER');
  74
  75 formatFigure;