From b857741d9d1b0131668eb6b321be7952d32d1a34 Mon Sep 17 00:00:00 2001 From: Ahmet Inan Date: Thu, 9 Feb 2023 23:50:43 +0100 Subject: [PATCH 01/22] experimenting with a scheme for Ribbit --- README.md | 4 +- decode.cc | 275 ++++++++++++------------------------------ encode.cc | 308 +++++++++--------------------------------------- freezer.cc | 4 +- polar_tables.hh | 4 +- 5 files changed, 136 insertions(+), 459 deletions(-) diff --git a/README.md b/README.md index 4de45f3..a02dd24 100644 --- a/README.md +++ b/README.md @@ -3,10 +3,10 @@ Quick start: -Create file ```uncoded.dat``` with ```1360``` bits of random data: +Create file ```uncoded.dat``` with ```2048``` bits of random data: ``` -dd if=/dev/urandom of=uncoded.dat bs=1 count=170 +dd if=/dev/urandom of=uncoded.dat bs=1 count=256 ``` Encode file ```uncoded.dat``` to ```encoded.wav``` [WAV](https://en.wikipedia.org/wiki/WAV) audio file with 8000 Hz sample rate, 16 bits and only 1 (real) channel: diff --git a/decode.cc b/decode.cc index d7e4aee..ebdec08 100644 --- a/decode.cc +++ b/decode.cc @@ -1,7 +1,7 @@ /* OFDM modem decoder -Copyright 2021 Ahmet Inan +Copyright 2023 Ahmet Inan */ #include @@ -10,7 +10,6 @@ Copyright 2021 Ahmet Inan #include namespace DSP { using std::abs; using std::min; using std::cos; using std::sin; } #include "bip_buffer.hh" -#include "theil_sen.hh" #include "xorshift.hh" #include "trigger.hh" #include "complex.hh" @@ -26,7 +25,6 @@ namespace DSP { using std::abs; using std::min; using std::cos; using std::sin; #include "fft.hh" #include "mls.hh" #include "crc.hh" -#include "osd.hh" #include "psk.hh" #include "polar_tables.hh" #include "polar_helper.hh" @@ -42,11 +40,12 @@ struct SchmidlCox DSP::FastFourierTransform fwd; DSP::FastFourierTransform bwd; DSP::SMA4 cor; - DSP::SMA4 pwr; + DSP::SMA4 pwr; DSP::SMA4 match; DSP::Delay delay; DSP::SchmittTrigger threshold; DSP::FallingEdgeTrigger falling; + DSP::Phasor fsh; cmplx tmp0[symbol_len], tmp1[symbol_len], tmp2[symbol_len]; cmplx seq[symbol_len], kern[symbol_len]; cmplx cmplx_shift = 0; @@ -72,8 +71,11 @@ public: value cfo_rad = 0; value frac_cfo = 0; - SchmidlCox(const cmplx *sequence) : threshold(value(0.17*match_len), value(0.19*match_len)) + SchmidlCox(const cmplx *sequence) : threshold(value(0.2*match_len), value(0.3*match_len)) { +#if 0 + fsh.omega(-1, symbol_len); +#endif for (int i = 0; i < symbol_len; ++i) seq[i] = sequence[i]; fwd(kern, sequence); @@ -82,8 +84,8 @@ public: } bool operator()(const cmplx *samples) { - cmplx P = cor(samples[search_pos+symbol_len] * conj(samples[search_pos+2*symbol_len])); - value R = value(0.5) * pwr(norm(samples[search_pos+2*symbol_len])); + cmplx P = cor(fsh() * samples[search_pos] * conj(samples[search_pos+symbol_len+guard_len])); + value R = value(0.5) * pwr(norm(samples[search_pos]) + norm(samples[search_pos+symbol_len+guard_len])); value min_R = 0.0001 * symbol_len; R = std::max(R, min_R); value timing = match(norm(P) / (R * R)); @@ -92,6 +94,8 @@ public: bool collect = threshold(timing); bool process = falling(collect); + //std::cout << timing << " " << phase << " " << process << std::endl; + if (!collect && !process) return false; @@ -101,6 +105,10 @@ public: index_max = match_del; } else if (index_max < symbol_len + guard_len + match_del) { ++index_max; + } else if (process) { + index_max = 0; + timing_max = 0; + return false; } if (!process) @@ -110,11 +118,11 @@ public: DSP::Phasor osc; osc.omega(frac_cfo); - symbol_pos = search_pos - index_max; + symbol_pos = search_pos - index_max + guard_len; index_max = 0; timing_max = 0; for (int i = 0; i < symbol_len; ++i) - tmp1[i] = samples[i+symbol_pos+symbol_len] * osc(); + tmp1[i] = samples[i+symbol_pos] * osc(); fwd(tmp0, tmp1); for (int i = 0; i < symbol_len; ++i) tmp1[i] = demod_or_erase(tmp0[i], tmp0[bin(i-1)]); @@ -151,15 +159,10 @@ public: } }; -void base37_decoder(char *str, long long int val, int len) -{ - for (int i = len-1; i >= 0; --i, val /= 37) - str[i] = " 0123456789ABCDEFGHIJKLMNOPQRSTUVWXYZ"[val%37]; -} - -template struct Decoder { + typedef float value; + typedef DSP::Complex cmplx; typedef int8_t code_type; #ifdef __AVX2__ typedef SIMD mesg_type; @@ -167,23 +170,19 @@ struct Decoder typedef SIMD mesg_type; #endif typedef DSP::Const Const; - static const int code_order = 11; + static const int sample_rate = 8000; + static const int code_order = 12; static const int mod_bits = 2; static const int code_len = 1 << code_order; - static const int symbol_len = (1280 * rate) / 8000; - static const int filter_len = (((21 * rate) / 8000) & ~3) | 1; + static const int symbol_len = 256; + static const int filter_len = 33; static const int guard_len = symbol_len / 8; - static const int max_bits = 1360 + 32; - static const int cons_cols = 256; - static const int cons_rows = 4; - static const int cons_total = cons_rows * cons_cols; - static const int code_off = - cons_cols / 2; - static const int mls0_len = 127; - static const int mls0_off = - mls0_len + 1; - static const int mls0_poly = 0b10001001; - static const int mls1_len = 255; - static const int mls1_off = - mls1_len / 2; - static const int mls1_poly = 0b100101011; + static const int data_bits = 2048; + static const int mesg_bits = data_bits + 32; + static const int subcarrier_count = 64; + static const int payload_symbols = 32; + static const int first_subcarrier = 16; + static const int cons_total = payload_symbols * subcarrier_count; static const int buffer_len = 6 * (symbol_len + guard_len); static const int search_pos = buffer_len - 4 * (symbol_len + guard_len); DSP::ReadPCM *pcm; @@ -191,29 +190,18 @@ struct Decoder DSP::BlockDC blockdc; DSP::Hilbert hilbert; DSP::BipBuffer input_hist; - DSP::TheilSenEstimator tse; - SchmidlCox correlator; - CODE::CRC crc0; - CODE::CRC crc1; - CODE::OrderedStatisticsDecoder<255, 71, 4> osddec; + SchmidlCox correlator; + CODE::CRC crc; CODE::PolarEncoder polarenc; CODE::PolarListDecoder polardec; - int8_t genmat[255*71]; - mesg_type mesg[max_bits], mess[code_len]; + mesg_type mesg[mesg_bits], mess[code_len]; code_type code[code_len]; - cmplx cons[cons_total], prev[cons_cols]; + cmplx cons[cons_total], prev[subcarrier_count]; cmplx fdom[symbol_len], tdom[symbol_len]; - value index[cons_cols], phase[cons_cols]; value cfo_rad, sfo_rad; const uint32_t *frozen_bits; int symbol_pos; - int oper_mode; - int crc_bits; - static int bin(int carrier) - { - return (carrier + symbol_len) % symbol_len; - } static int nrz(bool bit) { return 1 - 2 * bit; @@ -227,20 +215,20 @@ struct Decoder return 0; return cons; } - const cmplx *mls0_seq() + const cmplx *sync_seq() { - CODE::MLS seq0(mls0_poly); - for (int i = 0; i < symbol_len/2; ++i) + CODE::MLS seq(0b1100111); + for (int i = 0; i < symbol_len; ++i) fdom[i] = 0; - for (int i = 0; i < mls0_len; ++i) - fdom[(i+mls0_off/2+symbol_len/2)%(symbol_len/2)] = nrz(seq0()); + for (int i = first_subcarrier + 1; i < first_subcarrier + subcarrier_count; ++i) + fdom[i] = nrz(seq()); return fdom; } void systematic() { polarenc(mess, mesg, frozen_bits, code_order); int code_bits = 1 << code_order; - for (int i = 0, j = 0; i < code_bits && j < crc_bits; ++i) + for (int i = 0, j = 0; i < code_bits && j < mesg_bits; ++i) if (!((frozen_bits[i/32] >> (i%32)) & 1)) mesg[j++] = mess[i]; } @@ -258,25 +246,15 @@ struct Decoder } const cmplx *next_sample() { - cmplx tmp; - pcm->read(reinterpret_cast(&tmp), 1); - if (pcm->channels() == 1) - tmp = hilbert(blockdc(tmp.real())); - return input_hist(tmp); + value real; + pcm->read(&real, 1); + return input_hist(hilbert(blockdc(real))); } Decoder(uint8_t *out, int *len, DSP::ReadPCM *pcm, int skip_count) : - pcm(pcm), correlator(mls0_seq()), crc0(0xA8F4), crc1(0x8F6E37A0) + pcm(pcm), correlator(sync_seq()), crc(0x8F6E37A0) { - CODE::BoseChaudhuriHocquenghemGenerator<255, 71>::matrix(genmat, true, { - 0b100011101, 0b101110111, 0b111110011, 0b101101001, - 0b110111101, 0b111100111, 0b100101011, 0b111010111, - 0b000010011, 0b101100101, 0b110001011, 0b101100011, - 0b100011011, 0b100111111, 0b110001101, 0b100101101, - 0b101011111, 0b111111001, 0b111000011, 0b100111001, - 0b110101001, 0b000011111, 0b110000111, 0b110110001}); - - blockdc.samples(2*(symbol_len+guard_len)); - DSP::Phasor osc; + frozen_bits = frozen_4096_2080; + blockdc.samples(filter_len); const cmplx *buf; bool okay; do { @@ -290,70 +268,26 @@ struct Decoder symbol_pos = correlator.symbol_pos; cfo_rad = correlator.cfo_rad; std::cerr << "symbol pos: " << symbol_pos << std::endl; - std::cerr << "coarse cfo: " << cfo_rad * (rate / Const::TwoPi()) << " Hz " << std::endl; - - osc.omega(-cfo_rad); - for (int i = 0; i < symbol_len; ++i) - tdom[i] = buf[i+symbol_pos+(symbol_len+guard_len)] * osc(); - fwd(fdom, tdom); - CODE::MLS seq1(mls1_poly); - for (int i = 0; i < mls1_len; ++i) - fdom[bin(i+mls1_off)] *= nrz(seq1()); - int8_t soft[mls1_len]; - uint8_t data[(mls1_len+7)/8]; - for (int i = 0; i < mls1_len; ++i) - soft[i] = std::min(std::max( - std::nearbyint(127 * demod_or_erase( - fdom[bin(i+mls1_off)], fdom[bin(i-1+mls1_off)]).real()), - -128), 127); - bool unique = osddec(data, soft, genmat); - if (!unique) { - std::cerr << "OSD error." << std::endl; - continue; - } - uint64_t md = 0; - for (int i = 0; i < 55; ++i) - md |= (uint64_t)CODE::get_be_bit(data, i) << i; - uint16_t cs = 0; - for (int i = 0; i < 16; ++i) - cs |= (uint16_t)CODE::get_be_bit(data, i+55) << i; - crc0.reset(); - if (crc0(md<<9) != cs) { - std::cerr << "header CRC error." << std::endl; - continue; - } - oper_mode = md & 255; - if (oper_mode && (oper_mode < 14 || oper_mode > 16)) { - std::cerr << "operation mode " << oper_mode << " unsupported." << std::endl; - continue; - } - std::cerr << "oper mode: " << oper_mode << std::endl; - if ((md>>8) == 0 || (md>>8) >= 129961739795077L) { - std::cerr << "call sign unsupported." << std::endl; - continue; - } - char call_sign[10]; - base37_decoder(call_sign, md>>8, 9); - call_sign[9] = 0; - std::cerr << "call sign: " << call_sign << std::endl; + std::cerr << "coarse cfo: " << cfo_rad * (sample_rate / Const::TwoPi()) << " Hz " << std::endl; okay = true; } while (skip_count--); - *len = 0; - if (!okay || !oper_mode) + if (!okay) return; - for (int i = 0; i < symbol_pos+(symbol_len+guard_len); ++i) + DSP::Phasor osc; + osc.omega(-cfo_rad); + for (int i = 0; i < symbol_pos; ++i) buf = next_sample(); for (int i = 0; i < symbol_len; ++i) tdom[i] = buf[i] * osc(); for (int i = 0; i < guard_len; ++i) osc(); fwd(fdom, tdom); - for (int i = 0; i < cons_cols; ++i) - prev[i] = fdom[bin(i+code_off)]; + for (int i = 0; i < subcarrier_count; ++i) + prev[i] = fdom[first_subcarrier+i]; std::cerr << "demod "; - for (int j = 0; j < cons_rows; ++j) { + for (int j = 0; j < payload_symbols; ++j) { for (int i = 0; i < symbol_len+guard_len; ++i) buf = next_sample(); for (int i = 0; i < symbol_len; ++i) @@ -361,56 +295,30 @@ struct Decoder for (int i = 0; i < guard_len; ++i) osc(); fwd(fdom, tdom); - for (int i = 0; i < cons_cols; ++i) - cons[cons_cols*j+i] = demod_or_erase(fdom[bin(i+code_off)], prev[i]); - for (int i = 0; i < cons_cols; ++i) - prev[i] = fdom[bin(i+code_off)]; + for (int i = 0; i < subcarrier_count; ++i) + cons[subcarrier_count*j+i] = demod_or_erase(fdom[first_subcarrier+i], prev[i]); + for (int i = 0; i < subcarrier_count; ++i) + prev[i] = fdom[first_subcarrier+i]; std::cerr << "."; } std::cerr << " done" << std::endl; - if (1) { - value sum_slope = 0, sum_yint = 0; - for (int j = 0; j < cons_rows; ++j) { - for (int i = 0; i < cons_cols; ++i) { - code_type tmp[mod_bits]; - mod_hard(tmp, cons[cons_cols*j+i]); - index[i] = i + code_off; - phase[i] = arg(cons[cons_cols*j+i] * conj(mod_map(tmp))); - } - tse.compute(index, phase, cons_cols); - //std::cerr << "Theil-Sen slope = " << tse.slope() << std::endl; - //std::cerr << "Theil-Sen yint = " << tse.yint() << std::endl; - sum_slope += tse.slope(); - sum_yint += tse.yint(); - for (int i = 0; i < cons_cols; ++i) - cons[cons_cols*j+i] *= DSP::polar(1, -tse(i+code_off)); - } - value avg_slope = sum_slope / cons_rows; - value avg_yint = sum_yint / cons_rows; - //for (int i = 0; i < cons_cnt; ++i) - // cons[i] *= DSP::polar(1, -(avg_yint+avg_slope*((i%cons_cols)+code_off))); - sfo_rad -= avg_slope * symbol_len / value(symbol_len+guard_len); - cfo_rad += avg_yint / (symbol_len+guard_len); - std::cerr << "coarse sfo: " << 1000000 * sfo_rad / Const::TwoPi() << " ppm" << std::endl; - std::cerr << "finer cfo: " << cfo_rad * (rate / Const::TwoPi()) << " Hz " << std::endl; - } if (1) { std::cerr << "Es/N0 (dB):"; value sp = 0, np = 0; - for (int j = 0; j < cons_rows; ++j) { - for (int i = 0; i < cons_cols; ++i) { + for (int j = 0; j < payload_symbols; ++j) { + for (int i = 0; i < subcarrier_count; ++i) { code_type tmp[mod_bits]; - mod_hard(tmp, cons[cons_cols*j+i]); + mod_hard(tmp, cons[subcarrier_count*j+i]); cmplx hard = mod_map(tmp); - cmplx error = cons[cons_cols*j+i] - hard; + cmplx error = cons[subcarrier_count*j+i] - hard; sp += norm(hard); np += norm(error); } value precision = sp / np; value snr = DSP::decibel(precision); std::cerr << " " << snr; - for (int i = 0; i < cons_cols; ++i) - mod_soft(code+2*(cons_cols*j+i), cons[cons_cols*j+i], precision); + for (int i = 0; i < subcarrier_count; ++i) + mod_soft(code+mod_bits*(subcarrier_count*j+i), cons[subcarrier_count*j+i], precision); } std::cerr << std::endl; } else { @@ -418,25 +326,6 @@ struct Decoder for (int i = 0; i < cons_total; ++i) mod_soft(code+mod_bits*i, cons[i], precision); } - int data_bits = 0; - switch (oper_mode) { - case 14: - data_bits = 1360; - frozen_bits = frozen_2048_1392; - break; - case 15: - data_bits = 1024; - frozen_bits = frozen_2048_1056; - break; - case 16: - data_bits = 680; - frozen_bits = frozen_2048_712; - break; - default: - return; - } - *len = data_bits / 8; - crc_bits = data_bits + 32; CODE::PolarHelper::PATH metric[mesg_type::SIZE]; polardec(metric, mesg, code, frozen_bits, code_order); systematic(); @@ -446,19 +335,19 @@ struct Decoder std::sort(order, order+mesg_type::SIZE, [metric](int a, int b){ return metric[a] < metric[b]; }); int best = -1; for (int k = 0; k < mesg_type::SIZE; ++k) { - crc1.reset(); - for (int i = 0; i < crc_bits; ++i) - crc1(mesg[i].v[order[k]] < 0); - if (crc1() == 0) { + crc.reset(); + for (int i = 0; i < mesg_bits; ++i) + crc(mesg[i].v[order[k]] < 0); + if (crc() == 0) { best = order[k]; break; } } if (best < 0) { std::cerr << "payload decoding error." << std::endl; - *len = 0; return; } + *len = data_bits / 8; int flips = 0; for (int i = 0, j = 0; i < data_bits; ++i, ++j) { while ((frozen_bits[j / 32] >> (j % 32)) & 1) @@ -479,9 +368,6 @@ int main(int argc, char **argv) return 1; } - typedef float value; - typedef DSP::Complex cmplx; - const char *output_name = argv[1]; if (output_name[0] == '-' && output_name[1] == 0) output_name = "/dev/stdout"; @@ -489,10 +375,10 @@ int main(int argc, char **argv) if (input_name[0] == '-' && input_name[1] == 0) input_name = "/dev/stdin"; - DSP::ReadWAV input_file(input_name); + DSP::ReadWAV input_file(input_name); - if (input_file.channels() < 1 || input_file.channels() > 2) { - std::cerr << "Only real or analytic signal (one or two channels) supported." << std::endl; + if (input_file.channels() != 1) { + std::cerr << "Only real signal (one channel) supported." << std::endl; return 1; } @@ -500,28 +386,17 @@ int main(int argc, char **argv) if (argc > 3) skip_count = std::atoi(argv[3]); - const int data_max = 1360 / 8; + const int data_max = 2048 / 8; uint8_t *output_data = new uint8_t[data_max]; int data_len = 0; - switch (input_file.rate()) { - case 8000: - delete new Decoder(output_data, &data_len, &input_file, skip_count); - break; - case 16000: - delete new Decoder(output_data, &data_len, &input_file, skip_count); - break; - case 44100: - delete new Decoder(output_data, &data_len, &input_file, skip_count); - break; - case 48000: - delete new Decoder(output_data, &data_len, &input_file, skip_count); - break; - default: + if (input_file.rate() != 8000) { std::cerr << "Unsupported sample rate." << std::endl; return 1; } + delete new Decoder(output_data, &data_len, &input_file, skip_count); + std::ofstream output_file(output_name, std::ios::binary | std::ios::trunc); if (output_file.bad()) { std::cerr << "Couldn't open file \"" << output_name << "\" for writing." << std::endl; diff --git a/encode.cc b/encode.cc index 046be60..39dbdf4 100644 --- a/encode.cc +++ b/encode.cc @@ -1,7 +1,7 @@ /* OFDM modem encoder -Copyright 2021 Ahmet Inan +Copyright 2023 Ahmet Inan */ #include @@ -9,9 +9,9 @@ Copyright 2021 Ahmet Inan #include #include "xorshift.hh" #include "complex.hh" -#include "utils.hh" +#include "phasor.hh" #include "bitman.hh" -#include "decibel.hh" +#include "utils.hh" #include "fft.hh" #include "wav.hh" #include "pcm.hh" @@ -21,325 +21,131 @@ Copyright 2021 Ahmet Inan #include "polar_tables.hh" #include "polar_helper.hh" #include "polar_encoder.hh" -#include "bose_chaudhuri_hocquenghem_encoder.hh" -template struct Encoder { + typedef float value; + typedef DSP::Complex cmplx; typedef int8_t code_type; static const int mod_bits = 2; - static const int code_order = 11; + static const int code_order = 12; static const int code_len = 1 << code_order; - static const int symbol_len = (1280 * rate) / 8000; + static const int symbol_len = 256; static const int guard_len = symbol_len / 8; - static const int max_bits = 1360 + 32; - static const int cons_cols = 256; - static const int cons_rows = 4; - static const int mls0_len = 127; - static const int mls0_poly = 0b10001001; - static const int mls1_len = 255; - static const int mls1_poly = 0b100101011; - static const int mls2_poly = 0b100101010001; + static const int data_bits = 2048; + static const int mesg_bits = data_bits + 32; + static const int subcarrier_count = 64; + static const int payload_symbols = 32; + static const int first_subcarrier = 16; DSP::WritePCM *pcm; DSP::FastFourierTransform bwd; - DSP::FastFourierTransform<4*symbol_len, cmplx, -1> fwd4; - DSP::FastFourierTransform<4*symbol_len, cmplx, 1> bwd4; - CODE::CRC crc0; - CODE::CRC crc1; - CODE::BoseChaudhuriHocquenghemEncoder<255, 71> bchenc; + CODE::CRC crc; CODE::PolarSysEnc polarenc; - code_type code[code_len], mesg[max_bits]; - cmplx fdom[symbol_len], fdom4[4*symbol_len]; - cmplx tdom[symbol_len], tdom4[4*symbol_len]; - cmplx temp[symbol_len]; - cmplx guard[guard_len]; - cmplx papr_min, papr_max; - int code_off; - int mls0_off; - int mls1_off; + code_type code[code_len], mesg[mesg_bits]; + cmplx fdom[symbol_len], tdom[symbol_len], guard[guard_len]; - static int bin(int carrier) - { - return (carrier + symbol_len) % symbol_len; - } - static int bin4(int carrier) - { - return (carrier + 4*symbol_len) % (4*symbol_len); - } static int nrz(bool bit) { return 1 - 2 * bit; } - void improve_papr() + void symbol() { - for (int i = 0; i < 4*symbol_len; ++i) - fdom4[i] = 0; - for (int i = -symbol_len/2; i < symbol_len/2; ++i) - fdom4[bin4(i)] = fdom[bin(i)]; - bwd4(tdom4, fdom4); - for (int i = 0; i < 4*symbol_len; ++i) - tdom4[i] /= std::sqrt(value(4*symbol_len)); - for (int i = 0; i < 4*symbol_len; ++i) { - value amp = std::max(std::abs(tdom4[i].real()), std::abs(tdom4[i].imag())); - if (amp > value(1)) - tdom4[i] /= amp; - } - fwd4(fdom4, tdom4); - for (int i = -symbol_len/2; i < symbol_len/2; ++i) - if (norm(temp[bin(i)])) - temp[bin(i)] = fdom4[bin4(i)] / std::sqrt(value(4*symbol_len)); - else - temp[bin(i)] = 0; - } - void symbol(bool papr_reduction = true) - { - for (int i = 0; i < symbol_len; ++i) - temp[i] = fdom[i]; - if (papr_reduction) - improve_papr(); - bwd(tdom, temp); + bwd(tdom, fdom); for (int i = 0; i < symbol_len; ++i) tdom[i] /= std::sqrt(value(8*symbol_len)); for (int i = 0; i < guard_len; ++i) { value x = value(i) / value(guard_len - 1); - value ratio(0.5); - x = std::min(x, ratio) / ratio; x = value(0.5) * (value(1) - std::cos(DSP::Const::Pi() * x)); guard[i] = DSP::lerp(guard[i], tdom[i+symbol_len-guard_len], x); } - cmplx peak, mean; - for (int i = 0; i < symbol_len; ++i) { - cmplx power(tdom[i].real() * tdom[i].real(), tdom[i].imag() * tdom[i].imag()); - peak = cmplx(std::max(peak.real(), power.real()), std::max(peak.imag(), power.imag())); - mean += power; - } - if (mean.real() > 0 && mean.imag() > 0) { - cmplx papr(peak.real() / mean.real(), peak.imag() / mean.imag()); - papr *= value(symbol_len); - papr_min = cmplx(std::min(papr_min.real(), papr.real()), std::min(papr_min.imag(), papr.imag())); - papr_max = cmplx(std::max(papr_max.real(), papr.real()), std::max(papr_max.imag(), papr.imag())); - } pcm->write(reinterpret_cast(guard), guard_len, 2); pcm->write(reinterpret_cast(tdom), symbol_len, 2); for (int i = 0; i < guard_len; ++i) guard[i] = tdom[i]; } - void pilot_block() - { - CODE::MLS seq2(mls2_poly); - value code_fac = std::sqrt(value(symbol_len) / value(cons_cols)); - for (int i = 0; i < symbol_len; ++i) - fdom[i] = 0; - for (int i = code_off; i < code_off + cons_cols; ++i) - fdom[bin(i)] = code_fac * nrz(seq2()); - symbol(); - } void schmidl_cox() { - CODE::MLS seq0(mls0_poly); - value mls0_fac = std::sqrt(value(2 * symbol_len) / value(mls0_len)); + CODE::MLS seq(0b1100111); for (int i = 0; i < symbol_len; ++i) fdom[i] = 0; - fdom[bin(mls0_off-2)] = mls0_fac; - for (int i = 0; i < mls0_len; ++i) - fdom[bin(2*i+mls0_off)] = nrz(seq0()); - for (int i = 0; i < mls0_len; ++i) - fdom[bin(2*i+mls0_off)] *= fdom[bin(2*(i-1)+mls0_off)]; - symbol(false); - } - void meta_data(uint64_t md) - { - uint8_t data[9] = { 0 }, parity[23] = { 0 }; - for (int i = 0; i < 55; ++i) - CODE::set_be_bit(data, i, (md>>i)&1); - crc0.reset(); - uint16_t cs = crc0(md << 9); - for (int i = 0; i < 16; ++i) - CODE::set_be_bit(data, i+55, (cs>>i)&1); - bchenc(data, parity); - CODE::MLS seq4(mls1_poly); - value mls1_fac = std::sqrt(value(symbol_len) / value(mls1_len)); - for (int i = 0; i < symbol_len; ++i) - fdom[i] = 0; - fdom[bin(mls1_off-1)] = mls1_fac; - for (int i = 0; i < 71; ++i) - fdom[bin(i+mls1_off)] = nrz(CODE::get_be_bit(data, i)); - for (int i = 71; i < mls1_len; ++i) - fdom[bin(i+mls1_off)] = nrz(CODE::get_be_bit(parity, i-71)); - for (int i = 0; i < mls1_len; ++i) - fdom[bin(i+mls1_off)] *= fdom[bin(i-1+mls1_off)]; - for (int i = 0; i < mls1_len; ++i) - fdom[bin(i+mls1_off)] *= nrz(seq4()); +#if 0 + fdom[first_subcarrier+1] = std::sqrt(value(2 * symbol_len) / value(subcarrier_count)); + for (int i = first_subcarrier + 2; i < first_subcarrier + subcarrier_count + 1; ++i) + fdom[i] = fdom[i-1] * cmplx(nrz(seq())); symbol(); + for (int i = first_subcarrier; i < first_subcarrier + subcarrier_count; ++i) + fdom[i] = fdom[i+1]; + fdom[first_subcarrier+subcarrier_count] = 0; + symbol(); +#else + fdom[first_subcarrier] = std::sqrt(value(2 * symbol_len) / value(subcarrier_count)); + for (int i = first_subcarrier + 1; i < first_subcarrier + subcarrier_count; ++i) + fdom[i] = fdom[i-1] * cmplx(nrz(seq())); + symbol(); + symbol(); +#endif } cmplx mod_map(code_type *b) { return PhaseShiftKeying<4, cmplx, code_type>::map(b); } - Encoder(DSP::WritePCM *pcm, const uint8_t *inp, int freq_off, uint64_t call_sign, int oper_mode) : - pcm(pcm), crc0(0xA8F4), crc1(0x8F6E37A0), bchenc({ - 0b100011101, 0b101110111, 0b111110011, 0b101101001, - 0b110111101, 0b111100111, 0b100101011, 0b111010111, - 0b000010011, 0b101100101, 0b110001011, 0b101100011, - 0b100011011, 0b100111111, 0b110001101, 0b100101101, - 0b101011111, 0b111111001, 0b111000011, 0b100111001, - 0b110101001, 0b000011111, 0b110000111, 0b110110001}) + Encoder(DSP::WritePCM *pcm, const uint8_t *inp) : pcm(pcm), crc(0x8F6E37A0) { - int offset = (freq_off * symbol_len) / rate; - code_off = offset - cons_cols / 2; - mls0_off = offset - mls0_len + 1; - mls1_off = offset - mls1_len / 2; - papr_min = cmplx(1000, 1000), papr_max = cmplx(-1000, -1000); - pilot_block(); schmidl_cox(); - meta_data((call_sign << 8) | oper_mode); - if (oper_mode > 0) { - const uint32_t *frozen_bits = nullptr; - int data_bits = 0; - switch (oper_mode) { - case 14: - data_bits = 1360; - frozen_bits = frozen_2048_1392; - break; - case 15: - data_bits = 1024; - frozen_bits = frozen_2048_1056; - break; - case 16: - data_bits = 680; - frozen_bits = frozen_2048_712; - break; - default: - return; - } - for (int i = 0; i < data_bits; ++i) - mesg[i] = nrz(CODE::get_le_bit(inp, i)); - crc1.reset(); - for (int i = 0; i < data_bits / 8; ++i) - crc1(inp[i]); - for (int i = 0; i < 32; ++i) - mesg[i+data_bits] = nrz((crc1()>>i)&1); - polarenc(code, mesg, frozen_bits, code_order); - for (int j = 0; j < cons_rows; ++j) { - for (int i = 0; i < cons_cols; ++i) - fdom[bin(i+code_off)] *= - mod_map(code+mod_bits*(cons_cols*j+i)); - symbol(); - } + for (int i = 0; i < data_bits; ++i) + mesg[i] = nrz(CODE::get_le_bit(inp, i)); + crc.reset(); + for (int i = 0; i < data_bits / 8; ++i) + crc(inp[i]); + for (int i = 0; i < 32; ++i) + mesg[i+data_bits] = nrz((crc()>>i)&1); + polarenc(code, mesg, frozen_4096_2080, code_order); + for (int j = 0; j < payload_symbols; ++j) { + for (int i = 0; i < subcarrier_count; ++i) + fdom[first_subcarrier+i] *= + mod_map(code+mod_bits*(subcarrier_count*j+i)); + symbol(); } for (int i = 0; i < symbol_len; ++i) fdom[i] = 0; symbol(); - std::cerr << "real PAPR: " << DSP::decibel(papr_min.real()) << " .. " << DSP::decibel(papr_max.real()) << " dB" << std::endl; - if (pcm->channels() == 2) - std::cerr << "imag PAPR: " << DSP::decibel(papr_min.imag()) << " .. " << DSP::decibel(papr_max.imag()) << " dB" << std::endl; } }; -long long int base37_encoder(const char *str) -{ - long long int acc = 0; - for (char c = *str++; c; c = *str++) { - acc *= 37; - if (c >= '0' && c <= '9') - acc += c - '0' + 1; - else if (c >= 'a' && c <= 'z') - acc += c - 'a' + 11; - else if (c >= 'A' && c <= 'Z') - acc += c - 'A' + 11; - else if (c != ' ') - return -1; - } - return acc; -} - int main(int argc, char **argv) { - if (argc < 6 || argc > 8) { - std::cerr << "usage: " << argv[0] << " OUTPUT RATE BITS CHANNELS INPUT [OFFSET] [CALLSIGN]" << std::endl; + if (argc != 3) { + std::cerr << "usage: " << argv[0] << " OUTPUT INPUT" << std::endl; return 1; } const char *output_name = argv[1]; if (output_name[0] == '-' && output_name[1] == 0) output_name = "/dev/stdout"; - int output_rate = std::atoi(argv[2]); - int output_bits = std::atoi(argv[3]); - int output_chan = std::atoi(argv[4]); - const char *input_name = argv[5]; + const char *input_name = argv[2]; if (input_name[0] == '-' && input_name[1] == 0) input_name = "/dev/stdin"; - int freq_off = output_chan == 1 ? 2000 : 0; - if (argc >= 7) - freq_off = std::atoi(argv[6]); - - long long int call_sign = base37_encoder("ANONYMOUS"); - if (argc >= 8) - call_sign = base37_encoder(argv[7]); - - if (call_sign <= 0 || call_sign >= 129961739795077L) { - std::cerr << "Unsupported call sign." << std::endl; - return 1; - } - - const int band_width = 1600; - - if ((output_chan == 1 && freq_off < band_width / 2) || freq_off < band_width / 2 - output_rate / 2 || freq_off > output_rate / 2 - band_width / 2) { - std::cerr << "Unsupported frequency offset." << std::endl; - return 1; - } - - if (freq_off % 50) { - std::cerr << "Frequency offset must be divisible by 50." << std::endl; - return 1; - } - - typedef float value; - typedef DSP::Complex cmplx; - std::ifstream input_file(input_name, std::ios::binary); if (input_file.bad()) { std::cerr << "Couldn't open file \"" << input_name << "\" for reading." << std::endl; return 1; } - const int data_len = 1360 / 8; + const int data_len = 2048 / 8; uint8_t *input_data = new uint8_t[data_len]; for (int i = 0; i < data_len; ++i) input_data[i] = std::max(input_file.get(), 0); - int oper_mode = 0; - for (int i = 128; i < 170; ++i) - if (!oper_mode && input_data[i]) - oper_mode = 14; - for (int i = 85; i < 128; ++i) - if (!oper_mode && input_data[i]) - oper_mode = 15; - for (int i = 0; i < 85; ++i) - if (!oper_mode && input_data[i]) - oper_mode = 16; + CODE::Xorshift32 scrambler; for (int i = 0; i < data_len; ++i) input_data[i] ^= scrambler(); - DSP::WriteWAV output_file(output_name, output_rate, output_bits, output_chan); + const int output_rate = 8000; + const int output_bits = 16; + const int output_chan = 1; + DSP::WriteWAV output_file(output_name, output_rate, output_bits, output_chan); output_file.silence(output_rate); - switch (output_rate) { - case 8000: - delete new Encoder(&output_file, input_data, freq_off, call_sign, oper_mode); - break; - case 16000: - delete new Encoder(&output_file, input_data, freq_off, call_sign, oper_mode); - break; - case 44100: - delete new Encoder(&output_file, input_data, freq_off, call_sign, oper_mode); - break; - case 48000: - delete new Encoder(&output_file, input_data, freq_off, call_sign, oper_mode); - break; - default: - std::cerr << "Unsupported sample rate." << std::endl; - return 1; - } + delete new Encoder(&output_file, input_data); output_file.silence(output_rate); delete []input_data; diff --git a/freezer.cc b/freezer.cc index c60e62e..b398b18 100644 --- a/freezer.cc +++ b/freezer.cc @@ -33,8 +33,6 @@ void code(int N, int K) int main() { - code<11>(2048, 1360+32); - code<11>(2048, 1024+32); - code<11>(2048, 680+32); + code<12>(4096, 2048+32); return 0; } diff --git a/polar_tables.hh b/polar_tables.hh index 56c0b90..b477b47 100644 --- a/polar_tables.hh +++ b/polar_tables.hh @@ -1,3 +1 @@ -static const uint32_t frozen_2048_1392[64] = { 0xffffffff, 0xffffffff, 0xffffffff, 0xffffffff, 0xffffffff, 0xffffffff, 0x7fffffff, 0x11f7fff, 0xffffffff, 0x7fffffff, 0x17ffffff, 0x117177f, 0x177f7fff, 0x1037f, 0x1011f, 0x1, 0xffffffff, 0x177fffff, 0x77f7fff, 0x1011f, 0x1173fff, 0x10117, 0x10117, 0x0, 0x117177f, 0x17, 0x3, 0x0, 0x1, 0x0, 0x0, 0x0, 0x7fffffff, 0x11f7fff, 0x11717ff, 0x117, 0x17177f, 0x3, 0x1, 0x0, 0x1037f, 0x1, 0x1, 0x0, 0x1, 0x0, 0x0, 0x0, 0x1011f, 0x1, 0x1, 0x0, 0x1, 0x0, 0x0, 0x0, 0x1, 0x0, 0x0, 0x0, 0x0, 0x0, 0x0, 0x0, }; -static const uint32_t frozen_2048_1056[64] = { 0xffffffff, 0xffffffff, 0xffffffff, 0xffffffff, 0xffffffff, 0xffffffff, 0xffffffff, 0x7fffffff, 0xffffffff, 0xffffffff, 0xffffffff, 0x7fffffff, 0xffffffff, 0x177fffff, 0x177f7fff, 0x1017f, 0xffffffff, 0xffffffff, 0xffffffff, 0x177f7fff, 0x7fffffff, 0x13f7fff, 0x1171fff, 0x117, 0x3fffffff, 0x11717ff, 0x7177f, 0x1, 0x1017f, 0x1, 0x1, 0x0, 0xffffffff, 0x7fffffff, 0x7fffffff, 0x1171fff, 0x17ffffff, 0x7177f, 0x1037f, 0x1, 0x77f7fff, 0x1013f, 0x10117, 0x1, 0x10117, 0x0, 0x0, 0x0, 0x1173fff, 0x10117, 0x117, 0x0, 0x7, 0x0, 0x0, 0x0, 0x1, 0x0, 0x0, 0x0, 0x0, 0x0, 0x0, 0x0, }; -static const uint32_t frozen_2048_712[64] = { 0xffffffff, 0xffffffff, 0xffffffff, 0xffffffff, 0xffffffff, 0xffffffff, 0xffffffff, 0xffffffff, 0xffffffff, 0xffffffff, 0xffffffff, 0xffffffff, 0xffffffff, 0xffffffff, 0xffffffff, 0x177fffff, 0xffffffff, 0xffffffff, 0xffffffff, 0xffffffff, 0xffffffff, 0xffffffff, 0x7fffffff, 0x11f7fff, 0xffffffff, 0x7fffffff, 0x1fffffff, 0x17177f, 0x177fffff, 0x1037f, 0x1011f, 0x1, 0xffffffff, 0xffffffff, 0xffffffff, 0x7fffffff, 0xffffffff, 0x1fffffff, 0x177fffff, 0x1077f, 0xffffffff, 0x177f7fff, 0x13f7fff, 0x10117, 0x1171fff, 0x117, 0x7, 0x0, 0x7fffffff, 0x1173fff, 0x11717ff, 0x7, 0x3077f, 0x1, 0x1, 0x0, 0x1013f, 0x1, 0x1, 0x0, 0x1, 0x0, 0x0, 0x0, }; +static const uint32_t frozen_4096_2080[128] = { 0xffffffff, 0xffffffff, 0xffffffff, 0xffffffff, 0xffffffff, 0xffffffff, 0xffffffff, 0xffffffff, 0xffffffff, 0xffffffff, 0xffffffff, 0xffffffff, 0xffffffff, 0xffffffff, 0xffffffff, 0x7fffffff, 0xffffffff, 0xffffffff, 0xffffffff, 0xffffffff, 0xffffffff, 0xffffffff, 0xffffffff, 0x17ffffff, 0xffffffff, 0xffffffff, 0xffffffff, 0x77f7fff, 0x7fffffff, 0x1173fff, 0x11717ff, 0x117, 0xffffffff, 0xffffffff, 0xffffffff, 0xffffffff, 0xffffffff, 0xffffffff, 0x7fffffff, 0x11f7fff, 0xffffffff, 0x7fffffff, 0x1fffffff, 0x117177f, 0x177fffff, 0x1077f, 0x1013f, 0x1, 0xffffffff, 0x177fffff, 0x177f7fff, 0x1013f, 0x11f7fff, 0x10117, 0x10117, 0x1, 0x11717ff, 0x117, 0x7, 0x0, 0x1, 0x0, 0x0, 0x0, 0xffffffff, 0xffffffff, 0xffffffff, 0x7fffffff, 0xffffffff, 0x17ffffff, 0x177fffff, 0x1077f, 0xffffffff, 0x177f7fff, 0x11f7fff, 0x10117, 0x11717ff, 0x10117, 0x17, 0x0, 0x7fffffff, 0x1171fff, 0x117177f, 0x17, 0x3177f, 0x1, 0x1, 0x0, 0x1017f, 0x1, 0x1, 0x0, 0x1, 0x0, 0x0, 0x0, 0x177fffff, 0x7177f, 0x1037f, 0x1, 0x1011f, 0x1, 0x1, 0x0, 0x10117, 0x1, 0x1, 0x0, 0x0, 0x0, 0x0, 0x0, 0x10117, 0x0, 0x0, 0x0, 0x0, 0x0, 0x0, 0x0, 0x0, 0x0, 0x0, 0x0, 0x0, 0x0, 0x0, 0x0, }; From e57afffecf6e595f9039c8a1d66eb862cb984188 Mon Sep 17 00:00:00 2001 From: Ahmet Inan Date: Fri, 10 Feb 2023 09:15:52 +0100 Subject: [PATCH 02/22] need to align the rotation --- decode.cc | 5 +++-- encode.cc | 2 +- 2 files changed, 4 insertions(+), 3 deletions(-) diff --git a/decode.cc b/decode.cc index ebdec08..3d5bd49 100644 --- a/decode.cc +++ b/decode.cc @@ -73,7 +73,7 @@ public: SchmidlCox(const cmplx *sequence) : threshold(value(0.2*match_len), value(0.3*match_len)) { -#if 0 +#if 1 fsh.omega(-1, symbol_len); #endif for (int i = 0; i < symbol_len; ++i) @@ -84,7 +84,8 @@ public: } bool operator()(const cmplx *samples) { - cmplx P = cor(fsh() * samples[search_pos] * conj(samples[search_pos+symbol_len+guard_len])); + cmplx S = fsh(); + cmplx P = conj(S) * cor(S * samples[search_pos] * conj(samples[search_pos+symbol_len+guard_len])); value R = value(0.5) * pwr(norm(samples[search_pos]) + norm(samples[search_pos+symbol_len+guard_len])); value min_R = 0.0001 * symbol_len; R = std::max(R, min_R); diff --git a/encode.cc b/encode.cc index 39dbdf4..bef988e 100644 --- a/encode.cc +++ b/encode.cc @@ -68,7 +68,7 @@ struct Encoder CODE::MLS seq(0b1100111); for (int i = 0; i < symbol_len; ++i) fdom[i] = 0; -#if 0 +#if 1 fdom[first_subcarrier+1] = std::sqrt(value(2 * symbol_len) / value(subcarrier_count)); for (int i = first_subcarrier + 2; i < first_subcarrier + subcarrier_count + 1; ++i) fdom[i] = fdom[i-1] * cmplx(nrz(seq())); From 964c3f0ea3e3a710ae1eb6aed091ba8fac8a7859 Mon Sep 17 00:00:00 2001 From: Ahmet Inan Date: Fri, 10 Feb 2023 10:02:15 +0100 Subject: [PATCH 03/22] the distance needs to include the guard now --- decode.cc | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/decode.cc b/decode.cc index 3d5bd49..69dd89c 100644 --- a/decode.cc +++ b/decode.cc @@ -115,7 +115,7 @@ public: if (!process) return false; - frac_cfo = phase_max / value(symbol_len); + frac_cfo = phase_max / value(symbol_len+guard_len); DSP::Phasor osc; osc.omega(frac_cfo); From cca0199d6377de5f2487ecf538286e82795f409b Mon Sep 17 00:00:00 2001 From: Ahmet Inan Date: Fri, 10 Feb 2023 11:24:18 +0100 Subject: [PATCH 04/22] omit guard between synchronization symbols --- decode.cc | 8 ++++---- encode.cc | 9 +++++---- 2 files changed, 9 insertions(+), 8 deletions(-) diff --git a/decode.cc b/decode.cc index 69dd89c..470146e 100644 --- a/decode.cc +++ b/decode.cc @@ -85,8 +85,8 @@ public: bool operator()(const cmplx *samples) { cmplx S = fsh(); - cmplx P = conj(S) * cor(S * samples[search_pos] * conj(samples[search_pos+symbol_len+guard_len])); - value R = value(0.5) * pwr(norm(samples[search_pos]) + norm(samples[search_pos+symbol_len+guard_len])); + cmplx P = conj(S) * cor(S * samples[search_pos] * conj(samples[search_pos+symbol_len])); + value R = value(0.5) * pwr(norm(samples[search_pos]) + norm(samples[search_pos+symbol_len])); value min_R = 0.0001 * symbol_len; R = std::max(R, min_R); value timing = match(norm(P) / (R * R)); @@ -115,11 +115,11 @@ public: if (!process) return false; - frac_cfo = phase_max / value(symbol_len+guard_len); + frac_cfo = phase_max / value(symbol_len); DSP::Phasor osc; osc.omega(frac_cfo); - symbol_pos = search_pos - index_max + guard_len; + symbol_pos = search_pos - index_max; index_max = 0; timing_max = 0; for (int i = 0; i < symbol_len; ++i) diff --git a/encode.cc b/encode.cc index bef988e..f4fd0bd 100644 --- a/encode.cc +++ b/encode.cc @@ -48,7 +48,7 @@ struct Encoder { return 1 - 2 * bit; } - void symbol() + void symbol(bool output_guard = true) { bwd(tdom, fdom); for (int i = 0; i < symbol_len; ++i) @@ -58,7 +58,8 @@ struct Encoder x = value(0.5) * (value(1) - std::cos(DSP::Const::Pi() * x)); guard[i] = DSP::lerp(guard[i], tdom[i+symbol_len-guard_len], x); } - pcm->write(reinterpret_cast(guard), guard_len, 2); + if (output_guard) + pcm->write(reinterpret_cast(guard), guard_len, 2); pcm->write(reinterpret_cast(tdom), symbol_len, 2); for (int i = 0; i < guard_len; ++i) guard[i] = tdom[i]; @@ -76,13 +77,13 @@ struct Encoder for (int i = first_subcarrier; i < first_subcarrier + subcarrier_count; ++i) fdom[i] = fdom[i+1]; fdom[first_subcarrier+subcarrier_count] = 0; - symbol(); + symbol(false); #else fdom[first_subcarrier] = std::sqrt(value(2 * symbol_len) / value(subcarrier_count)); for (int i = first_subcarrier + 1; i < first_subcarrier + subcarrier_count; ++i) fdom[i] = fdom[i-1] * cmplx(nrz(seq())); symbol(); - symbol(); + symbol(false); #endif } cmplx mod_map(code_type *b) From b1d704513b3f988f0281cd76b7c968aac08c07a7 Mon Sep 17 00:00:00 2001 From: Ahmet Inan Date: Fri, 10 Feb 2023 11:28:59 +0100 Subject: [PATCH 05/22] frequency shifted sync symbol not worth the trouble --- decode.cc | 7 +------ encode.cc | 11 ----------- 2 files changed, 1 insertion(+), 17 deletions(-) diff --git a/decode.cc b/decode.cc index 470146e..0e62a86 100644 --- a/decode.cc +++ b/decode.cc @@ -45,7 +45,6 @@ struct SchmidlCox DSP::Delay delay; DSP::SchmittTrigger threshold; DSP::FallingEdgeTrigger falling; - DSP::Phasor fsh; cmplx tmp0[symbol_len], tmp1[symbol_len], tmp2[symbol_len]; cmplx seq[symbol_len], kern[symbol_len]; cmplx cmplx_shift = 0; @@ -73,9 +72,6 @@ public: SchmidlCox(const cmplx *sequence) : threshold(value(0.2*match_len), value(0.3*match_len)) { -#if 1 - fsh.omega(-1, symbol_len); -#endif for (int i = 0; i < symbol_len; ++i) seq[i] = sequence[i]; fwd(kern, sequence); @@ -84,8 +80,7 @@ public: } bool operator()(const cmplx *samples) { - cmplx S = fsh(); - cmplx P = conj(S) * cor(S * samples[search_pos] * conj(samples[search_pos+symbol_len])); + cmplx P = cor(samples[search_pos] * conj(samples[search_pos+symbol_len])); value R = value(0.5) * pwr(norm(samples[search_pos]) + norm(samples[search_pos+symbol_len])); value min_R = 0.0001 * symbol_len; R = std::max(R, min_R); diff --git a/encode.cc b/encode.cc index f4fd0bd..7b4d264 100644 --- a/encode.cc +++ b/encode.cc @@ -69,22 +69,11 @@ struct Encoder CODE::MLS seq(0b1100111); for (int i = 0; i < symbol_len; ++i) fdom[i] = 0; -#if 1 - fdom[first_subcarrier+1] = std::sqrt(value(2 * symbol_len) / value(subcarrier_count)); - for (int i = first_subcarrier + 2; i < first_subcarrier + subcarrier_count + 1; ++i) - fdom[i] = fdom[i-1] * cmplx(nrz(seq())); - symbol(); - for (int i = first_subcarrier; i < first_subcarrier + subcarrier_count; ++i) - fdom[i] = fdom[i+1]; - fdom[first_subcarrier+subcarrier_count] = 0; - symbol(false); -#else fdom[first_subcarrier] = std::sqrt(value(2 * symbol_len) / value(subcarrier_count)); for (int i = first_subcarrier + 1; i < first_subcarrier + subcarrier_count; ++i) fdom[i] = fdom[i-1] * cmplx(nrz(seq())); symbol(); symbol(false); -#endif } cmplx mod_map(code_type *b) { From a9d4249ea52862f7d4e6fe5ef0001c462d443ecc Mon Sep 17 00:00:00 2001 From: Ahmet Inan Date: Fri, 10 Feb 2023 11:43:04 +0100 Subject: [PATCH 06/22] removed redundant buffers --- decode.cc | 12 +++++------- 1 file changed, 5 insertions(+), 7 deletions(-) diff --git a/decode.cc b/decode.cc index 0e62a86..5451f70 100644 --- a/decode.cc +++ b/decode.cc @@ -45,8 +45,8 @@ struct SchmidlCox DSP::Delay delay; DSP::SchmittTrigger threshold; DSP::FallingEdgeTrigger falling; - cmplx tmp0[symbol_len], tmp1[symbol_len], tmp2[symbol_len]; - cmplx seq[symbol_len], kern[symbol_len]; + cmplx tmp0[symbol_len], tmp1[symbol_len]; + cmplx kern[symbol_len]; cmplx cmplx_shift = 0; value timing_max = 0; value phase_max = 0; @@ -72,8 +72,6 @@ public: SchmidlCox(const cmplx *sequence) : threshold(value(0.2*match_len), value(0.3*match_len)) { - for (int i = 0; i < symbol_len; ++i) - seq[i] = sequence[i]; fwd(kern, sequence); for (int i = 0; i < symbol_len; ++i) kern[i] = conj(kern[i]) / value(symbol_len); @@ -125,13 +123,13 @@ public: fwd(tmp0, tmp1); for (int i = 0; i < symbol_len; ++i) tmp0[i] *= kern[i]; - bwd(tmp2, tmp0); + bwd(tmp1, tmp0); int shift = 0; value peak = 0; value next = 0; for (int i = 0; i < symbol_len; ++i) { - value power = norm(tmp2[i]); + value power = norm(tmp1[i]); if (power > peak) { next = peak; peak = power; @@ -143,7 +141,7 @@ public: if (peak <= next * 4) return false; - int pos_err = std::nearbyint(arg(tmp2[shift]) * symbol_len / Const::TwoPi()); + int pos_err = std::nearbyint(arg(tmp1[shift]) * symbol_len / Const::TwoPi()); if (abs(pos_err) > guard_len / 2) return false; symbol_pos -= pos_err; From 19875626517b3d4347ba03d79f4649607fc5c0d2 Mon Sep 17 00:00:00 2001 From: Ahmet Inan Date: Fri, 10 Feb 2023 11:53:47 +0100 Subject: [PATCH 07/22] reduced size of buffering --- decode.cc | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/decode.cc b/decode.cc index 5451f70..43178cf 100644 --- a/decode.cc +++ b/decode.cc @@ -177,8 +177,8 @@ struct Decoder static const int payload_symbols = 32; static const int first_subcarrier = 16; static const int cons_total = payload_symbols * subcarrier_count; - static const int buffer_len = 6 * (symbol_len + guard_len); - static const int search_pos = buffer_len - 4 * (symbol_len + guard_len); + static const int buffer_len = 2 * symbol_len + guard_len; + static const int search_pos = symbol_len; DSP::ReadPCM *pcm; DSP::FastFourierTransform fwd; DSP::BlockDC blockdc; From 826992a2d023804e413487c34f879f0aaa7e15fd Mon Sep 17 00:00:00 2001 From: Ahmet Inan Date: Fri, 10 Feb 2023 12:08:01 +0100 Subject: [PATCH 08/22] use only 50% of the suffix for windowing --- encode.cc | 2 ++ 1 file changed, 2 insertions(+) diff --git a/encode.cc b/encode.cc index 7b4d264..b1a8ec9 100644 --- a/encode.cc +++ b/encode.cc @@ -55,6 +55,8 @@ struct Encoder tdom[i] /= std::sqrt(value(8*symbol_len)); for (int i = 0; i < guard_len; ++i) { value x = value(i) / value(guard_len - 1); + value ratio(0.5); + x = std::min(x, ratio) / ratio; x = value(0.5) * (value(1) - std::cos(DSP::Const::Pi() * x)); guard[i] = DSP::lerp(guard[i], tdom[i+symbol_len-guard_len], x); } From 2763a9c2e8762600813cb0e9752e3edb94f75b58 Mon Sep 17 00:00:00 2001 From: Ahmet Inan Date: Fri, 10 Feb 2023 12:32:38 +0100 Subject: [PATCH 09/22] renamed delay to align --- decode.cc | 6 ++---- 1 file changed, 2 insertions(+), 4 deletions(-) diff --git a/decode.cc b/decode.cc index 43178cf..40e2dd1 100644 --- a/decode.cc +++ b/decode.cc @@ -42,7 +42,7 @@ struct SchmidlCox DSP::SMA4 cor; DSP::SMA4 pwr; DSP::SMA4 match; - DSP::Delay delay; + DSP::Delay align; DSP::SchmittTrigger threshold; DSP::FallingEdgeTrigger falling; cmplx tmp0[symbol_len], tmp1[symbol_len]; @@ -83,13 +83,11 @@ public: value min_R = 0.0001 * symbol_len; R = std::max(R, min_R); value timing = match(norm(P) / (R * R)); - value phase = delay(arg(P)); + value phase = align(arg(P)); bool collect = threshold(timing); bool process = falling(collect); - //std::cout << timing << " " << phase << " " << process << std::endl; - if (!collect && !process) return false; From c305cc78626356c1073980784c3c74e1640a6309 Mon Sep 17 00:00:00 2001 From: Ahmet Inan Date: Sat, 11 Feb 2023 19:29:51 +0100 Subject: [PATCH 10/22] corrected simulation example --- README.md | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/README.md b/README.md index a02dd24..db6778e 100644 --- a/README.md +++ b/README.md @@ -46,7 +46,7 @@ Prerequisite: [disorders](https://github.com/aicodix/disorders) Encode ```uncoded.dat``` to [analytic](https://en.wikipedia.org/wiki/Analytic_signal) audio signal, add [multipath](https://en.wikipedia.org/wiki/Multipath_propagation), [CFO, SFO](https://en.wikipedia.org/wiki/Carrier_frequency_offset), [AWGN](https://en.wikipedia.org/wiki/Additive_white_Gaussian_noise), decode and compare: ``` -./encode - 8000 16 2 uncoded.dat 2000 | ../disorders/multipath - - ../disorders/multipath.txt 10 | ../disorders/cfo - - 234.567 | ../disorders/sfo - - 147 | ../disorders/awgn - - -30 | ./decode - - | diff -q -s uncoded.dat - +./encode - uncoded.dat | ../disorders/multipath - - ../disorders/multipath.txt 10 | ../disorders/cfo - - 234.567 | ../disorders/sfo - - 147 | ../disorders/awgn - - -30 | ./decode - - | diff -q -s uncoded.dat - ``` ### Reading From 13dd3481c77691151dafde0d9834825cbc9f8ff4 Mon Sep 17 00:00:00 2001 From: Ahmet Inan Date: Sat, 11 Feb 2023 21:18:14 +0100 Subject: [PATCH 11/22] added 8 bit operation mode symbol --- decode.cc | 44 +++++++++++++++++++++++++++++++++----------- encode.cc | 10 ++++++++++ 2 files changed, 43 insertions(+), 11 deletions(-) diff --git a/decode.cc b/decode.cc index 40e2dd1..ef901af 100644 --- a/decode.cc +++ b/decode.cc @@ -26,6 +26,7 @@ namespace DSP { using std::abs; using std::min; using std::cos; using std::sin; #include "mls.hh" #include "crc.hh" #include "psk.hh" +#include "hadamard_decoder.hh" #include "polar_tables.hh" #include "polar_helper.hh" #include "polar_encoder.hh" @@ -184,6 +185,7 @@ struct Decoder DSP::BipBuffer input_hist; SchmidlCox correlator; CODE::CRC crc; + CODE::HadamardDecoder<8> hadamard; CODE::PolarEncoder polarenc; CODE::PolarListDecoder polardec; mesg_type mesg[mesg_bits], mess[code_len]; @@ -247,6 +249,7 @@ struct Decoder { frozen_bits = frozen_4096_2080; blockdc.samples(filter_len); + DSP::Phasor osc; const cmplx *buf; bool okay; do { @@ -261,23 +264,42 @@ struct Decoder cfo_rad = correlator.cfo_rad; std::cerr << "symbol pos: " << symbol_pos << std::endl; std::cerr << "coarse cfo: " << cfo_rad * (sample_rate / Const::TwoPi()) << " Hz " << std::endl; + + osc.omega(-cfo_rad); + for (int i = 0; i < symbol_pos; ++i) + buf = next_sample(); + for (int i = 0; i < symbol_len; ++i) + tdom[i] = buf[i] * osc(); + for (int i = 0; i < guard_len; ++i) + osc(); + fwd(fdom, tdom); + for (int i = 0; i < subcarrier_count; ++i) + prev[i] = fdom[first_subcarrier+i]; + for (int i = 0; i < symbol_len+guard_len; ++i) + buf = next_sample(); + for (int i = 0; i < symbol_len; ++i) + tdom[i] = buf[i] * osc(); + for (int i = 0; i < guard_len; ++i) + osc(); + fwd(fdom, tdom); + for (int i = 0; i < subcarrier_count; ++i) + cons[i] = demod_or_erase(fdom[first_subcarrier+i], prev[i]); + for (int i = 0; i < subcarrier_count; ++i) + prev[i] = fdom[first_subcarrier+i]; + for (int i = 0; i < subcarrier_count; ++i) + mod_soft(code+mod_bits*i, cons[i], 8); + int oper_mode = hadamard(code); + if (oper_mode != 1) { + std::cerr << "operation mode " << oper_mode << " unsupported." << std::endl; + continue; + } + std::cerr << "oper mode: " << oper_mode << std::endl; okay = true; } while (skip_count--); if (!okay) return; - DSP::Phasor osc; - osc.omega(-cfo_rad); - for (int i = 0; i < symbol_pos; ++i) - buf = next_sample(); - for (int i = 0; i < symbol_len; ++i) - tdom[i] = buf[i] * osc(); - for (int i = 0; i < guard_len; ++i) - osc(); - fwd(fdom, tdom); - for (int i = 0; i < subcarrier_count; ++i) - prev[i] = fdom[first_subcarrier+i]; std::cerr << "demod "; for (int j = 0; j < payload_symbols; ++j) { for (int i = 0; i < symbol_len+guard_len; ++i) diff --git a/encode.cc b/encode.cc index b1a8ec9..0c76fc9 100644 --- a/encode.cc +++ b/encode.cc @@ -18,6 +18,7 @@ Copyright 2023 Ahmet Inan #include "mls.hh" #include "crc.hh" #include "psk.hh" +#include "hadamard_encoder.hh" #include "polar_tables.hh" #include "polar_helper.hh" #include "polar_encoder.hh" @@ -40,6 +41,7 @@ struct Encoder DSP::WritePCM *pcm; DSP::FastFourierTransform bwd; CODE::CRC crc; + CODE::HadamardEncoder<8> hadamard; CODE::PolarSysEnc polarenc; code_type code[code_len], mesg[mesg_bits]; cmplx fdom[symbol_len], tdom[symbol_len], guard[guard_len]; @@ -77,6 +79,13 @@ struct Encoder symbol(); symbol(false); } + void meta_data(int data) + { + hadamard(code, data); + for (int i = 0; i < subcarrier_count; ++i) + fdom[first_subcarrier+i] *= mod_map(code+mod_bits*i); + symbol(); + } cmplx mod_map(code_type *b) { return PhaseShiftKeying<4, cmplx, code_type>::map(b); @@ -84,6 +93,7 @@ struct Encoder Encoder(DSP::WritePCM *pcm, const uint8_t *inp) : pcm(pcm), crc(0x8F6E37A0) { schmidl_cox(); + meta_data(1); for (int i = 0; i < data_bits; ++i) mesg[i] = nrz(CODE::get_le_bit(inp, i)); crc.reset(); From ad1dc32e7fdfa703e2002840f1c4c819647d9c1f Mon Sep 17 00:00:00 2001 From: Ahmet Inan Date: Sun, 12 Feb 2023 11:10:33 +0100 Subject: [PATCH 12/22] lowered lowest allowed volume by an order of magnitude --- decode.cc | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/decode.cc b/decode.cc index ef901af..899fedb 100644 --- a/decode.cc +++ b/decode.cc @@ -81,7 +81,7 @@ public: { cmplx P = cor(samples[search_pos] * conj(samples[search_pos+symbol_len])); value R = value(0.5) * pwr(norm(samples[search_pos]) + norm(samples[search_pos+symbol_len])); - value min_R = 0.0001 * symbol_len; + value min_R = 0.00001 * symbol_len; R = std::max(R, min_R); value timing = match(norm(P) / (R * R)); value phase = align(arg(P)); From 0c2404519ddb5ec6ba47e1faec416783e8bb4aba Mon Sep 17 00:00:00 2001 From: Ahmet Inan Date: Sun, 12 Feb 2023 11:30:02 +0100 Subject: [PATCH 13/22] added half a second of leading noise --- encode.cc | 12 ++++++++++++ 1 file changed, 12 insertions(+) diff --git a/encode.cc b/encode.cc index 0c76fc9..8b792e4 100644 --- a/encode.cc +++ b/encode.cc @@ -68,6 +68,17 @@ struct Encoder for (int i = 0; i < guard_len; ++i) guard[i] = tdom[i]; } + void leading_noise() { + CODE::MLS seq(0b100101010001); + value factor = std::sqrt(value(symbol_len) / value(subcarrier_count)); + for (int i = 0; i < symbol_len; ++i) + fdom[i] = 0; + for (int j = 0; j < 14; ++j) { + for (int i = 0; i < subcarrier_count; ++i) + fdom[first_subcarrier+i] = factor * cmplx(nrz(seq()), nrz(seq())); + symbol(); + } + } void schmidl_cox() { CODE::MLS seq(0b1100111); @@ -92,6 +103,7 @@ struct Encoder } Encoder(DSP::WritePCM *pcm, const uint8_t *inp) : pcm(pcm), crc(0x8F6E37A0) { + leading_noise(); schmidl_cox(); meta_data(1); for (int i = 0; i < data_bits; ++i) From 0e6b41a9c557cf3f33c380c9cb2886fb0f27743a Mon Sep 17 00:00:00 2001 From: Ahmet Inan Date: Thu, 16 Feb 2023 10:29:51 +0100 Subject: [PATCH 14/22] added comment for debugging --- decode.cc | 3 +++ 1 file changed, 3 insertions(+) diff --git a/decode.cc b/decode.cc index 899fedb..70767d1 100644 --- a/decode.cc +++ b/decode.cc @@ -89,6 +89,9 @@ public: bool collect = threshold(timing); bool process = falling(collect); + // std::cout << timing << " " << process << " " << index_max << std::endl; + // plot "< arecord -r 8000 -c 1 -f S16_LE | ./decode /dev/null - 1 | tee /tmp/data.txt" u ($1/33) w l, "/tmp/data.txt" u ($0-$3-1+16):2 w i, 0.2, 0.3 + if (!collect && !process) return false; From 03da985dd08ea1ce37769e61f1eb3d72dd01b69c Mon Sep 17 00:00:00 2001 From: Ahmet Inan Date: Sun, 5 Mar 2023 08:32:13 +0100 Subject: [PATCH 15/22] scramble preamble --- decode.cc | 4 ++++ encode.cc | 4 ++++ 2 files changed, 8 insertions(+) diff --git a/decode.cc b/decode.cc index 70767d1..a98ffbd 100644 --- a/decode.cc +++ b/decode.cc @@ -170,6 +170,7 @@ struct Decoder static const int code_order = 12; static const int mod_bits = 2; static const int code_len = 1 << code_order; + static const int meta_len = 128; static const int symbol_len = 256; static const int filter_len = 33; static const int guard_len = symbol_len / 8; @@ -291,6 +292,9 @@ struct Decoder prev[i] = fdom[first_subcarrier+i]; for (int i = 0; i < subcarrier_count; ++i) mod_soft(code+mod_bits*i, cons[i], 8); + CODE::MLS seq(0b10000011); + for (int i = 0; i < meta_len; ++i) + code[i] *= nrz(seq()); int oper_mode = hadamard(code); if (oper_mode != 1) { std::cerr << "operation mode " << oper_mode << " unsupported." << std::endl; diff --git a/encode.cc b/encode.cc index 8b792e4..1c6cff2 100644 --- a/encode.cc +++ b/encode.cc @@ -31,6 +31,7 @@ struct Encoder static const int mod_bits = 2; static const int code_order = 12; static const int code_len = 1 << code_order; + static const int meta_len = 128; static const int symbol_len = 256; static const int guard_len = symbol_len / 8; static const int data_bits = 2048; @@ -93,6 +94,9 @@ struct Encoder void meta_data(int data) { hadamard(code, data); + CODE::MLS seq(0b10000011); + for (int i = 0; i < meta_len; ++i) + code[i] *= nrz(seq()); for (int i = 0; i < subcarrier_count; ++i) fdom[first_subcarrier+i] *= mod_map(code+mod_bits*i); symbol(); From c844a3a1817fcb71c5919fa6f6afc4c0004aaecc Mon Sep 17 00:00:00 2001 From: Ahmet Inan Date: Wed, 8 Mar 2023 08:24:47 +0100 Subject: [PATCH 16/22] center subcarriers at DC for decoding --- decode.cc | 20 ++++++++++++-------- 1 file changed, 12 insertions(+), 8 deletions(-) diff --git a/decode.cc b/decode.cc index a98ffbd..1b1efbf 100644 --- a/decode.cc +++ b/decode.cc @@ -178,7 +178,7 @@ struct Decoder static const int mesg_bits = data_bits + 32; static const int subcarrier_count = 64; static const int payload_symbols = 32; - static const int first_subcarrier = 16; + static const int first_subcarrier = -subcarrier_count / 2; static const int cons_total = payload_symbols * subcarrier_count; static const int buffer_len = 2 * symbol_len + guard_len; static const int search_pos = symbol_len; @@ -200,6 +200,10 @@ struct Decoder const uint32_t *frozen_bits; int symbol_pos; + static int bin(int carrier) + { + return (carrier + first_subcarrier + symbol_len) % symbol_len; + } static int nrz(bool bit) { return 1 - 2 * bit; @@ -218,8 +222,8 @@ struct Decoder CODE::MLS seq(0b1100111); for (int i = 0; i < symbol_len; ++i) fdom[i] = 0; - for (int i = first_subcarrier + 1; i < first_subcarrier + subcarrier_count; ++i) - fdom[i] = nrz(seq()); + for (int i = 1; i < subcarrier_count; ++i) + fdom[bin(i)] = nrz(seq()); return fdom; } void systematic() @@ -278,7 +282,7 @@ struct Decoder osc(); fwd(fdom, tdom); for (int i = 0; i < subcarrier_count; ++i) - prev[i] = fdom[first_subcarrier+i]; + prev[i] = fdom[bin(i)]; for (int i = 0; i < symbol_len+guard_len; ++i) buf = next_sample(); for (int i = 0; i < symbol_len; ++i) @@ -287,9 +291,9 @@ struct Decoder osc(); fwd(fdom, tdom); for (int i = 0; i < subcarrier_count; ++i) - cons[i] = demod_or_erase(fdom[first_subcarrier+i], prev[i]); + cons[i] = demod_or_erase(fdom[bin(i)], prev[i]); for (int i = 0; i < subcarrier_count; ++i) - prev[i] = fdom[first_subcarrier+i]; + prev[i] = fdom[bin(i)]; for (int i = 0; i < subcarrier_count; ++i) mod_soft(code+mod_bits*i, cons[i], 8); CODE::MLS seq(0b10000011); @@ -317,9 +321,9 @@ struct Decoder osc(); fwd(fdom, tdom); for (int i = 0; i < subcarrier_count; ++i) - cons[subcarrier_count*j+i] = demod_or_erase(fdom[first_subcarrier+i], prev[i]); + cons[subcarrier_count*j+i] = demod_or_erase(fdom[bin(i)], prev[i]); for (int i = 0; i < subcarrier_count; ++i) - prev[i] = fdom[first_subcarrier+i]; + prev[i] = fdom[bin(i)]; std::cerr << "."; } std::cerr << " done" << std::endl; From 752b155b7ac7918b21e19c5af48b6d753206ea6c Mon Sep 17 00:00:00 2001 From: Ahmet Inan Date: Thu, 11 May 2023 09:05:48 +0200 Subject: [PATCH 17/22] corrected test --- Makefile | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Makefile b/Makefile index afe04ed..244f271 100644 --- a/Makefile +++ b/Makefile @@ -11,7 +11,7 @@ CXX = clang++ -stdlib=libc++ -march=native all: encode decode test: encode decode - $(QEMU) ./encode encoded.wav 8000 8 1 /dev/urandom + $(QEMU) ./encode encoded.wav /dev/urandom $(QEMU) ./decode /dev/null encoded.wav encode: encode.cc From 9c2cebe202f3d098cd567345430e86eb72393bb1 Mon Sep 17 00:00:00 2001 From: Ahmet Inan Date: Thu, 11 May 2023 09:06:40 +0200 Subject: [PATCH 18/22] added aarch64 example --- Makefile | 3 +++ 1 file changed, 3 insertions(+) diff --git a/Makefile b/Makefile index 244f271..7b1379a 100644 --- a/Makefile +++ b/Makefile @@ -6,6 +6,9 @@ CXX = clang++ -stdlib=libc++ -march=native #CXX = armv7a-hardfloat-linux-gnueabi-g++ -static -mfpu=neon -march=armv7-a #QEMU = qemu-arm +#CXX = aarch64-unknown-linux-gnu-g++ -static -march=armv8-a+crc+simd -mtune=cortex-a72 +#QEMU = qemu-aarch64 + .PHONY: all all: encode decode From 5afd9fa8d6b5dccbaebaba61ab01e20ad8d52de4 Mon Sep 17 00:00:00 2001 From: Ahmet Inan Date: Thu, 11 May 2023 09:07:02 +0200 Subject: [PATCH 19/22] shuffle bits to improve decoding with burst errors --- decode.cc | 3 +++ encode.cc | 3 +++ 2 files changed, 6 insertions(+) diff --git a/decode.cc b/decode.cc index 1b1efbf..c0e26f6 100644 --- a/decode.cc +++ b/decode.cc @@ -13,6 +13,7 @@ namespace DSP { using std::abs; using std::min; using std::cos; using std::sin; #include "xorshift.hh" #include "trigger.hh" #include "complex.hh" +#include "permute.hh" #include "decibel.hh" #include "blockdc.hh" #include "hilbert.hh" @@ -192,6 +193,7 @@ struct Decoder CODE::HadamardDecoder<8> hadamard; CODE::PolarEncoder polarenc; CODE::PolarListDecoder polardec; + CODE::ReverseFisherYatesShuffle shuffle; mesg_type mesg[mesg_bits], mess[code_len]; code_type code[code_len]; cmplx cons[cons_total], prev[subcarrier_count]; @@ -352,6 +354,7 @@ struct Decoder mod_soft(code+mod_bits*i, cons[i], precision); } CODE::PolarHelper::PATH metric[mesg_type::SIZE]; + shuffle(code); polardec(metric, mesg, code, frozen_bits, code_order); systematic(); int order[mesg_type::SIZE]; diff --git a/encode.cc b/encode.cc index 1c6cff2..04d0ffb 100644 --- a/encode.cc +++ b/encode.cc @@ -9,6 +9,7 @@ Copyright 2023 Ahmet Inan #include #include "xorshift.hh" #include "complex.hh" +#include "permute.hh" #include "phasor.hh" #include "bitman.hh" #include "utils.hh" @@ -44,6 +45,7 @@ struct Encoder CODE::CRC crc; CODE::HadamardEncoder<8> hadamard; CODE::PolarSysEnc polarenc; + CODE::FisherYatesShuffle shuffle; code_type code[code_len], mesg[mesg_bits]; cmplx fdom[symbol_len], tdom[symbol_len], guard[guard_len]; @@ -118,6 +120,7 @@ struct Encoder for (int i = 0; i < 32; ++i) mesg[i+data_bits] = nrz((crc()>>i)&1); polarenc(code, mesg, frozen_4096_2080, code_order); + shuffle(code); for (int j = 0; j < payload_symbols; ++j) { for (int i = 0; i < subcarrier_count; ++i) fdom[first_subcarrier+i] *= From 0698e22d57c2e109d66ff1b7753db4be59d8dfac Mon Sep 17 00:00:00 2001 From: Ahmet Inan Date: Thu, 11 May 2023 13:06:45 +0200 Subject: [PATCH 20/22] use a horizontally modulated simplex code for the preamble --- decode.cc | 29 +++++++++++------------------ encode.cc | 15 +++++++-------- 2 files changed, 18 insertions(+), 26 deletions(-) diff --git a/decode.cc b/decode.cc index c0e26f6..b839692 100644 --- a/decode.cc +++ b/decode.cc @@ -27,7 +27,7 @@ namespace DSP { using std::abs; using std::min; using std::cos; using std::sin; #include "mls.hh" #include "crc.hh" #include "psk.hh" -#include "hadamard_decoder.hh" +#include "simplex_decoder.hh" #include "polar_tables.hh" #include "polar_helper.hh" #include "polar_encoder.hh" @@ -171,7 +171,7 @@ struct Decoder static const int code_order = 12; static const int mod_bits = 2; static const int code_len = 1 << code_order; - static const int meta_len = 128; + static const int meta_len = 63; static const int symbol_len = 256; static const int filter_len = 33; static const int guard_len = symbol_len / 8; @@ -190,7 +190,7 @@ struct Decoder DSP::BipBuffer input_hist; SchmidlCox correlator; CODE::CRC crc; - CODE::HadamardDecoder<8> hadamard; + CODE::SimplexDecoder<6> simplex; CODE::PolarEncoder polarenc; CODE::PolarListDecoder polardec; CODE::ReverseFisherYatesShuffle shuffle; @@ -275,16 +275,9 @@ struct Decoder std::cerr << "symbol pos: " << symbol_pos << std::endl; std::cerr << "coarse cfo: " << cfo_rad * (sample_rate / Const::TwoPi()) << " Hz " << std::endl; - osc.omega(-cfo_rad); for (int i = 0; i < symbol_pos; ++i) - buf = next_sample(); - for (int i = 0; i < symbol_len; ++i) - tdom[i] = buf[i] * osc(); - for (int i = 0; i < guard_len; ++i) - osc(); - fwd(fdom, tdom); - for (int i = 0; i < subcarrier_count; ++i) - prev[i] = fdom[bin(i)]; + next_sample(); + osc.omega(-cfo_rad); for (int i = 0; i < symbol_len+guard_len; ++i) buf = next_sample(); for (int i = 0; i < symbol_len; ++i) @@ -292,16 +285,16 @@ struct Decoder for (int i = 0; i < guard_len; ++i) osc(); fwd(fdom, tdom); - for (int i = 0; i < subcarrier_count; ++i) - cons[i] = demod_or_erase(fdom[bin(i)], prev[i]); + for (int i = 0; i < meta_len; ++i) + cons[i] = demod_or_erase(fdom[bin(i+1)], fdom[bin(i)]); for (int i = 0; i < subcarrier_count; ++i) prev[i] = fdom[bin(i)]; - for (int i = 0; i < subcarrier_count; ++i) - mod_soft(code+mod_bits*i, cons[i], 8); - CODE::MLS seq(0b10000011); + for (int i = 0; i < meta_len; ++i) + code[i] = DSP::clamp(16 * cons[i].real(), -127, 127); + CODE::MLS seq(0b1000011); for (int i = 0; i < meta_len; ++i) code[i] *= nrz(seq()); - int oper_mode = hadamard(code); + int oper_mode = simplex(code); if (oper_mode != 1) { std::cerr << "operation mode " << oper_mode << " unsupported." << std::endl; continue; diff --git a/encode.cc b/encode.cc index 04d0ffb..6a9e9c3 100644 --- a/encode.cc +++ b/encode.cc @@ -19,7 +19,7 @@ Copyright 2023 Ahmet Inan #include "mls.hh" #include "crc.hh" #include "psk.hh" -#include "hadamard_encoder.hh" +#include "simplex_encoder.hh" #include "polar_tables.hh" #include "polar_helper.hh" #include "polar_encoder.hh" @@ -32,7 +32,7 @@ struct Encoder static const int mod_bits = 2; static const int code_order = 12; static const int code_len = 1 << code_order; - static const int meta_len = 128; + static const int meta_len = 63; static const int symbol_len = 256; static const int guard_len = symbol_len / 8; static const int data_bits = 2048; @@ -43,7 +43,7 @@ struct Encoder DSP::WritePCM *pcm; DSP::FastFourierTransform bwd; CODE::CRC crc; - CODE::HadamardEncoder<8> hadamard; + CODE::SimplexEncoder<6> simplex; CODE::PolarSysEnc polarenc; CODE::FisherYatesShuffle shuffle; code_type code[code_len], mesg[mesg_bits]; @@ -95,12 +95,11 @@ struct Encoder } void meta_data(int data) { - hadamard(code, data); - CODE::MLS seq(0b10000011); + simplex(code, data); + CODE::MLS seq(0b1000011); + fdom[first_subcarrier] = std::sqrt(value(symbol_len) / value(subcarrier_count)); for (int i = 0; i < meta_len; ++i) - code[i] *= nrz(seq()); - for (int i = 0; i < subcarrier_count; ++i) - fdom[first_subcarrier+i] *= mod_map(code+mod_bits*i); + fdom[first_subcarrier+1+i] = fdom[first_subcarrier+i] * cmplx(code[i] * nrz(seq())); symbol(); } cmplx mod_map(code_type *b) From 532394d3596892e2e3c3bc465350f69a032ec257 Mon Sep 17 00:00:00 2001 From: Ahmet Inan Date: Fri, 12 May 2023 09:29:01 +0200 Subject: [PATCH 21/22] keep limits symmetric so it doesn't interfere with scrambling --- psk.hh | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/psk.hh b/psk.hh index a07e85f..f20e3bc 100644 --- a/psk.hh +++ b/psk.hh @@ -26,7 +26,7 @@ struct PhaseShiftKeying<2, TYPE, CODE> if (std::is_integral::value) value = std::nearbyint(value); if (std::is_same::value) - value = std::min(std::max(value, -128), 127); + value = std::min(std::max(value, -127), 127); return value; } @@ -65,7 +65,7 @@ struct PhaseShiftKeying<4, TYPE, CODE> if (std::is_integral::value) value = std::nearbyint(value); if (std::is_same::value) - value = std::min(std::max(value, -128), 127); + value = std::min(std::max(value, -127), 127); return value; } @@ -111,7 +111,7 @@ struct PhaseShiftKeying<8, TYPE, CODE> if (std::is_integral::value) value = std::nearbyint(value); if (std::is_same::value) - value = std::min(std::max(value, -128), 127); + value = std::min(std::max(value, -127), 127); return value; } From a9e221adf1d63144fc04e82d704acc92436bee02 Mon Sep 17 00:00:00 2001 From: Ahmet Inan Date: Fri, 12 May 2023 09:42:19 +0200 Subject: [PATCH 22/22] added bpsk and qpsk types to reduce code --- decode.cc | 27 ++++++++------------------- encode.cc | 8 ++------ 2 files changed, 10 insertions(+), 25 deletions(-) diff --git a/decode.cc b/decode.cc index b839692..4b19b9a 100644 --- a/decode.cc +++ b/decode.cc @@ -167,9 +167,10 @@ struct Decoder typedef SIMD mesg_type; #endif typedef DSP::Const Const; + typedef PhaseShiftKeying<2, cmplx, code_type> bpsk; + typedef PhaseShiftKeying<4, cmplx, code_type> qpsk; static const int sample_rate = 8000; static const int code_order = 12; - static const int mod_bits = 2; static const int code_len = 1 << code_order; static const int meta_len = 63; static const int symbol_len = 256; @@ -236,18 +237,6 @@ struct Decoder if (!((frozen_bits[i/32] >> (i%32)) & 1)) mesg[j++] = mess[i]; } - cmplx mod_map(code_type *b) - { - return PhaseShiftKeying<4, cmplx, code_type>::map(b); - } - void mod_hard(code_type *b, cmplx c) - { - PhaseShiftKeying<4, cmplx, code_type>::hard(b, c); - } - void mod_soft(code_type *b, cmplx c, value precision) - { - PhaseShiftKeying<4, cmplx, code_type>::soft(b, c, precision); - } const cmplx *next_sample() { value real; @@ -290,7 +279,7 @@ struct Decoder for (int i = 0; i < subcarrier_count; ++i) prev[i] = fdom[bin(i)]; for (int i = 0; i < meta_len; ++i) - code[i] = DSP::clamp(16 * cons[i].real(), -127, 127); + bpsk::soft(code+i, cons[i], 8); CODE::MLS seq(0b1000011); for (int i = 0; i < meta_len; ++i) code[i] *= nrz(seq()); @@ -327,9 +316,9 @@ struct Decoder value sp = 0, np = 0; for (int j = 0; j < payload_symbols; ++j) { for (int i = 0; i < subcarrier_count; ++i) { - code_type tmp[mod_bits]; - mod_hard(tmp, cons[subcarrier_count*j+i]); - cmplx hard = mod_map(tmp); + code_type tmp[2]; + qpsk::hard(tmp, cons[subcarrier_count*j+i]); + cmplx hard = qpsk::map(tmp); cmplx error = cons[subcarrier_count*j+i] - hard; sp += norm(hard); np += norm(error); @@ -338,13 +327,13 @@ struct Decoder value snr = DSP::decibel(precision); std::cerr << " " << snr; for (int i = 0; i < subcarrier_count; ++i) - mod_soft(code+mod_bits*(subcarrier_count*j+i), cons[subcarrier_count*j+i], precision); + qpsk::soft(code+2*(subcarrier_count*j+i), cons[subcarrier_count*j+i], precision); } std::cerr << std::endl; } else { value precision = 8; for (int i = 0; i < cons_total; ++i) - mod_soft(code+mod_bits*i, cons[i], precision); + qpsk::soft(code+2*i, cons[i], precision); } CODE::PolarHelper::PATH metric[mesg_type::SIZE]; shuffle(code); diff --git a/encode.cc b/encode.cc index 6a9e9c3..fed1426 100644 --- a/encode.cc +++ b/encode.cc @@ -29,7 +29,7 @@ struct Encoder typedef float value; typedef DSP::Complex cmplx; typedef int8_t code_type; - static const int mod_bits = 2; + typedef PhaseShiftKeying<4, cmplx, code_type> qpsk; static const int code_order = 12; static const int code_len = 1 << code_order; static const int meta_len = 63; @@ -102,10 +102,6 @@ struct Encoder fdom[first_subcarrier+1+i] = fdom[first_subcarrier+i] * cmplx(code[i] * nrz(seq())); symbol(); } - cmplx mod_map(code_type *b) - { - return PhaseShiftKeying<4, cmplx, code_type>::map(b); - } Encoder(DSP::WritePCM *pcm, const uint8_t *inp) : pcm(pcm), crc(0x8F6E37A0) { leading_noise(); @@ -123,7 +119,7 @@ struct Encoder for (int j = 0; j < payload_symbols; ++j) { for (int i = 0; i < subcarrier_count; ++i) fdom[first_subcarrier+i] *= - mod_map(code+mod_bits*(subcarrier_count*j+i)); + qpsk::map(code+2*(subcarrier_count*j+i)); symbol(); } for (int i = 0; i < symbol_len; ++i)