Discrete two-dimensional Fourier transform c. Discrete Fourier Transform (DFT)

This is one of the Fourier transforms widely used in digital signal processing algorithms (its modifications are used in audio compression in MP3, image compression in JPEG, etc.), as well as in other areas related to the analysis of frequencies in discrete (for example, digitized analog ) signal. The discrete Fourier transform requires as input discrete function. Such functions are often created by sampling (sampling values ​​from continuous functions). Discrete Fourier transforms help solve partial differential equations and perform operations such as convolutions. Discrete Fourier transforms are also actively used in statistics, in the analysis of time series. Transformations can be one-dimensional, two-dimensional and even three-dimensional.

Direct conversion:

Reverse conversion:

Designations:

§ N- the number of signal values ​​measured over a period, as well as the number of decomposition components;

§ - measured signal values ​​(at discrete time points with numbers that are input data for direct conversion and weekends for the return;

§ - N complex amplitudes of sinusoidal signals composing the original signal; are output data for direct conversion and input data for reverse conversion; since amplitudes are complex, it is possible to calculate both amplitude and phase from them;

§ is the usual (real) amplitude of the kth sinusoidal signal;

§ arg( Xk) - phase of the kth sinusoidal signal (argument of a complex number);

§ k- frequency of the kth signal, equal to , where T- the period of time during which the input data was taken.

From the latter it is clear that the transformation decomposes the signal into sinusoidal components (which are called harmonics) with frequencies from N oscillations per period to one oscillation per period. Since the sampling frequency itself is equal to N samples per period, high-frequency components cannot be displayed correctly - a moiré effect occurs. This leads to the fact that the second half of the N complex amplitudes is, in fact, a mirror image of the first and does not carry additional information.

Consider some periodic signal x(t) with a period equal to T. Let us expand it into a Fourier series:

Let us sample the signal so that there are N samples per period. Let's represent the discrete signal in the form of samples: x n = x(tn), where , then these readings through the Fourier series will be written as follows:

Using the relation: , we get:

Where

So we got inverse discrete Fourier transform.

Let us now scalarly multiply the expression for x n on and we get:

Here we use: a) an expression for the sum of a finite number of terms (exponents) geometric progression, and b) the expression of the Kronecker symbol as the limit of the ratio of the Euler functions for complex numbers. It follows that:

This formula describes direct discrete Fourier transform.

In the literature, it is customary to write the multiplier in the inverse transformation, and therefore the transformation formulas are usually written in the following form:

The discrete Fourier transform is linear transformation, which transforms a vector of time samples into a vector of spectral samples of the same length. So the transformation can be implemented as a multiplication square matrix to vector:

Fourier transform

When using Fourier transforms, the image is represented as a sum of complex exponential functions variables of amplitude, frequency and phase. The Fourier transform plays very important role in many areas of image processing, including enhancement, analysis, restoration and compression.

1. Basic definitions of the Fourier transform
2. Discrete Fourier transform, including fast Fourier transform
3. Applications of the Fourier transform (some examples of practical applications of the Fourier transform)

Basic definitions of the Fourier transform

If ƒ(m,n) is a function of two discrete spatial variables m and n, then two-dimensional transformation Fourier functions ƒ(m,n) can be represented by the following expression

The variables are angular frequencies. Thus, it represents a function ƒ(m,n) in the frequency domain. is a complex-valued function with corresponding frequencies. The frequencies are within the range , . Note that F(0,0) is represented as the sum of all variables ƒ(m,n). For this reason F(0,0) is often called the constant component of the Fourier transform.

The inverse two-dimensional Fourier transform is represented by the expression

Those. this expression represents ƒ(m,n) as a sum infinite number complex exponential functions (sine waves) with different frequencies. Amplitude and phase determine the contribution of frequencies to the representation.

Visualization of the Fourier transform

In illustrating the Fourier transform, let us assume that the function ƒ(m,n) equals 1 and is represented as a rectangle. To simplify the diagram, the function ƒ(m,n) will be represented as a continuous function of two discrete variables m And n.

Rectangular function

The figure below, using the mesh function, visualizes the amplitude values ​​obtained from the Fourier transform rectangular function shown in the previous figure. Amplitude visualization is also called Fourier transform visualization.

Rectangular function image amplitude

The peak of the function is in the center and displays the value F(0,0), which is the sum of all values ƒ(m,n). All other components represent the distribution of energy across vertical and horizontal frequencies.

Another way to visualize the Fourier transform is to display the values ​​as an image.

Logarithmic representation of the Fourier transform of a rectangular function

Let's look at examples of Fourier transforms of functions of various simple forms.

Examples of Fourier transforms of functions of various simple forms

Discrete cosine transform

Discrete cosine transforms represent an image as a sum of sinusoids with different amplitudes and frequencies. dct2 function in Image application Processing Toolbox implements two-dimensional discrete cosine transforms of images. One of the features of the discrete Fourier transform is that some local areas images can be characterized a small amount discrete Fourier transform coefficients. This property is very often used in the development of image compression methods. For example, the discrete cosine transform is the basis of an international standard used in the JPEG lossy image compression algorithm. The name of the format “JPEG” consists of the first letters of the name working group, which participated in the development of this standard (Joint Photographic Experts Group).

Two-dimensional discrete cosine matrix transform A with dimensions is implemented according to the following expression

Values Bpq are called the coefficients of the discrete cosine transformation of the matrix A.

(It should be noted that matrix indices in MATLAB always start at 1, not 0. Therefore, matrix elements that are represented in MATLAB as A(1,1) and B(1,1) will correspond to elements A 00 And B00 from the above formula.)

The inverse discrete cosine transform is implemented according to the expressions

The inverse discrete cosine transform expression can be interpreted as a matrix representation A with dimensions as the sum of the following functions

These functions are called the fundamental (basic) functions of the discrete cosine transform. Discrete Cosine Transform Coefficients Bpq can be considered as weights for each basic function. For example, for a matrix with element size there are 64 basic functions, which is shown in the image.

64 basic functions that are obtained for a matrix with element sizes

Horizontal frequencies increase from left to right, and vertical frequencies increase from top to bottom.

Discrete Cosine Transform Matrix

Application Image Processing Toolbox offers two different ways implementation of discrete cosine transformations. The first method is implemented in the dct2 function. The dct2 function uses fast conversion Fourier for speeding up calculations. The second method uses the discrete cosine transform matrix, which is returned by the dctmtx function. The transformation matrix T is formed according to the following expression

For matrix A with dimensions is a matrix with dimensions , where each column contains a one-dimensional discrete cosine transform A. Two-dimensional discrete cosine transform A calculated as B=T*A*T’. Inverse two-dimensional discrete cosine transform B calculated as T'*B*T.

Discrete Cosine Transforms and Image Compression

In the JPEG image compression algorithm, the original image is divided into blocks of dimensions or elements. Next, a two-dimensional discrete cosine transform is calculated for each block. The coefficients of discrete cosine transforms are quantized, encoded and transmitted. The JPEG receiver decodes the discrete cosine transform coefficients, calculates the inverse 2D discrete cosine transform in each block, and then stitches them together into a single image.

Let's consider an example of calculating two-dimensional discrete cosine transforms in blocks with the sizes of the elements of the original image. Further, when reconstructing the image, we will take into account only 10 coefficients from each block, the rest will be set to zero. When carrying out the described calculations, the transformation matrix will also be used.

I = imread("cameraman.tif"); I = im2double(I); T = dctmtx(8); B = blkproc(I,,"P1*x*P2",T,T"); mask = ; B2 = blkproc(B,,"P1.*x",mask); I2 = blkproc(B2,,"P1 *x*P2",T",T); imshow(I); figure, imshow(I2)

The figure shows two images - the original and the reconstructed one. Only 15% of discrete cosine transform coefficients were used in image reconstruction. However, it should be noted that the quality of the reconstructed image is quite acceptable. To view other properties of the discrete cosine transform, see the dctdemo function.

Radon transformations

The radon function in the Image Processing Toolbox calculates a matrix of image projections along given directions. The projection of a two-dimensional function f(x,y) is equal to the integral along the indicated line. The Radon function is the calculation of image projections on an axis, which are specified by angles in degrees relative to the horizontal counterclockwise. The figure shows the projection of a certain figure at a specified angle

Parallel beam projection with rotation angle theta

The figure below shows the horizontal and vertical projections for a simple two-dimensional function.

Horizontal and vertical projections of some simple function

Projections can be calculated along arbitrary angle theta. The radon function built into the Image Processing Toolbox calculates image projections along certain directions. The projection of a two-dimensional function f(x,y) onto the x’ axis is a linear integral

Thus, the x’ y’ axes are specified by rotating by an angle counterclockwise.

The image below illustrates the geometry of the Radon transform.

Radon Transform Geometry

Visualization of Radon Transforms

When performing Radon transformations, it is necessary to specify the original image and the vector of angles theta.

Radon(I,theta);

R is a matrix in which each column is the Radon transform for one of the angles contained in the vector theta. The vector xp contains the corresponding coordinates along the x axis. The central pixel I is determined according to the expression floor((size(I)+1)/2).

Let's look at how projections are calculated in Radon transforms. Let's consider projections at an angle of 0° and 45°.

I = zeros(100,100); I(25:75, 25:75) = 1; imshow(I)

Radon(I,); figure; plot(xp,R(:,1)); title("R_(0^o) (x\prime)")

Radon transformations at 0°

Figure; plot(xp,R(:,2)); title("R_(45^o) (x\prime)")

Radon transformations at 45°

The Radon transform at a large number of angles is often displayed as an image. In this example, the Radon transform is considered for an image in the form of a square with an angle range from 0° to 180° with a resolution of 1°.

Theta = 0:180;

= radon(I,theta); imagesc(theta,xp,R); title("R_(\theta) (X\prime)"); xlabel("\theta (degrees)"); ylabel("X\prime"); set(gca,"XTick",0:20:180); colormap(hot); colorbar

Radon transformations using 180 projections

Using Radon transforms when detecting lines

Radon transformations are similar to other well-known operations, which are known as Hoch transformations. The radon function can be used to detect straight lines. Let's look at the main stages of this process. R Largest peak in matrix corresponds to =1° and x´= -80. A line is drawn from the center of the original image at an angle at a distance x’. A straight line is drawn perpendicular to this line, which corresponds to the straight line on original image R. In addition, there are other lines in the image that are presented in the matrix

corresponding peaks.

Radon Transform Geometry for Straight Line Detection

Let us denote by

. (3.21)

a two-dimensional field (two-dimensional signal) describing a discrete image of the size of rows and columns. Outside the specified limits, this signal is not defined. Let us perform a periodic continuation of this finite signal by introducing a two-dimensional periodic signal

Rice. 3.4. Real (a) and periodically continued (b) images

Any periodic signal can be represented as a Fourier series, but, unlike one-dimensional signals, two-dimensional signals are described by a two-dimensional Fourier series, which has the form:

(3.23)

The basis functions of this two-dimensional representation are two-dimensional complex exponentials (sometimes called complex sinusoids)

having, like the signal, a rectangular periodicity with the same period. Here (,) is the two-dimensional number of the basis function, and the quantities have the meaning of spatial frequencies. Sometimes integer quantities and are called spatial frequencies. The Fourier coefficients of series (3.22) form a two-dimensional frequency spectrum

(3.24)

signal and are determined by the formula for the direct Fourier transform: Expression (3.22), which restores the signal from its spectrum, is the inverse Fourier transform. The validity of transformations (3.22) and (3.24), called two-dimensional DFT, can be verified by substituting (3.24) into (3.22) and bringing right side

Note that for an accurate representation of a discrete signal with a two-dimensional period of elements according to the FFT formulas, a finite number of basis functions (3.23) is sufficient - the series (3.22) is finite. This is understandable, since the represented signal itself contains in one period final number points, i.e. has a finite number of degrees of freedom. It is clear that the number of degrees of freedom in the spectrum cannot differ from the number of degrees of freedom in the signal itself.

Let us dwell on the most essential properties of the two-dimensional discrete Fourier spectrum. Let us calculate the spectral coefficients (3.24) at frequency points :

Since for any integer values ​​and the last factor in the resulting expression equal to one, then we have the equality:

,

denoting the rectangular periodicity of the two-dimensional DFT. Consequently, the picture of a two-dimensional DFT is similar to the picture of a two-dimensional periodically continuous signal, qualitatively shown in Fig. 3.4.b (if the spatial coordinates on it are replaced by frequency ones). However, it must be borne in mind that the spectral coefficients, as follows from (3.24), are complex numbers, including for a real signal. But then the question arises. The total number of spectral components is found to be . A complex number is equivalent to a pair of real numbers - the real and imaginary parts in algebraic notation or the modulus and phase in exponential notation. Therefore, the full spectrum is described real numbers, which is twice the dimension of the signal itself. At first glance, this contains a contradiction. It finds its clarification with further study of the properties of two-dimensional DFT.

Let us transform relation (3.25) as follows. First, let's substitute frequencies instead of frequencies. Secondly, we will perform a complex conjugation of both parts, which will not violate the equality. As a result, it is easy to obtain the expression:

,

which establishes an unambiguous connection between the spectral coefficients at two different points of the spectral rectangle. The resulting relationship removes the contradiction, since the number of independent spectral coefficients is reduced by half due to this spectral symmetry. According to the established property, the spectral coefficients belonging to the upper left and lower right corners of the rectangle are connected by a spectral conjugate dependence. Similarly, the Fourier coefficients from the upper right and lower left sections of the spectral rectangle are also related to each other.

At the end of this paragraph, we point out that when practical application two-dimensional DFT - both direct and inverse, there is no need to operate with periodic signals and spectra, as would seem to be assumed by transformations (3.22) and (3.24). Relations (3.22) and (3.24) themselves eliminate this need. In fact, the direct Fourier transform (3.24) contains on the right side the values ​​of the periodically continued signal only within one “main” rectangle. But within these limits, the original and periodically continued signals completely coincide, which makes it possible to use the original signal in formula (3.24). Similar explanations can be made regarding inverse conversion(3.22), from which it follows that practically in the process of calculations one should operate with the “main” part of the spectrum related to spectral region.

From the explanations made, which have only purely computational significance, one should not draw a conclusion about the artificiality and uselessness of the considered mathematical models periodic fields. When processing images, numerous problems arise, the correct interpretation and solution of which is only possible on the basis of these mathematical interpretations. One of these most important tasks is digital two-dimensional filtering in the spectral domain, the implementation of which is associated with the implementation of the so-called cyclic convolution.

Fourier transforms

It is convenient to analyze many signals by decomposing them into sinusoids (harmonics). There are several reasons for this. For example, the human ear works in a similar way. It decomposes sound into individual vibrations of different frequencies. Additionally, it can be shown that sinusoids are " own functions» linear systems (since they pass through linear systems, without changing the shape, but can only change the phase and amplitude). Another reason is that Kotelnikov's theorem is formulated in terms of the signal spectrum.

Fourier transform ) is the decomposition of functions into sinusoids (hereinafter we also call cosine functions sinusoids, since they differ from “real” sinusoids only in phase). There are several types of Fourier transform.

1. Non-periodic continuous signal can be expanded into a Fourier integral.

2. A periodic continuous signal can be expanded into an infinite Fourier series.

3. A non-periodic discrete signal can be expanded into a Fourier integral.

4. A periodic discrete signal can be expanded into a finite Fourier series.

A computer can only work with a limited amount of data, therefore, it can really only calculate the last type of Fourier transform. Let's take a closer look at it.

Real signal DFT

Let a discrete signal x have a period of N points. In this case, it can be represented as a finite series (i.e., a linear combination) of discrete sinusoids:

Equivalent notation (we decompose each cosine into sine and cosine, but now without a phase):

Rice. 6. Fourier series basis functions for an 8-point discrete signal. On the left are cosines, on the right are sines. Frequencies increase from top to bottom.

Basic sinusoids have multiple frequencies. The first term of the series (k =0) is a constant called constant component(DC offset) signal. The very first sinusoid (k = 1) has such a frequency that its period coincides with the period of the original signal itself. The highest frequency component (k =N /2) has such a frequency that its period is equal to two counts. CoefficientsA k and

B k is called the signal spectrum (spectrum). They show the amplitudes of the si-

nusoids that make up the signal. The frequency step between two adjacent sinusoids from the Fourier expansion is called frequency resolution spectrum

In Fig. Figure 6 shows the sinusoids used to decompose a discrete signal from 8 points. Each of the sinusoids consists of 8 points, that is, it is a regular discrete signal. Continuous sinusoids are shown in the figure for clarity.

convert the original signal by calculating the sum of the Fourier series at each point. Decomposing a signal into sinusoids (i.e. obtaining coefficients) is called direct Fourier transform. The reverse process - signal synthesis using sinusoids - is called inverse Fourier transform(inverse Fourier transform).

The algorithm for the inverse Fourier transform is obvious (it is contained in the formula of the Fourier series; to carry out the synthesis, you just need to substitute the coefficients into it). Let's consider the direct Fourier transform algorithm, i.e. finding the coefficients A k and B k .

tonal basis in the space of periodic discrete signals with period N. This means that to decompose any element of space (signal) into it, you need to calculate dot products this element with all the functions of the system, and the resulting coefficients are normalized. Then the basis expansion formula with coefficients A k and B k will be valid for the original signal.

So, the coefficients A k and B k are calculated as scalar products (in non-

in the discontinuous case – integrals of the product of functions, in the discrete case

– sums from the product of discrete signals):

Modern communication technology cannot be imagined without spectral analysis. Representation of signals in the frequency domain is necessary both for the analysis of their characteristics and for the analysis of blocks and units of transceivers of radio communication systems. To convert signals into the frequency domain, a direct Fourier transform is used. The generalized formula for the direct Fourier transform is written as follows:

As can be seen from this formula for frequency analysis, the calculation is made correlation dependence between a signal represented in the time domain and a complex exponential at a given frequency. In this case, according to Euler’s formula, the complex exponential is decomposed into a real and imaginary part:

A signal represented in the frequency domain can be converted back into a time domain using an inverse Fourier transform. The generalized formula for the inverse Fourier transform is written as follows:

The direct Fourier transform formula uses time integration from minus infinity to infinity. Naturally this is a mathematical abstraction. IN real conditions we can integrate from at this moment time, which we can denote as 0, before time T. The formula for the direct Fourier transform will be transformed to the following form:

As a result the properties of the Fourier transform change significantly. Signal spectrum instead continuous function becomes a discrete series of values. Now the minimum frequency and at the same time the step of frequency values ​​of the signal spectrum becomes:

Only functions sin and cos with frequencies k/T will be mutually orthogonal, and this is an indispensable condition for the Fourier transform. The set of the first functions of the Fourier series expansion is shown in Figure 1. In this case, the duration of the functions coincides with the duration of the analysis T.

Now the signal spectrum will look as shown in Figure 2.

IN in this case the formula for calculating the direct Fourier transform (4) is transformed to the following form:

The formula for the inverse Fourier transform for the case of determining the spectrum over a limited period of time will look like this:

In a similar way, you can determine the formula for the direct Fourier transform for digital signal samples. Considering that instead of a continuous signal its digital samples are used, in expression (6) the integral is replaced by a sum. In this case, the duration of the analyzed signal is determined by the number of digital samples N. The Fourier transform for digital signal samples is called the discrete Fourier transform and is written as follows:

Now let's look at how the properties of the discrete Fourier transform (DFT) have changed compared to the direct Fourier transform over a limited time interval. When we looked at sampling analog signal, we found that the spectrum of the input signal must be limited in frequency. This requirement limits the number of discrete components of the signal spectrum. Initially it may seem that we can limit the spectrum of the signal to the frequency f d/2, which corresponds to the number of frequency components K=N/2. However, it is not. Although the signal spectrum for real signal samples for positive frequencies and negative frequencies is symmetrical about 0, negative frequencies may be required for some spectrum algorithms, such as . The difference is even greater when performing a discrete Fourier transform on complex samples of the input signal. As a result for full description spectrum digital signal required N frequency samples ( k = 0, ..., N/2).