のどあめ

ググってもなかなか出てこなかった情報を垂れ流すブログ。

OpenCV 3.0で人物追跡したい!

最近社内でライトニングトークなるものがありました。 そのうち画像を使う機会がありそうで復習したかった半分、 OpenCV3.0使ってみたかった半分で単カメラの人物追跡をテーマに軽く解説をしてきました。

研究室時代ではいろんな論文の実装を組み込んで何とか実装していた単カメラ人物追跡ですが、 OpenCV3.0を使えば、100行くらいで人物追跡できます。素晴らしい限りです。

今回は、その実装について説明したいと思います。

ただし、以下の点には注意です。

  • 理論などはすべて省略します
  • 高速化は考えません(本当は大事)
  • 使っている技術は古いです(2005-2006くらいの技術)

用いた環境

  • C++11
  • OpenCV 3.0
  • boost 1.58 (ファイル読み込みのみ)

人物追跡

ここでいう人物追跡とは、単一のカメラ映像に写る人物の軌跡を取得することです。 人物検出をして、その結果を追跡することで実現できます。

人物検出

画像の中から人物の矩形を取得します。

OpenCVでは、cv::HOGDescriptor::getDefaultPeopleDetectorというズバリなものが用意されていますので、 今回はこれを利用します。

(人物)追跡

人物検出の結果を起点に、そのものが映像中に動いた奇跡を追跡します。 (人物以外でも使えます)

OpenCVでは、cv::Trackerというインタフェースが用意されていますので、 今回はこれを利用します。

実装

実装はこんな感じです。 ykicisk/tracking · GitHub

#include <iostream>
#include <opencv2/opencv.hpp>
#include <opencv2/tracking/tracker.hpp>
#include <boost/filesystem.hpp>

const cv::Size MAX_DETECT_SIZE = cv::Size(100, 200);
const int MAX_MISS_FRAME = 10;
const double MIN_NEW_DETECT_INTERSECTION_RATE = 0.5;

class MyTracker {
private:
    static int next_id;
    int id;
    int n_miss_frame = 0;
    cv::Rect2d rect;
    cv::Ptr<cv::Tracker> cv_tracker;
public:
    // フレーム画像と追跡対象(Rect)で初期化
    MyTracker(const cv::Mat& _frame, const cv::Rect2d& _rect) 
        : id(next_id++), rect(_rect) {
        cv_tracker = cv::Tracker::create("BOOSTING"); //  or "MIL"
        cv_tracker->init(_frame, _rect);
    }
    // 次フレームを入力にして、追跡対象の追跡(true)
    // MAX_MISS_FRAME以上検出が登録されていない場合は追跡失敗(false)
    bool update(const cv::Mat& _frame){
        n_miss_frame++;
        return cv_tracker->update(_frame, rect) && n_miss_frame < MAX_MISS_FRAME;
    }
    // 新しい検出(Rect)を登録。現在位置と近ければ受理してn_miss_frameをリセット(true)
    // そうでなければ(false)
    bool registerNewDetect(const cv::Rect2d& _new_detect){
        double intersection_rate = 1.0 * (_new_detect & rect).area() / (_new_detect | rect).area();
        bool is_registered = intersection_rate > MIN_NEW_DETECT_INTERSECTION_RATE;
        if (is_registered) n_miss_frame = 0;
        return is_registered;
    }
    // trackerの現在地を_imageに書き込む
    void draw(cv::Mat& _image) const{
        cv::rectangle(_image, rect, cv::Scalar(255, 0, 0), 2, 1);
        cv::putText(_image, cv::format("%03d", id), cv::Point(rect.x + 5, rect.y + 17), 
                cv::FONT_HERSHEY_SIMPLEX, 0.5, cv::Scalar(255,0,0), 1, CV_AA);
    }
};
int MyTracker::next_id = 0;


int main(int argc, char* argv[]){
    if(argc != 2){
        std::cout << "usage: " << argv[0] << " videodir" << std::endl;
        exit(1);
    }
    // フレーム画像列のパスを取得
    namespace fs = boost::filesystem;
    std::vector<std::string> frame_paths;
    for(auto it = fs::directory_iterator(argv[1]); it != fs::directory_iterator(); ++it){
        frame_paths.push_back(it->path().string());
    }
    // detector, trackerの宣言
    cv::HOGDescriptor detector;
    detector.setSVMDetector(cv::HOGDescriptor::getDefaultPeopleDetector());
    std::vector<MyTracker> trackers;
    // 1フレームずつループ
    for (auto& frame_path : frame_paths){
        std::cout << "frame : " << frame_path << std::endl;
        cv::Mat frame = cv::imread(frame_path);
        // 人物検出
        std::vector<cv::Rect> detections;
        detector.detectMultiScale(frame, detections);
        // trackerの更新(追跡に失敗したら削除)
        for (auto t_it = trackers.begin(); t_it != trackers.end();){
            t_it = (t_it->update(frame)) ? std::next(t_it) : trackers.erase(t_it);
        }
        // 新しい検出があればそれを起点にtrackerを作成。(既存Trackerに重なる検出は無視)
        for(auto& d_rect : detections){
            if (d_rect.size().area() > MAX_DETECT_SIZE.area()) continue;
            bool exists = std::any_of(trackers.begin(), trackers.end(), 
                    [&d_rect](MyTracker& t){return t.registerNewDetect(d_rect);});
            if(!exists) trackers.push_back(MyTracker(frame, d_rect));
        }
        // 人物追跡と人物検出の結果を表示
        cv::Mat image = frame.clone();
        for(auto& t : trackers) t.draw(image);
        for(auto& d_rect : detections) cv::rectangle(image, d_rect, cv::Scalar(0, 255, 0), 2, 1);
        cv::imshow("demo", image);
        cv::waitKey(1);
    }
    return 0;
}

動画のフレームごとに人物検出と追跡を繰り返して実現します。

そのままだとTrackerがどんどん増えていくので、 Trackerと重なる検出がでてこなくなったTrackerは消すようにしています。

実験

PETS2009 S2.L1のView001で使って動かしてみました。

緑が人物検出の結果、青が追跡結果です。

考察

思ったよりイマイチですね。

人物検出の誤検出・未検出について

OpenCV3.0のdetectMultiScaleを使った人物検出ではこれが限界だと思います。

もっとリッチな人物検出手法をつかったり、 シーンの知識(どの当たりにどんな大きさの人が写るなど)を利用すればもっと検出がよくなると思います。

追跡ミスについて

汎用追跡器を使ったこともありますが、オクルージョンにめちゃくちゃ弱いです。 真ん中のポールを通過した人はもれなくおかしなことになっていますね。

これに解決するためには、もっと人物に特化した追跡器を使ったり、 関連付けられた人物画像を一旦バラバラにして最適化で同一人物をくっつけなおしたりすると良いと思います。

まとめ

  • OpenCV3.0を使えば簡単に単カメラ人物追跡できる
  • 結果はいまいち。実際のカメラは更にいろいろあるので更におかしい事が起こる
  • 人物追跡にかぎらず、この辺りの分野は闇が深い

次回があれば、今回取得した人物画像列について、 複数カメラでの対応付けを行う方法について説明します。

でも次回はありません。

Web画像検索で画像を収集する

沢山の画像を利用して何やらしたい時があります.
たいていはImageNetやFlickrで収集すれば良いのですが, 時には普通のWebで画像検索した結果を利用したいこともあります.

以前に,Web画像検索した結果を収集したい時の最適解はなんだろうと検討したことがありましたので, その結果をメモついでに記事にしました.

画像収集に利用するサービスの検討

Google画像検索API

制限が厳しくなり,同じクエリに対して10枚×10ページの合計100枚までしか取得できないようです. 画像を収集するという意味では使いにくいですね.

基本的な使い方は

Bing Search API

同じクエリでも1000件以上収集できます(上限は確認していません).
無料枠でも月5000回Getできるので画像収集として使えますね.
今回はBing Search APIが最適解という結論に至りました.

API利用までの流れなどは以下の記事で詳しく紹介されています.
Bing Search APIを使ってWeb検索を行うには(Json編) | garicchi.com

画像収集スクリプト(Python)

今回はWindowsで利用したのでcp932とか書いてますが, 違う環境の場合は適宜修正してください.

# coding: utf-8
import sys, os
import argparse
import json
import urllib
import urllib2
import requests

#proxy setting
proxy_dict =  {"http":"your.proxy:8080"}
proxy = urllib2.ProxyHandler(proxy_dict)
opener = urllib2.build_opener(proxy)
urllib2.install_opener(opener)

#primary account key
api_key="*************************************"

def download_urllist(urllist,outpath):
    for url in urllist:
        try:
            print"url:"+url
            with open(outpath+"/"+os.path.basename(url),'wb') as f:
                img=urllib2.urlopen(url,timeout=5).read()
                f.write(img)
        except urllib2.URLError:
            print("URLError")
        except IOError:
            print("IOError")
        except UnicodeEncodeError:
            print("EncodeError")
        except OSError:
            print("OSError")

if __name__ == "__main__":
    parser = argparse.ArgumentParser()
    parser.add_argument("-q","--query", required=True, 
            help="query")
    parser.add_argument("-o","--outpath", required=True,
            help="output path")
    args = parser.parse_args()


    #Bing Image Search
    query = args.query.decode('cp932').encode('utf-8')
    query = urllib2.quote(query)

    step = 20
    num = 50
    
    url_param_dict={
            "Query":"'%s'"%(query),
            "Market":"'ja-JP'",
            "Adult":"'Off'",
            }
    url_param_base = urllib.urlencode(url_param_dict)
    url_param_base = url_param_base + "&$format=json&$top=%d&$skip="%(num)


    for skip in range(0,num*step,num):
        url_param = url_param_base + str(skip)
        url="https://api.datamarket.azure.com/Bing/Search/Image?"+url_param
        print url

        response = requests.get(url, 
                            auth=(api_key, api_key), 
                            headers={'User-Agent': 'My API Robot'},
                            proxies=proxy_dict)
        response=response.json()

        urllist=[item['MediaUrl'] for item in response['d']['results']]
        download_urllist(urllist,args.outpath)

画像収集スクリプトは以下の記事を参考にしました.
Bing Search API を使いたいと思ったのでPythonでラッパーを作ってみた - [[ともっくす alloc] init]

使い方は

python collect_images.py -q Query -o outputdir

です.

-qで指定したクエリについて,step*num件の画像をダウンロードします.
ダウンロードした結果にはゴミも多いので,いい感じにクリーニングして利用しましょう.

OpenCV 3.0 で cvSegmentFGMask

OpenCV 3.0.0 では,IplImage,legacy関数などの1.x系のサポートがなくなっています.
つまりcvHogehogeみたいな関数は使えません.
たいていの機能はcv::Hogehogeみたいな関数が実装されているのですが,一部の関数は代替関数がなくて残念なことになる場合があります.

私の環境でも,一部のプログラムで使っていたcvSegmentFGMaskが使えなくなり,困っていました.
cvSegmentFGMaskは,ノイズっぽいmask画像をいい感じにまとめてくれる関数です.

この度cvSegmentFGMaskをOpenCV 3.0で使えるように1.x系のソースを参考にして2.x系で実装しました.

cvSegmentFGMaskの2.0系実装

cvSegmentFGMask.cpp

#include <opencv2/opencv.hpp>

namespace cv {
  void SegmentFGMask(InputArray src, OutputArray dst, bool poly1Hull0, float perimScale, Point offset=Point(0,0))
  {
    Mat mask = src.getMat().clone();
    if(mask.type() != CV_8UC1) throw std::runtime_error("src.type() != CV_8UC1");

    Mat new_mask = Mat::zeros(mask.size(), mask.type());

    std::vector<std::vector<Point>> contours;
    findContours(mask, contours,RETR_EXTERNAL, CHAIN_APPROX_SIMPLE, offset);

    const Scalar white(255);
    for (int i = 0; i < contours.size(); i++){
      double len = cv::arcLength(contours[i], true);
      double q = (mask.rows + mask.cols) / perimScale; 

      if (len >= q){
        std::vector<std::vector<Point>> tmp_contours(1);
        if (poly1Hull0){
          cv::approxPolyDP(contours[i], tmp_contours[0], 2.0, true);
        } else{
          tmp_contours[0] = contours[i];
        }
        drawContours(new_mask, tmp_contours, 0, white, -1, 8,cv::noArray(),-1,Point(-offset.x,-offset.y));
      }
    }
    new_mask.copyTo(dst);
  }
}

int main(){
  cv::VideoCapture video("./WalkByShop1cor.mpg");
  cv::BackgroundSubtractorMOG bgs;

  while(1){
    cv::Mat image;
    video >> image;

    if(image.empty()) break;
    cv::Mat mask;
    bgs(image,mask);
    cv::Mat segmented_mask;
    cv::SegmentFGMask(mask,segmented_mask,true,16.0f);
    {
      cv::Mat tmp;
      cv::hconcat(mask,segmented_mask,tmp);
      cv::imshow("mask",tmp);
      char k = cv::waitKey(1);
      if(k == 'q') break;
    }
  }
  return 0;
}

結果

左が普通の背景差分.右がcv::SegmentFGMaskをかけたものです.
f:id:ykicisk:20150203000821p:plain

今回のまとめ

必要がなければOpenCV 2.x 系に戻しましょう.

VisualStudioの理不尽な問題(の一部)に対処する

私は普段のプログラミング(C++)でVisual Studioを良く使っています.
なんだかんだ言ってもIntelliSenseや高機能なデバッガは非常に便利です.今後Visual Studioオープンソース化するということなので,楽しみであります.
参考:Microsoft takes .NET open source and cross-platform, adds new development capabilities with Visual Studio 2015, .NET 2015 and Visual Studio Online | News Center

しかし,Visual Studioを使っていると時々理不尽な問題にぶち当たり,どうしていいかわからなくなる時があります.
最近こうした問題について後輩から何度か質問されたので,その辺りをまとめました.
今回対象とするのは最新版のVisual Studio2013のみですが,古いバージョンでも同じようなことが可能かもしれません.

Visual Studioの理不尽な問題

今回扱うVisual Studioの理不尽な問題は以下のとおりです.

以上3点です.これは酷い・・・
VisualStudioを使い始めでこのような問題にぶち当たったら投げ出す人が多いと思います.

いつの間にかIntelliSenseが効かなっている

特にActive Directoryを導入した環境でVisual Studioを使っている時に良く発生する気がします.こちらに関しては,フォールバック位置を変更することで対処する事ができる場合があります.

オプション→テキスト エディター→C/C++→詳細について,
常にフォールバック位置を使用を「True」
フォールバック位置をC直下などActive Directoryの管理下出ないフォルダに指定すればOKです.
指定したあとは,Visual Studioを再起動する必要があります.

また,最初からプロジェクトをローカルに作ってしまっても大丈夫だと思います.

ソースコードが元のバージョンと異なります」と言われてデバッグ時にブレークポイントで止まらない

こちら(Visual Studioでブレーク出来ない問題)のページによると,昔からある厄介なバグのようです.こちらはデバッグ情報の形式を変更することで対処することが出来る場合があります.

プロジェクトのプロパティ→C++→全般について
デバッグ情報の形式を「C7 互換 (/Z7)」に変更します.

理由はよくわかりませんが,これで大抵治ります.

ソースコードを編集したのにデバッグ時に反映されない

これに関しては,「最小リビルドを有効にする」というコンパイルオプションが悪さをしている場合があります.

プロジェクトのプロパティ→C++→コード生成について
最小リビルドを有効にするを「いいえ (/Gm-)」にします.

今回のまとめ

今回はVisual Studioの理不尽な問題(の一部)の解決法について紹介しました.
他にも,#includeをするときに<だと補完されるのにダブルクオーテーションだと補完されない等いろいろ細かい問題もありますが,それらについてはまた機会があれば紹介したいと思います.

文句をタラタラ言っていますが,最近はPythonを書くこともできますし,なんだかんだ言ってもVisual Studioは便利なツールです.
ぜひ使いこなして快適なVisual Studioライフを送りましょう!

OpenCV 3.0.0 beta で Dense samplingする

さて,先日の記事OpenCV 3.0.0 alpha で Dense samplingする - のどあめを書いたと思ったら,すぐにOpenCV 3.0.0 betaが出てしまいました.
(リリース予定くらいは確認すべきでした)

もちろん,研究室のうちのグループ(の私だけ)は,早速OpenCV 3.0.0 betaに切り替えました.しかし,先日実装したDenseSamplingのコードは,betaではコンパイルが通りません.
そこで,今回は先日の記事と同じようにOpenCV 3.0.0 betaでもDense Samplingの実装をしたいと思います.

Dense Samplingの実装(OpenCV 3.0.0 beta版)

前回のdetectImpl関数をdetectAndComputeに変更します.
detectなどはこの関数を内部で呼んでいるみたいです.
compute関数(DescriptorExtractor的な関数?)をoverrideしてassertを吐かせる方がいい気もしますが,今回は省略.

DenseFeatureDetector.h

#pragma once

#include <opencv2/opencv.hpp>
namespace cv {
	class DenseFeatureDetector : public cv::FeatureDetector {
	public:
		struct Param{
			int xStep = 4;
			int yStep = 4;

			float initFeatureScale = 1.0f;
			int scaleLevels = 1;
			float featureScaleMul = 1.2f;

			int imgBoundX = 0;
			int imgBoundY = 0;

			bool varyXyStepWithScale = true;
			bool varyImgBoundWithScale = false;
			Param(){}
		}param;

		DenseFeatureDetector(){}
		DenseFeatureDetector(const Param& p)
			:param(p)
		{
		}

		static cv::AlgorithmInfo& _info();
		virtual cv::AlgorithmInfo* info()const override;

		virtual void detectAndCompute( 
				cv::InputArray image, cv::InputArray mask,
				std::vector<cv::KeyPoint>& keypoints,
				cv::OutputArray descriptors,
				bool useProvidedKeypoints=false );

		static cv::Ptr<FeatureDetector> create(const Param& p=Param()){
			return cv::Ptr<FeatureDetector>(new DenseFeatureDetector(p));
		}
	};
}

DenseFeatureDetector.cpp

#include "DenseFeatureDetector.h"

namespace cv{
#define CV_INIT_ALGORITHM(classname, algname, memberinit) \
    static inline ::cv::Algorithm* create##classname##_hidden() \
    { \
        return new classname; \
    } \
    \
    static inline ::cv::Ptr< ::cv::Algorithm> create##classname##_ptr_hidden() \
    { \
        return ::cv::makePtr<classname>(); \
    } \
    \
    static inline ::cv::AlgorithmInfo& classname##_info() \
    { \
        static ::cv::AlgorithmInfo classname##_info_var(algname, create##classname##_hidden); \
        return classname##_info_var; \
    } \
    \
    static ::cv::AlgorithmInfo& classname##_info_auto = classname##_info(); \
    \
    ::cv::AlgorithmInfo* classname::info() const \
    { \
        static volatile bool initialized = false; \
        \
        if( !initialized ) \
        { \
            initialized = true; \
            classname obj; \
            memberinit; \
        } \
        return &classname##_info(); \
    }

	CV_INIT_ALGORITHM(DenseFeatureDetector, "Dense",);


	void DenseFeatureDetector::detectAndCompute( 
			cv::InputArray image, cv::InputArray mask,
			std::vector<cv::KeyPoint>& keypoints,
			cv::OutputArray descriptors,
			bool useProvidedKeypoints)
	{
		if (!useProvidedKeypoints) keypoints.clear();


		cv::Mat maskMat = mask.empty() ? cv::Mat(image.size(), CV_8UC1, cv::Scalar(255)) : mask.getMat();
		cv::Mat imageMat = image.getMat();

		int width = imageMat.cols;
		int height = imageMat.rows;

		int nowBoundX = param.imgBoundX;
		int nowBoundY = param.imgBoundY;

		int nowXStep = param.xStep;
		int nowYStep = param.yStep;

		float nowFeatureSize = param.initFeatureScale;

		for (int s = 0; s < param.scaleLevels; s++){
			for (int y = nowBoundY; y < height; y += nowYStep){
				unsigned char* maskline = maskMat.ptr<unsigned char>(y);

				for (int x = nowBoundX; x < width; x += nowYStep){
					if (maskline[x] > 128){
						cv::KeyPoint kp;
						kp.pt = cv::Point(x, y);
						kp.size = nowFeatureSize;
						kp.response = 100.0f;
						keypoints.push_back(kp);
					}
				}
			}

			if (param.varyImgBoundWithScale){
				nowBoundX *= param.featureScaleMul;
				nowBoundY *= param.featureScaleMul;
			}
			if (param.varyXyStepWithScale){
				nowXStep *= param.featureScaleMul;
				nowYStep *= param.featureScaleMul;
			}
			nowFeatureSize *= param.featureScaleMul;
		}
	}
}

基本的には関数名以外変える必要はありません.
これで,先日の記事のコードを動かすことで同様の結果が得られます.

今回のまとめ

  • 大きな仕様変更がなくてよかった.

Boost.Pythonを使ってC++からSparse Codingする.

人生のうちで,1度や2度,C++のプログラムでSparseCodingしたいときがあると思います.しかし,昨今機械学習の分野で優遇されているのはPythonばかりです.

Sparse CodingについてはPythonではscikit learnを使えば簡単に使えます.
sklearn.decomposition.DictionaryLearning — scikit-learn 0.15.2 documentation

しかし,C++だとMLPACK等のライブラリを別途入れる必要があります.
MLPACKは動作環境がwindowsだったりすると,導入は非常に困難です(私は諦めました).
また,人生のうち1度や2度しか使わないのでいちいち新しいライブラリを導入するのは正直面倒です.

こんなときは,Boost.Pythonを介してscikit learnを使うという荒業を使えば,
手軽にC++からSparseCodingをすることが可能です.
これはscikit learnやscipyに実装されているものならなんでもC++から使えてしまう方法なので,非常に便利だと思います.

今回は,Boost.Pythonを使ったSparse Codingを使って,手書き数字認識に挑戦してみたいと思います.

概要

手書き数字認識

今回は,データセットとしてMNISTを利用します.
訓練画像・テスト画像とも28x28の手書き数字(0〜9)です.
噂のDeep Learningを使えば99.7%以上の精度が出るみたいです.
すごいですね(小並感)

Sparse Coding

理論面は私もよくわかってないので省略します.
参考:scikit-learnでSparse CodingとDictionary Learning -理論編- - takminの書きっぱなし備忘録

Sparse Codingによる手書き文字認識

文献[1]のインスパイアで,訓練データをそのまま辞書として使ってスパースコーディングします.
テストデータと同じクラスの訓練データの基底の係数が非ゼロになることを期待しています.
もしかすると同じような考えのもとに誰かが手書き文字認識されているかもしれませんが,あまり気にしません.

Boost.Python

本来はPythonからC/C++のコードを呼び出すものですが,逆も可能です.
Python自体にC/C++用のインタフェースがありますが,Boost.Pythonを使うとより簡単に使えます.

実装

sparse_coding.py

# -*- coding: utf-8 -*-
import numpy as np
from sklearn.decomposition import SparseCoder

def sparse_coding(dict,data,transform_n_nonzero_coefs=10,transform_algorithm='lars',transform_alpha=1.0):
    coder = SparseCoder(dictionary=dict,transform_algorithm=transform_algorithm,transform_n_nonzero_coefs=transform_n_nonzero_coefs,transform_alpha=transform_alpha)
    return coder.transform(data)

こちらの関数をC++から呼び出します.

mnist.cpp

OpenCVとBoostを利用しています.

#define BOOST_PYTHON_STATIC_LIB
#define NPY_NO_DEPRECATED_API NPY_1_7_API_VERSION

#include <python.h>
#include <numpy/arrayobject.h>
#include <boost/python.hpp>
#include <boost/python/numeric.hpp>
#include <boost/python/suite/indexing/vector_indexing_suite.hpp>

#include <boost/progress.hpp>

#include <iostream>
#include <fstream>
#include <vector>
#include <opencv2/opencv.hpp>

typedef unsigned char BYTE;
const unsigned int NUM_TRAIN = 5000;
const unsigned int NUM_TEST = 500;

inline void endian_swap(unsigned int& x) {
    x = (x>>24) | ((x<<8) & 0x00FF0000) | ((x>>8) & 0x0000FF00) | (x<<24);
}

void loadImages(const std::string& path,cv::OutputArray image,unsigned int n_max=INT_MAX) {
  std::ifstream ifs(path,std::ios::binary);
  if(!ifs){
    std::cout << path <<" is not found." << std::endl;
    exit(1);
  }

  unsigned int magic_number,num_of_items,height,width;
  ifs.read((char*)&magic_number,sizeof(unsigned int));
  ifs.read((char*)&num_of_items,sizeof(unsigned int));
  ifs.read((char*)&height,sizeof(unsigned int));
  ifs.read((char*)&width,sizeof(unsigned int));

  //big-endianなのでlittle-endianに変更
  endian_swap(magic_number);
  endian_swap(num_of_items);
  endian_swap(height);
  endian_swap(width);

  num_of_items = std::min(num_of_items,n_max);
  unsigned int n_bytes=num_of_items*height*width;
  std::vector<BYTE> buf(n_bytes);
  ifs.read((char*)&buf[0],n_bytes);

  std::cout << "=== load images ===" << std::endl;
  std::cout << "path : "<< path  << std::endl;
  std::cout << "magic_number: "<< magic_number  << std::endl;
  std::cout << "num_of_items: "<< num_of_items  << std::endl;
  std::cout << "height: "<< height  << std::endl;
  std::cout << "width: "<< width  << std::endl;

  cv::Mat m(buf);

  m = m.reshape(1,num_of_items);
  m.convertTo(m,CV_64FC1);
  m.copyTo(image);
}

void loadLabels(const std::string& path,cv::OutputArray labels,unsigned int n_max=INT_MAX){
  std::ifstream ifs(path,std::ios::binary);
  if(!ifs){
    std::cout << path <<" is not found." << std::endl;
    exit(1);
  }

  unsigned int magic_number,num_of_items;
  ifs.read((char*)&magic_number,sizeof(unsigned int));
  ifs.read((char*)&num_of_items,sizeof(unsigned int));
  endian_swap(magic_number);
  endian_swap(num_of_items);
  num_of_items = std::min(num_of_items,n_max);

  unsigned int n_bytes=num_of_items;
  std::vector<BYTE> buf(n_bytes);
  ifs.read((char*)&buf[0],n_bytes);

  std::cout << "=== load labels ===" << std::endl;
  std::cout << "path : "<< path  << std::endl;
  std::cout << "magic_number: "<< magic_number  << std::endl;
  std::cout << "num_of_items: "<< num_of_items  << std::endl;

  cv::Mat m(buf);
  m = m.reshape(1,1);
  m.convertTo(m,CV_32SC1);
  m.copyTo(labels);
}


int main(int argc,char* argv[]){
  if(argc!=2){
    std::cout << "usage : " << argv[0] << " MNIST_path" <<std::endl;
    return -1;
  }

  std::string MINIST_path(argv[1]);
  std::string train_images_path = MINIST_path+"/train-images-idx3-ubyte";
  std::string train_labels_path = MINIST_path+"/train-labels-idx1-ubyte";
  std::string test_images_path  = MINIST_path+"/t10k-images-idx3-ubyte";
  std::string test_labels_path  = MINIST_path+"/t10k-labels-idx1-ubyte";

  //dataset load
  cv::Mat train_images, train_labels;
  loadImages(train_images_path,train_images,NUM_TRAIN);
  loadLabels(train_labels_path,train_labels,NUM_TRAIN);

  cv::Mat test_images, test_labels;
  loadImages(test_images_path,test_images,NUM_TEST);
  loadLabels(test_labels_path,test_labels,NUM_TEST);

  //predictのために工夫 : [1,2,3]を[1,0,0; 0,1,0; 0,0,1]に変換
  cv::Mat dot_mat;
  for(int i = 0;i < 10;i++){
    cv::Mat mask = (train_labels == i);
    mask.convertTo(mask,CV_64FC1);
    dot_mat.push_back(mask);
  }
  dot_mat/=255;

  std::cout << "=== test ===" << std::endl;
  {
    //sparce coding
    cv::Mat X;

    try{
      boost::progress_timer t;
      namespace bp = boost::python;

      Py_Initialize();

      // おまじない
      bp::numeric::array::set_module_and_type("numpy", "ndarray");
      bp::class_<std::vector<double>>("std::vector<double>")
        .def(bp::vector_indexing_suite<std::vector<double>>());

      // 利用するpyファイルの指定
      bp::object global_ns = bp::import("__main__").attr("__dict__");
      bp::exec_file("sparse_coding.py", global_ns, global_ns);

      // cv::Matからboost::python::numeric::arrayに変換
      std::vector<double>train_images_vec = train_images.reshape(1,1);
      bp::numeric::array dict(train_images_vec);
      dict.resize(train_images.rows,train_images.cols);

      std::vector<double>test_images_vec = test_images.reshape(1,1);
      bp::numeric::array data(test_images_vec);
      data.resize(test_images.rows,test_images.cols);

      // 呼び出す関数を指定
      bp::object py_func = global_ns["sparse_coding"];

      // pythonの関数を呼び出し
      bp::object py_val = py_func(dict,data,9); 
      
      // 戻り値をcv::Matに変換
      bp::numeric::array py_arr_val = bp::extract <bp::numeric::array>(py_val);

      PyArrayObject* py_arr_obj = (PyArrayObject*)(py_arr_val.ptr());

      double* data_ptr = (double*)(PyArray_DATA(py_arr_obj));

      npy_intp* dims = PyArray_DIMS(py_arr_obj);

      std::vector<double> vec_val(data_ptr, data_ptr + (dims[0]*dims[1]));

      std::cout << vec_val.size() << std::endl;
      cv::Mat m(vec_val);
      m = m.reshape(1,dims[0]);
      m.copyTo(X);// vec_valが開放されても大丈夫なようにdeep copy

      std::cout << "calc time : ";//boost::progres_timer
    }
    catch (...){
      PyErr_Print();
    }

    //predict
    int nT = 0, nF = 0;

    int* gtl=test_labels.ptr<int>(0);
    for(int i = 0;i < test_images.rows;i++){
      std::cout << "[";
      double max_p = 0.0;
      int max_d = -1;
      for(int d = 0;d < 10;d++){
        double p = X.row(i).dot(dot_mat.row(d));
        if(p > max_p){
          max_p = p;
          max_d = d;
        }
        std::cout << p << ",";
      }
      std::cout << "] predict("<< max_d <<")";
      std::cout << " GT(" << gtl[i] << ") " << std::endl;

      (gtl[i] == max_d)? nT++ : nF++;
    }

    std::cout << "Accuracy : " << 1.0*nT/(nT+nF) << std::endl;
  }
}

全データを使うと学習に時間がかかるので,学習データ5000,テストデータ500でやっています.
また,predict部分では係数の合計を認識スコアとして使っています.

ビルド&使い方

clang++ -std=c++11 -lboost_python $(shell python-config --cflags) $(shell pkg-config --libs opencv) mnist.cpp -o mnist.out
./mnist.out data

↑でdataはデータセットを入れたディレクトリへのパスです.

python-config
||
**出力
>||
(前略)
calc time : 86.88 s
[0,0,0,0,0,0,0,0.658375,0,0,] predict(7) GT(7)
[0,0,0.539294,0.0822132,0,0,0,0,0,0,] predict(2) GT(2)
[0,0.167486,0,0,0,0,0,0,0,0,] predict(1) GT(1)
[0.426151,0,0,0,0,0,0,0,0,0,] predict(0) GT(0)
[-0.0978972,0,0,0,0.235536,0,0,0,0,0.0391382,] predict(4) GT(4)
(中略)
[0,0,0,0,0,0,0,0,0,0.426833,] predict(9) GT(9)
[0,0,0,0,0.147261,0,0,0.0295174,0.275997,0,] predict(8) GT(4)
[0.394697,0,0,0,0,0,0,0,0.0479334,0,] predict(0) GT(0)
[0,0,0,0,0,0,0.19342,0,0.0134353,0,] predict(6) GT(6)
Accuracy : 0.902

精度は90.2%になりました.
うーん・・・学習データを全部使っていないとはいえイマイチですね.
やはりオーソドックスに辞書学習と組み合わせるほうが良さそうですね.

また,sparse_codingのtransform_n_nonzero_coefsの値によって精度が大きく変わりました.
大きくても小さくてもダメなようで,今回は一番精度の良かった9を使いました(これはひどい

今回のまとめ

Boost.Pythonを使うことでPythonで実装されている関数をC++から使うことができました.
今回はMacで動かしましたが,WindowsでもPyhtonをインストールすれば動かす事ができます.

しかし,体感ですがやはりC++Pythonのデータのやり取りにはそれなりに時間がかかるようです(大きい行列だと1分とか・・・).
何かをちょこっと試したいという場合に,ライブラリのインストールと格闘とする前に試してみるくらいには使えると思います.

でも,機械学習をしたい人は最初からPythonでプログラミングしましょう.

参考文献

[1] J. Wright, A. Y. Yang, A. Ganesh, S. S. Sastry, and Y. Ma: "Robust face recognition via sparse representation." (2009)

OpenCV 3.0.0 alpha で Dense samplingする

さて,OpenCV 3.0.0 alphaがリリースされて結構たちましたね.
研究室のうちのグループでは調子に乗って一足先にグループ内ライブラリにOpenCV 3.0.0-alphaを取り入れましたが,あれがないこれがない,などと色々と困ったことになったりしました.

さて,その中のうちの一つとして「OpenCV3.0.0 alphaからDense Samplingがなくなった」というものがあります.
参考:OpenCV3.0.0-alphaの特徴抽出・マッチングまとめ - whoopsidaisies's diary

これは,Bag-of-featuresの時代は終わったということなのでしょうか・・・

まあ,Dense Samplingぐらいすぐ実装できる,と思っていたのですが,
いざやってみると結構躓いたので記事にしました.

Dense Samplingの実装

OpenCVのcv::FeatureDetectorの子クラスとして,
Dense Samplingを実装しました.

DenseFeatureDetector.h

#pragma once

#include <opencv2/opencv.hpp>
namespace cv {
  class DenseFeatureDetector : public cv::FeatureDetector {
  public:
    struct Param{
      int xStep = 4;
      int yStep = 4;

      float initFeatureScale = 1.0f;
      int scaleLevels = 1;
      float featureScaleMul = 1.2f;

      int imgBoundX = 0;
      int imgBoundY = 0;

      bool varyXyStepWithScale = true;
      bool varyImgBoundWithScale = false;
      Param(){}
    }param;

    DenseFeatureDetector(){}
    DenseFeatureDetector(const Param& p)
      :param(p)
    {
    }

    static cv::AlgorithmInfo& _info();
    virtual cv::AlgorithmInfo* info()const override;

    virtual void detectImpl(cv::InputArray image, std::vector<cv::KeyPoint>& keypoints, cv::InputArray mask = cv::noArray()) const override;

    static cv::Ptr<FeatureDetector> create(const Param& p=Param()){
      return cv::Ptr<FeatureDetector>(new DenseFeatureDetector(p));
    }
  };
}

DenseFeatureDetector.cpp

#include "DenseFeatureDetector.h"

namespace cv{
#define CV_INIT_ALGORITHM(classname, algname, memberinit) \
  static inline ::cv::Algorithm* create##classname##_hidden() \
  { \
    return new classname; \
  } \
  \
  static inline ::cv::Ptr< ::cv::Algorithm> create##classname##_ptr_hidden() \
  { \
    return ::cv::makePtr<classname>(); \
  } \
  \
  static inline ::cv::AlgorithmInfo& classname##_info() \
  { \
    static ::cv::AlgorithmInfo classname##_info_var(algname, create##classname##_hidden); \
    return classname##_info_var; \
  } \
  \
  static ::cv::AlgorithmInfo& classname##_info_auto = classname##_info(); \
  \
  ::cv::AlgorithmInfo* classname::info() const \
  { \
    static volatile bool initialized = false; \
    \
    if( !initialized ) \
    { \
      initialized = true; \
      classname obj; \
      memberinit; \
    } \
    return &classname##_info(); \
  }

  CV_INIT_ALGORITHM(DenseFeatureDetector, "Dense",);

  void DenseFeatureDetector::detectImpl(cv::InputArray image, std::vector<cv::KeyPoint>& keypoints, cv::InputArray mask) const
  {
    cv::Mat maskMat = mask.empty() ? cv::Mat(image.size(), CV_8UC1, cv::Scalar(255)) : mask.getMat();
    cv::Mat imageMat = image.getMat();

    int width = imageMat.cols;
    int height = imageMat.rows;

    int nowBoundX = param.imgBoundX;
    int nowBoundY = param.imgBoundY;

    int nowXStep = param.xStep;
    int nowYStep = param.yStep;

    float nowFeatureSize = param.initFeatureScale;

    for (int s = 0; s < param.scaleLevels; s++){
      for (int y = nowBoundY; y < height; y += nowYStep){
        unsigned char* maskline = maskMat.ptr<unsigned char>(y);

        for (int x = nowBoundX; x < width; x += nowXStep){
          if (maskline[x] > 128){
            cv::KeyPoint kp;
            kp.pt = cv::Point(x, y);
            kp.size = nowFeatureSize;
            kp.response = 100.0f;
            keypoints.push_back(kp);
          }
        }
      }

      if (param.varyImgBoundWithScale){
        nowBoundX *= param.featureScaleMul;
        nowBoundY *= param.featureScaleMul;
      }
      if (param.varyXyStepWithScale){
        nowXStep *= param.featureScaleMul;
        nowYStep *= param.featureScaleMul;
      }
      nowFeatureSize *= param.featureScaleMul;
    }
  }
}

以上です.

使い方としては,

#include "DenseFeatureDetector.h"

int main(){
  cv::Mat img = cv::imread("lena.jpg");

  cv::DenseFeatureDetector::Param p;
  p.scaleLevels = 3;
  p.featureScaleMul = 1.5;
  p.xStep = 20;
  p.yStep = 20;
  cv::Ptr<cv::FeatureDetector> detector = cv::DenseFeatureDetector::create(p);

  std::vector<cv::KeyPoint> keypoints;
  detector->detect(img,keypoints);

  for(auto& k : keypoints){
    std::cout << "Pt : "<< k.pt << " Size :" << k.size << std::endl;
  }

  return 0;
}

出力

Pt : [0, 0] Size :1
Pt : [20, 0] Size :1
Pt : [40, 0] Size :1
Pt : [60, 0] Size :1
(中略)
Pt : [180, 270] Size :2.25
Pt : [225, 270] Size :2.25
Pt : [270, 270] Size :2.25
Pt : [315, 270] Size :2.25
Pt : [360, 270] Size :2.25

こんなかんじです.

cv::Algorithmを継承したクラスは基本的にCV_INIT_ALGORITHMを使う必要があるみたいですね.
(´・ω・`)知らんがな.
これに気づくまでに結構時間がかかりました.

本当はcv::FeatureDetector::create("Dense")みたいなことがしたいのですが,これはOpenCVを再コンパイルする必要があると思うので,今回はとりあえず使えればOKというスタンスで実装しました.

何にせよ,これを使えばOpenCV 3.0 alphaでもDense Samplingが出来るようになりました!めでたしめでたし.

なお,この実装はOpenCVのgitリポジトリがHEADだと動かないみたいなので,
3.0.0 alphaにチェックアウトしてから使ってください.

今回のまとめ

  • OpenCV 3.0 alpha をメインプログラムに使うのはやめよう
  • Dense Samplingしたい人はOpenCV 2.4.9を使おう

今回参考にしたWebページ


次回は,FisherVectorか,Sparse Codingについて記事を書きたいと思います.