CMS 3D CMS Logo

ZHistogrammer.cc
Go to the documentation of this file.
6 #include "TH1.h"
7 
8 class ZHistogrammer : public edm::EDAnalyzer {
9 public:
11 
12 private:
13  void analyze(const edm::Event &event, const edm::EventSetup &setup) override;
20  TH1F *h_invmMuMu_;
23  //TH1F *h_mZ2vs3MC_, *h_ptZ2vs3MC_, *h_phiZ2vs3MC_, *h_thetaZ2vs3MC_, *h_etaZ2vs3MC_, *h_rapidityZ2vs3MC_;
27 };
28 
36 #include <cmath>
37 #include <iostream>
38 
39 using namespace std;
40 using namespace reco;
41 using namespace edm;
42 
44  : zToken_(consumes<CandidateCollection>(pset.getParameter<InputTag>("z"))),
45  genToken_(consumes<CandidateCollection>(pset.getParameter<InputTag>("gen"))),
46  matchToken_(consumes<CandMatchMap>(pset.getParameter<InputTag>("match"))),
47  nbinsMass_(pset.getUntrackedParameter<unsigned int>("nbinsMass")),
48  nbinsPt_(pset.getUntrackedParameter<unsigned int>("nbinsPt")),
49  nbinsAng_(pset.getUntrackedParameter<unsigned int>("nbinsAng")),
50  nbinsMassRes_(pset.getUntrackedParameter<unsigned int>("nbinsMassRes")),
51  massMax_(pset.getUntrackedParameter<double>("massMax")),
52  ptMax_(pset.getUntrackedParameter<double>("ptMax")),
53  angMax_(pset.getUntrackedParameter<double>("angMax")),
54  massResMax_(pset.getUntrackedParameter<double>("massResMax")) {
55  cout << ">>> Z Histogrammer constructor" << endl;
57  TFileDirectory ZHisto = fs->mkdir("ZRecoHisto");
58  TFileDirectory ZMCHisto = fs->mkdir("ZMCHisto");
59  TFileDirectory ZResHisto = fs->mkdir("ZResolutionHisto");
60  //TFileDirectory Z2vs3MCHisto = fs->mkdir( "Z2vs3MCHisto" );
61  h_nZ_ = ZHisto.make<TH1F>("ZNumber", "number of Z particles", 11, -0.5, 10.5);
62  h_mZ_ = ZHisto.make<TH1F>("ZMass", "Z mass (GeV/c^{2})", nbinsMass_, 0, massMax_);
63  h_ptZ_ = ZHisto.make<TH1F>("ZPt", "Z p_{t} (GeV/c)", nbinsPt_, 0, ptMax_);
64  h_phiZ_ = ZHisto.make<TH1F>("ZPhi", "Z #phi", nbinsAng_, -angMax_, angMax_);
65  h_thetaZ_ = ZHisto.make<TH1F>("Ztheta", "Z #theta", nbinsAng_, 0, angMax_);
66  h_etaZ_ = ZHisto.make<TH1F>("ZEta", "Z #eta", nbinsAng_, -angMax_, angMax_);
67  h_rapidityZ_ = ZHisto.make<TH1F>("ZRapidity", "Z rapidity", nbinsAng_, -angMax_, angMax_);
68  h_invmMuMu_ = ZHisto.make<TH1F>("MuMuMass", "#mu #mu invariant mass", nbinsMass_, 0, massMax_);
69  h_nZMC_ = ZMCHisto.make<TH1F>("ZMCNumber", "number of Z MC particles", 11, -0.5, 10.5);
70  h_mZMC_ = ZMCHisto.make<TH1F>("ZMCMass", "Z MC mass (GeV/c^{2})", nbinsMass_, 0, massMax_);
71  h_ptZMC_ = ZMCHisto.make<TH1F>("ZMCPt", "Z MC p_{t} (GeV/c)", nbinsPt_, 0, ptMax_);
72  h_phiZMC_ = ZMCHisto.make<TH1F>("ZMCPhi", "Z MC #phi", nbinsAng_, -angMax_, angMax_);
73  h_thetaZMC_ = ZMCHisto.make<TH1F>("ZMCTheta", "Z MC #theta", nbinsAng_, 0, angMax_);
74  h_etaZMC_ = ZMCHisto.make<TH1F>("ZMCEta", "Z MC #eta", nbinsAng_, -angMax_, angMax_);
75  h_rapidityZMC_ = ZMCHisto.make<TH1F>("ZMCRapidity", "Z MC rapidity", nbinsAng_, -angMax_, angMax_);
76  h_invmMuMuMC_ = ZMCHisto.make<TH1F>("MuMuMCMass", "#mu #mu MC invariant mass", nbinsMass_, 0, massMax_);
77  /*
78  h_mZ2vs3MC_ = Z2vs3MCHisto.make<TH1F>("Z2vs3MCMass", "Z MC st 2 vs st 3 mass (GeV/c^{2})",
79  nbinsMassRes_, -massResMax_, massResMax_);
80  h_ptZ2vs3MC_ = Z2vs3MCHisto.make<TH1F>("Z2vs3MCPt", "Z MC st 2 vs st 3 p_{t} (GeV/c)",
81  nbinsPt_, -ptMax_, ptMax_);
82  h_phiZ2vs3MC_ = Z2vs3MCHisto.make<TH1F>("Z2vs3MCPhi", "Z MC st 2 vs st 3 #phi",
83  nbinsAng_, -angMax_, angMax_);
84  h_thetaZ2vs3MC_ = Z2vs3MCHisto.make<TH1F>("Z2vs3MCTheta", "Z MC st 2 vs st 3 #theta",
85  nbinsAng_, -angMax_, angMax_);
86  h_etaZ2vs3MC_ = Z2vs3MCHisto.make<TH1F>("Z2vs3MCEta", "Z MC st 2 vs st 3 #eta",
87  nbinsAng_, -angMax_, angMax_);
88  h_rapidityZ2vs3MC_ = Z2vs3MCHisto.make<TH1F>("Z2vs3MCRapidity", "Z MC st 2 vs st 3 rapidity",
89  nbinsAng_, -angMax_, angMax_);
90  */
91  h_mResZ_ = ZResHisto.make<TH1F>(
92  "ZMassResolution", "Z mass Resolution (GeV/c^{2})", nbinsMassRes_, -massResMax_, massResMax_);
93  h_ptResZ_ = ZResHisto.make<TH1F>("ZPtResolution", "Z p_{t} Resolution (GeV/c)", nbinsPt_, -ptMax_, ptMax_);
94  h_phiResZ_ = ZResHisto.make<TH1F>("ZPhiResolution", "Z #phi Resolution", nbinsAng_, -angMax_, angMax_);
95  h_thetaResZ_ = ZResHisto.make<TH1F>("ZThetaResolution", "Z #theta Resolution", nbinsAng_, -angMax_, angMax_);
96  h_etaResZ_ = ZResHisto.make<TH1F>("ZEtaResolution", "Z #eta Resolution", nbinsAng_, -angMax_, angMax_);
97  h_rapidityResZ_ = ZResHisto.make<TH1F>("ZRapidityResolution", "Z rapidity Resolution", nbinsAng_, -angMax_, angMax_);
98  h_mResZMuMu_ = ZResHisto.make<TH1F>("ZToMuMuRecoMassResolution",
99  "Z Reco vs matched final state #mu #mu mass Difference (GeV/c^{2})",
101  -massResMax_,
102  massResMax_);
104  ZResHisto.make<TH1F>("ZToMuMuRecoMassRatio", "Z Reco vs matched final state #mu #mu mass Ratio", 4000, 0, 2);
105  h_mResZMuMuMC_ = ZResHisto.make<TH1F>("ZToMuMuMCMassResolution",
106  "Z vs final state #mu #mu MC mass Difference (GeV/c^{2})",
107  nbinsMassRes_ / 2 + 1,
108  -2 * massResMax_ / nbinsMassRes_,
109  massResMax_);
111  ZResHisto.make<TH1F>("ZToMuMuMCMassRatio", "Z vs final state #mu #mu MC mass Ratio", 2002, 0.999, 2);
112 }
113 
115  cout << ">>> Z Histogrammer analyze" << endl;
119  event.getByToken(zToken_, z);
120  event.getByToken(genToken_, gen);
121  event.getByToken(matchToken_, match);
122  h_nZ_->Fill(z->size());
123  for (unsigned int i = 0; i < z->size(); ++i) {
124  const Candidate &zCand = (*z)[i];
125  h_mZ_->Fill(zCand.mass());
126  h_ptZ_->Fill(zCand.pt());
127  h_phiZ_->Fill(zCand.phi());
128  h_thetaZ_->Fill(zCand.theta());
129  h_etaZ_->Fill(zCand.eta());
130  h_rapidityZ_->Fill(zCand.rapidity());
131  CandidateRef zCandRef(z, i);
132  CandidateRef zMCMatch = (*match)[zCandRef];
133  if (zMCMatch.isNonnull() && zMCMatch->pdgId() == 23) {
134  h_mResZ_->Fill(zCandRef->mass() - zMCMatch->mass());
135  h_ptResZ_->Fill(zCandRef->pt() - zMCMatch->pt());
136  h_phiResZ_->Fill(zCandRef->phi() - zMCMatch->phi());
137  h_thetaResZ_->Fill(zCandRef->theta() - zMCMatch->theta());
138  h_etaResZ_->Fill(zCandRef->eta() - zMCMatch->eta());
139  h_rapidityResZ_->Fill(zCandRef->rapidity() - zMCMatch->rapidity());
140  const Candidate *dau0 = zMCMatch->daughter(0);
141  const Candidate *dau1 = zMCMatch->daughter(1);
142  for (unsigned int i0 = 0; i0 < dau0->numberOfDaughters(); ++i0) {
143  const Candidate *ddau0 = dau0->daughter(i0);
144  if (abs(ddau0->pdgId()) == 13 && ddau0->status() == 1) {
145  dau0 = ddau0;
146  break;
147  }
148  }
149  for (unsigned int i1 = 0; i1 < dau1->numberOfDaughters(); ++i1) {
150  const Candidate *ddau1 = dau1->daughter(i1);
151  if (abs(ddau1->pdgId()) == 13 && ddau1->status() == 1) {
152  dau1 = ddau1;
153  break;
154  }
155  }
156  assert(abs(dau0->pdgId()) == 13 && dau0->status() == 1);
157  assert(abs(dau1->pdgId()) == 13 && dau1->status() == 1);
158  double invMass = (dau0->p4() + dau1->p4()).mass();
159  h_invmMuMu_->Fill(invMass);
160  h_mResZMuMu_->Fill(zCand.mass() - invMass);
161  h_mRatioZMuMu_->Fill(zCand.mass() / invMass);
162  }
163  }
164  h_nZMC_->Fill(gen->size());
165  for (unsigned int i = 0; i < gen->size(); ++i) {
166  const Candidate &genCand = (*gen)[i];
167  if ((genCand.pdgId() == 23) && (genCand.status() == 2)) //this is an intermediate Z0
168  cout << ">>> intermediate Z0 found, with " << genCand.numberOfDaughters() << " daughters" << endl;
169  if ((genCand.pdgId() == 23) && (genCand.status() == 3)) { //this is a Z0
170  cout << ">>> Z0 found, with " << genCand.numberOfDaughters() << " daughters" << endl;
171  h_mZMC_->Fill(genCand.mass());
172  h_ptZMC_->Fill(genCand.pt());
173  h_phiZMC_->Fill(genCand.phi());
174  h_thetaZMC_->Fill(genCand.theta());
175  h_etaZMC_->Fill(genCand.eta());
176  h_rapidityZMC_->Fill(genCand.rapidity());
177  Particle::LorentzVector pZ(0, 0, 0, 0);
178  int nMu = 0;
179  for (unsigned int j = 0; j < genCand.numberOfDaughters(); ++j) {
180  const Candidate *dauGen = genCand.daughter(j);
181  /*
182  if((dauGen->pdgId() == 23) && (dauGen->status() == 2)) {
183  h_mZ2vs3MC_->Fill(genCand.mass() - dauGen->mass());
184  h_ptZ2vs3MC_->Fill(genCand.pt() - dauGen->pt());
185  h_phiZ2vs3MC_->Fill(genCand.phi() - dauGen->phi());
186  h_thetaZ2vs3MC_->Fill(genCand.theta() - dauGen->theta());
187  h_etaZ2vs3MC_->Fill(genCand.eta() - dauGen->eta());
188  h_rapidityZ2vs3MC_->Fill(genCand.rapidity() - dauGen->rapidity());
189  }
190  */
191  if ((abs(dauGen->pdgId()) == 13) && (dauGen->numberOfDaughters() != 0)) {
192  //we are looking for photons of final state radiation
193  cout << ">>> The muon " << j << " has " << dauGen->numberOfDaughters() << " daughters" << endl;
194  for (unsigned int k = 0; k < dauGen->numberOfDaughters(); ++k) {
195  const Candidate *dauMuGen = dauGen->daughter(k);
196  cout << ">>> Mu " << j << " daughter MC " << k << " PDG Id " << dauMuGen->pdgId() << ", status "
197  << dauMuGen->status() << ", charge " << dauMuGen->charge() << endl;
198  if (abs(dauMuGen->pdgId()) == 13 && dauMuGen->status() == 1) {
199  pZ += dauMuGen->p4();
200  nMu++;
201  }
202  }
203  }
204  }
205  assert(nMu == 2);
206  double mZ = pZ.mass();
207  h_invmMuMuMC_->Fill(mZ);
208  h_mResZMuMuMC_->Fill(genCand.mass() - mZ);
209  h_mRatioZMuMuMC_->Fill(genCand.mass() / mZ);
210  }
211  }
212 }
213 
215 
ZHistogrammer::zToken_
edm::EDGetTokenT< reco::CandidateCollection > zToken_
Definition: ZHistogrammer.cc:14
reco::Candidate::daughter
virtual const Candidate * daughter(size_type i) const =0
return daughter at a given position, i = 0, ... numberOfDaughters() - 1 (read only mode)
ZHistogrammer::h_ptZ_
TH1F * h_ptZ_
Definition: ZHistogrammer.cc:19
ZHistogrammer::h_rapidityResZ_
TH1F * h_rapidityResZ_
Definition: ZHistogrammer.cc:24
Handle.h
mps_fire.i
i
Definition: mps_fire.py:355
ZHistogrammer::h_phiZMC_
TH1F * h_phiZMC_
Definition: ZHistogrammer.cc:21
reco::Candidate::mass
virtual double mass() const =0
mass
reco::Candidate::eta
virtual double eta() const =0
momentum pseudorapidity
edm::EDGetTokenT
Definition: EDGetToken.h:33
ZHistogrammer::h_ptZMC_
TH1F * h_ptZMC_
Definition: ZHistogrammer.cc:21
TFileDirectory::make
T * make(const Args &... args) const
make new ROOT object
Definition: TFileDirectory.h:53
edm
HLT enums.
Definition: AlignableModifier.h:19
ZHistogrammer::matchToken_
edm::EDGetTokenT< reco::CandMatchMap > matchToken_
Definition: ZHistogrammer.cc:16
ZHistogrammer::h_rapidityZMC_
TH1F * h_rapidityZMC_
Definition: ZHistogrammer.cc:21
gather_cfg.cout
cout
Definition: gather_cfg.py:144
testProducerWithPsetDescEmpty_cfi.i1
i1
Definition: testProducerWithPsetDescEmpty_cfi.py:45
reco::Candidate::pt
virtual double pt() const =0
transverse momentum
cms::cuda::assert
assert(be >=bs)
ZHistogrammer::angMax_
double angMax_
Definition: ZHistogrammer.cc:18
ZHistogrammer::h_mRatioZMuMu_
TH1F * h_mRatioZMuMu_
Definition: ZHistogrammer.cc:25
reco::Candidate::status
virtual int status() const =0
status word
EDAnalyzer.h
TFileDirectory
Definition: TFileDirectory.h:24
reco
fixed size matrix
Definition: AlignmentAlgorithmBase.h:45
edm::Handle
Definition: AssociativeIterator.h:50
singleTopDQM_cfi.setup
setup
Definition: singleTopDQM_cfi.py:37
ZHistogrammer::ptMax_
double ptMax_
Definition: ZHistogrammer.cc:18
reco::Candidate::rapidity
virtual double rapidity() const =0
rapidity
reco::Candidate::theta
virtual double theta() const =0
momentum polar angle
CandMatchMap.h
edm::Ref
Definition: AssociativeIterator.h:58
edm::EDAnalyzer
Definition: EDAnalyzer.h:29
GenParticle.h
CandidateFwd.h
ZHistogrammer::massResMax_
double massResMax_
Definition: ZHistogrammer.cc:18
MakerMacros.h
ZHistogrammer::h_thetaZ_
TH1F * h_thetaZ_
Definition: ZHistogrammer.cc:19
DEFINE_FWK_MODULE
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:16
ZHistogrammer::h_mRatioZMuMuMC_
TH1F * h_mRatioZMuMuMC_
Definition: ZHistogrammer.cc:26
reco::Candidate::charge
virtual int charge() const =0
electric charge
Service.h
ZHistogrammer::h_mResZMuMuMC_
TH1F * h_mResZMuMuMC_
Definition: ZHistogrammer.cc:26
ZHistogrammer::nbinsPt_
unsigned int nbinsPt_
Definition: ZHistogrammer.cc:17
GenParticleFwd.h
ZHistogrammer::nbinsMass_
unsigned int nbinsMass_
Definition: ZHistogrammer.cc:17
DDAxes::z
ZHistogrammer::genToken_
edm::EDGetTokenT< reco::CandidateCollection > genToken_
Definition: ZHistogrammer.cc:15
ZHistogrammer::massMax_
double massMax_
Definition: ZHistogrammer.cc:18
ZHistogrammer::h_mZMC_
TH1F * h_mZMC_
Definition: ZHistogrammer.cc:21
reco::Candidate::numberOfDaughters
virtual size_type numberOfDaughters() const =0
number of daughters
dqmdumpme.k
k
Definition: dqmdumpme.py:60
ZHistogrammer::h_etaResZ_
TH1F * h_etaResZ_
Definition: ZHistogrammer.cc:24
gen
Definition: PythiaDecays.h:13
ZHistogrammer::h_nZMC_
TH1F * h_nZMC_
Definition: ZHistogrammer.cc:21
TFileService.h
ZHistogrammer::h_thetaResZ_
TH1F * h_thetaResZ_
Definition: ZHistogrammer.cc:24
TFileService::mkdir
TFileDirectory mkdir(const std::string &dir, const std::string &descr="")
create a new subdirectory
Definition: TFileService.h:69
edm::ParameterSet
Definition: ParameterSet.h:36
Event.h
edm::AssociationMap
Definition: AssociationMap.h:48
edm::Ref::isNonnull
bool isNonnull() const
Checks for non-null.
Definition: Ref.h:238
ZHistogrammer::analyze
void analyze(const edm::Event &event, const edm::EventSetup &setup) override
Definition: ZHistogrammer.cc:114
match
std::pair< typename Association::data_type::first_type, double > match(Reference key, Association association, bool bestMatchByMaxValue)
Generic matching function.
Definition: Utils.h:10
edm::Service< TFileService >
createfilelist.int
int
Definition: createfilelist.py:10
ZHistogrammer::h_rapidityZ_
TH1F * h_rapidityZ_
Definition: ZHistogrammer.cc:19
ZHistogrammer::h_nZ_
TH1F * h_nZ_
Definition: ZHistogrammer.cc:19
ZHistogrammer::nbinsAng_
unsigned int nbinsAng_
Definition: ZHistogrammer.cc:17
ZHistogrammer::ZHistogrammer
ZHistogrammer(const edm::ParameterSet &pset)
Definition: ZHistogrammer.cc:43
edm::EventSetup
Definition: EventSetup.h:57
ZHistogrammer::h_mResZMuMu_
TH1F * h_mResZMuMu_
Definition: ZHistogrammer.cc:25
reco::Candidate::pdgId
virtual int pdgId() const =0
PDG identifier.
ZHistogrammer::h_phiZ_
TH1F * h_phiZ_
Definition: ZHistogrammer.cc:19
InputTag.h
ZHistogrammer::h_mZ_
TH1F * h_mZ_
Definition: ZHistogrammer.cc:19
reco::Candidate
Definition: Candidate.h:27
reco::JetExtendedAssociation::LorentzVector
math::PtEtaPhiELorentzVectorF LorentzVector
Definition: JetExtendedAssociation.h:25
std
Definition: JetResolutionObject.h:76
ZHistogrammer::h_invmMuMuMC_
TH1F * h_invmMuMuMC_
Definition: ZHistogrammer.cc:22
EgHLTOffHistBins_cfi.mass
mass
Definition: EgHLTOffHistBins_cfi.py:34
ZHistogrammer::h_invmMuMu_
TH1F * h_invmMuMu_
Definition: ZHistogrammer.cc:20
relval_steps.gen
def gen(fragment, howMuch)
Production test section ####.
Definition: relval_steps.py:500
reco::Candidate::p4
virtual const LorentzVector & p4() const =0
four-momentum Lorentz vector
ZHistogrammer::h_mResZ_
TH1F * h_mResZ_
Definition: ZHistogrammer.cc:24
Candidate.h
ZHistogrammer::h_ptResZ_
TH1F * h_ptResZ_
Definition: ZHistogrammer.cc:24
funct::abs
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
ParameterSet.h
dqmiolumiharvest.j
j
Definition: dqmiolumiharvest.py:66
ZHistogrammer::nbinsMassRes_
unsigned int nbinsMassRes_
Definition: ZHistogrammer.cc:17
event
Definition: event.py:1
edm::Event
Definition: Event.h:73
ZHistogrammer::h_etaZMC_
TH1F * h_etaZMC_
Definition: ZHistogrammer.cc:21
reco::Candidate::phi
virtual double phi() const =0
momentum azimuthal angle
ZHistogrammer
Definition: ZHistogrammer.cc:8
edm::InputTag
Definition: InputTag.h:15
ZHistogrammer::h_phiResZ_
TH1F * h_phiResZ_
Definition: ZHistogrammer.cc:24
muonDTDigis_cfi.pset
pset
Definition: muonDTDigis_cfi.py:27
edm::OwnVector
Definition: OwnVector.h:24
ZHistogrammer::h_thetaZMC_
TH1F * h_thetaZMC_
Definition: ZHistogrammer.cc:21
ZHistogrammer::h_etaZ_
TH1F * h_etaZ_
Definition: ZHistogrammer.cc:19