diff --git a/MiniAODHelper/data/leptonSF/MuonID_Z_RunBCD_prompt80X_7p65.root b/MiniAODHelper/data/leptonSF/MuonID_Z_RunBCD_prompt80X_7p65.root new file mode 100644 index 0000000..47cc078 Binary files /dev/null and b/MiniAODHelper/data/leptonSF/MuonID_Z_RunBCD_prompt80X_7p65.root differ diff --git a/MiniAODHelper/data/leptonSF/MuonIso_Z_RunBCD_prompt80X_7p65.root b/MiniAODHelper/data/leptonSF/MuonIso_Z_RunBCD_prompt80X_7p65.root new file mode 100644 index 0000000..d2478eb Binary files /dev/null and b/MiniAODHelper/data/leptonSF/MuonIso_Z_RunBCD_prompt80X_7p65.root differ diff --git a/MiniAODHelper/data/leptonSF/SingleMuonTrigger_Z_RunBCD_prompt80X_7p65.root b/MiniAODHelper/data/leptonSF/SingleMuonTrigger_Z_RunBCD_prompt80X_7p65.root new file mode 100644 index 0000000..26dbe6d Binary files /dev/null and b/MiniAODHelper/data/leptonSF/SingleMuonTrigger_Z_RunBCD_prompt80X_7p65.root differ diff --git a/MiniAODHelper/interface/MiniAODHelper.h b/MiniAODHelper/interface/MiniAODHelper.h index 3803225..a97a32b 100644 --- a/MiniAODHelper/interface/MiniAODHelper.h +++ b/MiniAODHelper/interface/MiniAODHelper.h @@ -48,6 +48,8 @@ #include "DataFormats/PatCandidates/interface/Isolation.h" #include "DataFormats/PatCandidates/interface/Tau.h" #include "DataFormats/PatCandidates/interface/Particle.h" +#include "DataFormats/JetReco/interface/GenJetCollection.h" +#include "DataFormats/JetReco/interface/GenJet.h" #include "MiniAOD/BoostedObjects/interface/BoostedJet.h" #include "PhysicsTools/SelectorUtils/interface/JetIDSelectionFunctor.h" @@ -171,13 +173,13 @@ class MiniAODHelper{ std::vector GetSelectedJets(const std::vector&, const float, const float, const jetID::jetID, const char); std::vector GetUncorrectedJets(const std::vector&); std::vector GetUncorrectedJets(edm::Handle); - pat::Jet GetCorrectedJet(const pat::Jet&, const edm::Event&, const edm::EventSetup&, const sysType::sysType iSysType=sysType::NA, const bool& doJES=true, const bool& doJER=true, const float& corrFactor = 1, const float& uncFactor = 1); - float GetJetCorrectionFactor(const pat::Jet&, const edm::Event&, const edm::EventSetup&, const sysType::sysType iSysType=sysType::NA, const bool& doJES=true, const bool& doJER=true, const float& corrFactor = 1, const float& uncFactor = 1); - pat::Jet GetCorrectedAK8Jet(const pat::Jet&, const edm::Event&, const edm::EventSetup&, const sysType::sysType iSysType=sysType::NA, const bool& doJES=true, const bool& doJER=true, const float& corrFactor = 1, const float& uncFactor = 1); - float GetAK8JetCorrectionFactor(const pat::Jet&, const edm::Event&, const edm::EventSetup&, const sysType::sysType iSysType=sysType::NA, const bool& doJES=true, const bool& doJER=true, const float& corrFactor = 1, const float& uncFactor = 1); - std::vector GetCorrectedJets(const std::vector&, const edm::Event&, const edm::EventSetup&, const sysType::sysType iSysType=sysType::NA, const bool& doJES=true, const bool& doJER=true, const float& corrFactor = 1, const float& uncFactor = 1); + pat::Jet GetCorrectedJet(const pat::Jet&, const edm::Event&, const edm::EventSetup&, const edm::Handle&, const sysType::sysType iSysType=sysType::NA, const bool& doJES=true, const bool& doJER=true, const float& corrFactor = 1, const float& uncFactor = 1); + float GetJetCorrectionFactor(const pat::Jet&, const edm::Event&, const edm::EventSetup&, const edm::Handle&, const sysType::sysType iSysType=sysType::NA, const bool& doJES=true, const bool& doJER=true, const float& corrFactor = 1, const float& uncFactor = 1); + pat::Jet GetCorrectedAK8Jet(const pat::Jet&, const edm::Event&, const edm::EventSetup&, const edm::Handle&, const sysType::sysType iSysType=sysType::NA, const bool& doJES=true, const bool& doJER=true, const float& corrFactor = 1, const float& uncFactor = 1); + float GetAK8JetCorrectionFactor(const pat::Jet&, const edm::Event&, const edm::EventSetup&, const edm::Handle&, const sysType::sysType iSysType=sysType::NA, const bool& doJES=true, const bool& doJER=true, const float& corrFactor = 1, const float& uncFactor = 1); + std::vector GetCorrectedJets(const std::vector&, const edm::Event&, const edm::EventSetup&, const edm::Handle&, const sysType::sysType iSysType=sysType::NA, const bool& doJES=true, const bool& doJER=true, const float& corrFactor = 1, const float& uncFactor = 1); std::vector GetCorrectedJets(const std::vector&, const sysType::sysType iSysType=sysType::NA); - std::vector GetCorrectedBoostedJets(const std::vector& inputBoostedJets, const edm::Event&, const edm::EventSetup&, const sysType::sysType iSysType=sysType::NA, const bool& doJES=true, const bool& doJER=true, const float& corrFactor = 1, const float& uncFactor = 1); + std::vector GetCorrectedBoostedJets(const std::vector& inputBoostedJets, const edm::Event&, const edm::EventSetup&, const edm::Handle&, const sysType::sysType iSysType=sysType::NA, const bool& doJES=true, const bool& doJER=true, const float& corrFactor = 1, const float& uncFactor = 1); std::vector GetSelectedBoostedJets(const std::vector&, const float, const float, const float, const float, const jetID::jetID); bool passesMuonPOGIdTight(const pat::Muon&); bool isGoodMuon(const pat::Muon&, const float, const float, const muonID::muonID, const coneSize::coneSize, const corrType::corrType); @@ -222,7 +224,8 @@ class MiniAODHelper{ TLorentzVector * top =0 , TLorentzVector * antitop =0 ); - bool jetdPtMatched(const pat::Jet& inputJet); + bool GenJet_Match( const pat::Jet&, const edm::Handle&, reco::GenJet&, const double& ); + bool jetdPtMatched(const pat::Jet& inputJet, const reco::GenJet& genjet); double getJERfactor( const int, const double, const double, const double ); std::vector CorrectMET(const std::vector& oldJetsForMET, const std::vector& newJetsForMET, const std::vector& pfMETs); // Return weight factor dependent on number of true PU interactions diff --git a/MiniAODHelper/plugins/MiniAODAnalyzer.cc b/MiniAODHelper/plugins/MiniAODAnalyzer.cc index 1128fee..a18a5e4 100644 --- a/MiniAODHelper/plugins/MiniAODAnalyzer.cc +++ b/MiniAODHelper/plugins/MiniAODAnalyzer.cc @@ -97,6 +97,7 @@ class MiniAODAnalyzer : public edm::EDAnalyzer { edm::EDGetTokenT ak4jetToken; edm::EDGetTokenT > ak8jetToken; edm::EDGetTokenT beamspotToken; + edm::EDGetTokenT token_genjets; edm::EDGetTokenT ca12jetToken; edm::EDGetTokenT ca12filtjetToken; @@ -226,7 +227,8 @@ MiniAODAnalyzer::MiniAODAnalyzer(const edm::ParameterSet& iConfig) heptopfatjetToken = consumes (edm::InputTag(std::string("selectedPatJetsHEPTopFatPF"))); heptopsubjetToken = consumes (edm::InputTag(std::string("selectedPatJetsHEPTopSubPF"))); - + + token_genjets = consumes (edm::InputTag(std::string("slimmedGenJets"))); numEvents_ = 0; @@ -446,7 +448,8 @@ MiniAODAnalyzer::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup edm::Handle heptopsubpfjets; iEvent.getByToken(heptopsubjetToken,heptopsubpfjets); - + edm::Handle genjets; + iEvent.getByToken(token_genjets, genjets); edm::Handle > pfjetAK8Handle; iEvent.getByToken(ak8jetToken,pfjetAK8Handle); @@ -758,7 +761,7 @@ MiniAODAnalyzer::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup std::vector rawJets = miniAODhelper.GetUncorrectedJets(*pfjets); std::vector jetsNoMu = miniAODhelper.RemoveOverlaps(selectedMuons, rawJets); std::vector jetsNoEle = miniAODhelper.RemoveOverlaps(selectedElectrons, jetsNoMu); - std::vector correctedJets = miniAODhelper.GetCorrectedJets(jetsNoEle, iEvent, iSetup); + std::vector correctedJets = miniAODhelper.GetCorrectedJets(jetsNoEle, iEvent, iSetup, genjets); std::vector cleanSelectedJets = miniAODhelper.GetSelectedJets(correctedJets, 30., 2.4, jetID::jetLoose, '-' ); int nJet = 0; diff --git a/MiniAODHelper/src/LeptonSFHelper.cpp b/MiniAODHelper/src/LeptonSFHelper.cpp index bb03fdd..ae4507d 100644 --- a/MiniAODHelper/src/LeptonSFHelper.cpp +++ b/MiniAODHelper/src/LeptonSFHelper.cpp @@ -167,6 +167,11 @@ float LeptonSFHelper::GetElectronSF( float electronPt , float electronEta , int upval=nomval+error; downval=nomval-error; + + if(nomval==0) { + thisBin = h_ele_ID_abseta_pt_ratio->FindBin( searchEta , 20.001 ); + nomval=h_ele_ID_abseta_pt_ratio->GetBinContent( thisBin ); + } } else if ( type == "Trigger" ){ @@ -185,6 +190,11 @@ float LeptonSFHelper::GetElectronSF( float electronPt , float electronEta , int error=h_ele_ISO_abseta_pt_ratio->GetBinError( thisBin ); upval=nomval+error; downval=nomval-error; + + if(nomval==0) { + thisBin = h_ele_ISO_abseta_pt_ratio->FindBin( searchEta , 20.001 ); + nomval=h_ele_ISO_abseta_pt_ratio->GetBinContent( thisBin ); + } } else { @@ -225,13 +235,18 @@ float LeptonSFHelper::GetMuonSF( float muonPt , float muonEta , int syst , std: upval=upval*( 1.0+0.01 ); downval=downval*( 1.0-0.01 ); - + if(nomval==0) { + thisBin = h_mu_ID_abseta_pt_ratio->FindBin( searchEta , 20.001 ); + nomval=h_mu_ID_abseta_pt_ratio->GetBinContent( thisBin ); + } } else if ( type == "Trigger" ){ float mult4p2 = 0.2834; float mult4p3 = 0.7166; + //float mult4p2 = 0.0482; + //float mult4p3 = 0.9518; thisBin = h_mu_TRIGGER_abseta_pt_ratio4p3->FindBin(searchEta,searchPt); float nomval4p3=h_mu_TRIGGER_abseta_pt_ratio4p3->GetBinContent(thisBin); @@ -260,7 +275,11 @@ float LeptonSFHelper::GetMuonSF( float muonPt , float muonEta , int syst , std: downval=( nomval-error ); upval=upval*( 1.0+0.005 ); downval=downval*( 1.0-0.005 ); - + + if(nomval==0) { + thisBin = h_mu_ISO_abseta_pt_ratio->FindBin( searchEta , 20.001 ); + nomval=h_mu_ISO_abseta_pt_ratio->GetBinContent( thisBin ); + } } else { @@ -345,9 +364,10 @@ void LeptonSFHelper::SetElectronHistos( ){ void LeptonSFHelper::SetMuonHistos( ){ - std::string IDinputFile = std::string(getenv("CMSSW_BASE")) + "/src/MiniAOD/MiniAODHelper/data/leptonSF/" + "MuonID_Z_2016runB_2p6fb.root"; + std::string IDinputFile = std::string(getenv("CMSSW_BASE")) + "/src/MiniAOD/MiniAODHelper/data/leptonSF/" + "MuonID_Z_RunBCD_prompt80X_7p65.root"; std::string TRIGGERinputFile = std::string(getenv("CMSSW_BASE")) + "/src/MiniAOD/MiniAODHelper/data/leptonSF/" + "SingleMuonTrigger_Z_RunCD_Reco76X_Feb15.root"; - std::string ISOinputFile = std::string(getenv("CMSSW_BASE")) + "/src/MiniAOD/MiniAODHelper/data/leptonSF/" + "MuonISO_Z_2016runB_2p6fb.root"; + //std::string TRIGGERinputFile = std::string(getenv("CMSSW_BASE")) + "/src/MiniAOD/MiniAODHelper/data/leptonSF/" + "SingleMuonTrigger_Z_RunBCD_prompt80X_7p65.root"; + std::string ISOinputFile = std::string(getenv("CMSSW_BASE")) + "/src/MiniAOD/MiniAODHelper/data/leptonSF/" + "MuonIso_Z_RunBCD_prompt80X_7p65.root"; TFile *f_IDSF = new TFile(std::string(IDinputFile).c_str(),"READ"); TFile *f_TRIGGERSF = new TFile(std::string(TRIGGERinputFile).c_str(),"READ"); @@ -356,6 +376,8 @@ void LeptonSFHelper::SetMuonHistos( ){ h_mu_ID_abseta_pt_ratio = (TH2F*)f_IDSF->Get("MC_NUM_TightIDandIPCut_DEN_genTracks_PAR_pt_spliteta_bin1/abseta_pt_ratio"); h_mu_TRIGGER_abseta_pt_ratio4p3 = (TH2F*)f_TRIGGERSF->Get("runD_IsoMu20_OR_IsoTkMu20_HLTv4p3_PtEtaBins/abseta_pt_ratio"); h_mu_TRIGGER_abseta_pt_ratio4p2 = (TH2F*)f_TRIGGERSF->Get("runD_IsoMu20_OR_IsoTkMu20_HLTv4p2_PtEtaBins/abseta_pt_ratio"); + //h_mu_TRIGGER_abseta_pt_ratio4p3 = (TH2F*)f_TRIGGERSF->Get("IsoMu22_OR_IsoTkMu22_PtEtaBins_Run273158_to_274093/efficienciesDATA/abseta_pt_DATA"); + //h_mu_TRIGGER_abseta_pt_ratio4p2 = (TH2F*)f_TRIGGERSF->Get("IsoMu22_OR_IsoTkMu22_PtEtaBins_Run274094_to_276097/efficienciesDATA/abseta_pt_DATA"); h_mu_ISO_abseta_pt_ratio = (TH2F*)f_ISOSF->Get("MC_NUM_TightRelIso_DEN_TightID_PAR_pt_spliteta_bin1/abseta_pt_ratio"); } diff --git a/MiniAODHelper/src/MiniAODHelper.cc b/MiniAODHelper/src/MiniAODHelper.cc index 807e476..b3d9faa 100644 --- a/MiniAODHelper/src/MiniAODHelper.cc +++ b/MiniAODHelper/src/MiniAODHelper.cc @@ -362,7 +362,7 @@ std::vector MiniAODHelper::GetUncorrectedJets( pat::Jet -MiniAODHelper::GetCorrectedJet(const pat::Jet& inputJet, const edm::Event& event, const edm::EventSetup& setup, const sysType::sysType iSysType, const bool& doJES, const bool& doJER, const float& corrFactor, const float& uncFactor){ +MiniAODHelper::GetCorrectedJet(const pat::Jet& inputJet, const edm::Event& event, const edm::EventSetup& setup, const edm::Handle& genjets, const sysType::sysType iSysType, const bool& doJES, const bool& doJER, const float& corrFactor, const float& uncFactor){ if( !doJES && !doJER ) return inputJet; @@ -401,15 +401,17 @@ MiniAODHelper::GetCorrectedJet(const pat::Jet& inputJet, const edm::Event& event /// JER if( doJER){ double jerSF = 1.; - if( outputJet.genJet() && deltaR(outputJet,*outputJet.genJet())<0.4/2 && jetdPtMatched(outputJet)){ + reco::GenJet matched_genjet; + if ( GenJet_Match(outputJet, genjets, matched_genjet, 0.4) ) { + //if( outputJet.genJet() && deltaR(outputJet,*outputJet.genJet())<0.4/2 && jetdPtMatched(outputJet)){ if( iSysType == sysType::JERup ){ - jerSF = getJERfactor(uncFactor, fabs(outputJet.eta()), outputJet.genJet()->pt(), outputJet.pt()); + jerSF = getJERfactor(uncFactor, fabs(outputJet.eta()), matched_genjet.pt(), outputJet.pt()); } else if( iSysType == sysType::JERdown ){ - jerSF = getJERfactor(-uncFactor, fabs(outputJet.eta()), outputJet.genJet()->pt(), outputJet.pt()); + jerSF = getJERfactor(-uncFactor, fabs(outputJet.eta()), matched_genjet.pt(), outputJet.pt()); } else { - jerSF = getJERfactor(0, fabs(outputJet.eta()), outputJet.genJet()->pt(), outputJet.pt()); + jerSF = getJERfactor(0, fabs(outputJet.eta()), matched_genjet.pt(), outputJet.pt()); } // std::cout << "----->checking gen Jet pt " << jet.genJet()->pt() << ", jerSF is" << jerSF << std::endl; } @@ -423,7 +425,7 @@ MiniAODHelper::GetCorrectedJet(const pat::Jet& inputJet, const edm::Event& event float -MiniAODHelper::GetJetCorrectionFactor(const pat::Jet& inputJet, const edm::Event& event, const edm::EventSetup& setup, const sysType::sysType iSysType, const bool& doJES, const bool& doJER, const float& corrFactor, const float& uncFactor){ +MiniAODHelper::GetJetCorrectionFactor(const pat::Jet& inputJet, const edm::Event& event, const edm::EventSetup& setup, const edm::Handle& genjets, const sysType::sysType iSysType, const bool& doJES, const bool& doJER, const float& corrFactor, const float& uncFactor){ double factor = 1.; @@ -467,15 +469,17 @@ MiniAODHelper::GetJetCorrectionFactor(const pat::Jet& inputJet, const edm::Event /// JER if( doJER){ double jerSF = 1.; - if( outputJet.genJet() && deltaR(outputJet,*outputJet.genJet())<0.4/2 && jetdPtMatched(outputJet)){ + reco::GenJet matched_genjet; + if ( GenJet_Match(outputJet, genjets, matched_genjet, 0.4) ) { + //if( outputJet.genJet() && deltaR(outputJet,*outputJet.genJet())<0.4/2 && jetdPtMatched(outputJet)){ if( iSysType == sysType::JERup ){ - jerSF = getJERfactor(uncFactor, fabs(outputJet.eta()), outputJet.genJet()->pt(), outputJet.pt()); + jerSF = getJERfactor(uncFactor, fabs(outputJet.eta()), matched_genjet.pt(), outputJet.pt()); } else if( iSysType == sysType::JERdown ){ - jerSF = getJERfactor(-uncFactor, fabs(outputJet.eta()), outputJet.genJet()->pt(), outputJet.pt()); + jerSF = getJERfactor(-uncFactor, fabs(outputJet.eta()), matched_genjet.pt(), outputJet.pt()); } else { - jerSF = getJERfactor(0, fabs(outputJet.eta()), outputJet.genJet()->pt(), outputJet.pt()); + jerSF = getJERfactor(0, fabs(outputJet.eta()), matched_genjet.pt(), outputJet.pt()); } // std::cout << "----->checking gen Jet pt " << jet.genJet()->pt() << ", jerSF is" << jerSF << std::endl; } @@ -490,7 +494,7 @@ MiniAODHelper::GetJetCorrectionFactor(const pat::Jet& inputJet, const edm::Event pat::Jet -MiniAODHelper::GetCorrectedAK8Jet(const pat::Jet& inputJet, const edm::Event& event, const edm::EventSetup& setup, const sysType::sysType iSysType, const bool& doJES, const bool& doJER, const float& corrFactor, const float& uncFactor){ +MiniAODHelper::GetCorrectedAK8Jet(const pat::Jet& inputJet, const edm::Event& event, const edm::EventSetup& setup, const edm::Handle& genjets, const sysType::sysType iSysType, const bool& doJES, const bool& doJER, const float& corrFactor, const float& uncFactor){ if( !doJES && !doJER ) return inputJet; @@ -529,15 +533,17 @@ MiniAODHelper::GetCorrectedAK8Jet(const pat::Jet& inputJet, const edm::Event& ev /// JER if( doJER){ double jerSF = 1.; - if( outputJet.genJet() && deltaR(outputJet,*outputJet.genJet())<0.4/2 && jetdPtMatched(outputJet)){ + reco::GenJet matched_genjet; + if ( GenJet_Match(outputJet, genjets, matched_genjet, 0.8) ) { + //if( outputJet.genJet() && deltaR(outputJet,*outputJet.genJet())<0.4/2 && jetdPtMatched(outputJet)){ if( iSysType == sysType::JERup ){ - jerSF = getJERfactor(uncFactor, fabs(outputJet.eta()), outputJet.genJet()->pt(), outputJet.pt()); + jerSF = getJERfactor(uncFactor, fabs(outputJet.eta()), matched_genjet.pt(), outputJet.pt()); } else if( iSysType == sysType::JERdown ){ - jerSF = getJERfactor(-uncFactor, fabs(outputJet.eta()), outputJet.genJet()->pt(), outputJet.pt()); + jerSF = getJERfactor(-uncFactor, fabs(outputJet.eta()), matched_genjet.pt(), outputJet.pt()); } else { - jerSF = getJERfactor(0, fabs(outputJet.eta()), outputJet.genJet()->pt(), outputJet.pt()); + jerSF = getJERfactor(0, fabs(outputJet.eta()), matched_genjet.pt(), outputJet.pt()); } // std::cout << "----->checking gen Jet pt " << jet.genJet()->pt() << ", jerSF is" << jerSF << std::endl; } @@ -551,7 +557,7 @@ MiniAODHelper::GetCorrectedAK8Jet(const pat::Jet& inputJet, const edm::Event& ev float -MiniAODHelper::GetAK8JetCorrectionFactor(const pat::Jet& inputJet, const edm::Event& event, const edm::EventSetup& setup, const sysType::sysType iSysType, const bool& doJES, const bool& doJER, const float& corrFactor, const float& uncFactor){ +MiniAODHelper::GetAK8JetCorrectionFactor(const pat::Jet& inputJet, const edm::Event& event, const edm::EventSetup& setup, const edm::Handle& genjets, const sysType::sysType iSysType, const bool& doJES, const bool& doJER, const float& corrFactor, const float& uncFactor){ double factor = 1.; @@ -594,15 +600,17 @@ MiniAODHelper::GetAK8JetCorrectionFactor(const pat::Jet& inputJet, const edm::Ev /// JER if( doJER){ double jerSF = 1.; - if( outputJet.genJet() && deltaR(outputJet,*outputJet.genJet())<0.4/2 && jetdPtMatched(outputJet)){ + reco::GenJet matched_genjet; + if ( GenJet_Match(outputJet, genjets, matched_genjet, 0.8) ) { + //if( outputJet.genJet() && deltaR(outputJet,*outputJet.genJet())<0.4/2 && jetdPtMatched(outputJet)){ if( iSysType == sysType::JERup ){ - jerSF = getJERfactor(uncFactor, fabs(outputJet.eta()), outputJet.genJet()->pt(), outputJet.pt()); + jerSF = getJERfactor(uncFactor, fabs(outputJet.eta()), matched_genjet.pt(), outputJet.pt()); } else if( iSysType == sysType::JERdown ){ - jerSF = getJERfactor(-uncFactor, fabs(outputJet.eta()), outputJet.genJet()->pt(), outputJet.pt()); + jerSF = getJERfactor(-uncFactor, fabs(outputJet.eta()), matched_genjet.pt(), outputJet.pt()); } else { - jerSF = getJERfactor(0, fabs(outputJet.eta()), outputJet.genJet()->pt(), outputJet.pt()); + jerSF = getJERfactor(0, fabs(outputJet.eta()), matched_genjet.pt(), outputJet.pt()); } // std::cout << "----->checking gen Jet pt " << jet.genJet()->pt() << ", jerSF is" << jerSF << std::endl; } @@ -617,7 +625,7 @@ MiniAODHelper::GetAK8JetCorrectionFactor(const pat::Jet& inputJet, const edm::Ev std::vector -MiniAODHelper::GetCorrectedJets(const std::vector& inputJets, const edm::Event& event, const edm::EventSetup& setup, const sysType::sysType iSysType, const bool& doJES, const bool& doJER, const float& corrFactor, const float& uncFactor){ +MiniAODHelper::GetCorrectedJets(const std::vector& inputJets, const edm::Event& event, const edm::EventSetup& setup, const edm::Handle& genjets, const sysType::sysType iSysType, const bool& doJES, const bool& doJER, const float& corrFactor, const float& uncFactor){ if( !doJES && !doJER ) return inputJets; @@ -626,7 +634,7 @@ MiniAODHelper::GetCorrectedJets(const std::vector& inputJets, const ed std::vector outputJets; for( std::vector::const_iterator it = inputJets.begin(), ed = inputJets.end(); it != ed; ++it ){ - outputJets.push_back(GetCorrectedJet(*it,event,setup,iSysType,doJES,doJER,corrFactor,uncFactor)); + outputJets.push_back(GetCorrectedJet(*it,event,setup, genjets, iSysType,doJES,doJER,corrFactor,uncFactor)); } return outputJets; @@ -684,7 +692,7 @@ MiniAODHelper::GetCorrectedJets(const std::vector& inputJets, const sy std::vector -MiniAODHelper::GetCorrectedBoostedJets(const std::vector& inputBoostedJets, const edm::Event& event, const edm::EventSetup& setup, const sysType::sysType iSysType, const bool& doJES, const bool& doJER, const float& corrFactor, const float& uncFactor){ +MiniAODHelper::GetCorrectedBoostedJets(const std::vector& inputBoostedJets, const edm::Event& event, const edm::EventSetup& setup, const edm::Handle& genjets, const sysType::sysType iSysType, const bool& doJES, const bool& doJER, const float& corrFactor, const float& uncFactor){ if( !doJES && !doJER ) return inputBoostedJets; @@ -695,19 +703,19 @@ MiniAODHelper::GetCorrectedBoostedJets(const std::vector& i for( std::vector::const_iterator it = inputBoostedJets.begin(), ed = inputBoostedJets.end(); it != ed; ++it ){ boosted::BoostedJet outputBoostedJet = *it; - outputBoostedJet.fatjet = GetCorrectedAK8Jet(outputBoostedJet.fatjet,event,setup,iSysType,doJES,false,corrFactor,uncFactor); + outputBoostedJet.fatjet = GetCorrectedAK8Jet(outputBoostedJet.fatjet,event,setup,genjets,iSysType,doJES,false,corrFactor,uncFactor); - outputBoostedJet.nonW = GetCorrectedJet(outputBoostedJet.nonW,event,setup,iSysType,doJES,doJER,corrFactor,uncFactor); - outputBoostedJet.W1 = GetCorrectedJet(outputBoostedJet.W1 ,event,setup,iSysType,doJES,doJER,corrFactor,uncFactor); - outputBoostedJet.W2 = GetCorrectedJet(outputBoostedJet.W2 ,event,setup,iSysType,doJES,doJER,corrFactor,uncFactor); + outputBoostedJet.nonW = GetCorrectedJet(outputBoostedJet.nonW,event,setup,genjets,iSysType,doJES,doJER,corrFactor,uncFactor); + outputBoostedJet.W1 = GetCorrectedJet(outputBoostedJet.W1 ,event,setup,genjets,iSysType,doJES,doJER,corrFactor,uncFactor); + outputBoostedJet.W2 = GetCorrectedJet(outputBoostedJet.W2 ,event,setup,genjets,iSysType,doJES,doJER,corrFactor,uncFactor); - outputBoostedJet.subjets = GetCorrectedJets(outputBoostedJet.subjets ,event,setup,iSysType,doJES,doJER,corrFactor,uncFactor); - outputBoostedJet.filterjets = GetCorrectedJets(outputBoostedJet.filterjets,event,setup,iSysType,doJES,doJER,corrFactor,uncFactor); + outputBoostedJet.subjets = GetCorrectedJets(outputBoostedJet.subjets ,event,setup,genjets,iSysType,doJES,doJER,corrFactor,uncFactor); + outputBoostedJet.filterjets = GetCorrectedJets(outputBoostedJet.filterjets,event,setup,genjets,iSysType,doJES,doJER,corrFactor,uncFactor); outputBoostedJet.topjet.setP4(outputBoostedJet.nonW.p4()+outputBoostedJet.W1.p4()+outputBoostedJet.W2.p4()); // Correction of pruned mass - outputBoostedJet.prunedMass = outputBoostedJet.prunedMass * GetAK8JetCorrectionFactor(outputBoostedJet.fatjet,event,setup,iSysType,doJES,false,corrFactor,uncFactor); + outputBoostedJet.prunedMass = outputBoostedJet.prunedMass * GetAK8JetCorrectionFactor(outputBoostedJet.fatjet,event,setup,genjets,iSysType,doJES,false,corrFactor,uncFactor); // Recalculation of fRec double _mtmass = 172.3; @@ -2408,23 +2416,52 @@ std::vector MiniAODHelper::GetDeltaRCleanedJets( return outputJets; } - /// JER function -bool MiniAODHelper::jetdPtMatched(const pat::Jet& inputJet) { +bool MiniAODHelper::GenJet_Match(const pat::Jet& inputJet, const edm::Handle& genjets, reco::GenJet& matched_genjet, const double& Rcone) { + + double dpt_min=99999; + double dpt; + double dR; + int genjet_match = 0; + + // checking which genjets have dR < 0.2 and dpT < 3*sigma_mc for this particular jet + // if multiple genjets found satisfying this, then select the one with dpT minimum + + for (reco::GenJetCollection::const_iterator iter=genjets->begin();iter!=genjets->end();++iter){ + dpt = fabs(inputJet.pt()-iter->pt()); + dR = DeltaR( &inputJet , iter ); + if (dR<(Rcone/2)) { + if ( jetdPtMatched(inputJet,*iter) ) { + genjet_match = 1; + if (dpt <= dpt_min) { + matched_genjet = *(iter); + dpt_min = dpt; + } + } + } + } + + if (genjet_match==1) + return true; + else + return false; +} + +bool MiniAODHelper::jetdPtMatched(const pat::Jet& inputJet, const reco::GenJet& genjet) { const double jet_eta=inputJet.eta(); for( unsigned int i = 0 ; i < JER_etaMax.size() ; i ++){ - if(jet_eta < JER_etaMax[i] && jet_eta > JER_etaMin[i] && useRho < JER_rhoMax[i] && useRho > JER_rhoMin[i] ) { + if(jet_eta < JER_etaMax[i] && jet_eta >= JER_etaMin[i] && useRho < JER_rhoMax[i] && useRho >= JER_rhoMin[i] ) { double jet_pt=inputJet.pt(); if(jet_pt< JER_PtMin[i]){jet_pt=JER_PtMin[i];} if(jet_pt> JER_PtMax[i]){jet_pt=JER_PtMax[i];} const double sigma_mc=sqrt( JER_Par0[i]*fabs(JER_Par0[i]) / (jet_pt*jet_pt)+JER_Par1[i]*JER_Par1[i]*pow(jet_pt,JER_Par3[i])+JER_Par2[i]*JER_Par2[i]); - const double genjet_pt=inputJet.genJet()->pt(); + const double genjet_pt=genjet.pt(); if(fabs(jet_pt-genjet_pt)<3*sigma_mc*jet_pt) { //cout << "#############################################################" << endl; //cout << "JER dPt condition is satisfied" << endl; diff --git a/SkimDilep/plugins/SkimDilep.cc b/SkimDilep/plugins/SkimDilep.cc index 8268a2e..6b1d688 100644 --- a/SkimDilep/plugins/SkimDilep.cc +++ b/SkimDilep/plugins/SkimDilep.cc @@ -97,7 +97,7 @@ class SkimDilep : public edm::EDFilter { edm::EDGetTokenT electronToken; edm::EDGetTokenT muonToken; edm::EDGetTokenT ak4jetToken; - + edm::EDGetTokenT token_genjets; // HLTConfigProvider hlt_config_; // HLTConfigProvider filter_config_; @@ -133,6 +133,7 @@ SkimDilep::SkimDilep(const edm::ParameterSet& iConfig) electronToken = consumes (edm::InputTag(std::string("slimmedElectrons"))); muonToken = consumes (edm::InputTag(std::string("slimmedMuons"))); ak4jetToken = consumes (edm::InputTag(std::string("slimmedJets"))); + token_genjets = consumes (edm::InputTag(std::string("slimmedGenJets"))); // //Check to make sure we're not doing anything in consistent @@ -193,6 +194,9 @@ SkimDilep::filter(edm::Event& iEvent, const edm::EventSetup& iSetup) edm::Handle pfjets; iEvent.getByToken(ak4jetToken,pfjets); + edm::Handle genjets; + iEvent.getByToken(token_genjets, genjets); + // edm::Handle triggerResults; // iEvent.getByToken(triggerResultsToken, triggerResults); @@ -279,7 +283,7 @@ SkimDilep::filter(edm::Event& iEvent, const edm::EventSetup& iSetup) std::vector rawJets = miniAODhelper.GetUncorrectedJets(*pfjets); std::vector jetsNoMu = miniAODhelper.RemoveOverlaps(selectedMuonsLoose, rawJets); std::vector jetsNoEle = miniAODhelper.RemoveOverlaps(selectedElectronsLoose, jetsNoMu); - std::vector correctedJets = miniAODhelper.GetCorrectedJets(jetsNoEle, iEvent, iSetup); + std::vector correctedJets = miniAODhelper.GetCorrectedJets(jetsNoEle, iEvent, iSetup, genjets); std::vector cleanSelectedJets = miniAODhelper.GetSelectedJets(correctedJets, 15., 2.4, jetID::jetLoose, '-' ); // pt set to 15GeV // at least two jets