CMS 3D CMS Logo

List of all members | Public Member Functions | Protected Member Functions | Private Member Functions | Private Attributes
DQMExample_Step1 Class Reference

#include <DQMExample_Step1.h>

Inheritance diagram for DQMExample_Step1:
one::DQMEDAnalyzer< T > one::dqmimplementation::DQMBaseClass< T... >

Public Member Functions

 DQMExample_Step1 (const edm::ParameterSet &ps)
 
 ~DQMExample_Step1 () override
 
- Public Member Functions inherited from one::DQMEDAnalyzer< T >
 DQMEDAnalyzer ()=default
 
 DQMEDAnalyzer (DQMEDAnalyzer< T... > const &)=delete
 
 DQMEDAnalyzer (DQMEDAnalyzer< T... > &&)=delete
 
 ~DQMEDAnalyzer () override=default
 

Protected Member Functions

void analyze (edm::Event const &e, edm::EventSetup const &eSetup) override
 
void bookHistograms (DQMStore::IBooker &, edm::Run const &, edm::EventSetup const &) override
 
void dqmBeginRun (edm::Run const &, edm::EventSetup const &) override
 
void endRun (edm::Run const &run, edm::EventSetup const &eSetup) override
 

Private Member Functions

void bookHistos (DQMStore::IBooker &)
 
double calcDeltaPhi (double phi1, double phi2)
 
double Distance (const reco::Candidate &c1, const reco::Candidate &c2)
 
double DistancePhi (const reco::Candidate &c1, const reco::Candidate &c2)
 
bool MediumEle (const edm::Event &iEvent, const edm::EventSetup &iESetup, const reco::GsfElectron &electron)
 

Private Attributes

MonitorElementh_eEta_leading
 
MonitorElementh_eEta_leading_HLT
 
MonitorElementh_eEta_leading_HLT_matched
 
MonitorElementh_eEta_leading_matched
 
MonitorElementh_eMultiplicity
 
MonitorElementh_eMultiplicity_HLT
 
MonitorElementh_ePhi_leading
 
MonitorElementh_ePhi_leading_HLT
 
MonitorElementh_ePhi_leading_HLT_matched
 
MonitorElementh_ePhi_leading_matched
 
MonitorElementh_ePt_diff
 
MonitorElementh_ePt_leading
 
MonitorElementh_ePt_leading_HLT
 
MonitorElementh_ePt_leading_HLT_matched
 
MonitorElementh_ePt_leading_matched
 
MonitorElementh_jEta_leading
 
MonitorElementh_jMultiplicity
 
MonitorElementh_jPhi_leading
 
MonitorElementh_jPt_leading
 
MonitorElementh_pfMet
 
MonitorElementh_vertex_number
 
int nBJets
 
int nElectrons
 
double ptThrJet_
 
double ptThrL1_
 
double ptThrL2_
 
double ptThrMet_
 
math::XYZPoint PVPoint_
 
edm::EDGetTokenT< reco::BeamSpottheBSCollection_
 
edm::EDGetTokenT< reco::CaloJetCollectiontheCaloJetCollection_
 
edm::EDGetTokenT< reco::ConversionCollectiontheConversionCollection_
 
edm::EDGetTokenT< reco::GsfElectronCollectiontheElectronCollection_
 
edm::EDGetTokenT< reco::PFMETCollectionthePfMETCollection_
 
edm::EDGetTokenT< reco::VertexCollectionthePVCollection_
 
edm::EDGetTokenT< trigger::TriggerEventtriggerEvent_
 
edm::InputTag triggerFilter_
 
std::string triggerPath_
 
edm::EDGetTokenT< edm::TriggerResultstriggerResults_
 

Detailed Description

Definition at line 47 of file DQMExample_Step1.h.

Constructor & Destructor Documentation

DQMExample_Step1::DQMExample_Step1 ( const edm::ParameterSet ps)

Definition at line 21 of file DQMExample_Step1.cc.

References edm::ParameterSet::getParameter(), edm::ParameterSet::getUntrackedParameter(), ptThrJet_, ptThrL1_, ptThrL2_, ptThrMet_, AlCaHLTBitMon_QueryRunRegistry::string, theBSCollection_, theCaloJetCollection_, theConversionCollection_, theElectronCollection_, thePfMETCollection_, thePVCollection_, triggerEvent_, triggerFilter_, triggerPath_, and triggerResults_.

21  {
22  edm::LogInfo("DQMExample_Step1") << "Constructor DQMExample_Step1::DQMExample_Step1 " << std::endl;
23 
24  // Get parameters from configuration file
25  theElectronCollection_ = consumes<reco::GsfElectronCollection>(ps.getParameter<edm::InputTag>("electronCollection"));
26  theCaloJetCollection_ = consumes<reco::CaloJetCollection>(ps.getParameter<edm::InputTag>("caloJetCollection"));
27  thePfMETCollection_ = consumes<reco::PFMETCollection>(ps.getParameter<edm::InputTag>("pfMETCollection"));
29  consumes<reco::ConversionCollection>(ps.getParameter<edm::InputTag>("conversionsCollection"));
30  thePVCollection_ = consumes<reco::VertexCollection>(ps.getParameter<edm::InputTag>("PVCollection"));
31  theBSCollection_ = consumes<reco::BeamSpot>(ps.getParameter<edm::InputTag>("beamSpotCollection"));
32  triggerEvent_ = consumes<trigger::TriggerEvent>(ps.getParameter<edm::InputTag>("TriggerEvent"));
33  triggerResults_ = consumes<edm::TriggerResults>(ps.getParameter<edm::InputTag>("TriggerResults"));
34  triggerFilter_ = ps.getParameter<edm::InputTag>("TriggerFilter");
35  triggerPath_ = ps.getParameter<std::string>("TriggerPath");
36 
37  // cuts:
38  ptThrL1_ = ps.getUntrackedParameter<double>("PtThrL1");
39  ptThrL2_ = ps.getUntrackedParameter<double>("PtThrL2");
40  ptThrJet_ = ps.getUntrackedParameter<double>("PtThrJet");
41  ptThrMet_ = ps.getUntrackedParameter<double>("PtThrMet");
42 }
T getParameter(std::string const &) const
T getUntrackedParameter(std::string const &, T const &) const
edm::EDGetTokenT< reco::GsfElectronCollection > theElectronCollection_
edm::EDGetTokenT< reco::ConversionCollection > theConversionCollection_
edm::EDGetTokenT< trigger::TriggerEvent > triggerEvent_
edm::EDGetTokenT< reco::BeamSpot > theBSCollection_
edm::EDGetTokenT< edm::TriggerResults > triggerResults_
edm::EDGetTokenT< reco::PFMETCollection > thePfMETCollection_
edm::InputTag triggerFilter_
std::string triggerPath_
edm::EDGetTokenT< reco::VertexCollection > thePVCollection_
edm::EDGetTokenT< reco::CaloJetCollection > theCaloJetCollection_
DQMExample_Step1::~DQMExample_Step1 ( )
override

Definition at line 47 of file DQMExample_Step1.cc.

47  {
48  edm::LogInfo("DQMExample_Step1") << "Destructor DQMExample_Step1::~DQMExample_Step1 " << std::endl;
49 }

Member Function Documentation

void DQMExample_Step1::analyze ( edm::Event const &  e,
edm::EventSetup const &  eSetup 
)
overrideprotected

Definition at line 72 of file DQMExample_Step1.cc.

References funct::abs(), edm::HLTGlobalStatus::accept(), boostedElectronIsolation_cff::deltaR, Distance(), ElectronMcFakeValidator_cfi::electronCollection, reco::Particle::eta(), reco::LeafCandidate::eta(), edm::HandleBase::failedToGet(), MonitorElement::Fill(), trigger::TriggerEvent::filterIndex(), trigger::TriggerEvent::filterKeys(), edm::Event::getByToken(), trigger::TriggerEvent::getObjects(), h_eEta_leading, h_eEta_leading_HLT, h_eEta_leading_HLT_matched, h_eEta_leading_matched, h_eMultiplicity, h_eMultiplicity_HLT, h_ePhi_leading, h_ePhi_leading_HLT, h_ePhi_leading_HLT_matched, h_ePhi_leading_matched, h_ePt_diff, h_ePt_leading, h_ePt_leading_HLT, h_ePt_leading_HLT_matched, h_ePt_leading_matched, h_jEta_leading, h_jMultiplicity, h_jPhi_leading, h_jPt_leading, h_pfMet, h_vertex_number, edm::HandleBase::isValid(), relativeConstraints::keys, MediumEle(), trigger::TriggerObject::particle(), reco::Particle::pdgId(), reco::Particle::phi(), reco::LeafCandidate::phi(), reco::Particle::pt(), reco::LeafCandidate::pt(), ptThrJet_, ptThrL1_, ptThrL2_, PVPoint_, edm::TriggerNames::size(), trigger::TriggerEvent::sizeFilters(), theCaloJetCollection_, theElectronCollection_, thePfMETCollection_, thePVCollection_, triggerEvent_, triggerFilter_, edm::TriggerNames::triggerName(), edm::Event::triggerNames(), TriggerAnalyzer::triggerObjects, triggerPath_, triggerResults_, trigNames, findQualityFiles::v, and edm::HLTGlobalStatus::wasrun().

72  {
73  edm::LogInfo("DQMExample_Step1") << "DQMExample_Step1::analyze" << std::endl;
74 
75  //-------------------------------
76  //--- Vertex Info
77  //-------------------------------
79  e.getByToken(thePVCollection_, vertexHandle);
80  if (!vertexHandle.isValid()) {
81  edm::LogError("DQMClientExample") << "invalid collection: vertex"
82  << "\n";
83  return;
84  }
85 
86  int vertex_number = vertexHandle->size();
87  reco::VertexCollection::const_iterator v = vertexHandle->begin();
88 
89  math::XYZPoint PVPoint(-999, -999, -999);
90  if (vertex_number != 0)
91  PVPoint = math::XYZPoint(v->position().x(), v->position().y(), v->position().z());
92 
93  PVPoint_ = PVPoint;
94 
95  //-------------------------------
96  //--- MET
97  //-------------------------------
98  edm::Handle<reco::PFMETCollection> pfMETCollection;
99  e.getByToken(thePfMETCollection_, pfMETCollection);
100  if (!pfMETCollection.isValid()) {
101  edm::LogError("DQMClientExample") << "invalid collection: MET"
102  << "\n";
103  return;
104  }
105  //-------------------------------
106  //--- Electrons
107  //-------------------------------
109  e.getByToken(theElectronCollection_, electronCollection);
110  if (!electronCollection.isValid()) {
111  edm::LogError("DQMClientExample") << "invalid collection: electrons"
112  << "\n";
113  return;
114  }
115 
116  float nEle = 0;
117  int posEle = 0, negEle = 0;
118  const reco::GsfElectron *ele1 = nullptr;
119  const reco::GsfElectron *ele2 = nullptr;
120  for (reco::GsfElectronCollection::const_iterator recoElectron = electronCollection->begin();
121  recoElectron != electronCollection->end();
122  ++recoElectron) {
123  // decreasing pT
124  if (MediumEle(e, eSetup, *recoElectron)) {
125  if (!ele1 && recoElectron->pt() > ptThrL1_)
126  ele1 = &(*recoElectron);
127 
128  else if (!ele2 && recoElectron->pt() > ptThrL2_)
129  ele2 = &(*recoElectron);
130  }
131 
132  if (recoElectron->charge() == 1)
133  posEle++;
134  else if (recoElectron->charge() == -1)
135  negEle++;
136 
137  } // end of loop over electrons
138 
139  nEle = posEle + negEle;
140 
141  //-------------------------------
142  //--- Jets
143  //-------------------------------
144  edm::Handle<reco::CaloJetCollection> caloJetCollection;
145  e.getByToken(theCaloJetCollection_, caloJetCollection);
146  if (!caloJetCollection.isValid()) {
147  edm::LogError("DQMClientExample") << "invalid collection: jets"
148  << "\n";
149  return;
150  }
151 
152  int nJet = 0;
153  const reco::CaloJet *jet1 = nullptr;
154  const reco::CaloJet *jet2 = nullptr;
155 
156  for (reco::CaloJetCollection::const_iterator i_calojet = caloJetCollection->begin();
157  i_calojet != caloJetCollection->end();
158  ++i_calojet) {
159  // remove jet-ele matching
160  if (ele1)
161  if (Distance(*i_calojet, *ele1) < 0.3)
162  continue;
163 
164  if (ele2)
165  if (Distance(*i_calojet, *ele2) < 0.3)
166  continue;
167 
168  if (i_calojet->pt() < ptThrJet_)
169  continue;
170 
171  nJet++;
172 
173  if (!jet1)
174  jet1 = &(*i_calojet);
175 
176  else if (!jet2)
177  jet2 = &(*i_calojet);
178  }
179 
180  // ---------------------------
181  // ---- Analyze Trigger Event
182  // ---------------------------
183 
184  // check what is in the menu
186  e.getByToken(triggerResults_, hltresults);
187 
188  if (!hltresults.isValid()) {
189  edm::LogError("DQMClientExample") << "invalid collection: TriggerResults"
190  << "\n";
191  return;
192  }
193 
194  bool hasFired = false;
195  const edm::TriggerNames &trigNames = e.triggerNames(*hltresults);
196  unsigned int numTriggers = trigNames.size();
197 
198  for (unsigned int hltIndex = 0; hltIndex < numTriggers; ++hltIndex) {
199  if (trigNames.triggerName(hltIndex) == triggerPath_ && hltresults->wasrun(hltIndex) && hltresults->accept(hltIndex))
200  hasFired = true;
201  }
202 
203  // access the trigger event
205  e.getByToken(triggerEvent_, triggerEvent);
206  if (triggerEvent.failedToGet()) {
207  edm::LogError("DQMClientExample") << "invalid collection: TriggerEvent"
208  << "\n";
209  return;
210  }
211 
212  reco::Particle *ele1_HLT = nullptr;
213  int nEle_HLT = 0;
214 
215  size_t filterIndex = triggerEvent->filterIndex(triggerFilter_);
217  if (!(filterIndex >= triggerEvent->sizeFilters())) {
218  const trigger::Keys &keys = triggerEvent->filterKeys(filterIndex);
219  std::vector<reco::Particle> triggeredEle;
220 
221  for (size_t j = 0; j < keys.size(); ++j) {
222  trigger::TriggerObject foundObject = triggerObjects[keys[j]];
223  if (abs(foundObject.particle().pdgId()) != 11)
224  continue; // make sure that it is an electron
225 
226  triggeredEle.push_back(foundObject.particle());
227  ++nEle_HLT;
228  }
229 
230  if (!triggeredEle.empty())
231  ele1_HLT = &(triggeredEle.at(0));
232  }
233 
234  //-------------------------------
235  //--- Fill the histos
236  //-------------------------------
237 
238  // vertex
239  h_vertex_number->Fill(vertex_number);
240 
241  // met
242  h_pfMet->Fill(pfMETCollection->begin()->et());
243 
244  // multiplicities
245  h_eMultiplicity->Fill(nEle);
246  h_jMultiplicity->Fill(nJet);
247  h_eMultiplicity_HLT->Fill(nEle_HLT);
248 
249  // leading not matched
250  if (ele1) {
251  h_ePt_leading->Fill(ele1->pt());
252  h_eEta_leading->Fill(ele1->eta());
253  h_ePhi_leading->Fill(ele1->phi());
254  }
255  if (ele1_HLT) {
256  h_ePt_leading_HLT->Fill(ele1_HLT->pt());
257  h_eEta_leading_HLT->Fill(ele1_HLT->eta());
258  h_ePhi_leading_HLT->Fill(ele1_HLT->phi());
259  }
260  // leading Jet
261  if (jet1) {
262  h_jPt_leading->Fill(jet1->pt());
263  h_jEta_leading->Fill(jet1->eta());
264  h_jPhi_leading->Fill(jet1->phi());
265  }
266 
267  // fill only when the trigger candidate mathes with the reco one
268  if (ele1 && ele1_HLT && deltaR(*ele1_HLT, *ele1) < 0.3 && hasFired == true) {
269  h_ePt_leading_matched->Fill(ele1->pt());
270  h_eEta_leading_matched->Fill(ele1->eta());
271  h_ePhi_leading_matched->Fill(ele1->phi());
272 
273  h_ePt_leading_HLT_matched->Fill(ele1_HLT->pt());
274  h_eEta_leading_HLT_matched->Fill(ele1_HLT->eta());
275  h_ePhi_leading_HLT_matched->Fill(ele1_HLT->phi());
276 
277  h_ePt_diff->Fill(ele1->pt() - ele1_HLT->pt());
278  }
279 }
MonitorElement * h_ePhi_leading_matched
bool wasrun() const
Was at least one path run?
double pt() const
transverse momentum
Definition: Particle.h:104
MonitorElement * h_jMultiplicity
math::XYZPoint PVPoint_
bool MediumEle(const edm::Event &iEvent, const edm::EventSetup &iESetup, const reco::GsfElectron &electron)
double eta() const final
momentum pseudorapidity
edm::EDGetTokenT< reco::GsfElectronCollection > theElectronCollection_
Jets made from CaloTowers.
Definition: CaloJet.h:29
trigger::size_type sizeFilters() const
Definition: TriggerEvent.h:135
bool accept() const
Has at least one path accepted the event?
const Keys & filterKeys(trigger::size_type index) const
Definition: TriggerEvent.h:111
trigger::size_type filterIndex(const edm::InputTag &filterTag) const
find index of filter in data-member vector from filter tag
Definition: TriggerEvent.h:123
reco::Particle particle(reco::Particle::Charge q=0, const reco::Particle::Point &vertex=reco::Particle::Point(0, 0, 0), int status=0, bool integerCharge=true) const
Definition: TriggerObject.h:69
MonitorElement * h_ePhi_leading_HLT_matched
double pt() const final
transverse momentum
MonitorElement * h_ePhi_leading
MonitorElement * h_eEta_leading_HLT
MonitorElement * h_jPhi_leading
Strings::size_type size() const
Definition: TriggerNames.cc:31
edm::EDGetTokenT< trigger::TriggerEvent > triggerEvent_
double Distance(const reco::Candidate &c1, const reco::Candidate &c2)
int pdgId() const
PDG identifier.
Definition: Particle.h:134
void Fill(long long x)
MonitorElement * h_eEta_leading
MonitorElement * h_ePt_leading_HLT
Single trigger physics object (e.g., an isolated muon)
Definition: TriggerObject.h:22
MonitorElement * h_ePt_diff
MonitorElement * h_ePt_leading_HLT_matched
double phi() const
momentum azimuthal angle
Definition: Particle.h:106
MonitorElement * h_jEta_leading
const TriggerObjectCollection & getObjects() const
Definition: TriggerEvent.h:98
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
MonitorElement * h_eEta_leading_HLT_matched
MonitorElement * h_vertex_number
MonitorElement * h_ePt_leading
bool isValid() const
Definition: HandleBase.h:74
MonitorElement * h_eMultiplicity
MonitorElement * h_pfMet
std::vector< TriggerObject > TriggerObjectCollection
collection of trigger physics objects (e.g., all isolated muons)
Definition: TriggerObject.h:81
static const char *const trigNames[]
Definition: EcalDumpRaw.cc:74
bool failedToGet() const
Definition: HandleBase.h:78
edm::EDGetTokenT< edm::TriggerResults > triggerResults_
MonitorElement * h_jPt_leading
std::string const & triggerName(unsigned int index) const
Definition: TriggerNames.cc:22
XYZPointD XYZPoint
point in space with cartesian internal representation
Definition: Point3D.h:12
std::vector< size_type > Keys
MonitorElement * h_ePt_leading_matched
edm::EDGetTokenT< reco::PFMETCollection > thePfMETCollection_
edm::InputTag triggerFilter_
MonitorElement * h_ePhi_leading_HLT
std::string triggerPath_
MonitorElement * h_eMultiplicity_HLT
double eta() const
momentum pseudorapidity
Definition: Particle.h:110
double phi() const final
momentum azimuthal angle
MonitorElement * h_eEta_leading_matched
edm::EDGetTokenT< reco::VertexCollection > thePVCollection_
edm::EDGetTokenT< reco::CaloJetCollection > theCaloJetCollection_
void DQMExample_Step1::bookHistograms ( DQMStore::IBooker ibooker_,
edm::Run const &  ,
edm::EventSetup const &   
)
overrideprotected

Definition at line 62 of file DQMExample_Step1.cc.

References bookHistos().

62  {
63  edm::LogInfo("DQMExample_Step1") << "DQMExample_Step1::bookHistograms" << std::endl;
64 
65  // book at beginRun
66  bookHistos(ibooker_);
67 }
void bookHistos(DQMStore::IBooker &)
void DQMExample_Step1::bookHistos ( DQMStore::IBooker ibooker_)
private

Definition at line 293 of file DQMExample_Step1.cc.

References DQMStore::IBooker::book1D(), DQMStore::IBooker::cd(), h_eEta_leading, h_eEta_leading_HLT, h_eEta_leading_HLT_matched, h_eEta_leading_matched, h_eMultiplicity, h_eMultiplicity_HLT, h_ePhi_leading, h_ePhi_leading_HLT, h_ePhi_leading_HLT_matched, h_ePhi_leading_matched, h_ePt_diff, h_ePt_leading, h_ePt_leading_HLT, h_ePt_leading_HLT_matched, h_ePt_leading_matched, h_jEta_leading, h_jMultiplicity, h_jPhi_leading, h_jPt_leading, h_pfMet, h_vertex_number, and DQMStore::IBooker::setCurrentFolder().

Referenced by bookHistograms().

293  {
294  ibooker_.cd();
295  ibooker_.setCurrentFolder("Physics/TopTest");
296 
297  h_vertex_number = ibooker_.book1D("Vertex_number", "Number of event vertices in collection", 40, -0.5, 39.5);
298 
299  h_pfMet = ibooker_.book1D("pfMet", "Pf Missing E_{T}; GeV", 20, 0.0, 100);
300 
301  h_eMultiplicity = ibooker_.book1D("NElectrons", "# of electrons per event", 10, 0., 10.);
302  h_ePt_leading_matched = ibooker_.book1D("ElePt_leading_matched", "Pt of leading electron", 50, 0., 100.);
303  h_eEta_leading_matched = ibooker_.book1D("EleEta_leading_matched", "Eta of leading electron", 50, -5., 5.);
304  h_ePhi_leading_matched = ibooker_.book1D("ElePhi_leading_matched", "Phi of leading electron", 50, -3.5, 3.5);
305 
306  h_ePt_leading = ibooker_.book1D("ElePt_leading", "Pt of leading electron", 50, 0., 100.);
307  h_eEta_leading = ibooker_.book1D("EleEta_leading", "Eta of leading electron", 50, -5., 5.);
308  h_ePhi_leading = ibooker_.book1D("ElePhi_leading", "Phi of leading electron", 50, -3.5, 3.5);
309 
310  h_jMultiplicity = ibooker_.book1D("NJets", "# of electrons per event", 10, 0., 10.);
311  h_jPt_leading = ibooker_.book1D("JetPt_leading", "Pt of leading Jet", 150, 0., 300.);
312  h_jEta_leading = ibooker_.book1D("JetEta_leading", "Eta of leading Jet", 50, -5., 5.);
313  h_jPhi_leading = ibooker_.book1D("JetPhi_leading", "Phi of leading Jet", 50, -3.5, 3.5);
314 
315  h_eMultiplicity_HLT = ibooker_.book1D("NElectrons_HLT", "# of electrons per event @HLT", 10, 0., 10.);
316  h_ePt_leading_HLT = ibooker_.book1D("ElePt_leading_HLT", "Pt of leading electron @HLT", 50, 0., 100.);
317  h_eEta_leading_HLT = ibooker_.book1D("EleEta_leading_HLT", "Eta of leading electron @HLT", 50, -5., 5.);
318  h_ePhi_leading_HLT = ibooker_.book1D("ElePhi_leading_HLT", "Phi of leading electron @HLT", 50, -3.5, 3.5);
319 
320  h_ePt_leading_HLT_matched = ibooker_.book1D("ElePt_leading_HLT_matched", "Pt of leading electron @HLT", 50, 0., 100.);
322  ibooker_.book1D("EleEta_leading_HLT_matched", "Eta of leading electron @HLT", 50, -5., 5.);
324  ibooker_.book1D("ElePhi_leading_HLT_matched", "Phi of leading electron @HLT", 50, -3.5, 3.5);
325 
326  h_ePt_diff = ibooker_.book1D("ElePt_diff_matched", "pT(RECO) - pT(HLT) for mathed candidates", 100, -10, 10.);
327 
328  ibooker_.cd();
329 }
MonitorElement * h_ePhi_leading_matched
MonitorElement * h_jMultiplicity
MonitorElement * h_ePhi_leading_HLT_matched
MonitorElement * h_ePhi_leading
MonitorElement * h_eEta_leading_HLT
MonitorElement * h_jPhi_leading
MonitorElement * h_eEta_leading
MonitorElement * h_ePt_leading_HLT
MonitorElement * h_ePt_diff
MonitorElement * h_ePt_leading_HLT_matched
MonitorElement * h_jEta_leading
void setCurrentFolder(std::string const &fullpath)
Definition: DQMStore.cc:268
MonitorElement * book1D(Args &&...args)
Definition: DQMStore.h:106
MonitorElement * h_eEta_leading_HLT_matched
MonitorElement * h_vertex_number
MonitorElement * h_ePt_leading
MonitorElement * h_eMultiplicity
MonitorElement * h_pfMet
MonitorElement * h_jPt_leading
MonitorElement * h_ePt_leading_matched
MonitorElement * h_ePhi_leading_HLT
MonitorElement * h_eMultiplicity_HLT
MonitorElement * h_eEta_leading_matched
double DQMExample_Step1::calcDeltaPhi ( double  phi1,
double  phi2 
)
private

Definition at line 342 of file DQMExample_Step1.cc.

References hiPixelPairStep_cff::deltaPhi.

342  {
343  double deltaPhi = phi1 - phi2;
344  if (deltaPhi < 0)
345  deltaPhi = -deltaPhi;
346  if (deltaPhi > 3.1415926) {
347  deltaPhi = 2 * 3.1415926 - deltaPhi;
348  }
349  return deltaPhi;
350 }
double DQMExample_Step1::Distance ( const reco::Candidate c1,
const reco::Candidate c2 
)
private

Definition at line 335 of file DQMExample_Step1.cc.

References boostedElectronIsolation_cff::deltaR.

Referenced by analyze().

double DQMExample_Step1::DistancePhi ( const reco::Candidate c1,
const reco::Candidate c2 
)
private

Definition at line 337 of file DQMExample_Step1.cc.

References hiPixelPairStep_cff::deltaPhi, and reco::Candidate::p4().

337  {
338  return deltaPhi(c1.p4().phi(), c2.p4().phi());
339 }
virtual const LorentzVector & p4() const =0
four-momentum Lorentz vector
void DQMExample_Step1::dqmBeginRun ( edm::Run const &  ,
edm::EventSetup const &   
)
overrideprotected

Definition at line 55 of file DQMExample_Step1.cc.

55  {
56  edm::LogInfo("DQMExample_Step1") << "DQMExample_Step1::beginRun" << std::endl;
57 }
void DQMExample_Step1::endRun ( edm::Run const &  run,
edm::EventSetup const &  eSetup 
)
overrideprotected

Definition at line 285 of file DQMExample_Step1.cc.

285  {
286  edm::LogInfo("DQMExample_Step1") << "DQMExample_Step1::endRun" << std::endl;
287 }
bool DQMExample_Step1::MediumEle ( const edm::Event iEvent,
const edm::EventSetup iESetup,
const reco::GsfElectron electron 
)
private

Definition at line 356 of file DQMExample_Step1.cc.

References reco::GsfElectron::ambiguousGsfTracksSize(), reco::GsfElectron::deltaEtaSuperClusterTrackAtVtx(), reco::GsfElectron::deltaPhiSuperClusterTrackAtVtx(), PVValHelper::dxy, PVValHelper::dz, reco::GsfElectron::ecalEnergy(), reco::GsfElectron::eSuperClusterOverP(), PVValHelper::eta, reco::LeafCandidate::eta(), edm::Event::getByToken(), reco::GsfElectron::gsfTrack(), reco::GsfElectron::hadronicOverEm(), ConversionTools::hasMatchedConversion(), reco::GsfElectron::isEB(), reco::HitPattern::MISSING_INNER_HITS, reco::BeamSpot::position(), EnergyCorrector::pt, reco::LeafCandidate::pt(), PVPoint_, reco::GsfElectron::sigmaIetaIeta(), theBSCollection_, and theConversionCollection_.

Referenced by analyze().

358  {
359  //********* CONVERSION TOOLS
361  iEvent.getByToken(theConversionCollection_, conversions_h);
362 
363  bool isMediumEle = false;
364 
365  float pt = electron.pt();
366  float eta = electron.eta();
367 
368  int isEB = electron.isEB();
369  float sigmaIetaIeta = electron.sigmaIetaIeta();
370  float DetaIn = electron.deltaEtaSuperClusterTrackAtVtx();
371  float DphiIn = electron.deltaPhiSuperClusterTrackAtVtx();
372  float HOverE = electron.hadronicOverEm();
373  float ooemoop = (1.0 / electron.ecalEnergy() - electron.eSuperClusterOverP() / electron.ecalEnergy());
374 
375  int mishits = electron.gsfTrack()->hitPattern().numberOfLostHits(reco::HitPattern::MISSING_INNER_HITS);
376  int nAmbiguousGsfTracks = electron.ambiguousGsfTracksSize();
377 
378  reco::GsfTrackRef eleTrack = electron.gsfTrack();
379  float dxy = eleTrack->dxy(PVPoint_);
380  float dz = eleTrack->dz(PVPoint_);
381 
383  iEvent.getByToken(theBSCollection_, BSHandle);
384  const reco::BeamSpot BS = *BSHandle;
385 
386  bool isConverted = ConversionTools::hasMatchedConversion(electron, *conversions_h, BS.position());
387 
388  // default
389  if ((pt > 12.) && (fabs(eta) < 2.5) &&
390  (((isEB == 1) && (fabs(DetaIn) < 0.004)) || ((isEB == 0) && (fabs(DetaIn) < 0.007))) &&
391  (((isEB == 1) && (fabs(DphiIn) < 0.060)) || ((isEB == 0) && (fabs(DphiIn) < 0.030))) &&
392  (((isEB == 1) && (sigmaIetaIeta < 0.010)) || ((isEB == 0) && (sigmaIetaIeta < 0.030))) &&
393  (((isEB == 1) && (HOverE < 0.120)) || ((isEB == 0) && (HOverE < 0.100))) &&
394  (((isEB == 1) && (fabs(ooemoop) < 0.050)) || ((isEB == 0) && (fabs(ooemoop) < 0.050))) &&
395  (((isEB == 1) && (fabs(dxy) < 0.020)) || ((isEB == 0) && (fabs(dxy) < 0.020))) &&
396  (((isEB == 1) && (fabs(dz) < 0.100)) || ((isEB == 0) && (fabs(dz) < 0.100))) &&
397  (((isEB == 1) && (!isConverted)) || ((isEB == 0) && (!isConverted))) && (mishits == 0) &&
398  (nAmbiguousGsfTracks == 0))
399  isMediumEle = true;
400 
401  return isMediumEle;
402 }
GsfTrackRef gsfTrack() const override
reference to a GsfTrack
Definition: GsfElectron.h:186
math::XYZPoint PVPoint_
double eta() const final
momentum pseudorapidity
float eSuperClusterOverP() const
Definition: GsfElectron.h:249
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: Event.h:517
double pt() const final
transverse momentum
edm::EDGetTokenT< reco::ConversionCollection > theConversionCollection_
bool isEB() const
Definition: GsfElectron.h:356
edm::EDGetTokenT< reco::BeamSpot > theBSCollection_
float deltaEtaSuperClusterTrackAtVtx() const
Definition: GsfElectron.h:253
float sigmaIetaIeta() const
Definition: GsfElectron.h:440
float hadronicOverEm() const
Definition: GsfElectron.h:495
static bool hasMatchedConversion(const reco::GsfElectron &ele, const reco::ConversionCollection &convCol, const math::XYZPoint &beamspot, bool allowCkfMatch=true, float lxyMin=2.0, float probMin=1e-6, unsigned int nHitsBeforeVtxMax=0)
float deltaPhiSuperClusterTrackAtVtx() const
Definition: GsfElectron.h:256
GsfTrackRefVector::size_type ambiguousGsfTracksSize() const
Definition: GsfElectron.h:722
float ecalEnergy() const
Definition: GsfElectron.h:858
const Point & position() const
position
Definition: BeamSpot.h:62

Member Data Documentation

MonitorElement* DQMExample_Step1::h_eEta_leading
private

Definition at line 98 of file DQMExample_Step1.h.

Referenced by analyze(), and bookHistos().

MonitorElement* DQMExample_Step1::h_eEta_leading_HLT
private

Definition at line 106 of file DQMExample_Step1.h.

Referenced by analyze(), and bookHistos().

MonitorElement* DQMExample_Step1::h_eEta_leading_HLT_matched
private

Definition at line 109 of file DQMExample_Step1.h.

Referenced by analyze(), and bookHistos().

MonitorElement* DQMExample_Step1::h_eEta_leading_matched
private

Definition at line 101 of file DQMExample_Step1.h.

Referenced by analyze(), and bookHistos().

MonitorElement* DQMExample_Step1::h_eMultiplicity
private

Definition at line 96 of file DQMExample_Step1.h.

Referenced by analyze(), and bookHistos().

MonitorElement* DQMExample_Step1::h_eMultiplicity_HLT
private

Definition at line 104 of file DQMExample_Step1.h.

Referenced by analyze(), and bookHistos().

MonitorElement* DQMExample_Step1::h_ePhi_leading
private

Definition at line 99 of file DQMExample_Step1.h.

Referenced by analyze(), and bookHistos().

MonitorElement* DQMExample_Step1::h_ePhi_leading_HLT
private

Definition at line 107 of file DQMExample_Step1.h.

Referenced by analyze(), and bookHistos().

MonitorElement* DQMExample_Step1::h_ePhi_leading_HLT_matched
private

Definition at line 110 of file DQMExample_Step1.h.

Referenced by analyze(), and bookHistos().

MonitorElement* DQMExample_Step1::h_ePhi_leading_matched
private

Definition at line 102 of file DQMExample_Step1.h.

Referenced by analyze(), and bookHistos().

MonitorElement* DQMExample_Step1::h_ePt_diff
private

Definition at line 117 of file DQMExample_Step1.h.

Referenced by analyze(), and bookHistos().

MonitorElement* DQMExample_Step1::h_ePt_leading
private

Definition at line 97 of file DQMExample_Step1.h.

Referenced by analyze(), and bookHistos().

MonitorElement* DQMExample_Step1::h_ePt_leading_HLT
private

Definition at line 105 of file DQMExample_Step1.h.

Referenced by analyze(), and bookHistos().

MonitorElement* DQMExample_Step1::h_ePt_leading_HLT_matched
private

Definition at line 108 of file DQMExample_Step1.h.

Referenced by analyze(), and bookHistos().

MonitorElement* DQMExample_Step1::h_ePt_leading_matched
private

Definition at line 100 of file DQMExample_Step1.h.

Referenced by analyze(), and bookHistos().

MonitorElement* DQMExample_Step1::h_jEta_leading
private

Definition at line 114 of file DQMExample_Step1.h.

Referenced by analyze(), and bookHistos().

MonitorElement* DQMExample_Step1::h_jMultiplicity
private

Definition at line 112 of file DQMExample_Step1.h.

Referenced by analyze(), and bookHistos().

MonitorElement* DQMExample_Step1::h_jPhi_leading
private

Definition at line 115 of file DQMExample_Step1.h.

Referenced by analyze(), and bookHistos().

MonitorElement* DQMExample_Step1::h_jPt_leading
private

Definition at line 113 of file DQMExample_Step1.h.

Referenced by analyze(), and bookHistos().

MonitorElement* DQMExample_Step1::h_pfMet
private

Definition at line 94 of file DQMExample_Step1.h.

Referenced by analyze(), and bookHistos().

MonitorElement* DQMExample_Step1::h_vertex_number
private

Definition at line 92 of file DQMExample_Step1.h.

Referenced by analyze(), and bookHistos().

int DQMExample_Step1::nBJets
private

Definition at line 89 of file DQMExample_Step1.h.

int DQMExample_Step1::nElectrons
private

Definition at line 88 of file DQMExample_Step1.h.

double DQMExample_Step1::ptThrJet_
private

Definition at line 85 of file DQMExample_Step1.h.

Referenced by analyze(), and DQMExample_Step1().

double DQMExample_Step1::ptThrL1_
private

Definition at line 83 of file DQMExample_Step1.h.

Referenced by analyze(), and DQMExample_Step1().

double DQMExample_Step1::ptThrL2_
private

Definition at line 84 of file DQMExample_Step1.h.

Referenced by analyze(), and DQMExample_Step1().

double DQMExample_Step1::ptThrMet_
private

Definition at line 86 of file DQMExample_Step1.h.

Referenced by DQMExample_Step1().

math::XYZPoint DQMExample_Step1::PVPoint_
private

Definition at line 69 of file DQMExample_Step1.h.

Referenced by analyze(), and MediumEle().

edm::EDGetTokenT<reco::BeamSpot> DQMExample_Step1::theBSCollection_
private

Definition at line 77 of file DQMExample_Step1.h.

Referenced by DQMExample_Step1(), and MediumEle().

edm::EDGetTokenT<reco::CaloJetCollection> DQMExample_Step1::theCaloJetCollection_
private

Definition at line 74 of file DQMExample_Step1.h.

Referenced by analyze(), and DQMExample_Step1().

edm::EDGetTokenT<reco::ConversionCollection> DQMExample_Step1::theConversionCollection_
private

Definition at line 73 of file DQMExample_Step1.h.

Referenced by DQMExample_Step1(), and MediumEle().

edm::EDGetTokenT<reco::GsfElectronCollection> DQMExample_Step1::theElectronCollection_
private

Definition at line 72 of file DQMExample_Step1.h.

Referenced by analyze(), and DQMExample_Step1().

edm::EDGetTokenT<reco::PFMETCollection> DQMExample_Step1::thePfMETCollection_
private

Definition at line 75 of file DQMExample_Step1.h.

Referenced by analyze(), and DQMExample_Step1().

edm::EDGetTokenT<reco::VertexCollection> DQMExample_Step1::thePVCollection_
private

Definition at line 76 of file DQMExample_Step1.h.

Referenced by analyze(), and DQMExample_Step1().

edm::EDGetTokenT<trigger::TriggerEvent> DQMExample_Step1::triggerEvent_
private

Definition at line 78 of file DQMExample_Step1.h.

Referenced by analyze(), and DQMExample_Step1().

edm::InputTag DQMExample_Step1::triggerFilter_
private

Definition at line 80 of file DQMExample_Step1.h.

Referenced by analyze(), and DQMExample_Step1().

std::string DQMExample_Step1::triggerPath_
private

Definition at line 81 of file DQMExample_Step1.h.

Referenced by analyze(), and DQMExample_Step1().

edm::EDGetTokenT<edm::TriggerResults> DQMExample_Step1::triggerResults_
private

Definition at line 79 of file DQMExample_Step1.h.

Referenced by analyze(), and DQMExample_Step1().