1 /*<html><pre>  -<a                             href="index.htm#TOC"
2   >-------------------------------</a><a name="TOP">-</a>
3 
4    rboxlib.c
5      Generate input points
6 
7    notes:
8      For documentation, see prompt[] of rbox.c
9      50 points generated for 'rbox D4'
10 
11    WARNING:
12      incorrect range if qh_RANDOMmax is defined wrong (user.h)
13 */
14 
15 #include "random.h"
16 #include "libqhull.h"
17 
18 #include <ctype.h>
19 #include <limits.h>
20 #include <math.h>
21 #include <setjmp.h>
22 #include <string.h>
23 #include <time.h>
24 #include <stdio.h>
25 #include <stdlib.h>
26 
27 #ifdef _MSC_VER  /* Microsoft Visual C++ */
28 #pragma warning( disable : 4706)  /* assignment within conditional expression. */
29 #pragma warning( disable : 4996)  /* this function (strncat,sprintf,strcpy) or variable may be unsafe. */
30 #endif
31 
32 #define MAXdim 200
33 
34 /* ------------------------------ prototypes ----------------*/
35 int roundi( double a);
36 void out1( double a);
37 void out2n( double a, double b);
38 void out3n( double a, double b, double c);
39 
40 void    qh_fprintf_rbox(FILE *fp, int msgcode, const char *fmt, ... );
41 void    qh_free(void *mem);
42 void   *qh_malloc(size_t size);
43 int     qh_rand( void);
44 void    qh_srand( int seed);
45 
46 
47 /* ------------------------------ globals -------------------*/
48 
49 /* No state is carried between rbox requests */
50 typedef struct rboxT rboxT;
51 struct rboxT {
52   FILE *fout;
53   FILE *ferr;
54   int isinteger;
55   double out_offset;
56   jmp_buf errexit;        /* exit label for rboxpoints, defined by setjmp(), called by qh_errexit_rbox() */
57 };
58 
59 
60 int rbox_inuse= 0;
61 rboxT rbox;
62 
63 /*-<a                             href="qh-qhull.htm#TOC"
64   >-------------------------------</a><a name="rboxpoints">-</a>
65 
66   qh_rboxpoints( fout, ferr, rbox_command )
67     Generate points to fout according to rbox options
68     Report errors on ferr
69 
70   returns:
71     0 (qh_ERRnone) on success
72     1 (qh_ERRinput) on input error
73     4 (qh_ERRmem) on memory error
74     5 (qh_ERRqhull) on internal error
75 
76   notes:
77     To avoid stdio, redefine qh_malloc, qh_free, and qh_fprintf_rbox (user.c)
78     Rbox is not multithreaded.
79 
80   design:
81     Straight line code (consider defining a struct and functions):
82 
83     Parse arguments into variables
84     Determine the number of points
85     Generate the points
86 */
qh_rboxpoints(FILE * fout,FILE * ferr,char * rbox_command)87 int qh_rboxpoints(FILE* fout, FILE* ferr, char* rbox_command) {
88   int i,j,k;
89   int gendim;
90   int cubesize, diamondsize, seed=0, count, apex;
91   int dim=3 , numpoints= 0, totpoints, addpoints=0;
92   int issphere=0, isaxis=0,  iscdd= 0, islens= 0, isregular=0, iswidth=0, addcube=0;
93   int isgap=0, isspiral=0, NOcommand= 0, adddiamond=0;
94   int israndom=0, istime=0;
95   int isbox=0, issimplex=0, issimplex2=0, ismesh=0;
96   double width=0.0, gap=0.0, radius= 0.0;
97   double coord[MAXdim], offset, meshm=3.0, meshn=4.0, meshr=5.0;
98   double *simplex= NULL, *simplexp;
99   int nthroot, mult[MAXdim];
100   double norm, factor, randr, rangap, lensangle= 0, lensbase= 1;
101   double anglediff, angle, x, y, cube= 0.0, diamond= 0.0;
102   double box= qh_DEFAULTbox; /* scale all numbers before output */
103   double randmax= qh_RANDOMmax;
104   char command[200], seedbuf[200];
105   char *s= command, *t, *first_point= NULL;
106   time_t timedata = 0;
107   int exitcode;
108 
109   if (rbox_inuse) {
110     qh_fprintf_rbox(rbox.ferr, 6188, "rbox error: rbox in use by another process.  Please lock calls to rbox.\n");
111     return qh_ERRqhull;
112   }
113   rbox_inuse= True;
114   rbox.ferr= ferr;
115   rbox.fout= fout;
116 
117   exitcode= setjmp(rbox.errexit);
118   if (exitcode) {
119     /* same code for error exit and normal return */
120     if (simplex)
121         qh_free(simplex);
122     rbox_inuse= False;
123     return exitcode;
124   }
125 
126   *command= '\0';
127   strncat(command, rbox_command, sizeof(command)-strlen(command)-1);
128 
129   while (*s && !isspace(*s))  /* skip program name */
130     s++;
131   while (*s) {
132     while (*s && isspace(*s))
133       s++;
134     if (*s == '-')
135       s++;
136     if (!*s)
137       break;
138     if (isdigit(*s)) {
139       numpoints= qh_strtol(s, &s);
140       continue;
141     }
142     /* ============= read flags =============== */
143     switch (*s++) {
144     case 'c':
145       addcube= 1;
146       t= s;
147       while (isspace(*t))
148         t++;
149       if (*t == 'G')
150         cube= qh_strtod(++t, &s);
151       break;
152     case 'd':
153       adddiamond= 1;
154       t= s;
155       while (isspace(*t))
156         t++;
157       if (*t == 'G')
158         diamond= qh_strtod(++t, &s);
159       break;
160     case 'h':
161       iscdd= 1;
162       break;
163     case 'l':
164       isspiral= 1;
165       break;
166     case 'n':
167       NOcommand= 1;
168       break;
169     case 'r':
170       isregular= 1;
171       break;
172     case 's':
173       issphere= 1;
174       break;
175     case 't':
176       istime= 1;
177       if (isdigit(*s)) {
178         seed= qh_strtol(s, &s);
179         israndom= 0;
180       }else
181         israndom= 1;
182       break;
183     case 'x':
184       issimplex= 1;
185       break;
186     case 'y':
187       issimplex2= 1;
188       break;
189     case 'z':
190       rbox.isinteger= 1;
191       break;
192     case 'B':
193       box= qh_strtod(s, &s);
194       isbox= 1;
195       break;
196     case 'D':
197       dim= qh_strtol(s, &s);
198       if (dim < 1
199       || dim > MAXdim) {
200         qh_fprintf_rbox(rbox.ferr, 6189, "rbox error: dimension, D%d, out of bounds (>=%d or <=0)", dim, MAXdim);
201         qh_errexit_rbox(qh_ERRinput);
202       }
203       break;
204     case 'G':
205       if (isdigit(*s))
206         gap= qh_strtod(s, &s);
207       else
208         gap= 0.5;
209       isgap= 1;
210       break;
211     case 'L':
212       if (isdigit(*s))
213         radius= qh_strtod(s, &s);
214       else
215         radius= 10;
216       islens= 1;
217       break;
218     case 'M':
219       ismesh= 1;
220       if (*s)
221         meshn= qh_strtod(s, &s);
222       if (*s == ',') {
223         ++s;
224         meshm= qh_strtod(s, &s);
225       }else
226         meshm= 0.0;
227       if (*s == ',') {
228         ++s;
229         meshr= qh_strtod(s, &s);
230       }else
231         meshr= sqrt(meshn*meshn + meshm*meshm);
232       if (*s && !isspace(*s)) {
233         qh_fprintf_rbox(rbox.ferr, 7069, "rbox warning: assuming 'M3,4,5' since mesh args are not integers or reals\n");
234         meshn= 3.0, meshm=4.0, meshr=5.0;
235       }
236       break;
237     case 'O':
238       rbox.out_offset= qh_strtod(s, &s);
239       break;
240     case 'P':
241       if (!first_point)
242         first_point= s-1;
243       addpoints++;
244       while (*s && !isspace(*s))   /* read points later */
245         s++;
246       break;
247     case 'W':
248       width= qh_strtod(s, &s);
249       iswidth= 1;
250       break;
251     case 'Z':
252       if (isdigit(*s))
253         radius= qh_strtod(s, &s);
254       else
255         radius= 1.0;
256       isaxis= 1;
257       break;
258     default:
259       qh_fprintf_rbox(rbox.ferr, 7070, "rbox error: unknown flag at %s.\nExecute 'rbox' without arguments for documentation.\n", s);
260       qh_errexit_rbox(qh_ERRinput);
261     }
262     if (*s && !isspace(*s)) {
263       qh_fprintf_rbox(rbox.ferr, 7071, "rbox error: missing space between flags at %s.\n", s);
264       qh_errexit_rbox(qh_ERRinput);
265     }
266   }
267 
268   /* ============= defaults, constants, and sizes =============== */
269   if (rbox.isinteger && !isbox)
270     box= qh_DEFAULTzbox;
271   if (addcube) {
272     cubesize= (int)floor(ldexp(1.0,dim)+0.5);
273     if (cube == 0.0)
274       cube= box;
275   }else
276     cubesize= 0;
277   if (adddiamond) {
278     diamondsize= 2*dim;
279     if (diamond == 0.0)
280       diamond= box;
281   }else
282     diamondsize= 0;
283   if (islens) {
284     if (isaxis) {
285         qh_fprintf_rbox(rbox.ferr, 6190, "rbox error: can not combine 'Ln' with 'Zn'\n");
286         qh_errexit_rbox(qh_ERRinput);
287     }
288     if (radius <= 1.0) {
289         qh_fprintf_rbox(rbox.ferr, 6191, "rbox error: lens radius %.2g should be greater than 1.0\n",
290                radius);
291         qh_errexit_rbox(qh_ERRinput);
292     }
293     lensangle= asin(1.0/radius);
294     lensbase= radius * cos(lensangle);
295   }
296 
297   if (!numpoints) {
298     if (issimplex2)
299         ; /* ok */
300     else if (isregular + issimplex + islens + issphere + isaxis + isspiral + iswidth + ismesh) {
301         qh_fprintf_rbox(rbox.ferr, 6192, "rbox error: missing count\n");
302         qh_errexit_rbox(qh_ERRinput);
303     }else if (adddiamond + addcube + addpoints)
304         ; /* ok */
305     else {
306         numpoints= 50;  /* ./rbox D4 is the test case */
307         issphere= 1;
308     }
309   }
310   if ((issimplex + islens + isspiral + ismesh > 1)
311   || (issimplex + issphere + isspiral + ismesh > 1)) {
312     qh_fprintf_rbox(rbox.ferr, 6193, "rbox error: can only specify one of 'l', 's', 'x', 'Ln', or 'Mn,m,r' ('Ln s' is ok).\n");
313     qh_errexit_rbox(qh_ERRinput);
314   }
315 
316   /* ============= print header with total points =============== */
317   if (issimplex || ismesh)
318     totpoints= numpoints;
319   else if (issimplex2)
320     totpoints= numpoints+dim+1;
321   else if (isregular) {
322     totpoints= numpoints;
323     if (dim == 2) {
324         if (islens)
325           totpoints += numpoints - 2;
326     }else if (dim == 3) {
327         if (islens)
328           totpoints += 2 * numpoints;
329       else if (isgap)
330         totpoints += 1 + numpoints;
331       else
332         totpoints += 2;
333     }
334   }else
335     totpoints= numpoints + isaxis;
336   totpoints += cubesize + diamondsize + addpoints;
337 
338   /* ============= seed randoms =============== */
339   if (istime == 0) {
340     for (s=command; *s; s++) {
341       if (issimplex2 && *s == 'y') /* make 'y' same seed as 'x' */
342         i= 'x';
343       else
344         i= *s;
345       seed= 11*seed + i;
346     }
347   }else if (israndom) {
348     seed= (int)time(&timedata);
349     sprintf(seedbuf, " t%d", seed);  /* appends an extra t, not worth removing */
350     strncat(command, seedbuf, sizeof(command)-strlen(command)-1);
351     t= strstr(command, " t ");
352     if (t)
353       strcpy(t+1, t+3); /* remove " t " */
354   } /* else, seed explicitly set to n */
355   qh_RANDOMseed_(seed);
356 
357   /* ============= print header =============== */
358 
359   if (iscdd)
360       qh_fprintf_rbox(rbox.fout, 9391, "%s\nbegin\n        %d %d %s\n",
361       NOcommand ? "" : command,
362       totpoints, dim+1,
363       rbox.isinteger ? "integer" : "real");
364   else if (NOcommand)
365       qh_fprintf_rbox(rbox.fout, 9392, "%d\n%d\n", dim, totpoints);
366   else
367       qh_fprintf_rbox(rbox.fout, 9393, "%d %s\n%d\n", dim, command, totpoints);
368 
369   /* ============= explicit points =============== */
370   if ((s= first_point)) {
371     while (s && *s) { /* 'P' */
372       count= 0;
373       if (iscdd)
374         out1( 1.0);
375       while (*++s) {
376         out1( qh_strtod(s, &s));
377         count++;
378         if (isspace(*s) || !*s)
379           break;
380         if (*s != ',') {
381           qh_fprintf_rbox(rbox.ferr, 6194, "rbox error: missing comma after coordinate in %s\n\n", s);
382           qh_errexit_rbox(qh_ERRinput);
383         }
384       }
385       if (count < dim) {
386         for (k=dim-count; k--; )
387           out1( 0.0);
388       }else if (count > dim) {
389         qh_fprintf_rbox(rbox.ferr, 6195, "rbox error: %d coordinates instead of %d coordinates in %s\n\n",
390                   count, dim, s);
391         qh_errexit_rbox(qh_ERRinput);
392       }
393       qh_fprintf_rbox(rbox.fout, 9394, "\n");
394       while ((s= strchr(s, 'P'))) {
395         if (isspace(s[-1]))
396           break;
397       }
398     }
399   }
400 
401   /* ============= simplex distribution =============== */
402   if (issimplex+issimplex2) {
403     if (!(simplex= (double*)qh_malloc( dim * (dim+1) * sizeof(double)))) {
404       qh_fprintf_rbox(rbox.ferr, 6196, "rbox error: insufficient memory for simplex\n");
405       qh_errexit_rbox(qh_ERRmem); /* qh_ERRmem */
406     }
407     simplexp= simplex;
408     if (isregular) {
409       for (i=0; i<dim; i++) {
410         for (k=0; k<dim; k++)
411           *(simplexp++)= i==k ? 1.0 : 0.0;
412       }
413       for (k=0; k<dim; k++)
414         *(simplexp++)= -1.0;
415     }else {
416       for (i=0; i<dim+1; i++) {
417         for (k=0; k<dim; k++) {
418           randr= qh_RANDOMint;
419           *(simplexp++)= 2.0 * randr/randmax - 1.0;
420         }
421       }
422     }
423     if (issimplex2) {
424         simplexp= simplex;
425       for (i=0; i<dim+1; i++) {
426         if (iscdd)
427           out1( 1.0);
428         for (k=0; k<dim; k++)
429           out1( *(simplexp++) * box);
430         qh_fprintf_rbox(rbox.fout, 9395, "\n");
431       }
432     }
433     for (j=0; j<numpoints; j++) {
434       if (iswidth)
435         apex= qh_RANDOMint % (dim+1);
436       else
437         apex= -1;
438       for (k=0; k<dim; k++)
439         coord[k]= 0.0;
440       norm= 0.0;
441       for (i=0; i<dim+1; i++) {
442         randr= qh_RANDOMint;
443         factor= randr/randmax;
444         if (i == apex)
445           factor *= width;
446         norm += factor;
447         for (k=0; k<dim; k++) {
448           simplexp= simplex + i*dim + k;
449           coord[k] += factor * (*simplexp);
450         }
451       }
452       for (k=0; k<dim; k++)
453         coord[k] /= norm;
454       if (iscdd)
455         out1( 1.0);
456       for (k=0; k < dim; k++)
457         out1( coord[k] * box);
458       qh_fprintf_rbox(rbox.fout, 9396, "\n");
459     }
460     isregular= 0; /* continue with isbox */
461     numpoints= 0;
462   }
463 
464   /* ============= mesh distribution =============== */
465   if (ismesh) {
466     nthroot= (int)(pow((double)numpoints, 1.0/dim) + 0.99999);
467     for (k=dim; k--; )
468       mult[k]= 0;
469     for (i=0; i < numpoints; i++) {
470       for (k=0; k < dim; k++) {
471         if (k == 0)
472           out1( mult[0] * meshn + mult[1] * (-meshm));
473         else if (k == 1)
474           out1( mult[0] * meshm + mult[1] * meshn);
475         else
476           out1( mult[k] * meshr );
477       }
478       qh_fprintf_rbox(rbox.fout, 9397, "\n");
479       for (k=0; k < dim; k++) {
480         if (++mult[k] < nthroot)
481           break;
482         mult[k]= 0;
483       }
484     }
485   }
486   /* ============= regular points for 's' =============== */
487   else if (isregular && !islens) {
488     if (dim != 2 && dim != 3) {
489       qh_fprintf_rbox(rbox.ferr, 6197, "rbox error: regular points can be used only in 2-d and 3-d\n\n");
490       qh_errexit_rbox(qh_ERRinput);
491     }
492     if (!isaxis || radius == 0.0) {
493       isaxis= 1;
494       radius= 1.0;
495     }
496     if (dim == 3) {
497       if (iscdd)
498         out1( 1.0);
499       out3n( 0.0, 0.0, -box);
500       if (!isgap) {
501         if (iscdd)
502           out1( 1.0);
503         out3n( 0.0, 0.0, box);
504       }
505     }
506     angle= 0.0;
507     anglediff= 2.0 * M_PI/numpoints;
508     for (i=0; i < numpoints; i++) {
509       angle += anglediff;
510       x= radius * cos(angle);
511       y= radius * sin(angle);
512       if (dim == 2) {
513         if (iscdd)
514           out1( 1.0);
515         out2n( x*box, y*box);
516       }else {
517         norm= sqrt(1.0 + x*x + y*y);
518         if (iscdd)
519           out1( 1.0);
520         out3n( box*x/norm, box*y/norm, box/norm);
521         if (isgap) {
522           x *= 1-gap;
523           y *= 1-gap;
524           norm= sqrt(1.0 + x*x + y*y);
525           if (iscdd)
526             out1( 1.0);
527           out3n( box*x/norm, box*y/norm, box/norm);
528         }
529       }
530     }
531   }
532   /* ============= regular points for 'r Ln D2' =============== */
533   else if (isregular && islens && dim == 2) {
534     double cos_0;
535 
536     angle= lensangle;
537     anglediff= 2 * lensangle/(numpoints - 1);
538     cos_0= cos(lensangle);
539     for (i=0; i < numpoints; i++, angle -= anglediff) {
540       x= radius * sin(angle);
541       y= radius * (cos(angle) - cos_0);
542       if (iscdd)
543         out1( 1.0);
544       out2n( x*box, y*box);
545       if (i != 0 && i != numpoints - 1) {
546         if (iscdd)
547           out1( 1.0);
548         out2n( x*box, -y*box);
549       }
550     }
551   }
552   /* ============= regular points for 'r Ln D3' =============== */
553   else if (isregular && islens && dim != 2) {
554     if (dim != 3) {
555       qh_fprintf_rbox(rbox.ferr, 6198, "rbox error: regular points can be used only in 2-d and 3-d\n\n");
556       qh_errexit_rbox(qh_ERRinput);
557     }
558     angle= 0.0;
559     anglediff= 2* M_PI/numpoints;
560     if (!isgap) {
561       isgap= 1;
562       gap= 0.5;
563     }
564     offset= sqrt(radius * radius - (1-gap)*(1-gap)) - lensbase;
565     for (i=0; i < numpoints; i++, angle += anglediff) {
566       x= cos(angle);
567       y= sin(angle);
568       if (iscdd)
569         out1( 1.0);
570       out3n( box*x, box*y, 0.0);
571       x *= 1-gap;
572       y *= 1-gap;
573       if (iscdd)
574         out1( 1.0);
575       out3n( box*x, box*y, box * offset);
576       if (iscdd)
577         out1( 1.0);
578       out3n( box*x, box*y, -box * offset);
579     }
580   }
581   /* ============= apex of 'Zn' distribution + gendim =============== */
582   else {
583     if (isaxis) {
584       gendim= dim-1;
585       if (iscdd)
586         out1( 1.0);
587       for (j=0; j < gendim; j++)
588         out1( 0.0);
589       out1( -box);
590       qh_fprintf_rbox(rbox.fout, 9398, "\n");
591     }else if (islens)
592       gendim= dim-1;
593     else
594       gendim= dim;
595     /* ============= generate random point in unit cube =============== */
596     for (i=0; i < numpoints; i++) {
597       norm= 0.0;
598       for (j=0; j < gendim; j++) {
599         randr= qh_RANDOMint;
600         coord[j]= 2.0 * randr/randmax - 1.0;
601         norm += coord[j] * coord[j];
602       }
603       norm= sqrt(norm);
604       /* ============= dim-1 point of 'Zn' distribution ========== */
605       if (isaxis) {
606         if (!isgap) {
607           isgap= 1;
608           gap= 1.0;
609         }
610         randr= qh_RANDOMint;
611         rangap= 1.0 - gap * randr/randmax;
612         factor= radius * rangap / norm;
613         for (j=0; j<gendim; j++)
614           coord[j]= factor * coord[j];
615       /* ============= dim-1 point of 'Ln s' distribution =========== */
616       }else if (islens && issphere) {
617         if (!isgap) {
618           isgap= 1;
619           gap= 1.0;
620         }
621         randr= qh_RANDOMint;
622         rangap= 1.0 - gap * randr/randmax;
623         factor= rangap / norm;
624         for (j=0; j<gendim; j++)
625           coord[j]= factor * coord[j];
626       /* ============= dim-1 point of 'Ln' distribution ========== */
627       }else if (islens && !issphere) {
628         if (!isgap) {
629           isgap= 1;
630           gap= 1.0;
631         }
632         j= qh_RANDOMint % gendim;
633         if (coord[j] < 0)
634           coord[j]= -1.0 - coord[j] * gap;
635         else
636           coord[j]= 1.0 - coord[j] * gap;
637       /* ============= point of 'l' distribution =============== */
638       }else if (isspiral) {
639         if (dim != 3) {
640           qh_fprintf_rbox(rbox.ferr, 6199, "rbox error: spiral distribution is available only in 3d\n\n");
641           longjmp(rbox.errexit,qh_ERRinput);
642         }
643         coord[0]= cos(2*M_PI*i/(numpoints - 1));
644         coord[1]= sin(2*M_PI*i/(numpoints - 1));
645         coord[2]= 2.0*(double)i/(double)(numpoints-1) - 1.0;
646       /* ============= point of 's' distribution =============== */
647       }else if (issphere) {
648         factor= 1.0/norm;
649         if (iswidth) {
650           randr= qh_RANDOMint;
651           factor *= 1.0 - width * randr/randmax;
652         }
653         for (j=0; j<dim; j++)
654           coord[j]= factor * coord[j];
655       }
656       /* ============= project 'Zn s' point in to sphere =============== */
657       if (isaxis && issphere) {
658         coord[dim-1]= 1.0;
659         norm= 1.0;
660         for (j=0; j<gendim; j++)
661           norm += coord[j] * coord[j];
662         norm= sqrt(norm);
663         for (j=0; j<dim; j++)
664           coord[j]= coord[j] / norm;
665         if (iswidth) {
666           randr= qh_RANDOMint;
667           coord[dim-1] *= 1 - width * randr/randmax;
668         }
669       /* ============= project 'Zn' point onto cube =============== */
670       }else if (isaxis && !issphere) {  /* not very interesting */
671         randr= qh_RANDOMint;
672         coord[dim-1]= 2.0 * randr/randmax - 1.0;
673       /* ============= project 'Ln' point out to sphere =============== */
674       }else if (islens) {
675         coord[dim-1]= lensbase;
676         for (j=0, norm= 0; j<dim; j++)
677           norm += coord[j] * coord[j];
678         norm= sqrt(norm);
679         for (j=0; j<dim; j++)
680           coord[j]= coord[j] * radius/ norm;
681         coord[dim-1] -= lensbase;
682         if (iswidth) {
683           randr= qh_RANDOMint;
684           coord[dim-1] *= 1 - width * randr/randmax;
685         }
686         if (qh_RANDOMint > randmax/2)
687           coord[dim-1]= -coord[dim-1];
688       /* ============= project 'Wn' point toward boundary =============== */
689       }else if (iswidth && !issphere) {
690         j= qh_RANDOMint % gendim;
691         if (coord[j] < 0)
692           coord[j]= -1.0 - coord[j] * width;
693         else
694           coord[j]= 1.0 - coord[j] * width;
695       }
696       /* ============= write point =============== */
697       if (iscdd)
698         out1( 1.0);
699       for (k=0; k < dim; k++)
700         out1( coord[k] * box);
701       qh_fprintf_rbox(rbox.fout, 9399, "\n");
702     }
703   }
704 
705   /* ============= write cube vertices =============== */
706   if (addcube) {
707     for (j=0; j<cubesize; j++) {
708       if (iscdd)
709         out1( 1.0);
710       for (k=dim-1; k>=0; k--) {
711         if (j & ( 1 << k))
712           out1( cube);
713         else
714           out1( -cube);
715       }
716       qh_fprintf_rbox(rbox.fout, 9400, "\n");
717     }
718   }
719 
720   /* ============= write diamond vertices =============== */
721   if (adddiamond) {
722     for (j=0; j<diamondsize; j++) {
723       if (iscdd)
724         out1( 1.0);
725       for (k=dim-1; k>=0; k--) {
726         if (j/2 != k)
727           out1( 0.0);
728         else if (j & 0x1)
729           out1( diamond);
730         else
731           out1( -diamond);
732       }
733       qh_fprintf_rbox(rbox.fout, 9401, "\n");
734     }
735   }
736 
737   if (iscdd)
738     qh_fprintf_rbox(rbox.fout, 9402, "end\nhull\n");
739 
740   /* same code for error exit and normal return */
741   if (simplex)
742     qh_free(simplex);
743   rbox_inuse= False;
744   return qh_ERRnone;
745 } /* rboxpoints */
746 
747 /*------------------------------------------------
748 outxxx - output functions
749 */
roundi(double a)750 int roundi( double a) {
751   if (a < 0.0) {
752     if (a - 0.5 < INT_MIN) {
753       qh_fprintf_rbox(rbox.ferr, 6200, "rbox input error: negative coordinate %2.2g is too large.  Reduce 'Bn'\n", a);
754       qh_errexit_rbox(qh_ERRinput);
755     }
756     return (int)(a - 0.5);
757   }else {
758     if (a + 0.5 > INT_MAX) {
759       qh_fprintf_rbox(rbox.ferr, 6201, "rbox input error: coordinate %2.2g is too large.  Reduce 'Bn'\n", a);
760       qh_errexit_rbox(qh_ERRinput);
761     }
762     return (int)(a + 0.5);
763   }
764 } /* roundi */
765 
out1(double a)766 void out1(double a) {
767 
768   if (rbox.isinteger)
769     qh_fprintf_rbox(rbox.fout, 9403, "%d ", roundi( a+rbox.out_offset));
770   else
771     qh_fprintf_rbox(rbox.fout, 9404, qh_REAL_1, a+rbox.out_offset);
772 } /* out1 */
773 
out2n(double a,double b)774 void out2n( double a, double b) {
775 
776   if (rbox.isinteger)
777     qh_fprintf_rbox(rbox.fout, 9405, "%d %d\n", roundi(a+rbox.out_offset), roundi(b+rbox.out_offset));
778   else
779     qh_fprintf_rbox(rbox.fout, 9406, qh_REAL_2n, a+rbox.out_offset, b+rbox.out_offset);
780 } /* out2n */
781 
out3n(double a,double b,double c)782 void out3n( double a, double b, double c) {
783 
784   if (rbox.isinteger)
785     qh_fprintf_rbox(rbox.fout, 9407, "%d %d %d\n", roundi(a+rbox.out_offset), roundi(b+rbox.out_offset), roundi(c+rbox.out_offset));
786   else
787     qh_fprintf_rbox(rbox.fout, 9408, qh_REAL_3n, a+rbox.out_offset, b+rbox.out_offset, c+rbox.out_offset);
788 } /* out3n */
789 
qh_errexit_rbox(int exitcode)790 void qh_errexit_rbox(int exitcode)
791 {
792     longjmp(rbox.errexit, exitcode);
793 } /* rbox_errexit */
794