From 3834b66faa510239f2737f9ecfeedfc97a5f4fc3 Mon Sep 17 00:00:00 2001 From: Vadym Voitsekhovskyi Date: Tue, 11 Jul 2023 15:52:26 +0200 Subject: [PATCH 01/10] Update intensity_fitter.py Added variable to check if fit converged, and storing of this variable in resulting container --- ctapipe/image/muon/intensity_fitter.py | 4 ++++ 1 file changed, 4 insertions(+) diff --git a/ctapipe/image/muon/intensity_fitter.py b/ctapipe/image/muon/intensity_fitter.py index 1eb40309402..7f6fde70f95 100644 --- a/ctapipe/image/muon/intensity_fitter.py +++ b/ctapipe/image/muon/intensity_fitter.py @@ -536,6 +536,9 @@ def __call__(self, tel_id, center_x, center_y, radius, image, pedestal, mask=Non # Perform minimisation minuit.migrad() + # Check for convergence + validation = minuit.valid + # Get fitted values result = minuit.values @@ -545,4 +548,5 @@ def __call__(self, tel_id, center_x, center_y, radius, image, pedestal, mask=Non impact_y=result["impact_parameter"] * np.sin(result["phi"]) * u.m, width=u.Quantity(np.rad2deg(result["ring_width"]), u.deg), optical_efficiency=result["optical_efficiency_muon"], + fit_convergence=validation ) From ecc22c4c968fcc124b621572f0ee56f18be6cd3d Mon Sep 17 00:00:00 2001 From: Vadym Voitsekhovskyi Date: Tue, 11 Jul 2023 15:54:16 +0200 Subject: [PATCH 02/10] Update containers.py Added new field to store an information about convergence of the muon intensity fit --- ctapipe/containers.py | 1 + 1 file changed, 1 insertion(+) diff --git a/ctapipe/containers.py b/ctapipe/containers.py index bcccaef7f40..eb0dea53355 100644 --- a/ctapipe/containers.py +++ b/ctapipe/containers.py @@ -984,6 +984,7 @@ class MuonEfficiencyContainer(Container): impact_x = Field(nan * u.m, "impact parameter x position") impact_y = Field(nan * u.m, "impact parameter y position") optical_efficiency = Field(nan, "optical efficiency muon") + fit_convergence = Field(nan, "convergence of the fit") class MuonParametersContainer(Container): From 98bad8c0bca9d03defdab68f547b617f166cb744 Mon Sep 17 00:00:00 2001 From: Vadym Voitsekhovskyi Date: Tue, 11 Jul 2023 17:05:23 +0200 Subject: [PATCH 03/10] Update containers.py Implemented comments of @kosack and @maxnoe --- ctapipe/containers.py | 4 +++- 1 file changed, 3 insertions(+), 1 deletion(-) diff --git a/ctapipe/containers.py b/ctapipe/containers.py index eb0dea53355..b35f4f3c6b7 100644 --- a/ctapipe/containers.py +++ b/ctapipe/containers.py @@ -984,7 +984,9 @@ class MuonEfficiencyContainer(Container): impact_x = Field(nan * u.m, "impact parameter x position") impact_y = Field(nan * u.m, "impact parameter y position") optical_efficiency = Field(nan, "optical efficiency muon") - fit_convergence = Field(nan, "convergence of the fit") + is_valid = Field(False, "True if the fit converged successfully") + parameters_at_limit = Field(False, "True if any bounded parameter was fitted close to a bound") + likelihood_value = Field(nan, "cost function value at the minimum") class MuonParametersContainer(Container): From 26280e77ba3aa2944173a6735d3faa2e900be555 Mon Sep 17 00:00:00 2001 From: Vadym Voitsekhovskyi Date: Tue, 11 Jul 2023 17:46:36 +0200 Subject: [PATCH 04/10] Update intensity_fitter.py Added assignment of new fields --- ctapipe/image/muon/intensity_fitter.py | 4 +++- 1 file changed, 3 insertions(+), 1 deletion(-) diff --git a/ctapipe/image/muon/intensity_fitter.py b/ctapipe/image/muon/intensity_fitter.py index 7f6fde70f95..7ac332dee5a 100644 --- a/ctapipe/image/muon/intensity_fitter.py +++ b/ctapipe/image/muon/intensity_fitter.py @@ -548,5 +548,7 @@ def __call__(self, tel_id, center_x, center_y, radius, image, pedestal, mask=Non impact_y=result["impact_parameter"] * np.sin(result["phi"]) * u.m, width=u.Quantity(np.rad2deg(result["ring_width"]), u.deg), optical_efficiency=result["optical_efficiency_muon"], - fit_convergence=validation + is_valid=minuit.valid + parameters_at_limit=minuit.fmin.has_parameters_at_limit + likelihood_value=minuit.fval ) From 2e25d02112774b1bc4a39f465e107b535b160550 Mon Sep 17 00:00:00 2001 From: Vadym Voitsekhovskyi Date: Tue, 11 Jul 2023 17:47:16 +0200 Subject: [PATCH 05/10] Update intensity_fitter.py Deleted unnecessary variable --- ctapipe/image/muon/intensity_fitter.py | 3 --- 1 file changed, 3 deletions(-) diff --git a/ctapipe/image/muon/intensity_fitter.py b/ctapipe/image/muon/intensity_fitter.py index 7ac332dee5a..f3e3c7a1968 100644 --- a/ctapipe/image/muon/intensity_fitter.py +++ b/ctapipe/image/muon/intensity_fitter.py @@ -535,9 +535,6 @@ def __call__(self, tel_id, center_x, center_y, radius, image, pedestal, mask=Non # Perform minimisation minuit.migrad() - - # Check for convergence - validation = minuit.valid # Get fitted values result = minuit.values From 78507435a17fc06c139870de0f0bc4c579fcb5cd Mon Sep 17 00:00:00 2001 From: Vadym Voitsekhovskyi Date: Wed, 12 Jul 2023 10:28:13 +0200 Subject: [PATCH 06/10] Update intensity_fitter.py Added commas --- ctapipe/image/muon/intensity_fitter.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/ctapipe/image/muon/intensity_fitter.py b/ctapipe/image/muon/intensity_fitter.py index f3e3c7a1968..b4df7bd5910 100644 --- a/ctapipe/image/muon/intensity_fitter.py +++ b/ctapipe/image/muon/intensity_fitter.py @@ -545,7 +545,7 @@ def __call__(self, tel_id, center_x, center_y, radius, image, pedestal, mask=Non impact_y=result["impact_parameter"] * np.sin(result["phi"]) * u.m, width=u.Quantity(np.rad2deg(result["ring_width"]), u.deg), optical_efficiency=result["optical_efficiency_muon"], - is_valid=minuit.valid - parameters_at_limit=minuit.fmin.has_parameters_at_limit + is_valid=minuit.valid, + parameters_at_limit=minuit.fmin.has_parameters_at_limit, likelihood_value=minuit.fval ) From 89ca045b4a7da6e37b503489f06f0bf09abf6a20 Mon Sep 17 00:00:00 2001 From: Vadym Voitsekhovskyi Date: Wed, 12 Jul 2023 13:02:33 +0200 Subject: [PATCH 07/10] Added tests for the new fields --- ctapipe/image/muon/tests/test_intensity_fit.py | 3 +++ 1 file changed, 3 insertions(+) diff --git a/ctapipe/image/muon/tests/test_intensity_fit.py b/ctapipe/image/muon/tests/test_intensity_fit.py index a534da44aa0..25f8a98f620 100644 --- a/ctapipe/image/muon/tests/test_intensity_fit.py +++ b/ctapipe/image/muon/tests/test_intensity_fit.py @@ -89,6 +89,9 @@ def test_muon_efficiency_fit(prod5_lst): assert u.isclose(result.impact, impact_parameter, rtol=0.05) assert u.isclose(result.width, ring_width, rtol=0.05) assert u.isclose(result.optical_efficiency, efficiency, rtol=0.05) + assert result.is_valid + assert not result.parameters_at_limit + assert np.isfinite(result.likelihood_value) def test_scts(prod5_sst): From ba5c8159e3bb37ee9ec0f0f36d70b72e21998827 Mon Sep 17 00:00:00 2001 From: Vadym Voitsekhovskyi Date: Wed, 12 Jul 2023 14:39:39 +0200 Subject: [PATCH 08/10] Created file with description of channges --- docs/changes/2381.feature.rst | 1 + 1 file changed, 1 insertion(+) create mode 100644 docs/changes/2381.feature.rst diff --git a/docs/changes/2381.feature.rst b/docs/changes/2381.feature.rst new file mode 100644 index 00000000000..402e78f05a3 --- /dev/null +++ b/docs/changes/2381.feature.rst @@ -0,0 +1 @@ +Added new fields to the `MuonEfficiencyContainer` - `is_valid` to check if fit converged successfully, `parameters_at_limit` to check if parameters were fitted close to a bound and `likelihood_value` which represents cost function value atthe minimum. These fields were added to the output of the `MuonIntensityFitter` From 2cca40ab5c1ea3360b1776a59946952769323de5 Mon Sep 17 00:00:00 2001 From: Vadym Voitsekhovskyi Date: Wed, 12 Jul 2023 15:08:56 +0200 Subject: [PATCH 09/10] Replaced single ticks with double tickks --- docs/changes/2381.feature.rst | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/docs/changes/2381.feature.rst b/docs/changes/2381.feature.rst index 402e78f05a3..4d9e9177c9d 100644 --- a/docs/changes/2381.feature.rst +++ b/docs/changes/2381.feature.rst @@ -1 +1 @@ -Added new fields to the `MuonEfficiencyContainer` - `is_valid` to check if fit converged successfully, `parameters_at_limit` to check if parameters were fitted close to a bound and `likelihood_value` which represents cost function value atthe minimum. These fields were added to the output of the `MuonIntensityFitter` +Added new fields to the ``MuonEfficiencyContainer`` - ``is_valid`` to check if fit converged successfully, ``parameters_at_limit`` to check if parameters were fitted close to a bound and ``likelihood_value`` which represents cost function value atthe minimum. These fields were added to the output of the ``MuonIntensityFitter``. From 5ff3637d7c136adb7813ea80d2910d6b35566de1 Mon Sep 17 00:00:00 2001 From: Vadym Voitsekhovskyi Date: Wed, 12 Jul 2023 17:18:54 +0200 Subject: [PATCH 10/10] Resolved style issues --- ctapipe/containers.py | 4 +++- ctapipe/image/muon/intensity_fitter.py | 4 ++-- ctapipe/image/toymodel.py | 24 ++++++++++----------- ctapipe/visualization/bokeh.py | 30 ++++++++++++-------------- 4 files changed, 31 insertions(+), 31 deletions(-) diff --git a/ctapipe/containers.py b/ctapipe/containers.py index b35f4f3c6b7..060bbc903d6 100644 --- a/ctapipe/containers.py +++ b/ctapipe/containers.py @@ -985,7 +985,9 @@ class MuonEfficiencyContainer(Container): impact_y = Field(nan * u.m, "impact parameter y position") optical_efficiency = Field(nan, "optical efficiency muon") is_valid = Field(False, "True if the fit converged successfully") - parameters_at_limit = Field(False, "True if any bounded parameter was fitted close to a bound") + parameters_at_limit = Field( + False, "True if any bounded parameter was fitted close to a bound" + ) likelihood_value = Field(nan, "cost function value at the minimum") diff --git a/ctapipe/image/muon/intensity_fitter.py b/ctapipe/image/muon/intensity_fitter.py index b4df7bd5910..97046316607 100644 --- a/ctapipe/image/muon/intensity_fitter.py +++ b/ctapipe/image/muon/intensity_fitter.py @@ -535,7 +535,7 @@ def __call__(self, tel_id, center_x, center_y, radius, image, pedestal, mask=Non # Perform minimisation minuit.migrad() - + # Get fitted values result = minuit.values @@ -547,5 +547,5 @@ def __call__(self, tel_id, center_x, center_y, radius, image, pedestal, mask=Non optical_efficiency=result["optical_efficiency_muon"], is_valid=minuit.valid, parameters_at_limit=minuit.fmin.has_parameters_at_limit, - likelihood_value=minuit.fval + likelihood_value=minuit.fval, ) diff --git a/ctapipe/image/toymodel.py b/ctapipe/image/toymodel.py index bb9a9f5abf6..abad3a66199 100644 --- a/ctapipe/image/toymodel.py +++ b/ctapipe/image/toymodel.py @@ -18,15 +18,17 @@ >>> print(image.shape) (400,) """ -import numpy as np -from ctapipe.utils import linalg -from ctapipe.image.hillas import camera_to_shower_coordinates +from abc import ABCMeta, abstractmethod + import astropy.units as u +import numpy as np from astropy.coordinates import Angle -from scipy.stats import multivariate_normal, skewnorm, norm -from scipy.ndimage import convolve1d -from abc import ABCMeta, abstractmethod from numpy.random import default_rng +from scipy.ndimage import convolve1d +from scipy.stats import multivariate_normal, norm, skewnorm + +from ctapipe.image.hillas import camera_to_shower_coordinates +from ctapipe.utils import linalg __all__ = [ "WaveformModel", @@ -186,8 +188,7 @@ class ImageModel(metaclass=ABCMeta): @u.quantity_input(x=u.m, y=u.m) @abstractmethod def pdf(self, x, y): - """Probability density function. - """ + """Probability density function.""" def generate_image(self, camera, intensity=50, nsb_level_pe=20, rng=None): """Generate a randomized DL1 shower image. @@ -286,8 +287,7 @@ def pdf(self, x, y): class SkewedGaussian(ImageModel): - """A shower image that has a skewness along the major axis. - """ + """A shower image that has a skewness along the major axis.""" @u.quantity_input(x=u.m, y=u.m, length=u.m, width=u.m) def __init__(self, x, y, length, width, psi, skewness): @@ -327,8 +327,8 @@ def _moments_to_parameters(self): delta = np.sign(self.skewness) * np.sqrt( (np.pi / 2 * skew23) / (skew23 + (0.5 * (4 - np.pi)) ** (2 / 3)) ) - a = delta / np.sqrt(1 - delta ** 2) - scale = self.length.to_value(u.m) / np.sqrt(1 - 2 * delta ** 2 / np.pi) + a = delta / np.sqrt(1 - delta**2) + scale = self.length.to_value(u.m) / np.sqrt(1 - 2 * delta**2 / np.pi) loc = -scale * delta * np.sqrt(2 / np.pi) return a, loc, scale diff --git a/ctapipe/visualization/bokeh.py b/ctapipe/visualization/bokeh.py index c99b83aba39..117b7f66aab 100644 --- a/ctapipe/visualization/bokeh.py +++ b/ctapipe/visualization/bokeh.py @@ -1,32 +1,30 @@ import sys import tempfile from abc import ABCMeta -import matplotlib.pyplot as plt -from matplotlib.colors import to_hex +import astropy.units as u +import matplotlib.pyplot as plt import numpy as np - -from bokeh.io import output_notebook, push_notebook, show, output_file -from bokeh.plotting import figure +from bokeh.io import output_file, output_notebook, push_notebook, show from bokeh.models import ( - ColumnDataSource, - TapTool, + BoxZoomTool, + CategoricalColorMapper, ColorBar, - LinearColorMapper, - LogColorMapper, + ColumnDataSource, ContinuousColorMapper, - CategoricalColorMapper, - HoverTool, - BoxZoomTool, Ellipse, + HoverTool, Label, + LinearColorMapper, + LogColorMapper, + TapTool, ) -from bokeh.palettes import Viridis256, Magma256, Inferno256, Greys256, d3 -import astropy.units as u +from bokeh.palettes import Greys256, Inferno256, Magma256, Viridis256, d3 +from bokeh.plotting import figure +from matplotlib.colors import to_hex from ..instrument import CameraGeometry, PixelShape - PLOTARGS = dict(tools="", toolbar_location=None, outline_line_color="#595959") @@ -644,7 +642,7 @@ def _init_datasource(self, subarray, values, *, radius, frame, scale, alpha): for i, telescope_id in enumerate(telescope_ids): telescope = subarray.tel[telescope_id] tel_types.append(str(telescope)) - mirror_area = telescope.optics.mirror_area.to_value(u.m ** 2) + mirror_area = telescope.optics.mirror_area.to_value(u.m**2) mirror_radii[i] = np.sqrt(mirror_area) / np.pi if values is None: