/*
|||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||

    R E A L    R A N D O M

|||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
*/

/********************************************************************
 * 
 *   Park, S.K. and Miller, K.W., 1988, "Random Number Generators:
 *   Good Ones are Hard to Find", Communications of the ACM, Vol 31, 
 *   No 10, pp.1192-1201.
 * 
 *  Adapted by John Walker 
 *      Converted to Class Functions -> Octogre 1997
 *      ALife class version -> February 1998
 *      Copyright (c) John F. Walker 1997
 *
 *******************************************************************/
    /*  ---RAW CODE---
     *		extern double seed;
     *		
     *		#define A 16807.0
     *		#define M 2147483647.0
     *		#define Q 127773.0
     *		#define R 2836.0
     *		
     *		double random()
     *		{
     *		  int hi;
     *		  double rand;
     *		  double lo,test;
     *		  hi = seed/Q;
     *		  lo = seed-Q*hi;
     *		  test = A*lo - R*hi;
     *		  if(test>0.0)
     *		    seed = test;
     *		  else
     *		    seed = test + M;
     *		  rand = seed/M;
     *		  return(rand);
     *		}
     */


/********************************************************************
 * 
 *  I N S T R U C T I O N S    -->  realRandomJW.h
 * 
 *  If you chose not to read the instructions you are doomed:
 *    DELETE THIS FILE NOW AND FORGET YOU EVER SAW IT.
 * 
 ********************************************************************
 * 
 *  To use these functions, include this header file into 
 *    all the files containing functions which will be calling them.
 *    (I recommend including it everywhere to be sure...) and compiling
 *    the RealRandomJW.cxx file as part of your project.  (Yes, this 
 *    could easily be made into a library.)
 *
 *    Remember to call the SeedRealRandom or TimeSeedRealRandom 
 *    function before calling anything else.  This is NOT checked 
 *    for since it would hurt performance.
 * 
 *    The functions maintain an internal, consistent seed so any
 *    and all functions will not get repeating values.  Only the
 *    inital seed declaration is needed to maintain randomness.
 * 
 * 
 *  SeedRealRandom( double `seed', void | int )
 *      Provide the seed value for the randomization function.
 *
 *  TimeSeedRealRandom( void | int )
 *      Randomly generate the seed value for the randomization function.
 *      This returns the seed generated.
 *
 *     >>>  ONE OF THE ABOVE IS REQUIRED          <<<
 *     >>>    BEFORE CALLING ANY OTHER FUNCTION ! <<<
 *
 *  GetRealRandomSeed( double `seed' )
 *      Returns the (double) value of the current seed state.
 *
 *  RealRandom( void )
 *      Returns a double value between 0.0 and .9999999999999999
 *                                               (not a joke)
 *
 *  RealRandomDouble( double maximum )
 *      Returns a double value between 0.0 and 
 *                                     maximum*.9999999999999999
 *      `maximum' will never be returned, ever.
 *
 *  RealRandomFloat( float maximum )
 *      As above only returns a float
 *
 *  RealRandomInt( integer maximum )
 *      Returns an integer value from 0 to (maximum-1).
 *      Again, the value `maximum' will NEVER be returned.
 *
 *  GetM_CNST(void)
 *      This function returns the value of the `M constant' used
 *      by real_random.  This is provided primarily for creative
 *      seed generation and is not needed for normal use.
 *
 *
 *  private functions:
 *
 *      RealRandomCycle( void | int cycleNumber)
 *          Repeatedly calls the RealRandomJW function to cycle the
 *          generator and get it out of repeat cycles caused by 
 *          numerically similar seeds.  The defualt cycle is defined,
 *          but this can be changed with the `cycleNumber' parameter.
 *          If cycleNumber = 0, then the default value is used.
 * 
 ********************************************************************
 *
 *  This code is not intended for distribution.  If you chose to use
 *    it,  please give approprite credit to Park and Miller as
 *    the creators of the algorithm.
 *
 ********************************************************************/


#include "RealRandomJW.h"
    
#define A_CNST 16807.0
#define M_CNST 2147483647.0
#define Q_CNST 127773.0
#define R_CNST 2836.0

#define M_INV 4.6566128752457969E-10
#define Q_INV 7.8263795950631200E-06


/*
|||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||

	B E G I N        S O U R C E        C O D E

|||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
*/

//////////////////////////////////////////////////////////////////////
// Construction/Destruction
//////////////////////////////////////////////////////////////////////

CRealRandomJW::CRealRandomJW()
{
	jwRealRandomSeed = 0.;
}

CRealRandomJW::~CRealRandomJW()
{
	return;
}


/////////////////////////////////////////////////////
//
//  SeedRealRandom ->
//
/////////////////////////////////////////////////////

void CRealRandomJW::SeedRealRandom(double inputSeed)
{
    if (inputSeed < 0.)
		inputSeed += M_CNST;
    
    if (inputSeed > M_CNST)
		inputSeed -= M_CNST;
    
    jwRealRandomSeed = inputSeed;

	RealRandomCycle( 0 );
}


void CRealRandomJW::SeedRealRandom(double inputSeed, int genCycle)
{
    if (inputSeed < 0.)
		inputSeed += M_CNST;
    
    if (inputSeed > M_CNST)
		inputSeed -= M_CNST;
    
    jwRealRandomSeed = inputSeed;

	RealRandomCycle( genCycle );
}


/////////////////////////////////////////////////////
//
//  GetRealRandomSeed ->
//
/////////////////////////////////////////////////////

double CRealRandomJW::GetRealRandomSeed(void)
{
	return jwRealRandomSeed;
}

	
/////////////////////////////////////////////////////
//
//  TimeSeedRealRandom ->
//
/////////////////////////////////////////////////////

double CRealRandomJW::TimeSeedRealRandom(void)
{
    //srand( time( 0 ));			// JFW
    //r = drand48*M_CNST;			// JFW
    double r = time(0);
	SeedRealRandom(r);
	RealRandomCycle(0);
	return(r);
}

double CRealRandomJW::TimeSeedRealRandom(int genCycle)
{
    double r = time( 0 );
	SeedRealRandom( r );
	RealRandomCycle( genCycle );
	return(r);
}


/////////////////////////////////////////////////////
//
//  RealRandom ->
//
/////////////////////////////////////////////////////

double CRealRandomJW::RealRandom(void)
{
    int hi;
    double randomOut;
    double lo,test;
    
    hi = int(jwRealRandomSeed*Q_INV);
    lo = jwRealRandomSeed - Q_CNST*hi;
    test = A_CNST*lo - R_CNST*hi;
    
    if(test>0.0)
	jwRealRandomSeed = test;
    else
	jwRealRandomSeed = test + M_CNST;
    
    randomOut = jwRealRandomSeed*M_INV;
    return(randomOut);
}


/////////////////////////////////////////////////////
//
//  RealRandomDouble ->
//
/////////////////////////////////////////////////////

double CRealRandomJW::RealRandomDouble(double inputMultiplier)
{
    double output;
    
    output = RealRandom();
    if(output == 1.)
	return(0.);
    output *= inputMultiplier;
	return(output);
}


/////////////////////////////////////////////////////
//
//  RealRandomFloat ->
//
/////////////////////////////////////////////////////

float CRealRandomJW::RealRandomFloat(float inputMultiplier)
{
    double output;
    
    output = RealRandom();
    if(output == 1.)
	return(0.);
    output *= (float)inputMultiplier;
	return((float)output);
}


/////////////////////////////////////////////////////
//
//  RealRandomInt ->
//
/////////////////////////////////////////////////////

int CRealRandomJW::RealRandomInt(int inputMultiplier)
{
    double output;
    
    output = RealRandom();
    if(output == 1.)
	return(0);
    output *= (double)inputMultiplier;
    return((int)output);
}

/////////////////////////////////////////////////////
//
//  GetM_CONST ->
//
/////////////////////////////////////////////////////

double CRealRandomJW::GetM_CNST(void)
{
    return(M_CNST);
}


/////////////////////////////////////////////////////
//
//  RealRandomCycle ->
//
/////////////////////////////////////////////////////

void CRealRandomJW::RealRandomCycle(void)
{
	double junk = 0.;
	for(int i=0;i<20;i++)
		junk = junk + RealRandom();
	return;
}

void CRealRandomJW::RealRandomCycle(int cycleNum)
{
    double junk=0.;
    
    if (cycleNum < 1)
		for(int i=0;i<20;i++)
			junk = junk + RealRandom();
	else
		for(int i=0;i<cycleNum;i++)
			junk = junk + RealRandom();
	return;
}