Transformée en cosinus discrète

La transformée en cosinus discrète ou TCD est une transformation proche de la transformée de Fourier discrète. Le noyau de projection est un cosinus...



Catégories :

Transformée - Analyse harmonique discrète - Compression d'image - Transformée du signal

Page(s) en rapport avec ce sujet :

  • Transformée en Cosinus Discrète par Bloc 8x8 (DCT)... transposée de C. Par la suite la transformée par bloc 8x8 se réduit à deux multiplications matricielles :... (source : www-ljk.imag)
  • La transformée en cosinus discrète (DCT) 3.1.2. La méthode JPEG... C'est légèrement comme remplacer des additions par une multiplication.... (source : crteknologies.free)
  • ) de la transformation en cosinus discrète..... DCT de Chen et al. [19]. EFCT. N multiplications additions multiplications additions... (source : infoscience.epfl)

La transformée en cosinus discrète ou TCD (de l'anglais : DCT ou Discrete Cosine Transform) est une transformation proche de la transformée de Fourier discrète (DFT). Le noyau de projection est un cosinus et génère par conséquent des cœfficients réels, au contraire de la DFT, dont le noyau est une exponentielle complexe et qui génère par conséquent des cœfficients complexes. On peut cependant exprimer la DCT selon la DFT, qui est alors appliquée sur le signal symétrisé.

La variante la plus courante de la transformée en cosinus discret est la DCT type-II, fréquemment simplement nommée «la DCT». Son inverse, qui correspond au type-III est fréquemment simplement nommée «IDCT».

Une DCT 2D comparée à une DFT. Les zones claires représentent les cœfficients non-nuls

Applications

Définition

La DCT est une fonction linéaire inversible RNRN ou de manière équivalente une matrice carrée N × N inversible. Il existe plusieurs légères variantes de la DCT. Voici les quatre types les plus connus.

Le développement des algorithmes de calcul rapide des transformées DCT se basent sur la possibilité de décomposer la matrice de définition sous forme d'un produit de matrices dont le calcul est plus simple, et sert à diminuer le nombre de multiplications scalaires, en profitant des identités remarquables de périodicité et symétries des fonctions sinusoïdales. Ainsi, on peut décomposer toute transformée DCT de RN en transformées plus simples en décomposant N en produit de facteurs premiers, et en composant des sous-transformées dans Rnn est l'un de ces facteurs. Surtout, de nombreuses optimisations ont été développées lorsque N est une puissance de 2.

Cela revient à réécrire la matrice N × N sous forme de produit de sous-matrices semblables (disposées en pavage régulier et utilisant par conséquent des cœfficients réels communs ou différenciés seulement par leur signe) et de matrices à cœfficients unitaires ou nuls (-1, 0 ou 1), ces dernières ne nécessitant pas de multiplication.

DCT-I

X_k = \frac{1}{2} \left( x_0 + (-1)ˆk x_{N-1} \right) + \sum_{n=1}ˆ{N-2}{x_n \cos\left[ \frac{\pi}{N-1} n k \right]}

On peut rendre cette transformée orthogonale (à une constante multiplicative près) en multipliant x0 et xN-1 par √2 et réciproquement X0 et XN-1 par 1/√2. Cette normalisation casse cependant la correspondance avec une DFT.

On peut noter que la DCT-I n'est pas définie pour N \le 2, contrairement aux autres types qui sont définis pour tout N positif.

DCT-II

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

Cette variante DCT est la plus courante et la plus utilisée. Elle est le plus souvent simplement nommée «la DCT». De la même manière que pour la DCT-I, on peut rendre cette transformation orthogonale en multipliant X0 par 1/√2. Cette forme normalisée est particulièrement utilisée en pratique mais casse la correspondance avec la DFT.

Exemple pour N = 8

Un développement optimisé de cette transformée pour le cas N=8 (utilisé dans JPEG et MPEG) est obtenu en réécrivant la transformée sous forme matricielle et en factorisant la décomposition, pour diminuer le nombre de multiplications scalaires nécessaires. A titre d'exemple, la décomposition suivante est utilisée pour la factorisation par l'algorithme de Chen[1], ici orthogonalisée (voir remarque ci-dessus)  :

Cœfficients constants de calcul

\begin{pmatrix} C_1 \\ C_2 \\ C_3 \\ C_4 \\ C_5 \\ C_6 \\ C_7 \end{pmatrix}
= \sqrt{\frac{2}{N}} . \begin{pmatrix}
\cos{\frac{ \pi}{16}} \\
\cos{\frac{2\pi}{16}} \\
\cos{\frac{3\pi}{16}} \\
\cos{\frac{4\pi}{16}} \\
\cos{\frac{5\pi}{16}} \\
\cos{\frac{6\pi}{16}} \\
\cos{\frac{7\pi}{16}} \\
\end{pmatrix} \approx
\begin{pmatrix} 0,49039 \\ 0,46194 \\ 0,41573 \\ 0,35355 \\ 0,27779 \\ 0,19134 \\ 0,09755 \end{pmatrix}
DCT (8) (méthode de calcul rapide)

\begin{pmatrix} X_0 \\ X_2 \\ X_4 \\ X_6 \end{pmatrix} =
\begin{bmatrix}
C_4 &  C_4 &  C_4 &  C_4 \\
C_2 &  C_6 & -C_6 & -C_2 \\
C_4 & -C_4 & -C_4 &  C_4 \\
C_6 & -C_2 &  C_2 & -C_6
\end{bmatrix} .
\begin{pmatrix} x_0 + x_7 \\ x_1 + x_6 \\ x_2 + x_5 \\ x_3 + x_4 \end{pmatrix}

\begin{pmatrix} X_1 \\ X_3 \\ X_5 \\ X_7 \end{pmatrix} =
\begin{bmatrix}
C_1 &  C_3 &  C_5 &  C_7 \\
C_3 & -C_7 & -C_1 & -C_5 \\
C_5 & -C_1 &  C_7 &  C_3 \\
C_7 & -C_5 &  C_3 & -C_1
\end{bmatrix} .
\begin{pmatrix} x_0 - x_7 \\ x_1 - x_6 \\ x_2 - x_5 \\ x_3 - x_4 \end{pmatrix}

La formule optimisée pour une DCT unidimensionnelle est fréquemment utilisée telle quelle pour son utilisation dans l'espace bidimensionnel (par transposition et composition)  ; cette formule sert à diminuer de façon spectaculaire le calcul de 1024 multiplications (formule de base) à 256 multiplications uniquement dans le traitement d'un bloc image 8×8 (deux passes de 32 multiplications pour chaque ligne de 8 valeurs)  ; cependant, des optimisations sont toujours envisageable en optimisant la composition elle même des deux passes (horizontale et verticale) pour diminuer toujours de 256 à 91 multiplications uniquement (voire moins selon des recherches plus récentes).

On notera aussi que la première matrice ci-dessus permet aussi une réécriture de nombreuses multiplications communes (et par conséquent en fait la formule ci-dessus nécessite nettement moins que les 32 multiplications, en fait 16 si on regroupe les sous-expressions communes). On pourrait toujours décomposer aisément la première matrice car elle est elle-même une transformée DCT dans R4, décomposable en deux sous-matrices de R2.

De nombreuses études ont montré comment cette transformée peut être optimisée suivant les contraintes, surtout lorsque la transformée est utilisée pour la compression, car la transformée sert à concentrer la majeure partie de l'énergie dans les cœfficients obtenus xi d'indice faible, les autres concentrant peu d'énergie ont une contribution faible sur le signal spatial d'origine et sont réduits à zéro lors des étapes de quantification. Ainsi, la précision indispensable pour représenter les derniers cœfficients est plus faible ou alors nulle, et les cœfficients constants Ci utilisés pour le calcul des multiplications scalaires peuvent faire l'objet d'optimisation spécifique, en fixant leur précision, et en utilisant des techniques de multiplication par un nombre réduit d'additions-décalages sans avoir besoin d'utiliser une multiplication générique.

Néanmoins, cet algorithme de calcul (présenté tel quel, il calcule la DCT unidimensionnelle à 8 points avec 16 multiplications) est à la base de l'ensemble des optimisations suivantes par factorisation des sous-matrices. L'algorithme de Lœffler[2] est aujourd'hui le plus efficace ayant été publié (avec 11 multiplications pour la même DCT à 8 points au lieu de 16 avec l'algorithme de Chen, cependant certains cœfficients subissent deux multiplications et cela pourrait rendre l'algorithme moins stable). Il a même été demontré que le nombre minimum théorique de multiplications nécessaires pour la transformation DCT à 8 points ne peut être inférieur à 11, ce qui fait que les algorithmes à 11 multiplications scalaires sont optimums en termes de performance brute (ils se différencient uniquement en termes de stabilité suivant l'ordre dans lequel les multiplications sont réalisées, et la précision par conséquent indispensable pour les produits intermédiaires).

Cependant l'algorithme de Lœffler regroupe 8 des 11 multiplications scalaires sur les sorties, ce qui sert à regrouper ces multiplications avec l'étape suivante de quantification (ce qui en fait tout l'intérêt)  : pour une transformée 2D 8×8, il faut 8×11 multiplications pour la transformée des lignes, et uniquement 8×3 multiplications pour les colonnes, soit un total de 112 multiplications (au lieu de 256 avec l'algorithme de Chen) si les 64 dernières multiplications scalaires sont effectuées avec la quantification. Plus de détails sont disponibles dans les normes de compression JPEG et MPEG.

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]}

La DCT-III est la transformée inverse de la DCT-II. Elle est plus connue sous le nom de «DCT Inverse» et son acronyme (anglais) "IDCT".

De la même manière que pour la DCT-I, on peut rendre cette transformation orthogonale en multipliant x0 par √2. Cette forme normalisée est particulièrement utilisée en pratique mais casse la correspondance avec la DFT.

Exemple pour N = 8

En reprenant l'exemple ci-dessus, on obtient une décomposition inverse (ici orthogonalisée) aussi utilisée dans l'algorithme de Chen :

IDCT (8) (méthode de calcul rapide)

\begin{pmatrix} x_0 \\ x_1 \\ x_2 \\ x_3 \end{pmatrix}
=
\begin{bmatrix}
C_4 &  C_2 &  C_4 &  C_6 \\
C_4 &  C_6 & -C_4 & -C_2 \\
C_4 & -C_6 & -C_4 &  C_2 \\
C_4 & -C_2 &  C_4 & -C_6
\end{bmatrix}
.
\begin{pmatrix} X_0 \\ X_2 \\ X_4 \\ X_6 \end{pmatrix}
+
\begin{bmatrix}
C_1 &  C_3 &  C_5 &  C_7 \\
C_3 & -C_7 & -C_1 & -C_5 \\
C_5 & -C_1 &  C_7 &  C_3 \\
C_7 & -C_5 &  C_3 & -C_1
\end{bmatrix}
.
\begin{pmatrix} X_1 \\ X_3 \\ X_5 \\ X_7 \end{pmatrix}

\begin{pmatrix} x_7 \\ x_6 \\ x_5 \\ x_4 \end{pmatrix}
=
\begin{bmatrix}
C_4 &  C_2 &  C_4 &  C_6 \\
C_4 &  C_6 & -C_4 & -C_2 \\
C_4 & -C_6 & -C_4 &  C_2 \\
C_4 & -C_2 &  C_4 & -C_6
\end{bmatrix}
.
\begin{pmatrix} X_0 \\ X_2 \\ X_4 \\ X_6 \end{pmatrix}
-
\begin{bmatrix}
C_1 &  C_3 &  C_5 &  C_7 \\
C_3 & -C_7 & -C_1 & -C_5 \\
C_5 & -C_1 &  C_7 &  C_3 \\
C_7 & -C_5 &  C_3 & -C_1
\end{bmatrix}
.
\begin{pmatrix} X_1 \\ X_3 \\ X_5 \\ X_7 \end{pmatrix}

Là aussi, l'évaluation scalaire de ce produit matriciel contient de nombreuses sous-expressions communes permettant des réductions du nombre de multiplications scalaires nécessaires.

/*
 * Calcul d'une IDCT 2D en 2 IDCT 1D sur les lignes et les colonnes.
 * (Ici en flottant simple précision, 12 multiplications par IDCT 1D)
 * Factorisation utilisée ici (Normalisée avec c4=1): 
 *
 *   // Partie gauche
 *   f0 = (x0+x4);
 *   f1 = (x0-x4);
 *   f2 = c6*(x2+x6) + (c2-c6)*x2;
 *   f3 = c6*(x2+x6) + (-c2-c6)*x6;
 *
 *   // Partie droite
 *   e0 = c3*(x1+x3+x5+x7) + ( c5 - c3)*(x1+x5) + ( c1+c3-c5-c7)*x1 + ( c7 - c3)*(x1+x7);
 *   e1 = c3*(x1+x3+x5+x7) + (-c5 - c3)*(x3+x7) + ( c1+c3+c5-c7)*x3 + (-c1 - c3)*(x3+x5);
 *   e2 = c3*(x1+x3+x5+x7) + ( c5 - c3)*(x1+x5) + ( c1+c3-c5+c7)*x5 + (-c1 - c3)*(x3+x5);
 *   e3 = c3*(x1+x3+x5+x7) + (-c5 - c3)*(x3+x7) + (-c1+c3+c5-c7)*x7 + ( c7 - c3)*(x1+x7); 
 *
 *   out[0] = f0 + f2 + e0;
 *   out[1] = f1 + f3 + e1;
 *   out[2] = f1 - f3 + e2;
 *   out[3] = f0 - f2 + e3;
 *   out[4] = f0 - f2 - e3;
 *   out[5] = f1 - f3 - e2;
 *   out[6] = f1 + f3 - e1;
 *   out[7] = f0 + f2 - e0;
 */
 
#define CLAMP(i) if (i & 0xFF00) i = (((~i) >> 15) & 0xFF);
 
// Normalised to c4=1
#define c3   1.175875602f
#define c6   0.541196100f
 
#define k1   0.765366865f   //  c2-c6
#define k2  -1.847759065f   // -c2-c6
#define k3  -0.390180644f   //  c5-c3
#define k4  -1.961570561f   // -c5-c3
#define k5   1.501321110f   //  c1+c3-c5-c7
#define k6   2.053119869f   //  c1+c3-c5+c7
#define k7   3.072711027f   //  c1+c3+c5-c7
#define k8   0.298631336f   // -c1+c3+c5-c7
#define k9  -0.899976223f   //  c7-c3
#define k10 -2.562915448f   // -c1-c3
 
void idct(short *block,unsigned char *dest)
{
  float x0,x1,x2,x3,x4,x5,x6,x7;
  float e0,e1,e2,e3;
  float f0,f1,f2,f3;
  float x26,x1357,x15,x37,x17,x35;
  float out[64];
  short v;
 
  // Lignes
  for(int i=0;i<8;i++) {
 
    x0 = (float)block[0+i*8];    x1 = (float)block[1+i*8];
    x2 = (float)block[2+i*8];    x3 = (float)block[3+i*8];
    x4 = (float)block[4+i*8];    x5 = (float)block[5+i*8];
    x6 = (float)block[6+i*8];    x7 = (float)block[7+i*8];
 
    // Partie gauche
    f0 = (x0+x4);
    f1 = (x0-x4);
    x26 = c6*(x2+x6);
    f2 = x26 + k1*x2;
    f3 = x26 + k2*x6;
 
    // Partie droite
    x1357 = c3*(x1+x3+x5+x7);
    x15 = k3*(x1+x5);
    x37 = k4*(x3+x7);
    x17 = k9*(x1+x7);
    x35 = k10*(x3+x5);
 
    e0 = x1357 + x15 + k5*x1 + x17;
    e1 = x1357 + x37 + k7*x3 + x35;
    e2 = x1357 + x15 + k6*x5 + x35;
    e3 = x1357 + x37 + k8*x7 + x17;
 
    out[0+i*8] = f0 + f2 + e0;
    out[1+i*8] = f1 + f3 + e1;
    out[2+i*8] = f1 - f3 + e2;
    out[3+i*8] = f0 - f2 + e3;
    out[7+i*8] = f0 + f2 - e0;
    out[6+i*8] = f1 + f3 - e1;
    out[5+i*8] = f1 - f3 - e2;
    out[4+i*8] = f0 - f2 - e3;
 
  }
 
  // Colonnes
  for(int i=0;i<8;i++) {
 
    x0 = out[i+8*0];    x1 = out[i+8*1];
    x2 = out[i+8*2];    x3 = out[i+8*3];
    x4 = out[i+8*4];    x5 = out[i+8*5];
    x6 = out[i+8*6];    x7 = out[i+8*7];
 
    // Partie gauche
    f0 = (x0+x4);
    f1 = (x0-x4);
    x26 = c6*(x2+x6);
    f2 = x26 + k1*x2;
    f3 = x26 + k2*x6;
 
    // Partie droite
    x1357 = c3*(x1+x3+x5+x7);
    x15 = k3*(x1+x5);
    x37 = k4*(x3+x7);
    x17 = k9*(x1+x7);
    x35 = k10*(x3+x5);
 
    e0 = x1357 + x15 + k5*x1 + x17;
    e1 = x1357 + x37 + k7*x3 + x35;
    e2 = x1357 + x15 + k6*x5 + x35;
    e3 = x1357 + x37 + k8*x7 + x17;
 
    // Tronque dans l'intervalle [0,255]  
    // Note: La normalisation avec c4=1 implique une division par 8
    // Le +1028 sert à recentrer et à arrondir la valeur du pixel
    // v = (short)((.+.+.)/8.0 + 128.5f);
 
    v = (short)((f0 + f2 + e0) + 1028.0f) >> 3;
    CLAMP(v); dest[i+8*0] = (uchar)v;
    v = (short)((f1 + f3 + e1) + 1028.0f) >> 3;
    CLAMP(v); dest[i+8*1] = (uchar)v;
    v = (short)((f1 - f3 + e2) + 1028.0f) >> 3;
    CLAMP(v); dest[i+8*2] = (uchar)v;
    v = (short)((f0 - f2 + e3) + 1028.0f) >> 3;
    CLAMP(v); dest[i+8*3] = (uchar)v;
    v = (short)((f0 - f2 - e3) + 1028.0f) >> 3;
    CLAMP(v); dest[i+8*4] = (uchar)v;
    v = (short)((f1 - f3 - e2) + 1028.0f) >> 3;
    CLAMP(v); dest[i+8*5] = (uchar)v;
    v = (short)((f1 + f3 - e1) + 1028.0f) >> 3;
    CLAMP(v); dest[i+8*6] = (uchar)v;
    v = (short)((f0 + f2 - e0) + 1028.0f) >> 3;
    CLAMP(v); dest[i+8*7] = (uchar)v;
 
  }
 
}

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]}

La DCT-IV est une matrice orthogonale.

Références

  1. W. Chen, C. H. Smith, and S. C. Fralick, «A fast computational algorithm for the discrete cosine transform», IEEE Trans. Commun. , Vol. COM-25, p. 1004-1009, Sep. 1977.
  2. C. Lœffler, A. Ligtenberg et G. Moschytz, «Practical Fast 1D DCT Algorithms With 11 Multiplications», dans Proceedings of the International Conference on Acoustics, Speech, and Signal Processing, p. 988--991, 1989

Bibliographie

Voir aussi

Liens externes

Recherche sur Amazon (livres) :



Ce texte est issu de l'encyclopédie Wikipedia. Vous pouvez consulter sa version originale dans cette encyclopédie à l'adresse http://fr.wikipedia.org/wiki/Transform%C3%A9e_en_cosinus_discr%C3%A8te.
Voir la liste des contributeurs.
La version présentée ici à été extraite depuis cette source le 07/04/2010.
Ce texte est disponible sous les termes de la licence de documentation libre GNU (GFDL).
La liste des définitions proposées en tête de page est une sélection parmi les résultats obtenus à l'aide de la commande "define:" de Google.
Cette page fait partie du projet Wikibis.
Accueil Recherche Aller au contenuDébut page
ContactContact ImprimerImprimer liens d'évitement et raccourcis clavierAccessibilité
Aller au menu