Scale space implementation

By definition, the linear scale space representation of an N-dimensional signal is obtained by smoothing with the N-dimensional Gaussian kernel. Using the separability property of the Gaussian kernel, this operation can be decomposed into a set of separable smoothing steps with a one-dimensional kernel along each dimension.

$G(x,t)={\frac {1}{\sqrt {2{\pi }t}}}e^{-x^{2}/2t}$

When implementing this smoothing step in practice, there are several approaches that could be taken: The presumably simplest approach is to convolve the discrete signal with a sampled Gaussian kernel:

$G(n,t)={\frac {1}{\sqrt {2{\pi }t}}}e^{-n^{2}/2t}$

which in turn is truncated at the ends to give a filter with finite impulse response. A more refined approach is to convolve the original signal by the discrete Gaussian kernel .

$T(n,t)=e^{-t}I_{n}(t)$

where $I_n(t)$ denotes the modified Bessel functions of integer order. This filter can be truncated in the spatial domain or be implemented in terms of a closed-form expression for its discrete Fourier transform:

${\hat {T}}(\theta ,t)=e^{t(\cos \theta -1)}$

If computational efficiency is a primary criterion, implementation of the Gaussian kernel and Gaussian derivatives in terms of recursive filters is a common option . Please, note, however that at fine scales these separable approaches do not always give the best possible approximation to rotational symmetry, so non-separable implementations for two-dimensional images may be considered as an alternative.