#include "LHAPDFWrap.h"
// adapted for LHAPDFv4 by Mike Whalley
// adapted for LHAPDFv5 by Craig Group/Mike Whalley

LHAPDFWrap::LHAPDFWrap() {}

LHAPDFWrap::~LHAPDFWrap() {}

LHAPDFWrap::LHAPDFWrap(char *name, int memset) {
  finitpdfsetbyname_(*name);
  finitpdf_(memset);
}

LHAPDFWrap::LHAPDFWrap(char *name) {
  finitpdfsetbyname_(*name);
}


//For multiple set use
LHAPDFWrap::LHAPDFWrap(int nset, char *name, int memset) {
  finitpdfsetbynamem_(nset, *name);
  finitpdfm_(nset, memset);
}

LHAPDFWrap::LHAPDFWrap(int nset, char *name) {
  finitpdfsetbynamem_(nset, *name);
}

std::vector<double> LHAPDFWrap::xfx(const double &x, const double &Q) {  
  double f[13], mx = x, mQ = Q;
  fevolvepdf_(mx, mQ, f);
  std::vector<double> r;
  for (int i=0; i<13; i++) r.push_back(f[i]);
  return r;
}

double LHAPDFWrap::xfx(const double &x, const double &Q, int fl) {  
  double f[13], mx = x, mQ = Q;
  fevolvepdf_(mx, mQ, f);
  return f[fl+6];
}

std::vector<double> LHAPDFWrap::xfxp(const double &x, const double &Q, const double &P2, int ip) {  
  double f[13], mx = x, mQ = Q, mP2 = P2;
  int  mip = ip;
  fevolvepdfp_(mx, mQ, mP2, mip, f);
  std::vector<double> r;
  for (int i=0; i<13; i++) r.push_back(f[i]);
  return r;
}

double LHAPDFWrap::xfxp(const double &x, const double &Q, const double &P2, int ip, int fl) {  
  double f[13], mx = x, mQ = Q, mP2 = P2;
  int mip = ip;
  fevolvepdfp_(mx, mQ, mP2, mip, f);
  return f[fl+6];
}

std::vector<double> LHAPDFWrap::xfxa(const double &x, const double &Q, const double &a) {  
  double f[13], mx = x, mQ = Q, ma = a;
  fevolvepdfa_(mx, mQ, ma, f);
  std::vector<double> r;
  for (int i=0; i<13; i++) r.push_back(f[i]);
  return r;
}

double LHAPDFWrap::xfxa(const double &x, const double &Q, const double &a, int fl) {  
  double f[13], mx = x, mQ = Q, ma = a;
  fevolvepdfa_(mx, mQ, ma, f);
  return f[fl+6];
}

std::vector<double> LHAPDFWrap::xfxphoton(const double &x, const double &Q) {  
  double f[13], mx = x, mQ = Q, mphoton; 
  fevolvepdfphoton_(mx, mQ, f, mphoton );
  std::vector<double> r;
  for (int i=0; i<13; i++) r.push_back(f[i]);
  r.push_back(mphoton);
  return r;
}

double LHAPDFWrap::xfxphoton(const double &x, const double &Q, int fl) {  
  double f[13], fout[14], mx = x, mQ = Q, mphoton;
  fevolvepdfphoton_(mx, mQ, f, mphoton);
  for (int i=0; i<13; i++) fout[i] = f[i];
  fout[13] = mphoton;
  return fout[fl+6];
}

int LHAPDFWrap::numberPDF() {
  int N; 
  fnumberpdf_(N); 
  return N;
}

double LHAPDFWrap::alphasPDF(double Q) {
  double a;
  falphaspdf_(Q, a);
  return a;
}

int LHAPDFWrap::getOrderPDF() {
  int N;
  fgetorderpdf_(N);
  return N;
}

int LHAPDFWrap::getOrderAlphaS() {
  int N;
  fgetorderas_(N);
  return N;
}

double LHAPDFWrap::getQMass(int f) {
  double m;
  fgetqmass_(f, m);
  return m;
}

double LHAPDFWrap::getThreshold(int f) {
  double m;
  fgetthreshold_(f, m);
  return m;
}

int LHAPDFWrap::getNf() {
  int N; 
  fgetnf_(N); 
  return N;
}

double LHAPDFWrap::getLam4(int m){
  double l;
  fgetlam4_(m, l);
  return l;
}

double LHAPDFWrap::getLam5(int m){
  double l;
  fgetlam5_(m, l);
  return l;
}

double LHAPDFWrap::getXmin(int m){
  double xmin;
  fgetxmin_(m, xmin);
  return xmin;
}

double LHAPDFWrap::getXmax(int m){
  double xmax;
  fgetxmax_(m, xmax);
  return xmax;
}

double LHAPDFWrap::getQ2min(int m){
  double q2min;
  fgetq2min_(m, q2min);
  return q2min;
}

double LHAPDFWrap::getQ2max(int m){
  double q2max;
  fgetq2max_(m, q2max);
  return q2max;
}


//Functions below are the same as above but with M appended to their name.
//They are for multiple set use in LHAPDF.  They also take an extra paramater
//which is a lable for their memory location 


std::vector<double> LHAPDFWrap::xfxM(int nset, const double &x, const double &Q) {  
  double f[13], mx = x, mQ = Q;
  fevolvepdfm_(nset, mx, mQ, f);
  std::vector<double> r;
  for (int i=0; i<13; i++) r.push_back(f[i]);
  return r;
}

double LHAPDFWrap::xfxM(int nset, const double &x, const double &Q, int fl) {  
  double f[13], mx = x, mQ = Q;
  fevolvepdfm_(nset, mx, mQ, f);
  return f[fl+6];
}

std::vector<double> LHAPDFWrap::xfxpM(int nset, const double &x, const double &Q, const double &P2, int ip) {  
  double f[13], mx = x, mQ = Q, mP2 = P2;
  int  mip = ip;
  fevolvepdfpm_(nset, mx, mQ, mP2, mip, f);
  std::vector<double> r;
  for (int i=0; i<13; i++) r.push_back(f[i]);
  return r;
}

double LHAPDFWrap::xfxpM(int nset, const double &x, const double &Q, const double &P2, int ip, int fl) {  
  double f[13], mx = x, mQ = Q, mP2 = P2;
  int mip = ip;
  fevolvepdfpm_(nset, mx, mQ, mP2, mip, f);
  return f[fl+6];
}

std::vector<double> LHAPDFWrap::xfxaM(int nset, const double &x, const double &Q, const double &a) {  
  double f[13], mx = x, mQ = Q, ma = a;
  fevolvepdfam_(nset, mx, mQ, ma, f);
  std::vector<double> r;
  for (int i=0; i<13; i++) r.push_back(f[i]);
  return r;
}

double LHAPDFWrap::xfxaM(int nset, const double &x, const double &Q, const double &a, int fl) {  
  double f[13], mx = x, mQ = Q, ma = a;
  fevolvepdfam_(nset, mx, mQ, ma, f);
  return f[fl+6];
}

std::vector<double> LHAPDFWrap::xfxphotonM(int nset, const double &x, const double &Q) {  
  double f[13], mx = x, mQ = Q, mphoton;
  fevolvepdfphotonm_(nset, mx, mQ, f, mphoton );
  std::vector<double> r;
  for (int i=0; i<13; i++) r.push_back(f[i]);
  r.push_back(mphoton);
  return r;
}

double LHAPDFWrap::xfxphotonM(int nset, const double &x, const double &Q, int fl) {  
  double f[13], fout[14], mx = x, mQ = Q, mphoton;
  fevolvepdfphotonm_(nset, mx, mQ, f, mphoton);
  for (int i=0; i<13; i++) fout[1] = f[i];
  fout[14] = mphoton;
  return fout[fl+6];
}

int LHAPDFWrap::numberPDFM(int nset) {
  int N; 
  fnumberpdfm_(nset, N); 
  return N;
}

double LHAPDFWrap::alphasPDFM(int nset, double Q) {
  double a;
  falphaspdfm_(nset, Q, a);
  return a;
}

int LHAPDFWrap::getOrderPDFM(int nset) {
  int N;
  fgetorderpdfm_(nset, N);
  return N;
}

int LHAPDFWrap::getOrderAlphaSM(int nset) {
  int N;
  fgetorderasm_(nset, N);
  return N;
}

double LHAPDFWrap::getQMassM(int nset, int f) {
  double m;
  fgetqmassm_(nset, f, m);
  return m;
}

double LHAPDFWrap::getThresholdM(int nset, int f) {
  double m;
  fgetthresholdm_(nset, f, m);
  return m;
}

int LHAPDFWrap::getNfM(int nset) {
  int N; 
  fgetnfm_(nset, N); 
  return N;
}

double LHAPDFWrap::getLam4M(int nset, int m){
  double l;
  fgetlam4m_(nset, m, l);
  return l;
}

double LHAPDFWrap::getLam5M(int nset, int m){
  double l;
  fgetlam5m_(nset, m, l);
  return l;
}

double LHAPDFWrap::getXminM(int nset,int m){
  double xmin;
  fgetxminm_(nset, m, xmin);
  return xmin;
}

double LHAPDFWrap::getXmaxM(int nset,int m){
  double xmax;
  fgetxmaxm_(nset, m, xmax);
  return xmax;
}

double LHAPDFWrap::getQ2minM(int nset,int m){
  double q2min;
  fgetq2minm_(nset, m, q2min);
  return q2min;
}

double LHAPDFWrap::getQ2maxM(int nset,int m){
  double q2max;
  fgetq2maxm_(nset, m, q2max);
  return q2max;
}