Poisson Image Editing を画像ピラミッドで高速化する
概要
本稿では， Poisson Image Editing (Perez et al., 2003) を Python で実装してみます．ヤコビ法で高解像度画像を処理すると時間がかかるところ，あらかじめ画像ピラミッドを作成して低解像度側から順に処理することで処理時間を短縮できました．
考え方はマルチグリッド法と似ているところもあります（マルチグリッド法のほうが洗練されています）が，少しの書き換えだけで実装できるので，小難しいことを考えたくないときには便利です．
方法
Poisson Image Editing
Poisson Image Editing (Perez et al., 2003) （以降 PIE）は，画像の一部を別の画像で自然に置き換えたいときに使う手法です．この手法では，置き換えを行う領域に対し，(1) 継ぎ目の画素の値が元の画像と等しくなる， (2) 勾配は移植したい画像の勾配と等しくなる という条件を満たす解を求めます．このときにポアソン方程式を解くので，この名前がついています．
左から，合成先，合成元，マスク，合成結果
画像ピラミッド
画像ピラミッドとは，ある画像を繰り返し処理して得られる画像の組です．今回は，画像を繰り返しダウンサンプルして得られる最も一般的な画像ピラミッドを使いました．同じ処理を画像ピラミッドの各画像に対して適用することで，画像の異なる空間周波数成分に対して処理を適用することができます．
画像ピラミッドの中身
処理の流れ
まず，ピラミッドの中の最も低解像度の画像に対して PIE を行います．次に，得られた画像を1段階アップサンプルしたものを初期値として PIE を行います．これを目的の解像度まで繰り返します．
画像ピラミッドの作成
せっかく np.float64 で計算するので，ついでに簡易なガンマ補正処理を入れました．
def create_image_pyramid(image, mask, n_levels, gamma):
masks = []
images = []
for i in range(n_levels):
# mask
resampled_mask = np.copy(mask[::(2**i), ::(2**i)])
masks.append(resampled_mask)
# resample and degamma image
resampled_image = np.array(
image.resize((resampled_mask.shape[1], resampled_mask.shape[0]), Image.BOX)
).astype(np.float64) / 255
resampled_image_degammaed = np.power(resampled_image, 1/gamma)
images.append(resampled_image_degammaed)
return images, masks
PIE
PIE の中心部分は以下のように実装できます（参考にした講義動画）．
def apply_laplacian_filter(a):
return a[1:-1,:-2] + a[1:-1,2:] + a[2:,1:-1] + a[:-2,1:-1] - a[1:-1, 1:-1] * 4
def do_poisson_edit_per_channel(dst, src, mask, n, err_threshold):
max_err = np.nan
dst_ = np.copy(dst)
src_ = np.copy(src)
for i in range(n):
dst_new = (
apply_laplacian_filter(dst_) / 4
+ dst_[1:-1,1:-1]
- apply_laplacian_filter(src_)/4
)
max_err = np.max(
dst_[1:-1,1:-1][mask[1:-1,1:-1]]
- dst_new[mask[1:-1,1:-1]]
)
dst_new_padded = np.pad(dst_new, 1, mode="edge")
dst_[mask] = dst_new_padded[mask]
if (max_err < err_threshold):
break
return dst_, max_err, i
def format_info(result, channel):
return (
channel + " "
+ "Error: {:.5f}".format(result[1]) + " "
+ "Iter: " + str(result[2]).ljust(8)
)
def do_poisson_edit(dst, src, mask, n=10_000, err_threshold = 1e-3):
r_result = do_poisson_edit_per_channel(dst[:,:,0], src[:,:,0], mask, n, err_threshold)
print(format_info(r_result, "R"), end="")
g_result = do_poisson_edit_per_channel(dst[:,:,1], src[:,:,1], mask, n, err_threshold)
print(format_info(g_result, "G"), end="")
b_result = do_poisson_edit_per_channel(dst[:,:,2], src[:,:,2], mask, n, err_threshold)
print(format_info(b_result, "B"))
return np.dstack([r_result[0], g_result[0], b_result[0]])
画像ピラミッドとのデータのやりとり
max_iter = 10_000
err_threshold = 1e-2
gamma = 2.2
n_levels = 7
# initialise
dst_images, masks = create_image_pyramid(dst_image, mask, n_levels, gamma)
src_images, _ = create_image_pyramid(src_image, mask, n_levels, gamma)
results = [None for i in range(len(dst_images))]
# compute
for i in range(n_levels):
level = n_levels - 1 - i
s = src_images[level]
d = np.copy(dst_images[level])
m = masks[level]
if (i == 0):
results[level] = do_poisson_edit(d, s, m, n=max_iter, err_threshold=err_threshold)
else:
# upsample
r = Image.fromarray((results[level+1][:,:,0])).resize((d.shape[1], d.shape[0]), Image.BILINEAR)
g = Image.fromarray((results[level+1][:,:,1])).resize((d.shape[1], d.shape[0]), Image.BILINEAR)
b = Image.fromarray((results[level+1][:,:,2])).resize((d.shape[1], d.shape[0]), Image.BILINEAR)
# copy to next array
prev_d = np.dstack([np.array(r), np.array(g), np.array(b)]).astype(np.float64)
d[:,:,0][m] = prev_d[:,:,0][m]
d[:,:,1][m] = prev_d[:,:,1][m]
d[:,:,2][m] = prev_d[:,:,2][m]
# solve
results[level] = do_poisson_edit(d, s, m, n=max_iter, err_threshold=err_threshold)
結果
Google Colab を使った結果，上記パラメター（最大計算回数 10,000 回，画像ピラミッド7段，誤差の閾値 1e-2）で 900 x 600 の画像に対して 10 秒程度で計算を終えられました（1枚目）．画像ピラミッドなしだと6秒程度で計算が終了しましたが，自然な合成結果は得られませんでした（2枚目）．画像ピラミッドなしで結果をうまく収束させるためには，誤差を 1e-4 程度に抑える必要がありました（4枚目）．このとき，時間は5分以上かかっています．
このコードではマスクの部分だけ切り出して処理をするなどの工夫はせず，画像全体にラプラシアンフィルタを適用しています．マスク部分だけ切り出して処理をすれば，計算時間はもう少し短くできそうです．
画像は Unsplash から以下2点を使用しました．なお，鳥の画像はサイズ合わせの都合上，エッジの色で塗り足しました．
Discussion