バックナンバーはこちら。
https://www.simulationroom999.com/blog/compare-matlabpythonscilabjulia3-backnumber/
はじめに
非極大値抑制にHysteresis Thresholdを加えた、Canny法を実施して、2値化を行う。
今回はScilabで実施する。
登場人物
博識フクロウのフクさん
イラストACにて公開の「kino_k」さんのイラストを使用しています。
https://www.ac-illust.com/main/profile.php?id=iKciwKA9&area=1
エンジニア歴8年の太郎くん
イラストACにて公開の「しのみ」さんのイラストを使用しています。
https://www.ac-illust.com/main/profile.php?id=uCKphAW2&area=1
【再掲】Hysteresis Thresholdを実施するための手順
太郎くん
まずは手順を再掲。
- Sobelフィルタ等の微分フィルタで以下を推定
- x軸、y軸の濃淡変化量
- 変化強度(ノルム)
- 「x軸、y軸の濃淡変化量」から勾配方向角を推定
- arntan関数を利用
- 勾配方向を垂直(UD)、水平(LR)、斜め右上から右下(RULD)、斜め左上から右下の4パターンに丸め。
- 勾配方向角に応じて極大値評価をして非極大値だったら「変化強度(ノルム)」 を0値埋め
- Hysteresis Threshold
- High以上は白
- Low未満は黒
- High-Lowの間の場合は周辺を探査し、エッジが居れば白、いなければ黒
- 画像出力
フクさん
今回は、これをScilabで実現する。
Scilabコード
フクさん
Scilabコードは以下になる。
convolution2d.sci
function out = convolution2d(img, Kernel)
[m, n] = size(Kernel); // カーネルサイズ取得
// カーネル中心からみた幅
dy = int32((m-1)/2); // カーネル上下幅
dx = int32((n-1)/2); // カーネル左右幅
[h, w] = size(img); // イメージサイズ
out = zeros(h, w); // 出力用イメージ
// 畳み込み
for y = dy+1:(h - dy)
for x = dx+1:(w-dx)
out(y, x) = sum( double(img(y-dy:y+dy, x-dx:x+dx)).*Kernel );
end
end
endfunctionnon_maximum_suppression.sci
function out = non_maximum_suppression(G, theta)
[h, w] = size(G);
out = G;
// 勾配方向を4方向(LR,UD,RULD,LURD)に近似
theta( -22.5 <= theta & theta < 22.5) = 0; // LR ─
theta( 22.5 <= theta & theta < 67.5) = 45; // RULD ／
theta( 67.5 <= theta & theta < 112.5) = 90; // UD │
theta( 112.5 <= theta & theta < 157.5) = 135; // LURD ＼
theta( 157.5 <= theta & theta < 180.0) = 0; // LR ─
theta(-180.0 <= theta & theta < -157.5) = 0; // LR ─
theta(-157.5 <= theta & theta < -112.5) = 45; // RULD ／
theta(-112.5 <= theta & theta < -67.5) = 90; // UD │
theta( -67.5 <= theta & theta < -22.5) = 135; // LURD ＼
// 現画素の勾配方向に接する2つの画素値を比較し、現画素が極大値でなければ0にする。
for y = 2:(h - 1)
for x = 2:(w - 1)
if theta(y,x) == 0 // LR ─
if (G(y,x) < G(y,x+1)) | (G(y,x) < G(y,x-1))
out(y,x) = 0;
end
elseif theta(y,x) == 45 // RULD ／
if (G(y,x) < G(y-1,x+1)) | (G(y,x) < G(y+1,x-1))
out(y,x) = 0;
end
elseif theta(y,x) == 90 // UD │
if (G(y,x) < G(y+1,x)) | (G(y,x) < G(y-1,x))
out(y,x) = 0;
end
else // LURD ＼
if (G(y,x) < G(y+1,x+1)) | (G(y,x) < G(y-1,x-1))
out(y,x) = 0;
end
end
end
end
endfunction
hysteresis_threshold.sci
function out = hysteresis_threshold(img, low, high, r)
[h, w] = size(img);
out = img;
for y = (1+r):(h-r-1)
for x = (1+r):(w-r-1)
// 最大閾値より大きければ「エッジ」
if img(y,x) >= high
out(y,x) = 255;
// 最小閾値より小さければ「非エッジ」
elseif img(y,x) < low
out(y,x) = 0;
// 最小閾値と最大閾値の間で、半径rの範囲内に「エッジ」が1つでもあればエッジと判定
else
if max(img(y-r:y+r+1, x-r:x+r+1)) >= high
out(y,x) = 255;
else
out(y,x) = 0;
end
end
end
end
endfunction
canny_test.sci
function [] = canny_test()
// 入力画像の読み込み
img = imread('dog.jpg');
r = img(:,:,1);
g = img(:,:,2);
b = img(:,:,3);
// ガウシアンフィルタ用のkernel
kernel_gauss = [ 1/16 2/16 1/16;...
2/16 4/16 2/16;...
1/16 2/16 1/16];
// Sobelフィルタ用のKernel
kernel_sx = [-1 0 1;...
-2 0 2;...
-1 0 1];
kernel_sy = kernel_sx';
// SDTVグレースケール
gray_sdtv = uint8([0.2990 * double(r) ...
+ 0.5870 * double(g) + 0.1140 * double(b) ]);
// ガウシアンフィルタ
img_g = convolution2d(gray_sdtv, kernel_gauss);
// Sobelフィルタ
Gx = convolution2d(img_g, kernel_sx);
Gy = convolution2d(img_g, kernel_sy);
// 勾配強度と角度
G = sqrt( Gx.^2 + Gy.^2 );
theta = atan(Gy, Gx) * 180 / %pi;
// 極大値以外を除去(Non maximum Suppression)
G_nms = non_maximum_suppression(G, theta);
// Hysteresis Threshold(最小閾値除去、最大閾値以上を残し、且つそこと繋がっている最小閾値以上を残す)
G_canny = hysteresis_threshold(G_nms, 30, 65, 1); // 0 or 255に2値化
imwrite(uint8(G_canny),'dog_canny.jpg');
endfunction処理結果
フクさん
処理結果は以下。
メモリ不足問題
フクさん
ちなみに環境依存ではあるようだが、Scilabでこの処理をさせる場合、
メモリ不足になって処理が完了できないことがある。
この場合、Scilabのメニューの「編集」から「設定」選んで、Javaピープメモリを増加方向に調整してみると解決することが多い。
フクさん
画面上にも書いてあるが、再起動で反映される設定である点に注意だ。
太郎くん
なんかマニアックな問題を抱えてるな・・・。
考察
太郎くん
なんかメモリ不足問題とかあるようだけど、
処理としてはちゃんと2値化出来てるね。
ソースコードもMATLABとほぼ一緒だし。
フクさん
まぁ、今回も、コード自体はMATLABからコピペしたようなもんだしね。
まとめ
フクさん
まとめだよ。
- 非極大値抑制にHysteresis Thresholdを加えた、Canny法による2値化をScilabで実施。
- 基本的にはいままでのコードを再利用。
- 環境によってはメモリ不足に陥るのでJavaヒープメモリを調整する必要がある。
バックナンバーはこちら。
コメント