球

/********************************************************
*文件名:  Sphere.h
*文件描述:  创建一个用三次样条曲线光滑化的不规则球体
*创建人:  陈泽丹
*创建时间:  2010-04-03
*版本:   1.0
*备注:   用角度和R构成三次样条曲线,这样X可单调递增
    如果模型渲染时某些地方有闪烁现象,
    可能是给的区域有重叠的三角形.
    光照球时,不强烈但细微的差距可以肉眼识别(通过阴影)
    色彩球时,强烈但细微的差距看不出差别
*********************************************************/

#pragma once
#include "D3DPosition.h"

struct LightSphereVertex
{
 LightSphereVertex(){}

 LightSphereVertex(float x, float y, float z, float nx, float ny, float nz, D3DCOLOR dif, D3DCOLOR spe)
 {
  _x  = x;  _y  = y; _z  = z;
  _nx = nx; _ny = ny; _nz = nz;
  _dif = dif;
  _spe = spe;
 }
 float  _x,  _y,  _z;
 float _nx, _ny, _nz;
 D3DCOLOR _dif;
 D3DCOLOR _spe;
 static const DWORD FVF;
};

 

class SphereEx:public virtual D3DPosition
{
public:
 SphereEx(void);
 virtual ~SphereEx(void);

 //建立模型
 bool SetupSphere(IDirect3DDevice9* pDevice,
      const int in_iNewTS, const int in_iNewTE,
      const int in_iNewTStep,
      const int in_iNewPS, const int in_iNewPE,
      const int in_iNewPStep,
      const int in_iOldTS, const int in_iOldTE,
      const int in_iOldTStep,
      const int in_iOldPS, const int in_iOldPE,
      const int in_iOldPStep,
      const float** in_pfKnowData);

 //清空模型
 void CleanUp();

 //设置状态
 void SetState(IDirect3DDevice9* pDevice);

 //渲染
 void Rander(IDirect3DDevice9* pDevice);

 

protected:
 //设置数值
 void Transform(float& x, float& y, float& z,
  const int in_iTA, const int in_iPA, const float in_fValue);

 void AddTriangle_UpLeft(LightSphereVertex* v, int& in_ID,
  float*** fRet,int in_iRowCenter, int in_iLineCenter,
  int in_iRowNum, int in_iLineNum, bool in_bIsClockwise);
 void AddTriangle_UpRight(LightSphereVertex* v, int& in_ID,
  float*** fRet,int in_iRowCenter, int in_iLineCenter,
  int in_iRowNum, int in_iLineNum, bool in_bIsClockwise);
 void AddTriangle_DownRight(LightSphereVertex* v, int& in_ID,
  float*** fRet,int in_iRowCenter, int in_iLineCenter,
  int in_iRowNum, int in_iLineNum, bool in_bIsClockwise);
 void AddTriangle_DownLeft(LightSphereVertex* v, int& in_ID,
  float*** fRet,int in_iRowCenter, int in_iLineCenter,
  int in_iRowNum, int in_iLineNum, bool in_bIsClockwise);

 enum
 {
  LIGHT_DIRECT_ID = 0,
  XYZR = 4,
 };
 float m_fDirLightPos[3];
 float m_fRMAX;
 float m_fRMIN;

 int m_iTriangleCount;
 D3DLIGHT9 m_DirLight;
 IDirect3DVertexBuffer9* m_pVB;
};

 

////////////////////////////////////////////////////////////////////////////////////////

 

/****************************************
*文件名:  Sphere.cpp
*创建人:  陈泽丹
*创建时间:  2010-04-03
*版本:   1.0
*****************************************/

#include "stdafx.h"
#include "SphereEx.h"
#include "Math_D3D.h"
#include "Math_BiCubic.h"
#include <time.h>

const DWORD LightSphereVertex::FVF = D3DFVF_XYZ | D3DFVF_NORMAL | D3DFVF_DIFFUSE | D3DFVF_SPECULAR;

SphereEx::SphereEx(void):
m_fRMAX(50),
m_iTriangleCount(0),
m_pVB(0)
{
 m_fDirLightPos[0] = 0.00041f;
 m_fDirLightPos[1] = -0.04470f;
 m_fDirLightPos[2] = 0.998998f;
 srand(time(0));
}

SphereEx::~SphereEx(void)
{
}

 

//细腻型
void SphereEx::AddTriangle_UpRight(LightSphereVertex* v, int& in_ID,
          float*** fRet,int in_iRowCenter, int in_iLineCenter,
          int in_iRowNum, int in_iLineNum,bool in_bIsClockwise)
{
 D3DXVECTOR3 pNormal;
 int iY = in_iRowCenter, iX = in_iLineCenter;
 int iRX = iX+1, iLX = iX-1;
 if( iRX >= in_iLineNum) iRX = iRX-in_iLineNum;
 if( iLX < 0) iLX = iLX+in_iLineNum;
 int iUy = iY+1, iDy = iY-1;
 if( iUy >= in_iRowNum) iUy = iUy-in_iRowNum;
 if( iDy < 0) iDy = iDy+in_iRowNum;
 
 bool bIsClockwise = in_bIsClockwise;

 //up_left
 D3DXVECTOR3 p0(fRet[iY][iRX][0],fRet[iY][iRX][1],fRet[iY][iRX][2]);
 D3DXVECTOR3 p1(fRet[iY][iX][0],fRet[iY][iX][1],fRet[iY][iX][2]);
 D3DXVECTOR3 p2(fRet[iUy][iX][0],fRet[iUy][iX][1],fRet[iUy][iX][2]);

 pNormal = Math_D3D::ComputeCenterNormal_Sum((const float***)fRet,iY,iRX,bIsClockwise,
  in_iRowNum,in_iLineNum,1,1);
 v[in_ID++] = LightSphereVertex(p0.x,p0.y,p0.z, pNormal.x, pNormal.y, pNormal.z,
  Math_D3D::GetRGB(m_fRMIN,m_fRMAX,fRet[iY][iRX][3]),
  Math_D3D::GetRGB(m_fRMIN,m_fRMAX,fRet[iY][iRX][3]));

 pNormal = Math_D3D::ComputeCenterNormal_Sum((const float***)fRet,iY,iX,bIsClockwise,
  in_iRowNum,in_iLineNum,1,1);
 v[in_ID++] = LightSphereVertex(p1.x,p1.y,p1.z, pNormal.x, pNormal.y, pNormal.z,
  Math_D3D::GetRGB(m_fRMIN,m_fRMAX,fRet[iY][iX][3]),
  Math_D3D::GetRGB(m_fRMIN,m_fRMAX,fRet[iY][iX][3]));

 pNormal = Math_D3D::ComputeCenterNormal_Sum((const float***)fRet,iUy,iX,bIsClockwise,
  in_iRowNum,in_iLineNum,1,1);
 v[in_ID++] = LightSphereVertex(p2.x,p2.y,p2.z, pNormal.x, pNormal.y, pNormal.z,
  Math_D3D::GetRGB(m_fRMIN,m_fRMAX,fRet[iUy][iX][3]),
  Math_D3D::GetRGB(m_fRMIN,m_fRMAX,fRet[iUy][iX][3]));
}

//细腻型
void SphereEx::AddTriangle_DownRight(LightSphereVertex* v, int& in_ID,
          float*** fRet,int in_iRowCenter, int in_iLineCenter,
          int in_iRowNum, int in_iLineNum,bool in_bIsClockwise)
{
 D3DXVECTOR3 pNormal;
 int iY = in_iRowCenter, iX = in_iLineCenter;
 int iRX = iX+1, iLX = iX-1;
 if( iRX >= in_iLineNum) iRX = iRX-in_iLineNum;
 if( iLX < 0) iLX = iLX+in_iLineNum;
 int iUy = iY+1, iDy = iY-1;
 if( iUy >= in_iRowNum) iUy = iUy-in_iRowNum;
 if( iDy < 0) iDy = iDy+in_iRowNum;
 
 bool bIsClockwise = in_bIsClockwise;

 //up_left
 D3DXVECTOR3 p0(fRet[iY][iRX][0],fRet[iY][iRX][1],fRet[iY][iRX][2]);
 D3DXVECTOR3 p1(fRet[iDy][iX][0],fRet[iDy][iX][1],fRet[iDy][iX][2]);
 D3DXVECTOR3 p2(fRet[iY][iX][0],fRet[iY][iX][1],fRet[iY][iX][2]);

 pNormal = Math_D3D::ComputeCenterNormal_Sum((const float***)fRet,iY,iRX,bIsClockwise,
  in_iRowNum,in_iLineNum,1,1);
 v[in_ID++] = LightSphereVertex(p0.x,p0.y,p0.z, pNormal.x, pNormal.y, pNormal.z,
  Math_D3D::GetRGB(m_fRMIN,m_fRMAX,fRet[iY][iRX][3]),
  Math_D3D::GetRGB(m_fRMIN,m_fRMAX,fRet[iY][iRX][3]));

 pNormal = Math_D3D::ComputeCenterNormal_Sum((const float***)fRet,iDy,iX,bIsClockwise,
  in_iRowNum,in_iLineNum,1,1);
 v[in_ID++] = LightSphereVertex(p1.x,p1.y,p1.z, pNormal.x, pNormal.y, pNormal.z,
  Math_D3D::GetRGB(m_fRMIN,m_fRMAX,fRet[iDy][iX][3]),
  Math_D3D::GetRGB(m_fRMIN,m_fRMAX,fRet[iDy][iX][3]));

 pNormal = Math_D3D::ComputeCenterNormal_Sum((const float***)fRet,iY,iX,bIsClockwise,
  in_iRowNum,in_iLineNum,1,1);
 v[in_ID++] = LightSphereVertex(p2.x,p2.y,p2.z, pNormal.x, pNormal.y, pNormal.z,
  Math_D3D::GetRGB(m_fRMIN,m_fRMAX,fRet[iY][iX][3]),
  Math_D3D::GetRGB(m_fRMIN,m_fRMAX,fRet[iY][iX][3]));
}

//细腻型
void SphereEx::AddTriangle_DownLeft(LightSphereVertex* v, int& in_ID,
          float*** fRet,int in_iRowCenter, int in_iLineCenter,
          int in_iRowNum, int in_iLineNum,bool in_bIsClockwise)
{
 D3DXVECTOR3 pNormal;
 int iY = in_iRowCenter, iX = in_iLineCenter;
 int iRX = iX+1, iLX = iX-1;
 if( iRX >= in_iLineNum) iRX = iRX-in_iLineNum;
 if( iLX < 0) iLX = iLX+in_iLineNum;
 int iUy = iY+1, iDy = iY-1;
 if( iUy >= in_iRowNum) iUy = iUy-in_iRowNum;
 if( iDy < 0) iDy = iDy+in_iRowNum;
 
 bool bIsClockwise = in_bIsClockwise;

 //up_left
 D3DXVECTOR3 p0(fRet[iY][iX][0],fRet[iY][iX][1],fRet[iY][iX][2]);
 D3DXVECTOR3 p1(fRet[iDy][iX][0],fRet[iDy][iX][1],fRet[iDy][iX][2]);
 D3DXVECTOR3 p2(fRet[iY][iLX][0],fRet[iY][iLX][1],fRet[iY][iLX][2]);

 pNormal = Math_D3D::ComputeCenterNormal_Sum((const float***)fRet,iY,iX,bIsClockwise,
  in_iRowNum,in_iLineNum,1,1);
 v[in_ID++] = LightSphereVertex(p0.x,p0.y,p0.z, pNormal.x, pNormal.y, pNormal.z,
  Math_D3D::GetRGB(m_fRMIN,m_fRMAX,fRet[iY][iX][3]),
  Math_D3D::GetRGB(m_fRMIN,m_fRMAX,fRet[iY][iX][3]));

 pNormal = Math_D3D::ComputeCenterNormal_Sum((const float***)fRet,iDy,iX,bIsClockwise,
  in_iRowNum,in_iLineNum,1,1);
 v[in_ID++] = LightSphereVertex(p1.x,p1.y,p1.z, pNormal.x, pNormal.y, pNormal.z,
  Math_D3D::GetRGB(m_fRMIN,m_fRMAX,fRet[iDy][iX][3]),
  Math_D3D::GetRGB(m_fRMIN,m_fRMAX,fRet[iDy][iX][3]));

 pNormal = Math_D3D::ComputeCenterNormal_Sum((const float***)fRet,iY,iLX,bIsClockwise,
  in_iRowNum,in_iLineNum,1,1);
 v[in_ID++] = LightSphereVertex(p2.x,p2.y,p2.z, pNormal.x, pNormal.y, pNormal.z,
  Math_D3D::GetRGB(m_fRMIN,m_fRMAX,fRet[iY][iLX][3]),
  Math_D3D::GetRGB(m_fRMIN,m_fRMAX,fRet[iY][iLX][3]));
}

//细腻型
void SphereEx::AddTriangle_UpLeft(LightSphereVertex* v, int& in_ID,
          float*** fRet,int in_iRowCenter, int in_iLineCenter,
          int in_iRowNum, int in_iLineNum,bool in_bIsClockwise)
{
 D3DXVECTOR3 pNormal;
 int iY = in_iRowCenter, iX = in_iLineCenter;
 int iRX = iX+1, iLX = iX-1;
 if( iRX >= in_iLineNum) iRX = iRX-in_iLineNum;
 if( iLX < 0) iLX = iLX+in_iLineNum;
 int iUy = iY+1, iDy = iY-1;
 if( iUy >= in_iRowNum) iUy = iUy-in_iRowNum;
 if( iDy < 0) iDy = iDy+in_iRowNum;
 
 bool bIsClockwise = in_bIsClockwise;

 //up_left
 D3DXVECTOR3 p0(fRet[iUy][iX][0],fRet[iUy][iX][1],fRet[iUy][iX][2]);
 D3DXVECTOR3 p1(fRet[iY][iX][0],fRet[iY][iX][1],fRet[iY][iX][2]);
 D3DXVECTOR3 p2(fRet[iY][iLX][0],fRet[iY][iLX][1],fRet[iY][iLX][2]);

 pNormal = Math_D3D::ComputeCenterNormal_Sum((const float***)fRet,iUy,iX,bIsClockwise,
  in_iRowNum,in_iLineNum,1,1);
 v[in_ID++] = LightSphereVertex(p0.x,p0.y,p0.z, pNormal.x, pNormal.y, pNormal.z,
  Math_D3D::GetRGB(m_fRMIN,m_fRMAX,fRet[iUy][iX][3]),
  Math_D3D::GetRGB(m_fRMIN,m_fRMAX,fRet[iUy][iX][3]));

 pNormal = Math_D3D::ComputeCenterNormal_Sum((const float***)fRet,iY,iX,bIsClockwise,
  in_iRowNum,in_iLineNum,1,1);
 v[in_ID++] = LightSphereVertex(p1.x,p1.y,p1.z, pNormal.x, pNormal.y, pNormal.z,
  Math_D3D::GetRGB(m_fRMIN,m_fRMAX,fRet[iY][iX][3]),
  Math_D3D::GetRGB(m_fRMIN,m_fRMAX,fRet[iY][iX][3]));

 pNormal = Math_D3D::ComputeCenterNormal_Sum((const float***)fRet,iY,iLX,bIsClockwise,
  in_iRowNum,in_iLineNum,1,1);
 v[in_ID++] = LightSphereVertex(p2.x,p2.y,p2.z, pNormal.x, pNormal.y, pNormal.z,
  Math_D3D::GetRGB(m_fRMIN,m_fRMAX,fRet[iY][iLX][3]),
  Math_D3D::GetRGB(m_fRMIN,m_fRMAX,fRet[iY][iLX][3]));
}

//建立模型
bool SphereEx::SetupSphere(IDirect3DDevice9* pDevice,
      const int in_iNewTS, const int in_iNewTE,
      const int in_iNewTStep,
      const int in_iNewPS, const int in_iNewPE,
      const int in_iNewPStep,
      const int in_iOldTS, const int in_iOldTE,
      const int in_iOldTStep,
      const int in_iOldPS, const int in_iOldPE,
      const int in_iOldPStep,
      const float** in_pfKnowData)
{
 CleanUp();
 const int RET_ROWS = abs((int)(in_iNewPE-in_iNewPS)/in_iNewPStep)+1;
 const int RET_LINES = abs((int)(in_iNewTE-in_iNewTS)/in_iNewTStep)+1;
 const int KNOW_ROWS = abs((int)(in_iOldPE-in_iOldPS)/in_iOldPStep)+1;
 const int KNOW_LINES = abs((int)(in_iOldTE-in_iOldTS)/in_iOldTStep)+1;
 const int IMAGE_RET_ROWS = RET_ROWS+2;
 const int IMAGE_RET_LINES = RET_LINES+2;
 const int OLD_TSTEP =  in_iOldTStep;
 const int OLD_PSTEP =  in_iOldPStep;
 float*** fRet = NULL;
 Math_BiCubic::GetTable(fRet,RET_ROWS,RET_LINES,3);
 float*** fKnow = NULL;
 Math_BiCubic::GetTable(fKnow,KNOW_ROWS,KNOW_LINES,3);
 for(int i=0; i<RET_ROWS; i++)
 {
  for(int j=0; j<RET_LINES; j++)
  {
   fRet[i][j][0] = in_iNewTS + in_iNewTStep*j;
   fRet[i][j][1] = in_iNewPS + in_iNewPStep*i;
   fRet[i][j][2] = 0;
  }
 }
 for(int i=0; i<KNOW_ROWS; i++)
 {
  for(int j=0; j<KNOW_LINES; j++)
  {
   fKnow[i][j][0] = in_iOldTS + OLD_TSTEP*j;
   fKnow[i][j][1] = in_iOldPS + OLD_PSTEP*i;
   fKnow[i][j][2] = in_pfKnowData[i][j];
  }
 }
 Math_BiCubic::CubicSplineEx2_RowLine(fRet,RET_ROWS,RET_LINES,fKnow,KNOW_ROWS,KNOW_LINES);

 //平面图转为极化图
 m_fRMIN = 999999;
 m_fRMAX = -999999;
  float*** fSphereData = NULL;
 Math_BiCubic::GetTable(fSphereData,IMAGE_RET_ROWS,IMAGE_RET_LINES,4);
 for(int i=0; i<RET_ROWS; i++)
 {
  for(int j=0; j<RET_LINES; j++)
  {
   float x,y,z;
   Transform(x,y,z,fRet[i][j][0],fRet[i][j][1],fRet[i][j][2]);
   fSphereData[i+1][j+1][0] = x;
   fSphereData[i+1][j+1][1] = y;
   fSphereData[i+1][j+1][2] = z;
   fSphereData[i+1][j+1][3] = fRet[i][j][2];
   if( m_fRMIN > fRet[i][j][2])
    m_fRMIN = fRet[i][j][2];
   if( m_fRMAX < fRet[i][j][2])
    m_fRMAX = fRet[i][j][2];
  }
 }
 int iTri = (IMAGE_RET_ROWS)*(IMAGE_RET_LINES)*4;
 pDevice->CreateVertexBuffer(
  iTri * 3 * sizeof(LightSphereVertex),
  D3DUSAGE_WRITEONLY,
  LightSphereVertex::FVF,
  D3DPOOL_MANAGED,
  &m_pVB,
  0);
 LightSphereVertex* v;
 m_pVB->Lock(0, 0, (void**)&v, 0);
 int iID = 0;

 //测试数据
 for(int i=2; i<RET_ROWS; i++)
 {
  for(int j=2; j<RET_LINES; j++)
  {
   int iP = in_iNewPS + in_iNewPStep*i;
   bool bIsClockwise =  (iP > 180)? true: false;
   AddTriangle_DownLeft(v,iID,fSphereData,i,j,RET_ROWS, RET_LINES, bIsClockwise);
   AddTriangle_UpRight(v,iID,fSphereData,i-1,j-1,RET_ROWS, RET_LINES,bIsClockwise);
  }
 }
 //上边缘
 for(int i=1; i<RET_ROWS; i++)
 {
  int iP = in_iNewPS + in_iNewPStep*i;
  bool bIsClockwise =  (iP > 180)? true: false;
  AddTriangle_DownLeft(v,iID,fSphereData,i,1,RET_ROWS, RET_LINES, bIsClockwise);
  AddTriangle_UpRight(v,iID,fSphereData,i-1,0,RET_ROWS, RET_LINES,bIsClockwise);
 }
 //下边缘
 for(int i=1; i<RET_ROWS; i++)
 {
  int iP = in_iNewPS + in_iNewPStep*i;
  bool bIsClockwise =  (iP > 180)? true: false;
  AddTriangle_DownLeft(v,iID,fSphereData,i,RET_LINES+1,RET_ROWS, RET_LINES, bIsClockwise);
  AddTriangle_UpRight(v,iID,fSphereData,i-1,RET_LINES-1,RET_ROWS, RET_LINES,bIsClockwise);
 }
 //左边缘
 for(int i=1; i<RET_LINES; i++)
 {
  int iP = in_iNewPS;
  bool bIsClockwise = (iP > 180)? true: false;
  AddTriangle_DownLeft(v,iID,fSphereData,1,i,RET_ROWS, RET_LINES, bIsClockwise);
  AddTriangle_UpRight(v,iID,fSphereData,0,i-1,RET_ROWS, RET_LINES,bIsClockwise);
 }
 //右边缘
 for(int i=1; i<RET_LINES; i++)
 {
  int iP = in_iNewPS + in_iNewPStep*(RET_ROWS-1);
  bool bIsClockwise =  (iP > 180)? true: false;
  AddTriangle_DownLeft(v,iID,fSphereData,RET_ROWS+1,i,RET_ROWS, RET_LINES, bIsClockwise);
  AddTriangle_UpRight(v,iID,fSphereData,RET_ROWS-1,i-1,RET_ROWS, RET_LINES,bIsClockwise);
 }

 m_iTriangleCount = iID/3;
 m_pVB->Unlock();
 v = NULL;

 Math_BiCubic::DeleteTable(fKnow,KNOW_ROWS,KNOW_LINES);
 Math_BiCubic::DeleteTable(fRet,RET_ROWS,RET_LINES);
 Math_BiCubic::DeleteTable(fSphereData,IMAGE_RET_ROWS,IMAGE_RET_LINES);

 SetupDirectLight(pDevice,LIGHT_DIRECT_ID,m_DirLight);
 SetLightDir(m_DirLight,m_fDirLightPos[0],m_fDirLightPos[1], m_fDirLightPos[2]);
 return true;
}

//清空模型
void SphereEx::CleanUp()
{
 D3D::Release<IDirect3DVertexBuffer9*>(m_pVB);
}

//设置状态
void SphereEx::SetState(IDirect3DDevice9* pDevice)
{
 pDevice->LightEnable(LIGHT_DIRECT_ID, true);
 pDevice->SetRenderState(D3DRS_LIGHTING, true);
 pDevice->SetRenderState(D3DRS_SHADEMODE,D3DSHADE_GOURAUD);
 pDevice->SetRenderState(D3DRS_NORMALIZENORMALS, true);
 pDevice->SetRenderState(D3DRS_SPECULARENABLE, true);
 pDevice->SetRenderState(D3DRS_CULLMODE,D3DCULL_NONE);
}

//渲染
void SphereEx::Rander(IDirect3DDevice9* pDevice)
{
 SetState(pDevice);
 Transform_WORLD(pDevice);
 pDevice->SetStreamSource(0, m_pVB, 0, sizeof(LightSphereVertex));
 pDevice->SetFVF(LightSphereVertex::FVF);
 D3DMATERIAL9 mtrl;
 mtrl.Ambient = D3D::BLUE;
 mtrl.Diffuse  = D3D::WHITE;
 mtrl.Specular = D3D::WHITE;
 mtrl.Emissive = D3D::BLACK;
 mtrl.Power    = 120.0f;
 pDevice->SetMaterial(&mtrl);
 pDevice->DrawPrimitive(D3DPT_TRIANGLELIST, 0, m_iTriangleCount);
}

void SphereEx::Transform(float& x, float& y, float& z,
      const int in_iTA, const int in_iPA, const float in_fValue)
{
 float fValue = in_fValue;
 if( fValue < 0)
  fValue *= -1;
 Math_D3D::CumputerPosByTurnXY(x,y,z,0,fValue,0,in_iPA,in_iTA);
}

 

////////////////////////////////////////////////////

 

/****************************************
*文件名:  3D_Math.h
*文件描述:  建立3D开发里常用到的函数库
*创建人:  陈泽丹
*创建时间:  2010-04-03
*版本:   1.0
*****************************************/

#include <d3dx9.h>

 

namespace Math_D3D
{
 /*
 求in_fX,in_fY,in_fZ,绕Y轴转in_iA_XTAngle度[0,360],
 绕X轴转in_iA_XTAngle度[0,360]的位置out_fX, out_fY, out_fZ
 TurnXangle = A_Theta; TurnYangle = B_Phi;
 cosB  0   -sinB  0
 sinAsinB consA  sinAcosB 0
 cosAsinB -sinA  cosAcosB 0
 0   0   0   1
 */
 void CumputerPosByTurnXY(float& out_fX, float& out_fY, float& out_fZ,
        const float in_fX, const float in_fY, const float in_fZ,
        const int& in_iA_XTAngle, const int& in_iB_YTAngle);

 /*
 求p0,p1,p2三点组成的面法线
 */
 D3DXVECTOR3 ComputeNormal(const D3DXVECTOR3* p0, const D3DXVECTOR3* p1, const D3DXVECTOR3* p2);

 /*
 画矩型
 */
 void CreateRectangle(D3DXVECTOR3 vTrianglePoint[6],
       const float UpLeftX, const float UpLeftY, const float UpLeftZ,
       const float UpRightX, const float UpRightY, const float UpRightZ,
       const float DownLeftX, const float DownLeftY, const float DownLeftZ,
       const float DownRightX, const float DownRightY, const float DownRightZ);

 /*
 计算平面里以in_iCx,in_iCy为中心点,被四个三角形环绕共享时
 顺时针方向每个三角形三个点的ID
 */
 void ComputeCenterTriangle_UpRight(int out_iIDX[3], int out_iIDY[3],
            const int CenterX, const int CenterY,
            const int XPOINT_COUNT, const int YPOINT_COUNT,
            const int W, const int H);
 void ComputeCenterTriangle_DownRight(int out_iIDX[3], int out_iIDY[3],
            const int CenterX, const int CenterY,
            const int XPOINT_COUNT, const int YPOINT_COUNT,
            const int W, const int H);
 void ComputeCenterTriangle_DownLeft(int out_iIDX[3], int out_iIDY[3],
            const int CenterX, const int CenterY,
            const int XPOINT_COUNT, const int YPOINT_COUNT,
            const int W, const int H);
 void ComputeCenterTriangle_UpLef(int out_iIDX[3], int out_iIDY[3],
            const int CenterX, const int CenterY,
            const int XPOINT_COUNT, const int YPOINT_COUNT,
            const int W, const int H);

 /*
 计算以in_iCx,in_iCy为中心点,被四个三角形环绕共享时,
 中心点in_iCx,in_iCy的法线
 */
 D3DXVECTOR3 ComputeCenterNormal_Sum(const float*** in_fPoint,
         const int in_iCx, const int in_iCy,
         const bool in_bIsClockwise,
         const int XPOINT_COUNT, const int YPOINT_COUNT,
         const int W, const int H);

 //上色方案
 float GetRGB(const float in_fMin, const float in_fMax, const float in_fValue);
}

 

////////////////////////////////////////

 

/****************************************
*文件名:  3D_Math.cpp
*创建人:  陈泽丹
*创建时间:  2010-04-03
*版本:   1.0
*****************************************/

#include "stdafx.h"
#include "Math_D3D.h"
#include <math.h>

 

/*
计算以in_iCx,in_iCy为中心点,被四个三角形环绕共享时
顺时针方向每个三角形三个点的ID
*/
const int Math_D3D_UpRight = 0;
const int Math_D3D_DownRight = 1;
const int Math_D3D_DownLeft = 2;
const int Math_D3D_UpLeft = 3;
void ComputeCenterTriangle(int out_iIDX[3], int out_iIDY[3],
       const int TYPE, const int CenterX, const int CenterY,
       const int XPOINT_COUNT, const int YPOINT_COUNT,
       const int W, const int H)
{
 D3DXVECTOR3 vNormal;
 D3DXVECTOR3 vNormalUpRight,vNormaldownRight,vNormalDownLeft,vNormalUpLeft;

 int Cx = CenterX, Cy = CenterY;
 int Rx = Cx+W;
 while( Rx >= XPOINT_COUNT)
  Rx = Rx-XPOINT_COUNT;
 int Lx = Cx-W;
 while( Lx < 0)
  Lx = Lx+XPOINT_COUNT;
 int Uy = Cy+H;
 while( Uy >= YPOINT_COUNT)
  Uy = Uy-YPOINT_COUNT;
 int Dy = Cy-H;
 while( Dy < 0)
  Dy = Dy+YPOINT_COUNT;

 if( Math_D3D_UpRight == TYPE)
 {//up_right
  out_iIDX[0] = Cx;
  out_iIDY[0] = Cy;
  out_iIDX[1] = Cx;
  out_iIDY[1] = Uy;
  out_iIDX[2] = Rx;
  out_iIDY[2] = Cy;
 }
 else if( Math_D3D_DownRight == TYPE)
 {//down_right
  out_iIDX[0] = Cx;
  out_iIDY[0] = Cy;
  out_iIDX[1] = Rx;
  out_iIDY[1] = Cy;
  out_iIDX[2] = Cx;
  out_iIDY[2] = Dy;
 }
 else if( Math_D3D_DownLeft == TYPE)
 {//down_left
  out_iIDX[0] = Cx;
  out_iIDY[0] = Cy;
  out_iIDX[1] = Cx;
  out_iIDY[1] = Dy;
  out_iIDX[2] = Lx;
  out_iIDY[2] = Cy;
 }
 else if( Math_D3D_UpLeft == TYPE)
 {//up_left
  out_iIDX[0] = Cx;
  out_iIDY[0] = Cy;
  out_iIDX[1] = Lx;
  out_iIDY[1] = Cy;
  out_iIDX[2] = Cx;
  out_iIDY[2] = Uy;
 }
}

/*
求in_fX,in_fY,in_fZ,绕Y轴转in_iA_XTAngle度[0,360],
绕X轴转in_iA_XTAngle度[0,360]的位置out_fX, out_fY, out_fZ
TurnXangle = A_Theta; TurnYangle = B_Phi;
cosB  0   -sinB  0
sinAsinB consA  sinAcosB 0
cosAsinB -sinA  cosAcosB 0
0   0   0   1
*/
void Math_D3D::CumputerPosByTurnXY(float& out_fX, float& out_fY, float& out_fZ,
       const float in_fX, const float in_fY, const float in_fZ,
       const int& in_iA_XTAngle, const int& in_iB_YTAngle)
 {
 const float fA = (float)in_iA_XTAngle/180*3.14159;
 const float fB = (float)in_iB_YTAngle/180*3.14159;
 out_fX = in_fX*cos((float)fB)+in_fY*sin((float)fA)*sin((float)fB)+in_fZ*cos((float)fA)*sin((float)fB);
 out_fY = in_fY*cos((float)fA)-in_fZ*sin((float)fA);
 out_fZ = in_fX*(-1*sin((float)fB))+in_fY*sin((float)fA)*cos((float)fB)+in_fZ*cos((float)fA)*cos((float)fB);
}

/*
求p0,p1,p2三点组成的面法线
*/
D3DXVECTOR3 Math_D3D::ComputeNormal(const D3DXVECTOR3* p0, const D3DXVECTOR3* p1, const D3DXVECTOR3* p2)
{
 D3DXVECTOR3 v0ut;
 D3DXVECTOR3 u = *p1-*p0;
 D3DXVECTOR3 v = *p2-*p0;
 D3DXVec3Cross(&v0ut,&u,&v);
 D3DXVec3Normalize(&v0ut,&v0ut);
 return v0ut;
}

/*
画矩型
*/
void Math_D3D::CreateRectangle(D3DXVECTOR3 vTrianglePoint[6],
      const float UpLeftX, const float UpLeftY, const float UpLeftZ,
      const float UpRightX, const float UpRightY, const float UpRightZ,
      const float DownLeftX, const float DownLeftY, const float DownLeftZ,
      const float DownRightX, const float DownRightY, const float DownRightZ)
{
 vTrianglePoint[0].x = UpLeftX;
 vTrianglePoint[0].y = UpLeftY;
 vTrianglePoint[0].z = UpLeftZ;

 vTrianglePoint[1].x = DownRightX;
 vTrianglePoint[1].y = DownRightY;
 vTrianglePoint[1].z = DownRightZ;

 vTrianglePoint[2].x = DownLeftX;
 vTrianglePoint[2].y = DownLeftY;
 vTrianglePoint[2].z = DownLeftZ;

 vTrianglePoint[3].x = UpLeftX;
 vTrianglePoint[3].y = UpLeftY;
 vTrianglePoint[3].z = UpLeftZ;

 vTrianglePoint[4].x = UpRightX;
 vTrianglePoint[4].y = UpRightY;
 vTrianglePoint[4].z = UpRightZ;

 vTrianglePoint[5].x = DownRightX;
 vTrianglePoint[5].y = DownRightY;
 vTrianglePoint[5].z = DownRightZ;
}

void Math_D3D::ComputeCenterTriangle_UpRight(int out_iIDX[3], int out_iIDY[3],
           const int CenterX, const int CenterY,
           const int XPOINT_COUNT, const int YPOINT_COUNT,
           const int W, const int H)
{
 ComputeCenterTriangle(out_iIDX,out_iIDY,Math_D3D_UpRight,CenterX,CenterY,XPOINT_COUNT,YPOINT_COUNT,W,H);
}

void Math_D3D::ComputeCenterTriangle_DownRight(int out_iIDX[3], int out_iIDY[3],
           const int CenterX, const int CenterY,
           const int XPOINT_COUNT, const int YPOINT_COUNT,
           const int W, const int H)
{
 ComputeCenterTriangle(out_iIDX,out_iIDY,Math_D3D_DownRight,CenterX,CenterY,XPOINT_COUNT,YPOINT_COUNT,W,H);
}

void Math_D3D::ComputeCenterTriangle_DownLeft(int out_iIDX[3], int out_iIDY[3],
           const int CenterX, const int CenterY,
           const int XPOINT_COUNT, const int YPOINT_COUNT,
           const int W, const int H)
{
 ComputeCenterTriangle(out_iIDX,out_iIDY,Math_D3D_DownLeft,CenterX,CenterY,XPOINT_COUNT,YPOINT_COUNT,W,H);
}

void Math_D3D::ComputeCenterTriangle_UpLef(int out_iIDX[3], int out_iIDY[3],
           const int CenterX, const int CenterY,
           const int XPOINT_COUNT, const int YPOINT_COUNT,
           const int W, const int H)
{
 ComputeCenterTriangle(out_iIDX,out_iIDY,Math_D3D_UpLeft,CenterX,CenterY,XPOINT_COUNT,YPOINT_COUNT,W,H);
}

/*
计算以in_iCx,in_iCy为中心点,被四个三角形环绕共享时,
每个三角形中in_iCx,in_iCy的法线
*/
D3DXVECTOR3 ComputeCenterNormal(const float***  in_fPoint,
        const int TYPE,
        const bool in_bIsClockwise,
        const int CenterX, const int CenterY,
        const int XPOINT_COUNT, const int YPOINT_COUNT,
        const int W, const int H)
{
 D3DXVECTOR3 vNormal;
 D3DXVECTOR3 vNormalUpRight,vNormaldownRight,vNormalDownLeft,vNormalUpLeft;
 int IDX[3];
 int IDY[3];

 ComputeCenterTriangle(IDX,IDY, TYPE, CenterX,CenterY,XPOINT_COUNT,YPOINT_COUNT,W,H);
 D3DXVECTOR3 p0(in_fPoint[IDX[0]][IDY[0]][0],in_fPoint[IDX[0]][IDY[0]][1],in_fPoint[IDX[0]][IDY[0]][2]);
 D3DXVECTOR3 p1(in_fPoint[IDX[1]][IDY[1]][0],in_fPoint[IDX[1]][IDY[1]][1],in_fPoint[IDX[1]][IDY[1]][2]);
 D3DXVECTOR3 p2(in_fPoint[IDX[2]][IDY[2]][0],in_fPoint[IDX[2]][IDY[2]][1],in_fPoint[IDX[2]][IDY[2]][2]);
 if( in_bIsClockwise)
  vNormal = Math_D3D::ComputeNormal(&p0,&p1,&p2);
 else
  vNormal = Math_D3D::ComputeNormal(&p2,&p1,&p0);
 return vNormal;
}

/*
计算以in_iCx,in_iCy为中心点,被四个三角形环绕共享时,
中心点in_iCx,in_iCy的法线
*/
D3DXVECTOR3 Math_D3D::ComputeCenterNormal_Sum(const float*** in_fPoint,
         const int in_iCx, const int in_iCy,
         const bool in_bIsClockwise,
         const int XPOINT_COUNT, const int YPOINT_COUNT,
         const int W, const int H)
{
 D3DXVECTOR3 pNormal;
 D3DXVECTOR3 pNormalUpRight,pNormaldownRight,pNormalDownLeft,pNormalUpLeft;
 bool bIsClockwise = in_bIsClockwise;
 //up_right
 pNormalUpRight = ComputeCenterNormal(in_fPoint,
  Math_D3D_UpRight,
  bIsClockwise,
  in_iCx,in_iCy,
  XPOINT_COUNT,YPOINT_COUNT,
  W,H);

 //down_right
 pNormaldownRight = ComputeCenterNormal(in_fPoint,
  Math_D3D_DownRight,
  bIsClockwise,
  in_iCx,in_iCy,
  XPOINT_COUNT,YPOINT_COUNT,
  W,H);

 
 //down_left
 pNormalDownLeft = ComputeCenterNormal(in_fPoint,
  Math_D3D_DownLeft,
  bIsClockwise,
  in_iCx,in_iCy,
  XPOINT_COUNT,YPOINT_COUNT,
  W,H);

 //up_left
 pNormalUpLeft = ComputeCenterNormal(in_fPoint,
  Math_D3D_UpLeft,
  bIsClockwise,
  in_iCx,in_iCy,
  XPOINT_COUNT,YPOINT_COUNT,
  W,H);

 pNormal = (pNormalUpRight+pNormaldownRight+pNormalDownLeft+pNormalUpLeft)*0.25;
 return pNormal;
}

 //上色方案
float Math_D3D::GetRGB(const float in_fMin, const float in_fMax,
        const float in_fValue)
{
 float fRate = (float)(in_fValue-in_fMin)/(in_fMax-in_fMin);
 fRate = (fRate<1)? fRate: 1;
 fRate = (fRate>0)? fRate: 0;

 int R,G,B;
 if( fRate >= 0.75)
 {
  fRate -= 0.75;
  R = 255;
  G = 255-fRate/0.25*255;
  B = 0;
 }
 else if( fRate >= 0.5)
 {
  fRate -= 0.5;
  R = fRate/0.25*255;
  G = 255;
  B = 0;
 }
 else if( fRate >= 0.25)
 {
  fRate -= 0.25;
  R = 0;
  G = 255;
  B = 255-fRate/0.25*255;
 }
 else if( fRate >= 0)
 {
  R = 0;
  G = fRate/0.25*255;
  B = 255;
 }
 else
 {
  R = 0;
  G = 0;
  B = 0;
 }

 return D3DCOLOR_ARGB(255,R,G,B);
}

 

//////////////////////////////////////////////

 

#pragma once

namespace Math_BiCubic
{
 template<typename T>
 bool GetTable(T***& out_fTable, int in_iTH, int in_iTW, int in_iTK)
 {
  out_fTable = 0;
  out_fTable = new T**[in_iTH];
  if( 0 == out_fTable)
   return false;
  for(int i=0; i<in_iTH; i++)
  {
   out_fTable[i] = new T*[in_iTW];
   if( 0 == out_fTable[i])
    return false;
   for(int j=0; j<in_iTW; j++)
   {
    out_fTable[i][j] = new T[in_iTK];
    if( 0 == out_fTable[i][j])
     return false;
   }
  }
  for(int i=0; i<in_iTH; i++)
  {
   for(int j=0; j<in_iTW; j++)
   {
    for(int k=0; k<in_iTK; k++)
     out_fTable[i][j][k] = 0;
   }
  }
 }

 template<typename T>
 void DeleteTable(T*** in_dTable, int in_iTH, int in_iTW)
 {
  for(int i=0; i<in_iTH; i++)
  {
   for(int j=0; j<in_iTW; j++)
   {
    delete[] in_dTable[i][j];
    in_dTable[i][j] = 0;
   }
   delete[] in_dTable[i];
   in_dTable[i] = 0;
  }
  delete[] in_dTable;
  in_dTable = 0;
 }

 template<typename T>
 void Sort(T* in_tData, int in_iSize)
 {
  for(int i=0; i<in_iSize; i++)
  {
   int iRecored = i;
   for(int j=i+1; j<in_iSize; j++)
   {
    if( in_tData[j] < in_tData[iRecored])
     iRecored = j;
   }
   if( i != iRecored)
   {
    T tTemp = in_tData[i];
    in_tData[i] = in_tData[iRecored];
    in_tData[iRecored] = tTemp;
   }
  }
 }

 

 /* 三次样条曲线
 out_fVal   所求函数值
 out_fD1  所求插值点的一阶导数
 out_fD2 所求插值点的二阶导数
 in_iM 已知节点个数
 in_fPos    所求插值点的插值点的X轴坐标(无次序要求且可含已知点)
 in_fKnowPos,in_fKnowVal  已知点的节点坐标(要求升序排列)
 in_fKnowD1,  in_fKnowD2 已知点的节点的一阶导和二阶导(要求升序排列)
 in_iN 已知节点个数
 */
 void CubicSpline(float out_fVal[], float out_fD1[], float out_fD2[],
      float in_fPos[], int in_iM,
      float in_fKnowPos[], float in_fKnowVal[],
      float in_fKnowD1[], float in_fKnowD2[],
      int in_iN);

 //以估计的方法计算导数
 void CalculateDyDdy(float out_fD1[], float out_fD2[], 
      float in_fPos[], float in_fVal[], int in_iN,
      float in_fStartD1 = 0, float in_fEndD1 = 0,
      float in_fStartD2 = 0, float in_fEndD2 = 0);

 /* 三次样条曲线
 out_fVal   所求函数值
 in_fPos    所求插值点的插值点的X轴坐标(无次序要求且可含已知点)
 in_iM 已知节点个数
 in_fKnowPos  已知点的节点坐标(要求升序排列)
 in_fKnowVal  已知点的节点坐标值(要求升序排列)
 in_iN 已知节点个数
 */
 //三次样条曲线(1维(输出结果可含已知点))
 void CubicSplineEx(float* out_fVal, float* in_fPos, int in_iM,
  float* in_fKnowPos, float* in_fKnowVal, int in_iN);

 //三次样条曲线(2维)
 bool CubicSplineEx2_ROW(float*** inout_fVal,int in_iValRows,int in_iValLines,
  float*** in_fKnowVal,int in_iKnowRows,int in_iKnowLines);

 //三次样条曲线(2维)
 bool CubicSplineEx2_LINE(float*** inout_fVal,int in_iValRows,int in_iValLines,
  float*** in_fKnowVal,int in_iKnowRows,int in_iKnowLines);

 //三次样条曲线(2维)
 bool CubicSplineEx2_RowLine(float*** inout_fVal,int in_iValRows,int in_iValLines,
  float*** in_fKnowVal,int in_iKnowRows,int in_iKnowLines);

 //三次样条曲线(2维)
 bool CubicSplineEx2_LineRow(float*** inout_fVal,int in_iValRows,int in_iValLines,
  float*** in_fKnowVal,int in_iKnowRows,int in_iKnowLines);
}

 

#include "stdafx.h"
#include "Math_BiCubic.h"

 

//以估计的方法计算导数
void Math_BiCubic::CalculateDyDdy(float out_fD1[], float out_fD2[], 
     float in_fPos[], float in_fVal[], int in_iN,
     float in_fStartD1, float in_fEndD1,
     float in_fStartD2, float in_fEndD2)
{//以估计的方法计算导数
 out_fD1[0] = in_fStartD1;
 out_fD2[0] = in_fStartD2;
 out_fD1[in_iN-1] = in_fEndD1;
 out_fD2[in_iN-1] = in_fEndD2;
 for(int i=1; i<in_iN-1; i++)
 {
  float dx = in_fPos[i+1] - in_fPos[i];
  if( 0 != dx)
  {
   out_fD1[i] = (in_fVal[i+1] - in_fVal[i]) / (2*dx);
   out_fD2[i] = (in_fVal[i+1] + in_fVal[i-1] - 2*in_fVal[i]) / dx / dx;
  }
 }
}

/* 三次样条曲线
out_fVal   所求函数值
out_fD1  所求插值点的一阶导数
out_fD2 所求插值点的二阶导数
in_iM 已知节点个数
in_fPos    所求插值点的插值点的X轴坐标(无次序要求且可含已知点)
in_fKnowPos,in_fKnowVal  已知点的节点坐标(要求升序排列)
in_fKnowD1,  in_fKnowD2 已知点的节点的一阶导和二阶导(要求升序排列)
in_iN 已知节点个数
*/
void Math_BiCubic::CubicSpline(float out_fVal[], float out_fD1[], float out_fD2[],
     float in_fPos[], int in_iM,
     float in_fKnowPos[], float in_fKnowVal[],
     float in_fKnowD1[], float in_fKnowD2[],
     int in_iN)  
{
 //三次样条函数
 int   i,j;  
 float   h0,h1,alpha,beta,g;
 int SCOUN = (in_iM<400)? 400:(in_iM/3);
 float* s = new float[in_iM];

 in_fKnowD1[0]=-0.5;  
 h0=in_fKnowPos[1]-in_fKnowPos[0];  
 s[0]=3.0*(in_fKnowVal[1]-in_fKnowVal[0])/(2.0*h0)-in_fKnowD2[0]*h0/4.0;  
 for(j=1;j<=in_iN-2;j++)  
 {
   h1=in_fKnowPos[j+1]-in_fKnowPos[j];  
   alpha=h0/(h0+h1);  
   beta=(1.0-alpha)*(in_fKnowVal[j]-in_fKnowVal[j-1])/h0;  
   beta=3.0*(beta+alpha*(in_fKnowVal[j+1]-in_fKnowVal[j])/h1);
   in_fKnowD1[j]=-alpha/(2.0+(1.0-alpha)*in_fKnowD1[j-1]);  
   s[j]=(beta-(1.0-alpha)*s[j-1]);  
   s[j]=s[j]/(2.0+(1.0-alpha)*in_fKnowD1[j-1]);  
   h0=h1;  
 }

 in_fKnowD1[in_iN-1]=(3.0*(in_fKnowVal[in_iN-1]-in_fKnowVal[in_iN-2])/h1+in_fKnowD2[in_iN-1]*h1/2.0-s[in_iN-2])/(2.0+in_fKnowD1[in_iN-2]);  
 for(j=in_iN-2;j>=0;j--)
   in_fKnowD1[j]=in_fKnowD1[j]*in_fKnowD1[j+1]+s[j];  
 for(j=0;j<=in_iN-2;j++)  
   s[j]=in_fKnowPos[j+1]-in_fKnowPos[j];  
 for(j=0;j<=in_iN-2;j++)
 {
   h1=s[j]*s[j];
   in_fKnowD2[j]=6.0*(in_fKnowVal[j+1]-in_fKnowVal[j])/h1-2.0*(2.0*in_fKnowD1[j]+in_fKnowD1[j+1])/s[j];
 }  
 h1=s[in_iN-2]*s[in_iN-2];  
 in_fKnowD2[in_iN-1]=6.*(in_fKnowVal[in_iN-2]-in_fKnowVal[in_iN-1])/h1+2.*(2.*in_fKnowD1[in_iN-1]+in_fKnowD1[in_iN-2])/s[in_iN-2];  
 g=0.0;  
 for(i=0;i<=in_iN-2;i++)
 {
   h1=0.5*s[i]*(in_fKnowVal[i]+in_fKnowVal[i+1]);
   h1=h1-s[i]*s[i]*s[i]*(in_fKnowD2[i]+in_fKnowD2[i+1])/24.0;
   g=g+h1;
 }  

 for(j=0;j<=in_iM-1;j++)
 {
   if(in_fPos[j]>=in_fKnowPos[in_iN-1])  
    i=in_iN-2;
   else
   {
    i=0;
    while(in_fPos[j]>in_fKnowPos[i+1])  
     i=i+1;
   }  
   h1=(in_fKnowPos[i+1]-in_fPos[j])/s[i];  
   h0=h1*h1;  
   out_fVal[j]=(3.0*h0-2.0*h0*h1)*in_fKnowVal[i];  
   out_fVal[j]=out_fVal[j]+s[i]*(h0-h0*h1)*in_fKnowD1[i];  
   out_fD1[j]=6.0*(h0-h1)*in_fKnowVal[i]/s[i];  
   out_fD1[j]=out_fD1[j]+(3.0*h0-2.0*h1)*in_fKnowD1[i];  
   out_fD2[j]=(6.0-12.0*h1)*in_fKnowVal[i]/(s[i]*s[i]);  
   out_fD2[j]=out_fD2[j]+(2.0-6.0*h1)*in_fKnowD1[i]/s[i];  
   h1=(in_fPos[j]-in_fKnowPos[i])/s[i];  
   h0=h1*h1;  
   out_fVal[j]=out_fVal[j]+(3.0*h0-2.0*h0*h1)*in_fKnowVal[i+1];  
   out_fVal[j]=out_fVal[j]-s[i]*(h0-h0*h1)*in_fKnowD1[i+1];  
   out_fD1[j]=out_fD1[j]-6.0*(h0-h1)*in_fKnowVal[i+1]/s[i];  
   out_fD1[j]=out_fD1[j]+(3.0*h0-2.0*h1)*in_fKnowD1[i+1];  
   out_fD2[j]=out_fD2[j]+(6.0-12.0*h1)*in_fKnowVal[i+1]/(s[i]*s[i]);  
   out_fD2[j]=out_fD2[j]-(2.0-6.0*h1)*in_fKnowD1[i+1]/s[i];
 }

 if( NULL != s)
  delete s;
 s = NULL;
}

//三次样条曲线(1维(输出结果可含已知点))
void Math_BiCubic::CubicSplineEx(float* out_fVal, float* in_fPos, int in_iM, float* in_fKnowPos, float* in_fKnowVal, int in_iN)
{
 float* fKnowDy = new float[in_iN];
 float* fKnowDdy = new float[in_iN];
 float* fInsertDy = new float[in_iM];
 float* fInsertDdy = new float[in_iM];

 Math_BiCubic::CalculateDyDdy(fKnowDy, fKnowDdy, in_fKnowPos, in_fKnowVal, in_iN);
 Math_BiCubic::CubicSpline(out_fVal, fInsertDy, fInsertDdy, in_fPos, in_iM,
  in_fKnowPos, in_fKnowVal, fKnowDy, fKnowDdy, in_iN);

 for(int i=0; i<in_iM; i++)
 {
  int iRecored = i;
  for(int j=i+1; j<in_iM; j++)
  {
   if( in_fPos[j] < in_fPos[iRecored])
    iRecored = j;
  }
  if( i != iRecored)
  {
   float fTemp = in_fPos[i];
   in_fPos[i] = in_fPos[iRecored];
   in_fPos[iRecored] = fTemp;
   
   fTemp = out_fVal[i];
   out_fVal[i] = out_fVal[iRecored];
   out_fVal[iRecored] = fTemp;
  }
 }

 if( NULL != fInsertDdy)
  delete fInsertDdy;
 fInsertDdy = NULL;
 if( NULL != fInsertDy)
  delete fInsertDy;
 fInsertDy = NULL;
 if( NULL != fKnowDdy)
  delete fKnowDdy;
 fKnowDdy = NULL;
 if( NULL != fKnowDy)
  delete fKnowDy;
 fKnowDy = NULL;
}

//三次样条曲线(2维)
bool Math_BiCubic::CubicSplineEx2_ROW(float*** inout_fVal,int in_iValRows,int in_iValLines,
          float*** in_fKnowVal,int in_iKnowRows,int in_iKnowLines)
{
 if( in_iValRows < in_iKnowRows || in_iKnowRows <= 0)
  return false;

 // Rows
 const int ROWK = (double)(in_iValRows-1)/(in_iKnowRows-1)+0.5;
 float* fPos = new float[in_iValLines];
 float* fVal = new float[in_iValLines];
 float* fKnowPos = new float[in_iKnowLines];
 float* fKnowVal = new float[in_iKnowLines];
 for(int i=0; i<in_iValRows; i+=ROWK)
 {
  for(int j=0; j<in_iValLines; j++)
  {
   fPos[j] = inout_fVal[i][j][0];
  }
  for(int j=0; j<in_iKnowLines; j++)
  {
   fKnowPos[j] = in_fKnowVal[i/ROWK][j][0];
   fKnowVal[j] = in_fKnowVal[i/ROWK][j][2];
  }

  CubicSplineEx(fVal,fPos,in_iValLines,fKnowPos,fKnowVal,in_iKnowLines);

  for(int j=0; j<in_iValLines; j++)
  {
   inout_fVal[i][j][2] = fVal[j];
  }
 }
 if( NULL != fPos)
  delete fPos;
 fPos = NULL;
 if( NULL != fVal)
  delete fVal;
 fVal = NULL;
 if( NULL != fKnowPos)
  delete fKnowPos;
 fKnowPos = NULL;
 if( NULL != fKnowVal)
  delete fKnowVal;
 return true;
}

//三次样条曲线(2维)
bool Math_BiCubic::CubicSplineEx2_LINE(float*** inout_fVal,int in_iValRows,int in_iValLines,
          float*** in_fKnowVal,int in_iKnowRows,int in_iKnowLines)
{
 if( in_iValLines < in_iKnowLines || in_iKnowLines<=1 )
  return false;

 // Line
 const int LINIK = (in_iValLines-1)/(in_iKnowLines-1);
 float* fPos = new float[in_iValRows];
 float* fVal = new float[in_iValRows];
 float* fKnowPos = new float[in_iKnowRows];
 float* fKnowVal = new float[in_iKnowRows];
 for(int i=0; i<in_iValLines; i+=LINIK)
 {
  for(int j=0; j<in_iValRows; j++)
  {
   fPos[j] = inout_fVal[j][i][1];
  }
  for(int j=0; j<in_iKnowRows; j++)
  {
   fKnowPos[j] = in_fKnowVal[j][i/LINIK][1];
   fKnowVal[j] = in_fKnowVal[j][i/LINIK][2];
  }

  CubicSplineEx(fVal,fPos,in_iValRows,fKnowPos,fKnowVal,in_iKnowRows);

  for(int j=0; j<in_iValRows; j++)
  {
   inout_fVal[j][i][2] = fVal[j];
  }
 }
 if( NULL != fPos)
  delete fPos;
 fPos = NULL;
 if( NULL != fVal)
  delete fVal;
 fVal = NULL;
 if( NULL != fKnowPos)
  delete fKnowPos;
 fKnowPos = NULL;
 if( NULL != fKnowVal)
  delete fKnowVal;

 return true;
}

//三次样条曲线(2维)
bool Math_BiCubic::CubicSplineEx2_RowLine(float*** inout_fVal,int in_iValRows,int in_iValLines,
            float*** in_fKnowVal,int in_iKnowRows,int in_iKnowLines)
{
 CubicSplineEx2_ROW(inout_fVal,in_iValRows,in_iValLines,in_fKnowVal,in_iKnowRows,in_iKnowLines);
 const int ROWK = (in_iValRows-1)/(in_iKnowRows-1);
 float*** fKnow = NULL;
 GetTable(fKnow,in_iKnowRows,in_iValLines,3);
 for(int i=0; i<in_iKnowRows; i++)
 {
  for(int j=0; j<in_iValLines; j++)
  {
   fKnow[i][j][0] = inout_fVal[i*ROWK][j][0];
   fKnow[i][j][1] = inout_fVal[i*ROWK][j][1];
   fKnow[i][j][2] = inout_fVal[i*ROWK][j][2];
  }
 }
 CubicSplineEx2_LINE(inout_fVal,in_iValRows,in_iValLines,fKnow,in_iKnowRows,in_iValLines);
 DeleteTable(fKnow,in_iKnowRows,in_iValLines);
 return true;
}

//三次样条曲线(2维)
bool Math_BiCubic::CubicSplineEx2_LineRow(float*** inout_fVal,int in_iValRows,int in_iValLines,
            float*** in_fKnowVal,int in_iKnowRows,int in_iKnowLines)
{
 CubicSplineEx2_LINE(inout_fVal,in_iValRows,in_iValLines,in_fKnowVal,in_iKnowRows,in_iKnowLines);
 const int LINIK = (in_iValLines-1)/(in_iKnowLines-1);
 float*** fKnow = NULL;
 GetTable(fKnow,in_iValRows,in_iKnowLines,3);
 for(int i=0; i<in_iValRows; i++)
 {
  for(int j=0; j<in_iKnowLines; j++)
  {
   fKnow[i][j][0] = inout_fVal[i][j*LINIK][0];
   fKnow[i][j][1] = inout_fVal[i][j*LINIK][1];
   fKnow[i][j][2] = inout_fVal[i][j*LINIK][2];
  }
 }
 CubicSplineEx2_ROW(inout_fVal,in_iValRows,in_iValLines,fKnow,in_iValRows,in_iKnowLines);
 DeleteTable(fKnow,in_iValRows,in_iKnowLines);
 return true;
}