1gev min

src/HHbbVV/postprocessing/corrections.py

@@ -280,7 +280,8 @@ def postprocess_lpsfs(
                 td[key][jet_match] = events[key][jet_match].squeeze()
 
             # distortion uncertainty
-            dist_sf = np.clip(np.nan_to_num(events["lp_sf_dist"][jet_match][0], nan=1), 0.8, 1.2)
+            # dist_sf = np.clip(np.nan_to_num(events["lp_sf_dist"][jet_match][0], nan=1), 0.5, 2.0)
+            dist_sf = np.nan_to_num(events["lp_sf_dist"][jet_match][0], nan=1)
             td["lp_sf_dist_up"][jet_match] = td["lp_sf_nom"][jet_match] * dist_sf
             td["lp_sf_dist_down"][jet_match] = td["lp_sf_nom"][jet_match] / dist_sf
             # breakpoint()

src/HHbbVV/postprocessing/plotting.py

@@ -1272,3 +1272,21 @@ def cutsLinePlot(
         plt.show()
     else:
         plt.close()
+
+
+def plot_lund_plane(h: np.ndarray, ax=None, name: str = None, show: bool = False):
+    if ax is None:
+        fig, ax = plt.subplots(figsize=(12, 14))
+
+    extent = [-1, 8, -5, 7]
+    ax.imshow(h.T[::-1], extent=extent)
+    ax.colorbar()
+
+    if ax is None:
+        if len(name):
+            plt.savefig(name, bbox_inches="tight")
+
+        if show:
+            plt.show()
+        else:
+            plt.close()

src/HHbbVV/postprocessing/postprocessing.py

@@ -1345,7 +1345,7 @@ def lpsfs(
             for sig_key in sig_keys:
                 systs = rsysts[sig_key]
                 sf_table[sig_key] = {
-                    "SF": f"{systs['lp_sf']:.2f} +{systs['lp_sf'] * systs['lp_sf_unc_up']:.2f}-{systs['lp_sf'] * systs['lp_sf_unc_down']:.2f}-",
+                    "SF": f"{systs['lp_sf']:.2f} +{systs['lp_sf'] * systs['lp_sf_unc_up']:.2f}-{systs['lp_sf'] * systs['lp_sf_unc_down']:.2f}",
                     **systs["lp_sf_uncs_sym"],
                 }
                 for key in systs["lp_sf_uncs_asym"]["up"]:

src/HHbbVV/processors/corrections.py

@@ -754,7 +754,7 @@ def _get_lund_lookups(
         # 0s -> 1 in the ratio
         mc_sig_ratio = np.nan_to_num((mc_nom / mc_tot) / (sig_lp_hist.values() / sig_tot), nan=1.0)
         mc_sig_ratio[mc_sig_ratio == 0] = 1.0
-        mc_sig_ratio = np.clip(mc_sig_ratio, 0.5, 2.0)
+        # mc_sig_ratio = np.clip(mc_sig_ratio, 0.5, 2.0)
 
         ratio_dist = dense_lookup(mc_sig_ratio, ratio_edges)
 
@@ -849,6 +849,7 @@ def _get_lund_arrays(
     jec_fatjets: FatJetArray,
     fatjet_idx: int | ak.Array,
     num_prongs: int,
+    min_pt: float = 1.0,
 ):
     """
     Gets the ``num_prongs`` subjet pTs and Delta and kT per primary LP splitting of fatjets at
@@ -902,6 +903,8 @@ def _get_lund_arrays(
     kt_subjets_pt = kt_subjets_vec.pt * jec_correction
     # get constituents
     kt_subjet_consts = kt_clustering.exclusive_jets_constituents(num_prongs)
+    # breakpoint()
+    kt_subjet_consts = kt_subjet_consts[kt_subjet_consts.pt > min_pt]
 
     # then re-cluster with CA
     # won't need to flatten once https://github.com/scikit-hep/fastjet/pull/145 is released
@@ -942,8 +945,7 @@ def _calc_lund_SFs(
     pt_extrap_lookups: list[dense_lookup],
     max_pt_bin: int = MAX_PT_BIN,
     max_fparams: int = MAX_PT_FPARAMS,
-    clip_max: float = 10,
-    clip_min: float = 0.1,
+    CLIP: float = 5.0,
 ) -> np.ndarray:
     """
     Calculates scale factors for jets based on splittings in the primary Lund Plane.
@@ -985,12 +987,12 @@ def _calc_lund_SFs(
             # only recalculate if there are multiple pt param lookup tables
             if j == 0 or len(pt_extrap_lookups) > 1:
                 params = pt_extrap_lookup(hpt_logD, hpt_logkt)
-                pt_extrap_vals = np.maximum(
-                    np.minimum(np.sum(params * sj_pt_orders, axis=1), clip_max), clip_min
-                )
+                pt_extrap_vals = np.sum(params * sj_pt_orders, axis=1)
 
             ratio_vals[high_pt_sel] = pt_extrap_vals
 
+            ratio_vals = np.clip(ratio_vals, 1.0 / CLIP, CLIP)
+
             if len(ld_offsets) != 1:
                 # recover jagged event structure
                 reshaped_ratio_vals = ak.Array(
@@ -1131,6 +1133,8 @@ def get_lund_SFs(
             [pt_extrap_lookups_dict["params"]],
         )
 
+    print("lp sf sys")
+
     sfs["lp_sf_sys_down"] = _calc_lund_SFs(
         flat_logD,
         flat_logkt,
@@ -1153,6 +1157,7 @@ def get_lund_SFs(
 
     if ratio_dist is not None:
         # breakpoint()
+        print("lp sf dist")
         sfs["lp_sf_dist"] = _calc_lund_SFs(
             flat_logD,
             flat_logkt,

src/HHbbVV/processors/utils.py

@@ -37,6 +37,25 @@ def pad_val(
     return ret.to_numpy() if to_numpy else ret
 
 
+def ak_clip(array: ak.Array, min_val: float, max_val: float):
+    """
+    Clip the values in an Awkward Array to lie within the interval [min_val, max_val].
+
+    Parameters:
+    - array: ak.Array
+        The input Awkward Array.
+    - min_val: scalar
+        The minimum value to clip to.
+    - max_val: scalar
+        The maximum value to clip to.
+
+    Returns:
+    - ak.Array
+        A new Awkward Array with values clipped within the specified range.
+    """
+    return ak.where(array < min_val, min_val, ak.where(array > max_val, max_val, array))
+
+
 def add_selection(
     name: str,
     sel: np.ndarray,

src/condor/submit_configs/skimmer_24_07_24_signal_lp.yaml

@@ -16,8 +16,8 @@
               "VBF_HHTobbVV_CV_1_C2V_1_C3_0",
               "VBF_HHTobbVV_CV_1_C2V_1_C3_2",
             ],
-          "files_per_job": 10,
-          "max_files": 60,
+          "files_per_job": 8,
+          "max_files": 80,
           "chunksize": 80000,
           # "maxchunks": 20
         },