From 6e994838e4dea3588c50cf5282067fcf79e60772 Mon Sep 17 00:00:00 2001 From: giraffedata Date: Sat, 24 Nov 2012 20:14:44 +0000 Subject: Add -randomseed git-svn-id: http://svn.code.sf.net/p/netpbm/code/trunk@1780 9d0c8265-081b-0410-96cb-a4ca84ce46f8 --- generator/pgmcrater.c | 595 +++++++++++++++++++++++++++----------------------- 1 file changed, 318 insertions(+), 277 deletions(-) diff --git a/generator/pgmcrater.c b/generator/pgmcrater.c index ec592381..25b7b817 100644 --- a/generator/pgmcrater.c +++ b/generator/pgmcrater.c @@ -1,34 +1,34 @@ /* - Fractal cratering + Fractal cratering - Designed and implemented in November of 1989 by: + Designed and implemented in November of 1989 by: - John Walker - Autodesk SA - Avenue des Champs-Montants 14b - CH-2074 MARIN - Switzerland - Usenet: kelvin@Autodesk.com - Fax: 038/33 88 15 - Voice: 038/33 76 33 + John Walker + Autodesk SA + Avenue des Champs-Montants 14b + CH-2074 MARIN + Switzerland + Usenet: kelvin@Autodesk.com + Fax: 038/33 88 15 + Voice: 038/33 76 33 The algorithm used to determine crater size is as described on pages 31 and 32 of: - Peitgen, H.-O., and Saupe, D. eds., The Science Of Fractal - Images, New York: Springer Verlag, 1988. + Peitgen, H.-O., and Saupe, D. eds., The Science Of Fractal + Images, New York: Springer Verlag, 1988. - The mathematical technique used to calculate crater radii that + The mathematical technique used to calculate crater radii that obey the proper area law distribution from a uniformly distributed pseudorandom sequence was developed by Rudy Rucker. - Permission to use, copy, modify, and distribute this software and - its documentation for any purpose and without fee is hereby + Permission to use, copy, modify, and distribute this software and + its documentation for any purpose and without fee is hereby granted, without any conditions or restrictions. This software is provided "as is" without express or implied warranty. - PLUGWARE! + PLUGWARE! If you like this kind of stuff, you may also enjoy "James Gleick's Chaos--The Software" for MS-DOS, available for $59.95 from your @@ -36,10 +36,10 @@ Series, 2320 Marinship Way, Sausalito, CA 94965, USA. Telephone: (800) 688-2344 toll-free or, outside the U.S. (415) 332-2344 Ext 4886. Fax: (415) 289-4718. "Chaos--The Software" includes a more - comprehensive fractal forgery generator which creates + comprehensive fractal forgery generator which creates three-dimensional landscapes as well as clouds and planets, plus five more modules which explore other aspects of Chaos. The user - guide of more than 200 pages includes an introduction by James + guide of more than 200 pages includes an introduction by James Gleick and detailed explanations by Rudy Rucker of the mathematics and algorithms used by each program. @@ -57,75 +57,195 @@ #include "pm_c_util.h" #include "pgm.h" #include "mallocvar.h" +#include "shhopt.h" + + +struct CmdlineInfo { + /* All the information the user supplied in the command line, + in a form easy for the program to use. + */ + unsigned int number; + unsigned int height; + unsigned int width; + float gamma; + unsigned int randomseed; + unsigned int randomseedSpec; +}; + + + +static void +parseCommandLine(int argc, const char ** const argv, + struct CmdlineInfo * const cmdlineP) { +/*---------------------------------------------------------------------------- + Note that the file spec array we return is stored in the storage that + was passed to us as the argv array. +-----------------------------------------------------------------------------*/ + optEntry * option_def; + /* Instructions to OptParseOptions3 on how to parse our options. + */ + optStruct3 opt; + unsigned int option_def_index; + + unsigned int numberSpec, heightSpec, widthSpec, gammaSpec; + + MALLOCARRAY_NOFAIL(option_def, 100); + + option_def_index = 0; /* incremented by OPTENT3 */ + OPTENT3(0, "number", OPT_UINT, &cmdlineP->number, + &numberSpec, 0); + OPTENT3(0, "height", OPT_UINT, &cmdlineP->height, + &heightSpec, 0); + OPTENT3(0, "ysize", OPT_UINT, &cmdlineP->height, + &heightSpec, 0); + OPTENT3(0, "width", OPT_UINT, &cmdlineP->width, + &widthSpec, 0); + OPTENT3(0, "xsize", OPT_UINT, &cmdlineP->width, + &widthSpec, 0); + OPTENT3(0, "gamma", OPT_FLOAT, &cmdlineP->gamma, + &gammaSpec, 0); + OPTENT3(0, "randomseed", OPT_UINT, &cmdlineP->randomseed, + &cmdlineP->randomseedSpec, 0); + + opt.opt_table = option_def; + opt.short_allowed = FALSE; /* We have no short (old-fashioned) options */ + opt.allowNegNum = FALSE; /* We may have parms that are negative numbers */ + + pm_optParseOptions3(&argc, (char **)argv, opt, sizeof(opt), 0); + /* Uses and sets argc, argv, and some of *cmdlineP and others. */ + + if (argc-1 > 0) + pm_error("There are no non-option arguments. You specified %u", + argc-1); + + if (!numberSpec) + cmdlineP->number = 50000; + + if (cmdlineP->number == 0) + pm_error("-number must be positive"); + + if (!heightSpec) + cmdlineP->height = 256; + + if (cmdlineP->height == 0) + pm_error("-height must be positive"); + + if (!widthSpec) + cmdlineP->width = 256; + + if (cmdlineP->width == 0) + pm_error("-width must be positive"); + + if (!gammaSpec) + cmdlineP->gamma = 1.0; + + if (cmdlineP->gamma <= 0.0) + pm_error("gamma correction must be greater than 0"); + + free(option_def); +} + -static void gencraters ARGS((void)); -static void initseed ARGS((void)); /* Definitions for obtaining random numbers. */ -#define Cast(low, high) ((low)+((high)-(low)) * ((rand() & 0x7FFF) / arand)) +/* Display parameters */ + +#define SCRX screenxsize /* Screen width */ +#define SCRY screenysize /* Screen height */ +#define SCRGAMMA 1.0 /* Display gamma */ + +#define RGBQuant 255 -/* Data types */ -typedef int Boolean; -#define FALSE 0 -#define TRUE 1 +static double const ImageGamma = 0.5; /* Inherent gamma of mapped image */ -#define V (void) +static double const arand = 32767.0; /* Random number parameters */ -/* Display parameters */ +static double const CdepthPower = 1.5; /* Crater depth power factor */ -#define SCRX screenxsize /* Screen width */ -#define SCRY screenysize /* Screen height */ -#define SCRGAMMA 1.0 /* Display gamma */ +static double DepthBias; /* sqrt(.5) */ -/* Local variables */ +static int const slopemin = -52; +static int const slopemax = 52; -#define ImageGamma 0.5 /* Inherent gamma of mapped image */ -static int screenxsize = 256; /* Screen X size */ -static int screenysize = 256; /* Screen Y size */ -static double dgamma = SCRGAMMA; /* Display gamma */ -static double arand = 32767.0; /* Random number parameters */ -static long ncraters = 50000L; /* Number of craters to generate */ -static double CdepthPower = 1.5; /* Crater depth power factor */ -static double DepthBias = 0.707107; /* Depth bias */ +static double const +Cast(double const low, + double const high) { + + return low + (high - low) * ((rand() & 0x7FFF) / arand); +} + + + +static int +modulo(int const t, + int const n) { -static int modulo(int t, int n) -{ int m; - assert(n>0); + + assert(n > 0); + m = t % n; - while (m<0) { - m+=n; + + while (m < 0) { + m += n; } + return m; } -/* INITSEED -- Generate initial random seed, if needed. */ -static void initseed() -{ - unsigned int i; - srand(pm_randseed()); - for (i = 0; i < 7; ++i) - V rand(); +#define Auxadr(x, y) &aux[modulo(y, screenysize)*screenxsize+modulo(x, screenxsize)] + + + +static void +generateScreenImage(const unsigned short * const aux, + unsigned int const screenxsize, + unsigned int const screenysize, + unsigned char * const slopemap) { + + unsigned int row; + gray * pixrow; + + pgm_writepgminit(stdout, screenxsize, screenysize, RGBQuant, FALSE); + pixrow = pgm_allocrow(screenxsize); + + for (row = 0; row < screenysize; ++row) { + unsigned int col; + + for (col = 0; col < screenxsize; ++col) { + int j; + j = *Auxadr(col+1, row) - *Auxadr(col, row); + j = MIN(MAX(slopemin, j), slopemax); + pixrow[col] = slopemap[j - slopemin]; + } + pgm_writepgmrow(stdout, pixrow, screenxsize, RGBQuant, FALSE); + } + pm_close(stdout); + pgm_freerow(pixrow); + } -/* GENCRATERS -- Generate cratered terrain. */ -static void gencraters() -{ - int i, j, x, y; - long l; - unsigned short *aux; - int slopemin = -52, slopemax = 52; -#define RGBQuant 255 - unsigned char *slopemap; /* Slope to pixel map */ - gray *pixels; /* Pixel vector */ -#define Auxadr(x, y) &aux[modulo(y, SCRY)*SCRX+modulo(x, SCRX)] +static void +gencraters(struct CmdlineInfo const cmdline) { +/*---------------------------------------------------------------------------- + Generate cratered terrain +-----------------------------------------------------------------------------*/ + unsigned int const screenxsize = cmdline.width; /* screen X size */ + unsigned int const screenysize = cmdline.height; /* screen Y size */ + double const dgamma = cmdline.gamma; /* display gamma */ + unsigned int const ncraters = cmdline.number; /* num craters to gen */ + + int i, j; + unsigned int l; + unsigned short * aux; + unsigned char * slopemap; /* Slope to pixel map */ /* Acquire the elevation array and initialize it to mean surface elevation. */ @@ -138,121 +258,114 @@ static void gencraters() initial elevation. */ for (i = 0; i < SCRY; i++) { - unsigned short *zax = aux + (((long) SCRX) * i); + unsigned short *zax = aux + (((long) SCRX) * i); - for (j = 0; j < SCRX; j++) { - *zax++ = 32767; - } + for (j = 0; j < SCRX; j++) { + *zax++ = 32767; + } } /* Every time we go around this loop we plop another crater - on the surface. */ + on the surface. */ for (l = 0; l < ncraters; l++) { - double g; - int cx = Cast(0.0, ((double) SCRX - 1)), - cy = Cast(0.0, ((double) SCRY - 1)), - gx, gy, x, y; - unsigned long amptot = 0, axelev; - unsigned int npatch = 0; - - - /* Phase 1. Compute the mean elevation of the impact - area. We assume the impact area is a - fraction of the total crater size. */ - - /* Thanks, Rudy, for this equation that maps the uniformly - distributed numbers from Cast into an area-law - distribution as observed on cratered bodies. */ - - g = sqrt(1 / (M_PI * (1 - Cast(0, 0.9999)))); - - /* If the crater is tiny, handle it specially. */ - -#if 0 - fprintf(stderr, "crater=%lu ", (unsigned long)l); - fprintf(stderr, "cx=%d ", cx); - fprintf(stderr, "cy=%d ", cy); - fprintf(stderr, "size=%g\n", g); -#endif - - if (g < 3) { - - /* Set pixel to the average of its Moore neighbourhood. */ - - for (y = cy - 1; y <= cy + 1; y++) { - for (x = cx - 1; x <= cx + 1; x++) { - amptot += *Auxadr(x, y); - npatch++; - } - } - axelev = amptot / npatch; - - /* Perturb the mean elevation by a small random factor. */ - - x = (g >= 1) ? ((rand() >> 8) & 3) - 1 : 0; - *Auxadr(cx, cy) = axelev + x; - - /* Jam repaint sizes to correct patch. */ - - gx = 1; - gy = 0; - - } else { - - /* Regular crater. Generate an impact feature of the - correct size and shape. */ - - /* Determine mean elevation around the impact area. */ - - gx = MAX(2, (g / 3)); - gy = MAX(2, g / 3); - - for (y = cy - gy; y <= cy + gy; y++) { - for (x = cx-gx; x <= cx + gx; x++) { - amptot += *Auxadr(x,y); - npatch++; - } - } - axelev = amptot / npatch; - - gy = MAX(2, g); - g = gy; - gx = MAX(2, g); - - for (y = cy - gy; y <= cy + gy; y++) { - double dy = (cy - y) / (double) gy, - dysq = dy * dy; - - for (x = cx - gx; x <= cx + gx; x++) { - double dx = ((cx - x) / (double) gx), - cd = (dx * dx) + dysq, - cd2 = cd * 2.25, - tcz = DepthBias - sqrt(fabs(1 - cd2)), - cz = MAX((cd2 > 1) ? 0.0 : -10, tcz), - roll, iroll; - unsigned short av; - - cz *= pow(g, CdepthPower); - if (dy == 0 && dx == 0 && ((int) cz) == 0) { - cz = cz < 0 ? -1 : 1; - } - -#define rollmin 0.9 - roll = (((1 / (1 - MIN(rollmin, cd))) / - (1 / (1 - rollmin))) - (1 - rollmin)) / rollmin; - iroll = 1 - roll; - - av = (axelev + cz) * iroll + (*Auxadr(x,y) + cz) * roll; - av = MAX(1000, MIN(64000, av)); - *Auxadr(x,y) = av; - } - } - } - if ((l % 5000) == 4999) { - pm_message( "%ld craters generated of %ld (%ld%% done)", - l + 1, ncraters, ((l + 1) * 100) / ncraters); - } + double g; + int cx = Cast(0.0, ((double) SCRX - 1)), + cy = Cast(0.0, ((double) SCRY - 1)), + gx, gy, x, y; + unsigned int amptot = 0, axelev; + unsigned int npatch = 0; + + + /* Phase 1. Compute the mean elevation of the impact + area. We assume the impact area is a + fraction of the total crater size. */ + + /* Thanks, Rudy, for this equation that maps the uniformly + distributed numbers from Cast into an area-law + distribution as observed on cratered bodies. */ + + g = sqrt(1 / (M_PI * (1 - Cast(0, 0.9999)))); + + /* If the crater is tiny, handle it specially. */ + + if (g < 3) { + + /* Set pixel to the average of its Moore neighbourhood. */ + + for (y = cy - 1; y <= cy + 1; y++) { + for (x = cx - 1; x <= cx + 1; x++) { + amptot += *Auxadr(x, y); + npatch++; + } + } + axelev = amptot / npatch; + + /* Perturb the mean elevation by a small random factor. */ + + x = (g >= 1) ? ((rand() >> 8) & 3) - 1 : 0; + *Auxadr(cx, cy) = axelev + x; + + /* Jam repaint sizes to correct patch. */ + + gx = 1; + gy = 0; + + } else { + + /* Regular crater. Generate an impact feature of the + correct size and shape. */ + + /* Determine mean elevation around the impact area. */ + + gx = MAX(2, (g / 3)); + gy = MAX(2, g / 3); + + for (y = cy - gy; y <= cy + gy; y++) { + for (x = cx-gx; x <= cx + gx; x++) { + amptot += *Auxadr(x,y); + npatch++; + } + } + axelev = amptot / npatch; + + gy = MAX(2, g); + g = gy; + gx = MAX(2, g); + + for (y = cy - gy; y <= cy + gy; y++) { + double dy = (cy - y) / (double) gy, + dysq = dy * dy; + + for (x = cx - gx; x <= cx + gx; x++) { + double dx = ((cx - x) / (double) gx), + cd = (dx * dx) + dysq, + cd2 = cd * 2.25, + tcz = DepthBias - sqrt(fabs(1 - cd2)), + cz = MAX((cd2 > 1) ? 0.0 : -10, tcz), + roll, iroll; + unsigned short av; + + cz *= pow(g, CdepthPower); + if (dy == 0 && dx == 0 && ((int) cz) == 0) { + cz = cz < 0 ? -1 : 1; + } + +#define rollmin 0.9 + roll = (((1 / (1 - MIN(rollmin, cd))) / + (1 / (1 - rollmin))) - (1 - rollmin)) / rollmin; + iroll = 1 - roll; + + av = (axelev + cz) * iroll + (*Auxadr(x,y) + cz) * roll; + av = MAX(1000, MIN(64000, av)); + *Auxadr(x,y) = av; + } + } + } + if ((l % 5000) == 4999) { + pm_message( "%u craters generated of %u (%u%% done)", + l + 1, ncraters, ((l + 1) * 100) / ncraters); + } } i = MAX((slopemax - slopemin) + 1, 1); @@ -263,120 +376,48 @@ static void gencraters() for (i = slopemin; i <= slopemax; i++) { /* Confused? OK, we're using the left-to-right slope to - calculate a shade based on the sine of the angle with - respect to the vertical (light incident from the left). - Then, with one exponentiation, we account for both the - inherent gamma of the image (ad-hoc), and the - user-specified display gamma, using the identity: - - (x^y)^z = (x^(y*z)) */ - - slopemap[i - slopemin] = i > 0 ? - (128 + 127.0 * - pow(sin((M_PI / 2) * i / slopemax), - dgamma * ImageGamma)) : - (128 - 127.0 * - pow(sin((M_PI / 2) * i / slopemin), - dgamma * ImageGamma)); + calculate a shade based on the sine of the angle with + respect to the vertical (light incident from the left). + Then, with one exponentiation, we account for both the + inherent gamma of the image (ad-hoc), and the + user-specified display gamma, using the identity: + + (x^y)^z = (x^(y*z)) */ + + slopemap[i - slopemin] = i > 0 ? + (128 + 127.0 * + pow(sin((M_PI / 2) * i / slopemax), + dgamma * ImageGamma)) : + (128 - 127.0 * + pow(sin((M_PI / 2) * i / slopemin), + dgamma * ImageGamma)); } - /* Generate the screen image. */ - - pgm_writepgminit(stdout, SCRX, SCRY, RGBQuant, FALSE); - pixels = pgm_allocrow(SCRX); + generateScreenImage(aux, screenxsize, screenysize, slopemap); - for (y = 0; y < SCRY; y++) { - gray *pix = pixels; - - for (x = 0; x < SCRX; x++) { - int j = *Auxadr(x+1, y) - *Auxadr(x, y); - j = MIN(MAX(slopemin, j), slopemax); - *pix++ = slopemap[j - slopemin]; - } - pgm_writepgmrow(stdout, pixels, SCRX, RGBQuant, FALSE); - } - pm_close(stdout); - pgm_freerow(pixels); - -#undef Auxadr -#undef Scradr free((char *) slopemap); free((char *) aux); } -/* MAIN -- Main program. */ - -int main(argc, argv) - int argc; - char *argv[]; -{ - int i; - Boolean gammaspec = FALSE, numspec = FALSE, - widspec = FALSE, hgtspec = FALSE; - const char * const usage = "[-number ] [-width|-xsize ]\n\ - [-height|-ysize ] [-gamma ]"; - - DepthBias = sqrt(0.5); /* Get exact value for depth bias */ - - - pgm_init(&argc, argv); - - i = 1; - while ((i < argc) && (argv[i][0] == '-') && (argv[i][1] != '\0')) { - if (pm_keymatch(argv[i], "-gamma", 2)) { - if (gammaspec) { - pm_error("already specified gamma correction"); - } - i++; - if ((i == argc) || (sscanf(argv[i], "%lf", &dgamma) != 1)) - pm_usage(usage); - if (dgamma <= 0.0) { - pm_error("gamma correction must be greater than 0"); - } - gammaspec = TRUE; - } else if (pm_keymatch(argv[i], "-number", 2)) { - if (numspec) { - pm_error("already specified number of craters"); - } - i++; - if ((i == argc) || (sscanf(argv[i], "%ld", &ncraters) != 1)) - pm_usage(usage); - if (ncraters <= 0) { - pm_error("number of craters must be greater than 0!"); - } - numspec = TRUE; - } else if (pm_keymatch(argv[i], "-xsize", 2) || - pm_keymatch(argv[i], "-width", 2)) { - if (widspec) { - pm_error("already specified a width/xsize"); - } - i++; - if ((i == argc) || (sscanf(argv[i], "%d", &screenxsize) != 1)) - pm_usage(usage); - if (screenxsize <= 0) { - pm_error("screen width must be greater than 0"); - } - widspec = TRUE; - } else if (pm_keymatch(argv[i], "-ysize", 2) || - pm_keymatch(argv[i], "-height", 2)) { - if (hgtspec) { - pm_error("already specified a height/ysize"); - } - i++; - if ((i == argc) || (sscanf(argv[i], "%d", &screenysize) != 1)) - pm_usage(usage); - if (screenxsize <= 0) { - pm_error("screen height must be greater than 0"); - } - hgtspec = TRUE; - } else { - pm_usage(usage); - } - i++; - } - initseed(); - gencraters(); + +int +main(int argc, const char ** argv) { + + struct CmdlineInfo cmdline; + + pm_proginit(&argc, argv); + + parseCommandLine(argc, argv, &cmdline); + + srand(cmdline.randomseedSpec ? cmdline.randomseed : pm_randseed()); + + DepthBias = sqrt(0.5); + + gencraters(cmdline); exit(0); } + + + -- cgit 1.4.1