aboutsummaryrefslogtreecommitdiffstats
path: root/dpd
diff options
context:
space:
mode:
Diffstat (limited to 'dpd')
-rwxr-xr-xdpd/main.py81
-rw-r--r--dpd/src/Adapt.py30
-rw-r--r--dpd/src/Agc.py6
-rw-r--r--dpd/src/Dab_Util.py25
-rw-r--r--dpd/src/MER.py5
-rw-r--r--dpd/src/Model.py481
-rw-r--r--dpd/src/TX_Agc.py109
-rw-r--r--dpd/src/phase_align.py8
-rwxr-xr-xdpd/src/subsample_align.py6
9 files changed, 465 insertions, 286 deletions
diff --git a/dpd/main.py b/dpd/main.py
index cfd8126..5e67c90 100755
--- a/dpd/main.py
+++ b/dpd/main.py
@@ -12,10 +12,15 @@ predistortion module of ODR-DabMod."""
import datetime
import os
+import time
+
+import matplotlib
+matplotlib.use('GTKAgg')
import logging
+
dt = datetime.datetime.now().isoformat()
-logging_path = "/tmp/dpd_{}".format(dt).replace(".","_").replace(":","-")
+logging_path = "/tmp/dpd_{}".format(dt).replace(".", "_").replace(":", "-")
os.makedirs(logging_path)
logging.basicConfig(format='%(asctime)s - %(module)s - %(levelname)s - %(message)s',
datefmt='%Y-%m-%d %H:%M:%S',
@@ -33,27 +38,33 @@ console.setFormatter(formatter)
# add the handler to the root logger
logging.getLogger('').addHandler(console)
+import numpy as np
import traceback
import src.Measure as Measure
import src.Model as Model
import src.Adapt as Adapt
import src.Agc as Agc
+import src.TX_Agc as TX_Agc
+import src.Symbol_align
+import src.const
+import src.MER
import argparse
-parser = argparse.ArgumentParser(description="DPD Computation Engine for ODR-DabMod")
+parser = argparse.ArgumentParser(
+ description="DPD Computation Engine for ODR-DabMod")
parser.add_argument('--port', default='50055',
- help='port of DPD server to connect to (default: 50055)',
- required=False)
+ help='port of DPD server to connect to (default: 50055)',
+ required=False)
parser.add_argument('--rc-port', default='9400',
- help='port of ODR-DabMod ZMQ Remote Control to connect to (default: 9400)',
- required=False)
+ help='port of ODR-DabMod ZMQ Remote Control to connect to (default: 9400)',
+ required=False)
parser.add_argument('--samplerate', default='8192000',
- help='Sample rate',
- required=False)
+ help='Sample rate',
+ required=False)
parser.add_argument('--coefs', default='poly.coef',
- help='File with DPD coefficients, which will be read by ODR-DabMod',
- required=False)
-parser.add_argument('--txgain', default=65,
+ help='File with DPD coefficients, which will be read by ODR-DabMod',
+ required=False)
+parser.add_argument('--txgain', default=71,
help='TX Gain',
required=False,
type=int)
@@ -61,47 +72,62 @@ parser.add_argument('--rxgain', default=30,
help='TX Gain',
required=False,
type=int)
-parser.add_argument('--samps', default='20480',
- help='Number of samples to request from ODR-DabMod',
- required=False)
+parser.add_argument('--digital_gain', default=1,
+ help='Digital Gain',
+ required=False,
+ type=float)
+parser.add_argument('--samps', default='81920',
+ help='Number of samples to request from ODR-DabMod',
+ required=False)
parser.add_argument('-i', '--iterations', default='1',
- help='Number of iterations to run',
- required=False)
+ help='Number of iterations to run',
+ required=False)
parser.add_argument('-l', '--load-poly',
- help='Load existing polynomial',
- action="store_true")
+ help='Load existing polynomial',
+ action="store_true")
cli_args = parser.parse_args()
port = int(cli_args.port)
port_rc = int(cli_args.rc_port)
coef_path = cli_args.coefs
+digital_gain = cli_args.digital_gain
txgain = cli_args.txgain
rxgain = cli_args.rxgain
num_req = int(cli_args.samps)
samplerate = int(cli_args.samplerate)
num_iter = int(cli_args.iterations)
+SA = src.Symbol_align.Symbol_align(samplerate)
+MER = src.MER.MER(samplerate)
+c = src.const.const(samplerate)
+
meas = Measure.Measure(samplerate, port, num_req)
adapt = Adapt.Adapt(port_rc, coef_path)
coefs_am, coefs_pm = adapt.get_coefs()
if cli_args.load_poly:
- model = Model.Model(coefs_am, coefs_pm)
+ model = Model.Model(c, SA, MER, coefs_am, coefs_pm, plot=True)
else:
- model = Model.Model([1, 0, 0, 0, 0], [0, 0, 0, 0, 0])
+ model = Model.Model(c, SA, MER, [1.0, 0, 0, 0, 0], [0, 0, 0, 0, 0], plot=True)
+adapt.set_coefs(model.coefs_am, model.coefs_pm)
+adapt.set_digital_gain(digital_gain)
adapt.set_txgain(txgain)
adapt.set_rxgain(rxgain)
tx_gain = adapt.get_txgain()
rx_gain = adapt.get_rxgain()
+digital_gain = adapt.get_digital_gain()
dpd_coefs_am, dpd_coefs_pm = adapt.get_coefs()
logging.info(
- "TX gain {}, RX gain {}, dpd_coefs_am {}, dpd_coefs_pm {}".format(
- tx_gain, rx_gain, dpd_coefs_am, dpd_coefs_pm
+ "TX gain {}, RX gain {}, dpd_coefs_am {},"
+ " dpd_coefs_pm {}, digital_gain {}".format(
+ tx_gain, rx_gain, dpd_coefs_am, dpd_coefs_pm, digital_gain
)
)
+tx_agc = TX_Agc.TX_Agc(adapt)
+
# Automatic Gain Control
agc = Agc.Agc(meas, adapt)
agc.run()
@@ -110,8 +136,19 @@ for i in range(num_iter):
try:
txframe_aligned, tx_ts, rxframe_aligned, rx_ts, rx_median = meas.get_samples()
logging.debug("tx_ts {}, rx_ts {}".format(tx_ts, rx_ts))
+ assert tx_ts - rx_ts < 1e-5, "Time stamps do not match."
+
+ if tx_agc.adapt_if_necessary(txframe_aligned):
+ continue
+
coefs_am, coefs_pm = model.get_next_coefs(txframe_aligned, rxframe_aligned)
adapt.set_coefs(coefs_am, coefs_pm)
+
+ off = SA.calc_offset(txframe_aligned)
+ tx_mer = MER.calc_mer(txframe_aligned[off:off + c.T_U])
+ rx_mer = MER.calc_mer(rxframe_aligned[off:off + c.T_U], debug=True)
+ logging.info("MER with lag in it. {}: TX {}, RX {}".
+ format(i, tx_mer, rx_mer))
except Exception as e:
logging.warning("Iteration {} failed.".format(i))
logging.warning(traceback.format_exc())
diff --git a/dpd/src/Adapt.py b/dpd/src/Adapt.py
index 2fb596f..b4042d6 100644
--- a/dpd/src/Adapt.py
+++ b/dpd/src/Adapt.py
@@ -24,7 +24,7 @@ class Adapt:
"""
def __init__(self, port, coef_path):
- logging.info("Instantiate Adapt object")
+ logging.debug("Instantiate Adapt object")
self.port = port
self.coef_path = coef_path
self.host = "localhost"
@@ -60,7 +60,7 @@ class Adapt:
The message string that will be sent to the receiver.
"""
sock = self._connect()
- logging.info("Send message: %s" % message)
+ logging.debug("Send message: %s" % message)
msg_parts = message.split(" ")
for i, part in enumerate(msg_parts):
if i == len(msg_parts) - 1:
@@ -71,7 +71,7 @@ class Adapt:
sock.send(part.encode(), flags=f)
data = [el.decode() for el in sock.recv_multipart()]
- logging.info("Received message: %s" % message)
+ logging.debug("Received message: %s" % message)
return data
def set_txgain(self, gain):
@@ -85,12 +85,12 @@ class Adapt:
# TODO this is specific to the B200
if gain < 0 or gain > 89:
raise ValueError("Gain has to be in [0,89]")
- return self.send_receive("set uhd txgain %d" % gain)
+ return self.send_receive("set uhd txgain %.4f" % float(gain))
def get_txgain(self):
"""Get the txgain value in dB for the ODR-DabMod."""
# TODO handle failure
- return self.send_receive("get uhd txgain")
+ return float(self.send_receive("get uhd txgain")[0])
def set_rxgain(self, gain):
"""Set a new rxgain for the ODR-DabMod.
@@ -103,12 +103,28 @@ class Adapt:
# TODO this is specific to the B200
if gain < 0 or gain > 89:
raise ValueError("Gain has to be in [0,89]")
- return self.send_receive("set uhd rxgain %d" % gain)
+ return self.send_receive("set uhd rxgain %.4f" % float(gain))
def get_rxgain(self):
"""Get the rxgain value in dB for the ODR-DabMod."""
# TODO handle failure
- return self.send_receive("get uhd rxgain")
+ return float(self.send_receive("get uhd rxgain")[0])
+
+ def set_digital_gain(self, gain):
+ """Set a new rxgain for the ODR-DabMod.
+
+ Parameters
+ ----------
+ gain : int
+ new RX gain, in the same format as ODR-DabMod's config file
+ """
+ msg = "set gain digital %.5f" % gain
+ return self.send_receive(msg)
+
+ def get_digital_gain(self):
+ """Get the rxgain value in dB for the ODR-DabMod."""
+ # TODO handle failure
+ return float(self.send_receive("get gain digital")[0])
def _read_coef_file(self, path):
"""Load the coefficients from the file in the format given in the README,
diff --git a/dpd/src/Agc.py b/dpd/src/Agc.py
index 1fd11c8..978b607 100644
--- a/dpd/src/Agc.py
+++ b/dpd/src/Agc.py
@@ -47,7 +47,7 @@ class Agc:
def run(self):
self.adapt.set_rxgain(self.rxgain)
- for i in range(3):
+ for i in range(2):
# Measure
txframe_aligned, tx_ts, rxframe_aligned, rx_ts, rx_median= \
self.measure.get_samples()
@@ -65,7 +65,7 @@ class Agc:
))
self.adapt.set_rxgain(self.rxgain)
- time.sleep(1)
+ time.sleep(0.5)
def plot_estimates(self):
"""Plots the estimate of for Max, Median, Mean for different
@@ -74,7 +74,7 @@ class Agc:
time.sleep(1)
dt = datetime.datetime.now().isoformat()
- fig_path = logging_path + "/" + dt + "_agc.pdf"
+ fig_path = logging_path + "/" + dt + "_agc.svg"
fig, axs = plt.subplots(2, 2, figsize=(3*6,1*6))
axs = axs.ravel()
diff --git a/dpd/src/Dab_Util.py b/dpd/src/Dab_Util.py
index b0b3ce3..e3dbfe3 100644
--- a/dpd/src/Dab_Util.py
+++ b/dpd/src/Dab_Util.py
@@ -23,7 +23,7 @@ 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):
+ def __init__(self, sample_rate, plot=False):
"""
:param sample_rate: sample rate [sample/sec] to use for calculations
"""
@@ -31,6 +31,8 @@ class Dab_Util:
self.dab_bandwidth = 1536000 #Bandwidth of a dab signal
self.frame_ms = 96 #Duration of a Dab frame
+ self.plot=plot
+
def lag(self, sig_orig, sig_rec):
"""
Find lag between two signals
@@ -41,9 +43,9 @@ class Dab_Util:
off = sig_rec.shape[0]
c = np.abs(signal.correlate(sig_orig, sig_rec))
- if logging.getLogger().getEffectiveLevel() == logging.DEBUG:
+ if logging.getLogger().getEffectiveLevel() == logging.DEBUG and self.plot:
dt = datetime.datetime.now().isoformat()
- corr_path = (logging_path + "/" + dt + "_tx_rx_corr.pdf")
+ corr_path = (logging_path + "/" + dt + "_tx_rx_corr.svg")
plt.plot(c, label="corr")
plt.legend()
plt.savefig(corr_path)
@@ -95,9 +97,9 @@ class Dab_Util:
Returns an aligned version of sig_tx and sig_rx by cropping and subsample alignment
"""
- if logging.getLogger().getEffectiveLevel() == logging.DEBUG:
+ if logging.getLogger().getEffectiveLevel() == logging.DEBUG and self.plot:
dt = datetime.datetime.now().isoformat()
- fig_path = logging_path + "/" + dt + "_sync_raw.pdf"
+ fig_path = logging_path + "/" + dt + "_sync_raw.svg"
fig, axs = plt.subplots(2)
axs[0].plot(np.abs(sig_tx[:128]), label="TX Frame")
@@ -139,9 +141,9 @@ class Dab_Util:
sig_tx = sig_tx[:-1]
sig_rx = sig_rx[:-1]
- if logging.getLogger().getEffectiveLevel() == logging.DEBUG:
+ if logging.getLogger().getEffectiveLevel() == logging.DEBUG and self.plot:
dt = datetime.datetime.now().isoformat()
- fig_path = logging_path + "/" + dt + "_sync_sample_aligned.pdf"
+ fig_path = logging_path + "/" + dt + "_sync_sample_aligned.svg"
fig, axs = plt.subplots(2)
axs[0].plot(np.abs(sig_tx[:128]), label="TX Frame")
@@ -161,12 +163,11 @@ class Dab_Util:
fig.savefig(fig_path)
fig.clf()
-
sig_rx = sa.subsample_align(sig_rx, sig_tx)
- if logging.getLogger().getEffectiveLevel() == logging.DEBUG:
+ if logging.getLogger().getEffectiveLevel() == logging.DEBUG and self.plot:
dt = datetime.datetime.now().isoformat()
- fig_path = logging_path + "/" + dt + "_sync_subsample_aligned.pdf"
+ fig_path = logging_path + "/" + dt + "_sync_subsample_aligned.svg"
fig, axs = plt.subplots(2)
axs[0].plot(np.abs(sig_tx[:128]), label="TX Frame")
@@ -188,9 +189,9 @@ class Dab_Util:
sig_rx = pa.phase_align(sig_rx, sig_tx)
- if logging.getLogger().getEffectiveLevel() == logging.DEBUG:
+ if logging.getLogger().getEffectiveLevel() == logging.DEBUG and self.plot:
dt = datetime.datetime.now().isoformat()
- fig_path = logging_path + "/" + dt + "_sync_phase_aligned.pdf"
+ fig_path = logging_path + "/" + dt + "_sync_phase_aligned.svg"
fig, axs = plt.subplots(2)
axs[0].plot(np.abs(sig_tx[:128]), label="TX Frame")
diff --git a/dpd/src/MER.py b/dpd/src/MER.py
index 54ef6c6..00fcc23 100644
--- a/dpd/src/MER.py
+++ b/dpd/src/MER.py
@@ -75,8 +75,9 @@ class MER:
spectrum = self._calc_spectrum(tx)
- dt = datetime.datetime.now().isoformat()
- fig_path = logging_path + "/" + dt + "_MER.pdf"
+ if debug:
+ dt = datetime.datetime.now().isoformat()
+ fig_path = logging_path + "/" + dt + "_MER.svg"
MERs = []
for i, (x, y) in enumerate(self._split_in_carrier(
diff --git a/dpd/src/Model.py b/dpd/src/Model.py
index 2aa9feb..827027a 100644
--- a/dpd/src/Model.py
+++ b/dpd/src/Model.py
@@ -8,262 +8,221 @@
import datetime
import os
import logging
+
logging_path = os.path.dirname(logging.getLoggerClass().root.handlers[0].baseFilename)
import numpy as np
-import matplotlib
-matplotlib.use('agg')
import matplotlib.pyplot as plt
-from sklearn.linear_model import Ridge
+from sklearn import linear_model
class Model:
"""Calculates new coefficients using the measurement and the old
coefficients"""
- def __init__(self, coefs_am, coefs_pm):
+ def __init__(self,
+ c,
+ SA,
+ MER,
+ coefs_am,
+ coefs_pm,
+ learning_rate_am=1.,
+ learning_rate_pm=1.,
+ plot=False):
+ self.c = c
+ self.SA = SA
+ self.MER = MER
+
+ self.learning_rate_am = learning_rate_am
+ self.learning_rate_pm = learning_rate_pm
+
self.coefs_am = coefs_am
- self.coefs_history = [coefs_am, ]
- self.mses = [0, ]
- self.errs = [0, ]
+ self.coefs_am_history = [coefs_am, ]
+ self.mses_am = []
+ self.errs_am = []
+
+ self.tx_mers = []
+ self.rx_mers = []
self.coefs_pm = coefs_pm
self.coefs_pm_history = [coefs_pm, ]
- self.errs_phase = [0, ]
+ self.errs_pm = []
+
+ self.plot = plot
- def sample_uniformly(self, txframe_aligned, rxframe_aligned, n_bins=4):
+ def sample_uniformly(self, tx_dpd, rx_received, n_bins=5):
"""This function returns tx and rx samples in a way
that the tx amplitudes have an approximate uniform
- distribution with respect to the txframe_aligned amplitudes"""
- txframe_aligned_abs = np.abs(txframe_aligned)
+ distribution with respect to the tx_dpd amplitudes"""
+ mask = np.logical_and((np.abs(tx_dpd) > 0.01), (np.abs(rx_received) > 0.01))
+ tx_dpd = tx_dpd[mask]
+ rx_received = rx_received[mask]
+
+ txframe_aligned_abs = np.abs(tx_dpd)
ccdf_min = 0
ccdf_max = np.max(txframe_aligned_abs)
tx_hist, ccdf_edges = np.histogram(txframe_aligned_abs,
bins=n_bins,
range=(ccdf_min, ccdf_max))
- n_choise = np.min(tx_hist)
+ n_choise = np.min(tx_hist)
tx_choice = np.zeros(n_choise * n_bins, dtype=np.complex64)
rx_choice = np.zeros(n_choise * n_bins, dtype=np.complex64)
for idx, bin in enumerate(tx_hist):
indices = np.where((txframe_aligned_abs >= ccdf_edges[idx]) &
- (txframe_aligned_abs <= ccdf_edges[idx+1]))[0]
- indices_choise = np.random.choice(indices, n_choise, replace=False)
- rx_choice[idx*n_choise:(idx+1)*n_choise] = rxframe_aligned[indices_choise]
- tx_choice[idx*n_choise:(idx+1)*n_choise] = txframe_aligned[indices_choise]
+ (txframe_aligned_abs <= ccdf_edges[idx + 1]))[0]
+ indices_choise = np.random.choice(indices,
+ n_choise,
+ replace=False)
+ rx_choice[idx * n_choise:(idx + 1) * n_choise] = \
+ rx_received[indices_choise]
+ tx_choice[idx * n_choise:(idx + 1) * n_choise] = \
+ tx_dpd[indices_choise]
+
+ assert isinstance(rx_choice[0], np.complex64), \
+ "rx_choice is not complex64 but {}".format(rx_choice[0].dtype)
+ assert isinstance(tx_choice[0], np.complex64), \
+ "tx_choice is not complex64 but {}".format(tx_choice[0].dtype)
+
return tx_choice, rx_choice
- def get_next_coefs(self, txframe_aligned, rxframe_aligned):
- tx_choice, rx_choice = self.sample_uniformly(txframe_aligned, rxframe_aligned)
-
- # Calculate new coefficients for AM/AM correction
- rx_abs = np.abs(rx_choice)
- rx_A = np.vstack([rx_abs,
- rx_abs ** 3,
- rx_abs ** 5,
- rx_abs ** 7,
- rx_abs ** 9,
- ]).T
- rx_dpd = np.sum(rx_A * self.coefs_am, axis=1)
- rx_dpd = rx_dpd * (
- np.median(np.abs(tx_choice)) / np.median(np.abs(rx_dpd)))
-
- err = rx_dpd - np.abs(tx_choice)
- self.errs.append(np.mean(np.abs(err ** 2)))
-
- a_delta = np.linalg.lstsq(rx_A, err)[0]
- new_coefs = self.coefs_am - 0.1 * a_delta
- new_coefs = new_coefs * (self.coefs_am[0] / new_coefs[0])
- logging.debug("a_delta {}".format(a_delta))
- logging.debug("new coefs_am {}".format(new_coefs))
-
- # Calculate new coefficients for AM/PM correction
- phase_diff_rad = ((
- (np.angle(tx_choice) -
- np.angle(rx_choice) +
- np.pi) % (2 * np.pi)) -
- np.pi
- )
-
- tx_abs = np.abs(tx_choice)
- tx_abs_A = np.vstack([tx_abs,
- tx_abs ** 2,
- tx_abs ** 3,
- tx_abs ** 4,
- tx_abs ** 5,
+ def dpd_amplitude(self, sig, coefs=None):
+ if coefs is None:
+ coefs = self.coefs_am
+ assert isinstance(sig[0], np.complex64), "Sig is not complex64 but {}".format(sig[0].dtype)
+ sig_abs = np.abs(sig)
+ A_sig = np.vstack([np.ones(sig_abs.shape),
+ sig_abs ** 1,
+ sig_abs ** 2,
+ sig_abs ** 3,
+ sig_abs ** 4,
+ ]).T
+ sig_dpd = sig * np.sum(A_sig * coefs, axis=1)
+ return sig_dpd, A_sig
+
+ def dpd_phase(self, sig, coefs=None):
+ if coefs is None:
+ coefs = self.coefs_pm
+ assert isinstance(sig[0], np.complex64), "Sig is not complex64 but {}".format(sig[0].dtype)
+ sig_abs = np.abs(sig)
+ A_phase = np.vstack([np.ones(sig_abs.shape),
+ sig_abs ** 1,
+ sig_abs ** 2,
+ sig_abs ** 3,
+ sig_abs ** 4,
]).T
- phase_dpd = np.sum(tx_abs_A * self.coefs_pm, axis=1)
-
- err_phase = phase_dpd - phase_diff_rad
- self.errs_phase.append(np.mean(np.abs(err_phase ** 2)))
- a_delta = np.linalg.lstsq(tx_abs_A, err_phase)[0]
- new_coefs_pm = self.coefs_pm - 0.1 * a_delta
- logging.debug("a_delta {}".format(a_delta))
- logging.debug("new new_coefs_pm {}".format(new_coefs_pm))
-
- def dpd_phase(tx):
- tx_abs = np.abs(tx)
- tx_A_complex = np.vstack([tx,
- tx * tx_abs ** 1,
- tx * tx_abs ** 2,
- tx * tx_abs ** 3,
- tx * tx_abs ** 4,
- ]).T
- tx_dpd = np.sum(tx_A_complex * self.coefs_pm, axis=1)
- return tx_dpd
-
- tx_range = np.linspace(0, 2)
- phase_range_dpd = dpd_phase(tx_range)
-
- rx_A_complex = np.vstack([rx_choice,
- rx_choice * rx_abs ** 2,
- rx_choice * rx_abs ** 4,
- rx_choice * rx_abs ** 6,
- rx_choice * rx_abs ** 8,
- ]).T
- rx_post_distored = np.sum(rx_A_complex * self.coefs_am, axis=1)
- rx_post_distored = rx_post_distored * (
- np.median(np.abs(tx_choice)) /
- np.median(np.abs(rx_post_distored)))
- mse = np.mean(np.abs((tx_choice - rx_post_distored) ** 2))
- logging.debug("MSE: {}".format(mse))
- self.mses.append(mse)
-
- def dpd(tx):
- tx_abs = np.abs(tx)
- tx_A_complex = np.vstack([tx,
- tx * tx_abs ** 2,
- tx * tx_abs ** 4,
- tx * tx_abs ** 6,
- tx * tx_abs ** 8,
- ]).T
- tx_dpd = np.sum(tx_A_complex * self.coefs_am, axis=1)
- return tx_dpd
-
- rx_range = np.linspace(0, 1, num=100)
- rx_range_dpd = dpd(rx_range)
- rx_range = rx_range[(rx_range_dpd > 0) & (rx_range_dpd < 2)]
- rx_range_dpd = rx_range_dpd[(rx_range_dpd > 0) & (rx_range_dpd < 2)]
-
- if logging.getLogger().getEffectiveLevel() == logging.DEBUG:
- logging.debug("txframe: min %f, max %f, median %f" %
- (np.min(np.abs(txframe_aligned)),
- np.max(np.abs(txframe_aligned)),
- np.median(np.abs(txframe_aligned))
- ))
-
- logging.debug("rxframe: min %f, max %f, median %f" %
- (np.min(np.abs(rx_choice)),
- np.max(np.abs(rx_choice)),
- np.median(np.abs(rx_choice))
- ))
+ phase_diff_est = np.sum(A_phase * coefs, axis=1)
+ return phase_diff_est, A_phase
+ def _next_am_coefficent(self, tx_choice, rx_choice):
+ """Calculate new coefficients for AM/AM correction"""
+ rx_dpd, rx_A = self.dpd_amplitude(rx_choice)
+ rx_dpd = rx_dpd * (
+ np.median(np.abs(tx_choice)) /
+ np.median(np.abs(rx_dpd)))
+ err = np.abs(rx_dpd) - np.abs(tx_choice)
+ mse = np.mean(np.abs((rx_dpd - tx_choice) ** 2))
+ self.mses_am.append(mse)
+ self.errs_am.append(np.mean(err**2))
+
+ reg = linear_model.Ridge(alpha=0.00001)
+ reg.fit(rx_A, err)
+ a_delta = reg.coef_
+ new_coefs_am = self.coefs_am - self.learning_rate_am * a_delta
+ new_coefs_am = new_coefs_am * (self.coefs_am[0] / new_coefs_am[0])
+ return new_coefs_am
+
+ def _next_pm_coefficent(self, tx_choice, rx_choice):
+ """Calculate new coefficients for AM/PM correction
+ Assuming deviations smaller than pi/2"""
+ phase_diff_choice = np.angle(
+ (rx_choice * tx_choice.conjugate()) /
+ (np.abs(rx_choice) * np.abs(tx_choice))
+ )
+ plt.hist(phase_diff_choice)
+ plt.savefig('/tmp/hist_' + str(np.random.randint(0,1000)) + '.svg')
+ plt.clf()
+ phase_diff_est, phase_A = self.dpd_phase(rx_choice)
+ err_phase = phase_diff_est - phase_diff_choice
+ self.errs_pm.append(np.mean(np.abs(err_phase ** 2)))
+
+ reg = linear_model.Ridge(alpha=0.00001)
+ reg.fit(phase_A, err_phase)
+ p_delta = reg.coef_
+ new_coefs_pm = self.coefs_pm - self.learning_rate_pm * p_delta
+
+ return new_coefs_pm, phase_diff_choice
+
+ def get_next_coefs(self, tx_dpd, rx_received):
+ # Check data type
+ assert tx_dpd[0].dtype == np.complex64, \
+ "tx_dpd is not complex64 but {}".format(tx_dpd[0].dtype)
+ assert rx_received[0].dtype == np.complex64, \
+ "rx_received is not complex64 but {}".format(rx_received[0].dtype)
+ # Check if signals have same normalization
+ normalization_error = np.abs(np.median(np.abs(tx_dpd)) -
+ np.median(np.abs(rx_received))) / (
+ np.median(np.abs(tx_dpd)) + np.median(np.abs(rx_received)))
+ assert normalization_error < 0.01, "Non normalized signals"
+
+ tx_choice, rx_choice = self.sample_uniformly(tx_dpd, rx_received)
+ new_coefs_am = self._next_am_coefficent(tx_choice, rx_choice)
+ new_coefs_pm, phase_diff_choice = self._next_pm_coefficent(tx_choice, rx_choice)
+
+ logging.debug('txframe: min {:.2f}, max {:.2f}, ' \
+ 'median {:.2f}; rxframe: min {:.2f}, max {:.2f}, ' \
+ 'median {:.2f}; new coefs_am {};' \
+ 'new_coefs_pm {}'.format(
+ np.min(np.abs(tx_dpd)),
+ np.max(np.abs(tx_dpd)),
+ np.median(np.abs(tx_dpd)),
+ np.min(np.abs(rx_choice)),
+ np.max(np.abs(rx_choice)),
+ np.median(np.abs(rx_choice)),
+ new_coefs_am,
+ new_coefs_pm))
+
+ if logging.getLogger().getEffectiveLevel() == logging.DEBUG and self.plot:
+ off = self.SA.calc_offset(tx_dpd)
+ tx_mer = self.MER.calc_mer(tx_dpd[off:off + self.c.T_U])
+ rx_mer = self.MER.calc_mer(rx_received[off:off + self.c.T_U], debug=True)
+ self.tx_mers.append(tx_mer)
+ self.rx_mers.append(rx_mer)
+
+ if logging.getLogger().getEffectiveLevel() == logging.DEBUG and self.plot:
dt = datetime.datetime.now().isoformat()
- fig_path = logging_path + "/" + dt + "_Model.pdf"
+ fig_path = logging_path + "/" + dt + "_Model.svg"
- fig = plt.figure(figsize=(3*6, 1.5 * 6))
+ fig = plt.figure(figsize=(2 * 6, 2 * 6))
- ax = plt.subplot(3,3,1)
- ax.plot(np.abs(txframe_aligned[:128]),
- label="TX sent",
- linestyle=":")
- ax.plot(np.abs(rxframe_aligned[:128]),
- label="RX received",
- color="red")
- ax.set_title("Synchronized Signals of Iteration {}".format(len(self.coefs_history)))
- ax.set_xlabel("Samples")
- ax.set_ylabel("Amplitude")
- ax.text(0, 0, "TX (max {:01.3f}, mean {:01.3f}, median {:01.3f})".format(
- np.max(np.abs(txframe_aligned)),
- np.mean(np.abs(txframe_aligned)),
- np.median(np.abs(txframe_aligned))
- ), size = 8)
- ax.legend(loc=4)
+ i_sub = 1
- ax = plt.subplot(3,3,2)
- ax.plot(np.real(txframe_aligned[:128]),
+ ax = plt.subplot(4, 2, i_sub)
+ i_sub += 1
+ ax.plot(np.abs(tx_dpd[:128]),
label="TX sent",
linestyle=":")
- ax.plot(np.real(rxframe_aligned[:128]),
+ ax.plot(np.abs(rx_received[:128]),
label="RX received",
color="red")
- ax.set_title("Synchronized Signals")
- ax.set_xlabel("Samples")
- ax.set_ylabel("Real Part")
- ax.legend(loc=4)
-
- ax = plt.subplot(3,3,3)
- ax.plot(np.abs(txframe_aligned[:128]),
- label="TX Frame",
- linestyle=":",
- linewidth=0.5)
- ax.plot(np.abs(rxframe_aligned[:128]),
- label="RX Frame",
- linestyle="--",
- linewidth=0.5)
-
- rx_abs = np.abs(rxframe_aligned)
- rx_A = np.vstack([rx_abs,
- rx_abs ** 3,
- rx_abs ** 5,
- rx_abs ** 7,
- rx_abs ** 9,
- ]).T
- rx_dpd = np.sum(rx_A * self.coefs_am, axis=1)
- rx_dpd = rx_dpd * (
- np.median(np.abs(tx_choice)) / np.median(np.abs(rx_dpd)))
-
- ax.plot(np.abs(rx_dpd[:128]),
- label="RX DPD Frame",
- linestyle="-.",
- linewidth=0.5)
-
- tx_abs = np.abs(np.abs(txframe_aligned[:128]))
- tx_A = np.vstack([tx_abs,
- tx_abs ** 3,
- tx_abs ** 5,
- tx_abs ** 7,
- tx_abs ** 9,
- ]).T
- tx_dpd = np.sum(tx_A * new_coefs, axis=1)
- tx_dpd_norm = tx_dpd * (
- np.median(np.abs(tx_choice)) / np.median(np.abs(tx_dpd)))
-
- ax.plot(np.abs(tx_dpd_norm[:128]),
- label="TX DPD Frame Norm",
- linestyle="-.",
- linewidth=0.5)
- ax.legend(loc=4)
- ax.set_title("RX DPD")
+ ax.set_title("Synchronized Signals of Iteration {}"
+ .format(len(self.coefs_am_history)))
ax.set_xlabel("Samples")
ax.set_ylabel("Amplitude")
+ ax.text(0, 0, "TX (max {:01.3f}, mean {:01.3f}, " \
+ "median {:01.3f})".format(
+ np.max(np.abs(tx_dpd)),
+ np.mean(np.abs(tx_dpd)),
+ np.median(np.abs(tx_dpd))
+ ), size=8)
+ ax.legend(loc=4)
- ax = plt.subplot(3,3,4)
- ax.scatter(
- np.abs(tx_choice[:1024]),
- np.abs(rx_choice[:1024]),
- s=0.1)
- ax.plot(rx_range_dpd / self.coefs_am[0], rx_range, linewidth=0.25)
- ax.set_title("Amplifier Characteristic")
- ax.set_xlabel("TX Amplitude")
- ax.set_ylabel("RX Amplitude")
-
- ax = plt.subplot(3,3,5)
- ax.scatter(
- np.abs(tx_choice[:1024]),
- phase_diff_rad[:1024] * 180 / np.pi,
- s=0.1
- )
- ax.plot(tx_range, phase_range_dpd * 180 / np.pi, linewidth=0.25)
- ax.set_title("Amplifier Characteristic")
- ax.set_xlabel("TX Amplitude")
- ax.set_ylabel("Phase Difference [deg]")
-
- ax = plt.subplot(3,3,6)
+ ax = plt.subplot(4, 2, i_sub)
+ i_sub += 1
ccdf_min, ccdf_max = 0, 1
- tx_hist, ccdf_edges = np.histogram(np.abs(txframe_aligned),
- bins=60,
- range=(ccdf_min, ccdf_max))
- tx_hist_normalized = tx_hist.astype(float)/np.sum(tx_hist)
+ tx_hist, ccdf_edges = np.histogram(np.abs(tx_dpd),
+ bins=60,
+ range=(ccdf_min, ccdf_max))
+ tx_hist_normalized = tx_hist.astype(float) / np.sum(tx_hist)
ccdf = 1.0 - np.cumsum(tx_hist_normalized)
ax.semilogy(ccdf_edges[:-1], ccdf, label="CCDF")
ax.semilogy(ccdf_edges[:-1],
@@ -271,14 +230,51 @@ class Model:
label="Histogram",
drawstyle='steps')
ax.legend(loc=4)
- ax.set_ylim(1e-5,2)
+ ax.set_ylim(1e-5, 2)
ax.set_title("Complementary Cumulative Distribution Function")
ax.set_xlabel("TX Amplitude")
ax.set_ylabel("Ratio of Samples larger than x")
- ax = plt.subplot(3,3,7)
- coefs_history = np.array(self.coefs_history)
- for idx, coef_hist in enumerate(coefs_history.T):
+ ax = plt.subplot(4, 2, i_sub)
+ i_sub += 1
+ ax.semilogy(np.array(self.mses_am) + 1e-10, label="log(MSE)")
+ ax.semilogy(np.array(self.errs_am) + 1e-10, label="log(ERR)")
+ ax.legend(loc=4)
+ ax.set_title("MSE History")
+ ax.set_xlabel("Iterations")
+ ax.set_ylabel("MSE")
+
+ ax = plt.subplot(4, 2, i_sub)
+ i_sub += 1
+ ax.plot(self.tx_mers, label="TX MER")
+ ax.plot(self.rx_mers, label="RX MER")
+ ax.legend(loc=4)
+ ax.set_title("MER History")
+ ax.set_xlabel("Iterations")
+ ax.set_ylabel("MER")
+
+ ax = plt.subplot(4, 2, i_sub)
+ rx_range = np.linspace(0, 1, num=100, dtype=np.complex64)
+ rx_range_dpd = self.dpd_amplitude(rx_range)[0]
+ rx_range_dpd_new = self.dpd_amplitude(rx_range, new_coefs_am)[0]
+ i_sub += 1
+ ax.scatter(
+ np.abs(tx_choice),
+ np.abs(rx_choice),
+ s=0.1)
+ ax.plot(rx_range_dpd / self.coefs_am[0], rx_range, linewidth=0.25, label="current")
+ ax.plot(rx_range_dpd_new / self.coefs_am[0], rx_range, linewidth=0.25, label="next")
+ ax.set_ylim(0, 1)
+ ax.set_xlim(0, 1)
+ ax.legend()
+ ax.set_title("Amplifier Characteristic")
+ ax.set_xlabel("TX Amplitude")
+ ax.set_ylabel("RX Amplitude")
+
+ ax = plt.subplot(4, 2, i_sub)
+ i_sub += 1
+ coefs_am_history = np.array(self.coefs_am_history)
+ for idx, coef_hist in enumerate(coefs_am_history.T):
ax.plot(coef_hist,
label="Coef {}".format(idx),
linewidth=0.5)
@@ -287,10 +283,38 @@ class Model:
ax.set_xlabel("Iterations")
ax.set_ylabel("Coefficient Value")
- ax = plt.subplot(3,3,8)
- coefs_history = np.array(self.coefs_pm_history)
- for idx, coef_hist in enumerate(coefs_history.T):
- ax.plot(coef_hist,
+ phase_range = np.linspace(0, 1, num=100, dtype=np.complex64)
+ phase_range_dpd = self.dpd_phase(phase_range)[0]
+ phase_range_dpd_new = self.dpd_phase(phase_range,
+ coefs=new_coefs_pm)[0]
+ ax = plt.subplot(4, 2, i_sub)
+ i_sub += 1
+ ax.scatter(
+ np.abs(tx_choice),
+ np.rad2deg(phase_diff_choice),
+ s=0.1)
+ ax.plot(
+ np.abs(phase_range),
+ np.rad2deg(phase_range_dpd),
+ linewidth=0.25,
+ label="current")
+ ax.plot(
+ np.abs(phase_range),
+ np.rad2deg(phase_range_dpd_new),
+ linewidth=0.25,
+ label="next")
+ ax.set_ylim(-60, 60)
+ ax.set_xlim(0, 1)
+ ax.legend()
+ ax.set_title("Amplifier Characteristic")
+ ax.set_xlabel("TX Amplitude")
+ ax.set_ylabel("Phase Difference")
+
+ ax = plt.subplot(4, 2, i_sub)
+ i_sub += 1
+ coefs_pm_history = np.array(self.coefs_pm_history)
+ for idx, coef_phase_hist in enumerate(coefs_pm_history.T):
+ ax.plot(coef_phase_hist,
label="Coef {}".format(idx),
linewidth=0.5)
ax.legend(loc=4)
@@ -298,21 +322,12 @@ class Model:
ax.set_xlabel("Iterations")
ax.set_ylabel("Coefficient Value")
- ax = plt.subplot(3,3,9)
- coefs_history = np.array(self.coefs_history)
- ax.plot(self.mses, label="MSE")
- ax.plot(self.errs, label="ERR")
- ax.legend(loc=4)
- ax.set_title("MSE History")
- ax.set_xlabel("Iterations")
- ax.set_ylabel("MSE")
-
fig.tight_layout()
fig.savefig(fig_path)
fig.clf()
- self.coefs_am = new_coefs
- self.coefs_history.append(self.coefs_am)
+ self.coefs_am = new_coefs_am
+ self.coefs_am_history.append(self.coefs_am)
self.coefs_pm = new_coefs_pm
self.coefs_pm_history.append(self.coefs_pm)
return self.coefs_am, self.coefs_pm
diff --git a/dpd/src/TX_Agc.py b/dpd/src/TX_Agc.py
new file mode 100644
index 0000000..9274612
--- /dev/null
+++ b/dpd/src/TX_Agc.py
@@ -0,0 +1,109 @@
+# -*- coding: utf-8 -*-
+#
+# Automatic Gain Control
+#
+# http://www.opendigitalradio.org
+# Licence: The MIT License, see notice at the end of this file
+
+import datetime
+import os
+import logging
+import time
+
+logging_path = os.path.dirname(logging.getLoggerClass().root.handlers[0].baseFilename)
+
+import numpy as np
+import matplotlib
+
+matplotlib.use('agg')
+import matplotlib.pyplot as plt
+
+import src.Adapt as Adapt
+
+#TODO fix for float tx_gain
+class TX_Agc:
+ def __init__(self,
+ adapt,
+ max_txgain=89,
+ tx_median_target=0.1,
+ tx_median_threshold_max=0.12,
+ tx_median_threshold_min=0.08):
+ """
+ In order to avoid digital clipping, this class increases the
+ TX gain and reduces the digital gain. Digital clipping happens
+ when the digital analog converter receives values greater than
+ it's maximal output. This class solves that problem by adapting
+ the TX gain in a way that the peaks of the TX signal are in a
+ specified range. The TX gain is adapted accordingly. The TX peaks
+ are approximated by estimating it based on the signal median.
+
+ :param adapt: Instance of Adapt Class to update
+ txgain and coefficients
+ :param max_txgain: limit for TX gain
+ :param tx_median_threshold_max: if the median of TX is larger
+ than this value, then the digital gain is reduced
+ :param tx_median_threshold_min: if the median of TX is smaller
+ than this value, then the digital gain is increased
+ :param tx_median_target: The digital gain is reduced in a way that
+ the median TX value is expected to be lower than this value.
+ """
+ assert isinstance(adapt, Adapt.Adapt)
+ self.adapt = adapt
+ self.max_txgain = max_txgain
+ self.txgain = self.max_txgain
+
+ assert tx_median_threshold_max > tx_median_target,\
+ "The tolerated tx_median has to be larger then the goal tx_median"
+ self.tx_median_threshold_tolerate_max = tx_median_threshold_max
+ self.tx_median_threshold_tolerate_min = tx_median_threshold_min
+ self.tx_median_target = tx_median_target
+
+ def adapt_if_necessary(self, tx):
+ tx_median = np.median(np.abs(tx))
+ if tx_median > self.tx_median_threshold_tolerate_max or\
+ tx_median < self.tx_median_threshold_tolerate_min:
+ delta_db = 20 * np.log10(self.tx_median_target / tx_median)
+ new_txgain = self.adapt.get_txgain() - delta_db
+ assert new_txgain < self.max_txgain,\
+ "TX_Agc failed. New TX gain of {} is too large.".format(
+ new_txgain
+ )
+ self.adapt.set_txgain(new_txgain)
+ txgain = self.adapt.get_txgain()
+
+ digital_gain_factor = 10 ** (delta_db / 20.)
+ digital_gain = self.adapt.get_digital_gain() * digital_gain_factor
+ self.adapt.set_digital_gain(digital_gain)
+
+ logging.info(
+ "digital_gain = {}, txgain_new = {}, "\
+ "delta_db = {}, tx_median {}, "\
+ "digital_gain_factor = {}".
+ format(digital_gain, txgain, delta_db,
+ tx_median, digital_gain_factor))
+
+ time.sleep(1)
+ return True
+ return False
+
+# The MIT License (MIT)
+#
+# Copyright (c) 2017 Andreas Steger
+#
+# 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/dpd/src/phase_align.py b/dpd/src/phase_align.py
index 7045586..7f82392 100644
--- a/dpd/src/phase_align.py
+++ b/dpd/src/phase_align.py
@@ -15,7 +15,7 @@ import sys
import matplotlib.pyplot as plt
-def phase_align(sig, ref_sig):
+def phase_align(sig, ref_sig, plot=False):
"""Do phase alignment for sig relative to the reference signal
ref_sig.
@@ -26,9 +26,9 @@ def phase_align(sig, ref_sig):
real_diffs = np.cos(angle_diff)
imag_diffs = np.sin(angle_diff)
- if logging.getLogger().getEffectiveLevel() == logging.DEBUG:
+ if logging.getLogger().getEffectiveLevel() == logging.DEBUG and plot:
dt = datetime.datetime.now().isoformat()
- fig_path = logging_path + "/" + dt + "_phase_align.pdf"
+ fig_path = logging_path + "/" + dt + "_phase_align.svg"
plt.subplot(511)
plt.hist(angle_diff, bins=60, label="Angle Diff")
@@ -66,7 +66,7 @@ def phase_align(sig, ref_sig):
))
sig = sig * np.exp(1j * -angle)
- if logging.getLogger().getEffectiveLevel() == logging.DEBUG:
+ if logging.getLogger().getEffectiveLevel() == logging.DEBUG and plot:
plt.subplot(515)
plt.plot(np.angle(ref_sig[:128]), label="ref_sig")
plt.plot(np.angle(sig[:128]), label="sig")
diff --git a/dpd/src/subsample_align.py b/dpd/src/subsample_align.py
index d0f3dba..6d1cd2a 100755
--- a/dpd/src/subsample_align.py
+++ b/dpd/src/subsample_align.py
@@ -29,7 +29,7 @@ def gen_omega(length):
return omega
-def subsample_align(sig, ref_sig):
+def subsample_align(sig, ref_sig, plot=False):
"""Do subsample alignment for sig relative to the reference signal
ref_sig. The delay between the two must be less than sample
@@ -68,13 +68,13 @@ def subsample_align(sig, ref_sig):
if optim_result.success:
best_tau = optim_result.x
- if 1:
+ if plot:
corr = np.vectorize(correlate_for_delay)
ixs = np.linspace(-1, 1, 100)
taus = corr(ixs)
dt = datetime.datetime.now().isoformat()
- tau_path = (logging_path + "/" + dt + "_tau.pdf")
+ tau_path = (logging_path + "/" + dt + "_tau.svg")
plt.plot(ixs, taus)
plt.title("Subsample correlation, minimum is best: {}".format(best_tau))
plt.savefig(tau_path)