aboutsummaryrefslogtreecommitdiffstats
path: root/lib
diff options
context:
space:
mode:
authorpiotr <Piotr Krysik pkrysik@elka.pw.edu.pl>2014-02-04 17:57:25 +0100
committerpiotr <Piotr Krysik pkrysik@elka.pw.edu.pl>2014-02-04 17:57:25 +0100
commit437f5467a12ceddbc93a88f704697bd024378959 (patch)
tree68fb7a5ca2d4e72358c3444156072392a9fe627f /lib
Initial commit - gsm-receiver with removed quick hacks
Diffstat (limited to 'lib')
-rw-r--r--lib/CMakeLists.txt71
-rw-r--r--lib/assert.h67
-rw-r--r--lib/gsm_constants.h154
-rw-r--r--lib/receiver_config.cc84
-rw-r--r--lib/receiver_config.h164
-rw-r--r--lib/receiver_impl.cc770
-rw-r--r--lib/receiver_impl.h215
-rw-r--r--lib/sch.c333
-rw-r--r--lib/sch.h19
-rw-r--r--lib/system.h10
-rw-r--r--lib/test_gsm.cc47
-rw-r--r--lib/viterbi_detector.cc554
-rw-r--r--lib/viterbi_detector.h63
13 files changed, 2551 insertions, 0 deletions
diff --git a/lib/CMakeLists.txt b/lib/CMakeLists.txt
new file mode 100644
index 0000000..8d4c721
--- /dev/null
+++ b/lib/CMakeLists.txt
@@ -0,0 +1,71 @@
+# Copyright 2011,2012 Free Software Foundation, Inc.
+#
+# This file is part of GNU Radio
+#
+# 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.
+
+########################################################################
+# Setup library
+########################################################################
+include(GrPlatform) #define LIB_SUFFIX
+
+include_directories(${Boost_INCLUDE_DIR})
+link_directories(${Boost_LIBRARY_DIRS})
+
+list(APPEND gsm_sources
+ receiver_impl.cc
+ receiver_config.cc
+ viterbi_detector.cc
+ sch.c
+)
+
+add_library(gnuradio-gsm SHARED ${gsm_sources})
+target_link_libraries(gnuradio-gsm ${Boost_LIBRARIES} ${GNURADIO_RUNTIME_LIBRARIES})
+set_target_properties(gnuradio-gsm PROPERTIES DEFINE_SYMBOL "gnuradio_gsm_EXPORTS")
+
+########################################################################
+# Install built library files
+########################################################################
+install(TARGETS gnuradio-gsm
+ LIBRARY DESTINATION lib${LIB_SUFFIX} # .so/.dylib file
+ ARCHIVE DESTINATION lib${LIB_SUFFIX} # .lib file
+ RUNTIME DESTINATION bin # .dll file
+)
+
+########################################################################
+# Build and register unit test
+########################################################################
+include(GrTest)
+
+include_directories(${CPPUNIT_INCLUDE_DIRS})
+
+list(APPEND test_gsm_sources
+ ${CMAKE_CURRENT_SOURCE_DIR}/test_gsm.cc
+ ${CMAKE_CURRENT_SOURCE_DIR}/qa_gsm.cc
+ ${CMAKE_CURRENT_SOURCE_DIR}/qa_receiver.cc
+)
+
+add_executable(test-gsm ${test_gsm_sources})
+
+target_link_libraries(
+ test-gsm
+ ${GNURADIO_RUNTIME_LIBRARIES}
+ ${Boost_LIBRARIES}
+ ${CPPUNIT_LIBRARIES}
+ gnuradio-gsm
+)
+
+GR_ADD_TEST(test_gsm test-gsm)
diff --git a/lib/assert.h b/lib/assert.h
new file mode 100644
index 0000000..a80ae73
--- /dev/null
+++ b/lib/assert.h
@@ -0,0 +1,67 @@
+/*
+* Copyright 2007 Free Software Foundation, Inc.
+*
+* This software is distributed under the terms of the GNU Public License.
+* See the COPYING file in the main directory for details.
+
+ This program 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 of the License, or
+ (at your option) any later version.
+
+ This program 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 this program. If not, see <http://www.gnu.org/licenses/>.
+
+*/
+
+
+#ifndef ASSERT_H
+#define ASSERT_H
+
+#include "stdio.h"
+#include <iostream>
+
+//#define NDEBUG
+
+/**@name Macros for standard messages. */
+//@{
+#define COUT(text) { std::cout << text << std::endl; }
+#define CERR(text) { std::cerr << __FILE__ << ":" << __LINE__ << ": " << text; }
+#ifdef NDEBUG
+#define DCOUT(text) {}
+#define OBJDCOUT(text) {}
+#else
+#define DCOUT(text) { COUT(__FILE__ << ":" << __LINE__ << " " << text); }
+#define OBJDCOUT(text) { DCOUT(this << " " << text); }
+#endif
+//@}
+
+
+/** This is thrown by assert() so that gdb can catch it. */
+
+class assertion
+{
+
+ public:
+
+ assertion() {
+ fprintf( stderr,"Try setting a breakpoint @ %s:%u.\n",__FILE__,__LINE__ );
+ return;
+ }
+
+};
+
+#ifdef NDEBUG
+#define assert(EXPR) {};
+#else
+/** This replaces the regular assert() with a C++ exception. */
+#include "stdio.h"
+#define assert(E) { if (!(E)) { fprintf(stderr,"%s:%u failed assertion '%s'\n",__FILE__,__LINE__,#E); throw Assertion(); } }
+#endif
+
+#endif
diff --git a/lib/gsm_constants.h b/lib/gsm_constants.h
new file mode 100644
index 0000000..11913f1
--- /dev/null
+++ b/lib/gsm_constants.h
@@ -0,0 +1,154 @@
+#ifndef INCLUDED_GSM_CONSTANTS_H
+#define INCLUDED_GSM_CONSTANTS_H
+
+#define GSM_SYMBOL_RATE (1625000.0/6.0) //symbols per second
+#define GSM_SYMBOL_PERIOD (1.0/GSM_SYMBOL_RATE) //seconds per symbol
+
+//Burst timing
+#define TAIL_BITS 3
+#define GUARD_BITS 8
+#define GUARD_FRACTIONAL 0.25 //fractional part of guard period
+#define GUARD_PERIOD GUARD_BITS + GUARD_FRACTIONAL
+#define DATA_BITS 57 //size of 1 data block in normal burst
+#define STEALING_BIT 1
+#define N_TRAIN_BITS 26
+#define N_SYNC_BITS 64
+#define USEFUL_BITS 142 //(2*(DATA_BITS+STEALING_BIT) + N_TRAIN_BITS )
+#define FCCH_BITS USEFUL_BITS
+#define BURST_SIZE (USEFUL_BITS+2*TAIL_BITS)
+
+#define SCH_DATA_LEN 39
+#define TS_BITS (TAIL_BITS+USEFUL_BITS+TAIL_BITS+GUARD_BITS) //a full TS (156 bits)
+#define TS_PER_FRAME 8
+#define FRAME_BITS (TS_PER_FRAME * TS_BITS + 2) // 156.25 * 8
+#define FCCH_POS TAIL_BITS
+#define SYNC_POS 39
+#define TRAIN_POS ( TAIL_BITS + (DATA_BITS+STEALING_BIT) + 5) //first 5 bits of a training sequence
+ //aren't used for channel impulse response estimation
+#define TRAIN_BEGINNING 5
+#define SAFETY_MARGIN 6 //
+
+#define FCCH_HITS_NEEDED (USEFUL_BITS - 4)
+#define FCCH_MAX_MISSES 1
+#define FCCH_MAX_FREQ_OFFSET 100
+
+#define CHAN_IMP_RESP_LENGTH 5
+
+#define MAX_SCH_ERRORS 5 //maximum number of subsequent sch errors after which gsm receiver goes to find_next_fcch state
+
+typedef enum {empty, fcch_burst, sch_burst, normal_burst, rach_burst, dummy, dummy_or_normal} burst_type;
+typedef enum {unknown, multiframe_26, multiframe_51} multiframe_type;
+
+static const unsigned char SYNC_BITS[] = {
+ 1, 0, 1, 1, 1, 0, 0, 1, 0, 1, 1, 0, 0, 0, 1, 0,
+ 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1,
+ 0, 0, 1, 0, 1, 1, 0, 1, 0, 1, 0, 0, 0, 1, 0, 1,
+ 0, 1, 1, 1, 0, 1, 1, 0, 0, 0, 0, 1, 1, 0, 1, 1
+};
+
+const unsigned FCCH_FRAMES[] = {0, 10, 20, 30, 40};
+const unsigned SCH_FRAMES[] = {1, 11, 21, 31, 41};
+
+const unsigned BCCH_FRAMES[] = {2, 3, 4, 5}; //!!the receiver shouldn't care about logical
+ //!!channels so this will be removed from this header
+const unsigned TEST_CCH_FRAMES[] = {2, 3, 4, 5, 6, 7, 8, 9, 12, 13, 14, 15, 16, 17, 18, 19, 22, 23, 24, 25, 26, 27, 28, 29, 32, 33, 34, 35, 36, 37, 38, 39, 42, 43, 44, 45, 46, 47, 48, 49};
+const unsigned TRAFFIC_CHANNEL_F[] = {0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24};
+const unsigned TEST51[] = {0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33, 34, 35, 36, 37, 38, 39, 40, 41, 42, 43, 44, 45, 46, 47, 48, 49, 50};
+
+
+#define TSC0 0
+#define TSC1 1
+#define TSC2 2
+#define TSC3 3
+#define TSC4 4
+#define TSC5 5
+#define TSC6 6
+#define TSC7 7
+#define TS_DUMMY 8
+
+#define TRAIN_SEQ_NUM 9
+
+#define TIMESLOT0 0
+#define TIMESLOT1 1
+#define TIMESLOT2 2
+#define TIMESLOT3 3
+#define TIMESLOT4 4
+#define TIMESLOT5 5
+#define TIMESLOT6 6
+#define TIMESLOT7 7
+
+
+static const unsigned char train_seq[TRAIN_SEQ_NUM][N_TRAIN_BITS] = {
+ {0, 0, 1, 0, 0, 1, 0, 1, 1, 1, 0, 0, 0, 0, 1, 0, 0, 0, 1, 0, 0, 1, 0, 1, 1, 1},
+ {0, 0, 1, 0, 1, 1, 0, 1, 1, 1, 0, 1, 1, 1, 1, 0, 0, 0, 1, 0, 0, 1, 0, 1, 1, 1},
+ {0, 1, 0, 0, 0, 0, 1, 1, 1, 0, 1, 1, 1, 0, 1, 0, 0, 1, 0, 0, 0, 0, 1, 1, 1, 0},
+ {0, 1, 0, 0, 0, 1, 1, 1, 1, 0, 1, 1, 0, 1, 0, 0, 0, 1, 0, 0, 0, 1, 1, 1, 1, 0},
+ {0, 0, 0, 1, 1, 0, 1, 0, 1, 1, 1, 0, 0, 1, 0, 0, 0, 0, 0, 1, 1, 0, 1, 0, 1, 1},
+ {0, 1, 0, 0, 1, 1, 1, 0, 1, 0, 1, 1, 0, 0, 0, 0, 0, 1, 0, 0, 1, 1, 1, 0, 1, 0},
+ {1, 0, 1, 0, 0, 1, 1, 1, 1, 1, 0, 1, 1, 0, 0, 0, 1, 0, 1, 0, 0, 1, 1, 1, 1, 1},
+ {1, 1, 1, 0, 1, 1, 1, 1, 0, 0, 0, 1, 0, 0, 1, 0, 1, 1, 1, 0, 1, 1, 1, 1, 0, 0},
+ {0, 1, 1, 1, 0, 0, 0, 1, 0, 1, 1, 1, 0, 0, 0, 1, 0, 1, 1, 1, 0, 0, 0, 1, 0, 1} // DUMMY
+};
+
+
+//Dummy burst 0xFB 76 0A 4E 09 10 1F 1C 5C 5C 57 4A 33 39 E9 F1 2F A8
+static const unsigned char dummy_burst[] = {
+ 0, 0, 0,
+ 1, 1, 1, 1, 1, 0, 1, 1, 0, 1,
+ 1, 1, 0, 1, 1, 0, 0, 0, 0, 0,
+ 1, 0, 1, 0, 0, 1, 0, 0, 1, 1,
+ 1, 0, 0, 0, 0, 0, 1, 0, 0, 1,
+ 0, 0, 0, 1, 0, 0, 0, 0, 0, 0,
+ 0, 1, 1, 1, 1, 1, 0, 0,
+
+ 0, 1, 1, 1, 0, 0, 0, 1, 0, 1,
+ 1, 1, 0, 0, 0, 1, 0, 1, 1, 1,
+ 0, 0, 0, 1, 0, 1,
+
+ 0, 1, 1, 1, 0, 1, 0, 0, 1, 0,
+ 1, 0, 0, 0, 1, 1, 0, 0, 1, 1,
+ 0, 0, 1, 1, 1, 0, 0, 1, 1, 1,
+ 1, 0, 1, 0, 0, 1, 1, 1, 1, 1,
+ 0, 0, 0, 1, 0, 0, 1, 0, 1, 1,
+ 1, 1, 1, 0, 1, 0, 1, 0,
+ 0, 0, 0
+};
+
+
+/*
+ * The frequency correction burst is used for frequency synchronization
+ * of the mobile. This is broadcast in TS0 together with the SCH and
+ * BCCH.
+ *
+ * Modulating the bits below causes a spike at 62.5kHz above (below for
+ * COMPACT) the center frequency. One can use this spike with a narrow
+ * band filter to accurately determine the center of the channel.
+ */
+static const unsigned char fc_fb[] = {
+ 0, 0, 0, //I don't use this tables,
+ 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, //I copied this here from burst_types.h because
+ 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, //the description is very informative - p.krysik
+ 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
+ 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
+ 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
+ 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
+ 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
+ 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
+ 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
+ 0, 0, 0
+};
+
+static const unsigned char fc_compact_fb[] = {
+ 1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0,
+ 1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0,
+ 1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0,
+ 1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0,
+ 1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0,
+ 1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0,
+ 1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0,
+ 1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0,
+ 1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0
+};
+
+
+#endif /* INCLUDED_GSM_CONSTANTS_H */
diff --git a/lib/receiver_config.cc b/lib/receiver_config.cc
new file mode 100644
index 0000000..c7c996e
--- /dev/null
+++ b/lib/receiver_config.cc
@@ -0,0 +1,84 @@
+/* -*- c++ -*- */
+/*
+ * @file
+ * @author Piotr Krysik <pkrysik@stud.elka.pw.edu.pl>
+ * @section LICENSE
+ *
+ * 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.
+ *
+ * @section DESCRIPTION
+ * This file contains classes which define gsm_receiver configuration
+ * and the burst_counter which is used to store internal state of the receiver
+ * when it's synchronized
+ */
+#ifdef HAVE_CONFIG_H
+#include "config.h"
+#endif
+
+#include <receiver_config.h>
+
+burst_counter & burst_counter::operator++(int)
+{
+ d_timeslot_nr++;
+ if (d_timeslot_nr == TS_PER_FRAME) {
+ d_timeslot_nr = 0;
+
+ if ((d_t2 == 25) && (d_t3 == 50)) {
+ d_t1 = (d_t1 + 1) % (1 << 11);
+ }
+
+ d_t2 = (d_t2 + 1) % 26;
+ d_t3 = (d_t3 + 1) % 51;
+ }
+
+ //update offset - this is integer for d_OSR which is multiple of four
+ d_offset_fractional += GUARD_FRACTIONAL * d_OSR;
+ d_offset_integer = floor(d_offset_fractional);
+ d_offset_fractional = d_offset_fractional - d_offset_integer;
+ return (*this);
+}
+
+void burst_counter::set(uint32_t t1, uint32_t t2, uint32_t t3, uint32_t timeslot_nr)
+{
+ d_t1 = t1;
+ d_t2 = t2;
+ d_t3 = t3;
+ d_timeslot_nr = timeslot_nr;
+ double first_sample_position = (get_frame_nr() * 8 + timeslot_nr) * TS_BITS;
+ d_offset_fractional = first_sample_position - floor(first_sample_position);
+ d_offset_integer = 0;
+}
+
+burst_type channel_configuration::get_burst_type(burst_counter burst_nr)
+{
+ uint32_t timeslot_nr = burst_nr.get_timeslot_nr();
+ multiframe_type m_type = d_timeslots_descriptions[timeslot_nr].get_type();
+ uint32_t nr;
+
+ switch (m_type) {
+ case multiframe_26:
+ nr = burst_nr.get_t2();
+ break;
+ case multiframe_51:
+ nr = burst_nr.get_t3();
+ break;
+ default:
+ nr = 0;
+ break;
+ }
+
+ return d_timeslots_descriptions[timeslot_nr].get_burst_type(nr);
+}
diff --git a/lib/receiver_config.h b/lib/receiver_config.h
new file mode 100644
index 0000000..b7ba43a
--- /dev/null
+++ b/lib/receiver_config.h
@@ -0,0 +1,164 @@
+/* -*- c++ -*- */
+/*
+ * @file
+ * @author Piotr Krysik <pkrysik@stud.elka.pw.edu.pl>
+ * @section LICENSE
+ *
+ * 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.
+ *
+ * @section DESCRIPTION
+ * This file contains classes which define gsm_receiver configuration
+ * and the burst_counter which is used to store internal state of the receiver
+ * when it's synchronized
+ */
+#ifndef INCLUDED_GSM_RECEIVER_CONFIG_H
+#define INCLUDED_GSM_RECEIVER_CONFIG_H
+
+#include <vector>
+#include <algorithm>
+#include <math.h>
+#include <stdint.h>
+#include <gsm_constants.h>
+
+class multiframe_configuration
+{
+ private:
+ multiframe_type d_type;
+ std::vector<burst_type> d_burst_types;
+ public:
+ multiframe_configuration() {
+ d_type = unknown;
+ fill(d_burst_types.begin(), d_burst_types.end(), empty);
+ }
+
+ ~multiframe_configuration() {}
+
+ void set_type(multiframe_type type) {
+ if (type == multiframe_26) {
+ d_burst_types.resize(26);
+ } else {
+ d_burst_types.resize(51);
+ }
+
+ d_type = type;
+ }
+
+ void set_burst_type(int nr, burst_type type) {
+ d_burst_types[nr] = type;
+ }
+
+ multiframe_type get_type() {
+ return d_type;
+ }
+
+ burst_type get_burst_type(int nr) {
+ return d_burst_types[nr];
+ }
+};
+
+class burst_counter
+{
+ private:
+ const int d_OSR;
+ uint32_t d_t1, d_t2, d_t3, d_timeslot_nr;
+ double d_offset_fractional;
+ double d_offset_integer;
+ public:
+ burst_counter(int osr):
+ d_OSR(osr),
+ d_t1(0),
+ d_t2(0),
+ d_t3(0),
+ d_timeslot_nr(0),
+ d_offset_fractional(0.0),
+ d_offset_integer(0.0) {
+ }
+
+ burst_counter(int osr, uint32_t t1, uint32_t t2, uint32_t t3, uint32_t timeslot_nr):
+ d_OSR(osr),
+ d_t1(t1),
+ d_t2(t2),
+ d_t3(t3),
+ d_timeslot_nr(timeslot_nr),
+ d_offset_fractional(0.0),
+ d_offset_integer(0.0) {
+ double first_sample_position = (get_frame_nr() * 8 + timeslot_nr) * TS_BITS;
+ d_offset_integer = floor(first_sample_position);
+ d_offset_fractional = first_sample_position - floor(first_sample_position);
+ }
+
+ burst_counter & operator++(int);
+ void set(uint32_t t1, uint32_t t2, uint32_t t3, uint32_t timeslot_nr);
+
+ uint32_t get_t1() {
+ return d_t1;
+ }
+
+ uint32_t get_t2() {
+ return d_t2;
+ }
+
+ uint32_t get_t3() {
+ return d_t3;
+ }
+
+ uint32_t get_timeslot_nr() {
+ return d_timeslot_nr;
+ }
+
+ uint32_t get_frame_nr() {
+ return (51 * 26 * d_t1) + (51 * (((d_t3 + 26) - d_t2) % 26)) + d_t3;
+ }
+
+ uint32_t get_frame_nr_mod() {
+ return (d_t1 << 11) + (d_t3 << 5) + d_t2;
+ }
+
+ unsigned get_offset() {
+ return (unsigned)d_offset_integer;
+ }
+};
+
+class channel_configuration
+{
+ private:
+ multiframe_configuration d_timeslots_descriptions[TS_PER_FRAME];
+ public:
+ channel_configuration() {
+ for (int i = 0; i < TS_PER_FRAME; i++) {
+ d_timeslots_descriptions[i].set_type(unknown);
+ }
+ }
+
+ void set_multiframe_type(int timeslot_nr, multiframe_type type) {
+ d_timeslots_descriptions[timeslot_nr].set_type(type);
+ }
+
+ void set_burst_types(int timeslot_nr, const unsigned mapping[], unsigned mapping_size, burst_type b_type) {
+ unsigned i;
+ for (i = 0; i < mapping_size; i++) {
+ d_timeslots_descriptions[timeslot_nr].set_burst_type(mapping[i], b_type);
+ }
+ }
+
+ void set_single_burst_type(int timeslot_nr, int burst_nr, burst_type b_type) {
+ d_timeslots_descriptions[timeslot_nr].set_burst_type(burst_nr, b_type);
+ }
+
+ burst_type get_burst_type(burst_counter burst_nr);
+};
+
+#endif /* INCLUDED_GSM_RECEIVER_CONFIG_H */
diff --git a/lib/receiver_impl.cc b/lib/receiver_impl.cc
new file mode 100644
index 0000000..689079b
--- /dev/null
+++ b/lib/receiver_impl.cc
@@ -0,0 +1,770 @@
+/* -*- c++ -*- */
+/*
+ * Copyright 2014 <+YOU OR YOUR COMPANY+>.
+ *
+ * This 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.
+ *
+ * This software 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 this software; see the file COPYING. If not, write to
+ * the Free Software Foundation, Inc., 51 Franklin Street,
+ * Boston, MA 02110-1301, USA.
+ */
+
+#ifdef HAVE_CONFIG_H
+#include "config.h"
+#endif
+
+#include <gnuradio/io_signature.h>
+#include "receiver_impl.h"
+
+#include <gnuradio/io_signature.h>
+#include <gnuradio/math.h>
+#include <math.h>
+#include <boost/circular_buffer.hpp>
+#include <algorithm>
+#include <numeric>
+#include <viterbi_detector.h>
+#include <string.h>
+#include <sch.h>
+#include <iostream>
+#include <iomanip>
+
+#include <assert.h>
+
+#define SYNC_SEARCH_RANGE 30
+
+namespace gr {
+ namespace gsm {
+
+ typedef std::list<float> list_float;
+ typedef std::vector<float> vector_float;
+
+ typedef boost::circular_buffer<float> circular_buffer_float;
+
+ receiver::sptr
+ receiver::make(feval_dd * tuner, int osr)
+ {
+ return gnuradio::get_initial_sptr
+ (new receiver_impl(tuner, osr));
+ }
+
+ /*
+ * The private constructor
+ */
+ receiver_impl::receiver_impl(feval_dd * tuner, int osr)
+ : gr::block("receiver",
+ gr::io_signature::make(1, 1, sizeof(gr_complex)),
+ gr::io_signature::make(0, 1, 142 * sizeof(float))),
+ d_OSR(osr),
+ d_chan_imp_length(CHAN_IMP_RESP_LENGTH),
+ d_tuner(tuner),
+ d_counter(0),
+ d_fcch_start_pos(0),
+ d_freq_offset(0),
+ d_state(first_fcch_search),
+ d_burst_nr(osr),
+ d_failed_sch(0)
+ {
+ int i;
+ gmsk_mapper(SYNC_BITS, N_SYNC_BITS, d_sch_training_seq, gr_complex(0.0, -1.0));
+ for (i = 0; i < TRAIN_SEQ_NUM; i++) {
+ gr_complex startpoint;
+ if (i == 6 || i == 7) { //this is nasty hack
+ startpoint = gr_complex(-1.0, 0.0); //if I don't change it here all bits of normal bursts for BTSes with bcc=6 will have reversed values
+ } else {
+ startpoint = gr_complex(1.0, 0.0); //I've checked this hack for bcc==0,1,2,3,4,6
+ } //I don't know what about bcc==5 and 7 yet
+ //TODO:find source of this situation - this is purely mathematical problem I guess
+
+ gmsk_mapper(train_seq[i], N_TRAIN_BITS, d_norm_training_seq[i], startpoint);
+ }
+ }
+
+ /*
+ * Our virtual destructor.
+ */
+ receiver_impl::~receiver_impl()
+ {
+ }
+
+ void receiver_impl::forecast(int noutput_items, gr_vector_int &ninput_items_required)
+ {
+ ninput_items_required[0] = noutput_items * floor((TS_BITS + 2 * GUARD_PERIOD) * d_OSR);
+ }
+
+
+ int
+ receiver_impl::general_work(int noutput_items,
+ gr_vector_int &ninput_items,
+ gr_vector_const_void_star &input_items,
+ gr_vector_void_star &output_items)
+ {
+ const gr_complex *input = (const gr_complex *) input_items[0];
+ //float *out = (float *) output_items[0];
+ int produced_out = 0; //how many output elements were produced - this isn't used yet
+ //probably the gsm receiver will be changed into sink so this variable won't be necessary
+ switch (d_state) {
+ //bootstrapping
+ case first_fcch_search:
+ if (find_fcch_burst(input, ninput_items[0])) { //find frequency correction burst in the input buffer
+ set_frequency(d_freq_offset); //if fcch search is successful set frequency offset
+ //produced_out = 0;
+ d_state = next_fcch_search;
+ } else {
+ //produced_out = 0;
+ d_state = first_fcch_search;
+ }
+ break;
+
+ case next_fcch_search: { //this state is used because it takes some time (a bunch of buffered samples)
+ COUT("fcch");
+ float prev_freq_offset = d_freq_offset; //before previous set_frequqency cause change
+ if (find_fcch_burst(input, ninput_items[0])) {
+ if (abs(prev_freq_offset - d_freq_offset) > FCCH_MAX_FREQ_OFFSET) {
+ set_frequency(d_freq_offset); //call set_frequncy only frequency offset change is greater than some value
+ }
+ //produced_out = 0;
+ d_state = sch_search;
+ } else {
+ //produced_out = 0;
+ d_state = next_fcch_search;
+ }
+ break;
+ }
+
+
+ case sch_search: {
+ vector_complex channel_imp_resp(CHAN_IMP_RESP_LENGTH*d_OSR);
+ int t1, t2, t3;
+ int burst_start = 0;
+ unsigned char output_binary[BURST_SIZE];
+
+ if (reach_sch_burst(ninput_items[0])) { //wait for a SCH burst
+ burst_start = get_sch_chan_imp_resp(input, &channel_imp_resp[0]); //get channel impulse response from it
+ detect_burst(input, &channel_imp_resp[0], burst_start, output_binary); //detect bits using MLSE detection
+ if (decode_sch(&output_binary[3], &t1, &t2, &t3, &d_ncc, &d_bcc) == 0) { //decode SCH burst
+ COUT("sch burst_start: " << burst_start);
+ COUT("bcc: " << d_bcc << " ncc: " << d_ncc << " t1: " << t1 << " t2: " << t2 << " t3: " << t3);
+ d_burst_nr.set(t1, t2, t3, 0); //set counter of bursts value
+
+ //configure the receiver - tell him where to find which burst type
+ d_channel_conf.set_multiframe_type(TIMESLOT0, multiframe_51); //in the timeslot nr.0 bursts changes according to t3 counter
+ configure_receiver();//TODO: this shouldn't be here - remove it when gsm receiver's interface will be ready
+ d_channel_conf.set_burst_types(TIMESLOT0, FCCH_FRAMES, sizeof(FCCH_FRAMES) / sizeof(unsigned), fcch_burst); //tell where to find fcch bursts
+ d_channel_conf.set_burst_types(TIMESLOT0, SCH_FRAMES, sizeof(SCH_FRAMES) / sizeof(unsigned), sch_burst); //sch bursts
+ d_channel_conf.set_burst_types(TIMESLOT0, BCCH_FRAMES, sizeof(BCCH_FRAMES) / sizeof(unsigned), normal_burst);//!and maybe normal bursts of the BCCH logical channel
+ d_burst_nr++;
+
+ consume_each(burst_start + BURST_SIZE * d_OSR); //consume samples up to next guard period
+ d_state = synchronized;
+ } else {
+ d_state = next_fcch_search; //if there is error in the sch burst go back to fcch search phase
+ }
+ } else {
+ d_state = sch_search;
+ }
+ break;
+ }
+ //in this state receiver is synchronized and it processes bursts according to burst type for given burst number
+ case synchronized: {
+ vector_complex channel_imp_resp(CHAN_IMP_RESP_LENGTH*d_OSR);
+ int burst_start;
+ int offset = 0;
+ int to_consume = 0;
+ unsigned char output_binary[BURST_SIZE];
+
+ burst_type b_type = d_channel_conf.get_burst_type(d_burst_nr); //get burst type for given burst number
+
+ switch (b_type) {
+ case fcch_burst: { //if it's FCCH burst
+ const unsigned first_sample = ceil((GUARD_PERIOD + 2 * TAIL_BITS) * d_OSR) + 1;
+ const unsigned last_sample = first_sample + USEFUL_BITS * d_OSR - TAIL_BITS * d_OSR;
+ double freq_offset = compute_freq_offset(input, first_sample, last_sample); //extract frequency offset from it
+
+ d_freq_offset_vals.push_front(freq_offset);
+ //process_normal_burst(d_burst_nr, fc_fb);
+ if (d_freq_offset_vals.size() >= 10) {
+ double sum = std::accumulate(d_freq_offset_vals.begin(), d_freq_offset_vals.end(), 0);
+ double mean_offset = sum / d_freq_offset_vals.size(); //compute mean
+ d_freq_offset_vals.clear();
+ if (abs(mean_offset) > FCCH_MAX_FREQ_OFFSET) {
+ d_freq_offset -= mean_offset; //and adjust frequency if it have changed beyond
+ set_frequency(d_freq_offset); //some limit
+ DCOUT("mean_offset: " << mean_offset);
+ DCOUT("Adjusting frequency, new frequency offset: " << d_freq_offset << "\n");
+ }
+ }
+ }
+ break;
+ case sch_burst: { //if it's SCH burst
+ int t1, t2, t3, d_ncc, d_bcc;
+ burst_start = get_sch_chan_imp_resp(input, &channel_imp_resp[0]); //get channel impulse response
+ detect_burst(input, &channel_imp_resp[0], burst_start, output_binary); //MLSE detection of bits
+ //process_normal_burst(d_burst_nr, output_binary);
+ if (decode_sch(&output_binary[3], &t1, &t2, &t3, &d_ncc, &d_bcc) == 0) { //and decode SCH data
+ // d_burst_nr.set(t1, t2, t3, 0); //but only to check if burst_start value is correct
+ d_failed_sch = 0;
+ DCOUT("bcc: " << d_bcc << " ncc: " << d_ncc << " t1: " << t1 << " t2: " << t2 << " t3: " << t3);
+ offset = burst_start - floor((GUARD_PERIOD) * d_OSR); //compute offset from burst_start - burst should start after a guard period
+ DCOUT(offset);
+ to_consume += offset; //adjust with offset number of samples to be consumed
+ } else {
+ d_failed_sch++;
+ if (d_failed_sch >= MAX_SCH_ERRORS) {
+ // d_state = next_fcch_search; //TODO: this isn't good, the receiver is going wild when it goes back to next_fcch_search from here
+ // d_freq_offset_vals.clear();
+ DCOUT("many sch decoding errors");
+ }
+ }
+ }
+ break;
+
+ case normal_burst: //if it's normal burst
+ burst_start = get_norm_chan_imp_resp(input, &channel_imp_resp[0], d_bcc); //get channel impulse response for given training sequence number - d_bcc
+ detect_burst(input, &channel_imp_resp[0], burst_start, output_binary); //MLSE detection of bits
+ process_normal_burst(d_burst_nr, output_binary); //TODO: this shouldn't be here - remove it when gsm receiver's interface will be ready
+ break;
+
+ case dummy_or_normal: {
+ burst_start = get_norm_chan_imp_resp(input, &channel_imp_resp[0], TS_DUMMY);
+ detect_burst(input, &channel_imp_resp[0], burst_start, output_binary);
+
+ std::vector<unsigned char> v(20);
+ std::vector<unsigned char>::iterator it;
+ it = std::set_difference(output_binary + TRAIN_POS, output_binary + TRAIN_POS + 16, &train_seq[TS_DUMMY][5], &train_seq[TS_DUMMY][21], v.begin());
+ int different_bits = (it - v.begin());
+
+ if (different_bits > 2) {
+ burst_start = get_norm_chan_imp_resp(input, &channel_imp_resp[0], d_bcc);
+ detect_burst(input, &channel_imp_resp[0], burst_start, output_binary);
+ //if (!output_binary[0] && !output_binary[1] && !output_binary[2]) {
+ COUT("Normal burst");
+ process_normal_burst(d_burst_nr, output_binary); //TODO: this shouldn't be here - remove it when gsm receiver's interface will be ready
+ //}
+ } else {
+ //process_normal_burst(d_burst_nr, dummy_burst);
+ }
+ }
+ case rach_burst:
+ //implementation of this channel isn't possible in current gsm_receiver
+ //it would take some realtime processing, counter of samples from USRP to
+ //stay synchronized with this device and possibility to switch frequency from uplink
+ //to C0 (where sch is) back and forth
+
+ break;
+ case dummy: //if it's dummy
+ burst_start = get_norm_chan_imp_resp(input, &channel_imp_resp[0], TS_DUMMY); //read dummy
+ detect_burst(input, &channel_imp_resp[0], burst_start, output_binary); // but as far as I know it's pointless
+ break;
+ case empty: //if it's empty burst
+ break; //do nothing
+ }
+
+ d_burst_nr++; //go to next burst
+
+ to_consume += TS_BITS * d_OSR + d_burst_nr.get_offset(); //consume samples of the burst up to next guard period
+ //and add offset which is introduced by
+ //0.25 fractional part of a guard period
+ //burst_number computes this offset
+ //but choice of this class to do this was random
+ consume_each(to_consume);
+ }
+ break;
+ }
+
+ return produced_out;
+ }
+
+
+ bool receiver_impl::find_fcch_burst(const gr_complex *input, const int nitems)
+ {
+ circular_buffer_float phase_diff_buffer(FCCH_HITS_NEEDED * d_OSR); //circular buffer used to scan throug signal to find
+ //best match for FCCH burst
+ float phase_diff = 0;
+ gr_complex conjprod;
+ int start_pos = -1;
+ int hit_count = 0;
+ int miss_count = 0;
+ float min_phase_diff;
+ float max_phase_diff;
+ double best_sum = 0;
+ float lowest_max_min_diff = 99999;
+
+ int to_consume = 0;
+ int sample_number = 0;
+ bool end = false;
+ bool result = false;
+ circular_buffer_float::iterator buffer_iter;
+
+ /**@name Possible states of FCCH search algorithm*/
+ //@{
+ enum states {
+ init, ///< initialize variables
+ search, ///< search for positive samples
+ found_something, ///< search for FCCH and the best position of it
+ fcch_found, ///< when FCCH was found
+ search_fail ///< when there is no FCCH in the input vector
+ } fcch_search_state;
+ //@}
+
+ fcch_search_state = init;
+
+ while (!end) {
+ switch (fcch_search_state) {
+
+ case init: //initialize variables
+ hit_count = 0;
+ miss_count = 0;
+ start_pos = -1;
+ lowest_max_min_diff = 99999;
+ phase_diff_buffer.clear();
+ fcch_search_state = search;
+
+ break;
+
+ case search: // search for positive samples
+ sample_number++;
+
+ if (sample_number > nitems - FCCH_HITS_NEEDED * d_OSR) { //if it isn't possible to find FCCH because
+ //there's too few samples left to look into,
+ to_consume = sample_number; //don't do anything with those samples which are left
+ //and consume only those which were checked
+ fcch_search_state = search_fail;
+ } else {
+ phase_diff = compute_phase_diff(input[sample_number], input[sample_number-1]);
+
+ if (phase_diff > 0) { //if a positive phase difference was found
+ to_consume = sample_number;
+ fcch_search_state = found_something; //switch to state in which searches for FCCH
+ } else {
+ fcch_search_state = search;
+ }
+ }
+
+ break;
+
+ case found_something: {// search for FCCH and the best position of it
+ if (phase_diff > 0) {
+ hit_count++; //positive phase differencies increases hits_count
+ } else {
+ miss_count++; //negative increases miss_count
+ }
+
+ if ((miss_count >= FCCH_MAX_MISSES * d_OSR) && (hit_count <= FCCH_HITS_NEEDED * d_OSR)) {
+ //if miss_count exceeds limit before hit_count
+ fcch_search_state = init; //go to init
+ continue;
+ } else if (((miss_count >= FCCH_MAX_MISSES * d_OSR) && (hit_count > FCCH_HITS_NEEDED * d_OSR)) || (hit_count > 2 * FCCH_HITS_NEEDED * d_OSR)) {
+ //if hit_count and miss_count exceeds limit then FCCH was found
+ fcch_search_state = fcch_found;
+ continue;
+ } else if ((miss_count < FCCH_MAX_MISSES * d_OSR) && (hit_count > FCCH_HITS_NEEDED * d_OSR)) {
+ //find difference between minimal and maximal element in the buffer
+ //for FCCH this value should be low
+ //this part is searching for a region where this value is lowest
+ min_phase_diff = * (min_element(phase_diff_buffer.begin(), phase_diff_buffer.end()));
+ max_phase_diff = * (max_element(phase_diff_buffer.begin(), phase_diff_buffer.end()));
+
+ if (lowest_max_min_diff > max_phase_diff - min_phase_diff) {
+ lowest_max_min_diff = max_phase_diff - min_phase_diff;
+ start_pos = sample_number - FCCH_HITS_NEEDED * d_OSR - FCCH_MAX_MISSES * d_OSR; //store start pos
+ best_sum = 0;
+
+ for (buffer_iter = phase_diff_buffer.begin();
+ buffer_iter != (phase_diff_buffer.end());
+ buffer_iter++) {
+ best_sum += *buffer_iter - (M_PI / 2) / d_OSR; //store best value of phase offset sum
+ }
+ }
+ }
+
+ sample_number++;
+
+ if (sample_number >= nitems) { //if there's no single sample left to check
+ fcch_search_state = search_fail;//FCCH search failed
+ continue;
+ }
+
+ phase_diff = compute_phase_diff(input[sample_number], input[sample_number-1]);
+ phase_diff_buffer.push_back(phase_diff);
+ fcch_search_state = found_something;
+ }
+ break;
+
+ case fcch_found: {
+ DCOUT("fcch found on position: " << d_counter + start_pos);
+ to_consume = start_pos + FCCH_HITS_NEEDED * d_OSR + 1; //consume one FCCH burst
+
+ d_fcch_start_pos = d_counter + start_pos;
+
+ //compute frequency offset
+ double phase_offset = best_sum / FCCH_HITS_NEEDED;
+ double freq_offset = phase_offset * 1625000.0 / (12.0 * M_PI);
+ d_freq_offset -= freq_offset;
+ DCOUT("freq_offset: " << d_freq_offset);
+
+ end = true;
+ result = true;
+ break;
+ }
+
+ case search_fail:
+ end = true;
+ result = false;
+ break;
+ }
+ }
+
+ d_counter += to_consume;
+ consume_each(to_consume);
+
+ return result;
+ }
+
+
+ double receiver_impl::compute_freq_offset(const gr_complex * input, unsigned first_sample, unsigned last_sample)
+ {
+ double phase_sum = 0;
+ unsigned ii;
+
+ for (ii = first_sample; ii < last_sample; ii++) {
+ double phase_diff = compute_phase_diff(input[ii], input[ii-1]) - (M_PI / 2) / d_OSR;
+ phase_sum += phase_diff;
+ }
+
+ double phase_offset = phase_sum / (last_sample - first_sample);
+ double freq_offset = phase_offset * 1625000.0 / (12.0 * M_PI);
+ return freq_offset;
+ }
+
+ void receiver_impl::set_frequency(double freq_offset)
+ {
+ d_tuner->calleval(freq_offset);
+ }
+
+ inline float receiver_impl::compute_phase_diff(gr_complex val1, gr_complex val2)
+ {
+ gr_complex conjprod = val1 * conj(val2);
+ return fast_atan2f(imag(conjprod), real(conjprod));
+ }
+
+ bool receiver_impl::reach_sch_burst(const int nitems)
+ {
+ //it just consumes samples to get near to a SCH burst
+ int to_consume = 0;
+ bool result = false;
+ unsigned sample_nr_near_sch_start = d_fcch_start_pos + (FRAME_BITS - SAFETY_MARGIN) * d_OSR;
+
+ //consume samples until d_counter will be equal to sample_nr_near_sch_start
+ if (d_counter < sample_nr_near_sch_start) {
+ if (d_counter + nitems >= sample_nr_near_sch_start) {
+ to_consume = sample_nr_near_sch_start - d_counter;
+ } else {
+ to_consume = nitems;
+ }
+ result = false;
+ } else {
+ to_consume = 0;
+ result = true;
+ }
+
+ d_counter += to_consume;
+ consume_each(to_consume);
+ return result;
+ }
+
+ int receiver_impl::get_sch_chan_imp_resp(const gr_complex *input, gr_complex * chan_imp_resp)
+ {
+ vector_complex correlation_buffer;
+ vector_float power_buffer;
+ vector_float window_energy_buffer;
+
+ int strongest_window_nr;
+ int burst_start = 0;
+ int chan_imp_resp_center = 0;
+ float max_correlation = 0;
+ float energy = 0;
+
+ for (int ii = SYNC_POS * d_OSR; ii < (SYNC_POS + SYNC_SEARCH_RANGE) *d_OSR; ii++) {
+ gr_complex correlation = correlate_sequence(&d_sch_training_seq[5], N_SYNC_BITS - 10, &input[ii]);
+ correlation_buffer.push_back(correlation);
+ power_buffer.push_back(std::pow(abs(correlation), 2));
+ }
+
+ //compute window energies
+ vector_float::iterator iter = power_buffer.begin();
+ bool loop_end = false;
+ while (iter != power_buffer.end()) {
+ vector_float::iterator iter_ii = iter;
+ energy = 0;
+
+ for (int ii = 0; ii < (d_chan_imp_length) *d_OSR; ii++, iter_ii++) {
+ if (iter_ii == power_buffer.end()) {
+ loop_end = true;
+ break;
+ }
+ energy += (*iter_ii);
+ }
+ if (loop_end) {
+ break;
+ }
+ iter++;
+ window_energy_buffer.push_back(energy);
+ }
+
+ strongest_window_nr = max_element(window_energy_buffer.begin(), window_energy_buffer.end()) - window_energy_buffer.begin();
+ // d_channel_imp_resp.clear();
+
+ max_correlation = 0;
+ for (int ii = 0; ii < (d_chan_imp_length) *d_OSR; ii++) {
+ gr_complex correlation = correlation_buffer[strongest_window_nr + ii];
+ if (abs(correlation) > max_correlation) {
+ chan_imp_resp_center = ii;
+ max_correlation = abs(correlation);
+ }
+ // d_channel_imp_resp.push_back(correlation);
+ chan_imp_resp[ii] = correlation;
+ }
+
+ burst_start = strongest_window_nr + chan_imp_resp_center - 48 * d_OSR - 2 * d_OSR + 2 + SYNC_POS * d_OSR;
+ return burst_start;
+ }
+
+
+
+ void receiver_impl::detect_burst(const gr_complex * input, gr_complex * chan_imp_resp, int burst_start, unsigned char * output_binary)
+ {
+ float output[BURST_SIZE];
+ gr_complex rhh_temp[CHAN_IMP_RESP_LENGTH*d_OSR];
+ gr_complex rhh[CHAN_IMP_RESP_LENGTH];
+ gr_complex filtered_burst[BURST_SIZE];
+ int start_state = 3;
+ unsigned int stop_states[2] = {4, 12};
+
+ autocorrelation(chan_imp_resp, rhh_temp, d_chan_imp_length*d_OSR);
+ for (int ii = 0; ii < (d_chan_imp_length); ii++) {
+ rhh[ii] = conj(rhh_temp[ii*d_OSR]);
+ }
+
+ mafi(&input[burst_start], BURST_SIZE, chan_imp_resp, d_chan_imp_length*d_OSR, filtered_burst);
+
+ viterbi_detector(filtered_burst, BURST_SIZE, rhh, start_state, stop_states, 2, output);
+
+ for (int i = 0; i < BURST_SIZE ; i++) {
+ output_binary[i] = (output[i] > 0);
+ }
+ }
+
+ //TODO consider placing this funtion in a separate class for signal processing
+ void receiver_impl::gmsk_mapper(const unsigned char * input, int nitems, gr_complex * gmsk_output, gr_complex start_point)
+ {
+ gr_complex j = gr_complex(0.0, 1.0);
+
+ int current_symbol;
+ int encoded_symbol;
+ int previous_symbol = 2 * input[0] - 1;
+ gmsk_output[0] = start_point;
+
+ for (int i = 1; i < nitems; i++) {
+ //change bits representation to NRZ
+ current_symbol = 2 * input[i] - 1;
+ //differentially encode
+ encoded_symbol = current_symbol * previous_symbol;
+ //and do gmsk mapping
+ gmsk_output[i] = j * gr_complex(encoded_symbol, 0.0) * gmsk_output[i-1];
+ previous_symbol = current_symbol;
+ }
+ }
+
+ //TODO consider use of some generalized function for correlation and placing it in a separate class for signal processing
+ gr_complex receiver_impl::correlate_sequence(const gr_complex * sequence, int length, const gr_complex * input)
+ {
+ gr_complex result(0.0, 0.0);
+ int sample_number = 0;
+
+ for (int ii = 0; ii < length; ii++) {
+ sample_number = (ii * d_OSR) ;
+ result += sequence[ii] * conj(input[sample_number]);
+ }
+
+ result = result / gr_complex(length, 0);
+ return result;
+ }
+
+ //computes autocorrelation for positive arguments
+ //TODO consider placing this funtion in a separate class for signal processing
+ inline void receiver_impl::autocorrelation(const gr_complex * input, gr_complex * out, int nitems)
+ {
+ int i, k;
+ for (k = nitems - 1; k >= 0; k--) {
+ out[k] = gr_complex(0, 0);
+ for (i = k; i < nitems; i++) {
+ out[k] += input[i] * conj(input[i-k]);
+ }
+ }
+ }
+
+ //TODO consider use of some generalized function for filtering and placing it in a separate class for signal processing
+ inline void receiver_impl::mafi(const gr_complex * input, int nitems, gr_complex * filter, int filter_length, gr_complex * output)
+ {
+ int ii = 0, n, a;
+
+ for (n = 0; n < nitems; n++) {
+ a = n * d_OSR;
+ output[n] = 0;
+ ii = 0;
+
+ while (ii < filter_length) {
+ if ((a + ii) >= nitems*d_OSR)
+ break;
+ output[n] += input[a+ii] * filter[ii];
+ ii++;
+ }
+ }
+ }
+
+ //TODO: get_norm_chan_imp_resp is similar to get_sch_chan_imp_resp - consider joining this two functions
+ //TODO: this is place where most errors are introduced and can be corrected by improvements to this fuction
+ //especially computations of strongest_window_nr
+ int receiver_impl::get_norm_chan_imp_resp(const gr_complex *input, gr_complex * chan_imp_resp, int bcc)
+ {
+ vector_complex correlation_buffer;
+ vector_float power_buffer;
+ vector_float window_energy_buffer;
+
+ int strongest_window_nr;
+ int burst_start = 0;
+ int chan_imp_resp_center = 0;
+ float max_correlation = 0;
+ float energy = 0;
+
+ int search_center = (int)((TRAIN_POS + GUARD_PERIOD) * d_OSR);
+ int search_start_pos = search_center + 1;
+ // int search_start_pos = search_center - d_chan_imp_length * d_OSR;
+ int search_stop_pos = search_center + d_chan_imp_length * d_OSR + 2 * d_OSR;
+
+ for (int ii = search_start_pos; ii < search_stop_pos; ii++) {
+ gr_complex correlation = correlate_sequence(&d_norm_training_seq[bcc][TRAIN_BEGINNING], N_TRAIN_BITS - 10, &input[ii]);
+
+ correlation_buffer.push_back(correlation);
+ power_buffer.push_back(std::pow(abs(correlation), 2));
+ }
+
+ //compute window energies
+ vector_float::iterator iter = power_buffer.begin();
+ bool loop_end = false;
+ while (iter != power_buffer.end()) {
+ vector_float::iterator iter_ii = iter;
+ energy = 0;
+
+ for (int ii = 0; ii < (d_chan_imp_length - 2)*d_OSR; ii++, iter_ii++) {
+ // for (int ii = 0; ii < (d_chan_imp_length)*d_OSR; ii++, iter_ii++) {
+ if (iter_ii == power_buffer.end()) {
+ loop_end = true;
+ break;
+ }
+ energy += (*iter_ii);
+ }
+ if (loop_end) {
+ break;
+ }
+ iter++;
+
+ window_energy_buffer.push_back(energy);
+ }
+ //!why doesn't this work
+ int strongest_window_nr_new = max_element(window_energy_buffer.begin(), window_energy_buffer.end()) - window_energy_buffer.begin();
+ strongest_window_nr = 3; //! so I have to override it here
+
+ max_correlation = 0;
+ for (int ii = 0; ii < (d_chan_imp_length)*d_OSR; ii++) {
+ gr_complex correlation = correlation_buffer[strongest_window_nr + ii];
+ if (abs(correlation) > max_correlation) {
+ chan_imp_resp_center = ii;
+ max_correlation = abs(correlation);
+ }
+ // d_channel_imp_resp.push_back(correlation);
+ chan_imp_resp[ii] = correlation;
+ }
+ // We want to use the first sample of the impulseresponse, and the
+ // corresponding samples of the received signal.
+ // the variable sync_w should contain the beginning of the used part of
+ // training sequence, which is 3+57+1+6=67 bits into the burst. That is
+ // we have that sync_t16 equals first sample in bit number 67.
+
+ burst_start = search_start_pos + chan_imp_resp_center + strongest_window_nr - TRAIN_POS * d_OSR;
+
+ // GMSK modulator introduces ISI - each bit is expanded for 3*Tb
+ // and it's maximum value is in the last bit period, so burst starts
+ // 2*Tb earlier
+ burst_start -= 2 * d_OSR;
+ burst_start += 2;
+ COUT("Poczatek ###############################");
+ std::cout << " burst_start: " << burst_start << " center: " << ((float)(search_start_pos + strongest_window_nr + chan_imp_resp_center)) / d_OSR << " stronegest window nr: " << strongest_window_nr << "\n";
+ COUT("burst_start_new: " << (search_start_pos + strongest_window_nr_new - TRAIN_POS * d_OSR));
+ burst_start=(search_start_pos + strongest_window_nr_new - TRAIN_POS * d_OSR)
+ return burst_start;
+ }
+
+
+ void receiver_impl::process_normal_burst(burst_counter burst_nr, const unsigned char * burst_binary)
+ {
+ int ii;
+ //std::cout << "fn:" <<burst_nr.get_frame_nr() << " ts" << burst_nr.get_timeslot_nr() << " ";
+ for(ii=0;ii<148;ii++){
+ std::cout << std::setprecision(1) << static_cast<int>(burst_binary[ii]);
+ }
+ std::cout << std::endl;
+ }
+ //TODO: this shouldn't be here also - the same reason
+ void receiver_impl::configure_receiver()
+ {
+ d_channel_conf.set_multiframe_type(TSC0, multiframe_51);
+ d_channel_conf.set_burst_types(TIMESLOT0, TEST51, sizeof(TEST51) / sizeof(unsigned), dummy_or_normal);
+
+ d_channel_conf.set_burst_types(TSC0, TEST_CCH_FRAMES, sizeof(TEST_CCH_FRAMES) / sizeof(unsigned), dummy_or_normal);
+ d_channel_conf.set_burst_types(TSC0, FCCH_FRAMES, sizeof(FCCH_FRAMES) / sizeof(unsigned), fcch_burst);
+
+ // d_channel_conf.set_multiframe_type(TIMESLOT1, multiframe_26);
+ // d_channel_conf.set_burst_types(TIMESLOT1, TRAFFIC_CHANNEL_F, sizeof(TRAFFIC_CHANNEL_F) / sizeof(unsigned), dummy_or_normal);
+ // d_channel_conf.set_multiframe_type(TIMESLOT2, multiframe_26);
+ // d_channel_conf.set_burst_types(TIMESLOT2, TRAFFIC_CHANNEL_F, sizeof(TRAFFIC_CHANNEL_F) / sizeof(unsigned), dummy_or_normal);
+ // d_channel_conf.set_multiframe_type(TIMESLOT3, multiframe_26);
+ // d_channel_conf.set_burst_types(TIMESLOT3, TRAFFIC_CHANNEL_F, sizeof(TRAFFIC_CHANNEL_F) / sizeof(unsigned), dummy_or_normal);
+ // d_channel_conf.set_multiframe_type(TIMESLOT4, multiframe_26);
+ // d_channel_conf.set_burst_types(TIMESLOT4, TRAFFIC_CHANNEL_F, sizeof(TRAFFIC_CHANNEL_F) / sizeof(unsigned), dummy_or_normal);
+ // d_channel_conf.set_multiframe_type(TIMESLOT5, multiframe_26);
+ // d_channel_conf.set_burst_types(TIMESLOT5, TRAFFIC_CHANNEL_F, sizeof(TRAFFIC_CHANNEL_F) / sizeof(unsigned), dummy_or_normal);
+ // d_channel_conf.set_multiframe_type(TIMESLOT6, multiframe_26);
+ // d_channel_conf.set_burst_types(TIMESLOT6, TRAFFIC_CHANNEL_F, sizeof(TRAFFIC_CHANNEL_F) / sizeof(unsigned), dummy_or_normal);
+ // d_channel_conf.set_multiframe_type(TIMESLOT7, multiframe_26);
+ // d_channel_conf.set_burst_types(TIMESLOT7, TRAFFIC_CHANNEL_F, sizeof(TRAFFIC_CHANNEL_F) / sizeof(unsigned), dummy_or_normal);
+ d_channel_conf.set_multiframe_type(TIMESLOT1, multiframe_51);
+ d_channel_conf.set_burst_types(TIMESLOT1, TEST51, sizeof(TEST51) / sizeof(unsigned), dummy_or_normal);
+ d_channel_conf.set_multiframe_type(TIMESLOT2, multiframe_51);
+ d_channel_conf.set_burst_types(TIMESLOT2, TEST51, sizeof(TEST51) / sizeof(unsigned), dummy_or_normal);
+ d_channel_conf.set_multiframe_type(TIMESLOT3, multiframe_51);
+ d_channel_conf.set_burst_types(TIMESLOT3, TEST51, sizeof(TEST51) / sizeof(unsigned), dummy_or_normal);
+ d_channel_conf.set_multiframe_type(TIMESLOT4, multiframe_51);
+ d_channel_conf.set_burst_types(TIMESLOT4, TEST51, sizeof(TEST51) / sizeof(unsigned), dummy_or_normal);
+ d_channel_conf.set_multiframe_type(TIMESLOT5, multiframe_51);
+ d_channel_conf.set_burst_types(TIMESLOT5, TEST51, sizeof(TEST51) / sizeof(unsigned), dummy_or_normal);
+ d_channel_conf.set_multiframe_type(TIMESLOT6, multiframe_51);
+ d_channel_conf.set_burst_types(TIMESLOT6, TEST51, sizeof(TEST51) / sizeof(unsigned), dummy_or_normal);
+ d_channel_conf.set_multiframe_type(TIMESLOT7, multiframe_51);
+ d_channel_conf.set_burst_types(TIMESLOT7, TEST51, sizeof(TEST51) / sizeof(unsigned), dummy_or_normal);
+
+ }
+
+
+ } /* namespace gsm */
+} /* namespace gr */
+
diff --git a/lib/receiver_impl.h b/lib/receiver_impl.h
new file mode 100644
index 0000000..17dec98
--- /dev/null
+++ b/lib/receiver_impl.h
@@ -0,0 +1,215 @@
+/* -*- c++ -*- */
+/*
+ * Copyright 2014 <+YOU OR YOUR COMPANY+>.
+ *
+ * This 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.
+ *
+ * This software 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 this software; see the file COPYING. If not, write to
+ * the Free Software Foundation, Inc., 51 Franklin Street,
+ * Boston, MA 02110-1301, USA.
+ */
+
+#ifndef INCLUDED_GSM_RECEIVER_IMPL_H
+#define INCLUDED_GSM_RECEIVER_IMPL_H
+
+#include <gsm/receiver.h>
+#include <gsm_constants.h>
+#include <receiver_config.h>
+
+namespace gr {
+ namespace gsm {
+
+ typedef std::vector<gr_complex> vector_complex;
+
+ class receiver_impl : public receiver
+ {
+ private:
+ /**@name Configuration of the receiver */
+ //@{
+ const int d_OSR; ///< oversampling ratio
+ const int d_chan_imp_length; ///< channel impulse length
+ //@}
+
+ gr_complex d_sch_training_seq[N_SYNC_BITS]; ///<encoded training sequence of a SCH burst
+ gr_complex d_norm_training_seq[TRAIN_SEQ_NUM][N_TRAIN_BITS]; ///<encoded training sequences of a normal bursts and dummy bursts
+
+ feval_dd *d_tuner; ///<callback to a python object which is used for frequency tunning
+
+ /** Counts samples consumed by the receiver
+ *
+ * It is used in beetween find_fcch_burst and reach_sch_burst calls.
+ * My intention was to synchronize this counter with some internal sample
+ * counter of the USRP. Simple access to such USRP's counter isn't possible
+ * so this variable isn't used in the "synchronized" state of the receiver yet.
+ */
+ unsigned d_counter;
+
+ /**@name Variables used to store result of the find_fcch_burst fuction */
+ //@{
+ unsigned d_fcch_start_pos; ///< position of the first sample of the fcch burst
+ float d_freq_offset; ///< frequency offset of the received signal
+ //@}
+ std::list<double> d_freq_offset_vals;
+
+ /**@name Identifiers of the BTS extracted from the SCH burst */
+ //@{
+ int d_ncc; ///< network color code
+ int d_bcc; ///< base station color code
+ //@}
+
+ /**@name Internal state of the gsm receiver */
+ //@{
+ enum states {
+ first_fcch_search, next_fcch_search, sch_search, // synchronization search part
+ synchronized // receiver is synchronized in this state
+ } d_state;
+ //@}
+
+ /**@name Variables which make internal state in the "synchronized" state */
+ //@{
+ burst_counter d_burst_nr; ///< frame number and timeslot number
+ channel_configuration d_channel_conf; ///< mapping of burst_counter to burst_type
+ //@}
+
+ unsigned d_failed_sch; ///< number of subsequent erroneous SCH bursts
+
+ /** Function whis is used to search a FCCH burst and to compute frequency offset before
+ * "synchronized" state of the receiver
+ *
+ * TODO: Describe the FCCH search algorithm in the documentation
+ * @param input vector with input signal
+ * @param nitems number of samples in the input vector
+ * @return
+ */
+ bool find_fcch_burst(const gr_complex *input, const int nitems);
+
+ /** Computes frequency offset from FCCH burst samples
+ *
+ * @param input vector with input samples
+ * @param first_sample number of the first sample of the FCCH busrt
+ * @param last_sample number of the last sample of the FCCH busrt
+ * @return frequency offset
+ */
+ double compute_freq_offset(const gr_complex * input, unsigned first_sample, unsigned last_sample);
+
+ /** Calls d_tuner's method to set frequency offset from Python level
+ *
+ * @param freq_offset absolute frequency offset of the received signal
+ */
+ void set_frequency(double freq_offset);
+
+ /** Computes angle between two complex numbers
+ *
+ * @param val1 first complex number
+ * @param val2 second complex number
+ * @return
+ */
+ inline float compute_phase_diff(gr_complex val1, gr_complex val2);
+
+ /** Function whis is used to get near to SCH burst
+ *
+ * @param nitems number of samples in the gsm_receiver's buffer
+ * @return true if SCH burst is near, false otherwise
+ */
+ bool reach_sch_burst(const int nitems);
+
+ /** Extracts channel impulse response from a SCH burst and computes first sample number of this burst
+ *
+ * @param input vector with input samples
+ * @param chan_imp_resp complex vector where channel impulse response will be stored
+ * @return number of first sample of the burst
+ */
+ int get_sch_chan_imp_resp(const gr_complex *input, gr_complex * chan_imp_resp);
+
+ /** MLSE detection of a burst bits
+ *
+ * Detects bits of burst using viterbi algorithm.
+ * @param input vector with input samples
+ * @param chan_imp_resp vector with the channel impulse response
+ * @param burst_start number of the first sample of the burst
+ * @param output_binary vector with output bits
+ */
+ void detect_burst(const gr_complex * input, gr_complex * chan_imp_resp, int burst_start, unsigned char * output_binary);
+
+ /** Encodes differentially input bits and maps them into MSK states
+ *
+ * @param input vector with input bits
+ * @param nitems number of samples in the "input" vector
+ * @param gmsk_output bits mapped into MSK states
+ * @param start_point first state
+ */
+ void gmsk_mapper(const unsigned char * input, int nitems, gr_complex * gmsk_output, gr_complex start_point);
+
+ /** Correlates MSK mapped sequence with input signal
+ *
+ * @param sequence MKS mapped sequence
+ * @param length length of the sequence
+ * @param input_signal vector with input samples
+ * @return correlation value
+ */
+ gr_complex correlate_sequence(const gr_complex * sequence, int length, const gr_complex * input);
+
+ /** Computes autocorrelation of input vector for positive arguments
+ *
+ * @param input vector with input samples
+ * @param out output vector
+ * @param nitems length of the input vector
+ */
+ inline void autocorrelation(const gr_complex * input, gr_complex * out, int nitems);
+
+ /** Filters input signal through channel impulse response
+ *
+ * @param input vector with input samples
+ * @param nitems number of samples to pass through filter
+ * @param filter filter taps - channel impulse response
+ * @param filter_length nember of filter taps
+ * @param output vector with filtered samples
+ */
+ inline void mafi(const gr_complex * input, int nitems, gr_complex * filter, int filter_length, gr_complex * output);
+
+ /** Extracts channel impulse response from a normal burst and computes first sample number of this burst
+ *
+ * @param input vector with input samples
+ * @param chan_imp_resp complex vector where channel impulse response will be stored
+ * @param search_range possible absolute offset of a channel impulse response start
+ * @param bcc base station color code - number of a training sequence
+ * @return first sample number of normal burst
+ */
+ int get_norm_chan_imp_resp(const gr_complex * input, gr_complex * chan_imp_resp, int bcc);
+
+ /**
+ *
+ */
+ void process_normal_burst(burst_counter burst_nr, const unsigned char * burst_binary);
+
+ /**
+ *
+ */
+ void configure_receiver();
+
+ public:
+ receiver_impl(feval_dd * tuner, int osr);
+ ~receiver_impl();
+
+ void forecast(int noutput_items, gr_vector_int &ninput_items_required);
+
+ // Where all the action really happens
+ int general_work(int noutput_items,
+ gr_vector_int &ninput_items,
+ gr_vector_const_void_star &input_items,
+ gr_vector_void_star &output_items);
+ };
+ } // namespace gsm
+} // namespace gr
+
+#endif /* INCLUDED_GSM_RECEIVER_IMPL_H */
+
diff --git a/lib/sch.c b/lib/sch.c
new file mode 100644
index 0000000..6f141dd
--- /dev/null
+++ b/lib/sch.c
@@ -0,0 +1,333 @@
+#include "system.h"
+#include <stdio.h>
+#include <stdlib.h>
+#include <unistd.h>
+#include <string.h>
+#include "gsm_constants.h"
+
+/*
+ * Synchronization channel.
+ *
+ * Timeslot Repeat length Frame Number (mod repeat length)
+ * 0 51 1, 11, 21, 31, 41
+ */
+
+/*
+ * Parity (FIRE) for the GSM SCH.
+ *
+ * g(x) = x^10 + x^8 + x^6 + x^5 + x^4 + x^2 + 1
+ */
+#define DATA_BLOCK_SIZE 25
+#define PARITY_SIZE 10
+#define TAIL_BITS_SIZE 4
+#define PARITY_OUTPUT_SIZE (DATA_BLOCK_SIZE + PARITY_SIZE + TAIL_BITS_SIZE)
+
+static const unsigned char parity_polynomial[PARITY_SIZE + 1] = {
+ 1, 0, 1, 0, 1, 1, 1, 0, 1, 0, 1
+};
+
+static const unsigned char parity_remainder[PARITY_SIZE] = {
+ 1, 1, 1, 1, 1, 1, 1, 1, 1, 1
+};
+
+
+static void parity_encode(unsigned char *d, unsigned char *p)
+{
+
+ unsigned int i;
+ unsigned char buf[DATA_BLOCK_SIZE + PARITY_SIZE], *q;
+
+ memcpy(buf, d, DATA_BLOCK_SIZE);
+ memset(buf + DATA_BLOCK_SIZE, 0, PARITY_SIZE);
+
+ for (q = buf; q < buf + DATA_BLOCK_SIZE; q++)
+ if (*q)
+ for (i = 0; i < PARITY_SIZE + 1; i++)
+ q[i] ^= parity_polynomial[i];
+ for (i = 0; i < PARITY_SIZE; i++)
+ p[i] = !buf[DATA_BLOCK_SIZE + i];
+}
+
+
+static int parity_check(unsigned char *d)
+{
+
+ unsigned int i;
+ unsigned char buf[DATA_BLOCK_SIZE + PARITY_SIZE], *q;
+
+ memcpy(buf, d, DATA_BLOCK_SIZE + PARITY_SIZE);
+
+ for (q = buf; q < buf + DATA_BLOCK_SIZE; q++)
+ if (*q)
+ for (i = 0; i < PARITY_SIZE + 1; i++)
+ q[i] ^= parity_polynomial[i];
+ return memcmp(buf + DATA_BLOCK_SIZE, parity_remainder, PARITY_SIZE);
+}
+
+
+/*
+ * Convolutional encoding and Viterbi decoding for the GSM SCH.
+ * (Equivalent to the GSM SACCH.)
+ *
+ * G_0 = 1 + x^3 + x^4
+ * G_1 = 1 + x + x^3 + x^4
+ *
+ * i.e.,
+ *
+ * c_{2k} = u_k + u_{k - 3} + u_{k - 4}
+ * c_{2k + 1} = u_k + u_{k - 1} + u_{k - 3} + u_{k - 4}
+ */
+#define CONV_INPUT_SIZE PARITY_OUTPUT_SIZE
+#define CONV_SIZE (2 * CONV_INPUT_SIZE)
+#define K 5
+#define MAX_ERROR (2 * CONV_INPUT_SIZE + 1)
+
+
+/*
+ * Given the current state and input bit, what are the output bits?
+ *
+ * encode[current_state][input_bit]
+ */
+static const unsigned int encode[1 << (K - 1)][2] = {
+ {0, 3}, {3, 0}, {3, 0}, {0, 3},
+ {0, 3}, {3, 0}, {3, 0}, {0, 3},
+ {1, 2}, {2, 1}, {2, 1}, {1, 2},
+ {1, 2}, {2, 1}, {2, 1}, {1, 2}
+};
+
+
+/*
+ * Given the current state and input bit, what is the next state?
+ *
+ * next_state[current_state][input_bit]
+ */
+static const unsigned int next_state[1 << (K - 1)][2] = {
+ {0, 8}, {0, 8}, {1, 9}, {1, 9},
+ {2, 10}, {2, 10}, {3, 11}, {3, 11},
+ {4, 12}, {4, 12}, {5, 13}, {5, 13},
+ {6, 14}, {6, 14}, {7, 15}, {7, 15}
+};
+
+
+/*
+ * Given the previous state and the current state, what input bit caused
+ * the transition? If it is impossible to transition between the two
+ * states, the value is 2.
+ *
+ * prev_next_state[previous_state][current_state]
+ */
+static const unsigned int prev_next_state[1 << (K - 1)][1 << (K - 1)] = {
+ { 0, 2, 2, 2, 2, 2, 2, 2, 1, 2, 2, 2, 2, 2, 2, 2},
+ { 0, 2, 2, 2, 2, 2, 2, 2, 1, 2, 2, 2, 2, 2, 2, 2},
+ { 2, 0, 2, 2, 2, 2, 2, 2, 2, 1, 2, 2, 2, 2, 2, 2},
+ { 2, 0, 2, 2, 2, 2, 2, 2, 2, 1, 2, 2, 2, 2, 2, 2},
+ { 2, 2, 0, 2, 2, 2, 2, 2, 2, 2, 1, 2, 2, 2, 2, 2},
+ { 2, 2, 0, 2, 2, 2, 2, 2, 2, 2, 1, 2, 2, 2, 2, 2},
+ { 2, 2, 2, 0, 2, 2, 2, 2, 2, 2, 2, 1, 2, 2, 2, 2},
+ { 2, 2, 2, 0, 2, 2, 2, 2, 2, 2, 2, 1, 2, 2, 2, 2},
+ { 2, 2, 2, 2, 0, 2, 2, 2, 2, 2, 2, 2, 1, 2, 2, 2},
+ { 2, 2, 2, 2, 0, 2, 2, 2, 2, 2, 2, 2, 1, 2, 2, 2},
+ { 2, 2, 2, 2, 2, 0, 2, 2, 2, 2, 2, 2, 2, 1, 2, 2},
+ { 2, 2, 2, 2, 2, 0, 2, 2, 2, 2, 2, 2, 2, 1, 2, 2},
+ { 2, 2, 2, 2, 2, 2, 0, 2, 2, 2, 2, 2, 2, 2, 1, 2},
+ { 2, 2, 2, 2, 2, 2, 0, 2, 2, 2, 2, 2, 2, 2, 1, 2},
+ { 2, 2, 2, 2, 2, 2, 2, 0, 2, 2, 2, 2, 2, 2, 2, 1},
+ { 2, 2, 2, 2, 2, 2, 2, 0, 2, 2, 2, 2, 2, 2, 2, 1}
+};
+
+
+static inline unsigned int hamming_distance2(unsigned int w)
+{
+
+ return (w & 1) + !!(w & 2);
+}
+
+
+static void conv_encode(unsigned char *data, unsigned char *output)
+{
+
+ unsigned int i, state = 0, o;
+
+ // encode data
+ for (i = 0; i < CONV_INPUT_SIZE; i++) {
+ o = encode[state][data[i]];
+ state = next_state[state][data[i]];
+ *output++ = !!(o & 2);
+ *output++ = o & 1;
+ }
+}
+
+
+static int conv_decode(unsigned char *data, unsigned char *output)
+{
+
+ int i, t;
+ unsigned int rdata, state, nstate, b, o, distance, accumulated_error,
+ min_state, min_error, cur_state;
+
+ unsigned int ae[1 << (K - 1)];
+ unsigned int nae[1 << (K - 1)]; // next accumulated error
+ unsigned int state_history[1 << (K - 1)][CONV_INPUT_SIZE + 1];
+
+ // initialize accumulated error, assume starting state is 0
+ for (i = 0; i < (1 << (K - 1)); i++)
+ ae[i] = nae[i] = MAX_ERROR;
+ ae[0] = 0;
+
+ // build trellis
+ for (t = 0; t < CONV_INPUT_SIZE; t++) {
+
+ // get received data symbol
+ rdata = (data[2 * t] << 1) | data[2 * t + 1];
+
+ // for each state
+ for (state = 0; state < (1 << (K - 1)); state++) {
+
+ // make sure this state is possible
+ if (ae[state] >= MAX_ERROR)
+ continue;
+
+ // find all states we lead to
+ for (b = 0; b < 2; b++) {
+
+ // get next state given input bit b
+ nstate = next_state[state][b];
+
+ // find output for this transition
+ o = encode[state][b];
+
+ // calculate distance from received data
+ distance = hamming_distance2(rdata ^ o);
+
+ // choose surviving path
+ accumulated_error = ae[state] + distance;
+ if (accumulated_error < nae[nstate]) {
+
+ // save error for surviving state
+ nae[nstate] = accumulated_error;
+
+ // update state history
+ state_history[nstate][t + 1] = state;
+ }
+ }
+ }
+
+ // get accumulated error ready for next time slice
+ for (i = 0; i < (1 << (K - 1)); i++) {
+ ae[i] = nae[i];
+ nae[i] = MAX_ERROR;
+ }
+ }
+
+ // the final state is the state with the fewest errors
+ min_state = (unsigned int) - 1;
+ min_error = MAX_ERROR;
+ for (i = 0; i < (1 << (K - 1)); i++) {
+ if (ae[i] < min_error) {
+ min_state = i;
+ min_error = ae[i];
+ }
+ }
+
+ // trace the path
+ cur_state = min_state;
+ for (t = CONV_INPUT_SIZE; t >= 1; t--) {
+ min_state = cur_state;
+ cur_state = state_history[cur_state][t]; // get previous
+ output[t - 1] = prev_next_state[cur_state][min_state];
+ }
+
+ // return the number of errors detected (hard-decision)
+ return min_error;
+}
+
+
+int decode_sch(const unsigned char *buf, int * t1_o, int * t2_o, int * t3_o, int * ncc_o, int * bcc_o)
+{
+
+ int errors, t1, t2, t3p, t3, ncc, bcc;
+ unsigned char data[CONV_SIZE], decoded_data[PARITY_OUTPUT_SIZE];
+
+ // extract encoded data from synchronization burst
+ /* buf, 39 bit */
+ /* buf + 39 + 64 = 103, 39 */
+ memcpy(data, buf, SCH_DATA_LEN);
+ memcpy(data + SCH_DATA_LEN, buf + SCH_DATA_LEN + N_SYNC_BITS, SCH_DATA_LEN);
+
+ // Viterbi decode
+ if (errors = conv_decode(data, decoded_data)) {
+ // fprintf(stderr, "error: sch: conv_decode (%d)\n", errors);
+ DEBUGF("ERR: conv_decode %d\n", errors);
+ return errors;
+ }
+
+ // check parity
+ if (parity_check(decoded_data)) {
+ // fprintf(stderr, "error: sch: parity failed\n");
+ DEBUGF("ERR: parity_check failed\n");
+ return 1;
+ }
+
+ // Synchronization channel information, 44.018 page 171. (V7.2.0)
+ ncc =
+ (decoded_data[ 7] << 2) |
+ (decoded_data[ 6] << 1) |
+ (decoded_data[ 5] << 0);
+ bcc =
+ (decoded_data[ 4] << 2) |
+ (decoded_data[ 3] << 1) |
+ (decoded_data[ 2] << 0);
+ t1 =
+ (decoded_data[ 1] << 10) |
+ (decoded_data[ 0] << 9) |
+ (decoded_data[15] << 8) |
+ (decoded_data[14] << 7) |
+ (decoded_data[13] << 6) |
+ (decoded_data[12] << 5) |
+ (decoded_data[11] << 4) |
+ (decoded_data[10] << 3) |
+ (decoded_data[ 9] << 2) |
+ (decoded_data[ 8] << 1) |
+ (decoded_data[23] << 0);
+ t2 =
+ (decoded_data[22] << 4) |
+ (decoded_data[21] << 3) |
+ (decoded_data[20] << 2) |
+ (decoded_data[19] << 1) |
+ (decoded_data[18] << 0);
+ t3p =
+ (decoded_data[17] << 2) |
+ (decoded_data[16] << 1) |
+ (decoded_data[24] << 0);
+
+ t3 = 10 * t3p + 1;
+
+ // modulo arithmetic t3 - t2 mod 26
+// tt = ((t3 + 26) - t2) % 26;
+
+// fn = (51 * 26 * t1) + (51 * tt) + t3;
+
+ /*
+ * BSIC: Base Station Identification Code
+ * BCC: Base station Color Code
+ * NCC: Network Color Code
+ *
+ * FN: Frame Number
+ */
+
+// printf("bsic: %x (bcc: %u; ncc: %u)\tFN: %u\n", bsic, bsic & 7,
+// (bsic >> 3) & 7, fn);
+
+// if (fn_o)
+// *fn_o = fn;
+// if (bsic_o)
+ if (t1_o && t2_o && t3_o && ncc_o && bcc_o) {
+ *t1_o = t1;
+ *t2_o = t2;
+ *t3_o = t3;
+ *bcc_o = bcc;
+ *ncc_o = ncc;
+ }
+
+ return 0;
+}
diff --git a/lib/sch.h b/lib/sch.h
new file mode 100644
index 0000000..11a6d06
--- /dev/null
+++ b/lib/sch.h
@@ -0,0 +1,19 @@
+
+#ifndef __SCH_H__
+#define __SCH_H__ 1
+
+#include <gsm/api.h>
+
+#ifdef __cplusplus
+extern "C"
+{
+#endif
+
+ GSM_API int decode_sch(const unsigned char *buf, int * t1_o, int * t2_o, int * t3_o, int * ncc, int * bcc);
+
+#ifdef __cplusplus
+}
+#endif
+
+#endif
+
diff --git a/lib/system.h b/lib/system.h
new file mode 100644
index 0000000..407751e
--- /dev/null
+++ b/lib/system.h
@@ -0,0 +1,10 @@
+#ifndef __GSMTVOID_SYSTEM_H__
+#define __GSMTVOID_SYSTEM_H__ 1
+
+#define DEBUGF(a...) { \
+ fprintf(stderr, "%s:%d ", __FILE__, __LINE__); \
+ fprintf(stderr, a); \
+} while (0)
+
+#endif
+
diff --git a/lib/test_gsm.cc b/lib/test_gsm.cc
new file mode 100644
index 0000000..3fdb4ad
--- /dev/null
+++ b/lib/test_gsm.cc
@@ -0,0 +1,47 @@
+/* -*- c++ -*- */
+/*
+ * Copyright 2012 Free Software Foundation, Inc.
+ *
+ * This file is part of GNU Radio
+ *
+ * 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.
+ */
+
+#ifdef HAVE_CONFIG_H
+#include "config.h"
+#endif
+
+#include <cppunit/TextTestRunner.h>
+#include <cppunit/XmlOutputter.h>
+
+#include <gnuradio/unittests.h>
+#include "qa_gsm.h"
+#include <iostream>
+
+int
+main (int argc, char **argv)
+{
+ CppUnit::TextTestRunner runner;
+ std::ofstream xmlfile(get_unittest_path("gsm.xml").c_str());
+ CppUnit::XmlOutputter *xmlout = new CppUnit::XmlOutputter(&runner.result(), xmlfile);
+
+ runner.addTest(qa_gsm::suite());
+ runner.setOutputter(xmlout);
+
+ bool was_successful = runner.run("", false);
+
+ return was_successful ? 0 : 1;
+}
diff --git a/lib/viterbi_detector.cc b/lib/viterbi_detector.cc
new file mode 100644
index 0000000..f3445cf
--- /dev/null
+++ b/lib/viterbi_detector.cc
@@ -0,0 +1,554 @@
+/* -*- c++ -*- */
+/*
+ * @file
+ * @author Piotr Krysik <pkrysik@stud.elka.pw.edu.pl>
+ * @section LICENSE
+ *
+ * This program 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.
+ *
+ * This program 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 this program; see the file COPYING. If not, write to
+ * the Free Software Foundation, Inc., 51 Franklin Street,
+ * Boston, MA 02110-1301, USA.
+ */
+
+/*
+ * viterbi_detector:
+ * This part does the detection of received sequnece.
+ * Employed algorithm is viterbi Maximum Likehood Sequence Estimation.
+ * At this moment it gives hard decisions on the output, but
+ * it was designed with soft decisions in mind.
+ *
+ * SYNTAX: void viterbi_detector(
+ * const gr_complex * input,
+ * unsigned int samples_num,
+ * gr_complex * rhh,
+ * unsigned int start_state,
+ * const unsigned int * stop_states,
+ * unsigned int stops_num,
+ * float * output)
+ *
+ * INPUT: input: Complex received signal afted matched filtering.
+ * samples_num: Number of samples in the input table.
+ * rhh: The autocorrelation of the estimated channel
+ * impulse response.
+ * start_state: Number of the start point. In GSM each burst
+ * starts with sequence of three bits (0,0,0) which
+ * indicates start point of the algorithm.
+ * stop_states: Table with numbers of possible stop states.
+ * stops_num: Number of possible stop states
+ *
+ *
+ * OUTPUT: output: Differentially decoded hard output of the algorithm:
+ * -1 for logical "0" and 1 for logical "1"
+ *
+ * SUB_FUNC: none
+ *
+ * TEST(S): Tested with real world normal burst.
+ */
+
+#include <gnuradio/gr_complex.h>
+#include <gsm_constants.h>
+#define PATHS_NUM (1 << (CHAN_IMP_RESP_LENGTH-1))
+
+void viterbi_detector(const gr_complex * input, unsigned int samples_num, gr_complex * rhh, unsigned int start_state, const unsigned int * stop_states, unsigned int stops_num, float * output)
+{
+ float increment[8];
+ float path_metrics1[16];
+ float path_metrics2[16];
+ float * new_path_metrics;
+ float * old_path_metrics;
+ float * tmp;
+ float trans_table[BURST_SIZE][16];
+ float pm_candidate1, pm_candidate2;
+ bool real_imag;
+ float input_symbol_real, input_symbol_imag;
+ unsigned int i, sample_nr;
+
+/*
+* Setup first path metrics, so only state pointed by start_state is possible.
+* Start_state metric is equal to zero, the rest is written with some very low value,
+* which makes them practically impossible to occur.
+*/
+ for(i=0; i<PATHS_NUM; i++){
+ path_metrics1[i]=(-10e30);
+ }
+ path_metrics1[start_state]=0;
+
+/*
+* Compute Increment - a table of values which does not change for subsequent input samples.
+* Increment is table of reference levels for computation of branch metrics:
+* branch metric = (+/-)received_sample (+/-) reference_level
+*/
+ increment[0] = -rhh[1].imag() -rhh[2].real() -rhh[3].imag() +rhh[4].real();
+ increment[1] = rhh[1].imag() -rhh[2].real() -rhh[3].imag() +rhh[4].real();
+ increment[2] = -rhh[1].imag() +rhh[2].real() -rhh[3].imag() +rhh[4].real();
+ increment[3] = rhh[1].imag() +rhh[2].real() -rhh[3].imag() +rhh[4].real();
+ increment[4] = -rhh[1].imag() -rhh[2].real() +rhh[3].imag() +rhh[4].real();
+ increment[5] = rhh[1].imag() -rhh[2].real() +rhh[3].imag() +rhh[4].real();
+ increment[6] = -rhh[1].imag() +rhh[2].real() +rhh[3].imag() +rhh[4].real();
+ increment[7] = rhh[1].imag() +rhh[2].real() +rhh[3].imag() +rhh[4].real();
+
+
+/*
+* Computation of path metrics and decisions (Add-Compare-Select).
+* It's composed of two parts: one for odd input samples (imaginary numbers)
+* and one for even samples (real numbers).
+* Each part is composed of independent (parallelisable) statements like
+* this one:
+* pm_candidate1 = old_path_metrics[0] - input_symbol_real - increment[7];
+* pm_candidate2 = old_path_metrics[8] - input_symbol_real + increment[0];
+* if(pm_candidate1 > pm_candidate2){
+* new_path_metrics[0] = pm_candidate1;
+* trans_table[sample_nr][0] = -1.0;
+* }
+* else{
+* new_path_metrics[0] = pm_candidate2;
+* trans_table[sample_nr][0] = 1.0;
+* }
+* This is very good point for optimisations (SIMD or OpenMP) as it's most time
+* consuming part of this function.
+*/
+ sample_nr=0;
+ old_path_metrics=path_metrics1;
+ new_path_metrics=path_metrics2;
+ while(sample_nr<samples_num){
+ //Processing imag states
+ real_imag=1;
+ input_symbol_imag = input[sample_nr].imag();
+
+ pm_candidate1 = old_path_metrics[0] + input_symbol_imag - increment[2];
+ pm_candidate2 = old_path_metrics[8] + input_symbol_imag + increment[5];
+ if(pm_candidate1 > pm_candidate2){
+ new_path_metrics[0] = pm_candidate1;
+ trans_table[sample_nr][0] = -1.0;
+ }
+ else{
+ new_path_metrics[0] = pm_candidate2;
+ trans_table[sample_nr][0] = 1.0;
+ }
+
+ pm_candidate1 = old_path_metrics[0] - input_symbol_imag + increment[2];
+ pm_candidate2 = old_path_metrics[8] - input_symbol_imag - increment[5];
+ if(pm_candidate1 > pm_candidate2){
+ new_path_metrics[1] = pm_candidate1;
+ trans_table[sample_nr][1] = -1.0;
+ }
+ else{
+ new_path_metrics[1] = pm_candidate2;
+ trans_table[sample_nr][1] = 1.0;
+ }
+
+ pm_candidate1 = old_path_metrics[1] + input_symbol_imag - increment[3];
+ pm_candidate2 = old_path_metrics[9] + input_symbol_imag + increment[4];
+ if(pm_candidate1 > pm_candidate2){
+ new_path_metrics[2] = pm_candidate1;
+ trans_table[sample_nr][2] = -1.0;
+ }
+ else{
+ new_path_metrics[2] = pm_candidate2;
+ trans_table[sample_nr][2] = 1.0;
+ }
+
+ pm_candidate1 = old_path_metrics[1] - input_symbol_imag + increment[3];
+ pm_candidate2 = old_path_metrics[9] - input_symbol_imag - increment[4];
+ if(pm_candidate1 > pm_candidate2){
+ new_path_metrics[3] = pm_candidate1;
+ trans_table[sample_nr][3] = -1.0;
+ }
+ else{
+ new_path_metrics[3] = pm_candidate2;
+ trans_table[sample_nr][3] = 1.0;
+ }
+
+ pm_candidate1 = old_path_metrics[2] + input_symbol_imag - increment[0];
+ pm_candidate2 = old_path_metrics[10] + input_symbol_imag + increment[7];
+ if(pm_candidate1 > pm_candidate2){
+ new_path_metrics[4] = pm_candidate1;
+ trans_table[sample_nr][4] = -1.0;
+ }
+ else{
+ new_path_metrics[4] = pm_candidate2;
+ trans_table[sample_nr][4] = 1.0;
+ }
+
+ pm_candidate1 = old_path_metrics[2] - input_symbol_imag + increment[0];
+ pm_candidate2 = old_path_metrics[10] - input_symbol_imag - increment[7];
+ if(pm_candidate1 > pm_candidate2){
+ new_path_metrics[5] = pm_candidate1;
+ trans_table[sample_nr][5] = -1.0;
+ }
+ else{
+ new_path_metrics[5] = pm_candidate2;
+ trans_table[sample_nr][5] = 1.0;
+ }
+
+ pm_candidate1 = old_path_metrics[3] + input_symbol_imag - increment[1];
+ pm_candidate2 = old_path_metrics[11] + input_symbol_imag + increment[6];
+ if(pm_candidate1 > pm_candidate2){
+ new_path_metrics[6] = pm_candidate1;
+ trans_table[sample_nr][6] = -1.0;
+ }
+ else{
+ new_path_metrics[6] = pm_candidate2;
+ trans_table[sample_nr][6] = 1.0;
+ }
+
+ pm_candidate1 = old_path_metrics[3] - input_symbol_imag + increment[1];
+ pm_candidate2 = old_path_metrics[11] - input_symbol_imag - increment[6];
+ if(pm_candidate1 > pm_candidate2){
+ new_path_metrics[7] = pm_candidate1;
+ trans_table[sample_nr][7] = -1.0;
+ }
+ else{
+ new_path_metrics[7] = pm_candidate2;
+ trans_table[sample_nr][7] = 1.0;
+ }
+
+ pm_candidate1 = old_path_metrics[4] + input_symbol_imag - increment[6];
+ pm_candidate2 = old_path_metrics[12] + input_symbol_imag + increment[1];
+ if(pm_candidate1 > pm_candidate2){
+ new_path_metrics[8] = pm_candidate1;
+ trans_table[sample_nr][8] = -1.0;
+ }
+ else{
+ new_path_metrics[8] = pm_candidate2;
+ trans_table[sample_nr][8] = 1.0;
+ }
+
+ pm_candidate1 = old_path_metrics[4] - input_symbol_imag + increment[6];
+ pm_candidate2 = old_path_metrics[12] - input_symbol_imag - increment[1];
+ if(pm_candidate1 > pm_candidate2){
+ new_path_metrics[9] = pm_candidate1;
+ trans_table[sample_nr][9] = -1.0;
+ }
+ else{
+ new_path_metrics[9] = pm_candidate2;
+ trans_table[sample_nr][9] = 1.0;
+ }
+
+ pm_candidate1 = old_path_metrics[5] + input_symbol_imag - increment[7];
+ pm_candidate2 = old_path_metrics[13] + input_symbol_imag + increment[0];
+ if(pm_candidate1 > pm_candidate2){
+ new_path_metrics[10] = pm_candidate1;
+ trans_table[sample_nr][10] = -1.0;
+ }
+ else{
+ new_path_metrics[10] = pm_candidate2;
+ trans_table[sample_nr][10] = 1.0;
+ }
+
+ pm_candidate1 = old_path_metrics[5] - input_symbol_imag + increment[7];
+ pm_candidate2 = old_path_metrics[13] - input_symbol_imag - increment[0];
+ if(pm_candidate1 > pm_candidate2){
+ new_path_metrics[11] = pm_candidate1;
+ trans_table[sample_nr][11] = -1.0;
+ }
+ else{
+ new_path_metrics[11] = pm_candidate2;
+ trans_table[sample_nr][11] = 1.0;
+ }
+
+ pm_candidate1 = old_path_metrics[6] + input_symbol_imag - increment[4];
+ pm_candidate2 = old_path_metrics[14] + input_symbol_imag + increment[3];
+ if(pm_candidate1 > pm_candidate2){
+ new_path_metrics[12] = pm_candidate1;
+ trans_table[sample_nr][12] = -1.0;
+ }
+ else{
+ new_path_metrics[12] = pm_candidate2;
+ trans_table[sample_nr][12] = 1.0;
+ }
+
+ pm_candidate1 = old_path_metrics[6] - input_symbol_imag + increment[4];
+ pm_candidate2 = old_path_metrics[14] - input_symbol_imag - increment[3];
+ if(pm_candidate1 > pm_candidate2){
+ new_path_metrics[13] = pm_candidate1;
+ trans_table[sample_nr][13] = -1.0;
+ }
+ else{
+ new_path_metrics[13] = pm_candidate2;
+ trans_table[sample_nr][13] = 1.0;
+ }
+
+ pm_candidate1 = old_path_metrics[7] + input_symbol_imag - increment[5];
+ pm_candidate2 = old_path_metrics[15] + input_symbol_imag + increment[2];
+ if(pm_candidate1 > pm_candidate2){
+ new_path_metrics[14] = pm_candidate1;
+ trans_table[sample_nr][14] = -1.0;
+ }
+ else{
+ new_path_metrics[14] = pm_candidate2;
+ trans_table[sample_nr][14] = 1.0;
+ }
+
+ pm_candidate1 = old_path_metrics[7] - input_symbol_imag + increment[5];
+ pm_candidate2 = old_path_metrics[15] - input_symbol_imag - increment[2];
+ if(pm_candidate1 > pm_candidate2){
+ new_path_metrics[15] = pm_candidate1;
+ trans_table[sample_nr][15] = -1.0;
+ }
+ else{
+ new_path_metrics[15] = pm_candidate2;
+ trans_table[sample_nr][15] = 1.0;
+ }
+ tmp=old_path_metrics;
+ old_path_metrics=new_path_metrics;
+ new_path_metrics=tmp;
+
+ sample_nr++;
+ if(sample_nr==samples_num)
+ break;
+
+ //Processing real states
+ real_imag=0;
+ input_symbol_real = input[sample_nr].real();
+
+ pm_candidate1 = old_path_metrics[0] - input_symbol_real - increment[7];
+ pm_candidate2 = old_path_metrics[8] - input_symbol_real + increment[0];
+ if(pm_candidate1 > pm_candidate2){
+ new_path_metrics[0] = pm_candidate1;
+ trans_table[sample_nr][0] = -1.0;
+ }
+ else{
+ new_path_metrics[0] = pm_candidate2;
+ trans_table[sample_nr][0] = 1.0;
+ }
+
+ pm_candidate1 = old_path_metrics[0] + input_symbol_real + increment[7];
+ pm_candidate2 = old_path_metrics[8] + input_symbol_real - increment[0];
+ if(pm_candidate1 > pm_candidate2){
+ new_path_metrics[1] = pm_candidate1;
+ trans_table[sample_nr][1] = -1.0;
+ }
+ else{
+ new_path_metrics[1] = pm_candidate2;
+ trans_table[sample_nr][1] = 1.0;
+ }
+
+ pm_candidate1 = old_path_metrics[1] - input_symbol_real - increment[6];
+ pm_candidate2 = old_path_metrics[9] - input_symbol_real + increment[1];
+ if(pm_candidate1 > pm_candidate2){
+ new_path_metrics[2] = pm_candidate1;
+ trans_table[sample_nr][2] = -1.0;
+ }
+ else{
+ new_path_metrics[2] = pm_candidate2;
+ trans_table[sample_nr][2] = 1.0;
+ }
+
+ pm_candidate1 = old_path_metrics[1] + input_symbol_real + increment[6];
+ pm_candidate2 = old_path_metrics[9] + input_symbol_real - increment[1];
+ if(pm_candidate1 > pm_candidate2){
+ new_path_metrics[3] = pm_candidate1;
+ trans_table[sample_nr][3] = -1.0;
+ }
+ else{
+ new_path_metrics[3] = pm_candidate2;
+ trans_table[sample_nr][3] = 1.0;
+ }
+
+ pm_candidate1 = old_path_metrics[2] - input_symbol_real - increment[5];
+ pm_candidate2 = old_path_metrics[10] - input_symbol_real + increment[2];
+ if(pm_candidate1 > pm_candidate2){
+ new_path_metrics[4] = pm_candidate1;
+ trans_table[sample_nr][4] = -1.0;
+ }
+ else{
+ new_path_metrics[4] = pm_candidate2;
+ trans_table[sample_nr][4] = 1.0;
+ }
+
+ pm_candidate1 = old_path_metrics[2] + input_symbol_real + increment[5];
+ pm_candidate2 = old_path_metrics[10] + input_symbol_real - increment[2];
+ if(pm_candidate1 > pm_candidate2){
+ new_path_metrics[5] = pm_candidate1;
+ trans_table[sample_nr][5] = -1.0;
+ }
+ else{
+ new_path_metrics[5] = pm_candidate2;
+ trans_table[sample_nr][5] = 1.0;
+ }
+
+ pm_candidate1 = old_path_metrics[3] - input_symbol_real - increment[4];
+ pm_candidate2 = old_path_metrics[11] - input_symbol_real + increment[3];
+ if(pm_candidate1 > pm_candidate2){
+ new_path_metrics[6] = pm_candidate1;
+ trans_table[sample_nr][6] = -1.0;
+ }
+ else{
+ new_path_metrics[6] = pm_candidate2;
+ trans_table[sample_nr][6] = 1.0;
+ }
+
+ pm_candidate1 = old_path_metrics[3] + input_symbol_real + increment[4];
+ pm_candidate2 = old_path_metrics[11] + input_symbol_real - increment[3];
+ if(pm_candidate1 > pm_candidate2){
+ new_path_metrics[7] = pm_candidate1;
+ trans_table[sample_nr][7] = -1.0;
+ }
+ else{
+ new_path_metrics[7] = pm_candidate2;
+ trans_table[sample_nr][7] = 1.0;
+ }
+
+ pm_candidate1 = old_path_metrics[4] - input_symbol_real - increment[3];
+ pm_candidate2 = old_path_metrics[12] - input_symbol_real + increment[4];
+ if(pm_candidate1 > pm_candidate2){
+ new_path_metrics[8] = pm_candidate1;
+ trans_table[sample_nr][8] = -1.0;
+ }
+ else{
+ new_path_metrics[8] = pm_candidate2;
+ trans_table[sample_nr][8] = 1.0;
+ }
+
+ pm_candidate1 = old_path_metrics[4] + input_symbol_real + increment[3];
+ pm_candidate2 = old_path_metrics[12] + input_symbol_real - increment[4];
+ if(pm_candidate1 > pm_candidate2){
+ new_path_metrics[9] = pm_candidate1;
+ trans_table[sample_nr][9] = -1.0;
+ }
+ else{
+ new_path_metrics[9] = pm_candidate2;
+ trans_table[sample_nr][9] = 1.0;
+ }
+
+ pm_candidate1 = old_path_metrics[5] - input_symbol_real - increment[2];
+ pm_candidate2 = old_path_metrics[13] - input_symbol_real + increment[5];
+ if(pm_candidate1 > pm_candidate2){
+ new_path_metrics[10] = pm_candidate1;
+ trans_table[sample_nr][10] = -1.0;
+ }
+ else{
+ new_path_metrics[10] = pm_candidate2;
+ trans_table[sample_nr][10] = 1.0;
+ }
+
+ pm_candidate1 = old_path_metrics[5] + input_symbol_real + increment[2];
+ pm_candidate2 = old_path_metrics[13] + input_symbol_real - increment[5];
+ if(pm_candidate1 > pm_candidate2){
+ new_path_metrics[11] = pm_candidate1;
+ trans_table[sample_nr][11] = -1.0;
+ }
+ else{
+ new_path_metrics[11] = pm_candidate2;
+ trans_table[sample_nr][11] = 1.0;
+ }
+
+ pm_candidate1 = old_path_metrics[6] - input_symbol_real - increment[1];
+ pm_candidate2 = old_path_metrics[14] - input_symbol_real + increment[6];
+ if(pm_candidate1 > pm_candidate2){
+ new_path_metrics[12] = pm_candidate1;
+ trans_table[sample_nr][12] = -1.0;
+ }
+ else{
+ new_path_metrics[12] = pm_candidate2;
+ trans_table[sample_nr][12] = 1.0;
+ }
+
+ pm_candidate1 = old_path_metrics[6] + input_symbol_real + increment[1];
+ pm_candidate2 = old_path_metrics[14] + input_symbol_real - increment[6];
+ if(pm_candidate1 > pm_candidate2){
+ new_path_metrics[13] = pm_candidate1;
+ trans_table[sample_nr][13] = -1.0;
+ }
+ else{
+ new_path_metrics[13] = pm_candidate2;
+ trans_table[sample_nr][13] = 1.0;
+ }
+
+ pm_candidate1 = old_path_metrics[7] - input_symbol_real - increment[0];
+ pm_candidate2 = old_path_metrics[15] - input_symbol_real + increment[7];
+ if(pm_candidate1 > pm_candidate2){
+ new_path_metrics[14] = pm_candidate1;
+ trans_table[sample_nr][14] = -1.0;
+ }
+ else{
+ new_path_metrics[14] = pm_candidate2;
+ trans_table[sample_nr][14] = 1.0;
+ }
+
+ pm_candidate1 = old_path_metrics[7] + input_symbol_real + increment[0];
+ pm_candidate2 = old_path_metrics[15] + input_symbol_real - increment[7];
+ if(pm_candidate1 > pm_candidate2){
+ new_path_metrics[15] = pm_candidate1;
+ trans_table[sample_nr][15] = -1.0;
+ }
+ else{
+ new_path_metrics[15] = pm_candidate2;
+ trans_table[sample_nr][15] = 1.0;
+ }
+ tmp=old_path_metrics;
+ old_path_metrics=new_path_metrics;
+ new_path_metrics=tmp;
+
+ sample_nr++;
+ }
+
+/*
+* Find the best from the stop states by comparing their path metrics.
+* Not every stop state is always possible, so we are searching in
+* a subset of them.
+*/
+ unsigned int best_stop_state;
+ float stop_state_metric, max_stop_state_metric;
+ best_stop_state = stop_states[0];
+ max_stop_state_metric = old_path_metrics[best_stop_state];
+ for(i=1; i< stops_num; i++){
+ stop_state_metric = old_path_metrics[stop_states[i]];
+ if(stop_state_metric > max_stop_state_metric){
+ max_stop_state_metric = stop_state_metric;
+ best_stop_state = stop_states[i];
+ }
+ }
+
+/*
+* This table was generated with hope that it gives a litle speedup during
+* traceback stage.
+* Received bit is related to the number of state in the trellis.
+* I've numbered states so their parity (number of ones) is related
+* to a received bit.
+*/
+ static const unsigned int parity_table[PATHS_NUM] = { 0, 1, 1, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 1, 1, 0, };
+
+/*
+* Table of previous states in the trellis diagram.
+* For GMSK modulation every state has two previous states.
+* Example:
+* previous_state_nr1 = prev_table[current_state_nr][0]
+* previous_state_nr2 = prev_table[current_state_nr][1]
+*/
+ static const unsigned int prev_table[PATHS_NUM][2] = { {0,8}, {0,8}, {1,9}, {1,9}, {2,10}, {2,10}, {3,11}, {3,11}, {4,12}, {4,12}, {5,13}, {5,13}, {6,14}, {6,14}, {7,15}, {7,15}, };
+
+/*
+* Traceback and differential decoding of received sequence.
+* Decisions stored in trans_table are used to restore best path in the trellis.
+*/
+ sample_nr=samples_num;
+ unsigned int state_nr=best_stop_state;
+ unsigned int decision;
+ bool out_bit=0;
+
+ while(sample_nr>0){
+ sample_nr--;
+ decision = (trans_table[sample_nr][state_nr]>0);
+
+ if(decision != out_bit)
+ output[sample_nr]=-trans_table[sample_nr][state_nr];
+ else
+ output[sample_nr]=trans_table[sample_nr][state_nr];
+
+ out_bit = out_bit ^ real_imag ^ parity_table[state_nr];
+ state_nr = prev_table[state_nr][decision];
+ real_imag = !real_imag;
+ }
+}
diff --git a/lib/viterbi_detector.h b/lib/viterbi_detector.h
new file mode 100644
index 0000000..0360f83
--- /dev/null
+++ b/lib/viterbi_detector.h
@@ -0,0 +1,63 @@
+/* -*- c++ -*- */
+/*
+ * @file
+ * @author Piotr Krysik <pkrysik@stud.elka.pw.edu.pl>
+ * @section LICENSE
+ *
+ * This program 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.
+ *
+ * This program 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 this program; see the file COPYING. If not, write to
+ * the Free Software Foundation, Inc., 51 Franklin Street,
+ * Boston, MA 02110-1301, USA.
+ */
+
+/*
+ * viterbi_detector:
+ * This part does the detection of received sequnece.
+ * Employed algorithm is viterbi Maximum Likehood Sequence Estimation.
+ * At this moment it gives hard decisions on the output, but
+ * it was designed with soft decisions in mind.
+ *
+ * SYNTAX: void viterbi_detector(
+ * const gr_complex * input,
+ * unsigned int samples_num,
+ * gr_complex * rhh,
+ * unsigned int start_state,
+ * const unsigned int * stop_states,
+ * unsigned int stops_num,
+ * float * output)
+ *
+ * INPUT: input: Complex received signal afted matched filtering.
+ * samples_num: Number of samples in the input table.
+ * rhh: The autocorrelation of the estimated channel
+ * impulse response.
+ * start_state: Number of the start point. In GSM each burst
+ * starts with sequence of three bits (0,0,0) which
+ * indicates start point of the algorithm.
+ * stop_states: Table with numbers of possible stop states.
+ * stops_num: Number of possible stop states
+ *
+ *
+ * OUTPUT: output: Differentially decoded hard output of the algorithm:
+ * -1 for logical "0" and 1 for logical "1"
+ *
+ * SUB_FUNC: none
+ *
+ * TEST(S): Tested with real world normal burst.
+ */
+
+#ifndef INCLUDED_VITERBI_DETECTOR_H
+#define INCLUDED_VITERBI_DETECTOR_H
+
+void viterbi_detector(const gr_complex * input, unsigned int samples_num, gr_complex * rhh, unsigned int start_state, const unsigned int * stop_states, unsigned int stops_num, float * output);
+
+#endif /* INCLUDED_VITERBI_DETECTOR_H */