Das ist ja alles richtig und mir auch bekannt. Meine Info über die Phase bezog sich auch nicht auf die Berechnung deines Zeitsignals sondern lediglich darauf, dass er nicht diesen Phasenwert von
zum Schieben verwenden kann . Er bekommt zwar eine Phasenverschiebung, weiß dann aber nicht ob die Grundschwingung ein Sinus oder Cosinus war (mal vom Bild abgesehen).
Ich habe es aber gerade nochmal getestet und so geht es doch...egal ob Sinus oder Cosinus als Grundschwingung. Beim Cosinus gibt es dann aber eine Phasenverschiebung von 90 Grad. Aber die Grundwelle beginnt immer beim Nulldurchlauf.
ich verstehe weiterhin nicht was das Ziel mit dieser Phase sein soll...
Wenn ich den Code um die Phase wie von Dir vorgeschlagen erweitere ergeben sich verschobene Zeitverläufe, aber das Signal beginnt auch nicht bei null. Vielleicht kannst Du mal das gesamte Beispiel posten?
Die Phase gibt doch genau an, wie das Signal verschoben ist. Wenn ich sie kenne, kann ich das Signal auch so manipulieren, dass die Phasenverschiebung weg ist. Ob es Sinn macht, ist ja eine ganze andere Frage.
Es geht so aber doch nicht ganz. Ich habe es nur kurz getestet und mich auf Grund des Offsets des Testsignals täuschen lassen. Wenn man den Offset wie hier weglässt, treffe ich nicht exakt den Nulldurchgang. Obwohl meine Phasenverschiebung der Grundschwingung f1 60 Grad ist, beträgt sie mit dem Befehl angle() über 60 Grad Außerdem muss die Berechnung doch für Sinus und Cosinus als Testsignal unterschiedlich sein(wenn man die Phase manipulieren will). Diese Info hat man ja aber in der Regel nicht.
Code:
%% Signal
dt = 1/8000;
zeit = 0:dt:0.1;
f1 = 100;
f2 = 500;
Ampl1 = 10;
Ampl2 = 5;
offset = 0;
phase1 = (pi/3);
phase2 = 0;
signal = offset+(Ampl1*cos(2*pi*f1*zeit + phase1) + Ampl2*cos(2*pi*f2*zeit + phase2));
% FFT Berechnen
Ts = diff(zeit(1:2)); % Abtastzeit
N = length(zeit); % Länge des Datenvektors
f = [0:floor((N-1)/2)] / (N*Ts); % Frequenzvektor für Spektrum-Plot % Fouriertransformierte
SIGNAL = fft(signal); % Fouriertransformation
SIGNAL = SIGNAL / N; % Normierung
SIGNAL = [SIGNAL(1)2*SIGNAL(2:floor((N-1) / 2) + 1)]; % begrenzen auf f_max
% Fourierkoeffizienten a und b bestimmen for i=1:1:length(SIGNAL)
a(i) = real(SIGNAL(i));
b(i) = -imag(SIGNAL(i));
end;
% Frequenzvektor sortieren so dass später das Signal aus den Hauptfrequenzanteilen % wieder erzeugt bzw. approximiert werden kann [wert_f position_f] = sort(abs(SIGNAL));
[phase position_f2] = sort(angle(SIGNAL));
for i=1:1:length(f)
f_(i) = f(position_f(length(position_f)-i+1));
a_(i) = a(position_f(length(position_f)-i+1));
b_(i) = b(position_f(length(position_f)-i+1));
phase_(i) = phase(position_f2(length(position_f2)-i+1));
end
% Signal auf Grund der FFT-Werte erzeugen bis zur k-ten Harmonischen
k = 1;
a = a_;
b = b_;
phase = phase_;
% Signal-Berechnung auf Grund eines Filters
Anzahl_mw = 64;
b(1:Anzahl_mw) = 1/64;
signal_filt = filter(b,1,signal);
%% Plot der Berechnungen figure(1) subplot(311) plot(zeit,signal,'r','LineWidth',2);
xlabel('Zeit [s]');
ylabel('Signal');
grid on
subplot(312) plot(f,abs(SIGNAL),'rx','LineWidth',2,'Markersize',8);
xlabel('Frequenz [Hz]');
ylabel('Absolutwert');
grid on
subplot(313) plot(f,angledim(angle(SIGNAL),'radians','degrees'),'r.-','LineWidth',2,'Markersize',8);
xlabel('Frequenz [Hz]');
ylabel('Phase [Grad]');
grid on
figure(2) hold on
plot(zeit,signal,'r','LineWidth',2);
plot(zeit,signal_FFT,'b.-','LineWidth',2);
plot(zeit,signal_filt,'m','LineWidth',2);
hold off
legend('original','FFT','Filter',0) xlabel('Zeit [s]');
ylabel('Signal');
grid on
wenn Du den Zeitvektor auf 0,01 begrenzt erhältst Du "genau" eine Periode der Schwingung. Dann siehst Du in Deinem Phasen-Plot bei 100 Hz eine Phase von 1,07 (= pi/3 mit Unschärfe) und bei 500 Hz hast Du eine Phase von Null (wieder nicht exakt aber ziemlich gut).
Wenn Du aber nun das zeitliche Signal über einen längeren Zeitraum betrachtest kannst Du das so nicht mehr herauslesen. Dann verschmiert sich das Ergebnis genau so wie es bei Deinem Plot zu sehen ist.
aso, zur Info, Dein Plot ist nicht in Grad sondern in Radians
Das mit der Unschärfe wusste ich nicht...wieder was dazu gelernt. Wenn das Testsignal jetzt aber ein Sinus ist, geht es mit dem Addieren der Phase so nicht mehr. Dann habe ich keinen Nulldurchgang bei t=0.
Ich habe hier aber wohl ohnehin einen Denkfehler drin, da in den Koeff. a und b ja schon die Info über die Phase steckt. Also müsste man eigentlich a und b manipulieren. Aber es gibt ja auch folgende Beziehung...
Wenn mich nun nur die Amplitude und Frequenz der Grundschwinung interessiert, kann ich ja auch die Phase weglassen. Da er ja aber bei t=0 einen Nulldurchgang haben möchte, müsste ich die Phase = pi/2 setzen...da es ja ein cos ist.
Warum sollten das nicht Grad sein? Ich rechne doch mit angledim in Grad um, da angle mir den Winkel im Bogenmaß liefert.
Nun geht es bei cos und sin als Testsignal (auch mit längerer Messzeit). Ich hatte beim sortieren von Betrag und Phase noch einen Fehler. Jetzt kommt allerdings immer ein Sinus raus...
Code:
%% Signal clearall;
dt = 1/8000;
zeit = 0:dt:0.1;
f1 = 100;
f2 = 500;
Ampl1 = 10;
Ampl2 = 5;
offset = 0;
phase1 = (pi/3);
phase2 = 0;
signal = offset+(Ampl1*cos(2*pi*f1*zeit + phase1) + Ampl2*cos(2*pi*f2*zeit + phase2));
% FFT Berechnen
Ts = diff(zeit(1:2)); % Abtastzeit
N = length(zeit); % Länge des Datenvektors
f = [0:floor((N-1)/2)] / (N*Ts); % Frequenzvektor für Spektrum-Plot % Fouriertransformierte
SIGNAL = fft(signal); % Fouriertransformation
SIGNAL = SIGNAL / N; % Normierung
SIGNAL = [SIGNAL(1)2*SIGNAL(2:floor((N-1) / 2) + 1)]; % begrenzen auf f_max
% Fourierkoeffizienten a und b bestimmen for i=1:1:length(SIGNAL)
a(i) = real(SIGNAL(i));
b(i) = -imag(SIGNAL(i));
end;
% Frequenzvektor sortieren so dass später das Signal aus den Hauptfrequenzanteilen % wieder erzeugt bzw. approximiert werden kann [wert_f position_f] = sort(abs(SIGNAL));
[phase position_f2] = sort(angle(SIGNAL));
for i=1:1:length(f)
f_(i) = f(position_f(length(position_f)-i+1));
a_(i) = a(position_f(length(position_f)-i+1));
b_(i) = b(position_f(length(position_f)-i+1));
phase_(i) = phase(length(position_f2)-i+1);
betrag_(i) = wert_f(length(position_f)-i+1);
end
% Signal auf Grund der FFT-Werte erzeugen bis zur k-ten Harmonischen
k = 1;
a = a_;
b = b_;
phase = phase_;
betrag = betrag_;
% Signal-Berechnung auf Grund eines Filters
Anzahl_mw = 64;
b(1:Anzahl_mw) = 1/64;
%h_sinc = window_sinc_filter(100,400,8000,1) ;
signal_filt = filter(b,1,signal);
%% Plot der Berechnungen figure(1) subplot(311) plot(zeit,signal,'r','LineWidth',2);
xlabel('Zeit [s]');
ylabel('Signal');
grid on
subplot(312) plot(f,abs(SIGNAL),'rx','LineWidth',2,'Markersize',8);
xlabel('Frequenz [Hz]');
ylabel('Absolutwert');
grid on
subplot(313) plot(f,angledim(angle(SIGNAL),'radians','degrees'),'r.-','LineWidth',2,'Markersize',8);
xlabel('Frequenz [Hz]');
ylabel('Phase [Grad]');
grid on
figure(2) hold on
plot(zeit,signal,'r','LineWidth',2);
plot(zeit,signal_FFT,'b.-','LineWidth',2);
plot(zeit,signal_filt,'m','LineWidth',2);
%plot(zeit(1:length(signal)-32),signal_filt(1,33:length(signal)),'k--','LineWidth',2);
hold off
legend('original','FFT','Filter',0) xlabel('Zeit [s]');
ylabel('Signal');
grid on
pause closeall
Der von DSP gepostete Code scheint sehr gut zu funktionieren
Mittlerweile beschäftigt mich noch eine andere Fragestellung.
Hierfür habe ich euch mal das .mat-File des verrauschten, sinusförmigen Messsignals angehängt.
Es ist im Plot deutlich zu erkennen, dass das genannte Signal eine Phasenverschiebung besitzt.
Ich frage mich nun, ob diese Phasenverschiebung durch horizontales Verschieben des Signals auf der Zeitachse eliminiert werden kann.
Hier stößt man jedoch auf die folgenden Probleme:
1.) Wie ermittelt man eigentlich die vorhandene Phasenverschiebung des Signals, die es zu eliminieren gilt?
2.) Wie kann ich das Signal anschließend so auf der Zeitachse verschieben, dass ein Signal entsteht, welches:
a) bei x=0 seinen positiven Nulldurchgang hat und ansonsten unverändert bleibt?
b) dem Signal auch a) um 90° voreilt (Cosinus)?
Ich vermute, dass sich diese Aufgabe ebenfalls im Frequenzbereich lösen lässt indem man z.B. nach der fft jede Komponente mit exp(i*p) multipliziert (wobei p die in 1.) ermittelte Phasenverschiebung ist).
Die Implementierung will mir aber bisher nicht gelingen.
Habt ihr vielleicht eine gute Idee?
Ich hatte doch noch einen Fehler bei der Phase, der aber jetzt behoben ist. Ehrlich gesagt verstehe ich deine Fragen nicht. Ich habe nun dein beigefügtes Signal getestet und die Plots entsprechend angepasst. Die Phasenverschiebung der Grundwelle deines Zeisignals wird in der Legende angezeigt. Bei diesem Signal betrifft die Phasenverschiebung -140 Grad.
Da ich für die Änderungen aber die Phase weglasse erhältst du doch, ein Signal mit positivem Nulldurchgang bei t=0...was nur ein Sinus sein kann.
Warum willst du jetzt noch schieben? Es ist doch auch gegenüber dem Original in Amplitude und Frequenz identisch...nur die Phase ändere ich.
Die Phasenverschiebung berechnet man aus den Fourierkoeffizienten a und b = Real- und Imaginärteil so.
phi = arctan(b/a) -> phi in [rad]
in Matlab nutzt man einfach die Funtkion angle() oder um es wie hier angegeben zu berechnen, atan2(). Ich gebe dir im Command Window auch noch die Zeitverschiebung delta_t = (phi + 90 Grad)/freq. aus. Das Signal muss um diesen Wert nach links auf der Zeitachse verschoben werden, um keine Phasenverschiebung zu haben und somit deine Forderung erfüllt zu haben.
Teste doch einfach nochmal mit dem Testsignal rum, damit es klar wird. Der Cosinus eilt dem Sinus um 90 Grad voraus. Nehmen wir mein Psp. mit den pi/3 = 60 Grad Phasenverschiebung bei f1 an. Ist das Testsignal ein Sinus, erhältst du theoretisch -30 Grad Phase (real auf Grund von Ungenauigkeiten nicht ganz) aus der oben benannten Berechnung. Wenn du zu den -30 nun 90 Grad Verschiebung zum Cosinuns dazurechnest, erhältst du wieder die ursprunglichen 60 Grad deines Testsignals. Ist das Testsignal degegen ein Cosinus, wirst du gleich aus der Berechnung die 60 Grad erhalten.
%signal = offset+(Ampl1*sin(2*pi*f1*zeit + phase1) + Ampl2*sin(2*pi*f2*zeit + phase2));
% reales Signal
signal = y';
zeit = x;
% FFT Berechnen
Ts = diff(zeit(1:2)); % Abtastzeit
N = length(zeit); % Länge des Datenvektors
f = [0:floor((N-1)/2)] / (N*Ts); % Frequenzvektor für Spektrum-Plot % Fouriertransformierte
SIGNAL = fft(signal); % Fouriertransformation
SIGNAL = SIGNAL / N; % Normierung
SIGNAL = [SIGNAL(1)2*SIGNAL(2:floor((N-1) / 2) + 1)]; % begrenzen auf f_max
% Fourierkoeffizienten a und b bestimmen for i=1:1:length(SIGNAL)
a(i) = real(SIGNAL(i));
b(i) = -imag(SIGNAL(i));
end;
% Frequenzvektor sortieren so dass später das Signal aus den Hauptfrequenzanteilen % wieder erzeugt bzw. approximiert werden kann [wert_f position_f] = sort(abs(SIGNAL));
phase = angle(SIGNAL);
for i=1:1:length(f)
f_(i) = f(position_f(length(position_f)-i+1));
a_(i) = a(position_f(length(position_f)-i+1));
b_(i) = b(position_f(length(position_f)-i+1));
phase_(i) = phase(position_f(length(position_f)-i+1));
betrag_(i) = wert_f(length(position_f)-i+1);
end
% Signal auf Grund der FFT-Werte erzeugen bis zur k-ten Harmonischen
k = 1;
a = a_;
b = b_;
phase = phase_;
betrag = betrag_;
%% Plot der Berechnungen figure(1) subplot(311) plot(zeit,signal,'r','LineWidth',2);
xlabel('Zeit [s]');
ylabel('Signal y(t)');
title('Zeitsignal','FontSize',18);
grid on
subplot(312) plot(f,abs(SIGNAL),'rx','LineWidth',2,'Markersize',8);
xlabel('Frequenz [Hz]');
ylabel('Absolutwert');
title('Frequenzspektrum - Amplitude','FontSize',18);
grid on
subplot(313) plot(f,angledim(angle(SIGNAL),'radians','degrees'),'r.-','LineWidth',2,'Markersize',8);
xlabel('Frequenz [Hz]');
ylabel('Phase [Grad]');
title('Frequenzspektrum - Phase','FontSize',18);
grid on
figure(2) hold on
plot(zeit,signal,'r','LineWidth',2);
plot(zeit,signal_FFT,'b.-','LineWidth',1);
plot(zeit,signal_FFT_Phasenaenderung,'g--','LineWidth',4);
hold off
title('Ermittlung der Grundwelle des Signals','FontSize',18);
h=sprintf(['orig. Grundwelle mit Phase = ' num2str(angledim(phase(k),'radians','degrees'),'%3.1f') ' Grad'],'FontSize',14);
legend('original',h,'Grundwelle mit Phasenänderung');
xlabel('Zeit [s]');
ylabel('Signal y(t)');
grid on
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.