/*****************************************************************************
* Conversion routines from curves and surfaces to polygons and polylines.    *
*									     *
* Written by:  Gershon Elber				Ver 1.0, Apr 1992    *
*****************************************************************************/

#include "irit_sm.h"
#include "prsr_loc.h"
#include "allocate.h"
#include "ip_cnvrt.h"

#define MIN_VALID_CRV_PARAM_DOMAIN	1e-3

static RealType CubicBzrTol = 0.01;

static IPPolygonStruct *CagdPolylines2IritPolylines(CagdPolylineStruct *Polys);

/*****************************************************************************
* DESCRIPTION:                                                               *
* Routine to convert a curve into a polyline with SamplesPerCurve samples.   *
*                                                                            *
* PARAMETERS:                                                                *
*   Polys:     Polylines in cagd library format to convert.                  *
*                                                                            *
* RETURN VALUE:                                                              *
*   IPPolygonStruct *:   Same polylines in IRIT format.                      *
*****************************************************************************/
static IPPolygonStruct *CagdPolylines2IritPolylines(CagdPolylineStruct *Polys)
{
    int i, j, n;
    IPVertexStruct *V,
	*VHead = NULL,
	*VTail = NULL;
    IPPolygonStruct *P,
	*PHead = NULL;
    CagdPolylineStruct *CagdPoly;

    for (CagdPoly = Polys;
	 CagdPoly != NULL;
	 CagdPoly = CagdPoly -> Pnext) {
	n = CagdPoly -> Length;

	for (i = 0, VHead = NULL; i < n; i++) {	     /* Convert to vertices. */
	    V = IPAllocVertex(0, 0, NULL, NULL);

	    for (j = 0; j < 3; j++)	   	   /* Convert to our format. */
		V -> Coord[j] = CagdPoly -> Polyline[i].Pt[j];

	    if (VHead) {
		VTail -> Pnext = V;
		VTail = V;
	    }
	    else
		VHead = VTail = V;
	}

	P = IPAllocPolygon(0, 0, VHead, PHead);
	PHead = P;
    }

    CagdPolylineFreeList(Polys);

    return PHead;
}

/*****************************************************************************
* DESCRIPTION:                                                               M
* Routine to convert one curve into a polyline with SamplesPerCurve samples. M
*   If Optimal, attempt is made to optimally distribute the sampled points.  M
*                                                                            *
* PARAMETERS:                                                                M
*   Crv:               To approximate as a polyline .                        M
*   SamplesPerCurve:   Number of samples to compute on the polyline.         M
*   Optimal:           If FALSE samples are uniform in parametric space,     M
*                      otherwise attempt is made to sample optimally.        M
*                                                                            *
* RETURN VALUE:                                                              M
*   IPPolygonStruct *:   A polyline approximating Crv.                       M
*                                                                            *
* KEYWORDS:                                                                  M
*   IritCurve2Polylines, conversion, approximation                           M
*****************************************************************************/
IPPolygonStruct *IritCurve2Polylines(CagdCrvStruct *Crv,
				     int SamplesPerCurve,
				     int Optimal)
{
    CagdRType TMin, TMax;
    CagdPolylineStruct *CagdPoly;
    CagdBType
	NewCrv = FALSE;

    if (CAGD_IS_PERIODIC_CRV(Crv)) {
        NewCrv = TRUE;
        Crv = CnvrtPeriodic2FloatCrv(Crv);
    }

    CagdCrvDomain(Crv, &TMin, &TMax);
    if (TMax - TMin < MIN_VALID_CRV_PARAM_DOMAIN) {
	IPVertexStruct
	    *V2 = IPAllocVertex(0, 0, NULL, NULL),
	    *V1 = IPAllocVertex(0, 0, NULL, V2);
	IPPolygonStruct
	    *Pl = IPAllocPolygon(0, 0, V1, NULL);

	/* Construct a two point polyline eith two end points of curve. */
	CagdCoerceToE3(V1 -> Coord, Crv -> Points, 0, Crv -> PType);
	CagdCoerceToE3(V2 -> Coord, Crv -> Points, Crv -> Length - 1,
		       Crv -> PType);
	return Pl;
    }

    if (CAGD_IS_BSPLINE_CRV(Crv) && !BspCrvHasOpenEC(Crv)) {
	CagdCrvStruct
	    *TCrv = BspCrvOpenEnd(Crv);

	if (NewCrv)
	    CagdCrvFree(Crv);
	Crv = TCrv;
	NewCrv = TRUE;
    }

    if (SamplesPerCurve < 2)
	SamplesPerCurve = 2;

    CagdPoly = SymbCrv2Polyline(Crv, SamplesPerCurve, Optimal, TRUE);

    if (NewCrv)
	CagdCrvFree(Crv);

    return CagdPolylines2IritPolylines(CagdPoly);
}

/*****************************************************************************
* DESCRIPTION:                                                               M
* Routine to convert a single curve's control polygon into a polyline.	     M
*                                                                            *
* PARAMETERS:                                                                M
*   Crv:      To extract its control polygon as a polyline.                  M
*                                                                            *
* RETURN VALUE:                                                              M
*   IPPolygonStruct *:   A polyline representing Crv's control polygon.      M
*                                                                            *
* KEYWORDS:                                                                  M
*   IritCurve2CtlPoly, conversion                                            M
*****************************************************************************/
IPPolygonStruct *IritCurve2CtlPoly(CagdCrvStruct *Crv)
{
    CagdPolylineStruct
	*CagdPoly = CagdCrv2CtrlPoly(Crv);

    return CagdPolylines2IritPolylines(CagdPoly);
}

/*****************************************************************************
* DESCRIPTION:                                                               M
* Routine to convert a single surface into polylines with SamplesPerCurve    M
* samples, NumOfIsolines isolines into a polylines object.		     M
*   If NumOfIsolines has negative value, its absolute value is heuristically M
* used to derive a new NumOfIsolines number for it.			     M
*   If Optimal, attempt is made to optimally distribute the sampled points.  M
*                                                                            *
* PARAMETERS:                                                                M
*   Srf:               To approximate as a polyline .                        M
*   NumOfIsolines:     Number of isocurves to extract, in each direction.    M
*   SamplesPerCurve:   Number of samples to compute on the polyline.         M
*   Optimal:           If FALSE samples are uniform in parametric space,     M
*                      otherwise attempt is made to sample optimally.        M
*                                                                            *
* RETURN VALUE:                                                              M
*   IPPolygonStruct *:   A polylines object approximating Srf.               M
*                                                                            *
* KEYWORDS:                                                                  M
*   IritSurface2Polylines, conversion, approximation                         M
*****************************************************************************/
IPPolygonStruct *IritSurface2Polylines(CagdSrfStruct *Srf,
				       int NumOfIsolines[2],
				       int SamplesPerCurve,
				       int Optimal)
{
    CagdPolylineStruct *CagdPoly;
    CagdBType
	NewSrf = FALSE;

    if (CAGD_IS_PERIODIC_SRF(Srf)) {
        NewSrf = TRUE;
        Srf = CnvrtPeriodic2FloatSrf(Srf);
    }
    if (CAGD_IS_BSPLINE_SRF(Srf) && !BspSrfHasOpenEC(Srf)) {
	CagdSrfStruct
	    *TSrf = BspSrfOpenEnd(Srf);

	if (NewSrf)
	    CagdSrfFree(Srf);
	Srf = TSrf;
	NewSrf = TRUE;
    }

    if (NumOfIsolines[0] < 0) {
	if (Srf -> UOrder > 2)
	    NumOfIsolines[0] = (Srf -> ULength - NumOfIsolines[0]) / 2;
	else
	    NumOfIsolines[0] = -NumOfIsolines[0];
    }
    if (NumOfIsolines[0] < 2)
	NumOfIsolines[0] = 2;

    if (NumOfIsolines[1] < 0) {
	if (Srf -> VOrder > 2)
	    NumOfIsolines[1] = (Srf -> VLength - NumOfIsolines[1]) / 2;
	else
	    NumOfIsolines[1] = -NumOfIsolines[1];
    }
    if (NumOfIsolines[1] < 2)
	NumOfIsolines[1] = 2;

    if (SamplesPerCurve < 2)
	SamplesPerCurve = 2;

    CagdPoly = SymbSrf2Polylines(Srf, NumOfIsolines, SamplesPerCurve, Optimal);

    if (NewSrf)
	CagdSrfFree(Srf);

    return CagdPolylines2IritPolylines(CagdPoly);
}

/*****************************************************************************
* DESCRIPTION:                                                               M
* Routine to convert a single surface's control mesh into a polylines object.M
*                                                                            *
* PARAMETERS:                                                                M
*   Srf:      To extract its control mesh as a polylines.                    M
*                                                                            *
* RETURN VALUE:                                                              M
*   IPPolygonStruct *:   A polylines object representing Srf's control mesh. M
*                                                                            *
* KEYWORDS:                                                                  M
*   IritSurface2CtlMesh, conversion                                          M
*****************************************************************************/
IPPolygonStruct *IritSurface2CtlMesh(CagdSrfStruct *Srf)
{
    CagdPolylineStruct
	*CagdPoly = CagdSrf2CtrlMesh(Srf);

    return CagdPolylines2IritPolylines(CagdPoly);
}

/*****************************************************************************
* DESCRIPTION:                                                               M
* Routine to approximate a single surface by polygons.			     M
*   FourPerFlat cause four triangles from each flat patch, two otherwise.    M
*   FineNess is an estimate on number of polygons for each side of mesh if   M
* Optimal = 0, otherwise a FineNess curvature control on the maximal         M
* distance between the surface and its polygonal approximation.		     M
*   Optimal is a two digits number where the units hold the subdivision	     M
* strategy and the tens the termination criteria.			     M
*                                                                            *
* PARAMETERS:                                                                M
*   Srf:          To approximate using polygons.                             M
*   FourPerFlat:  If TRUE, four triangle per flat surface patch are created, M
*                 otherwise only two.					     M
*   FineNess:     Fineness control on polygonal approximation. The larger    M
*                 this number is the finer the approximation becomes. 10 is  M
*                 a good compromise when Optimal is FALSE.		     M
*   ComputeUV:    Do we want UV parameter values with the vertices of the    M
*                 triangles?						     M
*   Optimal:      If FALSE (0) then parametric space of Srf is sampled       M
*                 uniformely. Otherwise:				     M
*                 If the tenths digit of Optimal is equal to                 M
*                 1. Subdivide the surface alternatively in U and V.         M
*                 2. Subdivide the surface in the direction that minimizes   M
*                    the error of the approximation			     M
*                 If the units digit of Optimal is equal to		     M
*                 0. No real error analysis. The fastest way.		     M
*                 1. Use curvature surface analysis to decide where to       M
*                    subdivide. Much slower than 0.			     M
*                 2. Use a bilinear surface fit to estimate error. Somewhat  M
*                    slower than 0.					     M
*		  3. Combine 1 and 2 for a bilinar fit error measure with    M
*                    curvature analysis. Very slow.			     M
*                                                                            *
* RETURN VALUE:                                                              M
*   IPPolygonStruct *:   Resulting polygons that approximates Srf.           M
*                                                                            *
* KEYWORDS:                                                                  M
*   IritSurface2Polygons, approximation, conversion                          M
*****************************************************************************/
IPPolygonStruct *IritSurface2Polygons(CagdSrfStruct *Srf,
				      int FourPerFlat,
				      RealType FineNess,
				      int ComputeUV,
				      int Optimal)
{
    int i, j;
    CagdBType
	NewSrf = FALSE;
    char UVStr[LINE_LEN];
    VectorType Vin;
    SymbPlSubdivStrategyType Strategy;
    IPVertexStruct *V, *VHead,
	*VTail = NULL;
    IPPolygonStruct *P,
	*PHead = NULL;
    CagdPolygonStruct *CagdPolygon, *CagdPolygonHead;

    if (CAGD_IS_PERIODIC_SRF(Srf)) {
        NewSrf = TRUE;
        Srf = CnvrtPeriodic2FloatSrf(Srf);
    }
    if (CAGD_IS_BSPLINE_SRF(Srf) && !BspSrfHasOpenEC(Srf)) {
	CagdSrfStruct
	    *TSrf = BspSrfOpenEnd(Srf);

	if (NewSrf)
	    CagdSrfFree(Srf);
	Srf = TSrf;
	NewSrf = TRUE;
    }

    switch (Optimal % 10) {
	case 1:
	default:
	    Strategy = SYMB_SUBDIV_STRAT_ALTERNATE;
	    break;
	case 2:
	    Strategy = SYMB_SUBDIV_STRAT_MIN_MAX;
	    break;
	case 3:
	    Strategy = SYMB_SUBDIV_STRAT_MIN_MIN;
	    break;
    }
    Optimal /= 10;

    switch (Optimal) {
	case 0:
	    CagdPolygonHead = SymbSrf2Polygons(Srf, (int) FineNess, TRUE,
					       FourPerFlat, ComputeUV);
	    break;
	case 1:
	    SymbSrf2OptPolysCurvatureErrorPrep(Srf);
	    CagdPolygonHead = SymbSrf2OptimalPolygons(Srf, FineNess, Strategy,
						SymbSrf2OptPolysCurvatureError,
						TRUE, FourPerFlat, ComputeUV);
	    break;
	case 3:
	    /* Use bilinear fit (option 2) but with extreme curvature       */
	    /* locations as subdivision locations.			    */
	    SymbSrf2OptPolysIsoDirCurvatureErrorPrep(Srf);
	case 2:
	default:
	    CagdPolygonHead = SymbSrf2OptimalPolygons(Srf, FineNess, Strategy,
						SymbSrf2OptPolysBilinPolyError,
						TRUE, FourPerFlat, ComputeUV);
	    break;
    }

    for (CagdPolygon = CagdPolygonHead, VHead = NULL;
	 CagdPolygon != NULL;
	 CagdPolygon = CagdPolygon -> Pnext) {
	/* All polygons are triangles! */

	for (i = 0, VHead = NULL; i < 3; i++) {	     /* Convert to vertices. */
	    V = IPAllocVertex(0, 0, NULL, NULL);
	    IP_SET_NORMAL_VRTX(V);       	  /* This vertex has normal. */

	    for (j = 0; j < 3; j++)	     	   /* Convert to our format. */
		V -> Coord[j] = CagdPolygon -> Polygon[i].Pt[j];
	    for (j = 0; j < 3; j++)
		V -> Normal[j] = CagdPolygon -> Polygon[i].Nrml[j];

	    if (ComputeUV) {
		sprintf(UVStr, "%f %f",
			CagdPolygon -> Polygon[i].UV[0],
			CagdPolygon -> Polygon[i].UV[1]);
		AttrSetStrAttrib(&V -> Attrs, "uvvals", UVStr);
	    }

	    if (VHead) {
		VTail -> Pnext = V;
		VTail = V;
	    }
	    else
		VHead = VTail = V;
	}

	P = IPAllocPolygon(0, 0, VHead, PHead);
	IP_SET_CONVEX_POLY(P);

	VTail -> Pnext = VHead;
	PT_ADD(Vin, CagdPolygon -> Polygon[0].Pt,
	            CagdPolygon -> Polygon[0].Nrml);
	IritPrsrUpdatePolyPlane2(P, Vin);          /* Update plane equation. */

	if (!_IritPrsrPolyListCirc)
	    P -> PVertex -> Pnext -> Pnext -> Pnext = NULL;

	PHead = P;
    }

    CagdPolygonFreeList(CagdPolygonHead);

    if (NewSrf)
	CagdSrfFree(Srf);

    return PHead;
}

/*****************************************************************************
* DESCRIPTION:                                                               M
* Routine to convert a single trimmed surface into polylines with 	     M
* SamplesPerCurve samples, NumOfIsolines isolines into a polylines object.   M
*   If NumOfIsolines has negative value, its absolute value is heuristically M
* used to derive a new NumOfIsolines number for it.			     M
*   If Optimal, attempt is made to optimally distribute the sampled points.  M
*                                                                            *
* PARAMETERS:                                                                M
*   TrimSrf:           To approximate as a polyline .                        M
*   NumOfIsolines:     Number of isocurves to extract, in each direction.    M
*   SamplesPerCurve:   Number of samples to compute on the polyline.         M
*   Optimal:           If FALSE samples are uniform in parametric space,     M
*                      otherwise attempt is made to sample optimally.        M
*   TrimmingCurves:    Do we want the trimming curves as well.               M
*   IsoParamCurves:    Do we want trimmed isoparametric curves.              M
*                                                                            *
* RETURN VALUE:                                                              M
*   IPPolygonStruct *:   A polylines object approximating TrimSrf.           M
*                                                                            *
* KEYWORDS:                                                                  M
*   IritTrimSrf2Polylines, conversion, approximation                         M
*****************************************************************************/
IPPolygonStruct *IritTrimSrf2Polylines(TrimSrfStruct *TrimSrf,
				       int NumOfIsolines[2],
				       int SamplesPerCurve,
				       int Optimal,
				       int TrimmingCurves,
				       int IsoParamCurves)
{
    CagdPolylineStruct
	*CagdPoly1 = NULL,
	*CagdPoly2 = NULL;

    if (IsoParamCurves) {
	if (NumOfIsolines[0] < 0) {
	    if (TrimSrf -> Srf -> UOrder > 2)
		NumOfIsolines[0] = (TrimSrf -> Srf -> ULength -
				    NumOfIsolines[0]) / 2;
	    else
		NumOfIsolines[0] = -NumOfIsolines[0];
	}
	if (NumOfIsolines[0] < 2)
	    NumOfIsolines[0] = 2;

	if (NumOfIsolines[1] < 0) {
	    if (TrimSrf -> Srf -> VOrder > 2)
		NumOfIsolines[1] = (TrimSrf -> Srf -> VLength -
				    NumOfIsolines[1]) / 2;
	    else
		NumOfIsolines[1] = -NumOfIsolines[1];
	}
	if (NumOfIsolines[1] < 2)
	    NumOfIsolines[1] = 2;

	if (SamplesPerCurve < 2)
	    SamplesPerCurve = 2;

	CagdPoly1 = TrimSrf2Polylines(TrimSrf, NumOfIsolines,
				      SamplesPerCurve, Optimal);
    }

    if (TrimmingCurves) {
	CagdPoly2 = TrimCrvs2Polylines(TrimSrf, FALSE,
				       SamplesPerCurve, Optimal);
    }

    if (CagdPoly1 == NULL) {
	return CagdPolylines2IritPolylines(CagdPoly2);
    }
    else {
	CagdPolylineStruct
	    *CagdPolyLast = CagdListLast(CagdPoly1);

	CagdPolyLast -> Pnext = CagdPoly2;

	return CagdPolylines2IritPolylines(CagdPoly1);
    }
}

/*****************************************************************************
* DESCRIPTION:                                                               M
* Routine to convert a single trimmed surface's control mesh into a	     M
* polylines object.							     M
*                                                                            *
* PARAMETERS:                                                                M
*   TrimSrf:      To extract its control mesh as a polylines.                M
*                                                                            *
* RETURN VALUE:                                                              M
*   IPPolygonStruct *:   A polylines object representing TrimSrf's control   M
*		mesh.							     M
*                                                                            *
* KEYWORDS:                                                                  M
*   IritTrimSrf2CtlMesh, conversion                                          M
*****************************************************************************/
IPPolygonStruct *IritTrimSrf2CtlMesh(TrimSrfStruct *TrimSrf)
{
    CagdPolylineStruct
	*CagdPoly = CagdSrf2CtrlMesh(TrimSrf -> Srf);

    return CagdPolylines2IritPolylines(CagdPoly);
}

/*****************************************************************************
* DESCRIPTION:                                                               M
* Routine to approximate a single trivariate by polygons. Six faces of the   M
* trivariate are extracted as six surfaces that are displayed.		     M
*                                                                            *
* PARAMETERS:                                                                M
*   Trivar:       To approximate using polygons.                             M
*   FourPerFlat:  See IritSurface2Polygons.	                             M
*   FineNess:     See IritSurface2Polygons.	                             M
*   ComputeUV:    See IritSurface2Polygons.	                             M
*   Optimal:      See IritSurface2Polygons.	                             M
*                                                                            *
* RETURN VALUE:                                                              M
*   IPPolygonStruct *:   Resulting polygons that approximates Srf.           M
*                                                                            *
* KEYWORDS:                                                                  M
*   IritTrivar2Polygons, approximation, conversion                           M
*****************************************************************************/
IPPolygonStruct *IritTrivar2Polygons(TrivTVStruct *Trivar,
				     int FourPerFlat,
				     RealType FineNess,
				     int ComputeUV,
				     int Optimal)
{
    int i;
    CagdRType UMin, UMax, VMin, VMax, WMin, WMax;
    CagdSrfStruct *Srfs[6];
    IPPolygonStruct
	*Polys = NULL;

    TrivTVDomain(Trivar, &UMin, &UMax, &VMin, &VMax, &WMin, &WMax);

    Srfs[0] = TrivSrfFromTV(Trivar, UMin, TRIV_CONST_U_DIR);
    Srfs[1] = TrivSrfFromTV(Trivar, UMax, TRIV_CONST_U_DIR);
    Srfs[2] = TrivSrfFromTV(Trivar, VMin, TRIV_CONST_V_DIR);
    Srfs[3] = TrivSrfFromTV(Trivar, VMax, TRIV_CONST_V_DIR);
    Srfs[4] = TrivSrfFromTV(Trivar, WMin, TRIV_CONST_W_DIR);
    Srfs[5] = TrivSrfFromTV(Trivar, WMax, TRIV_CONST_W_DIR);

    for (i = 0; i < 6; i++) {
	IPPolygonStruct
	    *OnePolys = IritSurface2Polygons(Srfs[i], FourPerFlat, FineNess,
					     ComputeUV, Optimal),
	    *LastPoly = IritPrsrGetLastPoly(OnePolys);

	if (LastPoly != NULL) {
	    LastPoly -> Pnext = Polys;
	    Polys = OnePolys;
	}

	CagdSrfFree(Srfs[i]);
    }

    return Polys;
}

/*****************************************************************************
* DESCRIPTION:                                                               M
* Routine to convert a single trivariate function into polylines with	     M
* SamplesPerCurve samples, NumOfIsolines isolines into a polylines object.   M
*   If NumOfIsolines has negative value, its absolute value is heuristically M
* used to derive a new NumOfIsolines number for it.			     M
*   If Optimal, attempt is made to optimally distribute the sampled points.  M
*                                                                            *
* PARAMETERS:                                                                M
*   Trivar:            To approximate as a polyline .                        M
*   NumOfIsolines:     Number of isocurves to extract, in each direction.    M
*   SamplesPerCurve:   Number of samples to compute on the polyline.         M
*   Optimal:           If FALSE samples are uniform in parametric space,     M
*                      otherwise attempt is made to sample optimally.        M
*                                                                            *
* RETURN VALUE:                                                              M
*   IPPolygonStruct *:   A polylines object approximating Trivar.            M
*                                                                            *
* KEYWORDS:                                                                  M
*   IritTrivar2Polylines, conversion, approximation                          M
*****************************************************************************/
IPPolygonStruct *IritTrivar2Polylines(TrivTVStruct *Trivar,
				      int NumOfIsolines[3],
				      int SamplesPerCurve,
				      int Optimal)
{
    int Axis;
    CagdRType UMin, UMax, VMin, VMax, WMin, WMax;
    IPPolygonStruct
	*Polys = NULL;

    TrivTVDomain(Trivar, &UMin, &UMax, &VMin, &VMax, &WMin, &WMax);

    for (Axis = 0; Axis < 3; Axis++) {
	CagdRType Min, Max;
	TrivTVDirType Dir;
	int i, NumOfSrfIsolines[2];

	switch (Axis) {
	    case 0:
	        Dir = TRIV_CONST_U_DIR;
		Min = UMin;
		Max = UMax;
		NumOfSrfIsolines[0] = NumOfIsolines[1];
		NumOfSrfIsolines[1] = NumOfIsolines[2];
		break;
	    case 1:
	        Dir = TRIV_CONST_V_DIR;
		Min = VMin;
		Max = VMax;
		NumOfSrfIsolines[0] = NumOfIsolines[0];
		NumOfSrfIsolines[1] = NumOfIsolines[2];
		break;
	    default:
	    case 2:
	        Dir = TRIV_CONST_W_DIR;
		Min = WMin;
		Max = WMax;
		NumOfSrfIsolines[0] = NumOfIsolines[0];
		NumOfSrfIsolines[1] = NumOfIsolines[1];
		break;
	}

	for (i = 0; i < ABS(NumOfIsolines[Axis]); i++) {
	    CagdRType
		t = ((CagdRType) i) / (ABS(NumOfIsolines[Axis]) - 1);
	    CagdSrfStruct
		*Srf = TrivSrfFromTV(Trivar, Min * (1.0 - t) + Max * t, Dir);
	    IPPolygonStruct
		*Poly = IritSurface2Polylines(Srf, NumOfSrfIsolines,
					      SamplesPerCurve, Optimal);

	    if (Polys == NULL)
		Polys = Poly;
	    else {
		IPPolygonStruct
		    *PolyLast = IritPrsrGetLastPoly(Poly);

		PolyLast -> Pnext = Polys;
		Polys = Poly;
	    }

	    CagdSrfFree(Srf);
	}
    }

    return Polys;
}

/*****************************************************************************
* DESCRIPTION:                                                               M
* Routine to convert a single trivariate's control mesh into a polylines     M
* object.								     M
*                                                                            *
* PARAMETERS:                                                                M
*   Trivar:     To extract its control mesh as a polylines.                  M
*                                                                            *
* RETURN VALUE:                                                              M
*   IPPolygonStruct *:   A polylines object representing Trivar's   	     M
*		control mesh.						     M
*                                                                            *
* KEYWORDS:                                                                  M
*   IritTrivar2CtlMesh, conversion                                           M
*****************************************************************************/
IPPolygonStruct *IritTrivar2CtlMesh(TrivTVStruct *Trivar)
{
    CagdPolylineStruct
	*CagdPoly = TrivTV2CtrlMesh(Trivar);

    return CagdPolylines2IritPolylines(CagdPoly);
}

/*****************************************************************************
* DESCRIPTION:                                                               M
* Sets the tolerance that is used by the Bezier to Cubic Bezier 	     M
* conversion routines IritCurvesToCubicBzrCrvs and			     M
* IritSurfacesToCubicBzrCrvs						     M
*                                                                            *
* PARAMETERS:                                                                M
*   Tolerance:   Of approximation to use.                                    M
*                                                                            *
* RETURN VALUE:                                                              M
*   void                                                                     M
*                                                                            *
* KEYWORDS:                                                                  M
*   IritSetCurvesToCubicBzrTol, approximation                                M
*****************************************************************************/
void IritSetCurvesToCubicBzrTol(RealType Tolerance)
{
    CubicBzrTol = Tolerance;
}

/*****************************************************************************
* DESCRIPTION:                                                               M
* Approximates an arbitrary list of curves into cubic Beziers curves.	     M
*                                                                            *
* PARAMETERS:                                                                M
*   Crvs:          To approximate as cubic Bezier curves.                    M
*   CtlPolys:      If we want control polygons as well (DrawCtlPoly == TRUE) M
*                  they will be placed herein.				     M
*   DrawCurve:     Do we want to draw the curves?                            M
*   DrawCtlPoly:   Do we want to draw the control polygons?                  M
*   MaxArcLen:     Tolerance for cubic Bezier approximation. See function    M
*                  BzrApproxBzrCrvAsCubics.				     M
*                                                                            *
* RETURN VALUE:                                                              M
*   CagdCrvStruct *:   The cubic Bezier approximation, or NULL if DrawCurve  M
*                      is FALSE.                                             M
*                                                                            *
* KEYWORDS:                                                                  M
*   IritCurvesToCubicBzrCrvs, conversion, approximation                      M
*****************************************************************************/
CagdCrvStruct *IritCurvesToCubicBzrCrvs(CagdCrvStruct *Crvs,
					IPPolygonStruct **CtlPolys,
					CagdBType DrawCurve,
					CagdBType DrawCtlPoly,
					CagdRType MaxArcLen)
{
    CagdCrvStruct *BzrCrvs, *BzrCrv, *Crv, *CubicBzrCrvs, *TCrv,
	*AllCubicBzrCrvs = NULL;

    if (DrawCtlPoly)
	*CtlPolys = NULL;

    for (Crv = Crvs; Crv != NULL; Crv = Crv -> Pnext) {
	if (DrawCtlPoly) {
	    IPPolygonStruct
		*CagdPoly = CagdPolylines2IritPolylines(CagdCrv2CtrlPoly(Crv));

	    CagdPoly -> Pnext = *CtlPolys;
	    *CtlPolys = CagdPoly;
	}

	if (DrawCurve) {
	    if (CAGD_IS_BEZIER_CRV(Crv)) {
		CubicBzrCrvs = BzrApproxBzrCrvAsCubics(Crv, CubicBzrTol,
						       MaxArcLen, TRUE);
		for (TCrv = CubicBzrCrvs;
		     TCrv -> Pnext != NULL;
		     TCrv = TCrv -> Pnext);
		TCrv -> Pnext = AllCubicBzrCrvs;
		AllCubicBzrCrvs = CubicBzrCrvs;
	    }
	    else if (CAGD_IS_BSPLINE_CRV(Crv)) {
		BzrCrvs = CnvrtBspline2BezierCrv(Crv);
		for (BzrCrv = BzrCrvs;
		     BzrCrv != NULL;
		     BzrCrv = BzrCrv -> Pnext) {
		    CubicBzrCrvs = BzrApproxBzrCrvAsCubics(BzrCrv, CubicBzrTol,
							   MaxArcLen, TRUE);
		    for (TCrv = CubicBzrCrvs;
			 TCrv -> Pnext != NULL;
			 TCrv = TCrv -> Pnext);
		    TCrv -> Pnext = AllCubicBzrCrvs;
		    AllCubicBzrCrvs = CubicBzrCrvs;
		}
		CagdCrvFreeList(BzrCrvs);
	    }
	}
    }
    return AllCubicBzrCrvs;
}

/*****************************************************************************
* DESCRIPTION:                                                               M
* Approximates an arbitrary list of surfaces into cubic Beziers curves.	     M
*                                                                            *
* PARAMETERS:                                                                M
*   Srfs:          To approximate as cubic Bezier curves.                    M
*   CtlMeshes:     If we want control meshes as well (DrawMesh == TRUE)      M
*                  they will be placed herein.				     M
*   DrawSurface:   Do we want to draw the surfaces?                          M
*   DrawMesh:      Do we want to draw the control meshes?                    M
*   NumOfIsolines: Number of isocurves to extract, in each direction.        M
*   MaxArcLen:     Tolerance for cubic Bezier approximation. See function    M
*                  BzrApproxBzrCrvAsCubics.				     M
*                                                                            *
* RETURN VALUE:                                                              M
*   CagdCrvStruct *:  The cubic Bezier approximation, or NULL if DrawSurface M
*                     is FALSE.                                              M
*                                                                            *
* KEYWORDS:                                                                  M
*   IritSurfacesToCubicBzrCrvs, conversion, approximation                    M
*****************************************************************************/
CagdCrvStruct *IritSurfacesToCubicBzrCrvs(CagdSrfStruct *Srfs,
					  IPPolygonStruct **CtlMeshes,
					  CagdBType DrawSurface,
					  CagdBType DrawMesh,
					  int NumOfIsolines[2],
					  CagdRType MaxArcLen)
{
    CagdCrvStruct
	*AllCubicBzrCrvs = NULL;
    CagdSrfStruct *Srf;

    *CtlMeshes = NULL;

    for (Srf = Srfs; Srf != NULL; Srf = Srf -> Pnext) {
	if (DrawMesh) {
	    IPPolygonStruct
		*CagdPoly = CagdPolylines2IritPolylines(CagdSrf2CtrlMesh(Srf)),
		*CagdPolyTemp = CagdPoly;

	    if (CagdPolyTemp)
		while (CagdPolyTemp -> Pnext)
		    CagdPolyTemp = CagdPolyTemp -> Pnext;
	    CagdPolyTemp -> Pnext = *CtlMeshes;
	    *CtlMeshes = CagdPoly;
	}

	if (DrawSurface) {
	    CagdCrvStruct *CubicBzrCrvs, *TCrv,
		*Crvs = SymbSrf2Curves(Srf, NumOfIsolines);

	    CubicBzrCrvs = IritCurvesToCubicBzrCrvs(Crvs, NULL, TRUE, FALSE,
						    MaxArcLen);
	    CagdCrvFreeList(Crvs);
	    for (TCrv = CubicBzrCrvs;
		 TCrv -> Pnext != NULL;
		 TCrv = TCrv -> Pnext);
	    TCrv -> Pnext = AllCubicBzrCrvs;
	    AllCubicBzrCrvs = CubicBzrCrvs;
	}
    }
    return AllCubicBzrCrvs;
}

/*****************************************************************************
* DESCRIPTION:                                                               M
* Approximates an arbitrary list of trimmed surfaces into cubic Beziers      M
* curves. Only isoparametric curves are extarcted (no trimming curves).      M
*                                                                            *
* PARAMETERS:                                                                M
*   TrimSrfs:      To approximate as cubic Bezier curves.                    M
*   CtlMeshes:     If we want control meshes as well (DrawMesh == TRUE)      M
*                  they will be placed herein.				     M
*   DrawTrimSrf:   Do we want to draw the surfaces?                          M
*   DrawMesh:      Do we want to draw the control meshes?                    M
*   NumOfIsolines: Number of isocurves to extract, in each direction.        M
*   MaxArcLen:     Tolerance for cubic Bezier approximation. See function    M
*                  BzrApproxBzrCrvAsCubics.				     M
*                                                                            *
* RETURN VALUE:                                                              M
*   CagdCrvStruct *:  The cubic Bezier approximation, or NULL if DrawSurface M
*                     is FALSE.                                              M
*                                                                            *
* KEYWORDS:                                                                  M
*   IritTrimSrfsToCubicBzrCrvs, conversion, approximation                    M
*****************************************************************************/
CagdCrvStruct *IritTrimSrfsToCubicBzrCrvs(TrimSrfStruct *TrimSrfs,
					  IPPolygonStruct **CtlMeshes,
					  CagdBType DrawTrimSrf,
					  CagdBType DrawMesh,
					  int NumOfIsolines[2],
					  CagdRType MaxArcLen)
{
    CagdCrvStruct
	*AllCubicBzrCrvs = NULL;
    TrimSrfStruct *TrimSrf;

    if (DrawMesh)
	*CtlMeshes = NULL;

    for (TrimSrf = TrimSrfs; TrimSrf != NULL; TrimSrf = TrimSrf -> Pnext) {
	if (DrawMesh) {
	    IPPolygonStruct
		*CagdPoly = CagdPolylines2IritPolylines(
					   CagdSrf2CtrlMesh(TrimSrf -> Srf)),
		*CagdPolyTemp = CagdPoly;

	    if (CagdPolyTemp)
		while (CagdPolyTemp -> Pnext)
		    CagdPolyTemp = CagdPolyTemp -> Pnext;
	    CagdPolyTemp -> Pnext = *CtlMeshes;
	    *CtlMeshes = CagdPoly;
	}

	if (DrawTrimSrf) {
	    CagdCrvStruct *CubicBzrCrvs, *TCrv,
		*Crvs = TrimSrf2Curves(TrimSrf, NumOfIsolines);

	    CubicBzrCrvs = IritCurvesToCubicBzrCrvs(Crvs, NULL, TRUE, FALSE,
						    MaxArcLen);
	    CagdCrvFreeList(Crvs);
	    for (TCrv = CubicBzrCrvs;
		 TCrv -> Pnext != NULL;
		 TCrv = TCrv -> Pnext);
	    TCrv -> Pnext = AllCubicBzrCrvs;
	    AllCubicBzrCrvs = CubicBzrCrvs;
	}
    }
    return AllCubicBzrCrvs;
}

/*****************************************************************************
* DESCRIPTION:                                                               M
* Routine to convert a single trivariate function into polylines with	     M
* SamplesPerCurve samples, NumOfIsolines isolines into a polylines object.   M
*   If NumOfIsolines has negative value, its absolute value is heuristically M
* used to derive a new NumOfIsolines number for it.			     M
*   If Optimal, attempt is made to optimally distribute the sampled points.  M
*                                                                            *
* PARAMETERS:                                                                M
*   Trivar:        To approximate as a set of cubic bezier curves.	     M
*   CtlMeshes:     If we want control meshes as well (DrawMesh == TRUE)      M
*                  they will be placed herein.				     M
*   DrawTrivar:    Do we want to draw the trivariate function?               M
*   DrawMesh:      Do we want to draw the control meshes?                    M
*   NumOfIsolines: Number of isocurves to extract, in each direction.        M
*   MaxArcLen:     Tolerance for cubic Bezier approximation. See function    M
*                  BzrApproxBzrCrvAsCubics.				     M
*                                                                            *
* RETURN VALUE:                                                              M
*   CagdCrvStruct *:   A list of curves approximating Trivar.                M
*                                                                            *
* KEYWORDS:                                                                  M
*   IritTrivarToCubicBzrCrvs, conversion, approximation                      M
*****************************************************************************/
CagdCrvStruct *IritTrivarToCubicBzrCrvs(TrivTVStruct *Trivar,
					  IPPolygonStruct **CtlMeshes,
					  CagdBType DrawSurface,
					  CagdBType DrawMesh,
					  int NumOfIsolines[2],
					  CagdRType MaxArcLen)
{
    int Axis;
    CagdRType UMin, UMax, VMin, VMax, WMin, WMax;
    CagdCrvStruct
	*Crvs = NULL;

    TrivTVDomain(Trivar, &UMin, &UMax, &VMin, &VMax, &WMin, &WMax);

    *CtlMeshes = NULL;

    for (Axis = 0; Axis < 3; Axis++) {
	CagdRType Min, Max;
	TrivTVDirType Dir;
	int i, NumOfSrfIsolines[2];

	switch (Axis) {
	    case 0:
	        Dir = TRIV_CONST_U_DIR;
		Min = UMin;
		Max = UMax;
		NumOfSrfIsolines[0] = NumOfIsolines[1];
		NumOfSrfIsolines[1] = NumOfIsolines[2];
		break;
	    case 1:
	        Dir = TRIV_CONST_V_DIR;
		Min = VMin;
		Max = VMax;
		NumOfSrfIsolines[0] = NumOfIsolines[0];
		NumOfSrfIsolines[1] = NumOfIsolines[2];
		break;
	    default:
	    case 2:
	        Dir = TRIV_CONST_W_DIR;
		Min = WMin;
		Max = WMax;
		NumOfSrfIsolines[0] = NumOfIsolines[0];
		NumOfSrfIsolines[1] = NumOfIsolines[1];
		break;
	}

	for (i = 0; i < ABS(NumOfIsolines[Axis]); i++) {
	    CagdRType
		t = ((CagdRType) i) / (ABS(NumOfIsolines[Axis]) - 1);
	    CagdSrfStruct
		*Srf = TrivSrfFromTV(Trivar, Min * (1.0 - t) + Max * t, Dir);
	    IPPolygonStruct *CtlMesh;
	    CagdCrvStruct
		*Crv = IritSurfacesToCubicBzrCrvs(Srf, &CtlMesh,
						  DrawSurface,
						  DrawMesh,
						  NumOfSrfIsolines,
						  MaxArcLen);

	    if (Crv == NULL)
		Crvs = Crv;
	    else {
		CagdCrvStruct
		    *CrvLast = CagdListLast(Crv);

		if (CrvLast) {
		    CrvLast -> Pnext = Crvs;
		    Crvs = Crv;
		}
	    }

	    if (*CtlMeshes == NULL)
		*CtlMeshes = CtlMesh;
	    else {
		IPPolygonStruct
		    *MeshLast = IritPrsrGetLastPoly(CtlMesh);

		if (MeshLast) {
		    MeshLast -> Pnext = *CtlMeshes;
		    *CtlMeshes = CtlMesh;
		}
	    }

	    CagdSrfFree(Srf);
	}
    }

    return Crvs;
}
