7 #include "TStopwatch.h" 8 #include "TLorentzVector.h" 18 using namespace genie;
23 string cfg=
"MINOS-Near",
27 double enumin=0.0,
bool doaux=
false);
33 double enumin=0.0,
bool doaux=
false)
36 cout <<
"output file: " << fname <<
endl;
39 string cfg =
"MINOS-NearDet";
40 if ( det ==
"minos" ) cfg =
"MINOS-NearDet";
41 if ( det ==
"novaipnd" ) cfg =
"NOvA-IPNDShed";
42 if ( det ==
"uboone" ) cfg =
"microboone-numi";
44 cout <<
"det \"" << det <<
"\" ==> cfg " << cfg <<
endl;
47 if ( flux ==
"old" ) flxset = 0;
48 if ( flux ==
"old1" ) flxset = 10;
49 if ( flux ==
"old2" ) flxset = 20;
50 if ( flux ==
"lowcut" ) flxset = 1;
51 if ( flux ==
"lowcut1" ) flxset = 11;
52 if ( flux ==
"lowcut2" ) flxset = 21;
53 if ( flux ==
"flugg" ) flxset = 2;
54 if ( flux ==
"flugg1" ) flxset = 12;
55 if ( flux ==
"flugg2" ) flxset = 22;
56 if ( flux ==
"flugglowth" ) flxset = 3;
57 if ( flux ==
"flugglowth1" ) flxset = 13;
58 if ( flux ==
"flugglowth2" ) flxset = 23;
60 cout <<
"flux \"" << flux <<
"\" ==> flxset " << flxset <<
endl;
68 string fluxpatt =
"$FLUXPATH";
70 int flxver = flxset % 10;
73 fluxpatt +=
"/v19/fluka05_le010z185i/root/fluka05_le010z185i_";
76 fluxpatt +=
"/v19/fluka05_le010z185i_lowcut/root/fluka05_le010z185i_";
79 fluxpatt +=
"/flugg/flugg_le010z185i_run1/root/flugg_le010z185i_run1_";
82 fluxpatt +=
"/flugg_lowth/flugg_le010z185i_run1_lowth/root/flugg_le010z185i_run1_lowth_";
85 cout <<
" not a legal flux set " <<
endl;
89 if ( flxset/10 == 1 ) {
90 if ( flxver == 2 || flxver == 3 ) fluxpatt +=
"001";
92 }
else if ( flxset/10 == 2 ) {
93 if ( flxver == 2 || flxver == 3 ) fluxpatt +=
"00[12]";
94 else fluxpatt +=
"[12]";
95 }
else fluxpatt +=
"*";
98 string fluxfname(gSystem->ExpandPathName(fluxpatt.c_str()));
120 if (nentries == 0) nentries = 100000*20;
125 TFile*
file = TFile::Open(fname.c_str(),
"RECREATE");
126 TTree* fluxntp =
new TTree(
"flux",
"a simple flux n-tuple");
127 TTree* metantp =
new TTree(
"meta",
"metadata for flux n-tuple");
132 fluxntp->Branch(
"entry",&fentry);
133 fluxntp->Branch(
"numi",&fnumi);
134 if (doaux) fluxntp->Branch(
"aux",&faux);
135 metantp->Branch(
"meta",&fmeta);
137 TLorentzVector p4u, x4u;
143 double minwgt = +1.0e10;
144 double maxwgt = -1.0e10;
146 UInt_t metakey =
TString::Hash(fname.c_str(),strlen(fname.c_str()));
147 cout <<
"metakey " << metakey <<
endl;
149 metakey &= 0x7FFFFFFF;
150 cout <<
"metakey " << metakey <<
" after 0x7FFFFFFF" <<
endl;
151 cout <<
"=========================== Start " <<
endl;
155 int nok = 0, nwrite = 0;
156 while ( nok < nentries && gnumi->UsedPOTs() < pots ) {
168 fentry->
vtxx = x4u.X();
169 fentry->
vtxy = x4u.Y();
170 fentry->
vtxz = x4u.Z();
173 fentry->
px = p4u.Px();
174 fentry->
py = p4u.Py();
175 fentry->
pz = p4u.Pz();
201 if ( fentry->
E > enumin ) ++nok;
206 pdglist.insert(fentry->
pdg);
207 minwgt = TMath::Min(minwgt,fentry->
wgt);
208 maxwgt = TMath::Max(maxwgt,fentry->
wgt);
209 maxe = TMath::Max(maxe,fentry->
E);
212 cout <<
"=========================== Complete " <<
endl;
222 set<int>::const_iterator setitr = pdglist.begin();
223 for ( ; setitr != pdglist.end(); ++setitr) fmeta->
pdglist.push_back(*setitr);
231 TVector3 d1 = p1 -
p0;
232 TVector3 d2 = p2 -
p0;
248 string fname =
"flux_pattern=" + fluxfname;
249 fmeta->
infiles.push_back(fluxfname);
251 for (
size_t i = 0;
i < flist.size(); ++
i) {
252 fname =
"flux_infile=" + flist[
i];
253 fmeta->
infiles.push_back(fname);
256 fmeta->
seed = RandomGen::Instance()->GetSeed();
262 cout <<
"Generated " << nwrite <<
" (" << nok <<
" w/ E >" << enumin <<
" ) " 263 <<
" ( request " << nentries <<
" ) " 265 << gnumi->
UsedPOTs() <<
" POTs " <<
" ( request " << pots <<
" )" 269 <<
"Time to generate: " <<
endl;
272 cout <<
"===================================================" <<
endl;
274 cout <<
"Last GSimpleNtpEntry: " << endl
277 <<
"Last GSimpleNtpMeta: " << endl
298 double enumin=0.0,
bool doaux=
false)
Double_t E
energy in lab frame
void GetFluxWindow(TVector3 &p1, TVector3 &p2, TVector3 &p3) const
3 points define a plane in beam coordinate
void SetEntryReuse(long int nuse=1)
of times to use entry before moving to next
void fill_simple(string fname="fluxntgenie.root", string cfg="MINOS-Near", int flxset=0, long int nentries=0, double pots=1.0e30, double enumin=0.0, bool doaux=false)
THE MAIN GENIE PROJECT NAMESPACE
Double_t px
x momentum in lab frame
Long64_t GetEntryNumber()
index in chain
Double_t vtxy
y position in lab frame
Double_t vx
vertex position of hadron/muon decay
virtual const TLorentzVector & Position(void)=0
returns the flux neutrino 4-position (note: expect SI rather than physical units) ...
Int_t tptype
parent particle type at target exit
void gnumi2simple(string fname, string det, string flux, long int nentries=5000, double pots=1.0e30, double enumin=0.0, bool doaux=false)
std::vector< Double_t > auxdbl
additional doubles associated w/ entry
std::vector< Int_t > auxint
additional ints associated w/ entry
virtual void SetUpstreamZ(double z0)
A GENIE flux driver encapsulating the NuMI neutrino flux. It reads-in the official GNUMI neutrino flu...
double Weight(void)
returns the flux neutrino weight (if any)
virtual void LoadBeamSimData(const std::vector< std::string > &filenames, const std::string &det_loc)
std::vector< std::string > GetFileList()
list of files currently part of chain
const GNuMIFluxPassThroughInfo & PassThroughInfo(void)
GNuMIFluxPassThroughInfo.
Double_t vtxz
z position in lab frame
Int_t ppmedium
tracking medium where parent was produced
Double_t tpx
parent particle px at target exit
virtual bool GenerateNext(void)=0
generate the next flux neutrino (return false in err)
double GetDecayDist() const
dist (user units) from dk to current pos
std::size_t Hash(const std::vector< std::size_t > &vals)
Generate a unique hash from a sequence of integers.
Double_t vtxx
x position in lab frame
Double_t pz
z momentum in lab frame
virtual int PdgCode(void)=0
returns the flux neutrino pdg code
virtual const TLorentzVector & Momentum(void)=0
returns the flux neutrino 4-momentum
void make_simple_grid(string fname, string det, string flux, long int nentries=5000, double pots=1.0e30, double enumin=0.0, bool doaux=false)
void PrintConfig()
print the current configuration
Double_t dist
distance from hadron decay
virtual long int NFluxNeutrinos() const
of rays generated
UInt_t metakey
key to meta data
Double_t py
y momentum in lab frame
GENIE Interface for user-defined flux classes.
double UsedPOTs(void) const
of protons-on-target used