Discussion :

[Svart26]

Bonjour,

En complément du topic traitant de la création du .wav monofréquence, je passe actuellement à l'étape clé de la FFT. Grace à l'aide de oiffrig , la création du fichier .wav se fait de la sorte:

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
import math, numarray as NA, wave, numarray as NA2
 
freq = 200.
freqech = 44100
duree = 10.
omega = 2 * math.pi * freq
son = NA.arange(0., duree, 1./freqech, NA.Float64)
son = NA.sin(omega * son)
son *= 32768.
son = son.astype(NA.Int16)
 
 
wavfile = wave.open('AutreAmplitude200.wav','w')
wavfile.setparams((1, 2, freqech , len(son), 'NONE', 'not compressed'))
wavfile.writeframes(son.tostring())
wavfile.close()
 
 
Valeurs = son[0:222]

Dans la varable Valeur, je sauvegarde les valeurs prises par la fonction sinus (avec f = 200Hz) sur 1 période environ.
En appliquant la fonction fft issue du module numarray, je devrais donc obtenir quelque chose me disant qu'effectivement j'ai un signal de l'ordre de 200Hz.

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
from numarray import fft
print fft.fft(Valeurs,450).real

Or le résultat obtenu est le suivant:

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
3.51000000e+02   3.02301824e+06   2.29319478e+05  -1.96406314e+06
  -1.23762282e+04  -4.48999295e+05  -1.02429368e+04  -2.03945581e+05
  -9.59045653e+03  -1.15423303e+05  -9.26517011e+03  -7.30000882e+04
  -9.02794648e+03  -4.93465665e+04  -8.82407680e+03  -3.47937512e+04
  -8.63322297e+03  -2.52527407e+04  -8.44808818e+03  -1.86411046e+04
  -8.24596904e+03  -1.39443011e+04  -8.04285242e+03  -1.04709028e+04
  -7.82704600e+03  -7.87677142e+03  -7.59881077e+03  -5.88953325e+03
  -7.36427447e+03  -4.36747628e+03  -7.11766236e+03  -3.17780670e+03
ETC...

La première valeur correspond à la composante continue du signal, autrement sa valeur moyenne, OK
Mais Quid des autres valeurs? Amplitudes de quelles harmoniques, fréquences(cos, sin)? Pourquoi des valeurs négatives?

Si vous avez des suggestions....

[Svart26]

Bon, j'ai debusqué ça sur le net:

Results of FFT are in the form of complex vector of the same length as
data. Starting from position 0 it gives you an constant (DC factor), position
1 an amplitude and phase (remember - complex number gives you both amplitude
and phase) of wave of the length table/2 and so on. If you take real value of
this, you discard part of the information. AFAIR - taking real part gives you
sine component, taking imaginary part gives you cosine component, taking
absolute value gives you total amplitude and taking angle gives you phase

Autrement, pour récuperer l'amplitude des sinusoidale, il suffit de prendre le module du comple correspondant, mais quand à connaître la fréquence qui lui correspond...

[Svart26]

A titre de simplication:

la fft de [ 1. 0. 1. 0.] avec 8 samples donne:
[ 2.+0.j 1.-1.j 0.+0.j 1.+1.j 2.+0.j 1.-1.j 0.+0.j 1.+1.j]


la fft de [ 1. 0. 1. 0.] avec 4 samples donne:
[2 0 2 0]

Si quelqu'un comprend le résultat...

[parp1]

Et bien oui c'est simple:

Avec 8 samples tu as donc 8 echantillons
Avex 4 Samples tu en as donc 4 echantillons. (oui je sais que jusqu'ici je sais traduire.)

Mais lorsque tu regardes ton resultat avec huit sample et que tu prends un echantillon sur deux.

Ca te donne :

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

[2.+0.j 1.-1.j 0.+0.j 1.+1.j 2.+0.j 1.-1.j 0.+0.j 1.+1.j]

Ce qui donne:

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

2+0j,0+0j,2+0j,0+0j ou plus simplement 2,0,2,0

[Svart26]

Et bien oui c'est simple:

Avec 8 samples tu as donc 8 echantillons
Avex 4 Samples tu en as donc 4 echantillons. (oui je sais que jusqu'ici je sais traduire.)

Mais lorsque tu regardes ton resultat avec huit sample et que tu prends un echantillon sur deux.

Ca te donne :

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

[2.+0.j 1.-1.j 0.+0.j 1.+1.j 2.+0.j 1.-1.j 0.+0.j 1.+1.j]

Ce qui donne:

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

2+0j,0+0j,2+0j,0+0j ou plus simplement 2,0,2,0

Oui, tout à fait. Mais plus précisément, je cherche à savoir à quelles fréquences correspondent ces amplitudes

[fabrice-102]

Pour récupérer les fréquence, sous numpy, il suffit de faire :
freq = fftfreq(nbr,d=eps)


avec nbr le nombre d'échantillons et eps l'inverse du nombre d'échantillons récupérés par seconde. Pour vous : 1.0/44100 (ne pas oublie 1.0 pour travailler avec une précision suffisante !!)

Vous pouvez ensuite, par exemple, calculer la fréquence dont l'énergie est la plus importante (soit pour une sinusoïde, la fréquence de la sinusoïde) :
maxfreq = max(freq)

et vous pouvez ainsi calculer l'énergie correspondante :
fft = fft(Valeurs,len(Valeurs))
energie = fft[argmax(freq)]

[Svart26]

Pour récupérer les fréquence, sous numpy, il suffit de faire :
freq = fftfreq(nbr,d=eps)


avec nbr le nombre d'échantillons et eps l'inverse du nombre d'échantillons récupérés par seconde. Pour vous : 1.0/44100 (ne pas oublie 1.0 pour travailler avec une précision suffisante !!)

Vous pouvez ensuite, par exemple, calculer la fréquence dont l'énergie est la plus importante (soit pour une sinusoïde, la fréquence de la sinusoïde) :
maxfreq = max(freq)

et vous pouvez ainsi calculer l'énergie correspondante :
fft = fft(Valeurs,len(Valeurs))
energie = fft[argmax(freq)]

Merci beaucoup, tout paraît évident après lecture de votre réponse!

[Svart26]

Pour récupérer les fréquence, sous numpy, il suffit de faire :
freq = fftfreq(nbr,d=eps)


avec nbr le nombre d'échantillons et eps l'inverse du nombre d'échantillons récupérés par seconde. Pour vous : 1.0/44100 (ne pas oublie 1.0 pour travailler avec une précision suffisante !!)

Vous pouvez ensuite, par exemple, calculer la fréquence dont l'énergie est la plus importante (soit pour une sinusoïde, la fréquence de la sinusoïde) :
maxfreq = max(freq)

et vous pouvez ainsi calculer l'énergie correspondante :
fft = fft(Valeurs,len(Valeurs))
energie = fft[argmax(freq)]

Bien, donc petit bilan des essais que j'ai effectué aujourd'hui, histoire de donner un retour des choses:

Effectivement, fftfreq renvoie les fréquences, et je dois dire que cette fonction fait vraiment plaisir.

Mon besoin étant de déterminer la/les fréquences où l'amplitude est maximale, je réalise les étapes suivantes:

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
10
11
 
FFT = fft.fft(valeurs, 2048)
FFTabs = []
for i in range(50):
        FFTabs.append(abs(FFT[i])) #calcul de la valeur absolu des amplitudes (*) (pour les 50 premiers echantillons ici)
 
 
maxi = max(FFTabs)
indice = FFTabs.index(maxi)
freq = numpy.fft.fftfreq(2048,1.0/44100)
print 'frequence avec amplitude max: ', freq[indice]

(*)d'apres ce lien, l'amplitude totale de chacune des frequences obtenues lors de la decomposition est obtenue en prenant la valeur absolue du complexe . http://mail.python.org/pipermail/tut...st/040263.html

Bref, d'après mes premiers essais, tout est K, j'espère que ca ne changera pas...
Merci

[Matthieu Brucher]

A la place de numarray, commence par utiliser numpy complètement, ça résoudra peut-être quelques soucis !
Ensuite, si le retour est complexe, tu prends la valeur absolue de tes termes, et avec l'échelle qui t'a été donné, c'est bon.

Autre détail, la FFT fournie avec FFTW ne possède pas de terme correcteur, cf mon tuto : http://miles.developpez.com/tutoriels/algo/fft/

[Svart26]

A la place de numarray, commence par utiliser numpy complètement, ça résoudra peut-être quelques soucis !
Ensuite, si le retour est complexe, tu prends la valeur absolue de tes termes, et avec l'échelle qui t'a été donné, c'est bon.

Autre détail, la FFT fournie avec FFTW ne possède pas de terme correcteur, cf mon tuto : http://miles.developpez.com/tutoriels/algo/fft/

Sympathique tutoriel. Sinon oui en effet, avec l'amplitude "valeur absolue" du complexe et l'echelle de frequence c'est parfait, le tout complété par un graphique sous Tkinter est le résultat est ma foi fort sympathique. (même si les valeurs des amplitudes rendues par FFt() me tourmentent un peu).
Bref, demain je jetterai également un oeil à FFTW, Thanks.