Visual Studio (VC++) 応用編 （その7）
− USBカメラを用いた空間フィルタ法による速度計測のシミュレーション −
信州大学工学部 井澤裕司
(H18.3.18)
課 題 7
この課題では，「 空間フィルタ法 」を用いて「 USBカメラ 」の一次元的な「 動き 」を計測する「 シミュレータ 」を制作します．
このプログラムは，昨年の10月4日に配布した本研究会資料「 空間フィルタ法による速度計測 」の中で紹介したものです．
詳細な説明は，重複するので省略しますが，「 カメラ 」もしくは「 対象物 」を「 一定の速度 」で動かすことが難しいので，「 USBカメラ 」で撮影した1枚の画像を，「 仮想的に動かす 」ことにより，速度計測を行います．
なお，「 空間フィルタ出力 」の「 直流成分 」を除去するため，下の図にしめすような「 位相が180度 」シフトした「 2種類のフィルタ 」の「 差分出力 」を解析する構成となっています．
それでは，具体的なプログラミングの手順について，説明しましょう．
２．主な手順
「 USBカメラ 」 の画像をキャプチャして 「 画面上に表示 」 する手法については，十分習熟されたと思います．
本編でも， 「 詳細な手順 」 は省略し，フレーム間差分から動き量を計測する手法のプログラムの内容を重点的に説明します．
主な 「 手順 」 は，以下の通りです．
(1) リソースの 「 メニュー 」 の追加「 Visual C++ 」 の操作法について， 不明な点が残っているのであれば，再度 「 応用編の前半 」 に戻って復習して下さい．
(2) リソースの 「 ビットマップ 」 の追加
(3) ソースコード 「 Cspace_filterView() 」の追加・修正
(4) 「 ライブラリ (vfw32.lib) 」 のリンク設定
(5) 「 ビルド 」 の実行とデバッグ
(6) 完成した 「 プロジェクト 」 の実行
３．メニューの追加
これまでのように，新しいプロジェクトとして 「 space_filter 」 を生成し， 以下の表に示すように，「 新規メニュー 」 を追加して下さい．
次に，「 イベントハンドラ 」 を用いて，「 space_filterView() 」 の中に関数を追加します
メニュー サブメニュー ID 関数名 初期設定 解像度 ID_INIT OnInit() 画像 キャプチャ ID_CAPTURE OnCapture 表示 ID_DISP OnDisp シミュレーション 空間フィルタ ID_FILTER OnFilter DFT ID_DFT OnDft
４．ビットマップの追加
次の表のように，新しく 「 ビットマップ 」 を作成します．
ID 位置 (Height) 位置 (Width)
Colors IDB_BITMAP1 480 640 Trueカラー
５．ソースコードの追加・修正
Viewクラス 「 Cspace_filterView() 」 の最初に，ヘッダーファイルをインクルードし， 画像表示を行うための 「 変数群 」 を追加します．
‥（略）‥
#include "space_filterView.h"
#include "vfw.h" // 追加したヘッダーファイル
#include "winbase.h" // 追加したヘッダーファイル
#include "math.h" // 追加したヘッダーファイル
#define X_SIZE 640 // 水平画素数
#define Y_SIZE 480 // 垂直画素数
#define COLOR 4 // 色（RGB(A) 4 Byte）
#define PI 3.1415927 // 円周率
#define N_FILT 11 // 空間フィルタの数
#define N_LEN 256 // 空間フィルタの移動量
#define F_GAIN 10 // フィルタ描画の係数
#define P_GAIN 100 // パワースペクトル描画の係数
#define N_LINE 200 // 計測するライン
#define DISP_LIMIT 1100 // 描画の範囲
struct complex { // 複素数の構造体を定義
double re; // 実数部
double im; // 虚数部
};
static unsigned char org_img[Y_SIZE][X_SIZE][COLOR]; // キャプチャ画像
static unsigned char filt_img[Y_SIZE][X_SIZE][COLOR]; // フィルタ後の画像
static int filt_mask[[Y_SIZE][X_SIZE]; // フィルタのマスク
static int filt_out[N_FILT][X_SIZE]; // フィルタ出力
comlex dft_mtrx[N_LEN][N_LEN]; // DFTマトリクス
int h_band; // フィルタ形状（縦の間隔:40）
double power_spec[N_FILT][N_LEN]; // パワースペクトル
double power_tmp[N_FILT][DISP_LIMIT]; // パワースペクトル（補正）
double power_sum[DISP_LIMIT]; // パワースペクトル（総和）
CDC m_memDC; // メモリDC
HBITMAP m_hBitmap = NULL; // DIB Section
BITMAPINFOHEADER* m_lpImage = NULL;
static BYTE* m_lpData = NULL; // RGBデータ
HWND m_hCapWnd = NULL;
inｔ size; // データサイズ
BOOL dispFlag = FALSE;
BITMAPINFOHEADER* GetDIBHeader(){ // DIB(Device Independent BMP)のヘッダ
return m_lpImage;
}
BOOL CreateDIBHeader(int size){ // DIB(Device Independent BMP)のヘッダ作成
m_lpImage = (BITMAPINFOHEADER*)malloc(size);
return (m_lpImage != NULL);
}
以下にその編集画面を示します．
5-[b]
コントラクタ「 Cspace_filterView::Cspace_filterView() 」 の内部に，変換係数等の値を初期設定する前処理用コードを追加します．
Cspace_filterView::Cspace_filterView()
{
int ix, iy;
int jy;
int kx;
// filt_mask の設定
// 背景⇒1
for(iy = 0; iy < Y_SIZE; iy++){
for(ix = 0; ix < X_SIZE; ix++){
filt_mask[iy][ix] = 1;
}
}
h_band = Y_SIZE / (N_FILT + 1);
// 空間フィルタのマスキング部⇒0
for(iy = 0; iy < N_FILT; iy++){
for(jy = 0; jy < h_band - 2; jy++){
for(ix = -4 * (iy + 1) - 3; ix <= 4 * (iy + 1) + 2; ix++){
filt_mask[iy * h_band + h_band / 2 + jy][X_SIZE / 2 + ix] = 0;
}
}
// 空間フィルタの透過部(上下)⇒2
for(jy = 0; jy < h_band / 2 - 5; jy++){
for(ix = -3; ix <= 3; ix = ix + 2){
for(kx = 0; kx < iy + 1; kx++){
filt_mask[iy * h_band + h_band / 2 + 2 + jy][X_SIZE / 2 + ix * (iy + 1) + kx] = 2;
}
}
}
for(jy = 0; jy < h_band / 2 - 5; jy++){
for(ix = -4; ix <= 2; ix = ix + 2){
for(kx = 0; kx < iy + 1; kx++){
filt_mask[iy * h_band + h_band + jy][X_SIZE / 2 + ix * (iy + 1) + kx] = 2;
}
}
}
}
// DFT（離散フーリエ変換）マトリクスの初期設定
for(iy = 0; iy < N_LEN; iy++){
for(ix = 0; ix < N_LEN; ix++){
dft_mtrx[iy][ix].re = cos(2.0 * PI * ix * iy / N_LEN) / sqrt((double)N_LEN);
dft_mtrx[iy][ix].im = sin(2.0 * PI * ix * iy / N_LEN) / sqrt((double)N_LEN);
}
}
}
以下にその編集画面を示します．
5-[c]
デストラクタ「 Cspace_filterView::~Cspace_filterView() 」 の内部に，後処理用（リソースの解放）のコードを追加します．
Cspace_filterView::~Cspace_filterView()
{
capDriverDisconnect( m_hCapWnd );
::DestroyWindow( m_hCapWnd );
m_lpData = NULL;
free( m_lpImage );
m_lpImage = NULL;
DeleteDC( m_memDC );
DeleteObject( m_hBitmap );
以下にその編集画面を示します．
5-[d]
デストラクタ「 Cspace_filterView::~Cspace_filterView() 」 の後に，画像キャプチャ終了時に起動するコードを追加します．
static LRESULT FAR PASCAL CaptureImage(HWND hWnd，LPVIDEOHDR lpVHdr)
{
memcpy(m_lpData， (LPSTR)lpVHdr->lpData， m_lpImage->biSizeImage);
dispFlag = TRUE;
return (LRESULT)TRUE;
}
以下にその編集画面を示します．
5-[e]
「 Cspace_filterView::PreCreateWindow(CREATESTRUCT& cs) 」 の内部に，描画準備用のコードを追加します．
BOOL Cimage_fileView::PreCreateWindow(CREARESTRCT& cs)
{
int i;
dispFlag = FALSE;
m_hCapWnd = capCreateCaptureWindowA("space_filter",
WS_OVERLAPPEDWINDOW, 0， 0， X_SIZE， Y_SIZE， NULL， NULL);
capDriverConnect(m_hCapWnd， 0);
capSetCallbackOnFrame( m_hCapWnd， CaptureImage);
size = capGetVideoFormatSize(m_hCapWnd);
CreateDIBHeader(size);
m_lpImage = GetDIBHeader();
capGetVideoFormat(m_hCapWnd， m_lpImage， size);
// DIBsectionの確保
m_hBitmap = CreateDIBSection( m_memDC.GetSafeHdc()，
(BITMAPINFO*)m_lpImage，
DIB_RGB_COLORS， (void**)&m_lpData， NULL，0 );
for( i = 0; i < 5; i++ ){
capGrabFrame( m_hCapWnd ); // 立ち上げ時のキャプチャ命令無応答の一時的対策
}
return CView::PreCreateWindow(cs);
}
以下にその編集画面を示します．
5-[f]
メニューの「解像度」に対応する「 OnInit() 」 の内部に，以下に示す描画用のコードを追加します．
void Cspace_filterView::OnInit()
{
capDlgVideoFormat(m_hCapWnd);
size = capGetVideoFormatSize(m_hCapWnd);
CreateDIBHeader(size);
m_lpImage = GetDIBHeader();
capGetVideoFormat(m_hCapWnd， m_lpImage， size);
m_hBitmap = CreateDIBSection( m_memDC.GetSafeHdc()，
(BITMAPINFO*)m_lpImage,
DIB_RGB_COLORS，
(void**)&m_lpData， NULL， 0 );
}
以下にその編集画面を示します．
5-[g]
自動生成された「 OnCapture() 」 の内部に，以下に示す描画用のコードを追加します．
CClientDC dc(this);
CDC pM;
CBitmap pB， *pOld;
int ix， iy;
BYTE *m_lpData2;
unsigned char r, g, b;
unsigned char y; // 輝度信号
CPen bluePen(PS_SOLID, 1, RGB( 0, 0, 255));
CPen* pOldPen;
pM.CreateCompatibleDC( &dc );
pB.LoadBitmap( IDB_BITMAP1 );
// キャプチャ，配列の生成
dispFlag = FALSE;
capGrabFrame( m_hCapWnd ); // 1フレーム（1枚）を取得するコマンド
while( dispFlag == FALSE) { } // フレームの読み込みが終了するまで待つ
m_lpData2 = (BYTE*)m_lpData; // 画像メモリのポインタ
for( iy = 0; iy < Y_SIZE; iy++){
for( ix = 0; ix < X_SIZE; ix++){
b = *m_lpData2; m_lpData2++; // Blue （青）
g = *m_lpData2; m_lpData2++; // Green （緑）
r = *m_lpData2; m_lpData2++; // Red （赤）
y = (unsigned char)(0.299 * (double)r + 0.587 * (double)g + 0.114 * (double)b);
// Y（輝度）信号
org_img[Y_SIZE-iy-1][ix][0] = ｙ; // データは上下反転している
org_img[Y_SIZE-iy-1][ix][1] = y; // このため入れ替える
org_img[Y_SIZE-iy-1][ix][2] = y; // モノクローム画像
}
}
// 画面のクリア
dc.FillSolidRect(0, 0, 1100, 800, RGB(255, 255, 255));
// 画像の表示
SetBitmapBits( pB, X_SIZE * Y_SIZE * COLOR, &org_img);
pOld = pM.SelectObject(&pB);
dc.BitBlt( 50, 30, X_SIZE, Y_SIZE, &pM, 0, 0, SRCCOPY);
pM.SelectObject( pOld );
// 計測ラインの描画
pOldPen = dc.SelectObject(&bluePen);
dc.MoveTo(50 , 30 + N_LINE);
dc.LineTo(50 + X_SIZE, 30 + N_LINE);
以下にその編集画面を示します．
5-[h]
自動生成された「 OnDisp() 」 の内部に，以下に示す描画用のコードを追加します．
CClientDC dc(this);
CDC pM;
CBitmap pB， *pOld;
CPen bluePen(PS_SOLID, 1, RGB( 0, 0, 255));
CPen* pOldPen;
// 画面のクリア
dc.FillSolidRect(0, 0, 1100, 800, RGB(255, 255, 255));
// キャプチャせずに画像表示
pM.CreateCompatibleDC( &dc );
pB.LoadBitmap( IDB_BITMAP1 );
SetBitmapBits( pB, X_SIZE * Y_SIZE * COLOR, &org_img);
pOld = pM.SelectObject(&pB);
dc.BitBlt( 50, 30, X_SIZE, Y_SIZE, &pM, 0, 0, SRCCOPY);
pM.SelectObject( pOld );
// 計測ラインの描画
pOldPen = dc.SelectObject(&bluePen);
dc.MoveTo(50 , 30 + N_LINE);
dc.LineTo(50 + X_SIZE, 30 + N_LINE);
以下にその編集画面を示します．
自動生成された「 OnFilter() 」 の内部に，以下に示す描画用のコードを追加します．
CClientDC dc(this);
CDC pM;
CBitmap pB， *pOld;
int it;
int ix， iy;
int k;
unsigned char y; // 輝度信号
char buf[30]; // 文章表示用のバッファ
pM.CreateCompatibleDC( &dc );
pB.LoadBitmap( IDB_BITMAP1 );
// 空間フィルタ（2相）の処理
for( it = 0; it < N_LEN; it++){
for( iy = 0; iy < Y_SIZE; iy++){
// 選択したライン(N_LINE)の描画
for( ix = 0; ix < X_SIZE - it; ix++){
y = org_img[N_LINE][ix + it][0];
if(filt_mask[iy][ix] != 0){ // フィルタ透過部
filt_img[iy][ix][0] = y;
filt_img[iy][ix][1] = y;
filt_img[iy][ix][2] = y;
}else{ // マスキング部⇒黒
filt_img[iy][ix][0] = 0;
filt_img[iy][ix][1] = 0;
filt_img[iy][ix][2] = 0;
}
}
for(ix = X_SIZE - it; ix < X_SIZE; ix++){ // その他⇒白
filt_img[iy][ix][0] = 255;
filt_img[iy][ix][1] = 255;
filt_img[iy][ix][2] = 255;
}
}
// フィルタ処理
for(iy = 0; iy < N_FILT; iy++){
filt_out[iy][it] = 0;
for(k = 0; k < X_SIZE; k++){
if(filt_mask[(iy + 1) * h_band - 10] == 2){ // フィルタ上部(+)
filt_out[iy][it] = filt_out[iy][it] + org_img[N_LINE][k + it][0];
}
if(filt_mask[(iy + 1) * h_band + 10][k] == 2){ // フィルタ下部(-)
filt_out[iy][it] = filt_out[iy][it] - org_img[N_LINE][k + it][0];
}
}
}
// 画面のクリア
dc.FillSolidRect(0, 0, 1100, 800, RGB(255, 255, 255));
// 画像の表示
SetBitmapBits( pB, X_SIZE * Y_SIZE * COLOR, &filt_img);
pOld = pM.SelectObject(&pB);
dc.BitBlt( 50, 30, X_SIZE, Y_SIZE, &pM, 0, 0, SRCCOPY);
pM.SelectObject( pOld );
// フィルタ出力の描画
for( iy = 0; iy < N_FILT; iy++){
dc.MoveTo(700, 30 + (iy + 1) * h_band);
for( ix = 0; ix < it; ix++){
dc.LineTo(700 + ix, 30 + (iy + 1) * h_band - filt_out[iy][ix] / F_GAIN);
}
}
}
dc.TextOut(560, 20, (CString)"スリットピッチ");
dc.TextOut(780, 20, (CString)"スリット出力 波形");
// 座標軸
for( iy = 0; iy < N_FILT; iy++ ){
sprintf( buf, "%3d", iy + 1);
dc.TextOut( 590, 20 + (iy + 1) * h_band, (CString) buf);
dc.TextOut( 970, 20 + (iy + 1) * h_band, (CString) "t");
}
以下にその編集画面を示します．（画像サイズが大きいため，上下2つに分けています）．
5-[j]
「 OnDft() 」 で使用する以下の関数を追加します．
// 複素数の積
struct complex c_mult(struct complex a, struct complex b)
{
struct complex c;
c.re = a.re * b.re - a.im * b.im;
c.im = a.re * b.im + a.im * b.re;
return c;
}
// 複素数の和
struct complex c_add(struct complex a, struct complex b)
{
struct complex c;
c.re = a.re + b.re;
c.im = a.im + b.im;
return c;
}
// 複素数の絶対値の2乗
double c_sqr(struct complex a)
{
return a.re * a.re + a.im * a.im;
}
以下にその編集画面を示します．
自動生成された「 OnDft() 」 の内部に，以下に示す描画用のコードを追加します．
CClientDC dc(this);
int it;
int i， j, k;
complex sig[N_LEN];
complex tmp;
char buf[30]; // 文章表示用のバッファ
CPen redPen(PS_SOLID, 1, RGB(255, 0, 0));
CPen greenPen(PS_SOLID, 1, RGB( 0, 128, 0));
CPen magentaPen(PS_SOLID, 1, RGB(255, 0, 255));
CPen* pOldPen;
// 画面のクリア
dc.FillSolidRect(0, 0, 1000, 800, RGB(255, 255, 255));
// フィルタ出力の描画
for( it = 0; it < N_LEN; it++){
if ( it != 0 ){
for(k = 0; k < N_FILT; k++){
dc.MoveTo(50 + it , 70 + h_band * k - filt_out[k][it - 1] / F_GAIN);
dc.LineTo(50 + it + 1, 70 + h_band * k - filt_out[k][it] / F_GAIN);
}
}
}
dc.SetTextColor(RGB( 0, 0, 0));
dc.TextOut(100, 10, (CString)"スリット光 波形（差分）");
for(i = 0; i < N_FILT; i++){
dc.TextOut(320, 60 + h_band * i, (CString)"t");
}
dc.TextOut(340, 10, (CString)"スリットピッチ");
for(i = 0; i < N_FILT; i++){
sprintf(buf, "%2d", i + 1);
dc.TextOut(390, 60 + h_band * i, (CString)buf);
}
// DFT(離散フーリエ変換)
for(k = 0; k < N_FILT; k++){
for(i = 0; i < N_LEN; i++){
sig[i].re = (double)filt_out[k][i] / F_GAIN;
sig[i].im = 0.0;
}
for(i = 0; i < N_LEN; i++){
tmp.re = 0.0;
tmp.im = 0.0;
for(j = 0; j < N_LEN; j++){
tmp = c_add(tmp, c_mult(sig[j], dft_mtrx[i][j]));
}
power_spec[k][i] = c_sqr(tmp);
}
}
// パワースペクトルの描画
pOldPen = dc.SelectObject(&magentaPen);
for(k = 0; k < N_FILT; k++){
for(i = 1; i < N_LEN; i++){
dc.MoveTo(450 + i * 4 , 70 + h_band * k - (int)(power_spec[k][i - 1] / P_GAIN));
dc.LineTo(450 + (i + 1) * 4, 70 + h_band * k - (int)(power_spec[k][i] / P_GAIN));
}
}
dc.SetTextColor(RGB( 0, 0, 0));
dc.TextOut(650, 10, (CString)"パワースペクトル");
for(i = 0; i < N_FILT; i++){
dc.TextOut(950, 50 + h_band * i, (CString)"f");
}
// 正規化パワースペクトルの描画
for(k = 0; k < N_FILT; k++){
for(i = 0; i < DISP_LIMIT; i++){
power_tmp[k][i] = 0.0);
}
}
for(k = 0; k < N_FILT; k++)
for(i = 0; i < DISP_LIMIT; i++){
for(j = 0; j < k; j++){
if(i * k + j < DISP_LIMIT)
power_tmp[k][i * k + j] = power_spec[k][i];
}
}
}
pOldPen = dc.SelectObject(&greenPen);
for(k = 0; k < N_FILT; k++){
for(i = 1; i < DISP_LIMIT; i++){
if(i + 450 < DISP_LIMIT){
dc.MoveTo(450 + i , 520 - (int)(power_tmp[k][i - 1] / P_GAIN));
dc.LineTo(450 + i + 1, 520 - (int)(power_tmp[k][i] / P_GAIN));
}
}
}
dc.SetTextColor(RGB( 0, 0, 0));
dc.TextOut(650, 490, (CString)"推定速度（ピッチ×周波数）");
// 正規化パワースペクトルの重ね合せの描画
for(i = 0; i < DISP_LIMIT; i++){
power_sum[i] = 0.0;
for(k = 0; k < N_FILT; k++){
power_sum[i] = power_sum[i] + power_tmp[k][i];
}
}
pOldPen = dc.SelectObject(&redPen);
for(i = 1; i < DISP_LIMIT; i++){
if(i + 450 < DISP_LIMIT){
dc.MoveTo(450 + i , 620 - (int)(power_sum[i - 1] / P_GAIN));
dc.LineTo(450 + i + 1, 620 - (int)(power_sum[i] / P_GAIN));
}
}
dc.SetTextColor(RGB( 0, 0, 0));
dc.TextOut(650, 590, (CString)"推定速度（ピッチ×周波数）: 総和");
以下にその編集画面を示します．
画像サイズが大きいため，上下2つに分けています．
６．ライブラリ［vfw32.lib] のリンク設定
メニュー「 表示 」の「 プロパティマネージャ 」をクリックして，「 構成プロパティ 」を開き，「 リンカ 」，「 入力 」を選択します．
右側の「 追加の依存ファイル 」の欄に「 vfw32.lib 」と入力します．
７．ビルドとデバッグ
「 ビルド 」を実行して，エラーメッセージがなくなるまで，「 デバッグ 」を繰り返します．
８．プロジェクトの実行
メニュー「 画像 」の「 キャプチャ 」の実行例を以下に示します．
「 対象物 」には「 広帯域の信号成分 」が含まれている必要があるため，「 モルタルの壁 」を撮影しています．
なお，「 青い線 」は「 解析するライン 」を示しています．
メニュー「 シミュレーション 」の「 空間フィルタ 」の実行例を以下に示します．
なお，図に示すように，「 スリット 」は「 差動構成 」とし，ピッチが異なる11種類のフィルタを用いています．
「 差動構成 」により，「 フィルタ出力 」の直流成分が除去されていることがわかります．
メニュー「 シミュレーション 」の「 DFT 」の実行例を以下に示します．
図の右側に，フィルタ出力の「 パワースペクトル 」が表示されています．
それらの「 ピーク 」が，「 双曲線上 」に分布していることが分かります．
「 フィルタのピッチ 」と「 ピークとなる周波数 」の積が一定値になることは明らかです．
「 赤いグラフ は，これらを「 正規化 」して「 重ね合わせた 」ものです．
９．まとめ
本章では，「 空間フィルタ法 」を用いて「 対象物の動き量 」を計測するプログラムを作成しました．
このプログラムを用いて，「 様々な対象物 」を撮影して，その「 計測結果 」を評価して下さい．
このシミュレーションにより，対象物の明るさが，「 細かくランダム 」，すなわち「 白色雑音 」のように変化する必要があることが実感できると思います．
以上で，「 USBカメラ 」と「 VC++ 」を用いた「 速度計測手法 」の解説を終了します．
この他にも，「 色信号 」や「 ヒストグラム 」など，紹介できなかった重要な手法が数多くありますが，それらについては，「 画像処理 」や「 信号処理 」などの書籍で紹介されていますので，挑戦してみて下さい．
このコンテンツが，それらの力になることを願っています．
