#include "stdio.h"
#include "stdlib.h"
#include "math.h"

/*
	Create various data types for volume rendering experiments
*/

#define N 100

int main(int argc,char **argv)
{
	int i,j,k;
	int d,themin=1000,themax=-1000;
	double dx,dy,dz,r;
	unsigned char grid[N][N][N];
	int datatype;

	/* Check the arguments */
	if (argc < 2) {
		fprintf(stderr,"Usage: %s datatype\n",argv[0]);
		exit(-1);
	}
	if ((datatype = atoi(argv[1])) < 0 || datatype > 1) {
		fprintf(stderr,"Unexpected datatype\n");
		exit(-1);
	}

	/* Set everything to zero density */
	for (i=0;i<N;i++)
		for (j=0;j<N;j++)
			for (k=0;k<N;k++)
				grid[i][j][k] = 0;

	/* Create the data */
  	for (i=0;i<N;i++) {
     	for (j=0;j<N;j++) {
       	for (k=0;k<N;k++) {
				d = 0;
            dx = N/2 - i;
            dy = N/2 - j;
            dz = N/2 - k;
            r = sqrt(dx*dx + dy*dy + dz*dz);
				switch (datatype) {
				case 0: /* Decaying field */
					d = 200 * exp(-r/10);
					break;
            case 1: /* Constant field */
					if (r < N/2.5)
               	d = 10;
               break;
				}
				if (d < 0 || d > 255) {
					fprintf(stderr,"Out of range d of %d at (%d,%d,%d)\n",d,i,j,k);
					exit(-1);
				}
				grid[i][j][k] = d;
			}
		}
	}

	/* Calculate the range */
	for (i=0;i<N;i++) {
      for (j=0;j<N;j++) {
         for (k=0;k<N;k++) {
				if (grid[i][j][k] < themin)
					themin = grid[i][j][k];
            if (grid[i][j][k] > themax)
               themax = grid[i][j][k];
			}
		}
	}
	fprintf(stderr,"The range of densities is %d to %d\n",themin,themax);

	/* Write the data to standard out */
	for (i=0;i<N;i++) {
      for (j=0;j<N;j++) {
         for (k=0;k<N;k++) {
				putchar(grid[i][j][k]);
			}
		}
	}
}