Geant4-11
Public Member Functions | Static Public Member Functions | Static Protected Attributes | Static Private Member Functions | Private Attributes
CLHEP::RandBinomial Class Reference

#include <RandBinomial.h>

Inheritance diagram for CLHEP::RandBinomial:
CLHEP::HepRandom

Public Member Functions

HepRandomEngineengine ()
 
double fire ()
 
double fire (long n, double p)
 
void fireArray (const int size, double *vect)
 
void fireArray (const int size, double *vect, long n, double p)
 
double flat ()
 
double flat (HepRandomEngine *theNewEngine)
 
void flatArray (const int size, double *vect)
 
void flatArray (HepRandomEngine *theNewEngine, const int size, double *vect)
 
std::istream & get (std::istream &is)
 
std::string name () const
 
double operator() ()
 
double operator() (long n, double p)
 
std::ostream & put (std::ostream &os) const
 
 RandBinomial (HepRandomEngine &anEngine, long n=1, double p=0.5)
 
 RandBinomial (HepRandomEngine *anEngine, long n=1, double p=0.5)
 
virtual ~RandBinomial ()
 

Static Public Member Functions

static int createInstance ()
 
static std::string distributionName ()
 
static HepRandomEnginegetTheEngine ()
 
static HepRandomgetTheGenerator ()
 
static long getTheSeed ()
 
static const long * getTheSeeds ()
 
static void getTheTableSeeds (long *seeds, int index)
 
static std::istream & restoreDistState (std::istream &is)
 
static void restoreEngineStatus (const char filename[]="Config.conf")
 
static std::istream & restoreFullState (std::istream &is)
 
static std::istream & restoreStaticRandomStates (std::istream &is)
 
static std::ostream & saveDistState (std::ostream &os)
 
static void saveEngineStatus (const char filename[]="Config.conf")
 
static std::ostream & saveFullState (std::ostream &os)
 
static std::ostream & saveStaticRandomStates (std::ostream &os)
 
static void setTheEngine (HepRandomEngine *theNewEngine)
 
static void setTheSeed (long seed, int lxr=3)
 
static void setTheSeeds (const long *seeds, int aux=-1)
 
static double shoot ()
 
static double shoot (HepRandomEngine *anEngine)
 
static double shoot (HepRandomEngine *anEngine, long n, double p)
 
static double shoot (long n, double p)
 
static void shootArray (const int size, double *vect, long n=1, double p=0.5)
 
static void shootArray (HepRandomEngine *anEngine, const int size, double *vect, long n=1, double p=0.5)
 
static void showEngineStatus ()
 

Static Protected Attributes

static const long seedTable [215][2]
 

Static Private Member Functions

static double genBinomial (HepRandomEngine *anEngine, long n, double p)
 

Private Attributes

long defaultN
 
double defaultP
 
std::shared_ptr< HepRandomEnginelocalEngine
 

Detailed Description

Author

Definition at line 36 of file RandBinomial.h.

Constructor & Destructor Documentation

◆ RandBinomial() [1/2]

CLHEP::RandBinomial::RandBinomial ( HepRandomEngine anEngine,
long  n = 1,
double  p = 0.5 
)
inline

◆ RandBinomial() [2/2]

CLHEP::RandBinomial::RandBinomial ( HepRandomEngine anEngine,
long  n = 1,
double  p = 0.5 
)
inline

◆ ~RandBinomial()

CLHEP::RandBinomial::~RandBinomial ( )
virtual

Definition at line 31 of file RandBinomial.cc.

31 {
32}

Member Function Documentation

◆ createInstance()

int CLHEP::HepRandom::createInstance ( )
staticinherited

◆ distributionName()

static std::string CLHEP::RandBinomial::distributionName ( )
inlinestatic

Definition at line 98 of file RandBinomial.h.

98{return "RandBinomial";}

◆ engine()

HepRandomEngine & CLHEP::RandBinomial::engine ( )
virtual

Reimplemented from CLHEP::HepRandom.

Definition at line 29 of file RandBinomial.cc.

29{return *localEngine;}
std::shared_ptr< HepRandomEngine > localEngine
Definition: RandBinomial.h:105

References localEngine.

◆ fire() [1/2]

double CLHEP::RandBinomial::fire ( )
inline

Referenced by fireArray().

◆ fire() [2/2]

double CLHEP::RandBinomial::fire ( long  n,
double  p 
)

Definition at line 44 of file RandBinomial.cc.

44 {
45 return genBinomial( localEngine.get(), n, p );
46}
static double genBinomial(HepRandomEngine *anEngine, long n, double p)

References genBinomial(), localEngine, and CLHEP::detail::n.

◆ fireArray() [1/2]

void CLHEP::RandBinomial::fireArray ( const int  size,
double *  vect 
)

Definition at line 63 of file RandBinomial.cc.

64{
65 for( double* v = vect; v != vect+size; ++v )
67}

References defaultN, defaultP, and fire().

◆ fireArray() [2/2]

void CLHEP::RandBinomial::fireArray ( const int  size,
double *  vect,
long  n,
double  p 
)

Definition at line 69 of file RandBinomial.cc.

71{
72 for( double* v = vect; v != vect+size; ++v )
73 *v = fire(n,p);
74}

References fire(), and CLHEP::detail::n.

◆ flat() [1/2]

double CLHEP::HepRandom::flat ( )
inherited

◆ flat() [2/2]

double CLHEP::HepRandom::flat ( HepRandomEngine theNewEngine)
inlineinherited

◆ flatArray() [1/2]

void CLHEP::HepRandom::flatArray ( const int  size,
double *  vect 
)
inherited

◆ flatArray() [2/2]

void CLHEP::HepRandom::flatArray ( HepRandomEngine theNewEngine,
const int  size,
double *  vect 
)
inlineinherited

◆ genBinomial()

double CLHEP::RandBinomial::genBinomial ( HepRandomEngine anEngine,
long  n,
double  p 
)
staticprivate

Definition at line 127 of file RandBinomial.cc.

128{
129/******************************************************************
130 * *
131 * Binomial-Distribution - Acceptance Rejection/Inversion *
132 * *
133 ******************************************************************
134 * *
135 * Acceptance Rejection method combined with Inversion for *
136 * generating Binomial random numbers with parameters *
137 * n (number of trials) and p (probability of success). *
138 * For min(n*p,n*(1-p)) < 10 the Inversion method is applied: *
139 * The random numbers are generated via sequential search, *
140 * starting at the lowest index k=0. The cumulative probabilities *
141 * are avoided by using the technique of chop-down. *
142 * For min(n*p,n*(1-p)) >= 10 Acceptance Rejection is used: *
143 * The algorithm is based on a hat-function which is uniform in *
144 * the centre region and exponential in the tails. *
145 * A triangular immediate acceptance region in the centre speeds *
146 * up the generation of binomial variates. *
147 * If candidate k is near the mode, f(k) is computed recursively *
148 * starting at the mode m. *
149 * The acceptance test by Stirling's formula is modified *
150 * according to W. Hoermann (1992): The generation of binomial *
151 * random variates, to appear in J. Statist. Comput. Simul. *
152 * If p < .5 the algorithm is applied to parameters n, p. *
153 * Otherwise p is replaced by 1-p, and k is replaced by n - k. *
154 * *
155 ******************************************************************
156 * *
157 * FUNCTION: - btpec samples a random number from the binomial *
158 * distribution with parameters n and p and is *
159 * valid for n*min(p,1-p) > 0. *
160 * REFERENCE: - V. Kachitvichyanukul, B.W. Schmeiser (1988): *
161 * Binomial random variate generation, *
162 * Communications of the ACM 31, 216-222. *
163 * SUBPROGRAMS: - StirlingCorrection() *
164 * ... Correction term of the Stirling *
165 * approximation for log(k!) *
166 * (series in 1/k or table values *
167 * for small k) with long int k *
168 * - anEngine ... Pointer to a (0,1)-Uniform *
169 * engine *
170 * *
171 * Implemented by H. Zechner and P. Busswald, September 1992 *
172 ******************************************************************/
173
174#define C1_3 0.33333333333333333
175#define C5_8 0.62500000000000000
176#define C1_6 0.16666666666666667
177#define DMAX_KM 20L
178
179 static CLHEP_THREAD_LOCAL long int n_last = -1L, n_prev = -1L;
180 static CLHEP_THREAD_LOCAL double par,np,p0,q,p_last = -1.0, p_prev = -1.0;
181 static CLHEP_THREAD_LOCAL long b,m,nm;
182 static CLHEP_THREAD_LOCAL double pq, rc, ss, xm, xl, xr, ll, lr, c,
183 p1, p2, p3, p4, ch;
184
185 long bh,i, K, Km, nK;
186 double f, rm, U, V, X, T, E;
187
188 if (n != n_last || p != p_last) // set-up
189 {
190 n_last = n;
191 p_last = p;
192 par=std::min(p,1.0-p);
193 q=1.0-par;
194 np = n*par;
195
196// Check for invalid input values
197
198 if( np <= 0.0 ) return (-1.0);
199
200 rm = np + par;
201 m = (long int) rm; // mode, integer
202 if (np<10)
203 {
204 p0=std::exp(n*std::log(q)); // Chop-down
205 bh=(long int)(np+10.0*std::sqrt(np*q));
206 b=std::min(n,bh);
207 }
208 else
209 {
210 rc = (n + 1.0) * (pq = par / q); // recurr. relat.
211 ss = np * q; // variance
212 i = (long int) (2.195*std::sqrt(ss) - 4.6*q); // i = p1 - 0.5
213 xm = m + 0.5;
214 xl = (double) (m - i); // limit left
215 xr = (double) (m + i + 1L); // limit right
216 f = (rm - xl) / (rm - xl*par); ll = f * (1.0 + 0.5*f);
217 f = (xr - rm) / (xr * q); lr = f * (1.0 + 0.5*f);
218 c = 0.134 + 20.5/(15.3 + (double) m); // parallelogram
219 // height
220 p1 = i + 0.5;
221 p2 = p1 * (1.0 + c + c); // probabilities
222 p3 = p2 + c/ll; // of regions 1-4
223 p4 = p3 + c/lr;
224 }
225 }
226 if( np <= 0.0 ) return (-1.0);
227 if (np<10) //Inversion Chop-down
228 {
229 double pk;
230
231 K=0;
232 pk=p0;
233 U=anEngine->flat();
234 while (U>pk)
235 {
236 ++K;
237 if (K>b)
238 {
239 U=anEngine->flat();
240 K=0;
241 pk=p0;
242 }
243 else
244 {
245 U-=pk;
246 pk=(double)(((n-K+1)*par*pk)/(K*q));
247 }
248 }
249 return ((p>0.5) ? (double)(n-K):(double)K);
250 }
251
252 for (;;)
253 {
254 V = anEngine->flat();
255 if ((U = anEngine->flat() * p4) <= p1) // triangular region
256 {
257 K=(long int) (xm - U + p1*V);
258 return ((p>0.5) ? (double)(n-K):(double)K); // immediate accept
259 }
260 if (U <= p2) // parallelogram
261 {
262 X = xl + (U - p1)/c;
263 if ((V = V*c + 1.0 - std::fabs(xm - X)/p1) >= 1.0) continue;
264 K = (long int) X;
265 }
266 else if (U <= p3) // left tail
267 {
268 if ((X = xl + std::log(V)/ll) < 0.0) continue;
269 K = (long int) X;
270 V *= (U - p2) * ll;
271 }
272 else // right tail
273 {
274 if ((K = (long int) (xr - std::log(V)/lr)) > n) continue;
275 V *= (U - p3) * lr;
276 }
277
278 // acceptance test : two cases, depending on |K - m|
279 if ((Km = std::labs(K - m)) <= DMAX_KM || Km + Km + 2L >= ss)
280 {
281
282 // computation of p(K) via recurrence relationship from the mode
283 f = 1.0; // f(m)
284 if (m < K)
285 {
286 for (i = m; i < K; )
287 {
288 if ((f *= (rc / ++i - pq)) < V) break; // multiply f
289 }
290 }
291 else
292 {
293 for (i = K; i < m; )
294 {
295 if ((V *= (rc / ++i - pq)) > f) break; // multiply V
296 }
297 }
298 if (V <= f) break; // acceptance test
299 }
300 else
301 {
302
303 // lower and upper squeeze tests, based on lower bounds for log p(K)
304 V = std::log(V);
305 T = - Km * Km / (ss + ss);
306 E = (Km / ss) * ((Km * (Km * C1_3 + C5_8) + C1_6) / ss + 0.5);
307 if (V <= T - E) break;
308 if (V <= T + E)
309 {
310 if (n != n_prev || par != p_prev)
311 {
312 n_prev = n;
313 p_prev = par;
314
315 nm = n - m + 1L;
316 ch = xm * std::log((m + 1.0)/(pq * nm)) +
318 }
319 nK = n - K + 1L;
320
321 // computation of log f(K) via Stirling's formula
322 // final acceptance-rejection test
323 if (V <= ch + (n + 1.0)*std::log((double) nm / (double) nK) +
324 (K + 0.5)*std::log(nK * pq / (K + 1.0)) -
325 StirlingCorrection(K + 1L) - StirlingCorrection(nK)) break;
326 }
327 }
328 }
329 return ((p>0.5) ? (double)(n-K):(double)K);
330}
#define C5_8
#define C1_6
#define C1_3
#define DMAX_KM
static double StirlingCorrection(long int k)
Definition: RandBinomial.cc:94
static constexpr double m
static constexpr double nm
Definition: SystemOfUnits.h:93
T min(const T t1, const T t2)
brief Return the smallest of the two arguments
#define CLHEP_THREAD_LOCAL
Definition: thread_local.h:13

References C1_3, C1_6, C5_8, CLHEP_THREAD_LOCAL, DMAX_KM, CLHEP::HepRandomEngine::flat(), CLHEP::m, G4INCL::Math::min(), CLHEP::detail::n, CLHEP::nm, CLHEP::StirlingCorrection(), and anonymous_namespace{G4PionRadiativeDecayChannel.cc}::xl.

Referenced by fire(), and shoot().

◆ get()

std::istream & CLHEP::RandBinomial::get ( std::istream &  is)
virtual

Reimplemented from CLHEP::HepRandom.

Definition at line 343 of file RandBinomial.cc.

343 {
344 std::string inName;
345 is >> inName;
346 if (inName != name()) {
347 is.clear(std::ios::badbit | is.rdstate());
348 std::cerr << "Mismatch when expecting to read state of a "
349 << name() << " distribution\n"
350 << "Name found was " << inName
351 << "\nistream is left in the badbit state\n";
352 return is;
353 }
354 if (possibleKeywordInput(is, "Uvec", defaultN)) {
355 std::vector<unsigned long> t(2);
356 is >> defaultN >> defaultP;
357 is >> t[0] >> t[1]; defaultP = DoubConv::longs2double(t);
358 return is;
359 }
360 // is >> defaultN encompassed by possibleKeywordInput
361 is >> defaultP;
362 return is;
363}
static double longs2double(const std::vector< unsigned long > &v)
Definition: DoubConv.cc:110
std::string name() const
Definition: RandBinomial.cc:28
bool possibleKeywordInput(IS &is, const std::string &key, T &t)
Definition: RandomEngine.h:166

References defaultN, defaultP, CLHEP::DoubConv::longs2double(), name(), and CLHEP::possibleKeywordInput().

◆ getTheEngine()

HepRandomEngine * CLHEP::HepRandom::getTheEngine ( )
staticinherited

◆ getTheGenerator()

HepRandom * CLHEP::HepRandom::getTheGenerator ( )
staticinherited

◆ getTheSeed()

long CLHEP::HepRandom::getTheSeed ( )
staticinherited

◆ getTheSeeds()

const long * CLHEP::HepRandom::getTheSeeds ( )
staticinherited

◆ getTheTableSeeds()

void CLHEP::HepRandom::getTheTableSeeds ( long *  seeds,
int  index 
)
staticinherited

Definition at line 254 of file Random.cc.

255{
256 if ((index >= 0) && (index < 215)) {
257 seeds[0] = seedTable[index][0];
258 seeds[1] = seedTable[index][1];
259 }
260 else seeds = NULL;
261}
static const long seedTable[215][2]
Definition: Random.h:156

Referenced by CLHEP::HepJamesRandom::HepJamesRandom(), CLHEP::MTwistEngine::MTwistEngine(), CLHEP::RanecuEngine::RanecuEngine(), CLHEP::Ranlux64Engine::Ranlux64Engine(), CLHEP::RanluxEngine::RanluxEngine(), and CLHEP::RanecuEngine::setSeed().

◆ name()

std::string CLHEP::RandBinomial::name ( ) const
virtual

Reimplemented from CLHEP::HepRandom.

Definition at line 28 of file RandBinomial.cc.

28{return "RandBinomial";}

Referenced by source.g4viscp.G4Scene::create_scene(), get(), mcscore.MCParticle::printout(), put(), and source.g4viscp.G4Scene::update_scene().

◆ operator()() [1/2]

double CLHEP::RandBinomial::operator() ( )
inlinevirtual

Reimplemented from CLHEP::HepRandom.

◆ operator()() [2/2]

double CLHEP::RandBinomial::operator() ( long  n,
double  p 
)
inline

◆ put()

std::ostream & CLHEP::RandBinomial::put ( std::ostream &  os) const
virtual

Reimplemented from CLHEP::HepRandom.

Definition at line 332 of file RandBinomial.cc.

332 {
333 int pr=os.precision(20);
334 std::vector<unsigned long> t(2);
335 os << " " << name() << "\n";
336 os << "Uvec" << "\n";
338 os << defaultN << " " << defaultP << " " << t[0] << " " << t[1] << "\n";
339 os.precision(pr);
340 return os;
341}
static std::vector< unsigned long > dto2longs(double d)
Definition: DoubConv.cc:94

References defaultN, defaultP, CLHEP::DoubConv::dto2longs(), and name().

◆ restoreDistState()

static std::istream & CLHEP::HepRandom::restoreDistState ( std::istream &  is)
inlinestaticinherited

Definition at line 136 of file Random.h.

136{return is;}

◆ restoreEngineStatus()

void CLHEP::HepRandom::restoreEngineStatus ( const char  filename[] = "Config.conf")
staticinherited

Definition at line 283 of file Random.cc.

284{
285 theDefaults().theEngine->restoreStatus( filename );
286}

References CLHEP::anonymous_namespace{Random.cc}::theDefaults(), and CLHEP::anonymous_namespace{Random.cc}::defaults::theEngine.

◆ restoreFullState()

std::istream & CLHEP::HepRandom::restoreFullState ( std::istream &  is)
staticinherited

Definition at line 293 of file Random.cc.

293 {
294 is >> *getTheEngine();
295 return is;
296}
static HepRandomEngine * getTheEngine()
Definition: Random.cc:268

Referenced by CLHEP::RandFlat::restoreFullState(), and CLHEP::RandGauss::restoreFullState().

◆ restoreStaticRandomStates()

std::istream & CLHEP::HepRandom::restoreStaticRandomStates ( std::istream &  is)
staticinherited

Definition at line 302 of file Random.cc.

302 {
304}
static std::istream & restore(std::istream &is)

References CLHEP::StaticRandomStates::restore().

◆ saveDistState()

static std::ostream & CLHEP::HepRandom::saveDistState ( std::ostream &  os)
inlinestaticinherited

Definition at line 133 of file Random.h.

133{return os;}

◆ saveEngineStatus()

void CLHEP::HepRandom::saveEngineStatus ( const char  filename[] = "Config.conf")
staticinherited

◆ saveFullState()

std::ostream & CLHEP::HepRandom::saveFullState ( std::ostream &  os)
staticinherited

◆ saveStaticRandomStates()

std::ostream & CLHEP::HepRandom::saveStaticRandomStates ( std::ostream &  os)
staticinherited

Definition at line 298 of file Random.cc.

298 {
299 return StaticRandomStates::save(os);
300}
static std::ostream & save(std::ostream &os)

References CLHEP::StaticRandomStates::save().

◆ setTheEngine()

void CLHEP::HepRandom::setTheEngine ( HepRandomEngine theNewEngine)
staticinherited

Definition at line 273 of file Random.cc.

274{
275 theDefaults().theEngine.reset( theNewEngine, do_nothing_deleter() );
276}

References CLHEP::anonymous_namespace{Random.cc}::theDefaults(), and CLHEP::anonymous_namespace{Random.cc}::defaults::theEngine.

Referenced by CLHEP::StaticRandomStates::restore().

◆ setTheSeed()

void CLHEP::HepRandom::setTheSeed ( long  seed,
int  lxr = 3 
)
staticinherited

Definition at line 234 of file Random.cc.

235{
236 theDefaults().theEngine->setSeed(seed,lux);
237}
static constexpr double lux

References CLHEP::lux, CLHEP::anonymous_namespace{Random.cc}::theDefaults(), and CLHEP::anonymous_namespace{Random.cc}::defaults::theEngine.

◆ setTheSeeds()

void CLHEP::HepRandom::setTheSeeds ( const long *  seeds,
int  aux = -1 
)
staticinherited

◆ shoot() [1/4]

static double CLHEP::RandBinomial::shoot ( )
inlinestatic

Referenced by shootArray().

◆ shoot() [2/4]

static double CLHEP::RandBinomial::shoot ( HepRandomEngine anEngine)
inlinestatic

◆ shoot() [3/4]

double CLHEP::RandBinomial::shoot ( HepRandomEngine anEngine,
long  n,
double  p 
)
static

Definition at line 34 of file RandBinomial.cc.

35 {
36 return genBinomial( anEngine, n, p );
37}

References genBinomial(), and CLHEP::detail::n.

◆ shoot() [4/4]

double CLHEP::RandBinomial::shoot ( long  n,
double  p 
)
static

Definition at line 39 of file RandBinomial.cc.

39 {
40 HepRandomEngine *anEngine = HepRandom::getTheEngine();
41 return genBinomial( anEngine, n, p );
42}

References genBinomial(), CLHEP::HepRandom::getTheEngine(), and CLHEP::detail::n.

◆ shootArray() [1/2]

void CLHEP::RandBinomial::shootArray ( const int  size,
double *  vect,
long  n = 1,
double  p = 0.5 
)
static

Definition at line 48 of file RandBinomial.cc.

50{
51 for( double* v = vect; v != vect+size; ++v )
52 *v = shoot(n,p);
53}
static double shoot()

References CLHEP::detail::n, and shoot().

◆ shootArray() [2/2]

void CLHEP::RandBinomial::shootArray ( HepRandomEngine anEngine,
const int  size,
double *  vect,
long  n = 1,
double  p = 0.5 
)
static

Definition at line 55 of file RandBinomial.cc.

58{
59 for( double* v = vect; v != vect+size; ++v )
60 *v = shoot(anEngine,n,p);
61}

References CLHEP::detail::n, and shoot().

◆ showEngineStatus()

void CLHEP::HepRandom::showEngineStatus ( )
staticinherited

Field Documentation

◆ defaultN

long CLHEP::RandBinomial::defaultN
private

Definition at line 106 of file RandBinomial.h.

Referenced by fireArray(), get(), and put().

◆ defaultP

double CLHEP::RandBinomial::defaultP
private

Definition at line 107 of file RandBinomial.h.

Referenced by fireArray(), get(), and put().

◆ localEngine

std::shared_ptr<HepRandomEngine> CLHEP::RandBinomial::localEngine
private

Definition at line 105 of file RandBinomial.h.

Referenced by engine(), and fire().

◆ seedTable

const long CLHEP::HepRandom::seedTable
staticprotectedinherited

Definition at line 156 of file Random.h.


The documentation for this class was generated from the following files: