2次元画像を点広がり関数(PSF)でデコンボリューションしたいと思います。私は、1次元配列のために働くscipy.signal.deconvolve
関数があり、多次元配列を畳み込むのにscipy.signal.fftconvolve
があることを見ました。 scipyに2D配列のデコンボリューションを行うための特定の関数はありますか?2次元アレイに対応するscipy.signal.deconvolveはありますか?
Iはdivistionによってfftconvolveに製品を交換fftdeconvolve関数を定義した:
def fftdeconvolve(in1, in2, mode="full"):
"""Deconvolve two N-dimensional arrays using FFT. See convolve.
"""
s1 = np.array(in1.shape)
s2 = np.array(in2.shape)
complex_result = (np.issubdtype(in1.dtype, np.complex) or
np.issubdtype(in2.dtype, np.complex))
size = s1+s2-1
# Always use 2**n-sized FFT
fsize = 2**np.ceil(np.log2(size))
IN1 = fftpack.fftn(in1,fsize)
IN1 /= fftpack.fftn(in2,fsize)
fslice = tuple([slice(0, int(sz)) for sz in size])
ret = fftpack.ifftn(IN1)[fslice].copy()
del IN1
if not complex_result:
ret = ret.real
if mode == "full":
return ret
elif mode == "same":
if np.product(s1,axis=0) > np.product(s2,axis=0):
osize = s1
else:
osize = s2
return _centered(ret,osize)
elif mode == "valid":
return _centered(ret,abs(s2-s1)+1)
しかし、以下のコードは、畳み込みとデコンボリューション後に元の信号を回復しない。
sx, sy = 100, 100
X, Y = np.ogrid[0:sx, 0:sy]
star = stats.norm.pdf(np.sqrt((X - sx/2)**2 + (Y - sy/2)**2), 0, 4)
psf = stats.norm.pdf(np.sqrt((X - sx/2)**2 + (Y - sy/2)**2), 0, 10)
star_conv = fftconvolve(star, psf, mode="same")
star_deconv = fftdeconvolve(star_conv, psf, mode="same")
f, axes = plt.subplots(2,2)
axes[0,0].imshow(star)
axes[0,1].imshow(psf)
axes[1,0].imshow(star_conv)
axes[1,1].imshow(star_deconv)
plt.show()
得られた2D配列は、下の図の下の行に示されています。 FFTデコンボリューションを使用して元の信号を復元するにはどうすればよいですか? scipyのダウンロードのFFTPACKパッケージからfftn、関数ifftn場合、fftshiftとifftshiftを使用してこれらの機能が動作するように見える
水域を持つfドメインの分割について、 'spec_out = spec_signal/spec_to_deconv'ではなく、' spec_out = spec_signal /(spec_to_deconv + waterlevel) 'とはどういう意味ですか? 'spec_to_deconv'の小さな値の摂動を避け、' spec_out'に大きな変化をもたらします。 – Lee
これは基本的にウィーナーフィルタの背後にあるアイデアです。それは最高の激しい物理的正当性を有しており、閾値を設定することは、常に、雑音抑圧と導入される歪みとの間の妥協を伴うことになる。 –
RL繰り返しは1d時系列にも適用できますか? – Lee