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