gEvGen.cxx File Reference

#include <cstdlib>
#include <cassert>
#include <sstream>
#include <string>
#include <vector>
#include <map>
#include <TFile.h>
#include <TTree.h>
#include <TSystem.h>
#include <TVector3.h>
#include <TH1.h>
#include <TF1.h>
#include "Framework/Conventions/XmlParserStatus.h"
#include "Framework/Conventions/GBuild.h"
#include "Framework/Conventions/Controls.h"
#include "Framework/EventGen/EventRecord.h"
#include "Framework/EventGen/GFluxI.h"
#include "Framework/EventGen/GEVGDriver.h"
#include "Framework/EventGen/GMCJDriver.h"
#include "Framework/EventGen/GMCJMonitor.h"
#include "Framework/Interaction/Interaction.h"
#include "Framework/Messenger/Messenger.h"
#include "Framework/Ntuple/NtpWriter.h"
#include "Framework/Ntuple/NtpMCFormat.h"
#include "Framework/Numerical/RandomGen.h"
#include "Framework/Numerical/Spline.h"
#include "Framework/ParticleData/PDGCodes.h"
#include "Framework/ParticleData/PDGUtils.h"
#include "Framework/Utils/AppInit.h"
#include "Framework/Utils/RunOpt.h"
#include "Framework/Utils/XSecSplineList.h"
#include "Framework/Utils/StringUtils.h"
#include "Framework/Utils/PrintUtils.h"
#include "Framework/Utils/SystemUtils.h"
#include "Framework/Utils/CmdLnArgParser.h"

Include dependency graph for gEvGen.cxx:

Go to the source code of this file.

Functions
void	GetCommandLineArgs (int argc, char **argv)
void	Initialize (void)
void	PrintSyntax (void)
void	GenerateEventsAtFixedInitState (void)
int	main (int argc, char **argv)

Variables
int	kDefOptNevents = 0
NtpMCFormat_t	kDefOptNtpFormat = kNFGHEP
Long_t	kDefOptRunNu = 0
int	gOptNevents
double	gOptNuEnergy
double	gOptNuEnergyRange
int	gOptNuPdgCode
map< int, double >	gOptTgtMix
Long_t	gOptRunNu
string	gOptFlux
string	gOptFluxFactors
bool	gOptWeighted
bool	gOptForceInt
bool	gOptUsingFluxOrTgtMix = false
long int	gOptRanSeed
string	gOptInpXSecFile
string	gOptOutFileName
string	gOptStatFileName

Function Documentation

GenerateEventsAtFixedInitState()

void GenerateEventsAtFixedInitState

(

void

)

Definition at line 293 of file gEvGen.cxx.

{
  int neutrino = gOptNuPdgCode;
  int target   = gOptTgtMix.begin()->first;
  double Ev    = gOptNuEnergy;
  TLorentzVector nu_p4(0.,0.,Ev,Ev); // px,py,pz,E (GeV)
 
  // Create init state
  InitialState init_state(target, neutrino);
 
  // Create/config event generation driver
  GEVGDriver evg_driver;
  evg_driver.SetEventGeneratorList(RunOpt::Instance()->EventGeneratorList());
  evg_driver.SetUnphysEventMask(*RunOpt::Instance()->UnphysEventMask());
  evg_driver.Configure(init_state);
 
  // Initialize an Ntuple Writer
  NtpWriter ntpw(kDefOptNtpFormat, gOptRunNu, gOptRanSeed);
 
  // If an output file name has been specified... use it
  if (!gOptOutFileName.empty()){
    ntpw.CustomizeFilename(gOptOutFileName);
  }
  ntpw.Initialize();
 
 
  // Create an MC Job Monitor
  GMCJMonitor mcjmonitor(gOptRunNu);
  mcjmonitor.SetRefreshRate(RunOpt::Instance()->MCJobStatusRefreshRate());
 
  // If a status file name has been given... use it
  if (!gOptStatFileName.empty()){
    mcjmonitor.CustomizeFilename(gOptStatFileName);
  }
 
 
  LOG("gevgen", pNOTICE)
    << "\n ** Will generate " << gOptNevents << " events for \n"
    << init_state << " at Ev = " << Ev << " GeV";
 
  // Generate events / print the GHEP record / add it to the ntuple
  int ievent = 0;
  while (ievent < gOptNevents) {
     LOG("gevgen", pNOTICE)
        << " *** Generating event............ " << ievent;
 
     // generate a single event
     EventRecord * event = evg_driver.GenerateEvent(nu_p4);
 
     if(!event) {
        LOG("gevgen", pNOTICE)
          << "Last attempt failed. Re-trying....";
        continue;
     }
 
     LOG("gevgen", pNOTICE)
        << "Generated Event GHEP Record: " << *event;
 
     // add event at the output ntuple, refresh the mc job monitor & clean up
     ntpw.AddEventRecord(ievent, event);
     mcjmonitor.Update(ievent,event);
     ievent++;
     delete event;
  }
 
  // Save the generated MC events
  ntpw.Save();
}

References genie::NtpWriter::AddEventRecord(), genie::GEVGDriver::Configure(), genie::GMCJMonitor::CustomizeFilename(), genie::NtpWriter::CustomizeFilename(), genie::GEVGDriver::GenerateEvent(), gOptNevents, gOptNuEnergy, gOptNuPdgCode, gOptOutFileName, gOptRanSeed, gOptRunNu, gOptStatFileName, gOptTgtMix, genie::NtpWriter::Initialize(), genie::RunOpt::Instance(), kDefOptNtpFormat, LOG, pNOTICE, genie::NtpWriter::Save(), genie::GEVGDriver::SetEventGeneratorList(), genie::GMCJMonitor::SetRefreshRate(), genie::GEVGDriver::SetUnphysEventMask(), and genie::GMCJMonitor::Update().

Referenced by main().

GetCommandLineArgs()

void GetCommandLineArgs	(	int	argc,
		char **	argv )

Definition at line 644 of file gEvGen.cxx.

{
  LOG("gevgen", pINFO) << "Parsing command line arguments";
 
  // Common run options. Set defaults and read.
  RunOpt::Instance()->EnableBareXSecPreCalc(true);
  RunOpt::Instance()->ReadFromCommandLine(argc,argv);
 
  // Parse run options for this app
 
  CmdLnArgParser parser(argc,argv);
 
  // help?
  bool help = parser.OptionExists('h');
  if(help) {
      PrintSyntax();
      exit(0);
  }
 
  // number of events
  if( parser.OptionExists('n') ) {
    LOG("gevgen", pINFO) << "Reading number of events to generate";
    gOptNevents = parser.ArgAsInt('n');
  } else {
    LOG("gevgen", pINFO)
       << "Unspecified number of events to generate - Using default";
    gOptNevents = kDefOptNevents;
  }
 
  // run number
  if( parser.OptionExists('r') ) {
    LOG("gevgen", pINFO) << "Reading MC run number";
    gOptRunNu = parser.ArgAsLong('r');
  } else {
    LOG("gevgen", pINFO) << "Unspecified run number - Using default";
    gOptRunNu = kDefOptRunNu;
  }
 
  // Output file name
  if( parser.OptionExists('o') ) {
    LOG("gevgen", pINFO) << "Reading output file name";
    gOptOutFileName = parser.ArgAsString('o');
 
    gOptStatFileName = gOptOutFileName;
    // strip the output file format and replace with .status
    if (gOptOutFileName.find_last_of(".") != string::npos)
      gOptStatFileName =
        gOptStatFileName.substr(0, gOptOutFileName.find_last_of("."));
    gOptStatFileName .append(".status");
  }
 
  // flux functional form
  bool using_flux = false;
  if( parser.OptionExists('f') ) {
    LOG("gevgen", pINFO) << "Reading flux function";
    gOptFlux = parser.ArgAsString('f');
    using_flux = true;
  }
  if (parser.OptionExists('F') ) {
    LOG("gevgen", pINFO) << "Reading flux factors";
    gOptFluxFactors = parser.ArgAsString('F');
  }
 
  if(parser.OptionExists('s')) {
    LOG("gevgen", pWARN)
      << "-s option no longer available. Please read the revised code documentation";
    gAbortingInErr = true;
    exit(1);
  }
 
 
  // generate weighted events option (only relevant if using a flux)
  gOptWeighted = parser.OptionExists('w');
 
  // force interaction of all injected events (only relevant if using a flux)
  gOptForceInt = parser.OptionExists("force-flux-ray-interaction");
 
  // neutrino energy
  if( parser.OptionExists('e') ) {
    LOG("gevgen", pINFO) << "Reading neutrino energy";
    string nue = parser.ArgAsString('e');
 
    // is it just a value or a range (comma separated set of values)
    if(nue.find(",") != string::npos) {
       // split the comma separated list
       vector<string> nurange = utils::str::Split(nue, ",");
       assert(nurange.size() == 2);
       double emin = atof(nurange[0].c_str());
       double emax = atof(nurange[1].c_str());
       assert(emax>emin && emin>=0);
       gOptNuEnergy      = emin;
       gOptNuEnergyRange = emax-emin;
       if(!using_flux) {
          LOG("gevgen", pWARN)
             << "No flux was specified but an energy range was input!";
          LOG("gevgen", pWARN)
             << "Events will be generated at fixed E = " << gOptNuEnergy << " GeV";
          gOptNuEnergyRange = -1;
       }
    } else {
       gOptNuEnergy       = atof(nue.c_str());
       gOptNuEnergyRange = -1;
    }
  } else {
    LOG("gevgen", pFATAL) << "Unspecified neutrino energy - Exiting";
    PrintSyntax();
    exit(1);
  }
 
  // neutrino PDG code
  if( parser.OptionExists('p') ) {
    LOG("gevgen", pINFO) << "Reading neutrino PDG code";
    gOptNuPdgCode = parser.ArgAsInt('p');
  } else {
    LOG("gevgen", pFATAL) << "Unspecified neutrino PDG code - Exiting";
    PrintSyntax();
    exit(1);
  }
 
  // target mix (their PDG codes with their corresponding weights)
  bool using_tgtmix = false;
  if( parser.OptionExists('t') ) {
    LOG("gevgen", pINFO) << "Reading target mix";
    string stgtmix = parser.ArgAsString('t');
    gOptTgtMix.clear();
    vector<string> tgtmix = utils::str::Split(stgtmix,",");
    if(tgtmix.size()==1) {
         int    pdg = atoi(tgtmix[0].c_str());
         double wgt = 1.0;
         gOptTgtMix.insert(map<int, double>::value_type(pdg, wgt));
    } else {
      using_tgtmix = true;
      vector<string>::const_iterator tgtmix_iter = tgtmix.begin();
      for( ; tgtmix_iter != tgtmix.end(); ++tgtmix_iter) {
         string tgt_with_wgt = *tgtmix_iter;
         string::size_type open_bracket  = tgt_with_wgt.find("[");
         string::size_type close_bracket = tgt_with_wgt.find("]");
         string::size_type ibeg = 0;
         string::size_type iend = open_bracket;
         string::size_type jbeg = open_bracket+1;
         string::size_type jend = close_bracket-1;
         int    pdg = atoi(tgt_with_wgt.substr(ibeg,iend).c_str());
         double wgt = atof(tgt_with_wgt.substr(jbeg,jend).c_str());
         LOG("Main", pNOTICE)
            << "Adding to target mix: pdg = " << pdg << ", wgt = " << wgt;
         gOptTgtMix.insert(map<int, double>::value_type(pdg, wgt));
      }//tgtmix_iter
    }//>1
 
  } else {
    LOG("gevgen", pFATAL) << "Unspecified target PDG code - Exiting";
    PrintSyntax();
    exit(1);
  }
 
  gOptUsingFluxOrTgtMix = using_flux || using_tgtmix;
 
  // random number seed
  if( parser.OptionExists("seed") ) {
    LOG("gevgen", pINFO) << "Reading random number seed";
    gOptRanSeed = parser.ArgAsLong("seed");
  } else {
    LOG("gevgen", pINFO) << "Unspecified random number seed - Using default";
    gOptRanSeed = -1;
  }
 
  // input cross-section file
  if( parser.OptionExists("cross-sections") ) {
    LOG("gevgen", pINFO) << "Reading cross-section file";
    gOptInpXSecFile = parser.ArgAsString("cross-sections");
  } else {
    LOG("gevgen", pINFO) << "Unspecified cross-section file";
    gOptInpXSecFile = "";
  }
 
  //
  // print-out the command line options
  //
  LOG("gevgen", pNOTICE)
     << "\n"
     << utils::print::PrintFramedMesg("gevgen job configuration");
  LOG("gevgen", pNOTICE)
     << "MC Run Number: " << gOptRunNu;
  if(gOptRanSeed != -1) {
     LOG("gevgen", pNOTICE)
       << "Random number seed: " << gOptRanSeed;
  } else {
     LOG("gevgen", pNOTICE)
       << "Random number seed was not set, using default";
  }
  LOG("gevgen", pNOTICE)
       << "Number of events requested: " << gOptNevents;
  if(gOptInpXSecFile.size() > 0) {
     LOG("gevgen", pNOTICE)
       << "Using cross-section splines read from: " << gOptInpXSecFile;
  } else {
     LOG("gevgen", pNOTICE)
       << "No input cross-section spline file";
  }
  LOG("gevgen", pNOTICE)
       << "Flux: " << gOptFlux << " factors " << gOptFluxFactors;
  LOG("gevgen", pNOTICE)
       << "Generate weighted events? " << gOptWeighted;
  LOG("gevgen", pNOTICE)
       << "Force interaction of all flux rays? " << gOptForceInt;
  if(gOptNuEnergyRange>0) {
     LOG("gevgen", pNOTICE)
        << "Neutrino energy: ["
        << gOptNuEnergy << ", " << gOptNuEnergy+gOptNuEnergyRange << "]";
  } else {
     LOG("gevgen", pNOTICE)
        << "Neutrino energy: " << gOptNuEnergy;
  }
  LOG("gevgen", pNOTICE)
      << "Neutrino code (PDG): " << gOptNuPdgCode;
  LOG("gevgen", pNOTICE)
      << "Target code (PDG) & weight fraction (in case of multiple targets): ";
  map<int,double>::const_iterator iter;
  for(iter = gOptTgtMix.begin(); iter != gOptTgtMix.end(); ++iter) {
      int    tgtpdgc = iter->first;
      double wgt     = iter->second;
      LOG("gevgen", pNOTICE)
          << " >> " <<  tgtpdgc << " (weight fraction = " << wgt << ")";
  }
  LOG("gevgen", pNOTICE) << "\n";
 
  LOG("gevgen", pNOTICE) << *RunOpt::Instance();
 
}

References genie::CmdLnArgParser::ArgAsInt(), genie::CmdLnArgParser::ArgAsLong(), genie::CmdLnArgParser::ArgAsString(), genie::RunOpt::EnableBareXSecPreCalc(), genie::gAbortingInErr, gOptFlux, gOptFluxFactors, gOptForceInt, gOptInpXSecFile, gOptNevents, gOptNuEnergy, gOptNuEnergyRange, gOptNuPdgCode, gOptOutFileName, gOptRanSeed, gOptRunNu, gOptStatFileName, gOptTgtMix, gOptUsingFluxOrTgtMix, gOptWeighted, genie::RunOpt::Instance(), kDefOptNevents, kDefOptRunNu, LOG, genie::CmdLnArgParser::OptionExists(), pFATAL, pINFO, pNOTICE, genie::utils::print::PrintFramedMesg(), PrintSyntax(), pWARN, genie::RunOpt::ReadFromCommandLine(), and genie::utils::str::Split().

Referenced by main().

Initialize()

void Initialize

(

void

)

Definition at line 273 of file gEvGen.cxx.

{
 
  if ( ! RunOpt::Instance()->Tune() ) {
    LOG("gevgen", pFATAL) << " No TuneId in RunOption";
    exit(-1);
  }
  RunOpt::Instance()->BuildTune();
 
  // Initialization of random number generators, cross-section table,
  // messenger thresholds, cache file
  utils::app_init::MesgThresholds(RunOpt::Instance()->MesgThresholdFiles());
  utils::app_init::CacheFile(RunOpt::Instance()->CacheFile());
  utils::app_init::RandGen(gOptRanSeed);
  utils::app_init::XSecTable(gOptInpXSecFile, false);
 
  // Set GHEP print level
  GHepRecord::SetPrintLevel(RunOpt::Instance()->EventRecordPrintLevel());
}

References genie::RunOpt::BuildTune(), genie::utils::app_init::CacheFile(), gOptInpXSecFile, gOptRanSeed, genie::RunOpt::Instance(), LOG, genie::utils::app_init::MesgThresholds(), pFATAL, genie::utils::app_init::RandGen(), genie::GHepRecord::SetPrintLevel(), and genie::utils::app_init::XSecTable().

Referenced by main().

main()

int main	(	int	argc,
		char **	argv )

Definition at line 245 of file gEvGen.cxx.

{
  GetCommandLineArgs(argc,argv);
  Initialize();
 
  // throw on NaNs and Infs...
#if defined(HAVE_FENV_H) && defined(HAVE_FEENABLEEXCEPT)
  feenableexcept(FE_DIVBYZERO | FE_INVALID | FE_OVERFLOW);
#endif
  //
  // Generate neutrino events
  //
 
  if(gOptUsingFluxOrTgtMix) {
#ifdef __CAN_GENERATE_EVENTS_USING_A_FLUX_OR_TGTMIX__
        GenerateEventsUsingFluxOrTgtMix();
#else
  LOG("gevgen", pERROR)
    << "\n   To be able to generate neutrino events from a flux and/or a target mix"
    << "\n   you need to add the following config options at your GENIE installation:"
    << "\n   --enable-flux-drivers  --enable-geom-drivers \n" ;
#endif
  } else {
     GenerateEventsAtFixedInitState();
  }
  return 0;
}

References GenerateEventsAtFixedInitState(), GetCommandLineArgs(), gOptUsingFluxOrTgtMix, Initialize(), LOG, and pERROR.

PrintSyntax()

void PrintSyntax

(

void

)

Definition at line 874 of file gEvGen.cxx.

{
  LOG("gevgen", pNOTICE)
    << "\n\n" << "Syntax:" << "\n"
    << "\n      gevgen [-h]"
    << "\n              [-r run#]"
    << "\n               -n nev"
    << "\n               -e energy (or energy range) "
    << "\n               -p neutrino_pdg"
    << "\n               -t target_pdg "
    << "\n              [-f flux_description]"
    << "\n              [-o outfile_name]"
    << "\n              [-w]"
    << "\n              [--force-flux-ray-interaction]"
    << "\n              [--seed random_number_seed]"
    << "\n              [--cross-sections xml_file]"
    << RunOpt::RunOptSyntaxString(true)
    << "\n";
 
}

References LOG, pNOTICE, and genie::RunOpt::RunOptSyntaxString().

Referenced by GetCommandLineArgs().

Variable Documentation

gOptFlux

string gOptFlux

Definition at line 234 of file gEvGen.cxx.

Referenced by BuildSpectrum(), and GetCommandLineArgs().

gOptFluxFactors

string gOptFluxFactors

Definition at line 235 of file gEvGen.cxx.

Referenced by GetCommandLineArgs().

gOptForceInt

bool gOptForceInt

Definition at line 237 of file gEvGen.cxx.

Referenced by GetCommandLineArgs().

gOptInpXSecFile

string gOptInpXSecFile

Definition at line 240 of file gEvGen.cxx.

gOptNevents

int gOptNevents

Definition at line 228 of file gEvGen.cxx.

Referenced by GenerateEventsAtFixedInitState(), GetCommandLineArgs(), and main().

gOptNuEnergy

double gOptNuEnergy

Definition at line 229 of file gEvGen.cxx.

Referenced by GenerateEventsAtFixedInitState(), and GetCommandLineArgs().

gOptNuEnergyRange

double gOptNuEnergyRange

Definition at line 230 of file gEvGen.cxx.

Referenced by GetCommandLineArgs().

gOptNuPdgCode

int gOptNuPdgCode

Definition at line 231 of file gEvGen.cxx.

Referenced by GenerateEventsAtFixedInitState(), and GetCommandLineArgs().

gOptOutFileName

string gOptOutFileName

Definition at line 241 of file gEvGen.cxx.

Referenced by ConvertToGHad(), ConvertToGHepMock(), ConvertToGINuke(), ConvertToGRooTracker(), ConvertToGST(), ConvertToGTracker(), ConvertToGXML(), GenerateEventsAtFixedInitState(), GetCommandLineArgs(), and RunCherryPicker().

gOptRanSeed

long int gOptRanSeed

Definition at line 239 of file gEvGen.cxx.

gOptRunNu

Long_t gOptRunNu

Definition at line 233 of file gEvGen.cxx.

gOptStatFileName

string gOptStatFileName

Definition at line 242 of file gEvGen.cxx.

Referenced by GenerateEventsAtFixedInitState(), and GetCommandLineArgs().

gOptTgtMix

map<int,double> gOptTgtMix

Definition at line 232 of file gEvGen.cxx.

gOptUsingFluxOrTgtMix

bool gOptUsingFluxOrTgtMix = false

Definition at line 238 of file gEvGen.cxx.

Referenced by GetCommandLineArgs(), and main().

gOptWeighted

bool gOptWeighted

Definition at line 236 of file gEvGen.cxx.

Referenced by GetCommandLineArgs().

kDefOptNevents

int kDefOptNevents = 0

Definition at line 223 of file gEvGen.cxx.

Referenced by GetCommandLineArgs().

kDefOptNtpFormat

NtpMCFormat_t kDefOptNtpFormat = kNFGHEP

Definition at line 224 of file gEvGen.cxx.

kDefOptRunNu

Long_t kDefOptRunNu = 0

Definition at line 225 of file gEvGen.cxx.

Referenced by GetCommandLineArgs().