Line data Source code
1 : #include <iostream>
2 : #include "PhotosRandom.h"
3 : #include "Photos.h"
4 : #include "Log.h"
5 :
6 : namespace Photospp
7 : {
8 :
9 : bool PhotosRandom::init = false;
10 : int PhotosRandom::iseed[2]= { 1802, 9373 };
11 : int PhotosRandom::i97 = 96;
12 : int PhotosRandom::j97 = 32;
13 : double PhotosRandom::uran[97]= { 0.0 };
14 : double PhotosRandom::cran = 362436.0 /16777216.0;
15 : const double PhotosRandom::cdran = 7654321.0 /16777216.0;
16 : const double PhotosRandom::cmran = 16777213.0/16777216.0;
17 :
18 : void PhotosRandom::setSeed(int s1,int s2)
19 : {
20 0 : if(s1<0 || s1>31327) Log::Fatal("PhotosRandom::setSeed(): Seed(1) out of range [0,31327]",8);
21 0 : if(s2<0 || s2>30080) Log::Fatal("PhotosRandom::setSeed(): Seed(2) out of range [0,30080]",9);
22 0 : iseed[0]=s1;
23 0 : iseed[1]=s2;
24 0 : }
25 :
26 : /*******************************************************************************
27 : PHORIN: PHOton radiation in decays RANdom number generator init
28 :
29 : Purpose: Initialse PHORAN with the user specified seeds in the
30 : array iseed. For details see also: F. James CERN DD-
31 : Report November 1988.
32 :
33 : Author(s): B. van Eijk and F. James Created at: 27/09/89
34 : Last Update: 22/02/90
35 : Rewritten to C++: 18/10/10
36 : by T. Przedzinski (tprzedzi@cern.ch)
37 : *******************************************************************************/
38 : void PhotosRandom::initialize()
39 : {
40 : long IS1,IS2,IS3,IS4,IS5;
41 : double S,T;
42 :
43 : // Calculate Marsaglia and Zaman seeds (by F. James)
44 0 : IS1=(iseed[0]/177)%177+2;
45 0 : IS2= iseed[0]%177+2;
46 0 : IS3=(iseed[1]/169)%178+1;
47 0 : IS4= iseed[1]%169;
48 0 : for(int i=0;i<97;i++)
49 : {
50 : S=0.0;
51 : T=0.5;
52 0 : for(int j=0;j<24;j++)
53 : {
54 0 : IS5=( ((IS1*IS2)%179)*IS3 )%179;
55 : IS1=IS2;
56 : IS2=IS3;
57 : IS3=IS5;
58 0 : IS4=(53*IS4+1)%169;
59 0 : if( (IS4*IS5)%64>=32) S=S+T;
60 0 : T=0.5*T;
61 : }
62 0 : uran[i]=S;
63 : }
64 0 : init=true;
65 0 : Log::Debug(0)<<"PhotosRandom::inititalize(): seed: "<<iseed[0]<<", "<<iseed[1]<<std::endl;
66 0 : }
67 :
68 : /*******************************************************************************
69 : PHORAN: PHOton radiation in decays ret number generator based
70 : on Marsaglia Algorithm
71 :
72 : Purpose: Generate uniformly distributed random numbers between
73 : 0 and 1. Super long period: 2**144. See also:
74 : G. Marsaglia and A. Zaman, FSU-SCR-87-50, for seed mo-
75 : difications to this version see: F. James DD-Report,
76 : November 1988. The generator has to be initialized by
77 : a call to PHORIN ( C++ version: initialize() ).
78 :
79 : Author(s): B. van Eijk, G. Marsaglia and Created at: 27/09/89
80 : A. Zaman Last Update: 27/09/89
81 : Rewritten to C++: 18/10/10
82 : by T. Przedzinski (tprzedzi@cern.ch)
83 : *******************************************************************************/
84 : double PhotosRandom::randomReal()
85 : {
86 0 : if(!init) Log::Fatal("PhotosRandom::randomReal(): generator not initialized",1);
87 : double ret=0.0;
88 0 : while(true)
89 : {
90 0 : ret = uran[i97]-uran[j97];
91 0 : if(ret<0.0) ret+=1.;
92 0 : uran[i97]=ret;
93 0 : i97--;
94 0 : if(i97<0) i97=96;
95 0 : j97--;
96 0 : if(j97<0) j97=96;
97 0 : cran-=cdran;
98 0 : if(cran<0.0) cran+=cmran;
99 0 : ret-=cran;
100 0 : if(ret<0.0) ret+=1.0;
101 0 : if(ret>0.0) break;
102 : }
103 0 : return ret;
104 0 : }
105 :
106 : } // namespace Photospp
107 :
|
