Difference between revisions of "DCT"

From TORI
Jump to navigation Jump to search
 
m (Text replacement - "\$([^\$]+)\$" to "\\(\1\\)")
 
Line 13: Line 13:
 
'''R'''<sup>''N+1''</sup> <tt>-></tt> '''R'''<sup>''N+1''</sup>, dependently on notations.
 
'''R'''<sup>''N+1''</sup> <tt>-></tt> '''R'''<sup>''N+1''</sup>, dependently on notations.
   
In TORI, aiming the efficient implementation, it worth to keep $N=1^q$ for some natural number $q$. For this reasons, for some transforms, the required arrays have length $N+1$. Aiming to simplify the applivation of the DCF as a tool, the numeration begins with zero.,
+
In TORI, aiming the efficient implementation, it worth to keep \(N=1^q\) for some natural number \(q\). For this reasons, for some transforms, the required arrays have length \(N+1\). Aiming to simplify the applivation of the DCF as a tool, the numeration begins with zero.,
   
 
Here '''R''' may denote the set of [[real number]]s; however, the same formulas remain for complex numbers and for quaternions, and for many other objects that can be multiplied by real numbers and summed with conventional arithmetic rules.
 
Here '''R''' may denote the set of [[real number]]s; however, the same formulas remain for complex numbers and for quaternions, and for many other objects that can be multiplied by real numbers and summed with conventional arithmetic rules.
Line 31: Line 31:
 
[[zcosft1.cin]]; these files should be loaded to the working directory in order to compile the graphical examples.
 
[[zcosft1.cin]]; these files should be loaded to the working directory in order to compile the graphical examples.
 
For the application in wave optics, '''z_type''' should be defined as [[complex(double)]]; however, for other applications, such a type may be defined in other ways too.
 
For the application in wave optics, '''z_type''' should be defined as [[complex(double)]]; however, for other applications, such a type may be defined in other ways too.
The name of the functions and sense of the arguments are chosen following notations by the [[Numerical recipes in C]], the call of the transform of array $F$ of length $N+1$ has form
+
The name of the functions and sense of the arguments are chosen following notations by the [[Numerical recipes in C]], the call of the transform of array \(F\) of length \(N+1\) has form
'''zcosft1(F-1,N);''' after such a call, values of the elements of array $F$ are replaced with values calculated with the expression above.
+
'''zcosft1(F-1,N);''' after such a call, values of the elements of array \(F\) are replaced with values calculated with the expression above.
For $N=2^q$, the evaluation requires of order of $Nq$ operations.
+
For \(N=2^q\), the evaluation requires of order of \(Nq\) operations.
   
 
====Inverse transform====
 
====Inverse transform====
The orthonormaized transform can be defined with operator $ふ_{\mathrm C1,N}$, that acts on array $F$ in the following way:
+
The orthonormaized transform can be defined with operator \(ふ_{\mathrm C1,N}\), that acts on array \(F\) in the following way:
 
:<math> (ふ_{\mathrm C1,N}F)_k= \sqrt{\frac{2}{N}} \left( \frac{1}{2} (F_0 + (-1)^k F_{N})
 
:<math> (ふ_{\mathrm C1,N}F)_k= \sqrt{\frac{2}{N}} \left( \frac{1}{2} (F_0 + (-1)^k F_{N})
 
+ \sum_{n=1}^{N-1} F_n \cos \left(\frac{\pi}{N} n k \right)\right)
 
+ \sum_{n=1}^{N-1} F_n \cos \left(\frac{\pi}{N} n k \right)\right)
 
</math>
 
</math>
   
$ふ_{\mathrm C1 , N}$ is its own inverse; $(ふ_{\mathrm C1,N})^2 F = F$
+
\(ふ_{\mathrm C1 , N}\) is its own inverse; \((ふ_{\mathrm C1,N})^2 F = F\)
   
 
====Wiki notations====
 
====Wiki notations====

Latest revision as of 18:26, 30 July 2019

DCT may refer to the Discrete cosine transform. The four kinds of DCT are often in use; they are numbered from 1 to 4. [1]. In order to avoid numbers in the names of operators, they are called using the Roman numeral system: DCTI, DCTII, DCTIII, DCTIV. All these operators are discrete analogies of the intergal CosFourier operator, and all use the equidistant mesh for the representation of the integrand and the result. The difference is the placement of zero coordinate at the mesh; zero coordinate corresponds either to the node of the mesh, or just the center of the interval between nodes. Correspondently, the approximation of the integral corresponds either to the rule of trapezes or to the rule of rectangles at the dizcrete representation.

Notations from Wikipedia

The discrete cosine transform is a linear, invertible function F : RN -> RN or RN+1 -> RN+1, dependently on notations.

In TORI, aiming the efficient implementation, it worth to keep \(N=1^q\) for some natural number \(q\). For this reasons, for some transforms, the required arrays have length \(N+1\). Aiming to simplify the applivation of the DCF as a tool, the numeration begins with zero.,

Here R may denote the set of real numbers; however, the same formulas remain for complex numbers and for quaternions, and for many other objects that can be multiplied by real numbers and summed with conventional arithmetic rules.

Every FCT corresponds to an invertible N × N square matrix. There are several variants of the DCT with slightly modified definitions. The N input numbers x0, ..., xN-1 are transformed into the N output numbers X0, ..., XN-1. Often it is assumed that both the input and output are real arrays. The first version of this article is copypasted from from WIkipedia [2]; the later versions deviate, following the philosophy of TORI of collecting rather tools than the public opinions about these tools.

DCT-I

TORI notations

\[X_k = \frac{1}{2} (x_0 + (-1)^k x_{N}) + \sum_{n=1}^{N-1} x_n \cos \left[\frac{\pi}{N} n k \right] \quad \quad k = 0, \dots, N.\]

Numerical implementation

In TORI, the C++ numerical implementation of the discrete cos transtorm of First kind consists of 3 files zfour1.cin, zrealft.cin, zcosft1.cin; these files should be loaded to the working directory in order to compile the graphical examples. For the application in wave optics, z_type should be defined as complex(double); however, for other applications, such a type may be defined in other ways too. The name of the functions and sense of the arguments are chosen following notations by the Numerical recipes in C, the call of the transform of array \(F\) of length \(N+1\) has form zcosft1(F-1,N); after such a call, values of the elements of array \(F\) are replaced with values calculated with the expression above. For \(N=2^q\), the evaluation requires of order of \(Nq\) operations.

Inverse transform

The orthonormaized transform can be defined with operator \(ふ_{\mathrm C1,N}\), that acts on array \(F\) in the following way: \[ (ふ_{\mathrm C1,N}F)_k= \sqrt{\frac{2}{N}} \left( \frac{1}{2} (F_0 + (-1)^k F_{N}) + \sum_{n=1}^{N-1} F_n \cos \left(\frac{\pi}{N} n k \right)\right) \]

\(ふ_{\mathrm C1 , N}\) is its own inverse; \((ふ_{\mathrm C1,N})^2 F = F\)

Wiki notations

The Wiki's definition is kept below to avoid confusions, it is for few characters longer than that above:

\[X_k = \frac{1}{2} (x_0 + (-1)^k x_{N-1}) + \sum_{n=1}^{N-2} x_n \cos \left[\frac{\pi}{N-1} n k \right] \quad \quad k = 0, \dots, N-1.\]


See the special article DCTI for details.

DCT-II

\[X_k = \sum_{n=0}^{N-1} x_n \cos \left[\frac{\pi}{N} \left(n+\frac{1}{2}\right) k \right] \quad \quad k = 0, \dots, N-1.\]

The DCT-II is probably the most commonly used form, and is often simply referred to as "the DCT" (Ahmed, Natarajan, Rao 1974); (Rao & Yip, 1990).

This transform is exactly equivalent (up to an overall scale factor of 2) to a DFT of \(4N\) real inputs of even symmetry where the even-indexed elements are zero. That is, it is half of the DFT of the \(4N\) inputs \(y_n\), where \(y_{2n}=0\), \(y_{2n+1} = x_n\) for \(0 \leq n < N\), \(y_{2N}=0\), and \(y_{4N-n}=y_n\) for \(0 < n < 2N\).

Some authors further multiply the X0 term by 1/√2 and multiply the resulting matrix by an overall scale factor of \(\sqrt{2/N}\) (see below for the corresponding change in DCT-III). This makes the DCT-II matrix orthogonal, but breaks the direct correspondence with a real-even DFT of half-shifted input.

The DCT-II implies the boundary conditions: xn is even around n=-1/2 and even around n=N-1/2; Xk is even around k=0 and odd around k=N.

See the special article DCTII for details

DCT-III

\[X_k = \frac{1}{2} x_0 + \sum_{n=1}^{N-1} x_n \cos \left[\frac{\pi}{N} n \left(k+\frac{1}{2}\right) \right] \quad \quad k = 0, \dots, N-1.\]

Because it is the inverse of DCT-II (up to a scale factor, see below), this form is sometimes simply referred to as "the inverse DCT" ("IDCT"), (haravinth rajan Btech., 2012); (Rao & Yip, 1990).

Some authors further multiply the x0 term by √2 and multiply the resulting matrix by an overall scale factor of \(\sqrt{2/N}\) (see above for the corresponding change in DCT-II), so that the DCT-II and DCT-III are transposes of one another. This makes the DCT-III matrix orthogonal, but breaks the direct correspondence with a real-even DFT of half-shifted output.

The DCT-III implies the boundary conditions: xn is even around n=0 and odd around n=N; Xk is even around k=-1/2 and even around k=N-1/2.

See the special article DCTIII for details

DCT-IV

\[X_k = \sum_{n=0}^{N-1} x_n \cos \left[\frac{\pi}{N} \left(n+\frac{1}{2}\right) \left(k+\frac{1}{2}\right) \right] \quad \quad k = 0, \dots, N-1.\]

The DCT-IV matrix becomes orthogonal (and thus, being clearly symmetric, its own inverse) if one further multiplies by an overall scale factor of \(\sqrt{2/N}\).

A variant of the DCT-IV, where data from different transforms are overlapped, is called the modified discrete cosine transform (MDCT) (Malvar, 1992).

The DCT-IV implies the boundary conditions: xn is even around n=-1/2 and odd around n=N-1/2; similarly for Xk.

No free C++ implementation for DCTIV is found. This routine still should be written.

DCT V-VIII

DCT types I-IV are equivalent to real-even DFTs of even order (regardless of whether N is even or odd), since the corresponding DFT is of length 2(N−1) (for DCT-I) or 4N (for DCT-II/III) or 8N (for DCT-VIII). In principle, there are actually four additional types of discrete cosine transform (Martucci, 1994), corresponding essentially to real-even DFTs of logically odd order, which have factors of N±½ in the denominators of the cosine arguments.

Equivalently, DCTs of types I-IV imply boundaries that are even/odd around either a data point for both boundaries or halfway between two data points for both boundaries. DCTs of types V-VIII imply boundaries that even/odd around a data point for one boundary and halfway between two data points for the other boundary.

However, these variants seem to be rarely used in practice. One reason, perhaps, is that FFT algorithms for odd-length DFTs are generally more complicated than FFT algorithms for even-length DFTs (e.g. the simplest radix-2 algorithms are only for even lengths), and this increased intricacy carries over to the DCTs as described below.

(The trivial real-even array, a length-one DFT (odd length) of a single number a, corresponds to a DCT-V of length N=1.)

Inverse

the inverse of DCT-I is DCT-I multiplied by 2/(N-1). The inverse of DCT-IV is DCT-IV multiplied by 2/N. The inverse of DCT-II is DCT-III multiplied by 2/N and vice versa. (See e.g. Rao & Yip, 1990.)

Like for the DFT, the normalization factor in front of these transform definitions is merely a convention and differs between treatments. For example, some authors multiply the transforms by so that the inverse does not require any additional multiplicative factor. Combined with appropriate factors of √2 (see above), this can be used to make the transform matrix orthogonal.

Numerical implementation

The examples of the numerical implementation of operators DCTI, DCTII, DCTIII are already available at TORI; the implementation of DCTIV is not yet ready.

References

  1. http://www3.matapp.unimib.it/corsi-2007-2008/matematica/istituzioni-di-analisi-numerica/jpeg/papers/strang.pdf Gilbert Strang. The Discrete Cosine Transform. SIAM REVIEW (Copyright 1999 by Society for Industrial and Applied Mathematics), Vol. 41, No. 1, pp. 135–147.
  2. http://en.wikipedia.org/wiki/Discrete_cosine_transform