Package ivs :: Package sigproc :: Module filtering

Module filtering

Filter non-equidistant signals (spectra, timeseries...).

Currently implemented:

Gaussian
Sinc
Box
Iterative Nonlinear filter

Example usage:

Generate some data:

>>> x = np.linspace(0,150,10000)
>>> x_= np.linspace(0,150,100)
>>> y = np.sin(2*pi/20.*x)+np.random.normal(size=len(x))

Apply some filter:

>>> x1,y1,pnts = filter_signal(x,y,"gauss",sigma=1)
>>> x2,y2,pnts = filter_signal(x,y,"pijpers",delta=1)
>>> x3,y3,pnts = filter_signal(x,y,"box",window_width=10)
>>> x3_,y3_,pnts = filter_signal(x,y,"box",window_width=10,x_template=x_,threads=2)

>>> import pylab as pl
>>> p = pl.plot(x,y,'ko')
>>> p = pl.plot(x1,y1,'ro',mec='r',label='Gauss')
>>> p = pl.plot(x2,y2,'bo',mec='b',label='Pijpers')
>>> p = pl.plot(x3,y3,'go',mec='g',label='Box')
>>> p = pl.plot(x3_,y3_,'gs',label='Box with template')
>>> p = pl.legend()

Functions

[hide private]

test()

source code

Main filter program for convolution-based filters

tuple

filter_signal(x, y, ftype, f0=None, fn=None, step=1, x_template=None, **kwargs)
Filter a signal.

source code

Basic functions for convolution-based filters

(float,float)

_fixed_window(window_width=0)
Define general window

source code

Convolution-based filter kernels and Windows

(float,float)

gauss_window(sigma=1., limit=4.)
Define Gaussian_window

source code

gauss_kernel(times, signal, t, sigma=1., index=(0,-1), norm_weights=True)
Define Gaussian kernel.

source code

float

mad(a, c=0.6745)
Median Absolute Deviation of a 1D array:

source code

(float,float)

inl_window(**kwargs)
Defines the INL window

source code

inl_kernel(times_in, signal_in, t, c=0.6745, sig_level=3., tolerance=0.01, index=(0,-1), window_width=None)
Iterative Nonlinear Filter (Aigrain).

source code

(float,float)

pijpers_window(delta=1.)
Defines the window for the Pijpers filter.

source code

pijpers_kernel(times_in, signal_in, t, limit=0.000001, delta=1., sigma=1., gamma=0, r=0.0001, index0=0, indexn=-1, index=(0,-1), norm_weights=True)
Defines the Pijpers (2006) filter kernel.

source code

(float,float)

box_window(**kwargs)
Defines the window for the box-filter.

source code

box_kernel(times_in, signal_in, t, window_width=None, index=(0,-1), norm_weights=True)
Defines the Box filter kernel (moving average).

source code

Variables

[hide private]

logger = logging.getLogger("TS.FILTERING")

Function Details

[hide private]

filter_signal(x, y, ftype, f0=None, fn=None, step=1, x_template=None, **kwargs)

source code

Filter a signal.

Handles equidistant and non-equidistant data. x coordinate array must be monotonically increasing, but does not need to be strict monotonically increasing (e.g. their can be duplicate x coordinates).

It is possible to define a new template x coordinate array on which to filter to signal y.

For some filters, it is possible to set an adjustable window: that is, you can give an array to the window parameter, but it needs to be the same length as x_template.

Example usage:

Generate some data:

>>> import time
>>> x = np.linspace(0,150,10000)
>>> y = np.sin(2*np.pi/20.*x)+np.random.normal(size=len(x))

Apply some filter:

>>> c0 = time.time()
>>> x1,y1,pnts = filter_signal(x,y,"gauss",sigma=1)
>>> print time.time()-c0
>>> c0 = time.time()
>>> x1,y1,pnts = filter_signal(x,y,"gauss",sigma=1,threads=2)
>>> print time.time()-c0

Output parameter pnts shows you how many points of the original signal were used to compute the given point. E.g. for a box filter, it will give the number of points over which was averaged.

Parameters:

ftype (string) - one of 'gauss','pijpers','box','inl'

Returns: tuple

output from the used filter: typically (x, y, pnts)

Decorators:

@defaults_filtering
@parallel_pergram
@make_parallel

_fixed_window(window_width=0)

source code

Define general window

Returns: (float,float): window limits

gauss_window(sigma=1., limit=4.)

source code

Define Gaussian_window

Returns: (float,float): window limits

gauss_kernel(times, signal, t, sigma=1., index=(0,-1), norm_weights=True)

source code


Define Gaussian kernel.

#Import necessary modules:
    #>>> from pylab import plot,figure,title,subplot
    #>>> from sigproc.timeseries import tfreq

#Generate some data:
    #>>> x = linspace(0,150,10000)
    #>>> y = random.normal(size=len(x))

#Apply Gaussian filter twice, subtract and calculate periodogram:
    #>>> y1,pnts = filter_signal(x,y,"gauss",sigma=1)
    #>>> y2,pnts = filter_signal(x,y,"gauss",sigma=3)

#Plot the results:
    #>>> p=figure(figsize=(8,18))
    #>>> p=subplot(311);p=title("GAUSSIAN FILTER")
    #>>> p=plot(x,y,'ko')
    #>>> p=plot(x,y1,'r-',linewidth=2)
    #>>> p=plot(x,y2,'b-',linewidth=2)
    #>>> p=subplot(312)
    #>>> p=plot(x,y1-y2,'ko')
    #>>> p=subplot(313)
    #>>> p=plot(per1[0],per1[1],'k-')
    #>>> p=plot(per2[0],per2[1],'r-')

@rtype: (ndarray,)

mad(a, c=0.6745)

source code

Median Absolute Deviation of a 1D array:

median(abs(a - median(a))) / c

Returns: float

inl_window(**kwargs)

source code

Defines the INL window

Returns: (float,float): window limits

inl_kernel(times_in, signal_in, t, c=0.6745, sig_level=3., tolerance=0.01, index=(0,-1), window_width=None)

source code


Iterative Nonlinear Filter (Aigrain).

Example usage:

#Import necessary modules:
    #>>> from pylab import plot,figure,title,subplot
    #>>> from sigproc.timeseries import tfreq

#Generate some data:
    #>>> x = linspace(0,150,10000)
    #>>> y = sin(2*pi*x/20.*x)+random.normal(size=len(x))

#Apply INL filter twice, subtract and calculate periodogram:
    #>>> y1,pnts1,mads = filter_signal(x,y,"inl",window_width=1)
    #>>> y2,pnts2,mads = filter_signal(x,y,"inl",window_width=10)
    #>>> freq1,per1,stat1 = tfreq.scargle(x,y,norm='amplitude',fn=4)
    #>>> freq2,per2,stat2 = tfreq.scargle(x,y1-y2,norm='amplitude',fn=4)

#Plot the results:
    #>>> p=figure(figsize=(8,18))
    #>>> p=subplot(311);p=title("INL FILTER")
    #>>> p=plot(x,y,'ko')
    #>>> p=plot(x,y1,'r-',linewidth=2)
    #>>> p=plot(x,y2,'b-',linewidth=2)
    #>>> p=subplot(312)
    #>>> p=plot(x,y1-y2,'ko')
    #>>> p=subplot(313)
    #>>> p=plot(per1[0],per1[1],'k-')
    #>>> p=plot(per2[0],per2[1],'r-')

@rtype: (ndarray,ndarray)
@return: continuum,sigma

pijpers_window(delta=1.)

source code

Defines the window for the Pijpers filter.

Returns: (float,float): window limits

pijpers_kernel(times_in, signal_in, t, limit=0.000001, delta=1., sigma=1., gamma=0, r=0.0001, index0=0, indexn=-1, index=(0,-1), norm_weights=True)

source code


Defines the Pijpers (2006) filter kernel.

Equals box-multiplying in the frequency domain.

Example usage:

#Import necessary modules:
    #>>> from pylab import plot,figure,title,subplot
    #>>> from sigproc.timeseries import tfreq

#Generate some data:
    #>>> x = linspace(0,150,10000)
    #>>> y = random.normal(size=len(x))

#Apply Gaussian filter twice, subtract and calculate periodogram:
    #>>> y1,pnts1 = filter_signal(x,y,"pijpers",delta=1)
    #>>> y2,pnts2 = filter_signal(x,y,"pijpers",delta=3)
    #>>> freq1,per1,stat1 = tfreq.scargle(x,y,norm='amplitude',fn=4)
    #>>> freq2,per2,stat2 = tfreq.scargle(x,y1-y2,norm='amplitude',fn=4)

#Plot the results:
    #>>> p=figure(figsize=(8,18))
    #>>> p=subplot(311);p=title("PIJPERS FILTER")
    #>>> p=plot(x,y,'ko')
    #>>> p=plot(x,y1,'r-',linewidth=2)
    #>>> p=plot(x,y2,'b-',linewidth=2)
    #>>> p=subplot(312)
    #>>> p=plot(x,y1-y2,'ko')
    #>>> p=subplot(313)
    #>>> p=plot(per1[0],per1[1],'k-')
    #>>> p=plot(per2[0],per2[1],'r-')

@rtype: (ndarray,)
@return: convolved signal

box_window(**kwargs)

source code

Defines the window for the box-filter.

Returns: (float,float): window limits

box_kernel(times_in, signal_in, t, window_width=None, index=(0,-1), norm_weights=True)

source code


Defines the Box filter kernel (moving average).

Equals sinc-multiplication of frequency domain.

Example usage:

#Import necessary modules:
    #>>> from pylab import plot,figure,title,subplot
    #>>> from sigproc.timeseries import tfreq

#Generate some data:
    #>>> x = linspace(0,150,10000)
    #>>> y = random.normal(size=len(x))

#Apply Gaussian filter twice, subtract and calculate periodogram:
    #>>> y1,pnts = filter_signal(x,y,"box",window_width=1)
    #>>> y2,pnts = filter_signal(x,y,"box",window_width=3)

@rtype: (ndarray,)
@return: convolved signal