File indexing completed on 2024-04-06 12:01:09
0001 #ifndef CommonTools_RecoAlgos_PrimaryVertexSorter_
0002 #define CommonTools_RecoAlgos_PrimaryVertexSorter_
0003
0004
0005 #include <memory>
0006 #include <string>
0007
0008
0009 #include "FWCore/Framework/interface/Frameworkfwd.h"
0010 #include "FWCore/Framework/interface/stream/EDProducer.h"
0011 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0012
0013 #include "DataFormats/Common/interface/ValueMap.h"
0014 #include "DataFormats/Common/interface/Association.h"
0015
0016 #include "FWCore/Framework/interface/Event.h"
0017 #include "FWCore/Framework/interface/MakerMacros.h"
0018
0019 #include "DataFormats/VertexReco/interface/VertexFwd.h"
0020
0021 #include "CommonTools/RecoAlgos/interface/PrimaryVertexAssignment.h"
0022 #include "CommonTools/RecoAlgos/interface/PrimaryVertexSorting.h"
0023 #include "DataFormats/PatCandidates/interface/PackedCandidate.h"
0024
0025 #include "TrackingTools/TransientTrack/interface/TransientTrackBuilder.h"
0026 #include "TrackingTools/Records/interface/TransientTrackRecord.h"
0027
0028
0029
0030
0031
0032
0033 template <class ParticlesCollection>
0034
0035 class PrimaryVertexSorter : public edm::stream::EDProducer<> {
0036 public:
0037 typedef edm::Association<reco::VertexCollection> CandToVertex;
0038 typedef edm::ValueMap<int> CandToVertexQuality;
0039 typedef edm::ValueMap<float> VertexScore;
0040
0041 typedef ParticlesCollection PFCollection;
0042
0043 explicit PrimaryVertexSorter(const edm::ParameterSet&);
0044
0045 ~PrimaryVertexSorter() override {}
0046
0047 void produce(edm::Event&, const edm::EventSetup&) override;
0048
0049 static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0050
0051 private:
0052 PrimaryVertexAssignment assignmentAlgo_;
0053 PrimaryVertexSorting sortingAlgo_;
0054
0055
0056 edm::EDGetTokenT<PFCollection> tokenCandidates_;
0057
0058
0059 edm::EDGetTokenT<reco::VertexCollection> tokenVertices_;
0060 edm::EDGetTokenT<edm::View<reco::Candidate>> tokenJets_;
0061 edm::EDGetTokenT<edm::ValueMap<float>> tokenTrackTimeTag_;
0062 edm::EDGetTokenT<edm::ValueMap<float>> tokenTrackTimeResoTag_;
0063
0064 edm::ESGetToken<TransientTrackBuilder, TransientTrackRecord> tokenBuilder_;
0065
0066 bool produceOriginalMapping_;
0067 bool produceSortedVertices_;
0068 bool producePFPileUp_;
0069 bool producePFNoPileUp_;
0070 int qualityCut_;
0071 bool useMET_;
0072 bool useTiming_;
0073
0074 static void fillDescriptionsCommon(edm::ParameterSetDescription& descriptions);
0075
0076 void doConsumesForTiming(const edm::ParameterSet& iConfig);
0077 bool needsProductsForTiming();
0078 std::pair<int, PrimaryVertexAssignment::Quality> runAlgo(const reco::VertexCollection& vertices,
0079 const typename ParticlesCollection::value_type& pf,
0080 const edm::ValueMap<float>* trackTimeTag,
0081 const edm::ValueMap<float>* trackTimeResoTag,
0082 const edm::View<reco::Candidate>& jets,
0083 const TransientTrackBuilder& builder);
0084 };
0085
0086 #include "DataFormats/VertexReco/interface/Vertex.h"
0087 #include "FWCore/Framework/interface/ESHandle.h"
0088
0089
0090 #include "FWCore/Utilities/interface/Exception.h"
0091 #include "FWCore/Framework/interface/EventSetup.h"
0092
0093 template <class ParticlesCollection>
0094 PrimaryVertexSorter<ParticlesCollection>::PrimaryVertexSorter(const edm::ParameterSet& iConfig)
0095 : assignmentAlgo_(iConfig.getParameterSet("assignment")),
0096 sortingAlgo_(iConfig.getParameterSet("sorting")),
0097 tokenCandidates_(consumes<ParticlesCollection>(iConfig.getParameter<edm::InputTag>("particles"))),
0098 tokenVertices_(consumes<reco::VertexCollection>(iConfig.getParameter<edm::InputTag>("vertices"))),
0099 tokenJets_(consumes<edm::View<reco::Candidate>>(iConfig.getParameter<edm::InputTag>("jets"))),
0100 tokenBuilder_(esConsumes(edm::ESInputTag("", "TransientTrackBuilder"))),
0101 produceOriginalMapping_(iConfig.getParameter<bool>("produceAssociationToOriginalVertices")),
0102 produceSortedVertices_(iConfig.getParameter<bool>("produceSortedVertices")),
0103 producePFPileUp_(iConfig.getParameter<bool>("producePileUpCollection")),
0104 producePFNoPileUp_(iConfig.getParameter<bool>("produceNoPileUpCollection")),
0105 qualityCut_(iConfig.getParameter<int>("qualityForPrimary")),
0106 useMET_(iConfig.getParameter<bool>("usePVMET")),
0107 useTiming_(iConfig.getParameterSet("assignment").getParameter<bool>("useTiming")) {
0108 using namespace std;
0109 using namespace edm;
0110 using namespace reco;
0111
0112 if (produceOriginalMapping_) {
0113 produces<CandToVertex>("original");
0114 produces<CandToVertexQuality>("original");
0115 produces<VertexScore>("original");
0116 }
0117 if (produceSortedVertices_) {
0118 produces<reco::VertexCollection>();
0119 produces<CandToVertex>();
0120 produces<CandToVertexQuality>();
0121 produces<VertexScore>();
0122 }
0123
0124 if (producePFPileUp_) {
0125 if (produceOriginalMapping_)
0126 produces<PFCollection>("originalPileUp");
0127 if (produceSortedVertices_)
0128 produces<PFCollection>("PileUp");
0129 }
0130
0131 if (producePFNoPileUp_) {
0132 if (produceOriginalMapping_)
0133 produces<PFCollection>("originalNoPileUp");
0134 if (produceSortedVertices_)
0135 produces<PFCollection>("NoPileUp");
0136 }
0137
0138 if (useTiming_)
0139 doConsumesForTiming(iConfig);
0140 }
0141
0142 template <class ParticlesCollection>
0143 void PrimaryVertexSorter<ParticlesCollection>::produce(edm::Event& iEvent, const edm::EventSetup& iSetup) {
0144 using namespace std;
0145 using namespace edm;
0146 using namespace reco;
0147
0148 Handle<edm::View<reco::Candidate>> jets;
0149 iEvent.getByToken(tokenJets_, jets);
0150
0151 TransientTrackBuilder const& builder = iSetup.getData(tokenBuilder_);
0152
0153 Handle<VertexCollection> vertices;
0154 iEvent.getByToken(tokenVertices_, vertices);
0155
0156 Handle<ParticlesCollection> particlesHandle;
0157 iEvent.getByToken(tokenCandidates_, particlesHandle);
0158
0159 Handle<edm::ValueMap<float>> trackTimeTagHandle;
0160 Handle<edm::ValueMap<float>> trackTimeResoTagHandle;
0161
0162 const edm::ValueMap<float>* trackTimeTag = nullptr;
0163 const edm::ValueMap<float>* trackTimeResoTag = nullptr;
0164 if (useTiming_ && needsProductsForTiming()) {
0165 iEvent.getByToken(tokenTrackTimeTag_, trackTimeTagHandle);
0166 iEvent.getByToken(tokenTrackTimeResoTag_, trackTimeResoTagHandle);
0167
0168 trackTimeTag = trackTimeTagHandle.product();
0169 trackTimeResoTag = trackTimeResoTagHandle.product();
0170 }
0171
0172 ParticlesCollection particles = *particlesHandle.product();
0173 std::vector<int> pfToPVVector;
0174 std::vector<PrimaryVertexAssignment::Quality> pfToPVQualityVector;
0175
0176 std::vector<std::vector<int>> pvToPFVector(vertices->size());
0177 std::vector<std::vector<const reco::Candidate*>> pvToCandVector(vertices->size());
0178 std::vector<std::vector<PrimaryVertexAssignment::Quality>> pvToPFQualityVector(vertices->size());
0179 std::vector<float> vertexScoreOriginal(vertices->size());
0180 std::vector<float> vertexScore(vertices->size());
0181
0182 for (auto const& pf : particles) {
0183 std::pair<int, PrimaryVertexAssignment::Quality> vtxWithQuality =
0184 runAlgo(*vertices, pf, trackTimeTag, trackTimeResoTag, *jets, builder);
0185 pfToPVVector.push_back(vtxWithQuality.first);
0186 pfToPVQualityVector.push_back(vtxWithQuality.second);
0187 }
0188
0189
0190 for (size_t i = 0; i < pfToPVVector.size(); i++) {
0191 auto pv = pfToPVVector[i];
0192 auto qual = pfToPVQualityVector[i];
0193 if (pv >= 0 and qual >= qualityCut_) {
0194 pvToPFVector[pv].push_back(i);
0195
0196
0197
0198 pvToCandVector[pv].push_back(&particles[i]);
0199 pvToPFQualityVector[pv].push_back(qual);
0200 }
0201 }
0202
0203
0204 std::multimap<float, int> scores;
0205 for (unsigned int i = 0; i < vertices->size(); i++) {
0206 float s = sortingAlgo_.score((*vertices)[i], pvToCandVector[i], useMET_);
0207 vertexScoreOriginal[i] = s;
0208 scores.insert(std::pair<float, int>(-s, i));
0209 }
0210
0211
0212 std::vector<int> oldToNew(vertices->size()), newToOld(vertices->size());
0213 size_t newIdx = 0;
0214 for (auto const& idx : scores) {
0215
0216 vertexScore[newIdx] = -idx.first;
0217 oldToNew[idx.second] = newIdx;
0218 newToOld[newIdx] = idx.second;
0219 newIdx++;
0220 }
0221
0222 if (produceOriginalMapping_) {
0223 unique_ptr<CandToVertex> pfCandToOriginalVertexOutput(new CandToVertex(vertices));
0224 unique_ptr<CandToVertexQuality> pfCandToOriginalVertexQualityOutput(new CandToVertexQuality());
0225 CandToVertex::Filler cand2VertexFiller(*pfCandToOriginalVertexOutput);
0226 CandToVertexQuality::Filler cand2VertexQualityFiller(*pfCandToOriginalVertexQualityOutput);
0227
0228 cand2VertexFiller.insert(particlesHandle, pfToPVVector.begin(), pfToPVVector.end());
0229 cand2VertexQualityFiller.insert(particlesHandle, pfToPVQualityVector.begin(), pfToPVQualityVector.end());
0230
0231 cand2VertexFiller.fill();
0232 cand2VertexQualityFiller.fill();
0233 iEvent.put(std::move(pfCandToOriginalVertexOutput), "original");
0234 iEvent.put(std::move(pfCandToOriginalVertexQualityOutput), "original");
0235
0236 unique_ptr<VertexScore> vertexScoreOriginalOutput(new VertexScore);
0237 VertexScore::Filler vertexScoreOriginalFiller(*vertexScoreOriginalOutput);
0238 vertexScoreOriginalFiller.insert(vertices, vertexScoreOriginal.begin(), vertexScoreOriginal.end());
0239 vertexScoreOriginalFiller.fill();
0240 iEvent.put(std::move(vertexScoreOriginalOutput), "original");
0241 }
0242
0243 if (produceSortedVertices_) {
0244 std::vector<int> pfToSortedPVVector;
0245
0246 for (size_t i = 0; i < pfToPVVector.size(); i++) {
0247 pfToSortedPVVector.push_back(oldToNew[pfToPVVector[i]]);
0248
0249 }
0250
0251 unique_ptr<reco::VertexCollection> sortedVerticesOutput(new reco::VertexCollection);
0252 for (size_t i = 0; i < vertices->size(); i++) {
0253 sortedVerticesOutput->push_back((*vertices)[newToOld[i]]);
0254 }
0255 edm::OrphanHandle<reco::VertexCollection> oh = iEvent.put(std::move(sortedVerticesOutput));
0256 unique_ptr<CandToVertex> pfCandToVertexOutput(new CandToVertex(oh));
0257 unique_ptr<CandToVertexQuality> pfCandToVertexQualityOutput(new CandToVertexQuality());
0258 CandToVertex::Filler cand2VertexFiller(*pfCandToVertexOutput);
0259 CandToVertexQuality::Filler cand2VertexQualityFiller(*pfCandToVertexQualityOutput);
0260
0261 cand2VertexFiller.insert(particlesHandle, pfToSortedPVVector.begin(), pfToSortedPVVector.end());
0262 cand2VertexQualityFiller.insert(particlesHandle, pfToPVQualityVector.begin(), pfToPVQualityVector.end());
0263
0264 cand2VertexFiller.fill();
0265 cand2VertexQualityFiller.fill();
0266 iEvent.put(std::move(pfCandToVertexOutput));
0267 iEvent.put(std::move(pfCandToVertexQualityOutput));
0268
0269 unique_ptr<VertexScore> vertexScoreOutput(new VertexScore);
0270 VertexScore::Filler vertexScoreFiller(*vertexScoreOutput);
0271 vertexScoreFiller.insert(oh, vertexScore.begin(), vertexScore.end());
0272 vertexScoreFiller.fill();
0273 iEvent.put(std::move(vertexScoreOutput));
0274 }
0275
0276 unique_ptr<PFCollection> pfCollectionNOPUOriginalOutput(new PFCollection);
0277 unique_ptr<PFCollection> pfCollectionNOPUOutput(new PFCollection);
0278 unique_ptr<PFCollection> pfCollectionPUOriginalOutput(new PFCollection);
0279 unique_ptr<PFCollection> pfCollectionPUOutput(new PFCollection);
0280
0281 for (size_t i = 0; i < particles.size(); i++) {
0282 auto pv = pfToPVVector[i];
0283 auto qual = pfToPVQualityVector[i];
0284
0285 if (producePFNoPileUp_ && produceSortedVertices_)
0286 if (pv == newToOld[0] and qual >= qualityCut_)
0287 pfCollectionNOPUOutput->push_back(particles[i]);
0288
0289 if (producePFPileUp_ && produceSortedVertices_)
0290 if (pv != newToOld[0] and qual >= qualityCut_)
0291 pfCollectionPUOutput->push_back(particles[i]);
0292
0293 if (producePFNoPileUp_ && produceOriginalMapping_)
0294 if (pv == 0 and qual >= qualityCut_)
0295 pfCollectionNOPUOriginalOutput->push_back(particles[i]);
0296
0297 if (producePFPileUp_ && produceOriginalMapping_)
0298 if (pv != 0 and qual >= qualityCut_)
0299 pfCollectionPUOriginalOutput->push_back(particles[i]);
0300 }
0301 if (producePFNoPileUp_ && produceSortedVertices_)
0302 iEvent.put(std::move(pfCollectionNOPUOutput), "NoPileUp");
0303 if (producePFPileUp_ && produceSortedVertices_)
0304 iEvent.put(std::move(pfCollectionPUOutput), "PileUp");
0305 if (producePFNoPileUp_ && produceOriginalMapping_)
0306 iEvent.put(std::move(pfCollectionNOPUOriginalOutput), "originalNoPileUp");
0307 if (producePFPileUp_ && produceOriginalMapping_)
0308 iEvent.put(std::move(pfCollectionPUOriginalOutput), "originalPileUp");
0309 }
0310
0311 template <>
0312 inline void PrimaryVertexSorter<std::vector<reco::RecoChargedRefCandidate>>::fillDescriptions(
0313 edm::ConfigurationDescriptions& descriptions) {
0314 edm::ParameterSetDescription desc;
0315 fillDescriptionsCommon(desc);
0316 desc.add<edm::InputTag>("trackTimeTag", edm::InputTag(""));
0317 desc.add<edm::InputTag>("trackTimeResoTag", edm::InputTag(""));
0318 desc.add<edm::InputTag>("particles", edm::InputTag("trackRefsForJets"));
0319 desc.add<edm::InputTag>("vertices", edm::InputTag("offlinePrimaryVertices"));
0320 desc.add<edm::InputTag>("jets", edm::InputTag("ak4CaloJetsForTrk"));
0321 desc.add<bool>("produceAssociationToOriginalVertices", false);
0322 desc.add<bool>("produceSortedVertices", true);
0323 desc.add<bool>("producePileUpCollection", false);
0324 desc.add<bool>("produceNoPileUpCollection", false);
0325 descriptions.add("sortedPrimaryVertices", desc);
0326 }
0327
0328 template <>
0329 inline void PrimaryVertexSorter<std::vector<reco::PFCandidate>>::fillDescriptions(
0330 edm::ConfigurationDescriptions& descriptions) {
0331 edm::ParameterSetDescription desc;
0332 fillDescriptionsCommon(desc);
0333 desc.add<edm::InputTag>("particles", edm::InputTag("particleFlow"));
0334 desc.add<edm::InputTag>("vertices", edm::InputTag("offlinePrimaryVertices"));
0335 desc.add<edm::InputTag>("jets", edm::InputTag("ak4PFJets"));
0336 desc.add<bool>("produceAssociationToOriginalVertices", true);
0337 desc.add<bool>("produceSortedVertices", true);
0338 desc.add<bool>("producePileUpCollection", true);
0339 desc.add<bool>("produceNoPileUpCollection", true);
0340 descriptions.add("sortedPFPrimaryVertices", desc);
0341 }
0342
0343 template <>
0344 inline void PrimaryVertexSorter<std::vector<pat::PackedCandidate>>::fillDescriptions(
0345 edm::ConfigurationDescriptions& descriptions) {
0346 edm::ParameterSetDescription desc;
0347 fillDescriptionsCommon(desc);
0348 desc.add<edm::InputTag>("particles", edm::InputTag("packedPFCandidates"));
0349 desc.add<edm::InputTag>("vertices", edm::InputTag("offlineSlimmedPrimaryVertices"));
0350 desc.add<edm::InputTag>("jets", edm::InputTag("slimmedJets"));
0351 desc.add<bool>("produceAssociationToOriginalVertices", true);
0352 desc.add<bool>("produceSortedVertices", true);
0353 desc.add<bool>("producePileUpCollection", true);
0354 desc.add<bool>("produceNoPileUpCollection", true);
0355 descriptions.add("sortedPackedPrimaryVertices", desc);
0356 }
0357
0358 template <class ParticlesCollection>
0359 inline void PrimaryVertexSorter<ParticlesCollection>::fillDescriptionsCommon(edm::ParameterSetDescription& desc) {
0360 {
0361 edm::ParameterSetDescription psd0;
0362 desc.add<edm::ParameterSetDescription>("sorting", psd0);
0363 }
0364 {
0365 edm::ParameterSetDescription psd0;
0366 psd0.add<double>("maxDzSigForPrimaryAssignment", 5.0);
0367 psd0.add<double>("maxDzForPrimaryAssignment", 0.1);
0368 psd0.add<double>("maxDzErrorForPrimaryAssignment", 0.05);
0369 psd0.add<double>("maxDtSigForPrimaryAssignment", 3.0);
0370 psd0.add<double>("maxJetDeltaR", 0.5);
0371 psd0.add<double>("minJetPt", 25);
0372 psd0.add<double>("maxDistanceToJetAxis", 0.07);
0373 psd0.add<double>("maxDzForJetAxisAssigment", 0.1);
0374 psd0.add<double>("maxDxyForJetAxisAssigment", 0.1);
0375 psd0.add<double>("maxDxySigForNotReconstructedPrimary", 2);
0376 psd0.add<double>("maxDxyForNotReconstructedPrimary", 0.01);
0377 psd0.add<bool>("useTiming", false);
0378 psd0.add<bool>("useVertexFit", true);
0379 psd0.add<bool>("preferHighRanked", false);
0380 psd0.add<unsigned int>("NumOfPUVtxsForCharged", 0);
0381 psd0.add<double>("DzCutForChargedFromPUVtxs", 0.2);
0382 psd0.add<double>("PtMaxCharged", -1);
0383 psd0.add<double>("EtaMinUseDz", -1);
0384 psd0.add<bool>("OnlyUseFirstDz", false);
0385 desc.add<edm::ParameterSetDescription>("assignment", psd0);
0386 }
0387 desc.add<int>("qualityForPrimary", 3);
0388 desc.add<bool>("usePVMET", true);
0389 }
0390
0391 template <>
0392 inline void PrimaryVertexSorter<std::vector<reco::RecoChargedRefCandidate>>::doConsumesForTiming(
0393 const edm::ParameterSet& iConfig) {
0394 tokenTrackTimeTag_ = consumes<edm::ValueMap<float>>(iConfig.getParameter<edm::InputTag>("trackTimeTag"));
0395 tokenTrackTimeResoTag_ = consumes<edm::ValueMap<float>>(iConfig.getParameter<edm::InputTag>("trackTimeResoTag"));
0396 }
0397
0398 template <>
0399 inline void PrimaryVertexSorter<std::vector<reco::PFCandidate>>::doConsumesForTiming(const edm::ParameterSet& iConfig) {
0400 }
0401
0402 template <>
0403 inline void PrimaryVertexSorter<std::vector<pat::PackedCandidate>>::doConsumesForTiming(
0404 const edm::ParameterSet& iConfig) {}
0405
0406 template <>
0407 inline bool PrimaryVertexSorter<std::vector<reco::RecoChargedRefCandidate>>::needsProductsForTiming() {
0408 return true;
0409 }
0410
0411 template <>
0412 inline bool PrimaryVertexSorter<std::vector<reco::PFCandidate>>::needsProductsForTiming() {
0413 return false;
0414 }
0415
0416 template <>
0417 inline bool PrimaryVertexSorter<std::vector<pat::PackedCandidate>>::needsProductsForTiming() {
0418 return false;
0419 }
0420
0421 template <>
0422 inline std::pair<int, PrimaryVertexAssignment::Quality>
0423 PrimaryVertexSorter<std::vector<reco::RecoChargedRefCandidate>>::runAlgo(const reco::VertexCollection& vertices,
0424 const reco::RecoChargedRefCandidate& pf,
0425 const edm::ValueMap<float>* trackTimeTag,
0426 const edm::ValueMap<float>* trackTimeResoTag,
0427 const edm::View<reco::Candidate>& jets,
0428 const TransientTrackBuilder& builder) {
0429 return assignmentAlgo_.chargedHadronVertex(vertices, pf, trackTimeTag, trackTimeResoTag, jets, builder);
0430 }
0431
0432 template <>
0433 inline std::pair<int, PrimaryVertexAssignment::Quality> PrimaryVertexSorter<std::vector<reco::PFCandidate>>::runAlgo(
0434 const reco::VertexCollection& vertices,
0435 const reco::PFCandidate& pf,
0436 const edm::ValueMap<float>* trackTimeTag,
0437 const edm::ValueMap<float>* trackTimeResoTag,
0438 const edm::View<reco::Candidate>& jets,
0439 const TransientTrackBuilder& builder) {
0440 return assignmentAlgo_.chargedHadronVertex(vertices, pf, jets, builder);
0441 }
0442
0443 template <>
0444 inline std::pair<int, PrimaryVertexAssignment::Quality> PrimaryVertexSorter<std::vector<pat::PackedCandidate>>::runAlgo(
0445 const reco::VertexCollection& vertices,
0446 const pat::PackedCandidate& pf,
0447 const edm::ValueMap<float>* trackTimeTag,
0448 const edm::ValueMap<float>* trackTimeResoTag,
0449 const edm::View<reco::Candidate>& jets,
0450 const TransientTrackBuilder& builder) {
0451 return assignmentAlgo_.chargedHadronVertex(vertices, pf, jets, builder);
0452 }
0453
0454 #endif