Commit | Line | Data |
---|---|---|
5fae0077 | 1 | numSymbs = 2^16; |
1eeb62fb AIL |
2 | M = 4; |
3 | ||
4 | Rsym = 2.5e10; % symbol rate (sym/sec) | |
5 | ||
6 | rolloff = 0.25; | |
7 | span = 6; % filter span | |
5fae0077 | 8 | sps = 8; % samples per symbol |
1eeb62fb AIL |
9 | |
10 | fs = Rsym * sps; % sampling freq (Hz) | |
11 | Tsamp = 1 / fs; | |
12 | ||
5e9be3c4 | 13 | t = (0 : 1 / fs : numSymbs / Rsym + (1.5 * span * sps - 1) / fs).'; |
1eeb62fb | 14 | |
5fae0077 | 15 | EbN0_db = 0:0.5:14; |
1eeb62fb AIL |
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); | |
5e9be3c4 AIL |
28 | berNoComp = zeros(1, plotlen); |
29 | berAdapt = zeros(1, plotlen); | |
30 | berMatlabAdapt = zeros(1, plotlen); | |
1eeb62fb AIL |
31 | |
32 | data = randi([0 M - 1], numSymbs, 1); | |
5e9be3c4 | 33 | modData = pskmod(data, M, pi / M, 'gray'); |
1eeb62fb AIL |
34 | x = txFilter(modData, rolloff, span, sps); |
35 | ||
36 | %% Simulate chromatic dispersion | |
5e9be3c4 | 37 | D = 17; % ps / (nm km) |
1eeb62fb | 38 | lambda = 1550; % nm |
5fae0077 | 39 | z = 3000; % km |
1eeb62fb | 40 | |
5fae0077 AIL |
41 | |
42 | [xCD, xCDkstart] = chromaticDispersion(x, D, lambda, z, Tsamp); | |
43 | ||
44 | TsampOrig = Tsamp; | |
1eeb62fb | 45 | |
1eeb62fb | 46 | for i = 1:plotlen |
5fae0077 AIL |
47 | sps = 8; |
48 | ||
1eeb62fb | 49 | snr = EbN0_db(i) + 10 * log10(log2(M)) - 10 * log10(sps); |
1eeb62fb AIL |
50 | |
51 | y = awgn(xCD, snr, 'measured'); | |
52 | ||
f9a73e9e | 53 | r = rxFilter(y, rolloff, span, sps); |
5fae0077 AIL |
54 | |
55 | sps = 2; | |
56 | Tsamp = TsampOrig * 4; | |
57 | ||
58 | [rCDComp, CDCompkstart] = CDCompensation(r, D, lambda, z, Tsamp); | |
f9a73e9e | 59 | rCDComp = normalizeEnergy(rCDComp, numSymbs*sps, 1); |
1eeb62fb | 60 | |
5fae0077 AIL |
61 | rSampled = rCDComp(2:2:end); |
62 | rNoCompSampled = r(2:2:end); | |
5e9be3c4 AIL |
63 | |
64 | %% rotate rNoCompSampled to match original data | |
65 | theta = angle(-sum(rNoCompSampled .^ M)) / M; | |
66 | %% if theta approx +pi/M, wrap to -pi/M | |
67 | if abs(theta - pi / M) / (pi / M) < 0.1 | |
68 | theta = -pi / M; | |
69 | end | |
5e9be3c4 AIL |
70 | rNoCompSampled = rNoCompSampled .* exp(-j * theta); |
71 | ||
f9a73e9e AIL |
72 | |
73 | %% Not entirely sure why, but after using FFT instead of time-domain | |
74 | %% convolution for simulating CD, we now need to do the same rotation | |
75 | %% for rSampled as well, but this time with a positive rotation. | |
76 | theta = angle(-sum(rSampled .^ M)) / M; | |
77 | if abs(theta + pi / M) / (pi / M) < 0.1 | |
78 | theta = +pi / M; | |
79 | end | |
80 | rSampled = rSampled .* exp(-1j * theta); | |
81 | ||
82 | ||
83 | ||
5e9be3c4 AIL |
84 | %% adaptive filter |
85 | adaptFilterOut = adaptiveCMA(rSampled); | |
86 | ||
87 | demodData = pskdemod(rSampled, M, pi / M, 'gray'); | |
88 | demodNoComp = pskdemod(rNoCompSampled, M, pi / M, 'gray'); | |
89 | demodAdapt = pskdemod(adaptFilterOut, M, pi / M, 'gray'); | |
90 | %%demodMatlabAdapt = pskdemod(matlabEq, M, pi / M, 'gray'); | |
1eeb62fb AIL |
91 | |
92 | [bitErrors, ber(i)] = biterr(data, demodData); | |
5e9be3c4 AIL |
93 | [bitErrors, berNoComp(i)] = biterr(data, demodNoComp); |
94 | [~, berAdapt(i)] = biterr(data, demodAdapt); | |
95 | %%[~, berMatlabAdapt(i)] = biterr(data, demodMatlabAdapt); | |
96 | ||
f9a73e9e AIL |
97 | %{ |
98 | if EbN0_db(i) == 14 | |
5e9be3c4 | 99 | figure(1); |
f9a73e9e AIL |
100 | scatterplot(normalizeEnergy(rSampled, numSymbs, 1)); |
101 | formatFigure; | |
102 | title('Constellation after CD comp.', 'interpreter', 'latex'); | |
103 | xlabel('In-Phase', 'interpreter', 'latex'); | |
104 | ylabel('Quadrature', 'interpreter', 'latex'); | |
105 | set(gca, 'FontSize', 18); | |
5e9be3c4 AIL |
106 | %%scatterplot(modData); |
107 | %%title('Original constellation'); | |
f9a73e9e AIL |
108 | scatterplot(normalizeEnergy(rNoCompSampled, numSymbs, 1)); |
109 | formatFigure; | |
110 | title('Constellation without CD comp.', 'interpreter', 'latex'); | |
111 | xlabel('In-Phase', 'interpreter', 'latex'); | |
112 | ylabel('Quadrature', 'interpreter', 'latex'); | |
113 | set(gca, 'FontSize', 18); | |
114 | %scatterplot(adaptFilterOut); | |
115 | %title('Constellation with CD compensation and adaptive filter'); | |
5e9be3c4 AIL |
116 | %scatterplot(matlabEq); |
117 | %title('Matlab equalizer'); | |
118 | ber(i) | |
119 | %berNoComp(i) | |
f9a73e9e AIL |
120 | %berAdapt(i) |
121 | %berMatlabAdapt(i) | |
5e9be3c4 | 122 | end |
f9a73e9e | 123 | %} |
1eeb62fb AIL |
124 | end |
125 | ||
126 | figure(1); | |
127 | clf; | |
128 | ||
129 | %% Plot simulated results | |
130 | semilogy(EbN0_db, ber, 'r', 'LineWidth', 2); | |
131 | hold on; | |
f9a73e9e AIL |
132 | semilogy(EbN0_db, berNoComp, 'm', 'LineWidth', 2); |
133 | semilogy(EbN0_db, berAdapt, 'Color', [0, 0.6, 0], 'LineWidth', 2); | |
5e9be3c4 | 134 | %%%semilogy(EbN0_db, berMatlabAdapt, 'c', 'LineWidth', 1.4); |
1eeb62fb AIL |
135 | |
136 | theoreticalPSK(EbN0_db, M, 'b', 'LineWidth', 1); | |
f9a73e9e AIL |
137 | %%legend({'CD + AWGN + CD comp.', 'CD + AWGN + CD comp.~+ CMA', ... |
138 | %% 'Theoretical AWGN'}, 'Location', 'southwest'); | |
139 | %%legend({'CD + AWGN + CD comp.', 'CD + AWGN', 'Theoretical AWGN'}, ... | |
140 | %% 'Location', 'southwest'); | |
141 | legend({'CD + AWGN + CD comp.', 'CD + AWGN', ... | |
142 | 'CD + AWGN + CD comp.~+ CMA', 'Theoretical AWGN'}, 'Location', ... | |
143 | 'Southwest'); | |
144 | ||
145 | %%title(strcat(num2str(M), '-PSK with chromatic dispersion and compensation')); | |
146 | title({'QPSK with chromatic dispersion and compensation', ... | |
147 | strcat(['$D = 17$ ps/(nm km), $z = ', num2str(z), '$ km'])}); | |
1eeb62fb AIL |
148 | grid on; |
149 | xlabel('$E_b/N_0$ (dB)'); | |
150 | ylabel('BER'); | |
151 | ||
152 | formatFigure; |