CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
AcceptanceHistoProducer.cc
Go to the documentation of this file.
2 
4 
6 
8 
10  dqmDir_(cfg.getParameter<std::string>("dqmDir")),
11  srcGenParticles_(cfg.getParameter<edm::InputTag>("srcGenParticles")),
12  hPtPosPtNeg_(NULL), hEtaPosEtaNeg_(NULL), hPtPosEtaPos_(NULL), hPtNegEtaNeg_(NULL)
13 {
14 }
15 
17 {
18 }
19 
21 {
23 
24  dbe_->setCurrentFolder("MCEmbedding/ZmumuAcceptance/" + dqmDir_);
25  hPtPosPtNeg_ = dbe_->book2DD("PtPosPtNeg", "Positive muon transverse momentum vs. negative muon transverse momentum;p_{T}^{+};p_{T}^{-}", 500, 0.0, 500.0, 500, 0.0, 500.0);
26  hEtaPosEtaNeg_ = dbe_->book2DD("EtaPosEtaNeg", "Positive muon pseudorapdity vs. negative muon pseudorapidity;#eta^{+};#eta^{-}", 500, -2.5, 2.5, 500, -2.5, 2.5);
27  hPtPosEtaPos_ = dbe_->book2DD("PtPosEtaPos", "Positive muon transverse momentum vs. positive muon pseudorapidity;p_{T}^{+};#eta^{+}", 500, 0.0, 500.0, 500, -2.5, 2.5);
28  hPtNegEtaNeg_ = dbe_->book2DD("PtNegEtaNeg", "Negative muon transverse momentum vs. negative muon pseudorapidity;p_{T}^{-};#eta^{-}", 500, 0.0, 500.0, 500, -2.5, 2.5);
29 
30  hPtPosPtNeg_->setResetMe(true);
38 }
39 
41 {
42  assert(hPtPosPtNeg_ != NULL);
43  assert(hEtaPosEtaNeg_ != NULL);
44  assert(hPtPosEtaPos_ != NULL);
45  assert(hPtNegEtaNeg_ != NULL);
46 
48  evt.getByLabel(srcGenParticles_, genParticles);
49 
50  // Look for the two muons from the matrix element (status=3)
51  const reco::GenParticle* genPosMuonME = NULL;
52  const reco::GenParticle* genNegMuonME = NULL;
53  for(unsigned int i = 0; i < genParticles->size(); ++i)
54  {
55  const reco::GenParticle& part = (*genParticles)[i];
56  if(abs(part.pdgId()) == 13 && part.status() == 3)
57  {
58  if(part.charge() < 0)
59  if(!genNegMuonME)
60  genNegMuonME = &part;
61  if(part.charge() > 0)
62  if(!genPosMuonME)
63  genPosMuonME = &part;
64  }
65  }
66 
67  if(genPosMuonME && genNegMuonME)
68  {
69  const reco::GenParticle* genPosMuon = genPosMuonME;
70  const reco::GenParticle* genNegMuon = genNegMuonME;
71 
72  // Follow the decay chain to find the stable decay products (status=1)
73  while(genPosMuon->status() != 1 && genPosMuon->numberOfDaughters() > 0)
74  {
75  unsigned int i;
76  for(i = 0; i < genPosMuon->numberOfDaughters(); ++i)
77  {
78  const reco::GenParticle* daughter = dynamic_cast<const reco::GenParticle*>(genPosMuon->daughter(i));
79  if(abs(daughter->pdgId()) == 13)
80  { genPosMuon = daughter; break; }
81  }
82 
83  // No more muon daugthers? Maybe mu->e decay...
84  if(i == genPosMuon->numberOfDaughters())
85  break;
86  }
87 
88  while(genNegMuon->status() != 1 && genNegMuon->numberOfDaughters() > 0)
89  {
90  unsigned int i;
91  for(i = 0; i < genNegMuon->numberOfDaughters(); ++i)
92  {
93  const reco::GenParticle* daughter = dynamic_cast<const reco::GenParticle*>(genNegMuon->daughter(i));
94  if(abs(daughter->pdgId()) == 13)
95  { genNegMuon = daughter; break; }
96  }
97 
98  // No more muon daugthers? Maybe mu->e decay...
99  if(i == genNegMuon->numberOfDaughters())
100  break;
101  }
102 
103  // Loop might have broken earlier if there are no more daugthers, or if
104  // the decay chain does not end in a muon.
105  if(genPosMuon->status() == 1 && genNegMuon->status() == 1)
106  {
107  hPtPosPtNeg_->Fill(genPosMuon->pt(), genNegMuon->pt());
108  hEtaPosEtaNeg_->Fill(genPosMuon->eta(), genNegMuon->eta());
109  hPtPosEtaPos_->Fill(genPosMuon->pt(), genPosMuon->eta());
110  hPtNegEtaNeg_->Fill(genNegMuon->pt(), genNegMuon->eta());
111  }
112  }
113 }
114 
116 
int i
Definition: DBlmapReader.cc:9
virtual int pdgId() const
PDG identifier.
void analyze(const edm::Event &, const edm::EventSetup &)
virtual float pt() const
transverse momentum
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:17
virtual int status() const
status word
MonitorElement * book2DD(const char *name, const char *title, int nchX, double lowX, double highX, int nchY, double lowY, double highY)
Book 2D double histogram.
Definition: DQMStore.cc:1126
AcceptanceHistoProducer(const edm::ParameterSet &cfg)
#define NULL
Definition: scimark2.h:8
void Fill(long long x)
virtual float eta() const
momentum pseudorapidity
virtual size_t numberOfDaughters() const
number of daughters
virtual const Candidate * daughter(size_type) const
return daughter at a given position, i = 0, ... numberOfDaughters() - 1 (read only mode) ...
virtual int charge() const
electric charge
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
bool getByLabel(InputTag const &tag, Handle< PROD > &result) const
Definition: Event.h:390
part
Definition: HCALResponse.h:20
void setResetMe(bool)
void setLumiFlag(void)
this ME is meant to be stored for each luminosity section
void setCurrentFolder(const std::string &fullpath)
Definition: DQMStore.cc:667