[ VIGRA Homepage | Function Index | Class Index | Namespaces | File List | Main Page ]
vigra/random.hxx | ![]() |
00001 /************************************************************************/ 00002 /* */ 00003 /* Copyright 2008 by Ullrich Koethe */ 00004 /* */ 00005 /* This file is part of the VIGRA computer vision library. */ 00006 /* The VIGRA Website is */ 00007 /* http://hci.iwr.uni-heidelberg.de/vigra/ */ 00008 /* Please direct questions, bug reports, and contributions to */ 00009 /* ullrich.koethe@iwr.uni-heidelberg.de or */ 00010 /* vigra@informatik.uni-hamburg.de */ 00011 /* */ 00012 /* Permission is hereby granted, free of charge, to any person */ 00013 /* obtaining a copy of this software and associated documentation */ 00014 /* files (the "Software"), to deal in the Software without */ 00015 /* restriction, including without limitation the rights to use, */ 00016 /* copy, modify, merge, publish, distribute, sublicense, and/or */ 00017 /* sell copies of the Software, and to permit persons to whom the */ 00018 /* Software is furnished to do so, subject to the following */ 00019 /* conditions: */ 00020 /* */ 00021 /* The above copyright notice and this permission notice shall be */ 00022 /* included in all copies or substantial portions of the */ 00023 /* Software. */ 00024 /* */ 00025 /* THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND */ 00026 /* EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES */ 00027 /* OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND */ 00028 /* NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT */ 00029 /* HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, */ 00030 /* WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING */ 00031 /* FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR */ 00032 /* OTHER DEALINGS IN THE SOFTWARE. */ 00033 /* */ 00034 /************************************************************************/ 00035 00036 00037 #ifndef VIGRA_RANDOM_HXX 00038 #define VIGRA_RANDOM_HXX 00039 00040 #include "mathutil.hxx" 00041 #include "functortraits.hxx" 00042 00043 #include <ctime> 00044 00045 namespace vigra { 00046 00047 enum RandomSeedTag { RandomSeed }; 00048 00049 namespace detail { 00050 00051 enum RandomEngineTag { TT800, MT19937 }; 00052 00053 template<RandomEngineTag EngineTag> 00054 struct RandomState; 00055 00056 template <RandomEngineTag EngineTag> 00057 void seed(UInt32 theSeed, RandomState<EngineTag> & engine) 00058 { 00059 engine.state_[0] = theSeed; 00060 for(UInt32 i=1; i<RandomState<EngineTag>::N; ++i) 00061 { 00062 engine.state_[i] = 1812433253U * (engine.state_[i-1] ^ (engine.state_[i-1] >> 30)) + i; 00063 } 00064 } 00065 00066 template <class Iterator, RandomEngineTag EngineTag> 00067 void seed(Iterator init, UInt32 key_length, RandomState<EngineTag> & engine) 00068 { 00069 const UInt32 N = RandomState<EngineTag>::N; 00070 int k = (int)std::max(N, key_length); 00071 UInt32 i = 1, j = 0; 00072 Iterator data = init; 00073 for (; k; --k) 00074 { 00075 engine.state_[i] = (engine.state_[i] ^ ((engine.state_[i-1] ^ (engine.state_[i-1] >> 30)) * 1664525U)) 00076 + *data + j; /* non linear */ 00077 ++i; ++j; ++data; 00078 00079 if (i >= N) 00080 { 00081 engine.state_[0] = engine.state_[N-1]; 00082 i=1; 00083 } 00084 if (j>=key_length) 00085 { 00086 j=0; 00087 data = init; 00088 } 00089 } 00090 00091 for (k=N-1; k; --k) 00092 { 00093 engine.state_[i] = (engine.state_[i] ^ ((engine.state_[i-1] ^ (engine.state_[i-1] >> 30)) * 1566083941U)) 00094 - i; /* non linear */ 00095 ++i; 00096 if (i>=N) 00097 { 00098 engine.state_[0] = engine.state_[N-1]; 00099 i=1; 00100 } 00101 } 00102 00103 engine.state_[0] = 0x80000000U; /* MSB is 1; assuring non-zero initial array */ 00104 } 00105 00106 template <RandomEngineTag EngineTag> 00107 void seed(RandomSeedTag, RandomState<EngineTag> & engine) 00108 { 00109 static UInt32 globalCount = 0; 00110 UInt32 init[3] = { (UInt32)time(0), (UInt32)clock(), ++globalCount }; 00111 seed(init, 3, engine); 00112 } 00113 00114 00115 /* Tempered twister TT800 by M. Matsumoto */ 00116 template<> 00117 struct RandomState<TT800> 00118 { 00119 static const UInt32 N = 25, M = 7; 00120 00121 mutable UInt32 state_[N]; 00122 mutable UInt32 current_; 00123 00124 RandomState() 00125 : current_(0) 00126 { 00127 UInt32 seeds[N] = { 00128 0x95f24dab, 0x0b685215, 0xe76ccae7, 0xaf3ec239, 0x715fad23, 00129 0x24a590ad, 0x69e4b5ef, 0xbf456141, 0x96bc1b7b, 0xa7bdf825, 00130 0xc1de75b7, 0x8858a9c9, 0x2da87693, 0xb657f9dd, 0xffdc8a9f, 00131 0x8121da71, 0x8b823ecb, 0x885d05f5, 0x4e20cd47, 0x5a9ad5d9, 00132 0x512c0c03, 0xea857ccd, 0x4cc1d30f, 0x8891a8a1, 0xa6b7aadb 00133 }; 00134 00135 for(UInt32 i=0; i<N; ++i) 00136 state_[i] = seeds[i]; 00137 } 00138 00139 protected: 00140 00141 UInt32 get() const 00142 { 00143 if(current_ == N) 00144 generateNumbers<void>(); 00145 00146 UInt32 y = state_[current_++]; 00147 y ^= (y << 7) & 0x2b5b2500; 00148 y ^= (y << 15) & 0xdb8b0000; 00149 return y ^ (y >> 16); 00150 } 00151 00152 template <class DUMMY> 00153 void generateNumbers() const; 00154 00155 void seedImpl(RandomSeedTag) 00156 { 00157 seed(RandomSeed, *this); 00158 } 00159 00160 void seedImpl(UInt32 theSeed) 00161 { 00162 seed(theSeed, *this); 00163 } 00164 00165 template<class Iterator> 00166 void seedImpl(Iterator init, UInt32 length) 00167 { 00168 seed(init, length, *this); 00169 } 00170 }; 00171 00172 template <class DUMMY> 00173 void RandomState<TT800>::generateNumbers() const 00174 { 00175 UInt32 mag01[2]= { 0x0, 0x8ebfd028 }; 00176 00177 for(UInt32 i=0; i<N-M; ++i) 00178 { 00179 state_[i] = state_[i+M] ^ (state_[i] >> 1) ^ mag01[state_[i] % 2]; 00180 } 00181 for (UInt32 i=N-M; i<N; ++i) 00182 { 00183 state_[i] = state_[i+(M-N)] ^ (state_[i] >> 1) ^ mag01[state_[i] % 2]; 00184 } 00185 current_ = 0; 00186 } 00187 00188 /* Mersenne twister MT19937 by M. Matsumoto */ 00189 template<> 00190 struct RandomState<MT19937> 00191 { 00192 static const UInt32 N = 624, M = 397; 00193 00194 mutable UInt32 state_[N]; 00195 mutable UInt32 current_; 00196 00197 RandomState() 00198 : current_(0) 00199 { 00200 seed(19650218U, *this); 00201 } 00202 00203 protected: 00204 00205 UInt32 get() const 00206 { 00207 if(current_ == N) 00208 generateNumbers<void>(); 00209 00210 UInt32 x = state_[current_++]; 00211 x ^= (x >> 11); 00212 x ^= (x << 7) & 0x9D2C5680U; 00213 x ^= (x << 15) & 0xEFC60000U; 00214 return x ^ (x >> 18); 00215 } 00216 00217 template <class DUMMY> 00218 void generateNumbers() const; 00219 00220 static UInt32 twiddle(UInt32 u, UInt32 v) 00221 { 00222 return (((u & 0x80000000U) | (v & 0x7FFFFFFFU)) >> 1) 00223 ^ ((v & 1U) ? 0x9908B0DFU : 0x0U); 00224 } 00225 00226 void seedImpl(RandomSeedTag) 00227 { 00228 seed(RandomSeed, *this); 00229 generateNumbers<void>(); 00230 } 00231 00232 void seedImpl(UInt32 theSeed) 00233 { 00234 seed(theSeed, *this); 00235 generateNumbers<void>(); 00236 } 00237 00238 template<class Iterator> 00239 void seedImpl(Iterator init, UInt32 length) 00240 { 00241 seed(19650218U, *this); 00242 seed(init, length, *this); 00243 generateNumbers<void>(); 00244 } 00245 }; 00246 00247 template <class DUMMY> 00248 void RandomState<MT19937>::generateNumbers() const 00249 { 00250 for (unsigned int i = 0; i < (N - M); ++i) 00251 { 00252 state_[i] = state_[i + M] ^ twiddle(state_[i], state_[i + 1]); 00253 } 00254 for (unsigned int i = N - M; i < (N - 1); ++i) 00255 { 00256 state_[i] = state_[i + M - N] ^ twiddle(state_[i], state_[i + 1]); 00257 } 00258 state_[N - 1] = state_[M - 1] ^ twiddle(state_[N - 1], state_[0]); 00259 current_ = 0; 00260 } 00261 00262 } // namespace detail 00263 00264 00265 /** \addtogroup RandomNumberGeneration Random Number Generation 00266 00267 High-quality random number generators and related functors. 00268 */ 00269 //@{ 00270 00271 /** Generic random number generator. 00272 00273 The actual generator is passed in the template argument <tt>Engine</tt>. Two generators 00274 are currently available: 00275 <ul> 00276 <li> <tt>RandomMT19937</tt>: The state-of-the-art <a href="http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/emt.html">Mersenne Twister</a> with a state length of 2<sup>19937</sup> and very high statistical quality. 00277 <li> <tt>RandomTT800</tt>: (default) The Tempered Twister, a simpler predecessor of the Mersenne Twister with period length 2<sup>800</sup>. 00278 </ul> 00279 00280 Both generators have been designed by <a href="http://www.math.sci.hiroshima-u.ac.jp/~m-mat/eindex.html">Makoto Matsumoto</a>. 00281 00282 <b>Traits defined:</b> 00283 00284 <tt>FunctorTraits<RandomNumberGenerator<Engine> >::isInitializer</tt> is true (<tt>VigraTrueType</tt>). 00285 */ 00286 template <class Engine = detail::RandomState<detail::TT800> > 00287 class RandomNumberGenerator 00288 : public Engine 00289 { 00290 mutable double normalCached_; 00291 mutable bool normalCachedValid_; 00292 00293 public: 00294 00295 /** Create a new random generator object with standard seed. 00296 00297 Due to standard seeding, the random numbers generated will always be the same. 00298 This is useful for debugging. 00299 */ 00300 RandomNumberGenerator() 00301 : normalCached_(0.0), 00302 normalCachedValid_(false) 00303 {} 00304 00305 /** Create a new random generator object with a random seed. 00306 00307 The seed is obtained from the machines current <tt>clock()</tt> and <tt>time()</tt> 00308 values. 00309 00310 <b>Usage:</b> 00311 \code 00312 RandomNumberGenerator<> rnd = RandomNumberGenerator<>(RandomSeed); 00313 \endcode 00314 */ 00315 RandomNumberGenerator(RandomSeedTag) 00316 : normalCached_(0.0), 00317 normalCachedValid_(false) 00318 { 00319 this->seedImpl(RandomSeed); 00320 } 00321 00322 /** Create a new random generator object from the given seed. 00323 00324 The same seed will always produce identical random sequences. 00325 */ 00326 RandomNumberGenerator(UInt32 theSeed) 00327 : normalCached_(0.0), 00328 normalCachedValid_(false) 00329 { 00330 this->seedImpl(theSeed); 00331 } 00332 00333 /** Create a new random generator object from the given seed sequence. 00334 00335 Longer seed sequences lead to better initialization in the sense that the generator's 00336 state space is covered much better than is possible with 32-bit seeds alone. 00337 */ 00338 template<class Iterator> 00339 RandomNumberGenerator(Iterator init, UInt32 length) 00340 : normalCached_(0.0), 00341 normalCachedValid_(false) 00342 { 00343 this->seedImpl(init, length); 00344 } 00345 00346 00347 /** Re-initialize the random generator object with a random seed. 00348 00349 The seed is obtained from the machines current <tt>clock()</tt> and <tt>time()</tt> 00350 values. 00351 00352 <b>Usage:</b> 00353 \code 00354 RandomNumberGenerator<> rnd = ...; 00355 ... 00356 rnd.seed(RandomSeed); 00357 \endcode 00358 */ 00359 void seed(RandomSeedTag) 00360 { 00361 this->seedImpl(RandomSeed); 00362 normalCachedValid_ = false; 00363 } 00364 00365 /** Re-initialize the random generator object from the given seed. 00366 00367 The same seed will always produce identical random sequences. 00368 */ 00369 void seed(UInt32 theSeed) 00370 { 00371 this->seedImpl(theSeed); 00372 normalCachedValid_ = false; 00373 } 00374 00375 /** Re-initialize the random generator object from the given seed sequence. 00376 00377 Longer seed sequences lead to better initialization in the sense that the generator's 00378 state space is covered much better than is possible with 32-bit seeds alone. 00379 */ 00380 template<class Iterator> 00381 void seed(Iterator init, UInt32 length) 00382 { 00383 this->seedImpl(init, length); 00384 normalCachedValid_ = false; 00385 } 00386 00387 /** Return a uniformly distributed integer random number in [0, 2<sup>32</sup>). 00388 00389 That is, 0 <= i < 2<sup>32</sup>. 00390 */ 00391 UInt32 operator()() const 00392 { 00393 return this->get(); 00394 } 00395 00396 /** Return a uniformly distributed integer random number in [0, 2<sup>32</sup>). 00397 00398 That is, 0 <= i < 2<sup>32</sup>. 00399 */ 00400 UInt32 uniformInt() const 00401 { 00402 return this->get(); 00403 } 00404 00405 00406 #if 0 // difficult implementation necessary if low bits are not sufficiently random 00407 // in [0,beyond) 00408 UInt32 uniformInt(UInt32 beyond) const 00409 { 00410 if(beyond < 2) 00411 return 0; 00412 00413 UInt32 factor = factorForUniformInt(beyond); 00414 UInt32 res = this->get() / factor; 00415 00416 // Use rejection method to avoid quantization bias. 00417 // On average, we will need two raw random numbers to generate one. 00418 while(res >= beyond) 00419 res = this->get() / factor; 00420 return res; 00421 } 00422 #endif /* #if 0 */ 00423 00424 /** Return a uniformly distributed integer random number in [0, <tt>beyond</tt>). 00425 00426 That is, 0 <= i < <tt>beyond</tt>. 00427 */ 00428 UInt32 uniformInt(UInt32 beyond) const 00429 { 00430 // in [0,beyond) -- simple implementation when low bits are sufficiently random, which is 00431 // the case for TT800 and MT19937 00432 if(beyond < 2) 00433 return 0; 00434 00435 UInt32 remainder = (NumericTraits<UInt32>::max() - beyond + 1) % beyond; 00436 UInt32 lastSafeValue = NumericTraits<UInt32>::max() - remainder; 00437 UInt32 res = this->get(); 00438 00439 // Use rejection method to avoid quantization bias. 00440 // We will need two raw random numbers in amortized worst case. 00441 while(res > lastSafeValue) 00442 res = this->get(); 00443 return res % beyond; 00444 } 00445 00446 /** Return a uniformly distributed double-precision random number in [0.0, 1.0). 00447 00448 That is, 0.0 <= i < 1.0. All 53-bit bits of the mantissa are random (two 32-bit integers are used to 00449 create this number). 00450 */ 00451 double uniform53() const 00452 { 00453 // make full use of the entire 53-bit mantissa of a double, by Isaku Wada 00454 return ( (this->get() >> 5) * 67108864.0 + (this->get() >> 6)) * (1.0/9007199254740992.0); 00455 } 00456 00457 /** Return a uniformly distributed double-precision random number in [0.0, 1.0]. 00458 00459 That is, 0.0 <= i <= 1.0. This nuber is computed by <tt>uniformInt()</tt> / 2<sup>32</sup>, 00460 so it has effectively only 32 random bits. 00461 */ 00462 double uniform() const 00463 { 00464 return (double)this->get() / 4294967295.0; 00465 } 00466 00467 /** Return a uniformly distributed double-precision random number in [lower, upper]. 00468 00469 That is, <tt>lower</tt> <= i <= <tt>upper</tt>. This number is computed 00470 from <tt>uniform()</tt>, so it has effectively only 32 random bits. 00471 */ 00472 double uniform(double lower, double upper) const 00473 { 00474 vigra_precondition(lower < upper, 00475 "RandomNumberGenerator::uniform(): lower bound must be smaller than upper bound."); 00476 return uniform() * (upper-lower) + lower; 00477 } 00478 00479 /** Return a standard normal variate (Gaussian) random number. 00480 00481 Mean is zero, standard deviation is 1.0. It uses the polar form of the 00482 Box-Muller transform. 00483 */ 00484 double normal() const; 00485 00486 /** Return a normal variate (Gaussian) random number with the given mean and standard deviation. 00487 00488 It uses the polar form of the Box-Muller transform. 00489 */ 00490 double normal(double mean, double stddev) const 00491 { 00492 vigra_precondition(stddev > 0.0, 00493 "RandomNumberGenerator::normal(): standard deviation must be positive."); 00494 return normal()*stddev + mean; 00495 } 00496 00497 /** Access the global (program-wide) instance of the present random number generator. 00498 00499 Normally, you will create a local generator by one of the constructor calls. But sometimes 00500 it is useful to have all program parts access the same generator. 00501 */ 00502 static RandomNumberGenerator & global() 00503 { 00504 static RandomNumberGenerator generator; 00505 return generator; 00506 } 00507 00508 static UInt32 factorForUniformInt(UInt32 range) 00509 { 00510 return (range > 2147483648U || range == 0) 00511 ? 1 00512 : 2*(2147483648U / ceilPower2(range)); 00513 } 00514 }; 00515 00516 template <class Engine> 00517 double RandomNumberGenerator<Engine>::normal() const 00518 { 00519 if(normalCachedValid_) 00520 { 00521 normalCachedValid_ = false; 00522 return normalCached_; 00523 } 00524 else 00525 { 00526 double x1, x2, w; 00527 do 00528 { 00529 x1 = uniform(-1.0, 1.0); 00530 x2 = uniform(-1.0, 1.0); 00531 w = x1 * x1 + x2 * x2; 00532 } 00533 while ( w > 1.0 || w == 0.0); 00534 00535 w = std::sqrt( -2.0 * std::log( w ) / w ); 00536 00537 normalCached_ = x2 * w; 00538 normalCachedValid_ = true; 00539 00540 return x1 * w; 00541 } 00542 } 00543 00544 /** Shorthand for the TT800 random number generator class. 00545 */ 00546 typedef RandomNumberGenerator<> RandomTT800; 00547 00548 /** Shorthand for the TT800 random number generator class (same as RandomTT800). 00549 */ 00550 typedef RandomNumberGenerator<> TemperedTwister; 00551 00552 /** Shorthand for the MT19937 random number generator class. 00553 */ 00554 typedef RandomNumberGenerator<detail::RandomState<detail::MT19937> > RandomMT19937; 00555 00556 /** Shorthand for the MT19937 random number generator class (same as RandomMT19937). 00557 */ 00558 typedef RandomNumberGenerator<detail::RandomState<detail::MT19937> > MersenneTwister; 00559 00560 /** Access the global (program-wide) instance of the TT800 random number generator. 00561 */ 00562 inline RandomTT800 & randomTT800() { return RandomTT800::global(); } 00563 00564 /** Access the global (program-wide) instance of the MT19937 random number generator. 00565 */ 00566 inline RandomMT19937 & randomMT19937() { return RandomMT19937::global(); } 00567 00568 template <class Engine> 00569 class FunctorTraits<RandomNumberGenerator<Engine> > 00570 { 00571 public: 00572 typedef RandomNumberGenerator<Engine> type; 00573 00574 typedef VigraTrueType isInitializer; 00575 00576 typedef VigraFalseType isUnaryFunctor; 00577 typedef VigraFalseType isBinaryFunctor; 00578 typedef VigraFalseType isTernaryFunctor; 00579 00580 typedef VigraFalseType isUnaryAnalyser; 00581 typedef VigraFalseType isBinaryAnalyser; 00582 typedef VigraFalseType isTernaryAnalyser; 00583 }; 00584 00585 00586 /** Functor to create uniformely distributed integer random numbers. 00587 00588 This functor encapsulates the appropriate functions of the given random number 00589 <tt>Engine</tt> (usually <tt>RandomTT800</tt> or <tt>RandomMT19937</tt>) 00590 in an STL-compatible interface. 00591 00592 <b>Traits defined:</b> 00593 00594 <tt>FunctorTraits<UniformIntRandomFunctor<Engine> >::isInitializer</tt> and 00595 <tt>FunctorTraits<UniformIntRandomFunctor<Engine> >::isUnaryFunctor</tt> are true (<tt>VigraTrueType</tt>). 00596 */ 00597 template <class Engine = RandomTT800> 00598 class UniformIntRandomFunctor 00599 { 00600 UInt32 lower_, difference_, factor_; 00601 Engine const & generator_; 00602 bool useLowBits_; 00603 00604 public: 00605 00606 typedef UInt32 argument_type; ///< STL required functor argument type 00607 typedef UInt32 result_type; ///< STL required functor result type 00608 00609 /** Create functor for uniform random integers in the range [0, 2<sup>32</sup>) 00610 using the given engine. 00611 00612 That is, the generated numbers satisfy 0 <= i < 2<sup>32</sup>. 00613 */ 00614 explicit UniformIntRandomFunctor(Engine const & generator = Engine::global() ) 00615 : lower_(0), difference_(0xffffffff), factor_(1), 00616 generator_(generator), 00617 useLowBits_(true) 00618 {} 00619 00620 /** Create functor for uniform random integers in the range [<tt>lower</tt>, <tt>upper</tt>] 00621 using the given engine. 00622 00623 That is, the generated numbers satisfy <tt>lower</tt> <= i <= <tt>upper</tt>. 00624 \a useLowBits should be set to <tt>false</tt> when the engine generates 00625 random numbers whose low bits are significantly less random than the high 00626 bits. This does not apply to <tt>RandomTT800</tt> and <tt>RandomMT19937</tt>, 00627 but is necessary for simpler linear congruential generators. 00628 */ 00629 UniformIntRandomFunctor(UInt32 lower, UInt32 upper, 00630 Engine const & generator = Engine::global(), 00631 bool useLowBits = true) 00632 : lower_(lower), difference_(upper-lower), 00633 factor_(Engine::factorForUniformInt(difference_ + 1)), 00634 generator_(generator), 00635 useLowBits_(useLowBits) 00636 { 00637 vigra_precondition(lower < upper, 00638 "UniformIntRandomFunctor(): lower bound must be smaller than upper bound."); 00639 } 00640 00641 /** Return a random number as specified in the constructor. 00642 */ 00643 UInt32 operator()() const 00644 { 00645 if(difference_ == 0xffffffff) // lower_ is necessarily 0 00646 return generator_(); 00647 else if(useLowBits_) 00648 return generator_.uniformInt(difference_+1) + lower_; 00649 else 00650 { 00651 UInt32 res = generator_() / factor_; 00652 00653 // Use rejection method to avoid quantization bias. 00654 // On average, we will need two raw random numbers to generate one. 00655 while(res > difference_) 00656 res = generator_() / factor_; 00657 return res + lower_; 00658 } 00659 } 00660 00661 /** Return a uniformly distributed integer random number in the range [0, <tt>beyond</tt>). 00662 00663 That is, 0 <= i < <tt>beyond</tt>. This is a required interface for 00664 <tt>std::random_shuffle</tt>. It ignores the limits specified 00665 in the constructor and the flag <tt>useLowBits</tt>. 00666 */ 00667 UInt32 operator()(UInt32 beyond) const 00668 { 00669 if(beyond < 2) 00670 return 0; 00671 00672 return generator_.uniformInt(beyond); 00673 } 00674 }; 00675 00676 template <class Engine> 00677 class FunctorTraits<UniformIntRandomFunctor<Engine> > 00678 { 00679 public: 00680 typedef UniformIntRandomFunctor<Engine> type; 00681 00682 typedef VigraTrueType isInitializer; 00683 00684 typedef VigraTrueType isUnaryFunctor; 00685 typedef VigraFalseType isBinaryFunctor; 00686 typedef VigraFalseType isTernaryFunctor; 00687 00688 typedef VigraFalseType isUnaryAnalyser; 00689 typedef VigraFalseType isBinaryAnalyser; 00690 typedef VigraFalseType isTernaryAnalyser; 00691 }; 00692 00693 /** Functor to create uniformely distributed double-precision random numbers. 00694 00695 This functor encapsulates the function <tt>uniform()</tt> of the given random number 00696 <tt>Engine</tt> (usually <tt>RandomTT800</tt> or <tt>RandomMT19937</tt>) 00697 in an STL-compatible interface. 00698 00699 <b>Traits defined:</b> 00700 00701 <tt>FunctorTraits<UniformIntRandomFunctor<Engine> >::isInitializer</tt> is true (<tt>VigraTrueType</tt>). 00702 */ 00703 template <class Engine = RandomTT800> 00704 class UniformRandomFunctor 00705 { 00706 double offset_, scale_; 00707 Engine const & generator_; 00708 00709 public: 00710 00711 typedef double result_type; ///< STL required functor result type 00712 00713 /** Create functor for uniform random double-precision numbers in the range [0.0, 1.0] 00714 using the given engine. 00715 00716 That is, the generated numbers satisfy 0.0 <= i <= 1.0. 00717 */ 00718 UniformRandomFunctor(Engine const & generator = Engine::global()) 00719 : offset_(0.0), 00720 scale_(1.0), 00721 generator_(generator) 00722 {} 00723 00724 /** Create functor for uniform random double-precision numbers in the range [<tt>lower</tt>, <tt>upper</tt>] 00725 using the given engine. 00726 00727 That is, the generated numbers satisfy <tt>lower</tt> <= i <= <tt>upper</tt>. 00728 */ 00729 UniformRandomFunctor(double lower, double upper, 00730 Engine & generator = Engine::global()) 00731 : offset_(lower), 00732 scale_(upper - lower), 00733 generator_(generator) 00734 { 00735 vigra_precondition(lower < upper, 00736 "UniformRandomFunctor(): lower bound must be smaller than upper bound."); 00737 } 00738 00739 /** Return a random number as specified in the constructor. 00740 */ 00741 double operator()() const 00742 { 00743 return generator_.uniform() * scale_ + offset_; 00744 } 00745 }; 00746 00747 template <class Engine> 00748 class FunctorTraits<UniformRandomFunctor<Engine> > 00749 { 00750 public: 00751 typedef UniformRandomFunctor<Engine> type; 00752 00753 typedef VigraTrueType isInitializer; 00754 00755 typedef VigraFalseType isUnaryFunctor; 00756 typedef VigraFalseType isBinaryFunctor; 00757 typedef VigraFalseType isTernaryFunctor; 00758 00759 typedef VigraFalseType isUnaryAnalyser; 00760 typedef VigraFalseType isBinaryAnalyser; 00761 typedef VigraFalseType isTernaryAnalyser; 00762 }; 00763 00764 /** Functor to create normal variate random numbers. 00765 00766 This functor encapsulates the function <tt>normal()</tt> of the given random number 00767 <tt>Engine</tt> (usually <tt>RandomTT800</tt> or <tt>RandomMT19937</tt>) 00768 in an STL-compatible interface. 00769 00770 <b>Traits defined:</b> 00771 00772 <tt>FunctorTraits<UniformIntRandomFunctor<Engine> >::isInitializer</tt> is true (<tt>VigraTrueType</tt>). 00773 */ 00774 template <class Engine = RandomTT800> 00775 class NormalRandomFunctor 00776 { 00777 double mean_, stddev_; 00778 Engine const & generator_; 00779 00780 public: 00781 00782 typedef double result_type; ///< STL required functor result type 00783 00784 /** Create functor for standard normal random numbers 00785 using the given engine. 00786 00787 That is, mean is 0.0 and standard deviation is 1.0. 00788 */ 00789 NormalRandomFunctor(Engine const & generator = Engine::global()) 00790 : mean_(0.0), 00791 stddev_(1.0), 00792 generator_(generator) 00793 {} 00794 00795 /** Create functor for normal random numbers with goven mean and standard deviation 00796 using the given engine. 00797 */ 00798 NormalRandomFunctor(double mean, double stddev, 00799 Engine & generator = Engine::global()) 00800 : mean_(mean), 00801 stddev_(stddev), 00802 generator_(generator) 00803 { 00804 vigra_precondition(stddev > 0.0, 00805 "NormalRandomFunctor(): standard deviation must be positive."); 00806 } 00807 00808 /** Return a random number as specified in the constructor. 00809 */ 00810 double operator()() const 00811 { 00812 return generator_.normal() * stddev_ + mean_; 00813 } 00814 00815 }; 00816 00817 template <class Engine> 00818 class FunctorTraits<NormalRandomFunctor<Engine> > 00819 { 00820 public: 00821 typedef UniformRandomFunctor<Engine> type; 00822 00823 typedef VigraTrueType isInitializer; 00824 00825 typedef VigraFalseType isUnaryFunctor; 00826 typedef VigraFalseType isBinaryFunctor; 00827 typedef VigraFalseType isTernaryFunctor; 00828 00829 typedef VigraFalseType isUnaryAnalyser; 00830 typedef VigraFalseType isBinaryAnalyser; 00831 typedef VigraFalseType isTernaryAnalyser; 00832 }; 00833 00834 //@} 00835 00836 } // namespace vigra 00837 00838 #endif // VIGRA_RANDOM_HXX
© Ullrich Köthe (ullrich.koethe@iwr.uni-heidelberg.de) |
html generated using doxygen and Python
|