CVNEventDump_module.cc
Go to the documentation of this file.
1 ////////////////////////////////////////////////////////////////////////
2 // \file CVNEventDump_module.cc
3 // \brief Analyzer module for creating CVN PixelMap objects
4 // \author Dominick Rocco - rocco@physics.umn.edu
5 ////////////////////////////////////////////////////////////////////////
6 
7 // C/C++ includes
8 #include <iostream>
9 #include <sstream>
10 #include <vector>
11 
12 // ROOT includes
13 #include "TTree.h"
14 #include "TH2F.h"
15 
16 // Framework includes
23 #include "fhiclcpp/ParameterSet.h"
27 
28 
29 // NOvASoft includes
30 #include "MCCheater/BackTracker.h"
31 #include "RecoBase/Cluster.h"
32 #include "RecoBase/FilterList.h"
36 #include "Utilities/AssociationUtil.h"
37 
38 #include "CVN/func/AssignLabels.h"
39 #include "CVN/func/TrainingData.h"
44 
46 #include "CVN/func/ParticlesType.h"
47 
48 // Required for CAFCutter
51 
52 
53 namespace cvn {
54  class CVNEventDump : public art::EDAnalyzer {
55  public:
56  explicit CVNEventDump(fhicl::ParameterSet const& pset);
57  ~CVNEventDump();
58 
59  void analyze(const art::Event& evt);
60  void reconfigure(const fhicl::ParameterSet& pset);
61  void beginJob();
62  void endJob();
63 
64 
65 
66  private:
67 
72 
73  bool interactionTraining = true;
74  bool particlesTraining = false;
75 
78  TTree* fTrainTree;
79 
81 
82  /// Check rb::IsFiltered?
84  std::vector<std::string> fPreselectionLabels;
85 
86  /// Function to extract TH2 from PixelMap and write to TFile
87  void WriteMapTH2(const art::Event& evt, int slice, const PixelMap& pm);
88 
89  std::string fCAFLabel; ///< label for the process that made the standard records
90  bool fApplyCAFCuts; ///< should we apply the CAF-level cuts?
91  int fCutType; ///< what cuts to apply (see CAFCutter.h)
92 
93  recovalid::CAFCutter fCAFCutter; ///< the CAFCutter helper class
94  };
95 
96  //.......................................................................
98  : EDAnalyzer(pset)
99  {
100  this->reconfigure(pset);
101  }
102 
103  //......................................................................
105  {
106  }
107 
108  //......................................................................
110  {
111  fClusterLabel = pset.get<std::string> ("ClusterLabel");
112  fPixelMapInput = pset.get<std::string> ("PixelMapInput");
113  fChooseLabels = pset.get<std::string> ("ChooseLabels");
114  fWriteMapTH2 = pset.get<bool> ("WriteMapTH2");
115  fObeyPreselection = pset.get<bool> ("ObeyPreselection");
116  fPreselectionLabels = pset.get< std::vector<std::string> >("PreselectionLabels");
117  fCAFLabel = pset.get< std::string >("CAFLabel");
118  fApplyCAFCuts = pset.get< bool > ("ApplyCAFCuts");
119  fCutType = pset.get< int > ("CutType");
120  }
121 
122  //......................................................................
124  {
125 
127 
128  if(fChooseLabels == "particleType"){
129  particlesTraining = true;
130  interactionTraining = false;
131  }
132 
133  fTrainTree = tfs->make<TTree>("CVNTrainTree", "Training records");
135  fTrainTree->Branch("train", "cvn::TrainingData", &fTrain);
137  fTrainTree->Branch("train", "cvn::ParticlesTrainingData", &fParticlesTrain);
138 
139  }
140 
141  //......................................................................
143  {
144  }
145 
146  //......................................................................
148  {
149 
150  // get the slices
152  evt.getByLabel(fClusterLabel, clusterHandle);
153  const std::vector<rb::Cluster> & clusters = *clusterHandle;
155  for(unsigned int i = 0; i < clusterHandle->size(); ++i){
156  slices.push_back(art::Ptr<rb::Cluster>(clusterHandle, i));
157  }
158 
159  art::FindManyP<PixelMap> fmPixelMap(clusterHandle, evt, fPixelMapInput);
160 
161  // Get the Standard Records
162  art::FindManyP<caf::StandardRecord> recordFMP(slices, evt, fCAFLabel);
163 
164  std::vector<int> bestNuId;
165  std::vector<cheat::NeutrinoWithIndex> nus;
166  std::vector< std::vector<cheat::NeutrinoEffPur> > sEffPur;
167  if(fBT->HaveTruthInfo()){
168  nus = fBT->allMCTruth();
169  sEffPur = fBT->SlicesToMCTruthsTable(slices);
170  bestNuId = fBT->SliceToOrderedNuIdsByEnergy(nus, sEffPur);
171  }
172 
173  for(size_t iClust = 0; iClust < clusters.size(); ++iClust){
174 
175  if(fObeyPreselection && rb::IsFiltered(evt, clusterHandle, iClust,
176  fPreselectionLabels)) continue;
177  // get record associated with the slice
178  std::vector< art::Ptr<caf::StandardRecord> > records;
179  bool pass = true;
180  // Apply the CAF-level cuts
181  if(fApplyCAFCuts && recordFMP.isValid()){
182  records = recordFMP.at(iClust);
183  if(records.size() > 0) pass = fCAFCutter.passCuts(fCutType, records[0].get());
184  }
185 
186  if(!pass) continue;
187 
188  const rb::Cluster& cluster = clusters[iClust];
189  if(!fmPixelMap.isValid()) continue;
190 
191 
192  const std::vector< art::Ptr<PixelMap> > pixelMaps = fmPixelMap.at(iClust);
193 
194  if(pixelMaps.empty()) continue;
195  InteractionType interaction = kOther;
196  FinalStateType finalstate = kOtherFS;
197  FinalStateProngsType finalstateprongs = kOtherFSp;
198  ParentParticleType parentparticle = kUnknownPP;
199 
200  ParticlesType particles_type = kOther_PT;
201 
202  float nuEnergy = 0;
203  float lepEnergy = 0;
204  double vtxx = std::numeric_limits<double>::lowest();
205  double vtxy = std::numeric_limits<double>::lowest();
206  double vtxz = std::numeric_limits<double>::lowest();
208  if(fBT->HaveTruthInfo()){
209  if(bestNuId[iClust] != -1)
210  truth = sEffPur[iClust][bestNuId[iClust]].neutrinoInt;
211 
212  interaction = SliceClassify(truth, cluster, nuEnergy, lepEnergy,
213  vtxx, vtxy, vtxz);
214  particles_type = ParticlesSliceClassify(truth, cluster, nuEnergy, lepEnergy,
215  vtxx, vtxy, vtxz);
216  finalstate = GetFinalStateType(truth);
217  finalstateprongs = GetFinalStateProngsType(truth);
218  parentparticle = GetParentParticleType(truth, evt);
219 
220  }
221  else particles_type = kCosmic_PT; // ! Careful This will only work for cosmic data !
222 
224  TrainingData train(interaction, parentparticle, finalstate, finalstateprongs,
225  particles_type, nuEnergy, lepEnergy,
226  vtxx, vtxy, vtxz, *pixelMaps[0]);
227  if(fWriteMapTH2) WriteMapTH2(evt, (int)iClust, train.fPMap);
228  fTrain = &train;
229  fTrainTree->Fill();
230  }
231 
232  if(particlesTraining){
233  ParticlesTrainingData train(particles_type, nuEnergy, lepEnergy,
234  vtxx, vtxy, vtxz, *pixelMaps[0]);
235  if(fWriteMapTH2) WriteMapTH2(evt, (int)iClust, train.fPMap);
237  fTrainTree->Fill();
238  }
239  } // end loop over slices
240 
241  } // analyze
242 
243  //----------------------------------------------------------------------
244  void CVNEventDump::WriteMapTH2(const art::Event& evt, int slice, const PixelMap& pm)
245  {
246 
247  TString mapid = "_r"+std::to_string(evt.run())+"_s"+std::to_string(evt.subRun())+"_e"
248  +std::to_string(evt.event())+"_sl"+std::to_string(slice);
249  TString data = "PixelMap";
250  TString label = "TruthPixelMap";
251  TString object= "ObjectPixelMap";
252  TString time = "TimePixelMap";
253  TString nuId = "NuIdPixelMap";
254 
255  TH2F* hist = pm.ToTH2( data+mapid );
256  TH2F* histL = pm.ToLabTH2( label+mapid);
257  TH2F* histO = pm.ToObjTH2( object+mapid);
258  TH2F* histX = pm.SingleViewToTH2(0, data+"ViewX"+mapid);
259  TH2F* histY = pm.SingleViewToTH2(1, data+"ViewY"+mapid);
260  TH2F* histXL = pm.SingleViewToLabTH2(0, label+"ViewX"+mapid, false);
261  TH2F* histYL = pm.SingleViewToLabTH2(1, label+"ViewY"+mapid, false);
262  TH2F* histXT = pm.SingleViewToLabTH2(0, time+"ViewX"+mapid, true);
263  TH2F* histYT = pm.SingleViewToLabTH2(1, time+"ViewY"+mapid, true);
264  TH2F* histXO = pm.SingleViewToObjTH2(0, object+"ViewX"+mapid, false);
265  TH2F* histYO = pm.SingleViewToObjTH2(1, object+"ViewY"+mapid, false);
266  TH2F* histXN = pm.SingleViewToObjTH2(0, nuId+"ViewX"+mapid, true);
267  TH2F* histYN = pm.SingleViewToObjTH2(1, nuId+"ViewY"+mapid, true);
268 
270 
271  TH2F* histWrite = tfs->make<TH2F>(*hist);
272  histWrite->Write();
273  TH2F* histWriteL = tfs->make<TH2F>(*histL);
274  histWriteL->GetZaxis()->SetRangeUser(0,10);
275  histWriteL->Write();
276  TH2F* histWriteO = tfs->make<TH2F>(*histO);
277  //histWriteL->GetZaxis()->SetRangeUser(0,10);
278  histWriteO->Write();
279  TH2F* histWriteX = tfs->make<TH2F>(*histX);
280  histWriteX->Write();
281  TH2F* histWriteY = tfs->make<TH2F>(*histY);
282  histWriteY->Write();
283  TH2F* histWriteXL = tfs->make<TH2F>(*histXL);
284  histWriteXL->Write();
285  TH2F* histWriteYL = tfs->make<TH2F>(*histYL);
286  histWriteYL->Write();
287  TH2F* histWriteXO = tfs->make<TH2F>(*histXO);
288  histWriteXO->Write();
289  TH2F* histWriteYO = tfs->make<TH2F>(*histYO);
290  histWriteYO->Write();
291  TH2F* histWriteXT = tfs->make<TH2F>(*histXT);
292  histWriteXT->Write();
293  TH2F* histWriteYT = tfs->make<TH2F>(*histYT);
294  histWriteYT->Write();
295  TH2F* histWriteXN = tfs->make<TH2F>(*histXN);
296  histWriteXN->Write();
297  TH2F* histWriteYN = tfs->make<TH2F>(*histYN);
298  histWriteYN->Write();
299 
300  delete hist;
301  delete histWrite;
302  delete histL;
303  delete histWriteL;
304  delete histO;
305  delete histXO;
306  delete histYO;
307  delete histWriteO;
308  delete histWriteXO;
309  delete histWriteYO;
310  delete histXN;
311  delete histYN;
312  delete histWriteXN;
313  delete histWriteYN;
314  delete histX;
315  delete histWriteX;
316  delete histY;
317  delete histWriteY;
318  delete histXL;
319  delete histWriteXL;
320  delete histWriteXT;
321  delete histYL;
322  delete histWriteYL;
323  delete histWriteYT;
324 
325  } // WriteMapTH2
326 
328 } // end namespace cvn
329 ////////////////////////////////////////////////////////////////////////
SubRunNumber_t subRun() const
Definition: Event.h:72
back track the reconstruction to the simulation
InteractionType SliceClassify(const art::Ptr< simb::MCTruth > truth, const rb::Cluster &slice, float &nuEnergy, float &lepEnergy, double &vtxx, double &vtxy, double &vtxz)
The TrainingData objects contains a PixelMap and the output class type, and any other bit that goes i...
FinalStateType GetFinalStateType(const art::Ptr< simb::MCTruth > truth)
TH2F * SingleViewToTH2(const unsigned int &view, TString name) const
Definition: PixelMap.cxx:312
art::ServiceHandle< cheat::BackTracker > fBT
std::vector< NeutrinoWithIndex > allMCTruth() const
PixelMap fPMap
PixelMap for the event.
void WriteMapTH2(const art::Event &evt, int slice, const PixelMap &pm)
Function to extract TH2 from PixelMap and write to TFile.
A collection of associated CellHits.
Definition: Cluster.h:47
enum cvn::Interaction InteractionType
bool fObeyPreselection
Check rb::IsFiltered?
float vtxz
Definition: NusVarsTemp.cxx:32
ParticlesType ParticlesSliceClassify(const art::Ptr< simb::MCTruth > truth, const rb::Cluster &slice, float &nuEnergy, float &lepEnergy, double &vtxx, double &vtxy, double &vtxz)
DEFINE_ART_MODULE(TestTMapFile)
ParticlesTrainingData * fParticlesTrain
Defines an enumeration for prong classification.
Particle class.
std::vector< std::vector< cheat::NeutrinoEffPur > > SlicesToMCTruthsTable(const std::vector< const rb::Cluster * > &sliceList) const
Given ALL the slices in an event, including the noise slice, returns a vector of vector of structures...
Something else. Tau? Hopefully we don&#39;t use this.
const char * label
PixelMap fPMap
PixelMap for the event.
Definition: TrainingData.h:69
TH2F * ToLabTH2(TString name) const
Definition: PixelMap.cxx:258
const XML_Char const XML_Char * data
Definition: expat.h:268
std::vector< int > SliceToOrderedNuIdsByEnergy(const std::vector< cheat::NeutrinoWithIndex > &nusWithIdx, const std::vector< std::vector< cheat::NeutrinoEffPur >> &slTruthTable) const
TH2F * ToTH2(TString name) const
Definition: PixelMap.cxx:236
bool passCuts(int cut, const caf::StandardRecord *sr)
Definition: CAFCutter.cxx:37
ParentParticleType GetParentParticleType(const art::Ptr< simb::MCTruth > truth, const art::Event &evt)
void push_back(Ptr< U > const &p)
Definition: PtrVector.h:441
bool fApplyCAFCuts
should we apply the CAF-level cuts?
Whatever is left.
Something else. Tau? Hopefully we don&#39;t use this.
T get(std::string const &key) const
Definition: ParameterSet.h:231
FinalStateProngsType GetFinalStateProngsType(const art::Ptr< simb::MCTruth > truth)
bool IsFiltered(const art::Event &evt, art::Ptr< T > x, const std::vector< std::string > &labels)
Is this Ptr marked "filtered out"?
Definition: FilterList.h:96
EventNumber_t event() const
Definition: Event.h:67
TH2F * ToObjTH2(TString name) const
Definition: PixelMap.cxx:285
EDAnalyzer(Table< Config > const &config)
Definition: EDAnalyzer.h:100
enum cvn::FinalState FinalStateType
CVNEventDump(fhicl::ParameterSet const &pset)
TH2F * SingleViewToObjTH2(const unsigned int &view, TString name, bool nuID) const
Definition: PixelMap.cxx:359
std::vector< std::string > fPreselectionLabels
enum cvn::Particles ParticlesType
void reconfigure(const fhicl::ParameterSet &pset)
Something else. Tau? Hopefully we don&#39;t use this.
Catch all term, shouldn&#39;t be used.
::xsd::cxx::tree::string< char, simple_type > string
Definition: Database.h:154
int fCutType
what cuts to apply (see CAFCutter.h)
T * make(ARGS...args) const
recovalid::CAFCutter fCAFCutter
the CAFCutter helper class
float vtxy
Definition: NusVarsTemp.cxx:31
enum cvn::ParentParticle ParentParticleType
void train()
Definition: train.C:24
bool getByLabel(std::string const &label, std::string const &productInstanceName, Handle< PROD > &result) const
Definition: DataViewImpl.h:344
The TrainingData objects contains a PixelMap and the output class type, and any other bit that goes i...
Definition: TrainingData.h:23
PixelMap, basic input to CVN neural net.
Definition: PixelMap.h:23
bool HaveTruthInfo() const
Is this a file with truth info in? (Is BackTracker going to be any use to you?)
Definition: BackTracker.h:133
enum cvn::FinalStateProngs FinalStateProngsType
float vtxx
Definition: NusVarsTemp.cxx:30
TH2F * SingleViewToLabTH2(const unsigned int &view, TString name, bool time) const
Definition: PixelMap.cxx:331
Helper class for Reco Validation modules.
Definition: CAFCutter.h:58
std::string to_string(ModuleType mt)
Definition: ModuleType.h:32
void analyze(const art::Event &evt)
RunNumber_t run() const
Definition: Event.h:77
std::string fCAFLabel
label for the process that made the standard records