Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Friedmann Diaconis bins for SiPM threshold determination #115

Merged
merged 2 commits into from
Mar 12, 2025
Merged

Friedmann Diaconis bins for SiPM threshold determination #115

Changes from all commits
Commits
Show all changes
2 commits
Select commit Hold shift + click to select a range
93e10f1
friedman diaconis bins
Mar 12, 2025
1437d9c
raise error if number of wfs smaller than required value
Mar 12, 2025
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
19 changes: 13 additions & 6 deletions workflow/src/legenddataflow/scripts/par/spms/dsp/trigger_threshold.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -73,6 +73,12 @@ def par_spms_dsp_trg_thr() -> None:
)
log.warning(msg)

elif len(data) < settings.n_events:
msg = (
f"number of waveforms '{args.raw_table_name}/waveform_bit_drop' < {settings.n_events}"
"in {args.raw_file}, can't build histogram"
)
raise RuntimeError(msg)
else:
# run the DSP with the provided configuration
log.debug("running the DSP chain")
Expand All @@ -87,13 +93,14 @@ def par_spms_dsp_trg_thr() -> None:
# determine a cutoff for the histogram used to extract the FWHM
low_cutoff, high_cutoff = np.quantile(wf_current, [0.005, 0.995])

# make histogram of the curr values
h = (
hist.new.Regular(settings.n_baseline_bins, low_cutoff, high_cutoff)
.Double()
.fill(wf_current)
# determine hist edges with Friedmann Diaconis Estimator
bin_edges = np.histogram_bin_edges(
wf_current, bins="fd", range=(low_cutoff, high_cutoff)
)

# make histogram of the curr values
h = hist.new.Variable(bin_edges).Double().fill(wf_current)

# determine FWHM
counts = h.view()
idx_over_half = np.where(counts >= np.max(counts) / 2)[0]
Expand All @@ -102,7 +109,7 @@ def par_spms_dsp_trg_thr() -> None:
fwhm = edges[idx_over_half[-1]] - edges[idx_over_half[0]]

if fwhm <= 0:
msg = "determined FWHM of baseline derivative distribution is zero or negative"
msg = f"determined FWHM of baseline derivative distribution is so <= 0: {fwhm:.3f}"
raise RuntimeError(msg)

log.debug(f"writing out baseline_curr_fwhm = {fwhm}")
Expand Down