ChatGPTと壁打ちして、スマホにあるボケ強調機能であるポートレートモードをOpenCVで実装してみた。エッジをほとんど滲ませないで既にぼけている部分だけをぼかすためのEdge-Contained Pyramid Blurという新規手法を考案し、実装した。画像認識や空間認識の手法を使って前景と背景の区別をすることなしに、画像一枚に統計処理を適用するだけで実現でき、しかもアーティファクトがほとんど視認できない自然な結果が得られる。そのうえで、前景認識のGrabCutという手法と組み合わせて、より安定して動く実装を完成させた。
ポートレートモードとは
上述の写真は、ポートレートモードをかけた後のものだ。加工前の元画像は下のものだ。ぼけの感じが全く違うのがわかるだろう。この作例では、背景の葉っぱの乱雑さが明らかに注意散漫（distraction）の原因になっているが、ぼかすことでそれが解消され、一段上の作品になる。絞りを開いて被写界深度を狭めると花弁がぼけてしまうので、後処理で背景をぼかすことでのみこの結果に辿り着く。
人物写真のことをポートレートと言う。そこでは、典型的には、人物の顔（特に目）にピントを合わせて、そこに閲覧者の視線を誘導する。逆に、背景はよくぼけさせて、背景に気を取られないようにする。そのためには、絞りが大きいレンズや大型のセンサーのカメラなど、背景がよくボケる機材を使う必要がある。一方で、スマホは極小のレンズと極小のセンサーしか積んでいないので、背景をぼかしたポートレート的な表現をすることは難しい。そこで、スマホ各社はポートレートモードなる機能を提供して、撮影後の後処理で背景にぼかし加工を適用できるようにしている。今回はそれを作ってみようという話だ。
まずは、既にできている実装の使い方を説明したい。itb_stack.pyをインストールしてから、以下のようにコマンドを一発叩けば完了する。つまり、最初の引数で入力ファイルを指定し、-o（または--output）オプションで出力ファイルを指定し、--portraitオプションでぼかしの度合いを指定する。autoの場合は、画像サイズに合わせて適当に決めてくれる。
$ itb_stack.py input.jpg -o output.jpg --portrait=auto
こんな入力画像があったとしよう。ピントは手前のカレーの皿に合っていて背景のナンの皿はぼけている。この例でカレーに視線を集める効果を増すには、もっと背景をぼかしたくなる。
コマンドにかけると、こんな風に背景ボケを強調してくれる。カレーのシャープネスはほとんど落ちていないが、背景のナンがとろけて、ちょっと良いレンズで撮ったような効果が生まれている。
ポートレートモードでは、前景と背景のボケ量の違いを見て前景を判断するので、前景がはっきりしている写真だとうまく動く。言い換えると、ある程度は背景がボケた写真において、ボケ量を強調するのが仕事だ。
逆に言えば、もともとボケ量が少なくて、前景も背景もシャープネスが大して変わらない場合には、うまく動かない。この例だと、主要被写体であるバイクの一部も溶けてしまっている。
前景と背景が分離できていないからこそポートレートモードが使いたいのに、そのような場合にはうまく動かない。これは本質的な矛盾だが、統計的画像処理のアプローチではこれが限界だ。前景と背景がある程度分離できている状態で、その分離性を強調するツールとしてなら実用になる。
ガウシアンピラミッドとは
複雑な問題に取り組む際には、多くの場合、いきなり成功を目指すのは現実的ではない。素朴で分かりやすい方法を試してみて、敢えて失敗し、その限界を知ることから始めるのが良い。今回も、少しずつそれをやっていこう。まずは、入力画像を用意する。以下に示すように、白いノギスの手前側にピントが合っていて、奥のごちゃごちゃした置物は被写界深度から外れている。
これに何も考えずにガウシアンぼかしをかけると、当然ながら全体がぼける。ガウシアンぼかし半径5画素（直径11画素）の結果が以下だ。ガウシアンぼかしはカーネル関数なので、カーネルサイズの2乗の計算量がかかる。例えばカーネルサイズ11にすると、各画素に対して11^2=121画素の情報を見ることになる。カーネルサイズは21ぐらいまでが実用的な限界だろう。ここに大きな矛盾がある。大きな画像ほど大きなカーネルサイズにして大きくぼかしたいのに、大きな画像ほど画素数が多いから、カーネルサイズを大きくすると遅すぎて終わらなくなる。
ある程度以上大きな画像をぼかす際には、ガウシアンピラミッドを使う。ざっくり言うと、画像を縦横半分ずつ縮小してから縦横2倍に拡大すれば情報が1/4になってボケるという手法だ。それを1回やるならレベル1だ。折り畳みを2回やればレベル2だ。レベル1なら情報量は1/4になり、レベル2なら1/16になり、レベル3なら1/64になるという指数関数で強度が決まる。それでいてレベルの対数の計算量にすぎないので、どんな大きなぼけも作れる。レベル5にするとこんなにボケる。
ぼかしを強くかけると、画素を混ぜるので、コントラストが下がる。特にガウシアンピラミッドは一度に大量の画素を混ぜるので、コントラストの低下が激しい。よって、縮小した画像を復元する際に、コントラストを調整する実装を入れた。具体的には、同じ階層でぼかす前後の画像の輝度の標準偏差を取り、その比率を係数として線形補正を行う。そうすると、以下のようにコントラストが保たれる。
ガウシアンピラミッドはレベル指定だけだと強度の刻みが大雑把すぎるので、微調整する機構も作った。decayというパラメータで、各階層でぼかしがどの程度効くのかの重みを変えられるようにした。レベル5かつdecay=0.3だと、ぼけがこのくらい減る。レベルが同じであれば最大ぼけ半径は同じなのだが、decayを上げるとぼけの情報分布を中心に寄せるので、ぼけに芯が見られるようになる。
ガウシアンピラミッドは画像の縦と横の両方が2のレベル数乗の倍数でないとうまく動かない（実際には補完処理で動くのだが、各階層でサイズの整合性を取るのが面倒になる）。よって、入力画像がそのサイズに整合するように拡張してから、ガウシアンピラミッドによる処理を行って、その結果を元の画像のサイズにトリミングする。
以上をコードに落とすと、こんな感じになる。make_gaussian_pyramidはガウシアンピラミッドを作る関数で、blur_image_pyramid関数はそれを使ってぼかし演算を行う。blur_image_pyramidはECPBと直接の関係はないのだが、こいつの挙動を完璧に把握しておくことがこの先の理解に重要になる。
def make_gaussian_pyramid(image, levels): """ガウシアンピラミッドを作る""" assert image.dtype == np.float32 # ピラミッドの最下層は元画像 pyramid = [image] # レベル数分だけ繰り返す for _ in range(levels): # 今の階層の縮小画像を作ってリストに追加 image = cv2.pyrDown(image) pyramid.append(image) return pyramid def blur_image_pyramid(image, levels, decay=0.0, contrast=1.0): """ピラミッドブラーを画像に適用する""" assert image.dtype == np.float32 # 画像サイズを調べて、最大レベルを調整 h, w = image.shape[:2] levels = min(levels, int(math.log2(min(h, w))) - 1) # 画像の縦横が2^levelsの倍数になるように拡大 factor = 2 ** levels new_h = ((h + factor - 1) // factor) * factor new_w = ((w + factor - 1) // factor) * factor expanded = cv2.copyMakeBorder(image, 0, new_h - h, 0, new_w - w, cv2.BORDER_REPLICATE) # レベル毎にボケ半径が大きくなるので、それを相殺する重みを設定 bokehs = [2 ** i for i in range(0, levels)] alpha_weights = [b / bokehs[-1] for b in reversed(bokehs)] # ガウシアンピラミッドを作成 pyramid = make_gaussian_pyramid(expanded, levels) # ピラミッドの末尾からの復元作業のループ diffused = pyramid[-1] for i in range(levels - 1, -1, -1): # 現在の階層のサイズを取得し、縦横2倍ずつに展開 size = (pyramid[i].shape[1], pyramid[i].shape[0]) diffused = cv2.pyrUp(diffused, dstsize=size) # 加重平均の重みを計算 gamma_decay = decay ** 2 alpha = alpha_weights[i] * gamma_decay + (1 - gamma_decay) # 現在のガウシアンと展開してぼかした画像のコントラストを比較 std_ref = pyramid[i].std() std_diff = diffused.std() gain = std_ref / (std_diff + 1e-6) * contrast gain = np.clip(gain, 1.0, 2.0) # 展開してぼかした画像と現在のガウシアンを加重平均 diffused = alpha * diffused * gain + (1 - alpha) * pyramid[i] # サイズを元に戻して返す trimmed = diffused[:h, :w] return np.clip(trimmed, 0, 1)
ガウシアンピラミッドの肝であるcv2.pyrDownとcv2.pyrUpの実装がどうなっているかを知るのは有効だ。実際にはC++で最適化された実装がなされているが、ともかく畳み込みとガウシアンをやっているので、以下のコードで置き換えても動く。位置がずれないようにアフィン変換をしているのが肝であり、またその順序が、downではリサイズの前、upではリサイズの後にやるのが重要だ。
def pyramid_down(image): # 縦横0.5ピクセルずつずらす tx, ty = 0.5, 0.5 m = np.float32([[1, 0, tx], [0, 1, ty]]) image = cv2.warpAffine(image, m, (image.shape[1], image.shape[0]), flags=cv2.INTER_LINEAR, borderMode=cv2.BORDER_REFLECT101) # 縦横半分にした画像を作る resized = cv2.resize(image, None, fx=0.5, fy=0.5, interpolation=cv2.INTER_LINEAR) # ガウシアンブラーをかける resized = cv2.GaussianBlur(resized, (5, 5), sigmaX=0.45, sigmaY=0.45) return resized def pyramid_up(image): # 縦横2倍にした画像を作る image = cv2.resize(image, None, fx=2.0, fy=2.0, interpolation=cv2.INTER_LINEAR) # 縦横0.5ピクセルずつずらす tx, ty = -0.5, -0.5 m = np.float32([[1, 0, tx], [0, 1, ty]]) image = cv2.warpAffine(image, m, (image.shape[1], image.shape[0]), flags=cv2.INTER_LINEAR, borderMode=cv2.BORDER_REFLECT101) # ガウシアンブラーをかける image = cv2.GaussianBlur(image, (3, 3), sigmaX=0.3, sigmaY=0.3) return image
このピラミッド関数群の実装を知っておくと、いろいろ潰しが効いて便利だ。各階層でぼけた画像を元に戻すための差分情報であるラプラシアンピラミッドを作っておくと、元画像の完全復元もできる。ガウシアンブラーの代わりに別の関数を挿入すれば、広範囲の情報を拾ってゴニョゴニョする処理が簡単に書ける
シャープネスマップとは
大きなぼけが作れるようになったのは良いが、全体をぼかすのがポートレートモードの目的ではない。前景（主要被写体）はぼかさずに、背景だけをぼかしたいのだ。実際には主要被写体が奥にあって手前に前ボケがある構図もあるだろうが、ここでは便宜上、主要被写体のことを前景と呼ぶ。主要被写体にはピントが合っているものなので、言い換えれば、ピントが合っていてシャープな部分はぼかさずに、ピントが外れていて元々ぼけている部分はもっとぼけさせるということでもある。
ならば、画像の中でシャープな部分を検出できれば、その周辺をぼかすことでポートレートモードが成立しそうだ。シャープな部分はコントラストが高い部分であり、コントラストが高いということは輝度の差があるということであり、画素同士の輝度の差が大きい「エッジ」と呼ばれる構造を探せば良いことになる。エッジを数値化する手法はいろいろあるのだが、今回はラプラシアンとソーベルの二つの指標を足して使うことにした。意味も変域も異なる二つの指標を足すにあたっては、双方をZスコアとして平準化してから加算するという手法を取った。その後にさらにZスコアで平準化して、下位2%と上位2%をクリップして[0,1]の範囲に線形変換で収めたものを最終的なスコアとした。各画素のスコアを画像として表したものがシャープネスマップであり、以下のようになる。真っ白な部分が最もシャープで、真っ黒な部分は最もぼけているということになる。
シャープネスマップをうまく使いこなすには、その統計的な特徴を把握することが大事だ。ポートレート写真などの、被写界深度が浅くてぼけが大きい写真では、ピントが合った領域が小さいことが多い。言い換えると、画像の多くの領域はぼけていて、シャープな領域はごく一部に偏るというのが典型的ということだ。そうなると、シャープネスマップの統計値は、Zスコアに変換してもなお、正規分布には程遠い形状になる。平均は0に寄り、標準偏差は低くなり、歪度は正になり、尖度は高くなる。実際、上のシャープネスマップは、min=0.000, max=1.000, mean=0.117, stddev=0.185, skew=3.132, kurt=10.679 という統計値を持つ。これが意味するのは、シャープネスマップを使ってエッジを保護する機能は、そのままだと、画像のごく一部の、エッジのすぐ近くでしか働かないということだ。
シャープネスマップを作るコードを以下に示す。
def compute_sharpness(image, base_area=1000000, blur_radius=2): """Computes sharpness using normalized Laplacian and Sobel filters.""" assert image.dtype == np.float32 gray = cv2.cvtColor(image, cv2.COLOR_BGR2GRAY) h, w = gray.shape[:2] area = h * w is_scaled = False if area > base_area * 2: scale = (base_area / area) ** 0.5 gray = cv2.resize(gray, (math.ceil(w * scale), math.ceil(h * scale)), interpolation=cv2.INTER_AREA) is_scaled = True if blur_radius > 1: ksize = math.ceil(2 * blur_radius) + 1 gray = cv2.GaussianBlur(gray, (ksize, ksize), 0) laplacian = np.abs(cv2.Laplacian(gray, cv2.CV_32F, ksize=3)) laplacian = z_score_normalization(laplacian) sobel_x = np.abs(cv2.Sobel(gray, cv2.CV_32F, 1, 0, ksize=3)) sobel_y = np.abs(cv2.Sobel(gray, cv2.CV_32F, 0, 1, ksize=3)) sobel_e = np.sqrt(sobel_x**2 + sobel_y**2) sobel = z_score_normalization(sobel_e) sharpness = 0.5 * laplacian + 0.5 * sobel if is_scaled: sharpness = cv2.resize(sharpness, (w, h), interpolation=cv2.INTER_LANCZOS4) sharpness = z_score_normalization(sharpness) return np.clip(sharpness, -10, 10)
誰しも思い付くであろうことは、シャープネスマップの値を重み付けにして、元画像とぼけた画像を重み付き平均合成するというアイデアだ。コードにするとこんな感じ。
def merge_images_by_sharpness_weights(original, blurred): sharpness = compute_sharpness(original) weight_map_3d = np.repeat(sharpness[:, :, np.newaxis], 3, axis=2) merged = weight_map_3d * original + (1.0 - weight_map_3d) * blurred return np.clip(merged, 0, 1)
例えばレベル4のぼけと元画像を重み付き合成すると以下のようになる。
すぐ気づくように、この手法は駄目だ。確かにピントが合っている部分のノギスの目盛りは視認できるようになったが、ノギスの本体の白い部分が盛大にぼけて、ハロ（明るいモヤモヤ）として描写されてしまっている。このハロのように、不自然な模様が出てしまうことをアーティファクトという。写真を対象とした画像処理はアーティファクトの戦いなのだ。この失敗の理由は明白であり、既に説明した、シャープネスマップの局所性にある。
Edge-Contained Pyramid Blur
そこで、私とChatGPTが考えたのが、ECPB（Edge-Contained Pyramid Blur）という手法である。名が示す通り、エッジを閉じ込めてピラミッドぼかしをかけるものだ。そもそもなぜハロが発生するかというと、エッジの情報が周囲の画素に滲み出てしまうからだ。エントロピー増大の法則により、ぼけた情報を元に戻すということはできないので、そもそもエッジの情報がぼけないように配慮する必要がある。そのために各種のアルゴリズムをひたすら試した末にECPBが生まれたのだが、その過程について述べよう。まず、要件を整理する。
- ロバストネス。画像の内容や種類にあまり左右されずに安定した結果が得られる。
- フィードフォワード。画像自身の統計処理に終止し、AI等のモデルや外部のデータを参照しない。
- ベクトル演算。画素単位のループを自前では書かず、NumPyかOpenCVのベクトル演算で高速に処理する。
- ガウシアンピラミッドを使う。カーネル関数では大きいぼけを制御できない。
- 空間認識をしない。外部のモデルに頼らずにはできないし、できたとしても精度と速度が期待できない。
エッジの情報が滲んでしまうのを防ぐには、いかなるアルゴリズムを使うにせよ、エッジに由来する画素は畳み込みに使わないという哲学を持つ必要がある。この発想に至った時に、答えは半分出たようなものだ。これをガウシアンピラミッドの文脈で考えてみよう。単純化のために、一次元配列で考える。まず、以下の画像が合ったとする、
[0.7, 0.7, 0.7, 0.7, 1.0, 1.0, 1.0, 1.0, 0.0, 0.0, 0.0, 0.0, 0.3, 0.3, 0.3, 0.3]
これをナイーブなガウシアンピラミッドのpyrDown関数にかける。2画素毎に、カーネルサイズ5のガウシアン関数を取る。つまり、中心画素の周囲（おそらく5*5）の範囲の画素を読んで正規分布の離散近似での重み付けをした算術平均を取る。実際には2次元なので、縦横2*2の4画素ごとだ。画像の恥から溢れた分は端の画素が続いているとみなす。したがって、以下の画素が情報が作られる。この時点で、エッジが溶けてしまっている。
[0.79, 0.98, 0.31, 0.06, 0.09, 0.28, 0.30]
これをナイーブなガウシアンピラミッドのpyrUp関数にかける。ガウシアン関数で画素の値を周囲に拡散しながら復元するので、以下のようにさらにぼけた画像になる。元画像では1.0と0.0の間に明確なエッジがあったのだが、たった1階層のガウシアンピラミッドをかけただけで、明確なエッジは消滅した。
[0.70, 0.70, 0.73, 0.81, 0.88, 0.92, 0.85, 0.62, 0.35, 0.16, 0.07, 0.11, 0.19, 0.25, 0.28, 0.30]
画像サイズに依らずに一定以上のぼけ効果を得るには、ピラミッド方式の畳み込みは絶対に必要だ。よって、ガウシアンピラミッドの基盤として、エッジが溶け出さないような工夫をしよう。また元の画像を見て欲しい。
[0.7, 0.7, 0.7, 0.7, 1.0, 1.0, 1.0, 1.0, 0.0, 0.0, 0.0, 0.0, 0.3, 0.3, 0.3, 0.3]
こいつからシャープネスマップを作ると、1.0と0.0の間に明確なエッジを検出するだろうから、以下のような画像が得られるだろう。
[0.0, 0.0, 0.0, 0.1, 0.1, 0.0, 0.0, 0.5, 1.0, 1.0, 0.5, 0.0, 0.1, 0.01, 0.0, 0.0]
ガウシアンピラミッドのpyrDown関数を改造して、シャープネスマップの閾値が0.8以上であればエッジとみなすことにする。エッジとみなした画素はpyrDownで適用されるガウシアン関数での平均算出に参加させないようにする。それを実現するためには、エッジとみなした画素を0とみなしたうえで、エッジ以外の画素数で割れば良い。同様にして、エッジとみなした画素はpyrUp関数で適用されるガウシアン関数での平均算出にも参加させないようにする。この二重のガードで、エッジを超えたデータの滲みがほとんど起きなくなる。よって、pyrDownとpyrUpを行っても、以下のようなエッジ情報が残るはずだ。
[0.5, 0.5, 0.5, 0.5, 0.8, 0.9, 0.9, 1.0, 0.0, 0.1, 0.3, 0.3, 0.5, 0.5, 0.35 0.5]
この機構で、ピラミッドの第1階層は保護できることがわかった。では、第2階層より深い部分はどうするか。それに対しては、二つの工夫がある、一つは、浅い階層で相当する2*2の領域にエッジが2個ある画素は、下層でもエッジとみなすというものだ。2個のエッジがあれば壁があるという意味なので、その壁を超えて情報が滲まないようにするのだ。しかし、常にその条件が該当するわけじゃないので、壁を作り損ねることもある。そうすると滲みが広がってしまうので、壁を補強する措置が必要になる。そのために、各階層での画像を入力として改めてシャープネスマップを作り、その閾値が0.8を超えているものは改めてエッジとみなすのだ。下層からエッジ情報の伝搬と現在の層で新たな補強によって、滲みを抑止するための壁は、階層が深くなっても安定的に維持されることになる。ラプラシアンとソーベルというカーネル関数を使っているシャープネスマップは高周波のエッジしか取得しないので、折り畳んで低周波情報として使うのは、上述の壁バイナリ化をしても1階層か2階層くらいで効果が切れてしまう。したがって、各階層でシャープネスマップを再計算するのが必要になる。
エッジとみなした画素を平均値の算出に参加させないと言ったが、そうするとエッジの画素に相当する情報欠損をどうやって補うかが問題となる。そこで、ラプラシアンピラミッドを併用する。ラプラシアンピラミッドとは、ガウシアンピラミッドの各階層での元画像に相当する情報を復元するためのデータ構造だ。つまり、ガウシアンピラミッドとラプラシアンピラミッドを持っていれば、任意の階層（1/2^N解像度の画像）のスケールにおいて元画像の情報を提示できる。よって、エッジ部分は、ぼかしていない元画像のデータに置き換えれば良い。ここで重要なのは、エッジが滲まないように加工したガウシアンピラミッドとは別に普通のカーネル関数で作ったガウシアンピラミッドを作っておき、ラプラシアンピラミッドの作成と復元は普通のガウシアンピラミッドを使って行うということだ。そうしないと一貫性が崩れてしまう。ChatGPTにこれを分かってもらうまでにかなり苦労した。
さらに補足すると、エッジとみなした画素を平均値の算出に参加させないという処理をベクトル計算のままでやるのにちょっとした工夫を要した。ガウシアンピラミッドを処理するcv2.pyrDownやcv2.pyrUpといった関数を呼ぶ限りにおいては、ガウシアンピラミッドの処理は高速に行える。それらの関数の中身はC言語で実装されているだろうし、GPUを使ってもっと高速化する余地もある。しかし、規定のガウシアン関数と異なる振る舞いをさせるのは難しい。各画素毎に処理を行うループを自分で書くことは可能だが、Pythonの層でそれをやると遅すぎて実用にならない。そこで考えたのが、ゼロ置換という手法だ。エッジとみなした画素の値を0に置き換えたバイナリマスクを適用してからpyDownやpyUpを行って得られる値を「分子」とみなし、そのバイナリマスクにpyDownやpyUpを行って得られる値を「分母」とみなし、分子を分母で割れば、非エッジ画素の平均が手に入る。pyDownやpyUpの中でどんなカーネル関数が使われていたとしても、分子の計算と分母の計算で同じ関数が使われている限り、一貫性のある結果が得られる。この処理を始めとする複雑なものであっても、ECPBで使われる計算は全て画像横断のベクトル計算で完結するようになっているので、十分に速い実用的な速度で動作する。
エッジを加味するpyrDownだけを有効にすると、以下の画像が得られる。これはほとんど効果がないことがわかる。pyrDownは元画像から処理を開始するが、元画像にシャープネスマップをかけて得られるエッジはスパースだ。エッジの壁が得られたとしてもそれは細い。pyrDownのガウシアン関数のカーネルが参照する領域（おそらく5*5）がエッジの壁を超えた先を参照するので、滲みが抑えきれない。エッジを加味するpyrDownは、やらないよりはマシだが、単体で効果があるとは言い難い。
エッジを加味するpyrUpだけを有効にすると、以下の画像が得られる。これは劇的な効果がある事が分かる。pyrUpは最小の縮小画像から始めるが、そこまでエッジ抽出を重ねていくと、ガウシアン関数のカーネル関数が参照する領域（おそらく3*3）を覆う程度の太いエッジの壁ができているので、エッジを超えた情報伝達を抑えることができる。
エッジを加味するpyrDownとpyrUpの両方を有効にすると、以下の画像が得られる。ほぼpyUpだけの結果と変わらないが、pyrDownが加わったおかげで若干ながらアーティファクトが減っていることがわかる。
Edge-ContaindなpyrUpとpyrDownを行うだけだと、1ピクセル単位でのシャープネスを守り切ることはできない。そこで、仕上げに、元画像に適用したシャープネスマップをマスクに使い、明らかにエッジであると判定できる画素は元画像からそのままコピーしてくる処理を入れた。実際にはシャープネススコアを重みにした平均を取るので、この処理自体がアーティファクトを生むことは一切ない。この仕上げを適用すると、結果は以下のようになる。これがECPBの完成形だ。
完成形と言いつつ、画面の奥が不自然にぼけすぎだろと言われそうだ。その通りで、小さい画像にレベル5とかでかけるとオーバーシュートになってしまう。この大きさの画像ならばレベル3くらいが丁度よいだろう。
ECPBのソースを載せようかと思ったが、やたら長くなるので割愛する。詳しくはitb_stack.pyのソースコードを読んでたもれ。blur_image_portraitという関数だ。
Stacked ECPB
ECPBで強度を上げると不自然に感じてしまうのは、トランジションゾーンがバイナリ化するからだ。言い換えると、ボケ始めてから大きなぼけになるまでの移行が急激になり、そこに加工っぽい雰囲気が出てしまう。エッジかどうかの閾値をバイナリ判断しているのだから仕方がない。ならば、閾値を変えて複数実行したものを平均合成すれば良さそうだという発想が出てくる。レベル5を指定すると、レベル1から5までのECPBの結果を作って、それを幾何平均で合成する。レベルが異なることで、ボケの縁が分かりづらくなり、幾何平均を使うとこでハロの抑制が効果的にできる。また、エッジ判定の閾値をレベル毎に変えて、ボケ始めの位置を少しずつずらすことで、トランジションゾーンの滑らかさを演出する。
レベル4のStacked ECPBの結果は、劇的に美しい。加工したと言われなければ気づかない人も多いだろう。
レベル4では物足りないので、レベル5にしてみよう。結構強いぼかしを掛けている割には、単体ECPBの時に見られたトランジションの急激さがなく、自然な範疇の結果と言えるだろう。
前景抽出
ここで、正直に言おう。ECPBもStacked ECPBも単体ではポートレートモード使い物にならない。それらはエッジを残して美しくボケさせるための手法であって、エッジじゃないものは容赦なく溶かしてしまうのだ。しかし、ぼけさせたくない前景（＝主要被写体）がエッジばかりで構成されているとは限らない。コントラストが低いものが主要被写体であった場合、エッジ検出に基づくアルゴリズムでは、それを背景と区別しないので、前景も背景と同様に溶かしてしまう。それではポートレートモードとして成立しない。
ならば、前景を抽出して、そこにマスクをかけてボカシ処理から保護しようと考える。そのための手法も多数提案されているが、今回は定番のGrabCutという手法を採用した。入力画像を「確定背景」「多分背景」「多分前景」「確定背景」の四種類に分けたマスクを読ませると、「多分背景」と「多分前景」の双方を背景と前景に自動分類して、前景だけを示すバイナリマスクに変換してくれるというものだ。そこで当然、こう突っ込みたくなる。前景と背景が分からないからGrabCutを使いたいのに、どうやって入力データを作るんだと。
シャープネスマップを加工して作るしかないのだが、できるだけ高周波のエッジを拾うようにする。ピントが合っている場所にしか低周波の模様は出現しないという光学上の特性を利用するのだ。よって、高周波のエッジを多く拾うラプラシアンの重みを強くして、シャープネスマップを作るようにする。そうすると、ノイズに起因するエッジも拾ってしまうので、それを削る処理が必要になる。そのためには、画像を400個のタイルに分割して、個々のタイルでのラプラシアンの平均を取り、それをフロアノイズ（ホワイトノイズ）とみなして、全画素のラプラシアンの値から引く。
こうして得られた、主に高周波成分を拾ったシャープネスマップを使って前景の領域を推定する。その前準備として、シャープネスマップを100個のタイルに分割して、その個々のタイルの平均を取っておく。そのタイルの並びを使って、二つの方法で前景を推定する。一つは「合体矩形」方式だ。縦横最大5*5タイルの組み合わせの中で、その平均シャープネスが最も高いものを選ぶ。ChatGPTは単に平均の最大化を目指すアルゴリズムを提案してきたので、「いやいや、それだと1*1しか選ばれないじゃん、NLPで使うbrevity penaltyとかlength normalizationみたいに頑張ろう」と言ったら、ゴニョゴニョ言った挙げ句にタイル数の平方根で割ったスコアを使う方式に落ち着いた。これだと、前景が狭い範囲に集中している写真では1*1を選ぶし、そうでなければそれなりに広い領域を拾うようになる。ただし、合体矩形方式だと前景が複数箇所に分散している対処できないので、二の矢も考えた。「トップタイル」方式だ。これは矩形にはこだわらずに、平均値が高いタイルを上位から拾い上げていって、その平均シャープネスの最大化を目指す。平方根で割るとそれなりに落ち着いた数になるのだが、自由度が高い分、0.5乗ではなく0.7乗にすることでバランスを取った。上に挙げた作例で前景を推定した結果が以下のものだ。二つの方式が相互補完的であることが分かるだろう。
なお、合体矩形方式では、タイルの選択方法で選択位置が大きく変わるリスクがあるので、タイルで選択した矩形（赤線）の面積を2倍に拡大した枠（青線）を作ってから、その中のシャープネスマップの重心を取って、拡大枠と相似形で面積半分の枠（緑線）を引き直している。これによって、タイルの数をそんなに増やさなくても妥当な矩形が切り出せるようになる。トップタイル方式は、採用するタイルが適当に散らばるので、わざわざその補正は行っていない。というか、うまく位置を直す方法が思いつかなかった。
おそらく前景を含むであろう領域を選択することができたなら、それをGrabCutに渡して、前景だと思われる領域を矩形でない自然な形で推測させる。GrabCutは色情報を使って、画像を2種類に分類するアルゴリズムだ。ここまでの前景検出処理で持っている情報をGrabCutに入力するには、以下の方法がある。
- 合体矩形：オープン ＝ 選択した合体矩形を多分前景、それ以外を多分背景にする ＝ 合体矩形の外側も積極的に拾う
- 合体矩形：クローズド ＝ 選択肢した合体矩形を多分前景、それ以外を確定背景にする ＝ 合体矩形の内側だけを拾う
- トップタイル：オープン ＝ 選択肢したタイルを多分前景、それ以外を多分背景にする ＝ タイルの外側も積極的に拾う
- トップタイル：クローズド ＝ 選択したタイルを多分前景、それ以外を確定背景にする ＝ タイルの内側だけを拾う
それぞれのクローズド方式は保守的すぎるので、選択領域の外でシャープネスが高い領域は多分背景に格上げしている。
それぞれの方式をGrabCutにかけた結果が以下のものだ。赤が合体矩形の情報で、オープンとクローズどの両方で反応すれば明るい赤、どちらかだけなら暗い赤になる。青がトップタイルの情報で、オープンとクローズどの両方で反応すれば明るい青、どちらかだけなら暗い青になる。紫になっているのは合体矩形とトップタイルの両方に反応しているということだ。
これらの例を見る限り、GrabCutはかなりうまく働いている。各画素が空間的につながってうまく分けられるようにするのがGrabCutの考え方だが、その際にRGBの色情報を使う。画素間のRGB空間でのユークリッド距離を使ってクラスタリングを行うということだ。ゆえに、前景と背景の輝度や色味がはっきりと異なる例ではうまく動くが、似ていると失敗する。当然、入力するデータの精度によって結果は大きく左右される。４つの入力方式のどれが良いとは言えないので、とりあえずは全部使って、その結果のRMS（2乗平均平方根）を取ることにした。つまり、1つの方式のどれかでも拾ったら √(1^2 / 4) = 0.5になる。
GrabCutに入力するタイルの選び方には様々な工夫を重ねた。その最大のものは、タイル毎のシャープネスの分布のジニ係数を使った適応的なパラメータ調整だ。シャープネスマップを作る際には、ガウシアンぼかしの強度と、ラプラシアンとソーベルの混ぜ方という二つのパラメータがある。ノイズや模様で高周波成分が多すぎる画像はガウシアンぼかしを強めにかけたほうが前景と背景の分類がしやすい。また、前景が高周波で背景が低周波領域であるとくっきり区別できるなら、ラプラシアンを多く使い、ソーベルは引くくらいが丁度よい。そうでないならラプラシアンとソーベルを同じ重みで使った方が良い。そういったパラメータの組み合わせの中から、最善のものを選ぶのだ。そして、最善かどうかを決めるのにジニ係数を使う。タイル毎のスコアが全く同じならローレンツ曲線は直線になるが、スコアの偏りが大きいならば湾曲が増す。その度合いがジニ係数だ。経済学の指標をこんなところで使うとは、そしてそれをChatGPTが提案してくるとは、ちと驚いた。面白いことに、ガウシアンぼかしの強度とラプラシアン/ソーベル割合の双方で、ジニ係数は中庸に収束する。ぼかしが弱すぎても強すぎてもタイル同士が区別できないし、高周波を取りすぎても低周波を取りすぎてもタイル同士が区別できない。その塩梅を自動設定できることは意義が大きい。
もう一つ工夫がある。トップタイル方式では前景っぽい領域を全て選んでほしいが、スコアをタイル数の平方根で割る方式だと、選択されたタイルの数のみを見てバランスを決めることになる。つまり、構図を無視して面積だけでペナルティをかけている。しかし、写真の構図によっては、ごく一部のタイルのみを選んだほうが良いかもしれないし、半分以上のタイルを選んだ方がよいかもしれない。そこで、新たなスコア要素を導入した。選択肢したタイル毎の分散と選択肢していないタイル毎の分散の平均が画像全体の分散よりも下がっている方が良いとするものだ。それを導入したところ、かなり良い塩梅にタイルが選択できるようになった。
最終的には、GrabCutで作ったマスクを使って、Stacked ECPBの結果のボケ画像と、元画像をブレンドする。マスクの値を重みに使って算術平均を取るのだ。マスクの値が高い（前景っぽい）画素ほど多くの情報を元画像から使い、マスクの値が低い（背景っぽい）画素ほど多くの情報をボケ画像から引っ張る。結果として、前景はシャープになり、背景はボケて、どっちかよくわからない領域はうっすらボケる。4つの方式を混ぜたことにより、確信度が高い部分はしっかりしたエッジを残せるし、そうでもないトランジションゾーンがなだらかボケるのが美点だ。
GrabCutが賢いならECPBは要らないのでは？と思う人もいるかもしれないので、GrabCutと普通のガウシアンピラミッドの組み合わせの例も乗せておく。これで証明されるように、ECPBがエッジの滲みを抑制しないと、後でエッジを回復させても、周囲にアーティファクトが出てしまう。よって、GrabCutとECPBの組み合わせが必要だと言える。
ということで、ポートレートモードがここに完成し、itb_stack.pyから簡単に呼べるようになっている。--portrait=autoとすると、画像サイズに応じて適切なレベルのStacked ECPBを掛けてくれる。レベルを明示的に設定したい場合には--portrait=4とかやれば良い。--portrait=4:decay=0.1:edge_threshold=0.9 みたいな細かい設定もできるようになっている。
$ itb_stack.py input.jpg -o output.jpg --portrait=auto
ポートレートモードの前景抽出は自動的に行われ、前景だと思われる部分がぼかさないようにしてくれる。デフォルトでは画面内で最もシャープネスが高い領域を選ぶが、その位置を明示的に指定することもできる。例えば画面上のx軸（横方向左から）55%、y軸（縦方向上から）40%の位置に前景の中心があるなら、以下のようにする。重みのデフォルトは0.1だが、1とか2にすると、指定した位置をほぼ確実に選んでくれるようになる。
$ itb_stack.py input.jpg -o output.jpg --portrait=auto:attractor=0.55,0.4:attractor_weight=1
今後の課題。
Stacked ECPBと混合GrabCutマスクの併用により、各社のスマホのポートレートモードと同等またはそれ以上の精度で動作するポートレートモードができた。トランジションゾーンが自然なので、どんな被写体に使っても不自然さが少ない。AIとか使って前景と背景の区別していないのは欠点とも言えるが、空間認識の精度に依存しないで安定した結果が出るのは利点とも言える。
ポートレートモードの基本的な問題として、ガウシアンぼかしによる縁がじわっと溶けるボケの形と、均一濃度の錯乱円を重ねることによる光学的なボケの形が全く違うというのがある。つまり、ポートレートモードの強度をいくら強くしても、いくらエッジの認識が賢くなろうとも、光学的なボケと同じ描写にはできないのだ。錯乱円に似たようなボケを作るカーネル関数を実装すればソフトウェア的に実現できる可能性はあるが、OpenCVで高速に実現する方法についてまだ見つけられていない。
GrabCutによる前景選択には、まだまだ改良の余地がある。今のところ、最初のタイルを選択する際にラプラシアンによる高周波成分を積極的に拾うためにパラメータを動的に調整するアイデアと、そうして得られたエッジの情報を擬似色空間としてGrabCutに渡すアイデアを考えている。それについては追ってまとめる。
まとめ
エッジを残してぼけだけを強調するポートレートモードの効果と実装の詳細について説明した。単にぼけさせるだけでもかなりの数の様々なアルゴリズムを組み合わせていることが分かってもらえたと思う。開発途中で何度も課題にぶちあたってくじけそうになったが、なんとか実用になる機能として完成させることができた。画像処理の基本であるピラミッド処理とカーネル関数をふんだんに使っていし、画像処理やOpenCVの技術習得に丁度よい課題だった。タイル分割とGrabCutによる前景推定のアルゴリズムを考えるのも頭の体操になった。この開発作業を通じて私はかなり賢くなったし、その対話相手であるChatGPTも私の好みをどんどん覚えてくれて、いい感じに開発が進められるようになってきた。今後もこの修行を続けていこう。
