Public Member Functions | Private Member Functions | Private Attributes | List of all members
genie::GMCJDriver Class Reference

A GENIE `MC Job Driver'. Can be used for setting up complicated event generation cases involving detailed flux descriptions and detector geometry descriptions. More...

#include "/cvmfs/nova.opensciencegrid.org/externals/genie/v3_00_06_p01/Linux64bit+2.6-2.12-e17-debug/GENIE-Generator/src/Framework/EventGen/GMCJDriver.h"

Public Member Functions

 GMCJDriver ()
 
 ~GMCJDriver ()
 
void SetEventGeneratorList (string listname)
 
void SetUnphysEventMask (const TBits &mask)
 
void UseFluxDriver (GFluxI *flux)
 
void UseGeomAnalyzer (GeomAnalyzerI *geom)
 
void UseSplines (bool useLogE=true)
 
bool UseMaxPathLengths (string xml_filename)
 
void KeepOnThrowingFluxNeutrinos (bool keep_on)
 
void ForceSingleProbScale (void)
 
void PreSelectEvents (bool preselect=true)
 
bool PreCalcFluxProbabilities (void)
 
bool LoadFluxProbabilities (string filename)
 
void SaveFluxProbabilities (string outfilename)
 
void Configure (bool calc_prob_scales=true)
 
EventRecordGenerateEvent (void)
 
double GlobProbScale (void) const
 
long int NFluxNeutrinos (void) const
 
map< int, double > SumFluxIntProbs (void) const
 
const GFluxIFluxDriver (void) const
 
const GeomAnalyzerIGeomAnalyzer (void) const
 
GFluxIFluxDriverPtr (void) const
 
GeomAnalyzerIGeomAnalyzerPtr (void) const
 

Private Member Functions

void InitJob (void)
 
void InitEventGeneration (void)
 
void GetParticleLists (void)
 
void GetMaxPathLengthList (void)
 
void GetMaxFluxEnergy (void)
 
void PopulateEventGenDriverPool (void)
 
void BootstrapXSecSplines (void)
 
void BootstrapXSecSplineSummation (void)
 
void ComputeProbScales (void)
 
EventRecordGenerateEvent1Try (void)
 
bool GenerateFluxNeutrino (void)
 
bool ComputePathLengths (void)
 
double ComputeInteractionProbabilities (bool use_max_path_length)
 
int SelectTargetMaterial (double R)
 
void GenerateEventKinematics (void)
 
void GenerateVertexPosition (void)
 
void ComputeEventProbability (void)
 
double InteractionProbability (double xsec, double pl, int A)
 
double PreGenFluxInteractionProbability (void)
 

Private Attributes

GEVGPoolfGPool
 A pool of GEVGDrivers properly configured event generation drivers / one per init state. More...
 
GFluxIfFluxDriver
 [input] neutrino flux driver More...
 
GeomAnalyzerIfGeomAnalyzer
 [input] detector geometry analyzer More...
 
double fEmax
 [declared by the flux driver] maximum neutrino energy More...
 
PDGCodeList fNuList
 [declared by the flux driver] list of neutrino codes More...
 
PDGCodeList fTgtList
 [declared by the geom driver] list of target codes More...
 
PathLengthList fMaxPathLengths
 [declared by the geom driver] maximum path length list More...
 
PathLengthList fCurPathLengths
 [current] path length list for current flux neutrino More...
 
TLorentzVector fCurVtx
 [current] interaction vertex More...
 
EventRecordfCurEvt
 [current] generated event More...
 
int fSelTgtPdg
 [current] selected target material PDG code More...
 
map< int, double > fCurCumulProbMap
 [current] cummulative interaction probabilities More...
 
double fNFluxNeutrinos
 [current] number of flux nuetrinos fired by the flux driver so far More...
 
map< int, TH1D * > fPmax
 [computed at init] interaction probability scale /neutrino /energy for given geometry More...
 
double fGlobPmax
 [computed at init] global interaction probability scale for given flux & geometry More...
 
string fEventGenList
 [config] list of event generators loaded by this driver (what used to be the $GEVGL setting) More...
 
TBits * fUnphysEventMask
 [config] controls whether unphysical events are returned (what used to be the $GUNPHYSMASK setting) More...
 
string fMaxPlXmlFilename
 [config] input file with max density-weighted path lengths for all materials More...
 
bool fUseExtMaxPl
 [config] using external max path length estimate? More...
 
bool fUseSplines
 [config] compute all needed & not-loaded splines at init More...
 
bool fUseLogE
 [config] build splines = f(logE) (rather than f(E)) ? More...
 
bool fKeepThrowingFluxNu
 [config] keep firing flux neutrinos till one of them interacts More...
 
bool fGenerateUnweighted
 [config] force single probability scale? More...
 
bool fPreSelect
 [config] set whether to pre-select events using max interaction paths More...
 
TFile * fFluxIntProbFile
 [input] pre-generated flux interaction probability file More...
 
TTree * fFluxIntTree
 [computed-or-loaded] pre-computed flux interaction probabilities (expected tree name is "gFlxIntProbs") More...
 
double fBrFluxIntProb
 flux interaction probability (set to branch:"FluxIntProb") More...
 
int fBrFluxIndex
 corresponding entry in flux input tree (set to address of branch:"FluxEntry") More...
 
double fBrFluxEnu
 corresponding flux P4 (set to address of branch:"FluxP4") More...
 
double fBrFluxWeight
 corresponding flux weight (set to address of branch: "FluxWeight") More...
 
int fBrFluxPDG
 corresponding flux pdg code (set to address of branch: "FluxPDG") More...
 
string fFluxIntFileName
 whether to save pre-generated flux tree for use in later jobs More...
 
string fFluxIntTreeName
 name for tree holding flux probabilities More...
 
map< int, double > fSumFluxIntProbs
 map where the key is flux pdg code and the value is sum of fBrFluxWeight * fBrFluxIntProb for all these flux neutrinos More...
 

Detailed Description

A GENIE `MC Job Driver'. Can be used for setting up complicated event generation cases involving detailed flux descriptions and detector geometry descriptions.

Author
Costas Andreopoulos <costas.andreopoulos stfc.ac.uk> University of Liverpool & STFC Rutherford Appleton Lab

May 25, 2005

Copyright (c) 2003-2019, The GENIE Collaboration For the full text of the license visit http://copyright.genie-mc.org or see $GENIE/LICENSE

Definition at line 47 of file GMCJDriver.h.

Constructor & Destructor Documentation

GMCJDriver::GMCJDriver ( )

Definition at line 134 of file GMCJDriver.cxx.

135 {
136  this->InitJob();
137 }
void InitJob(void)
Definition: GMCJDriver.cxx:494
GMCJDriver::~GMCJDriver ( )

Definition at line 139 of file GMCJDriver.cxx.

References pmax.

140 {
142  if (fGPool) delete fGPool;
143 
144  map<int,TH1D*>::iterator pmax_iter = fPmax.begin();
145  for( ; pmax_iter != fPmax.end(); ++pmax_iter) {
146  TH1D * pmax = pmax_iter->second;
147  if(pmax) {
148  delete pmax; pmax = 0;
149  }
150  }
151  fPmax.clear();
152 
153  if(fFluxIntTree) delete fFluxIntTree;
155 }
TTree * fFluxIntTree
[computed-or-loaded] pre-computed flux interaction probabilities (expected tree name is "gFlxIntProbs...
Definition: GMCJDriver.h:131
const int pmax
Definition: cellShifts.C:13
TBits * fUnphysEventMask
[config] controls whether unphysical events are returned (what used to be the $GUNPHYSMASK setting) ...
Definition: GMCJDriver.h:122
GEVGPool * fGPool
A pool of GEVGDrivers properly configured event generation drivers / one per init state...
Definition: GMCJDriver.h:106
map< int, TH1D * > fPmax
[computed at init] interaction probability scale /neutrino /energy for given geometry ...
Definition: GMCJDriver.h:119
TFile * fFluxIntProbFile
[input] pre-generated flux interaction probability file
Definition: GMCJDriver.h:130

Member Function Documentation

void GMCJDriver::BootstrapXSecSplines ( void  )
private

Definition at line 642 of file GMCJDriver.cxx.

References genie::InitialState::AsString(), genie::GEVGDriver::CreateSplines(), LOG, pINFO, and pNOTICE.

Referenced by GeomAnalyzerPtr().

643 {
644 // Bootstrap cross section spline generation by the event generation drivers
645 // that handle each initial state.
646 
647  if(!fUseSplines) return;
648 
649  LOG("GMCJDriver", pNOTICE)
650  << "Asking event generation drivers to compute all needed xsec splines";
651 
652  PDGCodeList::const_iterator nuiter;
653  PDGCodeList::const_iterator tgtiter;
654  for(nuiter = fNuList.begin(); nuiter != fNuList.end(); ++nuiter){
655  for(tgtiter = fTgtList.begin(); tgtiter != fTgtList.end(); ++tgtiter) {
656  int target_pdgc = *tgtiter;
657  int neutrino_pdgc = *nuiter;
658  InitialState init_state(target_pdgc, neutrino_pdgc);
659  LOG("GMCJDriver", pINFO)
660  << "Computing all splines needed for init-state: "
661  << init_state.AsString();
662  GEVGDriver * evgdriver = fGPool->FindDriver(init_state);
663  evgdriver->CreateSplines(-1,-1,fUseLogE);
664  } // targets
665  } // neutrinos
666  LOG("GMCJDriver", pINFO) << "Finished creating cross section splines\n";
667 }
GEVGPool * fGPool
A pool of GEVGDrivers properly configured event generation drivers / one per init state...
Definition: GMCJDriver.h:106
PDGCodeList fTgtList
[declared by the geom driver] list of target codes
Definition: GMCJDriver.h:111
GEVGDriver * FindDriver(const InitialState &init) const
Definition: GEVGPool.cxx:50
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
bool fUseLogE
[config] build splines = f(logE) (rather than f(E)) ?
Definition: GMCJDriver.h:126
bool fUseSplines
[config] compute all needed & not-loaded splines at init
Definition: GMCJDriver.h:125
GENIE Event Generation Driver. A minimalist user interface object for generating neutrino interaction...
Definition: GEVGDriver.h:55
#define pINFO
Definition: Messenger.h:63
void CreateSplines(int nknots=-1, double emax=-1, bool inLogE=true)
Definition: GEVGDriver.cxx:578
#define pNOTICE
Definition: Messenger.h:62
PDGCodeList fNuList
[declared by the flux driver] list of neutrino codes
Definition: GMCJDriver.h:110
Initial State information.
Definition: InitialState.h:49
void GMCJDriver::BootstrapXSecSplineSummation ( void  )
private

Definition at line 669 of file GMCJDriver.cxx.

References ana::assert(), genie::GEVGDriver::CreateXSecSumSpline(), dE, LOG, cet::sqlite::max(), genie::Range1D_t::max, min(), genie::Range1D_t::min, pFATAL, pNOTICE, and genie::GEVGDriver::ValidEnergyRange().

Referenced by GeomAnalyzerPtr().

670 {
671 // Sum-up the cross section splines for all the interaction that can be
672 // simulated for each initial state
673 
674  LOG("GMCJDriver", pNOTICE)
675  << "Summing-up splines to get total cross section for each init state";
676 
677  GEVGPool::iterator diter;
678  for(diter = fGPool->begin(); diter != fGPool->end(); ++diter) {
679  string init_state = diter->first;
680  GEVGDriver * evgdriver = diter->second;
681  assert(evgdriver);
682  LOG("GMCJDriver", pNOTICE)
683  << "**** Summing xsec splines for init-state = " << init_state;
684 
685  Range1D_t rE = evgdriver->ValidEnergyRange();
686  if (fEmax>rE.max || fEmax<rE.min)
687  LOG("GMCJDriver",pFATAL)
688  << " rE (validEnergyRange) [" << rE.min << "," << rE.max << "] "
689  << " fEmax " << fEmax;
690  assert(fEmax<rE.max && fEmax>rE.min);
691 
692  // decide the energy range for the sum spline - extend the spline a little
693  // bit above the maximum beam energy (but below the maximum valid energy)
694  // to avoid the evaluation of the cubic spline around the viscinity of
695  // knots with zero y values (although the GENIE Spline object handles it)
696  double dE = fEmax/10.;
697  double min = rE.min;
698  double max = (fEmax+dE < rE.max) ? fEmax+dE : rE.max;
699  evgdriver->CreateXSecSumSpline(100,min,max,true);
700  }
701  LOG("GMCJDriver", pNOTICE)
702  << "Finished summing all interaction xsec splines per initial state";
703 }
void CreateXSecSumSpline(int nk, double Emin, double Emax, bool inlogE=true)
Definition: GEVGDriver.cxx:441
GEVGPool * fGPool
A pool of GEVGDrivers properly configured event generation drivers / one per init state...
Definition: GMCJDriver.h:106
A simple [min,max] interval for doubles.
Definition: Range1.h:43
#define pFATAL
Definition: Messenger.h:57
double dE
Range1D_t ValidEnergyRange(void) const
Definition: GEVGDriver.cxx:669
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
GENIE Event Generation Driver. A minimalist user interface object for generating neutrino interaction...
Definition: GEVGDriver.h:55
static float min(const float a, const float b, const float c)
Definition: absgeo.cxx:45
double max
Definition: Range1.h:54
assert(nhit_max >=nhit_nbins)
double min
Definition: Range1.h:53
#define pNOTICE
Definition: Messenger.h:62
T max(sqlite3 *const db, std::string const &table_name, std::string const &column_name)
Definition: statistics.h:68
double fEmax
[declared by the flux driver] maximum neutrino energy
Definition: GMCJDriver.h:109
void GMCJDriver::ComputeEventProbability ( void  )
private

Definition at line 1264 of file GMCJDriver.cxx.

References genie::units::A, ana::assert(), genie::pdg::IonPdgCodeToA(), P, genie::GHepParticle::P4(), genie::GHepParticle::Pdg(), pmax, ana::weight, and xsec.

Referenced by GeomAnalyzerPtr().

1265 {
1266 // Compute event probability for the given flux neutrino & detector geometry
1267 
1268  // get interaction cross section
1269  double xsec = fCurEvt->XSec();
1270 
1271  // get path length in detector along v direction for specified target material
1272  PathLengthList::const_iterator pliter = fCurPathLengths.find(fSelTgtPdg);
1273  double path_length = pliter->second;
1274 
1275  // get target material mass number
1277 
1278  // calculate interaction probability
1279  double P = this->InteractionProbability(xsec, path_length, A);
1280 
1281  //
1282  // get weight for selected event
1283  //
1284 
1285  GHepParticle * nu = fCurEvt->Probe();
1286  int nu_pdg = nu->Pdg();
1287  double Ev = nu->P4()->Energy();
1288 
1289  double weight = 1.0;
1290  if(!fGenerateUnweighted) {
1291  map<int,TH1D*>::const_iterator pmax_iter = fPmax.find(nu_pdg);
1292  assert(pmax_iter != fPmax.end());
1293  TH1D * pmax_hst = pmax_iter->second;
1294  assert(pmax_hst);
1295  double pmax = pmax_hst->GetBinContent(pmax_hst->FindBin(Ev));
1296  assert(pmax>0);
1297  weight = pmax/fGlobPmax;
1298  }
1299 
1300  // set probability & update weight
1301  fCurEvt->SetProbability(P);
1302  fCurEvt->SetWeight(weight * fCurEvt->Weight());
1303 }
const int pmax
Definition: cellShifts.C:13
virtual void SetWeight(double wght)
Definition: GHepRecord.h:131
const TLorentzVector * P4(void) const
Definition: GHepParticle.h:79
const Var weight
virtual void SetProbability(double prob)
Definition: GHepRecord.h:132
int IonPdgCodeToA(int pdgc)
Definition: PDGUtils.cxx:61
virtual double Weight(void) const
Definition: GHepRecord.h:125
virtual GHepParticle * Probe(void) const
Definition: GHepRecord.cxx:321
int Pdg(void) const
Definition: GHepParticle.h:64
PathLengthList fCurPathLengths
[current] path length list for current flux neutrino
Definition: GMCJDriver.h:113
#define P(a, b, c, d, e, x)
bool fGenerateUnweighted
[config] force single probability scale?
Definition: GMCJDriver.h:128
map< int, TH1D * > fPmax
[computed at init] interaction probability scale /neutrino /energy for given geometry ...
Definition: GMCJDriver.h:119
EventRecord * fCurEvt
[current] generated event
Definition: GMCJDriver.h:115
int fSelTgtPdg
[current] selected target material PDG code
Definition: GMCJDriver.h:116
double fGlobPmax
[computed at init] global interaction probability scale for given flux & geometry ...
Definition: GMCJDriver.h:120
Double_t xsec[nknots]
Definition: testXsec.C:47
double InteractionProbability(double xsec, double pl, int A)
static const double A
Definition: Units.h:82
assert(nhit_max >=nhit_nbins)
virtual double XSec(void) const
Definition: GHepRecord.h:127
STDHEP-like event record entry that can fit a particle or a nucleus.
Definition: GHepParticle.h:40
double GMCJDriver::ComputeInteractionProbabilities ( bool  use_max_path_length)
private

Definition at line 1108 of file GMCJDriver.cxx.

References genie::units::A, ana::assert(), genie::InitialState::AsString(), genie::units::cm2, genie::Spline::Evaluate(), exit(), genie::pdg::IonPdgCodeToA(), LOG, pDEBUG, pFATAL, std_candles::pl, pmax, pNOTICE, xsec, and genie::GEVGDriver::XSecSumSpline().

Referenced by GeomAnalyzerPtr().

1109 {
1110  LOG("GMCJDriver", pNOTICE)
1111  << "Computing relative interaction probabilities for each material";
1112 
1113  // current flux neutrino code & 4-p
1114  int nupdg = fFluxDriver->PdgCode();
1115  const TLorentzVector & nup4 = fFluxDriver->Momentum();
1116 
1117  fCurCumulProbMap.clear();
1118 
1119  const PathLengthList & path_length_list =
1120  (use_max_path_length) ? fMaxPathLengths : fCurPathLengths;
1121 
1122  double probsum=0;
1123  PathLengthList::const_iterator pliter;
1124 
1125  for(pliter = path_length_list.begin();
1126  pliter != path_length_list.end(); ++pliter) {
1127  int mpdg = pliter->first; // material PDG code
1128  double pl = pliter->second; // density x path-length
1129  int A = pdg::IonPdgCodeToA(mpdg);
1130  double xsec = 0.; // sum of xsecs for all modelled processes for given init state
1131  double prob = 0.; // interaction probability
1132  double probn = 0.; // normalized interaction probability
1133 
1134  // find the GEVGDriver object that is handling the current init state
1135  InitialState init_state(mpdg, nupdg);
1136  GEVGDriver * evgdriver = fGPool->FindDriver(init_state);
1137  if(!evgdriver) {
1138  LOG("GMCJDriver", pFATAL)
1139  << "\n * The MC Job driver isn't properly configured!"
1140  << "\n * No event generation driver could be found for init state: "
1141  << init_state.AsString();
1142  exit(1);
1143  }
1144  // compute the interaction xsec and probability (if path-length>0)
1145  if(pl>0.) {
1146  const Spline * totxsecspl = evgdriver->XSecSumSpline();
1147  if(!totxsecspl) {
1148  LOG("GMCJDriver", pFATAL)
1149  << "\n * The MC Job driver isn't properly configured!"
1150  << "\n * Couldn't retrieve total cross section spline for init state: "
1151  << init_state.AsString();
1152  exit(1);
1153  } else {
1154  xsec = totxsecspl->Evaluate( nup4.Energy() );
1155  }
1156  prob = this->InteractionProbability(xsec,pl,A);
1157  LOG("GMCJDriver", pDEBUG)
1158  << " (xsec, pl, A)=(" << xsec << "," << pl << "," << A << ")";
1159 
1160  // scale the interaction probability to the maximum one so as not
1161  // to have to throw few billions of flux neutrinos before getting
1162  // an interaction...
1163  double pmax = 0;
1164  if(fGenerateUnweighted) pmax = fGlobPmax;
1165  else {
1166  map<int,TH1D*>::const_iterator pmax_iter = fPmax.find(nupdg);
1167  assert(pmax_iter != fPmax.end());
1168  TH1D * pmax_hst = pmax_iter->second;
1169  assert(pmax_hst);
1170  int ie = pmax_hst->FindBin(nup4.Energy());
1171  pmax = pmax_hst->GetBinContent(ie);
1172  }
1173  assert(pmax>0);
1174  LOG("GMCJDriver", pDEBUG)
1175  << "Pmax=" << pmax;
1176  probn = prob/pmax;
1177  }
1178 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__
1179  LOG("GMCJDriver", pNOTICE)
1180  << "tgt: " << mpdg << " -> TotXSec = "
1181  << xsec/units::cm2 << " cm^2, Norm.Prob = " << 100*probn << "%";
1182 #endif
1183 
1184  probsum += probn;
1185  fCurCumulProbMap.insert(map<int,double>::value_type(mpdg,probsum));
1186  }
1187  return probsum;
1188 }
const int pmax
Definition: cellShifts.C:13
GEVGPool * fGPool
A pool of GEVGDrivers properly configured event generation drivers / one per init state...
Definition: GMCJDriver.h:106
map< int, double > fCurCumulProbMap
[current] cummulative interaction probabilities
Definition: GMCJDriver.h:117
A numeric analysis tool class for interpolating 1-D functions.
Definition: Spline.h:47
#define pFATAL
Definition: Messenger.h:57
int IonPdgCodeToA(int pdgc)
Definition: PDGUtils.cxx:61
GFluxI * fFluxDriver
[input] neutrino flux driver
Definition: GMCJDriver.h:107
const Spline * XSecSumSpline(void) const
Definition: GEVGDriver.h:88
GEVGDriver * FindDriver(const InitialState &init) const
Definition: GEVGPool.cxx:50
double Evaluate(double x) const
Definition: Spline.cxx:362
static const double cm2
Definition: Units.h:77
PathLengthList fCurPathLengths
[current] path length list for current flux neutrino
Definition: GMCJDriver.h:113
Object to be filled with the neutrino path-length, for all detector geometry materials, when starting from a position x and travelling along the direction of the neutrino 4-momentum.
bool fGenerateUnweighted
[config] force single probability scale?
Definition: GMCJDriver.h:128
map< int, TH1D * > fPmax
[computed at init] interaction probability scale /neutrino /energy for given geometry ...
Definition: GMCJDriver.h:119
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
GENIE Event Generation Driver. A minimalist user interface object for generating neutrino interaction...
Definition: GEVGDriver.h:55
double fGlobPmax
[computed at init] global interaction probability scale for given flux & geometry ...
Definition: GMCJDriver.h:120
Double_t xsec[nknots]
Definition: testXsec.C:47
double InteractionProbability(double xsec, double pl, int A)
static const double A
Definition: Units.h:82
virtual int PdgCode(void)=0
returns the flux neutrino pdg code
exit(0)
virtual const TLorentzVector & Momentum(void)=0
returns the flux neutrino 4-momentum
assert(nhit_max >=nhit_nbins)
#define pNOTICE
Definition: Messenger.h:62
PathLengthList fMaxPathLengths
[declared by the geom driver] maximum path length list
Definition: GMCJDriver.h:112
Initial State information.
Definition: InitialState.h:49
#define pDEBUG
Definition: Messenger.h:64
bool GMCJDriver::ComputePathLengths ( void  )
private

Definition at line 1079 of file GMCJDriver.cxx.

References LOG, pFATAL, and pNOTICE.

Referenced by GeomAnalyzerPtr().

1080 {
1081 // Ask the geometry driver to compute (pathLength x density x weight frac.)
1082 // for all detector materials for the neutrino generated by the flux driver
1083 // and make sure that things look ok...
1084 
1085  fCurPathLengths.clear();
1086 
1087  const TLorentzVector & nup4 = fFluxDriver -> Momentum ();
1088  const TLorentzVector & nux4 = fFluxDriver -> Position ();
1089 
1091 
1092  LOG("GMCJDriver", pNOTICE) << fCurPathLengths;
1093 
1094  if(fCurPathLengths.size() == 0) {
1095  LOG("GMCJDriver", pFATAL)
1096  << "\n *** Geometry driver error ***"
1097  << "\n Got an empty PathLengthList - No material found in geometry?";
1098  return false;
1099  }
1100 
1101  if(fCurPathLengths.AreAllZero()) {
1102  LOG("GMCJDriver", pNOTICE)
1103  << "current flux v doesn't cross any geometry material...";
1104  }
1105  return true;
1106 }
bool AreAllZero(void) const
#define pFATAL
Definition: Messenger.h:57
GFluxI * fFluxDriver
[input] neutrino flux driver
Definition: GMCJDriver.h:107
PathLengthList fCurPathLengths
[current] path length list for current flux neutrino
Definition: GMCJDriver.h:113
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
GeomAnalyzerI * fGeomAnalyzer
[input] detector geometry analyzer
Definition: GMCJDriver.h:108
virtual const PathLengthList & ComputePathLengths(const TLorentzVector &x, const TLorentzVector &p)=0
#define pNOTICE
Definition: Messenger.h:62
void GMCJDriver::ComputeProbScales ( void  )
private

Definition at line 705 of file GMCJDriver.cxx.

References genie::units::A, ana::assert(), genie::InitialState::AsString(), emax, emin, genie::Spline::Evaluate(), makeTrainCVSamples::int, genie::pdg::IonPdgCodeToA(), LOG, getGoodRuns4SAM::n, pDEBUG, pINFO, pmax, pNOTICE, pWARN, and genie::GEVGDriver::XSecSumSpline().

Referenced by GeomAnalyzerPtr().

706 {
707 // Computing interaction probability scales.
708 // To minimize the numbers or trials before choosing a neutrino+target init
709 // state for generating an event (note: each trial means selecting a flux
710 // neutrino, navigating it through the detector to compute path lengths,
711 // computing interaction probabilities for each material and so on...)
712 // a set of probability scales can be used for different neutrino species
713 // and at different energy bins.
714 // A global probability scale is also being constructed for keeping the correct
715 // proportions between differect flux neutrino species or flux neutrinos of
716 // different energies.
717 
718  LOG("GMCJDriver", pNOTICE)
719  << "Computing the max. interaction probability (probability scale)";
720 
721  // clean up global probability scale and maximum probabilties per neutrino
722  // type & energy bin
723  {
724  fGlobPmax = 0;
725  map<int,TH1D*>::iterator pmax_iter = fPmax.begin();
726  for( ; pmax_iter != fPmax.end(); ++pmax_iter) {
727  TH1D * pmax = pmax_iter->second;
728  if(pmax) {
729  delete pmax; pmax = 0;
730  }
731  }
732  fPmax.clear();
733  }
734 
735  // for maximum interaction probability vs E /for given geometry/ I will
736  // be using 300 bins up to the maximum energy for the input flux
737  // double de = fEmax/300.;//djk june 5, 2013
738  double de = fEmax/300.;//djk june 5, 2013
739  double emin = 0.0;
740  double emax = fEmax + de;
741  int n = 1 + (int) ((emax-emin)/de);
742 
743  PDGCodeList::const_iterator nuiter;
744  PDGCodeList::const_iterator tgtiter;
745 
746  // loop over all neutrino types generated by the flux driver
747  for(nuiter = fNuList.begin(); nuiter != fNuList.end(); ++nuiter) {
748  int neutrino_pdgc = *nuiter;
749  TH1D * pmax_hst = new TH1D("pmax_hst",
750  "max interaction probability vs E | geom",n,emin,emax);
751  pmax_hst->SetDirectory(0);
752 
753  // loop over energy bins
754  for(int ie = 1; ie <= pmax_hst->GetNbinsX(); ie++) {
755  double EvLow = pmax_hst->GetBinCenter(ie) - 0.5*pmax_hst->GetBinWidth(ie);
756  double EvHigh = pmax_hst->GetBinCenter(ie) + 0.5*pmax_hst->GetBinWidth(ie);
757  //double Ev = pmax_hst->GetBinCenter(ie);
758 
759  // loop over targets in input geometry, form initial state and compute
760  // the sum of maximum interaction probabilities at the current energy bin
761  //
762  for(tgtiter = fTgtList.begin(); tgtiter != fTgtList.end(); ++tgtiter) {
763  int target_pdgc = *tgtiter;
764 
765  InitialState init_state(target_pdgc, neutrino_pdgc);
766 
767  LOG("GMCJDriver", pDEBUG)
768  << "Computing Pmax for init-state: " << init_state.AsString() << " E from " << EvLow << "-" << EvHigh;
769 
770  // get the appropriate driver
771  GEVGDriver * evgdriver = fGPool->FindDriver(init_state);
772 
773  // get xsec sum over all modelled processes for given neutrino+target)
774  double sxsecLow = evgdriver->XSecSumSpline()->Evaluate(EvLow);
775  double sxsecHigh = evgdriver->XSecSumSpline()->Evaluate(EvHigh);
776 
777  // get max{path-length x density}
778  double plmax = fMaxPathLengths.PathLength(target_pdgc);
779 
780  // compute/store the max interaction probabiity (for given energy)
781  int A = pdg::IonPdgCodeToA(target_pdgc);
782  double pmaxLow = this->InteractionProbability(sxsecLow, plmax, A);
783  double pmaxHigh = this->InteractionProbability(sxsecHigh, plmax, A);
784 
785  double pmax = pmaxHigh;
786  if ( pmaxLow > pmaxHigh){
787  pmax = pmaxLow;
788  LOG("GMCJDriver", pWARN)
789  << "Lower energy neutrinos have a higher probability of interacting than those at higher energy."
790  << " pmaxLow(E=" << EvLow << ")=" << pmaxLow << " and " << " pmaxHigh(E=" << EvHigh << ")=" << pmaxHigh;
791  }
792 
793  pmax_hst->SetBinContent(ie, pmax_hst->GetBinContent(ie) + pmax);
794 
795  LOG("GMCJDriver", pDEBUG)
796  << "Pmax[" << init_state.AsString() << ", Ev from " << EvLow << "-" << EvHigh << "] = " << pmax;
797  } // targets
798 
799  pmax_hst->SetBinContent(ie, 1.2 * pmax_hst->GetBinContent(ie));
800 
801  LOG("GMCJDriver", pINFO)
802  << "Pmax[nu=" << neutrino_pdgc << ", Ev from " << EvLow << "-" << EvHigh << "] = "
803  << pmax_hst->GetBinContent(ie);
804  } // E
805 
806  fPmax.insert(map<int,TH1D*>::value_type(neutrino_pdgc,pmax_hst));
807  } // nu
808 
809  // Compute global probability scale
810  // Sum Probabilities {
811  // all neutrinos, all targets, @ max path length, @ max energy}
812  //
813  {
814  for(nuiter = fNuList.begin(); nuiter != fNuList.end(); ++nuiter) {
815  int neutrino_pdgc = *nuiter;
816  map<int,TH1D*>::const_iterator pmax_iter = fPmax.find(neutrino_pdgc);
817  assert(pmax_iter != fPmax.end());
818  TH1D * pmax_hst = pmax_iter->second;
819  assert(pmax_hst);
820 // double pmax = pmax_hst->GetBinContent(pmax_hst->FindBin(fEmax));
821  double pmax = pmax_hst->GetMaximum();
822  assert(pmax>0);
823 // fGlobPmax += pmax;
824  fGlobPmax = TMath::Max(pmax, fGlobPmax); // ?;
825  }
826  LOG("GMCJDriver", pNOTICE) << "*** Probability scale = " << fGlobPmax;
827  }
828 }
const int pmax
Definition: cellShifts.C:13
double PathLength(int pdgc) const
GEVGPool * fGPool
A pool of GEVGDrivers properly configured event generation drivers / one per init state...
Definition: GMCJDriver.h:106
PDGCodeList fTgtList
[declared by the geom driver] list of target codes
Definition: GMCJDriver.h:111
int IonPdgCodeToA(int pdgc)
Definition: PDGUtils.cxx:61
const Spline * XSecSumSpline(void) const
Definition: GEVGDriver.h:88
GEVGDriver * FindDriver(const InitialState &init) const
Definition: GEVGPool.cxx:50
double Evaluate(double x) const
Definition: Spline.cxx:362
const double emin
map< int, TH1D * > fPmax
[computed at init] interaction probability scale /neutrino /energy for given geometry ...
Definition: GMCJDriver.h:119
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
const double emax
GENIE Event Generation Driver. A minimalist user interface object for generating neutrino interaction...
Definition: GEVGDriver.h:55
#define pINFO
Definition: Messenger.h:63
double fGlobPmax
[computed at init] global interaction probability scale for given flux & geometry ...
Definition: GMCJDriver.h:120
#define pWARN
Definition: Messenger.h:61
double InteractionProbability(double xsec, double pl, int A)
static const double A
Definition: Units.h:82
assert(nhit_max >=nhit_nbins)
#define pNOTICE
Definition: Messenger.h:62
PathLengthList fMaxPathLengths
[declared by the geom driver] maximum path length list
Definition: GMCJDriver.h:112
PDGCodeList fNuList
[declared by the flux driver] list of neutrino codes
Definition: GMCJDriver.h:110
Initial State information.
Definition: InitialState.h:49
#define pDEBUG
Definition: Messenger.h:64
double fEmax
[declared by the flux driver] maximum neutrino energy
Definition: GMCJDriver.h:109
void GMCJDriver::Configure ( bool  calc_prob_scales = true)

Definition at line 446 of file GMCJDriver.cxx.

References LOG, pNOTICE, and genie::utils::print::PrintFramedMesg().

Referenced by GenerateEventsAtFixedInitState(), main(), and supernova::SnovaGen::SnovaGen().

447 {
448  LOG("GMCJDriver", pNOTICE)
449  << utils::print::PrintFramedMesg("Configuring GMCJDriver");
450 
451  // Get the list of neutrino types from the input flux driver and the list
452  // of target materials from the input geometry driver
453  this->GetParticleLists();
454 
455  // Ask the input GFluxI for the max. neutrino energy (to compute Pmax)
456  this->GetMaxFluxEnergy();
457 
458  // Create all possible initial states and for each one initialize,
459  // configure & store an GEVGDriver event generation driver object.
460  // Once an 'initial state' has been selected from the input flux / geom,
461  // the responsibility for generating the neutrino interaction will be
462  // delegated to one of these drivers.
464 
465  // If the user wants to use cross section splines in order to speed things
466  // up, then coordinate spline creation from all GEVGDriver objects pushed
467  // into GEVGPool. This will create all xsec splines needed for all (enabled)
468  // processes that can be simulated involving the particles in the input flux
469  // and geometry.
470  // Spline creation will be skipped for every spline that has been pre-loaded
471  // into the the XSecSplineList.
472  // Once more it is noted that computing cross section splines is a huge
473  // overhead. The user is encouraged to generate them in advance and load
474  // them into the XSecSplineList
475  this->BootstrapXSecSplines();
476 
477  // Create cross section splines describing the total interaction xsec
478  // for a given initial state (Create them by summing all xsec splines
479  // for each possible initial state)
481 
482  if(calc_prob_scales){
483  // Ask the input geometry driver to compute the max. path length for each
484  // material in the list of target materials (or load a precomputed list)
485  this->GetMaxPathLengthList();
486 
487  // Compute the max. interaction probability to scale all interaction
488  // probabilities to be computed by this driver
489  this->ComputeProbScales();
490  }
491  LOG("GMCJDriver", pNOTICE) << "Finished configuring GMCJDriver\n\n";
492 }
void PopulateEventGenDriverPool(void)
Definition: GMCJDriver.cxx:605
void GetMaxPathLengthList(void)
Definition: GMCJDriver.cxx:577
void GetMaxFluxEnergy(void)
Definition: GMCJDriver.cxx:595
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
void BootstrapXSecSplines(void)
Definition: GMCJDriver.cxx:642
void ComputeProbScales(void)
Definition: GMCJDriver.cxx:705
string PrintFramedMesg(string mesg, unsigned int nl=1, const char f='*')
Definition: PrintUtils.cxx:171
void BootstrapXSecSplineSummation(void)
Definition: GMCJDriver.cxx:669
#define pNOTICE
Definition: Messenger.h:62
void GetParticleLists(void)
Definition: GMCJDriver.cxx:560
const GFluxI& genie::GMCJDriver::FluxDriver ( void  ) const
inline

Definition at line 77 of file GMCJDriver.h.

References fFluxDriver.

77 { return *fFluxDriver; }
GFluxI * fFluxDriver
[input] neutrino flux driver
Definition: GMCJDriver.h:107
GFluxI* genie::GMCJDriver::FluxDriverPtr ( void  ) const
inline

Definition at line 79 of file GMCJDriver.h.

References fFluxDriver.

79 { return fFluxDriver; }
GFluxI * fFluxDriver
[input] neutrino flux driver
Definition: GMCJDriver.h:107
void GMCJDriver::ForceSingleProbScale ( void  )

Definition at line 219 of file GMCJDriver.cxx.

References LOG, and pNOTICE.

Referenced by GenerateEventsAtFixedInitState(), and main().

220 {
221 // Use a single probability scale. That generates unweighted events.
222 // (Note that generating unweighted event kinematics is a different thing)
223 //
224  fGenerateUnweighted = true;
225 
226  LOG("GMCJDriver", pNOTICE)
227  << "GMCJDriver will generate un-weighted events. "
228  << "Note: That does not force unweighted event kinematics!";
229 }
bool fGenerateUnweighted
[config] force single probability scale?
Definition: GMCJDriver.h:128
#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
EventRecord * GMCJDriver::GenerateEvent ( void  )

Definition at line 838 of file GMCJDriver.cxx.

References LOG, pINFO, and pNOTICE.

Referenced by GenerateEventsAtFixedInitState(), supernova::SnovaGen::GenerateNextRecord(), and main().

839 {
840  LOG("GMCJDriver", pNOTICE) << "Generating next event...";
841 
842  this->InitEventGeneration();
843 
844  while(1) {
845  bool flux_end = fFluxDriver->End();
846  if(flux_end) {
847  LOG("GMCJDriver", pNOTICE)
848  << "No more neutrinos can be thrown by the flux driver";
849  return 0;
850  }
851 
852  EventRecord * event = this->GenerateEvent1Try();
853  if(event) return event;
854 
855  if(fKeepThrowingFluxNu) {
856  LOG("GMCJDriver", pNOTICE)
857  << "Flux neutrino didn't interact - Trying the next one...";
858  continue;
859  }
860  break;
861  } // (w(1)
862 
863  LOG("GMCJDriver", pINFO) << "Returning NULL event!";
864  return 0;
865 }
void InitEventGeneration(void)
Definition: GMCJDriver.cxx:830
GFluxI * fFluxDriver
[input] neutrino flux driver
Definition: GMCJDriver.h:107
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
bool fKeepThrowingFluxNu
[config] keep firing flux neutrinos till one of them interacts
Definition: GMCJDriver.h:127
#define pINFO
Definition: Messenger.h:63
Generated Event Record. It is a GHepRecord object that can accept / be visited by EventRecordVisitorI...
Definition: EventRecord.h:38
EventRecord * GenerateEvent1Try(void)
Definition: GMCJDriver.cxx:867
virtual bool End(void)=0
true if no more flux nu&#39;s can be thrown (eg reaching end of beam sim ntuples)
#define pNOTICE
Definition: Messenger.h:62
EventRecord * GMCJDriver::GenerateEvent1Try ( void  )
private

Definition at line 867 of file GMCJDriver.cxx.

References exit(), genie::gAbortingInErr, genie::RandomGen::Instance(), genie::controls::kASmallNum, LOG, genie::utils::print::P4AsString(), pDEBUG, genie::utils::res::PdgCode(), pERROR, pFATAL, pNOTICE, pWARN, R, generate_hists::rnd, genie::RandomGen::RndEvg(), and genie::utils::print::X4AsString().

Referenced by GeomAnalyzerPtr().

868 {
869 // attempt generating a neutrino interaction by firing a single flux neutrino
870 //
872 
873  double Pno=0, Psum=0;
874  double R = rnd->RndEvg().Rndm();
875  LOG("GMCJDriver", pDEBUG) << "Rndm [0,1] = " << R;
876 
877  // Generate a neutrino using the input GFluxI & get current pdgc/p4/x4
878  bool flux_ok = this->GenerateFluxNeutrino();
879  if(!flux_ok) {
880  LOG("GMCJDriver", pERROR)
881  << "** Rejecting current flux neutrino (flux driver err)";
882  return 0;
883  }
884 
885  // Compute the interaction probabilities assuming max. path lengths
886  // and decide whether the neutrino would interact --
887  // Many flux neutrinos should be rejected here, drastically reducing
888  // the number of neutrinos that I need to propagate through the
889  // actual detector geometry (this is skipped when using
890  // pre-calculated flux interaction probabilities)
891  if(fPreSelect) {
892  LOG("GMCJDriver", pNOTICE)
893  << "Computing interaction probabilities for max. path lengths";
894 
895  Psum = this->ComputeInteractionProbabilities(true /* <- max PL*/);
896  Pno = 1-Psum;
897  LOG("GMCJDriver", pNOTICE)
898  << "The no-interaction probability (max. path lengths) is: "
899  << 100*Pno << " %";
900  if(Pno<0.) {
901  LOG("GMCJDriver", pFATAL)
902  << "Negative no-interaction probability! (P = " << 100*Pno << " %)"
903  << " Particle E=" << fFluxDriver->Momentum().E() << " type=" << fFluxDriver->PdgCode() << "Psum=" << Psum;
904  gAbortingInErr=true;
905  exit(1);
906  }
907  if(R>=1-Pno) {
908  LOG("GMCJDriver", pNOTICE)
909  << "** Rejecting current flux neutrino";
910  return 0;
911  }
912  } // preselect
913 
914  bool pl_ok = false;
915 
916 
917  // If possible use pre-generated flux neutrino interaction probabilities
918  if(fFluxIntTree){
919  Psum = this->PreGenFluxInteractionProbability();
920  }
921  // Else compute them in the usual manner
922  else {
923  // Compute (pathLength x density x weight fraction) for all materials
924  // in the input geometry, for the neutrino generated by the flux driver
925  pl_ok = this->ComputePathLengths();
926  if(!pl_ok) {
927  LOG("GMCJDriver", pERROR)
928  << "** Rejecting current flux neutrino (err computing path-lengths)";
929  return 0;
930  }
932  LOG("GMCJDriver", pNOTICE)
933  << "** Rejecting current flux neutrino (misses generation volume)";
934  return 0;
935  }
936  Psum = this->ComputeInteractionProbabilities(false /* <- actual PL */);
937  }
938 
939 
940  if(TMath::Abs(Psum) < controls::kASmallNum){
941  LOG("GMCJDriver", pNOTICE)
942  << "** Rejecting current flux neutrino (has null interaction probability)";
943  return 0;
944  }
945 
946  // Now decide whether the current neutrino interacts
947  Pno = 1-Psum;
948  LOG("GMCJDriver", pNOTICE)
949  << "The actual 'no interaction' probability is: " << 100*Pno << " %";
950  if(Pno<0.) {
951  LOG("GMCJDriver", pFATAL)
952  << "Negative no interactin probability! (P = " << 100*Pno << " %)";
953 
954  // print info about what caused the problem
955  int nupdg = fFluxDriver -> PdgCode ();
956  const TLorentzVector & nup4 = fFluxDriver -> Momentum ();
957  const TLorentzVector & nux4 = fFluxDriver -> Position ();
958 
959  LOG("GMCJDriver", pWARN)
960  << "\n [-] Problematic neutrino: "
961  << "\n |----o PDG-code : " << nupdg
962  << "\n |----o 4-momentum : " << utils::print::P4AsString(&nup4)
963  << "\n |----o 4-position : " << utils::print::X4AsString(&nux4)
964  << "\n Emax : " << fEmax;
965 
966  LOG("GMCJDriver", pWARN)
967  << "\n Problematic path lengths:" << fCurPathLengths;
968 
969  LOG("GMCJDriver", pWARN)
970  << "\n Maximum path lengths:" << fMaxPathLengths;
971 
972  exit(1);
973  }
974  if(R>=1-Pno) {
975  LOG("GMCJDriver", pNOTICE)
976  << "** Rejecting current flux neutrino";
977  return 0;
978  }
979 
980  //
981  // The flux neutrino interacts!
982  //
983 
984  // Calculate path lengths for first time and check potential mismatch if
985  // used pre-generated flux interaction probabilities
986  if(fFluxIntTree){
987  pl_ok = this->ComputePathLengths();
988  if(!pl_ok) {
989  LOG("GMCJDriver", pFATAL) << "** Cannot calculate path lenths!";
990  exit(1);
991  }
992  double Psum_curr = this->ComputeInteractionProbabilities(false /* <- actual PL */);
993  bool mismatch = TMath::Abs(Psum-Psum_curr) > controls::kASmallNum;
994  if(mismatch){
995  LOG("GMCJDriver", pFATAL) <<
996  "** Mismatch between pre-calculated and current interaction "<<
997  "probabilities!";
998  exit(1);
999  }
1000  }
1001 
1002  // Select a target material
1003  fSelTgtPdg = this->SelectTargetMaterial(R);
1004  if(fSelTgtPdg==0) {
1005  LOG("GMCJDriver", pERROR)
1006  << "** Rejecting current flux neutrino (failed to select tgt!)";
1007  return 0;
1008  }
1009 
1010  // Ask the GEVGDriver object to select and generate an interaction and
1011  // its kinematics for the selected initial state & neutrino 4-momentum
1012  this->GenerateEventKinematics();
1013  if(!fCurEvt) {
1014  LOG("GMCJDriver", pWARN)
1015  << "** Couldn't generate kinematics for selected interaction";
1016  return 0;
1017  }
1018 
1019  // Generate an 'interaction position' in the selected material (in the
1020  // detector coord system), along the direction of nup4 & set it
1021  this->GenerateVertexPosition();
1022 
1023  // Set the event probability (probability for this event to happen given
1024  // the detector setup & the selected flux neutrino)
1025  // Note for users:
1026  // The above probability is stored at GHepRecord::Probability()
1027  // For normalization purposes make sure that you take into account the
1028  // GHepRecord::Weight() -if event generation is weighted-, and
1029  // GFluxI::Weight() -if beam simulation is weighted-.
1030  this->ComputeEventProbability();
1031 
1032  return fCurEvt;
1033 }
TTree * fFluxIntTree
[computed-or-loaded] pre-computed flux interaction probabilities (expected tree name is "gFlxIntProbs...
Definition: GMCJDriver.h:131
bool AreAllZero(void) const
#define pERROR
Definition: Messenger.h:60
static RandomGen * Instance()
Access instance.
Definition: RandomGen.cxx:79
#define pFATAL
Definition: Messenger.h:57
bool fPreSelect
[config] set whether to pre-select events using max interaction paths
Definition: GMCJDriver.h:129
bool GenerateFluxNeutrino(void)
bool ComputePathLengths(void)
GFluxI * fFluxDriver
[input] neutrino flux driver
Definition: GMCJDriver.h:107
string P4AsString(const TLorentzVector *p)
Definition: PrintUtils.cxx:34
A singleton holding random number generator classes. All random number generation in GENIE should tak...
Definition: RandomGen.h:30
void ComputeEventProbability(void)
PathLengthList fCurPathLengths
[current] path length list for current flux neutrino
Definition: GMCJDriver.h:113
EventRecord * fCurEvt
[current] generated event
Definition: GMCJDriver.h:115
int PdgCode(Resonance_t res, int Q)
(resonance id, charge) -> PDG code
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
int fSelTgtPdg
[current] selected target material PDG code
Definition: GMCJDriver.h:116
double PreGenFluxInteractionProbability(void)
#define R(x)
static const double kASmallNum
Definition: Controls.h:40
TRandom3 & RndEvg(void) const
rnd number generator used by the event generation drivers
Definition: RandomGen.h:75
double ComputeInteractionProbabilities(bool use_max_path_length)
#define pWARN
Definition: Messenger.h:61
void GenerateEventKinematics(void)
int SelectTargetMaterial(double R)
void GenerateVertexPosition(void)
virtual int PdgCode(void)=0
returns the flux neutrino pdg code
exit(0)
virtual const TLorentzVector & Momentum(void)=0
returns the flux neutrino 4-momentum
#define pNOTICE
Definition: Messenger.h:62
string X4AsString(const TLorentzVector *x)
Definition: PrintUtils.cxx:64
PathLengthList fMaxPathLengths
[declared by the geom driver] maximum path length list
Definition: GMCJDriver.h:112
bool gAbortingInErr
Definition: Messenger.cxx:56
#define pDEBUG
Definition: Messenger.h:64
double fEmax
[declared by the flux driver] maximum neutrino energy
Definition: GMCJDriver.h:109
void GMCJDriver::GenerateEventKinematics ( void  )
private

Definition at line 1212 of file GMCJDriver.cxx.

References genie::InitialState::AsString(), exit(), genie::GEVGDriver::GenerateEvent(), LOG, pFATAL, pNOTICE, and genie::GEVGDriver::SetUnphysEventMask().

Referenced by GeomAnalyzerPtr().

1213 {
1214  int nupdg = fFluxDriver->PdgCode();
1215  const TLorentzVector & nup4 = fFluxDriver->Momentum();
1216 
1217  // Find the GEVGDriver object that generates interactions for the
1218  // given initial state (neutrino + target)
1219  InitialState init_state(fSelTgtPdg, nupdg);
1220  GEVGDriver * evgdriver = fGPool->FindDriver(init_state);
1221  if(!evgdriver) {
1222  LOG("GMCJDriver", pFATAL)
1223  << "No GEVGDriver object for init state: " << init_state.AsString();
1224  exit(1);
1225  }
1226 
1227  // propagate current unphysical event mask
1228  evgdriver->SetUnphysEventMask(*fUnphysEventMask);
1229 
1230  // Ask the GEVGDriver object to select and generate an interaction for
1231  // the selected initial state & neutrino 4-momentum
1232  LOG("GMCJDriver", pNOTICE)
1233  << "Asking the selected GEVGDriver object to generate an event";
1234  fCurEvt = evgdriver->GenerateEvent(nup4);
1235 }
TBits * fUnphysEventMask
[config] controls whether unphysical events are returned (what used to be the $GUNPHYSMASK setting) ...
Definition: GMCJDriver.h:122
GEVGPool * fGPool
A pool of GEVGDrivers properly configured event generation drivers / one per init state...
Definition: GMCJDriver.h:106
#define pFATAL
Definition: Messenger.h:57
GFluxI * fFluxDriver
[input] neutrino flux driver
Definition: GMCJDriver.h:107
GEVGDriver * FindDriver(const InitialState &init) const
Definition: GEVGPool.cxx:50
EventRecord * fCurEvt
[current] generated event
Definition: GMCJDriver.h:115
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
int fSelTgtPdg
[current] selected target material PDG code
Definition: GMCJDriver.h:116
GENIE Event Generation Driver. A minimalist user interface object for generating neutrino interaction...
Definition: GEVGDriver.h:55
virtual int PdgCode(void)=0
returns the flux neutrino pdg code
exit(0)
virtual const TLorentzVector & Momentum(void)=0
returns the flux neutrino 4-momentum
#define pNOTICE
Definition: Messenger.h:62
void SetUnphysEventMask(const TBits &mask)
Definition: GEVGDriver.cxx:220
Initial State information.
Definition: InitialState.h:49
EventRecord * GenerateEvent(const TLorentzVector &nu4p)
Definition: GEVGDriver.cxx:229
bool GMCJDriver::GenerateFluxNeutrino ( void  )
private

Definition at line 1035 of file GMCJDriver.cxx.

References LOG, genie::utils::print::P4AsString(), genie::utils::res::PdgCode(), pERROR, pFATAL, pNOTICE, and genie::utils::print::X4AsString().

Referenced by GeomAnalyzerPtr().

1036 {
1037 // Ask the neutrino flux driver to generate a flux neutrino and make sure
1038 // that things look ok...
1039 //
1040  LOG("GMCJDriver", pNOTICE) << "Generating a flux neutrino";
1041 
1042  bool ok = fFluxDriver->GenerateNext();
1043  if(!ok) {
1044  LOG("GMCJDriver", pERROR)
1045  << "*** The flux driver couldn't generate a flux neutrino!!";
1046  return false;
1047  }
1048 
1049  fNFluxNeutrinos++;
1050  int nupdg = fFluxDriver -> PdgCode ();
1051  const TLorentzVector & nup4 = fFluxDriver -> Momentum ();
1052  const TLorentzVector & nux4 = fFluxDriver -> Position ();
1053 
1054  LOG("GMCJDriver", pNOTICE)
1055  << "\n [-] Generated flux neutrino: "
1056  << "\n |----o PDG-code : " << nupdg
1057  << "\n |----o 4-momentum : " << utils::print::P4AsString(&nup4)
1058  << "\n |----o 4-position : " << utils::print::X4AsString(&nux4);
1059 
1060  if(nup4.Energy() > fEmax) {
1061  LOG("GMCJDriver", pFATAL)
1062  << "\n *** Flux driver error ***"
1063  << "\n Generated flux v with E = " << nup4.Energy() << " GeV"
1064  << "\n Max v energy (declared by flux driver) = " << fEmax << " GeV"
1065  << "\n My interaction probability scaling is invalidated!!";
1066  return false;
1067  }
1068  if(!fNuList.ExistsInPDGCodeList(nupdg)) {
1069  LOG("GMCJDriver", pFATAL)
1070  << "\n *** Flux driver error ***"
1071  << "\n Generated flux v with pdg = " << nupdg
1072  << "\n It does not belong to the declared list of flux neutrinos"
1073  << "\n I was not configured to handle this!!";
1074  return false;
1075  }
1076  return true;
1077 }
#define pERROR
Definition: Messenger.h:60
#define pFATAL
Definition: Messenger.h:57
GFluxI * fFluxDriver
[input] neutrino flux driver
Definition: GMCJDriver.h:107
string P4AsString(const TLorentzVector *p)
Definition: PrintUtils.cxx:34
bool ExistsInPDGCodeList(int pdg_code) const
int PdgCode(Resonance_t res, int Q)
(resonance id, charge) -> PDG code
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
double fNFluxNeutrinos
[current] number of flux nuetrinos fired by the flux driver so far
Definition: GMCJDriver.h:118
virtual bool GenerateNext(void)=0
generate the next flux neutrino (return false in err)
#define pNOTICE
Definition: Messenger.h:62
string X4AsString(const TLorentzVector *x)
Definition: PrintUtils.cxx:64
PDGCodeList fNuList
[declared by the flux driver] list of neutrino codes
Definition: GMCJDriver.h:110
double fEmax
[declared by the flux driver] maximum neutrino energy
Definition: GMCJDriver.h:109
void GMCJDriver::GenerateVertexPosition ( void  )
private

Definition at line 1237 of file GMCJDriver.cxx.

References plot_validation_datamc::c, CLHEP::dL, makeHTMLView::dt, genie::constants::kLightSpeed, LOG, genie::units::meter, make_associated_cosmic_defs::p4, pNOTICE, and genie::units::second.

Referenced by GeomAnalyzerPtr().

1238 {
1239  // Generate an 'interaction position' in the selected material, along
1240  // the direction of nup4
1241  LOG("GMCJDriver", pNOTICE)
1242  << "Asking the geometry analyzer to generate a vertex";
1243 
1244  const TLorentzVector & p4 = fFluxDriver->Momentum ();
1245  const TLorentzVector & x4 = fFluxDriver->Position ();
1246 
1247  const TVector3 & vtx = fGeomAnalyzer->GenerateVertex(x4, p4, fSelTgtPdg);
1248 
1249  TVector3 origin(x4.X(), x4.Y(), x4.Z());
1250  origin-=vtx; // computes vector dr = origin - vtx
1251 
1252  double dL = origin.Mag();
1254  double dt = dL/c;
1255 
1256  LOG("GMCJDriver", pNOTICE)
1257  << "|vtx - origin|: dL = " << dL << " m, dt = " << dt << " sec";
1258 
1259  fCurVtx.SetXYZT(vtx.x(), vtx.y(), vtx.z(), x4.T() + dt);
1260 
1262 }
static constexpr double dL
static const double second
Definition: Units.h:35
virtual const TLorentzVector & Position(void)=0
returns the flux neutrino 4-position (note: expect SI rather than physical units) ...
GFluxI * fFluxDriver
[input] neutrino flux driver
Definition: GMCJDriver.h:107
virtual const TVector3 & GenerateVertex(const TLorentzVector &x, const TLorentzVector &p, int tgtpdg)=0
Definition: Cand.cxx:23
EventRecord * fCurEvt
[current] generated event
Definition: GMCJDriver.h:115
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
int fSelTgtPdg
[current] selected target material PDG code
Definition: GMCJDriver.h:116
static const double meter
Definition: Units.h:33
static const double kLightSpeed
Definition: Constants.h:32
TLorentzVector fCurVtx
[current] interaction vertex
Definition: GMCJDriver.h:114
virtual void SetVertex(double x, double y, double z, double t)
Definition: GHepRecord.cxx:863
GeomAnalyzerI * fGeomAnalyzer
[input] detector geometry analyzer
Definition: GMCJDriver.h:108
virtual const TLorentzVector & Momentum(void)=0
returns the flux neutrino 4-momentum
#define pNOTICE
Definition: Messenger.h:62
const GeomAnalyzerI& genie::GMCJDriver::GeomAnalyzer ( void  ) const
inline

Definition at line 78 of file GMCJDriver.h.

References fGeomAnalyzer.

78 { return *fGeomAnalyzer; }
GeomAnalyzerI * fGeomAnalyzer
[input] detector geometry analyzer
Definition: GMCJDriver.h:108
GeomAnalyzerI* genie::GMCJDriver::GeomAnalyzerPtr ( void  ) const
inline
void GMCJDriver::GetMaxFluxEnergy ( void  )
private

Definition at line 595 of file GMCJDriver.cxx.

References LOG, and pNOTICE.

Referenced by GeomAnalyzerPtr().

596 {
597  LOG("GMCJDriver", pNOTICE)
598  << "Querying the flux driver for the maximum energy of flux neutrinos";
600 
601  LOG("GMCJDriver", pNOTICE)
602  << "Maximum flux neutrino energy = " << fEmax << " GeV";
603 }
virtual double MaxEnergy(void)=0
declare the max flux neutrino energy that can be generated (for init. purposes)
GFluxI * fFluxDriver
[input] neutrino flux driver
Definition: GMCJDriver.h:107
#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
double fEmax
[declared by the flux driver] maximum neutrino energy
Definition: GMCJDriver.h:109
void GMCJDriver::GetMaxPathLengthList ( void  )
private

Definition at line 577 of file GMCJDriver.cxx.

References LOG, and pNOTICE.

Referenced by GeomAnalyzerPtr().

578 {
579  if(fUseExtMaxPl) {
580  LOG("GMCJDriver", pNOTICE)
581  << "Loading external max path-length list for input geometry from "
584 
585  } else {
586  LOG("GMCJDriver", pNOTICE)
587  << "Querying the geometry driver to compute the max path-length list";
589  }
590  // Print maximum path lengths & neutrino energy
591  LOG("GMCJDriver", pNOTICE)
592  << "Maximum path length list: " << fMaxPathLengths;
593 }
bool fUseExtMaxPl
[config] using external max path length estimate?
Definition: GMCJDriver.h:124
virtual const PathLengthList & ComputeMaxPathLengths(void)=0
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
XmlParserStatus_t LoadFromXml(string filename)
string fMaxPlXmlFilename
[config] input file with max density-weighted path lengths for all materials
Definition: GMCJDriver.h:123
GeomAnalyzerI * fGeomAnalyzer
[input] detector geometry analyzer
Definition: GMCJDriver.h:108
#define pNOTICE
Definition: Messenger.h:62
PathLengthList fMaxPathLengths
[declared by the geom driver] maximum path length list
Definition: GMCJDriver.h:112
void GMCJDriver::GetParticleLists ( void  )
private

Definition at line 560 of file GMCJDriver.cxx.

References LOG, and pNOTICE.

Referenced by GeomAnalyzerPtr().

561 {
562  // Get the list of flux neutrinos from the flux driver
563  LOG("GMCJDriver", pNOTICE)
564  << "Asking the flux driver for its list of neutrinos";
566 
567  LOG("GMCJDriver", pNOTICE) << "Flux particles: " << fNuList;
568 
569  // Get the list of target materials from the geometry driver
570  LOG("GMCJDriver", pNOTICE)
571  << "Asking the geometry driver for its list of targets";
573 
574  LOG("GMCJDriver", pNOTICE) << "Target materials: " << fTgtList;
575 }
PDGCodeList fTgtList
[declared by the geom driver] list of target codes
Definition: GMCJDriver.h:111
GFluxI * fFluxDriver
[input] neutrino flux driver
Definition: GMCJDriver.h:107
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
virtual const PDGCodeList & FluxParticles(void)=0
declare list of flux neutrinos that can be generated (for init. purposes)
GeomAnalyzerI * fGeomAnalyzer
[input] detector geometry analyzer
Definition: GMCJDriver.h:108
#define pNOTICE
Definition: Messenger.h:62
PDGCodeList fNuList
[declared by the flux driver] list of neutrino codes
Definition: GMCJDriver.h:110
virtual const PDGCodeList & ListOfTargetNuclei(void)=0
double genie::GMCJDriver::GlobProbScale ( void  ) const
inline

Definition at line 72 of file GMCJDriver.h.

References fGlobPmax.

Referenced by main(), and supernova::SnovaGen::SnovaGen().

72 { return fGlobPmax; }
double fGlobPmax
[computed at init] global interaction probability scale for given flux & geometry ...
Definition: GMCJDriver.h:120
void GMCJDriver::InitEventGeneration ( void  )
private

Definition at line 830 of file GMCJDriver.cxx.

Referenced by GeomAnalyzerPtr().

831 {
832  fCurPathLengths.clear();
833  fCurEvt = 0;
834  fSelTgtPdg = 0;
835  fCurVtx.SetXYZT(0.,0.,0.,0.);
836 }
PathLengthList fCurPathLengths
[current] path length list for current flux neutrino
Definition: GMCJDriver.h:113
EventRecord * fCurEvt
[current] generated event
Definition: GMCJDriver.h:115
int fSelTgtPdg
[current] selected target material PDG code
Definition: GMCJDriver.h:116
TLorentzVector fCurVtx
[current] interaction vertex
Definition: GMCJDriver.h:114
void GMCJDriver::InitJob ( void  )
private

Definition at line 494 of file GMCJDriver.cxx.

References ana::assert(), MECModelEnuComparisons::i, genie::AlgConfigPool::Instance(), genie::Messenger::Instance(), and genie::GHepFlags::NFlags().

Referenced by GeomAnalyzerPtr().

495 {
496  fEventGenList = "Default"; // <-- set of event generators to be loaded by this driver
497 
498  fUnphysEventMask = new TBits(GHepFlags::NFlags()); //<-- unphysical event mask
499  //fUnphysEventMask->ResetAllBits(true);
500  for(unsigned int i = 0; i < GHepFlags::NFlags(); i++) {
501  fUnphysEventMask->SetBitNumber(i, true);
502  }
503 
504  fFluxDriver = 0; // <-- flux driver
505  fGeomAnalyzer = 0; // <-- geometry driver
506  fGPool = 0; // <-- pool of GEVGDriver event generation drivers
507  fEmax = 0; // <-- maximum neutrino energy
508  fMaxPlXmlFilename = ""; // <-- XML file with external path lengths
509  fUseExtMaxPl = false;
510  fUseSplines = false;
511  fNFluxNeutrinos = 0; // <-- number of flux neutrinos thrown so far
512 
513  fGlobPmax = 0; // <-- maximum interaction probability (global prob scale)
514  fPmax.clear(); // <-- maximum interaction probability per neutrino & per energy bin
515 
516  fGenerateUnweighted = false; // <-- default opt to generate weighted events
517  fPreSelect = true; // <-- default to use pre-selection based on maximum path lengths
518 
519  fSelTgtPdg = 0;
520  fCurEvt = 0;
521  fCurVtx.SetXYZT(0.,0.,0.,0.);
522 
523  fFluxIntProbFile = 0;
524  fFluxIntTreeName = "gFlxIntProb";
525  fFluxIntFileName = "";
526  fFluxIntTree = 0;
527  fBrFluxIntProb = -1.;
528  fBrFluxIndex = -1;
529  fBrFluxEnu = -1.;
530  fBrFluxWeight = -1.;
531  fBrFluxPDG = 0;
532  fSumFluxIntProbs.clear();
533 
534  // Throw as many flux neutrinos as necessary till one has interacted
535  // so that GenerateEvent() never returns NULL (except when in error)
536  this->KeepOnThrowingFluxNeutrinos(true);
537 
538  // Allow the selected GEVGDriver to go into recursive mode and regenerate
539  // an interaction that turns out to be unphysical.
540  //TBits unphysmask(GHepFlags::NFlags());
541  //unphysmask.ResetAllBits(false);
542  //this->FilterUnphysical(unphysmask);
543 
544  // Force early initialization of singleton objects that are typically
545  // would be initialized at their first use later on.
546  // This is purely cosmetic and I do it to send the banner and some prolific
547  // initialization printout at the top.
550 
551  // Clear the target and neutrino lists
552  fNuList.clear();
553  fTgtList.clear();
554 
555  // Clear the maximum path length list
556  fMaxPathLengths.clear();
557  fCurPathLengths.clear();
558 }
TTree * fFluxIntTree
[computed-or-loaded] pre-computed flux interaction probabilities (expected tree name is "gFlxIntProbs...
Definition: GMCJDriver.h:131
TBits * fUnphysEventMask
[config] controls whether unphysical events are returned (what used to be the $GUNPHYSMASK setting) ...
Definition: GMCJDriver.h:122
GEVGPool * fGPool
A pool of GEVGDrivers properly configured event generation drivers / one per init state...
Definition: GMCJDriver.h:106
double fBrFluxIntProb
flux interaction probability (set to branch:"FluxIntProb")
Definition: GMCJDriver.h:132
PDGCodeList fTgtList
[declared by the geom driver] list of target codes
Definition: GMCJDriver.h:111
bool fPreSelect
[config] set whether to pre-select events using max interaction paths
Definition: GMCJDriver.h:129
GFluxI * fFluxDriver
[input] neutrino flux driver
Definition: GMCJDriver.h:107
double fBrFluxEnu
corresponding flux P4 (set to address of branch:"FluxP4")
Definition: GMCJDriver.h:134
bool fUseExtMaxPl
[config] using external max path length estimate?
Definition: GMCJDriver.h:124
double fBrFluxWeight
corresponding flux weight (set to address of branch: "FluxWeight")
Definition: GMCJDriver.h:135
string fFluxIntFileName
whether to save pre-generated flux tree for use in later jobs
Definition: GMCJDriver.h:137
PathLengthList fCurPathLengths
[current] path length list for current flux neutrino
Definition: GMCJDriver.h:113
void KeepOnThrowingFluxNeutrinos(bool keep_on)
Definition: GMCJDriver.cxx:211
bool fGenerateUnweighted
[config] force single probability scale?
Definition: GMCJDriver.h:128
map< int, TH1D * > fPmax
[computed at init] interaction probability scale /neutrino /energy for given geometry ...
Definition: GMCJDriver.h:119
EventRecord * fCurEvt
[current] generated event
Definition: GMCJDriver.h:115
static unsigned int NFlags(void)
Definition: GHepFlags.h:77
TFile * fFluxIntProbFile
[input] pre-generated flux interaction probability file
Definition: GMCJDriver.h:130
int fSelTgtPdg
[current] selected target material PDG code
Definition: GMCJDriver.h:116
double fNFluxNeutrinos
[current] number of flux nuetrinos fired by the flux driver so far
Definition: GMCJDriver.h:118
int fBrFluxPDG
corresponding flux pdg code (set to address of branch: "FluxPDG")
Definition: GMCJDriver.h:136
bool fUseSplines
[config] compute all needed & not-loaded splines at init
Definition: GMCJDriver.h:125
static Messenger * Instance(void)
Definition: Messenger.cxx:71
int fBrFluxIndex
corresponding entry in flux input tree (set to address of branch:"FluxEntry")
Definition: GMCJDriver.h:133
string fFluxIntTreeName
name for tree holding flux probabilities
Definition: GMCJDriver.h:138
string fEventGenList
[config] list of event generators loaded by this driver (what used to be the $GEVGL setting) ...
Definition: GMCJDriver.h:121
double fGlobPmax
[computed at init] global interaction probability scale for given flux & geometry ...
Definition: GMCJDriver.h:120
TLorentzVector fCurVtx
[current] interaction vertex
Definition: GMCJDriver.h:114
string fMaxPlXmlFilename
[config] input file with max density-weighted path lengths for all materials
Definition: GMCJDriver.h:123
map< int, double > fSumFluxIntProbs
map where the key is flux pdg code and the value is sum of fBrFluxWeight * fBrFluxIntProb for all the...
Definition: GMCJDriver.h:139
GeomAnalyzerI * fGeomAnalyzer
[input] detector geometry analyzer
Definition: GMCJDriver.h:108
assert(nhit_max >=nhit_nbins)
PathLengthList fMaxPathLengths
[declared by the geom driver] maximum path length list
Definition: GMCJDriver.h:112
PDGCodeList fNuList
[declared by the flux driver] list of neutrino codes
Definition: GMCJDriver.h:110
static AlgConfigPool * Instance()
double fEmax
[declared by the flux driver] maximum neutrino energy
Definition: GMCJDriver.h:109
double GMCJDriver::InteractionProbability ( double  xsec,
double  pl,
int  A 
)
private

Definition at line 1305 of file GMCJDriver.cxx.

References genie::units::cm2, genie::units::gram, genie::units::kilogram, genie::constants::kNA, and genie::units::m2.

Referenced by GeomAnalyzerPtr().

1306 {
1307 // P = Na (Avogadro number, atoms/mole) *
1308 // 1/A (1/mass number, mole/gr) *
1309 // xsec (total interaction cross section, cm^2) *
1310 // pL (density-weighted path-length, gr/cm^2)
1311 //
1312  xsec = xsec / units::cm2;
1314 
1315  return kNA*(xsec*pL)/A;
1316 }
static const double cm2
Definition: Units.h:77
static const double kNA
Definition: Constants.h:50
static const double kilogram
Definition: Units.h:34
Double_t xsec[nknots]
Definition: testXsec.C:47
static const double gram
Definition: Units.h:141
static const double A
Definition: Units.h:82
static const double m2
Definition: Units.h:80
void GMCJDriver::KeepOnThrowingFluxNeutrinos ( bool  keep_on)

Definition at line 211 of file GMCJDriver.cxx.

References genie::utils::print::BoolAsYNString(), LOG, and pNOTICE.

212 {
213  LOG("GMCJDriver", pNOTICE)
214  << "Keep on throwing flux neutrinos till one interacts? : "
215  << utils::print::BoolAsYNString(keep_on);
216  fKeepThrowingFluxNu = keep_on;
217 }
string BoolAsYNString(bool b)
Definition: PrintUtils.cxx:115
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
bool fKeepThrowingFluxNu
[config] keep firing flux neutrinos till one of them interacts
Definition: GMCJDriver.h:127
#define pNOTICE
Definition: Messenger.h:62
bool GMCJDriver::LoadFluxProbabilities ( string  filename)

Definition at line 390 of file GMCJDriver.cxx.

References LOG, pERROR, pNOTICE, and pWARN.

Referenced by main().

391 {
392 // Load a pre-generated set of flux interaction probabilities from an external
393 // file. This is recommended when using large flux files (>100k entries) as
394 // for these the time to calculate the interaction probabilities can exceed
395 // ~20 minutes. After loading the input tree we call PreCalcFluxProbabilities
396 // to check that has successfully loaded
397 //
398  if(fFluxIntProbFile){
399  LOG("GMCJDriver", pWARN)
400  << "Can't load flux interaction prob file as one is already loaded";
401  return false;
402  }
403 
404  fFluxIntProbFile = new TFile(filename.c_str(), "OPEN");
405 
406  if(fFluxIntProbFile){
407  fFluxIntTree = dynamic_cast<TTree*>(fFluxIntProbFile->Get(fFluxIntTreeName.c_str()));
408  if(fFluxIntTree){
409  bool set_addresses =
410  fFluxIntTree->SetBranchAddress("FluxIntProb", &fBrFluxIntProb) >= 0 &&
411  fFluxIntTree->SetBranchAddress("FluxIndex", &fBrFluxIndex) >= 0 &&
412  fFluxIntTree->SetBranchAddress("FluxPDG", &fBrFluxPDG) >= 0 &&
413  fFluxIntTree->SetBranchAddress("FluxWeight", &fBrFluxWeight) >= 0 &&
414  fFluxIntTree->SetBranchAddress("FluxEnu", &fBrFluxEnu) >= 0;
415  if(set_addresses){
416  // Finally check that can use them
417  if(this->PreCalcFluxProbabilities()) {
418  LOG("GMCJDriver", pNOTICE)
419  << "Successfully loaded pre-generated flux interaction probabilities";
420  return true;
421  }
422  }
423  // If cannot load then delete tree
424  LOG("GMCJDriver", pERROR) <<
425  "Cannot find expected branches in input flux probability tree!";
426  delete fFluxIntTree; fFluxIntTree = 0;
427  }
428  else LOG("GMCJDriver", pERROR)
429  << "Cannot find tree: "<< fFluxIntTreeName.c_str();
430  }
431 
432  LOG("GMCJDriver", pWARN)
433  << "Unable to load flux interaction probabilities file";
434  return false;
435 }
TTree * fFluxIntTree
[computed-or-loaded] pre-computed flux interaction probabilities (expected tree name is "gFlxIntProbs...
Definition: GMCJDriver.h:131
bool PreCalcFluxProbabilities(void)
Definition: GMCJDriver.cxx:239
#define pERROR
Definition: Messenger.h:60
double fBrFluxIntProb
flux interaction probability (set to branch:"FluxIntProb")
Definition: GMCJDriver.h:132
string filename
Definition: shutoffs.py:106
double fBrFluxEnu
corresponding flux P4 (set to address of branch:"FluxP4")
Definition: GMCJDriver.h:134
double fBrFluxWeight
corresponding flux weight (set to address of branch: "FluxWeight")
Definition: GMCJDriver.h:135
TFile * fFluxIntProbFile
[input] pre-generated flux interaction probability file
Definition: GMCJDriver.h:130
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
int fBrFluxPDG
corresponding flux pdg code (set to address of branch: "FluxPDG")
Definition: GMCJDriver.h:136
int fBrFluxIndex
corresponding entry in flux input tree (set to address of branch:"FluxEntry")
Definition: GMCJDriver.h:133
string fFluxIntTreeName
name for tree holding flux probabilities
Definition: GMCJDriver.h:138
#define pWARN
Definition: Messenger.h:61
#define pNOTICE
Definition: Messenger.h:62
long int genie::GMCJDriver::NFluxNeutrinos ( void  ) const
inline

Definition at line 73 of file GMCJDriver.h.

References fNFluxNeutrinos.

73 { return (long int) fNFluxNeutrinos; }
double fNFluxNeutrinos
[current] number of flux nuetrinos fired by the flux driver so far
Definition: GMCJDriver.h:118
void GMCJDriver::PopulateEventGenDriverPool ( void  )
private

Definition at line 605 of file GMCJDriver.cxx.

References genie::InitialState::AsString(), genie::GEVGDriver::Configure(), LOG, pDEBUG, pNOTICE, genie::GEVGDriver::SetEventGeneratorList(), and genie::GEVGDriver::UseSplines().

Referenced by GeomAnalyzerPtr().

606 {
607  LOG("GMCJDriver", pDEBUG)
608  << "Creating GEVGPool & adding a GEVGDriver object per init-state";
609 
610  if (fGPool) delete fGPool;
611  fGPool = new GEVGPool;
612 
613  PDGCodeList::const_iterator nuiter;
614  PDGCodeList::const_iterator tgtiter;
615 
616  for(nuiter = fNuList.begin(); nuiter != fNuList.end(); ++nuiter) {
617  for(tgtiter = fTgtList.begin(); tgtiter != fTgtList.end(); ++tgtiter) {
618 
619  int target_pdgc = *tgtiter;
620  int neutrino_pdgc = *nuiter;
621 
622  InitialState init_state(target_pdgc, neutrino_pdgc);
623 
624  LOG("GMCJDriver", pNOTICE)
625  << "\n\n ---- Creating a GEVGDriver object configured for init-state: "
626  << init_state.AsString() << " ----\n\n";
627 
628  GEVGDriver * evgdriver = new GEVGDriver;
629  evgdriver->SetEventGeneratorList(fEventGenList); // specify list of generators
630  evgdriver->Configure(init_state);
631  evgdriver->UseSplines(); // check if all splines needed are loaded
632 
633  LOG("GMCJDriver", pDEBUG) << "Adding new GEVGDriver object to GEVGPool";
634  fGPool->insert( GEVGPool::value_type(init_state.AsString(), evgdriver) );
635  } // targets
636  } // neutrinos
637 
638  LOG("GMCJDriver", pNOTICE)
639  << "All necessary GEVGDriver object were pushed into GEVGPool\n";
640 }
GEVGPool * fGPool
A pool of GEVGDrivers properly configured event generation drivers / one per init state...
Definition: GMCJDriver.h:106
PDGCodeList fTgtList
[declared by the geom driver] list of target codes
Definition: GMCJDriver.h:111
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
void SetEventGeneratorList(string listname)
Definition: GEVGDriver.cxx:349
GENIE Event Generation Driver. A minimalist user interface object for generating neutrino interaction...
Definition: GEVGDriver.h:55
string fEventGenList
[config] list of event generators loaded by this driver (what used to be the $GEVGL setting) ...
Definition: GMCJDriver.h:121
void Configure(int nu_pdgc, int Z, int A)
Definition: GEVGDriver.cxx:138
void UseSplines(void)
Definition: GEVGDriver.cxx:509
#define pNOTICE
Definition: Messenger.h:62
PDGCodeList fNuList
[declared by the flux driver] list of neutrino codes
Definition: GMCJDriver.h:110
Initial State information.
Definition: InitialState.h:49
#define pDEBUG
Definition: Messenger.h:64
A pool of GEVGDriver objects with an initial state key.
Definition: GEVGPool.h:38
bool GMCJDriver::PreCalcFluxProbabilities ( void  )

Definition at line 239 of file GMCJDriver.cxx.

References ana::assert(), exit(), MECModelEnuComparisons::i, genie::controls::kASmallNum, LOG, pFATAL, pNOTICE, pWARN, and log::success().

Referenced by main().

240 {
241 // Loop over complete set of flux entries satisfying input config options
242 // (such as neutrino type) and save the interaction probability in a tree
243 // relating flux index (entry number in input flux tree) to interaction
244 // probability. If a pre-generated flux interaction probability tree has
245 // already been loaded then just returns true. Also save tree to a TFile
246 // for use in later jobs if flag is set
247 //
248  bool success = true;
249 
250  bool save_to_file = fFluxIntProbFile == 0 && fFluxIntFileName.size()>0;
251 
252  // Clear map storing sum(fBrFluxWeight*fBrFluxIntProb) for each neutrino pdg
253  fSumFluxIntProbs.clear();
254 
255  // check if already loaded flux interaction probs using LoadFluxProbTree
256  if(fFluxIntTree){
257  LOG("GMCJDriver", pNOTICE) <<
258  "Skipping pre-generation of flux interaction probabilities - "<<
259  "using pre-generated file";
260  success = true;
261  }
262  // otherwise create them on the fly now
263  else {
264 
265  if(save_to_file){
266  fFluxIntProbFile = new TFile(fFluxIntFileName.c_str(), "CREATE");
267  if(fFluxIntProbFile->IsZombie()){
268  LOG("GMCJDriver", pFATAL) << "Cannot overwrite an existing file. Exiting!";
269  exit(1);
270  }
271  }
272 
273  // Create the tree to store flux probs
274  fFluxIntTree = new TTree(fFluxIntTreeName.c_str(),
275  "Tree storing pre-calculated flux interaction probs");
276  fFluxIntTree->Branch("FluxIndex", &fBrFluxIndex, "FluxIndex/I");
277  fFluxIntTree->Branch("FluxIntProb", &fBrFluxIntProb, "FluxIntProb/D");
278  fFluxIntTree->Branch("FluxEnu", &fBrFluxEnu, "FluxEnu/D");
279  fFluxIntTree->Branch("FluxWeight", &fBrFluxWeight, "FluxWeight/D");
280  fFluxIntTree->Branch("FluxPDG", &fBrFluxPDG, "FluxPDG/I");
281  // Associate to file otherwise get std::bad_alloc when writing large trees
282  if(save_to_file) fFluxIntTree->SetDirectory(fFluxIntProbFile);
283 
285 
286  fGlobPmax = 1.0; // Force ComputeInteractionProbabilities to return absolute value
287 
288  // Loop over flux entries and calculate interaction probabilities
289  TStopwatch stopwatch;
290  stopwatch.Start();
291  long int first_index = -1;
292  bool first_loop = true;
293  // loop until at end of flux ntuple
294  while(fFluxDriver->End() == false){
295 
296  // get the next flux neutrino
297  bool gotnext = fFluxDriver->GenerateNext();
298  if(!gotnext){
299  LOG("GMCJDriver", pWARN) << "*** Couldn't generate next flux ray! ";
300  continue;
301  }
302 
303  // stop if completed a full cycle (this check is necessary as fluxdriver
304  // may be set to loop over more than one cycle before reaching end)
305  bool already_been_here = first_loop ? false : first_index == fFluxDriver->Index();
306  if(already_been_here) break;
307 
308  // compute the path lengths for current flux neutrino
309  if(this->ComputePathLengths() == false){ success = false; break;}
310 
311  // compute and store the interaction probability
312  double psum = this->ComputeInteractionProbabilities(false /*Based on actual PLs*/);
313  assert(psum+controls::kASmallNum > 0.);
314  fBrFluxIntProb = psum;
319  fFluxIntTree->Fill();
320 
321  // store the first index so know when have cycled exactly once
322  if(first_loop){
323  first_index = fFluxDriver->Index();
324  first_loop = false;
325  }
326  } // flux loop
327  stopwatch.Stop();
328  LOG("GMCJDriver", pNOTICE)
329  << "Finished pre-calculating flux interaction probabilities. "
330  << "Total CPU time to process "<< fFluxIntTree->GetEntries()
331  << " entries: "<< stopwatch.CpuTime();
332 
333  // reset the flux driver so can be used at next stage. N.B. This
334  // should also reset flux driver to throw de-weighted flux neutrinos
335  fFluxDriver->Clear("CycleHistory");
336  }
337 
338  // If successfully calculated/loaded interaction probabilities then set global
339  // probability scale and, if requested, save tree to output file
340  if(success){
341  fGlobPmax = 0.0;
342  double safety_factor = 1.01;
343  for(int i = 0; i< fFluxIntTree->GetEntries(); i++){
344  fFluxIntTree->GetEntry(i);
345  // Check have non-negative probabilities
348  // Update the global maximum
349  fGlobPmax = TMath::Max(fGlobPmax, fBrFluxIntProb*safety_factor);
350  // Update the sum of fBrFluxIntProb*fBrFluxWeight for different species
351  if(fSumFluxIntProbs.find(fBrFluxPDG) == fSumFluxIntProbs.end()){
353  }
355  }
356  LOG("GMCJDriver", pNOTICE) <<
357  "Updated global probability scale to fGlobPmax = "<< fGlobPmax;
358 
359  if(save_to_file){
360  LOG("GMCJDriver", pNOTICE) <<
361  "Saving pre-generated interaction probabilities to file: "<<
362  fFluxIntProbFile->GetName();
363  fFluxIntProbFile->cd();
364  fFluxIntTree->Write();
365  }
366 
367  // Also build index for use later
368  if(fFluxIntTree->BuildIndex("FluxIndex") != fFluxIntTree->GetEntries()){
369  LOG("GMCJDriver", pFATAL) <<
370  "Cannot build index using branch \"FluxIndex\" for flux prob tree!";
371  exit(1);
372  }
373 
374  // Now that have pre-generated flux probabilities need to trun off event
375  // preselection as this is only advantages when using max path lengths
376  this->PreSelectEvents(false);
377 
378  LOG("GMCJDriver", pNOTICE) << "Successfully generated/loaded pre-calculate flux interaction probabilities";
379  }
380  // Otherwise clean up
381  else if(fFluxIntTree){
382  delete fFluxIntTree;
383  fFluxIntTree = 0;
384  }
385 
386  // Return whether have successfully pre-calculated flux interaction probabilities
387  return success;
388 }
TTree * fFluxIntTree
[computed-or-loaded] pre-computed flux interaction probabilities (expected tree name is "gFlxIntProbs...
Definition: GMCJDriver.h:131
double fBrFluxIntProb
flux interaction probability (set to branch:"FluxIntProb")
Definition: GMCJDriver.h:132
#define pFATAL
Definition: Messenger.h:57
bool ComputePathLengths(void)
GFluxI * fFluxDriver
[input] neutrino flux driver
Definition: GMCJDriver.h:107
virtual void GenerateWeighted(bool gen_weighted)=0
set whether to generate weighted or unweighted neutrinos
double fBrFluxEnu
corresponding flux P4 (set to address of branch:"FluxP4")
Definition: GMCJDriver.h:134
double fBrFluxWeight
corresponding flux weight (set to address of branch: "FluxWeight")
Definition: GMCJDriver.h:135
string fFluxIntFileName
whether to save pre-generated flux tree for use in later jobs
Definition: GMCJDriver.h:137
virtual long int Index(void)=0
returns corresponding index for current flux neutrino (e.g. for a flux ntuple returns the current ent...
TFile * fFluxIntProbFile
[input] pre-generated flux interaction probability file
Definition: GMCJDriver.h:130
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
int fBrFluxPDG
corresponding flux pdg code (set to address of branch: "FluxPDG")
Definition: GMCJDriver.h:136
int fBrFluxIndex
corresponding entry in flux input tree (set to address of branch:"FluxEntry")
Definition: GMCJDriver.h:133
string fFluxIntTreeName
name for tree holding flux probabilities
Definition: GMCJDriver.h:138
static const double kASmallNum
Definition: Controls.h:40
def success(message)
Definition: log.py:5
double ComputeInteractionProbabilities(bool use_max_path_length)
double fGlobPmax
[computed at init] global interaction probability scale for given flux & geometry ...
Definition: GMCJDriver.h:120
#define pWARN
Definition: Messenger.h:61
virtual void Clear(Option_t *opt)=0
reset state variables based on opt
virtual bool GenerateNext(void)=0
generate the next flux neutrino (return false in err)
void PreSelectEvents(bool preselect=true)
Definition: GMCJDriver.cxx:231
map< int, double > fSumFluxIntProbs
map where the key is flux pdg code and the value is sum of fBrFluxWeight * fBrFluxIntProb for all the...
Definition: GMCJDriver.h:139
virtual int PdgCode(void)=0
returns the flux neutrino pdg code
exit(0)
virtual const TLorentzVector & Momentum(void)=0
returns the flux neutrino 4-momentum
assert(nhit_max >=nhit_nbins)
virtual bool End(void)=0
true if no more flux nu&#39;s can be thrown (eg reaching end of beam sim ntuples)
virtual double Weight(void)=0
returns the flux neutrino weight (if any)
#define pNOTICE
Definition: Messenger.h:62
double GMCJDriver::PreGenFluxInteractionProbability ( void  )
private

Definition at line 1318 of file GMCJDriver.cxx.

References ana::assert(), exit(), genie::controls::kASmallNum, LOG, pERROR, pFATAL, and log::success().

Referenced by GeomAnalyzerPtr().

1319 {
1320 // Return the pre-computed interaction probability for the current flux
1321 // neutrino index (entry number in flux file). Exit if not possible as
1322 // using meaningless interaction probability leads to incorrect physics
1323 //
1324  if(!fFluxIntTree){
1325  LOG("GMCJDriver", pERROR) <<
1326  "Cannot get pre-computed flux interaction probability as no tree!";
1327  exit(1);
1328  }
1329 
1330  assert(fFluxDriver->Index() >= 0); // Check trying to find meaningfull index
1331 
1332  // Check if can find relevant entry and no mismatch in energies -->
1333  // using correct pre-gen interaction prob file
1334  bool found_entry = fFluxIntTree->GetEntryWithIndex(fFluxDriver->Index()) > 0;
1335  bool enu_match = false;
1336  if(found_entry){
1337  double rel_err = fBrFluxEnu-fFluxDriver->Momentum().E();
1338  if(fBrFluxEnu > controls::kASmallNum) rel_err /= fBrFluxEnu;
1339  enu_match = TMath::Abs(rel_err)<controls::kASmallNum;
1340  if(enu_match == false){
1341  LOG("GMCJDriver", pERROR) <<
1342  "Mismatch between: Enu_curr = "<< fFluxDriver->Momentum().E() <<
1343  ", Enu_pre_gen = "<< fBrFluxEnu;
1344  }
1345  }
1346  else {
1347  LOG("GMCJDriver", pERROR) << "Cannot find flux entry in interaction prob tree!";
1348  }
1349 
1350  // Exit if not successful
1351  bool success = found_entry && enu_match;
1352  if(!success){
1353  LOG("GMCJDriver", pFATAL) <<
1354  "Cannot find pre-generated interaction probability! Check you "<<
1355  "are using the correct pre-generated interaction prob file " <<
1356  "generated using current flux input file with same input " <<
1357  "config (same geom TopVol, neutrino species list)";
1358  exit(1);
1359  }
1361  return fBrFluxIntProb/fGlobPmax;
1362 }
TTree * fFluxIntTree
[computed-or-loaded] pre-computed flux interaction probabilities (expected tree name is "gFlxIntProbs...
Definition: GMCJDriver.h:131
#define pERROR
Definition: Messenger.h:60
double fBrFluxIntProb
flux interaction probability (set to branch:"FluxIntProb")
Definition: GMCJDriver.h:132
#define pFATAL
Definition: Messenger.h:57
GFluxI * fFluxDriver
[input] neutrino flux driver
Definition: GMCJDriver.h:107
double fBrFluxEnu
corresponding flux P4 (set to address of branch:"FluxP4")
Definition: GMCJDriver.h:134
virtual long int Index(void)=0
returns corresponding index for current flux neutrino (e.g. for a flux ntuple returns the current ent...
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
static const double kASmallNum
Definition: Controls.h:40
def success(message)
Definition: log.py:5
double fGlobPmax
[computed at init] global interaction probability scale for given flux & geometry ...
Definition: GMCJDriver.h:120
exit(0)
virtual const TLorentzVector & Momentum(void)=0
returns the flux neutrino 4-momentum
assert(nhit_max >=nhit_nbins)
void GMCJDriver::PreSelectEvents ( bool  preselect = true)

Definition at line 231 of file GMCJDriver.cxx.

232 {
233 // Set whether to pre-select events based on a max-path lengths file. This
234 // should be turned off if using pre-generated interaction probabilities
235 // calculated from a given flux file.
236  fPreSelect = preselect;
237 }
bool fPreSelect
[config] set whether to pre-select events using max interaction paths
Definition: GMCJDriver.h:129
void GMCJDriver::SaveFluxProbabilities ( string  outfilename)

Definition at line 437 of file GMCJDriver.cxx.

References make_template_knob_config::outfilename.

Referenced by main().

438 {
439 // Configue the flux driver to save the calculated flux interaction
440 // probabilities to the specified output file name for use in later jobs. See
441 // the LoadFluxProbTree method for how they are fed into a later job.
442 //
444 }
string fFluxIntFileName
whether to save pre-generated flux tree for use in later jobs
Definition: GMCJDriver.h:137
string outfilename
knobs that need extra care
int GMCJDriver::SelectTargetMaterial ( double  R)
private

Definition at line 1190 of file GMCJDriver.cxx.

References LOG, pERROR, and pNOTICE.

Referenced by GeomAnalyzerPtr().

1191 {
1192 // Pick a target material using the pre-computed interaction probabilities
1193 // for a flux neutrino that has already been determined that interacts
1194 
1195  LOG("GMCJDriver", pNOTICE) << "Selecting target material";
1196  int tgtpdg = 0;
1197  map<int,double>::const_iterator probiter = fCurCumulProbMap.begin();
1198  for( ; probiter != fCurCumulProbMap.end(); ++probiter) {
1199  double prob = probiter->second;
1200  if(R<prob) {
1201  tgtpdg = probiter->first;
1202  LOG("GMCJDriver", pNOTICE)
1203  << "Selected target material = " << tgtpdg;
1204  return tgtpdg;
1205  }
1206  }
1207  LOG("GMCJDriver", pERROR)
1208  << "Could not select target material for an interacting neutrino";
1209  return 0;
1210 }
#define pERROR
Definition: Messenger.h:60
map< int, double > fCurCumulProbMap
[current] cummulative interaction probabilities
Definition: GMCJDriver.h:117
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
#define R(x)
#define pNOTICE
Definition: Messenger.h:62
void GMCJDriver::SetEventGeneratorList ( string  listname)

Definition at line 157 of file GMCJDriver.cxx.

References LOG, and pNOTICE.

Referenced by GenerateEventsAtFixedInitState(), main(), and supernova::SnovaGen::SnovaGen().

158 {
159  LOG("GMCJDriver", pNOTICE)
160  << "Setting event generator list: " << listname;
161 
162  fEventGenList = listname;
163 }
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
string fEventGenList
[config] list of event generators loaded by this driver (what used to be the $GEVGL setting) ...
Definition: GMCJDriver.h:121
#define pNOTICE
Definition: Messenger.h:62
void GMCJDriver::SetUnphysEventMask ( const TBits &  mask)

Definition at line 165 of file GMCJDriver.cxx.

References LOG, lem_server::mask, genie::GHepFlags::NFlags(), and pNOTICE.

Referenced by GenerateEventsAtFixedInitState().

166 {
168 
169  LOG("GMCJDriver", pNOTICE)
170  << "Setting unphysical event mask (bits: " << GHepFlags::NFlags() - 1
171  << " -> 0) : " << *fUnphysEventMask;
172 }
TBits * fUnphysEventMask
[config] controls whether unphysical events are returned (what used to be the $GUNPHYSMASK setting) ...
Definition: GMCJDriver.h:122
static unsigned int NFlags(void)
Definition: GHepFlags.h:77
#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
map<int, double> genie::GMCJDriver::SumFluxIntProbs ( void  ) const
inline

Definition at line 74 of file GMCJDriver.h.

References fSumFluxIntProbs.

Referenced by main().

74 { return fSumFluxIntProbs; }
map< int, double > fSumFluxIntProbs
map where the key is flux pdg code and the value is sum of fBrFluxWeight * fBrFluxIntProb for all the...
Definition: GMCJDriver.h:139
void GMCJDriver::UseFluxDriver ( GFluxI flux)

Definition at line 174 of file GMCJDriver.cxx.

Referenced by GenerateEventsAtFixedInitState(), main(), and supernova::SnovaGen::SnovaGen().

175 {
176  fFluxDriver = flux_driver;
177 }
GFluxI * fFluxDriver
[input] neutrino flux driver
Definition: GMCJDriver.h:107
void GMCJDriver::UseGeomAnalyzer ( GeomAnalyzerI geom)

Definition at line 179 of file GMCJDriver.cxx.

Referenced by GenerateEventsAtFixedInitState(), main(), and supernova::SnovaGen::SnovaGen().

180 {
181  fGeomAnalyzer = geom_analyzer;
182 }
GeomAnalyzerI * fGeomAnalyzer
[input] detector geometry analyzer
Definition: GMCJDriver.h:108
bool GMCJDriver::UseMaxPathLengths ( string  xml_filename)

Definition at line 190 of file GMCJDriver.cxx.

References LOG, and pWARN.

Referenced by main(), and supernova::SnovaGen::SnovaGen().

191 {
192 // If you supply the maximum path lengths for all materials in your geometry
193 // (eg for ROOT/GEANT geometries they can be computed running GENIE's gmxpl
194 // application, see $GENIE/src/stdapp/gMaxPathLengths.cxx ) you can speed up
195 // the driver init phase by quite a bit (especially for complex geometries).
196 
197  fMaxPlXmlFilename = xml_filename;
198 
199  bool is_accessible = !(gSystem->AccessPathName(fMaxPlXmlFilename.c_str()));
200 
201  if ( is_accessible ) fUseExtMaxPl = true;
202  else {
203  fUseExtMaxPl = false;
204  LOG("GMCJDriver", pWARN)
205  << "UseMaxPathLengths could not find file: \"" << xml_filename << "\"";
206  }
207  return fUseExtMaxPl;
208 
209 }
bool fUseExtMaxPl
[config] using external max path length estimate?
Definition: GMCJDriver.h:124
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
#define pWARN
Definition: Messenger.h:61
string fMaxPlXmlFilename
[config] input file with max density-weighted path lengths for all materials
Definition: GMCJDriver.h:123
void GMCJDriver::UseSplines ( bool  useLogE = true)

Definition at line 184 of file GMCJDriver.cxx.

Referenced by GenerateEventsAtFixedInitState(), and main().

185 {
186  fUseSplines = true;
187  fUseLogE = useLogE;
188 }
bool fUseLogE
[config] build splines = f(logE) (rather than f(E)) ?
Definition: GMCJDriver.h:126
bool fUseSplines
[config] compute all needed & not-loaded splines at init
Definition: GMCJDriver.h:125

Member Data Documentation

double genie::GMCJDriver::fBrFluxEnu
private

corresponding flux P4 (set to address of branch:"FluxP4")

Definition at line 134 of file GMCJDriver.h.

int genie::GMCJDriver::fBrFluxIndex
private

corresponding entry in flux input tree (set to address of branch:"FluxEntry")

Definition at line 133 of file GMCJDriver.h.

double genie::GMCJDriver::fBrFluxIntProb
private

flux interaction probability (set to branch:"FluxIntProb")

Definition at line 132 of file GMCJDriver.h.

int genie::GMCJDriver::fBrFluxPDG
private

corresponding flux pdg code (set to address of branch: "FluxPDG")

Definition at line 136 of file GMCJDriver.h.

double genie::GMCJDriver::fBrFluxWeight
private

corresponding flux weight (set to address of branch: "FluxWeight")

Definition at line 135 of file GMCJDriver.h.

map<int,double> genie::GMCJDriver::fCurCumulProbMap
private

[current] cummulative interaction probabilities

Definition at line 117 of file GMCJDriver.h.

EventRecord* genie::GMCJDriver::fCurEvt
private

[current] generated event

Definition at line 115 of file GMCJDriver.h.

PathLengthList genie::GMCJDriver::fCurPathLengths
private

[current] path length list for current flux neutrino

Definition at line 113 of file GMCJDriver.h.

TLorentzVector genie::GMCJDriver::fCurVtx
private

[current] interaction vertex

Definition at line 114 of file GMCJDriver.h.

double genie::GMCJDriver::fEmax
private

[declared by the flux driver] maximum neutrino energy

Definition at line 109 of file GMCJDriver.h.

string genie::GMCJDriver::fEventGenList
private

[config] list of event generators loaded by this driver (what used to be the $GEVGL setting)

Definition at line 121 of file GMCJDriver.h.

GFluxI* genie::GMCJDriver::fFluxDriver
private

[input] neutrino flux driver

Definition at line 107 of file GMCJDriver.h.

Referenced by FluxDriver(), and FluxDriverPtr().

string genie::GMCJDriver::fFluxIntFileName
private

whether to save pre-generated flux tree for use in later jobs

Definition at line 137 of file GMCJDriver.h.

TFile* genie::GMCJDriver::fFluxIntProbFile
private

[input] pre-generated flux interaction probability file

Definition at line 130 of file GMCJDriver.h.

TTree* genie::GMCJDriver::fFluxIntTree
private

[computed-or-loaded] pre-computed flux interaction probabilities (expected tree name is "gFlxIntProbs")

Definition at line 131 of file GMCJDriver.h.

string genie::GMCJDriver::fFluxIntTreeName
private

name for tree holding flux probabilities

Definition at line 138 of file GMCJDriver.h.

bool genie::GMCJDriver::fGenerateUnweighted
private

[config] force single probability scale?

Definition at line 128 of file GMCJDriver.h.

GeomAnalyzerI* genie::GMCJDriver::fGeomAnalyzer
private

[input] detector geometry analyzer

Definition at line 108 of file GMCJDriver.h.

Referenced by GeomAnalyzer(), and GeomAnalyzerPtr().

double genie::GMCJDriver::fGlobPmax
private

[computed at init] global interaction probability scale for given flux & geometry

Definition at line 120 of file GMCJDriver.h.

Referenced by GlobProbScale().

GEVGPool* genie::GMCJDriver::fGPool
private

A pool of GEVGDrivers properly configured event generation drivers / one per init state.

Definition at line 106 of file GMCJDriver.h.

bool genie::GMCJDriver::fKeepThrowingFluxNu
private

[config] keep firing flux neutrinos till one of them interacts

Definition at line 127 of file GMCJDriver.h.

PathLengthList genie::GMCJDriver::fMaxPathLengths
private

[declared by the geom driver] maximum path length list

Definition at line 112 of file GMCJDriver.h.

string genie::GMCJDriver::fMaxPlXmlFilename
private

[config] input file with max density-weighted path lengths for all materials

Definition at line 123 of file GMCJDriver.h.

double genie::GMCJDriver::fNFluxNeutrinos
private

[current] number of flux nuetrinos fired by the flux driver so far

Definition at line 118 of file GMCJDriver.h.

Referenced by NFluxNeutrinos().

PDGCodeList genie::GMCJDriver::fNuList
private

[declared by the flux driver] list of neutrino codes

Definition at line 110 of file GMCJDriver.h.

map<int,TH1D*> genie::GMCJDriver::fPmax
private

[computed at init] interaction probability scale /neutrino /energy for given geometry

Definition at line 119 of file GMCJDriver.h.

bool genie::GMCJDriver::fPreSelect
private

[config] set whether to pre-select events using max interaction paths

Definition at line 129 of file GMCJDriver.h.

int genie::GMCJDriver::fSelTgtPdg
private

[current] selected target material PDG code

Definition at line 116 of file GMCJDriver.h.

map<int, double> genie::GMCJDriver::fSumFluxIntProbs
private

map where the key is flux pdg code and the value is sum of fBrFluxWeight * fBrFluxIntProb for all these flux neutrinos

Definition at line 139 of file GMCJDriver.h.

Referenced by SumFluxIntProbs().

PDGCodeList genie::GMCJDriver::fTgtList
private

[declared by the geom driver] list of target codes

Definition at line 111 of file GMCJDriver.h.

TBits* genie::GMCJDriver::fUnphysEventMask
private

[config] controls whether unphysical events are returned (what used to be the $GUNPHYSMASK setting)

Definition at line 122 of file GMCJDriver.h.

bool genie::GMCJDriver::fUseExtMaxPl
private

[config] using external max path length estimate?

Definition at line 124 of file GMCJDriver.h.

bool genie::GMCJDriver::fUseLogE
private

[config] build splines = f(logE) (rather than f(E)) ?

Definition at line 126 of file GMCJDriver.h.

bool genie::GMCJDriver::fUseSplines
private

[config] compute all needed & not-loaded splines at init

Definition at line 125 of file GMCJDriver.h.


The documentation for this class was generated from the following files: