Add algorithm to compute image parameters via fitting the light distribution

ctapipe/containers.py

@@ -34,6 +34,7 @@
     "MorphologyContainer",
     "BaseHillasParametersContainer",
     "CameraHillasParametersContainer",
+    "ImageFitParametersContainer",
     "CameraTimingParametersContainer",
     "ParticleClassificationContainer",
     "PedestalContainer",
@@ -308,6 +309,37 @@ class HillasParametersContainer(BaseHillasParametersContainer):
     psi = Field(nan * u.deg, "rotation angle of ellipse", unit=u.deg)
 
 
+class ImageFitParametersContainer(HillasParametersContainer):
+    """
+    Hillas Parameters after fitting. The cog position
+    is given in meter from the camera center.
+    """
+
+    skewness_uncertainty = Field(nan, "measure of skewness uncertainty")
+    likelihood = Field(nan, "measure of likelihood")
+    goodness_of_fit = Field(
+        nan, "measure of goodness of fit, mean likelihood subtracted to the likelihood"
+    )
+    n_pixels = Field(-1, "number of pixels used in the fit")
+    free_parameters = Field(-1, "number of free parameters")
+    is_valid = Field(False, "True if the fit is valid")
+    is_accurate = Field(
+        False, "returns True if the fit is accurate. If False, the fit is not reliable."
+    )
+
+    fov_lon_uncertainty = Field(nan * u.deg, "centroid x uncertainty", unit=u.deg)
+    fov_lat_uncertainty = Field(nan * u.deg, "centroid y uncertainty", unit=u.deg)
+    r_uncertainty = Field(nan * u.deg, "centroid r uncertainty", unit=u.deg)
+    phi_uncertainty = Field(
+        nan * u.deg, "polar coordinate of centroid uncertainty", unit=u.deg
+    )
+    psi_uncertainty = Field(
+        nan * u.deg, "Uncertainty in rotation angle of ellipse", unit=u.deg
+    )
+    amplitude = Field(nan, "Amplitude of the fitted model")
+    amplitude_uncertainty = Field(nan, "error in amplitude from the fit")
+
+
 class LeakageContainer(Container):
     """
     Fraction of signal in 1 or 2-pixel width border from the edge of the

ctapipe/image/concentration.py

@@ -1,14 +1,14 @@
-import numpy as np
 import astropy.units as u
+import numpy as np
 
 from ..containers import (
+    CameraHillasParametersContainer,
     ConcentrationContainer,
     HillasParametersContainer,
-    CameraHillasParametersContainer,
 )
-from .hillas import camera_to_shower_coordinates
 from ..instrument import CameraGeometry
 from ..utils.quantities import all_to_value
+from .hillas import camera_to_shower_coordinates
 
 __all__ = ["concentration_parameters"]
 
@@ -53,15 +53,15 @@ def concentration_parameters(geom: CameraGeometry, image, hillas_parameters):
     delta_y = pix_y - y
 
     # take pixels within one pixel diameter from the cog
-    mask_cog = (delta_x ** 2 + delta_y ** 2) < pixel_width ** 2
+    mask_cog = (delta_x**2 + delta_y**2) < pixel_width**2
     conc_cog = np.sum(image[mask_cog]) / h.intensity
 
     if hillas_parameters.width.value != 0:
         # get all pixels inside the hillas ellipse
         longi, trans = camera_to_shower_coordinates(
             pix_x, pix_y, x, y, h.psi.to_value(u.rad)
         )
-        mask_core = (longi ** 2 / length ** 2) + (trans ** 2 / width ** 2) <= 1.0
+        mask_core = (longi**2 / length**2) + (trans**2 / width**2) <= 1.0
         conc_core = image[mask_core].sum() / h.intensity
     else:
         conc_core = 0.0