Jump to content
  • Advertisement
Sign in to follow this  

TCB spline interpolation

This topic is 3778 days old which is more than the 365 day threshold we allow for new replies. Please post a new topic.

If you intended to correct an error in the post then please contact us.

Recommended Posts

Okay, after googling for a while, and not finding anything useful, I decided to try my luck here. The title basically covers it all, I'm looking for the algorithm (or calculations) to handle TCB spline interpolation. The thing is, the data comes from Collada, which in turn, comes from 3D Studio Max. The difference between regular TCB and the one used in 3d studio max is that it also introduces an "ease-in" and "ease-out" variable. My guess is that that has something to do with the "s" variable of the hermite spline calculation. (my reference: http://www.cubic.org/docs/hermite.htm) EDIT: Code updated, works okay now as long as you don't use the easing :)
inline float TCBTangentEquationIncoming(float& pm1, float& p0, float& p1, float& t, float& c, float& b)
	return ((1-t)*(1-c)*(1+b))/2*(p0-pm1) + ((1-t)*(1+c)*(1-b))/2*(p1-p0);

inline float TCBTangentEquationOutgoing(float& pm1, float& p0, float& p1, float& t, float& c, float& b)
	return ((1-t)*(1+c)*(1+b))/2*(p0-pm1) + ((1-t)*(1-c)*(1-b))/2*(p1-p0);

float* CMatrix::ConcatFloat4Array(float* farr)
	float* V=new float[4];

	V[0] = farr[0]*Matrix[0][0] + farr[1]*Matrix[1][0] + farr[2]*Matrix[2][0] + farr[3]*Matrix[3][0];
	V[1] = farr[0]*Matrix[0][1] + farr[1]*Matrix[1][1] + farr[2]*Matrix[2][1] + farr[3]*Matrix[3][1];
	V[2] = farr[0]*Matrix[0][2] + farr[1]*Matrix[1][2] + farr[2]*Matrix[2][2] + farr[3]*Matrix[3][2];
	V[3] = farr[0]*Matrix[0][3] + farr[1]*Matrix[1][3] + farr[2]*Matrix[2][3] + farr[3]*Matrix[3][3];

	return V;

CMatrix& GetHermiteMatrix()
	static bool matrixInited=false;
	static CMatrix hm;
		return hm;

	hm.Matrix[0][0] = 2; hm.Matrix[1][0] = -2; hm.Matrix[2][0] = 1; hm.Matrix[3][0] = 1;
	hm.Matrix[0][1] = -3; hm.Matrix[1][1] = 3; hm.Matrix[2][1] = -2; hm.Matrix[3][1] = -1;
	hm.Matrix[0][2] = 0; hm.Matrix[1][2] = 0; hm.Matrix[2][2] = 1; hm.Matrix[3][2] = 0;
	hm.Matrix[0][3] = 1; hm.Matrix[1][3] = 0; hm.Matrix[2][3] = 0; hm.Matrix[3][3] = 0;

	matrixInited = true;
	return hm;

//! pm1 (y coordinate of the previous point)
//! p0 (y coordinate of the current point)
//! p1 (y coordinate of the next point)
//! p2 (y coordinate of the next+1 point)
//! s (interpolation value, use ApproximateCubicBezierParameter?...)
//! to0 (tangent out of (p0))
//! ti1 (tangent in of (p1))
float InterpolateTCB(float pm1, float p0, float p1, float p2, float s, float to0, float ti1)
	float S[4], C[4];
	S[3] = 1;
	S[2] = s; // s
	S[1] = s*s; // s^2
	S[0] = S[1]*s; // s^3

	C[0] = p0;
	C[1] = p1;
	C[2] = to0;
	C[3] = ti1;

        // EDIT: proper TCB interpolation code
	float* res=GetHermiteMatrix().ConcatFloat4Array(C);
	float ret=res[0]*S[0]+res[1]*S[1]+res[2]*S[2]+res[3]*S[3];

	delete res; 
	return ret;

float CFormula::Value(float X)
	unsigned int ps=PointCount();
	if(ps == 0)
		return 0.0f;
	else if(X < points[0].x)
		return points[0].y;
	else if(X >= points[ps - 1].x)
		return points[ps - 1].y;

	int segment=this->PointSegmentLocation(X);

	if(InterpolationType == INTERPOLATE_CUBICBEZIER)
		InterpolatePoint& p1=points[segment];
		InterpolatePoint& p2=points[segment+1];

		float s = ApproximateCubicBezierParameter( X, p1.x, p1.tangentOutX, p2.tangentInX, p2.x );

		float s1 = (1-s);
		float s2 = s1*s1;
		float ss2 = s*s;

		return p1.y*s2*s1+3*p1.tangentOutY*s*s2+3*p2.tangentInY*ss2*s1+p2.y*ss2*s;
	else if (InterpolationType == INTERPOLATE_TCB)
                // EDIT: it actually seems that this is also the method Collada/3DSM intended, as my test cases show identical movement at the beginning and end of the animation.
		InterpolatePoint& pm1=points[segment-1<0?segment:segment-1];
		InterpolatePoint& p0=points[segment];
		InterpolatePoint& p1=points[segment+1];
		InterpolatePoint& p2=points[segment+2>=PointCount()?segment+1:segment+2];

                // collada inverts these parameters for some reason, don't ask me why..
		float t= -p0.tangentInX; 
		float c= -p0.tangentInY;
		float b= -p0.tangentOutX;

                // ease variables.
                //float& easein= p0.tangentOutY;
                //float& easeout= p0.easeout;

		float to0x = TCBTangentEquationOutgoing(pm1.x, p0.x, p1.x, t, c, b);
		float to0y = TCBTangentEquationOutgoing(pm1.y, p0.y, p1.y, t, c, b);
		float ti1x = TCBTangentEquationIncoming(p0.x, p1.x, p2.x, t, c, b);
		float ti1y = TCBTangentEquationIncoming(p0.y, p1.y, p2.y, t, c, b);

                // EDIT: it is, and it just seems to be the Factor variable (linear interpolated range). Now it only needs to have easing implemented...
	        float Factor=(X-points[segment].x) / (points[segment+1].x-points[segment].x);

		return InterpolateTCB(pm1.y, p0.y, p1.y, p2.y, Factor, to0y, ti1y);

	float Factor=(X-points[segment].x) / (points[segment+1].x-points[segment].x);

		return InterpolateLinear(points[segment].y, points[segment+1].y, Factor);
		return InterpolateCosine(points[segment].y, points[segment+1].y, Factor);
		float y0 = points[segment-1<0?segment:segment-1].y, y3=points[segment+2>=PointCount()?segment+1:segment+2].y;
		return InterpolateCubic(y0, points[segment].y, points[segment+1].y, y3, Factor);
		float y0 = points[segment-1<0?segment:segment-1].y, y3=points[segment+2>=PointCount()?segment+1:segment+2].y;
		return InterpolateCatmullRom(y0, points[segment].y, points[segment+1].y, y3, Factor);
		float y0 = points[segment-1<0?segment:segment-1].y, y3=points[segment+2>=PointCount()?segment+1:segment+2].y;
		return InterpolateHermite(y0, points[segment].y, points[segment+1].y, y3, HermiteTension, HermiteBias, Factor);

//simply clamps a value between 0 .. 1

float ClampToZeroOne(float value) {
   if (value < .0f)
      return .0f;
   else if (value > 1.0f)
      return 1.0f;
      return value;

// the following function i got from the collada specification, but it's intended for the "BEZIER" type of interpolation, which works fine by the way..

 * Returns the approximated parameter of a parametric curve for the value X
 * @param atX At which value should the parameter be evaluated
 * @param P0_X The first interpolation point of a curve segment
 * @param C0_X The first control point of a curve segment
 * @param C1_X The second control point of a curve segment
 * @param P1_x The second interpolation point of a curve segment
 * @return The parametric argument that is used to retrieve atX using the parametric function representation of this curve

float ApproximateCubicBezierParameter (
         float atX, float P0_X, float C0_X, float C1_X, float P1_X ) {
   if (atX - P0_X < VERYSMALL)
      return 0.0;
   if (P1_X - atX < VERYSMALL) 
      return 1.0;
   long iterationStep = 0;
   float u = 0.0f; float v = 1.0f;
   //iteratively apply subdivision to approach value atX
   while (iterationStep < MAXIMUM_ITERATIONS) {
      // de Casteljau Subdivision.
      double a = (P0_X + C0_X)*0.5f;
      double b = (C0_X + C1_X)*0.5f;
      double c = (C1_X + P1_X)*0.5f;
      double d = (a + b)*0.5f;
      double e = (b + c)*0.5f;
      double f = (d + e)*0.5f; //this one is on the curve!
      //The curve point is close enough to our wanted atX
      if (fabs(f - atX) < APPROXIMATION_EPSILON)
         return ClampToZeroOne((u + v)*0.5f);
      if (f < atX) {
         P0_X = f;
         C0_X = e;
         C1_X = c;
         u = (u + v)*0.5f;
      } else {
         C0_X = a; C1_X = d; P1_X = f; v = (u + v)*0.5f;
   return ClampToZeroOne((u + v)*0.5f);

[Edited by - Seniltai on July 10, 2008 10:35:02 AM]

Share this post

Link to post
Share on other sites
Well, I think I figured it out. The above code is updated accordingly. The linear sampled test case and the original with TCB and BEZIER interpolation are now virtually indistinguishable afaik...

However, the only problem still residing are the easing parameters, as no source mentions anything about this. My guess is that it's probably not inherent to splines, and it's a separate equation.

This is what I managed to dig out using google:

Now the only question remaining, which equation is the one, if any... *sigh*
Also, it's got an easing factor, which determines the strength of the easing (and even worse, there are actually two, ease to, and ease from)

Help, please? :(

Share this post

Link to post
Share on other sites
Sign in to follow this  

  • Advertisement

Important Information

By using GameDev.net, you agree to our community Guidelines, Terms of Use, and Privacy Policy.

We are the game development community.

Whether you are an indie, hobbyist, AAA developer, or just trying to learn, GameDev.net is the place for you to learn, share, and connect with the games industry. Learn more About Us or sign up!

Sign me up!