Geant4 Cross Reference |
1 // -*- C++ -*- 1 // -*- C++ -*- 2 // 2 // 3 // ------------------------------------------- 3 // ----------------------------------------------------------------------- 4 // HEP Random 4 // HEP Random 5 // --- RandPoissonQ -- 5 // --- RandPoissonQ --- 6 // class header file 6 // class header file 7 // ------------------------------------------- 7 // ----------------------------------------------------------------------- 8 8 9 // Class defining RandPoissonQ, which is deriv 9 // Class defining RandPoissonQ, which is derived from RandPoison. 10 // The user interface is identical; but RandGa 10 // The user interface is identical; but RandGaussQ is much faster in all cases 11 // and a bit less accurate when mu > 100. 11 // and a bit less accurate when mu > 100. 12 12 13 // =========================================== 13 // ======================================================================= 14 // M. Fischler - Created: 4th Feb 2000 14 // M. Fischler - Created: 4th Feb 2000 15 // M Fischler - put and get to/from strea 15 // M Fischler - put and get to/from streams 12/10/04 16 // 16 // 17 // =========================================== 17 // ======================================================================= 18 18 19 #ifndef RandPoissonQ_h 19 #ifndef RandPoissonQ_h 20 #define RandPoissonQ_h 1 20 #define RandPoissonQ_h 1 21 21 22 #include "CLHEP/Random/Random.h" 22 #include "CLHEP/Random/Random.h" 23 #include "CLHEP/Random/RandPoisson.h" 23 #include "CLHEP/Random/RandPoisson.h" 24 24 25 namespace CLHEP { 25 namespace CLHEP { 26 26 27 /** 27 /** 28 * @author 28 * @author 29 * @ingroup random 29 * @ingroup random 30 */ 30 */ 31 class RandPoissonQ : public RandPoisson { 31 class RandPoissonQ : public RandPoisson { 32 32 33 public: 33 public: 34 34 35 inline RandPoissonQ ( HepRandomEngine& anEng 35 inline RandPoissonQ ( HepRandomEngine& anEngine, double b1=1.0 ); 36 inline RandPoissonQ ( HepRandomEngine* anEng 36 inline RandPoissonQ ( HepRandomEngine* anEngine, double b1=1.0 ); 37 // These constructors should be used to inst 37 // These constructors should be used to instantiate a RandPoissonQ 38 // distribution object defining a local engi 38 // distribution object defining a local engine for it. 39 // The static generator will be skipped usin 39 // The static generator will be skipped using the non-static methods 40 // defined below. 40 // defined below. 41 // If the engine is passed by pointer the co 41 // If the engine is passed by pointer the corresponding engine object 42 // will be deleted by the RandPoissonQ destr 42 // will be deleted by the RandPoissonQ destructor. 43 // If the engine is passed by reference the 43 // If the engine is passed by reference the corresponding engine object 44 // will not be deleted by the RandPoissonQ d 44 // will not be deleted by the RandPoissonQ destructor. 45 45 46 virtual ~RandPoissonQ(); 46 virtual ~RandPoissonQ(); 47 // Destructor 47 // Destructor 48 48 49 // Save and restore to/from streams 49 // Save and restore to/from streams 50 50 51 std::ostream & put ( std::ostream & os ) con 51 std::ostream & put ( std::ostream & os ) const; 52 std::istream & get ( std::istream & is ); 52 std::istream & get ( std::istream & is ); 53 53 54 // Methods to generate Poisson-distributed r 54 // Methods to generate Poisson-distributed random deviates. 55 55 56 // The method used for mu <= 100 is exact, 56 // The method used for mu <= 100 is exact, and 3-7 times faster than 57 // that used by RandPoisson. 57 // that used by RandPoisson. 58 // For mu > 100 then we use a corrected ve 58 // For mu > 100 then we use a corrected version of a 59 // (quick) Gaussian approximation. Naivel 59 // (quick) Gaussian approximation. Naively that would be: 60 // 60 // 61 // Poisson(mu) ~ floor( mu + .5 + Gaussian( 61 // Poisson(mu) ~ floor( mu + .5 + Gaussian(sqrt(mu)) ) 62 // 62 // 63 // but actually, that would give a slightl 63 // but actually, that would give a slightly incorrect sigma and a 64 // very different skew than a true Poisson 64 // very different skew than a true Poisson. Instead we return 65 // 65 // 66 // Poisson(mu) ~ floor( a0*mu + a1*g + a2*g 66 // Poisson(mu) ~ floor( a0*mu + a1*g + a2*g*g ) ) 67 // (with g a gaussian normal) 67 // (with g a gaussian normal) 68 // 68 // 69 // where a0, a1, a2 are chosen to give the 69 // where a0, a1, a2 are chosen to give the exctly correct mean, sigma, 70 // and skew for the Poisson distribution. 70 // and skew for the Poisson distribution. 71 71 72 // Static methods to shoot random values usi 72 // Static methods to shoot random values using the static generator 73 73 74 static long shoot( double mean=1.0 ); << 74 static long shoot( double m=1.0 ); 75 75 76 static void shootArray ( const int size, lo << 76 static void shootArray ( const int size, long* vect, double m=1.0 ); 77 77 78 // Static methods to shoot random values us 78 // Static methods to shoot random values using a given engine 79 // by-passing the static generator. 79 // by-passing the static generator. 80 80 81 static long shoot( HepRandomEngine* anEngin << 81 static long shoot( HepRandomEngine* anEngine, double m=1.0 ); 82 82 83 static void shootArray ( HepRandomEngine* a 83 static void shootArray ( HepRandomEngine* anEngine, 84 const int size, lo << 84 const int size, long* vect, double m=1.0 ); 85 85 86 // Methods using the localEngine to shoot r 86 // Methods using the localEngine to shoot random values, by-passing 87 // the static generator. 87 // the static generator. 88 88 89 long fire(); 89 long fire(); 90 long fire( double mean ); << 90 long fire( double m ); 91 91 92 void fireArray ( const int size, long* vect 92 void fireArray ( const int size, long* vect ); 93 void fireArray ( const int size, long* vect, << 93 void fireArray ( const int size, long* vect, double m); 94 94 95 double operator()(); 95 double operator()(); 96 double operator()( double mean ); << 96 double operator()( double m ); 97 97 98 std::string name() const; 98 std::string name() const; 99 HepRandomEngine & engine(); 99 HepRandomEngine & engine(); 100 100 101 static std::string distributionName() {retur 101 static std::string distributionName() {return "RandPoissonQ";} 102 // Provides the name of this distribution cl 102 // Provides the name of this distribution class 103 103 104 104 105 // static constants of possible interest to 105 // static constants of possible interest to users: 106 106 107 // RandPoisson will never return a deviate g 107 // RandPoisson will never return a deviate greater than this value: 108 static const double MAXIMUM_POISSON_DEVIATE; 108 static const double MAXIMUM_POISSON_DEVIATE; // Will be 2.0E9 109 109 110 static inline int tableBoundary(); 110 static inline int tableBoundary(); 111 111 112 private: 112 private: 113 113 114 // constructor helper 114 // constructor helper 115 void setupForDefaultMu(); 115 void setupForDefaultMu(); 116 116 117 // algorithm helper methods - all static sin 117 // algorithm helper methods - all static since the shoot methods mayneed them 118 static long poissonDeviateSmall ( HepRandomE 118 static long poissonDeviateSmall ( HepRandomEngine * e, double mean ); 119 static long poissonDeviateQuick ( HepRandomE 119 static long poissonDeviateQuick ( HepRandomEngine * e, double mean ); 120 static long poissonDeviateQuick ( HepRandomE 120 static long poissonDeviateQuick ( HepRandomEngine * e, 121 double A0, double A1, double A2, double si 121 double A0, double A1, double A2, double sig ); 122 122 123 // All the engine info, and the default mean 123 // All the engine info, and the default mean, are in the 124 // RandPoisson base class. 124 // RandPoisson base class. 125 125 126 // quantities for approximate Poisson by cor 126 // quantities for approximate Poisson by corrected Gaussian 127 double a0; 127 double a0; 128 double a1; 128 double a1; 129 double a2; 129 double a2; 130 double sigma; 130 double sigma; 131 131 132 // static data - constants only, so that sav 132 // static data - constants only, so that saveEngineStatus works properly! 133 133 134 // The following MUST MATCH the correspondin 134 // The following MUST MATCH the corresponding values used (in 135 // poissonTables.cc) when poissonTables.cdat 135 // poissonTables.cc) when poissonTables.cdat was created. 136 // poissonTables.cc gets these values by inc 136 // poissonTables.cc gets these values by including this header, 137 // but we must be careful not to change thes 137 // but we must be careful not to change these values, 138 // and rebuild RandPoissonQ before re-genera 138 // and rebuild RandPoissonQ before re-generating poissonTables.cdat. 139 139 140 // (These statics are given values near the 140 // (These statics are given values near the start of the .cc file) 141 141 142 static const double FIRST_MU; // lowest mu 142 static const double FIRST_MU; // lowest mu value in table 143 static const double LAST_MU; // highest mu 143 static const double LAST_MU; // highest mu value 144 static const double S; // Spacing be 144 static const double S; // Spacing between mu values 145 static const int BELOW; // Starting p 145 static const int BELOW; // Starting point for N is at mu - BELOW 146 static const int ENTRIES; // Number of 146 static const int ENTRIES; // Number of entries in each mu row 147 147 148 }; 148 }; 149 149 150 } // namespace CLHEP 150 } // namespace CLHEP 151 151 152 #include "CLHEP/Random/RandPoissonQ.icc" 152 #include "CLHEP/Random/RandPoissonQ.icc" 153 153 154 #endif 154 #endif 155 155