/* A C-program for MT19937: Real number version([0,1)-interval) */
/* (1999/10/28)                                                 */
/*   genrand() generates one pseudorandom real number (double)  */
/* which is uniformly distributed on [0,1)-interval, for each   */
/* call. sgenrand(seed) sets initial values to the working area */
/* of 624 words. Before genrand(), sgenrand(seed) must be       */
/* called once. (seed is any 32-bit integer.)                   */
/* Integer generator is obtained by modifying two lines.        */
/*   Coded by Takuji Nishimura, considering the suggestions by  */
/* Topher Cooper and Marc Rieffel in July-Aug. 1997.            */

/* This library is free software; you can redistribute it and/or   */
/* modify it under the terms of the GNU Library General Public     */
/* License as published by the Free Software Foundation; either    */
/* version 2 of the License, or (at your option) any later         */
/* version.                                                        */
/* This library 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 Library General Public License for more details.    */
/* You should have received a copy of the GNU Library General      */
/* Public License along with this library; if not, write to the    */
/* Free Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA   */
/* 02111-1307  USA                                                 */

/* Copyright (C) 1997, 1999 Makoto Matsumoto and Takuji Nishimura. */
/* Any feedback is very welcome. For any question, comments,       */
/* see http://www.math.keio.ac.jp/matumoto/emt.html or email       */
/* matumoto@math.keio.ac.jp                                        */

/* REFERENCE                                                       */
/* M. Matsumoto and T. Nishimura,                                  */
/* "Mersenne Twister: A 623-Dimensionally Equidistributed Uniform  */
/* Pseudo-Random Number Generator",                                */
/* ACM Transactions on Modeling and Computer Simulation,           */
/* Vol. 8, No. 1, January 1998, pp 3--30.                          */


/* Adapted by Martin Reinsprecht and Raphael Langerhorst */
/* for the G System */


#include "PseudoRNG.h"

namespace GBE
{

namespace Util
{

// Standard Constructor
PseudoRNG::PseudoRNG()
{}

// Standard Destructor
PseudoRNG::~PseudoRNG()
{}


double PseudoRNG::getNumberDouble()
{
  return this->genrandDouble();
}

unsigned long PseudoRNG::getNumberInt()
{
  return this->genrandInt();
}

void PseudoRNG::setNumber( unsigned long j )
{
  this->sgenrand( j );
}

/* Initializing the array with a seed */
void PseudoRNG::sgenrand( unsigned long seed )
{
  int i;

  for ( i = 0;i < N;i++ )
  {
    mt[ i ] = seed & 0xffff0000;
    seed = 69069 * seed + 1;
    mt[ i ] |= ( seed & 0xffff0000 ) >> 16;
    seed = 69069 * seed + 1;
  }
  mti = N;
}

/* Initialization by "sgenrand()" is an example. Theoretically,      */
/* there are 2^19937-1 possible states as an intial state.           */
/* This function allows to choose any of 2^19937-1 ones.             */
/* Essential bits in "seed_array[]" is following 19937 bits:         */
/*  (seed_array[0]&UPPER_MASK), seed_array[1], ..., seed_array[N-1]. */
/* (seed_array[0]&LOWER_MASK) is discarded.                          */
/* Theoretically,                                                    */
/*  (seed_array[0]&UPPER_MASK), seed_array[1], ..., seed_array[N-1]  */
/* can take any values except all zeros.                             */
void PseudoRNG::lsgenrand( unsigned long seed_array[] )
{
  /* the length of seed_array[] must be at least N */
  int i;

  for ( i = 0;i < N;i++ )
    mt[ i ] = seed_array[ i ];
  mti = N;
}

/* generating reals */
// alte Form
//double PseudoRNG::genrand(unsigned long y,unsigned long mag01[2])
double PseudoRNG::genrandDouble()
{
  return ( ( double ) genrandInt() * 2.3283064365386963e-10 ); /* reals: [0,1)-interval */
  /* return y; */ /* for integer generation */
}

/* unsigned long */
// alte Form
//double PseudoRNG::genrand(unsigned long y,unsigned long mag01[2])
unsigned long PseudoRNG::genrandInt()
{
  unsigned long mag01[ 2 ];
  unsigned long y;
  /* mag01[x] = x * MATRIX_A  for x=0,1 */
  mag01[ 0 ] = 0x0;
  mag01[ 1 ] = MATRIX_A;

  if ( mti >= N )
  { /* generate N words at one time */
    int kk;

    if ( mti == N + 1 )    /* if sgenrand() has not been called, */
      sgenrand( 4357 ); /* a default initial seed is used   */

    for ( kk = 0;kk < N - M;kk++ )
    {
      y = ( mt[ kk ] & UPPER_MASK ) | ( mt[ kk + 1 ] & LOWER_MASK );
      mt[ kk ] = mt[ kk + M ] ^ ( y >> 1 ) ^ mag01[ y & 0x1 ];
    }
    for ( ;kk < N - 1;kk++ )
    {
      y = ( mt[ kk ] & UPPER_MASK ) | ( mt[ kk + 1 ] & LOWER_MASK );
      mt[ kk ] = mt[ kk + ( M - N ) ] ^ ( y >> 1 ) ^ mag01[ y & 0x1 ];
    }
    y = ( mt[ N - 1 ] & UPPER_MASK ) | ( mt[ 0 ] & LOWER_MASK );
    mt[ N - 1 ] = mt[ M - 1 ] ^ ( y >> 1 ) ^ mag01[ y & 0x1 ];

    mti = 0;
  }

  y = mt[ mti++ ];
  y ^= TEMPERING_SHIFT_U( y );
  y ^= TEMPERING_SHIFT_S( y ) & TEMPERING_MASK_B;
  y ^= TEMPERING_SHIFT_T( y ) & TEMPERING_MASK_C;
  y ^= TEMPERING_SHIFT_L( y );

  return y;
}

/* This main() outputs first 1000 generated numbers.  */

// we need no main here :-)
/*
main()
{ 
    int i;
 
    sgenrand(4357);
    for (i=0; i<1000; i++) {
        printf("%10.8f ", genrand());
        if (i%5==4) printf("\n");
    }
}*/

}

}