diff --git a/CommonTools/RecoAlgos/plugins/SimHitRecHitAssocitionProducer.cc b/CommonTools/RecoAlgos/plugins/SimClusterRecHitAssocitionProducer.cc similarity index 69% rename from CommonTools/RecoAlgos/plugins/SimHitRecHitAssocitionProducer.cc rename to CommonTools/RecoAlgos/plugins/SimClusterRecHitAssocitionProducer.cc index 59ef9993233b8..6971ccd324bb5 100644 --- a/CommonTools/RecoAlgos/plugins/SimHitRecHitAssocitionProducer.cc +++ b/CommonTools/RecoAlgos/plugins/SimClusterRecHitAssocitionProducer.cc @@ -32,10 +32,10 @@ // typedef std::pair IdxAndFraction; -class SimHitRecHitAssociationProducer : public edm::stream::EDProducer<> { +class SimClusterRecHitAssociationProducer : public edm::stream::EDProducer<> { public: - explicit SimHitRecHitAssociationProducer(const edm::ParameterSet &); - ~SimHitRecHitAssociationProducer() override; + explicit SimClusterRecHitAssociationProducer(const edm::ParameterSet &); + ~SimClusterRecHitAssociationProducer() override; private: virtual void produce(edm::Event&, const edm::EventSetup&) override; @@ -47,7 +47,7 @@ class SimHitRecHitAssociationProducer : public edm::stream::EDProducer<> { edm::EDGetTokenT scCollectionToken_; }; -SimHitRecHitAssociationProducer::SimHitRecHitAssociationProducer(const edm::ParameterSet &pset) +SimClusterRecHitAssociationProducer::SimClusterRecHitAssociationProducer(const edm::ParameterSet &pset) : //caloSimhitCollectionTokens_(edm::vector_transform(pset.getParameter>("caloSimHits"), // [this](const edm::InputTag& tag) {return mayConsume(tag); })), //trackSimhitCollectionTokens_(edm::vector_transform(pset.getParameter("trackSimHits"), @@ -61,44 +61,68 @@ SimHitRecHitAssociationProducer::SimHitRecHitAssociationProducer(const edm::Para std::string label = tag.instance(); produces>(label+"ToSimClus"); } + produces>(); } -SimHitRecHitAssociationProducer::~SimHitRecHitAssociationProducer() {} +SimClusterRecHitAssociationProducer::~SimClusterRecHitAssociationProducer() {} // // member functions // // ------------ method called to produce the data ------------ -void SimHitRecHitAssociationProducer::produce(edm::Event &iEvent, const edm::EventSetup &iSetup) { +void SimClusterRecHitAssociationProducer::produce(edm::Event &iEvent, const edm::EventSetup &iSetup) { + auto simClusterToRecEnergy = std::make_unique>(); edm::Handle scCollection; iEvent.getByToken(scCollectionToken_, scCollection); std::unordered_map hitDetIdToIndex; + std::unordered_map hitDetIdToTotalSimFrac; for (size_t s = 0; s < scCollection->size(); s++) { const auto& sc = scCollection->at(s); + (*simClusterToRecEnergy)[s] = 0.; for (auto& hf : sc.hits_and_fractions()) { auto entry = hitDetIdToIndex.find(hf.first); // Update SimCluster assigment if detId has been found in no other SCs or if // SC has greater fraction of energy in DetId than the SC already found - if (entry == hitDetIdToIndex.end() || entry->second.second < hf.second) + if (entry == hitDetIdToIndex.end()) { + hitDetIdToTotalSimFrac[hf.first] = hf.second; hitDetIdToIndex[hf.first] = {s, hf.second}; + } + else { + hitDetIdToTotalSimFrac[hf.first] += hf.second; + if (entry->second.second < hf.second) + hitDetIdToIndex[hf.first] = {s, hf.second}; + } } } + std::unordered_map hitDetIdToTotalRecEnergy; + for (size_t i = 0; i < caloRechitCollectionTokens_.size(); i++) { std::string label = caloRechitTags_.at(i).instance(); std::vector rechitIndices; edm::Handle> caloRechitCollection; iEvent.getByToken(caloRechitCollectionTokens_.at(i), caloRechitCollection); - //edm::Handle caloSimhitCollection; - //iEvent.getByToken(caloSimhitCollectionTokens_.at(i), caloSimhitCollection); for (size_t h = 0; h < caloRechitCollection->size(); h++) { const CaloRecHit& caloRh = caloRechitCollection->at(h); size_t id = caloRh.detid().rawId(); + auto entry = hitDetIdToTotalRecEnergy.find(id); + float energy = caloRh.energy(); + if (entry == hitDetIdToTotalRecEnergy.end()) + hitDetIdToTotalRecEnergy[id] = energy; + else + hitDetIdToTotalRecEnergy.at(id) += energy; + int match = hitDetIdToIndex.find(id) == hitDetIdToIndex.end() ? -1 : hitDetIdToIndex.at(id).first; + float fraction = match != -1 ? hitDetIdToTotalSimFrac.at(id) : 1.; + if (simClusterToRecEnergy->find(match) == simClusterToRecEnergy->end()) + (*simClusterToRecEnergy)[match] = energy*fraction; + else + simClusterToRecEnergy->at(match) += energy*fraction; + rechitIndices.push_back(match); } @@ -108,8 +132,9 @@ void SimHitRecHitAssociationProducer::produce(edm::Event &iEvent, const edm::Eve filler.fill(); iEvent.put(std::move(assoc), label+"ToSimClus"); } + iEvent.put(std::move(simClusterToRecEnergy)); } // define this as a plug-in -DEFINE_FWK_MODULE(SimHitRecHitAssociationProducer); +DEFINE_FWK_MODULE(SimClusterRecHitAssociationProducer); diff --git a/IOMC/ParticleGuns/src/SealModule.cc b/IOMC/ParticleGuns/src/SealModule.cc index 020315ce98a68..d21e3d2b1f01f 100644 --- a/IOMC/ParticleGuns/src/SealModule.cc +++ b/IOMC/ParticleGuns/src/SealModule.cc @@ -23,7 +23,6 @@ #include "IOMC/ParticleGuns/interface/RandomtXiGunProducer.h" #include "IOMC/ParticleGuns/interface/RandomMultiParticlePGunProducer.h" #include "IOMC/ParticleGuns/interface/RandomXiThetaGunProducer.h" -#include "IOMC/ParticleGuns/interface/FlatEtaRangeGunProducer.h" // particle gun prototypes // @@ -68,5 +67,3 @@ using edm::RandomMultiParticlePGunProducer; DEFINE_FWK_MODULE(RandomMultiParticlePGunProducer); using edm::RandomXiThetaGunProducer; DEFINE_FWK_MODULE(RandomXiThetaGunProducer); -using edm::FlatEtaRangeGunProducer; -DEFINE_FWK_MODULE(FlatEtaRangeGunProducer); diff --git a/PhysicsTools/NanoAOD/plugins/ObjectIndexFromAssociationProducer.cc b/PhysicsTools/NanoAOD/plugins/ObjectIndexFromAssociationProducer.cc index 0ef18e2df143b..447fc98bedbfb 100644 --- a/PhysicsTools/NanoAOD/plugins/ObjectIndexFromAssociationProducer.cc +++ b/PhysicsTools/NanoAOD/plugins/ObjectIndexFromAssociationProducer.cc @@ -45,7 +45,7 @@ class ObjectIndexFromAssociationTableProducer : public edm::global::EDProducer<> std::vector keys; for (unsigned int i = 0; i < objs->size(); ++i) { - edm::Ref tk(objs, i); + edm::Ref tk(objs, i); if (cut_(*tk)) { edm::Ref match = (*assoc)[tk]; int key = match.isNonnull() ? match.key() : -1; diff --git a/PhysicsTools/NanoAOD/plugins/ObjectPropertyFromIndexMapProducer.cc b/PhysicsTools/NanoAOD/plugins/ObjectPropertyFromIndexMapProducer.cc new file mode 100644 index 0000000000000..2e80efe37a0fa --- /dev/null +++ b/PhysicsTools/NanoAOD/plugins/ObjectPropertyFromIndexMapProducer.cc @@ -0,0 +1,76 @@ +#include "FWCore/Framework/interface/global/EDProducer.h" +#include "FWCore/Framework/interface/Event.h" +#include "FWCore/ParameterSet/interface/ParameterSet.h" +#include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h" +#include "FWCore/ParameterSet/interface/ParameterSetDescription.h" +#include "DataFormats/NanoAOD/interface/FlatTable.h" +#include "DataFormats/Common/interface/View.h" +#include "SimDataFormats/CaloAnalysis/interface/SimCluster.h" +#include "SimDataFormats/CaloAnalysis/interface/SimClusterFwd.h" +#include "DataFormats/Common/interface/Association.h" +#include "CommonTools/Utils/interface/StringCutObjectSelector.h" + +#include +#include +#include + +template +class ObjectPropertyFromIndexMapTableProducer : public edm::global::EDProducer<> { +public: + ObjectPropertyFromIndexMapTableProducer(edm::ParameterSet const& params) + : objName_(params.getParameter("objName")), + branchName_(params.getParameter("branchName")), + doc_(params.getParameter("docString")), + src_(consumes(params.getParameter("src"))), + mapToken_(consumes>(params.getParameter("valueMap"))), + cut_(params.getParameter("cut"), true) { + produces(); + } + + ~ObjectPropertyFromIndexMapTableProducer() override {} + // + // Because I'm not sure if this can be templated, overload instead... + std::unique_ptr fillTable(const std::vector& values, const std::string& objName) const { + auto tab = std::make_unique(values.size(), objName, false, true); + tab->addColumn(branchName_, values, doc_); + return tab; + } + + // Because I'm not sure if this can be templated, overload instead... + std::unique_ptr fillTable(const std::vector& values, const std::string& objName) const { + auto tab = std::make_unique(values.size(), objName, false, true); + tab->addColumn(branchName_, values, doc_); + return tab; + } + + void produce(edm::StreamID id, edm::Event& iEvent, const edm::EventSetup& iSetup) const override { + edm::Handle objs; + iEvent.getByToken(src_, objs); + + edm::Handle> valueMap; + iEvent.getByToken(mapToken_, valueMap); + + std::vector values; + for (unsigned int i = 0; i < objs->size(); ++i) { + edm::Ref obj(objs, i); + if (cut_(*obj)) { + if (valueMap->find(i) == valueMap->end()) + throw cms::Exception("ObjectPropertyFromIndexMapTableProducer") << "No entry in value map for candidate " << i; + values.emplace_back(valueMap->at(i)); + } + } + + auto tab = fillTable(values, objName_); + iEvent.put(std::move(tab)); + } + +protected: + const std::string objName_, branchName_, doc_; + const edm::EDGetTokenT src_; + const edm::EDGetTokenT> mapToken_; + const StringCutObjectSelector cut_; +}; + +#include "FWCore/Framework/interface/MakerMacros.h" +typedef ObjectPropertyFromIndexMapTableProducer SimClusterRecEnergyTableProducer; +DEFINE_FWK_MODULE(SimClusterRecEnergyTableProducer); diff --git a/PhysicsTools/NanoAOD/plugins/SimpleFlatTableProducerPlugins.cc b/PhysicsTools/NanoAOD/plugins/SimpleFlatTableProducerPlugins.cc index d9a92164fc366..ed0107e48b74d 100644 --- a/PhysicsTools/NanoAOD/plugins/SimpleFlatTableProducerPlugins.cc +++ b/PhysicsTools/NanoAOD/plugins/SimpleFlatTableProducerPlugins.cc @@ -15,8 +15,6 @@ typedef SimpleFlatTableProducer SimpleCaloRecHitFlatTableProducer; #include "SimDataFormats/TrackingAnalysis/interface/TrackingParticle.h" typedef SimpleFlatTableProducer SimpleTrackingParticleFlatTableProducer; -#include "SimDataFormats/PFAnalysis/interface/PFParticle.h" -typedef SimpleFlatTableProducer SimplePFParticleFlatTableProducer; #include "SimDataFormats/CaloAnalysis/interface/CaloParticle.h" typedef SimpleFlatTableProducer SimpleCaloParticleFlatTableProducer; #include "DataFormats/ParticleFlowCandidate/interface/PFCandidate.h" @@ -50,7 +48,6 @@ DEFINE_FWK_MODULE(SimpleSimTrackFlatTableProducer); DEFINE_FWK_MODULE(SimpleSimClusterFlatTableProducer); DEFINE_FWK_MODULE(SimpleTrackingParticleFlatTableProducer); DEFINE_FWK_MODULE(SimpleTrackFlatTableProducer); -DEFINE_FWK_MODULE(SimplePFParticleFlatTableProducer); DEFINE_FWK_MODULE(SimplePFCandidateFlatTableProducer); DEFINE_FWK_MODULE(SimpleCaloParticleFlatTableProducer); DEFINE_FWK_MODULE(SimpleCandidateFlatTableProducer); diff --git a/PhysicsTools/NanoAOD/python/hgcRecHits_cff.py b/PhysicsTools/NanoAOD/python/hgcRecHits_cff.py index a010dbba7825b..a6b9f1ee10684 100644 --- a/PhysicsTools/NanoAOD/python/hgcRecHits_cff.py +++ b/PhysicsTools/NanoAOD/python/hgcRecHits_cff.py @@ -15,7 +15,7 @@ ) ) -hgcRecHitsToSimClusters = cms.EDProducer("SimHitRecHitAssociationProducer", +hgcRecHitsToSimClusters = cms.EDProducer("SimClusterRecHitAssociationProducer", caloRecHits = cms.VInputTag("HGCalRecHit:HGCEERecHits", "HGCalRecHit:HGCHEFRecHits", "HGCalRecHit:HGCHEBRecHits", ), @@ -31,6 +31,16 @@ docString = cms.string("SimCluster responsible for most sim energy in RecHit DetId") ) +# TODO: Pair with simclusters +simClusterRecEnergyTable = cms.EDProducer("SimClusterRecEnergyTableProducer", + src = cms.InputTag("mix:MergedCaloTruth"), + cut = cms.string(""), + objName = cms.string("SimCluster"), + branchName = cms.string("recEnergy"), + valueMap = cms.InputTag("hgcRecHitsToSimClusters"), + docString = cms.string("SimCluster deposited reconstructed energy associated to SimCluster") +) + hgcRecHitsToPFCands = cms.EDProducer("RecHitToPFCandAssociationProducer", caloRecHits = cms.VInputTag("HGCalRecHit:HGCEERecHits", "HGCalRecHit:HGCHEFRecHits", "HGCalRecHit:HGCHEBRecHits", @@ -92,6 +102,7 @@ hgcRecHitsSequence = cms.Sequence(hgcEERecHitsTable+hgcHEbackRecHitsTable+hgcHEfrontRecHitsTable +hgcRecHitsToSimClusters + +simClusterRecEnergyTable +hgcRecHitsToPFCands +hgcEERecHitsToPFCandTable+hgcHEfrontRecHitsToPFCandTable+hgcHEbackRecHitsToPFCandTable +hgcEERecHitsPositionTable diff --git a/PhysicsTools/NanoAOD/python/simClusters_cff.py b/PhysicsTools/NanoAOD/python/simClusters_cff.py index c26176afc5ee6..c62b31444ddae 100644 --- a/PhysicsTools/NanoAOD/python/simClusters_cff.py +++ b/PhysicsTools/NanoAOD/python/simClusters_cff.py @@ -6,17 +6,12 @@ cut = cms.string(""), name = cms.string("SimCluster"), doc = cms.string("SimCluster information"), - singleton = cms.bool(False), # the number of entries is variable - extension = cms.bool(False), # this is the main table for the muons + singleton = cms.bool(False), + extension = cms.bool(False), variables = cms.PSet(CandVars, lastPos_x = Var('g4Tracks.at(0).trackerSurfacePosition().x()', 'float', precision=14, doc='track x final position'), lastPos_y = Var('g4Tracks.at(0).trackerSurfacePosition().y()', 'float', precision=14, doc='track y final position'), lastPos_z = Var('g4Tracks.at(0).trackerSurfacePosition().z()', 'float', precision=14, doc='track z final position'), - # NOTE: This is the cms-pepr approach, which is needed for merged simclusters - #impactPoint_x = Var('impactPoint().x()', 'float', precision=14, doc='x position'), - #impactPoint_y = Var('impactPoint().y()', 'float', precision=14, doc='y position'), - #impactPoint_z = Var('impactPoint().z()', 'float', precision=14, doc='z position'), - #impactPoint_t = Var('impactPoint().t()', 'float', precision=14, doc='Impact time'), impactPoint_x = Var('g4Tracks().at(0).getPositionAtBoundary().x()', 'float', precision=14, doc='x position'), impactPoint_y = Var('g4Tracks().at(0).getPositionAtBoundary().y()', 'float', precision=14, doc='y position'), impactPoint_z = Var('g4Tracks().at(0).getPositionAtBoundary().z()', 'float', precision=14, doc='z position'), @@ -25,10 +20,10 @@ # are often referred to as rechits in the SimCluster class nHits = Var('numberOfRecHits', 'int', precision=-1, doc='number of simhits'), sumHitEnergy = Var('energy', 'float', precision=14, doc='total energy of simhits'), - #boundaryPmag = Var('impactMomentum.P()', 'float', precision=14, doc='magnitude of the boundary 3-momentum vector'), - #boundaryP4 = Var('impactMomentum.mag()', 'float', precision=14, doc='magnitude of four vector'), - #boundaryEnergy = Var('impactMomentum.energy()', 'float', precision=14, doc='magnitude of four vector'), - #boundaryPt = Var('impactMomentum.pt()', 'float', precision=14, doc='magnitude of four vector'), + boundaryPmag = Var('impactMomentum.P()', 'float', precision=14, doc='magnitude of the boundary 3-momentum vector'), + boundaryP4 = Var('impactMomentum.mag()', 'float', precision=14, doc='magnitude of four vector'), + boundaryEnergy = Var('impactMomentum.energy()', 'float', precision=14, doc='magnitude of four vector'), + boundaryPt = Var('impactMomentum.pt()', 'float', precision=14, doc='magnitude of four vector'), trackId = Var('g4Tracks().at(0).trackId()', 'int', precision=-1, doc='Geant track id'), trackIdAtBoundary = Var('g4Tracks().at(0).getIDAtBoundary()', 'int', precision=-1, doc='Track ID at boundary'), ) @@ -43,21 +38,4 @@ docString = cms.string("Index of CaloPart containing SimCluster") ) -hgcSimTruth = cms.EDProducer("HGCTruthProducer") - -simClusterToMergedSCTable = cms.EDProducer("SimClusterToSimClusterIndexTableProducer", - cut = simClusterTable.cut, - src = simClusterTable.src, - objName = simClusterTable.name, - branchName = cms.string("MergedSimCluster"), - objMap = cms.InputTag("hgcSimTruth"), - docString = cms.string("Index of Merged SimCluster containing SimCluster") -) - -mergedSimClusterTable = simClusterTable.clone() -mergedSimClusterTable.src = "hgcSimTruth" -mergedSimClusterTable.name = "MergedSimCluster" - simClusterTables = cms.Sequence(simClusterTable+simClusterToCaloPartTable) - -mergedSimClusterTables = cms.Sequence(hgcSimTruth+mergedSimClusterTable+simClusterToMergedSCTable) diff --git a/PhysicsTools/NanoAOD/python/trackSimHits_cff.py b/PhysicsTools/NanoAOD/python/trackSimHits_cff.py index 7efb49957b7b6..c1f55f132b517 100644 --- a/PhysicsTools/NanoAOD/python/trackSimHits_cff.py +++ b/PhysicsTools/NanoAOD/python/trackSimHits_cff.py @@ -1,8 +1,5 @@ import FWCore.ParameterSet.Config as cms from PhysicsTools.NanoAOD.common_cff import CandVars,Var -#TrackerHitsTECLowTof -#TrackerHitsPixelEndcapHighTof -#TrackerHitsPixelEndcapLowTof trackerHitsPixelEndcapLowTofTable = cms.EDProducer("SimplePSimHitFlatTableProducer", src = cms.InputTag("g4SimHits:TrackerHitsPixelEndcapLowTof"), diff --git a/SimDataFormats/Associations/src/classes_def.xml b/SimDataFormats/Associations/src/classes_def.xml index 3c66509430337..80d4a0f7009a6 100644 --- a/SimDataFormats/Associations/src/classes_def.xml +++ b/SimDataFormats/Associations/src/classes_def.xml @@ -38,6 +38,9 @@ + + +