user_eg2_r.c
Go to the documentation of this file.
1 /*<html><pre> -<a href="../libqhull_r/qh-user_r.htm"
2  >-------------------------------</a><a name="TOP">-</a>
3 
4  user_eg2_r.c
5 
6  sample code for calling qhull() from an application.
7 
8  See user_eg_r.c for a simpler method using qh_new_qhull().
9  The method used here and in unix_r.c gives you additional
10  control over Qhull.
11 
12  See user_eg3/user_eg3_r.cpp for a C++ example
13 
14  call with:
15 
16  user_eg2 "triangulated cube/diamond options" "delaunay options" "halfspace options"
17 
18  for example:
19 
20  user_eg2 # return summaries
21 
22  user_eg2 "n" "o" "Fp" # return normals, OFF, points
23 
24  user_eg2 "QR0 p" "QR0 v p" "QR0 Fp" # rotate input and return points
25  # 'v' returns Voronoi
26  # transform is rotated for halfspaces
27 
28  main() makes three runs of qhull.
29 
30  1) compute the convex hull of a cube, and incrementally add a diamond
31 
32  2a) compute the Delaunay triangulation of random points, and add points.
33 
34  2b) find the Delaunay triangle closest to a point.
35 
36  3) compute the halfspace intersection of a diamond, and add a cube
37 
38  notes:
39 
40  summaries are sent to stderr if other output formats are used
41 
42  derived from unix.c and compiled by 'make bin/user_eg2'
43 
44  see libqhull.h for data structures, macros, and user-callable functions.
45 
46  If you want to control all output to stdio and input to stdin,
47  set the #if below to "1" and delete all lines that contain "io.c".
48  This prevents the loading of io.o. Qhull will
49  still write to 'qh->ferr' (stderr) for error reporting and tracing.
50 
51  Defining #if 1, also prevents user.o from being loaded.
52 */
53 
54 #include "libqhull_r/qhull_ra.h"
55 
56 /*-------------------------------------------------
57 -internal function prototypes
58 */
59 void print_summary(qhT *qh);
60 void makecube(coordT *points, int numpoints, int dim);
61 void adddiamond(qhT *qh, coordT *points, int numpoints, int numnew, int dim);
62 void makeDelaunay(qhT *qh, coordT *points, int numpoints, int dim);
63 void addDelaunay(qhT *qh, coordT *points, int numpoints, int numnew, int dim);
64 void findDelaunay(qhT *qh, int dim);
65 void makehalf(coordT *points, int numpoints, int dim);
66 void addhalf(qhT *qh, coordT *points, int numpoints, int numnew, int dim, coordT *feasible);
67 
68 /*-------------------------------------------------
69 -print_summary(qh)
70 */
72  facetT *facet;
73  int k;
74 
75  printf("\n%d vertices and %d facets with normals:\n",
76  qh->num_vertices, qh->num_facets);
77  FORALLfacets {
78  for (k=0; k < qh->hull_dim; k++)
79  printf("%6.2g ", facet->normal[k]);
80  printf("\n");
81  }
82 }
83 
84 /*--------------------------------------------------
85 -makecube- set points to vertices of cube
86  points is numpoints X dim
87 */
88 void makecube(coordT *points, int numpoints, int dim) {
89  int j,k;
90  coordT *point;
91 
92  for (j=0; j<numpoints; j++) {
93  point= points + j*dim;
94  for (k=dim; k--; ) {
95  if (j & ( 1 << k))
96  point[k]= 1.0;
97  else
98  point[k]= -1.0;
99  }
100  }
101 } /*.makecube.*/
102 
103 /*--------------------------------------------------
104 -adddiamond- add diamond to convex hull
105  points is numpoints+numnew X dim.
106 
107 notes:
108  qh_addpoint() does not make a copy of the point coordinates.
109 
110  For inside points and some outside points, qh_findbestfacet performs
111  an exhaustive search for a visible facet. Algorithms that retain
112  previously constructed hulls should be faster for on-line construction
113  of the convex hull.
114 */
115 void adddiamond(qhT *qh, coordT *points, int numpoints, int numnew, int dim) {
116  int j,k;
117  coordT *point;
118  facetT *facet;
119  boolT isoutside;
120  realT bestdist;
121 
122  for (j= 0; j < numnew ; j++) {
123  point= points + (numpoints+j)*dim;
124  if (points == qh->first_point) /* in case of 'QRn' */
125  qh->num_points= numpoints+j+1;
126  /* qh.num_points sets the size of the points array. You may
127  allocate the points elsewhere. If so, qh_addpoint records
128  the point's address in qh->other_points
129  */
130  for (k=dim; k--; ) {
131  if (j/2 == k)
132  point[k]= (j & 1) ? 2.0 : -2.0;
133  else
134  point[k]= 0.0;
135  }
136  facet= qh_findbestfacet(qh, point, !qh_ALL, &bestdist, &isoutside);
137  if (isoutside) {
138  if (!qh_addpoint(qh, point, facet, False))
139  break; /* user requested an early exit with 'TVn' or 'TCn' */
140  }
141  printf("%d vertices and %d facets\n",
142  qh->num_vertices, qh->num_facets);
143  /* qh_produce_output(); */
144  }
145  if (qh->DOcheckmax)
146  qh_check_maxout(qh);
147  else if (qh->KEEPnearinside)
148  qh_nearcoplanar(qh);
149 } /*.adddiamond.*/
150 
151 /*--------------------------------------------------
152 -makeDelaunay- set points for dim-1 Delaunay triangulation of random points
153  points is numpoints X dim. Each point is projected to a paraboloid.
154 */
155 void makeDelaunay(qhT *qh, coordT *points, int numpoints, int dim) {
156  int j,k, seed;
157  coordT *point, realr;
158 
159  seed= (int)time(NULL); /* time_t to int */
160  printf("seed: %d\n", seed);
161  qh_RANDOMseed_(qh, seed);
162  for (j=0; j<numpoints; j++) {
163  point= points + j*dim;
164  for (k= 0; k < dim-1; k++) {
165  realr= qh_RANDOMint;
166  point[k]= 2.0 * realr/(qh_RANDOMmax+1) - 1.0;
167  }
168  }
169  qh_setdelaunay(qh, dim, numpoints, points);
170 } /*.makeDelaunay.*/
171 
172 /*--------------------------------------------------
173 -addDelaunay- add points to dim-1 Delaunay triangulation
174  points is numpoints+numnew X dim. Each point is projected to a paraboloid.
175 notes:
176  qh_addpoint() does not make a copy of the point coordinates.
177 
178  Since qh_addpoint() is not given a visible facet, it performs a directed
179  search of all facets. Algorithms that retain previously
180  constructed hulls may be faster.
181 */
182 void addDelaunay(qhT *qh, coordT *points, int numpoints, int numnew, int dim) {
183  int j,k;
184  coordT *point, realr;
185  facetT *facet;
186  realT bestdist;
187  boolT isoutside;
188 
189  for (j= 0; j < numnew ; j++) {
190  point= points + (numpoints+j)*dim;
191  if (points == qh->first_point) /* in case of 'QRn' */
192  qh->num_points= numpoints+j+1;
193  /* qh.num_points sets the size of the points array. You may
194  allocate the point elsewhere. If so, qh_addpoint records
195  the point's address in qh->other_points
196  */
197  for (k= 0; k < dim-1; k++) {
198  realr= qh_RANDOMint;
199  point[k]= 2.0 * realr/(qh_RANDOMmax+1) - 1.0;
200  }
201  qh_setdelaunay(qh, dim, 1, point);
202  facet= qh_findbestfacet(qh, point, !qh_ALL, &bestdist, &isoutside);
203  if (isoutside) {
204  if (!qh_addpoint(qh, point, facet, False))
205  break; /* user requested an early exit with 'TVn' or 'TCn' */
206  }
207  qh_printpoint(qh, stdout, "added point", point);
208  printf("%d points, %d extra points, %d vertices, and %d facets in total\n",
209  qh->num_points, qh_setsize(qh, qh->other_points),
210  qh->num_vertices, qh->num_facets);
211 
212  /* qh_produce_output(qh); */
213  }
214  if (qh->DOcheckmax)
215  qh_check_maxout(qh);
216  else if (qh->KEEPnearinside)
217  qh_nearcoplanar(qh);
218 } /*.addDelaunay.*/
219 
220 /*--------------------------------------------------
221 -findDelaunay- find Delaunay triangle for [0.5,0.5,...]
222  assumes dim < 100
223 notes:
224  calls qh_setdelaunay() to project the point to a parabaloid
225 warning:
226  This is not implemented for tricoplanar facets ('Qt'),
227  See <a href="../html/qh-code.htm#findfacet">locate a facet with qh_findbestfacet()</a>
228 */
229 void findDelaunay(qhT *qh, int dim) {
230  int k;
231  coordT point[ 100];
232  boolT isoutside;
233  realT bestdist;
234  facetT *facet;
235  vertexT *vertex, **vertexp;
236 
237  for (k= 0; k < dim-1; k++)
238  point[k]= 0.5;
239  qh_setdelaunay(qh, dim, 1, point);
240  facet= qh_findbestfacet(qh, point, qh_ALL, &bestdist, &isoutside);
241  if (facet->tricoplanar) {
242  fprintf(stderr, "findDelaunay: not implemented for triangulated, non-simplicial Delaunay regions (tricoplanar facet, f%d).\n",
243  facet->id);
244  qh_errexit(qh, qh_ERRqhull, facet, NULL);
245  }
246  FOREACHvertex_(facet->vertices) {
247  for (k=0; k < dim-1; k++)
248  printf("%5.2f ", vertex->point[k]);
249  printf("\n");
250  }
251 } /*.findDelaunay.*/
252 
253 /*--------------------------------------------------
254 -makehalf- set points to halfspaces for a (dim)-d diamond
255  points is numpoints X dim+1
256 
257  each halfspace consists of dim coefficients followed by an offset
258 */
259 void makehalf(coordT *points, int numpoints, int dim) {
260  int j,k;
261  coordT *point;
262 
263  for (j=0; j<numpoints; j++) {
264  point= points + j*(dim+1);
265  point[dim]= -1.0; /* offset */
266  for (k=dim; k--; ) {
267  if (j & ( 1 << k))
268  point[k]= 1.0;
269  else
270  point[k]= -1.0;
271  }
272  }
273 } /*.makehalf.*/
274 
275 /*--------------------------------------------------
276 -addhalf- add halfspaces for a (dim)-d cube to the intersection
277  points is numpoints+numnew X dim+1
278 notes:
279  assumes dim < 100.
280 
281  For makehalf(), points is the initial set of halfspaces with offsets.
282  It is transformed by qh_sethalfspace_all into a
283  (dim)-d set of newpoints. Qhull computed the convex hull of newpoints -
284  this is equivalent to the halfspace intersection of the
285  orginal halfspaces.
286 
287  For addhalf(), the remainder of points stores the transforms of
288  the added halfspaces. Qhull computes the convex hull of newpoints
289  and the added points. qh_addpoint() does not make a copy of these points.
290 
291  Since halfspace intersection is equivalent to a convex hull,
292  qh_findbestfacet may perform an exhaustive search
293  for a visible facet. Algorithms that retain previously constructed
294  intersections should be faster for on-line construction.
295 */
296 void addhalf(qhT *qh, coordT *points, int numpoints, int numnew, int dim, coordT *feasible) {
297  int j,k;
298  coordT *point, normal[100], offset, *next;
299  facetT *facet;
300  boolT isoutside;
301  realT bestdist;
302 
303  for (j= 0; j < numnew ; j++) {
304  offset= -1.0;
305  for (k=dim; k--; ) {
306  if (j/2 == k) {
307  normal[k]= sqrt((coordT)dim); /* to normalize as in makehalf */
308  if (j & 1)
309  normal[k]= -normal[k];
310  }else
311  normal[k]= 0.0;
312  }
313  point= points + (numpoints+j)* (dim+1); /* does not use point[dim] */
314  qh_sethalfspace(qh, dim, point, &next, normal, &offset, feasible);
315  facet= qh_findbestfacet(qh, point, !qh_ALL, &bestdist, &isoutside);
316  if (isoutside) {
317  if (!qh_addpoint(qh, point, facet, False))
318  break; /* user requested an early exit with 'TVn' or 'TCn' */
319  }
320  qh_printpoint(qh, stdout, "added offset -1 and normal", normal);
321  printf("%d points, %d extra points, %d vertices, and %d facets in total\n",
322  qh->num_points, qh_setsize(qh, qh->other_points),
323  qh->num_vertices, qh->num_facets);
324  /* qh_produce_output(qh); */
325  }
326  if (qh->DOcheckmax)
327  qh_check_maxout(qh);
328  else if (qh->KEEPnearinside)
329  qh_nearcoplanar(qh);
330 } /*.addhalf.*/
331 
332 #define DIM 3 /* dimension of points, must be < 31 for SIZEcube */
333 #define SIZEcube (1<<DIM)
334 #define SIZEdiamond (2*DIM)
335 #define TOTpoints (SIZEcube + SIZEdiamond)
336 
337 /*--------------------------------------------------
338 -main- Similar to unix.c, the main program for qhull
339 
340  see program header
341 
342  this contains three runs of Qhull for convex hull, Delaunay
343  triangulation or Voronoi vertices, and halfspace intersection
344 
345 */
346 int main(int argc, char *argv[]) {
347  boolT ismalloc;
348  int curlong, totlong, exitcode;
349  char options [2000];
350  qhT qh_qh;
351  qhT *qh= &qh_qh; /* Alternatively -- qhT *qh= (qhT*)malloc(sizeof(qhT)) */
352 
354 
355  printf("This is the output from user_eg2_r.c\n\n\
356 It shows how qhull() may be called from an application using qhull's\n\
357 static, reentrant library. user_eg2 is not part of qhull itself. If it\n\
358 appears accidently, please remove user_eg2_r.c from your project. If it fails\n\
359 immediately, user_eg2_r.c was incorrectly linked to the non-reentrant library.\n\
360 Also try 'user_eg2 T1 2>&1'\n\n");
361 
362  ismalloc= False; /* True if qh_freeqhull should 'free(array)' */
363  /*
364  Run 1: convex hull
365  */
366  qh_init_A(qh, stdin, stdout, stderr, 0, NULL);
367  exitcode= setjmp(qh->errexit);
368  if (!exitcode) {
369  coordT array[TOTpoints][DIM];
370 
371  qh->NOerrexit= False;
372  strcat(qh->rbox_command, "user_eg2 cube example");
373  sprintf(options, "qhull s Tcv Q11 %s ", argc >= 2 ? argv[1] : "");
374  qh_initflags(qh, options);
375  printf( "\ncompute triangulated convex hull of cube after rotating input\n");
376  makecube(array[0], SIZEcube, DIM);
377  qh_init_B(qh, array[0], SIZEcube, DIM, ismalloc);
378  qh_qhull(qh);
379  qh_check_output(qh);
380  qh_triangulate(qh); /* requires option 'Q11' if want to add points */
381  print_summary(qh);
382  if (qh->VERIFYoutput && !qh->STOPpoint && !qh->STOPcone)
383  qh_check_points(qh);
384  printf( "\nadd points in a diamond\n");
385  adddiamond(qh, array[0], SIZEcube, SIZEdiamond, DIM);
386  qh_check_output(qh);
387  print_summary(qh);
388  qh_produce_output(qh); /* delete this line to help avoid io.c */
389  if (qh->VERIFYoutput && !qh->STOPpoint && !qh->STOPcone)
390  qh_check_points(qh);
391  }
392  qh->NOerrexit= True;
393  qh_freeqhull(qh, !qh_ALL);
394  qh_memfreeshort(qh, &curlong, &totlong);
395  if (curlong || totlong)
396  fprintf(stderr, "qhull warning (user_eg2, run 1): did not free %d bytes of long memory (%d pieces)\n",
397  totlong, curlong);
398  /*
399  Run 2: Delaunay triangulation
400  */
401  qh_init_A(qh, stdin, stdout, stderr, 0, NULL);
402  exitcode= setjmp(qh->errexit);
403  if (!exitcode) {
404  coordT array[TOTpoints][DIM];
405 
406  qh->NOerrexit= False;
407  strcat(qh->rbox_command, "user_eg2 Delaunay example");
408  sprintf(options, "qhull s d Tcv %s", argc >= 3 ? argv[2] : "");
409  qh_initflags(qh, options);
410  printf( "\ncompute %d-d Delaunay triangulation\n", DIM-1);
411  makeDelaunay(qh, array[0], SIZEcube, DIM);
412  /* Instead of makeDelaunay with qh_setdelaunay, you may
413  produce a 2-d array of points, set DIM to 2, and set
414  qh->PROJECTdelaunay to True. qh_init_B will call
415  qh_projectinput to project the points to the paraboloid
416  and add a point "at-infinity".
417  */
418  qh_init_B(qh, array[0], SIZEcube, DIM, ismalloc);
419  qh_qhull(qh);
420  /* If you want Voronoi ('v') without qh_produce_output(), call
421  qh_setvoronoi_all() after qh_qhull() */
422  qh_check_output(qh);
423  print_summary(qh);
424  qh_produce_output(qh); /* delete this line to help avoid io.c */
425  if (qh->VERIFYoutput && !qh->STOPpoint && !qh->STOPcone)
426  qh_check_points(qh);
427  printf( "\nadd points to triangulation\n");
428  addDelaunay(qh, array[0], SIZEcube, SIZEdiamond, DIM);
429  qh_check_output(qh);
430  qh_produce_output(qh); /* delete this line to help avoid io.c */
431  if (qh->VERIFYoutput && !qh->STOPpoint && !qh->STOPcone)
432  qh_check_points(qh);
433  printf( "\nfind Delaunay triangle closest to [0.5, 0.5, ...]\n");
434  findDelaunay(qh, DIM);
435  }
436  qh->NOerrexit= True;
437  qh_freeqhull(qh, !qh_ALL);
438  qh_memfreeshort(qh, &curlong, &totlong);
439  if (curlong || totlong)
440  fprintf(stderr, "qhull warning (user_eg2, run 2): did not free %d bytes of long memory (%d pieces)\n",
441  totlong, curlong);
442  /*
443  Run 3: halfspace intersection
444  */
445  qh_init_A(qh, stdin, stdout, stderr, 0, NULL);
446  exitcode= setjmp(qh->errexit);
447  if (!exitcode) {
448  coordT array[TOTpoints][DIM+1]; /* +1 for halfspace offset */
449  pointT *points;
450 
451  qh->NOerrexit= False;
452  strcat(qh->rbox_command, "user_eg2 halfspace example");
453  sprintf(options, "qhull H0 s Tcv %s", argc >= 4 ? argv[3] : "");
454  qh_initflags(qh, options);
455  printf( "\ncompute halfspace intersection about the origin for a diamond\n");
456  makehalf(array[0], SIZEcube, DIM);
457  qh_setfeasible(qh, DIM); /* from io.c, sets qh->feasible_point from 'Hn,n' */
458  /* you may malloc and set qh->feasible_point directly. It is only used for
459  option 'Fp' */
460  points= qh_sethalfspace_all(qh, DIM+1, SIZEcube, array[0], qh->feasible_point);
461  qh_init_B(qh, points, SIZEcube, DIM, True); /* qh_freeqhull frees points */
462  qh_qhull(qh);
463  qh_check_output(qh);
464  qh_produce_output(qh); /* delete this line to help avoid io.c */
465  if (qh->VERIFYoutput && !qh->STOPpoint && !qh->STOPcone)
466  qh_check_points(qh);
467  printf( "\nadd halfspaces for cube to intersection\n");
468  addhalf(qh, array[0], SIZEcube, SIZEdiamond, DIM, qh->feasible_point);
469  qh_check_output(qh);
470  qh_produce_output(qh); /* delete this line to help avoid io.c */
471  if (qh->VERIFYoutput && !qh->STOPpoint && !qh->STOPcone)
472  qh_check_points(qh);
473  }
474  qh->NOerrexit= True;
475  qh->NOerrexit= True;
476  qh_freeqhull(qh, !qh_ALL);
477  qh_memfreeshort(qh, &curlong, &totlong);
478  if (curlong || totlong)
479  fprintf(stderr, "qhull warning (user_eg2, run 3): did not free %d bytes of long memory (%d pieces)\n",
480  totlong, curlong);
481  return exitcode;
482 } /* main */
483 
484 #if 1 /* use 1 to prevent loading of io.o and user.o */
485 /*-------------------------------------------
486 -errexit- return exitcode to system after an error
487  assumes exitcode non-zero
488  prints useful information
489  see qh_errexit2() in libqhull.c for 2 facets
490 */
491 void qh_errexit(qhT *qh, int exitcode, facetT *facet, ridgeT *ridge) {
492  QHULL_UNUSED(facet);
493  QHULL_UNUSED(ridge);
494 
495  if (qh->ERREXITcalled) {
496  fprintf(qh->ferr, "qhull error while processing previous error. Exit program\n");
497  exit(1);
498  }
499  qh->ERREXITcalled= True;
500  if (!qh->QHULLfinished)
501  qh->hulltime= (unsigned)clock() - qh->hulltime;
502  fprintf(qh->ferr, "\nWhile executing: %s | %s\n", qh->rbox_command, qh->qhull_command);
503  fprintf(qh->ferr, "Options selected:\n%s\n", qh->qhull_options);
504  if (qh->furthest_id >= 0) {
505  fprintf(qh->ferr, "\nLast point added to hull was p%d", qh->furthest_id);
506  if (zzval_(Ztotmerge))
507  fprintf(qh->ferr, " Last merge was #%d.", zzval_(Ztotmerge));
508  if (qh->QHULLfinished)
509  fprintf(qh->ferr, "\nQhull has finished constructing the hull.");
510  else if (qh->POSTmerging)
511  fprintf(qh->ferr, "\nQhull has started post-merging");
512  fprintf(qh->ferr, "\n\n");
513  }
514  if (qh->NOerrexit) {
515  fprintf(qh->ferr, "qhull error while ending program. Exit program\n");
516  exit(1);
517  }
518  if (!exitcode)
519  exitcode= qh_ERRqhull;
520  qh->NOerrexit= True;
521  longjmp(qh->errexit, exitcode);
522 } /* errexit */
523 
524 
525 /*-------------------------------------------
526 -errprint- prints out the information of the erroneous object
527  any parameter may be NULL, also prints neighbors and geomview output
528 */
529 void qh_errprint(qhT *qh, const char *string, facetT *atfacet, facetT *otherfacet, ridgeT *atridge, vertexT *atvertex) {
530 
531  fprintf(qh->ferr, "%s facets f%d f%d ridge r%d vertex v%d\n",
532  string, getid_(atfacet), getid_(otherfacet), getid_(atridge),
533  getid_(atvertex));
534 } /* errprint */
535 
536 
537 void qh_printfacetlist(qhT *qh, facetT *facetlist, setT *facets, boolT printall) {
538  facetT *facet, **facetp;
539 
540  /* remove these calls to help avoid io.c */
541  qh_printbegin(qh, qh->ferr, qh_PRINTfacets, facetlist, facets, printall);/*io.c*/
542  FORALLfacet_(facetlist) /*io.c*/
543  qh_printafacet(qh, qh->ferr, qh_PRINTfacets, facet, printall); /*io.c*/
544  FOREACHfacet_(facets) /*io.c*/
545  qh_printafacet(qh, qh->ferr, qh_PRINTfacets, facet, printall); /*io.c*/
546  qh_printend(qh, qh->ferr, qh_PRINTfacets, facetlist, facets, printall); /*io.c*/
547 
548  FORALLfacet_(facetlist)
549  fprintf( qh->ferr, "facet f%d\n", facet->id);
550 } /* printfacetlist */
551 
552 /* qh_printhelp_degenerate( fp )
553  prints descriptive message for precision error
554 
555  notes:
556  no message if qh_QUICKhelp
557 */
558 void qh_printhelp_degenerate(qhT *qh, FILE *fp) {
559 
560  if (qh->MERGEexact || qh->PREmerge || qh->JOGGLEmax < REALmax/2)
561  qh_fprintf(qh, fp, 9368, "\n\
562 A Qhull error has occurred. Qhull should have corrected the above\n\
563 precision error. Please send the input and all of the output to\n\
564 qhull_bug@qhull.org\n");
565  else if (!qh_QUICKhelp) {
566  qh_fprintf(qh, fp, 9369, "\n\
567 Precision problems were detected during construction of the convex hull.\n\
568 This occurs because convex hull algorithms assume that calculations are\n\
569 exact, but floating-point arithmetic has roundoff errors.\n\
570 \n\
571 To correct for precision problems, do not use 'Q0'. By default, Qhull\n\
572 selects 'C-0' or 'Qx' and merges non-convex facets. With option 'QJ',\n\
573 Qhull joggles the input to prevent precision problems. See \"Imprecision\n\
574 in Qhull\" (qh-impre.htm).\n\
575 \n\
576 If you use 'Q0', the output may include\n\
577 coplanar ridges, concave ridges, and flipped facets. In 4-d and higher,\n\
578 Qhull may produce a ridge with four neighbors or two facets with the same \n\
579 vertices. Qhull reports these events when they occur. It stops when a\n\
580 concave ridge, flipped facet, or duplicate facet occurs.\n");
581 #if REALfloat
582  qh_fprintf(qh, fp, 9370, "\
583 \n\
584 Qhull is currently using single precision arithmetic. The following\n\
585 will probably remove the precision problems:\n\
586  - recompile qhull for realT precision(#define REALfloat 0 in user.h).\n");
587 #endif
588  if (qh->DELAUNAY && !qh->SCALElast && qh->MAXabs_coord > 1e4)
589  qh_fprintf(qh, fp, 9371, "\
590 \n\
591 When computing the Delaunay triangulation of coordinates > 1.0,\n\
592  - use 'Qbb' to scale the last coordinate to [0,m] (max previous coordinate)\n");
593  if (qh->DELAUNAY && !qh->ATinfinity)
594  qh_fprintf(qh, fp, 9372, "\
595 When computing the Delaunay triangulation:\n\
596  - use 'Qz' to add a point at-infinity. This reduces precision problems.\n");
597 
598  qh_fprintf(qh, fp, 9373, "\
599 \n\
600 If you need triangular output:\n\
601  - use option 'Qt' to triangulate the output\n\
602  - use option 'QJ' to joggle the input points and remove precision errors\n\
603  - use option 'Ft'. It triangulates non-simplicial facets with added points.\n\
604 \n\
605 If you must use 'Q0',\n\
606 try one or more of the following options. They can not guarantee an output.\n\
607  - use 'QbB' to scale the input to a cube.\n\
608  - use 'Po' to produce output and prevent partitioning for flipped facets\n\
609  - use 'V0' to set min. distance to visible facet as 0 instead of roundoff\n\
610  - use 'En' to specify a maximum roundoff error less than %2.2g.\n\
611  - options 'Qf', 'Qbb', and 'QR0' may also help\n",
612  qh->DISTround);
613  qh_fprintf(qh, fp, 9374, "\
614 \n\
615 To guarantee simplicial output:\n\
616  - use option 'Qt' to triangulate the output\n\
617  - use option 'QJ' to joggle the input points and remove precision errors\n\
618  - use option 'Ft' to triangulate the output by adding points\n\
619  - use exact arithmetic (see \"Imprecision in Qhull\", qh-impre.htm)\n\
620 ");
621  }
622 } /* printhelp_degenerate */
623 
624 
625 /* qh_printhelp_narrowhull( minangle )
626  Warn about a narrow hull
627 
628  notes:
629  Alternatively, reduce qh_WARNnarrow in user.h
630 
631 */
632 void qh_printhelp_narrowhull(qhT *qh, FILE *fp, realT minangle) {
633 
634  qh_fprintf(qh, fp, 9375, "qhull precision warning: \n\
635 The initial hull is narrow (cosine of min. angle is %.16f).\n\
636 A coplanar point may lead to a wide facet. Options 'QbB' (scale to unit box)\n\
637 or 'Qbb' (scale last coordinate) may remove this warning. Use 'Pp' to skip\n\
638 this warning. See 'Limitations' in qh-impre.htm.\n",
639  -minangle); /* convert from angle between normals to angle between facets */
640 } /* printhelp_narrowhull */
641 
642 /* qh_printhelp_singular
643  prints descriptive message for singular input
644 */
645 void qh_printhelp_singular(qhT *qh, FILE *fp) {
646  facetT *facet;
647  vertexT *vertex, **vertexp;
648  realT min, max, *coord, dist;
649  int i,k;
650 
651  qh_fprintf(qh, fp, 9376, "\n\
652 The input to qhull appears to be less than %d dimensional, or a\n\
653 computation has overflowed.\n\n\
654 Qhull could not construct a clearly convex simplex from points:\n",
655  qh->hull_dim);
656  qh_printvertexlist(qh, fp, "", qh->facet_list, NULL, qh_ALL);
657  if (!qh_QUICKhelp)
658  qh_fprintf(qh, fp, 9377, "\n\
659 The center point is coplanar with a facet, or a vertex is coplanar\n\
660 with a neighboring facet. The maximum round off error for\n\
661 computing distances is %2.2g. The center point, facets and distances\n\
662 to the center point are as follows:\n\n", qh->DISTround);
663  qh_printpointid(qh, fp, "center point", qh->hull_dim, qh->interior_point, -1);
664  qh_fprintf(qh, fp, 9378, "\n");
665  FORALLfacets {
666  qh_fprintf(qh, fp, 9379, "facet");
667  FOREACHvertex_(facet->vertices)
668  qh_fprintf(qh, fp, 9380, " p%d", qh_pointid(qh, vertex->point));
669  zinc_(Zdistio);
670  qh_distplane(qh, qh->interior_point, facet, &dist);
671  qh_fprintf(qh, fp, 9381, " distance= %4.2g\n", dist);
672  }
673  if (!qh_QUICKhelp) {
674  if (qh->HALFspace)
675  qh_fprintf(qh, fp, 9382, "\n\
676 These points are the dual of the given halfspaces. They indicate that\n\
677 the intersection is degenerate.\n");
678  qh_fprintf(qh, fp, 9383,"\n\
679 These points either have a maximum or minimum x-coordinate, or\n\
680 they maximize the determinant for k coordinates. Trial points\n\
681 are first selected from points that maximize a coordinate.\n");
682  if (qh->hull_dim >= qh_INITIALmax)
683  qh_fprintf(qh, fp, 9384, "\n\
684 Because of the high dimension, the min x-coordinate and max-coordinate\n\
685 points are used if the determinant is non-zero. Option 'Qs' will\n\
686 do a better, though much slower, job. Instead of 'Qs', you can change\n\
687 the points by randomly rotating the input with 'QR0'.\n");
688  }
689  qh_fprintf(qh, fp, 9385, "\nThe min and max coordinates for each dimension are:\n");
690  for (k=0; k < qh->hull_dim; k++) {
691  min= REALmax;
692  max= -REALmin;
693  for (i=qh->num_points, coord= qh->first_point+k; i--; coord += qh->hull_dim) {
694  maximize_(max, *coord);
695  minimize_(min, *coord);
696  }
697  qh_fprintf(qh, fp, 9386, " %d: %8.4g %8.4g difference= %4.4g\n", k, min, max, max-min);
698  }
699  if (!qh_QUICKhelp) {
700  qh_fprintf(qh, fp, 9387, "\n\
701 If the input should be full dimensional, you have several options that\n\
702 may determine an initial simplex:\n\
703  - use 'QJ' to joggle the input and make it full dimensional\n\
704  - use 'QbB' to scale the points to the unit cube\n\
705  - use 'QR0' to randomly rotate the input for different maximum points\n\
706  - use 'Qs' to search all points for the initial simplex\n\
707  - use 'En' to specify a maximum roundoff error less than %2.2g.\n\
708  - trace execution with 'T3' to see the determinant for each point.\n",
709  qh->DISTround);
710 #if REALfloat
711  qh_fprintf(qh, fp, 9388, "\
712  - recompile qhull for realT precision(#define REALfloat 0 in libqhull.h).\n");
713 #endif
714  qh_fprintf(qh, fp, 9389, "\n\
715 If the input is lower dimensional:\n\
716  - use 'QJ' to joggle the input and make it full dimensional\n\
717  - use 'Qbk:0Bk:0' to delete coordinate k from the input. You should\n\
718  pick the coordinate with the least range. The hull will have the\n\
719  correct topology.\n\
720  - determine the flat containing the points, rotate the points\n\
721  into a coordinate plane, and delete the other coordinates.\n\
722  - add one or more points to make the input full dimensional.\n\
723 ");
724  if (qh->DELAUNAY && !qh->ATinfinity)
725  qh_fprintf(qh, fp, 9390, "\n\n\
726 This is a Delaunay triangulation and the input is co-circular or co-spherical:\n\
727  - use 'Qz' to add a point \"at infinity\" (i.e., above the paraboloid)\n\
728  - or use 'QJ' to joggle the input and avoid co-circular data\n");
729  }
730 } /* printhelp_singular */
731 
732 
733 /*-----------------------------------------
734 -user_memsizes- allocate up to 10 additional, quick allocation sizes
735 */
737 
738  QHULL_UNUSED(qh);
739  /* qh_memsize(qh, size); */
740 } /* user_memsizes */
741 
742 #endif
int hull_dim
Definition: libqhull.h:591
Definition: libqhull.h:465
void print_summary(qhT *qh)
Definition: user_eg2_r.c:71
facetT * qh_findbestfacet(pointT *point, boolT bestoutside, realT *bestdist, boolT *isoutside)
Definition: poly2.c:1239
void seed(unsigned int seed_value)
#define qh_ERRqhull
Definition: libqhull.h:198
#define FORALLfacets
Definition: libqhull.h:840
void addhalf(qhT *qh, coordT *points, int numpoints, int numnew, int dim, coordT *feasible)
Definition: user_eg2_r.c:296
#define zinc_(id)
Definition: stat.h:386
#define boolT
Definition: libqhull.h:121
#define qh_RANDOMmax
Definition: user.h:282
Definition: qset.h:83
int STOPpoint
Definition: libqhull.h:567
int main(int argc, char *argv[])
Definition: user_eg2_r.c:346
#define qh_INITIALmax
Definition: user.h:420
#define FOREACHvertex_(vertices)
Definition: libqhull.h:950
void qh_printpoint(FILE *fp, const char *string, pointT *point)
Definition: io.c:2852
jmp_buf errexit
Definition: libqhull.h:656
facetT * facet_list
Definition: libqhull.h:677
#define pointT
Definition: libqhull.h:96
realT MAXabs_coord
Definition: libqhull.h:631
void qh_printhelp_singular(qhT *qh, FILE *fp)
Definition: user_eg2_r.c:645
#define getid_(p)
Definition: libqhull.h:823
void qh_fprintf(FILE *fp, int msgcode, const char *fmt,...)
Definition: userprintf.c:42
boolT ATinfinity
Definition: libqhull.h:482
coordT * qh_sethalfspace_all(int dim, int count, coordT *halfspaces, pointT *feasible)
Definition: geom2.c:1918
void qh_check_output(void)
Definition: poly2.c:302
void qh_printend(FILE *fp, qh_PRINT format, facetT *facetlist, setT *facets, boolT printall)
Definition: io.c:1690
pointT * interior_point
Definition: libqhull.h:663
#define FORALLfacet_(facetlist)
Definition: poly.h:77
pointT * first_point
Definition: libqhull.h:594
qhT qh_qh
Definition: global.c:26
void adddiamond(qhT *qh, coordT *points, int numpoints, int numnew, int dim)
Definition: user_eg2_r.c:115
boolT DELAUNAY
Definition: libqhull.h:491
int num_facets
Definition: libqhull.h:694
int qh_pointid(pointT *point)
Definition: poly.c:1053
void addDelaunay(qhT *qh, coordT *points, int numpoints, int numnew, int dim)
Definition: user_eg2_r.c:182
boolT HALFspace
Definition: libqhull.h:501
#define qh_RANDOMint
Definition: user.h:283
void qh_nearcoplanar(void)
Definition: poly2.c:2198
void qh_init_A(FILE *infile, FILE *outfile, FILE *errfile, int argc, char *argv[])
Definition: global.c:487
setT * other_points
Definition: libqhull.h:762
#define REALmax
Definition: user.h:155
realT JOGGLEmax
Definition: libqhull.h:721
#define qh_QUICKhelp
Definition: user.h:643
FILE * ferr
Definition: libqhull.h:662
setT * vertices
Definition: libqhull.h:295
boolT QHULLfinished
Definition: libqhull.h:743
boolT qh_addpoint(pointT *furthest, facetT *facet, boolT checkdist)
Definition: libqhull.c:168
void qh_triangulate(void)
Definition: poly2.c:2767
boolT qh_sethalfspace(int dim, coordT *coords, coordT **nextp, coordT *normal, coordT *offset, coordT *feasible)
Definition: geom2.c:1844
void qh_printhelp_degenerate(qhT *qh, FILE *fp)
Definition: user_eg2_r.c:558
unsigned long hulltime
Definition: libqhull.h:712
#define coordT
Definition: libqhull.h:80
#define True
Definition: libqhull.h:129
realT DISTround
Definition: libqhull.h:630
void qh_check_points(void)
Definition: poly2.c:365
#define zzval_(id)
Definition: stat.h:413
#define qh
Definition: libqhull.h:457
int num_points
Definition: libqhull.h:593
void qh_setdelaunay(int dim, int count, pointT *points)
Definition: geom2.c:1806
void qh_printpointid(FILE *fp, const char *string, int dim, pointT *point, int id)
Definition: io.c:2858
boolT PREmerge
Definition: libqhull.h:526
pointT * point
Definition: libqhull.h:399
Definition: stat.h:106
void qh_errprint(qhT *qh, const char *string, facetT *atfacet, facetT *otherfacet, ridgeT *atridge, vertexT *atvertex)
Definition: user_eg2_r.c:529
void qh_produce_output(void)
Definition: io.c:39
void qh_freeqhull(boolT allmem)
Definition: global.c:425
#define maximize_(maxval, val)
Definition: geom.h:51
boolT ERREXITcalled
Definition: libqhull.h:791
boolT NOerrexit
Definition: libqhull.h:737
#define TOTpoints
Definition: user_eg2_r.c:335
boolT VERIFYoutput
Definition: libqhull.h:580
int STOPcone
Definition: libqhull.h:565
coordT * feasible_point
Definition: libqhull.h:588
void findDelaunay(qhT *qh, int dim)
Definition: user_eg2_r.c:229
int furthest_id
Definition: libqhull.h:717
void qh_init_B(coordT *points, int numpoints, int dim, boolT ismalloc)
Definition: global.c:534
void qh_setfeasible(int dim)
Definition: io.c:3967
#define FOREACHfacet_(facets)
Definition: libqhull.h:891
boolT KEEPnearinside
Definition: libqhull.h:590
#define QHULL_LIB_CHECK
Definition: libqhull.h:462
#define SIZEcube
Definition: user_eg2_r.c:333
void qh_printbegin(FILE *fp, qh_PRINT format, facetT *facetlist, setT *facets, boolT printall)
Definition: io.c:1309
#define DIM
Definition: user_eg2_r.c:332
int num_vertices
Definition: libqhull.h:696
void qh_printhelp_narrowhull(qhT *qh, FILE *fp, realT minangle)
Definition: user_eg2_r.c:632
void qh_memfreeshort(int *curlong, int *totlong)
Definition: mem.c:288
void qh_user_memsizes(qhT *qh)
Definition: user_eg2_r.c:736
boolT POSTmerging
Definition: libqhull.h:740
boolT MERGEexact
Definition: libqhull.h:511
#define REALmin
Definition: user.h:156
#define SIZEdiamond
Definition: user_eg2_r.c:334
#define QHULL_UNUSED(x)
Definition: qhull_a.h:109
void makehalf(coordT *points, int numpoints, int dim)
Definition: user_eg2_r.c:259
void makeDelaunay(qhT *qh, coordT *points, int numpoints, int dim)
Definition: user_eg2_r.c:155
void qh_printafacet(FILE *fp, qh_PRINT format, facetT *facet, boolT printall)
Definition: io.c:1113
void qh_check_maxout(void)
Definition: poly2.c:218
void qh_initflags(char *command)
Definition: global.c:615
void qh_errexit(qhT *qh, int exitcode, facetT *facet, ridgeT *ridge)
Definition: user_eg2_r.c:491
#define False
Definition: libqhull.h:128
#define qh_ALL
Definition: libqhull.h:180
boolT DOcheckmax
Definition: libqhull.h:586
void qh_distplane(pointT *point, facetT *facet, realT *dist)
Definition: geom.c:36
int qh_setsize(setT *set)
Definition: qset.c:1111
boolT SCALElast
Definition: libqhull.h:559
unsigned id
Definition: libqhull.h:311
#define qh_RANDOMseed_(seed)
Definition: user.h:284
char qhull_command[256]
Definition: libqhull.h:598
char qhull_options[512]
Definition: libqhull.h:601
char rbox_command[256]
Definition: libqhull.h:600
#define realT
Definition: user.h:154
void qh_printfacetlist(qhT *qh, facetT *facetlist, setT *facets, boolT printall)
Definition: user_eg2_r.c:537
void makecube(coordT *points, int numpoints, int dim)
Definition: user_eg2_r.c:88
coordT * normal
Definition: libqhull.h:274
void qh_qhull(void)
Definition: libqhull.c:60
void qh_printvertexlist(FILE *fp, const char *string, facetT *facetlist, setT *facets, boolT printall)
Definition: io.c:3266
#define minimize_(minval, val)
Definition: geom.h:59
flagT tricoplanar
Definition: libqhull.h:314


hpp-fcl
Author(s):
autogenerated on Fri Jun 2 2023 02:39:02