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

hpp-fcl
Author(s):
autogenerated on Fri Jan 26 2024 03:46:15