File indexing completed on 2022-06-29 02:25:58
0001 #ifndef L1Trigger_L1TCalorimeter_L1TStage2CaloLayer2Comp
0002 #define L1Trigger_L1TCalorimeter_L1TStage2CaloLayer2Comp
0003
0004 #include "FWCore/Framework/interface/MakerMacros.h"
0005
0006 #include "FWCore/Framework/interface/Event.h"
0007 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0008 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0009
0010 #include "DataFormats/L1Trigger/interface/EGamma.h"
0011 #include "DataFormats/L1Trigger/interface/Jet.h"
0012 #include "DataFormats/L1Trigger/interface/EtSum.h"
0013 #include "DataFormats/L1Trigger/interface/Tau.h"
0014
0015 #include "FWCore/Framework/interface/Event.h"
0016 #include "FWCore/Framework/interface/MakerMacros.h"
0017 #include "FWCore/Framework/interface/stream/EDProducer.h"
0018
0019 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0020
0021 #include "CondFormats/L1TObjects/interface/CaloParams.h"
0022 #include "CondFormats/DataRecord/interface/L1TCaloParamsRcd.h"
0023
0024 #include "DataFormats/L1TCalorimeter/interface/CaloTower.h"
0025 #include "DataFormats/L1TCalorimeter/interface/CaloCluster.h"
0026 #include "DataFormats/L1Trigger/interface/EGamma.h"
0027 #include "DataFormats/L1Trigger/interface/Tau.h"
0028 #include "DataFormats/L1Trigger/interface/Jet.h"
0029 #include "DataFormats/L1Trigger/interface/EtSum.h"
0030
0031
0032 #include "L1Trigger/L1TCalorimeter/interface/CaloTools.h"
0033
0034
0035 #include "L1Trigger/L1TCalorimeter/interface/Stage2Layer2JetAlgorithmFirmware.h"
0036 #include "L1Trigger/L1TCalorimeter/interface/AccumulatingSort.h"
0037 #include "L1Trigger/L1TCalorimeter/interface/BitonicSort.h"
0038 #include "DataFormats/Math/interface/LorentzVector.h"
0039
0040 namespace l1t {
0041 bool operator>(const l1t::Jet &a, l1t::Jet &b) { return a.hwPt() > b.hwPt(); }
0042 bool operator>(const l1t::EGamma &a, l1t::EGamma &b) { return a.hwPt() > b.hwPt(); }
0043 bool operator>(const l1t::Tau &a, l1t::Tau &b) { return a.hwPt() > b.hwPt(); }
0044 }
0045
0046
0047
0048
0049
0050
0051
0052 class L1TStage2CaloLayer2Comp : public edm::stream::EDProducer<> {
0053 public:
0054
0055
0056
0057
0058
0059
0060
0061
0062
0063
0064
0065
0066
0067 L1TStage2CaloLayer2Comp(const edm::ParameterSet &ps);
0068
0069 protected:
0070
0071
0072
0073
0074
0075
0076
0077
0078
0079
0080
0081
0082 void produce(edm::Event &, const edm::EventSetup &) override;
0083
0084 private:
0085
0086
0087
0088
0089
0090
0091
0092
0093
0094
0095
0096
0097
0098
0099
0100
0101
0102
0103
0104
0105
0106
0107 bool compareJets(const edm::Handle<l1t::JetBxCollection> &dataCol, const edm::Handle<l1t::JetBxCollection> &emulCol);
0108
0109
0110
0111
0112
0113
0114
0115
0116
0117
0118
0119
0120
0121
0122
0123
0124
0125
0126
0127
0128
0129
0130
0131
0132
0133 bool compareEGs(const edm::Handle<l1t::EGammaBxCollection> &dataCol,
0134 const edm::Handle<l1t::EGammaBxCollection> &emulCol);
0135
0136
0137
0138
0139
0140
0141
0142
0143
0144
0145
0146
0147
0148
0149
0150
0151
0152
0153
0154
0155
0156
0157
0158
0159
0160 bool compareTaus(const edm::Handle<l1t::TauBxCollection> &dataCol, const edm::Handle<l1t::TauBxCollection> &emulCol);
0161
0162
0163
0164
0165
0166
0167
0168
0169
0170
0171
0172
0173
0174
0175
0176
0177
0178
0179
0180
0181 bool compareSums(const edm::Handle<l1t::EtSumBxCollection> &dataCol,
0182 const edm::Handle<l1t::EtSumBxCollection> &emulCol);
0183
0184 void accuSort(std::vector<l1t::Jet> *jets);
0185 void accuSort(std::vector<l1t::Tau> *taus);
0186 void accuSort(std::vector<l1t::EGamma> *egs);
0187 void accuSort(std::vector<l1t::L1Candidate> &jets);
0188
0189 void dumpEventToFile();
0190 void dumpEventToEDM(edm::Event &e);
0191
0192
0193 edm::EDGetToken calol2JetCollectionData;
0194 edm::EDGetToken calol2JetCollectionEmul;
0195 edm::EDGetToken calol2EGammaCollectionData;
0196 edm::EDGetToken calol2EGammaCollectionEmul;
0197 edm::EDGetToken calol2TauCollectionData;
0198 edm::EDGetToken calol2TauCollectionEmul;
0199 edm::EDGetToken calol2EtSumCollectionData;
0200 edm::EDGetToken calol2EtSumCollectionEmul;
0201 edm::EDGetToken calol2CaloTowerCollectionData;
0202 edm::EDGetToken calol2CaloTowerCollectionEmul;
0203
0204
0205 edm::Handle<l1t::JetBxCollection> jetDataCol;
0206 edm::Handle<l1t::JetBxCollection> jetEmulCol;
0207 edm::Handle<l1t::EGammaBxCollection> egDataCol;
0208 edm::Handle<l1t::EGammaBxCollection> egEmulCol;
0209 edm::Handle<l1t::TauBxCollection> tauDataCol;
0210 edm::Handle<l1t::TauBxCollection> tauEmulCol;
0211 edm::Handle<l1t::EtSumBxCollection> sumDataCol;
0212 edm::Handle<l1t::EtSumBxCollection> sumEmulCol;
0213 edm::Handle<l1t::CaloTowerBxCollection> caloTowerDataCol;
0214 edm::Handle<l1t::CaloTowerBxCollection> caloTowerEmulCol;
0215
0216 const unsigned int currBx = 0;
0217
0218 bool dumpTowers = false;
0219 bool dumpWholeEvent = false;
0220 };
0221
0222 L1TStage2CaloLayer2Comp::L1TStage2CaloLayer2Comp(const edm::ParameterSet &ps)
0223 : calol2JetCollectionData(
0224 consumes<l1t::JetBxCollection>(ps.getParameter<edm::InputTag>("calol2JetCollectionData"))),
0225 calol2JetCollectionEmul(
0226 consumes<l1t::JetBxCollection>(ps.getParameter<edm::InputTag>("calol2JetCollectionEmul"))),
0227 calol2EGammaCollectionData(
0228 consumes<l1t::EGammaBxCollection>(ps.getParameter<edm::InputTag>("calol2EGammaCollectionData"))),
0229 calol2EGammaCollectionEmul(
0230 consumes<l1t::EGammaBxCollection>(ps.getParameter<edm::InputTag>("calol2EGammaCollectionEmul"))),
0231 calol2TauCollectionData(
0232 consumes<l1t::TauBxCollection>(ps.getParameter<edm::InputTag>("calol2TauCollectionData"))),
0233 calol2TauCollectionEmul(
0234 consumes<l1t::TauBxCollection>(ps.getParameter<edm::InputTag>("calol2TauCollectionEmul"))),
0235 calol2EtSumCollectionData(
0236 consumes<l1t::EtSumBxCollection>(ps.getParameter<edm::InputTag>("calol2EtSumCollectionData"))),
0237 calol2EtSumCollectionEmul(
0238 consumes<l1t::EtSumBxCollection>(ps.getParameter<edm::InputTag>("calol2EtSumCollectionEmul"))),
0239 calol2CaloTowerCollectionData(
0240 consumes<l1t::CaloTowerBxCollection>(ps.getParameter<edm::InputTag>("calol2CaloTowerCollectionData"))),
0241 calol2CaloTowerCollectionEmul(
0242 consumes<l1t::CaloTowerBxCollection>(ps.getParameter<edm::InputTag>("calol2CaloTowerCollectionEmul"))) {
0243 produces<l1t::JetBxCollection>("dataJet").setBranchAlias("dataJets");
0244 produces<l1t::JetBxCollection>("emulJet").setBranchAlias("emulJets");
0245 produces<l1t::EGammaBxCollection>("dataEg").setBranchAlias("dataEgs");
0246 produces<l1t::EGammaBxCollection>("emulEg").setBranchAlias("emulEgs");
0247 produces<l1t::TauBxCollection>("dataTau").setBranchAlias("dataTaus");
0248 produces<l1t::TauBxCollection>("emulTau").setBranchAlias("emulTaus");
0249 produces<l1t::EtSumBxCollection>("dataEtSum").setBranchAlias("dataEtSums");
0250 produces<l1t::EtSumBxCollection>("emulEtSum").setBranchAlias("emulEtSums");
0251 produces<l1t::CaloTowerBxCollection>("dataCaloTower").setBranchAlias("dataCaloTowers");
0252 produces<l1t::CaloTowerBxCollection>("emulCaloTower").setBranchAlias("emulCaloTowers");
0253
0254 dumpTowers = ps.getParameter<bool>("dumpTowers");
0255 dumpWholeEvent = ps.getParameter<bool>("dumpWholeEvent");
0256 }
0257
0258 void L1TStage2CaloLayer2Comp::produce(edm::Event &e, const edm::EventSetup &c) {
0259 bool eventGood = true;
0260
0261
0262 e.getByToken(calol2JetCollectionData, jetDataCol);
0263 e.getByToken(calol2JetCollectionEmul, jetEmulCol);
0264 e.getByToken(calol2EGammaCollectionData, egDataCol);
0265 e.getByToken(calol2EGammaCollectionEmul, egEmulCol);
0266 e.getByToken(calol2TauCollectionData, tauDataCol);
0267 e.getByToken(calol2TauCollectionEmul, tauEmulCol);
0268 e.getByToken(calol2EtSumCollectionData, sumDataCol);
0269 e.getByToken(calol2EtSumCollectionEmul, sumEmulCol);
0270 e.getByToken(calol2CaloTowerCollectionData, caloTowerDataCol);
0271 e.getByToken(calol2CaloTowerCollectionEmul, caloTowerEmulCol);
0272
0273 edm::LogProblem("l1tcalol2ebec") << "Processing event " << e.id() << std::endl;
0274
0275 if (caloTowerDataCol->size() == 0) {
0276 edm::LogProblem("l1tcalol2ebec") << "Empty caloTowers. Skipping event " << std::endl;
0277 return;
0278 }
0279
0280 if (!compareJets(jetDataCol, jetEmulCol)) {
0281 eventGood = false;
0282 }
0283
0284 if (!compareEGs(egDataCol, egEmulCol)) {
0285 eventGood = false;
0286 }
0287
0288 if (!compareTaus(tauDataCol, tauEmulCol)) {
0289 eventGood = false;
0290 }
0291
0292 if (!compareSums(sumDataCol, sumEmulCol)) {
0293 eventGood = false;
0294 }
0295
0296 if (!eventGood) {
0297 if (dumpWholeEvent) {
0298
0299 dumpEventToFile();
0300 }
0301
0302
0303 dumpEventToEDM(e);
0304
0305 edm::LogProblem("l1tcalol2ebec") << "Bad event! " << std::endl;
0306 }
0307 }
0308
0309
0310 bool L1TStage2CaloLayer2Comp::compareJets(const edm::Handle<l1t::JetBxCollection> &dataCol,
0311 const edm::Handle<l1t::JetBxCollection> &emulCol) {
0312 bool eventGood = true;
0313
0314 std::vector<l1t::Jet> *jets = new std::vector<l1t::Jet>;
0315 std::vector<l1t::Jet>::iterator dataIt;
0316 l1t::JetBxCollection::const_iterator dataBxIt = dataCol->begin(currBx);
0317 l1t::JetBxCollection::const_iterator emulBxIt = emulCol->begin(currBx);
0318
0319 if (dataCol->size(currBx) > 0) {
0320
0321 while (true) {
0322 jets->emplace_back(*dataBxIt);
0323
0324 dataBxIt++;
0325 if (dataBxIt == dataCol->end(currBx))
0326 break;
0327 }
0328
0329 accuSort(jets);
0330 dataIt = jets->begin();
0331 }
0332
0333
0334 if (jets->size() != emulCol->size(currBx)) {
0335 edm::LogProblem("l1tcalol2ebec") << "Jet collection size difference: "
0336 << "data = " << jets->size() << " emul = " << emulCol->size(currBx) << std::endl;
0337 return false;
0338 }
0339
0340 if (dataIt != jets->end() || emulBxIt != emulCol->end(currBx)) {
0341 int nPos = 0;
0342 int nNeg = 0;
0343
0344 while (true) {
0345 bool posGood = true;
0346 bool etGood = true;
0347
0348
0349 if (dataIt->hwPt() != emulBxIt->hwPt()) {
0350 etGood = false;
0351 eventGood = false;
0352 }
0353
0354
0355 if (dataIt->hwPhi() != emulBxIt->hwPhi()) {
0356 posGood = false;
0357 }
0358
0359
0360 if (dataIt->hwEta() != emulBxIt->hwEta()) {
0361 posGood = false;
0362 }
0363
0364
0365 if (etGood && !posGood) {
0366 l1t::JetBxCollection::const_iterator emulItCheckSort;
0367 l1t::JetBxCollection::const_iterator dataItCheckSort;
0368 for (emulItCheckSort = emulCol->begin(currBx); emulItCheckSort != emulCol->end(currBx); ++emulItCheckSort) {
0369 for (dataItCheckSort = jets->begin(); dataItCheckSort != jets->end(); ++dataItCheckSort) {
0370 if (dataItCheckSort->hwEta() > 0)
0371 ++nPos;
0372 if (dataItCheckSort->hwEta() < 0)
0373 ++nNeg;
0374
0375 if (dataIt->hwPt() == emulItCheckSort->hwPt() && dataIt->hwPhi() == emulItCheckSort->hwPhi() &&
0376 dataIt->hwEta() == emulItCheckSort->hwEta())
0377 posGood = true;
0378 }
0379 }
0380 }
0381
0382 if (etGood && dataIt->hwEta() > 0 && ((distance(dataIt, jets->end()) - nNeg) < 5))
0383 posGood = true;
0384 if (etGood && dataIt->hwEta() < 0 && (distance(dataIt, jets->end()) < 5))
0385 posGood = true;
0386
0387 if (!posGood)
0388 eventGood = false;
0389
0390
0391 if (!etGood || !posGood) {
0392 edm::LogProblem("l1tcalol2ebec") << "Jet Problem (data emul): "
0393 << "\tEt = " << dataIt->hwPt() << " " << emulBxIt->hwPt()
0394 << "\teta = " << dataIt->hwEta() << " " << emulBxIt->hwEta()
0395 << "\tphi = " << dataIt->hwPhi() << " " << emulBxIt->hwPhi() << std::endl;
0396 }
0397
0398
0399 ++dataIt;
0400 ++emulBxIt;
0401
0402 if (dataIt == jets->end() || emulBxIt == emulCol->end(currBx))
0403 break;
0404 }
0405 } else {
0406 if (!jets->empty() || emulCol->size(currBx) != 0)
0407 return false;
0408 }
0409
0410
0411
0412 return eventGood;
0413 }
0414
0415
0416 bool L1TStage2CaloLayer2Comp::compareEGs(const edm::Handle<l1t::EGammaBxCollection> &dataCol,
0417 const edm::Handle<l1t::EGammaBxCollection> &emulCol) {
0418 bool eventGood = true;
0419
0420 std::vector<l1t::EGamma> *egs = new std::vector<l1t::EGamma>;
0421 std::vector<l1t::EGamma>::iterator dataIt;
0422 l1t::EGammaBxCollection::const_iterator dataBxIt = dataCol->begin(currBx);
0423 l1t::EGammaBxCollection::const_iterator emulBxIt = emulCol->begin(currBx);
0424
0425 if (dataCol->size(currBx) > 0) {
0426
0427 while (true) {
0428 egs->emplace_back(*dataBxIt);
0429
0430 dataBxIt++;
0431 if (dataBxIt == dataCol->end(currBx))
0432 break;
0433 }
0434 accuSort(egs);
0435 dataIt = egs->begin();
0436 }
0437
0438
0439 if (egs->size() != emulCol->size(currBx)) {
0440 edm::LogProblem("l1tcalol2ebec") << "EG collection size difference: "
0441 << "data = " << egs->size() << " emul = " << emulCol->size(currBx) << std::endl;
0442 return false;
0443 }
0444
0445
0446 if (dataIt != egs->end() || emulBxIt != emulCol->end(currBx)) {
0447 int nPos = 0;
0448 int nNeg = 0;
0449
0450 while (true) {
0451 bool posGood = true;
0452 bool etGood = true;
0453
0454
0455 if (dataIt->hwPt() != emulBxIt->hwPt()) {
0456 etGood = false;
0457 eventGood = false;
0458 }
0459
0460
0461 if (dataIt->hwPhi() != emulBxIt->hwPhi()) {
0462 posGood = false;
0463 }
0464
0465
0466 if (dataIt->hwEta() != emulBxIt->hwEta()) {
0467 posGood = false;
0468 }
0469
0470
0471 if (etGood && !posGood) {
0472 l1t::EGammaBxCollection::const_iterator emulItCheckSort;
0473 l1t::EGammaBxCollection::const_iterator dataItCheckSort;
0474 for (emulItCheckSort = emulCol->begin(currBx); emulItCheckSort != emulCol->end(currBx); ++emulItCheckSort) {
0475 for (dataItCheckSort = egs->begin(); dataItCheckSort != egs->end(); ++dataItCheckSort) {
0476 if (dataItCheckSort->hwEta() > 0)
0477 ++nPos;
0478 if (dataItCheckSort->hwEta() < 0)
0479 ++nNeg;
0480
0481 if (dataIt->hwPt() == emulItCheckSort->hwPt() && dataIt->hwPhi() == emulItCheckSort->hwPhi() &&
0482 dataIt->hwEta() == emulItCheckSort->hwEta())
0483 posGood = true;
0484 }
0485 }
0486 }
0487
0488 if (etGood && dataIt->hwEta() > 0 && ((distance(dataIt, egs->end()) - nNeg) < 5))
0489 posGood = true;
0490 if (etGood && dataIt->hwEta() < 0 && (distance(dataIt, egs->end()) < 5))
0491 posGood = true;
0492
0493 if (!posGood)
0494 eventGood = false;
0495
0496
0497 if (!posGood || !etGood) {
0498 edm::LogProblem("l1tcalol2ebec") << "EG Problem (data emul): "
0499 << "\tEt = " << dataIt->hwPt() << " " << emulBxIt->hwPt()
0500 << "\teta = " << dataIt->hwEta() << " " << emulBxIt->hwEta()
0501 << "\tphi = " << dataIt->hwPhi() << " " << emulBxIt->hwPhi() << std::endl;
0502 }
0503
0504
0505 ++dataIt;
0506 ++emulBxIt;
0507
0508 if (dataIt == egs->end() || emulBxIt == emulCol->end(currBx))
0509 break;
0510 }
0511 } else {
0512 if (!egs->empty() || emulCol->size(currBx) != 0)
0513 return false;
0514 }
0515
0516
0517
0518 return eventGood;
0519 }
0520
0521
0522 bool L1TStage2CaloLayer2Comp::compareTaus(const edm::Handle<l1t::TauBxCollection> &dataCol,
0523 const edm::Handle<l1t::TauBxCollection> &emulCol) {
0524 bool eventGood = true;
0525
0526 std::vector<l1t::Tau> *taus = new std::vector<l1t::Tau>;
0527 std::vector<l1t::Tau>::iterator dataIt;
0528 l1t::TauBxCollection::const_iterator dataBxIt = dataCol->begin(currBx);
0529 l1t::TauBxCollection::const_iterator emulBxIt = emulCol->begin(currBx);
0530
0531 if (dataCol->size(currBx) > 0) {
0532
0533 while (true) {
0534 taus->emplace_back(*dataBxIt);
0535
0536 dataBxIt++;
0537 if (dataBxIt == dataCol->end(currBx))
0538 break;
0539 }
0540 accuSort(taus);
0541 dataIt = taus->begin();
0542 }
0543
0544
0545 if (taus->size() != emulCol->size(currBx)) {
0546 edm::LogProblem("l1tcalol2ebec") << "Tau collection size difference: "
0547 << "data = " << taus->size() << " emul = " << emulCol->size(currBx) << std::endl;
0548 return false;
0549 }
0550
0551
0552 if (dataIt != taus->end() || emulBxIt != emulCol->end(currBx)) {
0553 int nPos = 0;
0554 int nNeg = 0;
0555
0556 while (true) {
0557 bool posGood = true;
0558 bool etGood = true;
0559
0560
0561 if (dataIt->hwPt() != emulBxIt->hwPt()) {
0562 etGood = false;
0563 eventGood = false;
0564 }
0565
0566
0567 if (dataIt->hwPhi() != emulBxIt->hwPhi()) {
0568 posGood = false;
0569 }
0570
0571
0572 if (dataIt->hwEta() != emulBxIt->hwEta()) {
0573 posGood = false;
0574 }
0575
0576
0577 if (etGood && !posGood) {
0578 l1t::TauBxCollection::const_iterator emulItCheckSort;
0579 l1t::TauBxCollection::const_iterator dataItCheckSort;
0580 for (emulItCheckSort = emulCol->begin(currBx); emulItCheckSort != emulCol->end(currBx); ++emulItCheckSort) {
0581 for (dataItCheckSort = taus->begin(); dataItCheckSort != taus->end(); ++dataItCheckSort) {
0582 if (dataItCheckSort->hwEta() > 0)
0583 ++nPos;
0584 if (dataItCheckSort->hwEta() < 0)
0585 ++nNeg;
0586
0587 if (dataIt->hwPt() == emulItCheckSort->hwPt() && dataIt->hwPhi() == emulItCheckSort->hwPhi() &&
0588 dataIt->hwEta() == emulItCheckSort->hwEta())
0589 posGood = true;
0590 }
0591 }
0592 }
0593
0594 if (etGood && dataIt->hwEta() > 0 && ((distance(dataIt, taus->end()) - nNeg) < 5))
0595 posGood = true;
0596 if (etGood && dataIt->hwEta() < 0 && (distance(dataIt, taus->end()) < 5))
0597 posGood = true;
0598
0599 if (!posGood)
0600 eventGood = false;
0601
0602
0603 if (!posGood || !etGood) {
0604 edm::LogProblem("l1tcalol2ebec") << "Tau Problem (data emul): "
0605 << "\tEt = " << dataIt->hwPt() << " " << emulBxIt->hwPt()
0606 << "\teta = " << dataIt->hwEta() << " " << emulBxIt->hwEta()
0607 << "\tphi = " << dataIt->hwPhi() << " " << emulBxIt->hwPhi() << std::endl;
0608 }
0609
0610
0611 ++dataIt;
0612 ++emulBxIt;
0613
0614 if (dataIt == taus->end() || emulBxIt == emulCol->end(currBx))
0615 break;
0616 }
0617 } else {
0618 if (!taus->empty() || emulCol->size(currBx) != 0)
0619 return false;
0620 }
0621
0622
0623
0624 return eventGood;
0625 }
0626
0627
0628 bool L1TStage2CaloLayer2Comp::compareSums(const edm::Handle<l1t::EtSumBxCollection> &dataCol,
0629 const edm::Handle<l1t::EtSumBxCollection> &emulCol) {
0630 bool eventGood = true;
0631
0632 double dataEt = 0;
0633 double emulEt = 0;
0634
0635
0636
0637 if ((dataCol->isEmpty(currBx)) || (emulCol->isEmpty(currBx)) || (dataCol->size(currBx) != emulCol->size(currBx))) {
0638 edm::LogProblem("l1tcalol2ebec") << "EtSum collection size difference: "
0639 << "data = " << dataCol->size(currBx) << " emul = " << emulCol->size(currBx)
0640 << std::endl;
0641
0642 return false;
0643 }
0644
0645 for (unsigned int i = 0; i < emulCol->size(currBx); i++) {
0646 l1t::EtSum const &emulSum = emulCol->at(currBx, i);
0647 l1t::EtSum const &dataSum = dataCol->at(currBx, l1t::CaloTools::emul_to_data_sum_index_map[i]);
0648
0649 if (emulSum.getType() != dataSum.getType()) {
0650 edm::LogProblem("l1tcalol2ebec") << "EtSum type problem (data emul): " << dataSum.getType() << " "
0651 << emulSum.getType() << std::endl;
0652 }
0653
0654 dataEt = dataSum.hwPt();
0655 emulEt = emulSum.hwPt();
0656
0657 if (dataEt != emulEt) {
0658 eventGood = false;
0659 edm::LogProblem("l1tcalol2ebec") << "EtSum problem (data emul):\tType = " << emulSum.getType()
0660 << "\tEt = " << dataEt << " " << emulEt << std::endl;
0661 }
0662 }
0663
0664
0665
0666 return eventGood;
0667 }
0668
0669
0670 void L1TStage2CaloLayer2Comp::accuSort(std::vector<l1t::Jet> *jets) {
0671 math::PtEtaPhiMLorentzVector emptyP4;
0672 l1t::Jet tempJet(emptyP4, 0, 0, 0, 0);
0673 std::vector<std::vector<l1t::Jet> > jetEtaPos(41, std::vector<l1t::Jet>(18, tempJet));
0674 std::vector<std::vector<l1t::Jet> > jetEtaNeg(41, std::vector<l1t::Jet>(18, tempJet));
0675
0676 for (unsigned int iJet = 0; iJet < jets->size(); iJet++) {
0677 if (jets->at(iJet).hwEta() > 0)
0678 jetEtaPos.at(jets->at(iJet).hwEta() - 1).at((jets->at(iJet).hwPhi() - 1) / 4) = jets->at(iJet);
0679 else
0680 jetEtaNeg.at(-(jets->at(iJet).hwEta() + 1)).at((jets->at(iJet).hwPhi() - 1) / 4) = jets->at(iJet);
0681 }
0682
0683 AccumulatingSort<l1t::Jet> etaPosSorter(7);
0684 AccumulatingSort<l1t::Jet> etaNegSorter(7);
0685 std::vector<l1t::Jet> accumEtaPos;
0686 std::vector<l1t::Jet> accumEtaNeg;
0687
0688 for (int ieta = 0; ieta < 41; ++ieta) {
0689
0690 std::vector<l1t::Jet>::iterator start_, end_;
0691 start_ = jetEtaPos.at(ieta).begin();
0692 end_ = jetEtaPos.at(ieta).end();
0693 BitonicSort<l1t::Jet>(down, start_, end_);
0694 etaPosSorter.Merge(jetEtaPos.at(ieta), accumEtaPos);
0695
0696
0697 start_ = jetEtaNeg.at(ieta).begin();
0698 end_ = jetEtaNeg.at(ieta).end();
0699 BitonicSort<l1t::Jet>(down, start_, end_);
0700 etaNegSorter.Merge(jetEtaNeg.at(ieta), accumEtaNeg);
0701 }
0702
0703
0704 if (accumEtaPos.at(6).hwPt() == accumEtaPos.at(5).hwPt() && accumEtaPos.at(6).hwEta() == accumEtaPos.at(5).hwEta() &&
0705 accumEtaPos.at(6).hwPhi() > accumEtaPos.at(5).hwPhi()) {
0706 accumEtaPos.at(5) = accumEtaPos.at(6);
0707 }
0708 if (accumEtaNeg.at(6).hwPt() == accumEtaNeg.at(5).hwPt() && accumEtaNeg.at(6).hwEta() == accumEtaNeg.at(5).hwEta() &&
0709 accumEtaNeg.at(6).hwPhi() > accumEtaNeg.at(5).hwPhi()) {
0710 accumEtaNeg.at(5) = accumEtaNeg.at(6);
0711 }
0712
0713
0714 accumEtaPos.resize(6);
0715 accumEtaNeg.resize(6);
0716
0717
0718 jets->clear();
0719 for (const l1t::Jet &accjet : accumEtaPos) {
0720 if (accjet.hwPt() > 0)
0721 jets->push_back(accjet);
0722 }
0723 for (const l1t::Jet &accjet : accumEtaNeg) {
0724 if (accjet.hwPt() > 0)
0725 jets->push_back(accjet);
0726 }
0727 }
0728
0729
0730 void L1TStage2CaloLayer2Comp::accuSort(std::vector<l1t::EGamma> *egs) {
0731 math::PtEtaPhiMLorentzVector emptyP4;
0732 l1t::EGamma tempEGamma(emptyP4, 0, 0, 0, 0);
0733 std::vector<std::vector<l1t::EGamma> > jetEtaPos(41, std::vector<l1t::EGamma>(18, tempEGamma));
0734 std::vector<std::vector<l1t::EGamma> > jetEtaNeg(41, std::vector<l1t::EGamma>(18, tempEGamma));
0735
0736 for (unsigned int iEGamma = 0; iEGamma < egs->size(); iEGamma++) {
0737 if (egs->at(iEGamma).hwEta() > 0)
0738 jetEtaPos.at(egs->at(iEGamma).hwEta() - 1).at((egs->at(iEGamma).hwPhi() - 1) / 4) = egs->at(iEGamma);
0739 else
0740 jetEtaNeg.at(-(egs->at(iEGamma).hwEta() + 1)).at((egs->at(iEGamma).hwPhi() - 1) / 4) = egs->at(iEGamma);
0741 }
0742
0743 AccumulatingSort<l1t::EGamma> etaPosSorter(7);
0744 AccumulatingSort<l1t::EGamma> etaNegSorter(7);
0745 std::vector<l1t::EGamma> accumEtaPos;
0746 std::vector<l1t::EGamma> accumEtaNeg;
0747
0748 for (int ieta = 0; ieta < 41; ++ieta) {
0749
0750 std::vector<l1t::EGamma>::iterator start_, end_;
0751 start_ = jetEtaPos.at(ieta).begin();
0752 end_ = jetEtaPos.at(ieta).end();
0753 BitonicSort<l1t::EGamma>(down, start_, end_);
0754 etaPosSorter.Merge(jetEtaPos.at(ieta), accumEtaPos);
0755
0756
0757 start_ = jetEtaNeg.at(ieta).begin();
0758 end_ = jetEtaNeg.at(ieta).end();
0759 BitonicSort<l1t::EGamma>(down, start_, end_);
0760 etaNegSorter.Merge(jetEtaNeg.at(ieta), accumEtaNeg);
0761 }
0762
0763
0764 if (accumEtaPos.at(6).hwPt() == accumEtaPos.at(5).hwPt() && accumEtaPos.at(6).hwEta() == accumEtaPos.at(5).hwEta() &&
0765 accumEtaPos.at(6).hwPhi() > accumEtaPos.at(5).hwPhi()) {
0766 accumEtaPos.at(5) = accumEtaPos.at(6);
0767 }
0768 if (accumEtaNeg.at(6).hwPt() == accumEtaNeg.at(5).hwPt() && accumEtaNeg.at(6).hwEta() == accumEtaNeg.at(5).hwEta() &&
0769 accumEtaNeg.at(6).hwPhi() > accumEtaNeg.at(5).hwPhi()) {
0770 accumEtaNeg.at(5) = accumEtaNeg.at(6);
0771 }
0772
0773
0774 accumEtaPos.resize(6);
0775 accumEtaNeg.resize(6);
0776
0777
0778 egs->clear();
0779 for (const l1t::EGamma &accjet : accumEtaPos) {
0780 if (accjet.hwPt() > 0)
0781 egs->push_back(accjet);
0782 }
0783 for (const l1t::EGamma &accjet : accumEtaNeg) {
0784 if (accjet.hwPt() > 0)
0785 egs->push_back(accjet);
0786 }
0787 }
0788
0789
0790 void L1TStage2CaloLayer2Comp::accuSort(std::vector<l1t::Tau> *taus) {
0791 math::PtEtaPhiMLorentzVector emptyP4;
0792 l1t::Tau tempTau(emptyP4, 0, 0, 0, 0);
0793 std::vector<std::vector<l1t::Tau> > jetEtaPos(41, std::vector<l1t::Tau>(18, tempTau));
0794 std::vector<std::vector<l1t::Tau> > jetEtaNeg(41, std::vector<l1t::Tau>(18, tempTau));
0795
0796 for (unsigned int iTau = 0; iTau < taus->size(); iTau++) {
0797 if (taus->at(iTau).hwEta() > 0)
0798 jetEtaPos.at(taus->at(iTau).hwEta() - 1).at((taus->at(iTau).hwPhi() - 1) / 4) = taus->at(iTau);
0799 else
0800 jetEtaNeg.at(-(taus->at(iTau).hwEta() + 1)).at((taus->at(iTau).hwPhi() - 1) / 4) = taus->at(iTau);
0801 }
0802
0803 AccumulatingSort<l1t::Tau> etaPosSorter(7);
0804 AccumulatingSort<l1t::Tau> etaNegSorter(7);
0805 std::vector<l1t::Tau> accumEtaPos;
0806 std::vector<l1t::Tau> accumEtaNeg;
0807
0808 for (int ieta = 0; ieta < 41; ++ieta) {
0809
0810 std::vector<l1t::Tau>::iterator start_, end_;
0811 start_ = jetEtaPos.at(ieta).begin();
0812 end_ = jetEtaPos.at(ieta).end();
0813 BitonicSort<l1t::Tau>(down, start_, end_);
0814 etaPosSorter.Merge(jetEtaPos.at(ieta), accumEtaPos);
0815
0816
0817 start_ = jetEtaNeg.at(ieta).begin();
0818 end_ = jetEtaNeg.at(ieta).end();
0819 BitonicSort<l1t::Tau>(down, start_, end_);
0820 etaNegSorter.Merge(jetEtaNeg.at(ieta), accumEtaNeg);
0821 }
0822
0823
0824 if (accumEtaPos.at(6).hwPt() == accumEtaPos.at(5).hwPt() && accumEtaPos.at(6).hwEta() == accumEtaPos.at(5).hwEta() &&
0825 accumEtaPos.at(6).hwPhi() > accumEtaPos.at(5).hwPhi()) {
0826 accumEtaPos.at(5) = accumEtaPos.at(6);
0827 }
0828 if (accumEtaNeg.at(6).hwPt() == accumEtaNeg.at(5).hwPt() && accumEtaNeg.at(6).hwEta() == accumEtaNeg.at(5).hwEta() &&
0829 accumEtaNeg.at(6).hwPhi() > accumEtaNeg.at(5).hwPhi()) {
0830 accumEtaNeg.at(5) = accumEtaNeg.at(6);
0831 }
0832
0833
0834 accumEtaPos.resize(6);
0835 accumEtaNeg.resize(6);
0836
0837
0838 taus->clear();
0839 for (const l1t::Tau &accjet : accumEtaPos) {
0840 if (accjet.hwPt() > 0)
0841 taus->push_back(accjet);
0842 }
0843 for (const l1t::Tau &accjet : accumEtaNeg) {
0844 if (accjet.hwPt() > 0)
0845 taus->push_back(accjet);
0846 }
0847 }
0848
0849 void L1TStage2CaloLayer2Comp::dumpEventToFile() {
0850 edm::LogProblem("l1tcalol2ebec") << "==== Problems found, dumping full event contents ====" << std::endl;
0851
0852 edm::LogProblem("l1tcalol2ebec") << "==== Event contents in data: ====" << std::endl;
0853
0854 if (dumpTowers) {
0855 edm::LogProblem("l1tcalol2ebec") << "==== Towers: ====" << std::endl;
0856
0857 for (auto tower = caloTowerDataCol->begin(currBx); tower != caloTowerDataCol->end(currBx); ++tower)
0858 edm::LogProblem("l1tcalol2ebec") << "Tower: Et = " << tower->hwPt() << ", "
0859 << "eta = " << tower->hwEta() << ", "
0860 << "phi = " << tower->hwPhi() << std::endl;
0861 }
0862
0863 edm::LogProblem("l1tcalol2ebec") << "==== Jets: ====" << std::endl;
0864 for (auto jet = jetDataCol->begin(currBx); jet != jetDataCol->end(currBx); ++jet)
0865 edm::LogProblem("l1tcalol2ebec") << "Jet: Et = " << jet->hwPt() << ", "
0866 << "eta = " << jet->hwEta() << ", "
0867 << "phi = " << jet->hwPhi() << std::endl;
0868
0869 edm::LogProblem("l1tcalol2ebec") << "==== EGs: ====" << std::endl;
0870 for (auto eg = egDataCol->begin(currBx); eg != egDataCol->end(currBx); ++eg)
0871 edm::LogProblem("l1tcalol2ebec") << "EG: Et = " << eg->hwPt() << ", "
0872 << "eta = " << eg->hwEta() << ", "
0873 << "phi = " << eg->hwPhi() << std::endl;
0874
0875 edm::LogProblem("l1tcalol2ebec") << "==== Taus: ====" << std::endl;
0876 for (auto tau = tauDataCol->begin(currBx); tau != tauDataCol->end(currBx); ++tau)
0877 edm::LogProblem("l1tcalol2ebec") << "Tau: Et = " << tau->hwPt() << ", "
0878 << "eta = " << tau->hwEta() << ", "
0879 << "phi = " << tau->hwPhi() << std::endl;
0880
0881 edm::LogProblem("l1tcalol2ebec") << "==== Sums: ====" << std::endl;
0882 for (auto sum = sumDataCol->begin(currBx); sum != sumDataCol->end(currBx); ++sum)
0883 edm::LogProblem("l1tcalol2ebec") << "Sum: type = " << sum->getType() << " "
0884 << "Et = " << sum->hwPt() << ", "
0885 << "eta = " << sum->hwEta() << ", "
0886 << "phi = " << sum->hwPhi() << std::endl;
0887
0888 edm::LogProblem("l1tcalol2ebec") << "==== Event contents in emul: ====" << std::endl;
0889
0890 if (dumpTowers) {
0891 edm::LogProblem("l1tcalol2ebec") << "==== Towers: ====" << std::endl;
0892
0893 for (auto tower = caloTowerEmulCol->begin(currBx); tower != caloTowerEmulCol->end(currBx); ++tower)
0894 edm::LogProblem("l1tcalol2ebec") << "Tower: Et = " << tower->hwPt() << ", "
0895 << "eta = " << tower->hwEta() << ", "
0896 << "phi = " << tower->hwPhi() << std::endl;
0897 }
0898
0899 edm::LogProblem("l1tcalol2ebec") << "==== Jets: ====" << std::endl;
0900 for (auto jet = jetEmulCol->begin(currBx); jet != jetEmulCol->end(currBx); ++jet)
0901 edm::LogProblem("l1tcalol2ebec") << "Jet: Et = " << jet->hwPt() << ", "
0902 << "eta = " << jet->hwEta() << ", "
0903 << "phi = " << jet->hwPhi() << std::endl;
0904
0905 edm::LogProblem("l1tcalol2ebec") << "==== EGs: ====" << std::endl;
0906 for (auto eg = egEmulCol->begin(currBx); eg != egEmulCol->end(currBx); ++eg)
0907 edm::LogProblem("l1tcalol2ebec") << "EG: Et = " << eg->hwPt() << ", "
0908 << "eta = " << eg->hwEta() << ", "
0909 << "phi = " << eg->hwPhi() << std::endl;
0910
0911 edm::LogProblem("l1tcalol2ebec") << "==== Taus: ====" << std::endl;
0912 for (auto tau = tauEmulCol->begin(currBx); tau != tauEmulCol->end(currBx); ++tau)
0913 edm::LogProblem("l1tcalol2ebec") << "Tau: Et = " << tau->hwPt() << ", "
0914 << "eta = " << tau->hwEta() << ", "
0915 << "phi = " << tau->hwPhi() << std::endl;
0916
0917 edm::LogProblem("l1tcalol2ebec") << "==== Sums: ====" << std::endl;
0918 for (auto sum = sumEmulCol->begin(currBx); sum != sumEmulCol->end(currBx); ++sum)
0919 edm::LogProblem("l1tcalol2ebec") << "Sum: type = " << sum->getType() << " "
0920 << "Et = " << sum->hwPt() << ", "
0921 << "eta = " << sum->hwEta() << ", "
0922 << "phi = " << sum->hwPhi() << std::endl;
0923 }
0924
0925 void L1TStage2CaloLayer2Comp::dumpEventToEDM(edm::Event &e) {
0926
0927 std::unique_ptr<l1t::JetBxCollection> mpjets_data(new l1t::JetBxCollection(0, currBx, currBx));
0928 std::unique_ptr<l1t::JetBxCollection> mpjets_emul(new l1t::JetBxCollection(0, currBx, currBx));
0929
0930 for (auto jet = jetDataCol->begin(currBx); jet != jetDataCol->end(currBx); ++jet)
0931 mpjets_data->push_back(0, (*jet));
0932 for (auto jet = jetEmulCol->begin(currBx); jet != jetEmulCol->end(currBx); ++jet)
0933 mpjets_emul->push_back(0, (*jet));
0934
0935 e.put(std::move(mpjets_data), "dataJet");
0936 e.put(std::move(mpjets_emul), "emulJet");
0937
0938
0939 std::unique_ptr<l1t::EGammaBxCollection> mpEGammas_data(new l1t::EGammaBxCollection(0, currBx, currBx));
0940 std::unique_ptr<l1t::EGammaBxCollection> mpEGammas_emul(new l1t::EGammaBxCollection(0, currBx, currBx));
0941
0942 for (auto eg = egDataCol->begin(currBx); eg != egDataCol->end(currBx); ++eg)
0943 mpEGammas_data->push_back(0, (*eg));
0944 for (auto eg = egEmulCol->begin(currBx); eg != egEmulCol->end(currBx); ++eg)
0945 mpEGammas_emul->push_back(0, (*eg));
0946
0947 e.put(std::move(mpEGammas_data), "dataEg");
0948 e.put(std::move(mpEGammas_emul), "emulEg");
0949
0950
0951 std::unique_ptr<l1t::TauBxCollection> mptaus_data(new l1t::TauBxCollection(0, currBx, currBx));
0952 std::unique_ptr<l1t::TauBxCollection> mptaus_emul(new l1t::TauBxCollection(0, currBx, currBx));
0953
0954 for (auto tau = tauDataCol->begin(currBx); tau != tauDataCol->end(currBx); ++tau)
0955 mptaus_data->push_back(0, (*tau));
0956 for (auto tau = tauEmulCol->begin(currBx); tau != tauEmulCol->end(currBx); ++tau)
0957 mptaus_emul->push_back(0, (*tau));
0958
0959 e.put(std::move(mptaus_data), "dataTau");
0960 e.put(std::move(mptaus_emul), "emulTau");
0961
0962
0963 std::unique_ptr<l1t::EtSumBxCollection> mpsums_data(new l1t::EtSumBxCollection(0, currBx, currBx));
0964 std::unique_ptr<l1t::EtSumBxCollection> mpsums_emul(new l1t::EtSumBxCollection(0, currBx, currBx));
0965
0966 for (auto sum = sumDataCol->begin(currBx); sum != sumDataCol->end(currBx); ++sum)
0967 mpsums_data->push_back(0, (*sum));
0968 for (auto sum = sumEmulCol->begin(currBx); sum != sumEmulCol->end(currBx); ++sum)
0969 mpsums_emul->push_back(0, (*sum));
0970
0971 e.put(std::move(mpsums_data), "dataEtSum");
0972 e.put(std::move(mpsums_emul), "emulEtSum");
0973
0974
0975 std::unique_ptr<l1t::CaloTowerBxCollection> mptowers_data(new l1t::CaloTowerBxCollection(0, currBx, currBx));
0976 std::unique_ptr<l1t::CaloTowerBxCollection> mptowers_emul(new l1t::CaloTowerBxCollection(0, currBx, currBx));
0977
0978 for (auto tower = caloTowerDataCol->begin(currBx); tower != caloTowerDataCol->end(currBx); ++tower)
0979 mptowers_data->push_back(0, (*tower));
0980 for (auto tower = caloTowerEmulCol->begin(currBx); tower != caloTowerEmulCol->end(currBx); ++tower)
0981 mptowers_emul->push_back(0, (*tower));
0982
0983 e.put(std::move(mptowers_data), "dataCaloTower");
0984 e.put(std::move(mptowers_emul), "emulCaloTower");
0985 }
0986
0987 DEFINE_FWK_MODULE(L1TStage2CaloLayer2Comp);
0988
0989 #endif