TopDQMHelpers.h

Go to the documentation of this file.

#ifndef TOPDQMHELPERS
#define TOPDQMHELPERS

#include <string>
#include <vector>
#include <iostream>
#include "FWCore/Framework/interface/Event.h"
#include "FWCore/Common/interface/TriggerNames.h"
#include "DataFormats/Common/interface/TriggerResults.h"
#include "DataFormats/VertexReco/interface/Vertex.h"
#include "DataFormats/ParticleFlowCandidate/interface/PFCandidate.h"
#include <DataFormats/METReco/interface/PFMET.h>
#include "FWCore/Framework/interface/ConsumesCollector.h"
#include "FWCore/Framework/interface/EDConsumerBase.h"
using namespace std;
inline bool
accept(const edm::Event& event, const edm::TriggerResults& triggerTable, const std::string& triggerPath)
{
  bool passed=false;
  const edm::TriggerNames& triggerNames = event.triggerNames(triggerTable);
  for(unsigned int i=0; i<triggerNames.triggerNames().size(); ++i){
    if(triggerNames.triggerNames()[i] == triggerPath) {
      if(triggerTable.accept(i)){
    passed=true;
    break;
      }
    }
  }
  return passed;
}

inline bool
accept(const edm::Event& event, const edm::TriggerResults& triggerTable, const std::vector<std::string>& triggerPaths)
{
  bool passed=false;
  for(unsigned int j=0; j<triggerPaths.size(); ++j){
    if(accept(event, triggerTable, triggerPaths[j])){
      passed=true;
      break;
    }
  }
  return passed;
}


#include "DataFormats/JetReco/interface/Jet.h"
#include "FWCore/Framework/interface/EventSetup.h"

class Calculate {
 public:
  Calculate(int maxNJets, double wMass);
  ~Calculate(){};

  double massWBoson(const std::vector<reco::Jet>& jets);
  double massTopQuark(const std::vector<reco::Jet>& jets);
  //double massBTopQuark(const std::vector<reco::Jet>& jets, std::vector<bool> bjet);
  double massBTopQuark(const std::vector<reco::Jet>& jets, std::vector<double> VbtagWP, double btagWP_);

  double tmassWBoson(reco::RecoCandidate* lep, const reco::MET& met, const reco::Jet& b);
  double tmassTopQuark(reco::RecoCandidate* lep, const reco::MET& met, const reco::Jet& b);

 private:
  void operator()(const std::vector<reco::Jet>& jets);
  void operator2(const std::vector<reco::Jet>& , std::vector<double> , double);
  void operator()(const reco::Jet& bJet, reco::RecoCandidate* lepton, const reco::MET& met);
  

 private:
  bool failed_;
  int maxNJets_;
  double wMass_;
  double massWBoson_;
  double massTopQuark_;
  double massBTopQuark_;
  double tmassWBoson_;
  double tmassTopQuark_;
  


};


#include "DataFormats/JetReco/interface/JetID.h"
#include "DataFormats/JetReco/interface/PFJet.h"
#include "DataFormats/JetReco/interface/CaloJet.h"
#include "DataFormats/BTauReco/interface/JetTag.h"
#include "DataFormats/Common/interface/ValueMap.h"
#include "FWCore/ParameterSet/interface/ParameterSet.h"
#include "JetMETCorrections/Objects/interface/JetCorrector.h"
#include "DataFormats/EgammaCandidates/interface/GsfElectron.h"
#include "CommonTools/Utils/interface/StringCutObjectSelector.h"
#include "JetMETCorrections/Objects/interface/JetCorrectionsRecord.h"

template <typename Object>
class SelectionStep {
public:
  SelectionStep(const edm::ParameterSet& cfg, edm::ConsumesCollector && iC);
  ~SelectionStep(){};

  bool select(const edm::Event& event);
  bool select(const edm::Event& event, const std::string& type );
  bool select(const edm::Event& event, const edm::EventSetup& setup);
  bool selectVertex(const edm::Event& event);
private:
  edm::EDGetTokenT<edm::View<Object> > src_;
  int min_, max_;
  edm::EDGetTokenT<edm::ValueMap<float> > electronId_;
  int eidPattern_;
  std::string jetCorrector_;
  edm::EDGetTokenT<reco::JetTagCollection > btagLabel_;
  double btagWorkingPoint_;
  edm::EDGetTokenT<reco::JetIDValueMap> jetIDLabel_;

  edm::EDGetTokenT<edm::View<reco::Vertex> > pvs_; 
  edm::EDGetTokenT<edm::View<reco::GsfElectron> > gsfEs_; 
  StringCutObjectSelector<Object> select_;
  StringCutObjectSelector<reco::JetID>* jetIDSelect_;

};

template <typename Object> 
SelectionStep<Object>::SelectionStep(const edm::ParameterSet& cfg, edm::ConsumesCollector && iC): 
  select_(cfg.getParameter<std::string>("select")),
  jetIDSelect_( 0)
{
  //cout<<"In new constructor"<<endl;
  src_ = iC.consumes<edm::View<Object> >(cfg.getParameter<edm::InputTag>("src"));
  //cout<<"// construct min/max if the corresponding params"<<endl;
  // exist otherwise they are initialized with -1
  cfg.exists("min") ? min_= cfg.getParameter<int>("min") : min_= -1;
  cfg.exists("max") ? max_= cfg.getParameter<int>("max") : max_= -1;
  //cout<<"// read electron extras if they exist"<<endl;
  std::string mygSF = "gedGsfElectrons";
  gsfEs_ = iC.consumes<edm::View<reco::GsfElectron> >(cfg.getUntrackedParameter<edm::InputTag>("myGSF", mygSF));
  if(cfg.existsAs<edm::ParameterSet>("electronId")){ 
    edm::ParameterSet elecId=cfg.getParameter<edm::ParameterSet>("electronId");
    electronId_= iC.consumes<edm::ValueMap<float> >(elecId.getParameter<edm::InputTag>("src"));
    eidPattern_= elecId.getParameter<int>("pattern");
  }
  //cout<<"// read jet corrector label if it exists"<<endl;
  if(cfg.exists("jetCorrector")){ jetCorrector_= cfg.getParameter<std::string>("jetCorrector"); }
  //cout<<"// read btag information if it exists"<<endl;
  if(cfg.existsAs<edm::ParameterSet>("jetBTagger")){
    edm::ParameterSet jetBTagger=cfg.getParameter<edm::ParameterSet>("jetBTagger");
    btagLabel_=iC.consumes<reco::JetTagCollection>(jetBTagger.getParameter<edm::InputTag>("label"));
    btagWorkingPoint_=jetBTagger.getParameter<double>("workingPoint");
  }
  //cout<<"// read jetID information if it exists"<<endl;
  if(cfg.existsAs<edm::ParameterSet>("jetID")){
    edm::ParameterSet jetID=cfg.getParameter<edm::ParameterSet>("jetID");
    jetIDLabel_ =iC.consumes<reco::JetIDValueMap>(jetID.getParameter<edm::InputTag>("label"));
    jetIDSelect_= new StringCutObjectSelector<reco::JetID>(jetID.getParameter<std::string>("select"));
  }
  //cout<<"// end"<<endl;
}

template <typename Object>
bool SelectionStep<Object>::select(const edm::Event& event)
{
  // fetch input collection
  edm::Handle<edm::View<Object> > src; 
  if( !event.getByToken(src_, src) ) return false;

  // load electronId value map if configured such
  edm::Handle<edm::ValueMap<float> > electronId;
  if(!electronId_.isUninitialized()){
    if( !event.getByToken(electronId_, electronId) ) return false;
  }

  // determine multiplicity of selected objects
  int n=0;
  for(typename edm::View<Object>::const_iterator obj=src->begin(); obj!=src->end(); ++obj){
    // special treatment for electrons
    if(dynamic_cast<const reco::GsfElectron*>(&*obj)){
      unsigned int idx = obj-src->begin();
      if( electronId_.isUninitialized() ? true : ((int)(*electronId)[src->refAt(idx)] & eidPattern_) ){   
    if(select_(*obj))++n;
      }
    }
    // normal treatment
    else{
      if(select_(*obj))++n;
    }
  }
  bool accept=(min_>=0 ? n>=min_:true) && (max_>=0 ? n<=max_:true);
  return (min_<0 && max_<0) ? (n>0):accept;
}


template <typename Object>
bool SelectionStep<Object>::select(const edm::Event& event, const std::string& type)
{
  // fetch input collection
  edm::Handle<edm::View<Object> > src;
  if( !event.getByToken(src_, src) ) return false;
  
  // special for gsfElectron
  edm::Handle<edm::View<reco::GsfElectron> > elecs_gsf;
  
  // load electronId value map if configured such
  edm::Handle<edm::ValueMap<float> > electronId;
  if(!electronId_.isUninitialized()){
    if( !event.getByToken(electronId_, electronId) ) return false;
  }
  
  // determine multiplicity of selected objects
  int n=0;
  unsigned int idx_gsf = 0;
  for(typename edm::View<Object>::const_iterator obj=src->begin(); obj!=src->end(); ++obj){
    
    
    // special treatment for PF candidates
    if (dynamic_cast<const reco::PFCandidate*>(&*obj)){
      reco::PFCandidate objtmp = dynamic_cast<const reco::PFCandidate&>(*obj);
      
      if (objtmp.muonRef().isNonnull() && type == "muon") {

        if(select_(*obj)){++n;
        }
      }
      else if (objtmp.gsfElectronRef().isNonnull() && type == "electron") {
        if( !event.getByToken(gsfEs_, elecs_gsf) ) continue;
        if(select_(*obj)){
          if(elecs_gsf->refAt(idx_gsf).isNonnull()){
        int eID=0;
        if (!electronId_.isUninitialized()){
          eID = (int)(*electronId)[elecs_gsf->refAt(idx_gsf)];
        }
            if( electronId_.isUninitialized() ? true : ( (eID & eidPattern_)  && (eID >= 5) ) )
              ++n;
        
          }
        }
        idx_gsf++;
      }
    }
    
    // special treatment for electrons
    else if(dynamic_cast<const reco::GsfElectron*>(&*obj)){
      unsigned int idx = obj-src->begin();
      int eID = (int)(*electronId)[src->refAt(idx)];
      if( electronId_.isUninitialized() ? true : ( (eID & eidPattern_)  && (eID >= 5) ) ){
        if(select_(*obj))++n;
      }
    }
    
    
    // normal treatment
    else{
      if(select_(*obj))++n;
    }
  }
  bool accept=(min_>=0 ? n>=min_:true) && (max_>=0 ? n<=max_:true);
  return (min_<0 && max_<0) ? (n>0):accept;
}


template <typename Object> 
bool SelectionStep<Object>::selectVertex(const edm::Event& event)
{
  // fetch input collection
  edm::Handle<edm::View<Object> > src; 
  if( !event.getByToken(src_, src) ) return false;

  // load electronId value map if configured such
  edm::Handle<edm::ValueMap<float> > electronId;
  if(!electronId_.isUninitialized()) {
    if( !event.getByToken(electronId_, electronId) ) return false;
  }

  // determine multiplicity of selected objects
  int n=0;
  for(typename edm::View<Object>::const_iterator obj=src->begin(); obj!=src->end(); ++obj){

      if(select_(*obj))++n;
  }
  bool accept=(min_>=0 ? n>=min_:true) && (max_>=0 ? n<=max_:true);
  return (min_<0 && max_<0) ? (n>0):accept;
}

template <typename Object>
bool SelectionStep<Object>::select(const edm::Event& event, const edm::EventSetup& setup)
{
  // fetch input collection
  edm::Handle<edm::View<Object> > src; 
  if( !event.getByToken(src_, src) ) return false;

  // load btag collection if configured such
  // NOTE that the JetTagCollection needs an
  // edm::View to reco::Jets; we have to add
  // another Handle bjets for this purpose
  edm::Handle<edm::View<reco::Jet> > bjets;
  edm::Handle<reco::JetTagCollection> btagger;
  edm::Handle<edm::View<reco::Vertex> > pvertex; 
  if(!btagLabel_.isUninitialized()){ 
    if( !event.getByToken(src_, bjets) ) return false;
    if( !event.getByToken(btagLabel_, btagger) ) return false;
    if( !event.getByToken(pvs_, pvertex) ) return false;
  }

  // load jetID value map if configured such
  edm::Handle<reco::JetIDValueMap> jetID;
  if(jetIDSelect_){
    if( !event.getByToken(jetIDLabel_, jetID) ) return false;

  }

  // load jet corrector if configured such
  const JetCorrector* corrector=0;
  if(!jetCorrector_.empty()){
    // check whether a jet correcto is in the event setup or not
    if(setup.find( edm::eventsetup::EventSetupRecordKey::makeKey<JetCorrectionsRecord>() )){
      corrector = JetCorrector::getJetCorrector(jetCorrector_, setup);
    }
    else{
      edm::LogVerbatim( "TopDQMHelpers" )
        << "\n"
        << "------------------------------------------------------------------------------------- \n"
        << " No JetCorrectionsRecord available from EventSetup:                                   \n"
        << "  - Jets will not be corrected.                                                       \n"
        << "  - If you want to change this add the following lines to your cfg file               \n"
        << "                                                                                      \n"
        << "  ## load jet corrections                                                             \n"
        << "  process.load(\"JetMETCorrections.Configuration.JetCorrectionServicesAllAlgos_cff\") \n"
        << "  process.prefer(\"ak5CaloL2L3\")                                                     \n"
        << "                                                                                      \n"
        << "------------------------------------------------------------------------------------- \n";
    }
  }
  // determine multiplicity of selected objects
  int n=0;
  for(typename edm::View<Object>::const_iterator obj=src->begin(); obj!=src->end(); ++obj){
    // check for chosen btag discriminator to be above the
    // corresponding working point if configured such
    unsigned int idx = obj-src->begin();
    if( btagLabel_.isUninitialized() ? true : (*btagger)[bjets->refAt(idx)]>btagWorkingPoint_ ){   
      bool passedJetID=true;
      // check jetID for calo jets
      if( jetIDSelect_ && dynamic_cast<const reco::CaloJet*>(src->refAt(idx).get())){
    passedJetID=(*jetIDSelect_)((*jetID)[src->refAt(idx)]);
      }
      if(passedJetID){
    // scale jet energy if configured such
    Object jet=*obj; jet.scaleEnergy(corrector ? corrector->correction(*obj) : 1.);
    if(select_(jet))++n;
      }
    }
  }
  bool accept=(min_>=0 ? n>=min_:true) && (max_>=0 ? n<=max_:true);
  return (min_<0 && max_<0) ? (n>0):accept;
}








#endif

/* Local Variables: */
/* show-trailing-whitespace: t */
/* truncate-lines: t */
/* End: */