From ab07fe2126c51758e45ec0a5fed15e612c84eac8 Mon Sep 17 00:00:00 2001 From: Eric Teunis de Boone Date: Mon, 9 Jan 2023 11:51:03 +0100 Subject: [PATCH] ZH: mutliple k script saves results into hdf5 file + slightly improve power per antenna plot if tx is plotted --- .../ca_period_from_shower.py | 21 ++++++++++++++++++- 1 file changed, 20 insertions(+), 1 deletion(-) diff --git a/simulations/airshower_beacon_simulation/ca_period_from_shower.py b/simulations/airshower_beacon_simulation/ca_period_from_shower.py index 6ca54df..d532cc8 100755 --- a/simulations/airshower_beacon_simulation/ca_period_from_shower.py +++ b/simulations/airshower_beacon_simulation/ca_period_from_shower.py @@ -311,20 +311,27 @@ if __name__ == "__main__": axs.set_title(f"Grid Run {r}") axs.set_ylabel("vxvxB [km]") axs.set_xlabel(" vxB [km]") + axs.set_aspect('equal', 'datalim') sc = axs.scatter(xx/1e3, yy/1e3, c=maxima_per_loc, cmap='Spectral_r', alpha=0.6) fig.colorbar(sc, ax=axs) if fig_dir: + old_xlims = axs.get_xlim() + old_ylims = axs.get_ylim() fig.tight_layout() fig.savefig(path.join(fig_dir, __file__+f'.maxima.run{r}.pdf')) if True: axs.plot(tx.x/1e3, tx.y/1e3, marker='X', color='k') fig.tight_layout() fig.savefig(path.join(fig_dir, __file__+f'.maxima.run{r}.with_tx.pdf')) + axs.set_xlim(*old_xlims) + axs.set_ylim(*old_ylims) + fig.tight_layout() ## Save ks to file best_idx = np.argmax(maxima_per_loc) np.savetxt(fig_dir + __file__+f'.bestk.run{r}.txt', ks_per_loc[best_idx] ) + print("Max at location: ", locs[best_idx]) print('Best k:', ks_per_loc[best_idx]) # Abort if no improvement @@ -332,12 +339,24 @@ if __name__ == "__main__": print("No changes from previous grid, breaking") break + old_ks_per_loc = ks_per_loc[best_idx] + # Save best ks to hdf5 antenna file with h5py.File(antennas_fname, 'a') as fp: group = fp.require_group('antennas') for i, ant in enumerate(antennas): h5ant = group[ant.name] - h5ant.attrs['best_k'] = old_ks_per_loc[i] + + h5beacon_info = h5ant['beacon_info'] + + # find out freq_name + if freq_name is None: + freq_name = [ k for k in h5beacon_info.keys() if np.isclose(h5beacon_info[k].attrs['freq'], f_beacon)][0] + + h5attrs = h5beacon_info[freq_name].attrs + + h5attrs['best_k'] = old_ks_per_loc[i] + h5attrs['best_k_time'] = old_ks_per_loc[i]*dt/f_beacon plt.show()