diff options
Diffstat (limited to 'tii/analyse_tii_from_file.py')
| -rwxr-xr-x | tii/analyse_tii_from_file.py | 453 |
1 files changed, 453 insertions, 0 deletions
diff --git a/tii/analyse_tii_from_file.py b/tii/analyse_tii_from_file.py new file mode 100755 index 0000000..9d834ea --- /dev/null +++ b/tii/analyse_tii_from_file.py @@ -0,0 +1,453 @@ +#!/usr/bin/env python +# -*- coding: utf-8 -*- +# +# Analyse TII of iq file. complexf 2048000 +# +# It requires SciPy and matplotlib. +# +# http://www.opendigitalradio.org +# Licence: The MIT License, see notice at the end of this file + +import sys +import functools +import numpy as np +import matplotlib.pyplot as pp +import argparse +import subprocess +from scipy.misc import imsave +from tii import calculate_carrier_indices, calculate_reduced_carrier_indices + +all_tii_patterns = {} +cp_per_carrier = {} + +def prepare_tii_patterns(): + print("Building TII table...") + global all_tii_patterns + global cp_per_carrier + + for c in range(24): + for p in range(70): + carriers = calculate_carrier_indices(c, p) + all_tii_patterns[tuple(carriers)] = (c, p) + + carriers = calculate_reduced_carrier_indices(c, p) + for k in carriers: + if k not in cp_per_carrier: + cp_per_carrier[k] = set() + cp_per_carrier[k].add((c, p)) + +SIZEOF_SAMPLE = 8 # complex floats + +# Constants for TM 1 +NbSymbols = 76 +NbCarriers = 1536 +Spacing = 2048 +NullSize = 2656 +SymSize = 2552 +FicSizeOut = 288 + +CyclicPrefix = SymSize - Spacing + +prs_phases = np.array([0, -1, 1, -1, -1, -1, -1, 0, 0, 1, -1, -1, -1, 1, 1, 0, 0, -1, 1, + -1, -1, -1, -1, 0, 0, 1, -1, -1, -1, 1, 1, 0, 0, 1, 2, -1, 2, 1, 0, 0, + -1, -1, 0, -1, 2, -1, 2, 0, -1, 1, 2, -1, 2, 1, 0, 0, -1, -1, 0, -1, 2, + -1, 2, 0, -1, 1, 1, 1, -1, 1, -1, 2, 0, -1, -1, 1, -1, -1, 1, 2, 0, 1, + 1, 1, -1, 1, -1, 2, 0, -1, -1, 1, -1, -1, 1, 2, 0, 1, 0, -1, 0, 1, 2, + 0, 1, -1, 2, -1, 0, -1, 0, 0, 1, 1, 0, -1, 0, 1, 2, 0, 1, -1, 2, -1, 0, + -1, 0, 0, 1, 2, 0, 2, 2, 2, 2, -1, -1, 0, 2, 2, 2, 0, 0, -1, -1, 2, 0, + 2, 2, 2, 2, -1, -1, 0, 2, 2, 2, 0, 0, -1, -1, 2, -1, 0, -1, 2, 1, 1, 0, + 0, 1, 0, -1, 0, -1, 1, 0, 2, -1, 0, -1, 2, 1, 1, 0, 0, 1, 0, -1, 0, -1, + 1, 0, 1, 1, 1, -1, 1, -1, 2, 0, -1, -1, 1, -1, -1, 1, 2, 0, 1, 1, 1, + -1, 1, -1, 2, 0, -1, -1, 1, -1, -1, 1, 2, 0, 0, -1, 2, -1, 0, 1, -1, 0, + 2, 1, 2, -1, 2, -1, -1, 0, 0, -1, 2, -1, 0, 1, -1, 0, 2, 1, 2, -1, 2, + -1, -1, 0, 2, 0, 2, 2, 2, 2, -1, -1, 0, 2, 2, 2, 0, 0, -1, -1, 2, 0, 2, + 2, 2, 2, -1, -1, 0, 2, 2, 2, 0, 0, -1, -1, 2, -1, 0, -1, 2, 1, 1, 0, 0, + 1, 0, -1, 0, -1, 1, 0, 2, -1, 0, -1, 2, 1, 1, 0, 0, 1, 0, -1, 0, -1, 1, + 0, -1, -1, -1, 1, -1, 1, 0, 2, 1, 1, -1, 1, 1, -1, 0, 2, -1, -1, -1, 1, + -1, 1, 0, 2, 1, 1, -1, 1, 1, -1, 0, 2, -1, 2, 1, 2, -1, 0, 2, -1, 1, 0, + 1, 2, 1, 2, 2, -1, -1, 2, 1, 2, -1, 0, 2, -1, 1, 0, 1, 2, 1, 2, 2, -1, + 0, 2, 0, 0, 0, 0, 1, 1, 2, 0, 0, 0, 2, 2, 1, 1, 0, 2, 0, 0, 0, 0, 1, 1, + 2, 0, 0, 0, 2, 2, 1, 1, 2, -1, 0, -1, 2, 1, 1, 0, 0, 1, 0, -1, 0, -1, + 1, 0, 2, -1, 0, -1, 2, 1, 1, 0, 0, 1, 0, -1, 0, -1, 1, 0, 1, 1, 1, -1, + 1, -1, 2, 0, -1, -1, 1, -1, -1, 1, 2, 0, 1, 1, 1, -1, 1, -1, 2, 0, -1, + -1, 1, -1, -1, 1, 2, 0, -1, 2, 1, 2, -1, 0, 2, -1, 1, 0, 1, 2, 1, 2, 2, + -1, -1, 2, 1, 2, -1, 0, 2, -1, 1, 0, 1, 2, 1, 2, 2, -1, -1, 1, -1, -1, + -1, -1, 0, 0, 1, -1, -1, -1, 1, 1, 0, 0, -1, 1, -1, -1, -1, -1, 0, 0, + 1, -1, -1, -1, 1, 1, 0, 0, -1, 0, 1, 0, -1, 2, 2, 1, 1, 2, 1, 0, 1, 0, + 2, 1, -1, 0, 1, 0, -1, 2, 2, 1, 1, 2, 1, 0, 1, 0, 2, 1, -1, -1, -1, 1, + -1, 1, 0, 2, 1, 1, -1, 1, 1, -1, 0, 2, -1, -1, -1, 1, -1, 1, 0, 2, 1, + 1, -1, 1, 1, -1, 0, 2, 0, -1, 2, -1, 0, 1, -1, 0, 2, 1, 2, -1, 2, -1, + -1, 0, 0, -1, 2, -1, 0, 1, -1, 0, 2, 1, 2, -1, 2, -1, -1, 0, -1, 1, -1, + -1, -1, -1, 0, 0, 1, -1, -1, -1, 1, 1, 0, 0, -1, 1, -1, -1, -1, -1, 0, + 0, 1, -1, -1, -1, 1, 1, 0, 0, 0, 1, 2, 1, 0, -1, -1, 2, 2, -1, 2, 1, 2, + 1, -1, 2, 0, 1, 2, 1, 0, -1, -1, 2, 2, -1, 2, 1, 2, 1, -1, 2, 1, 1, 1, + -1, 1, -1, 2, 0, -1, -1, 1, -1, -1, 1, 2, 0, 1, 1, 1, -1, 1, -1, 2, 0, + -1, -1, 1, -1, -1, 1, 2, 0, 1, 0, -1, 0, 1, 2, 0, 1, -1, 2, -1, 0, -1, + 0, 0, 1, 1, 0, -1, 0, 1, 2, 0, 1, -1, 2, -1, 0, -1, 0, 0, 1, 0, 0, 0, + 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, + 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, + 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, + 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, + 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, + 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, + 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, + 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, + 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, + 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, + 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, + 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, + 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, + 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, + 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, +0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, +0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, +0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, +0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, +0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, +0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, -1, 1, 1, 1, 1, 2, 2, +-1, 1, 1, 1, -1, -1, 2, 2, 1, -1, 1, 1, 1, 1, 2, 2, -1, 1, 1, 1, -1, -1, 2, 2, +2, 1, 0, 1, 2, -1, 1, 2, 0, -1, 0, 1, 0, 1, 1, 2, 2, 1, 0, 1, 2, -1, 1, 2, 0, +-1, 0, 1, 0, 1, 1, 2, 0, 0, 0, 2, 0, 2, 1, -1, 2, 2, 0, 2, 2, 0, 1, -1, 0, 0, +0, 2, 0, 2, 1, -1, 2, 2, 0, 2, 2, 0, 1, -1, 1, 2, -1, 2, 1, 0, 0, -1, -1, 0, +-1, 2, -1, 2, 0, -1, 1, 2, -1, 2, 1, 0, 0, -1, -1, 0, -1, 2, -1, 2, 0, -1, -1, +1, -1, -1, -1, -1, 0, 0, 1, -1, -1, -1, 1, 1, 0, 0, -1, 1, -1, -1, -1, -1, 0, +0, 1, -1, -1, -1, 1, 1, 0, 0, 2, 1, 0, 1, 2, -1, 1, 2, 0, -1, 0, 1, 0, 1, 1, 2, +2, 1, 0, 1, 2, -1, 1, 2, 0, -1, 0, 1, 0, 1, 1, 2, 2, 2, 2, 0, 2, 0, -1, 1, 0, +0, 2, 0, 0, 2, -1, 1, 2, 2, 2, 0, 2, 0, -1, 1, 0, 0, 2, 0, 0, 2, -1, 1, -1, 0, +1, 0, -1, 2, 2, 1, 1, 2, 1, 0, 1, 0, 2, 1, -1, 0, 1, 0, -1, 2, 2, 1, 1, 2, 1, +0, 1, 0, 2, 1, 2, 0, 2, 2, 2, 2, -1, -1, 0, 2, 2, 2, 0, 0, -1, -1, 2, 0, 2, 2, +2, 2, -1, -1, 0, 2, 2, 2, 0, 0, -1, -1, 1, 0, -1, 0, 1, 2, 0, 1, -1, 2, -1, 0, +-1, 0, 0, 1, 1, 0, -1, 0, 1, 2, 0, 1, -1, 2, -1, 0, -1, 0, 0, 1, 2, 2, 2, 0, 2, +0, -1, 1, 0, 0, 2, 0, 0, 2, -1, 1, 2, 2, 2, 0, 2, 0, -1, 1, 0, 0, 2, 0, 0, 2, +-1, 1, -1, 0, 1, 0, -1, 2, 2, 1, 1, 2, 1, 0, 1, 0, 2, 1, -1, 0, 1, 0, -1, 2, 2, +1, 1, 2, 1, 0, 1, 0, 2, 1, 1, -1, 1, 1, 1, 1, 2, 2, -1, 1, 1, 1, -1, -1, 2, 2, +1, -1, 1, 1, 1, 1, 2, 2, -1, 1, 1, 1, -1, -1, 2, 2, 2, 1, 0, 1, 2, -1, 1, 2, 0, +-1, 0, 1, 0, 1, 1, 2, 2, 1, 0, 1, 2, -1, 1, 2, 0, -1, 0, 1, 0, 1, 1, 2, -1, -1, +-1, 1, -1, 1, 0, 2, 1, 1, -1, 1, 1, -1, 0, 2, -1, -1, -1, 1, -1, 1, 0, 2, 1, 1, +-1, 1, 1, -1, 0, 2, -1, 0, 1, 0, -1, 2, 2, 1, 1, 2, 1, 0, 1, 0, 2, 1, -1, 0, 1, +0, -1, 2, 2, 1, 1, 2, 1, 0, 1, 0, 2, 1, 2, 0, 2, 2, 2, 2, -1, -1, 0, 2, 2, 2, +0, 0, -1, -1, 2, 0, 2, 2, 2, 2, -1, -1, 0, 2, 2, 2, 0, 0, -1, -1, 2, 1, 0, 1, +2, -1, 1, 2, 0, -1, 0, 1, 0, 1, 1, 2, 2, 1, 0, 1, 2, -1, 1, 2, 0, -1, 0, 1, 0, +1, 1, 2, 2, 2, 2, 0, 2, 0, -1, 1, 0, 0, 2, 0, 0, 2, -1, 1, 2, 2, 2, 0, 2, 0, +-1, 1, 0, 0, 2, 0, 0, 2, -1, 1, 1, 2, -1, 2, 1, 0, 0, -1, -1, 0, -1, 2, -1, 2, +0, -1, 1, 2, -1, 2, 1, 0, 0, -1, -1, 0, -1, 2, -1, 2, 0, -1, 1, -1, 1, 1, 1, 1, +2, 2, -1, 1, 1, 1, -1, -1, 2, 2, 1, -1, 1, 1, 1, 1, 2, 2, -1, 1, 1, 1, -1, -1, +2, 2, -1, 2, 1, 2, -1, 0, 2, -1, 1, 0, 1, 2, 1, 2, 2, -1, -1, 2, 1, 2, -1, 0, +2, -1, 1, 0, 1, 2, 1, 2, 2, -1, 1, 1, 1, -1, 1, -1, 2, 0, -1, -1, 1, -1, -1, 1, +2, 0, 1, 1, 1, -1, 1, -1, 2, 0, -1, -1, 1, -1, -1, 1, 2, 0, 2, -1, 0, -1, 2, 1, +1, 0, 0, 1, 0, -1, 0, -1, 1, 0, 2, -1, 0, -1, 2, 1, 1, 0, 0, 1, 0, -1, 0, -1, +1, 0]) + +def main(): + parser = argparse.ArgumentParser(description="Plot TII") + parser.add_argument('--frame', default='0', help='Which transmission frame to analyse', + required=False) + parser.add_argument('--samplerate', default='2048000', help='Sample rate', + required=False) + parser.add_argument('--iq-file', default='ofdm.iq', help='File to read', + required=False) + parser.add_argument('--test-combs', action='store_const', const=True, help='Test all comb pattern values', + required=False) + parser.add_argument('--old-algo', action='store_const', const=True, help='Do an analysis with the old algorithm', + required=False) + + cli_args = parser.parse_args() + + if cli_args.test_combs: + prepare_tii_patterns() + prepare_comb_pattern_run(cli_args, 1, 11) + prepare_comb_pattern_run(cli_args, 1, 12) + elif cli_args.old_algo: + prepare_tii_patterns() + plot_tii_once(cli_args) + else: + prepare_tii_patterns() + algo1(cli_args) + +tii_ini_template = open("tii.ini.template", "r").read() + +def prepare_comb_pattern_run(options, comb, pattern): + fd = open("tii.ini", "w") + fd.write(tii_ini_template.format(comb=comb, pattern=pattern, oldvariant=0)) + fd.close() + + subprocess.call(['./odr-dabmod', 'tii.ini']) + + ret = plot_tii_once(options) + if ret[0] == "ok": + _, c, p, variant = ret + if c == comb and p == pattern: + print("== OK {}".format(ret)) + else: + print("== WRONG {}".format(ret)) + else: + print("== FAIL {}".format(ret)) + + +def plot_tii_once(options): + oversample = int(int(options.samplerate) / 2048000) + num_samples_per_tf = oversample * (NullSize + NbSymbols * SymSize) + print("{} samples per TF".format(num_samples_per_tf)) + + N = 8 + num_samples = num_samples_per_tf * (int(options.frame) + N) + + frames = np.fromfile(options.iq_file, np.complex64, count=num_samples) + + #skip to frame + frames = frames[num_samples_per_tf * int(options.frame):] + print("{} samples loaded".format(len(frames))) + + #fig = pp.figure() + #fig.suptitle("TII analysis") + #ax1 = fig.add_subplot(N, 1, i+1) + #ax1.set_title("Frame {}".format(i)) + #ax1.plot(np.abs(null_fft)) + + # analyse N successive frames + for i in range(N): + frame = frames[num_samples_per_tf * i:num_samples_per_tf * (i + 1)] + + if np.max(np.abs(frame[:NullSize])) != 0: + print("Frame {} has nonzero power in null symbol".format(i)) + + # Take the NULL symbol from that frame, but skip the cyclic prefix and truncate + null_skip = NullSize - Spacing + frame_null = frame[oversample*null_skip:oversample*(null_skip + Spacing)] + + null_fft = np.fft.fft(frame_null, oversample*Spacing) + null_fft_abs = np.abs(null_fft) + + # The phase reference symbol, skip the cyclic prefix + frame_prs = frame[oversample*(NullSize + CyclicPrefix):oversample * (NullSize + CyclicPrefix + Spacing)] + prs_fft = np.fft.fft(frame_prs, oversample*Spacing) + + threshold = 0.01 + carrier_indices = np.nonzero(null_fft_abs > threshold)[0] + print("Frame {} has {} carriers".format(i, len(carrier_indices))) + + if len(carrier_indices) != 32 and len(carrier_indices) != 0: + print("Invalid number of TII carriers") + return ('fail', ) + + def ix_to_k(ix): + "Convert a FFT index to a carrier index k" + # ix goes from 0 to 2047 + # 0..1024 positive, 1024..2047 negative, 0 is DC, 1..768 and 1280..2047 have power. + carrieroffset = int(oversample*Spacing) + if ix <= 1024: + return ix + else: + return ix - carrieroffset + + analyse_prs = False + if len(carrier_indices) and not analyse_prs: + print("k =" + " ".join("{:>4}".format(ix_to_k(ix)) for ix in carrier_indices)) + print("phases =" + " ".join("{:>4}".format(int(round(p))) for p in np.angle(prs_fft[null_fft_abs > threshold], deg=True))) + + phases_tii = np.angle(null_fft[null_fft_abs > threshold], deg=True) + phases_prs = np.angle(prs_fft[null_fft_abs > threshold], deg=True) + phase_error_wrong_impl = np.mod(np.asarray(np.around(phases_tii - phases_prs), dtype=np.int), 360) + + print("phases wrong =" + " ".join("{:>4}".format(p) for p in phase_error_wrong_impl)) + + pair_indices = carrier_indices.reshape(-1,2)[...,0] + correct_phase_indices = np.stack([pair_indices, pair_indices]).T.flatten() + phases_prs2 = np.angle(prs_fft[correct_phase_indices], deg=True) + phase_error_correct_impl = np.mod(np.asarray(np.around(phases_tii - phases_prs2), dtype=np.int), 360) + print("phases correct =" + " ".join("{:>4}".format(p) for p in phase_error_correct_impl)) + + variant = "unknown" + if np.all(phase_error_correct_impl == 0): + print("TII uses correct implementation variant") + variant = "correct" + elif np.all(phase_error_wrong_impl == 0): + print("TII uses old wrong implementation variant") + variant = "wrong" + else: + print("TII is wrong") + + carrier_indices_t = tuple(sorted(ix_to_k(c) for c in carrier_indices)) + print(carrier_indices_t) + if carrier_indices_t in all_tii_patterns: + c, p = all_tii_patterns[carrier_indices_t] + print("TII is using comb {} and pattern {}".format(c, p)) + return ("ok", c, p, variant) + else: + print("Unrecognised TII comb and pattern") + + return ("fail", variant) + + if analyse_prs: + delta = ((prs_phases[null_fft_abs > threshold] * 90) - np.angle(prs_fft[null_fft_abs > threshold], deg=True)) + print("prs=" + " ".join("{:>4}".format(int(round(p)) % 360) for p in delta)) + + #pp.show() + +def algo1(options): + oversample = 1 + if int(options.samplerate) != 2048000: + print("oversampling not supported for new algo") + sys.exit(1) + + num_samples_per_tf = NullSize + NbSymbols * SymSize + print("{} samples per TF".format(num_samples_per_tf)) + + def ix_to_k(ix): + "Convert a FFT index to a carrier index k" + # ix goes from 0 to 2047 + # 0..1024 positive, 1024..2047 negative, 0 is DC, 1..768 and 1280..2047 have power. + carrieroffset = int(Spacing) + if ix <= 1024: + return ix + else: + return ix - carrieroffset + + N = 8 + num_samples = num_samples_per_tf * (int(options.frame) + N) + + frames = np.fromfile(options.iq_file, np.complex64, count=num_samples) + + #skip to frame + frames = frames[num_samples_per_tf * int(options.frame):] + print("{} samples loaded".format(len(frames))) + + # analyse N successive frames + for i in range(N): + frame = frames[num_samples_per_tf * i:num_samples_per_tf * (i + 1)] + + # Take the NULL symbol from that frame, but skip the cyclic prefix and truncate + null_skip = NullSize - Spacing + frame_null = frame[null_skip:(null_skip + Spacing)] + + null_fft = np.fft.fft(frame_null, Spacing) + + # The phase reference symbol, skip the cyclic prefix + frame_prs = frame[NullSize + CyclicPrefix:NullSize + CyclicPrefix + Spacing] + prs_fft = np.fft.fft(frame_prs, Spacing) + + # In TM1, the carriers repeat four times: + # [-768, -384[ + # [-384, 0[ + # ]0, 384] + # ]384, 768] + # A consequence of the fact that the 0 bin is never used is that the + # first carrier of each pair is even for negative k, odd for positive k + + blocks = [null_fft[-769:-385], null_fft[-385:-1], null_fft[1:385], null_fft[385:769]] + blocks_multiplied = np.zeros(384//2, dtype=np.complex128) + for block in blocks: + even_odd = block.reshape(-1, 2) + b = even_odd[...,0] * np.conj(even_odd[...,1]) + blocks_multiplied += b + + threshold = 0.01 + indices_above_threshold = np.nonzero(np.abs(blocks_multiplied) > threshold)[0] + carrier_indices = indices_above_threshold * 2 + 1 + + if np.max(np.abs(frame[:NullSize])) != 0: + print("Frame {} has nonzero power in null symbol and {} TII carriers".format(i, len(carrier_indices))) + + if len(carrier_indices) >= 4: + carriers = [ix_to_k(ix) for ix in carrier_indices] + + cp_counts = {} + for k in carriers: + if k in cp_per_carrier: + cps = cp_per_carrier[k] + for cp in cps: + if cp not in cp_counts: + cp_counts[cp] = 0 + cp_counts[cp] += 1 + + if len(cp_counts) > 0: + for cp in cp_counts: + if cp_counts[cp] >= 4: + c, p = cp + valid, delay, err = analyse_phase(c, p, null_fft, prs_fft) + if valid: + print("TII likelihood {}: comb {} and pattern {}, delay {} samples, err {}".format( + cp_counts[cp], c, p, + delay, err)) + else: + print("TII likelihood {}: comb {} and pattern {}, invalid delay measurement, err {}".format( + cp_counts[cp], c, p, err)) + else: + print("Unrecognised TII comb and pattern") + +def convert_angles(angle): + """Convert 0 to 360 degree angles to 0 to 180, -179 to -1""" + if angle > 180: + return angle - 360 + elif angle < -180: + return angle + 360 + else: + return angle +phase_corrector = np.vectorize(convert_angles) + +def analyse_phase(c, p, null_fft, prs_fft): + carriers = np.array(calculate_carrier_indices(c, p)) + + #print("k =" + " ".join("{:>4}".format(k) for k in carriers)) + #print("phases =" + " ".join("{:>4}".format(int(round(p))) for p in np.angle(prs_fft[carriers], deg=True))) + + phases_tii = np.angle(null_fft[carriers], deg=True) + + pair_indices = carriers.reshape(-1,2)[...,0] + correct_phase_indices = np.stack([pair_indices, pair_indices]).T.flatten() + phases_prs2 = np.angle(prs_fft[correct_phase_indices], deg=True) + + fig = pp.figure() + fig.suptitle("Phase error C {} P {}".format(c, p)) + + phase_errors = {} + search_range = (-40, 40) + for err4 in range(*search_range): + err = err4 / 4 + rotate_vec = np.exp(2j * np.pi * err * carriers / 2048) + rotated = np.angle(null_fft[carriers] * rotate_vec, deg=True) + delta = np.mod(np.asarray(np.around(rotated - phases_prs2), dtype=np.int), 360) + delta = phase_corrector(delta) + #print("phases {:>3} =".format(err) + " ".join("{:>4}".format(p) for p in delta)) + phase_errors[err] = delta + + #ax1 = fig.add_subplot(7, 1, i+1) + #ax1.set_title("err {}".format(err)) + #ax1.plot(correct_phase_indices, delta, '.') + #ax1.plot(correct_phase_indices, slopepoly(correct_phase_indices)) + #pp.show() + + best_err = min(phase_errors.items(), key=(lambda e: np.sum(np.abs(e[1])))) + #print("Best err: {}".format(best_err[0])) + #print("phases =" + " ".join("{:>4}".format(p) for p in best_err[1])) + + meas_err = np.sum(np.abs(best_err[1])) + + if best_err[0] != search_range[0] and best_err[0] != search_range[1]: + return (True, best_err[0], meas_err) + else: + return (False, best_err[0], meas_err) + +main() + +# The MIT License (MIT) +# +# 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. |