/* @(#)ran3.c	17.1.1.1 (ES0-DMD) 01/25/02 17:34:55 */
/*===========================================================================
  Copyright (C) 1995 European Southern Observatory (ESO)
 
  This program is free software; you can redistribute it and/or 
  modify it under the terms of the GNU General Public License as 
  published by the Free Software Foundation; either version 2 of 
  the License, or (at your option) any later version.
 
  This program is distributed in the hope that it will be useful,
  but WITHOUT ANY WARRANTY; without even the implied warranty of
  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  GNU General Public License for more details.
 
  You should have received a copy of the GNU General Public 
  License along with this program; if not, write to the Free 
  Software Foundation, Inc., 675 Massachusetss Ave, Cambridge, 
  MA 02139, USA.
 
  Corresponding concerning ESO-MIDAS should be addressed as follows:
	Internet e-mail: midas@eso.org
	Postal address: European Southern Observatory
			Data Management Division 
			Karl-Schwarzschild-Strasse 2
			D 85748 Garching bei Muenchen 
			GERMANY
===========================================================================*/

/*++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
.COPYRIGHT  (c)  1990  European Southern Observatory
.IDENT      ran3.c
.LANGUAGE   C
.AUTHOR     P.Grosbol,  IPG/ESO
.COMMENT    Algorithm taken from 'Numerical Recipes in C' p.209
.KEYWORDS   uniformed distributions, random numbers
.VERSION    1.0  1990-Oct-28 : Creation, PJG
----------------------------------------------------------------------*/
#define MBIG    1000000000         /* Special constants for generator */
#define MSEED   161803398
#define MZ      0
#define FAC     (1.0/MBIG)

static  int     inext,inextp;
static  long    ma[55];
static  int     iff=0;


float randm(pdum)
/*++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
.PURPOSE   generate random numbers from a uniformed distribution in
           the range [0.0;1.0]. Initiatize sequence if *pdum<0.
.RETURN    random number from uniformed distribution [0.0;1.0].
----------------------------------------------------------------------*/
int *pdum;
{
  long mj, mk;
  int  i, ii, k;

  if (*pdum<0 || iff==0) {
     iff = 1;
     mj = MSEED - ((*pdum<0) ? -*pdum : *pdum);
     mj %= MBIG;
     ma[54] = mj;
     mk = 1;
     for (i=0; i<54; i++) {
        ii=(21*i)%55 - 1;
        ma[ii] = mk;
        mk = mj-mk;
        if (mk < MZ) mk += MBIG;
        mj = ma[ii];
      }
     for (k=0; k<4; k++)
        for (i=0; i<55; i++) {
           ma[i] -= ma[(i+30) % 55];
           if (ma[i] < MZ) ma[i] += MBIG;
         }
     inext = 0;
     inextp = 31;
     *pdum = 1;
   }
  if (++inext==55) inext=1;
  if (++inextp==55) inextp=1;
  mj = ma[inext] - ma[inextp];
  if (mj<MZ) mj += MBIG;
  ma[inext] = mj;

  return mj*FAC;
}

#undef MBIG
#undef MSEED
#undef MZ
#undef FAC