home *** CD-ROM | disk | FTP | other *** search
Text File | 1991-12-14 | 55.2 KB | 1,693 lines |
- Newsgroups: comp.sources.misc
- From: jef@well.sf.ca.us (Jef Poskanzer)
- Subject: v26i108: pbmplus - Extended Portable Bitmap Toolkit, Patch10dec91, Part03/05
- Message-ID: <1991Dec15.014640.20443@sparky.imd.sterling.com>
- X-Md4-Signature: 9a3f6109d610901d81f976509fa1525e
- Date: Sun, 15 Dec 1991 01:46:40 GMT
- Approved: kent@sparky.imd.sterling.com
-
- Submitted-by: jef@well.sf.ca.us (Jef Poskanzer)
- Posting-number: Volume 26, Issue 108
- Archive-name: pbmplus/patch10dec91/part03
- Environment: UNIX
- Patch-To: pbmplus: Volume 23, Issue 36-59
-
- #!/bin/sh
- # do not concatenate these parts, unpack them in order with /bin/sh
- # file pgm/pgmcrater.c continued
- #
- if test ! -r _shar_seq_.tmp; then
- echo 'Please unpack part 1 first!'
- exit 1
- fi
- (read Scheck
- if test "$Scheck" != 3; then
- echo Please unpack part "$Scheck" next!
- exit 1
- else
- exit 0
- fi
- ) < _shar_seq_.tmp || exit 1
- if test ! -f _shar_wnt_.tmp; then
- echo 'x - still skipping pgm/pgmcrater.c'
- else
- echo 'x - continuing file pgm/pgmcrater.c'
- sed 's/^X//' << 'SHAR_EOF' >> 'pgm/pgmcrater.c' &&
- X pages 31 and 32 of:
- X
- X Peitgen, H.-O., and Saupe, D. eds., The Science Of Fractal
- X Images, New York: Springer Verlag, 1988.
- X
- X The mathematical technique used to calculate crater radii that
- X obey the proper area law distribution from a uniformly distributed
- X pseudorandom sequence was developed by Rudy Rucker.
- X
- X Permission to use, copy, modify, and distribute this software and
- X its documentation for any purpose and without fee is hereby
- X granted, without any conditions or restrictions. This software is
- X provided "as is" without express or implied warranty.
- X
- X PLUGWARE!
- X
- X If you like this kind of stuff, you may also enjoy "James Gleick's
- X Chaos--The Software" for MS-DOS, available for $59.95 from your
- X local software store or directly from Autodesk, Inc., Attn: Science
- X Series, 2320 Marinship Way, Sausalito, CA 94965, USA. Telephone:
- X (800) 688-2344 toll-free or, outside the U.S. (415) 332-2344 Ext
- X 4886. Fax: (415) 289-4718. "Chaos--The Software" includes a more
- X comprehensive fractal forgery generator which creates
- X three-dimensional landscapes as well as clouds and planets, plus
- X five more modules which explore other aspects of Chaos. The user
- X guide of more than 200 pages includes an introduction by James
- X Gleick and detailed explanations by Rudy Rucker of the mathematics
- X and algorithms used by each program.
- X
- */
- X
- #include "pgm.h"
- #include <math.h>
- X
- /* Definitions for obtaining random numbers. */
- X
- #define Cast(low, high) ((low)+((high)-(low)) * ((random() & 0x7FFF) / arand))
- X
- /* Data types */
- X
- typedef int Boolean;
- #define FALSE 0
- #define TRUE 1
- X
- #define V (void)
- X
- /* Display parameters */
- X
- #define SCRX screenxsize /* Screen width */
- #define SCRY screenysize /* Screen height */
- #define SCRGAMMA 1.0 /* Display gamma */
- X
- /* Local variables */
- X
- #define ImageGamma 0.5 /* Inherent gamma of mapped image */
- X
- 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 */
- X
- /* INITSEED -- Generate initial random seed, if needed. */
- X
- static void initseed()
- {
- X int i;
- X
- X i = time((long *) 0) * 0xF37C;
- X srandom(i);
- X for (i = 0; i < 7; i++) {
- X V random();
- X }
- }
- X
- /* GENCRATERS -- Generate cratered terrain. */
- X
- static void gencraters()
- {
- X int i, j, x, y;
- X long l;
- X unsigned short *aux;
- X int slopemin = -52, slopemax = 52;
- #define RGBQuant 255
- X unsigned char *slopemap; /* Slope to pixel map */
- X gray *pixels; /* Pixel vector */
- X
- #define Auxadr(x, y) ((unsigned short *) (aux + ((((y)) * SCRX) + (x))))
- X
- X /* Acquire the elevation array and initialise it to mean
- X surface elevation. */
- X
- X aux = (unsigned short *) malloc(SCRX * SCRY * sizeof(short));
- X if (aux == (unsigned short *) 0) {
- X pm_error("out of memory allocating elevation array");
- X }
- X
- X /* Acquire the elevation buffer and initialise to mean
- X initial elevation. */
- X
- X for (i = 0; i < SCRY; i++) {
- X unsigned short *zax = aux + (((long) SCRX) * i);
- X
- X for (j = 0; j < SCRX; j++) {
- X *zax++ = 32767;
- X }
- X }
- X
- X /* Every time we go around this loop we plop another crater
- X on the surface. */
- X
- X for (l = 0; l < ncraters; l++) {
- X double g;
- X int cx = Cast(0.0, ((double) SCRX - 1)),
- X cy = Cast(0.0, ((double) SCRY - 1)),
- X gx, gy, x, y;
- X unsigned long amptot = 0, axelev;
- X unsigned int npatch = 0;
- X
- X
- X /* Phase 1. Compute the mean elevation of the impact
- X area. We assume the impact area is a
- X fraction of the total crater size. */
- X
- X /* Thanks, Rudy, for this equation that maps the uniformly
- X distributed numbers from Cast into an area-law
- X distribution as observed on cratered bodies. */
- X
- X g = sqrt(1 / (M_PI * (1 - Cast(0, 0.9999))));
- X
- X /* If the crater is tiny, handle it specially. */
- X
- X if (g < 3) {
- X
- X /* Set pixel to the average of its Moore neighbourhood. */
- X
- X for (y = max(0, cy - 1); y <= min(SCRY - 1, cy + 1); y++) {
- X int sx = max(0, cx - 1);
- X unsigned short *a = Auxadr(sx, y);
- X
- X for (x = sx; x <= min(SCRX - 1, cx + 1); x++) {
- X amptot += *a++;
- X npatch++;
- X }
- X }
- X axelev = amptot / npatch;
- X
- X /* Perturb the mean elevation by a small random factor. */
- X
- X x = (g >= 1) ? ((random() >> 8) & 3) - 1 : 0;
- X *Auxadr(cx, cy) = axelev + x;
- X
- X /* Jam repaint sizes to correct patch. */
- X
- X gx = 1;
- X gy = 0;
- X
- X } else {
- X
- X /* Regular crater. Generate an impact feature of the
- X correct size and shape. */
- X
- X /* Determine mean elevation around the impact area. */
- X
- X gx = max(2, (g / 3));
- X gy = max(2, g / 3);
- X
- X for (y = max(0, cy - gy); y <= min(SCRY - 1, cy + gy); y++) {
- X int sx = max(0, cx - gx);
- X unsigned short *a = Auxadr(sx, y);
- X
- X for (x = sx; x <= min(SCRX - 1, cx + gx); x++) {
- X amptot += *a++;
- X npatch++;
- X }
- X }
- X axelev = amptot / npatch;
- X
- X gy = max(2, g);
- X g = gy;
- X gx = max(2, g);
- X
- X for (y = max(0, cy - gy); y <= min(SCRY - 1, cy + gy); y++) {
- X int sx = max(0, cx - gx);
- X unsigned short *ax = Auxadr(sx, y);
- X double dy = (cy - y) / (double) gy,
- X dysq = dy * dy;
- X
- X for (x = sx; x <= min(SCRX - 1, cx + gx); x++) {
- X double dx = ((cx - x) / (double) gx),
- X cd = (dx * dx) + dysq,
- X cd2 = cd * 2.25,
- X tcz = DepthBias - sqrt(fabs(1 - cd2)),
- X cz = max((cd2 > 1) ? 0.0 : -10, tcz),
- X roll, iroll;
- X unsigned short av;
- X
- X cz *= pow(g, CdepthPower);
- X if (dy == 0 && dx == 0 && ((int) cz) == 0) {
- X cz = cz < 0 ? -1 : 1;
- X }
- X
- #define rollmin 0.9
- X roll = (((1 / (1 - min(rollmin, cd))) /
- X (1 / (1 - rollmin))) - (1 - rollmin)) / rollmin;
- X iroll = 1 - roll;
- X
- X av = (axelev + cz) * iroll + (*ax + cz) * roll;
- X av = max(1000, min(64000, av));
- X *ax++ = av;
- X }
- X }
- X }
- X if ((l % 5000) == 4999) {
- X pm_message( "%ld craters generated of %ld (%ld%% done)",
- X l + 1, ncraters, ((l + 1) * 100) / ncraters);
- X }
- X }
- X
- X i = max((slopemax - slopemin) + 1, 1);
- X slopemap = (unsigned char *) malloc(i * sizeof(unsigned char));
- X if (slopemap == (unsigned char *) 0) {
- X pm_error("out of memory allocating slope map");
- X }
- X for (i = slopemin; i <= slopemax; i++) {
- X
- X /* Confused? OK, we're using the left-to-right slope to
- X calculate a shade based on the sine of the angle with
- X respect to the vertical (light incident from the left).
- X Then, with one exponentiation, we account for both the
- X inherent gamma of the image (ad-hoc), and the
- X user-specified display gamma, using the identity:
- X
- X (x^y)^z = (x^(y*z)) */
- X
- X slopemap[i - slopemin] = i > 0 ?
- X (128 + 127.0 *
- X pow(sin((M_PI / 2) * i / slopemax),
- X dgamma * ImageGamma)) :
- X (128 - 127.0 *
- X pow(sin((M_PI / 2) * i / slopemin),
- X dgamma * ImageGamma));
- X }
- X
- X /* Generate the screen image. */
- X
- X pgm_writepgminit(stdout, SCRX, SCRY, RGBQuant, FALSE);
- X pixels = pgm_allocrow(SCRX);
- X
- X for (y = 0; y < SCRY; y++) {
- X unsigned short *ax = Auxadr(0, y);
- X gray *pix = pixels;
- X
- X for (x = 0; x < SCRX - 1; x++) {
- X int j = ax[1] - ax[0];
- X j = min(max(slopemin, j), slopemax);
- X *pix++ = slopemap[j - slopemin];
- X ax++;
- X }
- X pgm_writepgmrow(stdout, pixels, SCRX, RGBQuant, FALSE);
- X }
- X pm_close(stdout);
- X pgm_freerow(pixels);
- X
- #undef Auxadr
- #undef Scradr
- X free((char *) slopemap);
- X free((char *) aux);
- }
- X
- /* MAIN -- Main program. */
- X
- int main(argc, argv)
- X int argc;
- X char *argv[];
- {
- X int i;
- X Boolean gammaspec = FALSE, numspec = FALSE,
- X widspec = FALSE, hgtspec = FALSE;
- X char *usage = "[-number <n>] [-width|-xsize <w>]\n\
- X [-height|-ysize <h>] [-gamma <f>]";
- X
- X DepthBias = sqrt(0.5); /* Get exact value for depth bias */
- X
- X pgm_init(&argc, argv);
- X
- X i = 1;
- X while ((i < argc) && (argv[i][0] == '-') && (argv[i][1] != '\0')) {
- X if (pm_keymatch(argv[i], "-gamma", 2)) {
- X if (gammaspec) {
- X pm_error("already specified gamma correction");
- X }
- X i++;
- X if ((i == argc) || (sscanf(argv[i], "%lf", &dgamma) != 1))
- X pm_usage(usage);
- X if (dgamma <= 0.0) {
- X pm_error("gamma correction must be greater than 0");
- X }
- X gammaspec = TRUE;
- X } else if (pm_keymatch(argv[i], "-number", 2)) {
- X if (numspec) {
- X pm_error("already specified number of craters");
- X }
- X i++;
- X if ((i == argc) || (sscanf(argv[i], "%ld", &ncraters) != 1))
- X pm_usage(usage);
- X if (ncraters <= 0) {
- X pm_error("number of craters must be greater than 0!");
- X }
- X numspec = TRUE;
- X } else if (pm_keymatch(argv[i], "-xsize", 2) ||
- X pm_keymatch(argv[i], "-width", 2)) {
- X if (widspec) {
- X pm_error("already specified a width/xsize");
- X }
- X i++;
- X if ((i == argc) || (sscanf(argv[i], "%d", &screenxsize) != 1))
- X pm_usage(usage);
- X if (screenxsize <= 0) {
- X pm_error("screen width must be greater than 0");
- X }
- X widspec = TRUE;
- X } else if (pm_keymatch(argv[i], "-ysize", 2) ||
- X pm_keymatch(argv[i], "-height", 2)) {
- X if (hgtspec) {
- X pm_error("already specified a height/ysize");
- X }
- X i++;
- X if ((i == argc) || (sscanf(argv[i], "%d", &screenysize) != 1))
- X pm_usage(usage);
- X if (screenxsize <= 0) {
- X pm_error("screen height must be greater than 0");
- X }
- X hgtspec = TRUE;
- X } else {
- X pm_usage(usage);
- X }
- X i++;
- X }
- X
- X initseed();
- X gencraters();
- X
- X exit(0);
- }
- SHAR_EOF
- echo 'File pgm/pgmcrater.c is complete' &&
- chmod 0664 pgm/pgmcrater.c ||
- echo 'restore of pgm/pgmcrater.c failed'
- Wc_c="`wc -c < 'pgm/pgmcrater.c'`"
- test 10318 -eq "$Wc_c" ||
- echo 'pgm/pgmcrater.c: original size 10318, current size' "$Wc_c"
- rm -f _shar_wnt_.tmp
- fi
- # ============= pnm/Imakefile ==============
- if test ! -d 'pnm'; then
- echo 'x - creating directory pnm'
- mkdir 'pnm'
- fi
- if test -f 'pnm/Imakefile' -a X"$1" != X"-c"; then
- echo 'x - skipping pnm/Imakefile (File already exists)'
- rm -f _shar_wnt_.tmp
- else
- > _shar_wnt_.tmp
- echo 'x - extracting pnm/Imakefile (Text)'
- sed 's/^X//' << 'SHAR_EOF' > 'pnm/Imakefile' &&
- /* Imakefile for pnm tools
- X *
- X * Copyright (C) 1991 Rainer Klute
- X *
- X * Permission to use, copy, modify, distribute, and sell this software and
- X * its documentation for any purpose is hereby granted without fee, provided
- X * that the above copyright notice appear in all copies and that both that
- X * copyright notice and this permission notice appear in supporting
- X * documentation, and that the copyright holder's name not be used in
- X * advertising or publicity pertaining to distribution of the software
- X * without specific, written prior permission. The copyright holder makes
- X * no representations about the suitability of this software for any
- X * purpose. It is provided "as is" without express or implied warranty.
- X */
- X
- #define LibPnm libpnm.a
- #define DepLibPnm LibPnm
- #include <../Pbmplus.tmpl>
- X
- #if BuildLibTiff
- X CURRENTLIBS = $(LIBTIFF) $(LIBPNM) $(LIBPPM) $(LIBPGM) $(LIBPBM)
- CURRENTDEPLIBS = $(DEPLIBTIFF) $(DEPLIBPNM) $(DEPLIBPPM) $(DEPLIBPGM) $(DEPLIBPBM)
- X INCLUDES = -I.. -I$(PBMDIR) -I$(PGMDIR) -I$(PPMDIR) -I$(TIFFDIR)
- X DEFINES = -DLIBTIFF
- X MERGE = pnmmerge
- X TIFFMAN1 = tifftopnm.1 pnmtotiff.1
- X TIFFSRCS = tifftopnm.c pnmtotiff.c
- X TIFFOBJS = tifftopnm.o pnmtotiff.o
- X TIFFBINS = tifftopnm pnmtotiff
- #else
- X CURRENTLIBS = $(LIBPNM) $(LIBPPM) $(LIBPGM) $(LIBPBM)
- CURRENTDEPLIBS = $(DEPLIBPNM) $(DEPLIBPPM) $(DEPLIBPGM) $(DEPLIBPBM)
- X INCLUDES = -I.. -I$(PBMDIR) -I$(PGMDIR) -I$(PPMDIR)
- X MERGE = pnmmerge
- X TIFFMAN1 =
- X TIFFSRCS =
- X TIFFOBJS =
- X TIFFBINS =
- #endif
- X
- X MAN1 = pnmarith.1 pnmcat.1 pnmconvol.1 pnmcrop.1 pnmcut.1 \
- X pnmdepth.1 pnmenlarge.1 pnmfile.1 pnmflip.1 pnminvert.1 \
- X pnmnoraw.1 pnmpaste.1 pnmscale.1 pnmtile.1 pnmtops.1 \
- X pnmtorast.1 pnmtoxwd.1 rasttopnm.1 xwdtopnm.1 \
- X pnmgamma.1 pnmrotate.1 pnmshear.1 \
- X anytopnm.1 pnmindex.1 pnmmargin.1 pnmsmooth.1 \
- X $(TIFFMAN1)
- X MAN3 = libpnm.3
- X MAN5 = pnm.5
- X
- X SRCS = pnmarith.c pnmcat.c pnmconvol.c pnmcrop.c pnmcut.c \
- X pnmdepth.c pnmenlarge.c pnmfile.c pnmflip.c pnminvert.c \
- X pnmnoraw.c pnmpaste.c pnmscale.c pnmtile.c pnmtops.c \
- X pnmtorast.c pnmtoxwd.c rasttopnm.c xwdtopnm.c \
- X pnmgamma.c pnmrotate.c pnmshear.c \
- X $(TIFFSRCS)
- X
- X OBJS = pnmarith.o pnmcat.o pnmconvol.o pnmcrop.o pnmcut.o \
- X pnmdepth.o pnmenlarge.o pnmfile.o pnmflip.o pnminvert.o \
- X pnmnoraw.o pnmpaste.o pnmscale.o pnmtile.o pnmtops.o \
- X pnmtorast.o pnmtoxwd.o rasttopnm.o xwdtopnm.o \
- X pnmgamma.o pnmrotate.o pnmshear.o \
- X $(TIFFOBJS)
- X
- X BINS = pnmarith pnmcat pnmconvol pnmcrop pnmcut \
- X pnmdepth pnmenlarge pnmfile pnmflip pnminvert \
- X pnmnoraw pnmpaste pnmscale pnmtile pnmtops \
- X pnmtorast pnmtoxwd rasttopnm xwdtopnm \
- X pnmgamma pnmrotate pnmshear \
- X $(TIFFBINS)
- X
- includes:: anytopnm.script pnmindex.script pnmmargin.script pnmsmooth.script
- X
- anytopnm.script:
- X $(LN) anytopnm anytopnm.script
- X
- pnmindex.script:
- X $(LN) pnmindex pnmindex.script
- X
- pnmmargin.script:
- X $(LN) pnmmargin pnmmargin.script
- X
- pnmsmooth.script:
- X $(LN) pnmsmooth pnmsmooth.script
- X
- AllTarget($(LIBPNM) $(BINS))
- X
- DependTarget()
- X
- NormalPbmplusProgramTarget(pnmarith)
- NormalPbmplusProgramTarget(pnmcat)
- NormalPbmplusProgramTarget(pnmconvol)
- NormalPbmplusProgramTarget(pnmcrop)
- NormalPbmplusProgramTarget(pnmcut)
- NormalPbmplusProgramTarget(pnmdepth)
- NormalPbmplusProgramTarget(pnmenlarge)
- NormalPbmplusProgramTarget(pnmfile)
- NormalPbmplusProgramTarget(pnmflip)
- NormalPbmplusProgramTarget(pnminvert)
- NormalPbmplusProgramTarget(pnmnoraw)
- NormalPbmplusProgramTarget(pnmpaste)
- NormalPbmplusProgramTarget(pnmscale)
- NormalPbmplusProgramTarget(pnmtile)
- NormalPbmplusProgramTarget(pnmtops)
- NormalPbmplusProgramTarget(pnmtorast)
- NormalPbmplusProgramTarget(pnmtoxwd)
- NormalPbmplusProgramTarget(rasttopnm)
- NormalPbmplusProgramTarget(xwdtopnm)
- NormalPbmplusMathProgramTarget(pnmgamma)
- NormalPbmplusMathProgramTarget(pnmrotate)
- NormalPbmplusMathProgramTarget(pnmshear)
- #if BuildLibTiff
- NormalPbmplusProgramTarget(tifftopnm)
- NormalPbmplusProgramTarget(pnmtotiff)
- #endif
- X
- NormalLibraryObjectRule()
- NormalLibraryTarget(pnm,libpnm1.o libpnm2.o libpnm3.o libpnm4.o)
- X
- #if InstallMerged
- NormalProgramTarget($(MERGE),$(MERGE).o $(OBJS),$(CURRENTDEPLIBS),$(CURRENTLIBS),-lm)
- #if InstallBinaries
- InstallProgram($(MERGE),$(PBMPLUSDIR)$(PBMPLUSBINDIR))
- #endif
- #endif
- X
- #if InstallBinaries
- InstallPbmplusPrograms($(BINS),$(PBMPLUSDIR)$(PBMPLUSBINDIR),$(INSTPGMFLAGS))
- InstallScript(anytopnm,$(PBMPLUSDIR)$(PBMPLUSBINDIR))
- InstallScript(pnmindex,$(PBMPLUSDIR)$(PBMPLUSBINDIR))
- InstallScript(pnmmargin,$(PBMPLUSDIR)$(PBMPLUSBINDIR))
- InstallScript(pnmsmooth,$(PBMPLUSDIR)$(PBMPLUSBINDIR))
- #endif
- X
- #if InstallManuals
- InstallMultipleMan($(MAN1),$(PBMPLUSDIR)$(PBMPLUSMANDIR)/man1)
- InstallMultipleMan($(MAN3),$(PBMPLUSDIR)$(PBMPLUSMANDIR)/man3)
- InstallMultipleMan($(MAN5),$(PBMPLUSDIR)$(PBMPLUSMANDIR)/man5)
- #endif
- X
- #if InstallLibraries
- InstallLibrary(pnm,$(PBMPLUSDIR)$(PBMPLUSLIBDIR))
- #endif
- X
- #if InstallIncludes
- InstallMultipleFlags(pnm.h,$(PBMPLUSDIR)$(PBMPLUSINCDIR),$(INSTINCFLAGS))
- #endif
- SHAR_EOF
- chmod 0664 pnm/Imakefile ||
- echo 'restore of pnm/Imakefile failed'
- Wc_c="`wc -c < 'pnm/Imakefile'`"
- test 5158 -eq "$Wc_c" ||
- echo 'pnm/Imakefile: original size 5158, current size' "$Wc_c"
- rm -f _shar_wnt_.tmp
- fi
- # ============= ppm/Imakefile ==============
- if test ! -d 'ppm'; then
- echo 'x - creating directory ppm'
- mkdir 'ppm'
- fi
- if test -f 'ppm/Imakefile' -a X"$1" != X"-c"; then
- echo 'x - skipping ppm/Imakefile (File already exists)'
- rm -f _shar_wnt_.tmp
- else
- > _shar_wnt_.tmp
- echo 'x - extracting ppm/Imakefile (Text)'
- sed 's/^X//' << 'SHAR_EOF' > 'ppm/Imakefile' &&
- /* Imakefile for ppm tools
- X *
- X * Copyright (C) 1991 Rainer Klute
- X *
- X * Permission to use, copy, modify, distribute, and sell this software and
- X * its documentation for any purpose is hereby granted without fee, provided
- X * that the above copyright notice appear in all copies and that both that
- X * copyright notice and this permission notice appear in supporting
- X * documentation, and that the copyright holder's name not be used in
- X * advertising or publicity pertaining to distribution of the software
- X * without specific, written prior permission. The copyright holder makes
- X * no representations about the suitability of this software for any
- X * purpose. It is provided "as is" without express or implied warranty.
- X */
- X
- #define LibPpm libppm.a
- #define DepLibPpm LibPpm
- #include <../Pbmplus.tmpl>
- X
- X CURRENTLIBS = $(LIBPPM) $(LIBPGM) $(LIBPBM)
- CURRENTDEPLIBS = $(DEPLIBPPM) $(DEPLIBPGM) $(DEPLIBPBM)
- X INCLUDES = -I.. -I$(PBMDIR) -I$(PGMDIR)
- X DEFINES = -DRGB_DB=\"DefaultRGBDatabase\"
- X MERGE = ppmmerge
- X
- X MAN1 = giftoppm.1 gouldtoppm.1 ilbmtoppm.1 imgtoppm.1 mtvtoppm.1 \
- X pcxtoppm.1 pgmtoppm.1 pi1toppm.1 picttoppm.1 \
- X pjtoppm.1 ppmdither.1 ppmhist.1 ppmmake.1 \
- X ppmquant.1 ppmrelief.1 ppmtoacad.1 ppmtogif.1 ppmtoicr.1 \
- X ppmtoilbm.1 ppmtopcx.1 ppmtopgm.1 ppmtopi1.1 ppmtopict.1 \
- X ppmtopj.1 ppmtopuzz.1 ppmtorgb3.1 ppmtosixel.1 \
- X ppmtotga.1 ppmtouil.1 ppmtoxpm.1 ppmtoyuv.1 qrttoppm.1 \
- X rawtoppm.1 rgb3toppm.1 sldtoppm.1 spctoppm.1 sputoppm.1 \
- X tgatoppm.1 ximtoppm.1 xpmtoppm.1 yuvtoppm.1 \
- X ppmforge.1 ppmpat.1 \
- X ppmquantall.1
- X MAN3 = libppm.3
- X MAN5 = ppm.5
- X
- X SRCS = giftoppm.c gouldtoppm.c ilbmtoppm.c imgtoppm.c mtvtoppm.c \
- X pcxtoppm.c pgmtoppm.c pi1toppm.c picttoppm.c \
- X pjtoppm.c ppmdither.c ppmhist.c ppmmake.c \
- X ppmquant.c ppmrelief.c ppmtoacad.c ppmtogif.c ppmtoicr.c \
- X ppmtoilbm.c ppmtopcx.c ppmtopgm.c ppmtopi1.c ppmtopict.c \
- X ppmtopj.c ppmtopuzz.c ppmtorgb3.c ppmtosixel.c \
- X ppmtotga.c ppmtouil.c ppmtoxpm.c ppmtoyuv.c qrttoppm.c \
- X rawtoppm.c rgb3toppm.c sldtoppm.c spctoppm.c sputoppm.c \
- X tgatoppm.c ximtoppm.c xpmtoppm.c yuvtoppm.c \
- X ppmforge.c ppmpat.c
- X
- X OBJS = giftoppm.o gouldtoppm.o ilbmtoppm.o imgtoppm.o mtvtoppm.o \
- X pcxtoppm.o pgmtoppm.o pi1toppm.o picttoppm.o \
- X pjtoppm.o ppmdither.o ppmhist.o ppmmake.o \
- X ppmquant.o ppmrelief.o ppmtoacad.o ppmtogif.o ppmtoicr.o \
- X ppmtoilbm.o ppmtopcx.o ppmtopgm.o ppmtopi1.o ppmtopict.o \
- X ppmtopj.o ppmtopuzz.o ppmtorgb3.o ppmtosixel.o \
- X ppmtotga.o ppmtouil.o ppmtoxpm.o ppmtoyuv.o qrttoppm.o \
- X rawtoppm.o rgb3toppm.o sldtoppm.o spctoppm.o sputoppm.o \
- X tgatoppm.o ximtoppm.o xpmtoppm.o yuvtoppm.o \
- X ppmforge.o ppmpat.o
- X
- X BINS = giftoppm gouldtoppm ilbmtoppm imgtoppm mtvtoppm \
- X pcxtoppm pgmtoppm pi1toppm picttoppm \
- X pjtoppm ppmdither ppmhist ppmmake \
- X ppmquant ppmrelief ppmtoacad ppmtogif ppmtoicr \
- X ppmtoilbm ppmtopcx ppmtopgm ppmtopi1 ppmtopict \
- X ppmtopj ppmtopuzz ppmtorgb3 ppmtosixel \
- X ppmtotga ppmtouil ppmtoxpm ppmtoyuv qrttoppm \
- X rawtoppm rgb3toppm sldtoppm spctoppm sputoppm \
- X tgatoppm ximtoppm xpmtoppm yuvtoppm \
- X ppmforge ppmpat
- X
- includes:: ppmquantall.script
- X
- ppmquantall.script:
- X $(LN) ppmquantall ppmquantall.script
- X
- AllTarget($(LIBPPM) $(BINS))
- X
- DependTarget()
- X
- NormalPbmplusProgramTarget(giftoppm)
- NormalPbmplusProgramTarget(gouldtoppm)
- NormalPbmplusProgramTarget(ilbmtoppm)
- NormalPbmplusProgramTarget(imgtoppm)
- NormalPbmplusProgramTarget(mtvtoppm)
- NormalPbmplusProgramTarget(pcxtoppm)
- NormalPbmplusProgramTarget(pgmtoppm)
- NormalPbmplusProgramTarget(pi1toppm)
- NormalPbmplusProgramTarget(picttoppm)
- NormalPbmplusProgramTarget(pjtoppm)
- NormalPbmplusProgramTarget(ppmdither)
- NormalPbmplusProgramTarget(ppmhist)
- NormalPbmplusProgramTarget(ppmmake)
- NormalPbmplusProgramTarget(ppmquant)
- NormalPbmplusProgramTarget(ppmrelief)
- NormalPbmplusProgramTarget(ppmtoacad)
- NormalPbmplusProgramTarget(ppmtogif)
- NormalPbmplusProgramTarget(ppmtoicr)
- NormalPbmplusProgramTarget(ppmtoilbm)
- NormalPbmplusProgramTarget(ppmtopcx)
- NormalPbmplusProgramTarget(ppmtopgm)
- NormalPbmplusProgramTarget(ppmtopi1)
- NormalPbmplusProgramTarget(ppmtopict)
- NormalPbmplusProgramTarget(ppmtopj)
- NormalPbmplusProgramTarget(ppmtopuzz)
- NormalPbmplusProgramTarget(ppmtorgb3)
- NormalPbmplusProgramTarget(ppmtosixel)
- NormalPbmplusProgramTarget(ppmtotga)
- NormalPbmplusProgramTarget(ppmtouil)
- NormalPbmplusProgramTarget(ppmtoxpm)
- NormalPbmplusProgramTarget(ppmtoyuv)
- NormalPbmplusProgramTarget(qrttoppm)
- NormalPbmplusProgramTarget(rawtoppm)
- NormalPbmplusProgramTarget(rgb3toppm)
- NormalPbmplusProgramTarget(sldtoppm)
- NormalPbmplusProgramTarget(spctoppm)
- NormalPbmplusProgramTarget(sputoppm)
- NormalPbmplusProgramTarget(tgatoppm)
- NormalPbmplusProgramTarget(ximtoppm)
- NormalPbmplusProgramTarget(xpmtoppm)
- NormalPbmplusProgramTarget(yuvtoppm)
- NormalPbmplusMathProgramTarget(ppmforge)
- NormalPbmplusMathProgramTarget(ppmpat)
- X
- NormalLibraryObjectRule()
- NormalLibraryTarget(ppm,libppm1.o libppm2.o libppm3.o libppm4.o libppm5.o)
- X
- #if InstallMerged
- NormalProgramTarget($(MERGE),$(MERGE).o $(OBJS),$(CURRENTDEPLIBS),$(CURRENTLIBS),-lm)
- #if InstallBinaries
- InstallProgram($(MERGE),$(PBMPLUSDIR)$(PBMPLUSBINDIR))
- #endif
- #endif
- X
- #if InstallBinaries
- InstallPbmplusPrograms($(BINS),$(PBMPLUSDIR)$(PBMPLUSBINDIR),$(INSTPGMFLAGS))
- InstallScript(ppmquantall,$(PBMPLUSDIR)$(PBMPLUSBINDIR))
- #endif
- X
- #if InstallManuals
- InstallMultipleMan($(MAN1),$(PBMPLUSDIR)$(PBMPLUSMANDIR)/man1)
- InstallMultipleMan($(MAN3),$(PBMPLUSDIR)$(PBMPLUSMANDIR)/man3)
- InstallMultipleMan($(MAN5),$(PBMPLUSDIR)$(PBMPLUSMANDIR)/man5)
- #endif
- X
- #if InstallLibraries
- InstallLibrary(ppm,$(PBMPLUSDIR)$(PBMPLUSLIBDIR))
- #endif
- X
- #if InstallIncludes
- InstallMultipleFlags(ppm.h,$(PBMPLUSDIR)$(PBMPLUSINCDIR),$(INSTINCFLAGS))
- #endif
- SHAR_EOF
- chmod 0664 ppm/Imakefile ||
- echo 'restore of ppm/Imakefile failed'
- Wc_c="`wc -c < 'ppm/Imakefile'`"
- test 5811 -eq "$Wc_c" ||
- echo 'ppm/Imakefile: original size 5811, current size' "$Wc_c"
- rm -f _shar_wnt_.tmp
- fi
- # ============= ppm/ppmforge.1 ==============
- if test -f 'ppm/ppmforge.1' -a X"$1" != X"-c"; then
- echo 'x - skipping ppm/ppmforge.1 (File already exists)'
- rm -f _shar_wnt_.tmp
- else
- > _shar_wnt_.tmp
- echo 'x - extracting ppm/ppmforge.1 (Text)'
- sed 's/^X//' << 'SHAR_EOF' > 'ppm/ppmforge.1' &&
- .TH ppmforge 1 "25 October 1991"
- .IX ppmforge
- .IX fractals
- .IX clouds
- .IX planets
- .IX stars
- .SH NAME
- ppmforge - fractal forgeries of clouds, planets, and starry skies
- .SH SYNOPSIS
- .na
- .nh
- .B ppmforge
- .RB [ -clouds ]
- 'in +9n
- .RB [ -night ]
- .RB [ -dimension
- .IR dimen ]
- .RB [ -hour
- .IR hour ]
- .RB [ -inclination|-tilt
- .IR angle ]
- .RB [ -mesh
- .IR size ]
- .RB [ -power
- .IR factor ]
- .RB [ -glaciers
- .IR level ]
- .RB [ -ice
- .IR level ]
- .RB [ -saturation
- .IR sat ]
- .RB [ -seed
- .IR seed ]
- .RB [ -stars
- .IR fraction ]
- .RB [ -xsize|-width
- .IR width ]
- .RB [ -ysize|-height
- .IR height ]
- .in -9n
- .ad
- .hy
- .SH DESCRIPTION
- .B ppmforge
- generates three kinds of ``random fractal forgeries,'' the term coined
- by Richard F. Voss of the IBM Thomas J. Watson Research Center for
- seemingly realistic pictures of natural objects generated by simple
- algorithms embodying randomness and fractal self-similarity. The
- techniques used by
- .B ppmforge
- are essentially those
- given by Voss[1], particularly the technique of spectral synthesis
- explained in more detail by Dietmar Saupe[2].
- .PP
- The program generates two varieties of pictures: planets and clouds,
- which are just different renderings of data generated in an identical
- manner, illustrating the unity of the fractal structure of these very
- different objects. A third type of picture, a starry sky, is
- synthesised directly from pseudorandom numbers.
- .PP
- The generation of planets or clouds begins with the preparation of an
- array of random data in the frequency domain. The size of this
- array, the ``mesh size,'' can be set with the
- .B -mesh
- option; the larger the mesh the more realistic the pictures but the
- calculation time and memory requirement increases as the square of the
- mesh size. The fractal dimension, which you can specify with the
- .B -dimension
- option, determines the roughness of the terrain on the planet or the
- scale of detail in the clouds. As the fractal dimension is increased,
- more high frequency components are added into the random mesh.
- .PP
- Once the mesh is generated, an inverse two dimensional Fourier
- transform is performed upon it. This converts the original random
- frequency domain data into spatial amplitudes. We scale the real
- components that result from the Fourier transform into numbers from 0
- to 1 associated with each point on the mesh. You can further
- modify this number by applying a ``power law scale'' to it with the
- .B -power
- option. Unity scale
- leaves the numbers unmodified; a power scale of 0.5 takes the square
- root of the numbers in the mesh, while a power scale of 3 replaces the
- numbers in the mesh with their cubes. Power law scaling is best
- envisioned by thinking of the data as representing the elevation of
- terrain; powers less than 1 yield landscapes with vertical scarps that
- look like glacially-carved valleys; powers greater than one make
- fairy-castle spires (which require large mesh sizes and high
- resolution for best results).
- .PP
- After these calculations, we have a array of the specified size
- containing numbers that range from 0 to 1. The pixmaps are generated as
- follows:
- .TP 10
- .B Clouds
- A colour map is created that ranges from pure blue to white by
- increasing admixture (desaturation) of blue with white. Numbers less
- than 0.5 are coloured blue, numbers between 0.5 and 1.0 are coloured
- with corresponding levels of white, with 1.0 being pure white.
- .TP
- .B Planet
- The mesh is projected onto a sphere. Values less than 0.5 are treated
- as water and values between 0.5 and 1.0 as land. The water areas are
- coloured based upon the water depth, and land based on its elevation.
- The random depth data are used to create clouds over the oceans. An
- atmosphere approximately like the Earth's is simulated; its light
- absorption is calculated to create a blue cast around the limb of the
- planet. A function that rises from 0 to 1 based on latitude is
- modulated by the local elevation to generate polar ice caps--high
- altitude terrain carries glaciers farther from the pole. Based on the
- position of the star with respect to the observer, the apparent colour
- of each pixel of the planet is calculated by ray-tracing from the star
- to the planet to the observer and applying a lighting model that sums
- ambient light and diffuse reflection (for most planets ambient light
- is zero, as their primary star is the only source of illumination).
- Additional random data are used to generate stars around the planet.
- .TP
- .B Night
- A sequence of pseudorandom numbers is used to generate stars with a
- user specified density.
- .PP
- Cloud pictures always contain 256 or fewer colours and may be
- displayed on most colour mapped devices without further processing.
- Planet pictures often contain tens of thousands of colours which
- must be compressed with
- .B ppmquant
- or
- .B ppmdither
- before encoding in a colour mapped format. If the display resolution is
- high enough,
- .B ppmdither
- generally produces better looking planets.
- .B ppmquant
- tends to create discrete colour bands, particularly in the oceans,
- which are unrealistic and distracting. The number of colours in starry
- sky pictures generated with the
- .B -night
- option depends on the value specified for
- .BR -saturation .
- Small values limit the colour temperature distribution of the stars
- and reduce the number of colours in the image.
- If the
- .B -saturation
- is set to 0, none of the stars will be coloured and the resulting
- image will never contain more than 256 colours.
- Night sky pictures with many different star colours often look
- best when colour compressed by
- .B pnmdepth
- rather than
- .B ppmquant
- or
- .BR ppmdither .
- Try
- .I newmaxval
- settings of 63, 31, or 15 with
- .B pnmdepth
- to reduce the number of colours in the picture to 256 or fewer.
- .SH OPTIONS
- .TP 10
- .B -clouds
- Generate clouds. A pixmap of fractal clouds is generated. Selecting clouds
- sets the default for fractal dimension to 2.15 and power scale factor
- to 0.75.
- .TP
- .BI -dimension " dimen"
- Sets the fractal dimension to the specified
- .IR dimen ,
- which may be any floating point value between 0 and 3. Higher fractal
- dimensions create more ``chaotic'' images, which require higher
- resolution output and a larger FFT mesh size to look good. If no
- dimension is specified, 2.4 is used when generating planets and 2.15
- for clouds.
- .TP
- .BI -glaciers " level"
- The floating point
- .I level
- setting controls the extent to which terrain elevation causes ice to
- appear at lower latitudes. The default value of 0.75 makes the polar
- caps extend toward the equator across high terrain and forms glaciers
- in the highest mountains, as on Earth. Higher values make ice sheets
- that cover more and more of the land surface, simulating planets in the
- midst of an ice age. Lower values tend to be boring, resulting in
- unrealistic geometrically-precise ice cap boundaries.
- .TP
- .BI -hour " hour"
- When generating a planet,
- .I hour
- is used as the ``hour angle at the central meridian.'' If you specify
- .BR "-hour 12" ,
- for example, the planet will be fully illuminated, corresponding to
- high noon at the longitude at the centre of the screen. You can
- specify any floating point value between 0 and 24 for
- .IR hour ,
- but values which place most of the planet in darkness (0 to 4 and 20
- to 24) result in crescents which, while pretty, don't give you many
- illuminated pixels for the amount of computing that's required. If no
- .B -hour
- option is specified, a random hour angle is chosen, biased so that
- only 25% of the images generated will be crescents.
- .TP
- .BI -ice " level"
- Sets the extent of the polar ice caps to the given floating point
- .IR level .
- The default level of 0.4 produces ice caps similar to those of the Earth.
- Smaller values reduce the amount of ice, while larger
- .B -ice
- settings create more prominent ice caps. Sufficiently large values,
- such as 100 or more, in conjunction with small settings for
- .B -glaciers
- (try 0.1) create ``ice balls'' like Europa.
- .TP
- .BI -inclination|-tilt " angle"
- The inclination angle of the planet with regard to its primary star is
- set to
- .IR angle ,
- which can be any floating point value from -90 to 90. The inclination
- angle can be thought of as specifying, in degrees, the ``season'' the
- planet is presently experiencing or, more precisely, the latitude at
- which the star transits the zenith at local noon. If 0, the planet
- is at equinox; the star is directly overhead at the equator.
- Positive values represent summer in the northern hemisphere, negative
- values summer in the southern hemisphere. The Earth's inclination
- angle, for example, is about 23.5 at the June solstice, 0 at the
- equinoxes in March and September, and -23.5 at the December solstice.
- If no inclination angle is specified, a random value between -21.6 and
- 21.6 degrees is chosen.
- .TP
- .BI -mesh " size"
- A mesh of
- .IR size " by " size
- will be used for the fast Fourier transform (FFT). Note that memory
- requirements and computation speed increase as the square of
- .IR size ;
- if you double the mesh size, the program will use four times the
- memory and run four times as long. The default mesh is 256x256, which
- produces reasonably good looking pictures while using half a megabyte
- for the 256x256 array of single precision complex numbers
- required by the FFT. On machines with limited memory capacity, you
- may have to reduce the mesh size to avoid running out of RAM.
- Increasing the mesh size produces better looking pictures; the
- difference becomes particularly noticeable when generating high
- resolution images with relatively high fractal dimensions (between 2.2
- and 3).
- .TP
- .B -night
- A starry sky is generated. The stars are created by the same algorithm
- used for the stars that surround planet pictures, but the output
- consists exclusively of stars.
- .TP
- .BI -power " factor"
- Sets the ``power factor'' used to scale elevations synthesised from
- the FFT to
- .IR factor ,
- which can be any floating point number greater than zero. If no
- factor is specified a default of 1.2 is used if a planet is being
- generated, or 0.75 if clouds are selected by the
- .B -clouds
- option. The result of the FFT image synthesis is an array of elevation
- values between 0 and 1. A non-unity power factor exponentiates each
- of these elevations to the specified power. For example, a power
- factor of 2 squares each value, while a power factor of 0.5 replaces
- each with its square root. (Note that exponentiating values between 0
- and 1 yields values that remain within that range.) Power factors
- less than 1 emphasise large-scale elevation changes at the expense of
- small variations. Power factors greater than 1 increase the roughness
- of the terrain and, like high fractal dimensions, may require a larger
- FFT mesh size and/or higher screen resolution to look good.
- .TP
- .BI -saturation " sat"
- Controls the degree of colour saturation of the stars that surround planet
- pictures and fill starry skies created with the
- .B -night
- option. The default value of 125 creates stars which resemble the sky
- as seen by the human eye from Earth's surface. Stars are dim; only
- the brightest activate the cones in the human retina, causing colour
- to be perceived. Higher values of
- .I sat
- approximate the appearance of stars from Earth orbit, where better
- dark adaptation, absence of skyglow, and the concentration of light
- from a given star onto a smaller area of the retina thanks to the lack
- of atmospheric turbulence enhances the perception of colour. Values
- greater than 250 create ``science fiction'' skies that, while pretty,
- don't occur in this universe.
- .TP
- \
- Thanks to the inverse square law combined with Nature's love of
- mediocrity, there are many, many dim stars for every bright one.
- This population relationship is accurately reflected in the skies
- created by
- .BR ppmforge .
- Dim, low mass stars live much longer than bright massive stars,
- consequently there are many reddish stars for every blue giant. This
- relationship is preserved by
- .BR ppmforge .
- You can reverse the proportion, simulating the sky as seen in a starburst
- galaxy, by specifying a negative
- .I sat
- value.
- .TP
- .BI -seed " num"
- Sets the seed for the random number generator to the integer
- .IR num .
- The seed used to create each picture is displayed on standard output (unless
- suppressed with the
- .B -quiet
- option). Pictures generated with the same seed will be identical. If no
- .B -seed
- is specified, a random seed derived from the date and time will be
- chosen. Specifying an explicit seed allows you to re-render a picture
- you particularly like at a higher resolution or with different viewing
- parameters.
- .TP
- .BI -stars " fraction"
- Specifies the percentage of pixels, in tenths of a percent, which will
- appear as stars, either surrounding a planet or filling the entire
- frame if
- .B -night
- is specified. The default
- .I fraction
- is 100.
- .TP
- .BI -xsize|-width " width"
- Sets the width of the generated image to
- .I width
- pixels. The default width is 256 pixels. Images must be at least as
- wide as they are high; if a width less than the height is specified,
- it will be increased to equal the height. If you must have a long
- skinny pixmap, make a square one with
- .BR ppmforge ,
- then use
- .B pnmcut
- to extract a portion of the shape and size you require.
- .TP
- .BI -ysize|-height " height"
- Sets the height of the generated image to
- .I height
- pixels. The default height is 256 pixels. If the height specified
- exceeds the width, the width will be increased to equal the height.
- .PP
- All flags can be abbreviated to their shortest unique prefix.
- .SH BUGS
- .PP
- The algorithms require the output pixmap to be at least as wide as it
- is high, and the width to be an even number of pixels. These
- constraints are enforced by increasing the size of the requested
- pixmap if necessary.
- .PP
- You may have to reduce the FFT mesh size on machines with 16 bit
- integers and segmented pointer architectures.
- .SH "SEE ALSO"
- .PD
- .BR pnmcut (1),
- .BR pnmdepth (1),
- .BR ppmdither (1),
- .BR ppmquant (1),
- .BR ppm (5)
- .TP 5
- [1]
- Voss, Richard F., ``Random Fractal Forgeries,'' in Earnshaw
- et. al., Fundamental Algorithms for Computer Graphics, Berlin:
- Springer-Verlag, 1985.
- .TP
- [2]
- Peitgen, H.-O., and Saupe, D. eds., The Science Of Fractal Images,
- New York: Springer Verlag, 1988.
- .ne 10
- .SH AUTHOR
- .RS 5
- .nf
- John Walker
- Autodesk SA
- Avenue des Champs-Montants 14b
- CH-2074 MARIN
- Suisse/Schweiz/Svizzera/Svizra/Switzerland
- .PD 0
- .TP 9
- Usenet:
- kelvin@Autodesk.com
- .TP
- Fax:
- 038/33 88 15
- .TP
- Voice:
- 038/33 76 33
- .fi
- .RE
- .PD
- .PP
- 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.
- .PP
- .B 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
- local software store or directly from Autodesk, Inc., Attn: Science
- 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
- 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 Gleick and
- detailed explanations by Rudy Rucker of the mathematics and algorithms
- used by each program.
- SHAR_EOF
- chmod 0664 ppm/ppmforge.1 ||
- echo 'restore of ppm/ppmforge.1 failed'
- Wc_c="`wc -c < 'ppm/ppmforge.1'`"
- test 15385 -eq "$Wc_c" ||
- echo 'ppm/ppmforge.1: original size 15385, current size' "$Wc_c"
- rm -f _shar_wnt_.tmp
- fi
- # ============= ppm/ppmforge.c ==============
- if test -f 'ppm/ppmforge.c' -a X"$1" != X"-c"; then
- echo 'x - skipping ppm/ppmforge.c (File already exists)'
- rm -f _shar_wnt_.tmp
- else
- > _shar_wnt_.tmp
- echo 'x - extracting ppm/ppmforge.c (Text)'
- sed 's/^X//' << 'SHAR_EOF' > 'ppm/ppmforge.c' &&
- /*
- X
- X Fractal forgery generator for the PPM toolkit
- X
- X Originally designed and implemented in December of 1989 by
- X John Walker as a stand-alone program for the Sun and MS-DOS
- X under Turbo C. Adapted in September of 1991 for use with Jef
- X Poskanzer's raster toolkit.
- X
- X References cited in the comments are:
- X
- X Foley, J. D., and Van Dam, A., Fundamentals of Interactive
- X Computer Graphics, Reading, Massachusetts: Addison
- X Wesley, 1984.
- X
- X Peitgen, H.-O., and Saupe, D. eds., The Science Of Fractal
- X Images, New York: Springer Verlag, 1988.
- X
- X Press, W. H., Flannery, B. P., Teukolsky, S. A., Vetterling,
- X W. T., Numerical Recipes In C, New Rochelle: Cambridge
- X University Press, 1988.
- X
- X Author:
- X John Walker
- X Autodesk SA
- X Avenue des Champs-Montants 14b
- X CH-2074 MARIN
- X Switzerland
- X Usenet: kelvin@Autodesk.com
- X Fax: 038/33 88 15
- X Voice: 038/33 76 33
- X
- X Permission to use, copy, modify, and distribute this software and
- X its documentation for any purpose and without fee is hereby
- X granted, without any conditions or restrictions. This software is
- X provided "as is" without express or implied warranty.
- X
- X PLUGWARE!
- X
- X If you like this kind of stuff, you may also enjoy "James Gleick's
- X Chaos--The Software" for MS-DOS, available for $59.95 from your
- X local software store or directly from Autodesk, Inc., Attn: Science
- X Series, 2320 Marinship Way, Sausalito, CA 94965, USA. Telephone:
- X (800) 688-2344 toll-free or, outside the U.S. (415) 332-2344 Ext
- X 4886. Fax: (415) 289-4718. "Chaos--The Software" includes a more
- X comprehensive fractal forgery generator which creates
- X three-dimensional landscapes as well as clouds and planets, plus
- X five more modules which explore other aspects of Chaos. The user
- X guide of more than 200 pages includes an introduction by James
- X Gleick and detailed explanations by Rudy Rucker of the mathematics
- X and algorithms used by each program.
- X
- */
- X
- #include "ppm.h"
- #include <math.h>
- X
- #ifndef M_PI
- #define M_PI 3.14159265358979323846
- #endif
- #ifndef M_E
- #define M_E 2.7182818284590452354
- #endif
- X
- /* Definitions used to address real and imaginary parts in a two-dimensional
- X array of complex numbers as stored by fourn(). */
- X
- #define Real(v, x, y) v[1 + (((x) * meshsize) + (y)) * 2]
- #define Imag(v, x, y) v[2 + (((x) * meshsize) + (y)) * 2]
- X
- /* Co-ordinate indices within arrays. */
- X
- #define X 0
- #define Y 1
- #define Z 2
- X
- /* Definition for obtaining random numbers. */
- X
- #define nrand 4 /* Gauss() sample count */
- #define Cast(low, high) ((low)+(((high)-(low)) * ((random() & 0x7FFF) / arand)))
- X
- /* Utility definition to get an array's element count (at compile
- X time). For example:
- X
- X int arr[] = {1,2,3,4,5};
- X ...
- X printf("%d", ELEMENTS(arr));
- X
- X would print a five. ELEMENTS("abc") can also be used to tell how
- X many bytes are in a string constant INCLUDING THE TRAILING NULL. */
- X
- #define ELEMENTS(array) (sizeof(array)/sizeof((array)[0]))
- X
- /* Data types */
- X
- typedef int Boolean;
- #define FALSE 0
- #define TRUE 1
- X
- #define V (void)
- X
- /* Display parameters */
- X
- #define SCRSAT 255 /* Screen maximum intensity */
- X
- /* Local variables */
- X
- static double arand, gaussadd, gaussfac; /* Gaussian random parameters */
- static double fracdim; /* Fractal dimension */
- static double powscale; /* Power law scaling exponent */
- static int meshsize = 256; /* FFT mesh size */
- static unsigned int rseed; /* Current random seed */
- static Boolean seedspec = FALSE; /* Did the user specify a seed ? */
- static Boolean clouds = FALSE; /* Just generate clouds */
- static Boolean stars = FALSE; /* Just generate stars */
- static int screenxsize = 256; /* Screen X size */
- static int screenysize = 256; /* Screen Y size */
- static double inclangle, hourangle; /* Star position relative to planet */
- static Boolean inclspec = FALSE; /* No inclination specified yet */
- static Boolean hourspec = FALSE; /* No hour specified yet */
- static double icelevel; /* Ice cap theshold */
- static double glaciers; /* Glacier level */
- static int starfraction; /* Star fraction */
- static int starcolour; /* Star colour saturation */
- X
- /* FOURN -- Multi-dimensional fast Fourier transform
- X
- X Called with arguments:
- X
- X data A one-dimensional array of floats (NOTE!!! NOT
- X DOUBLES!!), indexed from one (NOTE!!! NOT ZERO!!),
- X containing pairs of numbers representing the complex
- X valued samples. The Fourier transformed results are
- X returned in the same array.
- X
- X nn An array specifying the edge size in each dimension.
- X THIS ARRAY IS INDEXED FROM ONE, AND ALL THE EDGE
- X SIZES MUST BE POWERS OF TWO!!!
- X
- X ndim Number of dimensions of FFT to perform. Set to 2 for
- X two dimensional FFT.
- X
- X isign If 1, a Fourier transform is done; if -1 the inverse
- X transformation is performed.
- X
- X This function is essentially as given in Press et al., "Numerical
- X Recipes In C", Section 12.11, pp. 467-470.
- */
- X
- static void fourn(data, nn, ndim, isign)
- X float data[];
- X int nn[], ndim, isign;
- {
- X register int i1, i2, i3;
- X int i2rev, i3rev, ip1, ip2, ip3, ifp1, ifp2;
- X int ibit, idim, k1, k2, n, nprev, nrem, ntot;
- X float tempi, tempr;
- X double theta, wi, wpi, wpr, wr, wtemp;
- X
- #define SWAP(a,b) tempr=(a); (a) = (b); (b) = tempr
- X
- X ntot = 1;
- X for (idim = 1; idim <= ndim; idim++)
- X ntot *= nn[idim];
- X nprev = 1;
- X for (idim = ndim; idim >= 1; idim--) {
- X n = nn[idim];
- X nrem = ntot / (n * nprev);
- X ip1 = nprev << 1;
- X ip2 = ip1 * n;
- X ip3 = ip2 * nrem;
- X i2rev = 1;
- X for (i2 = 1; i2 <= ip2; i2 += ip1) {
- X if (i2 < i2rev) {
- X for (i1 = i2; i1 <= i2 + ip1 - 2; i1 += 2) {
- X for (i3 = i1; i3 <= ip3; i3 += ip2) {
- X i3rev = i2rev + i3 - i2;
- X SWAP(data[i3], data[i3rev]);
- X SWAP(data[i3 + 1], data[i3rev + 1]);
- X }
- X }
- X }
- X ibit = ip2 >> 1;
- X while (ibit >= ip1 && i2rev > ibit) {
- X i2rev -= ibit;
- X ibit >>= 1;
- X }
- X i2rev += ibit;
- X }
- X ifp1 = ip1;
- X while (ifp1 < ip2) {
- X ifp2 = ifp1 << 1;
- X theta = isign * (M_PI * 2) / (ifp2 / ip1);
- X wtemp = sin(0.5 * theta);
- X wpr = -2.0 * wtemp * wtemp;
- X wpi = sin(theta);
- X wr = 1.0;
- X wi = 0.0;
- X for (i3 = 1; i3 <= ifp1; i3 += ip1) {
- X for (i1 = i3; i1 <= i3 + ip1 - 2; i1 += 2) {
- X for (i2 = i1; i2 <= ip3; i2 += ifp2) {
- X k1 = i2;
- X k2 = k1 + ifp1;
- X tempr = wr * data[k2] - wi * data[k2 + 1];
- X tempi = wr * data[k2 + 1] + wi * data[k2];
- X data[k2] = data[k1] - tempr;
- X data[k2 + 1] = data[k1 + 1] - tempi;
- X data[k1] += tempr;
- X data[k1 + 1] += tempi;
- X }
- X }
- X wr = (wtemp = wr) * wpr - wi * wpi + wr;
- X wi = wi * wpr + wtemp * wpi + wi;
- X }
- X ifp1 = ifp2;
- X }
- X nprev *= n;
- X }
- }
- #undef SWAP
- X
- /* INITGAUSS -- Initialise random number generators. As given in
- X Peitgen & Saupe, page 77. */
- X
- static void initgauss(seed)
- X unsigned int seed;
- {
- X /* Range of random generator */
- X arand = pow(2.0, 15.0) - 1.0;
- X gaussadd = sqrt(3.0 * nrand);
- X gaussfac = 2 * gaussadd / (nrand * arand);
- X srandom(seed);
- }
- X
- /* GAUSS -- Return a Gaussian random number. As given in Peitgen
- X & Saupe, page 77. */
- X
- static double gauss()
- {
- X int i;
- X double sum = 0.0;
- X
- X for (i = 1; i <= nrand; i++) {
- X sum += (random() & 0x7FFF);
- X }
- X return gaussfac * sum - gaussadd;
- }
- X
- /* SPECTRALSYNTH -- Spectrally synthesised fractal motion in two
- X dimensions. This algorithm is given under the
- X name SpectralSynthesisFM2D on page 108 of
- X Peitgen & Saupe. */
- X
- static void spectralsynth(x, n, h)
- X float **x;
- X unsigned int n;
- X double h;
- {
- X unsigned bl;
- X int i, j, i0, j0, nsize[3];
- X double rad, phase, rcos, rsin;
- X float *a;
- X
- X bl = ((((unsigned long) n) * n) + 1) * 2 * sizeof(float);
- X a = (float *) calloc(bl, 1);
- X if (a == (float *) 0) {
- X pm_error("Cannot allocate %d x %d result array (%ld bytes).",
- X n, n, bl);
- X }
- X *x = a;
- X
- X for (i = 0; i <= n / 2; i++) {
- X for (j = 0; j <= n / 2; j++) {
- X phase = 2 * M_PI * ((random() & 0x7FFF) / arand);
- X if (i != 0 || j != 0) {
- X rad = pow((double) (i * i + j * j), -(h + 1) / 2) * gauss();
- X } else {
- X rad = 0;
- X }
- X rcos = rad * cos(phase);
- X rsin = rad * sin(phase);
- X Real(a, i, j) = rcos;
- X Imag(a, i, j) = rsin;
- X i0 = (i == 0) ? 0 : n - i;
- X j0 = (j == 0) ? 0 : n - j;
- X Real(a, i0, j0) = rcos;
- X Imag(a, i0, j0) = - rsin;
- X }
- X }
- X Imag(a, n / 2, 0) = 0;
- X Imag(a, 0, n / 2) = 0;
- X Imag(a, n / 2, n / 2) = 0;
- X for (i = 1; i <= n / 2 - 1; i++) {
- X for (j = 1; j <= n / 2 - 1; j++) {
- X phase = 2 * M_PI * ((random() & 0x7FFF) / arand);
- X rad = pow((double) (i * i + j * j), -(h + 1) / 2) * gauss();
- X rcos = rad * cos(phase);
- X rsin = rad * sin(phase);
- X Real(a, i, n - j) = rcos;
- X Imag(a, i, n - j) = rsin;
- X Real(a, n - i, j) = rcos;
- X Imag(a, n - i, j) = - rsin;
- X }
- X }
- X
- X nsize[0] = 0;
- X nsize[1] = nsize[2] = n; /* Dimension of frequency domain array */
- X fourn(a, nsize, 2, -1); /* Take inverse 2D Fourier transform */
- }
- X
- /* INITSEED -- Generate initial random seed, if needed. */
- X
- static void initseed()
- {
- X int i = time((long *) 0) ^ 0xF37C;
- X V srandom(i);
- X for (i = 0; i < 7; i++)
- X V random();
- X rseed = random();
- }
- X
- /* TEMPRGB -- Calculate the relative R, G, and B components for a
- X black body emitting light at a given temperature.
- X The Planck radiation equation is solved directly for
- X the R, G, and B wavelengths defined for the CIE 1931
- X Standard Colorimetric Observer. The colour
- X temperature is specified in degrees Kelvin. */
- X
- static void temprgb(temp, r, g, b)
- X double temp;
- X double *r, *g, *b;
- {
- X double c1 = 3.7403e10,
- X c2 = 14384.0,
- X er, eg, eb, es;
- X
- /* Lambda is the wavelength in microns: 5500 angstroms is 0.55 microns. */
- X
- #define Planck(lambda) ((c1 * pow((double) lambda, -5.0)) / \
- X (pow(M_E, c2 / (lambda * temp)) - 1))
- X
- X er = Planck(0.7000);
- X eg = Planck(0.5461);
- X eb = Planck(0.4358);
- #undef Planck
- X
- X es = 1.0 / max(er, max(eg, eb));
- X
- X *r = er * es;
- X *g = eg * es;
- X *b = eb * es;
- }
- X
- /* ETOILE -- Set a pixel in the starry sky. */
- X
- static void etoile(pix)
- X pixel *pix;
- {
- X if ((random() % 1000) < starfraction) {
- #define StarQuality 0.5 /* Brightness distribution exponent */
- #define StarIntensity 8 /* Brightness scale factor */
- #define StarTintExp 0.5 /* Tint distribution exponent */
- X double v = StarIntensity * pow(1 / (1 - Cast(0, 0.9999)),
- X (double) StarQuality),
- X temp,
- X r, g, b;
- X
- X if (v > 255) {
- X v = 255;
- X }
- X
- X /* We make a special case for star colour of zero in order to
- X prevent floating point roundoff which would otherwise
- X result in more than 256 star colours. We can guarantee
- X that if you specify no star colour, you never get more than
- X 256 shades in the image. */
- X
- X if (starcolour == 0) {
- X int vi = v;
- X
- X PPM_ASSIGN(*pix, vi, vi, vi);
- X } else {
- X temp = 5500 + starcolour *
- X pow(1 / (1 - Cast(0, 0.9999)), StarTintExp) *
- X ((random() & 7) ? -1 : 1);
- X /* Constrain temperature to a reasonable value: >= 2600K
- X (S Cephei/R Andromedae), <= 28,000 (Spica). */
- X temp = max(2600, min(28000, temp));
- X temprgb(temp, &r, &g, &b);
- X PPM_ASSIGN(*pix, (int) (r * v + 0.499),
- X (int) (g * v + 0.499),
- X (int) (b * v + 0.499));
- X }
- X } else {
- X PPM_ASSIGN(*pix, 0, 0, 0);
- X }
- }
- X
- /* GENPLANET -- Generate planet from elevation array. */
- X
- static void genplanet(a, n)
- X float *a;
- X unsigned int n;
- {
- X int i, j;
- X unsigned char *cp, *ap;
- X double *u, *u1;
- X unsigned int *bxf, *bxc;
- X
- #define RGBQuant 255
- X pixel *pixels; /* Pixel vector */
- X pixel *rpix; /* Current pixel pointer */
- X
- #define Atthick 1.03 /* Atmosphere thickness as a percentage
- X of planet's diameter */
- X double athfac = sqrt(Atthick * Atthick - 1.0);
- X double sunvec[3];
- X
- X Boolean flipped = FALSE;
- X double shang, siang;
- X
- X ppm_writeppminit(stdout, screenxsize, screenysize,
- X (pixval) RGBQuant, FALSE);
- X
- X if (!stars) {
- X u = (double *) malloc((unsigned int) (screenxsize * sizeof(double)));
- X u1 = (double *) malloc((unsigned int) (screenxsize * sizeof(double)));
- X bxf = (unsigned int *) malloc((unsigned int) screenxsize *
- X sizeof(unsigned int));
- X bxc = (unsigned int *) malloc((unsigned int) screenxsize *
- X sizeof(unsigned int));
- X
- X if (u == (double *) 0 || u1 == (double *) 0 ||
- X bxf == (unsigned int *) 0 || bxc == (unsigned int *) 0) {
- X pm_error("Cannot allocate %d element interpolation tables.",
- X screenxsize);
- X }
- X
- X /* Compute incident light direction vector. */
- X
- X shang = hourspec ? hourangle : Cast(0, 2 * M_PI);
- X siang = inclspec ? inclangle : Cast(-M_PI * 0.12, M_PI * 0.12);
- X
- X sunvec[X] = sin(shang) * cos(siang);
- X sunvec[Y] = sin(siang);
- X sunvec[Z] = cos(shang) * cos(siang);
- X
- X /* Allow only 25% of random pictures to be crescents */
- X
- X if (!hourspec && ((random() % 100) < 75)) {
- X flipped = sunvec[Z] < 0 ? TRUE : FALSE;
- X sunvec[Z] = abs(sunvec[Z]);
- X }
- X pm_message("%s: -seed %d -dimension %.2f -power %.2f -mesh %d",
- X clouds ? "clouds" : "planet",
- X rseed, fracdim, powscale, meshsize);
- X if (!clouds) {
- X pm_message(
- X " -inclination %.0f -hour %d -ice %.2f -glaciers %.2f",
- X (siang * (180.0 / M_PI)),
- X (int) (((shang * (12.0 / M_PI)) + 12 +
- X (flipped ? 12 : 0)) + 0.5) % 24, icelevel, glaciers);
- X pm_message(" -stars %d -saturation %d.",
- X starfraction, starcolour);
- X }
- X
- X /* Prescale the grid points into intensities. */
- X
- X cp = (unsigned char *) malloc(n * n);
- X if (cp == (unsigned char *) 0)
- X return;
- X ap = cp;
- X for (i = 0; i < n; i++) {
- X for (j = 0; j < n; j++) {
- X *ap++ = (255.0 * (Real(a, i, j) + 1.0)) / 2.0;
- X }
- X }
- X
- X /* Fill the screen from the computed intensity grid by mapping
- X screen points onto the grid, then calculating each pixel value
- X by bilinear interpolation from the surrounding grid points.
- X (N.b. the pictures would undoubtedly look better when generated
- X with small grids if a more well-behaved interpolation were
- X used.)
- X
- X Before we get started, precompute the line-level interpolation
- X parameters and store them in an array so we don't have to do
- X this every time around the inner loop. */
- X
- #define UPRJ(a,size) ((a)/((size)-1.0))
- X
- X for (j = 0; j < screenxsize; j++) {
- X double bx = (n - 1) * UPRJ(j, screenxsize);
- X
- X bxf[j] = floor(bx);
- X bxc[j] = bxf[j] + 1;
- X u[j] = bx - bxf[j];
- X u1[j] = 1 - u[j];
- X }
- X } else {
- X pm_message("night: -seed %d -stars %d -saturation %d.",
- X rseed, starfraction, starcolour);
- X }
- X
- X pixels = ppm_allocrow(screenxsize);
- X for (i = 0; i < screenysize; i++) {
- X double t, t1, by, dy, dysq, sqomdysq, icet, svx, svy, svz,
- X azimuth;
- X int byf, byc, lcos;
- X
- X if (!stars) { /* Skip all this setup if just stars */
- X by = (n - 1) * UPRJ(i, screenysize);
- X dy = 2 * (((screenysize / 2) - i) / ((double) screenysize));
- X dysq = dy * dy;
- X sqomdysq = sqrt(1.0 - dysq);
- X svx = sunvec[X];
- X svy = sunvec[Y] * dy;
- X svz = sunvec[Z] * sqomdysq;
- X byf = floor(by) * n;
- SHAR_EOF
- true || echo 'restore of ppm/ppmforge.c failed'
- fi
- echo 'End of part 3'
- echo 'File ppm/ppmforge.c is continued in part 4'
- echo 4 > _shar_seq_.tmp
- exit 0
- exit 0 # Just in case...
- --
- Kent Landfield INTERNET: kent@sparky.IMD.Sterling.COM
- Sterling Software, IMD UUCP: uunet!sparky!kent
- Phone: (402) 291-8300 FAX: (402) 291-4362
- Please send comp.sources.misc-related mail to kent@uunet.uu.net.
-