Matlab-Low-Pass-filter mit fft

Ich implementiert einen einfachen low-pass-filter in matlab mit Hilfe einer vorwärts-und rückwärts-fft.
Es funktioniert im Prinzip, aber die minimum-und maximum-Werte vom original unterscheiden.

signal = data;
%% fourier spectrum
% number of elements in fft
NFFT = 1024;
% fft of data
Y = fft(signal,NFFT)/L;
% plot(freq_spectrum)

%% apply filter
fullw = zeros(1, numel(Y));
fullw( 1 : 20 ) = 1;
filteredData = Y.*fullw;

%% invers fft
iY = ifft(filteredData,NFFT);
% amplitude is in abs part
fY = abs(iY);
% use only the length of the original data
fY = fY(1:numel(signal));
filteredSignal = fY * NFFT; % correct maximum

clf; hold on;
plot(signal, 'g-')
plot(filteredSignal ,'b-')
hold off;

das resultierende Bild sieht aus wie diese

Matlab-Low-Pass-filter mit fft

Was mache ich falsch? Wenn ich Normalisiere beide Daten die gefilterte signal korrekt aussieht.

  • Ihr filter muss symmetrisch sein, wie das signal ist. warum erwarten Sie, dass min und max nicht zu ändern? es gibt keinen Grund zu.
  • Beachten Sie, dass man versucht, eine "ziegelwand" - filter, wie dies in der Frequenz-Domäne produzieren unschöne Artefakte - Sie brauchen, um eine glatte Funktion in den Frequenzbereich (in der Regel ein-Fenster-Funktion). Beachten Sie auch, dass Ihre filter-gain ist nicht normiert, und wie @thang Noten der filter sollte symmetrisch sein, sonst hat man komplexe time-domain-output-Daten.
InformationsquelleAutor Matthias Pospiech | 2015-03-02
  1. 17

    Nur um uns zu erinnern, wie MATLAB speichert die Frequenzen für Y = fft(y,N):

    • Y(1) ist der Konstante offset
    • Y(2:N/2 + 1) ist die Menge der positiven Frequenzen
    • Y(N/2 + 2:end) ist die Menge der negativen Frequenzen... (normalerweise würden wir dieses Grundstück Links der vertikalen Achse)

    Um eine echte low-pass-filter, die wir erhalten müssen sowohl die niedrigen positiven Frequenzen und die niedrigen negative Frequenzen.

    Hier ist ein Beispiel, dies zu tun mit einem multiplikativen Rechteck-filter in der Frequenz-Domäne, wie Sie getan haben: 

    % make our noisy function
t = linspace(1,10,1024);
x = -(t-5).^2  + 2;
y = awgn(x,0.5); 
Y = fft(y,1024);

r = 20; % range of frequencies we want to preserve

rectangle = zeros(size(Y));
rectangle(1:r+1) = 1;               % preserve low +ve frequencies
y_half = ifft(Y.*rectangle,1024);   % +ve low-pass filtered signal
rectangle(end-r+1:end) = 1;         % preserve low -ve frequencies
y_rect = ifft(Y.*rectangle,1024);   % full low-pass filtered signal

hold on;
plot(t,y,'g--'); plot(t,x,'k','LineWidth',2); plot(t,y_half,'b','LineWidth',2); plot(t,y_rect,'r','LineWidth',2);
legend('noisy signal','true signal','+ve low-pass','full low-pass','Location','southwest')

    Matlab-Low-Pass-filter mit fft

    Den vollen low-pass-fitler macht einen besseren job, aber Sie werden feststellen, dass der Wiederaufbau ist ein bisschen "wellig". Dies ist da die Multiplikation mit einem Rechteck-Funktion in der Frequenz-Domäne ist das gleiche wie ein Faltung mit einer sinc-Funktion im Zeitbereich. Faltung mit einer sinc-Funktion ersetzt jeden Punkt mit einer sehr ungleichen gewichteten Durchschnitt seiner Nachbarn, also die "wave" - Effekt.

    Einem Gauß-filter hat schönere low-pass-filter-Eigenschaften, da die fourier-Transformation einer Gauß-eine Gauß. Eine Gauß-Zerfälle auf null schön, so dass es nicht weit entfernten Nachbarn, die in den gewichteten Durchschnitt der während der Faltung. Hier ist ein Beispiel mit einem Gauß-filter Erhaltung der positiven und negativen Frequenzen:

    gauss = zeros(size(Y));
sigma = 8;                           % just a guess for a range of ~20
gauss(1:r+1) = exp(-(1:r+1).^ 2 / (2 * sigma ^ 2));  % +ve frequencies
gauss(end-r+1:end) = fliplr(gauss(2:r+1));           % -ve frequencies
y_gauss = ifft(Y.*gauss,1024);

hold on;
plot(t,x,'k','LineWidth',2); plot(t,y_rect,'r','LineWidth',2); plot(t,y_gauss,'c','LineWidth',2);
legend('true signal','full low-pass','gaussian','Location','southwest')

    Matlab-Low-Pass-filter mit fft

    Wie Sie sehen können, ist die Rekonstruktion wird auf diese Weise viel besser.

    • Sehr schöne Erklärung. Ich gebe zu, dass ich dazu Neige zu vergessen, dass die fourier-Raum ist gefaltet und filtern muss getan werden, auf beiden Seiten.
    InformationsquelleAutor eigenchris
Schreibe einen Kommentar