From 5173378320687cf941f1463685efabb909fefb66 Mon Sep 17 00:00:00 2001 From: Rob Riggs Date: Sat, 26 Jun 2021 10:59:38 -0500 Subject: [PATCH] Use same Viterbi code as m17-cxx-demod. --- TNC/Util.h | 72 +++++++++++++++++++++++++++++++++++++++++++++------ TNC/Viterbi.h | 53 +++++++++++++++++++++++-------------- 2 files changed, 98 insertions(+), 27 deletions(-) diff --git a/TNC/Util.h b/TNC/Util.h index 1b476fb..614f360 100644 --- a/TNC/Util.h +++ b/TNC/Util.h @@ -30,19 +30,47 @@ constexpr std::bitset make_bitset(std::index_sequence, Tup return result; } -template -constexpr std::array>, (((1 << (LLR - 1)) - 1) * 6 ) + 1> make_llr_map() +/** + * This is the max value for the LLR based on size N. + */ +template +constexpr size_t llr_limit() { - constexpr size_t size = (((1 << (LLR - 1)) - 1) * 6 ) + 1; + return (1 << (N - 1)) - 1; +} + +/** + * There are (2^(N-1)-1) elements (E) per segment (e.g. N=4, E=7; N=3, E=3). + * These contain the LLR values 1..E. There are 6 segments in the LLR map: + * 1. (-Inf,-2] + * 2. (-2, -1] + * 3. (-1, 0] + * 4. (0, 1] + * 5. (1, 2] + * 6. (2, Inf) + * + * Note the slight asymmetry. This is OK as we are dealing with floats and + * it only matters to an epsilon of the float type. + */ +template +constexpr size_t llr_size() +{ + return llr_limit() * 6 + 1; +} + +template +constexpr std::array>, llr_size()> make_llr_map() +{ + constexpr size_t size = llr_size(); std::array>, size> result; - constexpr int8_t limit = (1 << (LLR - 1)) - 1; + constexpr int8_t limit = llr_limit(); constexpr FloatType inc = 1.0 / FloatType(limit); int8_t i = limit; int8_t j = limit; // Output must be ordered by k, ascending. - FloatType k = -3.0; + FloatType k = -3.0 + inc; for (size_t index = 0; index != size; ++index) { auto& a = result[index]; @@ -50,19 +78,22 @@ constexpr std::array>, (((1 << std::get<0>(std::get<1>(a)) = i; std::get<1>(std::get<1>(a)) = j; - if (k + 1.0 < inc / -2.0) + if (k + 1.0 < 0) { j--; + if (j == 0) j = -1; if (j < -limit) j = -limit; } - else if (k - 1.0 < inc / -2.0) + else if (k - 1.0 < 0) { i--; + if (i == 0) i = -1; if (i < -limit) i = -limit; } else { j++; + if (j == 0) j = 1; if (j > limit) j = limit; } k += inc; @@ -96,8 +127,10 @@ template auto llr(FloatType sample) { static constexpr auto symbol_map = detail::make_llr_map(); + static constexpr FloatType MAX_VALUE = 3.0; + static constexpr FloatType MIN_VALUE = -3.0; - FloatType s = std::min(3.0, std::max(-3.0, sample)); + FloatType s = std::min(MAX_VALUE, std::max(MIN_VALUE, sample)); auto it = std::lower_bound(symbol_map.begin(), symbol_map.end(), s, [](std::tuple> const& e, FloatType s){ @@ -179,6 +212,7 @@ template constexpr bool get_bit_index(const std::array& input, size_t index) { auto byte_index = index >> 3; + assert(byte_index < N); auto bit_index = 7 - (index & 7); return (input[byte_index] & (1 << bit_index)) >> bit_index; @@ -188,6 +222,7 @@ template void set_bit_index(std::array& input, size_t index) { auto byte_index = index >> 3; + assert(byte_index < N); auto bit_index = 7 - (index & 7); input[byte_index] |= (1 << bit_index); } @@ -196,6 +231,7 @@ template void reset_bit_index(std::array& input, size_t index) { auto byte_index = index >> 3; + assert(byte_index < N); auto bit_index = 7 - (index & 7); input[byte_index] &= ~(1 << bit_index); } @@ -259,4 +295,24 @@ constexpr auto to_byte_array(std::array in) return out; } +template +constexpr void to_byte_array(std::array in, std::array& out) +{ + size_t i = 0; + size_t b = 0; + uint8_t tmp = 0; + for (auto c : in) + { + tmp |= (c << (7 - b)); + if (++b == 8) + { + out[i] = tmp; + tmp = 0; + ++i; + b = 0; + } + } + if (i < out.size()) out[i] = tmp; +} + } // mobilinkd diff --git a/TNC/Viterbi.h b/TNC/Viterbi.h index bf964cc..a92d83f 100644 --- a/TNC/Viterbi.h +++ b/TNC/Viterbi.h @@ -7,8 +7,11 @@ #include "Util.h" #include +#include +#include +#include +#include #include -#include namespace mobilinkd { @@ -111,7 +114,10 @@ struct Viterbi metrics_t prevMetrics, currMetrics; - std::array, 244> history_storage_; + // This is the maximum amount of storage needed for M17. If used for + // other modes, this may need to be increased. This will never overflow + // because of a static assertion in the decode() function. + std::array, 244> history_; Viterbi(Trellis_ trellis) : cost_(makeCost(trellis)) @@ -119,7 +125,6 @@ struct Viterbi , prevState_(makePrevState(trellis)) {} - [[gnu::noinline]] void calculate_path_metric( const std::array& cost0, const std::array& cost1, @@ -129,8 +134,8 @@ struct Viterbi auto& i0 = nextState_[j][0]; auto& i1 = nextState_[j][1]; - int16_t c0 = cost0[j]; - int16_t c1 = cost1[j]; + auto& c0 = cost0[j]; + auto& c1 = cost1[j]; auto& p0 = prevMetrics[j]; auto& p1 = prevMetrics[j + NumStates / 2]; @@ -157,14 +162,15 @@ struct Viterbi template size_t decode(std::array const& in, std::array& out) { - static_assert(sizeof(history_storage_) >= IN / 2); + static_assert(sizeof(history_) >= IN / 2); constexpr auto MAX_METRIC = std::numeric_limits::max() / 2; prevMetrics.fill(MAX_METRIC); prevMetrics[0] = 0; // Starting point. - std::span history(history_storage_.begin(), history_storage_.begin() + IN / 2); + auto hbegin = history_.begin(); + auto hend = history_.begin() + IN / 2; constexpr size_t BUTTERFLY_SIZE = NumStates / 2; @@ -172,31 +178,39 @@ struct Viterbi std::array cost0; std::array cost1; - for (size_t i = 0; i != IN; i += 2) + for (size_t i = 0; i != IN; i += 2, hindex += 1) { - auto& hist = history[hindex]; int16_t s0 = in[i]; int16_t s1 = in[i + 1]; + cost0.fill(0); + cost1.fill(0); for (size_t j = 0; j != BUTTERFLY_SIZE; ++j) { - int16_t c = std::abs(cost_[j][0] - s0) + std::abs(cost_[j][1] - s1); - cost0[j] = c; - cost1[j] = METRIC - c; + if (s0) // is not erased + { + cost0[j] = std::abs(cost_[j][0] - s0); + cost1[j] = std::abs(cost_[j][0] + s0); + } + if (s1) // is not erased + { + cost0[j] += std::abs(cost_[j][1] - s1); + cost1[j] += std::abs(cost_[j][1] + s1); + } } for (size_t j = 0; j != BUTTERFLY_SIZE; ++j) { - calculate_path_metric(cost0, cost1, hist, j); + calculate_path_metric(cost0, cost1, history_[hindex], j); } std::swap(currMetrics, prevMetrics); - hindex += 1; } // Find starting point. Should be 0 for properly flushed CCs. - // However, 0 may not be the path with the least errors. + // However, 0 may not be the path with the fewest errors. size_t min_element = 0; int32_t min_cost = prevMetrics[0]; + for (size_t i = 0; i != NumStates; ++i) { if (prevMetrics[i] < min_cost) @@ -206,21 +220,22 @@ struct Viterbi } } - size_t ber = min_cost / (METRIC >> 1); // Cost is at least equal to # of erasures. + size_t cost = std::round(min_cost / float(detail::llr_limit())); // Do chainback. auto oit = std::rbegin(out); - auto hit = std::rbegin(history); + auto hit = std::make_reverse_iterator(hend); // rbegin + auto hrend = std::make_reverse_iterator(hbegin); // rend size_t next_element = min_element; size_t index = IN / 2; - while (oit != std::rend(out) && hit != std::rend(history)) + while (oit != std::rend(out) && hit != hrend) { auto v = (*hit++)[next_element]; if (index-- <= OUT) *oit++ = next_element & 1; next_element = prevState_[next_element][v]; } - return ber; + return cost; } };