From be3041f4a23bddf3966ca16ead6ba9dccb85fce6 Mon Sep 17 00:00:00 2001
From: Eric Teunis de Boone <ericteunis@deboone.nl>
Date: Wed, 28 Sep 2022 15:10:59 +0200
Subject: [PATCH] ZH: modify clocks of all antennas

---
 .../airshower_beacon_simulation/README.md     | 14 ++++
 .../generate_beacon.py                        |  1 -
 .../modify_clocks.py                          | 77 +++++++++++++++++++
 3 files changed, 91 insertions(+), 1 deletion(-)
 create mode 100755 simulations/airshower_beacon_simulation/modify_clocks.py

diff --git a/simulations/airshower_beacon_simulation/README.md b/simulations/airshower_beacon_simulation/README.md
index c670c8c..98e8446 100644
--- a/simulations/airshower_beacon_simulation/README.md
+++ b/simulations/airshower_beacon_simulation/README.md
@@ -5,3 +5,17 @@ Using the beacon, the clocks can be calibrated.
 
 The ZHaires simulated airshower is stored in [./ZH_airshower](./ZH_airshower).
 The produced files can be read using [./earsim](./earsim).
+
+
+Steps:
+ 1. Setup
+	1. Beacon ([./generate_beacon.py])
+		1. Define tx position
+		2. Read in antennas
+		3. Sample beacon at each antenna
+		4. Add to relevant polarisation
+		5. Save antenna traces
+
+	2. Timeoffset ([./generate_timeoffset.pyt])
+		1. Generate timeoffsets for each antenna
+		2. Modify time samples
diff --git a/simulations/airshower_beacon_simulation/generate_beacon.py b/simulations/airshower_beacon_simulation/generate_beacon.py
index d9efbf1..ee0ae35 100755
--- a/simulations/airshower_beacon_simulation/generate_beacon.py
+++ b/simulations/airshower_beacon_simulation/generate_beacon.py
@@ -12,7 +12,6 @@ from earsim import REvent, Antenna
 import lib
 
 tx_fname = 'tx.json'
-clocks_fname = 'clocks.csv'
 antennas_fname = 'antennas.hdf5'
 
 def write_tx_file(fname, tx, f_beacon):
diff --git a/simulations/airshower_beacon_simulation/modify_clocks.py b/simulations/airshower_beacon_simulation/modify_clocks.py
new file mode 100755
index 0000000..915686d
--- /dev/null
+++ b/simulations/airshower_beacon_simulation/modify_clocks.py
@@ -0,0 +1,77 @@
+#!/usr/bin/env python3
+
+import numpy as np
+import json
+import h5py
+import os.path as path
+
+from copy import deepcopy as copy
+
+from earsim import REvent, Antenna
+
+import generate_beacon as beacon
+
+import lib
+
+clocks_fname = 'clocks.csv'
+
+
+if __name__ == "__main__":
+    from os import path
+    import sys
+
+    max_clock_offset = 100# ns
+    remake_clock_offsets = False
+
+    seed = 12345
+    rng = np.random.default_rng(seed)
+
+    fname = "ZH_airshower/mysim.sry"
+
+    ####
+    fname_dir = path.dirname(fname)
+    clocks_fname = path.join(fname_dir, clocks_fname)
+    antennas_fname = path.join(fname_dir, beacon.antennas_fname)
+
+
+    if path.isfile(clocks_fname) and not remake_clock_offsets:
+        print("Clock offsets previously generated")
+        sys.exit()
+
+    if not path.isfile(antennas_fname):
+        print("Antenna file cannot be found, did you try generating a beacon?")
+        sys.exit(1)
+    
+    # read in antennas
+    with h5py.File(antennas_fname, 'a') as fp:
+        if 'antennas' not in fp.keys():
+            print("Antenna file corrupted?")
+            sys.exit(1)
+
+        group = fp['antennas']
+
+        N_antennas = len(group.keys())
+
+        clock_offsets = max_clock_offset * (2*rng.uniform(size=N_antennas) - 1)
+
+        # modify time values of each antenna
+        for i, name in enumerate(group.keys()):
+            ant_group = group[name]
+            clk_offset = clock_offsets[i]
+
+            ant_attrs = ant_group.attrs
+            if 'clock_offset' in ant_attrs:
+                tmp_offset = ant_attrs['clock_offset']
+                if remake_clock_offsets:
+                    ant_group['traces'][0, :] -= tmp_offset 
+                else:
+                    clock_offsets[i] = tmp_offset
+                    continue
+
+            ant_attrs['clock_offset'] = clk_offset
+            ant_group['traces'][0, :] += clk_offset 
+
+    # save to simple csv
+    np.savetxt(clocks_fname, clock_offsets)
+
+    print("Antenna clocks modified in " + str(antennas_fname))