画像処理.js -エッジ検出 1
こんにちは，なめこです．
久しぶりの技術記事です．
大学での実験課題や趣味プログラミング，あとお仕事で画像処理に関する知識やら画像処理を行うjsプログラムやらが増えてきたので，せっかくなのでまとめていきたいと思います．
毎週(もしかしたら隔週)ちょっとずつ「画像処理.js」のタイトルで記事にしていくつもりです．もちろんjsで実装します．他の言語を使う予定なんてないです．処理速度?知らない子ですねぇ．
また基本的に画像処理のライブラリにある関数は使わずに泥臭くアルゴリズムそのものを実装をしていく予定です (画像処理ライブラリは入出力を簡単に行うためにあると思ってる)．
ライブラリ
画像処理用，というより画像の入出力用としてopencv4nodejsを使用します．この子は普通に優秀です．人のこと言えないですけどなんでjsで画像処理しようと思ったんですかね．不思議．
その他数学的な計算用にmath.jsを使用します．
opencv4nodejsでの基本的な操作をまとめておきます．
以降，画像の入出力やMat <-> Array変換の処理の記述は省略します(変数名の末尾がimgとなっているものは画像の配列)．
画像入出力
\\ カラー画像読み込み const colorImage = cv.imread(image_path, cv.CV_8UC3); \\ グレースケールで読み込み const grayImage = cv.imread(image_path, cv.CV_8UC3); \\ 画像表示 cv.imshow('image', colorImage); \\ キー入力があるまで待機 cv.waitKey(); \\ 画像保存 cv.imwrite(outImage_path, colorImage);
Mat <-> Arrayの相互変換
\\ Mat -> Array const img = image.getDataAsArray(); \\ Array -> Mat const image = new cv.Mat(img, cv.CV_8UC3);
画像に対する基本処理
畳み込み演算
画像処理系の基本である畳み込み演算の実装をします．
せっかくなので使い回しをしやすいようにモジュール化してみました．
--- 2月16日更新 ---
画像の領域外処理を0埋めとするか，画像の拡張を行うかの設定をできるようにしました．
エッジ検出を行う際は画像を拡張させて畳み込むことで，縁取りのようなものが出ないようになります．
const ZEROPADDING = 0; const EXPAND = 1; module.exports.ZEROPADDING = ZEROPADDING; module.exports.EXPAND = EXPAND; module.exports.conv = (img, imgSizes, filter, filterSize, {func=null, mode=ZEROPADDING} = {func:null, mode:ZEROPADDING}) => { const cols = imgSizes[1]; const rows = imgSizes[0]; const filterCols = filterSize[0]; const filterRows = filterSize[0]; const filterColsHalf = (filterCols-1)/2; const filterRowsHalf = (filterRows-1)/2; if(func == null){ if(mode==ZEROPADDING){ return Array.from(new Array(rows)).map((_, j) => Array.from(new Array(cols)).map((v, i) => filter.reduce((pre, cur, l) => { let row = j+l-filterRowsHalf; if(row < 0 || rows <= row)return pre; return pre+cur.reduce((pre, cur, k) => { let col = i+k-filterColsHalf; if(col < 0 || cols <= col)return pre; return pre+img[row][col]*cur; }, 0); }, 0) )); }else if(mode==EXPAND){ return Array.from(new Array(rows)).map((_, j) => Array.from(new Array(cols)).map((v, i) => filter.reduce((pre, cur, l) => { let row = j+l-filterRowsHalf; if(row < 0)row = 0; if(rows <= row)row = rows-1; return pre+cur.reduce((pre, cur, k) => { let col = i+k-filterColsHalf; if(col < 0)col = 0; if(cols <= col)col = cols-1; return pre+img[row][col]*cur; }, 0); }, 0) )); } }else{ if(mode==ZEROPADDING){ return Array.from(new Array(rows)).map((_, j) => Array.from(new Array(cols)).map((v, i) => func(filter.reduce((pre, cur, l) => { let row = j+l-filterRowsHalf; if(row < 0 || rows <= row)return pre; return pre+cur.reduce((pre, cur, k) => { let col = i+k-filterColsHalf; if(col < 0 || cols <= col)return pre; return pre+img[row][col]*cur; }, 0); }, 0)) )); }else if(mode==EXPAND){ return Array.from(new Array(rows)).map((_, j) => Array.from(new Array(cols)).map((v, i) => func(filter.reduce((pre, cur, l) => { let row = j+l-filterRowsHalf; if(row < 0)row = 0; if(rows <= row)row = rows-1; return pre+cur.reduce((pre, cur, k) => { let col = i+k-filterColsHalf; if(col < 0)col = 0; if(cols <= col)col = cols-1; return pre+img[row][col]*cur; }, 0); }, 0)) )); } } };
画像(2配列)とフィルタ，それぞれのサイズを引数として畳み込み後の画像(2次元配列)を返す関数として実装しました．
funcとして関数を与えると畳み込み後の各画素に対して任意の処理を追加することができる仕様です．
ちなみに領域外は0として処理してあります．
mapとかreduceとかをフル活用してるので可読性は低いですね．
ただまあ畳み込み後の配列を確保と計算を同時に実装した結果こうなりました．jsの配列に対するこれらの処理が好きなので変えるつもりはないです．
各画素への処理
毎回img.map.....って感じで書くのも面倒なのでこれもモジュール化します．
面倒と言ってもOne-Linerで済むのですが変なところでtypoするのも嫌ですし，何度も使うので簡単にできるものは簡単にしましょう．エンジニアは怠惰であるべき．
module.exports = (img, func) => img.map((v, j) => v.map((v, i) => func(v, i, j, img)));
はい，特段特殊なことはしてないです．全画素に対して行いたい処理を関数として渡せばいいだけですね．
必要に応じて添字とか画像の配列そのものを取得できます．
エッジ検出
準備が整ったので本題に入ります． 今回は一般的なエッジ検出ということでSobelフィルタとCanny法を実装していきます．
Sobelフィルタ
Sobelフィルタは基本的なエッジ検出で，その他の手法を用いる際の前段階として使われることが多いです．
処理の内容としては単純で，x軸方向，y軸方向それぞれで(離散的な)微分をとるだけです．
そうすると各軸方向の勾配の大きさを持つベクトルが得られるので，そのベクトルの大きさを計算すれば完了です．
各軸方向のエッジ(正確にはエッジに直交する)を得るためのフィルタは以下の通りです．
\begin{align} x軸方向 = \begin{bmatrix} -1 & 0 & 1 \\ -2 & 0 & 2 \\ -1 & 0 & 1 \end{bmatrix}, y軸方向 = \begin{bmatrix} -1 & -2 & -1 \\ 0 & 0 & 0 \\ 1 & 2 & 1 \end{bmatrix} \end{align}
このフィルタの適用後， \begin{align} g_{ij} = \sqrt{gx_{ij}^2 + gy_{ij}^2} \end{align}
とすれば各画素における濃度勾配が得られます．
実装
--- 2月16日更新 ---
畳み込み処理のモジュールの仕様変更に対応しました．
// Sobelフィルタ作成 const sobelX = [[-1, 0, 1], [-2, 0, 2], [-1, 0, 1]]; const sobelY = [[-1, -2, -1], [0, 0, 0], [1, 2, 1]]; const filterSizes = [sobelX.length, sobelX[0].length]; // Sobelフィルタ適用 const edgeXImg = conv.conv(img, image.sizes, sobelX, filterSizes, {func:val=>(math.abs(val)>255?255:math.abs(val)), mode:conv.EXPAND}); const edgeYImg = conv.conv(img, image.sizes, sobelY, filterSizes, {func:val=>(math.abs(val)>255?255:math.abs(val)), mode:conv.EXPAND}); const edgeImg = map(edgeXImg, (v, i, j) => { const r = math.sqrt(v*v+edgeYImg[j][i]*edgeYImg[j][i]); return r>255?255:r; });
絶対値とったり最終的にはuint8の形にしたいという理由で無理やり0~255の範囲に納めてます．
その他特に難しい処理はないので説明は省きます．
実行結果
結構綺麗にエッジが取れたんじゃないかと思います．
Canny法
今回初めてちゃんと調べたのですが，案外処理が多いんだなと思いました．
手順としては以下の通りです．
- Gaussianフィルタを適用
- Sobelフィルタを適用
- Non maximum supperession処理
- Hysteresis threshold処理
詳細を説明するのは面倒なので簡単にだけ説明すると，
- ノイズ除去
- 輪郭抽出(濃度勾配計算とついでに勾配方向の計算)
- 非極大抑制(輪郭の細線化)
- 濃度勾配が小さく，勾配が大きい画素に隣接していないエッジの除去
と言った感じです．
詳細を知りたい方は参考文献の方をご参照ください．
内容としては割と単純だったのでさくっと実装しましょう．
実装
--- 2月16日更新 ---
畳み込みモジュールの仕様変更への対応と一部処理の間違いの修正を行いました．
// 定数定義 const imgSize = image.sizes; let gaussianFilterSize = 5; let thresholdLow = 0.25; let thresholdHigh = 0.35; // Gaussianフィルタ作成 const gaussianFilterSizeHalf = (gaussianFilterSize-1)/2; const sigma = gaussianFilterSizeHalf/2; const gaussianFilter = Array.from(new Array(gaussianFilterSize)).map((_, j) => Array.from(new Array(gaussianFilterSize)).map((_, i) => { return math.exp(-(j*j+i*i)/(2/sigma/sigma))/(2*math.PI*sigma*sigma); })); const gaussianImg = conv.conv(img, imgSize, gaussianFilter, [gaussianFilterSize, gaussianFilterSize], {mode:conv.EXPAND}); // Sobelフィルタ作成 const sobelX = [[-1, 0, 1], [-2, 0, 2], [-1, 0, 1]]; const sobelY = [[-1, -2, -1], [0, 0, 0], [1, 2, 1]]; const filterSizes = [sobelX.length, sobelX[0].length]; // Sobelフィルタ適用 const sobelXImg = conv.conv(gaussianImg, imgSize, sobelX, filterSizes, {mode:conv.EXPAND}); const sobelYImg = conv.conv(gaussianImg, imgSize, sobelY, filterSizes, {mode:conv.EXPAND}); const sobelImg = map(sobelXImg, (v, i, j) => { const r = math.sqrt(v*v+sobelYImg[j][i]*sobelYImg[j][i]); return (r>255?255:r)/255; }); // エッジに直行する角度の計算 const thetaImg = map(sobelXImg, (v, i, j) => math.atan2(sobelYImg[j][i], v)); // non maximum supperession処理 const nmsImg = map(sobelImg, (v, i, j, img) => { const theta = thetaImg[j][i]; let di = 0; let dj = 0; const t = math.abs(theta); if(t < math.PI*3/8 || t > math.PI*5/8){ di = 1; } if(t > math.PI/8 && t < math.PI*7/8){ dj = (theta<0?theta+math.PI:theta)<(math.PI/2)?1:-1; } if((img[j+dj]?(img[j+dj][i+di] > v):true) || (img[j-dj]?(img[j-dj][i-di] > v):true)){ return 0; } return v; }); // hysteresis threshold処理 const cannyImg = map(nmsImg, (v, i, j, img) => { if(v < thresholdLow)return 0; if(v > thresholdHigh)return v; const theta = thetaImg[j][i]; let di = 0; let dj = 0; const t = math.abs(theta); if(t < math.PI*3/8 || t > math.PI*5/8){ dj = 1; } if(t > math.PI/8 && t < math.PI*7/8){ di = (theta<0?theta+math.PI:theta)<(math.PI/2)?-1:1; } if((img[j+dj]?(img[j+dj][i+di] > thresholdHigh):false) || (img[j-dj]?(img[j-dj][i-di] > thresholdHigh):false)){ return v; } return 0; });
はい，まあそこそこ長くなりましたね．
Gaussianフィルタは自作するところからしてます．ちなみにCannyを含む多くの場合，カーネルサイズは5とするようです．
Sobelフィルタの濃度勾配計算の際にさりげなく0~1の範囲に縮小してあったりします．
角度計算にはarctanが使われています．atan2は象限も考慮してくれる優れものです．
Non maximum supperessionとHysterisis thresholdは濃度勾配の角度に応じて参照する画素が変わるので，ちょっと処理が面倒でした．
x軸方向，y軸方向，2つの斜め方向の計4つの方向に分けて参照する画素が変わるような設計です．
本当は角度の値に応じて一つ一つ処理を書いていくほうが綺麗なコードになるのですが，きっともっと単純化できるはずだと中途半端に考えた結果妙なコードが生まれました．多分一応処理としては合ってます．違ってたらごめんなさい実際ちょっと間違えてました．ごめんなさい．(2月16日追記)
因みにthresholdLowとthresholdHighは天から降ってきた値です．その時々によって最適な閾値は変わるのでしかたないです．
また領域外処理を書くのが面倒でちょっとよろしくない書き方をしてあります．
jsの配列とundefinedには以下のような特徴があります．
const a = [1,2,3]; console.log(a[5]); // <- undefined console.log(a[-1]); // <- undefined console.log(undefined < 1); // <- false console.log(undefined ? 0 : 1); // <- 1
要するに領域外にアクセスしても即時エラーにはならず，またundefinedは評価の際は単体・比較に関わらずfalseとなる，ということです．
これを利用して「undefinedが出たら領域外ってことにしてね」という処理にしてあります．
いやほんとこれは良くないので良い子は真似しないでね．
実行結果
それぞれの処理の結果も入れておきます．
画像が1024x1024でそこそこ大きいせいかちょっと微妙な感じがありますね．パラメータチューニング次第で改善されるかもしれませんが，とりあえず実装はできてると思うのでこれでよしとしましょう．
--- 2月16日更新 ---
Canny法の処理の修正を行ったので結果が改善されました．
Sobelフィルタでの処理時には陰影にあたる部分も輪郭線として出ていましたが，Canny法ではそれらがなくなっていることがわかります．
その分境界がもともと曖昧だった輪郭は抽出されていません．
まとめ
今回は画像処理の基本的な処理を行うための関数の実装と一般的なエッジ検出を実装してみました．
エッジ検出をするときはSobelさんにはよくお世話になるので，Sobelフィルタくらいはすぐに実装できるようになっておくといいかもしれません．まあ画像処理系のライブラリを使えば普通Sobelフィルタの関数があるので自前で実装する機会はほとんどないと思いますが...
今回の実装も含め，jsによる画像処理の実装は以下のGithubにまとめてあります．
次回は最近見つけた「色情報を使ったエッジ検出」を実装してみたいと思います．
来週記事にできるといいなぁ．
参考文献
edge detection (Sobel filter/Prewitt filter)
