In [1]:
from __future__ import division
import os
import numpy as np
import matplotlib.pyplot as plt
import seaborn.apionly as sns
import comptools as comp
%matplotlib inline
In [2]:
config = 'IC86.2012'
num_groups = 2
comp_list = comp.get_comp_list(num_groups=num_groups)
In [3]:
energybins = comp.analysis.get_energybins()
In [4]:
df_sim = comp.load_sim(config=config, test_size=0, log_energy_min=6.0, log_energy_max=8.3)
In [5]:
reco_log_energy = df_sim.reco_log_energy.values
MC_log_energy = df_sim.MC_log_energy.values
In [6]:
reco_log_energy.min(), reco_log_energy.max()
Out[6]:
In [7]:
MC_log_energy.min(), MC_log_energy.max()
Out[7]:
In [8]:
h, xedges, yedges = np.histogram2d(MC_log_energy, reco_log_energy, bins=energybins.log_energy_bins)
extent = [xedges[0], xedges[-1], yedges[0], yedges[-1]]
h = h / h.sum(axis=0)
h = h.T
In [17]:
bins = energybins.log_energy_bins
bin_midpoints = energybins.log_energy_midpoints
fig, ax = plt.subplots(figsize=(10, 8))
im = ax.imshow(h, origin='lower', extent=extent)
# Current analysis range
# ax.axvline(6.4, marker='None', ls='-.', color='w')
# ax.axhline(6.4, marker='None', ls='-.', color='w')
# ax.axvline(7.9, marker='None', ls='-.', color='w')
# ax.axhline(7.9, marker='None', ls='-.', color='w')
# for i in bins:
# ax.axvline(i, marker='None', ls=':', color='w')
# ax.axhline(i, marker='None', ls=':', color='w')
ax.axvline(6.4, marker='None', ls='-.', color='C1')
ax.axhline(6.4, marker='None', ls='-.', color='C1')
# ax.axvline(7.8, marker='None', ls='-.', color='C1')
# ax.axhline(7.8, marker='None', ls='-.', color='C1')
# Add histogram value annotations
for i in range(len(xedges)-1):
for j in range(len(yedges)-1):
if h[j, i] != 0 and ~np.isnan(h)[j, i]:
ax.text(bin_midpoints[i], bin_midpoints[j], '{:0.2f}'.format(h[j, i]),
color='w', ha="center", va="center")
ax.set_xlabel('$\mathrm{\log_{10}(E_{true}/GeV)}$')
ax.set_ylabel('$\mathrm{\log_{10}(E_{reco}/GeV)}$')
ax.set_xlim(6.1, 8)
ax.set_ylim(6.1, 8)
# ax.text(6.375, 5.935, '6.4')
# ax.text(5.925, 6.375, '6.4')
plt.colorbar(im, label='$\mathrm{P(E_{true}|E_{reco})}$')
outfile = os.path.join(comp.paths.figures_dir, 'unfolding', 'analysis-bins.png')
plt.savefig(outfile)
plt.show()
In [ ]: