diff --git a/lib/adaptive_dfe_impl.cc b/lib/adaptive_dfe_impl.cc index 455bb37..1035e46 100644 --- a/lib/adaptive_dfe_impl.cc +++ b/lib/adaptive_dfe_impl.cc @@ -71,8 +71,6 @@ adaptive_dfe_impl::adaptive_dfe_impl(int sps, // samples per symbol , _mu(mu) , _alpha(alpha) , _use_symbol_taps(true) - // , _py_module_name(python_module_name) - // , _physicalLayer() , _taps_samples(nullptr) , _taps_symbols(nullptr) , _hist_symbols(nullptr) @@ -95,7 +93,7 @@ adaptive_dfe_impl::adaptive_dfe_impl(int sps, // samples per symbol , _num_samples_since_filter_update(0) , _rotated_samples() , _rotator() - , _control_loop(2*M_PI/100, 5e-5, -5e-5) + , _control_loop(2*M_PI/100, 5e-2, -5e-2) { GR_LOG_DECLARE_LOGPTR(d_logger); GR_LOG_ASSIGN_LOGPTR(d_logger, "adaptive_dfe"); @@ -136,7 +134,6 @@ adaptive_dfe_impl::general_work(int noutput_items, gr_vector_void_star &output_items) { gr::thread::scoped_lock lock(d_setlock); - //GR_LOG_DEBUG(d_logger, str(boost::format("work: %d") % noutput_items)); gr_complex const* in = (gr_complex const *)input_items[0]; gr_complex *out_symb = (gr_complex *)output_items[0]; gr_complex *out_taps = (gr_complex *)output_items[1]; @@ -149,7 +146,6 @@ adaptive_dfe_impl::general_work(int noutput_items, return 0; int const ninput = ninput_items[0] - _nGuard - _nF; - int nout = 0; // counter for produced output items switch (_state) { case WAIT_FOR_PREAMBLE: { @@ -159,13 +155,6 @@ adaptive_dfe_impl::general_work(int noutput_items, consume(0, ninput - history()+1); } else { tag_t const& tag = v.front(); - // uint64_t const offset = tag.offset - nitems_read(0) + history() - 1; - // std::cout << "========= offset= " << offset - // << " tag.offset= " << tag.offset - // << " nitems_read(0)= " << nitems_read(0) - // << " tag.offset-nitems_read(0)= " << tag.offset - nitems_read(0) << " ==========" << std::endl; - // for (int k=0; k= 0); - out_symb[nout] = filter(&_rotated_samples[0] + i - _nB, - &_rotated_samples[0] + i + _nF+1); + out_symb[nout] = filter(&_rotated_samples.front() + i - _nB, + &_rotated_samples.front() + i + _nF+1); std::memcpy(&out_taps[(_nB+_nF+1)*nout], _taps_samples, (_nB+_nF+1)*sizeof(gr_complex)); ++nout; } // next sample @@ -247,8 +244,6 @@ bool adaptive_dfe_impl::stop() gr_complex adaptive_dfe_impl::filter(gr_complex const* start, gr_complex const* end) { assert(end-start == _nB + _nF + 1); - _num_samples_since_filter_update += _sps; - // (1) run the filter filter gr_complex filter_output(0); // (1a) taps_samples @@ -259,9 +254,7 @@ gr_complex adaptive_dfe_impl::filter(gr_complex const* start, gr_complex const* // (1b) taps_symbols gr_complex dot_symbols(0); gr::digital::constellation_sptr constell = _constellations[_constellation_index]; - bool const update_taps = constell->bits_per_symbol() <= 3; - if (constell->bits_per_symbol() > 3) - _use_symbol_taps = false; + _use_symbol_taps = (constell->bits_per_symbol() <= 3); if (_use_symbol_taps) { for (int l=0; l<_nW; ++l) { assert(_hist_symbol_index+l < 2*_nW); @@ -273,6 +266,7 @@ gr_complex adaptive_dfe_impl::filter(gr_complex const* start, gr_complex const* gr_complex known_symbol = _symbols[_symbol_counter]; bool const is_known = std::abs(known_symbol) > 1e-5; + bool const update_taps = constell->bits_per_symbol() <= 3 || is_known; // (2) unknown symbols (=data): compute soft decisions if (not is_known) { gr_complex const descrambled_filter_output = std::conj(_scramble[_symbol_counter]) * filter_output; @@ -283,21 +277,18 @@ gr_complex adaptive_dfe_impl::filter(gr_complex const* start, gr_complex const* float const err = std::abs(descrambled_filter_output - descrambled_symbol); std::vector const soft_dec = constell->calc_soft_dec (descrambled_filter_output, _npwr[_constellation_index].filter(err)); - // std::cout << "SD: " << descrambled_filter_output << " : "; - for (int j=0, m=soft_dec.size(); j0.5) + if (std::abs(err)>0.7) std::cout << "err= " << std::abs(err) << std::endl; // taps_samples for (int j=0; j<_nB+_nF+1; ++j) { @@ -316,16 +307,23 @@ gr_complex adaptive_dfe_impl::filter(gr_complex const* start, gr_complex const* } } // (3b) control loop update for doppler correction using the adaptibve filter taps - if (_symbol_counter+1 == _symbols.size()) { - gr_complex acc(0); - for (int j=_nB+1-2*_sps; j<_nB+1+2*_sps+1; ++j) - acc += std::conj(_last_taps_samples[j]) * _taps_samples[j]; - float const frequency_err = gr::fast_atan2f(acc)/(1+0*_num_samples_since_filter_update); // frequency error (rad/sample) - GR_LOG_DEBUG(d_logger, str(boost::format("frequency_err= %f %d") % frequency_err % _num_samples_since_filter_update)); - _control_loop.advance_loop(frequency_err); - _control_loop.phase_wrap(); - _control_loop.frequency_limit(); - _rotator.set_phase_incr(gr_expj(_control_loop.get_frequency())); + if (update_taps) { + if (_symbol_counter != 0) { // a filter tap shift might have ocurred when _symbol_counter==0 + gr_complex acc(0); + for (int j=0; j<_nB+_nF+1; ++j) { + acc += std::conj(_last_taps_samples[j]) * _taps_samples[j]; + } + float const frequency_err = gr::fast_atan2f(acc)/(0+1*_num_samples_since_filter_update); // frequency error (rad/sample) + GR_LOG_DEBUG(d_logger, str(boost::format("frequency_err= %f %d") % frequency_err % _num_samples_since_filter_update)); + _control_loop.advance_loop(frequency_err); + _control_loop.phase_wrap(); + _control_loop.frequency_limit(); + _rotator.set_phase_incr(gr_expj(_control_loop.get_frequency())); + GR_LOG_DEBUG(d_logger, str(boost::format("frequency_err= %f %d %f") + % (frequency_err/(2*M_PI)*12000.0) + % _num_samples_since_filter_update + % _control_loop.get_frequency())); + } _num_samples_since_filter_update = 0; } @@ -336,29 +334,38 @@ gr_complex adaptive_dfe_impl::filter(gr_complex const* start, gr_complex const* int adaptive_dfe_impl::recenter_filter_taps() { - // get max(abs(taps)) - ssize_t const idx_max = std::distance(_taps_samples, +#if 0 + ssize_t const _idx_max = std::distance(_taps_samples, std::max_element(_taps_samples+_nB+1-3*_sps, _taps_samples+_nB+1+3*_sps, [](gr_complex a, gr_complex b) { return std::norm(a) < std::norm(b); })); +#else + float sum_w=0, sum_wi=0; + for (int i=0; i<_nB+_nF+1; ++i) { + float const w = std::norm(_taps_samples[i]); + sum_w += w; + sum_wi += w*i; + } + ssize_t const idx_max = ssize_t(0.5 + sum_wi/sum_w); +#endif // GR_LOG_DEBUG(d_logger, str(boost::format("idx_max=%2d abs(tap_max)=%f") % idx_max % std::abs(_taps_samples[idx_max]))); if (idx_max-_nB-1 > +2*_sps) { // maximum is right of the center tap // -> shift taps to the left left GR_LOG_DEBUG(d_logger, "shift left"); - std::copy(_taps_samples+2*_sps, _taps_samples+_nB+_nF+1, _taps_samples); - std::fill_n(_taps_samples+_nB+_nF+1-2*_sps, 2*_sps, gr_complex(0)); - return +2*_sps; + std::copy(_taps_samples+4*_sps, _taps_samples+_nB+_nF+1, _taps_samples); + std::fill_n(_taps_samples+_nB+_nF+1-4*_sps, 4*_sps, gr_complex(0)); + return +4*_sps; } if (idx_max-_nB-1 < -2*_sps) { // maximum is left of the center tap // -> shift taps to the right GR_LOG_DEBUG(d_logger, "shift right"); - std::copy_backward(_taps_samples, _taps_samples+_nB+_nF+1-2*_sps, + std::copy_backward(_taps_samples, _taps_samples+_nB+_nF+1-4*_sps, _taps_samples+_nB+_nF+1); - std::fill_n(_taps_samples, 2*_sps, gr_complex(0)); - return -2*_sps; + std::fill_n(_taps_samples, 4*_sps, gr_complex(0)); + return -4*_sps; } return 0; } @@ -378,10 +385,13 @@ void adaptive_dfe_impl::publish_frame_info() { pmt::pmt_t data = pmt::make_dict(); GR_LOG_DEBUG(d_logger, str(boost::format("publish_frame_info %d == %d") % _descrambled_symbols.size() % _symbols.size())); - data = pmt::dict_add(data, pmt::intern("symbols"), pmt::init_c32vector(_descrambled_symbols.size(), &_descrambled_symbols.front())); + data = pmt::dict_add(data, + pmt::intern("symbols"), + pmt::init_c32vector(_descrambled_symbols.size(), &_descrambled_symbols.front())); // for (int i=0; i<_vec_soft_decisions.size(); ++i) // _vec_soft_decisions[i] = std::max(-1.0f, std::min(1.0f, _vec_soft_decisions[i])); - data = pmt::dict_add(data, pmt::intern("soft_dec"), pmt::init_f32vector(_vec_soft_decisions.size(), &_vec_soft_decisions.front())); + data = pmt::dict_add(data, + pmt::intern("soft_dec"), pmt::init_f32vector(_vec_soft_decisions.size(), &_vec_soft_decisions.front())); message_port_pub(_msg_ports["frame_info"], data); _descrambled_symbols.clear(); } diff --git a/lib/doppler_correction_cc_impl.cc b/lib/doppler_correction_cc_impl.cc index 0bf0040..066dbec 100644 --- a/lib/doppler_correction_cc_impl.cc +++ b/lib/doppler_correction_cc_impl.cc @@ -150,10 +150,21 @@ doppler_correction_cc_impl::work(int noutput_items, } // CONSUME_AND_SKIP } // apply current doppler correction to all produced samples +#if 0 // rotateN is broken in some older VOLK versions _rotator.rotateN(out, in, nout); +#else + for (int i=0; i 2.0 + r['success'] = bool(np.median(tests) > 2.0) print('test:', np.abs(pks), tests) - if success: + if r['success']: print('doppler apks', np.abs(pks)) print('doppler ppks', np.angle(pks), np.diff(np.unwrap(np.angle(pks)))/m, np.mean(np.diff(np.unwrap(np.angle(pks)))/m)) - doppler = common.freq_est(pks)/m; - print('success=', success, 'doppler=', doppler) - return success,doppler + r['doppler'] = common.freq_est(pks)/m + print(r) + return r def set_mode(self, mode): pass + def get_mode(self): + return self._mode_description + def get_preamble_quality(self, symbols): print('get_preamble_quality', np.abs(np.mean(symbols[-32:])), symbols[-32:]) return np.abs(np.mean(symbols[-32:])) > 0.5 @@ -354,7 +359,6 @@ class PhysicalLayer(object): def is_HFXL(self): return self._mode_name == 'HFXL' - def decode_reinserted_preamble(self, symbols): ## decode D0,D1,D2 success = True @@ -385,13 +389,14 @@ class PhysicalLayer(object): if rate_info['baud'] == 'HFXL': self._mode_name = 'HFXL' + self._mode_description = '%s rate=%s intl=%s' % (self._mode_name, rate_info['baud'], intl_info['id']) + print('======== rate,interleaver:', rate_info, intl_info, self._mode_name) self._data_scramble_xor = np.zeros(256, dtype=np.uint8) self._data_scramble = np.ones (256, dtype=np.complex64) if self.is_12800bpsBurst(): self._scrp = ScrambleDataP() - self._constellation_index = MODE_BPSK# 64QAMp - ##self._data_scramble = QAM64p['points'][self._scrp.next() for _ in range(256)] + self._constellation_index = MODE_BPSK elif self.is_HFXL(): self._scramble.reset() num_bits = 3 @@ -511,23 +516,25 @@ class PhysicalLayer(object): soft_bits[2*i] = abs_soft_dec*(2*(b>>1)-1) soft_bits[2*i+1] = abs_soft_dec*(2*(b &1)-1) - return soft_bits>0 + return soft_bits>0,100.0 elif self.is_HFXL(): ## TODO - return [] + return np.zeros(0, dtype=np.float32),0.0 + elif self.is_plain_110C(): r = self._deintl_depunct.load(soft_dec) if r.shape[0] == 0: - return [] + return np.zeros(0, dtype=np.float32),0.0 self._viterbi_decoder.reset() decoded_bits = np.roll(self._viterbi_decoder.udpate(r), 7) print('bits=', decoded_bits[:100]) - print('quality={}% ({},{})'.format(120.0*self._viterbi_decoder.quality()/(2*len(decoded_bits)), + quality = 120.0*self._viterbi_decoder.quality()/(2*len(decoded_bits)) + print('quality={}% ({},{})'.format(quality, self._viterbi_decoder.quality(), len(decoded_bits))) - return decoded_bits + return decoded_bits,quality else: - return [] + return np.zeros(0, dtype=np.float32),0.0 @staticmethod def get_preamble(): diff --git a/python/physical_layer/STANAG_4285.py b/python/physical_layer/STANAG_4285.py index 93883af..fbb4fd7 100644 --- a/python/physical_layer/STANAG_4285.py +++ b/python/physical_layer/STANAG_4285.py @@ -43,7 +43,6 @@ class PhysicalLayer(object): def __init__(self, sps): """intialization""" self._sps = sps - ##self._mode = self.MODE_QPSK self._frame_counter = 0 self._is_first_frame = True self._constellations = [self.make_psk(2, [0,1]), @@ -52,16 +51,20 @@ class PhysicalLayer(object): self._preamble = self.get_preamble() self._data = self.get_data() self._viterbi_decoder = viterbi27(0x6d, 0x4f) + self._mode_description = None def set_mode(self, mode): - """set phase modultation type: 'BPS/S' or 'BPS/L'""" - print('set_mode', mode) + """set modulation and interleaver: 'BPS/S' or 'BPS/L'""" + self._mode_description = mode bps,intl = mode.split('/') self._mode = MODES[bps]['const'] self._deinterleaver = Deinterleaver(DEINTERLEAVER_INCR[intl] * MODES[bps]['deintl_multiple']) self._depuncturer = common.Depuncturer(repeat = MODES[bps]['repeat'], puncture_pattern = MODES[bps]['punct']) + def get_mode(self): + return self._mode_description + def get_constellations(self): return self._constellations @@ -71,31 +74,26 @@ class PhysicalLayer(object): [1] ... modulation type after descrambling [2] ... a boolean indicating if the processing should continue [3] ... a boolean indicating if the soft decision for the unknown symbols are saved""" - ## print('-------------------- get_frame --------------------', self._frame_counter, len(symbols)) if len(symbols) == 0: ## 1st preamble self._frame_counter = 0 success,frame_description = True,[] - if (self._frame_counter%2) == 0: + if (self._frame_counter%2) == 0: ## current frame is a data frame frame_description = [self._preamble,MODE_BPSK,success,False] - else: + else: ## current frame is a preamble frame idx = range(30,80) z = symbols[idx]*np.conj(self._preamble['symb'][idx]) - ## print('quality_preamble',np.sum(np.real(z)<0), symbols[idx]) - success = np.sum(np.real(z)<0) < 30 + success = bool(np.sum(np.real(z)<0) < 30) frame_description = [self._data,self._mode,success,True] self._frame_counter += 1 return frame_description def get_doppler(self, iq_samples): - """returns a tuple - [0] ... quality flag - [1] ... doppler estimate (rad/symbol) if available""" - ## print('-------------------- get_doppler --------------------', self._frame_counter,len(iq_samples)) - success,doppler = False,0 + r = {'success': False, ## -- quality flag + 'doppler': 0} ## -- doppler estimate (rad/symb) if len(iq_samples) == 0: - return success,doppler + return r sps = self._sps zp = np.array([x for x in self._preamble['symb'][9:40] @@ -104,20 +102,10 @@ class PhysicalLayer(object): imax = np.argmax(np.abs(cc[0:18*sps])) pks = cc[(imax,imax+31*sps),] tpks = cc[imax+15*sps:imax+16*sps] - ## print('doppler: ', np.abs(pks), np.abs(tpks)) - success = np.mean(np.abs(pks)) > 5*np.mean(np.abs(tpks)) - doppler = np.diff(np.unwrap(np.angle(pks)))[0]/31/self._sps if success else 0 - return success,doppler - def is_preamble(self): - return self._frame_counter == 0 - - def quality_data(self, s): - """quality check for the data frame""" - known_symbols = np.mod(range(176),48)>=32 - print('quality_data',np.sum(np.real(s[known_symbols])<0)) - success = np.sum(np.real(s[known_symbols])<0) < 20 - return success,0 ## no doppler estimate for data frames + r['success'] = bool(np.mean(np.abs(pks)) > 5*np.mean(np.abs(tpks))) + r['doppler'] = np.diff(np.unwrap(np.angle(pks)))[0]/31/self._sps if r['success'] else 0 + return r def get_preamble_z(self): """preamble symbols for preamble correlation""" @@ -132,9 +120,8 @@ class PhysicalLayer(object): r.extend(self._deinterleaver.fetch().tolist()) rd = self._depuncturer.process(np.array(r, dtype=np.float32)) decoded_bits = self._viterbi_decoder.udpate(rd) - print('bits=', decoded_bits) - print('quality={}%'.format(100.0*self._viterbi_decoder.quality()/(2*len(decoded_bits)))) - return decoded_bits + quality = 100.0*self._viterbi_decoder.quality()/(2*len(decoded_bits)) + return decoded_bits,quality @staticmethod def get_preamble(): @@ -143,8 +130,8 @@ class PhysicalLayer(object): taps = np.array([0,0,1,0,1], dtype=np.bool) p = np.zeros(80, dtype=np.uint8) for i in range(80): - p[i] = state[-1] - state = np.concatenate(([np.sum(state&taps)&1], state[0:-1])) + p[i] = state[-1] + state = np.concatenate(([np.sum(state&taps)&1], state[0:-1])) a = np.zeros(80, common.SYMB_SCRAMBLE_DTYPE) ## BPSK modulation constellation = PhysicalLayer.make_psk(2,range(2))['points'] @@ -156,7 +143,7 @@ class PhysicalLayer(object): def get_data(): """data symbols + scrambler; for unknown symbols 'symb'=0""" state = np.array([1,1,1,1,1,1,1,1,1], dtype=np.bool) - taps = np.array([0,0,0,0,1,0,0,0,1], dtype=np.bool) + taps = np.array([0,0,0,0,1,0,0,0,1], dtype=np.bool) p = np.zeros(176, dtype=np.uint8) for i in range(176): p[i] = np.sum(state[-3:]*[4,2,1]) diff --git a/python/physical_layer_driver.py b/python/physical_layer_driver.py index 09fa0cb..64f305c 100644 --- a/python/physical_layer_driver.py +++ b/python/physical_layer_driver.py @@ -103,3 +103,9 @@ class physical_layer_driver(gr.hier_block2): def set_mode(self, mode): self._physical_layer_driver_description.set_mode(mode) + + def get_quality(self): + return self._msg_proxy.get_quality() + + def get_mode(self): + return self._physical_layer_driver_description.get_mode()