Functions | Variables
gEvGenHadronNucleus.cxx File Reference
#include <cassert>
#include <cstdlib>
#include <TSystem.h>
#include <TFile.h>
#include <TTree.h>
#include <TH1D.h>
#include <TF1.h>
#include "Framework/Algorithm/AlgFactory.h"
#include "Framework/Conventions/Controls.h"
#include "Framework/EventGen/EventRecord.h"
#include "Framework/EventGen/GMCJMonitor.h"
#include "Framework/EventGen/EventRecordVisitorI.h"
#include "Framework/GHEP/GHepParticle.h"
#include "Framework/GHEP/GHepRecord.h"
#include "Framework/GHEP/GHepStatus.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/PDGLibrary.h"
#include "Framework/Utils/AppInit.h"
#include "Framework/Utils/StringUtils.h"
#include "Framework/Utils/PrintUtils.h"
#include "Framework/Utils/XSecSplineList.h"
#include "Framework/Utils/RunOpt.h"
#include "Framework/Utils/CmdLnArgParser.h"
#include "Physics/HadronTransport/INukeHadroFates.h"
#include "Physics/HadronTransport/INukeUtils.h"

Go to the source code of this file.

Functions

void GetCommandLineArgs (int argc, char **argv)
 
const EventRecordVisitorIGetIntranuke (void)
 
double GenProbeKineticEnergy (void)
 
EventRecordInitializeEvent (void)
 
void BuildSpectrum (void)
 
void PrintSyntax (void)
 
int main (int argc, char **argv)
 

Variables

int kDefOptNevents = 10000
 
Long_t kDefOptRunNu = 0
 
string kDefOptEvFilePrefix = "gntp.inuke"
 
string kDefOptMode = "hA"
 
string gOptMode
 
Long_t gOptRunNu
 
int gOptNevents
 
int gOptProbePdgCode
 
int gOptTgtPdgCode
 
double gOptProbeKE
 
double gOptProbeKEmin
 
double gOptProbeKEmax
 
string gOptFlux
 
string gOptEvFilePrefix
 
bool gOptUsingFlux =false
 
long int gOptRanSeed
 
TH1D * gSpectrum = 0
 

Function Documentation

void BuildSpectrum ( void  )

Definition at line 322 of file gEvGenHadronNucleus.cxx.

References ana::assert(), genie::Spline::Evaluate(), exit(), MakeMiniprodValidationCuts::f, genie::gAbortingInErr, gOptFlux, gOptProbeKEmax, gOptProbeKEmin, gOptUsingFlux, gSpectrum, MECModelEnuComparisons::i, genie::RandomGen::Instance(), genie::controls::kRjMaxIterations, LOG, getGoodRuns4SAM::n, pFATAL, pNOTICE, r(), genie::RandomGen::RndGen(), and ymax.

Referenced by ana::GenericSystematicDef< SRType >::GenericSystematicDef(), and main().

323 {
324 // Create kinetic energy spectrum from input function
325 //
326 
327  if(!gOptUsingFlux) return;
328 
329  if(gSpectrum) {
330  delete gSpectrum;
331  gSpectrum = 0;
332  }
333 
334  LOG("gevgen_hadron", pNOTICE)
335  << "Generating a flux histogram ... ";
336 
337  int flux_bins = 300;
338  int flux_entries = 1000000;
339  double ke_min = gOptProbeKEmin;
340  double ke_max = gOptProbeKEmax;
341  double dke = ke_max - ke_min;
342  assert(dke>0 && ke_min>=0.);
343 
344  // kinetic energy spectrum
345  gSpectrum = new TH1D(
346  "spectrum","hadron kinetic energy spectrum", flux_bins, ke_min, ke_max);
347 
348  // check whether the input flux is a file or a functional form
349  bool input_is_file = ! gSystem->AccessPathName(gOptFlux.c_str());
350  if(input_is_file) {
351  Spline * input_flux = new Spline(gOptFlux.c_str());
352 
353  // generate the flux hisogram from the input flux file
354  int n = 100;
355  double ke_step = (ke_max-ke_min)/(n-1);
356  double ymax = -1, ry = -1, gy = -1, ke = -1;
357  for(int i=0; i<n; i++) {
358  ke = ke_min + i*ke_step;
359  ymax = TMath::Max(ymax, input_flux->Evaluate(ke));
360  }
361  ymax *= 1.3;
362 
363  RandomGen * r = RandomGen::Instance();
364 
365  for(int ientry=0; ientry<flux_entries; ientry++) {
366  bool accept = false;
367  unsigned int iter=0;
368  while(!accept) {
369  iter++;
370  if(iter > kRjMaxIterations) {
371  LOG("gevgen_hadron", pFATAL)
372  << "Couldn't generate a flux histogram";
373  gAbortingInErr = true;
374  exit(1);
375  }
376  ke = ke_min + dke * r->RndGen().Rndm();
377  gy = ymax * r->RndGen().Rndm();
378  ry = input_flux->Evaluate(ke);
379  accept = gy < ry;
380  if(accept) gSpectrum->Fill(ke);
381  }
382  }
383  delete input_flux;
384 
385  } else {
386  // generate the flux histogram from the input functional form
387  TF1 * input_func = new TF1("input_func", gOptFlux.c_str(), ke_min, ke_max);
388  gSpectrum->FillRandom("input_func", flux_entries);
389  delete input_func;
390  }
391  TFile f("./input-hadron-flux.root","recreate");
392  gSpectrum->Write();
393  f.Close();
394 }
double gOptProbeKEmax
A numeric analysis tool class for interpolating 1-D functions.
Definition: Spline.h:47
#define pFATAL
Definition: Messenger.h:57
TH1D * gSpectrum
double Evaluate(double x) const
Definition: Spline.cxx:362
A singleton holding random number generator classes. All random number generation in GENIE should tak...
Definition: RandomGen.h:30
Double_t ymax
Definition: plot.C:25
bool gOptUsingFlux
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
string gOptFlux
TRandom3 & RndGen(void) const
rnd number generator for generic usage
Definition: RandomGen.h:81
exit(0)
assert(nhit_max >=nhit_nbins)
TRandom3 r(0)
double gOptProbeKEmin
#define pNOTICE
Definition: Messenger.h:62
static const unsigned int kRjMaxIterations
Definition: Controls.h:27
bool gAbortingInErr
Definition: Messenger.cxx:56
double GenProbeKineticEnergy ( void  )

Definition at line 316 of file gEvGenHadronNucleus.cxx.

References gOptProbeKE, gOptUsingFlux, and gSpectrum.

Referenced by InitializeEvent().

317 {
318  if(gOptUsingFlux) return gSpectrum->GetRandom(); // spectrum
319  else return gOptProbeKE; // mono-energetic
320 }
TH1D * gSpectrum
bool gOptUsingFlux
double gOptProbeKE
void GetCommandLineArgs ( int  argc,
char **  argv 
)

Definition at line 396 of file gEvGenHadronNucleus.cxx.

References genie::CmdLnArgParser::ArgAsInt(), genie::CmdLnArgParser::ArgAsLong(), genie::CmdLnArgParser::ArgAsString(), ana::assert(), exit(), genie::gAbortingInErr, gOptEvFilePrefix, gOptFlux, gOptMode, gOptNevents, gOptProbeKE, gOptProbeKEmax, gOptProbeKEmin, gOptProbePdgCode, gOptRanSeed, gOptRunNu, gOptTgtPdgCode, gOptUsingFlux, genie::RunOpt::Instance(), kDefOptEvFilePrefix, kDefOptMode, kDefOptNevents, kDefOptRunNu, kemax, kemin, LOG, genie::CmdLnArgParser::OptionExists(), plot_validation_datamc::parser, pDEBUG, pFATAL, pINFO, pNOTICE, genie::utils::print::PrintFramedMesg(), PrintSyntax(), genie::RunOpt::ReadFromCommandLine(), and genie::utils::str::Split().

Referenced by main().

397 {
398  LOG("gevgen_hadron", pINFO) << "Parsing command line arguments";
399 
400  // Common run options.
401  RunOpt::Instance()->ReadFromCommandLine(argc,argv);
402 
403  // Parse run options for this app
404 
405  CmdLnArgParser parser(argc,argv);
406 
407  // number of events
408  if( parser.OptionExists('n') ) {
409  LOG("gevgen_hadron", pINFO) << "Reading number of events to generate";
410  gOptNevents = parser.ArgAsInt('n');
411  } else {
412  LOG("gevgen_hadron", pINFO)
413  << "Unspecified number of events to generate - Using default";
415  }
416 
417  // run number
418  if( parser.OptionExists('r') ) {
419  LOG("gevgen_hadron", pINFO) << "Reading MC run number";
420  gOptRunNu = parser.ArgAsLong('r');
421  } else {
422  LOG("gevgen_hadron", pINFO) << "Unspecified run number - Using default";
424  }
425 
426  // incoming hadron PDG code
427  if( parser.OptionExists('p') ) {
428  LOG("gevgen_hadron", pINFO) << "Reading rescattering particle PDG code";
429  gOptProbePdgCode = parser.ArgAsInt('p');
430  } else {
431  LOG("gevgen_hadron", pFATAL) << "Unspecified PDG code - Exiting";
432  PrintSyntax();
433  gAbortingInErr = true;
434  exit(1);
435  }
436 
437  // target PDG code
438  if( parser.OptionExists('t') ) {
439  LOG("gevgen_hadron", pINFO) << "Reading target PDG code";
440  gOptTgtPdgCode = parser.ArgAsInt('t');
441  } else {
442  LOG("gevgen_hadron", pFATAL) << "Unspecified target PDG code - Exiting";
443  PrintSyntax();
444  gAbortingInErr = true;
445  exit(1);
446  }
447 
448  // target PDG code
449  if( parser.OptionExists('m') ) {
450  LOG("gevgen_hadron", pINFO) << "Reading mode";
451  gOptMode = parser.ArgAsString('m');
452  } else {
453  LOG("gevgen_hadron", pFATAL) << "Unspecified mode - Using default";
455  }
456 
457  // flux functional form or flux file
458  if( parser.OptionExists('f') ) {
459  LOG("gevgen_hadron", pINFO) << "Reading hadron's kinetic energy spectrum";
460  gOptFlux = parser.ArgAsString('f');
461  gOptUsingFlux = true;
462  }
463 
464  // incoming hadron kinetic energy (or kinetic energy range, if using flux)
465  if( parser.OptionExists('k') ) {
466  LOG("gevgen_hadron", pINFO) << "Reading probe kinetic energy";
467  string ke = parser.ArgAsString('k');
468  // is it just a value or a range (comma separated set of values)
469  if(ke.find(",") != string::npos) {
470  // split the comma separated list
471  vector<string> kerange = utils::str::Split(ke, ",");
472  assert(kerange.size() == 2);
473  double kemin = atof(kerange[0].c_str());
474  double kemax = atof(kerange[1].c_str());
475  assert(kemax>kemin && kemin>0);
476  gOptProbeKE = -1;
479  // if no flux was specified, generate uniformly within given range
480  if(!gOptUsingFlux) {
481  gOptFlux = "1";
482  gOptUsingFlux = true;
483  }
484  } else {
485  gOptProbeKE = atof(ke.c_str());
486  gOptProbeKEmin = -1;
487  gOptProbeKEmax = -1;
488  if(gOptUsingFlux) {
489  LOG("gevgen_hadron", pFATAL)
490  << "You specified an input flux without a kinetic energy range";
491  PrintSyntax();
492  gAbortingInErr = true;
493  exit(1);
494  }
495  }
496  } else {
497  LOG("gevgen_hadron", pFATAL) << "Unspecified kinetic energy - Exiting";
498  PrintSyntax();
499  gAbortingInErr = true;
500  exit(1);
501  }
502 
503  // event file prefix
504  if( parser.OptionExists('o') ) {
505  LOG("gevgen_hadron", pINFO) << "Reading the event filename prefix";
506  gOptEvFilePrefix = parser.ArgAsString('o');
507  } else {
508  LOG("gevgen_hadron", pDEBUG)
509  << "Will set the default event filename prefix";
511  } //-o
512 
513  // INTRANUKE mode
514  if( parser.OptionExists('m') ) {
515  LOG("gevgen_hadron", pINFO) << "Reading mode";
516  gOptMode = parser.ArgAsString('m');
517  } else {
518  LOG("gevgen_hadron", pDEBUG)
519  << "Unspecified mode - Using default";
521  }
522 
523  // random number seed
524  if( parser.OptionExists("seed") ) {
525  LOG("gevgen_hadron", pINFO) << "Reading random number seed";
526  gOptRanSeed = parser.ArgAsLong("seed");
527  } else {
528  LOG("gevgen_hadron", pINFO) << "Unspecified random number seed - Using default";
529  gOptRanSeed = -1;
530  }
531 
532 
533  LOG("gevgen_hadron", pNOTICE)
534  << "\n"
535  << utils::print::PrintFramedMesg("gevgen_hadron job configuration");
536 
537  LOG("gevgen_hadron", pNOTICE) << "MC Run Number = " << gOptRunNu;
538  LOG("gevgen_hadron", pNOTICE) << "Random number seed = " << gOptRanSeed;
539  LOG("gevgen_hadron", pNOTICE) << "Mode = " << gOptMode;
540  LOG("gevgen_hadron", pNOTICE) << "Number of events = " << gOptNevents;
541  LOG("gevgen_hadron", pNOTICE) << "Probe PDG code = " << gOptProbePdgCode;
542  LOG("gevgen_hadron", pNOTICE) << "Target PDG code = " << gOptTgtPdgCode;
543  if(gOptProbeKEmin<0 && gOptProbeKEmax<0) {
544  LOG("gevgen_hadron", pNOTICE)
545  << "Hadron input KE = " << gOptProbeKE;
546  } else {
547  LOG("gevgen_hadron", pNOTICE)
548  << "Hadron input KE range = ["
549  << gOptProbeKEmin << ", " << gOptProbeKEmax << "]";
550  }
551  if(gOptUsingFlux) {
552  LOG("gevgen_hadron", pNOTICE)
553  << "Input flux = "
554  << gOptFlux;
555  }
556 
557  LOG("gevgen_hadron", pNOTICE) << "\n";
558  LOG("gevgen_hadron", pNOTICE) << *RunOpt::Instance();
559 }
int gOptNevents
Long_t kDefOptRunNu
double gOptProbeKEmax
long int gOptRanSeed
Long_t gOptRunNu
#define pFATAL
Definition: Messenger.h:57
string kDefOptMode
const double kemax
int kDefOptNevents
bool gOptUsingFlux
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
void PrintSyntax(void)
#define pINFO
Definition: Messenger.h:63
int gOptProbePdgCode
string gOptFlux
string gOptMode
string gOptEvFilePrefix
vector< string > Split(string input, string delim)
Definition: StringUtils.cxx:42
string kDefOptEvFilePrefix
exit(0)
assert(nhit_max >=nhit_nbins)
string PrintFramedMesg(string mesg, unsigned int nl=1, const char f='*')
Definition: PrintUtils.cxx:171
double gOptProbeKEmin
Command line argument parser.
#define pNOTICE
Definition: Messenger.h:62
double gOptProbeKE
const double kemin
bool gAbortingInErr
Definition: Messenger.cxx:56
#define pDEBUG
Definition: Messenger.h:64
int gOptTgtPdgCode
const EventRecordVisitorI * GetIntranuke ( void  )

Definition at line 233 of file gEvGenHadronNucleus.cxx.

References ana::assert(), exit(), genie::gAbortingInErr, genie::AlgFactory::GetAlgorithm(), gOptMode, genie::AlgFactory::Instance(), LOG, and pFATAL.

Referenced by main().

234 {
235 // get the requested INTRANUKE module
236 
237  string sname = "";
238  string sconf = "";
239 
240  if(gOptMode.compare("hA")==0) {
241  sname = "genie::HAIntranuke";
242  sconf = "Default";
243  }
244  else
245  if(gOptMode.compare("hN")==0) {
246  sname = "genie::HNIntranuke";
247  sconf = "Default";
248  }
249  else
250  if(gOptMode.compare("hA2019")==0) {
251  sname = "genie::HAIntranuke2019";
252  sconf = "Default";
253  }
254  else
255  if(gOptMode.compare("hN2019")==0) {
256  sname = "genie::HNIntranuke2019";
257  sconf = "Default";
258  }
259  else
260  if(gOptMode.compare("hA2018")==0) {
261  sname = "genie::HAIntranuke2018";
262  sconf = "Default";
263  }
264  else
265  if(gOptMode.compare("hN2018")==0) {
266  sname = "genie::HNIntranuke2018";
267  sconf = "Default";
268  }
269  else {
270  LOG("gevgen_hadron", pFATAL) << "Invalid Intranuke mode - Exiting";
271  gAbortingInErr = true;
272  exit(1);
273  }
274 
275  AlgFactory * algf = AlgFactory::Instance();
276  const EventRecordVisitorI * intranuke =
277  dynamic_cast<const EventRecordVisitorI *> (algf->GetAlgorithm(sname,sconf));
278  assert(intranuke);
279 
280  return intranuke;
281 }
Defines the EventRecordVisitorI interface. Concrete implementations of this interface use the &#39;Visito...
#define pFATAL
Definition: Messenger.h:57
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
const Algorithm * GetAlgorithm(const AlgId &algid)
Definition: AlgFactory.cxx:86
string gOptMode
exit(0)
assert(nhit_max >=nhit_nbins)
The GENIE Algorithm Factory.
Definition: AlgFactory.h:40
bool gAbortingInErr
Definition: Messenger.cxx:56
EventRecord * InitializeEvent ( void  )

Definition at line 283 of file gEvGenHadronNucleus.cxx.

References genie::GHepRecord::AddParticle(), genie::GHepRecord::AttachSummary(), om::Icons::Find(), GenProbeKineticEnergy(), gOptProbePdgCode, gOptTgtPdgCode, genie::PDGLibrary::Instance(), genie::kIStInitialState, genie::utils::res::Mass(), and ana::Sqrt().

Referenced by main().

284 {
285 // Initialize event record. Inserting the probe and target particles.
286 
287  EventRecord * evrec = new EventRecord();
288  Interaction * interaction = new Interaction;
289  evrec->AttachSummary(interaction);
290 
291  // dummy vertex position
292  TLorentzVector x4null(0.,0.,0.,0.);
293 
294  // incident hadron & target nucleon masses
295  PDGLibrary * pdglib = PDGLibrary::Instance();
296  double mh = pdglib -> Find (gOptProbePdgCode) -> Mass();
297  double M = pdglib -> Find (gOptTgtPdgCode ) -> Mass();
298 
299  // incident hadron kinetic energy
300  double ke = GenProbeKineticEnergy();
301 
302  // form incident hadron and target 4-momenta
303  double Eh = mh + ke;
304  double pzh = TMath::Sqrt(TMath::Max(0.,Eh*Eh-mh*mh));
305  TLorentzVector p4h (0.,0.,pzh,Eh);
306  TLorentzVector p4tgt (0.,0.,0., M);
307 
308  // insert probe and target entries
310  evrec->AddParticle(gOptProbePdgCode, ist, -1,-1,-1,-1, p4h, x4null);
311  evrec->AddParticle(gOptTgtPdgCode, ist, -1,-1,-1,-1, p4tgt, x4null);
312 
313  return evrec;
314 }
enum genie::EGHepStatus GHepStatus_t
double GenProbeKineticEnergy(void)
double Mass(Resonance_t res)
resonance mass (GeV)
virtual void AttachSummary(Interaction *interaction)
Definition: GHepRecord.cxx:143
Summary information for an interaction.
Definition: Interaction.h:56
const char * Find(const char *fname)
Definition: Icons.cxx:12
int gOptProbePdgCode
Generated Event Record. It is a GHepRecord object that can accept / be visited by EventRecordVisitorI...
Definition: EventRecord.h:38
Singleton class to load & serve a TDatabasePDG.
Definition: PDGLibrary.h:30
Var Sqrt(const Var &v)
Use to take sqrt of a var.
Definition: Var.cxx:326
virtual void AddParticle(const GHepParticle &p)
Definition: GHepRecord.cxx:535
int gOptTgtPdgCode
int main ( int  argc,
char **  argv 
)

Definition at line 157 of file gEvGenHadronNucleus.cxx.

References genie::NtpWriter::AddEventRecord(), BuildSpectrum(), genie::RunOpt::BuildTune(), genie::NtpWriter::CustomizeFilenamePrefix(), exit(), GetCommandLineArgs(), GetIntranuke(), gOptEvFilePrefix, gOptNevents, gOptRanSeed, gOptRunNu, gSpectrum, genie::NtpWriter::Initialize(), InitializeEvent(), genie::RunOpt::Instance(), genie::kNFGHEP, LOG, genie::utils::app_init::MesgThresholds(), pFATAL, pNOTICE, genie::EventRecordVisitorI::ProcessEventRecord(), genie::utils::app_init::RandGen(), genie::NtpWriter::Save(), genie::GHepRecord::SetPrintLevel(), and genie::GMCJMonitor::Update().

158 {
159  // Parse command line arguments
160  GetCommandLineArgs(argc,argv);
161 
162  if ( ! RunOpt::Instance()->Tune() ) {
163  LOG("gmkspl", pFATAL) << " No TuneId in RunOption";
164  exit(-1);
165  }
166  RunOpt::Instance()->BuildTune();
167 
168  // Init random number generator generator with user-specified seed number,
169  // set user-specified mesg thresholds, set user-specified GHEP print-level
170  utils::app_init::MesgThresholds(RunOpt::Instance()->MesgThresholdFiles());
172  GHepRecord::SetPrintLevel(RunOpt::Instance()->EventRecordPrintLevel());
173 
174  // Build the incident hadron kinetic energy spectrum, if required
175  BuildSpectrum();
176 
177  LOG("gevgen_hadron",pNOTICE) << "finish setup";
178 
179  // Get the specified INTRANUKE model
180  const EventRecordVisitorI * intranuke = GetIntranuke();
181 
182  // Initialize an Ntuple Writer to save GHEP records into a ROOT tree
183  NtpWriter ntpw(kNFGHEP, gOptRunNu);
184  ntpw.CustomizeFilenamePrefix(gOptEvFilePrefix);
185  ntpw.Initialize();
186 
187  // Create an MC job monitor
188  GMCJMonitor mcjmonitor(gOptRunNu);
189 
190  LOG("gevgen_hadron",pNOTICE) << "ready to generate events";
191 
192  //
193  // Generate events
194  //
195 
196  int ievent = 0;
197  while (ievent < gOptNevents) {
198  LOG("gevgen_hadron", pNOTICE)
199  << " *** Generating event............ " << ievent;
200 
201  // initialize
202  EventRecord * evrec = InitializeEvent();
203 
204  // generate full h+A event
205  intranuke->ProcessEventRecord(evrec);
206 
207  // print generated events
208  LOG("gevgen_hadron", pNOTICE ) << *evrec;
209 
210  // add event at the output ntuple
211  ntpw.AddEventRecord(ievent, evrec);
212 
213  // refresh the mc job monitor
214  mcjmonitor.Update(ievent,evrec);
215 
216  ievent++;
217  delete evrec;
218 
219  } // end loop events
220 
221  // Save the generated MC events
222  ntpw.Save();
223 
224  // Clean-up
225  if(gSpectrum) {
226  delete gSpectrum;
227  gSpectrum = 0;
228  }
229 
230  return 0;
231 }
void RandGen(long int seed)
Definition: AppInit.cxx:31
int gOptNevents
void GetCommandLineArgs(int argc, char **argv)
virtual void ProcessEventRecord(GHepRecord *event_rec) const =0
Defines the EventRecordVisitorI interface. Concrete implementations of this interface use the &#39;Visito...
long int gOptRanSeed
Long_t gOptRunNu
#define pFATAL
Definition: Messenger.h:57
TH1D * gSpectrum
EventRecord * InitializeEvent(void)
void BuildSpectrum(void)
Simple class to create & update MC job status files and env. vars. This is used to be able to keep tr...
Definition: GMCJMonitor.h:30
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
Generated Event Record. It is a GHepRecord object that can accept / be visited by EventRecordVisitorI...
Definition: EventRecord.h:38
string gOptEvFilePrefix
A utility class to facilitate creating the GENIE MC Ntuple from the output GENIE GHEP event records...
Definition: NtpWriter.h:40
exit(0)
void MesgThresholds(string inpfile)
Definition: AppInit.cxx:100
#define pNOTICE
Definition: Messenger.h:62
const EventRecordVisitorI * GetIntranuke(void)
void PrintSyntax ( void  )

Definition at line 561 of file gEvGenHadronNucleus.cxx.

References LOG, and pNOTICE.

Referenced by GetCommandLineArgs().

562 {
563  LOG("gevgen_hadron", pNOTICE)
564  << "\n\n"
565  << "Syntax:" << "\n"
566  << " gevgen_hadron [-r run] [-n nev] -p hadron_pdg -t tgt_pdg -k KE [-m mode] "
567  << " [-f flux] "
568  << " [--seed random_number_seed]"
569  << " [--message-thresholds xml_file]"
570  << " [--event-record-print-level level]"
571  << " [--mc-job-status-refresh-rate rate]"
572  << "\n";
573 }
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
#define pNOTICE
Definition: Messenger.h:62

Variable Documentation

string gOptEvFilePrefix

Definition at line 150 of file gEvGenHadronNucleus.cxx.

Referenced by GetCommandLineArgs(), and main().

string gOptFlux

Definition at line 149 of file gEvGenHadronNucleus.cxx.

Referenced by BuildSpectrum(), and GetCommandLineArgs().

string gOptMode

Definition at line 141 of file gEvGenHadronNucleus.cxx.

Referenced by GetCommandLineArgs(), and GetIntranuke().

int gOptNevents

Definition at line 143 of file gEvGenHadronNucleus.cxx.

Referenced by GetCommandLineArgs(), and main().

double gOptProbeKE

Definition at line 146 of file gEvGenHadronNucleus.cxx.

Referenced by GenProbeKineticEnergy(), and GetCommandLineArgs().

double gOptProbeKEmax

Definition at line 148 of file gEvGenHadronNucleus.cxx.

Referenced by BuildSpectrum(), and GetCommandLineArgs().

double gOptProbeKEmin

Definition at line 147 of file gEvGenHadronNucleus.cxx.

Referenced by BuildSpectrum(), and GetCommandLineArgs().

int gOptProbePdgCode

Definition at line 144 of file gEvGenHadronNucleus.cxx.

Referenced by GetCommandLineArgs(), and InitializeEvent().

long int gOptRanSeed

Definition at line 152 of file gEvGenHadronNucleus.cxx.

Referenced by GetCommandLineArgs(), and main().

Long_t gOptRunNu

Definition at line 142 of file gEvGenHadronNucleus.cxx.

Referenced by GetCommandLineArgs(), and main().

int gOptTgtPdgCode

Definition at line 145 of file gEvGenHadronNucleus.cxx.

Referenced by GetCommandLineArgs(), and InitializeEvent().

bool gOptUsingFlux =false
TH1D* gSpectrum = 0

Definition at line 154 of file gEvGenHadronNucleus.cxx.

Referenced by BuildSpectrum(), GenProbeKineticEnergy(), and main().

string kDefOptEvFilePrefix = "gntp.inuke"

Definition at line 137 of file gEvGenHadronNucleus.cxx.

Referenced by GetCommandLineArgs().

string kDefOptMode = "hA"

Definition at line 138 of file gEvGenHadronNucleus.cxx.

Referenced by GetCommandLineArgs().

int kDefOptNevents = 10000

Definition at line 135 of file gEvGenHadronNucleus.cxx.

Referenced by GetCommandLineArgs().

Long_t kDefOptRunNu = 0

Definition at line 136 of file gEvGenHadronNucleus.cxx.

Referenced by GetCommandLineArgs().