2014-04-25 11 views
8

In numpy würde ich gerne die Punkte erkennen, bei denen das Signal von einer bestimmten Schwelle (früher) unterschritten wird und über einer bestimmten anderen Schwelle liegt. Dies ist für Dinge wie Entprellung oder genaue Nulldurchgänge in Gegenwart von Lärm usw.Wie finden Sie Nulldurchgänge mit Hysterese?

So:

import numpy 

# set up little test problem 
N = 1000 
values = numpy.sin(numpy.linspace(0, 20, N)) 
values += 0.4 * numpy.random.random(N) - 0.2 
v_high = 0.3 
v_low = -0.3 

# find transitions from below v_low to above v_high  
transitions = numpy.zeros_like(values, dtype=numpy.bool) 

state = "high" 

for i in range(N): 
    if values[i] > v_high: 
     # previous state was low, this is a low-to-high transition 
     if state == "low": 
      transitions[i] = True 
     state = "high" 
    if values[i] < v_low: 
     state = "low" 

ich einen Weg mag dies ausdrücklich über das Array ohne Schleifen zu tun: aber ich Ich kann mir keinen Weg vorstellen, da jeder Zustandswert vom vorherigen Zustand abhängt. Ist es möglich, auf eine Schleife zu verzichten?

def hyst(x, th_lo, th_hi, initial = False): 
    hi = x >= th_hi 
    lo_or_hi = (x <= th_lo) | hi 
    ind = np.nonzero(lo_or_hi)[0] 
    if not ind.size: # prevent index error if ind is empty 
     return np.zeros_like(x, dtype=bool) | initial 
    cnt = np.cumsum(lo_or_hi) # from 0 to len(x) 
    return np.where(cnt, hi[ind[cnt-1]], initial) 

Erläuterung::

Antwort

10

Dies kann etwa so erfolgen ind die Indizes aller Proben sind, wo das Signal unterhalb der unteren oder oberhalb der oberen Schwelle ist, und für die die Position der ‚Schalter "ist also klar definiert. Mit cumsum machen Sie eine Art Zähler, der auf den Index des letzten wohldefinierten Beispiels zeigt. Wenn der Anfang des Eingabevektors zwischen den beiden Schwellenwerten liegt, ist cnt 0, daher müssen Sie den entsprechenden Ausgang mit der Funktion where auf den Anfangswert setzen.

Guthaben: das ist ein Trick, den ich in einem old post auf irgendeinem Matlab Forum fand, das ich zu Numpy übersetzte. Dieser Code ist ein wenig schwer zu verstehen und muss auch verschiedene Zwischen-Arrays zuweisen. Es wäre besser, wenn Numpy eine dedizierte Funktion enthalten würde, die Ihrer einfachen for-Schleife ähnlich ist, aber in C für Geschwindigkeit implementiert wird.

Schnelltest:

x = np.linspace(0,20, 1000) 
y = np.sin(x) 
h1 = hyst(y, -0.5, 0.5) 
h2 = hyst(y, -0.5, 0.5, True) 
plt.plot(x, y, x, -0.5 + h1, x, -0.5 + h2) 
plt.legend(('input', 'output, start=0', 'output, start=1')) 
plt.title('Thresholding with hysteresis') 
plt.show() 

Ergebnis: enter image description here

+0

können Sie, dass 'hyst' Funktion C konvertieren? Vielen Dank. – SpaceDog

+0

@SpaceDog Sie könnten, aber wenn Sie C verwenden, ist es wahrscheinlich besser, eine einfache Schleife ähnlich wie in der ursprünglichen Frage zu schreiben. Der Trick in meiner Antwort ist in Python schneller, da er vektorisierten numpy-Code statt einer langsamen Python-Schleife verwendet. Der vektorisierte Code muss mehrmals über die Daten gehen, während eine einfache Schleife in C alles in einem Durchgang erledigen könnte. –

+0

Ich fragte das, weil ich versuche zu verstehen, was Ihre Funktion tut, also kann ich einen Code in C schreiben ... – SpaceDog

1

Modifikationen Ich hatte für meine Arbeit zu tun, die alle auf der Grundlage der Antwort oben durch Bas Swinckels, Detektion von Schwellendurchgang zu ermöglichen, wenn Standard verwendet sowie umgekehrt Schwellenwerte.

Ich bin nicht glücklich mit der harten Namensgebung, vielleicht ist es jetzt th_hi2lo und th_lo2hi statt th_lo und th_hi lesen sollte? Unter Verwendung der ursprünglichen Werte ist das Verhalten gleich schwer.

def hyst(x, th_lo, th_hi, initial = False): 
    """ 
    x : Numpy Array 
     Series to apply hysteresis to. 
    th_lo : float or int 
     Below this threshold the value of hyst will be False (0). 
    th_hi : float or int 
     Above this threshold the value of hyst will be True (1). 
    """   

    if th_lo > th_hi: # If thresholds are reversed, x must be reversed as well 
     x = x[::-1] 
     th_lo, th_hi = th_hi, th_lo 
     rev = True 
    else: 
     rev = False 

    hi = x >= th_hi 
    lo_or_hi = (x <= th_lo) | hi 

    ind = np.nonzero(lo_or_hi)[0] # Index für alle darunter oder darüber 
    if not ind.size: # prevent index error if ind is empty 
     x_hyst = np.zeros_like(x, dtype=bool) | initial 
    else: 
     cnt = np.cumsum(lo_or_hi) # from 0 to len(x) 
     x_hyst = np.where(cnt, hi[ind[cnt-1]], initial) 

    if rev: 
     x_hyst = x_hyst[::-1] 

    return x_hyst 

Und wie über einen Test des Codes, um zu sehen, was es tut:

x = np.linspace(0,20, 1000) 
y = np.sin(x) 
h1 = hyst(y, -0.2, 0.2) 
h2 = hyst(y, +0.5, -0.5) 
plt.plot(x, y, x, -0.2 + h1*0.4, x, -0.5 + h2) 
plt.legend(('input', 'output, classic, hyst(y, -0.2, +0.2)', 
      'output, reversed, hyst(y, +0.5, -0.5)')) 
plt.title('Thresholding with hysteresis') 
plt.show() 

Sine with two different settings for hysteresis.