1 /*<html><pre> -<a href="../libqhull_r/qh-qhull_r.htm"
2 >-------------------------------</a><a name="TOP">-</a>
3
4 qvoronoi_r.c
5 compute Voronoi diagrams and furthest-point Voronoi
6 diagrams using qhull
7
8 see unix_r.c for full interface
9
10 Copyright (c) 1993-2019, The Geometry Center
11 */
12
13 #include "libqhull_r/libqhull_r.h"
14
15 #include <stdio.h>
16 #include <stdlib.h>
17 #include <string.h>
18 #include <ctype.h>
19 #include <math.h>
20
21 #ifdef __cplusplus
22 extern "C" {
23 int isatty(int);
24 }
25
26 #elif defined(_MSC_VER)
27 #include <io.h>
28 #define isatty _isatty
29 /* int _isatty(int); */
30
31 #else
32 int isatty(int); /* returns 1 if stdin is a tty
33 if "Undefined symbol" this can be deleted along with call in main() */
34 #endif
35
36 /*-<a href="../libqhull_r/qh-qhull_r.htm#TOC"
37 >-------------------------------</a><a name="prompt">-</a>
38
39 qh_prompt
40 long prompt for qhull
41
42 notes:
43 restricted version of libqhull_r.c
44 same text as unix_r.c
45 see concise prompt below
46 limit maximum literal to 1800 characters
47 */
48
49 /* duplicated in qvoron_f.htm and qvoronoi.htm
50 QJ and Qt are deprecated, but allowed for backwards compatibility
51 */
52 char hidden_options[]=" d n m v H U Qb QB Qc Qf Qg Qi Qm Qr Qv Qx TR E V Fa FA FC FM Fp FS Ft FV Gt Pv Q0 Q1 Q2 Q3 Q4 Q5 Q6 Q7 Q8 Q9 Q10 Q11 Q15 ";
53
54 char qh_prompta[]= "\n\
55 qvoronoi -- compute the Voronoi diagram\n\
56 http://www.qhull.org %s\n\
57 \n\
58 input (stdin):\n\
59 first lines: dimension and number of points (or vice-versa).\n\
60 other lines: point coordinates, best if one point per line\n\
61 comments: start with a non-numeric character\n\
62 \n\
63 options:\n\
64 Qu - compute furthest-site Voronoi diagram\n\
65 \n\
66 Qhull control options:\n\
67 Qa - allow input with fewer or more points than coordinates\n\
68 QRn - random rotation (n=seed, n=0 time, n=-1 time/no rotate)\n\
69 Qs - search all points for the initial simplex\n\
70 Qz - add point-at-infinity to Voronoi diagram\n\
71 %s%s%s%s"; /* split up qh_prompt for Visual C++ */
72 char qh_promptb[]= "\
73 \n\
74 Qhull extra options:\n\
75 QGn - Voronoi vertices if visible from point n, -n if not\n\
76 QVn - Voronoi vertices for input point n, -n if not\n\
77 Qw - allow option warnings\n\
78 Q12 - allow wide facets and wide dupridge\n\
79 Q14 - merge pinched vertices that create a dupridge\n\
80 \n\
81 T options:\n\
82 TFn - report summary when n or more facets created\n\
83 TI file - input file, may be enclosed in single quotes\n\
84 TO file - output file, may be enclosed in single quotes\n\
85 Ts - statistics\n\
86 Tv - verify result: structure, convexity, and in-circle test\n\
87 Tz - send all output to stdout\n\
88 \n\
89 ";
90 char qh_promptc[]= "\
91 Trace options:\n\
92 T4 - trace at level n, 4=all, 5=mem/gauss, -1= events\n\
93 Ta - annotate output with message codes\n\
94 TAn - stop qhull after adding n vertices\n\
95 TCn - stop qhull after building cone for point n\n\
96 TVn - stop qhull after adding point n, -n for before\n\
97 Tc - check frequently during execution\n\
98 Tf - flush each qh_fprintf for debugging segfaults\n\
99 TPn - turn on tracing when point n added to hull\n\
100 TMn - turn on tracing at merge n\n\
101 TWn - trace merge facets when width > n\n\
102 \n\
103 Precision options:\n\
104 Cn - radius of centrum (roundoff added). Merge facets if non-convex\n\
105 An - cosine of maximum angle. Merge facets if cosine > n or non-convex\n\
106 C-0 roundoff, A-0.99/C-0.01 pre-merge, A0.99/C0.01 post-merge\n\
107 Rn - randomly perturb computations by a factor of [1-n,1+n]\n\
108 Wn - min facet width for non-coincident point (before roundoff)\n\
109 \n\
110 Output formats (may be combined; if none, summary to stdout):\n\
111 p - Voronoi vertices\n\
112 s - summary to stderr\n\
113 f - facet dump\n\
114 i - Delaunay regions (use 'Pp' to avoid warning)\n\
115 o - OFF format (dim, Voronoi vertices, and Voronoi regions)\n\
116 \n\
117 ";
118 char qh_promptd[]= "\
119 More formats:\n\
120 Fc - count plus coincident points (by Voronoi vertex)\n\
121 Fd - use cdd format for input (homogeneous with offset first)\n\
122 FD - use cdd format for output (offset first)\n\
123 FF - facet dump without ridges\n\
124 Fi - separating hyperplanes for bounded Voronoi regions\n\
125 FI - ID for each Voronoi vertex\n\
126 Fm - merge count for each Voronoi vertex (511 max)\n\
127 Fn - count plus neighboring Voronoi vertices for each Voronoi vertex\n\
128 FN - count and Voronoi vertices for each Voronoi region\n\
129 Fo - separating hyperplanes for unbounded Voronoi regions\n\
130 FO - options and precision constants\n\
131 FP - nearest point and distance for each coincident point\n\
132 FQ - command used for qvoronoi\n\
133 Fs - summary: #int (8), dimension, #points, tot vertices, tot facets,\n\
134 for output: #Voronoi regions, #Voronoi vertices,\n\
135 #coincident points, #non-simplicial regions\n\
136 #real (2), max outer plane and min vertex\n\
137 Fv - Voronoi diagram as Voronoi vertices between adjacent input sites\n\
138 Fx - extreme points of Delaunay triangulation (on convex hull)\n\
139 \n\
140 ";
141 char qh_prompte[]= "\
142 Geomview output (2-d only)\n\
143 Ga - all points as dots\n\
144 Gp - coplanar points and vertices as radii\n\
145 Gv - vertices as spheres\n\
146 Gc - centrums\n\
147 GDn - drop dimension n in 3-d and 4-d output\n\
148 Gh - hyperplane intersections\n\
149 Gi - inner planes only\n\
150 Gn - no planes\n\
151 Go - outer planes only\n\
152 Gr - ridges\n\
153 \n\
154 Print options:\n\
155 PAn - keep n largest Voronoi vertices by 'area'\n\
156 Pdk:n - drop facet if normal[k] <= n (default 0.0)\n\
157 PDk:n - drop facet if normal[k] >= n\n\
158 PFn - keep Voronoi vertices whose 'area' is at least n\n\
159 Pg - print good Voronoi vertices (needs 'QGn' or 'QVn')\n\
160 PG - print neighbors of good Voronoi vertices\n\
161 PMn - keep n Voronoi vertices with most merges\n\
162 Po - force output. If error, output neighborhood of facet\n\
163 Pp - do not report precision problems\n\
164 \n\
165 . - list of all options\n\
166 - - one line descriptions of all options\n\
167 -? - help with examples\n\
168 -V - version\n\
169 ";
170 /* for opts, don't assign 'e' or 'E' to a flag (already used for exponent) */
171
172 /*-<a href="../libqhull_r/qh-qhull_r.htm#TOC"
173 >-------------------------------</a><a name="prompt2">-</a>
174
175 qh_prompt2
176 synopsis for qhull
177 */
178 char qh_prompt2[]= "\n\
179 qvoronoi -- compute the Voronoi diagram. Qhull %s\n\
180 input (stdin): dimension, number of points, point coordinates\n\
181 comments start with a non-numeric character\n\
182 \n\
183 options (qvoronoi.htm):\n\
184 Qu - compute furthest-site Voronoi diagram\n\
185 Tv - verify result: structure, convexity, and in-circle test\n\
186 . - concise list of all options\n\
187 - - one-line description of all options\n\
188 -? - this message\n\
189 -V - version\n\
190 \n\
191 output options (subset):\n\
192 Fi - separating hyperplanes for bounded regions, 'Fo' for unbounded\n\
193 FN - count and Voronoi vertices for each Voronoi region\n\
194 Fv - Voronoi diagram as Voronoi vertices between adjacent input sites\n\
195 G - Geomview output (2-d only)\n\
196 o - OFF file format (dim, Voronoi vertices, and Voronoi regions)\n\
197 p - Voronoi vertices\n\
198 QVn - Voronoi vertices for input point n, -n if not\n\
199 s - summary of results (default)\n\
200 TI file - input file, may be enclosed in single quotes\n\
201 TO file - output file, may be enclosed in single quotes\n\
202 \n\
203 examples:\n\
204 rbox c P0 D2 | qvoronoi s o rbox c P0 D2 | qvoronoi Fi\n\
205 rbox c P0 D2 | qvoronoi Fo rbox c P0 D2 | qvoronoi Fv\n\
206 rbox c P0 D2 | qvoronoi s Qu Fv rbox c P0 D2 | qvoronoi Qu Fo\n\
207 rbox c G1 d D2 | qvoronoi s p rbox c P0 D2 | qvoronoi s Fv QV0\n\
208 \n\
209 ";
210 /* for opts, don't assign 'e' or 'E' to a flag (already used for exponent) */
211
212 /*-<a href="../libqhull_r/qh-qhull_r.htm#TOC"
213 >-------------------------------</a><a name="prompt3">-</a>
214
215 qh_prompt3
216 concise prompt for qhull
217 */
218 char qh_prompt3[]= "\n\
219 Qhull %s\n\
220 Except for 'F.' and 'PG', upper-case options take an argument.\n\
221 \n\
222 facet-dump Geomview i-delaunay off-format p-vertices\n\
223 summary\n\
224 \n\
225 Fcoincident Fd-cdd-in FD-cdd-out FF-dump-xridge Fi-bounded\n\
226 FIDs Fmerges Fneighbors FNeigh-region Fo-unbounded\n\
227 FOptions FPoint-near FQvoronoi Fsummary Fvoronoi\n\
228 Fxtremes\n\
229 \n\
230 Gall-points Gcentrums GDrop-dim Ghyperplanes Ginner\n\
231 Gno-planes Gouter Gpoints Gridges Gvertices\n\
232 \n\
233 PArea-keep Pdrop-d0:0D0 PFacet-area-keep Pgood PGood-neighbors\n\
234 PMerge-keep Poutput-forced Pprecision-not\n\
235 \n\
236 Qallow-short QG-vertex-good QRotate Qsearch-all Qupper-voronoi\n\
237 QV-point-good Qwarn-allow Qzinfinite Q12-allow-wide Q14-merge-pinched\n\
238 \n\
239 TFacet-log TInput-file TOutput-file Tstatistics Tverify\n\
240 Tz-stdout\n\
241 \n\
242 T4-trace Tannotate TAdd-stop Tcheck-often TCone-stop\n\
243 Tflush TMerge-trace TPoint-trace TVertex-stop TWide-trace\n\
244 \n\
245 Angle-max Centrum-size Random-dist Wide-outside\n\
246 ";
247
248 /*-<a href="../libqhull_r/qh-qhull_r.htm#TOC"
249 >-------------------------------</a><a name="main">-</a>
250
251 main( argc, argv )
252 processes the command line, calls qhull() to do the work, and exits
253
254 design:
255 initializes data structures
256 reads points
257 finishes initialization
258 computes convex hull and other structures
259 checks the result
260 writes the output
261 frees memory
262 */
main(int argc,char * argv[])263 int main(int argc, char *argv[]) {
264 int curlong, totlong; /* used !qh_NOmem */
265 int exitcode, numpoints, dim;
266 coordT *points;
267 boolT ismalloc;
268 qhT qh_qh;
269 qhT *qh= &qh_qh;
270
271 QHULL_LIB_CHECK /* Check for compatible library */
272
273 if ((argc == 1) && isatty( 0 /*stdin*/)) {
274 fprintf(stdout, qh_prompt2, qh_version);
275 exit(qh_ERRnone);
276 }
277 if (argc > 1 && *argv[1] == '-' && (*(argv[1] + 1) == '?' || *(argv[1] + 1) == '-')) { /* -? or --help */
278 fprintf(stdout, qh_prompt2, qh_version);
279 exit(qh_ERRnone);
280 }
281 if (argc > 1 && *argv[1] == '-' && !*(argv[1]+1)) {
282 fprintf(stdout, qh_prompta, qh_version,
283 qh_promptb, qh_promptc, qh_promptd, qh_prompte);
284 exit(qh_ERRnone);
285 }
286 if (argc > 1 && *argv[1] == '.' && !*(argv[1]+1)) {
287 fprintf(stdout, qh_prompt3, qh_version);
288 exit(qh_ERRnone);
289 }
290 if (argc > 1 && *argv[1] == '-' && *(argv[1]+1)=='V') {
291 fprintf(stdout, "%s\n", qh_version2);
292 exit(qh_ERRnone);
293 }
294 qh_init_A(qh, stdin, stdout, stderr, argc, argv); /* sets qh->qhull_command */
295 exitcode= setjmp(qh->errexit); /* simple statement for CRAY J916 */
296 if (!exitcode) {
297 qh->NOerrexit= False;
298 qh_option(qh, "voronoi _bbound-last _coplanar-keep", NULL, NULL);
299 qh->DELAUNAY= True; /* 'v' */
300 qh->VORONOI= True;
301 qh->SCALElast= True; /* 'Qbb' */
302 qh_checkflags(qh, qh->qhull_command, hidden_options);
303 qh_initflags(qh, qh->qhull_command);
304 points= qh_readpoints(qh, &numpoints, &dim, &ismalloc);
305 qh_init_B(qh, points, numpoints, dim, ismalloc);
306 qh_qhull(qh);
307 qh_check_output(qh);
308 qh_produce_output(qh);
309 if (qh->VERIFYoutput && !qh->FORCEoutput && !qh->STOPpoint && !qh->STOPcone)
310 qh_check_points(qh);
311 exitcode= qh_ERRnone;
312 }
313 qh->NOerrexit= True; /* no more setjmp */
314 #ifdef qh_NOmem
315 qh_freeqhull(qh, qh_ALL);
316 #else
317 qh_freeqhull(qh, !qh_ALL);
318 qh_memfreeshort(qh, &curlong, &totlong);
319 if (curlong || totlong)
320 qh_fprintf_stderr(7079, "qhull internal warning (main): did not free %d bytes of long memory(%d pieces)\n",
321 totlong, curlong);
322 #endif
323 return exitcode;
324 } /* main */
325
326