From eff4479059eda079b7f4d18110696a0b18efd227 Mon Sep 17 00:00:00 2001
From: Ian0sborne <ian.osborne@student.manchester.ac.uk>
Date: Wed, 25 Feb 2026 16:27:33 +0000
Subject: [PATCH 01/69] add sipm selection function using median and std

---
 invisible_cities/reco/wfm_functions.py | 22 ++++++++++++++++++++++
 1 file changed, 22 insertions(+)

diff --git a/invisible_cities/reco/wfm_functions.py b/invisible_cities/reco/wfm_functions.py
index 4ecced9e21..ff463d70f6 100644
--- a/invisible_cities/reco/wfm_functions.py
+++ b/invisible_cities/reco/wfm_functions.py
@@ -3,6 +3,8 @@
 authors: J.J. Gomez-Cadenas, G. Martinez
 """
 import numpy as np
+from typing import Optional
+from typing  import Callable
 
 from .. core.core_functions import define_window
 from .. calib               import calib_sensors_functions as csf
@@ -129,3 +131,23 @@ def compare_cwf_blr(cwf, pmtblr, event_list, window_size=500):
             DIFF.append(diff)
 
     return np.array(DIFF)
+
+
+def median_std_method(wfs    : np.ndarray,  
+                      nsigma : Optional[float] = 3.) -> np.ndarray:
+    """
+    Computes the median and standard deviation of the time summed SiPM 
+    waveforms and selects the SiPMs that are nsigma over the median.
+
+    Parameters
+    ----------
+    wfs    : 2D array of shape (n_sipms, n_time_bins) containing the waveforms of each SiPM.
+    nsigma : Number of standard deviations above the median, default 3.
+    
+    Returns
+    -------
+    Boolean numpy array of shape (n_sipms,) where True indicates that the SiPM is selected.
+    """
+    charges = np.sum(wfs, axis=1)
+    threshold = np.median(charges) + nsigma * np.std(charges)
+    return charges >= threshold
\ No newline at end of file

From 77c47ad2776f5a7113012f8e7833f6b829aad0e0 Mon Sep 17 00:00:00 2001
From: Ian0sborne <ian.osborne@student.manchester.ac.uk>
Date: Wed, 25 Feb 2026 16:31:38 +0000
Subject: [PATCH 02/69] add sipm selection method using a charge threshold

---
 invisible_cities/reco/wfm_functions.py | 20 +++++++++++++++++++-
 1 file changed, 19 insertions(+), 1 deletion(-)

diff --git a/invisible_cities/reco/wfm_functions.py b/invisible_cities/reco/wfm_functions.py
index ff463d70f6..9e0456bc54 100644
--- a/invisible_cities/reco/wfm_functions.py
+++ b/invisible_cities/reco/wfm_functions.py
@@ -150,4 +150,22 @@ def median_std_method(wfs    : np.ndarray,
     """
     charges = np.sum(wfs, axis=1)
     threshold = np.median(charges) + nsigma * np.std(charges)
-    return charges >= threshold
\ No newline at end of file
+    return charges >= threshold
+
+
+def charge_threshold_method(wfs       : np.ndarray,
+                            threshold : Optional[float] = 5.) -> np.ndarray:
+    """
+    Selects the SiPMs whose time summed waveforms are above a threshold.
+    
+    Parameters
+    ----------
+    wfs       : 2D array of shape (n_sipms, n_time_bins) containing the waveforms of each SiPM.
+    threshold : Charge threshold in PE, default 5.
+    
+    Returns
+    -------
+    Boolean numpy array of shape (n_sipms,) where True indicates that the SiPM is selected.
+    """
+    charges = np.sum(wfs, axis=1)
+    return charges >= threshold

From aa944154bc7cb3467575655d136c938356c27015 Mon Sep 17 00:00:00 2001
From: Ian0sborne <ian.osborne@student.manchester.ac.uk>
Date: Wed, 25 Feb 2026 16:32:32 +0000
Subject: [PATCH 03/69] add sipm selection function using top n most energetic
 sipms

---
 invisible_cities/reco/wfm_functions.py | 22 ++++++++++++++++++++++
 1 file changed, 22 insertions(+)

diff --git a/invisible_cities/reco/wfm_functions.py b/invisible_cities/reco/wfm_functions.py
index 9e0456bc54..62468b6ae8 100644
--- a/invisible_cities/reco/wfm_functions.py
+++ b/invisible_cities/reco/wfm_functions.py
@@ -169,3 +169,25 @@ def charge_threshold_method(wfs       : np.ndarray,
     """
     charges = np.sum(wfs, axis=1)
     return charges >= threshold
+
+
+def top_n_method(wfs : np.ndarray, 
+                 n   : Optional[int] = 10) -> np.ndarray:
+    """
+    Selects the SiPMs with the top n highest time summed waveforms.
+        
+    Parameters
+    ----------
+    wfs       : 2D array of shape (n_sipms, n_time_bins) containing the waveforms of each SiPM.
+    n         : Number of most energeticSiPMs to select, default 10.
+    
+    Returns
+    -------
+    Boolean numpy array of shape (n_sipms,) where True indicates that the SiPM is selected.
+    """
+    charges = np.sum(wfs, axis=1)
+    idx = np.argsort(charges)[-n:]
+
+    selected_ids = np.zeros_like(charges, dtype=bool)
+    selected_ids[idx] = True
+    return selected_ids

From b41d2cacd6b2be378a0f09d2f30a2801ce2f569e Mon Sep 17 00:00:00 2001
From: Ian0sborne <ian.osborne@student.manchester.ac.uk>
Date: Wed, 25 Feb 2026 16:46:58 +0000
Subject: [PATCH 04/69] add function that kills isolated sipms

---
 invisible_cities/reco/wfm_functions.py | 35 ++++++++++++++++++++++++++
 1 file changed, 35 insertions(+)

diff --git a/invisible_cities/reco/wfm_functions.py b/invisible_cities/reco/wfm_functions.py
index 62468b6ae8..625a337b3e 100644
--- a/invisible_cities/reco/wfm_functions.py
+++ b/invisible_cities/reco/wfm_functions.py
@@ -191,3 +191,38 @@ def top_n_method(wfs : np.ndarray,
     selected_ids = np.zeros_like(charges, dtype=bool)
     selected_ids[idx] = True
     return selected_ids
+
+
+def kill_isolated_sipms(selected_ids        : np.ndarray,
+                        sipm_x              : np.ndarray, 
+                        sipm_y              : np.ndarray, 
+                        proximity_threshold : float) -> np.ndarray:
+    """
+    For the SiPMs that have passed the previous selection, scans through the SiPMs to check if they 
+    have neighbouring SiPMs - i.e., within the proximity_threshold - that have also passed the selection. 
+    If no neighbours are found, the SiPMs are classed as isolated, and are removed.
+
+    Parameters
+    ----------
+    selected_ids        : Boolean array of shape (n_sipms,) indicating which SiPMs passed the previous selection.
+    sipm_x              : 1D array of shape (n_sipms,) containing the x positions of the SiPMs.
+    sipm_y              : 1D array of shape (n_sipms,) containing the y positions of the SiPMs.
+    proximity_threshold : Distance threshold in mm used to identify isolated SiPMs.
+
+    Returns
+    -------
+    selected_ids_no_isolated : Boolean array of shape (n_sipms,) where True indicates that the SiPM is selected.
+    """
+    selected_ids_no_isolated = selected_ids.copy()
+
+    for i in np.where(selected_ids)[0]:
+        x, y = sipm_x[i], sipm_y[i]
+
+        distances = np.sqrt((sipm_x - x)**2 + (sipm_y - y)**2)
+
+        n_neighbors = np.sum((distances < proximity_threshold) & selected_ids)
+
+        if n_neighbors <= 1:
+            selected_ids_no_isolated[i] = False
+            
+    return selected_ids_no_isolated
\ No newline at end of file

From d636192911410cfe76b98970cb072d56e018d814 Mon Sep 17 00:00:00 2001
From: Ian0sborne <ian.osborne@student.manchester.ac.uk>
Date: Wed, 25 Feb 2026 16:50:48 +0000
Subject: [PATCH 05/69] add function that creates circular padding around
 selected sipms

---
 invisible_cities/reco/wfm_functions.py | 31 +++++++++++++++++++++++++-
 1 file changed, 30 insertions(+), 1 deletion(-)

diff --git a/invisible_cities/reco/wfm_functions.py b/invisible_cities/reco/wfm_functions.py
index 625a337b3e..399bf1b0fc 100644
--- a/invisible_cities/reco/wfm_functions.py
+++ b/invisible_cities/reco/wfm_functions.py
@@ -225,4 +225,33 @@ def kill_isolated_sipms(selected_ids        : np.ndarray,
         if n_neighbors <= 1:
             selected_ids_no_isolated[i] = False
             
-    return selected_ids_no_isolated
\ No newline at end of file
+    return selected_ids_no_isolated
+
+
+def apply_circular_padding(selected_ids_no_isolated : np.ndarray,
+                           sipm_x                   : np.ndarray, 
+                           sipm_y                   : np.ndarray, 
+                           padding_radius           : float) -> np.ndarray:
+    """
+    For the SiPMs that pass the previous selection, creates circular padding of radius padding_radius, 
+    selecting all SiPMs within that radius. Stores the union of all selected SiPMs.
+
+    Parameters
+    ----------
+    selected_ids_no_isolated : Boolean array of shape (n_sipms,) indicating which SiPMs passed the previous selection.
+    sipm_x                   : 1D array of shape (n_sipms,) containing the x positions of the SiPMs.
+    sipm_y                   : 1D array of shape (n_sipms,) containing the y positions of the SiPMs.
+    padding_radius           : Distance threshold in mm used to create circular padding around selected SiPMs.
+
+    Returns
+    -------
+    sipm_ids_with_signal : Boolean array of shape (n_sipms,) where True indicates that the SiPM is selected.
+    """
+    sipm_ids_with_signal = np.zeros_like(selected_ids_no_isolated, dtype=bool)
+
+    for i in np.where(selected_ids_no_isolated)[0]:
+        x, y = sipm_x[i], sipm_y[i]
+        distances = np.sqrt((sipm_x - x)**2 + (sipm_y - y)**2)
+        sipm_ids_with_signal |= distances < padding_radius
+
+    return sipm_ids_with_signal
\ No newline at end of file

From 2d5679c32943c1c1134c30f17c900227d2431ce3 Mon Sep 17 00:00:00 2001
From: Ian0sborne <ian.osborne@student.manchester.ac.uk>
Date: Wed, 25 Feb 2026 16:51:50 +0000
Subject: [PATCH 06/69] add function that creates the overall sipm masks

---
 invisible_cities/reco/wfm_functions.py | 46 +++++++++++++++++++++++++-
 1 file changed, 45 insertions(+), 1 deletion(-)

diff --git a/invisible_cities/reco/wfm_functions.py b/invisible_cities/reco/wfm_functions.py
index 399bf1b0fc..a39ed2520e 100644
--- a/invisible_cities/reco/wfm_functions.py
+++ b/invisible_cities/reco/wfm_functions.py
@@ -254,4 +254,48 @@ def apply_circular_padding(selected_ids_no_isolated : np.ndarray,
         distances = np.sqrt((sipm_x - x)**2 + (sipm_y - y)**2)
         sipm_ids_with_signal |= distances < padding_radius
 
-    return sipm_ids_with_signal
\ No newline at end of file
+    return sipm_ids_with_signal
+
+
+def make_sipm_selection(wfs                 : np.ndarray, 
+                        selection_func      : Callable, 
+                        selection_kwargs    : dict, 
+                        proximity_threshold : float, 
+                        padding_radius      : float) -> np.ndarray: 
+    """
+    SiPM selection pipeline, applies SiPM cuts based on user input.
+    A first selection of SiPMs is made, isolated SiPMs are removed
+    and padding is added around the SiPMs that are left.
+
+    Parameters
+    ----------
+    wfs                 : 2D array of shape (n_sipms, n_time_bins) containing the waveforms of each SiPM.
+    selection_func      : Function returning a boolean mask of energetic SiPMs.
+    selection_kwargs    : Dictionary of arguments passed to selection_func.
+    proximity_threshold : Threshold used to identify isolated SiPMs.
+    padding_radius      : Radial padding added to each SiPM that passes the selections.
+
+    Returns
+    -------
+    sipm_ids_with_signal : Array of shape (n_sipms,) with boolean values indicating which SiPMs are selected.
+    """
+
+    sipm_x = np.array(detector_info.X)
+    sipm_y = np.array(detector_info.Y)
+
+    selected_ids = selection_func(wfs, **selection_kwargs)
+
+    selected_ids_no_isolated = kill_isolated_sipms(
+        selected_ids,
+        sipm_x,
+        sipm_y,
+        proximity_threshold
+    )
+
+    sipm_ids_with_signal = apply_circular_padding(
+        selected_ids_no_isolated,
+        sipm_x,
+        sipm_y,
+        padding_radius
+    )
+    return sipm_ids_with_signal

From 96cc2bec04480fefc7099a6611bf58eeb305d18f Mon Sep 17 00:00:00 2001
From: Ian0sborne <ian.osborne@student.manchester.ac.uk>
Date: Wed, 25 Feb 2026 16:54:08 +0000
Subject: [PATCH 07/69] load detector geometry in

---
 invisible_cities/reco/wfm_functions.py | 5 +++--
 1 file changed, 3 insertions(+), 2 deletions(-)

diff --git a/invisible_cities/reco/wfm_functions.py b/invisible_cities/reco/wfm_functions.py
index a39ed2520e..fb43889662 100644
--- a/invisible_cities/reco/wfm_functions.py
+++ b/invisible_cities/reco/wfm_functions.py
@@ -3,12 +3,13 @@
 authors: J.J. Gomez-Cadenas, G. Martinez
 """
 import numpy as np
-from typing import Optional
+from typing  import Optional
 from typing  import Callable
 
 from .. core.core_functions import define_window
 from .. calib               import calib_sensors_functions as csf
 from .. sierpe              import blr
+from .. database            import load_db
 
 def to_adc(wfs, adc_to_pes):
     """
@@ -279,7 +280,7 @@ def make_sipm_selection(wfs                 : np.ndarray,
     -------
     sipm_ids_with_signal : Array of shape (n_sipms,) with boolean values indicating which SiPMs are selected.
     """
-
+    detector_info = load_db.DataSiPM('next100', 0)
     sipm_x = np.array(detector_info.X)
     sipm_y = np.array(detector_info.Y)
 

From 4af4b0203dcdfd59fd9440b2e017111201f374e6 Mon Sep 17 00:00:00 2001
From: casper <johnwaiton@hotmail.co.uk>
Date: Wed, 25 Feb 2026 18:33:44 +0000
Subject: [PATCH 08/69] remove thresholds from `calibrate_sipms()`

---
 invisible_cities/calib/calib_sensors_functions.py | 10 +++++-----
 invisible_cities/cities/components.py             |  3 +--
 invisible_cities/cities/irene.py                  |  2 +-
 3 files changed, 7 insertions(+), 8 deletions(-)

diff --git a/invisible_cities/calib/calib_sensors_functions.py b/invisible_cities/calib/calib_sensors_functions.py
index f8db736c3b..d5387bdb08 100644
--- a/invisible_cities/calib/calib_sensors_functions.py
+++ b/invisible_cities/calib/calib_sensors_functions.py
@@ -136,15 +136,15 @@ def pmt_subtract_maw(cwfs, n_maw=100):
     return cwfs - maw
 
 
-def calibrate_sipms(sipm_wfs, adc_to_pes, thr, *, bls_mode=BlsMode.mode):
+def calibrate_sipms(sipm_wfs, adc_to_pes, *, bls_mode=BlsMode.mode):
     """
-    Subtracts the baseline, calibrates waveforms to pes
-    and suppresses values below `thr` (in pes).
+    Subtract baseline, and calibrates waveforms to pes.
     """
-    thr  = to_col_vector(np.full(sipm_wfs.shape[0], thr))
+    #thr  = to_col_vector(np.full(sipm_wfs.shape[0], thr))
     bls  = subtract_baseline(sipm_wfs, bls_mode=bls_mode)
     cwfs = calibrate_wfs(bls, adc_to_pes)
-    return np.where(cwfs > thr, cwfs, 0)
+    return cwfs
+    #return np.where(cwfs > thr, cwfs, 0)
 
 
 def subtract_mean  (wfs): return subtract_baseline(wfs, bls_mode=BlsMode.mean  )
diff --git a/invisible_cities/cities/components.py b/invisible_cities/cities/components.py
index 7445c505ae..f32d2fdc7c 100644
--- a/invisible_cities/cities/components.py
+++ b/invisible_cities/cities/components.py
@@ -802,7 +802,7 @@ def calibrate_pmts(cwf):# -> CCwfs:
     return calibrate_pmts
 
 
-def calibrate_sipms(dbfile, run_number, thr_sipm):
+def calibrate_sipms(dbfile, run_number):
     DataSiPM   = load_db.DataSiPM(dbfile, run_number)
     DataSiPM   = DataSiPM.loc[lambda df: df.Active==1]
     adc_to_pes = np.abs(DataSiPM.adc_to_pes.values)
@@ -810,7 +810,6 @@ def calibrate_sipms(dbfile, run_number, thr_sipm):
     def calibrate_sipms(rwf):
         return csf.calibrate_sipms(rwf,
                                    adc_to_pes = adc_to_pes,
-                                   thr        = thr_sipm,
                                    bls_mode   = BlsMode.mode)
 
     return calibrate_sipms
diff --git a/invisible_cities/cities/irene.py b/invisible_cities/cities/irene.py
index 472c0af1e4..49251e4bca 100644
--- a/invisible_cities/cities/irene.py
+++ b/invisible_cities/cities/irene.py
@@ -100,7 +100,7 @@ def irene( files_in        : OneOrManyFiles
                               out  = ("s1_indices", "s2_indices", "s2_energies"))
 
     # Remove baseline and calibrate SiPMs
-    sipm_rwf_to_cal  = fl.map(calibrate_sipms(detector_db, run_number, sipm_thr),
+    sipm_rwf_to_cal  = fl.map(calibrate_sipms(detector_db, run_number),
                               item = "sipm")
 
     event_count_in  = fl.spy_count()

From 416b618633b36d7993baba309aed756a5416082a Mon Sep 17 00:00:00 2001
From: casper <johnwaiton@hotmail.co.uk>
Date: Wed, 25 Feb 2026 18:34:09 +0000
Subject: [PATCH 09/69] create `CutAlgo()` class

---
 invisible_cities/types/symbols.py | 2 ++
 1 file changed, 2 insertions(+)

diff --git a/invisible_cities/types/symbols.py b/invisible_cities/types/symbols.py
index 793a2e076d..9b4ab17994 100644
--- a/invisible_cities/types/symbols.py
+++ b/invisible_cities/types/symbols.py
@@ -155,6 +155,8 @@ class XYReco(AutoNameEnumBase):
     barycenter = auto()
     corona     = auto()
 
+class CutAlgo(AutoNameEnumBase):
+    threshold  = auto()
 
 class WfType(AutoNameEnumBase):
     rwf  = auto()

From 5a5677a92e0937171e58263940cf1ae7279d6548 Mon Sep 17 00:00:00 2001
From: casper <johnwaiton@hotmail.co.uk>
Date: Wed, 25 Feb 2026 18:39:13 +0000
Subject: [PATCH 10/69] modify `charge_threshold_method()`

Done for backwards compatibility, this function now mirrors that
of `calibrate_sipms()` before thresholding was removed.
---
 invisible_cities/reco/wfm_functions.py | 11 +++++++----
 1 file changed, 7 insertions(+), 4 deletions(-)

diff --git a/invisible_cities/reco/wfm_functions.py b/invisible_cities/reco/wfm_functions.py
index fb43889662..16b2dbc06c 100644
--- a/invisible_cities/reco/wfm_functions.py
+++ b/invisible_cities/reco/wfm_functions.py
@@ -6,11 +6,12 @@
 from typing  import Optional
 from typing  import Callable
 
-from .. core.core_functions import define_window
+from .. core.core_functions import define_window, to_col_vector
 from .. calib               import calib_sensors_functions as csf
 from .. sierpe              import blr
 from .. database            import load_db
 
+
 def to_adc(wfs, adc_to_pes):
     """
     Convert waveform in pes to adc.
@@ -166,10 +167,12 @@ def charge_threshold_method(wfs       : np.ndarray,
     
     Returns
     -------
-    Boolean numpy array of shape (n_sipms,) where True indicates that the SiPM is selected.
+    2D array of shape (n_sipms, n_time_bins) containing the waveforms of
+    each SiPM with all values below threshold set to zero.
     """
-    charges = np.sum(wfs, axis=1)
-    return charges >= threshold
+
+    thr = to_col_vector(np.full(wfs.shape[0], threshold))
+    return np.where(wfs > thr, wfs, 0)
 
 
 def top_n_method(wfs : np.ndarray, 

From 742ed5f4e52abd1fb154a135ffaa7b9f65d04f9e Mon Sep 17 00:00:00 2001
From: casper <johnwaiton@hotmail.co.uk>
Date: Wed, 25 Feb 2026 18:41:59 +0000
Subject: [PATCH 11/69] Implement `threshold_sipm_selection()`

---
 invisible_cities/cities/components.py | 17 +++++++++++++++++
 1 file changed, 17 insertions(+)

diff --git a/invisible_cities/cities/components.py b/invisible_cities/cities/components.py
index f32d2fdc7c..00d6fe9f15 100644
--- a/invisible_cities/cities/components.py
+++ b/invisible_cities/cities/components.py
@@ -815,6 +815,23 @@ def calibrate_sipms(rwf):
     return calibrate_sipms
 
 
+
+
+def threshold_sipm_selection(thr_sipm_type, thr_sipm, detector_db, run_number):
+    '''
+    Function that applies thresholding to the sipms in standard irene manner
+    '''
+
+    # extract sipm threshold
+    sipm_thr = get_actual_sipm_thr(thr_sipm_type, thr_sipm, detector_db, run_number)
+
+    def threshold_sipm_selection(wfs):
+        return wfm.charge_threshold_method(wfs, threshold = sipm_thr)
+
+    return threshold_sipm_selection
+
+
+
 def calibrate_with_mean(dbfile, run_number):
     DataSiPM   = load_db.DataSiPM(dbfile, run_number)
     adc_to_pes = np.abs(DataSiPM.adc_to_pes.values)

From 6bb7dca3b358529b5a1d0fd28b604d002c9b8773 Mon Sep 17 00:00:00 2001
From: casper <johnwaiton@hotmail.co.uk>
Date: Wed, 25 Feb 2026 18:42:30 +0000
Subject: [PATCH 12/69] Include `apply_cutting_function()`

This is the main controlling function, that decides which type of cut
will be applied.
---
 invisible_cities/cities/components.py | 13 +++++++++++++
 1 file changed, 13 insertions(+)

diff --git a/invisible_cities/cities/components.py b/invisible_cities/cities/components.py
index 00d6fe9f15..75d62d9128 100644
--- a/invisible_cities/cities/components.py
+++ b/invisible_cities/cities/components.py
@@ -84,6 +84,7 @@
 from .. types  .ic_types          import        types_dict_summary
 from .. types  .ic_types          import         types_dict_tracks
 from .. types  .symbols           import                    WfType
+from .. types  .symbols           import                   CutAlgo
 from .. types  .symbols           import               RebinMethod
 from .. types  .symbols           import                SiPMCharge
 from .. types  .symbols           import                   BlsMode
@@ -815,6 +816,18 @@ def calibrate_sipms(rwf):
     return calibrate_sipms
 
 
+def apply_cutting_function(algo, **cutting_params):
+
+    if algo is CutAlgo.threshold:
+        func = threshold_sipm_selection(**cutting_params)
+    else:
+        # temporary solution, think of a nicer method
+        func = threshold_sipm_selection(**cutting_params)
+
+    def apply_cutting_function(wfs):
+        return func(wfs)
+
+    return apply_cutting_function
 
 
 def threshold_sipm_selection(thr_sipm_type, thr_sipm, detector_db, run_number):

From 3660704555f3816889c11ec783aa875616afbd3b Mon Sep 17 00:00:00 2001
From: casper <johnwaiton@hotmail.co.uk>
Date: Wed, 25 Feb 2026 18:44:47 +0000
Subject: [PATCH 13/69] Implement `apply_cut()` into irene

---
 invisible_cities/cities/irene.py | 70 ++++++++++++++++++++------------
 1 file changed, 45 insertions(+), 25 deletions(-)

diff --git a/invisible_cities/cities/irene.py b/invisible_cities/cities/irene.py
index 49251e4bca..d1d7701535 100644
--- a/invisible_cities/cities/irene.py
+++ b/invisible_cities/cities/irene.py
@@ -25,6 +25,7 @@
 from .. io   .trigger_io       import       trigger_writer
 from .. io   .dst_io           import            df_writer
 from .. types.symbols          import WfType
+from .. types.symbols          import CutAlgo
 from .. types.symbols          import SiPMThreshold
 
 from .. database.load_db import DataPMT
@@ -46,36 +47,51 @@
 from .  components import wf_from_files
 from .  components import get_number_of_pmts
 from .  components import compute_and_write_pmaps
-from .  components import get_actual_sipm_thr
+from .  components import apply_cutting_function
 from .  components import sensor_masker
 
+from typing import Dict
+from typing import Any
 
 @city
-def irene( files_in        : OneOrManyFiles
-         , file_out        : str
-         , compression     : str
-         , event_range     : EventRangeType
-         , print_mod       : int
-         , detector_db     : str
-         , run_number      : int
-         , n_baseline      : int
-         , n_maw           : int
-         , thr_maw         : float
-         , thr_sipm        : float
-         , thr_sipm_type   : SiPMThreshold
-         , s1_lmin         : int  , s1_lmax      : int
-         , s1_tmin         : float, s1_tmax      : float
-         , s1_rebin_stride : int  , s1_stride    : int
-         , thr_csum_s1     : float
-         , s2_lmin         : int  , s2_lmax      : int
-         , s2_tmin         : float, s2_tmax      : float
-         , s2_rebin_stride : int  , s2_stride    : int
-         , thr_csum_s2     : float, thr_sipm_s2  : float
-         , pmt_samp_wid    : float, sipm_samp_wid: float
-         , store_db        : bool = True
+def irene( files_in         : OneOrManyFiles
+         , file_out         : str
+         , compression      : str
+         , event_range      : EventRangeType
+         , print_mod        : int
+         , detector_db      : str
+         , run_number       : int
+         , n_baseline       : int
+         , n_maw            : int
+         , thr_maw          : float
+         , thr_sipm         : float
+         , thr_sipm_type    : SiPMThreshold
+         , s1_lmin          : int  , s1_lmax      : int
+         , s1_tmin          : float, s1_tmax      : float
+         , s1_rebin_stride  : int  , s1_stride    : int
+         , thr_csum_s1      : float
+         , s2_lmin          : int  , s2_lmax      : int
+         , s2_tmin          : float, s2_tmax      : float
+         , s2_rebin_stride  : int  , s2_stride    : int
+         , thr_csum_s2      : float, thr_sipm_s2  : float
+         , pmt_samp_wid     : float, sipm_samp_wid: float
+         , store_db         : bool = True
+         , cutting_function : CutAlgo
+         , cutting_params   : Dict[str, Any]
          ):
+    '''
+    `cutting_function` is defined within components.py, and can vary, resulting
+    in the need for `cutting_params`, which are defined as such to allow for
+    the prior function to run. Currently implemented are `threshold_sipm_selection`,
+    with differing selections methods added soon.
 
-    sipm_thr = get_actual_sipm_thr(thr_sipm_type, thr_sipm, detector_db, run_number)
+
+    params for `threshold_sipm_selection`:
+        thr_sipm_type : SiPMThreshold
+        thr_sipm      : float
+        detector_db   : str
+        run_number    : int
+    '''
 
     #### Define data transformations
 
@@ -99,9 +115,13 @@ def irene( files_in        : OneOrManyFiles
                               args = ("cwf_sum", "cwf_sum_maw"),
                               out  = ("s1_indices", "s2_indices", "s2_energies"))
 
+
     # Remove baseline and calibrate SiPMs
     sipm_rwf_to_cal  = fl.map(calibrate_sipms(detector_db, run_number),
                               item = "sipm")
+    # apply function depending on user input, from provided list of functions
+    apply_cut        = fl.map(apply_cutting_function(cutting_function, **cutting_params),
+                              item = "sipm")
 
     event_count_in  = fl.spy_count()
     event_count_out = fl.spy_count()
@@ -125,7 +145,7 @@ def irene( files_in        : OneOrManyFiles
                                          s1_lmax, s1_lmin, s1_rebin_stride, s1_stride, s1_tmax, s1_tmin,
                                          s2_lmax, s2_lmin, s2_rebin_stride, s2_stride, s2_tmax, s2_tmin,
                                          thr_sipm_s2,
-                                         h5out, sipm_rwf_to_cal)
+                                         h5out, apply_cut, sipm_rwf_to_cal)
 
         result = push(source = wf_from_files(files_in, WfType.rwf),
                       pipe   = pipe(fl.slice(*event_range, close_all=True),

From 6aa1c3b99f4a2b4b3f0b22fede1b21a3f9786a7b Mon Sep 17 00:00:00 2001
From: casper <johnwaiton@hotmail.co.uk>
Date: Wed, 25 Feb 2026 18:45:20 +0000
Subject: [PATCH 14/69] Implement `apply_cut()` into
 `compute_and_write_pmaps()`

This included the adjustment of the cuts, applying
`select_wfs_above_time_integrated_thr()` and the zeroing integration cut
before the waveform rebinning.
---
 invisible_cities/cities/components.py   | 10 ++++----
 invisible_cities/cities/irene.py        |  4 ++--
 invisible_cities/reco/peak_functions.py | 32 +++++++++++++++----------
 3 files changed, 27 insertions(+), 19 deletions(-)

diff --git a/invisible_cities/cities/components.py b/invisible_cities/cities/components.py
index 75d62d9128..962300f93e 100644
--- a/invisible_cities/cities/components.py
+++ b/invisible_cities/cities/components.py
@@ -763,7 +763,7 @@ def sensor_data(path, wf_type):
 
 def build_pmap(detector_db, run_number, pmt_samp_wid, sipm_samp_wid,
                s1_lmax, s1_lmin, s1_rebin_stride, s1_stride, s1_tmax, s1_tmin,
-               s2_lmax, s2_lmin, s2_rebin_stride, s2_stride, s2_tmax, s2_tmin, thr_sipm_s2):
+               s2_lmax, s2_lmin, s2_rebin_stride, s2_stride, s2_tmax, s2_tmin, apply_cut):
     s1_params = dict(time        = minmax(min = s1_tmin,
                                           max = s1_tmax),
                     length       = minmax(min = s1_lmin,
@@ -784,8 +784,8 @@ def build_pmap(detector_db, run_number, pmt_samp_wid, sipm_samp_wid,
     sipm_ids = np.argwhere(datasipm.Active.values==1).flatten()
     def build_pmap(ccwf, s1_indx, s2_indx, sipmzs): # -> PMap
         return pkf.get_pmap(ccwf, s1_indx, s2_indx, sipmzs,
-                            s1_params, s2_params, thr_sipm_s2, pmt_ids, sipm_ids,
-                            pmt_samp_wid, sipm_samp_wid)
+                            s1_params, s2_params, pmt_ids, sipm_ids,
+                            pmt_samp_wid, sipm_samp_wid, apply_cut)
 
     return build_pmap
 
@@ -1261,7 +1261,7 @@ def integrate_wfs(wfs):
 def compute_and_write_pmaps(detector_db, run_number, pmt_samp_wid, sipm_samp_wid,
                   s1_lmax, s1_lmin, s1_rebin_stride, s1_stride, s1_tmax, s1_tmin,
                   s2_lmax, s2_lmin, s2_rebin_stride, s2_stride, s2_tmax, s2_tmin, thr_sipm_s2,
-                  h5out, sipm_rwf_to_cal=None):
+                  h5out, apply_cut, sipm_rwf_to_cal=None):
 
     # Filter events without signal over threshold
     indices_pass    = fl.map(check_nonempty_indices,
@@ -1272,7 +1272,7 @@ def compute_and_write_pmaps(detector_db, run_number, pmt_samp_wid, sipm_samp_wid
     # Build the PMap
     compute_pmap     = fl.map(build_pmap(detector_db, run_number, pmt_samp_wid, sipm_samp_wid,
                                          s1_lmax, s1_lmin, s1_rebin_stride, s1_stride, s1_tmax, s1_tmin,
-                                         s2_lmax, s2_lmin, s2_rebin_stride, s2_stride, s2_tmax, s2_tmin, thr_sipm_s2),
+                                         s2_lmax, s2_lmin, s2_rebin_stride, s2_stride, s2_tmax, s2_tmin, apply_cut),
                               args = ("ccwfs", "s1_indices", "s2_indices", "sipm"),
                               out  = "pmap")
 
diff --git a/invisible_cities/cities/irene.py b/invisible_cities/cities/irene.py
index d1d7701535..f0e04539ae 100644
--- a/invisible_cities/cities/irene.py
+++ b/invisible_cities/cities/irene.py
@@ -120,8 +120,8 @@ def irene( files_in         : OneOrManyFiles
     sipm_rwf_to_cal  = fl.map(calibrate_sipms(detector_db, run_number),
                               item = "sipm")
     # apply function depending on user input, from provided list of functions
-    apply_cut        = fl.map(apply_cutting_function(cutting_function, **cutting_params),
-                              item = "sipm")
+    apply_cut        = apply_cutting_function(cutting_function, **cutting_params)
+
 
     event_count_in  = fl.spy_count()
     event_count_out = fl.spy_count()
diff --git a/invisible_cities/reco/peak_functions.py b/invisible_cities/reco/peak_functions.py
index 30ebaba508..830cc02b58 100644
--- a/invisible_cities/reco/peak_functions.py
+++ b/invisible_cities/reco/peak_functions.py
@@ -76,13 +76,20 @@ def build_pmt_responses(indices, times, widths, ccwf,
 
 
 def build_sipm_responses(indices, times, widths,
-                         sipm_wfs, sipm_ids, rebin_stride, thr_sipm_s2):
-    _, _, sipm_wfs_ = pick_slice_and_rebin(indices , times, widths,
+                         sipm_wfs, sipm_ids, rebin_stride, apply_cut):
+
+    if apply_cut is not None:
+        # apply cut before slicing and rebinning
+        (sipm_idx,
+         sipm_wfs)   = apply_cut(sipm_wfs)
+    else:
+        # give all sipm ids as index if no cut is applied
+        sipm_idx     = np.arange(sipm_wfs.shape[0])
+    # rebin
+    _, _, sipm_wfs = pick_slice_and_rebin(indices , times, widths,
                                            sipm_wfs, rebin_stride,
                                            pad_zeros = False)
-    (sipm_idx,
-     sipm_wfs)   = select_wfs_above_time_integrated_thr(sipm_wfs_,
-                                                        thr_sipm_s2)
+
     return SiPMResponses(sipm_ids[sipm_idx], sipm_wfs)
 
 
@@ -93,7 +100,7 @@ def build_peak(indices, times,
                pmt_samp_wid  = 25 * units.ns,
                sipm_samp_wid =  1 * units.mus,
                sipm_wfs      = None,
-               thr_sipm_s2   = 0):
+               apply_cut     = None):
     sipm_pmt_bin_ratio = int(sipm_samp_wid/pmt_samp_wid)
     (pk_times ,
      pk_widths,
@@ -107,7 +114,7 @@ def build_peak(indices, times,
                                       widths * sipm_pmt_bin_ratio,
                                       sipm_wfs, sipm_ids,
                                       rebin_stride // sipm_pmt_bin_ratio,
-                                      thr_sipm_s2)
+                                      apply_cut)
     else:
         sipm_r = SiPMResponses.build_empty_instance()
 
@@ -120,7 +127,8 @@ def find_peaks(ccwfs, index,
                Pk, pmt_ids, sipm_ids=None,
                pmt_samp_wid = 25*units.ns,
                sipm_samp_wid = 1*units.mus,
-               sipm_wfs=None, thr_sipm_s2=0):
+               sipm_wfs=None, apply_cut = None):
+
     ccwfs = np.array(ccwfs, ndmin=2)
 
     peaks           = []
@@ -136,20 +144,20 @@ def find_peaks(ccwfs, index,
                         rebin_stride,
                         with_sipms, Pk,
                         pmt_samp_wid, sipm_samp_wid,
-                        sipm_wfs, thr_sipm_s2)
+                        sipm_wfs, apply_cut)
         peaks.append(pk)
     return peaks
 
 
 def get_pmap(ccwf, s1_indx, s2_indx, sipm_zs_wf,
-             s1_params, s2_params, thr_sipm_s2, pmt_ids, sipm_ids,
-             pmt_samp_wid, sipm_samp_wid):
+             s1_params, s2_params, pmt_ids, sipm_ids,
+             pmt_samp_wid, sipm_samp_wid, apply_cut = None):
     return PMap(find_peaks(ccwf, s1_indx, Pk=S1, pmt_ids=pmt_ids,
                            pmt_samp_wid=pmt_samp_wid,
                            **s1_params),
                 find_peaks(ccwf, s2_indx, Pk=S2, pmt_ids=pmt_ids, sipm_ids=sipm_ids,
                            sipm_wfs      = sipm_zs_wf,
-                           thr_sipm_s2   = thr_sipm_s2,
+                           apply_cut     = apply_cut,
                            pmt_samp_wid  = pmt_samp_wid,
                            sipm_samp_wid = sipm_samp_wid,
                            **s2_params))

From 0b3710c46e9ac7a1b88db09b9ace7322b471d737 Mon Sep 17 00:00:00 2001
From: casper <johnwaiton@hotmail.co.uk>
Date: Wed, 25 Feb 2026 18:45:52 +0000
Subject: [PATCH 15/69] Update default irene config

---
 invisible_cities/config/irene.conf | 6 ++++++
 1 file changed, 6 insertions(+)

diff --git a/invisible_cities/config/irene.conf b/invisible_cities/config/irene.conf
index 9b536c7c1c..f62f962c4b 100644
--- a/invisible_cities/config/irene.conf
+++ b/invisible_cities/config/irene.conf
@@ -52,3 +52,9 @@ thr_sipm_s2 = 10 * pes  # Threshold for the full sipm waveform
 
 pmt_samp_wid  = 25 * ns
 sipm_samp_wid = 1 * mus
+
+cutting_function = threshold 
+cutting_params   = dict(  thr_sipm_type = thr_sipm_type 
+                        , thr_sipm      = thr_sipm 
+                        , detector_db   = detector_db 
+                        , run_number    = run_number)

From 05dd34a1edb561fb348060e2678b8313623598c4 Mon Sep 17 00:00:00 2001
From: casper <johnwaiton@hotmail.co.uk>
Date: Wed, 25 Feb 2026 18:46:10 +0000
Subject: [PATCH 16/69] Remove whitespaces, add newline

---
 invisible_cities/reco/wfm_functions.py | 43 +++++++++++++-------------
 1 file changed, 22 insertions(+), 21 deletions(-)

diff --git a/invisible_cities/reco/wfm_functions.py b/invisible_cities/reco/wfm_functions.py
index 16b2dbc06c..a66297d9c2 100644
--- a/invisible_cities/reco/wfm_functions.py
+++ b/invisible_cities/reco/wfm_functions.py
@@ -135,17 +135,17 @@ def compare_cwf_blr(cwf, pmtblr, event_list, window_size=500):
     return np.array(DIFF)
 
 
-def median_std_method(wfs    : np.ndarray,  
+def median_std_method(wfs    : np.ndarray,
                       nsigma : Optional[float] = 3.) -> np.ndarray:
     """
-    Computes the median and standard deviation of the time summed SiPM 
+    Computes the median and standard deviation of the time summed SiPM
     waveforms and selects the SiPMs that are nsigma over the median.
 
     Parameters
     ----------
     wfs    : 2D array of shape (n_sipms, n_time_bins) containing the waveforms of each SiPM.
     nsigma : Number of standard deviations above the median, default 3.
-    
+
     Returns
     -------
     Boolean numpy array of shape (n_sipms,) where True indicates that the SiPM is selected.
@@ -159,12 +159,12 @@ def charge_threshold_method(wfs       : np.ndarray,
                             threshold : Optional[float] = 5.) -> np.ndarray:
     """
     Selects the SiPMs whose time summed waveforms are above a threshold.
-    
+
     Parameters
     ----------
     wfs       : 2D array of shape (n_sipms, n_time_bins) containing the waveforms of each SiPM.
     threshold : Charge threshold in PE, default 5.
-    
+
     Returns
     -------
     2D array of shape (n_sipms, n_time_bins) containing the waveforms of
@@ -175,16 +175,16 @@ def charge_threshold_method(wfs       : np.ndarray,
     return np.where(wfs > thr, wfs, 0)
 
 
-def top_n_method(wfs : np.ndarray, 
+def top_n_method(wfs : np.ndarray,
                  n   : Optional[int] = 10) -> np.ndarray:
     """
     Selects the SiPMs with the top n highest time summed waveforms.
-        
+
     Parameters
     ----------
     wfs       : 2D array of shape (n_sipms, n_time_bins) containing the waveforms of each SiPM.
     n         : Number of most energeticSiPMs to select, default 10.
-    
+
     Returns
     -------
     Boolean numpy array of shape (n_sipms,) where True indicates that the SiPM is selected.
@@ -198,12 +198,12 @@ def top_n_method(wfs : np.ndarray,
 
 
 def kill_isolated_sipms(selected_ids        : np.ndarray,
-                        sipm_x              : np.ndarray, 
-                        sipm_y              : np.ndarray, 
+                        sipm_x              : np.ndarray,
+                        sipm_y              : np.ndarray,
                         proximity_threshold : float) -> np.ndarray:
     """
-    For the SiPMs that have passed the previous selection, scans through the SiPMs to check if they 
-    have neighbouring SiPMs - i.e., within the proximity_threshold - that have also passed the selection. 
+    For the SiPMs that have passed the previous selection, scans through the SiPMs to check if they
+    have neighbouring SiPMs - i.e., within the proximity_threshold - that have also passed the selection.
     If no neighbours are found, the SiPMs are classed as isolated, and are removed.
 
     Parameters
@@ -228,16 +228,16 @@ def kill_isolated_sipms(selected_ids        : np.ndarray,
 
         if n_neighbors <= 1:
             selected_ids_no_isolated[i] = False
-            
+
     return selected_ids_no_isolated
 
 
 def apply_circular_padding(selected_ids_no_isolated : np.ndarray,
-                           sipm_x                   : np.ndarray, 
-                           sipm_y                   : np.ndarray, 
+                           sipm_x                   : np.ndarray,
+                           sipm_y                   : np.ndarray,
                            padding_radius           : float) -> np.ndarray:
     """
-    For the SiPMs that pass the previous selection, creates circular padding of radius padding_radius, 
+    For the SiPMs that pass the previous selection, creates circular padding of radius padding_radius,
     selecting all SiPMs within that radius. Stores the union of all selected SiPMs.
 
     Parameters
@@ -261,11 +261,11 @@ def apply_circular_padding(selected_ids_no_isolated : np.ndarray,
     return sipm_ids_with_signal
 
 
-def make_sipm_selection(wfs                 : np.ndarray, 
-                        selection_func      : Callable, 
-                        selection_kwargs    : dict, 
-                        proximity_threshold : float, 
-                        padding_radius      : float) -> np.ndarray: 
+def make_sipm_selection(wfs                 : np.ndarray,
+                        selection_func      : Callable,
+                        selection_kwargs    : dict,
+                        proximity_threshold : float,
+                        padding_radius      : float) -> np.ndarray:
     """
     SiPM selection pipeline, applies SiPM cuts based on user input.
     A first selection of SiPMs is made, isolated SiPMs are removed
@@ -303,3 +303,4 @@ def make_sipm_selection(wfs                 : np.ndarray,
         padding_radius
     )
     return sipm_ids_with_signal
+

From 56d7922e9feeb9c38bdd1361c4a5f6a87c382147 Mon Sep 17 00:00:00 2001
From: jwaiton <john.waiton@postgrad.manchester.ac.uk>
Date: Fri, 27 Feb 2026 11:06:35 +0000
Subject: [PATCH 17/69] Allow `calibrate_sipms()` to apply threshold.

This is to avoid possible situations where `calibrate_sipms()` is used
outwith irene. `thr` is now a keyword argument to ensure that if a
threshold is applied it is done so with intent.
---
 invisible_cities/calib/calib_sensors_functions.py     | 11 +++++++----
 .../calib/calib_sensors_functions_test.py             |  6 +++---
 2 files changed, 10 insertions(+), 7 deletions(-)

diff --git a/invisible_cities/calib/calib_sensors_functions.py b/invisible_cities/calib/calib_sensors_functions.py
index d5387bdb08..ba2c8854eb 100644
--- a/invisible_cities/calib/calib_sensors_functions.py
+++ b/invisible_cities/calib/calib_sensors_functions.py
@@ -136,15 +136,18 @@ def pmt_subtract_maw(cwfs, n_maw=100):
     return cwfs - maw
 
 
-def calibrate_sipms(sipm_wfs, adc_to_pes, *, bls_mode=BlsMode.mode):
+def calibrate_sipms(sipm_wfs, adc_to_pes, *, thr = None, bls_mode=BlsMode.mode):
     """
     Subtract baseline, and calibrates waveforms to pes.
     """
-    #thr  = to_col_vector(np.full(sipm_wfs.shape[0], thr))
     bls  = subtract_baseline(sipm_wfs, bls_mode=bls_mode)
     cwfs = calibrate_wfs(bls, adc_to_pes)
-    return cwfs
-    #return np.where(cwfs > thr, cwfs, 0)
+    if thr is None:
+        return cwfs
+    else:
+        # if you apply a threshold here, apply as usual
+        thr  = to_col_vector(np.full(sipm_wfs.shape[0], thr))
+        return np.where(cwfs > thr, cwfs, 0)
 
 
 def subtract_mean  (wfs): return subtract_baseline(wfs, bls_mode=BlsMode.mean  )
diff --git a/invisible_cities/calib/calib_sensors_functions_test.py b/invisible_cities/calib/calib_sensors_functions_test.py
index 452d22d90a..bdbfc673eb 100644
--- a/invisible_cities/calib/calib_sensors_functions_test.py
+++ b/invisible_cities/calib/calib_sensors_functions_test.py
@@ -207,7 +207,7 @@ def test_calibrate_sipms_stat(oscillating_waveform_with_baseline,
      baseline )              = oscillating_waveform_with_baseline
     #n_maw                    = n_samples // 500
 
-    ccwfs = csf.calibrate_sipms(wfs, adc_to_pes, nsigma * noise_sigma, bls_mode=BlsMode.mode)
+    ccwfs = csf.calibrate_sipms(wfs, adc_to_pes, thr=nsigma * noise_sigma, bls_mode=BlsMode.mode)
 
     number_of_zeros = np.count_nonzero(ccwfs == 0)
     assert number_of_zeros > fraction * ccwfs.size
@@ -219,7 +219,7 @@ def test_calibrate_sipms_common_threshold(toy_sipm_signal):
      common_threshold, _) = toy_sipm_signal
 
     zs_wf = csf.calibrate_sipms(signal_adc, adc_to_pes,
-                                common_threshold, bls_mode=BlsMode.mode)
+                                thr=common_threshold, bls_mode=BlsMode.mode)
 
     for actual, expected in zip(zs_wf, signal_zs_common_threshold):
         assert actual == approx(expected)
@@ -232,7 +232,7 @@ def test_calibrate_sipms_individual_thresholds(toy_sipm_signal):
 
 
     zs_wf = csf.calibrate_sipms(signal_adc, adc_to_pes,
-                                individual_thresholds,
+                                thr=individual_thresholds,
                                 bls_mode=BlsMode.mode)
     for actual, expected in zip(zs_wf, signal_zs_individual_thresholds):
         assert actual == approx(expected)

From 55ca25fbda1dc0215be240ecdfc71e576bf9cf51 Mon Sep 17 00:00:00 2001
From: jwaiton <john.waiton@postgrad.manchester.ac.uk>
Date: Fri, 27 Feb 2026 11:20:10 +0000
Subject: [PATCH 18/69] Include `thr_sipm_s2` into `cutting_params`

---
 invisible_cities/config/irene.conf | 1 +
 1 file changed, 1 insertion(+)

diff --git a/invisible_cities/config/irene.conf b/invisible_cities/config/irene.conf
index f62f962c4b..1f2312a78d 100644
--- a/invisible_cities/config/irene.conf
+++ b/invisible_cities/config/irene.conf
@@ -56,5 +56,6 @@ sipm_samp_wid = 1 * mus
 cutting_function = threshold 
 cutting_params   = dict(  thr_sipm_type = thr_sipm_type 
                         , thr_sipm      = thr_sipm 
+                        , thr_sipm_s2   = thr_sipm_s2
                         , detector_db   = detector_db 
                         , run_number    = run_number)

From f51919b4d1de40eea414cad09e87da8b63de6b1b Mon Sep 17 00:00:00 2001
From: jwaiton <john.waiton@postgrad.manchester.ac.uk>
Date: Fri, 27 Feb 2026 11:33:08 +0000
Subject: [PATCH 19/69] Implement `thr_sipm_s2` in `charge_threshold_method()`

---
 invisible_cities/cities/components.py  | 22 +++++++++++++++-------
 invisible_cities/reco/wfm_functions.py | 15 +++++++++++----
 2 files changed, 26 insertions(+), 11 deletions(-)

diff --git a/invisible_cities/cities/components.py b/invisible_cities/cities/components.py
index 962300f93e..7c2c0a99ed 100644
--- a/invisible_cities/cities/components.py
+++ b/invisible_cities/cities/components.py
@@ -830,16 +830,24 @@ def apply_cutting_function(wfs):
     return apply_cutting_function
 
 
-def threshold_sipm_selection(thr_sipm_type, thr_sipm, detector_db, run_number):
+def threshold_sipm_selection(thr_sipm_type
+                            , thr_sipm
+                            , thr_sipm_s2
+                            , run_number
+                            , detector_db = None):
     '''
-    Function that applies thresholding to the sipms in standard irene manner
+    Function that applies thresholding to the sipms in standard irene manner,
+    by zeroing all waveform values below a threshold.
     '''
+    # assume that if the detector_db is None, you return sipm threshold as the number provided
+    if detector_db is None:
+        sipm_thr = thr_sipm
+    else:
+        # extract sipm threshold
+        sipm_thr = get_actual_sipm_thr(thr_sipm_type, thr_sipm, detector_db, run_number)
 
-    # extract sipm threshold
-    sipm_thr = get_actual_sipm_thr(thr_sipm_type, thr_sipm, detector_db, run_number)
-
-    def threshold_sipm_selection(wfs):
-        return wfm.charge_threshold_method(wfs, threshold = sipm_thr)
+    def threshold_sipm_selection(wfs, indices):
+        return wfm.charge_threshold_method(wfs, indices, zeroing_thr = sipm_thr, integration_thr=thr_sipm_s2)
 
     return threshold_sipm_selection
 
diff --git a/invisible_cities/reco/wfm_functions.py b/invisible_cities/reco/wfm_functions.py
index a66297d9c2..9426651ab6 100644
--- a/invisible_cities/reco/wfm_functions.py
+++ b/invisible_cities/reco/wfm_functions.py
@@ -155,10 +155,12 @@ def median_std_method(wfs    : np.ndarray,
     return charges >= threshold
 
 
-def charge_threshold_method(wfs       : np.ndarray,
-                            threshold : Optional[float] = 5.) -> np.ndarray:
+def charge_threshold_method(wfs             : np.ndarray,
+                            indices         : np.ndarray,
+                            zeroing_thr     : Optional[float] = 2.,
+                            integration_thr : Optional[float] = 5.) -> np.ndarray:
     """
-    Selects the SiPMs whose time summed waveforms are above a threshold.
+    Selects the SiPMs whose time summed waveforms within the s2 windows are above a threshold.
 
     Parameters
     ----------
@@ -168,8 +170,13 @@ def charge_threshold_method(wfs       : np.ndarray,
     Returns
     -------
     2D array of shape (n_sipms, n_time_bins) containing the waveforms of
-    each SiPM with all values below threshold set to zero.
+    each SiPM with all values below threshold set to zero, with a mask applied depending
+    on summed waveforms over certain indices.
     """
+    thr = to_col_vector(np.full(wfs.shape[0], zeroing_thr))
+
+    # zero entries below threshold
+    zwfs = np.where(wfs > thr, wfs, 0)
 
     thr = to_col_vector(np.full(wfs.shape[0], threshold))
     return np.where(wfs > thr, wfs, 0)

From e3e631d1a2b1e7f5503b758cb1120f8ba7e2a0be Mon Sep 17 00:00:00 2001
From: jwaiton <john.waiton@postgrad.manchester.ac.uk>
Date: Fri, 27 Feb 2026 11:35:39 +0000
Subject: [PATCH 20/69] Implement integration cut based on slices

Initially this was done over the entire waveforms due to the reordering
of the cuts when applied in `build_sipm_responses()`. As suggested by
@Ian0sborne, it would be wise to move these cuts out of
`compute_and_write_pmaps()`
---
 invisible_cities/cities/components.py   |  4 ++--
 invisible_cities/reco/peak_functions.py | 15 +++++++++++++++
 invisible_cities/reco/wfm_functions.py  |  6 +++---
 3 files changed, 20 insertions(+), 5 deletions(-)

diff --git a/invisible_cities/cities/components.py b/invisible_cities/cities/components.py
index 7c2c0a99ed..2c01dc00d9 100644
--- a/invisible_cities/cities/components.py
+++ b/invisible_cities/cities/components.py
@@ -824,8 +824,8 @@ def apply_cutting_function(algo, **cutting_params):
         # temporary solution, think of a nicer method
         func = threshold_sipm_selection(**cutting_params)
 
-    def apply_cutting_function(wfs):
-        return func(wfs)
+    def apply_cutting_function(wfs, indices):
+        return func(wfs, indices)
 
     return apply_cutting_function
 
diff --git a/invisible_cities/reco/peak_functions.py b/invisible_cities/reco/peak_functions.py
index 830cc02b58..4ff795140e 100644
--- a/invisible_cities/reco/peak_functions.py
+++ b/invisible_cities/reco/peak_functions.py
@@ -25,9 +25,24 @@ def indices_and_wf_above_threshold(wf, thr):
     return ZsWf(indices_above_thr, wf_above_thr)
 
 
+def select_wf_slices_above_time_integrated_thr(wfs, indices, thr):
+    '''
+    function that integrates over certain time slices, and passes the waveform
+    based on the aforementioned slices passing the threshold
+    '''
+    slice_ = slice(indices[0], indices[-1] + 1)
+    wfs_   = wfs[:, slice_]
+
+    selected_ids = np.where(np.sum(wfs_, axis = 1) >= thr)[0]
+    selected_wfs = wfs[selected_ids]
+
+    return selected_ids, selected_wfs
+
+
 def select_wfs_above_time_integrated_thr(wfs, thr):
     selected_ids = np.where(np.sum(wfs, axis=1) >= thr)[0]
     selected_wfs = wfs[selected_ids]
+
     return selected_ids, selected_wfs
 
 
diff --git a/invisible_cities/reco/wfm_functions.py b/invisible_cities/reco/wfm_functions.py
index 9426651ab6..82015d0b7b 100644
--- a/invisible_cities/reco/wfm_functions.py
+++ b/invisible_cities/reco/wfm_functions.py
@@ -10,7 +10,7 @@
 from .. calib               import calib_sensors_functions as csf
 from .. sierpe              import blr
 from .. database            import load_db
-
+from .. reco.peak_functions import select_wf_slices_above_time_integrated_thr
 
 def to_adc(wfs, adc_to_pes):
     """
@@ -178,8 +178,8 @@ def charge_threshold_method(wfs             : np.ndarray,
     # zero entries below threshold
     zwfs = np.where(wfs > thr, wfs, 0)
 
-    thr = to_col_vector(np.full(wfs.shape[0], threshold))
-    return np.where(wfs > thr, wfs, 0)
+    # returns selected ids and waveforms above integral
+    return select_wf_slices_above_time_integrated_thr(zwfs, indices, integration_thr)
 
 
 def top_n_method(wfs : np.ndarray,

From 6c7846dfcf5bf5687ae32b60a1151dd69ce162bb Mon Sep 17 00:00:00 2001
From: jwaiton <john.waiton@postgrad.manchester.ac.uk>
Date: Fri, 27 Feb 2026 11:38:15 +0000
Subject: [PATCH 21/69] Adjust tests to account for `cut_params`

---
 invisible_cities/reco/peak_functions_test.py | 21 ++++++++++++++------
 1 file changed, 15 insertions(+), 6 deletions(-)

diff --git a/invisible_cities/reco/peak_functions_test.py b/invisible_cities/reco/peak_functions_test.py
index dfd2116b91..f3acf49ff2 100644
--- a/invisible_cities/reco/peak_functions_test.py
+++ b/invisible_cities/reco/peak_functions_test.py
@@ -13,6 +13,7 @@
 from hypothesis.strategies    import integers
 from hypothesis.extra.numpy   import arrays
 
+from ..cities.components      import apply_cutting_function
 from ..core.testing_utils   import exactly
 from ..core.testing_utils   import previous_float
 from ..core.testing_utils   import assert_Peak_equality
@@ -26,6 +27,8 @@
 from ..evm .pmaps           import PMap
 from ..io  .pmaps_io        import load_pmaps
 from ..types.ic_types       import minmax
+from ..types.symbols        import SiPMThreshold
+from ..types.symbols        import CutAlgo
 from .                      import peak_functions as pf
 
 
@@ -350,11 +353,19 @@ def test_build_sipm_responses(wf_with_indices):
     wfs_slice       = wfs[:, indices]
     peak_integrals  = wfs_slice.sum(axis=1)
     below_thr_index = np.argmin (peak_integrals)
-    # next_float doesn't work here
     thr             = peak_integrals[below_thr_index] * 1.000001
     sipm_ids        = np.arange(len(wfs))
+
+    # next_float doesn't work here
+    cut_params = dict(detector_db   = 'None',
+                  thr_sipm_s2   = thr,
+                  thr_sipm      = 0,
+                  thr_sipm_type = SiPMThreshold.common,
+                  run_number    = 0)
+    apply_cut  = apply_cutting_function(CutAlgo.threshold, **cut_params)
+
     sipm_r          = pf.build_sipm_responses(indices, times, widths,
-                                              wfs, sipm_ids, 1, thr)
+                                              wfs, sipm_ids, 1, apply_cut)
 
     expected_ids = np.delete(      ids, below_thr_index)
     expected_wfs = np.delete(wfs_slice, below_thr_index, axis=0)
@@ -397,7 +408,7 @@ def test_build_peak_development(pmt_and_sipm_wfs_with_indices,
                          with_sipms   = with_sipms,
                          Pk           = Pk,
                          sipm_wfs     = sipm_wfs,
-                         thr_sipm_s2  = -1)
+                         )
 
     assert_Peak_equality(peak, expected_peak)
 
@@ -474,8 +485,7 @@ def test_find_peaks_s2_style(pmt_and_sipm_wfs_with_indices):
                           time_range, length_range,
                           stride, rebin_stride,
                           S2, pmt_ids, sipm_ids,
-                          sipm_wfs    = sipm_wfs,
-                          thr_sipm_s2 = -1)
+                          sipm_wfs    = sipm_wfs)
 
     (rebinned_times,
      rebinned_widths,
@@ -503,7 +513,6 @@ def test_get_pmap(s1_and_s2_with_indices):
     sipm_samp_wid = 1  * units.mus
     pmap = pf.get_pmap(pmt_wfs, s1_indx, s2_indx, sipm_wfs,
                        s1_params, s2_params,
-                       thr_sipm_s2 = -1,
                        pmt_ids     = pmt_ids, sipm_ids = sipm_ids,
                        pmt_samp_wid  = pmt_samp_wid ,
                        sipm_samp_wid = sipm_samp_wid)

From 705c38c7cc0ad9a53b8758735cdcba3d8cefb3c4 Mon Sep 17 00:00:00 2001
From: jwaiton <john.waiton@postgrad.manchester.ac.uk>
Date: Fri, 27 Feb 2026 11:38:35 +0000
Subject: [PATCH 22/69] Implement modifiable cuts in hypathia

---
 invisible_cities/cities/hypathia.py   | 61 ++++++++++++++++-----------
 invisible_cities/config/hypathia.conf |  7 +++
 2 files changed, 43 insertions(+), 25 deletions(-)

diff --git a/invisible_cities/cities/hypathia.py b/invisible_cities/cities/hypathia.py
index 91f2abaeb7..024e51555a 100644
--- a/invisible_cities/cities/hypathia.py
+++ b/invisible_cities/cities/hypathia.py
@@ -28,6 +28,7 @@
 from .. io   .run_and_event_io import run_and_event_writer
 from .. io   .      trigger_io import       trigger_writer
 from .. types.symbols          import WfType
+from .. types.symbols          import CutAlgo
 from .. types.symbols          import SiPMThreshold
 
 from .. dataflow            import dataflow as fl
@@ -48,32 +49,38 @@
 from .  components import calibrate_sipms
 from .  components import get_actual_sipm_thr
 from .  components import sensor_masker
+from .  components import apply_cutting_function
+
+from typing import Dict
+from typing import Any
 
 
 @city
-def hypathia( files_in        : OneOrManyFiles
-            , file_out        : str
-            , compression     : str
-            , event_range     : EventRangeType
-            , print_mod       : int
-            , detector_db     : str
-            , run_number      : int
-            , sipm_noise_cut  : float
-            , filter_padding  : int
-            , thr_sipm        : float
-            , thr_sipm_type   : SiPMThreshold
-            , pmt_wfs_rebin   : int
-            , pmt_pe_rms      : float
-            , s1_lmin         : int  , s1_lmax     : int
-            , s1_tmin         : float, s1_tmax     : float
-            , s1_rebin_stride : int  , s1_stride   : int
-            , thr_csum_s1     : float
-            , s2_lmin         : int  , s2_lmax     : int
-            , s2_tmin         : float, s2_tmax     : float
-            , s2_rebin_stride : int  , s2_stride   : int
-            , thr_csum_s2     : float, thr_sipm_s2 : float
-            , pmt_samp_wid    : float
-            , sipm_samp_wid   : float
+def hypathia( files_in         : OneOrManyFiles
+            , file_out         : str
+            , compression      : str
+            , event_range      : EventRangeType
+            , print_mod        : int
+            , detector_db      : str
+            , run_number       : int
+            , sipm_noise_cut   : float
+            , filter_padding   : int
+            , thr_sipm         : float
+            , thr_sipm_type    : SiPMThreshold
+            , pmt_wfs_rebin    : int
+            , pmt_pe_rms       : float
+            , s1_lmin          : int  , s1_lmax     : int
+            , s1_tmin          : float, s1_tmax     : float
+            , s1_rebin_stride  : int  , s1_stride   : int
+            , thr_csum_s1      : float
+            , s2_lmin          : int  , s2_lmax     : int
+            , s2_tmin          : float, s2_tmax     : float
+            , s2_rebin_stride  : int  , s2_stride   : int
+            , thr_csum_s2      : float, thr_sipm_s2 : float
+            , pmt_samp_wid     : float
+            , sipm_samp_wid    : float
+            , cutting_function : CutAlgo
+            , cutting_params   : Dict[str, Any]
             ):
 
     sipm_thr = get_actual_sipm_thr(thr_sipm_type, thr_sipm, detector_db, run_number)
@@ -115,9 +122,13 @@ def hypathia( files_in        : OneOrManyFiles
                               item="sipm")
 
     # SiPMs calibration
-    sipm_rwf_to_cal  = fl.map(calibrate_sipms(detector_db, run_number, sipm_thr),
+    sipm_rwf_to_cal  = fl.map(calibrate_sipms(detector_db, run_number),
                               item = "sipm")
 
+    # apply function depending on user input, from provided list of functions
+    apply_cut        = apply_cutting_function(cutting_function, **cutting_params)
+
+
     event_count_in  = fl.spy_count()
     event_count_out = fl.spy_count()
 
@@ -137,7 +148,7 @@ def hypathia( files_in        : OneOrManyFiles
                                              detector_db, run_number, pmt_samp_wid, sipm_samp_wid,
                                              s1_lmax, s1_lmin, s1_rebin_stride, s1_stride, s1_tmax, s1_tmin,
                                              s2_lmax, s2_lmin, s2_rebin_stride, s2_stride, s2_tmax, s2_tmin, thr_sipm_s2,
-                                             h5out, sipm_rwf_to_cal)
+                                             h5out, apply_cut, sipm_rwf_to_cal)
 
         result = push(source = wf_from_files(files_in, WfType.mcrd),
                       pipe   = pipe(fl.slice(*event_range, close_all=True),
diff --git a/invisible_cities/config/hypathia.conf b/invisible_cities/config/hypathia.conf
index 3aafcd8e41..d3029f25f9 100644
--- a/invisible_cities/config/hypathia.conf
+++ b/invisible_cities/config/hypathia.conf
@@ -52,3 +52,10 @@ thr_sipm_s2 = 5 * pes  # Threshold for the full sipm waveform
 
 pmt_samp_wid  = 25 * ns
 sipm_samp_wid = 1 * mus
+
+cutting_function = threshold 
+cutting_params   = dict(  thr_sipm_type = thr_sipm_type 
+                        , thr_sipm      = thr_sipm 
+                        , thr_sipm_s2   = thr_sipm_s2
+                        , detector_db   = detector_db 
+                        , run_number    = run_number)

From 3cc372adc98158e399fcae6307d5457285569f5f Mon Sep 17 00:00:00 2001
From: jwaiton <john.waiton@postgrad.manchester.ac.uk>
Date: Fri, 27 Feb 2026 11:47:05 +0000
Subject: [PATCH 23/69] Update docstrings, type annotation

---
 invisible_cities/reco/wfm_functions.py | 11 ++++++-----
 1 file changed, 6 insertions(+), 5 deletions(-)

diff --git a/invisible_cities/reco/wfm_functions.py b/invisible_cities/reco/wfm_functions.py
index 82015d0b7b..b9130311b2 100644
--- a/invisible_cities/reco/wfm_functions.py
+++ b/invisible_cities/reco/wfm_functions.py
@@ -5,6 +5,7 @@
 import numpy as np
 from typing  import Optional
 from typing  import Callable
+from typing  import Tuple
 
 from .. core.core_functions import define_window, to_col_vector
 from .. calib               import calib_sensors_functions as csf
@@ -158,9 +159,11 @@ def median_std_method(wfs    : np.ndarray,
 def charge_threshold_method(wfs             : np.ndarray,
                             indices         : np.ndarray,
                             zeroing_thr     : Optional[float] = 2.,
-                            integration_thr : Optional[float] = 5.) -> np.ndarray:
+                            integration_thr : Optional[float] = 5.) -> Tuple[np.ndarray, np.ndarray]:
     """
-    Selects the SiPMs whose time summed waveforms within the s2 windows are above a threshold.
+    Selects the SiPMs whose time summed waveforms within the s2 windows are above two thresholds:
+        - initial zero suprresion threshold (setting values in each waveform below a value to 0)
+        - threshold over each selected integrated slice of the waveforms
 
     Parameters
     ----------
@@ -169,9 +172,7 @@ def charge_threshold_method(wfs             : np.ndarray,
 
     Returns
     -------
-    2D array of shape (n_sipms, n_time_bins) containing the waveforms of
-    each SiPM with all values below threshold set to zero, with a mask applied depending
-    on summed waveforms over certain indices.
+    Tuple of np arrays including all passing sipm ids and the corresponding waveforms
     """
     thr = to_col_vector(np.full(wfs.shape[0], zeroing_thr))
 

From 9f3602d740f55be148b1ac0c16d0866e0fbd5758 Mon Sep 17 00:00:00 2001
From: Ian0sborne <ian.osborne@student.manchester.ac.uk>
Date: Tue, 3 Mar 2026 15:55:12 +0000
Subject: [PATCH 24/69] add pyrrha cutting method as an option in `CutAlgo` and
 `apply_cutting_function`

---
 invisible_cities/cities/components.py | 3 ++-
 invisible_cities/types/symbols.py     | 1 +
 2 files changed, 3 insertions(+), 1 deletion(-)

diff --git a/invisible_cities/cities/components.py b/invisible_cities/cities/components.py
index 2c01dc00d9..9a5cb7266b 100644
--- a/invisible_cities/cities/components.py
+++ b/invisible_cities/cities/components.py
@@ -820,6 +820,8 @@ def apply_cutting_function(algo, **cutting_params):
 
     if algo is CutAlgo.threshold:
         func = threshold_sipm_selection(**cutting_params)
+    elif algo is CutAlgo.pyrrha:
+        func = pyrrha_sipm_selection(**cutting_params)
     else:
         # temporary solution, think of a nicer method
         func = threshold_sipm_selection(**cutting_params)
@@ -852,7 +854,6 @@ def threshold_sipm_selection(wfs, indices):
     return threshold_sipm_selection
 
 
-
 def calibrate_with_mean(dbfile, run_number):
     DataSiPM   = load_db.DataSiPM(dbfile, run_number)
     adc_to_pes = np.abs(DataSiPM.adc_to_pes.values)
diff --git a/invisible_cities/types/symbols.py b/invisible_cities/types/symbols.py
index 9b4ab17994..4dab7ece07 100644
--- a/invisible_cities/types/symbols.py
+++ b/invisible_cities/types/symbols.py
@@ -157,6 +157,7 @@ class XYReco(AutoNameEnumBase):
 
 class CutAlgo(AutoNameEnumBase):
     threshold  = auto()
+    pyrrha     = auto()
 
 class WfType(AutoNameEnumBase):
     rwf  = auto()

From 6a0724d275e76b1cd1183bc20e6c002d2f7798e0 Mon Sep 17 00:00:00 2001
From: Ian0sborne <ian.osborne@student.manchester.ac.uk>
Date: Tue, 3 Mar 2026 16:30:24 +0000
Subject: [PATCH 25/69] pass `load_db` parameters into `make_sipm_selection`
 function

---
 invisible_cities/reco/wfm_functions.py | 6 ++++--
 1 file changed, 4 insertions(+), 2 deletions(-)

diff --git a/invisible_cities/reco/wfm_functions.py b/invisible_cities/reco/wfm_functions.py
index b9130311b2..9ad067549d 100644
--- a/invisible_cities/reco/wfm_functions.py
+++ b/invisible_cities/reco/wfm_functions.py
@@ -273,7 +273,9 @@ def make_sipm_selection(wfs                 : np.ndarray,
                         selection_func      : Callable,
                         selection_kwargs    : dict,
                         proximity_threshold : float,
-                        padding_radius      : float) -> np.ndarray:
+                        padding_radius      : float,
+                        run_number          : int,
+                        detector_db         : str) -> np.ndarray:
     """
     SiPM selection pipeline, applies SiPM cuts based on user input.
     A first selection of SiPMs is made, isolated SiPMs are removed
@@ -291,7 +293,7 @@ def make_sipm_selection(wfs                 : np.ndarray,
     -------
     sipm_ids_with_signal : Array of shape (n_sipms,) with boolean values indicating which SiPMs are selected.
     """
-    detector_info = load_db.DataSiPM('next100', 0)
+    detector_info = load_db.DataSiPM(detector_db, run_number)
     sipm_x = np.array(detector_info.X)
     sipm_y = np.array(detector_info.Y)
 

From 61127bf60e4e648d0348916d6fd35d2e8e75fcdd Mon Sep 17 00:00:00 2001
From: Ian0sborne <ian.osborne@student.manchester.ac.uk>
Date: Tue, 3 Mar 2026 16:33:13 +0000
Subject: [PATCH 26/69] add documentation to `make_sipm_selection`

---
 invisible_cities/reco/wfm_functions.py | 2 ++
 1 file changed, 2 insertions(+)

diff --git a/invisible_cities/reco/wfm_functions.py b/invisible_cities/reco/wfm_functions.py
index 9ad067549d..7bf3eaef7b 100644
--- a/invisible_cities/reco/wfm_functions.py
+++ b/invisible_cities/reco/wfm_functions.py
@@ -288,6 +288,8 @@ def make_sipm_selection(wfs                 : np.ndarray,
     selection_kwargs    : Dictionary of arguments passed to selection_func.
     proximity_threshold : Threshold used to identify isolated SiPMs.
     padding_radius      : Radial padding added to each SiPM that passes the selections.
+    run_number          : Run number used to load the detector database.
+    detector_db         : Database used to load the detector geometry.
 
     Returns
     -------

From 7eb86187e37687929edac286685333666e13875a Mon Sep 17 00:00:00 2001
From: Ian0sborne <ian.osborne@student.manchester.ac.uk>
Date: Tue, 3 Mar 2026 16:58:33 +0000
Subject: [PATCH 27/69] implement `spatial_selection_method` (previously named
 `make_sipm_selection`) into `pyrrha_sipm_selection`

---
 invisible_cities/cities/components.py  | 22 ++++++++++++++++++++++
 invisible_cities/reco/wfm_functions.py | 16 ++++++++--------
 2 files changed, 30 insertions(+), 8 deletions(-)

diff --git a/invisible_cities/cities/components.py b/invisible_cities/cities/components.py
index 9a5cb7266b..f45ae28703 100644
--- a/invisible_cities/cities/components.py
+++ b/invisible_cities/cities/components.py
@@ -854,6 +854,28 @@ def threshold_sipm_selection(wfs, indices):
     return threshold_sipm_selection
 
 
+def pyrrha_sipm_selection(selection_function : Callable
+                         , selection_kwargs    : dict
+                         , proximity_threshold : float
+                         , padding_radius      : float
+                         , run_number : int
+                         , detector_db : str):
+    '''
+    Function that applies a generic selection function to the sipms, which can be used to 
+    implement a spatial SiPM selection method (called Pyrrha).
+    '''
+    def pyrrha_sipm_selection(wfs, indices):
+        return wfm.spatial_selection_method(wfs, 
+                                            selection_function, 
+                                            selection_kwargs, 
+                                            proximity_threshold,
+                                            padding_radius,
+                                            run_number,
+                                            detector_db)
+
+    return pyrrha_sipm_selection
+
+
 def calibrate_with_mean(dbfile, run_number):
     DataSiPM   = load_db.DataSiPM(dbfile, run_number)
     adc_to_pes = np.abs(DataSiPM.adc_to_pes.values)
diff --git a/invisible_cities/reco/wfm_functions.py b/invisible_cities/reco/wfm_functions.py
index 7bf3eaef7b..a272b3b085 100644
--- a/invisible_cities/reco/wfm_functions.py
+++ b/invisible_cities/reco/wfm_functions.py
@@ -269,15 +269,15 @@ def apply_circular_padding(selected_ids_no_isolated : np.ndarray,
     return sipm_ids_with_signal
 
 
-def make_sipm_selection(wfs                 : np.ndarray,
-                        selection_func      : Callable,
-                        selection_kwargs    : dict,
-                        proximity_threshold : float,
-                        padding_radius      : float,
-                        run_number          : int,
-                        detector_db         : str) -> np.ndarray:
+def spatial_selection_method(wfs                 : np.ndarray,
+                             selection_func      : Callable,
+                             selection_kwargs    : dict,
+                             proximity_threshold : float,
+                             padding_radius      : float,
+                             run_number          : int,
+                             detector_db         : str) -> np.ndarray:
     """
-    SiPM selection pipeline, applies SiPM cuts based on user input.
+    SiPM selection function, applies SiPM cuts based on user input.
     A first selection of SiPMs is made, isolated SiPMs are removed
     and padding is added around the SiPMs that are left.
 

From c9db07e081dd9660fa87accdb825e676210c63d1 Mon Sep 17 00:00:00 2001
From: Ian0sborne <ian.osborne@student.manchester.ac.uk>
Date: Wed, 4 Mar 2026 15:35:19 +0000
Subject: [PATCH 28/69] modify `spatial_selection_method` to return the ids of
 the selected SiPMs and their waveforms

---
 invisible_cities/reco/wfm_functions.py | 9 +++++++--
 1 file changed, 7 insertions(+), 2 deletions(-)

diff --git a/invisible_cities/reco/wfm_functions.py b/invisible_cities/reco/wfm_functions.py
index a272b3b085..7e5f75ba9b 100644
--- a/invisible_cities/reco/wfm_functions.py
+++ b/invisible_cities/reco/wfm_functions.py
@@ -293,7 +293,8 @@ def spatial_selection_method(wfs                 : np.ndarray,
 
     Returns
     -------
-    sipm_ids_with_signal : Array of shape (n_sipms,) with boolean values indicating which SiPMs are selected.
+    selected_ids : Array of shape (n_sipms,) containing the indices of the selected SiPMs.
+    selected_wfs : 2D array of shape (n_selected_sipms, n_time_bins) with the waveforms of the selected SiPMs.
     """
     detector_info = load_db.DataSiPM(detector_db, run_number)
     sipm_x = np.array(detector_info.X)
@@ -314,5 +315,9 @@ def spatial_selection_method(wfs                 : np.ndarray,
         sipm_y,
         padding_radius
     )
-    return sipm_ids_with_signal
+
+    selected_ids = np.where(sipm_ids_with_signal)[0]
+    selected_wfs = wfs[selected_ids]
+
+    return selected_ids, selected_wfs
 

From 8ebeeb3c92c9d4e031e7e50a750030dfb1de5e4b Mon Sep 17 00:00:00 2001
From: Ian0sborne <ian.osborne@student.manchester.ac.uk>
Date: Wed, 4 Mar 2026 17:34:19 +0000
Subject: [PATCH 29/69] implement waveform index slicing in
 `spatial_selection_method`

---
 invisible_cities/cities/components.py  |  1 +
 invisible_cities/reco/wfm_functions.py | 18 +++++++++++-------
 2 files changed, 12 insertions(+), 7 deletions(-)

diff --git a/invisible_cities/cities/components.py b/invisible_cities/cities/components.py
index f45ae28703..6e1cef3856 100644
--- a/invisible_cities/cities/components.py
+++ b/invisible_cities/cities/components.py
@@ -866,6 +866,7 @@ def pyrrha_sipm_selection(selection_function : Callable
     '''
     def pyrrha_sipm_selection(wfs, indices):
         return wfm.spatial_selection_method(wfs, 
+                                            indices,
                                             selection_function, 
                                             selection_kwargs, 
                                             proximity_threshold,
diff --git a/invisible_cities/reco/wfm_functions.py b/invisible_cities/reco/wfm_functions.py
index 7e5f75ba9b..22e8de9c15 100644
--- a/invisible_cities/reco/wfm_functions.py
+++ b/invisible_cities/reco/wfm_functions.py
@@ -270,6 +270,7 @@ def apply_circular_padding(selected_ids_no_isolated : np.ndarray,
 
 
 def spatial_selection_method(wfs                 : np.ndarray,
+                             indices         : np.ndarray,
                              selection_func      : Callable,
                              selection_kwargs    : dict,
                              proximity_threshold : float,
@@ -300,24 +301,27 @@ def spatial_selection_method(wfs                 : np.ndarray,
     sipm_x = np.array(detector_info.X)
     sipm_y = np.array(detector_info.Y)
 
-    selected_ids = selection_func(wfs, **selection_kwargs)
+    slice_ = slice(indices[0], indices[-1] + 1)
+    wfs_   = wfs[:, slice_]
 
-    selected_ids_no_isolated = kill_isolated_sipms(
-        selected_ids,
+    starting_ids_ = selection_func(wfs_, **selection_kwargs)
+
+    selected_ids_no_isolated_ = kill_isolated_sipms(
+        starting_ids_,
         sipm_x,
         sipm_y,
         proximity_threshold
     )
 
-    sipm_ids_with_signal = apply_circular_padding(
-        selected_ids_no_isolated,
+    sipm_ids_with_signal_ = apply_circular_padding(
+        selected_ids_no_isolated_,
         sipm_x,
         sipm_y,
         padding_radius
     )
 
-    selected_ids = np.where(sipm_ids_with_signal)[0]
-    selected_wfs = wfs[selected_ids]
+    selected_ids = np.where(sipm_ids_with_signal_)[0]
+    selected_wfs = wfs_[selected_ids]
 
     return selected_ids, selected_wfs
 

From 2ad15f0fe88f3ed17df6c97dbd2e8bab55e317a0 Mon Sep 17 00:00:00 2001
From: Ian0sborne <ian.osborne@student.manchester.ac.uk>
Date: Wed, 4 Mar 2026 17:34:52 +0000
Subject: [PATCH 30/69] ammend function to return sliced waveforms

---
 invisible_cities/reco/peak_functions.py | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/invisible_cities/reco/peak_functions.py b/invisible_cities/reco/peak_functions.py
index 4ff795140e..40b73424e4 100644
--- a/invisible_cities/reco/peak_functions.py
+++ b/invisible_cities/reco/peak_functions.py
@@ -34,7 +34,7 @@ def select_wf_slices_above_time_integrated_thr(wfs, indices, thr):
     wfs_   = wfs[:, slice_]
 
     selected_ids = np.where(np.sum(wfs_, axis = 1) >= thr)[0]
-    selected_wfs = wfs[selected_ids]
+    selected_wfs = wfs_[selected_ids]
 
     return selected_ids, selected_wfs
 

From b27f4a17bfa312c8d931a5a39830b7c6075f21c9 Mon Sep 17 00:00:00 2001
From: Ian0sborne <ian.osborne@student.manchester.ac.uk>
Date: Wed, 4 Mar 2026 17:41:00 +0000
Subject: [PATCH 31/69] revert back to full waveform

---
 invisible_cities/reco/peak_functions.py | 2 +-
 invisible_cities/reco/wfm_functions.py  | 2 +-
 2 files changed, 2 insertions(+), 2 deletions(-)

diff --git a/invisible_cities/reco/peak_functions.py b/invisible_cities/reco/peak_functions.py
index 40b73424e4..4ff795140e 100644
--- a/invisible_cities/reco/peak_functions.py
+++ b/invisible_cities/reco/peak_functions.py
@@ -34,7 +34,7 @@ def select_wf_slices_above_time_integrated_thr(wfs, indices, thr):
     wfs_   = wfs[:, slice_]
 
     selected_ids = np.where(np.sum(wfs_, axis = 1) >= thr)[0]
-    selected_wfs = wfs_[selected_ids]
+    selected_wfs = wfs[selected_ids]
 
     return selected_ids, selected_wfs
 
diff --git a/invisible_cities/reco/wfm_functions.py b/invisible_cities/reco/wfm_functions.py
index 22e8de9c15..0e9e879a6b 100644
--- a/invisible_cities/reco/wfm_functions.py
+++ b/invisible_cities/reco/wfm_functions.py
@@ -321,7 +321,7 @@ def spatial_selection_method(wfs                 : np.ndarray,
     )
 
     selected_ids = np.where(sipm_ids_with_signal_)[0]
-    selected_wfs = wfs_[selected_ids]
+    selected_wfs = wfs[selected_ids]
 
     return selected_ids, selected_wfs
 

From 88afc488aaa63b4393834b6baf4070f8a04472a4 Mon Sep 17 00:00:00 2001
From: Ian0sborne <ian.osborne@student.manchester.ac.uk>
Date: Thu, 5 Mar 2026 15:53:31 +0000
Subject: [PATCH 32/69] implement `SiPMSelectionMethod` class to call specific
 selection methods for pyrrha

---
 invisible_cities/types/symbols.py | 3 +++
 1 file changed, 3 insertions(+)

diff --git a/invisible_cities/types/symbols.py b/invisible_cities/types/symbols.py
index 4dab7ece07..098d563800 100644
--- a/invisible_cities/types/symbols.py
+++ b/invisible_cities/types/symbols.py
@@ -159,6 +159,9 @@ class CutAlgo(AutoNameEnumBase):
     threshold  = auto()
     pyrrha     = auto()
 
+class SiPMSelectionMethod(AutoNameEnumBase):
+    median_std_method = auto()
+
 class WfType(AutoNameEnumBase):
     rwf  = auto()
     mcrd = auto()

From 866c532b44b95aadee8ae4d1e04e88e6d11430b1 Mon Sep 17 00:00:00 2001
From: Ian0sborne <ian.osborne@student.manchester.ac.uk>
Date: Thu, 5 Mar 2026 15:54:17 +0000
Subject: [PATCH 33/69] modify `pyrrha_sipm_selection` and
 `spatial_selection_method` to use SiPMSelectionMethod class

---
 invisible_cities/cities/components.py  |  9 +++++----
 invisible_cities/reco/wfm_functions.py | 16 +++++++++++-----
 2 files changed, 16 insertions(+), 9 deletions(-)

diff --git a/invisible_cities/cities/components.py b/invisible_cities/cities/components.py
index 6e1cef3856..614f3d07fb 100644
--- a/invisible_cities/cities/components.py
+++ b/invisible_cities/cities/components.py
@@ -85,6 +85,7 @@
 from .. types  .ic_types          import         types_dict_tracks
 from .. types  .symbols           import                    WfType
 from .. types  .symbols           import                   CutAlgo
+from .. types  .symbols           import       SiPMSelectionMethod
 from .. types  .symbols           import               RebinMethod
 from .. types  .symbols           import                SiPMCharge
 from .. types  .symbols           import                   BlsMode
@@ -854,12 +855,12 @@ def threshold_sipm_selection(wfs, indices):
     return threshold_sipm_selection
 
 
-def pyrrha_sipm_selection(selection_function : Callable
+def pyrrha_sipm_selection(selection_method     : SiPMSelectionMethod
                          , selection_kwargs    : dict
                          , proximity_threshold : float
                          , padding_radius      : float
-                         , run_number : int
-                         , detector_db : str):
+                         , run_number          : int
+                         , detector_db         : str):
     '''
     Function that applies a generic selection function to the sipms, which can be used to 
     implement a spatial SiPM selection method (called Pyrrha).
@@ -867,7 +868,7 @@ def pyrrha_sipm_selection(selection_function : Callable
     def pyrrha_sipm_selection(wfs, indices):
         return wfm.spatial_selection_method(wfs, 
                                             indices,
-                                            selection_function, 
+                                            selection_method, 
                                             selection_kwargs, 
                                             proximity_threshold,
                                             padding_radius,
diff --git a/invisible_cities/reco/wfm_functions.py b/invisible_cities/reco/wfm_functions.py
index 0e9e879a6b..123aa88546 100644
--- a/invisible_cities/reco/wfm_functions.py
+++ b/invisible_cities/reco/wfm_functions.py
@@ -13,6 +13,8 @@
 from .. database            import load_db
 from .. reco.peak_functions import select_wf_slices_above_time_integrated_thr
 
+from .. types  .symbols     import SiPMSelectionMethod
+
 def to_adc(wfs, adc_to_pes):
     """
     Convert waveform in pes to adc.
@@ -270,8 +272,8 @@ def apply_circular_padding(selected_ids_no_isolated : np.ndarray,
 
 
 def spatial_selection_method(wfs                 : np.ndarray,
-                             indices         : np.ndarray,
-                             selection_func      : Callable,
+                             indices             : np.ndarray,
+                             selection_method    : SiPMSelectionMethod,
                              selection_kwargs    : dict,
                              proximity_threshold : float,
                              padding_radius      : float,
@@ -285,8 +287,8 @@ def spatial_selection_method(wfs                 : np.ndarray,
     Parameters
     ----------
     wfs                 : 2D array of shape (n_sipms, n_time_bins) containing the waveforms of each SiPM.
-    selection_func      : Function returning a boolean mask of energetic SiPMs.
-    selection_kwargs    : Dictionary of arguments passed to selection_func.
+    selection_method    : Method used to select SiPMs.
+    selection_kwargs    : Dictionary of arguments passed to the selection function.
     proximity_threshold : Threshold used to identify isolated SiPMs.
     padding_radius      : Radial padding added to each SiPM that passes the selections.
     run_number          : Run number used to load the detector database.
@@ -304,7 +306,11 @@ def spatial_selection_method(wfs                 : np.ndarray,
     slice_ = slice(indices[0], indices[-1] + 1)
     wfs_   = wfs[:, slice_]
 
-    starting_ids_ = selection_func(wfs_, **selection_kwargs)
+    if selection_method is SiPMSelectionMethod.median_std_method:
+        starting_ids_ = median_std_method(wfs_, **selection_kwargs)
+    else:
+        # temporary solution, think of a nicer method
+        starting_ids_ = median_std_method(wfs_, **selection_kwargs)
 
     selected_ids_no_isolated_ = kill_isolated_sipms(
         starting_ids_,

From 7f8812ead264b5e64ce44fee7d4df0d169ced5e9 Mon Sep 17 00:00:00 2001
From: Ian0sborne <ian.osborne@student.manchester.ac.uk>
Date: Thu, 5 Mar 2026 15:58:12 +0000
Subject: [PATCH 34/69] include another method for initial SiPM selection and
 raise an error for incorrect input

---
 invisible_cities/reco/wfm_functions.py | 5 +++--
 1 file changed, 3 insertions(+), 2 deletions(-)

diff --git a/invisible_cities/reco/wfm_functions.py b/invisible_cities/reco/wfm_functions.py
index 123aa88546..4c9cb70d2b 100644
--- a/invisible_cities/reco/wfm_functions.py
+++ b/invisible_cities/reco/wfm_functions.py
@@ -308,9 +308,10 @@ def spatial_selection_method(wfs                 : np.ndarray,
 
     if selection_method is SiPMSelectionMethod.median_std_method:
         starting_ids_ = median_std_method(wfs_, **selection_kwargs)
+    elif selection_method is SiPMSelectionMethod.top_n_method:
+        starting_ids_ = top_n_method(wfs_, **selection_kwargs)
     else:
-        # temporary solution, think of a nicer method
-        starting_ids_ = median_std_method(wfs_, **selection_kwargs)
+        raise ValueError(f"Selection method {selection_method} not recognized.")
 
     selected_ids_no_isolated_ = kill_isolated_sipms(
         starting_ids_,

From 7770a7101729045609f021286e908b3507c53379 Mon Sep 17 00:00:00 2001
From: Ian0sborne <ian.osborne@student.manchester.ac.uk>
Date: Tue, 10 Mar 2026 17:33:52 +0000
Subject: [PATCH 35/69] remove time index dependency of SiPM cutting functions

---
 invisible_cities/cities/components.py   | 11 +++++------
 invisible_cities/reco/peak_functions.py |  8 +++-----
 invisible_cities/reco/wfm_functions.py  | 13 ++++---------
 3 files changed, 12 insertions(+), 20 deletions(-)

diff --git a/invisible_cities/cities/components.py b/invisible_cities/cities/components.py
index 614f3d07fb..620264ebdc 100644
--- a/invisible_cities/cities/components.py
+++ b/invisible_cities/cities/components.py
@@ -827,8 +827,8 @@ def apply_cutting_function(algo, **cutting_params):
         # temporary solution, think of a nicer method
         func = threshold_sipm_selection(**cutting_params)
 
-    def apply_cutting_function(wfs, indices):
-        return func(wfs, indices)
+    def apply_cutting_function(wfs):
+        return func(wfs)
 
     return apply_cutting_function
 
@@ -849,8 +849,8 @@ def threshold_sipm_selection(thr_sipm_type
         # extract sipm threshold
         sipm_thr = get_actual_sipm_thr(thr_sipm_type, thr_sipm, detector_db, run_number)
 
-    def threshold_sipm_selection(wfs, indices):
-        return wfm.charge_threshold_method(wfs, indices, zeroing_thr = sipm_thr, integration_thr=thr_sipm_s2)
+    def threshold_sipm_selection(wfs):
+        return wfm.charge_threshold_method(wfs, zeroing_thr = sipm_thr, integration_thr=thr_sipm_s2)
 
     return threshold_sipm_selection
 
@@ -865,9 +865,8 @@ def pyrrha_sipm_selection(selection_method     : SiPMSelectionMethod
     Function that applies a generic selection function to the sipms, which can be used to 
     implement a spatial SiPM selection method (called Pyrrha).
     '''
-    def pyrrha_sipm_selection(wfs, indices):
+    def pyrrha_sipm_selection(wfs):
         return wfm.spatial_selection_method(wfs, 
-                                            indices,
                                             selection_method, 
                                             selection_kwargs, 
                                             proximity_threshold,
diff --git a/invisible_cities/reco/peak_functions.py b/invisible_cities/reco/peak_functions.py
index 4ff795140e..a1db39c5c8 100644
--- a/invisible_cities/reco/peak_functions.py
+++ b/invisible_cities/reco/peak_functions.py
@@ -92,7 +92,9 @@ def build_pmt_responses(indices, times, widths, ccwf,
 
 def build_sipm_responses(indices, times, widths,
                          sipm_wfs, sipm_ids, rebin_stride, apply_cut):
-
+    _, _, sipm_wfs = pick_slice_and_rebin(indices , times, widths,
+                                           sipm_wfs, rebin_stride,
+                                           pad_zeros = False)
     if apply_cut is not None:
         # apply cut before slicing and rebinning
         (sipm_idx,
@@ -100,10 +102,6 @@ def build_sipm_responses(indices, times, widths,
     else:
         # give all sipm ids as index if no cut is applied
         sipm_idx     = np.arange(sipm_wfs.shape[0])
-    # rebin
-    _, _, sipm_wfs = pick_slice_and_rebin(indices , times, widths,
-                                           sipm_wfs, rebin_stride,
-                                           pad_zeros = False)
 
     return SiPMResponses(sipm_ids[sipm_idx], sipm_wfs)
 
diff --git a/invisible_cities/reco/wfm_functions.py b/invisible_cities/reco/wfm_functions.py
index 4c9cb70d2b..53e8725af7 100644
--- a/invisible_cities/reco/wfm_functions.py
+++ b/invisible_cities/reco/wfm_functions.py
@@ -11,7 +11,7 @@
 from .. calib               import calib_sensors_functions as csf
 from .. sierpe              import blr
 from .. database            import load_db
-from .. reco.peak_functions import select_wf_slices_above_time_integrated_thr
+from .. reco.peak_functions import select_wfs_above_time_integrated_thr
 
 from .. types  .symbols     import SiPMSelectionMethod
 
@@ -159,7 +159,6 @@ def median_std_method(wfs    : np.ndarray,
 
 
 def charge_threshold_method(wfs             : np.ndarray,
-                            indices         : np.ndarray,
                             zeroing_thr     : Optional[float] = 2.,
                             integration_thr : Optional[float] = 5.) -> Tuple[np.ndarray, np.ndarray]:
     """
@@ -182,7 +181,7 @@ def charge_threshold_method(wfs             : np.ndarray,
     zwfs = np.where(wfs > thr, wfs, 0)
 
     # returns selected ids and waveforms above integral
-    return select_wf_slices_above_time_integrated_thr(zwfs, indices, integration_thr)
+    return select_wfs_above_time_integrated_thr(zwfs, integration_thr)
 
 
 def top_n_method(wfs : np.ndarray,
@@ -272,7 +271,6 @@ def apply_circular_padding(selected_ids_no_isolated : np.ndarray,
 
 
 def spatial_selection_method(wfs                 : np.ndarray,
-                             indices             : np.ndarray,
                              selection_method    : SiPMSelectionMethod,
                              selection_kwargs    : dict,
                              proximity_threshold : float,
@@ -303,13 +301,10 @@ def spatial_selection_method(wfs                 : np.ndarray,
     sipm_x = np.array(detector_info.X)
     sipm_y = np.array(detector_info.Y)
 
-    slice_ = slice(indices[0], indices[-1] + 1)
-    wfs_   = wfs[:, slice_]
-
     if selection_method is SiPMSelectionMethod.median_std_method:
-        starting_ids_ = median_std_method(wfs_, **selection_kwargs)
+        starting_ids_ = median_std_method(wfs, **selection_kwargs)
     elif selection_method is SiPMSelectionMethod.top_n_method:
-        starting_ids_ = top_n_method(wfs_, **selection_kwargs)
+        starting_ids_ = top_n_method(wfs, **selection_kwargs)
     else:
         raise ValueError(f"Selection method {selection_method} not recognized.")
 

From 83abad785053063fd112f1a5750648e4b77a499b Mon Sep 17 00:00:00 2001
From: Ian0sborne <ian.osborne@student.manchester.ac.uk>
Date: Tue, 10 Mar 2026 17:34:34 +0000
Subject: [PATCH 36/69] remove `select_wf_slices_above_time_integrated_thr`

---
 invisible_cities/reco/peak_functions.py | 14 --------------
 1 file changed, 14 deletions(-)

diff --git a/invisible_cities/reco/peak_functions.py b/invisible_cities/reco/peak_functions.py
index a1db39c5c8..e55376c28c 100644
--- a/invisible_cities/reco/peak_functions.py
+++ b/invisible_cities/reco/peak_functions.py
@@ -25,20 +25,6 @@ def indices_and_wf_above_threshold(wf, thr):
     return ZsWf(indices_above_thr, wf_above_thr)
 
 
-def select_wf_slices_above_time_integrated_thr(wfs, indices, thr):
-    '''
-    function that integrates over certain time slices, and passes the waveform
-    based on the aforementioned slices passing the threshold
-    '''
-    slice_ = slice(indices[0], indices[-1] + 1)
-    wfs_   = wfs[:, slice_]
-
-    selected_ids = np.where(np.sum(wfs_, axis = 1) >= thr)[0]
-    selected_wfs = wfs[selected_ids]
-
-    return selected_ids, selected_wfs
-
-
 def select_wfs_above_time_integrated_thr(wfs, thr):
     selected_ids = np.where(np.sum(wfs, axis=1) >= thr)[0]
     selected_wfs = wfs[selected_ids]

From c862256709ab4feae18bb8254f76327562a22630 Mon Sep 17 00:00:00 2001
From: Ian0sborne <ian.osborne@student.manchester.ac.uk>
Date: Tue, 10 Mar 2026 17:35:14 +0000
Subject: [PATCH 37/69] update docstrings

---
 invisible_cities/reco/wfm_functions.py | 5 +++--
 1 file changed, 3 insertions(+), 2 deletions(-)

diff --git a/invisible_cities/reco/wfm_functions.py b/invisible_cities/reco/wfm_functions.py
index 53e8725af7..a5977aeab7 100644
--- a/invisible_cities/reco/wfm_functions.py
+++ b/invisible_cities/reco/wfm_functions.py
@@ -168,8 +168,9 @@ def charge_threshold_method(wfs             : np.ndarray,
 
     Parameters
     ----------
-    wfs       : 2D array of shape (n_sipms, n_time_bins) containing the waveforms of each SiPM.
-    threshold : Charge threshold in PE, default 5.
+    wfs             : 2D array of shape (n_sipms, n_time_bins) containing the waveforms of each SiPM.
+    zeroing_thr     : Charge threshold for zero suppression in PE, default 2.
+    integration_thr : Charge threshold for total SiPM waveform in PE, default 5.
 
     Returns
     -------

From 1b4dbd7436e7578287688d7a2b4424ab2d6d1c39 Mon Sep 17 00:00:00 2001
From: Ian0sborne <ian.osborne@student.manchester.ac.uk>
Date: Thu, 7 May 2026 16:22:42 +0100
Subject: [PATCH 38/69] remove unused variable from hypathia

---
 invisible_cities/cities/hypathia.py | 2 --
 1 file changed, 2 deletions(-)

diff --git a/invisible_cities/cities/hypathia.py b/invisible_cities/cities/hypathia.py
index 024e51555a..6593b84c6b 100644
--- a/invisible_cities/cities/hypathia.py
+++ b/invisible_cities/cities/hypathia.py
@@ -83,8 +83,6 @@ def hypathia( files_in         : OneOrManyFiles
             , cutting_params   : Dict[str, Any]
             ):
 
-    sipm_thr = get_actual_sipm_thr(thr_sipm_type, thr_sipm, detector_db, run_number)
-
     #### Define data transformations
     sd = sensor_data(files_in[0], WfType.mcrd)
 

From d9659d5ad711fedec7da01797524e38ed035ba0b Mon Sep 17 00:00:00 2001
From: Ian0sborne <ian.osborne@student.manchester.ac.uk>
Date: Thu, 7 May 2026 16:53:57 +0100
Subject: [PATCH 39/69] raise error for incorrect cutting function input

---
 invisible_cities/cities/components.py | 3 +--
 1 file changed, 1 insertion(+), 2 deletions(-)

diff --git a/invisible_cities/cities/components.py b/invisible_cities/cities/components.py
index 620264ebdc..37b6b70dab 100644
--- a/invisible_cities/cities/components.py
+++ b/invisible_cities/cities/components.py
@@ -824,8 +824,7 @@ def apply_cutting_function(algo, **cutting_params):
     elif algo is CutAlgo.pyrrha:
         func = pyrrha_sipm_selection(**cutting_params)
     else:
-        # temporary solution, think of a nicer method
-        func = threshold_sipm_selection(**cutting_params)
+        raise ValueError(f"Unsupported cutting algorithm: {algo!r}. Expected one of {list(CutAlgo)}")
 
     def apply_cutting_function(wfs):
         return func(wfs)

From 31452002894b810b6486d7eb42eb982bbe66d657 Mon Sep 17 00:00:00 2001
From: Ian0sborne <ian.osborne@student.manchester.ac.uk>
Date: Thu, 7 May 2026 18:26:43 +0100
Subject: [PATCH 40/69] add fixture and test for `charge_threshold_method()`

---
 invisible_cities/reco/wfm_functions_test.py | 33 +++++++++++++++++++++
 1 file changed, 33 insertions(+)

diff --git a/invisible_cities/reco/wfm_functions_test.py b/invisible_cities/reco/wfm_functions_test.py
index 60c8a9abf8..a94635001b 100644
--- a/invisible_cities/reco/wfm_functions_test.py
+++ b/invisible_cities/reco/wfm_functions_test.py
@@ -4,6 +4,7 @@
 import tables as tb
 
 from pytest import mark
+from pytest import fixture
 
 from .. database import load_db
 
@@ -48,3 +49,35 @@ def test_compare_cwf_blr(dbnew, ICDATADIR):
                                    event_list=range(NEVT), window_size=300)
 
     assert max(diff) < 0.15
+
+@fixture
+def sipm_wfs_for_sipm_energy_selection_testing():
+    """
+    2D array of SiPM waveforms. All but two SiPMs have 100 PEs integrated charge.
+    Two outliers:
+    - SiPM 2: 150 PEs, more than 3 sigma above the median
+    - SiPM 7: 120 PEs, between 1-2 sigma above the median
+    """
+    n_sipms     = 10
+    n_time_bins = 100
+
+    wfs = np.ones((n_sipms, n_time_bins), dtype=np.float32)
+
+    wfs[2, :] = 1.5
+    wfs[7, :] = 1.2
+    # median ~ 100, std ~ 15.5
+    return wfs, [2, 7], [2]
+
+
+def test_median_std_method(sipm_wfs_for_sipm_energy_selection_testing):
+    """
+    Test function median_std_method(). The test asserts that the function correctly 
+    identifies the outliers based on different standard deviation thresholds. 
+    """
+    wfs, expected_outliers_1sigma, expected_outliers_3sigma = sipm_wfs_for_sipm_energy_selection_testing
+
+    passing_sipms_1sigma = np.where(wfm.median_std_method(wfs, nsigma=1))[0].tolist()
+    passing_sipms_3sigma = np.where(wfm.median_std_method(wfs, nsigma=3))[0].tolist()
+
+    assert passing_sipms_1sigma == expected_outliers_1sigma
+    assert passing_sipms_3sigma == expected_outliers_3sigma
\ No newline at end of file

From d3448e3d16059049b7acef5cd1995e355c9b1c19 Mon Sep 17 00:00:00 2001
From: Ian0sborne <ian.osborne@student.manchester.ac.uk>
Date: Thu, 7 May 2026 18:32:27 +0100
Subject: [PATCH 41/69] add test for `charge_threshold_method()`

---
 invisible_cities/reco/wfm_functions_test.py | 16 +++++++++++++++-
 1 file changed, 15 insertions(+), 1 deletion(-)

diff --git a/invisible_cities/reco/wfm_functions_test.py b/invisible_cities/reco/wfm_functions_test.py
index a94635001b..9d7407850a 100644
--- a/invisible_cities/reco/wfm_functions_test.py
+++ b/invisible_cities/reco/wfm_functions_test.py
@@ -80,4 +80,18 @@ def test_median_std_method(sipm_wfs_for_sipm_energy_selection_testing):
     passing_sipms_3sigma = np.where(wfm.median_std_method(wfs, nsigma=3))[0].tolist()
 
     assert passing_sipms_1sigma == expected_outliers_1sigma
-    assert passing_sipms_3sigma == expected_outliers_3sigma
\ No newline at end of file
+    assert passing_sipms_3sigma == expected_outliers_3sigma
+
+
+def test_threshold_method(sipm_wfs_for_sipm_energy_selection_testing):
+    """
+    Test function threshold_method(). The test asserts that the function correctly
+    kills the SiPMs below a certain charge threshold. 
+    """
+    wfs, expected_outliers_110pes, expected_outliers_130pes = sipm_wfs_for_sipm_energy_selection_testing
+
+    passing_sipms_110pes, _ = wfm.charge_threshold_method(wfs, zeroing_thr=0, integration_thr=110)
+    passing_sipms_130pes, _ = wfm.charge_threshold_method(wfs, zeroing_thr=0, integration_thr=130)
+
+    assert passing_sipms_110pes.tolist() == expected_outliers_110pes
+    assert passing_sipms_130pes.tolist() == expected_outliers_130pes
\ No newline at end of file

From 78444128d017a8e325ad53c12a30571001bf7cac Mon Sep 17 00:00:00 2001
From: Ian0sborne <ian.osborne@student.manchester.ac.uk>
Date: Sat, 9 May 2026 02:31:18 +0100
Subject: [PATCH 42/69] simplify naming

---
 invisible_cities/reco/wfm_functions_test.py | 10 +++++-----
 1 file changed, 5 insertions(+), 5 deletions(-)

diff --git a/invisible_cities/reco/wfm_functions_test.py b/invisible_cities/reco/wfm_functions_test.py
index 9d7407850a..a455fcc36f 100644
--- a/invisible_cities/reco/wfm_functions_test.py
+++ b/invisible_cities/reco/wfm_functions_test.py
@@ -51,7 +51,7 @@ def test_compare_cwf_blr(dbnew, ICDATADIR):
     assert max(diff) < 0.15
 
 @fixture
-def sipm_wfs_for_sipm_energy_selection_testing():
+def sipm_wfs_for_sipm_selection_testing():
     """
     2D array of SiPM waveforms. All but two SiPMs have 100 PEs integrated charge.
     Two outliers:
@@ -69,12 +69,12 @@ def sipm_wfs_for_sipm_energy_selection_testing():
     return wfs, [2, 7], [2]
 
 
-def test_median_std_method(sipm_wfs_for_sipm_energy_selection_testing):
+def test_median_std_method(sipm_wfs_for_sipm_selection_testing):
     """
     Test function median_std_method(). The test asserts that the function correctly 
     identifies the outliers based on different standard deviation thresholds. 
     """
-    wfs, expected_outliers_1sigma, expected_outliers_3sigma = sipm_wfs_for_sipm_energy_selection_testing
+    wfs, expected_outliers_1sigma, expected_outliers_3sigma = sipm_wfs_for_sipm_selection_testing
 
     passing_sipms_1sigma = np.where(wfm.median_std_method(wfs, nsigma=1))[0].tolist()
     passing_sipms_3sigma = np.where(wfm.median_std_method(wfs, nsigma=3))[0].tolist()
@@ -83,12 +83,12 @@ def test_median_std_method(sipm_wfs_for_sipm_energy_selection_testing):
     assert passing_sipms_3sigma == expected_outliers_3sigma
 
 
-def test_threshold_method(sipm_wfs_for_sipm_energy_selection_testing):
+def test_threshold_method(sipm_wfs_for_sipm_selection_testing):
     """
     Test function threshold_method(). The test asserts that the function correctly
     kills the SiPMs below a certain charge threshold. 
     """
-    wfs, expected_outliers_110pes, expected_outliers_130pes = sipm_wfs_for_sipm_energy_selection_testing
+    wfs, expected_outliers_110pes, expected_outliers_130pes = sipm_wfs_for_sipm_selection_testing
 
     passing_sipms_110pes, _ = wfm.charge_threshold_method(wfs, zeroing_thr=0, integration_thr=110)
     passing_sipms_130pes, _ = wfm.charge_threshold_method(wfs, zeroing_thr=0, integration_thr=130)

From f3d837445184a575cdbf5db9e216a0f71d0c2ca3 Mon Sep 17 00:00:00 2001
From: Ian0sborne <ian.osborne@student.manchester.ac.uk>
Date: Sat, 9 May 2026 02:37:39 +0100
Subject: [PATCH 43/69] add test for `top_n_method`

---
 invisible_cities/reco/wfm_functions_test.py | 19 ++++++++++++++++++-
 1 file changed, 18 insertions(+), 1 deletion(-)

diff --git a/invisible_cities/reco/wfm_functions_test.py b/invisible_cities/reco/wfm_functions_test.py
index a455fcc36f..4eea01b234 100644
--- a/invisible_cities/reco/wfm_functions_test.py
+++ b/invisible_cities/reco/wfm_functions_test.py
@@ -94,4 +94,21 @@ def test_threshold_method(sipm_wfs_for_sipm_selection_testing):
     passing_sipms_130pes, _ = wfm.charge_threshold_method(wfs, zeroing_thr=0, integration_thr=130)
 
     assert passing_sipms_110pes.tolist() == expected_outliers_110pes
-    assert passing_sipms_130pes.tolist() == expected_outliers_130pes
\ No newline at end of file
+    assert passing_sipms_130pes.tolist() == expected_outliers_130pes
+
+
+def test_top_n_method(sipm_wfs_for_sipm_selection_testing):
+    """
+    Test function top_n_method(). The test asserts that the function  selects the correct
+    number of SiPMs andcorrectly identifies the top N SiPMs based on their integrated charge. 
+    """
+    wfs, expected_outliers_top2, expected_outliers_top1 = sipm_wfs_for_sipm_selection_testing
+
+    passing_sipms_top2 = np.where(wfm.top_n_method(wfs, n=2))[0].tolist()
+    passing_sipms_top1 = np.where(wfm.top_n_method(wfs, n=1))[0].tolist()
+
+    assert len(passing_sipms_top2) == 2
+    assert len(passing_sipms_top1) == 1
+
+    assert passing_sipms_top2 == expected_outliers_top2
+    assert passing_sipms_top1 == expected_outliers_top1
\ No newline at end of file

From 82092feb5f58e2914b3876883b9f00cc2f5a60e5 Mon Sep 17 00:00:00 2001
From: Ian0sborne <ian.osborne@student.manchester.ac.uk>
Date: Sun, 10 May 2026 01:09:37 +0100
Subject: [PATCH 44/69] add fixture and test for `kill_isolated_sipms()`

---
 invisible_cities/reco/wfm_functions_test.py | 60 ++++++++++++++++++++-
 1 file changed, 59 insertions(+), 1 deletion(-)

diff --git a/invisible_cities/reco/wfm_functions_test.py b/invisible_cities/reco/wfm_functions_test.py
index 4eea01b234..cbb3478c80 100644
--- a/invisible_cities/reco/wfm_functions_test.py
+++ b/invisible_cities/reco/wfm_functions_test.py
@@ -111,4 +111,62 @@ def test_top_n_method(sipm_wfs_for_sipm_selection_testing):
     assert len(passing_sipms_top1) == 1
 
     assert passing_sipms_top2 == expected_outliers_top2
-    assert passing_sipms_top1 == expected_outliers_top1
\ No newline at end of file
+    assert passing_sipms_top1 == expected_outliers_top1
+
+
+@fixture
+def sipm_grid_for_isolation_testing():
+    """
+    5x5 grid of with 10mm spacing containing 4 SiPMs:
+    - SiPM  7 (x=20, y=10) has nearest neighbours 12 and 13
+    - SiPM 12 (x=20, y=20) has nearest neighbours 7 and 13
+    - SiPM 13 (x=30, y=20) has nearestneighbours 12 and 7
+    - SiPM 24 (x=40, y=40) has no nearest neighbours
+    With proximity_threshold=15mm SiPMs 7, 12, 13 survive and SiPM 24 is killed.
+    With proximity_threshold=5mm, all SiPMs are killed.
+    """
+    spacing = 10
+    xs = np.arange(5) * spacing  # [0, 10, 20, 30, 40]
+    ys = np.arange(5) * spacing
+
+    grid_x, grid_y = np.meshgrid(xs, ys)
+    sipm_x = grid_x.flatten().astype(np.float32)  # shape (25,)
+    sipm_y = grid_y.flatten().astype(np.float32)
+
+    # visual representation of the grid with SiPMs plotted as their IDs:
+    # X  X  X  X  24
+    # X  X  X  X  X
+    # X  X  12 13 X
+    # X  X  7  X  X
+    # X  X  X  X  X
+    selected_ids = np.zeros(25, dtype=bool)
+    cluster_ids  = [7, 12, 13]
+    isolated_id  = [24]
+    for idx in cluster_ids + isolated_id:
+        selected_ids[idx] = True
+
+    # diagonal SiPM distance ~ 14mm
+    # only SiPMs with a nearest neighbour pass assuming proximity_threshold=15mm
+    # SiPM 24 should be killed
+    expected_survivors_15mm = np.zeros(25, dtype=bool)
+    for idx in cluster_ids:
+        expected_survivors_15mm[idx] = True
+
+    # given a SiPM distance of 10mm, there should be no survivors with proximity_threshold=5mm
+    expected_survivors_5mm = np.zeros(25, dtype=bool)
+
+    return sipm_x, sipm_y, selected_ids, expected_survivors_15mm, expected_survivors_5mm
+
+
+def test_kill_isolated_sipms(sipm_grid_for_isolation_testing):
+    """"
+    Test function kill_isolated_sipms(). The test asserts that the function correctly
+    identifies and removes isolated SiPMs based on their proximity to other SiPMs.
+    """
+    sipm_x, sipm_y, selected_ids, expected_survivors_15mm, expected_survivors_5mm = sipm_grid_for_isolation_testing
+    
+    surviving_sipms_15mm = wfm.kill_isolated_sipms(selected_ids, sipm_x, sipm_y, proximity_threshold=15.0)
+    surviving_sipms_5mm = wfm.kill_isolated_sipms(selected_ids, sipm_x, sipm_y, proximity_threshold=5.0)
+
+    assert surviving_sipms_15mm.tolist() == expected_survivors_15mm.tolist()
+    assert surviving_sipms_5mm.tolist() == expected_survivors_5mm.tolist()
\ No newline at end of file

From cc158daa400e7542b7075e7869ac0e77fbc9bd79 Mon Sep 17 00:00:00 2001
From: Ian0sborne <ian.osborne@student.manchester.ac.uk>
Date: Sun, 10 May 2026 11:49:38 +0100
Subject: [PATCH 45/69] amend `apply_circular_padding()` to not kill SiPMs when
 padding is set to 0

---
 invisible_cities/reco/wfm_functions.py | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/invisible_cities/reco/wfm_functions.py b/invisible_cities/reco/wfm_functions.py
index a5977aeab7..e0576c6e12 100644
--- a/invisible_cities/reco/wfm_functions.py
+++ b/invisible_cities/reco/wfm_functions.py
@@ -266,7 +266,7 @@ def apply_circular_padding(selected_ids_no_isolated : np.ndarray,
     for i in np.where(selected_ids_no_isolated)[0]:
         x, y = sipm_x[i], sipm_y[i]
         distances = np.sqrt((sipm_x - x)**2 + (sipm_y - y)**2)
-        sipm_ids_with_signal |= distances < padding_radius
+        sipm_ids_with_signal |= distances <= padding_radius
 
     return sipm_ids_with_signal
 

From b4c3057e40c09ee44560f32062def6b309975157 Mon Sep 17 00:00:00 2001
From: Ian0sborne <ian.osborne@student.manchester.ac.uk>
Date: Sun, 10 May 2026 11:53:14 +0100
Subject: [PATCH 46/69] add test for `apply_circular_padding()` and modify
 fixture to accomodate for it

---
 invisible_cities/reco/wfm_functions_test.py | 45 ++++++++++++++++++---
 1 file changed, 40 insertions(+), 5 deletions(-)

diff --git a/invisible_cities/reco/wfm_functions_test.py b/invisible_cities/reco/wfm_functions_test.py
index cbb3478c80..5864d07bc3 100644
--- a/invisible_cities/reco/wfm_functions_test.py
+++ b/invisible_cities/reco/wfm_functions_test.py
@@ -115,7 +115,7 @@ def test_top_n_method(sipm_wfs_for_sipm_selection_testing):
 
 
 @fixture
-def sipm_grid_for_isolation_testing():
+def sipm_grid_for_isolation_and_padding_testing():
     """
     5x5 grid of with 10mm spacing containing 4 SiPMs:
     - SiPM  7 (x=20, y=10) has nearest neighbours 12 and 13
@@ -124,6 +124,9 @@ def sipm_grid_for_isolation_testing():
     - SiPM 24 (x=40, y=40) has no nearest neighbours
     With proximity_threshold=15mm SiPMs 7, 12, 13 survive and SiPM 24 is killed.
     With proximity_threshold=5mm, all SiPMs are killed.
+
+    With a padding of 12mm around the cluster of SiPMs 7, 12, 13 you also
+    include SiPMs 2, 6, 8, 11, 14, 17, 18.
     """
     spacing = 10
     xs = np.arange(5) * spacing  # [0, 10, 20, 30, 40]
@@ -155,18 +158,50 @@ def sipm_grid_for_isolation_testing():
     # given a SiPM distance of 10mm, there should be no survivors with proximity_threshold=5mm
     expected_survivors_5mm = np.zeros(25, dtype=bool)
 
-    return sipm_x, sipm_y, selected_ids, expected_survivors_15mm, expected_survivors_5mm
+    # adding a padding of 12mm around the cluster of SiPMs 7, 12, 13 would include 
+    # SiPMs 2, 6, 8, 11, 14, 17, 18
+    # X  X  X  X  X
+    # X  X  17 18 X
+    # X  11 12 13 14
+    # X  6  7  8  X
+    # X  X  2  X  X
+    padded_cluster_ids = cluster_ids + [2, 6, 8, 11, 14, 17, 18]
+    expected_survivors_padding12 = np.zeros(25, dtype=bool)
+    for idx in padded_cluster_ids:
+        expected_survivors_padding12[idx] = True
+
+    # adding a padding of 0mm around the cluster of SiPMs 7, 12, 13 would include only the cluster itself
+    expected_survivors_padding0 = np.zeros(25, dtype=bool)
+    for idx in cluster_ids:
+        expected_survivors_padding0[idx] = True
 
+    return (sipm_x, sipm_y, selected_ids, expected_survivors_15mm, expected_survivors_5mm, 
+            expected_survivors_padding12, expected_survivors_padding0)
 
-def test_kill_isolated_sipms(sipm_grid_for_isolation_testing):
+
+def test_kill_isolated_sipms(sipm_grid_for_isolation_and_padding_testing):
     """"
     Test function kill_isolated_sipms(). The test asserts that the function correctly
     identifies and removes isolated SiPMs based on their proximity to other SiPMs.
     """
-    sipm_x, sipm_y, selected_ids, expected_survivors_15mm, expected_survivors_5mm = sipm_grid_for_isolation_testing
+    sipm_x, sipm_y, selected_ids, expected_survivors_15mm, expected_survivors_5mm, _, _ = sipm_grid_for_isolation_and_padding_testing
     
     surviving_sipms_15mm = wfm.kill_isolated_sipms(selected_ids, sipm_x, sipm_y, proximity_threshold=15.0)
     surviving_sipms_5mm = wfm.kill_isolated_sipms(selected_ids, sipm_x, sipm_y, proximity_threshold=5.0)
 
     assert surviving_sipms_15mm.tolist() == expected_survivors_15mm.tolist()
-    assert surviving_sipms_5mm.tolist() == expected_survivors_5mm.tolist()
\ No newline at end of file
+    assert surviving_sipms_5mm.tolist() == expected_survivors_5mm.tolist()
+
+
+def test_apply_circular_padding(sipm_grid_for_isolation_and_padding_testing):
+    """
+    Test function apply_circular_padding(). The test asserts that the function correctly applies 
+    a circular padding around selected SiPMs to include neighboring SiPMs within the specified radius.
+    """
+    sipm_x, sipm_y, _, selected_ids, _, expected_survivors_padding12, expected_survivors_padding0 = sipm_grid_for_isolation_and_padding_testing
+
+    padded_sipms_12mm = wfm.apply_circular_padding(selected_ids, sipm_x, sipm_y, padding_radius=12.0)
+    padded_sipms_0mm = wfm.apply_circular_padding(selected_ids, sipm_x, sipm_y, padding_radius=0.0)
+
+    assert padded_sipms_12mm.tolist() == expected_survivors_padding12.tolist()
+    assert padded_sipms_0mm.tolist() == expected_survivors_padding0.tolist()
\ No newline at end of file

From 549e3de8c26d07cf8e7d3d313ff2cbb5f3925867 Mon Sep 17 00:00:00 2001
From: Ian0sborne <ian.osborne@student.manchester.ac.uk>
Date: Mon, 11 May 2026 16:14:57 +0100
Subject: [PATCH 47/69] add `zero_wfs_below_threshold()` to zero wavefrom
 entries below a threshold outside of `charge_threshold_method()`

---
 invisible_cities/reco/wfm_functions.py | 22 +++++++++++++++++++---
 1 file changed, 19 insertions(+), 3 deletions(-)

diff --git a/invisible_cities/reco/wfm_functions.py b/invisible_cities/reco/wfm_functions.py
index e0576c6e12..0f2feba267 100644
--- a/invisible_cities/reco/wfm_functions.py
+++ b/invisible_cities/reco/wfm_functions.py
@@ -138,6 +138,24 @@ def compare_cwf_blr(cwf, pmtblr, event_list, window_size=500):
     return np.array(DIFF)
 
 
+def zero_wfs_below_threshold(wfs         : np.ndarray,
+                             zeroing_thr : Optional[float] = 2.) -> np.ndarray:
+    """
+    Zeroes the entries of the input waveforms that are below a given threshold.
+
+    Parameters
+    ----------
+    wfs         : 2D array of shape (n_sipms, n_time_bins) containing the waveforms of each SiPM.
+    zeroing_thr : Charge threshold for zero suppression in PE, default 2.
+
+    Returns
+    -------
+    2D array of shape (n_sipms, n_time_bins) containing the input waveforms with entries below threshold set to zero.
+    """
+    thr = to_col_vector(np.full(wfs.shape[0], zeroing_thr))
+    return np.where(wfs > thr, wfs, 0)
+
+
 def median_std_method(wfs    : np.ndarray,
                       nsigma : Optional[float] = 3.) -> np.ndarray:
     """
@@ -176,10 +194,8 @@ def charge_threshold_method(wfs             : np.ndarray,
     -------
     Tuple of np arrays including all passing sipm ids and the corresponding waveforms
     """
-    thr = to_col_vector(np.full(wfs.shape[0], zeroing_thr))
-
     # zero entries below threshold
-    zwfs = np.where(wfs > thr, wfs, 0)
+    zwfs = zero_wfs_below_threshold(wfs, zeroing_thr)
 
     # returns selected ids and waveforms above integral
     return select_wfs_above_time_integrated_thr(zwfs, integration_thr)

From c63d9d159bcd69735d7d55532f70d78618eaddb1 Mon Sep 17 00:00:00 2001
From: Ian0sborne <ian.osborne@student.manchester.ac.uk>
Date: Mon, 11 May 2026 16:33:47 +0100
Subject: [PATCH 48/69] update `calibrate_sipms()` and corresponding test it to
 rely on `zero_wfs_below_threshold()` for zeroing waveform entries

---
 .../calib/calib_sensors_functions.py          |  9 ++----
 .../calib/calib_sensors_functions_test.py     | 28 ++++++++++---------
 2 files changed, 17 insertions(+), 20 deletions(-)

diff --git a/invisible_cities/calib/calib_sensors_functions.py b/invisible_cities/calib/calib_sensors_functions.py
index ba2c8854eb..38d522f205 100644
--- a/invisible_cities/calib/calib_sensors_functions.py
+++ b/invisible_cities/calib/calib_sensors_functions.py
@@ -136,18 +136,13 @@ def pmt_subtract_maw(cwfs, n_maw=100):
     return cwfs - maw
 
 
-def calibrate_sipms(sipm_wfs, adc_to_pes, *, thr = None, bls_mode=BlsMode.mode):
+def calibrate_sipms(sipm_wfs, adc_to_pes, *, bls_mode=BlsMode.mode):
     """
     Subtract baseline, and calibrates waveforms to pes.
     """
     bls  = subtract_baseline(sipm_wfs, bls_mode=bls_mode)
     cwfs = calibrate_wfs(bls, adc_to_pes)
-    if thr is None:
-        return cwfs
-    else:
-        # if you apply a threshold here, apply as usual
-        thr  = to_col_vector(np.full(sipm_wfs.shape[0], thr))
-        return np.where(cwfs > thr, cwfs, 0)
+    return cwfs
 
 
 def subtract_mean  (wfs): return subtract_baseline(wfs, bls_mode=BlsMode.mean  )
diff --git a/invisible_cities/calib/calib_sensors_functions_test.py b/invisible_cities/calib/calib_sensors_functions_test.py
index bdbfc673eb..1b39485619 100644
--- a/invisible_cities/calib/calib_sensors_functions_test.py
+++ b/invisible_cities/calib/calib_sensors_functions_test.py
@@ -13,6 +13,7 @@
 
 from .          import calib_sensors_functions as csf
 from .. core    import core_functions          as cf
+from .. reco    import wfm_functions           as wfm
 from .. sierpe  import fee                     as FE
 from .. sierpe  import waveform_generator      as wfg
 
@@ -207,10 +208,11 @@ def test_calibrate_sipms_stat(oscillating_waveform_with_baseline,
      baseline )              = oscillating_waveform_with_baseline
     #n_maw                    = n_samples // 500
 
-    ccwfs = csf.calibrate_sipms(wfs, adc_to_pes, thr=nsigma * noise_sigma, bls_mode=BlsMode.mode)
+    ccwfs = csf.calibrate_sipms(wfs, adc_to_pes, bls_mode=BlsMode.mode)
+    zeroed_ccwfs = wfm.zero_wfs_below_threshold(ccwfs, zeroing_thr=nsigma * noise_sigma)
 
-    number_of_zeros = np.count_nonzero(ccwfs == 0)
-    assert number_of_zeros > fraction * ccwfs.size
+    number_of_zeros = np.count_nonzero(zeroed_ccwfs == 0)
+    assert number_of_zeros > fraction * zeroed_ccwfs.size
 
 
 def test_calibrate_sipms_common_threshold(toy_sipm_signal):
@@ -218,8 +220,8 @@ def test_calibrate_sipms_common_threshold(toy_sipm_signal):
      signal_zs_common_threshold, _,
      common_threshold, _) = toy_sipm_signal
 
-    zs_wf = csf.calibrate_sipms(signal_adc, adc_to_pes,
-                                thr=common_threshold, bls_mode=BlsMode.mode)
+    ccwf = csf.calibrate_sipms(signal_adc, adc_to_pes, bls_mode=BlsMode.mode)
+    zs_wf = wfm.zero_wfs_below_threshold(ccwf, zeroing_thr=common_threshold)
 
     for actual, expected in zip(zs_wf, signal_zs_common_threshold):
         assert actual == approx(expected)
@@ -230,10 +232,9 @@ def test_calibrate_sipms_individual_thresholds(toy_sipm_signal):
      _, signal_zs_individual_thresholds,
      _, individual_thresholds) = toy_sipm_signal
 
+    ccwf = csf.calibrate_sipms(signal_adc, adc_to_pes, bls_mode=BlsMode.mode)
+    zs_wf = wfm.zero_wfs_below_threshold(ccwf, zeroing_thr=individual_thresholds)
 
-    zs_wf = csf.calibrate_sipms(signal_adc, adc_to_pes,
-                                thr=individual_thresholds,
-                                bls_mode=BlsMode.mode)
     for actual, expected in zip(zs_wf, signal_zs_individual_thresholds):
         assert actual == approx(expected)
 
@@ -309,11 +310,12 @@ def test_area_of_sum_equals_sum_of_areas_pmts(square_pmt_and_sipm_waveforms):
 
 def test_area_of_sum_equals_sum_of_areas_sipms(square_pmt_and_sipm_waveforms):
     _, nsensors, _, _, _, sipms_wfm, _ = square_pmt_and_sipm_waveforms
-    adc_to_pes = np.full(nsensors, 100, dtype=float)
-    cwfs       = csf.calibrate_sipms(sipms_wfm, adc_to_pes, thr=10, bls_mode=BlsMode.mode)
-    stot       = np.sum(cwfs[0]) * nsensors
-    sums       = np.sum(cwfs, axis=1)
-    stot2      = reduce(add, sums)
+    adc_to_pes  = np.full(nsensors, 100, dtype=float)
+    cwfs        = csf.calibrate_sipms(sipms_wfm, adc_to_pes, bls_mode=BlsMode.mode)
+    zeroed_cwfs = wfm.zero_wfs_below_threshold(cwfs, zeroing_thr=10)
+    stot        = np.sum(zeroed_cwfs[0]) * nsensors
+    sums        = np.sum(zeroed_cwfs, axis=1)
+    stot2       = reduce(add, sums)
     assert stot == approx(stot2, rel=1e-3)
 
 

From 952fffdab3498665da65cafdd3b6ff486cb79aaf Mon Sep 17 00:00:00 2001
From: Ian0sborne <ian.osborne@student.manchester.ac.uk>
Date: Mon, 11 May 2026 16:49:45 +0100
Subject: [PATCH 49/69] add test for `zero_wfs_below_threshold()` with small
 tweak in fixture

---
 invisible_cities/reco/wfm_functions_test.py | 21 +++++++++++++++++----
 1 file changed, 17 insertions(+), 4 deletions(-)

diff --git a/invisible_cities/reco/wfm_functions_test.py b/invisible_cities/reco/wfm_functions_test.py
index 5864d07bc3..50ac35b497 100644
--- a/invisible_cities/reco/wfm_functions_test.py
+++ b/invisible_cities/reco/wfm_functions_test.py
@@ -66,7 +66,20 @@ def sipm_wfs_for_sipm_selection_testing():
     wfs[2, :] = 1.5
     wfs[7, :] = 1.2
     # median ~ 100, std ~ 15.5
-    return wfs, [2, 7], [2]
+    return wfs, [2, 7], [2], [0, 1, 3, 4, 5, 6, 8, 9]
+
+
+def test_zero_wfs_below_threshold(sipm_wfs_for_sipm_selection_testing):
+    """
+    Test function zero_wfs_below_threshold(). The test asserts that the function correctly
+    sets to zero the entries of the waveforms that are below a specified threshold.
+    """
+    wfs, passing_wfs_11_ids, _, zeroed_wfs_11_ids  = sipm_wfs_for_sipm_selection_testing
+
+    zeroed_wfs_11 = wfm.zero_wfs_below_threshold(wfs, zeroing_thr=1.1)
+
+    assert np.all(zeroed_wfs_11[zeroed_wfs_11_ids] == 0)
+    assert np.all(zeroed_wfs_11[passing_wfs_11_ids] == wfs[passing_wfs_11_ids])
 
 
 def test_median_std_method(sipm_wfs_for_sipm_selection_testing):
@@ -74,7 +87,7 @@ def test_median_std_method(sipm_wfs_for_sipm_selection_testing):
     Test function median_std_method(). The test asserts that the function correctly 
     identifies the outliers based on different standard deviation thresholds. 
     """
-    wfs, expected_outliers_1sigma, expected_outliers_3sigma = sipm_wfs_for_sipm_selection_testing
+    wfs, expected_outliers_1sigma, expected_outliers_3sigma, _ = sipm_wfs_for_sipm_selection_testing
 
     passing_sipms_1sigma = np.where(wfm.median_std_method(wfs, nsigma=1))[0].tolist()
     passing_sipms_3sigma = np.where(wfm.median_std_method(wfs, nsigma=3))[0].tolist()
@@ -88,7 +101,7 @@ def test_threshold_method(sipm_wfs_for_sipm_selection_testing):
     Test function threshold_method(). The test asserts that the function correctly
     kills the SiPMs below a certain charge threshold. 
     """
-    wfs, expected_outliers_110pes, expected_outliers_130pes = sipm_wfs_for_sipm_selection_testing
+    wfs, expected_outliers_110pes, expected_outliers_130pes, _ = sipm_wfs_for_sipm_selection_testing
 
     passing_sipms_110pes, _ = wfm.charge_threshold_method(wfs, zeroing_thr=0, integration_thr=110)
     passing_sipms_130pes, _ = wfm.charge_threshold_method(wfs, zeroing_thr=0, integration_thr=130)
@@ -102,7 +115,7 @@ def test_top_n_method(sipm_wfs_for_sipm_selection_testing):
     Test function top_n_method(). The test asserts that the function  selects the correct
     number of SiPMs andcorrectly identifies the top N SiPMs based on their integrated charge. 
     """
-    wfs, expected_outliers_top2, expected_outliers_top1 = sipm_wfs_for_sipm_selection_testing
+    wfs, expected_outliers_top2, expected_outliers_top1, _ = sipm_wfs_for_sipm_selection_testing
 
     passing_sipms_top2 = np.where(wfm.top_n_method(wfs, n=2))[0].tolist()
     passing_sipms_top1 = np.where(wfm.top_n_method(wfs, n=1))[0].tolist()

From 64939503ac6cc3e93b47962e970767df5276da72 Mon Sep 17 00:00:00 2001
From: Ian0sborne <ian.osborne@student.manchester.ac.uk>
Date: Fri, 5 Jun 2026 11:24:44 +0100
Subject: [PATCH 50/69] Change comment position

---
 invisible_cities/reco/peak_functions_test.py | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/invisible_cities/reco/peak_functions_test.py b/invisible_cities/reco/peak_functions_test.py
index f3acf49ff2..a9f2fbe4c4 100644
--- a/invisible_cities/reco/peak_functions_test.py
+++ b/invisible_cities/reco/peak_functions_test.py
@@ -353,10 +353,10 @@ def test_build_sipm_responses(wf_with_indices):
     wfs_slice       = wfs[:, indices]
     peak_integrals  = wfs_slice.sum(axis=1)
     below_thr_index = np.argmin (peak_integrals)
+    # next_float doesn't work here
     thr             = peak_integrals[below_thr_index] * 1.000001
     sipm_ids        = np.arange(len(wfs))
 
-    # next_float doesn't work here
     cut_params = dict(detector_db   = 'None',
                   thr_sipm_s2   = thr,
                   thr_sipm      = 0,

From f5a939a5dbecfa3550387956c1181f853fecf88e Mon Sep 17 00:00:00 2001
From: Ian0sborne <ian.osborne@student.manchester.ac.uk>
Date: Fri, 5 Jun 2026 11:30:04 +0100
Subject: [PATCH 51/69] Update `SiPMCalibMode` name to match new functionality

---
 invisible_cities/calib/calib_sensors_functions.py | 12 ++++++------
 invisible_cities/types/symbols.py                 | 12 ++++++------
 2 files changed, 12 insertions(+), 12 deletions(-)

diff --git a/invisible_cities/calib/calib_sensors_functions.py b/invisible_cities/calib/calib_sensors_functions.py
index 38d522f205..f98ea569f5 100644
--- a/invisible_cities/calib/calib_sensors_functions.py
+++ b/invisible_cities/calib/calib_sensors_functions.py
@@ -157,10 +157,10 @@ def sipm_subtract_median_and_calibrate(sipm_wfs, adc_to_pes): return calibrate_w
 
 # Dict of functions for SiPM processing
 sipm_processing = {
-    SiPMCalibMode.subtract_mode            :      subtract_mode                ,# For gain extraction
-    SiPMCalibMode.subtract_median          :      subtract_median              ,# For gain extraction
-    SiPMCalibMode.subtract_mode_calibrate  : sipm_subtract_mode_and_calibrate  ,# For PDF calculation
-    SiPMCalibMode.subtract_mean_calibrate  : sipm_subtract_mean_and_calibrate  ,# For PDF calculation
-    SiPMCalibMode.subtract_median_calibrate: sipm_subtract_median_and_calibrate,# For PDF calculation
-    SiPMCalibMode.subtract_mode_zs         : calibrate_sipms                    # For data processing
+    SiPMCalibMode.subtract_mode              :      subtract_mode                ,# For gain extraction
+    SiPMCalibMode.subtract_median            :      subtract_median              ,# For gain extraction
+    SiPMCalibMode.subtract_mode_calibrate    : sipm_subtract_mode_and_calibrate  ,# For PDF calculation
+    SiPMCalibMode.subtract_mean_calibrate    : sipm_subtract_mean_and_calibrate  ,# For PDF calculation
+    SiPMCalibMode.subtract_median_calibrate  : sipm_subtract_median_and_calibrate,# For PDF calculation
+    SiPMCalibMode.subtract_baseline_calibrate: calibrate_sipms                    # For data processing
 }
diff --git a/invisible_cities/types/symbols.py b/invisible_cities/types/symbols.py
index 098d563800..62fdfe6a16 100644
--- a/invisible_cities/types/symbols.py
+++ b/invisible_cities/types/symbols.py
@@ -132,12 +132,12 @@ class SensorType(AutoNameEnumBase):
 
 
 class SiPMCalibMode(AutoNameEnumBase):
-    subtract_mode             = auto()
-    subtract_median           = auto()
-    subtract_mode_calibrate   = auto()
-    subtract_mean_calibrate   = auto()
-    subtract_median_calibrate = auto()
-    subtract_mode_zs          = auto()
+    subtract_mode                        = auto()
+    subtract_median                      = auto()
+    subtract_mode_calibrate              = auto()
+    subtract_mean_calibrate              = auto()
+    subtract_median_calibrate            = auto()
+    subtract_baseline_calibrate          = auto()
 
 
 class SiPMCharge(AutoNameEnumBase):

From 4ebae1ea9359c5f93f7c335d333574f2df40772b Mon Sep 17 00:00:00 2001
From: Ian0sborne <ian.osborne@student.manchester.ac.uk>
Date: Fri, 5 Jun 2026 11:40:04 +0100
Subject: [PATCH 52/69] Remove default parameters in `wfm_functions.py`

---
 invisible_cities/reco/wfm_functions.py | 8 ++++----
 1 file changed, 4 insertions(+), 4 deletions(-)

diff --git a/invisible_cities/reco/wfm_functions.py b/invisible_cities/reco/wfm_functions.py
index 0f2feba267..a230b61c5b 100644
--- a/invisible_cities/reco/wfm_functions.py
+++ b/invisible_cities/reco/wfm_functions.py
@@ -139,7 +139,7 @@ def compare_cwf_blr(cwf, pmtblr, event_list, window_size=500):
 
 
 def zero_wfs_below_threshold(wfs         : np.ndarray,
-                             zeroing_thr : Optional[float] = 2.) -> np.ndarray:
+                             zeroing_thr : float) -> np.ndarray:
     """
     Zeroes the entries of the input waveforms that are below a given threshold.
 
@@ -177,8 +177,8 @@ def median_std_method(wfs    : np.ndarray,
 
 
 def charge_threshold_method(wfs             : np.ndarray,
-                            zeroing_thr     : Optional[float] = 2.,
-                            integration_thr : Optional[float] = 5.) -> Tuple[np.ndarray, np.ndarray]:
+                            zeroing_thr     : float,
+                            integration_thr : float) -> Tuple[np.ndarray, np.ndarray]:
     """
     Selects the SiPMs whose time summed waveforms within the s2 windows are above two thresholds:
         - initial zero suprresion threshold (setting values in each waveform below a value to 0)
@@ -202,7 +202,7 @@ def charge_threshold_method(wfs             : np.ndarray,
 
 
 def top_n_method(wfs : np.ndarray,
-                 n   : Optional[int] = 10) -> np.ndarray:
+                 n   : int) -> np.ndarray:
     """
     Selects the SiPMs with the top n highest time summed waveforms.
 

From b89561298ec69db76cfc9dc23786e44f28b7a061 Mon Sep 17 00:00:00 2001
From: Ian0sborne <ian.osborne@student.manchester.ac.uk>
Date: Fri, 5 Jun 2026 12:00:12 +0100
Subject: [PATCH 53/69] Update & improve docstrings

---
 invisible_cities/reco/wfm_functions.py | 28 ++++++++++++++++----------
 1 file changed, 17 insertions(+), 11 deletions(-)

diff --git a/invisible_cities/reco/wfm_functions.py b/invisible_cities/reco/wfm_functions.py
index a230b61c5b..2f09b9df5e 100644
--- a/invisible_cities/reco/wfm_functions.py
+++ b/invisible_cities/reco/wfm_functions.py
@@ -146,7 +146,7 @@ def zero_wfs_below_threshold(wfs         : np.ndarray,
     Parameters
     ----------
     wfs         : 2D array of shape (n_sipms, n_time_bins) containing the waveforms of each SiPM.
-    zeroing_thr : Charge threshold for zero suppression in PE, default 2.
+    zeroing_thr : Charge threshold for zero suppression in PE.
 
     Returns
     -------
@@ -187,8 +187,8 @@ def charge_threshold_method(wfs             : np.ndarray,
     Parameters
     ----------
     wfs             : 2D array of shape (n_sipms, n_time_bins) containing the waveforms of each SiPM.
-    zeroing_thr     : Charge threshold for zero suppression in PE, default 2.
-    integration_thr : Charge threshold for total SiPM waveform in PE, default 5.
+    zeroing_thr     : Charge threshold for zero suppression in PE.
+    integration_thr : Charge threshold for total SiPM waveform in PE.
 
     Returns
     -------
@@ -209,7 +209,7 @@ def top_n_method(wfs : np.ndarray,
     Parameters
     ----------
     wfs       : 2D array of shape (n_sipms, n_time_bins) containing the waveforms of each SiPM.
-    n         : Number of most energeticSiPMs to select, default 10.
+    n         : Number of most energeticSiPMs to select.
 
     Returns
     -------
@@ -228,13 +228,17 @@ def kill_isolated_sipms(selected_ids        : np.ndarray,
                         sipm_y              : np.ndarray,
                         proximity_threshold : float) -> np.ndarray:
     """
-    For the SiPMs that have passed the previous selection, scans through the SiPMs to check if they
-    have neighbouring SiPMs - i.e., within the proximity_threshold - that have also passed the selection.
-    If no neighbours are found, the SiPMs are classed as isolated, and are removed.
+    Receives a list of SiPM IDs corresponding to the SiPMs with the most signficant 
+    energy depositions. Scans through these SiPMs to check if they have neighbouring 
+    SiPMs - i.e., within the proximity_threshold - that are also in the initial list
+    of energetic SiPMs. If no neighbours are found, the SiPMs are classed as isolated, 
+    and are removed. Outputs a list of SiPM IDs representing the most energetic SiPMs 
+    belonging to a cluster, which generally occurs in the region where the event took 
+    place.
 
     Parameters
     ----------
-    selected_ids        : Boolean array of shape (n_sipms,) indicating which SiPMs passed the previous selection.
+    selected_ids        : Boolean array of shape (n_sipms,) corresponding to the most energetic SiPMs.
     sipm_x              : 1D array of shape (n_sipms,) containing the x positions of the SiPMs.
     sipm_y              : 1D array of shape (n_sipms,) containing the y positions of the SiPMs.
     proximity_threshold : Distance threshold in mm used to identify isolated SiPMs.
@@ -263,12 +267,14 @@ def apply_circular_padding(selected_ids_no_isolated : np.ndarray,
                            sipm_y                   : np.ndarray,
                            padding_radius           : float) -> np.ndarray:
     """
-    For the SiPMs that pass the previous selection, creates circular padding of radius padding_radius,
-    selecting all SiPMs within that radius. Stores the union of all selected SiPMs.
+    Receives a list of SiPM IDs corresponding to the most energetic SiPMs clustered 
+    near the event. For these SiPMs, creates circular padding of radius padding_radius, 
+    selecting all SiPMs within that radius. Stores the union of all selected SiPMs. 
+    Outputs the SiPM IDs which are relevant for a given event.
 
     Parameters
     ----------
-    selected_ids_no_isolated : Boolean array of shape (n_sipms,) indicating which SiPMs passed the previous selection.
+    selected_ids_no_isolated : Boolean array of shape (n_sipms,) corresponding to the "relevant" SiPMs.
     sipm_x                   : 1D array of shape (n_sipms,) containing the x positions of the SiPMs.
     sipm_y                   : 1D array of shape (n_sipms,) containing the y positions of the SiPMs.
     padding_radius           : Distance threshold in mm used to create circular padding around selected SiPMs.

From 0780c5c2cf91a98620914712c6628d087a8e2208 Mon Sep 17 00:00:00 2001
From: Ian0sborne <ian.osborne@student.manchester.ac.uk>
Date: Mon, 8 Jun 2026 12:01:36 +0100
Subject: [PATCH 54/69] Remove trailing underscores

---
 invisible_cities/reco/wfm_functions.py | 14 +++++++-------
 1 file changed, 7 insertions(+), 7 deletions(-)

diff --git a/invisible_cities/reco/wfm_functions.py b/invisible_cities/reco/wfm_functions.py
index 2f09b9df5e..47c1ce628a 100644
--- a/invisible_cities/reco/wfm_functions.py
+++ b/invisible_cities/reco/wfm_functions.py
@@ -325,27 +325,27 @@ def spatial_selection_method(wfs                 : np.ndarray,
     sipm_y = np.array(detector_info.Y)
 
     if selection_method is SiPMSelectionMethod.median_std_method:
-        starting_ids_ = median_std_method(wfs, **selection_kwargs)
+        starting_ids = median_std_method(wfs, **selection_kwargs)
     elif selection_method is SiPMSelectionMethod.top_n_method:
-        starting_ids_ = top_n_method(wfs, **selection_kwargs)
+        starting_ids = top_n_method(wfs, **selection_kwargs)
     else:
         raise ValueError(f"Selection method {selection_method} not recognized.")
 
-    selected_ids_no_isolated_ = kill_isolated_sipms(
-        starting_ids_,
+    selected_ids_no_isolated = kill_isolated_sipms(
+        starting_ids,
         sipm_x,
         sipm_y,
         proximity_threshold
     )
 
-    sipm_ids_with_signal_ = apply_circular_padding(
-        selected_ids_no_isolated_,
+    sipm_ids_with_signal = apply_circular_padding(
+        selected_ids_no_isolated,
         sipm_x,
         sipm_y,
         padding_radius
     )
 
-    selected_ids = np.where(sipm_ids_with_signal_)[0]
+    selected_ids = np.where(sipm_ids_with_signal)[0]
     selected_wfs = wfs[selected_ids]
 
     return selected_ids, selected_wfs

From 59fb6c125a7cd61929b71def3c35a8e8956eba59 Mon Sep 17 00:00:00 2001
From: Ian0sborne <ian.osborne@student.manchester.ac.uk>
Date: Mon, 8 Jun 2026 12:57:04 +0100
Subject: [PATCH 55/69] Implement `no_cut` option in cut algorithms

---
 invisible_cities/cities/components.py        | 12 ++++++++++++
 invisible_cities/reco/peak_functions.py      |  9 ++-------
 invisible_cities/reco/peak_functions_test.py |  9 +++++++--
 invisible_cities/types/symbols.py            |  1 +
 4 files changed, 22 insertions(+), 9 deletions(-)

diff --git a/invisible_cities/cities/components.py b/invisible_cities/cities/components.py
index 37b6b70dab..7734958458 100644
--- a/invisible_cities/cities/components.py
+++ b/invisible_cities/cities/components.py
@@ -823,6 +823,8 @@ def apply_cutting_function(algo, **cutting_params):
         func = threshold_sipm_selection(**cutting_params)
     elif algo is CutAlgo.pyrrha:
         func = pyrrha_sipm_selection(**cutting_params)
+    elif algo is CutAlgo.no_cut:
+        func = no_cut_sipm_selection()
     else:
         raise ValueError(f"Unsupported cutting algorithm: {algo!r}. Expected one of {list(CutAlgo)}")
 
@@ -876,6 +878,16 @@ def pyrrha_sipm_selection(wfs):
     return pyrrha_sipm_selection
 
 
+def no_cut_sipm_selection():
+    """"
+    Function that applies no cuts to the SiPM waveforms. 
+    """
+    def no_cut_sipm_selection(wfs):
+        sipm_ids = np.arange(wfs.shape[0])
+        return sipm_ids, wfs
+    return no_cut_sipm_selection
+
+
 def calibrate_with_mean(dbfile, run_number):
     DataSiPM   = load_db.DataSiPM(dbfile, run_number)
     adc_to_pes = np.abs(DataSiPM.adc_to_pes.values)
diff --git a/invisible_cities/reco/peak_functions.py b/invisible_cities/reco/peak_functions.py
index e55376c28c..ee41b9a9e2 100644
--- a/invisible_cities/reco/peak_functions.py
+++ b/invisible_cities/reco/peak_functions.py
@@ -81,13 +81,8 @@ def build_sipm_responses(indices, times, widths,
     _, _, sipm_wfs = pick_slice_and_rebin(indices , times, widths,
                                            sipm_wfs, rebin_stride,
                                            pad_zeros = False)
-    if apply_cut is not None:
-        # apply cut before slicing and rebinning
-        (sipm_idx,
-         sipm_wfs)   = apply_cut(sipm_wfs)
-    else:
-        # give all sipm ids as index if no cut is applied
-        sipm_idx     = np.arange(sipm_wfs.shape[0])
+    (sipm_idx,
+     sipm_wfs)   = apply_cut(sipm_wfs)
 
     return SiPMResponses(sipm_ids[sipm_idx], sipm_wfs)
 
diff --git a/invisible_cities/reco/peak_functions_test.py b/invisible_cities/reco/peak_functions_test.py
index a9f2fbe4c4..e14c12a47b 100644
--- a/invisible_cities/reco/peak_functions_test.py
+++ b/invisible_cities/reco/peak_functions_test.py
@@ -408,6 +408,7 @@ def test_build_peak_development(pmt_and_sipm_wfs_with_indices,
                          with_sipms   = with_sipms,
                          Pk           = Pk,
                          sipm_wfs     = sipm_wfs,
+                         apply_cut    = apply_cutting_function(CutAlgo.no_cut)
                          )
 
     assert_Peak_equality(peak, expected_peak)
@@ -485,7 +486,9 @@ def test_find_peaks_s2_style(pmt_and_sipm_wfs_with_indices):
                           time_range, length_range,
                           stride, rebin_stride,
                           S2, pmt_ids, sipm_ids,
-                          sipm_wfs    = sipm_wfs)
+                          sipm_wfs    = sipm_wfs,
+                          apply_cut    = apply_cutting_function(CutAlgo.no_cut)
+                          )
 
     (rebinned_times,
      rebinned_widths,
@@ -515,7 +518,9 @@ def test_get_pmap(s1_and_s2_with_indices):
                        s1_params, s2_params,
                        pmt_ids     = pmt_ids, sipm_ids = sipm_ids,
                        pmt_samp_wid  = pmt_samp_wid ,
-                       sipm_samp_wid = sipm_samp_wid)
+                       sipm_samp_wid = sipm_samp_wid,
+                       apply_cut    = apply_cutting_function(CutAlgo.no_cut)
+                       )
 
     (rebinned_times ,
      rebinned_widths,
diff --git a/invisible_cities/types/symbols.py b/invisible_cities/types/symbols.py
index 62fdfe6a16..4dbd642f73 100644
--- a/invisible_cities/types/symbols.py
+++ b/invisible_cities/types/symbols.py
@@ -158,6 +158,7 @@ class XYReco(AutoNameEnumBase):
 class CutAlgo(AutoNameEnumBase):
     threshold  = auto()
     pyrrha     = auto()
+    no_cut     = auto()
 
 class SiPMSelectionMethod(AutoNameEnumBase):
     median_std_method = auto()

From cc79f2b62c36420f39c0c13b4ea00a88bf938cd3 Mon Sep 17 00:00:00 2001
From: Ian0sborne <ian.osborne@student.manchester.ac.uk>
Date: Fri, 3 Jul 2026 16:18:13 +0100
Subject: [PATCH 56/69] Improve function name

---
 invisible_cities/cities/components.py        |  6 +++---
 invisible_cities/cities/hypathia.py          |  4 ++--
 invisible_cities/cities/irene.py             |  5 +++--
 invisible_cities/reco/peak_functions_test.py | 10 +++++-----
 4 files changed, 13 insertions(+), 12 deletions(-)

diff --git a/invisible_cities/cities/components.py b/invisible_cities/cities/components.py
index 7734958458..d20da203a3 100644
--- a/invisible_cities/cities/components.py
+++ b/invisible_cities/cities/components.py
@@ -817,7 +817,7 @@ def calibrate_sipms(rwf):
     return calibrate_sipms
 
 
-def apply_cutting_function(algo, **cutting_params):
+def select_cutting_algorithm(algo, **cutting_params):
 
     if algo is CutAlgo.threshold:
         func = threshold_sipm_selection(**cutting_params)
@@ -828,10 +828,10 @@ def apply_cutting_function(algo, **cutting_params):
     else:
         raise ValueError(f"Unsupported cutting algorithm: {algo!r}. Expected one of {list(CutAlgo)}")
 
-    def apply_cutting_function(wfs):
+    def select_cutting_algorithm(wfs):
         return func(wfs)
 
-    return apply_cutting_function
+    return select_cutting_algorithm
 
 
 def threshold_sipm_selection(thr_sipm_type
diff --git a/invisible_cities/cities/hypathia.py b/invisible_cities/cities/hypathia.py
index 6593b84c6b..4790ab636f 100644
--- a/invisible_cities/cities/hypathia.py
+++ b/invisible_cities/cities/hypathia.py
@@ -49,7 +49,7 @@
 from .  components import calibrate_sipms
 from .  components import get_actual_sipm_thr
 from .  components import sensor_masker
-from .  components import apply_cutting_function
+from .  components import select_cutting_algorithm
 
 from typing import Dict
 from typing import Any
@@ -124,7 +124,7 @@ def hypathia( files_in         : OneOrManyFiles
                               item = "sipm")
 
     # apply function depending on user input, from provided list of functions
-    apply_cut        = apply_cutting_function(cutting_function, **cutting_params)
+    apply_cut        = select_cutting_algorithm(cutting_function, **cutting_params)
 
 
     event_count_in  = fl.spy_count()
diff --git a/invisible_cities/cities/irene.py b/invisible_cities/cities/irene.py
index f0e04539ae..68822c9d66 100644
--- a/invisible_cities/cities/irene.py
+++ b/invisible_cities/cities/irene.py
@@ -47,8 +47,9 @@
 from .  components import wf_from_files
 from .  components import get_number_of_pmts
 from .  components import compute_and_write_pmaps
-from .  components import apply_cutting_function
 from .  components import sensor_masker
+from .  components import select_cutting_algorithm
+
 
 from typing import Dict
 from typing import Any
@@ -120,7 +121,7 @@ def irene( files_in         : OneOrManyFiles
     sipm_rwf_to_cal  = fl.map(calibrate_sipms(detector_db, run_number),
                               item = "sipm")
     # apply function depending on user input, from provided list of functions
-    apply_cut        = apply_cutting_function(cutting_function, **cutting_params)
+    apply_cut        = select_cutting_algorithm(cutting_function, **cutting_params)
 
 
     event_count_in  = fl.spy_count()
diff --git a/invisible_cities/reco/peak_functions_test.py b/invisible_cities/reco/peak_functions_test.py
index e14c12a47b..e4a48c85fa 100644
--- a/invisible_cities/reco/peak_functions_test.py
+++ b/invisible_cities/reco/peak_functions_test.py
@@ -13,7 +13,7 @@
 from hypothesis.strategies    import integers
 from hypothesis.extra.numpy   import arrays
 
-from ..cities.components      import apply_cutting_function
+from ..cities.components      import select_cutting_algorithm
 from ..core.testing_utils   import exactly
 from ..core.testing_utils   import previous_float
 from ..core.testing_utils   import assert_Peak_equality
@@ -362,7 +362,7 @@ def test_build_sipm_responses(wf_with_indices):
                   thr_sipm      = 0,
                   thr_sipm_type = SiPMThreshold.common,
                   run_number    = 0)
-    apply_cut  = apply_cutting_function(CutAlgo.threshold, **cut_params)
+    apply_cut  = select_cutting_algorithm(CutAlgo.threshold, **cut_params)
 
     sipm_r          = pf.build_sipm_responses(indices, times, widths,
                                               wfs, sipm_ids, 1, apply_cut)
@@ -408,7 +408,7 @@ def test_build_peak_development(pmt_and_sipm_wfs_with_indices,
                          with_sipms   = with_sipms,
                          Pk           = Pk,
                          sipm_wfs     = sipm_wfs,
-                         apply_cut    = apply_cutting_function(CutAlgo.no_cut)
+                         apply_cut    = select_cutting_algorithm(CutAlgo.no_cut)
                          )
 
     assert_Peak_equality(peak, expected_peak)
@@ -487,7 +487,7 @@ def test_find_peaks_s2_style(pmt_and_sipm_wfs_with_indices):
                           stride, rebin_stride,
                           S2, pmt_ids, sipm_ids,
                           sipm_wfs    = sipm_wfs,
-                          apply_cut    = apply_cutting_function(CutAlgo.no_cut)
+                          apply_cut    = select_cutting_algorithm(CutAlgo.no_cut)
                           )
 
     (rebinned_times,
@@ -519,7 +519,7 @@ def test_get_pmap(s1_and_s2_with_indices):
                        pmt_ids     = pmt_ids, sipm_ids = sipm_ids,
                        pmt_samp_wid  = pmt_samp_wid ,
                        sipm_samp_wid = sipm_samp_wid,
-                       apply_cut    = apply_cutting_function(CutAlgo.no_cut)
+                       apply_cut    = select_cutting_algorithm(CutAlgo.no_cut)
                        )
 
     (rebinned_times ,

From a941a426902a6231bcd96dab6be13c1822294dc8 Mon Sep 17 00:00:00 2001
From: Ian0sborne <ian.osborne@student.manchester.ac.uk>
Date: Wed, 8 Jul 2026 16:57:37 +0100
Subject: [PATCH 57/69] Remove unnecessary comments

---
 invisible_cities/cities/components.py  | 2 --
 invisible_cities/reco/wfm_functions.py | 2 --
 2 files changed, 4 deletions(-)

diff --git a/invisible_cities/cities/components.py b/invisible_cities/cities/components.py
index d20da203a3..72417ba25b 100644
--- a/invisible_cities/cities/components.py
+++ b/invisible_cities/cities/components.py
@@ -843,11 +843,9 @@ def threshold_sipm_selection(thr_sipm_type
     Function that applies thresholding to the sipms in standard irene manner,
     by zeroing all waveform values below a threshold.
     '''
-    # assume that if the detector_db is None, you return sipm threshold as the number provided
     if detector_db is None:
         sipm_thr = thr_sipm
     else:
-        # extract sipm threshold
         sipm_thr = get_actual_sipm_thr(thr_sipm_type, thr_sipm, detector_db, run_number)
 
     def threshold_sipm_selection(wfs):
diff --git a/invisible_cities/reco/wfm_functions.py b/invisible_cities/reco/wfm_functions.py
index 47c1ce628a..0bd0ae2c52 100644
--- a/invisible_cities/reco/wfm_functions.py
+++ b/invisible_cities/reco/wfm_functions.py
@@ -194,10 +194,8 @@ def charge_threshold_method(wfs             : np.ndarray,
     -------
     Tuple of np arrays including all passing sipm ids and the corresponding waveforms
     """
-    # zero entries below threshold
     zwfs = zero_wfs_below_threshold(wfs, zeroing_thr)
 
-    # returns selected ids and waveforms above integral
     return select_wfs_above_time_integrated_thr(zwfs, integration_thr)
 
 

From 1cc83678b4106afdfe336968a350c64a2bd35c6b Mon Sep 17 00:00:00 2001
From: Ian0sborne <ian.osborne@student.manchester.ac.uk>
Date: Wed, 8 Jul 2026 17:03:33 +0100
Subject: [PATCH 58/69] Clean up `select_cutting_algorithm()`

---
 invisible_cities/cities/components.py | 12 +++---------
 1 file changed, 3 insertions(+), 9 deletions(-)

diff --git a/invisible_cities/cities/components.py b/invisible_cities/cities/components.py
index 72417ba25b..94ef61c3ca 100644
--- a/invisible_cities/cities/components.py
+++ b/invisible_cities/cities/components.py
@@ -818,21 +818,15 @@ def calibrate_sipms(rwf):
 
 
 def select_cutting_algorithm(algo, **cutting_params):
-
     if algo is CutAlgo.threshold:
-        func = threshold_sipm_selection(**cutting_params)
+        return threshold_sipm_selection(**cutting_params)
     elif algo is CutAlgo.pyrrha:
-        func = pyrrha_sipm_selection(**cutting_params)
+        return pyrrha_sipm_selection(**cutting_params)
     elif algo is CutAlgo.no_cut:
-        func = no_cut_sipm_selection()
+        return no_cut_sipm_selection()
     else:
         raise ValueError(f"Unsupported cutting algorithm: {algo!r}. Expected one of {list(CutAlgo)}")
 
-    def select_cutting_algorithm(wfs):
-        return func(wfs)
-
-    return select_cutting_algorithm
-
 
 def threshold_sipm_selection(thr_sipm_type
                             , thr_sipm

From 37f5449c2eaec59400094ab88a5f3bca1ce68d94 Mon Sep 17 00:00:00 2001
From: Ian0sborne <ian.osborne@student.manchester.ac.uk>
Date: Thu, 9 Jul 2026 14:19:24 +0100
Subject: [PATCH 59/69] Make docstrings self contained

---
 invisible_cities/reco/wfm_functions.py | 39 ++++++++++++--------------
 1 file changed, 18 insertions(+), 21 deletions(-)

diff --git a/invisible_cities/reco/wfm_functions.py b/invisible_cities/reco/wfm_functions.py
index 0bd0ae2c52..3fe073a08c 100644
--- a/invisible_cities/reco/wfm_functions.py
+++ b/invisible_cities/reco/wfm_functions.py
@@ -226,20 +226,17 @@ def kill_isolated_sipms(selected_ids        : np.ndarray,
                         sipm_y              : np.ndarray,
                         proximity_threshold : float) -> np.ndarray:
     """
-    Receives a list of SiPM IDs corresponding to the SiPMs with the most signficant 
-    energy depositions. Scans through these SiPMs to check if they have neighbouring 
-    SiPMs - i.e., within the proximity_threshold - that are also in the initial list
-    of energetic SiPMs. If no neighbours are found, the SiPMs are classed as isolated, 
-    and are removed. Outputs a list of SiPM IDs representing the most energetic SiPMs 
-    belonging to a cluster, which generally occurs in the region where the event took 
-    place.
+    Removes isolated SiPMs from a boolean selection mask of SiPMs. A selected SiPM is considered 
+    isolated if none of the other selected SiPMs lie within the proximity_threshold of it. Isolated 
+    SiPMs are removed from the selection. The output keeps only SiPMs that belong to a cluster 
+    of two or more mutually nearby selected SiPMs.
 
     Parameters
     ----------
-    selected_ids        : Boolean array of shape (n_sipms,) corresponding to the most energetic SiPMs.
-    sipm_x              : 1D array of shape (n_sipms,) containing the x positions of the SiPMs.
-    sipm_y              : 1D array of shape (n_sipms,) containing the y positions of the SiPMs.
-    proximity_threshold : Distance threshold in mm used to identify isolated SiPMs.
+    selected_ids             : Boolean array of shape (n_sipms,) corresponding to the initial selection of SiPMs.
+    sipm_x                   : 1D array of shape (n_sipms,) containing the x positions of the SiPMs.
+    sipm_y                   : 1D array of shape (n_sipms,) containing the y positions of the SiPMs.
+    proximity_threshold      : Distance threshold in mm used to identify isolated SiPMs.
 
     Returns
     -------
@@ -265,21 +262,22 @@ def apply_circular_padding(selected_ids_no_isolated : np.ndarray,
                            sipm_y                   : np.ndarray,
                            padding_radius           : float) -> np.ndarray:
     """
-    Receives a list of SiPM IDs corresponding to the most energetic SiPMs clustered 
-    near the event. For these SiPMs, creates circular padding of radius padding_radius, 
-    selecting all SiPMs within that radius. Stores the union of all selected SiPMs. 
-    Outputs the SiPM IDs which are relevant for a given event.
+    Expands a boolean selection mask of SiPMs to include all SiPMs within a given region determined
+    by the padding_radius. For each selected SiPM, all SiPM within the padding_radius are added to
+    the selection; the result is the union of all such neighborhoods together with the original selection.
 
     Parameters
     ----------
-    selected_ids_no_isolated : Boolean array of shape (n_sipms,) corresponding to the "relevant" SiPMs.
+    selected_ids_no_isolated : Boolean array of shape (n_sipms,). True marks the SiPMs to pad around.
     sipm_x                   : 1D array of shape (n_sipms,) containing the x positions of the SiPMs.
     sipm_y                   : 1D array of shape (n_sipms,) containing the y positions of the SiPMs.
-    padding_radius           : Distance threshold in mm used to create circular padding around selected SiPMs.
+    padding_radius           : Distance threshold in mm used to include neighboring SiPMs around each 
+                               selected SiPM.
 
     Returns
     -------
-    sipm_ids_with_signal : Boolean array of shape (n_sipms,) where True indicates that the SiPM is selected.
+    sipm_ids_with_signal     : Boolean array of shape (n_sipms,). True for SiPMs that are either in the
+                               original selection or within the "padding_radius" of a selected SiPM.
     """
     sipm_ids_with_signal = np.zeros_like(selected_ids_no_isolated, dtype=bool)
 
@@ -299,9 +297,8 @@ def spatial_selection_method(wfs                 : np.ndarray,
                              run_number          : int,
                              detector_db         : str) -> np.ndarray:
     """
-    SiPM selection function, applies SiPM cuts based on user input.
-    A first selection of SiPMs is made, isolated SiPMs are removed
-    and padding is added around the SiPMs that are left.
+    SiPM selection function, applies SiPM cuts based on user input. A first selection of SiPMs is made, 
+    isolated SiPMs are removed and padding is added around the SiPMs that are left.
 
     Parameters
     ----------

From 3cdc97562ba9903ed990040c299c466dd2588193 Mon Sep 17 00:00:00 2001
From: Ian0sborne <ian.osborne@student.manchester.ac.uk>
Date: Thu, 9 Jul 2026 15:21:07 +0100
Subject: [PATCH 60/69] Implement `CutAlgo.no_cut` as default parameter

---
 invisible_cities/reco/peak_functions.py | 21 +++++++++++----------
 1 file changed, 11 insertions(+), 10 deletions(-)

diff --git a/invisible_cities/reco/peak_functions.py b/invisible_cities/reco/peak_functions.py
index ee41b9a9e2..2b176d9a6e 100644
--- a/invisible_cities/reco/peak_functions.py
+++ b/invisible_cities/reco/peak_functions.py
@@ -10,13 +10,14 @@
 
 import numpy        as np
 
-from .. core               import system_of_units as units
-from .. evm .ic_containers import ZsWf
-from .. evm .pmaps         import S1
-from .. evm .pmaps         import S2
-from .. evm .pmaps         import PMap
-from .. evm .pmaps         import PMTResponses
-from .. evm .pmaps         import SiPMResponses
+from .. core                 import system_of_units as units
+from .. evm   .ic_containers import ZsWf
+from .. evm   .pmaps         import S1
+from .. evm   .pmaps         import S2
+from .. evm   .pmaps         import PMap
+from .. evm   .pmaps         import PMTResponses
+from .. evm   .pmaps         import SiPMResponses
+from .. types .symbols       import CutAlgo
 
 
 def indices_and_wf_above_threshold(wf, thr):
@@ -94,7 +95,7 @@ def build_peak(indices, times,
                pmt_samp_wid  = 25 * units.ns,
                sipm_samp_wid =  1 * units.mus,
                sipm_wfs      = None,
-               apply_cut     = None):
+               apply_cut     = CutAlgo.no_cut):
     sipm_pmt_bin_ratio = int(sipm_samp_wid/pmt_samp_wid)
     (pk_times ,
      pk_widths,
@@ -121,7 +122,7 @@ def find_peaks(ccwfs, index,
                Pk, pmt_ids, sipm_ids=None,
                pmt_samp_wid = 25*units.ns,
                sipm_samp_wid = 1*units.mus,
-               sipm_wfs=None, apply_cut = None):
+               sipm_wfs=None, apply_cut = CutAlgo.no_cut):
 
     ccwfs = np.array(ccwfs, ndmin=2)
 
@@ -145,7 +146,7 @@ def find_peaks(ccwfs, index,
 
 def get_pmap(ccwf, s1_indx, s2_indx, sipm_zs_wf,
              s1_params, s2_params, pmt_ids, sipm_ids,
-             pmt_samp_wid, sipm_samp_wid, apply_cut = None):
+             pmt_samp_wid, sipm_samp_wid, apply_cut = CutAlgo.no_cut):
     return PMap(find_peaks(ccwf, s1_indx, Pk=S1, pmt_ids=pmt_ids,
                            pmt_samp_wid=pmt_samp_wid,
                            **s1_params),

From 476c5d7f77357872eba55b45fd01dbff8ce8bd09 Mon Sep 17 00:00:00 2001
From: Ian0sborne <ian.osborne@student.manchester.ac.uk>
Date: Thu, 9 Jul 2026 16:12:50 +0100
Subject: [PATCH 61/69] Remove  unused `thr_sipm*` params

---
 invisible_cities/cities/components.py |  2 +-
 invisible_cities/cities/hypathia.py   |  6 ++--
 invisible_cities/cities/irene.py      |  5 +---
 invisible_cities/cities/irene_test.py | 41 +++++++++++++++------------
 invisible_cities/config/hypathia.conf | 14 ++-------
 invisible_cities/config/irene.conf    | 13 ++-------
 6 files changed, 33 insertions(+), 48 deletions(-)

diff --git a/invisible_cities/cities/components.py b/invisible_cities/cities/components.py
index 94ef61c3ca..fbc01dc9c2 100644
--- a/invisible_cities/cities/components.py
+++ b/invisible_cities/cities/components.py
@@ -1295,7 +1295,7 @@ def integrate_wfs(wfs):
 # Compound components
 def compute_and_write_pmaps(detector_db, run_number, pmt_samp_wid, sipm_samp_wid,
                   s1_lmax, s1_lmin, s1_rebin_stride, s1_stride, s1_tmax, s1_tmin,
-                  s2_lmax, s2_lmin, s2_rebin_stride, s2_stride, s2_tmax, s2_tmin, thr_sipm_s2,
+                  s2_lmax, s2_lmin, s2_rebin_stride, s2_stride, s2_tmax, s2_tmin,
                   h5out, apply_cut, sipm_rwf_to_cal=None):
 
     # Filter events without signal over threshold
diff --git a/invisible_cities/cities/hypathia.py b/invisible_cities/cities/hypathia.py
index 4790ab636f..ef7ba0abbd 100644
--- a/invisible_cities/cities/hypathia.py
+++ b/invisible_cities/cities/hypathia.py
@@ -65,8 +65,6 @@ def hypathia( files_in         : OneOrManyFiles
             , run_number       : int
             , sipm_noise_cut   : float
             , filter_padding   : int
-            , thr_sipm         : float
-            , thr_sipm_type    : SiPMThreshold
             , pmt_wfs_rebin    : int
             , pmt_pe_rms       : float
             , s1_lmin          : int  , s1_lmax     : int
@@ -76,7 +74,7 @@ def hypathia( files_in         : OneOrManyFiles
             , s2_lmin          : int  , s2_lmax     : int
             , s2_tmin          : float, s2_tmax     : float
             , s2_rebin_stride  : int  , s2_stride   : int
-            , thr_csum_s2      : float, thr_sipm_s2 : float
+            , thr_csum_s2      : float
             , pmt_samp_wid     : float
             , sipm_samp_wid    : float
             , cutting_function : CutAlgo
@@ -145,7 +143,7 @@ def hypathia( files_in         : OneOrManyFiles
         compute_pmaps, empty_indices, empty_pmaps = compute_and_write_pmaps(
                                              detector_db, run_number, pmt_samp_wid, sipm_samp_wid,
                                              s1_lmax, s1_lmin, s1_rebin_stride, s1_stride, s1_tmax, s1_tmin,
-                                             s2_lmax, s2_lmin, s2_rebin_stride, s2_stride, s2_tmax, s2_tmin, thr_sipm_s2,
+                                             s2_lmax, s2_lmin, s2_rebin_stride, s2_stride, s2_tmax, s2_tmin,
                                              h5out, apply_cut, sipm_rwf_to_cal)
 
         result = push(source = wf_from_files(files_in, WfType.mcrd),
diff --git a/invisible_cities/cities/irene.py b/invisible_cities/cities/irene.py
index 68822c9d66..da6f091e4a 100644
--- a/invisible_cities/cities/irene.py
+++ b/invisible_cities/cities/irene.py
@@ -65,8 +65,6 @@ def irene( files_in         : OneOrManyFiles
          , n_baseline       : int
          , n_maw            : int
          , thr_maw          : float
-         , thr_sipm         : float
-         , thr_sipm_type    : SiPMThreshold
          , s1_lmin          : int  , s1_lmax      : int
          , s1_tmin          : float, s1_tmax      : float
          , s1_rebin_stride  : int  , s1_stride    : int
@@ -74,7 +72,7 @@ def irene( files_in         : OneOrManyFiles
          , s2_lmin          : int  , s2_lmax      : int
          , s2_tmin          : float, s2_tmax      : float
          , s2_rebin_stride  : int  , s2_stride    : int
-         , thr_csum_s2      : float, thr_sipm_s2  : float
+         , thr_csum_s2      : float
          , pmt_samp_wid     : float, sipm_samp_wid: float
          , store_db         : bool = True
          , cutting_function : CutAlgo
@@ -145,7 +143,6 @@ def irene( files_in         : OneOrManyFiles
                                          detector_db, run_number, pmt_samp_wid, sipm_samp_wid,
                                          s1_lmax, s1_lmin, s1_rebin_stride, s1_stride, s1_tmax, s1_tmin,
                                          s2_lmax, s2_lmin, s2_rebin_stride, s2_stride, s2_tmax, s2_tmin,
-                                         thr_sipm_s2,
                                          h5out, apply_cut, sipm_rwf_to_cal)
 
         result = push(source = wf_from_files(files_in, WfType.rwf),
diff --git a/invisible_cities/cities/irene_test.py b/invisible_cities/cities/irene_test.py
index d0fe3f4270..a3f7a41513 100644
--- a/invisible_cities/cities/irene_test.py
+++ b/invisible_cities/cities/irene_test.py
@@ -82,13 +82,15 @@ def test_irene_electrons_40keV(config_tmpdir, ICDATADIR, s12params,
     nrequired  = 2
 
     conf = configure('dummy invisible_cities/config/irene.conf'.split())
-    conf.update(dict(detector_db = DetDB.new,
-                     run_number    = 0,
+    cutting_params = conf['cutting_params'].copy()
+    cutting_params.update(thr_sipm_type = thr_sipm_type,
+                          thr_sipm      = thr_sipm_value)
+    conf.update(dict(detector_db    = DetDB.new,
+                     run_number     = 0,
                      files_in      = PATH_IN,
                      file_out      = PATH_OUT,
                      event_range   = (0, nrequired),
-                     thr_sipm_type = thr_sipm_type,
-                     thr_sipm      = thr_sipm_value,
+                     cutting_params = cutting_params,
                      **unpack_s12params(s12params)))
 
     cnt = irene(**conf)
@@ -414,20 +416,23 @@ def test_irene_sequential_times(config_tmpdir, ICDATADIR):
     PATH_OUT = os.path.join(config_tmpdir, 'test_pmaps.h5')
 
     conf = configure('dummy invisible_cities/config/irene.conf'.split())
-    conf.update(dict(files_in    = PATH_IN           ,
-                     file_out    = PATH_OUT          ,
-                     run_number  =   6351            ,
-                     n_baseline  =  48000            ,
-                     thr_sipm    =      1 * units.pes,
-                     s1_tmin     =      0 * units.mus,
-                     s1_tmax     =    640 * units.mus,
-                     s1_lmin     =      5            ,
-                     s1_lmax     =     30            ,
-                     s2_tmin     =    645 * units.mus,
-                     s2_tmax     =   1300 * units.mus,
-                     s2_lmin     =     80            ,
-                     s2_lmax     = 200000            ,
-                     thr_sipm_s2 =      5 * units.pes))
+    cutting_params = conf['cutting_params'].copy()
+    cutting_params.update(run_number    = 6351,
+                          thr_sipm_s2   = 5 * units.pes,
+                          thr_sipm      = 1 * units.pes)
+    conf.update(dict(files_in       = PATH_IN           ,
+                     file_out       = PATH_OUT          ,
+                     run_number     =   6351            ,
+                     n_baseline     =  48000            ,
+                     s1_tmin        =      0 * units.mus,
+                     s1_tmax        =    640 * units.mus,
+                     s1_lmin        =      5            ,
+                     s1_lmax        =     30            ,
+                     s2_tmin        =    645 * units.mus,
+                     s2_tmax        =   1300 * units.mus,
+                     s2_lmin        =     80            ,
+                     s2_lmax        = 200000            ,
+                     cutting_params = cutting_params))
 
     irene(**conf)
 
diff --git a/invisible_cities/config/hypathia.conf b/invisible_cities/config/hypathia.conf
index d3029f25f9..4d1838bb01 100644
--- a/invisible_cities/config/hypathia.conf
+++ b/invisible_cities/config/hypathia.conf
@@ -25,11 +25,6 @@ event_range =  0, 2
 thr_csum_s1 = 0.5 * pes
 thr_csum_s2 = 2.0 * pes
 
-# Set thresholds for SiPM
-thr_sipm      = 1.0 * pes
-thr_sipm_type = common
-
-
 # Set parameters to search for S1
 # Notice that in MC file S1 is in t=100 mus
 s1_tmin       =  99 * mus # position of S1 in MC files at 100 mus
@@ -47,15 +42,12 @@ s2_lmin     =     80       # 100 x 25 = 2.5 mus
 s2_lmax     = 100000       # maximum value of S2 width
 s2_rebin_stride = 40       # Rebin by default, 40 25 ns time bins to make one 1us time bin
 
-# Set S2Si parameters
-thr_sipm_s2 = 5 * pes  # Threshold for the full sipm waveform
-
 pmt_samp_wid  = 25 * ns
 sipm_samp_wid = 1 * mus
 
 cutting_function = threshold 
-cutting_params   = dict(  thr_sipm_type = thr_sipm_type 
-                        , thr_sipm      = thr_sipm 
-                        , thr_sipm_s2   = thr_sipm_s2
+cutting_params   = dict(  thr_sipm_type = common
+                        , thr_sipm      = 1.0 * pes # Threshold for each SiPM time bin
+                        , thr_sipm_s2   = 5 * pes   # Threshold for the full sipm waveform
                         , detector_db   = detector_db 
                         , run_number    = run_number)
diff --git a/invisible_cities/config/irene.conf b/invisible_cities/config/irene.conf
index 1f2312a78d..bc5bdff334 100644
--- a/invisible_cities/config/irene.conf
+++ b/invisible_cities/config/irene.conf
@@ -26,10 +26,6 @@ thr_maw =   3 * adc
 thr_csum_s1 = 0.5 * pes
 thr_csum_s2 = 1.0 * pes
 
-# Set thresholds for SiPM
-thr_sipm      = 3.5 * pes
-thr_sipm_type = common
-
 # Set parameters to search for S1
 # Notice that in MC file S1 is in t=100 mus
 s1_tmin       =  99 * mus # position of S1 in MC files at 100 mus
@@ -47,15 +43,12 @@ s2_lmin     =    100       # 100 x 25 = 2.5 mus
 s2_lmax     = 100000       # maximum value of S2 width
 s2_rebin_stride = 40       # Rebin by default, 40 25 ns time bins to make one 1us time bin
 
-# Set S2Si parameters
-thr_sipm_s2 = 10 * pes  # Threshold for the full sipm waveform
-
 pmt_samp_wid  = 25 * ns
 sipm_samp_wid = 1 * mus
 
 cutting_function = threshold 
-cutting_params   = dict(  thr_sipm_type = thr_sipm_type 
-                        , thr_sipm      = thr_sipm 
-                        , thr_sipm_s2   = thr_sipm_s2
+cutting_params   = dict(  thr_sipm_type = common 
+                        , thr_sipm      = 3.5 * pes # Threshold for each SiPM time bin
+                        , thr_sipm_s2   = 10 * pes  # Threshold for the full sipm waveform
                         , detector_db   = detector_db 
                         , run_number    = run_number)

From 9d752da98398d3ad16bb2b761ee51334a88e7a0a Mon Sep 17 00:00:00 2001
From: Ian0sborne <ian.osborne@student.manchester.ac.uk>
Date: Mon, 13 Jul 2026 13:34:11 +0100
Subject: [PATCH 62/69] Guard against zero-suppressed data

---
 .../calib/calib_sensors_functions.py           | 18 +++++++++++-------
 1 file changed, 11 insertions(+), 7 deletions(-)

diff --git a/invisible_cities/calib/calib_sensors_functions.py b/invisible_cities/calib/calib_sensors_functions.py
index f98ea569f5..e58d2b7f23 100644
--- a/invisible_cities/calib/calib_sensors_functions.py
+++ b/invisible_cities/calib/calib_sensors_functions.py
@@ -59,8 +59,9 @@ def modes  (wfs): return to_col_vector(mode  (wfs, axis=1))
 
 def subtract_baseline(wfs, *, bls_mode=BlsMode.mean):
     """
-    Subtract the baseline to all waveforms in the input
-    with a specific algorithm.
+    Subtract the baseline from all waveforms in the input
+    with a specific algorithm, leaving individual samples
+    that are already 0 untouched.
 
     Parameters
     ----------
@@ -75,16 +76,19 @@ def subtract_baseline(wfs, *, bls_mode=BlsMode.mean):
     Returns
     -------
     bls: np.ndarray with shape (n, m)
-        Baseline-subtracted waveforms.
+        Baseline-subtracted waveforms, with originally-zero
+        samples left as 0.
     """
 
-    if   bls_mode is BlsMode.mean     : return wfs - means     (wfs)
-    elif bls_mode is BlsMode.median   : return wfs - medians   (wfs)
-    elif bls_mode is BlsMode.mode     : return wfs - modes     (wfs)
-    elif bls_mode is BlsMode.scipymode: return wfs - scipy_mode(wfs, axis=1)
+    if   bls_mode is BlsMode.mean     : baseline = means     (wfs)
+    elif bls_mode is BlsMode.median   : baseline = medians   (wfs)
+    elif bls_mode is BlsMode.mode     : baseline = modes     (wfs)
+    elif bls_mode is BlsMode.scipymode: baseline = scipy_mode(wfs, axis=1)
     else:
         raise TypeError(f"Unrecognized baseline subtraction option: {bls_mode}")
 
+    return np.where(wfs == 0, 0, wfs - baseline) # subtracts the baseline to non-zero samples only
+
 
 def calibrate_wfs(wfs, adc_to_pes):
     """

From 7153a25a290231fc4874bf6f8249e6c486d66d8c Mon Sep 17 00:00:00 2001
From: Ian0sborne <ian.osborne@student.manchester.ac.uk>
Date: Tue, 14 Jul 2026 16:13:44 +0100
Subject: [PATCH 63/69] Improve `zero_wfs_below_threshold()`

---
 invisible_cities/reco/wfm_functions.py | 6 +++++-
 1 file changed, 5 insertions(+), 1 deletion(-)

diff --git a/invisible_cities/reco/wfm_functions.py b/invisible_cities/reco/wfm_functions.py
index 3fe073a08c..3cfa92dc7e 100644
--- a/invisible_cities/reco/wfm_functions.py
+++ b/invisible_cities/reco/wfm_functions.py
@@ -152,7 +152,11 @@ def zero_wfs_below_threshold(wfs         : np.ndarray,
     -------
     2D array of shape (n_sipms, n_time_bins) containing the input waveforms with entries below threshold set to zero.
     """
-    thr = to_col_vector(np.full(wfs.shape[0], zeroing_thr))
+    if isinstance(zeroing_thr, (int, float)):
+        thr = np.full((wfs.shape[0], 1), zeroing_thr)
+    else:
+        thr = np.reshape(zeroing_thr, (wfs.shape[0], 1))
+
     return np.where(wfs > thr, wfs, 0)
 
 

From 4c01f8e8698a4572dc99b7c0dd8b891450dafb42 Mon Sep 17 00:00:00 2001
From: Ian0sborne <ian.osborne@student.manchester.ac.uk>
Date: Thu, 16 Jul 2026 14:34:51 +0100
Subject: [PATCH 64/69] Add Irene config with Pyrrha settings

---
 invisible_cities/config/irene_pyrrha.conf | 55 +++++++++++++++++++++++
 1 file changed, 55 insertions(+)
 create mode 100644 invisible_cities/config/irene_pyrrha.conf

diff --git a/invisible_cities/config/irene_pyrrha.conf b/invisible_cities/config/irene_pyrrha.conf
new file mode 100644
index 0000000000..8fe2a9d230
--- /dev/null
+++ b/invisible_cities/config/irene_pyrrha.conf
@@ -0,0 +1,55 @@
+files_in = '$ICDIR/database/test_data/electrons_40keV_z25_RWF.h5'
+
+# REPLACE /tmp with your output directory
+file_out = '/tmp/electrons_40keV_z25_PMP.h5'
+
+# compression library
+compression = 'ZLIB4'
+
+# run number 0 is for MC
+run_number = 0
+detector_db = 'new'
+
+# How frequently to print events
+print_mod = 1
+
+# max number of events to run
+event_range =  1
+
+n_baseline =   28000 # for a window of 800 mus
+
+# Set MAW for calibrated sum
+n_maw   = 100
+thr_maw =   3 * adc
+
+# Set thresholds for calibrated sum
+thr_csum_s1 = 0.5 * pes
+thr_csum_s2 = 1.0 * pes
+
+# Set parameters to search for S1
+# Notice that in MC file S1 is in t=100 mus
+s1_tmin       =  99 * mus # position of S1 in MC files at 100 mus
+s1_tmax       = 101 * mus # change tmin and tmax if S1 not at 100 mus
+s1_stride     =   4       # minimum number of 25 ns bins in S1 searches
+s1_lmin       =   8       # 8 x 25 = 200 ns
+s1_lmax       =  20       # 20 x 25 = 500 ns
+s1_rebin_stride = 1       # Do not rebin S1 by default
+
+# Set parameters to search for S2
+s2_tmin     =    101 * mus # assumes S1 at 100 mus, change if S1 not at 100 mus
+s2_tmax     =   1199 * mus # end of the window
+s2_stride   =     40       #  40 x 25 = 1   mus
+s2_lmin     =    100       # 100 x 25 = 2.5 mus
+s2_lmax     = 100000       # maximum value of S2 width
+s2_rebin_stride = 40       # Rebin by default, 40 25 ns time bins to make one 1us time bin
+
+pmt_samp_wid  = 25 * ns
+sipm_samp_wid = 1 * mus
+
+cutting_function = pyrrha
+cutting_params   = dict(  selection_method    = median_std_method
+                        , selection_kwargs    = {'nsigma': 3} # method for selecting energetic SiPMs
+                        , proximity_threshold = 25  # energetic SiPMs with no neighbors within this distance (in mm) are discarded
+                        , padding_radius      = 50 # amount of padding around selected SiPMs (in mm)
+                        , run_number          = run_number
+                        , detector_db         = detector_db)

From 249da176408efe27e79f51f92427a653942fb1e9 Mon Sep 17 00:00:00 2001
From: Ian0sborne <ian.osborne@student.manchester.ac.uk>
Date: Thu, 16 Jul 2026 15:07:29 +0100
Subject: [PATCH 65/69] Change `apply_cut` to `sipm_selection_algo`

---
 invisible_cities/cities/components.py        | 10 +++--
 invisible_cities/cities/hypathia.py          |  4 +-
 invisible_cities/cities/irene.py             |  4 +-
 invisible_cities/reco/peak_functions.py      | 19 ++++-----
 invisible_cities/reco/peak_functions_test.py | 42 ++++++++++----------
 5 files changed, 41 insertions(+), 38 deletions(-)

diff --git a/invisible_cities/cities/components.py b/invisible_cities/cities/components.py
index fbc01dc9c2..e160b75593 100644
--- a/invisible_cities/cities/components.py
+++ b/invisible_cities/cities/components.py
@@ -764,7 +764,8 @@ def sensor_data(path, wf_type):
 
 def build_pmap(detector_db, run_number, pmt_samp_wid, sipm_samp_wid,
                s1_lmax, s1_lmin, s1_rebin_stride, s1_stride, s1_tmax, s1_tmin,
-               s2_lmax, s2_lmin, s2_rebin_stride, s2_stride, s2_tmax, s2_tmin, apply_cut):
+               s2_lmax, s2_lmin, s2_rebin_stride, s2_stride, s2_tmax, s2_tmin, 
+               sipm_selection_algo):
     s1_params = dict(time        = minmax(min = s1_tmin,
                                           max = s1_tmax),
                     length       = minmax(min = s1_lmin,
@@ -786,7 +787,7 @@ def build_pmap(detector_db, run_number, pmt_samp_wid, sipm_samp_wid,
     def build_pmap(ccwf, s1_indx, s2_indx, sipmzs): # -> PMap
         return pkf.get_pmap(ccwf, s1_indx, s2_indx, sipmzs,
                             s1_params, s2_params, pmt_ids, sipm_ids,
-                            pmt_samp_wid, sipm_samp_wid, apply_cut)
+                            pmt_samp_wid, sipm_samp_wid, sipm_selection_algo)
 
     return build_pmap
 
@@ -1296,7 +1297,7 @@ def integrate_wfs(wfs):
 def compute_and_write_pmaps(detector_db, run_number, pmt_samp_wid, sipm_samp_wid,
                   s1_lmax, s1_lmin, s1_rebin_stride, s1_stride, s1_tmax, s1_tmin,
                   s2_lmax, s2_lmin, s2_rebin_stride, s2_stride, s2_tmax, s2_tmin,
-                  h5out, apply_cut, sipm_rwf_to_cal=None):
+                  h5out, sipm_selection_algo, sipm_rwf_to_cal=None):
 
     # Filter events without signal over threshold
     indices_pass    = fl.map(check_nonempty_indices,
@@ -1307,7 +1308,8 @@ def compute_and_write_pmaps(detector_db, run_number, pmt_samp_wid, sipm_samp_wid
     # Build the PMap
     compute_pmap     = fl.map(build_pmap(detector_db, run_number, pmt_samp_wid, sipm_samp_wid,
                                          s1_lmax, s1_lmin, s1_rebin_stride, s1_stride, s1_tmax, s1_tmin,
-                                         s2_lmax, s2_lmin, s2_rebin_stride, s2_stride, s2_tmax, s2_tmin, apply_cut),
+                                         s2_lmax, s2_lmin, s2_rebin_stride, s2_stride, s2_tmax, s2_tmin, 
+                                         sipm_selection_algo),
                               args = ("ccwfs", "s1_indices", "s2_indices", "sipm"),
                               out  = "pmap")
 
diff --git a/invisible_cities/cities/hypathia.py b/invisible_cities/cities/hypathia.py
index ef7ba0abbd..d99458039e 100644
--- a/invisible_cities/cities/hypathia.py
+++ b/invisible_cities/cities/hypathia.py
@@ -122,7 +122,7 @@ def hypathia( files_in         : OneOrManyFiles
                               item = "sipm")
 
     # apply function depending on user input, from provided list of functions
-    apply_cut        = select_cutting_algorithm(cutting_function, **cutting_params)
+    sipm_selection_algo = select_cutting_algorithm(cutting_function, **cutting_params)
 
 
     event_count_in  = fl.spy_count()
@@ -144,7 +144,7 @@ def hypathia( files_in         : OneOrManyFiles
                                              detector_db, run_number, pmt_samp_wid, sipm_samp_wid,
                                              s1_lmax, s1_lmin, s1_rebin_stride, s1_stride, s1_tmax, s1_tmin,
                                              s2_lmax, s2_lmin, s2_rebin_stride, s2_stride, s2_tmax, s2_tmin,
-                                             h5out, apply_cut, sipm_rwf_to_cal)
+                                             h5out, sipm_selection_algo, sipm_rwf_to_cal)
 
         result = push(source = wf_from_files(files_in, WfType.mcrd),
                       pipe   = pipe(fl.slice(*event_range, close_all=True),
diff --git a/invisible_cities/cities/irene.py b/invisible_cities/cities/irene.py
index da6f091e4a..8856538030 100644
--- a/invisible_cities/cities/irene.py
+++ b/invisible_cities/cities/irene.py
@@ -119,7 +119,7 @@ def irene( files_in         : OneOrManyFiles
     sipm_rwf_to_cal  = fl.map(calibrate_sipms(detector_db, run_number),
                               item = "sipm")
     # apply function depending on user input, from provided list of functions
-    apply_cut        = select_cutting_algorithm(cutting_function, **cutting_params)
+    sipm_selection_algo = select_cutting_algorithm(cutting_function, **cutting_params)
 
 
     event_count_in  = fl.spy_count()
@@ -143,7 +143,7 @@ def irene( files_in         : OneOrManyFiles
                                          detector_db, run_number, pmt_samp_wid, sipm_samp_wid,
                                          s1_lmax, s1_lmin, s1_rebin_stride, s1_stride, s1_tmax, s1_tmin,
                                          s2_lmax, s2_lmin, s2_rebin_stride, s2_stride, s2_tmax, s2_tmin,
-                                         h5out, apply_cut, sipm_rwf_to_cal)
+                                         h5out, sipm_selection_algo, sipm_rwf_to_cal)
 
         result = push(source = wf_from_files(files_in, WfType.rwf),
                       pipe   = pipe(fl.slice(*event_range, close_all=True),
diff --git a/invisible_cities/reco/peak_functions.py b/invisible_cities/reco/peak_functions.py
index 2b176d9a6e..3f796176fc 100644
--- a/invisible_cities/reco/peak_functions.py
+++ b/invisible_cities/reco/peak_functions.py
@@ -77,13 +77,13 @@ def build_pmt_responses(indices, times, widths, ccwf,
     return pk_times, pk_widths, PMTResponses(pmt_ids, pmt_wfs)
 
 
-def build_sipm_responses(indices, times, widths,
-                         sipm_wfs, sipm_ids, rebin_stride, apply_cut):
+def build_sipm_responses(indices, times, widths, sipm_wfs, 
+                         sipm_ids, rebin_stride, sipm_selection_algo):
     _, _, sipm_wfs = pick_slice_and_rebin(indices , times, widths,
                                            sipm_wfs, rebin_stride,
                                            pad_zeros = False)
     (sipm_idx,
-     sipm_wfs)   = apply_cut(sipm_wfs)
+     sipm_wfs)   = sipm_selection_algo(sipm_wfs)
 
     return SiPMResponses(sipm_ids[sipm_idx], sipm_wfs)
 
@@ -95,7 +95,7 @@ def build_peak(indices, times,
                pmt_samp_wid  = 25 * units.ns,
                sipm_samp_wid =  1 * units.mus,
                sipm_wfs      = None,
-               apply_cut     = CutAlgo.no_cut):
+               sipm_selection_algo = CutAlgo.no_cut):
     sipm_pmt_bin_ratio = int(sipm_samp_wid/pmt_samp_wid)
     (pk_times ,
      pk_widths,
@@ -109,7 +109,7 @@ def build_peak(indices, times,
                                       widths * sipm_pmt_bin_ratio,
                                       sipm_wfs, sipm_ids,
                                       rebin_stride // sipm_pmt_bin_ratio,
-                                      apply_cut)
+                                      sipm_selection_algo)
     else:
         sipm_r = SiPMResponses.build_empty_instance()
 
@@ -122,7 +122,8 @@ def find_peaks(ccwfs, index,
                Pk, pmt_ids, sipm_ids=None,
                pmt_samp_wid = 25*units.ns,
                sipm_samp_wid = 1*units.mus,
-               sipm_wfs=None, apply_cut = CutAlgo.no_cut):
+               sipm_wfs=None, 
+               sipm_selection_algo = CutAlgo.no_cut):
 
     ccwfs = np.array(ccwfs, ndmin=2)
 
@@ -139,20 +140,20 @@ def find_peaks(ccwfs, index,
                         rebin_stride,
                         with_sipms, Pk,
                         pmt_samp_wid, sipm_samp_wid,
-                        sipm_wfs, apply_cut)
+                        sipm_wfs, sipm_selection_algo)
         peaks.append(pk)
     return peaks
 
 
 def get_pmap(ccwf, s1_indx, s2_indx, sipm_zs_wf,
              s1_params, s2_params, pmt_ids, sipm_ids,
-             pmt_samp_wid, sipm_samp_wid, apply_cut = CutAlgo.no_cut):
+             pmt_samp_wid, sipm_samp_wid, sipm_selection_algo = CutAlgo.no_cut):
     return PMap(find_peaks(ccwf, s1_indx, Pk=S1, pmt_ids=pmt_ids,
                            pmt_samp_wid=pmt_samp_wid,
                            **s1_params),
                 find_peaks(ccwf, s2_indx, Pk=S2, pmt_ids=pmt_ids, sipm_ids=sipm_ids,
                            sipm_wfs      = sipm_zs_wf,
-                           apply_cut     = apply_cut,
+                           sipm_selection_algo = sipm_selection_algo,
                            pmt_samp_wid  = pmt_samp_wid,
                            sipm_samp_wid = sipm_samp_wid,
                            **s2_params))
diff --git a/invisible_cities/reco/peak_functions_test.py b/invisible_cities/reco/peak_functions_test.py
index e4a48c85fa..4f54652e02 100644
--- a/invisible_cities/reco/peak_functions_test.py
+++ b/invisible_cities/reco/peak_functions_test.py
@@ -350,22 +350,22 @@ def test_build_pmt_responses(wf_with_indices):
 def test_build_sipm_responses(wf_with_indices):
     times, widths, wfs, indices = wf_with_indices
     ids = np.arange(wfs.shape[0])
-    wfs_slice       = wfs[:, indices]
-    peak_integrals  = wfs_slice.sum(axis=1)
-    below_thr_index = np.argmin (peak_integrals)
+    wfs_slice           = wfs[:, indices]
+    peak_integrals      = wfs_slice.sum(axis=1)
+    below_thr_index     = np.argmin (peak_integrals)
     # next_float doesn't work here
-    thr             = peak_integrals[below_thr_index] * 1.000001
-    sipm_ids        = np.arange(len(wfs))
+    thr                 = peak_integrals[below_thr_index] * 1.000001
+    sipm_ids            = np.arange(len(wfs))
 
-    cut_params = dict(detector_db   = 'None',
-                  thr_sipm_s2   = thr,
-                  thr_sipm      = 0,
-                  thr_sipm_type = SiPMThreshold.common,
-                  run_number    = 0)
-    apply_cut  = select_cutting_algorithm(CutAlgo.threshold, **cut_params)
+    cut_params          = dict(detector_db   = 'None',
+                               thr_sipm_s2   = thr,
+                               thr_sipm      = 0,
+                               thr_sipm_type = SiPMThreshold.common,
+                               run_number    = 0)
+    sipm_selection_algo  = select_cutting_algorithm(CutAlgo.threshold, **cut_params)
 
-    sipm_r          = pf.build_sipm_responses(indices, times, widths,
-                                              wfs, sipm_ids, 1, apply_cut)
+    sipm_r               = pf.build_sipm_responses(indices, times, widths,
+                                                   wfs, sipm_ids, 1, sipm_selection_algo)
 
     expected_ids = np.delete(      ids, below_thr_index)
     expected_wfs = np.delete(wfs_slice, below_thr_index, axis=0)
@@ -404,11 +404,11 @@ def test_build_peak_development(pmt_and_sipm_wfs_with_indices,
 
     peak = pf.build_peak(pmt_indices, times,
                          widths, pmt_wfs, pmt_ids, sipm_ids,
-                         rebin_stride = rebin,
-                         with_sipms   = with_sipms,
-                         Pk           = Pk,
-                         sipm_wfs     = sipm_wfs,
-                         apply_cut    = select_cutting_algorithm(CutAlgo.no_cut)
+                         rebin_stride        = rebin,
+                         with_sipms          = with_sipms,
+                         Pk                  = Pk,
+                         sipm_wfs            = sipm_wfs,
+                         sipm_selection_algo = select_cutting_algorithm(CutAlgo.no_cut)
                          )
 
     assert_Peak_equality(peak, expected_peak)
@@ -486,8 +486,8 @@ def test_find_peaks_s2_style(pmt_and_sipm_wfs_with_indices):
                           time_range, length_range,
                           stride, rebin_stride,
                           S2, pmt_ids, sipm_ids,
-                          sipm_wfs    = sipm_wfs,
-                          apply_cut    = select_cutting_algorithm(CutAlgo.no_cut)
+                          sipm_wfs = sipm_wfs,
+                          sipm_selection_algo = select_cutting_algorithm(CutAlgo.no_cut)
                           )
 
     (rebinned_times,
@@ -519,7 +519,7 @@ def test_get_pmap(s1_and_s2_with_indices):
                        pmt_ids     = pmt_ids, sipm_ids = sipm_ids,
                        pmt_samp_wid  = pmt_samp_wid ,
                        sipm_samp_wid = sipm_samp_wid,
-                       apply_cut    = select_cutting_algorithm(CutAlgo.no_cut)
+                       sipm_selection_algo = select_cutting_algorithm(CutAlgo.no_cut)
                        )
 
     (rebinned_times ,

From 3a1f8e7da9638f0bc7dfcf7b39b9bf1b7862a1d6 Mon Sep 17 00:00:00 2001
From: Ian0sborne <ian.osborne@student.manchester.ac.uk>
Date: Thu, 16 Jul 2026 16:53:10 +0100
Subject: [PATCH 66/69] Fix parameter ordering

---
 invisible_cities/cities/irene.py | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/invisible_cities/cities/irene.py b/invisible_cities/cities/irene.py
index 8856538030..e7776e741f 100644
--- a/invisible_cities/cities/irene.py
+++ b/invisible_cities/cities/irene.py
@@ -74,9 +74,9 @@ def irene( files_in         : OneOrManyFiles
          , s2_rebin_stride  : int  , s2_stride    : int
          , thr_csum_s2      : float
          , pmt_samp_wid     : float, sipm_samp_wid: float
-         , store_db         : bool = True
          , cutting_function : CutAlgo
          , cutting_params   : Dict[str, Any]
+         , store_db         : bool = True
          ):
     '''
     `cutting_function` is defined within components.py, and can vary, resulting

From 347192344b839c5c4d86dbcecc8c35c45a4ca011 Mon Sep 17 00:00:00 2001
From: Ian0sborne <ian.osborne@student.manchester.ac.uk>
Date: Thu, 23 Jul 2026 14:00:02 +0100
Subject: [PATCH 67/69] Fix indentation

---
 invisible_cities/cities/components.py |  4 ++--
 invisible_cities/types/symbols.py     | 12 ++++++------
 2 files changed, 8 insertions(+), 8 deletions(-)

diff --git a/invisible_cities/cities/components.py b/invisible_cities/cities/components.py
index e160b75593..2ba2109dcd 100644
--- a/invisible_cities/cities/components.py
+++ b/invisible_cities/cities/components.py
@@ -829,7 +829,7 @@ def select_cutting_algorithm(algo, **cutting_params):
         raise ValueError(f"Unsupported cutting algorithm: {algo!r}. Expected one of {list(CutAlgo)}")
 
 
-def threshold_sipm_selection(thr_sipm_type
+def threshold_sipm_selection( thr_sipm_type
                             , thr_sipm
                             , thr_sipm_s2
                             , run_number
@@ -849,7 +849,7 @@ def threshold_sipm_selection(wfs):
     return threshold_sipm_selection
 
 
-def pyrrha_sipm_selection(selection_method     : SiPMSelectionMethod
+def pyrrha_sipm_selection( selection_method     : SiPMSelectionMethod
                          , selection_kwargs    : dict
                          , proximity_threshold : float
                          , padding_radius      : float
diff --git a/invisible_cities/types/symbols.py b/invisible_cities/types/symbols.py
index 4dbd642f73..d31b63454e 100644
--- a/invisible_cities/types/symbols.py
+++ b/invisible_cities/types/symbols.py
@@ -132,12 +132,12 @@ class SensorType(AutoNameEnumBase):
 
 
 class SiPMCalibMode(AutoNameEnumBase):
-    subtract_mode                        = auto()
-    subtract_median                      = auto()
-    subtract_mode_calibrate              = auto()
-    subtract_mean_calibrate              = auto()
-    subtract_median_calibrate            = auto()
-    subtract_baseline_calibrate          = auto()
+    subtract_mode               = auto()
+    subtract_median             = auto()
+    subtract_mode_calibrate     = auto()
+    subtract_mean_calibrate     = auto()
+    subtract_median_calibrate   = auto()
+    subtract_baseline_calibrate = auto()
 
 
 class SiPMCharge(AutoNameEnumBase):

From 6116960a09a8c422c8b7f1541fa35162fca0702a Mon Sep 17 00:00:00 2001
From: Ian0sborne <ian.osborne@student.manchester.ac.uk>
Date: Thu, 23 Jul 2026 14:08:42 +0100
Subject: [PATCH 68/69] Improve docstrings

---
 invisible_cities/cities/components.py | 9 ++++++---
 1 file changed, 6 insertions(+), 3 deletions(-)

diff --git a/invisible_cities/cities/components.py b/invisible_cities/cities/components.py
index 2ba2109dcd..751b43b8a7 100644
--- a/invisible_cities/cities/components.py
+++ b/invisible_cities/cities/components.py
@@ -835,8 +835,11 @@ def threshold_sipm_selection( thr_sipm_type
                             , run_number
                             , detector_db = None):
     '''
-    Function that applies thresholding to the sipms in standard irene manner,
-    by zeroing all waveform values below a threshold.
+    Applies energy thresholds to SiPM S2 waveforms:
+    - thr_sipm:    applied per time bin. Waveform samples below this 
+                   threshold are set to zero.  
+    - thr_sipm_s2: applied to the integrated waveform charge. If the total
+                   charge is below this threshold, it is set to zero.
     '''
     if detector_db is None:
         sipm_thr = thr_sipm
@@ -856,7 +859,7 @@ def pyrrha_sipm_selection( selection_method     : SiPMSelectionMethod
                          , run_number          : int
                          , detector_db         : str):
     '''
-    Function that applies a generic selection function to the sipms, which can be used to 
+    Applies a generic selection function to the sipms, which can be used to 
     implement a spatial SiPM selection method (called Pyrrha).
     '''
     def pyrrha_sipm_selection(wfs):

From 814e1dcba7ca81bdca053e5c04d786be11be01d2 Mon Sep 17 00:00:00 2001
From: Ian0sborne <ian.osborne@student.manchester.ac.uk>
Date: Thu, 23 Jul 2026 14:14:35 +0100
Subject: [PATCH 69/69] Remove unused imports

---
 invisible_cities/reco/wfm_functions.py | 3 +--
 1 file changed, 1 insertion(+), 2 deletions(-)

diff --git a/invisible_cities/reco/wfm_functions.py b/invisible_cities/reco/wfm_functions.py
index 3cfa92dc7e..d13c15115c 100644
--- a/invisible_cities/reco/wfm_functions.py
+++ b/invisible_cities/reco/wfm_functions.py
@@ -4,10 +4,9 @@
 """
 import numpy as np
 from typing  import Optional
-from typing  import Callable
 from typing  import Tuple
 
-from .. core.core_functions import define_window, to_col_vector
+from .. core.core_functions import define_window
 from .. calib               import calib_sensors_functions as csf
 from .. sierpe              import blr
 from .. database            import load_db