From 1894b10b22c46a1d4700177119cd7f40757b57a3 Mon Sep 17 00:00:00 2001 From: Martin Holters Date: Mon, 2 Dec 2024 14:48:52 +0100 Subject: [PATCH 1/6] Support multi-column filtering for `DF2TFilter` Fixes #372. --- src/Filters/filt.jl | 66 +++++++++++++++++++++++++++------------------ src/dspbase.jl | 4 +-- test/filt.jl | 20 ++++++++++++++ 3 files changed, 62 insertions(+), 28 deletions(-) diff --git a/src/Filters/filt.jl b/src/Filters/filt.jl index 62551161e..80ab93a83 100644 --- a/src/Filters/filt.jl +++ b/src/Filters/filt.jl @@ -136,37 +136,39 @@ struct DF2TFilter{T<:FilterCoefficients{:z},S<:Array} coef::T state::S - function DF2TFilter{Ti,Si}(coef::PolynomialRatio{:z}, state::Vector) where {Ti,Si} - length(state) == max(length(coefa(coef)), length(coefb(coef)))-1 || + function DF2TFilter{Ti,Si}(coef::PolynomialRatio{:z}, state::Array) where {Ti,Si} + size(state, 1) == max(length(coefa(coef)), length(coefb(coef)))-1 || throw(ArgumentError("length of state vector must match filter order")) new{Ti,Si}(coef, state) end - function DF2TFilter{Ti,Si}(coef::SecondOrderSections{:z}, state::Matrix) where {Ti,Si} + function DF2TFilter{Ti,Si}(coef::SecondOrderSections{:z}, state::Array) where {Ti,Si} (size(state, 1) == 2 && size(state, 2) == length(coef.biquads)) || throw(ArgumentError("state must be 2 x nbiquads")) new{Ti,Si}(coef, state) end - function DF2TFilter{Ti,Si}(coef::Biquad{:z}, state::Vector) where {Ti,Si} - length(state) == 2 || throw(ArgumentError("length of state must be 2")) + function DF2TFilter{Ti,Si}(coef::Biquad{:z}, state::Array) where {Ti,Si} + size(state,1) == 2 || throw(ArgumentError("length of state must be 2")) new{Ti,Si}(coef, state) end end -DF2TFilter(coef::Union{PolynomialRatio{:z,T},Biquad{:z,T}}, state::Vector{S}) where {T,S} = - DF2TFilter{typeof(coef),Vector{S}}(coef, state) +DF2TFilter(coef::Union{PolynomialRatio{:z,T},Biquad{:z,T}}, state::S) where {T,S<:Array} = + DF2TFilter{typeof(coef),S}(coef, state) -DF2TFilter(coef::SecondOrderSections{:z,T,G}, state::Matrix{S}) where {T,G,S} = - DF2TFilter{SecondOrderSections{:z,T,G},Matrix{S}}(coef, state) +DF2TFilter(coef::SecondOrderSections{:z,T,G}, state::S) where {T,G,S<:Array} = + DF2TFilter{SecondOrderSections{:z,T,G},S}(coef, state) ## PolynomialRatio DF2TFilter(coef::PolynomialRatio{:z,T}, ::Type{V}=T) where {T,V} = DF2TFilter(coef, zeros(promote_type(T, V), max(length(coefa(coef)), length(coefb(coef))) - 1)) -function filt!(out::AbstractVector, f::DF2TFilter{<:PolynomialRatio,<:Vector}, x::AbstractVector) - length(x) != length(out) && throw(ArgumentError("out size must match x")) +function filt!(out::AbstractArray{<:Any, N}, f::DF2TFilter{<:PolynomialRatio,Array{T,N}} where T, x::AbstractArray{<:Any, N}) where N + size(x) != size(out) && throw(ArgumentError("out size must match x")) si = f.state - n = length(si) + 1 + size(x)[2:end] != size(si)[2:end] && throw(ArgumentError("state size must match x")) + + n = size(si, 1) + 1 b = coefb(f.coef) if n == 1 mul!(out, x, b[1]) @@ -180,12 +182,18 @@ function filt!(out::AbstractVector, f::DF2TFilter{<:PolynomialRatio,<:Vector}, x elseif bs < n append!(b, zero(eltype(b)) for _ in 1:(n-bs)) end - _filt_iir!(out, b, a, x, si, 1) + for col in CartesianIndices(axes(x)[2:end]) + _filt_iir!(out, b, a, x, view(si, :, col), col) + end elseif n <= SMALL_FILT_CUTOFF vtup = ntuple(j -> b[j], Val(length(b))) - si .= getfield.(_filt_fir!(out, vtup, x, si, Val(:DF2)), :value) + for col in CartesianIndices(axes(x)[2:end]) + si .= getfield.(_filt_fir!(out, vtup, x, view(si, :, col), col, Val(true)), :value) + end else - _filt_fir!(out, b, x, si, 1) + for col in CartesianIndices(axes(x)[2:end]) + _filt_fir!(out, b, x, view(si, :, col), col) + end end end return out @@ -195,25 +203,31 @@ end DF2TFilter(coef::SecondOrderSections{:z,T,G}, ::Type{V}=T) where {T,G,V} = DF2TFilter(coef, zeros(promote_type(T, G, V), 2, length(coef.biquads))) -function filt!(out::AbstractVector, f::DF2TFilter{<:SecondOrderSections,<:Matrix}, x::AbstractVector) - length(x) != length(out) && throw(ArgumentError("out size must match x")) - _filt!(out, f.state, f.coef, x, 1) - out +function filt!(out::AbstractArray{<:Any, N}, f::DF2TFilter{<:SecondOrderSections,<:Array}, x::AbstractArray{<:Any, N}) where N + size(x) != size(out) && throw(ArgumentError("out size must match x")) + size(x)[2:end] != size(f.state)[3:end] && throw(ArgumentError("state size must match x")) + for col in CartesianIndices(axes(x)[2:end]) + _filt!(out, view(f.state, :, :, col), f.coef, x, col) + end + return out end ## Biquad DF2TFilter(coef::Biquad{:z,T}, ::Type{V}=T) where {T,V} = DF2TFilter(coef, zeros(promote_type(T, V), 2)) -function filt!(out::AbstractVector, f::DF2TFilter{<:Biquad,<:Vector}, x::AbstractVector) - length(x) != length(out) && throw(ArgumentError("out size must match x")) +function filt!(out::AbstractArray{<:Any, N}, f::DF2TFilter{<:Biquad,Array{T, N}} where T, x::AbstractArray{<:Any, N}) where N + size(x) != size(out) && throw(ArgumentError("out size must match x")) si = f.state - (si[1], si[2]) = _filt!(out, si[1], si[2], f.coef, x, 1) - out + size(x)[2:end] != size(si)[2:end] && throw(ArgumentError("state size must match x")) + for col in CartesianIndices(axes(x)[2:end]) + (si[1, col], si[2, col]) = _filt!(out, si[1, col], si[2, col], f.coef, x, col) + end + return out end """ - filt(f::DF2TFilter{<:FilterCoefficients{:z},<:Array{T}}, x::AbstractVector{V}) where {T,V} + filt(f::DF2TFilter{<:FilterCoefficients{:z},<:Array{T}}, x::AbstractArray{V}) where {T,V} Apply the [stateful filter](@ref stateful-filter-objects) `f` on `x`. @@ -223,12 +237,12 @@ Apply the [stateful filter](@ref stateful-filter-objects) `f` on `x`. For more control over the output type, provide a preallocated output array `out` to `filt!(out, f, x)`. """ -filt(f::DF2TFilter{<:FilterCoefficients{:z},<:Array{T}}, x::AbstractVector{V}) where {T,V} = +filt(f::DF2TFilter{<:FilterCoefficients{:z},<:Array{T}}, x::AbstractArray{V}) where {T,V} = filt!(similar(x, promote_type(T, V)), f, x) # Fall back to SecondOrderSections DF2TFilter(coef::FilterCoefficients{:z}) = DF2TFilter(convert(SecondOrderSections, coef)) -DF2TFilter(coef::FilterCoefficients{:z}, arg::Union{Matrix,Type}) = +DF2TFilter(coef::FilterCoefficients{:z}, arg::Union{Array,Type}) = DF2TFilter(convert(SecondOrderSections, coef), arg) # diff --git a/src/dspbase.jl b/src/dspbase.jl index 49921bd7e..3cd4e094f 100644 --- a/src/dspbase.jl +++ b/src/dspbase.jl @@ -118,7 +118,7 @@ end const SMALL_FILT_VECT_CUTOFF = 19 # Transposed direct form II -@generated function _filt_fir!(out, b::NTuple{N,T}, x, siarr, colv) where {N,T} +@generated function _filt_fir!(out, b::NTuple{N,T}, x, siarr, colv, _::Val{AsVecElem}=Val(false)) where {N,T,AsVecElem} silen = N - 1 si_end = Symbol(:si_, silen) @@ -137,7 +137,7 @@ const SMALL_FILT_VECT_CUTOFF = 19 out[i, col] = val end end - if colv isa Val{:DF2} + if AsVecElem return Base.@ntuple $silen j -> VecElement(si_j) end end diff --git a/test/filt.jl b/test/filt.jl index 73e5cf1b6..91ba60385 100644 --- a/test/filt.jl +++ b/test/filt.jl @@ -89,6 +89,26 @@ end @test all(col -> col ≈ y_ref, eachslice(filt(PolynomialRatio(b, a), x); dims=slicedims)) end +@testset "multi-column DF2TFilter $D-D" for D in 1:4 + b = [0.1, 0.1] + a = [1.0, -0.8] + sz = (10, ntuple(n -> n+1, Val(D))...) + y_ref = filt(b, a, ones(2*sz[1])) + x = ones(sz) + + H = DF2TFilter(PolynomialRatio(b, a), zeros(1, sz[2:end]...)) + @test all(col -> col ≈ y_ref[1:sz[1]], eachslice(filt(H, x); dims=ntuple(n -> n+1, Val(D)))) + @test all(col -> col ≈ y_ref[sz[1]+1:end], eachslice(filt(H, x); dims=ntuple(n -> n+1, Val(D)))) + + H = DF2TFilter(SecondOrderSections(PolynomialRatio(b, a)), zeros(2, 1, sz[2:end]...)) + @test all(col -> col ≈ y_ref[1:sz[1]], eachslice(filt(H, x); dims=ntuple(n -> n+1, Val(D)))) + @test all(col -> col ≈ y_ref[sz[1]+1:end], eachslice(filt(H, x); dims=ntuple(n -> n+1, Val(D)))) + + H = DF2TFilter(Biquad(PolynomialRatio(b, a)), zeros(2, sz[2:end]...)) + @test all(col -> col ≈ y_ref[1:sz[1]], eachslice(filt(H, x); dims=ntuple(n -> n+1, Val(D)))) + @test all(col -> col ≈ y_ref[sz[1]+1:end], eachslice(filt(H, x); dims=ntuple(n -> n+1, Val(D)))) +end + # # filtfilt # From c9b44d71bc3e4cc3d721c916c8eeb85ee585daf2 Mon Sep 17 00:00:00 2001 From: Martin Holters Date: Mon, 2 Dec 2024 15:30:21 +0100 Subject: [PATCH 2/6] Support `coldims` in `DF2TFilter` constructor --- src/Filters/filt.jl | 25 ++++++++++++++----------- test/filt.jl | 6 +++--- 2 files changed, 17 insertions(+), 14 deletions(-) diff --git a/src/Filters/filt.jl b/src/Filters/filt.jl index 80ab93a83..9e1bcb039 100644 --- a/src/Filters/filt.jl +++ b/src/Filters/filt.jl @@ -117,14 +117,16 @@ filt!(out, f::FilterCoefficients{:z}, x) = filt!(out, convert(SecondOrderSection """ DF2TFilter(coef::FilterCoefficients{:z}[, si]) - DF2TFilter(coef::FilterCoefficients{:z}[, sitype::Type]) + DF2TFilter(coef::FilterCoefficients{:z}[, sitype::Type][, coldims::Tuple]) Construct a stateful direct form II transposed filter with coefficients `coef`. One can optionally specify as the second argument either - `si`, an array representing the initial filter state, or -- `sitype`, the eltype of a zeroed `si` +- `sitype`, the eltype of a zeroed `si` and/or `coldims`, the size of extra dimensions + of the input. E.g. to column-wise filter an input with dims `(L, N1, N2)`, set + `coldims` to `(N1, N2)`. The initial filter state defaults to zeros if called with one argument. @@ -159,8 +161,9 @@ DF2TFilter(coef::SecondOrderSections{:z,T,G}, state::S) where {T,G,S<:Array} = DF2TFilter{SecondOrderSections{:z,T,G},S}(coef, state) ## PolynomialRatio -DF2TFilter(coef::PolynomialRatio{:z,T}, ::Type{V}=T) where {T,V} = - DF2TFilter(coef, zeros(promote_type(T, V), max(length(coefa(coef)), length(coefb(coef))) - 1)) +DF2TFilter(coef::PolynomialRatio{:z,T}, coldims::Tuple{Vararg{Integer}}) where {T} = DF2TFilter(coef, T, coldims) +DF2TFilter(coef::PolynomialRatio{:z,T}, ::Type{V}=T, coldims::Tuple{Vararg{Integer}}=()) where {T,V} = + DF2TFilter(coef, zeros(promote_type(T, V), max(length(coefa(coef)), length(coefb(coef))) - 1, coldims...)) function filt!(out::AbstractArray{<:Any, N}, f::DF2TFilter{<:PolynomialRatio,Array{T,N}} where T, x::AbstractArray{<:Any, N}) where N size(x) != size(out) && throw(ArgumentError("out size must match x")) @@ -200,8 +203,9 @@ function filt!(out::AbstractArray{<:Any, N}, f::DF2TFilter{<:PolynomialRatio,Arr end ## SecondOrderSections -DF2TFilter(coef::SecondOrderSections{:z,T,G}, ::Type{V}=T) where {T,G,V} = - DF2TFilter(coef, zeros(promote_type(T, G, V), 2, length(coef.biquads))) +DF2TFilter(coef::SecondOrderSections{:z,T}, coldims::Tuple{Vararg{Integer}}) where {T} = DF2TFilter(coef, T, coldims) +DF2TFilter(coef::SecondOrderSections{:z,T,G}, ::Type{V}=T, coldims::Tuple{Vararg{Integer}}=()) where {T,G,V} = + DF2TFilter(coef, zeros(promote_type(T, G, V), 2, length(coef.biquads), coldims...)) function filt!(out::AbstractArray{<:Any, N}, f::DF2TFilter{<:SecondOrderSections,<:Array}, x::AbstractArray{<:Any, N}) where N size(x) != size(out) && throw(ArgumentError("out size must match x")) @@ -213,8 +217,9 @@ function filt!(out::AbstractArray{<:Any, N}, f::DF2TFilter{<:SecondOrderSections end ## Biquad -DF2TFilter(coef::Biquad{:z,T}, ::Type{V}=T) where {T,V} = - DF2TFilter(coef, zeros(promote_type(T, V), 2)) +DF2TFilter(coef::Biquad{:z,T}, coldims::Tuple{Vararg{Integer}}) where {T} = DF2TFilter(coef, T, coldims) +DF2TFilter(coef::Biquad{:z,T}, ::Type{V}=T, coldims::Tuple{Vararg{Integer}}=()) where {T,V} = + DF2TFilter(coef, zeros(promote_type(T, V), 2, coldims...)) function filt!(out::AbstractArray{<:Any, N}, f::DF2TFilter{<:Biquad,Array{T, N}} where T, x::AbstractArray{<:Any, N}) where N size(x) != size(out) && throw(ArgumentError("out size must match x")) @@ -241,9 +246,7 @@ filt(f::DF2TFilter{<:FilterCoefficients{:z},<:Array{T}}, x::AbstractArray{V}) wh filt!(similar(x, promote_type(T, V)), f, x) # Fall back to SecondOrderSections -DF2TFilter(coef::FilterCoefficients{:z}) = DF2TFilter(convert(SecondOrderSections, coef)) -DF2TFilter(coef::FilterCoefficients{:z}, arg::Union{Array,Type}) = - DF2TFilter(convert(SecondOrderSections, coef), arg) +DF2TFilter(coef::FilterCoefficients{:z}, args...) = DF2TFilter(convert(SecondOrderSections, coef), args...) # # filtfilt diff --git a/test/filt.jl b/test/filt.jl index 91ba60385..faba09c88 100644 --- a/test/filt.jl +++ b/test/filt.jl @@ -96,15 +96,15 @@ end y_ref = filt(b, a, ones(2*sz[1])) x = ones(sz) - H = DF2TFilter(PolynomialRatio(b, a), zeros(1, sz[2:end]...)) + H = DF2TFilter(PolynomialRatio(b, a),sz[2:end]) @test all(col -> col ≈ y_ref[1:sz[1]], eachslice(filt(H, x); dims=ntuple(n -> n+1, Val(D)))) @test all(col -> col ≈ y_ref[sz[1]+1:end], eachslice(filt(H, x); dims=ntuple(n -> n+1, Val(D)))) - H = DF2TFilter(SecondOrderSections(PolynomialRatio(b, a)), zeros(2, 1, sz[2:end]...)) + H = DF2TFilter(SecondOrderSections(PolynomialRatio(b, a)), sz[2:end]) @test all(col -> col ≈ y_ref[1:sz[1]], eachslice(filt(H, x); dims=ntuple(n -> n+1, Val(D)))) @test all(col -> col ≈ y_ref[sz[1]+1:end], eachslice(filt(H, x); dims=ntuple(n -> n+1, Val(D)))) - H = DF2TFilter(Biquad(PolynomialRatio(b, a)), zeros(2, sz[2:end]...)) + H = DF2TFilter(Biquad(PolynomialRatio(b, a)), sz[2:end]) @test all(col -> col ≈ y_ref[1:sz[1]], eachslice(filt(H, x); dims=ntuple(n -> n+1, Val(D)))) @test all(col -> col ≈ y_ref[sz[1]+1:end], eachslice(filt(H, x); dims=ntuple(n -> n+1, Val(D)))) end From 5880dc3aa164cd14a3aae21d0f2d9cfac5dde5ab Mon Sep 17 00:00:00 2001 From: wheeheee <104880306+wheeheee@users.noreply.github.com> Date: Tue, 3 Dec 2024 00:40:02 +0800 Subject: [PATCH 3/6] directly store inside `_filt_fir!` does not have performance repercussions for Vector and Matrix states --- src/Filters/filt.jl | 2 +- src/dspbase.jl | 8 ++++---- 2 files changed, 5 insertions(+), 5 deletions(-) diff --git a/src/Filters/filt.jl b/src/Filters/filt.jl index 9e1bcb039..0e63f392b 100644 --- a/src/Filters/filt.jl +++ b/src/Filters/filt.jl @@ -191,7 +191,7 @@ function filt!(out::AbstractArray{<:Any, N}, f::DF2TFilter{<:PolynomialRatio,Arr elseif n <= SMALL_FILT_CUTOFF vtup = ntuple(j -> b[j], Val(length(b))) for col in CartesianIndices(axes(x)[2:end]) - si .= getfield.(_filt_fir!(out, vtup, x, view(si, :, col), col, Val(true)), :value) + _filt_fir!(out, vtup, x, view(si, :, col), col, Val(true)) end else for col in CartesianIndices(axes(x)[2:end]) diff --git a/src/dspbase.jl b/src/dspbase.jl index 3cd4e094f..c2135e83a 100644 --- a/src/dspbase.jl +++ b/src/dspbase.jl @@ -118,12 +118,11 @@ end const SMALL_FILT_VECT_CUTOFF = 19 # Transposed direct form II -@generated function _filt_fir!(out, b::NTuple{N,T}, x, siarr, colv, _::Val{AsVecElem}=Val(false)) where {N,T,AsVecElem} +@generated function _filt_fir!(out, b::NTuple{N,T}, x, siarr, col, ::Val{StoreSI}=Val(false)) where {N,T,StoreSI} silen = N - 1 si_end = Symbol(:si_, silen) quote - col = colv isa Val{:DF2} ? CartesianIndex() : colv N <= SMALL_FILT_VECT_CUTOFF && checkbounds(siarr, $silen) Base.@nextract $silen si siarr for i in axes(x, 1) @@ -137,9 +136,10 @@ const SMALL_FILT_VECT_CUTOFF = 19 out[i, col] = val end end - if AsVecElem - return Base.@ntuple $silen j -> VecElement(si_j) + if StoreSI + Base.@nexprs $silen j -> siarr[j] = si_j end + return nothing end end From 5bfcfb0e801700a476b78a8e13c6b4b62109282a Mon Sep 17 00:00:00 2001 From: wheeheee <104880306+wheeheee@users.noreply.github.com> Date: Tue, 3 Dec 2024 12:54:36 +0800 Subject: [PATCH 4/6] prevent stack overflow --- src/Filters/filt.jl | 5 ++++- test/filt.jl | 3 +++ 2 files changed, 7 insertions(+), 1 deletion(-) diff --git a/src/Filters/filt.jl b/src/Filters/filt.jl index 0e63f392b..2ca81d535 100644 --- a/src/Filters/filt.jl +++ b/src/Filters/filt.jl @@ -246,7 +246,10 @@ filt(f::DF2TFilter{<:FilterCoefficients{:z},<:Array{T}}, x::AbstractArray{V}) wh filt!(similar(x, promote_type(T, V)), f, x) # Fall back to SecondOrderSections -DF2TFilter(coef::FilterCoefficients{:z}, args...) = DF2TFilter(convert(SecondOrderSections, coef), args...) +DF2TFilter(coef::FilterCoefficients{:z}, coldims::Tuple{Vararg{Integer}}=()) = + DF2TFilter(convert(SecondOrderSections, coef), coldims) +DF2TFilter(coef::FilterCoefficients{:z}, ::Type{V}, coldims::Tuple{Vararg{Integer}}=()) where {V} = + DF2TFilter(convert(SecondOrderSections, coef), V, coldims) # # filtfilt diff --git a/test/filt.jl b/test/filt.jl index faba09c88..79473b500 100644 --- a/test/filt.jl +++ b/test/filt.jl @@ -63,6 +63,9 @@ using DSP, Test, Random, FilterTestHelpers df = digitalfilter(Lowpass(0.25), Butterworth(4)) f = @test_nowarn DF2TFilter(df, ComplexF64) @test_nowarn filt(f, s) + + # DF2TFilter{ZPG/SOS} stackoverflow error + @test_throws MethodError DF2TFilter(ZeroPoleGain([1], [2], 3), :D, 'x', Ref(1)) end @testset "multi-column filt $D-D" for D in 1:4 From e3aba8204de22c920d267c9b1f15c51cf36f279f Mon Sep 17 00:00:00 2001 From: Martin Holters Date: Tue, 3 Dec 2024 11:14:15 +0100 Subject: [PATCH 5/6] Apply cosmetic changes suggested in code review Co-authored-by: wheeheee <104880306+wheeheee@users.noreply.github.com> --- src/Filters/filt.jl | 10 +++++----- 1 file changed, 5 insertions(+), 5 deletions(-) diff --git a/src/Filters/filt.jl b/src/Filters/filt.jl index 2ca81d535..58b38c742 100644 --- a/src/Filters/filt.jl +++ b/src/Filters/filt.jl @@ -144,12 +144,12 @@ struct DF2TFilter{T<:FilterCoefficients{:z},S<:Array} new{Ti,Si}(coef, state) end function DF2TFilter{Ti,Si}(coef::SecondOrderSections{:z}, state::Array) where {Ti,Si} - (size(state, 1) == 2 && size(state, 2) == length(coef.biquads)) || + size(state, 1) == 2 && size(state, 2) == length(coef.biquads) || throw(ArgumentError("state must be 2 x nbiquads")) new{Ti,Si}(coef, state) end function DF2TFilter{Ti,Si}(coef::Biquad{:z}, state::Array) where {Ti,Si} - size(state,1) == 2 || throw(ArgumentError("length of state must be 2")) + size(state, 1) == 2 || throw(ArgumentError("length of state must be 2")) new{Ti,Si}(coef, state) end end @@ -165,7 +165,7 @@ DF2TFilter(coef::PolynomialRatio{:z,T}, coldims::Tuple{Vararg{Integer}}) where { DF2TFilter(coef::PolynomialRatio{:z,T}, ::Type{V}=T, coldims::Tuple{Vararg{Integer}}=()) where {T,V} = DF2TFilter(coef, zeros(promote_type(T, V), max(length(coefa(coef)), length(coefb(coef))) - 1, coldims...)) -function filt!(out::AbstractArray{<:Any, N}, f::DF2TFilter{<:PolynomialRatio,Array{T,N}} where T, x::AbstractArray{<:Any, N}) where N +function filt!(out::AbstractArray{<:Any,N}, f::DF2TFilter{<:PolynomialRatio,Array{T,N}} where T, x::AbstractArray{<:Any,N}) where N size(x) != size(out) && throw(ArgumentError("out size must match x")) si = f.state @@ -207,7 +207,7 @@ DF2TFilter(coef::SecondOrderSections{:z,T}, coldims::Tuple{Vararg{Integer}}) whe DF2TFilter(coef::SecondOrderSections{:z,T,G}, ::Type{V}=T, coldims::Tuple{Vararg{Integer}}=()) where {T,G,V} = DF2TFilter(coef, zeros(promote_type(T, G, V), 2, length(coef.biquads), coldims...)) -function filt!(out::AbstractArray{<:Any, N}, f::DF2TFilter{<:SecondOrderSections,<:Array}, x::AbstractArray{<:Any, N}) where N +function filt!(out::AbstractArray{<:Any,N}, f::DF2TFilter{<:SecondOrderSections,<:Array}, x::AbstractArray{<:Any,N}) where N size(x) != size(out) && throw(ArgumentError("out size must match x")) size(x)[2:end] != size(f.state)[3:end] && throw(ArgumentError("state size must match x")) for col in CartesianIndices(axes(x)[2:end]) @@ -221,7 +221,7 @@ DF2TFilter(coef::Biquad{:z,T}, coldims::Tuple{Vararg{Integer}}) where {T} = DF2T DF2TFilter(coef::Biquad{:z,T}, ::Type{V}=T, coldims::Tuple{Vararg{Integer}}=()) where {T,V} = DF2TFilter(coef, zeros(promote_type(T, V), 2, coldims...)) -function filt!(out::AbstractArray{<:Any, N}, f::DF2TFilter{<:Biquad,Array{T, N}} where T, x::AbstractArray{<:Any, N}) where N +function filt!(out::AbstractArray{<:Any,N}, f::DF2TFilter{<:Biquad,Array{T,N}} where T, x::AbstractArray{<:Any,N}) where N size(x) != size(out) && throw(ArgumentError("out size must match x")) si = f.state size(x)[2:end] != size(si)[2:end] && throw(ArgumentError("state size must match x")) From 307428a6a74821cb1de7caa864dbe6deed691cdc Mon Sep 17 00:00:00 2001 From: Martin Holters Date: Tue, 3 Dec 2024 11:43:29 +0100 Subject: [PATCH 6/6] Rephrase `DF2TFilter` docs --- src/Filters/filt.jl | 21 +++++++++++++-------- 1 file changed, 13 insertions(+), 8 deletions(-) diff --git a/src/Filters/filt.jl b/src/Filters/filt.jl index 58b38c742..4ec76e003 100644 --- a/src/Filters/filt.jl +++ b/src/Filters/filt.jl @@ -116,19 +116,24 @@ filt(f::FilterCoefficients{:z}, x) = filt(convert(SecondOrderSections, f), x) filt!(out, f::FilterCoefficients{:z}, x) = filt!(out, convert(SecondOrderSections, f), x) """ - DF2TFilter(coef::FilterCoefficients{:z}[, si]) - DF2TFilter(coef::FilterCoefficients{:z}[, sitype::Type][, coldims::Tuple]) + DF2TFilter(coef::FilterCoefficients{:z}) + DF2TFilter(coef::FilterCoefficients{:z}, coldims::Tuple) + DF2TFilter(coef::FilterCoefficients{:z}, sitype::Type, coldims::Tuple = ()) + DF2TFilter(coef::FilterCoefficients{:z}, si) Construct a stateful direct form II transposed filter with coefficients `coef`. -One can optionally specify as the second argument either -- `si`, an array representing the initial filter state, or -- `sitype`, the eltype of a zeroed `si` and/or `coldims`, the size of extra dimensions - of the input. E.g. to column-wise filter an input with dims `(L, N1, N2)`, set - `coldims` to `(N1, N2)`. +The initial filter state defaults to zeros (of a type derived from `coef`) +suitable for vector input. Another element type of the state can be specified +with `sitype`. -The initial filter state defaults to zeros if called with one argument. +To allow column-wise filtering of higher-dimensional input, the size of the +extra dimensions have to be given in `coldims`. To e.g. column-wise filter an +input with size `(L, N1, N2)`, set `coldims` to `(N1, N2)`. + +Alternatively, an array representing the initial filter state can be passed +as `si`. If `coef` is a `PolynomialRatio`, `Biquad`, or `SecondOrderSections`, filtering is implemented directly. If `coef` is a `ZeroPoleGain`