Rational B-Spline

Vd4A	b_spline	(const Vd4& q0,const Vd4& q1,const Vd4& q2,const Vd4& q3,const long div_c)
{
	Vd4A	vd3a ;
	for (long index=0 ; index<=div_c ; index++) {
		double	t	= 1./div_c*index ;
		double	n0	= 1./6 * (1.-t)*(1.-t)*(1.-t) ;
		double	n1	= 1./2 * t*t*t - t*t + 2./3 ;
		double	n2	=-1./2 * t*t*t + 1./2 * t*t + 1./2 *t + 1./6 ;
		double	n3	= 1./6 * t*t*t ;
		double	px = ( (n0*q0.x)*q0.w + (n1*q1.x)*q1.w + (n2*q2.x)*q2.w + (n3*q3.x)*q3.w ) / ((n0)*q0.w + (n1)*q1.w + (n2)*q2.w + (n3)*q3.w) ;
		double	py = ( (n0*q0.y)*q0.w + (n1*q1.y)*q1.w + (n2*q2.y)*q2.w + (n3*q3.y)*q3.w ) / ((n0)*q0.w + (n1)*q1.w + (n2)*q2.w + (n3)*q3.w) ;
		double	pz = ( (n0*q0.z)*q0.w + (n1*q1.z)*q1.w + (n2*q2.z)*q2.w + (n3*q3.z)*q3.w ) / ((n0)*q0.w + (n1)*q1.w + (n2)*q2.w + (n3)*q3.w) ;
		double	pw = ( (n0*q0.w)*q0.w + (n1*q1.w)*q1.w + (n2*q2.w)*q2.w + (n3*q3.w)*q3.w ) / ((n0)*q0.w + (n1)*q1.w + (n2)*q2.w + (n3)*q3.w) ;
		vd3a.push_back(Vd4(px,py,pz,pw)) ;
		}
	return	vd3a ;
	}

pw は，1 とすべきなのか，これで良いのかは不明．

次の P の配列に対してのロジックは前回のものと同様で，Vd4 に変更．
  Vd4A	b_spline	(const   Vd4A&   pnts_,const long div_n=2)
v_Vd4A	b_spline	(const v_Vd4A& v_pnts_,const long div_n=2)


与えたデータは，
  10, 0, 0,  1
  10,10, 0,  1
   0,10, 0,  1
// x, y, z,  w
  10,10,0 の w の部分を 0.0 , 0.5 , 0.7 , 1 , 2 , 3 としたもの