Verfasst am: 30.01.2015, 12:44
Titel: Bestimmung der Grenzfrequenz aus FFT
Hallo zusammen,
ich habe folgendes Problem, ich habe ein Signal, Gierrate eines Fahrversuchs (siehe Anhang 'Plot_yaw_rate'), und würde dieses nun gern filtern.
Ich suche aktuell nach der Grenzfrequenz, dich ich beim Filter verwenden will.
Dazu hab ich mir das Betragsspektrum (wiederum siehe Anhang 'Frequenzbereich1' und 'Frequenzbereich2') mit Hilfe dieses Skriptes
http://www.gomatlab.de/signalverarb.....t,fftbetragsspektrum.html
ausgegeben.
Mein Input:
- signal = yaw_rate;
- nfft = 512;
- fa = 100; Aufzeichnungsrate des Messinstruments (ich hoffe, das ist die hier gemeinte Abtastrate??)
- scale = 0;
Wie kann ich jetzt hieraus feststellen, was mein Nutzsignal und was mein Störsignal (Rauschen der Sensoren o. ä.) ist?
Wie bestimme ich die Grenzfrequenz, die ich dann im Filter verwende?
Ich hoffe, ich hab das Ganze einigermaßen verständlich beschrieben.
Schon mal vielen Dank,
viele Grüße
Thomas
Wie kann ich jetzt hieraus feststellen, was mein Nutzsignal und was mein Störsignal (Rauschen der Sensoren o. ä.) ist?
Wie bestimme ich die Grenzfrequenz, die ich dann im Filter verwende?
Diese Entscheidung kann man am besten durch einen Vergleich der Spektren von Störsignal und Nutzsignal+Störsignal treffen. Die Messwertaufzeichnung beinhaltet in deinem Fall beides, weshalb ich die Messwerte getrennt habe. Wie du in figure 2 sehen kannst, erreicht das Störsignal das gleiche Amplitudenlevel wie das Nutzsignal bei ca. 8 Hz. Das Filter sollte also spätestens ab 8 Hz dämpfen. Da die verschiedenen Filtertypen (Butterworth, Chebyshev und Bessel) aber einen unterschiedlichen Verlauf der Dämpfung haben, ist es auch hier ratsam, sich das gefilterte Ergebnis im Zeit- und Frequenzbereich anzeigen und mit den ungefilterten Messwerten zu vergleichen.
Das Bessel Filter ist recht gut für Rechtecksignal-Formen und filtert die Störung hier am deutlichsten Weg. Die Frage ist allerdings ob dadurch evtl. auch schon Informationen des Nutzsignal mit weggefiltert werden. Aktuell sind im unteren Code die anderen Filtertypen auskommentiert. Nutze jeweils ein anderes Filter und schau dir die Ergebnisse in beiden Bereichen an.
Es gibt keine klare Regel wie man die Grenzfrequenz in deinem Fall festlegt. Ein Rechtecksignal hat nämlich über den gesamten Frequenzbereich Signalanteile vorliegen. Somit überlappen sich Stör- und Rechtecksignale immer. Je mehr ich vom hochfrequenten Bereich abschneide, desto mehr runden sich die Ecken des Rechtecksignals. Somit kommt es immer zu einem Informationsverlust. Hier muss man dann einfach einen Kompromiss zwischen Dämpfung der Störung und möglichsten wenig vom Nutzsignal abschneiden. Hier kommt es auf Erfahrung an und so ein Vergleich macht es deutlich einfacher.
Code:
%% get Data from figure file
h_data = openfig('Plot_yaw_rate.fig');
axesObjs = get(h_data, 'Children'); %axes handles
dataObjs = get(axesObjs, 'Children'); %handles to low-level graphics objects
xdata = get(dataObjs, 'xdata');
ydata = get(dataObjs, 'ydata');
%% Bessel IIR Tiefpass Filter
fc = 20; % Cut-off Frequenz in Hz
order = 5; % Ordnung [z,p,k] = besself(order,fc); % Bessel analog filter design [zd,pd,kd] = bilinear(z,p,k,fs); % Analog to digital mapping
sos = zp2sos(zd,pd,kd); % Convert to SOS form [b,a] = sos2tf(sos); % convert to num and den coefficents
%% Chebyshev Typ2 IIR Tiefpass Filter
fc = 5; % Cut-off Frequenz in Hz
order = 11; % Ordnung
damping = 40; % Dämpfung im Sperrbereich in dB
%[b,a] = cheby2(order,damping,2*fc/fs,'low');
%% Signal filtern -> bidirektional, damit das Filterdelay verschwindet und somit gefilteres und ungefilteres Signal übereinander liegen
filtsignal = filtfilt(b, a, nutzsignal); % filter with bessel IIR
%% Darstellung im Zeitbereich figure(1) plot(t,nutzsignal,'b',t,stoersignal,'r',t,filtsignal,'k');
xlabel('Zeit in s') ylabel('Amplitude') legend('Nutzsignal','Stoersignal','Nutzsignal gefiltert') grid on;
%% Darstellung im Frequenzbereich [~, Amplitude_dB, fv] = FFT_betragsspektrum( nutzsignal, N, fs, 1, 0);
figure(2) plot(fv,Amplitude_dB,'b','Linewidth',2);
xlabel('Frequenz in Hz') ylabel('Amplitude dB') grid on;
function[mag, mag_dB, fv] = FFT_betragsspektrum( signal, nfft, fa, window,scale) % Input: % Signal im Zeitbereich % nfft = Anzahl Messwerte für fft % wenn nfft > length(sig) -> fft(sig,nfft) führt Zeropadding durch % fa = Abtastfreq. % scale: 0 = keine Impulsantwort als Eingang, 1 = Impulsantwort % Output: % Magnitude des Spektrums linear und dB skaliert % Frequenzvektor fv in [Hz] von 0...fa/2
% un-,gerade Anzahl Messwerte? ifmod(nfft,2) == 0;
k = (nfft/2) + 1;
else
nfft = nfft + 1;
k = (nfft/2) + 1;
end
Du kannst Beiträge in dieses Forum schreiben. Du kannst auf Beiträge in diesem Forum antworten. Du kannst deine Beiträge in diesem Forum nicht bearbeiten. Du kannst deine Beiträge in diesem Forum nicht löschen. Du kannst an Umfragen in diesem Forum nicht mitmachen. Du kannst Dateien in diesem Forum posten Du kannst Dateien in diesem Forum herunterladen
MATLAB, Simulink, Stateflow, Handle Graphics, Real-Time Workshop, SimBiology, SimHydraulics, SimEvents, and xPC TargetBox are registered trademarks and The MathWorks, the L-shaped membrane logo, and Embedded MATLAB are trademarks of The MathWorks, Inc.
RSS