LegendSpecFits

Template for Julia packages.

BinnedPoissonProcess(intensity::AbstractMeasure, edges::AbstractVector{<:Real}) <: MeasureBase.AbstractMeasure

Binned view of a Poisson point process with the given intensity measure: its variates are vectors of non-negative integer bin counts for the given bin edges (logdensityof assumes valid counts, binned_poisson_likelihood validates them). The count in each bin is Poisson-distributed, with the intensity mass in the bin as expectation, computed via analytic MeasureBase.smf differences where available (all spectral shapes, distributed over weighted measures and superpositions) and via Simpson's rule otherwise. Distribution intensities are converted to measures automatically.

ExpDecay(μ, τ) <: SpectralShape

Exponential shape with density exp(-(x - μ) / τ), unit height at μ and decaying towards larger values. Infinite total mass over ℝ.

GaussPeak(μ, σ) <: SpectralShape

Gaussian peak shape with unit total mass.

GaussianStep(μ, σ) <: SpectralShape

Gaussian-smoothed step shape with unit plateau height below μ, erfc((x - μ) / (√2 σ)) / 2. Infinite total mass over ℝ.

HighETailPeak(μ, σ, θ) <: SpectralShape

Exponentially modified Gaussian peak shape with unit total mass and high-energy tail of scale θ, the mirror image of LowETailPeak. The shape underlying highEtail_peakshape.

LowETailPeak(μ, σ, θ) <: SpectralShape

Exponentially modified Gaussian peak shape with unit total mass and low-energy tail of scale θ (the LEGEND low-energy tail orientation). The shape underlying lowEtail_peakshape.

PolynomialShape(p0, p1, ...; x0 = 0) <: SpectralShape

Polynomial shape with density p0 + p1 (x - x0) + p2 (x - x0)^2 + .... The coefficients must keep the density non-negative on the fit region. An expansion point x0 near the region of interest improves the numerical conditioning of the coefficients.

abstract type SpectralShape <: BAT.BATMeasure

Supertype of the spectral shape measures.

_fit_fwhm_ft(e_grid::Matrix, e_grid_ft::StepRangeLen{Quantity{<:T}, Base.TwicePrecision{Quantity{<:T}}, Base.TwicePrecision{Quantity{<:T}}, Int}, rt::Unitful.RealOrRealQuantity, min_e::T, max_e::T, nbins::Int, rel_cut_fit::T; default_ft::Quantity{T}=3.0u"µs", ft_fwhm_tol::Unitful.Energy{<:Real} = 0.1u"keV") where {T <:Real}

Fit the FWHM values in e_grid for each FT in e_grid_ft with a Gamma Peakshape and return the optimal FT and the corresponding FWHM value. The cut values cut for each flat-top time a window for better histogramming.

Arguments

• e_grid: 2D array of energy values for each FT in e_grid_ft
• e_grid_ft: 1D array of FT values for which the FWHM values in e_grid are calculated
• rt: RT value for which the FWHM values in e_grid are calculated
• min_e: minimum energy value to consider for the fit
• max_e: maximum energy value to consider for the fit
• nbins: number of bins to use for the histogram of energy values
• rel_cut_fit: relative cut value to use for the fit
• ft_fwhm_tol: search for lowest "optimal" ft within minimum(fwhm) + ft_fwhm_tol to avoid artificially large ft

Returns

• ft: optimal FT value
• min_fwhm: corresponding FWHM value

_fit_fwhm_ft_ctc(e_grid::Matrix, e_grid_ft::StepRangeLen, qdrift::Vector{<:Real}, rt::Unitful.RealOrRealQuantity, min_e::T, max_e::T, nbins::Int, rel_cut_fit::T; default_ft::Quantity{T}=3.0u"µs", peak::Unitful.Energy{<:Real}=2614.5u"keV", window::Tuple{<:Unitful.Energy{<:Real}, <:Unitful.Energy{<:Real}}=(35.0u"keV", 25.0u"keV"), ft_fwhm_tol::Unitful.Energy{<:Real} = 0.1u"keV") where {T <:Real}

Fit the FWHM values in e_grid for each FT in e_grid_ft with a Gamma Peakshape and return the optimal FT and the corresponding FWHM value. The cut values cut for each flat-top time a window for better histogramming.

Arguments

• e_grid: 2D array of energy values for each FT in e_grid_ft
• e_grid_ft: 1D array of FT values for which the FWHM values in e_grid are calculated
• qdrift: drift time values for each energy value in e_grid
• rt: RT value for which the FWHM values in e_grid are calculated
• min_e: minimum energy value to consider for the fit
• max_e: maximum energy value to consider for the fit
• nbins: number of bins to use for the histogram of energy values
• rel_cut_fit: relative cut value to use for the fit
• ft_fwhm_tol: search for lowest "optimal" ft within minimum(fwhm) + ft_fwhm_tol to avoid artificially large ft

Returns

• ft: optimal FT value
• min_fwhm: corresponding FWHM value

_get_model_counts(f_fit::Base.Callable,v_ml::Union{NamedTuple, AbstractVector},bin_centers::StepRangeLen,bin_widths::StepRangeLen)

aux. function to get modelled peakshape based on histogram binning and best-fit parameter

_prepare_data(h::Histogram{<:Real,1})

aux. function to convert histogram data into bin edges, bin width and bin counts

aoe_compton_background_peakshape(
    x::Real, μ::Real, σ::Real,
    background::Real, δ::Real
)

Describes the background shape of a typical A/E Compton peak in a detector as a step like background for MSE events.

aoe_compton_peakshape(
    x::Real, μ::Real, σ::Real, n::Real,
    background::Real, δ::Real
)

Describes the shape of a typical A/E Compton peak in a detector as a gaussian SSE peak and a step like background for MSE events.

aoe_compton_peakshape_components(fit_func::Symbol; background_center::Real)

This function defines the components (signal, low/high-energy tail, backgrounds) of the fit function used in gamma specfits. These component functions are used in the fit-report and in plot receipes

aoe_compton_signal_peakshape(
    x::Real, μ::Real, σ::Real, n::Real
)

Describes the signal shape of a typical A/E Compton peak in a detector as a gaussian SSE peak.

array_to_tuple(a::AbstractArray, as_nt::NamedTuple)

Return a NamedTuple with the values of a in the order given by fieldnames(as_nt).

autocal_energy(e::AbstractArray{<:Real}, photon_lines::Vector{<:Unitful.RealOrRealQuantity}; min_e::Real=100, max_e_binning_quantile::Real=0.5, σ::Real = 2.0, threshold::Real = 50.0, min_n_peaks::Int = length(photon_lines), max_n_peaks::Int = 4 * length(photon_lines), α::Real = 0.01, rtol::Real = 5e-3)

Compute an energy calibration from raw reconstructed energy deposition values based on a given number of known photon lines which are contained in the spectrum

background_peakshape(
x::Real, μ::Real, σ::Real, 
step_amplitude::Real, background::Real; 
background_slope::Real = 0.0, background_exp = 0.0, background_center::Real = µ

)

Describes the background part of the shape of a typical gamma peak in a detector: components:

• step-function scaled with step_amplitude`
• energy-independent background: background`
• linear slope: background_slope (optional)
• exponential decay: background_exp (optional)

binned_poisson_likelihood(model, h::Histogram{<:Real,1})

Binned Poisson likelihood of the observed histogram h: the likelihood of BinnedPoissonProcess(model(p), h.edges) given the observed bin counts h.weights, model being a Markov kernel that maps parameters to an intensity measure (e.g. a @pf property function combining shape measures, see GaussPeak etc.).

A fixed intensity measure instead of model yields a parameter-independent likelihood.

The bin counts of h must be non-negative integers: a weighted histogram is not Poisson-distributed data and requires a different statistical treatment.

fit_chisq(x::AbstractVector{<:Real},y::AbstractVector{<:Real},yerr::AbstractVector{<:Real}, f_fit::Function;pull_t::Vector{<:NamedTuple} = fill(NamedTuple(), first(methods(f_fit)).nargs - 2), v_init::Vector = [])

Least square fit with chi2 minimization

Input:

• x : x-values
• y : y-values
• yerr : 1 sigma uncertainty on y
• ffit : fit/model function. e.g. for a linear function: flin(x,p1,p2) = p1 .* x .+ p2

The numer of fit parameter is determined with first(methods(f_fit)).nargs - 2. That's why it's important that ffit has the synthax f(x,arg1,arg2,arg3,...) pullt : pull term, a vector of NamedTuple with fields mean and std. A Gaussian pull term is added to the chi2 function to account for systematic uncertainties. If left blank, no pull term is used. vinit : initial value for fit parameter optimization. If left blank, the initial value is set to 1 or guessed roughly for all fit parameters correlated : if true, the returned parameters carry their full covariance, so that uncertainties of quantities derived from them, including `report.ffit, account for the parameter correlations. Iffalse`, the parameters are independent measurements with the same marginal uncertainties, which overestimates the uncertainty of derived quantities whenever the parameters are correlated.

Return:

• result : NamedTuple with the optimized fit parameter and the fit function
• report:

ctc_energy(e::Array{T}, qdrift::Array{T}, peak::T, window::T) where T<:Real

Correct for the drift time dependence of the energy by minimizing the ratio of the FWHM and the peak height of the peak around peak in e with a cut window of window. The drift time dependence is given by qdrift.

Returns

* `peak`: peak position
* `window`: window size
* `fct`: correction factor
* `fwhm_before`: FWHM before correction
* `fwhm_after`: FWHM after correction
* `func`: function to correct energy
* `func_generic`: generic function to correct energy

cut_single_peak(x::Array, min_x::Float64, max_x::Float64,; n_bins::Int=-1, relative_cut::Float64=0.5)

Cut out a single peak from the array x between min_x and max_x. The number of bins is the number of bins to use for the histogram. The relative cut is the fraction of the maximum counts to use for the cut.

Returns

* `max`: maximum position of the peak
* `low`: lower edge of the cut peak
* `high`: upper edge of the cut peak

double_gaussian(
    x::Real, μ1::Real, σ1::Real, n1::Real, 
    μ2::Real, σ2::Real, n2::Real
)

Evaluates the sum of two gaussians at x with parameters μ1, σ1, n1, μ2, σ2, n2.

estimate_combined_peak_stats(peakstats::StructArray,; calib_type::Symbol=:th228)

Estimate the peak position, FWHM, sigma, counts and background of a peak from a histogram.

estimate_fwhm(v::NamedTuple, v_err::NamedTuple)

Get the FWHM of a peak from the fit parameters.

Returns

* `fwhm`: the FWHM of the peak

estimate_fwhm_aoe_compton(v::NamedTuple, f_fit::Function)

Get the FWHM of a peak from the fit parameters.

Returns

* `fwhm`: the FWHM of the peak

estimate_single_peak_stats(h::Histogram, calib_type::Symbol=:th228)

Estimate statistics/parameters for a single peak in the given histogram h.

h must only contain a single peak. The peak should have a Gaussian-like shape. calib_type specifies the calibration type. Currently :th228, :psd and :simple is implemented..

Returns

NamedTuple with the fields * peak_pos: estimated position of the peak (in the middle of the peak) * peak_fwhm: full width at half maximum (FWHM) of the peak * peak_sigma: estimated standard deviation of the peak * peak_counts: estimated number of counts in the peak * mean_background: estimated mean background value

ex_gauss_pdf(x::Real, μ::Real, σ::Real, θ::Real)

The PDF of an Exponentially modified Gaussian distribution with Gaussian parameters μ, σ and exponential scale θ at x.

It is the PDF of the distribution that descibes the random process rand(Normal(μ, σ)) + rand(Exponential(θ)).

ex_step_gauss(x::Real, l::Real, k::Real, t::Real, d::Real)

Evaluates an extended step gauss model at x with parameters l, k, t and d.

expand_vars(v::NamedTuple)::StructArray

Expand all fields in v (scalars or arrays) to same array size and return a StructArray.

exponential_decay(x::Real, amplitude::Real, decay::Real, offset::Real)

Evaluates an exponential decay function at x with parameters amplitude, decay and offset.

fit_aoe_compton(peakhists::Array, peakstats::StructArray, compton_bands::Array{T}) where T<:Real

Fit the A/E Compton bands using the f_aoe_compton function consisting of a gaussian SSE peak and a step like background for MSE events.

Returns

* `result`: Dict of NamedTuples of the fit results containing values and errors for each compton band
* `report`: Dict of NamedTuples of the fit report which can be plotted for each compton band

fit_aoe_compton_combined(peakhists::Vector{<:Histogram}, peakstats::StructArray, compton_bands::Array{T}, result_corrections::NamedTuple; pars_aoe::NamedTuple{(:μ, :μ_err, :σ, :σ_err)}=NamedTuple{(:μ, :μ_err, :σ, :σ_err)}(nothing, nothing, nothing, nothing), uncertainty::Bool=false) where T<:Unitful.Energy{<:Real}

Performed a combined fit over all A/E Compton band using the f_aoe_compton function consisting of a gaussian SSE peak and a step like background for MSE events, assuming f_aoe_mu for μ and f_aoe_sigma for σ.

Returns

* `v_ml`: The fit result from the maximum-likelihood fit.
* `report_μ`: Report to plot the combined fit result for the enery-dependence of `μ`.
* `report_σ`: Report to plot the combined fit result for the enery-dependence of `σ`.
* `report_bands`: Dict of NamedTuples of the fit report which can be plotted for each compton band

fitaoecorrections(e::Array{<:Unitful.Energy{<:Real}}, μ::Array{<:Real}, σ::Array{<:Real})

Fit the corrections for the AoE value of the detector.

Returns

• e: Energy values
• μ: Mean values
• σ: Sigma values
• μ_scs: Fit result for the mean values
• f_μ_scs: Fit function for the mean values
• σ_scs: Fit result for the sigma values
• f_σ_scs: Fit function for the sigma values

fit_binned_double_gauss(h::Histogram, ps::NamedTuple{(:peak_pos, :peak_fwhm, :peak_sigma, :peak_counts, :mean_background, :μ, :σ), NTuple{7, T}}; uncertainty::Bool=true) where T<:Real

Perform a binned fit of the peakshape to the data in h using the initial values in ps while using the f_double_gauss function consisting of a double gaussian peak. The priors for the first gaussian peak are given by the ps tuple. For the priors of the second gaussian peak a wide window around the first peak is used.

Returns

* `result`: NamedTuple of the fit results containing values and errors
* `report`: NamedTuple of the fit report which can be plotted

fit_binned_half_trunc_gauss(x::Array, cuts::NamedTuple{(:low, :high, :max), Tuple{Float64, Float64, Float64}})

Fit a single truncated Gaussian to the data x between cut.low and cut.high. The peak center is fixed at μ and the peak is cut in half either in the left or right half.

Returns report and result with:

* `f_fit`: fitted function
* `μ`: mean of the Gaussian
* `σ`: standard deviation of the Gaussian

fit_binned_trunc_gauss(h::Histogram, ps::NamedTuple{(:peak_pos, :peak_fwhm, :peak_sigma, :peak_counts, :mean_background, :μ, :σ), NTuple{7, T}}; uncertainty::Bool=true) where T<:Real

Perform a binned fit of the peakshape to the data in h using the initial values in ps while using the f_gauss function consisting of a gaussian peak multiplied with an amplitude n.

Returns

* `result`: NamedTuple of the fit results containing values and errors
* `report`: NamedTuple of the fit report which can be plotted

fit_calibration(pol_order::Int, µ::AbstractVector{<:Union{Real,Measurement{<:Real}}}, peaks::AbstractVector{<:Union{Real,Measurement{<:Real}}}; pull_t::Vector{<:NamedTuple}=fill(NamedTuple(), pol_order+1), v_init::Vector = [], uncertainty::Bool=true )

Fit the calibration lines with polynomial function of polorder order polorder == 1 -> linear function pol_order == 2 -> quadratic function

Returns

* `result`: NamedTuple with the following fields
    * `par`: best-fit parameters
    * `gof`: godness of fit
    * `func`: the calibration curve as a LEGEND Julia expression of `e_expression`
    * `func_err`: the same expression with the uncertainty of the curve written out as `value ± error`,
      which includes the correlations between the fit parameters (see `chi2fit`)
* `report`:

fit_enc_sigmas(enc_grid::Matrix{T}, enc_grid_rt::StepRangeLen{Quantity{<:T}, Base.TwicePrecision{Quantity{<:T}}, Base.TwicePrecision{Quantity{<:T}}, Int}, min_enc::T, max_enc::T, nbins::Int, rel_cut_fit::T) where T<:Real

Fit the ENC values in enc_grid for each RT in enc_grid_rt with a Gaussian and return the optimal RT and the corresponding ENC value.

Arguments

• enc_grid: 2D array of ENC values for each RT in enc_grid_rt
• enc_grid_rt: 1D array of RT values for which the ENC values in enc_grid are calculated
• min_enc: minimum ENC value to consider for the fit
• max_enc: maximum ENC value to consider for the fit
• nbins: number of bins to use for the histogram of ENC values
• rel_cut_fit: relative cut value to use for the fit

Returns

• rt: optimal RT value
• min_enc: corresponding ENC value

fit_fwhm(pol_order::Int, peaks::Vector{<:Unitful.Energy}, fwhm::Vector{<:Unitful.Energy}; e_type_cal, e_expression, uncertainty, correlated, pull_fano)

Fit the resolution curve fwhm(E) = √(enc + fano·E [+ ct·E²]) to the FWHM of the peaks, with pol_order the degree of the polynomial under the square root. With pull_fano = true, the quadratic fit ties the linear term to Fano statistics of germanium (F = 0.112 ± 0.015) by a pull term, so that the E² term takes up the other contributions; by default, and in the linear fit, the linear term is free.

Returns

* `par`: `enc`, the width √enc in keV that the curve approaches at zero energy; `fano`, the Fano
  factor extracted from the linear term with the pair creation energy of germanium; and, for the
  quadratic fit, the dimensionless coefficient `ct` of E²
* `qbb`: the FWHM at 2039 keV
* `func`, `func_cal`: the resolution curve as LEGEND Julia expressions of `e_expression` and `e_type_cal`
* `func_err`, `func_cal_err`: the same expressions with the uncertainty of the curve written out as
  `value ± error`, which includes the correlations between the fit parameters (see `chi2fit`).
  Evaluating them at an energy that is itself a `Measurement` nests the two uncertainties.
* `f_fit` (report): the fitted function

fit_fwhm_ft(e_grid::Matrix, e_grid_ft::StepRangeLen, rt::Unitful.RealOrRealQuantity, min_e::T, max_e::T, rel_cut_fit::T, apply_ctc::Bool=true; kwargs...)
fit_fwhm_ft(e_grid::Matrix, e_grid_ft::StepRangeLen, rt::Unitful.RealOrRealQuantity, min_e, max_e, rel_cut_fit; kwargs...)
fit_fwhm_ft(e_grid::Matrix, e_grid_ft::StepRangeLen, qdrift::Vector{<:Real}, rt::Unitful.RealOrRealQuantity, min_e, max_e, rel_cut_fit; kwargs...)

Fit the FWHM values in e_grid for each FT in e_grid_ft with a Gamma Peakshape and return the optimal FT and the corresponding FWHM value. The cut values cut for each flat-top time a window for better histogramming. If the apply_ctc flag is set to true, the CTC correction is applied to the energy values. Othwise, if a qdrift vector is provided, the CTC correction is applied to the energy values.

Arguments

• e_grid: 2D array of energy values for each FT in e_grid_ft
• e_grid_ft: 1D array of FT values for which the FWHM values in e_grid are calculated
• rt: RT value for which the FWHM values in e_grid are calculated
• min_e: minimum energy value to consider for the fit
• max_e: maximum energy value to consider for the fit
• rel_cut_fit: relative cut value to use for the fit

Returns

• ft: optimal FT value
• min_fwhm: corresponding FWHM value

fit_half_centered_trunc_gauss(x::Array, cuts::NamedTuple{(:low, :high, :max), Tuple{Float64, Float64, Float64}})

Fit a single truncated Gaussian to the data x between cut.low and cut.high. The peak center is fixed at μ and the peak is cut in half either in the left or right half.

Returns report and result` with:

* `f_fit`: fitted function
* `μ`: mean of the Gaussian
* `σ`: standard deviation of the Gaussian

fit_half_trunc_gauss(x::Array, cuts::NamedTuple{(:low, :high, :max), Tuple{Float64, Float64, Float64}})

Fit a single truncated Gaussian to the data x between cut.low and cut.high. The peak center is fixed at μ and the peak is cut in half either in the left or right half.

Returns report and result with:

* `f_fit`: fitted function
* `μ`: mean of the Gaussian
* `σ`: standard deviation of the Gaussian

fit_lq_compton(peakhists::Array, peakstats::StructArray, compton_bands::Array{T}) where T<:Real

Fit the A/E Compton bands using the f_lq_compton function consisting of a gaussian SSE peak and a step like background for MSE events.

Returns

* `result`: Dict of NamedTuples of the fit results containing values and errors for each compton band
* `report`: Dict of NamedTuples of the fit report which can be plotted for each compton band

fit_peaks(peakhists::Array, peakstats::StructArray, th228_lines::Array,; calib_type::Symbol=:th228, uncertainty::Bool=true, low_e_tail::Bool=true)

Perform a fit of the peakshape to the data in peakhists using the initial values in peakstats to the calibration lines in th228_lines.

Returns

* `peak_fit_plots`: array of plots of the peak fits
* `return_vals`: dictionary of the fit results

fit_peaks_combined(peakhists::Array, peakstats::StructArray, th228_lines::Array{T},; calib_type::Symbol=:th228, uncertainty::Bool=true, fixed_position::Bool=false) where T<:Real

Fit the peaks of a histogram to a combined peakshape function while sharring parameters between peaks.

fit_sf_wl(dep_sep_data, a_grid_wl_sg, optimization_config)

Fit a A/E filter window length for the SEP data and return the optimal window length and the corresponding survival fraction.

Arguments

• dep_sep_data: NamedTuple with the DEP and SEP data
• a_grid_wl_sg: range of window lengths to sweep through
• optimization_config: configuration dictionary

Returns

• result: optimal window length and corresponding survival fraction
• report: report with all window lengths and survival fractions

fit_single_aoe_compton(h::Histogram, ps::NamedTuple{(:peak_pos, :peak_fwhm, :peak_sigma, :peak_counts, :mean_background, :μ, :σ), NTuple{7, T}}; uncertainty::Bool=true) where T<:Real

Perform a fit of the peakshape to the data in h using the initial values in ps while using the f_aoe_compton function consisting of a gaussian SSE peak and a step like background for MSE events.

Returns

* `result`: NamedTuple of the fit results containing values and errors
* `report`: NamedTuple of the fit report which can be plotted

fit_single_aoe_compton_with_fixed_μ_and_σ(h::Histogram, μ::Number, σ::Number, ps::NamedTuple; uncertainty::Bool=true)

Fit a single A/E Compton band using the f_aoe_compton function consisting of a gaussian SSE peak and a step like background for MSE events using fixed values for μ and σ.

Returns

* `neg_log_likelihood`: The negative log-likelihood of the likelihood fit
* `report`: Dict of NamedTuples of the fit report which can be plotted for each compton band

fit_single_lq_compton(h::Histogram, ps::NamedTuple{(:peak_pos, :peak_fwhm, :peak_sigma, :peak_counts, :mean_background, :μ, :σ), NTuple{7, T}}; uncertainty::Bool=true) where T<:Real

Perform a fit of the peakshape to the data in h using the initial values in ps while using the f_lq_compton function consisting of a gaussian SSE peak and a step like background for MSE events.

Returns

* `result`: NamedTuple of the fit results containing values and errors
* `report`: NamedTuple of the fit report which can be plotted

fit_single_peak_th228(h::Histogram, ps::NamedTuple{(:peak_pos, :peak_fwhm, :peak_sigma, :peak_counts, :mean_background), NTuple{5, T}};, uncertainty::Bool=true, fixed_position::Bool=false, low_e_tail::Bool=true) where T<:Real

Perform a fit of the peakshape to the data in h using the initial values in ps while using the gamma_peakshape with low-E tail. Also, FWHM is calculated from the fitted peakshape with MC error propagation. The peak position can be fixed to the value in ps by setting fixed_position=true. If the low-E tail should not be fitted, it can be disabled by setting low_e_tail=false.

Returns

* `result`: NamedTuple of the fit results containing values and errors
* `report`: NamedTuple of the fit report which can be plotted

fit_single_trunc_gauss(x::Array, cuts::NamedTuple{(:low, :high, :max), Tuple{Float64, Float64, Float64}})

Fit a single truncated Gaussian to the data x between min_x and max_x. Returns report and resultwith: *ffit: fitted function *μ: mean of the Gaussian *μerr: error of the mean *σ: standard deviation of the Gaussian *σ_err: error of the standard deviation *n`: number of counts in the peak

fit_sipm_spectrum(pe_cal::Vector{<:Real}, min_pe::Real=0.5, max_pe::Real=3.5; 
    n_mixtures::Int=ceil(Int, (max_pe - min_pe) * 4), nIter::Int=50, nInit::Int=50, 
    method::Symbol=:kmeans, kind=:diag, Δpe_peak_assignment::Real=0.3, f_uncal::Function=identity, uncertainty::Bool=true)

Fit a Gaussian Mixture Model to the given pe calibration data and return the fit parameters.

Arguments

• pe_cal::Vector{<:Real}: the pe calibration data
• min_pe::Real=0.5: the minimum pe to consider
• max_pe::Real=3.5: the maximum pe to consider
• n_mixtures::Int=ceil(Int, (max_pe - min_pe) * 4): the number of mixtures to fit
• nIter::Int=50: the number of iterations for the EM algorithm
• nInit::Int=50: the number of initializations for the EM algorithm
• method::Symbol=:kmeans: the method to use for initialization
• kind::Symbol=:diag: the kind of covariance matrix to use
• Δpe_peak_assignment::Real=0.3: the range to consider for peak assignment
• f_uncal::Function=identity: the function to use for uncalibration
• uncertainty::Bool=true: whether to calculate the uncertainty

Returns

• result: a tuple with the fit parameters
• report: a tuple with the fit report which can be plotted via a recipe

fit_sipm_threshold(thresholds::Vector{<:Real}, min_cut::Real=minimum(thresholds), max_cut::Real=maximum(thresholds); n_bins::Int=-1, relative_cut::Real=0.2, fit_thresholds::Bool=true, uncertainty::Bool=true)

Fit the SiPM threshold spectrum and return the optimal threshold.

Arguments

• thresholds: vector of thresholds
• min_cut: minimum threshold
• max_cut: maximum threshold
• n_bins: number of bins for histogram
• relative_cut: relative cut for threshold
• fit_thresholds: fit thresholds
• uncertainty: calculate uncertainty

Returns

• result: optimal threshold and corresponding gain, resolution and position of 1pe peak
• report: report with all thresholds and corresponding gains, resolutions and positions of 1pe peaks

fit_sipm_wl(trig_max_grid::VectorOfVectors{<:Real}, e_grid_wl::StepRangeLen)

Fit the SiPM spectrum for different window lengths and return the optimal window length.

Arguments

• trig_max_grid: grid of trigger maxima for different window lengths
• e_grid_wl: range of window lengths to sweep through

Returns

• result: optimal window length and corresponding gain, resolution and position of 1pe peak
• report: report with all window lengths and corresponding gains, resolutions and positions of 1pe peaks

fit_single_peak_th228(h::Histogram, ps::NamedTuple{(:peak_pos, :peak_fwhm, :peak_sigma, :peak_counts, :mean_background), NTuple{5, T}};, uncertainty::Bool=true, fixed_position::Bool=false, low_e_tail::Bool=true) where T<:Real

Perform a simultaneous fit of two peaks (h_survived and h_cut) that together would form a histogram h, from which the result h_result was already determined using fit_single_peak_th228. Also, FWHM is calculated from the fitted peakshape with MC error propagation. The peak position can be fixed to the value in ps by setting fixed_position=true. If the low-E tail should not be fitted, it can be disabled by setting low_e_tail=false.

Returns

* `result`: NamedTuple of the fit results containing values and errors, in particular the signal survival fraction `sf` and the background survival frachtion `bsf`.
* `report`: NamedTuple of the fit report which can be plotted

gamma_peakshape(
    x::Real, μ::Real, σ::Real, n::Real,
    step_amplitude::Real, skew_fraction::Real, skew_width::Real,
    background::Real;  
    skew_fraction_highE::Real = 0.0, skew_width_highE::Real= 0.0, 
    background_kwargs... 
)

Standard gamma peakshape: Describes the shape of a typical gamma peak in a detector. Components:

• Gaussian signal peak with μ, σ, n - skew_fraction - skew_fraction_highE
• low-energy tail: skew_fraction, skew_width
• high-energy tail: skew_fraction_highE, skew_width_highE (optional, default off)
• background:
  • energy-independent background
  • step-function scaled with step_amplitude from Compton scattered gammas
  • linear slope: background_slope (optional, default off)
  • exponential decay: background_exp (optional, default off)

LegendSpecFits.gauss_pdf(x::Real, μ::Real, σ::Real)

Equivalent to pdf(Normal(μ, σ), x)

generate_aoe_compton_bands(aoe::Vector{<:Real}, e::Vector{<:T}, compton_bands::Vector{<:T}, compton_window::T) where T<:Unitful.Energy{<:Real}

Generate histograms for the A/E Compton bands and estimate peak parameters. The compton bands are cutted out of the A/E spectrum and then binned using the Freedman-Diaconis Rule. For better performance the binning is only done in the area around the peak. The peak parameters are estimated using the estimate_single_peak_stats_psd function.

Returns

* `peakhists`: Array of histograms for each compton band
* `peakstats`: StructArray of peak parameters for each compton band
* `min_aoe`: Array of minimum A/E values for each compton band
* `max_aoe`: Array of maximum A/E values for each compton band
* `mean_peak_pos`: Mean peak position of all compton bands
* `std_peak_pos`: Standard deviation of the peak position of all compton bands
* `simple_pars_aoe_μ`: Simple curve fit parameters for the peak position energy depencence
* `simple_pars_error_aoe_μ`: Simple curve fit parameter errors for the peak position energy depencence
* `simple_pars_aoe_σ`: Simple curve fit parameters for the peak sigma energy depencence
* `simple_pars_error_aoe_σ`: Simple curve fit parameter errors for the peak sigma energy depencence

get_centered_gaussian_window_cut(x::Array, min_x::Float64, max_x::Float64, n_σ::Real, center::Float64=0.0, n_bins::Int=500, relative_cut::Float64=0.2, left::Bool=false)

Cut out a single peak from the array x between min_x and max_x by fitting a truncated one-sided Gaussian and extrapolating a window cut with n_σ standard deviations. The center and side of the fit can be specified with left and center variable.

Returns

* `low_cut`: lower edge of the cut peak
* `high_cut`: upper edge of the cut peak
* `µ`: center of the peak
* `σ`: standard deviation of the Gaussian
* `gof`: goodness of fit
* `low_cut_fit`: lower edge of the cut peak from the fit
* `high_cut_fit`: upper edge of the cut peak from the fit
* `max_cut_fit`: maximum position of the peak

get_continuum_survival_fraction(e, center, window, survival_flag; low_flag=survival_flag)
get_continuum_survival_fraction(cut_parameter, e, center, window; low_cut=-Inf, high_cut=Inf, selection=trues(length(cut_parameter)))

Count events in the continuum window before and after a selection. The interval overload forwards a Boolean flag to the counting method. The report retains the existing h_after_low and h_after_ds fields: the former applies only the lower bound, while the latter always represents the final selection, including high_cut and selection.

get_distribution_transform(d::Distribution, pprior::Prior)

Return a DistributionTransform for the given Distribution and Prior.

get_friedman_diaconis_bin_width(x::AbstractArray)

Return the bin width for the given data x using the Friedman-Diaconis rule.

get_low_aoe_cut(aoe::Vector{<:Unitful.RealOrRealQuantity}, e::Vector{<:T},; 
        dep::T=1592.53u"keV", window::Vector{<:T}=[12.0, 10.0]u"keV", dep_sf::Float64=0.9, rtol::Float64=0.001, maxiters::Int=300, sigma_high_sided::Float64=Inf,
        cut_search_interval::Tuple{<:Unitful.RealOrRealQuantity, <:Unitful.RealOrRealQuantity}=(-25.0*unit(first(aoe)), 1.0*unit(first(aoe))), 
        bin_width_window::T=3.0u"keV", max_e_plot::T=3000.0u"keV",  plot_window::Vector{<:T}=[12.0, 50.0]u"keV",
        fixed_position::Bool=true, fit_func::Symbol=:gamma_def, uncertainty::Bool=true) where T<:Unitful.Energy{<:Real}

Get the AoE cut value for a given dep and window size while performing a peak fit with fixed position. The AoE cut value is determined by finding the cut value for which the number of counts after the cut is equal to dep_sf times the number of counts before the cut. The algorhithm utilizes a root search algorithm to find the cut value with a relative tolerance of rtol.

Returns

• cut: AoE cut value
• n0: Number of counts before the cut
• nsf: Number of counts after the cut

get_mc_value_shapes(v::NamedTuple, v_err::NamedTuple, n::Integer)

Return a NamedTuple with the same fields as v and v_err but with Normal distributions for each field.

get_mc_value_shapes(v::NamedTuple, v_err::Matrix, n::Integer)

Generate n random samples of fit parameters using their respective best-fit values v and covariance matrix v_err

get_number_of_bins(x::AbstractArray,; method::Symbol=:sqrt)

Return the number of bins for the given data x using the given method.

get_peak_fwhm_aoe_compton(v_ml::NamedTuple, v_ml_err::NamedTuple)

Get the FWHM of a peak from the fit parameters while performing a MC error propagation.

Returns

* `fwhm`: the FWHM of the peak
* `fwhm_err`: the uncertainty of the FWHM of the peak

get_peak_fwhm_th228(v_ml::NamedTuple, v_ml_err::NamedTuple)

Get the FWHM of a peak from the fit parameters while performing a MC error propagation.

Returns

* `fwhm`: the FWHM of the peak
* `fwhm_err`: the uncertainty of the FWHM of the peak

get_peak_survival_fraction(h_before::Histogram{<:Real, 1}, h_survived::Histogram{<:Real, 1}, h_cut::Histogram{<:Real, 1}; uncertainty::Bool=true, fit_func::Symbol=:gamma_def)
get_peak_survival_fraction(e::AbstractVector{<:T}, survival_flag::AbstractVector{Bool}; uncertainty::Bool=true, fit_func::Symbol=:gamma_def) where T<:Unitful.Energy{<:Real}
get_peak_survival_fraction(e::AbstractVector{<:T}, peak::T, window::Tuple{T, T}, survival_flag::AbstractVector{Bool}; uncertainty::Bool=true, bin_width_window::T=2.0u"keV", fit_func::Symbol=:gamma_def) where T<:Unitful.Energy{<:Real}
get_peak_survival_fraction(cut_parameter::AbstractVector{<:Unitful.RealOrRealQuantity}, e::AbstractVector{<:T}, peak::T, window::Tuple{T, T}; low_cut::Unitful.RealOrRealQuantity=-Inf*unit(first(cut_parameter)), high_cut::Unitful.RealOrRealQuantity=Inf*unit(first(cut_parameter)), selection::AbstractVector{Bool}=trues(length(cut_parameter)), kwargs...) where T<:Unitful.Energy{<:Real}

Fit a peak before and after an arbitrary event selection. The histogram method accepts one-dimensional count histograms for all, surviving, and rejected events with the same unitless energy binning; the latter two must partition the first. It returns fit results and reports without a nominal peak label.

The two-vector method bins the supplied peak-region energies and returns the histogram method's outputs. The center/window method accepts a tuple of lower and upper energy widths and adds peak to both outputs. The cut-parameter method constructs the selection and forwards to the center/window method. Cut bounds are exclusive.

get_peaks_survival_fractions(e::AbstractVector{<:T}, peaks::AbstractVector{<:T}, peak_names::AbstractVector{Symbol}, windows::AbstractVector{<:Tuple{T, T}}, survival_flag::AbstractVector{Bool}; kwargs...) where T<:Unitful.Energy{<:Real}
get_peaks_survival_fractions(cut_parameter::AbstractVector{<:Unitful.RealOrRealQuantity}, e::AbstractVector{<:T}, peaks::AbstractVector{<:T}, peak_names::AbstractVector{Symbol}, windows::AbstractVector{<:Tuple{T, T}}; kwargs...) where T<:Unitful.Energy{<:Real}

Fit survival fractions for several peaks. Each window is a tuple of lower and upper energy widths. The overloads accepting separate vectors of widths construct those tuples and forward. The primary method accepts the final event selection directly. The interval-based method constructs that flag and forwards to the primary method. Either bound may be omitted; with both omitted, finite cut parameter values survive. Bounds are exclusive.

residuals(f_fit::Base.Callable, h::Histogram{<:Real,1},v_ml::Union{NamedTuple, AbstractVector})

Calculate bin-wise residuals and normalized residuals. Calcualte bin-wise p-value based on poisson distribution for each bin.

Input:

• f_fitfunction handle of fit function (peakshape)
• h histogram of data
• v_ml best-fit parameters

Returns:

• residuals difference: model - data (histogram bin count)
• residuals_norm normalized residuals: model - data / sqrt(model)
• p_value_binwise p-value for each bin based on poisson distribution
• bin_centers centers of the bins for which the residuals were determined

get_sf_after_aoe_cut(aoe_cut::Unitful.RealOrRealQuantity, aoe::Vector{<:Unitful.RealOrRealQuantity}, e::Vector{<:T}, peak::T, window::Vector{T}, bin_width::T, result_before::NamedTuple; uncertainty::Bool=true, fit_func::Symbol=:gamma_def) where T<:Unitful.Energy{<:Real}

Get the survival fraction after a AoE cut value aoe_cut for a given peak and window size from a combined fit to the survived and cut histograms.

Returns

• sf: Survival fraction after the cut

get_th228_fit_functions(; background_center::Union{Real,Nothing} = nothing)

This function defines the gamma peakshape fit functions used in the calibration specfits.

• gamma_def: "default" gamma peakshape with gaussian signal, low-energy tail, and background (flat + step)
• gamma_tails: default gamma peakshape + high-energy tail
• gamma_bckSlope: default gamma peakshape + linear background slope
• gamma_bckExp: default gamma peakshape + exponential background
• gamma_bckFlat: default gamma peakshape - step background (only flat component!)
• gammatailsbckFlat: default gamma peakshape + high-energy tail - step background (only flat component!)
• gamma_minimal: only Gaussian signal and flat background

highEtail_peakshape(
    x::Real, μ::Real, σ::Real, n::Real,
    skew_fraction::Real, skew_width::Real,
)

Describes the high-E signal tail part of the shape

hist_loglike(f_fit::Base.Callable, h::Histogram{<:Real,1})

Calculate the Poisson log-likelihood of a fit function f_fit(x) and a histogram h. f_fit must accept all values x on the horizontal axis of the histogram.

Currently uses a simple midpoint-rule integration of f_fit over the bins of h.

linear_function(x::Real, slope::Real, intercept::Real)

Evaluates a linear function at x with parameters slope and intercept.

lowEtail_peakshape(
    x::Real, μ::Real, σ::Real, n::Real,
    skew_fraction::Real, skew_width::Real,
)

Describes the low-E signal tail part of the shape of a typical gamma peak in a detector.

lq_ctc_correction(lq::Vector{<:AbstractFloat}, dt_eff::Vector{<:Unitful.RealOrRealQuantity}, e_cal::Vector{<:Unitful.Energy{<:Real}}, dep_µ::Unitful.AbstractQuantity, dep_σ::Unitful.AbstractQuantity; 
ctc_dep_edgesigma::Float64=3.0, ctc_lq_precut_relative_cut::Float64=0.25, lq_outlier_sigma::Float64 = 2.0, ctc_driftime_cutoff_method::Symbol=:percentile, dt_eff_outlier_sigma::Float64 = 2.0, lq_e_corr_expression::Union{String,Symbol}="lq / e", dt_eff_expression::Union{String,Symbol}="qdrift / e", ctc_dt_eff_low_quantile::Float64=0.15, ctc_dt_eff_high_quantile::Float64=0.95, pol_fit_order::Int=1, uncertainty::Bool=false)

Perform the drift time correction on the LQ data using the DEP peak. The function cuts outliers in lq and drift time, then performs a polynomial fit on the remaining data. The data is Corrected by subtracting the polynomial fit from the lq data.

Arguments

* `lq`: Energy corrected lq parameter
* `dt_eff`: Effective drift time
* `e_cal`: Energy
* `dep_µ`: Mean of the DEP peak
* `dep_σ`: Standard deviation of the DEP peak

Keywords

* `ctc_dep_edgesigma`: Number of standard deviations used to define the DEP edges
* `ctc_lq_precut_relative_cut`: Relative cut for cut_single_peak function
* `ctc_driftime_cutoff_method`: Method used to define the drift time cutoff
* `lq_outlier_sigma`: Number of standard deviations used to define the lq cutoff
* `dt_eff_outlier_sigma`: Number of standard deviations used to define the drift time cutoff
* `lq_e_corr_expression`: Expression for the energy corrected lq classifier 
* `dt_eff_expression`: Expression for the effective drift time 
* `ctc_dt_eff_low_quantile`: Lower quantile used to define the drift time cutoff
* `ctc_dt_eff_high_quantile`: Higher quantile used to define the drift time cutoff
* `pol_fit_order`: Order of the polynomial fit used for the drift time correction

Returns

* `result`: NamedTuple of the function used for the drift time correction, the polynomial fit result and the box constraints
* `report`: NamedTuple of the histograms used for the fit, the cutoff values and the DEP edges

lq_norm(dep_µ::Unitful.Energy, dep_σ::Unitful.Energy, e_cal::Vector{<:Unitful.Energy}, lq_classifier::Vector{<:AbstractFloat}; dep_sideband_sigma::Float64=4.5, cut_truncation_sigma::Float64=2.0,  uncertainty::Bool=true, lq_class_expression::Union{String,Symbol}="lq / e  - (slope * qdrift / e + y_inter)" )

Performs normalization of the charge-trapping-corrected LQ classifier using the double-escape peak (DEP) region of the Th calibration. It subtracts sidebands around the DEP, fits a truncated Gaussian to the resulting histogram, and generates a normalization expression so that a numerical value of one corresponds to one standard deviation (σ) of the fitted Gaussian.

Arguments

* `dep_µ`: Mean of the DEP peak
* `dep_σ`: Standard deviation of the DEP peak
* `e_cal`: Vector of Energy values
* `lq_classifier`: LQ classifier (typically charge-trapping-corrected)

Keywords

* `dep_sideband_sigma`: Number of standard deviations used to define the sideband edges
* `cut_truncation_sigma`: Number of standard deviations used for the precut of sideband subtracted histogram
* `uncertainty`: Boolean flag to include uncertainty in the fit (default: true)
* `lq_class_expression`: Expression for the used LQ classifier

Returns

* `result`: NamedTuple of the fit result and normalization function
* `report`: NamedTuple of the fit result, fit report and temporary histograms

nearestSPD(A::Matrix{<:Real})

Returns the nearest positive definite matrix to A Calculation is based on matrix factorization techniques described in https://www.sciencedirect.com/science/article/pii/0024379588902236

p_value(f_fit::Base.Callable, h::Histogram{<:Real,1},v_ml::Union{NamedTuple, AbstractVector})

calculate p-value based on least-squares, assuming gaussian uncertainty baseline method to get goodness-of-fit (gof)

input:

• f_fitfunction handle of fit function (peakshape)
• h histogram of data
• v_ml best-fit parameters

returns:

• pval p-value of chi2 test
• chi2 chi2 value
• dof degrees of freedom

p_value_MC(f_fit::Base.Callable, h::Histogram{<:Real,1},ps::NamedTuple{(:peak_pos, :peak_fwhm, :peak_sigma, :peak_counts, :mean_background)},v_ml::NamedTuple,;n_samples::Int64=1000)

alternative p-value calculation via Monte Carlo sampling. Warning: computational more expensive than pvaule() and pvalue_LogLikeRatio()

Input:

• f_fitfunction handle of fit function (peakshape)
• h histogram of data
• ps best-fit parameters
• v_ml best-fit parameters
• n_samples number of samples

Performed Steps:

• Create n_samples randomized histograms. For each bin, samples are drawn from a Poisson distribution with λ = model peak shape (best-fit parameter)
• Each sample histogram is fit using the model function f_fit
• For each sample fit, the max. loglikelihood fit is calculated

Returns

• % p value –> comparison of sample max. loglikelihood and max. loglikelihood of best-fit

p_value_poissonll(f_fit::Base.Callable, h::Histogram{<:Real,1},v_ml::Union{NamedTuple, AbstractVector})

p-value via poisson likelihood ratio: baseline for ML fits using Poisson statistics and bins with low number of counts

peak_centroid(v::NamedTuple)

calculate centroid of gamma peak from fit parameters

 peak_search_gamma(e_uncal::Vector{<:Real}, gamma_lines::Vector{<:Unitful.Energy{<:Real}}; peakfinder_σ::Real = 2.0, peakfinder_threshold::Real = 10.0, peak_quantile::ClosedInterval{<:Real} = 0.0..1.0, bin_quantile::ClosedInterval{<:Real} = peak_quantile, quantile_perc::Float64=NaN)

peak search in gamma-ray spectrum

inputs:

• e_uncal: uncalibrated energy array
• gamma_lines: array of gamma lines
• keyword arguments: * peakfinder_σ::Real=2.0: The expected sigma of a peak in the spectrum. In units of bins. (see RadiationSpectra.peakfinder) * peakfinder_threshold::Real=10.0: Threshold for being identified as a peak in the deconvoluted spectrum. A single bin is identified as an peak when its weight exceeds the threshold and the previous bin was not identified as an peak. (see RadiationSpectra.peakfinder) *peakquantile::ClosedInterval{<:Real}=0.5..1.0: quantiles that define energy window that is considered peak search (default: 0.0..1.0 == whole range). Allgammalineshave to be within this window. *binquantile::ClosedInterval{<:Real}=0.5..1.0`: quantiles that define energy window that is used to find optimal binning (default same as peakquantile) * quantile_perc::Float64=NaN: If NaN the standard peakfinder is used. If not NaN, the peak for simple calibration is set to given quantile energy. Useful for debugging/testing

output: result NamedTuple with the following fields:

• h_uncal: histogram of the uncalibrated energy array with the binning that is used in peak search
• c: simple calibration factor
• peak_guess: estimated peak energy in units of the uncalibrated energy array

peakhists_gamma(e::Vector{<:Unitful.Energy{<:Real}}, gamma_lines::Vector{<:Unitful.Energy{<:Real}}, window_sizes::Vector{<:Tuple{Unitful.Energy{<:Real}, Unitful.Energy{<:Real}}},; binning_peak_window::Unitful.Energy{<:Real}=10.0u"keV")

Create histograms around the calibration lines and return the histograms and the peak statistics.

input

* `e`: energy array
* `gamma_lines`: array of gamma lines
* `window_sizes`: array of tuples with left and right window sizes around the gamma lines

keyword Arguments:

* `binning_peak_window::Unitful.Energy{<:Real}`: energy window around each peak that is used to find optimal binning for peakhists

output

* `peakhists`: array of histograms around the calibration lines
* `peakstats`: array of statistics for the calibration line fits

peakshape_components(fit_func::Symbol; background_center::Real)

This function defines the components (signal, low/high-energy tail, backgrounds) of the fit function used in gamma specfits. These component functions are used in the fit-report and in plot receipes

prepare_sep_peakhist(e::Array{T}, dep::T,; relative_cut::T=0.5, n_bins_cut::Int=500) where T<:Real

Prepare an array of uncalibrated SEP energies for parameter extraction and calibration.

Returns

• result: Result of the initial fit
• report: Report of the initial fit

qc_window_cut(tab::Table, config::PropDict, cols::NTuple{<:Any, Symbol})

Perform quality control on the data using a window cut based on a centered Gaussian fit. The config is a PropDict containing the configuration for each column to be checked. Each column should have a min, max, and sigma value, along with optional keyword arguments for the cut.

Arguments

* `data`: a `Table` containing the data to be checked
* `config`: a `PropDict` containing the configuration for each column to be checked
* `cols`: a tuple of column names to be checked

Returns

* `result`: a `NamedTuple` containing the results of the quality control, with a `qc` field that is a string describing the cut condition for each column
* `report`: an `OrderedDict` containing the reports for each column which can be plotted, with fields like `f_fit`, `h`, `μ`, `σ`, `gof`, `low_cut`, and `high_cut`

set_memlimit(gig::Float64)
Set memory limit in GB for the whole current process ignoring if things could run in parallel

set_timelimit(sec::Float64)
Set time limit in seconds a single Optim.jl optimization step

signal_peakshape(
    x::Real, μ::Real, σ::Real, n::Real,
    skew_fraction::Real;  skew_fraction_highE::Real = 0.0
)

Describes the signal part of the shape of a typical gamma peak in a detector.

simple_calibration(e_uncal::Array, th228_lines::Array, window_size::Float64=25.0, n_bins::Int=15000, calib_type::String="th228")

Perform a simple calibration for the uncalibrated energy array e_uncal using the calibration type calib_type and the calibration lines th228_lines. The window size is the size of the window around the calibration line to use for the calibration. The number of bins is the number of bins to use for the histogram.

Returns * h_calsimple: histogram of the calibrated energy array * h_uncal: histogram of the uncalibrated energy array * c: calibration factor * fep_guess: estimated full energy peak (FEP) * peakhists: array of histograms around the calibration lines * peakstats: array of statistics for the calibration line fits

sipm_simple_calibration(pe_uncal::Vector{<:Real})
sipm_simple_calibration(pe_uncal_vov::VectorOfVectors{<:Real})

Simple SiPM calibration from the 1 PE and 2 PE peak positions found by a peakfinder. The noise/1PE cut is auto-detected at the centroid of the first valley after the noise peak.

n_fwhm_noise_cut caps the valley search window in units of the noise-peak half-width above its center. 0.0 skips valley detection and uses initial_min_amp. Negative values force the legacy formula cuts_1pe.max + n_fwhm_noise_cut * fwhm_noise (kept for backward-compatibility with configs that dip into the noise peak).

For VoV input, an optional second-stage QC (single_trigger_only=true) keeps only waveforms with exactly one trigger above the detected threshold.

kwargs: * initial_min_amp, initial_max_amp: histogram bounds for the noise/peak search * relative_cut_noise_cut: passed to cut_single_peak * n_fwhm_noise_cut: see above * single_trigger_only (VoV only) * min_pe_peak, max_pe_peak, peakfinder_*: peakfinder controls

Returns (result, report). result carries the calibration (f_simple_calib, c, offset, peakpos, noisepeakpos, noisepeakwidth, noise_threshold, noise_threshold_cal). report carries plotting data (peakpos*, h_uncal, h_calsimple, h_*_full for the unfiltered all-trigger spectra, noise_threshold*, valley_found).

LegendSpecFits.stabilize_dist(d::Distribution)

Return a numerically stable version of the distribution d.

Tries to ensure that the (log-)PDF of the result is finite over its whole support.

step_gauss(x::Real, μ::Real, σ::Real)

Evaluates the convulution of a Heaviside step function and the PDF of Normal(μ, σ) at x.

The result does not correspond to a PDF as it is not normalizable.

subhist(h::Histogram, r::Tuple{<:Real,<:Real})

Return a new Histogram with the bins in the range r.

tuple_to_array(nt::NamedTuple, fields::Vector{Symbol})

Return an array with the values of the fields in nt in the order given by fields.

weibull_from_mx(m::Real, x::Real, p_x::Real = 0.6827)::Weibull

Construct a Weibull distribution with a given median m and a given p_x-quantile x.

Useful to construct priors for positive quantities.