Peak-Erkennung: Wie sicher ist ein Treffer?

Ein Korrelationsausgang schlägt einmal aus — ist das ein echter Treffer oder nur Zufall? Standardabweichung, Sigma-Signifikanz, Umrechnung zu Dezibel und die praktische Umsetzung als Peak-Detektor.

Ein Korrelator vergleicht ein empfangenes Signal fortlaufend mit einem bekannten Muster und liefert dabei nicht nur an einer Stelle einen Wert, sondern über die gesamte Zeit hinweg einen schwankenden Ausgangsstrom. Irgendwo darin verbirgt sich — wenn das gesuchte Muster tatsächlich vorhanden ist — ein deutlicher Ausschlag. Die praktische Frage lautet dann: Ab wann ist ein Ausschlag ein echter Treffer und nicht bloß eine zufällig hohe Rauschspitze? Dieser Artikel beantwortet das allgemein, mit einem konkreten Beispiel aus PZF-Demodulation in der Praxis.

Die Standardabweichung: Wie stark schwankt das Rauschen?

Die Standardabweichung (Formelzeichen $\sigma$, „Sigma") beschreibt, wie stark eine Reihe von Messwerten üblicherweise um ihren Mittelwert streut.

  • Kleine Standardabweichung → die Werte liegen eng beieinander.
  • Große Standardabweichung → auch deutliche Ausreißer nach oben und unten sind normal.

Für eine Reihe von Werten $x_1, \dots, x_n$ mit Mittelwert $\bar{x}$ berechnet sie sich als:

$$\sigma = \sqrt{\frac{1}{n}\sum_{i=1}^{n}(x_i - \bar{x})^2}$$

In Worten: für jeden Wert den Abstand zum Mittelwert bilden, quadrieren (damit sich positive und negative Abweichungen nicht aufheben), diese Quadrate mitteln, und aus dem Ergebnis die Wurzel ziehen.

Praktische Ermittlung an einem Korrelationsausgang: Man nimmt den kompletten Ausgangsstrom des Korrelators über einen ausreichend langen Zeitraum (abzüglich einer eventuellen Einschwingzeit am Anfang) und berechnet Mittelwert und Standardabweichung über alle diese Werte. Da echte Treffer normalerweise nur kurze, seltene Ausschläge sind, verändern sie dieses Ergebnis kaum messbar — es entspricht damit in guter Näherung dem „Rauschteppich" zwischen den Treffern, auch wenn rechnerisch der gesamte Strom einfließt.

Vom Sigma-Wert zur Trefferwahrscheinlichkeit

Hat man Mittelwert $\mu$ und Standardabweichung $\sigma$ des Rauschens ermittelt, lässt sich für jeden einzelnen Ausschlag (Peak-Höhe $P$) berechnen, wie viele Standardabweichungen er über dem Rauschmittel liegt:

$$k = \frac{P - \mu}{\sigma}$$

Dieser Wert $k$ (in „Sigma" angegeben, z. B. „6,8σ") ist aussagekräftiger als die absolute Peak-Höhe, weil er direkt mit einer Wahrscheinlichkeit verknüpft ist. Rauschen verhält sich in guter Näherung wie eine Normalverteilung (Glockenkurve) um seinen Mittelwert. Für eine solche Verteilung gilt:

Sigma-Wert $k$ Wahrscheinlichkeit für einen zufälligen Ausreißer dieser Höhe (einseitig)
ca. 1 zu 741
ca. 1 zu 3,5 Millionen
ca. 1 zu 1 Milliarde
ca. 1 zu 390 Milliarden

Als Faustregel gilt: Ab etwa 5σ ist ein zufälliger Treffer praktisch ausgeschlossen. Dieselbe Sigma-Konvention wird z. B. auch in der Teilchenphysik verwendet, um die statistische Sicherheit einer Entdeckung anzugeben.

Sigma und Dezibel: zwei Sprachen für dieselbe Sache

In der Nachrichtentechnik wird Signalqualität meist als Signal-Rausch-Verhältnis (SNR) in Dezibel angegeben — ein Verhältnis von Signalleistung zu Rauschleistung, üblich z. B. bei einem durchgehend vorhandenen Träger. Bei der Peak-Erkennung geht es dagegen um eine andere Frage: nicht „wie gut ist der Empfang insgesamt", sondern „ist dieser eine Ausschlag echt". Beide Größen lassen sich aber ineinander umrechnen, wenn man den Peak als Signalamplitude und die Standardabweichung des Rauschens als Rauschamplitude auffasst:

$$\text{SNR}{\text{dB}} = 20 \cdot \log{10}(k) \qquad\Longleftrightarrow\qquad k = 10^{,\text{SNR}_{\text{dB}}/20}$$

(Der Faktor 20 statt 10 gilt, weil hier Amplituden verglichen werden, nicht Leistungen — Leistung ist proportional zum Quadrat der Amplitude, und $10\log_{10}(k^2) = 20\log_{10}(k)$.)

Sigma-Wert $k$ SNR in dB
9,5 dB
14,0 dB
6,8σ 16,7 dB
10σ 20,0 dB
Voraussetzung dieser Umrechnung: Sie gilt nur, wenn das Rauschen näherungsweise normalverteilt und der Peak eine einzelne, gegen dieses Rauschen abgehobene Amplitude ist — genau das Szenario eines Korrelators, der nach einem kurzen, bekannten Muster in einem sonst rauschenden Signal sucht. Für ein durchgehend moduliertes Signal (z. B. die Sprachqualität einer FM-Übertragung) ist die übliche Leistungs-SNR in dB die richtige, aber andere Größe.

Praktische Umsetzung: Schwellwert mit Hysterese

Um aus dem schwankenden Korrelationsausgang automatisch „hier ist ein Treffer" zu erkennen, reicht ein einzelner Schwellwert meist nicht aus: Liegt die Schwelle knapp unter dem Peak, kann das Signal beim Auf- und Abschwingen mehrfach kurz darüber und darunter pendeln und mehrere Treffer melden, obwohl nur einer stattfand. Üblich ist deshalb eine Hysterese mit zwei Schwellen:

  • eine obere Schwelle, die den Beginn eines Treffers markiert,
  • eine untere Schwelle (niedriger als die obere), deren Unterschreiten den Treffer als abgeschlossen zählt.

Zwischen beiden Schwellen wird währenddessen der höchste erreichte Wert mitgeführt — das ist die tatsächliche Peak-Höhe. Zusätzlich sperrt man die Erkennung für eine kurze Zeit nach jedem gemeldeten Treffer (eine Sperrzeit bzw. Refraktärzeit), damit unmittelbar folgendes Nachschwingen nicht als weiterer, eigenständiger Treffer gezählt wird.

Beide Schwellen werden anhand der gemessenen Verteilung festgelegt: oberhalb der üblichen Rauschausschläge, aber mit Abstand unterhalb der tatsächlichen Peak-Höhe.

Umsetzung als GNU-Radio-Block

In GNU Radio lässt sich das als Embedded Python Block (ein kleines, direkt im Flowgraph eingebettetes Stück Python-Code) umsetzen, der bei jedem erkannten Treffer eine Message (ein einzelnes Datenpaket, unabhängig vom fortlaufenden Signalstrom, auch PDU genannt) verschickt:

import numpy as np
from gnuradio import gr
import pmt

class blk(gr.sync_block):
    """Peak-Erkennung mit Hysterese und Sperrzeit, gibt eine PMT-Message
    pro erkanntem Peak aus (Wert + Sample-Index)."""

    def __init__(self, schwelle_hoch=0.026, schwelle_tief=0.016,
                 sperrzeit_samples=5000):
        gr.sync_block.__init__(
            self, name='Peak-Detektor (Message)',
            in_sig=[np.float32], out_sig=None
        )
        self.schwelle_hoch = schwelle_hoch
        self.schwelle_tief = schwelle_tief
        self.sperrzeit = sperrzeit_samples
        self.message_port_register_out(pmt.intern('peak'))
        self._ueber_schwelle = False
        self._max_wert = 0.0
        self._max_index_global = 0
        self._sample_zaehler = 0
        self._letzter_peak = -sperrzeit_samples

    def work(self, input_items, output_items):
        x = input_items[0]
        for i, wert in enumerate(x):
            idx = self._sample_zaehler + i
            if not self._ueber_schwelle:
                if wert > self.schwelle_hoch and (idx - self._letzter_peak) > self.sperrzeit:
                    self._ueber_schwelle = True
                    self._max_wert = wert
                    self._max_index_global = idx
            else:
                if wert > self._max_wert:
                    self._max_wert = wert
                    self._max_index_global = idx
                if wert < self.schwelle_tief:
                    meldung = pmt.make_dict()
                    meldung = pmt.dict_add(meldung, pmt.intern('value'), pmt.from_double(float(self._max_wert)))
                    meldung = pmt.dict_add(meldung, pmt.intern('sample'), pmt.from_uint64(int(self._max_index_global)))
                    self.message_port_pub(pmt.intern('peak'), pmt.cons(meldung, pmt.PMT_NIL))
                    self._letzter_peak = self._max_index_global
                    self._ueber_schwelle = False

        self._sample_zaehler += len(x)
        return len(x)
import numpy as np
from gnuradio import gr
import pmt

class blk(gr.sync_block):
    """Peak-Erkennung mit Hysterese und Sperrzeit, gibt eine PMT-Message
    pro erkanntem Peak aus (Wert + Sample-Index)."""

    def __init__(self, schwelle_hoch=0.026, schwelle_tief=0.016,
                 sperrzeit_samples=5000):
        gr.sync_block.__init__(
            self, name='Peak-Detektor (Message)',
            in_sig=[np.float32], out_sig=None
        )
        self.schwelle_hoch = schwelle_hoch
        self.schwelle_tief = schwelle_tief
        self.sperrzeit = sperrzeit_samples
        self.message_port_register_out(pmt.intern('peak'))
        self._ueber_schwelle = False
        self._max_wert = 0.0
        self._max_index_global = 0
        self._sample_zaehler = 0
        self._letzter_peak = -sperrzeit_samples

    def work(self, input_items, output_items):
        x = input_items[0]
        for i, wert in enumerate(x):
            idx = self._sample_zaehler + i
            if not self._ueber_schwelle:
                if wert > self.schwelle_hoch and (idx - self._letzter_peak) > self.sperrzeit:
                    self._ueber_schwelle = True
                    self._max_wert = wert
                    self._max_index_global = idx
            else:
                if wert > self._max_wert:
                    self._max_wert = wert
                    self._max_index_global = idx
                if wert < self.schwelle_tief:
                    meldung = pmt.make_dict()
                    meldung = pmt.dict_add(meldung, pmt.intern('value'), pmt.from_double(float(self._max_wert)))
                    meldung = pmt.dict_add(meldung, pmt.intern('sample'), pmt.from_uint64(int(self._max_index_global)))
                    self.message_port_pub(pmt.intern('peak'), pmt.cons(meldung, pmt.PMT_NIL))
                    self._letzter_peak = self._max_index_global
                    self._ueber_schwelle = False

        self._sample_zaehler += len(x)
        return len(x)

Der Block hat einen einzelnen Stream-Eingang (den Korrelationsausgang, float) und keinen Stream-Ausgang, dafür einen Message-Ausgang namens peak. Angeschlossen an einen Message Debug-Block lässt sich jede erkannte Meldung direkt anzeigen; ebenso lässt sich der Ausgang an eigene Folgelogik anschließen, etwa um aus dem erkannten Zeitpunkt eine Uhrzeit oder ein Steuersignal abzuleiten.

Praktisches Beispiel: In PZF-Demodulation in der Praxis wird genau dieser Block eingesetzt, um die einmal pro Sekunde auftretenden Korrelationspeaks der DCF77-Phasenmodulation zu erkennen. Die Schwellen wurden dort anhand der gemessenen Rausch- und Peak-Verteilung festgelegt (0,026 bzw. 0,016, bei Peaks zwischen 0,029 und 0,038 und Rauschspitzen bis etwa 0,02) und lieferten im Test exakt eine Meldung pro Sekunde, mit einer Abweichung von unter einer hundertstel Sekunde.

Zusammenfassung

  1. Die Standardabweichung des Korrelations-Rauschens lässt sich direkt aus dem laufenden Ausgangsstrom berechnen.
  2. Ein Peak wird nicht anhand seines Absolutwerts bewertet, sondern anhand seines Abstands zum Rauschen in Vielfachen der Standardabweichung — daraus folgt unmittelbar, wie unwahrscheinlich ein Zufallstreffer wäre.
  3. Dieser Sigma-Wert lässt sich in ein gewohntes SNR in Dezibel umrechnen ($\text{SNR}{\text{dB}} = 20\log{10}(k)$) und umgekehrt.
  4. Praktisch umgesetzt wird die Erkennung über zwei Schwellwerte (Hysterese) plus eine Sperrzeit, um Mehrfachauslösungen an einem einzelnen Peak zu vermeiden — in GNU Radio z. B. als Embedded Python Block mit Message-Ausgang.