diff --git a/Tracking/exampleConfigs/test_reduced_config.py b/Tracking/exampleConfigs/test_reduced_config.py index 03ab71118..96d87fbbf 100644 --- a/Tracking/exampleConfigs/test_reduced_config.py +++ b/Tracking/exampleConfigs/test_reduced_config.py @@ -26,11 +26,11 @@ p.termLogLevel = 0 -p.maxEvents = 1 +p.maxEvents = 100 p.run = 200 -p.histogramFile = f'hist_reducedAlgoV1Test_seederONLY.root' -p.outputFiles = [f'events_1_reducedAlgoV1Test_seederONLY.root'] +p.histogramFile = f'hist_reducedAlgoV1Test.root' +p.outputFiles = [f'events_1_reducedAlgoV1Test.root'] import LDMX.Ecal.EcalGeometry import LDMX.Ecal.ecal_hardcoded_conditions diff --git a/Tracking/include/Tracking/Event/ReducedTrack.h b/Tracking/include/Tracking/Event/ReducedTrack.h index 8654640c6..0a9fa860b 100644 --- a/Tracking/include/Tracking/Event/ReducedTrack.h +++ b/Tracking/include/Tracking/Event/ReducedTrack.h @@ -73,6 +73,12 @@ namespace ldmx { void setDistancetoEcalRecHit(double distance) {distance_to_Ecal_ = distance;} double getDistanceToRecHit() const { return distance_to_Ecal_; } + const std::vector>& getAllSensorPoints() const { return totalSensor_;}; + + void setAllSensorPoints(const std::vector>& sensorPoints) { + totalSensor_ = sensorPoints; + } + void setFirstSensorPosition(const std::array& firstSensor) { firstSensor_ = firstSensor; } @@ -131,6 +137,7 @@ namespace ldmx { double by_; double distance_to_Ecal_; + std::vector> totalSensor_; //! transient for ROOT std::array firstSensor_; std::array secondSensor_; std::array ecalRecHit_; diff --git a/Tracking/include/Tracking/Reco/ReducedSeedFinder.h b/Tracking/include/Tracking/Reco/ReducedSeedFinder.h index bad24611e..2059557d1 100644 --- a/Tracking/include/Tracking/Reco/ReducedSeedFinder.h +++ b/Tracking/include/Tracking/Reco/ReducedSeedFinder.h @@ -78,9 +78,9 @@ class ReducedSeedFinder : public TrackingGeometryUser { void produce(framework::Event& event) override; protected: - ldmx::ReducedTrack SeedTracker(const std::array recoilOne, const std::array recoilTwo, const std::array ecalOne); + ldmx::ReducedTrack SeedTracker(const std::array recoilOne, const std::array recoilTwo, const std::array ecalOne, const std::vector> allPoints); - std::pair>, std::vector>> combineMultiGlobalHits(const std::vector> &hitCollection); + std::tuple>, std::vector>, std::vector>> combineMultiGlobalHits(const std::vector> &hitCollection); std::vector> weightedAverage(const std::vector> &layer1, const std::vector> &layer2); std::tuple fit3DLine(const std::array &firstRecoil, const std::array &secondRecoil, const std::array &ECal); diff --git a/Tracking/python/eventDisplay.py b/Tracking/python/eventDisplay.py index b8841407f..ab6eeac55 100644 --- a/Tracking/python/eventDisplay.py +++ b/Tracking/python/eventDisplay.py @@ -20,8 +20,8 @@ def trackPlotter(tree, event_number, tag, save): recoilSimHits = addBranch(tree, 'SimTrackerHit', 'RecoilSimHits_{}'.format(tag)) digiRecoil = addBranch(tree, 'Measurement', 'DigiRecoilSimHits_{}'.format(tag)) ecalRecHit = addBranch(tree, 'EcalHit', 'EcalRecHits_{}'.format(tag)) - #recoil = addBranch(tree, 'ReducedTrack', 'ReducedTracks_{}'.format(tag)) - recoil = addBranch(tree, 'ReducedTrack', 'ReducedSeedTracks_{}'.format(tag)) + recoil = addBranch(tree, 'ReducedTrack', 'ReducedTracks_{}'.format(tag)) + #recoil = addBranch(tree, 'ReducedTrack', 'ReducedSeedTracks_{}'.format(tag)) tree.GetEntry(event_number) diff --git a/Tracking/src/Tracking/Reco/ReducedSeedFinder.cxx b/Tracking/src/Tracking/Reco/ReducedSeedFinder.cxx index a1c8d2e03..a51ec2121 100644 --- a/Tracking/src/Tracking/Reco/ReducedSeedFinder.cxx +++ b/Tracking/src/Tracking/Reco/ReducedSeedFinder.cxx @@ -63,12 +63,12 @@ void ReducedSeedFinder::produce(framework::Event& event) { return; } - auto [firstSensor, secondSensor] = combineMultiGlobalHits(digiPoints); + auto [firstSensor, secondSensor, allHits_noEDEP] = combineMultiGlobalHits(digiPoints); for (const auto& firstPoint : firstSensor) { for (const auto& secondPoint : secondSensor) { for (const auto& recHit : firstLayerEcalRecHits) { - ldmx::ReducedTrack seedTrack = SeedTracker(firstPoint, secondPoint, recHit); + ldmx::ReducedTrack seedTrack = SeedTracker(firstPoint, secondPoint, recHit, allHits_noEDEP); if (seedTrack.getChi2() > 0.0) { reduced_seed_tracks.push_back(seedTrack); } @@ -89,7 +89,7 @@ void ReducedSeedFinder::produce(framework::Event& event) { } //produce -ldmx::ReducedTrack ReducedSeedFinder::SeedTracker(const std::array recoilOne, const std::array recoilTwo, const std::array ecalOne) { +ldmx::ReducedTrack ReducedSeedFinder::SeedTracker(const std::array recoilOne, const std::array recoilTwo, const std::array ecalOne, std::vector> allPoints) { auto [ax, bx, ay, by] = fit3DLine(recoilOne, recoilTwo, ecalOne); std::array tempExtrapolatedPoint = {ecalOne[0], ax * ecalOne[0] + bx, ay * ecalOne[0] + by}; @@ -103,6 +103,7 @@ ldmx::ReducedTrack ReducedSeedFinder::SeedTracker(const std::array re trk.setAY(ay); trk.setBY(by); + trk.setAllSensorPoints(allPoints); trk.setFirstSensorPosition(recoilOne); trk.setSecondSensorPosition(recoilTwo); trk.setFirstLayerEcalRecHit(ecalOne); @@ -134,11 +135,16 @@ void ReducedSeedFinder::onProcessEnd() { //HAVE TO FIX THESE VALUES // ldmx_log(info) << " nfailthetacut=" << nfailtheta_; } //onProcessEnd -std::pair>, std::vector>> ReducedSeedFinder::combineMultiGlobalHits(const std::vector>& hitCollection) { +std::tuple>, + std::vector>, + std::vector>> +ReducedSeedFinder::combineMultiGlobalHits(const std::vector>& hitCollection) { std::vector> layer1, layer2, layer3, layer4; + std::vector> allPoints; // Split hits into layers based on z position for (const auto& point : hitCollection) { + allPoints.push_back({point[0], point[1], point[2]}); if (point[0] < 12) layer1.push_back(point); else if (point[0] < 20) layer2.push_back(point); else if (point[0] < 28) layer3.push_back(point); @@ -149,7 +155,7 @@ std::pair>, std::vector> auto firstSensorMergedHits = weightedAverage(layer1, layer2); auto secondSensorMergedHits = weightedAverage(layer3, layer4); - return {firstSensorMergedHits, secondSensorMergedHits}; + return {firstSensorMergedHits, secondSensorMergedHits, allPoints}; } //combineMultiGlobalHits std::vector> ReducedSeedFinder::weightedAverage(const std::vector>& layer1, const std::vector>& layer2) { diff --git a/Tracking/src/Tracking/Reco/ReducedTrackFinder.cxx b/Tracking/src/Tracking/Reco/ReducedTrackFinder.cxx index 1e8d014d5..a563d6c19 100644 --- a/Tracking/src/Tracking/Reco/ReducedTrackFinder.cxx +++ b/Tracking/src/Tracking/Reco/ReducedTrackFinder.cxx @@ -97,18 +97,15 @@ std::vector ReducedTrackFinder::findTracks(const std::vector ldmx::ReducedTrack bestSeed = *bestSeedIt; bestTrack.push_back(bestSeed); - // Remove any seeds that share Recoil sensor points with the best seed - auto sensorPoint1 = bestSeed.getFirstSensorPosition(); - auto sensorPoint2 = bestSeed.getSecondSensorPosition(); - seedsWithSameRecHit.erase( - std::remove_if(seedsWithSameRecHit.begin(), seedsWithSameRecHit.end(), - [&](const ldmx::ReducedTrack& seed) { - auto recoilPoint1 = seed.getFirstSensorPosition(); - auto recoilPoint2 = seed.getSecondSensorPosition(); - return (sensorPoint1 == recoilPoint1 || sensorPoint1 == recoilPoint2 || - sensorPoint2 == recoilPoint1 || sensorPoint2 == recoilPoint2); - }), - seedsWithSameRecHit.end()); + // Remove any seeds that share any sensor point with the best seed + auto sensorPoints = bestSeed.getAllSensorPoints(); + seedsWithSameRecHit.erase(std::remove_if(seedsWithSameRecHit.begin(), seedsWithSameRecHit.end(), [&](const ldmx::ReducedTrack& seed) { + auto seedPoints = seed.getAllSensorPoints(); + // Check if any sensor point in the seed overlaps with the best seed's points + return std::any_of(sensorPoints.begin(), sensorPoints.end(), [&](const auto& point) { + return std::find(seedPoints.begin(), seedPoints.end(), point) != seedPoints.end(); }); + }), + seedsWithSameRecHit.end()); } return bestTrack;