Details | Last modification | View Log | RSS feed
Rev | Author | Line No. | Line |
---|---|---|---|
2454 | craig | 1 | // |
2 | // Little cms - profiler construction set |
||
3 | // Copyright (C) 1998-2001 Marti Maria |
||
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., 59 Temple Place - Suite 330, Boston, MA 02111-1307, 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 | } |