diff --git a/simulations/11_pulsed_timing.py b/simulations/11_pulsed_timing.py index f762011..0876428 100755 --- a/simulations/11_pulsed_timing.py +++ b/simulations/11_pulsed_timing.py @@ -625,8 +625,8 @@ if __name__ == "__main__": # SNR time accuracy plot # if True: - threshold_markers = ['^', 'v', '8', 'o'] - mask_thresholds = [np.inf, N_residuals*0.5, N_residuals*0.1, 1, 0] + threshold_markers = ['^', 'v', '8', 'X'] # make sure to have filled markers here + mask_thresholds = np.array([np.inf, N_residuals*0.5, N_residuals*0.1, 1, 0]) fig, ax = plt.subplots() ax.set_title(f"Template matching SNR vs time accuracy") @@ -640,30 +640,66 @@ if __name__ == "__main__": f"antenna_dt={antenna_dt:0.1e}ns", ])) - if True: - ax.set_xscale('log') - ax.set_yscale('log') - # plot the values per template_dt slice template_dt_colors = [None]*len(template_dts) - for k, template_dt in enumerate(template_dts): + for a, template_dt in enumerate(template_dts): + for k, snr_sigma_factor in enumerate(snr_factors): + time_residuals, snrs, valid_mask = time_residuals_data[a][k] - # indicate masking values - for j, mask_threshold in enumerate(pairwise(mask_thresholds)): + valid_mask = np.array(valid_mask, dtype=bool) + + mean_residual = np.mean(time_residuals[valid_mask]) + time_accuracy = np.std(time_residuals[valid_mask]) + + residual_mean_deviation = np.sqrt( (time_residuals - mean_residual)**2 ) + + scatter_kwargs = dict( + ls='none', + marker='.', + alpha=0.3, + zorder=1.8, + ) + + y_values = residual_mean_deviation + + # snr_sigma_factor is a factor 2 too low + snr_sigma_factor *= 2 + + # plot all invalid datapoints + if True: + ax.plot(snrs[~valid_mask], y_values[~valid_mask], color='grey', **scatter_kwargs) + + # plot valid datapoints + if True: + if template_dt_colors[a] is not None: + scatter_kwargs['color'] = template_dt_colors[a] + + l = ax.plot(snrs[valid_mask], y_values[valid_mask], **scatter_kwargs) + + template_dt_colors[a] = l[0].get_color() + + masked_count = np.count_nonzero(~valid_mask) + + # plot accuracy indicating masking counts kwargs = dict( ls='none', - marker=threshold_markers[j], - color= None if template_dt_colors[k] is None else template_dt_colors[k] - ) - mask = mask_counts[k] >= mask_threshold[1] - mask &= mask_counts[k] < mask_threshold[0] + color= None if template_dt_colors[a] is None else template_dt_colors[a], + marker=threshold_markers[np.argmin( masked_count <= mask_thresholds)-1], + ms=10, + markeredgecolor='white', + markeredgewidth=1, + ) - l = ax.plot(snr_factors[mask], time_accuracies[k][mask], **kwargs) - template_dt_colors[k] = l[0].get_color() + #l = ax.plot(snr_sigma_factor, np.sqrt(np.mean(y_values[valid_mask])**2), **{**kwargs, **dict(ms=50)}) + l = ax.plot(snr_sigma_factor, np.std(time_residuals[valid_mask]), **kwargs) - # indicate threshold + # set color if not yet set + template_dt_colors[a] = l[0].get_color() + + + # indicate boxcar threshold if True: - ax.axhline(template_dt/np.sqrt(12), ls='--', alpha=0.7, color=template_dt_colors[k], label=f'Template dt:{template_dt:0.1e}ns') + ax.axhline(template_dt/np.sqrt(12), ls='--', alpha=0.7, color=template_dt_colors[a], label=f'Template dt:{template_dt:0.1e}ns') # Set horizontal line at 1 ns @@ -672,8 +708,32 @@ if __name__ == "__main__": ax.legend() - if True: # limit y-axis to 1e1 - ax.set_ylim([None, 1e1]) + fig.tight_layout() + fig.savefig(f"figures/11_time_res_vs_snr_full_linear.pdf") + + # logscaling + if True: + ax.set_xscale('log') + ax.set_yscale('log') + + # limit y-axis upper limit to 1e1 + if True: + this_lim = 1e1 + ax.set_ylim([None, this_lim]) + + # require y-axis lower limit to be at least 1e-1 + if True: + this_lim = 1e-1 + low_ylims = ax.get_ylim()[0] + if low_ylims >= this_lim: + ax.set_ylim([this_lim, None]) + + # .. but keep it above 1e-3 + if True: + this_lim = 1e-3 + low_ylims = ax.get_ylim()[0] + if low_ylims <= this_lim: + ax.set_ylim([this_lim, None]) if True: # require y-axis lower limit to be at least 1e-1 low_ylims = ax.get_ylim()[0]