From daf2209c739d8885ba16a4083ff70f3a8ab28fc0 Mon Sep 17 00:00:00 2001 From: Eric Teunis de Boone Date: Mon, 16 Jan 2023 19:42:21 +0100 Subject: [PATCH] ZH: new function read_clock_offsets This is moved out from the scripts that already have to modify timing --- .../aa_generate_beacon.py | 29 +++++++++++++++++++ .../ca_period_from_shower.py | 7 ++--- .../da_reconstruction.py | 8 ++--- 3 files changed, 35 insertions(+), 9 deletions(-) diff --git a/simulations/airshower_beacon_simulation/aa_generate_beacon.py b/simulations/airshower_beacon_simulation/aa_generate_beacon.py index d56d232..f3f3321 100755 --- a/simulations/airshower_beacon_simulation/aa_generate_beacon.py +++ b/simulations/airshower_beacon_simulation/aa_generate_beacon.py @@ -20,6 +20,35 @@ tx_fname = 'tx.json' antennas_fname = 'antennas.hdf5' c_light = lib.c_light +def read_antenna_clock_repair_offsets(antennas, mode='all', freq_name=None): + valid_modes = ['orig', 'ks', 'phases', 'all'] + + time_offsets = [] + for i, ant in enumerate(antennas): + _clock_delta = 0 + + # original timing + if mode == 'orig': + _clock_delta = ant.attrs['clock_offset'] + + # phase + if mode in ['all', 'phases']: + clock_phase = ant.beacon_info[freq_name]['sigma_phase_mean'] + f_beacon = ant.beacon_info[freq_name]['freq'] + clock_phase_time = clock_phase/(2*np.pi*f_beacon) + + _clock_delta += clock_phase_time + + # ks + if mode in ['all', 'ks']: + best_k_time = ant.beacon_info[freq_name]['best_k_time'] + + _clock_delta += best_k_time + + time_offsets.append(_clock_delta) + + return time_offsets + def write_tx_file(fname, tx, f_beacon, **kwargs): with open(fname, 'w') as fp: return json.dump( diff --git a/simulations/airshower_beacon_simulation/ca_period_from_shower.py b/simulations/airshower_beacon_simulation/ca_period_from_shower.py index 88654c0..1f0a3c7 100755 --- a/simulations/airshower_beacon_simulation/ca_period_from_shower.py +++ b/simulations/airshower_beacon_simulation/ca_period_from_shower.py @@ -201,14 +201,13 @@ if __name__ == "__main__": # Remove time due to true phase # and optionally remove the beacon + measured_repair_offsets = beacon.read_antenna_clock_repair_offsets(ev.antennas, mode='phases', freq_name=freq_name) for i, ant in enumerate(ev.antennas): - clock_phase = ant.beacon_info[freq_name]['sigma_phase_mean'] - clock_phase_time = clock_phase/(2*np.pi*f_beacon) - ev.antennas[i].orig_t = ev.antennas[i].t_AxB - ev.antennas[i].t_AxB += clock_phase_time + ev.antennas[i].t_AxB += measured_repair_offsets[i] if remove_beacon_from_trace: + clock_phase = measured_repair_offsets[i]*2*np.pi*f_beacon meas_phase = ant.beacon_info[freq_name]['phase'] f = ant.beacon_info[freq_name]['freq'] ampl = ant.beacon_info[freq_name]['amplitude'] diff --git a/simulations/airshower_beacon_simulation/da_reconstruction.py b/simulations/airshower_beacon_simulation/da_reconstruction.py index 6897ef0..3c70a3b 100755 --- a/simulations/airshower_beacon_simulation/da_reconstruction.py +++ b/simulations/airshower_beacon_simulation/da_reconstruction.py @@ -66,12 +66,10 @@ if __name__ == "__main__": f_beacon = ev.antennas[0].beacon_info[freq_name]['freq'] # Repair clock offsets with the measured offsets + measured_repair_offsets = beacon.read_antenna_clock_repair_offsets(ev.antennas, mode='all', freq_name=freq_name) for i, ant in enumerate(ev.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'] - - total_clock_time = best_k_time + clock_phase_time - ev.antennas[i].t += total_clock_time + # t_AxB will be set by the rit.set_pol_and_bp function + ev.antennas[i].t += measured_repair_offsets[i] N_X, Xlow, Xhigh = 23, 100, 1200 with joblib.parallel_backend("loky"):