/* 
 * Motif
 *
 * Copyright (c) 1987-2012, The Open Group. All rights reserved.
 *
 * These libraries and programs are free software; you can
 * redistribute them and/or modify them under the terms of the GNU
 * Lesser General Public License as published by the Free Software
 * Foundation; either version 2 of the License, or (at your option)
 * any later version.
 *
 * These libraries and programs are distributed in the hope that
 * they will be useful, but WITHOUT ANY WARRANTY; without even the
 * implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
 * PURPOSE. See the GNU Lesser General Public License for more
 * details.
 *
 * You should have received a copy of the GNU Lesser General Public
 * License along with these librararies and programs; if not, write
 * to the Free Software Foundation, Inc., 51 Franklin Street, Fifth
 * Floor, Boston, MA 02110-1301 USA
*/ 
/* 
 * 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);
}