From 8a784c714595170d79078bdde3d7ffe27c5510fe Mon Sep 17 00:00:00 2001
From: Pau Espin Pedrol <pespin@sysmocom.de>
Date: Tue, 30 Jul 2019 17:44:45 +0200
Subject: Move matlab files under utils/matlab

Change-Id: I15b687fbf436d662b264cb00f72b367ccd64b962
---
 utils/matlab/laurent.m          | 83 +++++++++++++++++++++++++++++++++++++++++
 utils/matlab/pulseApproximate.m | 15 ++++++++
 2 files changed, 98 insertions(+)
 create mode 100644 utils/matlab/laurent.m
 create mode 100644 utils/matlab/pulseApproximate.m

(limited to 'utils')

diff --git a/utils/matlab/laurent.m b/utils/matlab/laurent.m
new file mode 100644
index 0000000..ef15428
--- /dev/null
+++ b/utils/matlab/laurent.m
@@ -0,0 +1,83 @@
+%
+% Laurent decomposition of GMSK signals
+% Generates C0, C1, and C2 pulse shapes
+%
+% Pierre Laurent, "Exact and Approximate Construction of Digital Phase
+%   Modulations by Superposition of Amplitude Modulated Pulses", IEEE
+%   Transactions of Communications, Vol. 34, No. 2, Feb 1986.
+%
+% Author: Thomas Tsou <tom@tsou.cc>
+%
+
+% Modulation parameters
+oversamp = 16;
+L = 3;
+f = 270.83333e3;
+T = 1/f;
+h = 0.5;
+BT = 0.30;
+B = BT / T;
+
+% Generate sampling points for L symbol periods
+t = -(L*T/2):T/oversamp:(L*T/2);
+t = t(1:end-1) + (T/oversamp/2);
+
+% Generate Gaussian pulse
+g = qfunc(2*pi*B*(t - T/2)/(log(2)^.5)) - qfunc(2*pi*B*(t + T/2)/(log(2)^.5));
+g = g / sum(g) * pi/2;
+g = [0 g];
+
+% Integrate phase 
+q = 0;
+for i = 1:size(g,2);
+    q(i) = sum(g(1:i));
+end
+
+% Compute two sided "generalized phase pulse" function
+s = 0;
+for i = 1:size(g,2);
+    s(i) = sin(q(i)) / sin(pi*h);
+end
+for i = (size(g,2) + 1):(2 * size(g,2) - 1);
+    s(i) = sin(pi*h - q(i - (size(g,2) - 1))) / sin(pi*h);
+end
+
+% Compute C0 pulse: valid for all L values
+c0 = s(1:end-(oversamp*(L-1)));
+for i = 1:L-1;
+    c0 = c0 .* s((1 + i*oversamp):end-(oversamp*(L - 1 - i)));
+end
+
+% Compute C1 pulse: valid for L = 3 only!
+%   C1 = S0 * S4 * S2
+c1 = s(1:end-(oversamp*(4)));
+c1 = c1 .* s((1 + 4*oversamp):end-(oversamp*(4 - 1 - 3)));
+c1 = c1 .* s((1 + 2*oversamp):end-(oversamp*(4 - 1 - 1)));
+
+% Compute C2 pulse: valid for L = 3 only!
+%   C2 = S0 * S1 * S5
+c2 = s(1:end-(oversamp*(5)));
+c2 = c2 .* s((1 + 1*oversamp):end-(oversamp*(5 - 1 - 0)));
+c2 = c2 .* s((1 + 5*oversamp):end-(oversamp*(5 - 1 - 4)));
+
+% Plot C0, C1, C2 Laurent pulse series
+figure(1);
+hold off;
+plot((0:size(c0,2)-1)/oversamp - 2,c0, 'b');
+hold on;
+plot((0:size(c1,2)-1)/oversamp - 2,c1, 'r');
+plot((0:size(c2,2)-1)/oversamp - 2,c2, 'g');
+
+% Generate OpenBTS pulse
+numSamples = size(c0,2); 
+centerPoint = (numSamples - 1)/2;
+i = ((0:numSamples) - centerPoint) / oversamp;
+xP = .96*exp(-1.1380*i.^2 - 0.527*i.^4);
+xP = xP / max(xP) * max(c0);
+
+% Plot C0 pulse compared to OpenBTS pulse
+figure(2);
+hold off;
+plot((0:size(c0,2)-1)/oversamp, c0, 'b');
+hold on;
+plot((0:size(xP,2)-1)/oversamp, xP, 'r');
diff --git a/utils/matlab/pulseApproximate.m b/utils/matlab/pulseApproximate.m
new file mode 100644
index 0000000..2ff9234
--- /dev/null
+++ b/utils/matlab/pulseApproximate.m
@@ -0,0 +1,15 @@
+pp = [0 0 0.015 0.18 0.7 0.96 0.7 0.18 0.015 0 0];
+t = -2.5:0.5:2.5;
+
+v = -0.000:-0.001:-1.999;
+
+
+for ix1 = 1:length(v),
+	    disp(ix1);
+ for ix2 = 1:length(v),
+  p = exp(v(ix1)*t.^2+v(ix2)*t.^4);
+  r(ix1,ix2) = norm(p./max(abs(p)) - pp./max(abs(pp)));
+ end;
+end;
+
+
-- 
cgit v1.2.3