CMS 3D CMS Logo

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

#include <L1TStage2uGTCaloLayer2Comp.h>

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

Public Member Functions

 L1TStage2uGTCaloLayer2Comp (const edm::ParameterSet &ps)
 
- 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 (const edm::Event &, const edm::EventSetup &) override
 
void bookHistograms (DQMStore::IBooker &, const edm::Run &, const edm::EventSetup &) override
 

Private Types

enum  denumBins {
  EVENTS1 = 1, EVENTS2, EVENTS3, EVENTS4,
  EVENTS5, JETS1, JETS2, JETS3,
  EGS1, EGS2, EGS3, TAUS1,
  TAUS2, TAUS3, SUMS
}
 
enum  numeratorBins {
  EVENTBAD = 1, EVENTBADJETCOL, EVENTBADEGCOL, EVENTBADTAUCOL,
  EVENTBADSUMCOL, JETBADET, JETBADETA, JETBADPHI,
  EGBADET, EGBADETA, EGBADPHI, TAUBADET,
  TAUBADETA, TAUBADPHI, BADSUM
}
 

Private Member Functions

bool compareEGs (const edm::Handle< l1t::EGammaBxCollection > &col1, const edm::Handle< l1t::EGammaBxCollection > &col2)
 
bool compareJets (const edm::Handle< l1t::JetBxCollection > &col1, const edm::Handle< l1t::JetBxCollection > &col2)
 
bool compareSums (const edm::Handle< l1t::EtSumBxCollection > &col1, const edm::Handle< l1t::EtSumBxCollection > &col2)
 
bool compareTaus (const edm::Handle< l1t::TauBxCollection > &col1, const edm::Handle< l1t::TauBxCollection > &col2)
 

Private Attributes

std::string collection1Title
 
std::string collection2Title
 
MonitorElementcomparisonDenum
 
MonitorElementcomparisonNum
 
edm::EDGetTokenT< l1t::EGammaBxCollectionEGammaCollection1
 
edm::EDGetTokenT< l1t::EGammaBxCollectionEGammaCollection2
 
edm::EDGetTokenT< l1t::EtSumBxCollectionEtSumCollection1
 
edm::EDGetTokenT< l1t::EtSumBxCollectionEtSumCollection2
 
edm::EDGetTokenT< l1t::JetBxCollectionJetCollection1
 
edm::EDGetTokenT< l1t::JetBxCollectionJetCollection2
 
std::string monitorDir
 
edm::EDGetTokenT< l1t::TauBxCollectionTauCollection1
 
edm::EDGetTokenT< l1t::TauBxCollectionTauCollection2
 
bool verbose
 

Detailed Description

Class to perform event by event comparisons between CaloLayer2 ouputs and uGT inputs and produce summary of problems found

Module should be run as part of L1TStage2 online sequence and relies on collections of jets, e/g, tau and sum objects as they come out of CaloLayer2 and are unpacked by uGT. The purpose of the comparisions is to identify issues with the data transmission links or the unpacking process in the uGT FPGA firmare. The module is also used to compare the inputs of all uGT boards.

Summary differentiates between different types of errors and errors with different objects in attempt to identify issues with specific links. In the case of the sums, due to the large number of different types and their spread over many links, it was decided to not differentiate between the different types of objects and their properties for the current implementation.

Definition at line 34 of file L1TStage2uGTCaloLayer2Comp.h.

Member Enumeration Documentation

Enumerator
EVENTS1 
EVENTS2 
EVENTS3 
EVENTS4 
EVENTS5 
JETS1 
JETS2 
JETS3 
EGS1 
EGS2 
EGS3 
TAUS1 
TAUS2 
TAUS3 
SUMS 

Definition at line 207 of file L1TStage2uGTCaloLayer2Comp.h.

207  {
208  EVENTS1 = 1, // total no. events (used for taking a ratio) x5
209  EVENTS2,
210  EVENTS3,
211  EVENTS4,
212  EVENTS5,
213  JETS1, // total no. jets x3
214  JETS2,
215  JETS3,
216  EGS1, // total no. egs x3
217  EGS2,
218  EGS3,
219  TAUS1, // total no. taus x3
220  TAUS2,
221  TAUS3,
222  SUMS // total no. sums
223  };
Enumerator
EVENTBAD 
EVENTBADJETCOL 
EVENTBADEGCOL 
EVENTBADTAUCOL 
EVENTBADSUMCOL 
JETBADET 
JETBADETA 
JETBADPHI 
EGBADET 
EGBADETA 
EGBADPHI 
TAUBADET 
TAUBADETA 
TAUBADPHI 
BADSUM 

Definition at line 189 of file L1TStage2uGTCaloLayer2Comp.h.

189  {
190  EVENTBAD = 1, // number of (no.) bad events (where an error was found)
191  EVENTBADJETCOL, // no. events with a jet collection size difference
192  EVENTBADEGCOL, // no. events with a eg collection size difference
193  EVENTBADTAUCOL, // no. events with a tau collection size difference
194  EVENTBADSUMCOL, // no. events with a sum collection size difference
195  JETBADET, // no. jets with bad Et
196  JETBADETA, // no. jets with bad eta
197  JETBADPHI, // no. jets with bad phi
198  EGBADET, // no. egs with bad Et
199  EGBADETA, // no. egs with bad phi
200  EGBADPHI, // no. egs with bad eta
201  TAUBADET, // no. tau with bad Et
202  TAUBADETA, // no. tau with bad eta
203  TAUBADPHI, // no. tau with bad phi
204  BADSUM // no. sums with any disagreement
205  };

Constructor & Destructor Documentation

L1TStage2uGTCaloLayer2Comp::L1TStage2uGTCaloLayer2Comp ( const edm::ParameterSet ps)

Class constructor

Receives the set of parameters, specified by the python configuration file used to initialise the module as a part of a sequence. The contents of the set is used to configure the internal state of the objects of this class. Values from the parameter set are extracted and used to initialise bxcollections for jet, e/g, tau and sum objects reconstructed and unpacked by CaloLayer2 and uGT firmwares. These collections are the basis of the comparisons performed by this module.

Parameters
edm::ParamterSet& ps A pointer to the parameter set used

Definition at line 3 of file L1TStage2uGTCaloLayer2Comp.cc.

4  : monitorDir(ps.getUntrackedParameter<std::string>("monitorDir", "")),
5  collection1Title(ps.getUntrackedParameter<std::string>("collection1Title")),
6  collection2Title(ps.getUntrackedParameter<std::string>("collection2Title")),
7  JetCollection1(consumes <l1t::JetBxCollection>(ps.getParameter<edm::InputTag>("JetCollection1"))),
8  JetCollection2(consumes <l1t::JetBxCollection>(ps.getParameter<edm::InputTag>("JetCollection2"))),
9  EGammaCollection1(consumes <l1t::EGammaBxCollection>(ps.getParameter<edm::InputTag>("EGammaCollection1"))),
10  EGammaCollection2(consumes <l1t::EGammaBxCollection>(ps.getParameter<edm::InputTag>("EGammaCollection2"))),
11  TauCollection1(consumes <l1t::TauBxCollection>(ps.getParameter<edm::InputTag>("TauCollection1"))),
12  TauCollection2(consumes <l1t::TauBxCollection>(ps.getParameter<edm::InputTag>("TauCollection2"))),
13  EtSumCollection1(consumes <l1t::EtSumBxCollection>(ps.getParameter<edm::InputTag>("EtSumCollection1"))),
14  EtSumCollection2(consumes <l1t::EtSumBxCollection>(ps.getParameter<edm::InputTag>("EtSumCollection2"))),
15  verbose(ps.getUntrackedParameter<bool> ("verbose", false))
16 {}
T getParameter(std::string const &) const
T getUntrackedParameter(std::string const &, T const &) const
edm::EDGetTokenT< l1t::EtSumBxCollection > EtSumCollection1
edm::EDGetTokenT< l1t::EGammaBxCollection > EGammaCollection1
edm::EDGetTokenT< l1t::TauBxCollection > TauCollection1
edm::EDGetTokenT< l1t::TauBxCollection > TauCollection2
edm::EDGetTokenT< l1t::JetBxCollection > JetCollection2
edm::EDGetTokenT< l1t::EGammaBxCollection > EGammaCollection2
edm::EDGetTokenT< l1t::EtSumBxCollection > EtSumCollection2
edm::EDGetTokenT< l1t::JetBxCollection > JetCollection1

Member Function Documentation

void L1TStage2uGTCaloLayer2Comp::analyze ( const edm::Event e,
const edm::EventSetup c 
)
overrideprotected

Main method where the analysis code resides, executed once for each run

The main body of the module code is contained in this method. The different object collections are extracted from the run and are passed to the respective comparison methods for processing of each object type.

Parameters
edm::Eventconst & Reference to event object
edm::EventSetupconst & Reference to event configuration object
Returns
void

Definition at line 70 of file L1TStage2uGTCaloLayer2Comp.cc.

References compareEGs(), compareJets(), compareSums(), compareTaus(), comparisonDenum, comparisonNum, EGammaCollection1, EGammaCollection2, EtSumCollection1, EtSumCollection2, EVENTBAD, EVENTS1, EVENTS2, EVENTS3, EVENTS4, EVENTS5, MonitorElement::Fill(), edm::Event::getByToken(), JetCollection1, JetCollection2, TauCollection1, and TauCollection2.

72  {
73 
74  // define collections to hold lists of objects in event
83 
84  // map event contents to above collections
85  e.getByToken(JetCollection1, jetCol1);
86  e.getByToken(JetCollection2, jetCol2);
87  e.getByToken(EGammaCollection1, egCol1);
88  e.getByToken(EGammaCollection2, egCol2);
89  e.getByToken(TauCollection1, tauCol1);
90  e.getByToken(TauCollection2, tauCol2);
91  e.getByToken(EtSumCollection1, sumCol1);
92  e.getByToken(EtSumCollection2, sumCol2);
93 
94  bool eventGood = true;
95 
96  if (!compareJets(jetCol1, jetCol2)) {
97  eventGood = false;
98  }
99 
100  if (!compareEGs(egCol1, egCol2)) {
101  eventGood = false;
102  }
103 
104  if (!compareTaus(tauCol1, tauCol2)) {
105  eventGood = false;
106  }
107 
108  if (!compareSums(sumCol1, sumCol2)) {
109  eventGood = false;
110  }
111 
112  if (!eventGood) {
114  }
115 
121 }
bool compareTaus(const edm::Handle< l1t::TauBxCollection > &col1, const edm::Handle< l1t::TauBxCollection > &col2)
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: Event.h:517
edm::EDGetTokenT< l1t::EtSumBxCollection > EtSumCollection1
edm::EDGetTokenT< l1t::EGammaBxCollection > EGammaCollection1
void Fill(long long x)
bool compareSums(const edm::Handle< l1t::EtSumBxCollection > &col1, const edm::Handle< l1t::EtSumBxCollection > &col2)
edm::EDGetTokenT< l1t::TauBxCollection > TauCollection1
edm::EDGetTokenT< l1t::TauBxCollection > TauCollection2
edm::EDGetTokenT< l1t::JetBxCollection > JetCollection2
bool compareJets(const edm::Handle< l1t::JetBxCollection > &col1, const edm::Handle< l1t::JetBxCollection > &col2)
edm::EDGetTokenT< l1t::EGammaBxCollection > EGammaCollection2
edm::EDGetTokenT< l1t::EtSumBxCollection > EtSumCollection2
edm::EDGetTokenT< l1t::JetBxCollection > JetCollection1
bool compareEGs(const edm::Handle< l1t::EGammaBxCollection > &col1, const edm::Handle< l1t::EGammaBxCollection > &col2)
void L1TStage2uGTCaloLayer2Comp::bookHistograms ( DQMStore::IBooker ,
const edm::Run ,
const edm::EventSetup  
)
overrideprotected

Method to declare or "book" all histograms that will be part of module

Histograms that are to be visualised as part of the DQM module should be registered with the IBooker object any additional configuration such as title or axis labels and ranges. A good rule of thumb for the amount of configuration is that it should be possible to understnand the contents of the histogram using the configuration received from this method since the plots generated by this module would later be stored into ROOT files for transfer to the DQM system and it should be possible to extract useful information without the need for specific render plugins.

Parameters
DQMStore::IBooker&ibooker Object that handles the creation of plots
edm::Runconst & Reference to run object
edm::EventSetupconst & Reference to event configuration object
Returns
void

Definition at line 18 of file L1TStage2uGTCaloLayer2Comp.cc.

References BADSUM, DQMStore::IBooker::book1D(), collection1Title, collection2Title, comparisonDenum, comparisonNum, EGBADET, EGBADETA, EGBADPHI, EGS1, EGS2, EGS3, EVENTBAD, EVENTBADEGCOL, EVENTBADJETCOL, EVENTBADSUMCOL, EVENTBADTAUCOL, EVENTS1, EVENTS2, EVENTS3, EVENTS4, EVENTS5, MonitorElement::getTH1F(), JETBADET, JETBADETA, JETBADPHI, JETS1, JETS2, JETS3, monitorDir, MonitorElement::setBinLabel(), DQMStore::IBooker::setCurrentFolder(), SUMS, TAUBADET, TAUBADETA, TAUBADPHI, TAUS1, TAUS2, and TAUS3.

21  {
22 
23  ibooker.setCurrentFolder(monitorDir);
24 
25  // the index of the first bin in histogram should match value of first enum
26  comparisonNum = ibooker.book1D(
27  "errorSummaryNum",
28  collection1Title+" vs "+collection2Title+" Comparison - Numerator (# Disagreements)", 15, 1, 16);
29 
30  comparisonNum->setBinLabel(EVENTBAD, "# bad evts");
31  comparisonNum->setBinLabel(EVENTBADJETCOL, "# evts w/ bad jet col size");
32  comparisonNum->setBinLabel(EVENTBADEGCOL, "# evts w/ bad eg col size");
33  comparisonNum->setBinLabel(EVENTBADTAUCOL, "# evts w/ bad tau col size");
34  comparisonNum->setBinLabel(EVENTBADSUMCOL, "# evts w/ bad sum col size");
35  comparisonNum->setBinLabel(JETBADET, "# jets bad Et");
36  comparisonNum->setBinLabel(JETBADPHI, "# jets bad phi");
37  comparisonNum->setBinLabel(JETBADETA, "# jets bad eta");
38  comparisonNum->setBinLabel(EGBADET, "# egs bad Et");
39  comparisonNum->setBinLabel(EGBADPHI, "# egs bad phi");
40  comparisonNum->setBinLabel(EGBADETA, "# egs bad eta");
41  comparisonNum->setBinLabel(TAUBADET, "# taus bad Et");
42  comparisonNum->setBinLabel(TAUBADPHI, "# taus bad phi");
43  comparisonNum->setBinLabel(TAUBADETA, "# taus bad eta");
44  comparisonNum->setBinLabel(BADSUM, "# bad sums");
45 
46  comparisonDenum = ibooker.book1D(
47  "errorSummaryDen",
48  collection1Title+" vs "+collection2Title+" Comparison - Denominator (# Objects)", 15, 1, 16);
49 
55  comparisonDenum->setBinLabel(JETS1, "# jets");
56  comparisonDenum->setBinLabel(JETS2, "# jets");
57  comparisonDenum->setBinLabel(JETS3, "# jets");
58  comparisonDenum->setBinLabel(EGS1, "# egs");
59  comparisonDenum->setBinLabel(EGS2, "# egs");
60  comparisonDenum->setBinLabel(EGS3, "# egs");
61  comparisonDenum->setBinLabel(TAUS1, "# taus");
62  comparisonDenum->setBinLabel(TAUS2, "# taus");
63  comparisonDenum->setBinLabel(TAUS3, "# taus");
64  comparisonDenum->setBinLabel(SUMS, "# sums");
65  // Setting canExtend to false is needed to get the correct behaviour when running multithreaded.
66  // Otherwise, when merging the histgrams of the threads, TH1::Merge sums bins that have the same label in one bin.
67  // This needs to come after the calls to setBinLabel.
68  comparisonDenum->getTH1F()->GetXaxis()->SetCanExtend(false);
69 }
TH1F * getTH1F() const
void setBinLabel(int bin, const std::string &label, int axis=1)
set bin label for x, y or z axis (axis=1, 2, 3 respectively)
bool L1TStage2uGTCaloLayer2Comp::compareEGs ( const edm::Handle< l1t::EGammaBxCollection > &  col1,
const edm::Handle< l1t::EGammaBxCollection > &  col2 
)
private

Encapsulates the code required for performing a comparison of the e/gs contained in a given event.

Method is called once per each event with the e/g collections associated with the event being extracted for all bx. The implementation checks if the size of collections is the same and when so, compares the e/gs in the same positions within the calol2/ugt collections. The number and type of discrepancies are accumulated in different bins of a summary histogram.

Parameters
edm::Handle<l1t::EGammaBXCollection>&col1 Reference to e/gamma collection 1
edm::Handle<l1t::EGammaBXCollection>&col2 Reference to e/gamma collection 2
Returns
bool Flag of whether the agreement was perfect

Definition at line 190 of file L1TStage2uGTCaloLayer2Comp.cc.

References BXVector< T >::begin(), comparisonDenum, comparisonNum, EGBADET, EGBADETA, EGBADPHI, EGS1, EGS2, EGS3, BXVector< T >::end(), EVENTBADEGCOL, MonitorElement::Fill(), and BXVector< T >::size().

Referenced by analyze().

193 {
194  bool eventGood = true;
195 
198 
199  // check length of collections
200  if (col1->size() != col2->size()) {
202  return false;
203  }
204 
205  // processing continues only of length of object collections is the same
206  if (col1It != col1->end() ||
207  col2It != col2->end()) {
208 
209  while(true) {
210 
211  // object pt mismatch
212  if (col1It->hwPt() != col2It->hwPt()) {
214  eventGood = false;
215  }
216 
217  // object position mismatch (phi)
218  if (col1It->hwPhi() != col2It->hwPhi()) {
220  eventGood = false;
221  }
222 
223  // object position mismatch (eta)
224  if (col1It->hwEta() != col2It->hwEta()) {
226  eventGood = false;
227  }
228 
229  // keep track of number of objects
233 
234  // increment position of pointers
235  ++col1It;
236  ++col2It;
237 
238  if (col1It == col1->end() ||
239  col2It == col2->end())
240  break;
241  }
242  } else {
243  if (col1->size() != 0 || col2->size() != 0) {
245  return false;
246  }
247  }
248 
249  // return a boolean that states whether the eg data in the event is in
250  // agreement
251  return eventGood;
252 }
const_iterator end(int bx) const
unsigned size(int bx) const
void Fill(long long x)
const_iterator begin(int bx) const
std::vector< T >::const_iterator const_iterator
Definition: BXVector.h:20
bool L1TStage2uGTCaloLayer2Comp::compareJets ( const edm::Handle< l1t::JetBxCollection > &  col1,
const edm::Handle< l1t::JetBxCollection > &  col2 
)
private

Encapsulates the code required for performing a comparison of the jets contained in a given event.

Method is called once per each event with the jet collections associated with the event being extracted for all bx. The implementation checks if the size of collections is the same and when so, compares the jets in the same positions within the calol2/ugt collections. The number and type of discrepancies are accumulated in different bins of a summary histogram.

Parameters
edm::Handle<l1t::JetBXCollection>&col1 Reference to jet collection 1
edm::Handle<l1t::JetBXCollection>&col2 Reference to jet collection 2
Returns
bool Flag of whether the agreement was perfect

Definition at line 124 of file L1TStage2uGTCaloLayer2Comp.cc.

References BXVector< T >::begin(), comparisonDenum, comparisonNum, BXVector< T >::end(), EVENTBADJETCOL, MonitorElement::Fill(), JETBADET, JETBADETA, JETBADPHI, JETS1, JETS2, JETS3, and BXVector< T >::size().

Referenced by analyze().

127 {
128  bool eventGood = true;
129 
132 
133  // process jets
134  if (col1->size() != col2->size()) {
136  return false;
137  }
138 
139  int nJets = 0;
140  if (col1It != col1->end() ||
141  col2It != col2->end()) {
142  while(true) {
143 
144  ++nJets;
145 
146  // object pt mismatch
147  if (col1It->hwPt() != col2It->hwPt()) {
149  eventGood = false;
150  }
151 
152  // object position mismatch (phi)
153  if (col1It->hwPhi() != col2It->hwPhi()){
155  eventGood = false;
156  }
157 
158  // object position mismatch (eta)
159  if (col1It->hwEta() != col2It->hwEta()) {
161  eventGood = false;
162  }
163 
164  // keep track of jets
168 
169  // increment position of pointers
170  ++col1It;
171  ++col2It;
172 
173  if (col1It == col1->end() ||
174  col2It == col2->end())
175  break;
176  }
177  } else {
178  if (col1->size() != 0 || col2->size() != 0) {
180  return false;
181  }
182  }
183 
184  // return a boolean that states whether the jet data in the event is in
185  // agreement
186  return eventGood;
187 }
const_iterator end(int bx) const
unsigned size(int bx) const
void Fill(long long x)
const_iterator begin(int bx) const
std::vector< T >::const_iterator const_iterator
Definition: BXVector.h:20
bool L1TStage2uGTCaloLayer2Comp::compareSums ( const edm::Handle< l1t::EtSumBxCollection > &  col1,
const edm::Handle< l1t::EtSumBxCollection > &  col2 
)
private

Encapsulates the code required for performing a comparison of the taus contained in a given event.

Method is called once per each event with the sum collections associated with the event being extracted for all bx. The implementation loops over the collection and depending of the their type sums are compared separately but all sum errors are accumulated together.

Parameters
edm::Handle<l1t::TauBXCollection>&col1 Reference to sum collection 1
edm::Handle<l1t::TauBXCollection>&col2 Reference to sum collection 2
Returns
bool Flag of whether the agreement was perfect

Definition at line 319 of file L1TStage2uGTCaloLayer2Comp.cc.

References BADSUM, BXVector< T >::begin(), comparisonDenum, comparisonNum, BXVector< T >::end(), EVENTBADSUMCOL, MonitorElement::Fill(), L1Analysis::kMinBiasHFM0, L1Analysis::kMinBiasHFM1, L1Analysis::kMinBiasHFP0, L1Analysis::kMinBiasHFP1, L1Analysis::kMissingEt, L1Analysis::kMissingEtHF, L1Analysis::kMissingHt, L1Analysis::kMissingHtHF, L1Analysis::kTotalEt, L1Analysis::kTotalEtEm, L1Analysis::kTotalHt, L1Analysis::kTowerCount, BXVector< T >::size(), and SUMS.

Referenced by analyze().

322 {
323  bool eventGood = true;
324 
325  double col1Et = 0;
326  double col2Et = 0;
327  double col1Phi = 0;
328  double col2Phi = 0;
329 
330  // if the calol2 or ugt collections have different size, mark the event as
331  // bad (this should never occur in normal running)
332  if (col1->size() != col2->size()) {
334  return false;
335  }
336 
339 
340  while (col1It != col1->end() && col2It != col2->end()) {
341 
342  // ETT, ETTEM, HTT, TowCnt, MBHFP0, MBHFM0, MBHFP1 or MBHFM1
343  if ((l1t::EtSum::EtSumType::kTotalEt == col1It->getType()) || // ETT
344  (l1t::EtSum::EtSumType::kTotalEtEm == col1It->getType()) || // ETTEM
345  (l1t::EtSum::EtSumType::kTotalHt == col1It->getType()) || // HTT
346  (l1t::EtSum::EtSumType::kTowerCount == col1It->getType()) || // TowCnt
347  (l1t::EtSum::EtSumType::kMinBiasHFP0 == col1It->getType()) ||// MBHFP0
348  (l1t::EtSum::EtSumType::kMinBiasHFM0 == col1It->getType()) ||// MBHFM0
349  (l1t::EtSum::EtSumType::kMinBiasHFP1 == col1It->getType()) ||// MBHFP1
350  (l1t::EtSum::EtSumType::kMinBiasHFM1 == col1It->getType())) {// MBHFM1
351 
352  col1Et = col1It->hwPt();
353  col2Et = col2It->hwPt();
354 
355  if (col1Et != col2Et) {
356  eventGood = false;
358  }
359 
360  // update sum counters
362  }
363 
364  // MET, METHF, MHT or MHTHF
365  if ((l1t::EtSum::EtSumType::kMissingEt == col1It->getType()) || // MET
366  (l1t::EtSum::EtSumType::kMissingEtHF == col1It->getType()) || // METHF
367  (l1t::EtSum::EtSumType::kMissingHt == col1It->getType()) || // MHT
368  (l1t::EtSum::EtSumType::kMissingHtHF == col1It->getType())) { // MHTHF
369 
370  col1Et = col1It->hwPt();
371  col2Et = col2It->hwPt();
372 
373  col1Phi = col1It->hwPhi();
374  col2Phi = col2It->hwPhi();
375 
376  if ((col1Et != col2Et) || (col1Phi != col2Phi)) {
377  eventGood = false;
379  }
380 
381  // update sum counters
383  }
384 
385  ++col1It;
386  ++col2It;
387  }
388 
389  // return a boolean that states whether the sum data in the event is in
390  // agreement
391  return eventGood;
392 }
const_iterator end(int bx) const
unsigned size(int bx) const
void Fill(long long x)
const_iterator begin(int bx) const
std::vector< T >::const_iterator const_iterator
Definition: BXVector.h:20
bool L1TStage2uGTCaloLayer2Comp::compareTaus ( const edm::Handle< l1t::TauBxCollection > &  col1,
const edm::Handle< l1t::TauBxCollection > &  col2 
)
private

Encapsulates the code required for performing a comparison of the taus contained in a given event.

Method is called once per each event with the e/g collections associated with the event being extracted for all bx. The implementation checks if the size of collections is the same and when so, compares the taus in the same positions within the calol2/ugt collections. The number and type

Parameters
edm::Handle<l1t::TauBXCollection>&col1 Reference to tau collection 1
edm::Handle<l1t::TauBXCollection>&col2 Reference to tau collection 2
Returns
bool Flag of whether the agreement was perfect

Definition at line 255 of file L1TStage2uGTCaloLayer2Comp.cc.

References BXVector< T >::begin(), comparisonDenum, comparisonNum, BXVector< T >::end(), EVENTBADTAUCOL, MonitorElement::Fill(), BXVector< T >::size(), TAUBADET, TAUBADETA, TAUBADPHI, TAUS1, TAUS2, and TAUS3.

Referenced by analyze().

258 {
259  bool eventGood = true;
260 
263 
264  // check length of collections
265  if (col1->size() != col2->size()) {
267  return false;
268  }
269 
270  // processing continues only of length of object collections is the same
271  if (col1It != col1->end() ||
272  col2It != col2->end()) {
273 
274  while(true) {
275  // object Et mismatch
276  if (col1It->hwPt() != col2It->hwPt()) {
278  eventGood = false;
279  }
280 
281  // object position mismatch (phi)
282  if (col1It->hwPhi() != col2It->hwPhi()) {
284  eventGood = false;
285  }
286 
287  // object position mismatch (eta)
288  if (col1It->hwEta() != col2It->hwEta()) {
290  eventGood = false;
291  }
292 
293  // keep track of number of objects
297 
298  // increment position of pointers
299  ++col1It;
300  ++col2It;
301 
302  if (col1It == col1->end() ||
303  col2It == col2->end())
304  break;
305  }
306  } else {
307  if (col1->size() != 0 || col2->size() != 0) {
309  return false;
310  }
311  }
312 
313  // return a boolean that states whether the tau data in the event is in
314  // agreement
315  return eventGood;
316 }
const_iterator end(int bx) const
unsigned size(int bx) const
void Fill(long long x)
const_iterator begin(int bx) const
std::vector< T >::const_iterator const_iterator
Definition: BXVector.h:20

Member Data Documentation

std::string L1TStage2uGTCaloLayer2Comp::collection1Title
private

Definition at line 176 of file L1TStage2uGTCaloLayer2Comp.h.

Referenced by bookHistograms().

std::string L1TStage2uGTCaloLayer2Comp::collection2Title
private

Definition at line 177 of file L1TStage2uGTCaloLayer2Comp.h.

Referenced by bookHistograms().

MonitorElement* L1TStage2uGTCaloLayer2Comp::comparisonDenum
private
MonitorElement* L1TStage2uGTCaloLayer2Comp::comparisonNum
private
edm::EDGetTokenT<l1t::EGammaBxCollection> L1TStage2uGTCaloLayer2Comp::EGammaCollection1
private

Definition at line 182 of file L1TStage2uGTCaloLayer2Comp.h.

Referenced by analyze().

edm::EDGetTokenT<l1t::EGammaBxCollection> L1TStage2uGTCaloLayer2Comp::EGammaCollection2
private

Definition at line 183 of file L1TStage2uGTCaloLayer2Comp.h.

Referenced by analyze().

edm::EDGetTokenT<l1t::EtSumBxCollection> L1TStage2uGTCaloLayer2Comp::EtSumCollection1
private

Definition at line 186 of file L1TStage2uGTCaloLayer2Comp.h.

Referenced by analyze().

edm::EDGetTokenT<l1t::EtSumBxCollection> L1TStage2uGTCaloLayer2Comp::EtSumCollection2
private

Definition at line 187 of file L1TStage2uGTCaloLayer2Comp.h.

Referenced by analyze().

edm::EDGetTokenT<l1t::JetBxCollection> L1TStage2uGTCaloLayer2Comp::JetCollection1
private

Definition at line 180 of file L1TStage2uGTCaloLayer2Comp.h.

Referenced by analyze().

edm::EDGetTokenT<l1t::JetBxCollection> L1TStage2uGTCaloLayer2Comp::JetCollection2
private

Definition at line 181 of file L1TStage2uGTCaloLayer2Comp.h.

Referenced by analyze().

std::string L1TStage2uGTCaloLayer2Comp::monitorDir
private

Definition at line 173 of file L1TStage2uGTCaloLayer2Comp.h.

Referenced by bookHistograms().

edm::EDGetTokenT<l1t::TauBxCollection> L1TStage2uGTCaloLayer2Comp::TauCollection1
private

Definition at line 184 of file L1TStage2uGTCaloLayer2Comp.h.

Referenced by analyze().

edm::EDGetTokenT<l1t::TauBxCollection> L1TStage2uGTCaloLayer2Comp::TauCollection2
private

Definition at line 185 of file L1TStage2uGTCaloLayer2Comp.h.

Referenced by analyze().

bool L1TStage2uGTCaloLayer2Comp::verbose
private