From 59ff35e5b6a81150a87cc0b5a972a91bd64c3ab9 Mon Sep 17 00:00:00 2001 From: andreas128 Date: Mon, 29 May 2017 21:55:44 +0100 Subject: Add subsample_alignment and it's test --- calc_lag.py | 44 ++++++++++++++++++++++++++++++++++++++ grc/live_analyse_dab_poly.py | 51 ++++++++++++++++++-------------------------- src/dab_util.py | 40 ++++++++++++++++++++-------------- src/dab_util_test.py | 49 ++++++++++++++++++++++++++++++++++++++++++ src/subsample_align.py | 13 +++++------ 5 files changed, 145 insertions(+), 52 deletions(-) create mode 100644 calc_lag.py diff --git a/calc_lag.py b/calc_lag.py new file mode 100644 index 0000000..5f5cf70 --- /dev/null +++ b/calc_lag.py @@ -0,0 +1,44 @@ +from __future__ import absolute_import +from __future__ import division +from __future__ import print_function + +import argparse +import re +import sys +from tqdm import tqdm + +from glob import glob +from natsort import natsorted +import numpy as np +import pandas as pd +import matplotlib.pyplot as plt + +import src.dab_util as du + + + + +tx_paths = natsorted(glob(r'/home/andreas/dab/ODR-StaticPrecorrection/data/received/*_tx_record.iq')) +rx_paths = natsorted(glob(r'/home/andreas/dab/ODR-StaticPrecorrection/data/received/*_rx_record.iq')) + +res = [] + +for tx_path, rx_path in zip(tx_paths, rx_paths): + s1 = du.fromfile(tx_path) + s2 = du.fromfile(rx_path) + + res.append({ + '1':du.lag_upsampling(s2, s1, n_up=1), + '2':du.lag_upsampling(s2, s1, n_up=2), + '3':du.lag_upsampling(s2, s1, n_up=3), + '4':du.lag_upsampling(s2, s1, n_up=4), + '8':du.lag_upsampling(s2, s1, n_up=8), + '16':du.lag_upsampling(s2, s1, n_up=16), + '32':du.lag_upsampling(s2, s1, n_up=32), + }) + +df = pd.DataFrame(res) +df = df.reindex_axis(sorted(df.columns), axis=1) +print(df) +print(df.describe()) + diff --git a/grc/live_analyse_dab_poly.py b/grc/live_analyse_dab_poly.py index cc7672b..126c43d 100755 --- a/grc/live_analyse_dab_poly.py +++ b/grc/live_analyse_dab_poly.py @@ -3,7 +3,7 @@ ################################################## # GNU Radio Python Flow Graph # Title: Live Analyse Dab Poly -# Generated: Mon Apr 10 14:20:48 2017 +# Generated: Tue May 23 20:21:44 2017 ################################################## if __name__ == '__main__': @@ -60,13 +60,12 @@ class live_analyse_dab_poly(gr.top_block, Qt.QWidget): ################################################## # Variables ################################################## - self.samp_rate = samp_rate = 8e6 self.txgain = txgain = 80 self.shift_freq = shift_freq = 1 + self.samp_rate_2 = samp_rate_2 = 2048000 + self.samp_rate_1 = samp_rate_1 = 8192000 self.rxgain = rxgain = 10 self.freq = freq = 222e6 - self.f2 = f2 = samp_rate / 3.875 - self.f1 = f1 = samp_rate / 4 self.ampl = ampl = 0.4 self.a_8 = a_8 = 0 self.a_7 = a_7 = 0 @@ -122,7 +121,7 @@ class live_analyse_dab_poly(gr.top_block, Qt.QWidget): channels=range(1), ), ) - self.uhd_usrp_source_0.set_samp_rate(samp_rate) + self.uhd_usrp_source_0.set_samp_rate(samp_rate_2) self.uhd_usrp_source_0.set_center_freq(freq, 0) self.uhd_usrp_source_0.set_gain(rxgain, 0) self.uhd_usrp_sink_0 = uhd.usrp_sink( @@ -132,7 +131,7 @@ class live_analyse_dab_poly(gr.top_block, Qt.QWidget): channels=range(1), ), ) - self.uhd_usrp_sink_0.set_samp_rate(samp_rate) + self.uhd_usrp_sink_0.set_samp_rate(samp_rate_1) self.uhd_usrp_sink_0.set_center_freq(freq, 0) self.uhd_usrp_sink_0.set_gain(txgain, 0) self.uhd_amsg_source_0 = uhd.amsg_source(device_addr="", msgq=uhd_amsg_source_0_msgq_out) @@ -143,7 +142,7 @@ class live_analyse_dab_poly(gr.top_block, Qt.QWidget): 16000, #size firdes.WIN_BLACKMAN_hARRIS, #wintype 0, #fc - samp_rate, #bw + samp_rate_2, #bw "", #name 1 #number of inputs ) @@ -206,17 +205,6 @@ class live_analyse_dab_poly(gr.top_block, Qt.QWidget): event.accept() - def get_samp_rate(self): - return self.samp_rate - - def set_samp_rate(self, samp_rate): - self.samp_rate = samp_rate - self.set_f1(self.samp_rate / 4) - self.set_f2(self.samp_rate / 3.875) - self.qtgui_freq_sink_x_0_0.set_frequency_range(0, self.samp_rate) - self.uhd_usrp_sink_0.set_samp_rate(self.samp_rate) - self.uhd_usrp_source_0.set_samp_rate(self.samp_rate) - def get_txgain(self): return self.txgain @@ -231,6 +219,21 @@ class live_analyse_dab_poly(gr.top_block, Qt.QWidget): def set_shift_freq(self, shift_freq): self.shift_freq = shift_freq + def get_samp_rate_2(self): + return self.samp_rate_2 + + def set_samp_rate_2(self, samp_rate_2): + self.samp_rate_2 = samp_rate_2 + self.qtgui_freq_sink_x_0_0.set_frequency_range(0, self.samp_rate_2) + self.uhd_usrp_source_0.set_samp_rate(self.samp_rate_2) + + def get_samp_rate_1(self): + return self.samp_rate_1 + + def set_samp_rate_1(self, samp_rate_1): + self.samp_rate_1 = samp_rate_1 + self.uhd_usrp_sink_0.set_samp_rate(self.samp_rate_1) + def get_rxgain(self): return self.rxgain @@ -247,18 +250,6 @@ class live_analyse_dab_poly(gr.top_block, Qt.QWidget): self.uhd_usrp_sink_0.set_center_freq(self.freq, 0) self.uhd_usrp_source_0.set_center_freq(self.freq, 0) - def get_f2(self): - return self.f2 - - def set_f2(self, f2): - self.f2 = f2 - - def get_f1(self): - return self.f1 - - def set_f1(self, f1): - self.f1 = f1 - def get_ampl(self): return self.ampl diff --git a/src/dab_util.py b/src/dab_util.py index 617bd9a..3187036 100644 --- a/src/dab_util.py +++ b/src/dab_util.py @@ -2,6 +2,7 @@ import numpy as np import scipy import matplotlib.pyplot as plt import src.dabconst as dabconst +import src.subsample_align as sa from scipy import signal c = {} @@ -76,22 +77,26 @@ def lag_upsampling(sig_orig, sig_rec, n_up): l_orig = float(l) / n_up return l_orig -def fftlag(sig_orig, sig_rec, n_upsampling = 1): +def subsample_align(sig1, sig2): """ - Efficient way to find lag between two signals - Args: - sig_orig: The signal that has been sent - sig_rec: The signal that has been recored + Returns an aligned version of sig1 and sig2 by cropping and subsample alignment """ - #off = sig_rec.shape[0] - #fft1 = np.fft.fft(sig_orig, n=sig_orig.shape[0]) - #fft2 = np.fft.fft(np.flipud(sig_rec), n=sig_rec.shape[0]) - #fftc = fft1 * fft2 - #c = np.fft.ifft(fftc) - c = signal.convolve(sig_orig, np.flipud(sig_rec)) - #c = signal.correlate(sig_orig, sig_rec) - return c - return np.argmax(c) - off + 1 + off_meas = lag_upsampling(sig2, sig1, n_up=1) + off = int(abs(off_meas)) + + if off_meas > 0: + sig1 = sig1[:-off] + sig2 = sig2[off:] + elif off_meas < 0: + sig1 = sig1[off:] + sig2 = sig2[:-off] + + if off % 2 == 1: + sig1 = sig1[:-1] + sig2 = sig2[:-1] + + sig2 = sa.subsample_align(sig2, sig1) + return sig1, sig2 def get_amp_ratio(ampl_1, ampl_2, a_out_abs, a_in_abs): idxs = (a_in_abs > ampl_1) & (a_in_abs < ampl_2) @@ -108,5 +113,8 @@ def get_transmission_frame_indices(n_frames, offset, rate = 2048000): indices = [tm1.S_F * i + offset for i in range(n_frames)] return indices -def fromfile(filename, offset, length): - return np.memmap(filename, dtype=np.complex64, mode='r', offset=64/8*offset, shape=length) +def fromfile(filename, offset=0, length=None): + if length is None: + return np.memmap(filename, dtype=np.complex64, mode='r', offset=64/8*offset) + else: + return np.memmap(filename, dtype=np.complex64, mode='r', offset=64/8*offset, shape=length) diff --git a/src/dab_util_test.py b/src/dab_util_test.py index be36d53..3f9e941 100644 --- a/src/dab_util_test.py +++ b/src/dab_util_test.py @@ -1,5 +1,7 @@ from scipy import signal import numpy as np +import pandas as pd +from tqdm import tqdm import src.gen_source as gs import src.dab_util as du @@ -28,7 +30,54 @@ def test_phase_offset(lag_function, tol): res.append(np.abs(off-off_meas)