diff --git a/DPGAnalysis/HGCalNanoAOD/python/nanoHGCML_cff.py b/DPGAnalysis/HGCalNanoAOD/python/nanoHGCML_cff.py index 1320a72f05d38..bc17d73869416 100644 --- a/DPGAnalysis/HGCalNanoAOD/python/nanoHGCML_cff.py +++ b/DPGAnalysis/HGCalNanoAOD/python/nanoHGCML_cff.py @@ -26,7 +26,7 @@ genParticleTable.variables = cms.PSet(genParticleTable.variables, charge = CandVars.charge) -nanoHGCMLSequence = cms.Sequence(nanoMetadata+genVertexTables+genParticleTable+ +nanoHGCMLSequence = cms.Sequence(nanoMetadata+genVertexTable+genVertexT0Table+genParticleTable+ trackingParticleTable+caloParticleTable+ layerClusterTables+ simTrackTables+hgcSimHitsSequence+trackerSimHitTables+ diff --git a/DPGAnalysis/TrackNanoAOD/plugins/RecoTrackToSimClusterAssociation.cc b/DPGAnalysis/TrackNanoAOD/plugins/RecoTrackToSimClusterAssociation.cc new file mode 100644 index 0000000000000..f09eaa4817c3b --- /dev/null +++ b/DPGAnalysis/TrackNanoAOD/plugins/RecoTrackToSimClusterAssociation.cc @@ -0,0 +1,83 @@ +#include "FWCore/Framework/interface/stream/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 "DataFormats/Math/interface/deltaR.h" +#include "CommonTools/Utils/interface/StringCutObjectSelector.h" +#include "RecoHGCal/GraphReco/interface/HGCalTrackPropagator.h" +#include "FWCore/Framework/interface/MakerMacros.h" +#include "DataFormats/HGCalReco/interface/TICLLayerTile.h" +#include "SimDataFormats/CaloAnalysis/interface/SimCluster.h" +#include "SimDataFormats/CaloAnalysis/interface/SimClusterFwd.h" +#include <vector> + +class RecoTrackToSimClusterAssociation : public edm::stream::EDProducer<> { +public: + RecoTrackToSimClusterAssociation(edm::ParameterSet const& params) + : tracksToken_(consumes<reco::TrackCollection>(params.getParameter<edm::InputTag>("tracks"))), + simClusterToken_(consumes<SimClusterCollection>(params.getParameter<edm::InputTag>("simclusters"))) { + produces<nanoaod::FlatTable>(); + } + + ~RecoTrackToSimClusterAssociation() override {} + + void beginRun(const edm::Run&, const edm::EventSetup& iSetup) override { + trackprop_.getEventSetup(iSetup); + } + + void produce(edm::Event& iEvent, const edm::EventSetup& iSetup) override { + edm::Handle<reco::TrackCollection> tracks; + iEvent.getByToken(tracksToken_, tracks); + + edm::Handle<SimClusterCollection> simClusters; + iEvent.getByToken(simClusterToken_, simClusters); + + ticl::TICLLayerTile scTile; + + std::cout << "Track size is " << tracks->size() << std::endl; + std::cout << "simClusters " << simClusters->size() << std::endl; + + for (size_t i = 0; i < simClusters->size(); i++) { + SimClusterRef sc(simClusters, i); + scTile.fill(sc->impactPoint().eta(), sc->impactPoint().phi(), i); + std::cout << "Filling in bin " << scTile.globalBin(sc->impactPoint().eta(), sc->impactPoint().phi()) << std::endl; + } + + for (const auto& track : *tracks) { + auto prop = trackprop_.propagateObject(track); + float eta = prop.pos.eta(); + float phi = prop.pos.phi(); + if (std::abs(eta) < ticl::TileConstants::minEta || std::abs(eta) > ticl::TileConstants::maxEta) + continue; + + auto searchRange = scTile.searchBoxEtaPhi(eta-0.4, eta+0.4, phi-0.4, phi+0.4); + std::cout << "Looking in search range " << searchRange[0] << " " << searchRange[1] << " " << searchRange[1] << " " << searchRange[3] << std::endl; + for (int etab = searchRange[0]; etab < searchRange[1]; etab++) { + std::cout << "Looking at etab " << etab << std::endl; + for (int phib = searchRange[2]; phib < searchRange[3]; phib++) { + int bin = scTile.globalBin(etab, phib); + const std::vector<unsigned int>& matchIndices = scTile[bin]; + for (auto scidx : matchIndices) { + std::cout << "-->Found an overlap at SCIDX " << scidx << std::endl; + SimClusterRef sc(simClusters, scidx); + if (deltaR(prop.pos, sc->impactPoint()) < 0.4) { + std::cout << "These guys match in DR. The pt ratio is " << sc->impactMomentum().pt()/prop.obj->pt() << std::endl; + } + } + } + } + } + + } + +protected: + const edm::EDGetTokenT<reco::TrackCollection> tracksToken_; + const edm::EDGetTokenT<SimClusterCollection> simClusterToken_; + HGCalTrackPropagator trackprop_; +}; + +DEFINE_FWK_MODULE(RecoTrackToSimClusterAssociation); + diff --git a/DPGAnalysis/TrackNanoAOD/python/tracks_cff.py b/DPGAnalysis/TrackNanoAOD/python/tracks_cff.py index 819f0afab71c0..d5d37970f383a 100644 --- a/DPGAnalysis/TrackNanoAOD/python/tracks_cff.py +++ b/DPGAnalysis/TrackNanoAOD/python/tracks_cff.py @@ -51,4 +51,9 @@ ) ) -trackTables = cms.Sequence(generalTrackTable+generalTrackHGCPositionTable+trackConversionsTable+trackDisplacedTable) +trackSimClusterMatch = cms.EDProducer("RecoTrackToSimClusterAssociation", + tracks = cms.InputTag("generalTracks"), + simclusters = cms.InputTag("hgcSimTruth"), + ) + +trackTables = cms.Sequence(generalTrackTable+generalTrackHGCPositionTable+trackConversionsTable+trackDisplacedTable+trackSimClusterMatch) diff --git a/RecoParticleFlow/PFTruthProducer/plugins/PFTruthParticleProducer.cc b/RecoParticleFlow/PFTruthProducer/plugins/PFTruthParticleProducer.cc index 78ab62c939c5b..51d727df5ecd5 100644 --- a/RecoParticleFlow/PFTruthProducer/plugins/PFTruthParticleProducer.cc +++ b/RecoParticleFlow/PFTruthProducer/plugins/PFTruthParticleProducer.cc @@ -221,7 +221,7 @@ void PFTruthParticleProducer::produce(edm::StreamID, edm::Event &iEvent, const e matchedtptocp.push_back(i); tprefs.push_back(tpref); - std::cout << "added track " << tpref->pdgId() << " "<< tpref->momentum() << ", pt " << tpref->pt()<< std::endl;//DEBUG + //std::cout << "added track " << tpref->pdgId() << " "<< tpref->momentum() << ", pt " << tpref->pt()<< std::endl;//DEBUG } for(const auto g4t: tpref->g4Tracks()){ @@ -233,7 +233,7 @@ void PFTruthParticleProducer::produce(edm::StreamID, edm::Event &iEvent, const e } - std::cout << "PF particle with "<< tprefs.size() << " attached tracking particles " <<std::endl;//DEBUG + //std::cout << "PF particle with "<< tprefs.size() << " attached tracking particles " <<std::endl;//DEBUG PFTruthParticle pftp(tprefs,screfs); pftp.setP4(cp.p4()); pftp.setPdgId(cp.pdgId()); @@ -250,7 +250,7 @@ void PFTruthParticleProducer::produce(edm::StreamID, edm::Event &iEvent, const e idealPFTruth.push_back(pftp); } - std::cout << "ideal PF truth built, splitting" << std::endl;//DEBUG + //std::cout << "ideal PF truth built, splitting" << std::endl;//DEBUG // --- now we have the ideal PF truth, far from realistic. @@ -287,7 +287,7 @@ void PFTruthParticleProducer::produce(edm::StreamID, edm::Event &iEvent, const e scToPFpartIdx[sc.key()] = PFtruth->size(); PFtruth->push_back(npftp); } - std::cout << "Split PF Particle to "<< PFtruth->size()-pre << " neutrals " <<std::endl;//DEBUG + //std::cout << "Split PF Particle to "<< PFtruth->size()-pre << " neutrals " <<std::endl;//DEBUG } else { //remaining: the ones that have tracking particles associated @@ -321,13 +321,9 @@ void PFTruthParticleProducer::produce(edm::StreamID, edm::Event &iEvent, const e int scIdx = scRef.isNonnull() ? scRef.key() : -1; int pfIdx = scIdx >= 0 ? scToPFpartIdx[scIdx] : -1; if (scIdx >= 0) - std::cout << "SC Index is " << scIdx << " pfIdx is " << pfIdx << std::endl; rhToPFpartIdx.at(i) = pfIdx; } - - std::cout << "putting PFTruthParticles" << std::endl;//DEBUG - auto pfTruthHand = iEvent.put(std::move(PFtruth)); // Used to fill the associations in the ntuples