# [SciPy-dev] reimplementation of lfilter

Sturla Molden sturla@molden...
Tue Sep 22 20:01:53 CDT 2009

```Ralf Gommers skrev:
> Sure, I know. The docs I linked have improvements over svn. If you
> merge your changes, there'll be a conflict and those docs will not get
> merged soon - I was hoping to avoid that.

Oh, the docs, not the code... Would something like this be better?

The Cython module doesn't need to have docs in it.

from linear_filter import lfilter as _lfilter

def lfilter(b, a, x, axis=-1, zi=None, direction=1, inplace=-1):
"""Filter data along one-dimension with an IIR or FIR filter.

Description

Filter a data sequence, x, using a digital filter.  This works for many
fundamental data types (including Object type).  The filter is a direct
form II transposed implementation of the standard difference equation
(see "Algorithm").

Inputs:

b -- The numerator coefficient vector in a 1-D sequence.
a -- The denominator coefficient vector in a 1-D sequence. If a[0]
is not 1, then both a and b are normalized by a[0].
x -- An N-dimensional input array.
axis -- The axis of the input data array along which to apply the
linear filter. The filter is applied to each subarray along
this axis. If axis is -1, the filter is applied to x
flattened. (*Default* = -1)
zi -- Initial conditions for the filter delays.  It is a vector
(or array of vectors for an N-dimensional input) of length
max(len(a),len(b)). If zi=None or is not given then initial
rest is assumed. If axis is None, zi should be a vector
of length K = max(M,N), where M=len(b)-1 and N=len(a)-1. If
an axis is not None, zi should either be a vector of length
K, or an array with same shape as x, except for zi.shape[axis]
which should be K. If zi is a vector, the same initial
conditions are applied to each vector along the specified
axis. Default: zi=None (initial rest).
direction -- Filtered in the reverse direction if -1. (*Default* = 1)
inplace -- If possible, modification to x is done implace.
(*Default* = -1)

Outputs: (y, {zf})

y -- The output of the digital filter.
zf -- If zi is None, this is not returned, otherwise, zf holds the
final filter delay values.

Algorithm:

The filter function is implemented as a direct II transposed structure.
This means that the filter implements

a[0]*y[n] = b[0]*x[n] + b[1]*x[n-1] + ... + b[nb]*x[n-nb]
- a[1]*y[n-1] - ... - a[na]*y[n-na]

using the following difference equations:

y[m] = b[0]*x[m] + z[0,m-1]
z[0,m] = b[1]*x[m] + z[1,m-1] - a[1]*y[m]
...
z[n-3,m] = b[n-2]*x[m] + z[n-2,m-1] - a[n-2]*y[m]
z[n-2,m] = b[n-1]*x[m] - a[n-1]*y[m]

where m is the output sample number and n=max(len(a),len(b)) is the
model order.

The rational transfer function describing this filter in the
z-transform domain is
-1               -nb
b[0] + b[1]z  + ... + b[nb] z
Y(z) = ---------------------------------- X(z)
-1               -na
a[0] + a[1]z  + ... + a[na] z

"""
if isscalar(a):
a = [a]

if axis < 0: axis = None
if inplace < 0: inplace = None

return _lfilter(b, a, x, axis=axis, zi=zi,
inplace=inplace, direction=direction)

Regards,
S.M.

```