From 5155713e970946e2ce213c4f68ccc44e6737ead0 Mon Sep 17 00:00:00 2001 From: Martin Braun Date: Wed, 8 Dec 2010 19:32:06 +0100 Subject: CPM make checks --- gnuradio-core/src/lib/general/gr_cpm.cc | 211 ++++++++++++++++++++++++++++++++ 1 file changed, 211 insertions(+) create mode 100644 gnuradio-core/src/lib/general/gr_cpm.cc (limited to 'gnuradio-core/src/lib/general/gr_cpm.cc') diff --git a/gnuradio-core/src/lib/general/gr_cpm.cc b/gnuradio-core/src/lib/general/gr_cpm.cc new file mode 100644 index 000000000..a5d328edf --- /dev/null +++ b/gnuradio-core/src/lib/general/gr_cpm.cc @@ -0,0 +1,211 @@ +/* -*- c++ -*- */ +/* + * Copyright 2010 Free Software Foundation, Inc. + * + * GNU Radio is free software; you can redistribute it and/or modify + * it under the terms of the GNU General Public License as published by + * the Free Software Foundation; either version 3, or (at your option) + * any later version. + * + * GNU Radio is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with GNU Radio; see the file COPYING. If not, write to + * the Free Software Foundation, Inc., 51 Franklin Street, + * Boston, MA 02110-1301, USA. + */ + +// Calculate the taps for the CPM phase responses + +#ifdef HAVE_CONFIG_H +#include "config.h" +#endif + +#include +#include +#include + + +//! Normalised sinc function, sinc(x)=sin(pi*x)/pi*x +inline double +sinc(double x) +{ + if (x == 0) { + return 1.0; + } + + return sin(M_PI * x) / (M_PI * x); +} + + +//! Taps for L-RC CPM (Raised cosine of length L symbols) +std::vector +generate_cpm_lrc_taps(unsigned samples_per_sym, unsigned L) +{ + std::vector taps(samples_per_sym * L, 1.0/L/samples_per_sym); + for (unsigned i = 0; i < samples_per_sym * L; i++) { + taps[i] *= 1 - cos(M_TWOPI * i / L / samples_per_sym); + } + + return taps; +} + + +/*! Taps for L-SRC CPM (Spectral raised cosine of length L symbols). + * + * L-SRC has a time-continuous phase response function of + * + * g(t) = 1/LT * sinc(2t/LT) * cos(beta * 2pi t / LT) / (1 - (4beta / LT * t)^2) + * + * which is the Fourier transform of a cos-rolloff function with rolloff + * beta, and looks like a sinc-function, multiplied with a rolloff term. + * We return the main lobe of the sinc, i.e., everything between the + * zero crossings. + * The time-discrete IR is thus + * + * g(k) = 1/Ls * sinc(2k/Ls) * cos(beta * pi k / Ls) / (1 - (4beta / Ls * k)^2) + * where k = 0...Ls-1 + * and s = samples per symbol. + */ +std::vector +generate_cpm_lsrc_taps(unsigned samples_per_sym, unsigned L, double beta) +{ + double Ls = (double) L * samples_per_sym; + std::vector taps_d(L * samples_per_sym, 0.0); + std::vector taps(L * samples_per_sym, 0.0); + + double sum = 0; + for (unsigned i = 0; i < samples_per_sym * L; i++) { + double k = i - Ls/2; // Causal to acausal + + taps_d[i] = 1.0 / Ls * sinc(2.0 * k / Ls); + + // For k = +/-Ls/4*beta, the rolloff term's cos-function becomes zero + // and the whole thing converges to PI/4 (to prove this, use de + // l'hopital's rule). + if (fabs(fabs(k) - Ls/4/beta) < 2*DBL_EPSILON) { + taps_d[i] *= M_PI_4; + } else { + double tmp = 4.0 * beta * k / Ls; + taps_d[i] *= cos(beta * M_TWOPI * k / Ls) / (1 - tmp * tmp); + } + sum += taps_d[i]; + } + for (unsigned i = 0; i < samples_per_sym * L; i++) { + taps[i] = (float) taps_d[i] / sum; + } + + return taps; +} + + +//! Taps for L-REC CPM (Rectangular pulse shape of length L symbols) +std::vector +generate_cpm_lrec_taps(unsigned samples_per_sym, unsigned L) +{ + return std::vector(samples_per_sym * L, 1.0/L/samples_per_sym); +} + + +//! Helper function for TFM +double tfm_g0(double k, double sps) +{ + if (k < 2 * DBL_EPSILON) { + return 1.145393004159143; // 1 + pi^2/48 / sqrt(2) + } + + const double pi2_24 = 0.411233516712057; // pi^2/24 + double f = M_PI * k / sps; + return sinc(k/sps) - pi2_24 * (2 * sin(f) - 2*f*cos(f) - f*f*sin(f)) / (f*f*f); +} + + +//! Taps for TFM CPM (Tamed frequency modulation) +// +// See [2, Chapter 2.7.2]. +// +// [2]: Anderson, Aulin and Sundberg; Digital Phase Modulation +std::vector +generate_cpm_tfm_taps(unsigned sps, unsigned L) +{ + double causal_shift = (double) L * sps / 2; + std::vector taps_d(sps * L, 0.0); + std::vector taps(sps * L, 0.0); + + double sum = 0; + for (unsigned i = 0; i < sps * L; i++) { + double k = (double)i - causal_shift; // Causal to acausal + + taps_d[i] = tfm_g0(k - sps, sps) + + 2 * tfm_g0(k, sps) + + tfm_g0(k + sps, sps); + sum += taps_d[i]; + } + for (unsigned i = 0; i < sps * L; i++) { + taps[i] = (float) taps_d[i] / sum; + } + + return taps; +} + + +//! Taps for Gaussian CPM. Phase response is truncated after \p L symbols. +// \p bt sets the 3dB-time-bandwidth product. +// +// Note: for h = 0.5, this is the phase response for GMSK. +// +// This C99-compatible formula for the taps is taken straight +// from [1, Chapter 9.2.3]. +// A version in Q-notation can be found in [2, Chapter 2.7.2]. +// +// [1]: Karl-Dirk Kammeyer; Nachrichtenübertragung, 4th Edition. +// [2]: Anderson, Aulin and Sundberg; Digital Phase Modulation +// +std::vector +generate_cpm_gaussian_taps(unsigned samples_per_sym, unsigned L, double bt) +{ + double Ls = (double) L * samples_per_sym; + std::vector taps_d(L * samples_per_sym, 0.0); + std::vector taps(L * samples_per_sym, 0.0); + + // alpha = sqrt(2/ln(2)) * pi * BT + double alpha = 5.336446256636997 * bt; + for (unsigned i = 0; i < samples_per_sym * L; i++) { + double k = i - Ls/2; // Causal to acausal + taps_d[i] = (erf(alpha * (k / samples_per_sym + 0.5)) - + erf(alpha * (k / samples_per_sym - 0.5))) + * 0.5 / samples_per_sym; + taps[i] = (float) taps_d[i]; + } + + return taps; +} + + +std::vector +gr_cpm::phase_response(cpm_type type, unsigned samples_per_sym, unsigned L, double beta) +{ + switch (type) { + case LRC: + return generate_cpm_lrc_taps(samples_per_sym, L); + + case LSRC: + return generate_cpm_lsrc_taps(samples_per_sym, L, beta); + + case LREC: + return generate_cpm_lrec_taps(samples_per_sym, L); + + case TFM: + return generate_cpm_tfm_taps(samples_per_sym, L); + + case GAUSSIAN: + return generate_cpm_gaussian_taps(samples_per_sym, L, beta); + + default: + return generate_cpm_lrec_taps(samples_per_sym, 1); + } +} + -- cgit From dd35f129cc0ecbbf092a25a60414fef133e6dc04 Mon Sep 17 00:00:00 2001 From: Martin Braun Date: Thu, 9 Dec 2010 09:54:42 +0100 Subject: more elaborate checks, they also work now --- gnuradio-core/src/lib/general/gr_cpm.cc | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) (limited to 'gnuradio-core/src/lib/general/gr_cpm.cc') diff --git a/gnuradio-core/src/lib/general/gr_cpm.cc b/gnuradio-core/src/lib/general/gr_cpm.cc index a5d328edf..f2d032a22 100644 --- a/gnuradio-core/src/lib/general/gr_cpm.cc +++ b/gnuradio-core/src/lib/general/gr_cpm.cc @@ -113,7 +113,7 @@ generate_cpm_lrec_taps(unsigned samples_per_sym, unsigned L) //! Helper function for TFM double tfm_g0(double k, double sps) { - if (k < 2 * DBL_EPSILON) { + if (fabs(k) < 2 * DBL_EPSILON) { return 1.145393004159143; // 1 + pi^2/48 / sqrt(2) } @@ -131,13 +131,13 @@ double tfm_g0(double k, double sps) std::vector generate_cpm_tfm_taps(unsigned sps, unsigned L) { - double causal_shift = (double) L * sps / 2; + unsigned causal_shift = sps * L / 2; std::vector taps_d(sps * L, 0.0); std::vector taps(sps * L, 0.0); double sum = 0; for (unsigned i = 0; i < sps * L; i++) { - double k = (double)i - causal_shift; // Causal to acausal + double k = (double)(((int)i) - ((int)causal_shift)); // Causal to acausal taps_d[i] = tfm_g0(k - sps, sps) + 2 * tfm_g0(k, sps) + -- cgit From b12498643aa5c11a35a484925c565a7a9e746f75 Mon Sep 17 00:00:00 2001 From: Martin Braun Date: Fri, 10 Dec 2010 17:30:01 +0100 Subject: fixed: FM sensitivity and calling gr_cpm::phase_response() through SWIG --- gnuradio-core/src/lib/general/gr_cpm.cc | 3 +++ 1 file changed, 3 insertions(+) (limited to 'gnuradio-core/src/lib/general/gr_cpm.cc') diff --git a/gnuradio-core/src/lib/general/gr_cpm.cc b/gnuradio-core/src/lib/general/gr_cpm.cc index f2d032a22..a00526b52 100644 --- a/gnuradio-core/src/lib/general/gr_cpm.cc +++ b/gnuradio-core/src/lib/general/gr_cpm.cc @@ -28,6 +28,9 @@ #include #include +#ifndef M_TWOPI +# define M_TWOPI (2*M_PI) +#endif //! Normalised sinc function, sinc(x)=sin(pi*x)/pi*x inline double -- cgit