File indexing completed on 2021-10-06 02:54:13
0001
0002
0003
0004
0005
0006
0007
0008 #include "L1Trigger/L1TMuonOverlapPhase1/interface/Tools/PatternGenerator.h"
0009 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0010
0011 #include <boost/range/adaptor/reversed.hpp>
0012
0013 #include "TFile.h"
0014 #include "TDirectory.h"
0015
0016 PatternGenerator::PatternGenerator(const edm::ParameterSet& edmCfg,
0017 const OMTFConfiguration* omtfConfig,
0018 GoldenPatternVec<GoldenPatternWithStat>& gps)
0019 : PatternOptimizerBase(edmCfg, omtfConfig, gps), eventCntPerGp(gps.size(), 0) {
0020 edm::LogImportant("l1tOmtfEventPrint") << "constructing PatternGenerator, type: "
0021 << edmCfg.getParameter<string>("patternGenerator") << std::endl;
0022
0023 if (edmCfg.getParameter<string>("patternGenerator") == "patternGen")
0024 initPatternGen();
0025 }
0026
0027 PatternGenerator::~PatternGenerator() {}
0028
0029 void PatternGenerator::initPatternGen() {
0030
0031 unsigned int i = 0;
0032 for (auto& gp : goldenPatterns) {
0033 gp->setKeyNumber(i++);
0034
0035 if (gp->key().thePt == 0)
0036 continue;
0037
0038 gp->reset();
0039
0040 int statBinsCnt = 1024;
0041 gp->iniStatisitics(statBinsCnt, 1);
0042 }
0043
0044 edm::LogImportant("l1tOmtfEventPrint") << "PatternGenerator::initPatternGen():" << __LINE__
0045 << " goldenPatterns.size() " << goldenPatterns.size() << std::endl;
0046
0047
0048
0049
0050
0051
0052
0053
0054 for (auto& gp : goldenPatterns) {
0055 for (unsigned int iLayer = 0; iLayer < gp->getPdf().size(); ++iLayer) {
0056 for (unsigned int iRefLayer = 0; iRefLayer < gp->getPdf()[iLayer].size(); ++iRefLayer) {
0057 for (unsigned int iBin = 0; iBin < gp->getPdf()[iLayer][iRefLayer].size(); iBin++) {
0058 gp->pdfAllRef[iLayer][iRefLayer][iBin] = 1;
0059 }
0060 }
0061 }
0062 }
0063
0064
0065
0066
0067
0068
0069
0070
0071
0072
0073
0074
0075
0076
0077
0078
0079
0080
0081
0082
0083
0084
0085
0086
0087
0088
0089
0090
0091
0092
0093
0094
0095
0096
0097
0098 }
0099
0100 void PatternGenerator::updateStat() {
0101
0102 AlgoMuon* algoMuon = omtfCand.get();
0103 if (!algoMuon) {
0104 edm::LogImportant("l1tOmtfEventPrint") << ":" << __LINE__ << " algoMuon is null" << std::endl;
0105 throw runtime_error("algoMuon is null");
0106 }
0107
0108 double ptSim = simMuon->momentum().pt();
0109 int chargeSim = (abs(simMuon->type()) == 13) ? simMuon->type() / -13 : 0;
0110
0111 unsigned int exptPatNum = omtfConfig->getPatternNum(ptSim, chargeSim);
0112 GoldenPatternWithStat* exptCandGp = goldenPatterns.at(exptPatNum).get();
0113
0114 eventCntPerGp[exptPatNum]++;
0115
0116
0117
0118 int pdfMiddle = 1 << (omtfConfig->nPdfAddrBits() - 1);
0119
0120
0121 for (unsigned int iRefHit = 0; iRefHit < exptCandGp->getResults()[candProcIndx].size(); ++iRefHit) {
0122 auto& gpResult = exptCandGp->getResults()[candProcIndx][iRefHit];
0123
0124 unsigned int refLayer = gpResult.getRefLayer();
0125
0126 if (gpResult.getFiredLayerCnt() >= 3) {
0127 for (unsigned int iLayer = 0; iLayer < gpResult.getStubResults().size(); iLayer++) {
0128
0129
0130 bool fired = false;
0131 if (gpResult.getStubResults()[iLayer].getMuonStub()) {
0132 if (omtfConfig->isBendingLayer(iLayer)) {
0133 if (gpResult.getStubResults()[iLayer].getMuonStub()->qualityHw >= 4)
0134 fired = true;
0135 } else
0136 fired = true;
0137 }
0138
0139 if (fired) {
0140 int phiDist = gpResult.getStubResults()[iLayer].getPdfBin();
0141 phiDist += exptCandGp->meanDistPhiValue(iLayer, refLayer) - pdfMiddle;
0142
0143
0144
0145
0146
0147
0148
0149
0150
0151
0152 phiDist += exptCandGp->getStatistics()[iLayer][refLayer].size() / 2;
0153
0154
0155 if (phiDist > 0 && phiDist < (int)(exptCandGp->getStatistics()[iLayer][refLayer].size())) {
0156
0157
0158 exptCandGp->updateStat(iLayer, refLayer, phiDist, 0, 1);
0159 }
0160 } else {
0161 int phiDist = 0;
0162 exptCandGp->updateStat(iLayer, refLayer, phiDist, 0, 1);
0163 }
0164 }
0165 }
0166 }
0167 }
0168
0169 void PatternGenerator::observeEventEnd(const edm::Event& iEvent,
0170 std::unique_ptr<l1t::RegionalMuonCandBxCollection>& finalCandidates) {
0171 if (simMuon == nullptr || omtfCand->getGoldenPatern() == nullptr)
0172 return;
0173
0174 if (abs(simMuon->momentum().eta()) < 0.8 || abs(simMuon->momentum().eta()) > 1.24)
0175 return;
0176
0177 PatternOptimizerBase::observeEventEnd(iEvent, finalCandidates);
0178
0179 updateStat();
0180 }
0181
0182 void PatternGenerator::endJob() {
0183 if (edmCfg.getParameter<string>("patternGenerator") == "modifyClassProb")
0184 modifyClassProb(1);
0185 else if (edmCfg.getParameter<string>("patternGenerator") == "groupPatterns")
0186 groupPatterns();
0187 else if (edmCfg.getParameter<string>("patternGenerator") == "patternGen") {
0188 upadatePdfs();
0189 writeLayerStat = true;
0190 } else if (edmCfg.getParameter<string>("patternGenerator") == "patternGenFromStat") {
0191 std::string rootFileName = edmCfg.getParameter<edm::FileInPath>("patternsROOTFile").fullPath();
0192 edm::LogImportant("l1tOmtfEventPrint") << "PatternGenerator::endJob() rootFileName " << rootFileName << std::endl;
0193 TFile inFile(rootFileName.c_str());
0194 TDirectory* curDir = (TDirectory*)inFile.Get("layerStats");
0195
0196 ostringstream ostrName;
0197 for (auto& gp : goldenPatterns) {
0198 if (gp->key().thePt == 0)
0199 continue;
0200
0201 int statBinsCnt = 1024;
0202 gp->iniStatisitics(statBinsCnt, 1);
0203
0204 for (unsigned int iLayer = 0; iLayer < gp->getPdf().size(); ++iLayer) {
0205 for (unsigned int iRefLayer = 0; iRefLayer < gp->getPdf()[iLayer].size(); ++iRefLayer) {
0206 ostrName.str("");
0207 ostrName << "histLayerStat_PatNum_" << gp->key().theNumber << "_refLayer_" << iRefLayer << "_Layer_"
0208 << iLayer;
0209
0210 TH1I* histLayerStat = (TH1I*)curDir->Get(ostrName.str().c_str());
0211
0212 if (histLayerStat) {
0213 for (int iBin = 0; iBin < statBinsCnt; iBin++) {
0214 gp->updateStat(iLayer, iRefLayer, iBin, 0, histLayerStat->GetBinContent(iBin + 1));
0215 }
0216 } else {
0217 edm::LogImportant("l1tOmtfEventPrint")
0218 << "PatternGenerator::endJob() - reading histLayerStat: histogram not found " << ostrName.str()
0219 << std::endl;
0220 }
0221 }
0222 }
0223 }
0224
0225 TH1* simMuFoundByOmtfPt_fromFile = (TH1*)inFile.Get("simMuFoundByOmtfPt");
0226 for (unsigned int iGp = 0; iGp < eventCntPerGp.size(); iGp++) {
0227 eventCntPerGp[iGp] = simMuFoundByOmtfPt_fromFile->GetBinContent(simMuFoundByOmtfPt_fromFile->FindBin(iGp));
0228 edm::LogImportant("l1tOmtfEventPrint")
0229 << "PatternGenerator::endJob() - eventCntPerGp: iGp" << iGp << " - " << eventCntPerGp[iGp] << std::endl;
0230 }
0231
0232
0233 int group = 0;
0234 int indexInGroup = 0;
0235 for (auto& gp : goldenPatterns) {
0236 indexInGroup++;
0237 gp->key().setGroup(group);
0238 gp->key().setIndexInGroup(indexInGroup);
0239
0240
0241 edm::LogImportant("l1tOmtfEventPrint")
0242 << "setGroup(group): group " << group << " indexInGroup " << indexInGroup << std::endl;
0243
0244 if (gp->key().thePt <= 10 && indexInGroup == 2) {
0245 indexInGroup = 0;
0246 group++;
0247 }
0248
0249 if (gp->key().thePt > 10 && indexInGroup == 4) {
0250 indexInGroup = 0;
0251 group++;
0252 }
0253 }
0254
0255 upadatePdfs();
0256
0257 modifyClassProb(1);
0258
0259
0260
0261 reCalibratePt();
0262 this->writeLayerStat = true;
0263 }
0264
0265 PatternOptimizerBase::endJob();
0266 }
0267
0268 void PatternGenerator::upadatePdfs() {
0269
0270 for (auto& gp : goldenPatterns) {
0271 if (gp->key().thePt == 0)
0272 continue;
0273 for (unsigned int iLayer = 0; iLayer < gp->getPdf().size(); ++iLayer) {
0274 for (unsigned int iRefLayer = 0; iRefLayer < gp->getPdf()[iLayer].size(); ++iRefLayer) {
0275 if (gp->getDistPhiBitShift(iLayer, iRefLayer)) {
0276 throw runtime_error(
0277 string(__FUNCTION__) + ":" + to_string(__LINE__) +
0278 "gp->getDistPhiBitShift(iLayer, iRefLayer) != 0 - cannot change DistPhiBitShift then!!!!");
0279 }
0280
0281 if ((gp->key().thePt <= 10) && (iLayer == 1 || iLayer == 3 || iLayer == 5)) {
0282 gp->setDistPhiBitShift(1, iLayer, iRefLayer);
0283 } else
0284 gp->setDistPhiBitShift(0, iLayer, iRefLayer);
0285
0286
0287
0288
0289
0290
0291
0292
0293
0294
0295
0296
0297
0298
0299
0300
0301 }
0302 }
0303 }
0304
0305 double minHitCntThresh = 0.001;
0306
0307 for (auto& gp : goldenPatterns) {
0308 if (gp->key().thePt == 0)
0309 continue;
0310
0311 int minHitCnt = minHitCntThresh * eventCntPerGp[gp->key().number()];
0312 edm::LogImportant("l1tOmtfEventPrint")
0313 << "PatternGenerator::upadatePdfs() Calculating meanDistPhi " << gp->key() << " eventCnt "
0314 << eventCntPerGp[gp->key().number()] << " minHitCnt " << minHitCnt << std::endl;
0315
0316 for (unsigned int iLayer = 0; iLayer < gp->getPdf().size(); ++iLayer) {
0317 for (unsigned int iRefLayer = 0; iRefLayer < gp->getPdf()[iLayer].size(); ++iRefLayer) {
0318
0319 double meanDistPhi = 0;
0320 double count = 0;
0321 for (unsigned int iBin = 1; iBin < gp->getStatistics()[iLayer][iRefLayer].size(); iBin++) {
0322
0323 meanDistPhi += iBin * gp->getStatistics()[iLayer][iRefLayer][iBin][0];
0324 count += gp->getStatistics()[iLayer][iRefLayer][iBin][0];
0325 }
0326
0327 if (count != 0) {
0328 meanDistPhi /= count;
0329
0330 meanDistPhi -= (gp->getStatistics()[iLayer][iRefLayer].size() / 2);
0331
0332 if (count < minHitCnt)
0333 meanDistPhi = 0;
0334 else
0335 edm::LogImportant("l1tOmtfEventPrint")
0336 << __FUNCTION__ << ": " << __LINE__ << " " << gp->key() << " iLayer " << iLayer << " iRefLayer "
0337 << iRefLayer << " count " << count << " meanDistPhi " << meanDistPhi << endl;
0338 }
0339 gp->setMeanDistPhiValue(round(meanDistPhi), iLayer, iRefLayer);
0340 }
0341 }
0342 }
0343
0344 OMTFConfiguration::vector2D patternGroups = omtfConfig->getPatternGroups(goldenPatterns);
0345 edm::LogImportant("l1tOmtfEventPrint") << "patternGroups:" << std::endl;
0346 for (unsigned int iGroup = 0; iGroup < patternGroups.size(); iGroup++) {
0347 edm::LogImportant("l1tOmtfEventPrint") << "patternGroup " << std::setw(2) << iGroup << " ";
0348 for (unsigned int i = 0; i < patternGroups[iGroup].size(); i++) {
0349 edm::LogImportant("l1tOmtfEventPrint") << i << " patNum " << patternGroups[iGroup][i] << " ";
0350 }
0351 edm::LogImportant("l1tOmtfEventPrint") << std::endl;
0352 }
0353
0354
0355 for (unsigned int iLayer = 0; iLayer < goldenPatterns.at(0)->getPdf().size(); ++iLayer) {
0356 for (unsigned int iRefLayer = 0; iRefLayer < goldenPatterns.at(0)->getPdf()[iLayer].size(); ++iRefLayer) {
0357
0358
0359 {
0360 for (unsigned int iGroup = 0; iGroup < patternGroups.size(); iGroup++) {
0361 double meanDistPhi = 0;
0362 int mergedCnt = 0;
0363 for (unsigned int i = 0; i < patternGroups[iGroup].size(); i++) {
0364 auto gp = goldenPatterns.at(patternGroups[iGroup][i]).get();
0365 meanDistPhi += gp->meanDistPhiValue(iLayer, iRefLayer);
0366 if (gp->meanDistPhiValue(iLayer, iRefLayer) != 0)
0367 mergedCnt++;
0368 }
0369
0370 if (mergedCnt) {
0371 meanDistPhi /= mergedCnt;
0372
0373 for (unsigned int i = 0; i < patternGroups[iGroup].size(); i++) {
0374 auto gp = goldenPatterns.at(patternGroups[iGroup][i]).get();
0375 gp->setMeanDistPhiValue(round(meanDistPhi), iLayer, iRefLayer);
0376 edm::LogImportant("l1tOmtfEventPrint")
0377 << __FUNCTION__ << ": " << __LINE__ << " iGroup " << iGroup << " numInGroup " << i << " " << gp->key()
0378 << " iLayer " << iLayer << " iRefLayer " << iRefLayer << " meanDistPhi after averaging "
0379 << meanDistPhi << endl;
0380 }
0381 }
0382 }
0383 }
0384 }
0385 }
0386
0387
0388 for (auto& gp : goldenPatterns) {
0389 if (gp->key().thePt == 0)
0390 continue;
0391
0392
0393 int minHitCnt = 2 * minHitCntThresh * eventCntPerGp[gp->key().number()];
0394
0395 for (unsigned int iLayer = 0; iLayer < gp->getPdf().size(); ++iLayer) {
0396 for (unsigned int iRefLayer = 0; iRefLayer < gp->getPdf()[iLayer].size(); ++iRefLayer) {
0397 {
0398 double norm = 0;
0399 for (unsigned int iBin = 0; iBin < gp->getStatistics()[iLayer][iRefLayer].size();
0400 iBin++) {
0401 norm += gp->getStatistics()[iLayer][iRefLayer][iBin][0];
0402 }
0403
0404 int pdfMiddle = gp->getPdf()[iLayer][iRefLayer].size() / 2;
0405 int statBinGroupSize = 1 << gp->getDistPhiBitShift(iLayer, iRefLayer);
0406 for (unsigned int iBinPdf = 0; iBinPdf < gp->getPdf()[iLayer][iRefLayer].size(); iBinPdf++) {
0407 double pdfVal = 0;
0408 if (iBinPdf > 0) {
0409 int groupedBins = 0;
0410 for (int i = 0; i < statBinGroupSize; i++) {
0411 int iBinStat =
0412 statBinGroupSize * ((int)(iBinPdf)-pdfMiddle) + i + gp->meanDistPhiValue(iLayer, iRefLayer);
0413
0414 iBinStat += (gp->getStatistics()[iLayer][iRefLayer].size() / 2);
0415
0416 if (iBinStat >= 0 && iBinStat < (int)gp->getStatistics()[iLayer][iRefLayer].size()) {
0417 pdfVal += gp->getStatistics()[iLayer][iRefLayer][iBinStat][0];
0418 groupedBins++;
0419
0420 }
0421 }
0422 if (norm > minHitCnt) {
0423 pdfVal /= (norm * statBinGroupSize);
0424 } else
0425 pdfVal = 0;
0426
0427
0428
0429
0430 } else {
0431 int iBinStat = 0;
0432 if (norm > 0) {
0433 pdfVal = gp->getStatistics()[iLayer][iRefLayer][iBinStat][0] / norm;
0434 }
0435 edm::LogImportant("l1tOmtfEventPrint")
0436 << __FUNCTION__ << ": " << __LINE__ << " " << gp->key() << "calculating pdf: iLayer " << iLayer
0437 << " iRefLayer " << iRefLayer << " norm " << std::setw(5) << norm << " no hits cnt " << std::setw(5)
0438 << gp->getStatistics()[iLayer][iRefLayer][iBinStat][0] << " pdfVal " << pdfVal << endl;
0439 }
0440
0441 double minPdfValFactor = 1;
0442 const double minPlog = log(omtfConfig->minPdfVal() * minPdfValFactor);
0443 const double pdfMaxVal = omtfConfig->pdfMaxValue();
0444
0445 int digitisedVal = 0;
0446 if (pdfVal >= omtfConfig->minPdfVal() * minPdfValFactor) {
0447 digitisedVal = rint(pdfMaxVal - log(pdfVal) / minPlog * pdfMaxVal);
0448 }
0449
0450 gp->setPdfValue(digitisedVal, iLayer, iRefLayer, iBinPdf);
0451
0452 }
0453 }
0454 }
0455 }
0456 }
0457 }
0458
0459 void PatternGenerator::saveHists(TFile& outfile) {
0460 outfile.mkdir("ptDeltaPhiHists")->cd();
0461
0462
0463
0464
0465
0466
0467
0468
0469 }
0470
0471 void PatternGenerator::modifyClassProb(double step) {
0472 edm::LogImportant("l1tOmtfEventPrint") << __FUNCTION__ << ": " << __LINE__ << " Correcting P(C_k) " << std::endl;
0473 unsigned int iPdf = omtfConfig->nPdfBins() / 2;
0474 for (unsigned int iRefLayer = 0; iRefLayer < goldenPatterns[0]->getPdf()[0].size(); ++iRefLayer) {
0475 unsigned int refLayerLogicNumber = omtfConfig->getRefToLogicNumber()[iRefLayer];
0476 if (iRefLayer == 0 || iRefLayer == 2)
0477 step = 1.5;
0478 else if (iRefLayer == 5)
0479 step = 1.5;
0480 else if (iRefLayer == 1)
0481 step = 1.5;
0482 else if (iRefLayer == 3)
0483 step = 1.5;
0484 else if (iRefLayer == 5)
0485 step = 1.5;
0486 else if (iRefLayer == 6 || iRefLayer == 7)
0487 step = 1.5;
0488
0489 edm::LogImportant("l1tOmtfEventPrint")
0490 << __FUNCTION__ << ":" << __LINE__ << " RefLayer " << iRefLayer << " step " << step << std::endl;
0491 for (int sign = -1; sign <= 1; sign++) {
0492 for (auto& gp : boost::adaptors::reverse(goldenPatterns)) {
0493 if (gp->key().thePt == 0 || gp->key().theCharge != sign)
0494 continue;
0495
0496 double ptFrom = omtfConfig->getPatternPtRange(gp->key().theNumber).ptFrom;
0497 double ptTo = omtfConfig->getPatternPtRange(gp->key().theNumber).ptTo;
0498
0499 double ptRange = ptTo - ptFrom;
0500
0501 double minPdfValFactor = 0.1;
0502 double minPlog = log(omtfConfig->minPdfVal() * minPdfValFactor);
0503 double pdfMaxVal = omtfConfig->pdfMaxValue();
0504
0505 pdfMaxVal /= 3.;
0506 minPlog *= 2;
0507
0508
0509 if (ptRange > 800)
0510 ptRange = 800;
0511
0512 double norm = 0.001;
0513 double classProb = vxIntegMuRate(ptFrom, ptRange, 0.82, 1.24) * norm;
0514
0515 int digitisedVal = rint(pdfMaxVal - log(classProb) / minPlog * pdfMaxVal);
0516
0517 int newPdfVal = digitisedVal;
0518
0519 if (ptFrom == 0)
0520 newPdfVal += 15;
0521 if (ptFrom == 3.5)
0522 newPdfVal += 15;
0523 if (ptFrom == 4)
0524 newPdfVal += 12;
0525 if (ptFrom == 4.5)
0526 newPdfVal += 9;
0527 if (ptFrom == 5)
0528 newPdfVal += 7;
0529 if (ptFrom == 6)
0530 newPdfVal += 4;
0531 if (ptFrom == 7)
0532 newPdfVal += 2;
0533
0534 if (ptFrom == 100)
0535 newPdfVal = 16;
0536 if (ptFrom == 200)
0537 newPdfVal = 22;
0538
0539 gp->setPdfValue(newPdfVal, refLayerLogicNumber, iRefLayer, iPdf);
0540
0541 edm::LogImportant("l1tOmtfEventPrint")
0542 << gp->key() << " " << omtfConfig->getPatternPtRange(gp->key().theNumber).ptFrom << " - "
0543 << omtfConfig->getPatternPtRange(gp->key().theNumber).ptTo << " GeV"
0544 << " ptRange " << ptRange << " RefLayer " << iRefLayer << " newPdfVal " << newPdfVal << std::endl;
0545 }
0546 }
0547 }
0548 }
0549
0550 void PatternGenerator::reCalibratePt() {
0551 edm::LogImportant("l1tOmtfEventPrint") << __FUNCTION__ << ": " << __LINE__ << " reCalibratePt" << std::endl;
0552 std::map<int, float> ptMap;
0553
0554 ptMap[7] = 4.0;
0555 ptMap[8] = 4.5;
0556 ptMap[9] = 5.0;
0557 ptMap[10] = 5.5;
0558 ptMap[11] = 6.0;
0559 ptMap[13] = 7.0;
0560 ptMap[15] = 8.5;
0561 ptMap[17] = 10.0;
0562 ptMap[21] = 12.0;
0563 ptMap[25] = 14.0;
0564 ptMap[29] = 16.0;
0565 ptMap[33] = 18.5;
0566 ptMap[37] = 21.0;
0567 ptMap[41] = 23.0;
0568 ptMap[45] = 26.0;
0569 ptMap[49] = 28.0;
0570 ptMap[53] = 30.0;
0571 ptMap[57] = 32.0;
0572 ptMap[61] = 36.0;
0573 ptMap[71] = 40.0;
0574 ptMap[81] = 48.0;
0575 ptMap[91] = 54.0;
0576 ptMap[101] = 60.0;
0577 ptMap[121] = 70.0;
0578 ptMap[141] = 82.0;
0579 ptMap[161] = 96.0;
0580 ptMap[201] = 114.0;
0581 ptMap[401] = 200.0;
0582
0583 for (auto& gp : goldenPatterns) {
0584 if (gp->key().thePt == 0)
0585 continue;
0586
0587 int newPt = omtfConfig->ptGevToHw(ptMap[gp->key().thePt]);
0588 edm::LogImportant("l1tOmtfEventPrint") << gp->key().thePt << " -> " << newPt << std::endl;
0589
0590 gp->key().setPt(newPt);
0591 }
0592 }
0593
0594 void PatternGenerator::groupPatterns() {
0595 int group = 0;
0596 int indexInGroup = 0;
0597 for (auto& gp : goldenPatterns) {
0598 indexInGroup++;
0599 gp->key().setGroup(group);
0600 gp->key().setIndexInGroup(indexInGroup);
0601
0602
0603 edm::LogImportant("l1tOmtfEventPrint")
0604 << "setGroup(group): group " << group << " indexInGroup " << indexInGroup << std::endl;
0605
0606 if (gp->key().thePt <= 12 && indexInGroup == 2) {
0607 indexInGroup = 0;
0608 group++;
0609 }
0610
0611 if (gp->key().thePt > 12 && indexInGroup == 4) {
0612 indexInGroup = 0;
0613 group++;
0614 }
0615 }
0616
0617 OMTFConfiguration::vector2D patternGroups = omtfConfig->getPatternGroups(goldenPatterns);
0618 edm::LogImportant("l1tOmtfEventPrint") << "patternGroups:" << std::endl;
0619 for (unsigned int iGroup = 0; iGroup < patternGroups.size(); iGroup++) {
0620 edm::LogImportant("l1tOmtfEventPrint") << "patternGroup " << std::setw(2) << iGroup << " ";
0621 for (unsigned int i = 0; i < patternGroups[iGroup].size(); i++) {
0622 edm::LogImportant("l1tOmtfEventPrint") << i << " patNum " << patternGroups[iGroup][i] << " ";
0623 }
0624 edm::LogImportant("l1tOmtfEventPrint") << std::endl;
0625 }
0626
0627 int pdfBins = exp2(omtfConfig->nPdfAddrBits());
0628
0629 for (unsigned int iLayer = 0; iLayer < goldenPatterns.at(0)->getPdf().size(); ++iLayer) {
0630 for (unsigned int iRefLayer = 0; iRefLayer < goldenPatterns.at(0)->getPdf()[iLayer].size(); ++iRefLayer) {
0631
0632
0633 {
0634
0635 for (unsigned int iGroup = 0; iGroup < patternGroups.size(); iGroup++) {
0636 double meanDistPhi = 0;
0637 int mergedCnt = 0;
0638 for (unsigned int i = 0; i < patternGroups[iGroup].size(); i++) {
0639 auto gp = goldenPatterns.at(patternGroups[iGroup][i]).get();
0640 meanDistPhi += gp->meanDistPhiValue(iLayer, iRefLayer);
0641 if (gp->meanDistPhiValue(iLayer, iRefLayer) != 0)
0642 mergedCnt++;
0643 edm::LogImportant("l1tOmtfEventPrint")
0644 << __FUNCTION__ << ": " << __LINE__ << " iGroup " << iGroup << " numInGroup " << i << " " << gp->key()
0645 << " iLayer " << iLayer << " iRefLayer " << iRefLayer << " old meanDistPhiValue "
0646 << gp->meanDistPhiValue(iLayer, iRefLayer) << endl;
0647 }
0648
0649 if (mergedCnt) {
0650 meanDistPhi /= mergedCnt;
0651 meanDistPhi = (int)meanDistPhi;
0652
0653
0654 for (unsigned int i = 0; i < patternGroups[iGroup].size(); i++) {
0655 auto gp = goldenPatterns.at(patternGroups[iGroup][i]).get();
0656 unsigned int refLayerLogicNum = omtfConfig->getRefToLogicNumber()[iRefLayer];
0657 if (refLayerLogicNum != iLayer) {
0658 int shift = meanDistPhi - gp->meanDistPhiValue(iLayer, iRefLayer);
0659 edm::LogImportant("l1tOmtfEventPrint")
0660 << __FUNCTION__ << ": " << __LINE__ << " iGroup " << iGroup << " numInGroup " << i << " "
0661 << gp->key() << " iLayer " << iLayer << " iRefLayer " << iRefLayer
0662 << " new meanDistPhi after averaging " << meanDistPhi << " old meanDistPhiValue "
0663 << gp->meanDistPhiValue(iLayer, iRefLayer) << " shift " << shift << endl;
0664
0665 if (shift < 0) {
0666 for (int iBin = 1 - shift; iBin < pdfBins;
0667 iBin++) {
0668 auto pdfVal = gp->pdfValue(iLayer, iRefLayer, iBin);
0669 gp->setPdfValue(pdfVal, iLayer, iRefLayer, iBin + shift);
0670 edm::LogImportant("l1tOmtfEventPrint")
0671 << " iBin " << iBin << " iBin + shift " << iBin + shift << " pdfVal " << pdfVal << endl;
0672 }
0673 for (int iBin = pdfBins + shift; iBin < pdfBins; iBin++) {
0674 gp->setPdfValue(0, iLayer, iRefLayer, iBin);
0675 }
0676 } else if (shift > 0) {
0677 for (int iBin = pdfBins - 1 - shift; iBin > 0;
0678 iBin--) {
0679 auto pdfVal = gp->pdfValue(iLayer, iRefLayer, iBin);
0680 gp->setPdfValue(pdfVal, iLayer, iRefLayer, iBin + shift);
0681 edm::LogImportant("l1tOmtfEventPrint")
0682 << " iBin " << iBin << " iBin + shift " << iBin + shift << " pdfVal " << pdfVal << endl;
0683 }
0684 for (int iBin = shift; iBin > 0; iBin--) {
0685 gp->setPdfValue(0, iLayer, iRefLayer, iBin);
0686 }
0687 }
0688 }
0689
0690 gp->setMeanDistPhiValue(round(meanDistPhi), iLayer, iRefLayer);
0691 }
0692 }
0693 }
0694 }
0695 }
0696 }
0697 }