Hoernchen has uploaded this change for review. ( https://gerrit.osmocom.org/c/osmo-trx/+/30415 )
Change subject: vita demod by piotr krysik, modified ...................................................................... vita demod by piotr krysik, modified Had a few rounds of extensive cleanup. Uses gcc multiversioning for x86 targets. Change-Id: I5466c522cf4de984a4810ec46df43a10b52ed78f --- A Transceiver52M/grgsm_vitac/constants.h A Transceiver52M/grgsm_vitac/grgsm_vitac.cpp A Transceiver52M/grgsm_vitac/grgsm_vitac.h A Transceiver52M/grgsm_vitac/viterbi_detector.cc A Transceiver52M/grgsm_vitac/viterbi_detector.h 5 files changed, 959 insertions(+), 0 deletions(-) git pull ssh://gerrit.osmocom.org:29418/osmo-trx refs/changes/15/30415/1 diff --git a/Transceiver52M/grgsm_vitac/constants.h b/Transceiver52M/grgsm_vitac/constants.h new file mode 100644 index 0000000..0bcfd59 --- /dev/null +++ b/Transceiver52M/grgsm_vitac/constants.h @@ -0,0 +1,142 @@ +#pragma once +/* -*- c++ -*- */ +/* + * @file + * @author (C) 2009-2017 by Piotr Krysik <ptrkry...@gmail.com> + * @section LICENSE + * + * Gr-gsm 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. + * + * Gr-gsm 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 gr-gsm; see the file COPYING. If not, write to + * the Free Software Foundation, Inc., 51 Franklin Street, + * Boston, MA 02110-1301, USA. + */ + +#include <complex> + +#define gr_complex std::complex<float> + + +#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 ACCESS_BURST_SIZE 88 +#define PROCESSED_CHUNK BURST_SIZE+2*GUARD_PERIOD + +#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 (TAIL_BITS + 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 10 //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, normal_or_noise } 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 +}; diff --git a/Transceiver52M/grgsm_vitac/grgsm_vitac.cpp b/Transceiver52M/grgsm_vitac/grgsm_vitac.cpp new file mode 100644 index 0000000..aa2ab6f --- /dev/null +++ b/Transceiver52M/grgsm_vitac/grgsm_vitac.cpp @@ -0,0 +1,299 @@ +/* -*- c++ -*- */ +/* + * @file + * @author (C) 2009-2017 by Piotr Krysik <ptrkry...@gmail.com> + * @author Contributions by sysmocom - s.f.m.c. GmbH / Eric Wild <ew...@sysmocom.de> + * @section LICENSE + * + * Gr-gsm 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. + * + * Gr-gsm 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 gr-gsm; see the file COPYING. If not, write to + * the Free Software Foundation, Inc., 51 Franklin Street, + * Boston, MA 02110-1301, USA. + */ + +#include "constants.h" + +#ifdef HAVE_CONFIG_H +#include "config.h" +#endif +#include <complex> + + +#include <algorithm> +#include <string.h> +#include <iostream> +#include <numeric> +#include <vector> +#include <fstream> + +#include "viterbi_detector.h" +#include "grgsm_vitac.h" + +//signalVector mChanResp; +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 and dummy burst +const int d_chan_imp_length = CHAN_IMP_RESP_LENGTH; + +void initvita() { + + /** + * Prepare SCH sequence bits + * + * (TS_BITS + 2 * GUARD_PERIOD) + * Burst and two guard periods + * (one guard period is an arbitrary overlap) + */ + gmsk_mapper(SYNC_BITS, N_SYNC_BITS, + d_sch_training_seq, gr_complex(0.0, -1.0)); + for (auto &i : d_sch_training_seq) + i = conj(i); + + /* Prepare bits of training sequences */ + for (int i = 0; i < TRAIN_SEQ_NUM; i++) { + /** + * If first bit of the sequence is 0 + * => first symbol is 1, else -1 + */ + gr_complex startpoint = train_seq[i][0] == 0 ? + gr_complex(1.0, 0.0) : gr_complex(-1.0, 0.0); + gmsk_mapper(train_seq[i], N_TRAIN_BITS, + d_norm_training_seq[i], startpoint); + for (auto &i : d_norm_training_seq[i]) + i = conj(i); + } + +} + +MULTI_VER_TARGET_ATTR +void detect_burst(const gr_complex *input, gr_complex *chan_imp_resp, int burst_start, char *output_binary) +{ + std::vector<gr_complex> rhh_temp(CHAN_IMP_RESP_LENGTH * d_OSR); + unsigned int stop_states[2] = { 4, 12 }; + gr_complex filtered_burst[BURST_SIZE]; + gr_complex rhh[CHAN_IMP_RESP_LENGTH]; + float output[BURST_SIZE]; + int start_state = 3; + + // if(burst_start < 0 ||burst_start > 10) + // fprintf(stderr, "bo %d\n", burst_start); + + // burst_start = burst_start >= 0 ? burst_start : 0; + + autocorrelation(chan_imp_resp, &rhh_temp[0], 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] * -127; // pre flip bits! +} + +void +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); + gmsk_output[0] = start_point; + + int previous_symbol = 2 * input[0] - 1; + int current_symbol; + int encoded_symbol; + + 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; + } +} + +gr_complex +correlate_sequence(const gr_complex* sequence, + int length, const gr_complex* input) +{ + gr_complex result(0.0, 0.0); + + for (int ii = 0; ii < length; ii++) + result += sequence[ii] * input[ii * d_OSR]; + + return conj(result) / gr_complex(length, 0); +} + +/* Computes autocorrelation for positive arguments */ +inline void +autocorrelation(const gr_complex* input, + gr_complex* out, int nitems) +{ + for (int k = nitems - 1; k >= 0; k--) { + out[k] = gr_complex(0, 0); + for (int i = k; i < nitems; i++) + out[k] += input[i] * conj(input[i - k]); + } +} + +inline void +mafi(const gr_complex* input, int nitems, + gr_complex* filter, int filter_length, gr_complex* output) +{ + for (int n = 0; n < nitems; n++) { + int a = n * d_OSR; + output[n] = 0; + + for (int ii = 0; ii < filter_length; ii++) { + if ((a + ii) >= nitems * d_OSR) + break; + + output[n] += input[a + ii] * filter[ii]; + } + } +} + +int get_chan_imp_resp(const gr_complex *input, gr_complex *chan_imp_resp, int search_center, int search_start_pos, + int search_stop_pos, gr_complex *tseq, int tseqlen, float *corr_max) +{ + std::vector<gr_complex> correlation_buffer; + std::vector<float> window_energy_buffer; + std::vector<float> power_buffer; + + for (int ii = search_start_pos; ii < search_stop_pos; ii++) { + gr_complex correlation = correlate_sequence(tseq, tseqlen, &input[ii]); + correlation_buffer.push_back(correlation); + power_buffer.push_back(std::pow(abs(correlation), 2)); + } + + int strongest_corr_nr = max_element(power_buffer.begin(), power_buffer.end()) - power_buffer.begin(); + + /* Compute window energies */ + auto window_energy_start_offset = strongest_corr_nr - 6 * d_OSR; + window_energy_start_offset = window_energy_start_offset < 0 ? 0 : window_energy_start_offset; //can end up out of range.. + auto window_energy_end_offset = strongest_corr_nr + 6 * d_OSR + d_chan_imp_length * d_OSR; + auto iter = power_buffer.begin() + window_energy_start_offset; + auto iter_end = power_buffer.begin() + window_energy_end_offset; + while (iter != iter_end) { + std::vector<float>::iterator iter_ii = iter; + bool loop_end = false; + float energy = 0; + + int len = d_chan_imp_length * d_OSR; + for (int ii = 0; ii < len; ii++, iter_ii++) { + if (iter_ii == power_buffer.end()) { + loop_end = true; + break; + } + + energy += (*iter_ii); + } + + if (loop_end) + break; + + window_energy_buffer.push_back(energy); + iter++; + } + + /* Calculate the strongest window number */ + int strongest_window_nr = window_energy_start_offset + + max_element(window_energy_buffer.begin(), window_energy_buffer.end()) - + window_energy_buffer.begin(); + + // auto window_search_start = window_energy_buffer.begin() + strongest_corr_nr - 5* d_OSR; + // auto window_search_end = window_energy_buffer.begin() + strongest_corr_nr + 10* d_OSR; + // window_search_end = window_search_end >= window_energy_buffer.end() ? window_energy_buffer.end() : window_search_end; + + // /* Calculate the strongest window number */ + // int strongest_window_nr = max_element(window_search_start, window_search_end /* - d_chan_imp_length * d_OSR*/) - window_energy_buffer.begin(); + + // if (strongest_window_nr < 0) + // strongest_window_nr = 0; + + float 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) + max_correlation = abs(correlation); + chan_imp_resp[ii] = correlation; + } + + *corr_max = max_correlation; + + /** + * Compute first sample position, which corresponds + * to the first sample of the impulse response + */ + return search_start_pos + strongest_window_nr - search_center * d_OSR; +} + +/* + +3 + 57 + 1 + 26 + 1 + 57 + 3 + 8.25 + +search center = 3 + 57 + 1 + 5 (due to tsc 5+16+5 split) +this is +-5 samples around (+5 beginning) of truncated t16 tsc + +*/ +int get_norm_chan_imp_resp(const gr_complex *input, gr_complex *chan_imp_resp, float *corr_max, int bcc) +{ + const int search_center = TRAIN_POS; + const int search_start_pos = (search_center - 5) * d_OSR + 1; + const int search_stop_pos = (search_center + 5 + d_chan_imp_length) * d_OSR; + const auto tseq = &d_norm_training_seq[bcc][TRAIN_BEGINNING]; + const auto tseqlen = N_TRAIN_BITS - (2 * TRAIN_BEGINNING); + return get_chan_imp_resp(input, chan_imp_resp, search_center, search_start_pos, search_stop_pos, tseq, tseqlen, + corr_max); +} + +/* + +3 tail | 39 data | 64 tsc | 39 data | 3 tail | 8.25 guard +start 3+39 - 10 +end 3+39 + SYNC_SEARCH_RANGE + +*/ +int get_sch_chan_imp_resp(const gr_complex *input, gr_complex *chan_imp_resp) +{ + const int search_center = SYNC_POS + TRAIN_BEGINNING; + const int search_start_pos = (search_center - 10) * d_OSR; + const int search_stop_pos = (search_center + SYNC_SEARCH_RANGE) * d_OSR; + const auto tseq = &d_sch_training_seq[TRAIN_BEGINNING]; + const auto tseqlen = N_SYNC_BITS - (2 * TRAIN_BEGINNING); + + // strongest_window_nr + chan_imp_resp_center + SYNC_POS *d_OSR - 48 * d_OSR - 2 * d_OSR + 2 ; + float corr_max; + return get_chan_imp_resp(input, chan_imp_resp, search_center, search_start_pos, search_stop_pos, tseq, tseqlen, + &corr_max); +} + +int get_sch_buffer_chan_imp_resp(const gr_complex *input, gr_complex *chan_imp_resp, unsigned int len, float *corr_max) +{ + const auto tseqlen = N_SYNC_BITS - (2 * TRAIN_BEGINNING); + const int search_center = SYNC_POS + TRAIN_BEGINNING; + const int search_start_pos = 0; + // FIXME: proper end offset + const int search_stop_pos = len - (N_SYNC_BITS*8); + auto tseq = &d_sch_training_seq[TRAIN_BEGINNING]; + + return get_chan_imp_resp(input, chan_imp_resp, search_center, search_start_pos, search_stop_pos, tseq, tseqlen, + corr_max); +} \ No newline at end of file diff --git a/Transceiver52M/grgsm_vitac/grgsm_vitac.h b/Transceiver52M/grgsm_vitac/grgsm_vitac.h new file mode 100644 index 0000000..41fe4b6 --- /dev/null +++ b/Transceiver52M/grgsm_vitac/grgsm_vitac.h @@ -0,0 +1,62 @@ +#pragma once +/* -*- c++ -*- */ +/* + * @file + * @author (C) 2009-2017 by Piotr Krysik <ptrkry...@gmail.com> + * @section LICENSE + * + * Gr-gsm 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. + * + * Gr-gsm 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 gr-gsm; see the file COPYING. If not, write to + * the Free Software Foundation, Inc., 51 Franklin Street, + * Boston, MA 02110-1301, USA. + */ + +#include <vector> +#include "constants.h" + +#if defined(__has_attribute) +#if __has_attribute(target_clones) && defined(__x86_64) && true +#define MULTI_VER_TARGET_ATTR __attribute__((target_clones("avx", "sse4.2", "sse3", "sse2", "sse", "default"))) +#else +#define MULTI_VER_TARGET_ATTR +#endif +#endif + +#define SYNC_SEARCH_RANGE 30 +const int d_OSR(4); + +void initvita(); + +int process_vita_burst(gr_complex *input, int tsc, unsigned char *output_binary); +int process_vita_sc_burst(gr_complex *input, int tsc, unsigned char *output_binary, int *offset); + +MULTI_VER_TARGET_ATTR +void detect_burst(const gr_complex *input, gr_complex *chan_imp_resp, int burst_start, char *output_binary); +void gmsk_mapper(const unsigned char *input, int nitems, gr_complex *gmsk_output, gr_complex start_point); +gr_complex correlate_sequence(const gr_complex *sequence, int length, const gr_complex *input); +inline void autocorrelation(const gr_complex *input, gr_complex *out, int nitems); +inline void mafi(const gr_complex *input, int nitems, gr_complex *filter, int filter_length, gr_complex *output); +int get_sch_chan_imp_resp(const gr_complex *input, gr_complex *chan_imp_resp); +int get_norm_chan_imp_resp(const gr_complex *input, gr_complex *chan_imp_resp, float *corr_max, int bcc); +int get_sch_buffer_chan_imp_resp(const gr_complex *input, gr_complex *chan_imp_resp, unsigned int len, float *corr_max); + +enum class btype { NB, SCH }; +struct fdata { + btype t; + unsigned int fn; + int tn; + int bcc; + std::string fpath; + std::vector<gr_complex> data; + unsigned int data_start_offset; +}; \ No newline at end of file diff --git a/Transceiver52M/grgsm_vitac/viterbi_detector.cc b/Transceiver52M/grgsm_vitac/viterbi_detector.cc new file mode 100644 index 0000000..beee386 --- /dev/null +++ b/Transceiver52M/grgsm_vitac/viterbi_detector.cc @@ -0,0 +1,392 @@ +/* -*- c++ -*- */ +/* + * @file + * @author (C) 2009 by Piotr Krysik <ptrkry...@gmail.com> + * @section LICENSE + * + * Gr-gsm 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. + * + * Gr-gsm 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 gr-gsm; 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 "constants.h" +#include <cmath> + +#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 paths_difference; + 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_imag +increment[2]; +* pm_candidate2 = old_path_metrics[8] -input_symbol_imag -increment[5]; +* paths_difference=pm_candidate2-pm_candidate1; +* new_path_metrics[1]=(paths_difference<0) ? pm_candidate1 : pm_candidate2; +* trans_table[sample_nr][1] = paths_difference; +* 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]; + paths_difference=pm_candidate2-pm_candidate1; + new_path_metrics[0]=(paths_difference<0) ? pm_candidate1 : pm_candidate2; + trans_table[sample_nr][0] = paths_difference; + + pm_candidate1 = old_path_metrics[0] -input_symbol_imag +increment[2]; + pm_candidate2 = old_path_metrics[8] -input_symbol_imag -increment[5]; + paths_difference=pm_candidate2-pm_candidate1; + new_path_metrics[1]=(paths_difference<0) ? pm_candidate1 : pm_candidate2; + trans_table[sample_nr][1] = paths_difference; + + pm_candidate1 = old_path_metrics[1] +input_symbol_imag -increment[3]; + pm_candidate2 = old_path_metrics[9] +input_symbol_imag +increment[4]; + paths_difference=pm_candidate2-pm_candidate1; + new_path_metrics[2]=(paths_difference<0) ? pm_candidate1 : pm_candidate2; + trans_table[sample_nr][2] = paths_difference; + + pm_candidate1 = old_path_metrics[1] -input_symbol_imag +increment[3]; + pm_candidate2 = old_path_metrics[9] -input_symbol_imag -increment[4]; + paths_difference=pm_candidate2-pm_candidate1; + new_path_metrics[3]=(paths_difference<0) ? pm_candidate1 : pm_candidate2; + trans_table[sample_nr][3] = paths_difference; + + pm_candidate1 = old_path_metrics[2] +input_symbol_imag -increment[0]; + pm_candidate2 = old_path_metrics[10] +input_symbol_imag +increment[7]; + paths_difference=pm_candidate2-pm_candidate1; + new_path_metrics[4]=(paths_difference<0) ? pm_candidate1 : pm_candidate2; + trans_table[sample_nr][4] = paths_difference; + + pm_candidate1 = old_path_metrics[2] -input_symbol_imag +increment[0]; + pm_candidate2 = old_path_metrics[10] -input_symbol_imag -increment[7]; + paths_difference=pm_candidate2-pm_candidate1; + new_path_metrics[5]=(paths_difference<0) ? pm_candidate1 : pm_candidate2; + trans_table[sample_nr][5] = paths_difference; + + pm_candidate1 = old_path_metrics[3] +input_symbol_imag -increment[1]; + pm_candidate2 = old_path_metrics[11] +input_symbol_imag +increment[6]; + paths_difference=pm_candidate2-pm_candidate1; + new_path_metrics[6]=(paths_difference<0) ? pm_candidate1 : pm_candidate2; + trans_table[sample_nr][6] = paths_difference; + + pm_candidate1 = old_path_metrics[3] -input_symbol_imag +increment[1]; + pm_candidate2 = old_path_metrics[11] -input_symbol_imag -increment[6]; + paths_difference=pm_candidate2-pm_candidate1; + new_path_metrics[7]=(paths_difference<0) ? pm_candidate1 : pm_candidate2; + trans_table[sample_nr][7] = paths_difference; + + pm_candidate1 = old_path_metrics[4] +input_symbol_imag -increment[6]; + pm_candidate2 = old_path_metrics[12] +input_symbol_imag +increment[1]; + paths_difference=pm_candidate2-pm_candidate1; + new_path_metrics[8]=(paths_difference<0) ? pm_candidate1 : pm_candidate2; + trans_table[sample_nr][8] = paths_difference; + + pm_candidate1 = old_path_metrics[4] -input_symbol_imag +increment[6]; + pm_candidate2 = old_path_metrics[12] -input_symbol_imag -increment[1]; + paths_difference=pm_candidate2-pm_candidate1; + new_path_metrics[9]=(paths_difference<0) ? pm_candidate1 : pm_candidate2; + trans_table[sample_nr][9] = paths_difference; + + pm_candidate1 = old_path_metrics[5] +input_symbol_imag -increment[7]; + pm_candidate2 = old_path_metrics[13] +input_symbol_imag +increment[0]; + paths_difference=pm_candidate2-pm_candidate1; + new_path_metrics[10]=(paths_difference<0) ? pm_candidate1 : pm_candidate2; + trans_table[sample_nr][10] = paths_difference; + + pm_candidate1 = old_path_metrics[5] -input_symbol_imag +increment[7]; + pm_candidate2 = old_path_metrics[13] -input_symbol_imag -increment[0]; + paths_difference=pm_candidate2-pm_candidate1; + new_path_metrics[11]=(paths_difference<0) ? pm_candidate1 : pm_candidate2; + trans_table[sample_nr][11] = paths_difference; + + pm_candidate1 = old_path_metrics[6] +input_symbol_imag -increment[4]; + pm_candidate2 = old_path_metrics[14] +input_symbol_imag +increment[3]; + paths_difference=pm_candidate2-pm_candidate1; + new_path_metrics[12]=(paths_difference<0) ? pm_candidate1 : pm_candidate2; + trans_table[sample_nr][12] = paths_difference; + + pm_candidate1 = old_path_metrics[6] -input_symbol_imag +increment[4]; + pm_candidate2 = old_path_metrics[14] -input_symbol_imag -increment[3]; + paths_difference=pm_candidate2-pm_candidate1; + new_path_metrics[13]=(paths_difference<0) ? pm_candidate1 : pm_candidate2; + trans_table[sample_nr][13] = paths_difference; + + pm_candidate1 = old_path_metrics[7] +input_symbol_imag -increment[5]; + pm_candidate2 = old_path_metrics[15] +input_symbol_imag +increment[2]; + paths_difference=pm_candidate2-pm_candidate1; + new_path_metrics[14]=(paths_difference<0) ? pm_candidate1 : pm_candidate2; + trans_table[sample_nr][14] = paths_difference; + + pm_candidate1 = old_path_metrics[7] -input_symbol_imag +increment[5]; + pm_candidate2 = old_path_metrics[15] -input_symbol_imag -increment[2]; + paths_difference=pm_candidate2-pm_candidate1; + new_path_metrics[15]=(paths_difference<0) ? pm_candidate1 : pm_candidate2; + trans_table[sample_nr][15] = paths_difference; + 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]; + paths_difference=pm_candidate2-pm_candidate1; + new_path_metrics[0]=(paths_difference<0) ? pm_candidate1 : pm_candidate2; + trans_table[sample_nr][0] = paths_difference; + + pm_candidate1 = old_path_metrics[0] +input_symbol_real +increment[7]; + pm_candidate2 = old_path_metrics[8] +input_symbol_real -increment[0]; + paths_difference=pm_candidate2-pm_candidate1; + new_path_metrics[1]=(paths_difference<0) ? pm_candidate1 : pm_candidate2; + trans_table[sample_nr][1] = paths_difference; + + pm_candidate1 = old_path_metrics[1] -input_symbol_real -increment[6]; + pm_candidate2 = old_path_metrics[9] -input_symbol_real +increment[1]; + paths_difference=pm_candidate2-pm_candidate1; + new_path_metrics[2]=(paths_difference<0) ? pm_candidate1 : pm_candidate2; + trans_table[sample_nr][2] = paths_difference; + + pm_candidate1 = old_path_metrics[1] +input_symbol_real +increment[6]; + pm_candidate2 = old_path_metrics[9] +input_symbol_real -increment[1]; + paths_difference=pm_candidate2-pm_candidate1; + new_path_metrics[3]=(paths_difference<0) ? pm_candidate1 : pm_candidate2; + trans_table[sample_nr][3] = paths_difference; + + pm_candidate1 = old_path_metrics[2] -input_symbol_real -increment[5]; + pm_candidate2 = old_path_metrics[10] -input_symbol_real +increment[2]; + paths_difference=pm_candidate2-pm_candidate1; + new_path_metrics[4]=(paths_difference<0) ? pm_candidate1 : pm_candidate2; + trans_table[sample_nr][4] = paths_difference; + + pm_candidate1 = old_path_metrics[2] +input_symbol_real +increment[5]; + pm_candidate2 = old_path_metrics[10] +input_symbol_real -increment[2]; + paths_difference=pm_candidate2-pm_candidate1; + new_path_metrics[5]=(paths_difference<0) ? pm_candidate1 : pm_candidate2; + trans_table[sample_nr][5] = paths_difference; + + pm_candidate1 = old_path_metrics[3] -input_symbol_real -increment[4]; + pm_candidate2 = old_path_metrics[11] -input_symbol_real +increment[3]; + paths_difference=pm_candidate2-pm_candidate1; + new_path_metrics[6]=(paths_difference<0) ? pm_candidate1 : pm_candidate2; + trans_table[sample_nr][6] = paths_difference; + + pm_candidate1 = old_path_metrics[3] +input_symbol_real +increment[4]; + pm_candidate2 = old_path_metrics[11] +input_symbol_real -increment[3]; + paths_difference=pm_candidate2-pm_candidate1; + new_path_metrics[7]=(paths_difference<0) ? pm_candidate1 : pm_candidate2; + trans_table[sample_nr][7] = paths_difference; + + pm_candidate1 = old_path_metrics[4] -input_symbol_real -increment[3]; + pm_candidate2 = old_path_metrics[12] -input_symbol_real +increment[4]; + paths_difference=pm_candidate2-pm_candidate1; + new_path_metrics[8]=(paths_difference<0) ? pm_candidate1 : pm_candidate2; + trans_table[sample_nr][8] = paths_difference; + + pm_candidate1 = old_path_metrics[4] +input_symbol_real +increment[3]; + pm_candidate2 = old_path_metrics[12] +input_symbol_real -increment[4]; + paths_difference=pm_candidate2-pm_candidate1; + new_path_metrics[9]=(paths_difference<0) ? pm_candidate1 : pm_candidate2; + trans_table[sample_nr][9] = paths_difference; + + pm_candidate1 = old_path_metrics[5] -input_symbol_real -increment[2]; + pm_candidate2 = old_path_metrics[13] -input_symbol_real +increment[5]; + paths_difference=pm_candidate2-pm_candidate1; + new_path_metrics[10]=(paths_difference<0) ? pm_candidate1 : pm_candidate2; + trans_table[sample_nr][10] = paths_difference; + + pm_candidate1 = old_path_metrics[5] +input_symbol_real +increment[2]; + pm_candidate2 = old_path_metrics[13] +input_symbol_real -increment[5]; + paths_difference=pm_candidate2-pm_candidate1; + new_path_metrics[11]=(paths_difference<0) ? pm_candidate1 : pm_candidate2; + trans_table[sample_nr][11] = paths_difference; + + pm_candidate1 = old_path_metrics[6] -input_symbol_real -increment[1]; + pm_candidate2 = old_path_metrics[14] -input_symbol_real +increment[6]; + paths_difference=pm_candidate2-pm_candidate1; + new_path_metrics[12]=(paths_difference<0) ? pm_candidate1 : pm_candidate2; + trans_table[sample_nr][12] = paths_difference; + + pm_candidate1 = old_path_metrics[6] +input_symbol_real +increment[1]; + pm_candidate2 = old_path_metrics[14] +input_symbol_real -increment[6]; + paths_difference=pm_candidate2-pm_candidate1; + new_path_metrics[13]=(paths_difference<0) ? pm_candidate1 : pm_candidate2; + trans_table[sample_nr][13] = paths_difference; + + pm_candidate1 = old_path_metrics[7] -input_symbol_real -increment[0]; + pm_candidate2 = old_path_metrics[15] -input_symbol_real +increment[7]; + paths_difference=pm_candidate2-pm_candidate1; + new_path_metrics[14]=(paths_difference<0) ? pm_candidate1 : pm_candidate2; + trans_table[sample_nr][14] = paths_difference; + + pm_candidate1 = old_path_metrics[7] +input_symbol_real +increment[0]; + pm_candidate2 = old_path_metrics[15] +input_symbol_real -increment[7]; + paths_difference=pm_candidate2-pm_candidate1; + new_path_metrics[15]=(paths_difference<0) ? pm_candidate1 : pm_candidate2; + trans_table[sample_nr][15] = paths_difference; + + 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/Transceiver52M/grgsm_vitac/viterbi_detector.h b/Transceiver52M/grgsm_vitac/viterbi_detector.h new file mode 100644 index 0000000..2595343 --- /dev/null +++ b/Transceiver52M/grgsm_vitac/viterbi_detector.h @@ -0,0 +1,64 @@ +/* -*- c++ -*- */ +/* + * @file + * @author (C) 2009 Piotr Krysik <ptrkry...@gmail.com> + * @section LICENSE + * + * Gr-gsm 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. + * + * Gr-gsm 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 gr-gsm; 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 +#include "constants.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 */ -- To view, visit https://gerrit.osmocom.org/c/osmo-trx/+/30415 To unsubscribe, or for help writing mail filters, visit https://gerrit.osmocom.org/settings Gerrit-Project: osmo-trx Gerrit-Branch: master Gerrit-Change-Id: I5466c522cf4de984a4810ec46df43a10b52ed78f Gerrit-Change-Number: 30415 Gerrit-PatchSet: 1 Gerrit-Owner: Hoernchen <ew...@sysmocom.de> Gerrit-MessageType: newchange