source: Sophya/trunk/SophyaLib/BaseTools/randinterf.h

#ifndef RANDOMGENERATORINTERFACE_H_SEEN
#define RANDOMGENERATORINTERFACE_H_SEEN

// Classes RandomGeneratorInterface 
//
// R. Ansari          (C) UPS+LAL IN2P3/CNRS  2009
// C. Magneville      (C) DAPNIA/SPP  CEA     2009


#include "machdefs.h"
#include "anydataobj.h"
#include <complex>
#include <vector>

namespace SOPHYA {

//! enum definition for usual types of probability distributio finctions
enum RNDTypes {
  C_RND_Flat = 0,
  C_RND_Gaussian = 1,
  C_RND_Poisson = 2,
  C_RND_Exponential = 3
};
 
//! enum definition for the type of algorithm used for Gaussian distribution generation 
enum GaussianGenAlgo {
  C_Gaussian_BoxMuller = 0, 
  C_Gaussian_RandLibSNorm = 1,
  C_Gaussian_PolarBoxMuller = 2,
  C_Gaussian_RatioUnif = 3,
  C_Gaussian_LevaRatioUnif = 4,
  C_Gaussian_Ziggurat128 = 5
};

//! enum definition for the type of algorithm used for Poisson distribution generation 
enum PoissonGenAlgo {
  C_Poisson_Simple = 0, 
  C_Poisson_Ahrens = 1
};

//! enum definition for the type of algorithm used for exponential distribution generation 
enum ExponentialGenAlgo {
  C_Exponential_Simple = 0, 
  C_Exponential_Ahrens = 1
};

// Definition de l'interface pour les classes generateur aleatoire
class RandomGeneratorInterface : public AnyDataObj {
public:
  RandomGeneratorInterface();
  virtual ~RandomGeneratorInterface();

//! Select Gaussian generation algorithm 
  inline void SelectGaussianAlgo(GaussianGenAlgo typ=C_Gaussian_BoxMuller) 
    {  usegaussian_ = typ; }
//! Select Poisson generation algorithm 
  inline void SelectPoissonAlgo(PoissonGenAlgo typ=C_Poisson_Simple) 
    {  usepoisson_ = typ; }
//! Select Exponential generation algorithm 
  inline void SelectExponentialAlgo(ExponentialGenAlgo typ=C_Exponential_Simple) 
    {  useexpo_ = typ; }

//! Return the Gaussian generation algorithm type
  inline GaussianGenAlgo GetGaussianAlgo() 
    {  return usegaussian_; }
//! Return the Poisson generation algorithm type
  inline PoissonGenAlgo GetPoissonAlgo() 
    {  return usepoisson_; }
//! Return the Exponential generation algorithm type
//! Select Exponential generation algorithm 
  inline ExponentialGenAlgo GetExponentialAlgo() 
    {  return useexpo_; }

  virtual void GenerateSeedVector(int nseed,vector<uint_2>& seed,int lp=0);

  // --- Le tirage sur une distribution plate
  /*! \brief Return a random number \b r with a flat distribution  0<=r<1 */ 
  inline r_8 Flat01() {return Next();}
  /*! \brief Return a random number \b r with a flat distribution  -1<=r<1 */ 
  inline r_8 Flatpm1() {return 2.*Next()-1.;}

  // --- Le tirage sur une distribution gaussienne
  /*! \brief Return a random number following a gaussian distribution (sigma=1, mean=0)*/ 
  virtual r_8 Gaussian();

  virtual r_8 GaussianBoxMuller(); 
  virtual r_8 GaussianSNorm(); 
  virtual r_8 GaussianPolarBoxMuller(); 
  virtual r_8 GaussianRatioUnif(); 
  virtual r_8 GaussianLevaRatioUnif();
  virtual r_8 GaussianZiggurat128(); 
  /*! \brief Return a random number following a gaussian distribution "sigma", (mean=0)*/ 
  inline r_8 Gaussian(double sigma) {return sigma*Gaussian();}
 /*! \brief Return a random number following a gaussian distribution with mean=mu, and sigma */ 
  inline r_8 Gaussian(double sigma,double mu) {return sigma*Gaussian()+mu;}

  /*! \brief Return a random number following a gaussian tail distribution for x>sdev */ 
  virtual r_8 GaussianTail(double sdev,double slim=1.);

  // --- Le tirage sur une distribution de poisson
  /*! \brief Return a random number following a poisson distribution with mean mu */ 
  virtual uint_8 Poisson(double mu, double mumax=-1);

  virtual uint_8 PoissonSimple(double mu, double mumax=-1);
  virtual uint_8 PoissonAhrens(double mu);

  //! Return a random number with exponential distribution 
  virtual r_8 Exponential();

  virtual r_8 ExpoSimple();
  virtual r_8 ExpoAhrens();

  // --- Le tirage gaussien complexe (cf texte a la fin du .cc)
  inline complex< r_8 > ComplexGaussian(void)
         {return complex< r_8 >(Gaussian(),Gaussian());}
  inline complex< r_8 > ComplexGaussian(double sig)
         {return complex< r_8 >(sig*Gaussian(),sig*Gaussian());}
  /*! \brief Returns the module of a random complex number generated by ComplexGaussRan */
  inline double ModComplexGaussian(double sig=1.)
         {double r=-log(1.-Next()); return sig*sqrt(2.*r);}

  // --- Le tirage sur une distribution gaussienne a 2D
  /*! \brief Return 2 random numbers following a gaussian 2D distribution */ 
  virtual int  Gaussian2DRho(double &x,double &y,double mx,double my,double sx,double sy,double ro);
  virtual void Gaussian2DAng(double &x,double &y,double mx,double my,double sa,double sb,double teta);

  //! Return the pointer to the default global RandomGenerator object 
  static inline RandomGeneratorInterface* GetGlobalRandGenP()
    { return gl_rndgen_p; }

  //! Initialisation automatique
  virtual void AutoInit(int lp=0);

  //! Print
  virtual void ShowRandom();

  // Permet de definir l'instance global du generateur aleatoire
  static void SetGlobalRandGenP(RandomGeneratorInterface* rgp);

protected:
  /*! \brief Should return a random number in the range [0,1] with flat distribution 

     This pure virtual method should be implemented by inheriting classes  */ 
  virtual r_8 Next() = 0;

// Selection de type d'algo pour les aleatoires autres que plates
  GaussianGenAlgo usegaussian_;
  PoissonGenAlgo usepoisson_;
  ExponentialGenAlgo useexpo_;

  // Instance global d'un generateur par defaut 
  static RandomGeneratorInterface* gl_rndgen_p;
};

} /* namespace SOPHYA */

#endif