Merge branch 'more-makefile' into main

It allows to modify a config.mk to set variables across scripts.

simulations/airshower_beacon_simulation/.gitignore

@@ -0,0 +1,3 @@
+config.mk
+figures
+data

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
+	./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 | tee figures/ab.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}
+	./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 --no-show-plots --fig-dir=${FIG_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}
+	./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}
+	./bb_measure_clock_phase.py --no-show-plots --fig-dir=${FIG_DIR} --data-dir ${DATA_DIR}
-	./bc_baseline_phase_deltas.py --no-show-plots --fig-dir=${FIG_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}
+	./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}
+	./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}
+	./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}
+	./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}
+	./db_longitudinal_figure.py --no-show-plots --fig-dir=${FIG_DIR} --data-dir ${DATA_DIR}
-	./dc_grid_power_time_fixes.py --no-show-plots --fig-dir=${FIG_DIR}
+
 dist-clean:
-	rm -f ZH_airshower/antennas.hdf5
+	rm -vf ${DATA_DIR}/antennas.hdf5
-	rm -f ZH_airshower/clocks.csv
+	rm -vf ${DATA_DIR}/ca_breaked_run
-	rm -f ZH_airshower/tx.json
+	rm -vf ${DATA_DIR}/res.pkl
+	rm -vf ${DATA_DIR}/clocks.csv
+	rm -vf ${DATA_DIR}/tx.json

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('-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)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)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

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)
 

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 = (12,8)
-    figsize = (9,6)
 
     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)
 

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_only_AxB_trace = args.use_AxB_trace
-    use_beacon_trace = True # only applicable if AxB = False
+    use_only_beacon_trace = args.use_beacon_trace # only applicable if AxB = False
 
     show_plots = args.show_plots
 
-    fname = "ZH_airshower/mysim.sry"
+    figsize = (12,8)
-    figsize = (9,6)
 
-    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 use_only_AxB_trace:
-                if 'E_AxB' not in h5ant.keys():
+                traces_key = 'E_AxB'
-                    print(f"Antenna does not have 'E_AxB' in {name}")
+                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
+            # Only beacon
-                # use beacon but remove where E_AxB-Beacon != 0
+            elif use_only_beacon_trace:
-                if True:
+                traces_key = 'filtered_traces'
-                    N_pre, N_post = 250, 250 # TODO: make this configurable
+                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])
+                traces = h5ant[traces_key]
-
+                t_trace = traces[0]
-                    low_idx = max(0, max_idx-N_pre)
+                test_traces = [traces[4]]
-                    high_idx = min(len(t_trace), max_idx+N_post)
+                orients = ['B']
-
-                    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'
 
             # use separate polarisations
             else:
-                if 'traces' not in h5ant.keys():
+                traces_key = 'filtered_traces'
-                    print(f"Antenna file corrupted? no 'traces' in {name}")
+                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:
+            # Really only select the first component
-                    # only take the Beacon trace
+            if True:
-                    test_traces = [traces[4]]
+                test_traces = [test_traces[0]]
-                    orients = ['B']
+                orients = [orients[0]]
-                else:
-                    test_traces = traces[1:]
-                    orients = ['Ex', 'Ey', 'Ez', 'B']
 
-            # modify the length of the traces
+            # TODO: refine masking
-            if False:
+            # use beacon but remove where E_AxB-Beacon != 0
-                t_trace = t_trace[:len(t_trace)//2]
+            # Uses the first traces as reference
-                half_traces = []
+            if N_mask and orients[0] != 'B':
-                for trace in test_traces:
+                N_pre, N_post = N_mask//2, N_mask//2
-                    half_traces.append( trace[:len(trace)//2])
+
-                test_traces = half_traces
+                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)

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 = (12,8)
-    figsize = (9,6)
     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

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 = (12,8)
-    figsize = (9,6)
     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]")
 
 

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 = (12,8)
-    figsize = (9,6)
 
     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

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"
+    if not args.input_fname:
-    figsize = (9,6)
+        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
+    if args.quick_run: #quicky
-        x_coarse = np.linspace(-scale2d, scale2d, 4)
+        x_coarse = np.linspace(-scale2d, scale2d, 6)
-        y_coarse = np.linspace(-scale2d, scale2d, 4)
+        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)
 

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 = (12,8)
-    figsize = (9,6)
 
     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
 

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"
+    if not args.input_fname:
-    figsize = (9,6)
+        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
 

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')
 

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 = (12,8)
-    figsize = (9,6)
 
     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

simulations/airshower_beacon_simulation/figures/.gitignore

@@ -1,2 +0,0 @@
-*
-!.gitignore

simulations/airshower_beacon_simulation/figures/periods_from_shower_figures/shifts/.gitignore

@@ -1,2 +0,0 @@
-*
-!.gitignore

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')
+    group1 = parser.add_mutually_exclusive_group(required=False)
-    parser.add_argument('--no-show-plots', dest='show-plots', action='store_false')
+    group1.add_argument('--show-plots', action='store_true', default=default_show_plots, help='Default: %(default)s')
-    parser.set_defaults(show_plots=default_show_plots)
+    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

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 = (12,8)
-    figsize = (9,6)
 
     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)