1
/*<html><pre> -<a href="qh-qhull.htm"
2
>-------------------------------</a><a name="TOP">-</a>
5
user-level header file for using qhull.a library
7
see qh-qhull.htm, qhull_a.h
9
Copyright (c) 1993-2012 The Geometry Center.
10
$Id: //main/2011/qhull/src/libqhull/libqhull.h#7 $$Change: 1464 $
11
$DateTime: 2012/01/25 22:58:41 $$Author: bbarber $
13
NOTE: access to qh_qh is via the 'qh' macro. This allows
14
qh_qh to be either a pointer or a structure. An example
15
of using qh is "qh DROPdim" which accesses the DROPdim
16
field of qh_qh. Similarly, access to qh_qhstat is via
19
includes function prototypes for libqhull.c, geom.c, global.c, io.c, user.c
24
see unix.c for an example of using libqhull.h
26
recompile qhull if you change this file
30
#define qhDEFlibqhull 1
32
/*=========================== -included files ==============*/
34
#include "user.h" /* user definable constants (e.g., qh_QHpointer) */
41
#if __MWERKS__ && __POWERPC__
50
#error Neither __STDC__ nor __cplusplus is defined. Please use strict ANSI C or C++ to compile
51
#error Qhull. You may need to turn off compiler extensions in your project configuration. If
52
#error your compiler is a standard C compiler, you can delete this warning from libqhull.h
57
/*============ constants and basic types ====================*/
59
extern const char *qh_version; /* defined in global.c */
61
/*-<a href="qh-geom.htm#TOC"
62
>--------------------------------</a><a name="coordT">-</a>
65
coordinates and coefficients are stored as realT (i.e., double)
68
Qhull works well if realT is 'float'. If so joggle (QJ) is not effective.
70
Could use 'float' for data and 'double' for calculations (realT vs. coordT)
71
This requires many type casts, and adjusted error bounds.
72
Also C compilers may do expressions in double anyway.
76
/*-<a href="qh-geom.htm#TOC"
77
>--------------------------------</a><a name="pointT">-</a>
80
a point is an array of coordinates, usually qh.hull_dim
84
/*-<a href="qh-qhull.htm#TOC"
85
>--------------------------------</a><a name="flagT">-</a>
90
#define flagT unsigned int
92
/*-<a href="qh-qhull.htm#TOC"
93
>--------------------------------</a><a name="boolT">-</a>
96
boolean value, either True or False
99
needed for portability
100
Use qh_False/qh_True as synonyms
102
#define boolT unsigned int
114
/*-<a href="qh-qhull.htm#TOC"
115
>--------------------------------</a><a name="CENTERtype">-</a>
118
to distinguish facet->center
122
qh_ASnone = 0, qh_ASvoronoi, qh_AScentrum
126
/*-<a href="qh-qhull.htm#TOC"
127
>--------------------------------</a><a name="qh_PRINT">-</a>
130
output formats for printing (qh.PRINTout).
135
some of these names are similar to qh names. The similar names are only
136
used in switch statements in qh_printbegin() etc.
138
typedef enum {qh_PRINTnone= 0,
139
qh_PRINTarea, qh_PRINTaverage, /* 'Fa' 'FV' 'Fc' 'FC' */
140
qh_PRINTcoplanars, qh_PRINTcentrums,
141
qh_PRINTfacets, qh_PRINTfacets_xridge, /* 'f' 'FF' 'G' 'FI' 'Fi' 'Fn' */
142
qh_PRINTgeom, qh_PRINTids, qh_PRINTinner, qh_PRINTneighbors,
143
qh_PRINTnormals, qh_PRINTouter, qh_PRINTmaple, /* 'n' 'Fo' 'i' 'm' 'Fm' 'FM', 'o' */
144
qh_PRINTincidences, qh_PRINTmathematica, qh_PRINTmerges, qh_PRINToff,
145
qh_PRINToptions, qh_PRINTpointintersect, /* 'FO' 'Fp' 'FP' 'p' 'FQ' 'FS' */
146
qh_PRINTpointnearest, qh_PRINTpoints, qh_PRINTqhull, qh_PRINTsize,
147
qh_PRINTsummary, qh_PRINTtriangles, /* 'Fs' 'Ft' 'Fv' 'FN' 'Fx' */
148
qh_PRINTvertices, qh_PRINTvneighbors, qh_PRINTextremes,
149
qh_PRINTEND} qh_PRINT;
151
/*-<a href="qh-qhull.htm#TOC"
152
>--------------------------------</a><a name="qh_ALL">-</a>
155
argument flag for selecting everything
158
#define qh_NOupper True /* argument for qh_findbest */
159
#define qh_IScheckmax True /* argument for qh_findbesthorizon */
160
#define qh_ISnewfacets True /* argument for qh_findbest */
161
#define qh_RESETvisible True /* argument for qh_resetlists */
163
/*-<a href="qh-qhull.htm#TOC"
164
>--------------------------------</a><a name="qh_ERR">-</a>
167
Qhull exit codes, for indicating errors
168
See: MSG_ERROR and MSG_WARNING [user.h]
170
#define qh_ERRnone 0 /* no error occurred during qhull */
171
#define qh_ERRinput 1 /* input inconsistency */
172
#define qh_ERRsingular 2 /* singular input data */
173
#define qh_ERRprec 3 /* precision error */
174
#define qh_ERRmem 4 /* insufficient memory, matches mem.h */
175
#define qh_ERRqhull 5 /* internal error detected, matches mem.h */
177
/*-<a href="qh-qhull.htm#TOC"
178
>--------------------------------</a><a name="qh_FILEstderr">-</a>
181
Fake stderr to distinguish error output from normal output
182
For C++ interface. Must redefine qh_fprintf_qhull
184
#define qh_FILEstderr (FILE*)1
186
/* ============ -structures- ====================
187
each of the following structures is defined by a typedef
188
all realT and coordT fields occur at the beginning of a structure
189
(otherwise space may be wasted due to alignment)
190
define all flags together and pack into 32-bit number
193
typedef struct vertexT vertexT;
194
typedef struct ridgeT ridgeT;
195
typedef struct facetT facetT;
198
typedef struct setT setT; /* defined in qset.h */
203
typedef struct qhstatT qhstatT; /* defined in stat.h */
206
/*-<a href="qh-poly.htm#TOC"
207
>--------------------------------</a><a name="facetT">-</a>
213
qhull() generates the hull as a list of facets.
215
topological information:
216
f.previous,next doubly-linked list of facets
217
f.vertices set of vertices
218
f.ridges set of ridges
219
f.neighbors set of neighbors
220
f.toporient True if facet has top-orientation (else bottom)
222
geometric information:
223
f.offset,normal hyperplane equation
224
f.maxoutside offset to outer plane -- all points inside
225
f.center centrum for testing convexity
226
f.simplicial True if facet is simplicial
227
f.flipped True if facet does not include qh.interior_point
229
for constructing hull:
230
f.visible True if facet on list of visible facets (will be deleted)
231
f.newfacet True if facet on list of newly created facets
232
f.coplanarset set of points coplanar with this facet
233
(includes near-inside points for later testing)
234
f.outsideset set of points outside of this facet
235
f.furthestdist distance to furthest point of outside set
236
f.visitid marks visited facets during a loop
237
f.replace replacement facet for to-be-deleted, visible facets
238
f.samecycle,newcycle cycle of facets for merging into horizon facet
240
see below for other flags and fields
243
#if !qh_COMPUTEfurthest
244
coordT furthestdist;/* distance to furthest point of outsideset */
247
coordT maxoutside; /* max computed distance of point to facet
248
Before QHULLfinished this is an approximation
249
since maxdist not always set for mergefacet
250
Actual outer plane is +DISTround and
251
computed outer plane is +2*DISTround */
253
coordT offset; /* exact offset of hyperplane from origin */
254
coordT *normal; /* normal of hyperplane, hull_dim coefficients */
255
/* if tricoplanar, shared with a neighbor */
256
union { /* in order of testing */
257
realT area; /* area of facet, only in io.c if ->isarea */
258
facetT *replace; /* replacement facet if ->visible and NEWfacets
259
is NULL only if qh_mergedegen_redundant or interior */
260
facetT *samecycle; /* cycle of facets from the same visible/horizon intersection,
262
facetT *newcycle; /* in horizon facet, current samecycle of new facets */
263
facetT *trivisible; /* visible facet for ->tricoplanar facets during qh_triangulate() */
264
facetT *triowner; /* owner facet for ->tricoplanar, !isarea facets w/ ->keepcentrum */
266
coordT *center; /* centrum for convexity, qh CENTERtype == qh_AScentrum */
267
/* Voronoi center, qh CENTERtype == qh_ASvoronoi */
268
/* if tricoplanar, shared with a neighbor */
269
facetT *previous; /* previous facet in the facet_list */
270
facetT *next; /* next facet in the facet_list */
271
setT *vertices; /* vertices for this facet, inverse sorted by ID
272
if simplicial, 1st vertex was apex/furthest */
273
setT *ridges; /* explicit ridges for nonsimplicial facets.
274
for simplicial facets, neighbors define the ridges */
275
setT *neighbors; /* neighbors of the facet. If simplicial, the kth
276
neighbor is opposite the kth vertex, and the first
277
neighbor is the horizon facet for the first vertex*/
278
setT *outsideset; /* set of points outside this facet
279
if non-empty, last point is furthest
280
if NARROWhull, includes coplanars for partitioning*/
281
setT *coplanarset; /* set of points coplanar with this facet
282
> qh.min_vertex and <= facet->max_outside
283
a point is assigned to the furthest facet
284
if non-empty, last point is furthest away */
285
unsigned visitid; /* visit_id, for visiting all neighbors,
286
all uses are independent */
287
unsigned id; /* unique identifier from qh facet_id */
288
unsigned nummerge:9; /* number of merges */
289
#define qh_MAXnummerge 511 /* 2^9-1, 32 flags total, see "flags:" in io.c */
290
flagT tricoplanar:1; /* True if TRIangulate and simplicial and coplanar with a neighbor */
291
/* all tricoplanars share the same ->center, ->normal, ->offset, ->maxoutside */
292
/* all tricoplanars share the same apex */
293
/* if ->degenerate, does not span facet (one logical ridge) */
294
/* one tricoplanar has ->keepcentrum and ->coplanarset */
295
/* during qh_triangulate, f.trivisible points to original facet */
296
flagT newfacet:1; /* True if facet on qh newfacet_list (new or merged) */
297
flagT visible:1; /* True if visible facet (will be deleted) */
298
flagT toporient:1; /* True if created with top orientation
299
after merging, use ridge orientation */
300
flagT simplicial:1;/* True if simplicial facet, ->ridges may be implicit */
301
flagT seen:1; /* used to perform operations only once, like visitid */
302
flagT seen2:1; /* used to perform operations only once, like visitid */
303
flagT flipped:1; /* True if facet is flipped */
304
flagT upperdelaunay:1; /* True if facet is upper envelope of Delaunay triangulation */
305
flagT notfurthest:1; /* True if last point of outsideset is not furthest*/
307
/*-------- flags primarily for output ---------*/
308
flagT good:1; /* True if a facet marked good for output */
309
flagT isarea:1; /* True if facet->f.area is defined */
311
/*-------- flags for merging ------------------*/
312
flagT dupridge:1; /* True if duplicate ridge in facet */
313
flagT mergeridge:1; /* True if facet or neighbor contains a qh_MERGEridge
314
->normal defined (also defined for mergeridge2) */
315
flagT mergeridge2:1; /* True if neighbor contains a qh_MERGEridge (mark_dupridges */
316
flagT coplanar:1; /* True if horizon facet is coplanar at last use */
317
flagT mergehorizon:1; /* True if will merge into horizon (->coplanar) */
318
flagT cycledone:1;/* True if mergecycle_all already done */
319
flagT tested:1; /* True if facet convexity has been tested (false after merge */
320
flagT keepcentrum:1; /* True if keep old centrum after a merge, or marks owner for ->tricoplanar */
321
flagT newmerge:1; /* True if facet is newly merged for reducevertices */
322
flagT degenerate:1; /* True if facet is degenerate (degen_mergeset or ->tricoplanar) */
323
flagT redundant:1; /* True if facet is redundant (degen_mergeset) */
327
/*-<a href="qh-poly.htm#TOC"
328
>--------------------------------</a><a name="ridgeT">-</a>
334
a ridge is hull_dim-1 simplex between two neighboring facets. If the
335
facets are non-simplicial, there may be more than one ridge between
336
two facets. E.G. a 4-d hypercube has two triangles between each pair
337
of neighboring facets.
339
topological information:
340
vertices a set of vertices
341
top,bottom neighboring facets with orientation
343
geometric information:
344
tested True if ridge is clearly convex
345
nonconvex True if ridge is non-convex
348
setT *vertices; /* vertices belonging to this ridge, inverse sorted by ID
349
NULL if a degen ridge (matchsame) */
350
facetT *top; /* top facet this ridge is part of */
351
facetT *bottom; /* bottom facet this ridge is part of */
352
unsigned id:24; /* unique identifier, =>room for 8 flags, bit field matches qh.ridge_id */
353
flagT seen:1; /* used to perform operations only once */
354
flagT tested:1; /* True when ridge is tested for convexity */
355
flagT nonconvex:1; /* True if getmergeset detected a non-convex neighbor
356
only one ridge between neighbors may have nonconvex */
359
/*-<a href="qh-poly.htm#TOC"
360
>--------------------------------</a><a name="vertexT">-</a>
365
topological information:
366
next,previous doubly-linked list of all vertices
367
neighbors set of adjacent facets (only if qh.VERTEXneighbors)
369
geometric information:
370
point array of DIM3 coordinates
373
vertexT *next; /* next vertex in vertex_list */
374
vertexT *previous; /* previous vertex in vertex_list */
375
pointT *point; /* hull_dim coordinates (coordT) */
376
setT *neighbors; /* neighboring facets of vertex, qh_vertexneighbors()
377
inits in io.c or after first merge */
378
unsigned visitid:31; /* for use with qh vertex_visit, size must match */
379
flagT seen2:1; /* another seen flag */
380
unsigned id:24; /* unique identifier, bit field matches qh.vertex_id */
381
unsigned dim:4; /* dimension of point if non-zero, used by cpp */
382
/* =>room for 4 flags */
383
flagT seen:1; /* used to perform operations only once */
384
flagT delridge:1; /* vertex was part of a deleted ridge */
385
flagT deleted:1; /* true if vertex on qh del_vertices */
386
flagT newlist:1; /* true if vertex on qh newvertex_list */
389
#define MAX_vdim 15 /* Maximum size of vertex->dim */
391
/*======= -global variables -qh ============================*/
393
/*-<a href="qh-globa.htm#TOC"
394
>--------------------------------</a><a name="qh">-</a>
397
all global variables for qhull are in qh, qhmem, and qhstat
400
qhmem is defined in mem.h, qhstat is defined in stat.h, qhrbox is defined in rboxpoints.h
401
Access to qh_qh is via the "qh" macro. See qh_QHpointer in user.h
403
All global variables for qhull are in qh, qhmem, and qhstat
404
qh must be unique for each instance of qhull
405
qhstat may be shared between qhull instances.
406
qhmem may be shared across multiple instances of Qhull.
407
Rbox uses global variables rbox_inuse and rbox, but does not persist data across calls.
409
Qhull is not multithreaded. Global state could be stored in thread-local storage.
412
extern int qhull_inuse;
414
typedef struct qhT qhT;
415
#if qh_QHpointer_dllimport
417
__declspec(dllimport) extern qhT *qh_qh; /* allocated in global.c */
420
extern qhT *qh_qh; /* allocated in global.c */
423
__declspec(dllimport) extern qhT qh_qh; /* allocated in global.c */
431
/*-<a href="qh-globa.htm#TOC"
432
>--------------------------------</a><a name="qh-const">-</a>
435
configuration flags and constants for Qhull
438
The user configures Qhull by defining flags. They are
439
copied into qh by qh_setflags(). qh-quick.htm#options defines the flags.
441
boolT ALLpoints; /* true 'Qs' if search all points for initial simplex */
442
boolT ANGLEmerge; /* true 'Qa' if sort potential merges by angle */
443
boolT APPROXhull; /* true 'Wn' if MINoutside set */
444
realT MINoutside; /* 'Wn' min. distance for an outside point */
445
boolT ANNOTATEoutput; /* true 'Ta' if annotate output with message codes */
446
boolT ATinfinity; /* true 'Qz' if point num_points-1 is "at-infinity"
447
for improving precision in Delaunay triangulations */
448
boolT AVOIDold; /* true 'Q4' if avoid old->new merges */
449
boolT BESToutside; /* true 'Qf' if partition points into best outsideset */
450
boolT CDDinput; /* true 'Pc' if input uses CDD format (1.0/offset first) */
451
boolT CDDoutput; /* true 'PC' if print normals in CDD format (offset first) */
452
boolT CHECKfrequently; /* true 'Tc' if checking frequently */
453
realT premerge_cos; /* 'A-n' cos_max when pre merging */
454
realT postmerge_cos; /* 'An' cos_max when post merging */
455
boolT DELAUNAY; /* true 'd' if computing DELAUNAY triangulation */
456
boolT DOintersections; /* true 'Gh' if print hyperplane intersections */
457
int DROPdim; /* drops dim 'GDn' for 4-d -> 3-d output */
458
boolT FORCEoutput; /* true 'Po' if forcing output despite degeneracies */
459
int GOODpoint; /* 1+n for 'QGn', good facet if visible/not(-) from point n*/
460
pointT *GOODpointp; /* the actual point */
461
boolT GOODthreshold; /* true if qh lower_threshold/upper_threshold defined
462
false if qh SPLITthreshold */
463
int GOODvertex; /* 1+n, good facet if vertex for point n */
464
pointT *GOODvertexp; /* the actual point */
465
boolT HALFspace; /* true 'Hn,n,n' if halfspace intersection */
466
int IStracing; /* trace execution, 0=none, 1=least, 4=most, -1=events */
467
int KEEParea; /* 'PAn' number of largest facets to keep */
468
boolT KEEPcoplanar; /* true 'Qc' if keeping nearest facet for coplanar points */
469
boolT KEEPinside; /* true 'Qi' if keeping nearest facet for inside points
470
set automatically if 'd Qc' */
471
int KEEPmerge; /* 'PMn' number of facets to keep with most merges */
472
realT KEEPminArea; /* 'PFn' minimum facet area to keep */
473
realT MAXcoplanar; /* 'Un' max distance below a facet to be coplanar*/
474
boolT MERGEexact; /* true 'Qx' if exact merges (coplanar, degen, dupridge, flipped) */
475
boolT MERGEindependent; /* true 'Q2' if merging independent sets */
476
boolT MERGING; /* true if exact-, pre- or post-merging, with angle and centrum tests */
477
realT premerge_centrum; /* 'C-n' centrum_radius when pre merging. Default is round-off */
478
realT postmerge_centrum; /* 'Cn' centrum_radius when post merging. Default is round-off */
479
boolT MERGEvertices; /* true 'Q3' if merging redundant vertices */
480
realT MINvisible; /* 'Vn' min. distance for a facet to be visible */
481
boolT NOnarrow; /* true 'Q10' if no special processing for narrow distributions */
482
boolT NOnearinside; /* true 'Q8' if ignore near-inside points when partitioning */
483
boolT NOpremerge; /* true 'Q0' if no defaults for C-0 or Qx */
484
boolT ONLYgood; /* true 'Qg' if process points with good visible or horizon facets */
485
boolT ONLYmax; /* true 'Qm' if only process points that increase max_outside */
486
boolT PICKfurthest; /* true 'Q9' if process furthest of furthest points*/
487
boolT POSTmerge; /* true if merging after buildhull (Cn or An) */
488
boolT PREmerge; /* true if merging during buildhull (C-n or A-n) */
489
/* NOTE: some of these names are similar to qh_PRINT names */
490
boolT PRINTcentrums; /* true 'Gc' if printing centrums */
491
boolT PRINTcoplanar; /* true 'Gp' if printing coplanar points */
492
int PRINTdim; /* print dimension for Geomview output */
493
boolT PRINTdots; /* true 'Ga' if printing all points as dots */
494
boolT PRINTgood; /* true 'Pg' if printing good facets */
495
boolT PRINTinner; /* true 'Gi' if printing inner planes */
496
boolT PRINTneighbors; /* true 'PG' if printing neighbors of good facets */
497
boolT PRINTnoplanes; /* true 'Gn' if printing no planes */
498
boolT PRINToptions1st; /* true 'FO' if printing options to stderr */
499
boolT PRINTouter; /* true 'Go' if printing outer planes */
500
boolT PRINTprecision; /* false 'Pp' if not reporting precision problems */
501
qh_PRINT PRINTout[qh_PRINTEND]; /* list of output formats to print */
502
boolT PRINTridges; /* true 'Gr' if print ridges */
503
boolT PRINTspheres; /* true 'Gv' if print vertices as spheres */
504
boolT PRINTstatistics; /* true 'Ts' if printing statistics to stderr */
505
boolT PRINTsummary; /* true 's' if printing summary to stderr */
506
boolT PRINTtransparent; /* true 'Gt' if print transparent outer ridges */
507
boolT PROJECTdelaunay; /* true if DELAUNAY, no readpoints() and
508
need projectinput() for Delaunay in qh_init_B */
509
int PROJECTinput; /* number of projected dimensions 'bn:0Bn:0' */
510
boolT QUICKhelp; /* true if quick help message for degen input */
511
boolT RANDOMdist; /* true if randomly change distplane and setfacetplane */
512
realT RANDOMfactor; /* maximum random perturbation */
513
realT RANDOMa; /* qh_randomfactor is randr * RANDOMa + RANDOMb */
515
boolT RANDOMoutside; /* true if select a random outside point */
516
int REPORTfreq; /* buildtracing reports every n facets */
517
int REPORTfreq2; /* tracemerging reports every REPORTfreq/2 facets */
518
int RERUN; /* 'TRn' rerun qhull n times (qh.build_cnt) */
519
int ROTATErandom; /* 'QRn' seed, 0 time, >= rotate input */
520
boolT SCALEinput; /* true 'Qbk' if scaling input */
521
boolT SCALElast; /* true 'Qbb' if scale last coord to max prev coord */
522
boolT SETroundoff; /* true 'E' if qh DISTround is predefined */
523
boolT SKIPcheckmax; /* true 'Q5' if skip qh_check_maxout */
524
boolT SKIPconvex; /* true 'Q6' if skip convexity testing during pre-merge */
525
boolT SPLITthresholds; /* true if upper_/lower_threshold defines a region
526
used only for printing (!for qh ONLYgood) */
527
int STOPcone; /* 'TCn' 1+n for stopping after cone for point n */
528
/* also used by qh_build_withresart for err exit*/
529
int STOPpoint; /* 'TVn' 'TV-n' 1+n for stopping after/before(-)
531
int TESTpoints; /* 'QTn' num of test points after qh.num_points. Test points always coplanar. */
532
boolT TESTvneighbors; /* true 'Qv' if test vertex neighbors at end */
533
int TRACElevel; /* 'Tn' conditional IStracing level */
534
int TRACElastrun; /* qh.TRACElevel applies to last qh.RERUN */
535
int TRACEpoint; /* 'TPn' start tracing when point n is a vertex */
536
realT TRACEdist; /* 'TWn' start tracing when merge distance too big */
537
int TRACEmerge; /* 'TMn' start tracing before this merge */
538
boolT TRIangulate; /* true 'Qt' if triangulate non-simplicial facets */
539
boolT TRInormals; /* true 'Q11' if triangulate duplicates normals (sets Qt) */
540
boolT UPPERdelaunay; /* true 'Qu' if computing furthest-site Delaunay */
541
boolT USEstdout; /* true 'Tz' if using stdout instead of stderr */
542
boolT VERIFYoutput; /* true 'Tv' if verify output at end of qhull */
543
boolT VIRTUALmemory; /* true 'Q7' if depth-first processing in buildhull */
544
boolT VORONOI; /* true 'v' if computing Voronoi diagram */
546
/*--------input constants ---------*/
547
realT AREAfactor; /* 1/(hull_dim-1)! for converting det's to area */
548
boolT DOcheckmax; /* true if calling qh_check_maxout (qh_initqhull_globals) */
549
char *feasible_string; /* feasible point 'Hn,n,n' for halfspace intersection */
550
coordT *feasible_point; /* as coordinates, both malloc'd */
551
boolT GETarea; /* true 'Fa', 'FA', 'FS', 'PAn', 'PFn' if compute facet area/Voronoi volume in io.c */
552
boolT KEEPnearinside; /* true if near-inside points in coplanarset */
553
int hull_dim; /* dimension of hull, set by initbuffers */
554
int input_dim; /* dimension of input, set by initbuffers */
555
int num_points; /* number of input points */
556
pointT *first_point; /* array of input points, see POINTSmalloc */
557
boolT POINTSmalloc; /* true if qh first_point/num_points allocated */
558
pointT *input_points; /* copy of original qh.first_point for input points for qh_joggleinput */
559
boolT input_malloc; /* true if qh input_points malloc'd */
560
char qhull_command[256];/* command line that invoked this program */
561
int qhull_commandsiz2; /* size of qhull_command at qh_clear_outputflags */
562
char rbox_command[256]; /* command line that produced the input points */
563
char qhull_options[512];/* descriptive list of options */
564
int qhull_optionlen; /* length of last line */
565
int qhull_optionsiz; /* size of qhull_options at qh_build_withrestart */
566
int qhull_optionsiz2; /* size of qhull_options at qh_clear_outputflags */
567
int run_id; /* non-zero, random identifier for this instance of qhull */
568
boolT VERTEXneighbors; /* true if maintaining vertex neighbors */
569
boolT ZEROcentrum; /* true if 'C-0' or 'C-0 Qx'. sets ZEROall_ok */
570
realT *upper_threshold; /* don't print if facet->normal[k]>=upper_threshold[k]
571
must set either GOODthreshold or SPLITthreshold
572
if Delaunay, default is 0.0 for upper envelope */
573
realT *lower_threshold; /* don't print if facet->normal[k] <=lower_threshold[k] */
574
realT *upper_bound; /* scale point[k] to new upper bound */
575
realT *lower_bound; /* scale point[k] to new lower bound
576
project if both upper_ and lower_bound == 0 */
578
/*-<a href="qh-globa.htm#TOC"
579
>--------------------------------</a><a name="qh-prec">-</a>
581
qh precision constants
582
precision constants for Qhull
585
qh_detroundoff() computes the maximum roundoff error for distance
586
and other computations. It also sets default values for the
589
realT ANGLEround; /* max round off error for angles */
590
realT centrum_radius; /* max centrum radius for convexity (roundoff added) */
591
realT cos_max; /* max cosine for convexity (roundoff added) */
592
realT DISTround; /* max round off error for distances, 'E' overrides */
593
realT MAXabs_coord; /* max absolute coordinate */
594
realT MAXlastcoord; /* max last coordinate for qh_scalelast */
595
realT MAXsumcoord; /* max sum of coordinates */
596
realT MAXwidth; /* max rectilinear width of point coordinates */
597
realT MINdenom_1; /* min. abs. value for 1/x */
598
realT MINdenom; /* use divzero if denominator < MINdenom */
599
realT MINdenom_1_2; /* min. abs. val for 1/x that allows normalization */
600
realT MINdenom_2; /* use divzero if denominator < MINdenom_2 */
601
realT MINlastcoord; /* min. last coordinate for qh_scalelast */
602
boolT NARROWhull; /* set in qh_initialhull if angle < qh_MAXnarrow */
603
realT *NEARzero; /* hull_dim array for near zero in gausselim */
604
realT NEARinside; /* keep points for qh_check_maxout if close to facet */
605
realT ONEmerge; /* max distance for merging simplicial facets */
606
realT outside_err; /* application's epsilon for coplanar points
607
qh_check_bestdist() qh_check_points() reports error if point outside */
608
realT WIDEfacet; /* size of wide facet for skipping ridge in
609
area computation and locking centrum */
611
/*-<a href="qh-globa.htm#TOC"
612
>--------------------------------</a><a name="qh-codetern">-</a>
614
qh internal constants
615
internal constants for Qhull
617
char qhull[sizeof("qhull")]; /* "qhull" for checking ownership while debugging */
618
jmp_buf errexit; /* exit label for qh_errexit, defined by setjmp() */
619
char jmpXtra[40]; /* extra bytes in case jmp_buf is defined wrong by compiler */
620
jmp_buf restartexit; /* restart label for qh_errexit, defined by setjmp() */
621
char jmpXtra2[40]; /* extra bytes in case jmp_buf is defined wrong by compiler*/
622
FILE *fin; /* pointer to input file, init by qh_meminit */
623
FILE *fout; /* pointer to output file */
624
FILE *ferr; /* pointer to error file */
625
pointT *interior_point; /* center point of the initial simplex*/
626
int normal_size; /* size in bytes for facet normals and point coords*/
627
int center_size; /* size in bytes for Voronoi centers */
628
int TEMPsize; /* size for small, temporary sets (in quick mem) */
630
/*-<a href="qh-globa.htm#TOC"
631
>--------------------------------</a><a name="qh-lists">-</a>
633
qh facet and vertex lists
634
defines lists of facets, new facets, visible facets, vertices, and
635
new vertices. Includes counts, next ids, and trace ids.
639
facetT *facet_list; /* first facet */
640
facetT *facet_tail; /* end of facet_list (dummy facet) */
641
facetT *facet_next; /* next facet for buildhull()
642
previous facets do not have outside sets
643
NARROWhull: previous facets may have coplanar outside sets for qh_outcoplanar */
644
facetT *newfacet_list; /* list of new facets to end of facet_list */
645
facetT *visible_list; /* list of visible facets preceding newfacet_list,
646
facet->visible set */
647
int num_visible; /* current number of visible facets */
648
unsigned tracefacet_id; /* set at init, then can print whenever */
649
facetT *tracefacet; /* set in newfacet/mergefacet, undone in delfacet*/
650
unsigned tracevertex_id; /* set at buildtracing, can print whenever */
651
vertexT *tracevertex; /* set in newvertex, undone in delvertex*/
652
vertexT *vertex_list; /* list of all vertices, to vertex_tail */
653
vertexT *vertex_tail; /* end of vertex_list (dummy vertex) */
654
vertexT *newvertex_list; /* list of vertices in newfacet_list, to vertex_tail
655
all vertices have 'newlist' set */
656
int num_facets; /* number of facets in facet_list
657
includes visble faces (num_visible) */
658
int num_vertices; /* number of vertices in facet_list */
659
int num_outside; /* number of points in outsidesets (for tracing and RANDOMoutside)
660
includes coplanar outsideset points for NARROWhull/qh_outcoplanar() */
661
int num_good; /* number of good facets (after findgood_all) */
662
unsigned facet_id; /* ID of next, new facet from newfacet() */
663
unsigned ridge_id:24; /* ID of next, new ridge from newridge() */
664
unsigned vertex_id:24; /* ID of next, new vertex from newvertex() */
666
/*-<a href="qh-globa.htm#TOC"
667
>--------------------------------</a><a name="qh-var">-</a>
670
defines minimum and maximum distances, next visit ids, several flags,
671
and other global variables.
672
initialize in qh_initbuild or qh_maxmin if used in qh_buildhull
674
unsigned long hulltime; /* ignore time to set up input and randomize */
675
/* use unsigned to avoid wrap-around errors */
676
boolT ALLOWrestart; /* true if qh_precision can use qh.restartexit */
677
int build_cnt; /* number of calls to qh_initbuild */
678
qh_CENTER CENTERtype; /* current type of facet->center, qh_CENTER */
679
int furthest_id; /* pointid of furthest point, for tracing */
680
facetT *GOODclosest; /* closest facet to GOODthreshold in qh_findgood */
681
boolT hasAreaVolume; /* true if totarea, totvol was defined by qh_getarea */
682
boolT hasTriangulation; /* true if triangulation created by qh_triangulate */
683
realT JOGGLEmax; /* set 'QJn' if randomly joggle input */
684
boolT maxoutdone; /* set qh_check_maxout(), cleared by qh_addpoint() */
685
realT max_outside; /* maximum distance from a point to a facet,
686
before roundoff, not simplicial vertices
687
actual outer plane is +DISTround and
688
computed outer plane is +2*DISTround */
689
realT max_vertex; /* maximum distance (>0) from vertex to a facet,
690
before roundoff, due to a merge */
691
realT min_vertex; /* minimum distance (<0) from vertex to a facet,
692
before roundoff, due to a merge
693
if qh.JOGGLEmax, qh_makenewplanes sets it
694
recomputed if qh.DOcheckmax, default -qh.DISTround */
695
boolT NEWfacets; /* true while visible facets invalid due to new or merge
696
from makecone/attachnewfacets to deletevisible */
697
boolT findbestnew; /* true if partitioning calls qh_findbestnew */
698
boolT findbest_notsharp; /* true if new facets are at least 90 degrees */
699
boolT NOerrexit; /* true if qh.errexit is not available */
700
realT PRINTcradius; /* radius for printing centrums */
701
realT PRINTradius; /* radius for printing vertex spheres and points */
702
boolT POSTmerging; /* true when post merging */
703
int printoutvar; /* temporary variable for qh_printbegin, etc. */
704
int printoutnum; /* number of facets printed */
705
boolT QHULLfinished; /* True after qhull() is finished */
706
realT totarea; /* 'FA': total facet area computed by qh_getarea, hasAreaVolume */
707
realT totvol; /* 'FA': total volume computed by qh_getarea, hasAreaVolume */
708
unsigned int visit_id; /* unique ID for searching neighborhoods, */
709
unsigned int vertex_visit:31; /* unique ID for searching vertices, reset with qh_buildtracing */
710
boolT ZEROall_ok; /* True if qh_checkzero always succeeds */
711
boolT WAScoplanar; /* True if qh_partitioncoplanar (qh_check_maxout) */
713
/*-<a href="qh-globa.htm#TOC"
714
>--------------------------------</a><a name="qh-set">-</a>
717
defines sets for merging, initial simplex, hashing, extra input points,
720
setT *facet_mergeset; /* temporary set of merges to be done */
721
setT *degen_mergeset; /* temporary set of degenerate and redundant merges */
722
setT *hash_table; /* hash table for matching ridges in qh_matchfacets
724
setT *other_points; /* additional points */
725
setT *del_vertices; /* vertices to partition and delete with visible
726
facets. Have deleted set for checkfacet */
728
/*-<a href="qh-globa.htm#TOC"
729
>--------------------------------</a><a name="qh-buf">-</a>
732
defines buffers for maxtrix operations, input, and error messages
734
coordT *gm_matrix; /* (dim+1)Xdim matrix for geom.c */
735
coordT **gm_row; /* array of gm_matrix rows */
736
char* line; /* malloc'd input line of maxline+1 chars */
738
coordT *half_space; /* malloc'd input array for halfspace (qh normal_size+coordT) */
739
coordT *temp_malloc; /* malloc'd input array for points */
741
/*-<a href="qh-globa.htm#TOC"
742
>--------------------------------</a><a name="qh-static">-</a>
745
defines static variables for individual functions
748
do not use 'static' within a function. Multiple instances of qhull
751
do not assume zero initialization, 'QPn' may cause a restart
753
boolT ERREXITcalled; /* true during qh_errexit (prevents duplicate calls */
754
boolT firstcentrum; /* for qh_printcentrum */
755
boolT old_randomdist; /* save RANDOMdist flag during io, tracing, or statistics */
756
setT *coplanarfacetset; /* set of coplanar facets for searching qh_findbesthorizon() */
757
realT last_low; /* qh_scalelast parameters for qh_setdelaunay */
760
unsigned lastreport; /* for qh_buildtracing */
761
int mergereport; /* for qh_tracemerging */
762
qhstatT *old_qhstat; /* for saving qh_qhstat in save_qhull() and UsingLibQhull. Free with qh_free() */
763
setT *old_tempstack; /* for saving qhmem.tempstack in save_qhull */
764
int ridgeoutnum; /* number of ridges for 4OFF output (qh_printbegin,etc) */
767
/*=========== -macros- =========================*/
769
/*-<a href="qh-poly.htm#TOC"
770
>--------------------------------</a><a name="otherfacet_">-</a>
772
otherfacet_(ridge, facet)
773
return neighboring facet for a ridge in facet
775
#define otherfacet_(ridge, facet) \
776
(((ridge)->top == (facet)) ? (ridge)->bottom : (ridge)->top)
778
/*-<a href="qh-poly.htm#TOC"
779
>--------------------------------</a><a name="getid_">-</a>
782
return int ID for facet, ridge, or vertex
785
#define getid_(p) ((p) ? (int)((p)->id) : -1)
787
/*============== FORALL macros ===================*/
789
/*-<a href="qh-poly.htm#TOC"
790
>--------------------------------</a><a name="FORALLfacets">-</a>
793
assign 'facet' to each facet in qh.facet_list
796
uses 'facetT *facet;'
797
assumes last facet is a sentinel
800
FORALLfacet_( facetlist )
802
#define FORALLfacets for (facet=qh facet_list;facet && facet->next;facet=facet->next)
804
/*-<a href="qh-poly.htm#TOC"
805
>--------------------------------</a><a name="FORALLpoints">-</a>
808
assign 'point' to each point in qh.first_point, qh.num_points
811
coordT *point, *pointtemp;
813
#define FORALLpoints FORALLpoint_(qh first_point, qh num_points)
815
/*-<a href="qh-poly.htm#TOC"
816
>--------------------------------</a><a name="FORALLpoint_">-</a>
818
FORALLpoint_( points, num) { ... }
819
assign 'point' to each point in points array of num points
822
coordT *point, *pointtemp;
824
#define FORALLpoint_(points, num) for (point= (points), \
825
pointtemp= (points)+qh hull_dim*(num); point < pointtemp; point += qh hull_dim)
827
/*-<a href="qh-poly.htm#TOC"
828
>--------------------------------</a><a name="FORALLvertices">-</a>
830
FORALLvertices { ... }
831
assign 'vertex' to each vertex in qh.vertex_list
837
assumes qh.vertex_list terminated with a sentinel
839
#define FORALLvertices for (vertex=qh vertex_list;vertex && vertex->next;vertex= vertex->next)
841
/*-<a href="qh-poly.htm#TOC"
842
>--------------------------------</a><a name="FOREACHfacet_">-</a>
844
FOREACHfacet_( facets ) { ... }
845
assign 'facet' to each facet in facets
848
facetT *facet, **facetp;
851
<a href="qset.h#FOREACHsetelement_">FOREACHsetelement_</a>
853
#define FOREACHfacet_(facets) FOREACHsetelement_(facetT, facets, facet)
855
/*-<a href="qh-poly.htm#TOC"
856
>--------------------------------</a><a name="FOREACHneighbor_">-</a>
858
FOREACHneighbor_( facet ) { ... }
859
assign 'neighbor' to each neighbor in facet->neighbors
861
FOREACHneighbor_( vertex ) { ... }
862
assign 'neighbor' to each neighbor in vertex->neighbors
865
facetT *neighbor, **neighborp;
868
<a href="qset.h#FOREACHsetelement_">FOREACHsetelement_</a>
870
#define FOREACHneighbor_(facet) FOREACHsetelement_(facetT, facet->neighbors, neighbor)
872
/*-<a href="qh-poly.htm#TOC"
873
>--------------------------------</a><a name="FOREACHpoint_">-</a>
875
FOREACHpoint_( points ) { ... }
876
assign 'point' to each point in points set
879
pointT *point, **pointp;
882
<a href="qset.h#FOREACHsetelement_">FOREACHsetelement_</a>
884
#define FOREACHpoint_(points) FOREACHsetelement_(pointT, points, point)
886
/*-<a href="qh-poly.htm#TOC"
887
>--------------------------------</a><a name="FOREACHridge_">-</a>
889
FOREACHridge_( ridges ) { ... }
890
assign 'ridge' to each ridge in ridges set
893
ridgeT *ridge, **ridgep;
896
<a href="qset.h#FOREACHsetelement_">FOREACHsetelement_</a>
898
#define FOREACHridge_(ridges) FOREACHsetelement_(ridgeT, ridges, ridge)
900
/*-<a href="qh-poly.htm#TOC"
901
>--------------------------------</a><a name="FOREACHvertex_">-</a>
903
FOREACHvertex_( vertices ) { ... }
904
assign 'vertex' to each vertex in vertices set
907
vertexT *vertex, **vertexp;
910
<a href="qset.h#FOREACHsetelement_">FOREACHsetelement_</a>
912
#define FOREACHvertex_(vertices) FOREACHsetelement_(vertexT, vertices,vertex)
914
/*-<a href="qh-poly.htm#TOC"
915
>--------------------------------</a><a name="FOREACHfacet_i_">-</a>
917
FOREACHfacet_i_( facets ) { ... }
918
assign 'facet' and 'facet_i' for each facet in facets set
922
int facet_n, facet_i;
925
<a href="qset.h#FOREACHsetelement_i_">FOREACHsetelement_i_</a>
927
#define FOREACHfacet_i_(facets) FOREACHsetelement_i_(facetT, facets, facet)
929
/*-<a href="qh-poly.htm#TOC"
930
>--------------------------------</a><a name="FOREACHneighbor_i_">-</a>
932
FOREACHneighbor_i_( facet ) { ... }
933
assign 'neighbor' and 'neighbor_i' for each neighbor in facet->neighbors
935
FOREACHneighbor_i_( vertex ) { ... }
936
assign 'neighbor' and 'neighbor_i' for each neighbor in vertex->neighbors
940
int neighbor_n, neighbor_i;
943
<a href="qset.h#FOREACHsetelement_i_">FOREACHsetelement_i_</a>
945
#define FOREACHneighbor_i_(facet) FOREACHsetelement_i_(facetT, facet->neighbors, neighbor)
947
/*-<a href="qh-poly.htm#TOC"
948
>--------------------------------</a><a name="FOREACHpoint_i_">-</a>
950
FOREACHpoint_i_( points ) { ... }
951
assign 'point' and 'point_i' for each point in points set
955
int point_n, point_i;
958
<a href="qset.h#FOREACHsetelement_i_">FOREACHsetelement_i_</a>
960
#define FOREACHpoint_i_(points) FOREACHsetelement_i_(pointT, points, point)
962
/*-<a href="qh-poly.htm#TOC"
963
>--------------------------------</a><a name="FOREACHridge_i_">-</a>
965
FOREACHridge_i_( ridges ) { ... }
966
assign 'ridge' and 'ridge_i' for each ridge in ridges set
970
int ridge_n, ridge_i;
973
<a href="qset.h#FOREACHsetelement_i_">FOREACHsetelement_i_</a>
975
#define FOREACHridge_i_(ridges) FOREACHsetelement_i_(ridgeT, ridges, ridge)
977
/*-<a href="qh-poly.htm#TOC"
978
>--------------------------------</a><a name="FOREACHvertex_i_">-</a>
980
FOREACHvertex_i_( vertices ) { ... }
981
assign 'vertex' and 'vertex_i' for each vertex in vertices set
985
int vertex_n, vertex_i;
988
<a href="qset.h#FOREACHsetelement_i_">FOREACHsetelement_i_</a>
990
#define FOREACHvertex_i_(vertices) FOREACHsetelement_i_(vertexT, vertices,vertex)
992
/********* -libqhull.c prototypes (duplicated from qhull_a.h) **********************/
995
boolT qh_addpoint(pointT *furthest, facetT *facet, boolT checkdist);
996
void qh_printsummary(FILE *fp);
998
/********* -user.c prototypes (alphabetical) **********************/
1000
void qh_errexit(int exitcode, facetT *facet, ridgeT *ridge);
1001
void qh_errprint(const char* string, facetT *atfacet, facetT *otherfacet, ridgeT *atridge, vertexT *atvertex);
1002
int qh_new_qhull(int dim, int numpoints, coordT *points, boolT ismalloc,
1003
char *qhull_cmd, FILE *outfile, FILE *errfile);
1004
void qh_printfacetlist(facetT *facetlist, setT *facets, boolT printall);
1005
void qh_printhelp_degenerate(FILE *fp);
1006
void qh_printhelp_narrowhull(FILE *fp, realT minangle);
1007
void qh_printhelp_singular(FILE *fp);
1008
void qh_user_memsizes(void);
1010
/********* -usermem.c prototypes (alphabetical) **********************/
1011
void qh_exit(int exitcode);
1012
void qh_free(void *mem);
1013
void *qh_malloc(size_t size);
1015
/********* -userprintf.c and userprintf_rbox.c prototypes **********************/
1016
void qh_fprintf(FILE *fp, int msgcode, const char *fmt, ... );
1017
void qh_fprintf_rbox(FILE *fp, int msgcode, const char *fmt, ... );
1019
/***** -geom.c/geom2.c/random.c prototypes (duplicated from geom.h, random.h) ****************/
1021
facetT *qh_findbest(pointT *point, facetT *startfacet,
1022
boolT bestoutside, boolT newfacets, boolT noupper,
1023
realT *dist, boolT *isoutside, int *numpart);
1024
facetT *qh_findbestnew(pointT *point, facetT *startfacet,
1025
realT *dist, boolT bestoutside, boolT *isoutside, int *numpart);
1026
boolT qh_gram_schmidt(int dim, realT **rows);
1027
void qh_outerinner(facetT *facet, realT *outerplane, realT *innerplane);
1028
void qh_printsummary(FILE *fp);
1029
void qh_projectinput(void);
1030
void qh_randommatrix(realT *buffer, int dim, realT **row);
1031
void qh_rotateinput(realT **rows);
1032
void qh_scaleinput(void);
1033
void qh_setdelaunay(int dim, int count, pointT *points);
1034
coordT *qh_sethalfspace_all(int dim, int count, coordT *halfspaces, pointT *feasible);
1036
/***** -global.c prototypes (alphabetical) ***********************/
1038
unsigned long qh_clock(void);
1039
void qh_checkflags(char *command, char *hiddenflags);
1040
void qh_clear_outputflags(void);
1041
void qh_freebuffers(void);
1042
void qh_freeqhull(boolT allmem);
1043
void qh_freeqhull2(boolT allmem);
1044
void qh_init_A(FILE *infile, FILE *outfile, FILE *errfile, int argc, char *argv[]);
1045
void qh_init_B(coordT *points, int numpoints, int dim, boolT ismalloc);
1046
void qh_init_qhull_command(int argc, char *argv[]);
1047
void qh_initbuffers(coordT *points, int numpoints, int dim, boolT ismalloc);
1048
void qh_initflags(char *command);
1049
void qh_initqhull_buffers(void);
1050
void qh_initqhull_globals(coordT *points, int numpoints, int dim, boolT ismalloc);
1051
void qh_initqhull_mem(void);
1052
void qh_initqhull_outputflags(void);
1053
void qh_initqhull_start(FILE *infile, FILE *outfile, FILE *errfile);
1054
void qh_initqhull_start2(FILE *infile, FILE *outfile, FILE *errfile);
1055
void qh_initthresholds(char *command);
1056
void qh_option(const char *option, int *i, realT *r);
1058
void qh_restore_qhull(qhT **oldqh);
1059
qhT *qh_save_qhull(void);
1062
/***** -io.c prototypes (duplicated from io.h) ***********************/
1064
void dfacet( unsigned id);
1065
void dvertex( unsigned id);
1066
void qh_printneighborhood(FILE *fp, qh_PRINT format, facetT *facetA, facetT *facetB, boolT printall);
1067
void qh_produce_output(void);
1068
coordT *qh_readpoints(int *numpoints, int *dimension, boolT *ismalloc);
1071
/********* -mem.c prototypes (duplicated from mem.h) **********************/
1073
void qh_meminit(FILE *ferr);
1074
void qh_memfreeshort(int *curlong, int *totlong);
1076
/********* -poly.c/poly2.c prototypes (duplicated from poly.h) **********************/
1078
void qh_check_output(void);
1079
void qh_check_points(void);
1080
setT *qh_facetvertices(facetT *facetlist, setT *facets, boolT allfacets);
1081
facetT *qh_findbestfacet(pointT *point, boolT bestoutside,
1082
realT *bestdist, boolT *isoutside);
1083
vertexT *qh_nearvertex(facetT *facet, pointT *point, realT *bestdistp);
1084
pointT *qh_point(int id);
1085
setT *qh_pointfacet(void /*qh.facet_list*/);
1086
int qh_pointid(pointT *point);
1087
setT *qh_pointvertex(void /*qh.facet_list*/);
1088
void qh_setvoronoi_all(void);
1089
void qh_triangulate(void /*qh facet_list*/);
1091
/********* -rboxpoints.c prototypes **********************/
1092
int qh_rboxpoints(FILE* fout, FILE* ferr, char* rbox_command);
1093
void qh_errexit_rbox(int exitcode);
1095
/********* -stat.c prototypes (duplicated from stat.h) **********************/
1097
void qh_collectstatistics(void);
1098
void qh_printallstatistics(FILE *fp, const char *string);
1100
#endif /* qhDEFlibqhull */