gMakePhotonStrucFunc.cxx File Reference

#include "Framework/Messenger/Messenger.h"
#include "Framework/Conventions/Constants.h"
#include "Framework/ParticleData/PDGCodes.h"
#include "Framework/ParticleData/PDGUtils.h"
#include <TSystem.h>
#include <TMath.h>
#include <Math/IFunction.h>
#include <Math/Integrator.h>
#include <string>
#include <iostream>
#include <fstream>
#include "LHAPDF/LHAPDF.h"

Include dependency graph for gMakePhotonStrucFunc.cxx:

Go to the source code of this file.

Functions
int	main (int argc, char **argv)

Variables
int	fNucPdg = 0

Function Documentation

main()

int main	(	int	argc,
		char **	argv )

Definition at line 111 of file gMakePhotonStrucFunc.cxx.

{
 
  string basedir = "";
  if ( gSystem->Getenv("PHOTON_SF_DATA_PATH")==NULL ) basedir = string(gSystem->Getenv("GENIE")) + "/data/evgen/photon-sf";
  else                                                basedir = string(gSystem->Getenv("PHOTON_SF_DATA_PATH"));
  LOG("gmkphotonsf", pERROR) << "Base directory: " << basedir;
 
  if ( gSystem->AccessPathName( basedir.c_str(), kWritePermission ) ) {
      LOG("gmkphotonsf", pFATAL) << "Base directory doesnt exist or you dont have write permission.";
      LOG("gmkphotonsf", pFATAL) << "Remember!!!";
      LOG("gmkphotonsf", pFATAL) << "Path to base directory is defined with the enviroment variable PHOTON_SF_DATA_PATH.";
      LOG("gmkphotonsf", pFATAL) << "If not defined, default location is $GENIE/data/evgen/photon-sf";
      assert(0);
  }
 
  const int nx = 1000.;
 
  int nucs[2] = { kPdgProton, kPdgNeutron };
  int pdgs[6] = { kPdgNuE, kPdgAntiNuE, kPdgNuMu, kPdgAntiNuMu, kPdgNuTau, kPdgAntiNuTau };
  
#ifdef __GENIE_APFEL_ENABLED__
  for (int k=0; k<2; k++) {
    
    fNucPdg = nucs[k];
 
    // initialising APFEL framework
    LOG("gmkphotonsf", pINFO) << "Initialising APFEL..." ; 
    string pdfset;
    if      (pdg::IsProton (fNucPdg) ) pdfset = "NNPDF31_nnlo_as_0118_luxqed";
    else if (pdg::IsNeutron(fNucPdg) ) pdfset = "NNPDF31_nnlo_as_0118";
 
    const LHAPDF::PDFSet set(pdfset);
    pdf_cache = LHAPDF::mkPDF(pdfset,0);
 
    double xPDFmin,QPDFmin,QPDFmax,mc,mb,mt;
    stringstream( set.get_entry("MCharm") )  >> mc;
    stringstream( set.get_entry("MBottom") ) >> mb;
    stringstream( set.get_entry("MTop") )    >> mt;    
    stringstream( set.get_entry("QMin") )    >> QPDFmin;
    stringstream( set.get_entry("QMax") )    >> QPDFmax;  
    stringstream( set.get_entry("XMin") )    >> xPDFmin;
    LOG("gmkphotonsf", pINFO) << "xPDFmin = " << xPDFmin;
    LOG("gmkphotonsf", pINFO) << "QPDFmin = " << QPDFmin;
    LOG("gmkphotonsf", pINFO) << "QPDFmax = " << QPDFmax;
    LOG("gmkphotonsf", pINFO) << "mc = " << mc;
    LOG("gmkphotonsf", pINFO) << "mb = " << mb;
    LOG("gmkphotonsf", pINFO) << "mt = " << mt;
 
    APFEL::CleanUp();
 
    APFEL::SetPDFSet(pdfset);
    APFEL::SetReplica(0);
    APFEL::SetPerturbativeOrder(2);
    APFEL::SetQLimits(QPDFmin,QPDFmax);
    APFEL::SetMaxFlavourPDFs(6);
    APFEL::SetMaxFlavourAlpha(6);
    APFEL::SetNumberOfGrids(3);
    APFEL::SetGridParameters(1,100,5,1e-9);
    APFEL::SetGridParameters(2,40,3,1e-1);
    APFEL::SetGridParameters(3,60,3,8e-1);
    APFEL::SetGFermi(kGF);
    APFEL::SetPoleMasses(mc,mb,mt);
    APFEL::SetTheory("QUniD");
    APFEL::EnableLeptonEvolution(true);
    APFEL::SetFastEvolution(true);
    APFEL::SetPDFSet("leptexternal");
    APFEL::InitializeAPFEL();
 
    LOG("gmkphotonsf", pWARN) << "Init EvolveAPFEL";        
    APFEL::EvolveAPFEL(QPDFmin,kMw);
    LOG("gmkphotonsf", pWARN) << "End EvolveAPFEL";        
 
    // open file in which SF will be stored
 
    double x[nx];
    for ( int i=0; i<nx; i++ ) x[i] = TMath::Power( 10, TMath::Log10(xPDFmin) + i*(TMath::Log10(1.)-TMath::Log10(xPDFmin))/(1000.-1) );
 
    for(int j=0; j<6; j++) {
 
      string SFname = basedir + "/PhotonSF_hitnuc"+to_string(fNucPdg)+"_hitlep"+to_string(pdgs[j])+".dat";
      std::ofstream sf_stream(SFname);
      for ( int i=0; i<nx; i++ ) {
        double tmp = 0;
        if      ( pdg::IsNuE      (pdgs[j]) ) tmp = APFEL::xLepton( 1,x[i]);
        else if ( pdg::IsAntiNuE  (pdgs[j]) ) tmp = APFEL::xLepton(-1,x[i]);
        else if ( pdg::IsNuMu     (pdgs[j]) ) tmp = APFEL::xLepton( 2,x[i]);
        else if ( pdg::IsAntiNuMu (pdgs[j]) ) tmp = APFEL::xLepton(-2,x[i]);
        else if ( pdg::IsNuTau    (pdgs[j]) ) tmp = APFEL::xLepton( 3,x[i]);
        else if ( pdg::IsAntiNuTau(pdgs[j]) ) tmp = APFEL::xLepton(-3,x[i]);
        LOG("gmkphotonsf", pWARN) << "SF " << pdgs[j] << " [x=" << x[i] << "] = " << tmp;
        sf_stream << x[i] << " " << tmp << endl;
      }
      // Close file in which SF are stored
      sf_stream.close();
    }       
    
  }
#endif
 
}

References e, fNucPdg, genie::pdg::IsAntiNuE(), genie::pdg::IsAntiNuMu(), genie::pdg::IsAntiNuTau(), genie::pdg::IsNeutron(), genie::pdg::IsNuE(), genie::pdg::IsNuMu(), genie::pdg::IsNuTau(), genie::pdg::IsProton(), genie::constants::kGF, genie::constants::kMw, genie::kPdgAntiNuE, genie::kPdgAntiNuMu, genie::kPdgAntiNuTau, genie::kPdgNeutron, genie::kPdgNuE, genie::kPdgNuMu, genie::kPdgNuTau, genie::kPdgProton, LOG, pERROR, pFATAL, pINFO, pWARN, and xPDFmin.

Variable Documentation

fNucPdg

int fNucPdg = 0

Definition at line 50 of file gMakePhotonStrucFunc.cxx.

Referenced by main().