JetAnalyzers/src/CMSDAS11DijetTestAnalyzer.cc

0001 // CMSDAS11DijetTestAnalyzer.cc
0002 // Description: A basic dijet analyzer for the CMSDAS 2011
0003 // Author: John Paul Chou
0004 // Date: January 12, 2011
0005
0006 #include "RecoJets/JetAnalyzers/interface/CMSDAS11DijetTestAnalyzer.h"
0007
0008 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0009 #include "FWCore/Framework/interface/MakerMacros.h"
0010 #include "FWCore/ServiceRegistry/interface/Service.h"
0011 #include "CommonTools/UtilAlgos/interface/TFileService.h"
0012
0013 #include "DataFormats/VertexReco/interface/Vertex.h"
0014 #include "JetMETCorrections/Objects/interface/JetCorrector.h"
0015
0016 #include "SimDataFormats/GeneratorProducts/interface/GenEventInfoProduct.h"
0017 #include "SimDataFormats/GeneratorProducts/interface/GenRunInfoProduct.h"
0018
0019 #include <TH1D.h>
0020 #include <TH2D.h>
0021
0022 CMSDAS11DijetTestAnalyzer::CMSDAS11DijetTestAnalyzer(edm::ParameterSet const& params)
0023     : jetSrc(params.getParameter<edm::InputTag>("jetSrc")),
0024       vertexSrc(params.getParameter<edm::InputTag>("vertexSrc")),
0025       jetCorrections(params.getParameter<std::string>("jetCorrections")),
0026       innerDeltaEta(params.getParameter<double>("innerDeltaEta")),
0027       outerDeltaEta(params.getParameter<double>("outerDeltaEta")),
0028       JESbias(params.getParameter<double>("JESbias")) {
0029   // setup file service
0030   usesResource(TFileService::kSharedResource);
0031   edm::Service<TFileService> fs;
0032
0033   const int NBINS = 36;
0034   Double_t BOUNDARIES[NBINS] = {220,  244,  270,  296,  325,  354,  386,  419,  453,  489,  526,  565,
0035                                 606,  649,  693,  740,  788,  838,  890,  944,  1000, 1058, 1118, 1181,
0036                                 1246, 1313, 1383, 1455, 1530, 1607, 1687, 1770, 1856, 1945, 2037, 2132};
0037
0038   // setup histograms
0039   hVertexZ = fs->make<TH1D>("hVertexZ", "Z position of the Vertex", 50, -20, 20);
0040   hJetRawPt = fs->make<TH1D>("hJetRawPt", "Raw Jet Pt", 50, 0, 1000);
0041   hJetCorrPt = fs->make<TH1D>("hJetCorrPt", "Corrected Jet Pt", 50, 0, 1000);
0042   hJet1Pt = fs->make<TH1D>("hJet1Pt", "Corrected Jet1 Pt", 50, 0, 1000);
0043   hJet2Pt = fs->make<TH1D>("hJet2Pt", "Corrected Jet2 Pt", 50, 0, 1000);
0044
0045   hJetEta = fs->make<TH1D>("hJetEta", "Corrected Jet Eta", 50, -5, 5);
0046   hJet1Eta = fs->make<TH1D>("hJet1Eta", "Corrected Jet1 Eta", 50, -5, 5);
0047   hJet2Eta = fs->make<TH1D>("hJet2Eta", "Corrected Jet2 Eta", 50, -5, 5);
0048
0049   hJetPhi = fs->make<TH1D>("hJetPhi", "Corrected Jet Phi", 50, -3.1415, 3.1415);
0050   hJet1Phi = fs->make<TH1D>("hJet1Phi", "Corrected Jet1 Phi", 50, -3.1415, 3.1415);
0051   hJet2Phi = fs->make<TH1D>("hJet2Phi", "Corrected Jet2 Phi", 50, -3.1415, 3.1415);
0052
0053   hJetEMF = fs->make<TH1D>("hJetEMF", "EM Fraction of Jets", 50, 0, 1);
0054   hJet1EMF = fs->make<TH1D>("hJet1EMF", "EM Fraction of Jet1", 50, 0, 1);
0055   hJet2EMF = fs->make<TH1D>("hJet2EMF", "EM Fraction of Jet2", 50, 0, 1);
0056
0057   hCorDijetMass = fs->make<TH1D>("hCorDijetMass", "Corrected Dijet Mass", NBINS - 1, BOUNDARIES);
0058   hDijetDeltaPhi = fs->make<TH1D>("hDijetDeltaPhi", "Dijet |#Delta #phi|", 50, 0, 3.1415);
0059   hDijetDeltaEta = fs->make<TH1D>("hDijetDeltaEta", "Dijet |#Delta #eta|", 50, 0, 1.3);
0060   hDijetDeltaPhiNJets =
0061       fs->make<TH2D>("hDijetDeltaPhiNJets", "Dijet |#Delta #phi| v. the number of jets", 50, 0, 3.1415, 7, 0.5, 7.5);
0062   hDijetEta1Eta2 = fs->make<TH2D>("hDijetEta1Eta2", "Eta 1 versus Eta 2 of dijet events", 50, -5, 5, 50, -5, 5);
0063
0064   hInnerDijetMass = fs->make<TH1D>("hInnerDijetMass", "Corrected Inner Dijet Mass", NBINS - 1, BOUNDARIES);
0065   hOuterDijetMass = fs->make<TH1D>("hOuterDijetMass", "Corrected Outer Dijet Mass", NBINS - 1, BOUNDARIES);
0066 }
0067
0068 void CMSDAS11DijetTestAnalyzer::endJob(void) {}
0069
0070 void CMSDAS11DijetTestAnalyzer::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0071   ////////Get Weight, if this is MC//////////////
0072   double mWeight;
0073   edm::Handle<GenEventInfoProduct> hEventInfo;
0074   iEvent.getByLabel("generator", hEventInfo);
0075   if (hEventInfo.isValid()) {
0076     mWeight = hEventInfo->weight();
0077   } else
0078     mWeight = 1.;
0079
0080   ////////////////////////////////////////////
0081   // Get event ID information
0082   ////////////////////////////////////////////
0083   //int nrun=iEvent.id().run();
0084   //int nlumi=iEvent.luminosityBlock();
0085   //int nevent=iEvent.id().event();
0086
0087   ////////////////////////////////////////////
0088   // Get Primary Vertex Information
0089   ////////////////////////////////////////////
0090
0091   // magic to get the vertices from EDM
0092   edm::Handle<std::vector<reco::Vertex> > vertices_h;
0093   iEvent.getByLabel(vertexSrc, vertices_h);
0094   if (!vertices_h.isValid()) {
0095     std::cout << "Didja hear the one about the empty vertex collection?\n";
0096     return;
0097   }
0098
0099   // require in the event that there is at least one reconstructed vertex
0100   if (vertices_h->empty())
0101     return;
0102
0103   // pick the first (i.e. highest sum pt) verte
0104   const reco::Vertex* theVertex = &(vertices_h->front());
0105
0106   // require that the vertex meets certain criteria
0107   if (theVertex->ndof() < 5)
0108     return;
0109   if (fabs(theVertex->z()) > 24.0)
0110     return;
0111   if (fabs(theVertex->position().rho()) > 2.0)
0112     return;
0113
0114   ////////////////////////////////////////////
0115   // Get Jet Information
0116   ////////////////////////////////////////////
0117
0118   // magic to get the jets from EDM
0119   edm::Handle<reco::CaloJetCollection> jets_h;
0120   iEvent.getByLabel(jetSrc, jets_h);
0121   if (!jets_h.isValid()) {
0122     std::cout << "Didja hear the one about the empty jet collection?\n";
0123     return;
0124   }
0125
0126   // magic to get the jet energy corrections
0127   const JetCorrector* corrector = JetCorrector::getJetCorrector(jetCorrections, iSetup);
0128   corrector->vectorialCorrection();
0129   // collection of selected jets
0130   std::vector<reco::CaloJet> selectedJets;
0131
0132   // loop over the jet collection
0133   for (reco::CaloJetCollection::const_iterator j_it = jets_h->begin(); j_it != jets_h->end(); j_it++) {
0134     reco::CaloJet jet = *j_it;
0135
0136     // put the selected jets into a collection
0137     selectedJets.push_back(jet);
0138   }
0139
0140   hVertexZ->Fill(theVertex->z(), mWeight);
0141   return;
0142 }
0143
0144 DEFINE_FWK_MODULE(CMSDAS11DijetTestAnalyzer);