Eric | ac726b1 | 2022-11-28 18:57:58 +0100 | [diff] [blame^] | 1 | /* -*- c++ -*- */ |
| 2 | /* |
| 3 | * @file |
| 4 | * @author (C) 2009-2017 by Piotr Krysik <ptrkrysik@gmail.com> |
| 5 | * @author Contributions by sysmocom - s.f.m.c. GmbH / Eric Wild <ewild@sysmocom.de> |
| 6 | * @section LICENSE |
| 7 | * |
| 8 | * Gr-gsm is free software; you can redistribute it and/or modify |
| 9 | * it under the terms of the GNU General Public License as published by |
| 10 | * the Free Software Foundation; either version 3, or (at your option) |
| 11 | * any later version. |
| 12 | * |
| 13 | * Gr-gsm is distributed in the hope that it will be useful, |
| 14 | * but WITHOUT ANY WARRANTY; without even the implied warranty of |
| 15 | * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the |
| 16 | * GNU General Public License for more details. |
| 17 | * |
| 18 | * You should have received a copy of the GNU General Public License |
| 19 | * along with gr-gsm; see the file COPYING. If not, write to |
| 20 | * the Free Software Foundation, Inc., 51 Franklin Street, |
| 21 | * Boston, MA 02110-1301, USA. |
| 22 | */ |
| 23 | |
| 24 | #include "constants.h" |
| 25 | |
| 26 | #ifdef HAVE_CONFIG_H |
| 27 | #include "config.h" |
| 28 | #endif |
| 29 | #include <complex> |
| 30 | |
| 31 | |
| 32 | #include <algorithm> |
| 33 | #include <string.h> |
| 34 | #include <iostream> |
| 35 | #include <numeric> |
| 36 | #include <vector> |
| 37 | #include <fstream> |
| 38 | |
| 39 | #include "viterbi_detector.h" |
| 40 | #include "grgsm_vitac.h" |
| 41 | |
| 42 | //signalVector mChanResp; |
| 43 | gr_complex d_sch_training_seq[N_SYNC_BITS]; ///<encoded training sequence of a SCH burst |
| 44 | gr_complex d_norm_training_seq[TRAIN_SEQ_NUM][N_TRAIN_BITS]; ///<encoded training sequences of a normal and dummy burst |
| 45 | const int d_chan_imp_length = CHAN_IMP_RESP_LENGTH; |
| 46 | |
| 47 | void initvita() { |
| 48 | |
| 49 | /** |
| 50 | * Prepare SCH sequence bits |
| 51 | * |
| 52 | * (TS_BITS + 2 * GUARD_PERIOD) |
| 53 | * Burst and two guard periods |
| 54 | * (one guard period is an arbitrary overlap) |
| 55 | */ |
| 56 | gmsk_mapper(SYNC_BITS, N_SYNC_BITS, |
| 57 | d_sch_training_seq, gr_complex(0.0, -1.0)); |
| 58 | for (auto &i : d_sch_training_seq) |
| 59 | i = conj(i); |
| 60 | |
| 61 | /* Prepare bits of training sequences */ |
| 62 | for (int i = 0; i < TRAIN_SEQ_NUM; i++) { |
| 63 | /** |
| 64 | * If first bit of the sequence is 0 |
| 65 | * => first symbol is 1, else -1 |
| 66 | */ |
| 67 | gr_complex startpoint = train_seq[i][0] == 0 ? |
| 68 | gr_complex(1.0, 0.0) : gr_complex(-1.0, 0.0); |
| 69 | gmsk_mapper(train_seq[i], N_TRAIN_BITS, |
| 70 | d_norm_training_seq[i], startpoint); |
| 71 | for (auto &i : d_norm_training_seq[i]) |
| 72 | i = conj(i); |
| 73 | } |
| 74 | |
| 75 | } |
| 76 | |
| 77 | MULTI_VER_TARGET_ATTR |
| 78 | void detect_burst(const gr_complex *input, gr_complex *chan_imp_resp, int burst_start, char *output_binary) |
| 79 | { |
| 80 | std::vector<gr_complex> rhh_temp(CHAN_IMP_RESP_LENGTH * d_OSR); |
| 81 | unsigned int stop_states[2] = { 4, 12 }; |
| 82 | gr_complex filtered_burst[BURST_SIZE]; |
| 83 | gr_complex rhh[CHAN_IMP_RESP_LENGTH]; |
| 84 | float output[BURST_SIZE]; |
| 85 | int start_state = 3; |
| 86 | |
| 87 | // if(burst_start < 0 ||burst_start > 10) |
| 88 | // fprintf(stderr, "bo %d\n", burst_start); |
| 89 | |
| 90 | // burst_start = burst_start >= 0 ? burst_start : 0; |
| 91 | |
| 92 | autocorrelation(chan_imp_resp, &rhh_temp[0], d_chan_imp_length * d_OSR); |
| 93 | for (int ii = 0; ii < d_chan_imp_length; ii++) |
| 94 | rhh[ii] = conj(rhh_temp[ii * d_OSR]); |
| 95 | |
| 96 | mafi(&input[burst_start], BURST_SIZE, chan_imp_resp, |
| 97 | d_chan_imp_length * d_OSR, filtered_burst); |
| 98 | |
| 99 | viterbi_detector(filtered_burst, BURST_SIZE, rhh, |
| 100 | start_state, stop_states, 2, output); |
| 101 | |
| 102 | for (int i = 0; i < BURST_SIZE; i++) |
| 103 | output_binary[i] = output[i] * -127; // pre flip bits! |
| 104 | } |
| 105 | |
| 106 | void |
| 107 | gmsk_mapper(const unsigned char* input, |
| 108 | int nitems, gr_complex* gmsk_output, gr_complex start_point) |
| 109 | { |
| 110 | gr_complex j = gr_complex(0.0, 1.0); |
| 111 | gmsk_output[0] = start_point; |
| 112 | |
| 113 | int previous_symbol = 2 * input[0] - 1; |
| 114 | int current_symbol; |
| 115 | int encoded_symbol; |
| 116 | |
| 117 | for (int i = 1; i < nitems; i++) { |
| 118 | /* Change bits representation to NRZ */ |
| 119 | current_symbol = 2 * input[i] - 1; |
| 120 | |
| 121 | /* Differentially encode */ |
| 122 | encoded_symbol = current_symbol * previous_symbol; |
| 123 | |
| 124 | /* And do GMSK mapping */ |
| 125 | gmsk_output[i] = j * gr_complex(encoded_symbol, 0.0) |
| 126 | * gmsk_output[i - 1]; |
| 127 | |
| 128 | previous_symbol = current_symbol; |
| 129 | } |
| 130 | } |
| 131 | |
| 132 | gr_complex |
| 133 | correlate_sequence(const gr_complex* sequence, |
| 134 | int length, const gr_complex* input) |
| 135 | { |
| 136 | gr_complex result(0.0, 0.0); |
| 137 | |
| 138 | for (int ii = 0; ii < length; ii++) |
| 139 | result += sequence[ii] * input[ii * d_OSR]; |
| 140 | |
| 141 | return conj(result) / gr_complex(length, 0); |
| 142 | } |
| 143 | |
| 144 | /* Computes autocorrelation for positive arguments */ |
| 145 | inline void |
| 146 | autocorrelation(const gr_complex* input, |
| 147 | gr_complex* out, int nitems) |
| 148 | { |
| 149 | for (int k = nitems - 1; k >= 0; k--) { |
| 150 | out[k] = gr_complex(0, 0); |
| 151 | for (int i = k; i < nitems; i++) |
| 152 | out[k] += input[i] * conj(input[i - k]); |
| 153 | } |
| 154 | } |
| 155 | |
| 156 | inline void |
| 157 | mafi(const gr_complex* input, int nitems, |
| 158 | gr_complex* filter, int filter_length, gr_complex* output) |
| 159 | { |
| 160 | for (int n = 0; n < nitems; n++) { |
| 161 | int a = n * d_OSR; |
| 162 | output[n] = 0; |
| 163 | |
| 164 | for (int ii = 0; ii < filter_length; ii++) { |
| 165 | if ((a + ii) >= nitems * d_OSR) |
| 166 | break; |
| 167 | |
| 168 | output[n] += input[a + ii] * filter[ii]; |
| 169 | } |
| 170 | } |
| 171 | } |
| 172 | |
| 173 | int get_chan_imp_resp(const gr_complex *input, gr_complex *chan_imp_resp, int search_center, int search_start_pos, |
| 174 | int search_stop_pos, gr_complex *tseq, int tseqlen, float *corr_max) |
| 175 | { |
| 176 | std::vector<gr_complex> correlation_buffer; |
| 177 | std::vector<float> window_energy_buffer; |
| 178 | std::vector<float> power_buffer; |
| 179 | |
| 180 | for (int ii = search_start_pos; ii < search_stop_pos; ii++) { |
| 181 | gr_complex correlation = correlate_sequence(tseq, tseqlen, &input[ii]); |
| 182 | correlation_buffer.push_back(correlation); |
| 183 | power_buffer.push_back(std::pow(abs(correlation), 2)); |
| 184 | } |
| 185 | |
| 186 | int strongest_corr_nr = max_element(power_buffer.begin(), power_buffer.end()) - power_buffer.begin(); |
| 187 | |
| 188 | /* Compute window energies */ |
| 189 | auto window_energy_start_offset = strongest_corr_nr - 6 * d_OSR; |
| 190 | window_energy_start_offset = window_energy_start_offset < 0 ? 0 : window_energy_start_offset; //can end up out of range.. |
| 191 | auto window_energy_end_offset = strongest_corr_nr + 6 * d_OSR + d_chan_imp_length * d_OSR; |
| 192 | auto iter = power_buffer.begin() + window_energy_start_offset; |
| 193 | auto iter_end = power_buffer.begin() + window_energy_end_offset; |
| 194 | while (iter != iter_end) { |
| 195 | std::vector<float>::iterator iter_ii = iter; |
| 196 | bool loop_end = false; |
| 197 | float energy = 0; |
| 198 | |
| 199 | int len = d_chan_imp_length * d_OSR; |
| 200 | for (int ii = 0; ii < len; ii++, iter_ii++) { |
| 201 | if (iter_ii == power_buffer.end()) { |
| 202 | loop_end = true; |
| 203 | break; |
| 204 | } |
| 205 | |
| 206 | energy += (*iter_ii); |
| 207 | } |
| 208 | |
| 209 | if (loop_end) |
| 210 | break; |
| 211 | |
| 212 | window_energy_buffer.push_back(energy); |
| 213 | iter++; |
| 214 | } |
| 215 | |
| 216 | /* Calculate the strongest window number */ |
| 217 | int strongest_window_nr = window_energy_start_offset + |
| 218 | max_element(window_energy_buffer.begin(), window_energy_buffer.end()) - |
| 219 | window_energy_buffer.begin(); |
| 220 | |
| 221 | // auto window_search_start = window_energy_buffer.begin() + strongest_corr_nr - 5* d_OSR; |
| 222 | // auto window_search_end = window_energy_buffer.begin() + strongest_corr_nr + 10* d_OSR; |
| 223 | // window_search_end = window_search_end >= window_energy_buffer.end() ? window_energy_buffer.end() : window_search_end; |
| 224 | |
| 225 | // /* Calculate the strongest window number */ |
| 226 | // int strongest_window_nr = max_element(window_search_start, window_search_end /* - d_chan_imp_length * d_OSR*/) - window_energy_buffer.begin(); |
| 227 | |
| 228 | // if (strongest_window_nr < 0) |
| 229 | // strongest_window_nr = 0; |
| 230 | |
| 231 | float max_correlation = 0; |
| 232 | for (int ii = 0; ii < d_chan_imp_length * d_OSR; ii++) { |
| 233 | gr_complex correlation = correlation_buffer[strongest_window_nr + ii]; |
| 234 | if (abs(correlation) > max_correlation) |
| 235 | max_correlation = abs(correlation); |
| 236 | chan_imp_resp[ii] = correlation; |
| 237 | } |
| 238 | |
| 239 | *corr_max = max_correlation; |
| 240 | |
| 241 | /** |
| 242 | * Compute first sample position, which corresponds |
| 243 | * to the first sample of the impulse response |
| 244 | */ |
| 245 | return search_start_pos + strongest_window_nr - search_center * d_OSR; |
| 246 | } |
| 247 | |
| 248 | /* |
| 249 | |
| 250 | 3 + 57 + 1 + 26 + 1 + 57 + 3 + 8.25 |
| 251 | |
| 252 | search center = 3 + 57 + 1 + 5 (due to tsc 5+16+5 split) |
| 253 | this is +-5 samples around (+5 beginning) of truncated t16 tsc |
| 254 | |
| 255 | */ |
| 256 | int get_norm_chan_imp_resp(const gr_complex *input, gr_complex *chan_imp_resp, float *corr_max, int bcc) |
| 257 | { |
| 258 | const int search_center = TRAIN_POS; |
| 259 | const int search_start_pos = (search_center - 5) * d_OSR + 1; |
| 260 | const int search_stop_pos = (search_center + 5 + d_chan_imp_length) * d_OSR; |
| 261 | const auto tseq = &d_norm_training_seq[bcc][TRAIN_BEGINNING]; |
| 262 | const auto tseqlen = N_TRAIN_BITS - (2 * TRAIN_BEGINNING); |
| 263 | return get_chan_imp_resp(input, chan_imp_resp, search_center, search_start_pos, search_stop_pos, tseq, tseqlen, |
| 264 | corr_max); |
| 265 | } |
| 266 | |
| 267 | /* |
| 268 | |
| 269 | 3 tail | 39 data | 64 tsc | 39 data | 3 tail | 8.25 guard |
| 270 | start 3+39 - 10 |
| 271 | end 3+39 + SYNC_SEARCH_RANGE |
| 272 | |
| 273 | */ |
| 274 | int get_sch_chan_imp_resp(const gr_complex *input, gr_complex *chan_imp_resp) |
| 275 | { |
| 276 | const int search_center = SYNC_POS + TRAIN_BEGINNING; |
| 277 | const int search_start_pos = (search_center - 10) * d_OSR; |
| 278 | const int search_stop_pos = (search_center + SYNC_SEARCH_RANGE) * d_OSR; |
| 279 | const auto tseq = &d_sch_training_seq[TRAIN_BEGINNING]; |
| 280 | const auto tseqlen = N_SYNC_BITS - (2 * TRAIN_BEGINNING); |
| 281 | |
| 282 | // strongest_window_nr + chan_imp_resp_center + SYNC_POS *d_OSR - 48 * d_OSR - 2 * d_OSR + 2 ; |
| 283 | float corr_max; |
| 284 | return get_chan_imp_resp(input, chan_imp_resp, search_center, search_start_pos, search_stop_pos, tseq, tseqlen, |
| 285 | &corr_max); |
| 286 | } |
| 287 | |
| 288 | int get_sch_buffer_chan_imp_resp(const gr_complex *input, gr_complex *chan_imp_resp, unsigned int len, float *corr_max) |
| 289 | { |
| 290 | const auto tseqlen = N_SYNC_BITS - (2 * TRAIN_BEGINNING); |
| 291 | const int search_center = SYNC_POS + TRAIN_BEGINNING; |
| 292 | const int search_start_pos = 0; |
| 293 | // FIXME: proper end offset |
| 294 | const int search_stop_pos = len - (N_SYNC_BITS*8); |
| 295 | auto tseq = &d_sch_training_seq[TRAIN_BEGINNING]; |
| 296 | |
| 297 | return get_chan_imp_resp(input, chan_imp_resp, search_center, search_start_pos, search_stop_pos, tseq, tseqlen, |
| 298 | corr_max); |
| 299 | } |