61 unsigned int StopAfter)
const 64 std::vector<recob::Seed> ReturnVector;
67 std::vector<recob::SpacePoint> spts;
72 if (
int(spts.size()) <
fMinPointsInSeed) {
return std::vector<recob::Seed>(); }
80 std::vector<std::vector<std::vector<int>>> OrgHits(3);
81 for (
size_t n = 0;
n != 3; ++
n)
86 std::vector<std::vector<int>> SpacePointsPerHit(HitsFlat.
size(), std::vector<int>());
87 std::vector<std::vector<int>> HitsPerSpacePoint(spts.size(), std::vector<int>());
95 std::vector<char> HitStatus(HitsFlat.
size(), 0);
99 for (
size_t i = 0; i != HitsFlat.
size(); ++i) {
100 OrgHits[HitsFlat.
at(i)->View()][HitsFlat.
at(i)->Channel()].push_back(i);
105 for (
size_t iSP = 0; iSP != spts.size(); ++iSP) {
108 for (
size_t iH = 0; iH != HitsThisSP.
size(); ++iH) {
110 uint32_t ThisChannel = HitsThisSP.
at(iH)->Channel();
111 float ThisTime = HitsThisSP.
at(iH)->PeakTime();
114 for (
size_t iOrg = 0; iOrg != OrgHits[ThisView][ThisChannel].size(); ++iOrg) {
115 if (fabs(ThisTime - HitsFlat.
at(OrgHits[ThisView][ThisChannel][iOrg])->PeakTime()) <
117 SpacePointsPerHit.at(OrgHits[ThisView][ThisChannel][iOrg]).push_back(iSP);
118 HitsPerSpacePoint.at(iSP).push_back(OrgHits[ThisView][ThisChannel][iOrg]);
138 std::vector<char> PointStatus(spts.size(), 0);
140 std::vector<std::map<geo::View_t, std::vector<int>>> WhichHitsPerSeed;
142 bool KeepChopping =
true;
144 while (KeepChopping) {
146 std::vector<int> PointsUsed;
150 FindSeedAtEnd(detProp, spts, PointStatus, PointsUsed, HitsFlat, OrgHits);
157 for (
size_t iP = 0; iP != PointsUsed.size(); ++iP) {
158 for (
size_t iH = 0; iH != HitsPerSpacePoint.at(PointsUsed.at(iP)).
size(); ++iH) {
159 int UsedHitID = HitsPerSpacePoint.at(PointsUsed.at(iP)).at(iH);
160 HitStatus[UsedHitID] = 2;
163 PointStatus[PointsUsed.at(0)] = 1;
164 ConsolidateSeed(detProp, TheSeed, HitsFlat, HitStatus, OrgHits,
false);
169 std::vector<char> HitStatusGood;
171 for (
size_t r = 0;
r != (
unsigned int)
fRefits; ++
r) {
175 HitStatusGood = HitStatus;
177 std::vector<int> PresentHitList;
178 for (
size_t iH = 0; iH != HitStatus.size(); ++iH) {
179 if (HitStatus[iH] == 2) { PresentHitList.push_back(iH); }
181 double pt[3],
dir[3], err[3];
187 std::vector<double> ViewRMS;
188 std::vector<int> HitsPerView;
190 detProp, HitsFlat, PresentHitList, Center, Direction, ViewRMS, HitsPerView);
192 Direction = Direction.Unit() * TheSeed.
GetLength();
194 int nViewsWithHits(0);
195 for (
size_t n = 0;
n != 3; ++
n) {
197 dir[
n] = Direction[
n];
202 if (HitsPerView[
n] > 0) nViewsWithHits++;
205 if (nViewsWithHits < 2) TheSeed.
SetValidity(
false);
214 HitStatus = HitStatusGood;
227 WhichHitsPerSeed.push_back(std::map<
geo::View_t, std::vector<int>>());
230 for (
size_t iH = 0; iH != HitStatus.size(); ++iH) {
231 if (HitStatus.at(iH) == 2) {
232 WhichHitsPerSeed.at(WhichHitsPerSeed.size() - 1)[HitsFlat[iH]->View()].push_back(iH);
234 HitStatus.at(iH) = 1;
236 for (
size_t iSP = 0; iSP != SpacePointsPerHit.at(iH).size(); ++iSP) {
237 PointStatus[SpacePointsPerHit.at(iH).at(iSP)] = 1;
243 ReturnVector.push_back(TheSeed);
244 CataloguedHits.push_back(HitsWithThisSeed);
247 HitsWithThisSeed.
clear();
251 PointStatus.at(PointsUsed.at(0)) = 2;
254 int TotalSPsUsed = 0;
255 for (
size_t i = 0; i != PointStatus.size(); ++i) {
256 if (PointStatus[i] != 0) TotalSPsUsed++;
259 if ((
int(spts.size()) - TotalSPsUsed) <
fMinPointsInSeed) KeepChopping =
false;
261 if ((PointStatus[0] == 3) || (PointStatus.size() == 0)) KeepChopping =
false;
265 if ((ReturnVector.size() >= StopAfter) && (StopAfter > 0))
break;
273 if (ReturnVector.size() == 0) {
274 std::vector<int> ListAllHits;
275 for (
size_t i = 0; i != HitsFlat.
size(); ++i) {
276 ListAllHits.push_back(i);
279 TVector3 SeedCenter(0, 0, 0);
280 TVector3 SeedDirection(0, 0, 0);
282 std::vector<double> ViewRMS;
283 std::vector<int> HitsPerView;
285 std::vector<art::PtrVector<recob::Hit>> HitsInThisCollection(3);
288 detProp, HitsFlat, ListAllHits, SeedCenter, SeedDirection, ViewRMS, HitsPerView);
290 bool ThrowOutSeed =
false;
292 double PtArray[3], DirArray[3];
293 int nViewsWithHits(0);
294 for (
size_t n = 0;
n != 3; ++
n) {
295 PtArray[
n] = SeedCenter[
n];
296 DirArray[
n] = SeedDirection[
n];
297 if (HitsPerView[
n] > 0) nViewsWithHits++;
301 if (nViewsWithHits < 2 || (nViewsWithHits < 3 && !
fAllow2DSeeds)) ThrowOutSeed =
true;
304 ConsolidateSeed(detProp, TheSeed, HitsFlat, HitStatus, OrgHits,
false);
309 for (
size_t i = 0; i != HitStatus.size(); ++i) {
310 if (HitStatus.at(i) == 2) ListAllHits.push_back(i);
312 std::vector<int> HitsPerView;
314 detProp, HitsFlat, ListAllHits, SeedCenter, SeedDirection, ViewRMS, HitsPerView);
316 int nViewsWithHits(0);
317 for (
size_t n = 0;
n != 3; ++
n) {
318 PtArray[
n] = SeedCenter[
n];
319 DirArray[
n] = SeedDirection[
n];
321 if (HitsPerView[
n] > 0) nViewsWithHits++;
324 if (nViewsWithHits < 2 || (nViewsWithHits < 3 && !
fAllow2DSeeds)) ThrowOutSeed =
true;
330 if (
fMaxViewRMS.at(j) < ViewRMS.at(j)) { ThrowOutSeed =
true; }
335 if ((!ThrowOutSeed) && (TheSeed.
IsValid())) {
336 ReturnVector.push_back(TheSeed);
338 for (
size_t i = 0; i != ListAllHits.size(); ++i) {
339 HitsThisSeed.
push_back(HitsFlat.
at(ListAllHits.at(i)));
341 CataloguedHits.push_back(HitsThisSeed);
346 SpacePointsPerHit.clear();
347 HitsPerSpacePoint.clear();
352 for (
size_t i = 0; i != ReturnVector.size(); ++i) {
353 double CrazyValue = 1000000;
354 double Length = ReturnVector.at(i).GetLength();
355 if (!((Length >
fLengthCut) && (Length < CrazyValue))) {
356 ReturnVector.erase(ReturnVector.begin() + i);
357 CataloguedHits.erase(CataloguedHits.begin() + i);
372 std::vector<char>& HitStatus,
377 bool ThrowOutSeed =
false;
380 std::map<geo::View_t, std::map<uint32_t, std::vector<int>>> HitsInThisSeed;
382 int NHitsThisSeed = 0;
384 double MinS = 1000, MaxS = -1000;
385 for (
size_t i = 0; i != HitStatus.size(); ++i) {
386 if (HitStatus.at(i) == 2) {
396 if (s < MinS) MinS = s;
397 if (s > MaxS) MaxS = s;
398 HitsInThisSeed[HitsFlat.
at(i)->View()][HitsFlat.
at(i)->Channel()].push_back(i);
403 double LengthRescale = (MaxS - MinS) / 2.;
404 double PositionShift = (MaxS + MinS) / 2.;
406 double pt[3],
dir[3], err[3];
410 for (
size_t n = 0;
n != 3; ++
n) {
411 pt[
n] += dir[
n] * PositionShift;
412 dir[
n] *= LengthRescale;
419 for (
auto itP = HitsInThisSeed.begin(); itP != HitsInThisSeed.end(); ++itP) {
422 uint32_t LowestChan = itP->second.begin()->first;
423 uint32_t HighestChan = itP->second.rbegin()->first;
424 for (uint32_t c = LowestChan; c != HighestChan; ++c) {
425 for (
size_t h = 0; h != OrgHits[View][c].size(); ++h) {
426 if (HitStatus[OrgHits[View][c].at(h)] == 0) {
431 HitStatus[OrgHits[View][c].at(h)] = 2;
433 HitsInThisSeed[View][c].push_back(OrgHits[View][c].at(h));
436 HitStatus[OrgHits[View][c].at(h)] = 0;
442 if (NHitsThisSeed == 0) ThrowOutSeed =
true;
446 uint32_t LowestChanInSeed[3], HighestChanInSeed[3];
447 double Occupancy[] = {0., 0., 0.};
448 int nHitsPerView[] = {0, 0, 0};
450 for (
auto itP = HitsInThisSeed.begin(); itP != HitsInThisSeed.end(); ++itP) {
454 LowestChanInSeed[View] = itP->second.begin()->first;
455 HighestChanInSeed[View] = itP->second.rbegin()->first;
457 nHitsPerView[View]++;
459 int FilledChanCount = 0;
461 for (
size_t c = LowestChanInSeed[View]; c != HighestChanInSeed[View]; ++c) {
462 if (itP->second[c].size() > 0) ++FilledChanCount;
466 float(FilledChanCount) / float(HighestChanInSeed[View] - LowestChanInSeed[View]);
470 int nViewsWithHits(0);
471 for (
size_t n = 0;
n != 3; ++
n) {
473 if (nHitsPerView[
n] > 0) nViewsWithHits++;
480 if (nBelowCut > belowCut) ThrowOutSeed =
true;
482 if ((Extend) && (!ThrowOutSeed)) {
483 std::vector<std::vector<double>> ToAddNegativeS(3, std::vector<double>());
484 std::vector<std::vector<double>> ToAddPositiveS(3, std::vector<double>());
485 std::vector<std::vector<int>> ToAddNegativeH(3, std::vector<int>());
486 std::vector<std::vector<int>> ToAddPositiveH(3, std::vector<int>());
488 for (
auto itP = HitsInThisSeed.begin(); itP != HitsInThisSeed.end(); ++itP) {
493 if (LowestChanInSeed[View] > 0) {
494 for (uint32_t c = LowestChanInSeed[View] - 1; c != 0; --c) {
495 bool GotOneThisChannel =
false;
496 for (
size_t h = 0; h != OrgHits[View][c].size(); ++h) {
497 if (HitStatus[OrgHits[View][c][h]] == 0) {
500 GotOneThisChannel =
true;
502 ToAddNegativeS[View].push_back(s);
503 ToAddNegativeH[View].push_back(OrgHits[View][c].at(h));
506 ToAddPositiveS[View].push_back(s);
507 ToAddPositiveH[View].push_back(OrgHits[View][c].at(h));
512 if (GotOneThisChannel ==
false)
break;
517 for (uint32_t c = HighestChanInSeed[View] + 1; c !=
fNChannels; ++c) {
518 bool GotOneThisChannel =
false;
519 for (
size_t h = 0; h != OrgHits[View][c].size(); ++h) {
520 if (HitStatus[OrgHits[View][c][h]] == 0) {
523 GotOneThisChannel =
true;
526 ToAddNegativeS[View].push_back(s);
527 ToAddNegativeH[View].push_back(OrgHits[View][c].at(h));
530 ToAddPositiveS[View].push_back(s);
531 ToAddPositiveH[View].push_back(OrgHits[View][c].at(h));
536 if (GotOneThisChannel ==
false)
break;
540 double ExtendPositiveS = 0, ExtendNegativeS = 0;
542 if ((ToAddPositiveS[0].
size() > 0) && (ToAddPositiveS[1].
size() > 0) &&
543 (ToAddPositiveS[2].
size() > 0)) {
544 for (
size_t n = 0;
n != 3; ++
n) {
545 int n1 = (
n + 1) % 3;
546 int n2 = (
n + 2) % 3;
548 if ((ToAddPositiveS[
n].back() <= ToAddPositiveS[n1].back()) &&
549 (ToAddPositiveS[
n].back() <= ToAddPositiveS[n2].back())) {
550 ExtendPositiveS = ToAddPositiveS[
n].back();
555 if ((ToAddNegativeS[0].
size() > 0) && (ToAddNegativeS[1].
size() > 0) &&
556 (ToAddNegativeS[2].
size() > 0)) {
557 for (
size_t n = 0;
n != 3; ++
n) {
558 int n1 = (
n + 1) % 3;
559 int n2 = (
n + 2) % 3;
560 if ((ToAddNegativeS[
n].back() >= ToAddNegativeS[n1].back()) &&
561 (ToAddNegativeS[
n].back() >= ToAddNegativeS[n2].back())) {
562 ExtendNegativeS = ToAddNegativeS[
n].back();
567 if (fabs(ExtendNegativeS) < 1.) ExtendNegativeS = -1.;
568 if (fabs(ExtendPositiveS) < 1.) ExtendPositiveS = 1.;
570 LengthRescale = (ExtendPositiveS - ExtendNegativeS) / 2.;
571 PositionShift = (ExtendPositiveS + ExtendNegativeS) / 2.;
573 for (
size_t n = 0;
n != 3; ++
n) {
574 pt[
n] += dir[
n] * PositionShift;
575 dir[
n] *= LengthRescale;
577 for (
size_t i = 0; i != ToAddPositiveS[
n].size(); ++i) {
578 if (ToAddPositiveS[
n].at(i) < ExtendPositiveS)
579 HitStatus[ToAddPositiveH[
n].at(i)] = 2;
581 HitStatus[ToAddPositiveH[
n].at(i)] = 0;
584 for (
size_t i = 0; i != ToAddNegativeS[
n].size(); ++i) {
585 if (ToAddNegativeS[
n].at(i) > ExtendNegativeS)
586 HitStatus[ToAddNegativeH[
n].at(i)] = 2;
588 HitStatus[ToAddNegativeH[
n].at(i)] = 0;
612 auto const [start,
end] = wireReadoutGeom.WireEndPoints(
geo::WireID{planeID, 0});
618 double HitWidth = HitXHigh - HitXLow;
620 double pt[3],
dir[3], err[3];
625 TVector3 sPt(pt[0], pt[1], pt[2]);
626 TVector3 sDir(dir[0], dir[1], dir[2]);
627 TVector3 hPt(HitX, start.Y(), start.Z());
628 TVector3 hDir(0, start.Y() -
end.Y(), start.Z() -
end.Z());
630 s = (sPt - hPt).Dot(hDir * (hDir.Dot(sDir)) - sDir * (hDir.Dot(hDir))) /
631 (hDir.Dot(hDir) * sDir.Dot(sDir) - pow(hDir.Dot(sDir), 2));
633 disp = fabs((sPt - hPt).Dot(sDir.Cross(hDir)) / (sDir.Cross(hDir)).Mag()) / HitWidth;
642 std::vector<recob::SpacePoint>
const& Points,
643 std::vector<char>& PointStatus,
644 std::vector<int>& PointsInRange,
652 std::vector<recob::SpacePoint> PointsUsed;
655 PointsInRange.clear();
658 TVector3 HighestZPoint;
659 bool NoPointFound =
true;
660 int counter = Points.size() - 1;
661 while ((NoPointFound ==
true) && (counter >= 0)) {
662 if (PointStatus[counter] == 0) {
663 HighestZPoint = TVector3(
664 Points.at(counter).XYZ()[0], Points.at(counter).XYZ()[1], Points.at(counter).XYZ()[2]);
665 NoPointFound =
false;
682 for (
int index = Points.size() - 1; index != -1; --index) {
683 if (PointStatus[index] == 0) {
686 if ((HighestZPoint[2] - Points.at(index).XYZ()[2]) < TwiceLength) {
687 double DistanceToHighZ = pow(pow(HighestZPoint[1] - Points.at(index).XYZ()[1], 2) +
688 pow(HighestZPoint[2] - Points.at(index).XYZ()[2], 2),
690 if (DistanceToHighZ < TwiceLength) {
691 PointsInRange.push_back(index);
692 PointsUsed.push_back(Points.at(index));
700 TVector3 SeedCenter(0, 0, 0);
701 TVector3 SeedDirection(0, 0, 0);
705 int NPoints = PointsInRange.size();
709 std::map<int, bool>
HitMap;
712 for (
unsigned int i = 0; i != PointsInRange.size(); i++) {
713 std::vector<art::PtrVector<recob::Hit>> HitsInThisCollection(3);
718 itHit != HitsThisSP.
end();
720 uint32_t Channel = (*itHit)->Channel();
724 for (
size_t iH = 0; iH != OrgHits[View][Channel].size(); ++iH) {
725 if (fabs(HitsFlat[OrgHits[View][Channel][iH]]->PeakTime() - (*itHit)->PeakTime()) < eta) {
726 HitMap[OrgHits[View][Channel][iH]] =
true;
732 for (
auto itH = HitMap.begin(); itH != HitMap.end(); ++itH) {
733 HitList.push_back(itH->first);
736 std::vector<double> ViewRMS;
737 std::vector<int> HitsPerView;
740 detProp, HitsFlat, HitList, SeedCenter, SeedDirection, ViewRMS, HitsPerView);
747 bool ThrowOutSeed =
false;
749 double PtArray[3], DirArray[3];
751 pow(pow(SeedDirection.Y(), 2) + pow(SeedDirection.Z(), 2), 0.5) / SeedDirection.Mag();
753 int nViewsWithHits(0);
755 for (
size_t n = 0;
n != 3; ++
n) {
757 PtArray[
n] = SeedCenter[
n];
758 if (HitsPerView[
n] > 0) nViewsWithHits++;
761 if (nViewsWithHits < 2 || (nViewsWithHits < 3 && !
fAllow2DSeeds)) ThrowOutSeed =
true;
768 if (
fMaxViewRMS.at(j) < ViewRMS.at(j)) { ThrowOutSeed =
true; }
785 std::vector<int>& HitsToUse,
788 std::vector<double>& ViewRMS,
789 std::vector<int>& N)
const 793 std::map<uint32_t, bool> HitsClaimed;
797 std::vector<std::vector<double>> HitTimes(3);
798 std::vector<std::vector<double>> HitWires(3);
799 std::vector<std::vector<double>> HitWidths(3);
800 std::vector<double> MeanWireCoord(3, 0);
801 std::vector<double> MeanTimeCoord(3, 0);
805 std::vector<double>
x(3, 0),
y(3, 0),
xx(3, 0), xy(3, 0), yy(3, 0), sig(3, 0);
807 for (
size_t i = 0; i != HitsToUse.size(); ++i) {
808 auto itHit = HitsFlat.
begin() + HitsToUse[i];
812 auto const hitView = (*itHit)->View();
824 double WireCoord = (*itHit)->WireID().Wire *
fPitches.at(ViewIndex);
825 double TimeCoord = detProp.
ConvertTicksToX((*itHit)->PeakTime(), ViewIndex, 0, 0);
826 double TimeUpper = detProp.
ConvertTicksToX((*itHit)->PeakTimePlusRMS(), ViewIndex, 0, 0);
827 double TimeLower = detProp.
ConvertTicksToX((*itHit)->PeakTimeMinusRMS(), ViewIndex, 0, 0);
828 double Width = fabs(0.5 * (TimeUpper - TimeLower));
829 double Width2 = pow(Width, 2);
831 HitWires.at(ViewIndex).push_back(WireCoord);
832 HitTimes.at(ViewIndex).push_back(TimeCoord);
833 HitWidths.at(ViewIndex).push_back(fabs(0.5 * (TimeUpper - TimeLower)));
835 MeanWireCoord.at(ViewIndex) += WireCoord;
836 MeanTimeCoord.at(ViewIndex) += TimeCoord;
839 x.at(ViewIndex) += WireCoord / Width2;
840 y.at(ViewIndex) += TimeCoord / Width2;
841 xy.at(ViewIndex) += (TimeCoord * WireCoord) / Width2;
842 xx.at(ViewIndex) += (WireCoord * WireCoord) / Width2;
843 yy.at(ViewIndex) += (TimeCoord * TimeCoord) / Width2;
844 sig.at(ViewIndex) += 1. / Width2;
850 std::vector<double> ViewGrad(3);
851 std::vector<double> ViewOffset(3);
853 for (
size_t n = 0;
n != 3; ++
n) {
854 MeanWireCoord[
n] /= N[
n];
855 MeanTimeCoord[
n] /= N[
n];
857 double BigN = 1000000;
858 double SmallN = 1. / BigN;
861 double Numerator = (
y[
n] / sig[
n] - xy[
n] /
x[
n]);
862 double Denominator = (
x[
n] / sig[
n] -
xx[
n] /
x[
n]);
863 if (fabs(Denominator) > SmallN)
864 ViewGrad.at(
n) = Numerator / Denominator;
869 ViewGrad[
n] = xy[
n] /
xx[
n];
873 ViewOffset.at(
n) = (
y[
n] - ViewGrad[
n] *
x[
n]) / sig[
n];
874 ViewRMS.at(
n) = pow((yy[
n] + pow(ViewGrad[
n], 2) *
xx[n] + pow(ViewOffset[n], 2) * sig[n] -
875 2 * ViewGrad[n] * xy[n] - 2 * ViewOffset[n] *
y[n] +
876 2 * ViewGrad[n] * ViewOffset[n] *
x[n]) /
880 if (ViewGrad.at(n) != 0) ViewRMS[n] *= sin(atan(1. / ViewGrad.at(n)));
883 for (
size_t n = 0;
n != 3; ++
n) {
884 size_t n1 = (
n + 1) % 3;
885 size_t n2 = (
n + 2) % 3;
886 if ((N[
n] <= N[n1]) && (N[
n] <= N[n2])) {
888 if (N[n1] < N[n2]) { std::swap(n1, n2); }
889 if ((N[n1] == 0) || (N[n2] == 0))
continue;
905 double TimeCoord = 0.5 * (MeanTimeCoord[n1] + MeanTimeCoord[n2]);
906 double WireCoordIn1 = (TimeCoord - ViewOffset[n1]) / ViewGrad[n1] +
fWireZeroOffset[n1];
907 double WireCoordIn2 = (TimeCoord - ViewOffset[n2]) / ViewGrad[n2] +
fWireZeroOffset[n2];
913 ViewRMS[
n] = -fabs(ViewRMS[
n]);
914 ViewRMS[n1] = fabs(ViewRMS[n1]);
915 ViewRMS[n2] = fabs(ViewRMS[n2]);
933 fPitches.at(0) = fabs(wireReadoutGeom.Plane({0, 0, geo::kU}).WirePitch());
934 fPitches.at(1) = fabs(wireReadoutGeom.Plane({0, 0, geo::kV}).WirePitch());
935 fPitches.at(2) = fabs(wireReadoutGeom.Plane({0, 0, geo::kW}).WirePitch());
938 fXDir = TVector3(1, 0, 0);
939 fYDir = TVector3(0, 1, 0);
940 fZDir = TVector3(0, 0, 1);
946 double xyzStart1[3], xyzStart2[3];
947 double xyzEnd1[3], xyzEnd2[3];
950 for (
unsigned int n = 0;
n != 3u; ++
n) {
951 wireReadoutGeom.WireEndPoints(
geo::WireID{0, 0,
n, 0}, xyzStart1, xyzEnd1);
952 wireReadoutGeom.WireEndPoints(
geo::WireID{0, 0,
n, 1}, xyzStart2, xyzEnd2);
954 TVector3(xyzEnd1[0] - xyzStart1[0], xyzEnd1[1] - xyzStart1[1], xyzEnd1[2] - xyzStart1[2])
958 xyzEnd2[0] - xyzEnd1[0], xyzEnd2[1] - xyzEnd1[1], xyzEnd2[2] - xyzEnd1[2])) < 0)
973 unsigned int StopAfter)
const 975 return FindSeeds(clockData, detProp, Hits, HitCatalogue, StopAfter);
985 unsigned int StopAfter)
const 988 std::vector<std::vector<recob::Seed>> ReturnVec;
995 <<
"Warning: SpacePointAlg is does not have three views enabled. This may cause unexpected " 996 "behaviour in the bezier tracker.";
1000 SortedHits.at(
geo::kU).begin();
1001 itU != SortedHits.at(
geo::kU).end();
1004 SortedHits.at(
geo::kV).begin();
1005 itV != SortedHits.at(
geo::kV).end();
1008 SortedHits.at(
geo::kW).begin();
1009 itW != SortedHits.at(
geo::kW).end();
1014 for (
size_t i = 0; i != itU->size(); ++i)
1015 HitsThisComboFlat.
push_back(itU->at(i));
1018 for (
size_t i = 0; i != itV->size(); ++i)
1019 HitsThisComboFlat.
push_back(itV->at(i));
1022 for (
size_t i = 0; i != itW->size(); ++i)
1023 HitsThisComboFlat.
push_back(itW->at(i));
1025 std::vector<art::PtrVector<recob::Hit>> CataloguedHits;
1027 std::vector<recob::Seed> Seeds =
1028 FindSeeds(clockData, detProp, HitsThisComboFlat, CataloguedHits, StopAfter);
1031 HitsPerSeed.push_back(CataloguedHits);
1032 ReturnVec.push_back(Seeds);
1035 CataloguedHits.clear();
1041 <<
" bailed during hit map lookup - have you enabled all 3 planes?";
1042 ReturnVec.push_back(std::vector<recob::Seed>());
std::vector< double > fPitches
float Length(const PFPStruct &pfp)
SeedFinderAlgorithm(const fhicl::ParameterSet &pset)
bool enableW() const noexcept
bool enableV() const noexcept
static std::string ViewName(View_t view)
Returns the name of the specified view.
enum geo::_plane_proj View_t
Enumerate the possible plane projections.
Declaration of signal hit object.
void GetPoint(double *Pt, double *Err) const
void CalculateGeometricalElements()
The data type to uniquely identify a Plane.
void GetHitDistAndProj(detinfo::DetectorPropertiesData const &detProp, recob::Seed const &ASeed, art::Ptr< recob::Hit > const &AHit, double &disp, double &s) const
cout<< "Opened file "<< fin<< " ixs= "<< ixs<< endl;if(ixs==0) hhh=(TH1F *) fff-> Get("h1")
geo::WireID const & WireID() const
Initial tdc tick for hit.
recob::Seed FindSeedAtEnd(detinfo::DetectorPropertiesData const &detProp, std::vector< recob::SpacePoint > const &, std::vector< char > &, std::vector< int > &, art::PtrVector< recob::Hit > const &HitsFlat, std::vector< std::vector< std::vector< int >>> &OrgHits) const
const art::PtrVector< recob::Hit > & getAssociatedHits(const recob::SpacePoint &spt) const
void ConsolidateSeed(detinfo::DetectorPropertiesData const &detProp, recob::Seed &TheSeed, art::PtrVector< recob::Hit > const &, std::vector< char > &HitStatus, std::vector< std::vector< std::vector< int >>> &OrgHits, bool Extend) const
std::map< int, art::Ptr< recob::Hit > > HitMap
decltype(auto) constexpr end(T &&obj)
ADL-aware version of std::end.
decltype(auto) constexpr size(T &&obj)
ADL-aware version of std::size.
auto vector(Vector const &v)
Returns a manipulator which will print the specified array.
typename data_t::const_iterator const_iterator
auto counter(T begin, T end)
Returns an object to iterate values from begin to end in a range-for loop.
void SetPoint(double *Pt, double *Err)
void push_back(Ptr< U > const &p)
void makeSpacePoints(detinfo::DetectorClocksData const &clockData, detinfo::DetectorPropertiesData const &detProp, const art::PtrVector< recob::Hit > &hits, std::vector< recob::SpacePoint > &spts) const
T get(std::string const &key) const
void GetCenterAndDirection(detinfo::DetectorPropertiesData const &detProp, art::PtrVector< recob::Hit > const &HitsFlat, std::vector< int > &HitsToUse, TVector3 &Center, TVector3 &Direction, std::vector< double > &ViewRMS, std::vector< int > &HitsPerView) const
reference at(size_type n)
The data type to uniquely identify a TPC.
PlaneID_t Plane
Index of the plane within its TPC.
std::vector< recob::Seed > FindSeeds(detinfo::DetectorClocksData const &clockData, detinfo::DetectorPropertiesData const &detProp, art::PtrVector< recob::Hit > const &HitsFlat, std::vector< art::PtrVector< recob::Hit >> &CataloguedHits, unsigned int StopAfter) const
float PeakTimeMinusRMS(float sigmas=+1.) const
Returns a time sigmas RMS away from the peak time.
std::vector< recob::Seed > GetSeedsFromUnSortedHits(detinfo::DetectorClocksData const &clockData, detinfo::DetectorPropertiesData const &detProp, art::PtrVector< recob::Hit > const &, std::vector< art::PtrVector< recob::Hit >> &, unsigned int StopAfter=0) const
double ConvertTicksToX(double ticks, int p, int t, int c) const
cet::coded_exception< errors::ErrorCodes, ExceptionDetail::translate > Exception
float PeakTime() const
Time of the signal peak, in tick units.
Encapsulate the construction of a single detector plane .
double dist(const TReal *x, const TReal *y, const unsigned int dimension)
Contains all timing reference information for the detector.
bool enableU() const noexcept
std::vector< TVector3 > fWireDir
MaybeLogger_< ELseverityLevel::ELsev_warning, false > LogWarning
std::vector< TVector3 > fPitchDir
Planes which measure W (third view for Bo, MicroBooNE, etc).
float PeakTimePlusRMS(float sigmas=+1.) const
Returns a time sigmas RMS away from the peak time.
Algorithm for generating space points from hits.
std::set< art::Ptr< recob::Hit > > HitList
void SetValidity(bool Validity)
std::vector< std::vector< recob::Seed > > GetSeedsFromSortedHits(detinfo::DetectorClocksData const &clockData, detinfo::DetectorPropertiesData const &detProp, std::vector< std::vector< art::PtrVector< recob::Hit >>> const &SortedHits, std::vector< std::vector< art::PtrVector< recob::Hit >>> &HitsPerSeed, unsigned int StopAfter=0) const
void GetDirection(double *Dir, double *Err) const
std::vector< double > fWireZeroOffset
art framework interface to geometry description
void SetDirection(double *Dir, double *Err)
std::vector< double > fMaxViewRMS