70 art::ServiceHandle<cheat::ParticleInventoryService const> pi_serv;
71 art::ServiceHandle<cheat::BackTrackerService const> bt_serv;
72 art::ServiceHandle<geo::Geometry const> geo;
75 art::Handle<std::vector<recob::Cluster>> clustercol;
82 std::vector<art::Ptr<recob::Cluster>> clusters;
83 art::fill_ptr_vector(clusters, clustercol);
86 std::vector<art::Ptr<recob::Cluster>>::iterator itr = clusters.begin();
89 std::map<int, std::vector<std::pair<size_t, art::Ptr<recob::Cluster>>>> eveClusterMap;
92 for (
size_t c = 0; c < clusters.size(); ++c) {
94 std::pair<size_t, art::Ptr<recob::Cluster>> idxPtr(c, clusters[c]);
100 int eveID = floor((*itr)->ID() / 1000.);
102 eveClusterMap[eveID].push_back(idxPtr);
107 std::unique_ptr<std::vector<recob::Track>> trackcol(
new std::vector<recob::Track>);
108 std::unique_ptr<std::vector<recob::SpacePoint>> spcol(
new std::vector<recob::SpacePoint>);
109 std::unique_ptr<art::Assns<recob::Track, recob::SpacePoint>> tspassn(
110 new art::Assns<recob::Track, recob::SpacePoint>);
111 std::unique_ptr<art::Assns<recob::Track, recob::Cluster>> tcassn(
112 new art::Assns<recob::Track, recob::Cluster>);
113 std::unique_ptr<art::Assns<recob::Track, recob::Hit>> thassn(
114 new art::Assns<recob::Track, recob::Hit>);
116 auto const clockData = art::ServiceHandle<detinfo::DetectorClocksService const>()->DataFor(
evt);
118 for (
auto const& clusterMapItr : eveClusterMap) {
121 std::vector<std::pair<size_t, art::Ptr<recob::Cluster>>>
const& eveClusters =
122 clusterMapItr.second;
124 simb::MCParticle* part = pi_serv->ParticleList()[clusterMapItr.first];
128 if (
abs(part->PdgCode()) != 11 &&
abs(part->PdgCode()) != 22 &&
abs(part->PdgCode()) != 111) {
131 std::vector<TVector3> points;
132 std::vector<TVector3> moms;
134 mf::LogInfo(
"TrackCheater")
135 <<
"G4 id " << clusterMapItr.first <<
" is a track with pdg code " << part->PdgCode();
137 art::PtrVector<recob::Cluster> ptrvs;
138 std::vector<size_t> idxs;
139 for (
auto const& idxPtr : eveClusters) {
140 ptrvs.push_back(idxPtr.second);
141 idxs.push_back(idxPtr.first);
145 std::vector<art::Ptr<recob::Hit>> hits;
146 for (
size_t p = 0;
p < ptrvs.size(); ++
p) {
147 std::vector<art::Ptr<recob::Hit>> chits = fmh.at(idxs[
p]);
148 if (!chits.size())
continue;
150 hits.insert(hits.end(), chits.begin(), chits.end());
154 if (hits.size() < 2)
continue;
157 size_t spStart = spcol->size();
158 for (
size_t t = 0; t < hits.size(); ++t) {
159 std::vector<double> xyz = bt_serv->HitToXYZ(clockData, hits[t]);
160 TVector3 point(xyz[0], xyz[1], xyz[2]);
161 points.push_back(point);
163 std::vector<double> xyz1;
167 if (t < hits.size() - 1) { xyz1 = bt_serv->HitToXYZ(clockData, hits[t + 1]); }
169 xyz1 = bt_serv->HitToXYZ(clockData, hits[t - 1]);
174 dx = std::sqrt(std::pow(xyz1[0] - xyz[0], 2) + std::pow(xyz1[1] - xyz[1], 2) +
175 std::pow(xyz1[2] - xyz[2], 2));
179 double drmin = std::numeric_limits<double>::max();
180 for (
unsigned int itp = 0; itp < part->NumberTrajectoryPoints(); itp++) {
181 TVector3
p(part->Vx(itp), part->Vy(itp), part->Vz(itp));
182 double dr = (p - point).Mag();
191 moms.push_back(TVector3(mom * sign * (xyz1[0] - xyz[0]) / dx,
192 mom * sign * (xyz1[1] - xyz[1]) / dx,
193 mom * sign * (xyz1[2] - xyz[2]) / dx));
213 double xyzerr[6] = {1.e-3};
216 spcol->push_back(sp);
219 size_t spEnd = spcol->size();
233 clusterMapItr.first));
240 for (
size_t p = 0; p < ptrvs.size(); ++
p) {
241 hits = fmh.at(idxs[p]);
248 mf::LogInfo(
"TrackCheater") <<
"adding track: \n" << trackcol->back() <<
"\nto collection.";
254 evt.put(std::move(trackcol));
255 evt.put(std::move(spcol));
256 evt.put(std::move(tcassn));
257 evt.put(std::move(thassn));
258 evt.put(std::move(tspassn));
std::string const fCheatedClusterLabel
ROOT::Math::SMatrix< Double32_t, 5, 5, ROOT::Math::MatRepSym< Double32_t, 5 > > SMatrixSym55
std::vector< Vector_t > convertCollToVector(std::vector< Vector > const &coll)
TrackTrajectory::Flags_t Flags_t
A trajectory in space reconstructed from hits.
std::vector< Point_t > convertCollToPoint(std::vector< Point > const &coll)
bool CreateAssn(art::Event &evt, std::vector< T > const &a, art::Ptr< U > const &b, art::Assns< U, T > &assn, std::string a_instance, size_t index=UINT_MAX)
Creates a single one-to-one association.
Track from a non-cascading particle.A recob::Track consists of a recob::TrackTrajectory, plus additional members relevant for a "fitted" track:
Signal from collection planes.