onyx.signalprocessing.filter – Simple FIR and IIR signal processing

An unstable second-order filter

>>> f1 = type2stageDevelop((5, 6), (1, 2), 2)
>>> f1
type2stageDevelop((5, 6), (1, 2), 2)
>>> f1.process_one(100)
[200]
>>> f1.process_one(0)
[-800]
>>> f1.process_one(0)
[3200]

An FIR filter, using default coeff_k=1

>>> f2 = type2stageDevelop((), (2, 3))
>>> f2.process_some((1, 0, 0, 0, -1, -2, -3, 0, 0, 0))
[1, 2, 3, 0, -1, -4, -10, -12, -9, 0]

A first-order filter

>>> for y in type2stageDevelop((-0.75,), (), 30).process_some((1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0)):
...   print '-' * int(y), int(y)
------------------------------ 30
---------------------- 22
---------------- 16
------------ 12
--------- 9
------- 7
----- 5
---- 4
--- 3
-- 2
- 1
- 1

A stable second-order filter

>>> for y in type2stageDevelop((-1.5, 0.875), (), 10).process_some((1,) * 20):
...   print '-' * int(y), int(y)
---------- 10
------------------------- 25
-------------------------------------- 38
---------------------------------------------- 46
--------------------------------------------- 45
------------------------------------- 37
-------------------------- 26
----------------- 17
------------ 12
------------- 13
------------------- 19
--------------------------- 27
--------------------------------- 33
------------------------------------- 37
----------------------------------- 35
------------------------------- 31
------------------------- 25
--------------------- 21
------------------- 19
-------------------- 20

A first-order, bi-linear ‘pre-emphasis’ filter run at a low and a high frequency

>>> map(int, type2stageDevelop((0.5,), (-0.05,), 10).process_some((1, 0, -1, 0) * 4 + (0, 0, 0, 0) + (1, -1) * 8))
[10, -5, -7, 4, 7, -4, -7, 4, 7, -4, -7, 4, 7, -4, -7, 4, -2, 1, 0, 0, 9, -15, 18, -19, 20, -20, 20, -20, 20, -20, 20, -20, 20, -20, 20, -20]

>>> f4 = type2stageDevelop((-1.5, 0.875), (), 10)
>>> def yow(x): print 'yow', int(x)
>>> f4.get_recipient()
>>> f4.set_recipient(yow)
>>> f4.send_many([1] * 20)
yow 10
yow 25
yow 38
yow 46
yow 45
yow 37
yow 26
yow 17
yow 12
yow 13
yow 19
yow 27
yow 33
yow 37
yow 35
yow 31
yow 25
yow 21
yow 19
yow 20
>>> f4.get_recipient()  
<function yow at 0x...>
>>> f4.set_recipient(None)
>>> f4.send_many([0] * 20)
>>> f4.set_recipient(yow)
>>> f4.send_many([0] * 20)
yow 2
yow 0
yow -2
yow -4
yow -3
yow -1
yow 0
yow 2
yow 2
yow 2
yow 1
yow 0
yow -1
yow -2
yow -1
yow 0
yow 0
yow 1
yow 1
yow 1
>>> f4.send_many([1] * 5)
yow 10
yow 24
yow 37
yow 45
yow 44
>>> f4.send(0)
yow 27

Get the state, and then run it for a few

>>> s = f4.get_state()
>>> f4.send_many([1] * 5)
yow 12
yow 4
yow 5
yow 14
yow 27

Reset it and see that it’s dead

>>> f4.reset()
>>> f4.send_many([0] * 3)
yow 0
yow 0
yow 0

Restore the state and get the same output as we got prior to the reset

>>> s = f4.set_state(s)
>>> f4.send_many([1] * 5)
yow 12
yow 4
yow 5
yow 14
yow 27

>>> f4.reset()
>>> f4.send(0)
yow 0

>>> f4.send_many([1] * 5)
yow 10
yow 25
yow 38
yow 46
yow 45
>>> f4.send_signals(recipient=yow, reset=True)
>>> f4.send(0)
yow 0

The trivial case

>>> f3 = type2stageDevelop((), ())
>>> f3
type2stageDevelop((0, 0), (0, 0), 1)
>>> f3.process_some(xrange(20))
[0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19]

A degenerate case

>>> f4 = type2stageDevelop((), (), 0)
>>> f4
type2stageDevelop((0, 0), (0, 0), 0)
>>> f4.process_some(xrange(20))
[0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]

class onyx.signalprocessing.filter.HermanskyRasta(low_pass_coeff=0.984375, sendee=None, sending=True, bypass_types=())

Bases: onyx.dataflow.streamprocess.ProcessorBase

An implementation of the fourth-order RASTA filter found in RASTA-PLP Speech Analysis, by Hermansky, Morgan, et al., TR-910-069, December 1991.

Their transfer function gives this non-causal difference equation:

y[n-4] + -0.98 y[n-5] = 0.1 (2 x[n] + x[n-1] - x[n-3] - 2 x[n-4])

By adding four poles at the origin we can solve for causal y[n] without affecting transfer function amplitude:

y[n] = 0.98 y[n-1] + 0.1 (2 x[n] + x[n-1] - x[n-3] - 2 x[n-4])

The low_pass_coeff defaults to a low-precision value close to 0.98, and internally a low-precision approximation to 0.1 is used.

>>> hr = HermanskyRasta()
>>> hr
HermanskyRasta(0.984375,)

Look at impulse response; starts with the mixed FIR/IIR portion, then the one-pole IIR decay

>>> results = list()
>>> hr.set_sendee(results.append)
>>> hr.process(1)
>>> for i in xrange(20):
...   hr.process(0)

>>> for result in results: print int(round(10000 * result))
2031
3015
2968
1906
-155
-153
-150
-148
-146
-143
-141
-139
-137
-135
-132
-130
-128
-126
-124
-122
-120

Set IIR coefficient to 0; look at impulse response of FIR portion only

>>> results = list()
>>> hr2 = HermanskyRasta(0, sendee=results.append)
>>> hr2.process(1)
>>> for i in xrange(9):
...   hr2.process(0)

>>> for result in results: print int(round(10000 * result))
2031
1016
0
-1016
-2031
0
0
0
0
0

Works on numpy arrays

>>> results = list()
>>> hr3 = HermanskyRasta(sendee=results.append)
>>> hr3.process(np.array((100, -1000)))
>>> for _ in xrange(10):
...   hr3.process(np.array((0,0)))

>>> for result in results: print repr(result)
array([  20.3125, -203.125 ], dtype=float32)
array([  30.15136719, -301.51367188], dtype=float32)
array([  29.68025208, -296.80252075], dtype=float32)
array([  19.06024742, -190.60249329], dtype=float32)
array([ -1.5500679 ,  15.50067902], dtype=float32)
array([ -1.52584887,  15.25847626], dtype=float32)
array([ -1.50200641,  15.02006721], dtype=float32)
array([ -1.47853744,  14.7853775 ], dtype=float32)
array([ -1.45543575,  14.55435753], dtype=float32)
array([ -1.43269515,  14.32693291], dtype=float32)
array([ -1.41030943,  14.10307884], dtype=float32)

dc

A debug context for this processor. This attribute is an object returned by dcheck() in the onyx.util.debugprint module, and may be used for debug output. The tag for turning on such output is available as debug_tag

debug_tag

Activate this tag to get debugging information, see onyx.util.debugprint.DebugPrint

graph

Return a graph for this processor. By default this is just a single node whose label is the label of processor; derived classes may wish to override this property.

label

Return a label for this processor. By default this is just the name of the class; derived classes may wish to override this property by providing a different label to __init__().

process(event)

Process the numeric value, x[n], through the filter, sending the resulting y[n].

If the first value sent after initialization or reset() is a numpy array, then all value must have the same shape. In this case, each element of value is filtered independently.

reset()

Returns the filter to its initialized state.

send(result)

Internal function that pushes result into the sendee. Implementations of process() must call this to push results. To set up the sendee, (the target of the push), clients of the processor must either initialize the object with a sendee, or call set_sendee(). Processors created with a sendee of False will never send, but will not error if send is called.

sendee

The callable this processor will use to send events; see set_sendee()

sending

Whether this processor will currently send events at all; see set_sending()

set_sendee(sendee)

Clients call this to set up the callable where the processor will send its results.

set_sending(sending)

Clients call this to turn sending from a processor on or off.

static std_process_prologue(process_function)

Subclasses may use this decorater on their process function to implement the usual bypass and process semantics and to set up the debug context returned by dc().

class onyx.signalprocessing.filter.Type2Stage(a1, a2, b1, b2, k=1, sendee=None, sending=True, bypass_types=())

Bases: onyx.dataflow.streamprocess.ProcessorBase

A canonic Type II filter that implements the inhomogeneous second-order difference equation:

y[n] + a1 y[n-1] + a2 y[n-2] = k (x[n] + b1 x[n-1] + b2 x[n-2])

generating the solution y[n] given input sequence x[n]:

y[n] = -(a1 y[n-1] + a2 y[n-2]) + k (x[n] + b1 x[n-1] + b2 x[n-2])

The arguments a1, a2, b1, and b2 are used as above. The optional k argument is the overall scaling factor, (default 1).

An unstable second-order filter

>>> f1 = Type2Stage(5, 6, 1, 2, 2)
>>> f1
Type2Stage(5, 6, 1, 2, 2)

>>> results = list()
>>> f1.set_sendee(results.append)
>>> f1.process(100)
>>> f1.process(0)
>>> f1.process(0)
>>> for result in results: print result
200
-800
3200

Works on numpy vectors, in this case with a stable second-order filter

>>> results = list()
>>> f2 = Type2Stage(-1.5, 0.875, 0, 0, 10, sendee=results.append)
>>> f2.process(np.array((100, -1000)))
>>> for _ in xrange(7):
...   f2.process(np.array((0,0)))

>>> for result in results: print repr(result)
array([  1000., -10000.], dtype=float32)
array([  1500., -15000.], dtype=float32)
array([  1375., -13750.], dtype=float32)
array([  750., -7500.], dtype=float32)
array([ -78.125,  781.25 ], dtype=float32)
array([ -773.4375,  7734.375 ], dtype=float32)
array([ -1091.796875,  10917.96875 ], dtype=float32)
array([ -960.9375,  9609.375 ], dtype=float32)

dc

A debug context for this processor. This attribute is an object returned by dcheck() in the onyx.util.debugprint module, and may be used for debug output. The tag for turning on such output is available as debug_tag

debug_tag

Activate this tag to get debugging information, see onyx.util.debugprint.DebugPrint

graph

Return a graph for this processor. By default this is just a single node whose label is the label of processor; derived classes may wish to override this property.

label

Return a label for this processor. By default this is just the name of the class; derived classes may wish to override this property by providing a different label to __init__().

process(event)

Process the numeric value, x[n], through the filter, sending the resulting y[n].

If the first value sent after initialization or reset() is a numpy array, then all value must have the same shape. In this case, each element of value is filtered independently.

reset()

Returns the filter to its initialized state.

send(result)

Internal function that pushes result into the sendee. Implementations of process() must call this to push results. To set up the sendee, (the target of the push), clients of the processor must either initialize the object with a sendee, or call set_sendee(). Processors created with a sendee of False will never send, but will not error if send is called.

sendee

The callable this processor will use to send events; see set_sendee()

sending

Whether this processor will currently send events at all; see set_sending()

set_sendee(sendee)

Clients call this to set up the callable where the processor will send its results.

set_sending(sending)

Clients call this to turn sending from a processor on or off.

static std_process_prologue(process_function)

Subclasses may use this decorater on their process function to implement the usual bypass and process semantics and to set up the debug context returned by dc().

class onyx.signalprocessing.filter.nonlinear2pole(r, theta)

Bases: onyx.dataflow.processor

A time-variable, nonlinear, two-pole, discrete-time system, with zeros at +1 and -1. Argument r is the radius to each pole, and theta is the angle in radians. The overall gain of the system is set so that H|j theta| is approximately 1 for the value of r.

In the running system, you can vary r, via method set_r(), for time-varying or nonlinear behavior. Note that set_r() does not change the overall gain.

Demonstrate that abs(H(theta)) ~= 1 by plotting 100 y[n] for a system driven with a sine wave at theta

>>> r0, theta0 = 0.85, 0.5
>>> for y in nonlinear2pole(r0, theta0).process_some(math.sin(i*theta0) for i in xrange(36)):
...   p = int(50 * (y + 1))
...   print '-' * p, int(round(100 * y))
-------------------------------------------------- 0
----------------------------------------------------- 7
------------------------------------------------------------ 22
------------------------------------------------------------------- 35
-------------------------------------------------------------------- 37
-------------------------------------------------------------- 25
------------------------------------------------- -1
--------------------------------- -32
------------------- -60
------------ -75
-------------- -72
------------------------- -49
-------------------------------------------- -11
----------------------------------------------------------------- 31
----------------------------------------------------------------------------------- 68
---------------------------------------------------------------------------------------------- 89
---------------------------------------------------------------------------------------------- 88
---------------------------------------------------------------------------------- 65
-------------------------------------------------------------- 26
--------------------------------------- -21
------------------ -63
---- -91
-- -96
----------- -77
------------------------------ -40
----------------------------------------------------- 8
----------------------------------------------------------------------------- 54
--------------------------------------------------------------------------------------------- 87
--------------------------------------------------------------------------------------------------- 99
--------------------------------------------------------------------------------------------- 86
---------------------------------------------------------------------------- 53
----------------------------------------------------- 6
---------------------------- -42
--------- -80
 -99
--- -93

Linear ringing response to a step

>>> r1, theta1 = 0.95, 0.5
>>> h1 = nonlinear2pole(r1, theta1)
>>> for y in h1.process_some((7.5,) * 50):
...   p = int(35 * (y + 1))
...   print '-' * p, int(round(100 * y))
----------------------------------------------- 37
--------------------------------------------------------------------- 98
-------------------------------------------------------------------------------- 130
------------------------------------------------------------------------------- 128
-------------------------------------------------------------------- 97
-------------------------------------------------- 46
------------------------------- -11
------------- -60
--- -90
- -96
------- -78
------------------- -44
--------------------------------- -3
----------------------------------------------- 35
-------------------------------------------------------- 61
----------------------------------------------------------- 70
-------------------------------------------------------- 62
------------------------------------------------ 40
-------------------------------------- 11
---------------------------- -18
--------------------- -40
----------------- -50
------------------ -48
----------------------- -34
------------------------------ -14
------------------------------------- 7
------------------------------------------- 25
----------------------------------------------- 35
----------------------------------------------- 36
-------------------------------------------- 28
---------------------------------------- 15
---------------------------------- -1
----------------------------- -15
-------------------------- -24
------------------------- -27
--------------------------- -23
------------------------------ -14
---------------------------------- -2
------------------------------------- 8
---------------------------------------- 16
----------------------------------------- 19
----------------------------------------- 18
--------------------------------------- 12
------------------------------------ 4
--------------------------------- -4
------------------------------- -10
------------------------------ -14
------------------------------ -13
------------------------------- -10
--------------------------------- -5

Non-linear response to a step, where we move the poles outside the unit circle when the signal amplitude is small

>>> r2, theta2 = 0.98, 0.3
>>> h2 = nonlinear2pole(r2, theta2)
>>> for x in ((9.5,) * 120):
...   y, = h2.process_one(x)
...   h2.set_r(r2 if abs(y) >= 0.125 else 1.02)
...   p = int(35 * (y + 1))
...   print '-' * p, int(100 * y)
----------------------------------------- 18
----------------------------------------------------- 54
---------------------------------------------------------------- 83
----------------------------------------------------------------------- 103
--------------------------------------------------------------------------- 114
--------------------------------------------------------------------------- 114
----------------------------------------------------------------------- 104
----------------------------------------------------------------- 86
-------------------------------------------------------- 60
--------------------------------------------- 30
---------------------------------- 0
----------------------- -32
------------- -60
------ -80
-- -91
- -94
--- -89
-------- -76
-------------- -57
----------------------- -33
-------------------------------- -8
---------------------------------------- 16
------------------------------------------------ 39
------------------------------------------------------- 58
----------------------------------------------------------- 70
------------------------------------------------------------- 76
------------------------------------------------------------- 75
---------------------------------------------------------- 67
------------------------------------------------------ 54
----------------------------------------------- 36
---------------------------------------- 16
--------------------------------- -4
-------------------------- -25
------------------- -43
--------------- -55
------------- -62
------------ -63
-------------- -58
----------------- -48
---------------------- -35
---------------------------- -18
---------------------------------- -1
---------------------------------------- 14
--------------------------------------------- 29
------------------------------------------------- 40
--------------------------------------------------- 48
---------------------------------------------------- 51
---------------------------------------------------- 49
-------------------------------------------------- 43
---------------------------------------------- 34
------------------------------------------ 21
------------------------------------- 8
-------------------------------- -6
--------------------------- -21
----------------------- -33
-------------------- -40
------------------- -43
-------------------- -42
--------------------- -37
------------------------ -29
---------------------------- -19
-------------------------------- -8
----------------------------------- 2
--------------------------------------- 12
------------------------------------------ 22
--------------------------------------------- 29
---------------------------------------------- 33
----------------------------------------------- 34
---------------------------------------------- 32
-------------------------------------------- 28
------------------------------------------ 21
--------------------------------------- 12
----------------------------------- 2
-------------------------------- -8
---------------------------- -18
------------------------- -26
------------------------ -29
------------------------ -30
------------------------- -28
-------------------------- -23
---------------------------- -17
------------------------------- -9
---------------------------------- -2
------------------------------------ 4
-------------------------------------- 11
---------------------------------------- 16
------------------------------------------ 21
------------------------------------------- 23
------------------------------------------- 23
------------------------------------------ 22
----------------------------------------- 18
--------------------------------------- 13
------------------------------------- 6
---------------------------------- 0
-------------------------------- -8
----------------------------- -15
--------------------------- -20
--------------------------- -22
--------------------------- -21
---------------------------- -19
----------------------------- -15
------------------------------- -9
--------------------------------- -4
---------------------------------- 0
------------------------------------ 4
-------------------------------------- 8
--------------------------------------- 12
---------------------------------------- 15
---------------------------------------- 17
----------------------------------------- 17
---------------------------------------- 15
--------------------------------------- 13
-------------------------------------- 9
------------------------------------ 4
---------------------------------- -1
-------------------------------- -7
------------------------------ -13
----------------------------- -16
----------------------------- -16
----------------------------- -15

static approxHmax(r, theta)

Return the approximate largest magnitude of the transfer function for a two-pole discrete-time system with zeros at +1 and -1, and poles at r exp(+-j theta), by evaluating |H(z)| at z = exp(j theta).

configure(**kwargs)

process_one(input)

Process a single input sample. Returns a sequence with the result.

process_some(items)

Process each of the elements of ‘items’. Return a list of results.

send_many(inputs)

set_r(r)

Update the radius of the poles to be r while leaving other parameters of the system unchanged.

class onyx.signalprocessing.filter.type2stageDevelop(coeffs_a, coeffs_b, coeff_k=1)

Bases: onyx.dataflow.processor

A canonic Type II filter that implements the (inhomogeneous) second-order difference equation:

y[n] + a1 y[n-1] + a2 y[n-2] = k (x[n] + b1 x[n-1] + b2 x[n-2])

generating the solution y[n] given input sequence x[n]:

y[n] = -(a1 y[n-1] + a2 y[n-2]) + k (x[n] + b1 x[n-1] + b2 x[n-2])

The coeffs_a is a sequence of zero to two IIR coefficients, the a in the above. The coeffs_b is a sequence of zero to two FIR coefficients, the b in the above. In these two sequences, missing trailing coefficients are treated as having value 0. The optional coeff_k argument is the overall scaling factor, the k in the above, (default 1).

configure(**kwargs)

get_recipient()

get_state()

Returns an opaque object representing the internal state of the filter. This object can be used subsequently as an argument to set_state() to return the filter to its state when get_state() was called.

process_one(input)

process_some(items)

Process each of the elements of ‘items’. Return a list of results.

reset()

Resets the filter to its zero state.

send_many(inputs)

send_signals(**kwargs)

set_recipient(recipient)

set_state(state)

Set the state of the filter according to state, a value returned by the get_state() method.