bspline.cpp

建立b样条函数基函数的程序！ 这是从the nurbs book 翻译过来的版本！

#include <iostream>
#include <cmath>
using namespace std;

double OneBasicFuns(int p,int m,double U[],int i,double u);
int main()
{
	const int m=11;
	int p=2;
	int n=m-p-1;
	double u=2.5;
	double U[m]={0,0,0,1,2,3,4,4,5,5,5};
	
	double UU[m];
	for(int i=0;i<n;i=i+1)
	{
		UU[i]=OneBasicFuns(p,m,U,i,u);
	}
	return 0;
}


double OneBasicFuns(int p,int m,double U[],int i,double u)
{
	// 此函数算法见<<the nurbs books>>  P74
	double N[30];
	double saved,temp;
	int j,k;
	double Uleft,Uright;
	if( ((i==0) && (u==U[0])) || ((i==m-p-2) && (u==U[m-1])) )
	{
		return(1.0);
	}

	if(u<U[i] || u>=U[i+p+1])
	{
		return (0);
	}
	
	for(j=0;j<=p;j=j+1)
	{
		if(u>=U[i+j] && u<U[i+j+1])
		{
			N[j]=1.0;
		}
		else
		{
			N[j]=0.0;
		}
	}
	for(k=1;k<=p;k=k+1)
	{
		if(N[0]==0.0)
		{
			saved=0.0;
		}
		else
		{
			saved=(u-U[i])*N[0]/(U[i+k]-U[i]);
		}
		for(j=0;j<p-k+1;j=j+1)
		{
			Uleft=U[i+j+1];
			Uright=U[i+j+k+1];
			if(N[j+1]==0)
			{
				N[j]=saved;
				saved=0.0;
			}
			else
			{
				temp=N[j+1]/(Uright-Uleft);
				N[j]=saved+(Uright-u)*temp;
				saved=(u-Uleft)*temp;
			}
		}
	}
	return(N[0]);
}