スカスカな画像のエッジ抽出方法

Question

前提条件

例えば以下の画像のエッジを抽出する場合、OpenCVのCannyを利用する事で抽出出来ます。

img = cv2.imread("img1.png")
gray_img = cv2.cvtColor(img, cv2.COLOR_BGR2GRAY)
canny_img = cv2.Canny(gray_img, 100,200)

実現したいこと

スカスカの画像を入力した場合に、上記内容と同じようなエッジを抽出したいです。

人間の目で見ると境界線を引けそうに見えるのですが、どのような考え方で抽出すれば良いのかが思い付きません
なお、このスカスカな画像の点の間隔は、実際には一定ではありません
モルフォロジーは試したのですが、隙間が大きい箇所については膨張処理では補いきれませんでした。

用途

LIDARの点群を平面に透視投影した際の、色の境目からエッジを抽出したい
→ LIDARに映っている物の輪郭を画像として取得したい

良い方法などありましたら、ご教示いただければ幸いです。
よろしくお願いいたします。

Comments

[退会済みユーザー]

おもしろそうですね。これの「黒ではない部分」は規則正しく並んでいるのでしょうか（○○px単位で方眼状など）。

[fana]

LIDARの点群のデータなのであれば，
LIDARから見れば密（というか点同士の "隣接関係" がわかるデータ）なのでは？

（なので，その方向の隣接関係を用いて並べた絵をつくれば，隙間の無い絵になるのでは）

[yak_jun]

> これの「黒ではない部分」は規則正しく並んでいるのでしょうか（○○px単位で方眼状など）。
こちらに関しましては、並びは不定になります。
（サンプルの絵を作るのに雑なプログラムを書いた結果方眼状に並びました…）

> LIDARから見れば密（というか点同士の "隣接関係" がわかるデータ）なのでは？
確かにそうとも言えるかもしれません。
この観点でも考察を進めてみようと思います。
ありがとうございます。

[fana]

例えば，ある横方向の１ライン分のデータについて言えば，
「横方向にφ度単位でスキャンした」みたいなデータなんじゃないかな，と．
であれば，「作る画像の横方向の１pixelがφ度相当」という話で絵を作れば密なものになるであろう，と．
（それでも LIDAR だから，「データが無い」部分が存在し得るので，その点は別途何かしらうまく扱う必要はある）

Answer 1

アイデア自体はすぐ出てましたが、書き込む時間を待っていたらもうクローズされてしまってますね…。
参考まで。

---

完成形
以下に掲載するコードからいらない処理を削れば、恐らく見やすさを保ちながら全てこみこみ10-20行くらいまで短くできると思います。

---

考え方
基本的な考え方は、以下の通りです。

a. HSVならHを使って○○色っぽさの尺度で分離できる
b. ボカせば大体均せる
c. カッチリするなら3値化してやればいい

順を追って説明
0. [ダメな見本]読み込んだ画像をぼかした後にHSVで分割してエッジ検出
（hからそのままエッジ検出は筋が良くない）

ここから無理やりやるなら、閾値で白っぽいところを抜いてモルフォロジー＋scikit-imageでスケルトン検出とかでしょうか。

1. 読み込んだ画像をgaussianした後にHSVで分割して鮮やかに
（少し良くなりそう）

2. Kmeansで三値化（画像処理で言えば減色処理、統計的にはクラスタリング）する
（だいぶいい感じ）

3. hに対してLaplacianでエッジ検出
（BGRに対してエッジ検出をすると白だけでなくBやRの筋が残るのでhのみで処理が推奨です）
Cannyもいいかもしれませんね。

import numpy as np
import cv2

def kmeans_process(img):
    #######
    # Kmeans

    # http://labs.eecs.tottori-u.ac.jp/sd/Member/oyamada/OpenCV/html/py_tutorials/py_ml/py_kmeans/py_kmeans_opencv/py_kmeans_opencv.html
    Z = img_mod.reshape((-1,3))

    # convert to np.float32
    Z = np.float32(Z)

    # define criteria, number of clusters(K) and apply kmeans()
    criteria = (cv2.TERM_CRITERIA_EPS + cv2.TERM_CRITERIA_MAX_ITER, 10, 1.0)
    K = 3
    ret,label,center=cv2.kmeans(Z,K,None,criteria,10,cv2.KMEANS_RANDOM_CENTERS)

    # Now convert back into uint8, and make original image
    center = np.uint8(center)
    res = center[label.flatten()]
    img_kmeans = res.reshape((img_mod.shape))

    cv2.imshow('Kmeans',img_kmeans)
    return img_kmeans

###################################
# Minimal process
###################################

#######
# Read image
img = cv2.imread("img.png")
# cv2.imshow("img",img)

#######
# Gaussian process
img_mod = cv2.GaussianBlur(img, (15, 15), 3)
# cv2.imshow("img_gauss",img_mod)
# cv2.imwrite("img_gauss.png",img_mod.astype(np.uint8))
#######
# Split it!
h,s,v = cv2.split(cv2.cvtColor(img_mod,cv2.COLOR_BGR2HSV))

# # Not goot method
# img_not_good = cv2.Laplacian(h,cv2.CV_64F)
# cv2.imwrite("img_not_good.png",img_not_good.astype(np.uint8))

###################################
# Optional process
###################################
#######
# Multiply it!
mul = 10
s = np.clip(s*mul,0,255).astype(np.uint8)
v = np.clip(v*mul,0,255).astype(np.uint8)
img_mod = cv2.cvtColor(cv2.merge((h,s,v)),cv2.COLOR_HSV2BGR)
# cv2.imshow("img_vivid",img_mod)
# cv2.imwrite("img_vivid.png",img_mod.astype(np.uint8))

#######
# Kmeans
img_kmeans = kmeans_process(img_mod)
# cv2.imshow("img_kmeans",img_kmeans)
# cv2.imwrite("img_kmeans.png",img_kmeans.astype(np.uint8))


#######
# Split it!
h,s,v = cv2.split(cv2.cvtColor(img_kmeans,cv2.COLOR_BGR2HSV))

#######
# Edge detection
img_edge = cv2.Laplacian(h,cv2.CV_64F)

cv2.imshow('img_edge',img_edge)
# cv2.imwrite("img_edge.png",img_edge.astype(np.uint8))
cv2.waitKey(0)

Answer 2

（「質問への追記・修正」側にて，そもそも密な絵を作れるのではないか？という疑問を投げかけておりますが，それはそれとして）

提示されたようなスカスカなデータを相手にするとしたら，原始的な話をすれば，下図のようにエッジ検出に使うカーネルを適当なでかさにして，各マスの値を求む際には，そのマスの中にある「有効な」画素値の平均等を用いればよいでしょう（つまり黒の画素は無視して色がある画素だけを扱う）．

↑の図に描いたカーネルはかなり極端なでかさですが，まぁ相応のサイズにすればよいかと．

で，ほんとにでかいカーネルで計算するのは処理時間とかが嫌な感じでしょうから
【カーネルをでかくする ≒ 画像を縮小する】という方向で何かしら考えてやれば良いのではないかと．
（例えば，縮小結果の画素値を決める際に，元の画像の黒画素を無視するような）

---

とりあえず原始的なのを何の工夫も無くやってみた．
色を見るのは面倒なので，画像はグレースケールで扱った．
画像は左から

• グレースケールの原画
• でかいSobel (X方向) の結果をてきとーに可視化したもの
• 同 y方向
• {x方向, y方向}の２つから得たエッジ強度をてきとーに可視化したもの

カーネルがでかい分だけ，結果も太いが．

//作業関数．正方形領域内の，輝度値0でない画素群の輝度平均値を算出．
double AreaAve( const cv::Mat &Src8U, int left, int top, int S )
{
    unsigned int Sum=0, n=0;
    for( int dy=0; dy<S; ++dy )
    {
        const unsigned char *p = Src8U.ptr<unsigned char>( top + dy, left );
        for( int dx=0; dx<S; ++dx, ++p )
        {   if( *p ){   Sum += *p;  ++n;    }   }
    }
    return ( n ? double(Sum)/n : 0 );
}

//32bit float 画像をてきとーに表示するためのヘルパ．
void Show32F( const std::string &WndName, cv::Mat &Img )
{
    double Min,Max;
    cv::minMaxLoc( Img, &Min,&Max );
    cv::Mat ShowImg( Img.size(), CV_8U );
    double alpha = 255 / (Max-Min);
    double beta = -Min * alpha;
    Img.convertTo( ShowImg, ShowImg.type(), alpha,beta );
    cv::imshow( WndName, ShowImg );
}

//main
int main()
{
    cv::Mat SrcImg = cv::imread( "Dots.png", cv::IMREAD_GRAYSCALE );
    if( SrcImg.empty() )return 0;

    const int S = 7;    //1マスのサイズ
    const int KernelR = (S * 3) / 2;
    const int dx[8] = { -1,0,1, -1,1, -1,0,1 };
    const int dy[8] = { -1,-1,-1, 0,0, 1,1,1 };
    const int wx[8] = { -1,0,1, -2,2, -1,0,1 };
    const int wy[8] = { -1,-2,-1, 0,0, 1,2,1 };

    cv::Mat SX = cv::Mat::zeros( SrcImg.size(), CV_32F );
    cv::Mat SY = cv::Mat::zeros( SrcImg.size(), CV_32F );
    cv::Mat Mag = cv::Mat::zeros( SrcImg.size(), CV_32F );

    for( int cy=KernelR; cy+KernelR<SrcImg.rows; ++cy )
    {
        float *pSY = SY.ptr<float>( cy );
        float *pSX = SX.ptr<float>( cy );
        float *pMag = Mag.ptr<float>( cy );
        for( int cx=KernelR; cx+KernelR<SrcImg.cols; ++cx )
        {
            double Xdir=0, Ydir=0;
            for( int i=0; i<8; ++i )
            {
                double Ave = AreaAve( SrcImg,  cx + dx[i]*S - S/2, cy + dy[i]*S - S/2, S );
                Xdir += Ave * wx[i];
                Ydir += Ave * wy[i];
            }
            pSX[cx] = (float)(Xdir);
            pSY[cx] = (float)(Ydir);
            pMag[cx] = (float)sqrt( Xdir*Xdir + Ydir*Ydir );
        }
    }

    cv::imshow( "Src", SrcImg );
    Show32F( "SX", SX );
    Show32F( "SY", SY );
    Show32F( "Mag", Mag );
    cv::waitKey();
    return 0;
}

Comments

[yak_jun]

なるほど、結局疎であるという事は解像度低いのと変わらないと言うのはその通りと思っており、その頭で勝手にカーネルのストライドをカーネル幅と同一で考えていて、その頭で固定されておりました。

ストライドは1のままで動かせば線は太くぼやけた結果になりますが、こういう検知結果が得られるのは大きな前進となりそうです。(計算量の多さに関しては要検討ですが。。。)
ありがとうございました。

Answer 3

モルフォロジーは試したのですが、隙間が大きい箇所については膨張処理では補いきれませんでした。

kernelを大きくしてみるとか、iterationsを大きくしてみるとかのパラメータ調整で埋めてしまうこと自体は可能なはずです。

綺麗にエッジが取れるかはまた別の問題と言えるでしょう……。

Comments

[yak_jun]

回答ありがとうございます。
そうですね、埋めてしまう事自体は可能なのですが、どちらかと言うと隙間を埋めるというのは手段でしかなく、他のやり方も含めて色々と考えている状況であります。