/******************************************************************************
* ff_krnl.c - Cumputation of kernel for freeforms.			      *
*******************************************************************************
* (C) Gershon Elber, Technion, Israel Institute of Technology                 *
*******************************************************************************
* Written by Gershon Elber, Oct. 2002.					      *
******************************************************************************/

#include "irit_sm.h"
#include "allocate.h"
#include "attribut.h"
#include "iritprsr.h"
#include "cagd_lib.h"
#include "symb_lib.h"
#include "user_lib.h"
#include "geom_lib.h"
#include "bool_lib.h"
#include "ip_cnvrt.h"

STATIC_DATA jmp_buf LclLongJumpBuffer;	/* Used to trap fatal Boolean error. */

static void FFKrnlBoolFatalError(BoolFatalErrorType ErrID);
static IPObjectStruct *FFKrnlBooleanAND(IPObjectStruct *PObj1,
					IPObjectStruct *PObj2);

/*****************************************************************************
* DESCRIPTION:                                                               M
*   Computes the kernel of a freeform closed C^1 continuous surface.         M
*   The parabolic curves are computed and surface tangent planes swept along M
* these parabolic curves clip the volume, resulting with the kerenl.         M
*                                                                            *
* PARAMETERS:                                                                M
*   Srf:        To compute its kernel.					     M
*   Tolerance:  Accuracy of parabolic piecewise linear approximation.        M
*   SkipRate:   Step size over the parabolic points, 1 to process them all.  M
*                                                                            *
* RETURN VALUE:                                                              M
*   IPObjectStruct *:   A polyehdra approximating the kernel, or NULL if     M
*		empty set.						     M
*                                                                            *
* SEE ALSO:                                                                  M
*   SymbSrfGaussCurvature, UserCntrSrfWithPlane                              M
*                                                                            *
* KEYWORDS:                                                                  M
*   UserSrfKernel                                                            M
*****************************************************************************/
IPObjectStruct *UserSrfKernel(CagdSrfStruct *Srf,
			      CagdRType Tolerance,
			      int SkipRate)
{
    STATIC_DATA PlaneType
	GaussPlane = { 1.0, 0.0, 0.0, 1.050964e-12 };
    CagdBType
	OldInterpFlag = BspMultInterpFlag(FALSE);
    CagdRType MaxSize;
    CagdVType Sizes;
    CagdBBoxStruct BBox;
    CagdSrfStruct
	*GaussSrf = SymbSrfGaussCurvature(Srf, TRUE);
    IPPolygonStruct *Pl,
	*Plls = UserCntrSrfWithPlane(GaussSrf, GaussPlane, Tolerance);
    IPObjectStruct *PolyObj, *PlnObj;
    BoolFatalErrorFuncType
        OldBoolErrorFunc = BoolSetFatalErrorFunc(FFKrnlBoolFatalError);

    CagdSrfFree(GaussSrf);
    BspMultInterpFlag(OldInterpFlag);

    CagdSrfBBox(Srf, &BBox);
    VEC_SUB(Sizes, BBox.Max, BBox.Min);
    MaxSize = MAX(Sizes[0], MAX(Sizes[1], Sizes[2])) * 4;

    /* Build the plane to make the trimmings with. */
    {
        int Rvrsd;
        PointType Pt1, Pt2, Pt3, PtIn;

	Pt1[0] =  MaxSize * cos(M_PI / 6);
	Pt1[1] = -MaxSize * sin(M_PI / 6);
	Pt1[2] = 0;

	Pt2[0] = 0;
	Pt2[1] = MaxSize;
	Pt2[2] = 0;

	Pt3[0] = -MaxSize * cos(M_PI / 6);
	Pt3[1] = -MaxSize * sin(M_PI / 6);
	Pt3[2] = 0;

	PtIn[0] = 0;
	PtIn[1] = 0;
	PtIn[2] = 1;

	PlnObj = IPGenPOLYObject(PrimGenPolygon3Vrtx(Pt1, Pt2, Pt3,
						     PtIn, &Rvrsd, NULL));
    }

    /* Create a polygonal version of the surface. */
    PolyObj = IPGenPOLYObject(IPSurface2Polygons(Srf, FALSE, Tolerance,
						 FALSE, FALSE, FALSE));

    /* Trim the object with tangent planes at all parabolic points. */
    for (Pl = Plls; Pl != NULL; Pl = Pl -> Pnext) {
        int Skip = SkipRate;
	IPVertexStruct
	    *V = Pl -> PVertex;

	for ( ; V != NULL; V = V -> Pnext) {
	    if (--Skip == 0) {
	        CagdRType
		    *R = CagdSrfEval(Srf, V -> Coord[1], V -> Coord[2]);
		CagdPType Pos;
		CagdVecStruct
		    *Nrml = CagdSrfNormal(Srf, V -> Coord[1], V -> Coord[2],
					  TRUE);
		MatrixType Mat1, Mat2;
		IPObjectStruct *NewPolyObj, *TPlnObj;

		CagdCoerceToE3(Pos, &R, -1, Srf -> PType);

		/* Map the plane and use Booleans to perform the trimming. */
		MatGenMatTrans(Pos[0], Pos[1], Pos[2], Mat1);
		GMGenMatrixZ2Dir(Mat2, Nrml -> Vec);
		MatMultTwo4by4(Mat1, Mat2, Mat1);
		TPlnObj = GMTransformObject(PlnObj, Mat1);

		if ((NewPolyObj = FFKrnlBooleanAND(PolyObj,
						   TPlnObj)) != NULL) {
		    IPFreeObject(PolyObj);
		    PolyObj = NewPolyObj;
		}

		IPFreeObject(TPlnObj);

		Skip = SkipRate;
	    }
	}
    }

    IPFreePolygonList(Plls);
    IPFreeObject(PlnObj);

    BoolSetFatalErrorFunc(OldBoolErrorFunc);

    return PolyObj;
}

/*****************************************************************************
* DESCRIPTION:                                                               *
* Traps Bool_lib errors right here. Call back function of bool_lib.	     *
*                                                                            *
* PARAMETERS:                                                                *
*   ErrID:    Error number in bool_lib library.                              *
*                                                                            *
* RETURN VALUE:                                                              *
*   void                                                                     *
*****************************************************************************/
static void FFKrnlBoolFatalError(BoolFatalErrorType ErrID)
{
    longjmp(LclLongJumpBuffer, 1);
}

/*****************************************************************************
* DESCRIPTION:                                                               *
*   A Boolean between two objects while we capture here failures.            *
*   No errors in the Booleans are broadcasted beyond the returned NULL.      *
*                                                                            *
* PARAMETERS:                                                                *
*   PObj1, PObj2:   The two objects to compute their AND.                    *
*                                                                            *
* RETURN VALUE:                                                              *
*   IPObjectStruct *:    Result of computation or NULL if error.             *
*****************************************************************************/
static IPObjectStruct *FFKrnlBooleanAND(IPObjectStruct *PObj1,
					IPObjectStruct *PObj2)
{
    if (setjmp(LclLongJumpBuffer) == 0) {    /* Its the setjmp itself call! */
        return BooleanAND(PObj1, PObj2);
    }
    else {                   /* We gain control from fatal error long jump. */
        return NULL;
    }
}