/******************************************************************************
* Bzr2Poly.c - Bezier to polygon/polylines conversion routines.		      *
*******************************************************************************
* Written by Gershon Elber, Mar. 90.					      *
******************************************************************************/

#include "cagd_loc.h"

/*****************************************************************************
* Routine to convert a single bezier surface to set of triangles	     *
* approximating it. FineNess is a finess control on result and the bigger it *
* is more triangles may result. a value of 10 is a good start value.	     *
* NULL is returned in case of an error, otherwise list of CagdPolygonStruct. *
*****************************************************************************/
CagdPolygonStruct *BzrSrf2Polygons(CagdSrfStruct *Srf, int FineNess,
			 CagdBType ComputeNormals, CagdBType FourPerFlat)
{
    int i, j, FineNessU1, FineNessV1, FineNessU, FineNessV, BaseIndex;
    CagdRType *Pt;
    CagdPointType
	PType = Srf -> PType;
    CagdPtStruct PtCenter, *Pt1, *Pt2, *Pt3, *Pt4, *PtMesh, *PtMeshPtr;
    CagdVecStruct NlCenter, *Nl1, *Nl2, *Nl3, *Nl4, *PtNrml, *PtNrmlPtr;
    CagdPolygonStruct *Poly,
	*PolyHead = NULL;

    if (!CAGD_IS_BEZIER_SRF(Srf)) return NULL;

    /* Simple heuristic to estimate how many samples to compute. */
    FineNessU = Srf -> UOrder * FineNess / 10;
    FineNessV = Srf -> VOrder * FineNess / 10;

    if (FineNessU < 2) FineNessU = 2;
    if (FineNessV < 2) FineNessV = 2;

    switch (_CagdLin2Poly) {
	case CAGD_REG_POLY_PER_LIN:
	    break;
    	case CAGD_ONE_POLY_PER_LIN:
	    if (Srf -> UOrder == 2) FineNessU = 2;
	    if (Srf -> VOrder == 2) FineNessV = 2;
	    break;
    	case CAGD_ONE_POLY_PER_COLIN:
	    break;
    }

    FineNessU1 = FineNessU - 1;
    FineNessV1 = FineNessV - 1;

    /* Allocate a mesh to hold all vertices so common vertices need not be   */
    /* Evaluated twice, and evaluate the surface at these mesh points.	     */
    PtMeshPtr = PtMesh = (CagdPtStruct *) CagdMalloc(FineNessU * FineNessV *
							sizeof(CagdPtStruct));

    for (i = 0; i < FineNessU; i++)
	for (j = 0; j < FineNessV; j++) {
	    Pt = BzrSrfEvalAtParam(Srf, ((CagdRType) i) / FineNessU1,
				        ((CagdRType) j) / FineNessV1);
	    CagdCoerceToE3(PtMeshPtr -> Pt, &Pt, -1, PType);
	    PtMeshPtr++;
	}

    if (ComputeNormals) {
	PtNrmlPtr = PtNrml = (CagdVecStruct *) CagdMalloc(FineNessU * FineNessV *
						      sizeof(CagdVecStruct));
	for (i = 0; i < FineNessU; i++)
	    for (j = 0; j < FineNessV; j++) {
		Nl1 = BzrSrfNormal(Srf, ((CagdRType) i) / FineNessU1,
					((CagdRType) j) / FineNessV1);
		CAGD_COPY_VECTOR(*PtNrmlPtr, *Nl1);
		PtNrmlPtr++;
	    }
    }

    /* Now that we have the mesh, create the polygons. */
    for (i = 0; i < FineNessU1; i++)
	for (j = 0; j < FineNessV1; j++) {
	    BaseIndex = i * FineNessV + j;
	    Pt1 = &PtMesh[BaseIndex];        /* Cache the four flat corners. */
	    Pt2 = &PtMesh[BaseIndex + 1];
	    Pt3 = &PtMesh[BaseIndex + FineNessV + 1];
	    Pt4 = &PtMesh[BaseIndex + FineNessV];
	    if (ComputeNormals) {
		Nl1 = &PtNrml[BaseIndex];
		Nl2 = &PtNrml[BaseIndex + 1];
		Nl3 = &PtNrml[BaseIndex + FineNessV + 1];
		Nl4 = &PtNrml[BaseIndex + FineNessV];
	    }

	    if (FourPerFlat) {  /* Eval middle point and create 4 triangles. */
		Pt = BzrSrfEvalAtParam(Srf, (i + 0.5) / FineNessU1,
				            (j + 0.5) / FineNessV1);
		CagdCoerceToE3(PtCenter.Pt, &Pt, -1, PType);
		if (ComputeNormals) {
		    /* Average the four normals to find the middle one. */
		    CAGD_COPY_VECTOR(NlCenter, *Nl1);
		    CAGD_ADD_VECTOR(NlCenter, *Nl2);
		    CAGD_ADD_VECTOR(NlCenter, *Nl3);
		    CAGD_ADD_VECTOR(NlCenter, *Nl4);
		    CAGD_NORMALIZE_VECTOR(NlCenter);
		}

		Poly = _CagdMakePolygon(ComputeNormals, Pt1, Pt2, &PtCenter,
							Nl1, Nl2, &NlCenter);
		CAGD_LIST_PUSH(Poly, PolyHead);
		Poly = _CagdMakePolygon(ComputeNormals, Pt2, Pt3, &PtCenter,
						        Nl2, Nl3, &NlCenter);
		CAGD_LIST_PUSH(Poly, PolyHead);
		Poly = _CagdMakePolygon(ComputeNormals, Pt3, Pt4, &PtCenter,
						        Nl3, Nl4, &NlCenter);
		CAGD_LIST_PUSH(Poly, PolyHead);
		Poly = _CagdMakePolygon(ComputeNormals, Pt4, Pt1, &PtCenter,
						        Nl4, Nl1, &NlCenter);
		CAGD_LIST_PUSH(Poly, PolyHead);
	    }
	    else {			   /* Only two along the diagonal... */
		Poly = _CagdMakePolygon(ComputeNormals, Pt1, Pt2, Pt3,
						        Nl1, Nl2, Nl3);
		CAGD_LIST_PUSH(Poly, PolyHead);
		Poly = _CagdMakePolygon(ComputeNormals, Pt3, Pt4, Pt1,
						        Nl3, Nl4, Nl1);
		CAGD_LIST_PUSH(Poly, PolyHead);
	    }
	}

    CagdFree((VoidPtr) PtMesh);
    if (ComputeNormals) CagdFree((VoidPtr) PtNrml);

    return PolyHead;
}

/*****************************************************************************
* Routine to convert a single bezier surface to NumOfIsolines polylines list *
* in each param. direction with SamplesPerCurve in each isoparametric curve. *
* Polyline are always E3 of CagdPolylineStruct type.			     *
* Iso parametric curves are sampled equally spaced in parametric space.	     *
* NULL is returned in case of an error, otherwise list of CagdPolylineStruct.*
*****************************************************************************/
CagdPolylineStruct *BzrSrf2Polylines(CagdSrfStruct *Srf, int NumOfIsocurves,
							 int SamplesPerCurve)
{
    int i;
    CagdRType t;
    CagdCrvStruct *Crv;
    CagdPolylineStruct *PolyList = NULL, *Poly;

    if (!CAGD_IS_BEZIER_SRF(Srf)) return NULL;

    /* Make sure requested format is something reasonable. */
    if (SamplesPerCurve < 1) SamplesPerCurve = 1;
    if (SamplesPerCurve > CAGD_MAX_BEZIER_CACHE_ORDER)
	SamplesPerCurve = CAGD_MAX_BEZIER_CACHE_ORDER;
    if (NumOfIsocurves < 2) NumOfIsocurves = 2;

    for (i = 0; i < NumOfIsocurves; i++) {
	t = ((CagdRType) i) / (NumOfIsocurves - 1);
	if (t > 1.0) t = 1.0;	   	      /* In case of round off error. */

	Crv = BzrSrfCrvFromSrf(Srf, t, CAGD_CONST_U_DIR);
	Poly = BzrCrv2Polyline(Crv, SamplesPerCurve);
	Poly -> Pnext = PolyList;
	PolyList = Poly;
	CagdCrvFree(Crv);

	Crv = BzrSrfCrvFromSrf(Srf, t, CAGD_CONST_V_DIR);
	Poly = BzrCrv2Polyline(Crv, SamplesPerCurve);
	Poly -> Pnext = PolyList;
	PolyList = Poly;
	CagdCrvFree(Crv);
    }

    return PolyList;
}

/*****************************************************************************
* Routine to convert a single bezier curve to polyline with SamplesPerCurve  *
* samples. Polyline is always E3 of CagdPolylineStruct type.		     *
* Curve is sampled equally spaced in parametric space.			     *
* NULL is returned in case of an error, otherwise CagdPolylineStruct.	     *
*****************************************************************************/
CagdPolylineStruct *BzrCrv2Polyline(CagdCrvStruct *Crv, int SamplesPerCurve)
{
    int i, j,
	n = 1 << SamplesPerCurve,
	IsNotRational = !CAGD_IS_RATIONAL_CRV(Crv),
	MaxCoord = CAGD_NUM_OF_PT_COORD(Crv -> PType);
    CagdRType *Polyline[CAGD_MAX_PT_SIZE], Scaler;
    CagdPtStruct *NewPolyline;
    CagdPolylineStruct *P;

    if (!CAGD_IS_BEZIER_CRV(Crv)) return NULL;

    /* Make sure requested format is something reasonable. */
    if (SamplesPerCurve < 1) SamplesPerCurve = 1;
    if (SamplesPerCurve > CAGD_MAX_BEZIER_CACHE_ORDER)
	SamplesPerCurve = CAGD_MAX_BEZIER_CACHE_ORDER;

    P = CagdPolylineNew(n);
    NewPolyline = P -> Polyline;

    /* Allocate temporary memory to hold evaluated curve. */
    for (i = 0; i < CAGD_MAX_PT_SIZE; i++)
	Polyline[i] = (CagdRType *) CagdMalloc(sizeof(CagdRType) * n);

    if (MaxCoord > 3) MaxCoord = 3;

    BzrCrvEvalToPolyline(Crv, SamplesPerCurve, Polyline);

    for (i = n - 1; i >= 0; i--) {	          /* Convert to E3 polyline. */
	if (IsNotRational)
	    Scaler = 1.0;
	else
	    Scaler = Polyline[0][i];

	for (j = 0; j < MaxCoord; j++)
	   NewPolyline[i].Pt[j] = Polyline[j+1][i] / Scaler;
	for (j = MaxCoord; j < 3; j++)
	   NewPolyline[i].Pt[j] = 0.0;
    }
    for (i = 0; i < CAGD_MAX_PT_SIZE; i++) CagdFree((VoidPtr) Polyline[i]);

    return P;
}