GPGPUによる高速画像処理

GPGPU による高速画像処理 ~ リアルタイム画像処理への挑戦 ~ 名古屋大学大学院情報科学研究科 出口大輔

リアルタイム画像処理 2

3 発表の流れ GPGPU を始める前に GPGPU の基礎知識 CUDA の使い方 CUDA を使う前に プログラミングの予備知識 CUDA を使って Hello World GPGPU にチャレンジ 行列積の計算 テンプレートマッチング ガウシアンフィルタ SIFT 特徴量の計算

~GPGPU って何?~ 4

5 GPGPU って? GPGPU は何の略? General-Purpose computation on GPUs GPUを汎用計算に利用しようという試み 現在は GPUコンピューティング とも呼ばれる どうして GPGPU が注目されているのか?

6 CPU と GPU の性能比較 1600 [GFLOP/s] NVIDIA GPU 1400 Geforece GTX480 ピーク性能 1200 1000 800 600 400 200 Geforece 8800GTX Geforece GTX280 Core2 Duo 3.0GHz Quad Core Xeon 3.2GHz Intel CPU 0 1900 1900 1900 1900 1900 1900 1900 1900 1900 2003 2004 2005 2006 2007 2008 2009 2010 NVIDIA CUDA Programming Guide 4.0より引用

7 GPGPU って? GPGPU は何の略? General-Purpose computation on GPUs GPUを汎用計算に利用しようという試み 現在は GPUコンピューティング とも呼ばれる どうして GPGPU が注目されているのか? GPUの計算性能が飛躍的に向上 最新のGPUは 1.5 TFLOPS 以上の演算性能 (CPUの約 10 倍 ) GPUはCPUと比較して並列計算に優れている GeForce GTX580 では 512 コアによる並列計算が可能 手頃な価格で入手可能 GeForce GTX580 は約 5 万円で購入可能

8 GPGPU の活用例 動画像処理 CyberLink PowerDirector 7 ビデオエフェクトのレンダリングを高速化 TMPGEnc 4.0 XPress フィルター処理 デコード処理の高速化 画像処理 Adobe Photoshop CS4 / CS5,Adobe Premire CS5 各種フィルタ処理の高速化 OpenCV 各種画像処理の高速化 数値計算 MATLAB FFT の高速化

9 GPU の歴史 : 1999 年以前 1970 年 ~ 1990 年 初期開発の時代 ソフトウェアによるグラフィックス処理 プログラム可能な GPU に関する初期研究 Ikonas System [1], Pixel Machine [2] [1] J. N. England, A system for interactive modeling of physical curved surface objects, Proceedings of SIGGRAPH 78, pp.336-340. 1978 [2] M. Potmesil and E. M. Hoffert, The Pixel Machine: A Parallel Image Computer, Proceedings of SIGGRAPH89, pp.69-78, 1989 1990 年 ~ 1999 年 GPU 技術の黎明期 3D グラフィックス アクセラレータの開発 グラフィックス向けライブラリの開発 OpenGL(1992 年 ~), DirectX(1995 年 ~)

10 GPU の歴史 : 1999 年 ~ GPU の誕生 NVIDIA 社の GeForce 256 の登場 ハードウェア T&L をサポート CPU 負荷を大幅に削減 グラフィックスパイプライン ハードウェア固定の処理 自由表現な表現は不可 頂点処理 ( ハードウェア固定 ) ジオメトリ処理 ( クリッピング等 ) ラスタライズ処理 ピクセル処理 ( ハードウェア固定 ) Transform & Lighting 画面出力

11 GPU の歴史 : 2003 年 ~ プログラマブルシェーダの登場 Vertex Shader 頂点座標の変換処理 Pixel Shader 画素の輝度計算処理 シェーダ言語の進化 アセンブラから高級言語へ Cg, HLSL, GLSL 柔軟な映像表現が可能に グラフィックス以外への応用 GPGPU への関心が高まる Vertex Shader ( プログラマブル ) ジオメトリ処理 ( クリッピング等 ) ラスタライズ処理 Pixel Shader ( プログラマブル ) 画面出力

12 GPU の歴史 : 2007 年 ~ GPGPU の開発環境の整備 CUDA (NVIDIA) ATI Stream (AMD) OpenCL DirectCompute GPGPU 時代の到来 物理シミュレーション 数値計算 信号解析 画像処理 認識 統合型シェーダ Vertex Shader ( プログラマブル ) Geometry Shader ( プログラマブル ) ラスタライズ処理 Pixel Shader ( プログラマブル ) 画面出力

13 GPU の歴史 : 2009 年 ~ Fermi アーキテクチャの登場 ビデオカードによって サポート状況は異なる 汎用計算向けのアーキテクチャ L1/L2キャッシュの搭載 複数カーネルの同時実行のサポート 倍精度浮動小数点演算の高速化 ECCメモリのサポート アトミックなメモリ操作の高速化 C++ のフルサポート GPGPU 関連のライブラリ & ツールの充実 CUBLAS, CUFFT, Thrust, NPP Parallel Nsight, Visual Profiler

~CUDA を使う前に ~ 14

15 CUDA って何? CUDA(Compute Unified Device Architecture) 発音は クーダ もしくは キューダ NVIDIA 社が提供する GPU を利用するための統合開発環境 GeForce 8 以降のハードウェアで利用可能 グラフィックス処理 API の知識は不要 CPU でプログラムを実行する感覚 C/C++ を用いてプログラムの開発が可能 既存のアルゴリズムの移植も比較的容易

16 CUDA が利用可能な環境 CUDA 対応のグラフィックスカード GeForce GTX 580, GTX 480, GTX 260, 8800シリーズ, 他 Quadro Plex 2200 D2, FX 5800, FX 5600, 5000, 6000, 他 Tesla S2050, C2070, S1070, C1060, S870, D870, C870 TeslaはHPCに特化した製品であり映像出力を持たない OS の対応状況 (CUDA 4.0) Windows XP(32bit, 64bit) Windows Vista(32bit, 64bit) Windows Server 2008(32bit, 64bit) Linux(32bit, 64bit) Mac OS

17 CUDA を使うための準備 (Windows 編 ) CUDA ZONE から次をダウンロード NVIDIA Driver CUDA 対応のビデオドライバー CUDA Toolkit 4.0 コンパイラ (nvcc) CUBLASやCUFFT ライブラリ ドキュメント CUDA SDK 4.0 サンプル 2011 年 5 月 25 日にリリース http://www.nvidia.com/object/cuda_home.html

18 CUDA を動かしてみよう!! Volumetric Particle Shadows Image Denoising CUDA SDK 内のサンプルの実行結果

~ プログラミングの予備知識 ~ 19

20 CUDA 対応の GPU(G80, GT200) ストリーミング マルチプロセッサ (SM) 8 個のスカラープロセッサ (SP) 16KB の共有メモリ スレッド間の同期機構 マルチプロセッサ内でのみ可能 マルチプロセッサ間での同期には CPU 処理が必要 GeForce GTX280 の場合 30 基のマルチプロセッサを搭載 1GByte 以上のグローバルメモリ マルチプロセッサ #1 SP SP SP SP SP SP SP SP マルチプロセッサ #2 SP SP SP SP SP SP SP SP

21 CUDA 対応の GPU(Fermi) ストリーミング マルチプロセッサ (SM) 32 個のスカラープロセッサ (SP) 48KB の共有メモリ スレッド間の同期機構 マルチプロセッサ内でのみ可能 マルチプロセッサ間での同期には CPU 処理が必要 GeForce GTX480 の場合 15 基のマルチプロセッサを搭載 1GByte 以上のグローバルメモリ SM #1 SP SP SP SP SP SP SP SP SP SP SP SP SP SP SP SP SP SP SP SP SP SP SP SP SP SP SP SP SP SP SP SP SM #2 SP SP SP SP SP SP SP SP SP SP SP SP SP SP SP SP SP SP SP SP SP SP SP SP SP SP SP SP SP SP SP SP

22 CUDA のプログラミングモデル GPU の特徴 多数のスレッドが高い並列性をもって処理を実行 GPU のみでプログラムを実行できない CPUとの連携が不可欠 CUDA における計算の流れ GPUを並列演算可能なデバイスとして扱う 複数のスレッドを同時に実行できる外部 CPU 階層的にスレッドを管理 スレッドのまとまり = ブロック ブロックのまとまり = グリッド 問題を分割して計算する際に便利

23 CUDA におけるスレッド管理 スレッドを 3 次元的に配置 各スレッドの ID を X, Y, Z の3 要素で表現 スレッドのまとまりをブロックとして管理 グリッド ブロックブロック (0,0,0) (0,0,1) ブロックブロック (1,0,0) (1,0,1) スレッドスレッド (0,0,0) (0,0,1) スレッドスレッド (1,0,0) (0,1,0) (0,1,1) スレッドスレッド (2,0,0) (0,2,0) (0,2,1) スレッドスレッド (0,0,0) (0,0,1) スレッドスレッド (1,0,0) (0,1,0) (0,1,1) スレッドスレッド (2,0,0) (0,2,0) (0,2,1) スレッドスレッド (0,1,0) (1,0,0) (1,0,1) スレッドスレッド (1,1,0) (1,1,1) スレッドスレッド (2,1,0) (1,2,0) (1,2,1) スレッドスレッド (0,1,0) (1,0,0) (1,0,1) スレッドスレッド (1,1,0) (1,1,1) スレッドスレッド (2,1,0) (1,2,0) (1,2,1) ブロックブロック (0,1,0) (0,1,1) ブロックブロック (1,1,0) (1,1,1)

24 階層的スレッドの利用方法 一般的な画像処理で利用する場合 ブロック内のスレッド数を決定 16 16 = 256 スレッド 画像内にブロックを配置 ブロック内のスレッドが各画素を処理 計算範囲の求め方 ブロック ID スレッド ID ブロック内スレッド数 T 00 T 10 T 01 T 0M T N0 T NM ブロック #1 T 00 T 10 T N0 T 01 T 0M T NM

25 CUDA における計算の流れ グリッド 1 処理 1 ブロック (0,0,0) スレッド (0,0,0) スレッド (1,0,0) スレッド (2,0,0) ブロック (1,0,0) スレッド (0,0,0) スレッド (1,0,0) スレッド (2,0,0) スレッド (0,1,0) スレッド (1,1,0) スレッド (2,1,0) スレッド (0,1,0) スレッド (1,1,0) スレッド (2,1,0) ブロック (0,1,0) ブロック (1,1,0) 処理 2 グリッド 2 ブロック (0,0,0) ブロック (1,0,0)

26 CUDA のメモリモデル ブロック (0,0,0) 共有メモリ ブロック (1,0,0) 共有メモリ ローカルメモリ ローカルメモリ ローカルメモリ ローカルメモリ レジスタ レジスタ レジスタ レジスタ スレッド (0,0,0) スレッド (1,0,0) スレッド (0,0,0) スレッド (1,0,0) コンスタントメモリ テクスチャメモリ グローバルメモリ

27 CUDA で利用可能なメモリ グローバルメモリ 大量のメモリを利用可能 低速なメモリアクセス (400~600 クロック必要 ) Fermi アーキテクチャではキャッシュ機構 (L1/L2) を搭載 共有メモリ レジスタと同じ速度でアクセス可能 マルチプロセッサ 1 基あたり 16KB(Fermi は最大 48KB) Fermi の場合は一部を L1 キャッシュとして利用可能 テクスチャメモリ キャッシュ機構による高速なアクセスが可能 ハードウェア線形補間や正規化テクスチャ座標が利用可能 コンスタントメモリ キャッシュ機構による高速なアクセスが可能 マルチプロセッサ 1 基あたり 64KB

28 各 GPU で使用できる機能の違い 使用できる機能を Compute Capability で区別 現在 1.0 ~ 2.1 の GPU が存在 Compute Capability による機能の違い 1.0: 初期リリース GeForce 8800GTX, Quadro FX 5600, 他 1.1: アトミックなメモリ操作のサポート GeForce 9800 GTX, Quadro FX 3700, 他 1.3: 倍精度浮動小数点演算のサポート GeForce GTX 280, Quadro FX 5800, 他 2.x : グローバルメモリのキャッシュをサポート GeForce GTX480,Quadro 5000, 他

29 Compute Capability による機能の違い Compute Capability 1.0 1.1 1.2 1.3 2.x 3 次元グリッド 倍精度浮動小数 32 ビット整数のアトミック演算 64 ビット整数のアトミック演算 単精度浮動小数のアトミック演算

30 CUDA における制限事項 Compute Capability 1.0 1.1 1.2 1.3 2.x 1 ブロックあたりのスレッド数 1 マルチプロセッサあたりのスレッド数 1 マルチプロセッサあたりのレジスタ数 コンスタントメモリ 512 1024 768 1024 1536 8192 16384 32768 64KB 共有メモリ 16KB 48KB 2D テクスチャ 2 16 2 15 2 16 2 16

~ CUDA を使って Hello World ~ 31

32 CUDA を使って Hello World!! Hello World!! を表示するサンプル (main.cu) global void hello( char *data ) { char *text = "Hello World!! n"; data[ threadidx.x ] = text[ threadidx.x ]; } int main( int argc, char *argv[] ) { char *ddata, hdata[ 15 ]; cudamalloc( ( void ** )&ddata, sizeof( char ) * 15 ); dim3 nthreads( 15, 1 ); dim3 nblocks( 1, 1 ); hello<<< nblocks, nthreads >>>( ddata ); cudamemcpy( hdata, ddata, sizeof( char ) * 15, cudamemcpydevicetohost ); printf( "%s", hdata ); cudafree( ddata ); } return( 0 );

33 CUDA を使って Hello World!! Hello World!! を表示するサンプル (main.cu) global void hello( char *data ) { char *text = "Hello World!! n"; data[ threadidx.x ] = text[ threadidx.x ]; } int main( int argc, char *argv[] ) { char *ddata, hdata[ 15 ]; global, threadidx CUDA で拡張された部分 cudamalloc( ( void ** )&ddata, sizeof( char ) * 15 ); dim3 nthreads( 15, 1 ); dim3 nblocks( 1, 1 ); hello<<< nblocks, nthreads >>>( ddata ); dim3, <<< >>> cudamemcpy( hdata, ddata, sizeof( char ) * 15, cudamemcpydevicetohost ); printf( "%s", hdata ); cudafree( ddata ); } return( 0 );

34 CUDA を使って Hello World!! Hello World!! を表示するサンプル (main.cu) global void hello( char *data ) { char *text = "Hello World!! n"; data[ threadidx.x ] = text[ threadidx.x ]; } int main( int argc, char *argv[] ) { char *ddata, hdata[ 15 ]; cudamalloc( ( void ** )&ddata, sizeof( char ) * 15 ); dim3 nthreads( 15, 1 ); dim3 nblocks( 1, 1 ); hello<<< nblocks, nthreads >>>( ddata ); cudamemcpy( hdata, ddata, sizeof( char ) * 15, cudamemcpydevicetohost ); printf( "%s", hdata ); cudafree( ddata ); } return( 0 );

35 CUDA の言語拡張 (1) 言語拡張により追加された修飾子 関数に対する修飾子 global device host CPU から呼び出され,GPU で実行される関数 GPU から呼び出され,GPU で実行される関数 CPU から呼び出され,CPU で実行される関数 変数に対する修飾子 constant shared device GPU 上のコンスタントメモリに存在する変数 GPU 上の共有メモリに存在する変数 GPU 上のメモリに存在する変数

36 CUDA を使って Hello World!! Hello World!! を表示するサンプル (main.cu) global void hello( char *data ) { char *text = "Hello World!! n"; data[ threadidx.x ] = text[ threadidx.x ]; } int main( int argc, char *argv[] ) { char *ddata, hdata[ 15 ]; cudamalloc( ( void ** )&ddata, sizeof( char ) * 15 ); dim3 nthreads( 15, 1 ); dim3 nblocks( 1, 1 ); hello<<< nblocks, nthreads >>>( ddata ); cudamemcpy( hdata, ddata, sizeof( char ) * 15, cudamemcpydevicetohost ); printf( "%s", hdata ); cudafree( ddata ); } return( 0 );

37 CUDA の言語拡張 (2) CUDA で利用可能な組み込み型 dim3 uchar2, int2, float2, uchar3, int3, float3, uchar4, int4, float4, 整数 x, y, z からなる 3 次元ベクトル ( スレッド数やブロック数の指定に利用 ) x, y からなる 2 次元ベクトル x, y, z からなる 3 次元ベクトル x, y, z, w からなる 4 次元ベクトル GPU 内のスレッドを識別する組み込み変数 griddim blockidx blockdim threadidx グリッドの次数スレッドが属するブロックのインデックススレッドが属するブロックの次数ブロック内のスレッドのインデックス

38 スレッドを識別する組み込み変数 グリッド blockidx ブロックブロック (0,0,0) (0,0,1) スレッドスレッドスレッドスレッドスレッド (0,0,1) スレッドスレッド (0,1,1) (0,2,1) スレッドスレッド (0,0,0) (0,0,0) (1,0,0) (0,1,0) (2,0,0) (0,2,0) スレッドスレッドスレッドスレッド (1,0,1) (1,1,1) スレッド (1,2,1) スレッドスレッドスレッドスレッド (0,1,0) (1,0,0) (1,1,0) (1,1,0) (2,1,0) (1,2,0) ブロックブロック (1,0,0) (1,0,1) blockdim = 3 2 2 スレッドスレッドスレッド (0,0,1) (0,0,0) (0,0,0) スレッドスレッド (1,0,1) スレッド (0,1,0) (1,0,0) スレッドスレッドスレッド (0,1,1) (1,0,0) (0,1,0) スレッド (1,1,1) スレッドスレッド (1,1,0) (1,1,0) スレッド (0,2,1) スレッドスレッド (2,0,0) (0,2,0) スレッド (1,2,1) スレッドスレッド (2,1,0) (1,2,0) threadidx ブロックブロック (0,1,0) (0,1,1) ブロックブロック (1,1,0) (1,1,1) griddim = 2 2 2

39 C/C++ プログラミングとの相違点 CPU と GPU で実行するコードを明確に区別 CPUで実行する場合は host ( 省略可 ) を付与 GPUで実行する場合は global を付与 CPU と GPU で処理可能なメモリ空間の違い CPU と GPU 間でメモリ転送が必要 スレッド数を指定して GPU 上の関数を呼び出し 関数名 <<< ブロック数, スレッド数 >>>( ); ブロック数とスレッド数は実行時に指定可能 問題サイズに合わせて調整可能

40 CUDA を使って Hello World!! Hello World!! を表示するサンプル (main.cu) global void hello( char *data ) { char *text = "Hello World!! n"; data[ threadidx.x ] = text[ threadidx.x ]; } int main( int argc, char *argv[] ) { char *ddata, hdata[ 15 ]; cudamalloc( ( void ** )&ddata, sizeof( char ) * 15 ); dim3 nthreads( 15, 1 ); dim3 nblocks( 1, 1 ); hello<<< nblocks, nthreads >>>( ddata ); GPU 上で実行される関数 GPU 上にメモリを確保 cudamemcpy( hdata, ddata, sizeof( char ) * 15, cudamemcpydevicetohost ); ブロック内に配置するスレッド数 15 1 = 15 総ブロック数 1 1 = 1 printf( "%s", hdata ); cudafree( ddata ); CPU から GPU へメモリ転送 } return( 0 ); GPU 上のメモリを解放

41 スレッドのレジスタ数を調査 Visual Studio 2008 コマンドプロンプトを起動 コマンドプロンプト上で main.cu をコンパイル C: Your Path>nvcc main.cu --ptxas-options=-v --compile main.cu ptxas info : Compiling entry function '_Z5helloPc' for 'sm_10' ptxas info : Used 2 registers, 8+16 bytes smem, 15 bytes cmem[0] レジスタ数 合計レジスタ数をチェック 2 15 = 30 レジスタ数生成スレッド数合計

~Thrust ライブラリの使い方 ~ 42

43 Thrust ライブラリって? CUDA で利用できるテンプレートライブラリ CUDAとOpenMPに対応したC++ STLの並列版 URL: http://thrust.googlecode.com/ Thrust ライブラリの特徴 コンテナ CPUとGPUのデータをコンテナとして管理 アルゴリズム コンテナに対してアルゴリズムを適用可能 直感的なインターフェース CUDAの複雑なAPIに関する知識は不要 メモリの確保 / 解放 / 転送がとても簡単

44 Thrust ライブラリの使い方 GPU 上で 1,2,3,,100 の総和を計算 #include <thrust/host_vector.h> #include <thrust/device_vector.h> #include <thrust/sequence.h> #include <thrust/reduce.h> int main( int argc, char *argv[] ) { // CPU 上のメモリを確保 thrust::host_vector< int > hvec( 100 ); // データを初期化 1, 2, 3,... thrust::sequence( hvec.begin( ), hvec.end( ), 1 ); // GPU 上のメモリを確保 thrust::device_vector< int > dvec( 100 ); // CPU -> GPU のメモリ転送 dvec = hvec; int val = thrust::reduce( dvec.begin( ), dvec.end( ), 0, thrust::plus< int >( ) ); printf( "%d n", val ); } return( 0 );

45 コンテナの使い方 (1) メモリの確保 thrust::host_vector< int > hvec( 20 ); CPU 上のメモリに要素数 20の int 配列を確保 thrust::device_vector< float > dvec( 100 ); GPU 上のメモリに要素数 100の float 配列を確保 メモリの解放 コンテナオブジェクトの消滅時に自動解放 CPUの場合 :free( ) GPUの場合 :cudafree( ) 明示的なメモリ解放 :dvec.clear( )

46 コンテナの使い方 (2) メモリの転送 代入操作でメモリ転送が可能 thrust::host_vector< float > hvec( 20 ); thrust::device_vector< float > dvec( 20 ); dvec = hvec; コンテナの各要素へのアクセス 配列と同様にアクセス可能 hvec[ 5 ] = 100.0f; dvec[ 1 ] = 23.4f; 内部で cudamemcpy が呼ばれるので注意

47 イテレータ (1) コンテナの要素を指すポインタ STLのイテレータと同じように使用可能 thrust::device_vector< int > dvec( 4 ); thrust::device_vector< int >::iterator ite = dvec.begin( ); dvec.begin( ) dvec.end( ) dvec 0 1 2 3 dvec.begin( ) + 3

48 イテレータ (2) コンテナの各要素へのアクセス thrust::device_vector< int > dvec( 4 ); thrust::device_vector< int >::iterator ite = dvec.begin( ); *ite = 10; dvec[ 0 ] = 10 と等価 ++ite; *ite = 25; dvec[ 1 ] = 25 と等価

49 アルゴリズムと使用方法 CPU と GPU のデータに対して同じ様に適用可能 GPU の場合は GPU を使用して並列処理される アルゴリズムの使用例 コンテナのデータすべてに1を代入する thrust::fill( dvec.begin( ), dvec.end( ), 1 ); 1 1 1 1 1 1 1 1 コンテナのデータに 1, 2, 3, の値で初期化する thrust::sequence( dvec.begin( ), dvec.end( ) ); 1 2 3 4 5 6 7 N

50 Parallel Reduction(1) コンテナ内の要素を 1 つの値に集約する処理 GPUのスレッドが並列に処理を実行 例 : 総和計算の場合 GPU 上のデータ 1 2 3 4 5 6 7 8 3 7 11 15 10 26 36 36 CUDA スレッドの同期 CUDA スレッドの同期

51 Parallel Reduction(2) Thrust ライブラリを用いた集約処理 総和の計算 int val = thrust::reduce( dvec.begin( ), dvec.end( ) ); 1 2 3 4 5 6 7 10 55 最大値の計算 int val = thrust::reduce( dvec.begin( ), dvec.end( ), 1 10 3 8 5 9 7 2 6 10 0, thrust::maximum< int >( ) );

52 CUDA と Thrust を組み合わせる Hello World!! を表示するサンプル #include <thrust/host_vector.h> #include <thrust/device_vector.h> global void hello( char *data ) { char *text = "Hello World!! n"; data[ threadidx.x ] = text[ threadidx.x ]; } int main( int argc, char *argv[] ) { thrust::host_vector< char > hvec( 15 ); thrust::device_vector< char > dvec( 15 ); dim3 nthreads( 15, 1 ); dim3 nblocks( 1, 1 ); hello<<< nblocks, nthreads >>>( thrust::raw_pointer_cast( &dvec[ 0 ] ) ); hvec = dvec; printf( "%s", &hvec[ 0 ] ); return( 0 ); CPUからGPUへメモリ転送 } GPU 上のメモリを指すポインタを取得

53

54 GPGPU にチャレンジ 行列積の計算 CUDA における基本的な実装方法 共有メモリ テンプレートマッチング 2 次元テクスチャメモリ ハードウェア線形補間 正規化テクスチャ座標 Parallel Reduction アルゴリズム ガウシアンフィルタ 1 次元テクスチャメモリ コンスタントメモリ SIFT 特徴量の計算 1 2 次元テクスチャメモリ コンスタントメモリ

~ 行列積の計算 ~ 55

56 GPU による行列積の計算 行列積計算の特徴 各要素はすべて独立に計算可能 GPUによる並列計算が可能 同じメモリ領域へ頻繁にアクセス メモリアクセス速度がボトルネック GPU 上での実装方法 1. 行列積の計算式に忠実な実装 (GPU#1) 2. 共有メモリを利用した高速化 (GPU#2) ブロック単位 ( 部分行列 ) で行列積を計算共有メモリをキャッシュとして利用

57 行列積 C=A B の流れ 1. 行列 A, B, CのメモリをGPU 上に確保 2. 行列 A, BのデータをGPUに転送 3. 行列 Cの各要素 c mn を計算 N 次正方行列を仮定 ただし,N は16の倍数 c mn N k 1 mk b 4. 計算結果を CPU に転送 a kn m 行 a mk N 行列 A N n 列 行列 C b kn 行列 B c mn

58 処理の流れ CPU 側の処理 行列の初期化 GPU 側の処理 行列 C の各要素を計算 CPU GPU CPU GPU GPU GPU CPU 行列の初期化 メモリ確保 データ転送 行列積の計算 データ転送

59 行列積計算の CPU 処理 ( ソース ) int main( int argc, char *argv[] ) { int N = 512; // 行列 A, B, Cのサイズ float *ha, *hb, *hc; // CPU(host) 側で利用するメモリへのポインタ float *da, *db, *dc; // GPU(device) 側で利用するメモリへのポインタ /* CPU 側のメモリを確保 */ /* GPU 側のメモリを確保 */ /* CPU 側のメモリをGPU 側へ転送 */ /* 実行するGPUのスレッド数, ブロック数を設定 */ /* GPUのカーネルを実行し,C=A B の結果を dc に格納 */ /* GPUの計算結果をCPU 側へ転送 */ /* CPUとGPUそれぞれのメモリを解放 */ } return( 0 );

60 行列積計算の CPU 処理 (1) 行列積の計算で使用する変数 int N = 512; float *ha, *hb, *hc; float *da, *db, *dc; // 行列 A, B, Cのサイズ // CPU(host) 側で利用するメモリへのポインタ // GPU(device) 側で利用するメモリへのポインタ CPU と GPU でメモリを確保 CPU:malloc, new, cudamallochost でメモリを確保 GPU:cudaMalloc によりグローバルメモリを確保 /* CPU 側のメモリを確保 */ ha = ( float * )malloc( N * N * sizeof( float ) ); hb = ( float * )malloc( N * N * sizeof( float ) ); hc = ( float * )malloc( N * N * sizeof( float ) ); /* GPU 側のメモリを確保 */ cudamalloc( ( void ** )&da, N * N * sizeof( float ) ); cudamalloc( ( void ** )&db, N * N * sizeof( float ) ); cudamalloc( ( void ** )&dc, N * N * sizeof( float ) );

61 行列積計算の CPU 処理 (2) CPU と GPU 間のメモリ転送 CPU GPU 転送元 転送方向 cudamemcpy( dst, src, size, cudamemcpyhosttodevice ); 転送先 転送バイト数 GPU CPU cudamemcpy( dst, src, size, cudamemcpydevicetohost ); この部分を変更 CPU と GPU で確保したメモリの解放 CPU:free, delete, cudafreehost でメモリを解放 GPU:cudaFree でグローバルメモリを解放

62 行列積計算の CPU 処理 (3) GPU 上で処理を実行 計算範囲の設定 N 次正方行列を仮定 ただし N は 16 の倍数 行列 Cをサイズ16 16のブロックに分割 ブロック内の各要素を各スレッドが計算 合計 M M 個のブロックを配置 各スレッドが c mn を計算 n = threadidx.x + blockidx.x blockdim.x; m = threadidx.y + blockidx.y blockdim.y; GPU 上で関数を実行 multiply<<< ブロック数, スレッド数 >>>( ); 16 16 (threadidx.y, threadidx.x) 行列 C (blockidx.y, blockidx.x) 行方向のブロック数 (=16)

63 行列積計算の CPU 処理 ( まとめ ) int main( int argc, char *argv[] ) { int N = 512; // 行列 A, B, Cのサイズ float *ha, *hb, *hc; // CPU(host) 側で利用するメモリへのポインタ float *da, *db, *dc; // GPU(device) 側で利用するメモリへのポインタ /* CPU 側のメモリを確保 */ /* GPU 側のメモリを確保 */ /* CPU 側のメモリを GPU 側へ転送 */ GPU 上で処理を実行 /* 実行する GPU のスレッド数, ブロック数を設定 */ dim3 nthreads( 16, 16 ); dim3 nblocks( N / nthreads.x, N / nthreads.y ); /* GPU のカーネルを実行し,C=A B の結果を dc に格納 */ multiply<<< nblocks, nthreads >>>( da, db, dc, N ); /* GPU の計算結果を CPU 側へ転送 */ /* CPU と GPU それぞれのメモリを解放 */ GPU 上のメモリを指定 } return( 0 );

64 実装方法 (GPU#1) n 列 行列積 C=A B の各要素 c mn を計算 N b kn c mn N k 1 a mk b kn m 行 a mk N 行列 B c mn 要素 c mn を GPU 上で並列計算 行列 A 行列 C global void multiply1( float *da, float *db, float *dc, int N ) { int n = threadidx.x + blockidx.x * blockdim.x; int m = threadidx.y + blockidx.y * blockdim.y; float sum = 0.0f; for( int k = 0 ; k < N ; k++ ) { sum += da[ m + k * N ] * db[ k + n * N ]; } } dc[ m + n * N ] = sum;

65 実装方法 (GPU#2) 各要素 c mn を部分行列 ( ブロック ) の積で計算 ブロックサイズが 16 16 の場合 c mn N 16 t 1 16( t 1) k 16t ブロック内で同じメモリを参照 mk 高速な共有メモリを利用 a グローバルメモリへのアクセスを削減 2 N 16 2 回 16 分の1 2 N 16 回 b kn 16 16 N a mk N 行列 A b kn 16 c mn 16 行列 B 行列 C

66 実装方法 (GPU#2) global void multiply2( float *da, float *db, float *dc, int N ) { int n = threadidx.x + blockidx.x * blockdim.x; 共有メモリの宣言 int m = threadidx.y + blockidx.y * blockdim.y; float sum = 0.0f; for( int k = 0 ; k < N ; k += 16 ) { shared float ta[ 16 ][ 16 ]; shared float tb[ 16 ][ 16 ]; 各スレッドが独立してメモリアクセス ta[ threadidx.y ][ threadidx.x ] = da[ m + ( k + threadidx.x ) * N ]; tb[ threadidx.x ][ threadidx.y ] = db[ ( k + threadidx.y ) + n * N ]; syncthreads( ); ブロック内のスレッドを同期 16 行列 B for( int t = 0 ; t < 16 ; t++ ) { sum += ta[ threadidx.y ][ t ] * tb[ threadidx.x ][ t ]; } N b kn 16 } } syncthreads( ); dc[ m + n * N ] = sum; ブロック内のスレッドを同期 各スレッドが独立して行列積を計算 16 16 N a mk 行列 A c mn 行列 C

67 計算時間の比較 (1) CPU, GPU#1, GPU#2 の計算時間を比較 使用計算機 CPU: Intel Core i7-x980(3.33ghz) OpenMPを使用して6スレッドで並列計算 GPU: NVIDIA Geforce GTX480 マルチプロセッサ数 : 15(SP 数 = 15 32 = 480) Memory:1.5 GB OS:Windows 7 SP1

16 32 48 64 80 96 112 128 144 160 176 192 208 224 240 272 288 304 320 336 352 368 400 416 432 448 464 480 496 528 544 560 592 608 624 656 672 688 720 736 752 768 784 68 800 計算時間の比較 (2) 500 450 400 350 300 [msec.] 計算時間 (N = 560 の場合 ) CPU#1: 91.7 msec. GPU#1: 18.3 msec. ( 5.0) GPU#2: 06.7 msec. ( 13.7) 250 200 N = 560 CPU 150 100 GPU#1 50 0 行列のサイズ (N 次正方行列 ) GPU#2

69 まとめ ( 行列積 ) GPU 上での実装方法 1. 行列積の計算式に忠実な実装 (GPU#1) 2. 共有メモリを利用した高速化 (GPU#2) 大きさ16 16の部分行列積に分割共有メモリをキャッシュとして利用 CPU と GPU それぞれでの計算速度を比較 CPU < GPU#1 < GPU#2 5.0 13.7 メモリアクセスの工夫により約 34 倍の高速化

~ テンプレートマッチング ~ 70

71 テンプレートマッチング (1) 画像処理の分野で広く用いられている手法 基板の品質検査 画像中の特定物体 ( 人物など ) の検出 入力画像中からテンプレートに類似する部分を探索 入力画像中に窓を設定 窓を移動させながら類似度を評価 比較 テンプレート 窓

72 テンプレートマッチング (2) 基本的なテンプレートマッチングの戦略 入力画像中の部分画像とテンプレートの類似度を評価 窓を移動させながら類似度計算 位置 (X 軸方向, Y 軸方向 ) 拡大 / 縮小 ( スケール変化 ) 回転 類似度の例 SAD(Sum of Absolute Difference) SSD(Sum of Squared Difference) NCC(Normalized Cross Correlation) 窓の大きさに比例して計算コストが増加 膨大な回数の類似度評価が必要 計算コスト大

デモ ( テンプレートマッチング ) 73

74 テンプレートマッチングの実装方法 テンプレートマッチングの特徴 類似度は各位置で独立に計算可能 GPU による並列計算が可能 頻繁に画像メモリへアクセス テクスチャメモリのキャッシュ機能を利用 複数スケールでテンプレートと入力画像を比較 正規化テクスチャ座標とハードウェア線形補間の利用 解像度に依存しないメモリアクセス スケール

75 処理の流れ CPU 側の処理 入力画像とテンプレートの読み込み GPUスレッドの同期 GPU 側の処理 類似度計算 スケールの変更 CPU GPU CPU GPU GPU CPU GPU CPU 画像読み込み メモリ確保 データ転送 類似度計算 スレッド同期 データ転送

76 テクスチャメモリとは? GPU 上の読み取り専用の特殊なメモリ領域 キャッシュを利用した高速なアクセスが可能 2 次元アクセスに対して効率的 ハードウェアを利用した高速な線形補間が可能 テクスチャメモリの定義 texture<datatype, Type, ReadMode> texref; DataType: データの型 (int,float 等 ) Type: テクスチャの種類 (1 次元テクスチャなど ) cudatexturetype1d,cudatexturetype2d, ReadMode: 値の範囲 cudareadmodeelementtype: 各データ型の値を使用 cudareadmodenormalizedfloat: 0 ~ 1 に正規化 ( 符号付きは -1~1)

77 メモリ確保と転送 (1) 入力画像のメモリを GPU 上に確保 テクスチャの定義 ( 画素表現 :RGBA) texture< uchar4, cudatexturetype2d, cudareadmodeelementtype > imgtex; cudaarray を利用して 2 次元テクスチャを確保 cudaarray *imgarray; cudachannelformatdesc c1 = cudacreatechanneldesc< uchar4 >( ); cudamallocarray( &imgarray, &c1, width, height ); 画素表現の指定 CPU から GPU へメモリ転送 転送データ量 cudamemcpytoarray( imgarray, 0, 0, psrc, nbytes, cudamemcpyhosttodevice ); 転送元ポインタ テクスチャメモリへの対応付け cudabindtexturetoarray( imgtex, imgarray, c1 ); CPU から GPU へ転送

78 メモリ確保と転送 (2) テンプレートのメモリを GPU 上に確保 テクスチャの定義 ハードウェア線形補間の有効化 ハードウェア線形補間の利用 (1) texture< uchar4, cudatexturetype2d, cudareadmodenormalizedfloat > reftex; reftex.filtermode = cudafiltermodelinear; ハードウェア線形補間の利用 (2) 正規化テクスチャ座標の有効化 reftex.normalized = 1; (0, 0) (W-1, 0) (0, 0) (1, 0) (0, 0) (1, 0) (0, 0) (1, 0) (0, 1) (0, H-1) (W-1, H-1) (0, 1) (0, 1) (1, 1)

79 GPU 上で類似度計算 global void kernel( float *error, float *scale, int imgw, int imgh, int areaw, int areah, int maskw, int maskh, float s ) { int i = threadidx.x + blockdim.x * blockidx.x; int j = threadidx.y + blockdim.y * blockidx.y; float err = 0.0f; float _1_w = 1.0f / maskw; // テンプレートの幅に対するスケーリング係数 float _1_h = 1.0f / maskh; // テンプレートの高さに対するスケーリング係数 for( int n = 0 ; n < maskh ; n++ ) { for( int m = 0 ; m < maskw ; m++ ) { uchar4 p1 = tex2d( imgtex, i + m, j + n ); float4 p2 = tex2d( reftex, _1_w * m, _1_h * n ) * 255.0f; テクスチャからデータを読み取り } } } err += ( p1.x - p2.x ) * ( p1.x - p2.x ); err += ( p1.y - p2.y ) * ( p1.y - p2.y ); err += ( p1.z - p2.z ) * ( p1.z - p2.z ); err *= _1_w * _1_h; if( error[ i + j * imgw ] > err ) { error[ i + j * imgw ] = err; scale[ i + j * imgw ] = s; } SSD が最小のものを記録 SSD を計算

80 CPU による類似度計算の制御 // スケール & 誤差を保持する GPU 側のメモリ領域 thrust::device_vector< float > derror( img.size( ), 1.0e10f ); thrust::device_vector< float > dscale( img.size( ), 0 ); スケールを変更して再探索 float s = MIN_SCALE; while( s <= MSCALE ) { // 実行する GPU のスレッド数を指定 int W = ( int )img.width( ), H = ( int )img.height( ); int w = ( int )( ref.width( ) * s ), h = ( int )( ref.height( ) * s ); int threadnumx = 16, threadnumy = 16; int blocknumx = ( W - w ) / threadnumx + ( ( ( W - w ) % threadnumx ) == 0? 0 : 1 ); int blocknumy = ( H - h ) / threadnumy + ( ( ( H - h ) % threadnumy ) == 0? 0 : 1 ); dim3 nthreads( threadnumx, threadnumy, 1 ); dim3 nblocks( blocknumx, blocknumy, 1 ); // GPU 側で類似度を計算する kernel<<< nblocks, nthreads >>>( thrust::raw_pointer_cast( &derror[ 0 ] ), thrust::raw_pointer_cast( &dscale[ 0 ] ), W, H, W - w, H - h, w, h, s ); cudadevicesynchronize( ); // 処理が完了するまで待機する } s *= SCALE_FACTOR; // スケールを変更する // GPU の処理結果を転送する thrust::host_vector< float > herror = derror; thrust::host_vector< float > hscale = dscale; GPU を用いて類似度を計算

81 計算時間の比較 (1) CPUとGPUで計算時間を評価 使用計算機 CPU: Intel Core i7-x980(3.33ghz) OpenMPを使用して6スレッドで並列計算 GPU: NVIDIA GeForce GTX480 マルチプロセッサ数 : 15(SP 数 = 15 32 = 480) Memory:1.5 GB OS:Windows 7 SP1

82 計算時間の比較 (2) 実験パラメータ 入力画像 :800 600 画素 テンプレート :105 135 画素 スケール :0.3~1.9 倍 ( 拡大率 1.2) テンプレート 実験結果 CPU: 77.3 sec. GPU: 2.2 sec. 約 35 倍の高速化 結果画像

83 まとめ ( テンプレートマッチング ) テクスチャメモリを利用した高速化 キャッシュを利用した効率的なメモリアクセス 解像度に依存しないメモリアクセス ハードウェア線形補間 正規化テクスチャ座標 CPU と GPU で計算時間を比較 CPU: 77.3 sec. GPU: 2.2 sec. 約 35 倍の高速化を実現

~ テンプレートマッチングを高速化 ~ 84

85 処理の流れ 類似度計算の打ち切り処理を導入 各スケールで計算した類似度の最大値をGPUで計算 Parallel Reduction アルゴリズムを利用 スケールの変更 CPU GPU CPU GPU GPU CPU GPU GPU CPU 画像読み込み メモリ確保 データ転送 類似度計算 スレッド同期 最大類似度計算 データ転送

86 GPU 上で類似度計算 global void kernel( float *error, float *scale, int imgw, int imgh, int areaw, int areah, int maskw, int maskh, float s, float maxerr ) { int i = threadidx.x + blockdim.x * blockidx.x; int j = threadidx.y + blockdim.y * blockidx.y; float err = 0.0f; float _1_w = 1.0f / maskw; // テンプレートの幅に対するスケーリング係数 float _1_h = 1.0f / maskh; // テンプレートの高さに対するスケーリング係数 for( int n = 0 ; n < maskh ; n++ ) { for( int m = 0 ; m < maskw ; m++ ) { uchar4 p1 = tex2d( imgtex, i + m, j + n ); float4 p2 = tex2d( reftex, _1_w * m, _1_h * n ) * 255.0f; } err += ( p1.x - p2.x ) * ( p1.x - p2.x ); err += ( p1.y - p2.y ) * ( p1.y - p2.y ); err += ( p1.z - p2.z ) * ( p1.z - p2.z ); } } if( maxerr < err *_1_w * _1_h ) { break; }... 以下同じ... 前スケールの探索時における SSD が最小値より大きい場合は計算を打ち切り

87 Parallel Reduction による最小値探索 GPU 上のメモリから最小値を計算 GPUのスレッドが並列に値の比較処理を実行 GPU 上のデータ 11 4 27 25 5 13 6 9 20 12 14 7 2 19 3 15 4 25 5 6 12 7 2 3 4 5 7 2 4 2 2 2 同期同期同期

88 Parallel Reduction による最小値探索 GPU 上のメモリから最小値を計算 GPU のスレッドが並列に値の比較処理を実行 Thrust ライブラリの reduce アルゴリズムを利用 thrust::reduce( 先頭, 末尾, 初期値, thrust::minimum<t>( ) ); 末尾を指すイテレータ 先頭を指すイテレータ GPU 上のデータ 11 4 27 25 5 13 6 9 20 12 14 7 2 19 3 15 4 25 5 6 12 7 2 3 4 5 7 2 4 2 2 2 同期同期同期

89 CPU による類似度計算の制御 // スケール & 誤差を保持する GPU 側のメモリ領域 thrust::device_vector< float > derror( img.size( ), 1.0e10f ); thrust::device_vector< float > dscale( img.size( ), 0 ); float s = MIN_SCALE, maxerr = 1.0e10f; while( s <= MSCALE ) { // 実行する GPU のスレッド数を指定 int W = ( int )img.width( ), H = ( int )img.height( ); int w = ( int )( ref.width( ) * s ), h = ( int )( ref.height( ) * s ); int threadnumx = 16, threadnumy = 16; int blocknumx = ( W - w ) / threadnumx + ( ( ( W - w ) % threadnumx ) == 0? 0 : 1 ); int blocknumy = ( H - h ) / threadnumy + ( ( ( H - h ) % threadnumy ) == 0? 0 : 1 ); dim3 nthreads( threadnumx, threadnumy, 1 ); dim3 nblocks( blocknumx, blocknumy, 1 ); // GPU 側で類似度を計算する kernel<<< nblocks, nthreads >>>( thrust::raw_pointer_cast( &derror[ 0 ] ), thrust::raw_pointer_cast( &dscale[ 0 ] ), W, H, W - w, H - h, w, h, s, maxerr ); cudadevicesynchronize( ); // 処理が完了するまで待機する // 誤差の最大値を取得する maxerr = thrust::reduce( derror.begin( ), derror.end( ), 1.0e10f, thrust::minimum< float >( ) ); } s *= SCALE_FACTOR; // スケールを変更する // GPU の処理結果を転送する... Trust ライブラリを用いて GPU を活用して誤差の最小値を計算

90 計算時間の比較 (1) CPUとGPUで計算時間を評価 使用計算機 CPU: Intel Core i7-x980(3.33ghz) OpenMPを使用して6スレッドで並列計算 GPU: NVIDIA GeForce GTX480 マルチプロセッサ数 : 15(SP 数 = 15 32 = 480) Memory:1.5 GB OS:Windows 7 SP1

91 計算時間の比較 (2) 実験パラメータ 入力画像 :800 600 画素 テンプレート :105 135 画素 スケール :0.3~1.9 倍 ( 拡大率 1.2) テンプレート 実験結果 CPU: 8.700 sec. GPU: 0.258 sec. 約 33.7 倍の高速化 結果画像

92 まとめ ( テンプレートマッチング ) 類似度計算の打ち切り処理により高速化 Parallel Reduction アルゴリズムを利用 CPU と GPU で計算時間を比較 類似度の打ち切り なし あり CPU 77.3 秒 8.7 秒 GPU 2.2 秒 0.258 秒

~ ガウシアンフィルタ ~ 93

ガウシアンフィルタに挑戦 94

95 空間フィルタリング ( 線形フィルタ ) 線形フィルタとは? i, j h x, y f i x j y 画像処理の基本的な処理 ガウシアンフィルタ LoGフィルタ ガボールフィルタ 他 h w g, y h x w 出力画像フィルタ 入力画像 フィルタ h(i,j) 積和 出力 g(i,j) 周辺領域 注目画素 f(i,j) 入力画像

ガウシアンフィルタ フィルタ係数に 2 次元ガウス分布を利用 問題点 フィルタ半径が大きくなるにつれ計算コスト大 フィルタを 1 次元ガウス分布の積に分解 96 2 2 2 2 exp 2 1, y x y x h 2 2 2 2 2 2 2 2 exp 2 exp 2 1 2 exp 2 1, y x y x y x h h h y w w x y j x i f x y j i g, 2 exp 2 exp 2 1, 2 2 2 2 X 軸方向への 1 次元ガウシアンフィルタ Y 軸方向への 1 次元ガウシアンフィルタ

97 ガウシアンフィルタの特徴 フィルタ出力は各画素で独立に計算可能 GPU による並列計算が可能 フィルタ係数の算出 フィルタ適用前に事前計算が可能 フィルタ適用中は常に同じ値を参照 コンスタントメモリを利用 注目画素の周辺領域 (1 次元 ) にアクセス 1 次元テクスチャのキャッシュ機能を利用

98 コンスタントメモリとは? GPU 上に実装されている特殊なメモリ領域 読み取り専用 マルチプロセッサ1 基あたり 64 KB キャッシュを利用した高速なアクセスが可能 レジスタとほぼ同じ速度でアクセス可能 コンスタントメモリの定義 constant float coeff[ 512 ]; C 言語の配列のようにアクセス可能

99 処理の流れ GPU 上に 2 つの 1 次元配列を確保 X,Y 軸方向の処理時に入力と出力を入れ替え GPU 内でのメモリ転送コストを削減 CPU GPU CPU GPU GPU CPU GPU GPU CPU 画像読み込み メモリ確保 データ転送 ガウシ ( ア X ン軸フ ) ィルタ スレッド同期 ガウシ ( ア Y ン軸フ ) ィルタ データ転送 入力画像 メモリ 1 メモリ 2 メモリ 2

100 メモリ確保と転送 (1) 入力画像のメモリを GPU 上に確保 テクスチャの定義 texture< float4, 1 > imgtexx; texture< float4, 1 > imgtexy; GPU 上に 1 次元配列を確保 cudamalloc( ( void ** )&idata, nbytes ); cudamalloc( ( void ** )&odata, nbytes ); CPU から GPU へメモリ転送 cudamemcpy( odata, psrc, nbytes, cudamemcpyhosttodevice ); 1 次元テクスチャにマッピング cudabindtexture( 0, imgtexx, odata ); cudabindtexture( 0, imgtexy, idata );

101 メモリ確保と転送 (2) フィルタ係数をコンスタントメモリに確保 コンスタントメモリの定義 constant float coeff[ 512 ]; サイズ指定が必要コンスタントメモリを表す修飾子 CPU から GPU へメモリ転送 cudamemcpytosymbol( coeff, psrc, 512 * sizeof( float ) );

102 計算時間の比較 (1) CPU と GPU で計算時間を評価 使用計算機 CPU: Intel Core i7-x980(3.33ghz) OpenMP を使用して 6 スレッドで並列計算 GPU: NVIDIA Geforce GTX480 マルチプロセッサ数 : 15(SP 数 = 15 32 = 480) Memory:1.5 GB OS:Windows 7 SP1 実験パラメータ 画像サイズ :100 2 ~ 3000 2 画素 σ:5.0,10.0

103 計算時間の比較 (2) 400 350 300 250 [msec.] 計算時間 (N 2 = 2000 2 の場合 ) CPU(σ=5) : 105.2 msec. GPU(σ=5) : 20.9 msec. ( 5.0) CPU(σ=10) : 167.7 msec. GPU(σ=10) : 30.7 msec. ( 5.5) CPU(σ=10) CPU(σ=5) 200 150 100 GPU(σ=10) 50 0 2002 400 2 600 2 800 2 1000 2 1200 2 1400 2 1600 2 1800 2 2000 2 2200 2 2400 2 2600 2 2800 2 3000 2 画像サイズ N 2 GPU(σ=5)

104 まとめ ( ガウシアンフィルタ ) 2 次元ガウス分布を 1 次元ガウス分布の積で表現 ガウシアンフィルタを 1 次元フィルタに分解 コンスタントメモリを利用した高速化 キャッシュを利用したメモリアクセスの高速化 同じ値をスレッド間で共有する場合に有効 CPU と GPU で計算時間を比較 画像サイズが 2000 2 の場合 (σ=5) CPU: 105.2 msec. GPU: 20.9 msec. 約 5 倍の高速化を実現

~SIFT 特徴量の計算 ~ 105

106 SIFT(Scale Invariant Feature Transform) 回転 スケール変化等に頑健な特徴点の検出 画像間のマッチングや物体認識 検出に利用 1,2 1. David G. Lowe, Distinctive image features from scale-invariant keypoints, International Journal of Computer Vision, 60, 2, pp. 91-110, 2004. 2. 藤吉弘亘. "Gradient ベースの特徴抽出 - SIFT と HOG - ", 情報処理学会研究報告 CVIM 160, pp. 211-224, 2007.

107 SIFT のアルゴリズム (1) DoG 画像の作成 (Difference-of-Gaussian) 特徴点の検出 エッジ上の点を削除 サブピクセル位置推定 コントラストの小さい点を削除 3 k 2 k 異なるスケールの平滑化画像の差分 (DoG) を計算 k 周辺 26 画素に対して極値をとる位置を特徴点として検出 DoG3 DoG2 DoG1 DoG3 DoG2 DoG1

108 SIFT のアルゴリズム (2) DoG 画像の作成 (Difference-of-Gaussian) 特徴点の検出 エッジ上の点を削除 サブピクセル位置推定 コントラストの小さい点を削除

109 SIFT のアルゴリズム (3) オリエンテーションの算出 特徴ベクトルを算出 周辺領域の勾配方向と強度からオリエンテーションを算出

デモ ( SIFT ) 110

111 SIFT の実装方法 さまざまなスケールにおける DoG 計算 高速なガウシアンフィルタを利用 キーポイント検出 各画素で独立に判定可能 GPU による並列計算が可能 判定に 26 近傍の画素値が必要 3 次元テクスチャを利用できるか? テクスチャサイズの制限から利用は困難 2 次元テクスチャで代用 オリエンテーションの算出 各キーポイントで独立に計算可能

処理の流れ 112 画像読み込み CPU D o G 計算 GPU スレッド同期 CPU データ転送 CPU GPU データ転送 GPU CPU キーポイント検出 GPU キーポイント DoG 画像オリエンテーション入力画像キーポイントリスト作成 CPU オリエンテーション算出 GPU データ転送 GPU CPU

113 SIFT 計算におけるメモリ配置 複数スケールの DoG 画像をテクスチャ 1 枚に配置 ミップマップ を構築 DoG 画像の出力先をソフトウェア的に調整 1 1 1 2 1 4 1 8 CUDA は未サポート

114 計算時間の比較 (1) CPU と GPU で計算時間を評価 SIFT 特徴量としてオリエンテーションを計算 使用計算機 CPU: Intel Core2 Quad Q9550(2.83 GHz) 4 スレッドで並列計算 GPU: NVIDIA GeForce GTX280 マルチプロセッサ数 : 30(SP 数 = 30 8 = 240) Memory:1.0 GB OS:Windows Vista SP1 実験パラメータ 画像サイズ :50 2 ~ 1000 2 画素 σ 0 :1.6, 分割数 :3

115 計算時間の比較 (2) 2100 2000 1900 1800 1700 1600 1500 1400 1300 1200 1100 1000 900 800 700 600 500 400 300 200 100 [msec.] 計算時間 (N 2 = 600 2 の場合 ) CPU : 458.0 msec. GPU : 10.2 msec. ( 44.9) 0 50 2 100 2 150 2 200 2 250 2 300 2 350 2 400 2 450 2 500 2 550 2 600 2 650 2 700 2 750 2 800 2 850 2 900 2 950 2 10002 画像サイズ N 2 GPU CPU

デモ (SIFT) 116

117 まとめ (SIFT) テクスチャメモリとコンスタントメモリの利用 キャッシュを利用したメモリアクセスの高速化 ミップマップを構築 複数スケールの DoG 画像を 1 枚のテクスチャに格納 CPUとGPUで計算時間を比較 N 2 = 600 2 の場合 CPU: 458.0 msec. GPU: 10.2 msec. 約 44 倍の高速化を実現

~GPGPU による高速画像処理に挑戦して ~ 118

119 GPGPU への挑戦を終えて CUDA を利用することで容易に GPGPU が可能 スレッドプログラミングの経験があれば非常に簡単 既存プログラムの移植も比較的容易 Thrust ライブラリによる簡単な GPGPU GPU を意識せずにプログラミングが可能 GPU を使うと 10 倍以上の高速化が可能? 多くの画像処理アルゴリズムは高速化が可能 空間フィルタリング, 局所特徴量計算, 他 逐次型の画像処理アルゴリズムは高速化が困難 ラベリング, 細線化, 他

120 GPGPU の問題点と今後の展望 複数 GPU の利用 各 GPU 上で別々に処理を実行 CUDA 4.0 は単一 CPUスレッドから複数 GPUを利用可能 CPU との連携 GPUの苦手な処理をCPUで計算 CPUとGPUの役割分担が重要 高性能な GPU が登場 リアルタイム画像処理 ( 大規模計算 ) への挑戦 GPGPU 開発環境の標準化 (OpenCL)

121 参考文献 [1] M. J. Harris, G. Coombe, T. Scheuermann, and A. Lastra, Physically- Based Visual Simulation on Graphics Hardware, Proceedings of SIGGRAPH 2002 / Eurographics Workshop on Graphics Hardware 2002, pp.1-10, 2002.(GPGPU の起源が書かれている論文 ) [2] J. D. Owens, D. Luebke, N. Govindaraju, M. Harris, J. Krüger, A. E. Lefohn, and T. J. Purcell, A Survey of General-Purpose Computation on Graphics Hardware, Computer Graphics Forum, Vol.26, No.1, pp.80-113, 2007.( 最近の GPGPU が詳しく述べられている論文 ) [3] GPGPU, http://gpgpu.org/ [4] CUDA ZONE, http://www.nvidia.com/object/cuda_home.html [5] CUDA Programming Guide, http://www.nvidia.com/object/cuda_develop.html [6] OpenCL, http://www.khronos.org/opencl/

122 MIST(Media Integration Standard Toolkit) 複数メディアを扱うためのライブラリ 音声 画像処理のアルゴリズムを多数実装 C/C++ を用いた高速な処理を実現 C++ のテンプレートを用いた汎用的な実装 複数のプラットフォームで動作 充実した日本語チュートリアルを用意 オープンソースとして公開中 BSDスタイルのライセンス 商用の製品開発でも利用可能 http://mist.murase.m.is.nagoya-u.ac.jp/

123 ~Visual Studio 2010 の詳細設定 ~ CUDA Toolkit 4.0 以降

124 Visual Studio 2010 の簡易設定 (1).cu ファイルの簡易コンパイル設定 プロジェクトメニューのビルドのカスタマイズを表示 CUDA 4.0 を選択.cu ファイルの追加前に行う 1. これを表示 2. これを選択

125 Visual Studio 2010 の簡易設定 (2).cu ファイルのコンパイルオプション.cu ファイルのプロパティを表示 NVCCのコンパイルオプションをGUIで調整可能 これを表示

126 Visual Studio 2010 の簡易設定 (3) プロジェクト プロパティ を選択 構成プロパティ リンカー を選択 追加の依存ファイル cudart.lib を指定 ここに入力

127 Visual Studio 2010 の簡易設定 (3) プログラムのコンパイルと実行 Visual Studio の ビルド を実行 実行ファイルが作成されることを確認 デバッグなしで開始 コマンドラインに Hello World!! が表示される Hello World!! 続行するには何かキーを押してください... 本講演のプログラムを実行した場合

~ 補足資料 ~ 128

129 FLOPS FLoating point number Operations Per Second FLOPSやFLOP/s と表記される 1 秒あたりに実行可能な浮動小数点演算回数 スーパーコンピュータ等の性能を表す指標 代表的な CPU/GPU の FLOPS Core i7-965 51.20 GFLOPS GeforceGTX580 1.58 TFLOPS RadeonHD5870 2.72 TFLOPS GPU は積和演算の性能

130 CUDA におけるエラー処理 API 関数 ( メモリ確保, 他 ) の場合 各 API 関数の戻り値を評価 cudasuccess 実行に成功 GPU で実行する関数の場合 cudathreadsynchronize( ) により同期 cudagetlasterror( ) の戻り値を評価 cudasuccess 実行に成功

~OpenCV + CUDA~ 131

デモ ( 特徴点検出 & 対応付け ) 132

133 OpenCV をダウンロード Subversion 経由で最新版を入手 https://code.ros.org/svn/opencv/trunk/opencv Cmake をダウンロード http://www.cmake.org/

134 OpenCV ライブラリのビルド CUDA を有効にする Cmake の設定で WITH_CUDA にチェック Configure を実行 ビルド設定を反映 Generate を実行 ビルドファイルを生成 3. ビルドファイルの生成 1. ここにチェック 2. 設定を反映 Cmake の設定画面

135 CPU/GPU のコードの違い SURF を用いた特徴点検出 & マッチング (CPU) #include <opencv2/opencv.hpp>... cv::mat_< float > desc2; cv::surffeaturedetector detector2( th2 ); // 特徴点を検出 detector2.detect( frame_gray, keys2 ); // 特徴量を計算 cv::surfdescriptorextractor extractor; extractor.compute( frame_gray, keys2, desc2 ); // 特徴点を対応付け cv::bruteforcematcher< cv::l2< float > > matcher; matcher.match( desc1, desc2, matches );

136 CPU/GPU のコードの違い SURF を用いた特徴点検出 & マッチング (GPU) #include <opencv2/gpu/gpu.hpp>... cv::gpu::gpumat desc_gpu2; cv::gpu::surf_gpu detector2; detector2.hessianthreshold = th2; // データを GPU へ転送 cv::gpu::gpumat frame_gpu; frame_gpu.upload( frame_gray ); // 特徴点検出 & 特徴量計算 detector2( frame_gpu, cv::gpu::gpumat(), keys2, desc_gpu2 ); // 特徴点を対応付け cv::gpu::bruteforcematcher_gpu< cv::l2< float > > matcher; matcher.match( desc_gpu1, desc_gpu2, matches );

137 OpenCV で GPU を利用できる機能 行列演算 テンプレートマッチング 歩行者検出 (HOG) 特徴点検出 & 特徴量 (SURF) 特徴点対応付け 画像フィルタ ラプラシアン, ソーベル, ガウシアン, 他 カメラキャリブレーション