ImageFiltering.jl
ImageFiltering supports linear and nonlinear filtering operations on arrays, with an emphasis on the kinds of operations used in image processing. The core function is imfilter, and common kernels (filters) are organized in the Kernel and KernelFactors modules.
Demonstration
Let's start with a simple example of linear filtering:
julia> using ImageFiltering, TestImages
julia> img = testimage("mandrill");
julia> imgg = imfilter(img, Kernel.gaussian(3));
julia> imgl = imfilter(img, Kernel.Laplacian());When displayed, these three images look like this:
The most commonly used function for filtering is imfilter.
Linear filtering: noteworthy features
Correlation, not convolution
ImageFiltering uses the following formula to calculate the filtered image F from an input image A and kernel K:
\[F[I] = \sum_J A[I+J] K[J]\]
Consequently, the resulting image is the correlation, not convolution, of the input and the kernel. If you want the convolution, first call reflect on the kernel.
Kernel indices
ImageFiltering exploits a feature introduced into Julia 0.5, the ability to define arrays whose indices span an arbitrary range:
julia> Kernel.gaussian(1)
OffsetArray{Float64,2,Array{Float64,2}} with indices -2:2×-2:2:
0.00296902 0.0133062 0.0219382 0.0133062 0.00296902
0.0133062 0.0596343 0.0983203 0.0596343 0.0133062
0.0219382 0.0983203 0.162103 0.0983203 0.0219382
0.0133062 0.0596343 0.0983203 0.0596343 0.0133062
0.00296902 0.0133062 0.0219382 0.0133062 0.00296902The indices of this array span the range -2:2 along each axis, and the center of the gaussian is at position [0,0]. As a consequence, this filter "blurs" but does not "shift" the image; were the center instead at, say, [3,3], the filtered image would be shifted by 3 pixels downward and to the right compared to the original.
The centered function is a handy utility for converting an ordinary array to one that has coordinates [0,0,...] at its center position:
julia> centered([1 0 1; 0 1 0; 1 0 1])
OffsetArray{Int64,2,Array{Int64,2}} with indices -1:1×-1:1:
1 0 1
0 1 0
1 0 1See OffsetArrays for more information.
Factored kernels
A key feature of Gaussian kernels–-along with many other commonly-used kernels–-is that they are separable, meaning that K[j_1,j_2,...] can be written as $K_1[j_1] K_2[j_2] \cdots$. As a consequence, the correlation
\[F[i_1,i_2] = \sum_{j_1,j_2} A[i_1+j_1,i_2+j_2] K[j_1,j_2]\]
can be written
\[F[i_1,i_2] = \sum_{j_2} \left(\sum_{j_1} A[i_1+j_1,i_2+j_2] K_1[j_1]\right) K_2[j_2]\]
If the kernel is of size m×n, then the upper version line requires mn operations for each point of filtered, whereas the lower version requires m+n operations. Especially when m and n are larger, this can result in a substantial savings.
To enable efficient computation for separable kernels, imfilter accepts a tuple of kernels, filtering the image by each sequentially. You can either supply m×1 and 1×n filters directly, or (somewhat more efficiently) call kernelfactors on a tuple-of-vectors:
julia> kern1 = centered([1/3, 1/3, 1/3])
OffsetArray{Float64,1,Array{Float64,1}} with indices -1:1:
0.333333
0.333333
0.333333
julia> kernf = kernelfactors((kern1, kern1))
(ImageFiltering.KernelFactors.ReshapedOneD{Float64,2,0,OffsetArray{Float64,1,Array{Float64,1}}}([0.333333,0.333333,0.333333]),ImageFiltering.KernelFactors.ReshapedOneD{Float64,2,1,OffsetArray{Float64,1,Array{Float64,1}}}([0.333333,0.333333,0.333333]))
julia> kernp = broadcast(*, kernf...)
OffsetArray{Float64,2,Array{Float64,2}} with indices -1:1×-1:1:
0.111111 0.111111 0.111111
0.111111 0.111111 0.111111
0.111111 0.111111 0.111111
julia> imfilter(img, kernf) ≈ imfilter(img, kernp)
trueIf the kernel is a two dimensional array, imfilter will attempt to factor it; if successful, it will use the separable algorithm. You can prevent this automatic factorization by passing the kernel as a tuple, e.g., as (kernp,).
Popular kernels in Kernel and KernelFactors modules
The two modules Kernel and KernelFactors implement popular correlation kernels in "dense" and "factored" forms, respectively. Type ?Kernel or ?KernelFactors at the REPL to see which kernels are supported.
A common task in image processing and computer vision is computing image gradients (derivatives), for which there is the dedicated function imgradients.
Automatic choice of FIR or FFT
For linear filtering with a finite-impulse response filtering, one can either choose a direct algorithm or one based on the fast Fourier transform (FFT). By default, this choice is made based on kernel size. You can manually specify the algorithm using Algorithm.FFT() or Algorithm.FIR().
Multithreading
If you launch Julia with JULIA_NUM_THREADS=n (where n > 1), then FIR filtering will by default use multiple threads. You can control the algorithm by specifying a resource as defined by ComputationalResources. For example, imfilter(CPU1(Algorithm.FIR()), img, ...) would force the computation to be single-threaded.
Arbitrary operations over sliding windows
This package also exports mapwindow, which allows you to pass an arbitrary function to operate on the values within a sliding window.
mapwindow has optimized implementations for some functions (currently, extrema).