Table des matières Python

Autocorrélation d'un signal

1. Introduction

La fonction d'autocorrélation est un outil très important pour l'analyse des signaux, particulièrement pour les signaux aléatoires ([1]). Elle est notamment utilisée pour l'analyse des données dans les simulations de Monte-Carlo.

Ce document montre comment calculer une fonction d'autocorrélation pour un signal discret. On verra aussi comment elle peut être calculée à partir de la densité spectrale de puissance, en utilisant le théorème de Wiener-Khinchine.

2. Fonction d'autocorrélation

2.a. Définition

Soit x(t) un signal. La fonction d'autocorrélation temporelle ([1]) est définie par :

Il s'agit donc de la moyenne temporelle du produit du signal par lui-même décalé d'un temps τ. La fonction d'autocorrélation est paire; on peut donc l'étudier pour τ>0.

Les signaux réels sont limités dans le temps. L'intégrale définissant la fonction d'autocorrélation est alors calculée pour une durée T finie assez grande.

Soit xk un signal numérique obtenu par échantillonnage avec une période Te. Soit M le nombre de points utilisés pour calculer la moyenne (T=MTe). La fonction d'autocorrélation discrète est définie par :

La moyenne est ainsi calculée pour les points d'indices i à i+M-1. Soit N le nombre de points de la fonction d'autocorrélation. L'indice n varie de 0 à N-1.

2.b. Fonction de calcul

import numpy
                    

La fonction suivante effectue le calcul de l'autocorrélation pour un signal donné par un tableau numpy x. L'indice de départ i doit être supérieur à N.

def autocorrel(x,N,i,M):
    C = numpy.zeros(N)
    for k in range(i,i+M):
        for n in range(N):
            C[n] += x[k]*x[k-n]
    return C/M        
                    

2.c. Exemple

On considère comme exemple un signal aléatoire obtenu en ajoutant un bruit gaussien à un signal périodique :

import math
import random
from matplotlib.pyplot import *

def signal(t):
    return math.sin(2*math.pi*t)+0.5*math.sin(4*math.pi*t)+0.3*random.gauss(0,0.8)
Ns = 4000
T = 8.0
dt = T/Ns
t = numpy.zeros(Ns)
x = numpy.zeros(Ns)
for k in range(Ns):
    t[k] = k*dt
    x[k] = signal(t[k])
figure(figsize=(12,4))
plot(t,x)
xlabel('t')
ylabel('x')
grid()
    
                    
figAfigA.pdf

Voyons son autocorrélation calculée sur une durée 1000Te :

N = 1000
temps = numpy.arange(N)*dt
C = autocorrel(x,N,N,Ns-N)
figure(figsize=(12,4))
plot(temps,C)
xlabel("t")
ylabel("C")
                    
figBfigB.pdf

Comme on le voit sur cet exemple, la fonction d'autocorrélation permet d'obtenir la période d'un signal, même lorsqu'il est fortement perturbé par le bruit.

Par ailleurs, la valeur de l'autocorrélation pour τ=0 est :

print(C[0])
--> 0.6816976624914507

Pour un signal de valeur moyenne nulle, il s'agit de la fluctuation quadratique moyenne :

La fonction d'autocorrélation est très utile pour analyser un signal sonore :

import scipy.io.wavfile as wave
rate,data = wave.read('../../tfd/spectreson3/piano_la3.wav')
n = data.size
duree = 1.0*n/rate
dt = 1.0/rate
Ns = 4000
x1 = data[0:Ns]
t1 = numpy.arange(Ns)*dt
figure()
plot(t1,x1)
xlabel("t (s)")
                    
figCfigC.pdf
N=1000
C = autocorrel(x1,N,N,Ns-N)
figure(figsize=(12,4))
plot(numpy.arange(N)*dt,C)
xlabel("t")
ylabel("C")
                     
figDfigD.pdf

L'instant du premier maximum donne la période du signal.

3. Théorème de Wiener-Khinchine

3.a. Énoncé

La transformée de Fourier du signal x(t) est :

Le signal s'exprime en fonction de sa transformée de Fourier par la transformée de Fourier inverse :

La densité spectrale de puissance du signal est définie par :

Le théorème de Wiener-Khinchine ([2]) énonce que la densité spectrale de puissance est la transformée de Fourier de la fonction d'autocorrélation :

Inversement, la fonction d'autocorrélation est la transformée de Fourier inverse de la densité spectrale :

3.b. Exemple

On calcule la transformée de Fourier discrète du signal échantillonné précédent, puis son spectre de puissance.

import numpy.fft
tfd = numpy.fft.fft(x)
freq = numpy.arange(Ns)*1.0/T
S = numpy.square(numpy.absolute(tfd))/Ns
figure(figsize=(8,5))
plot(freq,S)
xlabel("f")
ylabel("S")
axis([0,20,0,S.max()])
                    
figEfigE.pdf

La transformée de Fourier inverse permet d'obtenir la fonction d'autocorrélation :

c = numpy.fft.ifft(S)
figure(figsize=(12,4))
plot(t,c)
xlabel("t")
ylabel("C")
                    
figFfigF.pdf
Références
[1]  E. Tisserand, J.F. Pautex, P. Schweitzer,  Analyse et traitement des signaux,  (Dunod, 2008)
[2]  F. Reif,  Statistical and thermal physics,  (McGraw-Hill, 1965)
Creative Commons LicenseTextes et figures sont mis à disposition sous contrat Creative Commons.