diff --git a/SimCalorimetry/HGCalAssociatorProducers/plugins/AllTracksterToSimTracksterAssociatorsByHitsProducer.cc b/SimCalorimetry/HGCalAssociatorProducers/plugins/AllTracksterToSimTracksterAssociatorsByHitsProducer.cc index 20d81fa1e19aa..90d4cffb866f5 100644 --- a/SimCalorimetry/HGCalAssociatorProducers/plugins/AllTracksterToSimTracksterAssociatorsByHitsProducer.cc +++ b/SimCalorimetry/HGCalAssociatorProducers/plugins/AllTracksterToSimTracksterAssociatorsByHitsProducer.cc @@ -54,7 +54,7 @@ AllTracksterToSimTracksterAssociatorsByHitsProducer::AllTracksterToSimTracksterA const auto& tracksterCollections = pset.getParameter>("tracksterCollections"); for (const auto& tag : tracksterCollections) { std::string label = tag.label(); - if (tag.instance() != "") { + if (!tag.instance().empty()) { label += tag.instance(); } tracksterCollectionTokens_.emplace_back(label, consumes>(tag)); @@ -69,7 +69,7 @@ AllTracksterToSimTracksterAssociatorsByHitsProducer::AllTracksterToSimTracksterA const auto& simTracksterCollections = pset.getParameter>("simTracksterCollections"); for (const auto& tag : simTracksterCollections) { std::string label = tag.label(); - if (tag.instance() != "") { + if (!tag.instance().empty()) { label += tag.instance(); } simTracksterCollectionTokens_.emplace_back(label, consumes>(tag)); @@ -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> 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> 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); diff --git a/SimCalorimetry/HGCalAssociatorProducers/plugins/AllTracksterToSimTracksterAssociatorsByLCsProducer.cc b/SimCalorimetry/HGCalAssociatorProducers/plugins/AllTracksterToSimTracksterAssociatorsByLCsProducer.cc index 02992b166ed43..358cb535fa601 100644 --- a/SimCalorimetry/HGCalAssociatorProducers/plugins/AllTracksterToSimTracksterAssociatorsByLCsProducer.cc +++ b/SimCalorimetry/HGCalAssociatorProducers/plugins/AllTracksterToSimTracksterAssociatorsByLCsProducer.cc @@ -43,7 +43,7 @@ AllTracksterToSimTracksterAssociatorsByLCsProducer::AllTracksterToSimTracksterAs const auto& tracksterCollections = pset.getParameter>("tracksterCollections"); for (const auto& tag : tracksterCollections) { std::string label = tag.label(); - if (tag.instance() != "") { + if (!tag.instance().empty()) { label += tag.instance(); } tracksterCollectionTokens_.emplace_back(label, consumes>(tag)); @@ -57,7 +57,7 @@ AllTracksterToSimTracksterAssociatorsByLCsProducer::AllTracksterToSimTracksterAs const auto& simTracksterCollections = pset.getParameter>("simTracksterCollections"); for (const auto& tag : simTracksterCollections) { std::string label = tag.label(); - if (tag.instance() != "") { + if (!tag.instance().empty()) { label += tag.instance(); } simTracksterCollectionTokens_.emplace_back(label, consumes>(tag)); @@ -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); diff --git a/Validation/HGCalValidation/interface/HGCalValidator.h b/Validation/HGCalValidation/interface/HGCalValidator.h index 56d83bd0916c0..08a6f094c2bb1 100644 --- a/Validation/HGCalValidation/interface/HGCalValidator.h +++ b/Validation/HGCalValidation/interface/HGCalValidator.h @@ -32,6 +32,8 @@ #include "SimDataFormats/Associations/interface/LayerClusterToSimClusterAssociator.h" #include "SimDataFormats/Associations/interface/TICLAssociationMap.h" +#include "CommonTools/RecoAlgos/interface/MultiVectorManager.h" + class PileupSummaryInfo; struct HGCalValidatorHistograms { @@ -65,7 +67,7 @@ class HGCalValidator : public DQMGlobalEDAnalyzer { std::vector& selected_cPeff, unsigned int layers, std::unordered_map const&, - std::vector const& hits) const; + MultiVectorManager const& hits) const; static void fillDescriptions(edm::ConfigurationDescriptions& descriptions); @@ -111,7 +113,7 @@ class HGCalValidator : public DQMGlobalEDAnalyzer { edm::EDGetTokenT associatorMapRtSim; std::unique_ptr histoProducerAlgo_; std::vector hits_label_; - std::vector> hits_token_; + std::vector> hits_tokens_; std::unique_ptr candidateVal_; std::vector> tracksterToTracksterAssociatorsTokens_; std::vector> tracksterToTracksterByHitsAssociatorsTokens_; diff --git a/Validation/HGCalValidation/interface/HGVHistoProducerAlgo.h b/Validation/HGCalValidation/interface/HGVHistoProducerAlgo.h index fbcc6934ef686..e1636f7ae3959 100644 --- a/Validation/HGCalValidation/interface/HGVHistoProducerAlgo.h +++ b/Validation/HGCalValidation/interface/HGVHistoProducerAlgo.h @@ -22,6 +22,8 @@ #include "DataFormats/HGCalReco/interface/Trackster.h" #include "RecoLocalCalo/HGCalRecAlgos/interface/RecHitTools.h" +#include "CommonTools/RecoAlgos/interface/MultiVectorManager.h" + #include "SimDataFormats/CaloAnalysis/interface/CaloParticle.h" #include "SimDataFormats/Vertex/interface/SimVertex.h" #include "RecoLocalCalo/HGCalRecProducers/interface/HGCalClusteringAlgoBase.h" @@ -280,7 +282,7 @@ class HGVHistoProducerAlgo { unsigned int layers, const ticl::RecoToSimCollection& recSimColl, const ticl::SimToRecoCollection& simRecColl, - std::vector const& hits) const; + MultiVectorManager const& hits) const; void layerClusters_to_SimClusters(const Histograms& histograms, const int count, edm::Handle clusterHandle, @@ -293,7 +295,7 @@ class HGVHistoProducerAlgo { unsigned int layers, const ticl::RecoToSimCollectionWithSimClusters& recSimColl, const ticl::SimToRecoCollectionWithSimClusters& simRecColl, - std::vector const& hits) const; + MultiVectorManager const& hits) const; void tracksters_to_SimTracksters_fp(const Histograms& histograms, const int count, @@ -312,7 +314,7 @@ class HGVHistoProducerAlgo { std::vector const& simVertices, unsigned int layers, std::unordered_map const&, - std::vector const& hits) const; + MultiVectorManager const& hits) const; void fill_generic_cluster_histos(const Histograms& histograms, const int count, edm::Handle clusterHandle, @@ -327,7 +329,7 @@ class HGVHistoProducerAlgo { std::vector thicknesses, const ticl::RecoToSimCollection& recSimColl, const ticl::SimToRecoCollection& simRecColl, - std::vector const& hits) const; + MultiVectorManager const& hits) const; void fill_simCluster_histos(const Histograms& histograms, std::vector const& simClusters, unsigned int layers, @@ -344,7 +346,7 @@ class HGVHistoProducerAlgo { unsigned int layers, const ticl::RecoToSimCollectionWithSimClusters& recSimColl, const ticl::SimToRecoCollectionWithSimClusters& simRecColl, - std::vector const& hits) const; + MultiVectorManager const& hits) const; void fill_cluster_histos(const Histograms& histograms, const int count, const reco::CaloCluster& cluster) const; void fill_trackster_histos(const Histograms& histograms, const int count, @@ -360,7 +362,7 @@ class HGVHistoProducerAlgo { std::vector const& cPSelectedIndices, std::unordered_map const& hitMap, unsigned int layers, - std::vector const& hits, + MultiVectorManager const& hits, bool mapsFound, const edm::Handle& trackstersToSimTrackstersByLCsMapH, const edm::Handle& simTrackstersToTrackstersByLCsMapH, @@ -378,7 +380,7 @@ class HGVHistoProducerAlgo { DetId findmaxhit(const reco::CaloCluster& cluster, std::unordered_map const&, - std::vector const& hits) const; + MultiVectorManager const& hits) const; struct detIdInfoInCluster { bool operator==(const detIdInfoInCluster& o) const { return clusterId == o.clusterId; }; diff --git a/Validation/HGCalValidation/plugins/HGCalValidator.cc b/Validation/HGCalValidation/plugins/HGCalValidator.cc index 0a45de66c27a3..4683112dd5818 100644 --- a/Validation/HGCalValidation/plugins/HGCalValidator.cc +++ b/Validation/HGCalValidation/plugins/HGCalValidator.cc @@ -1,3 +1,5 @@ +#include + #include "Validation/HGCalValidation/interface/HGCalValidator.h" #include "SimCalorimetry/HGCalAssociatorProducers/interface/AssociatorTools.h" @@ -24,7 +26,10 @@ namespace { const auto recoTrackstersProductId = tracksterHandle.id(); const auto simTrackstersProductId = simTracksterHandle.id(); const auto simTrackstersFromCPsProductId = simTracksterFromCPHandle.id(); - + if (recoTrackstersProductId == simTrackstersProductId or recoTrackstersProductId == simTrackstersFromCPsProductId) { + edm::LogInfo("MissingProduct") << "no SimTrackster to Simtrackster map available."; + return false; + } for (const auto& handle : tracksterToTracksterMapsHandles) { const auto& firstID = handle->getCollectionIDs().first.id(); const auto& secondID = handle->getCollectionIDs().second.id(); @@ -100,7 +105,7 @@ HGCalValidator::HGCalValidator(const edm::ParameterSet& pset) const edm::InputTag& label_cp_fake_tag = pset.getParameter("label_cp_fake"); for (auto& label : hits_label_) { - hits_token_.push_back(consumes(label)); + hits_tokens_.push_back(consumes(label)); } label_cp_effic = consumes>(label_cp_effic_tag); label_cp_fake = consumes>(label_cp_fake_tag); @@ -178,7 +183,7 @@ HGCalValidator::HGCalValidator(const edm::ParameterSet& pset) pset.getParameter("notConvertedOnlyCP"), pset.getParameter>("pdgIdCP")); - tools_.reset(new hgcal::RecHitTools()); + tools_ = std::make_shared(); particles_to_monitor_ = pset.getParameter>("pdgIdCP"); totallayers_to_monitor_ = pset.getParameter("totallayers_to_monitor"); @@ -339,7 +344,7 @@ void HGCalValidator::cpParametersAndSelection(const Histograms& histograms, std::vector& selected_cPeff, unsigned int layers, std::unordered_map const& hitMap, - std::vector const& hits) const { + MultiVectorManager const& hits) const { selected_cPeff.reserve(cPeff.size()); size_t j = 0; @@ -387,7 +392,7 @@ void HGCalValidator::dqmAnalyze(const edm::Event& event, edm::Handle>> simTrackstersMapHandle; event.getByToken(simTrackstersMap_, simTrackstersMapHandle); - const std::map> cpToSc_SimTrackstersMap = *simTrackstersMapHandle; + const std::map>& cpToSc_SimTrackstersMap = *simTrackstersMapHandle; edm::ESHandle geom = setup.getHandle(caloGeomToken_); tools_->setGeometry(*geom); @@ -395,20 +400,20 @@ void HGCalValidator::dqmAnalyze(const edm::Event& event, edm::Handle simtorecoCollectionH; event.getByToken(associatorMapStR, simtorecoCollectionH); - auto simRecColl = *simtorecoCollectionH; + const auto& simRecColl = *simtorecoCollectionH; edm::Handle recotosimCollectionH; event.getByToken(associatorMapRtS, recotosimCollectionH); - auto recSimColl = *recotosimCollectionH; + const auto& recSimColl = *recotosimCollectionH; edm::Handle> hitMapHandle; event.getByToken(hitMap_, hitMapHandle); - const std::unordered_map* hitMap = &*hitMapHandle; + const std::unordered_map& hitMap = *hitMapHandle; - std::vector hits; - for (auto& token : hits_token_) { - edm::Handle hitsHandle; + MultiVectorManager rechitManager; + for (const auto& token : hits_tokens_) { + Handle hitsHandle; event.getByToken(token, hitsHandle); - hits.insert(hits.end(), (*hitsHandle).begin(), (*hitsHandle).end()); + rechitManager.addVector(*hitsHandle); } //Some general info on layers etc. @@ -428,7 +433,7 @@ void HGCalValidator::dqmAnalyze(const edm::Event& event, LogTrace("HGCalValidator") << "\n# of CaloParticles: " << caloParticles.size() << "\n" << std::endl; std::vector selected_cPeff; cpParametersAndSelection( - histograms, caloParticles, simVertices, selected_cPeff, totallayers_to_monitor_, *hitMap, hits); + histograms, caloParticles, simVertices, selected_cPeff, totallayers_to_monitor_, hitMap, rechitManager); //get collections from the event //simClusters @@ -500,11 +505,11 @@ void HGCalValidator::dqmAnalyze(const edm::Event& event, simClusters, sCIndices, inputClusterMask, - *hitMap, + hitMap, totallayers_to_monitor_, recSimColl, simRecColl, - hits); + rechitManager); //General Info on simClusters LogTrace("HGCalValidator") << "\n# of SimClusters: " << nSimClusters @@ -525,13 +530,13 @@ void HGCalValidator::dqmAnalyze(const edm::Event& event, caloParticles, cPIndices, selected_cPeff, - *hitMap, + hitMap, cumulative_material_budget, totallayers_to_monitor_, thicknesses_to_monitor_, recSimColl, simRecColl, - hits); + rechitManager); for (unsigned int layerclusterIndex = 0; layerclusterIndex < clusters.size(); layerclusterIndex++) { histoProducerAlgo_->fill_cluster_histos(histograms.histoProducerAlgo, w, clusters[layerclusterIndex]); @@ -550,8 +555,6 @@ void HGCalValidator::dqmAnalyze(const edm::Event& event, edm::Handle tracksterHandle; event.getByToken(label_tstTokens[wml], tracksterHandle); const ticl::TracksterCollection& tracksters = *tracksterHandle; - if (tracksterHandle.id() == simTracksterHandle.id() or tracksterHandle.id() == simTracksterFromCPHandle.id()) - continue; edm::Handle trackstersToSimTrackstersMapH, simTrackstersToTrackstersMapH, trackstersToSimTrackstersFromCPsMapH, simTrackstersFromCPsToTrackstersMapH, trackstersToSimTrackstersByHitsMapH, simTrackstersToTrackstersByHitsMapH, @@ -587,9 +590,9 @@ void HGCalValidator::dqmAnalyze(const edm::Event& event, caloParticles, cPIndices, selected_cPeff, - *hitMap, + hitMap, totallayers_to_monitor_, - hits, + rechitManager, mapsFound, trackstersToSimTrackstersMapH, simTrackstersToTrackstersMapH, diff --git a/Validation/HGCalValidation/src/HGVHistoProducerAlgo.cc b/Validation/HGCalValidation/src/HGVHistoProducerAlgo.cc index 2917ae9b759f0..677cd9942ae76 100644 --- a/Validation/HGCalValidation/src/HGVHistoProducerAlgo.cc +++ b/Validation/HGCalValidation/src/HGVHistoProducerAlgo.cc @@ -1453,7 +1453,7 @@ void HGVHistoProducerAlgo::fill_caloparticle_histos(const Histograms& histograms std::vector const& simVertices, unsigned int layers, std::unordered_map const& hitMap, - std::vector const& hits) const { + MultiVectorManager const& hits) const { const auto eta = getEta(caloParticle.eta()); if (histograms.h_caloparticle_eta.count(pdgid)) { histograms.h_caloparticle_eta.at(pdgid)->Fill(eta); @@ -1715,7 +1715,7 @@ void HGVHistoProducerAlgo::HGVHistoProducerAlgo::fill_simClusterAssociation_hist unsigned int layers, const ticl::RecoToSimCollectionWithSimClusters& scsInLayerClusterMap, const ticl::SimToRecoCollectionWithSimClusters& lcsInSimClusterMap, - std::vector const& hits) const { + MultiVectorManager const& hits) const { //Each event to be treated as two events: an event in +ve endcap, //plus another event in -ve endcap. In this spirit there will be //a layer variable (layerid) that maps the layers in : @@ -1757,7 +1757,7 @@ void HGVHistoProducerAlgo::layerClusters_to_CaloParticles(const Histograms& hist unsigned int layers, const ticl::RecoToSimCollection& cpsInLayerClusterMap, const ticl::SimToRecoCollection& cPOnLayerMap, - std::vector const& hits) const { + MultiVectorManager const& hits) const { const auto nLayerClusters = clusters.size(); std::unordered_map> detIdToCaloParticleId_Map; @@ -2018,7 +2018,7 @@ void HGVHistoProducerAlgo::layerClusters_to_SimClusters( unsigned int layers, const ticl::RecoToSimCollectionWithSimClusters& scsInLayerClusterMap, const ticl::SimToRecoCollectionWithSimClusters& lcsInSimClusterMap, - std::vector const& hits) const { + MultiVectorManager const& hits) const { // Here fill the plots to compute the different metrics linked to // reco-level, namely fake-rate and merge-rate. In this loop should *not* // restrict only to the selected SimClusters. @@ -2199,7 +2199,7 @@ void HGVHistoProducerAlgo::fill_generic_cluster_histos(const Histograms& histogr std::vector thicknesses, const ticl::RecoToSimCollection& cpsInLayerClusterMap, const ticl::SimToRecoCollection& cPOnLayerMap, - std::vector const& hits) const { + MultiVectorManager const& hits) const { //Each event to be treated as two events: an event in +ve endcap, //plus another event in -ve endcap. In this spirit there will be //a layer variable (layerid) that maps the layers in : @@ -2684,7 +2684,7 @@ void HGVHistoProducerAlgo::fill_trackster_histos( std::vector const& cPSelectedIndices, std::unordered_map const& hitMap, unsigned int layers, - std::vector const& hits, + MultiVectorManager const& hits, bool mapsFound, const edm::Handle& trackstersToSimTrackstersByLCsMapH, const edm::Handle& simTrackstersToTrackstersByLCsMapH, @@ -2932,7 +2932,7 @@ void HGVHistoProducerAlgo::setRecHitTools(std::shared_ptr re DetId HGVHistoProducerAlgo::findmaxhit(const reco::CaloCluster& cluster, std::unordered_map const& hitMap, - std::vector const& hits) const { + MultiVectorManager const& hits) const { const auto& hits_and_fractions = cluster.hitsAndFractions(); DetId themaxid;