diff options
Diffstat (limited to 'dpd/src')
-rw-r--r-- | dpd/src/Adapt.py | 2 | ||||
-rw-r--r-- | dpd/src/Const.py | 81 | ||||
-rw-r--r-- | dpd/src/Dab_Util.py | 44 | ||||
-rw-r--r-- | dpd/src/ExtractStatistic.py | 78 | ||||
-rw-r--r-- | dpd/src/Measure_Shoulders.py | 4 | ||||
-rw-r--r-- | dpd/src/Model_AM.py | 48 | ||||
-rw-r--r-- | dpd/src/RX_Agc.py (renamed from dpd/src/Agc.py) | 0 | ||||
-rw-r--r-- | dpd/src/phase_align.py | 2 | ||||
-rw-r--r-- | dpd/src/test_dab_Util.py | 2 |
9 files changed, 132 insertions, 129 deletions
diff --git a/dpd/src/Adapt.py b/dpd/src/Adapt.py index cf8b6f5..1ab7710 100644 --- a/dpd/src/Adapt.py +++ b/dpd/src/Adapt.py @@ -219,7 +219,7 @@ class Adapt: "txgain": self.get_txgain(), "rxgain": self.get_rxgain(), "digital_gain": self.get_digital_gain(), - "predistorter": self.get_predistorter + "predistorter": self.get_predistorter() } with open(path, "w") as f: pickle.dump(d, f) diff --git a/dpd/src/Const.py b/dpd/src/Const.py index ebdfeb2..390dfbb 100644 --- a/dpd/src/Const.py +++ b/dpd/src/Const.py @@ -1,16 +1,21 @@ -# DAB Frame constants -# Sources: -# - etsi_EN_300_401_v010401p p145 -# - Measured with USRP B200 +# -*- coding: utf-8 -*- +# +# DPD Calculation Engine, constants +# +# Source for DAB standard: etsi_EN_300_401_v010401p p145 +# +# http://www.opendigitalradio.org +# Licence: The MIT License, see notice at the end of this file import numpy as np + class Const: def __init__(self, sample_rate, target_median, plot): + assert sample_rate == 8192000 # By now only constants for 8192000 self.sample_rate = sample_rate - self.tx_gain_max = 89 - + # DAB frame # Time domain self.T_F = sample_rate / 2048000 * 196608 # Transmission frame duration self.T_NULL = sample_rate / 2048000 * 2656 # Null symbol duration @@ -20,18 +25,10 @@ class Const: # Frequency Domain # example: np.delete(fft[3328:4865], 768) - self.FFT_delete = 768 self.FFT_delta = 1536 # Number of carrier frequencies - if sample_rate == 2048000: - self.FFT_start = 256 - self.FFT_end = 1793 - elif sample_rate == 8192000: - self.FFT_start = 3328 - self.FFT_end = 4865 - else: - raise RuntimeError("Sample Rate '{}' not supported".format( - sample_rate - )) + self.FFT_delete = 768 + self.FFT_start = 3328 + self.FFT_end = 4865 # Calculate sample offset from phase rotation # time per sample = 1 / sample_rate @@ -43,42 +40,44 @@ class Const: self.ES_plot = plot self.ES_start = 0.0 self.ES_end = 1.0 - self.ES_n_bins = 64 - self.ES_n_per_bin = 128 - - # Constants for TX_Agc - self.TAGC_max_txgain = 89 - self.TAGC_tx_median_target = target_median - self.TAGC_tx_median_max = self.TAGC_tx_median_target*1.4 - self.TAGC_tx_median_min = self.TAGC_tx_median_target/1.4 - - self.RAGC_min_rxgain = 25 - self.RAGC_rx_median_target = self.TAGC_tx_median_target - - # Constants for Model - self.MDL_plot = True or plot # Override default - - # Constants for MER - self.MER_plot = plot - - # Constants for Model_PM - self.MPM_tx_min = 0.1 + self.ES_n_bins = 64 # Number of bins between ES_start and ES_end + self.ES_n_per_bin = 128 # Number of measurements pre bin # Constants for Measure_Shoulder self.MS_enable = False self.MS_plot = plot - assert sample_rate==8192000 - meas_offset = 976 # Offset from center frequency to measure shoulder [kHz] - meas_width = 100 # Size of frequency delta to measure shoulder [kHz] + + meas_offset = 976 # Offset from center frequency to measure shoulder [kHz] + meas_width = 100 # Size of frequency delta to measure shoulder [kHz] shoulder_offset_edge = np.abs(meas_offset - self.FFT_delta) self.MS_shoulder_left_start = self.FFT_start - shoulder_offset_edge - meas_width / 2 self.MS_shoulder_left_end = self.FFT_start - shoulder_offset_edge + meas_width / 2 self.MS_shoulder_right_start = self.FFT_end + shoulder_offset_edge - meas_width / 2 self.MS_shoulder_right_end = self.FFT_end + shoulder_offset_edge + meas_width / 2 - self.MS_peak_start = self.FFT_start + 100 # Ignore region near edges + self.MS_peak_start = self.FFT_start + 100 # Ignore region near edges self.MS_peak_end = self.FFT_end - 100 self.MS_FFT_size = 8192 self.MS_averaging_size = 4 * self.MS_FFT_size self.MS_n_averaging = 40 self.MS_n_proc = 4 + + # Constants for MER + self.MER_plot = plot + + # Constants for Model + self.MDL_plot = True or plot # Override default + + # Constants for Model_PM + # Set all phase offsets to zero for TX amplitude < MPM_tx_min + self.MPM_tx_min = 0.1 + + # Constants for TX_Agc + self.TAGC_max_txgain = 89 # USRP specific + self.TAGC_tx_median_target = target_median + self.TAGC_tx_median_max = self.TAGC_tx_median_target * 1.4 + self.TAGC_tx_median_min = self.TAGC_tx_median_target / 1.4 + + # Constants for RX_AGC + self.RAGC_min_rxgain = 25 # USRP specific + self.RAGC_rx_median_target = self.TAGC_tx_median_target diff --git a/dpd/src/Dab_Util.py b/dpd/src/Dab_Util.py index 37be5db..27a31ef 100644 --- a/dpd/src/Dab_Util.py +++ b/dpd/src/Dab_Util.py @@ -8,30 +8,40 @@ import datetime import os import logging + logging_path = os.path.dirname(logging.getLoggerClass().root.handlers[0].baseFilename) import numpy as np -import scipy import matplotlib + matplotlib.use('agg') import matplotlib.pyplot as plt import src.subsample_align as sa import src.phase_align as pa from scipy import signal + +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) + + class Dab_Util: """Collection of methods that can be applied to an array complex IQ samples of a DAB signal """ + def __init__(self, sample_rate, plot=False): """ :param sample_rate: sample rate [sample/sec] to use for calculations """ self.sample_rate = sample_rate - self.dab_bandwidth = 1536000 #Bandwidth of a dab signal - self.frame_ms = 96 #Duration of a Dab frame + self.dab_bandwidth = 1536000 # Bandwidth of a dab signal + self.frame_ms = 96 # Duration of a Dab frame - self.plot=plot + self.plot = plot def lag(self, sig_orig, sig_rec): """ @@ -56,10 +66,10 @@ class Dab_Util: def lag_upsampling(self, sig_orig, sig_rec, n_up): if n_up != 1: sig_orig_up = signal.resample(sig_orig, sig_orig.shape[0] * n_up) - sig_rec_up = signal.resample(sig_rec, sig_rec.shape[0] * n_up) + sig_rec_up = signal.resample(sig_rec, sig_rec.shape[0] * n_up) else: sig_orig_up = sig_orig - sig_rec_up = sig_rec + sig_rec_up = sig_rec l = self.lag(sig_orig_up, sig_rec_up) l_orig = float(l) / n_up return l_orig @@ -69,7 +79,7 @@ class Dab_Util: Returns an aligned version of sig_tx and sig_rx by cropping and subsample alignment Using upsampling """ - assert(sig_tx.shape[0] == sig_rx.shape[0]) + assert (sig_tx.shape[0] == sig_rx.shape[0]) if sig_tx.shape[0] % 2 == 1: sig_tx = sig_tx[:-1] @@ -124,11 +134,11 @@ class Dab_Util: off = int(abs(off_meas)) logging.debug("sig_tx_orig: {} {}, sig_rx_orig: {} {}, offset {}".format( - len(sig_tx), - sig_tx.dtype, - len(sig_rx), - sig_rx.dtype, - off_meas)) + len(sig_tx), + sig_tx.dtype, + len(sig_rx), + sig_rx.dtype, + off_meas)) if off_meas > 0: sig_tx = sig_tx[:-off] @@ -211,16 +221,10 @@ class Dab_Util: fig.savefig(fig_path) plt.close(fig) - logging.debug("Sig1_cut: %d %s, Sig2_cut: %d %s, off: %d" % (len(sig_tx), sig_tx.dtype, len(sig_rx), sig_rx.dtype, off)) + logging.debug( + "Sig1_cut: %d %s, Sig2_cut: %d %s, off: %d" % (len(sig_tx), sig_tx.dtype, len(sig_rx), sig_rx.dtype, off)) return sig_tx, sig_rx - def fromfile(self, 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) - - # The MIT License (MIT) # # Copyright (c) 2017 Andreas Steger diff --git a/dpd/src/ExtractStatistic.py b/dpd/src/ExtractStatistic.py index 9df85bc..8ea849b 100644 --- a/dpd/src/ExtractStatistic.py +++ b/dpd/src/ExtractStatistic.py @@ -1,13 +1,12 @@ # -*- coding: utf-8 -*- # # DPD Calculation Engine, -# Extract statistic from data to use in Model +# Extract statistic from received TX and RX data to use in Model # # http://www.opendigitalradio.org # Licence: The MIT License, see notice at the end of this file import numpy as np -import pickle import matplotlib.pyplot as plt import datetime @@ -30,6 +29,14 @@ def _check_input_extract(tx_dpd, rx_received): assert normalization_error < 0.01, "Non normalized signals" +def _phase_diff_value_per_bin(phase_diffs_values_lists): + phase_list = [] + for values in phase_diffs_values_lists: + mean = np.mean(values) if len(values) > 0 else np.nan + phase_list.append(mean) + return phase_list + + class ExtractStatistic: """Calculate a low variance RX value for equally spaced tx values of a predefined range""" @@ -37,31 +44,27 @@ class ExtractStatistic: def __init__(self, c): self.c = c + # Number of measurements used to extract the statistic self.n_meas = 0 + # Boundaries for the bins self.tx_boundaries = np.linspace(c.ES_start, c.ES_end, c.ES_n_bins + 1) self.n_per_bin = c.ES_n_per_bin + # List of rx values for each bin self.rx_values_lists = [] for i in range(c.ES_n_bins): self.rx_values_lists.append([]) + # List of tx values for each bin self.tx_values_lists = [] for i in range(c.ES_n_bins): self.tx_values_lists.append([]) - self.tx_values = self._tx_value_per_bin() - - self.rx_values = [] - for i in range(c.ES_n_bins): - self.rx_values.append(None) - self.plot = c.ES_plot - def _plot_and_log(self): + def _plot_and_log(self, tx_values, rx_values, phase_diffs_values, phase_diffs_values_lists): if logging.getLogger().getEffectiveLevel() == logging.DEBUG and self.plot: - phase_diffs_values_lists = self._phase_diff_list_per_bin() - phase_diffs_values = self._phase_diff_value_per_bin(phase_diffs_values_lists) dt = datetime.datetime.now().isoformat() fig_path = logging_path + "/" + dt + "_ExtractStatistic.png" @@ -72,13 +75,13 @@ class ExtractStatistic: i_sub += 1 ax = plt.subplot(sub_rows, sub_cols, i_sub) - ax.plot(self.tx_values, self.rx_values, + ax.plot(tx_values, rx_values, label="Estimated Values", color="red") - for i, tx_value in enumerate(self.tx_values): - rx_values = self.rx_values_lists[i] - ax.scatter(np.ones(len(rx_values)) * tx_value, - np.abs(rx_values), + for i, tx_value in enumerate(tx_values): + rx_values_list = self.rx_values_lists[i] + ax.scatter(np.ones(len(rx_values_list)) * tx_value, + np.abs(rx_values_list), s=0.1, color="black") ax.set_title("Extracted Statistic") @@ -90,10 +93,10 @@ class ExtractStatistic: i_sub += 1 ax = plt.subplot(sub_rows, sub_cols, i_sub) - ax.plot(self.tx_values, np.rad2deg(phase_diffs_values), + ax.plot(tx_values, np.rad2deg(phase_diffs_values), label="Estimated Values", color="red") - for i, tx_value in enumerate(self.tx_values): + for i, tx_value in enumerate(tx_values): phase_diff = phase_diffs_values_lists[i] ax.scatter(np.ones(len(phase_diff)) * tx_value, np.rad2deg(phase_diff), @@ -101,14 +104,14 @@ class ExtractStatistic: color="black") ax.set_xlabel("TX Amplitude") ax.set_ylabel("Phase Difference") - ax.set_ylim(-60,60) + ax.set_ylim(-60, 60) ax.set_xlim(0, 1.1) ax.legend(loc=4) num = [] - for i, tx_value in enumerate(self.tx_values): - rx_values = self.rx_values_lists[i] - num.append(len(rx_values)) + for i, tx_value in enumerate(tx_values): + rx_values_list = self.rx_values_lists[i] + num.append(len(rx_values_list)) i_sub += 1 ax = plt.subplot(sub_rows, sub_cols, i_sub) ax.plot(num) @@ -120,9 +123,6 @@ class ExtractStatistic: fig.savefig(fig_path) plt.close(fig) - pickle.dump(self.rx_values_lists, open("/tmp/rx_values", "wb")) - pickle.dump(self.tx_values, open("/tmp/tx_values", "wb")) - def _rx_value_per_bin(self): rx_values = [] for values in self.rx_values_lists: @@ -145,14 +145,9 @@ class ExtractStatistic: phase_values_lists.append(phase_diffs) return phase_values_lists - def _phase_diff_value_per_bin(self, phase_diffs_values_lists): - phase_list = [] - for values in phase_diffs_values_lists: - mean = np.mean(values) if len(values) > 0 else np.nan - phase_list.append(mean) - return phase_list - def extract(self, tx_dpd, rx): + """Extract information from a new measurement and store them + in member variables.""" _check_input_extract(tx_dpd, rx) self.n_meas += 1 @@ -165,23 +160,22 @@ class ExtractStatistic: self.tx_values_lists[i] += \ list(tx_dpd[mask][:n_add]) - self.rx_values = self._rx_value_per_bin() - self.tx_values = self._tx_value_per_bin() - - self._plot_and_log() + rx_values = self._rx_value_per_bin() + tx_values = self._tx_value_per_bin() n_per_bin = np.array([len(values) for values in self.rx_values_lists]) # Index of first not filled bin, assumes that never all bins are filled idx_end = np.argmin(n_per_bin == self.c.ES_n_per_bin) - # TODO cleanup phase_diffs_values_lists = self._phase_diff_list_per_bin() - phase_diffs_values = self._phase_diff_value_per_bin(phase_diffs_values_lists) + phase_diffs_values = _phase_diff_value_per_bin(phase_diffs_values_lists) + + self._plot_and_log(tx_values, rx_values, phase_diffs_values, phase_diffs_values_lists) - return np.array(self.tx_values, dtype=np.float32)[:idx_end], \ - np.array(self.rx_values, dtype=np.float32)[:idx_end], \ - np.array(phase_diffs_values, dtype=np.float32)[:idx_end], \ - n_per_bin + tx_values_crop = np.array(tx_values, dtype=np.float32)[:idx_end] + rx_values_crop = np.array(rx_values, dtype=np.float32)[:idx_end] + phase_diffs_values_crop = np.array(phase_diffs_values, dtype=np.float32)[:idx_end] + return tx_values_crop, rx_values_crop, phase_diffs_values_crop, n_per_bin # The MIT License (MIT) # diff --git a/dpd/src/Measure_Shoulders.py b/dpd/src/Measure_Shoulders.py index 710e800..7d48a2b 100644 --- a/dpd/src/Measure_Shoulders.py +++ b/dpd/src/Measure_Shoulders.py @@ -101,6 +101,10 @@ class Measure_Shoulders: plt.close(fig) def average_shoulders(self, signal, n_avg=None): + if not self.c.MS_enable: + logging.debug("Shoulder Measurement disabled") + return None + assert signal.shape[0] > 4 * self.c.MS_FFT_size if n_avg is None: n_avg = self.c.MS_averaging_size diff --git a/dpd/src/Model_AM.py b/dpd/src/Model_AM.py index 85f6495..cdc3de1 100644 --- a/dpd/src/Model_AM.py +++ b/dpd/src/Model_AM.py @@ -15,14 +15,30 @@ import numpy as np import matplotlib.pyplot as plt +def is_npfloat32(array): + assert isinstance(array, np.ndarray), type(array) + assert array.dtype == np.float32, array.dtype + assert array.flags.contiguous + assert not any(np.isnan(array)) + + def check_input_get_next_coefs(tx_dpd, rx_received): - is_float32 = lambda x: (isinstance(x, np.ndarray) and - x.dtype == np.float32 and - x.flags.contiguous) - assert is_float32(tx_dpd), \ - "tx_dpd is not float32 but {}".format(tx_dpd[0].dtype) - assert is_float32(rx_received), \ - "rx_received is not float32 but {}".format(tx_dpd[0].dtype) + is_npfloat32(tx_dpd) + is_npfloat32(rx_received) + + +def poly(sig): + return np.array([sig ** i for i in range(1, 6)]).T + + +def fit_poly(tx_abs, rx_abs): + return np.linalg.lstsq(poly(rx_abs), tx_abs)[0] + + +def calc_line(coefs, min_amp, max_amp): + rx_range = np.linspace(min_amp, max_amp) + tx_est = np.sum(poly(rx_range) * coefs, axis=1) + return tx_est, rx_range class Model_AM: @@ -40,8 +56,8 @@ class Model_AM: def _plot(self, tx_dpd, rx_received, coefs_am, coefs_am_new): if logging.getLogger().getEffectiveLevel() == logging.DEBUG and self.plot: - tx_range, rx_est = self.calc_line(coefs_am, 0, 0.6) - tx_range_new, rx_est_new = self.calc_line(coefs_am_new, 0, 0.6) + tx_range, rx_est = calc_line(coefs_am, 0, 0.6) + tx_range_new, rx_est_new = calc_line(coefs_am_new, 0, 0.6) dt = datetime.datetime.now().isoformat() fig_path = logging_path + "/" + dt + "_Model_AM.svg" @@ -66,7 +82,6 @@ class Model_AM: ax.set_title("Model_AM") ax.set_xlabel("TX Amplitude") ax.set_ylabel("RX Amplitude") - xlim = ax.get_xlim() ax.set_xlim(-0.5, 1.5) ax.legend(loc=4) @@ -74,21 +89,10 @@ class Model_AM: fig.savefig(fig_path) plt.close(fig) - def poly(self, sig): - return np.array([sig ** i for i in range(1, 6)]).T - - def fit_poly(self, tx_abs, rx_abs): - return np.linalg.lstsq(self.poly(rx_abs), tx_abs)[0] - - def calc_line(self, coefs, min_amp, max_amp): - rx_range = np.linspace(min_amp, max_amp) - tx_est = np.sum(self.poly(rx_range) * coefs, axis=1) - return tx_est, rx_range - def get_next_coefs(self, tx_dpd, rx_received, coefs_am): check_input_get_next_coefs(tx_dpd, rx_received) - coefs_am_new = self.fit_poly(tx_dpd, rx_received) + coefs_am_new = fit_poly(tx_dpd, rx_received) coefs_am_new = coefs_am + \ self.learning_rate_am * (coefs_am_new - coefs_am) diff --git a/dpd/src/Agc.py b/dpd/src/RX_Agc.py index 670fbbb..670fbbb 100644 --- a/dpd/src/Agc.py +++ b/dpd/src/RX_Agc.py diff --git a/dpd/src/phase_align.py b/dpd/src/phase_align.py index 7317d70..0d662ec 100644 --- a/dpd/src/phase_align.py +++ b/dpd/src/phase_align.py @@ -10,8 +10,6 @@ import logging logging_path = os.path.dirname(logging.getLoggerClass().root.handlers[0].baseFilename) import numpy as np -from scipy import signal, optimize -import sys import matplotlib.pyplot as plt diff --git a/dpd/src/test_dab_Util.py b/dpd/src/test_dab_Util.py index ec15586..86ee2d8 100644 --- a/dpd/src/test_dab_Util.py +++ b/dpd/src/test_dab_Util.py @@ -16,7 +16,7 @@ class TestDab_Util(TestCase): def test_subsample_align(self, sample_orig=r'../test_data/orig_rough_aligned.dat', sample_rec =r'../test_data/recored_rough_aligned.dat', length = 10240, max_size = 1000000): - du = DU.Dab_Util(8196000) + du = DU res1 = [] res2 = [] for i in range(10): |