GMCJDriver.h
Go to the documentation of this file.
1 //____________________________________________________________________________
2 /*!
3 
4 \class genie::GMCJDriver
5 
6 \brief A GENIE `MC Job Driver'. Can be used for setting up complicated event
7  generation cases involving detailed flux descriptions and detector
8  geometry descriptions.
9 
10 \author Costas Andreopoulos <costas.andreopoulos \at stfc.ac.uk>
11  University of Liverpool & STFC Rutherford Appleton Lab
12 
13 \created May 25, 2005
14 
15 \cpright Copyright (c) 2003-2019, The GENIE Collaboration
16  For the full text of the license visit http://copyright.genie-mc.org
17  or see $GENIE/LICENSE
18 */
19 //____________________________________________________________________________
20 
21 #ifndef _GENIE_MC_JOB_DRIVER_H_
22 #define _GENIE_MC_JOB_DRIVER_H_
23 
24 #include <string>
25 #include <map>
26 
27 #include <TH1D.h>
28 #include <TLorentzVector.h>
29 #include <TFile.h>
30 #include <TTree.h>
31 #include <TBits.h>
32 
35 
36 using std::string;
37 using std::map;
38 
39 namespace genie {
40 
41 class EventRecord;
42 class GFluxI;
43 class GeomAnalyzerI;
44 class GENIE;
45 class GEVGPool;
46 
47 class GMCJDriver {
48 
49 public :
50  GMCJDriver();
51  ~GMCJDriver();
52 
53  // configure MC job
54  void SetEventGeneratorList (string listname);
55  void SetUnphysEventMask (const TBits & mask);
56  void UseFluxDriver (GFluxI * flux);
58  void UseSplines (bool useLogE = true);
59  bool UseMaxPathLengths (string xml_filename);
60  void KeepOnThrowingFluxNeutrinos (bool keep_on);
61  void ForceSingleProbScale (void);
62  void PreSelectEvents (bool preselect = true);
63  bool PreCalcFluxProbabilities (void);
64  bool LoadFluxProbabilities (string filename);
65  void SaveFluxProbabilities (string outfilename);
66  void Configure (bool calc_prob_scales = true);
67 
68  // generate single neutrino event for input flux & geometry
69  EventRecord * GenerateEvent (void);
70 
71  // info needed for computing the generated sample normalization
72  double GlobProbScale (void) const { return fGlobPmax; }
73  long int NFluxNeutrinos (void) const { return (long int) fNFluxNeutrinos; }
74  map<int, double> SumFluxIntProbs(void) const { return fSumFluxIntProbs; }
75 
76  // input flux and geometry drivers
77  const GFluxI & FluxDriver (void) const { return *fFluxDriver; }
78  const GeomAnalyzerI & GeomAnalyzer (void) const { return *fGeomAnalyzer; }
79  GFluxI * FluxDriverPtr (void) const { return fFluxDriver; }
80  GeomAnalyzerI * GeomAnalyzerPtr (void) const { return fGeomAnalyzer; }
81 
82 private:
83 
84  // private methods:
85  void InitJob (void);
86  void InitEventGeneration (void);
87  void GetParticleLists (void);
88  void GetMaxPathLengthList (void);
89  void GetMaxFluxEnergy (void);
90  void PopulateEventGenDriverPool (void);
91  void BootstrapXSecSplines (void);
92  void BootstrapXSecSplineSummation (void);
93  void ComputeProbScales (void);
95  bool GenerateFluxNeutrino (void);
96  bool ComputePathLengths (void);
97  double ComputeInteractionProbabilities (bool use_max_path_length);
98  int SelectTargetMaterial (double R);
99  void GenerateEventKinematics (void);
100  void GenerateVertexPosition (void);
101  void ComputeEventProbability (void);
102  double InteractionProbability (double xsec, double pl, int A);
104 
105  // private data members:
106  GEVGPool * fGPool; ///< A pool of GEVGDrivers properly configured event generation drivers / one per init state
107  GFluxI * fFluxDriver; ///< [input] neutrino flux driver
108  GeomAnalyzerI * fGeomAnalyzer; ///< [input] detector geometry analyzer
109  double fEmax; ///< [declared by the flux driver] maximum neutrino energy
110  PDGCodeList fNuList; ///< [declared by the flux driver] list of neutrino codes
111  PDGCodeList fTgtList; ///< [declared by the geom driver] list of target codes
112  PathLengthList fMaxPathLengths; ///< [declared by the geom driver] maximum path length list
113  PathLengthList fCurPathLengths; ///< [current] path length list for current flux neutrino
114  TLorentzVector fCurVtx; ///< [current] interaction vertex
115  EventRecord * fCurEvt; ///< [current] generated event
116  int fSelTgtPdg; ///< [current] selected target material PDG code
117  map<int,double> fCurCumulProbMap; ///< [current] cummulative interaction probabilities
118  double fNFluxNeutrinos; ///< [current] number of flux nuetrinos fired by the flux driver so far
119  map<int,TH1D*> fPmax; ///< [computed at init] interaction probability scale /neutrino /energy for given geometry
120  double fGlobPmax; ///< [computed at init] global interaction probability scale for given flux & geometry
121  string fEventGenList; ///< [config] list of event generators loaded by this driver (what used to be the $GEVGL setting)
122  TBits * fUnphysEventMask; ///< [config] controls whether unphysical events are returned (what used to be the $GUNPHYSMASK setting)
123  string fMaxPlXmlFilename; ///< [config] input file with max density-weighted path lengths for all materials
124  bool fUseExtMaxPl; ///< [config] using external max path length estimate?
125  bool fUseSplines; ///< [config] compute all needed & not-loaded splines at init
126  bool fUseLogE; ///< [config] build splines = f(logE) (rather than f(E)) ?
127  bool fKeepThrowingFluxNu; ///< [config] keep firing flux neutrinos till one of them interacts
128  bool fGenerateUnweighted; ///< [config] force single probability scale?
129  bool fPreSelect; ///< [config] set whether to pre-select events using max interaction paths
130  TFile* fFluxIntProbFile; ///< [input] pre-generated flux interaction probability file
131  TTree* fFluxIntTree; ///< [computed-or-loaded] pre-computed flux interaction probabilities (expected tree name is "gFlxIntProbs")
132  double fBrFluxIntProb; ///< flux interaction probability (set to branch:"FluxIntProb")
133  int fBrFluxIndex; ///< corresponding entry in flux input tree (set to address of branch:"FluxEntry")
134  double fBrFluxEnu; ///< corresponding flux P4 (set to address of branch:"FluxP4")
135  double fBrFluxWeight; ///< corresponding flux weight (set to address of branch: "FluxWeight")
136  int fBrFluxPDG; ///< corresponding flux pdg code (set to address of branch: "FluxPDG")
137  string fFluxIntFileName; ///< whether to save pre-generated flux tree for use in later jobs
138  string fFluxIntTreeName; ///< name for tree holding flux probabilities
139  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
140 };
141 
142 } // genie namespace
143 #endif // _GENIE_MC_JOB_DRIVER_H_
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
void SetUnphysEventMask(const TBits &mask)
Definition: GMCJDriver.cxx:165
GFluxI * FluxDriverPtr(void) const
Definition: GMCJDriver.h:79
void PopulateEventGenDriverPool(void)
Definition: GMCJDriver.cxx:605
bool PreCalcFluxProbabilities(void)
Definition: GMCJDriver.cxx:239
void InitJob(void)
Definition: GMCJDriver.cxx:494
map< int, double > SumFluxIntProbs(void) const
Definition: GMCJDriver.h:74
GEVGPool * fGPool
A pool of GEVGDrivers properly configured event generation drivers / one per init state...
Definition: GMCJDriver.h:106
THE MAIN GENIE PROJECT NAMESPACE
Definition: GeneratorBase.h:8
GeomAnalyzerI * GeomAnalyzerPtr(void) const
Definition: GMCJDriver.h:80
void GetMaxPathLengthList(void)
Definition: GMCJDriver.cxx:577
void SetEventGeneratorList(string listname)
Definition: GMCJDriver.cxx:157
map< int, double > fCurCumulProbMap
[current] cummulative interaction probabilities
Definition: GMCJDriver.h:117
void InitEventGeneration(void)
Definition: GMCJDriver.cxx:830
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
bool GenerateFluxNeutrino(void)
bool ComputePathLengths(void)
GFluxI * fFluxDriver
[input] neutrino flux driver
Definition: GMCJDriver.h:107
string filename
Definition: shutoffs.py:106
const GFluxI & FluxDriver(void) const
Definition: GMCJDriver.h:77
void ComputeEventProbability(void)
double fBrFluxEnu
corresponding flux P4 (set to address of branch:"FluxP4")
Definition: GMCJDriver.h:134
Loaders::FluxType flux
bool fUseExtMaxPl
[config] using external max path length estimate?
Definition: GMCJDriver.h:124
A list of PDG codes.
Definition: PDGCodeList.h:33
double fBrFluxWeight
corresponding flux weight (set to address of branch: "FluxWeight")
Definition: GMCJDriver.h:135
void UseFluxDriver(GFluxI *flux)
Definition: GMCJDriver.cxx:174
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
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.
void KeepOnThrowingFluxNeutrinos(bool keep_on)
Definition: GMCJDriver.cxx:211
void GetMaxFluxEnergy(void)
Definition: GMCJDriver.cxx:595
string outfilename
knobs that need extra care
bool fGenerateUnweighted
[config] force single probability scale?
Definition: GMCJDriver.h:128
A GENIE `MC Job Driver&#39;. Can be used for setting up complicated event generation cases involving deta...
Definition: GMCJDriver.h:47
double GlobProbScale(void) const
Definition: GMCJDriver.h:72
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
TFile * fFluxIntProbFile
[input] pre-generated flux interaction probability file
Definition: GMCJDriver.h:130
bool fKeepThrowingFluxNu
[config] keep firing flux neutrinos till one of them interacts
Definition: GMCJDriver.h:127
int fSelTgtPdg
[current] selected target material PDG code
Definition: GMCJDriver.h:116
void ForceSingleProbScale(void)
Definition: GMCJDriver.cxx:219
double fNFluxNeutrinos
[current] number of flux nuetrinos fired by the flux driver so far
Definition: GMCJDriver.h:118
double PreGenFluxInteractionProbability(void)
int fBrFluxPDG
corresponding flux pdg code (set to address of branch: "FluxPDG")
Definition: GMCJDriver.h:136
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
bool UseMaxPathLengths(string xml_filename)
Definition: GMCJDriver.cxx:190
void Configure(bool calc_prob_scales=true)
Definition: GMCJDriver.cxx:446
int fBrFluxIndex
corresponding entry in flux input tree (set to address of branch:"FluxEntry")
Definition: GMCJDriver.h:133
#define R(x)
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
const GeomAnalyzerI & GeomAnalyzer(void) const
Definition: GMCJDriver.h:78
void BootstrapXSecSplines(void)
Definition: GMCJDriver.cxx:642
double ComputeInteractionProbabilities(bool use_max_path_length)
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
TLorentzVector fCurVtx
[current] interaction vertex
Definition: GMCJDriver.h:114
EventRecord * GenerateEvent(void)
Definition: GMCJDriver.cxx:838
long int NFluxNeutrinos(void) const
Definition: GMCJDriver.h:73
string fMaxPlXmlFilename
[config] input file with max density-weighted path lengths for all materials
Definition: GMCJDriver.h:123
Generated Event Record. It is a GHepRecord object that can accept / be visited by EventRecordVisitorI...
Definition: EventRecord.h:38
double InteractionProbability(double xsec, double pl, int A)
::xsd::cxx::tree::string< char, simple_type > string
Definition: Database.h:154
static const double A
Definition: Units.h:82
void PreSelectEvents(bool preselect=true)
Definition: GMCJDriver.cxx:231
void GenerateEventKinematics(void)
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
int SelectTargetMaterial(double R)
void GenerateVertexPosition(void)
GeomAnalyzerI * fGeomAnalyzer
[input] detector geometry analyzer
Definition: GMCJDriver.h:108
void geom(int which=0)
Definition: geom.C:163
void ComputeProbScales(void)
Definition: GMCJDriver.cxx:705
bool LoadFluxProbabilities(string filename)
Definition: GMCJDriver.cxx:390
EventRecord * GenerateEvent1Try(void)
Definition: GMCJDriver.cxx:867
void BootstrapXSecSplineSummation(void)
Definition: GMCJDriver.cxx:669
Defines the GENIE Geometry Analyzer Interface.
Definition: GeomAnalyzerI.h:30
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
void UseGeomAnalyzer(GeomAnalyzerI *geom)
Definition: GMCJDriver.cxx:179
void UseSplines(bool useLogE=true)
Definition: GMCJDriver.cxx:184
void SaveFluxProbabilities(string outfilename)
Definition: GMCJDriver.cxx:437
void GetParticleLists(void)
Definition: GMCJDriver.cxx:560
GENIE Interface for user-defined flux classes.
Definition: GFluxI.h:37
A pool of GEVGDriver objects with an initial state key.
Definition: GEVGPool.h:38
double fEmax
[declared by the flux driver] maximum neutrino energy
Definition: GMCJDriver.h:109