Path: utzoo!utgpu!news-server.csri.toronto.edu!cs.utexas.edu!uunet!sparky!kent From: allen@viewlogic.com (Dave Allen) Newsgroups: comp.sources.misc Subject: REPOST: v18i003: planet - planet generation simulator, Part03/04 Message-ID: <1991Apr9.041754.8090@sparky.IMD.Sterling.COM> Date: 9 Apr 91 04:17:54 GMT Sender: kent@sparky.IMD.Sterling.COM (Kent Landfield) Organization: Sterling Software, IMD Lines: 1237 Approved: kent@sparky.imd.sterling.com X-Checksum-Snefru: 9e87f929 1fb6a493 125dafc0 59ca8a0d Submitted-by: Dave Allen Posting-number: Volume 18, Issue 3 Archive-name: planet/part03 Supersedes: tec: Volume 10, Issue 77-78 #! /bin/sh # into a shell via "sh file" or similar. To overwrite existing files, # type "sh file -c". # The tool that generated this appeared in the comp.sources.unix newsgroup; # send mail to comp-sources-unix@uunet.uu.net if you want that tool. # Contents: ./doc/about.tec ./src/clim.c ./src/heat.c ./src/pressure.c # ./src/tec1.c ./src/x.c # Wrapped by kent@sparky on Mon Apr 8 22:39:15 1991 PATH=/bin:/usr/bin:/usr/ucb ; export PATH echo If this archive is complete, you will see the following message: echo ' "shar: End of archive 3 (of 4)."' if test -f './doc/about.tec' -a "${1}" != "-c" ; then echo shar: Will not clobber existing file \"'./doc/about.tec'\" else echo shar: Extracting \"'./doc/about.tec'\" \(10300 characters\) sed "s/^X//" >'./doc/about.tec' <<'END_OF_FILE' XI. MOTIVATION AND APPROACH X XI wanted to write a program to generate maps of many imaginary worlds using Xa microcomputer. Fractal techniques produce unacceptably random-looking maps. XIn other related fields, people have used simplified models of physical Xprocesses and then tweaked the model parameters to give results that "looked Xright". Dole [1] simulated accretion of planets around stars. Hart [2] Xsimulated the evolution of a planetary atmosphere. Forrester [3] simulated Xpopulation dynamics for near-future Earth. My approach was to find a physical Xmodel for continent formation, and then write a program to implement that model. X XGeologists have proposed [4] that several times in history, the continents Xof Earth have drifted together to form a single "supercontinent", which later Xbroke apart. A recent paper [5] suggested a physical model for this X"supercontinent cycle", and I used that model for my program. X XTheir hypothesis runs like this. (Disclaimer - I don't know anything about Xgeology except what I got from reading a few papers.) The Earth generates Xheat by radioactive decay; this heat is conducted away more slowly by Xcontinental crust than by oceanic crust. Heat accumulates under large Xcontinents, causing them to dome upwards and eventually break apart. After Xa certain amount of heat has escaped through the new ocean basin, the Xcontinental fragments are drawn back together by the subsiding basin. XOver long periods of time, many oceans should open and close in the center Xof the land area, while the original superocean should remain about the same. X XAny map of the Earth shows evidence to support this hypothesis. The Atlantic Xand Mediterranean oceans are fairly recent, while the Pacific could be Xconsidered the superocean. The Himalayas, the Pyrenees, and the Ukraine Xmountain ranges could have been formed by collisions between continental Xfragments. The "Ring of Fire" around the Pacific represents subduction zones Xwhere new mountains, such as the Rockies, are forming as the continents drift Xfurther apart. The older, lower Appalachian mountains might have been created Xduring a collision between Europe and North America in the last supercontinent Xcycle. X X XII. REPRESENTATION X XOne data structure well-suited for storing topographic maps is a simple Xtwo-dimensional array, where each array element represents the height Xof the terrain at that point. However, the Earth is spherical and a sphere Xcan't be mapped onto a plane without some distortion. After trying many Xmapping methods, I came back to the simplest: I pretend Columbus was wrong, Xand the world really is flat. Since the continents in the physical model Xdrift toward and away from some central point, it is reasonable, if the map Xis large enough, to simply throw away the land area that falls off the Xedge of the map. I used a simple square, flat surface with no "wraparound." X XI could, and probably will at some point, write a couple of paragraphs about Xwhy doing plate movement on a sphere is so hard. All of the methods I Xtried required too much data storage for each point, or too much real Xnumber computation, or resulted in unacceptable distortion. The criterion XI used was to place a small continent, say a 5 x 5 square, anywhere on Xthe surface. It should be able to travel one great circle in any direction Xand remain undisturbed though the entire journey. The method should not Xrequire _any_ floating point computation once the initial velocity is Xcomputed since I want to do many iterations on a microcomputer. X XA second representation issue is how to cause plates to drift apart and Xback together. Maintaining the detailed structure of the ocean basin is Xcomplicated; consider a circular landmass which splits along a diameter. XThe ocean basin structure here is clear, since the rift remains in place and Xthe new basin extends out from it in both directions. Now suppose one of Xthe semicircular fragments splits in half perpendicular to the first split. XThe motion of the landmasses can be computed easily, but the new rift must Xdrift as well and the required motions in the basin seem self-contradictory. X XI chose to ignore the ocean basin and concentrate on the landmasses. Thus, Xrifts are not recorded in any way and the ocean floor does not move. Since Xthe basins are not simulated, some mechanism is needed to replace the physical Xmechanism of the sinking basin sucking the landmasses back together. This Xis done by a velocity scaling profile which is coded directly into the Xmovement routine. The profile says that after a short time, forward motion Xof a landmass slows and then stops; after this point the landmass begins to Xmove backwards towards its origin. X X XIII. ALGORITHM X XThe initial supercontinent is generated by a simple fractal algorithm; a Xfractal mountain is created by superimposing four two-dimensional fractal Xlines. The mountain is then turned into a binary blob, that is, everything Xabove a certain altitude is declared to be land, and everything else becomes Xwater. The blob is then "improved" by removing islands and internal oceans. XThis technique produces very irregular circular blobs of varying sizes. X XThe program then repeats a sequence of operations for each timestep. First, Xa rift may be generated, splitting some continent into pieces; new velocities Xare generated for the new landmasses. Second, the distance each landmass will Xmove this timestep is determined by applying the velocity profile described Xabove to the landmass velocity. Third, the landmasses are moved and mountains Xare built. Fourth, where the movement routine has caused landmasses to Xoverlap, new velocities are computed for the overlapping landmasses, and they Xmay be merged into a single landmass if they overlap sufficiently. Fifth, the Xentire topography is eroded; if a square is higher than the square next to Xit, the altitudes are made more equal. Finally, the screen is redrawn. X XIn step one, random coordinates are generated until a location is found Xwhich is at least five squares from any ocean. This location will be the Xcenter of a rift, which is represented as a randomly curving line. The Xdistance requirement is there to ensure that if a rift occurs, it occurs Xon a fairly large continent where internal heat would accumulate. When Xthe rift is drawn, it is allowed to bend randomly but is drawn using Xconventional digital line-drawing techniques (i.e. Bresenham's algorithm). XThe rift drawing routine terminates when the rift hits ocean. If the growing Xrift hits another landmass, which could occur if two landmasses were touching Xbut had not yet merged, the rift is aborted since the result would be a very Xunnatural looking U-shaped continent. X XOnce the rift is created, a segmentation algorithm determines how Xmany new landmasses were formed. There should only be two new landmasses, Xbut splinters can be formed if the rift runs along a narrow peninsula or Xtoo close to a coast. The algorithm detects and erases such splinters. XThe landmasses are given an initial velocity which makes them drift directly Xaway from the rift. The velocities are inversely dependent on the area Xof the landmass so that smaller masses move faster. After the initial Xvelocities have been determined, the rift is erased. X XThe second and third steps are concerned with moving the landmasses. This Xis done with a slightly simpler version of Bresenham's algorithm. The velocity Xprofile described in Section II is applied so that landmasses will drift Xapart at first, but will then slow down and come back together. Because a Xlandmass can split up again, the process does not simply reassemble the Xoriginal supercontinent when they drift back together. X XMountains are built under two conditions: where two landmasses collide, and Xwhere landmass subsumes ocean basin. Both conditions are detected and resolved Xin the movement routine. Movement is performed from a source array to a Xdestination array. A loop steps over every square in the source array. If Xthere is land in a square, its new position is computed based on what Xlandmass the square belongs to. If there is ocean basin there in the source Xarray, the square belongs to the leading edge of a landmass which is subsuming Xocean basin. A constant is added to the height of the land in that square. XIf there is land there in the destination array, then the two landmasses are Xcolliding. There is a collision array which records how many times each Xpair of landmasses collide in any given timestep, and the movement routine Xincrements the appropriate element of this array for each collision. XMountains are built in this case by adding half the height of the lower Xterrain to the height of the higher terrain. X XStep four uses the collision information provided by the movement routine Xto adjust the velocities of colliding landmasses and perhaps merge them. XFor each pair of colliding landmasses, the routine determines how "strong" Xthe collision is this timestep. The velocities are adjusted, again in Xinverse proportion to the landmass areas, so that the continents begin Xto come to rest with respect to each other. If the relative velocity of Xthe two landmasses, after adjustment, is small enough, then the masses Xare merged. This is the mechanism that reassembles the continental Xfragments into a supercontinent. X XErosion is the final, and most computation-intensive, step of the simulation. XEach pair of adjacent squares is considered exactly once. If either square Xis part of a landmass, erosion occurs. Some fraction of the difference between Xthe two squares' altitudes is subtracted from the taller and added to the Xshorter. This is a very simple algorithm since it ignores many factors, Xbut more detail would make the simulation even slower. X X XIV. REFERENCES X X1. Dole, Stephen H., "Computer Simulation of the Formation of Planetary X Systems", Icarus 13 (1970), pp 494-508. X2. Hart, Michael H., ""The Evolution of the Atmosphere of the Earth", Icarus X 33 (1978), pp 23-39. X3. Forrester, J.W., World Dynamics, Wright-Allen Press (1973). X4. Wilson, J. Tuzo, "Continental Drift", Scientific American 208, 4 (April X 1963), pp 86-100. X5. Nance, R.D., Worsley, T.R., and Moody, J.B., "The Supercontinent Cycle", X Scientific American (July 1988), pp 72-79. END_OF_FILE if test 10300 -ne `wc -c <'./doc/about.tec'`; then echo shar: \"'./doc/about.tec'\" unpacked with wrong size! fi # end of './doc/about.tec' fi if test -f './src/clim.c' -a "${1}" != "-c" ; then echo shar: Will not clobber existing file \"'./src/clim.c'\" else echo shar: Extracting \"'./src/clim.c'\" \(5481 characters\) sed "s/^X//" >'./src/clim.c' <<'END_OF_FILE' 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 all of the general-purpose routines which are not X machine specific: init, mainpar (called by fileio), onestep (called X by the machine-specific control routine), and the range-finding X function called by several computation routines. */ X X#include "const.h" X#include "clim.h" X X/* L holds the land values everybody starts with; lm is an edge map with X continents outlined in black and mountains outlined in white. */ Xunsigned char l[MAXX][MAXY], lm[MAXX][MAXY]; X X X/* These are the parameters used by all of the climate functions. They are X described in params.doc. */ Xint XSIZE = MAXX, YSIZE = MAXY, BSIZE = MAXB; Xint PRINTMODE = PRINTMODE_SHORT, MAXRANGE = 15; Xint MAXSTEP = 10000, ZCOAST = 0; X Xint change[] = { 1, 1, 1, 1, 1, 1 }, step = 0; Xextern int picktype; X X Xinit (s) char *s; { X if (s && *s) getparams (s); X fileinit (); X heatcomp (); presscomp (); windcomp (); X raincomp (); climcomp (); } X X Xonestep () { X switch (picktype) { X case M_HEAT: heatdraw (step % MAXB); break; X case M_PRESS: pressdraw (step % MAXB); break; X case M_WIND: windraw (step % MAXB); break; X case M_RAIN: raindraw (step % MAXB); break; X case M_CLIM: climdraw (); break; } } X X Xmainpar (s) char *s; { X /* This function is called by getparams() in fileio.c; it looks at each X parameter name read from the input file. If it recognizes the name, X the right function is called to set that parameter. The function X returns true if it recognizes the parameter name, and false otherwise. */ X X if (CMP ("LAND")) getland (); X else if (CMP ("BSIZE")) getlng (&BSIZE, M_MAIN); X else if (CMP ("MAXRANGE")) getlng (&MAXRANGE, M_MAIN); X else if (CMP ("PRINTMODE")) getlng (&PRINTMODE, M_MAIN); X else if (heatpar (s)) { } X else if (presspar (s)) { } X else if (windpar (s)) { } X else if (rainpar (s)) { } X else if (climpar (s)) { } X else return (0); X return (1); } X X Xgetland () { register int i, j, x; X /* This function is called by mainpar, above, when the LAND parameter X is encountered. It just calls getmat in fileio.c to do the work, but X then it also creates an edge map lm; this is used by a couple of the X drawing routines as background. */ X X getmat (&XSIZE, &YSIZE, M_MAIN, l); X X for (j=0; j 0) != (l[i-1][j] > 0)) x |= LINE_0V; X if (j) if ((l[i][j] > 0) != (l[i][j-1] > 0)) x |= LINE_0H; X X /* If both black and white lines are present, white wins */ X if ((x & LINE_0V) && (x & LINE_1V)) x &= (~LINE_0V); X if ((x & LINE_0H) && (x & LINE_1H)) x &= (~LINE_0H); X lm[i][j] = x; } } X X Xrange (rr) char rr[MAXX][MAXY]; { X /* This function is called by a number of climate routines. It takes an X input array with blobs of -1's on a background of 0's. The function winds X up replacing each 0 with the distance from that square to the nearest -1. X The function onerange() does all the work, but it will not compute ranges X greater than MAXRANGE. Therefore, after onerange() is called, any remaining X 0 values must be replaced with MAXRANGE, indicating that that square is X "very far" from any -1 value. */ X X register int i, j; X X onerange (rr); checkmouse (); X for (j=0; j'./src/heat.c' <<'END_OF_FILE' 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 the routines that compute local temperatures */ X X#include "const.h" X#include "clim.h" X X/* Some private defines. FREEZING is 0 degrees C in Kelvin; DEG2RAD is the X conversion factor from angular degrees to radians. */ X X#define FREEZING 273.0 X#define PI 3.14159 X#define DEG2RAD (PI / 180) X X/* The input array is l, from main.c; lm is used by heatdraw(). The output X array is ts, containing temperatures. Array t is an unscaled copy of the X temperatures; tmin and tmax are used for scaling, and tscale is the X computed scale factor. To convert the unscaled temperatures to degrees K, X divide by TEMPSCALE. */ X Xextern unsigned char l[MAXX][MAXY], lm[MAXX][MAXY]; Xunsigned char ts[MAXB][MAXX][MAXY]; Xint tt[MAXB][MAXX][MAXY], tmax, tmin; Xdouble tscale; int TEMPSCALE = 10; X X/* These parameters are defined in main.c */ X Xextern int BSIZE, XSIZE, YSIZE, PRINTMODE; X X/* These parameters are used here, and are described in params.doc */ X Xint PRINTEMP = 0; Xdouble TILT = 23.0, ECCENT = 0.0, ECCPHASE = 0.0; Xdouble LCOS = 45.0, LCONST = 275.0, LTILT = 1.0, LSMOOTH = 0.6, LDIV = 180.0; Xdouble OCOS = 30.0, OCONST = 275.0, OTILT = 0.2, OSMOOTH = 0.2, ODIV = 250.0; X X Xheatpar (s) char *s; { X /* This function is called by mainpar() in main.c; it simply tests input X parameter names to see if they are defined in this file. Each of the X above ints are defined in this file. If the input string matches here, X the function returns true. */ X if (CMP ("TILT")) getdbl (&TILT, M_HEAT); X else if (CMP ("ECCENT")) getdbl (&ECCENT, M_HEAT); X else if (CMP ("ECCPHASE")) getdbl (&ECCPHASE, M_HEAT); X X else if (CMP ("LCOS")) getdbl (&LCOS, M_HEAT); X else if (CMP ("LCONST")) getdbl (&LCONST, M_HEAT); X else if (CMP ("LTILT")) getdbl (<ILT, M_HEAT); X else if (CMP ("LSMOOTH")) getdbl (&LSMOOTH, M_HEAT); X else if (CMP ("LDIV")) getdbl (&LDIV, M_HEAT); X X else if (CMP ("OCOS")) getdbl (&OCOS, M_HEAT); X else if (CMP ("OCONST")) getdbl (&OCONST, M_HEAT); X else if (CMP ("OTILT")) getdbl (&OTILT, M_HEAT); X else if (CMP ("OSMOOTH")) getdbl (&OSMOOTH, M_HEAT); X else if (CMP ("ODIV")) getdbl (&ODIV, M_HEAT); X X else if (CMP ("PRINTEMP")) getlng (&PRINTEMP, M_HEAT); X else return (0); X return (1); } X X Xheatlut (x, s) int x; char *s; { X /* After the heatcomp routine is finished, a scale factor is computed X for converting degrees K to 0..255; this routine converts back. Functions X in the machine-dependent code can call here to print map keys. The caller X must provide a char buffer; a string containing degrees F is put there. */ X X double temp; X X temp = (((double) x / tscale) + (double) tmin) / TEMPSCALE; X temp = (temp - FREEZING) * 1.8 + 32; X sprintf (s, "%6.1f", temp); } X X Xheatcomp () { X /* This is the main routine for computing temperatures. After getheat() X is called to do all the work, this routine takes the ints from t and X finds the smallest and largest values. These are used to compute a scale X factor, tscale; the arrays ts are then filled with scaled values. Finally, X putmat() from main.c is called if needed to print results. */ X X register int i, j, buf; char s[20]; X X getheat (); tmin = 32000; tmax = 0; X X usermessage ("Scaling heat"); X /* Find minimum and maximum across all buffers */ X for (buf=0; buf tmax) tmax = tt[buf][i][j]; } } X X /* Compute scale; for every buffer, fill ts from t */ X tscale = 254.0 / ((double) (tmax - tmin)); X for (buf=0; buf 2.0 * PI) phase -= (2 * PI); X fact = (1.0 + ECCENT * cos (phase)) * TEMPSCALE; X x = (lat + cos (theta) * TILT * LTILT) * lscl; X tland[buf] = (LCONST + LCOS * cos (x)) * fact; X x = (lat + cos (theta) * TILT * OTILT) * sscl; X tsea[buf] = (OCONST + OCOS * cos (x)) * fact; } X for (i=0; i= YSIZE) jmax = YSIZE-1; X for (j1=jmin; j1<=jmax; j1++) for (i0=-5; i0<6; i0++) { X i1 = i0 + x; if (i1 < 0) i1 += XSIZE; X if (i1 >= XSIZE) i1 -= XSIZE; X sum += l[i1][j1]; } X return (sum); } END_OF_FILE if test 6893 -ne `wc -c <'./src/heat.c'`; then echo shar: \"'./src/heat.c'\" unpacked with wrong size! fi # end of './src/heat.c' fi if test -f './src/pressure.c' -a "${1}" != "-c" ; then echo shar: Will not clobber existing file \"'./src/pressure.c'\" else echo shar: Extracting \"'./src/pressure.c'\" \(9001 characters\) sed "s/^X//" >'./src/pressure.c' <<'END_OF_FILE' 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 the routines to compute global high and low X pressure areas. */ X X#include "const.h" X#include "clim.h" X X/* These are the data arrays required: l is an input array defining the X ocean, land and mountain areas; ts is the temperature array computed X by the functions in heat.c. Array pr is filled by ocean(), land() and X heateq(), below; PR_HIGH indicates a high pressure zone, PR_LOW indicates X a low, and PR_HEQ indicates the heat equator, a low zone. Array pm is X the output array for this file, and it contains an edge map which has X edges in color 1 surrounding lows and color 0 around highs. Array r is X temporary storage for local calls to range() in main.c. */ X Xextern unsigned char l[MAXX][MAXY], ts[MAXB][MAXX][MAXY]; Xstatic char r[MAXX][MAXY]; Xstatic unsigned char pm[MAXB][MAXX][MAXY]; Xunsigned char pr[MAXB][MAXX][MAXY]; X X X/* The externs below are declared in main.c; they represent global parameters. X The ints are described in params.doc. */ X Xextern int BSIZE, XSIZE, YSIZE, PRINTMODE; Xint OOTHRESH = 5, OLTHRESH = 1, LOTHRESH = 3, LLTHRESH = 7; Xint OHMIN = 130, OHMAX = 180, OLMIN = 40, OLMAX = 65; Xint LHMIN = 0, LHMAX = 20, LLMIN = 220, LLMAX = 255; Xint PRINTPR = 0; X X Xpresspar (s) char *s; { X /* This function is called by mainpar() in main.c; it simply tests input X parameters to see if they are defined in this file. Each of the above X ints is referred to in this function. If an input string matches here, X the function returns true. */ X if (CMP ("OOTHRESH")) getlng (&OOTHRESH, M_PRESS); X else if (CMP ("OLTHRESH")) getlng (&OLTHRESH, M_PRESS); X else if (CMP ("OHMIN")) getlng (&OHMIN, M_PRESS); X else if (CMP ("OHMAX")) getlng (&OHMAX, M_PRESS); X else if (CMP ("OLMIN")) getlng (&OLMIN, M_PRESS); X else if (CMP ("OLMAX")) getlng (&OLMAX, M_PRESS); X X else if (CMP ("LOTHRESH")) getlng (&LOTHRESH, M_PRESS); X else if (CMP ("LLTHRESH")) getlng (&LLTHRESH, M_PRESS); X else if (CMP ("LHMIN")) getlng (&LHMIN, M_PRESS); X else if (CMP ("LHMAX")) getlng (&LHMAX, M_PRESS); X else if (CMP ("LLMIN")) getlng (&LLMIN, M_PRESS); X else if (CMP ("LLMAX")) getlng (&LLMAX, M_PRESS); X X else if (CMP ("PRINTPR")) getlng (&PRINTPR, M_PRESS); X else return (0); X return (1); } X X Xpresscomp () { X /* The main routine for this file. It just calls the four routines X which do all the work. Ocean() finds pressure extremes on the ocean; X land() does the same for land; heateq() defines the heat equator, and X setpm() computes pm[][] from pr[][]. */ X X int buf; X X for (buf=0; buf OLTHRESH) ? -1 : 0; X range (r); X X /* For each array element, if it is at least OOTHRESH squares from the X nearest big piece of land, it might be the center of an ocean pressure X zone. The pressure zones are defined by temperature ranges; if the X temperature in ts is between OLMIN and OLMAX, a low is recorded, while X if the temperature is between OHMIN and OHMAX, a high is recorded. */ X for (j=0; j OOTHRESH) { X if ((x >= OLMIN) && (x <= OLMAX)) pr[buf][i][j] = PR_LOW; X if ((x >= OHMIN) && (x <= OHMAX)) pr[buf][i][j] = PR_HIGH; } } } X X Xland (buf) int buf; { X /* This function is simply the complement of ocean(): it finds land highs X and lows. A land high or low must occur over land, far from major oceans. X Two calls to range() are made to find the qualifying land areas; then X temperature criteria are used to select the actual pressure zones. */ X X register int i, j; int x; X X /* Set r to distance on water from coast. */ X for (j=0; j LOTHRESH) ? -1 : 0; X range (r); X X /* For each array element, if it is at least LLTHRESH squares from the X nearest large ocean, it might be the center of a land pressure zone. X The pressure zones are defined by temperature ranges; if the temperature X in ts is between LLMIN and LLMAX, a low is recorded, while if the X temperature is between LHMIN and LHMAX, a high is recorded. */ X for (j=0; j LLTHRESH) { X if ((x >= LLMIN) && (x <= LLMAX)) pr[buf][i][j] = PR_LOW; X if ((x >= LHMIN) && (x <= LHMAX)) pr[buf][i][j] = PR_HIGH; } } } X X Xfindheq (buf) int buf; { X /* This function finds the heat equator and marks it in pr. For each X vertical column of ts, the median position is found and marked. To X make the heat equator continuous, jlast is set to the position of the X heat equator in the previous column; a connection is made in the present X column to ensure continuity. */ X X register int i, j; int sum, jlast = 0, jnext; X X for (i=0; i>=1, j=0; j0; j++) sum -= ts[buf][i][j]; X X /* Mark this position and remember it with jnext */ X pr[buf][i][j] = PR_HEQ; jnext = j; X X /* For each column except the first (where i = 0), if the last heat X equator is above this one, move upwards to it, marking each square, X to ensure continuity; if below this one, move downwards to it. */ X X if (i && (j > jlast)) for (; j>=jlast; j--) pr[buf][i][j] = PR_HEQ; X else if (i && (j < jlast)) for (; j<=jlast; j++) pr[buf][i][j] = PR_HEQ; X X /* Remember this position for the next column. Note that no check is X done to ensure continuity at the wraparound point; this is bad. */ X jlast = jnext; } } X X Xsetpm (buf) int buf; { X /* Setpm() is called after the above three functions have filled pr with X the codes for high, low and heat equator. The purpose of this function X is to create an edge map surrounding lows with color 1 and highs with X color 0. */ X X register int i, j, k; int col; X X for (j=0; j'./src/tec1.c' <<'END_OF_FILE' 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; i0; 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 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 0) for (i=1; i 0) for (i=1; i'./src/x.c' <<'END_OF_FILE' 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 X#include X#include X#include X#include X#include X#include X#include 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