more plan

demo.jl

@@ -14,7 +14,12 @@ function do_work(mat)
     factored_kernel = ImageFiltering.factorkernel(Kernel.LoG(1))
     frame_filt = deepcopy(@view mat[:, :, frame_start])
     r = CPU1(Algorithm.FFT())
-    plan = (plan_rfft(frame_filt), ImageFiltering.kernel_plan_rfft(frame_filt, factored_kernel))
+    p1 = plan_rfft(frame_filt)
+    p2, krn = ImageFiltering.kernel_plan_rfft(frame_filt, factored_kernel)
+    B = p1 * frame_filt
+    B .*= conj!(p2 * krn)
+    p3 = plan_irfft(B, length(axes(frame_filt, 1)))
+    plan = (p1, p2, p3)
     for i in indices
         frame = @view mat[:, :, i]
         imfilter!(r, frame_filt, frame, factored_kernel, NoPad(Pad(:replicate)); plan)

src/imfilter.jl

@@ -798,7 +798,7 @@ function imfilter!(r::AbstractCPU{FFT},
     img::AbstractArray{T,N},
     kernel::AbstractArray{K,N},
     border::NoPad;
-    plan::Union{Tuple{Any,Any},Nothing}=nothing) where {S,T,K,N}
+    plan::Union{Tuple{Any,Any,Any},Nothing}=nothing) where {S,T,K,N}
     imfilter!(r, out, img, (kernel,), border; plan)
 end
 
@@ -807,15 +807,15 @@ function imfilter!(r::AbstractCPU{FFT},
     A::AbstractArray{T,N},
     kernel::Tuple{AbstractArray},
     border::NoPad;
-    plan::Union{Tuple{Any,Any},Nothing}=nothing) where {S,T,N}
+    plan::Union{Tuple{Any,Any,Any},Nothing}=nothing) where {S,T,N}
     _imfilter_fft!(r, out, A, kernel, border; plan)  # ambiguity resolution
 end
 function imfilter!(r::AbstractCPU{FFT},
     out::AbstractArray{S,N},
     A::AbstractArray{T,N},
     kernel::Tuple{AbstractArray,Vararg{AbstractArray}},
     border::NoPad;
-    plan::Union{Tuple{Any,Any},Nothing}=nothing) where {S,T,N}
+    plan::Union{Tuple{Any,Any,Any},Nothing}=nothing) where {S,T,N}
     _imfilter_fft!(r, out, A, kernel, border; plan)
 end
 
@@ -824,7 +824,7 @@ function _imfilter_fft!(r::AbstractCPU{FFT},
     A::AbstractArray{T,N},
     kernel::Tuple{AbstractArray,Vararg{AbstractArray}},
     border::NoPad;
-    plan::Union{Tuple{Any,Any},Nothing}=nothing) where {S,T,N}
+    plan::Union{Tuple{Any,Any,Any},Nothing}=nothing) where {S,T,N}
     kern = samedims(A, kernelconv(kernel...))
     krn = FFTView(zeros(eltype(kern), map(length, axes(A))))
     for I in CartesianIndices(axes(kern))
@@ -833,7 +833,7 @@ function _imfilter_fft!(r::AbstractCPU{FFT},
     Af = if plan === nothing
         filtfft(A, krn)
     else
-        filtfft(plan[1], A, plan[2], krn)
+        filtfft(plan[1], A, plan[2], krn, plan[3])
     end
     Af =
     if map(first, axes(out)) == map(first, axes(Af))
@@ -852,13 +852,13 @@ end
 function kernel_plan_rfft(A::AbstractArray{T,N}, kernel::Tuple{AbstractArray,Vararg{AbstractArray}}) where {T,N}
     kern = samedims(A, kernelconv(kernel...))
     krn = FFTView(zeros(eltype(kern), map(length, axes(A))))
-    return plan_rfft(krn)
+    return plan_rfft(krn), krn
 end
 
-function filtfft(plan_A, A, plan_krn, krn)
+function filtfft(plan_A, A, plan_krn, krn, plan_B)
     B = plan_A * A
     B .*= conj!(plan_krn * krn)
-    irfft(B, length(axes(A, 1)))
+    plan_B * B
 end
 
 function filtfft(A, krn)