www.mikrocontroller.net

Forum: Digitale Signalverarbeitung / DSP FIR Filter um Sinus zu detektieren


Important announcement: there is an English version of this forum on EmbDev.net. Posts you create there will be displayed on Mikrocontroller.net and EmbDev.net.
Autor: Techniker (Gast)
Datum:

Diesen Beitrag bewerten:
lesenswert
nicht lesenswert
Hallo,

ich möchte mit einem FIR Filter einen Sinus detektieren. Der Sinus liegt 
in einem definierten Zeitbereich. Könnte man hierfür einen FIR Bandpass 
einsetzen ?

Für weitere Ratschläge bin ich sehr dankbar!

Autor: Goerzel (Gast)
Datum:

Diesen Beitrag bewerten:
lesenswert
nicht lesenswert
Für die einfache detektion von Frequenzen (z.B. DTMF) wird meist der 
Goerzelalgorithmus verwendet.
Siehe z.B. hier:
http://netwerkt.wordpress.com/2011/08/25/goertzel-filter/

Autor: Techniker (Gast)
Datum:

Diesen Beitrag bewerten:
lesenswert
nicht lesenswert
Guten Morgen,

besten Dank für den Tipp. Gibt es in Matlab für den Goerzelalgorithmus 
eine Funktion ?

Autor: Techniker (Gast)
Datum:

Diesen Beitrag bewerten:
lesenswert
nicht lesenswert
Der Goertzelalgorithmus ist ja im prinzip ein IIR Filter. Ich möchte 
aber einen FIR-Filter hierfür verwenden. Würde es eigentlich Sinn machen 
einen Bandpass FIR-Filter zu verwenden ?

Autor: Techniker (Gast)
Datum:

Diesen Beitrag bewerten:
lesenswert
nicht lesenswert
In Matlab gibt es eine Funktion mit dem Namen goertzel. Leider besitze 
ich diese Funktion nicht. Nun habe ich mir gedacht ich programmiere den 
Goertzlalgorithmus mit der filter Funktion nach.

Siehe Übertragungsfunktion:
http://www.mathworks.de/help/toolbox/signal/ref/go...
%N = 256
for i=0:1:N-1
  f = filter([1 -exp(-i*2*pi*i/N)],[1 -2*cos(2*pi*i/N) 1],1:N);
end;

Müsste diese Implementierung so stimmen ?
Wie füge ich nun das Sinussignal (100kHz mit N = 256 Abtastwerte) hinzu 
?

Autor: H.B. (Gast)
Datum:

Diesen Beitrag bewerten:
lesenswert
nicht lesenswert
Sehe ich nicht, wie das gehen soll, was macht die Filterfunktion?

Musst Du nicht einfach nur für jeden Wert die jeweils beiden 
zurückliegenden Werte aufgreifen und mit dem COS-Faktor multiplizieren, 
wie es in WIki steht? Ich meine das englische Wiki!

Autor: Strubi (Gast)
Datum:

Diesen Beitrag bewerten:
lesenswert
nicht lesenswert
Tip: Etwas generischer geht's mit der Signal Processing toolbox so:
fgoer = tf([0.5], [ 1, -(r*kcos1), r*r ]);
filtered = filter(tf.num, tf.den, data);

tf plus entsprechende Argumente baut die Transfer-Funktion für den 
Goertzel.
r ist typischerweise fast 1 (0.9999...), kcos1 der cos(2*pi*k/N)-Term 
für die gewünschte Frequenz.

Autor: Strubi (Gast)
Datum:

Diesen Beitrag bewerten:
lesenswert
nicht lesenswert
Da hatte ich zu schnell geklickt: In der zweiten Zeile sollte tf = 
fgoer.

Sorry.

Autor: Techniker (Gast)
Datum:

Diesen Beitrag bewerten:
lesenswert
nicht lesenswert
Hallo Strubi,

danke für deine Hilfe. Ich nicht nachvollziehen wie du auf diese Zeile 
kommst:
>>fgoer = tf([0.5], [ 1, -(r*kcos1), r*r ]);

Mir fehlt eine genaue Herleitung des Goerzel Algorithmus.

Autor: Techniker (Gast)
Datum:

Diesen Beitrag bewerten:
lesenswert
nicht lesenswert
Ich habe den Filter so implementiert wie auf der Internetseite 
http://www.google.de/url?sa=t&rct=j&q=sa2011-lab3-... 
beschrieben (Seite 9).
f=filter([1],[1 -exp(-j*2*pi*k/N)],y)

y sind meine Eingangsdaten
N ist Konstant
und k ist der Laufindex von k=0... N-1

Nach der Berechnung mache ich dann eine Betragsbildung mit dem Befehl 
abs.
Stimmt meine Implemntierung in Matlab ?

Autor: Goerzel (Gast)
Datum:

Diesen Beitrag bewerten:
lesenswert
nicht lesenswert
Techniker schrieb:
> Stimmt meine Implemntierung in Matlab ?

Wie wäre es wenn Du Dir einfach ein Testsignal erzeugst, z.B. so was:
f=50;            %Signalfrequenz
fs=100*f;        %Abtastfrequenz
t=0:1/fs:1/f;    %Zeitvektor erstellen

y= sin(2*pi*f*t)+sin(2*pi*f*t)+sin(4*pi*f*t);

Dann lässt Du den Filter drüber laufen.
Anschließend ne FFT und Du weist es.

Grüße

Autor: Techniker (Gast)
Datum:

Diesen Beitrag bewerten:
lesenswert
nicht lesenswert
Irgendwie funktioniert das ganze nicht so richtig.
Was mache ich da noch falsch ?
f=50;            %Signalfrequenz
fs=100*f;        %Abtastfrequenz
t=0:1/fs:1/f;    %Zeitvektor erstellen

y= sin(2*pi*f*t);
max_y = max(abs(y))*1.1; 
subplot(2,1,1);
fig = figure(1); 
plot(t,y);
axis([0 t(end) -max_y max_y]);
grid on;

N = 128;      % Abtastwerte
k = 0:1:N-1;
g=filter([1],[1 -exp(-j*2*pi*k/N)],y);

H = fft(g,N); 
amplH = abs(H);
amplitudengang = fftshift(amplH/N);
subplot(2,1,2);
stem(amplitudengang);
grid on;

Autor: Mathegenie (Gast)
Datum:

Diesen Beitrag bewerten:
lesenswert
nicht lesenswert
Was funktioniert nicht richtig? Siehst du Deine Frequenz nicht, oder 
siehst Du gar nichts?

Autor: Jürgen Schuhmacher (engineer) Benutzerseite
Datum:
Angehängte Dateien:

Diesen Beitrag bewerten:
lesenswert
nicht lesenswert
Mathegenie schrieb:
> Was funktioniert nicht richtig? Siehst du Deine Frequenz nicht, oder
> siehst Du gar nichts?

Kann es sein, dass Du Deine Frequenz deshalb nicht findest, weil Du zu 
grob abtastest?

Die Schrittweite der Testfrequenzen mit denen Du rechnest, muss in einem 
sinnvollen Verhältnis zur Samplezahl stehen, weil mit zunehmender 
Samplezahl die DFT-Rechnung immer selektiver wird und Dir so eine 
Frequenz durchrutscht. Ich habe im Anhang ein altes Beispiel hinkopiert:

Die Schrittweite von 8 führt dazu, dass bei der blauen Kurve mit den 
Punkten die Frequenz einigermassen getroffen wird und sich ein peak 
bildet, bei 10 Schritten liegt man haarscharf daneben und der peak geht 
im Rauschen unter. Bei der gröber auflösenden, türkisen Kurve wurde die 
Berechung schon nach 2048 von 4096 samples abgebrochen.

Umgekehrt führt die Überlegung dahin, dass Du Deine Sanplezahl nicht zu 
gross wählen darfst, wenn Du Deine Frequenz nicht genau triffst, oder Du 
bekommst eben das Ergebnis, daß das Matching gering ist.


Techniker schrieb:
> Würde es eigentlich Sinn machen
> einen Bandpass FIR-Filter zu verwenden ?
Der Bandpass ist nicht anderes als eine Summe von Frequenzen als Maske, 
die gleichzeitig auf den Datenstrom losgelassen werden. Mathematisch ist 
der Bandpass ein Hochpass, von dem ein TP abgzogen wurde. In den 
Koeffizienten zeigt sich das durch ein Integral der oberen Grenzfrequenz 
hin zu 0, von der der Koeefizientensatz unterhalb der GF weggelassen 
wurde.

Wenn Du dir einander annäherst, hast Du am Ende eine Frequenz wie oben.

Wenn Du Deine Frequenz nicht genau kennst, ist der BP gfs besser.

Autor: Techniker (Gast)
Datum:

Diesen Beitrag bewerten:
lesenswert
nicht lesenswert
Hallo Jürgen,

danke für deinen Beitrag. Du meinst ich taste zu grob ab. Auf einem 
System kann habe ich nun mal pro Sinusperiode 128 Abtastwerte. Das 
müsste doch reichen. Ich habe die Vermutung das irgendwo noch ein Fehler 
drin ist. Bin leider noch nicht weitergekommen.

Autor: Techniker (Gast)
Datum:

Diesen Beitrag bewerten:
lesenswert
nicht lesenswert
Mit dem Goertzelalgorithmus möchte ich erkennen, ob der Sinus zum 
Beispiel mit 100Khz vorhanden ist oder nicht. Wenn Impulsstörungen oder 
auch Rauschen additiv mit dem Sinus überlagert sind, soll trotzdem der 
Sinus von 100Khz detektiert werden.

Autor: tsag (Gast)
Datum:

Diesen Beitrag bewerten:
lesenswert
nicht lesenswert
Techniker schrieb:
> vorhanden ist oder nicht
ist eine Frage der Definition. Du wirst beim Rauschen immer einen mehr 
oder weniger grossen Wert für die Korrelation erhalten. Du brauchst dann 
einen Grenzwert zum Entscheiden.

Autor: Jürgen Schuhmacher (engineer) Benutzerseite
Datum:
Angehängte Dateien:

Diesen Beitrag bewerten:
lesenswert
nicht lesenswert
Techniker schrieb:
> danke für deinen Beitrag. Du meinst ich taste zu grob ab. Auf einem
> System kann habe ich nun mal pro Sinusperiode 128 Abtastwerte.

Ich denke, wenn dann diesbezüglich "zu fein". Je mehr Punkte je Sinus, 
desto selektiver wird diese Maske, mit der Du den Datenstrom ansiehst. 
Wenn diese Punktezahl steigt, muss auch die "Maskenanzahl" steigen.

In der Grafik habe ich das veranschaulicht:

* Blau: Sehr genau abgetastet, man erkennt den Verlauf des Signals
* Rot: Mit einem 7-tel der Aufösung abtastet: Das Maximum rutscht durch
* Grün: Mit einem 7-tel der Auflösung bei 7-facher Breite
* Hellgrün genau abgetastet mit 7facher Breite

Der grüne ist also grob/breit genug, um das Signal noch zu erfassen und 
fein genug, um den Verlauf noch erahnen zu können.

Blau wäre etwa die obere Grenze. Grün etwa die untere Grenze.
.

Antwort schreiben

Die Angabe einer E-Mail-Adresse ist freiwillig. Wenn Sie automatisch per E-Mail über Antworten auf Ihren Beitrag informiert werden möchten, melden Sie sich bitte an.

Wichtige Regeln - erst lesen, dann posten!

  • Groß- und Kleinschreibung verwenden
  • Längeren Sourcecode nicht im Text einfügen, sondern als Dateianhang

Formatierung (mehr Informationen...)

  • [c]C-Code[/c]
  • [avrasm]AVR-Assembler-Code[/avrasm]
  • [code]Code in anderen Sprachen, ASCII-Zeichnungen[/code]
  • [math]Formel in LaTeX-Syntax[/math]
  • [[Titel]] - Link zu Artikel




Bitte das JPG-Format nur für Fotos und Scans verwenden!
Zeichnungen und Screenshots im PNG- oder GIF-Format hochladen.
Siehe Bildformate
Hinweis: der ursprüngliche Beitrag ist mehr als 6 Monate alt.
Bitte hier nur auf die ursprüngliche Frage antworten,
für neue Fragen einen neuen Beitrag erstellen.

Mit dem Abschicken erkennst du die Nutzungsbedingungen an.

webmaster@mikrocontroller.netImpressumNutzungsbedingungenWerbung auf Mikrocontroller.net