CLHEP VERSION Reference Documentation
   
CLHEP Home Page     CLHEP Documentation     CLHEP Bug Reports

Matrix/CLHEP/Random/RandGaussZiggurat.h
Go to the documentation of this file.
1 /*
2 Code adapted from:
3 http://www.jstatsoft.org/v05/i08/
4 
5 
6 Original disclaimer:
7 
8 The ziggurat method for RNOR and REXP
9 Combine the code below with the main program in which you want
10 normal or exponential variates. Then use of RNOR in any expression
11 will provide a standard normal variate with mean zero, variance 1,
12 while use of REXP in any expression will provide an exponential variate
13 with density exp(-x),x>0.
14 Before using RNOR or REXP in your main, insert a command such as
15 zigset(86947731 );
16 with your own choice of seed value>0, rather than 86947731.
17 (If you do not invoke zigset(...) you will get all zeros for RNOR and REXP.)
18 For details of the method, see Marsaglia and Tsang, "The ziggurat method
19 for generating random variables", Journ. Statistical Software.
20 
21 */
22 
23 
24 #ifndef RandGaussZiggurat_h
25 #define RandGaussZiggurat_h 1
26 
27 #include "cmath"
28 #include "CLHEP/Random/defs.h"
29 #include "CLHEP/Random/RandGauss.h"
30 
31 namespace CLHEP {
32 
37 class RandGaussZiggurat : public RandGauss {
38 
39 public:
40 
41  inline RandGaussZiggurat ( HepRandomEngine& anEngine, double mean=0.0, double stdDev=1.0 );
42  inline RandGaussZiggurat ( HepRandomEngine* anEngine, double mean=0.0, double stdDev=1.0 );
43 
44  // Destructor
45  virtual ~RandGaussZiggurat();
46 
47  // Static methods to shoot random values using the static generator
48 
49  static inline float shoot() {return ziggurat_RNOR(HepRandom::getTheEngine());};
50  static inline float shoot( float mean, float stdDev ) {return shoot()*stdDev + mean;};
51 
52  static void shootArray ( const int size, float* vect, float mean=0.0, float stdDev=1.0 );
53  static void shootArray ( const int size, double* vect, double mean=0.0, double stdDev=1.0 );
54 
55  // Static methods to shoot random values using a given engine
56  // by-passing the static generator.
57 
58  static inline float shoot( HepRandomEngine* anotherEngine ) {return ziggurat_RNOR(anotherEngine);};
59  static inline float shoot( HepRandomEngine* anotherEngine, float mean, float stdDev ) {return shoot(anotherEngine)*stdDev + mean;};
60 
61  static void shootArray ( HepRandomEngine* anotherEngine, const int size, float* vect, float mean=0.0, float stdDev=1.0 );
62  static void shootArray ( HepRandomEngine* anotherEngine, const int size, double* vect, double mean=0.0, double stdDev=1.0 );
63 
64  // Instance methods using the localEngine to instead of the static
65  // generator, and the default mean and stdDev established at construction
66 
67  inline float fire() {return ziggurat_RNOR(localEngine.get()) * defaultStdDev + defaultMean;};
68 
69  inline float fire ( float mean, float stdDev ) {return ziggurat_RNOR(localEngine.get()) * stdDev + mean;};
70 
71  void fireArray ( const int size, float* vect);
72  void fireArray ( const int size, double* vect);
73  void fireArray ( const int size, float* vect, float mean, float stdDev );
74  void fireArray ( const int size, double* vect, double mean, double stdDev );
75 
76  virtual double operator()();
77  virtual double operator()( double mean, double stdDev );
78 
79  // Save and restore to/from streams
80 
81  std::ostream & put ( std::ostream & os ) const;
82  std::istream & get ( std::istream & is );
83 
84  std::string name() const;
86 
87  static std::string distributionName() {return "RandGaussZiggurat";}
88  // Provides the name of this distribution class
89 
90  static bool ziggurat_init();
91 protected:
93  // Ziggurat Original code:
95  //static unsigned long jz,jsr=123456789;
96  //
97  //#define SHR3 (jz=jsr, jsr^=(jsr<<13), jsr^=(jsr>>17), jsr^=(jsr<<5),jz+jsr)
98  //#define UNI (.5 + (signed) SHR3*.2328306e-9)
99  //#define IUNI SHR3
100  //
101  //static long hz;
102  //static unsigned long iz, kn[128], ke[256];
103  //static float wn[128],fn[128], we[256],fe[256];
104  //
105  //#define RNOR (hz=SHR3, iz=hz&127, (fabs(hz)<kn[iz])? hz*wn[iz] : nfix())
106  //#define REXP (jz=SHR3, iz=jz&255, ( jz <ke[iz])? jz*we[iz] : efix())
107 
108  static unsigned long kn[128], ke[256];
109  static float wn[128],fn[128], we[256],fe[256];
110 
111  static bool ziggurat_is_init;
112 
113  static inline unsigned long ziggurat_SHR3(HepRandomEngine* anEngine) {return (unsigned int)(*anEngine);};
114  static inline float ziggurat_UNI(HepRandomEngine* anEngine) {return anEngine->flat();};
115  static inline float ziggurat_RNOR(HepRandomEngine* anEngine) {
116  long hz=(signed)ziggurat_SHR3(anEngine);
117  unsigned long iz=hz&127;
118  return ((unsigned long)abs(hz)<kn[iz]) ? hz*wn[iz] : ziggurat_nfix(hz,anEngine);
119  };
120  static float ziggurat_nfix(long hz,HepRandomEngine* anEngine);
121 
122 private:
123 
124  // Private copy constructor. Defining it here disallows use.
126 
127  // All the engine info, and the default mean and sigma, are in the RandGauss
128  // base class.
129 
130 };
131 
132 } // namespace CLHEP
133 
134 #ifdef ENABLE_BACKWARDS_COMPATIBILITY
135 // backwards compatibility will be enabled ONLY in CLHEP 1.9
136 using namespace CLHEP;
137 #endif
138 
139 namespace CLHEP {
140 
141 RandGaussZiggurat::RandGaussZiggurat(HepRandomEngine & anEngine, double mean,double stdDev ): RandGauss(anEngine, mean, stdDev)
142 {
144 }
145 
146 RandGaussZiggurat::RandGaussZiggurat(HepRandomEngine * anEngine, double mean,double stdDev ): RandGauss(anEngine, mean, stdDev)
147 {
149 }
150 
151 } // namespace CLHEP
152 
153 #endif
CLHEP::RandGaussZiggurat::distributionName
static std::string distributionName()
Definition: Matrix/CLHEP/Random/RandGaussZiggurat.h:87
CLHEP::RandGaussZiggurat::ziggurat_init
static bool ziggurat_init()
Definition: RandGaussZiggurat.cc:24
CLHEP::RandGaussZiggurat::get
std::istream & get(std::istream &is)
Definition: RandGaussZiggurat.cc:174
CLHEP::HepRandomEngine
Definition: Matrix/CLHEP/Random/RandomEngine.h:55
CLHEP::RandGaussZiggurat::fn
static float fn[128]
Definition: Matrix/CLHEP/Random/RandGaussZiggurat.h:109
CLHEP::RandGaussZiggurat::shoot
static float shoot(HepRandomEngine *anotherEngine, float mean, float stdDev)
Definition: Matrix/CLHEP/Random/RandGaussZiggurat.h:59
is
HepRotation and so forth isNear() norm2() rectify() static Rotation row1 row4(To avoid bloat in the code pulled in for programs which don 't use all these features, we split the implementation .cc files. Only isNear() goes into the original Rotation.cc) --------------------------------------- HepAxisAngle and HepEulerAngles classes --------------------------------------- These classes are very useful and simple structures for holding the result of a nice intuituve decomposition of a rotation there is no longer much content in the distinct ZOOM PhysicsVectors library The only content left in the library is the object files representing the various Exception objects When we build the CLHEP classes for the ZOOM we will set up so as to use ZOOM SpaceVector is(but we can disable namespace usage and most of our users do so at this point). What I do is leave Hep3Vector in the global namespace
CLHEP::RandGaussZiggurat::fe
static float fe[256]
Definition: Matrix/CLHEP/Random/RandGaussZiggurat.h:109
CLHEP::RandGaussZiggurat::fire
float fire()
Definition: Matrix/CLHEP/Random/RandGaussZiggurat.h:67
CLHEP::RandGaussZiggurat::ziggurat_RNOR
static float ziggurat_RNOR(HepRandomEngine *anEngine)
Definition: Matrix/CLHEP/Random/RandGaussZiggurat.h:115
size
user code seldom needs to call this function directly ZMerrno whether or not they are still recorded ZMerrno size() Return the(integer) number of ZMthrow 'n exceptions currently recorded. 5) ZMerrno.clear() Set an internal counter to zero. This counter is available(see next function) to user code to track ZMthrow 'n exceptions that have occurred during any arbitrary time interval. 6) ZMerrno.countSinceCleared() Return the(integer) number of ZMthrow 'n exceptions that have been recorded via ZMerrno.write()
CLHEP::RandGaussZiggurat::shootArray
static void shootArray(const int size, float *vect, float mean=0.0, float stdDev=1.0)
Definition: RandGaussZiggurat.cc:110
CLHEP::RandGauss::defaultStdDev
double defaultStdDev
Definition: Matrix/CLHEP/Random/RandGauss.h:153
CLHEP::RandGaussZiggurat::operator()
virtual double operator()()
Definition: RandGaussZiggurat.cc:102
CLHEP
Definition: ClhepVersion.h:13
CLHEP::RandGaussZiggurat
Definition: Matrix/CLHEP/Random/RandGaussZiggurat.h:37
CLHEP::RandGauss::defaultMean
double defaultMean
Definition: Matrix/CLHEP/Random/RandGauss.h:152
CLHEP::RandGaussZiggurat::ziggurat_nfix
static float ziggurat_nfix(long hz, HepRandomEngine *anEngine)
Definition: RandGaussZiggurat.cc:75
CLHEP::RandGaussZiggurat::shoot
static float shoot(float mean, float stdDev)
Definition: Matrix/CLHEP/Random/RandGaussZiggurat.h:50
CLHEP::RandGaussZiggurat::engine
HepRandomEngine & engine()
Definition: RandGaussZiggurat.cc:14
CLHEP::RandGaussZiggurat::shoot
static float shoot(HepRandomEngine *anotherEngine)
Definition: Matrix/CLHEP/Random/RandGaussZiggurat.h:58
CLHEP::RandGaussZiggurat::fireArray
void fireArray(const int size, float *vect)
Definition: RandGaussZiggurat.cc:138
CLHEP::RandGaussZiggurat::wn
static float wn[128]
Definition: Matrix/CLHEP/Random/RandGaussZiggurat.h:109
CLHEP::HepRandomEngine::flat
virtual double flat()=0
CLHEP::RandGaussZiggurat::RandGaussZiggurat
RandGaussZiggurat(HepRandomEngine &anEngine, double mean=0.0, double stdDev=1.0)
Definition: Matrix/CLHEP/Random/RandGaussZiggurat.h:141
CLHEP::RandGaussZiggurat::ziggurat_is_init
static bool ziggurat_is_init
Definition: Matrix/CLHEP/Random/RandGaussZiggurat.h:111
CLHEP::RandGaussZiggurat::name
std::string name() const
Definition: RandGaussZiggurat.cc:19
CLHEP::RandGaussZiggurat::put
std::ostream & put(std::ostream &os) const
Definition: RandGaussZiggurat.cc:166
CLHEP::RandGaussZiggurat::kn
static unsigned long kn[128]
Definition: Matrix/CLHEP/Random/RandGaussZiggurat.h:108
CLHEP::RandGauss
Definition: Matrix/CLHEP/Random/RandGauss.h:42
CLHEP::RandGaussZiggurat::~RandGaussZiggurat
virtual ~RandGaussZiggurat()
Definition: RandGaussZiggurat.cc:16
CLHEP::RandGaussZiggurat::we
static float we[256]
Definition: Matrix/CLHEP/Random/RandGaussZiggurat.h:109
CLHEP::RandGaussZiggurat::fire
float fire(float mean, float stdDev)
Definition: Matrix/CLHEP/Random/RandGaussZiggurat.h:69
CLHEP::RandGauss::localEngine
shared_ptr< HepRandomEngine > localEngine
Definition: Matrix/CLHEP/Random/RandGauss.h:155
CLHEP::RandGaussZiggurat::ziggurat_UNI
static float ziggurat_UNI(HepRandomEngine *anEngine)
Definition: Matrix/CLHEP/Random/RandGaussZiggurat.h:114
CLHEP::RandGaussZiggurat::ke
static unsigned long ke[256]
Definition: Matrix/CLHEP/Random/RandGaussZiggurat.h:108
CLHEP::HepRandom::getTheEngine
static HepRandomEngine * getTheEngine()
Definition: Random.cc:166
CLHEP::RandGaussZiggurat::shoot
static float shoot()
Definition: Matrix/CLHEP/Random/RandGaussZiggurat.h:49
CLHEP::RandGaussZiggurat::ziggurat_SHR3
static unsigned long ziggurat_SHR3(HepRandomEngine *anEngine)
Definition: Matrix/CLHEP/Random/RandGaussZiggurat.h:113