compensate using comb type pilots

This commit is contained in:
Ahmet Inan 2023-12-18 13:23:13 +01:00
commit b93079de81
2 changed files with 66 additions and 47 deletions

View file

@ -181,7 +181,11 @@ struct Decoder
static const int guard_len = symbol_len / 8;
static const int extended_len = symbol_len + guard_len;
static const int max_bits = 2720 + 32;
static const int cons_cols = 256;
static const int comb_cols = 8;
static const int code_cols = 256;
static const int cons_cols = code_cols + comb_cols;
static const int comb_dist = cons_cols / comb_cols;
static const int comb_off = comb_dist / 2;
static const int cons_rows = 4;
static const int cons_total = cons_rows * cons_cols;
static const int code_off = - cons_cols / 2;
@ -198,7 +202,7 @@ struct Decoder
DSP::BlockDC<value, value> blockdc;
DSP::Hilbert<cmplx, filter_len> hilbert;
DSP::BipBuffer<cmplx, buffer_len> input_hist;
DSP::TheilSenEstimator<value, cons_cols> tse;
DSP::TheilSenEstimator<value, comb_cols> tse;
SchmidlCox<value, cmplx, search_pos, symbol_len/2, guard_len> correlator;
CODE::CRC<uint16_t> crc0;
CODE::CRC<uint32_t> crc1;
@ -211,7 +215,7 @@ struct Decoder
code_type code[code_len];
cmplx cons[cons_total], prev[cons_cols];
cmplx fdom[symbol_len], tdom[symbol_len];
value index[cons_cols], phase[cons_cols];
value index[comb_cols], phase[comb_cols];
value cfo_rad, sfo_rad;
const uint32_t *frozen_bits;
int symbol_pos;
@ -361,6 +365,7 @@ struct Decoder
for (int i = 0; i < cons_cols; ++i)
prev[i] = fdom[bin(i+code_off)];
std::cerr << "demod ";
CODE::MLS seq0(mls0_poly);
for (int j = 0; j < cons_rows; ++j) {
for (int i = 0; i < extended_len; ++i)
buf = next_sample();
@ -371,48 +376,46 @@ struct Decoder
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]);
if (1) {
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;
for (int i = 0; i < cons_cols; ++i)
cons[cons_cols*j+i] *= DSP::polar<value>(1, -tse(i+code_off));
for (int i = 0; i < cons_cols; ++i)
prev[i] *= DSP::polar<value>(1, tse(i+code_off));
for (int i = 0; i < comb_cols; ++i)
cons[cons_cols*j+comb_dist*i+comb_off] *= nrz(seq0());
for (int i = 0; i < comb_cols; ++i) {
index[i] = code_off + comb_dist * i + comb_off;
phase[i] = arg(cons[cons_cols*j+comb_dist*i+comb_off]);
}
tse.compute(index, phase, comb_cols);
//std::cerr << "Theil-Sen slope = " << tse.slope() << std::endl;
//std::cerr << "Theil-Sen yint = " << tse.yint() << std::endl;
for (int i = 0; i < cons_cols; ++i)
cons[cons_cols*j+i] *= DSP::polar<value>(1, -tse(i+code_off));
for (int i = 0; i < cons_cols; ++i)
if (i % comb_dist == comb_off)
prev[i] = fdom[bin(i+code_off)];
else
prev[i] *= DSP::polar<value>(1, tse(i+code_off));
std::cerr << ".";
}
std::cerr << " done" << 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) {
code_type tmp[mod_bits];
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);
std::cerr << "Es/N0 (dB):";
value sp = 0, np = 0;
for (int j = 0, k = 0; j < cons_rows; ++j) {
for (int i = 0; i < comb_cols; ++i) {
cmplx hard(1, 0);
cmplx error = cons[cons_cols*j+comb_dist*i+comb_off] - hard;
sp += norm(hard);
np += norm(error);
}
value precision = sp / np;
// precision = 8;
value snr = DSP::decibel(precision);
std::cerr << " " << snr;
for (int i = 0; i < cons_cols; ++i) {
if (i % comb_dist == comb_off)
continue;
mod_soft(code+k, cons[cons_cols*j+i], precision);
k += mod_bits;
}
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);
}
std::cerr << std::endl;
int data_bits = 0;
switch (oper_mode) {
case 17:

View file

@ -35,7 +35,11 @@ struct Encoder
static const int symbol_len = (1280 * rate) / 8000;
static const int guard_len = symbol_len / 8;
static const int max_bits = 2720 + 32;
static const int cons_cols = 256;
static const int comb_cols = 8;
static const int code_cols = 256;
static const int cons_cols = code_cols + comb_cols;
static const int comb_dist = cons_cols / comb_cols;
static const int comb_off = comb_dist / 2;
static const int cons_rows = 4;
static const int mls0_len = 127;
static const int mls0_poly = 0b10001001;
@ -176,10 +180,10 @@ struct Encoder
CODE::set_be_bit(data, i+55, (cs>>i)&1);
bchenc(data, parity);
CODE::MLS seq1(mls1_poly);
value mls1_fac = std::sqrt(value(symbol_len) / value(mls1_len));
value cons_fac = std::sqrt(value(symbol_len) / value(cons_cols));
for (int i = 0; i < symbol_len; ++i)
fdom[i] = 0;
fdom[bin(mls1_off-1)] = mls1_fac;
fdom[bin(mls1_off-1)] = cons_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)
@ -188,6 +192,10 @@ struct Encoder
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(seq1());
for (int i = 0; i < comb_cols / 2; ++i) {
fdom[bin(i+code_off)] = cons_fac * nrz(seq1());
fdom[bin(i+mls1_off+mls1_len)] = cons_fac * nrz(seq1());
}
symbol();
}
cmplx mod_map(code_type *b)
@ -247,10 +255,17 @@ struct Encoder
shuffle(code);
for (int i = 0; i < cons_cols; ++i)
prev[i] = fdom[bin(i+code_off)];
for (int j = 0; j < cons_rows; ++j) {
for (int i = 0; i < cons_cols; ++i)
fdom[bin(i+code_off)] = prev[i] *
mod_map(code+mod_bits*(cons_cols*j+i));
CODE::MLS seq0(mls0_poly);
for (int j = 0, k = 0; j < cons_rows; ++j) {
for (int i = 0; i < cons_cols; ++i) {
if (i % comb_dist == comb_off) {
prev[i] *= nrz(seq0());
fdom[bin(i+code_off)] = prev[i];
} else {
fdom[bin(i+code_off)] = prev[i] * mod_map(code+k);
k += mod_bits;
}
}
symbol();
}
}
@ -308,8 +323,9 @@ int main(int argc, char **argv)
return 1;
}
const int comb_pilots = 8;
const int reserved_tones = 8;
const int band_width = 1600 + (25 * reserved_tones) / 4;
const int band_width = 1600 + (25 * (comb_pilots + reserved_tones)) / 4;
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;