B样条数据点反求控制点绘制曲线（源码）
作者：互联网
一、软件功能需求
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;
}
}
运行结果
致谢：强强（B样条原理）、瑶哥（矩阵求逆，递推公式绘图）
源码：https://download.csdn.net/download/qq_32059343/10927125
————————————————
版权声明：本文为CSDN博主「东辰叶落」的原创文章，遵循 CC 4.0 BY-SA 版权协议，转载请附上原文出处链接及本声明。
原文链接：https://blog.csdn.net/qq_32059343/article/details/86408359
Kernel-Tomposon
发布了129 篇原创文章 · 获赞 116 · 访问量 65万+
他的留言板
关注
标签：palette,endPoint,matrix,样条,反求,int,double,detal,源码 来源： https://blog.csdn.net/shixin_0125/article/details/104004008