MRST.h

// -*- C++ -*-
//
// MRST.h is a part of Herwig++ - A multi-purpose Monte Carlo event generator
// Copyright (C) 2002-2011 The Herwig Collaboration
//
// Herwig++ is licenced under version 2 of the GPL, see COPYING for details.
// Please respect the MCnet academic guidelines, see GUIDELINES for details.
//
#ifndef HERWIG_MRST_H
#define HERWIG_MRST_H

#include <ThePEG/PDF/PDFBase.h>
#include <iostream>
#include <fstream>
#include <string>
#include <cmath>

namespace Herwig {

using namespace ThePEG;

class MRST : public PDFBase {
  enum PDFFlavour { upValence = 1, dnValence, glu, upSea, chm, str, bot, dnSea };

  enum PDFType {Sea,Valence,Total};

public:

  MRST();
 
public:

  virtual bool canHandleParticle(tcPDPtr particle) const;

  virtual cPDVector partons(tcPDPtr p) const;

  virtual double xfx(tcPDPtr particle, tcPDPtr parton, Energy2 partonScale,
                     double x, double eps = 0.0,
                     Energy2 particleScale = ZERO) const;

  virtual double xfvx(tcPDPtr particle, tcPDPtr parton, Energy2 partonScale,
                      double x, double eps = 0.0,
                      Energy2 particleScale = ZERO) const;

  virtual double xfsx(tcPDPtr particle, tcPDPtr parton, Energy2 partonScale,
                      double x, double eps = 0.0,
                      Energy2 particleScale = ZERO) const;

public:

  void persistentOutput(PersistentOStream & os) const;

  void persistentInput(PersistentIStream & is, int version);

  static void Init();

protected:

  inline virtual IBPtr clone() const;

  inline virtual IBPtr fullclone() const;

protected:

  virtual void doinitrun();

private:

  static ClassDescription<MRST> initMRST;

  MRST & operator=(const MRST &);

private:

  inline double lookup(PDFFlavour i, int n, int m, double u, double t) const;


  double pdfValue(double x, Energy2 q2, 
                  tcPDPtr particle, tcPDPtr parton,PDFType type) const;

  inline int locate(double xx[],int n,double x) const;

  inline double polderivative(double x1, double x2, double x3,
                              double y1, double y2, double y3) const;

  virtual void readSetup(istream &);

  void initialize(bool reread = true);

private:

  static const int np=8;
  static const int nx=49;
  
  static const int nq=37;

  static const int nqc0=2; 
  
  static const int nqb0=11;

  static const int ntenth=23;
  
  static const double xmin;

  static const double xmax;
  
  static const Energy2 qsqmin;
  
  static const Energy2 qsqmax;
  
  static const Energy2 mc2;
  
  static const Energy2 mb2;

  unsigned _inter;

  double _xswitch;

  string _file;

  //  double data[np+1][nx+1][nq+1];
  vector<vector<vector<double> > > data;

  //  double data[np+1][nx+1][nq+1];
  vector<vector<vector<double> > > fdata;

  static double xx[nx+1];

  static double lxx[nx+1];

  static double lxxb[nx+1];

  static double qq[nq+1];

  static double lqq[nq+1];

  double c[np+1][nx][nq][5][5];

  static double n0[np+1];

  static bool initialized;
};

}

namespace ThePEG {

template <>
struct BaseClassTrait<Herwig::MRST,1> {
  typedef PDFBase NthBase;
};

template <>
struct ClassTraits<Herwig::MRST>: public ClassTraitsBase<Herwig::MRST> {
  static string className() { return "Herwig::MRST"; }
  static string library() { return "HwMRST.so"; }
};

}

#include "MRST.icc"

#endif