2009-05-31 15 views
2

Le traitement audio est assez nouveau pour moi. Et actuellement en utilisant Python Numpy pour le traitement des fichiers wave. Après avoir calculé la matrice FFT, j'obtiens des valeurs de puissance bruyantes pour des fréquences inexistantes. Je suis intéressé par la visualisation des données et la précision n'est pas une priorité élevée. Existe-t-il un moyen sûr de calculer la valeur d'écrêtage pour supprimer ces valeurs, ou devrais-je utiliser toutes les matrices FFT pour chaque ensemble d'échantillons pour trouver un nombre moyen?Clipping FFT Matrix

concernant

Edit:

from numpy import * 
    import wave 
    import pymedia.audio.sound as sound 
    import time, struct 
    from pylab import ion, plot, draw, show 

    fp = wave.open("500-200f.wav", "rb") 
    sample_rate = fp.getframerate() 
    total_num_samps = fp.getnframes() 
    fft_length = 2048. 
    num_fft = (total_num_samps/fft_length) - 2 
    temp = zeros((num_fft,fft_length), float) 

    for i in range(num_fft): 
     tempb = fp.readframes(fft_length); 
     data = struct.unpack("%dH"%(fft_length), tempb) 
     temp[i,:] = array(data, short) 
    pts = fft_length/2+1 
    data = (abs(fft.rfft(temp, fft_length))/(pts))[:pts] 

    x_axis = arange(pts)*sample_rate*.5/pts 
    spec_range = pts 
    plot(x_axis, data[0]) 
    show() 

Voici l'intrigue à l'échelle non logarithmique, pour le fichier d'onde synthétique contenant 500Hz (fondu) + 200Hz onde sinusoïdale créée en utilisant Goldwave.

+3

avez-vous vérifié vos sorties avec une bonne sortie FFT connue? (matlab ou fftw serait de bonnes sources). En outre, essayez d'entrer des sons purs, c'est-à-direune onde sinusoïdale à une fréquence connue et vérifiez votre sortie pour différentes tailles de fft. – basszero

Répondre

3

Les formes d'onde simulées ne doivent pas afficher les FFT comme votre figure, donc quelque chose est très faux, et probablement pas avec la FFT, mais avec la forme d'onde d'entrée. Le principal problème dans votre intrigue n'est pas les ondulations, mais les harmoniques autour de 1000 Hz, et le sous-harmonique à 500 Hz. Une forme d'onde simulée ne devrait pas montrer cela (par exemple, voir mon graphique ci-dessous). Tout d'abord, vous voudrez probablement essayer simplement de tracer la forme d'onde brute, ce qui indiquera probablement un problème évident. En outre, il semble étrange d'avoir une onde déballée en courts-circuits non signés, c'est-à-dire "H", et surtout après cela, de ne pas avoir une grande composante de fréquence nulle.

J'ai été en mesure d'obtenir une copie assez proche de votre FFT en appliquant un écrêtage à la forme d'onde, comme cela a été suggéré par les harmoniques subharmoniques et supérieures (et Trevor). Vous pourriez introduire un écrêtage dans la simulation ou le déballage. De toute façon, j'ai contourné cela en créant les formes d'onde en numpy pour commencer.

Voici à quoi devrait ressembler la FFT appropriée (c.-à-d.essentiellement parfait, à l'exception de l'élargissement des pics dus au fenêtrage)

alt text http://i43.tinypic.com/1rvsqx.png

est ici un à partir d'une forme d'onde qui a été coupé (et est très similaire à votre FFT, de la sous-harmonique au motif précis de la trois harmoniques supérieures environ 1000 Hz)

alt text http://i44.tinypic.com/mt4avd.png Voici le code que j'utilisé pour générer ces

from numpy import * 
from pylab import ion, plot, draw, show, xlabel, ylabel, figure 

sample_rate = 20000. 
times = arange(0, 10., 1./sample_rate) 
wfm0 = sin(2*pi*200.*times) 
wfm1 = sin(2*pi*500.*times) *(10.-times)/10. 
wfm = wfm0+wfm1 
# int test 
#wfm *= 2**8 
#wfm = wfm.astype(int16) 
#wfm = wfm.astype(float) 
# abs test 
#wfm = abs(wfm) 
# clip test 
#wfm = clip(wfm, -1.2, 1.2) 

fft_length = 5*2048. 
total_num_samps = len(times) 
num_fft = (total_num_samps/fft_length) - 2 
temp = zeros((num_fft,fft_length), float) 

for i in range(num_fft): 
    temp[i,:] = wfm[i*fft_length:(i+1)*fft_length] 
pts = fft_length/2+1 
data = (abs(fft.rfft(temp, fft_length))/(pts))[:pts] 

x_axis = arange(pts)*sample_rate*.5/pts 
spec_range = pts 
plot(x_axis, data[2], linewidth=3) 
xlabel("freq (Hz)") 
ylabel('abs(FFT)') 
show() 
+0

Vous soulevez un bon point que le test le plus rapide et le plus facile de savoir si la FFT est écrêtage, alias ou ajouter du bruit au signal est simplement de prendre la FFT du signal, puis inversez-le et voir si vous récupérez le signal original. Mon intuition est que la FFT fonctionne parfaitement et qu'il s'agit simplement d'une mauvaise interprétation de la sortie de la FFT. – las3rjock

2

FFT parce qu'ils sont fenêtré et sampled la cause aliasing et échantillonnage dans le domaine des fréquences aussi bien. Le filtrage dans le domaine temporel est simplement une multiplication dans le domaine fréquentiel. Vous pouvez donc simplement appliquer un filtre qui ne fait que multiplier chaque fréquence par une valeur pour la fonction du filtre que vous utilisez. Par exemple, multipliez par 1 dans la bande passante et par zéro tous les autres. Les valeurs inattendues sont probablement causées par un aliasing où les fréquences plus élevées sont repliées vers celles que vous voyez. Le original signal needs to be band limited to half your sampling rate ou vous obtiendrez aliasing. Plus inquiétant est le repliement qui déforme la zone d'intérêt car pour cette bande de fréquences vous voulez savoir que la fréquence est celle attendue.

L'autre chose à garder à l'esprit est que lorsque vous saisissez un morceau de données d'un fichier wave, vous le multipliez mathématiquement par une onde carrée. Cela provoque un sinx/x être convoluée avec la réponse en fréquence pour minimiser cela, vous pouvez multiplier le signal fenêtré d'origine avec quelque chose comme Hanning window.

1

Il convient de mentionner une FFT 1D que le le premier élément (index [0]) contient le terme courant continu (fréquence zéro), les éléments contiennent les fréquences positives et les éléments [N/2+1:N-1] contiennent les fréquences négatives. Puisque vous n'avez pas fourni d'échantillon de code ou d'informations supplémentaires sur la sortie de votre FFT, je ne peux pas exclure la possibilité que les "valeurs de puissance bruyantes à des fréquences inexistantes" ne soient pas seulement les fréquences négatives de votre spectre.


EDIT: Here est un exemple d'une FFT radix-2 mis en oeuvre dans le python pur avec une routine de test simple qui trouve la FFT d'une impulsion rectangulaire, [1.,1.,1.,1.,0.,0.,0.,0.]. Vous pouvez exécuter l'exemple sur codepad et voir que la FFT de cette séquence est

[0j,     Negative frequencies 
(1+0.414213562373j), ^
0j,      | 
(1+2.41421356237j),  | 
(4+0j),    <= DC term 
(1-2.41421356237j),  | 
0j,      v 
(1-0.414213562373j)] Positive frequencies 

Notez que les code affiche les coefficients de Fourier par ordre de fréquence ascendante, soit de la plus haute fréquence négative jusqu'à DC, puis jusqu'à la fréquence positive la plus élevée.

+0

Bon point. De même, en numpy, la fonction fftfreq renvoie les cases de fréquence dans l'ordre renvoyé par la fonction fft. – tom10

+0

Merci, semble fonctionner assez bien avec mes vagues calculées, mais j'ai encore des problèmes avec le traitement même des fichiers wave synthétiques, soit j'ai un problème avec ma routine de traitement de fichiers ou l'éditeur d'ondes lui-même. Même s'il n'y en a pas et que c'est la fuite que je ne peux pas contourner, je dois encore trouver un moyen sûr d'exclure des valeurs non significatives. Il semble pire avec l'échelle logarithmique bien sûr. –

1

Je ne sais pas assez de votre question pour réellement répondre à quelque chose de spécifique.

Mais voici quelques choses à essayer de mes propres TFR écrit d'expérience:

  • Assurez-vous que vous suivez la règle Nyquist
  • Si vous regardez la sortie linéaire de la FFT ... vous aura du mal à voir votre propre signal et pense que tout est cassé. Assurez-vous de regarder le dB de votre magnitude FFT. (c'est-à-dire "plot (10 * log10 (abs (fft (x))))")
  • Créez un unitTest pour votre fonction FFT() en alimentant les données générées comme une tonalité pure. Ensuite, alimentez les mêmes données générées à la FFT de Matlab(). Faites une diff valeur absolue entre les deux séries de données de sortie et assurez-vous que la différence de valeur absolue max est quelque chose comme 10^-6 (la seule différence est causée par de petites erreurs à virgule flottante)
  • Assurez-vous que vous êtes windowing your data

Si toutes ces trois choses fonctionnent, alors votre fft est bien. Et vos données d'entrée sont probablement le problème.

apparaît comme des images miroir du signal dans le domaine de fréquence à intervalles réguliers spécifiques avec moins d'amplitude.

+0

"Images en miroir" du signal? Une onde sinusoïdale d'écrêtage (comme illustré) produirait une distorsion harmonique. – endolith