Added baseband and passband simulations. Passband is buggy.

diff --git a/RRC_PSK_BER_SNR.m b/RRC_PSK_BER_SNR.m
deleted file mode 100644 (file)
index 73d67f3..0000000
--- a/RRC_PSK_BER_SNR.m
+++ /dev/null
@@ -1,94 +0,0 @@
-function RRC_PSK_BER_SNR(rolloff, M, numSymbs)
-  %% Set defaults for inputs
-  if nargin < 3
-    numSymbs = 1000;
-  end
-  if nargin < 2
-    M = 2;
-  end
-  if nargin < 1
-    rolloff = 0.5;
-  end
-
-  if isOctave()
-    pkg load communications
-  end
-
-  %% https://www.mathworks.com/help/signal/ref/rcosdesign.html
-  %% https://www.mathworks.com/help/comm/ug/pulse-shaping-using-a-raised-cosine-filter.html
-  span = 6; % filter span
-  sps = 4;
-
-  rrcFilter = rcosdesign(rolloff, span, sps, 'sqrt');
-
-  EbN0_db = 0:0.2:10;
-  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);
-
-  data = randi([0 M - 1], numSymbs, 1);
-  modData = pskmod(data, M, 0, 'gray');
-
-  txSig = upfirdn(modData, rrcFilter, sps);
-
-  for i = 1:plotlen
-    snr = EbN0_db(i) + 10 * log10(log2(M));% - 10 * log10(sps); % why sps?
-    rxSig = awgn(txSig, snr);
-
-    rxFilt = upfirdn(rxSig, rrcFilter, 1, sps);
-    rxFilt = rxFilt(span + 1 : end - span); % remove filter delay
-
-    demodData = pskdemod(rxFilt, M, 0, 'gray');
-
-    [bitErrors, ber(i)] = biterr(data, demodData);
-  end
-
-  fig1 = figure(1);
-  clf;
-
-  %% Plot simulated results
-  semilogy(EbN0_db, ber, 'r', 'LineWidth', 2);
-  hold on;
-
-  %% Plot theoretical curve
-  %% BPSK: bit error when noise Nr > sqrt(Eb)
-  %%   Pr(Nr > sqrt(Eb))
-  %% = Pr(Z > sqrt(Eb) / sqrt(N0/2))
-  %%
-  %% QPSK = 2 BPSKs, one real and one imaginary, each with one bit
-  %% so BER is the same as BPSK (assuming Gray code)
-  if M == 2 || M == 4
-    ber_th = qfunc(sqrt(2 * EbN0));
-    semilogy(EbN0_db, ber_th, 'b', 'LineWidth', 1);
-    legend('Simulated RRC', 'Discrete');
-  else
-    %% Approximation: J.G. Proakis and M. Salehi, 2000, Contemporary
-    %%                Communication Systems using MATLAB (Equations
-    %%                7.3.18 and 7.3.19), Brooks/Cole.
-    ber_ap = 2 * qfunc(sqrt(EbN0 * log2(M) * 2) * sin(pi / M)) / log2(M);
-    semilogy(EbN0_db, ber_ap, 'b', 'LineWidth', 1);
-    legend('Simulated RRC', 'Discrete');
-  end
-
-  title(strcat(num2str(M), '-PSK RRC with Gray code'));
-  grid on;
-  xlabel('$E_b/N_0$ (dB)');
-  ylabel('BER');
-
-  formatFigure;
-  %saveas(gcf, strcat('BER_SNR_', num2str(M), 'PSK_', num2str(numSymbs), ...
-  %                   '.svg'));
-
-  %scatterplot(rxFilt);
-  %eyediagram(rxFilt, sps);
-
-end

diff --git a/baseband.m b/baseband.m
new file mode 100644 (file)
index 0000000..faf37f9
--- /dev/null
+++ b/baseband.m
@@ -0,0 +1,115 @@
+function baseband(rolloff, M, numSymbs)
+  %% Set defaults for inputs
+  if nargin < 3
+    numSymbs = 1000;
+  end
+  if nargin < 2
+    M = 2;
+  end
+  if nargin < 1
+    rolloff = 0.5;
+  end
+
+
+  %% https://www.mathworks.com/help/comm/examples/passband-modulation-with-adjacent-channel-interference.html
+  Rsym = 1e6; % symbol rate (sym/sec)
+
+  span = 6; % filter span
+  sps = 4; % samples per symbol
+
+  txFilter = comm.RaisedCosineTransmitFilter...
+                 ('Shape', 'Square root', ...
+                  'RolloffFactor', rolloff, ...
+                  'FilterSpanInSymbols', span, ...
+                  'OutputSamplesPerSymbol', sps);
+  rxFilter = comm.RaisedCosineReceiveFilter...
+                 ('Shape', 'Square root', ...
+                  'RolloffFactor', rolloff, ...
+                  'FilterSpanInSymbols', span, ...
+                  'InputSamplesPerSymbol', sps, ...
+                  'DecimationFactor', 1);
+
+  fs = Rsym * sps; % sampling freq (Hz)
+
+  t = (0 : 1 / fs : numSymbs / Rsym + (1.5 * span * sps - 1) / fs)';
+
+
+  EbN0_db = 0:0.2:10;
+  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);
+
+  data = randi([0 M - 1], numSymbs, 1);
+  modData = pskmod(data, M, 0, 'gray');
+
+  xBaseband = txFilter([modData; zeros(span, 1)]);
+
+
+
+  for i = 1:plotlen
+    snr = EbN0_db(i) + 10 * log10(log2(M)) - 10 * log10(sps); % why sps?
+    noiseEnergy = 10 ^ (-snr / 10);
+
+    yBaseband = awgn(xBaseband, snr, 'measured');
+
+    rBaseband = rxFilter([yBaseband; zeros(span, 1)]);
+    %% truncate filter transients
+    rBaseband = rBaseband(span * sps / 2 + 1 : end);
+    %% normalize to unit energy
+    rBasebandEnergy = sum(abs(rBaseband) .^ 2) / numSymbs;
+    rBaseband = rBaseband .* sqrt((1 + noiseEnergy) / rBasebandEnergy);
+
+    rSampled = rBaseband(sps*span/2+1:sps:(numSymbs+span/2)*sps);
+
+    demodData = pskdemod(rSampled, M, 0, 'gray');
+
+    [bitErrors, ber(i)] = biterr(data, demodData);
+  end
+
+  fig1 = figure(1);
+  clf;
+
+  %% Plot simulated results
+  semilogy(EbN0_db, ber, 'r', 'LineWidth', 2);
+  hold on;
+
+  %% Plot theoretical curve
+  %% BPSK: bit error when noise Nr > sqrt(Eb)
+  %%   Pr(Nr > sqrt(Eb))
+  %% = Pr(Z > sqrt(Eb) / sqrt(N0/2))
+  %%
+  %% QPSK = 2 BPSKs, one real and one imaginary, each with one bit
+  %% so BER is the same as BPSK (assuming Gray code)
+  if M == 2 || M == 4
+    ber_th = qfunc(sqrt(2 * EbN0));
+    semilogy(EbN0_db, ber_th, 'b', 'LineWidth', 1);
+    legend('Simulated', 'Discrete');
+  else
+    %% Approximation: J.G. Proakis and M. Salehi, 2000, Contemporary
+    %%                Communication Systems using MATLAB (Equations
+    %%                7.3.18 and 7.3.19), Brooks/Cole.
+    ber_ap = 2 * qfunc(sqrt(EbN0 * log2(M) * 2) * sin(pi / M)) / log2(M);
+    semilogy(EbN0_db, ber_ap, 'b', 'LineWidth', 1);
+    legend('Simulated', 'Discrete');
+  end
+
+  title(strcat(num2str(M), '-PSK with Gray code'));
+  grid on;
+  xlabel('$E_b/N_0$ (dB)');
+  ylabel('BER');
+
+  formatFigure;
+  %saveas(gcf, strcat('BER_SNR_', num2str(M), 'PSK_', num2str(numSymbs), ...
+  %                   '.svg'));
+
+  %scatterplot(rxFilt);
+  %eyediagram(rxFilt, sps);
+
+end

diff --git a/passband.m b/passband.m
new file mode 100644 (file)
index 0000000..7abd749
--- /dev/null
+++ b/passband.m
@@ -0,0 +1,151 @@
+function passband(rolloff, M, numSymbs)
+  %% Set defaults for inputs
+  if nargin < 3
+    numSymbs = 1000;
+  end
+  if nargin < 2
+    M = 2;
+  end
+  if nargin < 1
+    rolloff = 0.5;
+  end
+
+
+  %% https://www.mathworks.com/help/comm/examples/passband-modulation-with-adjacent-channel-interference.html
+  Rsym = 1e6; % symbol rate (sym/sec)
+
+  span = 6; % filter span
+  sps = 4; % samples per symbol
+
+  txFilter = comm.RaisedCosineTransmitFilter...
+                 ('Shape', 'Square root', ...
+                  'RolloffFactor', rolloff, ...
+                  'FilterSpanInSymbols', span, ...
+                  'OutputSamplesPerSymbol', sps);
+  rxFilter = comm.RaisedCosineReceiveFilter...
+                 ('Shape', 'Square root', ...
+                  'RolloffFactor', rolloff, ...
+                  'FilterSpanInSymbols', span, ...
+                  'InputSamplesPerSymbol', sps, ...
+                  'DecimationFactor', 1);
+
+  fs = Rsym * sps; % sampling freq (Hz)
+
+  t = (0 : 1 / fs : numSymbs / Rsym + (1.5 * span * sps - 1) / fs)';
+
+
+
+
+  EbN0_db = 0:0.2:10;
+  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);
+
+
+
+
+  data = randi([0 M - 1], numSymbs, 1);
+  modData = pskmod(data, M, 0, 'gray');
+
+  xBaseband = txFilter([modData; zeros(span, 1)]);
+
+  %fc = 2.5e6; % Carrier freq (Hz)
+  %carrier = sqrt(2) * exp(j * 2 * pi * fc * t);
+
+  %xPassbandIdeal = normalizeEnergy...
+  %                   (real(xBaseband .* carrier(1:length(xBaseband))), numSymbs, 1);
+
+  txLOFreq = [2.49e6, 2.5e6, 2.51e6];
+  %%txLOEnergy = [0.05, 0.9, 0.05];
+  txLOEnergy = [0 1 0];
+
+  carrier = zeros(length(t), 1);
+  for i = 1 : length(txLOFreq)
+    carrier = carrier + ...
+              sqrt(2 * txLOEnergy(i)) * exp(j * 2 * pi * txLOFreq(i) * t);
+  end
+
+  xPassband = normalizeEnergy...
+                (real(xBaseband .* carrier(1:length(xBaseband))), numSymbs, 1);
+
+  sum(abs(xPassband) .^ 2) / numSymbs
+  input('pause')
+
+
+  for i = 1:plotlen
+    snr = EbN0_db(i) + 10 * log10(log2(M)) - 10 * log10(sps); % why sps?
+    noiseEnergy = 10 ^ (-snr / 10);
+
+
+    yPassband = awgn(xPassband, snr, 'measured');
+
+
+    rBaseband = rxFilter([yPassband .* carrier(1:length(yPassband)); zeros(span * sps, 1)]);
+    %% truncate filter transients
+    rBaseband = rBaseband(span * sps / 2 + 1 : end);
+    %% normalize energy
+    rBaseband = normalizeEnergy(rBaseband, numSymbs, 1 + noiseEnergy);
+
+
+    rSampled = rBaseband(sps*span/2+1:sps:(numSymbs + span/2) * sps);
+
+    demodData = pskdemod(rSampled, M, 0, 'gray');
+    [bitErrors, ber(i)] = biterr(data, demodData);
+
+  end
+
+  fig1 = figure(1);
+  clf;
+
+  %% Plot simulated results
+  semilogy(EbN0_db, ber, 'r', 'LineWidth', 2);
+  hold on;
+
+  %% Plot theoretical curve
+  %% BPSK: bit error when noise Nr > sqrt(Eb)
+  %%   Pr(Nr > sqrt(Eb))
+  %% = Pr(Z > sqrt(Eb) / sqrt(N0/2))
+  %%
+  %% QPSK = 2 BPSKs, one real and one imaginary, each with one bit
+  %% so BER is the same as BPSK (assuming Gray code)
+  if M == 2 || M == 4
+    ber_th = qfunc(sqrt(2 * EbN0));
+    semilogy(EbN0_db, ber_th, 'b', 'LineWidth', 1);
+    legend('Simulated RRC', 'Discrete');
+  else
+    %% Approximation: J.G. Proakis and M. Salehi, 2000, Contemporary
+    %%                Communication Systems using MATLAB (Equations
+    %%                7.3.18 and 7.3.19), Brooks/Cole.
+    ber_ap = 2 * qfunc(sqrt(EbN0 * log2(M) * 2) * sin(pi / M)) / log2(M);
+    semilogy(EbN0_db, ber_ap, 'b', 'LineWidth', 1);
+    legend('Simulated RRC', 'Discrete');
+  end
+
+  title(strcat(num2str(M), '-PSK RRC with Gray code'));
+  grid on;
+  xlabel('$E_b/N_0$ (dB)');
+  ylabel('BER');
+
+  formatFigure;
+  %saveas(gcf, strcat('BER_SNR_', num2str(M), 'PSK_', num2str(numSymbs), ...
+  %                   '.svg'));
+
+  %scatterplot(rxFilt);
+  %eyediagram(rxFilt, sps);
+
+end
+
+
+function y = normalizeEnergy(x, numSymbs, e)
+  energy = sum(abs(x) .^ 2) / numSymbs;
+  y = x * sqrt(e / energy);
+end