Merge branch 'master' into dl3_energy_dependent_alpha_cut

docs/conf.py

@@ -47,6 +47,7 @@
     "nbsphinx",
     "sphinx_automodapi.automodapi",
     "sphinxarg.ext",
+    "sphinxcontrib.mermaid",
 ]
 
 numpydoc_show_class_members = False

docs/examples/dl3_tool_config.json

@@ -1,22 +1,28 @@
 {
   "EventSelector": {
     "filters": {
-      "intensity": [100, Infinity],
+      "intensity": [0, Infinity],
       "width": [0, Infinity],
       "length": [0, Infinity],
       "r": [0, 1],
-      "wl": [0.1, 1],
-      "leakage_intensity_width_2": [0, 0.2],
+      "wl": [0.01, 1],
+      "leakage_intensity_width_2": [0, 1],
       "event_type": [32, 32]
     }
   },
   "DL3Cuts": {
+    "min_event_p_en_bin": 100,
     "global_gh_cut": 0.7,
     "gh_efficiency": 0.9,
+    "min_gh_cut": 0.1,
+    "max_gh_cut": 0.95,
     "global_alpha_cut": 10,
     "global_theta_cut": 0.2,
     "theta_containment": 0.68,
     "alpha_containment": 0.68,
+    "min_theta_cut": 0.05,
+    "max_theta_cut": 0.32,
+    "fill_theta_cut": 0.32,
     "allowed_tels": [1]
   }
 }

docs/examples/irf_tool_config.json

@@ -1,22 +1,28 @@
 {
   "EventSelector": {
     "filters": {
-      "intensity": [100, Infinity],
+      "intensity": [0, Infinity],
       "width": [0, Infinity],
       "length": [0, Infinity],
       "r": [0, 1],
-      "wl": [0.1, 1],
-      "leakage_intensity_width_2": [0, 0.2],
+      "wl": [0.01, 1],
+      "leakage_intensity_width_2": [0, 1],
       "event_type": [32, 32]
     }
   },
   "DL3Cuts": {
+    "min_event_p_en_bin": 100,
     "global_gh_cut": 0.7,
     "gh_efficiency": 0.9,
+    "min_gh_cut": 0.1,
+    "max_gh_cut": 0.95,
     "global_alpha_cut": 10,
     "global_theta_cut": 0.2,
     "theta_containment": 0.68,
     "alpha_containment": 0.68,
+    "min_theta_cut": 0.05,
+    "max_theta_cut": 0.32,
+    "fill_theta_cut": 0.32,
     "allowed_tels": [1]
   },
   "DataBinning": {

environment.yml

@@ -28,7 +28,8 @@ dependencies:
   - pandas
   - pymongo
   - seaborn
+  - jinja2=3.0 # pinned for bokeh 1.0 compatibility
   - pip:
-      - ctapipe_io_lst~=0.16.1
+      - ctapipe_io_lst~=0.18.0
       - ctaplot~=0.5.5
       - pyirf~=0.6.0

lstchain/data/lstchain_standard_config.json

@@ -35,6 +35,12 @@
     "r": [0, Infinity],
     "leakage_intensity_width_2": [0, Infinity]
   },
+  "n_training_events": {
+    "gamma_regressors": 1.0,
+    "gamma_tmp_regressors": 0.8,
+    "gamma_classifier": 0.2,
+    "proton_classifier": 1.0
+  },
 
   "tailcut": {
     "picture_thresh":8,

lstchain/high_level/significance_calculation.py

@@ -83,7 +83,7 @@ def analyze_wobble(config):
     for key, value in config['preselection'].items():
         LOGGER.debug('\nParameter: %s, range: %s, value type: %s', key, value, type(value))
 
-    selected_data = filter_events(data, config['preselection'])
+    selected_data = filter_events(data, filters=config['preselection'])
     # Add theta2 to selected data
     true_source_position = extract_source_position(selected_data, config['input']['observed_source'])
     plotting.plot_wobble(true_source_position, n_points, ax1)

lstchain/io/config.py

@@ -42,6 +42,7 @@ def get_standard_config():
     standard_config_file = os.path.join(os.path.dirname(__file__), "../data/lstchain_standard_config.json")
     return read_configuration_file(standard_config_file)
 
+
 def get_srcdep_config():
     """
     Load the config for source-dependent analysis from the file 'data/lstchain_src_dep_config.json'
@@ -53,6 +54,7 @@ def get_srcdep_config():
     srcdep_config_file = os.path.join(os.path.dirname(__file__), "../data/lstchain_src_dep_config.json")
     return read_configuration_file(srcdep_config_file)
 
+
 def replace_config(base_config, new_config):
     """
     Return a copy of the base_config with new configuration from new_config
@@ -92,3 +94,4 @@ def get_cleaning_parameters(config, clean_method_name):
     isolated_pixels = config[clean_method_name]['keep_isolated_pixels']
     min_n_picture_neighbors = config[clean_method_name]['min_number_picture_neighbors']
     return picture_th, boundary_th, isolated_pixels, min_n_picture_neighbors
+

lstchain/io/event_selection.py

@@ -52,6 +52,11 @@ class DL3Cuts(Component):
     Selection cuts for DL2 to DL3 conversion
     """
 
+    min_event_p_en_bin = Float(
+        help="Minimum events per energy bin, to evaluate percentile cuts",
+        default_value=100,
+    ).tag(config=True)
+
     global_gh_cut = Float(
         help="Global selection cut for gh_score (gammaness)",
         default_value=0.6,
@@ -62,13 +67,39 @@ class DL3Cuts(Component):
         default_value=0.95,
     ).tag(config=True)
 
+    min_gh_cut = Float(
+        help="Minimum gh_score (gammaness) cut in an energy bin",
+        default_value=0.1,
+    ).tag(config=True)
+
+    max_gh_cut = Float(
+        help="Maximum gh_score (gammaness) cut in an energy bin",
+        default_value=0.95,
+    ).tag(config=True)
+
+    min_theta_cut = Float(
+        help="Minimum theta cut (deg) in an energy bin",
+        default_value=0.05,
+    ).tag(config=True)
+
+    max_theta_cut = Float(
+        help="Maximum theta cut (deg) in an energy bin",
+        default_value=0.32,
+    ).tag(config=True)
+
+    fill_theta_cut = Float(
+        help="Fill value of theta cut (deg) in an energy bin with fewer " +
+            "than minimum number of events present",
+        default_value=0.32,
+    ).tag(config=True)
+
     theta_containment = Float(
         help="Percentage containment region for theta cuts",
         default_value=0.68,
     ).tag(config=True)
 
     global_theta_cut = Float(
-        help="Global selection cut for theta",
+        help="Global selection cut (deg) for theta",
         default_value=0.2,
     ).tag(config=True)
 
@@ -78,7 +109,7 @@ class DL3Cuts(Component):
     ).tag(config=True)
 
     global_alpha_cut = Float(
-        help="Global selection cut for alpha",
+        help="Global selection cut (deg) for alpha",
         default_value=20,
     ).tag(config=True)
 
@@ -95,8 +126,7 @@ def apply_global_gh_cut(self, data):
         return data[data["gh_score"] > self.global_gh_cut]
 
     def energy_dependent_gh_cuts(
-        self, data, energy_bins, min_value=0.1,
-        max_value=0.99, smoothing=None, min_events=10
+        self, data, energy_bins, smoothing=None
     ):
         """
         Evaluating energy-dependent gammaness cuts, in a given
@@ -108,12 +138,12 @@ def energy_dependent_gh_cuts(
             data["gh_score"],
             data["reco_energy"],
             bins=energy_bins,
-            min_value=min_value,
-            max_value=max_value,
+            min_value=self.min_gh_cut,
+            max_value=self.max_gh_cut,
             fill_value=data["gh_score"].max(),
             percentile=100 * (1 - self.gh_efficiency),
             smoothing=smoothing,
-            min_events=min_events,
+            min_events=self.min_event_p_en_bin,
         )
         return gh_cuts
 
@@ -128,7 +158,7 @@ def apply_energy_dependent_gh_cuts(self, data, gh_cuts):
             data["gh_score"],
             data["reco_energy"],
             gh_cuts,
-            operator.ge,
+            operator.gt,
         )
         return data[data["selected_gh"]]
 
@@ -139,9 +169,7 @@ def apply_global_theta_cut(self, data):
         return data[data["theta"].to_value(u.deg) < self.global_theta_cut]
 
     def energy_dependent_theta_cuts(
-        self, data, energy_bins, min_value=0.05 * u.deg,
-        fill_value=0.32 * u.deg, max_value=0.32 * u.deg,
-        smoothing=None, min_events=10
+        self, data, energy_bins, smoothing=None,
     ):
         """
         Evaluating an optimized energy-dependent theta cuts, in a given
@@ -155,12 +183,12 @@ def energy_dependent_theta_cuts(
             data["theta"],
             data["reco_energy"],
             bins=energy_bins,
-            min_value=min_value,
-            max_value=max_value,
-            fill_value=fill_value,
+            min_value=self.min_theta_cut * u.deg,
+            max_value=self.max_theta_cut * u.deg,
+            fill_value=self.fill_theta_cut * u.deg,
             percentile=100 * self.theta_containment,
             smoothing=smoothing,
-            min_events=min_events,
+            min_events=self.min_event_p_en_bin,
         )
         return theta_cuts
 
@@ -174,7 +202,7 @@ def apply_energy_dependent_theta_cuts(self, data, theta_cuts):
             data["theta"],
             data["reco_energy"],
             theta_cuts,
-            operator.le,
+            operator.lt,
         )
         return data[data["selected_theta"]]
 

lstchain/io/tests/test_eventselection.py

@@ -63,6 +63,7 @@ def test_dl3_energy_dependent_cuts():
 
     temp_cuts.gh_max_efficiency = 0.8
     temp_cuts.theta_containment = 0.68
+    temp_cuts.min_event_p_en_bin = 2
 
     temp_data = QTable({
         "gh_score": u.Quantity(np.tile(np.arange(0.35, 0.85, 0.05), 3)),
@@ -74,11 +75,11 @@ def test_dl3_energy_dependent_cuts():
     en_range = u.Quantity([0.01, 0.1, 1, 10, 100, np.inf], unit=u.TeV)
 
     theta_cut = temp_cuts.energy_dependent_theta_cuts(
-        temp_data, en_range, min_events=2
+        temp_data, en_range,
     )
 
     gh_cut = temp_cuts.energy_dependent_gh_cuts(
-        temp_data, en_range, min_events=2
+        temp_data, en_range,
     )
 
     data_th = temp_cuts.apply_energy_dependent_theta_cuts(temp_data, theta_cut)