unit s_MathPack;

//SANDMANのための演算パッケージ


// 参考資料
//
// AMD社      3DNow! テクノロジマニュアル
//            AMD K6-2 Processor Code Optimization
//            AMD Processor Recognition
//
//一部のルーチンの作成にあたって、AMD-3D SDKを参考にしましたm(__)m


//↓の一行をコメントアウトすれば、DirectX7仕様になります
{$DEFINE DX8}

//↓の一行をコメントアウトすれば、3DNow!を使わないコードになります
//現在は、3DNow!対応を見合わせていますので、必ずコメントアウトしてください(^^;)
//{$DEFINE USE_3DNOW}

interface

uses
  {$IFDEF DX8}
  s_3DNow,d3d8,Windows;
  {$ELSE}
  s_3DNow,d3dtypes,Windows;
  {$ENDIF}



//DirectX8用

{$IFDEF DX8}
const
  D3DFVF_VERTEX	= ( D3DFVF_XYZ or D3DFVF_NORMAL or D3DFVF_TEX1 );
  D3DFVF_LVERTEX	= ( D3DFVF_XYZ or D3DFVF_DIFFUSE or D3DFVF_SPECULAR or D3DFVF_TEX1 );
  D3DFVF_TLVERTEX	= ( D3DFVF_XYZRHW or D3DFVF_DIFFUSE or D3DFVF_SPECULAR or D3DFVF_TEX1);

type
  LPD3DTLVERTEX = ^D3DTLVERTEX;
  D3DTLVERTEX = record
	case Integer of
	0: (
	  sx	: Single;
	  sy	: Single;
	  sz	: Single;
	  rhw	: Single;
	  color	: D3DCOLOR;
	  specular	: D3DCOLOR;
	  tu	: Single;
	  tv	: Single;
	);
	1: (
	  dvSX	: Single;
	  dvSY	: Single;
	  dvSZ	: Single;
	  dvRHW	: Single;
	  dcColor	: D3DCOLOR;
	  dcSpecular	: D3DCOLOR;
	  dvTU	: Single;
	  dvTV	: Single;
	);
  end;

(*
 * Untransformed/lit vertices
 *)
  LPD3DLVERTEX	= ^D3DLVERTEX;
  D3DLVERTEX = record
	case Integer of
	0: (
	  x	: Single;
	  y	: Single;
	  z	: Single;
	  color	: D3DCOLOR;
	  specular	: D3DCOLOR;
	  tu	: Single;
	  tv	: Single;
	);
	1: (
	  dvX	: Single;
	  dvY	: Single;
	  dvZ	: Single;
	  dcColor	: D3DCOLOR;
	  dcSpecular	: D3DCOLOR;
	  dvTU	: Single;
	  dvTV	: Single;
	);
  end;

(*
 * Untransformed/unlit vertices
 *)

  LPD3DVERTEX	= ^D3DVERTEX;
  D3DVERTEX = record
	case Integer of
	0: (
	  x	: Single;
	  y	: Single;
	  z	: Single;
	  nx	: Single;
	  ny	: Single;
	  nz	: Single;
	  tu	: Single;
	  tv	: Single;
	);
	1: (
	  dvX	: Single;
	  dvY	: Single;
	  dvZ	: Single;
	  dvNX	: Single;
	  dvNY	: Single;
	  dvNZ	: Single;
	  dvTU	: Single;
	  dvTV	: Single;
	);
  end;

{$ENDIF}

type
  PSingle = ^Single;

//ベクトル(2D)
type TVector2D = packed record
  X,Y:Single;
end;
PVector2D = ^TVector2D;


//3次元のベクトルはD3DVector、ね。


//同次3次元のベクトル(同次なので4つ要素があります)
type THgVector = packed record
  X,Y,Z,W:Single;
end;
PHgVector = ^THgVector;


//クォータニオン(THgVectorとはあえて互換性を持たせない)
type TQuaternion = packed record
  case Integer of
    0: (X,Y,Z,W:Single);
    1: (
      V:D3DVector;
      );
end;
PQuaternion = ^TQuaternion;


//BOX(AABB)
type TBox3D = packed record
  Min,Max:D3DVector;
end;


//3x3行列
type TMatrix33 = packed record
  Case Integer of
  0: (
    _11,_12,_13:Single;
    _21,_22,_23:Single;
    _31,_32,_33:Single;
  );
  1: (
    A,B,C,D,E,F,G,H,I:Single;
  );
end;


//平面 (ax+by+cz+d = 0)
TPlane3D = packed record
  Case Integer of
  0:(A,B,C,D:Single;);
  1:(Normal:D3DVector;);
end;


//球
TSphere3D = packed record
  Case Integer of
  0:(X,Y,Z,R:Single);
  1:(
    Center:D3DVector;
    Radius:Single;
  );
end;

//線分…点A,Bを始点・終点とする線分
type TLineSeg3D = packed record
  Case Integer of
  0:(A,B:D3DVector;);
  1:(_Start,_End:D3DVector;);
end;

//直線、半直線…点Pを通り、傾きDirの直線
type TLine3D = packed record
  Case Integer of
  0:(P,Dir:D3DVector;);
  1:(Origin, Direction:D3DVector;);
end;

//円筒… Axisを軸にする、半径Radiusの円筒
type TCapsule3D = packed record
  Case Integer of
  0:(
    Axis:TLineSeg3D;
    Radius:Single;
  );
  1:(
    A,B:D3DVector;
    R:Single;
  );
end;

//表裏の有る三角形
type TFacedTriangle3D = packed record
  Case Integer of
  0:(
    A,B,C:D3Dvector;
    N:D3DVector;
  );
  1:(
    Vertices:Array[0..3] of D3Dvector;
    Normal:D3DVector;
  );
end;





//同次行列…は、D3DMATRIX

const
  //単位行列
  NowIdentityMatrix:D3DMATRIX = ( _11:1; _12:0; _13:0; _14:0;
                                  _21:0; _22:1; _23:0; _24:0;
                                  _31:0; _32:0; _33:1; _34:0;
                                  _41:0; _42:0; _43:0; _44:1);

  //ゼロ行列
  NowZeroMatrix:D3DMATRIX = (     _11:0; _12:0; _13:0; _14:0;
                                  _21:0; _22:0; _23:0; _24:0;
                                  _31:0; _32:0; _33:0; _34:0;
                                  _41:0; _42:0; _43:0; _44:0);

  //単位行列
  NowIdentityMatrix33:TMatrix33 = (   _11:1; _12:0; _13:0;
                                      _21:0; _22:1; _23:0;
                                      _31:0; _32:0; _33:1);

  //ゼロ行列
  NowZeroMatrix33:TMatrix33     = (   _11:0; _12:0; _13:0;
                                      _21:0; _22:0; _23:0;
                                      _31:0; _32:0; _33:0);

///////////////////////////////////////////////////////////////////////////////
// 3DNow!演算パッケージ
//
// Delphiにおける、引数の取り扱いについての注意
//
// EAX ... 第一引数
// EDX ... 第二引数
// ECX ... 第三引数または返り値
//
// var,constで引数が宣言されているならポインタ扱い。
// 返り値は、recordの場合、返り値へのポインタ
// 単純型の場合は、整数型はそのままでOK、実数型はST(0)に積む
//
// 4つ以上の引数がある場合は、スタックに積まれる。
// 尚、手続きに突入する際に、mov ebp,esp が行われるので、ebpを使って参照可
///////////////////////////////////////////////////////////////////////////////



//まずこれを呼ぶべし…Initialization部で呼ばれるけど
procedure NowInit;


//////////////////////コンストラクタ
function  Vector2D(_x,_y:Single):TVector2D;
function  Vector(_x,_y,_z:Single):D3DVECTOR;
function  HgVector(_x,_y,_z,_w:Single):THgVector;
function  Quaternion(_w,_x,_y,_z:Single):TQuaternion;
function  Sphere3D(_x,_y,_z,_r:Single):TSphere3D;
function  Plane3D(_a,_b,_c,_d:Single):TPlane3D;
function  Line3D(_p, _dir:D3DVector):TLine3D;
function  LineSeg3D(_start, _end:D3DVector):TLineSeg3D;

//////////////////////色調
function NowScaleColor(const col:D3DCOLOR; const scale:Single):D3DCOLOR;
function NowAddColor(const col1:D3DCOLOR; const col2:D3DCOLOR):D3DCOLOR;
function NowLerpColor(const col1:D3DCOLOR; const col2:D3DCOLOR; const t:Single):D3DCOLOR;
function NowHSBToRGB(h,s,b:Single):D3DCOLOR; 

//////////////////////変換

//Int <--> Single 変換(端数切捨て)
function  NowItoF2D(const a:TPoint):TVector2D;
function  NowFtoI2D(const a:TVector2D):TPoint;

//Int <--> Single 変換(端数四捨五入)
function  NowRound2D(const a:TVector2D):TPoint;

//同次ベクトルと普通のベクトルの変換
function  NowHomoginize(const a:D3DVector):THgVector;   //w = 1を付ける
function  NowHeteroginize(const a:THgVector; _div:Boolean = False):D3DVector; //_divがFalseなら、Wを捨てるだけ、divがTrueなら、wでx,y,zを割る

//LVからTLVへ
function  NowLVToTLV(const lv:D3DLVERTEX; const transMat:D3DMATRIX):D3DTLVERTEX;
procedure NowLVToTLVBatch(src:LPD3DLVERTEX; dest:LPD3DTLVERTEX; const transMat:D3DMATRIX; count:Integer);

//角度系の変換
//function NowDegree():Integer; //360度系から、4096度系へ

//TMatrix33からD3DNATRIXへ、
function NowMatrix33To44(m:TMatrix33):D3DMATRIX;
//44行列から左上3x3要素だけ抜いて33行列にする
function NowMatrix44To33(m:D3DMATRIX):TMatrix33;

//３頂点から平面(v1,v2,v3は法線に対して左回り、つまり、３角形が「見える」時、v1,v2,v3は視線に対して右回り)
function NowPlaneFromPoints(v1,v2,v3:D3DVector):TPlane3D;

//１頂点vと法線nから平面
function NowPlaneFromVertexNormal(v,n:D3DVector):TPlane3D;



//////////////////////補間

//スプライン補間(4点を与えられて、中間2点のブレンド比を与えて補間
function NowInterpolate(const p1,p2,p3,p4:Single; const rate:Single):Single;

//線形補間
function NowLerp(const p1,p2:Single; const rate:Single):Single;


//////////////////////その他スカラ演算

//Math入ってるArcTan2みたいなもの？ x=0でもOK
function NowArcTan2(Y, X: Extended): Extended;




//////////////////////2D演算

//ベクトルの足し算
function NowAdd2D(const a:TVector2D; const b:TVector2D):TVector2D;

//ベクトルの引き算
function NowSub2D(const a:TVector2D; const b:TVector2D):TVector2D;

//内積
function NowDotProduct2D(const a:TVector2D; const b:TVector2D):Single;

//2点の距離
function NowDistance2D(const a,b:TVector2D):Single;

//長さを算出
function NowLength2D(const a:TVector2D):Single;

//ベクトルの向きを4096度系で
function NowVecToAngle2D(const a:TVector2D):Single;

//回転(原点を右回りにangleだけ)
function NowRot2D(const a:TVector2D; const angle:Integer):TVector2D;

//正規化(長さを1に)
function NowNormalize2D(const a:TVector2D):TVector2D;

//スケーリング
function NowScale2D(const a:TVector2D; const scale:Single):TVector2D;

//ベクトルの長さをkにする
function NowChangeLength2D(const a:TVector2D; const k:Single):TVector2D;

//posからtargetに向いた長さkのベクトルを生成
function NowLookAt2D(const pos:TVector2D; const target:TVector2D; const k:Single):TVector2D;

//長さkのランダムなベクトルを作る
function NowRandom2D(const k:Single):TVector2D;

//線形補間
function NowLerp2D(const a:TVector2D; const b:TVector2D; const rate:Single):TVector2D;


//////////////////////3D演算(ベクトル)

//足し算
function NowAdd(const a:D3DVECTOR; const b:D3DVECTOR):D3DVECTOR;
function NowAddHg(const a:THgVector; const b:THgVector):THgVector;

//引き算
function NowSub(const a:D3DVECTOR; const b:D3DVECTOR):D3DVECTOR;
function NowSubHg(const a:THgVector; const b:THgVector):THgVector;

//内積
function NowDotProduct(const a:D3DVECTOR; const b:D3DVECTOR):Single;

//外積
function NowCrossProduct(const a:D3DVECTOR; const b:D3DVECTOR):D3DVECTOR;



//原点回りの回転(各軸を左回りに)
function NowRotX(const a:D3DVECTOR; const t:Integer):D3DVECTOR;
function NowRotY(const a:D3DVECTOR; const t:Integer):D3DVECTOR;
function NowRotZ(const a:D3DVECTOR; const t:Integer):D3DVECTOR;

//ベクトルの反転
function NowNeg(const a:D3DVECTOR):D3DVECTOR;

//正規化
function NowNormalize(const a:D3DVECTOR):D3DVECTOR;

//長さを得る
function NowLength(const a:D3DVECTOR):Single;

//スケーリング
function NowScale(const a:D3DVECTOR; const scale:Single):D3DVECTOR;

//長さを k に
function NowChangeLength(const a:D3DVECTOR; const k:Single):D3DVECTOR;

//線形補間
function NowLerp3D(const a:D3DVector; const b:D3DVector; const rate:Single):D3DVector;

//posからtargetに向いた長さkのベクトルを生成
function NowLookAt(const pos:D3DVECTOR; const target:D3DVECTOR; const k:Single):D3DVECTOR;

//正規化除算(x,y,z,w) → (x/w, y/w, z/w, 1/w)
function NowViewFrustumToScreen(const a:THgVector):THgVector;

//三つの頂点で張られる平面(三角形)の法線を返す。ただし、p1,p2,p3は法線に対して右回りと仮定する
function NowNormalFromPoints(const p1,p2,p3:D3DVector):D3DVector;

//ベクトルvを法線nなる平面にぶつけて、跳ね返った後のベクトルを得る(完全弾性衝突)
function NowReflection(const v:D3DVector; normal:D3Dvector):D3DVector; overload;

//ベクトルvを法線nなる平面にぶつけて、跳ね返った後のベクトルを得る(弾性係数e)
function NowReflection(const v:D3DVector; normal:D3Dvector; e:Single):D3DVector; overload;

//////////////////////3D演算(クォータニオン)

//足し算
function NowAddQuat(const a:TQuaternion; const b:TQuaternion):TQuaternion;

//引き算
function NowSubQuat(const a:TQuaternion; const b:TQuaternion):TQuaternion;

//乗算
function NowMulQuat(const a:TQuaternion; const b:TQuaternion):TQuaternion;

//内積
function NowDotQuat(const a:TQuaternion; const b:TQuaternion):Single;

//ノルム
function NowLengthQuat(const a:TQuaternion):Single;

//正規化
function NowNormalizeQuat(const a:TQuaternion):TQuaternion;

//axis回りの回転を表すクォータニオンを作る
function  NowQuatFromRotation(const axis:D3DVECTOR; const t:Integer):TQuaternion;

//クォータニオンから行列を作る
function  NowMatrixFromQuat(const a:TQuaternion):D3DMATRIX;

//axis回りの回転を表す行列を作る(QuatFromRotationとMatrixFromQuatの合成)
function NowMatrixFromRotation(const axis:D3DVECTOR; const t:Integer):D3DMATRIX;

//Slerp
function NowSlerpQuat(const a:TQuaternion; const b:TQuaternion; const rate:Single):TQuaternion;

//逆クォータニオン
function NowInvQuat( srcQuat: TQuaternion ): TQuaternion;

//普通の線形補間
function NowLerpQuat(const a:TQuaternion; const b:TQuaternion; const rate:Single):TQuaternion;

//////////////////////3D演算(マトリクス)

//ベクトルを変換(同次)
function NowTransform(const a:THgVector; const m:D3DMATRIX):THgVector; overload;
function NowTransform(const a:D3DVector; const m:D3DMATRIX):D3DVector; overload;

//一気にベクトルを変換
procedure NowTransformBatch(const src:Pointer; const m:D3DMATRIX; const num:Integer; const span:Integer);

//透視射影行列を作る(隠し関数)
function NowProjectFrustum(const halfWidth:Single; const halfHeight:Single; const nearZ:Single; const farZ:Single):D3DMATRIX;

//カメラとターゲットの位置から、ビュー行列を作る(隠し関数)
function NowViewMatrix(const from:D3DVECTOR; const at:D3DVECTOR; const world_up:D3DVECTOR; const roll:Integer):D3DMATRIX;

//姿勢・位置からビュー行列を作る(隠し関数)
function NowViewMatixFromOrientation(const from:D3DVECTOR; const dir:D3DVECTOR; const up:D3DVECTOR):D3DMATRIX;


//逆行列(姿勢・平行移動分を格納していると決めウチ、拡大・旋断は考慮しない)
function NowInvMatrix(const m:D3DMATRIX):D3DMATRIX;

//逆行列(回転だけを格納していると決めウチ)
function NowInvMatrixRot(const m:D3DMATRIX):D3DMATRIX;

//転置行列
function NowTMatrix(const m:D3DMATRIX):D3DMATRIX;

//回転部分だけ抜き出し
function NowExtractRotation(const m:D3DMATRIX):D3DMATRIX;

//平行移動部分だけ抜き出し
function NowExtractTranslation(const m:D3DMATRIX):D3DMATRIX;

//マトリクスの合成
function NowCompositeMatrix(const m:D3DMATRIX; const n:D3DMATRIX):D3DMATRIX;

//m.nの回転成分だけ抜き出して合成する
function NowCompositeRotMatrix(const m:D3DMATRIX; const n:D3DMATRIX):D3DMATRIX;




//回転行列
function NowRotXMatrix(const t:Integer):D3DMATRIX;
function NowRotYMatrix(const t:Integer):D3DMATRIX;
function NowRotZMatrix(const t:Integer):D3DMATRIX;

//スケーリング・平行移動行列
function NowScaleMatrix(const x,y,z:Single):D3DMATRIX;
function NowTranslateMatrix(const x,y,z:Single):D3DMATRIX;



//////////////////////3D演算(33マトリクス)

//乗算
function NowCompositeMatrix33(const a:TMatrix33; const b:TMatrix33):TMatrix33;


//逆行列
//valid には、逆行列が、ちゃんとあるか(det<>0)が入る
function NowInvMatrix33(const a:TMatrix33; var valid:Boolean):TMatrix33;

//一次変換
function NowTransform33(const a:D3DVector; const m:TMatrix33):D3DVector;

//行列式
function NowDeterminant33(const m:TMatrix33):Single;



//////////////////////3D演算(プリミティブ間距離の計算)

//2点の距離
function NowDistance(const a,b:D3DVECTOR):Single; overload;

//点と直線の距離
function NowDistance(const p:D3DVECTOR; const lin:TLine3D):Single; overload;

//点と線分の距離
function NowDistance(const p:D3DVECTOR; const lseg:TLineSeg3D):Single; overload;

//線分同士の距離
function NowDistance(const a,b:TLineSeg3D):Single; overload;

//点と平面の距離
function NowDistance(const a:D3Dvector; const plane:TPlane3D):Single; overload;

//////////////////////3D演算(包含・交差判定)

//pはboxに含まれるか？
function NowPointInBox(const p:D3DVector; const box:TBox3D):Boolean;

//aはbに包まれるか？
function NowBoxInBox(const a:TBox3D; const b:TBox3D):Boolean;

//pは 四面体 tri に含まれるか？
function NowPointInTetrahedron(const p:D3DVector; const t1,t2,t3,t4:D3DVector):Boolean;

//頂点が三角柱(triを底面とし、高さheight)にめり込んでいるか
function NowPointInPillar(const p:D3DVector; const tri:TFacedTriangle3D; const height:Single):Boolean;

//三角錐 tri と、boxは、交差ないし包含しあうか(ちょいと不正確)
//function NowIntersectTetrahedronBox(const tri:Array of D3DVector; const box:TBox3D):Boolean;

//線分abと平面sの交差判定 destに交差してる点を返す
function NowIntersectLinePlane(var dest:D3DVector; const seg:TLineSeg3D; const s:TPlane3D):Boolean;

//三角形と線分abの交差判定(あんまり自信ない(^^;))
function NowIntersectLineTriangle(const seg:TLineSeg3D; const t1:D3DVector; const t2:D3DVector; const t3:D3DVector):Boolean;

//ポリゴン同士の交差判定(あんまり自信ない(^^;))
function NowIntersectPolygons(const poly1:Array of D3DVector; const poly2:Array of D3DVector):Boolean;

//線分segとboxの交差判定
function NowIntersectLineBox(const seg:TLineSeg3D; const box:TBox3D):Boolean;

//四角錐とboxは交差しあうか？(四角錐を円錐で近似して計算、まだ不完全)
function NowIntersectPyramidBox(const top:D3DVector; const lt:D3DVector; const rt:D3DVector; const lb:D3DVector; const rb:D3DVector; const box:TBox3D):Boolean;

//boxは視界の中にあるか(viewは、boxの属する座標系にある頂点を、透視変換まで行う行列。clipの中に変換後のboxが入っていればTrue)…未実装
function NowIntersectViewBox(const view:D3DMatrix; const box:TBox3D; const clip:TRect):Boolean;

//カプセル同士は交差するか？
function NowIntersectCapsule(const a,b:TCapsule3D):Boolean;



//Sin,Cosテーブル(4096度系)
var
  NowSinT:packed Array[0..4095] of Single;
  NowCosT:packed Array[0..4095] of Single;
  NowSinCosT:packed Array[0..4095] of TVector2D;   //Sin Cos
  NowCosSinT:packed Array[0..4095] of TVector2D;   //Cos Sin
  NowCosNSinT:packed Array[0..4095] of TVector2D;  //Cos -Sin
  NowNSinCosT:packed Array[0..4095] of TVector2D;  //-Sin Cos


implementation

const
  PNowSinT:Pointer = @NowSinT;
  PNowCosT:Pointer = @NowCosT;
  PNowSinCosT:Pointer = @NowSinCosT;
  PNowCosSinT:Pointer = @NowCosSinT;
  PNowCosNSinT:Pointer = @NowCosNSinT;
  PNowNSinCosT:Pointer = @NowNSinCosT;



//テーブルや変数の初期化(Initializationで呼ばれる…けど)
procedure NowInit;
var
  s,c:Single;
  i:Integer;
  real_i:Single;
begin
  //WriteLn('テーブル作成ルーチンを通りました');
  //三角関数テーブルの作成
  For i:=0 to 4095 do begin
    real_i:=i;
    s:=Sin(real_i / 2048.0 * Pi);
    c:=Cos(real_i / 2048.0 * Pi);

    NowSinT[i]:=s;
    NowCosT[i]:=c;

    NowSinCosT[i]:=Vector2D(s,c);
    NowCosSinT[i]:=Vector2D(c,s);
    NowCosNSinT[i]:=Vector2D(c,-s);
    NowNSinCosT[i]:=Vector2D(-s,c);
  end;

end;



/////////////////////////////////////
// 各種構造体の初期化
/////////////////////////////////////
function Vector2D(_x,_y:Single):TVector2D;
begin
  result.X:=_x;
  result.Y:=_y;
end;

function Vector(_x,_y,_z:Single):D3DVECTOR;
begin
  result.X:=_x;
  result.Y:=_y;
  result.Z:=_z;
end;

function HgVector(_x,_y,_z,_w:Single):THgVector;
begin
  result.X:=_x;
  result.Y:=_y;
  result.Z:=_z;
  result.W:=_w;
end;


function Quaternion(_w,_x,_y,_z:Single):TQuaternion;
begin
  result.W:=_w;
  result.X:=_x;
  result.Y:=_y;
  result.Z:=_z;
end;


function  Sphere3D(_x,_y,_z,_r:Single):TSphere3D;
begin
  result.X:=_x;
  result.Y:=_y;
  result.Z:=_z;
  result.R:=_r;
end;

function  Plane3D(_a,_b,_c,_d:Single):TPlane3D;
var
  rlen:Single;
  v:D3DVector;
begin
  v:=Vector(_a,_b,_c);
  rlen:=1.0/NowLength(v);
  with result do begin
    Normal:=NowScale(v, rlen);
    D:=_d * rlen;
  end;
end;

function  Line3D(_p, _dir:D3DVector):TLine3D;
begin
  with result do begin
    P:=_p;
    Dir:=_dir;
  end;
end;

function  LineSeg3D(_start, _end:D3DVector):TLineSeg3D;
begin
  with result do begin
    A:=_start;
    B:=_end;
  end;
end;

///////////////////////////////////
// 色調
///////////////////////////////////

function NowScaleColor(const col:D3DCOLOR; const scale:Single):D3DCOLOR;
var
  scale2:DWord;
  r,g,b:Integer;
begin
  scale2:=Round(scale * 256.0);

  if scale <= 0.0 then begin
    Result:=0;
    exit;
  end;

  r:=col And $00FF0000 Shr 16 * scale2 Shr 8;
  if r < 0 then
    r:=0
  else if r > 255 then
    r:=255;

  g:=col And $0000FF00 Shr 8 * scale2 Shr 8;
  if g < 0 then
    g:=0
  else if g > 255 then
    g:=255;

  b:=col And $000000FF * scale2 Shr 8;
  if b < 0 then
    b:=0
  else if b > 255 then
    b:=255;


  result:=
    col And $FF000000 Or
    DWord(r) Shl 16 Or
    DWord(g) Shl 8 Or
    DWord(b);
end;

function NowAddColor(const col1:D3DCOLOR; const col2:D3DCOLOR):D3DCOLOR;
var
  c,d:DWord;
begin

  c:=col1 And $00FEFEFE + col2 And $00FEFEFE;
  d:=c;

  c:=c And $01010100;   //桁上がり検出マスク
  c:=c - c Shr 8;       //桁上がりをしたところは、0xffになる

  result := c Or d Or $FF000000;

end;

//色の線形補間 tはブレンド比率、0でcol1寄り、1でcol2寄り
function NowLerpColor(const col1:D3DCOLOR; const col2:D3DCOLOR; const t:Single):D3DCOLOR;
var
  scale1,scale2:DWord;
  a,r,g,b:Integer;
begin
  scale2:=Round(t * 256.0);
  scale1:=256 - scale2;

  a:=((col1 And $FF000000 Shr 24 * scale1)  +  (col2 And $FF000000 Shr 24 * scale2)) Shr 8;
  if a < 0 then
    a:=0
  else if a > 255 then
    a:=255;

  r:=((col1 And $00FF0000 Shr 16 * scale1)  +  (col2 And $00FF0000 Shr 16 * scale2)) Shr 8;
  if r < 0 then
    r:=0
  else if r > 255 then
    r:=255;

  g:=((col1 And $0000FF00 Shr 8 * scale1) + (col2 And $0000FF00 Shr 8 * scale2)) Shr 8;
  if g < 0 then
    g:=0
  else if g > 255 then
    g:=255;

  b:=((col1 And $000000FF * scale1) + (col2 And $000000FF * scale2)) Shr 8;
  if b < 0 then
    b:=0
  else if b > 255 then
    b:=255;

  result:=
    DWord(a) Shl 24 Or
    DWord(r) Shl 16 Or
    DWord(g) Shl 8 Or
    DWord(b);
end;

//HSBからRGB形式へ
function NowHSBToRGB(h,s,b:Single):D3DCOLOR;
var
  re,gr,bl:Single;
  ofs,phase,d1,d2,d3:Single;
begin

  //灰色
  if s = 0 then begin
    re:=b;
    gr:=b;
    bl:=b;
  end else begin

    if phase >= 0 then
      phase:=frac(h) * 6.0
    else
      phase:=(1.0 - frac(h)) * 6.0;

    ofs := frac(phase);
    d1 := b * ( 1.0 - s );
    d2 := b * ( 1.0 - s * ofs );
    d3 := b * ( 1.0 - s * (1.0 - ofs) );

    Case Trunc(phase) of

      0: begin
        re := b;
        gr := d3;
        bl := d1;
      end;

      1: begin
        re := d2;
        gr := b;
        bl := d1;
      end;

      2: begin
        re := d1;
        gr := b;
        bl := d3;
      end;

      3: begin
        re := d1;
        gr := d2;
        bl := b;
      end;

      4: begin
        re := d3;
        gr := d1;
        bl := b;
      end;

      5: begin
        re := b;
        gr := d1;
        bl := d2;
      end;

    end;
  end;

  result:=D3DCOLOR_COLORVALUE(re, gr, bl, 1.0);
end;


///////////////////////////////////
// 変換関数
///////////////////////////////////

//同次ベクトルと普通のベクトルの変換
function  NowHomoginize(const a:D3DVector):THgVector;
begin
  Move(a,result,Sizeof(D3DVector));
  result.W:=1.0;
end;

function  NowHeteroginize(const a:THgVector; _div:Boolean):D3DVector;
var
  d:Single;
begin
  if not _div then
    Move(a,result,Sizeof(D3DVector))
  else begin
    d:=1.0 / a.w;
    result.x := a.x * d;
    result.y := a.y * d;
    result.z := a.z * d;
  end;
end;




{$IFDEF USE_3DNOW}

/////////////////////////////////////
// Int <--> Float 変換
/////////////////////////////////////
function  NowItoF2D(const a:TPoint):TVector2D;
begin
  asm
    DW    FEMMS;

    //mov   eax,a;
    DW    MOVQ_ld;      DB MM0_mEAX;
    //mov   edx,@result;
    DW    _3DNow;       DB MM0_MM0;   DB PI2FD;
    DW    MOVQ_st;      DB MM0_mEDX;

    DW    FEMMS;
  end;
end;

function  NowFtoI2D(const a:TVector2D):TPoint;
begin
  asm
    DW  FEMMS;

    //mov   eax,a;
    DW    MOVQ_ld;      DB MM0_mEAX;
    //mov   edx,@result;
    DW    _3DNow;       DB MM0_MM0;   DB PF2ID;
    DW    MOVQ_st;      DB MM0_mEDX;

    DW  FEMMS;
  end;
end;


//四捨五入
function  NowRound2D(const a:TVector2D):TPoint;
begin
  asm
    DW  FEMMS;

    mov   ecx,$3F000000;  //0.5
    //mov   eax,a;
    DW    MOVQ_ld;      DB MM0_mEAX;
    DW    MOVD_ld;      DB MM1_ECX;
    DW    PUNPCKLDQ;    DB MM1_MM1;
    DW    _3DNow;       DB MM0_MM1;   DB PFADD;
    DW    _3DNow;       DB MM0_MM0;   DB PF2ID;
    DW    MOVQ_st;      DB MM0_mEDX;

    DW  FEMMS;
  end;
end;


//LVからTLVへ
function  NowLVToTLV(const lv:D3DLVERTEX; const transMat:D3DMATRIX):D3DTLVERTEX;
var
  hg:THgVector;  //頂点座標
begin
  Move(lv,hg,Sizeof(D3DVERTEX));
  hg.W:=1;

  hg:=NowViewFrustumToScreen(NowTransform(hg,transMat));
  Move(hg,result,sizeof(THgVector));
  Move(lv.color,result.color,sizeof(D3DTLVERTEX) - sizeof(THgVector));

end;

//LVからTLVへ、バッチ変換
procedure NowLVToTLVBatch(src:LPD3DLVERTEX; dest:LPD3DTLVERTEX; const transMat:D3DMATRIX; count:Integer);
begin
  asm
    DW      FEMMS;
    push    ebx
    push    edi


    mov     eax,src;
    mov     ebx,dest
    mov     ecx,count
    mov     edx,transMat
    mov     edi,$3F800000;    //1.0

    and     ecx,ecx;
    jz      @END



    @TOP:
    test  ecx,4;
    jnz   @NEXT;
    DW    PREFETCHW;  DB mm1_bEAX;  DB 128;
    DW    PREFETCHW;  DB mm1_bEBX;  DB 128;
    @NEXT:

    //まず、変換
    DW    movq_ld;    DB mm0_mEAX;		          // ay, ax	in mm0
    DW    movd_ld;    DB mm2_edi;
    DW    movd_ld;    DB mm1_bEAX;  DB 8;
    DW    punpckldq;  DB mm1_mm2;               //  1, az	in mm1

    //
    DW    movq_ld;    DB mm4_mEDX;      		    // m01, m00	in mm4
    DW    punpckhdq;	DB mm2_mm0;		            // ay, undef	in mm2
    //
    DW    movq_ld;    DB mm5_bEDX;  DB 16;	    // m11, m10	in mm5
    DW    punpckldq;  DB mm0_mm0;       		    // ax, ax	in mm2
    //
    DW    movq_ld;    DB mm6_bEDX;  DB 32;      // m21, m20	in mm6
    DW    _3DNow;     DB mm4_mm0;   DB pfmul;   // ax*m01, ax*m00	in mm4
    //
    DW    movq_ld;    DB mm7_bEDX;  DB 48;      // m31, m30	in mm7
    DW    punpckhdq;  DB mm2_mm2;           		// ay, ay	in mm2
    //
    DW    punpckhdq;  DB mm3_mm1;           		// aw, undef	in mm3
    DW    _3DNow;     DB mm5_mm2;   DB pfmul; 	// ay*m11, ay*m10
    //
    DW    punpckldq;  DB mm1_mm1;           		// az, az	in mm3
    DW    _3DNow;     DB mm0_bEDX;  DB 8;   DB pfmul;   	// ax*m03, ax*m02 in mm0
    //
    DW    punpckhdq;  DB mm3_mm3;            		// aw, aw	in mm3
    DW    _3DNow;     DB mm2_bEDX;  DB 24;  DB PFMUL;     // ay*13, ay*m12	in mm2
    //
    DW    _3DNow;     DB mm6_mm1;   DB pfmul; 	// az*m21, az*m20	in mm6
    DW    _3DNow;     DB mm5_mm4;   DB pfadd;
    //
    DW    _3DNow;     DB mm1_bEDX;  DB 40;    DB pfmul    // az*23, az*22	in mm1
    DW    _3DNow;     DB mm2_mm0;   DB pfadd;
    //
    DW    _3DNow;     DB mm7_mm3;   DB pfmul;   // aw*m31, aw*m30	in mm7
    DW    _3DNow;     DB mm6_mm5;   DB pfadd;
    //
    DW    _3DNow;     DB mm3_bEDX;  DB 56;    DB pfmul;   // aw*m33, aw*m32	in mm3
    DW    _3DNow;     DB mm2_mm1;   DB pfadd;
    //
    DW    _3DNow;     DB mm7_mm6;   DB pfadd    // ry, rx	in mm7
    DW    _3DNow;     DB mm3_mm2;   DB pfadd		// rw,rz		in mm3

    //正規化除算
    DW    MOVQ_ld;    DB mm4_mm3;
    DW    PSXXimq;    DB mm4_SRLim; DB 32;
    DW    _3DNow;     DB mm4_mm4;   DB pfrcp;   // 1/rw  | 1/rw
    DW    _3DNow;     DB mm7_mm4;   DB pfmul;   // ry/rw | rx/rw
    DW    _3DNow;     DB mm3_mm4;   DB pfmul;   // rw/rw | rz/rw
    DW    punpckldq;  DB mm3_mm4;               // 1/rw  | rz/rw

    //格納


    //転送
    DW    MOVQ_ld;    DB mm0_bEAX;  DB 16;
    DW    MOVQ_ld;    DB mm1_bEAX;  DB 24;
    DW    MOVQ_st;    DB mm7_mEBX;
    DW    MOVQ_st;    DB mm3_bEBX;  DB 8;
    DW    MOVQ_st;    DB mm0_bEBX;  DB 16;
    DW    MOVQ_st;    DB mm1_bEBX;  DB 24;

    add   eax,32;
    add   ebx,32;


    dec   ecx
    jnz   @TOP

    @END:

    pop     edi
    pop     ebx
    DW      FEMMS;
  end;
end;


/////////////////////////////////////
// 1D計算パッケージ
/////////////////////////////////////

//3次Bスプライン補間(4点)
//4点のうち、中間2点のブレンド比を与えて補間。rateは0〜1で、大きいほど出力はp3に近づく
procedure NowInterpolate(var dest:Single; const p1,p2,p3,p4:Single; const rate:Single);
begin
  asm
    DW    FEMMS;

    //まず、ブレンド比の計算
    mov   eax,rate;
    DW    MOVD_ld;    DB MM0_EAX;                   //MM0 = rate
    mov   ecx,$3F800000;
    DW    MOVD_ld;    DB MM7_ECX;                   //MM7 = 1.0;
    DW    MOVQ_ld;    DB MM3_MM7;
    DW    _3DNow;     DB MM3_MM0;   DB PFSUB;
    DW    MOVQ_ld;    DB MM4_MM3;
    DW    _3DNow;     DB MM3_MM4;   DB PFMUL;
    DW    _3DNow;     DB MM3_MM4;   DB PFMUL;       //mm3 = (1-rate)^3
    DW    MOVQ_ld;    DB MM1_MM0;
    DW    _3DNow;     DB MM1_MM1;   DB PFMUL;       //MM1 = rate^2
    DW    MOVQ_ld;    DB MM2_MM1;
    DW    _3DNow;     DB MM2_MM0;   DB PFMUL;       //MM2 = rate^3

    mov   ecx,$40400000;    //3.0
    DW    MOVQ_ld;    DB MM4_MM1;
    DW    MOVQ_ld;    DB MM5_MM2;
    DW    _3DNow;     DB MM4_MM4;   DB PFADD;
    mov   edx,$40800000     //4.0
    DW    _3DNow;     DB MM5_MM4;   DB PFSUB;
    DW    MOVD_ld;    DB MM4_ECX;
    DW    _3DNow;     DB MM4_MM5;   DB PFMUL;
    DW    MOVD_ld;    DB MM5_EDX;
    DW    _3DNow;     DB MM4_MM5;   DB PFADD;       //MM4 = 3t^3 - 6t^2 + 4

    DW    MOVQ_ld;    DB MM5_MM0;
    DW    _3DNow;     DB MM5_MM1;   DB PFADD;
    DW    _3DNow;     DB MM5_MM2;   DB PFSUB;
    DW    MOVD_ld;    DB MM6_ECX;
    DW    _3DNow;     DB MM5_MM6;   DB PFMUL;
    DW    _3DNow;     DB MM5_MM7;   DB PFADD;       //MM5 = 3(-t^3+t^2+t)+1

    DW    PUNPCKLDQ;  DB MM3_MM4;   // MM3 = r1 | r2
    DW    PUNPCKLDQ;  DB MM5_MM2;
    DW    MOVQ_ld;    DB MM0_MM5;   // MM0 = r3 | r4

    mov   eax,p1;
    mov   ecx,p2;
    DW    MOVD_ld;    DB MM2_EAX;
    DW    MOVD_ld;    DB MM4_ECX;
    DW    PUNPCKLDQ;  DB MM2_MM4;   // MM2 = p1 | p2

    mov   eax,p3;
    mov   ecx,p4;
    DW    MOVD_ld;    DB MM1_EAX;
    DW    MOVD_ld;    DB MM5_ECX;
    DW    PUNPCKLDQ;  DB MM1_MM5;   // MM1 = p3 | p4

    DW    _3DNow;     DB MM3_MM2;   DB PFMUL;
    DW    _3DNow;     DB MM0_MM1;   DB PFMUL;

    mov   ecx,$3E2AAAAB;    //1/6
    DW    _3DNow;     DB MM0_MM3;   DB PFADD;
    DW    MOVD_ld;    DB MM1_ECX;
    DW    _3DNow;     DB MM0_MM0;   DB PFACC;
    DW    _3DNow;     DB MM0_MM1;   DB PFMUL;

    mov   eax,dest;
    DW    MOVD_st;    DB MM0_mEAX;


    DW    FEMMS;
  end;

end;

//線形補間(rateは２点のブレンド比で0〜1、1に近いほどp2が濃くなる。
procedure NowLerp(var dest:Single; const p1,p2:Single; const rate:Single);
begin
  asm
    push  ebx;

    DW    FEMMS;

    mov   eax,dest;
    mov   edx,p1
    mov   ecx,p2

    DW    MOVD_ld;    DB MM0_EDX;   //p1
    DW    MOVD_ld;    DB MM1_ECX;   //p2;
    DW    PUNPCKLDQ;  DB MM0_MM1;   //MM0 = p1 | p2

    mov   ecx,$3F800000;  //1.0
    mov   ebx,rate;
    DW    MOVD_ld;    DB MM7_ECX;   //MM7 = 1.0 | 0.0
    DW    MOVD_ld;    DB MM6_EBX;   //MM6 = rate
    DW    MOVD_ld;    DB MM5_EBX;
    DW    _3DNow;     DB MM7_MM5;   DB PFSUB;   //MM7 = 1.0-rate | 0
    DW    PUNPCKLDQ;  DB MM7_MM6;   //MM7 = 1.0-rate | rate

    DW    _3DNow;     DB MM0_MM7;   DB PFMUL;
    DW    _3DNow;     DB MM0_MM0;   DB PFACC;

    DW    MOVD_st;    DB MM0_mEAX;

    DW    FEMMS;
    pop   ebx;
  end;
end;



//Math入ってるArcTan2みたいなもの？ x=0でもOK、値域は 0 <= result < Pi
function NowArcTan2(Y, X: Extended): Extended;
var
  rad:Extended;
begin
  if x <> 0 then begin
    rad := ArcTan(-y/x);

    if x < 0 then
      rad:=rad + Pi;
    else if y > 0 then
      rad:=rad + 2.0*Pi;

  end else begin
    if y > 0 then
      rad:=Pi * 0.5;
    else
      rad:=Pi * 1.5;
  end;
{
  if x <> 0 then begin
    rad := ArcTan(-y/x);

    if x < 0 then
      rad:=rad + Pi
    else if y > 0 then
      rad:=rad + 2.0*Pi;

  end else begin
    if y < 0 then
      rad:=Pi/4.0
    else
      rad:=Pi * (3.0/4.0);
  end;
}

  result:=rad;
end;




/////////////////////////////////////
// 2D計算パッケージ
/////////////////////////////////////

//ベクトルの足し算
function NowAdd2D(const a:TVector2D; const b:TVector2D):TVector2D;
begin
  asm
    DW    FEMMS;

    DW    MOVQ_ld;    DB MM0_mEAX;
    DW    MOVQ_ld;    DB MM1_mEDX;

    DW    _3DNow;     DB MM0_MM1;   DB PFADD;

    DW    MOVQ_st;    DB MM0_mECX;

    DW    FEMMS;
  end;

end;



//ベクトルの引き算
function NowSub2D(const a:TVector2D; const b:TVector2D):TVector2D;
begin
  asm
    DW    FEMMS;

    DW    MOVQ_ld;    DB MM0_mEAX;
    DW    MOVQ_ld;    DB MM1_mEDX;

    DW    _3DNow;     DB MM0_MM1;   DB PFSUB;

    DW    MOVQ_st;    DB MM0_mECX;

    DW    FEMMS;
  end;
end;


//内積
procedure NowDotProduct2D(var dest:Single; const a:TVector2D; const b:TVector2D);
begin

  asm
    DW    FEMMS;

    DW    MOVQ_ld;    DB MM0_mEDX;
    DW    MOVQ_ld;    DB MM1_mECX;

    DW    _3DNow;     DB MM0_MM1;   DB PFMUL;
    DW    _3DNow;     DB MM0_MM0;   DB PFACC;

    DW    MOVD_st;    DB MM0_mEAX;

    DW    FEMMS;
  end;

end;


//2点の距離
procedure NowDistance2D(var dest:Single; const a,b:TVector2D);
begin
  asm

    DW    FEMMS

    //読みこみ
    DW    MOVQ_ld;  DB    MM0_mEDX;
    DW    MOVQ_ld;  DB    MM1_mECX;

    //減算
    DW    _3DNow;   DB    MM0_MM1;    DB    PFSUB

    //自乗
    DW    _3DNow;   DB    MM0_MM0;    DB    PFMUL

    //上位と下位の加算
    DW    _3DNow;   DB    MM0_MM0;    DB    PFACC

    //平方根
    DW    _3DNow;   DB    MM1_MM0;    DB    PFRSQRT
    DW    _3DNow;   DB    MM0_MM1;    DB    PFMUL

    //結果の格納
    DW    MOVD_st;  DB    MM0_mEAX;


    DW    FEMMS
  end;
end;




//長さ
procedure NowLength2D(var dest:Single; const a:TVector2D);
begin
  asm
    DW    FEMMS;

    DW    MOVQ_ld;    DB MM0_mEDX;
    DW    _3DNow;     DB MM0_MM0;   DB PFMUL;
    DW    _3DNow;     DB MM0_MM0;   DB PFACC;
    DW    _3DNow;     DB MM1_MM0;   DB PFRSQRT;
    DW    _3DNow;     DB MM0_MM1;   DB PFMUL;

    DW    MOVD_st;    DB MM0_mEAX;

    DW    FEMMS;
  end;
end;


//スケーリング
function NowScale2D(const a:TVector2D; const scale:Single):TVector2D;
begin
  asm
    DW    FEMMS;

    mov   eax,a;
    mov   edx,[scale];
    DW    MOVD_ld;    DB MM1_EDX;
    DW    MOVQ_ld;    DB MM0_mEAX;

    DW    PUNPCKLDQ;  DB MM1_MM1;

    DW    _3DNow;     DB MM0_MM1;   DB  PFMUL;

    mov   ecx,@result;
    DW    MOVQ_st;    DB MM0_mECX;

    DW    FEMMS;
  end;
end;


//ベクトルの向きを4096度系で
procedure NowVecToAngle2D(var dest:Integer; const a:TVector2D);
var
  v:TVector2D;
  pv:PVector2D;
  sg:Integer;   //象限
  mid:Integer;  //2分探索用
begin
  v:=a;
  pv:=@v;

  //象限はどれ？
  {
  if (v.y > 0) then begin
    if v.x > 0 then
      sg:=1
    else begin
      sg:=2;
      v.x:=-v.x;
    end;
  end else if (v.y < 0) then begin
    if v.x >= 0 then
      sg:=4
    else begin
      sg:=3;
      v.x:=-v.x;
    end;
  end else begin
    if v.x >= 0 then
      sg:=1
    else begin
      sg:=3;
      v.x:=-v.x;
    end;
  end;
  }
  

  asm
    push  eax;
    push  ebx;
    push  edx;
    push  edi;
    DW    FEMMS;

    mov   edx,pv;

    //象限はどれ?
    DW    PXOR;       DB MM1_MM1;
    DW    MOVQ_ld;    DB MM0_mEDX;
    DW    MOVQ_ld;    DB MM7_MM0;

    DW    MOVQ_ld;    DB MM2_MM0;
    DW    _3DNow;     DB MM0_MM1;   DB PFCMPGE;   //v.x >= 0 ?
    DW    MOVD_st;    DB MM0_EAX;
    and   eax,eax;
    jnz   @PLUS_X;
    // x < 0
    mov    sg,2;
    DW    _3DNow;     DB MM7_MM1;   DB PFSUBR;    //X軸反転
    jmp   @JUDGE_Y;

    // x >= 0
    @PLUS_X:
    mov    sg,1;


    @JUDGE_Y:
    DW    PUNPCKHDQ;  DB MM2_MM2;
    DW    _3DNow;     DB MM2_MM1;   DB PFCMPGE;   //v.y >= 0 ?
    DW    MOVD_st;    DB MM2_EAX;
    and   eax,eax;
    jnz   @PLUS_Y;
    // y < 0
    add   sg,2;
    @PLUS_Y:



    //ノーマライズ

    DW    MOVQ_ld;  DB MM0_MM7;
    DW    MOVQ_ld;  DB MM2_MM0;
    DW    _3DNow;   DB MM0_MM0;   DB PFMUL;
    DW    _3DNow;   DB MM0_MM0;   DB PFACC;     //2乗和
    DW    MOVQ_ld;  DB MM4_MM0;
    DW    _3DNow;   DB MM0_MM0;   DB PFRSQRT;   //これで、1/kがでる

    DW    MOVQ_ld;  DB MM1_MM0;
    DW    _3DNow;   DB MM1_MM1;   DB PFMUL;
    DW    _3DNow;   DB MM4_MM1;   DB PFRSQIT1;
    DW    _3DNow;   DB MM4_MM0;   DB PFRCPIT2;

    DW    _3DNow;   DB MM4_MM2;   DB PFMUL;     //元の値と掛け算


    //Cosテーブル内を2分探索
    mov   ecx,1023;   //hi
    mov   ebx,0;      //lo
    mov   edx,511;    //mid


    @LOOP_TOP:

    DW    MOVQ_ld;  DB MM2_MM4;
    cmp   ecx,ebx;
    jle   @LOOP_END

    mov   eax,edx;
    shl   eax,2;
    add   eax,PNowCosT;

    DW    _3DNow;   DB MM2_mEAX;   DB PFCMPGE;   //v.x >= NowCosT[mid] ?
    DW    MOVD_st;  DB MM2_EDI;
    and   edi,edi
    jz    @PATH1 //false then ...

    mov   ecx,edx;  //hi = mid
    add   edx,ebx;
    shr   edx,1;    //mid = (hi + lo) / 2
    jmp   @LOOP_TOP;

    @PATH1:
    mov   ebx,edx;  //lo = mid + 1 , mid = (hi+lo)/2
    add   edx,ecx;
    shr   edx,1;
    inc   ebx;
    jmp   @LOOP_TOP;

    @LOOP_END:


    DW    FEMMS;

    mov   mid,edx;

    pop   edi;
    pop   edx;
    pop   ebx;
    pop   eax;


  end;

  //象限に応じて、結果を加工
  Case sg of
    1: dest:=mid;
    2: dest:=2047-mid;
    4: dest:=2048+mid;
    3: dest:=4095-mid;
  end;

end;

//回転
function NowRot2D(const a:TVector2D; const angle:Integer):TVector2D;
begin
  asm
    DW  FEMMS

    mov   eax,[angle];
    and   eax,4095;   //角度をノーマライズ
    mov   ecx,a;


    //読みこみ
    DW    MOVQ_ld;  DB MM2_mECX;
    DW    MOVQ_ld;  DB MM0_MM2;

    //MM2に、result.Xを求める
    mov   edx,PNowCosNSinT;
    shl   eax,3;
    add   edx,eax;
    DW    MOVQ_ld;  DB MM3_mEDX;
    DW    _3DNow;   DB MM2_MM3;   DB PFMUL;
    DW    _3DNow;   DB MM2_MM2;   DB PFACC;

    //MM0に、result.Yを求める
    mov   edx,PNowSinCosT;
    add   edx,eax;
    DW    MOVQ_ld;  DB MM3_mEDX;
    DW    _3DNow;   DB MM0_MM3;   DB PFMUL;
    DW    _3DNow;   DB MM0_MM0;   DB PFACC;

    //destに格納して終わり
    DW    PUNPCKHDQ;  DB MM2_MM0;
    mov   eax,@result;
    DW    MOVQ_st;    DB MM2_mEAX;

    DW  FEMMS

  end;
end;



//正規化
function NowNormalize2D(const a:TVector2D):TVector2D;
begin
  asm
    DW    FEMMS;

    DW    MOVQ_ld;  DB MM0_mEAX;
    DW    MOVQ_ld;  DB MM2_MM0;
    DW    _3DNow;   DB MM0_MM0;   DB PFMUL;
    DW    _3DNow;   DB MM0_MM0;   DB PFACC;     //2乗和
    DW    _3DNow;   DB MM0_MM0;   DB PFRSQRT;   //これで、1/kがでる
    DW    _3DNow;   DB MM0_MM2;   DB PFMUL;     //元の値と掛け算

    //結果の格納
    DW    MOVQ_st;  DB MM0_mEDX;

    DW    FEMMS;
  end;
end;


//ベクトルの長さをkにする
function NowChangeLength2D(const a:TVector2D; const k:Single):TVector2D;
begin
  asm
    DW    FEMMS;

    mov   ecx,a;
    mov   eax,[k];

    DW    MOVQ_ld;  DB MM0_mECX;
    DW    MOVQ_ld;  DB MM2_MM0;

    DW    _3DNow;   DB MM0_MM0;   DB PFMUL;
    DW    MOVD_ld;  DB MM1_EAX;

    DW    _3DNow;   DB MM0_MM0;   DB PFACC;     //2乗和
    DW    PUNPCKLDQ; DB MM1_MM1;

    DW    _3DNow;   DB MM0_MM0;   DB PFRSQRT;   //これで、1/長さ がでる
    DW    _3DNow;   DB MM0_MM2;   DB PFMUL;     //元の値と掛け算(正規化)
    DW    _3DNow;   DB MM0_MM1;   DB PFMUL;

    //結果の格納
    mov   eax,@result;
    DW    MOVQ_st;  DB MM0_mEAX;

    DW    FEMMS;

  end;
end;


//posからtargetに向いた長さkのベクトルを生成
function NowLookAt2D(const pos:TVector2D; const target:TVector2D; const k:Single):TVector2D;
begin
  asm
    DW    FEMMS;

    mov   eax,[k]
    mov   edx,pos;
    mov   ecx,target;
    DW    MOVQ_ld;    DB MM0_mEDX;
    DW    MOVQ_ld;    DB MM1_mECX;

    DW    _3DNow;     DB MM0_MM1;   DB PFSUBR;  //ターゲットを向いたベクトル

    DW    MOVQ_ld;    DB MM2_MM0;

    DW    _3DNow;     DB MM0_MM0;   DB PFMUL;
    DW    MOVD_ld;    DB MM1_EAX;

    DW    _3DNow;     DB MM0_MM0;   DB PFACC;     //2乗和
    DW    PUNPCKLDQ;  DB MM1_MM1;

    DW    _3DNow;     DB MM0_MM0;   DB PFRSQRT;   //これで、1/長さ がでる
    DW    _3DNow;     DB MM0_MM2;   DB PFMUL;     //元の値と掛け算(正規化)
    DW    _3DNow;     DB MM0_MM1;   DB PFMUL;

    //結果の格納
    mov   eax,@result;
    DW    MOVQ_st;    DB MM0_mEAX;

    DW    FEMMS;
  end;
end;


//長さkのランダムなベクトルを作る
function NowRandom2D(const k:Single):TVector2D;
begin
  result:=NowChangeLength2D(Vector2D(Random,Random),k);
end;


/////////////////////////////////////
// 3D計算パッケージ(ベクトル)
/////////////////////////////////////

//足し算
function NowAdd(const a:D3DVECTOR; const b:D3DVECTOR):D3DVECTOR;
begin
  asm
    DW    FEMMS;

    DW    MOVQ_ld;    DB MM0_mEAX;
    DW    MOVD_ld;    DB MM2_bEAX;  DB 8;

    DW    MOVQ_ld;    DB MM1_mEDX;
    DW    MOVD_ld;    DB MM3_bEDX;  DB 8;

    DW    _3DNow;     DB MM0_MM1;   DB PFADD;
    DW    _3DNow;     DB MM2_MM3;   DB PFADD;

    DW    MOVQ_st;    DB MM0_mECX;
    DW    MOVD_st;    DB MM2_bECX;  DB 8;

    DW    FEMMS;
  end;
end;

function NowAddHg(const a:THgVector; const b:THgVector):THgVector;
begin
  asm
    DW    FEMMS;

    DW    MOVQ_ld;    DB MM0_mEAX;
    DW    MOVQ_ld;    DB MM2_bEAX;  DB 8;

    DW    MOVQ_ld;    DB MM1_mEDX;
    DW    MOVQ_ld;    DB MM3_bEDX;  DB 8;

    DW    _3DNow;     DB MM0_MM1;   DB PFADD;
    DW    _3DNow;     DB MM2_MM3;   DB PFADD;

    DW    MOVQ_st;    DB MM0_mECX;
    DW    MOVQ_st;    DB MM2_bECX;  DB 8;

    DW    FEMMS;
  end;
end;


//引き算
function NowSub(const a:D3DVECTOR; const b:D3DVECTOR):D3DVECTOR;
begin
  asm
    DW    FEMMS;

    DW    MOVQ_ld;    DB MM0_mEAX;
    DW    MOVD_ld;    DB MM2_bEAX;  DB 8;

    DW    MOVQ_ld;    DB MM1_mEDX;
    DW    MOVD_ld;    DB MM3_bEDX;  DB 8;

    DW    _3DNow;     DB MM0_MM1;   DB PFSUB;
    DW    _3DNow;     DB MM2_MM3;   DB PFSUB;

    DW    MOVQ_st;    DB MM0_mECX;
    DW    MOVD_st;    DB MM2_bECX;  DB 8;

    DW    FEMMS;
  end;
end;

function NowSubHg(const a:THgVector; const b:THgVector):THgVector;
begin
  asm
    DW    FEMMS;

    DW    MOVQ_ld;    DB MM0_mEAX;
    DW    MOVQ_ld;    DB MM2_bEAX;  DB 8;

    DW    MOVQ_ld;    DB MM1_mEDX;
    DW    MOVQ_ld;    DB MM3_bEDX;  DB 8;

    DW    _3DNow;     DB MM0_MM1;   DB PFSUB;
    DW    _3DNow;     DB MM2_MM3;   DB PFSUB;

    DW    MOVQ_st;    DB MM0_mECX;
    DW    MOVQ_st;    DB MM2_bECX;  DB 8;

    DW    FEMMS;
  end;
end;


//内積
procedure NowDotProduct(var dest:Single; const a:D3DVECTOR; const b:D3DVECTOR);
begin
  asm
    DW    FEMMS;

    DW    MOVQ_ld;    DB MM0_mEDX;  //a
    DW    MOVQ_ld;    DB MM1_mECX;  //b
    DW    MOVD_ld;    DB MM2_bEDX;  DB 8;   //a.z
    DW    MOVD_ld;    DB MM3_bECX;  DB 8;   //b.z

    DW    _3DNow;     DB MM0_MM1;   DB PFMUL;
    DW    _3DNow;     DB MM0_MM0;   DB PFACC;
    DW    _3DNow;     DB MM2_MM3;   DB PFMUL;
    DW    _3DNow;     DB MM2_MM2;   DB PFACC;

    DW    _3DNow;     DB MM0_MM2;   DB PFADD;

    DW    MOVD_st;    DB MM0_mEAX;

    DW    FEMMS;
  end;
end;

//外積…AMDのSDKより引用
function NowCrossProduct(const a:D3DVECTOR; const b:D3DVECTOR):D3DVECTOR;
begin
  asm
    DW    FEMMS;

    DW    movd_ld;    DB mm0_bEDX;  DB 8;       //bz
    DW    movd_ld;    DB mm1_bEAX;  DB 8;       //az
    DW    movd_ld;    DB mm3_bEAX;  DB 4;       //ay
    DW    punpckldq;	DB MM0_MM0;               //b_z,b_z	in	mm0
    DW    _3DNow;     DB MM0_mEAX;  DB pfmul;   //ay*bz,ax*bz
    DW    punpckldq;  DB mm1_mm1;               //az,a_z	in	mm1
    DW    _3DNow;     DB mm1_mEDX;  DB pfmul;		//by*az,bx*az
    DW    movd_ld;    DB mm4_bEDX;  DB 4;     	//by
    DW    movq_ld;    DB mm2_mm0
    DW    _3DNow;     DB mm3_mEDX;  DB pfmul;		//undef,bx*ay
    DW    _3DNow;     DB mm0_MM1;   DB pfsub;		//x,undef
    DW    _3DNow;     DB mm4_mEAX;  DB pfmul;   //undef,ax*by
    DW    _3DNow;     DB mm1_mm2;   DB pfsub;		//undef,y
    DW    psXXimq;    DB mm0_SRLim; DB 32;
    DW    _3DNow;     DB mm4_mm3;   DB pfsub;
    DW    MOVD_st;    DB mm0_mECX;
    DW    MOVD_st;    Db mm1_bECX;  DB 4;
    DW    MOVD_st;    DB mm4_bECX;  DB 8;

    DW    FEMMS;
  end;
end;


//2点の距離
procedure NowDistance(var dest:Single;  const a,b:D3DVECTOR);
begin
  asm
    DW    FEMMS;

    DW    MOVQ_ld;    DB MM0_mEDX;
    DW    MOVD_ld;    DB MM1_bEDX;  DB 8;
    DW    MOVQ_ld;    DB MM2_mECX;
    DW    MOVD_ld;    DB MM3_bECX;  DB 8;

    DW    _3DNow;     DB MM0_MM2;   DB PFSUB;
    DW    _3DNow;     DB MM1_MM3;   DB PFSUB;
    DW    _3DNow;     DB MM0_MM0;   DB PFMUL;
    DW    _3DNow;     DB MM1_MM1;   DB PFMUL;

    DW    _3DNow;     DB MM0_MM0;   DB PFACC;
    DW    _3DNow;     DB MM1_MM1;   DB PFACC;
    DW    _3DNow;     DB MM0_MM1;   DB PFADD;

    DW    _3DNow;     DB MM1_MM0;   DB PFRSQRT;
    DW    _3Dnow;     DB MM0_MM1;   DB PFMUL;

    DW    MOVD_st;    DB MM0_mEAX;

    DW    FEMMS;
  end;
end;


//原点回りの回転(各軸を左回りに)
function NowRotX(const a:D3DVECTOR; const t:Integer):D3DVECTOR;
begin
  asm
    DW  FEMMS;

    // X' = X
    // Y'=  YcosT + ZsinT
    // Z' = YnsinT + ZcosT    //nsinは-sin
    and   edx,4095;
    shl   edx,3;

    DW    MOVD_ld;    DB MM7_mEAX;            //xを読みこむ = x'
    DW    MOVQ_ld;    DB MM0_bEAX;    DB 4;   //y,zを読み込む

    mov   eax,PNowCosSinT;
    add   eax,edx;
    DW    MOVQ_ld;    DB MM1_MM0;

    DW    _3DNow;     DB MM0_mEAX;    DB PFMUL;
    DW    _3DNow;     DB MM0_MM0;     DB PFACC;   //Y'

    mov   eax,PNowNSinCosT;
    add   eax,edx;
    DW    _3DNow;     DB MM1_mEAX;    DB PFMUL;
    DW    _3DNow;     DB MM1_MM1;     DB PFACC;   //Z'


    //結果を入れる
    DW    PUNPCKLDQ;  DB MM7_MM0;   //MM7  X' Y'
    DW    MOVQ_st;    DB MM7_mECX;
    DW    MOVD_st;    DB MM1_bECX;    DB 8; // Z'


    DW  FEMMS;
  end;
end;

function NowRotY(const a:D3DVECTOR; const t:Integer):D3DVECTOR;
begin
  asm
    DW  FEMMS;

    // X' = XcosT + ZsinT
    // Y'=  Y
    // Z' = XnsinT + ZcosT    //nsinは-sin
    and   edx,4095;
    shl   edx,3;

    DW    MOVD_ld;    DB MM7_bEAX;    DB 4;   //yを読みこむ = y'
    DW    MOVD_ld;    DB MM0_mEAX;
    DW    MOVD_ld;    DB MM1_bEAX;    DB 8;
    DW    PUNPCKLDQ;  DB MM0_MM1;             //x,zを読み込む

    mov   eax,PNowCosNSinT;
    add   eax,edx;
    DW    MOVQ_ld;    DB MM1_MM0;

    DW    _3DNow;     DB MM0_mEAX;    DB PFMUL;
    DW    _3DNow;     DB MM0_MM0;     DB PFACC;   //X'

    mov   eax,PNowSinCosT;
    add   eax,edx;
    DW    _3DNow;     DB MM1_mEAX;    DB PFMUL;
    DW    _3DNow;     DB MM1_MM1;     DB PFACC;   //Z'


    //結果を入れる
    DW    PUNPCKLDQ;  DB MM0_MM7;   //MM7  X' Y'
    DW    MOVQ_st;    DB MM0_mECX;
    DW    MOVD_st;    DB MM1_bECX;    DB 8;


    DW  FEMMS;
  end;
end;

function NowRotZ(const a:D3DVECTOR; const t:Integer):D3DVECTOR;
begin
  asm
    DW  FEMMS;

    // X' = XcosT + YsinT
    // Y' = XnsinT + YcosT
    // Z' = Z
    and   edx,4095;
    shl   edx,3;

    DW    MOVD_ld;    DB MM7_bEAX;    DB 8;   //zを読みこむ = z'
    DW    MOVQ_ld;    DB MM0_mEAX;            //x,yを読み込む

    mov   eax,PNowCosSinT;
    add   eax,edx;
    DW    MOVQ_ld;    DB MM1_MM0;

    DW    _3DNow;     DB MM0_mEAX;    DB PFMUL;
    DW    _3DNow;     DB MM0_MM0;     DB PFACC;   //X'

    mov   eax,PNowNSinCosT;
    add   eax,edx;
    DW    _3DNow;     DB MM1_mEAX;    DB PFMUL;
    DW    _3DNow;     DB MM1_MM1;     DB PFACC;   //Y'


    //結果を入れる
    DW    PUNPCKLDQ;  DB MM0_MM1;   //MM0  X' Y'
    DW    MOVQ_st;    DB MM0_mECX;
    DW    MOVD_st;    DB MM7_bECX;    DB 8;


    DW  FEMMS;
  end;
end;

//ベクトルの反転
function NowNeg(const a:D3DVECTOR):D3DVECTOR;
begin
  asm
    DW    FEMMS;

    DW    PXOR;       DB MM0_MM0; //0
    DW    PXOR;       DB MM1_MM1;

    DW    MOVD_ld;    DB MM2_bEAX;  DB 8;
    DW    _3DNow;     DB MM0_mEAX;  DB PFSUB;
    DW    _3DNow;     DB MM1_MM2;   DB PFSUB;

    DW    MOVQ_st;    DB MM0_mEDX;
    DW    MOVD_st;    DB MM1_bEDX;    DB 8;

    DW    FEMMS;
  end;
end;


//正規化
function NowNormalize(const a:D3DVECTOR):D3DVECTOR;
begin
  asm
    DW    FEMMS;

    DW    PXOR;       DB MM1_MM1;
    DW    MOVQ_ld;    DB MM0_mEAX;          //x y
    DW    MOVD_ld;    DB MM1_bEAX;  DB 8;   //z 0

    DW    MOVQ_ld;    DB MM2_MM0;
    DW    MOVQ_ld;    DB MM3_MM1;

    DW    _3DNow;     DB MM0_MM0;   DB PFMUL;
    DW    _3DNow;     DB MM1_MM1;   DB PFMUL;
    DW    _3DNow;     DB MM0_MM0;   DB PFACC;
    DW    _3DNow;     DB MM1_MM1;   DB PFACC;
    DW    _3DNow;     DB MM0_MM1;   DB PFADD;
    DW    _3DNow;     DB MM0_MM0;   DB PFRSQRT;   //長さの逆数

    DW    _3DNow;     DB MM2_MM0;   DB PFMUL;
    DW    _3DNow;     DB MM3_MM0;   DB PFMUL;     //元のベクトルを長さで割る

    DW    MOVQ_st;    DB MM2_mEDX;
    DW    MOVD_st;    DB MM3_bEDX;  DB 8;


    DW    FEMMS;
  end;
end;

//長さを得る
procedure NowLength(var dest:Single; const a:D3DVECTOR);
begin
  asm
    DW    FEMMS;

    DW    PXOR;       DB MM1_MM1;
    DW    MOVQ_ld;    DB MM0_mEDX;          //x y
    DW    MOVD_ld;    DB MM1_bEDX;  DB 8;   //z 0

    DW    _3DNow;     DB MM0_MM0;   DB PFMUL;
    DW    _3DNow;     DB MM1_MM1;   DB PFMUL;
    DW    _3DNow;     DB MM0_MM0;   DB PFACC;
    DW    _3DNow;     DB MM1_MM1;   DB PFACC;
    DW    _3DNow;     DB MM0_MM1;   DB PFADD;
    DW    _3DNow;     DB MM1_MM0;   DB PFRSQRT;
    DW    _3DNow;     DB MM0_MM1;   DB PFMUL;     //長さ

    DW    MOVD_st;    DB MM0_mEAX;

    DW    FEMMS;
  end;
end;

//スケーリング
function NowScale(const a:D3DVECTOR; const scale:Single):D3DVECTOR;
begin
  asm
    DW    FEMMS;

    mov   eax,a;
    mov   edx,scale;
    mov   ecx,@result;

    DW    MOVQ_ld;    DB MM0_mEAX;          // MM0 = x y
    DW    MOVQ_ld;    DB MM1_bEAX;  DB 8;   // MM1 = z w
    DW    MOVD_ld;    DB MM2_EDX;
    DW    PUNPCKLDQ;  DB MM2_MM2;           // MM2 = scale scale

    DW    _3DNow;     DB MM0_MM2;   DB PFMUL;
    DW    _3DNow;     DB MM1_MM2;   DB PFMUL;

    DW    MOVQ_st;    DB MM0_mECX;
    DW    MOVD_st;    DB MM1_bECX;  DB 8;

    DW    FEMMS;
  end;
end;

//長さを k に
function NowChangeLength(const a:D3DVECTOR; const k:Single):D3DVECTOR;
begin
  asm
    DW    FEMMS;

    mov   eax,a;

    DW    PXOR;       DB MM1_MM1;
    DW    MOVQ_ld;    DB MM0_mEAX;          //x y
    DW    MOVD_ld;    DB MM1_bEAX;  DB 8;   //z 0

    DW    MOVQ_ld;    DB MM2_MM0;
    DW    MOVQ_ld;    DB MM3_MM1;

    DW    _3DNow;     DB MM0_MM0;   DB PFMUL;
    DW    _3DNow;     DB MM1_MM1;   DB PFMUL;
    DW    _3DNow;     DB MM0_MM0;   DB PFACC;
    DW    _3DNow;     DB MM1_MM1;   DB PFACC;
    DW    _3DNow;     DB MM0_MM1;   DB PFADD;
    DW    _3DNow;     DB MM0_MM0;   DB PFRSQRT;   //長さの逆数

    mov   eax,k;
    DW    MOVD_ld;    DB MM5_EAX;
    DW    PUNPCKLDQ;  DB MM5_MM5;   //MM5 = k k
    DW    _3DNow;     DB MM0_MM5;   DB PFMUL;     //長さ * k

    DW    _3DNow;     DB MM2_MM0;   DB PFMUL;
    DW    _3DNow;     DB MM3_MM0;   DB PFMUL;     //元のベクトルを長さ*kで割る


    mov   edx,@result;

    DW    MOVQ_st;    DB MM2_mEDX;
    DW    MOVD_st;    DB MM3_bEDX;  DB 8;

    DW    FEMMS;
  end;
end;

//posからtargetに向いた長さkのベクトルを生成
function NowLookAt(const pos:D3DVECTOR; const target:D3DVECTOR; const k:Single):D3DVECTOR;
begin
  asm
    DW    FEMMS;

    mov   eax,pos;
    mov   edx,target;
    mov   ecx,k;

    DW    MOVQ_ld;    DB MM0_mEAX;          //x y (pos)
    DW    MOVD_ld;    DB MM1_bEAX;  DB 8;   //z 0

    mov   eax,@result;

    DW    MOVQ_ld;    DB MM2_mEDX;          //x y (target)
    DW    MOVD_ld;    DB MM3_bEDX;  DB 8;   //z 0

    DW    _3DNow;     DB MM0_MM2;  DB PFSUBR;
    DW    _3DNow;     DB MM1_MM3;  DB PFSUBR;


    DW    MOVQ_ld;    DB MM2_MM0;
    DW    MOVQ_ld;    DB MM3_MM1;

    DW    _3DNow;     DB MM0_MM0;   DB PFMUL;
    DW    _3DNow;     DB MM1_MM1;   DB PFMUL;
    DW    _3DNow;     DB MM0_MM0;   DB PFACC;
    DW    _3DNow;     DB MM1_MM1;   DB PFACC;
    DW    _3DNow;     DB MM0_MM1;   DB PFADD;
    DW    _3DNow;     DB MM0_MM0;   DB PFRSQRT;   //長さの逆数

    DW    MOVD_ld;    DB MM5_ECX;
    DW    PUNPCKLDQ;  DB MM5_MM5;   //MM5 = k k
    DW    _3DNow;     DB MM0_MM5;   DB PFMUL;     //長さ * k

    DW    _3DNow;     DB MM2_MM0;   DB PFMUL;
    DW    _3DNow;     DB MM3_MM0;   DB PFMUL;     //元のベクトルを長さ*kで割る



    DW    MOVQ_st;    DB MM2_mEAX;
    DW    MOVD_st;    DB MM3_bEAX;  DB 8;

    DW    FEMMS;
  end;
end;


//正規化除算(x,y,z,w) → (x/w, y/w, z/w, 1/w)
function NowViewFrustumToScreen(const a:THgVector):THgVector;
begin
  asm
    DW    FEMMS;
                                                  // L | H

    DW    MOVQ_ld;    DB MM0_mEAX;                // x | y
    DW    MOVQ_ld;    DB MM1_bEAX;  DB  8;        // z | w

    DW    MOVQ_ld;    DB MM2_MM1;
    DW    PUNPCKHDQ;  DB MM2_MM2;
    DW    _3DNow;     DB MM2_MM2;   DB PFRCP;     //1/w | 1/w

    DW    _3DNow;     DB MM0_MM2;   DB PFMUL;     // x/w | y/w
    DW    _3DNow;     DB MM1_MM2;   DB PFMUL;     // z/w | w/w

    DW    PUNPCKLDQ;  DB MM1_MM2;                 // z/w | 1/w

    DW    MOVQ_st;    DB MM0_mEDX;
    DW    MOVQ_st;    DB MM1_bEDX;  DB 8;

    DW    FEMMS;
  end;
end;


//////////////////////////////////////////////
//クォータニオン演算
//////////////////////////////////////////////


//足し算
function NowAddQuat(const a:TQuaternion; const b:TQuaternion):TQuaternion;
begin
  asm
    DW    FEMMS;

    DW    MOVQ_ld;    DB MM0_mEAX;
    DW    MOVQ_ld;    DB MM1_bEAX;  DB 8;
    DW    MOVQ_ld;    DB MM2_mEDX;
    DW    MOVQ_ld;    DB MM3_bEDX;  DB 8;

    DW    _3DNow;     DB MM0_MM2;   DB PFADD;
    DW    _3DNow;     DB MM1_MM3;   DB PFADD;

    DW    MOVQ_st;    DB MM0_mECX;
    DW    MOVQ_st;    DB MM1_bECX;  DB 8;

    DW    FEMMS;
  end;
end;


//引き算
function NowSubQuat(const a:TQuaternion; const b:TQuaternion):TQuaternion;
begin
  asm
    DW    FEMMS;

    DW    MOVQ_ld;    DB MM0_mEAX;
    DW    MOVQ_ld;    DB MM1_bEAX;  DB 8;
    DW    MOVQ_ld;    DB MM2_mEDX;
    DW    MOVQ_ld;    DB MM3_bEDX;  DB 8;

    DW    _3DNow;     DB MM0_MM2;   DB PFSUB;
    DW    _3DNow;     DB MM1_MM3;   DB PFSUB;

    DW    MOVQ_st;    DB MM0_mECX;
    DW    MOVQ_st;    DB MM1_bECX;  DB 8;

    DW    FEMMS;
  end;
end;


//乗算
function NowMulQuat(const a:TQuaternion; const b:TQuaternion):TQuaternion;
begin
  asm
    DW    FEMMS;

    //a1a1 a2a2　の格納

    DW    MOVD_ld;    DB MM5_bEAX;  DB 12;
    DW    PUNPCKLDQ;  DB MM5_MM5;               //a1 | a1
    DW    MOVQ_ld;    DB MM3_MM5;
    DW    _3DNow;     DB MM5_mEDX;  DB PFMUL;   //a1 * V2x | a1 *V2y
    DW    _3DNow;     DB MM3_bEDX;  DB 8;      DB PFMUL;    //a1 * V2z | a1 * V2w

    DW    MOVD_ld;    DB MM6_bEDX;  DB 12;
    DW    PUNPCKLDQ;  DB MM6_MM6;               //a2 | a2
    DW    MOVQ_ld;    DB MM4_MM6;
    DW    _3DNow;     DB MM6_mEAX;  DB PFMUL;   //a2 * V1x | a2 *V1y
    DW    _3DNow;     DB MM4_bEAX;  DB 8;      DB PFMUL;    //a2 * V1z | a2 * V1w

    DW    _3DNow;     DB MM6_MM5;   DB PFADD;
    DW    _3DNow;     DB MM4_MM3;   DB PFADD;
    DW    MOVQ_ld;    DB MM5_MM4;



    //外積の計算
                                                //        lo    hi
    DW    MOVQ_ld;    DB MM0_bEAX;  DB 4;       // MM0 = a.y | a.z
    DW    MOVD_ld;    DB MM2_mEDX;
    DW    MOVD_ld;    DB MM1_bEDX;  DB 8;
    DW    PUNPCKLDQ;  DB MM1_MM2;               // MM1 = b.z | b.x
    DW    _3DNow;     DB MM0_MM1;   DB  PFMUL;

    DW    MOVQ_ld;    DB MM2_bEDX;  DB 4        // MM2 = b.y | b.z
    DW    MOVD_ld;    DB MM4_mEAX;
    DW    MOVD_ld;    DB MM3_bEAX;  DB 8;
    DW    PUNPCKLDQ;  DB MM3_MM4;               // MM3 = a.z | b.x
    DW    _3DNow;     DB MM2_MM3;   DB PFMUL;

    DW    _3DNow;     DB MM0_MM2;   DB PFSUB;
    DW    _3DNow;     DB MM0_MM6;   DB PFADD;   //a1V2 + a2V1 を足す
    DW    MOVQ_st;    DB MM0_mECX;              //お先にx,yだけ格納


    DW    MOVD_ld;    DB MM0_mEAX;
    DW    MOVD_ld;    DB MM1_bEDX;  DB 4;
    DW    _3DNow;     DB MM0_MM1;   DB PFMUL;   // MM0 = a.x * b.y

    DW    MOVD_ld;    DB MM2_bEAX;  DB 4;
    DW    MOVD_ld;    DB MM3_mEDX;
    DW    _3DNow;     DB MM2_MM3;   DB PFMUL;   // MM2 = a.y * b.x

    DW    _3DNow;     DB MM0_MM2;   DB PFSUB;
    DW    _3DNow;     DB MM0_MM5;   DB PFADD;   //a1V2 + a2V1 を足す

    DW    MOVQ_ld;    DB MM7_MM0;               //MM7 = z | ?


    //内積の計算
    DW    MOVQ_ld;    DB MM0_mEAX;
    DW    MOVQ_ld;    DB MM1_mEDX;
    DW    _3DNow;     DB MM0_MM1;   DB PFMUL;
    DW    _3DNow;     DB MM0_MM0;   DB PFACC;

    DW    MOVQ_ld;    DB MM2_bEAX;  DB 8;
    DW    MOVD_ld;    DB MM3_bEDX;  DB 8;
    DW    _3DNow;     DB MM2_MM3;   DB PFMUL;
    DW    _3DNow;     DB MM2_MM2;   DB PFACC;

    DW    _3DNow;     DB MM0_MM2;   DB PFADD;   //MM0 = 内積

    DW    MOVD_ld;    DB MM1_bEAX;  DB 12;
    DW    MOVD_ld;    DB MM2_bEDX;  DB 12;
    DW    _3DNow;     DB MM1_MM2;   DB PFMUL;   //a1a2

    DW    _3DNow;     DB MM0_MM1;   DB PFSUBR;  //a1a2 - 内積


    //結果を格納
    DW    PUNPCKLDQ;  DB MM7_MM0;
    DW    MOVQ_st;    DB MM7_bECX;  DB 8;


    DW    FEMMS;
  end;
end;


//内積
procedure NowDotQuat(var dest:Single; const a:TQuaternion; const b:TQuaternion);
asm
  DW    FEMMS;

  DW    MOVQ_ld;    DB MM0_mEDX;  //a.x,y
  DW    MOVQ_ld;    DB MM1_mECX;  //b.x,y
  DW    MOVQ_ld;    DB MM2_bEDX;  DB 8;   //a.z,w
  DW    MOVQ_ld;    DB MM3_bECX;  DB 8;   //b.z,w

  DW    _3DNow;     DB MM0_MM1;   DB PFMUL;
  DW    _3DNow;     DB MM0_MM0;   DB PFACC;
  DW    _3DNow;     DB MM2_MM3;   DB PFMUL;
  DW    _3DNow;     DB MM2_MM2;   DB PFACC;

  DW    _3DNow;     DB MM0_MM2;   DB PFADD;

  DW    MOVD_st;    DB MM0_mEAX;  //dest

  DW    FEMMS;
end;

//ノルム
procedure NowLengthQuat(var dest:Single; const a:TQuaternion);
asm
  DW    FEMMS;

  DW    MOVQ_ld;    DB MM0_mEDX;          //a.x,y
  DW    MOVQ_ld;    DB MM2_bEDX;  DB 8;   //a.z,w

  DW    _3DNow;     DB MM0_MM0;   DB PFMUL;
  DW    _3DNow;     DB MM2_MM2;   DB PFMUL;
  DW    _3DNow;     DB MM0_MM0;   DB PFACC;
  DW    _3DNow;     DB MM2_MM2;   DB PFACC;

  DW    _3DNow;     DB MM0_MM2;   DB PFADD; //ここまでで、ノルムの２乗

  DW    MOVQ_ld;    DB MM1_MM0;               // 一旦退避
  DW    _3DNow;     DB MM0_MM1;   DB PFRSQRT; // 1/SQRT
  DW    _3DNow;     DB MM1_MM0;   DB PFMUL;   //掛けて出来上がり

  DW    MOVD_st;    DB MM1_mEAX;  //dest

  DW    FEMMS;
end;


//正規化
function NowNormalizeQuat(const a:TQuaternion):TQuaternion;
var
  len:Single;
begin
  NowLengthQuat(len, a);

  if len = 0.0 then begin
    result:=Quaternion(1, 0,0,0);
    exit;
  end;

  result.x:=a.X / len;
  result.y:=a.Y / len;
  result.z:=a.Z / len;
  result.w:=a.W / len;

end;



//回転を表すクォータニオンを作る
function  NowQuatFromRotation(const axis:D3DVECTOR; const t:Integer):TQuaternion;
begin
  asm
    DW    FEMMS;

    //q = ( cos T/2 , sinT/2 Vaxis )

    //まず、aを正規化
    DW    MOVQ_ld;    DB MM0_mEAX;          //x y
    DW    MOVD_ld;    DB MM1_bEAX;  DB 8;   //z 0

    DW    MOVQ_ld;    DB MM2_MM0;
    DW    MOVQ_ld;    DB MM3_MM1;

    DW    _3DNow;     DB MM0_MM0;   DB PFMUL;
    DW    _3DNow;     DB MM1_MM1;   DB PFMUL;
    DW    _3DNow;     DB MM0_MM0;   DB PFACC;
    DW    _3DNow;     DB MM1_MM1;   DB PFACC;
    DW    _3DNow;     DB MM0_MM1;   DB PFADD;
    DW    _3DNow;     DB MM0_MM0;   DB PFRSQRT;   //長さの逆数

    DW    _3DNow;     DB MM2_MM0;   DB PFMUL;
    DW    _3DNow;     DB MM3_MM0;   DB PFMUL;     //元のベクトルを長さで割る

    //MM2 = x y
    //MM3 = z 0


    //テーブルの参照とか
    sar   edx,1;
    mov   eax,PNowCosSinT;
    and   edx,4095;
    shl   edx,3;
    add   eax,edx;

    DW    MOVQ_ld;    DB MM7_mEAX;    //MM7(lo) = w;
    DW    MOVQ_ld;    DB MM6_MM7;
    DW    PUNPCKHDQ;  DB MM6_MM6;     //MM6 = sinT/2

    DW    _3DNow;     DB MM3_MM6;   DB PFMUL;
    DW    PUNPCKLDQ;  DB MM3_MM7;
    DW    _3DNow;     DB MM2_MM6;   DB PFMUL;


    //結果の格納
    DW    MOVQ_st;    DB MM2_mECX;
    DW    MOVQ_st;    DB MM3_bECX;  DB 8;


    DW    FEMMS;
  end;
end;


//正規化済みクォータニオンから行列を作る…ひでぇコード(^^;)
function  NowMatrixFromQuat(const a:TQuaternion):D3DMATRIX;
begin
  asm
    push  ebx;
    push  edi;

    DW    FEMMS;

    //読みこみ  MM0 = x   MM1 = y   MM2 = z   MM3 = w
    DW    MOVQ_ld;    DB MM0_mEAX;
    DW    MOVQ_ld;    DB MM1_MM0;
    DW    PUNPCKLDQ;  DB MM0_MM0;
    DW    PUNPCKHDQ;  DB MM1_MM1;

    DW    MOVQ_ld;    DB MM2_bEAX;   DB 8;
    DW    MOVQ_ld;    DB MM3_MM2;
    DW    PUNPCKLDQ;  DB MM2_MM2;
    DW    PUNPCKHDQ;  DB MM3_MM3;

    mov   eax,$3F800000;    //MM7 = 1.0
    DW    MOVD_ld;    DB MM7_EAX;


    // m11
    DW    MOVQ_ld;    DB MM4_MM1;
    DW    MOVQ_ld;    DB MM5_MM2;
    DW    PUNPCKLDQ;  DB MM4_MM5;   //MM4 = y | z
    DW    _3DNow;     DB MM4_MM4;   DB PFMUL;
    DW    _3DNow;     DB MM4_MM4;   DB PFACC;
    DW    _3DNow;     DB MM4_MM4;   DB PFACC;
    DW    _3DNow;     DB MM4_MM7;   DB PFSUBR;    //MM4 = 1-2(yy+zz)
    //m12
    DW    MOVQ_ld;    DB MM5_MM0;
    DW    MOVQ_ld;    DB MM6_MM1;
    DW    PUNPCKLDQ;  DB MM5_MM2;   //MM5 = x | z
    DW    PUNPCKLDQ;  DB MM6_MM3;   //MM6 = y | w
    DW    _3DNow;     DB MM5_MM6;   DB PFMUL;
    DW    _3Dnow;     DB MM5_MM5;   DB PFACC;
    DW    _3Dnow;     DB MM5_MM5;   DB PFACC;     //MM5 = 2(xy+wz)
    //m11,m12まとめる
    DW    PUNPCKLDQ;  DB MM4_MM5;
    DW    MOVQ_st;    DB MM4_mEDX;


    //m13 2(xz-wy)
    DW    MOVQ_ld;    DB MM4_MM0;
    DW    MOVQ_ld;    DB MM5_MM2;
    DW    PUNPCKLDQ;  DB MM4_MM3;   //MM4 = x | w
    DW    PUNPCKLDQ;  DB MM5_MM1;   //MM5 = z | y
    DW    _3DNow;     DB MM4_MM5;   DB PFMUL;
    DW    MOVQ_ld;    DB MM5_MM4;
    DW    PUNPCKHDQ;  DB MM5_MM5;
    DW    _3DNow;     DB MM4_MM5;   DB PFSUB;   //xz - wy
    DW    _3DNow;     DB MM4_MM4;   DB PFADD;   //2倍
    //m14 = 0
    DW    PXOR;       DB MM5_MM5;
    DW    PUNPCKLDQ;  DB MM4_MM5;
    DW    MOVQ_st;    DB MM4_bEDX;  DB 8;



    //m21 2(xy-wz)
    DW    MOVQ_ld;    DB MM4_MM0;
    DW    MOVQ_ld;    DB MM5_MM1;
    DW    PUNPCKLDQ;  DB MM4_MM3;   //MM4 = x | w
    DW    PUNPCKLDQ;  DB MM5_MM2;   //MM5 = y | z
    DW    _3DNow;     DB MM4_MM5;   DB PFMUL;
    DW    MOVQ_ld;    DB MM5_MM4;
    DW    PUNPCKHDQ;  DB MM5_MM5;
    DW    _3DNow;     DB MM4_MM5;   DB PFSUB;   //xy - wz
    DW    _3DNow;     DB MM4_MM4;   DB PFADD;   //2倍
    //m22 1-2(xx+zz)
    DW    MOVQ_ld;    DB MM5_MM0;
    DW    MOVQ_ld;    DB MM6_MM2;
    DW    PUNPCKLDQ;  DB MM5_MM6;   //MM4 = x | z
    DW    _3DNow;     DB MM5_MM5;   DB PFMUL;
    DW    _3DNow;     DB MM5_MM5;   DB PFACC;
    DW    _3DNow;     DB MM5_MM5;   DB PFACC;
    DW    _3DNow;     DB MM5_MM7;   DB PFSUBR;    //MM4 = 1-2(xx+zz)
    //m21,m22まとめる
    DW    PUNPCKLDQ;  DB MM4_MM5;
    DW    MOVQ_st;    DB MM4_bEDX;  DB 16

    //m23 2(yz+wx)
    DW    MOVQ_ld;    DB MM5_MM1;
    DW    MOVQ_ld;    DB MM6_MM2;
    DW    PUNPCKLDQ;  DB MM5_MM3;   //MM5 = y | w
    DW    PUNPCKLDQ;  DB MM6_MM0;   //MM6 = z | x
    DW    _3DNow;     DB MM5_MM6;   DB PFMUL;
    DW    _3Dnow;     DB MM5_MM5;   DB PFACC;
    DW    _3Dnow;     DB MM5_MM5;   DB PFACC;     //MM5 = 2(yz+wx)
    //m24 = 0
    DW    PXOR;       DB MM4_MM4;
    DW    PUNPCKLDQ;  DB MM5_MM4;
    DW    MOVQ_st;    DB MM5_bEDX;  DB 24;


    //m31 2(xz+wy)
    DW    MOVQ_ld;    DB MM4_MM0;
    DW    MOVQ_ld;    DB MM5_MM2;
    DW    PUNPCKLDQ;  DB MM4_MM3;   //MM4 = x | w
    DW    PUNPCKLDQ;  DB MM5_MM1;   //MM5 = z | y
    DW    _3DNow;     DB MM4_MM5;   DB PFMUL;
    DW    _3Dnow;     DB MM4_MM4;   DB PFACC;
    DW    _3Dnow;     DB MM4_MM4;   DB PFACC;     //MM5 = 2(xz+wy)
    //m32 2(yz-wx)
    DW    MOVQ_ld;    DB MM5_MM1;
    DW    MOVQ_ld;    DB MM6_MM2;
    DW    PUNPCKLDQ;  DB MM5_MM3;   //MM5 = y | w
    DW    PUNPCKLDQ;  DB MM6_MM0;   //MM6 = z | x
    DW    _3DNow;     DB MM5_MM6;   DB PFMUL;
    DW    MOVQ_ld;    DB MM6_MM5;
    DW    PUNPCKHDQ;  DB MM6_MM6;
    DW    _3DNow;     DB MM5_MM6;   DB PFSUB;   //xy - wz
    DW    _3DNow;     DB MM5_MM5;   DB PFADD;   //2倍
    //m31 & m32
    DW    PUNPCKLDQ;  DB MM4_MM5;
    DW    MOVQ_st;    DB MM4_bEDX;  DB 32;

    //m33 1-2(xx+yy)
    DW    MOVQ_ld;    DB MM4_MM0;
    DW    MOVQ_ld;    DB MM5_MM1;
    DW    PUNPCKLDQ;  DB MM4_MM5;   //MM4 = x | y
    DW    _3DNow;     DB MM4_MM4;   DB PFMUL;
    DW    _3DNow;     DB MM4_MM4;   DB PFACC;
    DW    _3DNow;     DB MM4_MM4;   DB PFACC;
    DW    _3DNow;     DB MM4_MM7;   DB PFSUBR;    //MM4 = 1-2(xx+yy)
    //m34 = 0とまとめる
    DW    PXOR;       DB MM5_MM5;
    DW    PUNPCKLDQ;  DB MM4_MM5;
    DW    MOVQ_st;    DB MM4_bEDX;  DB 40;

    //m41,m42 = 0
    DW    MOVQ_st;    DB MM5_bEDX;  DB 48;

    //m43=0 , m44=1;
    DW    PUNPCKLDQ;  DB MM5_MM7;
    DW    MOVQ_st;    DB MM5_bEDX;  DB 56;


    DW    FEMMS;

    pop   edi;
    pop   ebx;
  end;
end;

//axis回りの回転を表す行列を作る(QuadtFromRotationとMatrixFromQuatの合成)
function  NowMatrixFromRotation(const axis:D3DVECTOR; const t:Integer):D3DMATRIX;
begin
  result:=NowMatrixFromQuat(NowQuatFromRotation(axis,t));
end;


//Slerp rate=0.0ならaを100%、 rate=1.0ならbを100% の割合で球面線形補間
function  NowSlerpQuat(const a:TQuaternion; const b:TQuaternion; const rate:Single):TQuaternion;
var
  c,s:Single;     //a,bのなす角のCosとSin
  angle:Single; //a,bのなす角
  ra,rb:Single; //aとbの混合比

begin

  //Sin,Cosの使用が避けられないので、アセンブラ使っても大して高速化できないっしょ。

  NowDotQuat(c, a,b); //Cosθ

  //クォータニオンが正規化されてないっぽいので帰る
  if Abs(c) > 1.0 then begin
    result:=Quaternion(1, 0,0,0);
    exit;
  end;

  //ArcCosで、角度を求める
  angle:=ArcTan(sqrt (1-sqr (c)) /c);

  //角度が違わないなら帰る
  if angle = 0.0 then begin
    result:=a;
    exit;
  end;

  //Sinθを求めて
  s:=Sin(angle);

  //こんな感じで
  ra:=Sin( (1-rate) * angle ) / s;
  rb:=Sin( rate * angle ) / s;

  with result do begin
    x:=a.x*ra + b.x*rb;
    y:=a.y*ra + b.y*rb;
    z:=a.z*ra + b.z*rb;
    w:=a.w*ra + b.w*rb;
  end;

 result:=NowNormalizeQuat(result);
end;






//////////////////////3D演算(マトリクス)

//ベクトルを変換(同次…ジオメトリの変換に)
function NowTransform(const a:THgVector; const m:D3DMATRIX):THgVector;
begin
  asm
    DW    FEMMS;

    DW    movq_ld;    DB mm0_mEAX;		          // ay, ax	in mm0
    DW    movq_ld;    DB mm1_bEAX;  DB 8;		    // aw, az	in mm1
    //
    DW    movq_ld;    DB mm4_mEDX;      		    // m01, m00	in mm4
    DW    punpckhdq;	DB mm2_mm0;		            // ay, undef	in mm2
    //
    DW    movq_ld;       DB mm5_bEDX;  DB 16;	    // m11, m10	in mm5
    DW    punpckldq;  DB mm0_mm0;       		    // ax, ax	in mm2
    //
    DW    movq_ld;    DB mm6_bEDX;  DB 32;      // m21, m20	in mm6
    DW    _3DNow;     DB mm4_mm0;   DB pfmul;   // ax*m01, ax*m00	in mm4
    //
    DW    movq_ld;    DB mm7_bEDX;  DB 48;      // m31, m30	in mm7
    DW    punpckhdq;  DB mm2_mm2;           		// ay, ay	in mm2
    //
    DW    punpckhdq;  DB mm3_mm1;           		// aw, undef	in mm3
    DW    _3DNow;     DB mm5_mm2;   DB pfmul; 	// ay*m11, ay*m10
    //
    DW    punpckldq;  DB mm1_mm1;           		// az, az	in mm3
    DW    _3DNow;     DB mm0_bEDX;  DB 8;   DB pfmul;   	// ax*m03, ax*m02 in mm0
    //
    DW    punpckhdq;  DB mm3_mm3;            		// aw, aw	in mm3
    DW    _3DNow;     DB mm2_bEDX;  DB 24;  DB PFMUL;     // ay*13, ay*m12	in mm2
    //
    DW    _3DNow;     DB mm6_mm1;   DB pfmul; 	// az*m21, az*m20	in mm6
    DW    _3DNow;     DB mm5_mm4;   DB pfadd;
    //
    DW    _3DNow;     DB mm1_bEDX;  DB 40;    DB pfmul    // az*23, az*22	in mm1
    DW    _3DNow;     DB mm2_mm0;   DB pfadd;
    //
    DW    _3DNow;     DB mm7_mm3;   DB pfmul;   // aw*m31, aw*m30	in mm7
    DW    _3DNow;     DB mm6_mm5;   DB pfadd;
    //
    DW    _3DNow;     DB mm3_bEDX;  DB 56;    DB pfmul;   // aw*m33, aw*m32	in mm3
    DW    _3DNow;     DB mm2_mm1;   DB pfadd;
    //
    DW    _3DNow;     DB mm7_mm6;   DB pfadd    // ry, rx	in mm7
    DW    _3DNow;     DB mm3_mm2;   DB pfadd		// rw,rz		in mm3
    DW    movq_st;    DB mm7_mECX;
    DW    movq_st;    DB mm3_bECX;  DB 8;

    DW    FEMMS;

  end;
end;

//ベクトルを一気に変換
procedure NowTransformBatch(const src:Pointer; const m:D3DMATRIX; const num:Integer; const span:Integer);
begin
  asm
    push  ebx;

    DW    FEMMS;

    DW    PREFETCH;   DB MM0_bEDX;  DB 0;
    DW    PREFETCH;   DB MM0_bEDX;  DB 8;

    //まじめ計算

    DW    MOVQ_ld;    DB MM0_bEAX;  DB 0;     //MM0 = x
    DW    MOVQ_ld;    DB MM2_bEAX;  DB 8;     //MM2 = z
    DW    MOVQ_ld;    DB MM1_MM0;             //MM1 = y
    DW    MOVQ_ld;    DB MM3_MM2              //MM3 = w


    @TransformLoop:

    DW    PREFETCHW;  DB MM1_mEAX;            //この頂点と次の頂点をキャッシュに入れる

    DW    PUNPCKLDQ;  DB MM0_MM0;
    DW    MOVQ_ld;    DB MM4_mEDX;            //mm4 = m11 | m12
    DW    PUNPCKHDQ;  DB MM1_MM1;
    DW    MOVQ_ld;    DB MM5_bEDX;  DB 16;    //mm5 = m21 | m22
    DW    PUNPCKLDQ;  DB MM2_MM2;
    DW    MOVQ_ld;    DB MM6_bEDX;  DB 32     //mm4 = m31 | m32
    DW    PUNPCKHDQ;  DB MM3_MM3;
    DW    MOVQ_ld;    DB MM7_bEDX;  DB 48;    //mm5 = m41 | m42

    DW    _3DNow;     DB MM4_MM0;   DB PFMUL;
    DW    _3DNow;     DB MM5_MM1;   DB PFMUL;
    DW    _3DNow;     DB MM6_MM2;   DB PFMUL;
    DW    _3DNow;     DB MM7_MM3;   DB PFMUL;

    DW    _3DNow;     DB MM4_MM5;   DB PFADD;
    DW    MOVQ_ld;    DB MM5_bEDX;  DB 24 ;   //mm5 = m23 | m24
    DW    _3DNow;     DB MM6_MM7;   DB PFADD;
    DW    MOVQ_ld;    DB MM7_bEDX;  DB 56;    //mm5 = m43 | m44
    DW    _3DNow;     DB MM4_MM6;   DB PFADD;
    DW    MOVQ_ld;    DB MM6_bEDX;  DB 40     //mm4 = m33 | m34

    DW    MOVQ_st;    DB MM4_mEAX;            //x',y'を格納


    add   eax,span;
    DW    MOVQ_ld;    DB MM4_bEDX;  DB 8      //mm4 = m13 | m14


    DW    _3DNow;     DB MM4_MM0;   DB PFMUL;
    DW    MOVQ_ld;    DB MM0_mEAX;            //MM0 = x

    DW    _3DNow;     DB MM5_MM1;   DB PFMUL;
    DW    MOVQ_ld;    DB MM1_MM0;             //MM1 = y

    DW    _3DNow;     DB MM6_MM2;   DB PFMUL;
    DW    MOVQ_ld;    DB MM2_bEAX;  DB 8;     //MM2 = z

    DW    _3DNow;     DB MM7_MM3;   DB PFMUL;
    DW    MOVQ_ld;    DB MM3_MM2              //MM3 = w


    DW    _3DNow;     DB MM4_MM5;   DB PFADD;
    DW    _3DNow;     DB MM6_MM7;   DB PFADD;

    DW    _3DNow;     DB MM4_MM6;   DB PFADD;
    DW    MOVQ_st;    DB MM4_bEAX;  DB -8;     //z',w'を格納

    LOOP  @TransformLoop;


    {
    //不真面目ルーチン(w=1で決め打ち)

    DW    MOVD_ld;    DB MM6_bEDX;  DB 8
    DW    PUNPCKLDQ;  DB MM6_bEDX;  DB 24     //MM6 = m13 | m23
    DW    MOVD_ld;    DB MM7_bEDX;  DB 40;
    DW    PUNPCKLDQ;  DB MM7_bEDX;  DB 56;    //MM7 = m33 | m43
    mov   ebx,$3F800000;


    @TransformLoop:

    DW    PREFETCH;   DB MM0_mEAX;

    DW    MOVQ_ld;    DB MM0_mEAX;            //MM0 = x | x
    DW    MOVQ_ld;    DB MM2_bEAX;  DB 8;     //MM2 = z | z
    DW    MOVQ_ld;    DB MM1_MM0;             //MM1 = y | y
    DW    PUNPCKLDQ;  DB MM0_MM0;
    DW    PUNPCKHDQ;  DB MM1_MM1;
    DW    PUNPCKLDQ;  DB MM2_MM2;

    DW    MOVQ_ld;    DB MM3_MM0;
    DW    MOVQ_ld;    DB MM4_MM1;
    DW    MOVQ_ld;    DB MM5_MM2;
    DW    _3DNow;     DB MM3_mEDX;  DB PFMUL;
    DW    _3DNow;     DB MM4_bEDX;  DB 16;    DB PFMUL;
    DW    _3DNow;     DB MM5_bEDX;  DB 32;    DB PFMUL;
    DW    _3DNow;     DB MM3_MM4;   DB PFADD;
    DW    _3DNow;     DB MM3_MM5;   DB PFADD;
    DW    _3DNow;     DB MM3_bEDX;  DB 48;    DB PFADD;
    DW    MOVQ_st;    DB MM3_mEAX;            // x'とy'を格納

    DW    PUNPCKLDQ;  DB MM0_MM1;             //MM0 = x | y
    DW    _3DNow;     DB MM0_MM6;   DB PFMUL;
    DW    MOVD_ld;    DB MM4_EBX;
    DW    PUNPCKLDQ;  DB MM2_MM4;             //MM2 = z | 1

    DW    _3DNow;     DB MM2_MM7;   DB PFMUL;
    DW    _3DNow;     DB MM0_MM2;   DB PFADD;
    DW    _3DNow;     DB MM0_MM0;   DB PFACC;

    DW    PUNPCKLDQ;  DB MM0_MM4;

    DW    MOVQ_st;    DB MM0_bEAX;  DB 8;

    add   eax,16;

    LOOP  @TransformLoop;

    }

    DW    FEMMS;
    pop   ebx;
  end;
end;

//カメラとターゲットの位置から、ビュー行列を作る
//ワールド座標→視野座標
function NowViewMatrix(const from:D3DVECTOR; const at:D3DVECTOR; const world_up:D3DVECTOR; const roll:Integer):D3DMATRIX;
var
  up,right,direction:D3DVECTOR; //頭の中心→頭のてっぺん、頭の中心→右耳、視線の方向
begin
  direction:=NowNormalize(NowSub(at,from));
  right := NowCrossProduct(world_up, direction);
  up:=NowCrossProduct(direction,right);

  right:=NowNormalize(right);
  up:=NowNormalize(up);


  result:=NowIdentityMatrix;
  result._11 := right.x;
  result._21 := right.y;
  result._31 := right.z;
  result._12 := up.x;
  result._22 := up.y;
  result._32 := up.z;
  result._13 := direction.x;
  result._23 := direction.y;
  result._33 := direction.z;

  NowDotProduct(result._41, right, from);
  result._41:=-result._41;
  NowDotProduct(result._42, up, from);
  result._42:=-result._42;
  NowDotProduct(result._43, direction, from);
  result._43:=-result._43;

  result._44:=1;

  if roll <> 0 then begin
    result:=NowCompositeMatrix(NowRotZMatrix(roll),result);
  end;

end;


//姿勢・位置からビュー行列を作る
function NowViewMatixFromOrientation(const from:D3DVECTOR; const dir:D3DVECTOR; const up:D3DVECTOR):D3DMATRIX;
var
  vx,vy,vz:D3DVECTOR;
begin
  vy:=NowNormalize(up);
  vz:=NowNormalize(dir);
  vx:=NowNormalize(NowCrossProduct(up,dir));

  result:=NowIdentityMatrix;
  with result do begin
    _11:=vx.x;    _12:=vy.x;    _13:=vz.x;
    _21:=vx.y;    _22:=vy.y;    _23:=vz.y;
    _31:=vx.z;    _32:=vy.z;    _33:=vz.z;

    NowDotProduct(_41, vx,from);
    _41:=-_41;
    NowDotProduct(_42, vy,from);
    _42:=-_42;
    NowDotProduct(_43, vz,from);
    _43:=-_43;
  end;
end;


//透視射影行列(視野座標→除算前正規化座標)を作る
function NowProjectFrustum(const halfWidth:Single; const halfHeight:Single; const nearZ:Single; const farZ:Single):D3DMATRIX;
begin
  ZeroMemory(@result,sizeof(D3DMATRIX));

  with result do begin
    _11:=nearZ * halfHeight * 2;

    _22:=-_11;

    _31:=halfWidth;
    _32:=halfHeight;
    _33:=farZ / (farZ-nearZ);
    _34:=1.0;

    _43:=_33 * (-nearZ);
  end;
end;

function NowRotXMatrix(const t:Integer):D3DMATRIX;
var
  t2:Integer;
begin
  t2:=t And 4095;

  Move(NowIdentityMatrix,result,sizeof(D3DMATRIX));
  with result do begin
    _22:=NowCosT[t2];   _23:=NowSinT[t2];
    _32:=-NowSinT[t2];  _33:=NowCosT[t2];
  end;
end;

function NowRotYMatrix(const t:Integer):D3DMATRIX;
var
  t2:Integer;
begin
  t2:=t And 4095;

  Move(NowIdentityMatrix,result,sizeof(D3DMATRIX));
  with result do begin
    _11:=NowCosT[t2]; _13:=-NowSinT[t2];
    _31:=NowSinT[t2]; _33:=NowCosT[t2];   _34:=0;
  end;
end;

function NowRotZMatrix(const t:Integer):D3DMATRIX;
var
  t2:Integer;
begin
  t2:=t And 4095;

  Move(NowIdentityMatrix,result,sizeof(D3DMATRIX));
  with result do begin
    _11:=NowCosT[t2];   _12:=NowSinT[t2];
    _21:=-NowSinT[t2];  _22:=NowCosT[t2];
  end;
end;


//スケーリング行列
function NowScaleMatrix(const x,y,z:Single):D3DMATRIX;
begin
  with result do begin
    _11:=x;   _12:=0; _13:=0; _14:=0;
    _21:=0;   _22:=y; _23:=0; _24:=0;
    _31:=0;   _32:=0; _33:=z; _34:=0;
    _41:=0;   _42:=0; _43:=0; _44:=1;
  end;
end;

//平行移動行列
function NowTranslateMatrix(const x,y,z:Single):D3DMATRIX;
begin
  with result do begin
    _11:=0;   _12:=0; _13:=0; _14:=0;
    _21:=0;   _22:=0; _23:=0; _24:=0;
    _31:=0;   _32:=0; _33:=0; _34:=0;
    _41:=x;   _42:=y; _43:=z; _44:=1;
  end;
end;


//逆行列
function NowInvMatrix(const m:D3DMATRIX):D3DMATRIX;
begin
  with result do begin
    _11:=m._11;   _12:=m._21;   _13:=m._31;   _14:=0;
    _21:=m._12;   _22:=m._22;   _23:=m._32;   _24:=0;
    _31:=m._13;   _32:=m._23;   _33:=m._33;   _34:=0;

    _41:= - (m._11*m._41 + m._12*m._42 + m._13*m._43);
    _42:= - (m._21*m._41 + m._22*m._42 + m._23*m._43);
    _43:= - (m._31*m._41 + m._32*m._42 + m._33*m._43);
    _44:= 1;
  end;
end;

//転置行列
function NowTMatrix(const m:D3DMATRIX):D3DMATRIX;
begin
  with result do begin
    _11:=m._11;   _12:=m._21;   _13:=m._31;   _14:=m._41;
    _21:=m._12;   _22:=m._22;   _23:=m._32;   _24:=m._42;
    _31:=m._13;   _32:=m._23;   _33:=m._33;   _34:=m._43;
    _41:=m._14;   _42:=m._24;   _43:=m._34;   _44:=m._44;
  end;
end;


//回転部分だけ抜き出し
function NowExtractRotation(const m:D3DMATRIX):D3DMATRIX;
begin
  result:=NowIdentityMatrix;
  with result do begin
    _11:=m._11;   _12:=m._12;   _13:=m._13;
    _21:=m._21;   _22:=m._22;   _23:=m._23;
    _31:=m._31;   _32:=m._32;   _33:=m._33;
  end;
end;

//平行移動部分だけ抜き出し
function NowExtractTranslation(const m:D3DMATRIX):D3DMATRIX;
begin
  result:=NowIdentityMatrix;
  with result do begin
    _41:=m._41;   _42:=m._42;   _43:=m._43;
  end;
end;

//マトリクスの合成
function NowCompositeMatrix(const m:D3DMATRIX; const n:D3DMATRIX):D3DMATRIX;
begin
  asm
    push  ebx;

    DW    FEMMS;

    mov   ebx,4;


    @LOOP_1:;

    DW    PREFETCH;   DB MM0_bEAX;  DB 16;

    DW    movq_ld;    DB mm0_mEAX;		          // ay, ax	in mm0
    DW    movq_ld;    DB mm1_bEAX;  DB 8;		    // aw, az	in mm1
    //
    DW    movq_ld;    DB mm4_mEDX;      		    // m01, m00	in mm4
    DW    punpckhdq;	DB mm2_mm0;		            // ay, undef	in mm2
    //
    DW    movq_ld;       DB mm5_bEDX;  DB 16;	    // m11, m10	in mm5
    DW    punpckldq;  DB mm0_mm0;       		    // ax, ax	in mm2
    //
    DW    movq_ld;    DB mm6_bEDX;  DB 32;      // m21, m20	in mm6
    DW    _3DNow;     DB mm4_mm0;   DB pfmul;   // ax*m01, ax*m00	in mm4
    //
    DW    movq_ld;    DB mm7_bEDX;  DB 48;      // m31, m30	in mm7
    DW    punpckhdq;  DB mm2_mm2;           		// ay, ay	in mm2
    //
    DW    punpckhdq;  DB mm3_mm1;           		// aw, undef	in mm3
    DW    _3DNow;     DB mm5_mm2;   DB pfmul; 	// ay*m11, ay*m10
    //
    DW    punpckldq;  DB mm1_mm1;           		// az, az	in mm3
    DW    _3DNow;     DB mm0_bEDX;  DB 8;   DB pfmul;   	// ax*m03, ax*m02 in mm0
    //
    DW    punpckhdq;  DB mm3_mm3;            		// aw, aw	in mm3
    DW    _3DNow;     DB mm2_bEDX;  DB 24;  DB PFMUL;     // ay*13, ay*m12	in mm2
    //
    DW    _3DNow;     DB mm6_mm1;   DB pfmul; 	// az*m21, az*m20	in mm6
    DW    _3DNow;     DB mm5_mm4;   DB pfadd;
    //
    DW    _3DNow;     DB mm1_bEDX;  DB 40;    DB pfmul    // az*23, az*22	in mm1
    DW    _3DNow;     DB mm2_mm0;   DB pfadd;
    //
    DW    _3DNow;     DB mm7_mm3;   DB pfmul;   // aw*m31, aw*m30	in mm7
    DW    _3DNow;     DB mm6_mm5;   DB pfadd;
    //
    DW    _3DNow;     DB mm3_bEDX;  DB 56;    DB pfmul;   // aw*m33, aw*m32	in mm3
    DW    _3DNow;     DB mm2_mm1;   DB pfadd;
    //
    DW    _3DNow;     DB mm7_mm6;   DB pfadd    // ry, rx	in mm7
    DW    _3DNow;     DB mm3_mm2;   DB pfadd		// rw,rz		in mm3
    DW    movq_st;    DB mm7_mECX;
    DW    movq_st;    DB mm3_bECX;  DB 8;


    {
    DW    MOVQ_ld;    DB MM0_mEAX;            //MM0 = x
    DW    MOVQ_ld;    DB MM2_bEAX;    DB 8;   //MM2 = z
    DW    MOVQ_ld;    DB MM1_MM0;             //MM1 = y
    DW    MOVQ_ld;    DB MM3_MM2              //MM3 = w
    DW    PUNPCKLDQ;  DB MM0_MM0;
    DW    PUNPCKHDQ;  DB MM1_MM1;
    DW    PUNPCKLDQ;  DB MM2_MM2;
    DW    PUNPCKHDQ;  DB MM3_MM3;

    DW    MOVQ_ld;    DB MM4_mEDX;            //mm4 = m11 | m12
    DW    MOVQ_ld;    DB MM5_bEDX;  DB 16;    //mm5 = m21 | m22
    DW    _3DNow;     DB MM4_MM0;   DB PFMUL;
    DW    _3DNow;     DB MM5_MM1;   DB PFMUL;

    DW    MOVQ_ld;    DB MM6_bEDX;  DB 32     //mm4 = m31 | m32
    DW    MOVQ_ld;    DB MM7_bEDX;  DB 48;    //mm5 = m41 | m42
    DW    _3DNow;     DB MM6_MM2;   DB PFMUL;
    DW    _3DNow;     DB MM7_MM3;   DB PFMUL;

    DW    _3DNow;     DB MM4_MM5;   DB PFADD;
    DW    _3DNow;     DB MM6_MM7;   DB PFADD;
    DW    _3DNow;     DB MM4_MM6;   DB PFADD;

    DW    MOVQ_st;    DB MM4_mECX;            //x',y'を格納


    DW    MOVQ_ld;    DB MM4_bEDX;  DB 8      //mm4 = m13 | m14
    DW    MOVQ_ld;    DB MM5_bEDX;  DB 24 ;   //mm5 = m23 | m24
    DW    _3DNow;     DB MM4_MM0;   DB PFMUL;
    DW    _3DNow;     DB MM5_MM1;   DB PFMUL;

    DW    MOVQ_ld;    DB MM6_bEDX;  DB 40     //mm4 = m33 | m34
    DW    MOVQ_ld;    DB MM7_bEDX;  DB 56;    //mm5 = m43 | m44
    DW    _3DNow;     DB MM6_MM2;   DB PFMUL;
    DW    _3DNow;     DB MM7_MM3;   DB PFMUL;

    DW    _3DNow;     DB MM4_MM5;   DB PFADD;
    DW    _3DNow;     DB MM6_MM7;   DB PFADD;
    DW    _3DNow;     DB MM4_MM6;   DB PFADD;

    DW    MOVQ_st;    DB MM4_bECX;  DB 8;     //z',w'を格納
    }

    add   eax,16;
    add   ecx,16;

    dec   ebx;
    jnz   @LOOP_1;


    DW    FEMMS;
    pop   ebx;
  end;
end;


{$ELSE}

///////////////////////////////////////////////////////////////////////////////
// x87バージョン
///////////////////////////////////////////////////////////////////////////////


//////////////////////変換

function  NowItoF2D(const a:TPoint):TVector2D;
begin
  result.X:=a.X;
  result.Y:=a.Y;
end;

function  NowFtoI2D(const a:TVector2D):TPoint;
begin
  result.X:=Trunc(a.X);
  result.Y:=Trunc(a.Y);
end;

function  NowRound2D(const a:TVector2D):TPoint;
begin
  result.X:=Round(a.X);
  result.Y:=Round(a.Y);
end;

function  NowLVToTLV(const lv:D3DLVERTEX; const transMat:D3DMATRIX):D3DTLVERTEX;
var
  hg:THgVector;  //頂点座標
begin
  Move(lv,hg,Sizeof(D3DVECTOR));
  hg.W:=1;

  hg:=NowViewFrustumToScreen(NowTransform(hg,transMat));
  Move(hg,result,sizeof(THgVector));
  Move(lv.color,result.color,sizeof(D3DTLVERTEX) - sizeof(THgVector));

end;


procedure NowLVToTLVBatch(src:LPD3DLVERTEX; dest:LPD3DTLVERTEX; const transMat:D3DMATRIX; count:Integer);
var
  i:Integer;
begin
  for i:=0 to count-1 do begin
    dest^:=NowLVToTLV(src^,transMat);
    Inc(dest);
    Inc(src);
  end;
end;


//TMatrix33からD3DNATRIXへ、
function NowMatrix33To44(m:TMatrix33):D3DMATRIX;
begin
  result:=NowIdentityMatrix;
  with result do begin
    _11:=m._11; _12:=m._12; _13:=m._13;
    _21:=m._21; _22:=m._22; _23:=m._23;
    _31:=m._31; _32:=m._32; _33:=m._33;
  end;
end;

//44行列から左上3x3要素だけ抜いて33行列にする
function NowMatrix44To33(m:D3DMATRIX):TMatrix33;
begin
  with result do begin
    _11:=m._11; _12:=m._12; _13:=m._13;
    _21:=m._21; _22:=m._22; _23:=m._23;
    _31:=m._31; _32:=m._32; _33:=m._33;
  end;
end;

//３頂点から平面
function NowPlaneFromPoints(v1,v2,v3:D3DVector):TPlane3D;
begin
  result:=NowPlaneFromVertexNormal(v1, NowNormalFromPoints(v1,v2,v3));
end;

//１頂点vと法線nから平面
function NowPlaneFromVertexNormal(v,n:D3DVector):TPlane3D;
begin
  result.A:=n.x;
  result.B:=n.y;
  result.C:=n.z;
  result.D:=-NowDotProduct(v,n);
end;


//////////////////////1D演算

//スプライン補間(4点を与えられて、中間2点のブレンド比を与えて補間
function NowInterpolate(const p1,p2,p3,p4:Single; const rate:Single):Single;
var
  ir:Single;      //1-rate
  r2,r3:Single;   //rate^2　, rate^3
begin
  ir:=1-rate;
  r2:=rate*rate;
  r3:=r2*rate;

  result:= ( (p1*ir*ir*ir) + p2*(3*r3-6*r2+4) + p3*(-3*(r3-r2-rate)+1) + p4*r3)  / 6.0;
end;

//線形補間(rateは２点のブレンド比で0〜1、1に近いほどp2が濃くなる。
function NowLerp(const p1,p2:Single; const rate:Single):Single;
begin
  result:=rate*p2 + (1.0-rate)*p1;
end;

function SQR(x:Single):Single;
begin
  result:=x*x;
end;

function NowArcTan2(Y, X: Extended): Extended;
var
  rad:Extended;
begin
  if x <> 0 then begin
    rad := ArcTan(y/x);

    if x < 0 then
      rad:=rad + Pi
    else if y < 0 then
      rad:=rad + 2*Pi;

  end else begin
    if y > 0 then
      rad:=Pi*0.5
    else
      rad:=Pi*1.5;
  end;


  result:=rad;
end;


//////////////////////2D演算

//ベクトルの足し算
function NowAdd2D(const a:TVector2D; const b:TVector2D):TVector2D;
begin
  result.X:=a.X + b.X;
  result.Y:=a.Y + b.Y;
end;

//ベクトルの引き算
function NowSub2D(const a:TVector2D; const b:TVector2D):TVector2D;
begin
  result.X:=a.X - b.X;
  result.Y:=a.Y - b.Y;
end;

//内積
function NowDotProduct2D(const a:TVector2D; const b:TVector2D):Single;
begin
  result:=a.X * b.X + a.Y * b.Y;
end;

//2点の距離
function NowDistance2D(const a,b:TVector2D):Single;
begin
  result:=sqrt(SQR(a.X-b.X) + SQR(a.Y-b.Y));
end;

//長さを算出
function NowLength2D(const a:TVector2D):Single;
begin
  result:=sqrt(SQR(a.X) + SQR(a.Y));
end;

//ベクトルの向きを4096度系で
function NowVecToAngle2D(const a:TVector2D):Single;
var
  v:TVector2D;
  sg:Integer;   //象限
  hi,lo,mid:Integer;  //2分探索用
begin
  //まず、ベクトルをノーマライズする
  v:=NowNormalize2D(a);

  //象限はどれ？
  if (v.y > 0) then begin
    if v.x > 0 then
      sg:=1
    else begin
      sg:=2;
      v.x:=-v.x;
    end;
  end else if (v.y < 0) then begin
    if v.x >= 0 then
      sg:=4
    else begin
      sg:=3;
      v.x:=-v.x;
    end;
  end else begin
    if v.x >= 0 then
      sg:=1
    else begin
      sg:=3;
      v.x:=-v.x;
    end;
  end;


  hi:=1023;
  lo:=0;
  mid:=511;
  //Cosテーブル内を2分探索し、角度を求める
  while hi > lo do begin
    if NowCosT[mid] > v.x then
      lo:=mid+1
    else
      hi:=mid;
    mid:=(hi+lo) Shr 1;
  end;

  //象限に応じて、結果を加工
  Case sg of
    1: result:=mid;
    2: result:=2047-mid;
    3: result:=2048+mid;
    4: result:=4095-mid;
  end;

end;


//回転(原点を右回りにangleだけ)
function NowRot2D(const a:TVector2D; const angle:Integer):TVector2D;
var
  t:Integer;
begin
  t:=angle and 4095;
  result.X:=a.X * NowCosT[t] - a.Y * NowSinT[t];
  result.Y:=a.X * NowSinT[t] + a.Y * NowCosT[t];
end;

//正規化(長さを1に)
function NowNormalize2D(const a:TVector2D):TVector2D;
var
  len:Single;
begin
  len:=NowLength2D(a);
  if len = 0 then begin
    result:=Vector2D(0,0);
    Exit;
  end;

  result.X:=a.X / len;
  result.Y:=a.Y / len;
end;

//スケーリング
function NowScale2D(const a:TVector2D; const scale:Single):TVector2D;
begin
  result.X:=a.X * scale;
  result.Y:=a.Y * scale;
end;

//ベクトルの長さをkにする
function NowChangeLength2D(const a:TVector2D; const k:Single):TVector2D;
var
  len:Single;
begin
  len:=NowLength2D(a);
  if len = 0 then begin
    result:=Vector2D(0,0);
    Exit;
  end;

  result.X:=a.X * k / len;
  result.Y:=a.Y * k / len;
end;

//posからtargetに向いた長さkのベクトルを生成
function NowLookAt2D(const pos:TVector2D; const target:TVector2D; const k:Single):TVector2D;
var
  sub:TVector2D;
  len:Single;
begin
  sub:=NowSub2D(target,pos);
  len:=NowLength2D(sub);
  if len = 0 then begin
    result:=Vector2D(0,0);
    Exit;
  end;

  result.X:=sub.X * k / len;
  result.Y:=sub.Y * k / len;
end;

//長さkのランダムなベクトルを作る
function NowRandom2D(const k:Single):TVector2D;
begin
  result:=NowChangeLength2D(Vector2D(Random-0.5,Random-0.5),k);
end;

//線形補間
function NowLerp2D(const a:TVector2D; const b:TVector2D; const rate:Single):TVector2D;
begin
  result.x:=rate*b.x + (1.0-rate)*a.x;
  result.y:=rate*b.y + (1.0-rate)*a.y;
end;


//////////////////////3D演算(ベクトル)

//足し算
function NowAdd(const a:D3DVECTOR; const b:D3DVECTOR):D3DVECTOR;
begin
  result.X:=a.X + b.X;
  result.Y:=a.Y + b.Y;
  result.Z:=a.Z + b.Z;
end;

function NowAddHg(const a:THgVector; const b:THgVector):THgVector;
begin
  result.X:=a.X + b.X;
  result.Y:=a.Y + b.Y;
  result.Z:=a.Z + b.Z;
  result.W:=a.W + b.W;
end;


//引き算
function NowSub(const a:D3DVECTOR; const b:D3DVECTOR):D3DVECTOR;
begin
  result.X:=a.X - b.X;
  result.Y:=a.Y - b.Y;
  result.Z:=a.Z - b.Z;
end;

function NowSubHg(const a:THgVector; const b:THgVector):THgVector;
begin
  result.X:=a.X - b.X;
  result.Y:=a.Y - b.Y;
  result.Z:=a.Z - b.Z;
  result.W:=a.W - b.W;
end;


//内積
function NowDotProduct(const a:D3DVECTOR; const b:D3DVECTOR):Single;
begin
  result:=a.X*b.X + a.Y*b.Y + a.Z*b.Z;
end;

//外積
function NowCrossProduct(const a:D3DVECTOR; const b:D3DVECTOR):D3DVECTOR;
begin
  result.x := a.y*b.z - a.z*b.y;
  result.y := a.z*b.x - a.x*b.z;
  result.z := a.x*b.y - a.y*b.x;
end;


//2点の距離
function NowDistance(const a,b:D3DVECTOR):Single;
begin
  result:=sqrt(SQR(a.X-b.X) + SQR(a.Y-b.Y) + SQR(a.Z-b.Z));
end;


//点と直線の距離
function NowDistance(const p:D3DVECTOR; const lin:TLine3D):Single;
var
  h,dif:D3Dvector;
  t:Single;
begin
  dif := NowSub(p, lin.Origin);
  t := NowDotProduct(lin.Direction, dif);
  h := NowSub(dif, NowScale(lin.Direction, t));
  result:=NowLength(h);
end;

//点と線分の距離
function NowDistance(const p:D3DVECTOR; const lseg:TLineSeg3D):Single;
var
  dif:D3DVector;
  t:Single;
  sqrLen:Single;
  dseg:D3DVector;
begin
  dif := NowSub(p, lseg.A);
  dseg:= NowSub(lseg.B, lseg.A);
  sqrLen:=NowDotProduct(dseg,dseg);
  t := NowDotProduct(dif, dseg);

  if t >= 0.0 then begin
    if t >= sqrLen then begin
      dif:=NowSub(dif, dseg);
    end else begin
      t:=t / sqrLen;
      dif := NowSub(dif, NowScale(dseg, t));
    end;
  end;

  result:=NowLength(dif);
end;

//線分同士の距離
function NowDistance(const a,b:TLineSeg3D):Single;
var
  diff,dirA,dirB:D3DVector;
  a00,a01,a11,b0,c,det:Single;
  b1,s,t,sqrDist,tmp,invDet:Single;
begin
  diff := NowSub(a.A, b.A);
  dirA := NowSub(a.B, a.A);
  dirB := NowSub(b.B, b.A);

  a00 := NowDotProduct(dirA, dirA);
  a01 :=  - NowDotProduct(dirA, dirB);
  a11 := NowDotProduct(dirB,dirB);
  b0  := NowDotProduct(diff, dirA);
  c   := NowDotProduct(diff, diff);
  det := Abs(a00*a11 - a01*a01);

  if det >= 0 then begin
    //線分が並行ではない
    b1 := -NowDotProduct(diff,dirB);
    s  := a01*b1 - a11*b0;
    t  := a01*b0 - a00*b1;

    if s >= 0 then begin
      if s <= det then begin
        if t >= 0 then begin
          if t <= det then begin //region 0(interior) 
            invDet:=1.0 / det;
            s := s * invDet;
            t := t * invDet;
            sqrDist := s * (a00*s + a01*t + 2.0*b0) + t * (a01*s + a11*t + 2.0*b1) + c;
          end else begin //region 3 (side)
            t := 1.0;
            tmp := a01 + b0;
            if tmp >= 0 then begin
              s := 0;
              sqrDist := a11 + 2.0*b1 + c;
            end else if (-tmp >= a00) then begin
              s := 1.0;
              sqrDist := a00 + a11 + c + 2.0*(b1+tmp);
            end else begin
              s := -tmp / a00;
              sqrDist := tmp * s + a11 + 2.0*b1 + c;
            end;
          end;
        end else begin  //region 7 (side)
          t := 0;
          if b0 >= 0 then begin
            s := 0;
            sqrDist:=c;
          end else if -b0 >= a00 then begin
            s := 1.0;
            sqrDist := a00 + 2.0*b0 + c;
          end else begin
            s := -b0 / a00;
            sqrDist := b0 * s + c;
          end;
        end;

      end else begin

        if t >= 0 then begin
          if t < det then begin //region 1 (side)
            s := 1;
            tmp := a01 + b1;
            if tmp >= 0 then begin
              t := 0;
              sqrDist := a00 + 2.0*b0 + c;
            end else if -tmp >= a11 then begin
              t := 1.0;
              sqrDist := a00 + a11 + c + 2.0*(b0 + tmp);
            end else begin
              t := -tmp / a11;
              sqrDist := tmp * t + a00 + 2.0*b0 + c;
            end;
          end else begin  //region 2 (corner)
            tmp := a01 + b0;
            if -tmp <= a00 then begin
              t := 1.0;
              if tmp >= 0 then begin
                s := 0;
                sqrDist := a11 + 2.0*b1 + c;
              end else begin
                s := -tmp / a00;
                sqrDist := tmp * s + a11 + 2.0*b1 + c;
              end;
            end else begin
              s := 1;
              tmp := a01 + b1;
              if tmp >= 0 then begin
                t := 0;
                sqrDist := a00 + 2.0*b0 + c;
              end else if -tmp >= a11 then begin
                t := 1;
                sqrDist := a00 + a11 + c + 2.0*(b0 + tmp);
              end else begin
                t := -tmp / a11;
                sqrDist := tmp * t + a00 + 2.0*b0 + c;
              end;
            end;
          end;

        end else begin  // region 8 (corner)

          if -b0 < a00 then begin
            t := 0;
            if (b0 >= 0) then begin
              s:=0;
              sqrDist := c;
            end else begin
              s := -b0/a00;
              sqrDist := b0*s + c;
            end;
          end else begin
            s := 1.0;
            tmp := a01 + b1;
            if tmp >= 0 then begin
              t := 0;
              sqrDist := a00 + 2.0*b0 + c;
            end else if -tmp > a11 then begin
              t := 1;
              sqrDist := a00 + a11 + c + 2.0*(b0+tmp);
            end else begin
              t := -tmp / a11;
              sqrDist := tmp * t + a00 + 2.0*b0 + c;
            end;
          end;
        end;
      end;

    end else begin

      if t >= 0 then begin
        if t <= det then begin  // region 5 (side)
          s := 0;
          if b1 >= 0 then begin
            t := 0;
            sqrDist := c;
          end else if -b1 >= a11 then begin
            t := 1;
            sqrDist := a11 + 2.0*b1 + c;
          end else begin
            t := -b1 / a11;
            sqrDist := b1 * t + c;
          end;
        end else begin  // region 4 (corner)
          tmp := a01 + b0;
          if tmp < 0 then begin
            t := 1;
            if -tmp >= a00 then begin
              s := 1;
              sqrDist := a00 + a11 + c + 2.0*(b1 + tmp);
            end else begin
              s := -tmp / a00;
              sqrDist := tmp * s + a11 + 2.0*b1 + c;
            end;
          end else begin
            s := 0;
            if b1 >= 0 then begin
              t:=0;
              sqrDist := c;
            end else if -b1 >= a11 then begin
              t := 1;
              sqrDist := a11 + 2.0*b1 + c;
            end else begin
              t := -b1 / a11;
              sqrDist := b1 * t + c;
            end;
          end;
        end;

      end else begin  // region 6 (corner)
        if b0 < 0 then begin
          t := 0;
          if -b0 >= a00 then begin
            s := 1.0;
            sqrDist := a00 + 2.0*b0 +c;
          end else begin
            s := -b0 / a00;
            sqrDist := b0*s + c;
          end;
        end else begin
          s := 0;
          if b1 >= 0 then begin
            t := 0;
            sqrDist := c;
          end else if -b1 > a11 then begin
            t := 1;
            sqrDist := a11 + 2.0*b1 + c;
          end else begin
            t := -b1 / a11;
            sqrDist := b1*t + c;
          end;
        end;
      end;
    end;
  end else begin
    //平行なとき
    if a01 > 0 then begin
      if b0 > 0 then begin
        s := 0;
        t := 0;
        sqrDist := c;
      end else if -b0 <= a00 then begin
        s := -b0 / a00;
        t := 0;
        sqrDist := b0*s + c;
      end else begin
        b1 := - NowDotProduct(diff, dirB);
        s := 1;
        tmp := a00 + b0;
        if - tmp >= a01 then begin
          t := 1;
          sqrDist := a00 + a11 + c + 2.0*(a01+b0+b1);
        end else begin
          t := -tmp/a01;
          sqrDist := a00 + 2.0*b0 + c + t*(a11*t + 2.0*(a01+b1));
        end;
      end;
    end else begin
      //direction vectors from an acute angle
      if -b0 >= a00 then begin
        s := 1.0;
        t := 0;
        sqrDist := a00 + 2.0*b0 + c;
      end else if b0 <= 0 then begin
        s := -b0 / a00;
        t := 0;
        sqrDist := b0*s + c;
      end else begin
        b1 := -NowDotProduct(diff, dirB);
        s  := 0;
        if b0 >= -a01 then begin
          t := 1;
          sqrDist := a11 + 2.0*b1 + c;
        end else begin
          t := -b0 / a01;
          sqrDist := c + t * (2.0*b1+a11*t);
        end;
      end;
    end;
  end;


  result:=Sqrt(sqrDist);
end;

//点と平面の距離
function NowDistance(const a:D3Dvector; const plane:TPlane3D):Single;
begin
  result:=Abs(NowDotProduct(a, plane.Normal) + plane.D);
end;


//原点回りの回転(各軸を左回りに)
function NowRotX(const a:D3DVECTOR; const t:Integer):D3DVECTOR;
var
  r:Integer;
begin
  r:=t and 4095;

  result.X:=a.X;
  result.Y:=a.Y*NowCosT[r] + a.Z*NowSinT[r];
  result.Z:=-a.Y*NowSinT[r] + a.Z*NowCosT[r];
end;

function NowRotY(const a:D3DVECTOR; const t:Integer):D3DVECTOR;
var
  r:Integer;
begin
  r:=t and 4095;

  result.X:=a.X*NowCosT[r] - a.Z*NowSinT[r];
  result.Y:=a.Y;
  result.Z:=a.X*NowSinT[r] + a.Z*NowCosT[r];
end;

function NowRotZ(const a:D3DVECTOR; const t:Integer):D3DVECTOR;
var
  r:Integer;
begin
  r:=t and 4095;

  result.X:=a.X*NowCosT[r] + a.Y*NowSinT[r];
  result.Y:=-a.X*NowSinT[r] + a.Y*NowCosT[r];
  result.Z:=a.Z;
end;

//ベクトルの反転
function NowNeg(const a:D3DVECTOR):D3DVECTOR;
begin
  result.X:=-a.X;
  result.Y:=-a.Y;
  result.Z:=-a.Z;
end;

//正規化
function NowNormalize(const a:D3DVECTOR):D3DVECTOR;
var
  len:Single;
begin
  len:=NowLength(a);
  if len = 0 then begin
    result:=Vector(0,0,0);
    exit;
  end;

  result.X:=a.X / len;
  result.Y:=a.Y / len;
  result.Z:=a.Z / len;
end;

//長さを得る
function NowLength(const a:D3DVECTOR):Single;
begin
  result:=sqrt(SQR(a.X) + SQR(a.Y) + SQR(a.Z));
end;

//スケーリング
function NowScale(const a:D3DVECTOR; const scale:Single):D3DVECTOR;
begin
  result.X:=a.X * scale;
  result.Y:=a.Y * scale;
  result.Z:=a.Z * scale;
end;

//長さを k に
function NowChangeLength(const a:D3DVECTOR; const k:Single):D3DVECTOR;
var
  len:Single;
begin
  len:=NowLength(a);
  if len = 0 then begin
    result:=Vector(0,0,0);
    exit;
  end;

  result.X:=a.X * k / len;
  result.Y:=a.Y * k / len;
  result.Z:=a.Z * k / len;
end;

//posからtargetに向いた長さkのベクトルを生成
function NowLookAt(const pos:D3DVECTOR; const target:D3DVECTOR; const k:Single):D3DVECTOR;
var
  sub:D3DVECTOR;
  len:Single;
begin
  sub:=NowSub(target,pos);
  len:=NowLength(sub);
  if len = 0 then begin
    result:=Vector(0,0,0);
    exit;
  end;

  result.X:=sub.X * k / len;
  result.Y:=sub.Y * k / len;
  result.Z:=sub.Z * k / len;
end;


//ViewFrustum(x,y,z,w)→TLVertex座標(x/w,y/w,z/w,1/w)
function NowViewFrustumToScreen(const a:THgVector):THgVector;
begin
  with result do begin
    if a.W <> 0.0 then begin
      W:=1.0 / a.W;
      X:=a.X * W;
      Y:=a.Y * W;
      Z:=a.Z * W;
    end;
  end;
end;

//三つの頂点で張られる平面(三角形)の法線を返す。ただし、p1,p2,p3は法線に対して右回りと仮定する
function NowNormalFromPoints(const p1,p2,p3:D3DVector):D3DVector;
var
  a,b:D3DVector;
begin
  a:=NowSub(p1,p2);
  b:=NowSub(p3,p2);

  result:=NowNormalize(NowCrossProduct(b,a));
end;

//ベクトルvを法線nなる平面にぶつけて、跳ね返った後のベクトルを得る
function NowReflection(const v:D3DVector; normal:D3Dvector):D3DVector;
begin
  //反射ベクトルと入射ベクトルの対象軸になるのが法線 normal
  // (normal・v) normal で、対称軸が出来る

  result:=NowAdd(v, NowScale( normal, Abs(NowDotProduct(normal,v)) * 2.0 * NowLength(v)));
end;

//ベクトルvを法線nなる平面にぶつけて、跳ね返った後のベクトルを得る(弾性係数e)
function NowReflection(const v:D3DVector; normal:D3Dvector; e:Single):D3DVector;
begin
  //反射ベクトルと入射ベクトルの対象軸になるのが法線 normal
  // (normal・v) normal で、対称軸が出来る

  result:=NowAdd(v, NowScale( normal, Abs(NowDotProduct(normal,v)) * (1.0+e) * NowLength(v)));
end;


//////////////////////3D演算(クォータニオン)

//足し算
function NowAddQuat(const a:TQuaternion; const b:TQuaternion):TQuaternion;
begin
  result.X:=a.X + b.X;
  result.Y:=a.Y + b.Y;
  result.Z:=a.Z + b.Z;
  result.W:=a.W + b.W;
end;

//引き算
function NowSubQuat(const a:TQuaternion; const b:TQuaternion):TQuaternion;
begin
  result.X:=a.X - b.X;
  result.Y:=a.Y - b.Y;
  result.Z:=a.Z - b.Z;
  result.W:=a.W - b.W;
end;

//乗算
function NowMulQuat(const a:TQuaternion; const b:TQuaternion):TQuaternion;
begin
  result.X:=a.W*b.X + a.X*b.W + a.Y*b.Z - a.Z*b.Y;
  result.Y:=a.W*b.Y + b.W*a.Y + a.Z*b.X - a.X*b.Z;
  result.Z:=a.W*b.Z + b.W*a.Z + a.X*b.Y - a.Y*b.X;
  result.W:=a.W*b.W - a.X*b.X - a.Y*b.Y - a.Z*b.Z;
end;

//内積
function NowDotQuat(const a:TQuaternion; const b:TQuaternion):Single;
begin
  result:=a.w*b.w + a.x*b.x + a.y*b.y + a.z*b.z;
end;

//ノルム
function NowLengthQuat(const a:TQuaternion):Single;
begin
  result:=sqrt(a.w*a.w + a.x*a.x + a.y*a.y + a.z*a.z);
end;

//正規化
function NowNormalizeQuat(const a:TQuaternion):TQuaternion;
var
  len:Single;
begin
  len:=NowLengthQuat(a);

  if len = 0.0 then begin
    result:=Quaternion(1, 0,0,0);
    exit;
  end;

  if a.W < 0 then begin
    len:=-len;
  end;

  result.x:=a.X / len;
  result.y:=a.Y / len;
  result.z:=a.Z / len;
  result.w:=a.W / len;

end;


//線形補間
function NowLerp3D(const a:D3DVector; const b:D3DVector; const rate:Single):D3DVector;
begin
  result.x:=rate*b.x + (1.0-rate)*a.x;
  result.y:=rate*b.y + (1.0-rate)*a.y;
  result.z:=rate*b.z + (1.0-rate)*a.z;
end;


//回転を表すクォータニオンを作る
function  NowQuatFromRotation(const axis:D3DVECTOR; const t:Integer):TQuaternion;
var
  r:Integer;
  a:D3DVECTOR;
begin
  r:=(t div 2) And 4095;

  a:=NowNormalize(axis);
  result.X:=a.X * NowSinT[r];
  result.Y:=a.Y * NowSinT[r];
  result.Z:=a.Z * NowSinT[r];
  result.W:=NowCosT[r];
end;

//クォータニオンから行列を作る
function  NowMatrixFromQuat(const a:TQuaternion):D3DMATRIX;
begin
  result._11:=1.0-2.0*(a.y*a.y+a.z*a.z);
  result._12:=2.0*(a.x*a.y+a.w*a.z);
  result._13:=2.0*(a.x*a.z-a.w*a.y);
  result._14:=0;

  result._21:=2.0*(a.x*a.y-a.w*a.z);
  result._22:=1.0-2.0*(a.x*a.x+a.z*a.z);
  result._23:=2.0*(a.y*a.z+a.w*a.x);
  result._24:=0;

  result._31:=2.0*(a.x*a.z+a.w*a.y);
  result._32:=2.0*(a.y*a.z-a.w*a.x);
  result._33:=1.0-2.0*(a.x*a.x+a.y*a.y);
  result._34:=0;

  result._41:=0;
  result._42:=0;
  result._43:=0;
  result._44:=1;
end;

//axis回りの回転を表す行列を作る(QuadtFromRotationとMatrixFromQuatの合成)
function  NowMatrixFromRotation(const axis:D3DVECTOR; const t:Integer):D3DMATRIX;
begin
  result:=NowMatrixFromQuat(NowQuatFromRotation(axis,t));
end;


//Slerp rate=0.0ならaを100%、 rate=1.0ならbを100% の割合で球面線形補間
function  NowSlerpQuat(const a:TQuaternion; const b:TQuaternion; const rate:Single):TQuaternion;
var
  c,s:Single;     //a,bのなす角のCosとSin
  angle:Single; //a,bのなす角
  ra,rb:Single; //aとbの混合比

begin

  c:=NowDotQuat(a,b); //Cosθ

  //クォータニオンが正規化されてないっぽいけど、誤差が累積してる可能性があるのでc=±1として続行
  if c > 1.0 then
    c:=1.0
  else if c<-1.0 then
    c:=-1.0;

  //ArcCosで、角度を求める
  angle:=ArcTan(sqrt (1-sqr (c)) /c);

  //角度が違わないなら帰る
  if angle = 0.0 then begin
    result:=a;
    exit;
  end;

  //Sinθを求めて
  s:=Sin(angle);

  //こんな感じで
  ra:=Sin( (1-rate) * angle ) / s;
  rb:=Sin( rate * angle ) / s;

  with result do begin
    x:=a.x*ra + b.x*rb;
    y:=a.y*ra + b.y*rb;
    z:=a.z*ra + b.z*rb;
    w:=a.w*ra + b.w*rb;
  end;

 result:=NowNormalizeQuat(result);

end;

//逆クォータニオン
function NowInvQuat( srcQuat: TQuaternion ): TQuaternion;
var  n: Single;
begin
  with srcQuat do begin
    n:= ( w*w + x*x + y*y + z*z );
    if n <> 0.0 then
      Result:= Quaternion( w/n, -x/n, -y/n, -z/n )
    else
      result:=Quaternion(1,0,0,0);  //適当な値を返しちゃう
  end;
end;

//普通の線形補間
function NowLerpQuat(const a:TQuaternion; const b:TQuaternion; const rate:Single):TQuaternion;
var
  ra,rb:Single;
  dp:Single;
begin
  ra:=1-rate;

  dp:=NowDotQuat(a,b);
  if dp < 0.0 then
    rb:=-rate
  else
    rb:=rate;

  with result do begin
    w:=a.W*ra + b.W*rb;
    x:=a.X*ra + b.X*rb;
    y:=a.Y*ra + b.Y*rb;
    z:=a.Z*ra + b.Z*rb;
  end;

  result:=NowNormalizeQuat(result);
end;



//////////////////////3D演算(マトリクス)

//ベクトルを変換
function NowTransform(const a:THgVector; const m:D3DMATRIX):THgVector;
begin
  result.X:=a.X*m._11 + a.Y*m._21 + a.Z*m._31 + a.W*m._41;
  result.Y:=a.X*m._12 + a.Y*m._22 + a.Z*m._32 + a.W*m._42;
  result.Z:=a.X*m._13 + a.Y*m._23 + a.Z*m._33 + a.W*m._43;
  result.W:=a.X*m._14 + a.Y*m._24 + a.Z*m._34 + a.W*m._44;
end;

function NowTransform(const a:D3DVector; const m:D3DMATRIX):D3DVector;
begin
  result.X:=a.X*m._11 + a.Y*m._21 + a.Z*m._31 + m._41;
  result.Y:=a.X*m._12 + a.Y*m._22 + a.Z*m._32 + m._42;
  result.Z:=a.X*m._13 + a.Y*m._23 + a.Z*m._33 + m._43;
end;

//一気にベクトルを変換
procedure NowTransformBatch(const src:Pointer; const m:D3DMATRIX; const num:Integer; const span:Integer);
var
  i:Integer;
  a:THgVector;
  pSrc:PHgVector;
begin
  pSrc:=src;

  for i:=1 to num do begin
    a:=pSrc^;
    pSrc^.X:=a.X*m._11 + a.Y*m._21 + a.Z*m._31 + a.W*m._41;
    pSrc^.Y:=a.X*m._12 + a.Y*m._22 + a.Z*m._32 + a.W*m._42;
    pSrc^.Z:=a.X*m._13 + a.Y*m._23 + a.Z*m._33 + a.W*m._43;
    pSrc^.W:=a.X*m._14 + a.Y*m._24 + a.Z*m._34 + a.W*m._44;
    Inc(PBYTE(pSrc),span);
  end;
end;

//カメラとターゲットの位置から、ビュー行列を作る
function NowViewMatrix(const from:D3DVECTOR; const at:D3DVECTOR; const world_up:D3DVECTOR; const roll:Integer):D3DMATRIX;
var
  up,right,direction:D3DVECTOR; //頭の中心→頭のてっぺん、頭の中心→右耳、視線の方向
begin
  direction:=NowNormalize(NowSub(at,from));
  right := NowCrossProduct(world_up, direction);
  up:=NowCrossProduct(direction,right);

  right:=NowNormalize(right);
  up:=NowNormalize(up);


  result:=NowIdentityMatrix;
  result._11 := right.x;
  result._21 := right.y;
  result._31 := right.z;
  result._12 := up.x;
  result._22 := up.y;
  result._32 := up.z;
  result._13 := direction.x;
  result._23 := direction.y;
  result._33 := direction.z;

  result._41:=-NowDotProduct(right, from);
  result._42:=-NowDotProduct(up, from);
  result._43:=-NowDotProduct(direction, from);

  if roll <> 0 then begin
    result:=NowCompositeMatrix(NowRotZMatrix(roll),result);
  end;

end;


//姿勢・位置からビュー行列を作る
function NowViewMatixFromOrientation(const from:D3DVECTOR; const dir:D3DVECTOR; const up:D3DVECTOR):D3DMATRIX;
var
  vx,vy,vz:D3DVECTOR;
begin
  vx:=NowNormalize(NowCrossProduct(up,dir));
  vy:=NowNormalize(up);
  vz:=NowNormalize(dir);

  result:=NowIdentityMatrix;
  with result do begin
    _11:=vx.x;    _12:=vy.x;    _13:=vz.x;
    _21:=vx.y;    _22:=vy.y;    _23:=vz.y;
    _31:=vx.z;    _32:=vy.z;    _33:=vz.z;
    _41:=-from.x; _42:=-from.y; _43:=-from.z;
  end;
end;


//透視射影行列(視野座標→除算前正規化座標)を作る
function NowProjectFrustum(const halfWidth:Single; const halfHeight:Single; const nearZ:Single; const farZ:Single):D3DMATRIX;
begin
  ZeroMemory(@result,sizeof(D3DMATRIX));

  with result do begin
    _11:=(nearZ+1) * halfHeight;

    _22:=-_11;

    _31:=halfWidth;
    _32:=halfHeight;
    _33:=farZ / (farZ-nearZ);
    _34:=1.0;

    _43:=_33 * (-nearZ);
  end;
end;


function NowRotXMatrix(const t:Integer):D3DMATRIX;
var
  t2:Integer;
begin
  t2:=t And 4095;

  with result do begin
    _11:=1; _12:=0;             _13:=0;           _14:=0;
    _21:=0; _22:=NowCosT[t2];   _23:=NowSinT[t2]; _24:=0;
    _31:=0; _32:=-NowSinT[t2];  _33:=NowCosT[t2]; _34:=0;
    _41:=0; _42:=0;             _43:=0;           _44:=1;
  end;
end;

function NowRotYMatrix(const t:Integer):D3DMATRIX;
var
  t2:Integer;
begin
  t2:=t And 4095;

  with result do begin
    _11:=NowCosT[t2]; _12:=0; _13:=-NowSinT[t2];  _14:=0;
    _21:=0;           _22:=1; _23:=0;             _24:=0;
    _31:=NowSinT[t2]; _32:=0; _33:=NowCosT[t2];   _34:=0;
    _41:=0;           _42:=0; _43:=0;             _44:=1;
  end;
end;

function NowRotZMatrix(const t:Integer):D3DMATRIX;
var
  t2:Integer;
begin
  t2:=t And 4095;

  with result do begin
    _11:=NowCosT[t2];   _12:=NowSinT[t2]; _13:=0; _14:=0;
    _21:=-NowSinT[t2];  _22:=NowCosT[t2]; _23:=0; _24:=0;
    _31:=0;             _32:=0;           _33:=1; _34:=0;
    _41:=0;             _42:=0;           _43:=0; _44:=1;
  end;
end;

//スケーリング行列
function NowScaleMatrix(const x,y,z:Single):D3DMATRIX;
begin
  with result do begin
    _11:=x;   _12:=0; _13:=0; _14:=0;
    _21:=0;   _22:=y; _23:=0; _24:=0;
    _31:=0;   _32:=0; _33:=z; _34:=0;
    _41:=0;   _42:=0; _43:=0; _44:=1;
  end;
end;

//平行移動行列
function NowTranslateMatrix(const x,y,z:Single):D3DMATRIX;
begin
  with result do begin
    _11:=1;   _12:=0; _13:=0; _14:=0;
    _21:=0;   _22:=1; _23:=0; _24:=0;
    _31:=0;   _32:=0; _33:=1; _34:=0;
    _41:=x;   _42:=y; _43:=z; _44:=1;
  end;
end;


//逆行列(変換行列用の、手抜き計算)
function NowInvMatrix(const m:D3DMATRIX):D3DMATRIX;
begin
  with result do begin
    _11:=m._11;   _12:=m._21;   _13:=m._31;   _14:=0;
    _21:=m._12;   _22:=m._22;   _23:=m._32;   _24:=0;
    _31:=m._13;   _32:=m._23;   _33:=m._33;   _34:=0;

    _41:= - (m._11*m._41 + m._12*m._42 + m._13*m._43);
    _42:= - (m._21*m._41 + m._22*m._42 + m._23*m._43);
    _43:= - (m._31*m._41 + m._32*m._42 + m._33*m._43);
    _44:= 1;
  end;
end;

//逆行列(回転だけを格納していると決めウチ)
function NowInvMatrixRot(const m:D3DMATRIX):D3DMATRIX;
begin
  with result do begin
    _11:=m._11;   _12:=m._21;   _13:=m._31;   _14:=0;
    _21:=m._12;   _22:=m._22;   _23:=m._32;   _24:=0;
    _31:=m._13;   _32:=m._23;   _33:=m._33;   _34:=0;

    _41:=0;
    _42:=0;
    _43:=0;
    _44:=1;
  end;
end;


//転置行列
function NowTMatrix(const m:D3DMATRIX):D3DMATRIX;
begin
  with result do begin
    _11:=m._11;   _12:=m._21;   _13:=m._31;   _14:=m._41;
    _21:=m._12;   _22:=m._22;   _23:=m._32;   _24:=m._42;
    _31:=m._13;   _32:=m._23;   _33:=m._33;   _34:=m._43;
    _41:=m._14;   _42:=m._24;   _43:=m._34;   _44:=m._44;
  end;
end;

//回転部分だけ抜き出し
function NowExtractRotation(const m:D3DMATRIX):D3DMATRIX;
begin
  with result do begin
    _11:=m._11;   _12:=m._12;   _13:=m._13;   _14:=0;
    _21:=m._21;   _22:=m._22;   _23:=m._23;   _24:=0;
    _31:=m._31;   _32:=m._32;   _33:=m._33;   _34:=0;
    _41:=0;       _42:=0;       _43:=0;       _44:=1;
  end;
end;

//平行移動部分だけ抜き出し
function NowExtractTranslation(const m:D3DMATRIX):D3DMATRIX;
begin
  result:=NowIdentityMatrix;
  with result do begin
    _41:=m._41;   _42:=m._42;   _43:=m._43;
  end;
end;


//マトリクスの合成
function NowCompositeMatrix(const m:D3DMATRIX; const n:D3DMATRIX):D3DMATRIX;
begin

  result._11:= m._11*n._11 + m._12*n._21 + m._13*n._31 + m._14*n._41;
  result._12:= m._11*n._12 + m._12*n._22 + m._13*n._32 + m._14*n._42;
  result._13:= m._11*n._13 + m._12*n._23 + m._13*n._33 + m._14*n._43;
  result._14:= m._11*n._14 + m._12*n._24 + m._13*n._34 + m._14*n._44;

  result._21:= m._21*n._11 + m._22*n._21 + m._23*n._31 + m._24*n._41;
  result._22:= m._21*n._12 + m._22*n._22 + m._23*n._32 + m._24*n._42;
  result._23:= m._21*n._13 + m._22*n._23 + m._23*n._33 + m._24*n._43;
  result._24:= m._21*n._14 + m._22*n._24 + m._23*n._34 + m._24*n._44;

  result._31:= m._31*n._11 + m._32*n._21 + m._33*n._31 + m._34*n._41;
  result._32:= m._31*n._12 + m._32*n._22 + m._33*n._32 + m._34*n._42;
  result._33:= m._31*n._13 + m._32*n._23 + m._33*n._33 + m._34*n._43;
  result._34:= m._31*n._14 + m._32*n._24 + m._33*n._34 + m._34*n._44;

  result._41:= m._41*n._11 + m._42*n._21 + m._43*n._31 + m._44*n._41;
  result._42:= m._41*n._12 + m._42*n._22 + m._43*n._32 + m._44*n._42;
  result._43:= m._41*n._13 + m._42*n._23 + m._43*n._33 + m._44*n._43;
  result._44:= m._41*n._14 + m._42*n._24 + m._43*n._34 + m._44*n._44;

end;

//m.nの回転成分だけ抜き出して合成する
function NowCompositeRotMatrix(const m:D3DMATRIX; const n:D3DMATRIX):D3DMATRIX;
begin
  result._11:= m._11*n._11 + m._12*n._21 + m._13*n._31;
  result._12:= m._11*n._12 + m._12*n._22 + m._13*n._32;
  result._13:= m._11*n._13 + m._12*n._23 + m._13*n._33;
  result._14:= 0;

  result._21:= m._21*n._11 + m._22*n._21 + m._23*n._31;
  result._22:= m._21*n._12 + m._22*n._22 + m._23*n._32;
  result._23:= m._21*n._13 + m._22*n._23 + m._23*n._33;
  result._24:= 0;

  result._31:= m._31*n._11 + m._32*n._21 + m._33*n._31;
  result._32:= m._31*n._12 + m._32*n._22 + m._33*n._32;
  result._33:= m._31*n._13 + m._32*n._23 + m._33*n._33;
  result._34:= 0;

  result._41:= m._41*n._11 + m._42*n._21 + m._43*n._31;
  result._42:= m._41*n._12 + m._42*n._22 + m._43*n._32;
  result._43:= m._41*n._13 + m._42*n._23 + m._43*n._33;
  result._44:= 1;
end;


//乗算
function NowCompositeMatrix33(const a:TMatrix33; const b:TMatrix33):TMatrix33;
begin
  result._11:= a._11*b._11 + a._12*b._21 + a._13*b._31;
  result._12:= a._11*b._12 + a._12*b._22 + a._13*b._32;
  result._13:= a._11*b._13 + a._12*b._23 + a._13*b._33;

  result._21:= a._21*b._11 + a._22*b._21 + a._23*b._31;
  result._22:= a._21*b._12 + a._22*b._22 + a._23*b._32;
  result._23:= a._21*b._13 + a._22*b._23 + a._23*b._33;

  result._31:= a._31*b._11 + a._32*b._21 + a._33*b._31;
  result._32:= a._31*b._12 + a._32*b._22 + a._33*b._32;
  result._33:= a._31*b._13 + a._32*b._23 + a._33*b._33;
end;

//逆行列
function NowInvMatrix33(const a:TMatrix33; var valid:Boolean):TMatrix33;
var
  det,rdet:Single;
begin

  //まず、行列式の計算
  det:=NowDeterminant33(a);

  //行列式がゼロじゃ、ねぇ
  if det = 0.0 then begin
    valid:=False;
    exit;
  end;

  {
     -1     1     |   EI-FH  -(BI-HC)   BF-EC  |
    M   = ----- . | -(DI-FG)   AI-GC  -(AF-DC) |
          det M   |   DH-GE  -(AH-GB)   AE-BD  |
  }
  rdet:=1.0/det;
  with a do begin
    result._11 := (e*i-f*h) * rdet;
    result._12 := (-(b*i-h*c)) * rdet;
    result._13 := (b*f-e*c) * rdet;

    result._21 := (-(d*i-f*g)) * rdet;
    result._22 := (a*i-g*c) * rdet;
    result._23 := (-(a*f-d*c)) * rdet;

    result._31 := (d*h-g*e) * rdet;
    result._32 := (-(a*h-g*b)) * rdet;
    result._33 := (a*e-d*b) * rdet;
  end;

  valid:=True;
end;

//一次変換
function NowTransform33(const a:D3DVector; const m:TMatrix33):D3DVector;
begin
  result.x:=a.x*m._11 + a.y*m._21 + a.z*m._31;
  result.y:=a.x*m._12 + a.y*m._22 + a.z*m._32;
  result.z:=a.x*m._13 + a.y*m._32 + a.z*m._33;
end;

//行列式
function NowDeterminant33(const m:TMatrix33):Single;
begin
  with m do
    result:=a * (e*i - h*f) - b * (d*i - g*f) + c * (d*h - g*e);
end;


//pはboxに含まれるか？
function NowPointInBox(const p:D3DVector; const box:TBox3D):Boolean;
begin
  result:= (box.min.X <= p.X) And (p.X <= box.max.X) And
          (box.min.Y <= p.Y) And (p.Y <= box.max.Y) And
          (box.min.Z <= p.Y) And (p.Z <= box.max.Z);
end;


//aはbに包まれるか？
function NowBoxInBox(const a:TBox3D; const b:TBox3D):Boolean;
begin
  result:=
    NowPointInBox( a.min, b) And
    NowPointInBox( a.max, b);
end;

//pは 三角錐 t1〜t4 に含まれるか？
function NowPointInTetrahedron(const p:D3DVector; const t1,t2,t3,t4:D3DVector):Boolean;
var
  a,b,c,pp,t:D3DVector;
  m:TMatrix33;
  valid:Boolean;
begin
  //tri[0]を基準点として考える
  a:=NowSub(t2,t1);
  b:=NowSub(t3,t1);
  c:=NowSub(t4,t1);
  pp:=NowSub(p,t1);

  //pp =  t M → t = pp M^-1
  m.a:=a.x;   m.b:=a.y;   m.c:=a.z;
  m.d:=b.x;   m.e:=b.y;   m.f:=b.z;
  m.g:=c.x;   m.h:=c.y;   m.i:=c.z;

  //逆行列の計算
  m:=NowInvMatrix33(m,valid);

  //三角錐がひしゃげてるよ
  if not valid then begin
    result:=False;
    exit;
  end;

  t:=NowTransform33(pp, m);

  // t.x〜t.zが0〜1の範囲、かつ、t.x+t.y+t.zが0〜1の範囲なら、含まれてる

  result:= (0<=t.x) And (t.x<=1) And (0<=t.y) And (t.y<=1) And (0<=t.z) And (t.z<=1) And
          (0 <= t.x+t.y+t.z) And (t.x+t.y+t.z <= 1);

  //result:= (-0.01<=t.x) And (t.x<=1.01) And (-0.01<=t.y) And (t.y<=1.01) And (-0.01<=t.z) And (t.z<=1.01) And
  //        (-0.01 <= t.x+t.y+t.z) And (t.x+t.y+t.z <= 1.01);

end;

//三角錐 tri と、boxは、交差ないし包含しあうか
{
function NowIntersectTetrahedronBox(const tri:Array of D3DVector; const box:TBox3D):Boolean;
var
  i,j:Integer;
begin
  result:=False;

  //triの稜線がboxと一本でも交差する
  For i:=0 to 2 do begin
    For j:=i+1 to 3 do begin
      if NowIntersectLineBox(tri[i],tri[j], box) then begin
        result:=True;
        exit;
      end;
    end;
  end;


  //boxがtriの点を一つでも包んでいる
  for i:=0 to 3 do
  if NowPointInBox(tri[i],box) then begin
    result:=True;
    exit;
  end;

  //triがboxの点を一つでも包んでいる
  if NowPointInTetrahedron(Vector(box.min.x, box.min.y, box.min.z), tri) Or
    NowPointInTetrahedron(Vector(box.max.x, box.min.y, box.min.z), tri) Or
    NowPointInTetrahedron(Vector(box.min.x, box.max.y, box.min.z), tri) Or
    NowPointInTetrahedron(Vector(box.max.x, box.max.y, box.min.z), tri) Or
    NowPointInTetrahedron(Vector(box.min.x, box.min.y, box.max.z), tri) Or
    NowPointInTetrahedron(Vector(box.max.x, box.min.y, box.max.z), tri) Or
    NowPointInTetrahedron(Vector(box.min.x, box.max.y, box.max.z), tri) Or
    NowPointInTetrahedron(Vector(box.max.x, box.max.y, box.max.z), tri) then begin
    result:=True;
    exit;
  end;

end;
}

//頂点が三角柱(triを底面とし、高さheight)にめり込んでいるか
function NowPointInPillar(const p:D3DVector; const tri:TFacedTriangle3D; const height:Single):Boolean;
var
  dif:D3Dvector;  //
  prod:Single;
  h:D3DVector;    //pから「triを含む平面」への垂線の脚

  d1,d2,d3:Single;
  v1,v2,v3:D3DVector;
  prod1,prod2,prod3:D3DVector;
begin
  result:=False;

  dif:=NowSub(tri.A, p);
  prod:=NowDotProduct(dif, tri.N);

  //三角柱の上の面よりさらに上にpがある
  //三角柱の下の面よりさらに下にpがある
  if (prod < 0) or (prod > height) then
    exit;


  //三角形を含む平面と線分がヒットしてる？
  //hpがt1,t2,t3内部にあるか判定
  h:=NowAdd(p, NowScale(tri.N, prod));
  v1:=NowSub(tri.A,h);
  v2:=NowSub(tri.B,h);
  v3:=NowSub(tri.C,h);

  prod1:=NowCrossProduct(v1,v2);
  prod2:=NowCrossProduct(v2,v3);
  prod3:=NowCrossProduct(v3,v1);

  d1:=NowDotProduct(prod1, prod2);
  d2:=NowDotProduct(prod2, prod3);
  d3:=NowDotProduct(prod3, prod1);

  //d1〜3が全部同符号か、一個でもゼロが入ってたら、OK
  if d1*d2*d3 = 0.0 then begin
    result:=True;
  end else begin
    result:=( ((d1>0) And (d2>0) And (d3>0)) Or ((d1<0) And (d2<0) And (d3<0)));
  end;

end;


//線分abと平面sの交差判定
function NowIntersectLinePlane(var dest:D3DVector; const seg:TLineSeg3D; const s:TPlane3D):Boolean;
var
  n:D3Dvector;
  da,db:Single; //面の法線と、a,bとの内積
  det:Single;
  t:Single;
begin
  n:=NowNormalize(s.Normal);

  //面の裏と表、どっちに点がある？
  da:=NowDotProduct(seg.A, n);
  db:=NowDotProduct(seg.B, n);

  det:=(da + s.D)*(db + s.D);
  if det > 0 then begin
    //線分が面を跨いでないなら、交点なし
    result:=False;
    exit;
  end else if det = 0.0 then begin
    //aかbのどちらかが、面の上にある
    if da+s.D = 0.0 then
      dest:=seg.A
    else if db+s.D = 0.0 then
      dest:=seg.B;

    result:=True;
  end else begin
    //交点あり
    t := (-s.D - db) / (da - db);
    dest:=NowAdd( NowScale(seg.A, t), NowScale(seg.B, (1-t)));

    result:=True;
  end;

end;



//三角形と線分abの交差判定
function NowIntersectLineTriangle(const seg:TLineSeg3D; const t1:D3DVector; const t2:D3DVector; const t3:D3DVector):Boolean;
var
  n:D3DVector;
  hp:D3DVector; //線分と三角形を含む平面の交点
  s:TPlane3D;

  v1,v2,v3:D3DVector; //hpを基準とした、v1,v2,v3;
  prod1,prod2,prod3:D3DVector;

  d1,d2,d3:Single;    //内外判定用
begin

  //三角形を含む平面を作る
  s:=NowPlaneFromPoints(t1,t2,t3);

  //三角形を含む平面と線分がヒットしてる？
  result:=False;
  if NowIntersectLinePlane(hp, seg, s) then begin
    //hpがt1,t2,t3内部にあるか判定

    v1:=NowSub(t1,hp);
    v2:=NowSub(t2,hp);
    v3:=NowSub(t3,hp);

    prod1:=NowCrossProduct(v1,v2);
    prod2:=NowCrossProduct(v2,v3);
    prod3:=NowCrossProduct(v3,v1);

    d1:=NowDotProduct(prod1, prod2);
    d2:=NowDotProduct(prod2, prod3);
    d3:=NowDotProduct(prod3, prod1);

    //d1〜3が全部同符号か、一個でもゼロが入ってたら、OK
    if d1*d2*d3 = 0.0 then begin
      result:=True;
    end else begin
      result:=( ((d1>0) And (d2>0) And (d3>0)) Or ((d1<0) And (d2<0) And (d3<0)));
    end;

  end;

end;


//ポリゴン同士の当たり判定
function NowIntersectPolygons(const poly1:Array of D3DVector; const poly2:Array of D3DVector):Boolean;
var
  i,j:Integer;
begin

  //poly1の各辺と、poly2の各三角形が交わるか判定。どれかひとつでも交わってるなら、True
  For i:=0 to High(poly2) do begin
    For j:=0 to High(poly1)-1 do begin
      if NowIntersectLineTriangle(LineSeg3D(poly1[j],poly1[j+1]), poly2[i],poly2[(i+1) mod (High(Poly2)+1)],poly2[(i+2) mod (High(Poly2)+1)]) then begin
        result:=True;
        exit;
      end;
    end;
  end;

  //poly2の各辺と、poly1の各三角形が交わるか判定。どれかひとつでも交わってるなら、True
  For i:=0 to High(poly1) do begin
    For j:=0 to High(poly2)-1 do begin
      if NowIntersectLineTriangle(LineSeg3D(poly2[j],poly2[j+1]), poly1[i],poly1[(i+1) mod (High(Poly1)+1)],poly1[(i+2) mod (High(Poly1)+1)]) then begin
        result:=True;
        exit;
      end;
    end;
  end;

  result:=False;

end;


//線分abとboxの交差判定
function NowIntersectLineBox(const seg:TLineSeg3D; const box:TBox3D):Boolean;
var
  tmin,tmax:D3DVector;  //xyz各要素について、線分がboxに含まれるためにtの取りうる範囲
  amb:D3DVector;        //計算用
  btmin,btmax:Single;       //tの取りうる範囲

  //bmin〜bmaxを、min〜maxで狭める。これいじょう狭められないなら、result:=False;
  function AndBounds(var bmin:Single; var bmax:Single; min,max:Single):Boolean;
  begin

    if bmin<min then
      bmin:=min;
    if bmax>max then
      bmax:=max;

    result := (bmax>=bmin);
  end;
begin

  //abが線分じゃなくて点だった場合
  if ((seg.a.x) = (seg.b.x)) And ((seg.a.y) = (seg.b.y)) And ((seg.a.z) = (seg.b.z)) then begin
    result:=NowPointInBox(seg.a, box);
    exit;
  end;


  //(大文字はベクトルとして表記)
  //線分は、tA + (1-t)B で定義される。
  // 0<=t<=1の範囲で以上の関係を満たす点が箱の中に存在しうるとき、箱の中に線分abが存在する。

  amb:=NowSub(seg.a,seg.b); //A-B

  if amb.x > 0 then begin
    tmin.x := (box.min.x-seg.b.x) / (amb.x);
    tmax.x := (box.max.x-seg.b.x) / (amb.x);
  end else if amb.x = 0.0 then begin
    if (box.min.x <=seg.a.x) And (seg.a.x <= box.max.x) then begin
      tmin.x := 0;  //この要素が線分の両端で変化しない場合
      tmax.x := 1;
    end else begin
      result:=False;
      exit;
    end;
  end else begin
    tmin.x := (box.max.x-seg.b.x) / (amb.x);
    tmax.x := (box.min.x-seg.b.x) / (amb.x);
  end;

  if amb.y > 0 then begin
    tmin.y := (box.min.y-seg.b.y) / (amb.y);
    tmax.y := (box.max.y-seg.b.y) / (amb.y);
  end else if amb.y = 0.0 then begin
    if (box.min.y <=seg.a.y) And (seg.a.y <= box.max.y) then begin
      tmin.y := 0;  //この要素が線分の両端で変化しない場合
      tmax.y := 1;
    end else begin
      result:=False;
      exit;
    end;
  end else begin
    tmin.y := (box.max.y-seg.b.y) / (amb.y);
    tmax.y := (box.min.y-seg.b.y) / (amb.y);
  end;

  if amb.z > 0 then begin
    tmin.z := (box.min.z-seg.b.z) / (amb.z);
    tmax.z := (box.max.z-seg.b.z) / (amb.z);
  end else if amb.z = 0.0 then begin
    if (box.min.z <=seg.a.z) And (seg.a.z <= box.max.z) then begin
      tmin.z := 0;  //この要素が線分の両端で変化しない場合
      tmax.z := 1;
    end else begin
      result:=False;
      exit;
    end;
  end else begin
    tmin.z := (box.max.z-seg.b.z) / (amb.z);
    tmax.z := (box.min.z-seg.b.z) / (amb.z);
  end;


  //tmin,tmaxで示される範囲を、x,y,z各要素についてAndを取る
  btmin:=0;
  btmax:=1;

  result:=AndBounds(btmin,btmax, tmin.x,tmax.x);
  if not result then
    exit;

  result:=AndBounds(btmin,btmax, tmin.y,tmax.y);
  if not result then
    exit;

  result:=AndBounds(btmin,btmax, tmin.z,tmax.z);
  if not result then
    exit;


end;


//四角錐とboxは交差しあうか？(四角錐を円錐で近似して計算)
function NowIntersectPyramidBox(const top:D3DVector; const lt:D3DVector; const rt:D3DVector; const lb:D3DVector; const rb:D3DVector; const box:TBox3D):Boolean;
var
  dir:D3DVector;  //四角錐の頂点から底面への垂線
  thre_angle:Single;    //円錐の内外判定用
  thre_dist:Single;

  ridge:Array[0..3] of D3Dvector; //稜線
  bm,bn:D3DVector;  //ボックスの頂点(円錐の頂点からの相対座標)
  mbox:TBox3D;      //円錐の頂点からの相対座標で示したボックス

  //ptが円錐に入ってるか？
  function PointInCone(pt:D3DVector):Boolean;
  var
    d:Single;
    len:Single;
  begin
    d:=NowDotProduct(NowNormalize(pt), dir);

    if d >= thre_angle then begin
      //角度はOK、お次は距離
      len:=NowDotProduct(dir,pt);

      if len <= thre_dist then begin
        result:=True;
        exit;
      end;
    end;

    result:=False;
  end;
begin


  result:=False;

  {boxに四角錐の頂点うち一点でも含まれるか？}
  {
  if NowPointInBox(top, box) Or
     NowPointInBox(lt, box) Or
     NowPointInBox(rt, box) Or
     NowPointInBox(lb, box) Or
     NowPointInBox(rb, box) then begin
     result:=True;
     exit;
  end;
  }


  {円錐にboxの頂点のうち一点でも含まれるか？}

  //頂点を基準点にする
  ridge[0]:=NowSub(rt,top);
  ridge[1]:=NowSub(lt,top);
  ridge[2]:=NowSub(rb,top);
  ridge[3]:=NowSub(lb,top);

  //四角錐の向きを計算…頂点→底面の中心 なベクトル
  dir:=NowNormalize(NowScale(NowAdd(ridge[0],ridge[3]),0.5));

  //内外判定用の値を計算…四角錐を包むような円錐の広がりの角度をφとして、cosφ=thre_d
  thre_angle:=NowDotProduct(dir, NowNormalize(ridge[0]));
  thre_dist:=NowLength(NowScale(NowAdd(ridge[0],ridge[3]),0.5));

  //ボックスの各頂点をチェック
  bm:=NowSub(box.max, top);
  bn:=NowSub(box.min, top);

  if PointInCone(Vector(bm.x,bm.y,bm.z)) Or
      PointInCone(Vector(bn.x,bm.y,bm.z)) Or
      PointInCone(Vector(bm.x,bn.y,bm.z)) Or
      PointInCone(Vector(bn.x,bn.y,bm.z)) Or
      PointInCone(Vector(bm.x,bm.y,bn.z)) Or
      PointInCone(Vector(bn.x,bm.y,bn.z)) Or
      PointInCone(Vector(bm.x,bn.y,bn.z)) Or
      PointInCone(Vector(bn.x,bn.y,bn.z)) then begin
    result:=True;
    exit;
  end;
  

  {互いに貫きあっている場合が考えられるので、四角錐の稜線がboxと交わらないかチェック}
  mbox.min:=bn;
  mbox.max:=bm;
  if NowIntersectLineBox(LineSeg3D(Vector(0,0,0), ridge[0]), mbox) Or
    NowIntersectLineBox(LineSeg3D(Vector(0,0,0), ridge[1]), mbox) Or
    NowIntersectLineBox(LineSeg3D(Vector(0,0,0), ridge[2]), mbox) Or
    NowIntersectLineBox(LineSeg3D(Vector(0,0,0), ridge[3]), mbox) Or
    NowIntersectLineBox(LineSeg3D(ridge[0],ridge[1]), mbox) Or
    NowIntersectLineBox(LineSeg3D(ridge[1],ridge[3]), mbox) Or
    NowIntersectLineBox(LineSeg3D(ridge[0],ridge[2]), mbox) Or
    NowIntersectLineBox(LineSeg3D(ridge[2],ridge[3]), mbox) then begin
    result:=True;
    exit;
  end;


end;

//boxは視界の中にあるか(viewは、boxの属する座標系にある頂点を、透視変換まで行う行列。clipの中に変換後のboxが入っていればTrue)
function NowIntersectViewBox(const view:D3DMatrix; const box:TBox3D; const clip:TRect):Boolean;
var
  face:Array [0..5,0..3] of THgVector;  //boxの面
  p:Array [0..7] of THgVector;          //boxの頂点
  i:Integer;
begin

  p[0]:=HgVector(box.min.x, box.min.y, box.min.z, 1);
  p[1]:=HgVector(box.max.x, box.min.y, box.min.z, 1);
  p[2]:=HgVector(box.min.x, box.max.y, box.min.z, 1);
  p[3]:=HgVector(box.max.x, box.max.y, box.min.z, 1);
  p[4]:=HgVector(box.min.x, box.min.y, box.max.z, 1);
  p[5]:=HgVector(box.max.x, box.min.y, box.max.z, 1);
  p[6]:=HgVector(box.min.x, box.max.y, box.max.z, 1);
  p[7]:=HgVector(box.max.x, box.max.y, box.max.z, 1);

  //全部変換
  NowTransformBatch(@p[0], view, 8, Sizeof(THgVector));

  //一個でもW=0の頂点があるなら、視界に入ってるとみなす
  For i:=0 to 7 do begin
    if p[i].W = 0 then begin
      result:=True;
      exit;
    end;
  end;

  //正規化除算
  For i:=0 to 7 do
    p[i]:=NowViewFrustumToScreen(p[i]);

  //boxを面に分解する
  face[0,0]:=p[0];  face[0,1]:=p[2];  face[0,2]:=p[3];  face[0,3]:=p[1];
  face[1,0]:=p[3];  face[1,1]:=p[7];  face[1,2]:=p[5];  face[1,3]:=p[1];
  face[2,0]:=p[6];  face[2,1]:=p[2];  face[2,2]:=p[0];  face[2,3]:=p[4];
  face[3,0]:=p[7];  face[3,1]:=p[6];  face[3,2]:=p[4];  face[3,3]:=p[5];
  face[4,0]:=p[6];  face[4,1]:=p[7];  face[4,2]:=p[3];  face[4,3]:=p[2];
  face[5,0]:=p[5];  face[5,1]:=p[4];  face[5,2]:=p[0];  face[5,3]:=p[1];


  //各面がclip内に入っているかチェック
  For i:=0 to 5 do begin

  end;


end;


//カプセル同士は交差するか？
function NowIntersectCapsule(const a,b:TCapsule3D):Boolean;
begin
  result := (NowDistance(a.Axis, b.Axis) <= (a.Radius + B.Radius));
end;

{$ENDIF}






initialization

//なぜか一部の環境では、ここが呼ばれないので、自主的に呼び出すべし…
// コンパイラのバージョンも上がったし、いいかな
NowInit;

end.
