📄 rbox.c
字号:
/*<html><pre> -<a href="qh-man.htm"
>-------------------------------</a><a name="TOP">-</a>
rbox.c
Generate input points for qhull.
notes:
50 points generated for 'rbox D4'
This code needs a full rewrite. It needs separate procedures for each
distribution with common, helper procedures.
WARNING:
incorrect range if qh_RANDOMmax is defined wrong (user.h)
*/
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <ctype.h>
#include <math.h>
#include <limits.h>
#include <time.h>
#include "user.h"
#if __MWERKS__ && __POWERPC__
#include <SIOUX.h>
#include <Files.h>
#include <console.h>
#include <Desk.h>
#endif
#ifdef _MSC_VER /* Microsoft Visual C++ */
#pragma warning( disable : 4244) /* conversion from double to int */
#endif
#define MINVALUE 0.8
#define MAXdim 200
#define PI 3.1415926535897932384
#define DEFAULTzbox 1e6
char prompt[]= "\n\
-rbox- generate various point distributions. Default is random in cube.\n\
\n\
args (any order, space separated): Version: 1998/8/10\n\
3000 number of random points in cube, lens, spiral, sphere or grid\n\
Bn bounding box coordinates, default %2.2g\n\
D3 dimension 3-d\n\
l generate a regular 3-d spiral\n\
r generate a regular polygon, ('r s Z1 G0.1' makes a cone)\n\
s generate cospherical points\n\
x generate random points in simplex, may use 'r' or 'Wn'\n\
y same as 'x', plus simplex\n\
\n\
Ln lens distribution of radius n. Also 's', 'r', 'G', 'W'.\n\
Mn,m,r lattice (Mesh) rotated by [n,-m,0], [m,n,0], [0,0,r], ... \n\
'27 M1,0,1' is {0,1,2} x {0,1,2} x {0,1,2}. Try 'M3,4 z'.\n\
On offset coordinates by n\n\
Pn,m,r add point [n,m,r] first, pads with 0\n\
W0.1 random distribution within 0.1 of the cube's or sphere's surface\n\
Z0.5 s random points in a 0.5 disk projected to a sphere\n\
Z0.5 s G0.6 same as Z0.5 within a 0.6 gap\n\
\n\
c add a unit cube to the output ('c G2.0' sets size)\n\
d add a unit diamond to the output ('d G2.0' sets size)\n\
h output as homogeneous coordinates for cdd\n\
n remove command line from the first line of output\n\
t use time as the random number seed (default is command line)\n\
tn use n as the random number seed\n\
z print integer coordinates, default 'Bn' is %2.2g\n\
";
/* ------------------------------ prototypes ----------------*/
int roundi( double a);
void out1( double a);
void out2n( double a, double b);
void out3n( double a, double b, double c);
int qh_rand( void);
void qh_srand( int seed);
/* ------------------------------ globals -------------------*/
FILE *fp;
int isinteger= 0;
double out_offset= 0.0;
/*--------------------------------------------
-rbox- main procedure of rbox application
*/
int main(int argc, char **argv) {
int i,j,k;
int gendim;
int cubesize, diamondsize, seed=0, count, apex;
int dim=3 , numpoints= 0, totpoints, addpoints=0;
int issphere=0, isaxis=0, iscdd= 0, islens= 0, isregular=0, iswidth=0, addcube=0;
int isgap=0, isspiral=0, NOcommand= 0, adddiamond=0, istime=0;
int isbox=0, issimplex=0, issimplex2=0, ismesh=0;
double width=0.0, gap=0.0, radius= 0.0;
double coord[MAXdim], offset, meshm=3.0, meshn=4.0, meshr=5.0;
double *simplex, *simplexp;
int nthroot, mult[MAXdim];
double norm, factor, randr, rangap, lensangle= 0, lensbase= 1;
double anglediff, angle, x, y, cube= 0.0, diamond= 0.0;
double box= qh_DEFAULTbox; /* scale all numbers before output */
double randmax= qh_RANDOMmax;
char command[200], *s, seedbuf[200];
time_t timedata;
#if __MWERKS__ && __POWERPC__
char inBuf[BUFSIZ], outBuf[BUFSIZ], errBuf[BUFSIZ];
SIOUXSettings.showstatusline= False;
SIOUXSettings.tabspaces= 1;
SIOUXSettings.rows= 40;
if (setvbuf (stdin, inBuf, _IOFBF, sizeof(inBuf)) < 0 /* w/o, SIOUX I/O is slow*/
|| setvbuf (stdout, outBuf, _IOFBF, sizeof(outBuf)) < 0
|| (stdout != stderr && setvbuf (stderr, errBuf, _IOFBF, sizeof(errBuf)) < 0))
fprintf ( stderr, "qhull internal warning (main): could not change stdio to fully buffered.\n");
argc= ccommand(&argv);
#endif
if (argc == 1) {
printf (prompt, box, DEFAULTzbox);
exit(1);
}
if ((s = strrchr( argv[0], '\\'))) /* Borland gives full path */
strcpy (command, s+1);
else
strcpy (command, argv[0]);
if ((s= strstr (command, ".EXE"))
|| (s= strstr (command, ".exe")))
*s= '\0';
/* ============= read flags =============== */
for (i=1; i < argc; i++) {
if (strlen (command) + strlen(argv[i]) + 1 < sizeof(command) ) {
strcat (command, " ");
strcat (command, argv[i]);
}
if (isdigit (argv[i][0])) {
numpoints= atoi (argv[i]);
continue;
}
if (argv[i][0] == '-')
(argv[i])++;
switch (argv[i][0]) {
case 'c':
addcube= 1;
if (i+1 < argc && argv[i+1][0] == 'G')
cube= (double) atof (&argv[++i][1]);
break;
case 'd':
adddiamond= 1;
if (i+1 < argc && argv[i+1][0] == 'G')
diamond= (double) atof (&argv[++i][1]);
break;
case 'h':
iscdd= 1;
break;
case 'l':
isspiral= 1;
break;
case 'n':
NOcommand= 1;
break;
case 'r':
isregular= 1;
break;
case 's':
issphere= 1;
break;
case 't':
istime= 1;
if (isdigit (argv[i][1]))
seed= atoi (&argv[i][1]);
else {
seed= time (&timedata);
sprintf (seedbuf, "%d", seed);
strcat (command, seedbuf);
}
break;
case 'x':
issimplex= 1;
break;
case 'y':
issimplex2= 1;
break;
case 'z':
isinteger= 1;
break;
case 'B':
box= (double) atof (&argv[i][1]);
isbox= 1;
break;
case 'D':
dim= atoi (&argv[i][1]);
if (dim < 1
|| dim > MAXdim) {
fprintf (stderr, "rbox error: dim %d too large or too small\n", dim);
exit (1);
}
break;
case 'G':
if (argv[i][1])
gap= (double) atof (&argv[i][1]);
else
gap= 0.5;
isgap= 1;
break;
case 'L':
if (argv[i][1])
radius= (double) atof (&argv[i][1]);
else
radius= 10;
islens= 1;
break;
case 'M':
ismesh= 1;
s= argv[i]+1;
if (*s)
meshn= strtod (s, &s);
if (*s == ',')
meshm= strtod (++s, &s);
else
meshm= 0.0;
if (*s == ',')
meshr= strtod (++s, &s);
else
meshr= sqrt (meshn*meshn + meshm*meshm);
if (*s) {
fprintf (stderr, "rbox warning: assuming 'M3,4,5' since mesh args are not integers or reals\n");
meshn= 3.0, meshm=4.0, meshr=5.0;
}
break;
case 'O':
out_offset= (double) atof (&argv[i][1]);
break;
case 'P':
addpoints++;
break;
case 'W':
width= (double) atof (&argv[i][1]);
iswidth= 1;
break;
case 'Z':
if (argv[i][1])
radius= (double) atof (&argv[i][1]);
else
radius= 1.0;
isaxis= 1;
break;
default:
fprintf (stderr, "rbox warning: unknown flag %s.\nExecute 'rbox' without arguments for documentation.\n", argv[i]);
}
}
/* ============= defaults, constants, and sizes =============== */
if (isinteger && !isbox)
box= DEFAULTzbox;
if (addcube) {
cubesize= floor(ldexp(1.0,dim)+0.5);
if (cube == 0.0)
cube= box;
}else
cubesize= 0;
if (adddiamond) {
diamondsize= 2*dim;
if (diamond == 0.0)
diamond= box;
}else
diamondsize= 0;
if (islens) {
if (isaxis) {
fprintf (stderr, "rbox error: can not combine 'Ln' with 'Zn'\n");
exit(1);
}
if (radius <= 1.0) {
fprintf (stderr, "rbox error: lens radius %.2g should be greater than 1.0\n",
radius);
exit(1);
}
lensangle= asin (1.0/radius);
lensbase= radius * cos (lensangle);
}
if (!numpoints) {
if (issimplex2)
; /* ok */
else if (isregular + issimplex + islens + issphere + isaxis + isspiral + iswidth + ismesh) {
fprintf (stderr, "rbox error: missing count\n");
exit(1);
}else if (adddiamond + addcube + addpoints)
; /* ok */
else {
numpoints= 50; /* ./rbox D4 is the test case */
issphere= 1;
}
}
if ((issimplex + islens + isspiral + ismesh > 1)
|| (issimplex + issphere + isspiral + ismesh > 1)) {
fprintf (stderr, "rbox error: can only specify one of 'l', 's', 'x', 'Ln', or 'Mn,m,r' ('Ln s' is ok).\n");
exit(1);
}
fp= stdout;
/* ============= print header with total points =============== */
if (issimplex || ismesh)
totpoints= numpoints;
else if (issimplex2)
totpoints= numpoints+dim+1;
else if (isregular) {
totpoints= numpoints;
if (dim == 2) {
if (islens)
totpoints += numpoints - 2;
}else if (dim == 3) {
if (islens)
totpoints += 2 * numpoints;
else if (isgap)
totpoints += 1 + numpoints;
else
totpoints += 2;
}
}else
totpoints= numpoints + isaxis;
totpoints += cubesize + diamondsize + addpoints;
if (iscdd)
fprintf(fp, "%s\nbegin\n %d %d %s\n",
NOcommand ? "" : command,
totpoints, dim+1,
isinteger ? "integer" : "real");
else if (NOcommand)
fprintf(fp, "%d\n%d\n", dim, totpoints);
else
fprintf(fp, "%d %s\n%d\n", dim, command, totpoints);
/* ============= seed randoms =============== */
if (istime == 0) {
for (s=command; *s; s++) {
if (issimplex2 && *s == 'y') /* make 'y' same seed as 'x' */
i= 'x';
else
i= *s;
seed= 11*seed + i;
}
} /* else, seed explicitly set to n or to time */
qh_RANDOMseed_(seed);
/* ============= explicit points =============== */
for (i=1; i < argc; i++) {
if (argv[i][0] == 'P') {
s= argv[i]+1;
count= 0;
if (iscdd)
out1( 1.0);
while (*s) {
out1( strtod (s, &s));
count++;
if (*s) {
if (*s++ != ',') {
fprintf (stderr, "rbox error: missing comma after coordinate in %s\n\n", argv[i]);
exit (1);
}
}
}
if (count < dim) {
for (k= dim-count; k--; )
out1( 0.0);
}else if (count > dim) {
fprintf (stderr, "rbox error: %d coordinates instead of %d coordinates in %s\n\n",
count, dim, argv[i]);
exit (1);
}
fprintf (fp, "\n");
}
}
/* ============= simplex distribution =============== */
if (issimplex+issimplex2) {
if (!(simplex= malloc( dim * (dim+1) * sizeof(double)))) {
fprintf (stderr, "insufficient memory for simplex\n");
exit(0);
}
simplexp= simplex;
if (isregular) {
for (i= 0; i<dim; i++) {
for (k= 0; k<dim; k++)
*(simplexp++)= i==k ? 1.0 : 0.0;
}
for (k= 0; k<dim; k++)
*(simplexp++)= -1.0;
}else {
for (i= 0; i<dim+1; i++) {
for (k= 0; k<dim; k++) {
randr= qh_RANDOMint;
*(simplexp++)= 2.0 * randr/randmax - 1.0;
}
}
}
if (issimplex2) {
simplexp= simplex;
for (i= 0; i<dim+1; i++) {
if (iscdd)
out1( 1.0);
for (k= 0; k<dim; k++)
out1( *(simplexp++) * box);
fprintf (fp, "\n");
⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -