シリーズ６．マクロ言語を使った画像処理の応用編～ノイズ軽減①　空間フィルタ処理～ - IMACEL Academy -人工知能・画像解析の技術応用に向けて-| エルピクセル株式会社

【空間フィルタ処理】

(2n+1)×(2n+1)のカーネルと呼ばれる小さな画像を定義し、処理を施したい画像の各ピクセルに対応する輝度値と掛け合わせの計算を行い、計算結果を中心の新しい輝度値としてもとの輝度値と置き換える処理のことです。この計算処理を畳み込みといいます。

〇カーネル　※四角の中の数字はフィルタの種類によって変わります。

空間フィルタ処理には線形フィルタと非線形フィルタの2種類があります。
線形フィルタとは畳み込み演算を行う処理のことで、Mean（平均化フィルタ）とGaussian Blur（ガウシアンフィルタ）があります。

【Mean】

例　下記のような画像があり3x3のカーネルと赤枠で囲んだ領域の輝度値から計算された平均値で中央の輝度値128を置き換えることを考えます。

この場合の計算は以下のようになります。

上記のようになり、3x3の□の中の中央の輝度値128が111.777778に置き換えられます。
この計算処理を3x3の枠をずらしながらすべてのピクセルについて行い、平均値で輝度値を置き換えていく処理がMeanです。
画像の一番外側のピクセルについては3x3の枠を当てはめたときに計算に必要なピクセルが足りなくなってしまいますが、これは外側のピクセルを複製することによって解決することがあります。（詳しくは次回の記事で説明します。）

【Gaussian Blur】

ガウシアンフィルタは正規分布に従って中央のピクセルから近いピクセルの輝度値を優先的に用いるように重みづけをして画像の平滑化を行うフィルタです。
下記のような3x3のカーネルの場合、中央のピクセルは3、隣り合うピクセルは2、最も距離の遠い四隅のピクセルは1という重みをつけて計算します。

この計算処理を3x3の枠をずらしながらすべてのピクセルについて上記の計算を行い、輝度値を置き換えていく処理がガウシアンフィルタです。

※MeanとGaussian Blurは輝度値が小数点をとることがあるので処理を施したい画像を8bitから32bitに変換する必要があります。

【Median】

空間フィルタ処理のうち非線形フィルタに分類され、中央値フィルタとも呼ばれます。MeanやGaussian Blurのような畳み込みを行わないフィルタ処理で、指定した大きさの領域内のピクセルに対応する輝度値を大きさの順に並べたときの中央値を中心の輝度値と置き換える処理を行います。

赤枠内の輝度値を小さいほうから順に並べると、以下のようになります。
63, 72, 90, 91, 100, 112, 128, 150, 200
この場合の中央値は100なので赤枠内の中心の輝度値128は100に置き換えられます。

【ImageJにおけるフィルタ処理】

ImageJではカーネルを自分で指定することができます。
画像を開いてProcess→Filters→Convolveをクリックします。

3x3の平均値フィルタに使用するカーネルを設計するときは下記のように1を入力します。

※この画面での入力する数値と大きさを変えることによって最適なフィルタを探します。

【ImageJ　マクロ言語で記述してみましょう】

ここではMeanでフィルタ処理を行うスクリプトをマクロで記述することを考えてみます。画像の原点（左上）から右下に向かって3x3のカーネルを使って9ピクセルずつ輝度値を読み込みその平均値で中央の輝度値を置き換えればよいことになります。
　図のようにある画像に3x3のカーネルを設計しフィルタをかけるとします。本来は画像の端を複製するなどの対策をとる必要がありますが、ここでは上下左右の端1ピクセルについてはフィルタをかけないとします。最初は(0,0)、(0,1)、(0,2)、(1,0)、(1,1)、(1,2)、(2,0)、(2,1)、(2,2)の9個の座標に対応する輝度値を読み込みその平均値を計算します。そしてその平均値が図の青い四角で示した中央の輝度値と置き換えられます。この操作をx座標とy座標を1ずつ増やしていくことにより右下に向かってfor文を用いて繰り返します。図に示した位置から計算が始まるようにx=1、y=1でfor文を開始します。

最終的には下図にように青い領域がフィルタ処理されることになります。

＜スクリプト＞

width=getWidth();
height=getHeight();
N=2*1+1;//カーネルの大きさN=2n+1, nは任意の数値に変更可。この場合は3x3のカーネル
intensity=0;

run("32-bit");
for(y=1;y<height;y++){
        for(x=1; x<width;x++){ 
                        
intensity=(getPixel(x-1,y-1)
+getPixel(x,y-1)
+getPixel(x+1,y)
+getPixel(x-1,y)
+getPixel(x,y)
+getPixel(x+1,y)
+getPixel(x+1,y+1)
+getPixel(x, y+1)
+getPixel(x+1, y+1))
/(N*N);                
                        setPixel (x,y, intensity); 
                        
                }
        }

ノイズ軽減

左側が上記スクリプトを実行した画像です。右側のフィルタ処理をする前の原画像と比べるとフィルタがかかっていることがわかると思います。

ここで上下左右の端1ピクセルの輝度値が本当に変わっていないかピクセルの四角が見えるまで拡大して確認します。
赤枠の中はフィルタを書けたので原画像の輝度値とは変化しています。一方、赤枠の外にある端の輝度値は原画像と同じです。

したがって今回の処理では周囲の端1ピクセル分はフィルタ処理がかかっていないことが確認できました。

　ではここでもしx=0、y=0として上記のスクリプトを実行すると結果がどうなるのか考えてみましょう。x=0、y=0でfor文を書きなしてもスクリプトは動きます。

　しかし、(-1,-1)、（-1,0）、（0, -1）など負の値をとる座標が出てきます。試しに(-1,-1)に対応する輝度値をgetPixel(-1,-1)で取得し出力すると0と出てきます。
これはImageJが画像にはない負の値をとる座標に対応する輝度値を0として処理していることを意味します。つまり下図のように上下左右の周囲1ピクセルに対応する輝度値をすべて0として処理しているのでこれらのピクセルを含む3x3の9個のピクセルの平均輝度値が他の画像内のピクセルのみで計算された輝度値と比較して低くなります。

実際にx=0、y=0でのフィルタ処理の結果を見てみると、以下のようになります。

端を拡大するとさきほどのフィルタ処理後の画像に比べ端1ピクセルが暗くなっており、明らかに低い輝度値で周囲のピクセルが置き換えられていることが確かめられました。
これより、画像の外側1ピクセル分を輝度値0としてフィルタ処理をしてしまうのは不正確であると考えられます。
次回の記事で端をどのように処理すればよいかについて説明していきたいと思います。