Analyse temps fréquence

Le contenu d'un signal évolue en fonction du temps. La transformée de Fourier discrète de l'entièreté du signal permet de connaitre les fréquences contenues dans le signal mais ne permet pas de savoir si ces fréquences ont été présentes tout le temps.

L'analyse temps fréquence permet une analyse de Fourier du signal par intervalle et ainsi de voir l'évolution des fréquences du signal dans le temps.

Vous pouvez trouver plus d'information sur cette page.

Méthode

Le programme suivant permet d'analyser par morceau un signal.

1
"""
2
===========
3
Analyse temps-fréquence
4
Analyse fréquentielle d'un signal par morceau
5
===========
6
"""
7
import numpy as np
8
import matplotlib.pyplot as plt
9
from matplotlib.widgets import Slider, Button, RadioButtons
10
from scipy import signal
11
import soundfile as sf
12
import sounddevice as sd
13
14
import wx
15
16
def update_tobs(val):
17
    tobs = sl_tobs.val
18
    sl_orig.valmax = N / Fe - tobs
19
    sl_orig.ax.set_xlim(0,sl_orig.valmax)
20
    update(0)
21
22
def update_orig(val):
23
    orig = sl_orig.val
24
    sl_tobs.valmax = N / Fe - orig
25
    sl_tobs.ax.set_xlim(0,sl_tobs.valmax)
26
    update(0)
27
28
def update(val):
29
    global choix_fenetre
30
    der_tobs = sl_tobs.val
31
    der_orig = sl_orig.val
32
    k_orig = int(der_orig * Fe)
33
    nb_l = int (der_tobs * Fe)
34
    fenetre = np.zeros(N)
35
    fenetre[k_orig: k_orig + nb_l] = son.max()
36
    courbe_fenetre.set_ydata(fenetre)
37
    if len(son.shape) == 1:
38
        signal = son[k_orig:min(N, k_orig + nb_l)]
39
    else:
40
        signal = son[k_orig:min(N, k_orig + nb_l), :]
41
    freq = np.fft.fftfreq(nb_l)*Fe
42
    #ma_courbe_tempo.set_ydata(signal)
43
    S = np.fft.fft(signal)
44
    ma_courbe_freq.set_data(np.fft.fftshift(freq),np.fft.fftshift(np.abs(S).real/Fe))
45
    fig.canvas.draw_idle()
46
47
def reset(event):
48
    sl_orig.reset()
49
    sl_tobs.reset()
50
51
def change_fenetre(label):
52
    global choix_fenetre
53
    choix_fenetre = label
54
    update(0)
55
56
def annule_extreme(signal, freq_ech, duree_obs):
57
    nb_zero = int((signal.shape[0] / Fe - duree_obs)* Fe)
58
    if nb_zero == 0:
59
        return signal
60
    signal[0:nb_zero//2] = 0
61
    signal[-nb_zero//2:] = 0
62
    return signal
63
64
65
if __name__ == '__main__':
66
    my_app = wx.App()
67
    nom_fichier_son = wx.FileSelector("Fichier son",wildcard="*.wav")
68
    son , Fe = sf.read(nom_fichier_son)
69
    N = son.shape[0]
70
    if len(son.shape) == 2:
71
        son = son[:,0]
72
    print("Nombre d'échantillons : ", N)
73
    print("Fréquence d'échantillonnage : ", Fe)
74
    S = np.fft.fft(son, axis=0)
75
    if len(son.shape) == 1:
76
        nb_courbe = 1
77
    else:
78
        nb_courbe = son.shape[1]
79
    fenetre_active = { 'Rectangulaire':signal.boxcar,
80
                       'Hamming':signal.hamming,
81
                       'blackman':signal.blackman}
82
    fig, ax = plt.subplots(nrows=2, ncols=1)
83
    plt.subplots_adjust(left=0.25, bottom=0.25)
84
    te = np.arange(0,N)/Fe
85
    ma_courbe_tempo,  = ax[0].plot(te,son,label='x(t)', lw=2)
86
    axcolor = 'lightgoldenrodyellow'
87
    ax_orig = plt.axes([0.25, 0.1, 0.65, 0.03], facecolor=axcolor)
88
    ax_tobs = plt.axes([0.25, 0.15, 0.65, 0.03], facecolor=axcolor)
89
    der_orig = 0
90
    der_tobs = N / Fe
91
    sl_tobs = Slider(ax_tobs, 'T_Obs', 0, N/Fe, valinit=der_tobs, valstep=N/ Fe / 20)
92
    sl_orig = Slider(ax_orig, 'Origine', 0., N/Fe-der_tobs, valinit=der_orig, valstep=N/ Fe / 20)
93
    fenetre = np.zeros(N)
94
    fenetre[int(der_orig * Fe):int(der_orig * Fe)+int(der_tobs * Fe)] = son.max()
95
    courbe_fenetre, = ax[0].plot(te,fenetre,color='red',label='fenetre')
96
    ax[0].grid(True)
97
    delta_f = 1.0
98
    ax[0].margins(x=0)
99
    ax[0].set(xlabel='t (s)', ylabel='f(t)) (u.a.)')
100
    freq = np.fft.fftfreq(N)*Fe
101
    S = np.fft.fft(son)
102
    ma_courbe_freq, = ax[1].plot(np.fft.fftshift(freq),np.fft.fftshift(np.abs(S).real/Fe), label=r"$ |X(\nu)| $")
103
    ax[1].grid(True)
104
    ax[1].set(xlabel=r'$ \nu (Hz) $', ylabel='(u.a.)')
105
    ax[0].legend()
106
    ax[1].legend()
107
    sl_orig.on_changed(update_orig)
108
    sl_tobs.on_changed(update_tobs)
109
110
    reset_ax = plt.axes([0.8, 0.025, 0.1, 0.04])
111
    button = Button(reset_ax, 'Reset', color=axcolor, hovercolor='0.975')
112
    button.on_clicked(reset)
113
114
    plt.show()
115