From 9ccd9f2c3c5fb69c79661568e08e20d7a7582ab5 Mon Sep 17 00:00:00 2001 From: Thomas Tsou Date: Wed, 21 Aug 2013 13:59:52 -0400 Subject: Transceiver52M: Add 4 samples-per-symbol Laurent pulse shape When 4 samples-per-symbol operation is selected, replace the existing pulse approximation, which becomes inaccurate with non-unit oversampling, with the primary pulse, C0, from the Laurent linear pulse approximation. 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. Octave pulse generation code for the first three pulses of the linear approximation are included. Signed-off-by: Thomas Tsou --- Transceiver52M/laurent.m | 83 +++++++++++++++++++++++++++++++++++++++++++ Transceiver52M/sigProcLib.cpp | 49 ++++++++++++++++++++----- 2 files changed, 123 insertions(+), 9 deletions(-) create mode 100644 Transceiver52M/laurent.m diff --git a/Transceiver52M/laurent.m b/Transceiver52M/laurent.m new file mode 100644 index 0000000..ef15428 --- /dev/null +++ b/Transceiver52M/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 +% + +% 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/Transceiver52M/sigProcLib.cpp b/Transceiver52M/sigProcLib.cpp index b0ef21a..51584b7 100644 --- a/Transceiver52M/sigProcLib.cpp +++ b/Transceiver52M/sigProcLib.cpp @@ -422,11 +422,22 @@ void generateGSMPulse(int sps, int symbolLength) GSMPulse->empty->isRealOnly(true); *(GSMPulse->empty->begin()) = 1.0f; - len = sps * symbolLength; - if (len < 4) - len = 4; + /* + * For 4 samples-per-symbol use a precomputed single pulse Laurent + * approximation. This should yields below 2 degrees of phase error at + * the modulator output. Use the existing pulse approximation for all + * other oversampling factors. + */ + switch (sps) { + case 4: + len = 16; + break; + default: + len = sps * symbolLength; + if (len < 4) + len = 4; + } - /* GSM pulse approximation */ GSMPulse->buffer = convolve_h_alloc(len); GSMPulse->gaussian = new signalVector((complex *) GSMPulse->buffer, 0, len); @@ -435,12 +446,32 @@ void generateGSMPulse(int sps, int symbolLength) signalVector::iterator xP = GSMPulse->gaussian->begin(); - center = (float) (len - 1.0) / 2.0; + if (sps == 4) { + *xP++ = 4.46348606e-03; + *xP++ = 2.84385729e-02; + *xP++ = 1.03184855e-01; + *xP++ = 2.56065552e-01; + *xP++ = 4.76375085e-01; + *xP++ = 7.05961177e-01; + *xP++ = 8.71291644e-01; + *xP++ = 9.29453645e-01; + *xP++ = 8.71291644e-01; + *xP++ = 7.05961177e-01; + *xP++ = 4.76375085e-01; + *xP++ = 2.56065552e-01; + *xP++ = 1.03184855e-01; + *xP++ = 2.84385729e-02; + *xP++ = 4.46348606e-03; + *xP++ = 0.0; + } else { + center = (float) (len - 1.0) / 2.0; - for (int i = 0; i < len; i++) { - arg = ((float) i - center) / (float) sps; - *xP++ = 0.96 * exp(-1.1380 * arg * arg - - 0.527 * arg * arg * arg * arg); + /* GSM pulse approximation */ + for (int i = 0; i < len; i++) { + arg = ((float) i - center) / (float) sps; + *xP++ = 0.96 * exp(-1.1380 * arg * arg - + 0.527 * arg * arg * arg * arg); + } } float avgAbsval = sqrtf(vectorNorm2(*GSMPulse->gaussian)/sps); -- cgit v1.2.3