simplify likelihood calculation in chi2

1 files changed, 5 insertions, 8 deletions

diff --git a/utils/chi2 b/utils/chi2
index 4515b91..3b519e2 100755
--- a/utils/chi2
+++ b/utils/chi2
@@ -110,14 +110,11 @@ def make_nll(data, mc_hists):
 
         nll = 0
         for id in data_hists:
-            N = data_hists[id].sum()
-            nll -= poisson.logpmf(N,mc_hists[id].sum()*x[0])
-            if N > 0:
-                p = mc_hists[id]/mc_hists[id].sum()
-                # Fix a bug in scipy(). See https://github.com/scipy/scipy/issues/8235 (I think).
-                p += 1e-10
-                p /= p.sum()
-                nll -= multinomial.logpmf(data_hists[id],N,p)
+            oi = data_hists[id].sum()
+            ei = mc_hists[id] + EPSILON
+            N = oi.sum()
+            n = ei.sum()
+            nll -= -N - np.sum(gammaln(oi+1)) + np.sum(oi*np.log(ei))
         return nll - norm.logpdf(x[1],ENERGY_SCALE_MEAN,ENERGY_SCALE_UNCERTAINTY)
     return nll