/******************************************************************************
* Cagd_aux.c - auxiliary routine to interface to different free from types.   *
*******************************************************************************
* Written by Gershon Elber, July. 90.					      *
******************************************************************************/

#include "cagd_loc.h"

#define VEC_FIELD_TRIES	10
#define VEC_FIELD_START_STEP 1e-6


/*****************************************************************************
* DESCRIPTION:                                                               M
*   Returns the parametric domain of a curve.	     		             M
*                                                                            *
* PARAMETERS:                                                                M
*   Crv:       To get its parametric domain.                                 M
*   TMin:      Where to put the minimal domain's boundary.                   M
*   TMax:      Where to put the maximal domain's boundary.                   M
*                                                                            *
* RETURN VALUE:                                                              M
*   void                                                                     M
*                                                                            *
* KEYWORDS:                                                                  M
*   CagdCrvDomain, domain, parametric domain                                 M
*****************************************************************************/
void CagdCrvDomain(CagdCrvStruct *Crv, CagdRType *TMin, CagdRType *TMax)
{
    switch (Crv -> GType) {
	case CAGD_CBEZIER_TYPE:
	    *TMin = 0.0;
	    *TMax = 1.0;
	    break;
	case CAGD_CBSPLINE_TYPE:
	    BspCrvDomain(Crv, TMin, TMax);
	    break;
	case CAGD_CPOWER_TYPE:
	default:
	    CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_CRV);
	    break;
    }
}

/*****************************************************************************
* DESCRIPTION:                                                               M
* Given a curve and parameter value t, evaluate the curve at t.		     M
*                                                                            *
* PARAMETERS:                                                                M
*   Crv:      To evaluate at the given parametric location t.                M
*   t:        The parameter value at which the curve Crv is to be evaluated. M
*                                                                            *
* RETURN VALUE:                                                              M
*   CagdRType *:  A vector holding all the coefficients of all components    M
*                 of curve Crv's point type. If for example the curve's      M
*                 point type is P2, the W, X, and Y will be saved in the     M
*                 first three locations of the returned vector. The first    M
*                 location (index 0) of the returned vector is reserved for  M
*                 the rational coefficient W and XYZ always starts at second M
*                 location of the returned vector (index 1).                 M
*                                                                            *
* KEYWORDS:                                                                  M
*   CagdCrvEval, evaluation                                                  M
*****************************************************************************/
CagdRType *CagdCrvEval(CagdCrvStruct *Crv, CagdRType t)
{
    switch (Crv -> GType) {
	case CAGD_CBEZIER_TYPE:
	    return BzrCrvEvalAtParam(Crv, t);
	case CAGD_CBSPLINE_TYPE:
	    return BspCrvEvalAtParam(Crv, t);
	case CAGD_CPOWER_TYPE:
	    CAGD_FATAL_ERROR(CAGD_ERR_POWER_NO_SUPPORT);
	    return NULL;
	default:
	    CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_CRV);
	    return NULL;
    }
}

/*****************************************************************************
* DESCRIPTION:                                                               M
*   Returns the parametric domain of a surface.				     M
*                                                                            *
* PARAMETERS:                                                                M
*   Srf:       To get its parametric domain.                                 M
*   UMin:      Where to put the minimal U domain's boundary.                 M
*   UMax:      Where to put the maximal U domain's boundary.                 M
*   VMin:      Where to put the minimal V domain's boundary.                 M
*   VMax:      Where to put the maximal V domain's boundary.                 M
*                                                                            *
* RETURN VALUE:                                                              M
*   void                                                                     M
*                                                                            *
* KEYWORDS:                                                                  M
*   CagdSrfDomain, domain, parametric domain                                 M
*****************************************************************************/
void CagdSrfDomain(CagdSrfStruct *Srf,
		   CagdRType *UMin,
		   CagdRType *UMax,
		   CagdRType *VMin,
		   CagdRType *VMax)
{
    switch (Srf -> GType) {
	case CAGD_SBEZIER_TYPE:
	    *UMin = 0.0;
	    *UMax = 1.0;
	    *VMin = 0.0;
	    *VMax = 1.0;
	    break;
	case CAGD_SBSPLINE_TYPE:
	    BspSrfDomain(Srf, UMin, UMax, VMin, VMax);
	    break;
	case CAGD_SPOWER_TYPE:
	default:
	    CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_CRV);
	    break;
    }
}

/*****************************************************************************
* DESCRIPTION:                                                               M
* Given a surface and parameter values u, v, evaluate the surface at (u, v). M
*                                                                            *
* PARAMETERS:                                                                M
*   Srf:     To evaluate at the given parametric location (u, v).            M
*   u, v:    The parameter values at which the curve Crv is to be evaluated. M
*                                                                            *
* RETURN VALUE:                                                              M
*   CagdRType *:  A vector holding all the coefficients of all components    M
*                 of surface Srf's point type. If for example the surface's  M
*                 point type is P2, the W, X, and Y will be saved in the     M
*                 first three locations of the returned vector. The first    M
*                 location (index 0) of the returned vector is reserved for  M
*                 the rational coefficient W and XYZ always starts at second M
*                 location of the returned vector (index 1).                 M
*                                                                            *
* KEYWORDS:                                                                  M
*   CagdSrfEval, evaluation                                                  M
*****************************************************************************/
CagdRType *CagdSrfEval(CagdSrfStruct *Srf, CagdRType u, CagdRType v)
{
    switch (Srf -> GType) {
	case CAGD_SBEZIER_TYPE:
	    return BzrSrfEvalAtParam(Srf, u, v);
	case CAGD_SBSPLINE_TYPE:
	    return BspSrfEvalAtParam(Srf, u, v);
	case CAGD_SPOWER_TYPE:
	    CAGD_FATAL_ERROR(CAGD_ERR_POWER_NO_SUPPORT);
	    return NULL;
	default:
	    CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_SRF);
	    return NULL;
    }
}

/*****************************************************************************
* DESCRIPTION:                                                               M
* Evaluates a vector field surface to a unit size vector. If fails, moves a  M
* tad until success. Useful for normal field evaluations.		     M
*                                                                            *
* PARAMETERS:                                                                M
*   Vec:            Where resulting unit length vector is to be saved.       M
*   VecFieldSrf:    A surface representing a vector field.                   M
*   U, V:           Parameter locations.                                     M
*                                                                            *
* RETURN VALUE:                                                              M
*   void                                                                     M
*                                                                            *
* KEYWORDS:                                                                  M
*   CagdEvaluateSurfaceVecField, normal, vector field                        M
*****************************************************************************/
void CagdEvaluateSurfaceVecField(CagdVType Vec,
				 CagdSrfStruct *VecFieldSrf,
				 CagdRType U,
				 CagdRType V)
{
    CagdRType
	*R = CagdSrfEval(VecFieldSrf, U, V);

    CagdCoerceToE3(Vec, &R, -1, VecFieldSrf -> PType);

    if (PT_LENGTH(Vec) < IRIT_EPSILON) {
	int i = 0;
	CagdRType UMin, UMax, VMin, VMax, UMid, VMid,
	    Step = VEC_FIELD_START_STEP;

	CagdSrfDomain(VecFieldSrf, &UMin, &UMax, &VMin, &VMax);
	UMid = (UMin + UMax) / 2;
	VMid = (VMin + VMax) / 2;
	while (PT_LENGTH(Vec) < IRIT_EPSILON && i++ < VEC_FIELD_TRIES) {
	    U += U < UMid ? Step : -Step;
	    V += V < VMid ? Step : -Step;
	    Step *= 2.0;

	    R = CagdSrfEval(VecFieldSrf, U, V);
	    CagdCoerceToE3(Vec, &R, -1, VecFieldSrf -> PType);
	}
	if (i >= VEC_FIELD_TRIES)
	    CAGD_FATAL_ERROR(CAGD_ERR_CANNOT_COMP_VEC_FIELD);
    }

    PT_NORMALIZE(Vec);
}

/*****************************************************************************
* DESCRIPTION:                                                               M
* Given a curve, computes its derivative curve (Hodograph).                  M
*                                                                            *
* PARAMETERS:                                                                M
*   Crv:      To compute its Hodograph curve.                                M
*                                                                            *
* RETURN VALUE:                                                              M
*   CagdCrvStruct *:  Resulting hodograph.                                   M
*                                                                            *
* KEYWORDS:                                                                  M
*   CagdCrvDerive, derivatives, Hodograph                                    M
*****************************************************************************/
CagdCrvStruct *CagdCrvDerive(CagdCrvStruct *Crv)
{
    switch (Crv -> GType) {
	case CAGD_CBEZIER_TYPE:
	    return BzrCrvDerive(Crv);
	case CAGD_CBSPLINE_TYPE:
	    return BspCrvDerive(Crv);
	case CAGD_CPOWER_TYPE:
	    CAGD_FATAL_ERROR(CAGD_ERR_POWER_NO_SUPPORT);
	    return NULL;
	default:
	    CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_CRV);
	    return NULL;
    }
}

/*****************************************************************************
* DESCRIPTION:                                                               M
* Given a curve, compute its integral curve.                                 M
*                                                                            *
* PARAMETERS:                                                                M
*   Crv:      To compute its integral curve.                                 M
*                                                                            *
* RETURN VALUE:                                                              M
*   CagdCrvStruct *:  Resulting integral curve.                              M
*                                                                            *
* KEYWORDS:                                                                  M
*   CagdCrvIntegrate, integrals                                              M
*****************************************************************************/
CagdCrvStruct *CagdCrvIntegrate(CagdCrvStruct *Crv)
{
    switch (Crv -> GType) {
	case CAGD_CBEZIER_TYPE:
	    return BzrCrvIntegrate(Crv);
	case CAGD_CBSPLINE_TYPE:
	    return BspCrvIntegrate(Crv);
	case CAGD_CPOWER_TYPE:
	    CAGD_FATAL_ERROR(CAGD_ERR_POWER_NO_SUPPORT);
	    return NULL;
	default:
	    CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_CRV);
	    return NULL;
    }
}

/*****************************************************************************
* DESCRIPTION:                                                               M
* Given a surface, computes its partial derivative in the prescibed          M
* direction Dir.							     M
*                                                                            *
* PARAMETERS:                                                                M
*   Srf:      To compute its derivative surface in direction Dir.            M
*   Dir:      Direction of differentiation. Either U or V.                   M
*                                                                            *
* RETURN VALUE:                                                              M
*   CagdSrfStruct *:  Resulting partial derivative surface.                  M
*                                                                            *
* KEYWORDS:                                                                  M
*   CagdSrfDerive, derivatives, partial derivatives                          M
*****************************************************************************/
CagdSrfStruct *CagdSrfDerive(CagdSrfStruct *Srf, CagdSrfDirType Dir)
{
    switch (Srf -> GType) {
	case CAGD_SBEZIER_TYPE:
	    return BzrSrfDerive(Srf, Dir);
	case CAGD_SBSPLINE_TYPE:
	    return BspSrfDerive(Srf, Dir);
	case CAGD_SPOWER_TYPE:
	    CAGD_FATAL_ERROR(CAGD_ERR_POWER_NO_SUPPORT);
	    return NULL;
	default:
	    CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_SRF);
	    return NULL;
    }
}

/*****************************************************************************
* DESCRIPTION:                                                               M
* Given a curve - subdivides it into two curves at the given parameter       M
* value t.								     M
*    Returns pointer to first curve in a list of two subdivided curves.      M
*                                                                            *
* PARAMETERS:                                                                M
*   Crv:      To subdivide at the prescibed parameter value t.               M
*   t:        The parameter to subdivide the curve Crv at.                   M
*                                                                            *
* RETURN VALUE:                                                              M
*   CagdCrvStruct *:  A list of the two curves resulting from the process    M
*                     of subdivision.                                        M
*                                                                            *
* KEYWORDS:                                                                  M
*   CagdCrvSubdivAtParam, subdivision                                        M
*****************************************************************************/
CagdCrvStruct *CagdCrvSubdivAtParam(CagdCrvStruct *Crv, CagdRType t)
{
    switch (Crv -> GType) {
	case CAGD_CBEZIER_TYPE:
	    return BzrCrvSubdivAtParam(Crv, t);
	case CAGD_CBSPLINE_TYPE:
	    return BspCrvSubdivAtParam(Crv, t);
	case CAGD_CPOWER_TYPE:
	    CAGD_FATAL_ERROR(CAGD_ERR_POWER_NO_SUPPORT);
	    return NULL;
	default:
	    CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_CRV);
	    return NULL;
    }
}

/*****************************************************************************
* DESCRIPTION:                                                               M
* Given a curve - extracts a sub-region within the domain specified by t1    M
* and t2.                                                                    M
*                                                                            *
* PARAMETERS:                                                                M
*   Crv:       To extract a sub-region from.                                 M
*   t1, t2:    Parametric domain boundaries of sub-region.                   M
*                                                                            *
* RETURN VALUE:                                                              M
*   CagdCrvStruct *:  Sub-region extracted from Crv from t1 to t2.           M
*                                                                            *
* KEYWORDS:                                                                  M
*   CagdCrvRegionFromCrv, regions, subdivision                               M
*****************************************************************************/
CagdCrvStruct *CagdCrvRegionFromCrv(CagdCrvStruct *Crv,
				    CagdRType t1,
				    CagdRType t2)
{
    CagdRType TMin, TMax;
    CagdCrvStruct *Crvs;
    CagdBType
	BezCrv = FALSE,
	OpenEnd = TRUE,
	NewCrv = FALSE;

    switch (Crv -> GType) {
	case CAGD_CBEZIER_TYPE:
	    BezCrv = TRUE;
	    break;
	case CAGD_CBSPLINE_TYPE:
	    OpenEnd = BspCrvHasOpenEC(Crv);
	    break;
	case CAGD_CPOWER_TYPE:
	    CAGD_FATAL_ERROR(CAGD_ERR_POWER_NO_SUPPORT);
	    return NULL;
	default:
	    CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_CRV);
	    return NULL;
    }

    CagdCrvDomain(Crv, &TMin, &TMax);
    CAGD_DOMAIN_T_VERIFY(t1, TMin, TMax);
    CAGD_DOMAIN_T_VERIFY(t2, TMin, TMax);

    if (t1 > t2)
	SWAP(CagdRType, t1, t2);

    if (!APX_EQ(t1, TMin) || !OpenEnd) {
    	Crvs = CagdCrvSubdivAtParam(Crv, t1);
	Crv = Crvs -> Pnext;
	Crvs -> Pnext = NULL;
	CagdCrvFree(Crvs);			   /* Free the first region. */
	NewCrv = TRUE;
    }

    if (APX_EQ(t2, TMax) && OpenEnd)
	return NewCrv ? Crv : CagdCrvCopy(Crv);
    else {
	if (BezCrv)
	    t2 = (t2 - t1) / (TMax - t1);

	Crvs = CagdCrvSubdivAtParam(Crv, t2);

	if (NewCrv)
	    CagdCrvFree(Crv);

    	CagdCrvFree(Crvs -> Pnext);		  /* Free the second region. */
    	Crvs -> Pnext = NULL;
    	return Crvs;				/* Returns the first region. */
    }
}

/*****************************************************************************
* DESCRIPTION:                                                               M
* Given a curve - refines it at the given n knots as defined by vector t.    M
* If Replace is TRUE, the values in t replaces current knot vector.	     M
* Returns pointer to refined surface (Note a Bezier curve will be converted  M
* into a Bspline curve).                                                     M
*                                                                            *
* PARAMETERS:                                                                M
*   Crv:       To refine.                                                    M
*   Replace:   If TRUE, t holds knots in exactly the same length as the      M
*              length of the knot vector of Crv and t simply replaces the    M
*              knot vector.                                                  M
*   t:         Vector of knots with length of n.                             M
*   n:         Length of vector t.                                           M
*                                                                            *
* RETURN VALUE:                                                              M
*   CagdCrvStruct *:  A refined curve of Crv after insertion of all the      M
*                     knots as specified by vector t of length n.            M
*                                                                            *
* KEYWORDS:                                                                  M
*   CagdCrvRefineAtParams, refinement, subdivision                           M
*****************************************************************************/
CagdCrvStruct *CagdCrvRefineAtParams(CagdCrvStruct *Crv,
				     CagdBType Replace,
				     CagdRType *t,
				     int n)
{
    CagdCrvStruct *BspCrv, *TCrv;

    switch (Crv -> GType) {
	case CAGD_CBEZIER_TYPE:
    	    BspCrv = CnvrtBezier2BsplineCrv(Crv);
	    TCrv = BspCrvKnotInsertNDiff(BspCrv, Replace, t, n);
	    CagdCrvFree(BspCrv);
	    return TCrv;
	case CAGD_CBSPLINE_TYPE:
	    return BspCrvKnotInsertNDiff(Crv, Replace, t, n);
	case CAGD_CPOWER_TYPE:
	    CAGD_FATAL_ERROR(CAGD_ERR_POWER_NO_SUPPORT);
	    return NULL;
	default:
	    CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_CRV);
	    return NULL;
    }
}

/*****************************************************************************
* DESCRIPTION:                                                               M
* Returns a new curve that is the reversed curve of Crv by reversing the     M
* control polygon and the knot vector of Crv is a Bspline curve.             M
* See also BspKnotReverse.                                                   M
*                                                                            *
* PARAMETERS:                                                                M
*   Crv:       To be reversed.                                               M
*                                                                            *
* RETURN VALUE:                                                              M
*   CagdCrvStruct *:   Reversed curve of Crv.                                M
*                                                                            *
* KEYWORDS:                                                                  M
*   CagdCrvReverse, reverse                                                  M
*****************************************************************************/
CagdCrvStruct *CagdCrvReverse(CagdCrvStruct *Crv)
{
    CagdBType
	IsNotRational = !CAGD_IS_RATIONAL_CRV(Crv);
    int i, Len, Col,
	Length = Crv -> Length,
	MaxCoord = CAGD_NUM_OF_PT_COORD(Crv -> PType);
    CagdCrvStruct
        *ReversedCrv = CagdCrvCopy(Crv);
    CagdRType *KV,
	**Points = ReversedCrv -> Points;

    switch (Crv -> GType) {
	case CAGD_CBEZIER_TYPE:
	case CAGD_CBSPLINE_TYPE:
	    break;
	case CAGD_CPOWER_TYPE:
	    CAGD_FATAL_ERROR(CAGD_ERR_POWER_NO_SUPPORT);
	    return NULL;
	default:
	    CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_CRV);
	    return NULL;
    }

    /* Reverse the Ctl Polygon: */
    Len = Length / 2;
    for (Col = 0; Col < Len; Col++)
	for (i = IsNotRational; i <= MaxCoord; i++)
	    SWAP(CagdRType,
		 Points[i][Col],
		 Points[i][Length - Col - 1]);

    /* Reverse the knot vector if it exists: */
    if (Crv -> GType == CAGD_CBSPLINE_TYPE &&
	Crv -> KnotVector != NULL) {
	KV = BspKnotReverse(Crv -> KnotVector, Crv -> Order + Length);
	IritFree((VoidPtr) ReversedCrv -> KnotVector);
	ReversedCrv -> KnotVector = KV;
    }

    return ReversedCrv;
}

/*****************************************************************************
* DESCRIPTION:                                                               M
* Returns a new curve representing the same curve as Crv but with its degree M
* raised by one.							     M
*                                                                            *
* PARAMETERS:                                                                M
*   Crv:       To raise its degree.                                          M
*                                                                            *
* RETURN VALUE:                                                              M
*   CagdCrvStruct *:  A curve with same geometry as Crv but with one degree  M
*                     higher.                                                M
*                                                                            *
* KEYWORDS:                                                                  M
*   CagdCrvDegreeRaise, degree raising                                       M
*****************************************************************************/
CagdCrvStruct *CagdCrvDegreeRaise(CagdCrvStruct *Crv)
{
    switch (Crv -> GType) {
	case CAGD_CBEZIER_TYPE:
	    return BzrCrvDegreeRaise(Crv);
	case CAGD_CBSPLINE_TYPE:
	    return BspCrvDegreeRaise(Crv);
	case CAGD_CPOWER_TYPE:
	    CAGD_FATAL_ERROR(CAGD_ERR_POWER_NO_SUPPORT);
	    return NULL;
	default:
	    CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_CRV);
	    return NULL;
    }
}

/*****************************************************************************
* DESCRIPTION:                                                               M
* Returns a new curve representing the same curve as Crv but with its degree M
* raised to NewOrder							     M
*                                                                            *
* PARAMETERS:                                                                M
*   Crv:        To raise its degree.                                         M
*   NewOrder:   Expected new order of the raised curve.                      M
*                                                                            *
* RETURN VALUE:                                                              M
*   CagdCrvStruct *:  A curve with same geometry as Crv but with order that  M
*                     is equal to NewOrder.                                  M
*                                                                            *
* KEYWORDS:                                                                  M
*   CagdCrvDegreeRaiseN, degree raising                                      M
*****************************************************************************/
CagdCrvStruct *CagdCrvDegreeRaiseN(CagdCrvStruct *Crv, int NewOrder)
{
    switch (Crv -> GType) {
	case CAGD_CBEZIER_TYPE:
	    return BzrCrvDegreeRaiseN(Crv, NewOrder);
	case CAGD_CBSPLINE_TYPE:
	    return BspCrvDegreeRaiseN(Crv, NewOrder);
	case CAGD_CPOWER_TYPE:
	    CAGD_FATAL_ERROR(CAGD_ERR_POWER_NO_SUPPORT);
	    return NULL;
	default:
	    CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_CRV);
	    return NULL;
    }
}

/*****************************************************************************
* DESCRIPTION:                                                               M
* Returns a new surface representing the same surface as Srf but with its    M
* degree raised by one.							     M
*                                                                            *
* PARAMETERS:                                                                M
*   Srf:       To raise its degree.                                          M
*   Dir:       Direction of degree raising. Either U or V.		     M
*                                                                            *
* RETURN VALUE:                                                              M
*   CagdSrfStruct *:  A surface with same geometry as Srf but with one       M
*                     degree higher.                                         M
*                                                                            *
* KEYWORDS:                                                                  M
*   CagdSrfDegreeRaise, degree raising                                       M
*****************************************************************************/
CagdSrfStruct *CagdSrfDegreeRaise(CagdSrfStruct *Srf, CagdSrfDirType Dir)
{
    switch (Srf -> GType) {
	case CAGD_SBEZIER_TYPE:
	    return BzrSrfDegreeRaise(Srf, Dir);
	case CAGD_SBSPLINE_TYPE:
	    return BspSrfDegreeRaise(Srf, Dir);
	case CAGD_SPOWER_TYPE:
	    CAGD_FATAL_ERROR(CAGD_ERR_POWER_NO_SUPPORT);
	    return NULL;
	default:
	    CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_SRF);
	    return NULL;
    }
}

/*****************************************************************************
* DESCRIPTION:                                                               M
* Extracts an isoparametric curve from the surface Srf in direction Dir at   M
* the parameter value of t.                                                  M
*                                                                            *
* PARAMETERS:                                                                M
*   Srf:       To extract an isoparametric curve from.                       M
*   t:         Parameter value of extracted isoparametric curve.             M
*   Dir:       Direction of extracted isoparametric curve. Either U or V.    M
*                                                                            *
* RETURN VALUE:                                                              M
*   CagdCrvStruct *:   An isoparametric curve of Srf. This curve inherit the M
*                      order and continuity of surface Srf in direction Dir. M
*                                                                            *
* KEYWORDS:                                                                  M
*   CagdCrvFromSrf, isoparametric curves, curve from surface                 M
*****************************************************************************/
CagdCrvStruct *CagdCrvFromSrf(CagdSrfStruct *Srf,
			      CagdRType t,
			      CagdSrfDirType Dir)
{
    switch (Srf -> GType) {
	case CAGD_SBEZIER_TYPE:
	    return BzrSrfCrvFromSrf(Srf, t, Dir);
	case CAGD_SBSPLINE_TYPE:
	    return BspSrfCrvFromSrf(Srf, t, Dir);
	case CAGD_SPOWER_TYPE:
	    CAGD_FATAL_ERROR(CAGD_ERR_POWER_NO_SUPPORT);
	    return NULL;
	default:
	    CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_SRF);
	    return NULL;
    }
}

/*****************************************************************************
* DESCRIPTION:                                                               M
* Extracts a curve from the mesh of surface Srf in direction Dir at index    M
* Index.                                                                     M
*                                                                            *
* PARAMETERS:                                                                M
*   Srf:       To extract a curve from.  		                     M
*   Index:     Index along the mesh of Srf to extract the curve from.        M
*   Dir:       Direction of extracted curve. Either U or V.		     M
*                                                                            *
* RETURN VALUE:                                                              M
*   CagdCrvStruct *:   A curve from Srf. This curve inherit the  order and   M
*                      continuity of surface Srf in direction Dir. However,  M
*                      thiscurve is not on surface Srf, in general.          M
*                                                                            *
* KEYWORDS:                                                                  M
*   CagdCrvFromMesh, isoparametric curves, curve from mesh                   M
*****************************************************************************/
CagdCrvStruct *CagdCrvFromMesh(CagdSrfStruct *Srf,
			       int Index,
			       CagdSrfDirType Dir)
{
    switch (Srf -> GType) {
	case CAGD_SBEZIER_TYPE:
	    return BzrSrfCrvFromMesh(Srf, Index, Dir);
	case CAGD_SBSPLINE_TYPE:
	    return BspSrfCrvFromMesh(Srf, Index, Dir);
	case CAGD_SPOWER_TYPE:
	    CAGD_FATAL_ERROR(CAGD_ERR_POWER_NO_SUPPORT);
	    return NULL;
	default:
	    CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_SRF);
	    return NULL;
    }
}

/*****************************************************************************
* DESCRIPTION:                                                               M
* Substitutes a row/column of surface Srf from the given curve Crv at        M
* surface direction Dir and mesh index Index. Curve must have the same       M
* PtType/Length as the surface in the selected direction.		     M
*                                                                            *
* PARAMETERS:                                                                M
*   Crv:       To substitute into the surface Srf.                           M
*   Index:     Of mesh where the curve Crv should be substituted in.         M
*   Dir:       Either U or V.                                                M
*   Srf:       That a row or a column of should be replaced by Crv.          M
*                                                                            *
* RETURN VALUE:                                                              M
*   void                                                                     M
*                                                                            *
* KEYWORDS:                                                                  M
*   CagdCrvToMesh, curve from mesh                                           M
*****************************************************************************/
void CagdCrvToMesh(CagdCrvStruct *Crv,
		   int Index,
		   CagdSrfDirType Dir,
		   CagdSrfStruct *Srf)
{
    CagdBType
	IsNotRational = !CAGD_IS_RATIONAL_SRF(Srf);
    int i, j,
	Length = Crv -> Length,
	ULength = Srf -> ULength,
	VLength = Srf -> VLength,
	MaxCoord = CAGD_NUM_OF_PT_COORD(Srf -> PType);
    CagdRType *CrvP, *SrfP;

    if (Crv -> PType != Srf -> PType ||
	Length != (Dir == CAGD_CONST_U_DIR ? VLength : ULength))
	CAGD_FATAL_ERROR(CAGD_ERR_PT_OR_LEN_MISMATCH);

    switch (Dir) {
	case CAGD_CONST_U_DIR:
	    if (Index + 1 > ULength)
		CAGD_FATAL_ERROR(CAGD_ERR_INDEX_NOT_IN_MESH);

	    for (i = IsNotRational; i <= MaxCoord; i++) {
		CrvP = Crv -> Points[i];
		SrfP = Srf -> Points[i] + Index * CAGD_NEXT_U(Srf);
		for (j = 0; j < Length; j++) {
		    *SrfP = *CrvP++;
		    SrfP += CAGD_NEXT_V(Srf);
		}
	    }
	    break;
	case CAGD_CONST_V_DIR:
	    if (Index + 1 > VLength)
		CAGD_FATAL_ERROR(CAGD_ERR_INDEX_NOT_IN_MESH);

	    for (i = IsNotRational; i <= MaxCoord; i++) {
		CrvP = Crv -> Points[i];
		SrfP = Srf -> Points[i] + Index * CAGD_NEXT_V(Srf);
		for (j = 0; j < Length; j++) {
		    *SrfP = *CrvP++;
		    SrfP += CAGD_NEXT_U(Srf);
		}
	    }
	    break;
	default:
	    CAGD_FATAL_ERROR(CAGD_ERR_DIR_NOT_CONST_UV);
	    break;
    }
}

/*****************************************************************************
* DESCRIPTION:                                                               M
* Given a surface - subdivides it into two sub-surfaces at given parametric  M
* value t in the given direction Dir.                                        M
*   Returns pointer to first surface in a list of two subdivided surfaces.   M
*                                                                            *
* PARAMETERS:                                                                M
*   Srf:      To subdivide at the prescibed parameter value t.               M
*   t:        The parameter to subdivide the curve Crv at.                   M
*   Dir:      Direction of subdivision. Either U or V.                       M
*                                                                            *
* RETURN VALUE:                                                              M
*   CagdSrfStruct *:  A list of the two surfaces resulting from the process  M
*                     of subdivision.                                        M
*                                                                            *
* KEYWORDS:                                                                  M
*   CagdSrfSubdivAtParam, subdivision                                        M
*****************************************************************************/
CagdSrfStruct *CagdSrfSubdivAtParam(CagdSrfStruct *Srf,
				    CagdRType t,
				    CagdSrfDirType Dir)
{
    switch (Srf -> GType) {
	case CAGD_SBEZIER_TYPE:
	    return BzrSrfSubdivAtParam(Srf, t, Dir);
	case CAGD_SBSPLINE_TYPE:
	    return BspSrfSubdivAtParam(Srf, t, Dir);
	case CAGD_SPOWER_TYPE:
	    CAGD_FATAL_ERROR(CAGD_ERR_POWER_NO_SUPPORT);
	    return NULL;
	default:
	    CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_SRF);
	    return NULL;
    }
}
/*****************************************************************************
* DESCRIPTION:                                                               M
* Given a surface - extracts a sub-region within the domain specified by t1  M
* and t2, in the direction Dir.                                              M
*                                                                            *
* PARAMETERS:                                                                M
*   Srf:       To extract a sub-region from.                                 M
*   t1, t2:    Parametric domain boundaries of sub-region.                   M
*   Dir:       Direction of region extraction. Either U or V.                M
*                                                                            *
* RETURN VALUE:                                                              M
*   CagdSrfStruct *:  Sub-region extracted from Srf from t1 to t2.           M
*                                                                            *
* KEYWORDS:                                                                  M
*   CagdSrfRegionFromSrf, regions, subdivision                               M
*****************************************************************************/
CagdSrfStruct *CagdSrfRegionFromSrf(CagdSrfStruct *Srf,
				    CagdRType t1,
				    CagdRType t2,
				    CagdSrfDirType Dir)
{
    CagdRType TMin, TMax, R1, R2;
    CagdSrfStruct *Srfs;
    CagdBType
	OpenEnd = FALSE,
	NewSrf = FALSE;

    switch (Srf -> GType) {
	case CAGD_SBEZIER_TYPE:
	    break;
	case CAGD_SBSPLINE_TYPE:
	    OpenEnd = BspSrfHasOpenECDir(Srf, Dir);
	    break;
	case CAGD_SPOWER_TYPE:
	    CAGD_FATAL_ERROR(CAGD_ERR_POWER_NO_SUPPORT);
	    return NULL;
	default:
	    CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_SRF);
	    return NULL;
    }

    if (Dir == CAGD_CONST_U_DIR)
	CagdSrfDomain(Srf, &TMin, &TMax, &R1, &R2);
    else
	CagdSrfDomain(Srf, &R1, &R2, &TMin, &TMax);
    CAGD_DOMAIN_T_VERIFY(t1, TMin, TMax);
    CAGD_DOMAIN_T_VERIFY(t2, TMin, TMax);

    if (t1 > t2)
	SWAP(CagdRType, t1, t2);

    if (!APX_EQ(t1, TMin) || !OpenEnd) {
	Srfs = CagdSrfSubdivAtParam(Srf, t1, Dir);
	Srf = Srfs -> Pnext;
	Srfs -> Pnext = NULL;
	CagdSrfFree(Srfs);			   /* Free the first region. */
	NewSrf = TRUE;
    }

    if (APX_EQ(t2, TMax) && OpenEnd)
	return NewSrf ? Srf : CagdSrfCopy(Srf);
    else {
	Srfs = CagdSrfSubdivAtParam(Srf, t2, Dir);

	if (NewSrf)
	    CagdSrfFree(Srf);

    	CagdSrfFree(Srfs -> Pnext);		  /* Free the second region. */
    	Srfs -> Pnext = NULL;
	return Srfs;				/* Returns the first region. */
    }
}
/*****************************************************************************
* DESCRIPTION:                                                               M
* Given a surface - refines it at the given n knots as defined by vector t.  M
*   If Replace is TRUE, the values in t replaces current knot vector.	     M
*   Returns pointer to refined surface (Note a Bezier surface will be        M
* converted into a Bspline surface).                                         M
*                                                                            *
* PARAMETERS:                                                                M
*   Srf:       To refine.                                                    M
*   Dir:       Direction of refinement. Either U or V.                       M
*   Replace:   If TRUE, t holds knots in exactly the same length as the      M
*              length of the knot vector of Srf and t simply replaces the    M
*              knot vector.                                                  M
*   t:         Vector of knots with length of n.                             M
*   n:         Length of vector t.                                           M
*                                                                            *
* RETURN VALUE:                                                              M
*   CagdSrfStruct *:  A refined curve of Srf after insertion of all the      M
*                     knots as specified by vector t of length n.            M
*                                                                            *
* KEYWORDS:                                                                  M
*   CagdSrfRefineAtParams, refinement, subdivision                           M
*****************************************************************************/
CagdSrfStruct *CagdSrfRefineAtParams(CagdSrfStruct *Srf,
				     CagdSrfDirType Dir,
				     CagdBType Replace,
				     CagdRType *t,
				     int n)
{
    CagdSrfStruct *BspSrf, *TSrf;

    switch (Srf -> GType) {
	case CAGD_SBEZIER_TYPE:
    	    BspSrf = CnvrtBezier2BsplineSrf(Srf);
	    TSrf = BspSrfKnotInsertNDiff(BspSrf, Dir, Replace, t, n);
	    CagdSrfFree(BspSrf);
	    return TSrf;
	case CAGD_SBSPLINE_TYPE:
	    return BspSrfKnotInsertNDiff(Srf, Dir, Replace, t, n);
	case CAGD_SPOWER_TYPE:
	    CAGD_FATAL_ERROR(CAGD_ERR_POWER_NO_SUPPORT);
	    return NULL;
	default:
	    CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_SRF);
	    return NULL;
    }
}

/*****************************************************************************
* DESCRIPTION:                                                               M
* Given a curve Crv and a parameter value t, returns the unit tangent        M
* direction of Crv at t.                                                     M
*                                                                            *
* PARAMETERS:                                                                M
*   Crv:       To compute unit tangent vector for.                           M
*   t:         Location where to evaluate the tangent of Crv.                M
*                                                                            *
* RETURN VALUE:                                                              M
*   CagdVecStruct *:  A pointer to a static vector holding the unit tanegnt  M
*                     information.                                           M
*                                                                            *
* KEYWORDS:                                                                  M
*   CagdCrvTangent, tangent                                                  M
*****************************************************************************/
CagdVecStruct *CagdCrvTangent(CagdCrvStruct *Crv, CagdRType t)
{
    switch (Crv -> GType) {
	case CAGD_CBEZIER_TYPE:
	    return BzrCrvTangent(Crv, t);
	case CAGD_CBSPLINE_TYPE:
	    return BspCrvTangent(Crv, t);
	case CAGD_CPOWER_TYPE:
	    CAGD_FATAL_ERROR(CAGD_ERR_POWER_NO_SUPPORT);
	    return NULL;
	default:
	    CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_CRV);
	    return NULL;
    }
}

/*****************************************************************************
* DESCRIPTION:                                                               M
* Given a curve Crv and a parameter value t, returns the unit binormal       M
* direction of Crv at t.                                                     M
*                                                                            *
* PARAMETERS:                                                                M
*   Crv:       To compute unit binormal vector for.                          M
*   t:         Location where to evaluate the binormal of Crv.               M
*                                                                            *
* RETURN VALUE:                                                              M
*   CagdVecStruct *:  A pointer to a static vector holding the unit binormal M
*                     information.                                           M
*                                                                            *
* KEYWORDS:                                                                  M
*   CagdCrvBiNormal, binormal                                                M
*****************************************************************************/
CagdVecStruct *CagdCrvBiNormal(CagdCrvStruct *Crv, CagdRType t)
{
    switch (Crv -> GType) {
	case CAGD_CBEZIER_TYPE:
	    return BzrCrvBiNormal(Crv, t);
	case CAGD_CBSPLINE_TYPE:
	    return BspCrvBiNormal(Crv, t);
	case CAGD_CPOWER_TYPE:
	    CAGD_FATAL_ERROR(CAGD_ERR_POWER_NO_SUPPORT);
	    return NULL;
	default:
	    CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_CRV);
	    return NULL;
    }
}

/*****************************************************************************
* DESCRIPTION:                                                               M
* Given a curve Crv and a parameter value t, returns the unit normal         M
* direction of Crv at t.                                                     M
*                                                                            *
* PARAMETERS:                                                                M
*   Crv:       To compute unit normal vector for.                            M
*   t:         Location where to evaluate the normal of Crv.                 M
*                                                                            *
* RETURN VALUE:                                                              M
*   CagdVecStruct *:  A pointer to a static vector holding the unit normal   M
*                     information.                                           M
*                                                                            *
* KEYWORDS:                                                                  M
*   CagdCrvNormal, normal                                                    M
*****************************************************************************/
CagdVecStruct *CagdCrvNormal(CagdCrvStruct *Crv, CagdRType t)
{
    switch (Crv -> GType) {
	case CAGD_CBEZIER_TYPE:
	    return BzrCrvNormal(Crv, t);
	case CAGD_CBSPLINE_TYPE:
	    return BspCrvNormal(Crv, t);
	case CAGD_CPOWER_TYPE:
	    CAGD_FATAL_ERROR(CAGD_ERR_POWER_NO_SUPPORT);
	    return NULL;
	default:
	    CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_CRV);
	    return NULL;
    }
}

/*****************************************************************************
* DESCRIPTION:                                                               M
* Given a surface Srf and a parameter values u, v, returns the unit tangent  M
* vector of Srf in direction Dir.                                            M
*                                                                            *
* PARAMETERS:                                                                M
*   Srf:       To compute unit tangent vector for.                           M
*   u, v:      Location where to evaluate the tangent of Srf.                M
*   Dir:       Direction of tangent, Either U or V.			     *
*                                                                            *
* RETURN VALUE:                                                              M
*   CagdVecStruct *:  A pointer to a static vector holding the unit tangent  M
*                     information.                                           M
*                                                                            *
* KEYWORDS:                                                                  M
*   CagdSrfTangent, tangent                                                  M
*****************************************************************************/
CagdVecStruct *CagdSrfTangent(CagdSrfStruct *Srf,
			      CagdRType u,
			      CagdRType v,
			      CagdSrfDirType Dir)
{
    switch (Srf -> GType) {
	case CAGD_SBEZIER_TYPE:
	    return BzrSrfTangent(Srf, u, v, Dir);
	case CAGD_SBSPLINE_TYPE:
	    return BspSrfTangent(Srf, u, v, Dir);
	case CAGD_SPOWER_TYPE:
	    CAGD_FATAL_ERROR(CAGD_ERR_POWER_NO_SUPPORT);
	    return NULL;
	default:
	    CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_SRF);
	    return NULL;
    }
}

/*****************************************************************************
* DESCRIPTION:                                                               M
* Given a surface Srf and a parameter values u, v, returns the unit normal   M
* vector of Srf.                                                             M
*                                                                            *
* PARAMETERS:                                                                M
*   Srf:       To compute unit normal vector for.                            M
*   u, v:      Location where to evaluate the normal of Srf.                 M
*                                                                            *
* RETURN VALUE:                                                              M
*   CagdVecStruct *:  A pointer to a static vector holding the unit normal   M
*                     information.                                           M
*                                                                            *
* KEYWORDS:                                                                  M
*   CagdSrfNormal, normal                                                    M
*****************************************************************************/
CagdVecStruct *CagdSrfNormal(CagdSrfStruct *Srf, CagdRType u, CagdRType v)
{
    switch (Srf -> GType) {
	case CAGD_SBEZIER_TYPE:
	    return BzrSrfNormal(Srf, u, v);
	case CAGD_SBSPLINE_TYPE:
	    return BspSrfNormal(Srf, u, v);
	case CAGD_SPOWER_TYPE:
	    CAGD_FATAL_ERROR(CAGD_ERR_POWER_NO_SUPPORT);
	    return NULL;
	default:
	    CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_SRF);
	    return NULL;
    }
}

/*****************************************************************************
* DESCRIPTION:                                                               M
* Returns a new surface that is the reversed surface of Srf by reversing the M
* control mesh and the knot vector (if Bspline surface) of Srf in the U      M
* direction. See also BspKnotReverse.                                        M
*                                                                            *
* PARAMETERS:                                                                M
*   Srf:       To be reversed.                                               M
*                                                                            *
* RETURN VALUE:                                                              M
*   CagdSrfStruct *:   Reversed surface of Srf.                              M
*                                                                            *
* KEYWORDS:                                                                  M
*   CagdSrfReverse, reverse                                                  M
*****************************************************************************/
CagdSrfStruct *CagdSrfReverse(CagdSrfStruct *Srf)
{
    CagdBType
	IsNotRational = !CAGD_IS_RATIONAL_SRF(Srf);
    int i, Len, Row, Col,
	ULength = Srf -> ULength,
	VLength = Srf -> VLength,
	MaxCoord = CAGD_NUM_OF_PT_COORD(Srf -> PType);
    CagdSrfStruct
	*ReversedSrf = CagdSrfCopy(Srf);
    CagdRType *KV,
	**Points = ReversedSrf -> Points;

    switch (Srf -> GType) {
	case CAGD_SBEZIER_TYPE:
	case CAGD_SBSPLINE_TYPE:
	    break;
	case CAGD_SPOWER_TYPE:
	    CAGD_FATAL_ERROR(CAGD_ERR_POWER_NO_SUPPORT);
	    return NULL;
	default:
	    CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_SRF);
	    return NULL;
    }

    /* Reverse the Mesh: */
    Len = ULength / 2;
    for (Row = 0; Row < VLength; Row++)
	for (Col = 0; Col < Len; Col++)
	    for (i = IsNotRational; i <= MaxCoord; i++)
		SWAP(CagdRType,
		     Points[i][Row * ULength + Col],
		     Points[i][Row * ULength + ULength - Col - 1]);

    /* Reverse the U knot vector if it exists: */
    if (Srf -> GType == CAGD_SBSPLINE_TYPE &&
	Srf -> UKnotVector != NULL) {
	KV = BspKnotReverse(Srf -> UKnotVector, Srf -> UOrder + ULength);
	IritFree((VoidPtr) ReversedSrf -> UKnotVector);
	ReversedSrf -> UKnotVector = KV;
    }

    return ReversedSrf;
}

/*****************************************************************************
* DESCRIPTION:                                                               M
* Returns a new surface that is the reversed surface of Srf by flipping the  M
* U and the V directions of the surface.				     M
* See also BspKnotReverse.		                                     M
*                                                                            *
* PARAMETERS:                                                                M
*   Srf:       To be reversed.                                               M
*                                                                            *
* RETURN VALUE:                                                              M
*   CagdSrfStruct *:   Reversed surface of Srf.                              M
*                                                                            *
* KEYWORDS:                                                                  M
*   CagdSrfReverse, reverse                                                  M
*****************************************************************************/
CagdSrfStruct *CagdSrfReverse2(CagdSrfStruct *Srf)
{
    CagdBType
	IsNotRational = !CAGD_IS_RATIONAL_SRF(Srf);
    int i, Row, Col,
	ULength = Srf -> ULength,
	VLength = Srf -> VLength,
	MaxCoord = CAGD_NUM_OF_PT_COORD(Srf -> PType);
    CagdSrfStruct
	*ReversedSrf = CagdSrfCopy(Srf);
    CagdRType
	**Points = Srf -> Points,
	**RevPoints = ReversedSrf -> Points;

    switch (Srf -> GType) {
	case CAGD_SBEZIER_TYPE:
	case CAGD_SBSPLINE_TYPE:
	    break;
	case CAGD_SPOWER_TYPE:
	    CAGD_FATAL_ERROR(CAGD_ERR_POWER_NO_SUPPORT);
	    return NULL;
	default:
	    CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_SRF);
	    return NULL;
    }

    /* Reverse the Mesh: */
    for (Row = 0; Row < VLength; Row++)
	for (Col = 0; Col < ULength; Col++)
	    for (i = IsNotRational; i <= MaxCoord; i++)
		RevPoints[i][Col * VLength + Row] =
		    Points[i][Row * ULength + Col];

    /* Swap the U and the V knot vectors if the exists: */
    if (Srf -> GType == CAGD_SBSPLINE_TYPE) {
	SWAP(CagdRType *,
	     ReversedSrf -> UKnotVector, ReversedSrf -> VKnotVector);
    }

    /* And swap the orders and lengths. */
    SWAP(int, ReversedSrf -> UOrder, ReversedSrf -> VOrder);
    SWAP(int, ReversedSrf -> ULength, ReversedSrf -> VLength);

    return ReversedSrf;
}
