4 Rsym = 2.5e10; % symbol rate (sym/sec)
6 span = 6; % filter span
7 sps = 4; % samples per symbol
9 fs = Rsym * sps; % sampling freq (Hz)
12 t = (0 : 1 / fs : numSymbs / Rsym + (1.5 * span * sps - 1) / fs).';
16 EbN0 = 10 .^ (EbN0_db ./ 10);
22 EsN0 = EbN0 .* log2(M);
23 EsN0_db = 10 .* log10(EsN0);
26 data = randi([0 M - 1], numSymbs, 1);
27 modData = dpskmod(data, M, 0, 'gray');
29 x = txFilter(modData, rolloff, span, sps);
32 linewidthLO = 1e6; % Hz
34 [xPN, pTxLO] = phaseNoise(x, linewidthTx, linewidthLO, Tsamp);
36 snr = EbN0_db + 10 * log10(log2(M)) - 10 * log10(sps);
38 %%y = awgn(xPN, snr, 'measured');
41 r = rxFilter(y, rolloff, span, sps);
44 rSampled = r(sps*span/2+1:sps:(numSymbs+span/2)*sps);
45 [rSaPhEq, phiests] = phaseNoiseCorr(rSampled, M, 40);
47 demodData = dpskdemod(rSaPhEq, M, 0, 'gray');
49 [bitErrors, ber] = biterr(data, demodData.')
53 plot(repelem(-phiests, sps));
56 legend('estimate', 'actual');
57 title('Phase noise estimation');
61 plot(t(1:length(x)), real(normalizeEnergy(x, numSymbs*sps, 1)));
63 %%plot(t, real(rPhaseEq), 'r');
64 sampledTimes = t(sps*span/2+1:sps:(numSymbs+span/2)*sps);
65 %%plot(sampledTimes, real(rSampled), 'x');
66 plot(sampledTimes, real(normalizeEnergy(rSaPhEq, numSymbs, 1)), 'x');
67 legend('original signal', 'corrected received signal');
68 title('Phase noise correction, linewidth 1 MHz, $E_b/N_0=8$ dB');
69 ylabel('Real part of signals');
70 axis([t(1) t(300) -Inf +Inf]);