Fixing bug in calculation of sys shape effect for evaluation plots.

rdstar/offline_analysis/fitting/dedicated_fit_approach/plotting_tools.py

@@ -210,18 +210,20 @@ class SpecificShapePlotInfoContainer:
     latex_str: str
     subset_index: int
     normed_base_shape: np.ndarray
-    relative_shape_error: np.ndarray
-    relative_stat_error: np.ndarray
+    shape_error: np.ndarray
+    stat_error: np.ndarray
     observable: FitObservableInfo
     reco_ch_info: RecoChannelInfo
     component_info: ComponentInfo
     scale_factor: Optional[float] = None
+    add_statistical_uncertainty: bool = True
 
     def __post_init__(self) -> None:
         assert len(self.normed_base_shape.shape) == 1, (len(self.normed_base_shape.shape), self.normed_base_shape.shape)
-        assert self.normed_base_shape.shape == self.relative_shape_error.shape, (
+
+        assert self.normed_base_shape.shape == self.shape_error.shape, (
             self.normed_base_shape.shape,
-            self.relative_shape_error.shape,
+            self.shape_error.shape,
         )
 
 
@@ -273,6 +275,7 @@ class SystematicsShapePlotter:
         selected_relative_stat_error: np.ndarray,
         reco_ch_info: RecoChannelInfo,
         component_info: ComponentInfo,
+        add_statistical_uncertainty: bool,
     ) -> Generator[SpecificShapePlotInfoContainer, None, None]:
 
         assert len(selected_normed_base_shape.shape) == 1, selected_normed_base_shape.shape
@@ -315,7 +318,7 @@ class SystematicsShapePlotter:
         )
 
         # Calculate projections for uncertainties
-        projected_relative_shape_error_splits: Tuple[Tuple[np.ndarray, np.ndarray], ...] = tuple(
+        projected_shape_error_splits: Tuple[Tuple[np.ndarray, np.ndarray], ...] = tuple(
             (
                 (rel_err * _2d_normed_base_shape).sum(axis=1),
                 (rel_err * _2d_normed_base_shape).sum(axis=0),
@@ -324,14 +327,14 @@ class SystematicsShapePlotter:
         )
 
         sum_of_2d_relative_shape_error_splits: np.ndarray = np.sum(_2d_relative_shape_error_splits, axis=0)
-        projected_relative_shape_error_combination: Tuple[np.ndarray, np.ndarray] = (
+        projected_shape_error_combination: Tuple[np.ndarray, np.ndarray] = (
             (sum_of_2d_relative_shape_error_splits * _2d_normed_base_shape).sum(axis=1),
             (sum_of_2d_relative_shape_error_splits * _2d_normed_base_shape).sum(axis=0),
         )
 
         # Stat error handling
         _2d_relative_stat_error: np.ndarray = selected_relative_stat_error.reshape(binning_info)
-        projected_relative_stat_error: Tuple[np.ndarray, np.ndarray] = (
+        projected_stat_error: Tuple[np.ndarray, np.ndarray] = (
             np.sqrt((np.square(_2d_relative_stat_error * _2d_normed_base_shape)).sum(axis=1)),
             np.sqrt((np.square(_2d_relative_stat_error * _2d_normed_base_shape)).sum(axis=0)),
         )
@@ -352,12 +355,13 @@ class SystematicsShapePlotter:
                         latex_str=sys_shape_info.latex_str,
                         subset_index=subset_id,
                         normed_base_shape=projected_normed_base_shape,
-                        relative_shape_error=projected_relative_shape_error_splits[subset_id][obs_index],
-                        relative_stat_error=projected_relative_stat_error[obs_index],
+                        shape_error=projected_shape_error_splits[subset_id][obs_index],
+                        stat_error=projected_stat_error[obs_index],
                         observable=obs_info,
                         reco_ch_info=reco_ch_info,
                         component_info=component_info,
                         scale_factor=sys_shape_info.scale_factor,
+                        add_statistical_uncertainty=add_statistical_uncertainty,
                     )
 
             yield SpecificShapePlotInfoContainer(
@@ -365,17 +369,19 @@ class SystematicsShapePlotter:
                 latex_str=sys_shape_info.latex_str,
                 subset_index=-1,
                 normed_base_shape=projected_normed_base_shape,
-                relative_shape_error=projected_relative_shape_error_combination[obs_index],
-                relative_stat_error=projected_relative_stat_error[obs_index],
+                shape_error=projected_shape_error_combination[obs_index],
+                stat_error=projected_stat_error[obs_index],
                 observable=obs_info,
                 reco_ch_info=reco_ch_info,
                 component_info=component_info,
                 scale_factor=sys_shape_info.scale_factor,
+                add_statistical_uncertainty=add_statistical_uncertainty,
             )
 
     def create_systematics_shape_plots(
         self,
         sys_shape_info: ShapePlotInfoContainer,
+        add_statistical_uncertainties: bool = True,
     ) -> Generator[SpecificShapePlotInfoContainer, None, None]:
         target_dir_path: PathType = os.path.join(self.output_dir_path, sys_shape_info.name)
         os.makedirs(target_dir_path, exist_ok=True)
@@ -432,6 +438,7 @@ class SystematicsShapePlotter:
                     selected_relative_stat_error=_rel_shape_error_stat_reco_ch_splits[reco_ch_index],
                     reco_ch_info=reco_ch_info,
                     component_info=comp_info,
+                    add_statistical_uncertainty=add_statistical_uncertainties,
                 ):
                     self.plot_systematics_shape_effect_for(
                         spec_sys_shape_info=_this_shape_info,
@@ -485,12 +492,10 @@ class SystematicsShapePlotter:
         bin_widths: np.ndarray = bin_edges[1:] - bin_edges[:-1]
 
         norm_bin_count: np.ndarray = spec_sys_shape_info.normed_base_shape
-        u_bin_shape_uncert: np.ndarray = spec_sys_shape_info.normed_base_shape * spec_sys_shape_info.relative_shape_error
-        d_bin_shape_uncert: np.ndarray = spec_sys_shape_info.normed_base_shape * (
-            -1.0 * spec_sys_shape_info.relative_shape_error
-        )
+        u_bin_shape_uncert: np.ndarray = spec_sys_shape_info.shape_error
+        d_bin_shape_uncert: np.ndarray = -1.0 * spec_sys_shape_info.shape_error
 
-        bin_stat_uncert: np.ndarray = spec_sys_shape_info.relative_stat_error
+        bin_stat_uncert: np.ndarray = spec_sys_shape_info.stat_error
 
         assert np.all(u_bin_shape_uncert * d_bin_shape_uncert <= 0.0), (
             np.sum(u_bin_shape_uncert * d_bin_shape_uncert > 0.0),
@@ -554,17 +559,18 @@ class SystematicsShapePlotter:
         )
 
         # Pull Plot
-        ax2.bar(
-            x=bin_mids,
-            height=2 * bin_stat_uncert,
-            width=bin_widths,
-            bottom=-1.0 * bin_stat_uncert,
-            color=KITColors.grey,
-            alpha=0.4,
-            fill=True,
-            lw=0,
-            label="stat.",
-        )
+        if spec_sys_shape_info.add_statistical_uncertainty:
+            ax2.bar(
+                x=bin_mids,
+                height=2 * bin_stat_uncert,
+                width=bin_widths,
+                bottom=-1.0 * bin_stat_uncert,
+                color=KITColors.grey,
+                alpha=0.4,
+                fill=True,
+                lw=0,
+                label="stat.",
+            )
         ax2.hist(
             x=[bin_mids, bin_mids],
             bins=bin_edges,