From e3f722081b4e73b05588b242e3f61578c1742bbb Mon Sep 17 00:00:00 2001
From: Felix Hagemann <hagemann@mpp.mpg.de>
Date: Mon, 16 Dec 2024 14:46:23 +0100
Subject: [PATCH 1/3] Add `DerivativeFilter`

---
 src/LegendDSP.jl  |  4 +++-
 src/derivative.jl | 57 +++++++++++++++++++++++++++++++++++++++++++++++
 2 files changed, 60 insertions(+), 1 deletion(-)
 create mode 100644 src/derivative.jl

diff --git a/src/LegendDSP.jl b/src/LegendDSP.jl
index d3a2d0c..20dbe48 100644
--- a/src/LegendDSP.jl
+++ b/src/LegendDSP.jl
@@ -27,13 +27,15 @@ using LegendDataTypes
 
 import Adapt
 import RadiationDetectorDSP: fltinstance, rdfilt!, flt_output_length, 
-    flt_input_smpltype, flt_output_smpltype
+    flt_input_smpltype, flt_output_smpltype, _smpllen, smplinfo, SamplingInfo,
+    AbstractRadSigFilterInstance, LinearFiltering, _filterlen, _floattype
 
 using Unitful: RealOrRealQuantity as RealQuantity
 
 include("tailstats.jl")
 include("types.jl")
 include("utils.jl")
+include("derivative.jl")
 include("dsp_decaytime.jl")
 include("dsp_filter_optimization.jl")
 include("dsp_icpc.jl")
diff --git a/src/derivative.jl b/src/derivative.jl
new file mode 100644
index 0000000..6fa2888
--- /dev/null
+++ b/src/derivative.jl
@@ -0,0 +1,57 @@
+# This file is a part of LegendDSP.jl, licensed under the MIT License (MIT).
+
+"""
+    struct DerivativeFilter <: AbstractRadIIRFilter
+
+Working example:
+```julia
+using LegendDSP
+using Unitful
+
+signal = rand(100)
+t = range(0u"ms", 20u"ms", 100)
+
+wf = RDWaveform(t, signal)
+flt = DerivativeFilter()
+wf_new = flt(wf)
+```
+Constructors:
+
+* ```$(FUNCTIONNAME)(; fields...)```
+
+Fields:
+
+$(TYPEDFIELDS)
+"""
+Base.@kwdef struct DerivativeFilter{T <: RadiationDetectorDSP.RealQuantity} <: AbstractRadIIRFilter
+    "Filter gain"
+    gain::T = 1
+end
+export DerivativeFilter
+
+struct DerivativeFilterInstance{T <: RadiationDetectorDSP.RealQuantity, TG <: RadiationDetectorDSP.RealQuantity} <: AbstractRadSigFilterInstance{LinearFiltering}
+    gain::TG
+    n_input::Int
+end
+
+function RadiationDetectorDSP.fltinstance(flt::DerivativeFilter{TG}, si::SamplingInfo{T}) where {T, TG}
+    DerivativeFilterInstance{T, TG}(flt.gain, _smpllen(si))
+end
+
+
+RadiationDetectorDSP._filterlen(fi::DerivativeFilterInstance) = 1
+RadiationDetectorDSP.flt_output_smpltype(::DerivativeFilterInstance{T, TG})  where {T,TG} = promote_type(typeof(zero(T) * unit(TG)), TG)
+RadiationDetectorDSP.flt_input_smpltype(::DerivativeFilterInstance{T}) where {T} = T
+RadiationDetectorDSP.flt_output_length(fi::DerivativeFilterInstance) = fi.n_input
+RadiationDetectorDSP.flt_input_length(fi::DerivativeFilterInstance) = fi.n_input
+RadiationDetectorDSP.flt_output_time_axis(::DerivativeFilterInstance, time::AbstractVector) = time
+
+function rdfilt!(y::AbstractVector, fi::DerivativeFilterInstance, x::AbstractVector)
+
+    @assert firstindex(y) == firstindex(x)
+    @assert lastindex(y) == lastindex(x)
+    @inbounds @simd for i in eachindex(y)
+        y[i] = fi.gain * (x[max(i,begin+1)] - x[max(i-1,begin)])
+    end
+    y
+end
\ No newline at end of file

From 6a0d463e9a2b910b9577c2c0e37f6b1d64653516 Mon Sep 17 00:00:00 2001
From: Felix Hagemann <hagemann@mpp.mpg.de>
Date: Mon, 16 Dec 2024 14:46:34 +0100
Subject: [PATCH 2/3] Add tests for `DerivativeFilter`

---
 test/test_derivative.jl | 32 ++++++++++++++++++++++++++++++++
 1 file changed, 32 insertions(+)
 create mode 100644 test/test_derivative.jl

diff --git a/test/test_derivative.jl b/test/test_derivative.jl
new file mode 100644
index 0000000..a081537
--- /dev/null
+++ b/test/test_derivative.jl
@@ -0,0 +1,32 @@
+using Test
+using LegendDSP
+using RadiationDetectorSignals: RDWaveform
+using Unitful
+
+@testset "Test DerivativeFilter" begin
+    T = Float64
+    sig = rand(100)
+    wvf = RDWaveform(range(0u"ns", step = 16u"ns", length = length(sig)), sig)
+
+    flt = DerivativeFilter()
+    wvf_out = flt(wvf)
+    @test isapprox(wvf_out.signal, vcat(sig[2] - sig[1], diff(sig)...))
+
+    gain = rand()
+    flt = DerivativeFilter(gain)
+    wvf_out = flt(wvf)
+    @test isapprox(wvf_out.signal, gain * vcat(sig[2] - sig[1], diff(sig)...))
+
+    @testset "Test units and precision types" begin
+        for T1 in (Float16, Float32, Float64)
+            for T2 in (Float16, Float32, Float64)
+                sig = rand(T1, 100)
+                wvf = RDWaveform(range(0u"ns", step = 16u"ns", length = length(sig)), sig)
+                gain = rand(T1) * u"keV"
+                flt = DerivativeFilter(gain)
+                wvf_out = flt(wvf)
+                @test isapprox(wvf_out.signal, gain * vcat(sig[2] - sig[1], diff(sig)...))
+            end
+        end
+    end
+end
\ No newline at end of file

From 374131f96f363113b314a3efb11d2b33c184a9df Mon Sep 17 00:00:00 2001
From: Felix Hagemann <hagemann@mpp.mpg.de>
Date: Mon, 16 Dec 2024 15:52:27 +0100
Subject: [PATCH 3/3] Remove `_filterlen` for `DerivativeFilter`

---
 src/derivative.jl | 2 --
 1 file changed, 2 deletions(-)

diff --git a/src/derivative.jl b/src/derivative.jl
index 6fa2888..4af309b 100644
--- a/src/derivative.jl
+++ b/src/derivative.jl
@@ -38,8 +38,6 @@ function RadiationDetectorDSP.fltinstance(flt::DerivativeFilter{TG}, si::Samplin
     DerivativeFilterInstance{T, TG}(flt.gain, _smpllen(si))
 end
 
-
-RadiationDetectorDSP._filterlen(fi::DerivativeFilterInstance) = 1
 RadiationDetectorDSP.flt_output_smpltype(::DerivativeFilterInstance{T, TG})  where {T,TG} = promote_type(typeof(zero(T) * unit(TG)), TG)
 RadiationDetectorDSP.flt_input_smpltype(::DerivativeFilterInstance{T}) where {T} = T
 RadiationDetectorDSP.flt_output_length(fi::DerivativeFilterInstance) = fi.n_input