Go to the documentation of this file.00001 #include "DataFormats/Common/interface/Handle.h"
00002 #include "DataFormats/JetReco/interface/Jet.h"
00003 #include "DataFormats/PatCandidates/interface/Jet.h"
00004 #include "PhysicsTools/PatExamples/interface/AnalysisTasksAnalyzerJEC.h"
00005 #include "TH2.h"
00006 #include "TProfile.h"
00008 AnalysisTasksAnalyzerJEC::AnalysisTasksAnalyzerJEC(const edm::ParameterSet& cfg, TFileDirectory& fs):
00009 edm::BasicAnalyzer::BasicAnalyzer(cfg, fs),
00010 Jets_(cfg.getParameter<edm::InputTag>("Jets")),
00011 jecLevel_(cfg.getParameter<std::string>("jecLevel")),
00012 patJetCorrFactors_(cfg.getParameter<std::string>("patJetCorrFactors")),
00013 help_(cfg.getParameter<bool>("help"))
00014 {
00015
00016 hists_["Response" ] = fs.make<TH2F>("Response" , "response; #eta;p_{T}(reco)/p_{T}(gen)" , 5, -3.3, 3.3, 100, 0.4, 1.6);
00017 jetInEvents_=0;
00018 }
00020 AnalysisTasksAnalyzerJEC::~AnalysisTasksAnalyzerJEC()
00021 {
00022
00024 TFile* file=new TFile("myResponse"+TString(jecLevel_)+".root", "RECREATE");
00025 TProfile* prof = hists_["Response" ]->ProfileX();
00026 prof->Write();
00027 file->Write();
00028 file->Close();
00029 }
00031 void
00032 AnalysisTasksAnalyzerJEC::analyze(const edm::EventBase& event)
00033 {
00034
00035
00036 using pat::Jet;
00037
00038
00039 edm::Handle<std::vector<Jet> > Jets;
00040 event.getByLabel(Jets_, Jets);
00041
00042
00043 for(std::vector<Jet>::const_iterator jet_it=Jets->begin(); jet_it!=Jets->end(); ++jet_it){
00044
00045
00047 if(help_==true){
00048 if(jetInEvents_==0){
00049 std::cout<<"\n \n number of available JEC sets: "<< jet_it->availableJECSets().size() << std::endl;
00050 for(unsigned int k=0; k< jet_it->availableJECSets().size(); ++k){
00051 std::cout<<"\n \n available JEC Set: "<< jet_it->availableJECSets()[k] << std::endl;
00052 }
00053 std::cout<<"\n \n*** You found out which JEC Sets exist! Now correct it in your config file."<<std::endl;
00054 std::cout<<"\n \n number of available JEC levels "<< jet_it->availableJECLevels().size() << std::endl;
00055 for(unsigned int k=0; k< jet_it->availableJECLevels().size(); ++k){
00056 std::cout<<"\n \n available JEC: "<< jet_it->availableJECLevels(patJetCorrFactors_)[k] << std::endl;
00057 }
00058 std::cout<<"\n \n**** You did it correctly congratulations!!!! And you found out which JEC levels are saved within the jets. Correct this in your configuration file." <<std::endl;
00059 }
00060 }
00061
00062 if(jet_it->genParticlesSize()>0){
00063 hists_["Response" ]->Fill( jet_it->correctedJet(jecLevel_).eta(), jet_it->correctedJet(jecLevel_).pt()/ jet_it->genParticle(0)->pt());
00064 std::cout<<"\n \n**** You did it correctly congratulations!!!! "<< std::endl;
00065 }
00066 jetInEvents_+=1;
00067 }
00068 }