1
/*<html><pre> -<a href="qh-poly.htm"
2
>-------------------------------</a><a name="TOP">-</a>
5
implements polygons and simplices
7
see qh-poly.htm, poly.h and qhull.h
9
infrequent code is in poly2.c
10
(all but top 50 and their callers 12/3/95)
12
copyright (c) 1993-2003, The Geometry Center
17
/*======== functions in alphabetical order ==========*/
19
/*-<a href="qh-poly.htm#TOC"
20
>-------------------------------</a><a name="appendfacet">-</a>
22
qh_appendfacet( facet )
23
appends facet to end of qh.facet_list,
26
updates qh.newfacet_list, facet_next, facet_list
27
increments qh.numfacets
30
assumes qh.facet_list/facet_tail is defined (createsimplex)
36
void qh_appendfacet(facetT *facet) {
37
facetT *tail= qh facet_tail;
39
if (tail == qh newfacet_list)
40
qh newfacet_list= facet;
41
if (tail == qh facet_next)
43
facet->previous= tail->previous;
46
tail->previous->next= facet;
49
tail->previous= facet;
51
trace4((qh ferr, "qh_appendfacet: append f%d to facet_list\n", facet->id));
55
/*-<a href="qh-poly.htm#TOC"
56
>-------------------------------</a><a name="appendvertex">-</a>
58
qh_appendvertex( vertex )
59
appends vertex to end of qh.vertex_list,
63
updates qh.vertex_list, newvertex_list
64
increments qh.num_vertices
67
assumes qh.vertex_list/vertex_tail is defined (createsimplex)
70
void qh_appendvertex (vertexT *vertex) {
71
vertexT *tail= qh vertex_tail;
73
if (tail == qh newvertex_list)
74
qh newvertex_list= vertex;
75
vertex->newlist= True;
76
vertex->previous= tail->previous;
79
tail->previous->next= vertex;
81
qh vertex_list= vertex;
82
tail->previous= vertex;
84
trace4((qh ferr, "qh_appendvertex: append v%d to vertex_list\n", vertex->id));
88
/*-<a href="qh-poly.htm#TOC"
89
>-------------------------------</a><a name="attachnewfacets">-</a>
92
attach horizon facets to new facets in qh.newfacet_list
93
newfacets have neighbor and ridge links to horizon but not vice versa
94
only needed for qh.ONLYgood
98
horizon facets linked to new facets
99
ridges changed from visible facets to new facets
100
simplicial ridges deleted
101
qh.visible_list, no ridges valid
102
facet->f.replace is a newfacet (if any)
105
delete interior ridges and neighbor sets by
106
for each visible, non-simplicial facet
108
if last visit or if neighbor is simplicial
110
delete ridge for horizon's ridge set
113
attach horizon facets and new facets by
115
if corresponding horizon facet is simplicial
116
locate corresponding visible facet {may be more than one}
117
link visible facet to new facet
118
replace visible facet with new facet in horizon
119
else it's non-simplicial
120
for all visible neighbors of the horizon facet
121
link visible neighbor to new facet
122
delete visible neighbor from horizon facet
123
append new facet to horizon's neighbors
124
the first ridge of the new facet is the horizon ridge
125
link the new facet into the horizon ridge
127
void qh_attachnewfacets (void ) {
128
facetT *newfacet= NULL, *neighbor, **neighborp, *horizon, *visible;
129
ridgeT *ridge, **ridgep;
132
trace3((qh ferr, "qh_attachnewfacets: delete interior ridges\n"));
134
FORALLvisible_facets {
135
visible->visitid= qh visit_id;
136
if (visible->ridges) {
137
FOREACHridge_(visible->ridges) {
138
neighbor= otherfacet_(ridge, visible);
139
if (neighbor->visitid == qh visit_id
140
|| (!neighbor->visible && neighbor->simplicial)) {
141
if (!neighbor->visible) /* delete ridge for simplicial horizon */
142
qh_setdel (neighbor->ridges, ridge);
143
qh_setfree (&(ridge->vertices)); /* delete on 2nd visit */
144
qh_memfree (ridge, sizeof(ridgeT));
147
SETfirst_(visible->ridges)= NULL;
149
SETfirst_(visible->neighbors)= NULL;
151
trace1((qh ferr, "qh_attachnewfacets: attach horizon facets to new facets\n"));
153
horizon= SETfirstt_(newfacet->neighbors, facetT);
154
if (horizon->simplicial) {
156
FOREACHneighbor_(horizon) { /* may have more than one horizon ridge */
157
if (neighbor->visible) {
159
if (qh_setequal_skip (newfacet->vertices, 0, horizon->vertices,
160
SETindex_(horizon->neighbors, neighbor))) {
169
visible->f.replace= newfacet;
170
qh_setreplace (horizon->neighbors, visible, newfacet);
172
fprintf (qh ferr, "qhull internal error (qh_attachnewfacets): couldn't find visible facet for horizon f%d of newfacet f%d\n",
173
horizon->id, newfacet->id);
174
qh_errexit2 (qh_ERRqhull, horizon, newfacet);
176
}else { /* non-simplicial, with a ridge for newfacet */
177
FOREACHneighbor_(horizon) { /* may hold for many new facets */
178
if (neighbor->visible) {
179
neighbor->f.replace= newfacet;
180
qh_setdelnth (horizon->neighbors,
181
SETindex_(horizon->neighbors, neighbor));
182
neighborp--; /* repeat */
185
qh_setappend (&horizon->neighbors, newfacet);
186
ridge= SETfirstt_(newfacet->ridges, ridgeT);
187
if (ridge->top == horizon)
188
ridge->bottom= newfacet;
190
ridge->top= newfacet;
193
if (qh PRINTstatistics) {
194
FORALLvisible_facets {
195
if (!visible->f.replace)
196
zinc_(Zinsidevisible);
199
} /* attachnewfacets */
201
/*-<a href="qh-poly.htm#TOC"
202
>-------------------------------</a><a name="checkflipped">-</a>
204
qh_checkflipped( facet, dist, allerror )
205
checks facet orientation to interior point
208
tests against qh.DISTround
210
tests against 0 since tested against DISTround before
213
False if it flipped orientation (sets facet->flipped)
216
boolT qh_checkflipped (facetT *facet, realT *distp, boolT allerror) {
219
if (facet->flipped && !distp)
222
qh_distplane(qh interior_point, facet, &dist);
225
if ((allerror && dist > -qh DISTround)|| (!allerror && dist >= 0.0)) {
226
facet->flipped= True;
227
zzinc_(Zflippedfacets);
228
trace0((qh ferr, "qh_checkflipped: facet f%d is flipped, distance= %6.12g during p%d\n",
229
facet->id, dist, qh furthest_id));
230
qh_precision ("flipped facet");
236
/*-<a href="qh-poly.htm#TOC"
237
>-------------------------------</a><a name="delfacet">-</a>
240
removes facet from facet_list and frees up its memory
243
assumes vertices and ridges already freed
245
void qh_delfacet(facetT *facet) {
246
void **freelistp; /* used !qh_NOmem */
248
trace4((qh ferr, "qh_delfacet: delete f%d\n", facet->id));
249
if (facet == qh tracefacet)
251
if (facet == qh GOODclosest)
252
qh GOODclosest= NULL;
253
qh_removefacet(facet);
254
if (!facet->tricoplanar || facet->keepcentrum) {
255
qh_memfree_(facet->normal, qh normal_size, freelistp);
256
if (qh CENTERtype == qh_ASvoronoi) { /* uses macro calls */
257
qh_memfree_(facet->center, qh center_size, freelistp);
258
}else /* AScentrum */ {
259
qh_memfree_(facet->center, qh normal_size, freelistp);
262
qh_setfree(&(facet->neighbors));
264
qh_setfree(&(facet->ridges));
265
qh_setfree(&(facet->vertices));
266
if (facet->outsideset)
267
qh_setfree(&(facet->outsideset));
268
if (facet->coplanarset)
269
qh_setfree(&(facet->coplanarset));
270
qh_memfree_(facet, sizeof(facetT), freelistp);
274
/*-<a href="qh-poly.htm#TOC"
275
>-------------------------------</a><a name="deletevisible">-</a>
278
delete visible facets and vertices
281
deletes each facet and removes from facetlist
282
at exit, qh.visible_list empty (== qh.newfacet_list)
285
ridges already deleted
286
horizon facets do not reference facets on qh.visible_list
287
new facets in qh.newfacet_list
290
void qh_deletevisible (void /*qh visible_list*/) {
291
facetT *visible, *nextfacet;
292
vertexT *vertex, **vertexp;
293
int numvisible= 0, numdel= qh_setsize(qh del_vertices);
295
trace1((qh ferr, "qh_deletevisible: delete %d visible facets and %d vertices\n",
296
qh num_visible, numdel));
297
for (visible= qh visible_list; visible && visible->visible;
298
visible= nextfacet) { /* deleting current */
299
nextfacet= visible->next;
301
qh_delfacet(visible);
303
if (numvisible != qh num_visible) {
304
fprintf (qh ferr, "qhull internal error (qh_deletevisible): qh num_visible %d is not number of visible facets %d\n",
305
qh num_visible, numvisible);
306
qh_errexit (qh_ERRqhull, NULL, NULL);
309
zadd_(Zvisfacettot, numvisible);
310
zmax_(Zvisfacetmax, numvisible);
311
zzadd_(Zdelvertextot, numdel);
312
zmax_(Zdelvertexmax, numdel);
313
FOREACHvertex_(qh del_vertices)
314
qh_delvertex (vertex);
315
qh_settruncate (qh del_vertices, 0);
316
} /* deletevisible */
318
/*-<a href="qh-poly.htm#TOC"
319
>-------------------------------</a><a name="facetintersect">-</a>
321
qh_facetintersect( facetA, facetB, skipa, skipB, prepend )
322
return vertices for intersection of two simplicial facets
323
may include 1 prepended entry (if more, need to settemppush)
326
returns set of qh.hull_dim-1 + prepend vertices
327
returns skipped index for each test and checks for exactly one
330
does not need settemp since set in quick memory
333
qh_vertexintersect and qh_vertexintersect_new
334
use qh_setnew_delnthsorted to get nth ridge (no skip information)
337
locate skipped vertex by scanning facet A's neighbors
338
locate skipped vertex by scanning facet B's neighbors
339
intersect the vertex sets
341
setT *qh_facetintersect (facetT *facetA, facetT *facetB,
342
int *skipA,int *skipB, int prepend) {
344
int dim= qh hull_dim, i, j;
345
facetT **neighborsA, **neighborsB;
347
neighborsA= SETaddr_(facetA->neighbors, facetT);
348
neighborsB= SETaddr_(facetB->neighbors, facetT);
350
if (facetB == *neighborsA++)
352
else if (facetB == *neighborsA++)
354
else if (facetB == *neighborsA++)
357
for (i= 3; i < dim; i++) {
358
if (facetB == *neighborsA++) {
364
if (facetA == *neighborsB++)
366
else if (facetA == *neighborsB++)
368
else if (facetA == *neighborsB++)
371
for (j= 3; j < dim; j++) {
372
if (facetA == *neighborsB++) {
378
if (i >= dim || j >= dim) {
379
fprintf (qh ferr, "qhull internal error (qh_facetintersect): f%d or f%d not in others neighbors\n",
380
facetA->id, facetB->id);
381
qh_errexit2 (qh_ERRqhull, facetA, facetB);
383
intersect= qh_setnew_delnthsorted (facetA->vertices, qh hull_dim, *skipA, prepend);
384
trace4((qh ferr, "qh_facetintersect: f%d skip %d matches f%d skip %d\n",
385
facetA->id, *skipA, facetB->id, *skipB));
387
} /* facetintersect */
389
/*-<a href="qh-poly.htm#TOC"
390
>-------------------------------</a><a name="gethash">-</a>
392
qh_gethash( hashsize, set, size, firstindex, skipelem )
393
return hashvalue for a set with firstindex and skipelem
396
assumes at least firstindex+1 elements
397
assumes skipelem is NULL, in set, or part of hash
399
hashes memory addresses which may change over different runs of the same data
400
using sum for hash does badly in high d
402
unsigned qh_gethash (int hashsize, setT *set, int size, int firstindex, void *skipelem) {
403
void **elemp= SETelemaddr_(set, firstindex, void);
404
ptr_intT hash = 0, elem;
407
switch (size-firstindex) {
409
hash= (ptr_intT)(*elemp) - (ptr_intT) skipelem;
412
hash= (ptr_intT)(*elemp) + (ptr_intT)elemp[1] - (ptr_intT) skipelem;
415
hash= (ptr_intT)(*elemp) + (ptr_intT)elemp[1] + (ptr_intT)elemp[2]
416
- (ptr_intT) skipelem;
419
hash= (ptr_intT)(*elemp) + (ptr_intT)elemp[1] + (ptr_intT)elemp[2]
420
+ (ptr_intT)elemp[3] - (ptr_intT) skipelem;
423
hash= (ptr_intT)(*elemp) + (ptr_intT)elemp[1] + (ptr_intT)elemp[2]
424
+ (ptr_intT)elemp[3] + (ptr_intT)elemp[4] - (ptr_intT) skipelem;
427
hash= (ptr_intT)(*elemp) + (ptr_intT)elemp[1] + (ptr_intT)elemp[2]
428
+ (ptr_intT)elemp[3] + (ptr_intT)elemp[4]+ (ptr_intT)elemp[5]
429
- (ptr_intT) skipelem;
434
do { /* this is about 10% in 10-d */
435
if ((elem= (ptr_intT)*elemp++) != (ptr_intT)skipelem) {
436
hash ^= (elem << i) + (elem >> (32-i));
444
hash %= (ptr_intT) hashsize;
445
/* hash= 0; for debugging purposes */
449
/*-<a href="qh-poly.htm#TOC"
450
>-------------------------------</a><a name="makenewfacet">-</a>
452
qh_makenewfacet( vertices, toporient, horizon )
453
creates a toporient? facet from vertices
457
adds newfacet to qh.facet_list
458
newfacet->vertices= vertices
460
newfacet->neighbor= horizon, but not vice versa
461
newvertex_list updated with vertices
463
facetT *qh_makenewfacet(setT *vertices, boolT toporient,facetT *horizon) {
465
vertexT *vertex, **vertexp;
467
FOREACHvertex_(vertices) {
468
if (!vertex->newlist) {
469
qh_removevertex (vertex);
470
qh_appendvertex (vertex);
473
newfacet= qh_newfacet();
474
newfacet->vertices= vertices;
475
newfacet->toporient= toporient;
477
qh_setappend(&(newfacet->neighbors), horizon);
478
qh_appendfacet(newfacet);
483
/*-<a href="qh-poly.htm#TOC"
484
>-------------------------------</a><a name="makenewplanes">-</a>
487
make new hyperplanes for facets on qh.newfacet_list
490
all facets have hyperplanes or are marked for merging
491
doesn't create hyperplane if horizon is coplanar (will merge)
492
updates qh.min_vertex if qh.JOGGLEmax
495
facet->f.samecycle is defined for facet->mergehorizon facets
497
void qh_makenewplanes (void /* newfacet_list */) {
501
if (!newfacet->mergehorizon)
502
qh_setfacetplane (newfacet);
504
if (qh JOGGLEmax < REALmax/2)
505
minimize_(qh min_vertex, -wwval_(Wnewvertexmax));
506
} /* makenewplanes */
508
/*-<a href="qh-poly.htm#TOC"
509
>-------------------------------</a><a name="makenew_nonsimplicial">-</a>
511
qh_makenew_nonsimplicial( visible, apex, numnew )
512
make new facets for ridges of a visible facet
515
first newfacet, bumps numnew as needed
516
attaches new facets if !qh.ONLYgood
517
marks ridge neighbors for simplicial visible
519
ridges on newfacet, horizon, and visible
521
ridge and neighbors between newfacet and horizon
522
visible facet's ridges are deleted
525
qh.visit_id if visible has already been processed
526
sets neighbor->seen for building f.samecycle
527
assumes all 'seen' flags initially false
530
for each ridge of visible facet
531
get neighbor of visible facet
532
if neighbor was already processed
533
delete the ridge (will delete all visible facets later)
534
if neighbor is a horizon facet
537
adds newfacet to f.samecycle for later merging
539
updates neighbor's neighbor set
540
(checks for non-simplicial facet with multiple ridges to visible facet)
541
updates neighbor's ridge set
542
(checks for simplicial neighbor to non-simplicial visible facet)
543
(deletes ridge if neighbor is simplicial)
547
facetT *qh_makenew_nonsimplicial (facetT *visible, vertexT *apex, int *numnew) {
548
void **freelistp; /* used !qh_NOmem */
549
ridgeT *ridge, **ridgep;
550
facetT *neighbor, *newfacet= NULL, *samecycle;
555
FOREACHridge_(visible->ridges) {
557
neighbor= otherfacet_(ridge, visible);
558
if (neighbor->visible) {
560
if (neighbor->visitid == qh visit_id) {
561
qh_setfree (&(ridge->vertices)); /* delete on 2nd visit */
562
qh_memfree_(ridge, sizeof(ridgeT), freelistp);
565
}else { /* neighbor is an horizon facet */
566
toporient= (ridge->top == visible);
567
vertices= qh_setnew (qh hull_dim); /* makes sure this is quick */
568
qh_setappend (&vertices, apex);
569
qh_setappend_set (&vertices, ridge->vertices);
570
newfacet= qh_makenewfacet(vertices, toporient, neighbor);
572
if (neighbor->coplanar) {
573
newfacet->mergehorizon= True;
574
if (!neighbor->seen) {
575
newfacet->f.samecycle= newfacet;
576
neighbor->f.newcycle= newfacet;
578
samecycle= neighbor->f.newcycle;
579
newfacet->f.samecycle= samecycle->f.samecycle;
580
samecycle->f.samecycle= newfacet;
584
if (!neighbor->simplicial)
585
qh_setappend(&(newfacet->ridges), ridge);
586
}else { /* qh_attachnewfacets */
587
if (neighbor->seen) {
588
if (neighbor->simplicial) {
589
fprintf (qh ferr, "qhull internal error (qh_makenew_nonsimplicial): simplicial f%d sharing two ridges with f%d\n",
590
neighbor->id, visible->id);
591
qh_errexit2 (qh_ERRqhull, neighbor, visible);
593
qh_setappend (&(neighbor->neighbors), newfacet);
595
qh_setreplace (neighbor->neighbors, visible, newfacet);
596
if (neighbor->simplicial) {
597
qh_setdel (neighbor->ridges, ridge);
598
qh_setfree (&(ridge->vertices));
599
qh_memfree (ridge, sizeof(ridgeT));
601
qh_setappend(&(newfacet->ridges), ridge);
603
ridge->top= newfacet;
605
ridge->bottom= newfacet;
607
trace4((qh ferr, "qh_makenew_nonsimplicial: created facet f%d from v%d and r%d of horizon f%d\n",
608
newfacet->id, apex->id, ridgeid, neighbor->id));
611
neighbor->seen= True;
612
} /* for each ridge */
614
SETfirst_(visible->ridges)= NULL;
616
} /* makenew_nonsimplicial */
617
#else /* qh_NOmerge */
618
facetT *qh_makenew_nonsimplicial (facetT *visible, vertexT *apex, int *numnew) {
621
#endif /* qh_NOmerge */
623
/*-<a href="qh-poly.htm#TOC"
624
>-------------------------------</a><a name="makenew_simplicial">-</a>
626
qh_makenew_simplicial( visible, apex, numnew )
627
make new facets for simplicial visible facet and apex
630
attaches new facets if (!qh.ONLYgood)
631
neighbors between newfacet and horizon
634
nop if neighbor->seen or neighbor->visible (see qh_makenew_nonsimplicial)
637
locate neighboring horizon facet for visible facet
638
determine vertices and orientation
641
add new facet to f.samecycle
642
update horizon facet's neighbor list
644
facetT *qh_makenew_simplicial (facetT *visible, vertexT *apex, int *numnew) {
645
facetT *neighbor, **neighborp, *newfacet= NULL;
647
boolT flip, toporient;
648
int horizonskip, visibleskip;
650
FOREACHneighbor_(visible) {
651
if (!neighbor->seen && !neighbor->visible) {
652
vertices= qh_facetintersect(neighbor,visible, &horizonskip, &visibleskip, 1);
653
SETfirst_(vertices)= apex;
654
flip= ((horizonskip & 0x1) ^ (visibleskip & 0x1));
655
if (neighbor->toporient)
656
toporient= horizonskip & 0x1;
658
toporient= (horizonskip & 0x1) ^ 0x1;
659
newfacet= qh_makenewfacet(vertices, toporient, neighbor);
661
if (neighbor->coplanar && (qh PREmerge || qh MERGEexact)) {
663
newfacet->f.samecycle= newfacet;
664
newfacet->mergehorizon= True;
668
SETelem_(neighbor->neighbors, horizonskip)= newfacet;
669
trace4((qh ferr, "qh_makenew_simplicial: create facet f%d top %d from v%d and horizon f%d skip %d top %d and visible f%d skip %d, flip? %d\n",
670
newfacet->id, toporient, apex->id, neighbor->id, horizonskip,
671
neighbor->toporient, visible->id, visibleskip, flip));
675
} /* makenew_simplicial */
677
/*-<a href="qh-poly.htm#TOC"
678
>-------------------------------</a><a name="matchneighbor">-</a>
680
qh_matchneighbor( newfacet, newskip, hashsize, hashcount )
681
either match subridge of newfacet with neighbor or add to hash_table
684
duplicate ridges are unmatched and marked by qh_DUPLICATEridge
687
ridge is newfacet->vertices w/o newskip vertex
688
do not allocate memory (need to free hash_table cleanly)
689
uses linear hash chains
695
for each possible matching facet in qh.hash_table
697
set ismatch, if facets have opposite orientation
698
if ismatch and matching facet doesn't have a match
699
match the facets by updating their neighbor sets
701
indicate a duplicate ridge
702
set facet hyperplane for later testing
703
add facet to hashtable
704
unless the other facet was already a duplicate ridge
705
mark both facets with a duplicate ridge
706
add other facet (if defined) to hash table
708
void qh_matchneighbor (facetT *newfacet, int newskip, int hashsize, int *hashcount) {
709
boolT newfound= False; /* True, if new facet is already in hash chain */
712
facetT *facet, *matchfacet;
715
hash= (int)qh_gethash (hashsize, newfacet->vertices, qh hull_dim, 1,
716
SETelem_(newfacet->vertices, newskip));
717
trace4((qh ferr, "qh_matchneighbor: newfacet f%d skip %d hash %d hashcount %d\n",
718
newfacet->id, newskip, hash, *hashcount));
720
for (scan= hash; (facet= SETelemt_(qh hash_table, scan, facetT));
721
scan= (++scan >= hashsize ? 0 : scan)) {
722
if (facet == newfacet) {
727
if (qh_matchvertices (1, newfacet->vertices, newskip, facet->vertices, &skip, &same)) {
728
if (SETelem_(newfacet->vertices, newskip) ==
729
SETelem_(facet->vertices, skip)) {
730
qh_precision ("two facets with the same vertices");
731
fprintf (qh ferr, "qhull precision error: Vertex sets are the same for f%d and f%d. Can not force output.\n",
732
facet->id, newfacet->id);
733
qh_errexit2 (qh_ERRprec, facet, newfacet);
735
ismatch= (same == (newfacet->toporient ^ facet->toporient));
736
matchfacet= SETelemt_(facet->neighbors, skip, facetT);
737
if (ismatch && !matchfacet) {
738
SETelem_(facet->neighbors, skip)= newfacet;
739
SETelem_(newfacet->neighbors, newskip)= facet;
741
trace4((qh ferr, "qh_matchneighbor: f%d skip %d matched with new f%d skip %d\n",
742
facet->id, skip, newfacet->id, newskip));
745
if (!qh PREmerge && !qh MERGEexact) {
746
qh_precision ("a ridge with more than two neighbors");
747
fprintf (qh ferr, "qhull precision error: facets f%d, f%d and f%d meet at a ridge with more than 2 neighbors. Can not continue.\n",
748
facet->id, newfacet->id, getid_(matchfacet));
749
qh_errexit2 (qh_ERRprec, facet, newfacet);
751
SETelem_(newfacet->neighbors, newskip)= qh_DUPLICATEridge;
752
newfacet->dupridge= True;
753
if (!newfacet->normal)
754
qh_setfacetplane (newfacet);
755
qh_addhash (newfacet, qh hash_table, hashsize, hash);
758
qh_setfacetplane (facet);
759
if (matchfacet != qh_DUPLICATEridge) {
760
SETelem_(facet->neighbors, skip)= qh_DUPLICATEridge;
761
facet->dupridge= True;
763
qh_setfacetplane (facet);
765
matchskip= qh_setindex (matchfacet->neighbors, facet);
766
SETelem_(matchfacet->neighbors, matchskip)= qh_DUPLICATEridge;
767
matchfacet->dupridge= True;
768
if (!matchfacet->normal)
769
qh_setfacetplane (matchfacet);
770
qh_addhash (matchfacet, qh hash_table, hashsize, hash);
774
trace4((qh ferr, "qh_matchneighbor: new f%d skip %d duplicates ridge for f%d skip %d matching f%d ismatch %d at hash %d\n",
775
newfacet->id, newskip, facet->id, skip,
776
(matchfacet == qh_DUPLICATEridge ? -2 : getid_(matchfacet)),
778
return; /* end of duplicate ridge */
782
SETelem_(qh hash_table, scan)= newfacet; /* same as qh_addhash */
784
trace4((qh ferr, "qh_matchneighbor: no match for f%d skip %d at hash %d\n",
785
newfacet->id, newskip, hash));
786
} /* matchneighbor */
789
/*-<a href="qh-poly.htm#TOC"
790
>-------------------------------</a><a name="matchnewfacets">-</a>
793
match newfacets in qh.newfacet_list to their newfacet neighbors
796
qh.newfacet_list with full neighbor sets
797
get vertices with nth neighbor by deleting nth vertex
798
if qh.PREmerge/MERGEexact or qh.FORCEoutput
799
sets facet->flippped if flipped normal (also prevents point partitioning)
800
if duplicate ridges and qh.PREmerge/MERGEexact
802
missing neighbor links identifies extra ridges to be merging (qh_MERGEridge)
805
newfacets already have neighbor[0] (horizon facet)
806
assumes qh.hash_table is NULL
807
vertex->neighbors has not been updated yet
808
do not allocate memory after qh.hash_table (need to free it cleanly)
811
delete neighbor sets for all new facets
812
initialize a hash table
814
match facet with neighbors
815
if unmatched facets (due to duplicate ridges)
816
for each new facet with a duplicate ridge
817
match it with a facet
818
check for flipped facets
820
void qh_matchnewfacets (void /* qh newfacet_list */) {
821
int numnew=0, hashcount=0, newskip;
822
facetT *newfacet, *neighbor;
823
int dim= qh hull_dim, hashsize, neighbor_i, neighbor_n;
826
int facet_i, facet_n, numfree= 0;
830
trace1((qh ferr, "qh_matchnewfacets: match neighbors for new facets.\n"));
833
{ /* inline qh_setzero (newfacet->neighbors, 1, qh hull_dim); */
834
neighbors= newfacet->neighbors;
835
neighbors->e[neighbors->maxsize].i= dim+1; /*may be overwritten*/
836
memset ((char *)SETelemaddr_(neighbors, 1, void), 0, dim * SETelemsize);
839
qh_newhashtable (numnew*(qh hull_dim-1)); /* twice what is normally needed,
840
but every ridge could be DUPLICATEridge */
841
hashsize= qh_setsize (qh hash_table);
843
for (newskip=1; newskip<qh hull_dim; newskip++) /* furthest/horizon already matched */
844
qh_matchneighbor (newfacet, newskip, hashsize, &hashcount);
845
#if 0 /* use the following to trap hashcount errors */
848
facetT *facet, *neighbor;
851
FORALLfacet_(qh newfacet_list) { /* newfacet already in use */
852
for (k=1; k < qh hull_dim; k++) {
853
neighbor= SETelemt_(facet->neighbors, k, facetT);
854
if (!neighbor || neighbor == qh_DUPLICATEridge)
857
if (facet == newfacet)
860
if (count != hashcount) {
861
fprintf (qh ferr, "qh_matchnewfacets: after adding facet %d, hashcount %d != count %d\n",
862
newfacet->id, hashcount, count);
863
qh_errexit (qh_ERRqhull, newfacet, NULL);
866
#endif /* end of trap code */
870
if (newfacet->dupridge) {
871
FOREACHneighbor_i_(newfacet) {
872
if (neighbor == qh_DUPLICATEridge) {
873
qh_matchduplicates (newfacet, neighbor_i, hashsize, &hashcount);
874
/* this may report MERGEfacet */
881
fprintf (qh ferr, "qhull internal error (qh_matchnewfacets): %d neighbors did not match up\n",
883
qh_printhashtable (qh ferr);
884
qh_errexit (qh_ERRqhull, NULL, NULL);
887
if (qh IStracing >= 2) {
888
FOREACHfacet_i_(qh hash_table) {
892
fprintf (qh ferr, "qh_matchnewfacets: %d new facets, %d unused hash entries . hashsize %d\n",
893
numnew, numfree, qh_setsize (qh hash_table));
895
#endif /* !qh_NOtrace */
896
qh_setfree (&qh hash_table);
897
if (qh PREmerge || qh MERGEexact) {
898
if (qh IStracing >= 4)
899
qh_printfacetlist (qh newfacet_list, NULL, qh_ALL);
901
if (newfacet->normal)
902
qh_checkflipped (newfacet, NULL, qh_ALL);
904
}else if (qh FORCEoutput)
905
qh_checkflipped_all (qh newfacet_list); /* prints warnings for flipped */
906
} /* matchnewfacets */
909
/*-<a href="qh-poly.htm#TOC"
910
>-------------------------------</a><a name="matchvertices">-</a>
912
qh_matchvertices( firstindex, verticesA, skipA, verticesB, skipB, same )
913
tests whether vertices match with a single skip
914
starts match at firstindex since all new facets have a common vertex
917
true if matched vertices
918
skip index for each set
919
sets same iff vertices have the same orientation
922
assumes skipA is in A and both sets are the same size
926
scan both sets checking for a match
929
boolT qh_matchvertices (int firstindex, setT *verticesA, int skipA,
930
setT *verticesB, int *skipB, boolT *same) {
931
vertexT **elemAp, **elemBp, **skipBp=NULL, **skipAp;
933
elemAp= SETelemaddr_(verticesA, firstindex, vertexT);
934
elemBp= SETelemaddr_(verticesB, firstindex, vertexT);
935
skipAp= SETelemaddr_(verticesA, skipA, vertexT);
936
do if (elemAp != skipAp) {
937
while (*elemAp != *elemBp++) {
940
skipBp= elemBp; /* one extra like FOREACH */
945
*skipB= SETindex_(verticesB, skipB);
946
*same= !(((ptr_intT)skipA & 0x1) ^ ((ptr_intT)*skipB & 0x1));
947
trace4((qh ferr, "qh_matchvertices: matched by skip %d (v%d) and skip %d (v%d) same? %d\n",
948
skipA, (*skipAp)->id, *skipB, (*(skipBp-1))->id, *same));
950
} /* matchvertices */
952
/*-<a href="qh-poly.htm#TOC"
953
>-------------------------------</a><a name="newfacet">-</a>
959
all fields initialized or cleared (NULL)
960
preallocates neighbors set
962
facetT *qh_newfacet(void) {
964
void **freelistp; /* used !qh_NOmem */
966
qh_memalloc_(sizeof(facetT), freelistp, facet, facetT);
967
memset ((char *)facet, 0, sizeof(facetT));
968
if (qh facet_id == qh tracefacet_id)
969
qh tracefacet= facet;
970
facet->id= qh facet_id++;
971
facet->neighbors= qh_setnew(qh hull_dim);
972
#if !qh_COMPUTEfurthest
973
facet->furthestdist= 0.0;
976
if (qh FORCEoutput && qh APPROXhull)
977
facet->maxoutside= qh MINoutside;
979
facet->maxoutside= qh DISTround;
981
facet->simplicial= True;
983
facet->newfacet= True;
984
trace4((qh ferr, "qh_newfacet: created facet f%d\n", facet->id));
989
/*-<a href="qh-poly.htm#TOC"
990
>-------------------------------</a><a name="newridge">-</a>
995
ridgeT *qh_newridge(void) {
997
void **freelistp; /* used !qh_NOmem */
999
qh_memalloc_(sizeof(ridgeT), freelistp, ridge, ridgeT);
1000
memset ((char *)ridge, 0, sizeof(ridgeT));
1002
if (qh ridge_id == 0xFFFFFF) {
1004
qhull warning: more than %d ridges. ID field overflows and two ridges\n\
1005
may have the same identifier. Otherwise output ok.\n", 0xFFFFFF);
1007
ridge->id= qh ridge_id++;
1008
trace4((qh ferr, "qh_newridge: created ridge r%d\n", ridge->id));
1013
/*-<a href="qh-poly.htm#TOC"
1014
>-------------------------------</a><a name="pointid">-</a>
1017
return id for a point,
1018
returns -3 if null, -2 if interior, or -1 if not known
1022
id= ((unsigned long)point - (unsigned long)qh.first_point)/qh.normal_size;
1025
if point not in point array
1026
the code does a comparison of unrelated pointers.
1028
int qh_pointid (pointT *point) {
1033
else if (point == qh interior_point)
1035
else if (point >= qh first_point
1036
&& point < qh first_point + qh num_points * qh hull_dim) {
1037
offset= point - qh first_point;
1038
id= offset / qh hull_dim;
1039
}else if ((id= qh_setindex (qh other_points, point)) != -1)
1040
id += qh num_points;
1046
/*-<a href="qh-poly.htm#TOC"
1047
>-------------------------------</a><a name="removefacet">-</a>
1049
qh_removefacet( facet )
1050
unlinks facet from qh.facet_list,
1053
updates qh.facet_list .newfacet_list .facet_next visible_list
1054
decrements qh.num_facets
1059
void qh_removefacet(facetT *facet) {
1060
facetT *next= facet->next, *previous= facet->previous;
1062
if (facet == qh newfacet_list)
1063
qh newfacet_list= next;
1064
if (facet == qh facet_next)
1065
qh facet_next= next;
1066
if (facet == qh visible_list)
1067
qh visible_list= next;
1069
previous->next= next;
1070
next->previous= previous;
1071
}else { /* 1st facet in qh facet_list */
1072
qh facet_list= next;
1073
qh facet_list->previous= NULL;
1076
trace4((qh ferr, "qh_removefacet: remove f%d from facet_list\n", facet->id));
1080
/*-<a href="qh-poly.htm#TOC"
1081
>-------------------------------</a><a name="removevertex">-</a>
1083
qh_removevertex( vertex )
1084
unlinks vertex from qh.vertex_list,
1087
updates qh.vertex_list .newvertex_list
1088
decrements qh.num_vertices
1090
void qh_removevertex(vertexT *vertex) {
1091
vertexT *next= vertex->next, *previous= vertex->previous;
1093
if (vertex == qh newvertex_list)
1094
qh newvertex_list= next;
1096
previous->next= next;
1097
next->previous= previous;
1098
}else { /* 1st vertex in qh vertex_list */
1099
qh vertex_list= vertex->next;
1100
qh vertex_list->previous= NULL;
1103
trace4((qh ferr, "qh_removevertex: remove v%d from vertex_list\n", vertex->id));
1104
} /* removevertex */
1107
/*-<a href="qh-poly.htm#TOC"
1108
>-------------------------------</a><a name="updatevertices">-</a>
1111
update vertex neighbors and delete interior vertices
1114
if qh.VERTEXneighbors, updates neighbors for each vertex
1115
if qh.newvertex_list,
1116
removes visible neighbors from vertex neighbors
1118
adds new facets to vertex neighbors
1120
interior vertices added to qh.del_vertices for later partitioning
1123
if qh.VERTEXneighbors
1124
deletes references to visible facets from vertex neighbors
1125
appends new facets to the neighbor list for each vertex
1126
checks all vertices of visible facets
1127
removes visible facets from neighbor lists
1128
marks unused vertices for deletion
1130
void qh_updatevertices (void /*qh newvertex_list, newfacet_list, visible_list*/) {
1131
facetT *newfacet= NULL, *neighbor, **neighborp, *visible;
1132
vertexT *vertex, **vertexp;
1134
trace3((qh ferr, "qh_updatevertices: delete interior vertices and update vertex->neighbors\n"));
1135
if (qh VERTEXneighbors) {
1136
FORALLvertex_(qh newvertex_list) {
1137
FOREACHneighbor_(vertex) {
1138
if (neighbor->visible)
1139
SETref_(neighbor)= NULL;
1141
qh_setcompact (vertex->neighbors);
1144
FOREACHvertex_(newfacet->vertices)
1145
qh_setappend (&vertex->neighbors, newfacet);
1147
FORALLvisible_facets {
1148
FOREACHvertex_(visible->vertices) {
1149
if (!vertex->newlist && !vertex->deleted) {
1150
FOREACHneighbor_(vertex) { /* this can happen under merging */
1151
if (!neighbor->visible)
1155
qh_setdel (vertex->neighbors, visible);
1157
vertex->deleted= True;
1158
qh_setappend (&qh del_vertices, vertex);
1159
trace2((qh ferr, "qh_updatevertices: delete vertex p%d (v%d) in f%d\n",
1160
qh_pointid(vertex->point), vertex->id, visible->id));
1165
}else { /* !VERTEXneighbors */
1166
FORALLvisible_facets {
1167
FOREACHvertex_(visible->vertices) {
1168
if (!vertex->newlist && !vertex->deleted) {
1169
vertex->deleted= True;
1170
qh_setappend (&qh del_vertices, vertex);
1171
trace2((qh ferr, "qh_updatevertices: delete vertex p%d (v%d) in f%d\n",
1172
qh_pointid(vertex->point), vertex->id, visible->id));
1177
} /* updatevertices */