|
|
|
Offset bei der FFT-Zerlegung von Messwerten |
|
| Roachnbua |
Gast
|
 |
Beiträge: ---
|
|
|
|
Anmeldedatum: ---
|
|
|
|
Wohnort: ---
|
|
|
|
Version: ---
|
|
|
|
 |
|
Verfasst am: 13.02.2014, 10:13
Titel: Offset bei der FFT-Zerlegung von Messwerten
|
 |
| |
 |
 |
Hallo zusammen,
habe ein Programm geschrieben zur Zerlegung von Beschleunigungswerten in die verschiedenen Sinus bzw. Cosinusschwingungen, auch in die verschiedenen Ordnungen genannt. Beim Ergebnis ist die 0.Ordnung der Offset. Wenn ich diesen Wert, zusammen mit den anderen Ordnungen, nachher nehme und zur Kontrolle wieder die Funktion bilde, habe ich zwar den gleichen Verlauf der Funktion, nur der Offset stimmt nicht. Die generierte Funktion liegt immer oberhalb der Originalen. Lasse ich den Offsetwert ganz weg, liegt sie darunter. Teile ich den Offsetwert durch zwei, überdecken sich die Funktionen. Kann mir jemand erklären warum das so ist? Das Poblem ist eingerahmt und in fetter Schrift. Ich habe zur Kontrolle der Werte wieder eine Funktion gebildet, mit der Cosinusfunktion. Weiß jemand genau, ob die FFT-Zerlegung auf den Cosinus ausgelegt ist, mit dem Sinus und meinen Ergebnissen stimmt überhaupt nichts.
Danke und schöne Grüße, Kevin
P.S.: Damit ihr euch das besser vorstellen könnt ist im Anhang ein Bild.
%Mein Programm
clear
%-------------------------------------------------------------------------------------------
%PLOTTEN DER BESCHLEUNIGUNGSWERTE DES MITTELBOLZENS (SOLL)
dateipfad = 'Kinematische Berechnung.xlsx';
sheet = 3;
xlRange1 = 'H8:H152';
xlRange2 = 'K8:K152';
xlRange3 = 'B8:B152';
xlRange4 = 'N8:N152';
ax = xlsread(dateipfad, sheet, xlRange1);
ay = xlsread(dateipfad, sheet, xlRange2);
a = xlsread(dateipfad, sheet, xlRange4);
KW = xlsread(dateipfad, sheet, xlRange3);
fig = figure(1);
subplot(3,1,1);
plot(KW,ax);
title('Beschleunigung Ax des Mittelbolzens C SOLL');
xlabel('Kurbelwinkel (°KW)');
ylabel('Beschleunigung Ax (m/s²)');
axis([0 720 -500 1500])
subplot(3,1,2);
plot(KW,ay);
title('Beschleunigung Ay des Mittelbolzens C SOLL');
xlabel('Kurbelwinkel (°KW)');
ylabel('Beschleunigung Ay (m/s²)');
axis([0 720 -1000 1000])
subplot(3,1,3);
plot(KW,a);
title('Beschleunigung Ages des Mittelbolzens C SOLL');
xlabel('Kurbelwinkel (°KW)');
ylabel('Beschleunigung Ages (m/s²)');
axis([0 720 0 1500])
%----------------------------------------------------------------------------------------------
%ERMITTLUNG ALLER ORDNUNGEN MIT DER FOURIERZERLEGUNG
fa = 144; %Abtastfrequenz = 1 / Abtastperiode = 1 / (t(2)-t(1)) (Beispiel)
N1 = length(ax);
Fft_Analyse1 = fft(ax, N1)/N1;
X = 2*abs(Fft_Analyse1(1:N1/2+1));
N2 = length(ay);
Fft_Analyse2 = fft(ay, N2)/N2;
Y = 2*abs(Fft_Analyse2(1:N2/2+1));
fig = figure(fig + 1);
Ordnungen1 = 0:0.5:(size(X,1))/2-0.5;
subplot(2,1,1);
plot(Ordnungen1, X, 'rx', 'LineWidth', 2.0)
title('Anregungsordnungen')
xlabel('Ordnung')
ylabel('Amplitude')
axis([0 10 0 (max(X)+50)]) %in die Grafik reinzoomen, so dass nur die ersten 10 Ordnungen dargestellt werden
Ordnungen2 = 0:0.5:(size(Y,1))/2-0.5;
subplot(2,1,2);
plot(Ordnungen2, Y, 'rx', 'LineWidth', 2.0)
%title('Anregungsordnungen')
xlabel('Ordnung')
ylabel('Amplitude')
axis([0 10 0 (max(Y)+50)]) %in die Grafik reinzoomen, so dass nur die ersten 10 Ordnungen dargestellt werden
%--------------------------------------------------------------------------------------------------------------
%AUTOMATISCHES FINDEN DER 1. UND 2.ORDNUNG MIT OFFSET
%Amplitudenwerte von ax finden
Schwellwert = 0.050; %Schwellwert von 5%
MaxValueX = max(X); %maximale Amplitude im Spektrum suchen
SchwellwertX = Schwellwert * MaxValueX; %Schwellwert berechnen
% Maskierung der Werte mit Schwellwert -> nur noch Werte über dem Schwellwert vorhanden
maskeX = (X > SchwellwertX);
X_mask = maskeX.* X;
% Zusammenfassen von Frequenz (Ordnungen) und Amplitude in einem Vektor
X_result(:,1) = X_mask; %Amplituden in Spalte 1 speichern
X_result(:,2) = Ordnungen1; %Ordnugnen in Spalte 2 speichern eventuell ist hier noch eine Transponierung der Matrix nötig. Also X_result(:,2) = Ordnungen';
X_result(X_result(:,1)==0, =[];
%Amplitudenwerte von ay finden
MaxValueY = max(Y); %maximale Amplitude im Spektrum suchen
SchwellwertY = Schwellwert * MaxValueY; %Schwellwert berechnen
% Maskierung der Werte mit Schwellwert -> nur noch Werte über dem Schwellwert vorhanden
maskeY = (Y > SchwellwertY);
Y_mask = maskeY.* Y;
% Zusammenfassen von Frequenz (Ordnungen) und Amplitude in einem Vektor
Y_result(:,1) = Y_mask; %Amplituden in Spalte 1 speichern
Y_result(:,2) = Ordnungen2'; %Ordnugnen in Spalte 2 speichern eventuell ist hier noch eine Transponierung der Matrix nötig. Also X_result(:,2) = Ordnungen';
Y_result(Y_result(:,1)==0,=[];
%-------------------------------------------------------------------------------------------------------
%PHASENVERSCHIEBUNG ALLER ORDNUNGEN ERMITTELN
FFT_Analyse(:,1) = fft(ax);
Phasenwinkel(:,1) = abs(angle(FFT_Analyse(:,1)));
FFT_Analyse(:,2) = fft(ay);
Phasenwinkel(:,2) = abs(angle(FFT_Analyse(:,2)));
%----------------------------------------------------------------------------------------
%PLOTTEN DER ZERLEGTEN FUNKTIONEN
%Eintrag der Ergebnisse in die Result-Datei
Result(:,1) = Ordnungen1;
for i = [1:73]
Result(i,2) = X_mask(i);
end
for i = [1:73]
Result(i,3) = X(i);
end
for i = [1:73]
Result(i,4) = Phasenwinkel(i,1);
end
Result(:,5) = Ordnungen2;
for i = [1:73]
Result(i,6) = Y_mask(i);
end
for i = [1:73]
Result(i,7) = Y(i);
end
for i = [1:73]
Result(i, = Phasenwinkel(i,2);
end
%Abbild der FFT zerlegten Funktione zur Kontrolle
Omega = 157.0796; %Drehzahl der Kurbelwelle --> entspricht 1500 U/min
dateipfad = 'Kinematische Berechnung.xlsx';
sheet = 2;
xlRange1 = 'A8:A152';
t = xlsread(dateipfad, sheet, xlRange1);
%Korrekturfaktor der 0.Ordnung (Offset)
Ordnung_Ax = Result(1,3)/2;
Ordnung_Ay = Result(1,7)/2;
%Funktion ax mit nur 1. und 2. Ordnungen
for j = [2:73]
for k = [1:145]
Ax(j,k) = Result(j,2)*cos(Result(j,1)*Omega*t(k,1)+Result(j,4));
end
end
%Funktion ax mit allen Ordnungen
for j = [2:73]
for k = [1:145]
Ax_ges(j,k) = Result(j,3)*cos(Result(j,1)*Omega*t(k,1)+Result(j,4));
end
end
%Funktion ay mit nur 1. und 2. Ordnungen
for j = [2:73]
for k = [1:145]
Ay(j,k) = Result(j,6)*cos(Result(j,5)*Omega*t(k,1)+Result(j,);
end
end
%Funktion ay mit allen Ordnungen
for j = [2:73]
for k = [1:145]
Ay_ges(j,k) = Result(j,7)*cos(Result(j,5)*Omega*t(k,1)+Result(j,);
end
end
%-------------------------------------------------------------------------
%Hier liegt das Problem
Funktion_x(:,1) = Ordnung_Ax + sum(Ax(1:73,);
Funktion_x_ges(:,1) = sum(Ax_ges(1:73,);
Funktion_y(:,1) = Ordnung_Ay + sum(Ay(1:73,);
Funktion_y_ges(:,1) = sum(Ay_ges(1:73,);
%-------------------------------------------------------------------------
%Ploten der Funktionen ax und ay Vergleich
fig = figure( fig + 1);
subplot(2,1,1);
plot(KW,Funktion_x,'red',KW,Funktion_x_ges,'green',KW,ax,'black');
title('Beschleunigung Ax im Mittelgelenk C Vergleich');
xlabel('Kurbelwinkel (°KW)');
ylabel('Beschleunigung Ax [m/s²]');
axis([0 720 -500 1500]);
subplot(2,1,2);
plot(KW,Funktion_y,'red',KW,Funktion_y_ges,'green',KW,ay,'black');
title('Beschleunigung Ay im Mittelgelenk C Vergleich');
xlabel('Kurbelwinkel (°KW)');
ylabel('Beschleunigung Ay [m/s²]');
axis([0 720 -1000 1000]);
|
|
|
|
|
|
| DSP |
Forum-Meister
|
|
Beiträge: 2.117
|
|
|
|
Anmeldedatum: 28.02.11
|
|
|
|
Wohnort: ---
|
|
|
|
Version: R2014b
|
|
|
|
|
|
|
Verfasst am: 13.02.2014, 12:19
Titel:
|
|
Immer wieder gern gemachter Fehler bei der Skalierung.
Der Gleichsignalanteil Y(1) sowie bei der Nyquistfreq. Y(N/2+1) darf eben NICHT mit 2 multipliziert werden. Die FFT berechnet ein zweiseitiges Spektrum. Die beiden genannten Signalanteile kommen dort aber nur 1 mal vor...nicht wie der Rest doppelt.
http://www.gomatlab.de/viewtopic.ph.....9bbf0f67bdfd45810ed79fc3a
|
|
|
|
| Roachnbua |
Gast
|
|
Beiträge: ---
|
|
|
|
Anmeldedatum: ---
|
|
|
|
Wohnort: ---
|
|
|
|
Version: ---
|
|
|
|
|
|
|
Verfasst am: 14.02.2014, 11:13
Titel:
|
|
Schon mal danke für die Antwort. Die habe ich auch verstanden, nur die Essenz daraus noch nicht wirklich.
Ich habe mir auch den anderen Code angeschaut, aber irgendwie blicke ich da noch nicht ganz durch. Wenn ich die Multiplikation 2 weg lasse, stimmen meine Werte dazwischen drin nicht mehr und sonst die zwei Randwerte, also Y(1) und Y(N/2+1). Stimmt meine Berechnung wenn ich sozusagen einfach Y(1) durch zwei teile und aufsummiere? Kann ich meinen Code so belassen und sicher sein das er stimmt?
|
|
|
|
| DSP |
Forum-Meister
|
|
Beiträge: 2.117
|
|
|
|
Anmeldedatum: 28.02.11
|
|
|
|
Wohnort: ---
|
|
|
|
Version: R2014b
|
|
|
|
|
|
|
Verfasst am: 14.02.2014, 12:16
Titel:
|
|
So wäre die Skalierung der Amplituden richtig:
|
|
|
|
|
|
|
Einstellungen und Berechtigungen
|
|
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
|
|
Impressum
| Nutzungsbedingungen
| Datenschutz
| FAQ
|  RSS
Hosted by:
Copyright © 2007 - 2024
goMatlab.de | Dies ist keine offizielle Website der Firma The Mathworks
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.
|
|
