summaryrefslogtreecommitdiffstats
path: root/src/libs/lprof/cmshull.cpp
diff options
context:
space:
mode:
Diffstat (limited to 'src/libs/lprof/cmshull.cpp')
-rw-r--r--src/libs/lprof/cmshull.cpp1480
1 files changed, 1480 insertions, 0 deletions
diff --git a/src/libs/lprof/cmshull.cpp b/src/libs/lprof/cmshull.cpp
new file mode 100644
index 00000000..05b8dfa3
--- /dev/null
+++ b/src/libs/lprof/cmshull.cpp
@@ -0,0 +1,1480 @@
+/* */
+/* Little cms - profiler construction set */
+/* Copyright (C) 1998-2001 Marti Maria <marti@littlecms.com> */
+/* */
+/* THIS SOFTWARE IS PROVIDED "AS-IS" AND WITHOUT WARRANTY OF ANY KIND, */
+/* EXPRESS, IMPLIED OR OTHERWISE, INCLUDING WITHOUT LIMITATION, ANY */
+/* WARRANTY OF MERCHANTABILITY OR FITNESS FOR A PARTICULAR PURPOSE. */
+/* */
+/* IN NO EVENT SHALL MARTI MARIA BE LIABLE FOR ANY SPECIAL, INCIDENTAL, */
+/* INDIRECT OR CONSEQUENTIAL DAMAGES OF ANY KIND, */
+/* OR ANY DAMAGES WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR PROFITS, */
+/* WHETHER OR NOT ADVISED OF THE POSSIBILITY OF DAMAGE, AND ON ANY THEORY OF */
+/* LIABILITY, ARISING OUT OF OR IN CONNECTION WITH THE USE OR PERFORMANCE */
+/* OF THIS SOFTWARE. */
+/* */
+/* This file is free software; you can redistribute it and/or modify it */
+/* under the terms of the GNU General Public License as published by */
+/* the Free Software Foundation; either version 2 of the License, or */
+/* (at your option) any later version. */
+/* */
+/* This program is distributed in the hope that it will be useful, but */
+/* WITHOUT ANY WARRANTY; without even the implied warranty of */
+/* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU */
+/* General Public License for more details. */
+/* */
+/* You should have received a copy of the GNU General Public License */
+/* along with this program; if not, write to the Free Software */
+/* Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA. */
+/* */
+/* As a special exception to the GNU General Public License, if you */
+/* distribute this file as part of a program that contains a */
+/* configuration script generated by Autoconf, you may include it under */
+/* the same distribution terms that you use for the rest of that program. */
+/* */
+/* Version 1.09a */
+
+
+#include "lcmsprf.h"
+
+/* Convex hull management */
+
+LCMSHANDLE cdecl cmsxHullInit(void);
+void cdecl cmsxHullDone(LCMSHANDLE hHull);
+BOOL cdecl cmsxHullAddPoint(LCMSHANDLE hHull, int x, int y, int z);
+BOOL cdecl cmsxHullComputeHull(LCMSHANDLE hHull);
+char cdecl cmsxHullCheckpoint(LCMSHANDLE hHull, int x, int y, int z);
+BOOL cdecl cmsxHullDumpVRML(LCMSHANDLE hHull, const char* fname);
+
+/* --------------------------------------------------------------------- */
+
+
+
+/* This method is described in "Computational Geometry in C" Chapter 4. */
+/* */
+/* -------------------------------------------------------------------- */
+/* This code is Copyright 1998 by Joseph O'Rourke. It may be freely */
+/* redistributed in its entirety provided that this copyright notice is */
+/* not removed. */
+/* -------------------------------------------------------------------- */
+
+#define SWAP(t,x,y) { t = x; x = y; y = t; }
+
+#define XFREE(p)
+/* if (p) { free ((char *) p); p = NULL; } */
+
+
+#define ADD( head, p ) if ( head ) { \
+ p->Next = head; \
+ p->Prev = head->Prev; \
+ head->Prev = p; \
+ p->Prev->Next = p; \
+ } \
+ else { \
+ head = p; \
+ head->Next = head->Prev = p; \
+ }
+
+#define XDELETE( head, p ) if ( head ) { \
+ if ( head == head->Next ) \
+ head = NULL; \
+ else if ( p == head ) \
+ head = head->Next; \
+ p->Next->Prev = p->Prev; \
+ p->Prev->Next = p->Next; \
+ XFREE( p ); \
+ }
+
+/* Define Vertex indices. */
+#define X 0
+#define Y 1
+#define Z 2
+
+/* Define structures for vertices, edges and faces */
+
+typedef struct _vertex_struct VERTEX,FAR *LPVERTEX;
+typedef struct _edge_struct EDGE, FAR *LPEDGE;
+typedef struct _face_struct FACE, FAR *LPFACE;
+
+
+struct _edge_struct {
+
+ LPFACE AdjFace[2];
+ LPVERTEX EndPts[2];
+ LPFACE NewFace; /* pointer to incident cone face. */
+ BOOL DoDelete; /* T iff Edge should be delete. */
+
+ LPEDGE Next, Prev;
+};
+
+struct _face_struct {
+
+ LPEDGE Edge[3];
+ LPVERTEX Vertex[3];
+ BOOL Visible; /* T iff face Visible from new point. */
+
+ LPFACE Next, Prev;
+};
+
+struct _vertex_struct {
+
+ int v[3];
+ int vnum;
+ LPEDGE duplicate; /* pointer to incident cone Edge (or NULL) */
+ BOOL onhull; /* T iff point on hull. */
+ BOOL mark; /* T iff point already processed. */
+
+ LPVERTEX Next, Prev;
+};
+
+/* Define flags */
+
+#define ONHULL true
+#define REMOVED true
+#define VISIBLE true
+#define PROCESSED true
+#define SAFE 1000000 /* Range of safe coord values. */
+
+#define DIM 3 /* Dimension of points */
+typedef int VEC3I[DIM]; /* Type integer point */
+
+#define PMAX 10000 /* Max # of pts */
+
+typedef struct {
+
+ /* Global variable definitions */
+ LPVERTEX vertices;
+ LPEDGE edges;
+ LPFACE faces;
+
+ VEC3I Vertices[PMAX]; /* All the points */
+ VEC3I Faces[PMAX]; /* Each triangle face is 3 indices */
+ VEC3I Box[PMAX][2]; /* Box around each face */
+
+
+ VEC3I bmin, bmax;
+ int radius;
+ int vnumCounter;
+
+ int nfaces;
+ int nvertex;
+
+} HULL, FAR* LPHULL;
+
+/* static HULL Global; */
+
+
+/*---------------------------------------------------------------------
+MakeNullVertex: Makes a Vertex, nulls out fields.
+---------------------------------------------------------------------*/
+
+static
+LPVERTEX MakeNullVertex(LPHULL hull)
+{
+ LPVERTEX v;
+
+ v = (LPVERTEX) malloc(sizeof(VERTEX));
+ if (!v) return NULL;
+
+ v->duplicate = NULL;
+ v->onhull = !ONHULL;
+ v->mark = !PROCESSED;
+ ADD( hull->vertices, v );
+
+ return v;
+}
+
+
+
+/*---------------------------------------------------------------------
+MakeNullEdge creates a new cell and initializes all pointers to NULL
+and sets all flags to off. It returns a pointer to the empty cell.
+---------------------------------------------------------------------*/
+static
+LPEDGE MakeNullEdge(LPHULL hull)
+{
+ LPEDGE e;
+
+ e = (LPEDGE) malloc(sizeof(EDGE));
+ if (!e) return NULL;
+
+ e->AdjFace[0] = e->AdjFace[1] = e->NewFace = NULL;
+ e->EndPts[0] = e->EndPts[1] = NULL;
+ e->DoDelete = !REMOVED;
+ ADD( hull->edges, e );
+ return e;
+}
+
+/*--------------------------------------------------------------------
+MakeNullFace creates a new face structure and initializes all of its
+flags to NULL and sets all the flags to off. It returns a pointer
+to the empty cell.
+---------------------------------------------------------------------*/
+static
+LPFACE MakeNullFace(LPHULL hull)
+{
+ LPFACE f;
+ int i;
+
+ f = (LPFACE) malloc(sizeof(FACE));
+ if (!f) return NULL;
+
+ for ( i=0; i < 3; ++i ) {
+ f->Edge[i] = NULL;
+ f->Vertex[i] = NULL;
+ }
+ f->Visible = !VISIBLE;
+ ADD( hull->faces, f );
+ return f;
+}
+
+
+
+/*---------------------------------------------------------------------
+MakeFace creates a new face structure from three vertices (in ccw
+order). It returns a pointer to the face.
+---------------------------------------------------------------------*/
+static
+LPFACE MakeFace(LPHULL hull, LPVERTEX v0, LPVERTEX v1, LPVERTEX v2, LPFACE fold)
+{
+ LPFACE f;
+ LPEDGE e0, e1, e2;
+
+ /* Create edges of the initial triangle. */
+ if( !fold ) {
+ e0 = MakeNullEdge(hull);
+ e1 = MakeNullEdge(hull);
+ e2 = MakeNullEdge(hull);
+ }
+ else { /* Copy from fold, in reverse order. */
+ e0 = fold->Edge[2];
+ e1 = fold->Edge[1];
+ e2 = fold->Edge[0];
+ }
+ e0->EndPts[0] = v0; e0->EndPts[1] = v1;
+ e1->EndPts[0] = v1; e1->EndPts[1] = v2;
+ e2->EndPts[0] = v2; e2->EndPts[1] = v0;
+
+ /* Create face for triangle. */
+ f = MakeNullFace(hull);
+ f->Edge[0] = e0; f->Edge[1] = e1; f->Edge[2] = e2;
+ f->Vertex[0] = v0; f->Vertex[1] = v1; f->Vertex[2] = v2;
+
+ /* Link edges to face. */
+ e0->AdjFace[0] = e1->AdjFace[0] = e2->AdjFace[0] = f;
+
+ return f;
+}
+
+/*---------------------------------------------------------------------
+Collinear checks to see if the three points given are collinear,
+by checking to see if each element of the cross product is zero.
+---------------------------------------------------------------------*/
+static
+BOOL Collinear( LPVERTEX a, LPVERTEX b, LPVERTEX c )
+{
+ return
+ ( c->v[Z] - a->v[Z] ) * ( b->v[Y] - a->v[Y] ) -
+ ( b->v[Z] - a->v[Z] ) * ( c->v[Y] - a->v[Y] ) == 0
+ && ( b->v[Z] - a->v[Z] ) * ( c->v[X] - a->v[X] ) -
+ ( b->v[X] - a->v[X] ) * ( c->v[Z] - a->v[Z] ) == 0
+ && ( b->v[X] - a->v[X] ) * ( c->v[Y] - a->v[Y] ) -
+ ( b->v[Y] - a->v[Y] ) * ( c->v[X] - a->v[X] ) == 0 ;
+}
+
+/*---------------------------------------------------------------------
+VolumeSign returns the sign of the volume of the tetrahedron determined by f
+and p. VolumeSign is +1 iff p is on the negative side of f,
+where the positive side is determined by the rh-rule. So the volume
+is positive if the ccw normal to f points outside the tetrahedron.
+The final fewer-multiplications form is due to Bob Williamson.
+---------------------------------------------------------------------*/
+int VolumeSign( LPFACE f, LPVERTEX p )
+{
+ double vol;
+ double ax, ay, az, bx, by, bz, cx, cy, cz;
+
+ ax = f->Vertex[0]->v[X] - p->v[X];
+ ay = f->Vertex[0]->v[Y] - p->v[Y];
+ az = f->Vertex[0]->v[Z] - p->v[Z];
+ bx = f->Vertex[1]->v[X] - p->v[X];
+ by = f->Vertex[1]->v[Y] - p->v[Y];
+ bz = f->Vertex[1]->v[Z] - p->v[Z];
+ cx = f->Vertex[2]->v[X] - p->v[X];
+ cy = f->Vertex[2]->v[Y] - p->v[Y];
+ cz = f->Vertex[2]->v[Z] - p->v[Z];
+
+ vol = ax * (by*cz - bz*cy)
+ + ay * (bz*cx - bx*cz)
+ + az * (bx*cy - by*cx);
+
+
+ /* The volume should be an integer. */
+ if ( vol > 0.5 ) return 1;
+ else if ( vol < -0.5 ) return -1;
+ else return 0;
+}
+
+
+
+/*---------------------------------------------------------------------
+CleanEdges runs through the Edge list and cleans up the structure.
+If there is a NewFace then it will put that face in place of the
+Visible face and NULL out NewFace. It also deletes so marked edges.
+---------------------------------------------------------------------*/
+static
+void CleanEdges(LPHULL hull)
+{
+ LPEDGE e; /* Primary index into Edge list. */
+ LPEDGE t; /* Temporary Edge pointer. */
+
+ /* Integrate the NewFace's into the data structure. */
+ /* Check every Edge. */
+
+ e = hull ->edges;
+ do {
+ if ( e->NewFace ) {
+
+ if ( e->AdjFace[0]->Visible )
+ e->AdjFace[0] = e->NewFace;
+ else
+ e->AdjFace[1] = e->NewFace;
+
+ e->NewFace = NULL;
+ }
+
+ e = e->Next;
+
+ } while ( e != hull ->edges );
+
+ /* Delete any edges marked for deletion. */
+ while ( hull ->edges && hull ->edges->DoDelete ) {
+
+ e = hull ->edges;
+
+ XDELETE( hull ->edges, e );
+ }
+
+ e = hull ->edges->Next;
+
+ do {
+ if ( e->DoDelete ) {
+
+ t = e;
+ e = e->Next;
+ XDELETE( hull ->edges, t );
+ }
+ else e = e->Next;
+
+ } while ( e != hull ->edges );
+}
+
+/*---------------------------------------------------------------------
+CleanFaces runs through the face list and deletes any face marked Visible.
+---------------------------------------------------------------------*/
+static
+void CleanFaces(LPHULL hull)
+{
+ LPFACE f; /* Primary pointer into face list. */
+ LPFACE t; /* Temporary pointer, for deleting. */
+
+
+ while ( hull ->faces && hull ->faces->Visible ) {
+
+ f = hull ->faces;
+ XDELETE( hull ->faces, f );
+ }
+
+ f = hull ->faces->Next;
+
+ do {
+ if ( f->Visible ) {
+
+ t = f;
+ f = f->Next;
+ XDELETE( hull ->faces, t );
+ }
+ else f = f->Next;
+
+ } while ( f != hull ->faces );
+}
+
+
+
+/*---------------------------------------------------------------------
+CleanVertices runs through the Vertex list and deletes the
+vertices that are marked as processed but are not incident to any
+undeleted edges.
+---------------------------------------------------------------------*/
+static
+void CleanVertices(LPHULL hull)
+{
+ LPEDGE e;
+ LPVERTEX v, t;
+
+ /* Mark all vertices incident to some undeleted Edge as on the hull. */
+
+ e = hull ->edges;
+ do {
+ e->EndPts[0]->onhull = e->EndPts[1]->onhull = ONHULL;
+ e = e->Next;
+
+ } while (e != hull ->edges);
+
+
+ /* Delete all vertices that have been processed but
+ are not on the hull. */
+
+ while ( hull ->vertices && hull->vertices->mark && !hull ->vertices->onhull ) {
+
+ v = hull ->vertices;
+ XDELETE(hull ->vertices, v );
+ }
+
+
+ v = hull ->vertices->Next;
+ do {
+ if (v->mark && !v->onhull ) {
+ t = v;
+ v = v->Next;
+ XDELETE(hull ->vertices, t )
+ }
+ else
+ v = v->Next;
+
+ } while ( v != hull ->vertices );
+
+
+ /* Reset flags. */
+
+ v = hull ->vertices;
+ do {
+ v->duplicate = NULL;
+ v->onhull = !ONHULL;
+ v = v->Next;
+
+ } while (v != hull->vertices );
+
+}
+
+
+
+
+/*---------------------------------------------------------------------
+MakeCcw puts the vertices in the face structure in counterclock wise
+order. We want to store the vertices in the same
+order as in the Visible face. The third Vertex is always p.
+
+Although no specific ordering of the edges of a face are used
+by the code, the following condition is maintained for each face f:
+one of the two endpoints of f->Edge[i] matches f->Vertex[i].
+But note that this does not imply that f->Edge[i] is between
+f->Vertex[i] and f->Vertex[(i+1)%3]. (Thanks to Bob Williamson.)
+---------------------------------------------------------------------*/
+
+static
+void MakeCcw(LPFACE f, LPEDGE e, LPVERTEX p)
+{
+ LPFACE fv; /* The Visible face adjacent to e */
+ int i; /* Index of e->endpoint[0] in fv. */
+ LPEDGE s; /* Temporary, for swapping */
+
+ if (e->AdjFace[0]->Visible)
+
+ fv = e->AdjFace[0];
+ else
+ fv = e->AdjFace[1];
+
+ /* Set Vertex[0] & [1] of f to have the same orientation
+ as do the corresponding vertices of fv. */
+
+ for ( i=0; fv->Vertex[i] != e->EndPts[0]; ++i )
+ ;
+
+ /* Orient f the same as fv. */
+
+ if ( fv->Vertex[ (i+1) % 3 ] != e->EndPts[1] ) {
+
+ f->Vertex[0] = e->EndPts[1];
+ f->Vertex[1] = e->EndPts[0];
+ }
+ else {
+ f->Vertex[0] = e->EndPts[0];
+ f->Vertex[1] = e->EndPts[1];
+ SWAP( s, f->Edge[1], f->Edge[2] );
+ }
+
+ /* This swap is tricky. e is Edge[0]. Edge[1] is based on endpt[0],
+ Edge[2] on endpt[1]. So if e is oriented "forwards," we
+ need to move Edge[1] to follow [0], because it precedes. */
+
+ f->Vertex[2] = p;
+}
+
+/*---------------------------------------------------------------------
+MakeConeFace makes a new face and two new edges between the
+Edge and the point that are passed to it. It returns a pointer to
+the new face.
+---------------------------------------------------------------------*/
+
+static
+LPFACE MakeConeFace(LPHULL hull, LPEDGE e, LPVERTEX p)
+{
+ LPEDGE new_edge[2];
+ LPFACE new_face;
+ int i, j;
+
+ /* Make two new edges (if don't already exist). */
+
+ for ( i=0; i < 2; ++i )
+ /* If the Edge exists, copy it into new_edge. */
+ if ( !( new_edge[i] = e->EndPts[i]->duplicate) ) {
+
+ /* Otherwise (duplicate is NULL), MakeNullEdge. */
+ new_edge[i] = MakeNullEdge(hull);
+ new_edge[i]->EndPts[0] = e->EndPts[i];
+ new_edge[i]->EndPts[1] = p;
+ e->EndPts[i]->duplicate = new_edge[i];
+ }
+
+ /* Make the new face. */
+ new_face = MakeNullFace(hull);
+ new_face->Edge[0] = e;
+ new_face->Edge[1] = new_edge[0];
+ new_face->Edge[2] = new_edge[1];
+ MakeCcw( new_face, e, p );
+
+ /* Set the adjacent face pointers. */
+ for ( i=0; i < 2; ++i )
+ for ( j=0; j < 2; ++j )
+ /* Only one NULL link should be set to new_face. */
+ if ( !new_edge[i]->AdjFace[j] ) {
+ new_edge[i]->AdjFace[j] = new_face;
+ break;
+ }
+
+ return new_face;
+}
+
+
+/*---------------------------------------------------------------------
+AddOne is passed a Vertex. It first determines all faces Visible from
+that point. If none are Visible then the point is marked as not
+onhull. Next is a loop over edges. If both faces adjacent to an Edge
+are Visible, then the Edge is marked for deletion. If just one of the
+adjacent faces is Visible then a new face is constructed.
+---------------------------------------------------------------------*/
+static
+BOOL AddOne(LPHULL hull, LPVERTEX p)
+{
+ LPFACE f;
+ LPEDGE e, temp;
+ int vol;
+ BOOL vis = false;
+
+
+ /* Mark faces Visible from p. */
+ f = hull -> faces;
+
+ do {
+
+ vol = VolumeSign(f, p);
+
+ if ( vol < 0 ) {
+ f->Visible = VISIBLE;
+ vis = true;
+ }
+
+ f = f->Next;
+
+ } while ( f != hull ->faces );
+
+ /* If no faces are Visible from p, then p is inside the hull. */
+
+ if ( !vis ) {
+
+ p->onhull = !ONHULL;
+ return false;
+ }
+
+ /* Mark edges in interior of Visible region for deletion.
+ Erect a NewFace based on each border Edge. */
+
+ e = hull ->edges;
+
+ do {
+
+ temp = e->Next;
+
+ if ( e->AdjFace[0]->Visible && e->AdjFace[1]->Visible )
+ /* e interior: mark for deletion. */
+ e->DoDelete = REMOVED;
+
+ else
+ if ( e->AdjFace[0]->Visible || e->AdjFace[1]->Visible )
+ /* e border: make a new face. */
+ e->NewFace = MakeConeFace(hull, e, p );
+
+ e = temp;
+
+ } while ( e != hull ->edges );
+
+ return true;
+}
+
+
+/*---------------------------------------------------------------------
+ DoubleTriangle builds the initial double triangle. It first finds 3
+ noncollinear points and makes two faces out of them, in opposite order.
+ It then finds a fourth point that is not coplanar with that face. The
+ vertices are stored in the face structure in counterclockwise order so
+ that the volume between the face and the point is negative. Lastly, the
+ 3 newfaces to the fourth point are constructed and the data structures
+ are cleaned up.
+---------------------------------------------------------------------*/
+
+static
+BOOL DoubleTriangle(LPHULL hull)
+{
+ LPVERTEX v0, v1, v2, v3;
+ LPFACE f0, f1 = NULL;
+ int vol;
+
+ /* Find 3 noncollinear points. */
+ v0 = hull ->vertices;
+ while ( Collinear( v0, v0->Next, v0->Next->Next ) )
+ if ( ( v0 = v0->Next ) == hull->vertices )
+ return false; /* All points are Collinear! */
+
+ v1 = v0->Next;
+ v2 = v1->Next;
+
+ /* Mark the vertices as processed. */
+ v0->mark = PROCESSED;
+ v1->mark = PROCESSED;
+ v2->mark = PROCESSED;
+
+ /* Create the two "twin" faces. */
+ f0 = MakeFace(hull, v0, v1, v2, f1 );
+ f1 = MakeFace(hull, v2, v1, v0, f0 );
+
+ /* Link adjacent face fields. */
+ f0->Edge[0]->AdjFace[1] = f1;
+ f0->Edge[1]->AdjFace[1] = f1;
+ f0->Edge[2]->AdjFace[1] = f1;
+ f1->Edge[0]->AdjFace[1] = f0;
+ f1->Edge[1]->AdjFace[1] = f0;
+ f1->Edge[2]->AdjFace[1] = f0;
+
+ /* Find a fourth, noncoplanar point to form tetrahedron. */
+ v3 = v2->Next;
+ vol = VolumeSign( f0, v3 );
+
+ while ( !vol ) {
+
+ if ( ( v3 = v3->Next ) == v0 )
+ return false; /* All points are coplanar! */
+
+ vol = VolumeSign( f0, v3 );
+ }
+
+ /* Insure that v3 will be the first added. */
+ hull ->vertices = v3;
+ return true;
+}
+
+
+
+/*---------------------------------------------------------------------
+ConstructHull adds the vertices to the hull one at a time. The hull
+vertices are those in the list marked as onhull.
+---------------------------------------------------------------------*/
+static
+void ConstructHull(LPHULL hull)
+{
+ LPVERTEX v, vnext;
+ BOOL changed; /* T if addition changes hull; not used. */
+
+ v = hull->vertices;
+
+ do {
+ vnext = v->Next;
+
+ changed = false;
+
+ if (!v->mark ) {
+
+ v->mark = PROCESSED;
+ changed = AddOne(hull, v );
+
+ CleanEdges(hull);
+ CleanFaces(hull);
+ CleanVertices(hull);
+ }
+
+ v = vnext;
+
+ } while (v != hull->vertices );
+
+}
+
+
+
+/*-------------------------------------------------------------------*/
+
+
+static
+void AddVec( VEC3I q, VEC3I ray )
+{
+ int i;
+
+ for( i = 0; i < DIM; i++ )
+ ray[i] = q[i] + ray[i];
+}
+
+/*---------------------------------------------------------------------
+a - b ==> c.
+---------------------------------------------------------------------*/
+static
+void SubVec( VEC3I a, VEC3I b, VEC3I c )
+{
+ int i;
+
+ for( i = 0; i < DIM; i++ )
+ c[i] = a[i] - b[i];
+}
+
+
+/*---------------------------------------------------------------------
+Returns the dot product of the two input vectors.
+---------------------------------------------------------------------*/
+static
+double Dot( VEC3I a, LPVEC3 b )
+{
+ int i;
+ double sum = 0.0;
+
+ for( i = 0; i < DIM; i++ )
+ sum += a[i] * b->n[i];
+
+ return sum;
+}
+
+/*---------------------------------------------------------------------
+Compute the cross product of (b-a)x(c-a) and place into N.
+---------------------------------------------------------------------*/
+static
+void NormalVec( VEC3I a, VEC3I b, VEC3I c, LPVEC3 N )
+{
+ N->n[X] = ( c[Z] - a[Z] ) * ( b[Y] - a[Y] ) -
+ ( b[Z] - a[Z] ) * ( c[Y] - a[Y] );
+ N->n[Y] = ( b[Z] - a[Z] ) * ( c[X] - a[X] ) -
+ ( b[X] - a[X] ) * ( c[Z] - a[Z] );
+ N->n[Z] = ( b[X] - a[X] ) * ( c[Y] - a[Y] ) -
+ ( b[Y] - a[Y] ) * ( c[X] - a[X] );
+}
+
+
+
+
+static
+int InBox( VEC3I q, VEC3I bmin, VEC3I bmax )
+{
+
+ if( ( bmin[X] <= q[X] ) && ( q[X] <= bmax[X] ) &&
+ ( bmin[Y] <= q[Y] ) && ( q[Y] <= bmax[Y] ) &&
+ ( bmin[Z] <= q[Z] ) && ( q[Z] <= bmax[Z] ) )
+ return true;
+
+ return false;
+}
+
+
+
+/*
+ This function returns a char:
+ '0': the segment [ab] does not intersect (completely misses) the
+ bounding box surrounding the n-th triangle T. It lies
+ strictly to one side of one of the six supporting planes.
+ '?': status unknown: the segment may or may not intersect T.
+*/
+
+static
+char BoxTest(LPHULL hull, int n, VEC3I a, VEC3I b)
+{
+ int i; /* Coordinate index */
+ int w;
+
+ for ( i=0; i < DIM; i++ ) {
+
+ w = hull ->Box[ n ][0][i]; /* min: lower left */
+
+ if ( (a[i] < w) && (b[i] < w) ) return '0';
+
+ w = hull ->Box[ n ][1][i]; /* max: upper right */
+
+ if ( (a[i] > w) && (b[i] > w) ) return '0';
+ }
+
+ return '?';
+}
+
+
+
+/* Return a random ray endpoint */
+
+static
+void RandomRay( VEC3I ray, int radius )
+{
+ double x, y, z, w, t;
+
+ /* Generate a random point on a sphere of radius 1. */
+ /* the sphere is sliced at z, and a random point at angle t
+ generated on the circle of intersection. */
+
+ z = 2.0 * (double) rand() / RAND_MAX - 1.0;
+ t = 2.0 * M_PI * (double) rand() / RAND_MAX;
+ w = sqrt( 1 - z*z );
+ x = w * cos( t );
+ y = w * sin( t );
+
+ ray[X] = (int) ( radius * x );
+ ray[Y] = (int) ( radius * y );
+ ray[Z] = (int) ( radius * z );
+}
+
+
+
+static
+int ComputeBox(LPHULL hull, int F, VEC3I bmin, VEC3I bmax )
+{
+ int i, j;
+ double radius;
+
+ for( i = 0; i < F; i++ )
+ for( j = 0; j < DIM; j++ ) {
+
+ if( hull ->Vertices[i][j] < bmin[j] )
+ bmin[j] = hull ->Vertices[i][j];
+
+ if( hull ->Vertices[i][j] > bmax[j] )
+ bmax[j] = hull ->Vertices[i][j];
+ }
+
+ radius = sqrt( pow( (double)(bmax[X] - bmin[X]), 2.0 ) +
+ pow( (double)(bmax[Y] - bmin[Y]), 2.0 ) +
+ pow( (double)(bmax[Z] - bmin[Z]), 2.0 ) );
+
+ return (int)( radius +1 ) + 1;
+}
+
+
+/*---------------------------------------------------------------------
+Computes N & D and returns index m of largest component.
+---------------------------------------------------------------------*/
+static
+int PlaneCoeff(LPHULL hull, VEC3I T, LPVEC3 N, double *D )
+{
+ int i;
+ double t; /* Temp storage */
+ double biggest = 0.0; /* Largest component of normal vector. */
+ int m = 0; /* Index of largest component. */
+
+ NormalVec(hull ->Vertices[T[0]], hull ->Vertices[T[1]], hull ->Vertices[T[2]], N );
+ *D = Dot( hull ->Vertices[T[0]], N );
+
+ /* Find the largest component of N. */
+ for ( i = 0; i < DIM; i++ ) {
+ t = fabs( N->n[i] );
+ if ( t > biggest ) {
+ biggest = t;
+ m = i;
+ }
+ }
+ return m;
+}
+
+/*---------------------------------------------------------------------
+ 'p': The segment lies wholly within the plane.
+ 'q': The q endpoint is on the plane (but not 'p').
+ 'r': The r endpoint is on the plane (but not 'p').
+ '0': The segment lies strictly to one side or the other of the plane.
+ '1': The segement intersects the plane, and 'p' does not hold.
+---------------------------------------------------------------------*/
+static
+char SegPlaneInt(LPHULL hull, VEC3I T, VEC3I q, VEC3I r, LPVEC3 p, int *m)
+{
+ VEC3 N; double D;
+ VEC3I rq;
+ double num, denom, t;
+ int i;
+
+ *m = PlaneCoeff(hull, T, &N, &D );
+ num = D - Dot( q, &N );
+ SubVec( r, q, rq );
+ denom = Dot( rq, &N );
+
+ if ( denom == 0.0 ) { /* Segment is parallel to plane. */
+ if ( num == 0.0 ) /* q is on plane. */
+ return 'p';
+ else
+ return '0';
+ }
+ else
+ t = num / denom;
+
+ for( i = 0; i < DIM; i++ )
+ p->n[i] = q[i] + t * ( r[i] - q[i] );
+
+ if ( (0.0 < t) && (t < 1.0) )
+ return '1';
+ else if ( num == 0.0 ) /* t == 0 */
+ return 'q';
+ else if ( num == denom ) /* t == 1 */
+ return 'r';
+ else return '0';
+}
+
+
+
+static
+int AreaSign( VEC3I a, VEC3I b, VEC3I c )
+{
+ double area2;
+
+ area2 = ( b[0] - a[0] ) * (double)( c[1] - a[1] ) -
+ ( c[0] - a[0] ) * (double)( b[1] - a[1] );
+
+ /* The area should be an integer. */
+ if ( area2 > 0.5 ) return 1;
+ else if ( area2 < -0.5 ) return -1;
+ else return 0;
+}
+
+
+static
+char InTri2D( VEC3I Tp[3], VEC3I pp )
+{
+ int area0, area1, area2;
+
+ /* compute three AreaSign() values for pp w.r.t. each Edge of the face in 2D */
+ area0 = AreaSign( pp, Tp[0], Tp[1] );
+ area1 = AreaSign( pp, Tp[1], Tp[2] );
+ area2 = AreaSign( pp, Tp[2], Tp[0] );
+
+ if ( (( area0 == 0 ) && ( area1 > 0 ) && ( area2 > 0 )) ||
+ (( area1 == 0 ) && ( area0 > 0 ) && ( area2 > 0 )) ||
+ (( area2 == 0 ) && ( area0 > 0 ) && ( area1 > 0 )) )
+ return 'E';
+
+ if ( (( area0 == 0 ) && ( area1 < 0 ) && ( area2 < 0 )) ||
+ (( area1 == 0 ) && ( area0 < 0 ) && ( area2 < 0 )) ||
+ (( area2 == 0 ) && ( area0 < 0 ) && ( area1 < 0 )))
+ return 'E';
+
+ if ( (( area0 > 0 ) && ( area1 > 0 ) && ( area2 > 0 )) ||
+ (( area0 < 0 ) && ( area1 < 0 ) && ( area2 < 0 )))
+ return 'F';
+
+ if ( ( area0 == 0 ) && ( area1 == 0 ) && ( area2 == 0 ) )
+ return '?'; /* Error in InTriD */
+
+ if ( (( area0 == 0 ) && ( area1 == 0 )) ||
+ (( area0 == 0 ) && ( area2 == 0 )) ||
+ (( area1 == 0 ) && ( area2 == 0 )) )
+ return 'V';
+
+ else
+ return '0';
+}
+
+/* Assumption: p lies in the plane containing T.
+ Returns a char:
+ 'V': the query point p coincides with a Vertex of triangle T.
+ 'E': the query point p is in the relative interior of an Edge of triangle T.
+ 'F': the query point p is in the relative interior of a Face of triangle T.
+ '0': the query point p does not intersect (misses) triangle T.
+*/
+
+static
+char InTri3D(LPHULL hull, VEC3I T, int m, VEC3I p )
+{
+ int i; /* Index for X,Y,Z */
+ int j; /* Index for X,Y */
+ int k; /* Index for triangle Vertex */
+ VEC3I pp; /* projected p */
+ VEC3I Tp[3]; /* projected T: three new vertices */
+
+ /* Project out coordinate m in both p and the triangular face */
+ j = 0;
+ for ( i = 0; i < DIM; i++ ) {
+ if ( i != m ) { /* skip largest coordinate */
+ pp[j] = p[i];
+ for ( k = 0; k < 3; k++ )
+ Tp[k][j] = hull->Vertices[T[k]][i];
+ j++;
+ }
+ }
+ return( InTri2D( Tp, pp ) );
+}
+
+
+
+static
+int VolumeSign2( VEC3I a, VEC3I b, VEC3I c, VEC3I d )
+{
+ double vol;
+ double ax, ay, az, bx, by, bz, cx, cy, cz, dx, dy, dz;
+ double bxdx, bydy, bzdz, cxdx, cydy, czdz;
+
+ ax = a[X];
+ ay = a[Y];
+ az = a[Z];
+ bx = b[X];
+ by = b[Y];
+ bz = b[Z];
+ cx = c[X];
+ cy = c[Y];
+ cz = c[Z];
+ dx = d[X];
+ dy = d[Y];
+ dz = d[Z];
+
+ bxdx=bx-dx;
+ bydy=by-dy;
+ bzdz=bz-dz;
+ cxdx=cx-dx;
+ cydy=cy-dy;
+ czdz=cz-dz;
+ vol = (az-dz) * (bxdx*cydy - bydy*cxdx)
+ + (ay-dy) * (bzdz*cxdx - bxdx*czdz)
+ + (ax-dx) * (bydy*czdz - bzdz*cydy);
+
+
+ /* The volume should be an integer. */
+ if ( vol > 0.5 ) return 1;
+ else if ( vol < -0.5 ) return -1;
+ else return 0;
+}
+
+
+
+
+/*---------------------------------------------------------------------
+The signed volumes of three tetrahedra are computed, determined
+by the segment qr, and each Edge of the triangle.
+Returns a char:
+ 'v': the open segment includes a Vertex of T.
+ 'e': the open segment includes a point in the relative interior of an Edge
+ of T.
+ 'f': the open segment includes a point in the relative interior of a face
+ of T.
+ '0': the open segment does not intersect triangle T.
+---------------------------------------------------------------------*/
+
+static
+char SegTriCross(LPHULL hull, VEC3I T, VEC3I q, VEC3I r )
+{
+ int vol0, vol1, vol2;
+
+ vol0 = VolumeSign2( q, hull->Vertices[ T[0] ], hull->Vertices[ T[1] ], r );
+ vol1 = VolumeSign2( q, hull->Vertices[ T[1] ], hull->Vertices[ T[2] ], r );
+ vol2 = VolumeSign2( q, hull->Vertices[ T[2] ], hull->Vertices[ T[0] ], r );
+
+
+ /* Same sign: segment intersects interior of triangle. */
+ if ( ( ( vol0 > 0 ) && ( vol1 > 0 ) && ( vol2 > 0 ) ) ||
+ ( ( vol0 < 0 ) && ( vol1 < 0 ) && ( vol2 < 0 ) ) )
+ return 'f';
+
+ /* Opposite sign: no intersection between segment and triangle */
+ if ( ( ( vol0 > 0 ) || ( vol1 > 0 ) || ( vol2 > 0 ) ) &&
+ ( ( vol0 < 0 ) || ( vol1 < 0 ) || ( vol2 < 0 ) ) )
+ return '0';
+
+ else if ( ( vol0 == 0 ) && ( vol1 == 0 ) && ( vol2 == 0 ) )
+ return '?'; /* Error 1 in SegTriCross */
+
+ /* Two zeros: segment intersects Vertex. */
+ else if ( ( ( vol0 == 0 ) && ( vol1 == 0 ) ) ||
+ ( ( vol0 == 0 ) && ( vol2 == 0 ) ) ||
+ ( ( vol1 == 0 ) && ( vol2 == 0 ) ) )
+ return 'v';
+
+ /* One zero: segment intersects Edge. */
+ else if ( ( vol0 == 0 ) || ( vol1 == 0 ) || ( vol2 == 0 ) )
+ return 'e';
+
+ else
+ return '?'; /* Error 2 in SegTriCross */
+}
+
+
+
+static
+char SegTriInt(LPHULL hull, VEC3I T, VEC3I q, VEC3I r, LPVEC3 p )
+{
+ int code;
+ int m = -1;
+
+ code = SegPlaneInt(hull, T, q, r, p, &m );
+
+ if ( code == '0') return '0';
+ else if ( code == 'q') return InTri3D(hull, T, m, q );
+ else if ( code == 'r') return InTri3D(hull, T, m, r );
+ else if ( code == 'p') return 'p';
+ else if ( code == '1' ) return SegTriCross(hull, T, q, r );
+ else
+ return code; /* Error */
+}
+
+
+
+
+/*
+ This function returns a char:
+ 'i': the query point a is strictly interior to polyhedron P.
+ 'o': the query point a is strictly exterior to( or outside of) polyhedron P.
+*/
+char InPolyhedron(LPHULL hull, VEC3I q)
+{
+ int F = hull->nfaces;
+ VEC3I Ray; /* Ray endpoint. */
+ VEC3 p; /* Intersection point; not used. */
+ int f, k = 0, crossings = 0;
+ char code = '?';
+
+
+ /* If query point is outside bounding box, finished. */
+ if ( !InBox( q, hull->bmin, hull->bmax ) )
+ return 'o';
+
+ LOOP:
+ while( k++ < F ) {
+
+ crossings = 0;
+
+ RandomRay(Ray, hull->radius );
+
+ AddVec( q, Ray );
+
+
+ for ( f = 0; f < F; f++ ) { /* Begin check each face */
+
+ if ( BoxTest(hull, f, q, Ray ) == '0' ) {
+ code = '0';
+
+ }
+ else code = SegTriInt(hull, hull->Faces[f], q, Ray, &p );
+
+
+ /* If ray is degenerate, then goto outer while to generate another. */
+ if ( code == 'p' || code == 'v' || code == 'e' ) {
+
+ goto LOOP;
+ }
+
+ /* If ray hits face at interior point, increment crossings. */
+ else if ( code == 'f' ) {
+ crossings++;
+
+ }
+
+ /* If query endpoint q sits on a V/E/F, return inside. */
+ else if ( code == 'V' || code == 'E' || code == 'F' )
+ return code; /* 'i'; MM2 */
+
+ /* If ray misses triangle, do nothing. */
+ else if ( code == '0' )
+ ;
+
+ else
+ return '?'; /* Error */
+
+ }
+ /* No degeneracies encountered: ray is generic, so finished. */
+ break;
+
+ } /* End while loop */
+
+
+ /* q strictly interior to polyhedron iff an odd number of crossings. */
+ if( ( crossings % 2 ) == 1 )
+ return 'i';
+
+ else return 'o';
+}
+
+
+/*/ ---------------------------------------------------------------------------------- */
+
+
+
+
+static
+void StoreResults(LPHULL hull)
+{
+
+ int i, w;
+ LPVERTEX v;
+ LPFACE f;
+ int V = 0, F = 0;
+ int j, k;
+
+ /* Vertices */
+
+ v = hull ->vertices;
+ V = 0;
+ do {
+
+ v -> vnum = V;
+ hull ->Vertices[V][X] = v -> v[X];
+ hull ->Vertices[V][Y] = v -> v[Y];
+ hull ->Vertices[V][Z] = v -> v[Z];
+
+ v = v->Next;
+ V++;
+
+ } while ( v != hull ->vertices );
+
+ hull ->nvertex = V;
+
+ /* Faces */
+ f = hull ->faces;
+ F = 0;
+ do {
+
+ hull ->Faces[F][0] = f->Vertex[0]->vnum;
+ hull ->Faces[F][1] = f->Vertex[1]->vnum;
+ hull ->Faces[F][2] = f->Vertex[2]->vnum;
+
+ for ( j=0; j < 3; j++ ) {
+
+ hull ->Box[F][0][j] = hull ->Vertices[ hull ->Faces[F][0] ][j];
+ hull ->Box[F][1][j] = hull ->Vertices[ hull ->Faces[F][0] ][j];
+ }
+
+ /* Check k=1,2 vertices of face. */
+ for ( k=1; k < 3; k++ )
+ for ( j=0; j < 3; j++ ) {
+
+ w = hull ->Vertices[ hull ->Faces[F][k] ][j];
+ if ( w < hull ->Box[F][0][j] ) hull ->Box[F][0][j] = w;
+ if ( w > hull ->Box[F][1][j] ) hull ->Box[F][1][j] = w;
+ }
+
+
+ f = f->Next; F++;
+
+ } while ( f != hull ->faces );
+
+
+ hull ->nfaces = F;
+
+
+ /* Initialize the bounding box */
+ for ( i = 0; i < DIM; i++ )
+ hull ->bmin[i] = hull ->bmax[i] = hull ->Vertices[0][i];
+
+ hull ->radius = ComputeBox(hull, V, hull ->bmin, hull ->bmax );
+
+
+}
+
+
+LCMSHANDLE cmsxHullInit(void)
+{
+ LPHULL hull = (LPHULL) malloc(sizeof(HULL));
+
+ ZeroMemory(hull, sizeof(HULL));
+
+ hull->vnumCounter = 0;
+ hull->vertices = NULL;
+ hull->edges = NULL;
+ hull->faces = NULL;
+ hull->nfaces = 0;
+ hull->nvertex = 0;
+
+ return (LCMSHANDLE) (LPSTR) hull;
+}
+
+
+void cmsxHullDone(LCMSHANDLE hHull)
+{
+ LPHULL hull = (LPHULL) (LPSTR) hHull;
+
+ if (hull)
+ free((LPVOID) hull);
+}
+
+
+BOOL cmsxHullAddPoint(LCMSHANDLE hHull, int x, int y, int z)
+{
+ LPVERTEX v;
+ LPHULL hull = (LPHULL) (LPSTR) hHull;
+
+
+ v = MakeNullVertex(hull);
+ v->v[X] = x;
+ v->v[Y] = y;
+ v->v[Z] = z;
+ v->vnum = hull->vnumCounter++;
+
+ return true;
+}
+
+BOOL cmsxHullComputeHull(LCMSHANDLE hHull)
+{
+
+ LPHULL hull = (LPHULL) (LPSTR) hHull;
+
+ if (!DoubleTriangle(hull)) return false;
+
+ ConstructHull(hull);
+ StoreResults(hull);
+
+ return true;
+}
+
+
+char cmsxHullCheckpoint(LCMSHANDLE hHull, int x, int y, int z)
+{
+ VEC3I q;
+ LPHULL hull = (LPHULL) (LPSTR) hHull;
+
+ q[X] = x; q[Y] = y; q[Z] = z;
+
+ return InPolyhedron(hull, q ) ;
+}
+
+
+BOOL cmsxHullDumpVRML(LCMSHANDLE hHull, const char* fname)
+{
+ FILE* fp;
+ int i;
+ LPHULL hull = (LPHULL) (LPSTR) hHull;
+
+ fp = fopen (fname, "wt");
+ if (fp == NULL)
+ return false;
+
+ fprintf (fp, "#VRML V2.0 utf8\n");
+
+ /* set the viewing orientation and distance */
+ fprintf (fp, "DEF CamTest Group {\n");
+ fprintf (fp, "\tchildren [\n");
+ fprintf (fp, "\t\tDEF Cameras Group {\n");
+ fprintf (fp, "\t\t\tchildren [\n");
+ fprintf (fp, "\t\t\t\tDEF DefaultView Viewpoint {\n");
+ fprintf (fp, "\t\t\t\t\tposition 0 0 340\n");
+ fprintf (fp, "\t\t\t\t\torientation 0 0 1 0\n");
+ fprintf (fp, "\t\t\t\t\tdescription \"default view\"\n");
+ fprintf (fp, "\t\t\t\t}\n");
+ fprintf (fp, "\t\t\t]\n");
+ fprintf (fp, "\t\t},\n");
+ fprintf (fp, "\t]\n");
+ fprintf (fp, "}\n");
+
+ /* Output the background stuff */
+ fprintf (fp, "Background {\n");
+ fprintf (fp, "\tskyColor [\n");
+ fprintf (fp, "\t\t.5 .5 .5\n");
+ fprintf (fp, "\t]\n");
+ fprintf (fp, "}\n");
+
+ /* Output the shape stuff */
+ fprintf (fp, "Transform {\n");
+ fprintf (fp, "\tscale 8 8 8\n");
+ fprintf (fp, "\tchildren [\n");
+
+ /* Draw the axes as a shape: */
+ fprintf (fp, "\t\tShape {\n");
+ fprintf (fp, "\t\t\tappearance Appearance {\n");
+ fprintf (fp, "\t\t\t\tmaterial Material {\n");
+ fprintf (fp, "\t\t\t\t\tdiffuseColor 0 0.8 0\n");
+ fprintf (fp, "\t\t\t\t\temissiveColor 1.0 1.0 1.0\n");
+ fprintf (fp, "\t\t\t\t\tshininess 0.8\n");
+ fprintf (fp, "\t\t\t\t}\n");
+ fprintf (fp, "\t\t\t}\n");
+ fprintf (fp, "\t\t\tgeometry IndexedLineSet {\n");
+ fprintf (fp, "\t\t\t\tcoord Coordinate {\n");
+ fprintf (fp, "\t\t\t\t\tpoint [\n");
+ fprintf (fp, "\t\t\t\t\t0.0 0.0 0.0,\n");
+ fprintf (fp, "\t\t\t\t\t%f 0.0 0.0,\n", 255.0);
+ fprintf (fp, "\t\t\t\t\t0.0 %f 0.0,\n", 255.0);
+ fprintf (fp, "\t\t\t\t\t0.0 0.0 %f]\n", 255.0);
+ fprintf (fp, "\t\t\t\t}\n");
+ fprintf (fp, "\t\t\t\tcoordIndex [\n");
+ fprintf (fp, "\t\t\t\t\t0, 1, -1\n");
+ fprintf (fp, "\t\t\t\t\t0, 2, -1\n");
+ fprintf (fp, "\t\t\t\t\t0, 3, -1]\n");
+ fprintf (fp, "\t\t\t}\n");
+ fprintf (fp, "\t\t}\n");
+
+
+ /* Draw the triangles as a shape: */
+ fprintf (fp, "\t\tShape {\n");
+ fprintf (fp, "\t\t\tappearance Appearance {\n");
+ fprintf (fp, "\t\t\t\tmaterial Material {\n");
+ fprintf (fp, "\t\t\t\t\tdiffuseColor 0 0.8 0\n");
+ fprintf (fp, "\t\t\t\t\temissiveColor 0 0 0\n");
+ fprintf (fp, "\t\t\t\t\tshininess 0.8\n");
+ fprintf (fp, "\t\t\t\t}\n");
+ fprintf (fp, "\t\t\t}\n");
+ fprintf (fp, "\t\t\tgeometry IndexedFaceSet {\n");
+ fprintf (fp, "\t\t\t\tsolid false\n");
+
+ /* fill in the points here */
+ fprintf (fp, "\t\t\t\tcoord Coordinate {\n");
+ fprintf (fp, "\t\t\t\t\tpoint [\n");
+
+ for (i = 0; i < hull->nvertex; ++i)
+ {
+ fprintf (fp, "\t\t\t\t\t%g %g %g%c\n",
+ (double) hull->Vertices[i][X], (double) hull->Vertices[i][Y], (double) hull->Vertices[i][Z],
+ i == hull->nvertex-1? ']': ',');
+ }
+ fprintf (fp, "\t\t\t\t}\n");
+
+ /* fill in the Vertex indices (followed by -1) */
+
+
+ fprintf (fp, "\t\t\t\tcoordIndex [\n");
+ for (i = 0; i < hull->nfaces; ++i)
+ {
+ fprintf (fp, "\t\t\t\t\t%d, %d, %d, -1\n",
+ hull->Faces[i][0], hull->Faces[i][1], hull->Faces[i][2]);
+
+ }
+ fprintf (fp, "]\n");
+
+
+ /* fill in the face colors */
+ fprintf (fp, "\t\t\t\tcolor Color {\n");
+ fprintf (fp, "\t\t\t\t\tcolor [\n");
+ for (i = 0; i < hull->nfaces; ++i)
+ {
+ int vx, vy, vz;
+ double r, g, b;
+
+ vx = hull->Faces[i][0]; vy = hull->Faces[i][1]; vz = hull->Faces[i][2];
+ r = (double) (hull->Vertices[vx][X] + hull->Vertices[vy][X] + hull->Vertices[vz][X]) / (3* 255);
+ g = (double) (hull->Vertices[vx][Y] + hull->Vertices[vy][Y] + hull->Vertices[vz][Y]) / (3* 255);
+ b = (double) (hull->Vertices[vx][Z] + hull->Vertices[vy][Z] + hull->Vertices[vz][Z]) / (3* 255);
+
+ fprintf (fp, "\t\t\t\t\t%g %g %g%c\n", r, g, b,
+ i == hull->nfaces-1? ']': ',');
+ }
+ fprintf (fp, "\t\t\t}\n");
+
+ fprintf (fp, "\t\t\tcolorPerVertex false\n");
+
+ fprintf (fp, "\t\t\t}\n");
+ fprintf (fp, "\t\t}\n");
+ fprintf (fp, "\t]\n");
+ fprintf (fp, "}\n");
+
+ fclose (fp);
+
+ return true;
+}
generated by cgit v1.2.1 (git 2.18.0) at 2024-12-22 09:15:03 +0000