This program will produce a raster map layer of uniform deviates whose
magnitude can be expressed by the user. It is essentialy the same as
r.gauss.surf, but uses a linear random number generator instead. Both 
random number generators are taken from Press, Flannery, Teukolsky and 
Vetterling (1988) - Numerical Recipes in C.

[JWO 23.10.91]
/*
** Code Compiled by Jo Wood [JWO] 23rd October 1991
** Midlands Regional Research Laboratory (ASSIST)
**
**
*/

/* %W% %G% */
#include "gis.h"

main(argc,argv) 
    int argc;
    char *argv[];
{

    /****** INITIALISE ******/

    struct Option 	*out;		/* Structures required for the G_parser()	*/
   					/* call. These can be filled with the		*/
    struct Option	*min;		/* various defaults, mandatory paramters	*/
					/* etc. for the GRASS user interface.		*/
    struct Option	*max;

    G_gisinit (argv[0]);	/*	This GRASS library function MUST
					be called first to check for valid
					database and mapset. As usual argv[0]
					is the program name. This can be
					recalled using G_program_name(). */


    /****** SET PARSER OPTIONS ******/

    out   = G_define_option(); 	/*	Request pointer to memory for each option.	*/
    min   = G_define_option();	/* 	Minimum random value				*/
    max   = G_define_option();	/* 	Maximum random value				*/

    out->key		= "out";
    out->description	= "Name of the random surface to be produced";
    out->type		= TYPE_STRING;
    out->required	= YES;

    min->key		= "min";
    min->description	= "Minimum random value";
    min->type		= TYPE_INTEGER;
    min->answer		= "0";

    max->key		= "max";
    max->description	= "Maximum random value";
    max->type		= TYPE_INTEGER;
    max->answer		= "100";

    if (G_parser(argc,argv))
	exit(-1);		/*	Returns a 0 if sucessful		*/



    /****** CHECK THE CELL FILE (OUT) DOES NOT ALREADY EXIST******/

    if (G_legal_filename(out->answer)==NULL)
    {
	char err[256];
	sprintf(err,"Illegal file name. Please try another.");
	G_fatal_error(err);
    }
    else
    {
	if (G_find_cell(out->answer,"") !=NULL)
	{
	    char err[256];
	    sprintf(err,"Raster map [%s] already exists.\nPlease try another.",out);
	    G_fatal_error(err);
	}

    }

    /****** CREATE THE RANDOM CELL FILE  ******/

    randsurf(out->answer,atoi(min->answer),atoi(max->answer));


}
/*****************/
/** randsurf()	**/
/*****************/

#include "gis.h"
#include <math.h>

randsurf(out,min,max)

char	*out;				/* Name of cell files to be opened.	*/
int	min,max;			/* Minimum and maximum cell values.	*/

{

	int	nrows,ncols;		/* Number of cell rows and columns	*/

	

	CELL	*row_out;		/* Buffer just large enough to hold one	*/
    					/* row of the raster map layer.		*/

	int	fd_out;			/* File descriptor - used to identify 	*/
					/* open cell files.			*/
	

	int	row_count,col_count;

	float	rand1();
		
	/****** INITIALISE RANDOM NUMBER GENERATOR ******/

	rand1(-1* getpid());

	/****** OPEN CELL FILES AND GET CELL DETAILS ******/

	fd_out = G_open_cell_new(out);
	
	nrows = G_window_rows();
	ncols = G_window_cols();

	row_out = G_allocate_cell_buf();


	/****** PASS THROUGH EACH CELL ASSIGNING RANDOM VALUE ******/

	for (row_count=0;row_count<nrows;row_count++)
	{
	    for (col_count=0;col_count<ncols;col_count++)
		*(row_out+col_count) = rint(rand1(2742)*(max-min)+min);	
	   
	    /* Write contents row by row */
	    G_put_map_row(fd_out,row_out);
	}


	/****** CLOSE THE CELL FILE ******/

	G_close_cell(fd_out);


}



/****************************************************************/
/***    Random Number Generator (Uniform Deviates 0.0 -> 1.0) ***/
/***                                                          ***/
/*** Based on three linear congruential generators:           ***/
/***  One for most significant part,                          ***/
/***  One for the least significant,                          ***/
/***  One for a shuffling routine (Knuth, 1981).              ***/
/***                                     [Press et al, 1988]  ***/
/***                                                          ***/
/*** Coded  Oct 23 1991                                       ***/
/*** Version 1.0                                              ***/
/***                                                          ***/
/****************************************************************/

			/* Arbitrary constants */
#define M1 259200	
#define IA1 7141
#define IC1 54773
#define RM1 (1.0/M1)
#define M2 134456
#define IA2 8121
#define IC2 28411
#define RM2 (1.0/M2)
#define M3 243000
#define IA3 4561
#define IC3 51349


float rand1(seed)
int 	seed;


{
	static long 	ix1,ix2,ix3;
	float		temp;
	int		j;
	static float	r[98];
	static int	iff=0;
	
	if (seed < 0 || iff==0) /* Initialise if negative seed is given. */
				 /* This MUST be done on the first call.  */
	{
		iff=1;
		ix1=(IC1-seed) % M1;		/* Seed the 1st routine */
		ix1=(IA1*ix1+IC1) % M1;	/* and use it to seed the */
		ix2=ix1 % M2;			/* 2nd and 3rd routines */
		ix1=(IA1*ix1+IC1) % M1;
		ix3=ix1 % M3;
		
		for (j=1; j<=97; j++)
		{
			ix1=(IA1*ix1+IC1) % M1;  /* Fill table with sequen- */
			ix2=(IA2*ix2+IC2) % M2;  /* tial uniform deviates   */
			r[j]=(ix1+ix2*RM2)*RM1;  /* from the 1st 2 routines */
		}
		seed=1;
	}

	ix1=(IA1*ix1+IC1) % M1;	/* Generate the next number in the */
	ix2=(IA2*ix2+IC2) % M2;	/* sequence of routine		   */
	ix3=(IA3*ix3+IC3) % M3;

	j=1 + ((97*ix3) / M3);  /* Get random position ion array */
	temp = r[j];
	r[j]=(ix1+ix2*RM2)*RM1;        /* Refil array and return entry */
	return(temp);		        
}