一、軟體功能需求
1）所設計的軟體應具有圖形化使用者界面（GUI）；
2）使用者在軟體界面上可用随機數方式或手工方式輸入若幹曲線或曲面的資料點，例如起點、終點、清單型值點等，對于曲線，還可設定步長參數；對于曲面，還可設定步長與行距參數；曲線或曲面的類型不限。
3）具有“繪制理想圖形”按鈕，使用者完成資料點與參數輸入後，點選該按鈕，軟體可繪制出理想的曲線或曲面，若該曲線或曲面有特征多邊形，則還能通過型值點反算出控制頂點，并繪制出對應的特征多邊形；
4）具有“生成加工路徑”按鈕，點選該按鈕，可根據設定的步長與行距參數自動生成加工該曲線或曲面的刀具路徑；
5）具有“生成加工程式”按鈕，點選該按鈕，可自動生成一個實作該曲線或曲面加工的G代碼數控加工程式檔案。
二、B樣條反求控制點原理
參考書籍：計算機輔助幾何設計與非均勻有理B樣條（一般學校的電子圖書館可以線上閱讀）
P212頁：B樣條遞推定義
P237頁：齊次B樣條曲線的節點矢量
P257頁：反算三次B樣條插值曲線的控制頂點
P259頁：邊界條件：抛物線條件
二、軟體實作
軟體平台：Windows7 Qt5.9.3 語言：C++
由于預設界面框的原點在左上角，并且向下為Y+，向右為x+，是以在畫圖或者滑鼠捕捉時，要将其轉化到正常的直角坐标系，并繪制坐标系。
QPainter* painter = new QPainter(this);
painter->setWindow(-width()/2,height()/2,width(),-height()); //設定坐标原點在螢幕的正中心
// painter->drawLine(-296,0,296,0);
// painter->drawLine(0,296,0,-296);
//繪制x軸
QPointF startPoint(-295,0);
QPointF endPoint(295,0);
double x1, y1, x2, y2; //箭頭的兩點坐标
//求得箭頭兩點坐标
calcVertexes(startPoint.x(), startPoint.y(), endPoint.x(), endPoint.y(), x1, y1, x2, y2);
painter->drawLine(startPoint, endPoint);//繪制線段
painter->drawLine(endPoint.x(), endPoint.y(), x1, y1);//繪制箭頭一半
painter->drawLine(endPoint.x(), endPoint.y(), x2, y2);//繪制箭頭另一半
//繪制Y軸
startPoint.setX(0);
startPoint.setY(-295);
endPoint.setX(0);
endPoint.setY(295);
//求得箭頭兩點坐标
calcVertexes(startPoint.x(), startPoint.y(), endPoint.x(), endPoint.y(), x1, y1, x2, y2);
painter->drawLine(startPoint, endPoint);//繪制線段
painter->drawLine(endPoint.x(), endPoint.y(), x1, y1);//繪制箭頭一半
painter->drawLine(endPoint.x(), endPoint.y(), x2, y2);//繪制箭頭另一半
void Palette::calcVertexes(double start_x, double start_y, double end_x, double end_y, double& x1, double& y1, double& x2, double& y2)
{
double arrow_lenght_ = 10;//箭頭長度，一般固定
double arrow_degrees_ = 0.5;//箭頭角度，一般固定
double angle = atan2(end_y - start_y, end_x - start_x) + 3.1415926;//
x1 = end_x + arrow_lenght_ * cos(angle - arrow_degrees_);//求得箭頭點1坐标
y1 = end_y + arrow_lenght_ * sin(angle - arrow_degrees_);
x2 = end_x + arrow_lenght_ * cos(angle + arrow_degrees_);//求得箭頭點2坐标
y2 = end_y + arrow_lenght_ * sin(angle + arrow_degrees_);
} 采用随機數生成8個随機點，對第一組随機數進行排序，第二組随機數不排序，填充到QTableWidget。
void MainWindow::on_general_random_clicked()
{
QVector <int> a; //定義一個vector容器
srand((int)time(NULL)); //時間作為種子
for (int i=0;i<8;i++)
{
a.append(rand()%296-148); //向vector後面追加元素，區間根據顯示區域設定
}
qSort(a.begin(),a.end()); //對vector從小到大排序
for(int i=0;i<8;i++) //向X坐标寫入資料
{
ui->model->setItem(i,0,new QTableWidgetItem(QString("%1").arg(a[i])));
}
a.remove(0,8); //删去第一列
srand((int)time(NULL));
for (int i=0;i<8;i++)
{
a.append(rand()%296-148);
}
for(int i=0;i<8;i++)
{
ui->model->setItem(i,1,new QTableWidgetItem(QString("%1").arg(a[i])));
ui->model->setItem(i,2,new QTableWidgetItem(QString("%1").arg(0)));
}
} 擷取資料點，根據資料點反求控制點。
void MainWindow::on_general_curve_clicked()
{
QPoint point;
QPoint start,end;
QString QTW_data1,QTW_data2;
double a1,a2;
//提取資料表中資料
for(int i=0; i<8; i++)
{
QTW_data1=ui->model->item(i, 0)->text();
QTW_data2=ui->model->item(i, 1)->text();
a1=QTW_data1.toDouble();
a2=QTW_data2.toDouble();
point.setX(a1);
point.setY(a2);
if(i==0)
{
start.setX(a1);
start.setY(a2);
}
if(i==7)
{
end.setX(a1);
end.setY(a2);
}
palette->typepoints.append(point);
}
/*
此段程式将資料點反求出控制點,八個資料點
*/
double l[7]; //存放資料點之間的距離，l[0]代表第一個資料點到第二個資料點之間的距離
double L=0;//存放總長度
for(int i=0;i<7;i++)
{
l[i]=sqrt(pow((palette->typepoints[i+1].x()-palette->typepoints[i].x()),2)+
pow((palette->typepoints[i+1].y()-palette->typepoints[i].y()),2)); //兩點之間距離公式
L=L+l[i];
}
//節點矢量計算
for(int i=0;i<4;i++)
{
palette->u[i]=0;
}
for(int i=0;i<4;i++)
{
palette->u[i+10]=1;
}
palette->u[4]=l[0]/L;
palette->u[5]=(l[0]+l[1])/L;
palette->u[6]=(l[0]+l[1]+l[2])/L;
palette->u[7]=(l[0]+l[1]+l[2]+l[3])/L;
palette->u[8]=(l[0]+l[1]+l[2]+l[3]+l[4])/L;
palette->u[9]=(l[0]+l[1]+l[2]+l[3]+l[4]+l[5])/L;
for(int i=0;i<14;i++)
{
qDebug()<<"u[i]="<<palette->u[i]<<endl;
}
double detal[13];
for(int i=0;i<13;i++)
{
detal[i]=palette->u[i+1]-palette->u[i];
qDebug()<<"detal[i]="<<detal[i]<<endl;
}
//根據節點矢量，反算控制點，其中需要使用矩陣求逆
double a[8],b[8],c[8],e[8],f[8]; //角标0-8
a[0]=1-(detal[3]*detal[4])/pow((detal[3]+detal[4]),2);
b[0]=(detal[3]/(detal[3]+detal[4]))*(detal[4]/(detal[3]+detal[4])-detal[3]/(detal[3]+detal[4]+detal[5]));
c[0]=(pow((detal[3]),2))/((detal[3]+detal[4])*(detal[3]+detal[4]+detal[5]));
e[0]=(1.0/3)*(palette->typepoints[0].x()+2*palette->typepoints[1].x()); //X軸坐标
f[0]=(1.0/3)*(palette->typepoints[0].y()+2*palette->typepoints[1].y()); //Y軸坐标
a[7]=-(pow((detal[9]),2))/((detal[8]+detal[9])*(detal[8]+detal[8]+detal[9]));
b[7]=(detal[9]/(detal[8]+detal[9]))*(detal[9]/(detal[8]+detal[8]+detal[9])-detal[8]/(detal[8]+detal[9]));
c[7]=detal[8]*detal[9]/pow((detal[8]+detal[9]),2)-1;
e[7]=-(1.0/3)*(palette->typepoints[7].x()+2*palette->typepoints[6].x());
f[7]=-(1.0/3)*(palette->typepoints[7].y()+2*palette->typepoints[6].y());
for(int i=1;i<7;i++)
{
// e[8]=palette->typepoints[i].x();
a[i]=(pow((detal[i+3]),2))/(detal[i+1]+detal[i+2]+detal[i+3]);
b[i]=detal[i+3]*(detal[i+1]+detal[i+2])/(detal[i+1]+detal[i+2]+detal[i+3])+detal[i+2]*(detal[i+3]+detal[i+4])/(detal[i+2]+detal[i+3]+detal[i+4]);
c[i]=pow((detal[i+2]),2)/(detal[i+2]+detal[i+3]+detal[i+4]);
e[i]=(detal[i+2]+detal[i+3])*palette->typepoints[i].x();
f[i]=(detal[i+2]+detal[i+3])*palette->typepoints[i].y();
}
for(int i=0;i<8;i++)
{
qDebug()<<"a[i]="<<a[i]<<"b[i]="<<b[i]<<"c[i]="<<c[i]<<"e[i]"<<e[i]<<endl;
}
double matrix[8][8];
for(int i=0;i<8;i++)
{
for(int j=0;j<8;j++)
{
matrix[i][j]=0;
qDebug()<<matrix[i][j]<<" ";
}
qDebug()<<endl;
}
qDebug()<<"_________________________"<<endl;
/*對系數矩陣指派*/
matrix[0][0]=a[0];
matrix[0][1]=b[0];
matrix[0][2]=c[0];
for(int i=1;i<7;i++)
{
matrix[i][i-1]=a[i];
matrix[i][i]=b[i];
matrix[i][i+1]=c[i];
}
matrix[7][5]=a[7];
matrix[7][6]=b[7];
matrix[7][7]=c[7];
for(int i=0;i<8;i++)
{
for(int j=0;j<8;j++)
{
qDebug()<<matrix[i][j]<<" ";
}
}
qDebug()<<endl;
double matrix_after[N][N];
qDebug()<<"_________________"<<endl;
matrix_inv(matrix,matrix_after); //求逆矩陣
for(int i=0;i<8;i++)
{
for(int j=0;j<8;j++)
{
qDebug()<<matrix_after[i][j]<<" ";
}
}
qDebug()<<endl;
double sum_x[8]={},sum_y[8]={};
palette->ctrlPoints.append(start);
for(int i=0;i<8;i++)
{
for(int j=0;j<8;j++)
{
sum_x[i]=sum_x[i]+matrix_after[i][j]*e[j];
sum_y[i]=sum_y[i]+matrix_after[i][j]*f[j];
}
palette->ctrlPoints.append(QPoint(sum_x[i],sum_y[i]) );
qDebug()<<"sum_x[i]="<<sum_x[i]<<"sum_y[i]="<<sum_y[i]<<endl;
}
palette->ctrlPoints.append(end);
qDebug()<<palette->ctrlPoints<<endl;
palette->generateCurve(); //根據控制點繪制型值點
palette->update(); //更新畫闆
} 矩陣求逆C++實作：
void MainWindow::matrix_inv(double matrix_before[N][N],double matrix_after[N][N])
{
bool flag;//标志位，如果行列式為0，則結束程式
flag=GetMatrixInverse(matrix_before,N,matrix_after);
}
//按第一行展開計算|A|
double MainWindow::getA(double arcs[N][N],int n)
{
if(n==1)
{
return arcs[0][0];
}
double ans = 0;
double temp[N][N]={0.0};
int i,j,k;
for(i=0;i<n;i++)
{
for(j=0;j<n-1;j++)
{
for(k=0;k<n-1;k++)
{
temp[j][k] = arcs[j+1][(k>=i)?k+1:k];
}
}
double t = getA(temp,n-1);
if(i%2==0)
{
ans += arcs[0][i]*t;
}
else
{
ans -= arcs[0][i]*t;
}
}
return ans;
}
//計算每一行每一列的每個元素所對應的餘子式，組成A*
void MainWindow::getAStart(double arcs[N][N],int n,double ans[N][N])
{
if(n==1)
{
ans[0][0] = 1;
return;
}
int i,j,k,t;
double temp[N][N];
for(i=0;i<n;i++)
{
for(j=0;j<n;j++)
{
for(k=0;k<n-1;k++)
{
for(t=0;t<n-1;t++)
{
temp[k][t] = arcs[k>=i?k+1:k][t>=j?t+1:t];
}
}
ans[j][i] = getA(temp,n-1); //此處順便進行了轉置
if((i+j)%2 == 1)
{
ans[j][i] = - ans[j][i];
}
}
}
}
//得到給定矩陣src的逆矩陣儲存到des中。
bool MainWindow::GetMatrixInverse(double src[N][N],int n,double des[N][N])
{
double flag=getA(src,n);
double t[N][N];
if(0==flag)
{
// cout<< "原矩陣行列式為0，無法求逆。請重新運作" <<endl;
return false;//如果算出矩陣的行列式為0，則不往下進行
}
else
{
getAStart(src,n,t);
for(int i=0;i<n;i++)
{
for(int j=0;j<n;j++)
{
des[i][j]=t[i][j]/flag;
}
}
}
return true;
} 根據控制點繪制B樣條曲線：
void Palette::generateCurve()
{
curvePoints.clear();
int i=3;
for(double u1=0; u1<1; u1+=0.001)
{
QPointF tmp(0,0);
if(u1>u[i+1])
{i=i+1;}
for(int k=0; k<4;k++)
{
QPointF t = ctrlPoints[i-k];
t*=Nu(i-k,u1,3);
tmp+=t;
}
//qDebug()<<"tmp="<<tmp;
curvePoints.push_back(tmp);
}
}
double Palette::Nu(int i,double u2,int k)
{
if(k==0)
{
if((u2>=u[i])&&(u2<u[i+1]))
return 1;
else
return 0;
}
else{
double a1=(u2-u[i]),c1=(u[i+k]-u[i]);
double b1=(u[i+k+1]-u2),d1=(u[i+k+1]-u[i+1]);
if(d1<0.001)
{
b1=0;
d1=1;
}
if(c1<0.001)
{
a1=0;
c1=1;
}
double nu=(a1/c1)*Nu(i,u2,k-1)+(b1/d1)*Nu(i+1,u2,k-1);
return nu;
}
} 運作結果