From f4ca82137e850e30d31e7008b34800d8b2699e5d Mon Sep 17 00:00:00 2001 From: "Matthias P. Braendli" Date: Wed, 19 Dec 2018 16:11:58 +0100 Subject: DPD: Merge Model_PM and _AM into _Poly --- python/dpd/ExtractStatistic.py | 8 ++- python/dpd/GlobalConfig.py | 9 +-- python/dpd/Model_AM.py | 119 --------------------------------- python/dpd/Model_PM.py | 121 ---------------------------------- python/dpd/Model_Poly.py | 146 +++++++++++++++++++++++++++++++++++------ python/dpd/RX_Agc.py | 18 ++--- 6 files changed, 145 insertions(+), 276 deletions(-) delete mode 100644 python/dpd/Model_AM.py delete mode 100644 python/dpd/Model_PM.py (limited to 'python/dpd') diff --git a/python/dpd/ExtractStatistic.py b/python/dpd/ExtractStatistic.py index 1aa4391..a23fa1a 100644 --- a/python/dpd/ExtractStatistic.py +++ b/python/dpd/ExtractStatistic.py @@ -38,7 +38,7 @@ class ExtractStatistic: """Calculate a low variance RX value for equally spaced tx values of a predefined range""" - def __init__(self, c): + def __init__(self, c, peak_amplitude): self.c = c self._plot_data = None @@ -47,7 +47,7 @@ class ExtractStatistic: 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.tx_boundaries = np.linspace(0.0, peak_amplitude, c.ES_n_bins + 1) self.n_per_bin = c.ES_n_per_bin # List of rx values for each bin @@ -60,6 +60,10 @@ class ExtractStatistic: for i in range(c.ES_n_bins): self.tx_values_lists.append([]) + def get_bin_info(self): + return "Binning: {} bins used for amplitudes between {} and {}".format( + len(self.tx_boundaries), np.min(self.tx_boundaries), np.max(self.tx_boundaries)) + def plot(self, plot_path, title): if self._plot_data is not None: tx_values, rx_values, phase_diffs_values, phase_diffs_values_lists = self._plot_data diff --git a/python/dpd/GlobalConfig.py b/python/dpd/GlobalConfig.py index 99280f2..632a63b 100644 --- a/python/dpd/GlobalConfig.py +++ b/python/dpd/GlobalConfig.py @@ -26,6 +26,8 @@ class GlobalConfig: self.T_U = oversample * 2048 # Inverse of carrier spacing self.T_C = oversample * 504 # Duration of cyclic prefix + self.median_to_peak = 12 # Estimated value for a DAB OFDM signal + # Frequency Domain # example: np.delete(fft[3328:4865], 768) self.FFT_delta = 1536 # Number of carrier frequencies @@ -40,10 +42,8 @@ class GlobalConfig: self.phase_offset_per_sample = 1. / self.sample_rate * 2 * np.pi * 1000 # Constants for ExtractStatistic - self.ES_plot = plot - self.ES_start = 0.0 self.ES_end = 1.0 - self.ES_n_bins = 64 # Number of bins between ES_start and ES_end + self.ES_n_bins = 64 self.ES_n_per_bin = 128 # Number of measurements pre bin # Constants for Measure_Shoulder @@ -68,9 +68,6 @@ class GlobalConfig: # Constants for MER self.MER_plot = plot - # Constants for Model - self.MDL_plot = plot - # Constants for Model_PM # Set all phase offsets to zero for TX amplitude < MPM_tx_min self.MPM_tx_min = 0.1 diff --git a/python/dpd/Model_AM.py b/python/dpd/Model_AM.py deleted file mode 100644 index b07a5a5..0000000 --- a/python/dpd/Model_AM.py +++ /dev/null @@ -1,119 +0,0 @@ -# -*- coding: utf-8 -*- -# -# DPD Computation Engine, model implementation for Amplitude and not Phase -# -# http://www.opendigitalradio.org -# Licence: The MIT License, see notice at the end of this file - -import datetime -import os -import logging -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_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, rcond=None)[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: - """Calculates new coefficients using the measurement and the previous - coefficients""" - - def __init__(self, c, learning_rate_am=1): - self.c = c - self.learning_rate_am = learning_rate_am - self._plot_data = None - - def plot(self, plot_location, title): - if self._plot_data is not None: - tx_dpd, rx_received, coefs_am, coefs_am_new = self._plot_data - - 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) - - sub_rows = 1 - sub_cols = 1 - fig = plt.figure(figsize=(sub_cols * 6, sub_rows / 2. * 6)) - i_sub = 0 - - i_sub += 1 - ax = plt.subplot(sub_rows, sub_cols, i_sub) - ax.plot(tx_range, rx_est, - label="Estimated TX", - alpha=0.3, - color="gray") - ax.plot(tx_range_new, rx_est_new, - label="New Estimated TX", - color="red") - ax.scatter(tx_dpd, rx_received, - label="Binned Data", - color="blue", - s=1) - ax.set_title("Model_AM {}".format(title)) - ax.set_xlabel("TX Amplitude") - ax.set_ylabel("RX Amplitude") - ax.set_xlim(-0.5, 1.5) - ax.legend(loc=4) - - fig.tight_layout() - fig.savefig(plot_location) - plt.close(fig) - - def get_next_coefs(self, tx_dpd, rx_received, coefs_am): - """Calculate the next AM/AM coefficients using the extracted - statistic of TX and RX amplitude""" - check_input_get_next_coefs(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) - - self._plot_data = (tx_dpd, rx_received, coefs_am, coefs_am_new) - - return coefs_am_new - -# The MIT License (MIT) -# -# Copyright (c) 2017 Andreas Steger -# Copyright (c) 2018 Matthias P. Braendli -# -# Permission is hereby granted, free of charge, to any person obtaining a copy -# of this software and associated documentation files (the "Software"), to deal -# in the Software without restriction, including without limitation the rights -# to use, copy, modify, merge, publish, distribute, sublicense, and/or sell -# copies of the Software, and to permit persons to whom the Software is -# furnished to do so, subject to the following conditions: -# -# The above copyright notice and this permission notice shall be included in all -# copies or substantial portions of the Software. -# -# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR -# IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, -# FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE -# AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER -# LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, -# OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE -# SOFTWARE. diff --git a/python/dpd/Model_PM.py b/python/dpd/Model_PM.py deleted file mode 100644 index 40fa1d4..0000000 --- a/python/dpd/Model_PM.py +++ /dev/null @@ -1,121 +0,0 @@ -# -*- coding: utf-8 -*- -# -# DPD Computation Engine, model implementation for Amplitude and not Phase -# -# http://www.opendigitalradio.org -# Licence: The MIT License, see notice at the end of this file - -import datetime -import os -import logging -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, phase_diff): - is_npfloat32(tx_dpd) - is_npfloat32(phase_diff) - - -class Model_PM: - """Calculates new coefficients using the measurement and the previous - coefficients""" - - def __init__(self, c, learning_rate_pm=1): - self.c = c - self.learning_rate_pm = learning_rate_pm - self._plot_data = None - - def plot(self, plot_location, title): - if self._plot_data is not None: - tx_dpd, phase_diff, coefs_pm, coefs_pm_new = self._plot_data - - tx_range, phase_diff_est = self.calc_line(coefs_pm, 0, 0.6) - tx_range_new, phase_diff_est_new = self.calc_line(coefs_pm_new, 0, 0.6) - - sub_rows = 1 - sub_cols = 1 - fig = plt.figure(figsize=(sub_cols * 6, sub_rows / 2. * 6)) - i_sub = 0 - - i_sub += 1 - ax = plt.subplot(sub_rows, sub_cols, i_sub) - ax.plot(tx_range, phase_diff_est, - label="Estimated Phase Diff", - alpha=0.3, - color="gray") - ax.plot(tx_range_new, phase_diff_est_new, - label="New Estimated Phase Diff", - color="red") - ax.scatter(tx_dpd, phase_diff, - label="Binned Data", - color="blue", - s=1) - ax.set_title("Model_PM {}".format(title)) - ax.set_xlabel("TX Amplitude") - ax.set_ylabel("Phase DIff") - ax.legend(loc=4) - - fig.tight_layout() - fig.savefig(plot_location) - plt.close(fig) - - def _discard_small_values(self, tx_dpd, phase_diff): - """ Assumes that the phase for small tx amplitudes is zero""" - mask = tx_dpd < self.c.MPM_tx_min - phase_diff[mask] = 0 - return tx_dpd, phase_diff - - def poly(self, sig): - return np.array([sig ** i for i in range(0, 5)]).T - - def fit_poly(self, tx_abs, phase_diff): - return np.linalg.lstsq(self.poly(tx_abs), phase_diff, rcond=None)[0] - - def calc_line(self, coefs, min_amp, max_amp): - tx_range = np.linspace(min_amp, max_amp) - phase_diff = np.sum(self.poly(tx_range) * coefs, axis=1) - return tx_range, phase_diff - - def get_next_coefs(self, tx_dpd, phase_diff, coefs_pm): - """Calculate the next AM/PM coefficients using the extracted - statistic of TX amplitude and phase difference""" - tx_dpd, phase_diff = self._discard_small_values(tx_dpd, phase_diff) - check_input_get_next_coefs(tx_dpd, phase_diff) - - coefs_pm_new = self.fit_poly(tx_dpd, phase_diff) - - coefs_pm_new = coefs_pm + self.learning_rate_pm * (coefs_pm_new - coefs_pm) - self._plot_data = (tx_dpd, phase_diff, coefs_pm, coefs_pm_new) - - return coefs_pm_new - -# The MIT License (MIT) -# -# Copyright (c) 2017 Andreas Steger -# Copyright (c) 2018 Matthias P. Braendli -# -# Permission is hereby granted, free of charge, to any person obtaining a copy -# of this software and associated documentation files (the "Software"), to deal -# in the Software without restriction, including without limitation the rights -# to use, copy, modify, merge, publish, distribute, sublicense, and/or sell -# copies of the Software, and to permit persons to whom the Software is -# furnished to do so, subject to the following conditions: -# -# The above copyright notice and this permission notice shall be included in all -# copies or substantial portions of the Software. -# -# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR -# IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, -# FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE -# AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER -# LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, -# OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE -# SOFTWARE. diff --git a/python/dpd/Model_Poly.py b/python/dpd/Model_Poly.py index ca39492..5722531 100644 --- a/python/dpd/Model_Poly.py +++ b/python/dpd/Model_Poly.py @@ -8,15 +8,13 @@ import os import logging import numpy as np +import matplotlib.pyplot as plt -import dpd.Model_AM as Model_AM -import dpd.Model_PM as Model_PM - - -def assert_np_float32(x): - assert isinstance(x, np.ndarray) - assert x.dtype == np.float32 - assert x.flags.contiguous +def assert_np_float32(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_abs, rx_abs, phase_diff): @@ -44,12 +42,73 @@ class Poly: self.reset_coefs() - self.model_am = Model_AM.Model_AM(c) - self.model_pm = Model_PM.Model_PM(c) - def plot(self, am_plot_location, pm_plot_location, title): - self.model_am.plot(am_plot_location, title) - self.model_pm.plot(pm_plot_location, title) + if self._am_plot_data is not None: + tx_dpd, rx_received, coefs_am, coefs_am_new = self._am_plot_data + + tx_range, rx_est = self._am_calc_line(coefs_am, 0, 0.6) + tx_range_new, rx_est_new = self._am_calc_line(coefs_am_new, 0, 0.6) + + sub_rows = 1 + sub_cols = 1 + fig = plt.figure(figsize=(sub_cols * 6, sub_rows / 2. * 6)) + i_sub = 0 + + i_sub += 1 + ax = plt.subplot(sub_rows, sub_cols, i_sub) + ax.plot(tx_range, rx_est, + label="Estimated TX", + alpha=0.3, + color="gray") + ax.plot(tx_range_new, rx_est_new, + label="New Estimated TX", + color="red") + ax.scatter(tx_dpd, rx_received, + label="Binned Data", + color="blue", + s=1) + ax.set_title("Model AM {}".format(title)) + ax.set_xlabel("TX Amplitude") + ax.set_ylabel("RX Amplitude") + ax.set_xlim(-0.5, 1.5) + ax.legend(loc=4) + + fig.tight_layout() + fig.savefig(am_plot_location) + plt.close(fig) + + if self._pm_plot_data is not None: + tx_dpd, phase_diff, coefs_pm, coefs_pm_new = self._pm_plot_data + + tx_range, phase_diff_est = self._pm_calc_line(coefs_pm, 0, 0.6) + tx_range_new, phase_diff_est_new = self._pm_calc_line(coefs_pm_new, 0, 0.6) + + sub_rows = 1 + sub_cols = 1 + fig = plt.figure(figsize=(sub_cols * 6, sub_rows / 2. * 6)) + i_sub = 0 + + i_sub += 1 + ax = plt.subplot(sub_rows, sub_cols, i_sub) + ax.plot(tx_range, phase_diff_est, + label="Estimated Phase Diff", + alpha=0.3, + color="gray") + ax.plot(tx_range_new, phase_diff_est_new, + label="New Estimated Phase Diff", + color="red") + ax.scatter(tx_dpd, phase_diff, + label="Binned Data", + color="blue", + s=1) + ax.set_title("Model PM {}".format(title)) + ax.set_xlabel("TX Amplitude") + ax.set_ylabel("Phase DIff") + ax.legend(loc=4) + + fig.tight_layout() + fig.savefig(pm_plot_location) + plt.close(fig) def reset_coefs(self): self.coefs_am = np.zeros(5, dtype=np.float32) @@ -65,12 +124,8 @@ class Poly: """ _check_input_get_next_coefs(tx_abs, rx_abs, phase_diff) - if not lr is None: - self.model_am.learning_rate_am = lr - self.model_pm.learning_rate_pm = lr - - coefs_am_new = self.model_am.get_next_coefs(tx_abs, rx_abs, self.coefs_am) - coefs_pm_new = self.model_pm.get_next_coefs(tx_abs, phase_diff, self.coefs_pm) + coefs_am_new = self._am_get_next_coefs(tx_abs, rx_abs, self.coefs_am) + coefs_pm_new = self._pm_get_next_coefs(tx_abs, phase_diff, self.coefs_pm) self.coefs_am = self.coefs_am + (coefs_am_new - self.coefs_am) * self.learning_rate_am self.coefs_pm = self.coefs_pm + (coefs_pm_new - self.coefs_pm) * self.learning_rate_pm @@ -78,9 +133,62 @@ class Poly: def get_dpd_data(self): return "poly", self.coefs_am, self.coefs_pm + def _am_calc_line(self, coefs, min_amp, max_amp): + rx_range = np.linspace(min_amp, max_amp) + tx_est = np.sum(self._am_poly(rx_range) * coefs, axis=1) + return tx_est, rx_range + + def _am_poly(self, sig): + return np.array([sig ** i for i in range(1, 6)]).T + + def _am_fit_poly(self, tx_abs, rx_abs): + return np.linalg.lstsq(self._am_poly(rx_abs), tx_abs, rcond=None)[0] + + def _am_get_next_coefs(self, tx_dpd, rx_received, coefs_am): + """Calculate the next AM/AM coefficients using the extracted + statistic of TX and RX amplitude""" + + coefs_am_new = self._am_fit_poly(tx_dpd, rx_received) + coefs_am_new = coefs_am + \ + self.learning_rate_am * (coefs_am_new - coefs_am) + + self._am_plot_data = (tx_dpd, rx_received, coefs_am, coefs_am_new) + + return coefs_am_new + + def _pm_poly(self, sig): + return np.array([sig ** i for i in range(0, 5)]).T + + def _pm_calc_line(self, coefs, min_amp, max_amp): + tx_range = np.linspace(min_amp, max_amp) + phase_diff = np.sum(self._pm_poly(tx_range) * coefs, axis=1) + return tx_range, phase_diff + + def _discard_small_values(self, tx_dpd, phase_diff): + """ Assumes that the phase for small tx amplitudes is zero""" + mask = tx_dpd < self.c.MPM_tx_min + phase_diff[mask] = 0 + return tx_dpd, phase_diff + + def _pm_fit_poly(self, tx_abs, phase_diff): + return np.linalg.lstsq(self._pm_poly(tx_abs), phase_diff, rcond=None)[0] + + def _pm_get_next_coefs(self, tx_dpd, phase_diff, coefs_pm): + """Calculate the next AM/PM coefficients using the extracted + statistic of TX amplitude and phase difference""" + tx_dpd, phase_diff = self._discard_small_values(tx_dpd, phase_diff) + + coefs_pm_new = self._pm_fit_poly(tx_dpd, phase_diff) + + coefs_pm_new = coefs_pm + self.learning_rate_pm * (coefs_pm_new - coefs_pm) + self._pm_plot_data = (tx_dpd, phase_diff, coefs_pm, coefs_pm_new) + + return coefs_pm_new + # The MIT License (MIT) # # Copyright (c) 2017 Andreas Steger +# Copyright (c) 2018 Matthias P. Brandli # # Permission is hereby granted, free of charge, to any person obtaining a copy # of this software and associated documentation files (the "Software"), to deal diff --git a/python/dpd/RX_Agc.py b/python/dpd/RX_Agc.py index 4700e68..911f3c9 100644 --- a/python/dpd/RX_Agc.py +++ b/python/dpd/RX_Agc.py @@ -19,19 +19,19 @@ import dpd.Adapt as Adapt import dpd.Measure as Measure class Agc: - """The goal of the automatic gain control is to set the - RX gain to a value at which all received amplitudes can - be detected. This means that the maximum possible amplitude + """The goal of the automatic gain control is to set the + RX gain to a value at which all received amplitudes can + be detected. This means that the maximum possible amplitude should be quantized at the highest possible digital value. - A problem we have to face, is that the estimation of the - maximum amplitude by applying the max() function is very - unstable. This is due to the maximum’s rareness. Therefore - we estimate a far more robust value, such as the median, + A problem we have to face, is that the estimation of the + maximum amplitude by applying the max() function is very + unstable. This is due to the maximum’s rareness. Therefore + we estimate a far more robust value, such as the median, and then approximate the maximum amplitude from it. - Given this, we tune the RX gain in such a way, that the - received signal fulfills our desired property, of having + Given this, we tune the RX gain in such a way, that the + received signal fulfills our desired property, of having all amplitudes properly quantized.""" def __init__(self, measure, adapt, c): -- cgit v1.2.3