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

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