---
title: INFO911 (4c) Filtrage fréquentiel
type: slide
slideOptions:
transition: slide
progress: true
slideNumber: true
---
# Filtrage fréquentiel
## (Traitement et Analyse d'Image 4c)
> [name=Jacques-Olivier Lachaud][time=Novembre 2020][color=#907bf7]
> (Les images peuvent être soumises à des droits d'auteur. Elles sont utilisées ici exclusivement dans un but pédagogique)
###### tags: `info911`
Retour à [INFO911 (Main) Traitement et analyse d'image](https://codimd.math.cnrs.fr/s/UE_B59gMy)
---
# Fonctions périodiques, série de Fourier
==[Fourier]== Toute fonction **périodique** $f(x)$ de période $2\pi$ peut être approchée par une somme de fonctions $a_n \cos( n x)$ et $b_n\sin( n x)$, avec $n$ entier $\ge 0$.
| Approximations successives $f_N(x)$ | Décomposition sur harmoniques |
| -------- | -------- |
|  |  |
---
## Décomposition en fréquences
|  | $\begin{align*} c_n &:= \frac{1}{2\pi}\int_0^{2\pi} f(t)e^{-\mathbf{i}nt} dt, n \in \mathbb{Z} \quad \text{(complexe)} \\ a_n &:= c_n+c_{-n} \quad \text{(réel)} \\ b_n &:= \mathbf{i}(c_n-c_{-n}) \quad \text{(réel)} \\ f_N(s) &= \frac{a_0}{2} + \sum_{n=1}^N a_n cos(nx) + b_n sin(nx) \end{align*}$ |
| -------- | -------- |
| | |
:::info
$A_n$ amplitude, $\phi_n$ phase de l'harmonique $c_n$ => $\widehat{f}[n]=c_n=A_n e^{\mathbf{i}\phi_n}$ avec $\mathbf{i}=\sqrt{-1}$
$$
f_N(s) = \frac{a_0}{2} + \sum_{n=1}^N a_n cos(nx) + b_n sin(nx) ~ \Leftrightarrow ~ f_N(s) = \frac{A_0}{2} + \sum_{n=1}^N A_n cos(nx - \phi_n)
$$
:::
---
## Exemple
1. Observez $f(x)=cos(x)-sin(x)$
* vérifiez que $f$ est périodique
* Transformer $f$ en $Acos(x - \phi)$, avec $\phi$ un angle
* $A=\sqrt{1^2+(-1)^2}=\sqrt{2}$, et $\phi=\angle(0x, (1,-1)=-\pi/4$
* $f(x)=\sqrt{2}\cos(x + \pi/4))$
2. Même travail avec $3\cos(4x)+\sqrt{3}\sin(4x)$
---
## Test en sage `Cours/fourier.sage`
```python=
cm=sage.plot.colors.get_cmap('viridis')
var('x','n')
f(x)=1-(x/pi-1)^2 # Choix d'une fonction (cloche)
c0 = integral(f(x)/(2*pi),x,0,2*pi)
c(n)= integral(f(x)*exp(-i*n*x)/(2*pi),x,0,2*pi)
a0 =c0 +c0
a(n)=c(n)+c(-n)
b(n)=i*(c(n)-c(-n))
fN(x)=a0/2
F=[fN(x)]
C=[cm(0.)[0:3]]
for k in range(1,10):
fN(x) = fN(x) + a(k)*cos(k*x) + b(k)*sin(k*x)
F.append( fN(x) )
C.append( cm(k/5.0)[0:3] )
plot(f(x),(x,0,2*pi),color="red",thickness=8)+plot(F,(x,0,2*pi),color=C)
```
---
## $f(x)=1-(x/\pi-1)^2$
| | 0 | 1 | 2 | 3 | 4 | 5 | ... |
| ----- | ----- | ---------- | ---------- | ------------ | ------------ | ------------- | --- |
| $a_k$ | $4/3$ | $-4/\pi^2$ | $-1/\pi^2$ | $-4/9/\pi^2$ | $-1/4/\pi^2$ | $-4/25/\pi^2$ | ... |
| $b_k$ | | 0 | 0 | 0 | 0 | 0 | ... |
---
## $f(x) = 1-((x/\pi-1)\cos(2x^2/\pi))^2$
(par intégration numérique)
---
# Transformée de Fourier discrète
Si $f[j]$ **discret** et **périodique** de période $N$, alors sa transformée de Fourier est finie !
:::info
Transformée de Fourier discrète $\mathcal{F}$, $f:\{0,1,\ldots,N-1\} \rightarrow \mathbb{C}$
$$
f \stackrel{\mathcal{F}}{\longrightarrow}\widehat{f}\quad\text{avec}\quad\widehat{f}[k]:=\sum_{j=0}^{N-1}f[j]e^{-\mathbf{i}2\pi\frac{j k}{N}}=\sum_{j=0}^{N-1} f[j](\cos(2\pi\frac{j k}{N}) - \mathbf{i} \sin(2\pi\frac{j k}{N}))
$$
Transformée de Fourier discrète inverse $\mathcal{F}^{-1}$, $\widehat{f}:\{0,1,\ldots,N-1\} \rightarrow \mathbb{C}$
$$
\widehat{f} \stackrel{\mathcal{F}^{-1}}{\longrightarrow}f\quad\text{avec}\quad f[j]:=\frac{1}{N}\sum_{k=0}^{N-1}\widehat{f}[k]e^{\mathbf{i}2\pi\frac{jk}{N}}
$$
:::
---
# Transformée de Fourier discrète 2D
Si $f[j,k]$ **discret** et **périodique** de période $M,N$, alors sa transf. de Fourier est finie !
:::info
Transformée de Fourier discrète $\mathcal{F}$, $f:\{0,\ldots,M-1\} \times \{0,\ldots,N-1\}\rightarrow \mathbb{C}$
$$
f \stackrel{\mathcal{F}}{\longrightarrow}\widehat{f}\quad\text{avec}\quad\widehat{f}[m,n]:=\sum_{j=0}^{M-1}\sum_{k=0}^{N-1}f[j,k]e^{-\mathbf{i}2\pi(\frac{j m}{M}+\frac{k n}{M})}
$$
Transf. de Fourier discrète inverse $\mathcal{F}^{-1}$, $\widehat{f}: \{0,\ldots,M-1\}\times \{0,\ldots,N-1\}\rightarrow \mathbb{C}$
$$
\widehat{f} \stackrel{\mathcal{F}^{-1}}{\longrightarrow}f\quad\text{avec}\quad f[j,k]:=\frac{1}{MN}\sum_{m=0}^{M-1}\sum_{n=0}^{N-1}\widehat{f}[m,n]e^{\mathbf{i}2\pi(\frac{j m}{M}+\frac{k n}{M})}
$$
:::
---
## Exemples de transformée de Fourier discrète 2D
Magnitude seulement, recentrée, vidéo inverse
`video-fourier`
---
## Basses et hautes fréquences
|  |  |
| -------- | -------- |
| contour net : hautes fréquences (HF) | contour flou: moins de HF |
On peut éliminer les hautes fréquences par masque dans l'espace de Fourier
---
## Filtre passe-bas
==filtre passe-bas==: conserve les basses fréquences
| $\widehat{I}=\mathcal{F}(I)$ |  |
| -------- | -------- |
| $\mathcal{F}^{-1}(\widehat{I} \cdot \text{mask})$ |  |
---
## Filtre passe-haut
==filtre passe-haut==: conserve les hautes fréquences
| $\widehat{I}=\mathcal{F}(I)$ |  |
| -------- | -------- |
| $\mathcal{F}^{-1}(\widehat{I} \cdot \text{mask})$ | 
|
---
## Phases et amplitudes sont aussi importants
---
# Propriétés fondamentales
| | Domaine spatial | Domaine fréquentiel |
| ------------------- | -------------------- | ----------------------------------------- |
| notation | $f[j,k]$, $g[j,k]$ | $\widehat{f}[m,n]$, $\widehat{g}[m,n]$ |
| linéarité | $af[j,k]+bg[j,k]$ | $a\widehat{f}[m,n]+b\widehat{g}[m,n]$ |
| convolution/produit | $f[j,k] \ast g[j,k]$ | $\widehat{f}[m,n] \cdot \widehat{g}[m,n]$ |
| différence selon $j$ | $f[j,k] - f[j-1,k]$ | $(1-e^{-\mathbf{i}2\pi\frac{m}{M}})\widehat{f}[m,n]$ |
| différence selon $k$ | $f[j,k] - f[j,k-1]$ | $(1-e^{-\mathbf{i}2\pi\frac{k}{N}})\widehat{f}[m,n]$ |
---
# Algorithme de transformation (lent)
$$
f \stackrel{\mathcal{F}}{\longrightarrow}\widehat{f}\quad\text{avec}\quad\widehat{f}[k]:=\sum_{j=0}^{N-1}f[j]e^{-\mathbf{i}2\pi\frac{j k}{N}}
$$
```cpp=
// T is input array (if spatial domain)
// Y is output array (then frequency domain)
void ft( int n, complex<double> T[n], complex<double> Y[n] )
{ complex<double> i(0,1);
for ( int k = 0; k < n; k++ ) {
complex<double> v( 0, 0 );
for ( int j = 0; j < n; j++ )
v += T[ j ] * exp( -(i*2*pi*j*k)/n);
Y[ k ] = v;
}
}
```
---
# Algorithme de transformation rapide
```cpp=
void fft( int n, complex<double> T[n], complex<double> Y[n] )
{ complex<double> A[n/2], B[n/2], U[n/2], V[n/2];
complex<double> i(0,1), r = 1, w = 2*pi*i/n;
if ( n == 1 ) Y[ 0 ] = T[ 0 ];
else {
for ( auto j = 0; j < n/2; j++ )
A[ j ] = T[ 2*j], B[ j ] = T[ 2*j+1 ];
FFT( n/2, A, U ); FFT( n/2, B, V );
for ( auto j = 0; j < n/2; j++ ) {
Y[ i ] = U[ i ] + r*V[ i ];
Y[ i+n/2 ] = U[ i ] - r*V[ i ];
r *= w;
}
}
} // version 1D de FFT
```
---
# Algorithme de transformation rapide
Si t(n) est le temps d'exécution de FFT(n,T,Y), alors:
$\left\{ \begin{array}{l} t(1) = O(1) \\ t(n) = 2t(n/2) + O(n) \end{array} \right.$
D'où $t(n)=O(n \log n)$ par le "Master theorem"
==Fast Fourier Transform==: $\mathcal{F}$ et $\mathcal{F}^{-1}$ se calcule en $O(MN\log(MN))$
:::success
Il est plus efficace de calculer la convolution dans le domaine de Fourier.
:::
---
# Conclusion
* traitement bas-niveau efficace
* variantes utilisées dans compression JPEG, MPEG
* beaucoup d'évolutions: Transformée en Z, ondelettes
* ondelettes à la base de JPEG2000 et des nouveaux MPEG
* très utilisé en résolution d'Equations aux Dérivées Partielles
* équivalent au filtrage spatial linéaire et invariant par translation
* images couleur par séparation sur les canaux RGB (ou autres)
* donc aussi **limité**
---
# Exercices
$$
f \stackrel{\mathcal{F}}{\longrightarrow}\widehat{f}\quad\text{avec}\quad\widehat{f}[k]:=\sum_{j=0}^{N-1}f[j]e^{-\mathbf{i}2\pi\frac{j k}{N}}, \quad f[j]:=\frac{1}{N}\sum_{k=0}^{N-1}\widehat{f}[k]e^{\mathbf{i}2\pi\frac{jk}{N}}
$$
* Quel est la transformée de $f$ si $f[0]=1$ et $f[j]=0$ si $0<j<N$ ?
* Quel est la transformée de $f$ si $f[m]=1$ et $f[j]=0$ si $0\le j<N$ et $j \neq m$ ?
* Vérifiez que la transformée inverse donne bien le résultat initial
* Vérifiez la *formule de Parseval* sur ces deux exemples, i.e.
$$
\sum_{j=0}^{N-1} |\widehat{f}[j]|^2 = N \sum_{j=0}^{N-1} |f[j]|^2 \quad \text{(conservation de l'énergie)}
$$