Back to home page

Project CMSSW displayed by LXR
 
 

    

File indexing completed on 2023-10-25 10:01:22 

0001 #ifndef MuonIsolationProducers_MuIsolatorResultProducer_H
0002 #define MuonIsolationProducers_MuIsolatorResultProducer_H
0003 
0004 #include "FWCore/Framework/interface/stream/EDProducer.h"
0005 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0006 
0007 #include "RecoMuon/MuonIsolation/interface/MuIsoBaseIsolator.h"
0008 
0009 #include "FWCore/Framework/interface/Event.h"
0010 
0011 #include "DataFormats/Common/interface/RefToBase.h"
0012 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0013 
0014 #include "DataFormats/BeamSpot/interface/BeamSpot.h"
0015 #include "DataFormats/Common/interface/ValueMap.h"
0016 
0017 #include <string>
0018 
0019 namespace edm {
0020   class EventSetup;
0021 }
0022 
0023 struct muisorhelper {
0024   typedef muonisolation::MuIsoBaseIsolator Isolator;
0025   typedef Isolator::Result Result;
0026   typedef Isolator::ResultType ResultType;
0027   typedef std::vector<Result> Results;
0028   typedef Isolator::DepositContainer DepositContainer;
0029 
0030   template <typename BT>
0031   class CandMap {
0032   public:
0033     typedef typename edm::RefToBase<BT> key_type;
0034     typedef typename edm::Handle<edm::View<BT>> handle_type;
0035     typedef DepositContainer value_type;
0036     typedef std::pair<key_type, value_type> pair_type;
0037     typedef typename std::vector<pair_type> map_type;
0038     typedef typename map_type::iterator iterator;
0039 
0040     map_type& get() { return cMap_; }
0041     const map_type& get() const { return cMap_; }
0042     const handle_type& handle() const { return handle_; }
0043     void setHandle(const handle_type& rhs) { handle_ = rhs; }
0044 
0045   private:
0046     map_type cMap_;
0047     handle_type handle_;
0048   };
0049 };
0050 
0051 //! BT == base type
0052 template <typename BT = reco::Candidate>
0053 class MuIsolatorResultProducer : public edm::stream::EDProducer<> {
0054 public:
0055   MuIsolatorResultProducer(const edm::ParameterSet&);
0056 
0057   ~MuIsolatorResultProducer() override;
0058 
0059   void produce(edm::Event&, const edm::EventSetup&) override;
0060 
0061 private:
0062   typedef muisorhelper::Isolator Isolator;
0063   typedef muisorhelper::Result Result;
0064   typedef muisorhelper::ResultType ResultType;
0065   typedef muisorhelper::Results Results;
0066   typedef muisorhelper::DepositContainer DepositContainer;
0067 
0068   typedef muisorhelper::CandMap<BT> CandMap;
0069 
0070   struct DepositConf {
0071     edm::InputTag tag;
0072     double weight;
0073     double threshold;
0074   };
0075 
0076   struct VetoCuts {
0077     bool selectAll;
0078     double muAbsEtaMax;
0079     double muPtMin;
0080     double muAbsZMax;
0081     double muD0Max;
0082   };
0083 
0084   void callWhatProduces();
0085 
0086   unsigned int initAssociation(edm::Event& event, CandMap& candMapT) const;
0087 
0088   void initVetos(reco::TrackBase::Point const& theBeam,
0089                  std::vector<reco::IsoDeposit::Vetos*>& vetos,
0090                  CandMap& candMap) const;
0091 
0092   template <typename RT>
0093   void writeOutImpl(edm::Event& event, const CandMap& candMapT, const Results& results) const;
0094 
0095   void writeOut(edm::Event& event, const CandMap& candMap, const Results& results) const;
0096 
0097   edm::ParameterSet theConfig;
0098   std::vector<DepositConf> theDepositConfs;
0099 
0100   //!choose which muon vetos should be removed from all deposits
0101   bool theRemoveOtherVetos;
0102   VetoCuts theVetoCuts;
0103 
0104   //!the isolator
0105   Isolator* theIsolator;
0106   ResultType theResultType;
0107 
0108   //! beam spot
0109   std::string theBeamlineOption;
0110   edm::InputTag theBeamSpotLabel;
0111 };
0112 
0113 //! actually do the writing here
0114 template <typename BT>
0115 template <typename RT>
0116 inline void MuIsolatorResultProducer<BT>::writeOutImpl(edm::Event& event,
0117                                                        const CandMap& candMapT,
0118                                                        const Results& results) const {
0119   //! make an output vec of what's to be written with a concrete type
0120   std::vector<RT> resV(results.size());
0121   for (unsigned int i = 0; i < resV.size(); ++i)
0122     resV[i] = results[i].val<RT>();
0123   auto outMap = std::make_unique<edm::ValueMap<RT>>();
0124   typename edm::ValueMap<RT>::Filler filler(*outMap);
0125 
0126   //! fill/insert of non-empty values only
0127   if (!candMapT.get().empty()) {
0128     filler.insert(candMapT.handle(), resV.begin(), resV.end());
0129     filler.fill();
0130   }
0131 
0132   event.put(std::move(outMap));
0133 }
0134 
0135 //! choose which result type to write here
0136 template <typename BT>
0137 inline void MuIsolatorResultProducer<BT>::writeOut(edm::Event& event,
0138                                                    const CandMap& candMapT,
0139                                                    const Results& results) const {
0140   std::string metname = "RecoMuon|MuonIsolationProducers";
0141   LogDebug(metname) << "Before calling writeOutMap  with result type " << theIsolator->resultType();
0142 
0143   if (theResultType == Isolator::ISOL_INT_TYPE)
0144     writeOutImpl<int>(event, candMapT, results);
0145   if (theResultType == Isolator::ISOL_FLOAT_TYPE)
0146     writeOutImpl<float>(event, candMapT, results);
0147   if (theResultType == Isolator::ISOL_BOOL_TYPE)
0148     writeOutImpl<bool>(event, candMapT, results);
0149 }
0150 
0151 //! declare what's going to be produced
0152 template <typename BT>
0153 inline void MuIsolatorResultProducer<BT>::callWhatProduces() {
0154   if (theResultType == Isolator::ISOL_FLOAT_TYPE)
0155     produces<edm::ValueMap<float>>();
0156   if (theResultType == Isolator::ISOL_INT_TYPE)
0157     produces<edm::ValueMap<int>>();
0158   if (theResultType == Isolator::ISOL_BOOL_TYPE)
0159     produces<edm::ValueMap<bool>>();
0160 }
0161 
0162 // Framework
0163 #include "FWCore/Framework/interface/EventSetup.h"
0164 #include "DataFormats/Common/interface/Handle.h"
0165 
0166 #include "FWCore/Framework/interface/ESHandle.h"
0167 
0168 #include "DataFormats/MuonReco/interface/Muon.h"
0169 #include "DataFormats/RecoCandidate/interface/IsoDeposit.h"
0170 #include "DataFormats/RecoCandidate/interface/IsoDepositFwd.h"
0171 #include "DataFormats/Common/interface/RefToBase.h"
0172 #include "DataFormats/TrackReco/interface/Track.h"
0173 #include "DataFormats/TrackReco/interface/TrackFwd.h"
0174 
0175 #include "RecoMuon/MuonIsolation/interface/Range.h"
0176 #include "DataFormats/RecoCandidate/interface/IsoDepositDirection.h"
0177 
0178 #include "RecoMuon/MuonIsolation/interface/IsolatorByDeposit.h"
0179 #include "RecoMuon/MuonIsolation/interface/IsolatorByDepositCount.h"
0180 #include "RecoMuon/MuonIsolation/interface/IsolatorByNominalEfficiency.h"
0181 
0182 #include <string>
0183 
0184 //! constructor with config
0185 template <typename BT>
0186 MuIsolatorResultProducer<BT>::MuIsolatorResultProducer(const edm::ParameterSet& par)
0187     : theConfig(par), theRemoveOtherVetos(par.getParameter<bool>("RemoveOtherVetos")), theIsolator(nullptr) {
0188   LogDebug("RecoMuon|MuonIsolation") << " MuIsolatorResultProducer CTOR";
0189 
0190   //! read input config for deposit types and weights and thresholds to apply to them
0191   std::vector<edm::ParameterSet> depositInputs = par.getParameter<std::vector<edm::ParameterSet>>("InputMuIsoDeposits");
0192 
0193   std::vector<double> dWeights(depositInputs.size());
0194   std::vector<double> dThresholds(depositInputs.size());
0195 
0196   for (unsigned int iDep = 0; iDep < depositInputs.size(); ++iDep) {
0197     DepositConf dConf;
0198     dConf.tag = depositInputs[iDep].getParameter<edm::InputTag>("DepositTag");
0199     dConf.weight = depositInputs[iDep].getParameter<double>("DepositWeight");
0200     dConf.threshold = depositInputs[iDep].getParameter<double>("DepositThreshold");
0201 
0202     dWeights[iDep] = dConf.weight;
0203     dThresholds[iDep] = dConf.threshold;
0204 
0205     theDepositConfs.push_back(dConf);
0206   }
0207 
0208   edm::ParameterSet isoPset = par.getParameter<edm::ParameterSet>("IsolatorPSet");
0209   //! will switch to a factory at some point
0210   std::string isolatorType = isoPset.getParameter<std::string>("ComponentName");
0211   if (isolatorType == "IsolatorByDeposit") {
0212     std::string coneSizeType = isoPset.getParameter<std::string>("ConeSizeType");
0213     if (coneSizeType == "FixedConeSize") {
0214       float coneSize(isoPset.getParameter<double>("coneSize"));
0215 
0216       theIsolator = new muonisolation::IsolatorByDeposit(coneSize, dWeights, dThresholds);
0217 
0218       //      theIsolator = new IsolatorByDeposit(isoPset);
0219     } else if (coneSizeType == "CutsConeSize") {
0220       //! FIXME
0221       //       Cuts cuts(isoPset.getParameter<edm::ParameterSet>("CutsPSet"));
0222 
0223       //       theIsolator = new IsolatorByDeposit(coneSize, dWeights, dThresholds);
0224     }
0225   } else if (isolatorType == "IsolatorByNominalEfficiency") {
0226     //! FIXME: need to get the file name here
0227     theIsolator =
0228         new muonisolation::IsolatorByNominalEfficiency("noname", std::vector<std::string>(1, "8:0.97"), dWeights);
0229   } else if (isolatorType == "IsolatorByDepositCount") {
0230     std::string coneSizeType = isoPset.getParameter<std::string>("ConeSizeType");
0231     if (coneSizeType == "FixedConeSize") {
0232       float coneSize(isoPset.getParameter<double>("coneSize"));
0233 
0234       theIsolator = new muonisolation::IsolatorByDepositCount(coneSize, dThresholds);
0235 
0236       //      theIsolator = new IsolatorByDeposit(isoPset);
0237     } else if (coneSizeType == "CutsConeSize") {
0238       //       Cuts cuts(isoPset.getParameter<edm::ParameterSet>("CutsPSet"));
0239 
0240       //       theIsolator = new IsolatorByDeposit(coneSize, dWeights, dThresholds);
0241     }
0242   }
0243 
0244   if (theIsolator == nullptr) {
0245     edm::LogError("MuonIsolationProducers") << "Failed to initialize an isolator";
0246   }
0247   theResultType = theIsolator->resultType();
0248 
0249   callWhatProduces();
0250 
0251   if (theRemoveOtherVetos) {
0252     edm::ParameterSet vetoPSet = par.getParameter<edm::ParameterSet>("VetoPSet");
0253     theVetoCuts.selectAll = vetoPSet.getParameter<bool>("SelectAll");
0254 
0255     //! "other vetoes" is limited to the same collection now
0256     //! for non-trivial choice an external map with pre-made selection flags
0257     //! can be a better choice
0258     if (!theVetoCuts.selectAll) {
0259       theVetoCuts.muAbsEtaMax = vetoPSet.getParameter<double>("MuAbsEtaMax");
0260       theVetoCuts.muPtMin = vetoPSet.getParameter<double>("MuPtMin");
0261       theVetoCuts.muAbsZMax = vetoPSet.getParameter<double>("MuAbsZMax");
0262       theVetoCuts.muD0Max = vetoPSet.getParameter<double>("MuD0Max");
0263       theBeamlineOption = par.getParameter<std::string>("BeamlineOption");
0264       theBeamSpotLabel = par.getParameter<edm::InputTag>("BeamSpotLabel");
0265     }
0266   }
0267 }
0268 
0269 //! destructor
0270 template <typename BT>
0271 MuIsolatorResultProducer<BT>::~MuIsolatorResultProducer() {
0272   if (theIsolator)
0273     delete theIsolator;
0274   LogDebug("RecoMuon|MuIsolatorResultProducer") << " MuIsolatorResultProducer DTOR";
0275 }
0276 
0277 template <typename BT>
0278 void MuIsolatorResultProducer<BT>::produce(edm::Event& event, const edm::EventSetup& eventSetup) {
0279   std::string metname = "RecoMuon|MuonIsolationProducers";
0280   LogDebug(metname) << " Muon Deposit producing..."
0281                     << " BEGINING OF EVENT "
0282                     << "================================";
0283 
0284   reco::TrackBase::Point theBeam = reco::TrackBase::Point(0, 0, 0);
0285 
0286   //! do it only if needed
0287   if (theRemoveOtherVetos && !theVetoCuts.selectAll) {
0288     if (theBeamlineOption == "BeamSpotFromEvent") {
0289       //pick beamSpot
0290       reco::BeamSpot beamSpot;
0291       edm::Handle<reco::BeamSpot> beamSpotH;
0292 
0293       event.getByLabel(theBeamSpotLabel, beamSpotH);
0294 
0295       if (beamSpotH.isValid()) {
0296         theBeam = beamSpotH->position();
0297         LogTrace(metname) << "Extracted beam point at " << theBeam << std::endl;
0298       }
0299     }
0300   }
0301 
0302   //! "smart" container
0303   //! used to repackage deposits_type_candIndex into deposits_candIndex_type
0304   //! have to have it for veto removal (could do away without it otherwise)
0305   //! IMPORTANT: ALL THE REFERENCING BUSINESS IS DONE THROUGH POINTERS
0306   //! Access to the mapped values as reference type HAS TO BE AVAILABLE
0307   CandMap candMapT;
0308 
0309   unsigned int colSize = initAssociation(event, candMapT);
0310 
0311   //! isolator results will be here
0312   Results results(colSize);
0313 
0314   //! extra vetos will be filled here
0315   std::vector<reco::IsoDeposit::Vetos*> vetoDeps(theDepositConfs.size(), nullptr);
0316 
0317   if (colSize != 0) {
0318     if (theRemoveOtherVetos) {
0319       initVetos(theBeam, vetoDeps, candMapT);
0320     }
0321 
0322     //! call the isolator result, passing {[deposit,vetos]_type} set and the candidate
0323     for (unsigned int muI = 0; muI < colSize; ++muI) {
0324       results[muI] = theIsolator->result(candMapT.get()[muI].second, *(candMapT.get()[muI].first));
0325 
0326       if (results[muI].typeF() != theIsolator->resultType()) {
0327         edm::LogError("MuonIsolationProducers") << "Failed to get result from the isolator";
0328       }
0329     }
0330   }
0331 
0332   LogDebug(metname) << "Ready to write out results of size " << results.size();
0333   writeOut(event, candMapT, results);
0334 
0335   for (unsigned int iDep = 0; iDep < vetoDeps.size(); ++iDep) {
0336     //! do cleanup
0337     if (vetoDeps[iDep]) {
0338       delete vetoDeps[iDep];
0339       vetoDeps[iDep] = nullptr;
0340     }
0341   }
0342 }
0343 
0344 template <typename BT>
0345 unsigned int MuIsolatorResultProducer<BT>::initAssociation(edm::Event& event, CandMap& candMapT) const {
0346   std::string metname = "RecoMuon|MuonIsolationProducers";
0347 
0348   typedef reco::IsoDepositMap::container CT;
0349 
0350   for (unsigned int iMap = 0; iMap < theDepositConfs.size(); ++iMap) {
0351     edm::Handle<reco::IsoDepositMap> depH;
0352     event.getByLabel(theDepositConfs[iMap].tag, depH);
0353     LogDebug(metname) << "Got Deposits of size " << depH->size();
0354     if (depH->empty())
0355       continue;
0356 
0357     //! WARNING: the input ValueMaps are better be for a single key product ID
0358     //! no effort is done (FIXME) for more complex cases
0359     typename CandMap::handle_type keyH;
0360     event.get(depH->begin().id(), keyH);
0361     candMapT.setHandle(keyH);
0362     typename CT::const_iterator depHCI = depH->begin().begin();
0363     typename CT::const_iterator depEnd = depH->begin().end();
0364     unsigned int keyI = 0;
0365     for (; depHCI != depEnd; ++depHCI, ++keyI) {
0366       typename CandMap::key_type muPtr(keyH->refAt(keyI));
0367       //! init {muon, {[deposit,veto]_type}} container
0368       if (iMap == 0)
0369         candMapT.get().push_back(typename CandMap::pair_type(muPtr, DepositContainer(theDepositConfs.size())));
0370       typename CandMap::iterator muI = candMapT.get().begin();
0371       for (; muI != candMapT.get().end(); ++muI) {
0372         if (muI->first == muPtr)
0373           break;
0374       }
0375       if (muI->first != muPtr) {
0376         edm::LogError("MuonIsolationProducers") << "Failed to align muon map";
0377       }
0378       muI->second[iMap].dep = &*depHCI;
0379     }
0380   }
0381 
0382   LogDebug(metname) << "Picked and aligned nDeps = " << candMapT.get().size();
0383   return candMapT.get().size();
0384 }
0385 
0386 template <typename BT>
0387 void MuIsolatorResultProducer<BT>::initVetos(reco::TrackBase::Point const& theBeam,
0388                                              std::vector<reco::IsoDeposit::Vetos*>& vetos,
0389                                              CandMap& candMapT) const {
0390   std::string metname = "RecoMuon|MuonIsolationProducers";
0391 
0392   if (theRemoveOtherVetos) {
0393     LogDebug(metname) << "Start checking for vetos based on input/expected vetos.size of " << vetos.size()
0394                       << " passed at " << &vetos << " and an input map.size of " << candMapT.get().size();
0395 
0396     unsigned int muI = 0;
0397     for (; muI < candMapT.get().size(); ++muI) {
0398       typename CandMap::key_type mu = candMapT.get()[muI].first;
0399       double d0 = ((mu->vx() - theBeam.x()) * mu->py() - (mu->vy() - theBeam.y()) * mu->px()) / mu->pt();
0400       LogDebug(metname) << "Muon at index " << muI;
0401       if (theVetoCuts.selectAll || (fabs(mu->eta()) < theVetoCuts.muAbsEtaMax && mu->pt() > theVetoCuts.muPtMin &&
0402                                     fabs(mu->vz()) < theVetoCuts.muAbsZMax && fabs(d0) < theVetoCuts.muD0Max)) {
0403         LogDebug(metname) << "muon passes the cuts";
0404         for (unsigned int iDep = 0; iDep < candMapT.get()[muI].second.size(); ++iDep) {
0405           if (vetos[iDep] == nullptr)
0406             vetos[iDep] = new reco::IsoDeposit::Vetos();
0407 
0408           vetos[iDep]->push_back(candMapT.get()[muI].second[iDep].dep->veto());
0409         }
0410       }
0411     }
0412 
0413     LogDebug(metname) << "Assigning vetos";
0414     muI = 0;
0415     for (; muI < candMapT.get().size(); ++muI) {
0416       for (unsigned int iDep = 0; iDep < candMapT.get()[muI].second.size(); ++iDep) {
0417         candMapT.get()[muI].second[iDep].vetos = vetos[iDep];
0418       }
0419     }
0420     LogDebug(metname) << "Done with vetos";
0421   }
0422 }
0423 
0424 #endif
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated by the 2.2.1 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!