From ed1eb3f5cb02a826c061054bdf59d19b3bb4ecfc Mon Sep 17 00:00:00 2001
From: Eric Teunis de Boone <ericteunis@deboone.nl>
Date: Mon, 9 Jan 2023 11:48:58 +0100
Subject: [PATCH] ZH: initial reporting for measured total clock offset per
 antenna

---
 ...cb_report_measured_antenna_time_offsets.py | 80 +++++++++++++++++++
 1 file changed, 80 insertions(+)
 create mode 100755 simulations/airshower_beacon_simulation/cb_report_measured_antenna_time_offsets.py

diff --git a/simulations/airshower_beacon_simulation/cb_report_measured_antenna_time_offsets.py b/simulations/airshower_beacon_simulation/cb_report_measured_antenna_time_offsets.py
new file mode 100755
index 0000000..d4ea78b
--- /dev/null
+++ b/simulations/airshower_beacon_simulation/cb_report_measured_antenna_time_offsets.py
@@ -0,0 +1,80 @@
+#!/usr/bin/env python3
+# vim: fdm=indent ts=4
+
+"""
+Report best time offset per frequency for each antenna
+"""
+
+import h5py
+import matplotlib.pyplot as plt
+import numpy as np
+from os import path
+from scipy.interpolate import interp1d
+
+import aa_generate_beacon as beacon
+import lib
+
+if __name__ == "__main__":
+    import sys
+    import os
+    import matplotlib
+    if os.name == 'posix' and "DISPLAY" not in os.environ:
+        matplotlib.use('Agg')
+
+    fname = "ZH_airshower/mysim.sry"
+
+    fig_dir = "./figures/periods_from_shower_figures/"
+    fig_subdir = path.join(fig_dir, 'shifts/')
+    show_plots = False
+
+    ####
+    fname_dir = path.dirname(fname)
+    antennas_fname = path.join(fname_dir, beacon.antennas_fname)
+    time_diffs_fname = 'time_diffs.hdf5' if not True else antennas_fname
+
+    # create fig_dir
+    if fig_dir:
+        os.makedirs(fig_dir, exist_ok=True)
+        if fig_subdir:
+            os.makedirs(fig_subdir, exist_ok=True)
+
+    # Read in antennas from file
+    _, tx, antennas = beacon.read_beacon_hdf5(antennas_fname)
+
+    # For now only implement using one freq_name
+    freq_names = antennas[0].beacon_info.keys()
+    if len(freq_names) > 1:
+        raise NotImplementedError
+
+    freq_name = next(iter(freq_names))
+    f_beacon = antennas[0].beacon_info[freq_name]['freq']
+
+
+    for i, ant in enumerate(antennas):
+        clock_phase_time = ant.beacon_info[freq_name]['sigma_phase_mean']/(2*np.pi*f_beacon)
+
+        best_k_time = ant.beacon_info[freq_name]['best_k_time']
+        best_k = ant.beacon_info[freq_name]['best_k_time']
+
+        total_clock_time = best_k_time + clock_phase_time
+
+        actual_clock_time = ant.attrs['clock_offset']
+
+        # filter k == 0
+        if best_k == 0:
+            continue
+
+        # Report to stdout
+        print("Timing of antenna", ant.name)
+        print(" + k:", best_k, ";(-)", best_k_time, 'ns')
+        print(" + phase_time:", clock_phase_time, 'ns')
+        print("==============================")
+        print("Sum: (-) ", total_clock_time, 'ns')
+        print("Actual:", actual_clock_time, 'ns')
+        print("Diff (A - - S):", actual_clock_time - -total_clock_time, 'ns')
+        print()
+                
+
+
+
+