CVNAddTrainingData_module.cc
Go to the documentation of this file.
1 ////////////////////////////////////////////////////////////////////////
2 // \file CVNAddTrainingData_module.cc
3 // \brief Producer module for adding CVN TrainingData objects
4 // \author Grant Nikseresht - gnikseresht@hawk.iit.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
21 #include "fhiclcpp/ParameterSet.h"
25 
26 // NOvASoft includesOkay
27 #include "Geometry/Geometry.h"
28 #include "Calibrator/Calibrator.h"
29 #include "MCCheater/BackTracker.h"
30 #include "ShowerLID/ShowerLID.h"
31 #include "RecoBase/CellHit.h"
32 #include "RecoBase/Cluster.h"
33 #include "RecoBase/FilterList.h"
34 #include "RecoBase/Prong.h"
35 #include "RecoBase/Shower.h"
36 #include "RecoBase/Vertex.h"
40 #include "Utilities/AssociationUtil.h"
42 
43 #include "CVN/func/AssignLabels.h"
44 #include "CVN/func/TrainingData.h"
45 #include "CVN/func/ProngType.h"
52 
54 #include "CVN/func/ParticlesType.h"
55 
56 // Required for CAFCutter
59 
60 namespace cvn {
61 
63  public:
64  explicit CVNAddTrainingData(fhicl::ParameterSet const& pset);
66 
67  void produce(art::Event& evt);
68  void reconfigure(const fhicl::ParameterSet& pset);
69  void beginJob();
70  void endJob();
71 
72 
73 
74  private:
88  bool fTDProngs;
89 
90 
91 
92  bool interactionTraining = true;
93  bool particlesTraining = false;
94 
97 
99 
100  /// Check rb::IsFiltered?
102  std::vector<std::string> fPreselectionLabels;
103 
104  };
105 
107  {
108  this->reconfigure(pset);
109  produces< std::vector<cvn::TrainingData> >(fClusterTDLabel);
110  produces< std::vector<cvn::ProngTrainingData> >(fProngTDLabel);
111  produces< std::vector<cvn::TrainingData> >(fPMapTDLabel);
112  produces< std::vector<cvn::ProngTrainingData> >(fProngPMapTDLabel);
113  produces< art::Assns<cvn::TrainingData, rb::Cluster> >();
114  produces< art::Assns<cvn::ProngTrainingData, rb::Prong> >();
115  produces< art::Assns<cvn::TrainingData, cvn::PixelMap> >();
116  produces< art::Assns<cvn::ProngTrainingData, cvn::PixelMap> >();
117  }
118  //.......................................................................
120  {
121  fClusterLabel = pset.get<std::string> ("ClusterLabel");
122  fProngModLabel = pset.get<std::string> ("ProngModLabel");
123  fProng3DLabel = pset.get<std::string> ("Prong3DLabel");
124  fProng2DLabel = pset.get<std::string> ("Prong2DLabel");
125  fClusterTDLabel = pset.get<std::string> ("ClusterTDLabel");
126  fProngTDLabel = pset.get<std::string> ("ProngTDLabel");
127  fProngPMapTDLabel = pset.get<std::string> ("ProngPMapTDLabel");
128  fPMapTDLabel = pset.get<std::string> ("PMapTDLabel");
129  fPixelMapInput = pset.get<std::string> ("PixelMapInput");
130  fShowerLabel = pset.get<std::string> ("ShowerLabel");
131  fShowerLIDLabel = pset.get<std::string> ("ShowerLIDLabel");
132  fVertexLabel = pset.get<std::string> ("VertexLabel");
133  fChooseLabels = pset.get<std::string> ("ChooseLabels");
134  fObeyPreselection = pset.get<bool> ("ObeyPreselection");
135  fTDProngs = pset.get<bool> ("UseTDProngs");
136  fPreselectionLabels = pset.get< std::vector<std::string> >("PreselectionLabels");
137  }
138 
139  //......................................................................
141  {
142  //======================================================================
143  // Clean up any memory allocated by your module
144  //======================================================================
145  }
146 
147  //......................................................................
149  {
150  std::cout << "Beginning job" << '\n';
151  }
152 
153  //......................................................................
155  {
156  }
157 
158  //......................................................................
160  {
161  //std::cout << "Producing an art event in CVNAddTrainingData" << std::endl;
162 
163  // This section fills the labels needed for the TrainingData object
164  //----- Derived from CVNEventDump_module -------------//
165  // get the slices
166 
168 
170  evt.getByLabel(fClusterLabel, clusterHandle);
171  const std::vector<rb::Cluster> & clusters = *clusterHandle;
173  for(unsigned int i = 0; i < clusterHandle->size(); ++i){
174  slices.push_back(art::Ptr<rb::Cluster>(clusterHandle, i));
175  }
176 
177  std::unique_ptr<std::vector<TrainingData>>
178  trainData(new std::vector<TrainingData>);
179  std::unique_ptr<std::vector<ProngTrainingData>>
180  prongtrainData(new std::vector<ProngTrainingData>);
181  std::unique_ptr<std::vector<TrainingData>>
182  pmaptrainData(new std::vector<TrainingData>);
183  std::unique_ptr<std::vector<ProngTrainingData>>
184  pmapprongtrainData(new std::vector<ProngTrainingData>);
185  std::unique_ptr< art::Assns<TrainingData, rb::Cluster> >
186  slice_assoc(new art::Assns<TrainingData, rb::Cluster>);
187  std::unique_ptr< art::Assns<ProngTrainingData, rb::Prong> >
189  std::unique_ptr< art::Assns<TrainingData, PixelMap> >
190  pmap_assoc(new art::Assns<TrainingData, PixelMap>);
191  std::unique_ptr< art::Assns<ProngTrainingData, PixelMap> >
192  prongpmap_assoc(new art::Assns<ProngTrainingData, PixelMap>);
193  art::FindManyP<PixelMap> fmPixelMap(clusterHandle, evt, fPixelMapInput);
194 
195  std::vector<int> bestNuId;
196  std::vector<cheat::NeutrinoWithIndex> nus;
197  std::vector< std::vector<cheat::NeutrinoEffPur> > sEffPur;
198  //std::cout << "Has truth info: " << fBT->HaveTruthInfo() << std::endl;
199  if(fBT->HaveTruthInfo()){
200  nus = fBT->allMCTruth();
201  sEffPur = fBT->SlicesToMCTruthsTable(slices);
202  bestNuId = fBT->SliceToOrderedNuIdsByEnergy(nus, sEffPur);
203  }
204 
205  // for loop over slices
206  for(size_t iClust = 0; iClust < clusters.size(); ++iClust){
207 
208  if(fObeyPreselection && rb::IsFiltered(evt, clusterHandle, iClust,
209  fPreselectionLabels)) continue;
210 
211  std::vector< art::Ptr<caf::StandardRecord> > records;
212  bool pass = true;
213 
214  if(!pass) continue;
215 
216  const rb::Cluster& cluster = clusters[iClust];
217  if(!fmPixelMap.isValid()) continue;
218 
219  const std::vector< art::Ptr<PixelMap> > pixelMaps = fmPixelMap.at(iClust);
220 
221  if(pixelMaps.empty()) continue;
222 
223  InteractionType interaction = kOther;
224  FinalStateType finalstate = kOtherFS;
225  FinalStateProngsType finalstateprongs = kOtherFSp;
226  ParentParticleType parentparticle = kUnknownPP;
227 
228  ParticlesType particles_type = kOther_PT;
229 
230  float nuEnergy = 0;
231  float lepEnergy = 0;
232  double vtxx = std::numeric_limits<double>::lowest();
233  double vtxy = std::numeric_limits<double>::lowest();
234  double vtxz = std::numeric_limits<double>::lowest();
235 
237  if(fBT->HaveTruthInfo()){
238  if(bestNuId[iClust] != -1) {
239  truth = sEffPur[iClust][bestNuId[iClust]].neutrinoInt;
240 
241  interaction = SliceClassify(truth, cluster, nuEnergy, lepEnergy,
242  vtxx, vtxy, vtxz);
243  particles_type = ParticlesSliceClassify(truth, cluster, nuEnergy,
244  lepEnergy, vtxx, vtxy, vtxz);
245  finalstate = GetFinalStateType(truth);
246  finalstateprongs = GetFinalStateProngsType(truth);
247  parentparticle = GetParentParticleType(truth, evt);
248  }
249  }
250  else particles_type = kCosmic_PT;
251 
252  //----- End CVNEventDump_module duplication -------//
253  // Now the TrainingData object can be filled with volues
254  // set in the previous section
255 
256  // std::cout << "Interaction: " << interaction << '\n';
257  // std::cout << "ParentParticle: " << parentparticle << '\n';
258  // std::cout << "FinalState: " << finalstate << '\n';
259  // std::cout << "FinalStateProngs: " << finalstateprongs << '\n';
260  // std::cout << "ParticlesType: " << particles_type << '\n';
261  // std::cout << "nuEnergy: " << nuEnergy << '\n';
262  // std::cout << "lepEnergy: " << lepEnergy << '\n';
263  // std::cout << "vtxx: " << vtxx << '\n';
264  // std::cout << "vtxy: " << vtxy << '\n';
265  // std::cout << "vtxz: " << vtxz << '\n';
266 
267  TrainingData train(interaction, parentparticle, finalstate,
268  finalstateprongs, particles_type, nuEnergy,
269  lepEnergy, vtxx, vtxy, vtxz, *pixelMaps[0]);
270  trainData->push_back(train);
271  pmaptrainData->push_back(train);
272  util::CreateAssn(*this, evt, *(trainData.get()), slices[iClust], *(slice_assoc.get()), UINT_MAX, fClusterTDLabel);
273  util::CreateAssn(*this, evt, *(pmaptrainData.get()), pixelMaps[0], *(pmap_assoc.get()), UINT_MAX, fPMapTDLabel);
274  if(fTDProngs){
275 
276  art::FindManyP<rb::Vertex> fmv(clusterHandle, evt, fVertexLabel);
277  art::FindManyP<rb::Prong> fmProng3D(clusterHandle, evt,
279  art::FindManyP<rb::Prong> fmProng2D(clusterHandle, evt,
281 
282  std::vector<art::Ptr<rb::Prong>> prongs3D;
283  std::vector<art::Ptr<rb::Prong>> prongs2D;
284  std::vector< art::Ptr<rb::Vertex> > vert;
285 
286  if( fmProng2D.isValid() ) prongs2D = fmProng2D.at(iClust);
287  if( fmProng3D.isValid() ) prongs3D = fmProng3D.at(iClust);
288  if(fmv.isValid() ) vert = fmv.at(iClust);
289  else continue;
290 
291  art::FindManyP<PixelMap> fmPixelMap3D(prongs3D, evt, fPixelMapInput);
292  art::FindManyP<PixelMap> fmPixelMap2D(prongs2D, evt, fPixelMapInput);
293  art::FindOneP<rb::Shower> foSh(prongs3D, evt, fShowerLabel);
294 
295  TVector3 v;
296  if(!vert.empty() ) v = vert[0]->GetXYZ();
297  else v.SetXYZ(-999.,-999.,-999.);
298 
299  double length;
300  TVector3 prongdir;
301  TVector3 beamdir = geom->NuMIBeamDirection();
302  float angle;
303 
304 
305  for( unsigned int iProng = 0; iProng < prongs3D.size(); ++iProng ){
306 
307  length = -5.0;
308  angle = -5.0;
309 
310  if(!fmPixelMap3D.isValid()) continue;
311 
312  const std::vector< art::Ptr<PixelMap> > pixelMaps3D = fmPixelMap3D.at(iProng);
313  if(pixelMaps3D.empty()) continue;
314 
315  double gap = -5.0;
316  if(foSh.isValid()){
317  cet::maybe_ref< art::Ptr<rb::Shower> const > roSh(foSh.at(iProng));
318  art::Ptr<rb::Shower> shower = roSh.ref();
319 
320  std::vector< art::Ptr<rb::Shower> > showers;
321  showers.push_back(shower);
322 
323  art::FindOneP<slid::ShowerLID> foShLID(showers, evt,fShowerLIDLabel);
324  if (foShLID.isValid()){
325  cet::maybe_ref< art::Ptr<slid::ShowerLID> const > roShLID(foShLID.at(0));
326  art::Ptr<slid::ShowerLID> showerLID = roShLID.ref();
327  gap = showerLID->Gap();
328  }
329  }
330 
331  ProngTrainingData p_train;
332  ProngType ptype3D, ptypeX, ptypeY;
333  double purity3D, purityX, purityY, recE;
334  unsigned int ncellX, ncellY;
335  bool isprimary;
336 
337  ptype3D = ProngClassify(*prongs3D[iProng], &ptype3D, &ptypeX, &ptypeY,
338  &isprimary, &purity3D, &purityX, &purityY,
339  &recE, &ncellX, &ncellY);
340 
341  p_train = ProngTrainingData(ptype3D, ptypeX, ptypeY, isprimary,
342  purity3D, purityX, purityY, recE,
343  ncellX, ncellY, v.X(), v.Y(), v.Z(),
344  length, (double)angle, gap, *pixelMaps3D[0]);
345 
346  prongtrainData->push_back(p_train);
347  pmapprongtrainData->push_back(p_train);
348  util::CreateAssn(*this, evt, *(prongtrainData.get()),
349  prongs3D[iProng], *(prong_assoc.get()), UINT_MAX, fProngTDLabel);
350  util::CreateAssn(*this, evt, *(pmapprongtrainData.get()),
351  pixelMaps3D[0], *(prongpmap_assoc.get()), UINT_MAX, fProngPMapTDLabel);
352 
353  }
354  for( unsigned int iProng2 = 0; iProng2 < prongs2D.size(); ++iProng2 ){
355 
356  if(!fmPixelMap2D.isValid()) continue;
357 
358  const std::vector< art::Ptr<PixelMap> > pixelMaps2D = fmPixelMap2D.at(iProng2);
359  if(pixelMaps2D.empty()) continue;
360 
361  ProngTrainingData p2_train;
362  ProngType ptype2D, ptypeX2D, ptypeY2D;
363  double purity2D, purityX2D, purityY2D, recE;
364  unsigned int ncellX, ncellY;
365  bool isprimary;
366 
367  ptype2D = ProngClassify(*prongs2D[iProng2], &ptype2D, &ptypeX2D, &ptypeY2D,
368  &isprimary, &purity2D, &purityX2D, &purityY2D,
369  &recE, &ncellX, &ncellY);
370 
371  p2_train = ProngTrainingData(ptype2D, ptypeX2D, ptypeY2D, isprimary,
372  purity2D, purityX2D, purityY2D, recE,
373  ncellX, ncellY, v.X(), v.Y(), v.Z(),
374  0.0, 0.0, 0.0, *pixelMaps2D[0]);
375 
376  prongtrainData->push_back(p2_train);
377  pmapprongtrainData->push_back(p2_train);
378  util::CreateAssn(*this, evt, *(prongtrainData.get()),
379  prongs2D[iProng2], *(prong_assoc.get()), UINT_MAX, fProngTDLabel);
380  util::CreateAssn(*this, evt, *(pmapprongtrainData.get()),
381  pixelMaps2D[0], *(prongpmap_assoc.get()), UINT_MAX, fProngPMapTDLabel);
382 
383  }
384  } //fTDProngs
385  } // end for loop over slices
386  evt.put(std::move(trainData), fClusterTDLabel);
387  evt.put(std::move(prongtrainData), fProngTDLabel);
388  evt.put(std::move(pmaptrainData), fPMapTDLabel);
389  evt.put(std::move(pmapprongtrainData), fProngPMapTDLabel);
390  evt.put(std::move(slice_assoc));
391  evt.put(std::move(prong_assoc));
392  evt.put(std::move(pmap_assoc));
393  evt.put(std::move(prongpmap_assoc));
394  //std::cout << "Training data object created." << std::endl;
395  //----------------------------------------------------------------------
396  }
397 
399 } // end namespace cvn
400 ////////////////////////////////////////////////////////////////////////
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)
static bool CreateAssn(art::EDProducer const &prod, art::Event &evt, std::vector< T > &a, art::Ptr< U > b, art::Assns< T, U > &assn, size_t indx=UINT_MAX, std::string const &instance=std::string())
Create a 1 to 1 association between a new product and one already in the event.
The TrainingData objects contains a PixelMap and the output class type, and any other bit that goes i...
enum cvn::PType ProngType
Double_t angle
Definition: plot.C:86
FinalStateType GetFinalStateType(const art::Ptr< simb::MCTruth > truth)
std::vector< NeutrinoWithIndex > allMCTruth() const
float vtxx
A collection of associated CellHits.
Definition: Cluster.h:47
enum cvn::Interaction InteractionType
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)
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.
ProductID put(std::unique_ptr< PROD > &&product)
Definition: Event.h:102
std::vector< int > SliceToOrderedNuIdsByEnergy(const std::vector< cheat::NeutrinoWithIndex > &nusWithIdx, const std::vector< std::vector< cheat::NeutrinoEffPur >> &slTruthTable) const
void reconfigure(const fhicl::ParameterSet &pset)
length
Definition: demo0.py:21
ParentParticleType GetParentParticleType(const art::Ptr< simb::MCTruth > truth, const art::Event &evt)
void push_back(Ptr< U > const &p)
Definition: PtrVector.h:441
Whatever is left.
Something else. Tau? Hopefully we don&#39;t use this.
T get(std::string const &key) const
Definition: ParameterSet.h:231
int evt
float vtxy
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
enum cvn::FinalState FinalStateType
Vertex location in position and time.
enum cvn::Particles ParticlesType
The TrainingData objects contains a PixelMap and the output class type, and any other bit that goes i...
OStream cout
Definition: OStream.cxx:6
TVector3 NuMIBeamDirection() const
Direction of neutrinos from the NuMI beam (unit vector)
Something else. Tau? Hopefully we don&#39;t use this.
Catch all term, shouldn&#39;t be used.
float vtxz
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
ProngType ProngClassify(const rb::Prong &prong, ProngType *pType3D, ProngType *pTypeX, ProngType *pTypeY, bool *isprimary, double *purity3D, double *purityX, double *purityY, double *recE, unsigned int *ncellX, unsigned int *ncellY)
CVNAddTrainingData(fhicl::ParameterSet const &pset)
void geom(int which=0)
Definition: geom.C:163
bool HaveTruthInfo() const
Is this a file with truth info in? (Is BackTracker going to be any use to you?)
Definition: BackTracker.h:133
art::ServiceHandle< cheat::BackTracker > fBT
std::vector< std::string > fPreselectionLabels
enum cvn::FinalStateProngs FinalStateProngsType
bool fObeyPreselection
Check rb::IsFiltered?
ParticlesTrainingData * fParticlesTrain
Encapsulate the geometry of one entire detector (near, far, ndos)
enum BeamMode string