In signal processing, multidimensional discrete convolution refers to the mathematical operation between two functions f and g on an n-dimensional lattice that produces a third function, also of n-dimensions. Multidimensional discrete convolution is the discrete analog of the multidimensional convolution of functions on Euclidean space. It is also a special case of convolution on groups when the group is the group of n-tuples of integers.

Definition

Problem statement and basics

Similar to the one-dimensional case, an asterisk is used to represent the convolution operation. The number of dimensions in the given operation is reflected in the number of asterisks. For example, an M-dimensional convolution would be written with M asterisks. The following represents a M-dimensional convolution of discrete signals:

For discrete-valued signals, this convolution can be directly computed via the following:

The resulting output region of support of a discrete multidimensional convolution will be determined based on the size and regions of support of the two input signals.

Listed are several properties of the two-dimensional convolution operator. Note that these can also be extended for signals of -dimensions.

Commutative Property:

Associate Property:

Distributive Property:

These properties are seen in use in the figure below. Given some input that goes into a filter with impulse response and then another filter with impulse response , the output is given by . Assume that the output of the first filter is given by , this means that:

Further, that intermediate function is then convolved with the impulse response of the second filter, and thus the output can be represented by:

Using the associative property, this can be rewritten as follows:

meaning that the equivalent impulse response for a cascaded system is given by:

A similar analysis can be done on a set of parallel systems illustrated below.

In this case, it is clear that:

Using the distributive law, it is demonstrated that:

This means that in the case of a parallel system, the equivalent impulse response is provided by:

The equivalent impulse responses in both cascaded systems and parallel systems can be generalized to systems with -number of filters.

Motivation and applications

Convolution in one dimension was a powerful discovery that allowed the input and output of a linear shift-invariant (LSI) system (see LTI system theory) to be easily compared so long as the impulse response of the filter system was known. This notion carries over to multidimensional convolution as well, as simply knowing the impulse response of a multidimensional filter too allows for a direct comparison to be made between the input and output of a system. This is profound since several of the signals that are transferred in the digital world today are of multiple dimensions including images and videos. Similar to the one-dimensional convolution, the multidimensional convolution allows the computation of the output of an LSI system for a given input signal.

For example, consider an image that is sent over some wireless network subject to electro-optical noise. Possible noise sources include errors in channel transmission, the analog to digital converter, and the image sensor. Usually noise caused by the channel or sensor creates spatially-independent, high-frequency signal components that translates to arbitrary light and dark spots on the actual image. In order to rid the image data of the high-frequency spectral content, it can be multiplied by the frequency response of a low-pass filter, which based on the convolution theorem, is equivalent to convolving the signal in the time/spatial domain by the impulse response of the low-pass filter. Several impulse responses that do so are shown below.

In addition to filtering out spectral content, the multidimensional convolution can implement edge detection and smoothing. This once again is wholly dependent on the values of the impulse response that is used to convolve with the input image. Typical impulse responses for edge detection are illustrated below.

In addition to image processing, multidimensional convolution can be implemented to enable a variety of other applications. Since filters are widespread in digital communication systems, any system that must transmit multidimensional data is assisted by filtering techniques It is used in real-time video processing, neural network analysis, digital geophysical data analysis, and much more.

One typical distortion that occurs during image and video capture or transmission applications is blur that is caused by a low-pass filtering process. The introduced blur can be modeled using Gaussian low-pass filtering.

Row-column decomposition with separable signals

Separable signals

A signal is said to be separable if it can be written as the product of multiple one-dimensional signals. Mathematically, this is expressed as the following:

Some readily recognizable separable signals include the unit step function, and the dirac-delta impulse function.

(unit step function)

(dirac-delta impulse function)

Convolution is a linear operation. It then follows that the multidimensional convolution of separable signals can be expressed as the product of many one-dimensional convolutions. For example, consider the case where x and h are both separable functions.

By applying the properties of separability, this can then be rewritten as the following:

It is readily seen then that this reduces to the product of one-dimensional convolutions:

This conclusion can then be extended to the convolution of two separable M-dimensional signals as follows:

So, when the two signals are separable, the multidimensional convolution can be computed by computing one-dimensional convolutions.

Row-column decomposition The row-column method can be applied when one of the signals in the convolution is separable. The method exploits the properties of separability in order to achieve a method of calculating the convolution of two multidimensional signals that is more computationally efficient than direct computation of each sample (given that one of the signals are separable). The following shows the mathematical reasoning behind the row-column decomposition approach (typically is the separable signal):

The value of can now be re-used when evaluating other values with a shared value of :

Thus, the resulting convolution can be effectively calculated by first performing the convolution operation on all of the rows of , and then on all of its columns. This approach can be further optimized by taking into account how memory is accessed within a computer processor.

A processor will load in the signal data needed for the given operation. For modern processors, data will be loaded from memory into the processors cache, which has faster access times than memory. The cache itself is partitioned into lines. When a cache line is loaded from memory, multiple data operands are loaded at once. Consider the optimized case where a row of signal data can fit entirely within the processor's cache. This particular processor would be able to access the data row-wise efficiently, but not column-wise since different data operands in the same column would lie on different cache lines. In order to take advantage of the way in which memory is accessed, it is more efficient to transpose the data set and then access it row-wise rather than attempt to access it column-wise. The algorithm then becomes: Separate the separable two-dimensional signal into two one-dimensional signals and Perform row-wise convolution on the horizontal components of the signal using to obtain Transpose the vertical components of the signal resulting from Step 2. Perform row-wise convolution on the transposed vertical components of to get the desired output

Computational speedup from row-column decomposition Examine the case where an image of size is being passed through a separable filter of size . The image itself is not separable. If the result is calculated using the direct convolution approach without exploiting the separability of the filter, this will require approximately multiplications and additions. If the separability of the filter is taken into account, the filtering can be performed in two steps. The first step will have multiplications and additions and the second step will have , resulting in a total of or multiplications and additions. A comparison of the computational complexity between direct and separable convolution is given in the following image:

Circular convolution of discrete-valued multidimensional signals The premise behind the circular convolution approach on multidimensional signals is to develop a relation between the Convolution theorem and the Discrete Fourier transform (DFT) that can be used to calculate the convolution between two finite-extent, discrete-valued signals.

Convolution theorem in multiple dimensions

For one-dimensional signals, the Convolution Theorem states that the Fourier transform of the convolution between two signals is equal to the product of the Fourier Transforms of those two signals. Thus, convolution in the time domain is equal to multiplication in the frequency domain. Mathematically, this principle is expressed via the following:This principle is directly extendable to dealing with signals of multiple dimensions. This property is readily extended to the usage with the Discrete Fourier transform (DFT) as follows (note that linear convolution is replaced with circular convolution where is used to denote the circular convolution operation of size ):

When dealing with signals of multiple dimensions:The circular convolutions here will be of size .

Circular convolution approach

The motivation behind using the circular convolution approach is that it is based on the DFT. The premise behind circular convolution is to take the DFTs of the input signals, multiply them together, and then take the inverse DFT. Care must be taken such that a large enough DFT is used such that aliasing does not occur. The DFT is numerically computable when dealing with signals of finite-extent. One advantage this approach has is that since it requires taking the DFT and inverse DFT, it is possible to utilize efficient algorithms such as the Fast Fourier transform (FFT). Circular convolution can also be computed in the time/spatial domain and not only in the frequency domain.

Choosing DFT size to avoid aliasing

Consider the following case where two finite-extent signals x and h are taken. For both signals, there is a corresponding DFT as follows:

and

The region of support of is and and the region of support of is and .

The linear convolution of these two signals would be given as:Given the regions of support of and , the region of support of will then be given as the following:

Based on the regions of support of the two signals, a DFT of size must be used where and since the same size DFT must be used on both signals. In the event where a DFT size larger than the extent of a signal is needed, the signal is zero-padded until it reaches the required length. After multiplying the DFTs and taking the inverse DFT on the result, the resulting circular convolution is then given by:

for

The result will be that will be a spatially aliased version of the linear convolution result . This can be expressed as the following:

Then, in order to avoid aliasing between the spatially aliased replicas, and must be chosen to satisfy the following conditions:

If these conditions are satisfied, then the results of the circular convolution will equal that of the linear convolution (taking the main period of the circular convolution as the region of support). That is:

for

Summary of procedure using DFTs The Convolution theorem and circular convolution can thus be used in the following manner to achieve a result that is equal to performing the linear convolution: