-function phasenoise_AWGN(rolloff, M, numSymbs)
- %% Set defaults for inputs
- if nargin < 3
- numSymbs = 1000;
- end
- if nargin < 2
- M = 2;
- end
- if nargin < 1
- rolloff = 0.25;
- end
-
- plotted = 0;
-
- Rsym = 2.5e10; % symbol rate (sym/sec)
-
- span = 6; % filter span
- sps = 4; % samples per symbol
-
- fs = Rsym * sps; % sampling freq (Hz)
- Tsamp = 1 / fs;
-
- t = (0 : 1 / fs : numSymbs / Rsym + (1.5 * span * sps - 1) / fs)';
-
-
- EbN0_db = 0:0.2:14;
- EbN0 = 10 .^ (EbN0_db ./ 10);
-
- Es = 1;
- Eb = Es / log2(M);
- N0 = Eb ./ EbN0;
-
- EsN0 = EbN0 .* log2(M);
- EsN0_db = 10 .* log10(EsN0);
-
- plotlen = length(EbN0);
-
- ber = zeros(1, plotlen);
-
+numSymbs = 5e4;
+M = 4;
+rolloff = 0.5;
- data = randi([0 M - 1], numSymbs, 1);
- modData = pskmod(data, M, 0, 'gray');
+Rsym = 2.5e10; % symbol rate (sym/sec)
- x = txFilter(modData, rolloff, span, sps);
+span = 6; % filter span
+sps = 8; % samples per symbol
- linewidthTx = 0;%1e5; % Hz
- linewidthLO = 1e6; % Hz
- %%linewidthTx = Rsym * 1e-4; % Hz
- %%linewidthLO = Rsym * 1e-3; % Hz
+fs = Rsym * sps; % sampling freq (Hz)
+Tsamp = 1 / fs;
- totalIterations = 1;
+t = (0 : 1 / fs : numSymbs / Rsym + (1.5 * span * sps - 1) / fs).';
- for iter = 1:totalIterations
- [xPN, pTxLO] = phaseNoise(x, linewidthTx, linewidthLO, Tsamp);
- for i = 1:plotlen
- snr = EbN0_db(i) + 10 * log10(log2(M)) - 10 * log10(sps);
- noiseEnergy = 10 ^ (-snr / 10);
+EbN0_db = 0:0.2:14;
+EbN0 = 10 .^ (EbN0_db ./ 10);
- y = awgn(xPN, snr, 'measured');
+Es = 1;
+Eb = Es / log2(M);
+N0 = Eb ./ EbN0;
- r = rxFilter(y, rolloff, span, sps);
- %% normalize energy
- %r = normalizeEnergy(r, numSymbs, 1 + noiseEnergy);
+EsN0 = EbN0 .* log2(M);
+EsN0_db = 10 .* log10(EsN0);
- [rPhaseEq, phiests] = phaseNoiseCorr(r, M, 40 * sps);
+plotlen = length(EbN0);
- rSampled = rPhaseEq(sps*span/2+1:sps:(numSymbs + span/2) * sps);
- demodData = pskdemod(rSampled, M, 0, 'gray')';
+berPSK = zeros(1, plotlen);
+berDEPSK = zeros(1, plotlen);
+berDPSK = zeros(1, plotlen);
- %%[bitErrors, ber(i)] = biterr(data, demodData);
- [zzz, thisBER] = biterr(data, demodData);
- ber(i) = ber(i) + thisBER / totalIterations;
+data = randi([0 M - 1], numSymbs, 1);
+pskSym = pskmod(data, M, pi / M, 'gray');
+%% DEPSK: Part VII, M.G. Taylor (2009)
+depskSym = pskmod(data, M, 0, 'gray');
+for i = 2:numSymbs
+ depskSym(i) = depskSym(i) * depskSym(i-1);
+end
- if plotted == 0 && EbN0_db(i) >6 && ber(i) > 1e-1
- plotted=1
- figure(1234);
- plot(-phiests);
- hold on;
- plot(pTxLO);
- legend('estimate', 'actual');
- hold off;
-
- figure(100);
- %plot(t(1:length(x)), real(x));
- %%plot(t, real(x(1:length(t))));
- length(t)
- length(x)
- hold on
- %%plot(t(1:length(xPhaseNoise)), real(xPhaseNoise));
-
- plot(t, real(r), 'g')
-
- sampledTimes = t(sps*span/2+1:sps:(numSymbs+span/2)*sps);
-
- plot(sampledTimes, real(rSampled), 'x')
- hold off
-
- end
-
+dpskSym = dpskmod(data, M, pi / M, 'gray');
+
+xPSK = txFilter(pskSym, rolloff, span, sps);
+xDEPSK = txFilter(depskSym, rolloff, span, sps);
+xDPSK = txFilter(dpskSym, rolloff, span, sps);
+
+linewidthTx = 0; % Hz
+linewidthLO = 5e6; % Hz
+%linewidthLO = Rsym * 1e-3;
+
+iterations = 1;
+avgSa = 40;
+
+TsampOrig = Tsamp;
+
+for it = 1 : iterations
+ [xPSKpn, pTxLoPSK] = phaseNoise(xPSK, linewidthTx, linewidthLO, Tsamp);
+ [xDEPSKpn, pTxLoDEPSK] = phaseNoise(xDEPSK, linewidthTx, linewidthLO, Tsamp);
+ [xDPSKpn, pTxLoDPSK] = phaseNoise(xDPSK, linewidthTx, linewidthLO, Tsamp);
+
+ for i = 1:plotlen
+ Tsamp = TsampOrig;
+ sps = 8;
+
+ snr = EbN0_db(i) + 10 * log10(log2(M)) - 10 * log10(sps);
+ noiseEnergy = 10 ^ (-snr / 10);
+
+ yPSK = awgn(xPSKpn, snr, 'measured');
+ yDEPSK = awgn(xDEPSKpn, snr, 'measured');
+ yDPSK = awgn(xDPSKpn, snr, 'measured');
+
+ rPSK = rxFilter(yPSK, rolloff, span, sps);
+ rDEPSK = rxFilter(yDEPSK, rolloff, span, sps);
+ rDPSK = rxFilter(yDPSK, rolloff, span, sps);
+
+ sps = 2;
+ Tsamp = TsampOrig * 4;
+
+ rPSKSamp = rPSK(1:2:end);
+ rDEPSKSamp = rDEPSK(1:2:end);
+ rDPSKSamp = rDPSK(1:2:end);
+
+ [rPSKSampEq, phiestsPSK] = phaseNoiseCorr(rPSKSamp, M, pi/M, avgSa);
+ [rDEPSKSampEq, phiestsDEPSK] = phaseNoiseCorr(rDEPSKSamp, M, 0, avgSa);
+
+ demodPSK = pskdemod(rPSKSampEq, M, pi/M, 'gray').';
+ %% The decoding method described in Taylor (2009)
+ %% works on the complex symbols, i.e. after taking
+ %% the nearest symbol in the constellation, but before
+ %% converting them back to integers/bits.
+ %% MATLAB's pskdemod() does not provide this intermediate
+ %% result, so to be lazy, a pskmod() call is performed
+ %% to obtain the complex symbols.
+ demodDEPSK = pskdemod(rDEPSKSampEq, M, 0, 'gray').';
+ remodDEPSK = pskmod(demodDEPSK, M, 0, 'gray');
+ delayed = [1; remodDEPSK(1:end-1)];
+ demodDEPSK = pskdemod(remodDEPSK .* conj(delayed), M, 0, 'gray');
+
+ demodDPSK = dpskdemod(rDPSKSamp, M, pi/M, 'gray');
+
+ [~, ber] = biterr(data, demodPSK);
+ berPSK(i) = berPSK(i) + ber / iterations;
+ [~, ber] = biterr(data, demodDEPSK);
+ berDEPSK(i) = berDEPSK(i) + ber / iterations;
+ [~, ber] = biterr(data, demodDPSK);
+ berDPSK(i) = berDPSK(i) + ber / iterations;
+
+ if EbN0_db(i) == 8 && it == 1
+ figure(1234);
+ plot(repelem(-phiestsPSK, 8));
+ hold on;
+ plot(pTxLoPSK);
+ legend('estimate', 'actual');
+ hold off;
+
+ figure(1);
+ scatterplot(rPSKSampEq);
+ title('rPSKSampEq');
end
- end
-
-
- figure(1);
- clf;
-
- %% Plot simulated results
- semilogy(EbN0_db, ber, 'r', 'LineWidth', 2);
- hold on;
- theoreticalPSK(EbN0_db, M, 'b', 'LineWidth', 1);
- legend({'Simulated phase noise', 'Without phase noise'}, 'Location', 'southwest');
+ end
+end
- title(strcat(num2str(M), '-PSK with phase noise and correction'));
- grid on;
- xlabel('$E_b/N_0$ (dB)');
- ylabel('BER');
- formatFigure;
-end
+figure(1);
+clf;
+
+%% Plot simulated results
+semilogy(EbN0_db, berPSK, 'r', 'LineWidth', 1.5);
+hold on;
+semilogy(EbN0_db, berDEPSK, 'c', 'LineWidth', 2);
+semilogy(EbN0_db, berDPSK, 'Color', [0, 0.6, 0], 'LineWidth', 2.5);
+
+theoreticalPSK(EbN0_db, M, 'b', 'LineWidth', 1);
+DEPSKTheoretical = berawgn(EbN0_db, 'psk', M, 'diff');
+semilogy(EbN0_db, DEPSKTheoretical, 'Color', [1, 0.6, 0], 'LineWidth', 1);
+DPSKTheoretical = berawgn(EbN0_db, 'dpsk', M);
+semilogy(EbN0_db, DPSKTheoretical, 'm', 'LineWidth', 1);
+
+legend({'PSK with Viterbi-Viterbi', ...
+ 'DEPSK with Viterbi-Viterbi', ...
+ 'DPSK', ...
+ 'Theoretical PSK over AWGN', ...
+ 'Theoretical DEPSK over AWGN', ...
+ 'Theoretical DPSK over AWGN'}, ...
+ 'Location', 'southwest');
+
+title({'QPSK with phase nosie and correction', ...
+ strcat('$10^{', num2str(log10(numSymbs * log2(M))), ...
+ '}$~bits, LO~', ...
+ num2str(linewidthLO / 1e6), '~MHz, blocksize~', ...
+ num2str(avgSa), '~Sa')});
+grid on;
+xlabel('$E_b/N_0$ (dB)');
+ylabel('BER');
+
+formatFigure;
projects / 4yp.git / blobdiff