Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2023-03-17 11:16:43

0001 //
0002 
0003 /**
0004   \class    ObjectMultiplicityCounter"
0005   \brief    Matcher of number of reconstructed objects in the event to probe
0006 
0007   \author   Kalanand Mishra
0008 */
0009 
0010 #include "FWCore/Framework/interface/stream/EDProducer.h"
0011 #include "FWCore/Framework/interface/Event.h"
0012 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0013 #include "FWCore/Utilities/interface/InputTag.h"
0014 
0015 #include "DataFormats/Math/interface/deltaR.h"
0016 #include "DataFormats/Common/interface/ValueMap.h"
0017 #include "DataFormats/Common/interface/View.h"
0018 
0019 #include "DataFormats/Candidate/interface/CandidateFwd.h"
0020 #include "DataFormats/Candidate/interface/Candidate.h"
0021 #include "DataFormats/TrackReco/interface/Track.h"
0022 #include "DataFormats/VertexReco/interface/Vertex.h"
0023 #include "CommonTools/Utils/interface/StringCutObjectSelector.h"
0024 
0025 template <typename T>
0026 class ObjectMultiplicityCounter : public edm::stream::EDProducer<> {
0027 public:
0028   explicit ObjectMultiplicityCounter(const edm::ParameterSet& iConfig);
0029   ~ObjectMultiplicityCounter() override;
0030 
0031   void produce(edm::Event& iEvent, const edm::EventSetup& iSetup) override;
0032 
0033 private:
0034   edm::EDGetTokenT<edm::View<reco::Candidate>> probesToken_;
0035   edm::EDGetTokenT<edm::View<T>> objectsToken_;
0036   StringCutObjectSelector<T, true> objCut_;  // lazy parsing, to allow cutting on variables not in reco::Candidate class
0037 };
0038 
0039 template <typename T>
0040 ObjectMultiplicityCounter<T>::ObjectMultiplicityCounter(const edm::ParameterSet& iConfig)
0041     : probesToken_(consumes<edm::View<reco::Candidate>>(iConfig.getParameter<edm::InputTag>("probes"))),
0042       objectsToken_(consumes<edm::View<T>>(iConfig.getParameter<edm::InputTag>("objects"))),
0043       objCut_(
0044           iConfig.existsAs<std::string>("objectSelection") ? iConfig.getParameter<std::string>("objectSelection") : "",
0045           true) {
0046   produces<edm::ValueMap<float>>();
0047 }
0048 
0049 template <typename T>
0050 ObjectMultiplicityCounter<T>::~ObjectMultiplicityCounter() {}
0051 
0052 template <typename T>
0053 void ObjectMultiplicityCounter<T>::produce(edm::Event& iEvent, const edm::EventSetup& iSetup) {
0054   using namespace edm;
0055 
0056   // read input
0057   Handle<View<reco::Candidate>> probes;
0058   Handle<View<T>> objects;
0059   iEvent.getByToken(probesToken_, probes);
0060   iEvent.getByToken(objectsToken_, objects);
0061 
0062   // fill
0063   float count = 0.0;
0064   typename View<T>::const_iterator object, endobjects = objects->end();
0065   for (object = objects->begin(); object != endobjects; ++object) {
0066     if (!(objCut_(*object)))
0067       continue;
0068     count += 1.0;
0069   }
0070 
0071   // prepare vector for output
0072   std::vector<float> values(probes->size(), count);
0073 
0074   // convert into ValueMap and store
0075   auto valMap = std::make_unique<ValueMap<float>>();
0076   ValueMap<float>::Filler filler(*valMap);
0077   filler.insert(probes, values.begin(), values.end());
0078   filler.fill();
0079   iEvent.put(std::move(valMap));
0080 }
0081 
0082 typedef ObjectMultiplicityCounter<reco::Candidate> CandMultiplicityCounter;
0083 typedef ObjectMultiplicityCounter<reco::Track> TrackMultiplicityCounter;
0084 typedef ObjectMultiplicityCounter<reco::Vertex> VertexMultiplicityCounter;
0085 
0086 #include "FWCore/Framework/interface/MakerMacros.h"
0087 DEFINE_FWK_MODULE(CandMultiplicityCounter);
0088 DEFINE_FWK_MODULE(TrackMultiplicityCounter);
0089 DEFINE_FWK_MODULE(VertexMultiplicityCounter);