TRandom1.cxx

Go to the documentation of this file.

// @(#)root/mathcore:$Id: TRandom1.cxx 22866 2008-03-27 15:32:50Z rdm $
// Author: Rene Brun from CLHEP & CERNLIB  04/05/2006

//////////////////////////////////////////////////////////////////////////
//
// TRandom1
//
// The Ranlux Random number generator class
//
// The algorithm for this random engine has been taken from the original
// implementation in FORTRAN by Fred James as part of CLHEP.
//
// The initialisation is carried out using a Multiplicative Congruential
// generator using formula constants of L'Ecuyer as described in "F.James,
// Comp. Phys. Comm. 60 (1990) 329-344".
//
//////////////////////////////////////////////////////////////////////////

#include "TRandom1.h"
#include "TRandom3.h"
#include "TMath.h"
#include <stdlib.h>

// Number of instances with automatic seed selection
int TRandom1::fgNumEngines = 0;

// Maximum index into the seed table
int TRandom1::fgMaxIndex = 215;
#ifndef __CINT__
const UInt_t fgSeedTable[215][2] = {
                             {           9876, 54321            },
                             {     1299961164, 253987020        },
                             {      669708517, 2079157264       },
                             {      190904760, 417696270        },
                             {     1289741558, 1376336092       },
                             {     1803730167, 324952955        },
                             {      489854550, 582847132        },
                             {     1348037628, 1661577989       },
                             {      350557787, 1155446919       },
                             {      591502945, 634133404        },
                             {     1901084678, 862916278        },
                             {     1988640932, 1785523494       },
                             {     1873836227, 508007031        },
                             {     1146416592, 967585720        },
                             {     1837193353, 1522927634       },
                             {       38219936, 921609208        },
                             {      349152748, 112892610        },
                             {      744459040, 1735807920       },
                             {     1983990104, 728277902        },
                             {      309164507, 2126677523       },
                             {      362993787, 1897782044       },
                             {      556776976, 462072869        },
                             {     1584900822, 2019394912       },
                             {     1249892722, 791083656        },
                             {     1686600998, 1983731097       },
                             {     1127381380, 198976625        },
                             {     1999420861, 1810452455       },
                             {     1972906041, 664182577        },
                             {       84636481, 1291886301       },
                             {     1186362995, 954388413        },
                             {     2141621785, 61738584         },
                             {     1969581251, 1557880415       },
                             {     1150606439, 136325185        },
                             {       95187861, 1592224108       },
                             {      940517655, 1629971798       },
                             {      215350428, 922659102        },
                             {      786161212, 1121345074       },
                             {     1450830056, 1922787776       },
                             {     1696578057, 2025150487       },
                             {     1803414346, 1851324780       },
                             {     1017898585, 1452594263       },
                             {     1184497978, 82122239         },
                             {      633338765, 1829684974       },
                             {      430889421, 230039326        },
                             {      492544653, 76320266         },
                             {      389386975, 1314148944       },
                             {     1720322786, 709120323        },
                             {     1868768216, 1992898523       },
                             {      443210610, 811117710        },
                             {     1191938868, 1548484733       },
                             {      616890172, 159787986        },
                             {      935835339, 1231440405       },
                             {     1058009367, 1527613300       },
                             {     1463148129, 1970575097       },
                             {     1795336935, 434768675        },
                             {      274019517, 605098487        },
                             {      483689317, 217146977        },
                             {     2070804364, 340596558        },
                             {      930226308, 1602100969       },
                             {      989324440, 801809442        },
                             {      410606853, 1893139948       },
                             {     1583588576, 1219225407       },
                             {     2102034391, 1394921405       },
                             {     2005037790, 2031006861       },
                             {     1244218766, 923231061        },
                             {       49312790, 775496649        },
                             {      721012176, 321339902        },
                             {     1719909107, 1865748178       },
                             {     1156177430, 1257110891       },
                             {      307561322, 1918244397       },
                             {      906041433, 360476981        },
                             {     1591375755, 268492659        },
                             {      461522398, 227343256        },
                             {     2145930725, 2020665454       },
                             {     1938419274, 1331283701       },
                             {      174405412, 524140103        },
                             {      494343653,  18063908        },
                             {     1025534808, 181709577        },
                             {     2048959776, 1913665637       },
                             {      950636517, 794796256        },
                             {     1828843197, 1335757744       },
                             {      211109723, 983900607        },
                             {      825474095, 1046009991       },
                             {      374915657, 381856628        },
                             {     1241296328, 698149463        },
                             {     1260624655, 1024538273       },
                             {      900676210, 1628865823       },
                             {      697951025, 500570753        },
                             {     1007920268, 1708398558       },
                             {      264596520, 624727803        },
                             {     1977924811, 674673241        },
                             {     1440257718, 271184151        },
                             {     1928778847, 993535203        },
                             {     1307807366, 1801502463       },
                             {     1498732610, 300876954        },
                             {     1617712402, 1574250679       },
                             {     1261800762, 1556667280       },
                             {      949929273, 560721070        },
                             {     1766170474, 1953522912       },
                             {     1849939248, 19435166         },
                             {      887262858, 1219627824       },
                             {      483086133, 603728993        },
                             {     1330541052, 1582596025       },
                             {     1850591475, 723593133        },
                             {     1431775678, 1558439000       },
                             {      922493739, 1356554404       },
                             {     1058517206, 948567762        },
                             {      709067283, 1350890215       },
                             {     1044787723, 2144304941       },
                             {      999707003, 513837520        },
                             {     2140038663, 1850568788       },
                             {     1803100150, 127574047        },
                             {      867445693, 1149173981       },
                             {      408583729, 914837991        },
                             {     1166715497, 602315845        },
                             {      430738528, 1743308384       },
                             {     1388022681, 1760110496       },
                             {     1664028066, 654300326        },
                             {     1767741172, 1338181197       },
                             {     1625723550, 1742482745       },
                             {      464486085, 1507852127       },
                             {      754082421, 1187454014       },
                             {     1315342834, 425995190        },
                             {      960416608, 2004255418       },
                             {     1262630671, 671761697        },
                             {       59809238, 103525918        },
                             {     1205644919, 2107823293       },
                             {     1615183160, 1152411412       },
                             {     1024474681, 2118672937       },
                             {     1703877649, 1235091369       },
                             {     1821417852, 1098463802       },
                             {     1738806466, 1529062843       },
                             {      620780646, 1654833544       },
                             {     1070174101, 795158254        },
                             {      658537995, 1693620426       },
                             {     2055317555, 508053916        },
                             {     1647371686, 1282395762       },
                             {       29067379, 409683067        },
                             {     1763495989, 1917939635       },
                             {     1602690753, 810926582        },
                             {      885787576, 513818500        },
                             {     1853512561, 1195205756       },
                             {     1798585498, 1970460256       },
                             {     1819261032, 1306536501       },
                             {     1133245275, 37901            },
                             {      689459799, 1334389069       },
                             {     1730609912, 1854586207       },
                             {     1556832175, 1228729041       },
                             {      251375753, 683687209        },
                             {     2083946182, 1763106152       },
                             {     2142981854, 1365385561       },
                             {      763711891, 1735754548       },
                             {     1581256466, 173689858        },
                             {     2121337132, 1247108250       },
                             {     1004003636, 891894307        },
                             {      569816524, 358675254        },
                             {      626626425, 116062841        },
                             {      632086003, 861268491        },
                             {     1008211580, 779404957        },
                             {     1134217766, 1766838261       },
                             {     1423829292, 1706666192       },
                             {      942037869, 1549358884       },
                             {     1959429535, 480779114        },
                             {      778311037, 1940360875       },
                             {     1531372185, 2009078158       },
                             {      241935492, 1050047003       },
                             {      272453504, 1870883868       },
                             {      390441332, 1057903098       },
                             {     1230238834, 1548117688       },
                             {     1242956379, 1217296445       },
                             {      515648357, 1675011378       },
                             {      364477932, 355212934        },
                             {     2096008713, 1570161804       },
                             {     1409752526, 214033983        },
                             {     1288158292, 1760636178       },
                             {      407562666, 1265144848       },
                             {     1071056491, 1582316946       },
                             {     1014143949, 911406955        },
                             {      203080461, 809380052        },
                             {      125647866, 1705464126       },
                             {     2015685843, 599230667        },
                             {     1425476020, 668203729        },
                             {     1673735652, 567931803        },
                             {     1714199325, 181737617        },
                             {     1389137652, 678147926        },
                             {      288547803, 435433694        },
                             {      200159281, 654399753        },
                             {     1580828223, 1298308945       },
                             {     1832286107, 169991953        },
                             {      182557704, 1046541065       },
                             {     1688025575, 1248944426       },
                             {     1508287706, 1220577001       },
                             {       36721212, 1377275347       },
                             {     1968679856, 1675229747       },
                             {      279109231, 1835333261       },
                             {     1358617667, 1416978076       },
                             {      740626186, 2103913602       },
                             {     1882655908, 251341858        },
                             {      648016670, 1459615287       },
                             {      780255321, 154906988        },
                             {      857296483, 203375965        },
                             {     1631676846, 681204578        },
                             {     1906971307, 1623728832       },
                             {     1541899600, 1168449797       },
                             {     1267051693, 1020078717       },
                             {     1998673940, 1298394942       },
                             {     1914117058, 1381290704       },
                             {      426068513, 1381618498       },
                             {      139365577, 1598767734       },
                             {     2129910384, 952266588        },
                             {      661788054, 19661356         },
                             {     1104640222, 240506063        },
                             {      356133630, 1676634527       },
                             {      242242374, 1863206182       },
                             {      957935844, 1490681416       }};
#endif

ClassImp(TRandom1)

//______________________________________________________________________________
TRandom1::TRandom1(UInt_t seed, Int_t lux)
        : fIntModulus(0x1000000),
          fMantissaBit24( TMath::Power(0.5,24.) ),
          fMantissaBit12( TMath::Power(0.5,12.) )
{
// Luxury level is set in the same way as the original FORTRAN routine.
//  level 0  (p=24): equivalent to the original RCARRY of Marsaglia
//           and Zaman, very long period, but fails many tests.
//  level 1  (p=48): considerable improvement in quality over level 0,
//           now passes the gap test, but still fails spectral test.
//  level 2  (p=97): passes all known tests, but theoretically still
//           defective.
//  level 3  (p=223): DEFAULT VALUE.  Any theoretically possible
//           correlations have very small chance of being observed.
//  level 4  (p=389): highest possible luxury, all 24 bits chaotic.
   UInt_t seedlist[2]={0,0};

   fTheSeeds = &fSeed;
   fLuxury = lux;
   SetSeed2(seed, fLuxury);

   // setSeeds() wants a zero terminated array!
   seedlist[0]=fSeed;
   seedlist[1]=0;
   SetSeeds(seedlist, fLuxury);
}

//______________________________________________________________________________
TRandom1::TRandom1()
        : fIntModulus(0x1000000),
          fMantissaBit24( TMath::Power(0.5,24.) ),
          fMantissaBit12( TMath::Power(0.5,12.) )
{
   //default constructor
   fTheSeeds = &fSeed;
   UInt_t seed;
   UInt_t seedlist[2]={0,0};

   fLuxury = 3;
   int cycle = abs(int(fgNumEngines/fgMaxIndex));
   int curIndex = abs(int(fgNumEngines%fgMaxIndex));
   fgNumEngines +=1;
   UInt_t mask = ((cycle & 0x007fffff) << 8);
   GetTableSeeds( seedlist, curIndex );
   seed = seedlist[0]^mask;
   SetSeed2(seed, fLuxury);

   // setSeeds() wants a zero terminated array!
   seedlist[0]=fSeed; //<=============
   seedlist[1]=0;
   SetSeeds(seedlist, fLuxury);
}

//______________________________________________________________________________
TRandom1::TRandom1(int rowIndex, int colIndex, int lux)
        : fIntModulus(0x1000000),
          fMantissaBit24( TMath::Power(0.5,24.) ),
          fMantissaBit12( TMath::Power(0.5,12.) )
{
   //constructor
   fTheSeeds = &fSeed;
   UInt_t seed;
   UInt_t seedlist[2]={0,0};

   fLuxury = lux;
   int cycle = abs(int(rowIndex/fgMaxIndex));
   int row = abs(int(rowIndex%fgMaxIndex));
   int col = abs(int(colIndex%2));
   UInt_t mask = (( cycle & 0x000007ff ) << 20 );
   GetTableSeeds( seedlist, row );
   seed = ( seedlist[col] )^mask;
   SetSeed2(seed, fLuxury);

   // setSeeds() wants a zero terminated array!
   seedlist[0]=fSeed;
   seedlist[1]=0;
   SetSeeds(seedlist, fLuxury);
}

//______________________________________________________________________________
TRandom1::~TRandom1()
{
   //destructor
}

//______________________________________________________________________________
void TRandom1::GetTableSeeds(UInt_t* seeds, Int_t index)
{
   //static function returning the table of seeds
   if ((index >= 0) && (index < 215)) {
      seeds[0] = fgSeedTable[index][0];
      seeds[1] = fgSeedTable[index][1];
   }
   else seeds = 0;
}

//______________________________________________________________________________
Double_t TRandom1::Rndm(Int_t)
{
   //return a random number in ]0,1]
   float next_random;
   float uni;
   int i;

   uni = fFloatSeedTable[fJlag] - fFloatSeedTable[fIlag] - fCarry;
   if(uni < 0. ) {
      uni += 1.0;
      fCarry = fMantissaBit24;
   } else {
      fCarry = 0.;
   }

   fFloatSeedTable[fIlag] = uni;
   fIlag --;
   fJlag --;
   if(fIlag < 0) fIlag = 23;
   if(fJlag < 0) fJlag = 23;

   if( uni < fMantissaBit12 ){
      uni += fMantissaBit24 * fFloatSeedTable[fJlag];
      if( uni == 0) uni = fMantissaBit24 * fMantissaBit24;
   }
   next_random = uni;
   fCount24 ++;

// every 24th number generation, several random numbers are generated
// and wasted depending upon the fLuxury level.

   if(fCount24 == 24 ) {
      fCount24 = 0;
      for( i = 0; i != fNskip ; i++) {
         uni = fFloatSeedTable[fJlag] - fFloatSeedTable[fIlag] - fCarry;
         if(uni < 0. ) {
            uni += 1.0;
            fCarry = fMantissaBit24;
         } else {
            fCarry = 0.;
         }
         fFloatSeedTable[fIlag] = uni;
         fIlag --;
         fJlag --;
         if(fIlag < 0)fIlag = 23;
         if(fJlag < 0) fJlag = 23;
      }
   }
   return (double) next_random;
}

//______________________________________________________________________________
void TRandom1::RndmArray(const Int_t size, Float_t *vect)
{
   //return an array of random numbers in ]0,1]
   for (Int_t i=0;i<size;i++) vect[i] = Rndm();
}

//______________________________________________________________________________
void TRandom1::RndmArray(const Int_t size, Double_t *vect)
{
   //return an array of random numbers in ]0,1]
   float next_random;
   float uni;
   int i;
   int index;

   for (index=0; index<size; ++index) {
      uni = fFloatSeedTable[fJlag] - fFloatSeedTable[fIlag] - fCarry;
      if(uni < 0. ) {
         uni += 1.0;
         fCarry = fMantissaBit24;
      } else {
         fCarry = 0.;
      }

      fFloatSeedTable[fIlag] = uni;
      fIlag --;
      fJlag --;
      if(fIlag < 0) fIlag = 23;
      if(fJlag < 0) fJlag = 23;

      if( uni < fMantissaBit12 ){
         uni += fMantissaBit24 * fFloatSeedTable[fJlag];
         if( uni == 0) uni = fMantissaBit24 * fMantissaBit24;
      }
      next_random = uni;
      vect[index] = (double)next_random;
      fCount24 ++;

// every 24th number generation, several random numbers are generated
// and wasted depending upon the fLuxury level.

      if(fCount24 == 24 ) {
         fCount24 = 0;
         for( i = 0; i != fNskip ; i++) {
            uni = fFloatSeedTable[fJlag] - fFloatSeedTable[fIlag] - fCarry;
            if(uni < 0. ) {
               uni += 1.0;
               fCarry = fMantissaBit24;
            } else {
               fCarry = 0.;
            }
            fFloatSeedTable[fIlag] = uni;
            fIlag --;
            fJlag --;
            if(fIlag < 0)fIlag = 23;
            if(fJlag < 0) fJlag = 23;
         }
      }
   }
}


//______________________________________________________________________________
void TRandom1::SetSeeds(const UInt_t *seeds, int lux)
{
   //set seeds
   const int ecuyer_a = 53668;
   const int ecuyer_b = 40014;
   const int ecuyer_c = 12211;
   const int ecuyer_d = 2147483563;

   const int lux_levels[5] = {0,24,73,199,365};
   int i;
   UInt_t int_seed_table[24];
   Long64_t k_multiple,next_seed;
   const UInt_t *seedptr;

   fTheSeeds = seeds;
   seedptr   = seeds;

   if(seeds == 0) {
      SetSeed2(fSeed,lux);
      fTheSeeds = &fSeed;
      return;
   }

   fSeed = *seeds;

// number of additional random numbers that need to be 'thrown away'
// every 24 numbers is set using fLuxury level variable.

   if( (lux > 4)||(lux < 0) ) {
      if(lux >= 24) {
         fNskip = lux - 24;
      } else {
         fNskip = lux_levels[3]; // corresponds to default fLuxury level
      }
   } else {
      fLuxury = lux;
      fNskip  = lux_levels[fLuxury];
   }

   for( i = 0;(i != 24)&&(*seedptr != 0);i++) {
      int_seed_table[i] = *seedptr % fIntModulus;
      seedptr++;
   }

   if(i != 24){
      next_seed = int_seed_table[i-1];
      for(;i != 24;i++) {
         k_multiple = next_seed / ecuyer_a;
         next_seed = ecuyer_b * (next_seed - k_multiple * ecuyer_a)
         - k_multiple * ecuyer_c ;
         if(next_seed < 0)next_seed += ecuyer_d;
         int_seed_table[i] = next_seed % fIntModulus;
      }
   }

   for(i = 0;i != 24;i++)
      fFloatSeedTable[i] = int_seed_table[i] * fMantissaBit24;

   fIlag = 23;
   fJlag = 9;
   fCarry = 0. ;

   if( fFloatSeedTable[23] == 0. ) fCarry = fMantissaBit24;

   fCount24 = 0;
}

//______________________________________________________________________________
void TRandom1::SetSeed2(UInt_t seed, int lux)
{
// The initialisation is carried out using a Multiplicative
// Congruential generator using formula constants of L'Ecuyer
// as described in "A review of pseudorandom number generators"
// (Fred James) published in Computer Physics Communications 60 (1990)
// pages 329-344
//
// modified for the case of seed = 0. In that case a random 64 bits seed based on
// TUUID (using TRandom3(0) ) is generated in order to have a unique seed
//

   const int ecuyer_a = 53668;
   const int ecuyer_b = 40014;
   const int ecuyer_c = 12211;
   const int ecuyer_d = 2147483563;

   const int lux_levels[5] = {0,24,73,199,365};

   UInt_t int_seed_table[24];

   // case of seed == 0
   // use a random seed based on TRandom3(0) which is base don the UUID
   if (seed == 0) {
      TRandom3 r3(0);
      seed =  static_cast<UInt_t> (4294967296.*r3.Rndm());
   }


   Long64_t next_seed = seed;
   Long64_t k_multiple;
   int i;

   // number of additional random numbers that need to be 'thrown away'
   // every 24 numbers is set using fLuxury level variable.

   fSeed = seed;
   if( (lux > 4)||(lux < 0) ) {
      if(lux >= 24) {
         fNskip = lux - 24;
      } else {
         fNskip = lux_levels[3]; // corresponds to default fLuxury level
      }
   } else {
      fLuxury = lux;
      fNskip  = lux_levels[fLuxury];
   }


   for(i = 0;i != 24;i++) {
      k_multiple = next_seed / ecuyer_a;
      next_seed = ecuyer_b * (next_seed - k_multiple * ecuyer_a)
      - k_multiple * ecuyer_c ;
      if(next_seed < 0)next_seed += ecuyer_d;
      int_seed_table[i] = next_seed % fIntModulus;
   }

   for(i = 0;i != 24;i++)
      fFloatSeedTable[i] = int_seed_table[i] * fMantissaBit24;

   fIlag = 23;
   fJlag = 9;
   fCarry = 0. ;

   if( fFloatSeedTable[23] == 0. ) fCarry = fMantissaBit24;

   fCount24 = 0;
}

void TRandom1::SetSeed(UInt_t seed)
{
   // Set RanLux seed using default luxury level
   SetSeed2(seed);
}

---