1 // -*- C++ -*-
2 // $Id: DRand48Engine.cc,v 1.7 2010/07/29 16:50:34 garren Exp $
3 // -----------------------------------------------------------------------
4 // HEP Random
5 // --- DRand48Engine ---
6 // class implementation file
7 // -----------------------------------------------------------------------
8 // This file is part of Geant4 (simulation toolkit for HEP).
9
10 // =======================================================================
11 // G.Cosmo - Created: 5th September 1995
12 // - Minor corrections: 31st October 1996
13 // - Added methods for engine status: 19th November 1996
14 // - Added srand48(), seed48(), drand48() implementations
15 // for Windows/NT: 6th March 1997
16 // - Fixed bug in setSeeds(): 15th September 1997
17 // - Private copy constructor and operator=: 26th Feb 1998
18 // J.Marraffino - Added stream operators and related constructor.
19 // Added automatic seed selection from seed table and
20 // engine counter: 16th Feb 1998
21 // J.Marraffino - Remove dependence on hepString class 13 May 1999
22 // E.Tcherniaev - More accurate code for drand48() on NT base on
23 // a code extracted from GNU C Library 2.1.3: 8th Nov 2000
24 // M. Fischler - In restore, checkFile for file not found 03 Dec 2004
25 // M. Fischler - Methods for distrib. instacne save/restore 12/8/04
26 // M. Fischler - split get() into tag validation and
27 // getState() for anonymous restores 12/27/04
28 // M. Fischler - put/get for vectors of ulongs 3/8/05
29 // M. Fischler - State-saving using only ints, for portability 4/12/05
30 //
31 // =======================================================================
32
33 #include "CLHEP/Random/defs.h"
34 #include "CLHEP/Random/Random.h"
35 #include "CLHEP/Random/DRand48Engine.h"
36 #include "CLHEP/Random/RandomFunc.h"
37 #include "CLHEP/Random/engineIDulong.h"
38 #include <string.h> // for strcmp
39 #include <stdlib.h> // for std::abs(int)
40
41 //#define TRACE_IO
42
43 namespace CLHEP {
44
45 static const int MarkerLen = 64; // Enough room to hold a begin or end marker.
46 // Number of instances with automatic seed selection
47 int DRand48Engine::numEngines = 0;
48
name() const49 std::string DRand48Engine::name() const {return "DRand48Engine";}
50
51 // Maximum index into the seed table
52 const int DRand48Engine::maxIndex = 215;
53
DRand48Engine(long seed)54 DRand48Engine::DRand48Engine(long seed)
55 : HepRandomEngine()
56 {
57 setSeed(seed,0);
58 setSeeds(&theSeed,0);
59 }
60
DRand48Engine()61 DRand48Engine::DRand48Engine()
62 : HepRandomEngine()
63 {
64 long seeds[2];
65 long seed;
66
67 int cycle = abs(int(numEngines/maxIndex));
68 int curIndex = abs(int(numEngines%maxIndex));
69 ++numEngines;
70 long mask = ((cycle & 0x007fffff) << 8);
71 HepRandom::getTheTableSeeds( seeds, curIndex );
72 seed = seeds[0]^mask;
73 setSeed(seed,0);
74 setSeeds(&theSeed,0);
75 }
76
DRand48Engine(int rowIndex,int colIndex)77 DRand48Engine::DRand48Engine(int rowIndex, int colIndex)
78 : HepRandomEngine()
79 {
80 long seed;
81 long seeds[2];
82
83 int cycle = abs(int(rowIndex/maxIndex));
84 int row = abs(int(rowIndex%maxIndex));
85 int col = abs(int(colIndex%2));
86 long mask = ((cycle & 0x000007ff) << 20);
87 HepRandom::getTheTableSeeds( seeds, row );
88 seed = (seeds[col])^mask;
89 setSeed(seed,0);
90 setSeeds(&theSeed,0);
91 }
92
DRand48Engine(std::istream & is)93 DRand48Engine::DRand48Engine(std::istream& is)
94 : HepRandomEngine()
95 {
96 is >> *this;
97 }
98
~DRand48Engine()99 DRand48Engine::~DRand48Engine() {}
100
setSeed(long seed,int)101 void DRand48Engine::setSeed(long seed, int)
102 {
103 srand48( seed );
104 theSeed = seed;
105 }
106
setSeeds(const long * seeds,int)107 void DRand48Engine::setSeeds(const long* seeds, int)
108 {
109 setSeed(seeds ? *seeds : 19780503L, 0);
110 theSeeds = seeds;
111 }
112
saveStatus(const char filename[]) const113 void DRand48Engine::saveStatus( const char filename[] ) const
114 {
115 std::ofstream outFile( filename, std::ios::out ) ;
116
117 if (!outFile.bad()) {
118 outFile << "Uvec\n";
119 std::vector<unsigned long> v = put();
120 #ifdef TRACE_IO
121 std::cout << "Result of v = put() is:\n";
122 #endif
123 for (unsigned int i=0; i<v.size(); ++i) {
124 outFile << v[i] << "\n";
125 #ifdef TRACE_IO
126 std::cout << v[i] << " ";
127 if (i%6==0) std::cout << "\n";
128 #endif
129 }
130 #ifdef TRACE_IO
131 std::cout << "\n";
132 #endif
133 }
134
135 #ifdef REMOVED
136 unsigned short dummy[] = { 0, 0, 0 };
137 unsigned short* cseed = seed48(dummy);
138 if (!outFile.bad()) {
139 outFile << theSeed << std::endl;
140 for (int i=0; i<3; ++i) {
141 outFile << cseed[i] << std::endl;
142 dummy[i] = cseed[i];
143 }
144 seed48(dummy);
145 }
146 #endif
147 }
148
restoreStatus(const char filename[])149 void DRand48Engine::restoreStatus( const char filename[] )
150 {
151 std::ifstream inFile( filename, std::ios::in);
152 unsigned short cseed[3];
153
154 if (!checkFile ( inFile, filename, engineName(), "restoreStatus" )) {
155 std::cerr << " -- Engine state remains unchanged\n";
156 return;
157 }
158 if ( possibleKeywordInput ( inFile, "Uvec", theSeed ) ) {
159 std::vector<unsigned long> v;
160 unsigned long xin;
161 for (unsigned int ivec=0; ivec < VECTOR_STATE_SIZE; ++ivec) {
162 inFile >> xin;
163 #ifdef TRACE_IO
164 std::cout << "ivec = " << ivec << " xin = " << xin << " ";
165 if (ivec%3 == 0) std::cout << "\n";
166 #endif
167 if (!inFile) {
168 inFile.clear(std::ios::badbit | inFile.rdstate());
169 std::cerr << "\nDRand48Engine state (vector) description improper."
170 << "\nrestoreStatus has failed."
171 << "\nInput stream is probably mispositioned now." << std::endl;
172 return;
173 }
174 v.push_back(xin);
175 }
176 getState(v);
177 return;
178 }
179
180 if (!inFile.bad() && !inFile.eof()) {
181 inFile >> theSeed;
182 for (int i=0; i<3; ++i)
183 // inFile >> theSeed; removed -- encompased by possibleKeywordInput
184 seed48(cseed);
185 }
186 }
187
showStatus() const188 void DRand48Engine::showStatus() const
189 {
190 unsigned short dummy[] = { 0, 0, 0 };
191 unsigned short* cseed = seed48(dummy);
192 std::cout << std::endl;
193 std::cout << "-------- DRand48 engine status ---------" << std::endl;
194 std::cout << " Initial seed = " << theSeed << std::endl;
195 std::cout << " Current seeds = " << cseed[0] << ", ";
196 std::cout << cseed[1] << ", ";
197 std::cout << cseed[2] << std::endl;
198 std::cout << "----------------------------------------" << std::endl;
199 for (int i=0; i<3; ++i)
200 dummy[i] = cseed[i];
201 seed48(dummy);
202 }
203
flat()204 double DRand48Engine::flat()
205 {
206 double num = 0.;
207
208 while (num == 0.)
209 num = drand48();
210 return num;
211 }
212
flatArray(const int size,double * vect)213 void DRand48Engine::flatArray(const int size, double* vect)
214 {
215 int i;
216
217 for (i=0; i<size; ++i)
218 vect[i]=flat();
219 }
220
put(std::ostream & os) const221 std::ostream & DRand48Engine::put ( std::ostream& os ) const
222 {
223 char beginMarker[] = "DRand48Engine-begin";
224 os << beginMarker << "\nUvec\n";
225 std::vector<unsigned long> v = put();
226 for (unsigned int i=0; i<v.size(); ++i) {
227 os << v[i] << "\n";
228 }
229 return os;
230
231 #ifdef REMOVED
232 unsigned short dummy[] = { 0, 0, 0 };
233 unsigned short* cseed = seed48(dummy);
234 char endMarker[] = "DRand48Engine-end";
235 os << " " << beginMarker << " ";
236 os << theSeed << " ";
237 for (int i=0; i<3; ++i) {
238 dummy[i] = cseed[i];
239 os << cseed[i] << " ";
240 }
241 os << endMarker << " ";
242 seed48(dummy);
243 return os;
244 #endif
245 }
246
put() const247 std::vector<unsigned long> DRand48Engine::put () const {
248 std::vector<unsigned long> v;
249 v.push_back (engineIDulong<DRand48Engine>());
250 unsigned short dummy[] = { 0, 0, 0 };
251 unsigned short* cseed = seed48(dummy);
252 for (int i=0; i<3; ++i) {
253 dummy[i] = cseed[i];
254 v.push_back (static_cast<unsigned long>(cseed[i]));
255 }
256 seed48(dummy);
257 return v;
258 }
259
get(std::istream & is)260 std::istream & DRand48Engine::get ( std::istream& is )
261 {
262 char beginMarker [MarkerLen];
263 is >> std::ws;
264 is.width(MarkerLen); // causes the next read to the char* to be <=
265 // that many bytes, INCLUDING A TERMINATION \0
266 // (Stroustrup, section 21.3.2)
267 is >> beginMarker;
268 if (strcmp(beginMarker,"DRand48Engine-begin")) {
269 is.clear(std::ios::badbit | is.rdstate());
270 std::cerr << "\nInput stream mispositioned or"
271 << "\nDRand48Engine state description missing or"
272 << "\nwrong engine type found." << std::endl;
273 return is;
274 }
275 return getState(is);
276 }
277
beginTag()278 std::string DRand48Engine::beginTag ( ) {
279 return "DRand48Engine-begin";
280 }
281
getState(std::istream & is)282 std::istream & DRand48Engine::getState ( std::istream& is )
283 {
284 unsigned short cseed[3];
285 if ( possibleKeywordInput ( is, "Uvec", cseed[0] ) ) {
286 std::vector<unsigned long> v;
287 unsigned long uu;
288 #ifdef TRACE_IO
289 std::cout << "DRand48Engine::getState detected Uvec keyword\n";
290 uu = 999999;
291 #endif
292 for (unsigned int ivec=0; ivec < VECTOR_STATE_SIZE; ++ivec) {
293 uu = 999999;
294 is >> uu;
295 #ifdef TRACE_IO
296 std::cout << "ivec = " << ivec << " uu = " << uu << "\n";
297 #endif
298 if (!is) {
299 is.clear(std::ios::badbit | is.rdstate());
300 std::cerr << "\nDRand48Engine state (vector) description improper."
301 << "\ngetState() has failed."
302 << "\nInput stream is probably mispositioned now." << std::endl;
303 return is;
304 }
305 v.push_back(uu);
306 }
307 getState(v);
308 return (is);
309 }
310
311 // is >> cseed[0] was removed from loop, encompassed by possibleKeywordInput()
312
313 char endMarker [MarkerLen];
314 is >> theSeed;
315 for (int i=1; i<3; ++i) {
316 is >> cseed[i];
317 }
318 is >> std::ws;
319 is.width(MarkerLen);
320 is >> endMarker;
321 if (strcmp(endMarker,"DRand48Engine-end")) {
322 is.clear(std::ios::badbit | is.rdstate());
323 std::cerr << "\nDRand48Engine state description incomplete."
324 << "\nInput stream is probably mispositioned now." << std::endl;
325 return is;
326 }
327 seed48(cseed);
328 return is;
329 }
330
get(const std::vector<unsigned long> & v)331 bool DRand48Engine::get (const std::vector<unsigned long> & v) {
332 if ((v[0] & 0xffffffffUL) != engineIDulong<DRand48Engine>()) {
333 std::cerr <<
334 "\nDRand48Engine get:state vector has wrong ID word - state unchanged\n";
335 return false;
336 }
337 return getState(v);
338 }
339
getState(const std::vector<unsigned long> & v)340 bool DRand48Engine::getState (const std::vector<unsigned long> & v) {
341 if (v.size() != VECTOR_STATE_SIZE ) {
342 std::cerr <<
343 "\nDRand48Engine getState:state vector has wrong length - state unchanged\n";
344 return false;
345 }
346 unsigned short cseed[3];
347 for (int i=0; i<3; ++i) {
348 cseed[i] = static_cast<unsigned short>(v[i+1]);
349 }
350 seed48(cseed);
351 return true;
352 }
353
354 } // namespace CLHEP
355