Make tracksters Associators reproducible

Author: felicepantaleo

SimCalorimetry/HGCalAssociatorProducers/plugins/AllTracksterToSimTracksterAssociatorsByHitsProducer.cc

@@ -259,14 +259,15 @@ void AllTracksterToSimTracksterAssociatorsByHitsProducer::produce(edm::StreamID,
           const auto& recHit = rechitManager[hitIndex];
           float recoFraction = recoTracksterHitsAndFractions[i].fraction();
           float squaredRecoFraction = recoFraction * recoFraction;
-          float squaredRecHitEnergy = recHit.energy() * recHit.energy();
-          float recoSharedEnergy = recHit.energy() * recoFraction;
+          float rechitEnergy = recHit.energy();
+          float squaredRecHitEnergy = rechitEnergy * rechitEnergy;
+          float recoSharedEnergy = rechitEnergy * recoFraction;
           const auto& simTracksterVec = hitToAssociatedSimTracksterMap[i];
           for (const auto& simTracksterElement : simTracksterVec) {
             auto simTracksterIndex = simTracksterElement.index();
             auto simFraction = simTracksterElement.fraction();
             edm::Ref<std::vector<ticl::Trackster>> simTracksterRef(simTrackstersHandle, simTracksterIndex);
-            float sharedEnergy = std::min(simFraction * recHit.energy(), recoSharedEnergy);
+            float sharedEnergy = std::min(simFraction * rechitEnergy, recoSharedEnergy);
             float squaredFraction =
                 std::min(squaredRecoFraction, (recoFraction - simFraction) * (recoFraction - simFraction));
             float score = invDenominator * squaredFraction * squaredRecHitEnergy;
@@ -297,15 +298,16 @@ void AllTracksterToSimTracksterAssociatorsByHitsProducer::produce(edm::StreamID,
                         : std::find_if(hitToSimClusterMap[hitIndex].begin(),
                                        hitToSimClusterMap[hitIndex].end(),
                                        [simObjectIndex](const auto& pair) { return pair.index() == simObjectIndex; });
-          if (it != hitToCaloParticleMap[hitIndex].end() and it != hitToSimClusterMap[hitIndex].end()) {
+          if ((isSimTracksterFromCP and it != hitToCaloParticleMap[hitIndex].end()) or
+              (!isSimTracksterFromCP and it != hitToSimClusterMap[hitIndex].end())) {
             simFractions[i] = it->fraction();
           }
           float simFraction = simFractions[i];
           const auto& recHit = rechitManager[hitIndex];
+          float rechitEnergy = recHit.energy();
           float squaredSimFraction = simFraction * simFraction;
-          float squaredRecHitEnergy = recHit.energy() * recHit.energy();
+          float squaredRecHitEnergy = rechitEnergy * rechitEnergy;
           simToRecoScoresDenominator += squaredSimFraction * squaredRecHitEnergy;
-
           const auto& hitToRecoTracksterVec = hitToTracksterMap[hitIndex];
           for (const auto& recoTracksterElement : hitToRecoTracksterVec) {
             unsigned int recoTracksterIndex = recoTracksterElement.index();
@@ -333,22 +335,22 @@ void AllTracksterToSimTracksterAssociatorsByHitsProducer::produce(edm::StreamID,
           }
         }
 
+        assert(simToRecoScoresDenominator > 0.f);
         const float invDenominator = 1.f / simToRecoScoresDenominator;
-
         for (unsigned int i = 0; i < simTracksterHitsAndFractions.size(); ++i) {
           const auto& hitIndex = simTracksterHitsAndFractions[i].index();
           float simFraction = simFractions[i];
           const auto& recHit = rechitManager[hitIndex];
+          float rechitEnergy = recHit.energy();
           float squaredSimFraction = simFraction * simFraction;
-          float squaredRecHitEnergy = recHit.energy() * recHit.energy();
-          float simSharedEnergy = recHit.energy() * simFraction;
-
+          float squaredRecHitEnergy = rechitEnergy * rechitEnergy;
+          float simSharedEnergy = rechitEnergy * simFraction;
           const auto& hitToRecoTracksterVec = hitToAssociatedRecoTracksterMap[i];
           for (const auto& recoTracksterElement : hitToRecoTracksterVec) {
             auto recoTracksterIndex = recoTracksterElement.index();
             float recoFraction = recoTracksterElement.fraction();
             edm::Ref<std::vector<ticl::Trackster>> recoTracksterRef(recoTrackstersHandle, recoTracksterIndex);
-            float sharedEnergy = std::min(recoFraction * recHit.energy(), simSharedEnergy);
+            float sharedEnergy = std::min(recoFraction * rechitEnergy, simSharedEnergy);
             float squaredFraction =
                 std::min(squaredSimFraction, (recoFraction - simFraction) * (recoFraction - simFraction));
             float score = invDenominator * squaredFraction * squaredRecHitEnergy;
@@ -357,9 +359,16 @@ void AllTracksterToSimTracksterAssociatorsByHitsProducer::produce(edm::StreamID,
         }
       }
 
+      auto sortingFunc = [](const auto& a, const auto& b) {
+        if (a.score() != b.score())
+          return a.score() < b.score();
+        else
+          return a.index() < b.index();
+      };
+
       // Sort the maps by score in ascending order
-      tracksterToSimTracksterMap->sort([](const auto& a, const auto& b) { return a.score() < b.score(); });
-      simTracksterToTracksterMap->sort([](const auto& a, const auto& b) { return a.score() < b.score(); });
+      tracksterToSimTracksterMap->sort(sortingFunc);
+      simTracksterToTracksterMap->sort(sortingFunc);
 
       // After populating the maps, store them in the event
       iEvent.put(std::move(tracksterToSimTracksterMap), tracksterToken.first + "To" + simTracksterToken.first);

SimCalorimetry/HGCalAssociatorProducers/plugins/AllTracksterToSimTracksterAssociatorsByLCsProducer.cc

@@ -244,7 +244,7 @@ void AllTracksterToSimTracksterAssociatorsByLCsProducer::produce(edm::StreamID,
             }
           }
         }
-
+        assert(simToRecoScoresDenominator > 0.f);
         const float invDenominator = 1.f / simToRecoScoresDenominator;
 
         for (unsigned int i = 0; i < layerClustersIds.size(); ++i) {
@@ -269,9 +269,16 @@ void AllTracksterToSimTracksterAssociatorsByLCsProducer::produce(edm::StreamID,
           }
         }
       }
+      auto sortingFunc = [](const auto& a, const auto& b) {
+        if (a.score() != b.score())
+          return a.score() < b.score();
+        else
+          return a.index() < b.index();
+      };
+
       // Sort the maps by score in ascending order
-      tracksterToSimTracksterMap->sort([](const auto& a, const auto& b) { return a.score() < b.score(); });
-      simTracksterToTracksterMap->sort([](const auto& a, const auto& b) { return a.score() < b.score(); });
+      tracksterToSimTracksterMap->sort(sortingFunc);
+      simTracksterToTracksterMap->sort(sortingFunc);
 
       // After populating the maps, store them in the event
       iEvent.put(std::move(tracksterToSimTracksterMap), tracksterToken.first + "To" + simTracksterToken.first);