SlowMonopoleAna_module.cc
Go to the documentation of this file.
1 ////////////////////////////////////////////////////////////////////////
2 // Class: SlowMonopoleAna
3 // Module Type: analyzer
4 // File: SlowMonopoleAna_module.cc
5 //
6 // Generated at Mon Aug 11 11:18:01 2014 by Martin Frank using artmod
7 // from cetpkgsupport v1_06_02.
8 ////////////////////////////////////////////////////////////////////////
9 
19 #include "fhiclcpp/ParameterSet.h"
21 
26 
28 #include "Geometry/Geometry.h"
30 #include "MCCheater/BackTracker.h"
31 #include "RecoBase/CellHit.h"
32 #include "RawData/RawTrigger.h"
33 #include "Simulation/Particle.h"
34 
35 #include <TH1.h>
36 #include <TH2.h>
37 #include <TMath.h>
38 #include <TTree.h>
39 
40 #include <boost/algorithm/string.hpp>
41 
42 #include <math.h>
43 
44 #include <map>
45 #include <tuple>
46 
47 namespace mono
48 {
50  typedef std::tuple<unsigned, double, double> Binning1;
51  typedef std::tuple<unsigned, double, double,
52  unsigned, double, double> Binning2;
53 }
54 
56 {
57 public:
58  explicit SlowMonopoleAna(fhicl::ParameterSet const & p);
59 
60  SlowMonopoleAna(SlowMonopoleAna const &) = delete;
61  SlowMonopoleAna(SlowMonopoleAna &&) = delete;
62  SlowMonopoleAna & operator = (SlowMonopoleAna const &) = delete;
64  void analyze(art::Event const & e) override;
65  void beginJob() override;
66 
67 private:
68  template<class T>
69  void check(T & hists, std::string hist_name);
70  template<class T>
71  void tset(T & tree, std::string const branch_name, double const value);
72 
74  (double ddt_entry_time, double ddt_exit_time,
75  double mc_entry_time, double mc_exit_time) const;
77  void fill(std::string name, double value);
78  void fill(std::string name, std::string label, double weight);
79  void fill(std::string name, std::vector<double> const& values);
80  void fill(std::string name, double x_value, double y_value);
81  void fill_hit_tree(std::vector<art::Ptr<rb::CellHit> > const& hitptrs);
82  void fill_track_tree(std::vector<Track> const& tracks, art::Event const& e);
83  TVector3 get_position(sim::FLSHit const& hit) const;
84  void make_histo
85  (std::string const name, std::string const title, Binning1 const& bin);
86  void make_histo
87  (std::string const name, std::string const title, Binning2 const& bin);
88  double phi(sim::Particle const* particle) const;
89  void print_truth_info(art::Event const& e,
90  bool const& print_hit_details) const;
91  double theta(sim::Particle const* particle) const;
92 
94 
96 
98 
100  std::map<std::string, double> et_, tt_, ht_, rt_;
101  std::map<std::string, TH1*> h1_;
102  std::map<std::string, TH2*> h2_;
106 };
107 
108 
109 
110 template<class T>
112 {
113  if (hists.find(hist_name) == hists.end())
114  {
115  std::cerr << "\n\t Histogram " << hist_name << " does not exist! \n"
116  << std::endl;
117  assert(false);
118  }
119 }
120 
121 
122 
124  EDAnalyzer(p),
125  fill_track_tree_(p.get<bool>("fill_track_tree")),
126  fill_hit_tree_(p.get<bool>("fill_hit_tree")),
127  branch_default_value_(-9999),
128  to_degrees_(180.0 / TMath::Pi())
129 {}
130 
131 
132 
134 {
136  { "algorithm_time" };
137 
139  { "contains_gap_exceeding_cut", "plane_set_size", "cell_set_size",
140  "max_plane_gap", "max_cell_gap",
141  "counts_all", "counts_checked",
142  "counts_after_t_max_cell", "counts_after_t_max_plane" };
143 
144 
145  //
146  // Trees
147  //
148  // Event-Based Tree
149  std::vector<std::string> event_tree_names =
150  { "run_number", "subrun_number", "event_number", "event_length",
151  "n_hits", "n_tracks",
152  "ddt_trigger_decisions", "trigger_decisions",
153  "ddt_trigger_length", "trigger_length",
154  "is_mc", "mc_phi", "mc_theta", "mc_beta", "mc_n_particles",
155  "mc_entry_time_ns", "mc_exit_time_ns", "mc_n_hits", "mc_length",
156  "ddt_entry_time_ns", "ddt_exit_time_ns", "mc_ddt_missed_fraction",
157  "max_plane_gap_cut", "sparsification_factor" };
158  event_tree_names.insert(event_tree_names.end(),
159  smt_parameter_names_.begin(),
160  smt_parameter_names_.end());
161 
162  event_tree_ = tfs_->make<TTree>("Event", "Monopole Event Tree");
163  for (auto const& name : event_tree_names)
164  event_tree_->Branch(name.c_str(), &et_[name], (name + "/D").c_str());
165 
166  // Track-Based Tree
167  std::vector<std::string> track_tree_names =
168  { "run_number", "subrun_number", "event_number", "track_number",
169  "distance", "velocity", "vpro", "distance_in_cm", "beta",
170  "min_plane", "max_plane", "plane_difference", "max_plane_gap",
171  "min_cell", "max_cell", "is_good" };
172  track_tree_names.insert(track_tree_names.end(),
175 
176  track_tree_ = tfs_->make<TTree>("Track", "Monopole Track Tree");
177  for (auto const& name : track_tree_names)
178  track_tree_->Branch(name.c_str(), &tt_[name], (name + "/D").c_str());
179 
180  // Hit-Based Tree
181  std::vector<std::string> hit_tree_names =
182  { "view", "plane", "cell", "tns", "adc",
183  "mc_pdg", "mc_edep", "mc_path_length" };
184  hit_tree_ = tfs_->make<TTree>("Hit", "Monopole Hit Tree");
185  for (auto const& name : hit_tree_names)
186  hit_tree_->Branch(name.c_str(), &ht_[name], (name + "/D").c_str());
187 
188  // Rate Tree
189  std::vector<std::string> rate_tree_names =
190  { "trigger", "max_plane_gap_cut", "event_length" };
191  rate_tree_ = tfs_->make<TTree>("Rate", "Monopole Trigger Rate Tree");
192  for (auto const& name : rate_tree_names)
193  rate_tree_->Branch(name.c_str(), &rt_[name], (name + "/D").c_str());
194 
195 
196  //
197  // Histograms
198  //
199  make_histo("event_counter", "Event Counts", Binning1(10, 0, 10));
201 }
202 
203 
204 
207 {
208  h1_[name] = tfs_->make<TH1D>
209  (name.c_str(), title.c_str(),
210  std::get<0>(bin), std::get<1>(bin), std::get<2>(bin));
211 }
212 
213 
214 
216 (std::string const name, std::string const title, Binning2 const& bin)
217 {
218  h2_[name] = tfs_->make<TH2D>
219  (name.c_str(), title.c_str(),
220  std::get<0>(bin), std::get<1>(bin), std::get<2>(bin),
221  std::get<3>(bin), std::get<4>(bin), std::get<5>(bin));
222 }
223 
224 
225 
227 {
228  const std::vector<std::string> names =
229  { "all", "ddt triggered", "mf triggered", "empty",
230  "multiple ddt triggers", "multiple mf triggers" };
231 
232  for (auto const& name : names)
233  fill("event_counter", name, 0);
234 }
235 
236 
237 
239 (std::string name, double value)
240 {
241  check(h1_, name);
242  h1_[name]->Fill(value);
243 }
244 
245 
246 
249 {
250  check(h1_, name);
251  h1_[name]->Fill(label.c_str(), weight);
252 }
253 
254 
255 
257 (std::string name, std::vector<double> const& values)
258 {
259  check(h1_, name);
260 
261  for (auto const& value : values)
262  h1_[name]->Fill(value);
263 }
264 
265 
266 
268 (std::string name, double x_value, double y_value)
269 {
270  check(h2_, name);
271  h2_[name]->Fill(x_value, y_value);
272 }
273 
274 
275 
276 template<class T>
278 (T & tree, std::string const branch_name, double const value)
279 {
280  auto it = tree.find(branch_name);
281  if (it == tree.end())
282  {
283  std::cerr << "\n\t Branch " << branch_name << " does not exist! \n"
284  << std::endl;
285  assert(false);
286  }
287 
288  it->second = value;
289 }
290 
291 
292 
294 {
295  //
296  // Extract Objects from the Event
297  //
298  // High ADC Filtered Hits
299  art::Handle<novaddt::HitList> ddt_hits_handle;
300  e.getByLabel("smmhighadcfilter", "HEHits", ddt_hits_handle);
301 
302  // Raw Trigger Header
304  e.getByLabel("daq", raw_triggers);
305  assert(raw_triggers->size() == 1);
306  unsigned long long
307  event_time(raw_triggers->at(0).fTriggerTimingMarker_TimeStart);
308  uint32_t event_length(raw_triggers->at(0).fTriggerRange_TriggerLength);
309  // uint32_t extra_info(raw_triggers->front().fTriggerHeader_SourceSubID);
310  // std::bitset<novaddt::smt::Extra_Info::BITSET_SIZE> extra_bits(extra_info);
311 
312  // Offline Cell Hits
313  // art::Handle<std::vector<rb::CellHit> > hits;
314  // e.getByLabel("calhit", hits);
315  // std::vector<art::Ptr<rb::CellHit> > hitptrs;
316  // for (unsigned i = 0; i < hits->size(); ++i)
317  // hitptrs.emplace_back(hits, i);
318 
319  // MC Information
321  e.getByLabel("geantgen", mc_particles);
322 
323  // DDT Trigger Decisions
325  e.getByLabel("slowmonopoletrigger", ddt_decisions);
326 
327 
328 
329 
330  // parse extra_info and print to screen
331  // std::cout << "MF: extra_info = " << extra_info << std::endl;
332  // std::cout << "MF: extra_bits = " << extra_bits << std::endl;
333  // if (extra_bits.test(novaddt::smt::Extra_Info::INVALID_OLD_DEFAULT_BIT))
334  // std::cout << "MF: extra_info is invalid (i.e. set to old default)"
335  // << std::endl;
336  // else
337  // std::cout << "MF: (Check, Timeout, Invalid) = (" << std::boolalpha
338  // << extra_bits.test(novaddt::smt::Extra_Info::CHECK_BIT)
339  // << ", " << extra_bits.test(novaddt::smt::Extra_Info::TIMEOUT_BIT)
340  // << ", " << extra_bits.
341  // test(novaddt::smt::Extra_Info::INVALID_OLD_DEFAULT_BIT)
342  // << ")" << std::endl;
343 
344 
345 
346  //
347  // Run the Local Trigger
348  //
349  std::map<std::string, std::vector<novaddt::TriggerDecision> > decisions;
350  if (!ddt_decisions.failedToGet())
351  decisions["ddt"] = *ddt_decisions;
352 
353  Configuration global_config;
354  global_config.max_plane_gap_cut = 30;
355  global_config.sparsification_factor = 10;
356  // SlowMonopoleTrigger smt(*ddt_hits_handle, event_time, global_config);
357  // decisions["mf"] = smt.mf_trigger_decisions();
358 
359  std::map<int, unsigned> triggers;
360  // triggers[global_config.max_plane_gap_cut] = smt.trigger_decisions().size();
361 
362  // print_truth_info(e, true);
363 
364 
365 
366  //
367  // Fill Trees
368  //
369  // Event-Based Tree
370  for (auto & branch : et_)
371  branch.second = branch_default_value_;
372 
373  tset(et_, "run_number", e.run());
374  tset(et_, "subrun_number", e.subRun());
375  tset(et_, "event_number", e.event());
376  tset(et_, "event_length", event_length);
377 
378  tset(et_, "max_plane_gap_cut", global_config.max_plane_gap_cut);
379  tset(et_, "sparsification_factor", global_config.sparsification_factor);
380 
381  tset(et_, "n_hits", ddt_hits_handle->size());
382 
383  // for (auto const& p_name : smt_parameter_names_)
384  // if (smt.parameter_exists(p_name))
385  // tset(et_, p_name, smt.parameter(p_name));
386 
387  for (auto d : decisions)
388  {
389  tset(et_, d.first + "_trigger_decisions", d.second.size());
390 
391  if (d.second.size() > 0)
392  {
393  std::cout << "MF: " << d.first + "_trigger_length = "
394  << d.second.front().duration() << std::endl;
395  tset(et_, d.first + "_trigger_length", d.second.front().duration());
396  }
397  }
398 
399  float ddt_entry_time(0.0), ddt_exit_time(0.0);
400  if (!ddt_decisions.failedToGet() && !ddt_decisions->empty())
401  {
402  ddt_entry_time = ddt_decisions->front().start() - event_time;
403  ddt_exit_time = ddt_entry_time + ddt_decisions->front().duration();
404 
405  // convert to ns
406  float conversion_factor = 1e9 / 64e6;
407  ddt_entry_time *= conversion_factor;
408  ddt_exit_time *= conversion_factor;
409 
410  tset(et_, "ddt_entry_time_ns", ddt_entry_time);
411  tset(et_, "ddt_exit_time_ns", ddt_exit_time);
412  }
413 
414  if (bt_->HaveTruthInfo())
415  {
416  tset(et_, "is_mc", true);
417 
419  tset(et_, "mc_n_particles", pn.size());
420  if (pn.size() == 1)
421  {
422  sim::Particle const* particle = pn.Primary(0);
423  unsigned const trajectory = 0;
424 
425  tset(et_, "mc_phi", phi(particle));
426  tset(et_, "mc_theta", theta(particle));
427  tset(et_, "mc_beta", particle->P(trajectory) / particle->Mass());
428 
429  std::vector<sim::FLSHit> fls_hits =
430  bt_->ParticleToFLSHit(particle->TrackId());
431  tset(et_, "mc_n_hits", fls_hits.size());
432 
433  std::sort(fls_hits.begin(), fls_hits.end(),
434  [](sim::FLSHit const& lhs, sim::FLSHit const& rhs)
435  { return lhs.GetEntryT() < rhs.GetEntryT(); });
436 
437  auto entry_hit = fls_hits.front();
438  auto exit_hit = fls_hits.back();
439 
440 
441  // Monopole Distance
442  TVector3 entry_position = get_position(entry_hit);
443  std::cout << "MF: entry (x, y, z, plane, cell) = ("
444  << entry_position.x()
445  << ", " << entry_position.y()
446  << ", " << entry_position.z()
447  << ", " << entry_hit.GetPlaneID()
448  << ", " << entry_hit.GetCellID()
449  << ")" << std::endl;
450  TVector3 exit_position = get_position(exit_hit);
451  std::cout << "MF: exit (x, y, z, plane, cell) = ("
452  << exit_position.x()
453  << ", " << exit_position.y()
454  << ", " << exit_position.z()
455  << ", " << exit_hit.GetPlaneID()
456  << ", " << exit_hit.GetCellID()
457  << ")" << std::endl;
458  TVector3 track = exit_position - entry_position;
459  tset(et_, "mc_length", track.Mag());
460  std::cout << "MF: (dx, dy, dz, length) = ("
461  << track.x()
462  << ", " << track.y()
463  << ", " << track.z()
464  << ", " << track.Mag()
465  << ")" << std::endl;
466 
467 
468  // Monopole Timing
469  float mc_entry_time = entry_hit.GetEntryT();
470  float mc_exit_time = exit_hit.GetExitT();
471  tset(et_, "mc_entry_time_ns", mc_entry_time);
472  tset(et_, "mc_exit_time_ns", mc_exit_time);
473 
474  if (!ddt_decisions.failedToGet() && !ddt_decisions->empty())
475  {
476  double time_overlap = calculate_time_overlap
477  (ddt_entry_time, ddt_exit_time, mc_entry_time, mc_exit_time);
478 
479  double overlap_fraction = time_overlap / (mc_exit_time - mc_entry_time);
480  tset(et_, "mc_ddt_missed_fraction", 1 - overlap_fraction);
481  }
482  } else {
483  std::cerr << "\n\tParticle navigator contains "
484  << pn.size() << " particles.\n"
485  << std::endl;
486  }
487  } else {
488  tset(et_, "is_mc", false);
489  }
490 
491  event_tree_->Fill();
492 
493  // Track-Based Tree
494  // if (fill_track_tree_)
495  // fill_track_tree(tracks, e);
496 
497  // Hit-Based Tree
498  // if (fill_hit_tree_)
499  // fill_hit_tree(hitptrs);
500 
501  // Rate Tree
502  for (auto const& trig : triggers)
503  {
504  for (auto & branch : rt_)
505  branch.second = branch_default_value_;
506 
507  tset(rt_, "event_length", event_length);
508  tset(rt_, "max_plane_gap_cut", trig.first);
509  tset(rt_, "trigger", trig.second);
510 
511  std::cout << "MF: Ana: (Plane Cut, MF Trigger) = ("
512  << trig.first << ", " << trig.second << ")"
513  << std::endl;
514 
515  rate_tree_->Fill();
516  }
517 
518 
519 
520  //
521  // Fill Histograms
522  //
523  fill("event_counter", "all", 1);
524 
525  if (ddt_hits_handle->size() == 0)
526  fill("event_counter", "empty", 1);
527 
528  for (auto d : decisions)
529  {
530  if (d.second.size() >= 1)
531  fill("event_counter", d.first + " triggered", 1);
532 
533  if (d.second.size() > 1)
534  fill("event_counter", "multiple " + d.first + " triggered", 1);
535  }
536 }
537 
538 
539 
541 (double ddt_entry_time, double ddt_exit_time,
542  double mc_entry_time, double mc_exit_time) const
543 {
544  double end_time(0);
545  if (ddt_exit_time > mc_entry_time && mc_exit_time > ddt_exit_time)
546  end_time = ddt_exit_time;
547  else if (ddt_exit_time > mc_exit_time)
548  end_time = mc_exit_time;
549  else if (mc_entry_time > ddt_exit_time)
550  return 0.0;
551  else
553  << "Monopole/SlowMonopoleAna_module.cc\n"
554  << "There is a DDT/MC overlap configuration that we have not handled"
555  << std::endl;
556 
557  double start_time(0);
558  if (ddt_entry_time > mc_entry_time && mc_exit_time > ddt_entry_time)
559  start_time = ddt_entry_time;
560  else if (mc_entry_time > ddt_entry_time)
561  start_time = mc_entry_time;
562  else if (ddt_entry_time > mc_exit_time)
563  return 0.0;
564  else
566  << "Monopole/SlowMonopoleAna_module.cc\n"
567  << "There is a DDT/MC overlap configuration that we have not handled"
568  << std::endl;
569 
570  return end_time - start_time;
571 }
572 
573 
574 
576 {
577  TVector3 result;
578  geo_->Plane(hit.GetPlaneID())->Cell(hit.GetCellID())->
579  GetCenter(result, hit.GetZAverage());
580 
581  return result;
582 }
583 
584 
585 
586 double mono::SlowMonopoleAna::phi(sim::Particle const* particle) const
587 {
588  double Px(particle->Px()), Py(particle->Py());
589 
590  return std::atan2(Py, Px);
591 }
592 
593 
594 
595 double mono::SlowMonopoleAna::theta(sim::Particle const* particle) const
596 {
597  double Px(particle->Px()), Py(particle->Py()), Pz(particle->Pz());
598  double Pt = std::sqrt(Px * Px + Py * Py);
599 
600  return std::atan2(Pt, Pz);
601 }
602 
603 
604 
606 (std::vector<Track> const& tracks, art::Event const& e)
607 {
608  unsigned track_number = 0;
609  for (auto const& track : tracks)
610  {
611  for (auto & branch : tt_)
612  branch.second = branch_default_value_;
613 
614  tset(tt_, "run_number", e.run());
615  tset(tt_, "subrun_number", e.subRun());
616  tset(tt_, "event_number", e.event());
617  tset(tt_, "track_number", ++track_number);
618 
619  tset(tt_, "distance", track.distance());
620  tset(tt_, "velocity", track.velocity());
621  tset(tt_, "beta", track.beta());
622 
623  tset(tt_, "min_plane", track.min_plane());
624  tset(tt_, "max_plane", track.max_plane());
625  tset(tt_, "plane_difference", track.plane_difference());
626 
627  tset(tt_, "min_cell", track.min_cell());
628  tset(tt_, "max_cell", track.max_cell());
629 
630  tset(tt_, "is_good", track.is_good());
631 
632  for (auto const& name : smt_track_parameter_names_)
633  tset(tt_, name, track.get(name));
634 
635  track_tree_->Fill();
636  }
637 }
638 
639 
640 
642 (std::vector<art::Ptr<rb::CellHit> > const& hitptrs)
643 {
644  for (auto const& hitptr : hitptrs)
645  {
646  for (auto & branch : ht_)
647  branch.second = branch_default_value_;
648 
649  tset(ht_, "view", hitptr->View());
650  tset(ht_, "plane", hitptr->Plane());
651  tset(ht_, "cell", hitptr->Cell());
652  tset(ht_, "tns", hitptr->TNS());
653  tset(ht_, "adc", hitptr->ADC());
654 
655  if (hitptr->IsMC())
656  {
657  const std::vector<sim::FLSHit> fls_hits = bt_->HitToFLSHit(hitptr);
658 
659  if (fls_hits.size() >= 2)
660  for (auto const& hit: fls_hits)
661  std::cout << "MF: \t FLS Hit PDG = " << hit.GetPDG() << std::endl;
662 
663  if (fls_hits.size() >= 1)
664  {
665  auto fls_hit = fls_hits.front();
666  tset(ht_, "mc_pdg", fls_hit.GetPDG());
667  tset(ht_, "mc_edep", fls_hit.GetEdep());
668  tset(ht_, "mc_path_length", fls_hit.GetTotalPathLength());
669  }
670  }
671 
672  hit_tree_->Fill();
673  }
674 }
675 
676 
677 
679 (art::Event const& e, bool const& print_hit_details) const
680 {
681  if (bt_->HaveTruthInfo())
682  {
684  unsigned n_particle = 0;
685  for (auto p_handle : pn)
686  {
687  const sim::Particle* particle = p_handle.second;
688 
689  std::cout << "MF: MC Particle " << ++n_particle
690  << "\tEvent " << e.event()
691  << "\nMF: Mass = " << particle->Mass()
692  << "\nMF: N Trajectory Points = "
693  << particle->NumberTrajectoryPoints()
694  << std::endl;
695 
696  for (unsigned trajectory = 0;
697  trajectory != particle->NumberTrajectoryPoints(); ++trajectory)
698  {
699  double to_degrees = 180 / TMath::Pi();
700 
701  std::cout << "MF: Trajectory " << trajectory
702  << "\nMF:\t Energy = " << particle->E(trajectory)
703  << "\nMF:\t Momentum = " << particle->P(trajectory)
704  << "\nMF:\t Beta = "
705  << particle->P(trajectory) / particle->Mass()
706  << "\nMF:\t Time = " << particle->T(trajectory)
707  << "\nMF:\t Phi = " << to_degrees * phi(particle)
708  << "\nMF:\t Theta = " << to_degrees * theta(particle)
709  << "\nMF:\t [Vx, Vy, Vz] = ["
710  << particle->Vx(trajectory) << ", "
711  << particle->Vy(trajectory) << ", "
712  << particle->Vz(trajectory) << "]"
713  << "\nMF:\t [Px, Py, Pz] = ["
714  << particle->Px(trajectory) << ", "
715  << particle->Py(trajectory) << ", "
716  << particle->Pz(trajectory) << "]"
717  << std::endl;
718  }
719 
720  if (print_hit_details)
721  {
722  std::vector<sim::FLSHit> fls_hits =
723  bt_->ParticleToFLSHit(particle->TrackId());
724  std::cout << "MF: Hits:" << std::endl;
725 
726  std::sort(fls_hits.begin(), fls_hits.end(),
727  [](sim::FLSHit const& lhs, sim::FLSHit const& rhs)
728  { return lhs.GetPlaneID() < rhs.GetPlaneID(); });
729 
730  for (auto const& fls_hit : fls_hits)
731  {
732  float dt = fls_hit.GetExitT() - fls_hit.GetEntryT();
733  float dx = fls_hit.GetTotalPathLength();
734  float v = dx / dt * 1.0e7;
735  float beta = v / 3.0e8;
736  float dEdx = fls_hit.GetEdep() / dx;
737 
738  std::cout << "MF:\t [plane, cell, beta, dEdx, path, Edep, "
739  << "EntryT, ExitT] = ["
740  << fls_hit.GetPlaneID() << ", "
741  << fls_hit.GetCellID() << ", "
742  << beta << ", "
743  << dEdx << ", "
744  << fls_hit.GetTotalPathLength() << ", "
745  << fls_hit.GetEdep() << ", "
746  << fls_hit.GetEntryT() << ", "
747  << fls_hit.GetExitT()
748  << "]" << std::endl;
749 
750  for (int point = 0; point != fls_hit.GetNPoints(); ++point)
751  {
752  std::cout << "MF:\t\t [x, y, z] = ["
753  << fls_hit.GetX(point) << ", "
754  << fls_hit.GetY(point) << ", "
755  << fls_hit.GetZ(point) << "]" << std::endl;
756  }
757  }
758  }
759  }
760  }
761 }
762 
763 
764 
double E(const int i=0) const
Definition: MCParticle.h:232
const XML_Char * name
Definition: expat.h:151
unsigned int NumberTrajectoryPoints() const
Definition: MCParticle.h:217
SubRunNumber_t subRun() const
Definition: Event.h:72
back track the reconstruction to the simulation
void check(T &hists, std::string hist_name)
void print_truth_info(art::Event const &e, bool const &print_hit_details) const
double Py(const int i=0) const
Definition: MCParticle.h:230
set< int >::iterator it
int GetPlaneID() const
Plane ID.
Definition: FLSHit.h:37
Definition: Cluster.h:13
const sim::ParticleNavigator & ParticleNavigator() const
Get a reference to the ParticleNavigator.
Definition: BackTracker.h:744
const Var weight
size_type size() const
int GetCellID() const
Cell ID.
Definition: FLSHit.h:39
art::ServiceHandle< art::TFileService > tfs_
std::vector< std::string > smt_track_parameter_names_
const char * p
Definition: xmltok.h:285
std::map< std::string, double > ht_
T sqrt(T number)
Definition: d0nt_math.hpp:156
std::vector< sim::FLSHit > ParticleToFLSHit(const int &trackID) const
All the FLSHits that were created by the track id trackID, sorted from most to least light...
double Mass() const
Definition: MCParticle.h:238
double Px(const int i=0) const
Definition: MCParticle.h:229
Definition: event.h:19
OStream cerr
Definition: OStream.cxx:7
A single unit of energy deposition in the liquid scintillator.
Definition: FLSHit.h:19
const PlaneGeo * Plane(unsigned int i) const
Double_t beta
TString hists[nhists]
Definition: bdt_com.C:3
void make_histo(std::string const name, std::string const title, Binning1 const &bin)
DEFINE_ART_MODULE(TestTMapFile)
double phi(sim::Particle const *particle) const
TVector3 get_position(sim::FLSHit const &hit) const
int TrackId() const
Definition: MCParticle.h:209
const char * label
SlowMonopoleAna & operator=(SlowMonopoleAna const &)=delete
std::map< std::string, TH2 * > h2_
Definition: Cand.cxx:23
void fill_hit_tree(std::vector< art::Ptr< rb::CellHit > > const &hitptrs)
double dx[NP][NC]
const sim::Particle * Primary(const int) const
const XML_Char int const XML_Char * value
Definition: expat.h:331
SlowMonopoleAna(fhicl::ParameterSet const &p)
std::vector< sim::FLSHit > HitToFLSHit(const rb::CellHit &hit) const
All the FLSHits that contributed to this hit, sorted from most to least light.
double P(const int i=0) const
Definition: MCParticle.h:233
double T(const int i=0) const
Definition: MCParticle.h:223
void fill(std::string name, double value)
std::tuple< unsigned, double, double > Binning1
Float_t d
Definition: plot.C:236
std::tuple< unsigned, double, double, unsigned, double, double > Binning2
std::map< std::string, double > et_
std::map< std::string, TH1 * > h1_
EventNumber_t event() const
Definition: Event.h:67
void tset(T &tree, std::string const branch_name, double const value)
EDAnalyzer(Table< Config > const &config)
Definition: EDAnalyzer.h:100
double theta(sim::Particle const *particle) const
float bin[41]
Definition: plottest35.C:14
art::ServiceHandle< cheat::BackTracker > bt_
OStream cout
Definition: OStream.cxx:6
float GetEntryT() const
Definition: FLSHit.h:51
cet::coded_exception< errors::ErrorCodes, ExceptionDetail::translate > Exception
Definition: Exception.h:66
art::ServiceHandle< geo::Geometry > geo_
double Vx(const int i=0) const
Definition: MCParticle.h:220
::xsd::cxx::tree::string< char, simple_type > string
Definition: Database.h:154
T * make(ARGS...args) const
double calculate_time_overlap(double ddt_entry_time, double ddt_exit_time, double mc_entry_time, double mc_exit_time) const
std::map< std::string, double > tt_
bool getByLabel(std::string const &label, std::string const &productInstanceName, Handle< PROD > &result) const
Definition: DataViewImpl.h:344
Definition: structs.h:12
bool HaveTruthInfo() const
Is this a file with truth info in? (Is BackTracker going to be any use to you?)
Definition: BackTracker.h:133
double Pz(const int i=0) const
Definition: MCParticle.h:231
double Vz(const int i=0) const
Definition: MCParticle.h:222
assert(nhit_max >=nhit_nbins)
std::vector< std::string > smt_parameter_names_
double T
Definition: Xdiff_gwt.C:5
void analyze(art::Event const &e) override
float GetZAverage(const int step) const
Get Z-average for the step. This is in local coordinates.
Definition: FLSHit.h:80
Float_t e
Definition: plot.C:35
std::map< std::string, double > rt_
RunNumber_t run() const
Definition: Event.h:77
double Vy(const int i=0) const
Definition: MCParticle.h:221
Encapsulate the geometry of one entire detector (near, far, ndos)
void fill_track_tree(std::vector< Track > const &tracks, art::Event const &e)
bool failedToGet() const
Definition: Handle.h:196
T atan2(T number)
Definition: d0nt_math.hpp:72