1 <?xml version="1.0" encoding="UTF-8"?>
3 * Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
4 * Copyright (C) 1997 - INRIA
5 * Copyright (C) 2012 - Serge Steer - INRIA
7 * This file must be used under the terms of the CeCILL.
8 * This source file is licensed as described in the file COPYING, which
9 * you should have received as part of this distribution. The terms
10 * are also available at
11 * http://www.cecill.info/licences/Licence_CeCILL_V2-en.txt
14 <refentry xmlns="http://docbook.org/ns/docbook"
15 xmlns:xlink="http://www.w3.org/1999/xlink"
16 xmlns:svg="http://www.w3.org/2000/svg"
17 xmlns:mml="http://www.w3.org/1998/Math/MathML"
18 xmlns:db="http://docbook.org/ns/docbook"
19 xmlns:scilab="http://www.scilab.org" xml:lang="en" xml:id="dct">
21 <refname>dct</refname>
22 <refpurpose>Discrete cosine transform.</refpurpose>
24 <refnamediv xml:id="idct">
25 <refname>idct</refname>
26 <refpurpose>Inverse discrete cosine transform.</refpurpose>
29 <title>Calling Sequence</title>
30 <synopsis>X=dct(A [,sign] [,option])
31 X=dct(A,sign,selection [,option])
32 X=dct(A,sign,dims,incr [,option])
34 X=idct(A,selection [,option])
35 X=idct(A,dims,incr [,option])
39 <title>Arguments</title>
44 <para>a real or complex vector or real or complex array
45 (vector, matrix or N-D array.
52 a real or complex array with same shape as <literal>A</literal>.
58 an integer. with possible values <literal>1</literal> or
59 <literal>-1</literal>. Select direct or inverse
60 transform. The default value is <literal>-1</literal>
65 <term>selection</term>
67 a vector containing index on <literal>A</literal> array
68 dimensions. See the Description part for details.
75 a vector of positive numbers with integer values, or a
76 vector of positive integers. See the Description part for details.
78 Each element must be a divisor
79 of the total number of elements of <literal>A</literal>.
82 The product of the elements must be less than the total
83 number of elements of <literal>A</literal>.
90 a vector of positive numbers with integer values, or a
91 vector of positive integers. See the Description part for
94 <literal>incr</literal> must have the same number of
95 elements than <literal>dims</literal>.
98 Each element must be a divisor of the total number of
99 elements of <literal>A</literal>.
102 The <literal>incr</literal> elements must be in strictly
110 a character string. with possible values
111 <literal>"dct1"</literal>, <literal>"dct2"</literal>,
112 <literal>"dct4"</literal> or <literal>"dct"</literal> for
113 direct transform and <literal>"dct1"</literal>,
114 <literal>"dct3"</literal>, <literal>"dct4"</literal> or
115 <literal>"idct"</literal> for inverse transform. The
116 default value is <literal>"dct"</literal> for direct
117 transform and <literal>"idct"</literal> for inverse
118 transform. See the Description part for details.
124 <title>Description</title>
126 <title>Transform description</title>
128 This function realizes direct or
129 inverse 1-D or N-D Discrete Cosine Transforms with shift depending on the <literal>option</literal> parameter value. For a 1-D array <latex>$A$</latex> of length <latex>$n$</latex>:
134 For <literal>"dct1"</literal> the function computes the unnormalized DCT-I transform:
138 $X(k) = X(1)+(-1)^{k-1}X(n)+2\sum_{i=2}^{n-1} {A(i)
139 \cos\frac{\pi (i -1)(k-1)}{n-1}}, k=1\ldots n$
145 For <literal>"dct2"</literal> the function computes the unnormalized DCT-II transform:
149 $X(k) =2 \sum_{i=1}^{n} {A(i) \cos\frac{\pi (i
150 -1/2)(k-1)}{n}}, k=1\ldots n$
156 For <literal>"dct3"</literal> the function computes the unnormalized DCT-III transform:
160 $X(k) =X(1)+ 2 \sum_{i=2}^{n} {A(i) \cos\frac{\pi (i
161 -1)(k-1/2)}{n}}, k=1\ldots n$
169 For <literal>"dct4"</literal> the function computes the unnormalized DCT-IV transform:
173 $X(k) =2 \sum_{i=1}^{n} {A(i) \cos\frac{\pi (i
174 -1/2)(k-1/2)}{n}}, k=1\ldots n$
180 For <literal>"dct"</literal> the function computes the normalized DCT-II transform:
185 \sum_{i=1}^n {A(i) \cos\frac{\pi (i
186 -1/2)(k-1)}{n}}, k=1\ldots n \quad\text{with }
187 \omega(1)=\frac{1}{\sqrt{n}} \quad\text{and }
188 \omega(k)=\sqrt{\frac{2}{n}} , k>1$
195 For <literal>"idct"</literal> the function computes the normalized DCT-III transform:
199 $X(i) = \sum_{k=1}^n {\omega(k) A(k) \cos\frac{\pi (i
200 -1/2)(k-1)}{n}}, k=1\ldots n \quad\text{with }
201 \omega(1)=\frac{1}{\sqrt{n}} \quad\text{and }
202 \omega(k)=\sqrt{\frac{2}{n}} , k>1$
208 The multi-dimensional DCT transforms , in general, are the
209 separable product of the given 1d transform along each dimension
210 of the array. For unnormalized transforms , computing the
211 forward followed by the backward/inverse multi-dimensional
212 transform will result in the original array scaled by the
213 product of the dimension sizes.
216 The normalized multi-dimensional DCT transform of an array
217 <literal>A</literal> with dimensions <latex>$n_1, n_2, \ldots, n_p$</latex>
222 $\begin{array} \\X(k_1, \dots, k_p) =
223 \omega_1(k_1) \ldots \omega_p(k_p)
224 \sum_{i_1=1}^{n_1} \ldots \sum_{i_p=1}^{n_p}
225 {A(i_1,\ldots ,i_p) \cos\frac{\pi (2 i_1
226 -1)(k_1-1)}{2 n} \ldots \cos\frac{\pi (2 i_p
227 -1)(k_p-1)}{2 n}}, k_j=1\ldots n_j\\
229 \omega_j(1)=\frac{1}{\sqrt{n_j}}\\
230 \omega_j(k)=\sqrt{\frac{2}{n_j}} , k>1
235 The normalized multi-dimensional DCT inverse transform of an
236 array <literal>A</literal> with dimensions <latex>$n_1, n_2, \ldots, n_p$</latex>is given by
240 $\begin{array} \\X(i_1, \dots, i_p) = \sum_{k_1=1}^{n_1}
241 \ldots \sum_{k_p=1}^{n_p} {\omega_1(i_1) \ldots \omega_p(i_p)
242 A(k_1,\ldots ,k_p) \cos\frac{\pi (2 k_1 -1)(i_1-1)}{2 n}
243 \ldots \cos\frac{\pi (2 k_p -1)(i_p-1)}{2 n}}, i_j=1\ldots
246 \omega_j(1)=\frac{1}{\sqrt{n_j}}\\
247 \omega_j(k)=\sqrt{\frac{2}{n_j}} , k>1 \end{array}$
253 <title>Syntax description</title>
256 <term>Short syntax </term>
263 <literal>X=dct(A,-1 [,option])</literal> or
264 <literal>X=dct(A [,option])</literal> gives a direct
265 transform according to the option value. The default is normalized DCT-II direct transform.
268 If <literal>A</literal> is a vector (only one
269 dimension greater than 1) a 1-d transform is performed
270 and in the other cases a n-dimensional transform is
274 (the <literal>-1</literal> argument refers
275 to the sign of the exponent..., NOT to
284 <literal>X=dct(A,1 [,option])</literal> or
285 <literal>X=idct(A [,option])</literal>performs the inverse
289 If <literal>A</literal> is a vector (only one
290 dimension greater than 1) a 1-d transform is performed
291 and in the other cases a n-dimensional transform is
301 <term>Long syntax for DCT along specified dimensions</term>
306 <literal>X=dct(A,sign,selection [,option])</literal>
307 allows to perform efficiently all direct or inverse
308 dct of the "slices" of <literal>A</literal> along
312 For example, if <literal>A</literal> is a 3-D array
313 <literal>X=dct(A,-1,2)</literal> is equivalent to:
315 <programlisting role=""><![CDATA[
318 X(i1,:,i3)=dct(A(i1,:,i3),-1);
323 and <literal>X=dct(A,-1,[1 3])</literal> is equivalent to:
325 <programlisting role=""><![CDATA[
327 X(:,i2,:)=dct(A(:,i2,:),-1);
333 <literal>X=dct(A,sign,dims,incr)</literal> is
334 an old syntax that also allows to perform all direct or
335 inverse dct of the slices of <literal>A</literal> along
339 For example, if <literal>A</literal> is an array with
340 <literal>n1*n2*n3</literal> elements
341 <literal>X=dct(A,-1,n1,1)</literal> is equivalent to
342 <literal>X=dct(matrix(A,[n1,n2,n3]),-1,1)</literal>.
343 and <literal>X=dct(A,-1,[n1 n3],[1 n1*n2])</literal>
345 <literal>X=dct(matrix(A,[n1,n2,n3]),-1,[1,3])</literal>.
354 <title>Optimizing dct</title>
356 Remark: function automatically stores his last parameters in
357 memory to re-use it in a second time. This improves greatly
358 the time computation when consecutives calls (with same
359 parameters) are performed.
362 It is possible to go further in dct optimization using
363 <link linkend="get_fftw_wisdom">get_fftw_wisdom</link>, <link
364 linkend="set_fftw_wisdom">set_fftw_wisdom</link> functions.
369 <title>Algorithms</title>
371 This function uses the <ulink url="http://www.fftw.org/">fftw3</ulink> library.
375 <title>Examples</title>
377 <programlisting role="example"><![CDATA[
378 //Frequency components of a signal
379 //----------------------------------
380 // build a sampled at 1000hz containing pure frequencies
383 t = 0:1/sample_rate:0.6;
384 N=size(t,'*'); //number of samples
385 s=sin(2*%pi*50*t)+sin(2*%pi*70*t+%pi/4)+grand(1,N,'nor',0,1);
387 // zero low energy components
389 size(find(y1<>0),'*') //only 30 non zero coefficients out of 600
390 clf;plot(s,'b'),plot(dct(d,1),'r')
395 <programlisting role="example"><![CDATA[
397 A=eval3d(milk_drop,x,x);
402 clf();fig=gcf();fig.color_map=graycolormap(128);
403 subplot(121),grayplot(x,x,A)
404 subplot(122),grayplot(x,x,A1)
409 <refsection role="see also">
410 <title>See Also</title>
411 <simplelist type="inline">
413 <link linkend="fft">fft</link>
416 <link linkend="dst">dst</link>
419 <link linkend="fftw_flags">fftw_flags</link>
422 <link linkend="get_fftw_wisdom">get_fftw_wisdom</link>
425 <link linkend="set_fftw_wisdom">set_fftw_wisdom</link>
428 <link linkend="fftw_forget_wisdom">fftw_forget_wisdom</link>
433 <title>Bibliography</title>
435 Matteo Frigo and Steven G. Johnson, "FFTW Documentation" <ulink
436 url="http://www.fftw.org/#documentation">http://www.fftw.org/#documentation</ulink>