/* 
 *  @OPENGROUP_COPYRIGHT@
 *  COPYRIGHT NOTICE
 *  Copyright (c) 1990, 1991, 1992, 1993 Open Software Foundation, Inc.
 *  Copyright (c) 1996, 1997, 1998, 1999, 2000 The Open Group
 *  ALL RIGHTS RESERVED (MOTIF). See the file named COPYRIGHT.MOTIF for
 *  the full copyright text.
 *  
 *  This software is subject to an open license. It may only be
 *  used on, with or for operating systems which are themselves open
 *  source systems. You must contact The Open Group for a license
 *  allowing distribution and sublicensing of this software on, with,
 *  or for operating systems which are not Open Source programs.
 *  
 *  See http://www.opengroup.org/openmotif/license for full
 *  details of the license agreement. Any use, reproduction, or
 *  distribution of the program constitutes recipient's acceptance of
 *  this agreement.
 *  
 *  EXCEPT AS EXPRESSLY SET FORTH IN THIS AGREEMENT, THE PROGRAM IS
 *  PROVIDED ON AN "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY
 *  KIND, EITHER EXPRESS OR IMPLIED INCLUDING, WITHOUT LIMITATION, ANY
 *  WARRANTIES OR CONDITIONS OF TITLE, NON-INFRINGEMENT, MERCHANTABILITY
 *  OR FITNESS FOR A PARTICULAR PURPOSE
 *  
 *  EXCEPT AS EXPRESSLY SET FORTH IN THIS AGREEMENT, NEITHER RECIPIENT
 *  NOR ANY CONTRIBUTORS SHALL HAVE ANY LIABILITY FOR ANY DIRECT,
 *  INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
 *  DAMAGES (INCLUDING WITHOUT LIMITATION LOST PROFITS), HOWEVER CAUSED
 *  AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
 *  LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN
 *  ANY WAY OUT OF THE USE OR DISTRIBUTION OF THE PROGRAM OR THE
 *  EXERCISE OF ANY RIGHTS GRANTED HEREUNDER, EVEN IF ADVISED OF THE
 *  POSSIBILITY OF SUCH DAMAGES.
*/ 
/* 
 * HISTORY
*/ 
#ifdef REV_INFO
#ifndef lint
static char rcsid[] = "$XConsortium: Random.c /main/7 1995/07/14 11:36:33 drk $"
#endif
#endif
/***********************************************************************
  	@(#)Random.c	1.3.1.1	Date:1/22/91
  	Author: TAT
	History:
            05/24/90 SJS add to sccs
	Calls:

	Summary:
            This file contains routines that manipulate RANDOM NUMBERS
 
************************************************************************/
#include "xislib.h"
#include <math.h>

#define MAX_ERROR 1e-8
#define INV_MAX_ERROR 1e8
#define ROUND(x) ( floor((double)(x)*INV_MAX_ERROR + 0.5) * MAX_ERROR )


#define RANDOM_TABLE_SIZE 97

/*************************************************************
*           RANDOM INTEGER FROM 0.0 to 1.0                   *
**************************************************************
*
*       Returns a uniform random deviate between 0.0 and 1.0. 
*       Set seed to any negative value to reinitialize sequence.
*
**************************************************************/

static double xisPortableRandom(seed)           
long seed; 
{
    static first_pass = 1;
    static long m1=259200, ia1=7141, ic1=54773;
    static long m2=134456, ia2=8121, ic2=28411;
    static long m3=243000, ia3=4561, ic3=51349;
    static long glix1,glix2,glix3;
    static double glr[RANDOM_TABLE_SIZE];
    static double rm1,rm2;
    long j;
    double return_value;

    /** initialize table and glix's on first call **/

    if ((seed < 0)||first_pass) {  
        first_pass = 0;
        rm1 = (double)1.0/m1;
        rm2 = (double)1.0/m2;
        glix1 = (ic1 - seed) % m1;   /* Seed the first routine */                                               
        glix1 = (ia1*glix1 + ic1) % m1;
        glix2 = glix1 % m2;        /* and use it to seed the second */
        glix1 = (ia1*glix1 + ic1) % m1;
        glix3 = glix1 % m3;          /* and third routines. */

        /* Fill the table with sequential uniform deviates generated by the */
        /* first two routines.                                              */
        for (j=0; j < RANDOM_TABLE_SIZE; j++) {   
            glix1 = (ia1*glix1 + ic1) % m1;
            glix2 = (ia2*glix2 + ic2) % m2;
            glr[j] = (glix1 + glix2*rm2)*rm1; /* Combine low,high order pieces*/
        }
    }

    /*** Except when initializing, this is where we start. ****/

    /* Generate the next number for each sequence. */

    glix1 = (ia1*glix1 + ic1) % m1;  
    glix2 = (ia2*glix2 + ic2) % m2;  
    glix3 = (ia3*glix3 + ic3) % m3;  

    /* Use the third sequence to get an integer */

    j = (RANDOM_TABLE_SIZE*glix3) / m3;   

    if ((j >= RANDOM_TABLE_SIZE) || (j < 0))
        (*xisWarningMsg)("in xisPortableRandom() j=%d is out of range\n",j);
                                         
    return_value = glr[j];                   /* Return that table entry and */
    glr[j] = ROUND((glix1 + glix2*rm2)*rm1); /* fill it with a new table entry*/

    return (return_value);
}

/*************************************************************
*           INITIALIZE RANDOM
**************************************************************                         
*
*  Initializes the psuedo random number generator.
*
**************************************************************/
void xisInitRandom(num)
int num;
{
    if (num > 0)
        num = (-num);
    else if (num == 0)
        num = -5678;
    xisPortableRandom((long)(num));
}

/*************************************************************
*           SELECT AN ARBITRARY INTEGER
**************************************************************                         
*
*  Selects an arbitrary integer between lo and hi.
*
**************************************************************/
int xisArbitrary(lo,hi)
int lo,hi;
{
    long i;
    i = ROUND(((hi+0.9999)-lo)*xisPortableRandom(1) + lo);
    return(i);
}