[ VIGRA Homepage | Function Index | Class Index | Namespaces | File List | Main Page ]

random.hxx
1 /************************************************************************/
2 /* */
3 /* Copyright 2008 by Ullrich Koethe */
4 /* */
5 /* This file is part of the VIGRA computer vision library. */
6 /* The VIGRA Website is */
7 /* http://hci.iwr.uni-heidelberg.de/vigra/ */
8 /* Please direct questions, bug reports, and contributions to */
9 /* ullrich.koethe@iwr.uni-heidelberg.de or */
10 /* vigra@informatik.uni-hamburg.de */
11 /* */
12 /* Permission is hereby granted, free of charge, to any person */
13 /* obtaining a copy of this software and associated documentation */
14 /* files (the "Software"), to deal in the Software without */
15 /* restriction, including without limitation the rights to use, */
16 /* copy, modify, merge, publish, distribute, sublicense, and/or */
17 /* sell copies of the Software, and to permit persons to whom the */
18 /* Software is furnished to do so, subject to the following */
19 /* conditions: */
20 /* */
21 /* The above copyright notice and this permission notice shall be */
22 /* included in all copies or substantial portions of the */
23 /* Software. */
24 /* */
25 /* THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND */
26 /* EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES */
27 /* OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND */
28 /* NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT */
29 /* HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, */
30 /* WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING */
31 /* FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR */
32 /* OTHER DEALINGS IN THE SOFTWARE. */
33 /* */
34 /************************************************************************/
35 
36 
37 #ifndef VIGRA_RANDOM_HXX
38 #define VIGRA_RANDOM_HXX
39 
40 #include "mathutil.hxx"
41 #include "functortraits.hxx"
42 
43 #include <ctime>
44 
45 namespace vigra {
46 
47 enum RandomSeedTag { RandomSeed };
48 
49 namespace detail {
50 
51 enum RandomEngineTag { TT800, MT19937 };
52 
53 template<RandomEngineTag EngineTag>
54 struct RandomState;
55 
56 template <RandomEngineTag EngineTag>
57 void seed(UInt32 theSeed, RandomState<EngineTag> & engine)
58 {
59  engine.state_[0] = theSeed;
60  for(UInt32 i=1; i<RandomState<EngineTag>::N; ++i)
61  {
62  engine.state_[i] = 1812433253U * (engine.state_[i-1] ^ (engine.state_[i-1] >> 30)) + i;
63  }
64 }
65 
66 template <class Iterator, RandomEngineTag EngineTag>
67 void seed(Iterator init, UInt32 key_length, RandomState<EngineTag> & engine)
68 {
69  const UInt32 N = RandomState<EngineTag>::N;
70  int k = (int)std::max(N, key_length);
71  UInt32 i = 1, j = 0;
72  Iterator data = init;
73  for (; k; --k)
74  {
75  engine.state_[i] = (engine.state_[i] ^ ((engine.state_[i-1] ^ (engine.state_[i-1] >> 30)) * 1664525U))
76  + *data + j; /* non linear */
77  ++i; ++j; ++data;
78 
79  if (i >= N)
80  {
81  engine.state_[0] = engine.state_[N-1];
82  i=1;
83  }
84  if (j>=key_length)
85  {
86  j=0;
87  data = init;
88  }
89  }
90 
91  for (k=N-1; k; --k)
92  {
93  engine.state_[i] = (engine.state_[i] ^ ((engine.state_[i-1] ^ (engine.state_[i-1] >> 30)) * 1566083941U))
94  - i; /* non linear */
95  ++i;
96  if (i>=N)
97  {
98  engine.state_[0] = engine.state_[N-1];
99  i=1;
100  }
101  }
102 
103  engine.state_[0] = 0x80000000U; /* MSB is 1; assuring non-zero initial array */
104 }
105 
106 template <RandomEngineTag EngineTag>
107 void seed(RandomSeedTag, RandomState<EngineTag> & engine)
108 {
109  static UInt32 globalCount = 0;
110  UInt32 init[3] = { (UInt32)time(0), (UInt32)clock(), ++globalCount };
111  seed(init, 3, engine);
112 }
113 
114 
115  /* Tempered twister TT800 by M. Matsumoto */
116 template<>
117 struct RandomState<TT800>
118 {
119  static const UInt32 N = 25, M = 7;
120 
121  mutable UInt32 state_[N];
122  mutable UInt32 current_;
123 
124  RandomState()
125  : current_(0)
126  {
127  UInt32 seeds[N] = {
128  0x95f24dab, 0x0b685215, 0xe76ccae7, 0xaf3ec239, 0x715fad23,
129  0x24a590ad, 0x69e4b5ef, 0xbf456141, 0x96bc1b7b, 0xa7bdf825,
130  0xc1de75b7, 0x8858a9c9, 0x2da87693, 0xb657f9dd, 0xffdc8a9f,
131  0x8121da71, 0x8b823ecb, 0x885d05f5, 0x4e20cd47, 0x5a9ad5d9,
132  0x512c0c03, 0xea857ccd, 0x4cc1d30f, 0x8891a8a1, 0xa6b7aadb
133  };
134 
135  for(UInt32 i=0; i<N; ++i)
136  state_[i] = seeds[i];
137  }
138 
139  protected:
140 
141  UInt32 get() const
142  {
143  if(current_ == N)
144  generateNumbers<void>();
145 
146  UInt32 y = state_[current_++];
147  y ^= (y << 7) & 0x2b5b2500;
148  y ^= (y << 15) & 0xdb8b0000;
149  return y ^ (y >> 16);
150  }
151 
152  template <class DUMMY>
153  void generateNumbers() const;
154 
155  void seedImpl(RandomSeedTag)
156  {
157  seed(RandomSeed, *this);
158  }
159 
160  void seedImpl(UInt32 theSeed)
161  {
162  seed(theSeed, *this);
163  }
164 
165  template<class Iterator>
166  void seedImpl(Iterator init, UInt32 length)
167  {
168  seed(init, length, *this);
169  }
170 };
171 
172 template <class DUMMY>
173 void RandomState<TT800>::generateNumbers() const
174 {
175  UInt32 mag01[2]= { 0x0, 0x8ebfd028 };
176 
177  for(UInt32 i=0; i<N-M; ++i)
178  {
179  state_[i] = state_[i+M] ^ (state_[i] >> 1) ^ mag01[state_[i] % 2];
180  }
181  for (UInt32 i=N-M; i<N; ++i)
182  {
183  state_[i] = state_[i+(M-N)] ^ (state_[i] >> 1) ^ mag01[state_[i] % 2];
184  }
185  current_ = 0;
186 }
187 
188  /* Mersenne twister MT19937 by M. Matsumoto */
189 template<>
190 struct RandomState<MT19937>
191 {
192  static const UInt32 N = 624, M = 397;
193 
194  mutable UInt32 state_[N];
195  mutable UInt32 current_;
196 
197  RandomState()
198  : current_(0)
199  {
200  seed(19650218U, *this);
201  }
202 
203  protected:
204 
205  UInt32 get() const
206  {
207  if(current_ == N)
208  generateNumbers<void>();
209 
210  UInt32 x = state_[current_++];
211  x ^= (x >> 11);
212  x ^= (x << 7) & 0x9D2C5680U;
213  x ^= (x << 15) & 0xEFC60000U;
214  return x ^ (x >> 18);
215  }
216 
217  template <class DUMMY>
218  void generateNumbers() const;
219 
220  static UInt32 twiddle(UInt32 u, UInt32 v)
221  {
222  return (((u & 0x80000000U) | (v & 0x7FFFFFFFU)) >> 1)
223  ^ ((v & 1U) ? 0x9908B0DFU : 0x0U);
224  }
225 
226  void seedImpl(RandomSeedTag)
227  {
228  seed(RandomSeed, *this);
229  generateNumbers<void>();
230  }
231 
232  void seedImpl(UInt32 theSeed)
233  {
234  seed(theSeed, *this);
235  generateNumbers<void>();
236  }
237 
238  template<class Iterator>
239  void seedImpl(Iterator init, UInt32 length)
240  {
241  seed(19650218U, *this);
242  seed(init, length, *this);
243  generateNumbers<void>();
244  }
245 };
246 
247 template <class DUMMY>
248 void RandomState<MT19937>::generateNumbers() const
249 {
250  for (unsigned int i = 0; i < (N - M); ++i)
251  {
252  state_[i] = state_[i + M] ^ twiddle(state_[i], state_[i + 1]);
253  }
254  for (unsigned int i = N - M; i < (N - 1); ++i)
255  {
256  state_[i] = state_[i + M - N] ^ twiddle(state_[i], state_[i + 1]);
257  }
258  state_[N - 1] = state_[M - 1] ^ twiddle(state_[N - 1], state_[0]);
259  current_ = 0;
260 }
261 
262 } // namespace detail
263 
264 
265 /** \addtogroup RandomNumberGeneration Random Number Generation
266 
267  High-quality random number generators and related functors.
268 */
269 //@{
270 
271 /** Generic random number generator.
272 
273  The actual generator is passed in the template argument <tt>Engine</tt>. Two generators
274  are currently available:
275  <ul>
276  <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.
277  <li> <tt>RandomTT800</tt>: (default) The Tempered Twister, a simpler predecessor of the Mersenne Twister with period length 2<sup>800</sup>.
278  </ul>
279 
280  Both generators have been designed by <a href="http://www.math.sci.hiroshima-u.ac.jp/~m-mat/eindex.html">Makoto Matsumoto</a>.
281 
282  <b>Traits defined:</b>
283 
284  \verbatim FunctorTraits<RandomNumberGenerator<Engine> >::isInitializer \endverbatim
285  is true (<tt>VigraTrueType</tt>).
286 */
287 template <class Engine = detail::RandomState<detail::TT800> >
289 : public Engine
290 {
291  mutable double normalCached_;
292  mutable bool normalCachedValid_;
293 
294  public:
295 
296  /** Create a new random generator object with standard seed.
297 
298  Due to standard seeding, the random numbers generated will always be the same.
299  This is useful for debugging.
300  */
302  : normalCached_(0.0),
303  normalCachedValid_(false)
304  {}
305 
306  /** Create a new random generator object with a random seed.
307 
308  The seed is obtained from the machines current <tt>clock()</tt> and <tt>time()</tt>
309  values.
310 
311  <b>Usage:</b>
312  \code
313  RandomNumberGenerator<> rnd = RandomNumberGenerator<>(RandomSeed);
314  \endcode
315  */
316  RandomNumberGenerator(RandomSeedTag)
317  : normalCached_(0.0),
318  normalCachedValid_(false)
319  {
320  this->seedImpl(RandomSeed);
321  }
322 
323  /** Create a new random generator object from the given seed.
324 
325  The same seed will always produce identical random sequences.
326  */
328  : normalCached_(0.0),
329  normalCachedValid_(false)
330  {
331  this->seedImpl(theSeed);
332  }
333 
334  /** Create a new random generator object from the given seed sequence.
335 
336  Longer seed sequences lead to better initialization in the sense that the generator's
337  state space is covered much better than is possible with 32-bit seeds alone.
338  */
339  template<class Iterator>
340  RandomNumberGenerator(Iterator init, UInt32 length)
341  : normalCached_(0.0),
342  normalCachedValid_(false)
343  {
344  this->seedImpl(init, length);
345  }
346 
347 
348  /** Re-initialize the random generator object with a random seed.
349 
350  The seed is obtained from the machines current <tt>clock()</tt> and <tt>time()</tt>
351  values.
352 
353  <b>Usage:</b>
354  \code
355  RandomNumberGenerator<> rnd = ...;
356  ...
357  rnd.seed(RandomSeed);
358  \endcode
359  */
360  void seed(RandomSeedTag)
361  {
362  this->seedImpl(RandomSeed);
363  normalCachedValid_ = false;
364  }
365 
366  /** Re-initialize the random generator object from the given seed.
367 
368  The same seed will always produce identical random sequences.
369  */
370  void seed(UInt32 theSeed)
371  {
372  this->seedImpl(theSeed);
373  normalCachedValid_ = false;
374  }
375 
376  /** Re-initialize the random generator object from the given seed sequence.
377 
378  Longer seed sequences lead to better initialization in the sense that the generator's
379  state space is covered much better than is possible with 32-bit seeds alone.
380  */
381  template<class Iterator>
382  void seed(Iterator init, UInt32 length)
383  {
384  this->seedImpl(init, length);
385  normalCachedValid_ = false;
386  }
387 
388  /** Return a uniformly distributed integer random number in [0, 2<sup>32</sup>).
389 
390  That is, 0 &lt;= i &lt; 2<sup>32</sup>.
391  */
393  {
394  return this->get();
395  }
396 
397  /** Return a uniformly distributed integer random number in [0, 2<sup>32</sup>).
398 
399  That is, 0 &lt;= i &lt; 2<sup>32</sup>.
400  */
402  {
403  return this->get();
404  }
405 
406 
407 #if 0 // difficult implementation necessary if low bits are not sufficiently random
408  // in [0,beyond)
409  UInt32 uniformInt(UInt32 beyond) const
410  {
411  if(beyond < 2)
412  return 0;
413 
414  UInt32 factor = factorForUniformInt(beyond);
415  UInt32 res = this->get() / factor;
416 
417  // Use rejection method to avoid quantization bias.
418  // On average, we will need two raw random numbers to generate one.
419  while(res >= beyond)
420  res = this->get() / factor;
421  return res;
422  }
423 #endif /* #if 0 */
424 
425  /** Return a uniformly distributed integer random number in [0, <tt>beyond</tt>).
426 
427  That is, 0 &lt;= i &lt; <tt>beyond</tt>.
428  */
429  UInt32 uniformInt(UInt32 beyond) const
430  {
431  // in [0,beyond) -- simple implementation when low bits are sufficiently random, which is
432  // the case for TT800 and MT19937
433  if(beyond < 2)
434  return 0;
435 
436  UInt32 remainder = (NumericTraits<UInt32>::max() - beyond + 1) % beyond;
437  UInt32 lastSafeValue = NumericTraits<UInt32>::max() - remainder;
438  UInt32 res = this->get();
439 
440  // Use rejection method to avoid quantization bias.
441  // We will need two raw random numbers in amortized worst case.
442  while(res > lastSafeValue)
443  res = this->get();
444  return res % beyond;
445  }
446 
447  /** Return a uniformly distributed double-precision random number in [0.0, 1.0).
448 
449  That is, 0.0 &lt;= i &lt; 1.0. All 53-bit bits of the mantissa are random (two 32-bit integers are used to
450  create this number).
451  */
452  double uniform53() const
453  {
454  // make full use of the entire 53-bit mantissa of a double, by Isaku Wada
455  return ( (this->get() >> 5) * 67108864.0 + (this->get() >> 6)) * (1.0/9007199254740992.0);
456  }
457 
458  /** Return a uniformly distributed double-precision random number in [0.0, 1.0].
459 
460  That is, 0.0 &lt;= i &lt;= 1.0. This number is computed by <tt>uniformInt()</tt> / 2<sup>32</sup>,
461  so it has effectively only 32 random bits.
462  */
463  double uniform() const
464  {
465  return (double)this->get() / 4294967295.0;
466  }
467 
468  /** Return a uniformly distributed double-precision random number in [lower, upper].
469 
470  That is, <tt>lower</tt> &lt;= i &lt;= <tt>upper</tt>. This number is computed
471  from <tt>uniform()</tt>, so it has effectively only 32 random bits.
472  */
473  double uniform(double lower, double upper) const
474  {
475  vigra_precondition(lower < upper,
476  "RandomNumberGenerator::uniform(): lower bound must be smaller than upper bound.");
477  return uniform() * (upper-lower) + lower;
478  }
479 
480  /** Return a standard normal variate (Gaussian) random number.
481 
482  Mean is zero, standard deviation is 1.0. It uses the polar form of the
483  Box-Muller transform.
484  */
485  double normal() const;
486 
487  /** Return a normal variate (Gaussian) random number with the given mean and standard deviation.
488 
489  It uses the polar form of the Box-Muller transform.
490  */
491  double normal(double mean, double stddev) const
492  {
493  vigra_precondition(stddev > 0.0,
494  "RandomNumberGenerator::normal(): standard deviation must be positive.");
495  return normal()*stddev + mean;
496  }
497 
498  /** Access the global (program-wide) instance of the present random number generator.
499 
500  Normally, you will create a local generator by one of the constructor calls. But sometimes
501  it is useful to have all program parts access the same generator.
502  */
504  {
505  static RandomNumberGenerator generator;
506  return generator;
507  }
508 
509  static UInt32 factorForUniformInt(UInt32 range)
510  {
511  return (range > 2147483648U || range == 0)
512  ? 1
513  : 2*(2147483648U / ceilPower2(range));
514  }
515 };
516 
517 template <class Engine>
519 {
520  if(normalCachedValid_)
521  {
522  normalCachedValid_ = false;
523  return normalCached_;
524  }
525  else
526  {
527  double x1, x2, w;
528  do
529  {
530  x1 = uniform(-1.0, 1.0);
531  x2 = uniform(-1.0, 1.0);
532  w = x1 * x1 + x2 * x2;
533  }
534  while ( w > 1.0 || w == 0.0);
535 
536  w = std::sqrt( -2.0 * std::log( w ) / w );
537 
538  normalCached_ = x2 * w;
539  normalCachedValid_ = true;
540 
541  return x1 * w;
542  }
543 }
544 
545  /** Shorthand for the TT800 random number generator class.
546  */
548 
549  /** Shorthand for the TT800 random number generator class (same as RandomTT800).
550  */
552 
553  /** Shorthand for the MT19937 random number generator class.
554  */
556 
557  /** Shorthand for the MT19937 random number generator class (same as RandomMT19937).
558  */
560 
561  /** Access the global (program-wide) instance of the TT800 random number generator.
562  */
564 
565  /** Access the global (program-wide) instance of the MT19937 random number generator.
566  */
568 
569 template <class Engine>
570 class FunctorTraits<RandomNumberGenerator<Engine> >
571 {
572  public:
573  typedef RandomNumberGenerator<Engine> type;
574 
575  typedef VigraTrueType isInitializer;
576 
577  typedef VigraFalseType isUnaryFunctor;
578  typedef VigraFalseType isBinaryFunctor;
579  typedef VigraFalseType isTernaryFunctor;
580 
581  typedef VigraFalseType isUnaryAnalyser;
582  typedef VigraFalseType isBinaryAnalyser;
583  typedef VigraFalseType isTernaryAnalyser;
584 };
585 
586 
587 /** Functor to create uniformly distributed integer random numbers.
588 
589  This functor encapsulates the appropriate functions of the given random number
590  <tt>Engine</tt> (usually <tt>RandomTT800</tt> or <tt>RandomMT19937</tt>)
591  in an STL-compatible interface.
592 
593  <b>Traits defined:</b>
594 
595  \verbatim FunctorTraits<UniformIntRandomFunctor<Engine> >::isInitializer \endverbatim
596  and
597  \verbatim FunctorTraits<UniformIntRandomFunctor<Engine> >::isUnaryFunctor \endverbatim
598  are true (<tt>VigraTrueType</tt>).
599 */
600 template <class Engine = RandomTT800>
602 {
603  UInt32 lower_, difference_, factor_;
604  Engine const & generator_;
605  bool useLowBits_;
606 
607  public:
608 
609  typedef UInt32 argument_type; ///< STL required functor argument type
610  typedef UInt32 result_type; ///< STL required functor result type
611 
612  /** Create functor for uniform random integers in the range [0, 2<sup>32</sup>)
613  using the given engine.
614 
615  That is, the generated numbers satisfy 0 &lt;= i &lt; 2<sup>32</sup>.
616  */
617  explicit UniformIntRandomFunctor(Engine const & generator = Engine::global() )
618  : lower_(0), difference_(0xffffffff), factor_(1),
619  generator_(generator),
620  useLowBits_(true)
621  {}
622 
623  /** Create functor for uniform random integers in the range [<tt>lower</tt>, <tt>upper</tt>]
624  using the given engine.
625 
626  That is, the generated numbers satisfy <tt>lower</tt> &lt;= i &lt;= <tt>upper</tt>.
627  \a useLowBits should be set to <tt>false</tt> when the engine generates
628  random numbers whose low bits are significantly less random than the high
629  bits. This does not apply to <tt>RandomTT800</tt> and <tt>RandomMT19937</tt>,
630  but is necessary for simpler linear congruential generators.
631  */
633  Engine const & generator = Engine::global(),
634  bool useLowBits = true)
635  : lower_(lower), difference_(upper-lower),
636  factor_(Engine::factorForUniformInt(difference_ + 1)),
637  generator_(generator),
638  useLowBits_(useLowBits)
639  {
640  vigra_precondition(lower < upper,
641  "UniformIntRandomFunctor(): lower bound must be smaller than upper bound.");
642  }
643 
644  /** Return a random number as specified in the constructor.
645  */
647  {
648  if(difference_ == 0xffffffff) // lower_ is necessarily 0
649  return generator_();
650  else if(useLowBits_)
651  return generator_.uniformInt(difference_+1) + lower_;
652  else
653  {
654  UInt32 res = generator_() / factor_;
655 
656  // Use rejection method to avoid quantization bias.
657  // On average, we will need two raw random numbers to generate one.
658  while(res > difference_)
659  res = generator_() / factor_;
660  return res + lower_;
661  }
662  }
663 
664  /** Return a uniformly distributed integer random number in the range [0, <tt>beyond</tt>).
665 
666  That is, 0 &lt;= i &lt; <tt>beyond</tt>. This is a required interface for
667  <tt>std::random_shuffle</tt>. It ignores the limits specified
668  in the constructor and the flag <tt>useLowBits</tt>.
669  */
670  UInt32 operator()(UInt32 beyond) const
671  {
672  if(beyond < 2)
673  return 0;
674 
675  return generator_.uniformInt(beyond);
676  }
677 };
678 
679 template <class Engine>
680 class FunctorTraits<UniformIntRandomFunctor<Engine> >
681 {
682  public:
683  typedef UniformIntRandomFunctor<Engine> type;
684 
685  typedef VigraTrueType isInitializer;
686 
687  typedef VigraTrueType isUnaryFunctor;
688  typedef VigraFalseType isBinaryFunctor;
689  typedef VigraFalseType isTernaryFunctor;
690 
691  typedef VigraFalseType isUnaryAnalyser;
692  typedef VigraFalseType isBinaryAnalyser;
693  typedef VigraFalseType isTernaryAnalyser;
694 };
695 
696 /** Functor to create uniformly distributed double-precision random numbers.
697 
698  This functor encapsulates the function <tt>uniform()</tt> of the given random number
699  <tt>Engine</tt> (usually <tt>RandomTT800</tt> or <tt>RandomMT19937</tt>)
700  in an STL-compatible interface.
701 
702  <b>Traits defined:</b>
703 
704  \verbatim FunctorTraits<UniformIntRandomFunctor<Engine> >::isInitializer \endverbatim
705  is true (<tt>VigraTrueType</tt>).
706 */
707 template <class Engine = RandomTT800>
709 {
710  double offset_, scale_;
711  Engine const & generator_;
712 
713  public:
714 
715  typedef double result_type; ///< STL required functor result type
716 
717  /** Create functor for uniform random double-precision numbers in the range [0.0, 1.0]
718  using the given engine.
719 
720  That is, the generated numbers satisfy 0.0 &lt;= i &lt;= 1.0.
721  */
722  UniformRandomFunctor(Engine const & generator = Engine::global())
723  : offset_(0.0),
724  scale_(1.0),
725  generator_(generator)
726  {}
727 
728  /** Create functor for uniform random double-precision numbers in the range [<tt>lower</tt>, <tt>upper</tt>]
729  using the given engine.
730 
731  That is, the generated numbers satisfy <tt>lower</tt> &lt;= i &lt;= <tt>upper</tt>.
732  */
733  UniformRandomFunctor(double lower, double upper,
734  Engine & generator = Engine::global())
735  : offset_(lower),
736  scale_(upper - lower),
737  generator_(generator)
738  {
739  vigra_precondition(lower < upper,
740  "UniformRandomFunctor(): lower bound must be smaller than upper bound.");
741  }
742 
743  /** Return a random number as specified in the constructor.
744  */
745  double operator()() const
746  {
747  return generator_.uniform() * scale_ + offset_;
748  }
749 };
750 
751 template <class Engine>
752 class FunctorTraits<UniformRandomFunctor<Engine> >
753 {
754  public:
755  typedef UniformRandomFunctor<Engine> type;
756 
757  typedef VigraTrueType isInitializer;
758 
759  typedef VigraFalseType isUnaryFunctor;
760  typedef VigraFalseType isBinaryFunctor;
761  typedef VigraFalseType isTernaryFunctor;
762 
763  typedef VigraFalseType isUnaryAnalyser;
764  typedef VigraFalseType isBinaryAnalyser;
765  typedef VigraFalseType isTernaryAnalyser;
766 };
767 
768 /** Functor to create normal variate random numbers.
769 
770  This functor encapsulates the function <tt>normal()</tt> of the given random number
771  <tt>Engine</tt> (usually <tt>RandomTT800</tt> or <tt>RandomMT19937</tt>)
772  in an STL-compatible interface.
773 
774  <b>Traits defined:</b>
775 
776  \verbatim FunctorTraits<UniformIntRandomFunctor<Engine> >::isInitializer \endverbatim
777  is true (<tt>VigraTrueType</tt>).
778 */
779 template <class Engine = RandomTT800>
781 {
782  double mean_, stddev_;
783  Engine const & generator_;
784 
785  public:
786 
787  typedef double result_type; ///< STL required functor result type
788 
789  /** Create functor for standard normal random numbers
790  using the given engine.
791 
792  That is, mean is 0.0 and standard deviation is 1.0.
793  */
794  NormalRandomFunctor(Engine const & generator = Engine::global())
795  : mean_(0.0),
796  stddev_(1.0),
797  generator_(generator)
798  {}
799 
800  /** Create functor for normal random numbers with given mean and standard deviation
801  using the given engine.
802  */
803  NormalRandomFunctor(double mean, double stddev,
804  Engine & generator = Engine::global())
805  : mean_(mean),
806  stddev_(stddev),
807  generator_(generator)
808  {
809  vigra_precondition(stddev > 0.0,
810  "NormalRandomFunctor(): standard deviation must be positive.");
811  }
812 
813  /** Return a random number as specified in the constructor.
814  */
815  double operator()() const
816  {
817  return generator_.normal() * stddev_ + mean_;
818  }
819 
820 };
821 
822 template <class Engine>
823 class FunctorTraits<NormalRandomFunctor<Engine> >
824 {
825  public:
826  typedef UniformRandomFunctor<Engine> type;
827 
828  typedef VigraTrueType isInitializer;
829 
830  typedef VigraFalseType isUnaryFunctor;
831  typedef VigraFalseType isBinaryFunctor;
832  typedef VigraFalseType isTernaryFunctor;
833 
834  typedef VigraFalseType isUnaryAnalyser;
835  typedef VigraFalseType isBinaryAnalyser;
836  typedef VigraFalseType isTernaryAnalyser;
837 };
838 
839 //@}
840 
841 } // namespace vigra
842 
843 #endif // VIGRA_RANDOM_HXX

© Ullrich Köthe (ullrich.koethe@iwr.uni-heidelberg.de)
Heidelberg Collaboratory for Image Processing, University of Heidelberg, Germany

html generated using doxygen and Python
vigra 1.8.0 (Wed Sep 26 2012)