Forum: Digitale Signalverarbeitung / DSP / Machine Learning Fourieranalyse eines Signals


von So H. (tronix)



Lesenswert?

Hallo!

Ich habe ein Spannungssignal über einem Winkel vorliegen, welches sich 
2pi-periodisch wiederholt (signal_2pi.jpg) und habe davon 
Fourierkoeffizienten berechnet (fourier_signal_2pi.jpg), soweit so gut.

Nun möchte ich von diesem Signal nur ein Sechstel von 2pi betrachten 
(signal_2pi_6tel.jpg) und von diesem Teilstück ebenfalls die 
Fourierkoeff. berechnen. Ich hab das auch gemacht, aber es scheint etwas 
mit dem Spektrum nicht zu stimmen (fourier_signal_2pi_6tel.jpg).

Berechnet wurden die Fourierkoeffizienten folgendermaßen:

https://de.wikipedia.org/wiki/Fourierreihe#Darstellung_in_Sinus-Kosinus-Form

Fällt jemand ad hoc ein, was da nicht stimmen könnte?

: Bearbeitet durch User
von Bernd (Gast)


Lesenswert?

> Fällt jemand ad hoc ein, was da nicht stimmen könnte?

Ja - die Bilder fehlen.

von So H. (tronix)


Lesenswert?

Sind schon da, sorry!

von C Programmierer (Gast)


Lesenswert?

Ich würde mal behaupten, dass du einen Fehler in der 
Sprektrumsberechnung des pi/3 Signals hast.

von So H. (tronix)


Lesenswert?

Es wird natürlich jeweils die selbe Funktion zur Berechnung des 
Spektrums verwendet, und zwar folgende (Matlab):
1
function [a0 ck] = fouriercoe(os, x_werte, y_werte)
2
    
3
    signal_period = max(x_werte)-min(x_werte);
4
5
    a0 = ...
6
            2/signal_period ...                      
7
            * trapz(    x_werte, ...              
8
                        y_werte); 
9
    
10
    ck = [];
11
12
    for k = 1:os
13
    ak = ...
14
            2/signal_period ...
15
            * trapz(    x_werte, ...              
16
                        y_werte ...
17
                        .* cos(k * x_werte));
18
    bk = ...
19
            2/signal_period ...
20
            * trapz(    x_werte, ...              
21
                        y_werte ...
22
                        .* sin(k * x_werte));
23
24
    ck = [ck (ak+1i*bk)];
25
    end    
26
end

von C Programmierer (Gast)


Lesenswert?

Wofür steht denn jeweils das "..."? Wenn du schon Matlab benutzt, dann 
nimm doch einfach die DFT.

spektrum = fft(y_werte)

Wenn ich mich richtig erinnere, musst du das noch durch die Anzahl 
teilen.

spektrum = fft(y_werte) / length(y_werte)

von C Programmierer (Gast)


Lesenswert?

Wenn das ganze mit dem 2*pi Signal funktioniert, vermute ich, dass du 
beim pi/3 Signal falsche Parameter an die Funktion übergibst.

von anonymous (Gast)


Lesenswert?

C Programmierer schrieb:
> Wofür steht denn jeweils das "..."?

damit kann man in matlab lange befehle umbrechen in die nächste zeile

von Frank (Gast)


Lesenswert?


von C Programmierer (Gast)


Lesenswert?

@Frank: Daran liegts nicht. Hier wird eine Fourier Reihe von einem 
ungefähren Sinussignal über ziemlich genau die Periode des Sinus 
berechnet. Das ideale Spektrum müsste zumindest grob erkennbar sein. (So 
wie bei dem Spektrum mit der 2pi Periode). Hier ist es aber komplett 
anders. Ich denke, dass bei der Parameterübergabe etwas nicht stimmt.

von C Programmierer (Gast)


Lesenswert?

Für mich sieht das Spektrum sehr SI-Förmig aus. Das würde man wohl 
besser erkennen, wenn es nicht logarithmisch dargestellt wäre. Statt 
eines Sinus wurde wohl ausversehen ein Rechteck als y-Werte übergeben.

von Lalala (Gast)


Lesenswert?

Mal ein clear all eingeben. Oder das ganze Programm posten

von So H. (tronix)


Lesenswert?

Der Grund, warum ich nicht die DFT von Matlab nehmen kann ist, dass 
Spannung und Winkel über der Zeit aufgenommen werden, dargestellt und 
verarbeitet wird das Signal aber über den Winkel, wodurch die Samples 
nicht mehr ganz äquidistant sind...

Hier der Funktionsaufruf mit einem Signal aus der Retorte (in echt sind 
das  Messdaten), da jedoch das gleiche rauskommt sollte es aber auch 
passen (bzw. halt irgendwas nicht passen):
1
x_werte = 0:(2*pi/6)/514:2*pi/6;
2
y_werte = sin(x_werte*6);
3
[ao ck] = fouriercoeffs(60,x_werte,y_werte);

von Dumdi D. (dumdidum)


Lesenswert?

So H. schrieb:
> Hier der Funktionsaufruf mit einem Signal aus der Retorte (

ist das das funktionierende oder das nicht funktionierende Signal. bei 
den wenigen Infomationwn tippe ich mal, das Du die falschen k Werte 
nimmst (durch N teilen?) und so den Peak verpasst.

von Student (Gast)


Lesenswert?

Seit wann ist Fourieranalyse denn bei nicht-äquidistanter Abtastung 
möglich?

von F. F. (foldi)


Lesenswert?

Hallo ihr Fourier Spezies, gibt es ein Programm wo man nur die 
Frequenzen eingibt und das Programm dann entsprechende Berechnungen 
durchführt, die dann in C-Code übernommen werden können?

Ich habe leider nur Fachoberschulreife und das übersteigt meine 
mathematischen Fähigkeiten. Leider kann mein Sohn mir auch nicht mehr 
helfen.

Es geht eigentlich nur darum, bestimmte Frequenzen heraus zu rechnen, 
also um digitale Filter (Gehörschutz).

Will hier auch nicht den Thread kapern und bin dann gleich wieder weg.

Danke!

von So H. (tronix)


Lesenswert?

Dumdi D. schrieb:
> ist das das funktionierende oder das nicht funktionierende Signal. bei
> den wenigen Infomationwn tippe ich mal, das Du die falschen k Werte
> nimmst (durch N teilen?) und so den Peak verpasst.

Das nicht funktionierende, siehe x_werte, die ich nur bis 2*pi/6 laufen 
hab lassen.
Wie meinst du das genau, die k-Werte müssen doch ganzzahlig sein, oder?

Student schrieb:
> Seit wann ist Fourieranalyse denn bei nicht-äquidistanter Abtastung
> möglich?

http://de.mathworks.com/help/signal/examples/spectral-analysis-of-nonuniformly-sampled-signals.html

von Frank (Gast)


Lesenswert?

So H. schrieb:
> Der Grund, warum ich nicht die DFT von Matlab nehmen kann ist, dass
> Spannung und Winkel über der Zeit aufgenommen werden, dargestellt und
> verarbeitet wird das Signal aber über den Winkel, wodurch die Samples
> nicht mehr ganz äquidistant sind...

Du könntest dir Test halber mal äquidistante Punkte rechnen in damit die 
FFT füttern.

F. F. schrieb:
> Will hier auch nicht den Thread kapern und bin dann gleich wieder weg.

Dann mach es auch nicht ;-)
Mach einen anderen Artikel auf.
Nebenbei verbessere noch (im neuen Artikel) die Beschreibung.
Ich hab nicht kapiert was du willst.

von So H. (tronix)


Angehängte Dateien:

Lesenswert?

Frank schrieb:
> Du könntest dir Test halber mal äquidistante Punkte rechnen in damit die
> FFT füttern.

Anbei mein Testsignal, welches nicht funktioniert (eingangssignal.jpg) 
und daraus die berechneten Spektren, einmal mit der Funktion 
fouriercoeffs.m berechnet (spektrum_eigenbau_fourier.jpg) und einmal mit 
der Funktion fft() berechnet (spektrum_fft.jpg), sowie nachfolgend der 
Code:
1
x_werte = 0:(2*pi/6)/514:2*pi/6;
2
y_werte = sin(x_werte*6);
3
4
figure(1)
5
plot(x_werte,y_werte)
6
grid on;
7
8
n_fourier = length(y_werte)
9
10
[ao ck_cum] = fouriercoeffs(n_fourier,x_werte,y_werte);
11
12
figure(2)
13
plot(abs(ck_cum))
14
grid on;
15
16
NFFT=n_fourier;
17
figure(3)
18
y_fft = fft(y_werte);
19
plot(abs(y_fft(1:NFFT/2)*2/NFFT));
20
grid on

So richtig schlau werd ich daraus aber nicht wirklich...

von Frank (Gast)


Lesenswert?

Wieviele Punkte sind denn später jeweils in Figure2 und Figure3?

Für mich sieht es so aus als ob in Figure3 gefühlt 1000 mal weniger 
Punkte sind.

von So H. (tronix)


Lesenswert?

Die Spektren wurden jeweils mit 514 Punkten berechnet.

Es sieht irgendwie so aus, als ob die Funktion fouriercoeffs.m mit 
Ausschnitten die kürzer sind als 2*pi nicht umgehen kann, vielfache von 
2*pi sind keine Problem, d.h. ich kann 100 2*pi-Perioden reinschicken, 
und das Spektrum ändert sich (wie erwartet) nicht, aber bei Bruchteilen 
von 2*pi ists vorbei.

: Bearbeitet durch User
von Lalala (Gast)


Lesenswert?

Ersetzt mal k* xwerte durch
N=len (xwerte)
2*pi*k*(0:N-1)/N

Natuerlich mit der korrekten Matlab Syntax.

von So H. (tronix)


Lesenswert?

Lalala schrieb:
> Ersetzt mal k* xwerte durch
> N=len (xwerte)
> 2*pi*k*(0:N-1)/N
>
> Natuerlich mit der korrekten Matlab Syntax.

Das war des Rätsel Lösung, Vielen Dank!

Bitte melde dich an um einen Beitrag zu schreiben. Anmeldung ist kostenlos und dauert nur eine Minute.
Bestehender Account
Schon ein Account bei Google/GoogleMail? Keine Anmeldung erforderlich!
Mit Google-Account einloggen
Noch kein Account? Hier anmelden.