功能：根據參數u值和k（大小為階數值）與節點矢量，計算第i個k次B樣條基數
輸入參數： u—參數值；k—大小值為階數；i—第i個k次B樣條的支撐區間左端節點的下標；aNode為節點向量。
輸出參數：返回函數值。
double GetBaseFunVal(double u, int i, int k, vector <double> m_aNode)
{
double Val = 0.0;
double Val1 = 0.0;
double Val2 = 0.0;
if (k==0)
{
if (u < m_aNode[i] || u > m_aNode[i+1])
return Val;
else
{
Val = 1.0;
return Val;
}
}
if (k>0)
{
if (u < m_aNode[i] || u > m_aNode[i+k+1])
{
return Val;
}
else
{
double alpha = 0.0;
double beta = 0.0;
double dTemp = 0.0;
dTemp = m_aNode[i+k] - m_aNode[i];
if (dTemp == 0.0)
{
alpha = 0;
}
else
alpha = (u - m_aNode[i])/dTemp;
dTemp = m_aNode[i+k+1] - m_aNode[i+1];
if (dTemp == 0.0)
{
beta = 0.0;
}
else
beta = (m_aNode[i+k+1] - u)/dTemp;
Val1 = alpha * GetBaseFunVal(u, i, k-1, m_aNode);
Val2 = beta * GetBaseFunVal(u, i+1, k-1, m_aNode);
Val = Val1 + Val2;
}
}
return Val;
}
上述功能模塊摘自於計算機輔助幾何設計與非均勻有理B樣條。已知B樣條的n+1控制點坐標，以及相應的節點向量，可求得對應的曲線方程。
先計算各個控制點的基函數
各個基函數的求解可根據上述的功能模塊求出。
下面是我的C++實現：曲線是二維的，三維的情況，就Z坐標做同X,Y求解方式相同即可。在求解的過程中，我自己在CAD上畫了個樣條曲線，然后通過GetBaseFunVal(double u, int i, int k, vector <double> m_aNode)和頂點坐標，及節點向量求各個點的坐標。隨着u值的變化，計算各個X,Y,Z值。一個星期的摸爬滾打中，能輸出圖形，但是與原來的圖形對應不上。最終找到的原因在與基函數出問題了。在書本等相關資源中，基函數成員中的k表示的是次數，在我畫的樣條曲線中，階數顯示為3（為什么是3？CAD的標注里，實體塊中的 70 下一行，為3），所以我理所當然的寫為了2,。一直有問題。我將它改為3以后，竟然奇跡般的可以用了。而且跟原來圖形吻合。這個是我的相關經歷，希望對你們能有用。另外，哪位熱心人士可以說明下，為什么k改為階數大小，就可以呢？
#include <iostream>
#include <fstream>
#include <afxtempl.h>
using namespace std;
struct tPoint
{
double x;
double y;
double z;
};
double GetBaseFunVal(double u, int i, int k, vector <double> m_aNode)
{
double Val = 0.0;
double Val1 = 0.0;
double Val2 = 0.0;
if (k==0)
{
if (u < m_aNode[i] || u > m_aNode[i+1])
return Val;
else
{
Val = 1.0;
return Val;
}
}
if (k>0)
{
if (u < m_aNode[i] || u > m_aNode[i+k+1])
{
return Val;
}
else
{
double alpha = 0.0;
double beta = 0.0;
double dTemp = 0.0;
dTemp = m_aNode[i+k] - m_aNode[i];
if (dTemp == 0.0)
{
alpha = 0;
}
else
alpha = (u - m_aNode[i])/dTemp;
dTemp = m_aNode[i+k+1] - m_aNode[i+1];
if (dTemp == 0.0)
{
beta = 0.0;
}
else
beta = (m_aNode[i+k+1] - u)/dTemp;
Val1 = alpha * GetBaseFunVal(u, i, k-1, m_aNode);
Val2 = beta * GetBaseFunVal(u, i+1, k-1, m_aNode);
Val = Val1 + Val2;
}
}
return Val;
}
int main()
{
tPoint tData;
vector <tPoint> vtData;
vtData.clear();
vector <double> nodeVector;
nodeVector.push_back(0);
nodeVector.push_back(0);
nodeVector.push_back(0);
nodeVector.push_back(0);
nodeVector.push_back(1);
nodeVector.push_back(2);
nodeVector.push_back(3);
nodeVector.push_back(4);
nodeVector.push_back(5);
nodeVector.push_back(6);
nodeVector.push_back(6);
nodeVector.push_back(6);
nodeVector.push_back(6);
//節點向量nodeVector， 控制點坐標(0,3),(200,100), (750, 200), k=2
for (double u = 0; u < 6; u=u+0.01)
{
// 樣條的數據
tData.x = (GetBaseFunVal(u, 0, 3, nodeVector)*(-7585) + GetBaseFunVal(u, 1, 3, nodeVector)*(-3427.5) + GetBaseFunVal(u, 2, 3, nodeVector)*46087.5
+ GetBaseFunVal(u, 3, 3, nodeVector)*9220.0 + GetBaseFunVal(u, 4, 3, nodeVector)*(-14835.0) + GetBaseFunVal(u, 5, 3, nodeVector)*(-2002.5) + GetBaseFunVal(u, 6, 3, nodeVector)*71975
+ GetBaseFunVal(u, 7, 3, nodeVector)*45235 + GetBaseFunVal(u, 8, 3, nodeVector)*83150)/*/(GetBaseFunVal(u, 0, 3, nodeVector) + GetBaseFunVal(u, 1, 3, nodeVector) + GetBaseFunVal(u, 2, 3, nodeVector)
+ GetBaseFunVal(u, 3, 3, nodeVector) + GetBaseFunVal(u, 4, 3, nodeVector) + GetBaseFunVal(u, 5, 3, nodeVector))*/;
tData.y = (GetBaseFunVal(u, 0, 3, nodeVector)*(-3807.5) + GetBaseFunVal(u, 1, 3, nodeVector)*(19850.0) + GetBaseFunVal(u, 2, 3, nodeVector)*14335
+ GetBaseFunVal(u, 3, 3, nodeVector)*(-17582.5) + GetBaseFunVal(u, 4, 3, nodeVector)*(-5445.0) + GetBaseFunVal(u, 5, 3, nodeVector)*(-80735.0) + GetBaseFunVal(u, 6, 3, nodeVector)*(-23817.5)
+ GetBaseFunVal(u, 7, 3, nodeVector)*5037.5 + GetBaseFunVal(u, 8, 3, nodeVector)*(-9360))/*/(GetBaseFunVal(u, 0, 3, nodeVector) + GetBaseFunVal(u, 1, 3, nodeVector) + GetBaseFunVal(u, 2, 3, nodeVector)
+ GetBaseFunVal(u, 3, 3, nodeVector) + GetBaseFunVal(u, 4, 3, nodeVector) + GetBaseFunVal(u, 5, 3, nodeVector))*/;
tData.z = 0.0;
vtData.push_back(tData);
}
char *file = "C:/Users/Monkey/Desktop/新建文件夾 (2)/TEST/Last.txt";
ofstream out(file);
if (!out)
{
cout << "打開文件失敗！！！！" << endl;
}
for (int n = 0; n < vtData.size(); n++)
{
out << vtData[n].x << " " << vtData[n].y <<" " << vtData[n].z << endl;
}
out.close();
return 0;
}