前提
pythonのscipy.signalで、ある信号と機械特有のフィルターを使って畳み込み積分、逆畳み込み積分をしようとしています。
将来的に目標としているのは、測定された信号から、フィルター関数を逆畳み込みをして本来の信号を再現しようと試みています。
なのでテストとしてそれを行なっているのですが、うまくいきません。
実現したいこと
ある信号があったとします。
この信号をフィルター関数で畳み込み積分を行なって、フィルターにかけられた信号を得ます。
本来の測定に近づけるため、この畳み込まれた信号にノイズを含ませます。
ノイズを含めた信号をフィルターで逆畳み込みをして本来の信号に近づくか調べています。
しかし、ここで元の信号に近づくと思ったのですが、計算してみると発散してうまくいきません。
ノイズを含めない状態で逆畳み込みをすると完全に元に戻ることは確認済みです。
該当のソースコード
python
1import numpy as np 2import matplotlib.pyplot as plt 3import random 4from scipy import signal 5 6#元の信号の生成 7x = np.arange(-50,50,0.1) 8sigma = 0.7 9y_raw = [] 10for i in range(len(x)): 11 if x[i] <= -45 or 45 <= x[i]: 12 y_raw.append(0) 13 else: 14 y_raw.append(-0.7*np.exp(-x[i]**2 / 2*sigma**2) + 1 + random.uniform(-0.05,0.05)) 15 16#filter関数 17x_filter = np.arange(-6,6,1) 18y_pre_filter = [] 19for i in range(len(x_filter)): 20 y_pre_filter.append(np.exp(-x_filter[i]**2 / 2*0.1**2)) 21 22#正規化 23y_filter = [y_pre_filter[i]/sum(y_pre_filter) for i in range(len(y_pre_filter))] 24 25#畳み込み 26y_conv = np.convolve(y_raw,y_filter,mode='same') 27 28#ノイズ生成(ここがないと、うまく元の信号に戻る) 29random_list = np.random.normal(loc=0,scale=0.005,size=(1,len(y_conv))) 30F = [y_conv[i] * (1+random_list[0][i]) for i in range(len(y_conv))] 31y_conv = F 32 33#逆畳み込み 34y_deconv, _ = signal.deconvolve(y_conv, y_filter) 35 36signal = y_conv 37gauss = y_filter 38n = len(signal)-len(gauss)+1 39s = (len(signal)-n)//2 40 41deconv_res = np.zeros(len(signal)) 42deconv_res[s:len(signal)-s-1] = y_deconv 43y_deconv = deconv_res 44 45 46#プロット 47plt.plot(x,y_raw,label='obs')#本来の信号 48plt.plot(x,y_conv,label='conv')#フィルターで畳み込まれ、ノイズが入った信号 49plt.plot(x,y_deconv,label='deconv')#ノイズ信号を逆畳み込みした信号 50plt.legend() 51plt.show() 52
試したこと
いろんなサイトを調べてみて実装してみたのですが、どうもうまくいきません。
何かご存じでしたら、ご教授いただけないでしょうか。
補足情報
python 3.7.9, scipy 1.5.2, numpy 1.19.2
> ノイズを含めない状態で逆畳み込みをすると完全に元に戻ることは確認済み
ノイズ有りでの処理後のデータが、
https://jp.mathworks.com/help/images/deblurring-images-using-a-wiener-filter.html
の「動きによるブレとガウスノイズのシミュレーションと復元」の、「ノイズ推定を指定せずに deconvwnr を使用して...」と書かれてるところの次の、ノイズだらけの画像みたいになってるのではありませんか?
上記Webページのコードはpythonではなくmatlabのものなので、コードは気にせず文章と画像だけ見てください
おそらくその画像と似た状況になっています。
私の場合もノイズ含めたものの逆畳み込みは、元の信号の再現には到底及ばない状況です。
ノイズを含む信号のデコンボリューションは、単純な逆フィルターではノイズが増幅されるため、うまくいきません
Wiener Filter
Richardson-Lucy
Total Variation
あたりを調べてみてください

回答1件
あなたの回答
tips
プレビュー