29 #include "TStopwatch.h" 32 class CosmicTrackTagger;
55 auto const*
geo = lar::providerFrom<geo::Geometry>();
70 const double driftVelocity = detp.DriftVelocity(detp.Efield(), detp.Temperature());
73 2 *
geo->DetHalfWidth() / (driftVelocity * fSamplingRate / 1000);
74 fMinTickDrift = clock_data.Time2Tick(clock_data.TriggerTime());
77 produces<std::vector<anab::CosmicTag>>();
78 produces<art::Assns<recob::Track, anab::CosmicTag>>();
85 std::unique_ptr<std::vector<anab::CosmicTag>> cosmicTagTrackVector(
86 new std::vector<anab::CosmicTag>);
87 std::unique_ptr<art::Assns<recob::Track, anab::CosmicTag>> assnOutCosmicTagTrack(
94 std::vector<art::Ptr<recob::Track>> TrkVec;
103 for (
unsigned int iTrack = 0; iTrack < Trk_h->size(); iTrack++) {
109 std::vector<art::Ptr<recob::Hit>> HitVec = hitsSpill.at(iTrack);
111 if (iTrack != tTrack.
key()) { std::cout <<
"Mismatch in track index/key" << std::endl; }
114 auto tVector1 = tTrack->
Vertex();
115 auto tVector2 = tTrack->
End();
117 float trackEndPt1_X = tVector1.X();
118 float trackEndPt1_Y = tVector1.Y();
119 float trackEndPt1_Z = tVector1.Z();
120 float trackEndPt2_X = tVector2.X();
121 float trackEndPt2_Y = tVector2.Y();
122 float trackEndPt2_Z = tVector2.Z();
124 if (trackEndPt1_X != trackEndPt1_X || trackEndPt1_Y != trackEndPt1_Y ||
125 trackEndPt1_Z != trackEndPt1_Z || trackEndPt2_X != trackEndPt2_X ||
126 trackEndPt2_Y != trackEndPt2_Y || trackEndPt2_Z != trackEndPt2_Z) {
127 std::cerr <<
"!!! FOUND A PROBLEM... the length is: " << tTrack->
Length()
129 << tTrack << std::endl;
130 std::vector<float> tempPt1, tempPt2;
131 tempPt1.push_back(-999);
132 tempPt1.push_back(-999);
133 tempPt1.push_back(-999);
134 tempPt2.push_back(-999);
135 tempPt2.push_back(-999);
136 tempPt2.push_back(-999);
137 cosmicTagTrackVector->emplace_back(tempPt1, tempPt2, -999, tag_id);
149 for (
unsigned int p = 0; p < HitVec.size(); p++) {
151 std::cout <<
"###>> Hit key: " << HitVec[p].key()
152 <<
", peak - RMS: " << HitVec[p]->PeakTimeMinusRMS()
153 <<
", peak + RMS: " << HitVec[p]->PeakTimePlusRMS() << std::endl;
155 if (HitVec[p]->PeakTimeMinusRMS() < tick1) tick1 = HitVec[p]->PeakTimeMinusRMS();
156 if (HitVec[p]->PeakTimePlusRMS() > tick2) tick2 = HitVec[p]->PeakTimePlusRMS();
170 int nBdY = 0, nBdZ = 0;
189 if ((nBdY + nBdZ) > 1) {
198 else if ((nBdY + nBdZ) == 1) {
207 std::vector<float> endPt1;
208 std::vector<float> endPt2;
209 endPt1.push_back(trackEndPt1_X);
210 endPt1.push_back(trackEndPt1_Y);
211 endPt1.push_back(trackEndPt1_Z);
212 endPt2.push_back(trackEndPt2_X);
213 endPt2.push_back(trackEndPt2_Y);
214 endPt2.push_back(trackEndPt2_Z);
216 float cosmicScore = isCosmic > 0 ? 1 : 0;
217 if (isCosmic == 3) cosmicScore = 0.5;
228 cosmicTagTrackVector->emplace_back(endPt1, endPt2, cosmicScore, tag_id);
237 float dE = 0, dS = 0, temp = 0, IScore = 0;
238 unsigned int IndexE = 0, iTrk1 = 0, iTrk = 0;
241 for (iTrk = 0; iTrk < Trk_h->size(); iTrk++) {
243 if ((*cosmicTagTrackVector)[iTrk].CosmicScore() == 0) {
244 auto tStart = tTrk->
Vertex();
245 auto tEnd = tTrk->
End();
247 for (iTrk1 = 0; iTrk1 < Trk_h->size(); iTrk1++) {
249 float getScore = (*cosmicTagTrackVector)[iTrk1].CosmicScore();
250 if (getScore == 1 || getScore == 0.5) {
252 auto tStart1 = tTrk1->
Vertex();
253 auto tEnd1 = tTrk1->
End();
254 auto NumE = (tEnd - tStart1).Cross(tEnd - tEnd1);
255 auto DenE = tEnd1 - tStart1;
256 dE = NumE.R() / DenE.R();
273 auto tStartI = tTrkI->
Vertex();
274 auto tEndI = tTrkI->
End();
275 auto NumS = (tStart - tStartI).Cross(tStart - tEndI);
276 auto DenS = tEndI - tStartI;
277 dS = NumS.R() / DenS.R();
278 if (((dS < 5 && temp < 5) || (dS < temp && dS < 5)) && (tTrk->
Length() < 60)) {
279 (*cosmicTagTrackVector)[iTrk].CosmicScore() = IScore - 0.05;
280 (*cosmicTagTrackVector)[iTrk].CosmicType() = IType;
285 e.
put(std::move(cosmicTagTrackVector));
286 e.
put(std::move(assnOutCosmicTagTrack));
std::string getType(cet::LibraryManager const &lm, std::string const &fullSpec)
Utilities related to art service access.
enum anab::cosmic_tag_id CosmicTagID_t
Declaration of signal hit object.
EDProducer(fhicl::ParameterSet const &pset)
size_t NumberTrajectoryPoints() const
Various functions related to the presence and the number of (valid) points.
PutHandle< PROD > put(std::unique_ptr< PROD > &&edp, std::string const &instance={})
double Length(size_t p=0) const
Access to various track properties.
#define DEFINE_ART_MODULE(klass)
key_type key() const noexcept
Point_t const & Vertex() const
Access to track position at different points.
Provides recob::Track data product.
std::string fTrackModuleLabel
void produce(art::Event &e) override
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.
bool getByLabel(std::string const &label, std::string const &instance, Handle< PROD > &result) const
Utility object to perform functions of association.
Point_t const & End() const
Access to track position at different points.
CosmicTrackTagger(fhicl::ParameterSet const &p)
void fill_ptr_vector(std::vector< Ptr< T >> &ptrs, H const &h)
double sampling_rate(DetectorClocksData const &data)
Returns the period of the TPC readout electronics clock.
Namespace collecting geometry-related classes utilities.
art framework interface to geometry description