Skip to main content
The PDS (Photon Detection System) pipeline reconstructs scintillation light signals from the liquid argon detector. It is activated with the -pds flag and runs in workflow.py via pds_signal_proc and pds_reco.
PDS reconstruction is only configured for pdvd and cbbot. The pdhd detector explicitly does not support PDS reconstruction. Other detectors (cbtop, cb1top, cb1bot, dp, 50l) do not have PDS geometry configured in their settings files.

PDS Data Types

LARDON handles two types of PDS data simultaneously:
TypeVariableDescription
Streamingdc.data_stream_pdsContinuous readout data with timestamps relative to delay_pds_stream_time
Triggereddc.data_trig_pdsSelf-triggered data with timestamps relative to delay_pds_trig_time
Each type has its own channel arrays, sample counts (cf.n_pds_stream_sample, cf.n_pds_trig_sample), and mask arrays. Processing is applied identically to both. If either sample count is ≤ 0, PDS processing is skipped entirely.

PDS Signal Processing

pds_signal_proc()
1

First pedestal pass

ped.compute_pedestal_pds(first=True)
The first pass removes the waveform median and builds an initial signal mask:
s_med = bn.median(dc.data_stream_pds, axis=1)
dc.data_stream_pds -= s_med[:, None]
dc.mask_stream_pds = (data <= adc_thresh) | (data <= -abs(baseline) + 10)
The mask marks samples as baseline (True) if they are below raw_adc_thresh = 200 ADC or below a symmetric threshold around zero. This is a simple amplitude cut rather than an RMS-based cut.After the initial masking, the iterative pedestal computation from compute_pedestal_nb runs for n_iter = 3 iterations:
for i in range(n_iter):
    dc.mask_stream_pds = (data <= rms_thresh * rms) | (data <= -abs(baseline) + 10)
    s_mean, s_std = compute_pedestal_nb(dc.data_stream_pds, dc.mask_stream_pds, False)
    dc.data_stream_pds -= s_mean[:, None]
Each iteration refines the baseline estimate using rms_thresh = 3 as the mask threshold (samples within 3σ of zero are treated as baseline). The same procedure runs identically for trigger data.
"pds": {
  "pedestal": {
    "raw_adc_thresh": 200,
    "rms_thresh": 3,
    "n_iter": 3
  }
}
2

Median filter baseline flattening

noise.median_filter_pds()
A centered sliding median filter flattens slow baseline variations in both streaming and triggered data:
data_stream_masked = np.where(dc.mask_stream_pds, dc.data_stream_pds, np.nan)
baseline = centered_median_filter(data_stream_masked, window)
np.nan_to_num(baseline, copy=False, nan=0.0)
dc.data_stream_pds -= baseline
ROI samples (where the mask is False) are replaced with NaN before computing the median, so the filter only slides over baseline samples. The window size is flat_baseline_window = 4000 samples.This step is applied to streaming data if window ≤ n_pds_stream_sample, and to triggered data if window ≤ n_pds_trig_sample.
"noise": {
  "flat_baseline_window": 4000
}
3

Second pedestal pass

ped.compute_pedestal_pds(first=False)
After filtering, the pedestal refinement loop runs again (n_iter = 3 iterations) on both data types. The result is stored as noise_pds_filt on the event object and used as the noise baseline for PDS hit finding.

PDS Hit/Peak Finding

hf.find_all_pds_peak()
Peak finding is applied to both streaming and triggered data via find_pds_peak("stream") and find_pds_peak("trigger"). The algorithm is the same collection-type hit finder used for TPC channels (hit_search_collection_nb), applied to each PDS channel’s ROI:
  1. Contiguous ROI intervals (where the mask is False) are identified using np.diff.
  2. Intervals shorter than dt_min = 50 samples are discarded.
  3. Left padding of pad_left = 30 and right padding of pad_right = 50 samples are applied.
  4. The padded excerpt is passed to hit_search_collection_nb with thresholds:
    • thr1 = amp_sig[0] × RMS = 5 × RMS (lower threshold, starts the peak)
    • thr2 = amp_sig[1] × RMS = 8 × RMS (upper threshold, used to split overlapping peaks)
    • Absolute minimum: 0.5 ADC
"hit_finder": {
  "amp_sig": [5, 8],
  "dt_min": 50,
  "pad": { "left": 30, "right": 50 }
}
Each found peak is stored as a dc.pds_peak object with:
  • Global channel number (glob_ch)
  • Module assignment
  • Start, stop, and peak-time sample indices
  • Peak amplitude (max_adc)
  • Integrated charge (sum of ADC over the padded window)
  • Absolute timestamp: delta_time_ref + start / pds_sampling
The delta_time_ref is the time offset of the PDS readout window relative to the trigger, so timestamps are expressed in absolute detector time (μs).

Light Clustering

clu.light_clustering()
Light clustering groups PDS peaks from different channels that correspond to the same physical particle interaction. Step 1: Waveform alignment (align_waveforms) Peak timestamps are cross-correlated between channels to account for relative time offsets between PDS readout channels. For each channel, the best lag relative to channel 0 is found using sparse_crosscorr:
lag_scores = np.zeros(2 * max_lag + 1, dtype=int)
for i, lag in enumerate(lag_range):
    shifted = peaks_ch + lag / cf.pds_sampling
    lag_scores[i] = np.sum(np.isin(shifted, peaks_ref))
best_lag = lag_range[np.argmax(lag_scores)]
The search window is max_lag = 50 PDS ticks. Per-channel offsets are stored in dc.evt_list[-1].pds_time_offset. Step 2: Time-based clustering An R-tree is built over all PDS peaks indexed by (channel, corrected_timestamp). For each unclustered peak, all other peaks within time_tol = 5 μs are retrieved:
overlaps = list(rtree_idx.intersection(
    (0, i_start - time_tol, 9999, i_start + time_tol)
))
Peaks from the same channel as the query peak are excluded. Any group of two or more peaks from different channels is formed into a dc.pds_cluster object.
"cluster": {
  "max_lag": 50,
  "time_tol": 5
}
Each cluster stores the list of contributing peak IDs, global channel numbers, individual timestamps, charges, and an overall cluster timestamp (minimum of all peak timestamps).

TPC–PDS Matching

mat.matching_trk_pds()
When both -trk and -pds are active, reconstructed 3D tracks are matched to light clusters. This function is called from match_charge_and_pds in workflow.py.
This step requires both dc.tracks3D_list and dc.pds_cluster_list to be non-empty and runs after track timing has been computed.
Track classification 3D tracks are categorised by their crossing properties, and each category is tested against light clusters with a different time tolerance:
CategoryConditionTime tolerance
Anode-cathode crossersis_anode_crosser and is_cathode_crosser±5 μs (anode settings)
Anode crossersis_anode_crosser only±5 μs
Cathode crossersis_cathode_crosser only±5 μs
UnknownNeither flag set±5 μs
Default tolerances from the configuration:
"tpc_matching": {
  "min_cluster_size": 2,
  "trigger":         { "time_tol_bef": [5], "time_tol_aft": [5] },
  "anode_crosser":   { "time_tol_bef": [5], "time_tol_aft": [5] },
  "cathode_crosser": { "time_tol_bef": [5, 5], "time_tol_aft": [5, 5] },
  "unknown":         { "time_tol_bef": [5, 5], "time_tol_aft": [5, 5] }
}
The per-volume arrays (length 2 for cathode/unknown) select the tolerance based on which drift volume (trk_vol = module_ini / n_drift_volumes) contains the track. Matching procedure An R-tree is built over light clusters indexed by timestamp. For each track, a time window query retrieves candidate clusters:
pds_overlaps = list(rtree.intersection(
    (trk_start - time_tol_before[trk_vol], 0,
     trk_stop  + time_tol_after[trk_vol],  0)
))
Only clusters with size >= min_cluster_size = 2 and not yet matched to another track are considered. A match is recorded only when exactly one cluster falls within the time window. For cathode-crossing track pairs, both tracks in the pair are matched to the same cluster. Post-match timing correction When a match is found, the track’s absolute Z position is recalibrated using the light cluster timestamp:
t.set_times_from_light(clus.timestamp, v_drift[t.module_ini])
The distance between the track and each contributing PDS channel is also computed and stored on the cluster object for downstream analysis.

Build docs developers (and LLMs) love