バイラテラルフィルタ

今回の記事では，左の写真のようなノン・フォトリアリスティックな画像処理をしたいと思います．

「バイラテラルフィルタとは何か」ということに関しては，参考に載せてありますが，ガウシアンフィルタとバイラテラルフィルタというページに，詳しく解説されています．ガウシアンフィルタを発展させたフィルタで，単純ですが効果的な画像処理フィルタです．

実装する数式

一見難しい数式ですが，いきなり行きます．
とりあえず，フィルタの大きさは$5\times5$にします．

BF(i,j) = \frac{\sum_{n=-2}^2 \sum_{m=-2}^2 exp\bigl(-\alpha (n^2+m^2)\bigr)exp\Bigl(-\beta\bigl(f(i,j)-f(i+n, j+m)\bigr)^2\Bigr)f(i+n,j+m)}{\sum_{n=-2}^2 \sum_{m=-2}^2 exp\bigl(-\alpha (n^2+m^2)\bigr)exp\bigl(-\beta(f(i,j)-f(i+n, j+m))^2\bigr)}
\\
= \frac{\sum_{n=-2}^2 \sum_{m=-2}^2 exp\bigl(-\alpha \text{フィルタの中心からの距離}\bigr)exp\Bigl(-\beta\text{色の差分}\Bigr)f(i+n,j+m)}{\sum_{n=-2}^2 \sum_{m=-2}^2 exp\bigl(-\alpha \text{フィルタの中心からの距離}\bigr)exp\bigl(-\beta\text{色の差分}\bigr)}
\\
= \frac{\bigl(\text{セルの重み}\times\text{セルの画素値}\bigr)の和}{\text{セルの重みの和}}

$BF(i, j)$は画像:$f$の$i,j$にバイラテラルフィルタをかけた時の出力画素値とします．

ちょっと長めの数式で，怖気付きそうですが，少し耐えましょう！よくよくみてみると，分母と分子がほぼ同じなことがわかります．
思ったより簡単に解明するのが早そうです．
ひとまず，この数式の要素を全部解明していきたいと思います．

変数n,m

$\sum$で使われている変数です．
nはx軸方向，mはy軸方向のフィルタ内の中心からの変位です．

パラメータαとβ

$\alpha$: 値が大きいほど，フィルタの中心から離れたピクセルの影響を受けない．
$\beta$: 値が大きいほど，対象のピクセルとの色の差が大きいピクセルの影響を受けない．

それぞれ，値が大きいほど影響を受けないのは，$exp$の指数に$\text{マイナス}$が含まれているからです．

空間フィルタとは

さて，「フィルタの中心からの距離」だとか「対象のピクセル」というのはどういうものでしょうか．これを理解するためにも，まずは空間フィルタというものを理解しましょう．空間フィルタとは，「対象となるピクセルの周りのピクセルも用いて処理を行うフィルタ」のことです．

空間フィルタはこのような形をしていて，ど真ん中に対象ピクセルがやってきます．

画像処理とは，加工したい画像の「各ピクセルに処理を施す」ことを意味します．なので通常画像処理を行う際は，その加工したい画像の画素数分だけ，同じような処理が実行されています．それゆえ，多重のfor文や$\sum$がよく出てきて，初学者を混乱させてしまいます．

この記事内では「今処理をしたい各ピクセル」を「対象ピクセル」と呼んでいます．
$f$を加工したい画像（入力画像）として，$f(i, j)$は画像上の$(i, j)$のピクセルの情報を表します．
- 画像それ自体は数式による関数ではありませんが，処理フィルタを数式で表すためにも．
  $(i,j)$という座標を渡せば，その座標の画素データが帰ってくる関数$f$と解釈しています．
- コード上では，color pix_ij = pixels[j*width+i];としてあげれば良いでしょう．

空間とは領域のことであり，その領域の中心に，対象ピクセルが来るように空間フィルタをセットします．基本的にはフィルタの中心に対象画素がくるようにフィルタを配置し，それを全ピクセルに対してフィルタをかけます．

フィルタの各セルでは

\text{各セルの重み}\times\text{各セルの画素値}

を計算しています．これが分子の$\sum$の中身です．

Σとfor文

数式の$\sum$は，どんな意図で使われているか考えなくても，一応プログラムに落とし込むことができます．
$\sum_{i=0}^9$はfor(int i=0; i<10; i++) { ~ }になりますし，
$\sum_{i=a}-b$はfor(int i=a; i<=b; i++) { = }としてあげればプログラムに置き換えることができます．
そして，この{ ~~ }に$\sum$の中身を書いてあげれば良いのです．

式を読みにかかる

今回のバイラテラルフィルタでは


\text{各ピクセルの重み} = exp\bigl(-\alpha (n^2+m^2)\bigr)exp\Bigl(-\beta\bigl(f(i,j)-f(i+n, j+m)\bigr)^2\Bigr)

であり，式の大部分が，このピクセルの重みに使われているのがわかります．(というよりも，フィルタは重みの配列なので当たり前ですが...)



BF(i,j) = \frac{\sum_{n=-2}^2 \sum_{m=-2}^2 \text{重み}\times f(i+n,j+m)}{\sum_{n=-2}^2 \sum_{m=-2}^2 \text{重み}}
\\

そしてこの，ピクセルの重みはexpでできています．

そして，$exp\bigl(-\alpha (n^2+m^2)\bigr)$ と $exp\Bigl(-\beta\bigl(f(i,j)-f(i+n, j+m)\bigr)^2\Bigr)$のどちらについても，

$-\alpha (n^2+m^2) \le 0$
$-\beta\bigl(f(i,j)-f(i+n, j+m)\bigr)^2 \le 0$

となるので，$(n^2+m^2)$または$\bigl(f(i,j)-f(i+n, j+m)\bigr)^2$が大きくなればなるほど，そのセルの重みは小さくなることがわかります.

$n^2+m^2$はフィルタ内の中心からの距離を表しています．
$\bigl(f(i,j)-f(i+n, j+m)\bigr)^2$はフィルタ内の中心の画素値との画素差を示しています．

ガウス関数

最初に，バイラテラルフィルタはガウシアンフィルタの応用と説明しました．それゆえ，ガウス関数を使っているのですが，係数などを取り払ってしまい，ただの$exp$関数を使用しています．

実装する

このバイラテラルフィルタを，画像のRGB成分それぞれに対して適用します．

シンタックスハイライトのため，拡張子はpde→javaに変更

bilateral.java

// parameters
final float A = 0.001;
final float B = 0.01;
final int W = 10;

// resource
PImage koara;
color[] k_pixs;
color[] koara_bilateral;
int w, h;

void settings() {
  koara = loadImage("koara.jpg");
  w = koara.width;
  h = koara.height;
  size(w*2, h);
  k_pixs = koara.pixels;
  koara_bilateral = new color[w*h];
  for (int i=0; i<h*w; i++) { koara_bilateral[i]=color(100); }
}

void setup() {
}

void draw() {
  for (int y=W; y<h-W; y++) {
    for (int x=W; x<w-W; x++) {
      color c = bilateral(x, y);
      koara_bilateral[y*w+x] = c;
      set(x, y, koara_bilateral[y*w+x]);
    }
  }
  image(koara, w, 0);

  // 多重にかけるために，繰越処理をする．
  k_pixs=koara_bilateral.clone();
}

color bilateral(int i, int j) {
  float sum_r = 0.0;
  float W_3_r = 0.0;
  float sum_g = 0.0;
  float W_3_g = 0.0;
  float sum_b = 0.0;
  float W_3_b = 0.0;
  for (int n=-W; n<=W; n++) {
    for (int m=-W; m<=W; m++) {
      color c_ij = k_pixs[i+j*w];
      color c_nm = k_pixs[(i+n)+(j+m)*w];
      sum_r += gauss(n*n+m*m, pow(red(c_ij)-red(c_nm), 2), A, B)*red(c_nm);
      W_3_r += gauss(n*n+m*m, pow(red(c_ij)-red(c_nm), 2), A, B);
      sum_g += gauss(n*n+m*m, pow(green(c_ij)-green(c_nm), 2), A, B)*green(c_nm);
      W_3_g += gauss(n*n+m*m, pow(green(c_ij)-green(c_nm), 2), A,B);
      sum_b += gauss(n*n+m*m, pow(blue(c_ij)-blue(c_nm), 2), A,B)*blue(c_nm);
      W_3_b += gauss(n*n+m*m, pow(blue(c_ij)-blue(c_nm), 2), A, B);
    }
  }
  return color(sum_r/W_3_r, sum_g/W_3_g, sum_b/W_3_b);
}

float gauss(float x, float y, float a, float b) {
  return exp(-a*x-b*y);
}

バイラテラルフィルタを複数回かけてみた時の遷移画像です.
左上は1回で右下が6回です.
だいぶイラスト調な画像になっているのがわかります．

参考

ガウシアンフィルタとバイラテラルフィルタ