reboost.hpge package¶
Submodules¶
reboost.hpge.psd module¶
- reboost.hpge.psd._current_pulse_model(times, amax, mu, sigma, tail_fraction, tau, high_tail_fraction=0, high_tau=0)¶
Analytic model for the current pulse in a Germanium detector.
Consists of a Gaussian, a high side exponential tail and a low side tail:
\[\begin{split}\begin{align} A(t) = \; &A_\text{max} \times (1-p-p_h) \times \text{Gauss}(t;\mu,\sigma) \\ &+ A \times p \; \left(1 - \text{erf}\left(\frac{t-\mu}{\sigma_i}\right)\right) \times \frac{e^{t/\tau}}{2e^{\mu/\tau}} \\ &+ A \times p_h \; \left(1 - \text{erf}\left(-\frac{t-\mu}{\sigma_i}\right)\right) \times \frac{1}{2}e^{-t/\tau} \end{align}\end{split}\]- Parameters:
times (ArrayLike) – Array of times to compute current for.
amax (float) – Maximum current for the template
mu (float) – Time of the maximum current.
sigma (float) – Width of the current pulse
tail_fraction (float) – Fraction of the tail in the pulse.
tau (float) – Time constant of the low time tail.
high__tail_fraction – Fraction of the high tail in the pulse.
high_tau (float) – Time constant of the high time tail.
high_tail_fraction (float)
- Returns:
The predicted current waveform for this energy deposit.
- Return type:
NDArray
- reboost.hpge.psd._drift_time_heuristic_impl(dt, edep)¶
Low-level implementation of the HPGe drift-time-based pulse-shape heuristic.
Accepts Awkward arrays and uses Numba to speed up the computation.
For each hit (collection of steps), the drift times and corresponding energies are sorted in ascending order. The function finds the optimal split point \(m\) that maximizes the identification metric:
\[I = \frac{|T_1 - T_2|}{E_\text{s}(E_1, E_2)}\]where:
\[T_1 = \frac{\sum_{i < m} t_i E_i}{\sum_{i < m} E_i} \quad \text{and} \quad T_2 = \frac{\sum_{i \geq m} t_i E_i}{\sum_{i \geq m} E_i}\]are the energy-weighted mean drift times of the two groups.
\[E_\text{scale}(E_1, E_2) = \frac{1}{\sqrt{E_1 E_2}}\]is the scaling factor.
The function iterates over all possible values of \(m\) and selects the maximum I as the drift time heuristic value.
- reboost.hpge.psd._estimate_current_impl(edep, dt, dist_to_nplus, template, times, include_surface_effects, fccd, templates_surface, activeness_surface, surface_step_in_um)¶
Estimate the maximum current that would be measured in the HPGe detector.
This is based on extracting a waveform with
get_current_waveform()and finding the maxima of it.- Parameters:
edep (ak.Array) – Array of energies for each step.
dt (ak.Array) – Array of drift times for each step.
dist_to_nplus (ak.Array) – Array of distance to nplus contact for each step (can be None, in which case no surface effects are included.)
template (np.array) – array of the bulk pulse template
times (np.array) – time-stamps for the bulk pulse template
include_surface_effects (bool)
fccd (float)
templates_surface (np.array)
activeness_surface (np.array)
surface_step_in_um (float)
- Return type:
tuple[NDArray, NDArray, NDArray]
- reboost.hpge.psd._get_waveform_maximum_impl(t, e, dist, template, templates_surface, activeness_surface, tmin, tmax, start, fccd, n, time_step, surface_step_in_um, include_surface_effects)¶
Basic implementation to get the maximum of the waveform.
- Parameters:
t (ArrayLike) – drift time for each step.
e (ArrayLike) – energy for each step.
dist (ArrayLike) – distance to surface for each step.
template (ArrayLike)
templates_surface (ArrayLike)
activeness_surface (ArrayLike)
tmin (float)
tmax (float)
start (float)
fccd (float)
n (int)
time_step (int)
surface_step_in_um (float)
include_surface_effects (bool)
- reboost.hpge.psd._get_waveform_value(idx, edep, drift_time, template, start, dt)¶
Get the value of the waveform at a certain index.
- Parameters:
idx (int) – the index of the time array to find the waveform at.
edep (ak.Array) – Array of energies for each step
drift_time (ak.Array) – Array of drift times for each step
template (ArrayLike) – array of the template for the current waveforms
start (float) – first time value of the template
dt (float) – timestep (in ns) for the template.
- Returns:
Value of the current waveform
- Return type:
- reboost.hpge.psd._get_waveform_value_surface(idx, edep, drift_time, dist_to_nplus, bulk_template, templates_surface, activeness_surface, distance_step_in_um, fccd, start, dt)¶
Get the value of the waveform at a certain index.
- Parameters:
idx (int) – the index of the time array to find the waveform at.
edep (NDArray) – Array of energies for each step
drift_time (np.array) – Array of drift times for each step
template – array of the template for the current waveforms
templates_surface (ArrayLike) – The current templates from the surface.
activeness_surface (ArrayLike) – The total collected charge for each surface point.
dist_step_in_um – The binning in distance for the surface pulse library.
start (float) – first time value of the template
dt (float) – timestep (in ns) for the template.
dist_to_nplus (np.array)
bulk_template (ArrayLike)
distance_step_in_um (float)
fccd (float)
- Returns:
Value of the current waveform and the energy.
- Return type:
- reboost.hpge.psd._interpolate_pulse_model(template, time, start, end, dt, mu)¶
Interpolate to extract the pulse model given a particular mu.
- reboost.hpge.psd._njit_erf(x)¶
Error function that can take in a numpy array.
- Parameters:
x (ArrayLike)
- Return type:
NDArray
- reboost.hpge.psd.convolve_surface_response(surf_current, bulk_pulse)¶
Convolve the surface response pulse with the bulk current pulse.
This combines the current induced on the edge of the FCCD region with the bulk response on the p+ contact.
- reboost.hpge.psd.drift_time(xloc, yloc, zloc, dt_map, coord_offset=<Quantity([0 0 0], 'meter')>)¶
Calculates drift times for each step (cluster) in an HPGe detector.
- Parameters:
xloc (ArrayLike) – awkward array of x coordinate position.
yloc (ArrayLike) – awkward array of y coordinate position.
zloc (ArrayLike) – awkward array of z coordinate position.
dt_map (HPGeScalarRZField) – the drift time map.
coord_offset (Quantity | Position) – this (x, y, z) coordinates will be subtracted to (xloc, yloc, zloc)` before drift time computation. The length units must be the same as xloc, yloc and zloc.
- Return type:
- reboost.hpge.psd.drift_time_heuristic(drift_time, edep)¶
HPGe drift-time-based pulse-shape heuristic.
See
_drift_time_heuristic_impl()for a description of the algorithm.- Parameters:
drift_time (ArrayLike) – drift time of charges originating from steps/clusters. Can be calculated with
drift_time().edep (ArrayLike) – energy deposited in step/cluster (same shape as drift_time).
- Return type:
- reboost.hpge.psd.get_current_template(low=-1000, high=4000, step=1, mean_aoe=1, **kwargs)¶
Build the current template from the analytic model, defined by
_current_pulse_model().- Parameters:
low (float) – start of the template
high (float) – end of the template
step (float) – time-step, this should divide high-low
mean_aoe (float) – The mean AoE value for this detector (to normalise current pulses).
**kwargs – Other keyword arguments passed to
_current_pulse_model().
- Returns:
tuple of the (template,times)
- Return type:
tuple[NDArray, NDArray]
- reboost.hpge.psd.get_current_waveform(edep, drift_time, template, start, dt, range_t)¶
Estimate the current waveform.
Based on modelling the current as a sum over the current pulse model defined by the template.
\[A(t) = \sum_i E_i \times N f(t, dt_i, \vec{\theta})\]- Where:
\(f(t)\) is the template
:math`vec{theta}` are the parameters \((\sigma, p, \tau)\)
\(E_i\) and \(dt_i\) are the deposited energy and drift time.
\(N\) is a normalisation term
- Parameters:
edep (ak.Array) – Array of energies for each step
drift_time (ak.Array) – Array of drift times for each step
template (ArrayLike) – array of the template for the current waveforms
start (float) – first time value of the template
dt (float) – timestep (in ns) for the template.
range_t (tuple) – a range of times to search around
- Returns:
A tuple of the time and current for the current waveform for this event.
- Return type:
tuple(NDArray, NDArray)
- reboost.hpge.psd.make_convolved_surface_library(bulk_template, surface_library)¶
Make the convolved surface library out of the template.
This convolves every row of the surface_library with the template and reshapes the output to match the initial template. It returns a 2D array with one more row than the surface_library and each row the same length as the template. The final row is the bulk_template for easier interpolation.
- Parameters:
bulk_template (array) – The template for the bulk response
surface_library (array) – The 2D array of the surface library.
- Returns:
2D array of the surface library convolved with the bulk response.
- Return type:
NDArray
- reboost.hpge.psd.maximum_current(edep, drift_time, dist_to_nplus=None, *, template, times, fccd_in_um=0, templates_surface=None, activeness_surface=None, surface_step_in_um=10, return_mode='current')¶
Estimate the maximum current in the HPGe detector based on
_estimate_current_impl().- Parameters:
edep (ArrayLike) – Array of energies for each step.
drift_time (ArrayLike) – Array of drift times for each step.
dist_to_nplus (ArrayLike | None) – Distance to n-plus electrode, only needed if surface heuristics are enabled.
template (np.array) – array of the bulk pulse template
times (np.array) – time-stamps for the bulk pulse template
fccd – Value of the full-charge-collection depth, if None no surface corrections are performed.
surface_library – 2D array (distance, time) of the rate of charge arriving at the p-n junction. Each row should be an array of length 10000 giving the charge arriving at the p-n junction for each timestep (in ns). This is produced by
hpge.surface.get_surface_response()or other libraries.surface_step_in_um (float) – Distance step for the surface library.
return_mode (str) – either current, energy or max_time
fccd_in_um (float)
templates_surface (ArrayLike | None)
activeness_surface (ArrayLike | None)
- Returns:
An Array of the maximum current/ time / energy for each hit.
- Return type:
Array
reboost.hpge.surface module¶
- reboost.hpge.surface._advance_diffusion(charge, factor, recomb=0, recomb_depth=600, delta_x=10)¶
Make a step of diffusion using explicit Euler scheme.
- Parameters:
- Returns:
a tuple of the charge distribution at the next time step and the collected charge.
- reboost.hpge.surface._compute_diffusion_impl(init_charge, nsteps, factor, recomb=0, recomb_depth=600, delta_x=10)¶
Compute the charge collected as a function of time.
- reboost.hpge.surface.distance_to_surface(positions_x, positions_y, positions_z, hpge, det_pos, *, surface_type=None, unit='mm', distances_precompute=None, precompute_cutoff=None)¶
Computes the distance from each step to the detector surface.
The calculation can be performed for any surface type nplus, pplus, passive or None. In order to speed up the calculation we provide an option to only compute the distance for points within a certain distance of any surface (as computed by remage and stored in the “distances_precompute”) argument.
- Parameters:
positions_x (VectorOfVectors) – Global x positions for each step.
positions_y (VectorOfVectors) – Global y positions for each step.
positions_z (VectorOfVectors) – Global z positions for each step.
hpge (HPGe) – HPGe object.
det_pos (ArrayLike) – position of the detector origin, must be a 3 component array corresponding to (x,y,z).
surface_type (str | None) – string of which surface to use, can be nplus, pplus passive or None (in which case the distance to any surface is calculated).
unit (str) – unit for the hit tier positions table.
distances_precompute (VectorOfVectors | None) – VectorOfVectors of distance to any surface computed by remage.
precompute_cutoff (float | None) – cutoff on distances_precompute to not compute the distance for (in mm)
- Returns:
VectorOfVectors with the same shape as positions_x/y/z of the distance to the surface.
- Return type:
Note
positions_x/positions_y/positions_z must all have the same shape.
- reboost.hpge.surface.get_surface_library(fccd, dist_step_in_um, **kwargs)¶
Build the surface response library by calling reboost.hpge.surface.get_surface_response.
- Parameters:
- Returns:
2D array of the cumulative-charge arriving at the p-n junction as a function of time, for each distance.
- reboost.hpge.surface.get_surface_response(fccd, init=0, *, recomb_depth=500, recomb=0.002, init_size=0.0, factor=0.29, nsteps=10000, delta_x=10)¶
Extract the surface response current pulse based on diffusion.
This extracts the amount of charge arrived (cumulative) at the p-n junction as a function of time, based on diffusion. The final induced waveform on the p-n contact is obtained from convolution with the bulk pulse.
- Parameters:
fccd (float) – the full charge collection depth (in um)
recomb_depth (float) – the depth of the recombination region (in um)
init (float) – the initial position of the charge (in um)
recomb (float) – the recombination rate
init_size (float) – the initial size of the charge cloud (in um)
factor (float) – the factor for the explicit Euler scheme (the probability of charge diffusuion)
nsteps (int) – the number of time steps.
delta_x (float) – the width of each position bin.
reboost.hpge.utils module¶
- class reboost.hpge.utils.HPGeScalarRZField(φ, r_units, z_units, φ_units)¶
Bases:
NamedTupleA scalar field defined in the cylindrical-like (r, z) HPGe plane.
Create new instance of HPGeScalarRZField(φ, r_units, z_units, φ_units)
- _asdict()¶
Return a new dict which maps field names to their values.
- _field_defaults = {}¶
- _fields = ('φ', 'r_units', 'z_units', 'φ_units')¶
- classmethod _make(iterable)¶
Make a new HPGeScalarRZField object from a sequence or iterable
- _replace(**kwds)¶
Return a new HPGeScalarRZField object replacing specified fields with new values
- reboost.hpge.utils.get_hpge_scalar_rz_field(filename, obj, field, out_of_bounds_val=nan, **kwargs)¶
Create an interpolator for a gridded scalar HPGe field defined on (r, z).
Reads from disk the following data structure:
FILENAME/ └── OBJ · struct{r,z,FIELD} ├── r · array<1>{real} ── {'units': 'UNITS'} ├── z · array<1>{real} ── {'units': 'UNITS'} └── FIELD · array<2>{real} ── {'units': 'UNITS'}where
FILENAME,OBJandFIELDare provided as arguments to this function. obj is aStruct, r and z are one dimensional arrays specifying the radial and z coordinates of the rectangular grid — not the coordinates of each single grid point. In this coordinate system, the center of the p+ contact surface is at (0, 0), with the p+ contact facing downwards. field is instead a two-dimensional array specifying the field value at each grid point. The first and second dimensions are r and z, respectively. NaN values are interpreted as points outside the detector profile in the (r, z) plane.Before returning a
HPGeScalarRZField, the gridded field is fed toscipy.interpolate.RegularGridInterpolator.- Parameters:
- Return type: