2013-08-15 4 views
4

파이썬에서 적응 형 필터를 적용하고 싶습니다. 그러나 그러한 알고리즘을 구현하는 방법에 대한 온라인 설명서는 찾을 수 없습니다. 저는 scipy.signal 도구 상자를 사용하여 "정적"필터를 설계하는 것에 익숙하지만 적응 필터를 설계하는 방법을 모릅니다.Python에서 적응 필터를 적용하는 방법

명확히하기 : 잡음이 포함 된 녹음 신호 S가 있습니다. 이 녹음에는 액세스하고 싶은 "참"기능이 있는데 T이라고합니다. 예상치도 T입니다. 필터링 된 ST 사이의 오류가 최소화되도록 필터를 설계하고 싶습니다. 이 경우 정적 필터는 비 정적 신호를 필터링하려고하므로 유용하지 않습니다.

+2

적응 형 신호 처리를위한 두 개의 라이브러리가 있는데 첫 번째는 [adaptfilt] (https://pypi.python.org/pypi/adaptfilt/0.2)이고 두 번째 라이브러리는 [padasip] (https : // pypi. python.org/pypi/padasip) –

답변

7

여기 Numpy가있는 Python의 기본 LMS 적응 필터가 있습니다.
의견을 환영합니다. 대부분의 테스트 사례를 환영합니다.

enter image description here

""" lms.py: a simple python class for Least mean squares adaptive filter """ 

from __future__ import division 
import numpy as np 

__version__ = "2013-08-29 aug denis" 

#............................................................................... 
class LMS: 
    """ lms = LMS(Wt, damp=.5) Least mean squares adaptive filter 
    in: 
     Wt: initial weights, e.g. np.zeros(33) 
     damp: a damping factor for swings in Wt 

    # for t in range(1000): 

    yest = lms.est(X, y [verbose=]) 
    in: X: a vector of the same length as Wt 
     y: signal + noise, a scalar 
     optional verbose > 0: prints a line like "LMS: yest y c" 
    out: yest = Wt.dot(X) 
     lms.Wt updated 

    How it works: 
    on each call of est(X, y)/each timestep, 
    increment Wt with a multiple of this X: 
     Wt += c X 
    What c would give error 0 for *this* X, y ? 

     y = (Wt + c X) . X 
     => 
     c = (y - Wt . X) 
      -------------- 
       X . X 

    Swings in Wt are damped a bit with a damping factor a.k.a. mu in 0 .. 1: 
     Wt += damp * c * X 

    Notes: 
     X s are often cut from a long sequence of scalars, but can be anything: 
     samples at different time scales, seconds minutes hours, 
     or for images, cones in 2d or 3d x time. 

""" 

# See also: 
#  http://en.wikipedia.org/wiki/Least_mean_squares_filter 
#  Mahmood et al. Tuning-free step-size adaptation, 2012, 4p 
# todo: y vec, X (Wtlen,ylen) 

#............................................................................... 
    def __init__(self, Wt, damp=.5): 
     self.Wt = np.squeeze(getattr(Wt, "A", Wt)) # matrix -> array 
     self.damp = damp 

    def est(self, X, y, verbose=0): 
     X = np.squeeze(getattr(X, "A", X)) 
     yest = self.Wt.dot(X) 
     c = (y - yest)/X.dot(X) 
      # clip to cmax ? 
     self.Wt += self.damp * c * X 
     if verbose: 
      print "LMS: yest %-6.3g y %-6.3g err %-5.2g c %.2g" % (
       yest, y, yest - y, c) 
     return yest 

#............................................................................... 
if __name__ == "__main__": 
    import sys 

    filterlen = 10 
    damp = .1 
    nx = 500 
    f1 = 40 # chirp 
    noise = .05 * 2 # * swing 
    plot = 0 
    seed = 0 

    exec("\n".join(sys.argv[1:])) # run this.py n= ... from sh or ipython 
    np.set_printoptions(2, threshold=100, edgeitems=10, linewidth=80, suppress=True) 
    np.random.seed(seed) 

    def chirp(n, f0=2, f1=40, t1=1): # <-- your test function here 
     # from $scipy/signal/waveforms.py 
     t = np.arange(n + 0.)/n * t1 
     return np.sin(2*np.pi * f0 * (f1/f0)**t) 

    Xlong = chirp(nx, f1=f1) 
    # Xlong = np.cos(2*np.pi * freq * np.arange(nx)) 
    if noise: 
     Xlong += np.random.normal(scale=noise, size=nx) # laplace ... 
    Xlong *= 10 

    print 80 * "-" 
    title = "LMS chirp filterlen %d nx %d noise %.2g damp %.2g " % (
     filterlen, nx, noise, damp) 
    print title 
    ys = [] 
    yests = [] 

#............................................................................... 
    lms = LMS(np.zeros(filterlen), damp=damp) 
    for t in xrange(nx - filterlen): 
     X = Xlong[t:t+filterlen] 
     y = Xlong[t+filterlen] # predict 
     yest = lms.est(X, y, verbose = (t % 10 == 0)) 
     ys += [y] 
     yests += [yest] 

    y = np.array(ys) 
    yest = np.array(yests) 
    err = yest - y 
    averr = "av %.2g += %.2g" % (err.mean(), err.std()) 
    print "LMS yest - y:", averr 
    print "LMS weights:", lms.Wt 
    if plot: 
     import pylab as pl 
     fig, ax = pl.subplots(nrows=2) 
     fig.set_size_inches(12, 8) 
     fig.suptitle(title, fontsize=12) 
     ax[0].plot(y, color="orangered", label="y") 
     ax[0].plot(yest, label="yest") 
     ax[0].legend() 
     ax[1].plot(err, label=averr) 
     ax[1].legend() 
     if plot >= 2: 
      pl.savefig("tmp.png") 
     pl.show() 
+1

이 훌륭한 스크립트를 가져 주셔서 감사합니다. 나는 하나의 질문을 가지고있다 : 적응 필터에 원래 함수를 추정 할 수 있을까? 귀하의 예제 코드에서는 그냥 몇 가지 스칼라 변수 y를 전달 ...하지만 신호의 추정을 사용하여 함수에서 신호를 추출하는 적응 필터를 사용하려면 말하고 싶습니다. 스크립트로이 작업을 수행하는 방법이 분명하지 않습니다. – allhands

+1

@allhands, 플롯 및 테스트 코드의'yest'는 실제 신호 (여기에서는 삑 소리), 실제 잡음 (정상), X (이전 filterlen = 10 신호 + 잡음 입력)를 기반으로 한 신호 y의 추정치입니다. , 댐핑 팩터. 'yest'와'y'가 아주 가깝기 때문에 짹짹 소리 + 잡음에서 chirp를 추출하는 예제는 좋지 않습니다. 더 좋은 예가 필요하다. – denis

4

당신은 https://pypi.python.org/pypi/padasip

여러 적응 필터는이 라이브러리에서 구현되는 Padasip 라이브러리를 확인해야합니다. 당신은 당신의 사건에 대한 구현 적응 필터의 적용 가능성을 테스트하려면 당신은 직접 사용하거나 (소스 코드는 GitHub의 https://github.com/matousc89/padasip에) 필터가이 만들어지는 방법, 조사 할 수

,이 튜토리얼을 볼 Padasip Adaptive Filters Basics - Noise Cancelation, System Identification and Signal Prediction

관련 문제