LegendDSP

Template for Julia packages.

DSPConfig{T <: Real}

Configuration parameters for DSP algorithms.

Fields

• enc_pickoff::Quantity{<:T}: pick-off time for ENC noise calculations
• bl_window::ClosedInterval{Quantity{<:T}}: fit window for basline extraction
• tail_window::ClosedInterval{Quantity{<:T}}: fit window for decay time extraction
• inTraceCut_std_threshold::T: in-trace pile-up rejector threshold in standard deviations
• t0_threshold::T: ADC threshold for t0 determination
• e_grid_rt_trap::StepRangeLen{Quantity{<:T}}: rise time grid scan range for trapezoidal filter
• e_grid_ft_trap::StepRangeLen{Quantity{<:T}}: flat-top time grid scan range for trapezoidal filter
• e_grid_rt_zac::StepRangeLen{Quantity{<:T}}: rise time grid scan range for ZAC filter
• e_grid_ft_zac::StepRangeLen{Quantity{<:T}}: flat-top time grid scan range for ZAC filter
• e_grid_rt_cusp::StepRangeLen{Quantity{<:T}}: rise time grid scan range for CUSP filter
• e_grid_ft_cusp::StepRangeLen{Quantity{<:T}}: flat-top time grid scan range for CUSP filter
• a_grid_rt_sg::StepRangeLen{Quantity{<:T}}: window length grid scan range for SG filter

Examples

using LegendDSP
using LegendDataManagement

l200 = LegendData(:l200)
filekey = start_filekey(l200, (:p03, :r000, :cal))
dsp_config = DSPConfig(dataprod_config(l200).dsp(filekey).default)

struct IntersectMaximum <: Function

Finds the intersects of a Y with a threshold and picking the maximum in a given time window. Constructors:

• IntersectMaximum(; fields...)

Fields:

• mintot::Union{Unitful.Quantity{<:var"#s10"}, var"#s10"} where var"#s10"<:Real: minimum time-over-threshold

• maxtot::Union{Unitful.Quantity{<:var"#s10"}, var"#s10"} where var"#s10"<:Real: maximum time-over-threshold for max to appear

struct ModifiedSincFilter{T<:RealQuantity} <: AbstractRadFIRFilter{}

A Modified-Sinc filter Working example

using RadiationDetectorSignals
using Unitful

n = 600
noise = 1.
t = range(0u"μs", 20u"μs", 2*n)
signal = vcat(zeros(n), 10*ones(n)) + (noise*rand(2*n) .- noise/2)
wf = RDWaveform(t, signal)

# define filter parameters and filter
flt = ModifiedSincFilter(d=4, m=1u"μs")

# apply filter to signal
wf_new = flt(wf)

Constructors:

• ModifiedSincFilter(; fields...)

Fields:

• d::Int64: degree of the filter determining the number of extrema in the kernel

• m::Union{Real, Unitful.AbstractQuantity{<:Real}}: half-width of the kernel

struct MovingWindowFilter <: AbstractRadFIRFilter

apply left moving average window to signal. The exact computations are:
yₙ = yₙ₋₁ + (xₙ - x₁)/l, for n ∈ {1, …, l}
yₙ = yₙ₋₁ + (xₙ - xₙ₋ₗ)/l, for n ∈ {l+1, …, L}

Constructors:

• MovingWindowFilter(; fields...)

Fields:

• length::Union{Unitful.Quantity{<:var"#s10"}, var"#s10"} where var"#s10"<:Real: length of the moving average window

struct MovingWindowMultiFilter <: AbstractRadFIRFilter

apply left and right moving average windows to signal. Working example:

using RadiationDetectorSignals
using Unitful

signal = rand(100)
t = range(0u"ms", 20u"ms", 100)

wf = RDWaveform(t, signal)
flt = MovingWindowMultiFilter(1u"ms")
wf_new = flt(wf)

Constructors:

• MovingWindowMultiFilter(; fields...)

Fields:

• length::Union{Unitful.Quantity{<:var"#s10"}, var"#s10"} where var"#s10"<:Real: size of the moving average window

struct MultiIntersect <: Function Finds the x values at which a signal exceeds the provided threshold_ratios the first time. Constructors: * MultiIntersect(; fields...) Fields:

• threshold_ratios::Vector{Float64}: ratios which determine the thresholds

• mintot::Union{Unitful.Quantity{<:var"#s10"}, var"#s10"} where var"#s10"<:Real: minimum time-over-threshold

• n::Int64: half window length of polynomial fit

• d::Int64: degree of polynomial

• sampling_rate::Int64: upsampling rate

struct WeightedSavitzkyGolayFilter{T<:RealQuantity} <: AbstractRadFIRFilter

A Weighted-Savitzky-Golay filter. Working example:

using RadiationDetectorSignals
using Unitful

n = 600
noise = 1.
t = range(0u"μs", 20u"μs", 2*n)
signal = vcat(zeros(n), 10*ones(n)) + (noise*rand(2*n) .- noise/2)
wf = RDWaveform(t, signal)

# define filter parameters and filter
flt = WeightedSavitzkyGolayFilter(length=1u"μs", degree=3, weightType=2)

# apply filter to signal
wf_new = flt(wf)

Constructors:

• WeightedSavitzkyGolayFilter(; fields...)

Fields:

• length::Union{Real, Unitful.AbstractQuantity{<:Real}}: total filter length

• degree::Int64: Polynomial degree

• weightType::Int64: weight function to use

struct WhittakerHendersonFilter <: AbstractRadSigFilter{LinearFiltering}

A Whittaker-Henderson filter. Working example: ```julia using RadiationDetectorSignals using Unitful

n = 600 noise = 1. t = range(0u"μs", 20u"μs", 2n) signal = vcat(zeros(n), 10ones(n)) + (noiserand(2n) .- noise/2) wf = RDWaveform(t, signal)

define filter parameters and filter

flt = WhittakerHendersonFilter(p=3, λ=4)

apply filter to signal

wf_new = flt(wf)

_create_dsp_config(dsp_metadata::PropDicts.PropDict)

Create a DSPConfig from a PropDict of DSP metadata.

Arguments

• dsp_metadata::PropDicts.PropDict: DSP metadata

Returns

• dsp_config::DSPConfig: DSP configuration

banded_cholesky_solve!(x::AbstractVector, A::AbstractMatrix, y::AbstractVector)

solve the linear equation AA'x = y and store the result in x, where A is the cholesky decomposition of a banded centro symmetric matrix. A[1, i] is the diagonal, A[2, i] is the first subdiagonal and so on.

choleskyL!(b::AbstractMatrix{T})

inplace cholesky decomposition of a banded symmetric matrix, where b[1, i] contains the diagonal elements, b[2, i] the elements of the first subdiagonal and so on.

dsp_cusp_ft_optimization(wvfs::ArrayOfRDWaveforms, config::DSPConfig, τ::Quantity{T}, rt::Quantity{T}) where T<:Real

Get energy grid values for given CUSP grid rise times while varying the flat-top time.

Returns

- `e_grid`: Array energy values for the given CUSP rise time grid at a given flat-top time grid.

dsp_cusp_rt_optimization(wvfs::ArrayOfRDWaveforms, config::DSPConfig, τ::Quantity{T},; ft::Quantity{T}=4.0u"µs") where T<:Real

Get ENC noise grid values for given CUSP grid rise times.

Returns

- `enc_cusp_grid`: Array ENC noise values for the given CUSP rise time grid

dsp_decay_times(wvfs::ArrayOfRDWaveforms, bl_window::ClosedInterval{<:Unitful.Time{<:T}}, tail_window::ClosedInterval{<:Unitful.Time{<:T}})
dsp_decay_times(wvfs::ArrayOfRDWaveforms, config::DSPConfig)

Get statistics on the logarhithmic of the tail of the wvfs in the interval tail_window.

Returns

• τ: decay time in µs

dsp_icpc(data::Q, config::DSPConfig, τ::Quantity{T}, pars_filter::PropDict) where {Q <: Table, T<:Real}

DSP for ICPC detectors. It needs the decay time τ of the detector and the filter parameters pars_filter for the optimal filter parameters for the Trap, CUSP and ZAC filter.

Input data

The input data is a table with the following columns:

• waveform: waveform data
• baseline: baseline from FADC
• timestamp: timestamp from FADC
• eventnumber: event ID from FADC
• daqenergy: energy from FADC

Output data

The output data is a table with the following columns:

• blmean: baseline mean
• blsigma: baseline sigma
• blslope: baseline slope
• bloffset: baseline offset
• tailmean: tail mean after PZ correction
• tailsigma: tail sigma after PZ correction
• tailslope: tail slope after PZ correction
• tailoffset: tail offset after PZ correction
• t0: start time of waveform drift
• t10: timepoint of 10% of waveform maximum
• t50: timepoint of 50% of waveform maximum
• t80: timepoint of 80% of waveform maximum
• t90: timepoint of 90% of waveform maximum
• t99: timepoint of 99% of waveform maximum
• t50_current: timepoint of current rise to 50% of maximum
• tail_τ: tail decay time
• tail_mean: tail mean before PZ correction
• tail_sigma: tail sigma before PZ correction
• e_max: maximum of waveform
• e_min: minimum of waveform
• e_10410: energy of waveform with trapezoidal filter of 10µs rise time with 4µs flat-top
• e_313: energy of waveform with trapezoidal filter of 3µs rise time with 1µs flat-top
• e_10410_inv: maximum of inverted waveform with trapezoidal filter of 10µs rise time with 4µs flat-top
• e_313_inv: maximum of inverted waveform with trapezoidal filter of 3µs rise time with 1µs flat-top
• t0_inv: start time of inverted waveform drift
• e_trap: energy of waveform with trapezoidal filter of optimized rise and flat-top time
• e_cusp: energy of waveform with CUSP filter of optimized rise and flat-top time
• e_zac: energy of waveform with ZAC filter of optimized rise and flat-top time
• qdrift: Q-drift parameter
• lq: LQ parameter
• a: current maximum with optimal Savitzky-Golay filter length parameter
• blfc: baseline from FADC
• timestamp: timestamp from FADC
• eventID_fadc: event ID from FADC
• e_fc: energy from FADC
• pretrace_diff: difference between first sample and baseline mean
• drift_time: drift time between t0 and 90% of waveform maximum
• inTrace_intersect: position of in-trace pile-up
• inTrace_n: multiplicity of in-trace pile-up
• n_sat_low: number of samples the waveform is saturated at low of FADC range
• n_sat_high: number of samples the waveform is saturated at high of FADC range
• n_sat_low_cons: number of consecutive samples the waveform is saturated at low of FADC range
• n_sat_high_cons: number of consecutive samples the waveform is saturated at high of FADC range

dsp_icpc_compressed(data::Q, config::DSPConfig, τ::Quantity{T}, pars_filter::PropDict) where {Q <: Table, T<:Real}

DSP for ICPC detectors. It needs the decay time τ of the detector and the filter parameters pars_filter for the optimal filter parameters for the Trap, CUSP and ZAC filter.

Input data

The input data is a table with the following columns:

• waveform: waveform data
• baseline: baseline from FADC
• timestamp: timestamp from FADC
• eventnumber: event ID from FADC
• daqenergy: energy from FADC

Output data

The output data is a table with the following columns:

• blmean: baseline mean
• blsigma: baseline sigma
• blslope: baseline slope
• bloffset: baseline offset
• tailmean: tail mean after PZ correction
• tailsigma: tail sigma after PZ correction
• tailslope: tail slope after PZ correction
• tailoffset: tail offset after PZ correction
• t0: start time of waveform drift
• t10: timepoint of 10% of waveform maximum
• t50: timepoint of 50% of waveform maximum
• t80: timepoint of 80% of waveform maximum
• t90: timepoint of 90% of waveform maximum
• t99: timepoint of 99% of waveform maximum
• t50_current: timepoint of current rise to 50% of maximum
• tail_τ: tail decay time
• tail_mean: tail mean before PZ correction
• tail_sigma: tail sigma before PZ correction
• e_max: maximum of waveform
• e_min: minimum of waveform
• e_10410: energy of waveform with trapezoidal filter of 10µs rise time with 4µs flat-top
• e_313: energy of waveform with trapezoidal filter of 3µs rise time with 1µs flat-top
• e_10410_inv: maximum of inverted waveform with trapezoidal filter of 10µs rise time with 4µs flat-top
• e_313_inv: maximum of inverted waveform with trapezoidal filter of 3µs rise time with 1µs flat-top
• t0_inv: start time of inverted waveform drift
• e_trap: energy of waveform with trapezoidal filter of optimized rise and flat-top time
• e_cusp: energy of waveform with CUSP filter of optimized rise and flat-top time
• e_zac: energy of waveform with ZAC filter of optimized rise and flat-top time
• qdrift: Q-drift parameter
• lq: LQ parameter
• a: current maximum with optimal Savitzky-Golay filter length parameter
• blfc: baseline from FADC
• timestamp: timestamp from FADC
• eventID_fadc: event ID from FADC
• e_fc: energy from FADC
• pretrace_diff: difference between first sample and baseline mean
• drift_time: drift time between t0 and 90% of waveform maximum
• inTrace_intersect: position of in-trace pile-up
• inTrace_n: multiplicity of in-trace pile-up
• n_sat_low: number of samples the waveform is saturated at low of FADC range
• n_sat_high: number of samples the waveform is saturated at high of FADC range
• n_sat_low_cons: number of consecutive samples the waveform is saturated at low of FADC range
• n_sat_high_cons: number of consecutive samples the waveform is saturated at high of FADC range

dsp_puls(data::Q, config::DSPConfig) where {Q <: Table}

DSP function for pulser processing.

Input data

The input data is a table with the following columns:

• waveform: waveform data
• baseline: baseline data
• timestamp: timestamp data
• eventnumber: event number data
• daqenergy: energy data

Output data

The output data is a table with the following columns:

• blmean: baseline mean
• blsigma: baseline sigma
• blslope: baseline slope
• bloffset: baseline offset
• t50: timepoint of 50% of waveform maximum
• e_max: maximum of waveform
• e_10410: energy of waveform with trapezoidal filter of 10µs rise time with 4µs flat-top
• blfc: baseline from FADC
• timestamp: timestamp
• eventID_fadc: event number from FADC
• e_fc: energy from FADC

dsp_flt_optimization(wvfs::ArrayOfRDWaveforms, config::DSPConfig, τ::Unitful.Time{<:Real}, f_evaluate_qc::Function)

Get QC DSP for filter parameter optimization for a given waveform set.

dsp_qc_flt_optimization_compressed(wvfs::ArrayOfRDWaveforms, config::DSPConfig, τ::Unitful.Time{<:Real}, f_evaluate_qc::Function)

Get QC DSP for filter parameter optimization for a given waveform set.

dsp_qdrift_flt_optimization(wvfs::ArrayOfRDWaveforms, blmean::Vector{<:Real}, config::DSPConfig, τ::Unitful.Time{<:Real})

Get QDrift filter parameter optimization for a given waveform set.

dsp_sg_optimization(wvfs::ArrayOfRDWaveforms, config::DSPConfig, τ::Quantity{T}, pars_filter::PropDict) where T<:Real

Optimize the Savitzky-Golay filter parameters for a given waveform set.

Returns

- `aoe`: Array of efficiency values for the given Savitzky-Golay filter parameters
- `e`: Array of energy values for the given Savitzky-Golay filter parameters
- `blmean`: Baseline mean value
- `blslope`: Baseline slope value

dsp_sg_optimization_compressed(wvfs::ArrayOfRDWaveforms, config::DSPConfig, τ::Quantity{T}, pars_filter::PropDict) where T<:Real

Optimize the Savitzky-Golay filter parameters for a given waveform set.

Returns

- `aoe`: Array of efficiency values for the given Savitzky-Golay filter parameters
- `e`: Array of energy values for the given Savitzky-Golay filter parameters
- `blmean`: Baseline mean value
- `blslope`: Baseline slope value

dsp_sipm(data::Q, config::PropDict, pars_threshold::PropDict)

DSP routine for SiPM data. It needs the threshold parameters from the threshold scan for the given SiPM channel as well as the discharge threshold parameters.

Input data

The input data is a table with the following columns:

• waveform: waveform data
• baseline: baseline data
• timestamp: timestamp data
• eventnumber: event number data
• daqenergy: energy data

Output data

The output data is a table with the following columns:

• blfc: baseline from FADC
• timestamp: timestamp
• eventID_fadc: event number from FADC
• e_fc: energy from FADC
• trig_pos: trigger positions from DSP
• trig_max: trigger maxima from DSP
• trig_pos_DC: trigger positions of discharges
• trig_max_DC: trigger maxima of discharges

dsp_sipm_compressed(data::Q, config::PropDict, pars_threshold::PropDict)

DSP routine for SiPM data. It needs the threshold parameters from the threshold scan for the given SiPM channel as well as the discharge threshold parameters.

Input data

The input data is a table with the following columns:

• waveform: waveform data
• baseline: baseline data
• timestamp: timestamp data
• eventnumber: event number data
• daqenergy: energy data

Output data

The output data is a table with the following columns:

• blfc: baseline from FADC
• timestamp: timestamp
• eventID_fadc: event number from FADC
• e_fc: energy from FADC
• trig_pos: trigger positions from DSP
• trig_max: trigger maxima from DSP
• trig_pos_DC: trigger positions of discharges
• trig_max_DC: trigger maxima of discharges

dsp_sipm_thresholds_compressed(wvfs::ArrayOfRDWaveforms, config::PropDict)

This function calculates the baseline of the waveforms and the baseline of the waveforms with the sign flipped. The function is used to calculate the thresholds for the SiPMs.

Arguments

• wvfs::ArrayOfRDWaveforms: Array of RDWaveforms
• config::PropDict: Configuration parameters

Returns

• Table: Table with the baseline of the waveforms and the baseline of the waveforms with the sign flipped

dsp_trap_ft_optimization(wvfs::ArrayOfRDWaveforms, config::DSPConfig, τ::Quantity{T}, rt::Quantity{T}) where T<:Real

Get energy grid values for given trap grid rise times while varying the flat-top time.

Returns

- `e_grid`: Array energy values for the given trap rise time grid at a given flat-top time grid.

dsp_trap_rt_optimization(wvfs::ArrayOfRDWaveforms, config::DSPConfig, τ::Quantity{T},; ft::Quantity{T}=4.0u"µs") where T<:Real

Get ENC noise grid values for given trap grid rise times.

Returns

- `enc_trap_grid`: Array ENC noise values for the given trap rise time grid

dsp_zac_ft_optimization(wvfs::ArrayOfRDWaveforms, config::DSPConfig, τ::Quantity{T}, rt::Quantity{T}) where T<:Real

Get energy grid values for given ZAC grid rise times while varying the flat-top time.

Returns: - e_grid: Array energy values for the given ZAC rise time grid at a given flat-top time grid.

dsp_zac_rt_optimization(wvfs::ArrayOfRDWaveforms, config::DSPConfig, τ::Quantity{T},; ft::Quantity{T}=4.0u"µs") where T<:Real

Get ENC noise grid values for given ZAC grid rise times.

Returns

- `enc_zac_grid`: Array ENC noise values for the given ZAC rise time grid

extrema3points(y1::T,y2::T,y3::T)::T where {T<:AbstractFloat}

Calculate the extrema of a parabola defined by three points.

get_intracePileUp(wvfs::ArrayOfRDWaveforms, sigma_threshold::Real, bl_window::ClosedInterval{<:Unitful.Time{<:Real}}; mintot::Unitful.Time=100.0u"ns")

Get position and multiplicity of in-trace pile-up as intersect of reversed derivative signal with threshold as multiple of std. The wvfs have to be a current signal.

get_qc_classifier(wvfs::ArrayOfRDWaveforms, f_evaluate_qc::Function)

Get a classifier for quality cuts using the given waveforms and evaluation function which evulates a predefined SVM model.

get_qc_classifier_compressed(wvfs::ArrayOfRDWaveforms, f_evaluate_qc::Function)
get_qc_classifier_compressed(wvfs::ArrayOfRDWaveforms, f_evaluate_qc::Function, config::DSPConfig)

Get a classifier for quality cuts using the given waveforms and evaluation function which evulates a predefined SVM model.

get_qc_ml_func(dwts_norm::Matrix{<:Real}, dc_labels::Vector{<:Real}, hyperparams::PropDict)

Create a function that takes a signal and returns the prediction using the given hyperparameters and discrete wavelet transforms.

get_qdrift(wvfs::ArrayOfRDWaveforms, t_start::Array{Unitful.Time{T}}, Δt::UnitRange{Unitful.Time{T}}; pol_power::Int=3, sign_est_length::Unitful.Time=100u"ns")

Get the Q-drift parameter for each waveform in wvfs by integrating the waveform with gain = 1 and using a polynomial signal estimator of order pol_power and length sign_est_length to estimate the signal.

get_t0(wvfs_pz::ArrayOfRDWaveforms, t0_threshold::T) where T<:Real

Get t0 for each waveform in wvfs_pz by using a fast asymetric trapezoidal filter and a fixed threshold at t0_threshold. The filter is truncated to the range 0µs to 60µs where the Ge trigger is expected in FlashCam.

get_threshold(wvfs::ArrayOfRDWaveforms, threshold::Array{T}) where T<:Real

Get threshold for each waveform in wvfs by intersecting the waveform with a threshold per waveform or globally.

get_wvf_maximum(signal::AbstractSamples, start::Real, stop::Real)
get_wvf_maximum(signal::RDWaveform, start::RealQuantity, stop::RealQuantity)

Get the maximum of a signal in the interval (start,stop) by quadaratic interpolation.

makeFitWeights(d::Int, m::Int) where {U}

return the weights for the linear fit user for linear extrapolation at at the right boundary.

makeIPlusLambdaDprimeD(λ::T, p::Int, N::Int) where {T}

build the centro symmetric banded matrix which needs to be inverted later: I + λ*D'D D is the finite difference matrix of order p.

tailstats(signal::AbstractSamples, start::Real, stop::Real)
tailstats(signal::RDWaveform, start::RealQuantity, stop::RealQuantity)

Get statistics on the logarhithmic of the tail of a signal in the interval (start,stop).

tailstats(signal::AbstractSamples, start::Real, stop::Real)
tailstats(signal::RDWaveform, start::RealQuantity, stop::RealQuantity)

Get statistics on the logarhithmic of the tail of a signal in the interval (start,stop) by fitting a linear function.

weighted_linear_reg(w::AbstractVector, y::AbstractVector)

Do a weighted linear regression of the data y with weights w and return the offset and slope. its assumed that x is a range from 1 to length(y)