Back to home page

Project CMSSW displayed by LXR
 
 

    

File indexing completed on 2024-11-28 03:54:26 

0001 #ifndef TOPDQMHELPERS
0002 #define TOPDQMHELPERS
0003 
0004 #include <string>
0005 #include <vector>
0006 #include <iostream>
0007 #include "FWCore/Framework/interface/Event.h"
0008 #include "FWCore/Common/interface/TriggerNames.h"
0009 #include "DataFormats/Common/interface/TriggerResults.h"
0010 #include "DataFormats/VertexReco/interface/Vertex.h"
0011 #include "DataFormats/ParticleFlowCandidate/interface/PFCandidate.h"
0012 #include <DataFormats/METReco/interface/PFMET.h>
0013 #include "FWCore/Framework/interface/ConsumesCollector.h"
0014 #include "FWCore/Framework/interface/EDConsumerBase.h"
0015 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0016 
0017 #include "DataFormats/PatCandidates/interface/Muon.h"
0018 #include "DataFormats/PatCandidates/interface/Electron.h"
0019 #include "DataFormats/PatCandidates/interface/Jet.h"
0020 #include "DataFormats/PatCandidates/interface/MET.h"
0021 
0022 /**
0023    \fn      accept TopDQMHelpers.h "DQM/Physics/interface/TopDQMHelpers.h"
0024 
0025    \brief   Helper function to determine trigger accepts.
0026 
0027    Helper function to determine trigger accept for given TriggerResults and
0028    a given TriggerPath(s).
0029 */
0030 
0031 inline bool accept(const edm::Event& event, const edm::TriggerResults& triggerTable, const std::string& triggerPath) {
0032   bool passed = false;
0033   const edm::TriggerNames& triggerNames = event.triggerNames(triggerTable);
0034   for (unsigned int i = 0; i < triggerNames.triggerNames().size(); ++i) {
0035     if (triggerNames.triggerNames()[i] == triggerPath) {
0036       if (triggerTable.accept(i)) {
0037         passed = true;
0038         break;
0039       }
0040     }
0041   }
0042   return passed;
0043 }
0044 
0045 inline bool accept(const edm::Event& event,
0046                    const edm::TriggerResults& triggerTable,
0047                    const std::vector<std::string>& triggerPaths) {
0048   bool passed = false;
0049   for (unsigned int j = 0; j < triggerPaths.size(); ++j) {
0050     if (accept(event, triggerTable, triggerPaths[j])) {
0051       passed = true;
0052       break;
0053     }
0054   }
0055   return passed;
0056 }
0057 
0058 #include "DataFormats/JetReco/interface/Jet.h"
0059 #include "FWCore/Framework/interface/EventSetup.h"
0060 
0061 /**
0062    \class   Calculate TopDQMHelpers.h "DQM/Physics/interface/TopDQMHelpers.h"
0063 
0064    \brief   Helper class for the calculation of a top and a W boson mass estime.
0065 
0066    Helper class for the calculation of a top and a W boson mass estimate. The
0067    core implementation originates from the plugin TtSemiLepHypMaxSumPtWMass
0068    in TopQuarkAnalysis/TopJetCombination package. It may be extended to include
0069    b tag information.
0070 */
0071 class Calculate_miniAOD {
0072 public:
0073   /// default constructor
0074   Calculate_miniAOD(int maxNJets, double wMass);
0075   /// default destructor
0076   ~Calculate_miniAOD() {}
0077 
0078   /// calculate W boson mass estimate
0079   double massWBoson(const std::vector<pat::Jet>& jets);
0080   /// calculate t-quark mass estimate
0081   double massTopQuark(const std::vector<pat::Jet>& jets);
0082   /// calculate b-tagged t-quark mass estimate
0083   // double massBTopQuark(const std::vector<reco::Jet>& jets, std::vector<bool>
0084   // bjet);
0085   double massBTopQuark(const std::vector<pat::Jet>& jets, std::vector<double> VbtagWP, double btagWP_);
0086 
0087   /// calculate W boson transverse mass estimate
0088   double tmassWBoson(pat::Muon* lep, const pat::MET& met, const pat::Jet& b);
0089   /// calculate top quark transverse mass estimate
0090   double tmassTopQuark(pat::Muon* lep, const pat::MET& met, const pat::Jet& b);
0091 
0092   /// calculate W boson transverse mass estimate
0093   double tmassWBoson(pat::Electron* lep, const pat::MET& met, const pat::Jet& b);
0094   /// calculate top quark transverse mass estimate
0095   double tmassTopQuark(pat::Electron* lep, const pat::MET& met, const pat::Jet& b);
0096 
0097 private:
0098   /// do the calculation; this is called only once per event by the first
0099   /// function call to return a mass estimate. The once calculated values
0100   /// are cached afterwards
0101   void operator()(const std::vector<pat::Jet>& jets);
0102   /// do the calculation of the t-quark mass with one b-jet
0103   void operator2(const std::vector<pat::Jet>&, std::vector<double>, double);
0104   /// do the calculation of the transverse top and W masses
0105   void operator()(const pat::Jet& bJet, pat::Electron* lepton, const pat::MET& met);
0106 
0107   void operator()(const pat::Jet& bJet, pat::Muon* lepton, const pat::MET& met);
0108 
0109 private:
0110   /// indicate failed associations
0111   bool failed_;
0112   /// max. number of jets to be considered
0113   int maxNJets_;
0114   /// paramater of the w boson mass
0115   double wMass_;
0116   /// cache of w boson mass estimate
0117   double massWBoson_;
0118   /// cache of top quark mass estimate
0119   double massTopQuark_;
0120   /// cache of b-tagged top quark mass estimate
0121   double massBTopQuark_;
0122   /// cache of W boson transverse mass estimate
0123   double tmassWBoson_;
0124   /// cache of top quark transverse mass estimate
0125   double tmassTopQuark_;
0126 };
0127 
0128 class Calculate {
0129 public:
0130   /// default constructor
0131   Calculate(int maxNJets, double wMass);
0132   /// default destructor
0133   ~Calculate() {}
0134 
0135   /// calculate W boson mass estimate
0136   double massWBoson(const std::vector<reco::Jet>& jets);
0137   /// calculate t-quark mass estimate
0138   double massTopQuark(const std::vector<reco::Jet>& jets);
0139   /// calculate b-tagged t-quark mass estimate
0140   // double massBTopQuark(const std::vector<reco::Jet>& jets, std::vector<bool>
0141   // bjet);
0142   double massBTopQuark(const std::vector<reco::Jet>& jets, std::vector<double> VbtagWP, double btagWP_);
0143 
0144   /// calculate W boson transverse mass estimate
0145   double tmassWBoson(reco::RecoCandidate* lep, const reco::MET& met, const reco::Jet& b);
0146   /// calculate top quark transverse mass estimate
0147   double tmassTopQuark(reco::RecoCandidate* lep, const reco::MET& met, const reco::Jet& b);
0148 
0149 private:
0150   /// do the calculation; this is called only once per event by the first
0151   /// function call to return a mass estimate. The once calculated values
0152   /// are cached afterwards
0153   void operator()(const std::vector<reco::Jet>& jets);
0154   /// do the calculation of the t-quark mass with one b-jet
0155   void operator2(const std::vector<reco::Jet>&, std::vector<double>, double);
0156   /// do the calculation of the transverse top and W masses
0157   void operator()(const reco::Jet& bJet, reco::RecoCandidate* lepton, const reco::MET& met);
0158 
0159 private:
0160   /// indicate failed associations
0161   bool failed_;
0162   /// max. number of jets to be considered
0163   int maxNJets_;
0164   /// paramater of the w boson mass
0165   double wMass_;
0166   /// cache of w boson mass estimate
0167   double massWBoson_;
0168   /// cache of top quark mass estimate
0169   double massTopQuark_;
0170   /// cache of b-tagged top quark mass estimate
0171   double massBTopQuark_;
0172   /// cache of W boson transverse mass estimate
0173   double tmassWBoson_;
0174   /// cache of top quark transverse mass estimate
0175   double tmassTopQuark_;
0176 };
0177 
0178 #include "DataFormats/JetReco/interface/JetID.h"
0179 #include "DataFormats/JetReco/interface/PFJet.h"
0180 #include "DataFormats/JetReco/interface/CaloJet.h"
0181 #include "DataFormats/BTauReco/interface/JetTag.h"
0182 #include "DataFormats/Common/interface/Handle.h"
0183 #include "DataFormats/Common/interface/ValueMap.h"
0184 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0185 #include "JetMETCorrections/JetCorrector/interface/JetCorrector.h"
0186 #include "DataFormats/EgammaCandidates/interface/GsfElectron.h"
0187 #include "CommonTools/Utils/interface/StringCutObjectSelector.h"
0188 #include "FWCore/Utilities/interface/EDGetToken.h"
0189 #include "FWCore/Utilities/interface/InputTag.h"
0190 
0191 /**
0192    \class   SelectionStep TopDQMHelpers.h
0193    "DQM/Physics/interface/TopDQMHelpers.h"
0194 
0195    \brief   Templated helper class to allow a selection on a certain object
0196    collection.
0197 
0198    Templated helper class to allow a selection on a certain object collection,
0199    which may be monitored by a separate class afterwards. The class wraps and
0200    slightly extends the features of the StringCutParser to allow also to apply
0201    event based selections, according to a minimal or maximal number of elements
0202    in the collection after the object selection has been applied. It takes an
0203    edm::ParameterSet in the constructor, which should contain the following
0204    elements:
0205     - src          : the input collection (mandatory).
0206     - select       : the selection string (mandatory).
0207     - min          : whether there is a min value on which to reject the whole
0208    event after
0209                      the selection (optional).
0210     - max          : whether there is a max value on which to reject the whole
0211    event after
0212                      the selection (optional).
0213     - electronId   : input tag of an electronId association map and cut value
0214                      (optional).
0215     - jetCorrector : label of jet corrector (optional).
0216     - jetBTagger   : parameters defining the btag algorithm and working point of
0217    choice
0218                      (optional).
0219     - jetID        : parameters defining the jetID value map and selection
0220    (optional).
0221 
0222 
0223    The parameters _src_ and _select_ are mandatory. The parameters _min_ and
0224    _max_ are optional. The parameters _electronId_ and _jetCorrector_ are
0225    optional. They are added to keep the possibility to apply selections on id'ed
0226    electrons or on corrected jets.  They may be omitted in the PSet for
0227    simplification reasons if not needed at any time.  They are not effiective
0228    for other object collections but electrons or jets.  If none of the two
0229    parameters _min_ or _max_ is found in the event the select function returns
0230    true if at least one object fullfilled the requirements.
0231 
0232    The class has one template value, which is the object collection to apply the
0233    selection on. This has to be parsed to the StringCutParser class. The
0234    function select is overrided for jets to circumvent problems with the
0235    template specialisation. Note that for MET not type1 or muon corrections are
0236    supported on reco candidates.
0237 */
0238 
0239 template <typename Object>
0240 class SelectionStep {
0241 public:
0242   /// default constructor
0243   SelectionStep(const edm::ParameterSet& cfg, edm::ConsumesCollector&& iC);
0244   /// default destructor
0245   ~SelectionStep() {}
0246 
0247   /// apply selection
0248   bool select(const edm::Event& event);
0249   bool select(const edm::Event& event, const std::string& type);
0250   /// apply selection override for jets
0251   bool select(const edm::Event& event, const edm::EventSetup& setup);
0252   bool selectVertex(const edm::Event& event);
0253 
0254 private:
0255   /// input collection
0256   edm::EDGetTokenT<edm::View<Object> > src_;
0257   /// min/max for object multiplicity
0258   int min_, max_;
0259   /// electronId label as extra selection type
0260   edm::EDGetTokenT<edm::ValueMap<float> > electronId_;
0261   /// electronId pattern we expect the following pattern:
0262   ///  0: fails
0263   ///  1: passes electron ID only
0264   ///  2: passes electron Isolation only
0265   ///  3: passes electron ID and Isolation only
0266   ///  4: passes conversion rejection
0267   ///  5: passes conversion rejection and ID
0268   ///  6: passes conversion rejection and Isolation
0269   ///  7: passes the whole selection
0270   /// As described on
0271   /// https://twiki.cern.ch/twiki/bin/view/CMS/SimpleCutBasedEleID
0272   double eidCutValue_;
0273   /// jet corrector as extra selection type
0274   edm::EDGetTokenT<reco::JetCorrector> jetCorrector_;
0275   /// choice for b-tag as extra selection type
0276   edm::EDGetTokenT<reco::JetTagCollection> btagLabel_;
0277   /// choice of b-tag working point as extra selection type
0278   double btagWorkingPoint_;
0279   /// jetID as an extra selection type
0280   edm::EDGetTokenT<reco::JetIDValueMap> jetIDLabel_;
0281 
0282   edm::EDGetTokenT<edm::View<reco::Vertex> > pvs_;
0283   edm::EDGetTokenT<edm::View<reco::GsfElectron> > gsfEs_;
0284   /// string cut selector
0285   StringCutObjectSelector<Object> select_;
0286   /// selection string on the jetID
0287   StringCutObjectSelector<reco::JetID>* jetIDSelect_;
0288 };
0289 
0290 /// default constructor
0291 template <typename Object>
0292 SelectionStep<Object>::SelectionStep(const edm::ParameterSet& cfg, edm::ConsumesCollector&& iC)
0293     : select_(cfg.getParameter<std::string>("select")), jetIDSelect_(nullptr) {
0294   src_ = iC.consumes<edm::View<Object> >(cfg.getParameter<edm::InputTag>("src"));
0295   // exist otherwise they are initialized with -1
0296   cfg.exists("min") ? min_ = cfg.getParameter<int>("min") : min_ = -1;
0297   cfg.exists("max") ? max_ = cfg.getParameter<int>("max") : max_ = -1;
0298   std::string mygSF = "gedGsfElectrons";
0299   gsfEs_ = iC.consumes<edm::View<reco::GsfElectron> >(cfg.getUntrackedParameter<edm::InputTag>("myGSF", mygSF));
0300   if (cfg.existsAs<edm::ParameterSet>("electronId")) {
0301     edm::ParameterSet elecId = cfg.getParameter<edm::ParameterSet>("electronId");
0302     electronId_ = iC.consumes<edm::ValueMap<float> >(elecId.getParameter<edm::InputTag>("src"));
0303     eidCutValue_ = elecId.getParameter<double>("cutValue");
0304   }
0305   if (cfg.exists("jetCorrector")) {
0306     jetCorrector_ = iC.consumes<reco::JetCorrector>(edm::InputTag(cfg.getParameter<std::string>("jetCorrector")));
0307   }
0308   if (cfg.existsAs<edm::ParameterSet>("jetBTagger")) {
0309     edm::ParameterSet jetBTagger = cfg.getParameter<edm::ParameterSet>("jetBTagger");
0310     btagLabel_ = iC.consumes<reco::JetTagCollection>(jetBTagger.getParameter<edm::InputTag>("label"));
0311     btagWorkingPoint_ = jetBTagger.getParameter<double>("workingPoint");
0312   }
0313   if (cfg.existsAs<edm::ParameterSet>("jetID")) {
0314     edm::ParameterSet jetID = cfg.getParameter<edm::ParameterSet>("jetID");
0315     jetIDLabel_ = iC.consumes<reco::JetIDValueMap>(jetID.getParameter<edm::InputTag>("label"));
0316     jetIDSelect_ = new StringCutObjectSelector<reco::JetID>(jetID.getParameter<std::string>("select"));
0317   }
0318 }
0319 
0320 /// apply selection
0321 template <typename Object>
0322 bool SelectionStep<Object>::select(const edm::Event& event) {
0323   // fetch input collection
0324   edm::Handle<edm::View<Object> > src;
0325   if (!event.getByToken(src_, src))
0326     return false;
0327 
0328   // load electronId value map if configured such
0329   edm::Handle<edm::ValueMap<float> > electronId;
0330   if (!electronId_.isUninitialized()) {
0331     if (!event.getByToken(electronId_, electronId))
0332       return false;
0333   }
0334 
0335   // determine multiplicity of selected objects
0336   int n = 0;
0337   for (typename edm::View<Object>::const_iterator obj = src->begin(); obj != src->end(); ++obj) {
0338     // special treatment for electrons
0339     if (dynamic_cast<const reco::GsfElectron*>(&*obj)) {
0340       unsigned int idx = obj - src->begin();
0341       if (electronId_.isUninitialized() ? true : ((double)(*electronId)[src->refAt(idx)] >= eidCutValue_)) {
0342         if (select_(*obj))
0343           ++n;
0344       }
0345     }
0346     // normal treatment
0347     else {
0348       if (select_(*obj))
0349         ++n;
0350     }
0351   }
0352   bool accept = (min_ >= 0 ? n >= min_ : true) && (max_ >= 0 ? n <= max_ : true);
0353   return (min_ < 0 && max_ < 0) ? (n > 0) : accept;
0354 }
0355 
0356 /// apply selection with special treatment for PFCandidates
0357 template <typename Object>
0358 bool SelectionStep<Object>::select(const edm::Event& event, const std::string& type) {
0359   // fetch input collection
0360   edm::Handle<edm::View<Object> > src;
0361   if (!event.getByToken(src_, src))
0362     return false;
0363 
0364   // special for gsfElectron
0365   edm::Handle<edm::View<reco::GsfElectron> > elecs_gsf;
0366 
0367   // load electronId value map if configured such
0368   edm::Handle<edm::ValueMap<float> > electronId;
0369   if (!electronId_.isUninitialized()) {
0370     if (!event.getByToken(electronId_, electronId))
0371       return false;
0372   }
0373 
0374   // determine multiplicity of selected objects
0375   int n = 0;
0376   for (typename edm::View<Object>::const_iterator obj = src->begin(); obj != src->end(); ++obj) {
0377     // special treatment for PF candidates
0378     if (dynamic_cast<const reco::PFCandidate*>(&*obj)) {
0379       reco::PFCandidate objtmp = dynamic_cast<const reco::PFCandidate&>(*obj);
0380 
0381       if (type == "muon") {
0382         if (select_(*obj)) {
0383           ++n;
0384         }
0385       } else if (type == "electron") {
0386         if (select_(*obj)) {
0387           ++n;
0388         }
0389         //        idx_gsf++;
0390       }
0391     }
0392 
0393     // special treatment for electrons
0394     else if (dynamic_cast<const reco::GsfElectron*>(&*obj)) {
0395       unsigned int idx = obj - src->begin();
0396       if (electronId_.isUninitialized() ? true : ((double)(*electronId)[src->refAt(idx)] >= eidCutValue_)) {
0397         if (select_(*obj))
0398           ++n;
0399       }
0400     }
0401 
0402     // normal treatment
0403     else {
0404       if (select_(*obj))
0405         ++n;
0406     }
0407   }
0408   bool accept = (min_ >= 0 ? n >= min_ : true) && (max_ >= 0 ? n <= max_ : true);
0409   return (min_ < 0 && max_ < 0) ? (n > 0) : accept;
0410 }
0411 
0412 template <typename Object>
0413 bool SelectionStep<Object>::selectVertex(const edm::Event& event) {
0414   // fetch input collection
0415   edm::Handle<edm::View<Object> > src;
0416   if (!event.getByToken(src_, src))
0417     return false;
0418 
0419   // load electronId value map if configured such
0420   edm::Handle<edm::ValueMap<float> > electronId;
0421   if (!electronId_.isUninitialized()) {
0422     if (!event.getByToken(electronId_, electronId))
0423       return false;
0424   }
0425 
0426   // determine multiplicity of selected objects
0427   int n = 0;
0428   for (typename edm::View<Object>::const_iterator obj = src->begin(); obj != src->end(); ++obj) {
0429     if (select_(*obj))
0430       ++n;
0431   }
0432   bool accept = (min_ >= 0 ? n >= min_ : true) && (max_ >= 0 ? n <= max_ : true);
0433   return (min_ < 0 && max_ < 0) ? (n > 0) : accept;
0434 }
0435 
0436 /// apply selection (w/o using the template class Object), override for jets
0437 template <typename Object>
0438 bool SelectionStep<Object>::select(const edm::Event& event, const edm::EventSetup& setup) {
0439   // fetch input collection
0440   edm::Handle<edm::View<Object> > src;
0441   if (!event.getByToken(src_, src))
0442     return false;
0443 
0444   // load btag collection if configured such
0445   // NOTE that the JetTagCollection needs an
0446   // edm::View to reco::Jets; we have to add
0447   // another Handle bjets for this purpose
0448   edm::Handle<edm::View<reco::Jet> > bjets;
0449   edm::Handle<reco::JetTagCollection> btagger;
0450   edm::Handle<edm::View<reco::Vertex> > pvertex;
0451   if (!btagLabel_.isUninitialized()) {
0452     if (!event.getByToken(src_, bjets))
0453       return false;
0454     if (!event.getByToken(btagLabel_, btagger))
0455       return false;
0456     if (!event.getByToken(pvs_, pvertex))
0457       return false;
0458   }
0459 
0460   // load jetID value map if configured such
0461   edm::Handle<reco::JetIDValueMap> jetID;
0462   if (jetIDSelect_) {
0463     if (!event.getByToken(jetIDLabel_, jetID))
0464       return false;
0465   }
0466 
0467   // load jet corrector if configured such
0468   const reco::JetCorrector* corrector = nullptr;
0469   if (!jetCorrector_.isUninitialized()) {
0470     // check whether a jet corrector is in the event or not
0471     edm::Handle<reco::JetCorrector> correctorHandle = event.getHandle(jetCorrector_);
0472     if (correctorHandle.isValid()) {
0473       corrector = correctorHandle.product();
0474     } else {
0475       edm::LogVerbatim("TopDQMHelpers") << "\n"
0476                                         << "---------------------------------------------------------------\n"
0477                                         << " No reco::JetCorrrector available from Event:\n"
0478                                         << "  - Jets will not be corrected.\n"
0479                                         << "---------------------------------------------------------------"
0480                                            "\n";
0481     }
0482   }
0483   // determine multiplicity of selected objects
0484   int n = 0;
0485   for (typename edm::View<Object>::const_iterator obj = src->begin(); obj != src->end(); ++obj) {
0486     // check for chosen btag discriminator to be above the
0487     // corresponding working point if configured such
0488     unsigned int idx = obj - src->begin();
0489     if (btagLabel_.isUninitialized() ? true : (*btagger)[bjets->refAt(idx)] > btagWorkingPoint_) {
0490       bool passedJetID = true;
0491       // check jetID for calo jets
0492       if (jetIDSelect_ && dynamic_cast<const reco::CaloJet*>(src->refAt(idx).get())) {
0493         passedJetID = (*jetIDSelect_)((*jetID)[src->refAt(idx)]);
0494       }
0495       if (passedJetID) {
0496         // scale jet energy if configured such
0497         Object jet = *obj;
0498         jet.scaleEnergy(corrector ? corrector->correction(*obj) : 1.);
0499         if (select_(jet))
0500           ++n;
0501       }
0502     }
0503   }
0504   bool accept = (min_ >= 0 ? n >= min_ : true) && (max_ >= 0 ? n <= max_ : true);
0505   return (min_ < 0 && max_ < 0) ? (n > 0) : accept;
0506 }
0507 
0508 #endif
0509 
0510 /* Local Variables: */
0511 /* show-trailing-whitespace: t */
0512 /* truncate-lines: t */
0513 /* End: */
[ 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!