1
2 /*<html><pre> -<a href="qh-globa_r.htm"
3 >-------------------------------</a><a name="TOP">-</a>
4
5 global_r.c
6 initializes all the globals of the qhull application
7
8 see README
9
10 see libqhull_r.h for qh.globals and function prototypes
11
12 see qhull_ra.h for internal functions
13
14 Copyright (c) 1993-2020 The Geometry Center.
15 $Id: //main/2019/qhull/src/libqhull_r/global_r.c#19 $$Change: 3037 $
16 $DateTime: 2020/09/03 17:28:32 $$Author: bbarber $
17 */
18
19 #include "qhull_ra.h"
20
21 /*========= qh->definition -- globals defined in libqhull_r.h =======================*/
22
23 /*-<a href ="qh-globa_r.htm#TOC"
24 >--------------------------------</a><a name="version">-</a>
25
26 qh_version
27 version string by year and date
28 qh_version2 for Unix users and -V
29
30 the revision increases on code changes only
31
32 notes:
33 change date: Changes.txt, Announce.txt, index.htm, README.txt,
34 qhull-news.html, Eudora signatures, CMakeLists.txt
35 change version: README.txt, qh-get.htm, File_id.diz, Makefile.txt, CMakeLists.txt
36 check that CMakeLists.txt @version is the same as qh_version2
37 change year: Copying.txt
38 check download size
39 recompile user_eg_r.c, rbox_r.c, libqhull_r.c, qconvex_r.c, qdelaun_r.c qvoronoi_r.c, qhalf_r.c, testqset_r.c
40 */
41
42 const char qh_version[]= "2020.2.r 2020/08/31";
43 const char qh_version2[]= "qhull_r 8.0.2 (2020.2.r 2020/08/31)";
44
45 /*-<a href="qh-globa_r.htm#TOC"
46 >-------------------------------</a><a name="appendprint">-</a>
47
48 qh_appendprint(qh, printFormat )
49 append printFormat to qh.PRINTout unless already defined
50 */
qh_appendprint(qhT * qh,qh_PRINT format)51 void qh_appendprint(qhT *qh, qh_PRINT format) {
52 int i;
53
54 for (i=0; i < qh_PRINTEND; i++) {
55 if (qh->PRINTout[i] == format && format != qh_PRINTqhull)
56 break;
57 if (!qh->PRINTout[i]) {
58 qh->PRINTout[i]= format;
59 break;
60 }
61 }
62 } /* appendprint */
63
64 /*-<a href="qh-globa_r.htm#TOC"
65 >-------------------------------</a><a name="checkflags">-</a>
66
67 qh_checkflags(qh, commandStr, hiddenFlags )
68 errors if commandStr contains hiddenFlags
69 hiddenFlags starts and ends with a space and is space delimited (checked)
70
71 notes:
72 ignores first word (e.g., "qconvex i")
73 use qh_strtol/strtod since strtol/strtod may or may not skip trailing spaces
74
75 see:
76 qh_initflags() initializes Qhull according to commandStr
77 */
qh_checkflags(qhT * qh,char * command,char * hiddenflags)78 void qh_checkflags(qhT *qh, char *command, char *hiddenflags) {
79 char *s= command, *t, *chkerr; /* qh_skipfilename is non-const */
80 char key, opt, prevopt;
81 char chkkey[]= " "; /* check one character options ('s') */
82 char chkopt[]= " "; /* check two character options ('Ta') */
83 char chkopt2[]= " "; /* check three character options ('Q12') */
84 boolT waserr= False;
85
86 if (*hiddenflags != ' ' || hiddenflags[strlen(hiddenflags)-1] != ' ') {
87 qh_fprintf(qh, qh->ferr, 6026, "qhull internal error (qh_checkflags): hiddenflags must start and end with a space: \"%s\"\n", hiddenflags);
88 qh_errexit(qh, qh_ERRqhull, NULL, NULL);
89 }
90 if (strpbrk(hiddenflags, ",\n\r\t")) {
91 qh_fprintf(qh, qh->ferr, 6027, "qhull internal error (qh_checkflags): hiddenflags contains commas, newlines, or tabs: \"%s\"\n", hiddenflags);
92 qh_errexit(qh, qh_ERRqhull, NULL, NULL);
93 }
94 while (*s && !isspace(*s)) /* skip program name */
95 s++;
96 while (*s) {
97 while (*s && isspace(*s))
98 s++;
99 if (*s == '-')
100 s++;
101 if (!*s)
102 break;
103 key= *s++;
104 chkerr= NULL;
105 if (key == 'T' && (*s == 'I' || *s == 'O')) { /* TI or TO 'file name' */
106 s= qh_skipfilename(qh, ++s);
107 continue;
108 }
109 chkkey[1]= key;
110 if (strstr(hiddenflags, chkkey)) {
111 chkerr= chkkey;
112 }else if (isupper(key)) {
113 opt= ' ';
114 prevopt= ' ';
115 chkopt[1]= key;
116 chkopt2[1]= key;
117 while (!chkerr && *s && !isspace(*s)) {
118 opt= *s++;
119 if (isalpha(opt)) {
120 chkopt[2]= opt;
121 if (strstr(hiddenflags, chkopt))
122 chkerr= chkopt;
123 if (prevopt != ' ') {
124 chkopt2[2]= prevopt;
125 chkopt2[3]= opt;
126 if (strstr(hiddenflags, chkopt2))
127 chkerr= chkopt2;
128 }
129 }else if (key == 'Q' && isdigit(opt) && prevopt != 'b'
130 && (prevopt == ' ' || islower(prevopt))) {
131 if (isdigit(*s)) { /* Q12 */
132 chkopt2[2]= opt;
133 chkopt2[3]= *s++;
134 if (strstr(hiddenflags, chkopt2))
135 chkerr= chkopt2;
136 }else {
137 chkopt[2]= opt;
138 if (strstr(hiddenflags, chkopt))
139 chkerr= chkopt;
140 }
141 }else {
142 qh_strtod(s-1, &t);
143 if (s < t)
144 s= t;
145 }
146 prevopt= opt;
147 }
148 }
149 if (chkerr) {
150 *chkerr= '\'';
151 chkerr[strlen(chkerr)-1]= '\'';
152 qh_fprintf(qh, qh->ferr, 6029, "qhull option error: option %s is not used with this program.\n It may be used with qhull.\n", chkerr);
153 waserr= True;
154 }
155 }
156 if (waserr)
157 qh_errexit(qh, qh_ERRinput, NULL, NULL);
158 } /* checkflags */
159
160 /*-<a href="qh-globa_r.htm#TOC"
161 >-------------------------------</a><a name="clear_outputflags">-</a>
162
163 qh_clear_outputflags(qh)
164 Clear output flags for QhullPoints
165 */
qh_clear_outputflags(qhT * qh)166 void qh_clear_outputflags(qhT *qh) {
167 int i,k;
168
169 qh->ANNOTATEoutput= False;
170 qh->DOintersections= False;
171 qh->DROPdim= -1;
172 qh->FORCEoutput= False;
173 qh->GETarea= False;
174 qh->GOODpoint= 0;
175 qh->GOODpointp= NULL;
176 qh->GOODthreshold= False;
177 qh->GOODvertex= 0;
178 qh->GOODvertexp= NULL;
179 qh->IStracing= 0;
180 qh->KEEParea= False;
181 qh->KEEPmerge= False;
182 qh->KEEPminArea= REALmax;
183 qh->PRINTcentrums= False;
184 qh->PRINTcoplanar= False;
185 qh->PRINTdots= False;
186 qh->PRINTgood= False;
187 qh->PRINTinner= False;
188 qh->PRINTneighbors= False;
189 qh->PRINTnoplanes= False;
190 qh->PRINToptions1st= False;
191 qh->PRINTouter= False;
192 qh->PRINTprecision= True;
193 qh->PRINTridges= False;
194 qh->PRINTspheres= False;
195 qh->PRINTstatistics= False;
196 qh->PRINTsummary= False;
197 qh->PRINTtransparent= False;
198 qh->SPLITthresholds= False;
199 qh->TRACElevel= 0;
200 qh->TRInormals= False;
201 qh->USEstdout= False;
202 qh->VERIFYoutput= False;
203 for (k=qh->input_dim+1; k--; ) { /* duplicated in qh_initqhull_buffers and qh_clear_outputflags */
204 qh->lower_threshold[k]= -REALmax;
205 qh->upper_threshold[k]= REALmax;
206 qh->lower_bound[k]= -REALmax;
207 qh->upper_bound[k]= REALmax;
208 }
209
210 for (i=0; i < qh_PRINTEND; i++) {
211 qh->PRINTout[i]= qh_PRINTnone;
212 }
213
214 if (!qh->qhull_commandsiz2)
215 qh->qhull_commandsiz2= (int)strlen(qh->qhull_command); /* WARN64 */
216 else {
217 qh->qhull_command[qh->qhull_commandsiz2]= '\0';
218 }
219 if (!qh->qhull_optionsiz2)
220 qh->qhull_optionsiz2= (int)strlen(qh->qhull_options); /* WARN64 */
221 else {
222 qh->qhull_options[qh->qhull_optionsiz2]= '\0';
223 qh->qhull_optionlen= qh_OPTIONline; /* start a new line */
224 }
225 } /* clear_outputflags */
226
227 /*-<a href="qh-globa_r.htm#TOC"
228 >-------------------------------</a><a name="clock">-</a>
229
230 qh_clock()
231 return user CPU time in 100ths (qh_SECtick)
232 only defined for qh_CLOCKtype == 2
233
234 notes:
235 use first value to determine time 0
236 from Stevens '92 8.15
237 */
qh_clock(qhT * qh)238 unsigned long qh_clock(qhT *qh) {
239
240 #if (qh_CLOCKtype == 2)
241 struct tms time;
242 static long clktck; /* initialized first call and never updated */
243 double ratio, cpu;
244 unsigned long ticks;
245
246 if (!clktck) {
247 if ((clktck= sysconf(_SC_CLK_TCK)) < 0) {
248 qh_fprintf(qh, qh->ferr, 6030, "qhull internal error (qh_clock): sysconf() failed. Use qh_CLOCKtype 1 in user_r.h\n");
249 qh_errexit(qh, qh_ERRqhull, NULL, NULL);
250 }
251 }
252 if (times(&time) == -1) {
253 qh_fprintf(qh, qh->ferr, 6031, "qhull internal error (qh_clock): times() failed. Use qh_CLOCKtype 1 in user_r.h\n");
254 qh_errexit(qh, qh_ERRqhull, NULL, NULL);
255 }
256 ratio= qh_SECticks / (double)clktck;
257 ticks= time.tms_utime * ratio;
258 return ticks;
259 #else
260 qh_fprintf(qh, qh->ferr, 6032, "qhull internal error (qh_clock): use qh_CLOCKtype 2 in user_r.h\n");
261 qh_errexit(qh, qh_ERRqhull, NULL, NULL); /* never returns */
262 return 0;
263 #endif
264 } /* clock */
265
266 /*-<a href="qh-globa_r.htm#TOC"
267 >-------------------------------</a><a name="freebuffers">-</a>
268
269 qh_freebuffers()
270 free up global memory buffers
271
272 notes:
273 must match qh_initbuffers()
274 */
qh_freebuffers(qhT * qh)275 void qh_freebuffers(qhT *qh) {
276
277 trace5((qh, qh->ferr, 5001, "qh_freebuffers: freeing up global memory buffers\n"));
278 /* allocated by qh_initqhull_buffers */
279 qh_setfree(qh, &qh->other_points);
280 qh_setfree(qh, &qh->del_vertices);
281 qh_setfree(qh, &qh->coplanarfacetset);
282 qh_memfree(qh, qh->NEARzero, qh->hull_dim * (int)sizeof(realT));
283 qh_memfree(qh, qh->lower_threshold, (qh->input_dim+1) * (int)sizeof(realT));
284 qh_memfree(qh, qh->upper_threshold, (qh->input_dim+1) * (int)sizeof(realT));
285 qh_memfree(qh, qh->lower_bound, (qh->input_dim+1) * (int)sizeof(realT));
286 qh_memfree(qh, qh->upper_bound, (qh->input_dim+1) * (int)sizeof(realT));
287 qh_memfree(qh, qh->gm_matrix, (qh->hull_dim+1) * qh->hull_dim * (int)sizeof(coordT));
288 qh_memfree(qh, qh->gm_row, (qh->hull_dim+1) * (int)sizeof(coordT *));
289 qh->NEARzero= qh->lower_threshold= qh->upper_threshold= NULL;
290 qh->lower_bound= qh->upper_bound= NULL;
291 qh->gm_matrix= NULL;
292 qh->gm_row= NULL;
293
294 if (qh->line) /* allocated by qh_readinput, freed if no error */
295 qh_free(qh->line);
296 if (qh->half_space)
297 qh_free(qh->half_space);
298 if (qh->temp_malloc)
299 qh_free(qh->temp_malloc);
300 if (qh->feasible_point) /* allocated by qh_readfeasible */
301 qh_free(qh->feasible_point);
302 if (qh->feasible_string) /* allocated by qh_initflags */
303 qh_free(qh->feasible_string);
304 qh->line= qh->feasible_string= NULL;
305 qh->half_space= qh->feasible_point= qh->temp_malloc= NULL;
306 /* usually allocated by qh_readinput */
307 if (qh->first_point && qh->POINTSmalloc) {
308 qh_free(qh->first_point);
309 qh->first_point= NULL;
310 }
311 if (qh->input_points && qh->input_malloc) { /* set by qh_joggleinput */
312 qh_free(qh->input_points);
313 qh->input_points= NULL;
314 }
315 trace5((qh, qh->ferr, 5002, "qh_freebuffers: finished\n"));
316 } /* freebuffers */
317
318
319 /*-<a href="qh-globa_r.htm#TOC"
320 >-------------------------------</a><a name="freebuild">-</a>
321
322 qh_freebuild(qh, allmem )
323 free global memory used by qh_initbuild and qh_buildhull
324 if !allmem,
325 does not free short memory (e.g., facetT, freed by qh_memfreeshort)
326
327 design:
328 free centrums
329 free each vertex
330 for each facet
331 free ridges
332 free outside set, coplanar set, neighbor set, ridge set, vertex set
333 free facet
334 free hash table
335 free interior point
336 free merge sets
337 free temporary sets
338 */
qh_freebuild(qhT * qh,boolT allmem)339 void qh_freebuild(qhT *qh, boolT allmem) {
340 facetT *facet, *previousfacet= NULL;
341 vertexT *vertex, *previousvertex= NULL;
342 ridgeT *ridge, **ridgep, *previousridge= NULL;
343 mergeT *merge, **mergep;
344 int newsize;
345 boolT freeall;
346
347 /* free qhT global sets first, includes references from qh_buildhull */
348 trace5((qh, qh->ferr, 5004, "qh_freebuild: free global sets\n"));
349 FOREACHmerge_(qh->facet_mergeset) /* usually empty */
350 qh_memfree(qh, merge, (int)sizeof(mergeT));
351 FOREACHmerge_(qh->degen_mergeset) /* usually empty */
352 qh_memfree(qh, merge, (int)sizeof(mergeT));
353 FOREACHmerge_(qh->vertex_mergeset) /* usually empty */
354 qh_memfree(qh, merge, (int)sizeof(mergeT));
355 qh->facet_mergeset= NULL; /* temp set freed by qh_settempfree_all */
356 qh->degen_mergeset= NULL; /* temp set freed by qh_settempfree_all */
357 qh->vertex_mergeset= NULL; /* temp set freed by qh_settempfree_all */
358 qh_setfree(qh, &(qh->hash_table));
359 trace5((qh, qh->ferr, 5003, "qh_freebuild: free temporary sets (qh_settempfree_all)\n"));
360 qh_settempfree_all(qh);
361 trace1((qh, qh->ferr, 1005, "qh_freebuild: free memory from qh_inithull and qh_buildhull\n"));
362 if (qh->del_vertices)
363 qh_settruncate(qh, qh->del_vertices, 0);
364 if (allmem) {
365 while ((vertex= qh->vertex_list)) {
366 if (vertex->next)
367 qh_delvertex(qh, vertex);
368 else {
369 qh_memfree(qh, vertex, (int)sizeof(vertexT)); /* sentinel */
370 qh->newvertex_list= qh->vertex_list= NULL;
371 break;
372 }
373 previousvertex= vertex; /* in case of memory fault */
374 QHULL_UNUSED(previousvertex)
375 }
376 }else if (qh->VERTEXneighbors) {
377 FORALLvertices
378 qh_setfreelong(qh, &(vertex->neighbors));
379 }
380 qh->VERTEXneighbors= False;
381 qh->GOODclosest= NULL;
382 if (allmem) {
383 FORALLfacets {
384 FOREACHridge_(facet->ridges)
385 ridge->seen= False;
386 }
387 while ((facet= qh->facet_list)) {
388 if (!facet->newfacet || !qh->NEWtentative || qh_setsize(qh, facet->ridges) > 1) { /* skip tentative horizon ridges */
389 trace4((qh, qh->ferr, 4095, "qh_freebuild: delete the previously-seen ridges of f%d\n", facet->id));
390 FOREACHridge_(facet->ridges) {
391 if (ridge->seen)
392 qh_delridge(qh, ridge);
393 else
394 ridge->seen= True;
395 previousridge= ridge; /* in case of memory fault */
396 QHULL_UNUSED(previousridge)
397 }
398 }
399 qh_setfree(qh, &(facet->outsideset));
400 qh_setfree(qh, &(facet->coplanarset));
401 qh_setfree(qh, &(facet->neighbors));
402 qh_setfree(qh, &(facet->ridges));
403 qh_setfree(qh, &(facet->vertices));
404 if (facet->next)
405 qh_delfacet(qh, facet);
406 else {
407 qh_memfree(qh, facet, (int)sizeof(facetT));
408 qh->visible_list= qh->newfacet_list= qh->facet_list= NULL;
409 }
410 previousfacet= facet; /* in case of memory fault */
411 QHULL_UNUSED(previousfacet)
412 }
413 }else {
414 freeall= True;
415 if (qh_setlarger_quick(qh, qh->hull_dim + 1, &newsize))
416 freeall= False;
417 FORALLfacets {
418 qh_setfreelong(qh, &(facet->outsideset));
419 qh_setfreelong(qh, &(facet->coplanarset));
420 if (!facet->simplicial || freeall) {
421 qh_setfreelong(qh, &(facet->neighbors));
422 qh_setfreelong(qh, &(facet->ridges));
423 qh_setfreelong(qh, &(facet->vertices));
424 }
425 }
426 }
427 /* qh internal constants */
428 qh_memfree(qh, qh->interior_point, qh->normal_size);
429 qh->interior_point= NULL;
430 } /* freebuild */
431
432 /*-<a href="qh-globa_r.htm#TOC"
433 >-------------------------------</a><a name="freeqhull">-</a>
434
435 qh_freeqhull(qh, allmem )
436
437 free global memory and set qhT to 0
438 if !allmem,
439 does not free short memory (freed by qh_memfreeshort unless qh_NOmem)
440
441 notes:
442 sets qh.NOerrexit in case caller forgets to
443 Does not throw errors
444
445 see:
446 see qh_initqhull_start2()
447 For libqhull_r, qhstatT is part of qhT
448
449 design:
450 free global and temporary memory from qh_initbuild and qh_buildhull
451 free buffers
452 */
qh_freeqhull(qhT * qh,boolT allmem)453 void qh_freeqhull(qhT *qh, boolT allmem) {
454
455 qh->NOerrexit= True; /* no more setjmp since called at exit and ~QhullQh */
456 trace1((qh, qh->ferr, 1006, "qh_freeqhull: free global memory\n"));
457 qh_freebuild(qh, allmem);
458 qh_freebuffers(qh);
459 trace1((qh, qh->ferr, 1061, "qh_freeqhull: clear qhT except for qh.qhmem and qh.qhstat\n"));
460 /* memset is the same in qh_freeqhull() and qh_initqhull_start2() */
461 memset((char *)qh, 0, sizeof(qhT)-sizeof(qhmemT)-sizeof(qhstatT));
462 qh->NOerrexit= True;
463 } /* freeqhull */
464
465 /*-<a href="qh-globa_r.htm#TOC"
466 >-------------------------------</a><a name="init_A">-</a>
467
468 qh_init_A(qh, infile, outfile, errfile, argc, argv )
469 initialize memory and stdio files
470 convert input options to option string (qh.qhull_command)
471
472 notes:
473 infile may be NULL if qh_readpoints() is not called
474
475 errfile should always be defined. It is used for reporting
476 errors. outfile is used for output and format options.
477
478 argc/argv may be 0/NULL
479
480 called before error handling initialized
481 qh_errexit() may not be used
482 */
qh_init_A(qhT * qh,FILE * infile,FILE * outfile,FILE * errfile,int argc,char * argv[])483 void qh_init_A(qhT *qh, FILE *infile, FILE *outfile, FILE *errfile, int argc, char *argv[]) {
484 qh_meminit(qh, errfile);
485 qh_initqhull_start(qh, infile, outfile, errfile);
486 qh_init_qhull_command(qh, argc, argv);
487 } /* init_A */
488
489 /*-<a href="qh-globa_r.htm#TOC"
490 >-------------------------------</a><a name="init_B">-</a>
491
492 qh_init_B(qh, points, numpoints, dim, ismalloc )
493 initialize globals for points array
494
495 points has numpoints dim-dimensional points
496 points[0] is the first coordinate of the first point
497 points[1] is the second coordinate of the first point
498 points[dim] is the first coordinate of the second point
499
500 ismalloc=True
501 Qhull will call qh_free(points) on exit or input transformation
502 ismalloc=False
503 Qhull will allocate a new point array if needed for input transformation
504
505 qh.qhull_command
506 is the option string.
507 It is defined by qh_init_B(), qh_qhull_command(), or qh_initflags
508
509 returns:
510 if qh.PROJECTinput or (qh.DELAUNAY and qh.PROJECTdelaunay)
511 projects the input to a new point array
512
513 if qh.DELAUNAY,
514 qh.hull_dim is increased by one
515 if qh.ATinfinity,
516 qh_projectinput adds point-at-infinity for Delaunay tri.
517
518 if qh.SCALEinput
519 changes the upper and lower bounds of the input, see qh_scaleinput
520
521 if qh.ROTATEinput
522 rotates the input by a random rotation, see qh_rotateinput
523 if qh.DELAUNAY
524 rotates about the last coordinate
525
526 notes:
527 called after points are defined
528 qh_errexit() may be used
529 */
qh_init_B(qhT * qh,coordT * points,int numpoints,int dim,boolT ismalloc)530 void qh_init_B(qhT *qh, coordT *points, int numpoints, int dim, boolT ismalloc) {
531 qh_initqhull_globals(qh, points, numpoints, dim, ismalloc);
532 if (qh->qhmem.LASTsize == 0)
533 qh_initqhull_mem(qh);
534 /* mem_r.c and qset_r.c are initialized */
535 qh_initqhull_buffers(qh);
536 qh_initthresholds(qh, qh->qhull_command);
537 if (qh->PROJECTinput || (qh->DELAUNAY && qh->PROJECTdelaunay))
538 qh_projectinput(qh);
539 if (qh->SCALEinput)
540 qh_scaleinput(qh);
541 if (qh->ROTATErandom >= 0) {
542 qh_randommatrix(qh, qh->gm_matrix, qh->hull_dim, qh->gm_row);
543 if (qh->DELAUNAY) {
544 int k, lastk= qh->hull_dim-1;
545 for (k=0; k < lastk; k++) {
546 qh->gm_row[k][lastk]= 0.0;
547 qh->gm_row[lastk][k]= 0.0;
548 }
549 qh->gm_row[lastk][lastk]= 1.0;
550 }
551 qh_gram_schmidt(qh, qh->hull_dim, qh->gm_row);
552 qh_rotateinput(qh, qh->gm_row);
553 }
554 } /* init_B */
555
556 /*-<a href="qh-globa_r.htm#TOC"
557 >-------------------------------</a><a name="init_qhull_command">-</a>
558
559 qh_init_qhull_command(qh, argc, argv )
560 build qh.qhull_command from argc/argv
561 Calls qh_exit if qhull_command is too short
562
563 returns:
564 a space-delimited string of options (just as typed)
565
566 notes:
567 makes option string easy to input and output
568
569 argc/argv may be 0/NULL
570 */
qh_init_qhull_command(qhT * qh,int argc,char * argv[])571 void qh_init_qhull_command(qhT *qh, int argc, char *argv[]) {
572
573 if (!qh_argv_to_command(argc, argv, qh->qhull_command, (int)sizeof(qh->qhull_command))){
574 /* Assumes qh.ferr is defined. */
575 qh_fprintf(qh, qh->ferr, 6033, "qhull input error: more than %d characters in command line.\n",
576 (int)sizeof(qh->qhull_command));
577 qh_exit(qh_ERRinput); /* error reported, can not use qh_errexit */
578 }
579 } /* init_qhull_command */
580
581 /*-<a href="qh-globa_r.htm#TOC"
582 >-------------------------------</a><a name="initflags">-</a>
583
584 qh_initflags(qh, commandStr )
585 set flags and initialized constants from commandStr
586 calls qh_exit() if qh.NOerrexit
587
588 returns:
589 sets qh.qhull_command to command if needed
590
591 notes:
592 ignores first word (e.g., 'qhull' in "qhull d")
593 use qh_strtol/strtod since strtol/strtod may or may not skip trailing spaces
594
595 see:
596 qh_initthresholds() continues processing of 'Pdn' and 'PDn'
597 'prompt' in unix_r.c for documentation
598
599 design:
600 for each space-delimited option group
601 if top-level option
602 check syntax
603 append appropriate option to option string
604 set appropriate global variable or append printFormat to print options
605 else
606 for each sub-option
607 check syntax
608 append appropriate option to option string
609 set appropriate global variable or append printFormat to print options
610 */
qh_initflags(qhT * qh,char * command)611 void qh_initflags(qhT *qh, char *command) {
612 int k, i, lastproject;
613 char *s= command, *t, *prev_s, *start, key, *lastwarning= NULL;
614 boolT isgeom= False, wasproject;
615 realT r;
616
617 if(qh->NOerrexit){
618 qh_fprintf(qh, qh->ferr, 6245, "qhull internal error (qh_initflags): qh.NOerrexit was not cleared before calling qh_initflags(). It should be cleared after setjmp(). Exit qhull.\n");
619 qh_exit(qh_ERRqhull);
620 }
621 #ifdef qh_RANDOMdist
622 qh->RANDOMfactor= qh_RANDOMdist;
623 qh_option(qh, "Random-qh_RANDOMdist", NULL, &qh->RANDOMfactor);
624 qh->RANDOMdist= True;
625 #endif
626 if (command <= &qh->qhull_command[0] || command > &qh->qhull_command[0] + sizeof(qh->qhull_command)) {
627 if (command != &qh->qhull_command[0]) {
628 *qh->qhull_command= '\0';
629 strncat(qh->qhull_command, command, sizeof(qh->qhull_command)-strlen(qh->qhull_command)-1);
630 }
631 while (*s && !isspace(*s)) /* skip program name */
632 s++;
633 }
634 while (*s) {
635 while (*s && isspace(*s))
636 s++;
637 if (*s == '-')
638 s++;
639 if (!*s)
640 break;
641 prev_s= s;
642 switch (*s++) {
643 case 'd':
644 qh_option(qh, "delaunay", NULL, NULL);
645 qh->DELAUNAY= True;
646 break;
647 case 'f':
648 qh_option(qh, "facets", NULL, NULL);
649 qh_appendprint(qh, qh_PRINTfacets);
650 break;
651 case 'i':
652 qh_option(qh, "incidence", NULL, NULL);
653 qh_appendprint(qh, qh_PRINTincidences);
654 break;
655 case 'm':
656 qh_option(qh, "mathematica", NULL, NULL);
657 qh_appendprint(qh, qh_PRINTmathematica);
658 break;
659 case 'n':
660 qh_option(qh, "normals", NULL, NULL);
661 qh_appendprint(qh, qh_PRINTnormals);
662 break;
663 case 'o':
664 qh_option(qh, "offFile", NULL, NULL);
665 qh_appendprint(qh, qh_PRINToff);
666 break;
667 case 'p':
668 qh_option(qh, "points", NULL, NULL);
669 qh_appendprint(qh, qh_PRINTpoints);
670 break;
671 case 's':
672 qh_option(qh, "summary", NULL, NULL);
673 qh->PRINTsummary= True;
674 break;
675 case 'v':
676 qh_option(qh, "voronoi", NULL, NULL);
677 qh->VORONOI= True;
678 qh->DELAUNAY= True;
679 break;
680 case 'A':
681 if (!isdigit(*s) && *s != '.' && *s != '-') {
682 qh_fprintf(qh, qh->ferr, 7002, "qhull input warning: no maximum cosine angle given for option 'An'. A1.0 is coplanar\n");
683 lastwarning= s-1;
684 }else {
685 if (*s == '-') {
686 qh->premerge_cos= -qh_strtod(s, &s);
687 qh_option(qh, "Angle-premerge-", NULL, &qh->premerge_cos);
688 qh->PREmerge= True;
689 }else {
690 qh->postmerge_cos= qh_strtod(s, &s);
691 qh_option(qh, "Angle-postmerge", NULL, &qh->postmerge_cos);
692 qh->POSTmerge= True;
693 }
694 qh->MERGING= True;
695 }
696 break;
697 case 'C':
698 if (!isdigit(*s) && *s != '.' && *s != '-') {
699 qh_fprintf(qh, qh->ferr, 7003, "qhull input warning: no centrum radius given for option 'Cn'\n");
700 lastwarning= s-1;
701 }else {
702 if (*s == '-') {
703 qh->premerge_centrum= -qh_strtod(s, &s);
704 qh_option(qh, "Centrum-premerge-", NULL, &qh->premerge_centrum);
705 qh->PREmerge= True;
706 }else {
707 qh->postmerge_centrum= qh_strtod(s, &s);
708 qh_option(qh, "Centrum-postmerge", NULL, &qh->postmerge_centrum);
709 qh->POSTmerge= True;
710 }
711 qh->MERGING= True;
712 }
713 break;
714 case 'E':
715 if (*s == '-') {
716 qh_fprintf(qh, qh->ferr, 6363, "qhull option error: expecting a positive number for maximum roundoff 'En'. Got '%s'\n", s-1);
717 qh_errexit(qh, qh_ERRinput, NULL, NULL);
718 }else if (!isdigit(*s)) {
719 qh_fprintf(qh, qh->ferr, 7005, "qhull option warning: no maximum roundoff given for option 'En'\n");
720 lastwarning= s-1;
721 }else {
722 qh->DISTround= qh_strtod(s, &s);
723 qh_option(qh, "Distance-roundoff", NULL, &qh->DISTround);
724 qh->SETroundoff= True;
725 }
726 break;
727 case 'H':
728 start= s;
729 qh->HALFspace= True;
730 qh_strtod(s, &t);
731 while (t > s) {
732 if (*t && !isspace(*t)) {
733 if (*t == ',')
734 t++;
735 else {
736 qh_fprintf(qh, qh->ferr, 6364, "qhull option error: expecting 'Hn,n,n,...' for feasible point of halfspace intersection. Got '%s'\n", start-1);
737 qh_errexit(qh, qh_ERRinput, NULL, NULL);
738 }
739 }
740 s= t;
741 qh_strtod(s, &t);
742 }
743 if (start < t) {
744 if (!(qh->feasible_string= (char *)calloc((size_t)(t-start+1), (size_t)1))) {
745 qh_fprintf(qh, qh->ferr, 6034, "qhull error: insufficient memory for 'Hn,n,n'\n");
746 qh_errexit(qh, qh_ERRmem, NULL, NULL);
747 }
748 strncpy(qh->feasible_string, start, (size_t)(t-start));
749 qh_option(qh, "Halfspace-about", NULL, NULL);
750 qh_option(qh, qh->feasible_string, NULL, NULL);
751 }else
752 qh_option(qh, "Halfspace", NULL, NULL);
753 break;
754 case 'R':
755 if (!isdigit(*s)) {
756 qh_fprintf(qh, qh->ferr, 7007, "qhull option warning: missing random perturbation for option 'Rn'\n");
757 lastwarning= s-1;
758 }else {
759 qh->RANDOMfactor= qh_strtod(s, &s);
760 qh_option(qh, "Random-perturb", NULL, &qh->RANDOMfactor);
761 qh->RANDOMdist= True;
762 }
763 break;
764 case 'V':
765 if (!isdigit(*s) && *s != '-') {
766 qh_fprintf(qh, qh->ferr, 7008, "qhull option warning: missing visible distance for option 'Vn'\n");
767 lastwarning= s-1;
768 }else {
769 qh->MINvisible= qh_strtod(s, &s);
770 qh_option(qh, "Visible", NULL, &qh->MINvisible);
771 }
772 break;
773 case 'U':
774 if (!isdigit(*s) && *s != '-') {
775 qh_fprintf(qh, qh->ferr, 7009, "qhull option warning: missing coplanar distance for option 'Un'\n");
776 lastwarning= s-1;
777 }else {
778 qh->MAXcoplanar= qh_strtod(s, &s);
779 qh_option(qh, "U-coplanar", NULL, &qh->MAXcoplanar);
780 }
781 break;
782 case 'W':
783 if (*s == '-') {
784 qh_fprintf(qh, qh->ferr, 6365, "qhull option error: expecting a positive number for outside width 'Wn'. Got '%s'\n", s-1);
785 qh_errexit(qh, qh_ERRinput, NULL, NULL);
786 }else if (!isdigit(*s)) {
787 qh_fprintf(qh, qh->ferr, 7011, "qhull option warning: missing outside width for option 'Wn'\n");
788 lastwarning= s-1;
789 }else {
790 qh->MINoutside= qh_strtod(s, &s);
791 qh_option(qh, "W-outside", NULL, &qh->MINoutside);
792 qh->APPROXhull= True;
793 }
794 break;
795 /************ sub menus ***************/
796 case 'F':
797 while (*s && !isspace(*s)) {
798 switch (*s++) {
799 case 'a':
800 qh_option(qh, "Farea", NULL, NULL);
801 qh_appendprint(qh, qh_PRINTarea);
802 qh->GETarea= True;
803 break;
804 case 'A':
805 qh_option(qh, "FArea-total", NULL, NULL);
806 qh->GETarea= True;
807 break;
808 case 'c':
809 qh_option(qh, "Fcoplanars", NULL, NULL);
810 qh_appendprint(qh, qh_PRINTcoplanars);
811 break;
812 case 'C':
813 qh_option(qh, "FCentrums", NULL, NULL);
814 qh_appendprint(qh, qh_PRINTcentrums);
815 break;
816 case 'd':
817 qh_option(qh, "Fd-cdd-in", NULL, NULL);
818 qh->CDDinput= True;
819 break;
820 case 'D':
821 qh_option(qh, "FD-cdd-out", NULL, NULL);
822 qh->CDDoutput= True;
823 break;
824 case 'F':
825 qh_option(qh, "FFacets-xridge", NULL, NULL);
826 qh_appendprint(qh, qh_PRINTfacets_xridge);
827 break;
828 case 'i':
829 qh_option(qh, "Finner", NULL, NULL);
830 qh_appendprint(qh, qh_PRINTinner);
831 break;
832 case 'I':
833 qh_option(qh, "FIDs", NULL, NULL);
834 qh_appendprint(qh, qh_PRINTids);
835 break;
836 case 'm':
837 qh_option(qh, "Fmerges", NULL, NULL);
838 qh_appendprint(qh, qh_PRINTmerges);
839 break;
840 case 'M':
841 qh_option(qh, "FMaple", NULL, NULL);
842 qh_appendprint(qh, qh_PRINTmaple);
843 break;
844 case 'n':
845 qh_option(qh, "Fneighbors", NULL, NULL);
846 qh_appendprint(qh, qh_PRINTneighbors);
847 break;
848 case 'N':
849 qh_option(qh, "FNeighbors-vertex", NULL, NULL);
850 qh_appendprint(qh, qh_PRINTvneighbors);
851 break;
852 case 'o':
853 qh_option(qh, "Fouter", NULL, NULL);
854 qh_appendprint(qh, qh_PRINTouter);
855 break;
856 case 'O':
857 if (qh->PRINToptions1st) {
858 qh_option(qh, "FOptions", NULL, NULL);
859 qh_appendprint(qh, qh_PRINToptions);
860 }else
861 qh->PRINToptions1st= True;
862 break;
863 case 'p':
864 qh_option(qh, "Fpoint-intersect", NULL, NULL);
865 qh_appendprint(qh, qh_PRINTpointintersect);
866 break;
867 case 'P':
868 qh_option(qh, "FPoint-nearest", NULL, NULL);
869 qh_appendprint(qh, qh_PRINTpointnearest);
870 break;
871 case 'Q':
872 qh_option(qh, "FQhull", NULL, NULL);
873 qh_appendprint(qh, qh_PRINTqhull);
874 break;
875 case 's':
876 qh_option(qh, "Fsummary", NULL, NULL);
877 qh_appendprint(qh, qh_PRINTsummary);
878 break;
879 case 'S':
880 qh_option(qh, "FSize", NULL, NULL);
881 qh_appendprint(qh, qh_PRINTsize);
882 qh->GETarea= True;
883 break;
884 case 't':
885 qh_option(qh, "Ftriangles", NULL, NULL);
886 qh_appendprint(qh, qh_PRINTtriangles);
887 break;
888 case 'v':
889 /* option set in qh_initqhull_globals */
890 qh_appendprint(qh, qh_PRINTvertices);
891 break;
892 case 'V':
893 qh_option(qh, "FVertex-average", NULL, NULL);
894 qh_appendprint(qh, qh_PRINTaverage);
895 break;
896 case 'x':
897 qh_option(qh, "Fxtremes", NULL, NULL);
898 qh_appendprint(qh, qh_PRINTextremes);
899 break;
900 default:
901 s--;
902 qh_fprintf(qh, qh->ferr, 7012, "qhull option warning: unknown 'F' output option 'F%c', skip to next space\n", (int)s[0]);
903 lastwarning= s-1;
904 while (*++s && !isspace(*s));
905 break;
906 }
907 }
908 break;
909 case 'G':
910 isgeom= True;
911 qh_appendprint(qh, qh_PRINTgeom);
912 while (*s && !isspace(*s)) {
913 switch (*s++) {
914 case 'a':
915 qh_option(qh, "Gall-points", NULL, NULL);
916 qh->PRINTdots= True;
917 break;
918 case 'c':
919 qh_option(qh, "Gcentrums", NULL, NULL);
920 qh->PRINTcentrums= True;
921 break;
922 case 'h':
923 qh_option(qh, "Gintersections", NULL, NULL);
924 qh->DOintersections= True;
925 break;
926 case 'i':
927 qh_option(qh, "Ginner", NULL, NULL);
928 qh->PRINTinner= True;
929 break;
930 case 'n':
931 qh_option(qh, "Gno-planes", NULL, NULL);
932 qh->PRINTnoplanes= True;
933 break;
934 case 'o':
935 qh_option(qh, "Gouter", NULL, NULL);
936 qh->PRINTouter= True;
937 break;
938 case 'p':
939 qh_option(qh, "Gpoints", NULL, NULL);
940 qh->PRINTcoplanar= True;
941 break;
942 case 'r':
943 qh_option(qh, "Gridges", NULL, NULL);
944 qh->PRINTridges= True;
945 break;
946 case 't':
947 qh_option(qh, "Gtransparent", NULL, NULL);
948 qh->PRINTtransparent= True;
949 break;
950 case 'v':
951 qh_option(qh, "Gvertices", NULL, NULL);
952 qh->PRINTspheres= True;
953 break;
954 case 'D':
955 if (!isdigit(*s)) {
956 qh_fprintf(qh, qh->ferr, 7004, "qhull option warning: missing dimension for option 'GDn'\n");
957 lastwarning= s-2;
958 }else {
959 if (qh->DROPdim >= 0) {
960 qh_fprintf(qh, qh->ferr, 7013, "qhull option warning: can only drop one dimension. Previous 'GD%d' ignored\n",
961 qh->DROPdim);
962 lastwarning= s-2;
963 }
964 qh->DROPdim= qh_strtol(s, &s);
965 qh_option(qh, "GDrop-dim", &qh->DROPdim, NULL);
966 }
967 break;
968 default:
969 s--;
970 qh_fprintf(qh, qh->ferr, 7014, "qhull option warning: unknown 'G' geomview option 'G%c', skip to next space\n", (int)s[0]);
971 lastwarning= s-1;
972 while (*++s && !isspace(*s));
973 break;
974 }
975 }
976 break;
977 case 'P':
978 while (*s && !isspace(*s)) {
979 switch (*s++) {
980 case 'd': case 'D': /* see qh_initthresholds() */
981 key= s[-1];
982 i= qh_strtol(s, &s);
983 r= 0;
984 if (*s == ':') {
985 s++;
986 r= qh_strtod(s, &s);
987 }
988 if (key == 'd')
989 qh_option(qh, "Pdrop-facets-dim-less", &i, &r);
990 else
991 qh_option(qh, "PDrop-facets-dim-more", &i, &r);
992 break;
993 case 'g':
994 qh_option(qh, "Pgood-facets", NULL, NULL);
995 qh->PRINTgood= True;
996 break;
997 case 'G':
998 qh_option(qh, "PGood-facet-neighbors", NULL, NULL);
999 qh->PRINTneighbors= True;
1000 break;
1001 case 'o':
1002 qh_option(qh, "Poutput-forced", NULL, NULL);
1003 qh->FORCEoutput= True;
1004 break;
1005 case 'p':
1006 qh_option(qh, "Pprecision-ignore", NULL, NULL);
1007 qh->PRINTprecision= False;
1008 break;
1009 case 'A':
1010 if (!isdigit(*s)) {
1011 qh_fprintf(qh, qh->ferr, 7006, "qhull option warning: missing facet count for keep area option 'PAn'\n");
1012 lastwarning= s-2;
1013 }else {
1014 qh->KEEParea= qh_strtol(s, &s);
1015 qh_option(qh, "PArea-keep", &qh->KEEParea, NULL);
1016 qh->GETarea= True;
1017 }
1018 break;
1019 case 'F':
1020 if (!isdigit(*s)) {
1021 qh_fprintf(qh, qh->ferr, 7010, "qhull option warning: missing facet area for option 'PFn'\n");
1022 lastwarning= s-2;
1023 }else {
1024 qh->KEEPminArea= qh_strtod(s, &s);
1025 qh_option(qh, "PFacet-area-keep", NULL, &qh->KEEPminArea);
1026 qh->GETarea= True;
1027 }
1028 break;
1029 case 'M':
1030 if (!isdigit(*s)) {
1031 qh_fprintf(qh, qh->ferr, 7090, "qhull option warning: missing merge count for option 'PMn'\n");
1032 lastwarning= s-2;
1033 }else {
1034 qh->KEEPmerge= qh_strtol(s, &s);
1035 qh_option(qh, "PMerge-keep", &qh->KEEPmerge, NULL);
1036 }
1037 break;
1038 default:
1039 s--;
1040 qh_fprintf(qh, qh->ferr, 7015, "qhull option warning: unknown 'P' print option 'P%c', skip to next space\n", (int)s[0]);
1041 lastwarning= s-1;
1042 while (*++s && !isspace(*s));
1043 break;
1044 }
1045 }
1046 break;
1047 case 'Q':
1048 lastproject= -1;
1049 while (*s && !isspace(*s)) {
1050 switch (*s++) {
1051 case 'a':
1052 qh_option(qh, "Qallow-short", NULL, NULL);
1053 qh->ALLOWshort= True;
1054 break;
1055 case 'b': case 'B': /* handled by qh_initthresholds */
1056 key= s[-1];
1057 if (key == 'b' && *s == 'B') {
1058 s++;
1059 r= qh_DEFAULTbox;
1060 qh->SCALEinput= True;
1061 qh_option(qh, "QbBound-unit-box", NULL, &r);
1062 break;
1063 }
1064 if (key == 'b' && *s == 'b') {
1065 s++;
1066 qh->SCALElast= True;
1067 qh_option(qh, "Qbbound-last", NULL, NULL);
1068 break;
1069 }
1070 k= qh_strtol(s, &s);
1071 r= 0.0;
1072 wasproject= False;
1073 if (*s == ':') {
1074 s++;
1075 if ((r= qh_strtod(s, &s)) == 0.0) {
1076 t= s; /* need true dimension for memory allocation */
1077 while (*t && !isspace(*t)) {
1078 if (toupper(*t++) == 'B'
1079 && k == qh_strtol(t, &t)
1080 && *t++ == ':'
1081 && qh_strtod(t, &t) == 0.0) {
1082 qh->PROJECTinput++;
1083 trace2((qh, qh->ferr, 2004, "qh_initflags: project dimension %d\n", k));
1084 qh_option(qh, "Qb-project-dim", &k, NULL);
1085 wasproject= True;
1086 lastproject= k;
1087 break;
1088 }
1089 }
1090 }
1091 }
1092 if (!wasproject) {
1093 if (lastproject == k && r == 0.0)
1094 lastproject= -1; /* doesn't catch all possible sequences */
1095 else if (key == 'b') {
1096 qh->SCALEinput= True;
1097 if (r == 0.0)
1098 r= -qh_DEFAULTbox;
1099 qh_option(qh, "Qbound-dim-low", &k, &r);
1100 }else {
1101 qh->SCALEinput= True;
1102 if (r == 0.0)
1103 r= qh_DEFAULTbox;
1104 qh_option(qh, "QBound-dim-high", &k, &r);
1105 }
1106 }
1107 break;
1108 case 'c':
1109 qh_option(qh, "Qcoplanar-keep", NULL, NULL);
1110 qh->KEEPcoplanar= True;
1111 break;
1112 case 'f':
1113 qh_option(qh, "Qfurthest-outside", NULL, NULL);
1114 qh->BESToutside= True;
1115 break;
1116 case 'g':
1117 qh_option(qh, "Qgood-facets-only", NULL, NULL);
1118 qh->ONLYgood= True;
1119 break;
1120 case 'i':
1121 qh_option(qh, "Qinterior-keep", NULL, NULL);
1122 qh->KEEPinside= True;
1123 break;
1124 case 'm':
1125 qh_option(qh, "Qmax-outside-only", NULL, NULL);
1126 qh->ONLYmax= True;
1127 break;
1128 case 'r':
1129 qh_option(qh, "Qrandom-outside", NULL, NULL);
1130 qh->RANDOMoutside= True;
1131 break;
1132 case 's':
1133 qh_option(qh, "Qsearch-initial-simplex", NULL, NULL);
1134 qh->ALLpoints= True;
1135 break;
1136 case 't':
1137 qh_option(qh, "Qtriangulate", NULL, NULL);
1138 qh->TRIangulate= True;
1139 break;
1140 case 'T':
1141 qh_option(qh, "QTestPoints", NULL, NULL);
1142 if (!isdigit(*s)) {
1143 qh_fprintf(qh, qh->ferr, 7091, "qhull option warning: missing number of test points for option 'QTn'\n");
1144 lastwarning= s-2;
1145 }else {
1146 qh->TESTpoints= qh_strtol(s, &s);
1147 qh_option(qh, "QTestPoints", &qh->TESTpoints, NULL);
1148 }
1149 break;
1150 case 'u':
1151 qh_option(qh, "QupperDelaunay", NULL, NULL);
1152 qh->UPPERdelaunay= True;
1153 break;
1154 case 'v':
1155 qh_option(qh, "Qvertex-neighbors-convex", NULL, NULL);
1156 qh->TESTvneighbors= True;
1157 break;
1158 case 'x':
1159 qh_option(qh, "Qxact-merge", NULL, NULL);
1160 qh->MERGEexact= True;
1161 break;
1162 case 'z':
1163 qh_option(qh, "Qz-infinity-point", NULL, NULL);
1164 qh->ATinfinity= True;
1165 break;
1166 case '0':
1167 qh_option(qh, "Q0-no-premerge", NULL, NULL);
1168 qh->NOpremerge= True;
1169 break;
1170 case '1':
1171 if (!isdigit(*s)) {
1172 qh_option(qh, "Q1-angle-merge", NULL, NULL);
1173 qh->ANGLEmerge= True;
1174 break;
1175 }
1176 switch (*s++) {
1177 case '0':
1178 qh_option(qh, "Q10-no-narrow", NULL, NULL);
1179 qh->NOnarrow= True;
1180 break;
1181 case '1':
1182 qh_option(qh, "Q11-trinormals Qtriangulate", NULL, NULL);
1183 qh->TRInormals= True;
1184 qh->TRIangulate= True;
1185 break;
1186 case '2':
1187 qh_option(qh, "Q12-allow-wide", NULL, NULL);
1188 qh->ALLOWwide= True;
1189 break;
1190 case '4':
1191 #ifndef qh_NOmerge
1192 qh_option(qh, "Q14-merge-pinched-vertices", NULL, NULL);
1193 qh->MERGEpinched= True;
1194 #else
1195 /* ignore 'Q14' for q_benchmark testing of difficult cases for Qhull */
1196 qh_fprintf(qh, qh->ferr, 7099, "qhull option warning: option 'Q14-merge-pinched' disabled due to qh_NOmerge\n");
1197 #endif
1198 break;
1199 case '7':
1200 qh_option(qh, "Q15-check-duplicates", NULL, NULL);
1201 qh->CHECKduplicates= True;
1202 break;
1203 default:
1204 s--;
1205 qh_fprintf(qh, qh->ferr, 7016, "qhull option warning: unknown 'Q' qhull option 'Q1%c', skip to next space\n", (int)s[0]);
1206 lastwarning= s-1;
1207 while (*++s && !isspace(*s));
1208 break;
1209 }
1210 break;
1211 case '2':
1212 qh_option(qh, "Q2-no-merge-independent", NULL, NULL);
1213 qh->MERGEindependent= False;
1214 goto LABELcheckdigit;
1215 break; /* no gcc warnings */
1216 case '3':
1217 qh_option(qh, "Q3-no-merge-vertices", NULL, NULL);
1218 qh->MERGEvertices= False;
1219 LABELcheckdigit:
1220 if (isdigit(*s)) {
1221 qh_fprintf(qh, qh->ferr, 7017, "qhull option warning: can not follow '1', '2', or '3' with a digit. 'Q%c%c' skipped\n", *(s-1), *s);
1222 lastwarning= s-2;
1223 s++;
1224 }
1225 break;
1226 case '4':
1227 qh_option(qh, "Q4-avoid-old-into-new", NULL, NULL);
1228 qh->AVOIDold= True;
1229 break;
1230 case '5':
1231 qh_option(qh, "Q5-no-check-outer", NULL, NULL);
1232 qh->SKIPcheckmax= True;
1233 break;
1234 case '6':
1235 qh_option(qh, "Q6-no-concave-merge", NULL, NULL);
1236 qh->SKIPconvex= True;
1237 break;
1238 case '7':
1239 qh_option(qh, "Q7-no-breadth-first", NULL, NULL);
1240 qh->VIRTUALmemory= True;
1241 break;
1242 case '8':
1243 qh_option(qh, "Q8-no-near-inside", NULL, NULL);
1244 qh->NOnearinside= True;
1245 break;
1246 case '9':
1247 qh_option(qh, "Q9-pick-furthest", NULL, NULL);
1248 qh->PICKfurthest= True;
1249 break;
1250 case 'G':
1251 i= qh_strtol(s, &t);
1252 if (qh->GOODpoint) {
1253 qh_fprintf(qh, qh->ferr, 7018, "qhull option warning: good point already defined for option 'QGn'. Ignored\n");
1254 lastwarning= s-2;
1255 }else if (s == t) {
1256 qh_fprintf(qh, qh->ferr, 7019, "qhull option warning: missing good point id for option 'QGn'. Ignored\n");
1257 lastwarning= s-2;
1258 }else if (i < 0 || *s == '-') {
1259 qh->GOODpoint= i-1;
1260 qh_option(qh, "QGood-if-dont-see-point", &i, NULL);
1261 }else {
1262 qh->GOODpoint= i+1;
1263 qh_option(qh, "QGood-if-see-point", &i, NULL);
1264 }
1265 s= t;
1266 break;
1267 case 'J':
1268 if (!isdigit(*s) && *s != '-')
1269 qh->JOGGLEmax= 0.0;
1270 else {
1271 qh->JOGGLEmax= (realT) qh_strtod(s, &s);
1272 qh_option(qh, "QJoggle", NULL, &qh->JOGGLEmax);
1273 }
1274 break;
1275 case 'R':
1276 if (!isdigit(*s) && *s != '-') {
1277 qh_fprintf(qh, qh->ferr, 7020, "qhull option warning: missing random seed for option 'QRn'\n");
1278 lastwarning= s-2;
1279 }else {
1280 qh->ROTATErandom= i= qh_strtol(s, &s);
1281 if (i > 0)
1282 qh_option(qh, "QRotate-id", &i, NULL );
1283 else if (i < -1)
1284 qh_option(qh, "QRandom-seed", &i, NULL );
1285 }
1286 break;
1287 case 'V':
1288 i= qh_strtol(s, &t);
1289 if (qh->GOODvertex) {
1290 qh_fprintf(qh, qh->ferr, 7021, "qhull option warning: good vertex already defined for option 'QVn'. Ignored\n");
1291 lastwarning= s-2;
1292 }else if (s == t) {
1293 qh_fprintf(qh, qh->ferr, 7022, "qhull option warning: no good point id given for option 'QVn'. Ignored\n");
1294 lastwarning= s-2;
1295 }else if (i < 0) {
1296 qh->GOODvertex= i - 1;
1297 qh_option(qh, "QV-good-facets-not-point", &i, NULL);
1298 }else {
1299 qh_option(qh, "QV-good-facets-point", &i, NULL);
1300 qh->GOODvertex= i + 1;
1301 }
1302 s= t;
1303 break;
1304 case 'w':
1305 qh_option(qh, "Qwarn-allow", NULL, NULL);
1306 qh->ALLOWwarning= True;
1307 break;
1308 default:
1309 s--;
1310 qh_fprintf(qh, qh->ferr, 7023, "qhull option warning: unknown 'Q' qhull option 'Q%c', skip to next space\n", (int)s[0]);
1311 lastwarning= s-1;
1312 while (*++s && !isspace(*s));
1313 break;
1314 }
1315 }
1316 break;
1317 case 'T':
1318 while (*s && !isspace(*s)) {
1319 if (isdigit(*s) || *s == '-')
1320 qh->IStracing= qh_strtol(s, &s);
1321 else switch (*s++) {
1322 case 'a':
1323 qh_option(qh, "Tannotate-output", NULL, NULL);
1324 qh->ANNOTATEoutput= True;
1325 break;
1326 case 'c':
1327 qh_option(qh, "Tcheck-frequently", NULL, NULL);
1328 qh->CHECKfrequently= True;
1329 break;
1330 case 'f':
1331 qh_option(qh, "Tflush", NULL, NULL);
1332 qh->FLUSHprint= True;
1333 break;
1334 case 's':
1335 qh_option(qh, "Tstatistics", NULL, NULL);
1336 qh->PRINTstatistics= True;
1337 break;
1338 case 'v':
1339 qh_option(qh, "Tverify", NULL, NULL);
1340 qh->VERIFYoutput= True;
1341 break;
1342 case 'z':
1343 if (qh->ferr == qh_FILEstderr) {
1344 /* The C++ interface captures the output in qh_fprint_qhull() */
1345 qh_option(qh, "Tz-stdout", NULL, NULL);
1346 qh->USEstdout= True;
1347 }else if (!qh->fout) {
1348 qh_fprintf(qh, qh->ferr, 7024, "qhull option warning: output file undefined(stdout). Option 'Tz' ignored.\n");
1349 lastwarning= s-2;
1350 }else {
1351 qh_option(qh, "Tz-stdout", NULL, NULL);
1352 qh->USEstdout= True;
1353 qh->ferr= qh->fout;
1354 qh->qhmem.ferr= qh->fout;
1355 }
1356 break;
1357 case 'C':
1358 if (!isdigit(*s)) {
1359 qh_fprintf(qh, qh->ferr, 7025, "qhull option warning: missing point id for cone for trace option 'TCn'\n");
1360 lastwarning= s-2;
1361 }else {
1362 i= qh_strtol(s, &s);
1363 qh_option(qh, "TCone-stop", &i, NULL);
1364 qh->STOPcone= i + 1;
1365 }
1366 break;
1367 case 'F':
1368 if (!isdigit(*s)) {
1369 qh_fprintf(qh, qh->ferr, 7026, "qhull option warning: missing frequency count for trace option 'TFn'\n");
1370 lastwarning= s-2;
1371 }else {
1372 qh->REPORTfreq= qh_strtol(s, &s);
1373 qh_option(qh, "TFacet-log", &qh->REPORTfreq, NULL);
1374 qh->REPORTfreq2= qh->REPORTfreq/2; /* for tracemerging() */
1375 }
1376 break;
1377 case 'I':
1378 while (isspace(*s))
1379 s++;
1380 t= qh_skipfilename(qh, s);
1381 {
1382 char filename[qh_FILENAMElen];
1383
1384 qh_copyfilename(qh, filename, (int)sizeof(filename), s, (int)(t-s)); /* WARN64 */
1385 s= t;
1386 if (!freopen(filename, "r", stdin)) {
1387 qh_fprintf(qh, qh->ferr, 6041, "qhull option error: cannot open 'TI' file \"%s\"\n", filename);
1388 qh_errexit(qh, qh_ERRinput, NULL, NULL);
1389 }else {
1390 qh_option(qh, "TInput-file", NULL, NULL);
1391 qh_option(qh, filename, NULL, NULL);
1392 }
1393 }
1394 break;
1395 case 'O':
1396 while (isspace(*s))
1397 s++;
1398 t= qh_skipfilename(qh, s);
1399 {
1400 char filename[qh_FILENAMElen];
1401
1402 qh_copyfilename(qh, filename, (int)sizeof(filename), s, (int)(t-s)); /* WARN64 */
1403 if (!qh->fout) {
1404 qh_fprintf(qh, qh->ferr, 7092, "qhull option warning: qh.fout was not set by caller of qh_initflags. Cannot use option 'TO' to redirect output. Ignoring option 'TO'\n");
1405 lastwarning= s-2;
1406 }else if (!freopen(filename, "w", qh->fout)) {
1407 qh_fprintf(qh, qh->ferr, 6044, "qhull option error: cannot open file \"%s\" for writing as option 'TO'. It is already in use or read-only\n", filename);
1408 qh_errexit(qh, qh_ERRinput, NULL, NULL);
1409 }else {
1410 qh_option(qh, "TOutput-file", NULL, NULL);
1411 qh_option(qh, filename, NULL, NULL);
1412 }
1413 s= t;
1414 }
1415 break;
1416 case 'A':
1417 if (!isdigit(*s)) {
1418 qh_fprintf(qh, qh->ferr, 7093, "qhull option warning: missing count of added points for trace option 'TAn'\n");
1419 lastwarning= s-2;
1420 }else {
1421 i= qh_strtol(s, &t);
1422 qh->STOPadd= i + 1;
1423 qh_option(qh, "TA-stop-add", &i, NULL);
1424 }
1425 s= t;
1426 break;
1427 case 'P':
1428 if (*s == '-') {
1429 if (s[1] == '1' && !isdigit(s[2])) {
1430 s += 2;
1431 qh->TRACEpoint= qh_IDunknown; /* qh_buildhull done */
1432 qh_option(qh, "Trace-point", &qh->TRACEpoint, NULL);
1433 }else {
1434 qh_fprintf(qh, qh->ferr, 7100, "qhull option warning: negative point id for trace option 'TPn'. Expecting 'TP-1' for tracing after qh_buildhull and qh_postmerge\n");
1435 lastwarning= s-2;
1436 while (isdigit(*(++s)))
1437 ; /* skip digits */
1438 }
1439 }else if (!isdigit(*s)) {
1440 qh_fprintf(qh, qh->ferr, 7029, "qhull option warning: missing point id or -1 for trace option 'TPn'\n");
1441 lastwarning= s-2;
1442 }else {
1443 qh->TRACEpoint= qh_strtol(s, &s);
1444 qh_option(qh, "Trace-point", &qh->TRACEpoint, NULL);
1445 }
1446 break;
1447 case 'M':
1448 if (!isdigit(*s)) {
1449 qh_fprintf(qh, qh->ferr, 7030, "qhull option warning: missing merge id for trace option 'TMn'\n");
1450 lastwarning= s-2;
1451 }else {
1452 qh->TRACEmerge= qh_strtol(s, &s);
1453 qh_option(qh, "Trace-merge", &qh->TRACEmerge, NULL);
1454 }
1455 break;
1456 case 'R':
1457 if (!isdigit(*s)) {
1458 qh_fprintf(qh, qh->ferr, 7031, "qhull option warning: missing rerun count for trace option 'TRn'\n");
1459 lastwarning= s-2;
1460 }else {
1461 qh->RERUN= qh_strtol(s, &s);
1462 qh_option(qh, "TRerun", &qh->RERUN, NULL);
1463 }
1464 break;
1465 case 'V':
1466 i= qh_strtol(s, &t);
1467 if (s == t) {
1468 qh_fprintf(qh, qh->ferr, 7032, "qhull option warning: missing furthest point id for trace option 'TVn'\n");
1469 lastwarning= s-2;
1470 }else if (i < 0) {
1471 qh->STOPpoint= i - 1;
1472 qh_option(qh, "TV-stop-before-point", &i, NULL);
1473 }else {
1474 qh->STOPpoint= i + 1;
1475 qh_option(qh, "TV-stop-after-point", &i, NULL);
1476 }
1477 s= t;
1478 break;
1479 case 'W':
1480 if (!isdigit(*s)) {
1481 qh_fprintf(qh, qh->ferr, 7033, "qhull option warning: missing max width for trace option 'TWn'\n");
1482 lastwarning= s-2;
1483 }else {
1484 qh->TRACEdist= (realT) qh_strtod(s, &s);
1485 qh_option(qh, "TWide-trace", NULL, &qh->TRACEdist);
1486 }
1487 break;
1488 default:
1489 s--;
1490 qh_fprintf(qh, qh->ferr, 7034, "qhull option warning: unknown 'T' trace option 'T%c', skip to next space\n", (int)s[0]);
1491 lastwarning= s-2;
1492 while (*++s && !isspace(*s));
1493 break;
1494 }
1495 }
1496 break;
1497 default:
1498 qh_fprintf(qh, qh->ferr, 7094, "qhull option warning: unknown option '%c'(%x)\n",
1499 (int)s[-1], (int)s[-1]);
1500 lastwarning= s-2;
1501 break;
1502 }
1503 if (s-1 == prev_s && *s && !isspace(*s)) {
1504 qh_fprintf(qh, qh->ferr, 7036, "qhull option warning: missing space after option '%c'(%x), reserved for sub-options, ignoring '%c' options to next space\n",
1505 (int)*prev_s, (int)*prev_s, (int)*prev_s);
1506 lastwarning= s-1;
1507 while (*s && !isspace(*s))
1508 s++;
1509 }
1510 }
1511 if (qh->STOPcone && qh->JOGGLEmax < REALmax/2) {
1512 qh_fprintf(qh, qh->ferr, 7078, "qhull option warning: 'TCn' (stopCone) ignored when used with 'QJn' (joggle)\n");
1513 lastwarning= command;
1514 }
1515 if (isgeom && !qh->FORCEoutput && qh->PRINTout[1]) {
1516 qh_fprintf(qh, qh->ferr, 7037, "qhull option warning: additional output formats ('Fc',etc.) are not compatible with Geomview ('G'). Use option 'Po' to override\n");
1517 lastwarning= command;
1518 }
1519 if (lastwarning && !qh->ALLOWwarning) {
1520 qh_fprintf(qh, qh->ferr, 6035, "qhull option error: see previous warnings, use 'Qw' to override: '%s' (last offset %d)\n",
1521 command, (int)(lastwarning-command));
1522 qh_errexit(qh, qh_ERRinput, NULL, NULL);
1523 }
1524 trace4((qh, qh->ferr, 4093, "qh_initflags: option flags initialized\n"));
1525 /* set derived values in qh_initqhull_globals */
1526 } /* initflags */
1527
1528
1529 /*-<a href="qh-globa_r.htm#TOC"
1530 >-------------------------------</a><a name="initqhull_buffers">-</a>
1531
1532 qh_initqhull_buffers(qh)
1533 initialize global memory buffers
1534
1535 notes:
1536 must match qh_freebuffers()
1537 */
qh_initqhull_buffers(qhT * qh)1538 void qh_initqhull_buffers(qhT *qh) {
1539 int k;
1540
1541 qh->TEMPsize= (qh->qhmem.LASTsize - (int)sizeof(setT))/SETelemsize;
1542 if (qh->TEMPsize <= 0 || qh->TEMPsize > qh->qhmem.LASTsize)
1543 qh->TEMPsize= 8; /* e.g., if qh_NOmem */
1544 qh->other_points= qh_setnew(qh, qh->TEMPsize);
1545 qh->del_vertices= qh_setnew(qh, qh->TEMPsize);
1546 qh->coplanarfacetset= qh_setnew(qh, qh->TEMPsize);
1547 qh->NEARzero= (realT *)qh_memalloc(qh, qh->hull_dim * (int)sizeof(realT));
1548 qh->lower_threshold= (realT *)qh_memalloc(qh, (qh->input_dim+1) * (int)sizeof(realT));
1549 qh->upper_threshold= (realT *)qh_memalloc(qh, (qh->input_dim+1) * (int)sizeof(realT));
1550 qh->lower_bound= (realT *)qh_memalloc(qh, (qh->input_dim+1) * (int)sizeof(realT));
1551 qh->upper_bound= (realT *)qh_memalloc(qh, (qh->input_dim+1) * (int)sizeof(realT));
1552 for (k=qh->input_dim+1; k--; ) { /* duplicated in qh_initqhull_buffers and qh_clear_outputflags */
1553 qh->lower_threshold[k]= -REALmax;
1554 qh->upper_threshold[k]= REALmax;
1555 qh->lower_bound[k]= -REALmax;
1556 qh->upper_bound[k]= REALmax;
1557 }
1558 qh->gm_matrix= (coordT *)qh_memalloc(qh, (qh->hull_dim+1) * qh->hull_dim * (int)sizeof(coordT));
1559 qh->gm_row= (coordT **)qh_memalloc(qh, (qh->hull_dim+1) * (int)sizeof(coordT *));
1560 } /* initqhull_buffers */
1561
1562 /*-<a href="qh-globa_r.htm#TOC"
1563 >-------------------------------</a><a name="initqhull_globals">-</a>
1564
1565 qh_initqhull_globals(qh, points, numpoints, dim, ismalloc )
1566 initialize globals
1567 if ismalloc
1568 points were malloc'd and qhull should free at end
1569
1570 returns:
1571 sets qh.first_point, num_points, input_dim, hull_dim and others
1572 seeds random number generator (seed=1 if tracing)
1573 modifies qh.hull_dim if ((qh.DELAUNAY and qh.PROJECTdelaunay) or qh.PROJECTinput)
1574 adjust user flags as needed
1575 also checks DIM3 dependencies and constants
1576
1577 notes:
1578 do not use qh_point() since an input transformation may move them elsewhere
1579 qh_initqhull_start() sets default values for non-zero globals
1580 consider duplicate error checks in qh_readpoints. It is called before qh_initqhull_globals
1581
1582 design:
1583 initialize points array from input arguments
1584 test for qh.ZEROcentrum
1585 (i.e., use opposite vertex instead of cetrum for convexity testing)
1586 initialize qh.CENTERtype, qh.normal_size,
1587 qh.center_size, qh.TRACEpoint/level,
1588 initialize and test random numbers
1589 qh_initqhull_outputflags() -- adjust and test output flags
1590 */
qh_initqhull_globals(qhT * qh,coordT * points,int numpoints,int dim,boolT ismalloc)1591 void qh_initqhull_globals(qhT *qh, coordT *points, int numpoints, int dim, boolT ismalloc) {
1592 int seed, pointsneeded, extra= 0, i, randi, k;
1593 realT randr;
1594 realT factorial;
1595
1596 time_t timedata;
1597
1598 trace0((qh, qh->ferr, 13, "qh_initqhull_globals: for %s | %s\n", qh->rbox_command,
1599 qh->qhull_command));
1600 if (numpoints < 1 || numpoints > qh_POINTSmax) {
1601 qh_fprintf(qh, qh->ferr, 6412, "qhull input error (qh_initqhull_globals): expecting between 1 and %d points. Got %d %d-d points\n",
1602 qh_POINTSmax, numpoints, dim);
1603 qh_errexit(qh, qh_ERRinput, NULL, NULL);
1604 /* same error message in qh_readpoints */
1605 }
1606 qh->POINTSmalloc= ismalloc;
1607 qh->first_point= points;
1608 qh->num_points= numpoints;
1609 qh->hull_dim= qh->input_dim= dim;
1610 if (!qh->NOpremerge && !qh->MERGEexact && !qh->PREmerge && qh->JOGGLEmax > REALmax/2) {
1611 qh->MERGING= True;
1612 if (qh->hull_dim <= 4) {
1613 qh->PREmerge= True;
1614 qh_option(qh, "_pre-merge", NULL, NULL);
1615 }else {
1616 qh->MERGEexact= True;
1617 qh_option(qh, "Qxact-merge", NULL, NULL);
1618 }
1619 }else if (qh->MERGEexact)
1620 qh->MERGING= True;
1621 if (qh->NOpremerge && (qh->MERGEexact || qh->PREmerge))
1622 qh_fprintf(qh, qh->ferr, 7095, "qhull option warning: 'Q0-no-premerge' ignored due to exact merge ('Qx') or pre-merge ('C-n' or 'A-n')\n");
1623 if (!qh->NOpremerge && qh->JOGGLEmax > REALmax/2) {
1624 #ifdef qh_NOmerge
1625 qh->JOGGLEmax= 0.0;
1626 #endif
1627 }
1628 if (qh->TRIangulate && qh->JOGGLEmax < REALmax/2 && !qh->PREmerge && !qh->POSTmerge && qh->PRINTprecision)
1629 qh_fprintf(qh, qh->ferr, 7038, "qhull option warning: joggle ('QJ') produces simplicial output (i.e., triangles in 2-D). Unless merging is requested, option 'Qt' has no effect\n");
1630 if (qh->JOGGLEmax < REALmax/2 && qh->DELAUNAY && !qh->SCALEinput && !qh->SCALElast) {
1631 qh->SCALElast= True;
1632 qh_option(qh, "Qbbound-last-qj", NULL, NULL);
1633 }
1634 if (qh->MERGING && !qh->POSTmerge && qh->premerge_cos > REALmax/2
1635 && qh->premerge_centrum == 0.0) {
1636 qh->ZEROcentrum= True;
1637 qh->ZEROall_ok= True;
1638 qh_option(qh, "_zero-centrum", NULL, NULL);
1639 }
1640 if (qh->JOGGLEmax < REALmax/2 && REALepsilon > 2e-8 && qh->PRINTprecision)
1641 qh_fprintf(qh, qh->ferr, 7039, "qhull warning: real epsilon, %2.2g, is probably too large for joggle('QJn')\nRecompile with double precision reals(see user_r.h).\n",
1642 REALepsilon);
1643 #ifdef qh_NOmerge
1644 if (qh->MERGING) {
1645 qh_fprintf(qh, qh->ferr, 6045, "qhull option error: merging not installed (qh_NOmerge) for 'Qx', 'Cn' or 'An')\n");
1646 qh_errexit(qh, qh_ERRinput, NULL, NULL);
1647 }
1648 #endif
1649 if (qh->DELAUNAY && qh->KEEPcoplanar && !qh->KEEPinside) {
1650 qh->KEEPinside= True;
1651 qh_option(qh, "Qinterior-keep", NULL, NULL);
1652 }
1653 if (qh->VORONOI && !qh->DELAUNAY) {
1654 qh_fprintf(qh, qh->ferr, 6038, "qhull internal error (qh_initqhull_globals): if qh.VORONOI is set, qh.DELAUNAY must be set. Qhull constructs the Delaunay triangulation in order to compute the Voronoi diagram\n");
1655 qh_errexit(qh, qh_ERRqhull, NULL, NULL);
1656 }
1657 if (qh->DELAUNAY && qh->HALFspace) {
1658 qh_fprintf(qh, qh->ferr, 6046, "qhull option error: can not use Delaunay('d') or Voronoi('v') with halfspace intersection('H')\n");
1659 qh_errexit(qh, qh_ERRinput, NULL, NULL);
1660 /* same error message in qh_readpoints */
1661 }
1662 if (!qh->DELAUNAY && (qh->UPPERdelaunay || qh->ATinfinity)) {
1663 qh_fprintf(qh, qh->ferr, 6047, "qhull option error: use upper-Delaunay('Qu') or infinity-point('Qz') with Delaunay('d') or Voronoi('v')\n");
1664 qh_errexit(qh, qh_ERRinput, NULL, NULL);
1665 }
1666 if (qh->UPPERdelaunay && qh->ATinfinity) {
1667 qh_fprintf(qh, qh->ferr, 6048, "qhull option error: can not use infinity-point('Qz') with upper-Delaunay('Qu')\n");
1668 qh_errexit(qh, qh_ERRinput, NULL, NULL);
1669 }
1670 if (qh->MERGEpinched && qh->ONLYgood) {
1671 qh_fprintf(qh, qh->ferr, 6362, "qhull option error: can not use merge-pinched-vertices ('Q14') with good-facets-only ('Qg')\n");
1672 qh_errexit(qh, qh_ERRinput, NULL, NULL);
1673 }
1674 if (qh->MERGEpinched && qh->hull_dim == 2) {
1675 trace2((qh, qh->ferr, 2108, "qh_initqhull_globals: disable qh.MERGEpinched for 2-d. It has no effect"))
1676 qh->MERGEpinched= False;
1677 }
1678 if (qh->SCALElast && !qh->DELAUNAY && qh->PRINTprecision)
1679 qh_fprintf(qh, qh->ferr, 7040, "qhull option warning: option 'Qbb' (scale-last-coordinate) is normally used with 'd' or 'v'\n");
1680 qh->DOcheckmax= (!qh->SKIPcheckmax && (qh->MERGING || qh->APPROXhull));
1681 qh->KEEPnearinside= (qh->DOcheckmax && !(qh->KEEPinside && qh->KEEPcoplanar)
1682 && !qh->NOnearinside);
1683 if (qh->MERGING)
1684 qh->CENTERtype= qh_AScentrum;
1685 else if (qh->VORONOI)
1686 qh->CENTERtype= qh_ASvoronoi;
1687 if (qh->TESTvneighbors && !qh->MERGING) {
1688 qh_fprintf(qh, qh->ferr, 6049, "qhull option error: test vertex neighbors('Qv') needs a merge option\n");
1689 qh_errexit(qh, qh_ERRinput, NULL ,NULL);
1690 }
1691 if (qh->PROJECTinput || (qh->DELAUNAY && qh->PROJECTdelaunay)) {
1692 qh->hull_dim -= qh->PROJECTinput;
1693 if (qh->DELAUNAY) {
1694 qh->hull_dim++;
1695 if (qh->ATinfinity)
1696 extra= 1;
1697 }
1698 }
1699 if (qh->hull_dim <= 1) {
1700 qh_fprintf(qh, qh->ferr, 6050, "qhull error: dimension %d must be > 1\n", qh->hull_dim);
1701 qh_errexit(qh, qh_ERRinput, NULL, NULL);
1702 }
1703 for (k=2, factorial=1.0; k < qh->hull_dim; k++)
1704 factorial *= k;
1705 qh->AREAfactor= 1.0 / factorial;
1706 trace2((qh, qh->ferr, 2005, "qh_initqhull_globals: initialize globals. input_dim %d, numpoints %d, malloc? %d, projected %d to hull_dim %d\n",
1707 qh->input_dim, numpoints, ismalloc, qh->PROJECTinput, qh->hull_dim));
1708 qh->normal_size= qh->hull_dim * (int)sizeof(coordT);
1709 qh->center_size= qh->normal_size - (int)sizeof(coordT);
1710 pointsneeded= qh->hull_dim+1;
1711 if (qh->hull_dim > qh_DIMmergeVertex) {
1712 qh->MERGEvertices= False;
1713 qh_option(qh, "Q3-no-merge-vertices-dim-high", NULL, NULL);
1714 }
1715 if (qh->GOODpoint)
1716 pointsneeded++;
1717 #ifdef qh_NOtrace
1718 if (qh->IStracing || qh->TRACEmerge || qh->TRACEpoint != qh_IDnone || qh->TRACEdist < REALmax/2) {
1719 qh_fprintf(qh, qh->ferr, 6051, "qhull option error: tracing is not installed (qh_NOtrace in user_r.h). Trace options 'Tn', 'TMn', 'TPn' and 'TWn' mostly removed. Continue with 'Qw' (allow warning)\n");
1720 if (!qh->ALLOWwarning)
1721 qh_errexit(qh, qh_ERRinput, NULL, NULL);
1722 }
1723 #endif
1724 if (qh->RERUN > 1) {
1725 qh->TRACElastrun= qh->IStracing; /* qh_build_withrestart duplicates next conditional */
1726 if (qh->IStracing && qh->IStracing != -1) {
1727 qh_fprintf(qh, qh->ferr, 8162, "qh_initqhull_globals: trace last of TR%d runs at level %d\n", qh->RERUN, qh->IStracing);
1728 qh->IStracing= 0;
1729 }
1730 }else if (qh->TRACEpoint != qh_IDnone || qh->TRACEdist < REALmax/2 || qh->TRACEmerge) {
1731 qh->TRACElevel= (qh->IStracing ? qh->IStracing : 3);
1732 qh->IStracing= 0;
1733 }
1734 if (qh->ROTATErandom == 0 || qh->ROTATErandom == -1) {
1735 seed= (int)time(&timedata);
1736 if (qh->ROTATErandom == -1) {
1737 seed= -seed;
1738 qh_option(qh, "QRandom-seed", &seed, NULL );
1739 }else
1740 qh_option(qh, "QRotate-random", &seed, NULL);
1741 qh->ROTATErandom= seed;
1742 }
1743 seed= qh->ROTATErandom;
1744 if (seed == INT_MIN) /* default value */
1745 seed= 1;
1746 else if (seed < 0)
1747 seed= -seed;
1748 qh_RANDOMseed_(qh, seed);
1749 randr= 0.0;
1750 for (i=1000; i--; ) {
1751 randi= qh_RANDOMint;
1752 randr += randi;
1753 if (randi > qh_RANDOMmax) {
1754 qh_fprintf(qh, qh->ferr, 8036, "\
1755 qhull configuration error (qh_RANDOMmax in user_r.h): random integer %d > qh_RANDOMmax (%.8g)\n",
1756 randi, qh_RANDOMmax);
1757 qh_errexit(qh, qh_ERRinput, NULL, NULL);
1758 }
1759 }
1760 qh_RANDOMseed_(qh, seed);
1761 randr= randr/1000;
1762 if (randr < qh_RANDOMmax * 0.1
1763 || randr > qh_RANDOMmax * 0.9)
1764 qh_fprintf(qh, qh->ferr, 8037, "\
1765 qhull configuration warning (qh_RANDOMmax in user_r.h): average of 1000 random integers (%.2g) is much different than expected (%.2g). Is qh_RANDOMmax (%.2g) wrong?\n",
1766 randr, qh_RANDOMmax * 0.5, qh_RANDOMmax);
1767 qh->RANDOMa= 2.0 * qh->RANDOMfactor/qh_RANDOMmax;
1768 qh->RANDOMb= 1.0 - qh->RANDOMfactor;
1769 if (qh_HASHfactor < 1.1) {
1770 qh_fprintf(qh, qh->ferr, 6052, "qhull internal error (qh_initqhull_globals): qh_HASHfactor %d must be at least 1.1. Qhull uses linear hash probing\n",
1771 qh_HASHfactor);
1772 qh_errexit(qh, qh_ERRqhull, NULL, NULL);
1773 }
1774 if (numpoints+extra < pointsneeded) {
1775 qh_fprintf(qh, qh->ferr, 6214, "qhull input error: not enough points(%d) to construct initial simplex (need %d)\n",
1776 numpoints, pointsneeded);
1777 qh_errexit(qh, qh_ERRinput, NULL, NULL);
1778 }
1779 qh_initqhull_outputflags(qh);
1780 } /* initqhull_globals */
1781
1782 /*-<a href="qh-globa_r.htm#TOC"
1783 >-------------------------------</a><a name="initqhull_mem">-</a>
1784
1785 qh_initqhull_mem(qh )
1786 initialize mem_r.c for qhull
1787 qh.hull_dim and qh.normal_size determine some of the allocation sizes
1788 if qh.MERGING,
1789 includes ridgeT
1790 calls qh_user_memsizes (user_r.c) to add up to 10 additional sizes for quick allocation
1791 (see numsizes below)
1792
1793 returns:
1794 mem_r.c already for qh_memalloc/qh_memfree (errors if called beforehand)
1795
1796 notes:
1797 qh_produceoutput() prints memsizes
1798
1799 */
qh_initqhull_mem(qhT * qh)1800 void qh_initqhull_mem(qhT *qh) {
1801 int numsizes;
1802 int i;
1803
1804 numsizes= 8+10;
1805 qh_meminitbuffers(qh, qh->IStracing, qh_MEMalign, numsizes,
1806 qh_MEMbufsize, qh_MEMinitbuf);
1807 qh_memsize(qh, (int)sizeof(vertexT));
1808 if (qh->MERGING) {
1809 qh_memsize(qh, (int)sizeof(ridgeT));
1810 qh_memsize(qh, (int)sizeof(mergeT));
1811 }
1812 qh_memsize(qh, (int)sizeof(facetT));
1813 i= (int)sizeof(setT) + (qh->hull_dim - 1) * SETelemsize; /* ridge.vertices */
1814 qh_memsize(qh, i);
1815 qh_memsize(qh, qh->normal_size); /* normal */
1816 i += SETelemsize; /* facet.vertices, .ridges, .neighbors */
1817 qh_memsize(qh, i);
1818 qh_user_memsizes(qh);
1819 qh_memsetup(qh);
1820 } /* initqhull_mem */
1821
1822 /*-<a href="qh-globa_r.htm#TOC"
1823 >-------------------------------</a><a name="initqhull_outputflags">-</a>
1824
1825 qh_initqhull_outputflags
1826 initialize flags concerned with output
1827
1828 returns:
1829 adjust user flags as needed
1830
1831 see:
1832 qh_clear_outputflags() resets the flags
1833
1834 design:
1835 test for qh.PRINTgood (i.e., only print 'good' facets)
1836 check for conflicting print output options
1837 */
qh_initqhull_outputflags(qhT * qh)1838 void qh_initqhull_outputflags(qhT *qh) {
1839 boolT printgeom= False, printmath= False, printcoplanar= False;
1840 int i;
1841
1842 trace3((qh, qh->ferr, 3024, "qh_initqhull_outputflags: %s\n", qh->qhull_command));
1843 if (!(qh->PRINTgood || qh->PRINTneighbors)) {
1844 if (qh->DELAUNAY || qh->KEEParea || qh->KEEPminArea < REALmax/2 || qh->KEEPmerge
1845 || (!qh->ONLYgood && (qh->GOODvertex || qh->GOODpoint))) {
1846 qh->PRINTgood= True;
1847 qh_option(qh, "Pgood", NULL, NULL);
1848 }
1849 }
1850 if (qh->PRINTtransparent) {
1851 if (qh->hull_dim != 4 || !qh->DELAUNAY || qh->VORONOI || qh->DROPdim >= 0) {
1852 qh_fprintf(qh, qh->ferr, 6215, "qhull option error: transparent Delaunay('Gt') needs 3-d Delaunay('d') w/o 'GDn'\n");
1853 qh_errexit(qh, qh_ERRinput, NULL, NULL);
1854 }
1855 qh->DROPdim= 3;
1856 qh->PRINTridges= True;
1857 }
1858 for (i=qh_PRINTEND; i--; ) {
1859 if (qh->PRINTout[i] == qh_PRINTgeom)
1860 printgeom= True;
1861 else if (qh->PRINTout[i] == qh_PRINTmathematica || qh->PRINTout[i] == qh_PRINTmaple)
1862 printmath= True;
1863 else if (qh->PRINTout[i] == qh_PRINTcoplanars)
1864 printcoplanar= True;
1865 else if (qh->PRINTout[i] == qh_PRINTpointnearest)
1866 printcoplanar= True;
1867 else if (qh->PRINTout[i] == qh_PRINTpointintersect && !qh->HALFspace) {
1868 qh_fprintf(qh, qh->ferr, 6053, "qhull option error: option 'Fp' is only used for \nhalfspace intersection('Hn,n,n').\n");
1869 qh_errexit(qh, qh_ERRinput, NULL, NULL);
1870 }else if (qh->PRINTout[i] == qh_PRINTtriangles && (qh->HALFspace || qh->VORONOI)) {
1871 qh_fprintf(qh, qh->ferr, 6054, "qhull option error: option 'Ft' is not available for Voronoi vertices ('v') or halfspace intersection ('H')\n");
1872 qh_errexit(qh, qh_ERRinput, NULL, NULL);
1873 }else if (qh->PRINTout[i] == qh_PRINTcentrums && qh->VORONOI) {
1874 qh_fprintf(qh, qh->ferr, 6055, "qhull option error: option 'FC' is not available for Voronoi vertices('v')\n");
1875 qh_errexit(qh, qh_ERRinput, NULL, NULL);
1876 }else if (qh->PRINTout[i] == qh_PRINTvertices) {
1877 if (qh->VORONOI)
1878 qh_option(qh, "Fvoronoi", NULL, NULL);
1879 else
1880 qh_option(qh, "Fvertices", NULL, NULL);
1881 }
1882 }
1883 if (printcoplanar && qh->DELAUNAY && qh->JOGGLEmax < REALmax/2) {
1884 if (qh->PRINTprecision)
1885 qh_fprintf(qh, qh->ferr, 7041, "qhull option warning: 'QJ' (joggle) will usually prevent coincident input sites for options 'Fc' and 'FP'\n");
1886 }
1887 if (printmath && (qh->hull_dim > 3 || qh->VORONOI)) {
1888 qh_fprintf(qh, qh->ferr, 6056, "qhull option error: Mathematica and Maple output is only available for 2-d and 3-d convex hulls and 2-d Delaunay triangulations\n");
1889 qh_errexit(qh, qh_ERRinput, NULL, NULL);
1890 }
1891 if (printgeom) {
1892 if (qh->hull_dim > 4) {
1893 qh_fprintf(qh, qh->ferr, 6057, "qhull option error: Geomview output is only available for 2-d, 3-d and 4-d\n");
1894 qh_errexit(qh, qh_ERRinput, NULL, NULL);
1895 }
1896 if (qh->PRINTnoplanes && !(qh->PRINTcoplanar + qh->PRINTcentrums
1897 + qh->PRINTdots + qh->PRINTspheres + qh->DOintersections + qh->PRINTridges)) {
1898 qh_fprintf(qh, qh->ferr, 6058, "qhull option error: no output specified for Geomview\n");
1899 qh_errexit(qh, qh_ERRinput, NULL, NULL);
1900 }
1901 if (qh->VORONOI && (qh->hull_dim > 3 || qh->DROPdim >= 0)) {
1902 qh_fprintf(qh, qh->ferr, 6059, "qhull option error: Geomview output for Voronoi diagrams only for 2-d\n");
1903 qh_errexit(qh, qh_ERRinput, NULL, NULL);
1904 }
1905 /* can not warn about furthest-site Geomview output: no lower_threshold */
1906 if (qh->hull_dim == 4 && qh->DROPdim == -1 &&
1907 (qh->PRINTcoplanar || qh->PRINTspheres || qh->PRINTcentrums)) {
1908 qh_fprintf(qh, qh->ferr, 7042, "qhull option warning: coplanars, vertices, and centrums output not available for 4-d output(ignored). Could use 'GDn' instead.\n");
1909 qh->PRINTcoplanar= qh->PRINTspheres= qh->PRINTcentrums= False;
1910 }
1911 }
1912 if (!qh->KEEPcoplanar && !qh->KEEPinside && !qh->ONLYgood) {
1913 if ((qh->PRINTcoplanar && qh->PRINTspheres) || printcoplanar) {
1914 if (qh->QHULLfinished) {
1915 qh_fprintf(qh, qh->ferr, 7072, "qhull output warning: ignoring coplanar points, option 'Qc' was not set for the first run of qhull.\n");
1916 }else {
1917 qh->KEEPcoplanar= True;
1918 qh_option(qh, "Qcoplanar", NULL, NULL);
1919 }
1920 }
1921 }
1922 qh->PRINTdim= qh->hull_dim;
1923 if (qh->DROPdim >=0) { /* after Geomview checks */
1924 if (qh->DROPdim < qh->hull_dim) {
1925 qh->PRINTdim--;
1926 if (!printgeom || qh->hull_dim < 3)
1927 qh_fprintf(qh, qh->ferr, 7043, "qhull option warning: drop dimension 'GD%d' is only available for 3-d/4-d Geomview\n", qh->DROPdim);
1928 }else
1929 qh->DROPdim= -1;
1930 }else if (qh->VORONOI) {
1931 qh->DROPdim= qh->hull_dim-1;
1932 qh->PRINTdim= qh->hull_dim-1;
1933 }
1934 } /* qh_initqhull_outputflags */
1935
1936 /*-<a href="qh-globa_r.htm#TOC"
1937 >-------------------------------</a><a name="initqhull_start">-</a>
1938
1939 qh_initqhull_start(qh, infile, outfile, errfile )
1940 allocate memory if needed and call qh_initqhull_start2()
1941 */
qh_initqhull_start(qhT * qh,FILE * infile,FILE * outfile,FILE * errfile)1942 void qh_initqhull_start(qhT *qh, FILE *infile, FILE *outfile, FILE *errfile) {
1943
1944 qh_initstatistics(qh);
1945 qh_initqhull_start2(qh, infile, outfile, errfile);
1946 } /* initqhull_start */
1947
1948 /*-<a href="qh-globa_r.htm#TOC"
1949 >-------------------------------</a><a name="initqhull_start2">-</a>
1950
1951 qh_initqhull_start2(qh, infile, outfile, errfile )
1952 start initialization of qhull
1953 initialize statistics, stdio, default values for global variables
1954 assumes qh is allocated
1955 notes:
1956 report errors elsewhere, error handling and g_qhull_output [Qhull.cpp, QhullQh()] not in initialized
1957 see:
1958 qh_maxmin() determines the precision constants
1959 qh_freeqhull()
1960 */
qh_initqhull_start2(qhT * qh,FILE * infile,FILE * outfile,FILE * errfile)1961 void qh_initqhull_start2(qhT *qh, FILE *infile, FILE *outfile, FILE *errfile) {
1962 time_t timedata;
1963 int seed;
1964
1965 qh_CPUclock; /* start the clock(for qh_clock). One-shot. */
1966 /* memset is the same in qh_freeqhull() and qh_initqhull_start2() */
1967 memset((char *)qh, 0, sizeof(qhT)-sizeof(qhmemT)-sizeof(qhstatT)); /* every field is 0, FALSE, NULL */
1968 qh->NOerrexit= True;
1969 qh->DROPdim= -1;
1970 qh->ferr= errfile;
1971 qh->fin= infile;
1972 qh->fout= outfile;
1973 qh->furthest_id= qh_IDunknown;
1974 #ifndef qh_NOmerge
1975 qh->JOGGLEmax= REALmax;
1976 #else
1977 qh->JOGGLEmax= 0.0; /* Joggle ('QJ') if qh_NOmerge */
1978 #endif
1979 qh->KEEPminArea= REALmax;
1980 qh->last_low= REALmax;
1981 qh->last_high= REALmax;
1982 qh->last_newhigh= REALmax;
1983 qh->last_random= 1; /* reentrant only */
1984 qh->lastcpu= 0.0;
1985 qh->max_outside= 0.0;
1986 qh->max_vertex= 0.0;
1987 qh->MAXabs_coord= 0.0;
1988 qh->MAXsumcoord= 0.0;
1989 qh->MAXwidth= -REALmax;
1990 qh->MERGEindependent= True;
1991 qh->MINdenom_1= fmax_(1.0/REALmax, REALmin); /* used by qh_scalepoints */
1992 qh->MINoutside= 0.0;
1993 qh->MINvisible= REALmax;
1994 qh->MAXcoplanar= REALmax;
1995 qh->outside_err= REALmax;
1996 qh->premerge_centrum= 0.0;
1997 qh->premerge_cos= REALmax;
1998 qh->PRINTprecision= True;
1999 qh->PRINTradius= 0.0;
2000 qh->postmerge_cos= REALmax;
2001 qh->postmerge_centrum= 0.0;
2002 qh->ROTATErandom= INT_MIN;
2003 qh->MERGEvertices= True;
2004 qh->totarea= 0.0;
2005 qh->totvol= 0.0;
2006 qh->TRACEdist= REALmax;
2007 qh->TRACEpoint= qh_IDnone; /* recompile to trace a point, or use 'TPn' */
2008 qh->tracefacet_id= UINT_MAX; /* recompile to trace a facet, set to UINT_MAX when done, see userprintf_r.c/qh_fprintf */
2009 qh->traceridge_id= UINT_MAX; /* recompile to trace a ridge, set to UINT_MAX when done, see userprintf_r.c/qh_fprintf */
2010 qh->tracevertex_id= UINT_MAX; /* recompile to trace a vertex, set to UINT_MAX when done, see userprintf_r.c/qh_fprintf */
2011 seed= (int)time(&timedata);
2012 qh_RANDOMseed_(qh, seed);
2013 qh->run_id= qh_RANDOMint;
2014 if(!qh->run_id)
2015 qh->run_id++; /* guarantee non-zero */
2016 qh_option(qh, "run-id", &qh->run_id, NULL);
2017 strcat(qh->qhull, "qhull");
2018 } /* initqhull_start2 */
2019
2020 /*-<a href="qh-globa_r.htm#TOC"
2021 >-------------------------------</a><a name="initthresholds">-</a>
2022
2023 qh_initthresholds(qh, commandString )
2024 set thresholds for printing and scaling from commandString
2025
2026 returns:
2027 sets qh.GOODthreshold or qh.SPLITthreshold if 'Pd0D1' used
2028
2029 see:
2030 qh_initflags(), 'Qbk' 'QBk' 'Pdk' and 'PDk'
2031 qh_inthresholds()
2032
2033 design:
2034 for each 'Pdn' or 'PDn' option
2035 check syntax
2036 set qh.lower_threshold or qh.upper_threshold
2037 set qh.GOODthreshold if an unbounded threshold is used
2038 set qh.SPLITthreshold if a bounded threshold is used
2039 */
qh_initthresholds(qhT * qh,char * command)2040 void qh_initthresholds(qhT *qh, char *command) {
2041 realT value;
2042 int idx, maxdim, k;
2043 char *s= command; /* non-const due to strtol */
2044 char *lastoption, *lastwarning= NULL;
2045 char key;
2046
2047 maxdim= qh->input_dim;
2048 if (qh->DELAUNAY && (qh->PROJECTdelaunay || qh->PROJECTinput))
2049 maxdim++;
2050 while (*s) {
2051 if (*s == '-')
2052 s++;
2053 if (*s == 'P') {
2054 lastoption= s++;
2055 while (*s && !isspace(key= *s++)) {
2056 if (key == 'd' || key == 'D') {
2057 if (!isdigit(*s)) {
2058 qh_fprintf(qh, qh->ferr, 7044, "qhull option warning: no dimension given for Print option 'P%c' at: %s. Ignored\n",
2059 key, s-1);
2060 lastwarning= lastoption;
2061 continue;
2062 }
2063 idx= qh_strtol(s, &s);
2064 if (idx >= qh->hull_dim) {
2065 qh_fprintf(qh, qh->ferr, 7045, "qhull option warning: dimension %d for Print option 'P%c' is >= %d. Ignored\n",
2066 idx, key, qh->hull_dim);
2067 lastwarning= lastoption;
2068 continue;
2069 }
2070 if (*s == ':') {
2071 s++;
2072 value= qh_strtod(s, &s);
2073 if (fabs((double)value) > 1.0) {
2074 qh_fprintf(qh, qh->ferr, 7046, "qhull option warning: value %2.4g for Print option 'P%c' is > +1 or < -1. Ignored\n",
2075 value, key);
2076 lastwarning= lastoption;
2077 continue;
2078 }
2079 }else
2080 value= 0.0;
2081 if (key == 'd')
2082 qh->lower_threshold[idx]= value;
2083 else
2084 qh->upper_threshold[idx]= value;
2085 }
2086 }
2087 }else if (*s == 'Q') {
2088 lastoption= s++;
2089 while (*s && !isspace(key= *s++)) {
2090 if (key == 'b' && *s == 'B') {
2091 s++;
2092 for (k=maxdim; k--; ) {
2093 qh->lower_bound[k]= -qh_DEFAULTbox;
2094 qh->upper_bound[k]= qh_DEFAULTbox;
2095 }
2096 }else if (key == 'b' && *s == 'b')
2097 s++;
2098 else if (key == 'b' || key == 'B') {
2099 if (!isdigit(*s)) {
2100 qh_fprintf(qh, qh->ferr, 7047, "qhull option warning: no dimension given for Qhull option 'Q%c'\n",
2101 key);
2102 lastwarning= lastoption;
2103 continue;
2104 }
2105 idx= qh_strtol(s, &s);
2106 if (idx >= maxdim) {
2107 qh_fprintf(qh, qh->ferr, 7048, "qhull option warning: dimension %d for Qhull option 'Q%c' is >= %d. Ignored\n",
2108 idx, key, maxdim);
2109 lastwarning= lastoption;
2110 continue;
2111 }
2112 if (*s == ':') {
2113 s++;
2114 value= qh_strtod(s, &s);
2115 }else if (key == 'b')
2116 value= -qh_DEFAULTbox;
2117 else
2118 value= qh_DEFAULTbox;
2119 if (key == 'b')
2120 qh->lower_bound[idx]= value;
2121 else
2122 qh->upper_bound[idx]= value;
2123 }
2124 }
2125 }else {
2126 while (*s && !isspace(*s))
2127 s++;
2128 }
2129 while (isspace(*s))
2130 s++;
2131 }
2132 for (k=qh->hull_dim; k--; ) {
2133 if (qh->lower_threshold[k] > -REALmax/2) {
2134 qh->GOODthreshold= True;
2135 if (qh->upper_threshold[k] < REALmax/2) {
2136 qh->SPLITthresholds= True;
2137 qh->GOODthreshold= False;
2138 break;
2139 }
2140 }else if (qh->upper_threshold[k] < REALmax/2)
2141 qh->GOODthreshold= True;
2142 }
2143 if (lastwarning && !qh->ALLOWwarning) {
2144 qh_fprintf(qh, qh->ferr, 6036, "qhull option error: see previous warnings, use 'Qw' to override: '%s' (last offset %d)\n",
2145 command, (int)(lastwarning-command));
2146 qh_errexit(qh, qh_ERRinput, NULL, NULL);
2147 }
2148 } /* initthresholds */
2149
2150 /*-<a href="qh-globa_r.htm#TOC"
2151 >-------------------------------</a><a name="lib_check">-</a>
2152
2153 qh_lib_check( qhullLibraryType, qhTsize, vertexTsize, ridgeTsize, facetTsize, setTsize, qhmemTsize )
2154 Report error if library does not agree with caller
2155
2156 notes:
2157 NOerrors -- qh_lib_check can not call qh_errexit()
2158 */
qh_lib_check(int qhullLibraryType,int qhTsize,int vertexTsize,int ridgeTsize,int facetTsize,int setTsize,int qhmemTsize)2159 void qh_lib_check(int qhullLibraryType, int qhTsize, int vertexTsize, int ridgeTsize, int facetTsize, int setTsize, int qhmemTsize) {
2160 int last_errcode= qh_ERRnone;
2161
2162 #if defined(_MSC_VER) && defined(_DEBUG) && defined(QHULL_CRTDBG) /* user_r.h */
2163 /*_CrtSetBreakAlloc(744);*/ /* Break at memalloc {744}, or 'watch' _crtBreakAlloc */
2164 _CrtSetDbgFlag( _CRTDBG_ALLOC_MEM_DF | _CRTDBG_DELAY_FREE_MEM_DF | _CRTDBG_LEAK_CHECK_DF | _CrtSetDbgFlag(_CRTDBG_REPORT_FLAG) );
2165 _CrtSetReportMode( _CRT_ERROR, _CRTDBG_MODE_FILE | _CRTDBG_MODE_DEBUG );
2166 _CrtSetReportFile( _CRT_ERROR, _CRTDBG_FILE_STDERR );
2167 _CrtSetReportMode( _CRT_WARN, _CRTDBG_MODE_FILE | _CRTDBG_MODE_DEBUG );
2168 _CrtSetReportFile( _CRT_WARN, _CRTDBG_FILE_STDERR );
2169 _CrtSetReportMode( _CRT_ASSERT, _CRTDBG_MODE_FILE | _CRTDBG_MODE_DEBUG );
2170 _CrtSetReportFile( _CRT_ASSERT, _CRTDBG_FILE_STDERR );
2171 #endif
2172
2173 if (qhullLibraryType==QHULL_NON_REENTRANT) { /* 0 */
2174 qh_fprintf_stderr(6257, "qh_lib_check: Incorrect qhull library called. Caller uses non-reentrant Qhull with a static qhT. Qhull library is reentrant.\n");
2175 last_errcode= 6257;
2176 }else if (qhullLibraryType==QHULL_QH_POINTER) { /* 1 */
2177 qh_fprintf_stderr(6258, "qh_lib_check: Incorrect qhull library called. Caller uses non-reentrant Qhull with a dynamic qhT via qh_QHpointer. Qhull library is reentrant.\n");
2178 last_errcode= 6258;
2179 }else if (qhullLibraryType != QHULL_REENTRANT) { /* 2 */
2180 qh_fprintf_stderr(6262, "qh_lib_check: Expecting qhullLibraryType QHULL_NON_REENTRANT(0), QHULL_QH_POINTER(1), or QHULL_REENTRANT(2). Got %d\n", qhullLibraryType);
2181 last_errcode= 6262;
2182 }
2183 if (qhTsize != (int)sizeof(qhT)) {
2184 qh_fprintf_stderr(6249, "qh_lib_check: Incorrect qhull library called. Size of qhT for caller is %d, but for qhull library is %d.\n", qhTsize, (int)sizeof(qhT));
2185 last_errcode= 6249;
2186 }
2187 if (vertexTsize != (int)sizeof(vertexT)) {
2188 qh_fprintf_stderr(6250, "qh_lib_check: Incorrect qhull library called. Size of vertexT for caller is %d, but for qhull library is %d.\n", vertexTsize, (int)sizeof(vertexT));
2189 last_errcode= 6250;
2190 }
2191 if (ridgeTsize != (int)sizeof(ridgeT)) {
2192 qh_fprintf_stderr(6251, "qh_lib_check: Incorrect qhull library called. Size of ridgeT for caller is %d, but for qhull library is %d.\n", ridgeTsize, (int)sizeof(ridgeT));
2193 last_errcode= 6251;
2194 }
2195 if (facetTsize != (int)sizeof(facetT)) {
2196 qh_fprintf_stderr(6252, "qh_lib_check: Incorrect qhull library called. Size of facetT for caller is %d, but for qhull library is %d.\n", facetTsize, (int)sizeof(facetT));
2197 last_errcode= 6252;
2198 }
2199 if (setTsize && setTsize != (int)sizeof(setT)) {
2200 qh_fprintf_stderr(6253, "qh_lib_check: Incorrect qhull library called. Size of setT for caller is %d, but for qhull library is %d.\n", setTsize, (int)sizeof(setT));
2201 last_errcode= 6253;
2202 }
2203 if (qhmemTsize && qhmemTsize != sizeof(qhmemT)) {
2204 qh_fprintf_stderr(6254, "qh_lib_check: Incorrect qhull library called. Size of qhmemT for caller is %d, but for qhull library is %d.\n", qhmemTsize, sizeof(qhmemT));
2205 last_errcode= 6254;
2206 }
2207 if (last_errcode) {
2208 qh_fprintf_stderr(6259, "qhull internal error (qh_lib_check): Cannot continue due to QH%d. '%s' is not reentrant (e.g., qhull.so) or out-of-date. Exit with %d\n",
2209 last_errcode, qh_version2, last_errcode - 6200);
2210 qh_exit(last_errcode - 6200); /* can not use qh_errexit(), must be less than 255 */
2211 }
2212 } /* lib_check */
2213
2214 /*-<a href="qh-globa_r.htm#TOC"
2215 >-------------------------------</a><a name="option">-</a>
2216
2217 qh_option(qh, option, intVal, realVal )
2218 add an option description to qh.qhull_options
2219
2220 notes:
2221 NOerrors -- qh_option can not call qh_errexit() [qh_initqhull_start2]
2222 will be printed with statistics ('Ts') and errors
2223 strlen(option) < 40
2224 */
qh_option(qhT * qh,const char * option,int * i,realT * r)2225 void qh_option(qhT *qh, const char *option, int *i, realT *r) {
2226 char buf[200];
2227 int buflen, remainder;
2228
2229 if (strlen(option) > sizeof(buf)-30-30) {
2230 qh_fprintf(qh, qh->ferr, 6408, "qhull internal error (qh_option): option (%d chars) has more than %d chars. May overflow temporary buffer. Option '%s'\n",
2231 (int)strlen(option), (int)sizeof(buf)-30-30, option);
2232 qh_errexit(qh, qh_ERRqhull, NULL, NULL);
2233 }
2234 sprintf(buf, " %s", option);
2235 if (i)
2236 sprintf(buf+strlen(buf), " %d", *i);
2237 if (r)
2238 sprintf(buf+strlen(buf), " %2.2g", *r);
2239 buflen= (int)strlen(buf); /* WARN64 */
2240 qh->qhull_optionlen += buflen;
2241 remainder= (int)(sizeof(qh->qhull_options) - strlen(qh->qhull_options)) - 1; /* WARN64 */
2242 maximize_(remainder, 0);
2243 if (qh->qhull_optionlen >= qh_OPTIONline && remainder > 0) {
2244 strncat(qh->qhull_options, "\n", (unsigned int)remainder);
2245 --remainder;
2246 qh->qhull_optionlen= buflen;
2247 }
2248 if (buflen > remainder) {
2249 trace1((qh, qh->ferr, 1058, "qh_option: option would overflow qh.qhull_options. Truncated '%s'\n", buf));
2250 }
2251 strncat(qh->qhull_options, buf, (unsigned int)remainder);
2252 } /* option */
2253
2254 /*-<a href="qh-globa_r.htm#TOC"
2255 >-------------------------------</a><a name="zero">-</a>
2256
2257 qh_zero( qh, errfile )
2258 Initialize and zero Qhull's memory for qh_new_qhull()
2259
2260 notes:
2261 Not needed in global_r.c because static variables are initialized to zero
2262 */
qh_zero(qhT * qh,FILE * errfile)2263 void qh_zero(qhT *qh, FILE *errfile) {
2264 memset((char *)qh, 0, sizeof(qhT)); /* every field is 0, FALSE, NULL */
2265 qh->NOerrexit= True;
2266 qh_meminit(qh, errfile);
2267 } /* zero */
2268
2269