Transcript PPT

ＦＰＧＡを用いた
レイトレーシングの実装
福永研究室
西村建郎
中村真魚
工藤健史
発表の流れ
•
•
•
•
•
•
•
レイトレーシングの概要（西村）
ＰＯＣＯの概要（中村）
ＩＥＥＥ754（浮動小数点規格）について（中村）
レイトレーシングを行うための改良点（工藤）
アセンブラの解説（中村）
コンパイラの解説（西村）
ＦＰＧＡへの実装（工藤）
レイトレーシングの概要
• レイトレーシングとは、3次元の物体の座標
データや視点・光源の位置などのデータを計
算して、画像を描画（レンダリング）する手法
の一つである。
（例）球に左から照射した場合
レイトレーシングの概要
• 物体の表面の反射率や透明度・屈折率など
を細かく反映させることが出来るのが特徴。
• 欠点としては、計算量が多い。
レイトレーシングの概要
• 視点からスクリーンの画素に向かう視線
（Ray）を延長（追跡、trace）し、物体と交差す
れば、物体の色や、光の反射等の計算をし、
色を決定します。物体に当たらなければ、そ
こには何もないので、背景の色となる。
ピクセル
視点
・
・
・
・
・
・
光源
・
・
・
スクリーン
レイトレーシングの概要
• 物体をどのように表すか。
• また、どのように物体との交差を判定するか。
レイトレーシングの概要
• 球の場合
P:点（変数）, Pc:球の中心, R:半径
Pv:視点, V:視線ベクトル, t>1:任意の実数 とおくと
（視線ベクトルは、[スクリーンの各座標]-[視点]で求める）
球の方程式：（P－Pc）²＝R²
視線：P＝Pv＋V×t
ピクセル
視点
・
・
・
・
・
・
・
・
・
スクリーン
レイトレーシングの概要
球の方程式：（P－Pc）²＝R²
視線の方程式：P＝Pv＋V×t
視線の式を球の方程式に代入して、
tについて解く
実数解が存在する⇒球と交わる⇒色を計算
実数解が存在しない⇒交わらない⇒背景色
（判別式を用いる）
レイトレーシングの概要
• 球以外の形状の場合
処理が容易な三角形の集合としてモデリングする。
すべての三角形の各頂点を座標で与える。
（例）平面
レイトレーシングの概要
• 視線と三角形の交差判定
三角形を含む平面の方程式を求める
平面と視線の交点を求める
交点が三角形の中にあるかを判定する
レイトレーシングの概要
• 三角形モデリングの欠点
モデルが複雑になると、視線との交差判定を
すべての三角形と行わなくてはならないため、
計算量が膨大になる。
今回の研究では、球のみを扱うことにする。
レイトレーシングの概要
今までの判定法だけを用いて球を描画すると
まだ光についてなにも考慮されていない。
明るさを決める要素
①環境光
面の傾きに関係なく明るさが決まる
（光源からの光が直接当たっていなくても、真っ
暗にはならない）
明るさを決める要素
②拡散反射
光源から入ってくる光がすべての方向に均一
に反射する。
反射する光の量は、入射する光の量に応じ、
入射する光の量は、光の当たる角度によるの
で、物体に影がつく。
反射光の強さ＝光源の明るさ×反射率×cosα
入射光と法線の成す角
明るさを決める要素
• 拡散反射のみの画像
明るさを決める要素
③鏡面反射
鏡のように、一方向からの光が、入射した角度
のちょうど反対側に強く反射する。反射した光の先
に視点があるとき、一点周辺だけ極端に光って見
える。
n:鏡面度（10～20）
反射光の強さ＝光源の明るさ×反射率×cosⁿγ
反射光と視線の成す角
明るさの決め方
・視線と物体が点Ｐで交差した場合
点Ｐでの拡散反射光と鏡面反射光の明るさを
計算し、環境光・拡散反射光・鏡面反射光のう
ち最も明るいものを点Ｐの明るさとする。
点Ｐの明るさ
＝ max {環境光, 拡散反射光, 鏡面反射光}
レイトレーシング結果
・C言語で実行
POCOの概要
• ゼミではコンピュータアーキテクチャの勉強を行った（4月～7月）。
– 教材：「作りながら学ぶコンピュータアーキテクチャ」天野英晴,西村克信
– 題材：ハードウェア記述言語「Verilog HDL」を使ったディジタル回路の設計
– POCOとは
•
教本中で開発されたRISC (Reduced Instruction Set Computer)。
演算装置（ALU）、データメモリ、命令メモリ、レジスタ、
プログラムカウンタで構成され、
整数の加算および減算と、論理演算が可能。
– 使用したソフトウェア
•
ISE Simulator( ISIM ) / Xilinx
コンピュータの基本形
CPU（中央処理装置）
• 制御装置
• 演算装置
• レジスタ
etc…
I/O
Memory
• モニタ
• マウス
• キーボード
etc…
アキュムレータマシンの基本形
ＡＬＵ
ＡＬＵ
• Arithmetic Logic Unit の略
• 加算と減算、論理演算を行う回路。
– 論理演算の例
• SL,SR
• AND,OR
accum ← X = 01010
SL … accm ← X<<1 = 10100
SR … accm ← X >>1 = 00101
AND 10101 … accum ← X AND {10101} = 00000
OR 10101 … accum ← X OR {10101} = 11111
アキュムレータマシンの基本形
プログラムカウンタ
命令メモリ
アキュムレータ
データメモリ
アキュムレータマシンの動作
1. pc の指し示す命令を取ってくる（命令フェッ
チ）。同時にpcに１を加える。
2. 命令を解読する（命令デコード）。
3. 命令を実行する。その後１に戻る
POCOの概要
POCOの概要
• 主な変更点
1. アキュムレータマシンから汎用レジスタマシンに
一時的なデータの格納場所が増え、用途に応じたレジスタの使い分
けがし易くなった（データメモリのポインタの実現）。
2. 命令デコード部分の複雑化
命令の種類として、
1.
2.
3.
ALUを使う命令
レジスタが対象の命令
pcが対象の命令
に分類できる。
3. pcのためのジャンプ命令の追加
プログラム中に繰り返し使われる命令群（サブルーチン）を呼び出す
ことが可能になる。
POCOで可能な命令
ビット幅 : 16bit
機械語命令長 : 16bit
オペコード 5bit
オペランド 3bit
レジスタ : 8個
メモリ : 16bit×アドレス2^16 = アドレス空間64K
テストベンチ
POCO
ALU
と階層的になっており、
各々部品のように分かれている。
これをモジュールと呼ぶ。
•
•
•
POCOの構造
テストベンチ
作成した回路が動くかシミュレーションを
行う際に、信号を与えたり観察を行う役
割がある。ここではPOCOにクロック信号
等を与えている。
クロック
CPUの動作基準となる時間の単位。
動作周波数
1秒間に刻むクロックの数を表す。単位は
Hz(ヘルツ）。
1秒間に100回のクロックを刻むと100Hz。
POCOでは1クロック=10ns。
1秒間に100000000回クロックを刻むため、
動作周波数は100MHzとなる。
IEEE754の概要
• IEEE(Institute of Electrical and Electronics Engineers)
電子部品や通信方式の研究や標準化を行う機関。
「米国電気電子学会」とも呼ぶ。
• IEEE754
IEEEが定めた浮動小数点数の演算の標準規格。
基本形式
有限数は符号s,仮数f,指数eの三つの整数で表
現でき
(-1)^s × f × b^e
と表せる(bは基数)。
例) -3.25 = (-1)^1 × 325 × 10^(-2)
(11.01)2 = (-1)^0 × (1101)2 × 2^(-2)
固定小数点数形式の二進数は桁が上がると重みが2倍になる
1
1
.
0
1
=
2^1
2^0
2^-1
2^-2
3.25
基本形式
• 単精度(32bit)
Sign bit
Exponent bit
8bit
• 倍精度(64bit)
– Sign bit … 1bit
– Exponent bit … 11bit
– Fraction bit … 52bit
Fraction bit
23bit
基本形式
例) (5.25)10
= (101.01)2
= (1.0101)2 × 2^2
符号 sign = 0 (正)
指数 exponent = 2 + 127(バイアス) = 129
仮数 fraction = 0101 0000 0000 … 0000
0 1000 0001 0101 0000 … 0000
バイアスを足すことで指数を
-127 ~ 128 の範囲で表せる。（ゲタばき表現）
表現可能な数
指数部(ゲタばき)
仮数部
意味
0
0
0
0
≠0
±非正規化数
1 ~ 254
任意
浮動小数点数（正規化数）
255
0
±無限大
255
≠0
Not a number(NaN) … 非数
• 例外処理
– NaN … 無効な演算の結果。（0/0や∞-∞など）
– 無限大 … 0以外の数を0で割った演算の結果(0除算)。
また正規化数で表せないほど大きな値（オーバーフロー）。
– 非正規化数 … 0より大きく、正規化数より小さい値。
非正規化数の範囲（絶対値）
0.0000 0000 0000 0000 0000 001 × 2^(-126)
~
0.1111 1111 1111 1111 1111 111 × 2^(-126)
正規化数の範囲（絶対値）
1.0000 0000 0000 0000 0000 000 × 2^(-126)
~
1.1111 1111 1111 1111 1111 111 × 2^(127)
丸め
• 最近接丸め（偶数）
最も近くの値に丸める。表現できる２つの値の中間の時（末
尾より下位の桁が末尾×1/2と等しい）、末尾が0になる方に
丸める
• 最近接丸め（0から遠い方）
最も近くの値に丸める。表現できる２つの値の中間の時（末尾より下位の
桁が末尾×1/2と等しい）、正なら切り上げ、負なら切り下げ
• 切り捨て
0に近い側に丸める
• 切り上げ
正の無限大に近い方に丸める
• 切り下げ
負の無限大に近い方に丸める
最近接丸め
1 2 3
・ ・ ・ ・
23
1.
仮数部の演算結果
Guard bit
Round bit
Round bit
• 丸めのために仮数部の下位３bitを使用する
– ガードビット … 仮数部の最下位ビット（ulp）の1/2の重さを持つ
– ラウンドビット …仮数部の最下位ビット（ulp）の1/4の重さを持つ
– スティッキービット … 仮数部の下位3ビット目以降のORを取った値。
最近接丸め
仮数部の最下位
ビット
Guard ビット
Round ビット
Sticky ビット
丸めの結果
0
0
0
1
0
0
最近偶数への丸め
（切り捨て）
0
1
0
1
切り上げ
0
1
1
1
切り上げ
1
0
1
1
0
0
最近偶数への丸め
（切り上げ）
1
1
0
1
切り上げ
1
1
1
1
切り上げ
切り捨て
切り捨て
計算機概要
変更点
– レジスタの拡張 8 → 32×2
（
16bit： 整数用 32bit：浮動小数点用)
– Moduleの追加 「fpu」及び「frfile」の実装
– Instructionの拡張 16bit→25bit
s-reg d-regに加えてもう一つ指定するregisterを追加
（計算結果を3番目に指定するレジスタに入れる）
– 命令の追加
浮動小数点での計算・分岐命令
F_ADD
F_SUB
1. 仮数同士を計算させる前に指数、符号を確認
2. 指数によって小数点の位置がずれているので、指数が小さい
（大きい）方を差分右（左）にシフトする （桁合わせ）
3. 数の符号、大小によって計算する命令を変える
4-(-5)=4+5
4 + ( - 6 ) = - (6 - 4)
4. 計算
5. 繰り上がりや繰り下がりの有無を確認し修正(正規化）
F_MUL
• 符号
• 指数
両数の指数を足して、バイアスを引く（筆算）
• 仮数の計算 boothアルゴリズム
• 一般に仮数の計算は筆算と同じ要領で行われる
乗数にある「1」 の分だけ計算が行われる（計算量が多い）
3
0
1
1
6
1
1
0
0
+6
0
1
+ 12
0
1
1
18
1
0
0
1
1
0
3*6=18
Boothアルゴリズム
3
被乗数を二桁ずつ見ていく手法
被乗数の先頭と最後に0を拡張
後ろから２bitずつ見ていく
「00」「11」なら、何も加算しない
「10」であれば乗数を「引く」
「01」であれば乗数を「足す」
7
1
1
１
3=11であるが
0 1 1 0と拡張する
+3
1
+6
1
+12
考え方として、右の表のような計算で
（3*4)+(3*2)+(3*1)=21でなく
(3*8)-(3*1)=21といった計算方法をとる
1
1
21
1
1
1
0
1
- 3
+24
1
1
21
1
0
1
1
1
0
1
1
1
0
1
0及び1が連続すれば計算を省略できるが、交互に出ると計算量は多
くなる
F_DIV
• 引き放し法
要領は筆算と同じ
被除数が大きければ１、小さければ０．商として
書き込んだら被除数を左シフトする
を繰り返す
• Newton-Raphson法
商を近似させていく方式
被除数をa、除数をbとすると a/b=a*(1/b)
下記の式はx=1/bとして、xを近似させていく
F(x)=(1/x)-bとして、F(x)=0となったら近似完了
F_SQRT
• Newton-Raphson法による近似によって生成
• 最終的にx(i)がaの平方根の逆数となるように
近似していく
• x(0)を初期値として、ある程度aの値によって
近似された初期値をテーブル化しておく
• ３回のループ
追加命令 その他
• IECONV
浮動小数点レジスタにおける値を2進数の値に変換
ppmファイルに書き込む際に変換する必要がある
• 浮動小数点レジスタの値を参照する分岐命令の
追加
• F_LD命令
参照するレジスタを整数用のものとした
→ 整数用のレジスタをポインタとして使用
• 丸め
最近接偶数丸め （5以下切り捨て）
下３桁のみをチェック
レイトレーシングのアセンブリ
• レイトレーシング法とは
– 視点からスクリーンを1画素ずつ見て、そのときの
視線ベクトルと球が交差するとき、光源からの光
線の拡散反射による光の強度を計算する。
• Lambertの余弦法則とは
– 物体の表面で反射する光の強度は、光線の入射
ベクトルLと表面に対する法線ベクトルNとが成す
角度の余弦に比例する。
レイトレーシングのアセンブリ
全体の流れ
1. 交差判定
視線ベクトルPと球の方程式
P=Pv+vt
(P-Pc)^2 = R^2
からtに関する二次方程式を作る。
実数解が存在するとき視線が球と交わる。
2. 拡散反射の計算
光線ベクトルvLと交点の法線nが成す角αに対し、
反射光の強さは
(Idcosα)×KD
となり、
A・B = |A||B|cosα
から
cosα = ((-vL)・n) / (|vL||n|)
= -vL・n / √(vL・vL)√(n・n)
で求められる。
3. スクリーンの画素を１つ進めて、１に戻る
レイトレーシングのアセンブリ
交差判定
C言語
A=InnerProduct(v,v);
B2=InnerProduct(Pvc,v);
C=InnerProduct(Pvc,Pvc)-R*R;
D=B2*B2-A*C;
アセンブリ
//==============
// A = v・v
//==============
FMV 25 8
FMV 26 9
FMV 27 10
FMV 28 8
FMV 29 9
FMV 30 10
JAL InnerProduct
FMV 11 25
//==============
//
B2 = Pvc・v
//==============
FMV 25 5
FMV 26 6
FMV 27 7
JAL InnerProduct
FMV 12 25
//====================
//C = Pvc・Pvc - R*R
//====================
FMV 25 5
FMV 26 6
FMV 27 7
FMV 28 5
FMV 29 6
FMV 30 7
JAL InnerProduct
FMV 13 25
FLD 25 16
FMUL 25 25 25
FSUB 13 25 13
//===================
// D = B2*B2 - A*C
//===================
FMUL 12 12 14
FMUL 11 13 25
FSUB 14 25 14
内積のサブルーチン
Ｃ言語
double InnerProduct( double a[3], double b[3] )
{
return a[0]*b[0] + a[1]*b[1] + a[2]*b[2];
}
アセンブリ
//######################
//#
InnerProduct
//######################
FMUL 25 28 25
FMUL 26 29 26
FMUL 27 30 27
FADD 25 26 25
FADD 25 27 25
JR 31
レイトレーシングのアセンブリ
拡散反射
C言語
アセンブリ
t0 = (-B2-sqrt(D))/A;
テキスト参照
for(k=0;k<3;k++) {
P[k] = Pv[k]+v[k]*t0;
N[k] = P[k]-Pc[k];
vL[k] = P[k]-PL[k];
}
cosA = -InnerProduct(vL,N);
cosA /=
sqrt(InnerProduct(vL,vL)*InnerProduct(N,N));
計算結果
• 視点座標
(0.0,0.0,2000.0)
• 球の中心座標
(0.0,0.0,-6000.0)
• 球の半径
R = 200
• 球の色
r = 255, g = 0, b = 0
• 背景の色
r = 0, g = 0, b = 255
• 光源座標
(20000.0,60000.0,90000.0)
• スクリーンの大きさ
640 × 640
コンパイラの概要
• コンパイラの役割
アセンブリ表記をPOCO用の機械語に変換する
（人間の手で機械語に直すのは大変すぎる）
（例）LD 2 3 → 0000000010000110000001001
・ 使用したプログラム言語
Haskell（純粋関数型言語）
Haskellの特徴
• Haskellでの関数とは、ある一つ以上の引数が与
えられると、ある一つの結果を返す変換器であ
る。
• 関数は同じ引数で実行されると、必ず同じ結果
を返すことが保証されている。
⇨ 変数が存在しない
• 関数を組み合わせて新しい関数を作ることが、
Haskellのプログラミングである。
Haskellの特徴
• リストを多用する（for文の代用）
• 記述が数学的
• 遅延評価
（例） sum [1..100] ⇨ 5050
[(x,y)|x<-[1..5], y<-[1..5], x+y==5]
⇨[(1,4),(2,3),(3,2),(4,1)]
zip[1..] [“a”, “b”, “c”]
⇨[(1,”a”),(2,”b”),(3,”c”)]
コンパイラの大まかな流れ
① コメントと空白行を消去
（プログラムから抜粋）
comment::[String]->[String] ⇐関数の型宣言
comment xs = filter comment2 xs
comment2::String->Bool
comment2 "" =False
comment2 ('/':('/':_)) =False
comment2 _ = True
コンパイラの大まかな流れ
②ラベルを数値に変換
First: LD 1 2
ADD 3 4
ST 5 6
JMP First
LD
ADD
ST
JMP
1 2
3 4
5 6
-4
各行に、行番号を与えることにより計算
コンパイラの大まかな流れ
• 各命令、数値をそれぞれ対応するバイナリー
（２進数）へ変換する
ADD →00110
LD →01001
2 →00010
5 →00110
コンパイラの実行例
FPGAにおける実行
FPGAにUSBが付いていないので実行されたデ
ータを読み込めない
320×320pixel一つ一つのデータを見るのは
厳しい（時間が掛かる）

Chipscopeによって、FPGAを流れる波形
データの要点だけを確認
Hz ヘルツ
• PCの速度を表す単位
• Hzは1秒あたりの電磁放射の振動の数を指す
1キロヘルツ
kHz
103Hz
1 000Hz
1メガヘルツ
MHz
106Hz
1 000 000Hz
1ギガヘルツ
GHz
109Hz
1 000 000 000Hz
1テラヘルツ
THz
1012Hz
1 000 000 000 000Hz
1ペタヘルツ
PHz
1015Hz
1 000 000 000 000 000Hz
1エクサヘルツ
EHz
1018Hz
1 000 000 000 000 000 000Hz
背景が青 球の色が赤となる実行結果を使用
赤となる値を書き込む時に「count」を1加算していく
1. 0から1になった時のx,y座標の値
2. 命令が終了した際の「count」の値を確認
3. fwrite_opの後にcount及びpcが遷移しているか
の確認
この３つの整合性が取れれば、FPGAでもシュミレ
ーションと同じ計算ができているはず
LEDランプ
常に光っている
クロックと同期して点灯
スタンバイ状態
STARTを押すと「3」が点灯する
稼働中
countが0から1になるとランプ「4」に移動
4 countがとりうる最大値になるとランプ「3」に移動
5 終了すると点灯
0
1
2
3
スイッチ
0 RESET
1 START
1. Isimにおける「count」が０から１への変化した
Ｘ，Ｙ座標
• Ｘ：c0c00000
• Ｙ：42140000
2. 最終的な描画数 4077
3. fend_op 命令の時、pcが止まる（
imem=5‘b11111)
CPI=４ ４クロックで１つの命令
35MHzまでの動作を確認
• １クロックで除算や平方根の近似を行っている為
と考えられる
FPGAにおける実行
• CPI=8
• 5MHzで実行している
35MHzだとランプの動作が確認できない為
課題
• パイプライン化による計算の高速化
– パイプライン化したmoduleを作成したが、NOPを2
つ分岐命令毎に挿入する必要がある
• div,sqrtによる計算をクロックを跨いでの命令に
変更
– 40MHz以上で実行できない理由は恐らくこの2つ
– もしくはアルゴリズム自体の変更