[comp.sources.misc] v18i004: planet - planet generation simulator, Part04/04

allen@viewlogic.com (Dave Allen) (04/08/91)

Submitted-by: Dave Allen <allen@viewlogic.com>
Posting-number: Volume 18, Issue 4
Archive-name: planet/part04
Supersedes: tec: Volume 10, Issue 77-78

#! /bin/sh
# This is a shell archive.  Remove anything before this line, then unpack
# it by saving it into a file and typing "sh file".  To overwrite existing
# files, type "sh file -c".  You can also feed this as standard input via
# unshar, or by typing "sh <file", e.g..  If this archive is complete, you
# will see the following message at the end:
#		"End of archive 4 (of 4)."
# Contents:  src/const.h src/tec.h src/unix.c src/x.c src/tec1.c
#   src/tec2.c src/tec3.c src/globe.c
# Wrapped by allen@baja on Sat Mar 23 16:18:38 1991
PATH=/bin:/usr/bin:/usr/ucb ; export PATH
if test -f src/const.h -a "${1}" != "-c" ; then 
  echo shar: Will not over-write existing file \"src/const.h\"
else
echo shar: Extracting \"src/const.h\" \(467 characters\)
sed "s/^X//" >src/const.h <<'END_OF_src/const.h'
X#define MAXX 100
X#define MAXY 100
X#define PRINTMODE_NONE  0
X#define PRINTMODE_LONG  1
X#define PRINTMODE_SCALE 2
X#define PRINTMODE_SHORT 3
X#define PRINTMODE_GREY  4
X#define PRINTMODE_CLIM  5
X#define DRAW_GREY  1
X#define DRAW_LAND  2
X#define DRAW_CLIM  3
X#define DRAW_TEC   4
X#define LINE_DIAG  1
X#define LINE_CORN  2
X#define LINE_NONE  3
X
X#define ABS(xx) ((xx < 0) ? -xx: xx)
X#define CMP(x) !strcmp (s, x)
X
Xextern double cos(), sin(), sqrt(), atof(), asin(), atan2();
END_OF_src/const.h
if test 467 -ne `wc -c <src/const.h`; then
    echo shar: \"src/const.h\" unpacked with wrong size!
fi
# end of overwriting check
fi
if test -f src/tec.h -a "${1}" != "-c" ; then 
  echo shar: Will not over-write existing file \"src/tec.h\"
else
echo shar: Extracting \"src/tec.h\" \(531 characters\)
sed "s/^X//" >src/tec.h <<'END_OF_src/tec.h'
X#define MAXSPLINTER 10
X#define MAXFRAG 100
X#define MAXPLATE 40
X#define MAXCHANGE 10
X#define REALSCALE 100
X
Xstruct plate { int dx, dy, odx, ody, rx, ry, age, area, id, next; };
X
Xextern int XSIZE,	YSIZE,		MAXSTEP;
Xextern int MAXLIFE,	MAXBUMP,	BUMPTOL;
Xextern int DRAWEVERY,	HYDROPCT,	PRINTMODE;
Xextern int ZINIT,	ZSUBSUME,	ZCOAST;
Xextern int ZSHELF,	ZMOUNTAIN,	ZERODE;
Xextern int RIFTPCT,	DOERODE,	ERODERND;
Xextern int MAXCTRTRY,	RIFTDIST,	BENDEVERY;
Xextern int BENDBY,	SPEEDBASE,	SPEEDRNG;
Xextern int UNDERSCAN;
Xextern double MR[];
END_OF_src/tec.h
if test 531 -ne `wc -c <src/tec.h`; then
    echo shar: \"src/tec.h\" unpacked with wrong size!
fi
# end of overwriting check
fi
if test -f src/unix.c -a "${1}" != "-c" ; then 
  echo shar: Will not over-write existing file \"src/unix.c\"
else
echo shar: Extracting \"src/unix.c\" \(862 characters\)
sed "s/^X//" >src/unix.c <<'END_OF_src/unix.c'
X/* This program is Copyright (c) 1991 David Allen.  It may be freely
X   distributed as long as you leave my name and copyright notice on it.
X   I'd really like your comments and feedback; send e-mail to
X   allen@viewlogic.com, or send us-mail to David Allen, 10 O'Moore Ave,
X   Maynard, MA 01754. */
X
X#include "const.h"
Xextern int XSIZE, YSIZE, MAXSTEP, step;
Xint picktype;
X
Xmain (argc,argv) int argc; char **argv; {
X   unsigned short seed16v[3];
X
X   /* Initialize random number generator */
X   srand (time ((long *) 0));	/* initialize rand() */
X   seed16v[0] = rand(); seed16v[1] = rand(); seed16v[2] = rand();
X   seed48 (seed16v);	/* initialize lrand48() */
X
X   init (*++argv);
X   for (step=0; step<MAXSTEP; step++) onestep();
X   return (0); }
X
Xcheckmouse () { }
X 
Xpanic (s) char *s; { printf ("PANIC: %s\n", s); exit (1); }
X
Xdraw (ctype, ltype, cra, lra) { }
END_OF_src/unix.c
if test 862 -ne `wc -c <src/unix.c`; then
    echo shar: \"src/unix.c\" unpacked with wrong size!
fi
# end of overwriting check
fi
if test -f src/x.c -a "${1}" != "-c" ; then 
  echo shar: Will not over-write existing file \"src/x.c\"
else
echo shar: Extracting \"src/x.c\" \(7119 characters\)
sed "s/^X//" >src/x.c <<'END_OF_src/x.c'
X/* This program is Copyright (c) 1991 David Allen.  It may be freely
X   distributed as long as you leave my name and copyright notice on it.
X   I'd really like your comments and feedback; send e-mail to
X   allen@viewlogic.com, or send us-mail to David Allen, 10 O'Moore Ave,
X   Maynard, MA 01754.
X
X   Based on code written by by George Ferguson (ferguson@cs.rochester.edu),
X   28 Apr 1990.  That code was intended for an X version of "tec", never
X   released. */
X
X#include <math.h>
X#include <X11/StringDefs.h>
X#include <X11/Intrinsic.h>
X#include <X11/keysym.h>
X#include <X11/Xaw/Box.h>
X#include <X11/Xaw/Command.h>
X#include <X11/Xaw/Viewport.h>
X#include <X11/Shell.h>
X#include "const.h"
X#include "clim.h"
X
X#define SIZE 6
X#define WINDOW_HEIGHT   (SIZE*MAXX)
X#define WINDOW_WIDTH    (SIZE*MAXY)
X
Xextern int MAXSTEP, XSIZE, YSIZE, step;
Xint picktype = M_HEAT;
X
Xstatic XtAppContext appc;
Xstatic Display *Dis;
Xstatic Window Win;
Xstatic Pixmap Pix;
Xstatic Pixmap Pixsave;
Xstatic GC gc_def;
Xstatic GC Gc[32];
X
X
Xstatic Arg simple_args[] = {
X   {XtNfromVert, (XtArgVal) NULL},
X   {XtNheight, (XtArgVal) WINDOW_HEIGHT},
X   {XtNwidth, (XtArgVal) WINDOW_WIDTH} };
X
X
Xunsigned char landcols [] = { 11, 28, 3 }, greycols[256], tecols[256];
Xunsigned char climcols [] = { 9, 20, 18, 2, 3, 18, 22, 26, 14, 17, 12 };
Xstatic char *color_names[] = {
X   "black","white","grey75", "grey50",
X   "MediumPurple1","purple4","purple3","purple2","purple1",
X   "blue4","blue3","blue2","blue1",
X   "DarkGreen","green4","green3","green2","green1",
X   "DarkGoldenrod4","yellow4","yellow3","yellow2","yellow1",
X   "orange4","orange3","orange2","orange1",
X   "brown4","red4","red3","red2","red1" };
X
X
Xmain (argc,argv) int argc; char **argv; {
X   Widget toplevel, top_form, surface;
X   XWindowAttributes wattr; ushort seed16v[3];
X
X   XtToolkitInitialize ();
X   appc = XtCreateApplicationContext ();
X
X   Dis = XtOpenDisplay (appc, NULL, NULL, "Demo", NULL, 0, &argc, argv);
X   toplevel = XtAppCreateShell (NULL, NULL, applicationShellWidgetClass,
X      Dis, NULL, 0);
X   top_form = XtCreateManagedWidget ("top form", formWidgetClass,
X      toplevel, NULL, 0);
X   surface = XtCreateManagedWidget ("drawing surface", simpleWidgetClass,
X      top_form, simple_args, XtNumber(simple_args));
X   XtRealizeWidget (toplevel);
X
X   Win = XtWindow (surface);
X   XGetWindowAttributes (Dis, Win, &wattr);
X   Pix = XCreatePixmap (Dis, Win, WINDOW_WIDTH, WINDOW_HEIGHT, wattr.depth);
X   Pixsave = XCreatePixmap (Dis, Win, WINDOW_WIDTH, WINDOW_HEIGHT, wattr.depth);
X   XSelectInput (Dis, Win,
X      ExposureMask | ButtonPressMask | ButtonReleaseMask |
X      KeyPressMask | StructureNotifyMask );
X
X   gc_def = DefaultGC (Dis, DefaultScreen (Dis));
X   assign_colors ();
X   XFlush (Dis);
X
X   /* Initialize random number generator */
X   srand (time ((long *) 0));	/* initialize rand() */
X   seed16v[0] = rand(); seed16v[1] = rand(); seed16v[2] = rand();
X   seed48 (seed16v);	/* initialize lrand48() */
X
X   XFillRectangle (Dis, Pix, Gc[0], 0, 0, MAXX*SIZE, MAXY*SIZE); redraw ();
X   init (*++argv);
X   for (step=0; step<MAXSTEP; step++) {
X      onestep();
X      redraw();
X      checkmouse (); }
X   XtDestroyApplicationContext (appc);
X   return (0); }
X
X
Xassign_colors () {
X   Colormap colormap; XColor screen_in_out,exact; int i;
X
X   for (i=0; i<256; i++) greycols[i] = 4.0 + (double) i * 27.0 / 254.0;
X   for (i=0; i<16; i++)  tecols[i] = 0;
X   for (; i<80; i++)     tecols[i] = 4.0 + (double) (i-16) * 27.0 / 62.0;
X   for (; i<256; i++)    tecols[i] = 1;
X
X   colormap = DefaultColormap (Dis, DefaultScreen(Dis));
X   for (i=0; i < 32; i++) {
X      XAllocNamedColor (Dis, colormap, color_names[i], &screen_in_out, &exact);
X      Gc[i] = XCreateGC (Dis, Win, 0, (XGCValues *) NULL);
X      XSetForeground (Dis, Gc[i], screen_in_out.pixel);
X      if (i == 0) XSetForeground (Dis, gc_def, screen_in_out.pixel); } }
X
X
X/* Each call to this routine dispatches any pending events.  We are only
X   interested in expose events which require a redraw, and the letter 'q'
X   being pressed to quit.  */
Xcheckmouse () { XEvent e; char c; KeySym sym; XComposeStatus status;
X 
X   while (XtAppPending (appc)) { /* while there are events */
X      XtAppNextEvent (appc, &e);    /* get one */
X      switch (e.type) {             /* and process it */
X         case Expose:
X            if (e.xexpose.window != Win) XtDispatchEvent(&e);
X            else if (e.xexpose.count == 0) redraw ();
X            break;
X         case KeyPress:
X            XLookupString(&e,&c,1,&sym,&status);
X            switch (sym) {
X               case XK_q: step     = MAXSTEP; break;
X               case XK_h: picktype = M_HEAT;  break;
X               case XK_p: picktype = M_PRESS; break;
X               case XK_w: picktype = M_WIND;  break;
X               case XK_r: picktype = M_RAIN;  break;
X               case XK_c: picktype = M_CLIM;  break; }
X            break;
X         default: XtDispatchEvent(&e); } } }
X
X
Xredraw () {
X   XCopyArea (Dis, Pix, Win, Gc[0], 0, 0, WINDOW_WIDTH, WINDOW_HEIGHT, 0, 0);
X   XFlush (Dis); }
X
X
Xrnd (top) short top; { return (lrand48 () % top); }
X
X
Xpanic (s) char *s; { printf ("PANIC: %s\n", s); exit (1); }
X
X
Xdraw (ctype, ltype, cra, lra)
X   int ctype, ltype;
X   unsigned char cra[MAXX][MAXY], lra[MAXX][MAXY]; {
X   register short i, j, k, x; unsigned char *lut;
X
X   switch (ctype) {
X      case DRAW_GREY: lut = greycols; break;
X      case DRAW_LAND: lut = landcols; break;
X      case DRAW_CLIM: lut = climcols; break;
X      case DRAW_TEC:  lut = tecols;   break; }
X   for (j=0; j < YSIZE; j++) for (i=0; i < XSIZE-1; i++) {
X      x = lut[cra[i][j]]; k = i+1;
X      while ((lut[cra[k][j]] == x) && (k < XSIZE-1)) k++;
X      XFillRectangle (Dis, Pix, Gc[x], i*SIZE, j*SIZE, (k-i)*SIZE, SIZE);
X      i = k-1; }
X   switch (ltype) {
X      case LINE_DIAG: diagonal (lra); break;
X      case LINE_CORN: corner   (lra); break;
X      case LINE_NONE: break; } }
X
X
X#define LEFT(i,j,x) \
X   XDrawLine (Dis, Pix, Gc[x], i*SIZE, j*SIZE, i*SIZE, (j+1)*SIZE)
X#define ABOVE(i,j,x) \
X   XDrawLine (Dis, Pix, Gc[x], i*SIZE, j*SIZE, (i+1)*SIZE, j*SIZE)
X#define UPLEFT(i,j,x) \
X   XDrawLine (Dis, Pix, Gc[x], i*SIZE, j*SIZE, (i+1)*SIZE, (j+1)*SIZE)
X#define UPRIGHT(i,j,x) \
X   XDrawLine (Dis, Pix, Gc[x], (i+1)*SIZE, j*SIZE, i*SIZE, (j+1)*SIZE)
X
X
Xdiagonal (ra) unsigned char ra[MAXX][MAXY]; { register short i, j;
X   for (j=0; j<YSIZE; j++) for (i=0; i<XSIZE; i++)
X      switch (ra[i][j]) {
X         case 0:   break;
X         case N:   LEFT    (i, j, 0); break;
X         case S:   LEFT    (i, j, 1); break;
X         case E:   ABOVE   (i, j, 1); break;
X         case W:   ABOVE   (i, j, 0); break;
X         case N|E: UPRIGHT (i, j, 0); break;
X         case N|W: UPLEFT  (i, j, 0); break;
X         case S|E: UPLEFT  (i, j, 1); break;
X         case S|W: UPRIGHT (i, j, 1); break; } }
X
X
Xcorner (ra) unsigned char ra[MAXX][MAXY]; { register short i, j, x;
X   for (j=0; j<YSIZE; j++) for (i=0; i<XSIZE; i++) {
X      x = ra[i][j];
X      if      (x & LINE_0V) LEFT  (i, j, 0);
X      else if (x & LINE_1V) LEFT  (i, j, 1);
X      if      (x & LINE_0H) ABOVE (i, j, 0);
X      else if (x & LINE_1H) ABOVE (i, j, 1); } }
END_OF_src/x.c
if test 7119 -ne `wc -c <src/x.c`; then
    echo shar: \"src/x.c\" unpacked with wrong size!
fi
# end of overwriting check
fi
if test -f src/tec1.c -a "${1}" != "-c" ; then 
  echo shar: Will not over-write existing file \"src/tec1.c\"
else
echo shar: Extracting \"src/tec1.c\" \(10750 characters\)
sed "s/^X//" >src/tec1.c <<'END_OF_src/tec1.c'
X/* This program is Copyright (c) 1991 David Allen.  It may be freely
X   distributed as long as you leave my name and copyright notice on it.
X   I'd really like your comments and feedback; send e-mail to
X   allen@viewlogic.com, or send us-mail to David Allen, 10 O'Moore Ave,
X   Maynard, MA 01754. */
X
X/* This is the file containing all of the important functions except for
X   trysplit (), which splits a continent into pieces.  Also, all of the main
X   arrays are declared here, even a couple that are only used by functions in
X   tec2.c.  The array declarations are first, followed by the sequencing
X   function onestep () and some miscellaneous routines including the text
X   output routines; initialization routines and the routines that do all
X   the interesting stuff are last. */
X
X#include "const.h"
X#include "tec.h"
X
X/* These are the parameters and their default values; each is defined in
X   params.doc */ 
Xint XSIZE = 100,	YSIZE = 100,	MAXSTEP = 100;
Xint MAXBUMP = 50,	BUMPTOL = 50,	UNDERSCAN = 0;
Xint HYDROPCT = 70,	DRAWEVERY = 0,	PRINTMODE = PRINTMODE_NONE;
Xint ZINIT = 22,		ZSUBSUME = 16,	ZCOAST = 16;
Xint ZSHELF = 8,		ZMOUNTAIN = 48, ZERODE = 16;
Xint RIFTPCT = 40,	DOERODE = 1,	ERODERND = 4;
Xint MAXCTRTRY = 50,	RIFTDIST = 5,	BENDEVERY = 6;
Xint BENDBY = 50,	SPEEDRNG = 300,	SPEEDBASE = 200;
X
Xdouble MR [] = { 1.0,1.0,1.0,0.7,0.4,0.1,-0.2,-0.5,-0.8,-1.0 };
Xint MAXLIFE = 10; /* Length of MR vector */
Xint change[1]; /* needed for fileio */
X
X/* The following arrays are global and most are used by functions in both
X   source files.  The two main ones are m and t.  Each is set up to be two
X   2-d arrays, where each array is the size of the whole world.  M is the
X   map; elements in m are indices of plates, showing which squares are
X   covered by which plate.  T is the topography; elements in t are altitudes. */
X 
Xchar m[2][MAXX][MAXY]; unsigned char t[2][MAXX][MAXY];
X 
X/* Several arrays are used by the binary blob segmenter, segment() in tec2.c.
X   These include r, which is used to store fragment indices; many fragments
X   make up one region during a segmentation.  Kid is a lookup table; fragment
X   k belongs to region kid[k] after a segmentation is finished.  Karea[k]
X   is the area of fragment k. */
X 
Xchar r[MAXX][MAXY], kid[MAXFRAG]; int karea [MAXFRAG];
X 
X/* The merge routine gets information from the move routine; when the move
X   routine puts a square of one plate on top of another plate, that information
X   is recorded in the merge matrix mm. */
X 
Xchar mm[MAXPLATE][MAXPLATE];
X 
X/* The erosion routine needs an array to store delta information; during an
X   erosion, the increases or decreases in elevation are summed in e and then
X   applied all at once to the topography. */
X 
Xchar e[MAXX][MAXY];
X 
X/* Several routines need temporary storage for areas and plate identifiers. */
X 
Xint tarea[MAXPLATE]; int ids[MAXPLATE];
X 
X/* The plates in use are stored in this data structure.  Dx,dy are the
X   values to move by THIS STEP ONLY; odx,ody are the permanent move
X   values; rx,ry are the remainder x and y values used by newdxy() to
X   determine dx,dy; age is the age of the plate, in steps; area is the
X   area of the plate, in squares; id is the value in the m array which
X   corresponds to this plate; next is a pointer to the next occupied
X   element of the plate array. */
X 
Xstruct plate p [MAXPLATE];
X 
X/* The linked list header for available plates and used plates are global,
X   as is the step counter.  */
X 
Xint pavail, phead, step;
X 
X 
Xonestep () {
X   /* This is the sequencing routine called by main once per step.
X   It just calls the important subfunctions in order:
X   - trysplit   finds a plate to break up, and computes new velocities
X   - newdxy     computes the deltas to move each plate this step
X   - move       moves the plates
X   - merge      determines results when plates rub together
X   - erode      erodes the terrain, adding or subtracting altitude
X   - draw       draw the resulting array once every DRAWEVERY steps
X   The m and t arrays are double-buffered in the sense that operations go
X   from m[0] to m[1] or vice-versa; src and dest determine which is which. */
X 
X   int src, dest;
X 
X   src = step % 2; dest = 1 - src;
X   if (rnd (100) < RIFTPCT) trysplit (src);
X   newdxy ();
X   move (src, dest);
X   merge (dest);
X   if (DOERODE) erode (dest);
X   draw (DRAW_TEC, LINE_NONE, t[dest], 0);
X
X   if (DRAWEVERY) if (step && !(step % DRAWEVERY)) tecst (dest);
X   if (!DRAWEVERY && (step == MAXSTEP - 1)) tecst (dest); }
X 
X 
Xpalloc () {
X   /* Allocate a plate from the array and return its index.  All the fields
X   of the plate are initialized to 0, except `next'.  That field is used to
X   link together the plate structures in use.  */
X 
X   int x;
X
X   if (!pavail) panic ("No more objects");
X   x = pavail; pavail = p[x].next;
X   p[x].next = phead; phead = x;
X   p[x].area = 0; p[x].age = 0;
X   p[x].rx = 0; p[x].ry = 0;
X   p[x].odx = 0; p[x].ody = 0;
X   p[x].dx = 0; p[x].dy = 0;
X   return (x); }
X 
X 
Xpfree (n) int n; {
X   /* Return a plate array element to the pool of available elements.
X   To check for infinite loops, the variable guard is incremented
X   at each operation; if the number of operations exceeds the maximum
X   possible number, the program panics. */
X 
X   int i, guard = 0;
X 
X   if (phead == n) phead = p[n].next;
X   else {
X      for (i=phead; p[i].next!=n; i=p[i].next)
X         if (++guard > MAXPLATE) panic ("Infinite loop in pfree");
X      p[i].next = p[n].next; }
X   p[n].next = pavail; pavail = n; }
X 
X
Xmainpar (s) char *s; { 
X   if      (CMP ("XSIZE"))     getdim  (&XSIZE, 0);
X   else if (CMP ("YSIZE"))     getdim  (&YSIZE, 0);
X   else if (CMP ("MOVERATE"))  getdvec (&MAXLIFE, MR, 0);
X   else if (CMP ("MAXSTEP"))   getlng  (&MAXSTEP, 0);
X   else if (CMP ("MAXBUMP"))   getlng  (&MAXBUMP, 0);
X   else if (CMP ("BUMPTOL"))   getlng  (&BUMPTOL, 0);
X   else if (CMP ("DRAWEVERY")) getlng  (&DRAWEVERY, 0);
X   else if (CMP ("PRINTMODE")) getlng  (&PRINTMODE, 0);
X   else if (CMP ("HYDROPCT"))  getlng  (&HYDROPCT, 0);
X   else if (CMP ("ZINIT"))     getlng  (&ZINIT, 0);
X   else if (CMP ("ZSUBSUME"))  getlng  (&ZSUBSUME, 0);
X   else if (CMP ("ZCOAST"))    getlng  (&ZCOAST, 0);
X   else if (CMP ("ZSHELF"))    getlng  (&ZSHELF, 0);
X   else if (CMP ("ZMOUNTAIN")) getlng  (&ZMOUNTAIN, 0);
X   else if (CMP ("ZERODE"))    getlng  (&ZERODE, 0);
X   else if (CMP ("RIFTPCT"))   getlng  (&RIFTPCT, 0);
X   else if (CMP ("DOERODE"))   getlng  (&DOERODE, 0);
X   else if (CMP ("ERODERND"))  getlng  (&ERODERND, 0);
X   else if (CMP ("MAXCTRTRY")) getlng  (&MAXCTRTRY, 0);
X   else if (CMP ("RIFTDIST"))  getlng  (&RIFTDIST, 0);
X   else if (CMP ("BENDEVERY")) getlng  (&BENDEVERY, 0);
X   else if (CMP ("BENDBY"))    getlng  (&BENDBY, 0);
X   else if (CMP ("SPEEDBASE")) getlng  (&SPEEDBASE, 0);
X   else if (CMP ("SPEEDRNG"))  getlng  (&SPEEDRNG, 0);
X   else if (CMP ("UNDERSCAN")) getlng  (&UNDERSCAN, 0);
X   else return (0);
X   return (1); }
X
X
Xtecst (src) int src; {
X   /* This function is called whenever map output is called for.  It looks
X   at the parameter `printmode' to decide between long text, simple text,
X   and PostScript output formats.  Note that the default for this
X   function is no output at all, corresponding to PRINTMODE_NONE.  If only
X   one output map is desired, then move the coastline up or down to meet the
X   desired hydrographic percentage. */
X
X   register int i, j, zcoast; int hist[256], goal; unsigned char sk[MAXX][MAXY];
X
X   if (!PRINTMODE) return (0);
X   if (!DRAWEVERY) {
X      /* Create a histogram of the output array */
X      for (i=0; i<256; i++) hist[i] = 0;
X      for (i=0; i<XSIZE; i++)
X         for (j=0; j<YSIZE; j++) hist[t[src][i][j]]++;
X
X      /* Starting from the highest altitude, move down until number of */
X      /* squares above water is slightly greater than the exact goal */
X      goal = XSIZE * YSIZE;
X      goal = (goal * (100 - HYDROPCT)) / 100;
X      for (zcoast=255, i=0; zcoast>0; zcoast--)
X         if ((i += hist[zcoast]) > goal) break;
X
X      /* If the new coast level is zero, then there wasn't enough land */
X      /* to meet the goal, even going right down to the ocean floor.  The */
X      /* only possible result is to panic since the goal can't be met. */
X      if (!zcoast) panic ("Scaled till oceans dried up");
X      ZCOAST = zcoast; }
X
X   if (PRINTMODE != PRINTMODE_SHORT) putmat ("LAND", -1, PRINTMODE, t[src], 0);
X   else {
X      for (i=0; i<XSIZE; i++) for (j=0; j<YSIZE; j++) {
X         if (t[src][i][j] < ZCOAST) sk[i][j] = 0;
X         else if (t[src][i][j] > ZMOUNTAIN) sk[i][j] = 2;
X         else sk[i][j] = 1; }
X      putmat ("LAND", -1, PRINTMODE, sk, 0); }
X   return (0); }
X
X
Xdouble greyscale (x) int x; {
X   /* Called by the PostScript print routine, this function simply computes
X   the intensity from 0-1 corresponding to the altitude 0-255 */
X   if (x < ZCOAST) return ((float) 0);
X   return (1.0 - ((x > 128) ? 128 : x) / 128.0); }
X
X
Xinit (s) char *s; {
X   /* This is the catchall function that initializes everything.  First,
X   it calls getparams() in fileio.c to allow the user to set parameters.  Next,
X   it links together the plates onto the free list and starts the used list
X   at empty.  The first plate is created by a fractal technique and then
X   improved.  Finally, the fractal is copied to the data array and drawn.
X   There are two kinds of improvement done here.  First, islands are
X   eliminated by segmenting the blob and erasing all the regions except
X   for the biggest.  Second, oceans inside the blob (holes) are eliminated
X   by segmenting the _ocean_ and filling in all regions except the biggest. */
X 
X   int besti, x; register int i, j;
X 
X   if (s) if (*s) getparams (s); fileinit ();
X 
X   for (i=1; i<MAXPLATE; i++) p[i].next = i + 1;
X   p[MAXPLATE-1].next = 0;
X   pavail = 1; phead = 0;
X 
X   /* Allocate a plate structure for the first plate and make a blob */
X   x = palloc (); makefrac (0, x);
X
X   /* Segment m[0] looking for x, set besti to the largest region, */
X   /* and zero out all the other regions.  This eliminates islands. */
X   besti = singlefy (0, x);
X   if (besti > 0) for (i=1; i<XSIZE; i++) for (j=1; j<YSIZE; j++)
X      if (kid[r[i][j]] != besti) m[0][i][j] = 0;
X 
X   /* Segment m[0] looking for 0 (ocean), set besti to the largest region, */
X   /* and fill in all the other regions.  This eliminates holes in the blob. */
X   besti = singlefy (0, 0);
X   if (besti > 0) for (i=1; i<XSIZE; i++) for (j=1; j<YSIZE; j++)
X      if (kid[r[i][j]] != besti) m[0][i][j] = x;
X 
X   /* Fill the topo structure with the blob shape while finding its area */
X   for (i=0; i<XSIZE; i++) for (j=0; j<YSIZE; j++)
X      if (m[0][i][j]) { t[0][i][j] = ZINIT; p[x].area++; }
X 
X   /* Draw the blob */
X   if (DRAWEVERY) draw (DRAW_TEC, LINE_NONE, t[0], 0); }
END_OF_src/tec1.c
if test 10750 -ne `wc -c <src/tec1.c`; then
    echo shar: \"src/tec1.c\" unpacked with wrong size!
fi
# end of overwriting check
fi
if test -f src/tec2.c -a "${1}" != "-c" ; then 
  echo shar: Will not over-write existing file \"src/tec2.c\"
else
echo shar: Extracting \"src/tec2.c\" \(14644 characters\)
sed "s/^X//" >src/tec2.c <<'END_OF_src/tec2.c'
X/* This program is Copyright (c) 1991 David Allen.  It may be freely
X   distributed as long as you leave my name and copyright notice on it.
X   I'd really like your comments and feedback; send e-mail to
X   allen@viewlogic.com, or send us-mail to David Allen, 10 O'Moore Ave,
X   Maynard, MA 01754. */
X
X/* This file contains the function trysplit(), which is called from
X   onestep() in tec1.c, and its subfunctions.  One of its subfunctions,
X   segment(), is also called from init() in tec1.c.   Trysplit is the
X   function in charge of splitting one plate into smaller plates. */
X
X
X#include "const.h"
X#include "tec.h"
X
X#define PI 3.14159
X#define TWOPI 6.28318
X#define TWOMILLIPI 0.00628318
X
X/* RIFTMARK is the temporary indicator placed in the arrays to indicate
X   the squares a rift has just appeared.  The function stoprift() puts
X   them in, and trysplit() takes them out before anybody can see them. */
X#define RIFTMARK -1
X
X/* These are all defined in tec1.c */
Xextern char m[2][MAXX][MAXY], r[MAXX][MAXY], kid[MAXFRAG];
Xextern unsigned char t[2][MAXX][MAXY];
Xextern int karea[MAXFRAG], tarea[MAXPLATE], ids[MAXPLATE], step;
Xextern struct plate p [MAXPLATE];
X
Xtrysplit (src) int src; {
X   /* Trysplit is called at most once per step in only 40% of the steps.
X   It first draws a rift on one of the plates, then it segments the result
X   into some number of new plates and some splinters.  If exactly two new
X   non-splinter plates are found, new plate structures are allocated, new
X   dx and dy values are computed, and the old plate is freed.  If anything
X   goes wrong, the rift is erased from the array, returning the array to its
X   previous state.  The functions newrift, segment and newplates do most
X   of the work. */
X
X   register int i, j, a; int count, old, frag, reg;
X
X   if (newrift (src, &old)) if (segment (src, old, &frag, &reg)) if (reg > 1) {
X
X         /* Set tarea[i] to areas of the final segmented regions */
X         for (i=0; i<MAXPLATE; i++) tarea[i] = 0;
X         for (i=1; i<=frag; i++) tarea[kid[i]] += karea[i];
X
X         /* Give up unless exactly two regions are large enough */
X         for (i=1, count=0; i<=reg; i++) if (tarea[i] > MAXSPLINTER) count++; 
X         if (count == 2) {
X
X            /* Compute new dx,dy; update m with the ids of the new plates */
X            newplates (src, old);
X            for (i=0, count=0; i<XSIZE; i++) for (j=0; j<YSIZE; j++) {
X               if (a = r[i][j]) m[src][i][j] = ids[kid[a]];
X               if (m[src][i][j] == RIFTMARK) {
X                  m[src][i][j] = 0; t[src][i][j] = 0; } }
X            pfree (old); return (0); } }
X
X   /* If execution reaches here, the split operation failed; remove rift */
X   if (old) for (i=0; i<XSIZE; i++) for (j=0; j<YSIZE; j++)
X      if (m[src][i][j] == RIFTMARK) {
X         m[src][i][j] = old; p[old].area++; }
X   return (0); }
X
X
Xnewrift (src, old) int src, *old; {
X   /* This function randomly picks a center for a new rift, and draws in
X   a curving line until the line hits either the coast or another plate.
X   If another plate is hit, the rift is invalid and the function returns 0.
X   To find a center, the function generates random x,y values until it
X   finds one that is at least RIFTDIST squares from any ocean square.  If a
X   center is found, a random angle is generated; the rift will pass through
X   the center at that angle.  Next, halfrift() is called twice.  Each call
X   generates the rift leaving the center in one direction.  If everything
X   works out, the function returns the id of the plate the rift is on. */
X
X   int x, y, lx, rx, ty, by, i, j, tries = 0, which; double tt;
X
X   /* Generate a random x, y value */
X   getctr: if (tries > MAXCTRTRY) { *old = 0; return (0); }
X   x = rnd (XSIZE); y = rnd (YSIZE);
X
X   /* If the location is ocean, try again */
X   if (!m[src][x][y]) { tries++; goto getctr; }
X
X   /* Set lx,rx,ty,by to the coordinate values of a box 2*RIFTDIST on a side */
X   /* centered on the center.  Clip the values to make sure they are on the */
X   /* array.  Loop through the box; if a point is ocean, try another center. */
X   lx = (x < RIFTDIST) ? 0 : x - RIFTDIST;
X   rx = (x > XSIZE - RIFTDIST - 1) ? XSIZE - 1 : x + RIFTDIST;
X   ty = (y < RIFTDIST) ? 0 : y - RIFTDIST;
X   by = (y > YSIZE - RIFTDIST - 1) ? YSIZE - 1 : y + RIFTDIST;
X   for (i=lx; i<rx; i++) for (j=ty; j<by; j++)
X      if (!m[src][i][j]) { tries++; goto getctr; }
X
X   /* Found a good center, on plate `which'.  Put a rift indicator in the */
X   /* center.  Generate a random angle t, which is really an integer in the */
X   /* range 0-499 multiplied by 2 PI / 1000.  Call halfrift once for each */
X   /* half of the rift; t is the initial angle for the first call, and */
X   /* t + PI is the initial angle for the second call.  If halfrift() */
X   /* returns zero, abort and return 0; otherwise, return the plate id. */
X   which = m[src][x][y]; m[src][x][y] = RIFTMARK;
X   tt = rnd (500) * TWOMILLIPI;
X   if (!halfrift (src, x, y, which, tt)) { *old = which; return (0); }
X   tt += PI; if (tt > TWOPI) tt -= TWOPI;
X   if (!halfrift (src, x, y, which, tt)) { *old = which; return (0); }
X   *old = which; return (1); }
X
X
Xhalfrift (src, cx, cy, which, tt) int src, cx, cy, which; double tt; {
X   /* Draw a rift from cx,cy on plate `which' at angle t.  At the beginning,
X   digitize the angle using Bresenham's algorithm; once in a while thereafter,
X   modify the angle randomly and digitize it again.  For each square travelled,
X   call stoprift() to see if the rift has left the plate. */
X
X   int ddx, ddy, rdx, rdy, draw, i, a; double dx, dy, adx, ady;
X
X   checkmouse ();
X   /* For-loop against SIZE to guard against infinite loops */
X   for (i=0; i<XSIZE; i++) {
X
X      /* If first square or 1/6 chance at each step, digitize */
X      if (!i || !rnd (BENDEVERY)) {
X
X         /* If not first step, modify angle a little */
X         if (i) tt = tt + (rnd (BENDBY<<1) * TWOMILLIPI) - (BENDBY*TWOMILLIPI);
X         if (tt > TWOPI) tt -= TWOPI; if (tt < 0) tt += TWOPI;
X
X         /* Compute dx and dy, scaled so that larger is exactly +1.0 or -1.0 */
X         dy = sin (tt); dx = cos (tt); adx = ABS(dx); ady = ABS(dy);
X         if (adx > ady) { dy = dy / adx; dx = (dx < 0) ? -1.0: 1.0; }
X         else { dx = dx / ady; dy = (dy < 0) ? -1.0: 1.0; }
X
X         /* Convert to integer value and initialize remainder */
X         /* for each coordinate to half value */
X         ddx = REALSCALE * dx; ddy = REALSCALE * dy;
X         rdx = ddx >> 1; rdy = ddy >> 1; }
X
X      /* Main part of loop, draws one square along line.  The basic idea */
X      /* of Bresenham's algorithm is that if the slope of the line is less */
X      /* than 45 degrees, each time you step one square in X and maybe step */
X      /* one square in Y.  If the slope is greater than 45, step one square */
X      /* in Y and maybe one square in X.  Here, if the slope is less than 45 */
X      /* then ddx == REALSCALE (or -REALSCALE) and the first call to */
X      /* stoprift() is guaranteed.  If stoprift returns <0, all is ok; */
X      /* if zero, the rift ran into the ocean, so stop now; if positive, the */
X      /* rift ran into another plate, which is a perverse condition and the */
X      /* rift must be abandoned.  */
X      rdx += ddx; rdy += ddy;
X      if (rdx >=  REALSCALE) { cx++; rdx -= REALSCALE; draw = 1; }
X      if (rdx <= -REALSCALE) { cx--; rdx += REALSCALE; draw = 1; }
X      if (draw == 1) {
X         a = stoprift (src, cx, cy, which); if (a >= 0) return (!a); }
X      if (rdy >=  REALSCALE) { cy++; rdy -= REALSCALE; draw = 2; }
X      if (rdy <= -REALSCALE) { cy--; rdy += REALSCALE; draw = 2; }
X      if (draw == 2) {
X         a = stoprift (src, cx, cy, which); if (a >= 0) return (!a); } }
X   return (1); }
X
X
Xstoprift (src, x, y, which) int src, x, y, which; {
X   /* This function is called once for each square the rift enters.  It
X   puts a rift marker into m[src] and decides whether the rift can go on.
X   It looks at all four adjacent squares.  If one of them contains ocean
X   or another plate, return immediately so that the rift stops (if ocean)
X   or aborts (if another plate).  If none of them do, then check to make
X   sure the point is within underscan boundaries.  If so, return ok. */
X
X   register int w, a;
X
X   w = which; p[w].area--; m[src][x][y] = RIFTMARK;
X   a = m[src][x][y+1]; if ((a != w) && (a!= RIFTMARK)) return (a != 0);
X   a = m[src][x][y-1]; if ((a != w) && (a!= RIFTMARK)) return (a != 0);
X   a = m[src][x+1][y]; if ((a != w) && (a!= RIFTMARK)) return (a != 0);
X   a = m[src][x-1][y]; if ((a != w) && (a!= RIFTMARK)) return (a != 0);
X   if ((x < UNDERSCAN) || (x > XSIZE - UNDERSCAN)) return (1);
X   if ((y < UNDERSCAN) || (y > YSIZE - UNDERSCAN)) return (1);
X   return (-1); }
X
X
Xsegment (src, match, frag, reg) int src, match, *frag, *reg; {
X   /* This routine implements a standard binary-blob segmentation.  It looks
X   at the array m[src]; match is the value of the blob, and everything else
X   is background.  The result is placed into array r and vectors kid and karea.
X   One 8-connected region can be made up of many fragments; each fragment is
X   assigned a unique index.  Array r contains the frag indices k, while kid[k]
X   is the region frag k belongs to and karea[k] is the area of frag k.
X   Variables frag and reg are set on output to the number of fragments and
X   regions found during the segmentation.  The private vector kk provides one
X   level of indirection for merging fragments; fragment k is merged with
X   fragment kk[k] where kk[k] is the smallest frag index in the region. */
X
X   register int i, j, k, k1, k2, k3, l;
X   char kk [MAXFRAG];
X
X   /* Initialize all frag areas to zero and every frag to merge with itself */
X   for (k=0; k<MAXFRAG; k++) { kk[k] = k; karea[k] = 0; }
X   checkmouse ();
X
X   /* Look at every point in the array */
X   for (k=0, i=0; i<XSIZE; i++) for (j=0; j<YSIZE; j++) {
X
X      /* If too many fragments, give up */
X      if (k == MAXFRAG) return (0);
X
X      /* If this square isn't part of the blob, try the next square */
X      if (m[src][i][j] != match) { r[i][j] = 0; goto bottom; }
X
X      /* It is part of the blob.  Set k1 to the frag id of the square to */
X      /* its left, and set k2 to the frag id of the square above it.  Note */
X      /* that because of the for-loop direction, both of these squares have */
X      /* already been processed. */
X      k1 = i ? kk [r [i-1] [j]] : 0; k2 = j ? kk [r [i] [j-1]] : 0;
X
X      /* If k1 and k2 are both background, start a new fragment */
X      if (!k1 && !k2) { r[i][j] = ++k; karea[k]++; goto bottom; }
X
X      /* If k1 and k2 are part of the same frag, add this square to it */
X      if (k1 && (k1 == k2)) { r[i][j] = k1; karea[k1]++; goto bottom; }
X
X      /* If k1 and k2 belong to different frags, merge them by finding */
X      /* all the frags merged with max(k1,k2) and merging them instead */
X      /* with min(k1,k2).  Add k to that fragment as well. */
X      if (k1 && k2) {
X         if (k2 < k1) { k3 = k1; k1 = k2; k2 = k3; }
X         for (l=1; l<=k; l++) if (kk[l] == k2) kk[l] = k1;
X         r[i][j] = k1; karea[k1]++; goto bottom; }
X
X      /* Default case is that one of k1,k2 is a fragment and the other is */
X      /* background.  Add k to the fragment. */
X      k3 = (k1) ? k1 : k2; r[i][j] = k3; karea[k3]++;
X   bottom: continue; }
X
X   /* Set up vector kid to map from fragments to regions by using i to count */
X   /* unique groups of fragments.  A unique group of fragments is */
X   /* characterized by kk[k] == k; otherwise, frag k is merged with some */
X   /* other fragment. */
X   for (i=0, j=1; j<=k; j++) {
X      if (j == kk[j]) kid[j] = ++i;
X      else kid[j] = kid [kk [j]]; }
X
X   /* Make sure the id of the background is zero; set up return values */
X   kid[0] = 0; *frag = k; *reg = i; return (1); }
X
X
X
Xnewplates (src, old) int src, old; {
X   /* Compute new dx and dy values for plates right after fragmentation.  This
X   function looks at the rift markers in m[src]; variable old is the index of
X   the plate from which the new plates were created.  For each plate adjacent
X   to the rift, this function subtracts the number of plate squares to the left
X   of the rift from the number to the right; this gives some indication of
X   whether the plate should move left or right, and how fast.  The same is done
X   for squares above and below the rift.  The results are put into dx[] and
X   dy[].  At this point some unscaled movement vector is available for both of
X   the new plates.  The vectors are then scaled by the relative sizes of the
X   plates.  The idea is that if one plate is much larger than the other, the
X   small one should move faster.  New plate structures are allocated for the
X   new plates, and the computed dx and dy values are put in them. */
X
X   int dx[MAXPLATE], dy[MAXPLATE];
X   register int i, j, a; int totarea=0, maxmag=0; double scale, b;
X
X   for (i=1; i<MAXPLATE; i++) { dx[i] = 0; dy[i] = 0; ids[i] = 0; }
X   checkmouse ();
X
X   /* For every point in the array, set a to the region id (kid is the */
X   /* lookup table and r contains frag indices); if a is nonzero and */
X   /* the rift is adjacent, adjust counters appropriately */
X   for (i=0; i<XSIZE; i++) for (j=0; j<YSIZE; j++) if (a = kid[r[i][j]]) {
X      if ((i-1 > -1)    && (m[src][i-1][j] == RIFTMARK)) (dx[a])++;
X      if ((i+1 < XSIZE) && (m[src][i+1][j] == RIFTMARK)) (dx[a])--;
X      if ((j-1 > -1)    && (m[src][i][j-1] == RIFTMARK)) (dy[a])++;
X      if ((j+1 < XSIZE) && (m[src][i][j+1] == RIFTMARK)) (dy[a])--; }
X
X   /* For those regions larger than splinters (tarea is set up in trysplit), */
X   /* allocate a new plate structure and initialize its area; compute the */
X   /* magnitude of the dx dy vector and remember the maximum magnitude; also */
X   /* record the total area of new regions */
X   for (i=1; i<MAXPLATE; i++) if (tarea[i] > MAXSPLINTER) {
X      ids[i] = palloc (); p[ids[i]].area = tarea[i];
X      totarea += tarea[i];
X      a =sqrt ((double) ((dx[i]*dx[i]) + (dy[i]*dy[i])));
X      if (a > maxmag) maxmag = a; }
X
X   /* Generate a random speed and predivide so that all speeds computed */
X   /* below are less than the random speed. */
X   scale = (double) (rnd (SPEEDRNG) + SPEEDBASE) / (maxmag * totarea);
X
X   /* Compute the dx and dy for each new plate; note that the speed the */
X   /* plate was moving at before splitting is given by p[old].odx,ody */
X   /* but those must be multiplied by MR to get the actual values */
X   for (i=1; i<MAXPLATE; i++) if (ids[i]) {
X      b = scale * (totarea - tarea[i]);
X      p[ids[i]].odx = p[old].odx * MR [p[old].age] + dx[i] * b;
X      p[ids[i]].ody = p[old].ody * MR [p[old].age] + dy[i] * b; } }
END_OF_src/tec2.c
if test 14644 -ne `wc -c <src/tec2.c`; then
    echo shar: \"src/tec2.c\" unpacked with wrong size!
fi
# end of overwriting check
fi
if test -f src/tec3.c -a "${1}" != "-c" ; then 
  echo shar: Will not over-write existing file \"src/tec3.c\"
else
echo shar: Extracting \"src/tec3.c\" \(14070 characters\)
sed "s/^X//" >src/tec3.c <<'END_OF_src/tec3.c'
X/* This program is Copyright (c) 1991 David Allen.  It may be freely
X   distributed as long as you leave my name and copyright notice on it.
X   I'd really like your comments and feedback; send e-mail to
X   allen@viewlogic.com, or send us-mail to David Allen, 10 O'Moore Ave,
X   Maynard, MA 01754. */
X
X#include "const.h"
X#include "tec.h"
X
X/* These are all defined in tec1.c */
Xextern char m[2][MAXX][MAXY], r[MAXX][MAXY], e[MAXX][MAXY], kid[MAXFRAG],
X   mm[MAXPLATE][MAXPLATE];
Xextern unsigned char t[2][MAXX][MAXY];
Xextern int step, phead, tarea[MAXPLATE], karea[MAXFRAG];
Xextern struct plate p [MAXPLATE];
X
X
Xmakefrac (src, match) int src, match; {
X   /* This function uses a very simple fractal technique to draw a blob.  Four
X   one-dimensional fractals are created and stored in array x, then the
X   fractals are superimposed on a square array, summed and thresholded to
X   produce a binary blob.  Squares in the blob are set to the value in `match'.
X   A one-dimensional fractal of length n is computed like this.  First,
X   set x [n/2] to some height and set the endpoints (x[0] and x[n]) to 0.
X   Then do log-n iterations.  The first iteration computes 2 more values:
X   x[n/4] = average of x[0] and x[n/2], plus some random number, and
X   x[3n/4] = average of x[n/2] and x[n], plus some random number.  The second
X   iteration computes 4 more values (x[n/8], x[3n/8], ...) and so on.  The
X   random number gets smaller by a factor of two each time also.
X 
X   Anyway, you wind up with a number sequence that looks like the cross-section
X   of a mountain.  If you sum two fractals, one horizontal and one vertical,
X   you get a 3-d mountain; but it looks too symmetric.  If you sum four,
X   including the two 45 degree diagonals, you get much better results. */
X 
X   register int xy, dist, n, inc, i, j, k; char x[4][65];
X   int xofs, yofs, xmin, ymin, xmax, ymax, bloblevel, hist[256];
X 
X   /* Compute offsets to put center of blob in center of the world, and
X      compute array limits to clip blob to world size */
X   xofs = (XSIZE - 64) >> 1; yofs = (YSIZE - 64) >> 1;
X   if (xofs < 0) { xmin = -xofs; xmax = 64 + xofs; }
X   else { xmin = 0; xmax = 64; }
X   if (yofs < 0) { ymin = -yofs; ymax = 64 + yofs; }
X   else { ymin = 0; ymax = 64; }
X 
X   for (xy=0; xy<4; xy++) {
X      /* Initialize loop values and fractal endpoints */
X      x [xy] [0] = 0; x [xy] [64] = 0; dist = 32;
X      x [xy] [32] = 24 + rnd (16); n = 2; inc = 16;
X 
X      /* Loop log-n times, each time halving distance and doubling iterations */
X      for (i=0; i<5; i++, dist>>=1, n<<=1, inc>>=1)
X         for (j=0, k=0; j<n; j++, k+=dist)
X            x[xy][k+inc] = ((x[xy][k]+x[xy][k+dist])>>1) + rnd (dist) - inc; }
X 
X   /* Superimpose fractals into the output array.  x[0] is horizontal, x[1] */
X   /* vertical, x[2] diagonal from top left to bottom right, x[3] diagonal */
X   /* from TR to BL.  While superimposing, create a histogram of the values.*/
X   for (i=0; i<256; i++) hist[i] = 0;
X   for (i=xmin; i<xmax; i++) for (j=ymin; j<ymax; j++) {
X      k = x[0][i] + x[1][j] + x[2][(i - j + 64) >> 1] + x[3][(i + j) >> 1];
X      if (k < 0) k = 0; if (k > 255) k = 255;
X      hist[k]++; m[src][i+xofs][j+yofs] = k; }
X
X   /* Pick a threshhold to get as close to the goal number of squares as */
X   /* possible, then go back through the array and adjust it */
X   bloblevel = XSIZE * YSIZE * (100 - HYDROPCT) / 100;
X   for (k=255, i=0; k>=0; k--) if ((i += hist[k]) > bloblevel) break;
X   for (i=xmin; i<xmax; i++) for (j=ymin; j<ymax; j++)
X      m[src][i+xofs][j+yofs] = (m[src][i+xofs][j+yofs] > k) ? match : 0; }
X 
X 
Xsinglefy (src, match) int src, match; {
X   /* This is a subfunction of init() which is called twice to improve the
X   fractal blob.  It calls segment() and then interprets the result.  If
X   only one region was found, no improvement is needed; otherwise, the
X   area of each region must be computed by summing the areas of all its
X   fragments, and the index of the largest region is returned. */
X 
X   int i, reg, frag, besti, besta;
X 
X   segment (src, match, &frag, &reg);
X   if (reg == 1) return (-1); /* No improvement needed */
X 
X   /* Initialize the areas to zero, then sum frag areas */
X   for (i=1; i<=reg; i++) tarea[i] = 0;
X   for (i=1; i<=frag; i++) tarea [kid[i]] += karea [i];
X 
X   /* Pick largest area of all regions and return it */
X   for (i=1, besta=0, besti=0; i<=reg; i++)
X      if (besta < tarea[i]) { besti = i; besta = tarea[i]; }
X   return (besti); }
X 
X 
Xnewdxy () {
X   /* For each plate, compute how many squares it should move this step.
X   Multiply the plate's basic movement vector odx,ody by the age modifier
X   MR[], then add the remainders rx,ry from the last move to get some large
X   integers.  Then turn the large integers into small ones by dividing by
X   REALSCALE and putting the remainders back into rx,ry.  This function also
X   increases the age of each plate, but doesn't let the age of any plate
X   go above MAXLIFE.  This is done to make sure that MR[] does not need to
X   be a long vector. */
X 
X   register int i, a;
X 
X   /* If there is only a single supercontinent, anchor it */
X   for (i=phead, a=0; i; i=p[i].next, a++);
X   if (a == 1) if (p[phead].odx || p[phead].ody) {
X      p[phead].odx = 0; p[phead].ody = 0; }
X
X   for (i=phead; i; i=p[i].next) {
X      a = (p[i].odx * MR[p[i].age]) + p[i].rx;
X      p[i].dx = a / REALSCALE; p[i].rx = a % REALSCALE;
X      a = (p[i].ody * MR[p[i].age]) + p[i].ry;
X      p[i].dy = a / REALSCALE; p[i].ry = a % REALSCALE;
X      if (p[i].age < MAXLIFE-1) (p[i].age)++; } }
X 
X 
Xmove (src, dest) int src, dest; {
X   /* This function moves all the plates that are drifting.  The amount to
X   move by is determined in newdxy().  The function simply steps through
X   every square in the array; if there's a plate in a square, its new location
X   is found and the topography is moved there.  Overlaps between plates are
X   detected and recorded so that merge() can resolve the collision; mountain
X   growing is performed.  If two land squares wind up on top of each other,
X   folded mountains are produced.  If a land square winds up where ocean was
X   previously, that square is the leading edge of a continent and grows a
X   mountain by subsuming the ocean basin. */
X 
X   register int i, j; int a, b, c, x, y;
X 
X   /* Clear out the merge matrix and the destination arrays */
X   for (i=1; i<MAXPLATE; i++) for (j=1; j<MAXPLATE; j++) mm[i][j] = 0;
X   for (i=0; i<XSIZE; i++) for (j=0; j<YSIZE; j++) {
X      m[dest][i][j] = 0; t[dest][i][j] = 0; }
X 
X   checkmouse ();
X   /* Look at every square which belongs to a plate */
X   for (i=0; i<XSIZE; i++) for (j=0; j<YSIZE; j++) if ((a = m[src][i][j]) > 0) {
X 
X      /* Add the plate's dx,dy to the position to get the square's new */
X      /* location; if it is off the map, throw it away */
X      x = p[a].dx + i; y = p[a].dy + j;
X      if ((x >= XSIZE) || (x < 0) || (y >= YSIZE) || (y < 0)) p[a].area--;
X      else { /* It IS on the map */
X
X         /* If the destination is occupied, remove the other guy but */
X         /* remember that the two plates overlapped; set the new height */
X         /* to the larger height plus half the smaller. */
X         if (c = m[dest][x][y]) {
X            (mm[a][c])++; (p[c].area)--;
X            b = t[src][i][j]; c = t[dest][x][y];
X            t[dest][x][y] = (b > c) ? b + (c>>1) : c + (b>>1); }
X
X         /* The destination isn't occupied.  Just copy the height. */
X         else t[dest][x][y] = t[src][i][j];
X
X         /* If this square is over ocean, increase its height. */
X         if (t[src][x][y] < ZCOAST) t[dest][x][y] += ZSUBSUME;
X
X         /* Plate A now owns this square */
X         m[dest][x][y] = a; } } }
X 
X 
Xmerge (dest) int dest; {
X   /* Since move has set up the merge matrix, most of the work is done.  This
X   function calls bump once for each pair of plates which are rubbing; note
X   that a and b below loop through the lower diagonal part of the matrix.
X   One subtle feature is that a plate can bump with several other plates in
X   a step; suppose that the plate is merged with the first plate it bumped.
X   The loop will try to bump the vanished plate with the other plates, which
X   would be wrong.  To avoid this, the lookup table lut is used to provide
X   a level of indirection.  When a plate is merged with another, its lut
X   entry is changed to indicate that future merges with the vanished plate
X   should be applied to the plate it has just been merged with. */
X 
X   char lut[MAXPLATE]; int a, aa, b, bb, i;
X 
X   for (a=1; a<MAXPLATE; a++) lut[a] = a;
X   for (a=2; a<MAXPLATE; a++) for (b=1; b<a; b++) if (mm[a][b] || mm[b][a]) {
X      aa = lut [a]; bb = lut[b];
X      if (aa != bb) if (bump (dest, aa, bb)) {
X         lut[aa] = bb;
X         for (i=1; i<MAXPLATE; i++) if (lut[i] == aa) lut[i] = bb; } } }
X 
X 
Xbump (dest, a, b) int dest, a, b; {
X   /* Plates a and b have been moved on top of each other by some amount;
X   alter their movement rates for a slow collision, possibly merging them.
X   The collision "strength" is a ratio of the area overlap (provided by
X   move ()) to the total area of the plates involved.  A fraction of each
X   plate's current movement vector is subtracted from the movement vector
X   of the other plate.  If the two vectors are now within some tolerance
X   of each other, they are almost at rest so merge them with each other. */
X 
X   double maa, mab, ta, tb, rat, area; register int i, j, x;
X 
X   checkmouse();
X   /* Find a ratio describing how strong the collision is */
X   x = mm[a][b] + mm[b][a]; area = p[a].area + p[b].area;
X   rat = x / (MAXBUMP + (area / 20)); if (rat > 1.0) rat = 1.0;
X 
X   /* Do some math to update the move vectors.  This looks complicated */
X   /* because a plate's actual movement vector must be multiplied by */
X   /* MR[age], and because I have rewritten the equations to maximize */
X   /* use of common factors.  Trust me, it's just inelastic collision. */
X   maa = p[a].area * MR[p[a].age]; mab = p[b].area * MR[p[b].age];
X   ta = MR[p[a].age] * area;
X   p[a].odx = (p[a].odx * maa + p[b].odx * mab * rat) / ta;
X   p[a].ody = (p[a].ody * maa + p[b].ody * mab * rat) / ta;
X   tb = MR[p[b].age] * area;
X   p[b].odx = (p[b].odx * mab + p[a].odx * maa * rat) / tb;
X   p[b].ody = (p[b].ody * mab + p[a].ody * maa * rat) / tb;
X 
X   /* For each axis, compute the remaining relative velocity.  If it is */
X   /* too large, return without merging the plates */
X   if (ABS (p[a].odx*MR[p[a].age] - p[b].odx*MR[p[b].age]) > BUMPTOL) return(0);
X   if (ABS (p[a].ody*MR[p[a].age] - p[b].ody*MR[p[b].age]) > BUMPTOL) return(0);
X 
X   /* The relative velocity is small enough, so merge the plates.  Replace */
X   /* all references to a with b, free a, and tell merge() a was freed. */
X   for (i=0; i<XSIZE; i++) for (j=0; j<YSIZE; j++)
X      if (m[dest][i][j] == a) m[dest][i][j] = b;
X   p[b].area += p[a].area; pfree (a);
X   return (a); }
X 
X 
Xerode (dest) int dest; {
X   /* This function takes the topography in t[dest] and smooths it, lowering
X   mountains and raising lowlands and continental margins.  It does this by
X   stepping across the entire array and doing a computation once for each
X   pair of 8-connected pixels.  The computation is done by onerode(), below.
X   The computation result for a pair is a small delta for each square, which
X   is summed in the array e.  When the computation is finished, the delta
X   is applied; if this pushes an ocean square high enough, it is added to
X   an adjacent plate if one can be found.  Also, if a land square is eroded
X   low enough, it is turned into ocean and removed from its plate. */
X 
X   register int i, j, x, z, xx;
X 
X   /* Zero out the array for the deltas first */
X   for (i=0; i<XSIZE; i++) for (j=0; j<YSIZE; j++) e[i][j] = 0;
X   checkmouse();
X
X   /* Step across the entire array; each pixel is adjacent to 8 others, and */
X   /* it turns out that if four pairs are considered for each pixel, each */
X   /* pair is considered exactly once.  This is important for even erosion */
X   for (i=1; i<XSIZE; i++) for (j=1; j<YSIZE; j++) {
X      onerode (dest, i, j,  0, -1);
X      onerode (dest, i, j, -1, -1);
X      onerode (dest, i, j, -1,  0);
X      if (i < XSIZE-1) onerode (dest, i, j, -1, 1); }
X 
X   /* Now go back across the array, applying the delta values from e[][] */
X   for (i=0; i<XSIZE; i++) for (j=0; j<YSIZE; j++) {
X      z = t[dest][i][j] + e[i][j]; if (z < 0) z = 0; if (z > 255) z = 255;
X 
X      /* If the square just rose above shelf level, look at the four */
X      /* adjacent squares.  If one is a plate, add the square to that plate */
X      if ((z >= ZSHELF) && (t[dest][i][j] < ZSHELF)) { xx = 0;
X         if (i > 1)       if (x = m[dest][i-1][j]) xx = x;
X         if (i < XSIZE-1) if (x = m[dest][i-1][j]) xx = x;
X         if (j > 1)       if (x = m[dest][i][j-1]) xx = x;
X         if (j < YSIZE-1) if (x = m[dest][i][j+1]) xx = x;
X         if (xx) { p[xx].area++; m[dest][i][j] = xx; } }
X 
X      /* Add the increment to the old value; if the square is lower than */
X      /* shelf level but belongs to a plate, remove it from the plate */
X      t[dest][i][j] = z;
X      if ((z < ZSHELF) && (x = m[dest][i][j])) {
X         p[x].area--; m[dest][i][j] = 0; } } }
X
X
Xonerode (dest, i, j, di, dj) int dest, i, j, di, dj; {
X   /* This function is called once per pair of squares in the array.  The
X   amount each square loses to an adjacent but lower square in each step is
X   one-eighth the difference in altitude.  This is coded as a shift right 3
X   bits, but since -1 >> 3 is still -1, the code must be repeated to avoid
X   shifting negative numbers.  Also, since erosion is reduced below the
X   waterline, if an altitude is lower than ZERODE, ZERODE is used instead. */
X
X   register int z, t1, t2;
X
X   t1 = t[dest][i][j]; t2 = t[dest][i+di][j+dj]; z = 0;
X   if ((t1 > t2) && (t1 > ZERODE)) {
X      if (t2 < ZERODE) t2 = ZERODE;
X      z = ((t1 - t2 + ERODERND) >> 3); }
X   else if ((t2 > t1) && (t2 > ZERODE)) {
X      if (t1 < ZERODE) t1 = ZERODE;
X      z = -((t2 - t1 + ERODERND) >> 3); }
X   if (z) { e[i][j] -= z; e[i+di][j+dj] += z; } }
END_OF_src/tec3.c
if test 14070 -ne `wc -c <src/tec3.c`; then
    echo shar: \"src/tec3.c\" unpacked with wrong size!
fi
# end of overwriting check
fi
if test -f src/globe.c -a "${1}" != "-c" ; then 
  echo shar: Will not over-write existing file \"src/globe.c\"
else
echo shar: Extracting \"src/globe.c\" \(3176 characters\)
sed "s/^X//" >src/globe.c <<'END_OF_src/globe.c'
X/* This program is Copyright (c) 1991 David Allen.  It may be freely
X   distributed as long as you leave my name and copyright notice on it.
X   I'd really like your comments and feedback; send e-mail to
X   allen@viewlogic.com, or send us-mail to David Allen, 10 O'Moore
X   Avenue, Maynard, MA 01754. */
X
X/* This file contains a spinning globe animation */
X
X#include "const.h"
X#define N ((MAXX / 2) - 1)
X#define AREA (MAXX * MAXX)
X#define RADIUS 45
X#define PI 3.14159
X#define TWOPI (2.0 * PI)
X#define HALFPI (PI / 2.0)
X#define SET(n,x,y,i,j) m[n][(x)+N][(y)+N] = m[0][(i)+N][(j)+N]
X
Xunsigned char m[21][MAXX][MAXY];
Xint rx[] = {0, 0, 0, 0, 0}, ry[] = {0, N, 0, N, 0};
Xint xofs, prime, xflip, sflip, sxflip, change[1], step = 0;
Xint MAXSTEP = 10000, PRINTMODE = PRINTMODE_SCALE, XSIZE = MAXX, YSIZE = MAXX;
Xint ZCOAST = 0;
X
X
Xinit (s) char *s; { int i;
X   fileinit ();
X   if (!s || !(*s)) panic ("no valid filename"); else getparams (s);
X   for (i=1, xofs=45; i<6; xofs-=10, i++) {
X      prime = 16-i; xflip = i+5;  sflip = 15+i; sxflip = 6-i; compute (); } }
X
X
Xonestep () { static int i = 0;
X   draw (DRAW_TEC, LINE_NONE, m[i], 0); i++; step++;
X   if (i == 5) i = 6; if (i == 16) i = 17; if (i == 21) i = 1; }
X
X
Xmainpar (s) char *s; { int vec [AREA], index = AREA; register int i, j, k;
X
X   if (CMP ("LAND")) {
X      getlvec (&index, vec, 0);
X      for (k=0, j=0; j<MAXX; j++) for (i=0; i<MAXX; i++, k++)
X         m[0][i][j] = (vec[k] > 16) ? vec[k] : 16;
X      return (1); }
X   return (0); }
X
X
Xcompute () {
X   register int i, j, top; int x, y;
X
X   rx[0] = xofs + N; rx[4] = rx[0]; rx[2] = xofs - N;
X   for (j=0; j<RADIUS; j++) {
X      checkmouse ();
X      top = sqrt ((double) (RADIUS * RADIUS - j * j + 0.5));
X      for (i=1-top; i<top; i++) {
X         point (i, j, &x, &y);
X         SET (prime,   i,  j,  x,  y); SET (prime,   i, -j,  x, -y);
X         SET (xflip,  -i,  j, -x,  y); SET (xflip,  -i, -j, -x, -y);
X         permute (x, y, &x, &y);
X         SET (sxflip,  i, -j, -x,  y); SET (sxflip,  i,  j, -x, -y);
X         SET (sflip,  -i, -j,  x,  y); SET (sflip,  -i,  j,  x, -y); } } }
X
X
Xpoint (x, y, i, j) int x, y, *i, *j; {
X   double r, theta, topx, topy, botx, boty; int bot, top;
X
X   r = asin (sqrt ((double)(x * x + y * y)) / RADIUS) * 2.0 / PI;
X   if (x || y) theta = atan2 ((double)y, (double)x);
X   else theta = 0;
X   if (theta <= 0) theta += TWOPI;
X   bot = 0; while (theta > HALFPI) { bot++; theta -= HALFPI; }
X   top = bot + 1; theta = theta / HALFPI;
X   mapr (bot, r, &botx, &boty); mapr (top, r, &topx, &topy);
X   *i = botx + theta * (topx - botx) + 0.5;
X   *j = boty + theta * (topy - boty) + 0.5;
X   if (*i >= N) *i = N-1; if (*j >= N) *j = N-1; }
X
X
Xmapr (n, r, x, y) int n; double r, *x, *y; {
X   *x = r * (rx[n] - xofs) + xofs;
X   *y = r * ry[n];
X   if (*x >= N) { *y = *x - N; *x = N; if (n == 4) *y = -(*y); } }
X
X
Xpermute (x, y, i, j) int x, y, *i, *j; {
X   if ((x >= 0) && (y >= 0)) { *i =  N - y - 1; *j =  x - N + 0; }
X   if ((x <  0) && (y >= 0)) { *i =  y - N + 0; *j = -N - x - 1; }
X   if ((x >= 0) && (y <  0)) { *i =  y + N + 0; *j =  N - x - 1; }
X   if ((x <  0) && (y <  0)) { *i = -N - y - 1; *j =  x + N + 0; } }
X
X
Xgreyscale (x) int x; { }
END_OF_src/globe.c
if test 3176 -ne `wc -c <src/globe.c`; then
    echo shar: \"src/globe.c\" unpacked with wrong size!
fi
# end of overwriting check
fi
echo shar: End of archive 4 \(of 4\).
cp /dev/null ark4isdone
MISSING=""
for I in 1 2 3 4 ; do
    if test ! -f ark${I}isdone ; then
	MISSING="${MISSING} ${I}"
    fi
done
if test "${MISSING}" = "" ; then
    echo You have unpacked all 4 archives.
    rm -f ark[1-9]isdone
else
    echo You still need to unpack the following archives:
    echo "        " ${MISSING}
fi
##  End of shell archive.
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.