2012/11/14

HOG特徴量の計算

概要

HOG(Histogram of Oriented Gradients)を知らないゴリラはいない。
中部大学 藤吉研究室 http://www.vision.cs.chubu.ac.jp/joint_hog/pdf/HOG+Boosting_LN.pdf
東工大 画像解析論 http://www.isl.titech.ac.jp/~nagahashilab/member/longb/imageanalysis/LectureNotes/ImageAnalysis07.pdf 
以下では、ある点の近傍を1ブロックとして、そのHOG特徴量のみを計算する。
一般的なHOGの計算では画像領域を格子状のブロックに分け、各ブロックの位置(格子点)で特徴量を計算してそれらを連結する。
しかし、これを素朴に実装すると問題が発生する。

たとえば人物検出で探索窓をスライドさせるとき、計算に用いるセルは数ピクセル単位で移動することになる。
また、ランダムに分布する点に対して、その近傍をブロックとしてHOGを求めたいこともある。
そんなとき、
  1. ブロック内を複数のセルに分割、
  2. そのセル内の勾配ヒストグラムを計算、 ← ここ!!
  3. それらを連結、
  4. 正規化する
これをそれぞれのブロックに対して素朴に計算すると、
「ここ!!」のせいでなかなかの計算コストがかかっちまうってわけだ。
そこで
Integral histogram for fast HoG feature calculation http://avresearch.wordpress.com/2011/08/05/integral-histogram-for-fast-hog-feature-calculation/
を参考に、インテグラルヒストグラムを用いて高速に計算しよう。


インテグラル??という人はググってもらうとして、
積分画像をどうやって使うんだ??と考えた君へ。

答え:
ヒストグラムのビンごとに積分画像を作る。
実にクレバーな方法だ。(優れたゴリラの仕業だろう)

実装

まず積分画像の生成部分。
// ヒストグラムのビン数
#define N_BIN 9
// 何度ずつに分けて投票するか(分解能)
#define THETA (180 / N_BIN)

// 積分画像生成
vector<Mat> calculateIntegralHOG(const Mat& image) {
    // X, Y方向に微分
    Mat xsobel, ysobel;
    Sobel(image, xsobel, CV_32F, 1, 0);
    Sobel(image, ysobel, CV_32F, 0, 1);

    // 角度別の画像を生成しておく
    vector<Mat> bins(N_BIN);
    for (int i = 0; i < N_BIN; i++)
        bins[i] = Mat::zeros(image.size(), CV_32F);

    // X, Y微分画像を勾配方向と強度に変換
    Mat Imag, Iang;
    cartToPolar(xsobel, ysobel, Imag, Iang, true);
    // 勾配方向を[0, 180)にする
    add(Iang, Scalar(180), Iang, Iang < 0);
    add(Iang, Scalar(-180), Iang, Iang >= 180);
    // 勾配方向を[0, 1, ..., 8]にする準備(まだfloat)
    Iang /= THETA;

    // 勾配方向を強度で重みをつけて、角度別に投票する
    for (int y = 0; y < image.rows; y++) {
        for (int x = 0; x < image.cols; x++) {
            int ind = Iang.at<float>(y, x);
            bins[ind].at<float>(y, x) += Imag.at<float>(y, x);
        }
    }

    // 角度別に積分画像生成
    vector<Mat> integrals(N_BIN);
    for (int i = 0; i < N_BIN; i++) {
        // 積分画像をつくる、OpenCVの関数がある
        integral(bins[i], integrals[i]);
    }

    return integrals;
}

// ある矩形領域の勾配ヒストグラムを求める
// ここでいう矩形はHOG特徴量のセルに該当
void calculateHOGInCell(Mat& hogCell, Rect roi, const vector<Mat>& integrals) {
    int x0 = roi.x, y0 = roi.y;
    int x1 = x0 + roi.width, y1 = y0 + roi.height;
    for (int i = 0; i < N_BIN; i++) {
        Mat integral = integrals[i];
        float a = integral.at<double>(y0, x0);
        float b = integral.at<double>(y1, x1);
        float c = integral.at<double>(y0, x1);
        float d = integral.at<double>(y1, x0);
        hogCell.at<float>(0, i) = (a + b) - (c + d);
    }
}
以上で積分画像に関する部分は終了だ。
以下はその呼び出し部分だ。
// セルの大きさ(ピクセル数)
#define CELL_SIZE 20
// ブロックの大きさ(セル数)奇数
#define BLOCK_SIZE 3
// ブロックの大きさの半分(ピクセル数)
#define R (CELL_SIZE*(BLOCK_SIZE)*0.5)

// HOG特徴量を計算する
// pt: ブロックの中心点
Mat getHOG(Point pt, const vector<Mat>& integrals) {
    // ブロックが画像からはみ出していないか確認
    if (pt.x - R < 0 ||
        pt.y - R < 0 ||
        pt.x + R >= integrals[0].cols ||
        pt.y + R >= integrals[0].rows
    ) {
            return Mat();
    }

    // 与点を中心としたブロックで、
    // セルごとに勾配ヒストグラムを求めて連結
    Mat hist(Size(N_BIN*BLOCK_SIZE*BLOCK_SIZE, 1), CV_32F);
    Point tl(0, pt.y - R);
    int c = 0;
    for (int i = 0; i < BLOCK_SIZE; i++) {
        tl.x = pt.x - R;
        for (int j = 0; j < BLOCK_SIZE; j++) {
            calculateHOGInCell(hist.colRange(c, c+N_BIN),
                               Rect(tl, tl+Point(CELL_SIZE, CELL_SIZE)),
                               integrals);
            tl.x += CELL_SIZE;
            c += N_BIN;
        }
        tl.y += CELL_SIZE;
    }
    // L2ノルムで正規化
    normalize(hist, hist, 1, 0, NORM_L2);
    return hist;
}

int main() {
    // 画像をグレイスケールで読み込む
    string fileName = "images\\lena.jpg";
    Mat originalImage = imread(fileName, CV_LOAD_IMAGE_GRAYSCALE);

    // 積分画像生成
    vector<Mat> integrals = calculateIntegralHOG(originalImage);
    // ある点(x, y)のHOG特徴量を求めるには
    // Mat hist = getHOG(Point(x, y), integrals);
    // とする。histはSize(81, 1) CV_32FのMat


    /* ****************** *
     * 以下、表示のための処理
     * ****************** */

    // 表示用画像を用意(半分の輝度に)
    Mat image = originalImage.clone();
    image *= 0.5;
    
    // 格子点でHOG計算
    Mat meanHOGInBlock(Size(N_BIN, 1), CV_32F);
    for (int y = CELL_SIZE/2; y < image.rows; y += CELL_SIZE) {
        for (int x = CELL_SIZE/2; x < image.cols; x += CELL_SIZE) {
            // (x, y)でのHOGを取得
            Mat hist = getHOG(Point(x, y), integrals);
            // ブロックが画像からはみ出ていたら continue
            if (hist.empty()) continue;

            // ブロックごとに勾配方向ヒストグラム生成
            meanHOGInBlock = Scalar(0);
            for (int i = 0; i < N_BIN; i++) {
                for (int j = 0; j < BLOCK_SIZE*BLOCK_SIZE; j++) {
                    meanHOGInBlock.at<float>(0, i) += hist.at<float>(0, i+j*N_BIN);
                }
            }
            // L2ノルムで正規化(強い方向が強調される)
            normalize(meanHOGInBlock, meanHOGInBlock, 1, 0, CV_L2);

            // 角度ごとに線を描画
            Point center(x, y);
            for (int i = 0; i < N_BIN; i++) {
                double theta = (i * THETA + 90.0 ) * CV_PI / 180.0;
                Point rd(CELL_SIZE*0.5*cos(theta), CELL_SIZE*0.5*sin(theta));
                Point rp = center -   rd;
                Point lp = center -  -rd;
                line(image, rp, lp, Scalar(255*meanHOGInBlock.at<float>(0, i), 255, 255));
            }
        }
    }

    // 表示
    imshow("out", image);
    waitKey(0);

    return 0;
}

結果

格子点でHOGを計算しない場合を想定して、積分画像を用いた。
でも出力を確認する上では、格子点で計算するとわかりやすい。

以下はHOG特徴量を計算した結果である。
といってもHOGの結果それ自体ではない。
なぜならHOG特徴量は次元数が大きいからだ。
20度ずつで投票する場合、特徴量はセル数×9次元であり、そのままでは可視化できない。
(君が100次元程度の世界に住んでいれば話は別だが)

ということで本プログラムでは、
  • セルごとに計算したヒストグラムを(セル数×9次元)、
  • ブロックごとに足し合わせてしまい(9次元)、
  • ブロックごとに強い勾配を求めて表示している。
これはセルが20x20ピクセル、ブロックが3x3セルでHOGを計算した結果である。

セル20x20、ブロック3x3

比較対象として、ブロックの大きさを1x1としたものを見てみよう。
セル20x20、ブロック1x1
これは、単純にブロックの範囲で輝度勾配を調べていることになる。
(もはやHOGでもなんでもない)
表示するだけではこれのほうがおもしろい、なんていわないで。


セル30x30、ブロック3x3
これはセルのサイズを大きくしたものだ。


1 件のコメント:

Sebastian Rasilla さんのコメント...
このコメントは投稿者によって削除されました。
Related Posts Plugin for WordPress, Blogger...