2012-04-04 15 views
5

私は、Weierstrass変換を使用して与えられた信号を平滑化するためのPythonコードを作っています。これは、基本的に正規化ガウス関数と信号の畳み込みです。次のようにscipy/numpy fftでゼロ埋め込みのために発生する境界効果を削除するにはどうすればよいですか?

コードは次のとおりです。


#Importing relevant libraries 
from __future__ import division 
from scipy.signal import fftconvolve 
import numpy as np 

def smooth_func(sig, x, t= 0.002): 
    N = len(x) 
    x1 = x[-1] 
    x0 = x[0]  


# defining a new array y which is symmetric around zero, to make the gaussian symmetric. 
    y = np.linspace(-(x1-x0)/2, (x1-x0)/2, N) 
    #gaussian centered around zero. 
    gaus = np.exp(-y**(2)/t)  

#using fftconvolve to speed up the convolution; gaus.sum() is the normalization constant. 
    return fftconvolve(sig, gaus/gaus.sum(), mode='same') 

私はステップ関数と言うためにこのコードを実行する場合は、それがコーナーを平滑化するが、境界で、それは別のコーナーを解釈していることを平滑化しますその結果、境界で不必要な動作が発生します。私はこれを以下のリンクの図で説明します。
Boundary effects

この問題は、畳み込みを見つけるために直接統合すると発生しません。したがって、問題はWeierstrass変換ではないため、問題はscipyのfftconvolve関数にあります。

なぜこの問題が発生するのか理解するためには、最初にscipyでfftconvolveの働きを理解する必要があります。
fftconvolve関数は基本的に畳み込み定理を使用して計算を高速化します。

畳み込み(int1、int2)= ifft(fft(int1)* fft(int2))
この定理を直接適用すると、結果が得られません。目的の結果を得るには、fftをmax(int1、int2)の2倍の大きさにする必要があります。しかし、これは望ましくない境界効果をもたらす。これは、fftコードでsize(int)がsize(fftを取るべき大きさ)より大きい場合、入力をゼロにしてからfftを取るためです。このゼロ埋め込みは、望ましくない境界効果の原因となるものです。

この境界効果を削除する方法を提案できますか?

私は単純なトリックで削除しようとしました。関数を平滑化した後、平滑化された信号の値と境界付近の元の信号とを比較し、一致しなければ平滑化された関数の値をその点の入力信号と置き換えます。
以下の通りである:


i = 0 
eps=1e-3 
while abs(smooth[i]-sig[i])> eps: #compairing the signals on the left boundary 
    smooth[i] = sig[i] 
    i = i + 1 
j = -1 

while abs(smooth[j]-sig[j])> eps: # compairing on the right boundary. 
    smooth[j] = sig[j] 
    j = j - 1 

以下に示すように、なぜなら平滑化関数の小さなジャンプが存在するイプシロンを使用する、この方法に問題がある:
jumps in the smooth func

この境界問題を解決するために上記の方法を変更することはできますか?

+0

デュエのhttp://math.stackexchange.com/q/127875/2206 – endolith

答えて

3

最後に対称フィルタカーネルが生成するのは、データが最後を超えていると仮定しているかどうかによって決まります。

現在の結果の外観が気に入らない場合は、両端を超えてゼロが仮定されます。データの反映や多項式回帰継続などの別の仮定でデータを拡張してみてください。 (拡張子が0の場合を除き、非円形の畳み込みに必要な既存のゼロパディングではなくなります)、フィルタカーネルの長さの半分以上でデータを両端で拡張します。その後、フィルタリング後に追加されたエンド拡張を削除し、あなたの前提の外観が好きかどうかを確認してください。そうでない場合は、別の前提をお試しください。または、実際のデータを使用している場合は、端を超えて実際のデータを使用することをお勧めします。

+0

お返事ありがとうございます。ご返信いただければ、さまざまな可能性が開かれます。 – Omkar

6

ベストアプローチはmode = 'valid'を使用することが考えられます:

The output consists only of those elements that do not rely on the zero-padding.

あなたの信号をラップすることができ、または処理される信号は大きな信号(場合からの抜粋である場合を除き:プロセスフル信号次に関心のある領域を切り抜く)、畳み込みを行うときに常にエッジ効果を持つことになります。どのように対処するかを選択する必要があります。 mode = validを使用すると、それらは切り取られます。これはかなり良い解決策です。信号が常に「ステップ状」であることがわかっている場合は、処理された信号のフロントとエンドを適切に延長することができます。

+0

ありがとうございました。 – Omkar

関連する問題