CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
List of all members | Public Member Functions | Static Public Member Functions | Private Attributes
CMSDAS11DijetAnalyzer Class Reference

#include <CMSDAS11DijetAnalyzer.h>

Inheritance diagram for CMSDAS11DijetAnalyzer:
edm::EDAnalyzer edm::EDConsumerBase

Public Member Functions

void analyze (const edm::Event &, const edm::EventSetup &)
 
virtual void beginJob ()
 
 CMSDAS11DijetAnalyzer (const edm::ParameterSet &)
 
virtual void endJob (void)
 
virtual ~CMSDAS11DijetAnalyzer ()
 
- Public Member Functions inherited from edm::EDAnalyzer
void callWhenNewProductsRegistered (std::function< void(BranchDescription const &)> const &func)
 
 EDAnalyzer ()
 
ModuleDescription const & moduleDescription () const
 
std::string workerType () const
 
virtual ~EDAnalyzer ()
 
- Public Member Functions inherited from edm::EDConsumerBase
std::vector< ConsumesInfoconsumesInfo () const
 
 EDConsumerBase ()
 
ProductHolderIndexAndSkipBit indexFrom (EDGetToken, BranchType, TypeID const &) const
 
void itemsMayGet (BranchType, std::vector< ProductHolderIndexAndSkipBit > &) const
 
void itemsToGet (BranchType, std::vector< ProductHolderIndexAndSkipBit > &) const
 
std::vector
< ProductHolderIndexAndSkipBit >
const & 
itemsToGetFromEvent () const
 
void labelsForToken (EDGetToken iToken, Labels &oLabels) const
 
void modulesDependentUpon (const std::string &iProcessName, std::vector< const char * > &oModuleLabels) const
 
void modulesWhoseProductsAreConsumed (std::vector< ModuleDescription const * > &modules, ProductRegistry const &preg, std::map< std::string, ModuleDescription const * > const &labelsToDesc, std::string const &processName) const
 
bool registeredToConsume (ProductHolderIndex, bool, BranchType) const
 
bool registeredToConsumeMany (TypeID const &, BranchType) const
 
void updateLookup (BranchType iBranchType, ProductHolderIndexHelper const &)
 
virtual ~EDConsumerBase ()
 

Static Public Member Functions

static bool compare_JetPt (const reco::CaloJet &jet1, const reco::CaloJet &jet2)
 
- Static Public Member Functions inherited from edm::EDAnalyzer
static const std::string & baseType ()
 
static void fillDescriptions (ConfigurationDescriptions &descriptions)
 
static void prevalidate (ConfigurationDescriptions &)
 

Private Attributes

TH1D * hCorDijetMass
 
TH1D * hCorDijetXsec
 
TH1D * hDijetDeltaEta
 
TH1D * hDijetDeltaPhi
 
TH1D * hInnerDijetMass
 
TH1D * hJet1EMF
 
TH1D * hJet1Eta
 
TH1D * hJet1Phi
 
TH1D * hJet1Pt
 
TH1D * hJet2EMF
 
TH1D * hJet2Eta
 
TH1D * hJet2Phi
 
TH1D * hJet2Pt
 
TH1D * hJetCorrPt
 
TH1D * hJetEMF
 
TH1D * hJetEta
 
TH1D * hJetPhi
 
TH1D * hJetRawPt
 
TH1D * hOuterDijetMass
 
TH1D * hRawDijetMass
 
TH1D * hVertexZ
 
double innerDeltaEta
 
double JESbias
 
std::string jetCorrections
 
edm::InputTag jetSrc
 
double outerDeltaEta
 
edm::InputTag vertexSrc
 

Additional Inherited Members

- Public Types inherited from edm::EDAnalyzer
typedef EDAnalyzer ModuleType
 
- Protected Member Functions inherited from edm::EDConsumerBase
template<typename ProductType , BranchType B = InEvent>
EDGetTokenT< ProductType > consumes (edm::InputTag const &tag)
 
EDGetToken consumes (const TypeToGet &id, edm::InputTag const &tag)
 
template<BranchType B>
EDGetToken consumes (TypeToGet const &id, edm::InputTag const &tag)
 
ConsumesCollector consumesCollector ()
 Use a ConsumesCollector to gather consumes information from helper functions. More...
 
template<typename ProductType , BranchType B = InEvent>
void consumesMany ()
 
void consumesMany (const TypeToGet &id)
 
template<BranchType B>
void consumesMany (const TypeToGet &id)
 
template<typename ProductType , BranchType B = InEvent>
EDGetTokenT< ProductType > mayConsume (edm::InputTag const &tag)
 
EDGetToken mayConsume (const TypeToGet &id, edm::InputTag const &tag)
 
template<BranchType B>
EDGetToken mayConsume (const TypeToGet &id, edm::InputTag const &tag)
 

Detailed Description

Definition at line 20 of file CMSDAS11DijetAnalyzer.h.

Constructor & Destructor Documentation

CMSDAS11DijetAnalyzer::CMSDAS11DijetAnalyzer ( const edm::ParameterSet params)

Definition at line 20 of file CMSDAS11DijetAnalyzer.cc.

References hCorDijetMass, hDijetDeltaEta, hDijetDeltaPhi, hInnerDijetMass, hJet1EMF, hJet1Eta, hJet1Phi, hJet1Pt, hJet2EMF, hJet2Eta, hJet2Phi, hJet2Pt, hJetCorrPt, hJetEMF, hJetEta, hJetPhi, hJetRawPt, hOuterDijetMass, hVertexZ, TFileService::make(), and NBINS.

20  :
22  jetSrc(params.getParameter<edm::InputTag>("jetSrc")),
23  vertexSrc(params.getParameter<edm::InputTag>("vertexSrc")),
24  jetCorrections(params.getParameter<std::string>("jetCorrections")),
25  innerDeltaEta(params.getParameter<double>("innerDeltaEta")),
26  outerDeltaEta(params.getParameter<double>("outerDeltaEta")),
27  JESbias(params.getParameter<double>("JESbias"))
28 {
29  // setup file service
31 
32  const int NBINS=36;
33  Double_t BOUNDARIES[NBINS] = { 220, 244, 270, 296, 325, 354, 386, 419, 453,
34  489, 526, 565, 606, 649, 693, 740, 788, 838,
35  890, 944, 1000, 1058, 1118, 1181, 1246, 1313, 1383,
36  1455, 1530, 1607, 1687, 1770, 1856, 1945, 2037, 2132 };
37 
38  // setup histograms
39  hVertexZ = fs->make<TH1D>("hVertexZ", "Z position of the Vertex",50,-20,20);
40  hJetRawPt = fs->make<TH1D>("hJetRawPt","Raw Jet Pt",50,0,1000);
41  hJetCorrPt = fs->make<TH1D>("hJetCorrPt","Corrected Jet Pt",50,0,1000);
42  hJet1Pt = fs->make<TH1D>("hJet1Pt","Corrected Jet1 Pt",50,0,1000);
43  hJet2Pt = fs->make<TH1D>("hJet2Pt","Corrected Jet2 Pt",50,0,1000);
44 
45  hJetEta = fs->make<TH1D>("hJetEta","Corrected Jet Eta", 10,-5,5);
46  hJet1Eta = fs->make<TH1D>("hJet1Eta","Corrected Jet1 Eta",10,-5,5);
47  hJet2Eta = fs->make<TH1D>("hJet2Eta","Corrected Jet2 Eta",10,-5,5);
48 
49  hJetPhi = fs->make<TH1D>("hJetPhi","Corrected Jet Phi", 10,-3.1415,3.1415);
50  hJet1Phi = fs->make<TH1D>("hJet1Phi","Corrected Jet1 Phi",10,-3.1415,3.1415);
51  hJet2Phi = fs->make<TH1D>("hJet2Phi","Corrected Jet2 Phi",10,-3.1415,3.1415);
52 
53  hJetEMF = fs->make<TH1D>("hJetEMF","EM Fraction of Jets",50,0,1);
54  hJet1EMF = fs->make<TH1D>("hJet1EMF","EM Fraction of Jet1",50,0,1);
55  hJet2EMF = fs->make<TH1D>("hJet2EMF","EM Fraction of Jet2",50,0,1);
56 
57  hCorDijetMass = fs->make<TH1D>("hCorDijetMass","Corrected Dijet Mass",NBINS-1,BOUNDARIES);
58  hDijetDeltaPhi= fs->make<TH1D>("hDijetDeltaPhi","Dijet |#Delta #phi|",50,0,3.1415);
59  hDijetDeltaEta= fs->make<TH1D>("hDijetDeltaEta","Dijet |#Delta #eta|",50,0,6);
60 
61  hInnerDijetMass = fs->make<TH1D>("hInnerDijetMass","Corrected Inner Dijet Mass",NBINS-1,BOUNDARIES);
62  hOuterDijetMass = fs->make<TH1D>("hOuterDijetMass","Corrected Outer Dijet Mass",NBINS-1,BOUNDARIES);
63 }
T getParameter(std::string const &) const
T * make(const Args &...args) const
make new ROOT object
Definition: TFileService.h:64
const int NBINS
virtual CMSDAS11DijetAnalyzer::~CMSDAS11DijetAnalyzer ( )
inlinevirtual

Definition at line 24 of file CMSDAS11DijetAnalyzer.h.

24 {}

Member Function Documentation

void CMSDAS11DijetAnalyzer::analyze ( const edm::Event iEvent,
const edm::EventSetup iSetup 
)
virtual

Implements edm::EDAnalyzer.

Definition at line 69 of file CMSDAS11DijetAnalyzer.cc.

References compare_JetPt(), JetCorrector::correction(), mvaPFMET_cff::corrector, gather_cfg::cout, HLTFastRecoForTau_cff::deltaEta, reco::CaloJet::emEnergyFraction(), reco::LeafCandidate::eta(), eta(), edm::Event::getByLabel(), JetCorrector::getJetCorrector(), hCorDijetMass, hDijetDeltaEta, hDijetDeltaPhi, hInnerDijetMass, hJet1EMF, hJet1Eta, hJet1Phi, hJet1Pt, hJet2EMF, hJet2Eta, hJet2Phi, hJet2Pt, hJetCorrPt, hJetEMF, hJetEta, hJetPhi, hJetRawPt, hOuterDijetMass, hVertexZ, innerDeltaEta, edm::HandleBase::isValid(), JESbias, metsig::jet, jetCorrections, jetSrc, reco::Vertex::ndof(), outerDeltaEta, reco::LeafCandidate::p4(), phi, reco::LeafCandidate::phi(), reco::Vertex::position(), EnergyCorrector::pt, reco::LeafCandidate::pt(), pileupReCalc_HLTpaths::scale, reco::Jet::scaleEnergy(), patRefSel_refMuJets_cfi::selectedJets, python.multivaluedict::sort(), vertexSrc, and reco::Vertex::z().

70 {
71 
73  double mWeight;
75  iEvent.getByLabel("generator", hEventInfo);
76  if (hEventInfo.isValid())
77  mWeight = hEventInfo->weight();
78  else
79  mWeight = 1.;
80 
82  // Get event ID information
84  //int nrun=iEvent.id().run();
85  //int nlumi=iEvent.luminosityBlock();
86  //int nevent=iEvent.id().event();
87 
89  // Get Primary Vertex Information
91 
92  // magic to get the vertices from EDM
94  iEvent.getByLabel(vertexSrc, vertices_h);
95  if (!vertices_h.isValid()) {
96  std::cout<<"Didja hear the one about the empty vertex collection?\n";
97  return;
98  }
99 
100  // require in the event that there is at least one reconstructed vertex
101  if(vertices_h->size()<=0) return;
102 
103  // pick the first (i.e. highest sum pt) verte
104  const reco::Vertex* theVertex=&(vertices_h->front());
105 
106  // require that the vertex meets certain criteria
107  if(theVertex->ndof()<5) return;
108  if(fabs(theVertex->z())>24.0) return;
109  if(fabs(theVertex->position().rho())>2.0) return;
110 
112  // Get Jet Information
114 
115  // magic to get the jets from EDM
117  iEvent.getByLabel(jetSrc, jets_h);
118  if (!jets_h.isValid()) {
119  std::cout<<"Didja hear the one about the empty jet collection?\n";
120  return;
121  }
122 
123  // magic to get the jet energy corrections
125 
126  // collection of selected jets
127  std::vector<reco::CaloJet> selectedJets;
128 
129  // loop over the jet collection
130  for(reco::CaloJetCollection::const_iterator j_it = jets_h->begin(); j_it!=jets_h->end(); j_it++) {
131  reco::CaloJet jet = *j_it;
132 
133  // calculate and apply the correction
134  double scale = corrector->correction(jet.p4());
135 
136  // Introduce a purposeful bias to the correction, to show what happens
137  scale *= JESbias;
138 
139  // fill the histograms
140  hJetRawPt->Fill(jet.pt(),mWeight);
141  hJetEta->Fill(jet.eta(),mWeight);
142  hJetPhi->Fill(jet.phi(),mWeight);
143  hJetEMF->Fill(jet.emEnergyFraction(),mWeight);
144 
145  // correct the jet energy
146  jet.scaleEnergy(scale);
147 
148  // now fill the correct jet pt after the correction
149  hJetCorrPt->Fill(jet.pt(),mWeight);
150 
151  // put the selected jets into a collection
152  selectedJets.push_back(jet);
153  }
154 
155  // require at least two jets to continue
156  if(selectedJets.size()<2) return;
157 
158  //sort by corrected pt (not the same order as raw pt, sometimes)
159  sort(selectedJets.begin(), selectedJets.end(), compare_JetPt);
160 
161  // select high pt, central, non-noise-like jets
162  if (selectedJets[0].pt()<50.0) return;
163  if (fabs(selectedJets[0].eta())>2.5) return;
164  if (selectedJets[0].emEnergyFraction()<0.01) return;
165  if (selectedJets[1].pt()<50.0) return;
166  if (fabs(selectedJets[1].eta())>2.5) return;
167  if (selectedJets[1].emEnergyFraction()<0.01) return;
168 
169 
170  // fill histograms for the jets in our dijets, only
171  hJet1Pt ->Fill(selectedJets[0].pt(),mWeight);
172  hJet1Eta->Fill(selectedJets[0].eta(),mWeight);
173  hJet1Phi->Fill(selectedJets[0].phi(),mWeight);
174  hJet1EMF->Fill(selectedJets[0].emEnergyFraction(),mWeight);
175  hJet2Pt ->Fill(selectedJets[1].pt(),mWeight);
176  hJet2Eta->Fill(selectedJets[1].eta(),mWeight);
177  hJet2Phi->Fill(selectedJets[1].phi(),mWeight);
178  hJet2EMF->Fill(selectedJets[1].emEnergyFraction(),mWeight);
179 
180  //Get the mass of the two leading jets. Needs their 4-vectors...
181  double corMass = (selectedJets[0].p4()+selectedJets[1].p4()).M();
182  double deltaEta = fabs(selectedJets[0].eta()-selectedJets[1].eta());
183  if (corMass < 489) return;
184  if (deltaEta > 1.3) return;
185  hCorDijetMass->Fill(corMass,mWeight);
186  hDijetDeltaPhi->Fill(fabs(selectedJets[0].phi()-selectedJets[1].phi()) ,mWeight);
187  hDijetDeltaEta->Fill(deltaEta ,mWeight);
188 
189  //Fill the inner and outer dijet mass spectra, to make the ratio from
190  if (deltaEta < innerDeltaEta) hInnerDijetMass->Fill(corMass,mWeight);
191  else if (deltaEta < outerDeltaEta) hOuterDijetMass->Fill(corMass,mWeight);
192 
193  hVertexZ->Fill(theVertex->z(),mWeight);
194  return;
195 }
Jets made from CaloTowers.
Definition: CaloJet.h:29
virtual void scaleEnergy(double fScale)
scale energy of the jet
static bool compare_JetPt(const reco::CaloJet &jet1, const reco::CaloJet &jet2)
virtual double correction(const LorentzVector &fJet) const =0
get correction using Jet information only
T eta() const
const Point & position() const
position
Definition: Vertex.h:106
virtual double eta() const
momentum pseudorapidity
virtual double pt() const
transverse momentum
tuple corrector
Definition: mvaPFMET_cff.py:89
double z() const
y coordinate
Definition: Vertex.h:112
bool isValid() const
Definition: HandleBase.h:75
bool getByLabel(InputTag const &tag, Handle< PROD > &result) const
Definition: Event.h:405
double ndof() const
Definition: Vertex.h:102
static const JetCorrector * getJetCorrector(const std::string &fName, const edm::EventSetup &fSetup)
retrieve corrector from the event setup. troughs exception if something is missing ...
Definition: JetCorrector.cc:50
tuple cout
Definition: gather_cfg.py:121
virtual double phi() const
momentum azimuthal angle
virtual const LorentzVector & p4() const
four-momentum Lorentz vector
Definition: LeafCandidate.h:99
float emEnergyFraction() const
Definition: CaloJet.h:102
Definition: DDAxes.h:10
virtual void CMSDAS11DijetAnalyzer::beginJob ( void  )
inlinevirtual

Reimplemented from edm::EDAnalyzer.

Definition at line 25 of file CMSDAS11DijetAnalyzer.h.

25 {}
static bool CMSDAS11DijetAnalyzer::compare_JetPt ( const reco::CaloJet jet1,
const reco::CaloJet jet2 
)
inlinestatic

Definition at line 28 of file CMSDAS11DijetAnalyzer.h.

References reco::LeafCandidate::pt().

Referenced by analyze().

28  {
29  return (jet1.pt() > jet2.pt() );
30  }
virtual double pt() const
transverse momentum
void CMSDAS11DijetAnalyzer::endJob ( void  )
virtual

Reimplemented from edm::EDAnalyzer.

Definition at line 65 of file CMSDAS11DijetAnalyzer.cc.

65  {
66 }

Member Data Documentation

TH1D* CMSDAS11DijetAnalyzer::hCorDijetMass
private

Definition at line 50 of file CMSDAS11DijetAnalyzer.h.

Referenced by analyze(), and CMSDAS11DijetAnalyzer().

TH1D* CMSDAS11DijetAnalyzer::hCorDijetXsec
private

Definition at line 51 of file CMSDAS11DijetAnalyzer.h.

TH1D* CMSDAS11DijetAnalyzer::hDijetDeltaEta
private

Definition at line 61 of file CMSDAS11DijetAnalyzer.h.

Referenced by analyze(), and CMSDAS11DijetAnalyzer().

TH1D* CMSDAS11DijetAnalyzer::hDijetDeltaPhi
private

Definition at line 60 of file CMSDAS11DijetAnalyzer.h.

Referenced by analyze(), and CMSDAS11DijetAnalyzer().

TH1D* CMSDAS11DijetAnalyzer::hInnerDijetMass
private

Definition at line 63 of file CMSDAS11DijetAnalyzer.h.

Referenced by analyze(), and CMSDAS11DijetAnalyzer().

TH1D* CMSDAS11DijetAnalyzer::hJet1EMF
private

Definition at line 55 of file CMSDAS11DijetAnalyzer.h.

Referenced by analyze(), and CMSDAS11DijetAnalyzer().

TH1D* CMSDAS11DijetAnalyzer::hJet1Eta
private

Definition at line 53 of file CMSDAS11DijetAnalyzer.h.

Referenced by analyze(), and CMSDAS11DijetAnalyzer().

TH1D* CMSDAS11DijetAnalyzer::hJet1Phi
private

Definition at line 54 of file CMSDAS11DijetAnalyzer.h.

Referenced by analyze(), and CMSDAS11DijetAnalyzer().

TH1D* CMSDAS11DijetAnalyzer::hJet1Pt
private

Definition at line 52 of file CMSDAS11DijetAnalyzer.h.

Referenced by analyze(), and CMSDAS11DijetAnalyzer().

TH1D* CMSDAS11DijetAnalyzer::hJet2EMF
private

Definition at line 59 of file CMSDAS11DijetAnalyzer.h.

Referenced by analyze(), and CMSDAS11DijetAnalyzer().

TH1D* CMSDAS11DijetAnalyzer::hJet2Eta
private

Definition at line 57 of file CMSDAS11DijetAnalyzer.h.

Referenced by analyze(), and CMSDAS11DijetAnalyzer().

TH1D* CMSDAS11DijetAnalyzer::hJet2Phi
private

Definition at line 58 of file CMSDAS11DijetAnalyzer.h.

Referenced by analyze(), and CMSDAS11DijetAnalyzer().

TH1D* CMSDAS11DijetAnalyzer::hJet2Pt
private

Definition at line 56 of file CMSDAS11DijetAnalyzer.h.

Referenced by analyze(), and CMSDAS11DijetAnalyzer().

TH1D* CMSDAS11DijetAnalyzer::hJetCorrPt
private

Definition at line 43 of file CMSDAS11DijetAnalyzer.h.

Referenced by analyze(), and CMSDAS11DijetAnalyzer().

TH1D* CMSDAS11DijetAnalyzer::hJetEMF
private

Definition at line 47 of file CMSDAS11DijetAnalyzer.h.

Referenced by analyze(), and CMSDAS11DijetAnalyzer().

TH1D* CMSDAS11DijetAnalyzer::hJetEta
private

Definition at line 45 of file CMSDAS11DijetAnalyzer.h.

Referenced by analyze(), and CMSDAS11DijetAnalyzer().

TH1D* CMSDAS11DijetAnalyzer::hJetPhi
private

Definition at line 46 of file CMSDAS11DijetAnalyzer.h.

Referenced by analyze(), and CMSDAS11DijetAnalyzer().

TH1D* CMSDAS11DijetAnalyzer::hJetRawPt
private

Definition at line 44 of file CMSDAS11DijetAnalyzer.h.

Referenced by analyze(), and CMSDAS11DijetAnalyzer().

TH1D* CMSDAS11DijetAnalyzer::hOuterDijetMass
private

Definition at line 64 of file CMSDAS11DijetAnalyzer.h.

Referenced by analyze(), and CMSDAS11DijetAnalyzer().

TH1D* CMSDAS11DijetAnalyzer::hRawDijetMass
private

Definition at line 49 of file CMSDAS11DijetAnalyzer.h.

TH1D* CMSDAS11DijetAnalyzer::hVertexZ
private

Definition at line 42 of file CMSDAS11DijetAnalyzer.h.

Referenced by analyze(), and CMSDAS11DijetAnalyzer().

double CMSDAS11DijetAnalyzer::innerDeltaEta
private

Definition at line 37 of file CMSDAS11DijetAnalyzer.h.

Referenced by analyze().

double CMSDAS11DijetAnalyzer::JESbias
private

Definition at line 39 of file CMSDAS11DijetAnalyzer.h.

Referenced by analyze().

std::string CMSDAS11DijetAnalyzer::jetCorrections
private

Definition at line 36 of file CMSDAS11DijetAnalyzer.h.

Referenced by analyze().

edm::InputTag CMSDAS11DijetAnalyzer::jetSrc
private

Definition at line 34 of file CMSDAS11DijetAnalyzer.h.

Referenced by analyze().

double CMSDAS11DijetAnalyzer::outerDeltaEta
private

Definition at line 38 of file CMSDAS11DijetAnalyzer.h.

Referenced by analyze().

edm::InputTag CMSDAS11DijetAnalyzer::vertexSrc
private

Definition at line 35 of file CMSDAS11DijetAnalyzer.h.

Referenced by analyze().