/*
								+----------------------------------+
								|                                  |
								| *** Kochanek-Bartels splines *** |
								|                                  |
								|   Copyright © -tHE SWINe- 2014   |
								|                                  |
								|        KochanekBartels.h         |
								|                                  |
								+----------------------------------+
*/

#pragma once
#ifndef __KOCHANEK_BARTELS_SPLINE_INCLUDED
#define __KOCHANEK_BARTELS_SPLINE_INCLUDED

/**
 *	@file lml/KochanekBartels.h
 *	@date 2014
 *	@author -tHE SWINe-
 *	@brief Kochanek-Bartels splines
 */

#include "Spline.h"

/**
 *	@brief simple Kochanek-Bartels spline implementation
 */
template <class CPointClass = Vector3f>
class CKochanekBartelsSplineImpl {
public:
	typedef CPointClass _TyPoint; /**< @brief type of points the spline is defined on */
	typedef typename CPointClass::_TyScalar _TyScalar; /**< @brief type of scalar */

	/**
	 *	@brief spline traits, stored as enum
	 */
	enum {
		b_has_explicit_time = true, /**< @brief the spline keypoints do have f_time member */
		b_has_rotation = false /**< @brief the spline keypoints do not have t_rotation member */
	};

	/**
	 *	@brief spline keypoint type
	 */
	struct _TyKeyPoint {
		_TyPoint v_position; /**< @brief interpolated quantity */
		_TyScalar f_time; /**< @brief explicit keypoint time */
		_TyScalar f_tension; /**< @brief tension */
		_TyScalar f_continuity; /**< @brief continuity */
		_TyScalar f_bias; /**< @brief bias */

		/**
		 *	@brief default constructor; has no effect
		 */
		inline _TyKeyPoint()
		{}

		/**
		 *	@brief constructor
		 *
		 *	@param[in] r_v_position is interpolated quantity
		 *	@param[in] _f_time is keypoint time
		 *	@param[in] _f_tension is tension
		 *	@param[in] _f_continuity is continuity
		 *	@param[in] _f_bias is bias
		 */
		inline _TyKeyPoint(const _TyPoint &r_v_position, _TyScalar _f_time,
			_TyScalar _f_tension = 0, _TyScalar _f_continuity = 0, _TyScalar _f_bias = 0)
			:v_position(r_v_position), f_time(_f_time), f_tension(_f_tension),
			f_continuity(_f_continuity), f_bias(_f_bias)
		{}
	};

protected:
	std::vector<_TyKeyPoint> m_point_list;

public:
	/**
	 *	@brief gets number of spline arcs
	 *	@return Returns number of spline arcs.
	 */
	inline size_t n_Arc_Num() const
	{
		return (m_point_list.size() > 3)? m_point_list.size() - 3 : 0;
	}

	/**
	 *	@brief calculates time length of a spline arc
	 *	@param[in] n_arc is zero-based index of spline arc
	 *	@return Returns time length of the selected arc.
	 */
	inline _TyScalar f_Arc_Time(size_t n_arc) const
	{
		_ASSERTE(n_arc < n_Arc_Num());
		const _TyKeyPoint &r_pt1 = m_point_list[n_arc + 1];
		const _TyKeyPoint &r_pt2 = m_point_list[n_arc + 2];
		return r_pt2.f_time - r_pt1.f_time;
	}

	/**
	 *	@brief calculates time length of this spline
	 *	@return Returns time length of the whole spline.
	 */
	inline _TyScalar f_Spline_Time() const
	{
		if(!n_Arc_Num())
			return 0;
		const _TyKeyPoint &r_pt1 = m_point_list[1];
		const _TyKeyPoint &r_pt2 = m_point_list[m_point_list.size() - 2];
		return r_pt2.f_time - r_pt1.f_time;
	}

	/**
	 *	@brief calculates position on a spline arc
	 *
	 *	@param[in] f_t is position on a spline arc in range [0, 1]
	 *	@param[in] n_arc is zero-based index of spline arc
	 *
	 *	@return Returns point at the given position along the selected arc.
	 */
	_TyPoint v_Arc_Position(_TyScalar f_t, size_t n_arc) const
	{
		_ASSERTE(n_arc < n_Arc_Num());

		const _TyKeyPoint &r_pt0 = m_point_list[n_arc];
		const _TyKeyPoint &r_pt1 = m_point_list[n_arc + 1];
		const _TyKeyPoint &r_pt2 = m_point_list[n_arc + 2];
		const _TyKeyPoint &r_pt3 = m_point_list[n_arc + 3];
		// get four points

		_TyPoint v_tangent0, v_tangent1;
		Get_Tangents(v_tangent0, v_tangent1, r_pt0, r_pt1, r_pt2, r_pt3);
		// calculate the tangents at r_pt1 and r_pt2

		//_ASSERTE(f_t >= 0 && f_t <= 1); // don't, may be slightly out of bounds
		_TyScalar f_t2 = f_t * f_t;
		_TyScalar f_t3 = f_t2 * f_t;
		// calculate t^2, t^3 (the times are already in [0, 1] range)

		_TyScalar f_weight0 = 2 * f_t3 - 3 * f_t2 + 1;	// [ 2 -3  0  1]
		_TyScalar f_weight1 = -2 * f_t3 + 3 * f_t2;		// [-2  3  0  0] 
		_TyScalar f_weight2 = f_t3 - 2 * f_t2 + f_t;	// [ 1 -2  1  0]
		_TyScalar f_weight3 = f_t3 - f_t2;				// [ 1 -1  0  0]
		// calculate polynomials

		return r_pt1.v_position * f_weight0 + r_pt2.v_position * f_weight1 +
			v_tangent0 * f_weight2 + v_tangent1 * f_weight3;
		// interpolate
	}

	/**
	 *	@brief calculates tangent vector on a spline arc
	 *
	 *	@param[in] f_t is position on a spline arc in range [0, 1]
	 *	@param[in] n_arc is zero-based index of spline arc
	 *
	 *	@return Returns tangent at the given position along the selected arc.
	 */
	_TyPoint v_Arc_Derivative(_TyScalar f_t, size_t n_arc) const
	{
		_ASSERTE(n_arc < n_Arc_Num());

		const _TyKeyPoint &r_pt0 = m_point_list[n_arc];
		const _TyKeyPoint &r_pt1 = m_point_list[n_arc + 1];
		const _TyKeyPoint &r_pt2 = m_point_list[n_arc + 2];
		const _TyKeyPoint &r_pt3 = m_point_list[n_arc + 3];
		// get four points

		_TyPoint v_tangent0, v_tangent1;
		Get_Tangents(v_tangent0, v_tangent1, r_pt0, r_pt1, r_pt2, r_pt3);
		// calculate the tangents at r_pt1 and r_pt2

		//_ASSERTE(f_t >= 0 && f_t <= 1); // don't, may be slightly out of bounds
		_TyScalar f_t2 = f_t * f_t;
		//_TyScalar f_t3 = f_t2 * f_t; // unused
		// calculate t^2, t^3 (the times are already in [0, 1] range)

		_TyScalar f_weight0 = 6 * f_t2 - 6 * f_t;
		_TyScalar f_weight1 = -6 * f_t2 + 6 * f_t;
		_TyScalar f_weight2 = 3 * f_t2 - 4 * f_t + 1;
		_TyScalar f_weight3 = 3 * f_t2 - 2 * f_t;
		// calculate polynomials

		return r_pt1.v_position * f_weight0 + r_pt2.v_position * f_weight1 +
			v_tangent0 * f_weight2 + v_tangent1 * f_weight3;
		// interpolate
	}

protected:
	static void Get_Tangents(_TyPoint &r_v_tangent0, _TyPoint &r_v_tangent1, const _TyKeyPoint &r_pt0,
		const _TyKeyPoint &r_pt1, const _TyKeyPoint &r_pt2, const _TyKeyPoint &r_pt3)
	{
		_TyPoint v_diff_21 = r_pt2.v_position - r_pt1.v_position,
			v_diff_10 = r_pt1.v_position - r_pt0.v_position,
			v_diff_32 = r_pt3.v_position - r_pt2.v_position;
		{
			const _TyScalar t = r_pt1.f_tension, c = r_pt1.f_continuity, b = r_pt1.f_bias;
			_TyScalar f_dt_i1 = r_pt2.f_time - r_pt1.f_time, f_dt_i_plus_dt_i1 = r_pt2.f_time - r_pt0.f_time;
			_TyScalar f_t_corr = 2 * f_dt_i1 / f_dt_i_plus_dt_i1; // adjust for nonuniform time
			_TyScalar f_weight0 = f_t_corr * ((1 - t) * (1 + c) * (1 - b)) / 2;
			_TyScalar f_weight1 = f_t_corr * ((1 - t) * (1 - c) * (1 + b)) / 2;
			r_v_tangent0 = v_diff_21 * f_weight0 + v_diff_10 * f_weight1;
		}
		{
			const _TyScalar t = r_pt2.f_tension, c = r_pt2.f_continuity, b = r_pt2.f_bias;
			_TyScalar f_dt_i = r_pt2.f_time - r_pt1.f_time, f_dt_i_plus_dt_i1 = r_pt3.f_time - r_pt1.f_time;
			_TyScalar f_t_corr = 2 * f_dt_i / f_dt_i_plus_dt_i1; // adjust for nonuniform time
			_TyScalar f_weight0 = f_t_corr * ((1 - t) * (1 - c) * (1 - b)) / 2;
			_TyScalar f_weight1 = f_t_corr * ((1 - t) * (1 + c) * (1 + b)) / 2;
			r_v_tangent1 = v_diff_32 * f_weight0 + v_diff_21 * f_weight1;
		}
		// calculate the tangents at r_pt1 and r_pt2
	}
};

/**
 *	@brief implementation of Kochanek-Bartels spline with rotation
 */
template <class CPointClass = Vector3f>
class CKochanekBartelsQuatSplineImpl {
public:
	typedef CPointClass _TyPosition; /**< @brief type of position the spline is defined on */
	typedef typename CPointClass::_TyScalar _TyScalar; /**< @brief type of scalar */
	typedef Quat<_TyScalar> _TyRotation; /**< @brief quaternion type for rotation representation */

	/**
	 *	@brief spline traits, stored as enum
	 */
	enum {
		b_has_explicit_time = true, /**< @brief the spline keypoints do have f_time member */
		b_has_rotation = true /**< @brief the spline keypoints do not have t_rotation member */
	};

	/**
	 *	@brief type of poses the spline is defined on
	 */
	struct _TyPoint {
		typedef typename CPointClass::_TyScalar _TyScalar; /**< @brief type of scalar */

		_TyPosition v_position; /**< @brief interpolated position */
		_TyRotation t_rotation; /**< @brief interpolated rotation */

		/**
		 *	@brief default constructor; has no effect
		 */
		inline _TyPoint()
		{}

		/**
		 *	@brief constructor
		 *
		 *	@param[in] r_v_position is interpolated position
		 *	@param[in] r_t_rotation is interpolated rotation
		 */
		inline _TyPoint(const _TyPosition &r_v_position, const _TyRotation &r_t_rotation)
			:v_position(r_v_position), t_rotation(r_t_rotation)
		{}

		/**
		 *	@brief gets length of the positional part of the pose
		 *	@return Returns Euclidean length of the positional part of the pose.
		 *	@note This is required for calculation of spline lenght.
		 */
		inline _TyScalar f_Length() const
		{
			return v_position.f_Length();
		}
	};

	/**
	 *	@brief spline keypoint type
	 */
	struct _TyKeyPoint : public _TyPoint {
		//_TyPosition v_position; /**< @brief interpolated position */
		//_TyRotation t_rotation; /**< @brief interpolated rotation */
		_TyScalar f_time; /**< @brief explicit keypoint time */
		_TyScalar f_tension; /**< @brief tension */
		_TyScalar f_continuity; /**< @brief continuity */
		_TyScalar f_bias; /**< @brief bias */

		/**
		 *	@brief default constructor; has no effect
		 */
		inline _TyKeyPoint()
		{}

		/**
		 *	@brief constructor
		 *
		 *	@param[in] r_v_position is interpolated position
		 *	@param[in] r_t_rotation is interpolated rotation
		 *	@param[in] _f_time is keypoint time
		 *	@param[in] _f_tension is tension
		 *	@param[in] _f_continuity is continuity
		 *	@param[in] _f_bias is bias
		 */
		inline _TyKeyPoint(const _TyPosition &r_v_position, const _TyRotation &r_t_rotation,
			_TyScalar _f_time, _TyScalar _f_tension = 0, _TyScalar _f_continuity = 0,
			_TyScalar _f_bias = 0)
			:_TyPoint(r_v_position, r_t_rotation),/*v_position(r_v_position), t_rotation(r_t_rotation),*/
			f_time(_f_time), f_tension(_f_tension), f_continuity(_f_continuity), f_bias(_f_bias)
		{}
	};

protected:
	std::vector<_TyKeyPoint> m_point_list;

public:
	/**
	 *	@brief gets number of spline arcs
	 *	@return Returns number of spline arcs.
	 */
	inline size_t n_Arc_Num() const
	{
		return (m_point_list.size() > 3)? m_point_list.size() - 3 : 0;
	}

	/**
	 *	@brief calculates time length of a spline arc
	 *	@param[in] n_arc is zero-based index of spline arc
	 *	@return Returns time length of the selected arc.
	 */
	inline _TyScalar f_Arc_Time(size_t n_arc) const
	{
		_ASSERTE(n_arc < n_Arc_Num());
		const _TyKeyPoint &r_pt1 = m_point_list[n_arc + 1];
		const _TyKeyPoint &r_pt2 = m_point_list[n_arc + 2];
		return r_pt2.f_time - r_pt1.f_time;
	}

	/**
	 *	@brief calculates time length of this spline
	 *	@return Returns time length of the whole spline.
	 */
	inline _TyScalar f_Spline_Time() const
	{
		if(!n_Arc_Num())
			return 0;
		const _TyKeyPoint &r_pt1 = m_point_list[1];
		const _TyKeyPoint &r_pt2 = m_point_list[m_point_list.size() - 2];
		return r_pt2.f_time - r_pt1.f_time;
	}

	/**
	 *	@brief calculates position on a spline arc
	 *
	 *	@param[in] f_t is position on a spline arc in range [0, 1]
	 *	@param[in] n_arc is zero-based index of spline arc
	 *
	 *	@return Returns point at the given position along the selected arc.
	 */
	_TyPoint v_Arc_Position(_TyScalar f_t, size_t n_arc) const
	{
		_ASSERTE(n_arc < n_Arc_Num());

		const _TyKeyPoint &r_pt0 = m_point_list[n_arc];
		const _TyKeyPoint &r_pt1 = m_point_list[n_arc + 1];
		const _TyKeyPoint &r_pt2 = m_point_list[n_arc + 2];
		const _TyKeyPoint &r_pt3 = m_point_list[n_arc + 3];
		// get four points

		_TyPoint t_tangent0, t_tangent1;
		Get_Tangents(t_tangent0, t_tangent1, r_pt0, r_pt1, r_pt2, r_pt3);
		// calculate the tangents at r_pt1 and r_pt2

		//_ASSERTE(f_t >= 0 && f_t <= 1); // don't, may be slightly out of bounds
		_TyScalar f_t2 = f_t * f_t;
		_TyScalar f_t3 = f_t2 * f_t;
		// calculate t^2, t^3 (the times are already in [0, 1] range)

		_TyScalar f_weight0 = 2 * f_t3 - 3 * f_t2 + 1;	// [ 2 -3  0  1]
		_TyScalar f_weight1 = -2 * f_t3 + 3 * f_t2;		// [-2  3  0  0] 
		_TyScalar f_weight2 = f_t3 - 2 * f_t2 + f_t;	// [ 1 -2  1  0]
		_TyScalar f_weight3 = f_t3 - f_t2;				// [ 1 -1  0  0]
		// calculate polynomials

		_TyPosition v_position = r_pt1.v_position * f_weight0 + r_pt2.v_position * f_weight1 +
			t_tangent0.v_position * f_weight2 + t_tangent1.v_position * f_weight3;

		/*_TyRotation t_pq = r_pt1.t_rotation.t_Slerp(f_t, (r_pt1.t_rotation.f_Dot(r_pt2.t_rotation) < 0)?
			-r_pt2.t_rotation : r_pt2.t_rotation);
		_TyRotation t_ab = t_tangent0.t_rotation.t_Slerp(f_t, t_tangent1.t_rotation);
		_TyRotation t_rotation = t_pq.t_Slerp(2 * f_t * (1 - f_t), (t_pq.f_Dot(t_ab) < 0)? -t_ab : t_ab);*/
		/*_TyRotation t_pq = r_pt1.t_rotation.t_SlerpShort(f_t, r_pt2.t_rotation);
		_TyRotation t_ab = t_tangent0.t_rotation.t_Slerp(f_t, t_tangent1.t_rotation); // uses the "long" slerp
		_TyRotation t_rotation = t_pq.t_SlerpShort(2 * f_t * (1 - f_t), t_ab);*/
		// squad with quaternion sign fixup
		// not needed anymore, t_Slerp() and t_Squad() fixed now

		_TyRotation t_rotation = _TyRotation::t_Squad(f_t, r_pt1.t_rotation,
		     t_tangent0.t_rotation, t_tangent1.t_rotation, r_pt2.t_rotation);

		return _TyPoint(v_position, t_rotation);
		// interpolate
	}

	/**
	 *	@brief calculates tangent vector on a spline arc
	 *
	 *	@param[in] f_t is position on a spline arc in range [0, 1]
	 *	@param[in] n_arc is zero-based index of spline arc
	 *
	 *	@return Returns tangent at the given position along the selected arc.
	 *
	 *	@note The derivative of the rotation may not be an unit quaternion.
	 *	@todo The derivative of the rotation is untested.
	 */
	_TyPoint v_Arc_Derivative(_TyScalar f_t, size_t n_arc) const
	{
		_ASSERTE(n_arc < n_Arc_Num());

		const _TyKeyPoint &r_pt0 = m_point_list[n_arc];
		const _TyKeyPoint &r_pt1 = m_point_list[n_arc + 1];
		const _TyKeyPoint &r_pt2 = m_point_list[n_arc + 2];
		const _TyKeyPoint &r_pt3 = m_point_list[n_arc + 3];
		// get four points

		_TyPoint t_tangent0, t_tangent1;
		Get_Tangents(t_tangent0, t_tangent1, r_pt0, r_pt1, r_pt2, r_pt3);
		// calculate the tangents at r_pt1 and r_pt2

		//_ASSERTE(f_t >= 0 && f_t <= 1); // don't, may be slightly out of bounds
		_TyScalar f_t2 = f_t * f_t;
		//_TyScalar f_t3 = f_t2 * f_t; // unused
		// calculate t^2, t^3 (the times are already in [0, 1] range)

		_TyScalar f_weight0 = 6 * f_t2 - 6 * f_t;
		_TyScalar f_weight1 = -6 * f_t2 + 6 * f_t;
		_TyScalar f_weight2 = 3 * f_t2 - 4 * f_t + 1;
		_TyScalar f_weight3 = 3 * f_t2 - 2 * f_t;
		// calculate polynomials

		_TyPosition v_position = r_pt1.v_position * f_weight0 + r_pt2.v_position * f_weight1 +
			t_tangent0.v_position * f_weight2 + t_tangent1.v_position * f_weight3;

		_TyRotation r_pt2_t_rotation = (r_pt1.t_rotation.f_Dot(r_pt2.t_rotation) < 0)?
			-r_pt2.t_rotation : r_pt2.t_rotation;
		_TyRotation t_pq = r_pt1.t_rotation.t_Slerp(f_t, r_pt2_t_rotation);
		_TyRotation t_ab = t_tangent0.t_rotation.t_Slerp(f_t, t_tangent1.t_rotation);
		/*if(t_pq.f_Dot(t_ab) < 0)
			t_ab = -t_ab;*/ // unused and likely just plain wrong
		_TyRotation t_pq_prime = t_pq * (r_pt1.t_rotation.t_Inverse() * r_pt2_t_rotation).t_Log();
		_TyRotation t_ab_prime = t_ab * (t_tangent0.t_rotation.t_Inverse() * t_tangent1.t_rotation).t_Log();
		_TyRotation t_w = t_pq.t_Inverse() * t_ab;
		_TyRotation t_w_prime = t_pq.t_Inverse() * t_ab_prime - t_pq.t_Power(-2) * t_pq_prime * t_ab;
		_TyScalar f_t_lerp2 = 2 * f_t * (1 - f_t);
		_TyRotation t_rotation = t_pq * (t_w.t_Power(f_t_lerp2) * t_w.t_Log() * (2 - 4 * f_t) +
			t_w.t_Power(f_t_lerp2 - 1) * t_w_prime * f_t_lerp2) + t_pq_prime * t_w.t_Power(f_t_lerp2);
		// squad derivative with quaternion sign fixup

		return _TyPoint(v_position, t_rotation);
		// interpolate
	}

protected:
	static void Get_Tangents(_TyPoint &r_v_tangent0, _TyPoint &r_v_tangent1, const _TyKeyPoint &r_pt0,
		const _TyKeyPoint &r_pt1, const _TyKeyPoint &r_pt2, const _TyKeyPoint &r_pt3)
	{
		_TyPosition v_diff_32 = r_pt3.v_position - r_pt2.v_position,
			v_diff_21 = r_pt2.v_position - r_pt1.v_position,
			v_diff_10 = r_pt1.v_position - r_pt0.v_position;

		_TyRotation r_pt1_t_rotation = (r_pt0.t_rotation.f_Dot(r_pt1.t_rotation) < 0)?
			-r_pt1.t_rotation : r_pt1.t_rotation;
		_TyRotation r_pt2_t_rotation = (r_pt1_t_rotation.f_Dot(r_pt2.t_rotation) < 0)?
			-r_pt2.t_rotation : r_pt2.t_rotation;
		_TyRotation r_pt3_t_rotation = (r_pt2_t_rotation.f_Dot(r_pt3.t_rotation) < 0)?
			-r_pt3.t_rotation : r_pt3.t_rotation;
		// fixup quaternion signs locally before calculating the tangents

		_TyRotation t_diff_10 = r_pt0.t_rotation.t_Conjugate() * r_pt1_t_rotation,
			t_diff_21 = r_pt1_t_rotation.t_Conjugate() * r_pt2_t_rotation,
			t_diff_32 = r_pt2_t_rotation.t_Conjugate() * r_pt3_t_rotation;
		t_diff_32.Log();
		t_diff_21.Log();
		t_diff_10.Log();
		// get the differences

		{
			const _TyScalar t = r_pt1.f_tension, c = r_pt1.f_continuity, b = r_pt1.f_bias;
			_TyScalar f_dt_i1 = r_pt2.f_time - r_pt1.f_time, f_dt_i_plus_dt_i1 = r_pt2.f_time - r_pt0.f_time;
			_TyScalar f_t_corr = 2 * f_dt_i1 / f_dt_i_plus_dt_i1; // adjust for nonuniform time
			_TyScalar f_weight0 = f_t_corr * ((1 - t) * (1 + c) * (1 - b)) / 2;
			_TyScalar f_weight1 = f_t_corr * ((1 - t) * (1 - c) * (1 + b)) / 2;
			// calculate weights

			r_v_tangent0.v_position = v_diff_21 * f_weight0 + v_diff_10 * f_weight1;
			// interpolate position

			_TyRotation t_tangent = t_diff_21 * f_weight0 + t_diff_10 * f_weight1;
			r_v_tangent0.t_rotation = r_pt1_t_rotation * ((t_tangent - t_diff_21) * _TyScalar(.5)).t_Exp();
			// interpolate rotation (note that the quaternions need to be normalized)
		}
		{
			const _TyScalar t = r_pt2.f_tension, c = r_pt2.f_continuity, b = r_pt2.f_bias;
			_TyScalar f_dt_i = r_pt2.f_time - r_pt1.f_time, f_dt_i_plus_dt_i1 = r_pt3.f_time - r_pt1.f_time;
			_TyScalar f_t_corr = 2 * f_dt_i / f_dt_i_plus_dt_i1; // adjust for nonuniform time
			_TyScalar f_weight0 = f_t_corr * ((1 - t) * (1 - c) * (1 - b)) / 2;
			_TyScalar f_weight1 = f_t_corr * ((1 - t) * (1 + c) * (1 + b)) / 2;
			// calculate weights

			r_v_tangent1.v_position = v_diff_32 * f_weight0 + v_diff_21 * f_weight1;
			// interpolate position

			_TyRotation t_tangent = t_diff_32 * f_weight0 + t_diff_21 * f_weight1;
			r_v_tangent1.t_rotation = r_pt2_t_rotation * ((t_diff_21 - t_tangent) * _TyScalar(.5)).t_Exp();
			//if(r_v_tangent0.t_rotation.f_Dot(r_v_tangent1.t_rotation) < 0)
			//	r_v_tangent1.t_rotation = -r_v_tangent1.t_rotation; // fixup sign for t_Slerp() // this seem to be wrong: unable to make smooth rotations by angles close to pi
			// interpolate rotation (note that the quaternions need to be normalized)
		}
		// calculate the tangents at r_pt1 and r_pt2 (note that the quaternion parts are not
		// tangents, rather they are the segment end rotations multiplied by the tangents)
	}
};

/**
 *	@brief simple Kochanek - Bartels spline template, parametrized by point type
 */
template <class CPointClass>
class CKochanekBartelsSpline : public CSplineIface<CKochanekBartelsSplineImpl<CPointClass> > {
public:
	typedef CSplineIface<CKochanekBartelsSplineImpl<CPointClass> > _TyIface; /**< @brief base class */
	typedef typename _TyIface::_TyKeyPoint _TyKeyPoint; /**< @brief keypoint type */

public:
	/**
	 *	@copydoc CSplineBase::CSplineBase()
	 */
	inline CKochanekBartelsSpline()
		:_TyIface()
	{}

	/**
	 *	@copydoc CSplineBase::CSplineBase(const std::vector<_TyKeyPoint>&)
	 */
	inline CKochanekBartelsSpline(const std::vector<_TyKeyPoint> &r_point_list) // throw(std::bad_alloc)
		:_TyIface(r_point_list)
	{}

	/**
	 *	@copydoc CSplineBase::CSplineBase(CPointIter,CPointIter)
	 */
	template <class CPointIter>
	inline CKochanekBartelsSpline(CPointIter p_point_begin, CPointIter p_point_end) // throw(std::bad_alloc)
		:_TyIface(p_point_begin, p_point_end)
	{}

	/**
	 *	@copydoc CSplineBase::CSplineBase(const CSplineBase&)
	 */
	inline CKochanekBartelsSpline(const CKochanekBartelsSpline &r_other) // throw(std::bad_alloc)
		:_TyIface(r_other)
	{}

	/**
	 *	@copydoc CSplineBase::operator=()
	 */
	inline CKochanekBartelsSpline &operator =(const CKochanekBartelsSpline &r_other) // throw(std::bad_alloc)
	{
		*(_TyIface*)this = r_other;
		return *this;
	}
};

/**
 *	@brief simple Kochanek - Bartels spline template, parametrized by point type
 *	@note Integration of the spline derivative is costy (the derivative function also
 *		calculates derivative of the quaternion, which is not used).
 */
template <class CPointClass>
class CKochanekBartelsQuatSpline : public CSplineIface<CKochanekBartelsQuatSplineImpl<CPointClass> > {
public:
	typedef CSplineIface<CKochanekBartelsQuatSplineImpl<CPointClass> > _TyIface; /**< @brief base class */
	typedef typename _TyIface::_TyKeyPoint _TyKeyPoint; /**< @brief keypoint type */

public:
	/**
	 *	@copydoc CSplineBase::CSplineBase()
	 */
	inline CKochanekBartelsQuatSpline()
		:_TyIface()
	{}

	/**
	 *	@copydoc CSplineBase::CSplineBase(const std::vector<_TyKeyPoint>&)
	 */
	inline CKochanekBartelsQuatSpline(const std::vector<_TyKeyPoint> &r_point_list) // throw(std::bad_alloc)
		:_TyIface(r_point_list)
	{}

	/**
	 *	@copydoc CSplineBase::CSplineBase(CPointIter,CPointIter)
	 */
	template <class CPointIter>
	inline CKochanekBartelsQuatSpline(CPointIter p_point_begin, CPointIter p_point_end) // throw(std::bad_alloc)
		:_TyIface(p_point_begin, p_point_end)
	{}

	/**
	 *	@copydoc CSplineBase::CSplineBase(const CSplineBase&)
	 */
	inline CKochanekBartelsQuatSpline(const CKochanekBartelsQuatSpline &r_other) // throw(std::bad_alloc)
		:_TyIface(r_other)
	{}

	/**
	 *	@copydoc CSplineBase::operator=()
	 */
	inline CKochanekBartelsQuatSpline &operator =(const CKochanekBartelsQuatSpline &r_other) // throw(std::bad_alloc)
	{
		*(_TyIface*)this = r_other;
		return *this;
	}
};

#endif // !__KOCHANEK_BARTELS_SPLINE_INCLUDED
