WICHTIG: Der Betrieb von goMatlab.de wird privat finanziert fortgesetzt. - Mehr Infos...

Mein MATLAB Forum - goMatlab.de

Mein MATLAB Forum

  
Gast > Registrieren       Autologin?   

Partner:




Forum
      Option
[Erweitert]
  • Diese Seite per Mail weiterempfehlen
     


Gehe zu:  
Neues Thema eröffnen Neue Antwort erstellen

Probleme bei der FFT eines Rechtecksignals
 

Phtagen

Gast

Beiträge: ---
Anmeldedatum: ---
Wohnort: ---
Version: ---
     Beitrag Verfasst am: 04.03.2016, 10:58     Titel: Probleme bei der FFT eines Rechtecksignals
  Antworten mit Zitat      
Hallo liebes Forum,

ich beschäftige mich gerade mit den Grundlagen den Signalverarbeitung und bin auf ein Problem gestoßen. Zwar benutze ich nicht Matlab sondern Octave 3.8.2, aber ich denke ich bin hier trotzdem richtig aufgehoben.

Ich erzeuge mir mit rectpuls(…) ein Rechtecksignal im Intervall [-Pi, Pi] und möchte anschließend hiervon eine FFT durchführen (Details siehe Code).

Das Amplitudenspektrum welches ich erzeuge jedoch, wirft einige Fragen bei mir auf.

1) Ich erwarte als FFT eines Rechtecksignals eigentlich die sinc-Funktion. Das Spektrum sieht jedoch anders aus.
2) Die Höhe der Amplitude kann ich mir nicht erklären. Ich hätte bei einem Einheitsrechteckimpuls auch eine Amplitude von 1 erwartet. Habe ich hier evtl. falsch normiert?
3) Wie verfahre ich in diesem konkreten Fall mit dem Gleichanteil bei omega = 0?
4) Die Frage ist wieder eher Octave spezifisch und von untergeordnetem Interesse. Beim Plotten mit Gnuplot wird die Achsenskalierung am linken Bildrand abgeschnitten. Was kann ich dagegen tun?

Ok, das wären soweit meine drei (vier) Problemchen. Wäre echt froh, wenn mir da jemand helfen könnte. Ach, nebenbei noch: Bitte nicht wundern, ich habe den Frequenzvektor als Kreisfrequenz formuliert, da das in dem Beispiel was ich gerade nachvollziehe auch so gemacht wird und ich die Vergleichbarkeit der Ergebnisse gewärhrleisten wollte.

Ok, hier noch der Code. Viele Grüße

Phtagen

Code:

octave:1> t=linspace(-pi,pi,10000);    %Samplelänge
octave:2> sp=length(t);                %Anzahl Samplepoints
octave:3> T=2*pi/sp;                   %Periodendauer
octave:4> fs=1/T;                      %Abtastfrequenz
octave:6> w=0.4;                       %Fensterbreite
octave:8> y=rectpuls(t,w);             %Rechteckpuls
octave:9> Y=fftshift(fft(y));          %FFT durchführen, zurechtrücken
octave:10> P2=2*abs(Y)/sp;              %Amplitudenspektrum normiert
octave:12> f2=linspace(-fs/2*2*pi, fs/2*2*pi, sp);   %Frequenzvektor
octave:16> plot(f2(4800:5200),P2(4800:5200));      %Ergebnis plotten
 


RechteckSpek.pdf
 Beschreibung:

Download
 Dateiname:  RechteckSpek.pdf
 Dateigröße:  15.31 KB
 Heruntergeladen:  462 mal
RechteckSpek.pdf
 Beschreibung:

Download
 Dateiname:  RechteckSpek.pdf
 Dateigröße:  15.31 KB
 Heruntergeladen:  430 mal


DSP
Forum-Meister

Forum-Meister


Beiträge: 2.117
Anmeldedatum: 28.02.11
Wohnort: ---
Version: R2014b
     Beitrag Verfasst am: 04.03.2016, 11:40     Titel:
  Antworten mit Zitat      
Hallo,

ich habe zwar kein Octave...wenn aber die Umsetzung der FFT wie in Matalb erfolgt, ist deine Skalierung auf jeden Fall falsch. Da du ein zweiseitiges Spektrum erstellst, darf der Betrag nicht mit 2 multipliziert werden. Das gilt nur für ein einseitges Spektrum.

Siehe: http://www.gomatlab.de/fft-umfassendes-beispiel-t777.html

Code:

N = ...; % Anzahl Messwerte

% Berechnung der FFT
% ------------------
H = fft(y, N);
% Berechnung des Amplitudengangs aus dem komplexen Frequenzvektor H:
amplH = abs(H);
% Amplitudenskalierung (Normierung auf N) und verschieben der Elemente des
% Amplitudenvektors, so dass die Darstellung des Amplitudengangs von -fn...0...fn
% erfolgen kann:
amplitudengang = fftshift(amplH/N);
 


Der Verlauf ist für mich nach einem kurzen Blick mal abgesehen von der Normierung einleuchtend.

Gruß DSP
Private Nachricht senden Benutzer-Profile anzeigen
 
Phtagen

Gast

Beiträge: ---
Anmeldedatum: ---
Wohnort: ---
Version: ---
     Beitrag Verfasst am: 04.03.2016, 15:45     Titel:
  Antworten mit Zitat      
Hallo DSP,

danke für deine schnelle Antwort.

Bist du dir mit der Skalierung sicher? Warum mache ich das nur bei dem einseitigen Spektrum. Ich hätte gedacht, dass das einseitige Spektrum (in der Matlab-Hilfe mit P1 bezeichnet) tatsächlich nur die erste Hälfte des P2-Spektrums ist: http://de.mathworks.com/help/matlab/ref/fft.html?!?!?!? Nichts desto trotz, meine Amplitude beträgt ja nur irgendwas um 0.14. Würde ich das noch halbieren wäre ich bei ca. 0.7. Das ist von meiner erwarteten Amplitude von 1 immer noch sehr weit weg. Oder habe ich einfach falsche Erwartungen Shocked ?

Und zum Verlauf. In dieser Internetquelle https://www.akustik.tu-berlin.de/fi.....alverarbeitung_gesamt.pdf auf S. 47 findet sich einmal die Darstellung des Amplitudenspektrums eines Rechteckimpulses als sinc-Funktion und darunter ein ähnlicher Verlauf wie ich ihn erhalten habe. Kann mir jemand erklären, worin der Unterschied besteht und wieso mein Verlauf so aussieht wie er aussieht?

Vielen Dank und viele Grüße

Phtagen
 
DSP
Forum-Meister

Forum-Meister


Beiträge: 2.117
Anmeldedatum: 28.02.11
Wohnort: ---
Version: R2014b
     Beitrag Verfasst am: 05.03.2016, 13:24     Titel:
  Antworten mit Zitat      
Hallo,

bei der Skalierung bin ich mir 100% sicher. Die FFT Doku ist hier schon seit Jahren nicht ganz korrekt. Erstens fehlt eine Unterscheidung, ob das Signal ein Impuls ist oder nicht. Der fouriertransformierte Impuls hat über alle Frequenzen den Betrag = 1. Diese Ergebnis erhält man aus dem Matlab Beispiel aber nicht. Zur Darstellung des positiven Freq.-bereichs ist die Doku auch nicht korrekt, da der Betrag von Gleichsignalanteil und Nyquistfrequenz ebenfalls falsch skaliert wird.

Wie du die Amplituden richtig skalierst, findest du in dem Link aus meinem vorherigen Post.

Nun zu der Frage der Amplitudehöhe. Das Beispiel in dem Buch beschreibt ein Rechteckfenster (über alle t -> y = 1) und kein Rechtecksignal wie in deinem Beispiel. Für das Rechteckfenster mit der Amplitude 1 ergibt sich auch im Betragsspektrum eine 1 für die Amplitude. Teile ich nun das Rechtecksignal gleichmäßig in High und Low auf, ergibt sich nur noch 0.5 als Amplitude im Betragsspektrum. Die Abweichung wird höher, je kleiner der High Anteil ist.

Gruß DSP
Private Nachricht senden Benutzer-Profile anzeigen
 
Phtagen

Gast

Beiträge: ---
Anmeldedatum: ---
Wohnort: ---
Version: ---
     Beitrag Verfasst am: 06.03.2016, 16:28     Titel:
  Antworten mit Zitat      
Hallo DSP,

nochmals vielen dank für deine Erläuterungen. Das mit der Amplitude klingt ziemlich einleuchtend.

Dann hätte ich aber noch mal eine Frage zum einseitigen Spektrum. Wie gehe ich mit dem Gleichanteil und dem von dir erwähnten Betrag der Nyquistfrequenz um? Und auch auf die Gefahr hin, dass meine (wiederholte) Frage ziemlich rum nervt, aber der Verlauf meines Ergebnisses ist keine sinc-Funktion. Woran liegt das?

Ich habe mir einen Rechteckimpuls gebastelt von 0 - 1 welcher von 0 - 0.5 gleich 1 ist und im weiteren Verlauf 0. Ich erwarte eine Amplitude von 0.5 im Spektrum. Für das zweiseitige Spektrum passt das dann auch aber das einseitige Spektrum muss ich doch noch mit 2 multiplizieren. Dann habe ich aber eine Amplitude von 1 im einseitigen Spektrum (siehe Pdf und Code).

Sorry im Übrigen für mein ständiges Nachhaken, aber der ganze Vorgang mit der FFT-Umsetzung in Matlab ist mir doch noch nicht 100% klar und ich habe festgestellt, dass man zwar sehr leicht irgend was fouriertransformiert bekommt, das aber nicht gleichbedeutend ist mit dem Generieren eines sinnvollen Ergebnisses Confused Wenn du verstehst *g

Vielen Dank auf jeden Fall für deine Hilfe.

Gruß Phtagen

Code:

octave:135> fs=8000;              %Abtastfrequenz
octave:136> fn=fs/2;              %Nyquistfrequenz
octave:137> x=0:1/fs:1-1/fs;      %X-Vektor
octave:139> sp=length(x);         %Signallänge
octave:140> df=fs/sp;             %Frequenzauflösung
octave:141> y=zeros(1,sp);        %Nullvektor
octave:142> y(1:sp/2)=1;          %Rechteckimpuls
octave:143> Y=fft(y);             %DFT
octave:144> P2=abs(Y)/sp;         %Normierung auf N
octave:145> P1=2*P2(1:sp/2);      %einseitiges Spektrum
octave:146> f1=0:df:fs/2-df;      %Frequenzvektor
octave:147> plot(f1(1:100),P1(1:100))
 


SpekRechteck1seitig.pdf
 Beschreibung:

Download
 Dateiname:  SpekRechteck1seitig.pdf
 Dateigröße:  9.98 KB
 Heruntergeladen:  445 mal
 
DSP
Forum-Meister

Forum-Meister


Beiträge: 2.117
Anmeldedatum: 28.02.11
Wohnort: ---
Version: R2014b
     Beitrag Verfasst am: 06.03.2016, 17:01     Titel:
  Antworten mit Zitat      
Das verlinkte Skript zeigt doch genau wie ein einseitiges Spektrum zu skalieren ist.

Code:

N = ... % Anzahl Messwerte

% Berechnung der FFT
% ------------------
H = fft(y, N);
% Berechnung des Amplitudengangs aus dem komplexen Frequenzvektor H:
amplH = abs(H);


% Darstellung des interessierenden Frequenzbereichs des
% Amplitudengangs (0...fn-df) und
% daran angepasste Amplitudenskalierung (Normierung auf N/2):
amplitudengang = [amplH(1)/N amplH(2:N/2)/(N/2)]; % DC-Bin auf N normieren!
 


Die Nyquistfreq. (= fs/2) wird hier wie auch in deinem Beispiel nicht dargestellt. Ansonsten müsste sie auch durch N geteilt werden. Die fft berechnet ein zweiseitiges Spektrum,welches N+1 Messwerte hat. Gleichsignalanteil und Nyquistfreq.anteil kommen nur einmal in dem Ergebnisvektor der fft vor. Die anderen Signalfreq. kommen doppelt vor (konjugiert komplex)...eben für negativen und positiven Bereich. Daher muss unterschiedlich skaliert werden.

Zitat:
Und auch auf die Gefahr hin, dass meine (wiederholte) Frage ziemlich rum nervt, aber der Verlauf meines Ergebnisses ist keine sinc-Funktion. Woran liegt das?


Der Verlauf ist für mich stimmig. Was erwartest du denn?

Gruß DSP
Private Nachricht senden Benutzer-Profile anzeigen
 
Phtagen

Gast

Beiträge: ---
Anmeldedatum: ---
Wohnort: ---
Version: ---
     Beitrag Verfasst am: 07.03.2016, 11:30     Titel:
  Antworten mit Zitat      
Hallo DSP,

die Zeile hab ich doch glatt überlesen Shocked

Ok, dann hab ich das ganze noch mal durch exerziert und den Gleichanteil entsprechend normiert.

Code:

octave:88> fs=8000;              %Abtastfrequenz
octave:89> fn=fs/2;              %Nyquistfrequenz
octave:90> x=0:1/fs:1-1/fs;      %X-Vektor
octave:91> sp=length(x);         %Signallänge
octave:92> df=fs/sp;             %Frequenzauflösung
octave:93> y=zeros(1,sp);        %Nullvektor
octave:94> y(1:sp/2)=1;          %Rechteckimpuls
octave:95> Y=fft(y);             %DFT
octave:96> P2=abs(Y)/sp;         %Normierung auf N
octave:97> P1=2*P2(1:sp/2);      %einseitiges Spektrum
octave:98> P1(1)=0.5*P1(1);      %Gleichanteil
octave:99> f1=0:df:fs/2-df;      %Frequenzvektor
octave:100> plot(f1(1:100),P1(1:100))
 


Entsprechend erhalte ich angehängte erste Grafik Rechteck1seitig). Dass der Gleichanteil betragsmäßig kleiner ist als der Peak direkt daneben verwundert mich halt. Ist das korrekt?

Ich hätte laut div. Literaturangaben einen Verlauf des Amplitudenspektrums erwartet wie er in der zweiten Abb. angedeutet ist (hat mit dem Spektrum des Rechteckimpulses nichts zu tun, dient nur der Veranschaulichung dessen, was ich meine), halt der Verlauf der sinc-Funktion (oder Si-Funktion). Stattdessen sieht das von mir erhaltene Spektrum anders aus. Es ist halt keine glatte Si-Funktion.

Hoffe du verstehst was ich meine. Vielen Dank auf jeden Fall für deine Erklärungen.

Gruß Phtagen

erwarteter Verlauf.pdf
 Beschreibung:

Download
 Dateiname:  erwarteter Verlauf.pdf
 Dateigröße:  49.13 KB
 Heruntergeladen:  418 mal
Rechteck1seitig.pdf
 Beschreibung:

Download
 Dateiname:  Rechteck1seitig.pdf
 Dateigröße:  10.86 KB
 Heruntergeladen:  461 mal
 
DSP
Forum-Meister

Forum-Meister


Beiträge: 2.117
Anmeldedatum: 28.02.11
Wohnort: ---
Version: R2014b
     Beitrag Verfasst am: 07.03.2016, 12:01     Titel:
  Antworten mit Zitat      
Hallo,

der erwartete Verlauf ist schon mal wegen der Betragsbildung unmöglich. Für mich ist es nicht nachvollziehbar wie hier negative y-Werte vorkommen sollen. Schon in dem von dir verlinkten Buch, wird angeblich ein lineares Spektrum mit ähnlichem Verlauf gezeigt. Eine Beschriftung der Achsen ist dort auch nicht vorhanden. Daher ist für mich dieser lineare Verlauf des Spektrums nicht nachvollziehbar.

Zitat:
Dass der Gleichanteil betragsmäßig kleiner ist als der Peak direkt daneben verwundert mich halt. Ist das korrekt?


Ja, das ist korrekt. Der Gleichsignalteil ist nichts weiter als eine Mittelwertbildung aller Messwerte. Bei einer Amplitude von 1 und gleichem Anteil high/low muss da 0.5 rauskommen. Für die eigentliche Signalamplitude müsste ja eigentlich 1 herauskommen. Das hatte ich aber schon bereits erklärt wieso nicht.

Wenn du dir mal die Messpunkte im Plot darstellen lässt, wird auch klar, warum es einen zackigen Verlauf gibt. Es sind einfach zu wenige Punkte um einen runden Verlauf zu erzeugen.

Code:

plot(f1,P1, 'Marker','*','MarkerSize',3)
xlim([0 100]) % einfacher als die array Begrenzung f1(1:100)
 


Gruß DSP
Private Nachricht senden Benutzer-Profile anzeigen
 
Ptagen

Gast

Beiträge: ---
Anmeldedatum: ---
Wohnort: ---
Version: ---
     Beitrag Verfasst am: 07.03.2016, 12:23     Titel:
  Antworten mit Zitat      
Hallo DSP,

mhmmmpfff… Das mit der Betragsbildung ist wohl… ziemlich offensichtlich. Manchmal sehe ich halt den Wald vor lauter Bäumen nicht. Ich hab mir auch noch mal

y(t) = |sinc(2*pi*t)|

geplottet. Die Ähnlichkeit zu meiner aller ersten Grafik ist schon frappierend Very Happy .

Ok, dann habe ich zumindest etwas mehr Sicherheit bei der FFT/DFT bekommen.

Ich danke dir auf jeden Fall für die Diskussion. Hat mir sehr weitergeholfen. Dann werde ich mich mal mit dem nächsten Thema auf meiner Agenda beschäftigen und behaupte jetzt einfach mal, dass ich mich sicherlich demnächst wieder hier melden werde.

Viele Grüße

Phtagen
 
Neues Thema eröffnen Neue Antwort erstellen



Einstellungen und Berechtigungen
Beiträge der letzten Zeit anzeigen:

 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 | goMatlab RSS Button 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.