From ef95114fe81e5f9c53e80a35fc96b7cb2747d3ee Mon Sep 17 00:00:00 2001
From: Pim Nelissen <pim.nelissen@protonmail.com>
Date: Mon, 27 Oct 2025 15:40:12 +0100
Subject: [PATCH] update histogram style

---
 I.py | 21 ++++++++++++++-------
 1 file changed, 14 insertions(+), 7 deletions(-)

diff --git a/I.py b/I.py
index c9bda0d..c0f152c 100644
--- a/I.py
+++ b/I.py
@@ -20,7 +20,8 @@ else:
         lengths[i] = band.length/band.a
 
 # histogram
-values, bin_edges = np.histogram(lengths, bins=int(2*N), range=(-N, N))
+bin_edges = np.arange(-N-1, N+1, 2)
+values, _ = np.histogram(lengths, bins=bin_edges)
 p_mc = values / values.sum() # normalise histogram
 
 bin_centers = bin_edges[:-1] + np.diff(bin_edges)/2
@@ -33,14 +34,21 @@ def P(N, L_over_a):
     """
     return scp.special.binom(N, (N+L_over_a)/2) / 2**N
 
-p_true = P(N, bin_centers-.5)
+# Calculate p_true, scale with NUM_BANDS,
+# mask low statistics bins,
+# calculate chi^2/ndf
 
-chi2 = np.sum((p_mc - p_true)**2 / p_true)
+p_true = P(N, bin_centers)
+p_true_scaled = p_true*np.sum(values)
+values_masked = values[values > 5]
+p_true_scaled_masked = p_true_scaled[values > 5]
+chi2 = np.sum((values_masked - p_true_scaled_masked)**2 / p_true_scaled_masked)
+chi2_over_ndf =  chi2 / (len(values_masked)-1)
 
 # plot hist vs true dist.
-plt.plot(bin_centers-.5, p_true, label="$P(L)$")
-plt.step(bin_centers, p_mc, color='r', alpha = 0.5, label="$\hat{P}(L)$")
-plt.text(-45, 0.07, f"$\chi^2 = {round(chi2, 3)}$\n({len(lengths):.0E} samples)")
+plt.plot(bin_centers, p_true, label="$P(L)$")
+plt.bar(bin_centers, p_mc, color='r', alpha = 0.5, label="$\hat{P}(L)$")
+plt.text(-45, 0.07, f"$\chi^2/ndf = {round(chi2_over_ndf, 3)}$\n({len(lengths):.0E} samples)")
 plt.legend()
 plt.xlabel("$L/a$")
 plt.ylabel("P($L/a$)")
@@ -53,7 +61,6 @@ plt.figure()
 plt.scatter(bin_centers, p_ratio, marker='x')
 plt.hlines(1, *xlim, color='r', linestyle='--')
 plt.xlim(xlim)
-plt.ylim(0.5, 1.5)
 plt.xlabel("$L/a$")
 plt.ylabel("$\hat{P}(L)/P(L)$")
 plt.show()