Stabilize muon intensity fit

ctapipe/image/muon/intensity_fitter.py

@@ -413,7 +413,11 @@ def negative_log_likelihood(
         # scale prediction by optical efficiency of the telescope
         prediction *= optical_efficiency_muon
 
-        return neg_log_likelihood_approx(image, prediction, spe_width, pedestal)
+        # neg_log_likelihood_approx is build choosen to reseble a chi2,
+        # use chi2 per dof for numerical stability in minimization later
+        return neg_log_likelihood_approx(image, prediction, spe_width, pedestal) / len(
+            image
+        )
 
     return negative_log_likelihood
 

ctapipe/image/pixel_likelihood.py

@@ -43,6 +43,7 @@
 ]
 
 EPSILON = 5.0e-324
+LOG_2_PI = np.log(2 * np.pi)
 
 
 class PixelLikelihoodError(RuntimeError):
@@ -72,8 +73,8 @@ def neg_log_likelihood_approx(image, prediction, spe_width, pedestal):
 
         - \\ln{P} = \\frac{\\ln{2 π} + \\ln{θ}}{2} + \\frac{(s - μ)^2}{2 θ}
 
-        We keep the constants in this because the actual value of the likelihood
-        can be used to calculate a goodness-of-fit value
+    We keep the constants in this because the actual value of the likelihood
+    can be used to calculate a goodness-of-fit value
 
 
     Parameters
@@ -89,15 +90,16 @@ def neg_log_likelihood_approx(image, prediction, spe_width, pedestal):
 
     Returns
     -------
-    float
+    result: float
+        Negative log likelihood value.
     """
     theta = pedestal**2 + prediction * (1 + spe_width**2)
 
-    neg_log_l = 0.5 * (
-        np.log(2 * np.pi * theta + EPSILON) + (image - prediction) ** 2 / theta
-    )
+    neg_log_l = np.log(theta + EPSILON) + (image - prediction) ** 2 / theta
 
-    return np.sum(neg_log_l)
+    # neg_log_l provides the variable term, add constants here to only
+    # compute them once
+    return 0.5 * (len(neg_log_l) * LOG_2_PI + np.sum(neg_log_l))
 
 
 def neg_log_likelihood_numeric(

docs/changes/2425.bugfix.rst

@@ -0,0 +1,4 @@
+Fixes a bug (see #2415), where the test for the ``MuonIntensityFitter`` would fail as the ``Minuit`` result
+was considered as not valid due to not satisfiying the ``is_above_max_edm`` criterion.
+Using a log-likelihood reduced by the number of dof lowered the likelihood-values by several 
+orders of magnitude, resulting in a valid fit-result.