IT++ Logo
random.h
Go to the documentation of this file.
00001 
00029 #ifndef RANDOM_H
00030 #define RANDOM_H
00031 
00032 #include <itpp/base/random_dsfmt.h>
00033 #include <itpp/base/operators.h>
00034 
00035 
00036 namespace itpp
00037 {
00038 
00040 
00050 typedef DSFMT_19937_RNG Random_Generator;
00051 
00054 
00056 void RNG_reset(unsigned int seed);
00058 void RNG_reset();
00060 void RNG_randomize();
00062 void RNG_get_state(ivec &state);
00064 void RNG_set_state(const ivec &state);
00066 
00071 class Bernoulli_RNG
00072 {
00073 public:
00075   Bernoulli_RNG(double prob) { setup(prob); }
00077   Bernoulli_RNG() { p = 0.5; }
00079   void setup(double prob) {
00080     it_assert(prob >= 0.0 && prob <= 1.0, "The Bernoulli source probability "
00081               "must be between 0 and 1");
00082     p = prob;
00083   }
00085   double get_setup() const { return p; }
00087   bin operator()() { return sample(); }
00089   bvec operator()(int n) { bvec temp(n); sample_vector(n, temp); return temp; }
00091   bmat operator()(int h, int w) { bmat temp(h, w); sample_matrix(h, w, temp); return temp; }
00093   bin sample() { return RNG.genrand_close_open() < p ? bin(1) : bin(0); }
00095   void sample_vector(int size, bvec &out) {
00096     out.set_size(size, false);
00097     for (int i = 0; i < size; i++) out(i) = sample();
00098   }
00100   void sample_matrix(int rows, int cols, bmat &out) {
00101     out.set_size(rows, cols, false);
00102     for (int i = 0; i < rows*cols; i++) out(i) = sample();
00103   }
00104 protected:
00105 private:
00107   double p;
00109   Random_Generator RNG;
00110 };
00111 
00129 class I_Uniform_RNG
00130 {
00131 public:
00133   I_Uniform_RNG(int min = 0, int max = 1);
00135   void setup(int min, int max);
00137   void get_setup(int &min, int &max) const;
00139   int operator()() { return sample(); }
00141   ivec operator()(int n);
00143   imat operator()(int h, int w);
00145   int sample() {
00146     return floor_i(RNG.genrand_close_open() * (hi - lo + 1)) + lo;
00147   }
00148 private:
00150   int lo;
00152   int hi;
00154   Random_Generator RNG;
00155 };
00156 
00161 class Uniform_RNG
00162 {
00163 public:
00165   Uniform_RNG(double min = 0, double max = 1.0);
00167   void setup(double min, double max);
00169   void get_setup(double &min, double &max) const;
00171   double operator()() { return (sample() * (hi_bound - lo_bound) + lo_bound); }
00173   vec operator()(int n) {
00174     vec temp(n);
00175     sample_vector(n, temp);
00176     temp *= hi_bound - lo_bound;
00177     temp += lo_bound;
00178     return temp;
00179   }
00181   mat operator()(int h, int w) {
00182     mat temp(h, w);
00183     sample_matrix(h, w, temp);
00184     temp *= hi_bound - lo_bound;
00185     temp += lo_bound;
00186     return temp;
00187   }
00189   double sample() {  return RNG.genrand_close_open(); }
00191   void sample_vector(int size, vec &out) {
00192     out.set_size(size, false);
00193     for (int i = 0; i < size; i++) out(i) = sample();
00194   }
00196   void sample_matrix(int rows, int cols, mat &out) {
00197     out.set_size(rows, cols, false);
00198     for (int i = 0; i < rows*cols; i++) out(i) = sample();
00199   }
00200 protected:
00201 private:
00203   double lo_bound, hi_bound;
00205   Random_Generator RNG;
00206 };
00207 
00212 class Exponential_RNG
00213 {
00214 public:
00216   Exponential_RNG(double lambda = 1.0);
00218   void setup(double lambda) { l = lambda; }
00220   double get_setup() const;
00222   double operator()() { return sample(); }
00224   vec operator()(int n);
00226   mat operator()(int h, int w);
00227 private:
00228   double sample() { return (-std::log(RNG.genrand_open_close()) / l); }
00229   double l;
00230   Random_Generator RNG;
00231 };
00232 
00247 class Normal_RNG
00248 {
00249 public:
00251   Normal_RNG(double meanval, double variance):
00252       mean(meanval), sigma(std::sqrt(variance)) {}
00254   Normal_RNG(): mean(0.0), sigma(1.0) {}
00256   void setup(double meanval, double variance)
00257   { mean = meanval; sigma = std::sqrt(variance); }
00259   void get_setup(double &meanval, double &variance) const;
00261   double operator()() { return (sigma*sample() + mean); }
00263   vec operator()(int n) {
00264     vec temp(n);
00265     sample_vector(n, temp);
00266     temp *= sigma;
00267     temp += mean;
00268     return temp;
00269   }
00271   mat operator()(int h, int w) {
00272     mat temp(h, w);
00273     sample_matrix(h, w, temp);
00274     temp *= sigma;
00275     temp += mean;
00276     return temp;
00277   }
00279   double sample();
00280 
00282   void sample_vector(int size, vec &out) {
00283     out.set_size(size, false);
00284     for (int i = 0; i < size; i++) out(i) = sample();
00285   }
00286 
00288   void sample_matrix(int rows, int cols, mat &out) {
00289     out.set_size(rows, cols, false);
00290     for (int i = 0; i < rows*cols; i++) out(i) = sample();
00291   }
00292 private:
00293   double mean, sigma;
00294   static const double ytab[128];
00295   static const unsigned int ktab[128];
00296   static const double wtab[128];
00297   static const double PARAM_R;
00298   Random_Generator RNG;
00299 };
00300 
00317 class Gamma_RNG
00318 {
00319 public:
00321   Gamma_RNG(double a = 1.0, double b = 1.0): alpha(a), beta(b) {}
00323   void setup(double a, double b) { alpha = a; beta = b; }
00325   double operator()() { return sample(); }
00327   vec operator()(int n);
00329   mat operator()(int r, int c);
00331   double sample();
00332 private:
00333   double alpha;
00334   double beta;
00335   Random_Generator RNG;
00336   Normal_RNG NRNG;
00337 };
00338 
00343 class Laplace_RNG
00344 {
00345 public:
00347   Laplace_RNG(double meanval = 0.0, double variance = 1.0);
00349   void setup(double meanval, double variance);
00351   void get_setup(double &meanval, double &variance) const;
00353   double operator()() { return sample(); }
00355   vec operator()(int n);
00357   mat operator()(int h, int w);
00359   double sample() {
00360     double u = RNG.genrand_open_open();
00361     double l = sqrt_12var;
00362     if (u < 0.5)
00363       l *= std::log(2.0 * u);
00364     else
00365       l *= -std::log(2.0 * (1 - u));
00366     return (l + mean);
00367   }
00368 private:
00369   double mean, var, sqrt_12var;
00370   Random_Generator RNG;
00371 };
00372 
00377 class Complex_Normal_RNG
00378 {
00379 public:
00381   Complex_Normal_RNG(std::complex<double> mean, double variance):
00382       norm_factor(1.0 / std::sqrt(2.0)) {
00383     setup(mean, variance);
00384   }
00386   Complex_Normal_RNG(): m(0.0), sigma(1.0), norm_factor(1.0 / std::sqrt(2.0)) {}
00388   void setup(std::complex<double> mean, double variance) {
00389     m = mean;
00390     sigma = std::sqrt(variance);
00391   }
00393   void get_setup(std::complex<double> &mean, double &variance) {
00394     mean = m;
00395     variance = sigma * sigma;
00396   }
00398   std::complex<double> operator()() { return sigma*sample() + m; }
00400   cvec operator()(int n) {
00401     cvec temp(n);
00402     sample_vector(n, temp);
00403     temp *= sigma;
00404     temp += m;
00405     return temp;
00406   }
00408   cmat operator()(int h, int w) {
00409     cmat temp(h, w);
00410     sample_matrix(h, w, temp);
00411     temp *= sigma;
00412     temp += m;
00413     return temp;
00414   }
00416   std::complex<double> sample() {
00417     double a = nRNG.sample() * norm_factor;
00418     double b = nRNG.sample() * norm_factor;
00419     return std::complex<double>(a, b);
00420   }
00421 
00423   void sample_vector(int size, cvec &out) {
00424     out.set_size(size, false);
00425     for (int i = 0; i < size; i++) out(i) = sample();
00426   }
00427 
00429   void sample_matrix(int rows, int cols, cmat &out) {
00430     out.set_size(rows, cols, false);
00431     for (int i = 0; i < rows*cols; i++) out(i) = sample();
00432   }
00433 
00435   Complex_Normal_RNG & operator=(const Complex_Normal_RNG&) { return *this; }
00436 
00437 private:
00438   std::complex<double> m;
00439   double sigma;
00440   const double norm_factor;
00441   Normal_RNG nRNG;
00442 };
00443 
00448 class AR1_Normal_RNG
00449 {
00450 public:
00452   AR1_Normal_RNG(double meanval = 0.0, double variance = 1.0,
00453                  double rho = 0.0);
00455   void setup(double meanval, double variance, double rho);
00457   void get_setup(double &meanval, double &variance, double &rho) const;
00459   void reset();
00461   double operator()() { return sample(); }
00463   vec operator()(int n);
00465   mat operator()(int h, int w);
00466 private:
00467   double sample() {
00468     mem *= r;
00469     if (odd) {
00470       r1 = m_2pi * RNG.genrand_open_close();
00471       r2 = std::sqrt(factr * std::log(RNG.genrand_open_close()));
00472       mem += r2 * std::cos(r1);
00473     }
00474     else {
00475       mem += r2 * std::sin(r1);
00476     }
00477     odd = !odd;
00478     return (mem + mean);
00479   }
00480   double mem, r, factr, mean, var, r1, r2;
00481   bool odd;
00482   Random_Generator RNG;
00483 };
00484 
00489 typedef Normal_RNG Gauss_RNG;
00490 
00495 typedef AR1_Normal_RNG AR1_Gauss_RNG;
00496 
00501 class Weibull_RNG
00502 {
00503 public:
00505   Weibull_RNG(double lambda = 1.0, double beta = 1.0);
00507   void setup(double lambda, double beta);
00509   void get_setup(double &lambda, double &beta) { lambda = l; beta = b; }
00511   double operator()() { return sample(); }
00513   vec operator()(int n);
00515   mat operator()(int h, int w);
00516 private:
00517   double sample() {
00518     return (std::pow(-std::log(RNG.genrand_open_close()), 1.0 / b) / l);
00519   }
00520   double l, b;
00521   double mean, var;
00522   Random_Generator RNG;
00523 };
00524 
00529 class Rayleigh_RNG
00530 {
00531 public:
00533   Rayleigh_RNG(double sigma = 1.0);
00535   void setup(double sigma) { sig = sigma; }
00537   double get_setup() { return sig; }
00539   double operator()() { return sample(); }
00541   vec operator()(int n);
00543   mat operator()(int h, int w);
00544 private:
00545   double sample() {
00546     double s1 = nRNG.sample();
00547     double s2 = nRNG.sample();
00548     // s1 and s2 are N(0,1) and independent
00549     return (sig * std::sqrt(s1*s1 + s2*s2));
00550   }
00551   double sig;
00552   Normal_RNG nRNG;
00553 };
00554 
00559 class Rice_RNG
00560 {
00561 public:
00563   Rice_RNG(double sigma = 1.0, double v = 1.0);
00565   void setup(double sigma, double v) { sig = sigma; s = v; }
00567   void get_setup(double &sigma, double &v) { sigma = sig; v = s; }
00569   double operator()() { return sample(); }
00571   vec operator()(int n);
00573   mat operator()(int h, int w);
00574 private:
00575   double sample() {
00576     double s1 = nRNG.sample() + s;
00577     double s2 = nRNG.sample();
00578     // s1 and s2 are N(0,1) and independent
00579     return (sig * std::sqrt(s1*s1 + s2*s2));
00580   }
00581   double sig, s;
00582   Normal_RNG nRNG;
00583 };
00584 
00587 
00589 inline bin randb(void) { Bernoulli_RNG src; return src.sample(); }
00591 inline void randb(int size, bvec &out) { Bernoulli_RNG src; src.sample_vector(size, out); }
00593 inline bvec randb(int size) { bvec temp; randb(size, temp); return temp; }
00595 inline void randb(int rows, int cols, bmat &out) { Bernoulli_RNG src; src.sample_matrix(rows, cols, out); }
00597 inline bmat randb(int rows, int cols) { bmat temp; randb(rows, cols, temp); return temp; }
00598 
00600 inline double randu(void) { Uniform_RNG src; return src.sample(); }
00602 inline void randu(int size, vec &out) { Uniform_RNG src; src.sample_vector(size, out); }
00604 inline vec randu(int size) { vec temp; randu(size, temp); return temp; }
00606 inline void randu(int rows, int cols, mat &out) { Uniform_RNG src; src.sample_matrix(rows, cols, out); }
00608 inline mat randu(int rows, int cols) { mat temp; randu(rows, cols, temp); return temp; }
00609 
00611 inline int randi(int low, int high) { I_Uniform_RNG src; src.setup(low, high); return src(); }
00613 inline ivec randi(int size, int low, int high) { I_Uniform_RNG src; src.setup(low, high); return src(size); }
00615 inline imat randi(int rows, int cols, int low, int high) { I_Uniform_RNG src; src.setup(low, high); return src(rows, cols); }
00616 
00618 inline vec randray(int size, double sigma = 1.0) { Rayleigh_RNG src; src.setup(sigma); return src(size); }
00619 
00621 inline vec randrice(int size, double sigma = 1.0, double s = 1.0) { Rice_RNG src; src.setup(sigma, s); return src(size); }
00622 
00624 inline vec randexp(int size, double lambda = 1.0) { Exponential_RNG src; src.setup(lambda); return src(size); }
00625 
00627 inline double randn(void) { Normal_RNG src; return src.sample(); }
00629 inline void randn(int size, vec &out) { Normal_RNG src; src.sample_vector(size, out); }
00631 inline vec randn(int size) { vec temp; randn(size, temp); return temp; }
00633 inline void randn(int rows, int cols, mat &out)  { Normal_RNG src; src.sample_matrix(rows, cols, out); }
00635 inline mat randn(int rows, int cols) { mat temp; randn(rows, cols, temp); return temp; }
00636 
00641 inline std::complex<double> randn_c(void) { Complex_Normal_RNG src; return src.sample(); }
00643 inline void randn_c(int size, cvec &out)  { Complex_Normal_RNG src; src.sample_vector(size, out); }
00645 inline cvec randn_c(int size) { cvec temp; randn_c(size, temp); return temp; }
00647 inline void randn_c(int rows, int cols, cmat &out) { Complex_Normal_RNG src; src.sample_matrix(rows, cols, out); }
00649 inline cmat randn_c(int rows, int cols) { cmat temp; randn_c(rows, cols, temp); return temp; }
00650 
00652 
00653 } // namespace itpp
00654 
00655 #endif // #ifndef RANDOM_H
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines
SourceForge Logo

Generated on Wed Jul 27 2011 16:27:04 for IT++ by Doxygen 1.7.4