← Back to branch summary

~ubuntu-branches/ubuntu/saucy/digikam/saucy

« back to all changes in this revision

Viewing changes to libs/3rdparty/lprof/cmshull.cpp

  • Committer: Package Import Robot
  • Author(s): Felix Geyer, Rohan Garg, Philip Muškovac, Felix Geyer
  • Date: 2011-09-23 18:18:55 UTC
  • mfrom: (1.2.36 upstream)
  • Revision ID: package-import@ubuntu.com-20110923181855-ifs67wxkugshev9k
Tags: 2:2.1.1-0ubuntu1
https://launchpad.net/bugs/834190
https://launchpad.net/bugs/804894
[ Rohan Garg ]
* New upstream release (LP: #834190)
  - debian/control
    + Build with libqtwebkit-dev
 - debian/kipi-plugins-common
    + Install libkvkontakte required by kipi-plugins
 - debian/digikam
    + Install panoramagui

[ Philip Muškovac ]
* New upstream release
  - debian/control:
    + Add libcv-dev, libcvaux-dev, libhighgui-dev, libboost-graph1.46-dev,
      libksane-dev, libxml2-dev, libxslt-dev, libqt4-opengl-dev, libqjson-dev,
      libgpod-dev and libqca2-dev to build-deps
    + Add packages for kipi-plugins, libmediawiki, libkface, libkgeomap and
      libkvkontakte
  - debian/rules:
    + Don't build with gphoto2 since it doesn't build with it.
  - Add kubuntu_fix_test_linking.diff to fix linking of the dngconverter test
  - update install files
  - update kubuntu_01_mysqld_executable_name.diff for new cmake layout
    and rename to kubuntu_mysqld_executable_name.diff
* Fix typo in digikam-data description (LP: #804894)
* Fix Vcs links

[ Felix Geyer ]
* Move library data files to the new packages libkface-data, libkgeomap-data
  and libkvkontakte-data.
* Override version of the embedded library packages to 1.0~digikam<version>.
* Exclude the library packages from digikam-dbg to prevent file conflicts in
  the future.
* Call dh_install with --list-missing.

Show diffs side-by-side

added added

removed removed

Lines of Context:
libs/3rdparty/lprof/cmshull.cpp
1
 
/* */
2
 
/*  Little cms - profiler construction set */
3
 
/*  Copyright (C) 1998-2001 Marti Maria <marti@littlecms.com> */
4
 
/* */
5
 
/* THIS SOFTWARE IS PROVIDED "AS-IS" AND WITHOUT WARRANTY OF ANY KIND, */
6
 
/* EXPRESS, IMPLIED OR OTHERWISE, INCLUDING WITHOUT LIMITATION, ANY */
7
 
/* WARRANTY OF MERCHANTABILITY OR FITNESS FOR A PARTICULAR PURPOSE. */
8
 
/* */
9
 
/* IN NO EVENT SHALL MARTI MARIA BE LIABLE FOR ANY SPECIAL, INCIDENTAL, */
10
 
/* INDIRECT OR CONSEQUENTIAL DAMAGES OF ANY KIND, */
11
 
/* OR ANY DAMAGES WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR PROFITS, */
12
 
/* WHETHER OR NOT ADVISED OF THE POSSIBILITY OF DAMAGE, AND ON ANY THEORY OF */
13
 
/* LIABILITY, ARISING OUT OF OR IN CONNECTION WITH THE USE OR PERFORMANCE */
14
 
/* OF THIS SOFTWARE. */
15
 
/* */
16
 
/* This file is free software; you can redistribute it and/or modify it */
17
 
/* under the terms of the GNU General Public License as published by */
18
 
/* the Free Software Foundation; either version 2 of the License, or */
19
 
/* (at your option) any later version. */
20
 
/* */
21
 
/* This program is distributed in the hope that it will be useful, but */
22
 
/* WITHOUT ANY WARRANTY; without even the implied warranty of */
23
 
/* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU */
24
 
/* General Public License for more details. */
25
 
/* */
26
 
/* You should have received a copy of the GNU General Public License */
27
 
/* along with this program; if not, write to the Free Software */
28
 
/* Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA. */
29
 
/* */
30
 
/* As a special exception to the GNU General Public License, if you */
31
 
/* distribute this file as part of a program that contains a */
32
 
/* configuration script generated by Autoconf, you may include it under */
33
 
/* the same distribution terms that you use for the rest of that program. */
34
 
/* */
35
 
/* Version 1.09a */
36
 
 
37
 
 
38
 
#include "lcmsprf.h"
39
 
 
40
 
/* Convex hull management */
41
 
 
42
 
LCMSHANDLE cdecl cmsxHullInit(void);
43
 
void       cdecl cmsxHullDone(LCMSHANDLE hHull);
44
 
BOOL       cdecl cmsxHullAddPoint(LCMSHANDLE hHull, int x, int y, int z);
45
 
BOOL       cdecl cmsxHullComputeHull(LCMSHANDLE hHull);
46
 
char       cdecl cmsxHullCheckpoint(LCMSHANDLE hHull, int x, int y, int z);
47
 
BOOL       cdecl cmsxHullDumpVRML(LCMSHANDLE hHull, const char* fname);
48
 
 
49
 
/* --------------------------------------------------------------------- */
50
 
 
51
 
 
52
 
 
53
 
/* This method is described in "Computational Geometry in C" Chapter 4.   */
54
 
/*  */
55
 
/* -------------------------------------------------------------------- */
56
 
/* This code is Copyright 1998 by Joseph O'Rourke.  It may be freely  */
57
 
/* redistributed in its entirety provided that this copyright notice is */
58
 
/* not removed. */
59
 
/* -------------------------------------------------------------------- */
60
 
 
61
 
#define SWAP(t,x,y)    { t = x; x = y; y = t; }
62
 
 
63
 
#define XFREE(p)
64
 
/* if (p) { free ((char *) p); p = NULL; } */
65
 
 
66
 
 
67
 
#define ADD( head, p )  if ( head )  { \
68
 
                p->Next = head; \
69
 
                p->Prev = head->Prev; \
70
 
                head->Prev = p; \
71
 
                p->Prev->Next = p; \
72
 
            } \
73
 
            else { \
74
 
                head = p; \
75
 
                head->Next = head->Prev = p; \
76
 
            }
77
 
 
78
 
#define XDELETE( head, p ) if ( head )  { \
79
 
                if ( head == head->Next ) \
80
 
                    head = NULL;  \
81
 
                else if ( p == head ) \
82
 
                    head = head->Next; \
83
 
                p->Next->Prev = p->Prev;  \
84
 
                p->Prev->Next = p->Next;  \
85
 
                XFREE( p ); \
86
 
            }
87
 
 
88
 
/* Define Vertex indices.  */
89
 
#define X   0
90
 
#define Y   1
91
 
#define Z   2
92
 
 
93
 
/* Define structures for vertices, edges and faces  */
94
 
 
95
 
typedef struct _vertex_struct VERTEX,FAR *LPVERTEX;
96
 
typedef struct _edge_struct   EDGE, FAR *LPEDGE;
97
 
typedef struct _face_struct   FACE, FAR *LPFACE;
98
 
 
99
 
 
100
 
struct _edge_struct {
101
 
 
102
 
   LPFACE    AdjFace[2];
103
 
   LPVERTEX  EndPts[2];
104
 
   LPFACE    NewFace;                /* pointer to incident cone face. */
105
 
   BOOL      DoDelete;               /* T iff Edge should be delete. */
106
 
 
107
 
   LPEDGE    Next, Prev;
108
 
};
109
 
 
110
 
struct _face_struct {
111
 
 
112
 
   LPEDGE    Edge[3];
113
 
   LPVERTEX  Vertex[3];
114
 
   BOOL      Visible;             /* T iff face Visible from new point. */
115
 
 
116
 
   LPFACE    Next, Prev;
117
 
};
118
 
 
119
 
struct _vertex_struct {
120
 
 
121
 
   int      v[3];
122
 
   int      vnum;
123
 
   LPEDGE   duplicate;            /* pointer to incident cone Edge (or NULL) */
124
 
   BOOL     onhull;               /* T iff point on hull. */
125
 
   BOOL     mark;                 /* T iff point already processed. */
126
 
 
127
 
   LPVERTEX  Next, Prev;
128
 
};
129
 
 
130
 
/* Define flags */
131
 
 
132
 
#define ONHULL       true
133
 
#define REMOVED      true
134
 
#define VISIBLE      true
135
 
#define PROCESSED    true
136
 
#define SAFE        1000000        /* Range of safe coord values. */
137
 
 
138
 
#define DIM 3                  /* Dimension of points */
139
 
typedef int    VEC3I[DIM];   /* Type integer point */
140
 
 
141
 
#define PMAX 10000     /* Max # of pts */
142
 
 
143
 
typedef struct {
144
 
 
145
 
    /* Global variable definitions */
146
 
    LPVERTEX vertices;
147
 
    LPEDGE edges;
148
 
    LPFACE faces;
149
 
 
150
 
    VEC3I Vertices[PMAX];        /* All the points */
151
 
    VEC3I Faces[PMAX];           /* Each triangle face is 3 indices */
152
 
    VEC3I Box[PMAX][2];          /* Box around each face */
153
 
 
154
 
 
155
 
    VEC3I bmin, bmax;
156
 
    int radius;
157
 
    int vnumCounter;
158
 
 
159
 
    int nfaces;
160
 
    int nvertex;
161
 
 
162
 
} HULL, FAR* LPHULL;
163
 
 
164
 
/* static HULL Global; */
165
 
 
166
 
 
167
 
/*---------------------------------------------------------------------
168
 
MakeNullVertex: Makes a Vertex, nulls out fields.
169
 
---------------------------------------------------------------------*/
170
 
 
171
 
static
172
 
LPVERTEX MakeNullVertex(LPHULL hull)
173
 
{
174
 
   LPVERTEX  v;
175
 
 
176
 
   v = (LPVERTEX) malloc(sizeof(VERTEX));
177
 
   if (!v) return NULL;
178
 
 
179
 
   v->duplicate = NULL;
180
 
   v->onhull = !ONHULL;
181
 
   v->mark = !PROCESSED;
182
 
   ADD( hull->vertices, v );
183
 
 
184
 
   return v;
185
 
}
186
 
 
187
 
 
188
 
 
189
 
/*---------------------------------------------------------------------
190
 
MakeNullEdge creates a new cell and initializes all pointers to NULL
191
 
and sets all flags to off.  It returns a pointer to the empty cell.
192
 
---------------------------------------------------------------------*/
193
 
static
194
 
LPEDGE MakeNullEdge(LPHULL hull)
195
 
{
196
 
   LPEDGE  e;
197
 
 
198
 
   e = (LPEDGE) malloc(sizeof(EDGE));
199
 
   if (!e) return NULL;
200
 
 
201
 
   e->AdjFace[0] = e->AdjFace[1] = e->NewFace = NULL;
202
 
   e->EndPts[0] = e->EndPts[1] = NULL;
203
 
   e->DoDelete = !REMOVED;
204
 
   ADD( hull->edges, e );
205
 
   return e;
206
 
}
207
 
 
208
 
/*--------------------------------------------------------------------
209
 
MakeNullFace creates a new face structure and initializes all of its
210
 
flags to NULL and sets all the flags to off.  It returns a pointer
211
 
to the empty cell.
212
 
---------------------------------------------------------------------*/
213
 
static
214
 
LPFACE MakeNullFace(LPHULL hull)
215
 
{
216
 
   LPFACE  f;
217
 
   int    i;
218
 
 
219
 
   f = (LPFACE) malloc(sizeof(FACE));
220
 
   if (!f) return NULL;
221
 
 
222
 
   for ( i=0; i < 3; ++i ) {
223
 
      f->Edge[i] = NULL;
224
 
      f->Vertex[i] = NULL;
225
 
   }
226
 
   f->Visible = !VISIBLE;
227
 
   ADD( hull->faces, f );
228
 
   return f;
229
 
}
230
 
 
231
 
 
232
 
 
233
 
/*---------------------------------------------------------------------
234
 
MakeFace creates a new face structure from three vertices (in ccw
235
 
order).  It returns a pointer to the face.
236
 
---------------------------------------------------------------------*/
237
 
static
238
 
LPFACE MakeFace(LPHULL hull, LPVERTEX v0, LPVERTEX v1, LPVERTEX v2, LPFACE fold)
239
 
{
240
 
   LPFACE  f;
241
 
   LPEDGE  e0, e1, e2;
242
 
 
243
 
   /* Create edges of the initial triangle. */
244
 
   if( !fold ) {
245
 
     e0 = MakeNullEdge(hull);
246
 
     e1 = MakeNullEdge(hull);
247
 
     e2 = MakeNullEdge(hull);
248
 
   }
249
 
   else { /* Copy from fold, in reverse order. */
250
 
     e0 = fold->Edge[2];
251
 
     e1 = fold->Edge[1];
252
 
     e2 = fold->Edge[0];
253
 
   }
254
 
   e0->EndPts[0] = v0; e0->EndPts[1] = v1;
255
 
   e1->EndPts[0] = v1; e1->EndPts[1] = v2;
256
 
   e2->EndPts[0] = v2; e2->EndPts[1] = v0;
257
 
 
258
 
   /* Create face for triangle. */
259
 
   f = MakeNullFace(hull);
260
 
   f->Edge[0]   = e0;  f->Edge[1]   = e1; f->Edge[2]   = e2;
261
 
   f->Vertex[0] = v0;  f->Vertex[1] = v1; f->Vertex[2] = v2;
262
 
 
263
 
   /* Link edges to face. */
264
 
   e0->AdjFace[0] = e1->AdjFace[0] = e2->AdjFace[0] = f;
265
 
 
266
 
   return f;
267
 
}
268
 
 
269
 
/*---------------------------------------------------------------------
270
 
Collinear checks to see if the three points given are collinear,
271
 
by checking to see if each element of the cross product is zero.
272
 
---------------------------------------------------------------------*/
273
 
static
274
 
BOOL Collinear( LPVERTEX a, LPVERTEX b, LPVERTEX c )
275
 
{
276
 
   return
277
 
         ( c->v[Z] - a->v[Z] ) * ( b->v[Y] - a->v[Y] ) -
278
 
         ( b->v[Z] - a->v[Z] ) * ( c->v[Y] - a->v[Y] ) == 0
279
 
      && ( b->v[Z] - a->v[Z] ) * ( c->v[X] - a->v[X] ) -
280
 
         ( b->v[X] - a->v[X] ) * ( c->v[Z] - a->v[Z] ) == 0
281
 
      && ( b->v[X] - a->v[X] ) * ( c->v[Y] - a->v[Y] ) -
282
 
         ( b->v[Y] - a->v[Y] ) * ( c->v[X] - a->v[X] ) == 0  ;
283
 
}
284
 
 
285
 
/*---------------------------------------------------------------------
286
 
VolumeSign returns the sign of the volume of the tetrahedron determined by f
287
 
and p.  VolumeSign is +1 iff p is on the negative side of f,
288
 
where the positive side is determined by the rh-rule.  So the volume
289
 
is positive if the ccw normal to f points outside the tetrahedron.
290
 
The final fewer-multiplications form is due to Bob Williamson.
291
 
---------------------------------------------------------------------*/
292
 
int  VolumeSign( LPFACE f, LPVERTEX p )
293
 
{
294
 
   double  vol;
295
 
   double  ax, ay, az, bx, by, bz, cx, cy, cz;
296
 
 
297
 
   ax = f->Vertex[0]->v[X] - p->v[X];
298
 
   ay = f->Vertex[0]->v[Y] - p->v[Y];
299
 
   az = f->Vertex[0]->v[Z] - p->v[Z];
300
 
   bx = f->Vertex[1]->v[X] - p->v[X];
301
 
   by = f->Vertex[1]->v[Y] - p->v[Y];
302
 
   bz = f->Vertex[1]->v[Z] - p->v[Z];
303
 
   cx = f->Vertex[2]->v[X] - p->v[X];
304
 
   cy = f->Vertex[2]->v[Y] - p->v[Y];
305
 
   cz = f->Vertex[2]->v[Z] - p->v[Z];
306
 
 
307
 
   vol =   ax * (by*cz - bz*cy)
308
 
         + ay * (bz*cx - bx*cz)
309
 
         + az * (bx*cy - by*cx);
310
 
 
311
 
 
312
 
   /* The volume should be an integer. */
313
 
   if      ( vol >  0.5 )  return  1;
314
 
   else if ( vol < -0.5 )  return -1;
315
 
   else                    return  0;
316
 
}
317
 
 
318
 
 
319
 
 
320
 
/*---------------------------------------------------------------------
321
 
CleanEdges runs through the Edge list and cleans up the structure.
322
 
If there is a NewFace then it will put that face in place of the
323
 
Visible face and NULL out NewFace. It also deletes so marked edges.
324
 
---------------------------------------------------------------------*/
325
 
static
326
 
void CleanEdges(LPHULL hull)
327
 
{
328
 
   LPEDGE  e;    /* Primary index into Edge list. */
329
 
   LPEDGE  t;    /* Temporary Edge pointer. */
330
 
 
331
 
   /* Integrate the NewFace's into the data structure. */
332
 
   /* Check every Edge. */
333
 
 
334
 
   e = hull ->edges;
335
 
   do {
336
 
        if ( e->NewFace ) {
337
 
 
338
 
            if ( e->AdjFace[0]->Visible )
339
 
                    e->AdjFace[0] = e->NewFace;
340
 
            else
341
 
                    e->AdjFace[1] = e->NewFace;
342
 
 
343
 
            e->NewFace = NULL;
344
 
        }
345
 
 
346
 
      e = e->Next;
347
 
 
348
 
   } while ( e != hull ->edges );
349
 
 
350
 
   /* Delete any edges marked for deletion. */
351
 
   while ( hull ->edges && hull ->edges->DoDelete ) {
352
 
 
353
 
      e = hull ->edges;
354
 
 
355
 
      XDELETE( hull ->edges, e );
356
 
   }
357
 
 
358
 
   e = hull ->edges->Next;
359
 
 
360
 
   do {
361
 
        if ( e->DoDelete ) {
362
 
 
363
 
            t = e;
364
 
            e = e->Next;
365
 
            XDELETE( hull ->edges, t );
366
 
        }
367
 
        else e = e->Next;
368
 
 
369
 
   } while ( e != hull ->edges );
370
 
}
371
 
 
372
 
/*---------------------------------------------------------------------
373
 
CleanFaces runs through the face list and deletes any face marked Visible.
374
 
---------------------------------------------------------------------*/
375
 
static
376
 
void CleanFaces(LPHULL hull)
377
 
{
378
 
   LPFACE  f;    /* Primary pointer into face list. */
379
 
   LPFACE  t;    /* Temporary pointer, for deleting. */
380
 
 
381
 
 
382
 
   while ( hull ->faces && hull ->faces->Visible ) {
383
 
 
384
 
        f = hull ->faces;
385
 
        XDELETE( hull ->faces, f );
386
 
   }
387
 
 
388
 
   f = hull ->faces->Next;
389
 
 
390
 
   do {
391
 
      if ( f->Visible ) {
392
 
 
393
 
            t = f;
394
 
            f = f->Next;
395
 
            XDELETE( hull ->faces, t );
396
 
      }
397
 
      else f = f->Next;
398
 
 
399
 
   } while ( f != hull ->faces );
400
 
}
401
 
 
402
 
 
403
 
 
404
 
/*---------------------------------------------------------------------
405
 
CleanVertices runs through the Vertex list and deletes the
406
 
vertices that are marked as processed but are not incident to any
407
 
undeleted edges.
408
 
---------------------------------------------------------------------*/
409
 
static
410
 
void CleanVertices(LPHULL hull)
411
 
{
412
 
   LPEDGE    e;
413
 
   LPVERTEX  v, t;
414
 
 
415
 
   /* Mark all vertices incident to some undeleted Edge as on the hull. */
416
 
 
417
 
   e = hull ->edges;
418
 
   do {
419
 
        e->EndPts[0]->onhull = e->EndPts[1]->onhull = ONHULL;
420
 
        e = e->Next;
421
 
 
422
 
   } while (e != hull ->edges);
423
 
 
424
 
 
425
 
   /* Delete all vertices that have been processed but
426
 
      are not on the hull. */
427
 
 
428
 
   while ( hull ->vertices && hull->vertices->mark && !hull ->vertices->onhull ) {
429
 
 
430
 
        v = hull ->vertices;
431
 
        XDELETE(hull ->vertices, v );
432
 
   }
433
 
 
434
 
 
435
 
   v = hull ->vertices->Next;
436
 
   do {
437
 
        if (v->mark && !v->onhull ) {
438
 
                t = v;
439
 
                v = v->Next;
440
 
                XDELETE(hull ->vertices, t )
441
 
        }
442
 
        else
443
 
                v = v->Next;
444
 
 
445
 
   } while ( v != hull ->vertices );
446
 
 
447
 
 
448
 
   /* Reset flags. */
449
 
 
450
 
   v = hull ->vertices;
451
 
   do {
452
 
            v->duplicate = NULL;
453
 
            v->onhull = !ONHULL;
454
 
            v = v->Next;
455
 
 
456
 
   } while (v != hull->vertices );
457
 
 
458
 
}
459
 
 
460
 
 
461
 
 
462
 
 
463
 
/*---------------------------------------------------------------------
464
 
MakeCcw puts the vertices in the face structure in counterclock wise
465
 
order.  We want to store the vertices in the same
466
 
order as in the Visible face.  The third Vertex is always p.
467
 
 
468
 
Although no specific ordering of the edges of a face are used
469
 
by the code, the following condition is maintained for each face f:
470
 
one of the two endpoints of f->Edge[i] matches f->Vertex[i].
471
 
But note that this does not imply that f->Edge[i] is between
472
 
f->Vertex[i] and f->Vertex[(i+1)%3].  (Thanks to Bob Williamson.)
473
 
---------------------------------------------------------------------*/
474
 
 
475
 
static
476
 
void MakeCcw(LPFACE f, LPEDGE e, LPVERTEX p)
477
 
{
478
 
   LPFACE  fv;   /* The Visible face adjacent to e */
479
 
   int    i;    /* Index of e->endpoint[0] in fv. */
480
 
   LPEDGE  s;    /* Temporary, for swapping */
481
 
 
482
 
   if  (e->AdjFace[0]->Visible)
483
 
 
484
 
        fv = e->AdjFace[0];
485
 
   else
486
 
        fv = e->AdjFace[1];
487
 
 
488
 
   /* Set Vertex[0] & [1] of f to have the same orientation
489
 
      as do the corresponding vertices of fv. */
490
 
 
491
 
   for ( i=0; fv->Vertex[i] != e->EndPts[0]; ++i )
492
 
      ;
493
 
 
494
 
   /* Orient f the same as fv. */
495
 
 
496
 
   if ( fv->Vertex[ (i+1) % 3 ] != e->EndPts[1] ) {
497
 
 
498
 
        f->Vertex[0] = e->EndPts[1];
499
 
        f->Vertex[1] = e->EndPts[0];
500
 
   }
501
 
   else {
502
 
        f->Vertex[0] = e->EndPts[0];
503
 
        f->Vertex[1] = e->EndPts[1];
504
 
        SWAP( s, f->Edge[1], f->Edge[2] );
505
 
   }
506
 
 
507
 
   /* This swap is tricky. e is Edge[0]. Edge[1] is based on endpt[0],
508
 
      Edge[2] on endpt[1].  So if e is oriented "forwards," we
509
 
      need to move Edge[1] to follow [0], because it precedes. */
510
 
 
511
 
   f->Vertex[2] = p;
512
 
}
513
 
 
514
 
/*---------------------------------------------------------------------
515
 
MakeConeFace makes a new face and two new edges between the
516
 
Edge and the point that are passed to it. It returns a pointer to
517
 
the new face.
518
 
---------------------------------------------------------------------*/
519
 
 
520
 
static
521
 
LPFACE MakeConeFace(LPHULL hull, LPEDGE e, LPVERTEX p)
522
 
{
523
 
   LPEDGE  new_edge[2];
524
 
   LPFACE  new_face;
525
 
   int       i, j;
526
 
 
527
 
   /* Make two new edges (if don't already exist). */
528
 
 
529
 
   for ( i=0; i < 2; ++i )
530
 
      /* If the Edge exists, copy it into new_edge. */
531
 
      if ( !( new_edge[i] = e->EndPts[i]->duplicate) ) {
532
 
 
533
 
        /* Otherwise (duplicate is NULL), MakeNullEdge. */
534
 
        new_edge[i] = MakeNullEdge(hull);
535
 
        new_edge[i]->EndPts[0] = e->EndPts[i];
536
 
        new_edge[i]->EndPts[1] = p;
537
 
        e->EndPts[i]->duplicate = new_edge[i];
538
 
      }
539
 
 
540
 
    /* Make the new face. */
541
 
    new_face = MakeNullFace(hull);
542
 
    new_face->Edge[0] = e;
543
 
    new_face->Edge[1] = new_edge[0];
544
 
    new_face->Edge[2] = new_edge[1];
545
 
    MakeCcw( new_face, e, p );
546
 
 
547
 
   /* Set the adjacent face pointers. */
548
 
    for ( i=0; i < 2; ++i )
549
 
        for ( j=0; j < 2; ++j )
550
 
            /* Only one NULL link should be set to new_face. */
551
 
            if ( !new_edge[i]->AdjFace[j] ) {
552
 
                    new_edge[i]->AdjFace[j] = new_face;
553
 
                    break;
554
 
            }
555
 
 
556
 
   return new_face;
557
 
}
558
 
 
559
 
 
560
 
/*---------------------------------------------------------------------
561
 
AddOne is passed a Vertex.  It first determines all faces Visible from
562
 
that point.  If none are Visible then the point is marked as not
563
 
onhull.  Next is a loop over edges.  If both faces adjacent to an Edge
564
 
are Visible, then the Edge is marked for deletion.  If just one of the
565
 
adjacent faces is Visible then a new face is constructed.
566
 
---------------------------------------------------------------------*/
567
 
static
568
 
BOOL AddOne(LPHULL hull, LPVERTEX p)
569
 
{
570
 
   LPFACE  f;
571
 
   LPEDGE  e, temp;
572
 
   int       vol;
573
 
   BOOL      vis = false;
574
 
 
575
 
 
576
 
   /* Mark faces Visible from p. */
577
 
   f = hull -> faces;
578
 
 
579
 
   do {
580
 
 
581
 
       vol = VolumeSign(f, p);
582
 
 
583
 
      if ( vol < 0 ) {
584
 
            f->Visible = VISIBLE;
585
 
            vis = true;
586
 
      }
587
 
 
588
 
      f = f->Next;
589
 
 
590
 
   } while ( f != hull ->faces );
591
 
 
592
 
   /* If no faces are Visible from p, then p is inside the hull. */
593
 
 
594
 
   if ( !vis ) {
595
 
 
596
 
            p->onhull = !ONHULL;
597
 
            return false;
598
 
   }
599
 
 
600
 
   /* Mark edges in interior of Visible region for deletion.
601
 
      Erect a NewFace based on each border Edge. */
602
 
 
603
 
   e = hull ->edges;
604
 
 
605
 
   do {
606
 
 
607
 
      temp = e->Next;
608
 
 
609
 
      if ( e->AdjFace[0]->Visible && e->AdjFace[1]->Visible )
610
 
                /* e interior: mark for deletion. */
611
 
                e->DoDelete = REMOVED;
612
 
 
613
 
      else
614
 
          if ( e->AdjFace[0]->Visible || e->AdjFace[1]->Visible )
615
 
                /* e border: make a new face. */
616
 
                e->NewFace = MakeConeFace(hull, e, p );
617
 
 
618
 
        e = temp;
619
 
 
620
 
   } while ( e != hull ->edges );
621
 
 
622
 
   return true;
623
 
}
624
 
 
625
 
 
626
 
/*---------------------------------------------------------------------
627
 
 DoubleTriangle builds the initial double triangle.  It first finds 3
628
 
 noncollinear points and makes two faces out of them, in opposite order.
629
 
 It then finds a fourth point that is not coplanar with that face.  The
630
 
 vertices are stored in the face structure in counterclockwise order so
631
 
 that the volume between the face and the point is negative. Lastly, the
632
 
 3 newfaces to the fourth point are constructed and the data structures
633
 
 are cleaned up.
634
 
---------------------------------------------------------------------*/
635
 
 
636
 
static
637
 
BOOL DoubleTriangle(LPHULL hull)
638
 
{
639
 
   LPVERTEX  v0, v1, v2, v3;
640
 
   LPFACE    f0, f1 = NULL;
641
 
   int      vol;
642
 
 
643
 
   /* Find 3 noncollinear points. */
644
 
   v0 = hull ->vertices;
645
 
   while ( Collinear( v0, v0->Next, v0->Next->Next ) )
646
 
      if ( ( v0 = v0->Next ) == hull->vertices )
647
 
                return false; /* All points are Collinear! */
648
 
 
649
 
   v1 = v0->Next;
650
 
   v2 = v1->Next;
651
 
 
652
 
   /* Mark the vertices as processed. */
653
 
   v0->mark = PROCESSED;
654
 
   v1->mark = PROCESSED;
655
 
   v2->mark = PROCESSED;
656
 
 
657
 
   /* Create the two "twin" faces. */
658
 
   f0 = MakeFace(hull, v0, v1, v2, f1 );
659
 
   f1 = MakeFace(hull, v2, v1, v0, f0 );
660
 
 
661
 
   /* Link adjacent face fields. */
662
 
   f0->Edge[0]->AdjFace[1] = f1;
663
 
   f0->Edge[1]->AdjFace[1] = f1;
664
 
   f0->Edge[2]->AdjFace[1] = f1;
665
 
   f1->Edge[0]->AdjFace[1] = f0;
666
 
   f1->Edge[1]->AdjFace[1] = f0;
667
 
   f1->Edge[2]->AdjFace[1] = f0;
668
 
 
669
 
   /* Find a fourth, noncoplanar point to form tetrahedron. */
670
 
   v3 = v2->Next;
671
 
   vol = VolumeSign( f0, v3 );
672
 
 
673
 
   while ( !vol )   {
674
 
 
675
 
      if ( ( v3 = v3->Next ) == v0 )
676
 
            return false; /* All points are coplanar! */
677
 
 
678
 
      vol = VolumeSign( f0, v3 );
679
 
   }
680
 
 
681
 
   /* Insure that v3 will be the first added. */
682
 
   hull ->vertices = v3;
683
 
   return true;
684
 
}
685
 
 
686
 
 
687
 
 
688
 
/*---------------------------------------------------------------------
689
 
ConstructHull adds the vertices to the hull one at a time.  The hull
690
 
vertices are those in the list marked as onhull.
691
 
---------------------------------------------------------------------*/
692
 
static
693
 
void ConstructHull(LPHULL hull)
694
 
{
695
 
 LPVERTEX  v, vnext;
696
 
 BOOL     changed;    /* T if addition changes hull; not used. */
697
 
 
698
 
 v = hull->vertices;
699
 
 
700
 
 do {
701
 
        vnext = v->Next;
702
 
 
703
 
        changed = false;
704
 
 
705
 
        if (!v->mark ) {
706
 
 
707
 
                v->mark = PROCESSED;
708
 
                changed = AddOne(hull,  v );
709
 
 
710
 
                CleanEdges(hull);
711
 
                CleanFaces(hull);
712
 
                CleanVertices(hull);
713
 
    }
714
 
 
715
 
    v = vnext;
716
 
 
717
 
 } while (v != hull->vertices );
718
 
 
719
 
}
720
 
 
721
 
 
722
 
 
723
 
/*-------------------------------------------------------------------*/
724
 
 
725
 
 
726
 
static
727
 
void AddVec( VEC3I q, VEC3I ray )
728
 
{
729
 
  int i;
730
 
 
731
 
  for( i = 0; i < DIM; i++ )
732
 
    ray[i] = q[i] + ray[i];
733
 
}
734
 
 
735
 
/*---------------------------------------------------------------------
736
 
a - b ==> c.
737
 
---------------------------------------------------------------------*/
738
 
static
739
 
void  SubVec( VEC3I a, VEC3I b, VEC3I c )
740
 
{
741
 
   int i;
742
 
 
743
 
   for( i = 0; i < DIM; i++ )
744
 
      c[i] = a[i] - b[i];
745
 
}
746
 
 
747
 
 
748
 
/*---------------------------------------------------------------------
749
 
Returns the dot product of the two input vectors.
750
 
---------------------------------------------------------------------*/
751
 
static
752
 
double    Dot( VEC3I a, LPVEC3 b )
753
 
{
754
 
    int i;
755
 
    double sum = 0.0;
756
 
 
757
 
    for( i = 0; i < DIM; i++ )
758
 
       sum += a[i] * b->n[i];
759
 
 
760
 
    return  sum;
761
 
}
762
 
 
763
 
/*---------------------------------------------------------------------
764
 
Compute the cross product of (b-a)x(c-a) and place into N.
765
 
---------------------------------------------------------------------*/
766
 
static
767
 
void    NormalVec( VEC3I a, VEC3I b, VEC3I c, LPVEC3 N )
768
 
{
769
 
    N->n[X] = ( c[Z] - a[Z] ) * ( b[Y] - a[Y] ) -
770
 
           ( b[Z] - a[Z] ) * ( c[Y] - a[Y] );
771
 
    N->n[Y] = ( b[Z] - a[Z] ) * ( c[X] - a[X] ) -
772
 
           ( b[X] - a[X] ) * ( c[Z] - a[Z] );
773
 
    N->n[Z] = ( b[X] - a[X] ) * ( c[Y] - a[Y] ) -
774
 
           ( b[Y] - a[Y] ) * ( c[X] - a[X] );
775
 
}
776
 
 
777
 
 
778
 
 
779
 
 
780
 
static
781
 
int InBox( VEC3I q, VEC3I bmin, VEC3I bmax )
782
 
{
783
 
 
784
 
  if( ( bmin[X] <= q[X] ) && ( q[X] <= bmax[X] ) &&
785
 
      ( bmin[Y] <= q[Y] ) && ( q[Y] <= bmax[Y] ) &&
786
 
      ( bmin[Z] <= q[Z] ) && ( q[Z] <= bmax[Z] ) )
787
 
    return true;
788
 
 
789
 
  return false;
790
 
}
791
 
 
792
 
 
793
 
 
794
 
/*
795
 
  This function returns a char:
796
 
    '0': the segment [ab] does not intersect (completely misses) the
797
 
         bounding box surrounding the n-th triangle T.  It lies
798
 
         strictly to one side of one of the six supporting planes.
799
 
    '?': status unknown: the segment may or may not intersect T.
800
 
*/
801
 
 
802
 
static
803
 
char BoxTest(LPHULL hull, int n, VEC3I a, VEC3I b)
804
 
{
805
 
   int i; /* Coordinate index */
806
 
   int w;
807
 
 
808
 
   for ( i=0; i < DIM; i++ ) {
809
 
 
810
 
       w = hull ->Box[ n ][0][i]; /* min: lower left */
811
 
 
812
 
       if ( (a[i] < w) && (b[i] < w) ) return '0';
813
 
 
814
 
       w = hull ->Box[ n ][1][i]; /* max: upper right */
815
 
 
816
 
       if ( (a[i] > w) && (b[i] > w) ) return '0';
817
 
   }
818
 
 
819
 
   return '?';
820
 
}
821
 
 
822
 
 
823
 
 
824
 
/* Return a random ray endpoint */
825
 
 
826
 
static
827
 
void RandomRay( VEC3I ray, int radius )
828
 
{
829
 
  double x, y, z, w, t;
830
 
 
831
 
  /* Generate a random point on a sphere of radius 1. */
832
 
  /* the sphere is sliced at z, and a random point at angle t
833
 
     generated on the circle of intersection. */
834
 
 
835
 
  z = 2.0 * (double) rand() / RAND_MAX - 1.0;
836
 
  t = 2.0 * M_PI * (double) rand() / RAND_MAX;
837
 
  w = sqrt( 1 - z*z );
838
 
  x = w * cos( t );
839
 
  y = w * sin( t );
840
 
 
841
 
  ray[X] = (int) ( radius * x );
842
 
  ray[Y] = (int) ( radius * y );
843
 
  ray[Z] = (int) ( radius * z );
844
 
}
845
 
 
846
 
 
847
 
 
848
 
static
849
 
int ComputeBox(LPHULL hull, int F, VEC3I bmin, VEC3I bmax )
850
 
{
851
 
  int i, j;
852
 
  double radius;
853
 
 
854
 
  for( i = 0; i < F; i++ )
855
 
    for( j = 0; j < DIM; j++ ) {
856
 
 
857
 
      if( hull ->Vertices[i][j] < bmin[j] )
858
 
            bmin[j] = hull ->Vertices[i][j];
859
 
 
860
 
      if( hull ->Vertices[i][j] > bmax[j] )
861
 
                bmax[j] = hull ->Vertices[i][j];
862
 
    }
863
 
 
864
 
  radius = sqrt( pow( (double)(bmax[X] - bmin[X]), 2.0 ) +
865
 
                 pow( (double)(bmax[Y] - bmin[Y]), 2.0 ) +
866
 
                 pow( (double)(bmax[Z] - bmin[Z]), 2.0 ) );
867
 
 
868
 
  return (int)( radius +1 ) + 1;
869
 
}
870
 
 
871
 
 
872
 
/*---------------------------------------------------------------------
873
 
Computes N & D and returns index m of largest component.
874
 
---------------------------------------------------------------------*/
875
 
static
876
 
int    PlaneCoeff(LPHULL hull,  VEC3I T, LPVEC3 N, double *D )
877
 
{
878
 
    int i;
879
 
    double t;              /* Temp storage */
880
 
    double biggest = 0.0;  /* Largest component of normal vector. */
881
 
    int m = 0;             /* Index of largest component. */
882
 
 
883
 
    NormalVec(hull ->Vertices[T[0]], hull ->Vertices[T[1]], hull ->Vertices[T[2]], N );
884
 
    *D = Dot( hull ->Vertices[T[0]], N );
885
 
 
886
 
    /* Find the largest component of N. */
887
 
    for ( i = 0; i < DIM; i++ ) {
888
 
      t = fabs( N->n[i] );
889
 
      if ( t > biggest ) {
890
 
        biggest = t;
891
 
        m = i;
892
 
      }
893
 
    }
894
 
    return m;
895
 
}
896
 
 
897
 
/*---------------------------------------------------------------------
898
 
    'p': The segment lies wholly within the plane.
899
 
    'q': The q endpoint is on the plane (but not 'p').
900
 
    'r': The r endpoint is on the plane (but not 'p').
901
 
    '0': The segment lies strictly to one side or the other of the plane.
902
 
    '1': The segement intersects the plane, and 'p' does not hold.
903
 
---------------------------------------------------------------------*/
904
 
static
905
 
char SegPlaneInt(LPHULL hull, VEC3I T, VEC3I q, VEC3I r, LPVEC3 p, int *m)
906
 
{
907
 
    VEC3 N; double D;
908
 
    VEC3I rq;
909
 
    double num, denom, t;
910
 
    int i;
911
 
 
912
 
    *m = PlaneCoeff(hull, T, &N, &D );
913
 
    num = D - Dot( q, &N );
914
 
    SubVec( r, q, rq );
915
 
    denom = Dot( rq, &N );
916
 
 
917
 
    if ( denom == 0.0 ) {  /* Segment is parallel to plane. */
918
 
       if ( num == 0.0 )   /* q is on plane. */
919
 
           return 'p';
920
 
       else
921
 
           return '0';
922
 
    }
923
 
    else
924
 
       t = num / denom;
925
 
 
926
 
    for( i = 0; i < DIM; i++ )
927
 
       p->n[i] = q[i] + t * ( r[i] - q[i] );
928
 
 
929
 
    if ( (0.0 < t) && (t < 1.0) )
930
 
         return '1';
931
 
    else if ( num == 0.0 )   /* t == 0 */
932
 
         return 'q';
933
 
    else if ( num == denom ) /* t == 1 */
934
 
         return 'r';
935
 
    else return '0';
936
 
}
937
 
 
938
 
 
939
 
 
940
 
static
941
 
int  AreaSign( VEC3I a, VEC3I b, VEC3I c )
942
 
{
943
 
    double area2;
944
 
 
945
 
    area2 = ( b[0] - a[0] ) * (double)( c[1] - a[1] ) -
946
 
            ( c[0] - a[0] ) * (double)( b[1] - a[1] );
947
 
 
948
 
    /* The area should be an integer. */
949
 
    if      ( area2 >  0.5 ) return  1;
950
 
    else if ( area2 < -0.5 ) return -1;
951
 
    else                     return  0;
952
 
}
953
 
 
954
 
 
955
 
static
956
 
char     InTri2D( VEC3I Tp[3], VEC3I pp )
957
 
{
958
 
   int area0, area1, area2;
959
 
 
960
 
   /* compute three AreaSign() values for pp w.r.t. each Edge of the face in 2D */
961
 
   area0 = AreaSign( pp, Tp[0], Tp[1] );
962
 
   area1 = AreaSign( pp, Tp[1], Tp[2] );
963
 
   area2 = AreaSign( pp, Tp[2], Tp[0] );
964
 
 
965
 
   if ( (( area0 == 0 ) && ( area1 > 0 ) && ( area2 > 0 )) ||
966
 
        (( area1 == 0 ) && ( area0 > 0 ) && ( area2 > 0 )) ||
967
 
        (( area2 == 0 ) && ( area0 > 0 ) && ( area1 > 0 )) )
968
 
     return 'E';
969
 
 
970
 
   if ( (( area0 == 0 ) && ( area1 < 0 ) && ( area2 < 0 )) ||
971
 
        (( area1 == 0 ) && ( area0 < 0 ) && ( area2 < 0 )) ||
972
 
        (( area2 == 0 ) && ( area0 < 0 ) && ( area1 < 0 )))
973
 
     return 'E';
974
 
 
975
 
   if ( (( area0 >  0 ) && ( area1 > 0 ) && ( area2 > 0 )) ||
976
 
        (( area0 <  0 ) && ( area1 < 0 ) && ( area2 < 0 )))
977
 
     return 'F';
978
 
 
979
 
   if ( ( area0 == 0 ) && ( area1 == 0 ) && ( area2 == 0 ) )
980
 
     return '?'; /* Error in InTriD */
981
 
 
982
 
   if ( (( area0 == 0 ) && ( area1 == 0 )) ||
983
 
        (( area0 == 0 ) && ( area2 == 0 )) ||
984
 
        (( area1 == 0 ) && ( area2 == 0 )) )
985
 
     return 'V';
986
 
 
987
 
   else
988
 
     return '0';
989
 
}
990
 
 
991
 
/* Assumption: p lies in the plane containing T.
992
 
    Returns a char:
993
 
     'V': the query point p coincides with a Vertex of triangle T.
994
 
     'E': the query point p is in the relative interior of an Edge of triangle T.
995
 
     'F': the query point p is in the relative interior of a Face of triangle T.
996
 
     '0': the query point p does not intersect (misses) triangle T.
997
 
*/
998
 
 
999
 
static
1000
 
char     InTri3D(LPHULL hull,  VEC3I T, int m, VEC3I p )
1001
 
{
1002
 
   int i;           /* Index for X,Y,Z           */
1003
 
   int j;           /* Index for X,Y             */
1004
 
   int k;           /* Index for triangle Vertex */
1005
 
   VEC3I pp;      /* projected p */
1006
 
   VEC3I Tp[3];   /* projected T: three new vertices */
1007
 
 
1008
 
   /* Project out coordinate m in both p and the triangular face */
1009
 
   j = 0;
1010
 
   for ( i = 0; i < DIM; i++ ) {
1011
 
     if ( i != m ) {    /* skip largest coordinate */
1012
 
       pp[j] = p[i];
1013
 
       for ( k = 0; k < 3; k++ )
1014
 
            Tp[k][j] = hull->Vertices[T[k]][i];
1015
 
       j++;
1016
 
     }
1017
 
   }
1018
 
   return( InTri2D( Tp, pp ) );
1019
 
}
1020
 
 
1021
 
 
1022
 
 
1023
 
static
1024
 
int     VolumeSign2( VEC3I a, VEC3I b, VEC3I c, VEC3I d )
1025
 
{
1026
 
   double vol;
1027
 
   double ax, ay, az, bx, by, bz, cx, cy, cz, dx, dy, dz;
1028
 
   double bxdx, bydy, bzdz, cxdx, cydy, czdz;
1029
 
 
1030
 
   ax = a[X];
1031
 
   ay = a[Y];
1032
 
   az = a[Z];
1033
 
   bx = b[X];
1034
 
   by = b[Y];
1035
 
   bz = b[Z];
1036
 
   cx = c[X];
1037
 
   cy = c[Y];
1038
 
   cz = c[Z];
1039
 
   dx = d[X];
1040
 
   dy = d[Y];
1041
 
   dz = d[Z];
1042
 
 
1043
 
   bxdx=bx-dx;
1044
 
   bydy=by-dy;
1045
 
   bzdz=bz-dz;
1046
 
   cxdx=cx-dx;
1047
 
   cydy=cy-dy;
1048
 
   czdz=cz-dz;
1049
 
   vol =   (az-dz) * (bxdx*cydy - bydy*cxdx)
1050
 
         + (ay-dy) * (bzdz*cxdx - bxdx*czdz)
1051
 
         + (ax-dx) * (bydy*czdz - bzdz*cydy);
1052
 
 
1053
 
 
1054
 
   /* The volume should be an integer. */
1055
 
   if      ( vol > 0.5 )   return  1;
1056
 
   else if ( vol < -0.5 )  return -1;
1057
 
   else                    return  0;
1058
 
}
1059
 
 
1060
 
 
1061
 
 
1062
 
 
1063
 
/*---------------------------------------------------------------------
1064
 
The signed volumes of three tetrahedra are computed, determined
1065
 
by the segment qr, and each Edge of the triangle.
1066
 
Returns a char:
1067
 
   'v': the open segment includes a Vertex of T.
1068
 
   'e': the open segment includes a point in the relative interior of an Edge
1069
 
   of T.
1070
 
   'f': the open segment includes a point in the relative interior of a face
1071
 
   of T.
1072
 
   '0': the open segment does not intersect triangle T.
1073
 
---------------------------------------------------------------------*/
1074
 
 
1075
 
static
1076
 
char SegTriCross(LPHULL hull,  VEC3I T, VEC3I q, VEC3I r )
1077
 
{
1078
 
   int vol0, vol1, vol2;
1079
 
 
1080
 
   vol0 = VolumeSign2( q, hull->Vertices[ T[0] ], hull->Vertices[ T[1] ], r );
1081
 
   vol1 = VolumeSign2( q, hull->Vertices[ T[1] ], hull->Vertices[ T[2] ], r );
1082
 
   vol2 = VolumeSign2( q, hull->Vertices[ T[2] ], hull->Vertices[ T[0] ], r );
1083
 
 
1084
 
 
1085
 
   /* Same sign: segment intersects interior of triangle. */
1086
 
   if ( ( ( vol0 > 0 ) && ( vol1 > 0 ) && ( vol2 > 0 ) ) ||
1087
 
        ( ( vol0 < 0 ) && ( vol1 < 0 ) && ( vol2 < 0 ) ) )
1088
 
      return 'f';
1089
 
 
1090
 
   /* Opposite sign: no intersection between segment and triangle */
1091
 
   if ( ( ( vol0 > 0 ) || ( vol1 > 0 ) || ( vol2 > 0 ) ) &&
1092
 
        ( ( vol0 < 0 ) || ( vol1 < 0 ) || ( vol2 < 0 ) ) )
1093
 
      return '0';
1094
 
 
1095
 
   else if ( ( vol0 == 0 ) && ( vol1 == 0 ) && ( vol2 == 0 ) )
1096
 
     return '?'; /* Error 1 in SegTriCross */
1097
 
 
1098
 
   /* Two zeros: segment intersects Vertex. */
1099
 
   else if ( ( ( vol0 == 0 ) && ( vol1 == 0 ) ) ||
1100
 
             ( ( vol0 == 0 ) && ( vol2 == 0 ) ) ||
1101
 
             ( ( vol1 == 0 ) && ( vol2 == 0 ) ) )
1102
 
      return 'v';
1103
 
 
1104
 
   /* One zero: segment intersects Edge. */
1105
 
   else if ( ( vol0 == 0 ) || ( vol1 == 0 ) || ( vol2 == 0 ) )
1106
 
      return 'e';
1107
 
 
1108
 
   else
1109
 
     return '?'; /* Error 2 in SegTriCross */
1110
 
}
1111
 
 
1112
 
 
1113
 
 
1114
 
static
1115
 
char    SegTriInt(LPHULL hull, VEC3I T, VEC3I q, VEC3I r, LPVEC3 p )
1116
 
{
1117
 
    int code;
1118
 
    int m = -1;
1119
 
 
1120
 
    code = SegPlaneInt(hull, T, q, r, p, &m );
1121
 
 
1122
 
    if ( code == '0')        return '0';
1123
 
    else if ( code == 'q')   return InTri3D(hull, T, m, q );
1124
 
    else if ( code == 'r')   return InTri3D(hull,  T, m, r );
1125
 
    else if ( code == 'p')   return 'p';
1126
 
    else if ( code == '1' )  return SegTriCross(hull, T, q, r );
1127
 
    else
1128
 
       return code;        /* Error */
1129
 
}
1130
 
 
1131
 
 
1132
 
 
1133
 
 
1134
 
/*
1135
 
  This function returns a char:
1136
 
    'i': the query point a is strictly interior to polyhedron P.
1137
 
    'o': the query point a is strictly exterior to( or outside of) polyhedron P.
1138
 
*/
1139
 
char InPolyhedron(LPHULL hull, VEC3I q)
1140
 
{
1141
 
   int F = hull->nfaces;
1142
 
   VEC3I Ray;  /* Ray endpoint. */
1143
 
   VEC3 p;  /* Intersection point; not used. */
1144
 
   int f, k = 0, crossings = 0;
1145
 
   char code = '?';
1146
 
 
1147
 
 
1148
 
   /* If query point is outside bounding box, finished. */
1149
 
   if ( !InBox( q, hull->bmin, hull->bmax ) )
1150
 
      return 'o';
1151
 
 
1152
 
   LOOP:
1153
 
   while( k++ < F ) {
1154
 
 
1155
 
      crossings = 0;
1156
 
 
1157
 
      RandomRay(Ray, hull->radius );
1158
 
 
1159
 
      AddVec( q, Ray );
1160
 
 
1161
 
 
1162
 
      for ( f = 0; f < F; f++ ) {  /* Begin check each face */
1163
 
 
1164
 
         if ( BoxTest(hull,  f, q, Ray ) == '0' ) {
1165
 
              code = '0';
1166
 
 
1167
 
         }
1168
 
         else code = SegTriInt(hull, hull->Faces[f], q, Ray, &p );
1169
 
 
1170
 
 
1171
 
         /* If ray is degenerate, then goto outer while to generate another. */
1172
 
         if ( code == 'p' || code == 'v' || code == 'e' ) {
1173
 
 
1174
 
            goto LOOP;
1175
 
         }
1176
 
 
1177
 
         /* If ray hits face at interior point, increment crossings. */
1178
 
         else if ( code == 'f' ) {
1179
 
            crossings++;
1180
 
 
1181
 
         }
1182
 
 
1183
 
         /* If query endpoint q sits on a V/E/F, return inside. */
1184
 
         else if ( code == 'V' || code == 'E' || code == 'F' )
1185
 
            return code; /* 'i'; MM2 */
1186
 
 
1187
 
         /* If ray misses triangle, do nothing. */
1188
 
         else if ( code == '0' )
1189
 
            ;
1190
 
 
1191
 
         else
1192
 
            return '?'; /* Error */
1193
 
 
1194
 
      }
1195
 
      /* No degeneracies encountered: ray is generic, so finished. */
1196
 
      break;
1197
 
 
1198
 
   } /* End while loop */
1199
 
 
1200
 
 
1201
 
   /* q strictly interior to polyhedron iff an odd number of crossings. */
1202
 
   if( ( crossings % 2 ) == 1 )
1203
 
      return   'i';
1204
 
 
1205
 
   else return 'o';
1206
 
}
1207
 
 
1208
 
 
1209
 
/*/ ---------------------------------------------------------------------------------- */
1210
 
 
1211
 
 
1212
 
 
1213
 
 
1214
 
static
1215
 
void StoreResults(LPHULL hull)
1216
 
{
1217
 
 
1218
 
   int       i, w;
1219
 
   LPVERTEX  v;
1220
 
   LPFACE    f;
1221
 
   int       V = 0, F = 0;
1222
 
   int       j, k;
1223
 
 
1224
 
   /* Vertices */
1225
 
 
1226
 
   v = hull ->vertices;
1227
 
   V = 0;
1228
 
   do {
1229
 
 
1230
 
      v -> vnum = V;
1231
 
      hull ->Vertices[V][X] = v -> v[X];
1232
 
      hull ->Vertices[V][Y] = v -> v[Y];
1233
 
      hull ->Vertices[V][Z] = v -> v[Z];
1234
 
 
1235
 
      v = v->Next;
1236
 
      V++;
1237
 
 
1238
 
   } while ( v != hull ->vertices );
1239
 
 
1240
 
   hull ->nvertex = V;
1241
 
 
1242
 
   /* Faces */
1243
 
   f = hull ->faces;
1244
 
   F = 0;
1245
 
   do {
1246
 
 
1247
 
      hull ->Faces[F][0] = f->Vertex[0]->vnum;
1248
 
      hull ->Faces[F][1] = f->Vertex[1]->vnum;
1249
 
      hull ->Faces[F][2] = f->Vertex[2]->vnum;
1250
 
 
1251
 
      for ( j=0; j < 3; j++ ) {
1252
 
 
1253
 
       hull ->Box[F][0][j] = hull ->Vertices[ hull ->Faces[F][0] ][j];
1254
 
       hull ->Box[F][1][j] = hull ->Vertices[ hull ->Faces[F][0] ][j];
1255
 
      }
1256
 
 
1257
 
      /* Check k=1,2 vertices of face. */
1258
 
      for ( k=1; k < 3; k++ )
1259
 
        for ( j=0; j < 3; j++ ) {
1260
 
 
1261
 
            w = hull ->Vertices[ hull ->Faces[F][k] ][j];
1262
 
            if ( w < hull ->Box[F][0][j] ) hull ->Box[F][0][j] = w;
1263
 
            if ( w > hull ->Box[F][1][j] ) hull ->Box[F][1][j] = w;
1264
 
        }
1265
 
 
1266
 
 
1267
 
      f = f->Next; F++;
1268
 
 
1269
 
   } while ( f != hull ->faces );
1270
 
 
1271
 
 
1272
 
  hull ->nfaces = F;
1273
 
 
1274
 
 
1275
 
    /* Initialize the bounding box */
1276
 
  for ( i = 0; i < DIM; i++ )
1277
 
    hull ->bmin[i] = hull ->bmax[i] = hull ->Vertices[0][i];
1278
 
 
1279
 
  hull ->radius = ComputeBox(hull, V, hull ->bmin, hull ->bmax );
1280
 
 
1281
 
 
1282
 
}
1283
 
 
1284
 
 
1285
 
LCMSHANDLE cmsxHullInit(void)
1286
 
{
1287
 
    LPHULL hull = (LPHULL) malloc(sizeof(HULL));
1288
 
 
1289
 
    ZeroMemory(hull, sizeof(HULL));
1290
 
 
1291
 
    hull->vnumCounter  = 0;
1292
 
    hull->vertices     = NULL;
1293
 
    hull->edges        = NULL;
1294
 
    hull->faces        = NULL;
1295
 
    hull->nfaces       = 0;
1296
 
    hull->nvertex      = 0;
1297
 
 
1298
 
    return (LCMSHANDLE) (LPSTR) hull;
1299
 
}
1300
 
 
1301
 
 
1302
 
void cmsxHullDone(LCMSHANDLE hHull)
1303
 
{
1304
 
    LPHULL hull = (LPHULL) (LPSTR) hHull;
1305
 
 
1306
 
    if (hull)
1307
 
        free((LPVOID) hull);
1308
 
}
1309
 
 
1310
 
 
1311
 
BOOL cmsxHullAddPoint(LCMSHANDLE hHull, int x, int y, int z)
1312
 
{
1313
 
     LPVERTEX  v;
1314
 
     LPHULL hull = (LPHULL) (LPSTR) hHull;
1315
 
 
1316
 
 
1317
 
      v = MakeNullVertex(hull);
1318
 
      v->v[X] = x;
1319
 
      v->v[Y] = y;
1320
 
      v->v[Z] = z;
1321
 
      v->vnum = hull->vnumCounter++;
1322
 
 
1323
 
      return true;
1324
 
}
1325
 
 
1326
 
BOOL cmsxHullComputeHull(LCMSHANDLE hHull)
1327
 
{
1328
 
 
1329
 
  LPHULL hull = (LPHULL) (LPSTR) hHull;
1330
 
 
1331
 
  if (!DoubleTriangle(hull)) return false;
1332
 
 
1333
 
  ConstructHull(hull);
1334
 
  StoreResults(hull);
1335
 
 
1336
 
  return true;
1337
 
}
1338
 
 
1339
 
 
1340
 
char cmsxHullCheckpoint(LCMSHANDLE hHull, int x, int y, int z)
1341
 
{
1342
 
     VEC3I q;
1343
 
     LPHULL hull = (LPHULL) (LPSTR) hHull;
1344
 
 
1345
 
     q[X] = x; q[Y] = y; q[Z] = z;
1346
 
 
1347
 
     return InPolyhedron(hull, q ) ;
1348
 
}
1349
 
 
1350
 
 
1351
 
BOOL cmsxHullDumpVRML(LCMSHANDLE hHull, const char* fname)
1352
 
{
1353
 
    FILE*           fp;
1354
 
    int             i;
1355
 
    LPHULL hull = (LPHULL) (LPSTR) hHull;
1356
 
 
1357
 
    fp = fopen (fname, "wt");
1358
 
    if (fp == NULL)
1359
 
        return false;
1360
 
 
1361
 
    fprintf (fp, "#VRML V2.0 utf8\n");
1362
 
 
1363
 
    /* set the viewing orientation and distance */
1364
 
    fprintf (fp, "DEF CamTest Group {\n");
1365
 
    fprintf (fp, "\tchildren [\n");
1366
 
    fprintf (fp, "\t\tDEF Cameras Group {\n");
1367
 
    fprintf (fp, "\t\t\tchildren [\n");
1368
 
    fprintf (fp, "\t\t\t\tDEF DefaultView Viewpoint {\n");
1369
 
    fprintf (fp, "\t\t\t\t\tposition 0 0 340\n");
1370
 
    fprintf (fp, "\t\t\t\t\torientation 0 0 1 0\n");
1371
 
    fprintf (fp, "\t\t\t\t\tdescription \"default view\"\n");
1372
 
    fprintf (fp, "\t\t\t\t}\n");
1373
 
    fprintf (fp, "\t\t\t]\n");
1374
 
    fprintf (fp, "\t\t},\n");
1375
 
    fprintf (fp, "\t]\n");
1376
 
    fprintf (fp, "}\n");
1377
 
 
1378
 
    /* Output the background stuff */
1379
 
    fprintf (fp, "Background {\n");
1380
 
    fprintf (fp, "\tskyColor [\n");
1381
 
    fprintf (fp, "\t\t.5 .5 .5\n");
1382
 
    fprintf (fp, "\t]\n");
1383
 
    fprintf (fp, "}\n");
1384
 
 
1385
 
    /* Output the shape stuff */
1386
 
    fprintf (fp, "Transform {\n");
1387
 
    fprintf (fp, "\tscale 8 8 8\n");
1388
 
    fprintf (fp, "\tchildren [\n");
1389
 
 
1390
 
    /* Draw the axes as a shape: */
1391
 
    fprintf (fp, "\t\tShape {\n");
1392
 
    fprintf (fp, "\t\t\tappearance Appearance {\n");
1393
 
    fprintf (fp, "\t\t\t\tmaterial Material {\n");
1394
 
    fprintf (fp, "\t\t\t\t\tdiffuseColor 0 0.8 0\n");
1395
 
    fprintf (fp, "\t\t\t\t\temissiveColor 1.0 1.0 1.0\n");
1396
 
    fprintf (fp, "\t\t\t\t\tshininess 0.8\n");
1397
 
    fprintf (fp, "\t\t\t\t}\n");
1398
 
    fprintf (fp, "\t\t\t}\n");
1399
 
    fprintf (fp, "\t\t\tgeometry IndexedLineSet {\n");
1400
 
    fprintf (fp, "\t\t\t\tcoord Coordinate {\n");
1401
 
    fprintf (fp, "\t\t\t\t\tpoint [\n");
1402
 
    fprintf (fp, "\t\t\t\t\t0.0 0.0 0.0,\n");
1403
 
    fprintf (fp, "\t\t\t\t\t%f 0.0 0.0,\n",  255.0);
1404
 
    fprintf (fp, "\t\t\t\t\t0.0 %f 0.0,\n",  255.0);
1405
 
    fprintf (fp, "\t\t\t\t\t0.0 0.0 %f]\n",  255.0);
1406
 
    fprintf (fp, "\t\t\t\t}\n");
1407
 
    fprintf (fp, "\t\t\t\tcoordIndex [\n");
1408
 
    fprintf (fp, "\t\t\t\t\t0, 1, -1\n");
1409
 
    fprintf (fp, "\t\t\t\t\t0, 2, -1\n");
1410
 
    fprintf (fp, "\t\t\t\t\t0, 3, -1]\n");
1411
 
    fprintf (fp, "\t\t\t}\n");
1412
 
    fprintf (fp, "\t\t}\n");
1413
 
 
1414
 
 
1415
 
    /* Draw the triangles as a shape: */
1416
 
    fprintf (fp, "\t\tShape {\n");
1417
 
    fprintf (fp, "\t\t\tappearance Appearance {\n");
1418
 
    fprintf (fp, "\t\t\t\tmaterial Material {\n");
1419
 
    fprintf (fp, "\t\t\t\t\tdiffuseColor 0 0.8 0\n");
1420
 
    fprintf (fp, "\t\t\t\t\temissiveColor 0 0 0\n");
1421
 
    fprintf (fp, "\t\t\t\t\tshininess 0.8\n");
1422
 
    fprintf (fp, "\t\t\t\t}\n");
1423
 
    fprintf (fp, "\t\t\t}\n");
1424
 
    fprintf (fp, "\t\t\tgeometry IndexedFaceSet {\n");
1425
 
    fprintf (fp, "\t\t\t\tsolid false\n");
1426
 
 
1427
 
    /* fill in the points here */
1428
 
    fprintf (fp, "\t\t\t\tcoord Coordinate {\n");
1429
 
    fprintf (fp, "\t\t\t\t\tpoint [\n");
1430
 
 
1431
 
    for (i = 0; i < hull->nvertex; ++i)
1432
 
    {
1433
 
        fprintf (fp, "\t\t\t\t\t%g %g %g%c\n",
1434
 
            (double) hull->Vertices[i][X], (double) hull->Vertices[i][Y], (double) hull->Vertices[i][Z],
1435
 
            i == hull->nvertex-1? ']': ',');
1436
 
    }
1437
 
    fprintf (fp, "\t\t\t\t}\n");
1438
 
 
1439
 
    /* fill in the Vertex indices (followed by -1) */
1440
 
 
1441
 
 
1442
 
    fprintf (fp, "\t\t\t\tcoordIndex [\n");
1443
 
    for (i = 0; i < hull->nfaces; ++i)
1444
 
    {
1445
 
        fprintf (fp, "\t\t\t\t\t%d, %d, %d, -1\n",
1446
 
            hull->Faces[i][0], hull->Faces[i][1], hull->Faces[i][2]);
1447
 
 
1448
 
    }
1449
 
    fprintf (fp, "]\n");
1450
 
 
1451
 
 
1452
 
    /* fill in the face colors */
1453
 
    fprintf (fp, "\t\t\t\tcolor Color {\n");
1454
 
    fprintf (fp, "\t\t\t\t\tcolor [\n");
1455
 
    for (i = 0; i < hull->nfaces; ++i)
1456
 
    {
1457
 
        int vx, vy, vz;
1458
 
        double r, g, b;
1459
 
 
1460
 
        vx = hull->Faces[i][0]; vy = hull->Faces[i][1]; vz = hull->Faces[i][2];
1461
 
        r = (double) (hull->Vertices[vx][X] + hull->Vertices[vy][X] + hull->Vertices[vz][X]) / (3* 255);
1462
 
        g = (double) (hull->Vertices[vx][Y] + hull->Vertices[vy][Y] + hull->Vertices[vz][Y]) / (3* 255);
1463
 
        b = (double) (hull->Vertices[vx][Z] + hull->Vertices[vy][Z] + hull->Vertices[vz][Z]) / (3* 255);
1464
 
 
1465
 
        fprintf (fp, "\t\t\t\t\t%g %g %g%c\n", r, g, b,
1466
 
            i == hull->nfaces-1? ']': ',');
1467
 
    }
1468
 
    fprintf (fp, "\t\t\t}\n");
1469
 
 
1470
 
    fprintf (fp, "\t\t\tcolorPerVertex false\n");
1471
 
 
1472
 
    fprintf (fp, "\t\t\t}\n");
1473
 
    fprintf (fp, "\t\t}\n");
1474
 
    fprintf (fp, "\t]\n");
1475
 
    fprintf (fp, "}\n");
1476
 
 
1477
 
    fclose (fp);
1478
 
 
1479
 
    return true;
1480
 
}