diff --git a/Makefile b/Makefile index afe04ed..67da107 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 8000 8 1 2000 6 ANONYMOUS /dev/urandom $(QEMU) ./decode /dev/null encoded.wav encode: encode.cc diff --git a/README.md b/README.md index c513c31..094e587 100644 --- a/README.md +++ b/README.md @@ -12,7 +12,7 @@ dd if=/dev/urandom of=uncoded.dat bs=1 count=5380 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: ``` -./encode encoded.wav 8000 16 1 uncoded.dat +./encode encoded.wav 8000 16 1 2000 6 CALLSIGN uncoded.dat ``` Start recording to ```recorded.wav``` audio file and stop after 20 seconds: @@ -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 -s uncoded.dat - +./encode - 8000 16 2 2000 6 CALLSIGN uncoded.dat | ../disorders/multipath - - ../disorders/multipath.txt 10 | ../disorders/cfo - - 234.567 | ../disorders/sfo - - 147 | ../disorders/awgn - - -30 | ./decode - - | diff -s uncoded.dat - ``` ### Reading diff --git a/decode.cc b/decode.cc index c69a97f..cce016b 100644 --- a/decode.cc +++ b/decode.cc @@ -7,6 +7,7 @@ Copyright 2021 Ahmet Inan #include #include #include +#include #include namespace DSP { using std::abs; using std::min; using std::cos; using std::sin; } #include "bip_buffer.hh" @@ -501,41 +502,40 @@ struct Decoder std::cerr << "coarse sfo: " << 1000000 * sfo_rad / Const::TwoPi() << " ppm" << std::endl; std::cerr << "finer cfo: " << cfo_rad * (rate / Const::TwoPi()) << " Hz " << std::endl; } - value precision = 16; if (1) { + std::cerr << "Es/N0 (dB):"; value sp = 0, np = 0; - for (int i = 0; i < cons_cnt; ++i) { - code_type tmp[mod_max]; - mod_hard(tmp, cons[i]); - cmplx hard = mod_map(tmp); - cmplx error = cons[i] - hard; - sp += norm(hard); - np += norm(error); + for (int j = 0; j < cons_rows; ++j) { + for (int i = 0; i < cons_cols; ++i) { + code_type tmp[mod_max]; + mod_hard(tmp, cons[cons_cols*j+i]); + cmplx hard = mod_map(tmp); + cmplx error = cons[cons_cols*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+mod_bits*(cons_cols*j+i), cons[cons_cols*j+i], precision); } - value snr = DSP::decibel(sp / np); - std::cerr << "init Es/N0: " << snr << " dB" << std::endl; - // $LLR=log(\frac{p(x=+1|y)}{p(x=-1|y)})$ - // $p(x|\mu,\sigma)=\frac{1}{\sqrt{2\pi}\sigma}}e^{-\frac{(x-\mu)^2}{2\sigma^2}}$ - value sigma = std::sqrt(np / (2 * sp)); - precision = 1 / (sigma * sigma); + std::cerr << std::endl; + } else { + value precision = 8; + for (int i = 0; i < cons_cnt; ++i) + mod_soft(code+mod_bits*i, cons[i], precision); } - for (int i = 0; i < cons_cnt; ++i) - mod_soft(code+mod_bits*i, cons[i], precision); lengthen(); - CODE::PolarHelper::PATH metric[mesg_type::SIZE]; - polardec(metric, mesg, code, frozen_bits, code_order); + polardec(nullptr, mesg, code, frozen_bits, code_order); systematic(); - int order[mesg_type::SIZE]; - for (int k = 0; k < mesg_type::SIZE; ++k) - order[k] = k; - 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); + crc1(mesg[i].v[k] < 0); if (crc1() == 0) { - best = order[k]; + best = k; break; } } @@ -543,8 +543,16 @@ struct Decoder std::cerr << "payload decoding error." << std::endl; return; } - for (int i = 0; i < data_bits; ++i) - CODE::set_le_bit(out, i, mesg[i].v[best] < 0); + int flips = 0; + for (int i = 0, j = 0; i < data_bits; ++i, ++j) { + while ((frozen_bits[j / 32] >> (j % 32)) & 1) + ++j; + bool received = code[j] < 0; + bool decoded = mesg[i].v[best] < 0; + flips += received != decoded; + CODE::set_le_bit(out, i, decoded); + } + std::cerr << "bit flips: " << flips << std::endl; } }; diff --git a/encode.cc b/encode.cc index 5ed87fa..f09bd7a 100644 --- a/encode.cc +++ b/encode.cc @@ -6,6 +6,7 @@ Copyright 2021 Ahmet Inan #include #include +#include #include #include "xorshift.hh" #include "complex.hh" @@ -30,6 +31,7 @@ struct Encoder static const int symbol_len = (1280 * rate) / 8000; static const int guard_len = symbol_len / 8; static const int data_bits = 43040; + static const int data_bytes = data_bits / 8; static const int crc_bits = data_bits + 32; static const int mls0_len = 127; static const int mls0_poly = 0b10001001; @@ -266,7 +268,7 @@ struct Encoder cons_rows = cons_cnt / cons_cols; return true; } - Encoder(DSP::WritePCM *pcm, const uint8_t *inp, int freq_off, uint64_t call_sign, int oper_mode) : + Encoder(DSP::WritePCM *pcm, const uint8_t *inp, int count, int freq_off, uint64_t call_sign, int oper_mode) : pcm(pcm), crc0(0xA8F4), crc1(0xD419CC15), bchenc({ 0b100011101, 0b101110111, 0b111110011, 0b101101001, 0b110111101, 0b111100111, 0b100101011, 0b111010111, @@ -284,25 +286,27 @@ struct Encoder 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); - pilot_block(); - 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); - for (int i = crc_bits; i < mesg_bits; ++i) - mesg[i] = 1; - polarenc(code, mesg, frozen_bits, code_order); - shorten(); - 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 k = 0; k < count; ++k) { + schmidl_cox(); + meta_data((call_sign << 8) | oper_mode); + pilot_block(); + for (int i = 0; i < data_bits; ++i) + mesg[i] = nrz(CODE::get_le_bit(inp+k*data_bytes, i)); + crc1.reset(); + for (int i = 0; i < data_bytes; ++i) + crc1(inp[i+k*data_bytes]); + for (int i = 0; i < 32; ++i) + mesg[i+data_bits] = nrz((crc1()>>i)&1); + for (int i = crc_bits; i < mesg_bits; ++i) + mesg[i] = 1; + polarenc(code, mesg, frozen_bits, code_order); + shorten(); + 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 < symbol_len; ++i) fdom[i] = 0; @@ -332,8 +336,8 @@ long long int base37_encoder(const char *str) int main(int argc, char **argv) { - if (argc < 6 || argc > 9) { - std::cerr << "usage: " << argv[0] << " OUTPUT RATE BITS CHANNELS INPUT [OFFSET] [CALLSIGN] [MODE]" << std::endl; + if (argc < 9) { + std::cerr << "usage: " << argv[0] << " OUTPUT RATE BITS CHANNELS OFFSET MODE CALLSIGN INPUT.." << std::endl; return 1; } @@ -343,30 +347,18 @@ int main(int argc, char **argv) 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]; - 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; - } - - int oper_mode = 6; - if (argc >= 9) - oper_mode = std::atoi(argv[8]); + int freq_off = std::atoi(argv[5]); + int oper_mode = std::atoi(argv[6]); if (oper_mode < 6 || oper_mode > 13) { std::cerr << "Unsupported operation mode." << std::endl; return 1; } + long long int call_sign = base37_encoder(argv[7]); + if (call_sign <= 0 || call_sign >= 129961739795077L) { + std::cerr << "Unsupported call sign." << std::endl; + return 1; + } int band_width; switch (oper_mode) { @@ -407,33 +399,40 @@ int main(int argc, char **argv) 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 input_count = argc - 8; const int data_len = 43040 / 8; - uint8_t *input_data = new uint8_t[data_len]; - for (int i = 0; i < data_len; ++i) - input_data[i] = input_file.get(); - CODE::Xorshift32 scrambler; - for (int i = 0; i < data_len; ++i) - input_data[i] ^= scrambler(); + uint8_t *input_data = new uint8_t[data_len*input_count]; + for (int j = 0; j < input_count; ++j) { + const char *input_name = argv[j+8]; + if (argc == 9 && input_name[0] == '-' && input_name[1] == 0) + input_name = "/dev/stdin"; + 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; + } + for (int i = 0; i < data_len; ++i) + input_data[j*data_len+i] = input_file.get(); + CODE::Xorshift32 scrambler; + for (int i = 0; i < data_len; ++i) + input_data[j*data_len+i] ^= scrambler(); + } 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); + delete new Encoder(&output_file, input_data, input_count, freq_off, call_sign, oper_mode); break; case 16000: - delete new Encoder(&output_file, input_data, freq_off, call_sign, oper_mode); + delete new Encoder(&output_file, input_data, input_count, freq_off, call_sign, oper_mode); break; case 44100: - delete new Encoder(&output_file, input_data, freq_off, call_sign, oper_mode); + delete new Encoder(&output_file, input_data, input_count, freq_off, call_sign, oper_mode); break; case 48000: - delete new Encoder(&output_file, input_data, freq_off, call_sign, oper_mode); + delete new Encoder(&output_file, input_data, input_count, freq_off, call_sign, oper_mode); break; default: std::cerr << "Unsupported sample rate." << std::endl;