Geant4 Cross Reference

Cross-Referencing   Geant4
Geant4/externals/clhep/src/JamesRandom.cc

Version: [ ReleaseNotes ] [ 1.0 ] [ 1.1 ] [ 2.0 ] [ 3.0 ] [ 3.1 ] [ 3.2 ] [ 4.0 ] [ 4.0.p1 ] [ 4.0.p2 ] [ 4.1 ] [ 4.1.p1 ] [ 5.0 ] [ 5.0.p1 ] [ 5.1 ] [ 5.1.p1 ] [ 5.2 ] [ 5.2.p1 ] [ 5.2.p2 ] [ 6.0 ] [ 6.0.p1 ] [ 6.1 ] [ 6.2 ] [ 6.2.p1 ] [ 6.2.p2 ] [ 7.0 ] [ 7.0.p1 ] [ 7.1 ] [ 7.1.p1 ] [ 8.0 ] [ 8.0.p1 ] [ 8.1 ] [ 8.1.p1 ] [ 8.1.p2 ] [ 8.2 ] [ 8.2.p1 ] [ 8.3 ] [ 8.3.p1 ] [ 8.3.p2 ] [ 9.0 ] [ 9.0.p1 ] [ 9.0.p2 ] [ 9.1 ] [ 9.1.p1 ] [ 9.1.p2 ] [ 9.1.p3 ] [ 9.2 ] [ 9.2.p1 ] [ 9.2.p2 ] [ 9.2.p3 ] [ 9.2.p4 ] [ 9.3 ] [ 9.3.p1 ] [ 9.3.p2 ] [ 9.4 ] [ 9.4.p1 ] [ 9.4.p2 ] [ 9.4.p3 ] [ 9.4.p4 ] [ 9.5 ] [ 9.5.p1 ] [ 9.5.p2 ] [ 9.6 ] [ 9.6.p1 ] [ 9.6.p2 ] [ 9.6.p3 ] [ 9.6.p4 ] [ 10.0 ] [ 10.0.p1 ] [ 10.0.p2 ] [ 10.0.p3 ] [ 10.0.p4 ] [ 10.1 ] [ 10.1.p1 ] [ 10.1.p2 ] [ 10.1.p3 ] [ 10.2 ] [ 10.2.p1 ] [ 10.2.p2 ] [ 10.2.p3 ] [ 10.3 ] [ 10.3.p1 ] [ 10.3.p2 ] [ 10.3.p3 ] [ 10.4 ] [ 10.4.p1 ] [ 10.4.p2 ] [ 10.4.p3 ] [ 10.5 ] [ 10.5.p1 ] [ 10.6 ] [ 10.6.p1 ] [ 10.6.p2 ] [ 10.6.p3 ] [ 10.7 ] [ 10.7.p1 ] [ 10.7.p2 ] [ 10.7.p3 ] [ 10.7.p4 ] [ 11.0 ] [ 11.0.p1 ] [ 11.0.p2 ] [ 11.0.p3, ] [ 11.0.p4 ] [ 11.1 ] [ 11.1.1 ] [ 11.1.2 ] [ 11.1.3 ] [ 11.2 ] [ 11.2.1 ] [ 11.2.2 ] [ 11.3.0 ]

Diff markup

Differences between /externals/clhep/src/JamesRandom.cc (Version 11.3.0) and /externals/clhep/src/JamesRandom.cc (Version 9.5.p2)


                                                   >>   1 // $Id:$
  1 // -*- C++ -*-                                      2 // -*- C++ -*-
  2 //                                                  3 //
  3 // -------------------------------------------      4 // -----------------------------------------------------------------------
  4 //                             HEP Random           5 //                             HEP Random
  5 //                       --- HepJamesRandom --      6 //                       --- HepJamesRandom ---
  6 //                      class implementation f      7 //                      class implementation file
  7 // -------------------------------------------      8 // -----------------------------------------------------------------------
  8 // This file is part of Geant4 (simulation too      9 // This file is part of Geant4 (simulation toolkit for HEP).
  9 //                                                 10 //
 10 // This algorithm implements the original univ     11 // This algorithm implements the original universal random number generator
 11 // as proposed by Marsaglia & Zaman in report      12 // as proposed by Marsaglia & Zaman in report FSU-SCRI-87-50 and coded
 12 // in FORTRAN77 by Fred James as the RANMAR ge     13 // in FORTRAN77 by Fred James as the RANMAR generator, part of the MATHLIB
 13 // HEP library.                                    14 // HEP library.
 14                                                    15 
 15 // ===========================================     16 // =======================================================================
 16 // Gabriele Cosmo - Created: 5th September 199     17 // Gabriele Cosmo - Created: 5th September 1995
 17 //                - Fixed a bug in setSeed():      18 //                - Fixed a bug in setSeed(): 26th February 1996
 18 //                - Minor corrections: 31st Oc     19 //                - Minor corrections: 31st October 1996
 19 //                - Added methods for engine s     20 //                - Added methods for engine status: 19th November 1996
 20 //                - Fixed bug in setSeeds(): 1     21 //                - Fixed bug in setSeeds(): 15th September 1997
 21 // J.Marraffino   - Added stream operators and     22 // J.Marraffino   - Added stream operators and related constructor.
 22 //                  Added automatic seed selec     23 //                  Added automatic seed selection from seed table and
 23 //                  engine counter: 16th Feb 1     24 //                  engine counter: 16th Feb 1998
 24 // Ken Smith      - Added conversion operators     25 // Ken Smith      - Added conversion operators:  6th Aug 1998
 25 // J. Marraffino  - Remove dependence on hepSt     26 // J. Marraffino  - Remove dependence on hepString class  13 May 1999
 26 // V. Innocente   - changed pointers to indice     27 // V. Innocente   - changed pointers to indices     3 may 2000
 27 // M. Fischler    - In restore, checkFile for      28 // M. Fischler    - In restore, checkFile for file not found    03 Dec 2004
 28 // M. Fischler    - Methods for distrib. insta     29 // M. Fischler    - Methods for distrib. instacne save/restore  12/8/04    
 29 // M. Fischler    - split get() into tag valid     30 // M. Fischler    - split get() into tag validation and 
 30 //                  getState() for anonymous r     31 //                  getState() for anonymous restores           12/27/04    
 31 // M. Fischler    - Enforcement that seeds be      32 // M. Fischler    - Enforcement that seeds be non-negative
 32 //        (lest the sequence be non-random)        33 //        (lest the sequence be non-random)          2/14/05    
 33 // M. Fischler    - put/get for vectors of ulo     34 // M. Fischler    - put/get for vectors of ulongs   3/14/05
 34 // M. Fischler    - State-saving using only in     35 // M. Fischler    - State-saving using only ints, for portability 4/12/05
 35 //                                                 36 //        
 36 // ===========================================     37 // =======================================================================
 37                                                    38 
 38 #include "CLHEP/Random/Random.h"                   39 #include "CLHEP/Random/Random.h"
 39 #include "CLHEP/Random/JamesRandom.h"              40 #include "CLHEP/Random/JamesRandom.h"
 40 #include "CLHEP/Random/engineIDulong.h"            41 #include "CLHEP/Random/engineIDulong.h"
 41 #include "CLHEP/Random/DoubConv.h"                 42 #include "CLHEP/Random/DoubConv.h"
 42 #include "CLHEP/Utility/atomic_int.h"          <<  43 #include <string.h> // for strcmp
 43                                                << 
 44 #include <atomic>                              << 
 45 #include <cmath>                                   44 #include <cmath>
 46 #include <cstdlib>                                 45 #include <cstdlib>
 47 #include <iostream>                            << 
 48 #include <string.h> // for strcmp              << 
 49 #include <string>                              << 
 50 #include <vector>                              << 
 51                                                    46 
 52 namespace CLHEP {                                  47 namespace CLHEP {
 53                                                    48 
 54 namespace {                                    << 
 55   // Number of instances with automatic seed s << 
 56   CLHEP_ATOMIC_INT_TYPE numberOfEngines(0);    << 
 57                                                << 
 58   // Maximum index into the seed table         << 
 59   const int maxIndex = 215;                    << 
 60 }                                              << 
 61                                                << 
 62 static const int MarkerLen = 64; // Enough roo     49 static const int MarkerLen = 64; // Enough room to hold a begin or end marker. 
 63                                                    50 
 64 std::string HepJamesRandom::name() const {retu     51 std::string HepJamesRandom::name() const {return "HepJamesRandom";}
 65                                                    52 
                                                   >>  53 // Number of instances with automatic seed selection
                                                   >>  54 int HepJamesRandom::numEngines = 0;
                                                   >>  55 
                                                   >>  56 // Maximum index into the seed table
                                                   >>  57 int HepJamesRandom::maxIndex = 215;
                                                   >>  58 
 66 HepJamesRandom::HepJamesRandom(long seed)          59 HepJamesRandom::HepJamesRandom(long seed)
 67 : HepRandomEngine()                                60 : HepRandomEngine()
 68 {                                                  61 {
 69   setSeed(seed,0);                                 62   setSeed(seed,0);
 70   setSeeds(&theSeed,0);                            63   setSeeds(&theSeed,0);
 71 }                                                  64 }
 72                                                    65 
 73 HepJamesRandom::HepJamesRandom()      // 15 Fe     66 HepJamesRandom::HepJamesRandom()      // 15 Feb. 1998  JMM
 74 : HepRandomEngine()                                67 : HepRandomEngine()
 75 {                                                  68 {
 76   long seeds[2];                                   69   long seeds[2];
 77   long seed;                                       70   long seed;
 78                                                    71 
 79   int numEngines = numberOfEngines++;          << 
 80   int cycle = std::abs(int(numEngines/maxIndex     72   int cycle = std::abs(int(numEngines/maxIndex));
 81   int curIndex = std::abs(int(numEngines%maxIn     73   int curIndex = std::abs(int(numEngines%maxIndex));
 82                                                <<  74   ++numEngines;
 83   long mask = ((cycle & 0x007fffff) << 8);         75   long mask = ((cycle & 0x007fffff) << 8);
 84   HepRandom::getTheTableSeeds( seeds, curIndex     76   HepRandom::getTheTableSeeds( seeds, curIndex );
 85   seed = seeds[0]^mask;                            77   seed = seeds[0]^mask;
 86   setSeed(seed,0);                                 78   setSeed(seed,0);
 87   setSeeds(&theSeed,0);                            79   setSeeds(&theSeed,0);
 88 }                                                  80 }
 89                                                    81 
 90 HepJamesRandom::HepJamesRandom(int rowIndex, i     82 HepJamesRandom::HepJamesRandom(int rowIndex, int colIndex) // 15 Feb. 1998  JMM
 91 : HepRandomEngine()                                83 : HepRandomEngine()
 92 {                                                  84 {
 93   long seed;                                       85   long seed;
 94    long seeds[2];                                  86    long seeds[2];
 95                                                    87 
 96   int cycle = std::abs(int(rowIndex/maxIndex))     88   int cycle = std::abs(int(rowIndex/maxIndex));
 97   int row = std::abs(int(rowIndex%maxIndex));      89   int row = std::abs(int(rowIndex%maxIndex));
 98   int col = std::abs(int(colIndex%2));             90   int col = std::abs(int(colIndex%2));
 99   long mask = ((cycle & 0x000007ff) << 20);        91   long mask = ((cycle & 0x000007ff) << 20);
100   HepRandom::getTheTableSeeds( seeds, row );       92   HepRandom::getTheTableSeeds( seeds, row );
101   seed = (seeds[col])^mask;                        93   seed = (seeds[col])^mask;
102   setSeed(seed,0);                                 94   setSeed(seed,0);
103   setSeeds(&theSeed,0);                            95   setSeeds(&theSeed,0);
104 }                                                  96 }
105                                                    97 
106 HepJamesRandom::HepJamesRandom(std::istream& i     98 HepJamesRandom::HepJamesRandom(std::istream& is)
107 : HepRandomEngine()                                99 : HepRandomEngine()
108 {                                                 100 {
109   is >> *this;                                    101   is >> *this;
110 }                                                 102 }
111                                                   103 
112 HepJamesRandom::~HepJamesRandom() {}              104 HepJamesRandom::~HepJamesRandom() {}
113                                                   105 
114 void HepJamesRandom::saveStatus( const char fi    106 void HepJamesRandom::saveStatus( const char filename[] ) const
115 {                                                 107 {
116   std::ofstream outFile( filename, std::ios::o    108   std::ofstream outFile( filename, std::ios::out ) ;
117                                                   109 
118   if (!outFile.bad()) {                           110   if (!outFile.bad()) {
119     outFile << "Uvec\n";                          111     outFile << "Uvec\n";
120     std::vector<unsigned long> v = put();         112     std::vector<unsigned long> v = put();
121     for (unsigned int i=0; i<v.size(); ++i) {     113     for (unsigned int i=0; i<v.size(); ++i) {
122       outFile << v[i] << "\n";                    114       outFile << v[i] << "\n";
123     }                                             115     }
124   }                                               116   }
125 }                                                 117 }
126                                                   118 
127 void HepJamesRandom::restoreStatus( const char    119 void HepJamesRandom::restoreStatus( const char filename[] )
128 {                                                 120 {
129    int ipos, jpos;                                121    int ipos, jpos;
130    std::ifstream inFile( filename, std::ios::i    122    std::ifstream inFile( filename, std::ios::in);
131    if (!checkFile ( inFile, filename, engineNa    123    if (!checkFile ( inFile, filename, engineName(), "restoreStatus" )) {
132      std::cerr << "  -- Engine state remains u    124      std::cerr << "  -- Engine state remains unchanged\n";
133      return;                                      125      return;
134    }                                              126    }
135   if ( possibleKeywordInput ( inFile, "Uvec",     127   if ( possibleKeywordInput ( inFile, "Uvec", theSeed ) ) {
136     std::vector<unsigned long> v;                 128     std::vector<unsigned long> v;
137     unsigned long xin;                            129     unsigned long xin;
138     for (unsigned int ivec=0; ivec < VECTOR_ST    130     for (unsigned int ivec=0; ivec < VECTOR_STATE_SIZE; ++ivec) {
139       inFile >> xin;                              131       inFile >> xin;
140       if (!inFile) {                              132       if (!inFile) {
141         inFile.clear(std::ios::badbit | inFile    133         inFile.clear(std::ios::badbit | inFile.rdstate());
142         std::cerr << "\nJamesRandom state (vec    134         std::cerr << "\nJamesRandom state (vector) description improper."
143          << "\nrestoreStatus has failed."         135          << "\nrestoreStatus has failed."
144          << "\nInput stream is probably mispos    136          << "\nInput stream is probably mispositioned now." << std::endl;
145         return;                                   137         return;
146       }                                           138       }
147       v.push_back(xin);                           139       v.push_back(xin);
148     }                                             140     }
149     getState(v);                                  141     getState(v);
150     return;                                       142     return;
151   }                                               143   }
152                                                   144 
153    if (!inFile.bad() && !inFile.eof()) {          145    if (!inFile.bad() && !inFile.eof()) {
154 //     inFile >> theSeed;  removed -- encompas    146 //     inFile >> theSeed;  removed -- encompased by possibleKeywordInput
155      for (int i=0; i<97; ++i)                     147      for (int i=0; i<97; ++i)
156        inFile >> u[i];                            148        inFile >> u[i];
157      inFile >> c; inFile >> cd; inFile >> cm;     149      inFile >> c; inFile >> cd; inFile >> cm;
158      inFile >> jpos;                              150      inFile >> jpos;
159      ipos = (64+jpos)%97;                         151      ipos = (64+jpos)%97;
160      i97 = ipos;                                  152      i97 = ipos;
161      j97 = jpos;                                  153      j97 = jpos;
162    }                                              154    }
163 }                                                 155 }
164                                                   156 
165 void HepJamesRandom::showStatus() const           157 void HepJamesRandom::showStatus() const
166 {                                                 158 {
167    std::cout << std::endl;                        159    std::cout << std::endl;
168    std::cout << "----- HepJamesRandom engine s    160    std::cout << "----- HepJamesRandom engine status -----" << std::endl;
169    std::cout << " Initial seed = " << theSeed     161    std::cout << " Initial seed = " << theSeed << std::endl;
170    std::cout << " u[] = ";                        162    std::cout << " u[] = ";
171    for (int i=0; i<97; ++i)                       163    for (int i=0; i<97; ++i)
172      std::cout << u[i] << " ";                    164      std::cout << u[i] << " ";
173    std::cout << std::endl;                        165    std::cout << std::endl;
174    std::cout << " c = " << c << ", cd = " << c    166    std::cout << " c = " << c << ", cd = " << cd << ", cm = " << cm
175        << std::endl;                              167        << std::endl;
176    std::cout << " i97 = " << i97 << ", u[i97]     168    std::cout << " i97 = " << i97 << ", u[i97] = " << u[i97] << std::endl;
177    std::cout << " j97 = " << j97 << ", u[j97]     169    std::cout << " j97 = " << j97 << ", u[j97] = " << u[j97] << std::endl;
178    std::cout << "-----------------------------    170    std::cout << "----------------------------------------" << std::endl;
179 }                                                 171 }
180                                                   172 
181 void HepJamesRandom::setSeed(long seed, int)      173 void HepJamesRandom::setSeed(long seed, int)
182 {                                                 174 {
183   // The input value for "seed" should be with    175   // The input value for "seed" should be within the range [0,900000000]
184   //                                              176   //
185   // Negative seeds result in serious flaws in    177   // Negative seeds result in serious flaws in the randomness;
186   // seeds above 900000000 are OK because of t    178   // seeds above 900000000 are OK because of the %177 in the expression for i,
187   // but may have the same effect as other see    179   // but may have the same effect as other seeds below 900000000.
188                                                   180 
189   int m, n;                                       181   int m, n;
190   float s, t;                                     182   float s, t;
191   long mm;                                        183   long mm;
192                                                   184 
193   if (seed < 0) {                                 185   if (seed < 0) {
194     std::cout << "Seed for HepJamesRandom must    186     std::cout << "Seed for HepJamesRandom must be non-negative\n" 
195       << "Seed value supplied was " << seed       187       << "Seed value supplied was " << seed  
196   << "\nUsing its absolute value instead\n";      188   << "\nUsing its absolute value instead\n";
197     seed = -seed;                                 189     seed = -seed;
198   }                                               190   }
199                                                   191   
200   long ij = seed/30082;                           192   long ij = seed/30082;
201   long kl = seed - 30082*ij;                      193   long kl = seed - 30082*ij;
202   long i = (ij/177) % 177 + 2;                    194   long i = (ij/177) % 177 + 2;
203   long j = ij % 177 + 2;                          195   long j = ij % 177 + 2;
204   long k = (kl/169) % 178 + 1;                    196   long k = (kl/169) % 178 + 1;
205   long l = kl % 169;                              197   long l = kl % 169;
206                                                   198 
207   theSeed = seed;                                 199   theSeed = seed;
208                                                   200 
209   for ( n = 1 ; n < 98 ; n++ ) {                  201   for ( n = 1 ; n < 98 ; n++ ) {
210     s = 0.0;                                      202     s = 0.0;
211     t = 0.5;                                      203     t = 0.5;
212     for ( m = 1 ; m < 25 ; m++) {                 204     for ( m = 1 ; m < 25 ; m++) {
213       mm = ( ( (i*j) % 179 ) * k ) % 179;         205       mm = ( ( (i*j) % 179 ) * k ) % 179;
214       i = j;                                      206       i = j;
215       j = k;                                      207       j = k;
216       k = mm;                                     208       k = mm;
217       l = ( 53 * l + 1 ) % 169;                   209       l = ( 53 * l + 1 ) % 169;
218       if ( (l*mm % 64 ) >= 32 )                   210       if ( (l*mm % 64 ) >= 32 )
219         s += t;                                   211         s += t;
220       t *= 0.5;                                   212       t *= 0.5;
221     }                                             213     }
222     u[n-1] = s;                                   214     u[n-1] = s;
223   }                                               215   }
224   c = 362436.0 / 16777216.0;                      216   c = 362436.0 / 16777216.0;
225   cd = 7654321.0 / 16777216.0;                    217   cd = 7654321.0 / 16777216.0;
226   cm = 16777213.0 / 16777216.0;                   218   cm = 16777213.0 / 16777216.0;
227                                                   219 
228   i97 = 96;                                       220   i97 = 96;
229   j97 = 32;                                       221   j97 = 32;
230                                                   222 
231 }                                                 223 }
232                                                   224 
233 void HepJamesRandom::setSeeds(const long* seed    225 void HepJamesRandom::setSeeds(const long* seeds, int)
234 {                                                 226 {
235   setSeed(seeds ? *seeds : 19780503L, 0);         227   setSeed(seeds ? *seeds : 19780503L, 0);
236   theSeeds = seeds;                               228   theSeeds = seeds;
237 }                                                 229 }
238                                                   230 
239 double HepJamesRandom::flat()                     231 double HepJamesRandom::flat()
240 {                                                 232 {
241    double uni;                                    233    double uni;
242                                                   234 
243    do {                                           235    do {
244       uni = u[i97] - u[j97];                      236       uni = u[i97] - u[j97];
245       if ( uni < 0.0 ) uni++;                     237       if ( uni < 0.0 ) uni++;
246       u[i97] = uni;                               238       u[i97] = uni;
247                                                   239       
248       if (i97 == 0) i97 = 96;                     240       if (i97 == 0) i97 = 96;
249       else i97--;                                 241       else i97--;
250                                                   242       
251       if (j97 == 0) j97 = 96;                     243       if (j97 == 0) j97 = 96;
252       else j97--;                                 244       else j97--;
253                                                   245       
254       c -= cd;                                    246       c -= cd;
255       if (c < 0.0) c += cm;                       247       if (c < 0.0) c += cm;
256                                                   248       
257       uni -= c;                                   249       uni -= c;
258       if (uni < 0.0) uni += 1.0;                  250       if (uni < 0.0) uni += 1.0;
259    } while ( uni <= 0.0 || uni >= 1.0 );          251    } while ( uni <= 0.0 || uni >= 1.0 );
260                                                   252    
261    return uni;                                    253    return uni;
262 }                                                 254 }
263                                                   255 
264 void HepJamesRandom::flatArray(const int size,    256 void HepJamesRandom::flatArray(const int size, double* vect)
265 {                                                 257 {
266 //   double uni;                                  258 //   double uni;
267    int i;                                         259    int i;
268                                                   260 
269    for (i=0; i<size; ++i) {                       261    for (i=0; i<size; ++i) {
270      vect[i] = flat();                            262      vect[i] = flat();
271    }                                              263    }   
272 }                                                 264 }
273                                                   265 
274 HepJamesRandom::operator double() {            << 
275    return flat();                              << 
276 }                                              << 
277                                                << 
278 HepJamesRandom::operator float() {             << 
279    return float( flat() );                     << 
280 }                                              << 
281                                                << 
282 HepJamesRandom::operator unsigned int() {         266 HepJamesRandom::operator unsigned int() {
283    return ((unsigned int)(flat() * exponent_bi    267    return ((unsigned int)(flat() * exponent_bit_32()) & 0xffffffff )  |
284          (((unsigned int)( u[i97] * exponent_b    268          (((unsigned int)( u[i97] * exponent_bit_32())>>16)  & 0xff);
285 }                                                 269 }
286                                                   270 
287 std::ostream & HepJamesRandom::put ( std::ostr    271 std::ostream & HepJamesRandom::put ( std::ostream& os ) const {
288   char beginMarker[] = "JamesRandom-begin";       272   char beginMarker[] = "JamesRandom-begin";
289   os << beginMarker << "\nUvec\n";                273   os << beginMarker << "\nUvec\n";
290   std::vector<unsigned long> v = put();           274   std::vector<unsigned long> v = put();
291   for (unsigned int i=0; i<v.size(); ++i) {       275   for (unsigned int i=0; i<v.size(); ++i) {
292      os <<  v[i] <<  "\n";                        276      os <<  v[i] <<  "\n";
293   }                                               277   }
294   return os;                                      278   return os;  
295 }                                                 279 }
296                                                   280 
297 std::vector<unsigned long> HepJamesRandom::put    281 std::vector<unsigned long> HepJamesRandom::put () const {
298   std::vector<unsigned long> v;                   282   std::vector<unsigned long> v;
299   v.push_back (engineIDulong<HepJamesRandom>()    283   v.push_back (engineIDulong<HepJamesRandom>());
300   std::vector<unsigned long> t;                   284   std::vector<unsigned long> t;
301   for (int i=0; i<97; ++i) {                      285   for (int i=0; i<97; ++i) {
302     t = DoubConv::dto2longs(u[i]);                286     t = DoubConv::dto2longs(u[i]);
303     v.push_back(t[0]); v.push_back(t[1]);         287     v.push_back(t[0]); v.push_back(t[1]);
304   }                                               288   }
305   t = DoubConv::dto2longs(c);                     289   t = DoubConv::dto2longs(c);
306   v.push_back(t[0]); v.push_back(t[1]);           290   v.push_back(t[0]); v.push_back(t[1]);
307   t = DoubConv::dto2longs(cd);                    291   t = DoubConv::dto2longs(cd);
308   v.push_back(t[0]); v.push_back(t[1]);           292   v.push_back(t[0]); v.push_back(t[1]);
309   t = DoubConv::dto2longs(cm);                    293   t = DoubConv::dto2longs(cm);
310   v.push_back(t[0]); v.push_back(t[1]);           294   v.push_back(t[0]); v.push_back(t[1]);
311   v.push_back(static_cast<unsigned long>(j97))    295   v.push_back(static_cast<unsigned long>(j97));
312   return v;                                       296   return v;
313 }                                                 297 }
314                                                   298 
315                                                   299 
316 std::istream & HepJamesRandom::get  ( std::ist    300 std::istream & HepJamesRandom::get  ( std::istream& is) {
317   char beginMarker [MarkerLen];                   301   char beginMarker [MarkerLen];
318   is >> std::ws;                                  302   is >> std::ws;
319   is.width(MarkerLen);  // causes the next rea    303   is.width(MarkerLen);  // causes the next read to the char* to be <=
320       // that many bytes, INCLUDING A TERMINAT    304       // that many bytes, INCLUDING A TERMINATION \0 
321       // (Stroustrup, section 21.3.2)             305       // (Stroustrup, section 21.3.2)
322   is >> beginMarker;                              306   is >> beginMarker;
323   if (strcmp(beginMarker,"JamesRandom-begin"))    307   if (strcmp(beginMarker,"JamesRandom-begin")) {
324      is.clear(std::ios::badbit | is.rdstate())    308      is.clear(std::ios::badbit | is.rdstate());
325      std::cerr << "\nInput stream mispositione    309      std::cerr << "\nInput stream mispositioned or"
326          << "\nJamesRandom state description m    310          << "\nJamesRandom state description missing or"
327          << "\nwrong engine type found." << st    311          << "\nwrong engine type found." << std::endl;
328      return is;                                   312      return is;
329   }                                               313   }
330   return getState(is);                            314   return getState(is);
331 }                                                 315 }
332                                                   316 
333 std::string HepJamesRandom::beginTag ( )  {       317 std::string HepJamesRandom::beginTag ( )  { 
334   return "JamesRandom-begin";                     318   return "JamesRandom-begin"; 
335 }                                                 319 }
336                                                   320 
337 std::istream & HepJamesRandom::getState  ( std    321 std::istream & HepJamesRandom::getState  ( std::istream& is) {
338   if ( possibleKeywordInput ( is, "Uvec", theS    322   if ( possibleKeywordInput ( is, "Uvec", theSeed ) ) {
339     std::vector<unsigned long> v;                 323     std::vector<unsigned long> v;
340     unsigned long uu;                             324     unsigned long uu;
341     for (unsigned int ivec=0; ivec < VECTOR_ST    325     for (unsigned int ivec=0; ivec < VECTOR_STATE_SIZE; ++ivec) {
342       is >> uu;                                   326       is >> uu;
343       if (!is) {                                  327       if (!is) {
344         is.clear(std::ios::badbit | is.rdstate    328         is.clear(std::ios::badbit | is.rdstate());
345         std::cerr << "\nJamesRandom state (vec    329         std::cerr << "\nJamesRandom state (vector) description improper."
346     << "\ngetState() has failed."                 330     << "\ngetState() has failed."
347          << "\nInput stream is probably mispos    331          << "\nInput stream is probably mispositioned now." << std::endl;
348         return is;                                332         return is;
349       }                                           333       }
350       v.push_back(uu);                            334       v.push_back(uu);
351     }                                             335     }
352     getState(v);                                  336     getState(v);
353     return (is);                                  337     return (is);
354   }                                               338   }
355                                                   339 
356 //  is >> theSeed;  Removed, encompassed by po    340 //  is >> theSeed;  Removed, encompassed by possibleKeywordInput()
357                                                   341 
358   int ipos, jpos;                                 342   int ipos, jpos;
359   char   endMarker [MarkerLen];                   343   char   endMarker [MarkerLen];
360   for (int i=0; i<97; ++i) {                      344   for (int i=0; i<97; ++i) {
361      is >> u[i];                                  345      is >> u[i];
362   }                                               346   }
363   is >> c; is >> cd; is >> cm;                    347   is >> c; is >> cd; is >> cm;
364   is >> jpos;                                     348   is >> jpos;
365   is >> std::ws;                                  349   is >> std::ws;
366   is.width(MarkerLen);                            350   is.width(MarkerLen);
367   is >> endMarker;                                351   is >> endMarker;
368   if(strcmp(endMarker,"JamesRandom-end")) {       352   if(strcmp(endMarker,"JamesRandom-end")) {
369      is.clear(std::ios::badbit | is.rdstate())    353      is.clear(std::ios::badbit | is.rdstate());
370      std::cerr << "\nJamesRandom state descrip    354      std::cerr << "\nJamesRandom state description incomplete."
371          << "\nInput stream is probably mispos    355          << "\nInput stream is probably mispositioned now." << std::endl;
372      return is;                                   356      return is;
373   }                                               357   }
374                                                   358 
375   ipos = (64+jpos)%97;                            359   ipos = (64+jpos)%97;
376   i97 = ipos;                                     360   i97 = ipos;
377   j97 = jpos;                                     361   j97 = jpos;
378   return is;                                      362   return is;
379 }                                                 363 }
380                                                   364 
381 bool HepJamesRandom::get (const std::vector<un    365 bool HepJamesRandom::get (const std::vector<unsigned long> & v) {
382   if ( (v[0] & 0xffffffffUL) != engineIDulong<    366   if ( (v[0] & 0xffffffffUL) != engineIDulong<HepJamesRandom>()) {
383     std::cerr <<                                  367     std::cerr << 
384       "\nHepJamesRandom get:state vector has w    368       "\nHepJamesRandom get:state vector has wrong ID word - state unchanged\n";
385     return false;                                 369     return false;
386   }                                               370   }
387   return getState(v);                             371   return getState(v);
388 }                                                 372 }
389                                                   373 
390 bool HepJamesRandom::getState (const std::vect    374 bool HepJamesRandom::getState (const std::vector<unsigned long> & v) {
391   if (v.size() != VECTOR_STATE_SIZE ) {           375   if (v.size() != VECTOR_STATE_SIZE ) {
392     std::cerr <<                                  376     std::cerr << 
393       "\nHepJamesRandom get:state vector has w    377       "\nHepJamesRandom get:state vector has wrong length - state unchanged\n";
394     return false;                                 378     return false;
395   }                                               379   }
396   std::vector<unsigned long> t(2);                380   std::vector<unsigned long> t(2);
397   for (int i=0; i<97; ++i) {                      381   for (int i=0; i<97; ++i) {
398     t[0] = v[2*i+1]; t[1] = v[2*i+2];             382     t[0] = v[2*i+1]; t[1] = v[2*i+2];
399     u[i] = DoubConv::longs2double(t);             383     u[i] = DoubConv::longs2double(t);
400   }                                               384   }
401   t[0] = v[195]; t[1] = v[196]; c  = DoubConv:    385   t[0] = v[195]; t[1] = v[196]; c  = DoubConv::longs2double(t);
402   t[0] = v[197]; t[1] = v[198]; cd = DoubConv:    386   t[0] = v[197]; t[1] = v[198]; cd = DoubConv::longs2double(t);
403   t[0] = v[199]; t[1] = v[200]; cm = DoubConv:    387   t[0] = v[199]; t[1] = v[200]; cm = DoubConv::longs2double(t);
404   j97  = (int)v[201];                          << 388   j97  = v[201];
405   i97  = (64+j97)%97;                             389   i97  = (64+j97)%97; 
406   return true;                                    390   return true;
407 }                                                 391 }
408                                                   392 
409 }  // namespace CLHEP                             393 }  // namespace CLHEP
410                                                   394