ShowerMuonCoincidence_module.cc
Go to the documentation of this file.
1 ////////////////////////////////////////////////////////////////////////
2 // Class: ShowerMuonCoincidence
3 // Plugin Type: analyzer (art v2_12_01)
4 // File: ShowerMuonCoincidence_module.cc
5 //
6 // Generated at Mon Jul 15 11:08:30 2019 by Dung Phan using cetskelgen
7 // from cetlib version v3_06_00.
8 ////////////////////////////////////////////////////////////////////////
9 
23 #include "fhiclcpp/ParameterSet.h"
25 
26 #include "HoughTransform.h"
27 #include "NNbarUtilities.h"
28 #include "RecoBase/CellHit.h"
29 #include "RecoBase/Cluster.h"
30 #include "RecoBase/Prong.h"
31 #include "RecoBase/RecoHit.h"
32 #include "RecoBase/Vertex.h"
33 #include "Utilities/AssociationUtil.h"
37 
38 #include <cmath>
39 #include <iterator>
40 #include <math.h>
41 #include <string>
42 #include <vector>
43 
44 #include <TCanvas.h>
45 #include <TColor.h>
46 #include <TGraph.h>
47 #include <TH1.h>
48 #include <TH2.h>
49 #include <TLatex.h>
50 #include <TLegend.h>
51 #include <TMath.h>
52 #include <TROOT.h>
53 #include <TStyle.h>
54 #include <TTree.h>
55 
56 namespace nnbar {
57 class ShowerMuonCoincidence;
58 }
59 
61  public:
67 
68  // Required functions.
69  void analyze(art::Event const &e) override;
70  double distancePointToLine(TVector3 point, TVector3 pointOnLineA, TVector3 pointOnLineB);
71 
72  // Selected optional functions.
73  void beginJob() override;
74  void beginRun(art::Run const &r) override;
75  void beginSubRun(art::SubRun const &sr) override;
76  void endJob() override;
77  void endRun(art::Run const &r) override;
78  void endSubRun(art::SubRun const &sr) override;
79 
80  private:
84 
86 };
87 
89  fhicl::ParameterSet const &p)
90  : EDAnalyzer(p) {
91  _ShowerClusterLabel = "shower:showers";
92  _MuonClusterLabel = "shower:cosmics";
93  _CosmicTrackLabel = "windowtrack";
94 }
95 
98  distFromShowerToCoincidenceMuon = tfs->make<TH1D>("distFromShowerToCoincidenceMuon", ";Distance (cm);", 300, 0, 6000);
99 }
100 
102 
104 
106  art::Handle<std::vector<rb::Cluster> > muonClusterHandle;
107  e.getByLabel(_MuonClusterLabel, muonClusterHandle);
108 
109  art::Handle<std::vector<rb::Cluster> > showerClusterHandle;
110  e.getByLabel(_ShowerClusterLabel, showerClusterHandle);
111 
112  art::FindManyP<rb::Track> fmtrack(muonClusterHandle, e, "windowtrack");
113 
114  for (unsigned int iShower = 0; iShower < showerClusterHandle->size(); iShower++) {
115  rb::Cluster showerSlice = (*showerClusterHandle)[iShower];
116  auto hitTimeMin = showerSlice.MinTNS();
117  auto hitTimeMax = showerSlice.MaxTNS();
118 
119  auto coincidenceWindowMin = hitTimeMin - 500;
120  auto coincidenceWindowMax = hitTimeMax + 500;
121  unsigned int countSlicers = 0;
122  std::vector<art::Ptr<rb::Track>> tracks;
123  for (unsigned int iMuon = 0; iMuon < muonClusterHandle->size(); iMuon++){
124  rb::Cluster muonSlice = (*muonClusterHandle)[iMuon];
125  auto muonHitTimeMin = muonSlice.MinTNS();
126  auto muonHitTimeMax = muonSlice.MaxTNS();
127  bool muonBeforeShower = coincidenceWindowMin <= muonHitTimeMax && muonHitTimeMax <= coincidenceWindowMax;
128  bool muonAfterShower = coincidenceWindowMin <= muonHitTimeMin && muonHitTimeMin <= coincidenceWindowMax;
129  if (muonBeforeShower && muonAfterShower) {
130  countSlicers++;
131  if (fmtrack.isValid()){
132  tracks = fmtrack.at(iMuon);
133  }
134  }
135  }
136 
137  if (countSlicers == 1) {
138  // find vertex
139  double vertexXinCM;
140  double vertexYinCM;
141  double vertexZinCM;
142  nnbar::NNbarUtilities::energyBalancingVertex(showerSlice, vertexXinCM, vertexYinCM, vertexZinCM);
143  TVector3 point(vertexXinCM, vertexYinCM, vertexZinCM);
144 
145  // check associated track
146  if (tracks.size() == 1) {
147  if (tracks.at(0)->Is3D()) {
148  auto startPoint = tracks.at(0)->Start();
149  auto endPoint = tracks.at(0)->Stop();
150  double dist = distancePointToLine(point, startPoint, endPoint);
151  std::cout << "dist: " << dist << std::endl;
153  }
154  }
155  }
156  }
157 }
158 
159 double nnbar::ShowerMuonCoincidence::distancePointToLine(TVector3 point0, TVector3 point1, TVector3 point2) {
160  // See: http://mathworld.wolfram.com/Point-LineDistance3-Dimensional.html
161 
162  double distSq21 = (point2 - point1).Mag2();
163  double distSq10 = (point1 - point0).Mag2();
164  double scalarProd = (point1 - point0).Dot(point2 - point1);
165  double numerator = distSq10 * distSq21 - TMath::Power(scalarProd, 2);
166  return TMath::Sqrt(numerator / distSq21);
167 }
168 
170 }
171 
173 
175 
void beginSubRun(art::SubRun const &sr) override
float Dot(const Proxy &v) const
const char * p
Definition: xmltok.h:285
A collection of associated CellHits.
Definition: Cluster.h:47
static void energyBalancingVertex(rb::Cluster cluster, double &vertexX, double &vertexY, double &vertexZ)
DEFINE_ART_MODULE(TestTMapFile)
Particle class.
Definition: Run.h:31
double dist
Definition: runWimpSim.h:113
void beginRun(art::Run const &r) override
double MinTNS() const
Definition: Cluster.cxx:482
ShowerMuonCoincidence(fhicl::ParameterSet const &p)
void analyze(art::Event const &e) override
void endSubRun(art::SubRun const &sr) override
caf::StandardRecord * sr
ShowerMuonCoincidence & operator=(ShowerMuonCoincidence const &)=delete
EDAnalyzer(Table< Config > const &config)
Definition: EDAnalyzer.h:100
Vertex location in position and time.
double MaxTNS() const
Definition: Cluster.cxx:528
OStream cout
Definition: OStream.cxx:6
float Mag2() const
double distancePointToLine(TVector3 point, TVector3 pointOnLineA, TVector3 pointOnLineB)
T * make(ARGS...args) const
bool getByLabel(std::string const &label, std::string const &productInstanceName, Handle< PROD > &result) const
Definition: DataViewImpl.h:344
Var Sqrt(const Var &v)
Use to take sqrt of a var.
Definition: Var.cxx:326
TRandom3 r(0)
Float_t e
Definition: plot.C:35
void endRun(art::Run const &r) override