C++ API

The public phepex:: kernels, grouped by header. Include everything via the umbrella header <phepex/phepex.hpp>, or a single header for one kernel. Each function’s description below is rendered from its own header documentation.

Deconvolution & upsampling

Declared in <phepex/deconvolve.hpp>.

void phepex::deconvolve_upsample(const float *waveforms, int n_ch, int n_pix, int n_samples, int upsampling, const float *pole_zero, const float *baseline, const float *scale, float *out)

Pole-zero deconvolution + upsampling of PMT/SiPM waveforms.

Applies, independently to every (channel, pixel) waveform: upsampling by upsampling, baseline subtraction and scaling, a single-pole-decay (pole-zero) correction, and two upsampling-wide moving-average smoothings (a forward+backward boxcar pass, so the smoothing introduces no net time shift).

pole_zero, baseline and scale are per-(channel, pixel) arrays of n_ch*n_pix floats, indexed row-major like waveforms. Each may be null, in which case a constant default is used for every pixel: pole_zero = 0 (no pole-zero correction), baseline = 0, scale = 1.

The smoothing produces 2*(upsampling-1) invalid samples at each end of every output waveform; where pole_zero != 0, the first 3*upsampling-2 samples are invalid (see deconvolve_valid_range() below).

Parameters:

  • waveforms – input, row-major (n_ch, n_pix, n_samples), float32

  • pole_zero – per-(channel, pixel) pole-zero factors, or null for 0 (no correction)

  • baseline – per-(channel, pixel) baselines subtracted before scaling, or null for 0

  • scale – per-(channel, pixel) scale factors, or null for 1

  • out – caller-allocated output of n_ch*n_pix*n_samples*upsampling floats, row-major (n_ch, n_pix, n_samples*upsampling)

struct SampleRange

Half-open sample range [lo, hi).

Public Members

int lo

lower limit (inclusive)

int hi

upper limit (exclusive)

SampleRange phepex::deconvolve_valid_range(int upsampling, int n_samples, double pole_zero)

Trustworthy (non-edge) output-sample range of deconvolve_upsample(): the upsample+filtfilt boxcar contaminates 2*(upsampling-1) samples at each end, and a non-zero pole_zero adds upsampling more invalid samples at the start.

Returns the range into the length upsampling*n_samples deconvolution output.

Waveform preprocessing

Declared in <phepex/preprocess.hpp>. Single-waveform upsampling + pole-zero deconvolution with optional Deriche (1992) Gaussian smoothing, plus its DVR-convention valid range.

struct SmoothingCoefficients

Coefficients for the delay-compensated second-order IIR (Deriche 1992) Gaussian-like smoothing applied by preprocess_waveform().

Public Members

double fwhm

target full width at half maximum of the impulse response

std::array<double, 2> n

numerator coefficients for the forward filter step

std::array<double, 2> m

numerator coefficients for the backward filter step

std::array<double, 2> d

denominator coefficients

SmoothingCoefficients phepex::calculate_smoothing_coefficients(double smoothing_fwhm)

Compute the smoothing coefficients for a given impulse-response FWHM (in samples).

void phepex::preprocess_waveform(const std::uint16_t *src, int n_samples, int upsampling, float pole_zero, const SmoothingCoefficients *smoothing, float offset, float scale, float *out, float *scratch = nullptr)

Upsampling + pole-zero deconvolution of a single waveform, with optional Gaussian smoothing.

Repeats each of the n_samples inputs upsampling times, subtracts offset and applies scale, corrects a single-pole decay pole_zero, and smooths with two upsampling-wide moving averages; when smoothing != nullptr an additional Deriche IIR pass is applied. upsampling == 1 degenerates to scale*(src - offset).

Writes n_samples*upsampling floats to out (caller-allocated). All arithmetic is float32; offset/scale/pole_zero are float on purpose (bit-exactness).

Parameters:

  • src – input waveform of n_samples samples

  • smoothing – optional Deriche coefficients; nullptr disables smoothing

  • out – caller-allocated output of n_samples*upsampling floats

  • scratch – optional caller-provided workspace of at least n_samples*upsampling floats, used only when smoothing != nullptr; if null a buffer is allocated internally for the duration of the call

void phepex::preprocess_waveform(const float *src, int n_samples, int upsampling, float pole_zero, const SmoothingCoefficients *smoothing, float offset, float scale, float *out, float *scratch = nullptr)
SampleRange phepex::preprocess_valid_range(int upsampling, float pole_zero, const SmoothingCoefficients *smoothing, int num_samples)

Trustworthy (non-edge) sample range of preprocess_waveform(): the boxcar smoothing contaminates 2*(upsampling-1) samples at each end, a non-zero pole_zero adds upsampling more invalid samples at the start, and smoothing widens both margins by floor(fwhm).

The upper bound follows the DVR convention num_samples - right (using the RAW, pre-upsample num_samples) — this differs from deconvolve_valid_range(). Returns {0, 0} if the margins exceed num_samples.

Clipping

Declared in <phepex/clip.hpp>.

void phepex::pos_soft_clip(const float *waveforms, int n_ch, int n_pix, int n_up, float scale, int sample_lo, int sample_hi, float *out)

Positive soft clip: max(y/(1+|y|), 0) with y = x/scale, applied over samples [sample_lo, sample_hi) of each (channel, pixel) waveform; samples outside that range are set to 0.

Passing sample_lo == sample_hi == 0 means the full trace. The soft clip already bounds the result to (-1,1), so only negatives are clamped (to 0).

Parameters:

  • waveforms – input (n_ch, n_pix, n_up) float32

  • out – caller-allocated n_ch*n_pix*n_up floats (fully written; 0 outside range)

Neighbour peak finding

Declared in <phepex/neighbor.hpp>.

void phepex::neighbor_peak_indices(const float *waveforms, int n_ch, int n_pix, int n_up, const std::int32_t *indptr, const std::int32_t *indices, int local_weight, const bool *broken_pixels, int sample_lo, int sample_hi, std::int64_t *peak_out, std::int32_t *neighbor_count = nullptr, float *scratch = nullptr)

Per-pixel peak index (over samples in [sample_lo, sample_hi)) of the neighbour-summed waveform: average = waveforms[pixel]*local_weight + sum over non-broken neighbours.

First-max tie-break, float32 accumulation, no normalisation. Passing sample_lo == sample_hi == 0 means the full trace.

Neighbours are given as a CSR adjacency (compressed sparse row: row pointers + column indices): pixel p’s neighbours are indices[indptr[p] .. indptr[p+1]).

Parameters:

  • waveforms – input (n_ch, n_pix, n_up) float32

  • indptr – CSR row pointers, length n_pix+1

  • indices – CSR column indices (neighbour pixel ids), length indptr[n_pix]

  • broken_pixels – (n_ch, n_pix) bool; broken neighbours are skipped

  • peak_out – caller-allocated n_ch*n_pix int64 (argmax sample index per pixel)

  • neighbor_count – optional n_ch*n_pix int32 (caller-allocated); if non-null, gets the number of non-broken neighbours summed per pixel

  • scratch – optional caller-provided workspace of at least n_up floats; if null, a buffer is allocated internally for the duration of the call

Charge & timing extraction

Declared in <phepex/extract.hpp>.

void phepex::extract_around_peak(const float *waveforms, int n_ch, int n_pix, int n_up, const std::int64_t *peak_index, int width, int shift, double sampling_rate_ghz, float *charge, float *peak_time)

Window integration + amplitude-weighted peak time, per (channel, pixel).

Sums the waveform over [peak-shift, peak-shift+width) (clamped to the trace) and computes the amplitude-weighted centroid over strictly-positive samples, in ns (divided by sampling_rate_ghz). float64 accumulation; float32 outputs.

Parameters:

  • waveforms – input (n_ch, n_pix, n_up) float32

  • peak_index – (n_ch, n_pix) int64 peak sample index per pixel

  • charge – caller-allocated n_ch*n_pix float32 (integrated charge)

  • peak_time – caller-allocated n_ch*n_pix float32 (peak time, ns)

void phepex::adaptive_centroid(const float *waveforms, int n_ch, int n_pix, int n_up, const std::int64_t *peak_index, double rel_descend_limit, float *centroids)

Leading-edge weighted centroid (in sample units), per (channel, pixel): walks left then right from peak while samples exceed rel_descend_limit*peak_amplitude, accumulating an amplitude-weighted index centroid.

float64 accumulation; float32 output.

Parameters:

  • waveforms – input (n_ch, n_pix, n_up) float32

  • peak_index – (n_ch, n_pix) int64 peak sample index per pixel

  • centroids – caller-allocated n_ch*n_pix float32 (centroid, sample units)

Waveform generation

Declared in <phepex/generate.hpp>.

void phepex::generate_waveforms(const double *charge, const double *time_ns, int n_events, int n_pix, const double *reference_pulse, int n_ref, double ref_sample_width_ns, double sample_width_ns, int n_samples, int upsampling, double nsb_rate_ghz, std::uint64_t seed, float *out)

Synthesize sampled waveforms from per-pixel (charge, peak_time): for the signal component, each pixel’s charge is deposited at its sub-sample peak time and convolved with the reference pulse shape; a temporally-distributed Poisson night-sky-background (NSB) process at nsb_rate_ghz is added on top.

Intended as a fast generator for testing/benchmarking extractors.

Parameters:

  • charge – (n_events, n_pix) float64 photo-electron charge per pixel

  • time_ns – (n_events, n_pix) float64 pulse peak time (ns) per pixel

  • reference_pulse – (n_ref) float64 single-channel reference pulse shape

  • ref_sample_width_ns – sample width of the reference pulse (ns)

  • sample_width_ns – readout sample width (ns) = 1 / sampling_rate

  • upsampling – internal upsampling used to place sub-sample deposits (e.g. 10)

  • nsb_rate_ghz – NSB rate in GHz (0 disables NSB)

  • seed – RNG seed (per-event deterministic)

  • out – caller-allocated n_events*n_pix*n_samples float32, row-major (n_events, n_pix, n_samples)