From dae6dd73a02020c09fd05b347b2044bba366ba24 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Esteban=20Le=C3=B3n?= <54425957+esleon97@users.noreply.github.com> Date: Fri, 10 Nov 2023 13:31:51 -0500 Subject: [PATCH 01/42] Update min_max.py Added min/max normalization for ML data cleaning. --- src/dspeed/processors/min_max.py | 53 ++++++++++++++++++++++++++++++++ 1 file changed, 53 insertions(+) diff --git a/src/dspeed/processors/min_max.py b/src/dspeed/processors/min_max.py index 01521d6..cad858a 100644 --- a/src/dspeed/processors/min_max.py +++ b/src/dspeed/processors/min_max.py @@ -70,3 +70,56 @@ def min_max( a_max[0] = w_in[max_index] t_min[0] = float(min_index) t_max[0] = float(max_index) + +@guvectorize( + [ + "void(float32[:], float32[:], float32[:], float32[:])", + "void(float64[:], float64[:], float64[:], float64[:])", + ], + "(n),(),()->(n)", + **nb_kwargs, +) +def min_max_norm( + w_in: np.ndarray, a_min: float, a_max: float, w_out: np.ndarray +) -> None: + """Normalize waveform by minimum or maximum value, whichever is + greater in absolute value. + + + Parameters + ---------- + w_in + the input waveform + a_min + the minimum value + a_max + the maximum value + w_out + the normalized output waveform + + JSON Configuration Example + -------------------------- + + .. code-block :: json + + "wf_norm": { + "function": "min_max_norm", + "module": "dspeed.processors", + "args": ["wf_blsub", "wf_min", "wf_max", "wf_norm"], + "unit": ["ADC"] + } + """ + + if np.isnan(w_in).any(): + return + + w_out[:] = np.nan + + if abs(a_max[0]) == 0 or abs(a_min[0]) == 0: + w_out[:] = w_in[:] + + elif abs(a_max[0]) >= abs(a_min[0]): + w_out[:] = w_in[:] / abs(a_max[0]) + + elif abs(a_max[0]) < abs(a_min[0]): + w_out[:] = w_in[:] / abs(a_min[0]) From a27b3f2abae147d55a02dc217650c141c3c71478 Mon Sep 17 00:00:00 2001 From: Esteban Leon Date: Fri, 10 Nov 2023 15:21:50 -0800 Subject: [PATCH 02/42] Added min_max_norm to processors/__init__.py --- src/dspeed/processors/__init__.py | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/src/dspeed/processors/__init__.py b/src/dspeed/processors/__init__.py index 78256c9..53d0639 100644 --- a/src/dspeed/processors/__init__.py +++ b/src/dspeed/processors/__init__.py @@ -70,7 +70,7 @@ from .histogram import histogram, histogram_stats from .linear_slope_fit import linear_slope_diff, linear_slope_fit from .log_check import log_check -from .min_max import min_max +from .min_max import min_max, min_max_norm from .moving_windows import ( avg_current, moving_window_left, @@ -112,6 +112,7 @@ "linear_slope_diff", "log_check", "min_max", + "min_max_norm", "avg_current", "moving_window_left", "moving_window_multi", From 30b21927e0162df97108dcd433c75472772ddaad Mon Sep 17 00:00:00 2001 From: ggmarshall Date: Mon, 4 Dec 2023 23:02:39 +0100 Subject: [PATCH 03/42] new general convolution functions moved kernel generation to other files --- src/dspeed/processors/convolutions.py | 461 ++++---------------------- 1 file changed, 67 insertions(+), 394 deletions(-) diff --git a/src/dspeed/processors/convolutions.py b/src/dspeed/processors/convolutions.py index 4c43336..f5ce486 100644 --- a/src/dspeed/processors/convolutions.py +++ b/src/dspeed/processors/convolutions.py @@ -1,422 +1,95 @@ from __future__ import annotations -from typing import Callable - import numpy as np from numba import guvectorize +from scipy.signal import fftconvolve from ..errors import DSPFatal from ..utils import numba_defaults_kwargs as nb_kwargs - -def cusp_filter(length: int, sigma: float, flat: int, decay: int) -> Callable: - """Apply a CUSP filter to the waveform. - - Note - ---- - This processor is composed of a factory function that is called using the - `init_args` argument. The input and output waveforms are passed using - `args`. - - Parameters - ---------- - length - the length of the filter to be convolved. - sigma - the curvature of the rising and falling part of the kernel. - flat - the length of the flat section. - decay - the decay constant of the exponential to be convolved. - - JSON Configuration Example - -------------------------- - - .. code-block :: json - - "wf_cusp": { - "function": "cusp_filter", - "module": "dspeed.processors", - "args": ["wf_bl", "wf_cusp(101,f)"], - "unit": "ADC", - "init_args": ["len(wf_bl)-100", "40*us", "3*us", "45*us"] - } - """ - if length <= 0: - raise DSPFatal("The length of the filter must be positive") - - if np.floor(length) != length: - raise DSPFatal("The length of the filter must be an integer") - - if sigma < 0: - raise DSPFatal("The curvature parameter must be positive") - - if flat < 0: - raise DSPFatal("The length of the flat section must be positive") - - if np.floor(flat) != flat: - raise DSPFatal("The length of the flat section must be an integer") - - if decay < 0: - raise DSPFatal("The decay constant must be positive") - - lt = int((length - flat) / 2) - flat_int = int(flat) - cusp = np.zeros(length) - for ind in range(0, lt, 1): - cusp[ind] = float(np.sinh(ind / sigma) / np.sinh(lt / sigma)) - for ind in range(lt, lt + flat_int + 1, 1): - cusp[ind] = 1 - for ind in range(lt + flat_int + 1, length, 1): - cusp[ind] = float(np.sinh((length - ind) / sigma) / np.sinh(lt / sigma)) - - den = [1, -np.exp(-1 / decay)] - cuspd = np.convolve(cusp, den, "same") - - @guvectorize( - ["void(float32[:], float32[:])", "void(float64[:], float64[:])"], - "(n),(m)", +@guvectorize( + ["void(float32[:], float32[:], char, float32[:])", + "void(float64[:], float64[:], char, float64[:])"], + "(n),(m),(),(p)", **nb_kwargs( cache=False, forceobj=True, ), ) - def cusp_out(w_in: np.ndarray, w_out: np.ndarray) -> None: - """ - Parameters - ---------- - w_in - the input waveform. - w_out - the filtered waveform. - """ - w_out[:] = np.nan - - if np.isnan(w_in).any(): - return - - if len(cuspd) > len(w_in): - raise DSPFatal("The filter is longer than the input waveform") - - w_out[:] = np.convolve(w_in, cuspd, mode="valid") - - return cusp_out - - -def zac_filter(length: int, sigma: float, flat: int, decay: int) -> Callable: - """Apply a ZAC (Zero Area CUSP) filter to the waveform. - - Note - ---- - This processor is composed of a factory function that is called using the - `init_args` argument. The input and output waveforms are passed using - `args`. - - Parameters - ---------- - length - the length of the filter to be convolved. - sigma - the curvature of the rising and falling part of the kernel. - flat - the length of the flat section. - decay - the decay constant of the exponential to be convolved. - - JSON Configuration Example - -------------------------- - - .. code-block :: json - - "wf_zac": { - "function": "zac_filter", - "module": "dspeed.processors", - "args": ["wf_bl", "wf_zac(101,f)"], - "unit": "ADC", - "init_args": ["len(wf_bl)-100", "40*us", "3*us", "45*us"] - } +def convolve_wf(w_in: np.ndarray, kernel:np.array, mode_in: np.int8, w_out: np.ndarray) -> None: # """ - if length <= 0: - raise DSPFatal("The length of the filter must be positive") - - if np.floor(length) != length: - raise DSPFatal("The length of the filter must be an integer") - - if sigma < 0: - raise DSPFatal("The curvature parameter must be positive") - - if flat < 0: - raise DSPFatal("The length of the flat section must be positive") - - if np.floor(flat) != flat: - raise DSPFatal("The length of the flat section must be an integer") - - if decay < 0: - raise DSPFatal("The decay constant must be positive") - - lt = int((length - flat) / 2) - flat_int = int(flat) - - # calculate cusp filter and negative parables - cusp = np.zeros(length) - par = np.zeros(length) - for ind in range(0, lt, 1): - cusp[ind] = float(np.sinh(ind / sigma) / np.sinh(lt / sigma)) - par[ind] = np.power(ind - lt / 2, 2) - np.power(lt / 2, 2) - for ind in range(lt, lt + flat_int + 1, 1): - cusp[ind] = 1 - for ind in range(lt + flat_int + 1, length, 1): - cusp[ind] = float(np.sinh((length - ind) / sigma) / np.sinh(lt / sigma)) - par[ind] = np.power(length - ind - lt / 2, 2) - np.power(lt / 2, 2) - - # calculate area of cusp and parables - areapar, areacusp = 0, 0 - for i in range(0, length, 1): - areapar += par[i] - areacusp += cusp[i] - - # normalize parables area - par = -par / areapar * areacusp - - # create zac filter - zac = cusp + par - - # deconvolve zac filter - den = [1, -np.exp(-1 / decay)] - zacd = np.convolve(zac, den, "same") - - @guvectorize( - ["void(float32[:], float32[:])", "void(float64[:], float64[:])"], - "(n),(m)", - **nb_kwargs( - cache=False, - forceobj=True, - ), - ) - def zac_out(w_in: np.ndarray, w_out: np.ndarray) -> None: - """ - Parameters - ---------- - w_in - the input waveform. - w_out - the filtered waveform. - """ - w_out[:] = np.nan - - if np.isnan(w_in).any(): - return - - if len(zacd) > len(w_in): - raise DSPFatal("The filter is longer than the input waveform") - - w_out[:] = np.convolve(w_in, zacd, mode="valid") - - return zac_out - - -def t0_filter(rise: int, fall: int) -> Callable: - """Apply a modified, asymmetric trapezoidal filter to the waveform. - - Note - ---- - This processor is composed of a factory function that is called using the - `init_args` argument. The input and output waveforms are passed using - `args`. - Parameters ---------- - rise - the length of the rise section. This is the linearly increasing - section of the filter that performs a weighted average. - fall - the length of the fall section. This is the simple averaging part - of the filter. - - JSON Configuration Example - -------------------------- - - .. code-block :: json - - "wf_t0_filter": { - "function": "t0_filter", - "module": "dspeed.processors", - "args": ["wf_pz", "wf_t0_filter(3748,f)"], - "unit": "ADC", - "init_args": ["128*ns", "2*us"] - } + w_in + the input waveform. + kernel + the kernel to convolve with + mode + mode of convolution options are f : full, v : valid or s : same, + explained here: https://numpy.org/doc/stable/reference/generated/numpy.convolve.html + w_out + the filtered waveform. """ - if rise < 0: - raise DSPFatal("The length of the rise section must be positive") - - if fall < 0: - raise DSPFatal("The length of the fall section must be positive") - - t0_kern = np.zeros(int(rise) + int(fall)) - for i in range(int(rise)): - t0_kern[i] = 2 * (int(rise) - i) / (rise**2) - for i in range(int(rise), len(t0_kern), 1): - t0_kern[i] = -1 / fall - - @guvectorize( - ["void(float32[:], float32[:])", "void(float64[:], float64[:])"], - "(n),(m)", + w_out[:] = np.nan + + if np.isnan(w_in).any(): + return + + if np.isnan(kernel).any(): + return + + if len(kernel) > len(w_in): + raise DSPFatal("The filter is longer than the input waveform") + + if chr(mode_in) not in ["f", "v", "s"]: + raise DSPFatal("Invalid mode") + + w_out[:] = np.convolve(w_in, kernel, mode=chr(mode_in)) + + +@guvectorize( + ["void(float32[:], float32[:], char, float32[:])", + "void(float64[:], float64[:], char, float64[:])"], + "(n),(m),(),(p)", **nb_kwargs( cache=False, forceobj=True, ), ) - def t0_filter_out(w_in: np.ndarray, w_out: np.ndarray) -> None: - """ - Parameters - ---------- - w_in - The input waveform - w_out - The filtered waveform - """ - w_out[:] = np.nan - - if np.isnan(w_in).any(): - return - - if len(t0_kern) > len(w_in): - raise DSPFatal("The filter is longer than the input waveform") - - w_out[:] = np.convolve(w_in, t0_kern)[: len(w_in)] - - return t0_filter_out - - -def moving_slope(length): - """Calculates the linear slope of a waveform in sections of length - - Note - ---- - This processor is composed of a factory function that is called using the - `init_args` argument. The input and output waveforms are passed using - `args`. - - Parameters - ---------- - length - the length of the section to calculate slope - - JSON Configuration Example - -------------------------- - - .. code-block :: json - - "wf_slopes": { - "function": "moving_slope", - "module": "dspeed.processors", - "args": ["wf_pz", "wf_slopes(len(wf_pz)-11,f)"], - "unit": "ADC", - "init_args": ["12"] - } +def fft_convolve_wf(w_in: np.ndarray, kernel:np.array, mode_in: np.int8, w_out: np.ndarray) -> None: # """ - - if length <= 0: - raise DSPFatal("The length of the filter must be positive") - - if np.floor(length) != length: - raise DSPFatal("The length of the filter must be an integer") - - sum_x = length * (length + 1) / 2 - sum_x2 = length * (length + 1) * (2 * length + 1) / 6 - - kernel = (np.arange(1, length + 1, 1) * length) - (np.ones(length) * sum_x) - kernel /= length * sum_x2 - sum_x * sum_x - kernel = kernel[::-1] - - @guvectorize( - ["void(float32[:], float32[:])", "void(float64[:], float64[:])"], - "(n),(m)", - **nb_kwargs( - cache=False, - forceobj=True, - ), - ) - def moving_slope_out(w_in, w_out): - w_out[:] = np.nan - - if np.isnan(w_in).any(): - return - - if len(kernel) > len(w_in): - raise DSPFatal("The filter is longer than the input waveform") - - w_out[:] = np.convolve(w_in, kernel, mode="valid") - - return moving_slope_out - - -def step(length: int) -> Callable: - """Process waveforms with a step function. - - Note - ---- - This processor is composed of a factory function that is called using the - `init_args` argument. The input and output waveforms are passed using - `args`. - Parameters ---------- - length - length of the step function. - weight_pos - relative weight of positive step side. - - JSON Configuration Example - -------------------------- - - .. code-block :: json - - "wf_step": { - "function": "step", - "module": "dspeed.processors", - "args": ["waveform", "wf_step(len(waveform)-15,f)"], - "unit": "ADC", - "init_args": ["16"] - } + w_in + the input waveform. + kernel + the kernel to convolve with + mode + mode of convolution options are f : full, v : valid or s : same, + explained here: https://numpy.org/doc/stable/reference/generated/numpy.convolve.html + w_out + the filtered waveform. """ + w_out[:] = np.nan + + if np.isnan(w_in).any(): + return + + if np.isnan(kernel).any(): + return + + if len(kernel) > len(w_in): + raise DSPFatal("The filter is longer than the input waveform") + + if chr(mode_in) == "f": + mode = "full" + elif chr(mode_in) == "v": + mode = "valid" + elif chr(mode_in) == "s": + mode = "same" + else: + raise DSPFatal("Invalid mode") + + w_out[:] = fftconvolve(w_in, kernel, mode=mode) - x = np.arange(length) - y = np.piecewise( - x, - [ - ((x >= 0) & (x < length / 4)), - ((x >= length / 4) & (x < 3 * length / 4)), - ((x >= 3 * length / 4) & (x < length)), - ], - [-1, 1, -1], - ) - - @guvectorize( - ["void(float32[:], float32[:])", "void(float64[:], float64[:])"], - "(n),(m)", - **nb_kwargs( - cache=False, - forceobj=True, - ), - ) - def step_out(w_in: np.ndarray, w_out: np.ndarray) -> None: - """ - Parameters - ---------- - w_in - the input waveform. - w_out - the filtered waveform. - """ - - w_out[:] = np.nan - - if np.isnan(w_in).any(): - return - - if len(y) > len(w_in): - raise DSPFatal("The filter is longer than the input waveform") - w_out[:] = np.convolve(w_in, y, mode="valid") - return step_out From 08574fcf017cfd1bb7ba4b4b7f9305401b4e1c66 Mon Sep 17 00:00:00 2001 From: ggmarshall Date: Mon, 4 Dec 2023 23:03:11 +0100 Subject: [PATCH 04/42] added kernel generation functions moved dplms into energy_kernels --- src/dspeed/processors/__init__.py | 7 +- src/dspeed/processors/dplms.py | 170 --------------- src/dspeed/processors/energy_kernels.py | 271 ++++++++++++++++++++++++ src/dspeed/processors/kernels.py | 137 ++++++++++++ 4 files changed, 413 insertions(+), 172 deletions(-) delete mode 100644 src/dspeed/processors/dplms.py create mode 100644 src/dspeed/processors/energy_kernels.py create mode 100644 src/dspeed/processors/kernels.py diff --git a/src/dspeed/processors/__init__.py b/src/dspeed/processors/__init__.py index 515af9d..4757214 100644 --- a/src/dspeed/processors/__init__.py +++ b/src/dspeed/processors/__init__.py @@ -60,15 +60,16 @@ """ from .bl_subtract import bl_subtract -from .convolutions import cusp_filter, moving_slope, step, t0_filter, zac_filter -from .dplms import dplms +from .convolutions import convolve_wf, fft_convolve_wf from .dwt import discrete_wavelet_transform +from .energy_kernels import cusp_filter, zac_filter, dplms from .fftw import dft, inv_dft, psd from .fixed_time_pickoff import fixed_time_pickoff from .gaussian_filter1d import gaussian_filter1d from .get_multi_local_extrema import get_multi_local_extrema from .get_wf_centroid import get_wf_centroid from .histogram import histogram, histogram_stats +from .kernels import t0_filter, moving_slope, step from .linear_slope_fit import linear_slope_diff, linear_slope_fit from .log_check import log_check from .min_max import min_max @@ -100,6 +101,8 @@ __all__ = [ "bl_subtract", + "convolve_wf", + "fft_convolve_wf", "cusp_filter", "t0_filter", "zac_filter", diff --git a/src/dspeed/processors/dplms.py b/src/dspeed/processors/dplms.py deleted file mode 100644 index 8ba62f4..0000000 --- a/src/dspeed/processors/dplms.py +++ /dev/null @@ -1,170 +0,0 @@ -from __future__ import annotations - -from typing import Callable - -import numpy as np -from numba import guvectorize - -from ..errors import DSPFatal -from ..utils import numba_defaults_kwargs as nb_kwargs - - -def dplms( - length: int, - noise_mat: list, - reference: list, - a1: float, - a2: float, - a3: float, - ff: int, - coefficients: list, -) -> Callable: - """Calculate and apply an optimum DPLMS filter to the waveform. - - The processor takes the noise matrix and the reference signal as input and - calculates the optimum filter according to the provided length and - penalized coefficients [DPLMS]_. - - .. [DPLMS] V. D'Andrea et al. “Optimum Filter Synthesis with DPLMS - Method for Energy Reconstruction” Eur. Phys. J. C 83, 149 (2023). - https://doi.org/10.1140/epjc/s10052-023-11299-z - - Note - ---- - This processor is composed of a factory function that is called using the - `init_args` argument. The input and output waveforms are passed using - `args`. - - Parameters - ---------- - noise_mat - noise matrix - reference - reference signal - length - length of the calculated filter - a1 - penalized coefficient for the noise matrix - a2 - penalized coefficient for the reference matrix - a3 - penalized coefficient for the zero area matrix - ff - flat top length for the reference signal - coefficients - filter coefficients (if given avoid filter calculation) - - - JSON Configuration Example - -------------------------- - - .. code-block :: json - - "wf_dplms": { - "function": "dplms", - "module": "dspeed.processors", - "args": ["wf_diff", "wf_dplms(len(wf_diff)-49, 'f')"], - "unit": "ADC", - "init_args": [ - "db.dplms.noise_matrix", - "db.dplms.reference", - "50", "0.1", "1", "0", "1"] - } - """ - - if length <= 0: - raise DSPFatal("The length of the filter must be positive") - - if len(coefficients) > 0: - if length != len(coefficients): - raise DSPFatal("The length and the provided filter are not matching") - - x = np.array(coefficients) - else: - noise_mat = np.array(noise_mat) - reference = np.array(reference) - - if length != noise_mat.shape[0]: - raise DSPFatal( - "The length of the filter is not consistent with the noise matrix" - ) - - if len(reference) <= 0: - raise DSPFatal("The length of the reference signal must be positive") - - if a1 <= 0: - raise DSPFatal("The penalized coefficient for the noise must be positive") - - if a2 <= 0: - raise DSPFatal( - "The penalized coefficient for the reference must be positive" - ) - - if a3 <= 0: - raise DSPFatal( - "The penalized coefficient for the zero area must be positive" - ) - - if ff <= 0: - raise DSPFatal( - "The penalized coefficient for the ref matrix must be positive" - ) - - # reference matrix - ssize = len(reference) - flo = int(ssize / 2 - length / 2) - fhi = int(ssize / 2 + length / 2) - ref_mat = np.zeros([length, length]) - ref_sig = np.zeros([length]) - if ff == 0: - ff = [0] - elif ff == 1: - ff = [-1, 0, 1] - else: - raise DSPFatal( - "The penalized coefficient for the ref matrix must be 0 or 1" - ) - for i in ff: - ref_mat += np.outer( - reference[flo + i : fhi + i], reference[flo + i : fhi + i] - ) - ref_sig += reference[flo + i : fhi + i] - ref_mat /= len(ff) - - # filter calculation - mat = a1 * noise_mat + a2 * ref_mat + a3 * np.ones([length, length]) - x = np.linalg.solve(mat, ref_sig) - y = np.convolve(reference, np.flip(x), mode="valid") - maxy = np.amax(y) - x /= maxy - y /= maxy - - @guvectorize( - ["void(float32[:], float32[:])", "void(float64[:], float64[:])"], - "(n),(m)", - **nb_kwargs( - cache=False, - forceobj=True, - ), - ) - def dplms_out(w_in: np.ndarray, w_out: np.ndarray) -> None: - """ - Parameters - ---------- - w_in - the input waveform. - w_out - the filtered waveform. - """ - - w_out[:] = np.nan - - if np.isnan(w_in).any(): - return - - if len(x) > len(w_in): - raise DSPFatal("The filter is longer than the input waveform") - - w_out[:] = np.convolve(w_in, np.flip(x), mode="valid") - - return dplms_out diff --git a/src/dspeed/processors/energy_kernels.py b/src/dspeed/processors/energy_kernels.py new file mode 100644 index 0000000..586cf20 --- /dev/null +++ b/src/dspeed/processors/energy_kernels.py @@ -0,0 +1,271 @@ +from __future__ import annotations + +import numpy as np +from numba import guvectorize + +from ..errors import DSPFatal +from ..utils import numba_defaults_kwargs as nb_kwargs + +@guvectorize( + ["void(float32, float32, float32, float32[:])", + "void(float64, float64, float64, float64[:])"], + "(),(),(),(n)", + **nb_kwargs( + cache=False, + forceobj=True, + ), + ) + +def cusp_filter(sigma: float, flat: int, decay: int, kernel: np.array) -> None: + """Calculates CUSP kernel. + + Parameters + ---------- + sigma + the curvature of the rising and falling part of the kernel. + flat + the length of the flat section. + decay + the decay constant of the exponential to be convolved. + kernel + the calculated kernel + + JSON Configuration Example + -------------------------- + + .. code-block :: json + + "kern_cusp": { + "function": "cusp_filter", + "module": "dspeed.processors", + "args": ["10*us", "3*us", "400*us", "kern_cusp"], + "unit": "ADC" + } + """ + + if sigma < 0: + raise DSPFatal("The curvature parameter must be positive") + + if flat < 0: + raise DSPFatal("The length of the flat section must be positive") + + if np.floor(flat) != flat: + raise DSPFatal("The length of the flat section must be an integer") + + if decay < 0: + raise DSPFatal("The decay constant must be positive") + + lt = int((len(kernel) - flat) / 2) + flat_int = int(flat) + for ind in range(0, lt, 1): + kernel[ind] = float(np.sinh(ind / sigma) / np.sinh(lt / sigma)) + for ind in range(lt, lt + flat_int + 1, 1): + kernel[ind] = 1 + for ind in range(lt + flat_int + 1, len(kernel), 1): + kernel[ind] = float(np.sinh((len(kernel) - ind) / sigma) / np.sinh(lt / sigma)) + + den = [1, -np.exp(-1 / decay)] + kernel[:] = np.convolve(kernel, den, "same") + + +@guvectorize( + ["void(float32, float32, float32, float32[:])", + "void(float64, float64, float64, float64[:])"], + "(),(),(),(n)", + **nb_kwargs( + cache=False, + forceobj=True, + ), + ) + +def zac_filter(sigma: float, flat: int, decay: int, kernel: np.array) -> None: + """Calculates ZAC (Zero Area CUSP) kernel. + + Parameters + ---------- + sigma + the curvature of the rising and falling part of the kernel. + flat + the length of the flat section. + decay + the decay constant of the exponential to be convolved. + kernel + the calculated kernel + + JSON Configuration Example + -------------------------- + + .. code-block :: json + + "kern_zac": { + "function": "zac_filter", + "module": "dspeed.processors", + "args": ["10*us", "3*us", "400*us", "kern_zac"], + "unit": "ADC" + } + """ + + if sigma < 0: + raise DSPFatal("The curvature parameter must be positive") + + if flat < 0: + raise DSPFatal("The length of the flat section must be positive") + + if np.floor(flat) != flat: + raise DSPFatal("The length of the flat section must be an integer") + + if decay < 0: + raise DSPFatal("The decay constant must be positive") + + lt = int((len(kernel) - flat) / 2) + flat_int = int(flat) + + # calculate cusp filter and negative parables + par = np.zeros(len(kernel)) + for ind in range(0, lt, 1): + kernel[ind] = float(np.sinh(ind / sigma) / np.sinh(lt / sigma)) + kernel[ind] = np.power(ind - lt / 2, 2) - np.power(lt / 2, 2) + for ind in range(lt, lt + flat_int + 1, 1): + kernel[ind] = 1 + for ind in range(lt + flat_int + 1, len(kernel), 1): + kernel[ind] = float(np.sinh((len(kernel) - ind) / sigma) / np.sinh(lt / sigma)) + par[ind] = np.power(len(kernel) - ind - lt / 2, 2) - np.power(lt / 2, 2) + + # calculate area of cusp and parables + areapar, areacusp = 0, 0 + for i in range(0, len(kernel), 1): + areapar += par[i] + areacusp += kernel[i] + + # normalize parables area + par = -par / areapar * areacusp + + # create zac filter + kernel += par + + # deconvolve zac filter + den = [1, -np.exp(-1 / decay)] + kernel = np.convolve(kernel, den, "same") + + +@guvectorize( + ["void(float32[:,:], float32[:], float32, float32, float32, float32, float32[:])", + "void(float64[:,:], float64[:], float64, float64, float64, float64, float64[:])"], + "(n,n),(m),(),(),(),(),(n)", + **nb_kwargs( + cache=False, + forceobj=True, + ), + ) + +def dplms( + noise_mat: list, + reference: list, + a1: float, + a2: float, + a3: float, + ff: int, + kernel:np.array +) -> None: + """Calculate and apply an optimum DPLMS filter to the waveform. + + The processor takes the noise matrix and the reference signal as input and + calculates the optimum filter according to the provided length and + penalized coefficients [DPLMS]_. + + .. [DPLMS] V. D'Andrea et al. “Optimum Filter Synthesis with DPLMS + Method for Energy Reconstruction” Eur. Phys. J. C 83, 149 (2023). + https://doi.org/10.1140/epjc/s10052-023-11299-z + + + Parameters + ---------- + noise_mat + noise matrix + reference + reference signal + a1 + penalized coefficient for the noise matrix + a2 + penalized coefficient for the reference matrix + a3 + penalized coefficient for the zero area matrix + ff + flat top length for the reference signal + kernel + output kernel + + + JSON Configuration Example + -------------------------- + + .. code-block :: json + + "kern_dplms": { + "function": "dplms", + "module": "dspeed.processors", + "args": ["db.dplms.noise_matrix", + "db.dplms.reference", + "50", "0.1", "1", "1" + "kern_dplms"], + "unit": "ADC", + } + """ + + noise_mat = np.array(noise_mat) + reference = np.array(reference) + + if len(kernel) != noise_mat.shape[0]: + raise DSPFatal( + "The length of the filter is not consistent with the noise matrix" + ) + + if len(reference) <= 0: + raise DSPFatal("The length of the reference signal must be positive") + + if a1 <= 0: + raise DSPFatal("The penalized coefficient for the noise must be positive") + + if a2 <= 0: + raise DSPFatal( + "The penalized coefficient for the reference must be positive" + ) + + if a3 <= 0: + raise DSPFatal( + "The penalized coefficient for the zero area must be positive" + ) + + if ff <= 0: + raise DSPFatal( + "The penalized coefficient for the ref matrix must be positive" + ) + + # reference matrix + length = len(kernel) + ssize = len(reference) + flo = int(ssize / 2 - length / 2) + fhi = int(ssize / 2 + length / 2) + ref_mat = np.zeros([length, length]) + ref_sig = np.zeros([length]) + if ff == 0: + ff = [0] + elif ff == 1: + ff = [-1, 0, 1] + else: + raise DSPFatal( + "The penalized coefficient for the ref matrix must be 0 or 1" + ) + for i in ff: + ref_mat += np.outer( + reference[flo + i : fhi + i], reference[flo + i : fhi + i] + ) + ref_sig += reference[flo + i : fhi + i] + ref_mat /= len(ff) + + # filter calculation + mat = a1 * noise_mat + a2 * ref_mat + a3 * np.ones([length, length]) + kernel[:] = np.flip(np.linalg.solve(mat, ref_sig)) + y = np.convolve(reference,kernel, mode="valid") + maxy = np.amax(y) + kernel[:] /= maxy \ No newline at end of file diff --git a/src/dspeed/processors/kernels.py b/src/dspeed/processors/kernels.py new file mode 100644 index 0000000..42c9268 --- /dev/null +++ b/src/dspeed/processors/kernels.py @@ -0,0 +1,137 @@ +from __future__ import annotations + +import numpy as np +from numba import guvectorize + +from ..errors import DSPFatal +from ..utils import numba_defaults_kwargs as nb_kwargs + +@guvectorize( + ["void(float32, float32, float32[:])", + "void(float64, float64, float64[:])"], + "(),(),(n)", + **nb_kwargs( + cache=False, + forceobj=True, + ), +) +def t0_filter(rise: int, fall: int, kernel:np.array) -> None: + """Apply a modified, asymmetric trapezoidal filter to the waveform. + + Parameters + ---------- + rise + the length of the rise section. This is the linearly increasing + section of the filter that performs a weighted average. + fall + the length of the fall section. This is the simple averaging part + of the filter. + kernel + the output kernel + + JSON Configuration Example + -------------------------- + + .. code-block :: json + + "t0_filter": { + "function": "t0_filter", + "module": "dspeed.processors", + "args": ["128*ns", "2*us", "t0_filter"], + "unit": "ADC", + "init_args": ["128*ns", "2*us"] + } + """ + if rise < 0: + raise DSPFatal("The length of the rise section must be positive") + + if fall < 0: + raise DSPFatal("The length of the fall section must be positive") + + for i in range(int(rise)): + kernel[i] = 2 * (int(rise) - i) / (rise**2) + for i in range(int(rise), len(kernel), 1): + kernel[i] = -1 / fall + +@guvectorize( + ["void(float32, float32[:])", + "void(float64, float64[:])"], + "(),(n)", + **nb_kwargs( + cache=False, + forceobj=True, + ), +) +def moving_slope(length, kernel): + """Calculates the linear slope of kernel + + Parameters + ---------- + length + the length of the section to calculate slope + kernel + the output kernel + + JSON Configuration Example + -------------------------- + + .. code-block :: json + + "kern_slopes": { + "function": "moving_slope", + "module": "dspeed.processors", + "args": ["12", "kern_slopes"], + "unit": "ADC" + } + """ + length = len(kernel) + + sum_x = length * (length + 1) / 2 + sum_x2 = length * (length + 1) * (2 * length + 1) / 6 + + kernel[:] = (np.arange(1, length + 1, 1) * length) - (np.ones(length) * sum_x) + kernel[:] /= length * sum_x2 - sum_x * sum_x + kernel[:] = kernel[::-1] + +@guvectorize( + ["void(float32, float32[:])", + "void(float64, float64[:])"], + "(),(n)", + **nb_kwargs( + cache=False, + forceobj=True, + ), + ) +def step(weight_pos: int, kernel:np.array) -> None: + """Process waveforms with a step function. + + Parameters + ---------- + weight_pos + relative weight of positive step side. + kernel + output kernel + + JSON Configuration Example + -------------------------- + + .. code-block :: json + + "kern_step": { + "function": "step", + "module": "dspeed.processors", + "args": ["16", "kern_step"], + "unit": "ADC" + } + """ + + x = np.arange(len(kernel)) + kernel[:] = np.piecewise( + x, + [ + ((x >= 0) & (x < len(kernel) / 4)), + ((x >= len(kernel) / 4) & (x < 3 * len(kernel) / 4)), + ((x >= 3 * len(kernel) / 4) & (x < len(kernel))), + ], + [-1, 1, -1], + ) From b6d818862166d4cd46e4e546faa5a5664bce6026 Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Mon, 4 Dec 2023 22:04:16 +0000 Subject: [PATCH 05/42] style: pre-commit fixes --- src/dspeed/processors/__init__.py | 6 +- src/dspeed/processors/convolutions.py | 67 +++++++++++--------- src/dspeed/processors/energy_kernels.py | 84 ++++++++++++------------- src/dspeed/processors/kernels.py | 28 ++++----- 4 files changed, 93 insertions(+), 92 deletions(-) diff --git a/src/dspeed/processors/__init__.py b/src/dspeed/processors/__init__.py index 4757214..ace9276 100644 --- a/src/dspeed/processors/__init__.py +++ b/src/dspeed/processors/__init__.py @@ -62,14 +62,14 @@ from .bl_subtract import bl_subtract from .convolutions import convolve_wf, fft_convolve_wf from .dwt import discrete_wavelet_transform -from .energy_kernels import cusp_filter, zac_filter, dplms +from .energy_kernels import cusp_filter, dplms, zac_filter from .fftw import dft, inv_dft, psd from .fixed_time_pickoff import fixed_time_pickoff from .gaussian_filter1d import gaussian_filter1d from .get_multi_local_extrema import get_multi_local_extrema from .get_wf_centroid import get_wf_centroid from .histogram import histogram, histogram_stats -from .kernels import t0_filter, moving_slope, step +from .kernels import moving_slope, step, t0_filter from .linear_slope_fit import linear_slope_diff, linear_slope_fit from .log_check import log_check from .min_max import min_max @@ -101,7 +101,7 @@ __all__ = [ "bl_subtract", - "convolve_wf", + "convolve_wf", "fft_convolve_wf", "cusp_filter", "t0_filter", diff --git a/src/dspeed/processors/convolutions.py b/src/dspeed/processors/convolutions.py index f5ce486..f49f893 100644 --- a/src/dspeed/processors/convolutions.py +++ b/src/dspeed/processors/convolutions.py @@ -7,16 +7,21 @@ from ..errors import DSPFatal from ..utils import numba_defaults_kwargs as nb_kwargs + @guvectorize( - ["void(float32[:], float32[:], char, float32[:])", - "void(float64[:], float64[:], char, float64[:])"], - "(n),(m),(),(p)", - **nb_kwargs( - cache=False, - forceobj=True, - ), - ) -def convolve_wf(w_in: np.ndarray, kernel:np.array, mode_in: np.int8, w_out: np.ndarray) -> None: # + [ + "void(float32[:], float32[:], char, float32[:])", + "void(float64[:], float64[:], char, float64[:])", + ], + "(n),(m),(),(p)", + **nb_kwargs( + cache=False, + forceobj=True, + ), +) +def convolve_wf( + w_in: np.ndarray, kernel: np.array, mode_in: np.int8, w_out: np.ndarray +) -> None: # """ Parameters ---------- @@ -24,8 +29,8 @@ def convolve_wf(w_in: np.ndarray, kernel:np.array, mode_in: np.int8, w_out: np.n the input waveform. kernel the kernel to convolve with - mode - mode of convolution options are f : full, v : valid or s : same, + mode + mode of convolution options are f : full, v : valid or s : same, explained here: https://numpy.org/doc/stable/reference/generated/numpy.convolve.html w_out the filtered waveform. @@ -34,29 +39,33 @@ def convolve_wf(w_in: np.ndarray, kernel:np.array, mode_in: np.int8, w_out: np.n if np.isnan(w_in).any(): return - + if np.isnan(kernel).any(): return if len(kernel) > len(w_in): raise DSPFatal("The filter is longer than the input waveform") - + if chr(mode_in) not in ["f", "v", "s"]: raise DSPFatal("Invalid mode") - + w_out[:] = np.convolve(w_in, kernel, mode=chr(mode_in)) @guvectorize( - ["void(float32[:], float32[:], char, float32[:])", - "void(float64[:], float64[:], char, float64[:])"], - "(n),(m),(),(p)", - **nb_kwargs( - cache=False, - forceobj=True, - ), - ) -def fft_convolve_wf(w_in: np.ndarray, kernel:np.array, mode_in: np.int8, w_out: np.ndarray) -> None: # + [ + "void(float32[:], float32[:], char, float32[:])", + "void(float64[:], float64[:], char, float64[:])", + ], + "(n),(m),(),(p)", + **nb_kwargs( + cache=False, + forceobj=True, + ), +) +def fft_convolve_wf( + w_in: np.ndarray, kernel: np.array, mode_in: np.int8, w_out: np.ndarray +) -> None: # """ Parameters ---------- @@ -64,8 +73,8 @@ def fft_convolve_wf(w_in: np.ndarray, kernel:np.array, mode_in: np.int8, w_out: the input waveform. kernel the kernel to convolve with - mode - mode of convolution options are f : full, v : valid or s : same, + mode + mode of convolution options are f : full, v : valid or s : same, explained here: https://numpy.org/doc/stable/reference/generated/numpy.convolve.html w_out the filtered waveform. @@ -74,13 +83,13 @@ def fft_convolve_wf(w_in: np.ndarray, kernel:np.array, mode_in: np.int8, w_out: if np.isnan(w_in).any(): return - + if np.isnan(kernel).any(): return if len(kernel) > len(w_in): raise DSPFatal("The filter is longer than the input waveform") - + if chr(mode_in) == "f": mode = "full" elif chr(mode_in) == "v": @@ -89,7 +98,5 @@ def fft_convolve_wf(w_in: np.ndarray, kernel:np.array, mode_in: np.int8, w_out: mode = "same" else: raise DSPFatal("Invalid mode") - - w_out[:] = fftconvolve(w_in, kernel, mode=mode) - + w_out[:] = fftconvolve(w_in, kernel, mode=mode) diff --git a/src/dspeed/processors/energy_kernels.py b/src/dspeed/processors/energy_kernels.py index 586cf20..6d2e8ec 100644 --- a/src/dspeed/processors/energy_kernels.py +++ b/src/dspeed/processors/energy_kernels.py @@ -6,16 +6,18 @@ from ..errors import DSPFatal from ..utils import numba_defaults_kwargs as nb_kwargs -@guvectorize( - ["void(float32, float32, float32, float32[:])", - "void(float64, float64, float64, float64[:])"], - "(),(),(),(n)", - **nb_kwargs( - cache=False, - forceobj=True, - ), - ) +@guvectorize( + [ + "void(float32, float32, float32, float32[:])", + "void(float64, float64, float64, float64[:])", + ], + "(),(),(),(n)", + **nb_kwargs( + cache=False, + forceobj=True, + ), +) def cusp_filter(sigma: float, flat: int, decay: int, kernel: np.array) -> None: """Calculates CUSP kernel. @@ -69,15 +71,16 @@ def cusp_filter(sigma: float, flat: int, decay: int, kernel: np.array) -> None: @guvectorize( - ["void(float32, float32, float32, float32[:])", - "void(float64, float64, float64, float64[:])"], - "(),(),(),(n)", - **nb_kwargs( - cache=False, - forceobj=True, - ), - ) - + [ + "void(float32, float32, float32, float32[:])", + "void(float64, float64, float64, float64[:])", + ], + "(),(),(),(n)", + **nb_kwargs( + cache=False, + forceobj=True, + ), +) def zac_filter(sigma: float, flat: int, decay: int, kernel: np.array) -> None: """Calculates ZAC (Zero Area CUSP) kernel. @@ -149,15 +152,16 @@ def zac_filter(sigma: float, flat: int, decay: int, kernel: np.array) -> None: @guvectorize( - ["void(float32[:,:], float32[:], float32, float32, float32, float32, float32[:])", - "void(float64[:,:], float64[:], float64, float64, float64, float64, float64[:])"], - "(n,n),(m),(),(),(),(),(n)", - **nb_kwargs( - cache=False, - forceobj=True, - ), - ) - + [ + "void(float32[:,:], float32[:], float32, float32, float32, float32, float32[:])", + "void(float64[:,:], float64[:], float64, float64, float64, float64, float64[:])", + ], + "(n,n),(m),(),(),(),(),(n)", + **nb_kwargs( + cache=False, + forceobj=True, + ), +) def dplms( noise_mat: list, reference: list, @@ -165,7 +169,7 @@ def dplms( a2: float, a3: float, ff: int, - kernel:np.array + kernel: np.array, ) -> None: """Calculate and apply an optimum DPLMS filter to the waveform. @@ -227,19 +231,13 @@ def dplms( raise DSPFatal("The penalized coefficient for the noise must be positive") if a2 <= 0: - raise DSPFatal( - "The penalized coefficient for the reference must be positive" - ) + raise DSPFatal("The penalized coefficient for the reference must be positive") if a3 <= 0: - raise DSPFatal( - "The penalized coefficient for the zero area must be positive" - ) + raise DSPFatal("The penalized coefficient for the zero area must be positive") if ff <= 0: - raise DSPFatal( - "The penalized coefficient for the ref matrix must be positive" - ) + raise DSPFatal("The penalized coefficient for the ref matrix must be positive") # reference matrix length = len(kernel) @@ -253,19 +251,15 @@ def dplms( elif ff == 1: ff = [-1, 0, 1] else: - raise DSPFatal( - "The penalized coefficient for the ref matrix must be 0 or 1" - ) + raise DSPFatal("The penalized coefficient for the ref matrix must be 0 or 1") for i in ff: - ref_mat += np.outer( - reference[flo + i : fhi + i], reference[flo + i : fhi + i] - ) + ref_mat += np.outer(reference[flo + i : fhi + i], reference[flo + i : fhi + i]) ref_sig += reference[flo + i : fhi + i] ref_mat /= len(ff) # filter calculation mat = a1 * noise_mat + a2 * ref_mat + a3 * np.ones([length, length]) kernel[:] = np.flip(np.linalg.solve(mat, ref_sig)) - y = np.convolve(reference,kernel, mode="valid") + y = np.convolve(reference, kernel, mode="valid") maxy = np.amax(y) - kernel[:] /= maxy \ No newline at end of file + kernel[:] /= maxy diff --git a/src/dspeed/processors/kernels.py b/src/dspeed/processors/kernels.py index 42c9268..2a1174c 100644 --- a/src/dspeed/processors/kernels.py +++ b/src/dspeed/processors/kernels.py @@ -6,16 +6,16 @@ from ..errors import DSPFatal from ..utils import numba_defaults_kwargs as nb_kwargs + @guvectorize( - ["void(float32, float32, float32[:])", - "void(float64, float64, float64[:])"], + ["void(float32, float32, float32[:])", "void(float64, float64, float64[:])"], "(),(),(n)", **nb_kwargs( cache=False, forceobj=True, ), ) -def t0_filter(rise: int, fall: int, kernel:np.array) -> None: +def t0_filter(rise: int, fall: int, kernel: np.array) -> None: """Apply a modified, asymmetric trapezoidal filter to the waveform. Parameters @@ -53,9 +53,9 @@ def t0_filter(rise: int, fall: int, kernel:np.array) -> None: for i in range(int(rise), len(kernel), 1): kernel[i] = -1 / fall + @guvectorize( - ["void(float32, float32[:])", - "void(float64, float64[:])"], + ["void(float32, float32[:])", "void(float64, float64[:])"], "(),(n)", **nb_kwargs( cache=False, @@ -93,16 +93,16 @@ def moving_slope(length, kernel): kernel[:] /= length * sum_x2 - sum_x * sum_x kernel[:] = kernel[::-1] + @guvectorize( - ["void(float32, float32[:])", - "void(float64, float64[:])"], - "(),(n)", - **nb_kwargs( - cache=False, - forceobj=True, - ), - ) -def step(weight_pos: int, kernel:np.array) -> None: + ["void(float32, float32[:])", "void(float64, float64[:])"], + "(),(n)", + **nb_kwargs( + cache=False, + forceobj=True, + ), +) +def step(weight_pos: int, kernel: np.array) -> None: """Process waveforms with a step function. Parameters From f702b92060816ad1f6e00b57cacc817f177d6025 Mon Sep 17 00:00:00 2001 From: iguinn Date: Tue, 5 Dec 2023 09:58:26 -0800 Subject: [PATCH 06/42] lh5_groups option in build_dsp now works with wildcards --- src/dspeed/build_dsp.py | 4 +++- 1 file changed, 3 insertions(+), 1 deletion(-) diff --git a/src/dspeed/build_dsp.py b/src/dspeed/build_dsp.py index aba0581..d834d98 100644 --- a/src/dspeed/build_dsp.py +++ b/src/dspeed/build_dsp.py @@ -107,7 +107,9 @@ def build_dsp( if lh5_tables is None: lh5_tables = lh5.ls(f_raw) elif isinstance(lh5_tables, str): - lh5_tables = [lh5_tables] + lh5_tables = lh5.ls(f_raw, lh5_tables) + elif isinstance(lh5_tables, list): + lh5_tables = [tab for tab_wc in lh5_tables for tab in lh5.ls(f_raw, tab_wc)] elif not ( hasattr(lh5_tables, "__iter__") and all(isinstance(el, str) for el in lh5_tables) From 10083c0001736c3e100a77982955827da7761f30 Mon Sep 17 00:00:00 2001 From: iguinn Date: Tue, 5 Dec 2023 09:58:57 -0800 Subject: [PATCH 07/42] Fixed stripped units warning --- src/dspeed/processing_chain.py | 9 ++++++--- 1 file changed, 6 insertions(+), 3 deletions(-) diff --git a/src/dspeed/processing_chain.py b/src/dspeed/processing_chain.py index db139b9..6cb37fd 100644 --- a/src/dspeed/processing_chain.py +++ b/src/dspeed/processing_chain.py @@ -693,7 +693,10 @@ def _parse_expr( op, op_form = ast_ops_dict[type(node.op)] if not (isinstance(lhs, ProcChainVar) or isinstance(rhs, ProcChainVar)): - return op(lhs, rhs) + ret = op(lhs, rhs) + if isinstance(ret, Quantity) and ureg.is_compatible_with(ret.u, ureg.dimensionless): + ret = ret.to(ureg.dimensionless).magnitude + return ret name = "(" + op_form.format(str(lhs), str(rhs)) + ")" if isinstance(lhs, ProcChainVar) and isinstance(rhs, ProcChainVar): @@ -1236,7 +1239,7 @@ def __init__( # Convert scalar to right type, including units if isinstance(param, (Quantity, Unit)): if ureg.is_compatible_with(ureg.dimensionless, param): - param = float(param) + param = param.to(ureg.dimensionless).magnitude elif not isinstance( grid, CoordinateGrid ) or not ureg.is_compatible_with(grid.period, param): @@ -1245,7 +1248,7 @@ def __init__( f"CoordinateGrid is {grid}" ) else: - param = float(param / grid.period) + param = (param / grid.period).to(ureg.dimensionless).magnitude if np.issubdtype(dtype, np.integer): param = dtype.type(round(param)) else: From b006c7857113a98e2b876714720d7c92bb7ee383 Mon Sep 17 00:00:00 2001 From: ggmarshall Date: Wed, 6 Dec 2023 01:10:08 +0100 Subject: [PATCH 08/42] updated concolve to pass full string as mode, removed uneeded length from moving slope --- src/dspeed/processors/convolutions.py | 10 ++++++++-- src/dspeed/processors/kernels.py | 2 +- 2 files changed, 9 insertions(+), 3 deletions(-) diff --git a/src/dspeed/processors/convolutions.py b/src/dspeed/processors/convolutions.py index f5ce486..b32055d 100644 --- a/src/dspeed/processors/convolutions.py +++ b/src/dspeed/processors/convolutions.py @@ -41,10 +41,16 @@ def convolve_wf(w_in: np.ndarray, kernel:np.array, mode_in: np.int8, w_out: np.n if len(kernel) > len(w_in): raise DSPFatal("The filter is longer than the input waveform") - if chr(mode_in) not in ["f", "v", "s"]: + if chr(mode_in) == "f": + mode = "full" + elif chr(mode_in) == "v": + mode = "valid" + elif chr(mode_in) == "s": + mode = "same" + else: raise DSPFatal("Invalid mode") - w_out[:] = np.convolve(w_in, kernel, mode=chr(mode_in)) + w_out[:] = np.convolve(w_in, kernel, mode=mode) @guvectorize( diff --git a/src/dspeed/processors/kernels.py b/src/dspeed/processors/kernels.py index 42c9268..b0c1072 100644 --- a/src/dspeed/processors/kernels.py +++ b/src/dspeed/processors/kernels.py @@ -62,7 +62,7 @@ def t0_filter(rise: int, fall: int, kernel:np.array) -> None: forceobj=True, ), ) -def moving_slope(length, kernel): +def moving_slope(kernel): """Calculates the linear slope of kernel Parameters From 24658b5405997352a8265c2579f2188f32802bab Mon Sep 17 00:00:00 2001 From: ggmarshall Date: Wed, 6 Dec 2023 01:11:20 +0100 Subject: [PATCH 09/42] merge --- src/dspeed/processors/__init__.py | 6 +- src/dspeed/processors/convolutions.py | 63 ++++++++++--------- src/dspeed/processors/energy_kernels.py | 84 ++++++++++++------------- src/dspeed/processors/kernels.py | 28 ++++----- 4 files changed, 91 insertions(+), 90 deletions(-) diff --git a/src/dspeed/processors/__init__.py b/src/dspeed/processors/__init__.py index 4757214..ace9276 100644 --- a/src/dspeed/processors/__init__.py +++ b/src/dspeed/processors/__init__.py @@ -62,14 +62,14 @@ from .bl_subtract import bl_subtract from .convolutions import convolve_wf, fft_convolve_wf from .dwt import discrete_wavelet_transform -from .energy_kernels import cusp_filter, zac_filter, dplms +from .energy_kernels import cusp_filter, dplms, zac_filter from .fftw import dft, inv_dft, psd from .fixed_time_pickoff import fixed_time_pickoff from .gaussian_filter1d import gaussian_filter1d from .get_multi_local_extrema import get_multi_local_extrema from .get_wf_centroid import get_wf_centroid from .histogram import histogram, histogram_stats -from .kernels import t0_filter, moving_slope, step +from .kernels import moving_slope, step, t0_filter from .linear_slope_fit import linear_slope_diff, linear_slope_fit from .log_check import log_check from .min_max import min_max @@ -101,7 +101,7 @@ __all__ = [ "bl_subtract", - "convolve_wf", + "convolve_wf", "fft_convolve_wf", "cusp_filter", "t0_filter", diff --git a/src/dspeed/processors/convolutions.py b/src/dspeed/processors/convolutions.py index b32055d..66f0bdd 100644 --- a/src/dspeed/processors/convolutions.py +++ b/src/dspeed/processors/convolutions.py @@ -7,16 +7,21 @@ from ..errors import DSPFatal from ..utils import numba_defaults_kwargs as nb_kwargs + @guvectorize( - ["void(float32[:], float32[:], char, float32[:])", - "void(float64[:], float64[:], char, float64[:])"], - "(n),(m),(),(p)", - **nb_kwargs( - cache=False, - forceobj=True, - ), - ) -def convolve_wf(w_in: np.ndarray, kernel:np.array, mode_in: np.int8, w_out: np.ndarray) -> None: # + [ + "void(float32[:], float32[:], char, float32[:])", + "void(float64[:], float64[:], char, float64[:])", + ], + "(n),(m),(),(p)", + **nb_kwargs( + cache=False, + forceobj=True, + ), +) +def convolve_wf( + w_in: np.ndarray, kernel: np.array, mode_in: np.int8, w_out: np.ndarray +) -> None: # """ Parameters ---------- @@ -24,8 +29,8 @@ def convolve_wf(w_in: np.ndarray, kernel:np.array, mode_in: np.int8, w_out: np.n the input waveform. kernel the kernel to convolve with - mode - mode of convolution options are f : full, v : valid or s : same, + mode + mode of convolution options are f : full, v : valid or s : same, explained here: https://numpy.org/doc/stable/reference/generated/numpy.convolve.html w_out the filtered waveform. @@ -34,7 +39,7 @@ def convolve_wf(w_in: np.ndarray, kernel:np.array, mode_in: np.int8, w_out: np.n if np.isnan(w_in).any(): return - + if np.isnan(kernel).any(): return @@ -54,15 +59,19 @@ def convolve_wf(w_in: np.ndarray, kernel:np.array, mode_in: np.int8, w_out: np.n @guvectorize( - ["void(float32[:], float32[:], char, float32[:])", - "void(float64[:], float64[:], char, float64[:])"], - "(n),(m),(),(p)", - **nb_kwargs( - cache=False, - forceobj=True, - ), - ) -def fft_convolve_wf(w_in: np.ndarray, kernel:np.array, mode_in: np.int8, w_out: np.ndarray) -> None: # + [ + "void(float32[:], float32[:], char, float32[:])", + "void(float64[:], float64[:], char, float64[:])", + ], + "(n),(m),(),(p)", + **nb_kwargs( + cache=False, + forceobj=True, + ), +) +def fft_convolve_wf( + w_in: np.ndarray, kernel: np.array, mode_in: np.int8, w_out: np.ndarray +) -> None: # """ Parameters ---------- @@ -70,8 +79,8 @@ def fft_convolve_wf(w_in: np.ndarray, kernel:np.array, mode_in: np.int8, w_out: the input waveform. kernel the kernel to convolve with - mode - mode of convolution options are f : full, v : valid or s : same, + mode + mode of convolution options are f : full, v : valid or s : same, explained here: https://numpy.org/doc/stable/reference/generated/numpy.convolve.html w_out the filtered waveform. @@ -80,13 +89,13 @@ def fft_convolve_wf(w_in: np.ndarray, kernel:np.array, mode_in: np.int8, w_out: if np.isnan(w_in).any(): return - + if np.isnan(kernel).any(): return if len(kernel) > len(w_in): raise DSPFatal("The filter is longer than the input waveform") - + if chr(mode_in) == "f": mode = "full" elif chr(mode_in) == "v": @@ -95,7 +104,5 @@ def fft_convolve_wf(w_in: np.ndarray, kernel:np.array, mode_in: np.int8, w_out: mode = "same" else: raise DSPFatal("Invalid mode") - - w_out[:] = fftconvolve(w_in, kernel, mode=mode) - + w_out[:] = fftconvolve(w_in, kernel, mode=mode) diff --git a/src/dspeed/processors/energy_kernels.py b/src/dspeed/processors/energy_kernels.py index 586cf20..6d2e8ec 100644 --- a/src/dspeed/processors/energy_kernels.py +++ b/src/dspeed/processors/energy_kernels.py @@ -6,16 +6,18 @@ from ..errors import DSPFatal from ..utils import numba_defaults_kwargs as nb_kwargs -@guvectorize( - ["void(float32, float32, float32, float32[:])", - "void(float64, float64, float64, float64[:])"], - "(),(),(),(n)", - **nb_kwargs( - cache=False, - forceobj=True, - ), - ) +@guvectorize( + [ + "void(float32, float32, float32, float32[:])", + "void(float64, float64, float64, float64[:])", + ], + "(),(),(),(n)", + **nb_kwargs( + cache=False, + forceobj=True, + ), +) def cusp_filter(sigma: float, flat: int, decay: int, kernel: np.array) -> None: """Calculates CUSP kernel. @@ -69,15 +71,16 @@ def cusp_filter(sigma: float, flat: int, decay: int, kernel: np.array) -> None: @guvectorize( - ["void(float32, float32, float32, float32[:])", - "void(float64, float64, float64, float64[:])"], - "(),(),(),(n)", - **nb_kwargs( - cache=False, - forceobj=True, - ), - ) - + [ + "void(float32, float32, float32, float32[:])", + "void(float64, float64, float64, float64[:])", + ], + "(),(),(),(n)", + **nb_kwargs( + cache=False, + forceobj=True, + ), +) def zac_filter(sigma: float, flat: int, decay: int, kernel: np.array) -> None: """Calculates ZAC (Zero Area CUSP) kernel. @@ -149,15 +152,16 @@ def zac_filter(sigma: float, flat: int, decay: int, kernel: np.array) -> None: @guvectorize( - ["void(float32[:,:], float32[:], float32, float32, float32, float32, float32[:])", - "void(float64[:,:], float64[:], float64, float64, float64, float64, float64[:])"], - "(n,n),(m),(),(),(),(),(n)", - **nb_kwargs( - cache=False, - forceobj=True, - ), - ) - + [ + "void(float32[:,:], float32[:], float32, float32, float32, float32, float32[:])", + "void(float64[:,:], float64[:], float64, float64, float64, float64, float64[:])", + ], + "(n,n),(m),(),(),(),(),(n)", + **nb_kwargs( + cache=False, + forceobj=True, + ), +) def dplms( noise_mat: list, reference: list, @@ -165,7 +169,7 @@ def dplms( a2: float, a3: float, ff: int, - kernel:np.array + kernel: np.array, ) -> None: """Calculate and apply an optimum DPLMS filter to the waveform. @@ -227,19 +231,13 @@ def dplms( raise DSPFatal("The penalized coefficient for the noise must be positive") if a2 <= 0: - raise DSPFatal( - "The penalized coefficient for the reference must be positive" - ) + raise DSPFatal("The penalized coefficient for the reference must be positive") if a3 <= 0: - raise DSPFatal( - "The penalized coefficient for the zero area must be positive" - ) + raise DSPFatal("The penalized coefficient for the zero area must be positive") if ff <= 0: - raise DSPFatal( - "The penalized coefficient for the ref matrix must be positive" - ) + raise DSPFatal("The penalized coefficient for the ref matrix must be positive") # reference matrix length = len(kernel) @@ -253,19 +251,15 @@ def dplms( elif ff == 1: ff = [-1, 0, 1] else: - raise DSPFatal( - "The penalized coefficient for the ref matrix must be 0 or 1" - ) + raise DSPFatal("The penalized coefficient for the ref matrix must be 0 or 1") for i in ff: - ref_mat += np.outer( - reference[flo + i : fhi + i], reference[flo + i : fhi + i] - ) + ref_mat += np.outer(reference[flo + i : fhi + i], reference[flo + i : fhi + i]) ref_sig += reference[flo + i : fhi + i] ref_mat /= len(ff) # filter calculation mat = a1 * noise_mat + a2 * ref_mat + a3 * np.ones([length, length]) kernel[:] = np.flip(np.linalg.solve(mat, ref_sig)) - y = np.convolve(reference,kernel, mode="valid") + y = np.convolve(reference, kernel, mode="valid") maxy = np.amax(y) - kernel[:] /= maxy \ No newline at end of file + kernel[:] /= maxy diff --git a/src/dspeed/processors/kernels.py b/src/dspeed/processors/kernels.py index b0c1072..5cb8e76 100644 --- a/src/dspeed/processors/kernels.py +++ b/src/dspeed/processors/kernels.py @@ -6,16 +6,16 @@ from ..errors import DSPFatal from ..utils import numba_defaults_kwargs as nb_kwargs + @guvectorize( - ["void(float32, float32, float32[:])", - "void(float64, float64, float64[:])"], + ["void(float32, float32, float32[:])", "void(float64, float64, float64[:])"], "(),(),(n)", **nb_kwargs( cache=False, forceobj=True, ), ) -def t0_filter(rise: int, fall: int, kernel:np.array) -> None: +def t0_filter(rise: int, fall: int, kernel: np.array) -> None: """Apply a modified, asymmetric trapezoidal filter to the waveform. Parameters @@ -53,9 +53,9 @@ def t0_filter(rise: int, fall: int, kernel:np.array) -> None: for i in range(int(rise), len(kernel), 1): kernel[i] = -1 / fall + @guvectorize( - ["void(float32, float32[:])", - "void(float64, float64[:])"], + ["void(float32, float32[:])", "void(float64, float64[:])"], "(),(n)", **nb_kwargs( cache=False, @@ -93,16 +93,16 @@ def moving_slope(kernel): kernel[:] /= length * sum_x2 - sum_x * sum_x kernel[:] = kernel[::-1] + @guvectorize( - ["void(float32, float32[:])", - "void(float64, float64[:])"], - "(),(n)", - **nb_kwargs( - cache=False, - forceobj=True, - ), - ) -def step(weight_pos: int, kernel:np.array) -> None: + ["void(float32, float32[:])", "void(float64, float64[:])"], + "(),(n)", + **nb_kwargs( + cache=False, + forceobj=True, + ), +) +def step(weight_pos: int, kernel: np.array) -> None: """Process waveforms with a step function. Parameters From 19d3c70a6a4a7fa8153653f3f360f84aaef31363 Mon Sep 17 00:00:00 2001 From: ggmarshall Date: Wed, 6 Dec 2023 01:18:42 +0100 Subject: [PATCH 10/42] fix tests --- src/dspeed/processors/convolutions.py | 8 --- tests/configs/icpc-dsp-config.json | 91 ++++++++++++++++++--------- tests/configs/sipm-dplms-config.json | 33 ++++++---- tests/processors/test_dplms.py | 33 +++------- 4 files changed, 91 insertions(+), 74 deletions(-) diff --git a/src/dspeed/processors/convolutions.py b/src/dspeed/processors/convolutions.py index cca64cc..66f0bdd 100644 --- a/src/dspeed/processors/convolutions.py +++ b/src/dspeed/processors/convolutions.py @@ -45,7 +45,6 @@ def convolve_wf( if len(kernel) > len(w_in): raise DSPFatal("The filter is longer than the input waveform") -<<<<<<< HEAD if chr(mode_in) == "f": mode = "full" @@ -57,13 +56,6 @@ def convolve_wf( raise DSPFatal("Invalid mode") w_out[:] = np.convolve(w_in, kernel, mode=mode) -======= - - if chr(mode_in) not in ["f", "v", "s"]: - raise DSPFatal("Invalid mode") - - w_out[:] = np.convolve(w_in, kernel, mode=chr(mode_in)) ->>>>>>> b6d818862166d4cd46e4e546faa5a5664bce6026 @guvectorize( diff --git a/tests/configs/icpc-dsp-config.json b/tests/configs/icpc-dsp-config.json index fefd532..8f98101 100644 --- a/tests/configs/icpc-dsp-config.json +++ b/tests/configs/icpc-dsp-config.json @@ -73,11 +73,24 @@ "args": ["wf_pz[1500:]", "pz_mean", "pz_std", "pz_slope", "pz_intercept"], "unit": ["ADC", "ADC", "ADC", "ADC"] }, - "wf_t0_filter": { + "t0_kernel":{ "function": "t0_filter", "module": "dspeed.processors", - "args": ["wf_pz", "wf_t0_filter(len(wf_pz), 'f', grid=wf_pz.grid)"], - "init_args": ["128*ns/wf_pz.period", "2*us/wf_pz.period"], + "args": [ + "128*ns/wf_pz.period", "2*us/wf_pz.period", + "t0_kernel(round((128*ns+2*us)/wf_pz.period), 'f')" + ], + "unit": "ADC" + }, + "wf_t0_filter": { + "function": "convolve_wf", + "module": "dspeed.processors", + "args": [ + "wf_pz", + "t0_kernel", + "'s'", + "wf_t0_filter(len(wf_pz), 'f', grid=wf_pz.grid)" + ], "unit": "ADC" }, "wf_atrap": { @@ -148,23 +161,33 @@ "db.etrap.sample": "0.8" } }, - "wf_cusp": { + "cusp_kernel":{ "function": "cusp_filter", "module": "dspeed.processors", - "args": ["wf_blsub", "wf_cusp(101, 'f')"], - "init_args": [ - "len(wf_blsub)-100", - "db.cusp.sigma/wf_blsub.period", - "round(db.cusp.flat/wf_blsub.period)", - "db.pz.tau" + "args":[ + "db.cusp.sigma/wf_blsub.period", + "round(db.cusp.flat/wf_blsub.period)", + "db.pz.tau/wf_blsub.period", + "cusp_kernel(round(len(wf_blsub)-(33.6*us/wf_blsub.period)-(4.8*us/wf_blsub.period)), 'f')" ], "defaults": { - "db.cusp.sigma": "20*us", - "db.cusp.flat": "3*us", - "db.pz.tau": "27460.5" - }, - "unit": "ADC" - }, + "db.cusp.sigma": "20*us", + "db.cusp.flat": "3*us", + "db.pz.tau": "450*us" + }, + "unit": "ADC" + }, + "wf_cusp": { + "function": "fft_convolve_wf", + "module": "dspeed.processors", + "args": [ + "wf_blsub[:round(len(wf_blsub)-(33.6*us/wf_blsub.period))]", + "cusp_kernel", + "'v'", + "wf_cusp(round((4.8*us/wf_blsub.period)+1), 'f')" + ], + "unit": "ADC" + }, "cuspEmax": { "function": "amax", "module": "numpy", @@ -179,23 +202,33 @@ "unit": "ADC", "defaults": { "db.cusp.sample": "50" } }, - "wf_zac": { + "zac_kernel":{ "function": "zac_filter", "module": "dspeed.processors", - "args": ["wf_blsub", "wf_zac(101, 'f')"], - "init_args": [ - "len(wf_blsub)-100", - "db.zac.sigma/wf_blsub.period", - "round(db.zac.flat/wf_blsub.period)", - "db.pz.tau" + "args":[ + "db.zac.sigma/wf_blsub.period", + "round(db.zac.flat/wf_blsub.period)", + "db.pz.tau/wf_blsub.period", + "zac_kernel(round(len(wf_blsub)-(33.6*us/wf_blsub.period)-(4.8*us/wf_blsub.period)), 'f')" ], "defaults": { - "db.zac.sigma": "20*us", - "db.zac.flat": "3*us", - "db.pz.tau": "27460.5" - }, - "unit": "ADC" - }, + "db.zac.sigma": "20*us", + "db.zac.flat": "3*us", + "db.pz.tau": "450*us" + }, + "unit": "ADC" + }, + "wf_zac": { + "function": "fft_convolve_wf", + "module": "dspeed.processors", + "args": [ + "wf_blsub[:round(len(wf_blsub)-(33.6*us/wf_blsub.period))]", + "zac_kernel", + "'v'", + "wf_zac(round((4.8*us/wf_blsub.period)+1), 'f')" + ], + "unit": "ADC" + }, "zacEmax": { "function": "amax", "module": "numpy", diff --git a/tests/configs/sipm-dplms-config.json b/tests/configs/sipm-dplms-config.json index c588c8d..3b1fbad 100644 --- a/tests/configs/sipm-dplms-config.json +++ b/tests/configs/sipm-dplms-config.json @@ -94,20 +94,29 @@ "args": ["waveform", 1, "wf_diff(len(waveform)-1)"], "unit": "ADC" }, - "wf_dplms": { + "dplms_kernel": { "function": "dplms_filter", "module": "dspeed.processors", - "args": ["wf_diff", "wf_dplms(len(wf_diff)-49, 'f')"], - "unit": "ADC", - "init_args": [ - "db.dplms.noise_matrix", - "db.dplms.reference", - "50", - "0.01", - "1", - "0", - "0" - ] + "args": ["db.dplms.noise_matrix", + "db.dplms.reference", + "0.01", + "1", + "0", + "0", + "dplms_kernel(50, 'f')"], + "unit": "ADC" + }, + "wf_dplms": { + "description": "convolve optimised cusp filter", + "function": "convolve_wf", + "module": "dspeed.processors", + "args": [ + "wf_diff", + "dplms_kernel", + "'s'", + "wf_dplms(len(wf_diff)-49, 'f')" + ], + "unit": "ADC" }, "h_weights , h_borders": { "function": "histogram", diff --git a/tests/processors/test_dplms.py b/tests/processors/test_dplms.py index 631e4b9..54eb433 100644 --- a/tests/processors/test_dplms.py +++ b/tests/processors/test_dplms.py @@ -11,42 +11,25 @@ def test_dplms(compare_numba_vs_python): with open(Path(__file__).parent / "dplms_noise_mat.dat") as f: nmat = [[float(num) for num in line.split(" ")] for line in f] - length = 50 + kernel = np.zeros(50) len_wf = 100 ref = np.zeros(len_wf) ref[int(len_wf / 2 - 1) : int(len_wf / 2)] = 1 # ensure the DSPFatal is raised for a negative length with pytest.raises(DSPFatal): - dplms(-1, nmat, ref, 1, 1, 1, 1, []) - - # ensure the DSPFatal is raised for length not equal to noise matrix shape - with pytest.raises(DSPFatal): - dplms(10, nmat, ref, 1, 1, 1, 1, []) + dplms(nmat, [], 1, 1, 1, 1, kernel) # ensure the DSPFatal is raised for negative coefficients with pytest.raises(DSPFatal): - dplms(length, nmat, ref, -1, 1, 1, 1, []) + dplms(nmat, ref, -1, 1, 1, 1, kernel) with pytest.raises(DSPFatal): - dplms(length, nmat, ref, 1, -1, 1, 1, []) + dplms(nmat, ref, 1, -1, 1, 1, kernel) with pytest.raises(DSPFatal): - dplms(length, nmat, ref, 1, 1, -1, 1, []) + dplms(nmat, ref, 1, 1, -1, 1, kernel) with pytest.raises(DSPFatal): - dplms(length, nmat, ref, 1, 1, 1, -1, []) + dplms(nmat, ref, 1, 1, 1, -1, kernel) with pytest.raises(DSPFatal): - dplms(length, nmat, ref, 1, 1, 1, 2, []) - - dplms_func = dplms(length, nmat, ref, 1, 1, 1, 1, []) - - # ensure to have a valid output - w_in = np.ones(len_wf) - w_out = np.empty(len_wf - length + 1) - - assert np.all(compare_numba_vs_python(dplms_func, w_in, w_out)) - - # test if there is a nan in w_in - w_in = np.ones(len_wf) - w_in[4] = np.nan - w_out = np.empty(len_wf - length + 1) + dplms(nmat, ref, 1, 1, 1, 2, kernel) - assert np.all(compare_numba_vs_python(dplms_func, w_in, w_out)) + assert np.all(compare_numba_vs_python(dplms, nmat, ref, 1, 1, 1, 1, kernel)) From cde2654b3f3cfeac5ecaa62f1d9110747b7af9a9 Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Wed, 6 Dec 2023 00:19:07 +0000 Subject: [PATCH 11/42] style: pre-commit fixes --- src/dspeed/processors/convolutions.py | 4 +- tests/configs/icpc-dsp-config.json | 95 ++++++++++++++------------- tests/configs/sipm-dplms-config.json | 22 ++++--- 3 files changed, 62 insertions(+), 59 deletions(-) diff --git a/src/dspeed/processors/convolutions.py b/src/dspeed/processors/convolutions.py index 66f0bdd..89de408 100644 --- a/src/dspeed/processors/convolutions.py +++ b/src/dspeed/processors/convolutions.py @@ -45,7 +45,7 @@ def convolve_wf( if len(kernel) > len(w_in): raise DSPFatal("The filter is longer than the input waveform") - + if chr(mode_in) == "f": mode = "full" elif chr(mode_in) == "v": @@ -54,7 +54,7 @@ def convolve_wf( mode = "same" else: raise DSPFatal("Invalid mode") - + w_out[:] = np.convolve(w_in, kernel, mode=mode) diff --git a/tests/configs/icpc-dsp-config.json b/tests/configs/icpc-dsp-config.json index 8f98101..28af292 100644 --- a/tests/configs/icpc-dsp-config.json +++ b/tests/configs/icpc-dsp-config.json @@ -73,12 +73,13 @@ "args": ["wf_pz[1500:]", "pz_mean", "pz_std", "pz_slope", "pz_intercept"], "unit": ["ADC", "ADC", "ADC", "ADC"] }, - "t0_kernel":{ + "t0_kernel": { "function": "t0_filter", "module": "dspeed.processors", "args": [ - "128*ns/wf_pz.period", "2*us/wf_pz.period", - "t0_kernel(round((128*ns+2*us)/wf_pz.period), 'f')" + "128*ns/wf_pz.period", + "2*us/wf_pz.period", + "t0_kernel(round((128*ns+2*us)/wf_pz.period), 'f')" ], "unit": "ADC" }, @@ -87,8 +88,8 @@ "module": "dspeed.processors", "args": [ "wf_pz", - "t0_kernel", - "'s'", + "t0_kernel", + "'s'", "wf_t0_filter(len(wf_pz), 'f', grid=wf_pz.grid)" ], "unit": "ADC" @@ -161,33 +162,33 @@ "db.etrap.sample": "0.8" } }, - "cusp_kernel":{ + "cusp_kernel": { "function": "cusp_filter", "module": "dspeed.processors", - "args":[ - "db.cusp.sigma/wf_blsub.period", - "round(db.cusp.flat/wf_blsub.period)", - "db.pz.tau/wf_blsub.period", - "cusp_kernel(round(len(wf_blsub)-(33.6*us/wf_blsub.period)-(4.8*us/wf_blsub.period)), 'f')" + "args": [ + "db.cusp.sigma/wf_blsub.period", + "round(db.cusp.flat/wf_blsub.period)", + "db.pz.tau/wf_blsub.period", + "cusp_kernel(round(len(wf_blsub)-(33.6*us/wf_blsub.period)-(4.8*us/wf_blsub.period)), 'f')" ], "defaults": { - "db.cusp.sigma": "20*us", - "db.cusp.flat": "3*us", - "db.pz.tau": "450*us" - }, - "unit": "ADC" - }, - "wf_cusp": { - "function": "fft_convolve_wf", - "module": "dspeed.processors", - "args": [ - "wf_blsub[:round(len(wf_blsub)-(33.6*us/wf_blsub.period))]", + "db.cusp.sigma": "20*us", + "db.cusp.flat": "3*us", + "db.pz.tau": "450*us" + }, + "unit": "ADC" + }, + "wf_cusp": { + "function": "fft_convolve_wf", + "module": "dspeed.processors", + "args": [ + "wf_blsub[:round(len(wf_blsub)-(33.6*us/wf_blsub.period))]", "cusp_kernel", "'v'", - "wf_cusp(round((4.8*us/wf_blsub.period)+1), 'f')" - ], - "unit": "ADC" - }, + "wf_cusp(round((4.8*us/wf_blsub.period)+1), 'f')" + ], + "unit": "ADC" + }, "cuspEmax": { "function": "amax", "module": "numpy", @@ -202,33 +203,33 @@ "unit": "ADC", "defaults": { "db.cusp.sample": "50" } }, - "zac_kernel":{ + "zac_kernel": { "function": "zac_filter", "module": "dspeed.processors", - "args":[ - "db.zac.sigma/wf_blsub.period", - "round(db.zac.flat/wf_blsub.period)", - "db.pz.tau/wf_blsub.period", - "zac_kernel(round(len(wf_blsub)-(33.6*us/wf_blsub.period)-(4.8*us/wf_blsub.period)), 'f')" + "args": [ + "db.zac.sigma/wf_blsub.period", + "round(db.zac.flat/wf_blsub.period)", + "db.pz.tau/wf_blsub.period", + "zac_kernel(round(len(wf_blsub)-(33.6*us/wf_blsub.period)-(4.8*us/wf_blsub.period)), 'f')" ], "defaults": { - "db.zac.sigma": "20*us", - "db.zac.flat": "3*us", - "db.pz.tau": "450*us" - }, - "unit": "ADC" - }, - "wf_zac": { - "function": "fft_convolve_wf", - "module": "dspeed.processors", - "args": [ - "wf_blsub[:round(len(wf_blsub)-(33.6*us/wf_blsub.period))]", + "db.zac.sigma": "20*us", + "db.zac.flat": "3*us", + "db.pz.tau": "450*us" + }, + "unit": "ADC" + }, + "wf_zac": { + "function": "fft_convolve_wf", + "module": "dspeed.processors", + "args": [ + "wf_blsub[:round(len(wf_blsub)-(33.6*us/wf_blsub.period))]", "zac_kernel", "'v'", - "wf_zac(round((4.8*us/wf_blsub.period)+1), 'f')" - ], - "unit": "ADC" - }, + "wf_zac(round((4.8*us/wf_blsub.period)+1), 'f')" + ], + "unit": "ADC" + }, "zacEmax": { "function": "amax", "module": "numpy", diff --git a/tests/configs/sipm-dplms-config.json b/tests/configs/sipm-dplms-config.json index 3b1fbad..dd69bac 100644 --- a/tests/configs/sipm-dplms-config.json +++ b/tests/configs/sipm-dplms-config.json @@ -97,13 +97,15 @@ "dplms_kernel": { "function": "dplms_filter", "module": "dspeed.processors", - "args": ["db.dplms.noise_matrix", - "db.dplms.reference", - "0.01", - "1", - "0", - "0", - "dplms_kernel(50, 'f')"], + "args": [ + "db.dplms.noise_matrix", + "db.dplms.reference", + "0.01", + "1", + "0", + "0", + "dplms_kernel(50, 'f')" + ], "unit": "ADC" }, "wf_dplms": { @@ -112,9 +114,9 @@ "module": "dspeed.processors", "args": [ "wf_diff", - "dplms_kernel", - "'s'", - "wf_dplms(len(wf_diff)-49, 'f')" + "dplms_kernel", + "'s'", + "wf_dplms(len(wf_diff)-49, 'f')" ], "unit": "ADC" }, From c04fe73682a9ce1c809c048eccf7ca46bfec2fef Mon Sep 17 00:00:00 2001 From: ggmarshall Date: Wed, 6 Dec 2023 01:22:29 +0100 Subject: [PATCH 12/42] fix moving slope args --- src/dspeed/processors/kernels.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/dspeed/processors/kernels.py b/src/dspeed/processors/kernels.py index 5cb8e76..7baa552 100644 --- a/src/dspeed/processors/kernels.py +++ b/src/dspeed/processors/kernels.py @@ -55,8 +55,8 @@ def t0_filter(rise: int, fall: int, kernel: np.array) -> None: @guvectorize( - ["void(float32, float32[:])", "void(float64, float64[:])"], - "(),(n)", + ["void(float32[:])", "void(float64[:])"], + "(n)", **nb_kwargs( cache=False, forceobj=True, From e3fcee92f15c8f4b28ca6e77546a0870967960ce Mon Sep 17 00:00:00 2001 From: ggmarshall Date: Wed, 6 Dec 2023 15:10:02 +0100 Subject: [PATCH 13/42] fix zac kernel --- src/dspeed/processors/energy_kernels.py | 24 +++++++++++++----------- 1 file changed, 13 insertions(+), 11 deletions(-) diff --git a/src/dspeed/processors/energy_kernels.py b/src/dspeed/processors/energy_kernels.py index 6d2e8ec..495e31d 100644 --- a/src/dspeed/processors/energy_kernels.py +++ b/src/dspeed/processors/energy_kernels.py @@ -122,33 +122,35 @@ def zac_filter(sigma: float, flat: int, decay: int, kernel: np.array) -> None: lt = int((len(kernel) - flat) / 2) flat_int = int(flat) + length = len(kernel) # calculate cusp filter and negative parables - par = np.zeros(len(kernel)) + cusp = np.zeros(length) + par = np.zeros(length) for ind in range(0, lt, 1): - kernel[ind] = float(np.sinh(ind / sigma) / np.sinh(lt / sigma)) - kernel[ind] = np.power(ind - lt / 2, 2) - np.power(lt / 2, 2) + cusp[ind] = float(np.sinh(ind / sigma) / np.sinh(lt / sigma)) + par[ind] = np.power(ind - lt / 2, 2) - np.power(lt / 2, 2) for ind in range(lt, lt + flat_int + 1, 1): - kernel[ind] = 1 - for ind in range(lt + flat_int + 1, len(kernel), 1): - kernel[ind] = float(np.sinh((len(kernel) - ind) / sigma) / np.sinh(lt / sigma)) - par[ind] = np.power(len(kernel) - ind - lt / 2, 2) - np.power(lt / 2, 2) + cusp[ind] = 1 + for ind in range(lt + flat_int + 1, length, 1): + cusp[ind] = float(np.sinh((length - ind) / sigma) / np.sinh(lt / sigma)) + par[ind] = np.power(length - ind - lt / 2, 2) - np.power(lt / 2, 2) # calculate area of cusp and parables areapar, areacusp = 0, 0 - for i in range(0, len(kernel), 1): + for i in range(0, length, 1): areapar += par[i] - areacusp += kernel[i] + areacusp += cusp[i] # normalize parables area par = -par / areapar * areacusp # create zac filter - kernel += par + zac = cusp + par # deconvolve zac filter den = [1, -np.exp(-1 / decay)] - kernel = np.convolve(kernel, den, "same") + kernel[:] = np.convolve(zac, den, "same") @guvectorize( From 84f74d1c92aa7f349a9f6b2c4d4cc38135251a6e Mon Sep 17 00:00:00 2001 From: iguinn Date: Wed, 6 Dec 2023 06:43:37 -0800 Subject: [PATCH 14/42] Added ability to set const variables --- src/dspeed/processing_chain.py | 64 ++++++++++++++++++++++++++++++---- 1 file changed, 58 insertions(+), 6 deletions(-) diff --git a/src/dspeed/processing_chain.py b/src/dspeed/processing_chain.py index 9ef4594..dee1f9d 100644 --- a/src/dspeed/processing_chain.py +++ b/src/dspeed/processing_chain.py @@ -133,6 +133,7 @@ def __init__( grid: CoordinateGrid = auto, unit: str | Unit = auto, is_coord: bool = auto, + is_const: bool = False ) -> None: """ Parameters @@ -154,6 +155,10 @@ def __init__( is_coord If ``True``, variable represents an array index and can be converted into a unitted number using grid. + is_const + If ``True``, variable is a constant. Variable will be set before + executing, and will not be recomputed. Does not have outer + dimension of size _block_width """ assert isinstance(proc_chain, ProcessingChain) and isinstance(name, str) self.proc_chain = proc_chain @@ -168,6 +173,7 @@ def __init__( self.grid = grid self.unit = unit self.is_coord = is_coord + self.is_const = is_const log.debug(f"added variable: {self.description()}") @@ -203,18 +209,23 @@ def __setattr__(self, name: str, value: Any) -> None: super().__setattr__(name, value) + def _make_buffer(self) -> np.ndarray: + shape = self.shape if self.is_const else (self.proc_chain._block_width,) + self.shape + len = np.product(shape) + # Flattened array, with padding to allow memory alignment + buf = np.zeros(len + 64//self.dtype.itemsize, dtype=self.dtype) + # offset to ensure memory alignment + offset = (64 - buf.ctypes.data)%64//self.dtype.itemsize + return buf[offset:offset+len].reshape(shape) + def get_buffer(self, unit: str | Unit = None) -> np.ndarray: # If buffer needs to be created, do so now if self._buffer is None: if self.shape is auto: raise ProcessingChainError(f"cannot deduce shape of {self.name}") if self.dtype is auto: - raise ProcessingChainError(f"cannot deduce shape of {self.name}") - - # create the buffer so that the array start is aligned in memory on a multiple of 64 bytes - self._buffer = np.zeros( - shape=(self.proc_chain._block_width,) + self.shape, dtype=self.dtype - ) + raise ProcessingChainError(f"cannot deduce dtype of {self.name}") + self._buffer = self._make_buffer() if isinstance(self._buffer, np.ndarray): if self.is_coord is True: @@ -418,6 +429,47 @@ def add_variable( self._vars_dict[name] = var return var + def set_constant( + self, + varname: str, + val: np.ndarray | int | float | Quantity, + dtype: DType = None, + unit: str | Unit | Quantity = None, + ) -> ProcChainVar: + """Make a variable act as a constant and set it to val. + + Parameters + ---------- + varname + name of internal variable to set. If it does not exist, create + it; otherwise, set existing variable to be constant + val + value of constant + dtype + dtype of constant + unit + unit of constant + """ + + param = get_variable(varname) + assert(param.is_constant or param._buffer is None) + param.is_constant = True + + if isinstance(val, Quantity): + unit = val.unit + val = val.magnitude + + val = np.array(val, dtype=dtype) + + param.update_auto( + shape=val.shape, + dtype=val.dtype, + unit=unit, + ) + np.copyto(var.get_buffer(), val, casting='unsafe') + log.debug(f"set constant: {self.description()} = {val}") + return param + def link_input_buffer( self, varname: str, buff: np.ndarray | LGDO = None ) -> np.ndarray | LGDO: From c1b77891448681c1ac845396795bdfb023edc8e4 Mon Sep 17 00:00:00 2001 From: iguinn Date: Wed, 6 Dec 2023 06:44:03 -0800 Subject: [PATCH 15/42] build_processing_chain will automatically make variables that have no dependencies constant --- src/dspeed/processing_chain.py | 36 +++++++++++++++++++++++++++++++--- 1 file changed, 33 insertions(+), 3 deletions(-) diff --git a/src/dspeed/processing_chain.py b/src/dspeed/processing_chain.py index dee1f9d..f901895 100644 --- a/src/dspeed/processing_chain.py +++ b/src/dspeed/processing_chain.py @@ -902,7 +902,10 @@ def get_index(slice_value): return attr # Otherwise this is probably a ProcChainVar - val = self._parse_expr(node.value, expr, dry_run, var_name_list) + # Note that we are excluding this variable from the vars list + # because it does not strictly need to be computed before as a + # prerequisite before accessing its attributes + val = self._parse_expr(node.value, expr, dry_run, []) if val is None: return None return getattr(val, node.attr) @@ -2032,10 +2035,10 @@ def resolve_dependencies( module = importlib.import_module(recipe["module"]) func = getattr(module, recipe["function"]) args = recipe["args"] + new_vars = [k for k in re.split(",| ", proc_par) if k != ""] # Initialize the new variables, if needed if "unit" in recipe: - new_vars = [k for k in re.split(",| ", proc_par) if k != ""] for i, name in enumerate(new_vars): unit = recipe.get("unit", auto) if isinstance(unit, list): @@ -2093,7 +2096,34 @@ def resolve_dependencies( except KeyError: pass - proc_chain.add_processor(func, *args, **kwargs) + # Check if new variables should be treated as constants + if not recipe["prereqs"]: + arg_params = [] + kwarg_params = {} + out_is_arg = False + for arg in args: + if isinstance(arg, str): + arg = proc_chain.get_variable(arg) + if isinstance(arg, dict): + kwarg_params.update(arg) + arg = arg.values()[0] + else: + arg_params.append(arg) + if isinstance(arg, ProcChainVar) and arg.name in new_vars: + out_is_arg = True + arg.is_const = True + #arg = arg.get_buffer() + + if out_is_arg: + proc_man = ProcessorManager(proc_chain, func, arg_params, kwarg_params, kwargs.get('signature', None), kwargs.get('types', None)) + proc_man.execute() + else: + const_val = func(*arg_params, **kwarg_params) + if len(new_vars)==1: const_val = (const_val) + for var, val in zip(new_vars, const_val): + proc_chain.set_constant(var, val) + else: + proc_chain.add_processor(func, *args, **kwargs) except Exception as e: raise ProcessingChainError( "Exception raised while attempting to add processor:\n" From 88d6801aaa28812a8dbfb08a313649a1884a8c82 Mon Sep 17 00:00:00 2001 From: iguinn Date: Wed, 6 Dec 2023 06:55:22 -0800 Subject: [PATCH 16/42] Fixed failed test --- src/dspeed/processing_chain.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/dspeed/processing_chain.py b/src/dspeed/processing_chain.py index f901895..b51f779 100644 --- a/src/dspeed/processing_chain.py +++ b/src/dspeed/processing_chain.py @@ -2106,7 +2106,7 @@ def resolve_dependencies( arg = proc_chain.get_variable(arg) if isinstance(arg, dict): kwarg_params.update(arg) - arg = arg.values()[0] + arg = list(arg.values())[0] else: arg_params.append(arg) if isinstance(arg, ProcChainVar) and arg.name in new_vars: From 8fa8292998abf514f3ffb446c85185054650f2d8 Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Wed, 6 Dec 2023 14:59:42 +0000 Subject: [PATCH 17/42] style: pre-commit fixes --- src/dspeed/processing_chain.py | 44 ++++++++++++++++++++++------------ 1 file changed, 29 insertions(+), 15 deletions(-) diff --git a/src/dspeed/processing_chain.py b/src/dspeed/processing_chain.py index b51f779..1b7d57b 100644 --- a/src/dspeed/processing_chain.py +++ b/src/dspeed/processing_chain.py @@ -133,7 +133,7 @@ def __init__( grid: CoordinateGrid = auto, unit: str | Unit = auto, is_coord: bool = auto, - is_const: bool = False + is_const: bool = False, ) -> None: """ Parameters @@ -210,13 +210,17 @@ def __setattr__(self, name: str, value: Any) -> None: super().__setattr__(name, value) def _make_buffer(self) -> np.ndarray: - shape = self.shape if self.is_const else (self.proc_chain._block_width,) + self.shape + shape = ( + self.shape + if self.is_const + else (self.proc_chain._block_width,) + self.shape + ) len = np.product(shape) # Flattened array, with padding to allow memory alignment - buf = np.zeros(len + 64//self.dtype.itemsize, dtype=self.dtype) + buf = np.zeros(len + 64 // self.dtype.itemsize, dtype=self.dtype) # offset to ensure memory alignment - offset = (64 - buf.ctypes.data)%64//self.dtype.itemsize - return buf[offset:offset+len].reshape(shape) + offset = (64 - buf.ctypes.data) % 64 // self.dtype.itemsize + return buf[offset : offset + len].reshape(shape) def get_buffer(self, unit: str | Unit = None) -> np.ndarray: # If buffer needs to be created, do so now @@ -437,7 +441,7 @@ def set_constant( unit: str | Unit | Quantity = None, ) -> ProcChainVar: """Make a variable act as a constant and set it to val. - + Parameters ---------- varname @@ -452,7 +456,7 @@ def set_constant( """ param = get_variable(varname) - assert(param.is_constant or param._buffer is None) + assert param.is_constant or param._buffer is None param.is_constant = True if isinstance(val, Quantity): @@ -462,11 +466,11 @@ def set_constant( val = np.array(val, dtype=dtype) param.update_auto( - shape=val.shape, - dtype=val.dtype, - unit=unit, + shape=val.shape, + dtype=val.dtype, + unit=unit, ) - np.copyto(var.get_buffer(), val, casting='unsafe') + np.copyto(var.get_buffer(), val, casting="unsafe") log.debug(f"set constant: {self.description()} = {val}") return param @@ -751,7 +755,9 @@ def _parse_expr( if not (isinstance(lhs, ProcChainVar) or isinstance(rhs, ProcChainVar)): ret = op(lhs, rhs) - if isinstance(ret, Quantity) and ureg.is_compatible_with(ret.u, ureg.dimensionless): + if isinstance(ret, Quantity) and ureg.is_compatible_with( + ret.u, ureg.dimensionless + ): ret = ret.to(ureg.dimensionless).magnitude return ret @@ -2112,14 +2118,22 @@ def resolve_dependencies( if isinstance(arg, ProcChainVar) and arg.name in new_vars: out_is_arg = True arg.is_const = True - #arg = arg.get_buffer() + # arg = arg.get_buffer() if out_is_arg: - proc_man = ProcessorManager(proc_chain, func, arg_params, kwarg_params, kwargs.get('signature', None), kwargs.get('types', None)) + proc_man = ProcessorManager( + proc_chain, + func, + arg_params, + kwarg_params, + kwargs.get("signature", None), + kwargs.get("types", None), + ) proc_man.execute() else: const_val = func(*arg_params, **kwarg_params) - if len(new_vars)==1: const_val = (const_val) + if len(new_vars) == 1: + const_val = const_val for var, val in zip(new_vars, const_val): proc_chain.set_constant(var, val) else: From 688ee0b6794cc6362534d3fd9ecceaf53d05cc5a Mon Sep 17 00:00:00 2001 From: iguinn Date: Wed, 6 Dec 2023 07:03:53 -0800 Subject: [PATCH 18/42] Fixed things for precommit bot --- src/dspeed/processing_chain.py | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/src/dspeed/processing_chain.py b/src/dspeed/processing_chain.py index 1b7d57b..c556069 100644 --- a/src/dspeed/processing_chain.py +++ b/src/dspeed/processing_chain.py @@ -437,7 +437,7 @@ def set_constant( self, varname: str, val: np.ndarray | int | float | Quantity, - dtype: DType = None, + dtype: str | np.dtype = None, unit: str | Unit | Quantity = None, ) -> ProcChainVar: """Make a variable act as a constant and set it to val. @@ -455,7 +455,7 @@ def set_constant( unit of constant """ - param = get_variable(varname) + param = self.get_variable(varname) assert param.is_constant or param._buffer is None param.is_constant = True @@ -470,7 +470,7 @@ def set_constant( dtype=val.dtype, unit=unit, ) - np.copyto(var.get_buffer(), val, casting="unsafe") + np.copyto(param.get_buffer(), val, casting="unsafe") log.debug(f"set constant: {self.description()} = {val}") return param From 6fba8998ef11df8eaa7b60d4c3f99a38e5740bc8 Mon Sep 17 00:00:00 2001 From: iguinn Date: Wed, 6 Dec 2023 18:08:38 -0800 Subject: [PATCH 19/42] Improved logging in processing chain --- src/dspeed/processing_chain.py | 43 ++++++++++++++++++---------------- 1 file changed, 23 insertions(+), 20 deletions(-) diff --git a/src/dspeed/processing_chain.py b/src/dspeed/processing_chain.py index c556069..8e39176 100644 --- a/src/dspeed/processing_chain.py +++ b/src/dspeed/processing_chain.py @@ -264,6 +264,7 @@ def get_buffer(self, unit: str | Unit = None) -> np.ndarray: conversion_manager = UnitConversionManager(self, unit) self._buffer.append((conversion_manager.out_buffer, unit)) self.proc_chain._proc_managers.append(conversion_manager) + log.debug(f"added conversion: {conversion_manager}") return conversion_manager.out_buffer @property @@ -471,7 +472,7 @@ def set_constant( unit=unit, ) np.copyto(param.get_buffer(), val, casting="unsafe") - log.debug(f"set constant: {self.description()} = {val}") + log.debug(f"set constant: {param.description()} = {val}") return param def link_input_buffer( @@ -632,6 +633,7 @@ def add_processor( proc_man = ProcessorManager(self, func, params, kw_params, signature, types) self._proc_managers.append(proc_man) + log.debug(f"added processor: {proc_man}") def execute(self, start: int = 0, stop: int = None) -> None: """Execute the dsp chain on the entire input/output buffers.""" @@ -786,7 +788,9 @@ def _parse_expr( is_coord=rhs.is_coord, ) - self._proc_managers.append(ProcessorManager(self, op, [lhs, rhs, out])) + proc_man = ProcessorManager(self, op, [lhs, rhs, out]) + self._proc_managers.append(proc_man) + log.debug(f"added processor: {proc_man}") return out # define unary operators (-) @@ -807,7 +811,9 @@ def _parse_expr( operand.unit, operand.is_coord, ) - self._proc_managers.append(ProcessorManager(self, op, [operand, out])) + proc_man = ProcessorManager(self, op, [operand, out]) + self._proc_managers.append(proc_man) + log.debug(f"added processor: {proc_man}") else: out = op(operand) @@ -869,9 +875,9 @@ def get_index(slice_value): new_off = ProcChainVar( self, name=f"({str(off)}+{str(start)})", is_coord=True ) - self._proc_managers.append( - ProcessorManager(self, np.add, [off, start, new_off]) - ) + proc_man = ProcessorManager(self, np.add, [off, start, new_off]) + self._proc_managers.append(proc_man) + log.debug(f"added processor: {proc_man}") off = new_off else: off += start @@ -1047,6 +1053,7 @@ def _round( conversion_manager = UnitConversionManager(var, grid, round=True) out._buffer = conversion_manager.out_buffer var.proc_chain._proc_managers.append(conversion_manager) + log.debug(f"added conversion: {conversion_manager}") else: out = ProcChainVar( var.proc_chain, @@ -1080,16 +1087,16 @@ def _astype(var: ProcChainVar, dtype: str) -> ProcChainVar: # noqa: N805 var.unit, var.is_coord, ) - var.proc_chain._proc_managers.append( - ProcessorManager( - var.proc_chain, - np.copyto, - [out, var], - kw_params={"casting": "'unsafe'"}, - signature="(),(),()", - types=f"{dtype.char}{var.dtype.char}", - ) + proc_man = ProcessorManager( + var.proc_chain, + np.copyto, + [out, var], + kw_params={"casting": "'unsafe'"}, + signature="(),(),()", + types=f"{dtype.char}{var.dtype.char}", ) + var.proc_chain._proc_managers.append(proc_man) + log.debug(f"added processor: {proc_man}") return out def _loadlh5(path_to_file, path_in_file: str) -> np.array: # noqa: N805 @@ -1403,8 +1410,6 @@ def __init__( else: self.kwargs[arg_name] = param - log.debug(f"added processor: {self}") - def execute(self) -> None: self.processor(*self.args, **self.kwargs) @@ -1511,8 +1516,6 @@ def __init__( ] self.kwargs = {} - log.debug(f"added conversion: {self}") - class IOManager(metaclass=ABCMeta): r"""Base class. @@ -2152,7 +2155,7 @@ def resolve_dependencies( buf_in = lh5_in.get(copy_par) if buf_in is None: log.warning( - f"I don't know what to do with {copy_par}. Building output without it!" + f"Did not find {copy_par} in either input file or parameter list. Building output without it!" ) else: lh5_out.add_field(copy_par, buf_in) From c7b5b1e368b71f23adf62329e989965ebbe8fcb5 Mon Sep 17 00:00:00 2001 From: iguinn Date: Wed, 6 Dec 2023 18:09:31 -0800 Subject: [PATCH 20/42] Changed how build_processing_chain figures out what is a constant --- src/dspeed/processing_chain.py | 67 +++++++++++++++++++--------------- 1 file changed, 38 insertions(+), 29 deletions(-) diff --git a/src/dspeed/processing_chain.py b/src/dspeed/processing_chain.py index 8e39176..821b4f2 100644 --- a/src/dspeed/processing_chain.py +++ b/src/dspeed/processing_chain.py @@ -914,10 +914,7 @@ def get_index(slice_value): return attr # Otherwise this is probably a ProcChainVar - # Note that we are excluding this variable from the vars list - # because it does not strictly need to be computed before as a - # prerequisite before accessing its attributes - val = self._parse_expr(node.value, expr, dry_run, []) + val = self._parse_expr(node.value, expr, dry_run, var_name_list) if val is None: return None return getattr(val, node.attr) @@ -1904,8 +1901,8 @@ def build_processing_chain( for db_var in db_parser.findall(arg): try: db_node = db_dict - for key in db_var[3:].split("."): - db_node = db_node[key] + for db_key in db_var[3:].split("."): + db_node = db_node[db_key] log.debug(f"database lookup: found {db_node} for {db_var}") except (KeyError, TypeError): try: @@ -2106,41 +2103,53 @@ def resolve_dependencies( pass # Check if new variables should be treated as constants - if not recipe["prereqs"]: - arg_params = [] - kwarg_params = {} - out_is_arg = False - for arg in args: - if isinstance(arg, str): - arg = proc_chain.get_variable(arg) - if isinstance(arg, dict): - kwarg_params.update(arg) - arg = list(arg.values())[0] - else: - arg_params.append(arg) - if isinstance(arg, ProcChainVar) and arg.name in new_vars: - out_is_arg = True - arg.is_const = True - # arg = arg.get_buffer() - - if out_is_arg: + params = [] + kw_params = {} + out_params = [] + is_const=True + for param in args: + if isinstance(param, str): + param = proc_chain.get_variable(param) + if isinstance(param, dict): + kw_params.update(param) + param = list(param.values())[0] + elif isinstance(param, str): + params.append(f"'{param}'") + else: + params.append(param) + + if isinstance(param, ProcChainVar): + if param.name in new_vars: + out_params.append(param) + elif not param.is_const: + is_const = False + + if is_const: + if out_params: + for param in out_params: + param.is_const = True proc_man = ProcessorManager( proc_chain, func, - arg_params, - kwarg_params, + params, + kw_params, kwargs.get("signature", None), kwargs.get("types", None), ) proc_man.execute() + for param in out_params: + log.debug(f"set constant: {param.description()} = {param.get_buffer()}") + else: - const_val = func(*arg_params, **kwarg_params) + const_val = func(*params, **kw_params) if len(new_vars) == 1: - const_val = const_val + const_val = [const_val] for var, val in zip(new_vars, const_val): proc_chain.set_constant(var, val) + else: - proc_chain.add_processor(func, *args, **kwargs) + proc_chain.add_processor(func, *params, kw_params, **kwargs) + except Exception as e: raise ProcessingChainError( "Exception raised while attempting to add processor:\n" From a25b7284cd7e3a5d53047a1695d0b9113f09fb40 Mon Sep 17 00:00:00 2001 From: iguinn Date: Wed, 6 Dec 2023 18:10:04 -0800 Subject: [PATCH 21/42] Added errors for wrong sized arrays in convolutions --- src/dspeed/processors/convolutions.py | 6 ++++++ 1 file changed, 6 insertions(+) diff --git a/src/dspeed/processors/convolutions.py b/src/dspeed/processors/convolutions.py index 89de408..0eac256 100644 --- a/src/dspeed/processors/convolutions.py +++ b/src/dspeed/processors/convolutions.py @@ -48,10 +48,16 @@ def convolve_wf( if chr(mode_in) == "f": mode = "full" + if len(w_out) != len(w_in) + len(kernel) - 1: + raise DSPFatal(f"Output waveform has length {len(w_out)}; expect {len(w_in) + len(kernel) - 1}") elif chr(mode_in) == "v": mode = "valid" + if len(w_out) != abs(len(w_in) - len(kernel)) + 1: + raise DSPFatal(f"Output waveform has length {len(w_out)}; expect {abs(len(w_in) - len(kernel)) + 1}") elif chr(mode_in) == "s": mode = "same" + if len(w_out) != max(len(w_in), len(kernel)): + raise DSPFatal("Output waveform has length {len(w_out)}; expect {max(len(w_in), len(kernel))}") else: raise DSPFatal("Invalid mode") From ee9190fcd64b6712059445315e5c7ff76ff8378c Mon Sep 17 00:00:00 2001 From: iguinn Date: Wed, 6 Dec 2023 18:16:17 -0800 Subject: [PATCH 22/42] Cache convolution functions --- src/dspeed/processors/convolutions.py | 2 -- src/dspeed/processors/energy_kernels.py | 3 --- src/dspeed/processors/gaussian_filter1d.py | 1 - src/dspeed/processors/kernels.py | 3 --- 4 files changed, 9 deletions(-) diff --git a/src/dspeed/processors/convolutions.py b/src/dspeed/processors/convolutions.py index 0eac256..fa28c26 100644 --- a/src/dspeed/processors/convolutions.py +++ b/src/dspeed/processors/convolutions.py @@ -15,7 +15,6 @@ ], "(n),(m),(),(p)", **nb_kwargs( - cache=False, forceobj=True, ), ) @@ -71,7 +70,6 @@ def convolve_wf( ], "(n),(m),(),(p)", **nb_kwargs( - cache=False, forceobj=True, ), ) diff --git a/src/dspeed/processors/energy_kernels.py b/src/dspeed/processors/energy_kernels.py index 495e31d..7b0b65a 100644 --- a/src/dspeed/processors/energy_kernels.py +++ b/src/dspeed/processors/energy_kernels.py @@ -14,7 +14,6 @@ ], "(),(),(),(n)", **nb_kwargs( - cache=False, forceobj=True, ), ) @@ -77,7 +76,6 @@ def cusp_filter(sigma: float, flat: int, decay: int, kernel: np.array) -> None: ], "(),(),(),(n)", **nb_kwargs( - cache=False, forceobj=True, ), ) @@ -160,7 +158,6 @@ def zac_filter(sigma: float, flat: int, decay: int, kernel: np.array) -> None: ], "(n,n),(m),(),(),(),(),(n)", **nb_kwargs( - cache=False, forceobj=True, ), ) diff --git a/src/dspeed/processors/gaussian_filter1d.py b/src/dspeed/processors/gaussian_filter1d.py index 944d23a..d160f8e 100644 --- a/src/dspeed/processors/gaussian_filter1d.py +++ b/src/dspeed/processors/gaussian_filter1d.py @@ -92,7 +92,6 @@ def _gaussian_kernel1d(sigma, radius): ], "(n),(m)", **nb_kwargs( - cache=False, forceobj=True, ), ) diff --git a/src/dspeed/processors/kernels.py b/src/dspeed/processors/kernels.py index 7baa552..925aeee 100644 --- a/src/dspeed/processors/kernels.py +++ b/src/dspeed/processors/kernels.py @@ -11,7 +11,6 @@ ["void(float32, float32, float32[:])", "void(float64, float64, float64[:])"], "(),(),(n)", **nb_kwargs( - cache=False, forceobj=True, ), ) @@ -58,7 +57,6 @@ def t0_filter(rise: int, fall: int, kernel: np.array) -> None: ["void(float32[:])", "void(float64[:])"], "(n)", **nb_kwargs( - cache=False, forceobj=True, ), ) @@ -98,7 +96,6 @@ def moving_slope(kernel): ["void(float32, float32[:])", "void(float64, float64[:])"], "(),(n)", **nb_kwargs( - cache=False, forceobj=True, ), ) From b84a0290d4b8bbe38168a546ae237b39c84e0cc4 Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Thu, 7 Dec 2023 02:19:35 +0000 Subject: [PATCH 23/42] style: pre-commit fixes --- src/dspeed/processing_chain.py | 10 +++++++--- src/dspeed/processors/convolutions.py | 12 +++++++++--- 2 files changed, 16 insertions(+), 6 deletions(-) diff --git a/src/dspeed/processing_chain.py b/src/dspeed/processing_chain.py index 821b4f2..23cbdf7 100644 --- a/src/dspeed/processing_chain.py +++ b/src/dspeed/processing_chain.py @@ -875,7 +875,9 @@ def get_index(slice_value): new_off = ProcChainVar( self, name=f"({str(off)}+{str(start)})", is_coord=True ) - proc_man = ProcessorManager(self, np.add, [off, start, new_off]) + proc_man = ProcessorManager( + self, np.add, [off, start, new_off] + ) self._proc_managers.append(proc_man) log.debug(f"added processor: {proc_man}") off = new_off @@ -2106,7 +2108,7 @@ def resolve_dependencies( params = [] kw_params = {} out_params = [] - is_const=True + is_const = True for param in args: if isinstance(param, str): param = proc_chain.get_variable(param) @@ -2138,7 +2140,9 @@ def resolve_dependencies( ) proc_man.execute() for param in out_params: - log.debug(f"set constant: {param.description()} = {param.get_buffer()}") + log.debug( + f"set constant: {param.description()} = {param.get_buffer()}" + ) else: const_val = func(*params, **kw_params) diff --git a/src/dspeed/processors/convolutions.py b/src/dspeed/processors/convolutions.py index fa28c26..cf1aa48 100644 --- a/src/dspeed/processors/convolutions.py +++ b/src/dspeed/processors/convolutions.py @@ -48,15 +48,21 @@ def convolve_wf( if chr(mode_in) == "f": mode = "full" if len(w_out) != len(w_in) + len(kernel) - 1: - raise DSPFatal(f"Output waveform has length {len(w_out)}; expect {len(w_in) + len(kernel) - 1}") + raise DSPFatal( + f"Output waveform has length {len(w_out)}; expect {len(w_in) + len(kernel) - 1}" + ) elif chr(mode_in) == "v": mode = "valid" if len(w_out) != abs(len(w_in) - len(kernel)) + 1: - raise DSPFatal(f"Output waveform has length {len(w_out)}; expect {abs(len(w_in) - len(kernel)) + 1}") + raise DSPFatal( + f"Output waveform has length {len(w_out)}; expect {abs(len(w_in) - len(kernel)) + 1}" + ) elif chr(mode_in) == "s": mode = "same" if len(w_out) != max(len(w_in), len(kernel)): - raise DSPFatal("Output waveform has length {len(w_out)}; expect {max(len(w_in), len(kernel))}") + raise DSPFatal( + "Output waveform has length {len(w_out)}; expect {max(len(w_in), len(kernel))}" + ) else: raise DSPFatal("Invalid mode") From 07362bcd10897c56dfea44001395262644816cd1 Mon Sep 17 00:00:00 2001 From: iguinn Date: Wed, 6 Dec 2023 18:31:50 -0800 Subject: [PATCH 24/42] Replaced np.product with np.prod --- src/dspeed/processing_chain.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/dspeed/processing_chain.py b/src/dspeed/processing_chain.py index 23cbdf7..d35e5dc 100644 --- a/src/dspeed/processing_chain.py +++ b/src/dspeed/processing_chain.py @@ -215,7 +215,7 @@ def _make_buffer(self) -> np.ndarray: if self.is_const else (self.proc_chain._block_width,) + self.shape ) - len = np.product(shape) + len = np.prod(shape) # Flattened array, with padding to allow memory alignment buf = np.zeros(len + 64 // self.dtype.itemsize, dtype=self.dtype) # offset to ensure memory alignment From a9c8bdb59edd8f5b8111097870c6812cb4bb3f9e Mon Sep 17 00:00:00 2001 From: Luigi Pertoldi Date: Fri, 8 Dec 2023 11:06:18 +0100 Subject: [PATCH 25/42] [ci] remove ignore github actions in dependabot config --- .github/dependabot.yml | 2 -- 1 file changed, 2 deletions(-) diff --git a/.github/dependabot.yml b/.github/dependabot.yml index 368d295..f9ecf57 100644 --- a/.github/dependabot.yml +++ b/.github/dependabot.yml @@ -5,5 +5,3 @@ updates: directory: "/" schedule: interval: "monthly" - ignore: - - dependency-name: "actions/*" From f20b0c41ae34be6357e8d5d185a1095960bc90a6 Mon Sep 17 00:00:00 2001 From: "dependabot[bot]" <49699333+dependabot[bot]@users.noreply.github.com> Date: Fri, 8 Dec 2023 10:06:39 +0000 Subject: [PATCH 26/42] Bump actions/setup-python from 2 to 5 Bumps [actions/setup-python](https://github.com/actions/setup-python) from 2 to 5. - [Release notes](https://github.com/actions/setup-python/releases) - [Commits](https://github.com/actions/setup-python/compare/v2...v5) --- updated-dependencies: - dependency-name: actions/setup-python dependency-type: direct:production update-type: version-update:semver-major ... Signed-off-by: dependabot[bot] --- .github/workflows/main.yml | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/.github/workflows/main.yml b/.github/workflows/main.yml index d2e5439..25b4e17 100644 --- a/.github/workflows/main.yml +++ b/.github/workflows/main.yml @@ -26,7 +26,7 @@ jobs: steps: - uses: actions/checkout@v2 - name: Set up Python ${{ matrix.python-version }} - uses: actions/setup-python@v2 + uses: actions/setup-python@v5 with: python-version: ${{ matrix.python-version }} - name: Get dependencies and install dspeed @@ -44,7 +44,7 @@ jobs: - uses: actions/checkout@v2 with: fetch-depth: 2 - - uses: actions/setup-python@v2 + - uses: actions/setup-python@v5 with: python-version: '3.10' @@ -63,7 +63,7 @@ jobs: - uses: actions/checkout@v2 with: fetch-depth: 0 - - uses: actions/setup-python@v2 + - uses: actions/setup-python@v5 with: python-version: '3.10' - name: Setup build environment From c003fe3865072f031ad780237e17b00967512029 Mon Sep 17 00:00:00 2001 From: "dependabot[bot]" <49699333+dependabot[bot]@users.noreply.github.com> Date: Fri, 8 Dec 2023 10:06:42 +0000 Subject: [PATCH 27/42] Bump actions/checkout from 2 to 4 Bumps [actions/checkout](https://github.com/actions/checkout) from 2 to 4. - [Release notes](https://github.com/actions/checkout/releases) - [Changelog](https://github.com/actions/checkout/blob/main/CHANGELOG.md) - [Commits](https://github.com/actions/checkout/compare/v2...v4) --- updated-dependencies: - dependency-name: actions/checkout dependency-type: direct:production update-type: version-update:semver-major ... Signed-off-by: dependabot[bot] --- .github/workflows/distribution.yml | 2 +- .github/workflows/main.yml | 6 +++--- 2 files changed, 4 insertions(+), 4 deletions(-) diff --git a/.github/workflows/distribution.yml b/.github/workflows/distribution.yml index add54cb..8a702f3 100644 --- a/.github/workflows/distribution.yml +++ b/.github/workflows/distribution.yml @@ -13,7 +13,7 @@ jobs: dist: runs-on: ubuntu-latest steps: - - uses: actions/checkout@v3 + - uses: actions/checkout@v4 with: fetch-depth: 0 diff --git a/.github/workflows/main.yml b/.github/workflows/main.yml index d2e5439..3e07337 100644 --- a/.github/workflows/main.yml +++ b/.github/workflows/main.yml @@ -24,7 +24,7 @@ jobs: os: [ubuntu-latest, macOS-latest] steps: - - uses: actions/checkout@v2 + - uses: actions/checkout@v4 - name: Set up Python ${{ matrix.python-version }} uses: actions/setup-python@v2 with: @@ -41,7 +41,7 @@ jobs: name: Calculate and upload test coverage runs-on: ubuntu-latest steps: - - uses: actions/checkout@v2 + - uses: actions/checkout@v4 with: fetch-depth: 2 - uses: actions/setup-python@v2 @@ -60,7 +60,7 @@ jobs: name: Build documentation runs-on: ubuntu-latest steps: - - uses: actions/checkout@v2 + - uses: actions/checkout@v4 with: fetch-depth: 0 - uses: actions/setup-python@v2 From 34f190f8a856e8c363a089677e241597e30efbd1 Mon Sep 17 00:00:00 2001 From: Luigi Pertoldi Date: Fri, 8 Dec 2023 16:24:03 +0100 Subject: [PATCH 28/42] [docs] import dev guide from pygama --- docs/source/developer.rst | 323 ++++++++++++++++++++++++++++++++++++++ docs/source/index.rst | 17 ++ 2 files changed, 340 insertions(+) create mode 100644 docs/source/developer.rst diff --git a/docs/source/developer.rst b/docs/source/developer.rst new file mode 100644 index 0000000..ac50543 --- /dev/null +++ b/docs/source/developer.rst @@ -0,0 +1,323 @@ +Developer's guide +================= + +.. note:: + + The https://learn.scientific-python.org webpages are an extremely valuable + learning resource for Python software developer. The reader is referred to + that for any detail not covered in the following guide. + +The following rules and conventions have been established for the package +development and are enforced throughout the entire code base. Merge requests +that do not comply to the following directives will be rejected. + +To start developing :mod:`dspeed`, fork the remote repository to your personal +GitHub account (see `About Forks +`_). +If you have not set up your ssh keys on the computer you will be working on, +please follow `GitHub's instructions +`_. +Once you have your own fork, you can clone it via +(replace "yourusername" with your GitHub username): + +.. code-block:: console + + $ git clone git@github.com:yourusername/dspeed.git + +All extra tools needed to develop :mod:`dspeed` are listed as optional +dependencies and can be installed via pip by running: + +.. code-block:: console + + $ cd dspeed + $ pip install -e '.[all]' # single quotes are not needed on bash + +.. important:: + + Pip's ``--editable | -e`` flag let's you install the package in "developer + mode", meaning that any change to the source code will be directly + propagated to the installed package and importable in scripts. + +.. tip:: + + It is strongly recommended to work inside a virtual environment, which + guarantees reproductibility and isolation. For more details, see + `learn.scientific-python.org + `_. + +Code style +---------- + +* All functions and methods (arguments and return types) must be + `type-annotated `_. Type + annotations for variables like class attributes are also highly appreciated. +* Messaging to the user is managed through the :mod:`logging` module. Do not + add :func:`print` statements. To make a logging object available in a module, + add this: + + .. code-block:: python + + import logging + log = logging.getLogger(__name__) + + at the top. In general, try to keep the number of :func:`logging.debug` calls + low and use informative messages. :func:`logging.info` calls should be + reserved for messages from high-level routines (like + :func:`dspeed.dsp.build_dsp`) and very sporadic. Good code is never too + verbose. +* If an error condition leading to undefined behavior occurs, raise an + exception. try to find the most suitable between the `built-in exceptions + `_, otherwise ``raise + RuntimeError("message")``. Do not raise ``Warning``\ s, use + :func:`logging.warning` for that and don't abort the execution. +* Warning messages (emitted when a problem is encountered that does not lead to + undefined behavior) must be emitted through :func:`logging.warning` calls. + +A set of `pre-commit `_ hooks is configured to make +sure that :mod:`dspeed` coherently follows standard coding style conventions. +The pre-commit tool is able to identify common style problems and automatically +fix them, wherever possible. Configured hooks are listed in the +``.pre-commit-config.yaml`` file at the project root folder. They are run +remotely on the GitHub repository through the `pre-commit bot +`_, but should also be run locally before submitting a +pull request: + +.. code-block:: console + + $ cd dspeed + $ pip install '.[test]' + $ pre-commit run --all-files # analyse the source code and fix it wherever possible + $ pre-commit install # install a Git pre-commit hook (strongly recommended) + +For a more comprehensive guide, check out the `learn.scientific-python.org +documentation about code style +`_. + +Testing +------- + +* The :mod:`dspeed` test suite is available below ``tests/``. We use `pytest + `_ to run tests and analyze their output. As + a starting point to learn how to write good tests, reading of `the + relevant learn.scientific-python.org webpage + `_ is + recommended. + +* Unit tests are automatically run for every push event and pull request to the + remote Git repository on a remote server (currently handled by GitHub + actions). Every pull request must pass all tests before being approved for + merging. Running the test suite is simple: + + .. code-block:: console + + $ cd dspeed + $ pip install '.[test]' + $ pytest + +* Additionally, pull request authors are required to provide tests with + sufficient code coverage for every proposed change or addition. If necessary, + high-level functional tests should be updated. We currently rely on + `codecov.io `_ to keep track of + test coverage. A local report, which must be inspected before submitting pull + requests, can be generated by running: + + .. code-block:: console + + $ pytest --cov=dspeed + +Testing Numba-Wrapped Functions +^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ + +When using Numba to vectorize Python functions, the Python version of the function +does not, by default, get directly tested, but the Numba version instead. In +this case, we need to unwrap the Numba function and test the pure Python version. +With various processors in :mod:`dspeed.processors`, this means that testing +and triggering the code coverage requires this unwrapping. + +Within the testing suite, we use the :func:`@pytest.fixture()` +decorator to include a helper function called ``compare_numba_vs_python`` that +can be used in any test. This function runs both the Numba and pure Python versions +of a function, asserts that they are equal up to floating precision, and returns the +output value. + +As an example, we show a snippet from the test for +:func:`dspeed.processors.fixed_time_pickoff`, a processor which uses the +:func:`@numba.guvectorize()` decorator. + +.. code-block:: python + + def test_fixed_time_pickoff(compare_numba_vs_python): + """Testing function for the fixed_time_pickoff processor.""" + + len_wf = 20 + + # test for nan if w_in has a nan + w_in = np.ones(len_wf) + w_in[4] = np.nan + assert np.isnan(compare_numba_vs_python(fixed_time_pickoff, w_in, 1, ord("i"))) + +In the assertion that the output is what we expect, we use +``compare_numba_vs_python(fixed_time_pickoff, w_in, 1, ord("i"))`` in place of +``fixed_time_pickoff(w_in, 1, ord("i"))``. In general, the replacement to make is +``func(*inputs)`` becomes ``compare_numba_vs_python(func, *inputs)``. + +Note, that in cases of testing for the raising of errors, it is recommended +to instead run the function twice: once with the Numba version, and once using the +:func:`inspect.unwrap` function. We again show a snippet from the test for +:func:`dspeed.processors.fixed_time_pickoff` below. We include the various +required imports in the snippet for verbosity. + +.. code-block:: python + + import inspect + + import numpy as np + import pytest + + from dspeed.dsp.errors import DSPFatal + from dspeed.dsp.processors import fixed_time_pickoff + + def test_fixed_time_pickoff(compare_numba_vs_python): + "skipping parts of function..." + # test for DSPFatal errors being raised + # noninteger t_in with integer interpolation + with pytest.raises(DSPFatal): + w_in = np.ones(len_wf) + fixed_time_pickoff(w_in, 1.5, ord("i")) + + with pytest.raises(DSPFatal): + a_out = np.empty(len_wf) + inspect.unwrap(fixed_time_pickoff)(w_in, 1.5, ord("i"), a_out) + +In this case, the general idea is to use :func:`pytest.raises` twice, once with +``func(*inputs)``, and again with ``inspect.unwrap(func)(*inputs)``. + +Testing Factory Functions that Return Numba-Wrapped Functions +^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ + +As in the previous section, we also have processors that are first initialized +with a factory function, which then returns a callable Numba-wrapped function. +In this case, there is a slightly different way of testing the function to ensure +full code coverage when using ``compare_numba_vs_python``, as the function +signature is generally different. + +As an example, we show a snippet from the test for +:func:`dspeed.processors.dwt.discrete_wavelet_transform`, a processor which uses +a factory function to return a function wrapped by the +:func:`@numba.guvectorize()` decorator. + +.. code-block:: python + + import numpy as np + import pytest + + from dspeed.dsp.errors import DSPFatal + from dspeed.dsp.processors import discrete_wavelet_transform + + def test_discrete_wavelet_transform(compare_numba_vs_python): + """Testing function for the discrete_wavelet_transform processor.""" + + # set up values to use for each test case + len_wf_in = 16 + wave_type = 'haar' + level = 2 + len_wf_out = 4 + + # ensure the DSPFatal is raised for a negative level + with pytest.raises(DSPFatal): + discrete_wavelet_transform(wave_type, -1) + + # ensure that a valid input gives the expected output + w_in = np.ones(len_wf_in) + w_out = np.empty(len_wf_out) + w_out_expected = np.ones(len_wf_out) * 2**(level / 2) + + dwt_func = discrete_wavelet_transform(wave_type, level) + assert np.allclose( + compare_numba_vs_python(dwt_func, w_in, w_out), + w_out_expected, + ) + ## rest of test function is truncated in this example + +In this case, the error is raised outside of the Numba-wrapped function, and +we only need to test for the error once. For the comparison of the calculated +values to expectation, we must initialize the output array and pass it to the +list of inputs that should be used in the comparison. This is different than +the previous section, where we are instead now updating the outputted values +in place. + +Documentation +------------- + +We adopt best practices in writing and maintaining :mod:`dspeed`'s +documentation. When contributing to the project, make sure to implement the +following: + +* Documentation should be exclusively available on the Project website + `dspeed.readthedocs.io `_. No READMEs, + GitHub/LEGEND wiki pages should be written. +* Pull request authors are required to provide sufficient documentation for + every proposed change or addition. +* Documentation for functions, classes, modules and packages should be provided + as `Docstrings `_ along with the respective + source code. Docstrings are automatically converted to HTML as part of the + :mod:`dspeed` package API documentation. +* General guides, comprehensive tutorials or other high-level documentation + (e.g. referring to how separate parts of the code interact between each + other) must be provided as separate pages in ``docs/source/`` and linked in + the table of contents. +* Jupyter notebooks should be added to the main Git repository below + ``docs/source/notebooks``. +* Before submitting a pull request, contributors are required to build the + documentation locally and resolve and warnings or errors. + +Writing documentation +^^^^^^^^^^^^^^^^^^^^^ + +We adopt the following guidelines for writing documentation: + +* Documentation source files must formatted in reStructuredText (reST). A + reference format specification is available on the `Sphinx reST usage guide + `_. + Usage of `Cross-referencing syntax + `_ + in general and `for Python objects + `_ + in particular is recommended. We also support cross-referencing external + documentation via `sphinx.ext.intersphinx + `_, + when referring for example to :class:`pandas.DataFrame`. +* To document Python objects, we also adopt the `NumPy Docstring style + `_. Examples are + available `here + `_. +* We support also the Markdown format through the `MyST-Parser + `_. +* Jupyter notebooks placed below ``docs/source/notebooks`` are automatically + rendered to HTML pages by the `nbsphinx `_ + extension. + +Building documentation +^^^^^^^^^^^^^^^^^^^^^^ + +Scripts and tools to build documentation are located below ``docs/``. To build +documentation, ``sphinx`` and a couple of additional Python packages are +required. You can get all the needed dependencies by running: + +.. code-block:: console + + $ cd dspeed + $ pip install '.[docs]' + +`Pandoc `_ is also required to render +Jupyter notebooks. To build documentation, run the following commands: + +.. code-block:: console + + $ cd docs + $ make clean + $ make + +Documentation can be then displayed by opening ``build/html/index.html`` with a +web browser. Documentation for the :mod:`dspeed` website is built and deployed by +`Read the Docs `_. diff --git a/docs/source/index.rst b/docs/source/index.rst index a4ff76b..50d0aaa 100644 --- a/docs/source/index.rst +++ b/docs/source/index.rst @@ -51,6 +51,23 @@ Table of Contents tutorials Package API reference +.. toctree:: + :maxdepth: 1 + :caption: Related projects + + LEGEND Data Objects + Decoding Digitizer Data + Data processing and analysis + +.. toctree:: + :maxdepth: 1 + :caption: Development + + Source Code + License + Changelog + developer + .. _2: https://pypi.org/project/dspeed .. _3: https://pip.pypa.io/en/stable/getting-started .. |dspeed| replace:: *dspeed* From b6f1a3d6beee20cfa708370893ecdca811a20470 Mon Sep 17 00:00:00 2001 From: Luigi Pertoldi Date: Fri, 8 Dec 2023 16:49:13 +0100 Subject: [PATCH 29/42] [docs] add versioning section to dev's guide --- docs/source/developer.rst | 40 +++++++++++++++++++++++++++++++++++++++ 1 file changed, 40 insertions(+) diff --git a/docs/source/developer.rst b/docs/source/developer.rst index ac50543..d9a2d3b 100644 --- a/docs/source/developer.rst +++ b/docs/source/developer.rst @@ -321,3 +321,43 @@ Jupyter notebooks. To build documentation, run the following commands: Documentation can be then displayed by opening ``build/html/index.html`` with a web browser. Documentation for the :mod:`dspeed` website is built and deployed by `Read the Docs `_. + +Versioning +---------- + +Collaborators with push access to the GitHub repository that wish to release a +new project version must implement the following procedures: + +* `Semantic versioning `_ is adopted. The version string + uses the ``MAJOR.MINOR.PATCH`` format. +* To release a new **minor** or **major version**, the following procedure + should be followed: + + 1. A new branch with name ``releases/vMAJOR.MINOR`` (note the ``v``) containing + the code at the intended stage is created + 2. The commit is tagged with a descriptive message: ``git tag vMAJOR.MINOR.0 + -m 'short descriptive message here'`` (note the ``v``) + 3. Changes are pushed to the remote: + + .. code-block:: console + + $ git push origin releases/vMAJOR.MINOR + $ git push origin refs/tags/vMAJOR.MINOR.0 + +* To release a new **patch version**, the following procedure should be followed: + + 1. A commit with the patch is created on the relevant release branch + ``releases/vMAJOR.MINOR`` + 2. The commit is tagged: ``git tag vMAJOR.MINOR.PATCH`` (note the ``v``) + 3. Changes are pushed to the remote: + + .. code-block:: console + + $ git push origin releases/vMAJOR.MINOR + $ git push origin refs/tags/vMAJOR.MINOR.PATCH + +* To upload the release to the `Python Package Index + `_, a new release must be created through + `the GitHub interface `_, + associated to the just created tag. Usage of the "Generate release notes" + option is recommended. From f27e95d02b74f943420a1231e33b3e3da0b1b6a0 Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Wed, 13 Dec 2023 16:46:42 +0000 Subject: [PATCH 30/42] style: pre-commit fixes --- src/dspeed/processors/min_max.py | 1 + 1 file changed, 1 insertion(+) diff --git a/src/dspeed/processors/min_max.py b/src/dspeed/processors/min_max.py index cad858a..92cdd2f 100644 --- a/src/dspeed/processors/min_max.py +++ b/src/dspeed/processors/min_max.py @@ -71,6 +71,7 @@ def min_max( t_min[0] = float(min_index) t_max[0] = float(max_index) + @guvectorize( [ "void(float32[:], float32[:], float32[:], float32[:])", From d84a33d45e800c880386d3c95122e39ae111cd3e Mon Sep 17 00:00:00 2001 From: Esteban Leon Date: Thu, 14 Dec 2023 16:30:49 -0800 Subject: [PATCH 31/42] Added svm processor, added min_max_norm processor test, refactored dwt processor without use of factory function --- src/dspeed/processors/__init__.py | 2 + src/dspeed/processors/dwt.py | 77 ++++++++++++++------------- src/dspeed/processors/svm.py | 74 +++++++++++++++++++++++++ tests/processors/test_dwt.py | 2 +- tests/processors/test_min_max_norm.py | 45 ++++++++++++++++ 5 files changed, 161 insertions(+), 39 deletions(-) create mode 100644 src/dspeed/processors/svm.py create mode 100644 tests/processors/test_min_max_norm.py diff --git a/src/dspeed/processors/__init__.py b/src/dspeed/processors/__init__.py index d545d50..12bc5d9 100644 --- a/src/dspeed/processors/__init__.py +++ b/src/dspeed/processors/__init__.py @@ -90,6 +90,7 @@ from .round_to_nearest import round_to_nearest from .saturation import saturation from .soft_pileup_corr import soft_pileup_corr, soft_pileup_corr_bl +from .svm import svm_predict from .time_over_threshold import time_over_threshold from .time_point_thresh import interpolated_time_point_thresh, time_point_thresh from .transfer_function_convolver import transfer_function_convolver @@ -139,6 +140,7 @@ "peak_snr_threshold", "soft_pileup_corr", "soft_pileup_corr_bl", + "svm_predict", "time_point_thresh", "interpolated_time_point_thresh", "asym_trap_filter", diff --git a/src/dspeed/processors/dwt.py b/src/dspeed/processors/dwt.py index 137e4bc..514113c 100644 --- a/src/dspeed/processors/dwt.py +++ b/src/dspeed/processors/dwt.py @@ -1,7 +1,5 @@ from __future__ import annotations -from typing import Callable - import numpy as np from numba import guvectorize from pywt import downcoef @@ -10,67 +8,70 @@ from ..utils import numba_defaults_kwargs as nb_kwargs -def discrete_wavelet_transform(wave_type: str, level: int, coeff: str) -> Callable: +@guvectorize( + [ + "void(float32[:], int32, char, char, float32[:])", + "void(float64[:], int64, char, char, float64[:])", + ], + "(n),(),(),(),(m)", + **nb_kwargs( + forceobj=True, + ), +) + +def discrete_wavelet_transform( + w_in: np.ndarray, level: int, wave_type: np.int8, coeff: np.int8, w_out: np.ndarray +) -> None: """ Apply a discrete wavelet transform to the waveform and return only - the approximate coefficients. - - Note - ---- - This processor is composed of a factory function that is called using the - ``init_args`` argument. The input and output waveforms are passed using - ``args``. The output waveform dimension must be specified. + the detailed or approximate coefficients. Parameters ---------- - wave_type - The wavelet type for discrete convolution ``('haar', 'db1', ...)`` + + w_in + The input waveform level The level of decompositions to be performed ``(1, 2, ...)`` + wave_type + The wavelet type for discrete convolution ``('h' = 'haar', 'd' = 'db1')``. coeff The coefficients to be saved ``('a', 'd')`` + w_out + The approximate coefficients. The dimension of this array can be calculated + by ``out_dim = len(w_in)/(filter_length^level)``, where ``filter_length`` + can be obtained via ``pywt.Wavelet(wave_type).dec_len`` JSON Configuration Example -------------------------- .. code-block :: json - "dwt":{ + "dwt_haar":{ "function": "discrete_wavelet_transform", "module": "dspeed.processors", - "args": ["wf_blsub", "dwt(250)"], + "args": ["wf_blsub", 5, "'h'", "'a'", "dwt_haar(256, 'f')"], "unit": "ADC", "prereqs": ["wf_blsub"], - "init_args": ["'haar'", "3", "a"] } """ if level <= 0: raise DSPFatal("The level must be a positive integer") - - @guvectorize( - ["void(float32[:], float32[:])", "void(float64[:], float64[:])"], - "(n),(m)", - **nb_kwargs, - forceobj=True, - ) - def dwt_out(w_in: np.ndarray, w_out: np.ndarray) -> None: - """ - Parameters - ---------- - w_in - The input waveform - w_out - The approximate coefficients. The dimension of this array can be calculated - by ``out_dim = wf_length/(filter_length^level)``, where ``filter_length`` - can be obtained via ``pywt.Wavelet(wave_type).dec_len``. - """ - w_out[:] = np.nan - - if np.isnan(w_in).any(): + + if np.isnan(w_in).any(): return - w_out[:] = downcoef(coeff, w_in, wave_type, level=level) + w_out[:] = np.nan + + coeff = chr(coeff) + + if chr(wave_type) == 'h': + wave_type = 'haar' + + elif chr(wave_type) == 'd': + wave_type = 'db1' + + w_out[:] = downcoef(coeff, w_in, wave_type, level=level) - return dwt_out diff --git a/src/dspeed/processors/svm.py b/src/dspeed/processors/svm.py new file mode 100644 index 0000000..92590df --- /dev/null +++ b/src/dspeed/processors/svm.py @@ -0,0 +1,74 @@ +from __future__ import annotations + +from typing import Callable + +import pickle +import numpy as np +from numba import guvectorize + +from ..errors import DSPFatal +from ..utils import numba_defaults_kwargs as nb_kwargs + + +def svm_predict(svm_file: str) -> Callable: + """ + Apply a Support Vector Machine (SVM) to an input waveform to + predict a data cleaning label. + + Note + ---- + This processor is composed of a factory function that is called + using the ``init_args`` argument. The input waveform and output + label are passed using ``args``. + + + Parameters + ---------- + svm_file + The name of the file with the trained SVM ``svm_p*_r*_T***Z.sav`` + + + JSON Configuration Example + -------------------------- + .. code-block :: json + + "svm_label":{ + "function": "svm_predict", + "module": "dspeed.processors", + "args": ["dwt_norm", "svm_label"], + "unit": "", + "prereqs": ["dwt_norm"], + "init_args": ["'svm_p*_r*_T***Z.sav'"] + } + """ + + with open(svm_file, 'rb') as f: + svm = pickle.load(f) + + @guvectorize( + [ + "void(float32[:], float32[:])", + "void(float64[:], float64[:])", + ], + "(n),()", + **nb_kwargs( + forceobj=True, + ), + ) + def svm_out(w_in: np.ndarray, label_out: float) -> None: + """ + Parameters + ---------- + w_in + The input waveform (has to be a max_min normalized discrete wavelet transform) + label_out + The predicted label by the trained SVM for the input waveform. + """ + label_out[0] = np.nan + + if w_in.ndim == 1: + label_out[0] = svm.predict(w_in.reshape(1,-1)) + else: + label_out[0] = svm.predict(w_in) + + return svm_out diff --git a/tests/processors/test_dwt.py b/tests/processors/test_dwt.py index 54421f5..e269aca 100644 --- a/tests/processors/test_dwt.py +++ b/tests/processors/test_dwt.py @@ -10,7 +10,7 @@ def test_discrete_wavelet_transform(compare_numba_vs_python): # set up values to use for each test case len_wf_in = 16 - wave_type = "haar" + wave_type = "h" level = 2 coeff = "a" len_wf_out = 4 diff --git a/tests/processors/test_min_max_norm.py b/tests/processors/test_min_max_norm.py new file mode 100644 index 0000000..d9dede6 --- /dev/null +++ b/tests/processors/test_min_max_norm.py @@ -0,0 +1,45 @@ +import numpy as np +import pytest + +from dspeed.errors import DSPFatal +from dspeed.processors import min_max_norm + + +def test_min_max_norm(compare_numba_vs_python): + """Testing function for the max_min_norm processor.""" + + # set up values to use for each test case + len_wf = 10 + + # test for nan if w_in has a nan + w_in = np.ones(len_wf) + w_in[4] = np.nan + w_out = np.empty(len_wf) + assert np.isnan(compare_numba_vs_python(min_max_norm, w_in, -1, 1, w_out)) + + # test for division by 0 + w_in = np.ones(len_wf) + w_out = np.ones(len_wf) + assert np.all(compare_numba_vs_python(min_max_norm, w_in, 0, 0, w_out)) + + # test for abs(a_max) > abs(a_min) + w_in = np.ones(len_wf) + a_max = 2 + a_min = -1 + w_out = np.ones(len_wf) + w_out_expected = np.ones(len_wf)/abs(a_max) + assert np.allclose( + compare_numba_vs_python(min_max_norm, w_in, a_min, a_max, w_out), + w_out_expected, + ) + + # test for abs(a_max) < abs(a_min) + w_in = np.ones(len_wf) + a_max = 1 + a_min = -2 + w_out = np.ones(len_wf) + w_out_expected = np.ones(len_wf)/abs(a_min) + assert np.allclose( + compare_numba_vs_python(min_max_norm, w_in, a_min, a_max, w_out), + w_out_expected, + ) \ No newline at end of file From 5356dc4f7576f2e86fc4ef74281c4a63dde7ebf5 Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Fri, 15 Dec 2023 00:36:10 +0000 Subject: [PATCH 32/42] style: pre-commit fixes --- src/dspeed/processors/dwt.py | 22 ++++++++++------------ src/dspeed/processors/svm.py | 17 ++++++++--------- tests/processors/test_min_max_norm.py | 14 ++++++-------- 3 files changed, 24 insertions(+), 29 deletions(-) diff --git a/src/dspeed/processors/dwt.py b/src/dspeed/processors/dwt.py index 514113c..b33a34c 100644 --- a/src/dspeed/processors/dwt.py +++ b/src/dspeed/processors/dwt.py @@ -18,7 +18,6 @@ forceobj=True, ), ) - def discrete_wavelet_transform( w_in: np.ndarray, level: int, wave_type: np.int8, coeff: np.int8, w_out: np.ndarray ) -> None: @@ -59,19 +58,18 @@ def discrete_wavelet_transform( if level <= 0: raise DSPFatal("The level must be a positive integer") - + if np.isnan(w_in).any(): - return + return w_out[:] = np.nan - + coeff = chr(coeff) - - if chr(wave_type) == 'h': - wave_type = 'haar' - - elif chr(wave_type) == 'd': - wave_type = 'db1' - - w_out[:] = downcoef(coeff, w_in, wave_type, level=level) + if chr(wave_type) == "h": + wave_type = "haar" + + elif chr(wave_type) == "d": + wave_type = "db1" + + w_out[:] = downcoef(coeff, w_in, wave_type, level=level) diff --git a/src/dspeed/processors/svm.py b/src/dspeed/processors/svm.py index 92590df..46c329d 100644 --- a/src/dspeed/processors/svm.py +++ b/src/dspeed/processors/svm.py @@ -1,24 +1,23 @@ from __future__ import annotations +import pickle from typing import Callable -import pickle import numpy as np from numba import guvectorize -from ..errors import DSPFatal from ..utils import numba_defaults_kwargs as nb_kwargs def svm_predict(svm_file: str) -> Callable: """ - Apply a Support Vector Machine (SVM) to an input waveform to + Apply a Support Vector Machine (SVM) to an input waveform to predict a data cleaning label. Note ---- - This processor is composed of a factory function that is called - using the ``init_args`` argument. The input waveform and output + This processor is composed of a factory function that is called + using the ``init_args`` argument. The input waveform and output label are passed using ``args``. @@ -41,8 +40,8 @@ def svm_predict(svm_file: str) -> Callable: "init_args": ["'svm_p*_r*_T***Z.sav'"] } """ - - with open(svm_file, 'rb') as f: + + with open(svm_file, "rb") as f: svm = pickle.load(f) @guvectorize( @@ -67,8 +66,8 @@ def svm_out(w_in: np.ndarray, label_out: float) -> None: label_out[0] = np.nan if w_in.ndim == 1: - label_out[0] = svm.predict(w_in.reshape(1,-1)) + label_out[0] = svm.predict(w_in.reshape(1, -1)) else: - label_out[0] = svm.predict(w_in) + label_out[0] = svm.predict(w_in) return svm_out diff --git a/tests/processors/test_min_max_norm.py b/tests/processors/test_min_max_norm.py index d9dede6..692bd1a 100644 --- a/tests/processors/test_min_max_norm.py +++ b/tests/processors/test_min_max_norm.py @@ -1,7 +1,5 @@ import numpy as np -import pytest -from dspeed.errors import DSPFatal from dspeed.processors import min_max_norm @@ -16,18 +14,18 @@ def test_min_max_norm(compare_numba_vs_python): w_in[4] = np.nan w_out = np.empty(len_wf) assert np.isnan(compare_numba_vs_python(min_max_norm, w_in, -1, 1, w_out)) - - # test for division by 0 + + # test for division by 0 w_in = np.ones(len_wf) w_out = np.ones(len_wf) assert np.all(compare_numba_vs_python(min_max_norm, w_in, 0, 0, w_out)) - + # test for abs(a_max) > abs(a_min) w_in = np.ones(len_wf) a_max = 2 a_min = -1 w_out = np.ones(len_wf) - w_out_expected = np.ones(len_wf)/abs(a_max) + w_out_expected = np.ones(len_wf) / abs(a_max) assert np.allclose( compare_numba_vs_python(min_max_norm, w_in, a_min, a_max, w_out), w_out_expected, @@ -38,8 +36,8 @@ def test_min_max_norm(compare_numba_vs_python): a_max = 1 a_min = -2 w_out = np.ones(len_wf) - w_out_expected = np.ones(len_wf)/abs(a_min) + w_out_expected = np.ones(len_wf) / abs(a_min) assert np.allclose( compare_numba_vs_python(min_max_norm, w_in, a_min, a_max, w_out), w_out_expected, - ) \ No newline at end of file + ) From ab13113f0912c889f0a6de07ec68ea8efe3d39bc Mon Sep 17 00:00:00 2001 From: Esteban Leon Date: Mon, 18 Dec 2023 12:42:50 -0800 Subject: [PATCH 33/42] Updated dwt and min_max_norm tests --- src/dspeed/processors/dwt.py | 12 ++++---- src/dspeed/processors/min_max.py | 4 +-- tests/processors/test_dwt.py | 25 +++++++-------- tests/processors/test_min_max_norm.py | 44 +++++++++++++++------------ 4 files changed, 45 insertions(+), 40 deletions(-) diff --git a/src/dspeed/processors/dwt.py b/src/dspeed/processors/dwt.py index 514113c..8c70f89 100644 --- a/src/dspeed/processors/dwt.py +++ b/src/dspeed/processors/dwt.py @@ -13,14 +13,14 @@ "void(float32[:], int32, char, char, float32[:])", "void(float64[:], int64, char, char, float64[:])", ], - "(n),(),(),(),(m)", + "(n),(),(),()->(m)", **nb_kwargs( forceobj=True, ), ) def discrete_wavelet_transform( - w_in: np.ndarray, level: int, wave_type: np.int8, coeff: np.int8, w_out: np.ndarray + w_in: np.ndarray, level: int, wave_type: int, coeff: int, w_out: np.ndarray ) -> None: """ Apply a discrete wavelet transform to the waveform and return only @@ -56,14 +56,14 @@ def discrete_wavelet_transform( "prereqs": ["wf_blsub"], } """ - + + w_out[:] = np.nan + if level <= 0: raise DSPFatal("The level must be a positive integer") if np.isnan(w_in).any(): - return - - w_out[:] = np.nan + return coeff = chr(coeff) diff --git a/src/dspeed/processors/min_max.py b/src/dspeed/processors/min_max.py index 92cdd2f..45b5283 100644 --- a/src/dspeed/processors/min_max.py +++ b/src/dspeed/processors/min_max.py @@ -111,11 +111,11 @@ def min_max_norm( } """ + w_out[:] = np.nan + if np.isnan(w_in).any(): return - w_out[:] = np.nan - if abs(a_max[0]) == 0 or abs(a_min[0]) == 0: w_out[:] = w_in[:] diff --git a/tests/processors/test_dwt.py b/tests/processors/test_dwt.py index e269aca..bcdb546 100644 --- a/tests/processors/test_dwt.py +++ b/tests/processors/test_dwt.py @@ -10,30 +10,31 @@ def test_discrete_wavelet_transform(compare_numba_vs_python): # set up values to use for each test case len_wf_in = 16 - wave_type = "h" + wave_type = ord('h') level = 2 - coeff = "a" + coeff = ord('a') len_wf_out = 4 # ensure the DSPFatal is raised for a negative level + w_in = np.ones(len_wf_in) + w_out = np.empty(len_wf_out) with pytest.raises(DSPFatal): - discrete_wavelet_transform(wave_type, -1, coeff) + discrete_wavelet_transform(w_in, -1, wave_type, coeff, w_out) # ensure that a valid input gives the expected output - w_in = np.ones(len_wf_in) - w_out = np.empty(len_wf_out) w_out_expected = np.ones(len_wf_out) * 2 ** (level / 2) - - dwt_func = discrete_wavelet_transform(wave_type, level, coeff) assert np.allclose( - compare_numba_vs_python(dwt_func, w_in, w_out), + compare_numba_vs_python(discrete_wavelet_transform, w_in, level, wave_type, coeff, w_out), w_out_expected, ) # ensure that if there is a nan in w_in, all nans are outputted w_in = np.ones(len_wf_in) w_in[4] = np.nan - w_out = np.empty(len_wf_out) - - dwt_func = discrete_wavelet_transform(wave_type, level, coeff) - assert np.all(np.isnan(compare_numba_vs_python(dwt_func, w_in, w_out))) + assert np.all( + np.isnan( + compare_numba_vs_python( + discrete_wavelet_transform, w_in, level, wave_type, coeff + ) + ) + ) \ No newline at end of file diff --git a/tests/processors/test_min_max_norm.py b/tests/processors/test_min_max_norm.py index d9dede6..2ef9c74 100644 --- a/tests/processors/test_min_max_norm.py +++ b/tests/processors/test_min_max_norm.py @@ -11,35 +11,39 @@ def test_min_max_norm(compare_numba_vs_python): # set up values to use for each test case len_wf = 10 - # test for nan if w_in has a nan + # ensure that if there is a nan in w_in, all nans are outputted w_in = np.ones(len_wf) w_in[4] = np.nan - w_out = np.empty(len_wf) - assert np.isnan(compare_numba_vs_python(min_max_norm, w_in, -1, 1, w_out)) + a = np.array([1]) + assert np.all( + np.isnan( + compare_numba_vs_python(min_max_norm, w_in, a, a) + ) + ) - # test for division by 0 + # test for division by 0 w_in = np.ones(len_wf) - w_out = np.ones(len_wf) - assert np.all(compare_numba_vs_python(min_max_norm, w_in, 0, 0, w_out)) + a = np.array([0]) + w_out_expected = np.ones(len_wf) + assert np.allclose( + compare_numba_vs_python(min_max_norm, w_in, a, a), + w_out_expected + ) # test for abs(a_max) > abs(a_min) - w_in = np.ones(len_wf) - a_max = 2 - a_min = -1 - w_out = np.ones(len_wf) - w_out_expected = np.ones(len_wf)/abs(a_max) + a_max = np.array([2]) + a_min = np.array([-1]) + w_out_expected = np.ones(len_wf)/abs(a_max[0]) assert np.allclose( - compare_numba_vs_python(min_max_norm, w_in, a_min, a_max, w_out), - w_out_expected, + compare_numba_vs_python(min_max_norm, w_in, a_min, a_max), + w_out_expected ) # test for abs(a_max) < abs(a_min) - w_in = np.ones(len_wf) - a_max = 1 - a_min = -2 - w_out = np.ones(len_wf) - w_out_expected = np.ones(len_wf)/abs(a_min) + a_max = np.array([1]) + a_min = np.array([-2]) + w_out_expected = np.ones(len_wf)/abs(a_min[0]) assert np.allclose( - compare_numba_vs_python(min_max_norm, w_in, a_min, a_max, w_out), - w_out_expected, + compare_numba_vs_python(min_max_norm, w_in, a_min, a_max), + w_out_expected ) \ No newline at end of file From 9e65d9d0e6cb783faf02c791cbee2b0382b2e060 Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Mon, 18 Dec 2023 22:10:24 +0000 Subject: [PATCH 34/42] style: pre-commit fixes --- src/dspeed/processors/dwt.py | 22 +++++++++++----------- src/dspeed/processors/min_max.py | 2 +- tests/processors/test_dwt.py | 10 ++++++---- tests/processors/test_min_max_norm.py | 23 ++++++++--------------- 4 files changed, 26 insertions(+), 31 deletions(-) diff --git a/src/dspeed/processors/dwt.py b/src/dspeed/processors/dwt.py index 399a14a..081526d 100644 --- a/src/dspeed/processors/dwt.py +++ b/src/dspeed/processors/dwt.py @@ -55,21 +55,21 @@ def discrete_wavelet_transform( "prereqs": ["wf_blsub"], } """ - + w_out[:] = np.nan - + if level <= 0: raise DSPFatal("The level must be a positive integer") if np.isnan(w_in).any(): - return - + return + coeff = chr(coeff) - - if chr(wave_type) == 'h': - wave_type = 'haar' - - elif chr(wave_type) == 'd': - wave_type = 'db1' - + + if chr(wave_type) == "h": + wave_type = "haar" + + elif chr(wave_type) == "d": + wave_type = "db1" + w_out[:] = downcoef(coeff, w_in, wave_type, level=level) diff --git a/src/dspeed/processors/min_max.py b/src/dspeed/processors/min_max.py index 45b5283..50462ac 100644 --- a/src/dspeed/processors/min_max.py +++ b/src/dspeed/processors/min_max.py @@ -112,7 +112,7 @@ def min_max_norm( """ w_out[:] = np.nan - + if np.isnan(w_in).any(): return diff --git a/tests/processors/test_dwt.py b/tests/processors/test_dwt.py index 515de9d..1a23eaf 100644 --- a/tests/processors/test_dwt.py +++ b/tests/processors/test_dwt.py @@ -10,9 +10,9 @@ def test_discrete_wavelet_transform(compare_numba_vs_python): # set up values to use for each test case len_wf_in = 16 - wave_type = ord('h') + wave_type = ord("h") level = 2 - coeff = ord('a') + coeff = ord("a") len_wf_out = 4 # ensure the DSPFatal is raised for a negative level @@ -24,7 +24,9 @@ def test_discrete_wavelet_transform(compare_numba_vs_python): # ensure that a valid input gives the expected output w_out_expected = np.ones(len_wf_out) * 2 ** (level / 2) assert np.allclose( - compare_numba_vs_python(discrete_wavelet_transform, w_in, level, wave_type, coeff, w_out), + compare_numba_vs_python( + discrete_wavelet_transform, w_in, level, wave_type, coeff, w_out + ), w_out_expected, ) @@ -37,4 +39,4 @@ def test_discrete_wavelet_transform(compare_numba_vs_python): discrete_wavelet_transform, w_in, level, wave_type, coeff, w_out ) ) - ) \ No newline at end of file + ) diff --git a/tests/processors/test_min_max_norm.py b/tests/processors/test_min_max_norm.py index 31d09ce..efda86c 100644 --- a/tests/processors/test_min_max_norm.py +++ b/tests/processors/test_min_max_norm.py @@ -13,35 +13,28 @@ def test_min_max_norm(compare_numba_vs_python): w_in = np.ones(len_wf) w_in[4] = np.nan a = np.array([1]) - assert np.all( - np.isnan( - compare_numba_vs_python(min_max_norm, w_in, a, a) - ) - ) - + assert np.all(np.isnan(compare_numba_vs_python(min_max_norm, w_in, a, a))) + # test for division by 0 w_in = np.ones(len_wf) a = np.array([0]) w_out_expected = np.ones(len_wf) assert np.allclose( - compare_numba_vs_python(min_max_norm, w_in, a, a), - w_out_expected + compare_numba_vs_python(min_max_norm, w_in, a, a), w_out_expected ) - + # test for abs(a_max) > abs(a_min) a_max = np.array([2]) a_min = np.array([-1]) - w_out_expected = np.ones(len_wf)/abs(a_max[0]) + w_out_expected = np.ones(len_wf) / abs(a_max[0]) assert np.allclose( - compare_numba_vs_python(min_max_norm, w_in, a_min, a_max), - w_out_expected + compare_numba_vs_python(min_max_norm, w_in, a_min, a_max), w_out_expected ) # test for abs(a_max) < abs(a_min) a_max = np.array([1]) a_min = np.array([-2]) - w_out_expected = np.ones(len_wf)/abs(a_min[0]) + w_out_expected = np.ones(len_wf) / abs(a_min[0]) assert np.allclose( - compare_numba_vs_python(min_max_norm, w_in, a_min, a_max), - w_out_expected + compare_numba_vs_python(min_max_norm, w_in, a_min, a_max), w_out_expected ) From fa6191c6e8c5641bdcc7b28d7c1c3cd39dfd6a30 Mon Sep 17 00:00:00 2001 From: Esteban Leon Date: Mon, 18 Dec 2023 18:57:19 -0800 Subject: [PATCH 35/42] Modified dwt processor output and test --- src/dspeed/processors/dwt.py | 2 +- tests/processors/test_dwt.py | 1 + 2 files changed, 2 insertions(+), 1 deletion(-) diff --git a/src/dspeed/processors/dwt.py b/src/dspeed/processors/dwt.py index 081526d..96098eb 100644 --- a/src/dspeed/processors/dwt.py +++ b/src/dspeed/processors/dwt.py @@ -13,7 +13,7 @@ "void(float32[:], int32, char, char, float32[:])", "void(float64[:], int64, char, char, float64[:])", ], - "(n),(),(),()->(m)", + "(n),(),(),(),(m)", **nb_kwargs( forceobj=True, ), diff --git a/tests/processors/test_dwt.py b/tests/processors/test_dwt.py index 1a23eaf..0a2cb2a 100644 --- a/tests/processors/test_dwt.py +++ b/tests/processors/test_dwt.py @@ -33,6 +33,7 @@ def test_discrete_wavelet_transform(compare_numba_vs_python): # ensure that if there is a nan in w_in, all nans are outputted w_in = np.ones(len_wf_in) w_in[4] = np.nan + w_out = np.empty(len_wf_out) assert np.all( np.isnan( compare_numba_vs_python( From 4b431e29abcca150f93a82208f43b3e3b93d6b63 Mon Sep 17 00:00:00 2001 From: "dependabot[bot]" <49699333+dependabot[bot]@users.noreply.github.com> Date: Mon, 1 Jan 2024 15:11:33 +0000 Subject: [PATCH 36/42] Bump actions/upload-artifact from 3 to 4 Bumps [actions/upload-artifact](https://github.com/actions/upload-artifact) from 3 to 4. - [Release notes](https://github.com/actions/upload-artifact/releases) - [Commits](https://github.com/actions/upload-artifact/compare/v3...v4) --- updated-dependencies: - dependency-name: actions/upload-artifact dependency-type: direct:production update-type: version-update:semver-major ... Signed-off-by: dependabot[bot] --- .github/workflows/distribution.yml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/.github/workflows/distribution.yml b/.github/workflows/distribution.yml index 8a702f3..57ab315 100644 --- a/.github/workflows/distribution.yml +++ b/.github/workflows/distribution.yml @@ -20,7 +20,7 @@ jobs: - name: Build SDist and wheel run: pipx run build - - uses: actions/upload-artifact@v3 + - uses: actions/upload-artifact@v4 with: path: dist/* From e7bf362d051c443f16c31f54dc4b0a808332efeb Mon Sep 17 00:00:00 2001 From: "dependabot[bot]" <49699333+dependabot[bot]@users.noreply.github.com> Date: Mon, 1 Jan 2024 15:11:37 +0000 Subject: [PATCH 37/42] Bump actions/download-artifact from 3 to 4 Bumps [actions/download-artifact](https://github.com/actions/download-artifact) from 3 to 4. - [Release notes](https://github.com/actions/download-artifact/releases) - [Commits](https://github.com/actions/download-artifact/compare/v3...v4) --- updated-dependencies: - dependency-name: actions/download-artifact dependency-type: direct:production update-type: version-update:semver-major ... Signed-off-by: dependabot[bot] --- .github/workflows/distribution.yml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/.github/workflows/distribution.yml b/.github/workflows/distribution.yml index 8a702f3..9800e2c 100644 --- a/.github/workflows/distribution.yml +++ b/.github/workflows/distribution.yml @@ -33,7 +33,7 @@ jobs: if: github.event_name == 'release' && github.event.action == 'published' steps: - - uses: actions/download-artifact@v3 + - uses: actions/download-artifact@v4 with: name: artifact path: dist From 6ba29afda629cbc0fd8b1595702bae70a090e788 Mon Sep 17 00:00:00 2001 From: Luigi Pertoldi Date: Mon, 1 Jan 2024 18:52:24 +0100 Subject: [PATCH 38/42] Update code to latest lgdo alpha version --- setup.cfg | 2 +- src/dspeed/build_dsp.py | 7 +++---- src/dspeed/processing_chain.py | 8 +++----- src/dspeed/processors/wiener_filter.py | 6 +++--- src/dspeed/vis/waveform_browser.py | 17 +++++++++-------- tests/conftest.py | 6 +++--- tests/processors/test_histogram.py | 4 ++-- tests/test_build_dsp.py | 4 ++-- tests/test_list_parsing.py | 4 ++-- tests/test_numpy_constants_parsing.py | 6 +++--- 10 files changed, 31 insertions(+), 33 deletions(-) diff --git a/setup.cfg b/setup.cfg index 0f3af95..9d653ee 100644 --- a/setup.cfg +++ b/setup.cfg @@ -35,7 +35,7 @@ install_requires = colorlog h5py>=3.2 iminuit - legend-pydataobj~=1.0 + legend-pydataobj>=1.5.0a1 matplotlib numba!=0.53.*,!=0.54.* parse diff --git a/src/dspeed/build_dsp.py b/src/dspeed/build_dsp.py index d834d98..e02b9f4 100644 --- a/src/dspeed/build_dsp.py +++ b/src/dspeed/build_dsp.py @@ -9,9 +9,8 @@ import os import h5py -import lgdo.lh5_store as lh5 import numpy as np -from lgdo.lgdo_utils import expand_path +from lgdo import lh5 from tqdm.auto import tqdm from .errors import DSPFatal @@ -133,7 +132,7 @@ def build_dsp( # get the database parameters. For now, this will just be a dict in a json # file, but eventually we will want to interface with the metadata repo if isinstance(database, str): - with open(expand_path(database)) as db_file: + with open(lh5.utils.expand_path(database)) as db_file: database = json.load(db_file) if database and not isinstance(database, dict): @@ -190,7 +189,7 @@ def build_dsp( e.wf_range = f"{e.wf_range[0]+start_row}-{e.wf_range[1]+start_row}" raise e - raw_store.write_object( + raw_store.write( obj=tb_out, name=tb_name, lh5_file=f_dsp, diff --git a/src/dspeed/processing_chain.py b/src/dspeed/processing_chain.py index d35e5dc..65d7098 100644 --- a/src/dspeed/processing_chain.py +++ b/src/dspeed/processing_chain.py @@ -16,10 +16,8 @@ from typing import Any import lgdo -import lgdo.lh5_store as lh5 import numpy as np -from lgdo import LGDO -from lgdo.lgdo_utils import expand_path +from lgdo import LGDO, lh5 from numba import vectorize from pint import Quantity, Unit @@ -1111,7 +1109,7 @@ def _loadlh5(path_to_file, path_in_file: str) -> np.array: # noqa: N805 """ try: - loaded_data = sto.read_object(path_in_file, path_to_file)[0].nda + loaded_data = sto.read(path_in_file, path_to_file)[0].nda except ValueError: raise ProcessingChainError(f"LH5 file not found: {path_to_file}") @@ -1863,7 +1861,7 @@ def build_processing_chain( proc_chain = ProcessingChain(block_width, lh5_in.size) if isinstance(dsp_config, str): - with open(expand_path(dsp_config)) as f: + with open(lh5.utils.expand_path(dsp_config)) as f: dsp_config = json.load(f) elif dsp_config is None: dsp_config = {"outputs": [], "processors": {}} diff --git a/src/dspeed/processors/wiener_filter.py b/src/dspeed/processors/wiener_filter.py index 060d669..868f5a1 100644 --- a/src/dspeed/processors/wiener_filter.py +++ b/src/dspeed/processors/wiener_filter.py @@ -1,7 +1,7 @@ from __future__ import annotations -import lgdo.lh5_store as lh5 import numpy as np +from lgdo import lh5 from numba import guvectorize from ..errors import DSPFatal @@ -65,10 +65,10 @@ def wiener_filter(file_name_array: list[str]) -> np.ndarray: # Read in the data - superpulse, _ = sto.read_object("spms/processed/superpulse", file_name) + superpulse, _ = sto.read("spms/processed/superpulse", file_name) superpulse = superpulse.nda - noise_wf, _ = sto.read_object("spms/processed/noise_wf", file_name) + noise_wf, _ = sto.read("spms/processed/noise_wf", file_name) noise_wf = noise_wf.nda # Now check that the data are valid diff --git a/src/dspeed/vis/waveform_browser.py b/src/dspeed/vis/waveform_browser.py index b1ce3a0..b16c874 100644 --- a/src/dspeed/vis/waveform_browser.py +++ b/src/dspeed/vis/waveform_browser.py @@ -5,12 +5,13 @@ import string import sys -import lgdo.lh5_store as lh5 +import lgdo import matplotlib.pyplot as plt import numpy as np import pandas import pint from cycler import cycler +from lgdo import lh5 from matplotlib.lines import Line2D from ..processing_chain import build_processing_chain @@ -28,7 +29,7 @@ class WaveformBrowser: def __init__( self, - files_in: str | list[str] | lgdo.LH5Iterator, # noqa: F821 + files_in: str | list[str] | lh5.LH5Iterator, # noqa: F821 lh5_group: str | list[str] = "", base_path: str = "", entry_list: list[int] | list[list[int]] = None, @@ -293,7 +294,7 @@ def __init__( # If we still have no x_unit get it from the first waveform we can find if self.x_unit is None: for wf in self.lh5_out.values(): - if not isinstance(wf, lh5.WaveformTable): + if not isinstance(wf, lgdo.WaveformTable): continue self.x_unit = ureg(wf.dt_units) @@ -390,7 +391,7 @@ def find_entry( ref_time = 0 elif isinstance(self.align_par, str): data = self.lh5_out.get(self.align_par, None) - if isinstance(data, lh5.Array): + if isinstance(data, lgdo.Array): ref_time = data.nda[i_tb] unit = data.attrs.get("units", None) if unit and unit in ureg and ureg.is_compatible_with(unit, self.x_unit): @@ -405,7 +406,7 @@ def find_entry( for name, lines in self.lines.items(): # Get the data; note this is implicitly copying it! data = self.lh5_out.get(name, None) - if isinstance(data, lh5.WaveformTable): + if isinstance(data, lgdo.WaveformTable): y = data.values.nda[i_tb, :] / norm dt = data.dt.nda[i_tb] * float(ureg(data.dt_units) / self.x_unit) t0 = ( @@ -416,13 +417,13 @@ def find_entry( lines.append(Line2D(x, y)) self._update_auto_limit(x, y) - elif isinstance(data, lh5.ArrayOfEqualSizedArrays): + elif isinstance(data, lgdo.ArrayOfEqualSizedArrays): y = data.nda[i_tb, :] / norm x = np.arange(len(y), dtype="float") lines.append(Line2D(x, y)) self._update_auto_limit(x, y) - elif isinstance(data, lh5.Array): + elif isinstance(data, lgdo.Array): val = data.nda[i_tb] unit = data.attrs.get("units", None) if unit and unit in ureg and ureg.is_compatible_with(unit, self.x_unit): @@ -453,7 +454,7 @@ def find_entry( if not data: data = ureg.Quantity(self.aux_vals[name][entry]) - elif isinstance(data, lh5.Array): + elif isinstance(data, lgdo.Array): unit = data.attrs.get("units", None) if unit and unit in ureg: data = data.nda[i_tb] * ureg(unit) diff --git a/tests/conftest.py b/tests/conftest.py index 8c515b5..71dc0a8 100644 --- a/tests/conftest.py +++ b/tests/conftest.py @@ -11,7 +11,7 @@ import numpy as np import pytest from legendtestdata import LegendTestData -from lgdo import LH5Store +from lgdo.lh5 import LH5Store import dspeed.processors # noqa: F401 @@ -42,7 +42,7 @@ def lgnd_test_data(): @pytest.fixture(scope="session") def geds_raw_tbl(lgnd_test_data): store = LH5Store() - obj, _ = store.read_object( + obj, _ = store.read( "/geds/raw", lgnd_test_data.get_path("lh5/LDQTA_r117_20200110T105115Z_cal_geds_raw.lh5"), n_rows=10, @@ -53,7 +53,7 @@ def geds_raw_tbl(lgnd_test_data): @pytest.fixture(scope="session") def spms_raw_tbl(lgnd_test_data): store = LH5Store() - obj, _ = store.read_object( + obj, _ = store.read( "/ch0/raw", lgnd_test_data.get_path("lh5/L200-comm-20211130-phy-spms.lh5"), n_rows=10, diff --git a/tests/processors/test_histogram.py b/tests/processors/test_histogram.py index 77aaeb9..af130b3 100644 --- a/tests/processors/test_histogram.py +++ b/tests/processors/test_histogram.py @@ -1,6 +1,6 @@ import os -import lgdo.lh5_store as store +from lgdo import lh5 from dspeed import build_dsp @@ -28,7 +28,7 @@ def test_histogram_fixed_width(lgnd_test_data, tmptestdir): ) assert os.path.exists(dsp_file) - df = store.load_nda(dsp_file, ["hist_weights", "hist_borders"], "geds/dsp/") + df = lh5.load_nda(dsp_file, ["hist_weights", "hist_borders"], "geds/dsp/") assert len(df["hist_weights"][0]) + 1 == len(df["hist_borders"][0]) for i in range(2, len(df["hist_borders"][0])): diff --git a/tests/test_build_dsp.py b/tests/test_build_dsp.py index 8ff600c..3248fe0 100644 --- a/tests/test_build_dsp.py +++ b/tests/test_build_dsp.py @@ -3,7 +3,7 @@ import lgdo import pytest -from lgdo.lh5_store import LH5Store, ls +from lgdo.lh5 import LH5Store, ls from dspeed import build_dsp @@ -72,7 +72,7 @@ def test_build_dsp_spms_channelwise(dsp_test_file_spm): ] store = LH5Store() - lh5_obj, n_rows = store.read_object("/ch0/dsp/energies", dsp_test_file_spm) + lh5_obj, n_rows = store.read("/ch0/dsp/energies", dsp_test_file_spm) assert isinstance(lh5_obj, lgdo.ArrayOfEqualSizedArrays) assert len(lh5_obj) == 5 assert len(lh5_obj.nda[0]) == 20 diff --git a/tests/test_list_parsing.py b/tests/test_list_parsing.py index b99d558..a0b11dd 100644 --- a/tests/test_list_parsing.py +++ b/tests/test_list_parsing.py @@ -1,8 +1,8 @@ import os from pathlib import Path -import lgdo.lh5_store as store import numpy as np +from lgdo import lh5 from dspeed import build_dsp @@ -33,6 +33,6 @@ def test_list_parisng(lgnd_test_data, tmptestdir): ) assert os.path.exists(dsp_file) - df = store.load_nda(dsp_file, ["wf_out"], "geds/dsp/") + df = lh5.load_nda(dsp_file, ["wf_out"], "geds/dsp/") assert np.all(df["wf_out"][:] == np.array([7, 9, 11, 13, 15])) diff --git a/tests/test_numpy_constants_parsing.py b/tests/test_numpy_constants_parsing.py index 09d3c4a..593e8b0 100644 --- a/tests/test_numpy_constants_parsing.py +++ b/tests/test_numpy_constants_parsing.py @@ -1,8 +1,8 @@ import os from pathlib import Path -import lgdo.lh5_store as store import numpy as np +from lgdo import lh5 from dspeed import build_dsp @@ -24,7 +24,7 @@ def test_build_dsp(lgnd_test_data, tmptestdir): def test_numpy_math_constants_dsp(tmptestdir): dsp_file = f"{tmptestdir}/LDQTA_r117_20200110T105115Z_cal_geds__numpy_test_dsp.lh5" - df = store.load_nda(dsp_file, ["timestamp", "calc1", "calc2", "calc3"], "geds/dsp/") + df = lh5.load_nda(dsp_file, ["timestamp", "calc1", "calc2", "calc3"], "geds/dsp/") a1 = df["timestamp"] - df["timestamp"] - np.pi * df["timestamp"] a2 = df["timestamp"] - df["timestamp"] - np.pi @@ -41,7 +41,7 @@ def test_numpy_math_constants_dsp(tmptestdir): def test_numpy_infinity_and_nan_dsp(tmptestdir): dsp_file = f"{tmptestdir}/LDQTA_r117_20200110T105115Z_cal_geds__numpy_test_dsp.lh5" - df = store.load_nda(dsp_file, ["calc4", "calc5", "calc6"], "geds/dsp/") + df = lh5.load_nda(dsp_file, ["calc4", "calc5", "calc6"], "geds/dsp/") assert (np.isnan(df["calc4"])).all() assert (np.isneginf(df["calc5"])).all() From ca7f2ff1bb3a5f88d3964adb2e04302dc2365871 Mon Sep 17 00:00:00 2001 From: Luigi Pertoldi Date: Mon, 1 Jan 2024 18:58:05 +0100 Subject: [PATCH 39/42] [docs] update Jupyter notebooks --- docs/source/notebooks/DSPTutorial.ipynb | 6 +++--- docs/source/notebooks/IntroToDSP.ipynb | 4 ++-- docs/source/notebooks/WaveformBrowser.ipynb | 4 ++-- 3 files changed, 7 insertions(+), 7 deletions(-) diff --git a/docs/source/notebooks/DSPTutorial.ipynb b/docs/source/notebooks/DSPTutorial.ipynb index 3d0dd90..546f90c 100644 --- a/docs/source/notebooks/DSPTutorial.ipynb +++ b/docs/source/notebooks/DSPTutorial.ipynb @@ -46,14 +46,14 @@ "import matplotlib.pyplot as plt\n", "import numpy as np\n", "\n", - "import lgdo.lh5_store as lh5\n", + "from lgdo import lh5\n", "from legendtestdata import LegendTestData\n", "\n", "# Get some sample waveforms from LEGEND test data\n", "ldata = LegendTestData()\n", "raw_file = ldata.get_path(\"lh5/LDQTA_r117_20200110T105115Z_cal_geds_raw.lh5\")\n", "st = lh5.LH5Store()\n", - "tab = st.read_object(\"geds/raw\", raw_file)[0]\n", + "tab = st.read(\"geds/raw\", raw_file)[0]\n", "wfs = tab[\"waveform\"].values.nda.astype(\"float32\")\n", "t = tab[\"waveform\"].dt.nda.reshape((100, 1)) * np.arange(wfs.shape[1])\n", "baselines = tab[\"baseline\"].nda\n", @@ -731,7 +731,7 @@ "\n", "lh5.show(\"test_dsp.lh5\")\n", "print()\n", - "print(st.read_object(\"geds/dsp\", \"test_dsp.lh5\")[0].get_dataframe())" + "print(st.read(\"geds/dsp\", \"test_dsp.lh5\")[0].get_dataframe())" ] }, { diff --git a/docs/source/notebooks/IntroToDSP.ipynb b/docs/source/notebooks/IntroToDSP.ipynb index 6794dcc..4b1fcfe 100644 --- a/docs/source/notebooks/IntroToDSP.ipynb +++ b/docs/source/notebooks/IntroToDSP.ipynb @@ -21,7 +21,7 @@ "metadata": {}, "outputs": [], "source": [ - "from lgdo import LH5Store, ls, show, load_dfs\n", + "from lgdo.lh5 import LH5Store, ls, show, load_dfs\n", "from dspeed.vis import WaveformBrowser\n", "from dspeed import units\n", "from dspeed import build_dsp\n", @@ -510,7 +510,7 @@ "outputs": [], "source": [ "store = LH5Store()\n", - "obj, nrows = store.read_object(\"geds/dsp\", \"./example_dsp_file.lh5\")\n", + "obj, nrows = store.read(\"geds/dsp\", \"./example_dsp_file.lh5\")\n", "print(\n", " f\"We have saved the following list of outputs from our DSP routine: {list(obj.keys())}\"\n", ")\n", diff --git a/docs/source/notebooks/WaveformBrowser.ipynb b/docs/source/notebooks/WaveformBrowser.ipynb index 57530d5..4a64099 100644 --- a/docs/source/notebooks/WaveformBrowser.ipynb +++ b/docs/source/notebooks/WaveformBrowser.ipynb @@ -29,7 +29,7 @@ "\n", "u = pint.get_application_registry()\n", "\n", - "import lgdo.lh5_store as lh5\n", + "from lgdo import lh5\n", "from dspeed.vis.waveform_browser import WaveformBrowser\n", "from legendtestdata import LegendTestData\n", "\n", @@ -132,7 +132,7 @@ "metadata": {}, "outputs": [], "source": [ - "from lgdo import load_dfs\n", + "from lgdo.lh5 import load_dfs\n", "\n", "df = load_dfs(dsp_file, [\"trapEmax\", \"AoE\"], \"geds/dsp\")" ] From e57ced3e4d7790db5b89d1d32dc91a53a24d9cf9 Mon Sep 17 00:00:00 2001 From: Luigi Pertoldi Date: Mon, 1 Jan 2024 19:18:00 +0100 Subject: [PATCH 40/42] [ci skip] forgot to update some docstrings --- src/dspeed/processing_chain.py | 4 ++-- src/dspeed/vis/waveform_browser.py | 2 +- 2 files changed, 3 insertions(+), 3 deletions(-) diff --git a/src/dspeed/processing_chain.py b/src/dspeed/processing_chain.py index 65d7098..14ab488 100644 --- a/src/dspeed/processing_chain.py +++ b/src/dspeed/processing_chain.py @@ -1778,8 +1778,8 @@ def build_processing_chain( block_width: int = 16, ) -> tuple[ProcessingChain, list[str], lgdo.Table]: """Produces a :class:`ProcessingChain` object and an LH5 - :class:`~lgdo.table.Table` for output parameters from an input LH5 - :class:`~lgdo.table.Table` and a JSON recipe. + :class:`~lgdo.types.table.Table` for output parameters from an input LH5 + :class:`~lgdo.types.table.Table` and a JSON recipe. Parameters ---------- diff --git a/src/dspeed/vis/waveform_browser.py b/src/dspeed/vis/waveform_browser.py index b16c874..abbb3e6 100644 --- a/src/dspeed/vis/waveform_browser.py +++ b/src/dspeed/vis/waveform_browser.py @@ -62,7 +62,7 @@ def __init__( and group, they must be the same size. If a file is wild-carded, the same group will be assigned to each file found base_path - base path for file. See :class:`~lgdo.lh5_store.LH5Store`. + base path for file. See :class:`~lgdo.lh5.LH5Store`. entry_list list of event indices to draw. If it is a nested list, use local From fddbcb4c49d15265519112f82b9e0fac329b7bd8 Mon Sep 17 00:00:00 2001 From: Luigi Pertoldi Date: Tue, 2 Jan 2024 14:35:01 +0100 Subject: [PATCH 41/42] Add missing SciPy dependency to setup.cfg --- setup.cfg | 1 + 1 file changed, 1 insertion(+) diff --git a/setup.cfg b/setup.cfg index 9d653ee..c65dedd 100644 --- a/setup.cfg +++ b/setup.cfg @@ -41,6 +41,7 @@ install_requires = parse pint pyfftw + scipy tqdm>=4.27 python_requires = >=3.9 include_package_data = True From 3dfed10d62f7f36f7b34c84fc9ca83b16e2f1e6a Mon Sep 17 00:00:00 2001 From: Luigi Pertoldi Date: Tue, 2 Jan 2024 14:43:51 +0100 Subject: [PATCH 42/42] Add stupid test that just imports all processors --- tests/processors/test_import.py | 5 +++++ 1 file changed, 5 insertions(+) create mode 100644 tests/processors/test_import.py diff --git a/tests/processors/test_import.py b/tests/processors/test_import.py new file mode 100644 index 0000000..0edf400 --- /dev/null +++ b/tests/processors/test_import.py @@ -0,0 +1,5 @@ +from dspeed.processors import * # noqa: F403, F401 + + +def test_import(): + pass