#ifndef STSRAND_H_SEEN
#define STSRAND_H_SEEN

// Classe RandomGenerator 
// Generateur aleatoire compatible multi-thread
//
// R. Ansari          (C) UPS+LAL IN2P3/CNRS  2007
// C. Magneville      (C) DAPNIA/SPP  CEA     2007


#include "machdefs.h"
#include "objfio.h"   // Pour rendre la classe PPersist

#include "ndatablock.h"

namespace SOPHYA {

class STSRandGen : public AnyDataObj {

 public:
  STSRandGen(size_t n = 1024, bool tsafe=true);
  STSRandGen(STSRandGen const & rg);

  void SetBuffSize(size_t n);

  virtual ~STSRandGen();

  /*! \brief Automatic initialization using the current time */
  static  void AutoInit(int lp=0);
  /*! \brief Initialization using the specified seed  */
  static  void Init(long seed_val, int lp=0);
  /*! \brief Initialization using the specified seed  */
  static  void Init(unsigned short seed_16v[3], int lp=0);
  /*! \brief Getting the current seed state */
  static  void GetSeed(unsigned short seed_16v[3], int lp=0);

  /*! \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.;}

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

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

  /*! \brief alias for Gaussian() */ 
  inline r_8 NorRand() {return Gaussian();}

  /*! \brief Return a random number following a poisson distribution with mean mu */ 
  virtual  uint_8 Poisson(double mu, double mumax);
  
  /*! \brief Return a random number following a poisson distribution with mean mu */ 
  inline uint_8 Poisson(double mu) {return Poisson(mu, -1.);}

   /*! \brief Return a random number following a poisson distribution with mean mu */ 
  inline uint_8  PoissRandLimit(double mu,double mumax=10.) {return Poisson(mu, mumax);}

  /*! \brief Returns a random complex number such that real and imaginary parts are gaussians with variance sig^2 */
  inline complex< r_8 > ComplexGaussRan(void)
         {return complex< r_8 >(NorRand(),NorRand());}
  inline complex< r_8 > ComplexGaussRan(double sig)
         {return complex< r_8 >(sig*NorRand(),sig*NorRand());}
  /*! \brief Returns the module of a random complex number generated by ComplexGaussRan */
  inline double ModComplexGaussRan(double sig=1.)
         {double r=-log(1.-Flat01()); return sig*sqrt(2.*r);}

 
  //  Pour la gestion de persistance PPF
  friend class ObjFileIO<STSRandGen> ;

 protected:
  // ---- protected data members
  NDataBlock<r_8>  rseq_;  
  size_t idx_;
  bool fg_nothrsafe;  // if true --> NOT thread-safe
  // ---- protected methods 
  void GenSeq();
  inline r_8 Next() 
    {
      if (rseq_.Size() == 0)  return drand48(); 
      else { 
	if(idx_==rseq_.Size()) GenSeq(); 
	return(rseq_(idx_++));
      }
    }
  // Non thread-safe version of Init() and GetSeed()
  static  void Init_P(unsigned short seed_16v[3]);
  static  void GetSeed_P(unsigned short seed_16v[3]);

};  // Fin de la classe STSRandGen

// Classe pour la gestion de persistance PPF :  ObjFileIO<STSRandGen>

/*! Writes the random generator object state in the POutPersist stream \b os */
inline POutPersist& operator << (POutPersist& os, STSRandGen & obj)
{ ObjFileIO<STSRandGen> fio(&obj);  fio.Write(os);  return(os); }
/*! Reads the random generator object state from the PInPersist stream \b is */
inline PInPersist& operator >> (PInPersist& is, STSRandGen & obj)
{ ObjFileIO<STSRandGen> fio(&obj); is.SkipToNextObject(); fio.Read(is); return(is); }

} /* namespace SOPHYA */

#endif 
