NCAna_module.cc
Go to the documentation of this file.
1 ///
2 /// \brief
3 /// \author Xuebing Bu - xbbu@fnal.gov
4 ///
5 //#include "RecoJMShower/JMShower.h"
6 //#include "RecoJMShower/EID.h"
7 #include "RecVarPID/RVP.h"
8 #include "LEM/PIDDetails.h"
9 #include "ShowerLID/EventLID.h"
10 #include "ShowerLID/ShowerLID.h"
11 #include "Calibrator/Calibrator.h"
12 
13 #include "Geometry/Geometry.h"
14 #include "GeometryObjects/CellGeo.h"
15 
16 // ROOT includes
17 #include "TVector3.h"
18 #include "TH1F.h"
19 #include "TMath.h"
20 #include "TTree.h"
21 #include "TProfile2D.h"
22 #include "time.h"
23 
24 #include "nusimdata/SimulationBase/MCFlux.h"
25 #include "nusimdata/SimulationBase/MCTruth.h"
26 #include "Simulation/GENIEReweightTable.h"
27 #include "Simulation/Particle.h"
28 #include "Simulation/ParticleNavigator.h"
29 #include "SummaryData/POTSum.h"
30 #include "SummaryData/SpillData.h"
31 #include "RecoBase/CellHit.h"
32 #include "RecoBase/Cluster.h"
33 #include "RecoBase/RecoHit.h"
34 #include "RecoBase/Track.h"
35 #include "RecoBase/Vertex.h"
36 #include "RecoBase/Prong.h"
37 #include "RecoBase/HoughResult.h"
38 #include "RecoBase/Shower.h"
39 
40 #include "Utilities/AssociationUtil.h"
41 #include "MCCheater/BackTracker.h"
42 #include "Utilities/func/MathUtil.h"
43 
44 // Framework includes
45 #include "art/Framework/Core/EDAnalyzer.h"
46 #include "art/Framework/Core/ModuleMacros.h"
47 #include "canvas/Persistency/Common/PtrVector.h"
48 #include "canvas/Persistency/Common/Assns.h"
49 #include "canvas/Persistency/Common/FindManyP.h"
50 #include "canvas/Persistency/Common/FindMany.h"
51 #include "canvas/Persistency/Common/FindOneP.h"
52 #include "canvas/Persistency/Common/FindOne.h"
53 #include "art/Framework/Principal/Event.h"
54 #include "art/Framework/Principal/SubRun.h"
55 #include "fhiclcpp/ParameterSet.h"
56 #include "art/Framework/Principal/Handle.h"
57 #include "art/Framework/Services/Registry/ServiceHandle.h"
58 #include "art/Framework/Services/Optional/TFileService.h"
59 #include "art/Framework/Services/Optional/TFileDirectory.h"
60 #include "messagefacility/MessageLogger/MessageLogger.h"
61 
62 #include <cmath>
63 #include <vector>
64 #include <string>
65 #include <vector>
66 
67 //include TMVA
68 #include "TMVA/Tools.h"
69 #include "TMVA/Reader.h"
70 #include "TMVA/MethodCuts.h"
71 
72 class TVector3;
73 class TH1F;
74 class TTree;
75 
76 namespace rb{class Cluster; class Track;}
77 namespace sim{class Particle;}
78 namespace simb{class MCFlux; class MCTruth; class MCNeutrino;}
79 
80 namespace ncs {
81  class NCAna : public art::EDAnalyzer {
82  public:
83  explicit NCAna(fhicl::ParameterSet const& pset);
84  virtual ~NCAna();
85 
86  void beginJob();
87  void analyze(const art::Event& evt);
88  void reconfigure(const fhicl::ParameterSet& p);
89  void endSubRun(const art::SubRun& sr);
90  void endJob();
91  double SimpleOscProb(const simb::MCFlux& flux,
92  const simb::MCNeutrino& nu) const;
93  double fTotalPOT;
94  int fTotalSpill;
95  TTree *_btree;
96 
97  private:
98 
99  const rb::RecoHit MakeRecoHit(art::Ptr<rb::CellHit> const& chit,
100  double const& w);
101 
102  //folders to get objects from
103  std::string fSliceLabel; ///label for slices
104  //std::string fDTrackLabel; ///label for discrete tracks
105  std::string fKTrackLabel; ///label for kalman tracks
106  std::string fProngLabel; ///label for prong
107  std::string fInstLabel3D; ///instance label for prongs 3D
108  std::string fInstLabel2D; ///instance label for prongs 2D
109  std::string fMCFluxLabel; ///label for generator information
110  std::string fVertexLabel; ///label for vertex
111  std::string fCellHitLabel; ///label for cell hits
112  //std::string fShowerLabel; ///label for jmshower
113  std::string frvpLabel; ///label for rvp
114  std::string flidLabel;
115  std::string flemLabel;
116  std::string fDataSpillLabel; //label for data spill
117  std::string fweightLabel;//genie weights
118 
119  int isMC; //mode for MC
120  int isROCK; //mode for ROCK muon
121  int is_nue;//switch to select unoscillated nue only events
122 
123  //save event info to ntuples
124  int runnum;
125  int subrunnum;
126  int evtnum;
127  int Nslices;
128  double spillPOT;
129  int goodSpill;
130  double horncurrent;
131  double posX;
132  double posY;
133  double widthX;
134  double widthY;
135  int deltaspilltimensec;
136 
137  int evt_year;
138  int evt_month;
139  int evt_day;
140  int evt_hour;
141  int evt_min;
142  int evt_sec;
143 
144  int Ndeaddcm;
145  std::string deaddcm_location[14];
146  int Ndb1apd_24cells;
147  int Ndb2apd_24cells;
148  int Ndb3apd_24cells;
149  int Nmcapd_24cells;
150  int Ndb1apd_28cells;
151  int Ndb2apd_28cells;
152  int Ndb3apd_28cells;
153  int Nmcapd_28cells;
154  int Nouttimehits;
155  int Nintimehits;
156  float Fouttimehits_noisyapd;
157 
158  TTree *_otree;
159  //TTree *_evttree;
160  int sl_runnum;
161  int sl_subrunnum;
162  int sl_evtnum;
163  int sl_Nslices;
164  int sl_goodSpill;
165  double sl_spillPOT;
166  double sl_horncurrent;
167  double sl_posX;
168  double sl_posY;
169  double sl_widthX;
170  double sl_widthY;
171  int sl_deltaspilltimensec;
172 
173  int sl_Ndeaddcm;
174  std::string sl_deaddcm_location[14];
175  int sl_Ndb1apd_24cells;
176  int sl_Ndb2apd_24cells;
177  int sl_Ndb3apd_24cells;
178  int sl_Nmcapd_24cells;
179  int sl_Ndb1apd_28cells;
180  int sl_Ndb2apd_28cells;
181  int sl_Ndb3apd_28cells;
182  int sl_Nmcapd_28cells;
183  int sl_Nouttimehits;
184  int sl_Nintimehits;
185  float sl_Fouttimehits_noisyapd;
186 
187  /*** Vertex ***/
188  float Vx;//true vertex X
189  float Vy;//true vertex Y
190  float Vz;//true vertex Z
191  float recoVx;//reco vertex X
192  float recoVy;//reco vertex Y
193  float recoVz;//reco vertex Z
194 
195  /*** true info ***/
196  float truePt;//true Pt
197  float trueE;//true E
198  float trueTheta;//true theta
199  int ccnc;
200  int PDG;
201  int origPDG;
202  int mode;
203  bool isQE;
204  int intType;
205  int hitNucl;
206  float HadX;
207  float HadY;
208  float HadW;
209  float qsqr;
210  float OscWt;
211  int lepPDG;
212  float lepPx;
213  float lepPy;
214  float lepPz;
215  float lepE;
216  int mesonPDG;
217  float mesonPx;
218  float mesonPy;
219  float mesonPz;
220 
221  //GENIE weights
222  int Nweights;
223  float plus1sigma[68];
224  float plus2sigma[68];
225  float minus1sigma[68];
226  float minus2sigma[68];
227 
228  //for pi0->2gamma from most energetic pi0
229  int HasP2GG;//has pi0->2gamma
230  float Epi0;
231  float gam1E;
232  float gam1Px;
233  float gam1Py;
234  float gam1Pz;
235  float gam2E;
236  float gam2Px;
237  float gam2Py;
238  float gam2Pz;
239 
240  //all pi0, first 20 if possible ?
241  int Npi0;
242  float pi0Px[20];
243  float pi0Py[20];
244  float pi0Pz[20];
245  float pi0E[20];
246  float pi0Vx[20];
247  float pi0Vy[20];
248  float pi0Vz[20];
249 
250  /*** slicer information ***/
251  int Ncells;//# of associate cells
252  int Nplanes;//# of extended planes
253  float recoE;//reco energy
254  float PlaneEnergy[192];//set the maximal plane # to 192
255  int MinPlane;//min plane index
256  int PlaneNcells[192];//# of cells in each plane
257  int sliceIndex;
258  float slice_meantns;
259  float slice_mintns;
260  float slice_maxtns;
261 
262  /*** longest kalman track ***/
263  float KtrackLength;//length of longest track
264  int KtrackNcells;//# of cells associated to longest track
265  float KtrackEnergy;//energy associated to longest track
266  float KtrackNcells_ratio;//# of cells associated to longest track over Ncells
267  float KtrackEnergy_ratio;//energy associated to longest track over recoE
268  float KtrackNplanes;
269  float KtrackTheta;
270  float KtrackPhi;
271  float KtrackStartX;
272  float KtrackStartY;
273  float KtrackStartZ;
274  float KtrackStopX;
275  float KtrackStopY;
276  float KtrackStopZ;
277  int NKtracks;//# of tracks
278  int Ktrack3D;
279 
280  int trackcellNcell;
281  float trackcellX[500];
282  float trackcellY[500];
283  float trackcellZ[500];
284  float trackcellE[500];
285  int trackcellADC[500];
286  float trackcellTNS[500];
287  int trackcellPlane[500];
288  int trackcellCell[500];
289  int trackcellView[500];
290  float trackcellPE[500];
291  float trackcellPECorr[500];
292  float trackcellW[500];
293 
294  /*** mip information ***/
295  float Nmip;//# of mip cells
296  float Fmip;//# of mip cells over Ncells
297 
298  /*** cell information ***/
299  int slcellNcell;//# of good cells
300  float slcellX[500];
301  float slcellY[500];
302  float slcellZ[500];
303  float slcellE[500];
304  int slcellADC[500];
305  int slcellPlane[500];
306  int slcellCell[500];
307  float slcellPECorr[500];
308  float slcellW[500];
309  float slcellTNS[500];
310 
311  /*** prong information ***/
312  int Nprongs;
313  float Ebalance;
314  float Dphi;
315 
316  float prongEnergy[20];
317  float prongPhi[20];
318  float prongTheta[20];
319  float prongStartX[20];
320  float prongStartY[20];
321  float prongStartZ[20];
322  float prongLength[20];
323  float prongStopX[20];
324  float prongStopY[20];
325  float prongStopZ[20];
326  float prongEnergyXView[20];//energy for XView cells
327  float prongEnergyYView[20];//energy for YView cells
328  int prongNCells[20];//# of cells
329  int prongNCellsXView[20];//# of cells for XView
330  int prongNCellsYView[20];//# of cells for YView
331  int prongNplanes[20];
332 
333  int allprongcellNcell;
334  int allprongIndex[500];
335  float allprongcellX[500];
336  float allprongcellY[500];
337  float allprongcellZ[500];
338  float allprongcellE[500];
339  int allprongcellADC[500];
340  int allprongcellPlane[500];
341  int allprongcellCell[500];
342  float allprongcellPECorr[500];
343  float allprongcellW[500];
344 
345 /*
346  int q1prongcellNcell;
347  float q1prongcellX[500];
348  float q1prongcellY[500];
349  float q1prongcellZ[500];
350  float q1prongcellE[500];
351  int q1prongcellADC[500];
352  int q1prongcellPlane[500];
353  int q1prongcellCell[500];
354  float q1prongcellPECorr[500];
355  float q1prongcellW[500];
356 
357  int q2prongcellNcell;
358  float q2prongcellX[500];
359  float q2prongcellY[500];
360  float q2prongcellZ[500];
361  float q2prongcellE[500];
362  int q2prongcellADC[500];
363  int q2prongcellPlane[500];
364  int q2prongcellCell[500];
365  float q2prongcellPECorr[500];
366  float q2prongcellW[500];
367 */
368  //leading prong
369  int prong1cellNcell;
370  float prong1Fmip;
371  float prong1cellX[500];
372  float prong1cellY[500];
373  float prong1cellZ[500];
374  float prong1cellE[500];
375  int prong1cellADC[500];
376  int prong1cellPlane[500];
377  int prong1cellCell[500];
378  float prong1cellPECorr[500];
379  float prong1cellW[500];
380 
381 /*
382  //sub-leading prong
383  int prong2cellNcell;
384  float prong2cellX[500];
385  float prong2cellY[500];
386  float prong2cellZ[500];
387  float prong2cellE[500];
388  int prong2cellADC[500];
389  float prong2cellTNS[500];
390  int prong2cellPlane[500];
391  int prong2cellCell[500];
392  int prong2cellView[500];
393  float prong2cellPE[500];
394  float prong2cellPECorr[500];
395  float prong2cellW[500];
396 */
397 
398  //2D prongs
399  int Nprongs2D;
400  int prong2DInXView[20];//prong in X View or not
401  float prongEnergy2D[20];//energy
402  int prongNCells2D[20];//# of cells
403  float prongPhi2D[20];
404  float prongTheta2D[20];
405  float prongStartX2D[20];
406  float prongStartY2D[20];
407  float prongStartZ2D[20];
408  float prongLength2D[20];
409  float prongStopX2D[20];
410  float prongStopY2D[20];
411  float prongStopZ2D[20];
412  float Ebalance2D;
413 
414  int allprong2DcellNcell;
415  int allprong2DIndex[500];
416  float allprong2DcellX[500];
417  float allprong2DcellY[500];
418  float allprong2DcellZ[500];
419  float allprong2DcellE[500];
420  int allprong2DcellADC[500];
421  int allprong2DcellPlane[500];
422  int allprong2DcellCell[500];
423  float allprong2DcellPECorr[500];
424 
425  float rvp;
426  float lid;
427  float lide;
428  float lem;
429  int Nshowers;
430  float shwLID[20];
431  int shw3D[20];
432  int shwNCells[20];
433  int shwNCellsXview[20];
434  int shwNCellsYview[20];
435  float shwLength[20];
436  float shwStartX[20];
437  float shwStartY[20];
438  float shwStartZ[20];
439  float shwStopX[20];
440  float shwStopY[20];
441  float shwStopZ[20];
442  float shwTheta[20];
443  float shwEcells[20];
444  float shwEnu[20];
445 
446  int shwcellNcell;
447  int shwIndex[500];
448  float shwcellX[500];
449  float shwcellY[500];
450  float shwcellZ[500];
451  float shwcellE[500];
452  int shwcellADC[500];
453  int shwcellPlane[500];
454  int shwcellCell[500];
455  float shwcellPECorr[500];
456  float shwcellW[500];
457  };
458 }
459 
460 namespace ncs{
461 
462  /*
463  //BDT reader
464  TMVA::Reader *_bdt_ncana=new TMVA::Reader();
465  //input variable values
466  float inputvars_ncana[12]={0.};
467  //input variable names
468  TString inputs_ncana[12]={"_ncells","_energy","_ltrack","_ltrack_ncells_ratio",
469  "_fmip","_nmip","_efrac_20planes","_efrac_slide2",
470  "_efrac_slide6","_efrac_2sigma","_nprongs","_ebalance"};
471  */
472 
473 
474  //......................................................................
475  NCAna::NCAna(fhicl::ParameterSet const& pset)
476  : EDAnalyzer(pset)
477  {
478  mf::LogInfo("NCAna")<<__PRETTY_FUNCTION__<<"\n";
479  reconfigure(pset);
480  }
481 
482  //......................................................................
483  NCAna::~NCAna()
484  {
485  }
486 
487  //......................................................................
488  void NCAna::reconfigure(const fhicl::ParameterSet& pset)
489  {
490  //inputPDG = pset.get<int>("pdg");
491  fSliceLabel = pset.get< std::string >("SliceLabel");
492  //fDTrackLabel = pset.get< std::string >("DTrackLabel");
493  fKTrackLabel = pset.get< std::string >("KTrackLabel");
494  fProngLabel = pset.get< std::string >("ProngLabel");
495  fInstLabel3D = pset.get< std::string >("InstLabel3D");
496  fInstLabel2D = pset.get< std::string >("InstLabel2D");
497  fMCFluxLabel = pset.get< std::string >("MCFluxLabel");
498  fCellHitLabel = pset.get< std::string >("CellHitLabel");
499  fVertexLabel = pset.get< std::string >("VertexLabel");
500  //fShowerLabel = pset.get< std::string >("ShowerLabel");
501  fweightLabel = pset.get< std::string >("weightLabel");
502  //
503  frvpLabel = pset.get< std::string >("rvpLabel");
504  flidLabel = pset.get< std::string >("lidLabel");
505  flemLabel = pset.get< std::string >("lemLabel");
506 
507  isMC = pset.get<int>("MCevts");
508  isROCK = pset.get<int>("ROCKevts");
509  is_nue = pset.get<int>("unoscNue");
510  fDataSpillLabel = pset.get<std::string> ("DataSpillLabel");
511  }
512 
513  //......................................................................
514  void NCAna::beginJob()
515  {
516  art::ServiceHandle<art::TFileService> tfs;
517 
518  fTotalPOT=0.;
519  fTotalSpill=0;
520 
521  //for( int i=0; i<12; ++i ){
522  // _bdt_ncana->AddVariable(inputs_ncana[i],&inputvars_ncana[i]);
523  // }
524 
525  //_bdt_ncana->BookMVA("BDT method","BDT_12vars_20121102.xml");
526 
527  //histo
528  //Ereco=tfs->make<TH1F>("Ereco","Ereco;Ereco;# of events",100,0.,5.);
529  //personal tuples
530  _otree = tfs->make<TTree>("NCAna","NCAna particle");
531 
532  _otree->Branch("sl_runnum",&sl_runnum,"sl_runnum/I");
533  _otree->Branch("sl_subrunnum",&sl_subrunnum,"sl_subrunnum/I");
534  _otree->Branch("sl_evtnum",&sl_evtnum,"sl_evtnum/I");
535  _otree->Branch("sl_Nslices",&sl_Nslices,"sl_Nslices/I");
536  _otree->Branch("sl_goodSpill",&sl_goodSpill,"sl_goodSpill/I");
537  _otree->Branch("sl_spillPOT",&sl_spillPOT,"sl_spillPOT/D");
538  _otree->Branch("sl_orncurrent",&sl_horncurrent,"sl_horncurrent/D");
539  _otree->Branch("sl_posX",&sl_posX,"sl_posX/D");
540  _otree->Branch("sl_posY",&sl_posY,"sl_posY/D");
541  _otree->Branch("sl_widthX",&sl_widthX,"sl_widthX/D");
542  _otree->Branch("sl_widthY",&sl_widthY,"sl_widthY/D");
543  _otree->Branch("sl_deltaspilltimensec",&sl_deltaspilltimensec,"sl_deltaspilltimensec/I");
544 
545  _otree->Branch("sl_Ndeaddcm",&sl_Ndeaddcm,"sl_Ndeaddcm/I");
546  _otree->Branch("sl_deaddcm_location",&sl_deaddcm_location,"sl_deaddcm_location[sl_Ndeaddcm]/S");
547  _otree->Branch("sl_Ndb1apd_24cells",&sl_Ndb1apd_24cells,"sl_Ndb1apd_24cells/I");
548  _otree->Branch("sl_Ndb2apd_24cells",&sl_Ndb2apd_24cells,"sl_Ndb2apd_24cells/I");
549  _otree->Branch("sl_Ndb3apd_24cells",&sl_Ndb3apd_24cells,"sl_Ndb3apd_24cells/I");
550  _otree->Branch("sl_Nmcapd_24cells",&sl_Nmcapd_24cells,"sl_Nmcapd_24cells/I");
551  _otree->Branch("sl_Ndb1apd_28cells",&sl_Ndb1apd_28cells,"sl_Ndb1apd_28cells/I");
552  _otree->Branch("sl_Ndb2apd_28cells",&sl_Ndb2apd_28cells,"sl_Ndb2apd_28cells/I");
553  _otree->Branch("sl_Ndb3apd_28cells",&sl_Ndb3apd_28cells,"sl_Ndb3apd_28cells/I");
554  _otree->Branch("sl_Nmcapd_28cells",&sl_Nmcapd_28cells,"sl_Nmcapd_28cells/I");
555  _otree->Branch("sl_Nintimehits",&sl_Nintimehits,"sl_Nintimehits/I");
556  _otree->Branch("sl_Nouttimehits",&sl_Nouttimehits,"sl_Nouttimehits/I");
557  _otree->Branch("sl_Fouttimehits_noisyapd",&sl_Fouttimehits_noisyapd,"sl_Fouttimehits_noisyapd/F");
558 
559  _otree->Branch("Vx",&Vx,"Vx/F");
560  _otree->Branch("Vy",&Vy,"Vy/F");
561  _otree->Branch("Vz",&Vz,"Vz/F");
562  _otree->Branch("recoVx",&recoVx,"recoVx/F");
563  _otree->Branch("recoVy",&recoVy,"recoVy/F");
564  _otree->Branch("recoVz",&recoVz,"recoVz/F");
565 
566  _otree->Branch("truePt",&truePt,"truePt/F");
567  _otree->Branch("trueE",&trueE,"trueE/F");
568  _otree->Branch("trueTheta",&trueTheta,"trueTheta/F");
569  _otree->Branch("ccnc",&ccnc,"ccnc/I");
570  _otree->Branch("PDG",&PDG,"PDG/I");
571  _otree->Branch("origPDG",&origPDG,"origPDG/I");
572  _otree->Branch("mode",&mode,"mode/I");
573  _otree->Branch("isQE",&isQE,"isQE/B");
574  _otree->Branch("intType",&intType,"intType/I");
575  _otree->Branch("hitNucl",&hitNucl,"hitNucl/I");
576  _otree->Branch("HadX",&HadX,"HadX/F");
577  _otree->Branch("HadY",&HadY,"HadY/F");
578  _otree->Branch("HadW",&HadW,"HadW/F");
579  _otree->Branch("qsqr",&qsqr,"qsqr/F");
580  _otree->Branch("OscWt",&OscWt,"OscWt/F");
581  _otree->Branch("lepE",&lepE,"lepE/F");
582  _otree->Branch("lepPx",&lepPx,"lepPx/F");
583  _otree->Branch("lepPy",&lepPy,"lepPy/F");
584  _otree->Branch("lepPz",&lepPz,"lepPz/F");
585  _otree->Branch("lepPDG",&lepPDG,"lepPDG/I");
586  _otree->Branch("mesonPDG",&mesonPDG,"mesonPDG/I");
587  _otree->Branch("mesonPx",&mesonPx,"mesonPx/F");
588  _otree->Branch("mesonPy",&mesonPy,"mesonPy/F");
589  _otree->Branch("mesonPz",&mesonPz,"mesonPz/F");
590 
591  _otree->Branch("Nweights",&Nweights,"Nweights/I");
592  _otree->Branch("plus1sigma",&plus1sigma,"plus1sigma[Nweights]/F");
593  _otree->Branch("plus2sigma",&plus2sigma,"plus2sigma[Nweights]/F");
594  _otree->Branch("minus1sigma",&minus1sigma,"minus1sigma[Nweights]/F");
595  _otree->Branch("minus2sigma",&minus2sigma,"minus2sigma[Nweights]/F");
596 
597  _otree->Branch("Epi0",&Epi0,"Epi0/F");
598  _otree->Branch("gam1E",&gam1E,"gam1E/F");
599  _otree->Branch("gam1Px",&gam1Px,"gam1Px/F");
600  _otree->Branch("gam1Py",&gam1Py,"gam1Py/F");
601  _otree->Branch("gam1Pz",&gam1Pz,"gam1Pz/F");
602  _otree->Branch("HasP2GG",&HasP2GG,"HasP2GG/I");
603 
604  _otree->Branch("gam2E",&gam2E,"gam2E/F");
605  _otree->Branch("gam2Px",&gam2Px,"gam2Px/F");
606  _otree->Branch("gam2Py",&gam2Py,"gam2Py/F");
607  _otree->Branch("gam2Pz",&gam2Pz,"gam2Pz/F");
608 
609  _otree->Branch("Npi0",&Npi0,"Npi0/I");
610  _otree->Branch("pi0Px",pi0Px,"pi0Px[Npi0]/F");
611  _otree->Branch("pi0Py",pi0Py,"pi0Py[Npi0]/F");
612  _otree->Branch("pi0Pz",pi0Pz,"pi0Pz[Npi0]/F");
613  _otree->Branch("pi0E",pi0E,"pi0E[Npi0]/F");
614  _otree->Branch("pi0Vx",pi0Vx,"pi0Vx[Npi0]/F");
615  _otree->Branch("pi0Vy",pi0Vy,"pi0Vy[Npi0]/F");
616  _otree->Branch("pi0Vz",pi0Vz,"pi0Vz[Npi0]/F");
617 
618  _otree->Branch("slcellNcell",&slcellNcell,"slcellNcell/I");
619  _otree->Branch("slcellX",&slcellX,"slcellX[slcellNcell]/F");
620  _otree->Branch("slcellY",&slcellY,"slcellY[slcellNcell]/F");
621  _otree->Branch("slcellZ",&slcellZ,"slcellZ[slcellNcell]/F");
622  _otree->Branch("slcellE",&slcellE,"slcellE[slcellNcell]/F");
623  _otree->Branch("slcellADC",&slcellADC,"slcellADC[slcellNcell]/I");
624  _otree->Branch("slcellPlane",&slcellPlane,"slcellPlane[slcellNcell]/I");
625  _otree->Branch("slcellCell",&slcellCell,"slcellCell[slcellNcell]/I");
626  _otree->Branch("slcellW",&slcellW,"slcellW[slcellNcell]/F");
627  _otree->Branch("slcellTNS",&slcellTNS,"slcellTNS[slcellNcell]/F");
628  _otree->Branch("slcellPECorr",&slcellPECorr,"slcellPECorr[slcellNcell]/F");
629 
630  _otree->Branch("NKtracks",&NKtracks,"NKtracks/I");
631  _otree->Branch("KtrackLength",&KtrackLength,"KtrackLength/F");
632  _otree->Branch("KtrackNcells",&KtrackNcells,"KtrackNcells/I");
633  _otree->Branch("KtrackEnergy",&KtrackEnergy,"KtrackEnergy/F");
634  _otree->Branch("KtrackNcells_ratio",&KtrackNcells_ratio,"KtrackNcells_ratio/F");
635  _otree->Branch("KtrackEnergy_ratio",&KtrackEnergy_ratio,"KtrackEnergy_ratio/F");
636  _otree->Branch("KtrackNplanes",&KtrackNplanes,"KtrackNplanes/F");
637  _otree->Branch("KtrackPhi",&KtrackPhi,"KtrackPhi/F");
638  _otree->Branch("KtrackTheta",&KtrackTheta,"KtrackTheta/F");
639  _otree->Branch("KtrackStartX",&KtrackStartX,"KtrackStartX/F");
640  _otree->Branch("KtrackStartY",&KtrackStartY,"KtrackStartY/F");
641  _otree->Branch("KtrackStartZ",&KtrackStartZ,"KtrackStartZ/F");
642  _otree->Branch("KtrackStopX",&KtrackStopX,"KtrackStopX/F");
643  _otree->Branch("KtrackStopY",&KtrackStopY,"KtrackStopY/F");
644  _otree->Branch("KtrackStopZ",&KtrackStopZ,"KtrackStopZ/F");
645  _otree->Branch("Ktrack3D",&Ktrack3D,"Ktrack3D/I");
646 
647  _otree->Branch("trackcellNcell",&trackcellNcell,"trackcellNcell/I");
648  _otree->Branch("trackcellX",&trackcellX,"trackcellX[trackcellNcell]/F");
649  _otree->Branch("trackcellY",&trackcellY,"trackcellY[trackcellNcell]/F");
650  _otree->Branch("trackcellZ",&trackcellZ,"trackcellZ[trackcellNcell]/F");
651  _otree->Branch("trackcellE",&trackcellE,"trackcellE[trackcellNcell]/F");
652  _otree->Branch("trackcellADC",&trackcellADC,"trackcellADC[trackcellNcell]/I");
653  _otree->Branch("trackcellPlane",&trackcellPlane,"trackcellPlane[trackcellNcell]/I");
654  _otree->Branch("trackcellCell",&trackcellCell,"trackcellCell[trackcellNcell]/I");
655  _otree->Branch("trackcellW",&trackcellW,"trackcellW[trackcellNcell]/F");
656  _otree->Branch("trackcellPECorr",&trackcellPECorr,"trackcellPECorr[trackcellNcell]/F");
657 
658  _otree->Branch("Ncells",&Ncells,"Ncells/I");
659  _otree->Branch("Nplanes",&Nplanes,"Nplanes/I");
660  _otree->Branch("recoE",&recoE,"recoE/F");
661  _otree->Branch("MinPlane",&MinPlane,"MinPlane/I");
662  _otree->Branch("sliceIndex",&sliceIndex,"sliceIndex/I");
663  _otree->Branch("slice_meantns",&slice_meantns,"slice_meantns/F");
664  _otree->Branch("slice_mintns",&slice_mintns,"slice_mintns/F");
665  _otree->Branch("slice_maxtns",&slice_maxtns,"slice_maxtns/F");
666 
667  _otree->Branch("Nmip",&Nmip,"Nmip/F");
668  _otree->Branch("Fmip",&Fmip,"Fmip/F");
669 
670  _otree->Branch("Nprongs",&Nprongs,"Nprongs/I");
671  _otree->Branch("Ebalance",&Ebalance,"Ebalance/F");
672  _otree->Branch("Dphi",&Dphi,"Dphi/F");
673 
674  _otree->Branch("prongLength",&prongLength,"prongLength[Nprongs]/F");
675  _otree->Branch("prongStopX",&prongStopX,"prongStopX[Nprongs]/F");
676  _otree->Branch("prongStopY",&prongStopY,"prongStopY[Nprongs]/F");
677  _otree->Branch("prongStopZ",&prongStopZ,"prongStopZ[Nprongs]/F");
678 
679  _otree->Branch("prongEnergy",&prongEnergy,"prongEnergy[Nprongs]/F");
680  _otree->Branch("prongPhi",&prongPhi,"prongPhi[Nprongs]/F");
681  _otree->Branch("prongTheta",&prongTheta,"prongTheta[Nprongs]/F");
682  _otree->Branch("prongStartX",&prongStartX,"prongStartX[Nprongs]/F");
683  _otree->Branch("prongStartY",&prongStartY,"prongStartY[Nprongs]/F");
684  _otree->Branch("prongStartZ",&prongStartZ,"prongStartZ[Nprongs]/F");
685  _otree->Branch("prongEnergyXView",&prongEnergyXView,"prongEnergyXView[Nprongs]/F");
686  _otree->Branch("prongEnergyYView",&prongEnergyYView,"prongEnergyYView[Nprongs]/F");
687  _otree->Branch("prongNCells",&prongNCells,"prongNCells[Nprongs]/I");
688  _otree->Branch("prongNCellsXView",&prongNCellsXView,"prongNCellsXView[Nprongs]/I");
689  _otree->Branch("prongNCellsYView",&prongNCellsYView,"prongNCellsYView[Nprongs]/I");
690  _otree->Branch("prongNplanes",&prongNplanes,"prongNplanes[Nprongs]/I");
691 
692  _otree->Branch("prong1cellNcell",&prong1cellNcell,"prong1cellNcell/I");
693  _otree->Branch("prong1Fmip",&prong1Fmip,"prong1Fmip/F");
694  _otree->Branch("prong1cellX",&prong1cellX,"prong1cellX[prong1cellNcell]/F");
695  _otree->Branch("prong1cellY",&prong1cellY,"prong1cellY[prong1cellNcell]/F");
696  _otree->Branch("prong1cellZ",&prong1cellZ,"prong1cellZ[prong1cellNcell]/F");
697  _otree->Branch("prong1cellE",&prong1cellE,"prong1cellE[prong1cellNcell]/F");
698  _otree->Branch("prong1cellADC",&prong1cellADC,"prong1cellADC[prong1cellNcell]/I");
699  _otree->Branch("prong1cellPlane",&prong1cellPlane,"prong1cellPlane[prong1cellNcell]/I");
700  _otree->Branch("prong1cellCell",&prong1cellCell,"prong1cellCell[prong1cellNcell]/I");
701  _otree->Branch("prong1cellW",&prong1cellW,"prong1cellW[prong1cellNcell]/F");
702  _otree->Branch("prong1cellPECorr",&prong1cellPECorr,"prong1cellPECorr[prong1cellNcell]/F");
703 
704  _otree->Branch("allprongcellNcell",&allprongcellNcell,"allprongcellNcell/I");
705  _otree->Branch("allprongIndex",&allprongIndex,"allprongIndex[allprongcellNcell]/I");
706  _otree->Branch("allprongcellX",&allprongcellX,"allprongcellX[allprongcellNcell]/F");
707  _otree->Branch("allprongcellY",&allprongcellY,"allprongcellY[allprongcellNcell]/F");
708  _otree->Branch("allprongcellZ",&allprongcellZ,"allprongcellZ[allprongcellNcell]/F");
709  _otree->Branch("allprongcellE",&allprongcellE,"allprongcellE[allprongcellNcell]/F");
710  _otree->Branch("allprongcellADC",&allprongcellADC,"allprongcellADC[allprongcellNcell]/I");
711  _otree->Branch("allprongcellPlane",&allprongcellPlane,"allprongcellPlane[allprongcellNcell]/I");
712  _otree->Branch("allprongcellCell",&allprongcellCell,"allprongcellCell[allprongcellNcell]/I");
713  _otree->Branch("allprongcellW",&allprongcellW,"allprongcellW[allprongcellNcell]/F");
714  _otree->Branch("allprongcellPECorr",&allprongcellPECorr,"allprongcellPECorr[allprongcellNcell]/F");
715 
716  _otree->Branch("allprong2DcellNcell",&allprong2DcellNcell,"allprong2DcellNcell/I");
717  _otree->Branch("allprong2DIndex",&allprong2DIndex,"allprong2DIndex[allprong2DcellNcell]/I");
718  _otree->Branch("allprong2DcellX",&allprong2DcellX,"allprong2DcellX[allprong2DcellNcell]/F");
719  _otree->Branch("allprong2DcellY",&allprong2DcellY,"allprong2DcellY[allprong2DcellNcell]/F");
720  _otree->Branch("allprong2DcellZ",&allprong2DcellZ,"allprong2DcellZ[allprong2DcellNcell]/F");
721  _otree->Branch("allprong2DcellE",&allprong2DcellE,"allprong2DcellE[allprong2DcellNcell]/F");
722  _otree->Branch("allprong2DcellADC",&allprong2DcellADC,"allprong2DcellADC[allprong2DcellNcell]/I");
723  _otree->Branch("allprong2DcellPlane",&allprong2DcellPlane,"allprong2DcellPlane[allprong2DcellNcell]/I");
724  _otree->Branch("allprong2DcellCell",&allprong2DcellCell,"allprong2DcellCell[allprong2DcellNcell]/I");
725  _otree->Branch("allprong2DcellPECorr",&allprong2DcellPECorr,"allprong2DcellPECorr[allprong2DcellNcell]/F");
726 
727  _otree->Branch("shwcellNcell",&shwcellNcell,"shwcellNcell/I");
728  _otree->Branch("shwIndex",&shwIndex,"shwIndex[shwcellNcell]/I");
729  _otree->Branch("shwcellX",&shwcellX,"shwcellX[shwcellNcell]/F");
730  _otree->Branch("shwcellY",&shwcellY,"shwcellY[shwcellNcell]/F");
731  _otree->Branch("shwcellZ",&shwcellZ,"shwcellZ[shwcellNcell]/F");
732  _otree->Branch("shwcellE",&shwcellE,"shwcellE[shwcellNcell]/F");
733  _otree->Branch("shwcellADC",&shwcellADC,"shwcellADC[shwcellNcell]/I");
734  _otree->Branch("shwcellPlane",&shwcellPlane,"shwcellPlane[shwcellNcell]/I");
735  _otree->Branch("shwcellCell",&shwcellCell,"shwcellCell[shwcellNcell]/I");
736  _otree->Branch("shwcellW",&shwcellW,"shwcellW[shwcellNcell]/F");
737  _otree->Branch("shwcellPECorr",&shwcellPECorr,"shwcellPECorr[shwcellNcell]/F");
738 /*
739  _otree->Branch("q2prongcellNcell",&q2prongcellNcell,"q2prongcellNcell/I");
740  _otree->Branch("q2prongcellX",&q2prongcellX,"q2prongcellX[q2prongcellNcell]/F");
741  _otree->Branch("q2prongcellY",&q2prongcellY,"q2prongcellY[q2prongcellNcell]/F");
742  _otree->Branch("q2prongcellZ",&q2prongcellZ,"q2prongcellZ[q2prongcellNcell]/F");
743  _otree->Branch("q2prongcellE",&q2prongcellE,"q2prongcellE[q2prongcellNcell]/F");
744  _otree->Branch("q2prongcellADC",&q2prongcellADC,"q2prongcellADC[q2prongcellNcell]/I");
745  _otree->Branch("q2prongcellPlane",&q2prongcellPlane,"q2prongcellPlane[q2prongcellNcell]/I");
746  _otree->Branch("q2prongcellCell",&q2prongcellCell,"q2prongcellCell[q2prongcellNcell]/I");
747  _otree->Branch("q2prongcellW",&q2prongcellW,"q2prongcellW[q2prongcellNcell]/F");
748  _otree->Branch("q2prongcellPECorr",&q2prongcellPECorr,"q2prongcellPECorr[q2prongcellNcell]/F");
749 
750  _otree->Branch("q1prongcellNcell",&q1prongcellNcell,"q1prongcellNcell/I");
751  _otree->Branch("q1prongcellX",&q1prongcellX,"q1prongcellX[q1prongcellNcell]/F");
752  _otree->Branch("q1prongcellY",&q1prongcellY,"q1prongcellY[q1prongcellNcell]/F");
753  _otree->Branch("q1prongcellZ",&q1prongcellZ,"q1prongcellZ[q1prongcellNcell]/F");
754  _otree->Branch("q1prongcellE",&q1prongcellE,"q1prongcellE[q1prongcellNcell]/F");
755  _otree->Branch("q1prongcellADC",&q1prongcellADC,"q1prongcellADC[q1prongcellNcell]/I");
756  _otree->Branch("q1prongcellPlane",&q1prongcellPlane,"q1prongcellPlane[q1prongcellNcell]/I");
757  _otree->Branch("q1prongcellCell",&q1prongcellCell,"q1prongcellCell[q1prongcellNcell]/I");
758  _otree->Branch("q1prongcellW",&q1prongcellW,"q1prongcellW[q1prongcellNcell]/F");
759  _otree->Branch("q1prongcellPECorr",&q1prongcellPECorr,"q1prongcellPECorr[q1prongcellNcell]/F");
760 
761  _otree->Branch("prong2cellNcell",&prong2cellNcell,"prong2cellNcell/I");
762  _otree->Branch("prong2cellX",&prong2cellX,"prong2cellX[prong2cellNcell]/F");
763  _otree->Branch("prong2cellY",&prong2cellY,"prong2cellY[prong2cellNcell]/F");
764  _otree->Branch("prong2cellZ",&prong2cellZ,"prong2cellZ[prong2cellNcell]/F");
765  _otree->Branch("prong2cellE",&prong2cellE,"prong2cellE[prong2cellNcell]/F");
766  _otree->Branch("prong2cellADC",&prong2cellADC,"prong2cellADC[prong2cellNcell]/I");
767  _otree->Branch("prong2cellPlane",&prong2cellPlane,"prong2cellPlane[prong2cellNcell]/I");
768  _otree->Branch("prong2cellCell",&prong2cellCell,"prong2cellCell[prong2cellNcell]/I");
769  _otree->Branch("prong2cellW",&prong2cellW,"prong2cellW[prong2cellNcell]/F");
770  _otree->Branch("prong2cellPECorr",&prong2cellPECorr,"prong2cellPECorr[prong2cellNcell]/F");
771 */
772  _otree->Branch("Nprongs2D",&Nprongs2D,"Nprongs2D/I");
773  _otree->Branch("prong2DInXView",&prong2DInXView,"prong2DInXView[Nprongs2D]/I");
774  _otree->Branch("prongEnergy2D",&prongEnergy2D,"prongEnergy2D[Nprongs2D]/F");
775  _otree->Branch("prongNCells2D",&prongNCells2D,"prongNCells2D[Nprongs2D]/I");
776  _otree->Branch("prongLength2D",&prongLength2D,"prongLength2D[Nprongs2D]/F");
777  _otree->Branch("prongStopX2D",&prongStopX2D,"prongStopX2D[Nprongs2D]/F");
778  _otree->Branch("prongStopY2D",&prongStopY2D,"prongStopY2D[Nprongs2D]/F");
779  _otree->Branch("prongStopZ2D",&prongStopZ2D,"prongStopZ2D[Nprongs2D]/F");
780  _otree->Branch("prongPhi2D",&prongPhi2D,"prongPhi2D[Nprongs2D]/F");
781  _otree->Branch("prongTheta2D",&prongTheta2D,"prongTheta2D[Nprongs2D]/F");
782  _otree->Branch("prongStartX2D",&prongStartX2D,"prongStartX2D[Nprongs2D]/F");
783  _otree->Branch("prongStartY2D",&prongStartY2D,"prongStartY2D[Nprongs2D]/F");
784  _otree->Branch("prongStartZ2D",&prongStartZ2D,"prongStartZ2D[Nprongs2D]/F");
785  _otree->Branch("Ebalance2D",&Ebalance2D,"Ebalance2D/F");
786 
787  _otree->Branch("rvp",&rvp,"rvp/F");
788  _otree->Branch("lid",&lid,"lid/F");
789  _otree->Branch("lide",&lide,"lide/F");
790  _otree->Branch("lem",&lem,"lem/F");
791  _otree->Branch("Nshowers",&Nshowers,"Nshowers/I");
792  _otree->Branch("shwLID",&shwLID,"shwLID[Nshowers]/F");
793  _otree->Branch("shw3D",&shw3D,"shw3D[Nshowers]/I");
794  _otree->Branch("shwNCells",&shwNCells,"shwNCells[Nshowers]/I");
795  _otree->Branch("shwNCellsXview",&shwNCellsXview,"shwNcellsXview[Nshowers]/I");
796  _otree->Branch("shwNCellsYview",&shwNCellsYview,"shwNcellsYview[Nshowers]/I");
797  _otree->Branch("shwLength",&shwLength,"shwLength[Nshowers]/F");
798  _otree->Branch("shwStartX",&shwStartX,"shwStartX[Nshowers]/F");
799  _otree->Branch("shwStartY",&shwStartY,"shwStartY[Nshowers]/F");
800  _otree->Branch("shwStartZ",&shwStartZ,"shwStartZ[Nshowers]/F");
801  _otree->Branch("shwStopX",&shwStopX,"shwStopX[Nshowers]/F");
802  _otree->Branch("shwStopY",&shwStopY,"shwStopY[Nshowers]/F");
803  _otree->Branch("shwStopZ",&shwStopZ,"shwStopZ[Nshowers]/F");
804  _otree->Branch("shwTheta",&shwTheta,"shwTheta[Nshowers]/F");
805  _otree->Branch("shwEcells",&shwEcells,"shwEcells[Nshowers]/F");
806  _otree->Branch("shwEnu",&shwEnu,"shwEnu[Nshowers]/F");
807 
808  _btree=new TTree("beamInfo","beam info");
809  _btree->Branch("runnum",&runnum,"runnum/I");
810  _btree->Branch("subrunnum",&subrunnum,"subrunnum/I");
811  _btree->Branch("evtnum",&evtnum,"evtnum/I");
812  _btree->Branch("Nslices",&Nslices,"Nslices/I");
813  _btree->Branch("spillPOT",&spillPOT,"spillPOT/D");
814  _btree->Branch("goodSpill",&goodSpill,"goodSpill/I");
815  _btree->Branch("horncurrent",&horncurrent,"horncurrent/D");
816  _btree->Branch("posX",&posX,"posX/D");
817  _btree->Branch("posY",&posY,"posY/D");
818  _btree->Branch("widthX",&widthX,"widthX/D");
819  _btree->Branch("widthY",&widthY,"widthY/D");
820  _btree->Branch("deltaspilltimensec",&deltaspilltimensec,"deltaspilltimensec/I");
821 
822  _btree->Branch("Ndeaddcm",&Ndeaddcm,"Ndeaddcm/I");
823  _btree->Branch("deaddcm_location",&deaddcm_location,"deaddcm_location[Ndeaddcm]/S");
824  _btree->Branch("Ndb1apd_24cells",&Ndb1apd_24cells,"Ndb1apd_24cells/I");
825  _btree->Branch("Ndb2apd_24cells",&Ndb2apd_24cells,"Ndb2apd_24cells/I");
826  _btree->Branch("Ndb3apd_24cells",&Ndb3apd_24cells,"Ndb3apd_24cells/I");
827  _btree->Branch("Nmcapd_24cells",&Nmcapd_24cells,"Nmcapd_24cells/I");
828  _btree->Branch("Ndb1apd_28cells",&Ndb1apd_28cells,"Ndb1apd_28cells/I");
829  _btree->Branch("Ndb2apd_28cells",&Ndb2apd_28cells,"Ndb2apd_28cells/I");
830  _btree->Branch("Ndb3apd_28cells",&Ndb3apd_28cells,"Ndb3apd_28cells/I");
831  _btree->Branch("Nmcapd_28cells",&Nmcapd_28cells,"Nmcapd_28cells/I");
832  _btree->Branch("Nintimehits",&Nintimehits,"Nintimehits/I");
833  _btree->Branch("Nouttimehits",&Nouttimehits,"Nouttimehits/I");
834  _btree->Branch("Fouttimehits_noisyapd",&Fouttimehits_noisyapd,"Fouttimehits_noisyapd/F");
835 
836  _btree->Branch("evt_year",&evt_year,"evt_year/I");
837  _btree->Branch("evt_month",&evt_month,"evt_month/I");
838  _btree->Branch("evt_day",&evt_day,"evt_day/I");
839  _btree->Branch("evt_hour",&evt_hour,"evt_hour/I");
840  _btree->Branch("evt_min",&evt_min,"evt_min/I");
841  _btree->Branch("evt_sec",&evt_sec,"evt_sec/I");
842 
843  }
844 
845  //......................................................................
846  void NCAna::endSubRun(const art::SubRun& sr)
847  {
848  art::Handle<sumdata::POTSum> pot;
849  if( isMC ) sr.getByLabel(fMCFluxLabel, pot);
850  else sr.getByLabel(fDataSpillLabel, pot);
851 
852  if( !pot.failedToGet() ){
853  fTotalPOT += pot->totgoodpot;
854  fTotalSpill += pot->goodspills;
855 
856  std::cout << "POT in subrun " << sr.subRun()
857  << " = " << pot->totgoodpot <<" with "<<fTotalSpill<<" spills"<< std::endl;
858  }
859  }
860 
861  //......................................................................
862  const rb::RecoHit NCAna::MakeRecoHit(art::Ptr<rb::CellHit> const& chit,
863  double const& w)
864  {
865  art::ServiceHandle<calib::Calibrator> cal;
866  const rb::RecoHit rhit(cal->MakeRecoHit(*chit, w));
867 
868  return rhit;
869  }
870 
871  //......................................................................
872  void NCAna::analyze(const art::Event& evt)
873  {
874  sl_runnum=-1;
875  sl_subrunnum=-1;
876  sl_evtnum=-1;
877  sl_Nslices=-1;
878  sl_goodSpill=-1;
879  sl_spillPOT=-1.;
880  sl_horncurrent=-1.;
881  sl_posX=-1.;
882  sl_posY=-1.;
883  sl_widthX=-1.;
884  sl_widthY=-1.;
885  sl_deltaspilltimensec=-1;
886 
887  sl_Ndeaddcm=-1;
888  sl_Ndb1apd_24cells=-1;
889  sl_Ndb2apd_24cells=-1;
890  sl_Ndb3apd_24cells=-1;
891  sl_Nmcapd_24cells=-1;
892  sl_Ndb1apd_28cells=-1;
893  sl_Ndb2apd_28cells=-1;
894  sl_Ndb3apd_28cells=-1;
895  sl_Nmcapd_28cells=-1;
896  sl_Nintimehits=-1;
897  sl_Nouttimehits=-1;
898  for( int i=0; i<14; ++i )
899  sl_deaddcm_location[i]="NULL";
900  sl_Fouttimehits_noisyapd=9999.;
901 
902  runnum=-1;
903  subrunnum=-1;
904  evtnum=-1;
905  Nslices=-1;
906  spillPOT=0.;
907  goodSpill=-1;
908  horncurrent=-1.;
909  posX=-1.;
910  posY=-1.;
911  widthX=-1.;
912  widthY=-1.;
913  deltaspilltimensec=-1;
914 
915  Ndeaddcm=-1;
916  Ndb1apd_24cells=-1;
917  Ndb2apd_24cells=-1;
918  Ndb3apd_24cells=-1;
919  Nmcapd_24cells=-1;
920  Ndb1apd_28cells=-1;
921  Ndb2apd_28cells=-1;
922  Ndb3apd_28cells=-1;
923  Nmcapd_28cells=-1;
924  Nintimehits=-1;
925  Nouttimehits=-1;
926  for( int i=0; i<14; ++i )
927  deaddcm_location[i]="NULL";
928  Fouttimehits_noisyapd=9999.;
929 
930  evt_year=0;
931  evt_month=0;
932  evt_day=-1;
933  evt_hour=-1;
934  evt_min=-1;
935  evt_sec=-1;
936 
937  rvp=-99.;
938  lid=-99.;
939  lide=-99.;
940  lem=-99.;
941 
942  Nshowers=-1;
943  for( int i=0; i<20; ++i ){
944  shwLID[i] = -99.;
945  shw3D[i]=-1;
946  shwNCells[i]=-1;
947  shwNCellsXview[i]=-1;
948  shwNCellsYview[i]=-1;
949  shwLength[i]=-1.;
950  shwStartX[i]=-999.;
951  shwStartY[i]=-999.;
952  shwStartZ[i]=-999.;
953  shwStopX[i]=-999.;
954  shwStopY[i]=-999.;
955  shwStopZ[i]=-999.;
956  shwTheta[i]=-1.;
957  shwEcells[i]=-1.;
958  shwEnu[i]=-1.;
959  }
960 
961  Vx=-999.;
962  Vy=-999.;
963  Vz=-999.;
964  recoVx=-999.;
965  recoVy=-999.;
966  recoVz=-999.;
967 
968  truePt=-1.;
969  trueE=-1.;
970  ccnc=-1;
971  PDG=-1;
972  origPDG=-1;
973  mode=-1;
974  isQE=false;
975  intType=-1;
976  hitNucl=-1;
977  HadX=-999.;
978  HadY=-999.;
979  HadW=-999.;
980  qsqr=-999.;
981  OscWt=-999.;
982 
983  lepPDG=-1;
984  lepPx=-1.;
985  lepPy=-1.;
986  lepPz=-1.;
987  mesonPDG=-1;
988  mesonPx=-1.;
989  mesonPy=-1.;
990  mesonPz=-1.;
991 
992  Nweights=0;
993  for( int i=0; i<68; ++i ){
994  plus1sigma[i]=1.;
995  plus2sigma[i]=1.;
996  minus1sigma[i]=1.;
997  minus2sigma[i]=1.;
998  }
999 
1000  HasP2GG=-1;
1001  Epi0=-1.;
1002  gam1E=-1.;
1003  gam1Px=-1.;
1004  gam1Py=-1.;
1005  gam1Pz=-1.;
1006  gam2E=-1.;
1007  gam2Px=-1.;
1008  gam2Py=-1.;
1009  gam2Pz=-1.;
1010 
1011  Npi0=0;
1012  for( int i=0; i<20; ++i ){
1013  pi0Px[i]=-1.;
1014  pi0Py[i]=-1.;
1015  pi0Pz[i]=-1.;
1016  pi0E[i]=-1.;
1017  pi0Vx[i]=-1.;
1018  pi0Vy[i]=-1.;
1019  pi0Vz[i]=-1.;
1020  }
1021 
1022  Ncells=0;
1023  Nplanes=0;
1024  recoE=-1.;
1025  MinPlane=-1;
1026  sliceIndex=-1;
1027  slice_meantns=-1.;
1028  slice_mintns=-1.;
1029  slice_maxtns=-1.;
1030  for( int i=0; i<192; ++i ){
1031  PlaneEnergy[i] = 0.;
1032  PlaneNcells[i] = 0;
1033  }
1034 
1035  KtrackLength=0.;
1036  KtrackNcells=0;
1037  KtrackEnergy=0.;
1038  KtrackEnergy_ratio=0.;
1039  KtrackNcells_ratio=0.;
1040  KtrackNplanes=0;
1041  KtrackTheta=-1.;
1042  KtrackPhi=-1.;
1043  KtrackStartX=-999.;
1044  KtrackStartY=-999.;
1045  KtrackStartZ=-999.;
1046  KtrackStopX=-999.;
1047  KtrackStopY=-999.;
1048  KtrackStopZ=-999.;
1049  NKtracks=0.;
1050 
1051  Nmip=-1.;
1052  Fmip=-1.;
1053 
1054  Nprongs=0;
1055  Ebalance=0.;
1056  Dphi=-1.;
1057 
1058  for( int i=0; i<20; ++i ){
1059  prongEnergy[i]=0.;
1060  prongPhi[i]=-1.;
1061  prongTheta[i]=-1.;
1062  prongStartX[i]=-999.;
1063  prongStartY[i]=-999.;
1064  prongStartZ[i]=-999.;
1065  prongStopX[i]=-999.;
1066  prongStopY[i]=-999.;
1067  prongStopZ[i]=-999.;
1068  prongLength[i]=-1.;
1069  prongEnergyXView[i]=0.;
1070  prongEnergyYView[i]=0.;
1071  prongNCells[i]=0;
1072  prongNCellsXView[i]=0;
1073  prongNCellsYView[i]=0;
1074  prongNplanes[i]=-1;
1075 
1076  prong2DInXView[i]=0;
1077  prongEnergy2D[i]=0.;
1078  prongNCells2D[i]=0;
1079  prongLength2D[i]=-1.;
1080  prongPhi2D[i]=-1.;
1081  prongTheta2D[i]=-1.;
1082  prongStartX2D[i]=-999.;
1083  prongStartY2D[i]=-999.;
1084  prongStartZ2D[i]=-999.;
1085  prongStopX2D[i]=-999.;
1086  prongStopY2D[i]=-999.;
1087  prongStopZ2D[i]=-999.;
1088  }
1089  Nprongs2D=0;
1090  Ebalance2D=0.;
1091 
1092  slcellNcell=0;
1093  prong1cellNcell=0;
1094  prong1Fmip=-1.;
1095  //prong2cellNcell=0;
1096  trackcellNcell=0;
1097  allprongcellNcell=0;
1098  allprong2DcellNcell=0;
1099  shwcellNcell=0;
1100  //q1prongcellNcell=0;
1101  //q2prongcellNcell=0;
1102  for( int i=0; i<500; ++i ){
1103  slcellX[i]=-999.;
1104  slcellY[i]=-999.;
1105  slcellZ[i]=-999.;
1106  slcellE[i]=-999.;
1107  slcellADC[i]=0;
1108  slcellPlane[i]=-1;
1109  slcellCell[i]=-1;
1110  slcellW[i]=-1.;
1111  slcellTNS[i]=-1.;
1112  slcellPECorr[i]=-1.;
1113 
1114  prong1cellX[i]=-999.;
1115  prong1cellY[i]=-999.;
1116  prong1cellZ[i]=-999.;
1117  prong1cellE[i]=-999.;
1118  prong1cellADC[i]=0;
1119  prong1cellPlane[i]=-1;
1120  prong1cellCell[i]=-1;
1121  prong1cellW[i]=-1.;
1122  prong1cellPECorr[i]=-1.;
1123 
1124  allprongIndex[i]=-1;
1125  allprongcellX[i]=-999.;
1126  allprongcellY[i]=-999.;
1127  allprongcellZ[i]=-999.;
1128  allprongcellE[i]=-999.;
1129  allprongcellADC[i]=0;
1130  allprongcellPlane[i]=-1;
1131  allprongcellCell[i]=-1;
1132  allprongcellW[i]=-1.;
1133  allprongcellPECorr[i]=-1.;
1134 
1135  allprong2DIndex[i]=-1;
1136  allprong2DcellX[i]=-999.;
1137  allprong2DcellY[i]=-999.;
1138  allprong2DcellZ[i]=-999.;
1139  allprong2DcellE[i]=-999.;
1140  allprong2DcellADC[i]=0;
1141  allprong2DcellPlane[i]=-1;
1142  allprong2DcellCell[i]=-1;
1143  allprong2DcellPECorr[i]=-1.;
1144 
1145  shwIndex[i]=-1;
1146  shwcellX[i]=-999.;
1147  shwcellY[i]=-999.;
1148  shwcellZ[i]=-999.;
1149  shwcellE[i]=-999.;
1150  shwcellADC[i]=0;
1151  shwcellPlane[i]=-1;
1152  shwcellCell[i]=-1;
1153  shwcellW[i]=-1.;
1154  shwcellPECorr[i]=-1.;
1155 
1156 /*
1157  q1prongcellX[i]=-999.;
1158  q1prongcellY[i]=-999.;
1159  q1prongcellZ[i]=-999.;
1160  q1prongcellE[i]=-999.;
1161  q1prongcellADC[i]=0;
1162  q1prongcellPlane[i]=-1;
1163  q1prongcellCell[i]=-1;
1164  q1prongcellW[i]=-1.;
1165  q1prongcellPECorr[i]=-1.;
1166 
1167  q2prongcellX[i]=-999.;
1168  q2prongcellY[i]=-999.;
1169  q2prongcellZ[i]=-999.;
1170  q2prongcellE[i]=-999.;
1171  q2prongcellADC[i]=0;
1172  q2prongcellPlane[i]=-1;
1173  q2prongcellCell[i]=-1;
1174  q2prongcellW[i]=-1.;
1175  q2prongcellPECorr[i]=-1.;
1176 
1177  prong2cellX[i]=-999.;
1178  prong2cellY[i]=-999.;
1179  prong2cellZ[i]=-999.;
1180  prong2cellE[i]=-999.;
1181  prong2cellADC[i]=0;
1182  prong2cellPlane[i]=-1;
1183  prong2cellCell[i]=-1;
1184  prong2cellW[i]=-1.;
1185  prong2cellPECorr[i]=-1.;
1186 */
1187  trackcellX[i]=-999.;
1188  trackcellY[i]=-999.;
1189  trackcellZ[i]=-999.;
1190  trackcellE[i]=-999.;
1191  trackcellADC[i]=0;
1192  trackcellPlane[i]=-1;
1193  trackcellCell[i]=-1;
1194  trackcellW[i]=-1.;
1195  trackcellPECorr[i]=-1.;
1196 
1197  }
1198 
1199  std::string dcmName[14]={"DB1-DCM1","DB1-DCM2","DB1-DCM3","DB1-DCM4",
1200  "DB2-DCM1","DB2-DCM2","DB2-DCM3","DB2-DCM4",
1201  "DB3-DCM1","DB3-DCM2","DB3-DCM3","DB3-DCM4",
1202  "MC-DCM1","MC-DCM2"};
1203 
1204  art::ServiceHandle<calib::Calibrator> cal;
1205  art::ServiceHandle<geo::Geometry> geom;
1206 
1207  //get data spill info
1208  art::Handle<sumdata::SpillData> spillPot;
1209  if( isMC ) evt.getByLabel(fMCFluxLabel, spillPot);
1210  else evt.getByLabel(fDataSpillLabel, spillPot);
1211 
1212  if( spillPot.failedToGet() ) return;
1213 
1214  //std::cout<<"POT info "<<spillPot->spillpot/1e+13<<std::endl;
1215  spillPOT = spillPot->spillpot;
1216  if( isMC ) spillPOT=spillPOT/1.0e+12;
1217  if( !isMC && !isROCK ){
1218  if( spillPot->goodbeam ) goodSpill=1;
1219  horncurrent = spillPot->hornI;
1220  posX = spillPot->posx;
1221  posY = spillPot->posy;
1222  widthX = spillPot->widthx;
1223  widthY = spillPot->widthy;
1224  deltaspilltimensec = spillPot->deltaspilltimensec;
1225  }
1226 
1227  //event info
1228  runnum = evt.run();
1229  subrunnum = evt.subRun();
1230  evtnum = evt.id().event();
1231 
1232  //get all hits info
1233  art::Handle<std::vector<rb::CellHit> > chits;
1234  evt.getByLabel(fCellHitLabel, chits);
1235 
1236  std::vector<int> sq_plane;
1237  sq_plane.clear();
1238  std::vector<int> sq_xview;
1239  sq_xview.clear();
1240  std::vector<float> sq_tns;
1241  sq_tns.clear();
1242  std::vector<int> sq_cell;
1243  sq_cell.clear();
1244 
1245  int sq_ncells=chits->size();
1246  for(int i = 0; i < sq_ncells; ++i){
1247  art::Ptr< rb::CellHit> hit(chits, i);
1248  double hittime=hit->TNS()/1000.;
1249  sq_plane.push_back(hit->Plane());
1250  sq_tns.push_back(hittime);
1251  sq_cell.push_back(hit->Cell());
1252  if( (hit->View())==geo::kX ) sq_xview.push_back(1);
1253  else sq_xview.push_back(0);
1254  }//loop all hits
1255 
1256  int ncells_intime=0;
1257  int ncells_outtime=0;
1258  //find dead dcms
1259  int Ndbdcm_intime[14]={0};
1260  int Ndbdcm_outtime[14]={0};
1261  //find dead and noisy apds
1262  int Ndb1apd[192]={0};
1263  int Ndb2apd[192]={0};
1264  int Ndb3apd[192]={0};
1265  int Nmcapd[55]={0};
1266 
1267  for( int ic=0; ic<sq_ncells; ++ic ){
1268  if( sq_tns[ic]>218. && sq_tns[ic]<228. ){
1269  ncells_intime++;
1270  }//in beam window
1271  else {
1272  ncells_outtime++;
1273 
1274  if( sq_plane[ic]<64 ){
1275  for( int ip=0; ip<64; ++ip ){
1276  if( ip == sq_plane[ic] ){
1277  if( sq_cell[ic]<32 ) Ndb1apd[0+ip*3]++;
1278  else if( sq_cell[ic]<64 ) Ndb1apd[1+ip*3]++;
1279  else Ndb1apd[2+ip*3]++;
1280  }//for each certain plane
1281  }//loop 64 planes
1282  }//db1
1283  else if( sq_plane[ic]<128 ){
1284  for( int ip=0; ip<64; ++ip ){
1285  if( (ip+64) == sq_plane[ic] ){
1286  if( sq_cell[ic]<32 ) Ndb2apd[0+ip*3]++;
1287  else if( sq_cell[ic]<64 ) Ndb2apd[1+ip*3]++;
1288  else Ndb2apd[2+ip*3]++;
1289  }//for each certain plane
1290  }//loop 64 planes
1291  }//db2
1292  else if( sq_plane[ic]<192 ){
1293  for( int ip=0; ip<64; ++ip ){
1294  if( (ip+128) == sq_plane[ic] ){
1295  if( sq_cell[ic]<32 ) Ndb3apd[0+ip*3]++;
1296  else if( sq_cell[ic]<64 ) Ndb3apd[1+ip*3]++;
1297  else Ndb3apd[2+ip*3]++;
1298  }//for each certain plane
1299  }//loop 64 planes
1300  }//db3
1301  else {
1302  for( int ip=0; ip<22; ++ip ){
1303  if( (ip+192) == sq_plane[ic] ){
1304  if( (ip+192)%2 ==0 ){
1305  if( sq_cell[ic]<32 ) Nmcapd[0+5*ip/2]++;
1306  else Nmcapd[1+5*ip/2]++;
1307  }//YZ view - 22 modules
1308  else{
1309  if( sq_cell[ic]<32 ) Nmcapd[2+5*(ip-1)/2]++;
1310  else if( sq_cell[ic]<64 ) Nmcapd[3+5*(ip-1)/2]++;
1311  else Nmcapd[4+5*(ip-1)/2]++;
1312  }//XZ view - 33 modules
1313  }//for each certain plane
1314  }//loop 22 planes
1315  }//MC
1316  }//out-beam-window
1317  if( sq_plane[ic]<64 ){
1318  if( sq_xview[ic]==1 ){
1319  if( sq_cell[ic]>63 ){
1320  if( sq_tns[ic]>218. && sq_tns[ic]<228. ) Ndbdcm_intime[0]++;
1321  else Ndbdcm_outtime[0]++;
1322  }//dcm 1
1323  else{
1324  if( sq_tns[ic]>218. && sq_tns[ic]<228. ) Ndbdcm_intime[1]++;
1325  else Ndbdcm_outtime[1]++;
1326  }//dcm 2
1327  }//XZ view
1328  else {
1329  if( sq_cell[ic]>31 ){
1330  if( sq_tns[ic]>218. && sq_tns[ic]<228. ) Ndbdcm_intime[2]++;
1331  else Ndbdcm_outtime[2]++;
1332  }//dcm 3
1333  else {
1334  if( sq_tns[ic]>218. && sq_tns[ic]<228. ) Ndbdcm_intime[3]++;
1335  else Ndbdcm_outtime[3]++;
1336  }//dcm 4
1337  }//YZ view
1338  }//DB1
1339  else if( sq_plane[ic]<128 ){
1340  if( sq_xview[ic]==1 ){
1341  if( sq_cell[ic]>63 ){
1342  if( sq_tns[ic]>218. && sq_tns[ic]<228. ) Ndbdcm_intime[4]++;
1343  else Ndbdcm_outtime[4]++;
1344  }//dcm 1
1345  else{
1346  if( sq_tns[ic]>218. && sq_tns[ic]<228. ) Ndbdcm_intime[5]++;
1347  else Ndbdcm_outtime[5]++;
1348  }//dcm 2
1349  }//XZ view
1350  else {
1351  if( sq_cell[ic]>31 ){
1352  if( sq_tns[ic]>218. && sq_tns[ic]<228. ) Ndbdcm_intime[6]++;
1353  else Ndbdcm_outtime[6]++;
1354  }//dcm 3
1355  else {
1356  if( sq_tns[ic]>218. && sq_tns[ic]<228. ) Ndbdcm_intime[7]++;
1357  else Ndbdcm_outtime[7]++;
1358  }//dcm 4
1359  }//YZ view
1360  }//DB2
1361  else if( sq_plane[ic]<192 ){
1362  if( sq_xview[ic]==1 ){
1363  if( sq_cell[ic]>63 ){
1364  if( sq_tns[ic]>218. && sq_tns[ic]<228. ) Ndbdcm_intime[8]++;
1365  else Ndbdcm_outtime[8]++;
1366  }//dcm 1
1367  else{
1368  if( sq_tns[ic]>218. && sq_tns[ic]<228. ) Ndbdcm_intime[9]++;
1369  else Ndbdcm_outtime[9]++;
1370  }//dcm 2
1371  }//XZ view
1372  else {
1373  if( sq_cell[ic]>31 ){
1374  if( sq_tns[ic]>218. && sq_tns[ic]<228. ) Ndbdcm_intime[10]++;
1375  else Ndbdcm_outtime[10]++;
1376  }//dcm 3
1377  else {
1378  if( sq_tns[ic]>218. && sq_tns[ic]<228. ) Ndbdcm_intime[11]++;
1379  else Ndbdcm_outtime[11]++;
1380  }//dcm 4
1381  }//YZ view
1382  }//DB3
1383  else {
1384  if( sq_xview[ic]==1 ){
1385  if( sq_tns[ic]>218. && sq_tns[ic]<228. ) Ndbdcm_intime[12]++;
1386  else Ndbdcm_outtime[12]++;
1387  }//XZ view
1388  else {
1389  if( sq_tns[ic]>218. && sq_tns[ic]<228. ) Ndbdcm_intime[13]++;
1390  else Ndbdcm_outtime[13]++;
1391  }//YZ view
1392  }//MC
1393  }//loop all cells
1394 
1395  int ndeaddcms=0;
1396  std::vector<std::string> dcmlocation;
1397  dcmlocation.clear();
1398  for( int id=0; id<14; ++id ){
1399  if( Ndbdcm_outtime[id]==0 && Ndbdcm_intime[id]==0 ){
1400  ndeaddcms++;
1401  dcmlocation.push_back(dcmName[id]);
1402  }
1403  }//loop all 14 dcms associated cells
1404 
1405  int Napd_24cells[4]={0};
1406  int Napd_28cells[4]={0};
1407  for( int ia=0; ia<192; ++ia ){
1408  //DB1
1409  if( Ndb1apd[ia]>=24 ) Napd_24cells[0]+=1.;
1410  if( Ndb1apd[ia]>=28 ) Napd_28cells[0]+=1.;
1411  //DB2
1412  if( Ndb2apd[ia]>=24 ) Napd_24cells[1]+=1.;
1413  if( Ndb2apd[ia]>=28 ) Napd_28cells[1]+=1.;
1414  //DB3
1415  if( Ndb3apd[ia]>=24 ) Napd_24cells[2]+=1.;
1416  if( Ndb3apd[ia]>=28 ) Napd_28cells[2]+=1.;
1417  //MC
1418  if( ia<55 ){
1419  if( Nmcapd[ia]>=24 ) Napd_24cells[3]+=1.;
1420  if( Nmcapd[ia]>=28 ) Napd_28cells[3]+=1.;
1421  }
1422  }//end loop of db planes
1423 
1424  int Napd_24cells_tot=0;
1425  int Napd_28cells_tot=0;
1426  for( int id=0; id<4; ++id ){
1427  Napd_24cells_tot += Napd_24cells[id];
1428  Napd_28cells_tot += Napd_28cells[id];
1429  }
1430 
1431  float ratio_outtimehits_deadapds=9999.;
1432  if( Napd_28cells_tot>0. ) ratio_outtimehits_deadapds=ncells_outtime*1.0/(Napd_28cells_tot*1.);
1433 
1434  Ndeaddcm=ndeaddcms;
1435  if( ndeaddcms>0 ){
1436  for( int id=0; id<ndeaddcms; ++id )
1437  deaddcm_location[id]=dcmlocation[id];
1438  }
1439  Nouttimehits=ncells_outtime;
1440  Nintimehits=ncells_intime;
1441  Ndb1apd_24cells=Napd_24cells[0];
1442  Ndb2apd_24cells=Napd_24cells[1];
1443  Ndb3apd_24cells=Napd_24cells[2];
1444  Nmcapd_24cells=Napd_24cells[3];
1445  Ndb1apd_28cells=Napd_28cells[0];
1446  Ndb2apd_28cells=Napd_28cells[1];
1447  Ndb3apd_28cells=Napd_28cells[2];
1448  Nmcapd_28cells=Napd_28cells[3];
1449  Fouttimehits_noisyapd=ratio_outtimehits_deadapds;
1450 
1451  art::Timestamp ts = evt.time();
1452  const time_t timeSec = ts.timeHigh();
1453  struct tm* timeStruct = localtime(&timeSec);
1454  evt_year=timeStruct->tm_year+1900;
1455  evt_month=timeStruct->tm_mon+1;
1456  evt_day=timeStruct->tm_mday;
1457  evt_hour=timeStruct->tm_hour+1;
1458  evt_min=timeStruct->tm_min;
1459  evt_sec=timeStruct->tm_sec;
1460 
1461 
1462  //get reco slicer info
1463  art::Handle< std::vector<rb::Cluster> > slices;
1464  evt.getByLabel(fSliceLabel, slices);
1465 
1466  //get tracks, verticies and showers associated with the slices
1467  //art::FindManyP<rb::Track> fmDTrack(slices, evt, fDTrackLabel);
1468  art::FindManyP<rb::Track> fmKTrack(slices, evt, fKTrackLabel);
1469  art::FindManyP<rb::Vertex> fmVertex(slices, evt, fVertexLabel);
1470  //art::FindManyP<jmshower::JMShower> fmShower(slices, evt, fShowerLabel);
1471 
1472  //std::cout<<evt.run()<<" "<<evt.subRun()<<" "<<evt.id().event()<<std::endl;
1473  //std::cout<<"# of slices: "<<slices->size()<<std::endl;
1474 
1475 
1476  //art::FindMany<rb::HoughResult> fmHough(slices, evt, "multihough");
1477 
1478  int nslice=0;
1479  for(unsigned int sliceIdx = 0; sliceIdx < slices->size(); ++sliceIdx){
1480 
1481  const rb::Cluster& slice = (*slices)[sliceIdx];
1482  if(slice.IsNoise()) continue;
1483  ++nslice;
1484 
1485  /*
1486  const rb::HoughResult* hrx = 0;
1487  const rb::HoughResult* hry = 0;
1488  int nX=0;
1489  int nY=0;
1490  std::vector<const rb::HoughResult*> hr = fmHough.at(sliceIdx);
1491  for (int j=0; j<hr.size(); ++j) {
1492  if (hr[j]->fView==geo::kX) {
1493  hrx = hr[j];
1494  ++nX;
1495  }
1496  if (hr[j]->fView==geo::kY) {
1497  hry = hr[j];
1498  ++nY;
1499  }
1500  }
1501 
1502  std::cout<<"# of hough lines "<<hr.size()<<std::endl;
1503  std::cout<<"in XZ view: "<<nX<<std::endl;
1504  std::cout<<"in YZ view: "<<nY<<std::endl;
1505  if (hrx==0 || hry==0) continue;
1506  std::cout<<"XZ peak: "<<hrx->fPeak.size()<<std::endl;
1507  std::cout<<"YZ peak: "<<hry->fPeak.size()<<std::endl;
1508  */
1509 
1510  }//loop slices
1511 
1512  Nslices=nslice;
1513  //std::cout<<"# of slice: "<<nslice<<std::endl;
1514 
1515  _btree->Fill();
1516 
1517  for(unsigned int sliceIdx = 0; sliceIdx < slices->size(); ++sliceIdx){
1518 
1519  const rb::Cluster& slice = (*slices)[sliceIdx];
1520  if(slice.IsNoise()) continue;
1521 
1522  //std::cout<<"*** # of slices: "<<slices->size()<<std::endl;
1523  //reco veretx
1524  //art::Handle<std::vector<rb::Vertex> > vert;
1525  //evt.getByLabel(fVertexLabel, vert);
1526  //now getting the list of verticies associated with the best slice
1527  //NOTE: right now it is not possible for a given slice to have more then 1
1528  //vertex so the size of vert should never be larger then 1,
1529  //so it should be safe to always use the first index
1530  if (!(fmVertex.isValid())) continue;
1531  std::vector<art::Ptr<rb::Vertex>> vert = fmVertex.at(sliceIdx);
1532 
1533  //std::cout<<"# of vertex: "<<vert.size()<<std::endl;
1534 
1535  if (vert.size() != 1) continue;
1536  //art::Ptr<rb::Vertex> v(vert, 0);
1537  //std::cout<<"Vertex Position: "<<vert[0]->GetX()<<" "<<vert[0]->GetY()<<" "<<vert[0]->GetZ()<<std::endl;
1538 
1539  //get MC truth info
1540 
1541  if( isMC ){
1542 
1543  art::ServiceHandle<cheat::BackTracker> bt;
1544  //art::Handle<std::vector<rb::CellHit> > chits;
1545  //evt.getByLabel(fCellHitLabel, chits);
1546  art::PtrVector<rb::CellHit> allhits;
1547  for(size_t h = 0; h < chits->size(); ++h)
1548  allhits.push_back(art::Ptr<rb::CellHit>(chits, h));
1549  const art::Ptr<rb::Cluster> allslice(slices,sliceIdx);
1550  std::vector<cheat::NeutrinoEffPur> iacts = bt->SliceToNeutrinoInteractions(allslice, allhits);
1551 
1552  if( iacts.empty() ){
1553  std::cout<<"no matching between slice and nu interactions!"<<std::endl;
1554  return;
1555  }
1556 
1557  // Pick the best-matching one
1558  art::Ptr<simb::MCTruth>& truth = iacts[0].neutrinoInt;
1559  // Assn interface is needlessly confusing. Have to pack stuff up like this
1560  art::PtrVector<simb::MCTruth> pv;
1561  pv.push_back(truth);
1562  // Get the flux associated with this truth
1563  //art::PtrVector<simb::MCFlux> fluxs = util::FindManyP<simb::MCFlux>(pv, evt, "geantgen", 0);
1564  art::FindManyP<simb::MCFlux> fmf(pv, evt, fMCFluxLabel);
1565  std::vector< art::Ptr<simb::MCFlux> > fluxs = fmf.at(0);
1566  simb::MCFlux flux;
1567  if(fluxs.size() == 1){
1568  flux = *fluxs[0];
1569  }
1570  else{
1571  std::cout<<"There is no assoicated neutrino intereaction !"<<std::endl;
1572  }
1573  const simb::MCNeutrino& nu = truth->GetNeutrino();
1574 
1575  //select unoscillated nue only events
1576  if( is_nue ){
1577  //std::cout<<"Select nue only events!"<<std::endl;
1578  if( !(abs(flux.fntype)==12 && abs(nu.Nu().PdgCode())==12 ) ) continue;
1579  if( !(nu.CCNC()==0) ) continue;
1580  }
1581 
1582 /*
1583  //get genie reweight table
1584  art::FindManyP<rwgt::GENIEReweightTable> gweights(pv, evt, fweightLabel);
1585  const std::vector<art::Ptr<rwgt::GENIEReweightTable>> prods = gweights.at(0);
1586  if( !prods.empty()){
1587  const unsigned int N = prods[0]->NShifts();
1588  //std::cout<<"*** # of weights "<<N<<std::endl;
1589  if( N>0 ){
1590  Nweights=N;
1591  if( Nweights>68 ) Nweights=68;//make sure index does blow-off the arrays...
1592  for(unsigned int i = 0; i < N; ++i){
1593  //std::cout<<i<<" "<<prods[0]->Plus1Sigma(i)<<" "<<prods[0]->Plus2Sigma(i)<<"; "<<prods[0]->Minus1Sigma(i)<<" "<<prods[0]->Minus2Sigma(i)<<std::endl;
1594  plus1sigma[i]=prods[0]->Plus1Sigma(i);
1595  plus2sigma[i]=prods[0]->Plus2Sigma(i);
1596  minus1sigma[i]=prods[0]->Minus1Sigma(i);
1597  minus2sigma[i]=prods[0]->Minus2Sigma(i);
1598  }
1599  }//number of weights>0
1600  }//having genie weights
1601 */
1602 
1603  trueE = nu.Nu().E();
1604  truePt = nu.Pt();
1605  trueTheta = nu.Theta();
1606  ccnc = nu.CCNC();
1607  PDG = nu.Nu().PdgCode();
1608  origPDG = flux.fntype;
1609  mode = nu.Mode();
1610  isQE = nu.InteractionType() == simb::kQE ||
1611  nu.InteractionType() == simb::kCCQE;
1612  intType = nu.InteractionType();
1613  hitNucl = nu.HitNuc();
1614  HadX = nu.X();
1615  HadY = nu.Y();
1616  HadW = nu.W();
1617  qsqr = nu.QSqr();
1618  OscWt = SimpleOscProb(flux, nu);
1619 
1620  Vx=nu.Nu().Vx();
1621  Vy=nu.Nu().Vy();
1622  Vz=nu.Nu().Vz();
1623 
1624  lepPDG = nu.Lepton().PdgCode();
1625  lepE = nu.Lepton().E();
1626  lepPx = nu.Lepton().Px();
1627  lepPy = nu.Lepton().Py();
1628  lepPz = nu.Lepton().Pz();
1629 
1630  //std::cout<<"Truth info: "<<nu.CCNC()<<" "<<nu.Nu().PdgCode()<<" "<<flux.fntype<<std::endl;
1631  //}
1632 
1633  art::Handle<std::vector<simb::MCTruth> > mclist;
1634  evt.getByLabel(fMCFluxLabel,mclist);
1635  art::Handle< std::vector<simb::MCFlux> > fllist;
1636  evt.getByLabel(fMCFluxLabel, fllist);
1637  if( mclist->size()>0 && fllist->size()>0 ){
1638  for( unsigned int i_intx=0; i_intx<mclist->size(); ++i_intx){
1639  art::Ptr<simb::MCTruth> mctruth(mclist,i_intx);
1640  const simb::MCNeutrino& mcnu = mctruth->GetNeutrino();
1641  TLorentzVector Tslnu(nu.Nu().Px(),nu.Nu().Py(),
1642  nu.Nu().Pz(),nu.Nu().E());
1643  TLorentzVector Tmcnu(mcnu.Nu().Px(),mcnu.Nu().Py(),
1644  mcnu.Nu().Pz(),mcnu.Nu().E());
1645  if(mcnu.Nu().PdgCode()==nu.Nu().PdgCode() &&
1646  mcnu.CCNC()==nu.CCNC() && mcnu.Mode()==nu.Mode() &&
1647  mcnu.Lepton().PdgCode()==nu.Lepton().PdgCode() &&
1648  Tslnu==Tmcnu ){
1649  art::Ptr<simb::MCFlux> flux(fllist,i_intx);
1650  mesonPDG = flux->ftptype;
1651  mesonPx = flux->ftpx;
1652  mesonPy = flux->ftpy;
1653  mesonPz = flux->ftpz;
1654  }
1655  }
1656  }
1657 
1658  //select the most energetic true pi0
1659  int pi0Idx=-1;
1660  double epi0=-99.;
1661  const sim::ParticleNavigator& pn = bt->ParticleNavigator();
1662  std::vector<const sim::Particle*> AllPi0;
1663  AllPi0.clear();
1664  for(int ip = 0; ip < pn.NumberOfPrimaries(); ++ip){
1665  const sim::Particle* genpar = pn.Primary(ip);
1666  if(genpar->PdgCode() != 111) continue;
1667  if( genpar->Vz()<0. || genpar->Vz()>1300. ) continue;//ND only, excluding the MC
1668  if( genpar->Vx()>200. || genpar->Vx()<-200. ) continue;
1669  if( genpar->Vy()>200. || genpar->Vy()<-200. ) continue;
1670  if( genpar->NumberDaughters() != 2 ) continue;
1671  const sim::Particle* gam1 = pn[genpar->Daughter(0)];
1672  const sim::Particle* gam2 = pn[genpar->Daughter(1)];
1673  if( gam1->PdgCode() != 22 ) continue;
1674  if( gam2->PdgCode() != 22 ) continue;
1675 
1676  AllPi0.push_back(pn.Primary(ip));
1677 
1678  if(genpar->E() > epi0){
1679  pi0Idx = ip;
1680  epi0 = genpar->E();
1681  }
1682  }//loop true particles
1683 
1684  if( AllPi0.size()>0 ){
1685  int npi0=20;
1686  if( AllPi0.size()<20 ) npi0=AllPi0.size();
1687  Npi0=npi0;
1688  for( int i=0; i<npi0; ++i ){
1689  pi0Px[i]=AllPi0[i]->Px();
1690  pi0Py[i]=AllPi0[i]->Py();
1691  pi0Pz[i]=AllPi0[i]->Pz();
1692  pi0E[i]=AllPi0[i]->E();
1693  pi0Vx[i]=AllPi0[i]->Vx();
1694  pi0Vy[i]=AllPi0[i]->Vy();
1695  pi0Vz[i]=AllPi0[i]->Vz();
1696  }
1697  }
1698 
1699  if( pi0Idx>=0 ){
1700  const sim::Particle* selectedpi0 = pn.Primary(pi0Idx);
1701  const sim::Particle* gam1 = pn[selectedpi0->Daughter(0)];
1702  const sim::Particle* gam2 = pn[selectedpi0->Daughter(1)];
1703  HasP2GG=1;
1704 
1705  Epi0=selectedpi0->E();
1706  gam1E=gam1->E();
1707  gam1Px=gam1->Px();
1708  gam1Py=gam1->Py();
1709  gam1Pz=gam1->Pz();
1710  gam2E=gam2->E();
1711  gam2Px=gam2->Px();
1712  gam2Py=gam2->Py();
1713  gam2Pz=gam2->Pz();
1714  }//has pi0
1715 
1716  }//only for MC events
1717 
1718  //select the longest kalman track
1719  std::vector<art::Ptr<rb::Track>> Ktracks = fmKTrack.at(sliceIdx);
1720 
1721  //for rock events only
1722  if( isROCK ){
1723  if( !fmKTrack.isValid() ) continue;
1724  if( Ktracks.size() != 1 ) continue;
1725  if( Ktracks[0]->TotalLength()<1300. ) continue;
1726  TVector3 dir = Ktracks[0]->Dir();
1727  double track_theta=dir.Theta();
1728  if( cos(track_theta)<0.9 ) continue;//track point forward
1729  TVector3 start = Ktracks[0]->Start();
1730  double track_startz = start.Z();
1731  if( track_startz>30. ) continue;//track starts within first 4 planes
1732  }
1733 
1734  double longestKTrack = 0;
1735  int longestKTrack_ncells=0;
1736  double longestKTrack_energy=0.;
1737  int longestKTrack_nplanes=0;
1738  float longestKTrack_phi=0.;
1739  float longestKTrack_theta=0.;
1740  float longestKTrack_startx=0.;
1741  float longestKTrack_starty=0.;
1742  float longestKTrack_startz=0.;
1743  float longestKTrack_stopx=0.;
1744  float longestKTrack_stopy=0.;
1745  float longestKTrack_stopz=0.;
1746  int longestKTrack_is3D=0;
1747 
1748  if( fmKTrack.isValid() ){
1749  int trkIndex=-1;
1750  for(unsigned int n = 0; n < Ktracks.size(); ++n){
1751  //if( !(Ktracks[n]->Is3D()) ) continue;
1752  //tracks3d.push_back(Ktracks[n]);
1753 
1754  const double L = Ktracks[n]->TotalLength();
1755  if(L > longestKTrack){
1756  longestKTrack = L;
1757  longestKTrack_ncells=Ktracks[n]->NCell();
1758  longestKTrack_energy=Ktracks[n]->TotalGeV();
1759  longestKTrack_nplanes=Ktracks[n]->ExtentPlane();
1760  TVector3 dir = Ktracks[n]->Dir();
1761  longestKTrack_phi=dir.Phi();
1762  longestKTrack_theta=dir.Theta();
1763  TVector3 start = Ktracks[n]->Start();
1764  TVector3 stop = Ktracks[n]->Stop();
1765  longestKTrack_startx = start.X();
1766  longestKTrack_starty = start.Y();
1767  longestKTrack_startz = start.Z();
1768  longestKTrack_stopx = stop.X();
1769  longestKTrack_stopy = stop.Y();
1770  longestKTrack_stopz = stop.Z();
1771  if( Ktracks[n]->Is3D() ) longestKTrack_is3D=1;
1772  trkIndex=n;
1773  }//for longest track
1774  }//loop all tracks
1775 
1776  if( trkIndex>-1 ){
1777  int ncells_trk=0;
1778 
1779  std::vector<float> tcellX;
1780  std::vector<float> tcellY;
1781  std::vector<float> tcellZ;
1782  std::vector<float> tcellE;
1783  std::vector<int> tcellPlane;
1784  std::vector<int> tcellCell;
1785  std::vector<int> tcellADC;
1786  std::vector<float> tcellPECorr;
1787  std::vector<float> tcellTNS;
1788  std::vector<float> tcellW;
1789 
1790 
1791  for( unsigned int icell=0; icell<Ktracks[trkIndex]->NCell(); ++ icell){
1792  const art::Ptr<rb::CellHit>& chit = Ktracks[trkIndex]->Cell(icell);
1793  const rb::RecoHit rhit(cal->MakeRecoHit(*chit, Ktracks[trkIndex]->W(chit.get())));
1794  if( !rhit.IsCalibrated() ) continue;
1795 
1796  ++ncells_trk;
1797 
1798  tcellX.push_back(rhit.X());
1799  tcellY.push_back(rhit.Y());
1800  tcellZ.push_back(rhit.Z());
1801  tcellE.push_back(rhit.GeV());
1802  tcellADC.push_back(chit->ADC());
1803  tcellPlane.push_back(chit->Plane());
1804  tcellCell.push_back(chit->Cell());
1805 
1806  tcellW.push_back(Ktracks[trkIndex]->W(chit.get()));
1807  tcellPECorr.push_back(rhit.PECorr());
1808  }//loop track associated cells
1809 
1810  trackcellNcell = ncells_trk;
1811  for( int ic=0; ic<ncells_trk; ++ic ){
1812  trackcellX[ic]=tcellX[ic];
1813  trackcellY[ic]=tcellY[ic];
1814  trackcellZ[ic]=tcellZ[ic];
1815 
1816  trackcellE[ic]=tcellE[ic];
1817  trackcellADC[ic]=tcellADC[ic];
1818  trackcellPlane[ic]=tcellPlane[ic];
1819  trackcellCell[ic]=tcellCell[ic];
1820  trackcellPECorr[ic]=tcellPECorr[ic];
1821  trackcellW[ic]=tcellW[ic];
1822  }//loop selected cells
1823  }//longest track
1824  }//has tracks
1825 
1826  //std::cout<<"# of tracks: "<<Ktracks.size()<<std::endl;
1827  //std::cout<<"# of cells: "<<slice.NCell()<<std::endl;
1828  //mip cells and reco energy for slicer
1829  double nCells = slice.NCell();
1830  int lengthPlanes = slice.ExtentPlane();
1831  if( slice.ExtentPlane()>192 ) lengthPlanes=192;
1832 
1833  double recoEnergy=0.;
1834  double nmip=0.;
1835  double Ntotal=0.;
1836  for( unsigned int icell=0; icell<slice.NCell(); ++ icell){
1837  const art::Ptr<rb::CellHit>& chit = slice.Cell(icell);
1838  const rb::RecoHit rhit(cal->MakeRecoHit(*chit, slice.W(chit.get())));
1839 
1840  if( !rhit.IsCalibrated() ) continue;
1841  recoEnergy += rhit.GeV();
1842 
1843  if( rhit.PECorr()>100. && rhit.PECorr()<245. ) nmip++;
1844  Ntotal++;
1845  }
1846  //preselection
1847  if( nCells<20 || nCells>500 ) continue;
1848 
1849  if( !isROCK ){
1850  if( fabs(vert[0]->GetX())>140. ) continue;
1851  if( fabs(vert[0]->GetY())>140. ) continue;
1852  if( vert[0]->GetZ()<100. ) continue;
1853  if( vert[0]->GetZ()>700. ) continue;
1854  }//fiducial cuts for non-rock events
1855 
1856  sl_runnum = evt.run();
1857  sl_subrunnum = evt.subRun();
1858  sl_evtnum = evt.id().event();
1859  sl_Nslices=nslice;
1860 
1861  sl_Ndeaddcm=ndeaddcms;
1862  if( ndeaddcms>0 ){
1863  for( int id=0; id<ndeaddcms; ++id )
1864  sl_deaddcm_location[id]=dcmlocation[id];
1865  }
1866  sl_Nouttimehits=ncells_outtime;
1867  sl_Nintimehits=ncells_intime;
1868  sl_Ndb1apd_24cells=Napd_24cells[0];
1869  sl_Ndb2apd_24cells=Napd_24cells[1];
1870  sl_Ndb3apd_24cells=Napd_24cells[2];
1871  sl_Nmcapd_24cells=Napd_24cells[3];
1872  sl_Ndb1apd_28cells=Napd_28cells[0];
1873  sl_Ndb2apd_28cells=Napd_28cells[1];
1874  sl_Ndb3apd_28cells=Napd_28cells[2];
1875  sl_Nmcapd_28cells=Napd_28cells[3];
1876  sl_Fouttimehits_noisyapd=ratio_outtimehits_deadapds;
1877 
1878  if( !isMC && !isROCK){
1879  if( spillPot->goodbeam ) sl_goodSpill=1;
1880  else sl_goodSpill=0;
1881  sl_spillPOT = spillPot->spillpot;
1882  sl_horncurrent = spillPot->hornI;
1883  sl_posX = spillPot->posx;
1884  sl_posY = spillPot->posy;
1885  sl_widthX = spillPot->widthx;
1886  sl_widthY = spillPot->widthy;
1887  sl_deltaspilltimensec = spillPot->deltaspilltimensec;
1888  }
1889 
1890  NKtracks = Ktracks.size();
1891  KtrackLength = longestKTrack;
1892  KtrackEnergy = longestKTrack_energy;
1893  KtrackNcells = longestKTrack_ncells;
1894  KtrackEnergy_ratio = longestKTrack_energy/recoEnergy;
1895  KtrackNcells_ratio = longestKTrack_ncells/nCells;
1896  KtrackNplanes = longestKTrack_nplanes;
1897  KtrackPhi = longestKTrack_phi;
1898  KtrackTheta = longestKTrack_theta;
1899  KtrackStartX = longestKTrack_startx;
1900  KtrackStartY = longestKTrack_starty;
1901  KtrackStartZ = longestKTrack_startz;
1902  KtrackStopX = longestKTrack_stopx;
1903  KtrackStopY = longestKTrack_stopy;
1904  KtrackStopZ = longestKTrack_stopz;
1905  Ktrack3D = longestKTrack_is3D;
1906 
1907  sliceIndex = sliceIdx;
1908  recoE = recoEnergy;
1909  Ncells = nCells;
1910  Nplanes = lengthPlanes;
1911  MinPlane = slice.MinPlane();
1912  slice_meantns=slice.MeanTNS()/1000.;
1913  slice_mintns=slice.MinTNS()/1000.;
1914  slice_maxtns=slice.MaxTNS()/1000.;
1915  Nmip=nmip;
1916  Fmip=nmip/Ntotal;
1917 
1918  recoVx=vert[0]->GetX();
1919  recoVy=vert[0]->GetY();
1920  recoVz=vert[0]->GetZ();
1921 
1922 
1923  if( !isROCK ){
1924 
1925 /*
1926  //RVP
1927  art::FindManyP<rvp::RVP> fmRvp(slices, evt, frvpLabel);
1928  std::vector<art::Ptr<rvp::RVP>> rvps = fmRvp.at(sliceIdx);
1929  if( rvps.size()>0 ){
1930  rvp = rvps[0]->Value();
1931  }
1932 */
1933  art::FindOneP<slid::EventLID> fslLID(slices, evt, flidLabel);
1934  cet::maybe_ref<art::Ptr<slid::EventLID> const> rlid(fslLID.at(sliceIdx));
1935  art::Ptr<slid::EventLID> elid = rlid.ref();
1936  if( elid ){
1937  lid = elid->Value();
1938  lide = elid->ValueWithE();
1939  //std::cout<<"*** LID is "<<lid<<std::endl;
1940  }
1941 
1942  art::FindMany<slid::ShowerLID> fmShwLID(slices, evt, flidLabel);
1943  //std::vector<art::Ptr<slid::ShowerLID> lids = fmShwLID.at(sliceIdx);
1944  std::vector<const slid::ShowerLID* > lids = fmShwLID.at(sliceIdx);
1945  sort(lids.begin(),lids.end(),slid::CompareByEnergy);
1946  art::FindOneP<rb::Shower> fos(lids, evt, flidLabel);
1947  int nshowers=lids.size();
1948  if( nshowers>20 ) nshowers=20;
1949  Nshowers=nshowers;
1950  if( lids.size()>0 ){
1951 
1952  int shw_ncell=0;//all cells assoicated to prong
1953  std::vector<int> shw_index;
1954  std::vector<float> shwX;
1955  std::vector<float> shwY;
1956  std::vector<float> shwZ;
1957  std::vector<float> shwE;
1958  std::vector<int> shwPlane;
1959  std::vector<int> shwCell;
1960  std::vector<int> shwADC;
1961  std::vector<float> shwPECorr;
1962 
1963  for( int is=0; is<nshowers; ++is ){
1964  shwLID[is] = lids[is]->Value();
1965  //std::cout<<"shower "<<is<<" lid is "<<lids[is]->Value()<<std::endl;
1966  cet::maybe_ref<art::Ptr<rb::Shower> const> rshw(fos.at(is));
1967  art::Ptr<rb::Shower> shw = rshw.ref();
1968  if( shw->Is3D() ) shw3D[is]=1;
1969  else shw3D[is]=0;
1970  shwNCells[is] = shw->NCell();
1971  shwNCellsXview[is] = shw->NXCell();
1972  shwNCellsYview[is] = shw->NYCell();
1973  shwLength[is] = shw->TotalLength();
1974  shwStartX[is] = shw->Start().X();
1975  shwStartY[is] = shw->Start().Y();
1976  shwStartX[is] = shw->Start().Z();
1977  shwStopX[is] = shw->Stop().X();
1978  shwStopY[is] = shw->Stop().Y();
1979  shwStopZ[is] = shw->Stop().Z();
1980 
1981 /*
1982  TVector3 Sstart = shw->Start();
1983  TVector3 Sdir = shw->Dir();
1984  TVector3 Sstop = Sstart + (shw->TotalLength())*Sdir;
1985  std::cout<<"stop x "<<shw->Stop().X()<<" vs "<<Sstop.X()<<std::endl;
1986  std::cout<<"stop y "<<shw->Stop().Y()<<" vs "<<Sstop.Y()<<std::endl;
1987  std::cout<<"stop z "<<shw->Stop().Z()<<" vs "<<Sstop.Z()<<std::endl;
1988  std::cout<<"********************************************"<<std::endl;
1989 */
1990 
1991  shwTheta[is] = shw->Dir().Theta();
1992  shwEcells[is] = shw->TotalGeV();//simple sum of the estimated GeV of all the hits
1993  shwEnu[is] = shw->CalorimetricEnergy();//simple estimate of neutrino energy
1994 
1995  for( unsigned int icell=0; icell<shw->NCell(); ++icell ){
1996  const art::Ptr<rb::CellHit>& chit = shw->Cell(icell);
1997  const rb::RecoHit rhit(cal->MakeRecoHit(*chit, shw->W(chit.get())));
1998  if( !rhit.IsCalibrated() ) continue;
1999  shw_ncell++;
2000 
2001  shw_index.push_back(is);
2002  shwX.push_back(rhit.X());
2003  shwY.push_back(rhit.Y());
2004  shwZ.push_back(rhit.Z());
2005  shwE.push_back(rhit.GeV());
2006  shwPlane.push_back(chit->Plane());
2007  shwCell.push_back(chit->Cell());
2008  shwADC.push_back(chit->ADC());
2009  shwPECorr.push_back(rhit.PECorr());
2010  }//loop shower associated cells
2011  }//loop all showers
2012 
2013  if( shw_ncell>500 ) shw_ncell=500;
2014  shwcellNcell = shw_ncell;
2015  for( int ic=0; ic<shw_ncell; ++ic ){
2016  shwIndex[ic] = shw_index[ic];
2017  shwcellX[ic] = shwX[ic];
2018  shwcellY[ic] = shwY[ic];
2019  shwcellZ[ic] = shwZ[ic];
2020  shwcellE[ic] = shwE[ic];
2021  shwcellADC[ic] = shwADC[ic];
2022  shwcellPECorr[ic] = shwPECorr[ic];
2023  shwcellPlane[ic] = shwPlane[ic];
2024  shwcellCell[ic] = shwCell[ic];
2025  }//loop all save showers' hits
2026  }//have shower ids
2027 
2028 
2029 /*
2030  art::FindManyP<lem::PIDDetails> fmLem(slices, evt, flemLabel);
2031  std::vector<art::Ptr<lem::PIDDetails>> lems = fmLem.at(sliceIdx);
2032  if( lems.size()>0 ){
2033  lem = lems[0]->Value();
2034  }
2035 */
2036  //std::cout<<rvps.size()<<" "<<lems.size()<<" "<<lids.size()<<std::endl;
2037  //std::cout<<"PID: "<<rvp<<" "<<lem<<" "<<lid<<std::endl;
2038  }//for non-rockmuon events
2039 
2040  //reco prongs
2041  art::FindManyP<rb::Prong> fmProng3D(vert, evt, art::InputTag(fProngLabel,fInstLabel3D));
2042  std::vector<art::Ptr<rb::Prong>> SelectedProng = fmProng3D.at(0);
2043 
2044  //std::cout<<"number of prongs: "<<SelectedProng.size()<<std::endl;
2045 
2046  int nprongs=SelectedProng.size();
2047  if( nprongs>20 ) nprongs=20;
2048  Nprongs=nprongs;
2049  std::vector<double> ProngEnergy;
2050  std::vector<double> ProngEnergyXView;
2051  std::vector<double> ProngEnergyYView;
2052  std::vector<double> ProngPhi;
2053 
2054  if( SelectedProng.size()<1 ) continue;
2055 
2056  if( SelectedProng.size()>0 ){
2057 
2058  int allprong_ncell=0;//all cells assoicated to prong
2059  std::vector<int> prongIndex;
2060  std::vector<float> pcellX;
2061  std::vector<float> pcellY;
2062  std::vector<float> pcellZ;
2063  std::vector<float> pcellE;
2064  std::vector<int> pcellPlane;
2065  std::vector<int> pcellCell;
2066  std::vector<int> pcellADC;
2067  std::vector<float> pcellPECorr;
2068  std::vector<float> pcellW;
2069 
2070  for(int ip=0; ip<nprongs; ++ip ){
2071  //reco direction
2072  TVector3 dir1 = SelectedProng[ip]->Dir();
2073  //recoEta[ip]=dir1.Eta();
2074  prongPhi[ip]=dir1.Phi();
2075  prongTheta[ip]=dir1.Theta();
2076 
2077  TVector3 start1 = SelectedProng[ip]->Start();
2078  prongStartX[ip] = start1.X();
2079  prongStartY[ip] = start1.Y();
2080  prongStartZ[ip] = start1.Z();
2081 
2082  prongLength[ip] = SelectedProng[ip]->TotalLength();
2083 
2084  TVector3 Pstop = start1 + (SelectedProng[ip]->TotalLength())*dir1;
2085  prongStopX[ip] = Pstop.X();
2086  prongStopY[ip] = Pstop.Y();
2087  prongStopZ[ip] = Pstop.Z();
2088 
2089  prongNplanes[ip] = SelectedProng[ip]->ExtentPlane();
2090 
2091  //reco momentum
2092  double Eprong=0.;
2093  double EprongX=0.;
2094  double EprongY=0.;
2095 
2096  int npcells=0;;
2097  int nxpcells=0;
2098  int nypcells=0;
2099 
2100  for( unsigned int icell=0; icell<SelectedProng[ip]->NCell(); ++ icell){
2101  const art::Ptr<rb::CellHit>& chit = SelectedProng[ip]->Cell(icell);
2102  const rb::RecoHit rhit(cal->MakeRecoHit(*chit, SelectedProng[ip]->W(chit.get())));
2103  double wx=SelectedProng[ip]->W(chit.get());
2104  double wy=wx;
2105 
2106  if( !rhit.IsCalibrated() ) continue;
2107  Eprong += rhit.GeV();
2108  npcells++;
2109 
2110  if( chit->View() == geo::kX ){
2111  EprongX += rhit.GeV();
2112  ++nxpcells;
2113  pcellW.push_back(wx);
2114  }
2115  if( chit->View() == geo::kY ){
2116  EprongY += rhit.GeV();
2117  ++nypcells;
2118  pcellW.push_back(wy);
2119  }
2120 
2121  ++allprong_ncell;
2122  prongIndex.push_back(ip);
2123  pcellX.push_back(rhit.X());
2124  pcellY.push_back(rhit.Y());
2125  pcellZ.push_back(rhit.Z());
2126  pcellE.push_back(rhit.GeV());
2127 
2128  pcellADC.push_back(chit->ADC());
2129  pcellPlane.push_back(chit->Plane());
2130  pcellCell.push_back(chit->Cell());
2131  pcellPECorr.push_back(rhit.PECorr());
2132 
2133  }//loop prong associated cells
2134  ProngEnergy.push_back(Eprong);
2135  ProngEnergyXView.push_back(EprongX);
2136  ProngEnergyYView.push_back(EprongY);
2137  ProngPhi.push_back(dir1.Phi());
2138  prongEnergy[ip]=Eprong;
2139  prongEnergyXView[ip]=EprongX;
2140  prongEnergyYView[ip]=EprongY;
2141  prongNCells[ip]=npcells;
2142  prongNCellsXView[ip]=nxpcells;
2143  prongNCellsYView[ip]=nypcells;
2144  }//loop prongs
2145 
2146  if( allprong_ncell>500 ) allprong_ncell=500;
2147  allprongcellNcell=allprong_ncell;
2148  for( int ic=0; ic<allprong_ncell; ++ic ){
2149  allprongIndex[ic]=prongIndex[ic];
2150  allprongcellX[ic]=pcellX[ic];
2151  allprongcellY[ic]=pcellY[ic];
2152  allprongcellZ[ic]=pcellZ[ic];
2153  allprongcellE[ic]=pcellE[ic];
2154 
2155  allprongcellADC[ic]=pcellADC[ic];
2156  allprongcellPlane[ic]=pcellPlane[ic];
2157  allprongcellCell[ic]=pcellCell[ic];
2158  allprongcellPECorr[ic]=pcellPECorr[ic];
2159  allprongcellW[ic]=pcellW[ic];
2160  }//all prongs
2161 
2162  double ratio=1.;
2163  double dphi=0.;
2164  int Index1=-1;
2165  int Index2=-1;
2166  double MaxEnergy1=-999.;
2167  double MaxEnergy2=-999.;
2168  if( nprongs>0 ){
2169  for( int i=0;i<nprongs; ++i ){
2170  if( MaxEnergy1<ProngEnergy[i] ){
2171  MaxEnergy1=ProngEnergy[i];
2172  Index1=i;
2173  }
2174  }//leading
2175  }
2176 
2177  //if( prongLength[Index1]<50. ) continue;
2178  //if( prongLength[Index1]>600. ) continue;
2179  //if( fabs(prongStopX[Index1])>180. ) continue;
2180  //if( fabs(prongStopY[Index1])>180. ) continue;
2181  //if( prongStopZ[Index1]<50. ) continue;
2182  //if( prongStopZ[Index1]>1250. ) continue;
2183 
2184  if( nprongs>1 ){
2185  for( int i=0;i<nprongs; ++i ){
2186  if( i==Index1 ) continue;
2187  if( MaxEnergy2<ProngEnergy[i] ){
2188  MaxEnergy2=ProngEnergy[i];
2189  Index2=i;
2190  }
2191  }//sub-leading
2192  double Eden=ProngEnergy[Index1]+ProngEnergy[Index2];
2193  double Enum=ProngEnergy[Index1]-ProngEnergy[Index2];
2194  if( Eden>0.0 ) ratio=Enum/Eden;
2195  dphi=fabs(ProngPhi[Index1]-ProngPhi[Index2]);
2196  if( dphi>2.*TMath::Pi() ) dphi=dphi-2.*TMath::Pi();
2197  if( dphi>TMath::Pi() ) dphi=2.*TMath::Pi()-dphi;
2198  }
2199 
2200  Ebalance=ratio;
2201  Dphi=dphi;
2202 
2203  if( Index1>=0 ){
2204 
2205  int p1_ncells=0;
2206  int p1_mipcells=0;
2207  std::vector<float> p1cellX;
2208  std::vector<float> p1cellY;
2209  std::vector<float> p1cellZ;
2210  std::vector<float> p1cellE;
2211  std::vector<int> p1cellPlane;
2212  std::vector<int> p1cellCell;
2213  std::vector<int> p1cellADC;
2214  std::vector<float> p1cellPE;
2215  std::vector<float> p1cellPECorr;
2216  std::vector<int> p1cellXView;
2217  std::vector<float> p1cellTNS;
2218  std::vector<float> p1cellW;
2219 
2220  for( unsigned int icell=0; icell<SelectedProng[Index1]->NCell(); ++ icell){
2221  const art::Ptr<rb::CellHit>& chit = SelectedProng[Index1]->Cell(icell);
2222  const rb::RecoHit rhit(cal->MakeRecoHit(*chit, SelectedProng[Index1]->W(chit.get())));
2223  double wx=SelectedProng[Index1]->W(chit.get());
2224  double wy=wx;
2225  if( !rhit.IsCalibrated() ) continue;
2226  ++p1_ncells;
2227  if( rhit.PECorr()>100. && rhit.PECorr()<240. ) ++p1_mipcells;
2228 
2229  if( chit->View()==geo::kX ){
2230  p1cellW.push_back(wx);
2231  }//XZ view
2232  else{
2233  p1cellW.push_back(wy);
2234  }//YZ view
2235 
2236  p1cellX.push_back(rhit.X());
2237  p1cellY.push_back(rhit.Y());
2238  p1cellZ.push_back(rhit.Z());
2239  p1cellE.push_back(rhit.GeV());
2240 
2241  p1cellADC.push_back(chit->ADC());
2242  p1cellPlane.push_back(chit->Plane());
2243  p1cellCell.push_back(chit->Cell());
2244  p1cellPECorr.push_back(rhit.PECorr());
2245 
2246  }//loop energetic prong associated cells
2247 
2248  if( p1_ncells>500 ) p1_ncells=500;
2249  prong1cellNcell=p1_ncells;
2250  if(p1_ncells>0 ) prong1Fmip = (float)p1_mipcells/p1_ncells;
2251  for( int ic=0; ic<p1_ncells; ++ic ){
2252  prong1cellX[ic]=p1cellX[ic];
2253  prong1cellY[ic]=p1cellY[ic];
2254  prong1cellZ[ic]=p1cellZ[ic];
2255  prong1cellE[ic]=p1cellE[ic];
2256 
2257  prong1cellADC[ic]=p1cellADC[ic];
2258  prong1cellPlane[ic]=p1cellPlane[ic];
2259  prong1cellCell[ic]=p1cellCell[ic];
2260  prong1cellPECorr[ic]=p1cellPECorr[ic];
2261  prong1cellW[ic]=p1cellW[ic];
2262  }//loop selected cells
2263  }//for most energecti prong
2264 
2265  }//has at least one prong
2266 
2267  //slice associated cell information
2268  int ncells_slice=0.;
2269  std::vector<float> slicecellX;
2270  std::vector<float> slicecellY;
2271  std::vector<float> slicecellZ;
2272  std::vector<float> slicecellE;
2273  std::vector<int> slicecellPlane;
2274  std::vector<int> slicecellCell;
2275  std::vector<int> slicecellADC;
2276  std::vector<float> slicecellPE;
2277  std::vector<float> slicecellPECorr;
2278  std::vector<int> slicecellXView;
2279  std::vector<float> slicecellTNS;
2280  std::vector<float> slicecellW;
2281  std::vector<int> slicecellNoise;
2282 
2283  for( unsigned int icell=0; icell<slice.NCell(); ++ icell){
2284  const art::Ptr<rb::CellHit>& chit = slice.Cell(icell);
2285  const rb::RecoHit rhit(cal->MakeRecoHit(*chit, slice.W(chit.get())));
2286 
2287 
2288  if( !rhit.IsCalibrated() ) continue;
2289 
2290  ncells_slice++;
2291 
2292  slicecellW.push_back(slice.W(chit.get()));
2293 
2294  slicecellX.push_back(rhit.X());
2295  slicecellY.push_back(rhit.Y());
2296  slicecellZ.push_back(rhit.Z());
2297  slicecellE.push_back(rhit.GeV());
2298 
2299  slicecellADC.push_back(chit->ADC());
2300  slicecellPlane.push_back(chit->Plane());
2301  slicecellCell.push_back(chit->Cell());
2302  slicecellPECorr.push_back(rhit.PECorr());
2303  slicecellTNS.push_back(chit->TNS()/1000.);
2304  }
2305 
2306  if( ncells_slice>500 ) ncells_slice=500;
2307  slcellNcell=ncells_slice;
2308 
2309  for( int ic=0; ic<ncells_slice; ++ic ){
2310 
2311  slcellX[ic]=slicecellX[ic];
2312  slcellY[ic]=slicecellY[ic];
2313  slcellZ[ic]=slicecellZ[ic];
2314  slcellE[ic]=slicecellE[ic];
2315 
2316  slcellADC[ic]=slicecellADC[ic];
2317  slcellPlane[ic]=slicecellPlane[ic];
2318  slcellCell[ic]=slicecellCell[ic];
2319  slcellPECorr[ic]=slicecellPECorr[ic];
2320  slcellW[ic]=slicecellW[ic];
2321  slcellTNS[ic]=slicecellTNS[ic];
2322  }//loop selected cells
2323 
2324  //reco prong View(): 0 for XZ, 1 for YZ, 2 for 3D
2325  //study 2D prongs
2326  art::FindManyP<rb::Prong> fmProng2D(vert, evt, art::InputTag(fProngLabel,fInstLabel2D));
2327  std::vector<art::Ptr<rb::Prong>> SelectedProng2D = fmProng2D.at(0);
2328 
2329  //std::cout<<"# of 2D prongs: "<<SelectedProng2D.size()<<std::endl;
2330 
2331  std::vector<int> Prong2DInXView;
2332  std::vector<double> ProngEnergy2D;
2333  int nprongs2D=SelectedProng2D.size();
2334  if( nprongs2D>20 ) nprongs2D=20;
2335  Nprongs2D = nprongs2D;
2336 
2337  if( nprongs2D>0 ){
2338  int allprong2d_ncell=0;//all cells assoicated to prong
2339  std::vector<int> prong2dIndex;
2340  std::vector<float> p2dcellX;
2341  std::vector<float> p2dcellY;
2342  std::vector<float> p2dcellZ;
2343  std::vector<float> p2dcellE;
2344  std::vector<int> p2dcellPlane;
2345  std::vector<int> p2dcellCell;
2346  std::vector<int> p2dcellADC;
2347  std::vector<float> p2dcellPECorr;
2348  std::vector<float> p2dcellW;
2349 
2350  for(int ip=0; ip<nprongs2D; ++ip ){
2351  prongNCells2D[ip]=SelectedProng2D[ip]->NCell();
2352  int prongview=SelectedProng2D[ip]->View();
2353 
2354  prong2DInXView[ip]=prongview;//0 for XZ view
2355  prongLength2D[ip]=SelectedProng2D[ip]->TotalLength();
2356 
2357  TVector3 dir1 = SelectedProng2D[ip]->Dir();
2358  prongPhi2D[ip]=dir1.Phi();
2359  prongTheta2D[ip]=dir1.Theta();
2360 
2361  TVector3 start1 = SelectedProng2D[ip]->Start();
2362  prongStartX2D[ip] = start1.X();
2363  prongStartY2D[ip] = start1.Y();
2364  prongStartZ2D[ip] = start1.Z();
2365 
2366  TVector3 Pstop = start1 + (SelectedProng2D[ip]->TotalLength())*dir1;
2367  prongStopX2D[ip] = Pstop.X();
2368  prongStopY2D[ip] = Pstop.Y();
2369  prongStopZ2D[ip] = Pstop.Z();
2370 
2371  Prong2DInXView.push_back(prongview);
2372  //reco momentum
2373  double Eprong=0.;
2374 
2375  for( unsigned int icell=0; icell<SelectedProng2D[ip]->NCell(); ++ icell){
2376  const art::Ptr<rb::CellHit>& chit = SelectedProng2D[ip]->Cell(icell);
2377  const TVector3 mean = SelectedProng2D[ip]->MeanXYZ();
2378  const rb::RecoHit rhit(cal->MakeRecoHit(*chit, SelectedProng2D[ip]->W(chit.get())));
2379 
2380  if( !rhit.IsCalibrated() ) continue;
2381  Eprong += rhit.GeV();
2382 
2383  ++allprong2d_ncell;
2384 
2385  prong2dIndex.push_back(ip);
2386  p2dcellX.push_back(rhit.X());
2387  p2dcellY.push_back(rhit.Y());
2388  p2dcellZ.push_back(rhit.Z());
2389  p2dcellE.push_back(rhit.GeV());
2390 
2391  p2dcellADC.push_back(chit->ADC());
2392  p2dcellPlane.push_back(chit->Plane());
2393  p2dcellCell.push_back(chit->Cell());
2394  p2dcellPECorr.push_back(rhit.PECorr());
2395  }
2396 
2397  prongEnergy2D[ip]=Eprong;
2398  ProngEnergy2D.push_back(Eprong);
2399  }//loop 2D prongs
2400 
2401  if( allprong2d_ncell>500 ) allprong2d_ncell=500;
2402  allprong2DcellNcell=allprong2d_ncell;
2403  for( int ic=0; ic<allprong2d_ncell; ++ic ){
2404  allprong2DIndex[ic]=prong2dIndex[ic];
2405  allprong2DcellX[ic]=p2dcellX[ic];
2406  allprong2DcellY[ic]=p2dcellY[ic];
2407  allprong2DcellZ[ic]=p2dcellZ[ic];
2408  allprong2DcellE[ic]=p2dcellE[ic];
2409 
2410  allprong2DcellADC[ic]=p2dcellADC[ic];
2411  allprong2DcellPlane[ic]=p2dcellPlane[ic];
2412  allprong2DcellCell[ic]=p2dcellCell[ic];
2413  allprong2DcellPECorr[ic]=p2dcellPECorr[ic];
2414  }//all prongs
2415 
2416  int index3D=-1;
2417  double max3D=-99.;
2418  for( int i=0; i<nprongs; ++i ){
2419  if( max3D<ProngEnergy[i] ){
2420  max3D=ProngEnergy[i];
2421  index3D=i;
2422  }
2423  }//loop 3D prongs
2424  int index2D=-1;
2425  double max2D=-99.;
2426  for( int i=0;i<nprongs2D; ++i ){
2427  if( max2D<ProngEnergy2D[i] ){
2428  max2D=ProngEnergy2D[i];
2429  index2D=i;
2430  }
2431  }//loop 2D prongs
2432  double ebalance2D=1.;
2433  if( nprongs>0 && nprongs2D>0 ){
2434  double Eden=0.;
2435  double Enum=0.;
2436  if( Prong2DInXView[index2D]==0 ){//XZ view
2437  Eden=ProngEnergy2D[index2D]+ProngEnergyXView[index3D];
2438  Enum=ProngEnergy2D[index2D]-ProngEnergyXView[index3D];
2439  }
2440  else {//YZ view
2441  Eden=ProngEnergy2D[index2D]+ProngEnergyYView[index3D];
2442  Enum=ProngEnergy2D[index2D]-ProngEnergyYView[index3D];
2443  }
2444  if( Eden>0. ) ebalance2D=fabs(Enum)/Eden;
2445  }
2446  Ebalance2D=ebalance2D;
2447  }//having 2D prongs
2448 
2449  _otree->Fill();
2450 
2451  //_evttree->Fill();
2452  }//end loop of slice
2453 
2454  //_evttree->Fill();
2455  return;
2456  }//end analyze
2457 
2458  void NCAna::endJob()
2459  {
2460  art::ServiceHandle<art::TFileService> tfs;
2461  TH1* hPOT = tfs->make<TH1F>("hTotalPOT", ";; POT", 1, 0, 1);
2462  hPOT->Fill(.5, fTotalPOT);
2463  TH1* hSPILL = tfs->make<TH1F>("hTotalSpill",";; SPILL", 1, 0, 1);
2464  hSPILL->Fill(.5, fTotalSpill);
2465  }
2466 
2467  double NCAna::SimpleOscProb(const simb::MCFlux& flux,
2468  const simb::MCNeutrino& nu) const
2469  {
2470  const double L = 810.; // km
2471  const double ldm = 2.35e-3; // large delta m^2, m23
2472  const double sin22t13 = 0.1; // sin^2(2theta_13)
2473  const double sin22t23 = 1.; // maximal sin^2(2theta_23)
2474  const double ssth23 = 0.5; // sin^2(theta_23) corresponding to above
2475 
2476  // Signal
2477  if(nu.CCNC() == 0 && abs(nu.Nu().PdgCode()) == 12 && abs(flux.fntype) == 14){
2478  // Nue appearance
2479  return sin22t13*ssth23*util::sqr(sin(1.267*ldm*L/nu.Nu().E()));
2480  }
2481  if(nu.CCNC() == 0 && abs(nu.Nu().PdgCode()) == 14 && abs(flux.fntype) == 14){
2482  // CC mu disappearance
2483  return 1-sin22t23*util::sqr(sin(1.267*ldm*L/nu.Nu().E()));
2484  }
2485 
2486  // Background
2487  if(nu.CCNC() == 0 && abs(nu.Nu().PdgCode()) == 12 && abs(flux.fntype) == 12){
2488  // Beam nue
2489  return 1-sin22t13*util::sqr(sin(1.267*ldm*L/nu.Nu().E()));
2490  }
2491  if(nu.CCNC() == 1 && nu.Nu().PdgCode() == flux.fntype){
2492  // NC
2493  return 1;
2494  }
2495  if(nu.CCNC() == 1 && nu.Nu().PdgCode() != flux.fntype){
2496  // Swapped NC. We play clever reweighting games so this is OK
2497  return 1;
2498  }
2499  if(nu.CCNC() == 0 && abs(nu.Nu().PdgCode()) == 14 && abs(flux.fntype) == 12){
2500  // CC mu appearance
2501  return sin22t13*ssth23*util::sqr(sin(1.267*ldm*L/nu.Nu().E()));
2502  }
2503  // Don't know what this is
2504  return 0;
2505  }
2506 }//end namespace
2507 
2508 namespace ncs{
2509 
2510  DEFINE_ART_MODULE(NCAna)
2511 
2512 }
simb::MCParticle::E
double E(const int i=0) const
Definition: MCParticle.h:232
xml_schema::id
::xsd::cxx::tree::id< char, ncname > id
Definition: Database.h:165
ncs::NCAna::slcellX
float slcellX[500]
Definition: NCAna_module.cc:300
rb::CellHit::TNS
float TNS() const
Definition: CellHit.h:46
ncs::NCAna::shwcellCell
int shwcellCell[500]
Definition: NCAna_module.cc:454
ncs::NCAna::allprong2DcellCell
int allprong2DcellCell[500]
Definition: NCAna_module.cc:422
ncs::NCAna::prongStartY
float prongStartY[20]
Definition: NCAna_module.cc:320
ncs::NCAna::allprongcellPECorr
float allprongcellPECorr[500]
Definition: NCAna_module.cc:342
art::Event::subRun
SubRunNumber_t subRun() const
Definition: Event.h:72
BackTracker.h
back track the reconstruction to the simulation
ncs::NCAna::sl_Ndb1apd_24cells
int sl_Ndb1apd_24cells
Definition: NCAna_module.cc:175
ncs::NCAna::KtrackNplanes
float KtrackNplanes
Definition: NCAna_module.cc:268
rb::Cluster::W
virtual double W(const rb::CellHit *chit) const
Estimate the unmeasured coordinate of chit.
Definition: Cluster.cxx:121
ncs::NCAna::fInstLabel3D
std::string fInstLabel3D
label for prong
Definition: NCAna_module.cc:107
simb::MCParticle::PdgCode
int PdgCode() const
Definition: MCParticle.h:211
simb::MCNeutrino::CCNC
int CCNC() const
Definition: MCNeutrino.h:148
ncs::NCAna::slcellE
float slcellE[500]
Definition: NCAna_module.cc:303
ncs::NCAna::allprong2DcellE
float allprong2DcellE[500]
Definition: NCAna_module.cc:419
ncs::NCAna::prongEnergy
float prongEnergy[20]
Definition: NCAna_module.cc:316
simb::MCNeutrino::QSqr
double QSqr() const
Definition: MCNeutrino.h:157
ncs::NCAna::allprong2DcellY
float allprong2DcellY[500]
Definition: NCAna_module.cc:417
ncs::NCAna::sl_horncurrent
double sl_horncurrent
Definition: NCAna_module.cc:166
simb::MCNeutrino::Theta
double Theta() const
angle between incoming and outgoing leptons, in radians
Definition: MCNeutrino.cxx:63
rb::Cluster::NCell
unsigned int NCell(geo::View_t view) const
Number of cells in view view.
Definition: Cluster.cxx:134
simb::MCTruth::GetNeutrino
const simb::MCNeutrino & GetNeutrino() const
Definition: MCTruth.h:76
ncs::NCAna::prongEnergy2D
float prongEnergy2D[20]
Definition: NCAna_module.cc:401
simb::MCParticle::Py
double Py(const int i=0) const
Definition: MCParticle.h:230
ncs::NCAna::PlaneNcells
int PlaneNcells[192]
Definition: NCAna_module.cc:256
ncs::NCAna::trackcellCell
int trackcellCell[500]
Definition: NCAna_module.cc:288
ncs::NCAna::KtrackStopZ
float KtrackStopZ
Definition: NCAna_module.cc:276
ssth23
double ssth23
Definition: make_nue_pidtuning.C:51
ncs::NCAna::pi0Pz
float pi0Pz[20]
Definition: NCAna_module.cc:244
ncs::NCAna::prongTheta
float prongTheta[20]
Definition: NCAna_module.cc:318
simb::MCFlux::ftpx
double ftpx
Definition: MCFlux.h:79
ncs::NCAna::KtrackStartX
float KtrackStartX
Definition: NCAna_module.cc:271
stan::math::fabs
fvar< T > fabs(const fvar< T > &x)
Definition: fabs.hpp:15
art::SubRun::subRun
SubRunNumber_t subRun() const
Definition: SubRun.h:44
ncs::NCAna::trackcellADC
int trackcellADC[500]
Definition: NCAna_module.cc:285
ncs::NCAna::trackcellNcell
int trackcellNcell
Definition: NCAna_module.cc:280
ncs::NCAna::Nmcapd_24cells
int Nmcapd_24cells
Definition: NCAna_module.cc:149
nprongs
int nprongs
Definition: FillTrainingTrees.h:26
cheat::BackTracker::SliceToNeutrinoInteractions
std::vector< NeutrinoEffPur > SliceToNeutrinoInteractions(const std::vector< const rb::CellHit * > &sliceHits, const std::vector< const rb::CellHit * > &allHits, bool sortPur=false) const
Given a collection of hits (often a slice), returns vector of structures of neutrino interactions...
Definition: BackTracker_service.cc:1062
ncs::NCAna::SimpleOscProb
double SimpleOscProb(const simb::MCFlux &flux, const simb::MCNeutrino &nu) const
Definition: NCAna_module.cc:2467
cheat::BackTracker::ParticleNavigator
const sim::ParticleNavigator & ParticleNavigator() const
Get a reference to the ParticleNavigator.
Definition: BackTracker.h:744
ncs::NCAna::endSubRun
void endSubRun(const art::SubRun &sr)
Definition: NCAna_module.cc:846
ncs::NCAna::prongPhi
float prongPhi[20]
Definition: NCAna_module.cc:317
ncs::NCAna::Ndb2apd_24cells
int Ndb2apd_24cells
Definition: NCAna_module.cc:147
ncs::NCAna::shwStartZ
float shwStartZ[20]
Definition: NCAna_module.cc:438
ncs::NCAna::fKTrackLabel
std::string fKTrackLabel
label for slices
Definition: NCAna_module.cc:105
mf::LogInfo
MaybeLogger_< ELseverityLevel::ELsev_info, false > LogInfo
Definition: MessageLogger.h:214
rb::CellHit::Plane
unsigned short Plane() const
Definition: CellHit.h:39
ncs::NCAna::pi0E
float pi0E[20]
Definition: NCAna_module.cc:245
ncs::NCAna::KtrackStopY
float KtrackStopY
Definition: NCAna_module.cc:275
rb::CellHit::View
geo::View_t View() const
Definition: CellHit.h:41
simb::MCFlux::ftptype
int ftptype
Definition: MCFlux.h:82
ncs::NCAna::trackcellE
float trackcellE[500]
Definition: NCAna_module.cc:284
p
const char * p
Definition: xmltok.h:285
simb::MCNeutrino::Nu
const simb::MCParticle & Nu() const
Definition: MCNeutrino.h:146
ncs::NCAna::shwLID
float shwLID[20]
Definition: NCAna_module.cc:430
ncs::NCAna::shwcellX
float shwcellX[500]
Definition: NCAna_module.cc:448
art::Timestamp::timeHigh
constexpr std::uint32_t timeHigh() const
Definition: Timestamp.h:35
ncs::NCAna::prong1cellNcell
int prong1cellNcell
Definition: NCAna_module.cc:369
ncs::NCAna::slcellADC
int slcellADC[500]
Definition: NCAna_module.cc:304
ncs::NCAna::shwNCells
int shwNCells[20]
Definition: NCAna_module.cc:432
geo::kX
Vertical planes which measure X.
Definition: PlaneGeo.h:28
ncs::NCAna::prong1cellW
float prong1cellW[500]
Definition: NCAna_module.cc:379
simb::MCNeutrino::Pt
double Pt() const
transverse momentum of interaction, in GeV/c
Definition: MCNeutrino.cxx:74
ncs::NCAna::allprongcellZ
float allprongcellZ[500]
Definition: NCAna_module.cc:337
simb::MCParticle::Px
double Px(const int i=0) const
Definition: MCParticle.h:229
ncs::NCAna::sl_Nmcapd_28cells
int sl_Nmcapd_28cells
Definition: NCAna_module.cc:182
ncs::NCAna::shwIndex
int shwIndex[500]
Definition: NCAna_module.cc:447
calib::Calibrator::MakeRecoHit
rb::RecoHit MakeRecoHit(rb::CellHit const &cellhit, double w)
Definition: Calibrator_service.cc:575
ncs::NCAna::widthY
double widthY
Definition: NCAna_module.cc:134
simb::MCNeutrino::HitNuc
int HitNuc() const
Definition: MCNeutrino.h:152
ncs::NCAna::allprong2DcellX
float allprong2DcellX[500]
Definition: NCAna_module.cc:416
ncs::NCAna::sl_Nmcapd_24cells
int sl_Nmcapd_24cells
Definition: NCAna_module.cc:178
rb::Cluster
A collection of associated CellHits.
Definition: Cluster.h:47
ncs::NCAna::prong1cellCell
int prong1cellCell[500]
Definition: NCAna_module.cc:377
ncs::NCAna::sl_Fouttimehits_noisyapd
float sl_Fouttimehits_noisyapd
Definition: NCAna_module.cc:185
ncs::NCAna::shwStartY
float shwStartY[20]
Definition: NCAna_module.cc:437
ncs::NCAna::deltaspilltimensec
int deltaspilltimensec
Definition: NCAna_module.cc:135
ncs::NCAna::~NCAna
virtual ~NCAna()
Definition: NCAna_module.cc:483
abs
void abs(TH1 *hist)
Definition: wrong_sign_uncertainty.C:716
ncs::NCAna::trackcellPlane
int trackcellPlane[500]
Definition: NCAna_module.cc:287
ncs::NCAna::fProngLabel
std::string fProngLabel
label for kalman tracks
Definition: NCAna_module.cc:106
ncs::NCAna::shwStopY
float shwStopY[20]
Definition: NCAna_module.cc:440
ncs::NCAna::KtrackStopX
float KtrackStopX
Definition: NCAna_module.cc:274
ip
TString ip
Definition: loadincs.C:5
ncs::NCAna::spillPOT
double spillPOT
Definition: NCAna_module.cc:128
ncs::NCAna::pi0Vz
float pi0Vz[20]
Definition: NCAna_module.cc:248
ncs::NCAna::shwcellZ
float shwcellZ[500]
Definition: NCAna_module.cc:450
simb::kCCQE
charged current quasi-elastic
Definition: MCNeutrino.h:96
ncs::NCAna::prongStartX2D
float prongStartX2D[20]
Definition: NCAna_module.cc:405
PandAna.reco_validation.prod5_pid_validation.ratio
def ratio(spec1, spec2, nbins, pot, binrange=(0, 1))
Definition: prod5_pid_validation.py:554
ncs::NCAna::sl_Nouttimehits
int sl_Nouttimehits
Definition: NCAna_module.cc:183
ncs::NCAna::shwcellE
float shwcellE[500]
Definition: NCAna_module.cc:451
ncs::NCAna::allprong2DcellPlane
int allprong2DcellPlane[500]
Definition: NCAna_module.cc:421
ncs::NCAna::prong2DInXView
int prong2DInXView[20]
Definition: NCAna_module.cc:400
ncs::NCAna::slcellW
float slcellW[500]
Definition: NCAna_module.cc:308
ncs::NCAna::shwEnu
float shwEnu[20]
Definition: NCAna_module.cc:444
ncs::NCAna::shwStopZ
float shwStopZ[20]
Definition: NCAna_module.cc:441
ncs::NCAna::prong1Fmip
float prong1Fmip
Definition: NCAna_module.cc:370
ncs::NCAna::prongStopY
float prongStopY[20]
Definition: NCAna_module.cc:324
ncs::NCAna::fVertexLabel
std::string fVertexLabel
label for generator information
Definition: NCAna_module.cc:110
util::sqr
T sqr(T x)
More efficient square function than pow(x,2)
Definition: MathUtil.h:23
ncs::NCAna::Fouttimehits_noisyapd
float Fouttimehits_noisyapd
Definition: NCAna_module.cc:156
ncs::NCAna::shwEcells
float shwEcells[20]
Definition: NCAna_module.cc:443
ncs::NCAna::prongEnergyXView
float prongEnergyXView[20]
Definition: NCAna_module.cc:326
ncs::NCAna::Ndb3apd_28cells
int Ndb3apd_28cells
Definition: NCAna_module.cc:152
flux
Loaders::FluxType flux
Definition: plot_prediction_systs.C:48
slid::CompareByEnergy
bool CompareByEnergy(const slid::ShowerLID *a, const slid::ShowerLID *b)
Definition: ShowerLID.cxx:51
simb::MCParticle::NumberDaughters
int NumberDaughters() const
Definition: MCParticle.h:216
ncs::NCAna::prong1cellX
float prong1cellX[500]
Definition: NCAna_module.cc:371
geo::kY
Horizontal planes which measure Y.
Definition: PlaneGeo.h:29
ncs::NCAna::sl_widthX
double sl_widthX
Definition: NCAna_module.cc:169
MCFlux.h
object containing MC flux information
ncs::NCAna::pi0Vy
float pi0Vy[20]
Definition: NCAna_module.cc:247
ncs::NCAna::prongPhi2D
float prongPhi2D[20]
Definition: NCAna_module.cc:403
ncs::NCAna::allprongcellADC
int allprongcellADC[500]
Definition: NCAna_module.cc:339
simb::MCParticle::Daughter
int Daughter(const int i) const
Definition: MCParticle.cxx:112
rb::RecoHit
Calibrated quantities relying on position in the orthogonal view. To generate a rb::CellHit from a rb...
Definition: RecoHit.h:19
febshutoff_auto.start
start
Definition: febshutoff_auto.py:53
ncs::NCAna::prong1cellZ
float prong1cellZ[500]
Definition: NCAna_module.cc:373
ncs::NCAna::fCellHitLabel
std::string fCellHitLabel
label for vertex
Definition: NCAna_module.cc:111
ncs::NCAna::KtrackStartZ
float KtrackStartZ
Definition: NCAna_module.cc:273
ncs::NCAna::sl_Nintimehits
int sl_Nintimehits
Definition: NCAna_module.cc:184
ncs::NCAna::allprongcellE
float allprongcellE[500]
Definition: NCAna_module.cc:338
lem
PID
Definition: FillPIDs.h:14
ncs::NCAna::_otree
TTree * _otree
Definition: NCAna_module.cc:158
ncs::NCAna::prongStopZ
float prongStopZ[20]
Definition: NCAna_module.cc:325
ncs::NCAna::sl_Ndb3apd_24cells
int sl_Ndb3apd_24cells
Definition: NCAna_module.cc:177
rb::CellHit::Cell
unsigned short Cell() const
Definition: CellHit.h:40
ncs::NCAna::fTotalPOT
double fTotalPOT
Definition: NCAna_module.cc:93
ncs::NCAna::allprong2DcellADC
int allprong2DcellADC[500]
Definition: NCAna_module.cc:420
ncs::NCAna::_btree
TTree * _btree
Definition: NCAna_module.cc:95
rb::Cluster::MinTNS
double MinTNS() const
Definition: Cluster.cxx:482
simb::MCNeutrino::InteractionType
int InteractionType() const
Definition: MCNeutrino.h:150
ncs::NCAna::trackcellW
float trackcellW[500]
Definition: NCAna_module.cc:292
ncs::NCAna::KtrackStartY
float KtrackStartY
Definition: NCAna_module.cc:272
ncs::NCAna::shwTheta
float shwTheta[20]
Definition: NCAna_module.cc:442
ncs::NCAna::minus2sigma
float minus2sigma[68]
Definition: NCAna_module.cc:226
ncs::DEFINE_ART_MODULE
DEFINE_ART_MODULE(ROCKMRE)
ncs::NCAna::Ndb3apd_24cells
int Ndb3apd_24cells
Definition: NCAna_module.cc:148
check_time_usage.float
float
Definition: check_time_usage.py:7
simb::MCNeutrino::Lepton
const simb::MCParticle & Lepton() const
Definition: MCNeutrino.h:147
sim::ParticleNavigator::Primary
const sim::Particle * Primary(const int) const
Definition: ParticleNavigator.cxx:155
ncs::NCAna::shwLength
float shwLength[20]
Definition: NCAna_module.cc:435
simb::MCNeutrino::W
double W() const
Definition: MCNeutrino.h:154
ncs::NCAna::Nmcapd_28cells
int Nmcapd_28cells
Definition: NCAna_module.cc:153
ncs::NCAna::fweightLabel
std::string fweightLabel
Definition: NCAna_module.cc:117
nshowers
int nshowers
Definition: FillTrainingTrees.h:42
sumdata::SpillData::deltaspilltimensec
signed long long int deltaspilltimensec
Definition: SpillData.h:24
breakpoints::beginJob
void beginJob()
ncs::NCAna::prongNCells
int prongNCells[20]
Definition: NCAna_module.cc:328
art::PtrVector::push_back
void push_back(Ptr< U > const &p)
Definition: PtrVector.h:441
ncs::NCAna::allprongcellY
float allprongcellY[500]
Definition: NCAna_module.cc:336
CLHEP::L
static constexpr double L
Definition: SystemOfUnits.h:104
simb::MCNeutrino::Y
double Y() const
Definition: MCNeutrino.h:156
ncs::NCAna::allprongIndex
int allprongIndex[500]
Definition: NCAna_module.cc:334
ncs::NCAna::trackcellZ
float trackcellZ[500]
Definition: NCAna_module.cc:283
ncs::NCAna::prongNCellsYView
int prongNCellsYView[20]
Definition: NCAna_module.cc:330
ncs::NCAna::KtrackNcells_ratio
float KtrackNcells_ratio
Definition: NCAna_module.cc:266
ncs::NCAna::prongStartY2D
float prongStartY2D[20]
Definition: NCAna_module.cc:406
ncs::NCAna::prong1cellY
float prong1cellY[500]
Definition: NCAna_module.cc:372
ncs::NCAna::prongStopZ2D
float prongStopZ2D[20]
Definition: NCAna_module.cc:411
ncs::NCAna::prongNCellsXView
int prongNCellsXView[20]
Definition: NCAna_module.cc:329
make_syst_table_plots.h
h
Definition: make_syst_table_plots.py:240
fhicl::ParameterSet::get
T get(std::string const &key) const
Definition: ParameterSet.h:231
ncs::NCAna::prongStopY2D
float prongStopY2D[20]
Definition: NCAna_module.cc:410
ncs::NCAna::flidLabel
std::string flidLabel
label for rvp
Definition: NCAna_module.cc:114
sumdata::SpillData::widthx
double widthx
mm
Definition: SpillData.h:36
evt
int evt
Definition: FillTrainingTrees.h:22
ncs::NCAna::Ebalance
float Ebalance
Definition: NCAna_module.cc:313
ncs::NCAna::KtrackLength
float KtrackLength
Definition: NCAna_module.cc:263
ncs::NCAna::prongNCells2D
int prongNCells2D[20]
Definition: NCAna_module.cc:402
ncs::NCAna::slice_mintns
float slice_mintns
Definition: NCAna_module.cc:259
ncs::NCAna::Ndb2apd_28cells
int Ndb2apd_28cells
Definition: NCAna_module.cc:151
ncs::NCAna::sl_Ndb3apd_28cells
int sl_Ndb3apd_28cells
Definition: NCAna_module.cc:181
pot
#define pot
Definition: template_basic.C:49
ncs::NCAna::sl_Ndb2apd_24cells
int sl_Ndb2apd_24cells
Definition: NCAna_module.cc:176
ncs::NCAna::trueTheta
float trueTheta
Definition: NCAna_module.cc:198
sr
caf::StandardRecord * sr
Definition: FillTrainingTrees.h:18
ncs::NCAna::allprong2DcellNcell
int allprong2DcellNcell
Definition: NCAna_module.cc:414
simb::MCNeutrino::X
double X() const
Definition: MCNeutrino.h:155
ncs::NCAna::prongStopX
float prongStopX[20]
Definition: NCAna_module.cc:323
sumdata::SpillData::hornI
double hornI
kA
Definition: SpillData.h:27
ncs::NCAna::allprongcellX
float allprongcellX[500]
Definition: NCAna_module.cc:335
GetY
double GetY(int dcm=1, int feb=0, int pix=0)
Definition: PlotOnMon.C:119
simb::MCFlux::ftpz
double ftpz
Definition: MCFlux.h:81
ncs::NCAna::prongTheta2D
float prongTheta2D[20]
Definition: NCAna_module.cc:404
ncs::NCAna::slcellCell
int slcellCell[500]
Definition: NCAna_module.cc:306
ncs::NCAna::allprongcellNcell
int allprongcellNcell
Definition: NCAna_module.cc:333
Vertex.h
Vertex location in position and time.
ncs::NCAna::slcellY
float slcellY[500]
Definition: NCAna_module.cc:301
rb
Perform a "2 point" Hough transform on a collection of hits.
Definition: dEdxCalculator.h:20
ncs::NCAna::fInstLabel2D
std::string fInstLabel2D
instance label for prongs 3D
Definition: NCAna_module.cc:108
sim
Example routine for calculating the "ultimate e-m mother" of a particle in a simulated event...
Definition: FillTruth.h:16
sumdata::SpillData::posy
double posy
mm
Definition: SpillData.h:35
ncs::NCAna::beginJob
void beginJob()
Definition: NCAna_module.cc:514
ncs::NCAna::shwStartX
float shwStartX[20]
Definition: NCAna_module.cc:436
simb
This class describes a particle created in the detector Monte Carlo simulation.
Definition: TrackEva_module.cc:41
ncs::NCAna::shwNCellsYview
int shwNCellsYview[20]
Definition: NCAna_module.cc:434
rb::Cluster::MaxTNS
double MaxTNS() const
Definition: Cluster.cxx:528
ncs::NCAna::trackcellPECorr
float trackcellPECorr[500]
Definition: NCAna_module.cc:291
ncs::NCAna::allprongcellPlane
int allprongcellPlane[500]
Definition: NCAna_module.cc:340
ncs::NCAna::KtrackEnergy_ratio
float KtrackEnergy_ratio
Definition: NCAna_module.cc:267
ncs::NCAna::Ndb1apd_28cells
int Ndb1apd_28cells
Definition: NCAna_module.cc:150
om::cout
OStream cout
Definition: OStream.cxx:6
ncs::NCAna::prong1cellPECorr
float prong1cellPECorr[500]
Definition: NCAna_module.cc:378
ncs::NCAna::sl_spillPOT
double sl_spillPOT
Definition: NCAna_module.cc:165
ncs::NCAna::allprong2DIndex
int allprong2DIndex[500]
Definition: NCAna_module.cc:415
rb::Cluster::MinPlane
unsigned int MinPlane(geo::View_t view=geo::kXorY) const
Definition: Cluster.cxx:462
Event.h
ncs::NCAna::shwcellPECorr
float shwcellPECorr[500]
Definition: NCAna_module.cc:455
art::EventID::event
EventNumber_t event() const
Definition: EventID.h:116
ncs::NCAna::horncurrent
double horncurrent
Definition: NCAna_module.cc:130
rb::Cluster::Cell
art::Ptr< rb::CellHit > Cell(geo::View_t view, unsigned int viewIdx) const
Get the ith cell from view view.
Definition: Cluster.cxx:145
ncs::NCAna::fDataSpillLabel
std::string fDataSpillLabel
Definition: NCAna_module.cc:116
ncs::NCAna::allprong2DcellPECorr
float allprong2DcellPECorr[500]
Definition: NCAna_module.cc:423
ncs::NCAna::allprongcellW
float allprongcellW[500]
Definition: NCAna_module.cc:343
ncs::NCAna::prongLength2D
float prongLength2D[20]
Definition: NCAna_module.cc:408
simb::MCParticle::Vx
double Vx(const int i=0) const
Definition: MCParticle.h:220
ncs::NCAna::shwcellY
float shwcellY[500]
Definition: NCAna_module.cc:449
ncs::NCAna::pi0Py
float pi0Py[20]
Definition: NCAna_module.cc:243
ncs::NCAna::KtrackPhi
float KtrackPhi
Definition: NCAna_module.cc:270
ncs::NCAna::trackcellY
float trackcellY[500]
Definition: NCAna_module.cc:282
ncs::NCAna::prong1cellE
float prong1cellE[500]
Definition: NCAna_module.cc:374
art::TFileDirectory::make
T * make(ARGS...args) const
Definition: TFileDirectory.h:82
ncs::NCAna::analyze
void analyze(const art::Event &evt)
Definition: NCAna_module.cc:872
rb::RecoHit::GeV
float GeV() const
Definition: RecoHit.cxx:69
std::sin
T sin(T number)
Definition: d0nt_math.hpp:132
ncs::NCAna::MakeRecoHit
const rb::RecoHit MakeRecoHit(art::Ptr< rb::CellHit > const &chit, double const &w)
Definition: NCAna_module.cc:862
rawdata::RawDigit::ADC
int16_t ADC(uint32_t i) const
Definition: RawDigit.cxx:58
HoughResult.h
Data resulting from a Hough transform on the cell hit positions.
ncs::NCAna::shwcellADC
int shwcellADC[500]
Definition: NCAna_module.cc:452
art::Ptr::get
T const * get() const
Definition: Ptr.h:321
dir
TDirectory * dir
Definition: macro.C:5
simb::MCFlux::fntype
int fntype
Definition: MCFlux.h:51
art::DataViewImpl::getByLabel
bool getByLabel(std::string const &label, std::string const &productInstanceName, Handle< PROD > &result) const
Definition: DataViewImpl.h:344
ncs::NCAna::prongNplanes
int prongNplanes[20]
Definition: NCAna_module.cc:331
ncs::NCAna::KtrackTheta
float KtrackTheta
Definition: NCAna_module.cc:269
ncs::NCAna::prongStartZ2D
float prongStartZ2D[20]
Definition: NCAna_module.cc:407
sumdata::POTSum::goodspills
int goodspills
Definition: POTSum.h:31
ncs::NCAna::widthX
double widthX
Definition: NCAna_module.cc:133
hit
Definition: structs.h:12
ncs::NCAna::PlaneEnergy
float PlaneEnergy[192]
Definition: NCAna_module.cc:254
rb::Cluster::MeanTNS
double MeanTNS(rb::AveragingScheme scheme=kDefaultScheme) const
Definition: Cluster.cxx:554
ncs::NCAna::plus2sigma
float plus2sigma[68]
Definition: NCAna_module.cc:224
geom
void geom(int which=0)
Definition: geom.C:163
ncs::NCAna::slcellZ
float slcellZ[500]
Definition: NCAna_module.cc:302
simb::MCParticle::Pz
double Pz(const int i=0) const
Definition: MCParticle.h:231
ncs::NCAna::slice_meantns
float slice_meantns
Definition: NCAna_module.cc:258
rb::Cluster::ExtentPlane
unsigned int ExtentPlane(geo::View_t view=geo::kXorY) const
Definition: Cluster.h:250
ncs::NCAna::prongStartX
float prongStartX[20]
Definition: NCAna_module.cc:319
simb::MCParticle::Vz
double Vz(const int i=0) const
Definition: MCParticle.h:222
std::cos
T cos(T number)
Definition: d0nt_math.hpp:78
ncs::NCAna::deaddcm_location
std::string deaddcm_location[14]
Definition: NCAna_module.cc:145
art::Event::time
Timestamp time() const
Definition: Event.h:61
ncs::NCAna::shwcellPlane
int shwcellPlane[500]
Definition: NCAna_module.cc:453
ncs::NCAna::flemLabel
std::string flemLabel
Definition: NCAna_module.cc:115
ncs::NCAna::shwcellW
float shwcellW[500]
Definition: NCAna_module.cc:456
ncs::NCAna::prongLength
float prongLength[20]
Definition: NCAna_module.cc:322
ncs::NCAna::shwNCellsXview
int shwNCellsXview[20]
Definition: NCAna_module.cc:433
ncs::NCAna::fSliceLabel
std::string fSliceLabel
Definition: NCAna_module.cc:103
ncs::NCAna::prongStopX2D
float prongStopX2D[20]
Definition: NCAna_module.cc:409
simb::MCFlux::ftpy
double ftpy
Definition: MCFlux.h:80
ncs::NCAna::sl_Ndb2apd_28cells
int sl_Ndb2apd_28cells
Definition: NCAna_module.cc:180
ncs::NCAna::shw3D
int shw3D[20]
Definition: NCAna_module.cc:431
ncs::NCAna::reconfigure
void reconfigure(const fhicl::ParameterSet &p)
Definition: NCAna_module.cc:488
ncs::NCAna::prongStartZ
float prongStartZ[20]
Definition: NCAna_module.cc:321
ncs::NCAna::sl_deaddcm_location
std::string sl_deaddcm_location[14]
Definition: NCAna_module.cc:174
ncs::NCAna::sl_widthY
double sl_widthY
Definition: NCAna_module.cc:170
ncs::NCAna::Ebalance2D
float Ebalance2D
Definition: NCAna_module.cc:412
rb::Cluster::IsNoise
bool IsNoise() const
Is the noise flag set?
Definition: Cluster.h:163
ncs::NCAna::pi0Px
float pi0Px[20]
Definition: NCAna_module.cc:242
ncs::NCAna::shwStopX
float shwStopX[20]
Definition: NCAna_module.cc:439
ncs::NCAna::prongEnergyYView
float prongEnergyYView[20]
Definition: NCAna_module.cc:327
sumdata::POTSum::totgoodpot
double totgoodpot
normalized by 10^12 POT
Definition: POTSum.h:28
ncs::NCAna::frvpLabel
std::string frvpLabel
label for cell hits
Definition: NCAna_module.cc:113
ncs::NCAna::sl_posY
double sl_posY
Definition: NCAna_module.cc:168
art::Event::run
RunNumber_t run() const
Definition: Event.h:77
ncs::NCAna::prong1cellPlane
int prong1cellPlane[500]
Definition: NCAna_module.cc:376
ncs::NCAna::fMCFluxLabel
std::string fMCFluxLabel
instance label for prongs 2D
Definition: NCAna_module.cc:109
simb::MCNeutrino
Event generator information.
Definition: MCNeutrino.h:18
ncs::NCAna::KtrackEnergy
float KtrackEnergy
Definition: NCAna_module.cc:265
w
Float_t w
Definition: plot.C:20
simb::MCNeutrino::Mode
int Mode() const
Definition: MCNeutrino.h:149
sumdata::SpillData::posx
double posx
mm
Definition: SpillData.h:34
W
#define W(x)
ncs::NCAna::slice_maxtns
float slice_maxtns
Definition: NCAna_module.cc:260
ncs::NCAna::slcellPECorr
float slcellPECorr[500]
Definition: NCAna_module.cc:307
ncs::NCAna::sl_Ndb1apd_28cells
int sl_Ndb1apd_28cells
Definition: NCAna_module.cc:179
ncs::NCAna::pi0Vx
float pi0Vx[20]
Definition: NCAna_module.cc:246
ncs::NCAna::minus1sigma
float minus1sigma[68]
Definition: NCAna_module.cc:225
ncs::NCAna::slcellTNS
float slcellTNS[500]
Definition: NCAna_module.cc:309
sumdata::SpillData::widthy
double widthy
mm
Definition: SpillData.h:37
ncs::NCAna::prong1cellADC
int prong1cellADC[500]
Definition: NCAna_module.cc:375
ncs::NCAna::sl_posX
double sl_posX
Definition: NCAna_module.cc:167
sumdata::SpillData::spillpot
double spillpot
POT for spill normalized by 10^12.
Definition: SpillData.h:26
rvp
Definition: FillPIDs.h:15
art::Event::id
EventID id() const
Definition: Event.h:56
simb::MCParticle::Vy
double Vy(const int i=0) const
Definition: MCParticle.h:221
Geometry.h
Encapsulate the geometry of one entire detector (near, far, ndos)
ncs::NCAna::slcellPlane
int slcellPlane[500]
Definition: NCAna_module.cc:305
GetX
double GetX(int ndb=14, int db=1, int feb=0, int pix=0)
Definition: PlotOnMon.C:111
art::Handle::failedToGet
bool failedToGet() const
Definition: Handle.h:196
cet::maybe_ref
Definition: maybe_ref.h:35
ncs::NCAna::Ndb1apd_24cells
int Ndb1apd_24cells
Definition: NCAna_module.cc:146
ncs::NCAna::allprongcellCell
int allprongcellCell[500]
Definition: NCAna_module.cc:341
ncs::NCAna::plus1sigma
float plus1sigma[68]
Definition: NCAna_module.cc:223
ncs::NCAna::trackcellX
float trackcellX[500]
Definition: NCAna_module.cc:281
ncs::NCAna::sl_deltaspilltimensec
int sl_deltaspilltimensec
Definition: NCAna_module.cc:171
ncs::NCAna::allprong2DcellZ
float allprong2DcellZ[500]
Definition: NCAna_module.cc:418
string
enum BeamMode string
Generated on Wed Feb 24 2021 14:13:08 for NOvA by   doxygen 1.8.11