diff --git a/simulations/airshower_beacon_simulation/.gitignore b/simulations/airshower_beacon_simulation/.gitignore new file mode 100644 index 0000000..5791952 --- /dev/null +++ b/simulations/airshower_beacon_simulation/.gitignore @@ -0,0 +1,3 @@ +config.mk +figures +data diff --git a/simulations/airshower_beacon_simulation/Makefile b/simulations/airshower_beacon_simulation/Makefile index d7e3fb1..c74642f 100644 --- a/simulations/airshower_beacon_simulation/Makefile +++ b/simulations/airshower_beacon_simulation/Makefile @@ -1,35 +1,61 @@ .PHONY: all - FIG_DIR := ./figures -all: beacon clocks phases findks reconstruct +-include config.mk +MAX_K ?= 3 +CLK_DEV ?= 20 +TRACE_N ?= 4096 +TRACE_PRE ?= 500 +NOISE_SIGMA ?= 1e4 + +BEAC_AMP ?= 1e3 0 0 +BEAC_F ?= 51.53e-3 +BEAC_DECAY ?= 0 + +PB_LOW ?= 0.03 +PB_HIGH ?= 0.08 + +REF_ANTS ?= + +N_MASK ?= 500 + +DATA_DIR ?= ./data +INPUT_DIR ?= ./ZH_airshower/ +SEED ?= 12345 + +all: beacon clocks phases findks vary-fixes reconstruct beacon: - ./aa_generate_beacon.py | tee figures/aa.log - ./ab_modify_clocks.py 0 | tee figures/ab.log - ./ac_show_signal_to_noise.py --no-show-plots --fig-dir=${FIG_DIR} - ./view_beaconed_antenna.py 72 -p x -p y -p z --no-show-plots --fig-dir=${FIG_DIR} + ./aa_generate_beacon.py -N ${TRACE_N} -P ${TRACE_PRE} -n ${NOISE_SIGMA} -a ${BEAC_AMP} -f ${BEAC_F} -l ${PB_LOW} -u ${PB_HIGH} -d ${BEAC_DECAY} --data-dir ${DATA_DIR} --input-fname ${INPUT_DIR} | tee ${FIG_DIR}/aa.log + ./ab_modify_clocks.py 0 --data-dir ${DATA_DIR} | tee ${FIG_DIR}/ab.log + ./ac_show_signal_to_noise.py --no-show-plots --fig-dir=${FIG_DIR} --data-dir ${DATA_DIR} + ./view_beaconed_antenna.py 72 -p x -p y -p z -p n -p b --ft --no-show-plots --fig-dir=${FIG_DIR} --data-dir ${DATA_DIR} clocks: - ./ab_modify_clocks.py 15 --gaussian + ./ab_modify_clocks.py ${CLK_DEV} --seed ${SEED} --gaussian --data-dir ${DATA_DIR} phases: - ./ba_measure_beacon_phase.py --no-show-plots --fig-dir=${FIG_DIR} - ./bb_measure_clock_phase.py --no-show-plots --fig-dir=${FIG_DIR} - ./bc_baseline_phase_deltas.py --no-show-plots --fig-dir=${FIG_DIR} - ./bd_antenna_phase_deltas.py --no-show-plots --fig-dir=${FIG_DIR} + ./ba_measure_beacon_phase.py --N-mask ${N_MASK} --no-show-plots --fig-dir=${FIG_DIR} --data-dir ${DATA_DIR} + ./bb_measure_clock_phase.py --no-show-plots --fig-dir=${FIG_DIR} --data-dir ${DATA_DIR} + ./bc_baseline_phase_deltas.py ${REF_ANTS} --no-show-plots --fig-dir=${FIG_DIR} --data-dir ${DATA_DIR} + ./bd_antenna_phase_deltas.py --no-show-plots --fig-dir=${FIG_DIR} --data-dir ${DATA_DIR} findks: - ./ca_period_from_shower.py --no-show-plots --fig-dir=${FIG_DIR} - ./cb_report_measured_antenna_time_offsets.py --no-show-plots --fig-dir=${FIG_DIR} + ./ca_period_from_shower.py --input-fname ${INPUT_DIR} --max-k ${MAX_K} --no-show-plots --fig-dir=${FIG_DIR} --data-dir ${DATA_DIR} -l ${PB_LOW} -u ${PB_HIGH} + ./cb_report_measured_antenna_time_offsets.py --no-show-plots --fig-dir=${FIG_DIR} --data-dir ${DATA_DIR} + +vary-fixes: + ./dc_grid_power_time_fixes.py --no-show-plots --fig-dir=${FIG_DIR} --data-dir ${DATA_DIR} --input-fname ${INPUT_DIR} reconstruct: - ./da_reconstruction.py --no-show-plots --fig-dir=${FIG_DIR} - ./db_longitudinal_figure.py --no-show-plots --fig-dir=${FIG_DIR} - ./dc_grid_power_time_fixes.py --no-show-plots --fig-dir=${FIG_DIR} + ./da_reconstruction.py --no-show-plots --fig-dir=${FIG_DIR} --data-dir ${DATA_DIR} --input-fname ${INPUT_DIR} + ./db_longitudinal_figure.py --no-show-plots --fig-dir=${FIG_DIR} --data-dir ${DATA_DIR} + dist-clean: - rm -f ZH_airshower/antennas.hdf5 - rm -f ZH_airshower/clocks.csv - rm -f ZH_airshower/tx.json + rm -vf ${DATA_DIR}/antennas.hdf5 + rm -vf ${DATA_DIR}/ca_breaked_run + rm -vf ${DATA_DIR}/res.pkl + rm -vf ${DATA_DIR}/clocks.csv + rm -vf ${DATA_DIR}/tx.json diff --git a/simulations/airshower_beacon_simulation/aa_generate_beacon.py b/simulations/airshower_beacon_simulation/aa_generate_beacon.py index 0f77fd2..8fb41a6 100755 --- a/simulations/airshower_beacon_simulation/aa_generate_beacon.py +++ b/simulations/airshower_beacon_simulation/aa_generate_beacon.py @@ -282,23 +282,31 @@ if __name__ == "__main__": # Bandpass parser.add_argument('-p', '--use-passband', type=bool, default=True, help='(Default: %(default)d)') - parser.add_argument('-l', '--passband-low', type=float, default=30e-3, help='Lower frequency [GHz] of the passband filter. (set -1 for np.inf) (Default: %(default)d)') - parser.add_argument('-u', '--passband-high', type=float, default=80e-3, help='Upper frequency [GHz] of the passband filter. (set -1 for np.inf) (Default: %(default)d)') + parser.add_argument('-l', '--passband-low', type=float, default=30e-3, help='Lower frequency [GHz] of the passband filter. (set -1 for np.inf) (Default: %(default)g)') + parser.add_argument('-u', '--passband-high', type=float, default=80e-3, help='Upper frequency [GHz] of the passband filter. (set -1 for np.inf) (Default: %(default)g)') # Trace length modification parser.add_argument('-N', '--new-trace-length', type=float, help='resize airshower trace (Default: %(default)d)', default=1e4) parser.add_argument('-P', '--pre-trace-length', type=float, help='amount of trace to prepend the airshower when resizing (Default: %(default)d)', default=2e3) + # Input directory + parser.add_argument('--input-fname', type=str, default=None, help='Path to mysim.sry, either directory or path. If empty it takes DATA_DIR and appends mysim.sry. (Default: %(default)s)') + parser.add_argument('--data-dir', type=str, default="./ZH_airshower", help='Path to Data Directory. (Default: %(default)s)') + args = parser.parse_args() + if not args.input_fname: + args.input_fname = args.data_dir + + if path.isdir(args.input_fname): + args.input_fname = path.join(args.input_fname, "mysim.sry") + ## ## End of ArgumentParsing ## rng = np.random.default_rng() - fname = "ZH_airshower/mysim.sry" - # Noise properties noise_sigma = args.noise_sigma # mu V/m set to None to ignore @@ -337,7 +345,7 @@ if __name__ == "__main__": block_filter = lambda x, dt, low, high: x #### - fname_dir = path.dirname(fname) + fname_dir = args.data_dir tx_fname = path.join(fname_dir, tx_fname) antennas_fname = path.join(fname_dir, antennas_fname) @@ -354,7 +362,7 @@ if __name__ == "__main__": print("Noise sigma [muV/m]:", noise_sigma) # read in antennas - ev = REvent(fname) + ev = REvent(args.input_fname) N_antennas = len(ev.antennas) # initialize hdf5 file diff --git a/simulations/airshower_beacon_simulation/ab_modify_clocks.py b/simulations/airshower_beacon_simulation/ab_modify_clocks.py index 0c33300..1532cdd 100755 --- a/simulations/airshower_beacon_simulation/ab_modify_clocks.py +++ b/simulations/airshower_beacon_simulation/ab_modify_clocks.py @@ -23,22 +23,23 @@ if __name__ == "__main__": parser = ArgumentParser() parser.add_argument('max_clock_offset', nargs='?', type=float, default=25, help='(Default: %(default)d)') + parser.add_argument('-s', '--seed', type=int, nargs='?', default=None, help='Fix seed if supplied.') parser.add_argument('--uniform', action='store_const', const='uniform', dest='dist_type') parser.add_argument('--gaussian', action='store_const', const='gauss', dest='dist_type') parser.set_defaults(dist_type='gauss') + parser.add_argument('--data-dir', type=str, default="./data", help='Path to Data Directory. (Default: %(default)s)') + args = parser.parse_args() max_clock_offset = args.max_clock_offset # ns remake_clock_offsets = True - seed = 12345 + seed = args.seed rng = np.random.default_rng(seed) - fname = "ZH_airshower/mysim.sry" - #### - fname_dir = path.dirname(fname) + fname_dir = args.data_dir clocks_fname = path.join(fname_dir, clocks_fname) antennas_fname = path.join(fname_dir, beacon.antennas_fname) diff --git a/simulations/airshower_beacon_simulation/ac_show_signal_to_noise.py b/simulations/airshower_beacon_simulation/ac_show_signal_to_noise.py index cda7241..1c8fad9 100755 --- a/simulations/airshower_beacon_simulation/ac_show_signal_to_noise.py +++ b/simulations/airshower_beacon_simulation/ac_show_signal_to_noise.py @@ -33,14 +33,13 @@ if __name__ == "__main__": args = parser.parse_args() - fname = "ZH_airshower/mysim.sry" - figsize = (9,6) + figsize = (12,8) fig_dir = args.fig_dir show_plots = args.show_plots #### - fname_dir = path.dirname(fname) + fname_dir = args.data_dir antennas_fname = path.join(fname_dir, beacon.antennas_fname) tx_fname = path.join(fname_dir, beacon.tx_fname) diff --git a/simulations/airshower_beacon_simulation/ba_measure_beacon_phase.py b/simulations/airshower_beacon_simulation/ba_measure_beacon_phase.py index 613acda..235de64 100755 --- a/simulations/airshower_beacon_simulation/ba_measure_beacon_phase.py +++ b/simulations/airshower_beacon_simulation/ba_measure_beacon_phase.py @@ -23,25 +23,31 @@ if __name__ == "__main__": from scriptlib import MyArgumentParser parser = MyArgumentParser() - parser.add_argument('--no-use-AxB', dest='use_AxB_trace', action='store_false') + + group1 = parser.add_mutually_exclusive_group() + group1.add_argument('--AxB', dest='use_AxB_trace', action='store_true', help='Only use AxB trace, if both AxB and beacon are not used, we use the antenna polarisations.') + group1.add_argument('--beacon', dest='use_beacon_trace', action='store_true', help='Only use the beacon trace') + + parser.add_argument('--N-mask', type=float, default=500, help='Mask N_MASK samples around the absolute maximum of the trace. (Default: %(default)d)') + args = parser.parse_args() f_beacon_band = (49e-3,55e-3) #GHz allow_frequency_fitting = False read_frequency_from_file = True + N_mask = int(args.N_mask) - use_AxB_trace = args.use_AxB_trace - use_beacon_trace = True # only applicable if AxB = False + use_only_AxB_trace = args.use_AxB_trace + use_only_beacon_trace = args.use_beacon_trace # only applicable if AxB = False show_plots = args.show_plots - fname = "ZH_airshower/mysim.sry" - figsize = (9,6) + figsize = (12,8) - print("use_AxB_trace:", use_AxB_trace, "use_beacon_trace:",use_beacon_trace) + print("use_only_AxB_trace:", use_only_AxB_trace, "use_only_beacon_trace:", use_only_beacon_trace) #### - fname_dir = path.dirname(fname) + fname_dir = args.data_dir antennas_fname = path.join(fname_dir, beacon.antennas_fname) fig_dir = args.fig_dir # set None to disable saving @@ -84,59 +90,64 @@ if __name__ == "__main__": h5ant = group[name] # use E_AxB only instead of polarisations - if use_AxB_trace: - if 'E_AxB' not in h5ant.keys(): - print(f"Antenna does not have 'E_AxB' in {name}") + if use_only_AxB_trace: + traces_key = 'E_AxB' + if traces_key not in h5ant.keys(): + print(f"Antenna does not have '{traces_key}' in {name}") sys.exit(1) - traces = h5ant['E_AxB'] - + traces = h5ant[traces_key] t_trace = traces[0] test_traces = [ traces[1] ] orients = ['E_AxB'] - # TODO: refine masking - # use beacon but remove where E_AxB-Beacon != 0 - if True: - N_pre, N_post = 250, 250 # TODO: make this configurable + # Only beacon + elif use_only_beacon_trace: + traces_key = 'filtered_traces' + if traces_key not in h5ant.keys(): + print(f"Antenna file corrupted? no '{traces_key}' in {name}") + sys.exit(1) - max_idx = np.argmax(test_traces[0]) - - low_idx = max(0, max_idx-N_pre) - high_idx = min(len(t_trace), max_idx+N_post) - - t_mask = np.ones(len(t_trace), dtype=bool) - t_mask[low_idx:high_idx] = False - - t_trace = t_trace[t_mask] - for j, t in enumerate(test_traces): - test_traces[j] = t[t_mask] - orients[j] = orients[j] + ' masked' + traces = h5ant[traces_key] + t_trace = traces[0] + test_traces = [traces[4]] + orients = ['B'] # use separate polarisations else: - if 'traces' not in h5ant.keys(): - print(f"Antenna file corrupted? no 'traces' in {name}") + traces_key = 'filtered_traces' + if traces_key not in h5ant.keys(): + print(f"Antenna file corrupted? no '{traces_key}' in {name}") sys.exit(1) - traces = h5ant['traces'] + traces = h5ant[traces_key] t_trace = traces[0] + test_traces = [traces[i] for i in range(1,4)] + orients = ['Ex', 'Ey', 'Ez'] - if use_beacon_trace: - # only take the Beacon trace - test_traces = [traces[4]] - orients = ['B'] - else: - test_traces = traces[1:] - orients = ['Ex', 'Ey', 'Ez', 'B'] + # Really only select the first component + if True: + test_traces = [test_traces[0]] + orients = [orients[0]] - # modify the length of the traces - if False: - t_trace = t_trace[:len(t_trace)//2] - half_traces = [] - for trace in test_traces: - half_traces.append( trace[:len(trace)//2]) - test_traces = half_traces + # TODO: refine masking + # use beacon but remove where E_AxB-Beacon != 0 + # Uses the first traces as reference + if N_mask and orients[0] != 'B': + N_pre, N_post = N_mask//2, N_mask//2 + + max_idx = np.argmax(test_traces[0]) + + low_idx = max(0, max_idx-N_pre) + high_idx = min(len(t_trace), max_idx+N_post) + + t_mask = np.ones(len(t_trace), dtype=bool) + t_mask[low_idx:high_idx] = False + + t_trace = t_trace[t_mask] + for j, t in enumerate(test_traces): + test_traces[j] = t[t_mask] + orients[j] = orients[j] + ' masked' # Do Fourier Transforms # to find phases and amplitudes @@ -156,7 +167,7 @@ if __name__ == "__main__": amps = [ 3e-7 ] # choose highest amp - idx = 0 + idx = np.argmax(amps) if False and len(beacon_phases) > 1: #idx = np.argmax(amplitudes, axis=-1) raise NotImplementedError @@ -176,7 +187,7 @@ if __name__ == "__main__": # for reporting using plots found_data[i] = frequency, beacon_phase, amplitude - if (show_plots or fig_dir) and (i == 0 or i == 72 or i == 70): + if (show_plots or fig_dir) and (i == 0 or i == 72): p2t = lambda phase: phase/(2*np.pi*f_beacon) fig, ax = plt.subplots(figsize=figsize) diff --git a/simulations/airshower_beacon_simulation/bb_measure_clock_phase.py b/simulations/airshower_beacon_simulation/bb_measure_clock_phase.py index a16fa76..a83f717 100755 --- a/simulations/airshower_beacon_simulation/bb_measure_clock_phase.py +++ b/simulations/airshower_beacon_simulation/bb_measure_clock_phase.py @@ -22,8 +22,7 @@ if __name__ == "__main__": parser = MyArgumentParser() args = parser.parse_args() - fname = "ZH_airshower/mysim.sry" - figsize = (9,6) + figsize = (12,8) c_light = lib.c_light show_plots = args.show_plots @@ -31,7 +30,7 @@ if __name__ == "__main__": remove_absolute_phase_offset_minimum = True #### - fname_dir = path.dirname(fname) + fname_dir = args.data_dir antennas_fname = path.join(fname_dir, beacon.antennas_fname) fig_dir = args.fig_dir # set None to disable saving diff --git a/simulations/airshower_beacon_simulation/bc_baseline_phase_deltas.py b/simulations/airshower_beacon_simulation/bc_baseline_phase_deltas.py index 887cd0e..a3c067f 100755 --- a/simulations/airshower_beacon_simulation/bc_baseline_phase_deltas.py +++ b/simulations/airshower_beacon_simulation/bc_baseline_phase_deltas.py @@ -21,16 +21,16 @@ if __name__ == "__main__": from scriptlib import MyArgumentParser parser = MyArgumentParser() + parser.add_argument('ref_ant_idx', default=None, nargs='*', type=int, help='Reference Antenna Indices for Baselines(ref_ant_idx, *). Leave empty to use all available antennas. (Default: %(default)s) ') args = parser.parse_args() - fname = "ZH_airshower/mysim.sry" - figsize = (9,6) + figsize = (12,8) c_light = 3e8*1e-9 show_plots = args.show_plots - ref_ant_id = None if True else [50] # leave None for all baselines + ref_ant_id = args.ref_ant_idx # leave None for all baselines #### - fname_dir = path.dirname(fname) + fname_dir = args.data_dir antennas_fname = path.join(fname_dir, beacon.antennas_fname) time_diffs_fname = 'time_diffs.hdf5' if False else antennas_fname @@ -44,14 +44,15 @@ if __name__ == "__main__": f_beacon, tx, antennas = beacon.read_beacon_hdf5(antennas_fname) # run over all baselines - if ref_ant_id is None: + if not ref_ant_id: print("Doing all baselines") - baselines = list(combinations(antennas,2)) + baselines = combinations(antennas,2) + # use ref_ant else: ref_ants = [antennas[i] for i in ref_ant_id] print("Doing all baselines with {}".format([int(a.name) for a in ref_ants])) - baselines = list(product(ref_ants, antennas)) + baselines = product(ref_ants, antennas) # For now, only one beacon_frequency is supported freq_names = antennas[0].beacon_info.keys() @@ -60,7 +61,11 @@ if __name__ == "__main__": freq_name = next(iter(freq_names)) + # Collect baselines from optional generators + baselines = list(baselines) + # Get phase difference per baseline + phase_diffs = np.empty( (len(baselines), 2) ) for i, base in enumerate(baselines): if i%1000==0: @@ -116,10 +121,10 @@ if __name__ == "__main__": if plot_residuals: phase_residuals = lib.phase_mod(phase_diffs[:,1] - my_phase_diffs) - fig.suptitle("Difference between Measured and Actual phase difference\n for Baselines (i,j" + (')' if ref_ant_id is None else '='+str([ int(a.name) for a in ref_ants])+')')) + fig.suptitle("Difference between Measured and Actual phase difference\n for Baselines (i,j" + (')' if not ref_ant_id else '='+str([ int(a.name) for a in ref_ants])+')')) axs[-1].set_xlabel("Baseline Phase Residual $\\Delta\\varphi_{ij_{meas}} - \\Delta\\varphi_{ij_{true}}$ [rad]") else: - fig.suptitle("Comparison Measured and Actual phase difference\n for Baselines (i,j" + (')' if ref_ant_id is None else '='+str([ int(a.name) for a in ref_ants])+')')) + fig.suptitle("Comparison Measured and Actual phase difference\n for Baselines (i,j" + (')' if not ref_ant_id else '='+str([ int(a.name) for a in ref_ants])+')')) axs[-1].set_xlabel("Baseline Phase $\\Delta\\varphi_{ij}$ [rad]") diff --git a/simulations/airshower_beacon_simulation/bd_antenna_phase_deltas.py b/simulations/airshower_beacon_simulation/bd_antenna_phase_deltas.py index 746ca9b..6abd475 100755 --- a/simulations/airshower_beacon_simulation/bd_antenna_phase_deltas.py +++ b/simulations/airshower_beacon_simulation/bd_antenna_phase_deltas.py @@ -25,14 +25,13 @@ if __name__ == "__main__": parser = MyArgumentParser() args = parser.parse_args() - fname = "ZH_airshower/mysim.sry" - figsize = (9,6) + figsize = (12,8) show_plots = args.show_plots ref_ant_id = None # leave None for all baselines #### - fname_dir = path.dirname(fname) + fname_dir = args.data_dir antennas_fname = path.join(fname_dir, beacon.antennas_fname) time_diffs_fname = 'time_diffs.hdf5' if False else antennas_fname fig_dir = args.fig_dir # set None to disable saving diff --git a/simulations/airshower_beacon_simulation/ca_period_from_shower.py b/simulations/airshower_beacon_simulation/ca_period_from_shower.py index 6f9e85f..30435c5 100755 --- a/simulations/airshower_beacon_simulation/ca_period_from_shower.py +++ b/simulations/airshower_beacon_simulation/ca_period_from_shower.py @@ -151,11 +151,24 @@ if __name__ == "__main__": from scriptlib import MyArgumentParser parser = MyArgumentParser(default_fig_dir="./figures/periods_from_shower_figures/") + parser.add_argument('--quick_run', action='store_true', help='Use a very coarse grid (6x6)') + parser.add_argument('--max-k', type=float, default=2, help='Maximum abs(k) allowed to be shifted. (Default: %(default)d)') + parser.add_argument('-N', '--N_runs', type=int, default=5, help='Maximum amount of iterations to grid search. (Default: %(default)d)') + + parser.add_argument('-l', '--passband-low', type=float, default=30e-3, help='Lower frequency [GHz] of the passband filter. (set -1 for np.inf) (Default: %(default)g)') + parser.add_argument('-u', '--passband-high', type=float, default=80e-3, help='Upper frequency [GHz] of the passband filter. (set -1 for np.inf) (Default: %(default)g)') + + parser.add_argument('--input-fname', type=str, default=None, help='Path to mysim.sry, either directory or path. If empty it takes DATA_DIR and appends mysim.sry. (Default: %(default)s)') args = parser.parse_args() - fname = "ZH_airshower/mysim.sry" - figsize = (9,6) + if not args.input_fname: + args.input_fname = args.data_dir + + if path.isdir(args.input_fname): + args.input_fname = path.join(args.input_fname, "mysim.sry") + + figsize = (12,8) fig_dir = args.fig_dir fig_subdir = path.join(fig_dir, 'shifts/') @@ -165,12 +178,15 @@ if __name__ == "__main__": allowed_ks = np.arange(-max_k, max_k+1, dtype=int) Xref = 400 - N_runs = 3 + N_runs = args.N_runs remove_beacon_from_trace = True apply_signal_window_from_max = True + low_bp = args.passband_low if args.passband_low >= 0 else np.inf # GHz + high_bp = args.passband_high if args.passband_high >= 0 else np.inf # GHz + #### - fname_dir = path.dirname(fname) + fname_dir = args.data_dir antennas_fname = path.join(fname_dir, beacon.antennas_fname) time_diffs_fname = 'time_diffs.hdf5' if not True else antennas_fname @@ -188,7 +204,7 @@ if __name__ == "__main__": # Read in antennas from file _, tx, antennas = beacon.read_beacon_hdf5(antennas_fname) # Read original REvent - ev = REvent(fname) + ev = REvent(args.input_fname) # .. patch in our antennas ev.antennas = antennas @@ -201,7 +217,7 @@ if __name__ == "__main__": f_beacon = ev.antennas[0].beacon_info[freq_name]['freq'] # Prepare polarisation and passbands - rit.set_pol_and_bp(ev, low=0.03, high=0.08) + rit.set_pol_and_bp(ev, low=low_bp, high=high_bp) ## ## Manipulate time and traces of each antenna @@ -297,14 +313,13 @@ if __name__ == "__main__": scale2d = dXref*np.tan(np.deg2rad(2.)) scale4d = dXref*np.tan(np.deg2rad(4.)) - if False: #quicky - x_coarse = np.linspace(-scale2d, scale2d, 4) - y_coarse = np.linspace(-scale2d, scale2d, 4) + if args.quick_run: #quicky + x_coarse = np.linspace(-scale2d, scale2d, 6) + y_coarse = np.linspace(-scale2d, scale2d, 6) x_fine = x_coarse/4 y_fine = y_coarse/4 else: # long - N_runs = 5 x_coarse = np.linspace(-scale4d, scale4d, 40) y_coarse = np.linspace(-scale4d, scale4d, 40) 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 index 6a00681..547216f 100755 --- a/simulations/airshower_beacon_simulation/cb_report_measured_antenna_time_offsets.py +++ b/simulations/airshower_beacon_simulation/cb_report_measured_antenna_time_offsets.py @@ -22,14 +22,13 @@ if __name__ == "__main__": parser = MyArgumentParser() args = parser.parse_args() - fname = "ZH_airshower/mysim.sry" - figsize = (9,6) + figsize = (12,8) fig_dir = args.fig_dir # set None to disable saving show_plots = args.show_plots #### - fname_dir = path.dirname(fname) + fname_dir = args.data_dir antennas_fname = path.join(fname_dir, beacon.antennas_fname) time_diffs_fname = 'time_diffs.hdf5' if not True else antennas_fname diff --git a/simulations/airshower_beacon_simulation/da_reconstruction.py b/simulations/airshower_beacon_simulation/da_reconstruction.py index 6bc4291..d85d764 100755 --- a/simulations/airshower_beacon_simulation/da_reconstruction.py +++ b/simulations/airshower_beacon_simulation/da_reconstruction.py @@ -31,10 +31,16 @@ if __name__ == "__main__": from scriptlib import MyArgumentParser parser = MyArgumentParser() + parser.add_argument('--input-fname', type=str, default=None, help='Path to mysim.sry, either directory or path. If empty it takes DATA_DIR and appends mysim.sry. (Default: %(default)s)') args = parser.parse_args() - fname = "ZH_airshower/mysim.sry" - figsize = (9,6) + if not args.input_fname: + args.input_fname = args.data_dir + + if path.isdir(args.input_fname): + args.input_fname = path.join(args.input_fname, "mysim.sry") + + figsize = (12,8) fig_dir = args.fig_dir fig_subdir = path.join(fig_dir, 'reconstruction') @@ -43,7 +49,7 @@ if __name__ == "__main__": remove_beacon_from_traces = True #### - fname_dir = path.dirname(fname) + fname_dir = args.data_dir antennas_fname = path.join(fname_dir, beacon.antennas_fname) pickle_fname = path.join(fname_dir, 'res.pkl') tx_fname = path.join(fname_dir, beacon.tx_fname) @@ -58,7 +64,7 @@ if __name__ == "__main__": _, tx, antennas = beacon.read_beacon_hdf5(antennas_fname) _, __, txdata = beacon.read_tx_file(tx_fname) # Read original REvent - ev = REvent(fname) + ev = REvent(args.input_fname) # .. patch in our antennas ev.antennas = antennas diff --git a/simulations/airshower_beacon_simulation/db_longitudinal_figure.py b/simulations/airshower_beacon_simulation/db_longitudinal_figure.py index 63f1da2..c5da800 100755 --- a/simulations/airshower_beacon_simulation/db_longitudinal_figure.py +++ b/simulations/airshower_beacon_simulation/db_longitudinal_figure.py @@ -28,13 +28,11 @@ if __name__ == "__main__": parser = MyArgumentParser() args = parser.parse_args() - fname = "ZH_airshower/mysim.sry" - fig_dir = args.fig_dir show_plots = args.show_plots #### - fname_dir = path.dirname(fname) + fname_dir = args.data_dir antennas_fname = path.join(fname_dir, beacon.antennas_fname) pickle_fname = path.join(fname_dir, 'res.pkl') diff --git a/simulations/airshower_beacon_simulation/dc_grid_power_time_fixes.py b/simulations/airshower_beacon_simulation/dc_grid_power_time_fixes.py index 585f080..210e158 100755 --- a/simulations/airshower_beacon_simulation/dc_grid_power_time_fixes.py +++ b/simulations/airshower_beacon_simulation/dc_grid_power_time_fixes.py @@ -32,6 +32,7 @@ if __name__ == "__main__": from scriptlib import MyArgumentParser parser = MyArgumentParser() + parser.add_argument('--input-fname', type=str, default=None, help='Path to mysim.sry, either directory or path. If empty it takes DATA_DIR and appends mysim.sry. (Default: %(default)s)') group = parser.add_argument_group('figures') for case in valid_cases: @@ -39,13 +40,18 @@ if __name__ == "__main__": args = parser.parse_args() + if not args.input_fname: + args.input_fname = args.data_dir + + if path.isdir(args.input_fname): + args.input_fname = path.join(args.input_fname, "mysim.sry") + wanted_cases = args.figures if not wanted_cases or 'all' in wanted_cases: wanted_cases = valid_cases - fname = "ZH_airshower/mysim.sry" - figsize = (9,6) + figsize = (12,8) fig_dir = args.fig_dir show_plots = args.show_plots @@ -54,7 +60,7 @@ if __name__ == "__main__": apply_signal_window_from_max = True #### - fname_dir = path.dirname(fname) + fname_dir = args.data_dir antennas_fname = path.join(fname_dir, beacon.antennas_fname) pickle_fname = path.join(fname_dir, 'res.pkl') tx_fname = path.join(fname_dir, beacon.tx_fname) @@ -67,7 +73,7 @@ if __name__ == "__main__": _, tx, antennas = beacon.read_beacon_hdf5(antennas_fname) _, __, txdata = beacon.read_tx_file(tx_fname) # Read original REvent - ev = REvent(fname) + ev = REvent(args.input_fname) bak_ants = ev.antennas # .. patch in our antennas ev.antennas = antennas diff --git a/simulations/airshower_beacon_simulation/figures/.gitignore b/simulations/airshower_beacon_simulation/figures/.gitignore deleted file mode 100644 index d6b7ef3..0000000 --- a/simulations/airshower_beacon_simulation/figures/.gitignore +++ /dev/null @@ -1,2 +0,0 @@ -* -!.gitignore diff --git a/simulations/airshower_beacon_simulation/figures/periods_from_shower_figures/shifts/.gitignore b/simulations/airshower_beacon_simulation/figures/periods_from_shower_figures/shifts/.gitignore deleted file mode 100644 index d6b7ef3..0000000 --- a/simulations/airshower_beacon_simulation/figures/periods_from_shower_figures/shifts/.gitignore +++ /dev/null @@ -1,2 +0,0 @@ -* -!.gitignore diff --git a/simulations/airshower_beacon_simulation/scriptlib.py b/simulations/airshower_beacon_simulation/scriptlib.py index 165f32f..90f0046 100644 --- a/simulations/airshower_beacon_simulation/scriptlib.py +++ b/simulations/airshower_beacon_simulation/scriptlib.py @@ -4,9 +4,14 @@ Some preconfigured ArgumentParser from argparse import ArgumentParser -def MyArgumentParser(default_fig_dir='./figures', default_show_plots=False, **kwargs): +def MyArgumentParser( + default_fig_dir='./figures', + default_show_plots=False, + default_data_dir='./ZH_airshower', + **kwargs): """ - A somewhat preconfigured ArgumentParser + A somewhat preconfigured ArgumentParser to be shared across + multiple scripts. Set show_plots=True to by default enable showing plots. Likewise, set fig_dir=None to by default disable saving figures. @@ -14,11 +19,14 @@ def MyArgumentParser(default_fig_dir='./figures', default_show_plots=False, **kw parser = ArgumentParser(**kwargs) # Whether to show plots - parser.add_argument('--show-plots', action='store_true') - parser.add_argument('--no-show-plots', dest='show-plots', action='store_false') - parser.set_defaults(show_plots=default_show_plots) + group1 = parser.add_mutually_exclusive_group(required=False) + group1.add_argument('--show-plots', action='store_true', default=default_show_plots, help='Default: %(default)s') + group1.add_argument('--no-show-plots', dest='show-plots', action='store_false') + + # Data directory + parser.add_argument('--data-dir', type=str, default=default_data_dir, help='Path to Data Directory. (Default: %(default)s)') # Figures directory - parser.add_argument('--fig-dir', type=str, default=default_fig_dir) + parser.add_argument('--fig-dir', type=str, default=default_fig_dir, help='Set None to disable saving figures. (Default: %(default)s)') return parser diff --git a/simulations/airshower_beacon_simulation/view_beaconed_antenna.py b/simulations/airshower_beacon_simulation/view_beaconed_antenna.py index 5af931e..22ec5e5 100755 --- a/simulations/airshower_beacon_simulation/view_beaconed_antenna.py +++ b/simulations/airshower_beacon_simulation/view_beaconed_antenna.py @@ -27,8 +27,7 @@ if __name__ == "__main__": args = parser.parse_args() - fname = "ZH_airshower/mysim.sry" - figsize = (9,6) + figsize = (12,8) plot_ft_amplitude = args.ft plot_geometry = args.geom @@ -40,7 +39,7 @@ if __name__ == "__main__": args.polarisations = ['x','y', 'z'] #### - fname_dir = path.dirname(fname) + fname_dir = args.data_dir antennas_fname = path.join(fname_dir, beacon.antennas_fname) tx_fname = path.join(fname_dir, beacon.tx_fname)