//--------------------------------------------------------------------------
//
// Description:
// 	class ClassEngine :
// encapsulation of class calls
//
//
// Author List:
//	Stephane Plaszczynski (plaszczy@lal.in2p3.fr)
//
// History (add to end):
//	creation:   ven. nov. 4 11:02:20 CET 2011 
//
//-----------------------------------------------------------------------

#ifndef ClassEngine_hh
#define ClassEngine_hh

//CLASS
extern  "C"{
#include"class.h"
}
#include"Engine.hh"
//STD
#include<string>
#include<vector>
#include<utility>
#include<ostream>
#include<sstream>

using std::string;

//general utility to convert safely numerical types to string
template<typename T> std::string str(const T &x);
//specialisations
template<> std::string str (const float &x);
template<> std::string str (const double &x);
template<> std::string str (const bool &x); //converts to "yes" or "no"
template<> std::string str (const std::string &x);

std::string str(const char* x);
//////////////////////////////////////////////////////////////////////////
//class to encapsulate CLASS parameters from any type (numerical or string)
class ClassParams{
public:

  ClassParams(){};
  ClassParams( const ClassParams& o):pars(o.pars){};

  //add class parameter with single value
  template<typename T> unsigned inline 
  add(const string& key,const T& val){
  pars.push_back(make_pair(key,str(val)));
  return pars.size();
  }

  //returs a class string from a list of values (assumed not empty)
  template<typename T> inline string classList(const std::vector<T>& vals){
    string valname=str(vals[0]);
    for (size_t i=1;i<vals.size();i++){
      valname+=",";
      valname+=str(vals[i]);
    }
    return valname;
  }

  //add a class parameter with a list of values
  template<typename T> inline unsigned addList(const string& key,const std::vector<T>& vals){
    string s=classList(vals);
    pars.push_back(make_pair(key,s));
    return pars.size();
  }


  //accesors
  size_t size() const;
  size_t findIndex(const string& keyname);
  string key(const unsigned& i) const;
  string value(const unsigned& i) const;
  //modifier
  void updateVal(const unsigned& i,const string& newval);


private:
  //ordered list of key-value pairs
  std::vector<std::pair<string,string> > pars;
};

///////////////////////////////////////////////////////////////////////////
class ClassEngine : public Engine
{

  friend class ClassParams;

public:
  //constructors
  ClassEngine(const ClassParams& pars);
  //with a class .pre file
  ClassEngine(const ClassParams& pars,const string & precision_file);
  

  // destructor
  ~ClassEngine();

  //obsolote: to be removed
  bool updateParValues(const std::vector<double>& par);

  //update from ClassParams class (should be same size than in constructor)
  bool updateParValues(const ClassParams& par);

  //services
  virtual bool has_Background() const {return true;}
//   //virtual bool has_CMB_Cl() const {return (_output.find("tCl")!=std::string::npos);}
//   virtual bool has_CMB_ClPol() const {return (_output.find("pCl")!=std::string::npos);}
//   virtual bool has_CMB_Lensing() const {return (_output.find("lCl")!=std::string::npos);}
//   virtual bool has_Pklin() const {return (_output.find("mPk")!=std::string::npos);}
//   virtual bool has_PkNL() const {return _nonlin;}

  //class based
  inline bool has_CMB_Cl() const {return (pt.has_cl_cmb_temperature==_TRUE_);}
  inline bool has_CMB_ClPol() const {return (pt.has_cl_cmb_polarization==_TRUE_);}
  inline bool has_CMB_Lensing() const {return (pt.has_cl_lensing_potential==_TRUE_);}
  inline bool has_Pklin() const {return (pt.has_pk_matter==_TRUE_);}
  inline bool has_PkNL() const {return (pt.has_nl_corrections_based_on_delta_m==_TRUE_);}


  //get value at l ( 2<l<lmax): in units = (micro-K)^2
  //don't call if FAILURE returned previously
  //throws std::execption if pb
  virtual std::string name() const {return "ClassEngine";}

  double getCl(Engine::cltype t,const long &l);  
  void getCls(const std::vector<unsigned>& lVec, //input 
	      std::vector<double>& cltt, 
	      std::vector<double>& clte, 
	      std::vector<double>& clee, 
	      std::vector<double>& clbb);

  
  bool getLensing(const std::vector<unsigned>& lVec, //input 
	      std::vector<double>& clphiphi, 
	      std::vector<double>& cltphi, 
	      std::vector<double>& clephi);


  //recombination
  inline double z_rec() const {return th.z_rec;}
  inline double rs_rec() const {return th.rs_rec;} 

 //baryon drag 
  inline double z_drag() const {return th.z_d;}
  inline double rs_drag() const {return th.rs_d;} 

  double get_Dv(double z);
  double get_Da(double z);
  double get_sigma8(double z);
  double get_f(double z);

  double get_Fz(double z);
  double get_H(double z);
  double get_DMod(double z);


  //PK related stuuff
  double get_Pklin(double k, double z);
  double get_PkNL(double k, double z);

  //nodes used by class to perfrom spline interpolation (at z=0)
  //knodes/pknodes resized
  bool get_PklinNodes(std::vector<double>& knodes,std::vector<double>& Pknodes);
  bool get_PkNLNodes(std::vector<double>& knodes,std::vector<double>& Pknodes);
  // the followingdepends on switch (prefer previous methods)
  inline bool get_PkNodes(std::vector<double>& knodes,std::vector<double>& Pknodes) 
  {return (has_PkNL()? get_PkNLNodes(knodes,Pknodes) : get_PklinNodes(knodes,Pknodes));} 

#ifdef FASTPK  
  //vectorized access (if macro defined): with switch according to has_NL 
  // outout vector pks must be properly sized (see Engine.hh)
  bool get_Pkvec(const std::vector<double>& kin, const std::vector<double>& z,std::vector<double>& pks);  
  inline bool get_Pkvec(const std::vector<double>& kvec,const double& z, std::vector<double>& pks){
    std::vector<double> zvec(1,z);
    return  this->get_Pkvec(kvec,zvec,pks);
  }
#endif



  //combile distance
  double com_distance(double z);

  //linear growth density term
  double get_growthD(double z);

  //expansion
  inline double get_H0() const {return ba.H0*_c_/1e3;}
  inline double get_h() const {return get_H0()/100.;}
  inline double theta() const {return th.rs_rec/th.ra_rec;}

  inline double Omega_m() const {return ba.Omega0_b+ba.Omega0_cdm+ba.Omega0_ncdm_tot;}
  //reio
  inline double get_tau_reio() const {return th.tau_reio;}
  double get_z_reio() const {return th.z_reio;}


  inline double YHe() const { return th.YHe;}

  //may need that
  inline int numCls() const {return sp.ct_size;};
  inline double Tcmb() const {return ba.T_cmb;}

  inline int l_max_scalars() const {return _lmax;}

  //print content of file_content
  void printFC();


  //access to class content:dangerous to use (but sometimes mandatory...)
  struct file_content fc;

  //structures class en commun
  struct precision pr;        /* for precision parameters */
  struct background ba;       /* for cosmological background */
  struct thermo th;           /* for thermodynamics */
  struct perturbs pt;         /* for source functions */
  struct transfers tr;        /* for transfer functions */
  struct primordial pm;       /* for primordial spectra */
  struct spectra sp;          /* for output spectra */
  struct nonlinear nl;        /* for non-linear spectra */
  struct lensing le;          /* for lensed spectra */
  struct output op;           /* for output files */

private:

  ErrorMsg _errmsg;            /* for error messages */
  double * cl;

  //helpers
  bool dofree;
  int freeStructs();

  //call once /model
  int computeCls();

  int class_main(
		 struct file_content *pfc,
		 struct precision * ppr,
		 struct background * pba,
		 struct thermo * pth,
		 struct perturbs * ppt,
		 struct transfers * ptr,
		 struct primordial * ppm,
		 struct spectra * psp,
		 struct nonlinear * pnl,
		 struct lensing * ple,
		 struct output * pop,
		 ErrorMsg errmsg);
  //parnames
  std::vector<std::string> parNames;


protected:
 
  std::string _output;
  double* pvecback; // helper to avoid too much allocations


};

#endif