convert_g4lbne.C
Go to the documentation of this file.
1 //
2 // Convert a g4lbne ntuple to the new format
3 //
4 // this script can be run using:
5 // $ root $(DK2NU)/snippets/load_dk2nu.C convert_g4lbne.C("file.root")
6 //
7 // rhatcher@fnal.gov 2012-11-06
8 //====================================================================
9 
10 #include <iostream>
11 #include <iomanip>
12 #include <string>
13 using namespace std;
14 
15 #include "TMath.h"
16 
17 #include "convert/common_convert.C"
18 #include "convert/g4lbne/g4lbne.C"
19 
20 void copy_g4lbne_to_dk2nu(const g4lbne& g4lbneObj);
21 void g4lbneCrossChecks(const g4lbne& g4lbneObj);
22 
23 void convert_g4lbne(string ifname="../fluxfiles/generic_g4lbne.root",
24  int jobnum=42,
25  string locfile="${DK2NU}/etc/LBNElocations.txt",
26  Long64_t maxentries=-1,
27  Long64_t moddump=-1) // modulo for dump
28 {
29  // set globals
30  myjob = jobnum; // allow override because perhaps g4lbne files forgot to set this
31  pots = 0;
32  // allowance in location energy/weight cross-check
33  frac_diff_tolerance = 2.5e-4;
34 
35  cout << "locations from " << locfile << endl;
36 
37 
38  int highest_potnum = 0;
39 
40  // open input ntuple
41  TFile* fin = TFile::Open(ifname.c_str());
42  if ( ! fin ) {
43  cout << "couldn't open input file: " << ifname << endl;
44  return;
45  }
46  const char* objName = "nudata";
47  TTree* tin = 0;
48  fin->GetObject(objName,tin);
49 
50  if ( ! tin ) {
51  cout << "couldn't find " << objName << endl;
52  return;
53  }
54  cout << endl << "Input file: " << ifname << endl;
55 
56  // use makeclass interface for g4lbne
57  g4lbne g4lbneObj(tin);
58 
59  // construct output filename
60  string ofname = construct_outfilename(ifname);
61  cout << "Output file: " << ofname << endl << endl;
62 
63  ConvertInit(locfile);
64  ConvertBookNtuple(ofname);
65 
66  ///
67  /// loop over entries
68  ///
69  Long64_t nentries = g4lbneObj.fChain->GetEntriesFast();
70  if ( maxentries > 0 ) {
71  cout << "limit # of entries to " << nentries
72  << " / " << maxentries << endl;
73  nentries = TMath::Min(nentries,maxentries);
74  }
75  cout << endl;
76 
77  Long64_t nbytes = 0, nb = 0;
78  for (Long64_t jentry = 0; jentry < nentries; ++jentry) {
79  Long64_t ientry = g4lbneObj.LoadTree(jentry);
80  if (ientry < 0) break;
81  nb = g4lbneObj.fChain->GetEntry(jentry);
82  nbytes += nb;
83 
84  // always clear the dk2nu object
85  dk2nu->clear(); // !!! important !!! always do this
86 
87  // fill the dk2nu object from the g4lbne entry
88  copy_g4lbne_to_dk2nu(g4lbneObj);
89 
90  // fill location specific p3, energy and weights
91  // locations to fill are in the metadata
92  // assumes that prior copying filled the first entry w/ random decay
94 
95  // keep track of potnum
96  if ( dk2nu->potnum > highest_potnum )
97  highest_potnum = dk2nu->potnum;
98 
99  // push entry out to tree
100  dk2nuTree->Fill();
101 
102  // just for fun print every n entries
103  if ( moddump > 0 && jentry%moddump == 0 ) cout << endl << *dk2nu << endl;
104 
105  g4lbneCrossChecks(g4lbneObj);
106 
107  }
108  cout << endl;
109 
110  pots = estimate_pots(highest_potnum);
111  cout << "estimated pots " << pots << " (" << highest_potnum << ")" << endl;
112  cout << "maximum observed fractional difference in cross checks: "
113  << obs_frac_diff_max << endl
114  << "print ## message when exceeded " << frac_diff_tolerance << endl;
115 
116  dkmeta->job = myjob;
117  dkmeta->pots = pots;
118  dkmeta->beamsim = "convert_g4lbne.C";
119  dkmeta->physics = "bogus";
120 
121  dkmeta->tgtcfg = TString::Format("Z=%6.1fcm",g4lbneObj.nuTarZ);
122  dkmeta->horncfg = TString::Format("I=%6.1fkA",g4lbneObj.hornCurrent);
123 
124  dkmeta->beam0x = g4lbneObj.beamX;
125  dkmeta->beam0y = g4lbneObj.beamY;
126  dkmeta->beamhwidth = g4lbneObj.beamHWidth;
127  dkmeta->beamvwidth = g4lbneObj.beamVWidth;
128 
129  // print last entries
130  cout << endl << "========== last entries" << endl << endl;
131  cout << *dk2nu << endl << endl;
132  cout << *dkmeta << endl;
133 
134  // finish off ntuples (including meta data)
135  ConvertFinish();
136 
137  // close input file
138  fin->Close();
139 }
140 
141 void copy_g4lbne_to_dk2nu(const g4lbne& g4lbneObj)
142 {
143  // fill the global dk2nu object
144 
145  dk2nu->job = myjob;
146  dk2nu->potnum = g4lbneObj.evtno;
147 
148  // calcLocationWeights needs random decay entries filled first
149  // but don't copy any other (i.e. near/far) values as they'll
150  // be recalculated for the whole set of locations
151  double pzrndm = g4lbneObj.Npz;
152  double pxrndm = g4lbneObj.Ndxdz * pzrndm;
153  double pyrndm = g4lbneObj.Ndydz * pzrndm;
154  bsim::NuRay nuray(pxrndm,pyrndm,pzrndm,g4lbneObj.Nenergy,1.0);
155  dk2nu->nuray.push_back(nuray);
156 
157  // calcLocationWeights needs these filled if it isn't going assert()
158  // really need to fill the other bits at this point as well:
159  // ntype, ptype, vx, vy, vz, pdpx, pdpy, pdpz, necm,
160  // ppenergy, ppdxdz, ppdydz, pppz,
161  // muparpx, muparpy, muparpz, mupare
162 
163  dk2nu->decay.norig = g4lbneObj.Norig;
164  dk2nu->decay.ndecay = g4lbneObj.Ndecay;
165  dk2nu->decay.ntype = Convert5xToPdg(g4lbneObj.Ntype);
166  dk2nu->decay.vx = g4lbneObj.Vx;
167  dk2nu->decay.vy = g4lbneObj.Vy;
168  dk2nu->decay.vz = g4lbneObj.Vz;
169  dk2nu->decay.pdpx = g4lbneObj.pdPx;
170  dk2nu->decay.pdpy = g4lbneObj.pdPy;
171  dk2nu->decay.pdpz = g4lbneObj.pdPz;
172  dk2nu->decay.ppdxdz = g4lbneObj.ppdxdz;
173  dk2nu->decay.ppdydz = g4lbneObj.ppdydz;
174  dk2nu->decay.pppz = g4lbneObj.pppz;
175  dk2nu->decay.ppenergy = g4lbneObj.ppenergy;
176  dk2nu->decay.ppmedium = g4lbneObj.ppmedium;
177  dk2nu->decay.ptype = ConvertGeantToPdg(g4lbneObj.ptype,"ptype");
178  dk2nu->decay.muparpx = g4lbneObj.muparpx;
179  dk2nu->decay.muparpy = g4lbneObj.muparpy;
180  dk2nu->decay.muparpz = g4lbneObj.muparpz;
181  dk2nu->decay.mupare = g4lbneObj.mupare;
182 
183  dk2nu->decay.necm = g4lbneObj.Necm;
184  dk2nu->decay.nimpwt = g4lbneObj.Nimpwt;
185 
186  dk2nu->ppvx = g4lbneObj.ppvx;
187  dk2nu->ppvy = g4lbneObj.ppvy;
188  dk2nu->ppvz = g4lbneObj.ppvz;
189 
190  //not-in-new//dk2nu->xpoint = g4lbneObj.xpoint;
191  //not-in-new//dk2nu->ypoint = g4lbneObj.ypoint;
192  //not-in-new//dk2nu->zpoint = g4lbneObj.zpoint;
193 
194  dk2nu->tgtexit.tvx = g4lbneObj.tvx;
195  dk2nu->tgtexit.tvy = g4lbneObj.tvy;
196  dk2nu->tgtexit.tvz = g4lbneObj.tvz;
197  dk2nu->tgtexit.tpx = g4lbneObj.tpx;
198  dk2nu->tgtexit.tpy = g4lbneObj.tpy;
199  dk2nu->tgtexit.tpz = g4lbneObj.tpz;
200  dk2nu->tgtexit.tptype = ConvertGeantToPdg(g4lbneObj.tptype,"tptype");
201  dk2nu->tgtexit.tgen = g4lbneObj.tgen;
202 
203  //no-equiv// dk2nu->tgptype = ConvertGeantToPdg(g4lbneObj.tgptype,"tgptype");
204  //no-equiv// dk2nu->tgppx = g4lbneObj.tgppx;
205  //no-equiv// dk2nu->tgppy = g4lbneObj.tgppy;
206  //no-equiv// dk2nu->tgppz = g4lbneObj.tgppz;
207  //no-equiv// dk2nu->tprivx = g4lbneObj.tprivx;
208  //no-equiv// dk2nu->tprivy = g4lbneObj.tprivy;
209  //no-equiv// dk2nu->tprivz = g4lbneObj.tprivz;
210 
211  //not-in-new//dk2nu->beamx = g4lbneObj.protonX;
212  //not-in-new//dk2nu->beamy = g4lbneObj.protonY;
213  //not-in-new//dk2nu->beamz = g4lbneObj.protonZ;
214  //not-in-new//dk2nu->beampx = g4lbneObj.protonPx;
215  //not-in-new//dk2nu->beampy = g4lbneObj.protonPy;
216  //not-in-new//dk2nu->beampz = g4lbneObj.protonPz;
217 
218  //no-equiv// dk2nu->tgptype = g4lbneObj.tgptype;
219 
220 #ifdef G4LBNE_HAS_ANCESTOR_HISTORY
221  // now copy ancestor history
222  if ( g4lbneObj.overflow ) dk2nu->flagbits |= bsim::kFlgOverflow; // flag overflow
223  int nmx = TMath::Min(g4lbneObj.ntrajectory,10);
224  for (int ian=0; ian < nmx; ++ian) {
225 
226  const double mm2cm = 0.1;
227  const double mev2gev = 0.001;
228 
230  /*
231  dk2nu->apdg.push_back(g4lbneObj.pdg[ian]);
232  dk2nu->trackid.push_back(g4lbneObj.trackId[ian]);
233  dk2nu->parentid.push_back(g4lbneObj.parentId[ian]);
234 
235  dk2nu->startx.push_back(g4lbneObj.startx[ian]*mm2cm);
236  dk2nu->starty.push_back(g4lbneObj.starty[ian]*mm2cm);
237  dk2nu->startz.push_back(g4lbneObj.startz[ian]*mm2cm);
238  dk2nu->stopx.push_back(g4lbneObj.stopx[ian]*mm2cm);
239  dk2nu->stopy.push_back(g4lbneObj.stopy[ian]*mm2cm);
240  dk2nu->stopz.push_back(g4lbneObj.stopz[ian]*mm2cm);
241 
242  dk2nu->startpx.push_back(g4lbneObj.startpx[ian]*mev2gev);
243  dk2nu->startpy.push_back(g4lbneObj.startpy[ian]*mev2gev);
244  dk2nu->startpz.push_back(g4lbneObj.startpz[ian]*mev2gev);
245  dk2nu->stoppx.push_back(g4lbneObj.stoppx[ian]*mev2gev);
246  dk2nu->stoppy.push_back(g4lbneObj.stoppy[ian]*mev2gev);
247  dk2nu->stoppz.push_back(g4lbneObj.stoppz[ian]*mev2gev);
248 
249  dk2nu->pprodpx.push_back(g4lbneObj.pprodpx[ian]*mev2gev);
250  dk2nu->pprodpy.push_back(g4lbneObj.pprodpy[ian]*mev2gev);
251  dk2nu->pprodpz.push_back(g4lbneObj.pprodpz[ian]*mev2gev);
252 
253  // TString -> std::string via char*
254  dk2nu->proc.push_back(g4lbneObj.proc[ian].Data());
255  dk2nu->ivol.push_back(g4lbneObj.ivol[ian].Data());
256  dk2nu->fvol.push_back(g4lbneObj.fvol[ian].Data());
257  */
258 
259  ancestor.pdg = g4lbneObj.pdg[ian];
260  ancestor.SetStartXYZT(g4lbneObj.startx[ian]*mm2cm,
261  g4lbneObj.starty[ian]*mm2cm,
262  g4lbneObj.startz[ian]*mm2cm,
263  0); // don't have a time
264  ancestor.SetStartP(g4lbneObj.startpx[ian]*mev2gev,
265  g4lbneObj.startpy[ian]*mev2gev,
266  g4lbneObj.startpz[ian]*mev2gev);
267  ancestor.SetStopP(g4lbneObj.stoppx[ian]*mev2gev,
268  g4lbneObj.stoppy[ian]*mev2gev,
269  g4lbneObj.stoppz[ian]*mev2gev);
270  ancestor.SetPProdP(g4lbneObj.pprodpx[ian]*mev2gev,
271  g4lbneObj.pprodpy[ian]*mev2gev,
272  g4lbneObj.pprodpz[ian]*mev2gev);
273 
274  // ancestor.polx = ancestor.poly = ancestor.polz = ?
275  // ancestor.nucleau = ?
276  ancestor.proc = g4lbneObj.proc[ian].Data();
277  ancestor.ivol = g4lbneObj.ivol[ian].Data();
278  //not-lbne//ancestor.imat = g4lbneObj.imat[ian].Data();
279 
280  dk2nu->ancestor.push_back(ancestor);
281  }
282 #endif
283 
284 }
285 
286 void g4lbneCrossChecks(const g4lbne& g4lbneObj)
287 {
288 
289  static bool first = true;
290  static size_t indxMinosN = -1, indxMinosF = -1, indxNovaF = -1, indxMini = -1;
291  static size_t ioldMinosN = 0, ioldMinosF = 0, ioldNovaF = 1, ioldMini = 8;
292  if ( first ) {
293  first = false;
294  /// find indices for new weights so we can compare later
295  indxMinosN = find_loc_index("MINOS NearDet");
296  indxMinosF = find_loc_index("MINOS FarDet");
297  //indxNovaF = find_loc_index("NOvA FarDet");
298  indxNovaF = find_loc_index("Ash River per Lbne");
299  indxMini = find_loc_index("MiniBooNE");
300  cout << "indx MINOS Near " << indxMinosN << " Far " << indxMinosF << endl;
301  cout << "indx NOvA Far " << indxNovaF << endl;
302  cout << "indx MiniBooNE " << indxMini << endl;
303 
304  /**** lbne g4numi/src/NumiDataInput.cc
305  * RWH: these appear to be in meters ... not sure about what
306  * the "NovaX" are meant to represent
307  *
308  //Near & Far Detector location
309  nNear=11;//was 9 without the different energy for the ND positions.
310  nFar=2;
311  G4double xdetNear[] = { 0 , 0. , 7. , 11. ,
312  14. , 14. , 14. , 0. ,
313  25.84 , 4.8/2. , -4.8/2. };
314  G4double ydetNear[] = { 0 , -3. , -5. , -5. ,
315  -6. , -3. , 0. , 71. ,
316  78.42 , 3.8/2. , -3.8/2. };
317  G4double zdetNear[] = { 1040 , 1010. , 975. , 958. ,
318  940. , 840. , 740. , 940. ,
319  745.25 , 1040+16.6/2. , 1040-16.6/2. };
320  G4String detNameNear[] = { "Near", "Nova1a","Nova1b","Nova1c",
321  "Nova2a", "Nova2b", "Nova3","MSB",
322  "MiniBooNE","Near +x +y +z","Near -x -y -z"};
323  G4double xdetFar[] = {0 , 28.81258 };
324  G4double ydetFar[] = {0 , 81.39258 };
325  G4double zdetFar[] = {735000, 811400 };
326  G4String detNameFar[] = {"Far" , "Ash River"};
327  *
328  */
329  }
330  // cross check location energy/weights
331  //
332  histCompare(g4lbneObj.NenergyN[ioldMinosN],dk2nu->nuray[indxMinosN].E, "MINOS near energy");
333  histCompare(g4lbneObj.NWtNear[ioldMinosN], dk2nu->nuray[indxMinosN].wgt, "MINOS near wgt");
334  histCompare(g4lbneObj.NenergyF[ioldMinosF],dk2nu->nuray[indxMinosF].E, "MINOS far energy");
335  histCompare(g4lbneObj.NWtFar[ioldMinosF], dk2nu->nuray[indxMinosF].wgt, "MINOS far wgt");
336  if ( indxNovaF > 0 ) {
337  histCompare(g4lbneObj.NenergyF[ioldNovaF], dk2nu->nuray[indxNovaF].E, "NOvA far energy");
338  histCompare(g4lbneObj.NWtFar[ioldNovaF], dk2nu->nuray[indxNovaF].wgt, "NOvA far wgt");
339  } else {
340  static int nmsgnovaf = 10;
341  if ( nmsgnovaf-- > 0 ) {
342  cout << "## no relevant index for NOvA FarDet that will match lbne version of location" << endl;
343  }
344  }
345  histCompare(g4lbneObj.NenergyN[ioldMini], dk2nu->nuray[indxMini].E, "MiniBooNE energy");
346  histCompare(g4lbneObj.NWtNear[ioldMini], dk2nu->nuray[indxMini].wgt, "MiniBooNE wgt");
347 
348  static int nmsg = 0;
349  double a = g4lbneObj.NWtNear[ioldMini];
350  double b = dk2nu->nuray[indxMini].wgt;
351  double myfdiff = ( (a+b) > 1.0e-30 ) ? TMath::Abs(2.0*(a-b)/(a+b)) : 0;
352  if ( (nmsg < 1) && (myfdiff > 1.8) ) {
353  ++nmsg;
354  cout << "#### extreme case: mini wtg " << a << " " << b << " fdiff " << myfdiff << endl;
355  cout << "## g4 MinosN=" << ioldMinosN << " MinosF=" << ioldMinosF << " NovaF=" << ioldNovaF << " MiniBooNE=" << ioldMini << endl;
356  cout << "g4lbne NenergyN ";
357  for (int i=0; i<11; ++i) cout << ((i!=0&&i%4==0)?"\n ":" ")
358  << std::setw(12) << g4lbneObj.NenergyN[i];
359  cout << endl;
360  cout << "g4lbne NWtNear ";
361  for (int i=0; i<11; ++i) cout << ((i!=0&&i%4==0)?"\n ":" ")
362  << std::setw(12) << g4lbneObj.NWtNear[i];
363  cout << endl;
364  cout << "g4lbne NenergyF ";
365  for (int i=0; i<2; ++i) cout << " " << std::setw(12) << g4lbneObj.NenergyF[i];
366  cout << endl;
367  cout << "g4lbne NWtFar ";
368  for (int i=0; i<2; ++i) cout << " " << std::setw(12) << g4lbneObj.NWtFar[i];
369  cout << endl;
370  cout << "## dk MinosN=" << indxMinosN << " MinosF=" << indxMinosF << " NovaF=" << indxNovaF << " MiniBooNE=" << indxMini << endl;
371  cout << *dk2nu << endl;
372  }
373 
374 #ifdef G4LBNE_HAS_ANCESTOR_HISTORY
375  int nmx = TMath::Min(g4lbneObj.ntrajectory,10);
376  for (int ian=0; ian < nmx; ++ian) {
377  if ( ian <= 0 ) continue;
378  // except for first particle start[i] should == stop[i-1]
379  bool match = true;
380  match = match && ( histCompare(g4lbneObj.startx[ian],g4lbneObj.stopx[ian-1],"startx") );
381  match = match && ( histCompare(g4lbneObj.starty[ian],g4lbneObj.stopy[ian-1],"starty") );
382  match = match && ( histCompare(g4lbneObj.startz[ian],g4lbneObj.stopz[ian-1],"startz") );
383  match = match && ( g4lbneObj.ivol[ian] == g4lbneObj.fvol[ian-1] );
384  if ( ! match ) {
385  cout << "## ancestor " << ian << " didn't start where previous entry stopped" << endl;
386  }
387  }
388 #endif
389 
390 }
391 
TString fin
Definition: Style.C:24
Int_t Ndecay
Definition: g4lbne.h:54
Double_t ppdydz
% direction of nu parent at its production point
Definition: dk2nu.h:139
Double_t muparpx
%
Definition: dk2nu.h:160
Int_t job
identifying job #
Definition: dk2nu.h:323
Float_t Ndxdz
Definition: g4lbne.h:41
Int_t Ntype
Definition: g4lbne.h:55
Double_t pppz
% z momentum of nu parent at its production point
Definition: dk2nu.h:140
void SetStopP(Double_t px, Double_t py, Double_t pz)
TTree * dk2nuTree
std::string ivol
name of the volume where the particle starts
Definition: dk2nu.h:225
Double_t pdpx
% px momentum of nu parent at (vx,vy,vz)
Definition: dk2nu.h:133
void clear(const std::string &opt="")
reset everything
std::vector< bsim::NuRay > nuray
rays through detector fixed points
Definition: dk2nu.h:328
Float_t pppz
Definition: g4lbne.h:64
Int_t potnum
proton # processed by simulation
Definition: dk2nu.h:324
std::string proc
name of the process that creates this particle
Definition: dk2nu.h:224
int pots
Float_t Necm
Definition: g4lbne.h:76
virtual Long64_t LoadTree(Long64_t entry)
Double_t beamvwidth
vertical width of beam
Definition: dkmeta.h:110
Int_t tgen
Definition: dk2nu.h:276
Double_t NWtFar[3]
Definition: g4lbne.h:52
Double_t nimpwt
% cumulative importance weight prod to decay
Definition: dk2nu.h:166
Double_t beamhwidth
horizontal width of beam
Definition: dkmeta.h:109
Float_t ppvz
Definition: g4lbne.h:71
Float_t NenergyN[5]
Definition: g4lbne.h:47
Float_t tpx
Definition: g4lbne.h:84
Double_t necm
% nu energy in center-of-mass frame
Definition: dk2nu.h:165
double estimate_pots(int highest_potnum)
Double_t Nimpwt
Definition: g4lbne.h:77
bsim::Decay decay
basic decay information
Definition: dk2nu.h:327
int myjob
Double_t ppvz
production vertex z of nu parent
Definition: dk2nu.h:336
int ConvertGeantToPdg(int geant_code, std::string tag="?")
Float_t Vz
Definition: g4lbne.h:58
Double_t tvy
y position of nu ancestor as it exits target
Definition: dk2nu.h:270
string construct_outfilename(string infilename)
Int_t tgen
Definition: g4lbne.h:88
Int_t job
identifying job # (keep files distinct)
Definition: dkmeta.h:86
void SetStartXYZT(Double_t x, Double_t y, Double_t z, Double_t t)
set triplets
bsim::Dk2Nu * dk2nu
Double_t ppdxdz
% direction of nu parent at its production point
Definition: dk2nu.h:138
Double_t ppvy
production vertex y of nu parent
Definition: dk2nu.h:335
Double_t vz
% neutrino production vertex z
Definition: dk2nu.h:132
Float_t ppdxdz
Definition: g4lbne.h:62
bsim::TgtExit tgtexit
info about leaving the target
Definition: dk2nu.h:338
Float_t tvz
Definition: g4lbne.h:83
TTree * fChain
Definition: g4lbne.h:21
Float_t muparpz
Definition: g4lbne.h:74
Double_t beam0y
y of beam center at start
Definition: dkmeta.h:107
Double_t pots
protons-on-target
Definition: dkmeta.h:87
std::string tgtcfg
target config e.g. "minos/epoch3/-10cm"
Definition: dkmeta.h:98
Int_t tptype
species of ancestor exiting the target
Definition: dk2nu.h:275
Int_t pdg
ancestor species
Definition: dk2nu.h:188
Float_t hornCurrent
Definition: g4lbne.h:40
void ConvertFinish()
double frac_diff_tolerance
Double_t pdpz
% pz momentum of nu parent at (vx,vy,vz)
Definition: dk2nu.h:135
Double_t muparpy
%
Definition: dk2nu.h:161
Double_t ppvx
production vertex x of nu parent
Definition: dk2nu.h:334
Double_t NWtNear[5]
Definition: g4lbne.h:48
Double_t muparpz
%
Definition: dk2nu.h:162
Float_t pdPy
Definition: g4lbne.h:60
Double_t tpz
z momentum of nu ancestor as it exits target
Definition: dk2nu.h:274
Long64_t nentries
Float_t ppdydz
Definition: g4lbne.h:63
void SetStartP(Double_t px, Double_t py, Double_t pz)
std::string physics
e.g. "fluka08", "g4.9.3p01"
Definition: dkmeta.h:96
Double_t ppenergy
% energy of nu parent at its production point
Definition: dk2nu.h:141
Float_t muparpx
Definition: g4lbne.h:72
const double a
std::string beamsim
e.g. "flugg" or "g4numi/<tag>"
Definition: dkmeta.h:95
Double_t tpy
y momentum of nu ancestor as it exits target
Definition: dk2nu.h:273
Int_t ntype
% neutrino flavor (PDG? code)
Definition: dk2nu.h:128
Float_t beamVWidth
Definition: g4lbne.h:30
void ConvertBookNtuple(std::string ofname="test_dk2nu.root")
Float_t pdPx
Definition: g4lbne.h:59
double obs_frac_diff_max
Float_t Npz
Definition: g4lbne.h:43
Float_t mupare
Definition: g4lbne.h:75
void g4lbneCrossChecks(const g4lbne &g4lbneObj)
Float_t NenergyF[3]
Definition: g4lbne.h:51
OStream cout
Definition: OStream.cxx:6
Float_t Vx
Definition: g4lbne.h:56
Float_t ppmedium
Definition: g4lbne.h:66
void SetPProdP(Double_t px, Double_t py, Double_t pz)
Int_t tptype
Definition: g4lbne.h:87
std::string horncfg
horn config e.g. "FHC/185A/LE/h1xoff=1mm"
Definition: dkmeta.h:99
Double_t tvz
z position of nu ancestor as it exits target
Definition: dk2nu.h:271
Float_t beamHWidth
Definition: g4lbne.h:29
bool histCompare(double a, double b, string s)
bsim::DkMeta * dkmeta
Float_t ppenergy
Definition: g4lbne.h:65
Double_t beam0x
x of beam center at start
Definition: dkmeta.h:106
Float_t nuTarZ
Definition: g4lbne.h:39
Float_t Vy
Definition: g4lbne.h:57
Long64_t nbytes
Int_t ptype
Definition: g4lbne.h:67
Float_t tvx
Definition: g4lbne.h:81
Definition: g4lbne.h:19
int Convert5xToPdg(int old_ntype)
Float_t tpz
Definition: g4lbne.h:86
size_t find_loc_index(string match)
Double_t tpx
x momentum of nu ancestor as it exits target
Definition: dk2nu.h:272
const hit & b
Definition: hits.cxx:21
Float_t Ndydz
Definition: g4lbne.h:42
Int_t Norig
Definition: g4lbne.h:53
Double_t vx
% neutrino production vertex x
Definition: dk2nu.h:130
void copy_g4lbne_to_dk2nu(const g4lbne &g4lbneObj)
void calcLocationWeights(const bsim::DkMeta *dkmeta, bsim::Dk2Nu *dk2nu)
user interface
Float_t tpy
Definition: g4lbne.h:85
Float_t beamY
Definition: g4lbne.h:32
Float_t muparpy
Definition: g4lbne.h:73
Int_t ndecay
decay process (see dkproc_t)
Definition: dk2nu.h:127
void Format(TGraph *gr, int lcol, int lsty, int lwid, int mcol, int msty, double msiz)
Definition: Style.cxx:154
Float_t tvy
Definition: g4lbne.h:82
Double_t mupare
% energy of nu grandparent
Definition: dk2nu.h:163
Double_t vy
% neutrino production vertex y
Definition: dk2nu.h:131
static const double nb
Definition: Units.h:89
Int_t flagbits
Definition: dk2nu.h:345
Float_t e
Definition: plot.C:35
Float_t ppvy
Definition: g4lbne.h:70
Int_t norig
not used?
Definition: dk2nu.h:126
Int_t ppmedium
material nu parent was produced in
Definition: dk2nu.h:143
Float_t pdPz
Definition: g4lbne.h:61
Double_t tvx
x position of nu ancestor as it exits target
Definition: dk2nu.h:269
Int_t evtno
Definition: g4lbne.h:27
std::vector< bsim::Ancestor > ancestor
chain from proton to neutrino
Definition: dk2nu.h:329
Double_t pdpy
% py momentum of nu parent at (vx,vy,vz)
Definition: dk2nu.h:134
Float_t beamX
Definition: g4lbne.h:31
Float_t Nenergy
Definition: g4lbne.h:44
Definition: nuray.py:1
void convert_g4lbne(string ifname="../fluxfiles/generic_g4lbne.root", int jobnum=42, string locfile="${DK2NU}/etc/LBNElocations.txt", Long64_t maxentries=-1, Long64_t moddump=-1)
Float_t ppvx
Definition: g4lbne.h:69
Int_t ptype
% nu parent species (PDG? code)
Definition: dk2nu.h:144
void ConvertInit(std::string locfilename="$(DK2NU)/etc/locations.txt")