Pandora
Pandora source code navigator
Loading...
Searching...
No Matches
ShowerTrackFinderCheater_tool.cc
Go to the documentation of this file.
1//############################################################################
2//### Name: ShowerTrackFinderCheater ###
3//### Author: Ed Tyley ###
4//### Date: 16.07.19 ###
5//### Description: Cheating tool using truth for shower direction ###
6//############################################################################
7
8//Framework Includes
9#include "art/Utilities/ToolMacros.h"
10
11//LArSoft Includes
12#include "lardata/DetectorInfoServices/DetectorClocksService.h"
13#include "lardata/DetectorInfoServices/DetectorPropertiesService.h"
14#include "lardataobj/RecoBase/Cluster.h"
15#include "lardataobj/RecoBase/PFParticle.h"
18
19namespace ShowerRecoTools {
20
22
23 public:
24 ShowerTrackFinderCheater(const fhicl::ParameterSet& pset);
25
26 //Generic Direction Finder
27 int CalculateElement(const art::Ptr<recob::PFParticle>& pfparticle,
28 art::Event& Event,
29 reco::shower::ShowerElementHolder& ShowerEleHolder) override;
30
31 private:
32 //Algorithm functions
34
35 //fcl
36 const bool fDebugEVD;
37 const art::InputTag fPFParticleLabel;
38 const art::InputTag fHitModuleLabel;
39
40 const std::string fTrueParticleIntputLabel;
42 const std::string fShowerDirectionInputTag;
45 };
46
48 : IShowerTool(pset.get<fhicl::ParameterSet>("BaseTools"))
49 , fLArPandoraShowerCheatingAlg(pset.get<fhicl::ParameterSet>("LArPandoraShowerCheatingAlg"))
50 , fDebugEVD(pset.get<bool>("DebugEVD"))
51 , fPFParticleLabel(pset.get<art::InputTag>("PFParticleLabel"))
52 , fHitModuleLabel(pset.get<art::InputTag>("HitModuleLabel"))
53 , fTrueParticleIntputLabel(pset.get<std::string>("TrueParticleIntputLabel"))
54 , fShowerStartPositionInputTag(pset.get<std::string>("ShowerStartPositionInputTag"))
55 , fShowerDirectionInputTag(pset.get<std::string>("ShowerDirectionInputTag"))
56 , fInitialTrackHitsOutputLabel(pset.get<std::string>("InitialTrackHitsOutputLabel"))
57 , fInitialTrackSpacePointsOutputLabel(
58 pset.get<std::string>("InitialTrackSpacePointsOutputLabel"))
59 {}
60
61 int ShowerTrackFinderCheater::CalculateElement(const art::Ptr<recob::PFParticle>& pfparticle,
62 art::Event& Event,
63 reco::shower::ShowerElementHolder& ShowerEleHolder)
64 {
65
66 const simb::MCParticle* trueParticle;
67
68 //Get the hits from the shower:
69 auto const pfpHandle = Event.getValidHandle<std::vector<recob::PFParticle>>(fPFParticleLabel);
70 auto const clockData =
71 art::ServiceHandle<detinfo::DetectorClocksService const>()->DataFor(Event);
72
73 //Get the clusters
74 auto const clusHandle = Event.getValidHandle<std::vector<recob::Cluster>>(fPFParticleLabel);
75 art::FindManyP<recob::Cluster> fmc(pfpHandle, Event, fPFParticleLabel);
76 std::vector<art::Ptr<recob::Cluster>> clusters = fmc.at(pfparticle.key());
77
78 //Get the hit association
79 art::FindManyP<recob::Hit> fmhc(clusHandle, Event, fPFParticleLabel);
80
81 std::vector<art::Ptr<recob::Hit>> showerHits;
82 for (auto const& cluster : clusters) {
83
84 //Get the hits
85 std::vector<art::Ptr<recob::Hit>> hits = fmhc.at(cluster.key());
86 showerHits.insert(showerHits.end(), hits.begin(), hits.end());
87 }
88
89 if (ShowerEleHolder.CheckElement(fTrueParticleIntputLabel)) {
90 ShowerEleHolder.GetElement(fTrueParticleIntputLabel, trueParticle);
91 }
92 else {
93
94 //Could store these in the shower element holder and just calculate once?
95 std::map<int, const simb::MCParticle*> trueParticles =
97 std::map<int, std::vector<int>> showersMothers =
99
100 //Get the true particle from the shower
101 std::pair<int, double> ShowerTrackInfo =
103 clockData, showersMothers, showerHits, 2);
104
105 if (ShowerTrackInfo.first == -99999) {
106 mf::LogError("ShowerStartPosition") << "True Shower Not Found";
107 return 1;
108 }
109 trueParticle = trueParticles[ShowerTrackInfo.first];
110 ShowerEleHolder.SetElement(trueParticle, fTrueParticleIntputLabel);
111 }
112
113 if (!trueParticle) {
114 mf::LogError("ShowerDirectionCheater") << "True shower not found, returning";
115 return 1;
116 }
117
118 //This is all based on the shower vertex being known. If it is not lets not do the track
119 if (!ShowerEleHolder.CheckElement(fShowerStartPositionInputTag)) {
120 mf::LogError("ShowerTrackFinderCheater") << "Start position not set, returning " << std::endl;
121 return 1;
122 }
123 if (!ShowerEleHolder.CheckElement(fShowerDirectionInputTag)) {
124 mf::LogError("ShowerTrackFinderCheater") << "Direction not set, returning " << std::endl;
125 return 1;
126 }
127
128 geo::Point_t ShowerStartPosition = {-999, -999, -999};
129 ShowerEleHolder.GetElement(fShowerStartPositionInputTag, ShowerStartPosition);
130
131 geo::Vector_t ShowerDirection = {-999, -999, -999};
132 ShowerEleHolder.GetElement(fShowerDirectionInputTag, ShowerDirection);
133
134 auto const hitHandle = Event.getValidHandle<std::vector<recob::Hit>>(fHitModuleLabel);
135
136 // Get the hits associated with the space points
137 art::FindManyP<recob::SpacePoint> fmsph(hitHandle, Event, fPFParticleLabel);
138 if (!fmsph.isValid()) {
139 throw cet::exception("ShowerTrackFinderCheater")
140 << "Spacepoint and hit association not valid. Stopping.";
141 }
142
143 std::vector<int> trueParticleIdVec;
144
145 // If we have an electron, take the hits from the primary only
146 // This will also cover the cases Pandora misclassify a track as a shower
147 if (trueParticle->PdgCode() != 22) { trueParticleIdVec.push_back(trueParticle->TrackId()); }
148 else {
149 // To check if we are rolling up showers, check the number of daughters the photon has
150 const int nDaughters = trueParticle->NumberDaughters();
151 if (nDaughters == 0) {
152 // If we roll up showers, we have no choice but to take all of the hits from the photon
153 trueParticleIdVec.push_back(-trueParticle->TrackId());
154 }
155 else {
156 // If we do not roll up the showers, take all of the primary daughters
157 for (int i = 0; i < nDaughters; i++) {
158 trueParticleIdVec.push_back(trueParticle->Daughter(i));
159 }
160 }
161 }
162
163 std::vector<art::Ptr<recob::Hit>> trackHits;
164 std::vector<art::Ptr<recob::SpacePoint>> trackSpacePoints;
165
166 //Get the hits from the true particle
167 for (auto hit : showerHits) {
168 int trueHitId = fLArPandoraShowerCheatingAlg.TrueParticleID(clockData, hit);
169 if (std::find(trueParticleIdVec.cbegin(), trueParticleIdVec.cend(), trueHitId) !=
170 trueParticleIdVec.cend()) {
171 trackHits.push_back(hit);
172 std::vector<art::Ptr<recob::SpacePoint>> sps = fmsph.at(hit.key());
173 if (sps.size() == 1) { trackSpacePoints.push_back(sps.front()); }
174 }
175 }
176
177 if (trackHits.empty() || trackSpacePoints.empty())
178 mf::LogWarning("ShowerTrackFinderCheater")
179 << "Creating intial track with " << trackHits.size() << " hits and "
180 << trackSpacePoints.size() << " spacepoints" << std::endl;
181
182 ShowerEleHolder.SetElement(trackHits, fInitialTrackHitsOutputLabel);
183 ShowerEleHolder.SetElement(trackSpacePoints, fInitialTrackSpacePointsOutputLabel);
184
185 if (fDebugEVD) {
187 clockData, trueParticle, Event, ShowerEleHolder, pfparticle);
188 }
189
190 return 0;
191 }
192}
193
int CalculateElement(const art::Ptr< recob::PFParticle > &pfparticle, art::Event &Event, reco::shower::ShowerElementHolder &ShowerEleHolder) override
shower::LArPandoraShowerCheatingAlg fLArPandoraShowerCheatingAlg
ShowerTrackFinderCheater(const fhicl::ParameterSet &pset)
int GetElement(const std::string &Name, T &Element) const
void SetElement(T &dataproduct, const std::string &Name, bool checktag=false)
bool CheckElement(const std::string &Name) const
std::pair< int, double > TrueParticleIDFromTrueChain(detinfo::DetectorClocksData const &clockData, std::map< int, std::vector< int > > const &ShowersMothers, std::vector< art::Ptr< recob::Hit > > const &hits, int planeid) const
void CheatDebugEVD(detinfo::DetectorClocksData const &clockData, const simb::MCParticle *trueParticle, art::Event const &Event, reco::shower::ShowerElementHolder &ShowerEleHolder, const art::Ptr< recob::PFParticle > &pfparticle) const
std::map< int, const simb::MCParticle * > GetTrueParticleMap() const
std::map< int, std::vector< int > > GetTrueChain(std::map< int, const simb::MCParticle * > &trueParticles) const
int TrueParticleID(detinfo::DetectorClocksData const &clockData, const art::Ptr< recob::Hit > &hit) const