5 Transformées et Compression

Dans les chapitres précédents, toutes les opérations ont été réalisées dans le domaine spatial, où les algorithmes agissent directement sur les valeurs d’intensité des pixels.

Ce chapitre présente une approche complémentaire : le domaine fréquentiel, dans lequel l’image est représentée par les variations spatiales d’intensité, et non seulement par les valeurs individuelles des pixels.

Le concept de fréquence spatiale décrit la rapidité avec laquelle l’intensité varie le long de l’image. Les variations lentes correspondent aux basses fréquences, tandis que les contours, les détails fins et les bruits correspondent aux hautes fréquences.

Cette représentation repose sur le fait que toute image numérique discrète peut être décomposée en une combinaison de fonctions orthogonales. La Transformée de Fourier utilise une base d’exponentielles complexes bidimensionnelles (équivalentes à des sinusoïdes avec une orientation et une fréquence spécifiques). D’autres transformées, comme la Transformée en Cosinus (DCT) et la Transformée Wavelet (DWT), utilisent différentes familles de fonctions de base — des cosinus bidimensionnels dans le cas de la DCT, et des fonctions à support compact dans le cas des wavelets.

Parmi les principales applications de cette représentation, on distingue :

Le filtrage dans le domaine fréquentiel, pour atténuer ou rehausser certaines bandes de fréquences ;
L’analyse multirésolution au moyen de transformées wavelet, qui représente les structures à différentes échelles ;
La compression d’images, par la réduction du nombre de coefficients nécessaires pour représenter l’image.

5.1 Objectifs

À la fin de ce chapitre, vous serez capable de :

Interpréter le spectre de Fourier d’une image, en distinguant l’amplitude, la phase et les composantes de fréquence ;
Appliquer le théorème de convolution pour effectuer un filtrage dans le domaine fréquentiel à l’aide de la transformée de Fourier rapide (FFT) ;
Concevoir et analyser des filtres dans le domaine fréquentiel, en comprenant le fonctionnement des filtres passe-bas, passe-haut et coupe-bande ;
Comprendre l’analyse multirésolution par transformées en ondelettes et son application à la représentation hiérarchique des images ;
Décrire le processus de compression d’images, y compris la transformée en cosinus discrète (DCT) et la quantification des coefficients ;
Choisir des formats de stockage d’images, tels que JPEG, PNG et WebP, en fonction des exigences de l’application.

5.2 Configuration de l’environnement

import os, urllib.request

os.makedirs("tmp/state", exist_ok=True)  # artefacts de build du parcours C++

url = "https://raw.githubusercontent.com/fzampirolli/pdi-vc/master/morph/config.py"
if not os.path.exists("config.py"):
    urllib.request.urlretrieve(url, "config.py")

# Noyau Python même dans le parcours C++. cpp=True télécharge morph.hpp + stb ; les
# cellules %%writefile *.cpp de ce chapitre compilent AVEC OpenCV
# (-DMM_USE_OPENCV + pkg-config opencv4).
import config
config.setup(cpp=True)
from morph import mm
import numpy as np

✅ Environnement prêt. Morph : 1.1.9 | OpenCV : 5.0.0

5.3 Transformada de Fourier Discreta 2D

L’analyse de Fourier repose sur le principe selon lequel tout signal périodique peut être représenté comme une somme de fonctions sinusoïdales de différentes fréquences, amplitudes et phases. Ce concept s’applique également aux images numériques, permettant de les représenter dans le domaine fréquentiel plutôt que dans le domaine spatial.

La Figure 5.1 illustre cette décomposition pour un signal unidimensionnel. Dans le cas d’une image, la Transformée de Fourier Discrète (TFD) convertit la matrice d’intensités \(f(x,y)\) en un ensemble de coefficients décrivant la contribution des différentes fréquences spatiales présentes dans l’image.

%%writefile tmp/fig_decomposicao_1d.cpp
#define MM_OUT "tmp/fig_decomposicao_1d.png"
//| label: fig-decomposicao-1d
//| fig-cap: "Décomposition de Fourier 1D: une onde carrée (ligne pointillée) est approchée par la somme des premières sinusoïdes (lignes colorées). Plus il y a de termes, meilleure est l'approximation."
//| echo: false
//| output: true

#include <iostream>
#include <vector>
#include <cmath>
#include "morph.hpp"
#include <filesystem>

int main() {
    // Créer l'axe x de 0 à 2π avec 400 points
    std::vector<double> x(400);
    for (int i = 0; i < 400; ++i) {
        x[i] = (2.0 * M_PI) * i / 399.0;
    }

    // Onde carrée : 1 pour x < π, -1 sinon
    std::vector<double> square(400);
    for (int i = 0; i < 400; ++i) {
        square[i] = (x[i] < M_PI) ? 1.0 : -1.0;
    }

    // Somme des harmoniques
    std::vector<double> soma(400, 0.0);
    std::vector<std::vector<double>> ys;
    std::vector<double> h1(400), h2(400), h3(400);

    for (int n = 1; n <= 3; ++n) {
        double coeff = (4.0 / M_PI) * (1.0 / (2 * n - 1));
        std::vector<double> h(400);
        for (int i = 0; i < 400; ++i) {
            h[i] = coeff * std::sin((2 * n - 1) * x[i]);
            soma[i] += h[i];
        }
        if (n == 1) h1 = h;
        else if (n == 2) h2 = h;
        else h3 = h;
    }

    // Préparer les courbes
    ys = {square, h1, h2, h3, soma};

    std::vector<std::string> labels = {
        "Onde carrée idéale", "1ère harmonique", "3ème harmonique", 
        "5ème harmonique", "Somme (3 premières)"
    };
    std::vector<cv::Scalar> colors = {
        cv::Scalar(40, 40, 40), cv::Scalar(60, 160, 80), 
        cv::Scalar(180, 120, 60), cv::Scalar(150, 80, 160), 
        cv::Scalar(60, 60, 220)
    };

    // Créer le graphique
    mm::Image chart = mm::lineChart(
        x, ys, labels, colors,
        "Synthèse de Fourier : des sinus à une onde carrée",
        "Position", "Intensité"
    );

    mm::show(std::vector<mm::Image>{chart}, MM_OUT, 
             std::vector<std::string>{"Décomposition de Fourier 1D"}, 1);

    
// [pdi:panel-io] auto-generated — do not edit by hand
std::filesystem::create_directories("tmp");
mm::write(chart, "tmp/fig_decomposicao_1d_0.png");
// [pdi:panel-io:end]
return 0;
}

Overwriting tmp/fig_decomposicao_1d.cpp

!g++ -I. -std=c++17 -DMM_USE_OPENCV -I/usr/include/opencv4 tmp/fig_decomposicao_1d.cpp -o tmp/fig_decomposicao_1d -lopencv_stitching -lopencv_alphamat -lopencv_aruco -lopencv_barcode -lopencv_bgsegm -lopencv_bioinspired -lopencv_ccalib -lopencv_dnn_objdetect -lopencv_dnn_superres -lopencv_dpm -lopencv_face -lopencv_freetype -lopencv_fuzzy -lopencv_hdf -lopencv_hfs -lopencv_img_hash -lopencv_intensity_transform -lopencv_line_descriptor -lopencv_mcc -lopencv_quality -lopencv_rapid -lopencv_reg -lopencv_rgbd -lopencv_saliency -lopencv_shape -lopencv_stereo -lopencv_structured_light -lopencv_phase_unwrapping -lopencv_superres -lopencv_optflow -lopencv_surface_matching -lopencv_tracking -lopencv_highgui -lopencv_datasets -lopencv_text -lopencv_plot -lopencv_ml -lopencv_videostab -lopencv_videoio -lopencv_viz -lopencv_wechat_qrcode -lopencv_ximgproc -lopencv_video -lopencv_xobjdetect -lopencv_objdetect -lopencv_calib3d -lopencv_imgcodecs -lopencv_features2d -lopencv_dnn -lopencv_flann -lopencv_xphoto -lopencv_photo -lopencv_imgproc -lopencv_core \
  && ./tmp/fig_decomposicao_1d \
  && test -f "tmp/fig_decomposicao_1d.png" \
  || echo "⚠ mm::show não gravou tmp/fig_decomposicao_1d.png"

[1] Décomposition de Fourier 1D

try:
    mm.show(
        [
            mm.read("tmp/fig_decomposicao_1d_0.png"),
        ],
        titles=[
            'Decomposicao de Fourier 1D',
        ],
        cols=1,
    )
except Exception as _e:
    print("figura indisponivel nesta trilha (C++): " + repr(_e) + " tmp/fig_decomposicao_1d_0.png (ver a versao Python)")

Figure 5.1: Decomposição de Fourier 1D: uma onda quadrada (linha tracejada) é aproximada pela soma das primeiras senoides (linhas coloridas). Quanto mais termos, melhor a aproximação.

5.3.1 Simulateur : Reconstruire des signaux avec des sinusoïdes

Avant d’étudier les images bidimensionnelles, le simulateur de la Figure 5.2 illustre le principe de l’analyse de Fourier pour des signaux unidimensionnels : une forme d’onde peut être approchée par la somme de sinusoïdes de différentes fréquences et amplitudes.

À mesure que de nouveaux termes sont ajoutés, la somme des sinusoïdes (courbe noire) se rapproche de la forme d’onde de référence (en pointillés). Le graphique inférieur présente le spectre d’amplitudes, indiquant la contribution de chaque fréquence à la reconstruction du signal.

Activité

Explorez le simulateur et répondez :

Combien de termes sont nécessaires pour obtenir une bonne approximation de l’onde carrée ?
Laquelle des trois formes d’onde converge le plus rapidement ? Justifiez votre réponse.
Comment le spectre d’amplitudes se modifie-t-il en remplaçant l’onde carrée par l’onde triangulaire ?

Réponses

1. Combien de termes sont nécessaires pour une bonne approximation de l’onde carrée ?

Avec environ 15 à 20 termes, la forme de l’onde se rapproche déjà bien de la référence. Cependant, à proximité des discontinuités subsiste une petite oscillation, connue sous le nom de phénomène de Gibbs, qui ne disparaît pas même avec l’ajout de termes supplémentaires.

2. Quelle forme converge le plus rapidement ? Pourquoi ?

L’onde triangulaire converge plus rapidement, car les amplitudes de ses harmoniques décroissent plus vite que celles de l’onde carrée et de l’onde en dents de scie. Par conséquent, peu de termes suffisent déjà à produire une bonne approximation.

3. Comment le spectre change-t-il entre l’onde carrée et l’onde triangulaire ?

Toutes deux ne possèdent que des harmoniques impairs, mais, dans l’onde triangulaire, les amplitudes diminuent beaucoup plus rapidement. Ainsi, peu d’harmoniques suffisent pour reconstruire le signal avec une bonne précision.

∿ Simulateur : Décomposition de Fourier 1D somme de sinusoïdes

Termes

Erreur RMS

–

Forme cible

carrée

Forme cible

Nombre de termes

Affichage

Figure 5.2: Simulateur interactif de la décomposition de Fourier 1D : visualisation de la somme de sinusoïdes avec différentes fréquences, amplitudes et phases. Ajoutez des termes et observez la convergence vers des formes d’onde arbitraires.

5.3.2 Interprétation du spectre de fréquences

En appliquant la Transformée de Fourier Discrète (TFD) à une image et en visualisant le module de ses coefficients (voir Figure 5.5), on obtient le spectre de magnitude, qui montre la distribution des fréquences spatiales présentes dans l’image.

Le coefficient situé à l’origine de la TFD, appelé composante continue (Direct Current), correspond à la fréquence nulle et représente l’intensité moyenne de l’image. Par convention, ce coefficient est stocké dans le coin supérieur gauche du spectre. Pour faciliter son interprétation, on applique l’opération FFT Shift, qui déplace la composante continue vers le centre de l’image. Après ce décalage, les basses fréquences se concentrent dans la région centrale, tandis que les hautes fréquences se trouvent près des bords, comme le résume la Table 5.1.

Table 5.1: Correspondance entre les régions du spectre de magnitude après application du FFT Shift.

Région du spectre	Composantes prédominantes	Exemples dans l’image
Centre (basses fréquences)	Variations spatiales lentes	Éclairage, régions homogènes et formes globales
Région intermédiaire (fréquences moyennes)	Variations à échelle intermédiaire	Textures et motifs répétitifs
Bords (hautes fréquences)	Variations spatiales rapides	Contours, détails fins et bruit

Cette organisation facilite l’interprétation du spectre et la conception de filtres. L’atténuation des basses fréquences réduit les variations globales d’intensité, tandis que l’atténuation des hautes fréquences lisse l’image en réduisant les détails fins et une partie du bruit.

5.3.3 L’Expérience de la Grille : Construire une Image à partir d’un Unique Coefficient

Avant de présenter la formulation mathématique de la Transformée de Fourier Discrète (TFD), il est utile d’analyser son inverse, appelée Transformée de Fourier Discrète Inverse (TFDI). Considérons un spectre où tous les coefficients sont nuls, à l’exception d’un seul. Un exemple de cette construction est présenté dans le code de la Figure 5.3 et peut être exploré de manière interactive dans le simulateur de la Figure 5.4..

L’image reconstruite est une sinusoïde bidimensionnelle. La position du coefficient dans le spectre détermine son orientation et sa fréquence spatiale, tandis que sa magnitude et sa phase définissent, respectivement, son amplitude et son déplacement spatial. Ainsi, chaque coefficient de la TFD représente une composante sinusoïdale, et l’image originale peut être reconstruite par la somme de toutes ces composantes.

%%writefile tmp/fig_05_grade_2d.cpp
#define MM_OUT "tmp/fig_05_grade_2d.png"
//| label: fig-05-grade-2d
//| fig-cap: "Toute fréquence dans le spectre (point isolé) correspond à une onde sinusoïdale 2D rotatée dans le domaine spatial."
//| echo: true
//| output: true

#include <opencv2/opencv.hpp>
#include <vector>
#include <string>
#include <cmath>
#include "morph.hpp"
#include <filesystem>

// Helper: inverse fftshift (swap quadrants)
cv::Mat ifftshift(const cv::Mat& input) {
    int cx = input.cols / 2;
    int cy = input.rows / 2;
    cv::Mat output = input.clone();

    // Top-left quadrant
    cv::Mat q0(input, cv::Rect(0, 0, cx, cy));
    // Top-right quadrant
    cv::Mat q1(input, cv::Rect(cx, 0, input.cols - cx, cy));
    // Bottom-left quadrant
    cv::Mat q2(input, cv::Rect(0, cy, cx, input.rows - cy));
    // Bottom-right quadrant
    cv::Mat q3(input, cv::Rect(cx, cy, input.cols - cx, input.rows - cy));

    // Swap quadrants: TL->BR, TR->BL, BL->TR, BR->TL
    cv::Mat tmp;
    q0.copyTo(tmp);
    q3.copyTo(q0);
    tmp.copyTo(q3);

    q1.copyTo(tmp);
    q2.copyTo(q1);
    tmp.copyTo(q2);

    return output;
}

int main() {
    int N_grid = 100;
    cv::Mat espectro_vazio = cv::Mat::zeros(N_grid, N_grid, CV_64FC2);

    // Allumant un seul point (fréquence) hors du centre
    int u0 = 10, v0 = 5;
    espectro_vazio.at<cv::Vec2d>(N_grid/2 - v0, N_grid/2 - u0) = cv::Vec2d(1000, 0);

    // Retour au domaine spatial (IDFT)
    cv::Mat espectro_shifted = ifftshift(espectro_vazio);
    cv::Mat onda_2d_complex;
    cv::idft(espectro_shifted, onda_2d_complex, cv::DFT_SCALE | cv::DFT_COMPLEX_OUTPUT);

    // Extraire la partie réelle
    cv::Mat onda_2d_parts[2];
    cv::split(onda_2d_complex, onda_2d_parts);
    cv::Mat onda_2d = onda_2d_parts[0];

    // Normalisation pour visualisation
    cv::Mat onda_vis;
    cv::normalize(onda_2d, onda_vis, 0, 255, cv::NORM_MINMAX, CV_8U);

    // Magnitude du spectre
    cv::Mat espectro_mag[2];
    cv::split(espectro_vazio, espectro_mag);
    cv::Mat espectro_abs;
    cv::magnitude(espectro_mag[0], espectro_mag[1], espectro_abs);

    cv::Mat espectro_vis;
    cv::normalize(espectro_abs, espectro_vis, 0, 255, cv::NORM_MINMAX, CV_8U);

    // Accent visuel du point
    cv::Mat espectro_color;
    cv::cvtColor(espectro_vis, espectro_color, cv::COLOR_GRAY2BGR);
    cv::circle(espectro_color, cv::Point(N_grid/2 - u0, N_grid/2 - v0), 2, cv::Scalar(0, 0, 255), -1);

    mm::show(std::vector<mm::Image>{espectro_color, onda_vis},
             MM_OUT,
             std::vector<std::string>{"Spectre (1 point actif)", "Onde 2D Résultante (IDFT)"},
             2);

    
// [pdi:panel-io] auto-generated — do not edit by hand
std::filesystem::create_directories("tmp");
mm::write(espectro_color, "tmp/fig_05_grade_2d_0.png");
mm::write(onda_vis, "tmp/fig_05_grade_2d_1.png");
// [pdi:panel-io:end]
return 0;
}

Overwriting tmp/fig_05_grade_2d.cpp

!g++ -I. -std=c++17 -DMM_USE_OPENCV -I/usr/include/opencv4 tmp/fig_05_grade_2d.cpp -o tmp/fig_05_grade_2d -lopencv_stitching -lopencv_alphamat -lopencv_aruco -lopencv_barcode -lopencv_bgsegm -lopencv_bioinspired -lopencv_ccalib -lopencv_dnn_objdetect -lopencv_dnn_superres -lopencv_dpm -lopencv_face -lopencv_freetype -lopencv_fuzzy -lopencv_hdf -lopencv_hfs -lopencv_img_hash -lopencv_intensity_transform -lopencv_line_descriptor -lopencv_mcc -lopencv_quality -lopencv_rapid -lopencv_reg -lopencv_rgbd -lopencv_saliency -lopencv_shape -lopencv_stereo -lopencv_structured_light -lopencv_phase_unwrapping -lopencv_superres -lopencv_optflow -lopencv_surface_matching -lopencv_tracking -lopencv_highgui -lopencv_datasets -lopencv_text -lopencv_plot -lopencv_ml -lopencv_videostab -lopencv_videoio -lopencv_viz -lopencv_wechat_qrcode -lopencv_ximgproc -lopencv_video -lopencv_xobjdetect -lopencv_objdetect -lopencv_calib3d -lopencv_imgcodecs -lopencv_features2d -lopencv_dnn -lopencv_flann -lopencv_xphoto -lopencv_photo -lopencv_imgproc -lopencv_core \
  && ./tmp/fig_05_grade_2d \
  && test -f "tmp/fig_05_grade_2d.png" \
  || echo "⚠ mm::show não gravou tmp/fig_05_grade_2d.png"

[1] Spectre (1 point actif)
[2] Onde 2D Résultante (IDFT)

try:
    mm.show(
        [
            mm.read("tmp/fig_05_grade_2d_0.png"),
            mm.read("tmp/fig_05_grade_2d_1.png"),
        ],
        titles=[
            'Espectro (1 ponto ativo)',
            'Onda 2D Resultante (IDFT)',
        ],
        cols=2,
        figsize=(10, 4),
    )
except Exception as _e:
    print("figura indisponivel nesta trilha (C++): " + repr(_e) + " tmp/fig_05_grade_2d_0.png (ver a versao Python)")

Figure 5.3: Toda frequência no espectro (ponto isolado) corresponde a uma onda senoidal 2D rotacionada no domínio espacial.

∿ Simulateur : Synthèse de fréquence 2D (IDFT) Espace de Fourier

Fréquence u

Fréquence v

Distance R

11.18

Angle θ

26.6°

Spectre (Cliquez pour déplacer le point)

➔

Onde 2D résultante (domaine spatial)

Déplacement u (X) : 10

Déplacement v (Y) : 5

Figure 5.4: Simulateur interactif de la synthèse de Fourier 2D. Modifiez la position horizontale (\(u\)) et verticale (\(v\)) du coefficient dans le spectre de fréquences centré et observez comment la distance par rapport au centre dicte la fréquence spatiale (épaisseur) et l’angle dicte l’orientation de l’onde sinusoïdale générée.

5.3.4 Définition mathématique

Considérons une image \(f(x,y)\) de dimensions \(M \times N\). Sa Transformée de Fourier discrète 2D (TFD) est définie par :

\[ F(u,v) = \sum_{x=0}^{M-1} \sum_{y=0}^{N-1} f(x,y)\, e^{-j2\pi\left(\frac{ux}{M}+\frac{vy}{N}\right)} \tag{5.1}\]

où \(u = 0, 1, \ldots, M-1\) et \(v = 0, 1, \ldots, N-1\) représentent les fréquences discrètes dans les directions horizontale et verticale, respectivement. Le terme exponentiel correspond à une sinusoïde bidimensionnelle, dont la fréquence et l’orientation sont déterminées par les indices \((u,v)\).

La Transformée de Fourier discrète inverse 2D (TFDI) reconstruit l’image originale à partir de ses coefficients :

\[ f(x,y) = \frac{1}{MN} \sum_{u=0}^{M-1} \sum_{v=0}^{N-1} F(u,v)\, e^{j2\pi\left(\frac{ux}{M}+\frac{vy}{N}\right)} \tag{5.2}\]

Les équations Équation 5.1 et Équation 5.2 montrent que la TFD et la TFDI forment une paire de transformations : la première convertit l’image dans le domaine fréquentiel, tandis que la seconde reconstruit exactement l’image originale à partir de ses coefficients.

À propos du symbole \(j\)

Le terme \(j\) désigne l’unité imaginaire, définie par \(j^2 = -1\). En ingénierie et en traitement du signal, on adopte \(j\) plutôt que \(i\) afin d’éviter toute confusion avec la notation du courant électrique. Son utilisation dans l’exponentielle complexe, régie par la formule d’Euler (\(e^{j\theta} = \cos\theta + j\sin\theta\)), permet de représenter de manière compacte l’amplitude et la phase de chaque fréquence spatiale présente dans l’image.

Qu’est-ce que la composante DC ?

Le coefficient \(F(0,0)\), appelé composante DC (courant continu), est égal à la somme des intensités de tous les pixels de l’image (voir Figure 5.5) :

\[ F(0,0)=MN\,\bar{f}, \]

où \(\bar{f}\) est l’intensité moyenne de l’image. Ainsi, la composante DC représente le niveau moyen d’intensité et, dans la plupart des images naturelles, possède la plus grande amplitude du spectre.

Les autres coefficients représentent les variations autour de cette moyenne. Après application du décalage de FFT, la composante DC est déplacée vers le centre du spectre, concentrant les basses fréquences dans la région centrale et les hautes fréquences sur les bords.

Anatomie du spectre de Fourier 2D (après fftshift)

Figure 5.5: Diagramme conceptuel du spectre de Fourier 2D centré.

5.3.5 Magnitude et Phase

Chaque coefficient de la Transformée de Fourier Discrète (TFD) est un nombre complexe et peut s’écrire sous la forme

\[ F(u,v)=R(u,v)+j\,I(u,v), \]

où \(R(u,v)\) et \(I(u,v)\) correspondent respectivement aux parties réelle et imaginaire du coefficient. À partir de l’équation Équation 5.1, on obtient

\[ R(u,v)= \sum_{x=0}^{M-1}\sum_{y=0}^{N-1} f(x,y) \cos\!\left( 2\pi\left(\frac{ux}{M}+\frac{vy}{N}\right) \right), \]

\[ I(u,v)= - \sum_{x=0}^{M-1}\sum_{y=0}^{N-1} f(x,y) \sin\!\left( 2\pi\left(\frac{ux}{M}+\frac{vy}{N}\right) \right). \]

À partir de cette représentation, on définit deux grandeurs fondamentales :

Magnitude, qui indique l’intensité de la composante fréquentielle,

\[ |F(u,v)|=\sqrt{R(u,v)^2+I(u,v)^2}; \]

Phase, qui détermine l’alignement (ou le décalage) spatial de la composante,

\[ \phi(u,v)=\operatorname{atan2}\!\left(I(u,v),\,R(u,v)\right). \]

Ainsi, chaque coefficient peut également être écrit sous sa forme polaire,

\[ F(u,v)=|F(u,v)|\,e^{j\phi(u,v)}. \]

Le spectre de Fourier peut donc être visualisé au moyen de deux images distinctes : le spectre de magnitude, généralement utilisé pour analyser la distribution des fréquences, et le spectre de phase, qui décrit l’organisation spatiale des composantes sinusoïdales.

Bien que le spectre de magnitude soit le plus utilisé pour l’inspection visuelle, la phase contient une grande partie des informations structurelles de l’image. La combinaison de la magnitude et de la phase permet de reconstruire exactement l’image originale au moyen de la TFD inverse.

5.3.6 O que a Magnitude e a Fase carregam?

5.3.7 Ce que portent la magnitude et la phase ?

Une démonstration classique consiste à combiner la magnitude d’une image avec la phase d’une autre et à reconstruire le résultat. Cette expérience met en évidence que :

La phase préserve la structure spatiale de l’image, y compris la position des objets, leurs contours et leur géométrie. De petites modifications de la phase peuvent provoquer de grands changements visuels.
La magnitude contrôle la manière dont l’énergie est distribuée entre les fréquences spatiales, influençant principalement le contraste et la texture.

Lorsqu’une image est reconstruite avec la magnitude de A et la phase de B, le résultat tend à ressembler davantage à B qu’à A, ce qui montre que la phase est le principal composant responsable de l’organisation spatiale de la scène. Cependant, la magnitude reste importante, car elle module le contraste des structures reconstruites. Ainsi, une reconstruction fidèle dépend de la combinaison cohérente entre magnitude et phase.

Un exemple de ce comportement est présenté dans la Figure 5.6.

Analogie avec l’audio : limites et précautions

La phase d’un signal joue des rôles distincts dans l’audio et les images :

Audio stéréo ou multicanal : la phase relative entre les canaux est fondamentale pour la perception de la position des sources sonores, via les différences interaurales de temps (ITD, Interaural Time Differences).
Audio monaural : la phase absolue exerce peu d’influence perceptuelle directe.
Images (TFD) : la phase est le principal facteur responsable de l’organisation spatiale de la scène, tandis que la magnitude module le contraste et la distribution de l’énergie entre les fréquences.

Dans les deux domaines, la magnitude est liée à l’intensité des composantes de fréquence : en audio, elle influence le timbre et l’intensité perçue ; en images, elle influence le contraste et la texture.

%%writefile tmp/fig_05_dft_intro.cpp
#define MM_OUT "tmp/fig_05_dft_intro.png"
//| label: fig-05-dft-intro
//| fig-cap: "Experimento de troca de fase: Imagem A (moedas) e Imagem B (padrão geométrico) reconstruídas com magnitudes e fases trocadas. O resultado mostra que a estrutura visual é **muito mais sensível à fase** do que à magnitude: quando a fase de B é mantida, a imagem resultante preserva a organização espacial de B, mesmo com a magnitude de A. A magnitude, por sua vez, influencia principalmente o contraste e a textura. Observe que a qualidade da reconstrução não é perfeita — há artefatos visíveis —, evidenciando a interdependência entre fase e magnitude para uma representação fiel da imagem."
//| echo: true
//| output: true

#include <opencv2/opencv.hpp>
#include <opencv2/core.hpp>
#include <opencv2/imgproc.hpp>
#include <opencv2/highgui.hpp>
#include <iostream>
#include <string>
#include <vector>
#include <cmath>
#include <complex>
#include <filesystem>
#include <algorithm>
#include "morph.hpp"

// Helper para trocar quadrantes (equivalente a np.fft.fftshift em 2D)
static cv::Mat fftShift(const cv::Mat& mat) {
    int cx = mat.cols / 2;
    int cy = mat.rows / 2;
    cv::Mat output;
    cv::Mat q0(mat, cv::Rect(0, 0, cx, cy));
    cv::Mat q1(mat, cv::Rect(cx, 0, mat.cols - cx, cy));
    cv::Mat q2(mat, cv::Rect(0, cy, cx, mat.rows - cy));
    cv::Mat q3(mat, cv::Rect(cx, cy, mat.cols - cx, mat.rows - cy));
    cv::Mat tmp;
    q0.copyTo(tmp);
    q3.copyTo(q0);
    tmp.copyTo(q3);
    q1.copyTo(tmp);
    q2.copyTo(q1);
    tmp.copyTo(q2);
    mat.copyTo(output);
    return output;
}

// Helper para construir imagem a partir de magnitude e fase
static cv::Mat buildFromMagPhase(const cv::Mat& mag, const cv::Mat& phase) {
    cv::Mat magF, phaseF;
    mag.convertTo(magF, CV_64F);
    phase.convertTo(phaseF, CV_64F);

    // Calcular cos e sin da fase
    cv::Mat cosPhase, sinPhase;
    cv::Mat phaseArr[] = {phaseF, phaseF};
    cv::Mat cosArr[] = {cosPhase, sinPhase};

    // Aplicar cos e sin elemento a elemento via loop
    cosPhase = cv::Mat(phaseF.size(), CV_64F);
    sinPhase = cv::Mat(phaseF.size(), CV_64F);
    for (int i = 0; i < phaseF.rows; i++) {
        for (int j = 0; j < phaseF.cols; j++) {
            double p = phaseF.at<double>(i, j);
            cosPhase.at<double>(i, j) = std::cos(p);
            sinPhase.at<double>(i, j) = std::sin(p);
        }
    }

    // real = mag * cos(phase), imag = mag * sin(phase)
    cv::Mat realPart = magF.mul(cosPhase);
    cv::Mat imagPart = magF.mul(sinPhase);

    // Construir complexa: real + i*imag
    cv::Mat complexParts[] = {realPart, imagPart};
    cv::Mat complexMat;
    cv::merge(complexParts, 2, complexMat);

    // IDFT
    cv::Mat result;
    cv::idft(complexMat, result, cv::DFT_SCALE | cv::DFT_REAL_OUTPUT);

    // Converter para 8-bit
    cv::Mat result8u;
    result.convertTo(result8u, CV_8U);

    return result8u;
}

int main() {
    // ── Experimento: A Importância da Fase ───────────────────────────────────────
    // ── Carregamento da imagem ────────────────────────────────────────────────────
    std::string url     = "https://upload.wikimedia.org/wikipedia/commons/2/25/GAZI.MD.AHAD_11.jpg";
    std::string caminho = "imagens/coins.jpg";

    mm::Image img_obj;
    if (!std::filesystem::exists(caminho)) {
        std::filesystem::create_directories("imagens");
        img_obj = mm::read(url);
        mm::write(img_obj, caminho);
    } else {
        img_obj = mm::read(caminho);
    }

    cv::Mat img_color = img_obj;
    if (img_color.channels() == 3) {
        cv::Mat temp;
        cv::cvtColor(img_color, temp, cv::COLOR_BGR2GRAY);
        img_color = temp;
    }
    mm::Image img_gray_mm = mm::gray(img_obj);
    cv::Mat img_gray = img_gray_mm;

    cv::Mat img_a;
    cv::resize(img_gray, img_a, cv::Size(400, 400));

    // Criar uma imagem B sintética (padrão geométrico)
    cv::Mat img_b = cv::Mat::zeros(400, 400, CV_8UC1);
    cv::rectangle(img_b, cv::Point(100, 100), cv::Point(300, 300), cv::Scalar(255), -1);
    cv::circle(img_b, cv::Point(200, 200), 150, cv::Scalar(128), 10);

    // FFT de A e B
    cv::Mat img_a_f;
    img_a.convertTo(img_a_f, CV_64F);
    cv::Mat img_b_f;
    img_b.convertTo(img_b_f, CV_64F);

    cv::Mat FA, FB;
    cv::dft(img_a_f, FA, cv::DFT_COMPLEX_OUTPUT);
    cv::dft(img_b_f, FB, cv::DFT_COMPLEX_OUTPUT);

    // Separar magnitude e fase de FA
    cv::Mat FA_parts[2];
    cv::split(FA, FA_parts);
    cv::Mat magA, phaseA;
    cv::magnitude(FA_parts[0], FA_parts[1], magA);
    cv::phase(FA_parts[0], FA_parts[1], phaseA);

    // Separar magnitude e fase de FB
    cv::Mat FB_parts[2];
    cv::split(FB, FB_parts);
    cv::Mat magB, phaseB;
    cv::magnitude(FB_parts[0], FB_parts[1], magB);
    cv::phase(FB_parts[0], FB_parts[1], phaseB);

    // Troca de Fase: rec_A_mag_B_fase = |FA| * exp(i * phase(FB))
    cv::Mat rec_A_mag_B_fase = buildFromMagPhase(magA, phaseB);

    // Troca de Fase: rec_B_mag_A_fase = |FB| * exp(i * phase(FA))
    cv::Mat rec_B_mag_A_fase = buildFromMagPhase(magB, phaseA);

    // Converter para mm::Image para exibição
    mm::Image img_a_mm = img_a;
    mm::Image img_b_mm = img_b;
    mm::Image rec_A_mm = rec_A_mag_B_fase;
    mm::Image rec_B_mm = rec_B_mag_A_fase;

    // Manter a variável img_gray como mm::Image no final (para persistência)
    img_gray = img_gray_mm;

    // Exibir resultados
    mm::show(
        std::vector<mm::Image>{img_a_mm, img_b_mm, rec_A_mm, rec_B_mm},
        MM_OUT,
        std::vector<std::string>{"Imagem A", "Imagem B", "Mag(A) + Fase(B)", "Mag(B) + Fase(A)"},
        4
    );

    // Exibir mensagem
    std::cout << "💡 A fase preserva bordas e contornos; a magnitude controla contraste e" << std::endl;
    std::cout << "textura. Em áudio estéreo, a fase afeta a localização espacial; em" << std::endl;
    std::cout << "imagens, determina a organização da cena." << std::endl;

    
// [pdi:state-io] auto-generated — do not edit by hand
std::filesystem::create_directories("tmp/state");
mm::write(img_gray, "tmp/state/img_gray_20.png");
// [pdi:state-io:end]

// [pdi:panel-io] auto-generated — do not edit by hand
std::filesystem::create_directories("tmp");
mm::write(img_a, "tmp/fig_05_dft_intro_0.png");
mm::write(img_b, "tmp/fig_05_dft_intro_1.png");
mm::write(rec_A_mag_B_fase, "tmp/fig_05_dft_intro_2.png");
mm::write(rec_B_mag_A_fase, "tmp/fig_05_dft_intro_3.png");
// [pdi:panel-io:end]
return 0;
}

Overwriting tmp/fig_05_dft_intro.cpp

!g++ -I. -std=c++17 -DMM_USE_OPENCV -I/usr/include/opencv4 tmp/fig_05_dft_intro.cpp -o tmp/fig_05_dft_intro -lopencv_stitching -lopencv_alphamat -lopencv_aruco -lopencv_barcode -lopencv_bgsegm -lopencv_bioinspired -lopencv_ccalib -lopencv_dnn_objdetect -lopencv_dnn_superres -lopencv_dpm -lopencv_face -lopencv_freetype -lopencv_fuzzy -lopencv_hdf -lopencv_hfs -lopencv_img_hash -lopencv_intensity_transform -lopencv_line_descriptor -lopencv_mcc -lopencv_quality -lopencv_rapid -lopencv_reg -lopencv_rgbd -lopencv_saliency -lopencv_shape -lopencv_stereo -lopencv_structured_light -lopencv_phase_unwrapping -lopencv_superres -lopencv_optflow -lopencv_surface_matching -lopencv_tracking -lopencv_highgui -lopencv_datasets -lopencv_text -lopencv_plot -lopencv_ml -lopencv_videostab -lopencv_videoio -lopencv_viz -lopencv_wechat_qrcode -lopencv_ximgproc -lopencv_video -lopencv_xobjdetect -lopencv_objdetect -lopencv_calib3d -lopencv_imgcodecs -lopencv_features2d -lopencv_dnn -lopencv_flann -lopencv_xphoto -lopencv_photo -lopencv_imgproc -lopencv_core \
  && ./tmp/fig_05_dft_intro \
  && test -f "tmp/fig_05_dft_intro.png" \
  || echo "⚠ mm::show não gravou tmp/fig_05_dft_intro.png"

[1] Imagem A
[2] Imagem B
[3] Mag(A) + Fase(B)
[4] Mag(B) + Fase(A)
💡 A fase preserva bordas e contornos; a magnitude controla contraste e
textura. Em áudio estéreo, a fase afeta a localização espacial; em
imagens, determina a organização da cena.

try:
    mm.show(
        [
            mm.read("tmp/fig_05_dft_intro_0.png"),
            mm.read("tmp/fig_05_dft_intro_1.png"),
            mm.read("tmp/fig_05_dft_intro_2.png"),
            mm.read("tmp/fig_05_dft_intro_3.png"),
        ],
        titles=[
            'Imagem A',
            'Imagem B',
            'Mag(A) + Fase(B)',
            'Mag(B) + Fase(A)',
        ],
        cols=4,
        figsize=(16, 4),
    )
except Exception as _e:
    print("figura indisponivel nesta trilha (C++): " + repr(_e) + " tmp/fig_05_dft_intro_0.png (ver a versao Python)")

Figure 5.6: Experimento de troca de fase: Imagem A (moedas) e Imagem B (padrão geométrico) reconstruídas com magnitudes e fases trocadas. O resultado mostra que a estrutura visual é **muito mais sensível à fase** do que à magnitude: quando a fase de B é mantida, a imagem resultante preserva a organização espacial de B, mesmo com a magnitude de A. A magnitude, por sua vez, influencia principalmente o contraste e a textura. Observe que a qualidade da reconstrução não é perfeita — há artefatos visíveis —, evidenciando a interdependência entre fase e magnitude para uma representação fiel da imagem.

5.4 Théorème de convolution et stratégies de filtrage

Le théorème de convolution établit une relation fondamentale entre les domaines spatial et fréquentiel :

\[ f(x,y) \circledast h(x,y) \;\overset{\mathcal{F}}{\longleftrightarrow}\; F(u,v)\,H(u,v) \tag{5.3}\]

où \(\circledast\) représente la convolution circulaire discrète. Ainsi, la convolution entre une image \(f(x,y)\) et un filtre \(h(x,y)\) peut être remplacée par la multiplication de leurs spectres.

En pratique, pour obtenir le même résultat que la convolution linéaire effectuée dans le domaine spatial, on applique un remplissage par zéros (zero-padding) avant la Transformée de Fourier rapide (FFT), afin d’éviter les artefacts sur les bords de l’image.

Cependant, le filtrage dans le domaine fréquentiel n’est pas toujours l’alternative la plus efficace. Pour des filtres tels que le filtre gaussien et le filtre moyenneur (Box Filter), la propriété de séparabilité permet de réduire considérablement le coût computationnel de la convolution dans le domaine spatial.

5.4.1 Noyau séparable vs. non séparable

Un noyau séparable peut être écrit comme le produit externe de deux vecteurs unidimensionnels,

\[ H = v\,h^T, \]

ce qui permet de remplacer la convolution bidimensionnelle par deux convolutions unidimensionnelles successives : l’une dans la direction horizontale et l’autre dans la direction verticale.

En revanche, un noyau non séparable n’admet pas cette décomposition et, par conséquent, sa convolution doit être réalisée directement sur le voisinage bidimensionnel.

En pratique, pour un noyau de dimension \(K \times K\), la convolution directe exige \(K^2\) multiplications par pixel, tandis qu’un noyau séparable ne requiert que \(2K\) multiplications, réduisant ainsi considérablement le coût computationnel.

5.4.2 Analyse de l’efficacité computationnelle

Considérons une image de dimensions \(M \times N\) et un filtre carré de taille \(K \times K\). La Table 5.2 compare la complexité des principales stratégies de filtrage.

Table 5.2: Comparaison de la complexité de la convolution directe, séparable et via la Transformée de Fourier rapide (FFT).

Méthode de filtrage	Complexité asymptotique	Dépendance de \(K\)	Application typique
Spatiale non séparable	\(\mathcal{O}(MNK^2)\)	Quadratique	Kernels petits et non séparables
Spatiale séparable	\(\mathcal{O}(MNK)\)	Linéaire	Filtres gaussien et de moyenne
Via FFT	\(\mathcal{O}(MN\log(MN))\)	Indépendante de \(K\)	Kernels grands

Pour de petits kernels, la convolution spatiale, surtout lorsque le filtre est séparable, s’avère généralement plus efficace en raison du faible coût des opérations. À mesure que la taille du kernel augmente, le filtrage via FFT devient plus avantageux, car son coût ne dépend pratiquement pas de la dimension du filtre.

5.4.3 Discussion des résultats expérimentaux

Le graphique obtenu lors de l’essai avec l’image des pièces (\(2560 \times 1920\)), présenté à la Figure 5.7, confirme le comportement prédit par l’analyse de complexité computationnelle.

Convolution non séparable (\(\mathcal{O}(MNK^2)\))
La convolution directe présente une croissance quadratique avec la taille du kernel. Pour de petites valeurs de \(K\), le coût est faible, mais il augmente rapidement à mesure que le kernel grandit, devenant irréalisable pour des applications en temps réel.
Filtrage par FFT (\(\mathcal{O}(MN \log(MN))\))
Le coût de la FFT ne dépend que de la taille de l’image, étant indépendant de \(K\). Par conséquent, ses performances restent approximativement constantes lorsque le kernel varie, ce qui la rend avantageuse pour les filtres de grande taille ou non séparables.
Convolution séparable (\(\mathcal{O}(MNK)\))
La décomposition du kernel en deux filtres unidimensionnels réduit significativement le coût computationnel. En pratique, cette approche tend à être la plus efficace pour les filtres séparables, en particulier dans les implémentations optimisées.

En général, le choix de la méthode dépend de la taille et de la structure du kernel. Les filtres séparables sont plus efficaces dans le domaine spatial, tandis que la FFT devient plus avantageuse pour les kernels de grande taille ou pour de multiples convolutions dans le domaine fréquentiel.

\[ g = \mathcal{F}^{-1}\bigl[\mathcal{F}(f)\cdot \mathcal{F}(h)\bigr] \quad \text{(FFT)} \qquad g = f \circledast h \quad \text{(convolution directe)} \qquad g = (f \circledast v) \circledast h^T \quad \text{(séparable)} \tag{5.4}\]

où :

\(f(x,y)\) représente l’image d’entrée ;
\(h(x,y)\) est le kernel bidimensionnel du filtre ;
\(v\) et \(h^T\) sont, respectivement, les vecteurs vertical et horizontal qui composent le kernel séparable.

%%writefile tmp/fig_05_conv_eficiencia.cpp
#define MM_OUT "tmp/fig_05_conv_eficiencia.png"
//| label: fig-05-conv-eficiencia
//| fig-cap: "Comparação de eficiência: Convolução Não Separável (Espacial 2D), Separável (Espacial 1D) e via FFT."
//| echo: false
//| output: true

#include <iostream>
#include <vector>
#include <string>
#include <cmath>
#include "morph.hpp"
#include <filesystem>

int main() {
    // Trilha C++: curvas de CUSTO relativo (contagem de operacoes) — a forma das
    // curvas (K^2 vs 2K vs log N) e o que importa; a trilha py mede tempos reais.
    std::vector<double> K = {3, 7, 11, 15, 21, 31, 41, 51};
    double N2 = 256.0 * 256.0;

    // t_nao_sep = N2 * K^2 / 1e6
    std::vector<double> t_nao_sep(K.size());
    // t_sep = N2 * 2 * K / 1e6
    std::vector<double> t_sep(K.size());
    // t_fft = constante (N2 * log2(N2) / 1e6)
    double fft_val = N2 * std::log2(N2) / 1e6;
    std::vector<double> t_fft(K.size(), fft_val);

    for (size_t i = 0; i < K.size(); ++i) {
        t_nao_sep[i] = N2 * K[i] * K[i] / 1e6;
        t_sep[i]     = N2 * 2.0 * K[i] / 1e6;
    }

    // Construção do gráfico de linhas com mm::lineChart (aceita vetores de vetores)
    std::vector<std::vector<double>> xs = {K, K, K};
    std::vector<std::vector<double>> ys = {t_nao_sep, t_sep, t_fft};
    std::vector<std::string> labels = {
        "Nao Separavel  O(N^2 K^2)", 
        "Separavel  O(N^2 . 2K)", 
        "Via FFT  O(N^2 log N)"
    };

    mm::Image chart = mm::lineChart(
        xs, ys, labels,
        {}, // colors vazio — usa o padrão
        "Custo relativo: Espacial vs Frequencia",
        "Tamanho do kernel  K", 
        "Operacoes  (x10^6)"
    );

    mm::show(std::vector<mm::Image>{chart}, MM_OUT, {"Comparacao de complexidade"}, 1);

    
// [pdi:panel-io] auto-generated — do not edit by hand
std::filesystem::create_directories("tmp");
mm::write(chart, "tmp/fig_05_conv_eficiencia_0.png");
// [pdi:panel-io:end]
return 0;
}

Overwriting tmp/fig_05_conv_eficiencia.cpp

!g++ -I. -std=c++17 -DMM_USE_OPENCV -I/usr/include/opencv4 tmp/fig_05_conv_eficiencia.cpp -o tmp/fig_05_conv_eficiencia -lopencv_stitching -lopencv_alphamat -lopencv_aruco -lopencv_barcode -lopencv_bgsegm -lopencv_bioinspired -lopencv_ccalib -lopencv_dnn_objdetect -lopencv_dnn_superres -lopencv_dpm -lopencv_face -lopencv_freetype -lopencv_fuzzy -lopencv_hdf -lopencv_hfs -lopencv_img_hash -lopencv_intensity_transform -lopencv_line_descriptor -lopencv_mcc -lopencv_quality -lopencv_rapid -lopencv_reg -lopencv_rgbd -lopencv_saliency -lopencv_shape -lopencv_stereo -lopencv_structured_light -lopencv_phase_unwrapping -lopencv_superres -lopencv_optflow -lopencv_surface_matching -lopencv_tracking -lopencv_highgui -lopencv_datasets -lopencv_text -lopencv_plot -lopencv_ml -lopencv_videostab -lopencv_videoio -lopencv_viz -lopencv_wechat_qrcode -lopencv_ximgproc -lopencv_video -lopencv_xobjdetect -lopencv_objdetect -lopencv_calib3d -lopencv_imgcodecs -lopencv_features2d -lopencv_dnn -lopencv_flann -lopencv_xphoto -lopencv_photo -lopencv_imgproc -lopencv_core \
  && ./tmp/fig_05_conv_eficiencia \
  && test -f "tmp/fig_05_conv_eficiencia.png" \
  || echo "⚠ mm::show não gravou tmp/fig_05_conv_eficiencia.png"

[1] Comparacao de complexidade

try:
    mm.show(
        [
            mm.read("tmp/fig_05_conv_eficiencia_0.png"),
        ],
        titles=[
            'Comparacao de complexidade',
        ],
        cols=1,
    )
except Exception as _e:
    print("figura indisponivel nesta trilha (C++): " + repr(_e) + " tmp/fig_05_conv_eficiencia_0.png (ver a versao Python)")

Figure 5.7: Comparação de eficiência: Convolução Não Separável (Espacial 2D), Separável (Espacial 1D) e via FFT.

Le problème de la convolution circulaire (wrap-around)

La Transformée de Fourier Discrète (TFD) suppose que l’image est périodiquement étendue dans l’espace, c’est-à-dire que ses bords se répètent indéfiniment.

Dans cette condition, la multiplication dans le domaine fréquentiel correspond à une convolution circulaire dans le domaine spatial. Par conséquent, des régions opposées de l’image (haut et bas, gauche et droite) interagissent artificiellement, comme illustré dans la Figure 5.8.

L’application d’un zero-padding avant la FFT réduit cet effet en étendant l’image avec des valeurs nulles sur les bords, rapprochant ainsi le résultat de la convolution linéaire. Ce comportement peut être interprété à la lumière du Théorème de Convolution, présenté dans la Figure 5.9.

%%writefile tmp/fig_05_padding_error.cpp
#define MM_OUT "tmp/fig_05_padding_error.png"
//| label: fig-05-padding-error
//| fig-cap: "Sem *padding*, um deslocamento severo faz a imagem vazar para o lado oposto (convolução circular)."
//| echo: true
//| output: true

#include <opencv2/opencv.hpp>
#include <complex>
#include <vector>
#include <string>
#include <cmath>
#include "morph.hpp"
#include <filesystem>

int main() {
// [pdi:state-io] auto-generated — do not edit by hand
mm::Image img_gray = mm::_read_state("tmp/state/img_gray_20.png");
// [pdi:state-io:end]

    // img_gray is provided (mm::Image)
    int M = img_gray.h;
    int N = img_gray.w;

    // Converter para CV_64F para cálculos complexos
    cv::Mat img_gray_mat(img_gray.h, img_gray.w, CV_8UC1, img_gray.data.data());
    cv::Mat img_gray_float;
    img_gray_mat.convertTo(img_gray_float, CV_64F);

    // Simulação de um filtro de deslocamento brutal
    cv::Mat H_shift(M, N, CV_64FC2, cv::Scalar(0, 0));
    for (int u = 0; u < M; ++u) {
        for (int v = 0; v < N; ++v) {
            double angle = -2.0 * M_PI * (u * 120.0 / M + v * 120.0 / N);
            H_shift.at<cv::Vec2d>(u, v) = cv::Vec2d(std::cos(angle), std::sin(angle));
        }
    }

    // Filtragem SEM padding (causa o wrap-around)
    cv::Mat F_img;
    cv::dft(img_gray_float, F_img, cv::DFT_COMPLEX_OUTPUT);

    // Multiplicação no domínio da frequência
    cv::Mat produto;
    cv::mulSpectrums(F_img, H_shift, produto, 0);

    cv::Mat img_vazada;
    cv::idft(produto, img_vazada, cv::DFT_SCALE | cv::DFT_REAL_OUTPUT);

    // Normalizar para visualização
    cv::Mat img_vazada_norm;
    cv::normalize(img_vazada, img_vazada_norm, 0, 255, cv::NORM_MINMAX, CV_8U);

    // Converter para mm::Image
    mm::Image img_vazada_vis(img_vazada_norm);

    // Mostrar resultados
    mm::show(std::vector<mm::Image>{img_gray, img_vazada_vis},
             MM_OUT,
             std::vector<std::string>{"Original", "Filtragem s/ Padding (Vazamento)"},
             2);

    
// [pdi:panel-io] auto-generated — do not edit by hand
std::filesystem::create_directories("tmp");
mm::write(img_gray, "tmp/fig_05_padding_error_0.png");
mm::write(img_vazada_vis, "tmp/fig_05_padding_error_1.png");
// [pdi:panel-io:end]
return 0;
}

Overwriting tmp/fig_05_padding_error.cpp

!g++ -I. -std=c++17 -DMM_USE_OPENCV -I/usr/include/opencv4 tmp/fig_05_padding_error.cpp -o tmp/fig_05_padding_error -lopencv_stitching -lopencv_alphamat -lopencv_aruco -lopencv_barcode -lopencv_bgsegm -lopencv_bioinspired -lopencv_ccalib -lopencv_dnn_objdetect -lopencv_dnn_superres -lopencv_dpm -lopencv_face -lopencv_freetype -lopencv_fuzzy -lopencv_hdf -lopencv_hfs -lopencv_img_hash -lopencv_intensity_transform -lopencv_line_descriptor -lopencv_mcc -lopencv_quality -lopencv_rapid -lopencv_reg -lopencv_rgbd -lopencv_saliency -lopencv_shape -lopencv_stereo -lopencv_structured_light -lopencv_phase_unwrapping -lopencv_superres -lopencv_optflow -lopencv_surface_matching -lopencv_tracking -lopencv_highgui -lopencv_datasets -lopencv_text -lopencv_plot -lopencv_ml -lopencv_videostab -lopencv_videoio -lopencv_viz -lopencv_wechat_qrcode -lopencv_ximgproc -lopencv_video -lopencv_xobjdetect -lopencv_objdetect -lopencv_calib3d -lopencv_imgcodecs -lopencv_features2d -lopencv_dnn -lopencv_flann -lopencv_xphoto -lopencv_photo -lopencv_imgproc -lopencv_core \
  && ./tmp/fig_05_padding_error \
  && test -f "tmp/fig_05_padding_error.png" \
  || echo "⚠ mm::show não gravou tmp/fig_05_padding_error.png"

[1] Original
[2] Filtragem s/ Padding (Vazamento)

try:
    mm.show(
        [
            mm.read("tmp/fig_05_padding_error_0.png"),
            mm.read("tmp/fig_05_padding_error_1.png"),
        ],
        titles=[
            'Original',
            'Filtragem s/ Padding (Vazamento)',
        ],
        cols=2,
        figsize=(10, 4),
    )
except Exception as _e:
    print("figura indisponivel nesta trilha (C++): " + repr(_e) + " tmp/fig_05_padding_error_0.png (ver a versao Python)")

Figure 5.8: Sem *padding*, um deslocamento severo faz a imagem vazar para o lado oposto (convolução circular).

%%writefile tmp/fig_05_conv_teorema.cpp
#define MM_OUT "tmp/fig_05_conv_teorema.png"
//| label: fig-05-conv-teorema
//| fig-cap: "Teorema da Convolution : filtrer dans le domaine spatial (Gaussienne) équivaut à multiplier le spectre par H(u,v) en fréquence. Les sorties coïncident visuellement."
//| echo: true
//| output: true

#include <opencv2/opencv.hpp>
#include <vector>
#include <string>
#include "morph.hpp"
#include <filesystem>

int main() {
// [pdi:state-io] auto-generated — do not edit by hand
mm::Image img_gray = mm::_read_state("tmp/state/img_gray_20.png");
// [pdi:state-io:end]

    // Même lissage par deux chemins : espace vs. fréquence.
    int M = img_gray.h;
    int N = img_gray.w;
    cv::Mat H = mm::gaussFilter(M, N, 30);      // H(u,v) : transfert passe-bas gaussien
    mm::Image f_freq = mm::freqFilter(img_gray, H);  // convolution via multiplication en fréquence
    mm::Image f_esp = mm::gaussian(img_gray, 31, 5); // même lissage réalisé dans l'espace

    mm::show(
        {img_gray, f_esp, f_freq},
        MM_OUT,
        {"Original", "Espace: Gaussienne", "Frequence: H(u,v).F(u,v)"},
        3
    );

    
// [pdi:panel-io] auto-generated — do not edit by hand
std::filesystem::create_directories("tmp");
mm::write(img_gray, "tmp/fig_05_conv_teorema_0.png");
mm::write(f_esp, "tmp/fig_05_conv_teorema_1.png");
mm::write(f_freq, "tmp/fig_05_conv_teorema_2.png");
// [pdi:panel-io:end]
return 0;
}

Overwriting tmp/fig_05_conv_teorema.cpp

!g++ -I. -std=c++17 -DMM_USE_OPENCV -I/usr/include/opencv4 tmp/fig_05_conv_teorema.cpp -o tmp/fig_05_conv_teorema -lopencv_stitching -lopencv_alphamat -lopencv_aruco -lopencv_barcode -lopencv_bgsegm -lopencv_bioinspired -lopencv_ccalib -lopencv_dnn_objdetect -lopencv_dnn_superres -lopencv_dpm -lopencv_face -lopencv_freetype -lopencv_fuzzy -lopencv_hdf -lopencv_hfs -lopencv_img_hash -lopencv_intensity_transform -lopencv_line_descriptor -lopencv_mcc -lopencv_quality -lopencv_rapid -lopencv_reg -lopencv_rgbd -lopencv_saliency -lopencv_shape -lopencv_stereo -lopencv_structured_light -lopencv_phase_unwrapping -lopencv_superres -lopencv_optflow -lopencv_surface_matching -lopencv_tracking -lopencv_highgui -lopencv_datasets -lopencv_text -lopencv_plot -lopencv_ml -lopencv_videostab -lopencv_videoio -lopencv_viz -lopencv_wechat_qrcode -lopencv_ximgproc -lopencv_video -lopencv_xobjdetect -lopencv_objdetect -lopencv_calib3d -lopencv_imgcodecs -lopencv_features2d -lopencv_dnn -lopencv_flann -lopencv_xphoto -lopencv_photo -lopencv_imgproc -lopencv_core \
  && ./tmp/fig_05_conv_teorema \
  && test -f "tmp/fig_05_conv_teorema.png" \
  || echo "⚠ mm::show não gravou tmp/fig_05_conv_teorema.png"

[1] Original
[2] Espace: Gaussienne
[3] Frequence: H(u,v).F(u,v)

try:
    mm.show(
        [
            mm.read("tmp/fig_05_conv_teorema_0.png"),
            mm.read("tmp/fig_05_conv_teorema_1.png"),
            mm.read("tmp/fig_05_conv_teorema_2.png"),
        ],
        titles=[
            'Original',
            'Espaco: Gaussiana',
            'Frequencia: H(u,v).F(u,v)',
        ],
        cols=3,
    )
except Exception as _e:
    print("figura indisponivel nesta trilha (C++): " + repr(_e) + " tmp/fig_05_conv_teorema_0.png (ver a versao Python)")

Figure 5.9: Teorema da Convolução: filtrar no domínio do espaço (Gaussiana) equivale a multiplicar o espectro por \(H(u,v)\) na frequência. As saídas coincidem visualmente.

À propos de la différence numérique

La différence résiduelle de l’ordre de \(10^{-13}\) ne viole pas le Théorème de Convolution, mais reflète des limitations computationnelles inhérentes à l’arithmétique en virgule flottante (double précision, ~\(10^{-16}\)) et à l’ordre des opérations entre les deux méthodes :

Convolution spatiale : somme pondérée des voisins avec des arrondis successifs.
Convolution fréquentielle : implique trois transformées FFT et une multiplication complexe, sujette à des erreurs de troncature et de quantification.

Par conséquent, l’égalité théorique est exacte, mais l’implémentation numérique produit une différence pratiquement nulle (erreur relative < \(10^{-12}\)), confirmant le théorème dans la précision de la machine.

5.5 Filtres dans le Domaine Fréquentiel

Un filtre dans le domaine fréquentiel peut être interprété comme une fonction de transfert appliquée au spectre de l’image. Dans cette représentation, chaque coefficient de fréquence est multiplié par une valeur comprise entre 0 et 1, qui détermine son atténuation ou sa préservation. La forme de cette fonction définit l’effet visuel du filtre.

Coupure brutale et ringing. Les filtres idéaux avec une transition instantanée à une fréquence de coupure \(D_0\) produisent des discontinuités dans le domaine fréquentiel. Cette discontinuité se reflète dans le domaine spatial sous forme d’oscillations près des bords, connues sous le nom de ringing. Cet effet est associé à la convolution avec des fonctions à support infini dans l’espace, comme la fonction sinc, comme illustré dans la Figure 5.10.

Filtres à transition douce. Des alternatives telles que les filtres gaussien et de Butterworth lissent la transition entre les régions de passage et de rejet, réduisant ainsi le ringing. En contrepartie, ce lissage implique une frontière de séparation moins nette entre les fréquences préservées et atténuées.

%%writefile tmp/fig_05_conv_teorema_zoom.cpp
#define MM_OUT "tmp/fig_05_conv_teorema_zoom.png"
//| label: fig-05-conv-teorema-zoom
//| fig-cap: "A Dualidade Perigosa: o corte abrupto na Frequência (cilindro Ideal) vira obrigatoriamente uma *sinc* no espaço. Suas ondulações causam o *ringing* fantasma nas bordas da imagem."
//| echo: true
//| output: true
#include <opencv2/opencv.hpp>
#include "morph.hpp"
#include <vector>
#include <string>
#include <filesystem>

int main() {
    // Filtro Ideal na frequência (cilindro) e sua resposta espacial (sinc 2D).
    int N = 128;
    cv::Mat H_freq = mm::idealFilter(N, N, 20);  // 1 dentro do raio 20, 0 fora
    mm::Image h_space = mm::spatialKernel(H_freq); // ifft2(H) -> ondulações da sinc (ringing)

    mm::Image H_freq_img = H_freq; // convertendo para exibição
    mm::show(
        std::vector<mm::Image>{H_freq_img, h_space},
        MM_OUT,
        std::vector<std::string>{
            "Frequencia: filtro Ideal (cilindro)",
            "Espaco: ondulacoes da sinc (causa do ringing)",
        },
        2
    );

    
// [pdi:panel-io] auto-generated — do not edit by hand
std::filesystem::create_directories("tmp");
mm::write(H_freq, "tmp/fig_05_conv_teorema_zoom_0.png");
mm::write(h_space, "tmp/fig_05_conv_teorema_zoom_1.png");
// [pdi:panel-io:end]
return 0;
}

Overwriting tmp/fig_05_conv_teorema_zoom.cpp

!g++ -I. -std=c++17 -DMM_USE_OPENCV -I/usr/include/opencv4 tmp/fig_05_conv_teorema_zoom.cpp -o tmp/fig_05_conv_teorema_zoom -lopencv_stitching -lopencv_alphamat -lopencv_aruco -lopencv_barcode -lopencv_bgsegm -lopencv_bioinspired -lopencv_ccalib -lopencv_dnn_objdetect -lopencv_dnn_superres -lopencv_dpm -lopencv_face -lopencv_freetype -lopencv_fuzzy -lopencv_hdf -lopencv_hfs -lopencv_img_hash -lopencv_intensity_transform -lopencv_line_descriptor -lopencv_mcc -lopencv_quality -lopencv_rapid -lopencv_reg -lopencv_rgbd -lopencv_saliency -lopencv_shape -lopencv_stereo -lopencv_structured_light -lopencv_phase_unwrapping -lopencv_superres -lopencv_optflow -lopencv_surface_matching -lopencv_tracking -lopencv_highgui -lopencv_datasets -lopencv_text -lopencv_plot -lopencv_ml -lopencv_videostab -lopencv_videoio -lopencv_viz -lopencv_wechat_qrcode -lopencv_ximgproc -lopencv_video -lopencv_xobjdetect -lopencv_objdetect -lopencv_calib3d -lopencv_imgcodecs -lopencv_features2d -lopencv_dnn -lopencv_flann -lopencv_xphoto -lopencv_photo -lopencv_imgproc -lopencv_core \
  && ./tmp/fig_05_conv_teorema_zoom \
  && test -f "tmp/fig_05_conv_teorema_zoom.png" \
  || echo "⚠ mm::show não gravou tmp/fig_05_conv_teorema_zoom.png"

[1] Frequencia: filtro Ideal (cilindro)
[2] Espaco: ondulacoes da sinc (causa do ringing)

try:
    mm.show(
        [
            mm.read("tmp/fig_05_conv_teorema_zoom_0.png"),
            mm.read("tmp/fig_05_conv_teorema_zoom_1.png"),
        ],
        titles=[
            'Frequencia: filtro Ideal (cilindro)',
            'Espaco: ondulacoes da sinc (causa do ringing)',
        ],
        cols=2,
    )
except Exception as _e:
    print("figura indisponivel nesta trilha (C++): " + repr(_e) + " tmp/fig_05_conv_teorema_zoom_0.png (ver a versao Python)")

Figure 5.10: A Dualidade Perigosa: o corte abrupto na Frequência (cilindro Ideal) vira obrigatoriamente uma *sinc* no espaço. Suas ondulações causam o *ringing* fantasma nas bordas da imagem.

5.5.1 Filtres Passe-Bas

Les filtres passe-bas atténuent les composantes haute fréquence, ce qui entraîne un lissage de l’image et une réduction du bruit. Après la centralisation du spectre (FFT Shift), la distance de chaque point au centre est donnée par :

\[ D(u,v) = \sqrt{\left(u - \tfrac{M}{2}\right)^2 + \left(v - \tfrac{N}{2}\right)^2} \tag{5.5}\]

Filtre idéal (LPFI) : \[ H_{\text{ideal}}(u,v) = \begin{cases} 1, & D(u,v) \leq D_0 \\ 0, & D(u,v) > D_0 \end{cases} \tag{5.6}\]

La coupure brutale à \(D_0\) introduit des discontinuités dans le domaine fréquentiel, entraînant des oscillations dans le domaine spatial, connues sous le nom de ringing. Cet effet est associé à la convolution avec des fonctions à support infini.

Filtre gaussien (LPFG) : \[ H_{\text{gauss}}(u,v) = e^{-D^2(u,v)/(2\sigma^2)} \tag{5.7}\]

La régularité de la fonction gaussienne dans le domaine fréquentiel évite les discontinuités, ce qui élimine le ringing et produit une transition progressive entre les fréquences conservées et atténuées.

Filtre de Butterworth (LPFB) d’ordre \(n\) : \[ H_{\text{BW}}(u,v) = \frac{1}{1 + \left[D(u,v)/D_0\right]^{2n}} \tag{5.8}\]

Le paramètre \(n\) contrôle la régularité de la transition entre la transmission et la rejection des fréquences. Les petites valeurs produisent des transitions douces, tandis que les grandes valeurs rapprochent le comportement du filtre idéal, avec un risque accru de ringing. Un exemple comparatif est présenté à la Figure 5.11..

Profils des filtres passe-bas — comparaison visuelle (D₀ = 30)

À mesure que l’ordre du Butterworth augmente, le profil se rapproche du filtre idéal — et le ringing augmente.

Figure 5.11: Filtros passe-bas.

5.5.2 Filtres Passe-Haut et Passe-Bande

Les filtres passe-haut peuvent être obtenus à partir d’un filtre passe-bas complémentaire, défini comme :

\[ H_{\text{HP}}(u,v) = 1 - H_{\text{LP}}(u,v) \]

Ce type de filtre préserve les composantes haute fréquence, en mettant en évidence les contours et les détails, tout en atténuant les régions de variation lente.

Les filtres passe-bande préservent uniquement une bande intermédiaire de fréquences, limitée par deux rayons \(D_L\) et \(D_H\) :

\[ H_{\text{BP}}(u,v) = H_{\text{LP}}^{(D_H)}(u,v)\cdot \left[1 - H_{\text{LP}}^{(D_L)}(u,v)\right] \]

Ce type de filtrage est utile lorsque l’on souhaite supprimer simultanément les composantes basse et haute fréquence, en ne préservant que les structures d’échelle intermédiaire.

Une application importante est la suppression du bruit périodique, dans lequel des motifs réguliers apparaissent comme des pics localisés dans le spectre de magnitude. Ces pics peuvent être atténués à l’aide de filtres notch (réjecteur de bande), positionnés spécifiquement sur les fréquences indésirables.

Des exemples de filtres dans le domaine fréquentiel sont présentés dans le simulateur de la Figure 5.12, Figure 5.13 et Figure 5.14..

🎛️ Simulateur : Filtres dans le domaine fréquentiel Passe-bas / Passe-haut

Fréquence de coupure D₀ 30 Type de filtre

passe-bas passe-haut

Réponse H(D)

Spectre filtré |F · H|

Signal 1D — original vs filtré

Énergie retenue par bande (%)

Figure 5.12: Simulateur interactif de filtres dans le domaine fréquentiel.

%%writefile tmp/fig_05_filtros_freq.cpp
#define MM_OUT "tmp/fig_05_filtros_freq.png"
//| label: fig-05-filtros-freq
//| fig-cap: "Comparação entre filtros passa-baixa: Ideal (D₀=30), Gaussiano (D₀=30) e Butterworth (D₀=30, n=2). Perfis de H(u,v) ao longo de uma linha central e imagens filtradas correspondentes."
//| echo: true
//| output: true

#include <opencv2/opencv.hpp>
#include <vector>
#include <string>
#include <cmath>
#include "morph.hpp"

// Fonction pour visualiser H (normalisation 0-255)
static cv::Mat H_vis(const cv::Mat& H) {
    cv::Mat h8;
    H.convertTo(h8, CV_8UC1, 255.0);
    cv::normalize(h8, h8, 0, 255, cv::NORM_MINMAX);
    return h8;
}

int main() {
// [pdi:state-io] auto-generated — do not edit by hand
mm::Image img_gray = mm::_read_state("tmp/state/img_gray_20.png");
// [pdi:state-io:end]

    int M = img_gray.h;
    int N = img_gray.w;
    int D0 = 30;
    int n_bw = 2;

    // ── Fonctions de transfert (H centré) via les constructeurs de morph.hpp ──
    //   Ideal : 1 si D<=D0, sinon 0
    //   Gaussien : exp(-D^2/(2 D0^2))
    //   Butterworth : 1 / (1 + (D/D0)^(2n))
    cv::Mat H_ideal = mm::idealFilter(M, N, D0);
    cv::Mat H_gauss = mm::gaussFilter(M, N, D0);
    cv::Mat H_bw    = mm::butterFilter(M, N, D0, n_bw);

    // ── Images filtrées via FFT ─────────────────────────────────────────────
    mm::Image img_ideal = mm::freqFilter(img_gray, H_ideal);
    mm::Image img_gauss = mm::freqFilter(img_gray, H_gauss);
    mm::Image img_bw    = mm::freqFilter(img_gray, H_bw);

    // ── Profils de H(u,v) sur la ligne centrale, dessinés avec cv::line ─────
    int linha = M / 2;
    int Wp = 512, Hp = 288;
    cv::Mat perfil(Hp, Wp, CV_8UC3, cv::Scalar(255, 255, 255));

    auto traca = [&](const cv::Mat& row, cv::Scalar cor) {
        cv::Point ant(-1, -1);
        bool first = true;
        for (int x = 0; x < N; ++x) {
            int px = (int)std::round(x * (Wp - 1.0) / (N - 1.0));
            int py = (int)std::round(10 + (1.0 - row.at<double>(0, x)) * (Hp - 20));
            if (!first) {
                cv::line(perfil, ant, cv::Point(px, py), cor, 2, cv::LINE_AA);
            }
            ant = cv::Point(px, py);
            first = false;
        }
    };

    // H_ideal[linha] est une ligne de la matrice CV_64F → extraire comme cv::Mat 1×N
    cv::Mat H_ideal_row = H_ideal.row(linha);
    cv::Mat H_gauss_row = H_gauss.row(linha);
    cv::Mat H_bw_row    = H_bw.row(linha);

    traca(H_ideal_row, cv::Scalar(216, 90, 48));
    traca(H_gauss_row, cv::Scalar(29, 158, 117));
    traca(H_bw_row,    cv::Scalar(83, 74, 183));

    mm::Image img_gray_c = img_gray;
    mm::Image img_ideal_c = img_ideal;
    mm::Image img_gauss_c = img_gauss;
    mm::Image img_bw_c = img_bw;
    mm::Image h_ideal_v = mm::Image(H_vis(H_ideal));
    mm::Image h_gauss_v = mm::Image(H_vis(H_gauss));
    mm::Image h_bw_v = mm::Image(H_vis(H_bw));
    mm::Image perfil_img = mm::Image(perfil);

    mm::show(
        std::vector<mm::Image>{img_gray_c, img_ideal_c, img_gauss_c, img_bw_c,
                                h_ideal_v, h_gauss_v, h_bw_v, perfil_img},
        MM_OUT,
        std::vector<std::string>{
            "Original", "LPF Ideal", "LPF Gaussien", "LPF Butterworth (n=2)",
            "H Ideal", "H Gaussien", "H Butterworth", "Profils H(u,v)"
        },
        4
    );

    return 0;
}

Overwriting tmp/fig_05_filtros_freq.cpp

!g++ -I. -std=c++17 -DMM_USE_OPENCV -I/usr/include/opencv4 tmp/fig_05_filtros_freq.cpp -o tmp/fig_05_filtros_freq -lopencv_stitching -lopencv_alphamat -lopencv_aruco -lopencv_barcode -lopencv_bgsegm -lopencv_bioinspired -lopencv_ccalib -lopencv_dnn_objdetect -lopencv_dnn_superres -lopencv_dpm -lopencv_face -lopencv_freetype -lopencv_fuzzy -lopencv_hdf -lopencv_hfs -lopencv_img_hash -lopencv_intensity_transform -lopencv_line_descriptor -lopencv_mcc -lopencv_quality -lopencv_rapid -lopencv_reg -lopencv_rgbd -lopencv_saliency -lopencv_shape -lopencv_stereo -lopencv_structured_light -lopencv_phase_unwrapping -lopencv_superres -lopencv_optflow -lopencv_surface_matching -lopencv_tracking -lopencv_highgui -lopencv_datasets -lopencv_text -lopencv_plot -lopencv_ml -lopencv_videostab -lopencv_videoio -lopencv_viz -lopencv_wechat_qrcode -lopencv_ximgproc -lopencv_video -lopencv_xobjdetect -lopencv_objdetect -lopencv_calib3d -lopencv_imgcodecs -lopencv_features2d -lopencv_dnn -lopencv_flann -lopencv_xphoto -lopencv_photo -lopencv_imgproc -lopencv_core \
  && ./tmp/fig_05_filtros_freq \
  && test -f "tmp/fig_05_filtros_freq.png" \
  || echo "⚠ mm::show não gravou tmp/fig_05_filtros_freq.png"

[1] Original
[2] LPF Ideal
[3] LPF Gaussien
[4] LPF Butterworth (n=2)
[5] H Ideal
[6] H Gaussien
[7] H Butterworth
[8] Profils H(u,v)

try:
    mm.show(mm.read("tmp/fig_05_filtros_freq.png"), figsize=(16, 9))
except Exception as _e:
    print("figura indisponivel nesta trilha (C++): " + repr(_e) + " tmp/fig_05_filtros_freq.png (ver a versao Python)")

Figure 5.13: Comparação entre filtros passa-baixa: Ideal (D₀=30), Gaussiano (D₀=30) e Butterworth (D₀=30, n=2). Perfis de H(u,v) ao longo de uma linha central e imagens filtradas correspondentes.

%%writefile tmp/fig_05_filtros_passa_alta.cpp
#define MM_OUT "tmp/fig_05_filtros_passa_alta.png"
#include <opencv2/opencv.hpp>
#include <vector>
#include <string>
#include "morph.hpp"
#include <filesystem>

int main() {
// [pdi:state-io] auto-generated — do not edit by hand
mm::Image img_gray = mm::_read_state("tmp/state/img_gray_20.png");
// [pdi:state-io:end]

    //| label: fig-05-filtros-passa-alta
    //| fig-cap: "Filtro passa-alta Gaussiano. (a) Original; (b) Filtro passa-alta (D₀=30) - as bordas das moedas e fundo texturizado são realçados."

    // Filtro passa-alta = complemento do passa-baixa Gaussiano (1 - H_gauss),
    // construído com mm.gaussFilter(..., highpass=true) e aplicado via FFT.
    int M = img_gray.h;
    int N = img_gray.w;
    double D0 = 30;
    cv::Mat H_alta = mm::gaussFilter(M, N, D0, true);
    mm::Image img_alta = mm::freqFilter(img_gray, H_alta);

    mm::show(
        std::vector<mm::Image>{img_gray, img_alta},
        MM_OUT,
        {"Original", "Passa-alta Gaussiano (D0=30)"},
        2
    );

    
// [pdi:panel-io] auto-generated — do not edit by hand
std::filesystem::create_directories("tmp");
mm::write(img_gray, "tmp/fig_05_filtros_passa_alta_0.png");
mm::write(img_alta, "tmp/fig_05_filtros_passa_alta_1.png");
// [pdi:panel-io:end]
return 0;
}

Overwriting tmp/fig_05_filtros_passa_alta.cpp

!g++ -I. -std=c++17 -DMM_USE_OPENCV -I/usr/include/opencv4 tmp/fig_05_filtros_passa_alta.cpp -o tmp/fig_05_filtros_passa_alta -lopencv_stitching -lopencv_alphamat -lopencv_aruco -lopencv_barcode -lopencv_bgsegm -lopencv_bioinspired -lopencv_ccalib -lopencv_dnn_objdetect -lopencv_dnn_superres -lopencv_dpm -lopencv_face -lopencv_freetype -lopencv_fuzzy -lopencv_hdf -lopencv_hfs -lopencv_img_hash -lopencv_intensity_transform -lopencv_line_descriptor -lopencv_mcc -lopencv_quality -lopencv_rapid -lopencv_reg -lopencv_rgbd -lopencv_saliency -lopencv_shape -lopencv_stereo -lopencv_structured_light -lopencv_phase_unwrapping -lopencv_superres -lopencv_optflow -lopencv_surface_matching -lopencv_tracking -lopencv_highgui -lopencv_datasets -lopencv_text -lopencv_plot -lopencv_ml -lopencv_videostab -lopencv_videoio -lopencv_viz -lopencv_wechat_qrcode -lopencv_ximgproc -lopencv_video -lopencv_xobjdetect -lopencv_objdetect -lopencv_calib3d -lopencv_imgcodecs -lopencv_features2d -lopencv_dnn -lopencv_flann -lopencv_xphoto -lopencv_photo -lopencv_imgproc -lopencv_core \
  && ./tmp/fig_05_filtros_passa_alta \
  && test -f "tmp/fig_05_filtros_passa_alta.png" \
  || echo "⚠ mm::show não gravou tmp/fig_05_filtros_passa_alta.png"

[1] Original
[2] Passa-alta Gaussiano (D0=30)

try:
    mm.show(
        [
            mm.read("tmp/fig_05_filtros_passa_alta_0.png"),
            mm.read("tmp/fig_05_filtros_passa_alta_1.png"),
        ],
        titles=[
            'Original',
            'Passa-alta Gaussiano ($D_0=30$)',
        ],
        cols=2,
    )
except Exception as _e:
    print("figura indisponivel nesta trilha (C++): " + repr(_e) + " tmp/fig_05_filtros_passa_alta_0.png (ver a versao Python)")

Figure 5.14: Filtro passa-alta Gaussiano. (a) Original; (b) Filtro passa-alta (D₀=30) - as bordas das moedas e fundo texturizado são realçados.

5.5.3 Suppression du bruit périodique

Le bruit périodique — associé à des interférences électriques, à des motifs réguliers de capteurs ou à des artefacts de numérisation — apparaît dans le spectre de Fourier sous forme de pics ponctuels symétriques autour du centre.

Le filtre rejette-bande (notch) atténue sélectivement ces fréquences, tout en préservant les autres composantes de l’image. Un exemple d’application est présenté dans la Figure 5.15..

%%writefile tmp/fig_05_ruido_periodico.cpp
#define MM_OUT "tmp/fig_05_ruido_periodico.png"
//| label: fig-05-ruido-periodico
//| fig-cap: "Remoção de ruído periódico via filtro *notch* no domínio da frequência: (a) imagem com ruído senoidal, (b) espectro mostrando os picos do ruído, (c) máscara *notch* centrada nos picos, (d) imagem restaurada."
//| echo: true
//| output: true

// Ruído senoidal 2D -> espectro -> máscara notch nos 4 picos simétricos -> restauração.
#include <opencv2/opencv.hpp>
#include <cmath>
#include <iostream>
#include <vector>
#include <string>
#include "morph.hpp"

int main() {
// [pdi:state-io] auto-generated — do not edit by hand
mm::Image img_gray = mm::_read_state("tmp/state/img_gray_20.png");
// [pdi:state-io:end]

    int h_img = img_gray.h;
    int w_img = img_gray.w;

    // Criar meshgrid X, Y (equivalente a np.meshgrid)
    std::vector<std::vector<double>> X(h_img, std::vector<double>(w_img));
    std::vector<std::vector<double>> Y(h_img, std::vector<double>(w_img));
    for (int i = 0; i < h_img; ++i) {
        for (int j = 0; j < w_img; ++j) {
            X[i][j] = j;
            Y[i][j] = i;
        }
    }

    int u0 = 20, v0 = 20;  // frequências exatas do ruído

    // Converter img_gray para cv::Mat e depois para float64
    cv::Mat img_gray_mat(h_img, w_img, CV_8UC1, img_gray.data.data());
    cv::Mat img_gray_float;
    img_gray_mat.convertTo(img_gray_float, CV_64F);

    // Calcular ruído: ruido = 40 * sin(2 * pi * (u0*X/w_img + v0*Y/h_img))
    cv::Mat ruido(h_img, w_img, CV_64F);
    for (int i = 0; i < h_img; ++i) {
        for (int j = 0; j < w_img; ++j) {
            double arg = 2.0 * M_PI * (u0 * j / static_cast<double>(w_img) + v0 * i / static_cast<double>(h_img));
            ruido.at<double>(i, j) = 40.0 * std::sin(arg);
        }
    }

    // img_ruidosa = clip(img_gray + ruido, 0, 255).astype(uint8)
    cv::Mat img_ruidosa_float = img_gray_float + ruido;
    cv::Mat img_ruidosa_uint;
    img_ruidosa_float.convertTo(img_ruidosa_uint, CV_8UC1, 1.0, 0.0);
    cv::Mat img_ruidosa;
    cv::normalize(img_ruidosa_uint, img_ruidosa, 0, 255, cv::NORM_MINMAX, CV_8UC1);

    // Espectro de magnitude (log) da imagem ruidosa
    mm::Image img_ruidosa_mm(img_ruidosa);
    mm::Image mag_vis = mm::spectrumMag(img_ruidosa_mm);

    // Máscara notch: 1.0 em todo o plano, disco de 0.0 em cada um dos 4 picos
    cv::Mat mascara = cv::Mat::ones(h_img, w_img, CV_64F);
    int r_notch = 8;
    int cy = h_img / 2, cx = w_img / 2;

    // Coordenadas relativas para os 4 picos simétricos
    std::vector<std::pair<int, int>> picos = {{v0, u0}, {-v0, -u0}, {v0, -u0}, {-v0, u0}};

    for (const auto& pico : picos) {
        int dy = pico.first;
        int dx = pico.second;
        int py = cy + dy;
        int px = cx + dx;

        for (int i = 0; i < h_img; ++i) {
            for (int j = 0; j < w_img; ++j) {
                double dist = std::sqrt(std::pow(i - py, 2) + std::pow(j - px, 2));
                if (dist <= r_notch) {
                    mascara.at<double>(i, j) = 0.0;
                }
            }
        }
    }

    // Versão visual da máscara (0-255)
    cv::Mat mascara_vis_float;
    cv::normalize(mascara, mascara_vis_float, 0, 255, cv::NORM_MINMAX, CV_8UC1);
    mm::Image mascara_vis(mascara_vis_float);

    // Filtragem: aplica a máscara centrada como filtro no domínio da frequência
    mm::Image img_rest_vis = mm::freqFilter(img_ruidosa_mm, mascara);

    // Calcular PSNR
    cv::Mat img_rest_cv = img_rest_vis;
    double psnr = cv::PSNR(img_gray_mat, img_rest_cv);
    std::cout << "PSNR (original vs restaurada): " << psnr << " dB" << std::endl;

    // Mostrar resultados
    std::vector<mm::Image> imagens = {img_ruidosa_mm, mag_vis, mascara_vis, img_rest_vis};
    std::vector<std::string> titulos = {
        "Com ruído periódico",
        "Espectro (log)",
        "Máscara notch",
        "Restaurada (PSNR=" + std::to_string(psnr).substr(0, 4) + " dB)"
    };
    mm::show(imagens, MM_OUT, titulos, 4);

    return 0;
}

Overwriting tmp/fig_05_ruido_periodico.cpp

!g++ -I. -std=c++17 -DMM_USE_OPENCV -I/usr/include/opencv4 tmp/fig_05_ruido_periodico.cpp -o tmp/fig_05_ruido_periodico -lopencv_stitching -lopencv_alphamat -lopencv_aruco -lopencv_barcode -lopencv_bgsegm -lopencv_bioinspired -lopencv_ccalib -lopencv_dnn_objdetect -lopencv_dnn_superres -lopencv_dpm -lopencv_face -lopencv_freetype -lopencv_fuzzy -lopencv_hdf -lopencv_hfs -lopencv_img_hash -lopencv_intensity_transform -lopencv_line_descriptor -lopencv_mcc -lopencv_quality -lopencv_rapid -lopencv_reg -lopencv_rgbd -lopencv_saliency -lopencv_shape -lopencv_stereo -lopencv_structured_light -lopencv_phase_unwrapping -lopencv_superres -lopencv_optflow -lopencv_surface_matching -lopencv_tracking -lopencv_highgui -lopencv_datasets -lopencv_text -lopencv_plot -lopencv_ml -lopencv_videostab -lopencv_videoio -lopencv_viz -lopencv_wechat_qrcode -lopencv_ximgproc -lopencv_video -lopencv_xobjdetect -lopencv_objdetect -lopencv_calib3d -lopencv_imgcodecs -lopencv_features2d -lopencv_dnn -lopencv_flann -lopencv_xphoto -lopencv_photo -lopencv_imgproc -lopencv_core \
  && ./tmp/fig_05_ruido_periodico \
  && test -f "tmp/fig_05_ruido_periodico.png" \
  || echo "⚠ mm::show não gravou tmp/fig_05_ruido_periodico.png"

PSNR (original vs restaurada): 33.0718 dB
[1] Com ruído periódico
[2] Espectro (log)
[3] Máscara notch
[4] Restaurada (PSNR=33.0 dB)

try:
    mm.show(mm.read("tmp/fig_05_ruido_periodico.png"), figsize=(16, 4))
except Exception as _e:
    print("figura indisponivel nesta trilha (C++): " + repr(_e) + " tmp/fig_05_ruido_periodico.png (ver a versao Python)")

Figure 5.15: Remoção de ruído periódico via filtro *notch* no domínio da frequência: (a) imagem com ruído senoidal, (b) espectro mostrando os picos do ruído, (c) máscara *notch* centrada nos picos, (d) imagem restaurada.

📌 Synthèse — Filtres spectraux

Filtre	Effet visuel	Artefact	Utilisation
Passe-bas idéal	Lissage intense	Ring	Illustratif
Passe-bas gaussien	Lissage doux	N’présente pas de ring	Lissage général
Passe-bas de Butterworth	Lissage contrôlé	Ring (ordres élevés)	Compromis entre lissage et sélectivité
Passe-haut	Renforcement des contours	Amplification du bruit	Détection des contours
Notch	Suppression sélective de fréquences	Distorsions locales possibles	Suppression du bruit périodique

La conception de filtres dans le domaine fréquentiel consiste à définir des masques spectraux. Cependant, des effets dans le domaine spatial, tels que le ring et le flou, émergent directement de ces choix dans le spectre.

5.6 Ondelettes et Multirésolution

La Transformée de Fourier décompose le signal en fréquences globales : chaque coefficient \(F(u,v)\) reçoit des contributions de toute l’image, sans information explicite sur la localisation spatiale de ces fréquences. Ainsi, les structures localisées, comme les bords, sont représentées de manière distribuée dans le spectre.

Les ondelettes (ondelettes) surmontent cette limitation en utilisant des fonctions de base localisées dans l’espace, qui peuvent être déplacées et mises à l’échelle. Ces fonctions possèdent un support compact, c’est-à-dire qu’elles sont différentes de zéro uniquement dans une région finie du domaine, permettant une représentation simultanée en termes de fréquence et de localisation spatiale.

5.6.1 La Limite de la Transformée de Fourier : localisation spatiale

La Transformée de Fourier décrit avec précision quelles fréquences sont présentes dans un signal, mais ne représente pas explicitement où ces fréquences se produisent dans l’espace.

Dans l’expérience présentée dans la Figure 5.16, deux images avec des structures localisées à des positions différentes produisent des spectres de magnitude pratiquement identiques. Cela est dû au fait que la représentation de Fourier est globale : chaque coefficient reçoit une contribution de l’ensemble de l’image.

Par conséquent, le spectre de magnitude ne représente pas explicitement la localisation des contours ou d’autres structures, mais uniquement la distribution des fréquences présentes. Cette limitation a motivé le développement de représentations multirésolution, telles que la Transformée Wavelet Discrète (DWT), capables de décrire simultanément la fréquence et la localisation spatiale des structures de l’image.

%%writefile tmp/fig_05_fracasso_fourier.cpp
#define MM_OUT "tmp/fig_05_fracasso_fourier.png"
//| label: fig-05-fracasso-fourier
//| fig-cap: "Fourier global é cego para a posição. Os espectros não dizem onde as bordas estão."
//| echo: true
//| output: true

#include <opencv2/opencv.hpp>
#include <vector>
#include <string>
#include "morph.hpp"
#include <filesystem>

cv::Mat fftshift(const cv::Mat& input) {
    // Desloca o espectro para o centro
    int cx = input.cols / 2;
    int cy = input.rows / 2;
    cv::Mat output = input.clone();
    // Trocar os quadrantes
    cv::Mat q0(input, cv::Rect(0, 0, cx, cy));       // superior esquerdo
    cv::Mat q1(input, cv::Rect(cx, 0, input.cols - cx, cy)); // superior direito
    cv::Mat q2(input, cv::Rect(0, cy, cx, input.rows - cy)); // inferior esquerdo
    cv::Mat q3(input, cv::Rect(cx, cy, input.cols - cx, input.rows - cy)); // inferior direito

    cv::Mat temp;
    q0.copyTo(temp); q3.copyTo(q0); temp.copyTo(q3); // trocar q0 <-> q3
    q1.copyTo(temp); q2.copyTo(q1); temp.copyTo(q2); // trocar q1 <-> q2
    return output;
}

int main() {
    // Criar os sinais
    cv::Mat img_sinal1 = cv::Mat::zeros(128, 128, CV_64F);
    img_sinal1.colRange(20, 25) = 1.0;
    img_sinal1.rowRange(100, 105) = 1.0;

    cv::Mat img_sinal2 = cv::Mat::zeros(128, 128, CV_64F);
    img_sinal2.colRange(90, 95) = 1.0;
    img_sinal2.rowRange(30, 35) = 1.0;

    // Calcular os espectros com FFT e log1p
    cv::Mat fft_mat1, mag1, fft_mat2, mag2;
    cv::dft(img_sinal1, fft_mat1, cv::DFT_COMPLEX_OUTPUT);
    cv::dft(img_sinal2, fft_mat2, cv::DFT_COMPLEX_OUTPUT);

    // Dividir canais complexos
    std::vector<cv::Mat> planes1, planes2;
    cv::split(fft_mat1, planes1);
    cv::split(fft_mat2, planes2);

    // Magnitude
    cv::Mat mag1_shift, mag2_shift;
    cv::magnitude(planes1[0], planes1[1], mag1);
    cv::magnitude(planes2[0], planes2[1], mag2);

    // Aplicar fftshift e log1p (log(1+|F|))
    mag1 = fftshift(mag1);
    mag2 = fftshift(mag2);
    cv::log(1.0 + mag1, mag1_shift);
    cv::log(1.0 + mag2, mag2_shift);

    // Normalizar para visualização
    cv::normalize(mag1_shift, mag1_shift, 0, 255, cv::NORM_MINMAX, CV_8U);
    cv::normalize(mag2_shift, mag2_shift, 0, 255, cv::NORM_MINMAX, CV_8U);

    // Converter os sinais para 8-bit para exibição
    cv::Mat img1_8u, img2_8u;
    cv::normalize(img_sinal1, img1_8u, 0, 255, cv::NORM_MINMAX, CV_8U);
    cv::normalize(img_sinal2, img2_8u, 0, 255, cv::NORM_MINMAX, CV_8U);

    // Mostrar os resultados
    std::vector<mm::Image> images = {mm::Image(img1_8u), mm::Image(mag1_shift), 
                                     mm::Image(img2_8u), mm::Image(mag2_shift)};
    std::vector<std::string> titles = {"Sinal A", "Espectro A", "Sinal B (Deslocado)", "Espectro B"};
    mm::show(images, MM_OUT, titles, 4);

    
// [pdi:panel-io] auto-generated — do not edit by hand
std::filesystem::create_directories("tmp");
mm::write(img_sinal1, "tmp/fig_05_fracasso_fourier_0.png");
mm::write(mag1, "tmp/fig_05_fracasso_fourier_1.png");
mm::write(img_sinal2, "tmp/fig_05_fracasso_fourier_2.png");
mm::write(mag2, "tmp/fig_05_fracasso_fourier_3.png");
// [pdi:panel-io:end]
return 0;
}

Overwriting tmp/fig_05_fracasso_fourier.cpp

!g++ -I. -std=c++17 -DMM_USE_OPENCV -I/usr/include/opencv4 tmp/fig_05_fracasso_fourier.cpp -o tmp/fig_05_fracasso_fourier -lopencv_stitching -lopencv_alphamat -lopencv_aruco -lopencv_barcode -lopencv_bgsegm -lopencv_bioinspired -lopencv_ccalib -lopencv_dnn_objdetect -lopencv_dnn_superres -lopencv_dpm -lopencv_face -lopencv_freetype -lopencv_fuzzy -lopencv_hdf -lopencv_hfs -lopencv_img_hash -lopencv_intensity_transform -lopencv_line_descriptor -lopencv_mcc -lopencv_quality -lopencv_rapid -lopencv_reg -lopencv_rgbd -lopencv_saliency -lopencv_shape -lopencv_stereo -lopencv_structured_light -lopencv_phase_unwrapping -lopencv_superres -lopencv_optflow -lopencv_surface_matching -lopencv_tracking -lopencv_highgui -lopencv_datasets -lopencv_text -lopencv_plot -lopencv_ml -lopencv_videostab -lopencv_videoio -lopencv_viz -lopencv_wechat_qrcode -lopencv_ximgproc -lopencv_video -lopencv_xobjdetect -lopencv_objdetect -lopencv_calib3d -lopencv_imgcodecs -lopencv_features2d -lopencv_dnn -lopencv_flann -lopencv_xphoto -lopencv_photo -lopencv_imgproc -lopencv_core \
  && ./tmp/fig_05_fracasso_fourier \
  && test -f "tmp/fig_05_fracasso_fourier.png" \
  || echo "⚠ mm::show não gravou tmp/fig_05_fracasso_fourier.png"

[1] Sinal A
[2] Espectro A
[3] Sinal B (Deslocado)
[4] Espectro B

try:
    mm.show(
        [
            mm.read("tmp/fig_05_fracasso_fourier_0.png"),
            mm.read("tmp/fig_05_fracasso_fourier_1.png"),
            mm.read("tmp/fig_05_fracasso_fourier_2.png"),
            mm.read("tmp/fig_05_fracasso_fourier_3.png"),
        ],
        titles=[
            'Sinal A',
            'Espectro A',
            'Sinal B (Deslocado)',
            'Espectro B',
        ],
        cols=4,
        figsize=(14, 4),
    )
except Exception as _e:
    print("figura indisponivel nesta trilha (C++): " + repr(_e) + " tmp/fig_05_fracasso_fourier_0.png (ver a versao Python)")

Figure 5.16: Fourier global é cego para a posição. Os espectros não dizem onde as bordas estão.

5.6.2 Transformée Wavelet Discrète 2D

La Transformée Wavelet Discrète (DWT) applique, séparément dans les directions horizontale et verticale, deux filtres complémentaires : un passe-bas \(h\) (approximation) et un passe-haut \(g\) (détails), suivis d’un sous-échantillonnage par un facteur de 2 dans chaque dimension. Ce processus produit quatre sous-bandes, dont les noms indiquent la combinaison des filtres appliqués dans chaque direction (L = Low-pass, passe-bas ; H = High-pass, passe-haut). Les caractéristiques de chaque sous-bande sont résumées dans la Table 5.3.

\[ \text{DWT}(f)=\{\underbrace{\text{LL}}_{\text{approx.}},\; \underbrace{\text{LH}}_{\text{détails horizontaux}},\; \underbrace{\text{HL}}_{\text{détails verticaux}},\; \underbrace{\text{HH}}_{\text{détails diagonaux}}\}. \]

Table 5.3: Sous-bandes produites par la Transformée Wavelet Discrète 2D (DWT), indiquant les filtres appliqués dans chaque direction et le contenu prédominant de chaque composante.

Sous-bande	Filtres appliqués	Contenu visuel
LL	bas × bas	Approximation de l’image (version lissée et réduite)
LH	bas × haut	Bords horizontaux et variations verticales
HL	haut × bas	Bords verticaux et variations horizontales
HH	haut × haut	Détails diagonaux et textures

La décomposition peut être appliquée récursivement sur la sous-bande LL, générant une représentation multi-résolution. Après \(J\) niveaux, on obtient une structure avec \(3J+1\) sous-bandes, où chaque nouveau niveau réduit la résolution de la composante d’approximation.

Lien avec les CNN

La décomposition multi-résolution des wavelets possède une relation conceptuelle avec les représentations hiérarchiques utilisées dans les réseaux de neurones convolutifs (CNN). Dans les deux cas, des étapes successives de filtrage et de réduction de résolution produisent des descriptions de plus en plus abstraites de l’image. Cependant, les wavelets utilisent des filtres mathématiquement définis et reconstructibles, tandis que les CNN apprennent leurs filtres pendant l’entraînement.

5.6.3 Familles d’ondelettes

Différentes familles d’ondelettes présentent des compromis distincts entre support spatial, régularité et capacité de compression. Le support correspond à l’étendue de la fonction ondelette dans le domaine spatial : plus le support est petit, plus la fonction est localisée ; plus il est grand, plus sa représentation tend à être régulière, mais avec un coût de calcul plus élevé. La Table 5.4 compare certaines des familles les plus utilisées.

Table 5.4: Comparaison entre familles d’ondelettes, mettant en évidence la longueur du support, le nombre de moments nuls, la symétrie et les applications typiques.

Ondelette	Longueur du support	Moments nuls	Symétrie	Usage typique
Haar	2	1	Asymétrique	Introduction et analyse de base
Daubechies db4	8	4	Asymétrique	Compression et analyse générale
Symlet sym4	8	4	Quasi symétrique	Reconstruction de signaux
Biorthogonale 5/3	5/3	2/2	Symétrique	JPEG 2000 sans perte
Biorthogonale 9/7	9/7	4/4	Symétrique	JPEG 2000 avec perte

Les moments nuls mesurent la capacité de l’ondelette à représenter les régions lisses de l’image avec peu de coefficients non nuls. Une ondelette à \(p\) moments nuls annule exactement les polynômes de degré jusqu’à \(p-1\). Par conséquent, plus le nombre de moments nuls est élevé, plus l’efficacité de compression tend à être grande dans les régions homogènes, bien que cela implique généralement des fonctions à support plus long.

La Figure 5.17 présente les fonctions de base (ondelettes) \(\psi(t)\) dans le domaine spatial. Ces fonctions possèdent un support compact, c’est-à-dire qu’elles sont non nulles uniquement sur une région finie du domaine, contrairement aux sinusoïdes de la Transformée de Fourier, qui s’étendent sur tout le domaine.

%%writefile tmp/fig_05_wavelet_functions.cpp
#define MM_OUT "tmp/fig_05_wavelet_functions.png"
#include <opencv2/opencv.hpp>
#include "morph.hpp"
#include <vector>
#include <string>
#include <filesystem>

//| label: fig-05-wavelet-functions
//| fig-cap: "Funções da *Wavelet* (ψ). Note como elas rapidamente decaem para zero (suporte compacto), ao contrário das senoides infinitas de Fourier."
//| echo: true
//| output: true

int main() {
    // psi via mm.wavefun (algoritmo em cascata) — suporte compacto: decai a zero.
    std::vector<double> x_h, phi_h, psi_h;
    mm::wavefun("haar", 6, x_h, phi_h, psi_h);
    std::vector<double> x_d, phi_d, psi_d;
    mm::wavefun("db4", 6, x_d, phi_d, psi_d);

    std::vector<std::vector<double>> xs_h = {x_h};
    std::vector<std::vector<double>> ys_h = {psi_h};
    std::vector<std::vector<double>> xs_d = {x_d};
    std::vector<std::vector<double>> ys_d = {psi_d};

    mm::Image chart_h = mm::lineChart(xs_h, ys_h, {"psi Haar"},
                                      {cv::Scalar(180, 60, 40)},
                                      "Ondaleta Haar (psi)", "t");
    mm::Image chart_d = mm::lineChart(xs_d, ys_d, {"psi Daubechies 4"},
                                      {cv::Scalar(60, 140, 40)},
                                      "Ondaleta Daubechies 4 (psi)", "t");

    mm::show(std::vector<mm::Image>{chart_h, chart_d}, MM_OUT,
             std::vector<std::string>{"Haar", "Daubechies 4"}, 2);

    
// [pdi:panel-io] auto-generated — do not edit by hand
std::filesystem::create_directories("tmp");
mm::write(chart_h, "tmp/fig_05_wavelet_functions_0.png");
mm::write(chart_d, "tmp/fig_05_wavelet_functions_1.png");
// [pdi:panel-io:end]
return 0;
}

Overwriting tmp/fig_05_wavelet_functions.cpp

!g++ -I. -std=c++17 -DMM_USE_OPENCV -I/usr/include/opencv4 tmp/fig_05_wavelet_functions.cpp -o tmp/fig_05_wavelet_functions -lopencv_stitching -lopencv_alphamat -lopencv_aruco -lopencv_barcode -lopencv_bgsegm -lopencv_bioinspired -lopencv_ccalib -lopencv_dnn_objdetect -lopencv_dnn_superres -lopencv_dpm -lopencv_face -lopencv_freetype -lopencv_fuzzy -lopencv_hdf -lopencv_hfs -lopencv_img_hash -lopencv_intensity_transform -lopencv_line_descriptor -lopencv_mcc -lopencv_quality -lopencv_rapid -lopencv_reg -lopencv_rgbd -lopencv_saliency -lopencv_shape -lopencv_stereo -lopencv_structured_light -lopencv_phase_unwrapping -lopencv_superres -lopencv_optflow -lopencv_surface_matching -lopencv_tracking -lopencv_highgui -lopencv_datasets -lopencv_text -lopencv_plot -lopencv_ml -lopencv_videostab -lopencv_videoio -lopencv_viz -lopencv_wechat_qrcode -lopencv_ximgproc -lopencv_video -lopencv_xobjdetect -lopencv_objdetect -lopencv_calib3d -lopencv_imgcodecs -lopencv_features2d -lopencv_dnn -lopencv_flann -lopencv_xphoto -lopencv_photo -lopencv_imgproc -lopencv_core \
  && ./tmp/fig_05_wavelet_functions \
  && test -f "tmp/fig_05_wavelet_functions.png" \
  || echo "⚠ mm::show não gravou tmp/fig_05_wavelet_functions.png"

[1] Haar
[2] Daubechies 4

try:
    mm.show(
        [
            mm.read("tmp/fig_05_wavelet_functions_0.png"),
            mm.read("tmp/fig_05_wavelet_functions_1.png"),
        ],
        titles=[
            'Haar',
            'Daubechies 4',
        ],
        cols=2,
    )
except Exception as _e:
    print("figura indisponivel nesta trilha (C++): " + repr(_e) + " tmp/fig_05_wavelet_functions_0.png (ver a versao Python)")

Figure 5.17: Funções da *Wavelet* (ψ). Note como elas rapidamente decaem para zero (suporte compacto), ao contrário das senoides infinitas de Fourier.

Le diagramme de la Figure 5.18 illustre l’analyse multirésolution effectuée par la DWT, dans laquelle la sous-bande d’approximation (LL) est successivement décomposée, formant une représentation hiérarchique à deux niveaux.

Décomposition en ondelettes 2D — Structure multirésolution (2 niveaux)

Figure 5.18: Schéma de la décomposition wavelet 2D en deux niveaux.

Le simulateur de la Figure 5.19 permet d’explorer de manière interactive la Transformée Wavelet Discrète 2D (TWD) à l’aide de la wavelet de Haar. La décomposition en sous-bandes met en évidence la séparation entre la composante d’approximation et les composantes de détail de l’image.

Les différents motifs d’entrée permettent d’observer le comportement directionnel des filtres. Dans les images comportant des bords horizontaux et verticaux, les sous-bandes LH et HL mettent en évidence, respectivement, les variations verticales et horizontales de l’intensité. Dans les régions à variation douce, la majeure partie de l’énergie se concentre dans la sous-bande d’approximation LL, tandis que les sous-bandes de détail présentent des coefficients proches de zéro.

L’analyse multirésolution peut également être observée en augmentant le nombre de niveaux de décomposition. Dans ce cas, seule la sous-bande \(\text{LL}_1\) est à nouveau décomposée, donnant naissance aux sous-bandes \(\text{LL}_2\), \(\text{LH}_2\), \(\text{HL}_2\) et \(\text{HH}_2\), qui forment le deuxième niveau de la représentation hiérarchique.

Dans les motifs constitués de régions homogènes de grande étendue, comme un dégradé doux ou un damier composé de grands blocs, l’énergie demeure principalement concentrée dans la sous-bande LL. Dans le dégradé, cela s’explique par le fait que les différences entre pixels voisins sont faibles. Dans le damier, en revanche, les pixels possèdent pratiquement la même intensité à l’intérieur de chaque bloc, de sorte que seules les frontières entre les blocs produisent des coefficients non nuls dans les sous-bandes de détail. Comme ces frontières n’occupent qu’une petite fraction de l’image, leur contribution à l’énergie totale reste réduite.

Pour permettre l’analyse visuelle de ces variations subtiles, le simulateur intègre un contrôle de gain de contraste des détails (variant de 1 à 8). Ce paramètre agit comme un facteur d’amplification linéaire appliqué exclusivement aux coefficients des sous-bandes de détail (LH, HL et HH) avant leur rendu à l’écran. Dans les scénarios de transition douce (comme le dégradé) ou d’uniformité locale (comme l’intérieur des blocs du damier), les différences numériques calculées par le filtre passe-haut de Haar donnent des coefficients très proches de zéro, ce qui rendrait les quadrants correspondants sombres et imperceptibles à l’œil nu. En multipliant ces valeurs par le gain, le simulateur fait ressortir visuellement les structures de haute fréquence cachées et met en évidence l’orientation des bords restants.

Le graphique de l’énergie par sous-bande quantifie cette répartition entre la composante d’approximation et les composantes de détail, démontrant que le gain visuel ne modifie pas la métrique originale de l’énergie. Dans les images naturelles, la majeure partie de l’énergie se concentre dans la sous-bande LL, tandis que les sous-bandes LH, HL et HH représentent principalement les bords, les textures et d’autres variations locales de l’intensité.

🌊 Simulateur : Décomposition Wavelet 2D Transformée de Haar

Motif Synthétique

Niveaux de Décomposition

1 niveau 2 niveaux

Gain de Contraste : 3.0

Image Originale

Décomposition Wavelet (Mosaïque)

Énergie par Sous-bande (%) — Somme Préservée (Parseval)

LL — Approximation

Version lissée et réduite de l'image

LH — Détail Horizontal

Met en évidence les bords horizontaux (variation verticale)

HL — Détail Vertical

Met en évidence les bords verticaux (variation horizontale)

HH — Détail Diagonal

Textures et coins (variation dans les deux directions)

Figure 5.19: Simulation de la décomposition wavelet 2D.

5.6.4 Analyse multirésolution avec la DWT 2D

La transformée wavelet discrète 2D (DWT) décompose une image en composantes d’approximation et de détail, organisées de manière hiérarchique à différentes échelles et orientations. Comme les sous-bandes de détail dans les images naturelles présentent souvent des coefficients de faible contraste, les exemples pratiques suivants utilisent un motif géométrique synthétique généré en Python. Cette approche reproduit le comportement du simulateur Figure 5.19, rendant visuellement explicites les effets du filtrage spatial et de la décomposition multirésolution.

5.6.4.1 Décomposition en Mosaïque à Multiples Niveaux

La Figure 5.20 illustre la structure hiérarchique de la DWT à deux niveaux en utilisant la wavelet de Haar. Le processus repose sur l’application combinée de filtres passe-bas et passe-haut dans les directions horizontale et verticale, suivis d’un sous-échantillonnage par un facteur de 2.

Au premier niveau, l’image originale donne naissance à la sous-bande d’approximation (\(LL_1\)) ainsi qu’aux composantes de détail horizontale (\(LH_1\)), verticale (\(HL_1\)) et diagonale (\(HH_1\)). Dans l’analyse multirésolution, la sous-bande \(LL_1\) est à nouveau filtrée et sous-échantillonnée, générant le deuxième niveau de décomposition (\(LL_2\), \(LH_2\), \(HL_2\) et \(HH_2\)).

Pour faciliter l’interprétation visuelle des composantes de détail, le code extrait la valeur absolue de leurs coefficients et applique une normalisation linéaire (min-max) afin d’occuper toute la plage dynamique des niveaux de gris [0, 255]. Cette opération transforme les régions homogènes (coefficients nuls) en noir et met en évidence en blanc les contours et textures extraits à chaque échelle et orientation.

%%writefile tmp/fig_05_dwt_subbandas.cpp
#define MM_OUT "tmp/fig_05_dwt_subbandas.png"
//| label: fig-05-dwt-subbandas
//| fig-cap: "Decomposição *wavelet* 2D de 2 níveis com *wavelet* Haar: subbandas LL, LH, HL, HH em cada nível. As subbandas de detalhe revelam estruturas orientadas em diferentes escalas utilizando um padrão sintético."
//| echo: true
//| output: true

#include <opencv2/opencv.hpp>
#include "morph.hpp"
#include <iostream>
#include <vector>
#include <string>
#include <cmath>

// Geração da Imagem Sintética (Mesmo padrão 'combined' do simulador)
cv::Mat gerar_imagem_sintetica(int N = 256) {
    cv::Mat img = cv::Mat::zeros(N, N, CV_64F);
    for (int y = 0; y < N; y++) {
        for (int x = 0; x < N; x++) {
            double v = 55 + 35 * ((double)x / N) + 15 * std::sin(y / 24.0);
            // Quadrado
            if (24 < x && x < 100 && 24 < y && y < 100) {
                v = 225;
            }
            // Círculo
            int cx = 190, cy = 76, r = 34;
            if ((x - cx) * (x - cx) + (y - cy) * (y - cy) < r * r) {
                v = 205;
            }
            // Textura periódica (inferior)
            if (y > 164 && y < 244) {
                int p = 12;
                v = (((x / p + y / p) % 2) == 0) ? 185 : 65;
            }
            // Linha diagonal
            if (std::abs(x - y) < 4) {
                v = 240;
            }
            img.at<double>(y, x) = std::max(0.0, std::min(v, 255.0));
        }
    }
    cv::Mat out;
    img.convertTo(out, CV_8U);
    return out;
}

// Função para visualizar subbanda
cv::Mat sb_vis(const cv::Mat& sb) {
    cv::Mat abs_sb;
    cv::absdiff(sb, cv::Scalar(0), abs_sb);
    cv::Mat normalized;
    cv::normalize(abs_sb, normalized, 0, 255, cv::NORM_MINMAX, CV_8U);
    return normalized;
}

int main() {
    // Substitui a imagem escura de moedas pelo padrão sintético claro
    cv::Mat img_gray = gerar_imagem_sintetica(256);

    // Decomposição wavelet 2 níveis
    std::string wavelet = "haar";
    cv::Mat img_float;
    img_gray.convertTo(img_float, CV_64F);

    // Nível 1
    mm::Subbands coefs1 = mm::dwt2(img_float, wavelet);
    cv::Mat LL1 = coefs1.LL;
    cv::Mat LH1 = coefs1.LH;
    cv::Mat HL1 = coefs1.HL;
    cv::Mat HH1 = coefs1.HH;

    // Nível 2 (aplicado sobre LL1)
    mm::Subbands coefs2 = mm::dwt2(LL1, wavelet);
    cv::Mat LL2 = coefs2.LL;
    cv::Mat LH2 = coefs2.LH;
    cv::Mat HL2 = coefs2.HL;
    cv::Mat HH2 = coefs2.HH;

    std::cout << "Forma original     : " << img_gray.rows << "x" << img_gray.cols << std::endl;
    std::cout << "LL1 (nível 1)      : " << LL1.rows << "x" << LL1.cols << "  |  LH1/HL1/HH1: " << LH1.rows << "x" << LH1.cols << std::endl;
    std::cout << "LL2 (nível 2)      : " << LL2.rows << "x" << LL2.cols << "    |  LH2/HL2/HH2: " << LH2.rows << "x" << LH2.cols << std::endl;

    std::vector<mm::Image> imgs_dwt = {mm::Image(img_gray), mm::Image(sb_vis(LL1)), mm::Image(sb_vis(LH1)), mm::Image(sb_vis(HL1)), mm::Image(sb_vis(HH1)),
                                       mm::Image(sb_vis(LL2)), mm::Image(sb_vis(LH2)), mm::Image(sb_vis(HL2)), mm::Image(sb_vis(HH2))};
    std::vector<std::string> titles_dwt = {"Original",
                                           "LL₁ (aprox.)", "LH₁ (horiz.)", "HL₁ (vert.)", "HH₁ (diag.)",
                                           "LL₂ (aprox.)", "LH₂ (horiz.)", "HL₂ (vert.)", "HH₂ (diag.)"};

    mm::show(imgs_dwt, MM_OUT, titles_dwt, 5);

    return 0;
}

Overwriting tmp/fig_05_dwt_subbandas.cpp

!g++ -I. -std=c++17 -DMM_USE_OPENCV -I/usr/include/opencv4 tmp/fig_05_dwt_subbandas.cpp -o tmp/fig_05_dwt_subbandas -lopencv_stitching -lopencv_alphamat -lopencv_aruco -lopencv_barcode -lopencv_bgsegm -lopencv_bioinspired -lopencv_ccalib -lopencv_dnn_objdetect -lopencv_dnn_superres -lopencv_dpm -lopencv_face -lopencv_freetype -lopencv_fuzzy -lopencv_hdf -lopencv_hfs -lopencv_img_hash -lopencv_intensity_transform -lopencv_line_descriptor -lopencv_mcc -lopencv_quality -lopencv_rapid -lopencv_reg -lopencv_rgbd -lopencv_saliency -lopencv_shape -lopencv_stereo -lopencv_structured_light -lopencv_phase_unwrapping -lopencv_superres -lopencv_optflow -lopencv_surface_matching -lopencv_tracking -lopencv_highgui -lopencv_datasets -lopencv_text -lopencv_plot -lopencv_ml -lopencv_videostab -lopencv_videoio -lopencv_viz -lopencv_wechat_qrcode -lopencv_ximgproc -lopencv_video -lopencv_xobjdetect -lopencv_objdetect -lopencv_calib3d -lopencv_imgcodecs -lopencv_features2d -lopencv_dnn -lopencv_flann -lopencv_xphoto -lopencv_photo -lopencv_imgproc -lopencv_core \
  && ./tmp/fig_05_dwt_subbandas \
  && test -f "tmp/fig_05_dwt_subbandas.png" \
  || echo "⚠ mm::show não gravou tmp/fig_05_dwt_subbandas.png"

Forma original     : 256x256
LL1 (nível 1)      : 128x128  |  LH1/HL1/HH1: 128x128
LL2 (nível 2)      : 64x64    |  LH2/HL2/HH2: 64x64
[1] Original
[2] LL₁ (aprox.)
[3] LH₁ (horiz.)
[4] HL₁ (vert.)
[5] HH₁ (diag.)
[6] LL₂ (aprox.)
[7] LH₂ (horiz.)
[8] HL₂ (vert.)
[9] HH₂ (diag.)

try:
    mm.show(mm.read("tmp/fig_05_dwt_subbandas.png"), figsize=(16, 7))
except Exception as _e:
    print("figura indisponivel nesta trilha (C++): " + repr(_e) + " tmp/fig_05_dwt_subbandas.png (ver a versao Python)")

Figure 5.20: Decomposição *wavelet* 2D de 2 níveis com *wavelet* Haar: subbandas LL, LH, HL, HH em cada nível. As subbandas de detalhe revelam estruturas orientadas em diferentes escalas utilizando um padrão sintético.

5.6.4.2 Le Compromis entre Localisation et Régularité

Le choix de la fonction de base (ondelette) influence directement la manière dont les caractéristiques de l’image sont distribuées et codées par les coefficients de la DWT. La Figure 5.21 compare les résultats pratiques obtenus en appliquant quatre familles distinctes sur le motif géométrique synthétique : haar, db4, sym4 et bior2.2.

En raison de son support court et de sa forme en fonction échelon, l’ondelette de Haar produit des coefficients hautement localisés au niveau des discontinuités spatiales, générant des bords fins et nets dans les sous-bandes de détail. En revanche, des familles comme Daubechies (db4) et Symlets (sym4), qui présentent un support plus étendu (filtres plus longs) et un plus grand nombre de moments nuls, génèrent des réponses plus lisses et plus distribuées autour des transitions, ce qui peut introduire de légères oscillations ou un adoucissement aux frontières abruptes.

Ce comportement met en évidence le compromis classique (trade-off) de l’analyse multirésolution : des supports plus petits favorisent la localisation spatiale exacte des bords, tandis que des supports plus grands et un nombre plus élevé de moments nuls tendent à produire des représentations plus parcimonieuses et plus régulières. Cette régularité et cette capacité d’atténuation des hautes fréquences garantissent une plus grande efficacité dans la compaction de l’énergie, des caractéristiques fondamentales pour les applications de compression de données et de débruitage (denoising).

%%writefile tmp/fig_05_dwt_wavelets.cpp
#define MM_OUT "tmp/fig_05_dwt_wavelets.png"
//| label: fig-05-dwt-wavelets
//| fig-cap: "Comparação entre famílias de *wavelets*: Haar, db4, sym4 e bior2.2. Subbanda LL₁ (aproximação) e HH₁ (diagonal) para cada escolha, ilustrando o compromisso entre compactação e suavidade com base no padrão sintético."
//| echo: true
//| output: true

#include <opencv2/opencv.hpp>
#include <vector>
#include <string>
#include <cmath>
#include "morph.hpp"

// Função auxiliar para visualizar subbandas (normaliza para 8 bits)
mm::Image sb_vis(const cv::Mat& subband) {
    cv::Mat normalized;
    cv::normalize(subband, normalized, 0, 255, cv::NORM_MINMAX, CV_8U);
    return mm::Image(normalized);
}

int main() {
    // Garante que img_gray e img_float utilizem o mesmo padrão sintético claro
    // (fallback: gera imagem sintética se não disponível de outra célula)
    cv::Mat img_gray;
    bool synthetic_available = false;

    // Nota: como não temos a função gerar_imagem_sintetica definida em outra célula,
    // geramos diretamente aqui.
    int N = 256;
    img_gray = cv::Mat(N, N, CV_8UC1);
    for (int y = 0; y < N; y++) {
        for (int x = 0; x < N; x++) {
            double v = 55 + 35 * (double(x) / N) + 15 * std::sin(y / 24.0);
            if (x > 24 && x < 100 && y > 24 && y < 100) v = 225;
            double cx = 190, cy = 76, r = 34;
            if (std::pow(x - cx, 2) + std::pow(y - cy, 2) < r*r) v = 205;
            if (y > 164 && y < 244) {
                int p = 12;
                v = (((x / p + y / p) % 2) == 0) ? 185 : 65;
            }
            if (std::abs(x - y) < 4) v = 240;
            v = std::max(0.0, std::min(v, 255.0));
            img_gray.at<uchar>(y, x) = static_cast<uchar>(std::round(v));
        }
    }

    cv::Mat img_float;
    img_gray.convertTo(img_float, CV_64F);

    std::vector<std::string> wavelets_comp = {"haar", "db4", "sym4", "bior2.2"};
    std::vector<mm::Image> imgs_comp;
    std::vector<std::string> titles_comp;

    for (const auto& wname : wavelets_comp) {
        mm::Subbands s = mm::dwt2(img_float, wname);
        imgs_comp.push_back(sb_vis(s.LL));
        imgs_comp.push_back(sb_vis(s.HH));
        titles_comp.push_back(wname + " — LL₁");
        titles_comp.push_back(wname + " — HH₁");
    }

    mm::show(imgs_comp, MM_OUT, titles_comp, 4);
    return 0;
}

Overwriting tmp/fig_05_dwt_wavelets.cpp

!g++ -I. -std=c++17 -DMM_USE_OPENCV -I/usr/include/opencv4 tmp/fig_05_dwt_wavelets.cpp -o tmp/fig_05_dwt_wavelets -lopencv_stitching -lopencv_alphamat -lopencv_aruco -lopencv_barcode -lopencv_bgsegm -lopencv_bioinspired -lopencv_ccalib -lopencv_dnn_objdetect -lopencv_dnn_superres -lopencv_dpm -lopencv_face -lopencv_freetype -lopencv_fuzzy -lopencv_hdf -lopencv_hfs -lopencv_img_hash -lopencv_intensity_transform -lopencv_line_descriptor -lopencv_mcc -lopencv_quality -lopencv_rapid -lopencv_reg -lopencv_rgbd -lopencv_saliency -lopencv_shape -lopencv_stereo -lopencv_structured_light -lopencv_phase_unwrapping -lopencv_superres -lopencv_optflow -lopencv_surface_matching -lopencv_tracking -lopencv_highgui -lopencv_datasets -lopencv_text -lopencv_plot -lopencv_ml -lopencv_videostab -lopencv_videoio -lopencv_viz -lopencv_wechat_qrcode -lopencv_ximgproc -lopencv_video -lopencv_xobjdetect -lopencv_objdetect -lopencv_calib3d -lopencv_imgcodecs -lopencv_features2d -lopencv_dnn -lopencv_flann -lopencv_xphoto -lopencv_photo -lopencv_imgproc -lopencv_core \
  && ./tmp/fig_05_dwt_wavelets \
  && test -f "tmp/fig_05_dwt_wavelets.png" \
  || echo "⚠ mm::show não gravou tmp/fig_05_dwt_wavelets.png"

[1] haar — LL₁
[2] haar — HH₁
[3] db4 — LL₁
[4] db4 — HH₁
[5] sym4 — LL₁
[6] sym4 — HH₁
[7] bior2.2 — LL₁
[8] bior2.2 — HH₁

try:
    mm.show(mm.read("tmp/fig_05_dwt_wavelets.png"), figsize=(14, 8))
except Exception as _e:
    print("figura indisponivel nesta trilha (C++): " + repr(_e) + " tmp/fig_05_dwt_wavelets.png (ver a versao Python)")

Figure 5.21: Comparação entre famílias de *wavelets*: Haar, db4, sym4 e bior2.2. Subbanda LL₁ (aproximação) e HH₁ (diagonal) para cada escolha, ilustrando o compromisso entre compactação e suavidade com base no padrão sintético.

5.6.4.3 Limiarização de Coeficientes e Compressão

Uma das principais aplicações da Transformada Wavelet Discreta (DWT) é a compressão de dados, impulsionada pela capacidade de representação esparsa dos coeficientes. A Figure 5.22 ilustra o efeito da limiarização abrupta (hard thresholding), técnica na qual coeficientes de detalhe com magnitude inferior a um limiar \(T\) são integralmente anulados antes do processo de síntese realizado pela Transformada Wavelet Discreta Inversa (IDWT).

À medida que o limiar \(T\) é elevado, um volume crescente de coeficientes de alta frequência é zerado. Por concentrarem menor energia, a remoção dessas componentes reduz consideravelmente a quantidade de informação necessária para representar a imagem, mantendo a componente de aproximação global (a subbanda \(LL\) mais profunda) intacta para preservar a estrutura macro. Visualmente, esse descarte de coeficientes manifesta-se através do desaparecimento progressivo de texturas finas e da suavização de transições abruptas de intensidade.

A fidelidade da imagem reconstruída frente à original é quantificada pela métrica de Pico da Relação Sinal-Ruído (PSNR, Peak Signal-to-Noise Ratio), expressa em decibéis (dB). Valores mais altos de PSNR indicam menor distorção e maior proximidade matemática com o sinal original. O experimento prático evidencia o decaimento gradual do PSNR conforme a agressividade da limiarização aumenta, permitindo avaliar numericamente o limiar ótimo para o balanço entre compressão e degradação visual.

5.6.4.4 Seuil des Coefficients et Compression

L’une des principales applications de la Transformée Wavelet Discrète (DWT) est la compression de données, portée par la capacité de représentation parcimonieuse des coefficients. La Figure 5.22 illustre l’effet du seuillage abrupt (hard thresholding), technique dans laquelle les coefficients de détail dont la magnitude est inférieure à un seuil \(T\) sont intégralement annulés avant le processus de synthèse réalisé par la Transformée Wavelet Discrète Inverse (IDWT).

À mesure que le seuil \(T\) est augmenté, un volume croissant de coefficients haute fréquence est mis à zéro. Comme ils concentrent moins d’énergie, la suppression de ces composantes réduit considérablement la quantité d’information nécessaire pour représenter l’image, tout en maintenant intacte la composante d’approximation globale (la sous-bande \(LL\) la plus profonde) afin de préserver la structure macroscopique. Visuellement, ce rejet de coefficients se manifeste par la disparition progressive des textures fines et par l’adoucissement des transitions abruptes d’intensité.

La fidélité de l’image reconstruite par rapport à l’originale est quantifiée par la métrique du Pic du Rapport Signal sur Bruit (PSNR, Peak Signal-to-Noise Ratio), exprimée en décibels (dB). Des valeurs plus élevées de PSNR indiquent une distorsion moindre et une plus grande proximité mathématique avec le signal original. L’expérience pratique met en évidence la décroissance progressive du PSNR à mesure que l’agressivité du seuillage augmente, permettant d’évaluer numériquement le seuil optimal pour l’équilibre entre compression et dégradation visuelle.

%%writefile tmp/fig_05_dwt_reconstrucao.cpp
#define MM_OUT "tmp/fig_05_dwt_reconstrucao.png"
//| label: fig-05-dwt-reconstrucao
//| fig-cap: "Reconstrução *wavelet* com limiarização de coeficientes (*hard thresholding*): à medida que o limiar aumenta, mais detalhes são zerados, produzindo imagens progressivamente mais suaves. Métrica PSNR quantifica a perda de qualidade sobre o padrão sintético."
//| echo: true
//| output: true

#include <opencv2/opencv.hpp>
#include <string>
#include <vector>
#include <cmath>
#include <algorithm>
#include "morph.hpp"

// Garante que img_gray utilize o mesmo padrão sintético claro
static cv::Mat gerar_imagem_sintetica(int N = 256) {
    cv::Mat img(N, N, CV_8UC1);
    for (int y = 0; y < N; ++y) {
        for (int x = 0; x < N; ++x) {
            double v = 55.0 + 35.0 * (static_cast<double>(x) / N) + 15.0 * std::sin(static_cast<double>(y) / 24.0);
            if (24 < x && x < 100 && 24 < y && y < 100) { v = 225.0; }
            int cx = 190, cy = 76, r = 34;
            if (std::pow(x - cx, 2) + std::pow(y - cy, 2) < r * r) { v = 205.0; }
            if (y > 164 && y < 244) {
                int p = 12;
                v = (((x / p + y / p) % 2) == 0) ? 185.0 : 65.0;
            }
            if (std::abs(x - y) < 4) { v = 240.0; }
            img.at<unsigned char>(y, x) = static_cast<unsigned char>(std::max(0.0, std::min(255.0, v)));
        }
    }
    return img;
}

static cv::Mat dwt_threshold_reconstruct(const cv::Mat& img, const std::string& wavelet = "db4", int nivel = 2, double threshold = 0.0) {
    // Decompõe, aplica limiar e reconstrói via IDWT
    cv::Mat img64f;
    img.convertTo(img64f, CV_64F);
    mm::WaveDec2 c = mm::wavedec2(img64f, wavelet, nivel);

    // Copia e aplica hard thresholding em todos os detalhes
    std::vector<std::vector<cv::Mat>> coefs_t;
    for (size_t j = 0; j < c.detail.size(); ++j) {
        std::vector<cv::Mat> sb(3);
        cv::Mat LH_t = mm::wave_threshold(c.detail[j][0], threshold, "hard");
        cv::Mat HL_t = mm::wave_threshold(c.detail[j][1], threshold, "hard");
        cv::Mat HH_t = mm::wave_threshold(c.detail[j][2], threshold, "hard");
        sb[0] = LH_t;
        sb[1] = HL_t;
        sb[2] = HH_t;
        coefs_t.push_back(sb);
    }

    // Reconstrução
    mm::WaveDec2 c_t;
    c_t.LL = c.LL;
    c_t.detail = coefs_t;
    cv::Mat rec = mm::waverec2(c_t, wavelet);

    // Recorte para dimensão original
    cv::Mat rec_cut = rec(cv::Rect(0, 0, img.cols, img.rows)).clone();
    cv::Mat rec_8u;
    cv::normalize(rec_cut, rec_8u, 0, 255, cv::NORM_MINMAX, CV_8U);
    return rec_8u;
}

int main() {
    // Cria imagem sintética de teste
    cv::Mat img_gray = gerar_imagem_sintetica(256);

    std::vector<int> thresholds = {0, 10, 30, 60, 100};
    std::vector<cv::Mat> imgs_thr;
    imgs_thr.push_back(img_gray);
    std::vector<std::string> titles_thr = {"Original"};

    for (double t : thresholds) {
        cv::Mat rec = dwt_threshold_reconstruct(img_gray, "db4", 2, t);
        double psnr = cv::PSNR(img_gray, rec);
        imgs_thr.push_back(rec);
        titles_thr.push_back("T=" + std::to_string(static_cast<int>(t)) + "  PSNR=" + 
                             std::to_string(psnr).substr(0, std::to_string(psnr).find(".") + 2) + " dB");
    }

    // Converte para mm::Image para exibição
    std::vector<mm::Image> imgs_mm;
    for (size_t i = 0; i < imgs_thr.size(); ++i) {
        imgs_mm.push_back(mm::Image(imgs_thr[i]));
    }
    mm::show(imgs_mm, MM_OUT, titles_thr, 3);

    return 0;
}

Overwriting tmp/fig_05_dwt_reconstrucao.cpp

!g++ -I. -std=c++17 -DMM_USE_OPENCV -I/usr/include/opencv4 tmp/fig_05_dwt_reconstrucao.cpp -o tmp/fig_05_dwt_reconstrucao -lopencv_stitching -lopencv_alphamat -lopencv_aruco -lopencv_barcode -lopencv_bgsegm -lopencv_bioinspired -lopencv_ccalib -lopencv_dnn_objdetect -lopencv_dnn_superres -lopencv_dpm -lopencv_face -lopencv_freetype -lopencv_fuzzy -lopencv_hdf -lopencv_hfs -lopencv_img_hash -lopencv_intensity_transform -lopencv_line_descriptor -lopencv_mcc -lopencv_quality -lopencv_rapid -lopencv_reg -lopencv_rgbd -lopencv_saliency -lopencv_shape -lopencv_stereo -lopencv_structured_light -lopencv_phase_unwrapping -lopencv_superres -lopencv_optflow -lopencv_surface_matching -lopencv_tracking -lopencv_highgui -lopencv_datasets -lopencv_text -lopencv_plot -lopencv_ml -lopencv_videostab -lopencv_videoio -lopencv_viz -lopencv_wechat_qrcode -lopencv_ximgproc -lopencv_video -lopencv_xobjdetect -lopencv_objdetect -lopencv_calib3d -lopencv_imgcodecs -lopencv_features2d -lopencv_dnn -lopencv_flann -lopencv_xphoto -lopencv_photo -lopencv_imgproc -lopencv_core \
  && ./tmp/fig_05_dwt_reconstrucao \
  && test -f "tmp/fig_05_dwt_reconstrucao.png" \
  || echo "⚠ mm::show não gravou tmp/fig_05_dwt_reconstrucao.png"

[1] Original
[2] T=0  PSNR=19.5 dB
[3] T=10  PSNR=19.6 dB
[4] T=30  PSNR=22.9 dB
[5] T=60  PSNR=21.5 dB
[6] T=100  PSNR=18.5 dB

try:
    mm.show(mm.read("tmp/fig_05_dwt_reconstrucao.png"), figsize=(14, 10))
except Exception as _e:
    print("figura indisponivel nesta trilha (C++): " + repr(_e) + " tmp/fig_05_dwt_reconstrucao.png (ver a versao Python)")

Figure 5.22: Reconstrução *wavelet* com limiarização de coeficientes (*hard thresholding*): à medida que o limiar aumenta, mais detalhes são zerados, produzindo imagens progressivamente mais suaves. Métrica PSNR quantifica a perda de qualidade sobre o padrão sintético.

Synthèse — Fourier vs. wavelets : quand utiliser chaque approche ?

La Table 5.5 synthétise les principales différences structurelles et opérationnelles entre la transformée de Fourier discrète (TFD) et la transformée en ondelettes discrète (TOD).

Table 5.5: Comparaison entre la transformée de Fourier discrète (TFD) et la transformée en ondelettes discrète (TOD), mettant en évidence leurs principales caractéristiques et applications.

Critère	Fourier (TFD)	Ondelettes (TOD)
Fonctions de base	Sinusoïdes à support infini	Fonctions à support compact
Localisation spatiale	Non explicite (globale)	Explicite (locale)
Filtrage spectral	Excellent pour un contrôle fin des fréquences	Basé sur des sous-bandes (échelles)
Compression d’images	Base de la TCD (JPEG traditionnel)	Base de la TOD (JPEG 2000)
Analyse multi-échelle	Non	Oui
Suppression du bruit périodique	Très efficace	Peu indiquée
Signaux non stationnaires	Limitée	Très efficace

En termes pratiques, la TFD s’impose comme l’outil idéal pour l’analyse spectrale pure, la conception de filtres sélectifs dans le domaine fréquentiel et l’atténuation des bruits périodiques et harmoniques. En revanche, la TOD excelle dans les scénarios exigeant une préservation rigoureuse de la localisation spatiale des caractéristiques associée à leur contenu fréquentiel, se distinguant dans la compression de données, l’analyse multirésolution et le traitement des transitions abruptes. Ainsi, les deux transformées doivent être comprises comme des techniques parfaitement complémentaires, traçant des voies distinctes et spécifiques pour la résolution de problèmes en TNI-VC.

Analogies avec l’audio : limites et précautions

Lorsqu’on établit des analogies entre le traitement d’images et l’audio, il est important de prendre en compte les différences fondamentales :

Dans les systèmes audio stéréo/multicanaux, la phase entre les canaux est cruciale pour la perception de la localisation spatiale (différences interaurales de phase et de temps).
Dans les systèmes monauraux, la phase a une influence perceptuelle limitée — l’oreille humaine est relativement insensible à la phase absolue des composantes sinusoïdales isolées.
Dans les images, la phase de la TFD est toujours fondamentale pour la localisation spatiale des structures, qu’il s’agisse d’une image monochrome ou couleur.

L’analogie entre la phase en audio et la phase en images doit être utilisée avec prudence, en soulignant que, bien que toutes deux portent des informations sur l’organisation spatiale/temporelle du signal, les mécanismes perceptuels sont fondamentalement différents.

5.7 Compression d’images

Alors que les wavelets établissent la fondation théorique du standard JPEG 2000, le standard JPEG traditionnel repose sur la Transformée en Cosinus Discrète (DCT, Discrete Cosine Transform). Malgré les différences structurelles, les deux approches partagent le même principe fondamental : compacter l’énergie de l’image en un nombre réduit de coefficients et éliminer les composantes de moindre importance avec un impact visuel minimal.

L’objectif central de la compression est de réduire le volume de données nécessaire au stockage ou à la transmission d’une image. Ce processus est rendu possible par l’identification et l’élimination des redondances structurelles et perceptuelles.

5.7.1 Taxonomie des Redondances

Le développement des algorithmes de compression repose sur l’identification et l’élimination de trois catégories principales de redondance, synthétisées dans la Table 5.6.

Table 5.6: Catégories de redondance dans les images numériques et leurs mécanismes d’exploitation respectifs.

Type	Définition	Approche d’exploitation
Spatiale (interpixel)	Forte corrélation et dépendance statistique entre pixels voisins.	DCT, DWT et codage prédictif.
Spectrale (intercanal)	Corrélation statistique entre les canaux de couleur d’une même image.	Transformations d’espace colorimétrique (ex : RGB vers \(YC_bC_r\)).
Psychovisuelle	Insensibilité du système visuel humain (SVH) aux variations de haute fréquence et de faible contraste.	Processus de quantification sélective des coefficients.

Selon la préservation de l’information originale après le processus de décodage, les méthodes de compression se divisent en deux classes fondamentales :

Sans perte (lossless) : Garantit une reconstruction bit à bit identique à l’image originale. Elle est employée dans des scénarios où l’intégrité des données est strictement critique, comme dans l’imagerie médicale, les diagnostics par imagerie et le stockage de documents textuels.
Avec perte (lossy) : Admet l’introduction d’une distorsion contrôlée du signal en échange de taux de compression substantiellement plus élevés. C’est l’approche standard pour les photographies grand public et le streaming vidéo, écosystèmes dans lesquels le SVH tolère de légères atténuations haute fréquence sans perception de dégradation de la qualité visuelle.

5.7.2 Transformée en Cosinus Discrète (DCT-II 2D)

La Transformée en Cosinus Discrète (DCT) constitue l’opération centrale du standard JPEG. Contrairement à la TFD, qui utilise une base complexe, la DCT repose sur des fonctions trigonométriques purement réelles. Pour un bloc d’image \(f(x,y)\) de dimensions \(N \times N\), la DCT-II 2D mappe le signal spatial vers le domaine des fréquences spatiales, générant la matrice de coefficients \(C(u,v)\) au moyen de :

\[ C(u,v) = \alpha(u)\,\alpha(v) \sum_{x=0}^{N-1}\sum_{y=0}^{N-1} f(x,y)\, \cos\!\left[\frac{\pi(2x+1)u}{2N}\right] \cos\!\left[\frac{\pi(2y+1)v}{2N}\right] \tag{5.9}\]

où les facteurs de normalisation orthogonale sont donnés par \(\alpha(0) = \sqrt{1/N}\) et \(\alpha(k) = \sqrt{2/N}\) pour \(k > 0\).

Chaque coefficient \(C(u,v)\) quantifie la contribution — ou le « poids » — d’une fréquence spatiale spécifique au sein de ce bloc. Le terme \(C(0,0)\) est appelé composante DC et représente l’intensité moyenne du bloc (fréquence nulle). Les autres coefficients, désignés sous le nom de composantes AC (Alternating Current), correspondent aux fréquences spatiales progressivement plus élevées.

5.7.3 Les fonctions de base de la DCT

D’un point de vue géométrique, la Équation 5.9 réalise la projection du bloc de pixels sur un ensemble de fonctions orthogonales. Pour le cas standard du JPEG (\(N=8\)), le bloc spatial est décomposé en une combinaison linéaire de 64 fonctions de base bidimensionnelles, notées \(B_{u,v}(x,y)\) et générées par le produit de fonctions cosinusoïdales :

\[B_{u,v}(x,y) = \cos\left[ \frac{\pi (2x+1)u}{16} \right] \cos\left[ \frac{\pi (2y+1)v}{16} \right]\]

Ainsi, l’opération inverse peut être interprétée comme la reconstruction exacte du bloc original par la somme pondérée de ces 64 matrices de base, où chaque coefficient \(C(u,v)\) agit comme le poids analytique de sa composante harmonique respective.

La fréquence spatiale indiquée par les indices \((u,v)\) détermine le nombre de cycles d’oscillation le long des dimensions horizontales et verticales du bloc. Comme illustré dans la Figure 5.23 — dont le code isole chaque base en appliquant la transformation inverse sur des impulsions unitaires —, ces 64 fonctions sont organisées en une matrice \(8 \times 8\). Le coin supérieur gauche (\(u=0, v=0\)) présente le motif uniforme de fréquence nulle (DC), tandis que la progression vers la droite (axe \(u\)) ou vers le bas (axe \(v\)) représente des variations harmoniques progressivement plus importantes, traduisant des transitions rapides, des contours et des textures dans les orientations horizontales, verticales et diagonales.

DCT vs DFT : avantage de la compaction d’énergie

La DCT et la DFT mappent toutes deux un bloc spatial \(N \times N\) en une matrice de coefficients de même dimension. Cependant, pour les images naturelles, la DCT présente une plus grande efficacité en matière de compaction d’énergie dans les basses fréquences. Cela s’explique par le fait que la DCT suppose implicitement une symétrie paire du signal aux frontières du bloc, ce qui équivaut à une extension périodique continue, minimisant ainsi l’effet de diffusion spectrale (ringing). Par conséquent, la plupart des coefficients AC décroissent rapidement vers des valeurs proches de zéro, optimisant le pipeline de compression sans introduire de dégradation visuelle perceptible.

%%writefile tmp/fig_05_dct_basis.cpp
#define MM_OUT "tmp/fig_05_dct_basis.png"
//| label: fig-05-dct-basis
//| fig-cap: "O Alfabeto Visual do JPEG: As 64 funções de base da DCT-II. O coeficiente DC fica no topo esquerdo (suave). Ao descer e avançar à direita, a oscilação espacial aumenta drasticamente."
//| echo: true
//| output: true

#include <opencv2/opencv.hpp>
#include <vector>
#include <string>
#include "morph.hpp"
#include <filesystem>

int main() {
    // Cada base é a IDCT de um único coeficiente unitário — montadas num mosaico 8×8.
    int tile = 32;
    cv::Mat montagem = cv::Mat::zeros(8 * tile, 8 * tile, CV_8U);
    for (int i = 0; i < 8; i++) {
        for (int j = 0; j < 8; j++) {
            cv::Mat coef = cv::Mat::zeros(8, 8, CV_64F);
            coef.at<double>(i, j) = 1.0;
            cv::Mat base = mm::idct2(coef);
            cv::Mat base_normalized;
            cv::normalize(base, base_normalized, 0, 255, cv::NORM_MINMAX, CV_8U);
            cv::Mat base_resized;
            cv::resize(base_normalized, base_resized, cv::Size(tile, tile), 0, 0, cv::INTER_NEAREST);
            base_resized.copyTo(montagem(cv::Rect(j * tile, i * tile, tile, tile)));
        }
    }

    std::string titulo = "As 64 bases da DCT-II 8×8 (DC no topo-esquerdo)";
    mm::show(std::vector<mm::Image>{mm::Image(montagem)}, MM_OUT, {titulo}, 1);

    
// [pdi:panel-io] auto-generated — do not edit by hand
std::filesystem::create_directories("tmp");
mm::write(montagem, "tmp/fig_05_dct_basis_0.png");
// [pdi:panel-io:end]
return 0;
}

Overwriting tmp/fig_05_dct_basis.cpp

!g++ -I. -std=c++17 -DMM_USE_OPENCV -I/usr/include/opencv4 tmp/fig_05_dct_basis.cpp -o tmp/fig_05_dct_basis -lopencv_stitching -lopencv_alphamat -lopencv_aruco -lopencv_barcode -lopencv_bgsegm -lopencv_bioinspired -lopencv_ccalib -lopencv_dnn_objdetect -lopencv_dnn_superres -lopencv_dpm -lopencv_face -lopencv_freetype -lopencv_fuzzy -lopencv_hdf -lopencv_hfs -lopencv_img_hash -lopencv_intensity_transform -lopencv_line_descriptor -lopencv_mcc -lopencv_quality -lopencv_rapid -lopencv_reg -lopencv_rgbd -lopencv_saliency -lopencv_shape -lopencv_stereo -lopencv_structured_light -lopencv_phase_unwrapping -lopencv_superres -lopencv_optflow -lopencv_surface_matching -lopencv_tracking -lopencv_highgui -lopencv_datasets -lopencv_text -lopencv_plot -lopencv_ml -lopencv_videostab -lopencv_videoio -lopencv_viz -lopencv_wechat_qrcode -lopencv_ximgproc -lopencv_video -lopencv_xobjdetect -lopencv_objdetect -lopencv_calib3d -lopencv_imgcodecs -lopencv_features2d -lopencv_dnn -lopencv_flann -lopencv_xphoto -lopencv_photo -lopencv_imgproc -lopencv_core \
  && ./tmp/fig_05_dct_basis \
  && test -f "tmp/fig_05_dct_basis.png" \
  || echo "⚠ mm::show não gravou tmp/fig_05_dct_basis.png"

[1] As 64 bases da DCT-II 8×8 (DC no topo-esquerdo)

try:
    mm.show(
        [
            mm.read("tmp/fig_05_dct_basis_0.png"),
        ],
        titles=[
            'As 64 bases da DCT-II 8×8 (DC no topo-esquerdo)',
        ],
        cols=1,
    )
except Exception as _e:
    print("figura indisponivel nesta trilha (C++): " + repr(_e) + " tmp/fig_05_dct_basis_0.png (ver a versao Python)")

Figure 5.23: O Alfabeto Visual do JPEG: As 64 funções de base da DCT-II. O coeficiente DC fica no topo esquerdo (suave). Ao descer e avançar à direita, a oscilação espacial aumenta drasticamente.

5.7.4 Concentration d’Énergie et Reconstruction Progressive

Avant l’application de la DCT, les pixels du bloc d’intensité sont systématiquement translatés (en soustrayant \(128\) pour les images 8 bits) afin de centrer le signal autour de zéro, éliminant ainsi les composantes continues superflues. Lors du calcul de la DCT sur le bloc résultant, la propriété de compaction d’énergie devient évidente : la quasi-totalité de la variance et de l’information de l’image originale se concentre dans le coefficient DC (\(C(0,0)\)) et dans les premiers harmoniques AC de basse fréquence.

La Figure 5.24 illustre ce phénomène par une reconstruction progressive par troncature abrupte. Au lieu d’utiliser l’ensemble des 64 coefficients, l’algorithme ne conserve que les \(k\) premiers composants — sélectionnés sur la base d’un balayage qui privilégie les basses fréquences spatiales — et annule les autres.

La synthèse inverse (IDCT) réalisée avec seulement une fraction des coefficients (comme 15 % ou 30 %) est déjà capable de récupérer les structures et l’éclairage macro du bloc de pixels original. À mesure que les harmoniques de fréquences plus élevées sont progressivement réincorporés, les détails fins et les transitions rapides sont restaurés. Ce comportement valide le principe de la compression perceptuelle : les hautes fréquences écartées possèdent peu d’énergie et leur absence, en conditions normales, génère un impact visuel secondaire sur la perception de l’observateur.

%%writefile tmp/fig_05_dct_bloco.cpp
#define MM_OUT "tmp/fig_05_dct_bloco.png"
//| label: fig-05-dct-bloco
//| fig-cap: "DCT 2D en bloc 8×8 : coefficients et reconstruction progressive."
//| echo: true
//| output: true

#include <opencv2/opencv.hpp>
#include <vector>
#include <string>
#include <algorithm>
#include <cmath>
#include <iostream>
#include <iomanip>
#include "morph.hpp"

int main() {
// [pdi:state-io] auto-generated — do not edit by hand
mm::Image img_gray = mm::_read_state("tmp/state/img_gray_20.png");
// [pdi:state-io:end]

    // img_gray est fourni automatiquement (mm::Image)

    // ── Bloc 8×8 centré de l'image ─────────────────────────────────────────
    int cy = img_gray.h / 2;
    int cx = img_gray.w / 2;

    // Extraction du bloc 8×8 et conversion en CV_64F avec soustraction de 128
    cv::Mat img_mat = img_gray;
    cv::Mat bloco_f;
    img_mat(cv::Rect(cx, cy, 8, 8)).convertTo(bloco_f, CV_64F);
    bloco_f -= 128.0;

    // DCT 2D du bloc
    cv::Mat C = mm::dct2(bloco_f);

    std::cout << "Coefficients DCT du bloc 8×8:" << std::endl;
    for(int i = 0; i < 8; i++) {
        for(int j = 0; j < 8; j++) {
            std::cout << std::setw(5) << std::round(C.at<double>(i,j)) << " ";
        }
        std::cout << std::endl;
    }

    // Calcul des énergies
    double energia_dc = C.at<double>(0,0) * C.at<double>(0,0);
    double energia_total = 0.0;
    for(int i = 0; i < 8; i++) {
        for(int j = 0; j < 8; j++) {
            energia_total += C.at<double>(i,j) * C.at<double>(i,j);
        }
    }

    std::cout << "\nÉnergie DC     : " << std::fixed << std::setprecision(1) << energia_dc << std::endl;
    std::cout << "Énergie totale : " << energia_total << std::endl;
    std::cout << "Fraction en DC : " << std::fixed << std::setprecision(0) 
              << (energia_dc / energia_total * 100.0) << "% ← concentration d'énergie" << std::endl;

    // ── Reconstruction progressive ──────────────────────────────────────────
    cv::Mat bloco_orig;
    bloco_f.copyTo(bloco_orig);
    bloco_orig += 128.0;
    cv::Mat bloco_orig_8u;
    bloco_orig.convertTo(bloco_orig_8u, CV_8U);

    std::vector<mm::Image> imgs_rec;
    imgs_rec.push_back(mm::Image(bloco_orig_8u));

    std::vector<std::string> titles_rec;
    titles_rec.push_back("Bloc original\n(8×8 pixels)");

    // Liste des paires (u,v) triées par somme u+v
    std::vector<std::pair<int,int>> indices;
    for(int u = 0; u < 8; u++) {
        for(int v = 0; v < 8; v++) {
            indices.push_back({u, v});
        }
    }
    std::sort(indices.begin(), indices.end(), 
              [](const std::pair<int,int>& a, const std::pair<int,int>& b) {
                  return (a.first + a.second) < (b.first + b.second);
              });

    std::vector<int> keeps = {1, 4, 10, 20, 40, 64};

    for(int keep : keeps) {
        // Création de la matrice tronquée
        cv::Mat C_trunc = cv::Mat::zeros(8, 8, CV_64F);
        for(int k = 0; k < keep; k++) {
            int u = indices[k].first;
            int v = indices[k].second;
            C_trunc.at<double>(u,v) = C.at<double>(u,v);
        }

        // Reconstruction IDCT
        cv::Mat rec_f = mm::idct2(C_trunc);
        rec_f += 128.0;

        cv::Mat rec_8u;
        rec_f.convertTo(rec_8u, CV_8U);
        cv::Mat rec_clipped;
        cv::threshold(rec_8u, rec_clipped, 255, 255, cv::THRESH_TRUNC);
        cv::threshold(rec_clipped, rec_clipped, 0, 0, cv::THRESH_TOZERO);

        imgs_rec.push_back(mm::Image(rec_clipped));

        double pct = (keep / 64.0) * 100.0;
        titles_rec.push_back(std::to_string(keep) + " coef.\n(" + 
                           std::to_string((int)pct) + "% du total)");
    }

    mm::show(imgs_rec, MM_OUT, titles_rec, 4);

    return 0;
}

Overwriting tmp/fig_05_dct_bloco.cpp

!g++ -I. -std=c++17 -DMM_USE_OPENCV -I/usr/include/opencv4 tmp/fig_05_dct_bloco.cpp -o tmp/fig_05_dct_bloco -lopencv_stitching -lopencv_alphamat -lopencv_aruco -lopencv_barcode -lopencv_bgsegm -lopencv_bioinspired -lopencv_ccalib -lopencv_dnn_objdetect -lopencv_dnn_superres -lopencv_dpm -lopencv_face -lopencv_freetype -lopencv_fuzzy -lopencv_hdf -lopencv_hfs -lopencv_img_hash -lopencv_intensity_transform -lopencv_line_descriptor -lopencv_mcc -lopencv_quality -lopencv_rapid -lopencv_reg -lopencv_rgbd -lopencv_saliency -lopencv_shape -lopencv_stereo -lopencv_structured_light -lopencv_phase_unwrapping -lopencv_superres -lopencv_optflow -lopencv_surface_matching -lopencv_tracking -lopencv_highgui -lopencv_datasets -lopencv_text -lopencv_plot -lopencv_ml -lopencv_videostab -lopencv_videoio -lopencv_viz -lopencv_wechat_qrcode -lopencv_ximgproc -lopencv_video -lopencv_xobjdetect -lopencv_objdetect -lopencv_calib3d -lopencv_imgcodecs -lopencv_features2d -lopencv_dnn -lopencv_flann -lopencv_xphoto -lopencv_photo -lopencv_imgproc -lopencv_core \
  && ./tmp/fig_05_dct_bloco \
  && test -f "tmp/fig_05_dct_bloco.png" \
  || echo "⚠ mm::show não gravou tmp/fig_05_dct_bloco.png"

Coefficients DCT du bloc 8×8:
  192   -48     1     8     0    -1    -0    -0 
  -96    54     7   -10    -0     0     0     0 
   14   -14     8    -0     0    -0    -0     0 
   -0     0   -11     0     0     0    -0     0 
    9   -11    -0     0     1    -0     0    -0 
   -0     0    -0    -0    -0    -0    -0     0 
   -0     0     0    -0    -1     0     0    -0 
    0    -0    -0    -0     0     0     0     0 

Énergie DC     : 36768.1
Énergie totale : 52234.0
Fraction en DC : 70% ← concentration d'énergie
[1] Bloc original
(8×8 pixels)
[2] 1 coef.
(1% du total)
[3] 4 coef.
(6% du total)
[4] 10 coef.
(15% du total)
[5] 20 coef.
(31% du total)
[6] 40 coef.
(62% du total)
[7] 64 coef.
(100% du total)

try:
    mm.show(mm.read("tmp/fig_05_dct_bloco.png"), figsize=(12, 7))
except Exception as _e:
    print("figura indisponivel nesta trilha (C++): " + repr(_e) + " tmp/fig_05_dct_bloco.png (ver a versao Python)")

Figure 5.24: DCT 2D em bloco 8×8: coeficientes e reconstrução progressiva.

5.7.5 Le Pipeline de Compression JPEG

La norme JPEG opère en divisant l’image en blocs disjoints de \(8 \times 8\) pixels, traités par une séquence de transformations spatiales, perceptuelles et statistiques. Le pipeline complet de codage est structuré en six étapes principales :

\[ \text{RGB} \xrightarrow{\text{(1) } YC_bC_r} \xrightarrow{\text{(2) Sous-échantillonnage}} \xrightarrow{\text{(3) Blocs } 8 \times 8} \xrightarrow{\text{(4) DCT}} \xrightarrow{\text{(5) Quantification}} \xrightarrow{\text{(6) Codage entropique}} \]

La Table 5.7 détaille la fonction analytique et le fondement perceptuel qui justifient chacune de ces étapes.

Table 5.7: Étapes du pipeline de compression JPEG et leurs fondements de conception respectifs.

Étape	Opération	Fondement perceptuel et statistique
1	Conversion \(RGB \rightarrow YC_bC_r\)	Sépare la luminance (\(Y\)) de la chrominance (\(C_b, C_r\)). Le système visuel humain (SVH) présente une plus grande sensibilité aux variations de luminosité qu’aux variations de couleur.
2	Sous-échantillonnage de la chrominance (ex. : 4:2:0)	Réduit la résolution spatiale des canaux de couleur de moitié, en éliminant des données redondantes avec un impact visuel négligeable.
3–4	Centrage et application de la DCT \(8 \times 8\)	Translates les pixels dans l’intervalle \([-128, 127]\) et compresse l’énergie spectrale du bloc dans les coefficients de basse fréquence.
5	Quantification linéaire sélective	Divise chaque coefficient \(C(u,v)\) par l’élément correspondant de la matrice \(Q(u,v)\), avec arrondi entier. Constitue la principale source de compression avec perte.
6	Balayage en zigzag et codage	Ordonne les coefficients quantifiés pour maximiser les séquences nulles consécutives, optimisant le codage par longueur de plage (RLE) et le codage de Huffman.

La matrice de quantification \(Q(u,v)\) est le mécanisme central de contrôle du compromis entre taux de compression et qualité visuelle. Dans l’algorithme pratique de la Figure 5.25, le facteur de qualité spécifié par l’utilisateur (échelle de 1 à 100) est converti en un scalaire qui paramètre la sévérité de la matrice \(Q\). Des valeurs de qualité réduites élargissent les diviseurs de \(Q(u,v)\), forçant le tronquement en masse des coefficients AC à zéro. Lorsque cette élimination est excessive, la discontinuité aux frontières des blocs adjacents n’est pas atténuée lors de la reconstruction, générant ce que l’on appelle les artefacts de bloc (blocking artifacts).

La Logique du Balayage en Zigzag

L’efficacité du codeur entropique subséquent à la quantification dépend directement de l’ordonnancement des données. Comme la DCT concentre l’énergie vitale dans le coin supérieur gauche de la matrice (basses fréquences) et repousse les coefficients nuls vers les extrémités opposées, la lecture linéaire par lignes ou par colonnes fragmenterait les séquences de zéros.

L’ordonnancement en zigzag résout cette limitation en parcourant la matrice en diagonale, dans l’ordre croissant de la fréquence spatiale. Ce mappage regroupe les coefficients significatifs au début du vecteur et concentre les coefficients nuls en une seule séquence continue à la fin du tableau, permettant à l’algorithme RLE de coder de grands blocs de données de manière compacte et efficace.

Qu’est-ce que le RLE ?

RLE (Run-Length Encoding) est une technique de compression sans perte qui code des séquences consécutives de valeurs identiques — en particulier des zéros — comme une paire (compteur, valeur). Dans le JPEG, après le balayage en zigzag, les coefficients quantifiés sont organisés de sorte que les zéros se concentrent à la fin du vecteur. Le RLE compresse ensuite cette longue course de zéros avec une extrême efficacité, optimisant le stockage et la transmission de l’image compressée.

%%writefile tmp/fig_05_jpeg_pipeline.cpp
#define MM_OUT "tmp/fig_05_jpeg_pipeline.png"
#include <opencv2/opencv.hpp>
#include "morph.hpp"
#include <string>
#include <vector>

int main() {
    //| label: fig-05-jpeg-pipeline
    //| fig-cap: "*Pipeline* JPEG simplifié appliqué à l'image classique du *Cameraman* : DCT en blocs 8×8, quantification avec différents facteurs de qualité et reconstruction via IDCT. Les artefacts de bloc (*blocking artifacts*) deviennent visuellement évidents pour des facteurs de qualité réduits ($Q=10$ et $Q=25$)."
    //| echo: true
    //| output: true

    // Pipeline JPEG (DCT 8x8 -> quantification -> IDCT) via mm.jpegCompress,
    // sur l'image classique du Cameraman (asset du chapitre) réduite à 256x256.
    mm::Image img_orig = mm::read("imagens/cameraman.png");
    cv::Mat img_resized;
    cv::resize(cv::Mat(img_orig), img_resized, cv::Size(256, 256));
    mm::Image img_src = mm::gray(mm::Image(img_resized));

    std::vector<mm::Image> imgs;
    std::vector<std::string> titles;
    imgs.push_back(img_src);
    titles.push_back("Original (Cameraman)");

    std::vector<int> qualities = {10, 25, 50, 75, 90};
    for (int q : qualities) {
        mm::Image rec = mm::jpegCompress(img_src, q);
        imgs.push_back(rec);
        double psnr_val = mm::psnr(img_src, rec);
        titles.push_back("Q=" + std::to_string(q) + " (PSNR=" + 
                         std::to_string(psnr_val).substr(0, std::to_string(psnr_val).find(".") + 2) + " dB)");
    }

    mm::show(imgs, MM_OUT, titles, 3);

    return 0;
}

Overwriting tmp/fig_05_jpeg_pipeline.cpp

!g++ -I. -std=c++17 -DMM_USE_OPENCV -I/usr/include/opencv4 tmp/fig_05_jpeg_pipeline.cpp -o tmp/fig_05_jpeg_pipeline -lopencv_stitching -lopencv_alphamat -lopencv_aruco -lopencv_barcode -lopencv_bgsegm -lopencv_bioinspired -lopencv_ccalib -lopencv_dnn_objdetect -lopencv_dnn_superres -lopencv_dpm -lopencv_face -lopencv_freetype -lopencv_fuzzy -lopencv_hdf -lopencv_hfs -lopencv_img_hash -lopencv_intensity_transform -lopencv_line_descriptor -lopencv_mcc -lopencv_quality -lopencv_rapid -lopencv_reg -lopencv_rgbd -lopencv_saliency -lopencv_shape -lopencv_stereo -lopencv_structured_light -lopencv_phase_unwrapping -lopencv_superres -lopencv_optflow -lopencv_surface_matching -lopencv_tracking -lopencv_highgui -lopencv_datasets -lopencv_text -lopencv_plot -lopencv_ml -lopencv_videostab -lopencv_videoio -lopencv_viz -lopencv_wechat_qrcode -lopencv_ximgproc -lopencv_video -lopencv_xobjdetect -lopencv_objdetect -lopencv_calib3d -lopencv_imgcodecs -lopencv_features2d -lopencv_dnn -lopencv_flann -lopencv_xphoto -lopencv_photo -lopencv_imgproc -lopencv_core \
  && ./tmp/fig_05_jpeg_pipeline \
  && test -f "tmp/fig_05_jpeg_pipeline.png" \
  || echo "⚠ mm::show não gravou tmp/fig_05_jpeg_pipeline.png"

[1] Original (Cameraman)
[2] Q=10 (PSNR=28.0 dB)
[3] Q=25 (PSNR=30.6 dB)
[4] Q=50 (PSNR=32.8 dB)
[5] Q=75 (PSNR=35.1 dB)
[6] Q=90 (PSNR=40.0 dB)

try:
    mm.show(mm.read("tmp/fig_05_jpeg_pipeline.png"), figsize=(14, 10))
except Exception as _e:
    print("figura indisponivel nesta trilha (C++): " + repr(_e) + " tmp/fig_05_jpeg_pipeline.png (ver a versao Python)")

Figure 5.25: *Pipeline* JPEG simplificado aplicado à imagem clássica do *Cameraman*: DCT em blocos 8×8, quantização com diferentes fatores de qualidade e reconstrução via IDCT. Os artefatos de bloco (*blocking artifacts*) tornam-se visualmente evidentes em fatores de qualidade reduzidos (\(Q=10\) e \(Q=25\)).

5.7.6 Simulateur interactif : Quantification DCT

Le simulateur de la Figure 5.26 permet d’explorer l’impact du processus de quantification sur un bloc \(8 \times 8\) extrait d’une image réelle, en synthétisant en temps réel les composantes suivantes :

Bloc original et reconstruit : Représentation directe des pixels dans le domaine spatial en niveaux de gris [0, 255].
Coefficients DCT : Distribution de l’énergie mappée de manière logarithmique sur un dégradé chromatique, mettant en évidence la concentration de l’intensité dans le coin supérieur gauche (basses fréquences).
Coefficients quantifiés : Affichage des valeurs entières résultant de la division par la matrice \(Q(u,v)\), rendant visuellement explicite l’apparition massive de coefficients nuls (en tons sombres) à mesure que le facteur de qualité est réduit.
Métriques de compression : Panneau de surveillance qui quantifie l’Erreur Quadratique Moyenne (MSE), le nombre de coefficients préservés et le volume de zéros générés pour le codage entropique.

⊞ Simulateur : Quantification DCT-JPEG (bloc 8×8) blocs 8×8

Qualité

Coef. ≠ 0

–

Zéros

–

Erreur MSE

–

Bloc original (8×8)

Coef. DCT (abs, log)

Coef. quantifiés

Bloc reconstruit

Qualité JPEG : 50

Figure 5.26: Simulateur interactif de compression DCT-JPEG : ajustez le facteur de qualité et visualisez en temps réel les coefficients nuls, le bloc reconstruit et l’erreur de quantification.

5.8 Comparação de Formatos de Imagem

Le choix d’un format de stockage numérique impacte directement le compromis entre qualité visuelle, taille de fichier et coût computationnel de décodage. Les trois formats les plus pertinents pour les architectures web et les systèmes de calcul visuel sont le JPEG, le PNG et le WebP.

5.8.1 Caractéristiques des Formats

La Table 5.8 synthétise les propriétés structurelles des principaux formats d’image matriciels.

Table 5.8: Comparaison structurelle entre les principaux formats d’image matriciels.

Caractéristique	JPEG	PNG	WebP
Compression	Avec perte	Sans perte	Avec et sans perte.
Transparence (canal alpha)	Non	Oui	Oui.
Prise en charge de l’animation	Non	Limitée (APNG)	Oui.
Algorithme de base	DCT + Huffman	DEFLATE (LZ77 + Huffman)	VP8 / VP8L.
Idéal pour	Photographie	Graphiques, texte et icônes	Usage universel en environnement Web.
Moins adapté pour	Texte et contours nets	Images photographiques complexes	Compatibilité héritée.

5.8.2 Métriques d’Évaluation de la Qualité

Deux métriques objectives sont largement adoptées pour quantifier la distorsion introduite par les processus de compression :

Pic du Rapport Signal-à-Bruit (PSNR, Peak Signal-to-Noise Ratio) : \[ \text{PSNR} = 10\,\log_{10}\!\left(\frac{L^2}{\text{MSE}}\right) \quad [\text{dB}] \tag{5.10}\]

où \(L = 255\) pour les images quantifiées sur 8 bits et \(\text{MSE}\) représente l’Erreur Quadratique Moyenne (Mean Squared Error). Des valeurs de PSNR supérieures à 40 dB indiquent une excellente fidélité ; entre 30 dB et 40 dB, elles représentent une bonne qualité ; et des valeurs inférieures à 30 dB correspondent à des dégradations visuelles facilement perceptibles.

Indice de Similarité Structurale (SSIM, Structural Similarity Index) : \[ \text{SSIM}(f,g) = \frac{(2\mu_f\mu_g + c_1)(2\sigma_{fg} + c_2)}{(\mu_f^2+\mu_g^2+c_1)(\sigma_f^2+\sigma_g^2+c_2)} \tag{5.11}\]

Le SSIM évalue des fenêtres locales de l’image sur la base de trois composantes complémentaires : la luminance (\(\mu_f, \mu_g\)), le contraste (\(\sigma_f, \sigma_g\)) et la structure (\(\sigma_{fg}\)), pondérées par des constantes de stabilité \(c_1\) et \(c_2\). L’indice varie dans l’intervalle \([-1, 1]\), où l’unité représente l’identité parfaite. Contrairement au PSNR, le SSIM prend en compte l’organisation spatiale des erreurs, s’alignant ainsi sur la perception du système visuel humain (SVH).

PSNR vs SSIM : Application de Métriques Perceptuelles

Le PSNR possède une formulation mathématique simple et un faible coût computationnel ; toutefois, il tend à surestimer la qualité des images présentant des distorsions localisées ou à la sous-estimer en cas de variations globales de luminosité tolérées par l’observateur. Le SSIM modélise plus fidèlement la perception biologique, mais exige un effort de traitement plus important. Pour des analyses rigoureuses des codecs, il est recommandé de rapporter les deux métriques statistiques à titre complémentaire.

5.8.3 Inspection Visuelle : Nature des Artefacts de Compression

La nature mathématique du codec détermine le type de dégradation introduit à des débits binaires réduits. Comme illustré dans la Figure 5.27, la compression agressive via DCT dans le standard JPEG segmente l’image en mailles rigides, générant les artefacts de bloc (blocking artifacts). En revanche, les algorithmes basés sur la codification prédictive ou les représentations soumises à des transformées spatiales avancées (comme le WebP et le JPEG 2000) éliminent les discontinuités de bloc, mais introduisent une perte de texture fine et des flous caractéristiques autour des contours à fort contraste.

%%writefile tmp/fig_05_zoom_artefatos.cpp
#define MM_OUT "tmp/fig_05_zoom_artefatos.png"
//| label: fig-05-zoom-artefatos
//| fig-cap: "Análise comparativa de artefatos de compressão sob fator de qualidade reduzido ($Q=10$). À esquerda, observa-se o artefato de bloco característico da discretização por DCT no JPEG. À direita, evidencia-se o efeito de atenuação e suavização de bordas intrínseco ao padrão WebP."
//| echo: true
//| output: true
#include <opencv2/opencv.hpp>
#include <vector>
#include <string>
#include <filesystem>
#include "morph.hpp"

int main() {
    // Lê a imagem, converte para escala de cinza e redimensiona
    cv::Mat img_src_mat = cv::imread("imagens/cameraman.png", cv::IMREAD_GRAYSCALE);
    cv::Mat img_src_resized;
    cv::resize(img_src_mat, img_src_resized, cv::Size(256, 256));
    mm::Image img_src = img_src_resized;

    // Salva com qualidade Q=10 em JPEG e WebP
    cv::imwrite("tmp/zoom_q10.jpg", cv::Mat(img_src), {cv::IMWRITE_JPEG_QUALITY, 10});
    cv::imwrite("tmp/zoom_q10.webp", cv::Mat(img_src), {cv::IMWRITE_WEBP_QUALITY, 10});

    // Função de zoom com interpolação vizinho mais próximo
    auto zoom = [](const cv::Mat& img) {
        cv::Mat crop = img(cv::Rect(150, 120, 80, 80));
        cv::Mat enlarged;
        cv::resize(crop, enlarged, cv::Size(320, 320), 0, 0, cv::INTER_NEAREST);
        return enlarged;
    };

    // Aplica zoom nas três imagens
    cv::Mat z_orig = zoom(cv::Mat(img_src));
    cv::Mat z_jpeg = zoom(cv::imread("tmp/zoom_q10.jpg", cv::IMREAD_GRAYSCALE));
    cv::Mat z_webp = zoom(cv::imread("tmp/zoom_q10.webp", cv::IMREAD_GRAYSCALE));

    // Exibe as três imagens lado a lado
    std::vector<mm::Image> images = {
        mm::Image(z_orig),
        mm::Image(z_jpeg),
        mm::Image(z_webp)
    };
    std::vector<std::string> titles = {
        "Zoom Original",
        "JPEG Q=10 (Artefato de Bloco)",
        "WebP Q=10 (Suavizacao)"
    };
    mm::show(images, MM_OUT, titles, 3);

    
// [pdi:panel-io] auto-generated — do not edit by hand
std::filesystem::create_directories("tmp");
mm::write(z_orig, "tmp/fig_05_zoom_artefatos_0.png");
mm::write(z_jpeg, "tmp/fig_05_zoom_artefatos_1.png");
mm::write(z_webp, "tmp/fig_05_zoom_artefatos_2.png");
// [pdi:panel-io:end]
return 0;
}

Overwriting tmp/fig_05_zoom_artefatos.cpp

!g++ -I. -std=c++17 -DMM_USE_OPENCV -I/usr/include/opencv4 tmp/fig_05_zoom_artefatos.cpp -o tmp/fig_05_zoom_artefatos -lopencv_stitching -lopencv_alphamat -lopencv_aruco -lopencv_barcode -lopencv_bgsegm -lopencv_bioinspired -lopencv_ccalib -lopencv_dnn_objdetect -lopencv_dnn_superres -lopencv_dpm -lopencv_face -lopencv_freetype -lopencv_fuzzy -lopencv_hdf -lopencv_hfs -lopencv_img_hash -lopencv_intensity_transform -lopencv_line_descriptor -lopencv_mcc -lopencv_quality -lopencv_rapid -lopencv_reg -lopencv_rgbd -lopencv_saliency -lopencv_shape -lopencv_stereo -lopencv_structured_light -lopencv_phase_unwrapping -lopencv_superres -lopencv_optflow -lopencv_surface_matching -lopencv_tracking -lopencv_highgui -lopencv_datasets -lopencv_text -lopencv_plot -lopencv_ml -lopencv_videostab -lopencv_videoio -lopencv_viz -lopencv_wechat_qrcode -lopencv_ximgproc -lopencv_video -lopencv_xobjdetect -lopencv_objdetect -lopencv_calib3d -lopencv_imgcodecs -lopencv_features2d -lopencv_dnn -lopencv_flann -lopencv_xphoto -lopencv_photo -lopencv_imgproc -lopencv_core \
  && ./tmp/fig_05_zoom_artefatos \
  && test -f "tmp/fig_05_zoom_artefatos.png" \
  || echo "⚠ mm::show não gravou tmp/fig_05_zoom_artefatos.png"

[1] Zoom Original
[2] JPEG Q=10 (Artefato de Bloco)
[3] WebP Q=10 (Suavizacao)

try:
    mm.show(
        [
            mm.read("tmp/fig_05_zoom_artefatos_0.png"),
            mm.read("tmp/fig_05_zoom_artefatos_1.png"),
            mm.read("tmp/fig_05_zoom_artefatos_2.png"),
        ],
        titles=[
            'Zoom Original',
            'JPEG Q=10 (Artefato de Bloco)',
            'WebP Q=10 (Suavizacao)',
        ],
        cols=3,
        figsize=(14, 5),
    )
except Exception as _e:
    print("figura indisponivel nesta trilha (C++): " + repr(_e) + " tmp/fig_05_zoom_artefatos_0.png (ver a versao Python)")

Figure 5.27: Análise comparativa de artefatos de compressão sob fator de qualidade reduzido (\(Q=10\)). À esquerda, observa-se o artefato de bloco característico da discretização por DCT no JPEG. À direita, evidencia-se o efeito de atenuação e suavização de bordas intrínseco ao padrão WebP.

5.8.4 Évaluation Quantitative et Spatiale de la Compression

La validation des algorithmes de compression avec perte exige une analyse qui corrèle le coût de stockage à la fidélité du signal reconstruit. Cette évaluation est réalisée de manière complémentaire à travers des courbes de performance globale et par la cartographie locale des distorsions induites par les codeurs.

5.8.4.1 Courbes Débit-Distorsion

La Figure 5.28 présente l’évaluation empirique du pipeline JPEG et WebP au moyen de courbes débit-distorsion, qui surveillent le gain de compression (taille du fichier en Ko) en fonction du PSNR. Le format PNG sert de ligne de base idéale (\(\text{PSNR} = \infty\)), car sa nature lossless empêche toute dégradation, bien qu’il exige un volume de données substantiellement plus important.

L’analyse des courbes démontre la supériorité et l’efficacité du standard WebP par rapport au JPEG traditionnel : pour atteindre un même niveau de fidélité mathématique (comme la plage d’excellente qualité, où \(\text{PSNR} > 40\text{ dB}\)), le codeur WebP génère des fichiers significativement plus petits. Ce comportement traduit l’impact pratique de l’évolution des algorithmes sur l’optimisation des systèmes de transmission et de stockage numérique.

%%writefile tmp/fig_05_formatos_comparacao.cpp
#define MM_OUT "tmp/fig_05_formatos_comparacao.png"
//| label: fig-05-formatos-comparacao
//| fig-cap: "Curva taxa-distorção: PSNR vs tamanho de arquivo para JPEG, WebP e PNG aplicada à imagem do *Cameraman*."
//| echo: true
//| output: true

#include <opencv2/opencv.hpp>
#include <vector>
#include <string>
#include <filesystem>
#include <cstdio>
#include "morph.hpp"

int main() {
    // Lecture et redimensionnement de l'image Cameraman
    cv::Mat src_full = cv::imread("imagens/cameraman.png", cv::IMREAD_GRAYSCALE);
    cv::Mat src_resized;
    cv::resize(src_full, src_resized, cv::Size(256, 256));
    mm::Image src = src_resized;

    // Tableaux pour courbe JPEG
    std::vector<double> jpeg_kb, jpeg_psnr;
    std::vector<int> jpeg_qualities = {10, 20, 30, 40, 50, 60, 70, 80, 90, 95};
    for (int q : jpeg_qualities) {
        std::string p = "tmp/fmt_q" + std::to_string(q) + ".jpg";
        cv::imwrite(p, cv::Mat(src), {cv::IMWRITE_JPEG_QUALITY, q});
        cv::Mat rec = cv::imread(p, cv::IMREAD_GRAYSCALE);
        jpeg_kb.push_back((double)std::filesystem::file_size(p) / 1024.0);
        jpeg_psnr.push_back(mm::psnr(src, mm::Image(rec)));
    }

    // Tableaux pour courbe WebP
    std::vector<double> webp_kb, webp_psnr;
    std::vector<int> webp_qualities = {30, 50, 70, 85, 95};
    for (int q : webp_qualities) {
        std::string p = "tmp/fmt_w" + std::to_string(q) + ".webp";
        cv::imwrite(p, cv::Mat(src), {cv::IMWRITE_WEBP_QUALITY, q});
        cv::Mat rec = cv::imread(p, cv::IMREAD_GRAYSCALE);
        webp_kb.push_back((double)std::filesystem::file_size(p) / 1024.0);
        webp_psnr.push_back(mm::psnr(src, mm::Image(rec)));
    }

    // Compression PNG sans perte
    std::string p_png = "tmp/fmt.png";
    cv::imwrite(p_png, cv::Mat(src), {cv::IMWRITE_PNG_COMPRESSION, 9});
    double png_kb = (double)std::filesystem::file_size(p_png) / 1024.0;

    // Créer le graphique comparatif
    std::vector<std::vector<double>> xs = {jpeg_kb, webp_kb};
    std::vector<std::vector<double>> ys = {jpeg_psnr, webp_psnr};
    std::vector<std::string> labels = {"JPEG", "WebP"};
    char title[128];
    snprintf(title, sizeof(title), "Curva Taxa-Distorcao (PNG sem perda: %.1f KB)", png_kb);

    mm::Image chart = mm::lineChart(
        xs, ys, labels, {}, title,
        "Tamanho do arquivo (KB)", "PSNR (dB)"
    );

    // Affichage du graphique
    mm::show(std::vector<mm::Image>{chart}, MM_OUT, {"JPEG x WebP x PNG"}, 1);
    
// [pdi:panel-io] auto-generated — do not edit by hand
std::filesystem::create_directories("tmp");
mm::write(chart, "tmp/fig_05_formatos_comparacao_0.png");
// [pdi:panel-io:end]
return 0;
}

Overwriting tmp/fig_05_formatos_comparacao.cpp

!g++ -I. -std=c++17 -DMM_USE_OPENCV -I/usr/include/opencv4 tmp/fig_05_formatos_comparacao.cpp -o tmp/fig_05_formatos_comparacao -lopencv_stitching -lopencv_alphamat -lopencv_aruco -lopencv_barcode -lopencv_bgsegm -lopencv_bioinspired -lopencv_ccalib -lopencv_dnn_objdetect -lopencv_dnn_superres -lopencv_dpm -lopencv_face -lopencv_freetype -lopencv_fuzzy -lopencv_hdf -lopencv_hfs -lopencv_img_hash -lopencv_intensity_transform -lopencv_line_descriptor -lopencv_mcc -lopencv_quality -lopencv_rapid -lopencv_reg -lopencv_rgbd -lopencv_saliency -lopencv_shape -lopencv_stereo -lopencv_structured_light -lopencv_phase_unwrapping -lopencv_superres -lopencv_optflow -lopencv_surface_matching -lopencv_tracking -lopencv_highgui -lopencv_datasets -lopencv_text -lopencv_plot -lopencv_ml -lopencv_videostab -lopencv_videoio -lopencv_viz -lopencv_wechat_qrcode -lopencv_ximgproc -lopencv_video -lopencv_xobjdetect -lopencv_objdetect -lopencv_calib3d -lopencv_imgcodecs -lopencv_features2d -lopencv_dnn -lopencv_flann -lopencv_xphoto -lopencv_photo -lopencv_imgproc -lopencv_core \
  && ./tmp/fig_05_formatos_comparacao \
  && test -f "tmp/fig_05_formatos_comparacao.png" \
  || echo "⚠ mm::show não gravou tmp/fig_05_formatos_comparacao.png"

[1] JPEG x WebP x PNG

try:
    mm.show(
        [
            mm.read("tmp/fig_05_formatos_comparacao_0.png"),
        ],
        titles=[
            'JPEG x WebP x PNG',
        ],
        cols=1,
    )
except Exception as _e:
    print("figura indisponivel nesta trilha (C++): " + repr(_e) + " tmp/fig_05_formatos_comparacao_0.png (ver a versao Python)")

Figure 5.28: Curva taxa-distorção: PSNR vs tamanho de arquivo para JPEG, WebP e PNG aplicada à imagem do *Cameraman*.

Taille originale de l’image

L’image Cameraman (\(256 \times 256\) pixels en niveaux de gris) occupe 64 Ko en format brut (sans compression). À titre de référence, le PNG lossless compresse ce volume à 36,2 Ko — mettant en évidence que la compression sans perte réduit déjà significativement le stockage pour les images comportant des régions homogènes. En revanche, les formats avec perte (JPEG et WebP) atteignent des tailles encore plus réduites : le JPEG avec une qualité de 95 occupe 22,3 Ko (PSNR ≈ 45 dB), tandis que le WebP avec une qualité de 90 atteint 12,5 Ko avec un PSNR équivalent, démontrant sa supériorité en matière d’efficacité de compression.

5.8.4.2 Cartographie Spatiale des Erreurs et Corrélation Perceptuelle

Bien que le PSNR offre un indicateur numérique rapide, les métriques globales ne parviennent pas à distinguer comment la perte d’informations se répartit géométriquement sur l’image. La Figure 5.29 résout cette limitation en associant les reconstructions à différentes qualités à leurs cartes d’erreur absolue respectives et au SSIM.

Les cartes résiduelles — obtenues par la différence absolue normalisée entre l’image originale et l’image compressée — révèlent la signature spatiale intrinsèque de chaque architecture de codage :

À des qualités élevées (\(Q=95\) à \(Q=75\)) : Les distorsions se concentrent principalement autour des transitions abruptes d’intensité (bords), résultant du repliement spectral dû à l’élimination des hautes fréquences. L’indice SSIM reste proche de l’unité, attestant de l’intégrité des structures originales.
À des qualités agressives (\(Q=50\) à \(Q=25\)) : L’erreur adopte une structure de maillage orthogonal régularisé. Ce motif géométrique met en évidence l’apparition des artefacts de bloc (blocking artifacts), indiquant que la quantification sévère a corrompu la corrélation spatiale entre les blocs adjacents de \(8 \times 8\) pixels.

Le SSIM capture cette dégradation morphologique de manière beaucoup plus sensible que le PSNR, pénalisant le score final à mesure que l’organisation structurelle et les textures fines — auxquelles le système visuel humain est hautement réactif — sont éliminées par le codeur.

%%writefile tmp/fig_05_ssim_artefatos.cpp
#define MM_OUT "tmp/fig_05_ssim_artefatos.png"
//| label: fig-05-ssim-artefatos
//| fig-cap: "Analyse spatiale de la dégradation : images reconstruites et cartes d'erreur absolue normalisées pour différents facteurs de qualité JPEG."
//| echo: true
//| output: true

#include <opencv2/opencv.hpp>
#include <vector>
#include <string>
#include <cmath>
#include "morph.hpp"

// Fonction SSIM (Wang et al., 2004) — version fenêtrée avec filtre gaussien
// Retourne la valeur moyenne SSIM et la carte SSIM complète.
static double compute_ssim_full(const cv::Mat& a, const cv::Mat& b, cv::Mat& ssim_map) {
    // Convertir en flottant 64 bits
    cv::Mat A, B;
    a.convertTo(A, CV_64F);
    b.convertTo(B, CV_64F);

    const double C1 = 6.5025;  // (0.01 * 255)^2
    const double C2 = 58.5225; // (0.03 * 255)^2

    // Filtre gaussien 11x11, sigma = 1.5
    cv::Mat mu1, mu2, mu1_sq, mu2_sq, mu1_mu2, sigma1_sq, sigma2_sq, sigma12;
    cv::GaussianBlur(A, mu1, cv::Size(11, 11), 1.5);
    cv::GaussianBlur(B, mu2, cv::Size(11, 11), 1.5);
    cv::multiply(mu1, mu1, mu1_sq);
    cv::multiply(mu2, mu2, mu2_sq);
    cv::multiply(mu1, mu2, mu1_mu2);

    cv::Mat A_sq, B_sq, A_B;
    cv::multiply(A, A, A_sq);
    cv::multiply(B, B, B_sq);
    cv::multiply(A, B, A_B);

    cv::Mat sigma1_sq_tmp, sigma2_sq_tmp, sigma12_tmp;
    cv::GaussianBlur(A_sq, sigma1_sq_tmp, cv::Size(11, 11), 1.5);
    cv::GaussianBlur(B_sq, sigma2_sq_tmp, cv::Size(11, 11), 1.5);
    cv::GaussianBlur(A_B, sigma12_tmp, cv::Size(11, 11), 1.5);

    cv::subtract(sigma1_sq_tmp, mu1_sq, sigma1_sq);
    cv::subtract(sigma2_sq_tmp, mu2_sq, sigma2_sq);
    cv::subtract(sigma12_tmp, mu1_mu2, sigma12);

    cv::Mat cs_map = (2 * sigma12 + C2) / (sigma1_sq + sigma2_sq + C2);
    cv::Mat ssim_map_tmp = ((2 * mu1_mu2 + C1) / (mu1_sq + mu2_sq + C1));
    cv::multiply(ssim_map_tmp, cs_map, ssim_map);

    // Valeur moyenne
    cv::Scalar mssim = cv::mean(ssim_map);
    return mssim[0];
}

int main() {
    // Cameraman (asset du chapitre), 256x256 — même image que la piste py.
    std::string img_path = "imagens/cameraman.png";
    cv::Mat img_full = cv::imread(img_path, cv::IMREAD_GRAYSCALE);
    if (img_full.empty()) {
        // Fallback : lire via mm::read si le chemin direct échoue
        mm::Image tmp = mm::read(img_path);
        img_full = cv::Mat(tmp.h, tmp.w, CV_8UC1, tmp.data.data());
    }
    cv::Mat img_resized;
    cv::resize(img_full, img_resized, cv::Size(256, 256));
    mm::Image img_gray = img_resized;

    std::vector<mm::Image> imgs_ssim;
    std::vector<std::string> titles_ssim;
    imgs_ssim.push_back(img_gray);
    titles_ssim.push_back("Original");

    for (int q : {25, 50, 75, 95}) {
        // DCT 8x8 -> quantisation -> IDCT
        mm::Image rec = mm::jpegCompress(img_gray, q);
        double psnr_v = mm::psnr(img_gray, rec);

        // Convertir en cv::Mat pour le calcul SSIM
        cv::Mat img_cv = img_gray;   // conversion implicite
        cv::Mat rec_cv = rec;        // conversion implicite

        cv::Mat ssim_map;
        double ssim_v = compute_ssim_full(img_cv, rec_cv, ssim_map);

        // Carte d'erreur absolue normalisée (0..255)
        cv::Mat img_f, rec_f;
        img_cv.convertTo(img_f, CV_64F);
        rec_cv.convertTo(rec_f, CV_64F);
        cv::Mat diff = cv::abs(img_f - rec_f);
        cv::Mat diff_norm;
        cv::normalize(diff, diff_norm, 0, 255, cv::NORM_MINMAX);
        cv::Mat diff_vis;
        diff_norm.convertTo(diff_vis, CV_8U);

        imgs_ssim.push_back(rec);
        imgs_ssim.push_back(mm::Image(diff_vis));
        titles_ssim.push_back("Q=" + std::to_string(q) + " (PSNR=" + 
                             std::to_string(psnr_v).substr(0, 3) + "dB | SSIM=" + 
                             std::to_string(ssim_v).substr(0, 4) + ")");
        titles_ssim.push_back("Carte d'erreur (Q=" + std::to_string(q) + 
                             ") - bords et blocage");
    }

    mm::show(imgs_ssim, MM_OUT, titles_ssim, 3);
    return 0;
}

Overwriting tmp/fig_05_ssim_artefatos.cpp

!g++ -I. -std=c++17 -DMM_USE_OPENCV -I/usr/include/opencv4 tmp/fig_05_ssim_artefatos.cpp -o tmp/fig_05_ssim_artefatos -lopencv_stitching -lopencv_alphamat -lopencv_aruco -lopencv_barcode -lopencv_bgsegm -lopencv_bioinspired -lopencv_ccalib -lopencv_dnn_objdetect -lopencv_dnn_superres -lopencv_dpm -lopencv_face -lopencv_freetype -lopencv_fuzzy -lopencv_hdf -lopencv_hfs -lopencv_img_hash -lopencv_intensity_transform -lopencv_line_descriptor -lopencv_mcc -lopencv_quality -lopencv_rapid -lopencv_reg -lopencv_rgbd -lopencv_saliency -lopencv_shape -lopencv_stereo -lopencv_structured_light -lopencv_phase_unwrapping -lopencv_superres -lopencv_optflow -lopencv_surface_matching -lopencv_tracking -lopencv_highgui -lopencv_datasets -lopencv_text -lopencv_plot -lopencv_ml -lopencv_videostab -lopencv_videoio -lopencv_viz -lopencv_wechat_qrcode -lopencv_ximgproc -lopencv_video -lopencv_xobjdetect -lopencv_objdetect -lopencv_calib3d -lopencv_imgcodecs -lopencv_features2d -lopencv_dnn -lopencv_flann -lopencv_xphoto -lopencv_photo -lopencv_imgproc -lopencv_core \
  && ./tmp/fig_05_ssim_artefatos \
  && test -f "tmp/fig_05_ssim_artefatos.png" \
  || echo "⚠ mm::show não gravou tmp/fig_05_ssim_artefatos.png"

[1] Original
[2] Q=25 (PSNR=30.dB | SSIM=0.86)
[3] Carte d'erreur (Q=25) - bords et blocage
[4] Q=50 (PSNR=32.dB | SSIM=0.90)
[5] Carte d'erreur (Q=50) - bords et blocage
[6] Q=75 (PSNR=35.dB | SSIM=0.93)
[7] Carte d'erreur (Q=75) - bords et blocage
[8] Q=95 (PSNR=44.dB | SSIM=0.98)
[9] Carte d'erreur (Q=95) - bords et blocage

try:
    mm.show(mm.read("tmp/fig_05_ssim_artefatos.png"), figsize=(14, 14))
except Exception as _e:
    print("figura indisponivel nesta trilha (C++): " + repr(_e) + " tmp/fig_05_ssim_artefatos.png (ver a versao Python)")

Figure 5.29: Análise espacial de degradação: imagens reconstruídas e respectivos mapas de erro absoluto normalizados para diferentes fatores de qualidade JPEG.

Interprétation des cartes d’erreur

Les cartes d’erreur présentées ont été normalisées individuellement (cv2.NORM_MINMAX) afin de maximiser le contraste visuel et de révéler la structure spatiale des distorsions. Cela signifie que :

Pour Q=95, l’erreur absolue est de l’ordre de 0,5 à 1,5 niveaux de gris (imperceptible visuellement), mais la normalisation l’amplifie en noir et blanc pour mettre en évidence sa localisation au niveau des bords et des transitions.
Pour Q=25, l’erreur absolue est 10 à 20 fois plus grande (5 à 15 niveaux de gris), mais la normalisation l’amène également dans la même plage [0, 255].

Par conséquent, l’intensité du blanc dans les cartes n’est PAS comparable entre différentes qualités — les cartes servent uniquement à révéler la signature spatiale de l’erreur (bords vs blocs), et non son ampleur. L’ampleur correcte est fournie par les valeurs PSNR et SSIM, qui montrent clairement que Q=95 présente une erreur bien inférieure à celle de Q=25.

Synthèse — Compression JPEG

Le processus de compression dans le standard JPEG repose sur l’application combinée de transformations spatiales, perceptuelles et statistiques pour réduire les redondances d’une image. La Table 5.9 résume le rôle de chaque étape dans le pipeline et son impact respectif sur la réduction des données.

Table 5.9: Synthèse des étapes du pipeline de compression JPEG et de leurs impacts respectifs.

Étape	Opération Analytique	Mécanisme de Gain / Compression
Conversion \(YC_bC_r\)	Isolation des canaux de luminance et de chrominance.	Modélise la perception du système visuel humain (SVH), permettant de traiter la couleur et la luminosité de manière indépendante.
Sous-échantillonnage 4:2:0	Réduction de la résolution spatiale des canaux de couleur (\(C_b\) et \(C_r\)).	Élimine environ 50 % des données brutes avec un impact visuel minimal.
DCT \(8 \times 8\)	Mappage du domaine spatial vers le domaine des fréquences spatiales.	Compactage de l’énergie, concentrant l’information essentielle dans les premiers coefficients.
Quantification Linéaire	Division entière des coefficients par une matrice de pondération \(Q(u,v)\).	Principale source de compression avec perte ; élimine les hautes fréquences imperceptibles.
Codage Entropique	Application d’algorithmes RLE et de codage de Huffman.	Compression statistique sans perte, optimisée par les longues séquences de coefficients nuls.

Artefacts de dégradation caractéristiques

L’application de taux de compression excessivement agressifs (facteurs de qualité réduits) introduit des distorsions prévisibles dans l’image reconstruite, découlant des limitations mathématiques du modèle :

Artefacts de bloc (blocking artifacts) : Discontinuités géométriques visibles aux frontières des blocs de \(8 \times 8\) pixels, causées par la perte de corrélation spatiale après la quantification sévère des composantes AC.
Effet de sonnerie (ringing) : Oscillations fantômes ou distorsions de « fumée » autour des bords nets et à fort contraste, provoquées par l’élimination brutale des harmoniques de haute fréquence nécessaires pour reconstruire les fonctions échelon.
Perte de texture fine : Atténuation des détails à haute fréquence et à faible contraste (comme les pelouses, les tissus ou la porosité), donnant aux régions initialement texturées un aspect excessivement lisse ou homogénéisé.

5.9 Application Pratique : Débruitage par Filtrage Hybride

En réunissant les techniques consolidées tout au long de ce chapitre, on présente un pipeline complet de restauration d’images qui combine l’analyse spectrale dans le domaine fréquentiel avec le filtrage adaptatif dans le domaine spatial. L’objectif est d’atténuer un bruit mixte (composé d’une dégradation gaussienne et d’une interférence périodique) tout en préservant au maximum les détails structurels de l’image originale.

\[ \text{Image Bruitée} \xrightarrow{\text{FFT2}} \xrightarrow{\text{Filtre Notch Gaussien}} \xrightarrow{\text{IFFT2}} \xrightarrow{\text{Filtre Bilatéral}} \text{Image Restaurée} \]

Évaluation Complémentaire : PSNR vs. SSIM

La paire de métriques statistiques PSNR et SSIM fournit une évaluation qualitative et morphologique complémentaire du processus de restauration :

PSNR : Pénalise uniformément l’écart quadratique moyen pixel par pixel.
SSIM : Évalue la préservation des structures locales perceptuellement pertinentes (luminance, contraste et contours).

En pratique, il existe un compromis analytique (trade-off) entre réduction du bruit et préservation des détails : des filtres spatiaux excessivement agressifs atténuent bien le bruit haute fréquence, mais dégradent les textures fines et lissent les contours nets — ce qui réduit simultanément à la fois le PSNR et le SSIM par rapport à l’image originale. Le défi de la conception de filtres est de trouver le point d’équilibre qui maximise les deux métriques, garantissant une restauration fidèle et visuellement agréable.

5.9.1 Analyse des performances et conclusion du chapitre

Les résultats numériques et visuels générés par Figure 5.30 démontrent la pertinence pratique d’associer différents domaines de traitement. L’insertion simultanée de bruit périodique et stochastique corrompt les propriétés morphologiques du signal, réduisant sévèrement les indices de similarité et le rapport signal-bruit de l’image de référence.

L’isolation et la suppression des pics harmoniques dans le domaine fréquentiel au moyen du masque notch éliminent les franges d’interférence sinusoïdales réparties sur l’espace bidimensionnel. Comme le montrent les données imprimées de Figure 5.30, ce filtrage chirurgical induit un saut immédiat et substantiel de la métrique PSNR. Toutefois, le bruit gaussien haute fréquence reste actif de manière homogène dans le spectre, exigeant une approche complémentaire.

La restauration finale est consolidée dans le domaine spatial avec l’introduction du filtre bilatéral. Contrairement aux opérateurs passe-bas classiques (tels que le filtre gaussien ou le filtre moyenneur), qui lisseraient indifféremment le bruit et les contours structurels, le filtrage bilatéral calcule des poids pondérés par la proximité géométrique et par la différence d’intensité radiométrique. Ce comportement adaptatif atténue les fluctuations stochastiques résiduelles dans les régions de transition douce et préserve la netteté des bords spatiaux.

La convergence des deux approches aboutit à une amélioration substantielle et simultanée du PSNR et du SSIM par rapport à l’image bruitée — bien que les valeurs finales demeurent inférieures à celles de l’image originale (PSNR = \(\infty\), SSIM = 1,0), en raison de la perte inévitable d’informations spectrales et texturales lors des processus de filtrage. L’atténuation douce (gaussienne) des pics dans le spectre évite les artefacts de ringing, tandis que le filtre bilatéral élimine le bruit stochastique résiduel sans compromettre la netteté des bords. Les résultats prouvent l’efficacité et la complémentarité pratique des outils d’analyse fréquentielle présentés dans ce chapitre, démontrant que le filtrage hybride (fréquence + spatial) est supérieur à toute approche isolée pour la restauration d’images dégradées par un bruit mixte.

%%writefile tmp/fig_05_pipeline_denoising.cpp
#define MM_OUT "tmp/fig_05_pipeline_denoising.png"
#include <opencv2/opencv.hpp>
#include <cmath>
#include <random>
#include <vector>
#include <string>
#include <iostream>
#include "morph.hpp"

int main() {
    // Cameraman (asset do capítulo), 256x256 — mesma imagem da trilha py.
    cv::Mat cam = cv::imread("imagens/cameraman.png", cv::IMREAD_COLOR);
    cv::resize(cam, cam, cv::Size(256, 256));
    mm::Image img_gray = mm::gray(cam);
    int h_img = img_gray.h;
    int w_img = img_gray.w;

    // ── 1. Ruído misto: gaussiano + periódico ───────────────────────────────
    std::mt19937 rng(42);
    std::normal_distribution<double> dist(0.0, 15.0);

    // u0, v0 frequências da interferência periódica
    int u0 = 15, v0 = 10;

    // X2, Y2 = malhas de coordenadas
    cv::Mat X2(h_img, w_img, CV_64F);
    cv::Mat Y2(h_img, w_img, CV_64F);
    for (int y = 0; y < h_img; ++y) {
        for (int x = 0; x < w_img; ++x) {
            X2.at<double>(y, x) = (double)x;
            Y2.at<double>(y, x) = (double)y;
        }
    }

    // constrói imagem ruidosa
    cv::Mat img_noisy_cv(h_img, w_img, CV_8U);
    for (int y = 0; y < h_img; ++y) {
        for (int x = 0; x < w_img; ++x) {
            double g = img_gray.data[y * w_img + x];
            double ruido_gauss = dist(rng);
            double ruido_period = 30.0 * std::sin(2.0 * M_PI *
                (u0 * X2.at<double>(y, x) / w_img + v0 * Y2.at<double>(y, x) / h_img));
            double v = g + ruido_gauss + ruido_period;
            v = std::max(0.0, std::min(255.0, v));
            img_noisy_cv.at<unsigned char>(y, x) = (unsigned char)std::lround(v);
        }
    }
    mm::Image img_noisy(img_noisy_cv);

    // ── 2. Espectro (log-magnitude) da imagem ruidosa ───────────────────────
    mm::Image mag_n = mm::spectrumMag(img_noisy);

    // ── 3. Máscara notch gaussiana nos 4 picos periódicos ───────────────────
    auto suprimir_pico_gaussiano = [](cv::Mat& mask, int cy, int cx, double sigma) {
        int M = mask.rows, N = mask.cols;
        for (int yy = 0; yy < M; ++yy) {
            for (int xx = 0; xx < N; ++xx) {
                double notch = std::exp(-((double)((yy - cy) * (yy - cy)) +
                                          (double)((xx - cx) * (xx - cx))) /
                                        (2.0 * sigma * sigma));
                mask.at<double>(yy, xx) *= (1.0 - notch);
            }
        }
    };

    int cy0 = h_img / 2, cx0 = w_img / 2;
    cv::Mat mascara_notch(h_img, w_img, CV_64F, cv::Scalar(1.0));
    int offsets[4][2] = {{v0, u0}, {-v0, -u0}, {v0, -u0}, {-v0, u0}};
    for (auto& off : offsets) {
        suprimir_pico_gaussiano(mascara_notch, cy0 + off[0], cx0 + off[1], 3.0);
    }

    // ── 4. Filtragem: notch no domínio da frequência + bilateral ────────────
    mm::Image img_notch = mm::freqFilter(img_noisy, mascara_notch);
    cv::Mat img_notch_cv = img_notch;
    cv::Mat img_den_cv;
    cv::bilateralFilter(img_notch_cv, img_den_cv, 7, 25, 7);
    mm::Image img_den(img_den_cv);

    cv::Mat mascara_vis_cv;
    mascara_notch.convertTo(mascara_vis_cv, CV_8U, 255.0);
    mm::Image mascara_vis(mascara_vis_cv);

    double psnr_n  = mm::psnr(img_gray, img_noisy);
    double psnr_no = mm::psnr(img_gray, img_notch);
    double psnr_d  = mm::psnr(img_gray, img_den);
    std::cout << "Ruidosa: PSNR=" << psnr_n << " dB | Apos notch: " << psnr_no
              << " dB | Notch+bilateral: " << psnr_d << " dB" << std::endl;

    char buf[128];
    std::snprintf(buf, sizeof(buf), "Ruidosa (PSNR=%.1f dB)", psnr_n);
    std::string t1 = buf;
    std::snprintf(buf, sizeof(buf), "Apos notch (PSNR=%.1f dB)", psnr_no);
    std::string t2 = buf;
    std::snprintf(buf, sizeof(buf), "Notch + bilateral (PSNR=%.1f dB)", psnr_d);
    std::string t3 = buf;

    mm::show(
        {img_gray, img_noisy, mag_n, mascara_vis, img_notch, img_den},
        MM_OUT,
        {
            "Original",
            t1,
            "Espectro (picos visiveis)",
            "Mascara notch (gaussiana)",
            t2,
            t3,
        },
        6
    );

    return 0;
}

Overwriting tmp/fig_05_pipeline_denoising.cpp

!g++ -I. -std=c++17 -DMM_USE_OPENCV -I/usr/include/opencv4 tmp/fig_05_pipeline_denoising.cpp -o tmp/fig_05_pipeline_denoising -lopencv_stitching -lopencv_alphamat -lopencv_aruco -lopencv_barcode -lopencv_bgsegm -lopencv_bioinspired -lopencv_ccalib -lopencv_dnn_objdetect -lopencv_dnn_superres -lopencv_dpm -lopencv_face -lopencv_freetype -lopencv_fuzzy -lopencv_hdf -lopencv_hfs -lopencv_img_hash -lopencv_intensity_transform -lopencv_line_descriptor -lopencv_mcc -lopencv_quality -lopencv_rapid -lopencv_reg -lopencv_rgbd -lopencv_saliency -lopencv_shape -lopencv_stereo -lopencv_structured_light -lopencv_phase_unwrapping -lopencv_superres -lopencv_optflow -lopencv_surface_matching -lopencv_tracking -lopencv_highgui -lopencv_datasets -lopencv_text -lopencv_plot -lopencv_ml -lopencv_videostab -lopencv_videoio -lopencv_viz -lopencv_wechat_qrcode -lopencv_ximgproc -lopencv_video -lopencv_xobjdetect -lopencv_objdetect -lopencv_calib3d -lopencv_imgcodecs -lopencv_features2d -lopencv_dnn -lopencv_flann -lopencv_xphoto -lopencv_photo -lopencv_imgproc -lopencv_core \
  && ./tmp/fig_05_pipeline_denoising \
  && test -f "tmp/fig_05_pipeline_denoising.png" \
  || echo "⚠ mm::show não gravou tmp/fig_05_pipeline_denoising.png"

Ruidosa: PSNR=20.1968 dB | Apos notch: 22.3552 dB | Notch+bilateral: 23.5442 dB
[1] Original
[2] Ruidosa (PSNR=20.2 dB)
[3] Espectro (picos visiveis)
[4] Mascara notch (gaussiana)
[5] Apos notch (PSNR=22.4 dB)
[6] Notch + bilateral (PSNR=23.5 dB)

try:
    mm.show(mm.read("tmp/fig_05_pipeline_denoising.png"), figsize=(20, 4))
except Exception as _e:
    print("figura indisponivel nesta trilha (C++): " + repr(_e) + " tmp/fig_05_pipeline_denoising.png (ver a versao Python)")

Figure 5.30: *Pipeline* completo de remoção de ruído misto: (1) adição de ruído gaussiano e periódico; (2) identificação de picos de interferência no espectro de frequências; (3) aplicação de máscara *notch* com atenuação gaussiana suave; (4) pós-processamento via filtro bilateral para eliminação do ruído estocástico residual.

5.10 Résumé du Chapitre

La transition du domaine spatial vers le domaine fréquentiel révèle la distribution spectrale de l’énergie de l’image, établissant ainsi la base analytique pour le filtrage avancé, la restauration et la compression des données. L’articulation structurelle de ces concepts est synthétisée dans la carte conceptuelle de la Figure 5.31.

Figure 5.31: Carte conceptuelle des transformations et propriétés dans le domaine fréquentiel.

Fundamentos Essenciais

TFD et Perception Visuelle : Le spectre décompose l’image en composantes harmoniques. La phase conserve l’intelligibilité géométrique de la scène et la localisation des contours, tandis que la magnitude détermine la distribution du contraste et les amplitudes globales.
Efficacité Algorithmique : Le Théorème de Convolution permet le traitement de masques à grande échelle dans le domaine fréquentiel via la FFT, réduisant la complexité computationnelle asymptotique de \(O(N^2 K^2)\) dans l’espace à \(O(N^2 \log N)\).
Phénomène de Ringing : Les coupures abruptes dans le spectre (filtres idéaux) génèrent des oscillations spatiales indésirables (phénomène de Gibbs). L’atténuation douce par des filtres de Butterworth ou gaussiens élimine ces discontinuités.
Analyse Multirésolution via Wavelets : Dépassant le caractère purement global de Fourier, la DWT capture simultanément la fréquence et la localisation spatiale, fondant le standard JPEG 2000 et appuyant des représentations hiérarchiques analogues aux extractions de caractéristiques dans les Réseaux de Neurones Convolutifs (CNN).
Compression Perceptuelle (DCT) : Le pipeline JPEG exploite les limitations de contraste du système visuel humain aux hautes fréquences spatiales. La DCT isole l’énergie de blocs \(8 \times 8\), permettant à la quantification d’éliminer les coefficients AC des détails fins sans préjudice perceptuel sévère.

Prochaines Étapes : Le Chapitre 6 inaugure la Partie II de l’ouvrage, appliquant les outils de traitement d’images à la résolution de problèmes réels d’inspection industrielle. Seront explorées des techniques de segmentation et d’analyse de formes pour la détection automatique de défauts sur les lignes de production — depuis l’identification de défauts superficiels sur les pièces jusqu’à la lecture QRCode sur les épreuves, consolidant le pont entre la théorie présentée dans la Partie I et les exigences pratiques de la vision par ordinateur.

5.11 🤖 Utilisation de Gemini Notebook comme Tuteur Complémentaire

Dans cette édition, l’utilisation de la plateforme Gemini Notebook est encouragée comme outil d’apprentissage complémentaire — et non comme substitut de la lecture attentive, de la résolution d’exercices ou de l’expérimentation pratique. Fondé sur des architectures d’intelligence artificielle, le système utilise exclusivement le matériel pédagogique et les documents fournis par l’auteur comme base de connaissances, garantissant que les réponses générées soient conceptuellement alignées sur le contenu programmatique et l’approche pédagogique adoptée tout au long de cet ouvrage.

Accès au Tuteur Intelligent

🚀 ACCÉDER À Gemini Notebook : CHAPITRE 05

🌐 Langue et Langage de Programmation

Le projet de ce chapitre dans Gemini Notebook a été construit uniquement avec le texte en portugais et les exemples de code en Python. Si vous étudiez à partir de l’édition en anglais ou en français, ou si vous suivez le parcours en C++, les réponses du tuteur peuvent ne pas correspondre exactement à la version que vous lisez.

Directives concernant le Contenu Généré par Intelligence Artificielle

Bien que les outils d’intelligence artificielle constituent des alliés efficaces dans le processus d’apprentissage et de révision, le contenu généré est sujet à des incohérences ou à des imprécisions techniques. Par conséquent, la consultation systématique de manuels, d’articles scientifiques et de sources académiques indexées est indispensable pour une validation rigoureuse des informations. Il est vivement recommandé d’exécuter et de modifier les exemples pratiques en Python fournis dans ce chapitre comme méthode principale de vérification expérimentale des résultats.

5.12 Liste d’exercices

(10 %) Implémentation directe de la TFD 2D : Implémentez analytiquement la Transformée de Fourier Discrète 2D (TFD) sans l’aide de fonctions natives de bibliothèques (comme np.fft.fft2), en utilisant strictement la formulation mathématique définie dans la Équation 5.1 pour une matrice de dimensions \(16 \times 16\). Effectuez la validation numérique en comparant les coefficients générés avec les résultats de la fonction np.fft.fft2, en vous assurant que l’écart absolu maximal soit inférieur à \(10^{-8}\). Mesurez les temps d’exécution des deux méthodes et présentez une justification théorique pour la disparité observée en termes de complexité asymptotique.
(15 %) Suppression du bruit périodique : Ajoutez des interférences sinusoïdales avec des fréquences spatiales \((u_0, v_0) \in \{(5,10), (20,5), (30,30)\}\) à l’image de test du Cameraman. Pour chaque scénario de dégradation, concevez un masque de filtrage notch spécifique dans le domaine fréquentiel afin d’isoler et d’atténuer les pics harmoniques indésirables. Évaluez quantitativement l’efficacité du processus de restauration par le calcul des métriques PSNR et SSIM. Discutez analytiquement du compromis entre l’atténuation du bruit sinusoïdal et l’atténuation indésirable des caractéristiques structurelles légitimes de l’image.
(15 %) Analyse comparative des opérateurs passe-bas : Réalisez une étude comparative entre les filtres passe-bas idéal, gaussien et de Butterworth (avec des ordres harmoniques \(n = 1, 2, 4\)), paramétrés avec des fréquences de coupure \(D_0 = 20, 40, 60\) pixels. Pour chaque combinaison structurelle, calculez les indices PSNR et SSIM de l’image résultante par rapport au signal original de référence. Organisez les données quantitatives dans un tableau structuré et tracez les graphiques unidimensionnels des fonctions de transfert correspondantes le long du profil horizontal \(H(u, 0)\).
(15 %) Banque de filtres multirésolution de Haar : Développez un script pour exécuter manuellement la décomposition wavelet discrète 2D de premier niveau en utilisant la famille de Haar. L’algorithme doit calculer les coefficients des filtres correspondants passe-bas (\(h\)) et passe-haut (\(g\)), en les appliquant de manière séparable sur les lignes et les colonnes de la matrice, suivis de l’opération de décimation (sous-échantillonnage spatial par un facteur de 2). Validez numériquement l’exactitude de votre implémentation en confrontant les sous-bandes obtenues avec la sortie de la fonction pywt.dwt2(img, 'haar').
(15 %) Compression éparse par seuillage wavelet : Appliquez la technique de filtrage par seuillage abrupt (hard thresholding) sur les coefficients de détail de la décomposition wavelet, en adoptant les seuils numériques \(T \in \{5, 10, 20, 40, 80\}\) pour les familles de Haar, Daubechies (db4) et Symlets (sym4). Après avoir réalisé le processus de synthèse au moyen de la transformée inverse (pywt.waverec2), calculez les valeurs de PSNR et SSIM de chaque image reconstruite. Identifiez et justifiez quelle combinaison de famille wavelet et de seuil \(T\) maximise la similarité structurelle.
(15 %) Construction d’un encodeur JPEG simplifié : Implémentez le pipeline complet de compression de données simulant la norme JPEG. Le flux doit englober : la conversion spatiale \(RGB \rightarrow YC_bC_r\), le sous-échantillonnage chromatique dans la proportion 4:2:0, la segmentation de la luminance en blocs disjoints de \(8 \times 8\) pixels, l’application de la DCT-II 2D orthogonale, et la quantification linéaire basée sur la matrice normalisée de luminance mise à l’échelle par des facteurs de qualité souhaités. Réalisez le décodage inverse et comparez quantitativement les reconstructions avec les fichiers générés par la fonction cv2.imencode pour les facteurs de qualité de 20, 50 et 80.
(15 %) Analyse perceptuelle sur des contenus hétérogènes : Développez une image synthétique composée de trois régions distinctes et aux caractéristiques spectrales contrastées : une texture photographique complexe (représentant de hautes fréquences stochastiques), une zone de texte vectorisé avec des bords nets (représentant des transitions en échelon pures) et un gradient linéaire continu (représentant de basses fréquences homogènes). Soumettez cette image mixte aux processus de compression sous les formats JPEG, PNG et WebP. Évaluez et interprétez les résultats en corrélant la taille finale du fichier sur disque avec les métriques PSNR et SSIM obtenues, en justifiant quel format présente les meilleures performances pour des signaux de nature hétérogène et pourquoi cet avantage survient en termes de compactage d’énergie et de préservation perceptuelle.

Références du chapitre

Le fondement théorique et le développement analytique des concepts abordés dans ce chapitre s’appuient sur les ouvrages de référence suivants :

Gonzalez (2018) — Formulations classiques des Transformées de Fourier Discrètes 2D (TFD), conception de filtres analytiques dans le domaine fréquentiel, Transformée en Cosinus Discrète (TCD) et principes fondamentaux des systèmes de compression d’images.
Oppenheim (2010) — Théorie formelle des signaux et des systèmes appliqués dans le domaine discret, couvrant les propriétés mathématiques de la TFD et la modélisation analytique du Théorème de Convolution.
Mallat (1999) — Fondements mathématiques de la théorie des ondelettes, formalisation de l’analyse multirésolution (AMR) et architecture des bancs de filtres dyadiques.
Wallace (1991) — Spécification originale et aspects techniques du standard de compression ISO/CEI JPEG, avec un accent sur les critères psychovisuels pour la conception des matrices de quantification TCD.
Szeliski (2022) — Modélisation computationnelle et caractérisation des métriques modernes de fidélité et de qualité perceptuelle (PSNR et SSIM), ainsi que l’analyse comparative des formats d’images tramées à hautes performances.

5.13 💻 Partie pratique avec exercices de programmation

La présente liste d’exercices de programmation (EP) consolide les formulations théoriques présentées tout au long du Chapitre 5 — Transformées et Compression — au moyen d’un parcours pratique appliqué. Les exercices sont structurés à partir de matrices de dimensions réduites, permettant la validation analytique et l’inspection manuelle de chaque coefficient, tout en maintenant la cohérence méthodologique adoptée dans les chapitres précédents.

L’enchaînement des exercices reproduit rigoureusement le flux conceptuel du chapitre : on commence par l’implémentation explicite de la Transformée de Fourier Discrète (TFD) à partir de sa définition mathématique fondamentale ; on progresse vers la conception de filtres passe-bas et de masques notch dans le domaine fréquentiel ; on applique la quantification des coefficients (noyau de la compression avec perte) ; et l’on conclut par l’intégration de ces étapes dans la construction d’un pipeline de compression JPEG simplifié et dans l’analyse perceptuelle des formats d’image.

Directives pour la résolution des exercices de programmation

Dans tous les exercices de ce chapitre, les coordonnées du centre du spectre (origine des fréquences spatiales après application du décalage fftshift) doivent être déterminées par division entière. Pour une matrice de \(L\) lignes et \(C\) colonnes, la composante de fréquence nulle se situe à la position :

\[ (c_y, c_x) = \left( \left\lfloor \frac{L}{2} \right\rfloor, \left\lfloor \frac{C}{2} \right\rfloor \right) \]

Cette convention est rigoureusement identique à celle adoptée par la fonction np.fft.fftshift. De plus, à toutes les étapes exigeant une discrétisation ou un arrondi numérique (que ce soit dans la quantification des coefficients AC ou dans la reconstruction finale des pixels), on doit employer l’arrondi standard au plus proche entier (round half away from zero), atténuant les ambiguïtés sur les valeurs dont la fraction est exactement égale à \(0.5\).

🎯 Objectif de ce cahier

Ce cahier permet de développer, valider, organiser et tester des solutions d’Exercices de Programmation (EPs) dans des environnements interactifs, comme Colab, avec les mêmes cas de test que Moodle, en les y copiant uniquement au moment d’enregistrer la note officielle.

Téléchargement

Téléchargez morph.py et testsuite.py en exécutant la cellule ci-dessous :

import os, urllib.request

os.makedirs("tmp/state", exist_ok=True)

url = "https://raw.githubusercontent.com/fzampirolli/pdi-vc/master/morph/config.py"
if not os.path.exists("config.py"):
    urllib.request.urlretrieve(url, "config.py")

import config
config.setup(testsuite=True, cpp=True)
from morph import mm
from testsuite import TestSuite

✅ Environnement prêt. Morph : 1.1.9 | OpenCV : 5.0.0 | TestSuite : 1.1.2

Exécution des tests

Pour évaluer les tests, exécutez TestSuite("EP05_01.extensão").run() dans une nouvelle cellule, en remplaçant l’extension par celle du langage utilisé (.py, .java, .c, .cpp, .js ou .r). Le système télécharge les cas de test depuis GitHub, exécute le programme et calcule automatiquement la note.

Pour tester directement du code Python, sans sauvegarder de fichier, utilisez run_code(codigo) en passant le code sous forme de chaîne de caractères dans une variable codigo :

codigo = """
from morph import mm
# ... votre code ici ...
"""
TestSuite("EP05_01").run_code(codigo)

5.13.1 EP05_01 🟢 Filtre Passe-Bas Idéal par Distance dans le Spectre

Dans un scanner de documents anciens, le capteur capte le papier froissé et la texture des fibres en même temps que le texte — un bruit haute fréquence qui « pollue » le spectre sur les bords. Le technicien de maintenance n’a pas accès à l’image originale, seulement au spectre de magnitude déjà calculé par le logiciel du scanner. Son travail est simple et chirurgical : ne conserver que le cercle central des basses fréquences (la structure globale du document) et effacer tout ce qui se trouve hors du rayon \(D_0\), éliminant la texture fine sans même avoir à toucher à l’image spatiale.

C’est le filtre passe-bas idéal (LPFI) : l’opération spectrale la plus directe du chapitre, mais aussi celle qui révèle le mieux l’anatomie d’un spectre centré.

5.13.1.1 📋 Directives d’implémentation

Dimensions : Lire les entiers \(L\) (lignes) et \(C\) (colonnes) du spectre de magnitude — déjà fourni centré (équivalent à la sortie de np.fft.fftshift).
Fréquence de coupure : Lire l’entier \(D_0\).
Données : Lire les valeurs entières de la matrice de magnitude, ligne par ligne.
Centre du spectre : Calculer \((c_y, c_x) = (L \mathbin{//} 2,\; C \mathbin{//} 2)\).
Distance : Pour chaque position \((u,v)\), calculer \[ D(u,v) = \sqrt{(u-c_y)^2 + (v-c_x)^2} \]
Masque idéal : Appliquer \[ H(u,v) = \begin{cases} 1, & D(u,v) \le D_0 \\ 0, & D(u,v) > D_0 \end{cases} \]
Filtrage : La valeur de sortie est \(\text{mag}'(u,v) = \text{mag}(u,v) \cdot H(u,v)\).
Sortie : Afficher la matrice filtrée avec les dimensions \(L \times C\).

5.13.1.2 📌 Contraintes computationnelles

Comparaison non stricte : le critère utilise \(D(u,v) \le D_0\) (la frontière appartient au filtre, c’est-à-dire qu’elle est conservée).
Type : toutes les valeurs d’entrée et de sortie sont des entiers ; la distance est calculée en virgule flottante uniquement en interne.
Pas d’arrondi de magnitude : comme l’entrée est déjà entière et que le masque est binaire (0 ou 1), la sortie n’a jamais besoin d’arrondi.

5.13.1.3 🧠 Fondement théorique

Région	Distance au centre	Effet du filtre
Centre (\(D \le D_0\))	Basses fréquences	Préservées — structure globale conservée
Bords (\(D > D_0\))	Hautes fréquences	Mises à zéro — texture et bruit supprimés
\(D_0\) petit	—	L’image reconstruite serait très floue
\(D_0\) grand	—	Peu de filtrage ; presque toute l’énergie est préservée

5.13.1.4 📦 Spécification d’entrée et de sortie (VPL)

Entrée :

Ligne 1 : Entier \(L\).
Ligne 2 : Entier \(C\).
Ligne 3 : Entier \(D_0\).
Lignes suivantes : Éléments entiers de la matrice de magnitude (centrée).

Sortie :

Matrice filtrée en \(L\) lignes et \(C\) colonnes, séparés par des espaces.

5.13.1.5 📌 Exemples

Entrée	Sortie	Observation
3 3 1 10 20 30 40 50 60 70 80 90	0 20 0 40 50 60 0 80 0	Centre \((1,1)\). Les coins ont \(D=\sqrt{2}\approx1.41 > 1\), donc ils sont mis à zéro ; les voisins orthogonaux ont \(D=1 \le 1\) et sont conservés.
1 3 0 5 9 7	0 9 0	\(L=1, C=3\) : centre en \((0,1)\). Seule la position centrale elle-même (\(D=0\)) survit à \(D_0=0\).

🎮 Simulateur EP05_01 : Filtre Passe-Bas Idéal H = (D ≤ D₀) ? 1 : 0

Rayon de coupure (D₀) : 1

Ajustez D₀ et observez quelles positions du spectre 5×5 survivent au filtre.

Spectre Original (Magnitude)

Résultat Filtré

–

Figure 5.32: Simulateur EP05_01 : Filtre passe-bas idéal dans le spectre

%%writefile EP05_01.cpp
// your solution

Overwriting EP05_01.cpp

TestSuite("EP05_01.cpp").run()

✔️ EP05_01.cases existe déjà dans casos/
📋 5 cas chargé(s) depuis casos/EP05_01.cases

🔍 Test de C++ : EP05_01.cpp
⚠️ EP05_01.cpp : fichier vide (moins de 3 lignes). Tests ignorés.

5.13.2 EP05_02 🟡 Filtre Notch : Suppression des pics périodiques

Une caméra d’inspection industrielle capture des images de circuits imprimés, mais l’alimentation électrique de la ligne de production introduit une interférence électrique périodique — un motif de stries quasi imperceptible à l’œil nu, mais qui apparaît dans le spectre de Fourier comme des paires de pics brillants symétriquement positionnés autour du centre. L’équipe de vision par ordinateur ne peut pas recapturer l’image : elle doit localiser et effacer chirurgicalement ces paires de pics dans le spectre, tout en préservant le reste de l’information utile de l’image.

C’est le rôle du filtre coupe-bande notch : contrairement au passe-bas (qui affecte une région continue), il cible des points spécifiques et leurs symétriques, laissant le reste du spectre intact.

5.13.2.1 📋 Directives d’implémentation

Dimensions : Lire les entiers \(L\) (lignes) et \(C\) (colonnes) du spectre de magnitude centré.
Données : Lire les valeurs entières de la matrice de magnitude, ligne par ligne.
Pics : Lire l’entier \(K\) (nombre de paires de pics à supprimer).
Pour chacun des \(K\) pics : lire trois entiers \(\Delta v\), \(\Delta u\), \(r\) — déplacement vertical, déplacement horizontal et rayon du notch.
Centre du spectre : \((c_y, c_x) = (L \mathbin{//} 2,\; C \mathbin{//} 2)\).
Suppression symétrique : pour chaque pic, mettre à zéro toutes les positions \((u,v)\) telles que la distance au point \((c_y+\Delta v,\, c_x+\Delta u)\) soit \(\le r\), et également toutes les positions à distance \(\le r\) du point symétrique \((c_y-\Delta v,\, c_x-\Delta u)\).
Sortie : Afficher la matrice résultante avec les dimensions \(L \times C\).

5.13.2.2 📌 Contraintes computationnelles

Symétrie obligatoire : chaque pic fourni génère deux disques mis à zéro (le point et son symétrique par rapport au centre) — oublier le symétrique est l’erreur la plus courante.
Chevauchement : si deux disques se chevauchent, la position reste à zéro (ni « addition » ni restauration).
Comparaison non stricte : une position est mise à zéro si \(\text{distance} \le r\).
Ordre de lecture : les \(K\) pics doivent être traités dans l’ordre où ils apparaissent en entrée, mais le résultat final est indépendant de l’ordre (les opérations de mise à zéro sont commutatives).

5.13.2.3 🧠 Fondement théorique

Concept	Rôle dans le filtre notch
Pic en \((\Delta v, \Delta u)\)	Fréquence de l’interférence périodique détectée visuellement dans le spectre
Point symétrique \((-\Delta v,-\Delta u)\)	Toute DFT d’un signal réel est hermitienne : les pics apparaissent toujours en paires symétriques par rapport au centre
Rayon \(r\)	Contrôle la « largeur » de la réjection — un \(r\) grand élimine davantage d’énergie autour du pic, mais aussi davantage d’information utile

5.13.2.4 📦 Spécification d’entrée et de sortie (VPL)

Entrée :

Ligne 1 : Entier \(L\).
Ligne 2 : Entier \(C\).
Lignes suivantes : Éléments entiers de la matrice de magnitude (centrée), \(L\) lignes.
Ligne suivante : Entier \(K\).
\(K\) lignes suivantes : trois entiers \(\Delta v\), \(\Delta u\), \(r\) (séparés par des espaces).

Sortie :

Matrice résultante sur \(L\) lignes et \(C\) colonnes, séparés par des espaces.

5.13.2.5 📌 Exemples

Entrée	Sortie	Observation
5 5 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 1 1 1 0	1 2 3 4 5 6 0 8 9 10 11 12 13 14 15 16 17 18 0 20 21 22 23 24 25	Centre \((c_y, c_x) = (2, 2)\). Pic fourni \((\Delta v, \Delta u) = (1, 1)\) génère le point \((3, 3)\) (valeur 19) et son symétrique \((1, 1)\) (valeur 7), tous deux mis à zéro avec \(r=0\) (seuls les points exacts).

🎮 Simulateur EP05_02 : Filtre Notch Paire symétrique

Δv 1

Δu 1

Rayon (r) 0

Déplacez Δv e Δu pour choisir le pic — observez que la paire symétrique est également filtrée.

Spectre 5×5 (Rouge = Supprimé par le filtre)

–

Figure 5.33: Simulateur EP05_02: Filtre Notch

%%writefile EP05_02.cpp
// your solution

Overwriting EP05_02.cpp

TestSuite("EP05_02.cpp").run()

✔️ EP05_02.cases existe déjà dans casos/
📋 5 cas chargé(s) depuis casos/EP05_02.cases

🔍 Test de C++ : EP05_02.cpp
⚠️ EP05_02.cpp : fichier vide (moins de 3 lignes). Tests ignorés.

5.13.3 EP05_03 🟠 Quantification DCT : la véritable source de compression

Une application de galerie de photos doit réduire la taille de milliers d’images avant de les téléverser vers le cloud, sans tout recoder de zéro. L’ingénieur responsable dispose déjà des coefficients DCT de chaque bloc \(4\times4\) calculés (l’étape coûteuse en calcul est déjà faite) — il ne reste plus qu’à appliquer la table de quantification, l’étape qui élimine réellement de l’information et génère la compression. Les coefficients haute fréquence, moins perceptibles à l’œil humain, reçoivent de grands diviseurs et tendent à devenir zéro ; les coefficients basse fréquence, plus perceptibles, reçoivent de petits diviseurs et survivent presque intacts.

Vous allez implémenter exactement cette étape : quantifier et déquantifier (diviser, arrondir, multiplier en retour) — le cœur de la compression lossy du JPEG.

5.13.3.1 📋 Directives d’implémentation

Dimension du bloc : Lire l’entier \(N\) (bloc \(N \times N\)).
Coefficients : Lire la matrice \(C\) des coefficients DCT, \(N\) lignes avec \(N\) entiers chacune (ils peuvent être négatifs).
Table de quantification : Lire la matrice \(Q\), \(N\) lignes avec \(N\) entiers positifs chacune.
Quantification : Pour chaque position \((u,v)\), calculer l’indice quantifié \[ \tilde{C}(u,v) = \text{round}\!\left(\frac{C(u,v)}{Q(u,v)}\right) \] en utilisant l’arrondi standard à l’entier le plus proche (les valeurs intermédiaires .5 ne se produisent jamais dans les cas de test).
Déquantification (reconstruction) : Calculer \[ C'(u,v) = \tilde{C}(u,v) \times Q(u,v) \]
Sortie : Afficher la matrice reconstruite \(C'\), \(N \times N\), entiers.

5.13.3.2 📌 Contraintes informatiques

Aller-retour complet : la sortie est le coefficient reconstruit (\(\tilde{C} \times Q\)), pas l’indice quantifié isolé.
Division en virgule flottante : la division \(C(u,v)/Q(u,v)\) doit être effectuée en virgule flottante avant l’arrondi — une division entière tronquée produirait un résultat incorrect.
Signe préservé : les coefficients négatifs conservent leur signe après quantification et reconstruction.
\(Q(u,v) > 0\) toujours : aucune gestion de division par zéro n’est nécessaire.

5.13.3.3 🧠 Fondements théoriques

Coefficient	Fréquence	Valeur typique de \(Q\)	Effet de la quantification
\(C(0,0)\)	DC (moyenne du bloc)	Petite	Survit presque toujours — domine l’énergie
\(C(u,v)\) faible \(u+v\)	Basse fréquence	Petite/moyenne	Partiellement préservé
\(C(u,v)\) élevé \(u+v\)	Haute fréquence	Grande	Devient souvent zéro — source de la compression

5.13.3.4 📦 Spécification d’entrée et de sortie (VPL)

Entrée :

Ligne 1 : Entier \(N\).
\(N\) lignes suivantes : matrice \(C\) (coefficients DCT, entiers, peuvent être négatifs).
\(N\) lignes suivantes : matrice \(Q\) (table de quantification, entiers positifs).

Sortie :

Matrice reconstruite \(C'\), \(N \times N\), entiers séparés par des espaces.

5.13.3.5 📌 Exemples

Entrée	Sortie	Observation
4 50 10 -5 0 8 -3 2 1 0 1 0 0 2 0 0 -1 2 5 7 8 4 7 8 11 6 8 11 12 9 11 12 14	50 10 -7 0 8 0 0 0 0 0 0 0 0 0 0 0	\(C(0,0)=50/2=25 \to 25\times2=50\) (préservé). \(C(0,2)=-5/7\approx-0.71\to-1\to-1\times7=-7\). Quant à \(C(1,1)=-3/7\approx-0.43\to0\) : mis à zéro par la quantification — la majeure partie du bloc devient zéro, illustrant la compaction de l’énergie dans le coin supérieur gauche.

🎮 Simulateur EP05_03 : Quantification DCT round(C / Q) × Q

Échelle de Q (Agressivité) : 1.00×

Ajustez l'échelle de Q et voyez combien de coefficients survivent (non nuls) après l'aller-retour.

Coefficients DCT (C)

Reconstruit (round(C / Q) · Q)

–

Figure 5.34: Simulateur EP05_03 : Quantification DCT (round-trip)

%%writefile EP05_03.cpp
// your solution

Overwriting EP05_03.cpp

TestSuite("EP05_03.cpp").run()

✔️ EP05_03.cases existe déjà dans casos/
📋 5 cas chargé(s) depuis casos/EP05_03.cases

🔍 Test de C++ : EP05_03.cpp
⚠️ EP05_03.cpp : fichier vide (moins de 3 lignes). Tests ignorés.

5.13.4 EP05_04 🔴 Implémentation de la TFD 2D à partir de la définition

Un laboratoire de recherche en astronomie computationnelle a reçu, d’une mission ancienne, un petit capteur expérimental dont les données brutes ne peuvent pas être traitées par les bibliothèques modernes de FFT — l’environnement de validation est isolé et ne permet que des opérations arithmétiques de base. L’équipe doit réimplémenter la Transformée de Fourier Discrète 2D à partir de la définition mathématique elle-même, cellule par cellule, pour ensuite comparer bit à bit avec np.fft.fft2 dans un autre environnement.

C’est l’exercice le plus conceptuel de la liste : il n’y a pas de raccourcis. Vous allez implémenter la double sommation de la Équation 5.1 directement, en mettant en évidence pourquoi la FFT existe — et le coût computationnel qu’elle évite.

5.13.4.1 📋 Directives d’implémentation

Dimensions : Lire les entiers \(M\) (lignes) et \(N\) (colonnes) de l’image \(f(x,y)\).
Données : Lire les valeurs entières de \(f(x,y)\), ligne par ligne.
TFD 2D : Pour chaque paire de fréquences \((u,v)\) avec \(u=0,\ldots,M-1\) et \(v=0,\ldots,N-1\), calculer \[ F(u,v) = \sum_{x=0}^{M-1}\sum_{y=0}^{N-1} f(x,y)\, e^{-j2\pi\left(\frac{ux}{M}+\frac{vy}{N}\right)} \] en utilisant l’identité d’Euler \(e^{-j\theta} = \cos(\theta) - j\sin(\theta)\) pour séparer les parties réelle et imaginaire — n’utilisez aucune fonction de FFT prête à l’emploi.
Magnitude : Calculer \(|F(u,v)| = \sqrt{\text{Re}(F)^2 + \text{Im}(F)^2}\) et arrondir à l’entier le plus proche.
Sortie : Afficher la matrice des magnitudes arrondies, \(M \times N\), dans le même ordre (sans fftshift — la composante DC reste en \((0,0)\)).

5.13.4.2 📌 Contraintes computationnelles

Interdiction d’utiliser des bibliothèques de FFT : l’implémentation doit calculer les double sommations explicitement (boucles imbriquées), même si c’est plus lent.
Sans fftshift : la sortie conserve la convention brute de la TFD, avec la composante DC en \(F(0,0)\) (coin supérieur gauche).
Arrondi : la magnitude finale doit être arrondie à l’entier le plus proche ; dans les cas de test, il n’y a pas d’ambiguïté .5.
Précision : de petites erreurs de virgule flottante (de l’ordre de \(10^{-6}\)) avant l’arrondi sont attendues et n’affectent pas le résultat entier final.

5.13.4.3 🧠 Fondement théorique

Élément	Signification
\(F(0,0)\)	Composante DC — somme de tous les pixels, \(F(0,0) = \sum f(x,y)\)
Partie réelle \(\text{Re}(F)\)	Projection du signal sur les cosinus
Partie imaginaire \(\text{Im}(F)\)	Projection du signal sur les sinus
Complexité de cette implémentation	\(\mathcal{O}((MN)^2)\) — c’est pourquoi la FFT, avec \(\mathcal{O}(MN\log(MN))\), est indispensable pour les images réelles

5.13.4.4 📦 Spécification d’entrée et de sortie (VPL)

Entrée :

Ligne 1 : Entier \(M\).
Ligne 2 : Entier \(N\).
Lignes suivantes : Éléments entiers de \(f(x,y)\), \(M\) lignes.

Sortie :

Matrice des magnitudes \(|F(u,v)|\) arrondies, \(M \times N\), séparées par des espaces.

5.13.4.5 📌 Exemples

Entrée	Sortie	Remarque
2 2 1 2 3 4	10 2 4 0	\(F(0,0)=1+2+3+4=10\) (DC = somme totale). \(F(0,1)=(1-2)+(3-4)=-2 \to \|F\|=2\). \(F(1,0)=(1+2)-(3+4)=-4\to\|F\|=4\). \(F(1,1)=(1-2)-(3-4)=0\).

🎮 Simulateur EP05_04 : DFT 2D — Définition directe ΣΣ f(x,y) e^-j2π(…)

Cliquez sur les cellules de f(x,y) pour modifier les valeurs (incrément +1 ; Maj + clic décrément -1) et voyez |F(u,v)| recalculé en direct.

f(x,y) — Domaine spatial

|F(u,v)| — Magnitude (sans décalage)

–

Figure 5.35: Simulateur EP05_04: DFT 2D manuel

%%writefile EP05_04.cpp
// your solution

Overwriting EP05_04.cpp

TestSuite("EP05_04.cpp").run()

✔️ EP05_04.cases existe déjà dans casos/
📋 5 cas chargé(s) depuis casos/EP05_04.cases

🔍 Test de C++ : EP05_04.cpp
⚠️ EP05_04.cpp : fichier vide (moins de 3 lignes). Tests ignorés.

5.13.5 EP05_05 🏆 Pipeline JPEG complet : DCT, quantification et reconstruction

Vous avez été chargé de créer, à partir de zéro, un codec JPEG didactique dans un environnement embarqué, sans aucune bibliothèque d’image disponible — uniquement des opérations mathématiques de base. Le client veut comprendre exactement où la qualité est perdue et où elle est récupérée, bloc par bloc. C’est le défi final du chapitre : intégrer tout ce qui a été étudié — la DCT-II orthonormale, la quantification perceptuelle et la reconstruction via IDCT — dans un seul pipeline de bout en bout, traitant un bloc \(N \times N\) du début à la fin, exactement comme le fait le standard JPEG en interne, \(8\times8\) pixels à la fois.

5.13.5.1 📋 Directives d’implémentation

Dimension du bloc : Lire l’entier \(N\).
Bloc original : Lire la matrice de pixels \(f(x,y)\), \(N\) lignes avec \(N\) entiers dans \([0,255]\).
Table de quantification : Lire la matrice \(Q\), \(N \times N\) entiers positifs.
Centrage : Soustraire 128 de chaque pixel : \(g(x,y) = f(x,y) - 128\).
DCT-II 2D orthonormale : Calculer \[ C(u,v) = \alpha(u)\,\alpha(v)\sum_{x=0}^{N-1}\sum_{y=0}^{N-1} g(x,y)\,\cos\!\left[\frac{\pi(2x+1)u}{2N}\right]\cos\!\left[\frac{\pi(2y+1)v}{2N}\right] \] avec \(\alpha(0)=\sqrt{1/N}\) et \(\alpha(k)=\sqrt{2/N}\) pour \(k>0\).
Quantification : \(\tilde{C}(u,v) = \text{round}(C(u,v)/Q(u,v))\).
Déquantification : \(C'(u,v) = \tilde{C}(u,v)\times Q(u,v)\).
IDCT-II 2D (inverse orthonormale) : Calculer \(g'(x,y)\) à partir de \(C'(u,v)\) en utilisant la transformée inverse correspondante (même base, somme sur \(u,v\)).
Inversion du centrage et arrondi : \(f'(x,y) = \text{round}(g'(x,y) + 128)\), restreint à l’intervalle \([0,255]\) (clipping).
Sortie : Afficher le bloc reconstruit \(f'\), \(N \times N\), entiers.

5.13.5.2 📌 Contraintes computationnelles

Pipeline complet obligatoire : toutes les six étapes (centrer, DCT, quantifier, déquantifier, IDCT, inverser) doivent être implémentées — sauter la quantification ne réussit pas les tests, car le résultat serait identique à l’original.
Clipping : les valeurs reconstruites hors de \([0,255]\) doivent être tronquées (0 si négatif, 255 si supérieur à 255).
Arrondi : à la fois dans la quantification et dans la reconstruction finale des pixels, utilisez un arrondi standard ; les cas de test évitent toute ambiguïté .5.
Base orthonormale : la normalisation \(\alpha(u)\) et \(\alpha(v)\) doit être appliquée exactement comme spécifié — sans elle, l’IDCT ne reconstruit pas correctement.

5.13.5.3 🧠 Fondement théorique

Étape	Analogue dans le standard JPEG réel	Où la qualité est perdue
Centrage	Identique — la DCT suppose un signal centré sur zéro	Aucune perte
DCT-II	Étapes 3–4 du pipeline (Table 5.7)	Aucune perte (transformation exacte et réversible)
Quantification	Étape 5 — division par \(Q(u,v)\)	Principale source de perte — les coefficients de haute fréquence deviennent zéro
IDCT	Reconstruction finale	Reconstruit exactement les coefficients quantifiés, pas les originaux

5.13.5.4 📦 Spécification d’entrée et de sortie (VPL)

Entrée :

Ligne 1 : Entier \(N\).
\(N\) lignes suivantes : bloc original \(f(x,y)\), entiers dans \([0,255]\).
\(N\) lignes suivantes : table de quantification \(Q\), entiers positifs.

Sortie :

Bloc reconstruit \(f'(x,y)\), \(N \times N\), entiers dans \([0,255]\), séparés par des espaces.

5.13.5.5 📌 Exemples

Entrée	Sortie	Observation
4 120 130 125 128 115 140 135 122 118 150 160 130 110 120 145 138 4 6 8 10 6 8 10 12 8 10 12 16 10 12 16 20	118 126 119 131 114 143 140 119 117 149 159 130 107 121 146 139	Après DCT, quantification agressive des hautes fréquences (grandes valeurs de \(Q\) en bas à droite) et reconstruction via IDCT, le bloc reste proche de l’original, mais pas identique — la différence est le coût de la compression lossy.

5.13.5.6 💡 Conseil de débogage

Si le résultat ne correspond pas, vérifiez dans cet ordre : (1) les coefficients DCT bruts (avant quantification) — ils doivent reconstruire l’original exactement via IDCT si vous sautez les étapes 6–7 ; (2) la table \(\alpha(u)\) — erreur courante : appliquer \(\sqrt{2/N}\) aussi pour \(u=0\) ; (3) l’arrondi de la quantification, qui doit se produire avant de multiplier à nouveau par \(Q\).

🎮 Simulateur EP05_05 : Pipeline JPEG (Bloc 4×4) DCT → Q → IDCT

Échelle de Q (1 = Table de base, Plus grand = Plus de perte) : 1.00×

Ajustez le facteur d'échelle de quantisation et observez le bloc reconstruit s'éloigner (ou se rapprocher) de l'original.

Bloc original

Reconstruit (DCT → Q → IDCT)

–

Figure 5.36: Simulateur EP05_05 : Pipeline JPEG complet par bloc

%%writefile EP05_05.cpp
// your solution

Overwriting EP05_05.cpp

TestSuite("EP05_05.cpp").run()

✔️ EP05_05.cases existe déjà dans casos/
📋 5 cas chargé(s) depuis casos/EP05_05.cases

🔍 Test de C++ : EP05_05.cpp
⚠️ EP05_05.cpp : fichier vide (moins de 3 lignes). Tests ignorés.