Back to home page

Project CMSSW displayed by LXR
 
 

    

File indexing completed on 2024-05-15 04:21:52 

0001 /*
0002  * PatternGenerator.cc
0003  *
0004  *  Created on: Nov 8, 2019
0005  *      Author: kbunkow
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                                    CandidateSimMuonMatcher* candidateSimMuonMatcher)
0020     : PatternOptimizerBase(edmCfg, omtfConfig, gps),
0021       updateStatFunction([this]() { updateStat(); }),
0022       candidateSimMuonMatcher(candidateSimMuonMatcher),
0023       eventCntPerGp(gps.size(), 0) {
0024   edm::LogImportant("l1tOmtfEventPrint") << "constructing PatternGenerator, type: "
0025                                          << edmCfg.getParameter<string>("patternGenerator") << std::endl;
0026 
0027   if (edmCfg.getParameter<string>("patternGenerator") == "patternGen" ||
0028       edmCfg.getParameter<string>("patternGenerator") == "2DHists" ||
0029       edmCfg.getParameter<string>("patternGenerator") == "deltaPhiVsPhiRef")
0030     initPatternGen();
0031 
0032   //2DHists are done for the displaced muons, then using the propagation for the matching is needed
0033   if (edmCfg.getParameter<string>("patternGenerator") == "2DHists" ||
0034       edmCfg.getParameter<string>("patternGenerator") == "deltaPhiVsPhiRef")
0035     updateStatFunction = [this]() { updateStatUsingMatcher2(); };
0036 
0037   if (edmCfg.exists("simTracksTag") == false)
0038     throw cms::Exception("PatternGenerator::PatternGenerator(): no simTracksTag !!!!!!!!!!!!!!!!!");
0039 
0040   if (!candidateSimMuonMatcher) {
0041     edm::LogImportant("l1tOmtfEventPrint") << "PatternGenerator: candidateSimMuonMatcher is null!!!!!!" << std::endl;
0042   }
0043 }
0044 
0045 PatternGenerator::~PatternGenerator() {}
0046 
0047 void PatternGenerator::initPatternGen() {
0048   //reseting the golden patterns
0049   unsigned int i = 0;
0050   for (auto& gp : goldenPatterns) {
0051     gp->setKeyNumber(i++);  //needed  if patterns were added
0052 
0053     if (gp->key().thePt == 0)
0054       continue;
0055 
0056     gp->reset();
0057 
0058     //1024 x 2048 is the maximum size that fits into 64GB of memory!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
0059     int statBinsCnt1 = 1024 * 2;  //TODO should be big enough to comprise the pdf tails
0060 
0061     int statBinsCnt2 = 1;  //for normal pattern generation
0062 
0063     if (edmCfg.getParameter<string>("patternGenerator") == "2DHists")
0064       statBinsCnt2 = 1024 * 2;
0065     //for 2D distribution, phiB vs phiDist, but if done for 8 ref layers, consumes too much memory
0066     else if (edmCfg.getParameter<string>("patternGenerator") == "deltaPhiVsPhiRef")
0067       statBinsCnt2 = omtfConfig->nPhiBins() / omtfConfig->nProcessors();
0068     //for 2D distribution, phiB vs phiDist, but if done for 8 ref layers, consumes too much memory
0069     //if(statBinsCnt2 > 10 && omtfConfig->nRefLayers() > 2)
0070     //  throw cms::Exception("PatternGenerator::initPatternGen(): statBinsCnt2 and omtfConfig->nRefLayers() too big, will consume too much memory");
0071 
0072     gp->iniStatisitics(statBinsCnt1, statBinsCnt2);
0073 
0074     if (statBinsCnt2 < 10 && sizeof(gp->getStatistics()[0][0][0][0]) < 4) {
0075       edm::LogImportant("l1tOmtfEventPrint")
0076           << "PatternGenerator::initPatternGen():" << __LINE__ << "sizeof gp statistics "
0077           << sizeof(gp->getStatistics()[0][0][0][0]) << std::endl;
0078       throw cms::Exception("PatternGenerator::initPatternGen(): getStatistics type is short!!!!");
0079     }
0080   }
0081 
0082   edm::LogImportant("l1tOmtfEventPrint") << "PatternGenerator::initPatternGen():" << __LINE__
0083                                          << " goldenPatterns.size() " << goldenPatterns.size() << std::endl;
0084 
0085   //setting all pdf to 1, this will cause that  when the OmtfProcessor process the input, the result will be based only on the number of fired layers,
0086   //and then the omtfCand will come from the processor that has the biggest number of fired layers
0087   //however, if the GoldenPatternResult::finalise3() is used - which just count the number of muonStubs (but do not check if it is valid, i.e. fired the pdf)
0088   // - the below does not matter
0089   for (auto& gp : goldenPatterns) {
0090     for (unsigned int iLayer = 0; iLayer < gp->getPdf().size(); ++iLayer) {
0091       for (unsigned int iRefLayer = 0; iRefLayer < gp->getPdf()[iLayer].size(); ++iRefLayer) {
0092         for (unsigned int iBin = 0; iBin < gp->getPdf()[iLayer][iRefLayer].size(); iBin++) {
0093           gp->pdfAllRef[iLayer][iRefLayer][iBin] = 1;
0094         }
0095       }
0096     }
0097   }
0098 
0099   //TODO uncomment if filling the ptDeltaPhiHist is needed
0100   /*  ptDeltaPhiHists.resize(2);
0101   for(unsigned int iCharge = 0; iCharge <= 1; iCharge++) {
0102     for(unsigned int iLayer = 0; iLayer < omtfConfig->nLayers(); ++iLayer) { //for the moment filing only ref layer, remove whe
0103       if(iLayer == 0 || iLayer == 2 || iLayer == 4 || iLayer == 6 || iLayer == 7 || iLayer == 10 || iLayer == 11  || iLayer == 16 || //refLayars
0104          iLayer == 1 || iLayer == 3 || iLayer == 5  ) //bending layers
0105       {
0106         ostringstream name;
0107         name<<"ptDeltaPhiHist_ch_"<<iCharge<<"_Layer_"<<iLayer;
0108         int phiFrom = -10;
0109         int phiTo   = 300; //TODO
0110         int phiBins = phiTo - phiFrom;
0111 
0112         if(iCharge == 1) {
0113           phiFrom = -300; //TODO
0114           phiTo = 10;
0115         }
0116 
0117         TH2I* ptDeltaPhiHist = new TH2I(name.str().c_str(), name.str().c_str(), 400, 0, 200, phiBins, phiFrom -0.5, phiTo -0.5);
0118         //cout<<"BinLowEdge "<<ptDeltaPhiHist->GetYaxis()->GetBinLowEdge(100)<<" BinUpEdge "<<ptDeltaPhiHist->GetYaxis()->GetBinUpEdge(100);
0119         ptDeltaPhiHists[iCharge].push_back(ptDeltaPhiHist);
0120       }
0121       else
0122         ptDeltaPhiHists[iCharge].push_back(nullptr);
0123     }
0124   }*/
0125 }
0126 
0127 void PatternGenerator::updateStat() {
0128   //cout<<__FUNCTION__<<":"<<__LINE__<<" omtfCand "<<*omtfCand<<std::endl;;
0129   AlgoMuon* algoMuon = omtfCand.get();
0130   if (!algoMuon) {
0131     edm::LogImportant("l1tOmtfEventPrint") << ":" << __LINE__ << " algoMuon is null" << std::endl;
0132     throw runtime_error("algoMuon is null");
0133   }
0134 
0135   simMuEta->Fill(simMuon->momentum().eta());
0136   candEta->Fill(omtfConfig->hwEtaToEta(regionalMuonCand.hwEta()));
0137 
0138   double ptSim = simMuon->momentum().pt();
0139   int chargeSim = (abs(simMuon->type()) == 13) ? simMuon->type() / -13 : 0;
0140 
0141   unsigned int exptPatNum = omtfConfig->getPatternNum(ptSim, chargeSim);
0142   GoldenPatternWithStat* exptCandGp = goldenPatterns.at(exptPatNum).get();  // expected pattern
0143 
0144   eventCntPerGp[exptPatNum]++;
0145 
0146   //edm::LogImportant("l1tOmtfEventPrint")<<"\n" <<__FUNCTION__<<": "<<__LINE__<<" exptCandGp "<<exptCandGp->key()<<" candProcIndx "<<candProcIndx<<" ptSim "<<ptSim<<" chargeSim "<<chargeSim<<std::endl;
0147 
0148   int pdfMiddle = 1 << (omtfConfig->nPdfAddrBits() - 1);
0149 
0150   //iRefHit is the index of the hit
0151   for (unsigned int iRefHit = 0; iRefHit < exptCandGp->getResults()[candProcIndx].size(); ++iRefHit) {
0152     auto& gpResult = exptCandGp->getResults()[candProcIndx][iRefHit];
0153 
0154     unsigned int refLayer = gpResult.getRefLayer();
0155 
0156     if (gpResult.getFiredLayerCnt() >= 3) {
0157       for (unsigned int iLayer = 0; iLayer < gpResult.getStubResults().size(); iLayer++) {
0158         //updating statistic for the gp which should have fired
0159 
0160         bool fired = false;
0161         if (gpResult.getStubResults()[iLayer].getMuonStub()) {
0162           if (omtfConfig->isBendingLayer(iLayer)) {
0163             if (gpResult.getStubResults()[iLayer].getMuonStub()->qualityHw >= 4)  //TODO change quality cut if needed
0164               fired = true;
0165           } else
0166             fired = true;
0167         }
0168 
0169         if (fired) {  //the result is not empty
0170           int phiDist = gpResult.getStubResults()[iLayer].getPdfBin();
0171           phiDist += exptCandGp->meanDistPhiValue(iLayer, refLayer) - pdfMiddle;
0172           //removing the shift applied in the GoldenPatternBase::process1Layer1RefLayer
0173 
0174           //TODO uncomment if filling ptDeltaPhiHists is needed
0175           /*
0176           unsigned int refLayerLogicNum = omtfConfig->getRefToLogicNumber()[iRefHit];
0177           if(ptDeltaPhiHists[iCharge][iLayer] != nullptr &&
0178               (iLayer == refLayerLogicNum || omtfConfig->getLogicToLogic().at(iLayer) == (int)refLayerLogicNum) )
0179             ptDeltaPhiHists[iCharge][iLayer]->Fill(ttAlgoMuon->getPt(), phiDist); //TODO correct
0180            */
0181 
0182           phiDist += exptCandGp->getStatistics()[iLayer][refLayer].size() / 2;
0183 
0184           //edm::LogImportant("l1tOmtfEventPrint")<<__FUNCTION__<<":"<<__LINE__<<" refLayer "<<refLayer<<" iLayer "<<iLayer<<" phiDist "<<phiDist<<" getPdfBin "<<gpResult.getStubResults()[iLayer].getPdfBin()<<std::endl;
0185           if (phiDist > 0 && phiDist < (int)(exptCandGp->getStatistics()[iLayer][refLayer].size())) {
0186             //updating statistic for the gp which found the candidate
0187             //edm::LogImportant("l1tOmtfEventPrint")<<__FUNCTION__<<":"<<__LINE__<<" updating statistic "<<std::endl;
0188             exptCandGp->updateStat(iLayer, refLayer, phiDist, 0, 1);
0189           }
0190         } else {  //if there is no hit at all in a given layer, the bin = 0 is filled
0191           int phiDist = 0;
0192           exptCandGp->updateStat(iLayer, refLayer, phiDist, 0, 1);
0193         }
0194       }
0195     }
0196   }
0197 }
0198 
0199 void PatternGenerator::updateStatUsingMatcher2() {
0200   //cout<<__FUNCTION__<<":"<<__LINE__<<" omtfCand "<<*omtfCand<<std::endl;;
0201 
0202   std::vector<MatchingResult> matchingResults = candidateSimMuonMatcher->getMatchingResults();
0203   LogTrace("l1tOmtfEventPrint") << "matchingResults.size() " << matchingResults.size() << std::endl;
0204 
0205   //candidateSimMuonMatcher should use the  trackingParticles, because the simTracks are not stored for the pile-up events
0206   for (auto& matchingResult : matchingResults) {
0207     if (matchingResult.muonCand && matchingResult.simTrack) {
0208       //&& matchingResult.muonCand->hwQual() >= 12 &&
0209       //matchingResult.muonCand->hwPt() > 38
0210 
0211       AlgoMuon* algoMuon = matchingResult.procMuon.get();
0212       if (!algoMuon) {
0213         edm::LogImportant("l1tOmtfEventPrint") << ":" << __LINE__ << " algoMuon is null" << std::endl;
0214         throw runtime_error("algoMuon is null");
0215       }
0216 
0217       double ptSim = matchingResult.simTrack->momentum().pt();
0218       int chargeSim = (abs(matchingResult.simTrack->type()) == 13) ? matchingResult.simTrack->type() / -13 : 0;
0219 
0220       double muDxy = (-1 * matchingResult.simVertex->position().x() * matchingResult.simTrack->momentum().py() +
0221                       matchingResult.simVertex->position().y() * matchingResult.simTrack->momentum().px()) /
0222                      matchingResult.simTrack->momentum().pt();
0223 
0224       simMuEta->Fill(matchingResult.simTrack->momentum().eta());
0225 
0226       simMuPtVsDispl->Fill(matchingResult.simTrack->momentum().pt(), muDxy);
0227       simMuPtVsRho->Fill(matchingResult.simTrack->momentum().pt(), matchingResult.simVertex->position().rho());
0228 
0229       unsigned int exptPatNum = omtfConfig->getPatternNum(ptSim, chargeSim);
0230       GoldenPatternWithStat* exptCandGp = goldenPatterns.at(exptPatNum).get();  // expected pattern
0231 
0232       eventCntPerGp[exptPatNum]++;
0233 
0234       candProcIndx =
0235           omtfConfig->getProcIndx(matchingResult.muonCand->processor(), matchingResult.muonCand->trackFinderType());
0236 
0237       //edm::LogImportant("l1tOmtfEventPrint")<<"\n" <<__FUNCTION__<<": "<<__LINE__<<" exptCandGp "<<exptCandGp->key()<<" candProcIndx "<<candProcIndx<<" ptSim "<<ptSim<<" chargeSim "<<chargeSim<<std::endl;
0238 
0239       int pdfMiddle = 1 << (omtfConfig->nPdfAddrBits() - 1);
0240       LogTrace("l1tOmtfEventPrint") << "updateStatUsingMatcher2 " << __LINE__ << std::endl;
0241       //iRefHit is the index of the hit
0242       for (unsigned int iRefHit = 0; iRefHit < exptCandGp->getResults()[candProcIndx].size(); ++iRefHit) {
0243         auto& gpResult = exptCandGp->getResults()[candProcIndx][iRefHit];
0244 
0245         if (gpResult.getFiredLayerCnt() >= 3) {
0246           int refLayer = gpResult.getRefLayer();
0247 
0248           if (refLayer < 0 || !gpResult.isValid())
0249             LogTrace("l1tOmtfEventPrint") << "updateStatUsingMatcher2 " << __LINE__ << " refLayer " << refLayer
0250                                           << " gpResult.isValid() " << gpResult.isValid() << std::endl;
0251 
0252           int refLayerLogicNumber = omtfConfig->getRefToLogicNumber()[refLayer];
0253 
0254           LogTrace("l1tOmtfEventPrint") << __FUNCTION__ << ":" << __LINE__ << " updating statistic: candProcIndx "
0255                                         << candProcIndx << " iRefHit " << iRefHit << " refLayer " << refLayer
0256                                         << " exptPatNum " << exptPatNum << " ptSim " << ptSim << " chargeSim "
0257                                         << chargeSim << " muDxy " << muDxy << " muRho "
0258                                         << matchingResult.simVertex->position().rho() << " x "
0259                                         << matchingResult.simVertex->position().x() << " y "
0260                                         << matchingResult.simVertex->position().y() << " z "
0261                                         << matchingResult.simVertex->position().z() << std::endl;
0262 
0263           int refPhiB = 0;
0264 
0265           if (omtfConfig->isBendingLayer(refLayerLogicNumber + 1))
0266             //if(refLayerLogicNumber < 5)
0267             refPhiB = gpResult.getStubResults()[refLayerLogicNumber].getMuonStub()->phiBHw;
0268 
0269           int refPhiBShifted = refPhiB + exptCandGp->getStatistics()[0][refLayer][0].size() / 2;
0270 
0271           if (edmCfg.getParameter<string>("patternGenerator") == "deltaPhiVsPhiRef") {
0272             refPhiBShifted = gpResult.getStubResults()[refLayerLogicNumber].getMuonStub()->phiHw;
0273           }
0274 
0275           if (refPhiBShifted < 0 || refPhiBShifted >= (int)exptCandGp->getStatistics()[0][refLayer][0].size()) {
0276             edm::LogImportant("l1tOmtfEventPrint") << "\n"
0277                                                    << __FUNCTION__ << ": " << __LINE__ << " wrong refPhiB " << refPhiB
0278                                                    << " refPhiBShifted " << refPhiBShifted;
0279             continue;
0280           }
0281 
0282           for (unsigned int iLayer = 0; iLayer < gpResult.getStubResults().size(); iLayer++) {
0283             //updating statistic for the gp which should have fired
0284 
0285             bool fired = false;
0286             if (gpResult.getStubResults()[iLayer].getMuonStub()) {
0287               fired = true;
0288             }
0289 
0290             if (fired) {            //the result is not empty
0291               int meanDistPhi = 0;  //exptCandGp->meanDistPhiValue(iLayer, refLayer, refPhiB);  //should be 0 here
0292 
0293               int phiDist = gpResult.getStubResults()[iLayer].getPdfBin() + meanDistPhi - pdfMiddle;
0294               //removing the shift applied in the GoldenPatternBase::process1Layer1RefLayer
0295 
0296               int lutMiddle = exptCandGp->getStatistics()[iLayer][refLayer].size() / 2;
0297 
0298               LogTrace("l1tOmtfEventPrint")
0299                   << __FUNCTION__ << ":" << __LINE__ << " refLayer " << refLayer << " iLayer " << iLayer << " phiDist "
0300                   << phiDist << " getPdfBin " << gpResult.getStubResults()[iLayer].getPdfBin() << " phiMean "
0301                   << meanDistPhi << " refPhiB " << refPhiB << std::endl;
0302 
0303               //updating statistic for the gp which found the candidate
0304               //edm::LogImportant("l1tOmtfEventPrint")<<__FUNCTION__<<":"<<__LINE__<<" updating statistic "<<std::endl;
0305 
0306               int phiDistCorr = phiDist + lutMiddle;
0307 
0308               if (phiDistCorr > 0 && phiDistCorr < (int)(exptCandGp->getStatistics()[iLayer][refLayer].size())) {
0309                 LogTrace("l1tOmtfEventPrint") << __FUNCTION__ << ":" << __LINE__ << " phiDistCorr "
0310                                               << phiDistCorr
0311                                               //<< " phiDistCorr + lutMiddle "<< phiDistCorr + lutMiddle
0312                                               << " refPhiBShifted " << refPhiBShifted << std::endl;
0313                 exptCandGp->updateStat(iLayer, refLayer, phiDistCorr, refPhiBShifted, 1);
0314               }
0315             } else {  //if there is no hit at all in a given layer, the bin = 0 is filled
0316               int phiDist = 0;
0317               exptCandGp->updateStat(iLayer, refLayer, phiDist, refPhiBShifted, 1);
0318             }
0319           }
0320         }
0321       }
0322     }
0323   }
0324 }
0325 
0326 void PatternGenerator::observeEventEnd(const edm::Event& iEvent,
0327                                        std::unique_ptr<l1t::RegionalMuonCandBxCollection>& finalCandidates) {
0328   if (simMuon == nullptr || omtfCand->getGoldenPatern() == nullptr)  //no sim muon or empty candidate
0329     return;
0330 
0331   if (abs(simMuon->momentum().eta()) < 0.8 || abs(simMuon->momentum().eta()) > 1.24)
0332     return;
0333 
0334   PatternOptimizerBase::observeEventEnd(iEvent, finalCandidates);
0335 
0336   //updateStat();
0337   //updateStatUsingMatcher2();
0338   updateStatFunction();
0339 }
0340 
0341 void PatternGenerator::endJob() {
0342   if (edmCfg.getParameter<string>("patternGenerator") == "modifyClassProb")
0343     modifyClassProb(1);
0344   else if (edmCfg.getParameter<string>("patternGenerator") == "groupPatterns")
0345     groupPatterns();
0346   else if (edmCfg.getParameter<string>("patternGenerator") == "patternGen") {
0347     upadatePdfs();
0348     writeLayerStat = true;
0349   } else if (edmCfg.getParameter<string>("patternGenerator") == "2DHists") {
0350     upadatePdfs();
0351     writeLayerStat = true;
0352   } else if (edmCfg.getParameter<string>("patternGenerator") == "deltaPhiVsPhiRef") {
0353     upadatePdfs();
0354     writeLayerStat = true;
0355   } else if (edmCfg.getParameter<string>("patternGenerator") == "patternGenFromStat") {
0356     std::string rootFileName = edmCfg.getParameter<edm::FileInPath>("patternsROOTFile").fullPath();
0357     edm::LogImportant("l1tOmtfEventPrint") << "PatternGenerator::endJob() rootFileName " << rootFileName << std::endl;
0358     TFile inFile(rootFileName.c_str());
0359     TDirectory* curDir = (TDirectory*)inFile.Get("layerStats");
0360 
0361     ostringstream ostrName;
0362     for (auto& gp : goldenPatterns) {
0363       if (gp->key().thePt == 0)
0364         continue;
0365 
0366       for (unsigned int iLayer = 0; iLayer < gp->getPdf().size(); ++iLayer) {
0367         for (unsigned int iRefLayer = 0; iRefLayer < gp->getPdf()[iLayer].size(); ++iRefLayer) {
0368           ostrName.str("");
0369           ostrName << "histLayerStat_PatNum_" << gp->key().theNumber << "_refLayer_" << iRefLayer << "_Layer_"
0370                    << iLayer;
0371 
0372           TH1I* histLayerStat = (TH1I*)curDir->Get(ostrName.str().c_str());
0373 
0374           if (histLayerStat) {
0375             int statBinsCnt = 1024 * 2;
0376             if ((int)(gp->getStatistics()[iLayer][iRefLayer].size()) != histLayerStat->GetNbinsX()) {
0377               statBinsCnt = histLayerStat->GetNbinsX();
0378               gp->iniStatisitics(statBinsCnt, 1);  //TODO
0379               edm::LogImportant("l1tOmtfEventPrint")
0380                   << "PatternGenerator::endJob() - " << ostrName.str() << "statBinsCnt = " << statBinsCnt << std::endl;
0381             }
0382 
0383             for (int iBin = 0; iBin < statBinsCnt; iBin++) {
0384               gp->updateStat(iLayer, iRefLayer, iBin, 0, histLayerStat->GetBinContent(iBin + 1));
0385             }
0386           } else {
0387             edm::LogImportant("l1tOmtfEventPrint")
0388                 << "PatternGenerator::endJob() - reading histLayerStat: histogram not found " << ostrName.str()
0389                 << std::endl;
0390           }
0391         }
0392       }
0393     }
0394 
0395     TH1* simMuFoundByOmtfPt_fromFile = (TH1*)inFile.Get("simMuFoundByOmtfPt");
0396     for (unsigned int iGp = 0; iGp < eventCntPerGp.size(); iGp++) {
0397       eventCntPerGp[iGp] = simMuFoundByOmtfPt_fromFile->GetBinContent(simMuFoundByOmtfPt_fromFile->FindBin(iGp));
0398       edm::LogImportant("l1tOmtfEventPrint")
0399           << "PatternGenerator::endJob() - eventCntPerGp: iGp" << iGp << " - " << eventCntPerGp[iGp] << std::endl;
0400     }
0401 
0402     //TODO chose the desired grouping ///////////////
0403     int group = 0;
0404     int indexInGroup = 0;
0405     for (auto& gp : goldenPatterns) {
0406       indexInGroup++;
0407       gp->key().setGroup(group);
0408       gp->key().setIndexInGroup(indexInGroup);
0409       //indexInGroup is counted from 1
0410 
0411       edm::LogImportant("l1tOmtfEventPrint")
0412           << "setGroup(group): group " << group << " indexInGroup " << indexInGroup << std::endl;
0413 
0414       if (gp->key().thePt <= 10 && indexInGroup == 2) {  //TODO
0415         indexInGroup = 0;
0416         group++;
0417       }
0418 
0419       if (gp->key().thePt > 10 && indexInGroup == 4) {  //TODO
0420         indexInGroup = 0;
0421         group++;
0422       }
0423     }  /////////////////////////////////////////////
0424 
0425     upadatePdfs();
0426 
0427     modifyClassProb(1);
0428 
0429     //groupPatterns(); IMPORTANT don't call grouping here, just set the groups above!!!!
0430 
0431     reCalibratePt();
0432     this->writeLayerStat = true;
0433   }
0434 
0435   PatternOptimizerBase::endJob();
0436 }
0437 
0438 void PatternGenerator::upadatePdfs() {
0439   //TODO setting the DistPhiBitShift i.e. grouping of the pdfBins
0440   for (auto& gp : goldenPatterns) {
0441     if (gp->key().thePt == 0)
0442       continue;
0443     for (unsigned int iLayer = 0; iLayer < gp->getPdf().size(); ++iLayer) {
0444       for (unsigned int iRefLayer = 0; iRefLayer < gp->getPdf()[iLayer].size(); ++iRefLayer) {
0445         if (gp->getDistPhiBitShift(iLayer, iRefLayer)) {
0446           throw runtime_error(
0447               string(__FUNCTION__) + ":" + to_string(__LINE__) +
0448               "gp->getDistPhiBitShift(iLayer, iRefLayer) != 0 -  cannot change DistPhiBitShift then!!!!");
0449         }
0450 
0451         //watch out - the pt here is the hardware pt before the recalibration
0452         //the shift for given pattern and layer should be the same same for all refLayers
0453         //otherwise the firmware does not compile - at least the phase-1
0454         //for phase2
0455         /*if ((gp->key().thePt <= 10) &&
0456             (iLayer == 1 || iLayer == 3 || iLayer == 5)) {  //iRefLayer: MB2, iLayer: MB1 and MB2 phiB
0457           gp->setDistPhiBitShift(2, iLayer, iRefLayer);
0458         } else if ((gp->key().thePt <= 10) && (iLayer == 10)) {  //iRefLayer: MB2, iLayer: RB1_in
0459           gp->setDistPhiBitShift(1, iLayer, iRefLayer);
0460         } else if ((gp->key().thePt >= 11 && gp->key().thePt <= 17) && (iLayer == 1)) {  //MB1 phiB
0461           //due to grouping the patterns 4-7, the pdfs for the layer 1 in the pattern go outside of the range
0462           //so the shift must be increased (or the group should be divided into to 2 groups, but it will increase fw occupancy
0463           gp->setDistPhiBitShift(2, iLayer, iRefLayer);
0464         } else if ((gp->key().thePt >= 11 && gp->key().thePt <= 17) && (iLayer == 3 || iLayer == 5)) {  //MB1 phiB
0465           //due to grouping the patterns 4-7, the pdfs for the layer 1 in the pattern go outside of the range
0466           //so the shift must be increased (or the group should be divided into to 2 groups, but it will increase fw occupancy
0467           gp->setDistPhiBitShift(1, iLayer, iRefLayer);
0468         } else
0469           gp->setDistPhiBitShift(0, iLayer, iRefLayer);*/
0470 
0471         //for phase1
0472         if ((gp->key().thePt <= 8) &&
0473             (iLayer == 1 || iLayer == 3 || iLayer == 5)) {  //iRefLayer: MB2, iLayer: MB1 and MB2 phiB
0474           gp->setDistPhiBitShift(2, iLayer, iRefLayer);
0475         } else if ((gp->key().thePt <= 10) && (iLayer == 10)) {  //iRefLayer: MB2, iLayer: RB1_in
0476           gp->setDistPhiBitShift(1, iLayer, iRefLayer);
0477         } else if ((gp->key().thePt <= 10) &&
0478                    (iLayer == 1 || iLayer == 3 || iLayer == 5)) {  //iRefLayer: MB2, iLayer: MB1 and MB2 phiB
0479           gp->setDistPhiBitShift(1, iLayer, iRefLayer);
0480         } else if ((gp->key().thePt <= 17) && (iLayer == 1)) {  //MB1 phiB
0481           //due to grouping the patterns 4-7, the pdfs for the layer 1 in the pattern go outside of the range
0482           //so the shift must be increased (or the group should be divided into to 2 groups, but it will increase fw occupancy
0483           gp->setDistPhiBitShift(1, iLayer, iRefLayer);
0484         } else if ((gp->key().thePt >= 11 && gp->key().thePt <= 17) && (iLayer == 3 || iLayer == 5)) {  //MB1 phiB
0485           //due to grouping the patterns 4-7, the pdfs for the layer 1 in the pattern go outside of the range
0486           //so the shift must be increased (or the group should be divided into to 2 groups, but it will increase fw occupancy
0487           gp->setDistPhiBitShift(0, iLayer, iRefLayer);
0488         } else
0489           gp->setDistPhiBitShift(0, iLayer, iRefLayer);
0490       }
0491     }
0492   }
0493 
0494   double minHitCntThresh = 0.001;
0495   //Calculating meanDistPhi
0496   for (auto& gp : goldenPatterns) {
0497     if (gp->key().thePt == 0)
0498       continue;
0499 
0500     int minHitCnt = minHitCntThresh * eventCntPerGp[gp->key().number()];  // //TODO tune threshold <<<<<<<<<<<<<<<<<<
0501     edm::LogImportant("l1tOmtfEventPrint")
0502         << "PatternGenerator::upadatePdfs() Calculating meanDistPhi " << gp->key() << " eventCnt "
0503         << eventCntPerGp[gp->key().number()] << " minHitCnt " << minHitCnt << std::endl;
0504 
0505     for (unsigned int iLayer = 0; iLayer < gp->getPdf().size(); ++iLayer) {
0506       for (unsigned int iRefLayer = 0; iRefLayer < gp->getPdf()[iLayer].size(); ++iRefLayer) {
0507         //calculate meanDistPhi
0508         double meanDistPhi = 0;
0509         double count = 0;
0510         for (unsigned int iBin = 1; iBin < gp->getStatistics()[iLayer][iRefLayer].size(); iBin++) {
0511           //iBin = 0 is reserved for the no hit
0512           meanDistPhi += iBin * gp->getStatistics()[iLayer][iRefLayer][iBin][0];
0513           count += gp->getStatistics()[iLayer][iRefLayer][iBin][0];
0514         }
0515 
0516         if (count != 0) {
0517           meanDistPhi /= count;
0518 
0519           meanDistPhi -= (gp->getStatistics()[iLayer][iRefLayer].size() / 2);
0520 
0521           if (count < minHitCnt)
0522             meanDistPhi = 0;
0523           else
0524             edm::LogImportant("l1tOmtfEventPrint")
0525                 << __FUNCTION__ << ": " << __LINE__ << " " << gp->key() << " iLayer " << iLayer << " iRefLayer "
0526                 << iRefLayer << " count " << count << " meanDistPhi " << meanDistPhi << endl;
0527         }
0528         gp->setMeanDistPhiValue(round(meanDistPhi), iLayer, iRefLayer);
0529       }
0530     }
0531   }
0532 
0533   OMTFConfiguration::vector2D patternGroups = omtfConfig->getPatternGroups(goldenPatterns);
0534   edm::LogImportant("l1tOmtfEventPrint") << "patternGroups:" << std::endl;
0535   for (unsigned int iGroup = 0; iGroup < patternGroups.size(); iGroup++) {
0536     edm::LogImportant("l1tOmtfEventPrint") << "patternGroup " << std::setw(2) << iGroup << " ";
0537     for (unsigned int i = 0; i < patternGroups[iGroup].size(); i++) {
0538       edm::LogImportant("l1tOmtfEventPrint") << i << " patNum " << patternGroups[iGroup][i] << " ";
0539     }
0540     edm::LogImportant("l1tOmtfEventPrint") << std::endl;
0541   }
0542 
0543   //averaging the meanDistPhi for the gp belonging to the same group
0544   for (unsigned int iLayer = 0; iLayer < goldenPatterns.at(0)->getPdf().size(); ++iLayer) {
0545     for (unsigned int iRefLayer = 0; iRefLayer < goldenPatterns.at(0)->getPdf()[iLayer].size(); ++iRefLayer) {
0546       //unsigned int refLayerLogicNum = omtfConfig->getRefToLogicNumber()[iRefLayer];
0547       //if(refLayerLogicNum == iLayer)
0548       {
0549         for (unsigned int iGroup = 0; iGroup < patternGroups.size(); iGroup++) {
0550           double meanDistPhi = 0;
0551           int mergedCnt = 0;
0552           double norm = 0;
0553           for (unsigned int i = 0; i < patternGroups[iGroup].size(); i++) {
0554             auto gp = goldenPatterns.at(patternGroups[iGroup][i]).get();
0555             if (gp->meanDistPhiValue(iLayer, iRefLayer) != 0) {
0556               double weight = 1. / gp->key().thePt;
0557               meanDistPhi += weight * gp->meanDistPhiValue(iLayer, iRefLayer);
0558               mergedCnt++;
0559               norm += weight;
0560             }
0561           }
0562 
0563           if (mergedCnt) {
0564             //because for some gps the statistics can be too low, and then the meanDistPhiValue is 0, so it should not contribute to meanDistPhi, therefore it is divide by mergedCnt
0565             //meanDistPhi /= mergedCnt;
0566             ///weighted average, weight is 1/pt
0567             //for low pT patterns it shifts the pdf of the pattern with bigger width (i.e. lower pt) towards the center of LUT
0568             //then higher value of shift can be avoided (sometimes). So this is just a simple trick
0569             meanDistPhi /= norm;
0570 
0571             /*
0572             //setting the meanDistPhi to 0 if it is already small - this should save logic in FPGA - but seems it does not
0573             if (iLayer == 2) {
0574               //the meanDistPhi for the iLayer == 2 i.e. MB2 is used to calculate the algoMuon output phi
0575               //therefore it is not zero-ed, as it will affect this output phi, phi and thus e.g. ghostbusting
0576             } else if (abs(round(meanDistPhi)) <= 3)
0577               meanDistPhi = 0;
0578             else if (goldenPatterns.at(patternGroups[iGroup][0]).get()->key().thePt >= 13) {
0579               //RPC layers, one strip is 4.7 units, the minimal possible spacing between two RPC hits is 2 strips
0580               if (iLayer >= 10 && abs(round(meanDistPhi)) <= 8)
0581                 meanDistPhi = 0;
0582               else if (abs(round(meanDistPhi)) <= 5)
0583                 meanDistPhi = 0;
0584             } */
0585 
0586             for (unsigned int i = 0; i < patternGroups[iGroup].size(); i++) {
0587               auto gp = goldenPatterns.at(patternGroups[iGroup][i]).get();
0588               gp->setMeanDistPhiValue(round(meanDistPhi), iLayer, iRefLayer);
0589               edm::LogImportant("l1tOmtfEventPrint")
0590                   << __FUNCTION__ << ": " << __LINE__ << " iGroup " << iGroup << " numInGroup " << i << " " << gp->key()
0591                   << " iLayer " << iLayer << " iRefLayer " << iRefLayer << " meanDistPhi after averaging "
0592                   << meanDistPhi << endl;
0593             }
0594           }
0595         }
0596       }
0597     }
0598   }
0599 
0600   //calculating the pdfs
0601   for (auto& gp : goldenPatterns) {
0602     if (gp->key().thePt == 0)
0603       continue;
0604 
0605     //TODO tune threshold <<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
0606     int minHitCnt = 2 * minHitCntThresh * eventCntPerGp[gp->key().number()];
0607 
0608     for (unsigned int iLayer = 0; iLayer < gp->getPdf().size(); ++iLayer) {
0609       for (unsigned int iRefLayer = 0; iRefLayer < gp->getPdf()[iLayer].size(); ++iRefLayer) {
0610         {
0611           double norm = 0;
0612           for (unsigned int iBin = 0; iBin < gp->getStatistics()[iLayer][iRefLayer].size();
0613                iBin++) {  //iBin = 0 i.e. no hit is included here, to have the proper norm
0614             norm += gp->getStatistics()[iLayer][iRefLayer][iBin][0];
0615           }
0616 
0617           int pdfMiddle = gp->getPdf()[iLayer][iRefLayer].size() / 2;
0618           int statBinGroupSize = 1 << gp->getDistPhiBitShift(iLayer, iRefLayer);
0619           for (unsigned int iBinPdf = 0; iBinPdf < gp->getPdf()[iLayer][iRefLayer].size(); iBinPdf++) {
0620             double pdfVal = 0;
0621             if (iBinPdf > 0) {
0622               for (int i = 0; i < statBinGroupSize; i++) {
0623                 int iBinStat =
0624                     statBinGroupSize * ((int)(iBinPdf)-pdfMiddle) + i + gp->meanDistPhiValue(iLayer, iRefLayer);
0625 
0626                 iBinStat += (gp->getStatistics()[iLayer][iRefLayer].size() / 2);
0627 
0628                 if (iBinStat >= 0 && iBinStat < (int)gp->getStatistics()[iLayer][iRefLayer].size()) {
0629                   pdfVal += gp->getStatistics()[iLayer][iRefLayer][iBinStat][0];
0630                 }
0631               }
0632               if (norm > minHitCnt) {
0633                 pdfVal /= (norm * statBinGroupSize);
0634               } else
0635                 pdfVal = 0;
0636             } else {  //iBinPdf == 0 i.e. no hit
0637               int iBinStat = 0;
0638               if (norm > 0) {
0639                 pdfVal = gp->getStatistics()[iLayer][iRefLayer][iBinStat][0] / norm;
0640               }
0641               edm::LogImportant("l1tOmtfEventPrint")
0642                   << __FUNCTION__ << ": " << __LINE__ << " " << gp->key() << "calculating pdf: iLayer " << iLayer
0643                   << " iRefLayer " << iRefLayer << " norm " << std::setw(5) << norm << " no hits cnt " << std::setw(5)
0644                   << gp->getStatistics()[iLayer][iRefLayer][iBinStat][0] << " pdfVal " << pdfVal << endl;
0645             }
0646 
0647             double minPdfValFactor = 1;
0648             const double minPlog = log(omtfConfig->minPdfVal() * minPdfValFactor);
0649             const double pdfMaxVal = omtfConfig->pdfMaxValue();
0650 
0651             int digitisedVal = 0;
0652             if (pdfVal >= omtfConfig->minPdfVal() * minPdfValFactor) {
0653               digitisedVal = rint(pdfMaxVal - log(pdfVal) / minPlog * pdfMaxVal);
0654             }
0655 
0656             gp->setPdfValue(digitisedVal, iLayer, iRefLayer, iBinPdf);
0657           }
0658         }
0659       }
0660     }
0661   }
0662 }
0663 
0664 void PatternGenerator::saveHists(TFile& outfile) {
0665   outfile.mkdir("ptDeltaPhiHists")->cd();
0666   //TODO uncomment if ptDeltaPhiHists are needed
0667   /*  for(unsigned int iCharge = 0; iCharge <= 1; iCharge++) {
0668     for(unsigned int iLayer = 0; iLayer < omtfConfig->nLayers(); ++iLayer) { //for the moment filing only ref layer, remove whe
0669       if(ptDeltaPhiHists[iCharge][iLayer]) {
0670         ptDeltaPhiHists[iCharge][iLayer]->Write();
0671       }
0672     }
0673   }*/
0674 }
0675 
0676 void PatternGenerator::modifyClassProb(double step) {
0677   edm::LogImportant("l1tOmtfEventPrint") << __FUNCTION__ << ": " << __LINE__ << " Correcting P(C_k) " << std::endl;
0678   unsigned int iPdf = omtfConfig->nPdfBins() / 2;  // <<(omtfConfig->nPdfAddrBits()-1);
0679   for (unsigned int iRefLayer = 0; iRefLayer < goldenPatterns[0]->getPdf()[0].size(); ++iRefLayer) {
0680     unsigned int refLayerLogicNumber = omtfConfig->getRefToLogicNumber()[iRefLayer];
0681     if (iRefLayer == 0 || iRefLayer == 2)  //DT
0682       step = 1.5;
0683     else if (iRefLayer == 5)  //DT
0684       step = 1.5;
0685     else if (iRefLayer == 1)  //CSC
0686       step = 1.5;
0687     else if (iRefLayer == 3)  //CSC
0688       step = 1.5;
0689     else if (iRefLayer == 5)  //RE2/3
0690       step = 1.5;
0691     else if (iRefLayer == 6 || iRefLayer == 7)  //bRPC
0692       step = 1.5;
0693 
0694     edm::LogImportant("l1tOmtfEventPrint")
0695         << __FUNCTION__ << ":" << __LINE__ << " RefLayer " << iRefLayer << " step " << step << std::endl;
0696     for (int sign = -1; sign <= 1; sign++) {
0697       for (auto& gp : boost::adaptors::reverse(goldenPatterns)) {
0698         if (gp->key().thePt == 0 || gp->key().theCharge != sign)
0699           continue;
0700 
0701         double ptFrom = omtfConfig->getPatternPtRange(gp->key().theNumber).ptFrom;
0702         double ptTo = omtfConfig->getPatternPtRange(gp->key().theNumber).ptTo;
0703 
0704         double ptRange = ptTo - ptFrom;
0705 
0706         double minPdfValFactor = 0.1;
0707         double minPlog = log(omtfConfig->minPdfVal() * minPdfValFactor);
0708         double pdfMaxVal = omtfConfig->pdfMaxValue();
0709 
0710         pdfMaxVal /= 3.;
0711         minPlog *= 2;
0712 
0713         //last bin of the ptRange goes to 10000, so here we change it to 1000
0714         if (ptRange > 800)
0715           ptRange = 800;
0716 
0717         double norm = 0.001;
0718         double classProb = vxIntegMuRate(ptFrom, ptRange, 0.82, 1.24) * norm;
0719 
0720         int digitisedVal = rint(pdfMaxVal - log(classProb) / minPlog * pdfMaxVal);
0721 
0722         int newPdfVal = digitisedVal;  //gp->getPdf()[refLayerLogicNumber][iRefLayer][iPdf]
0723         //watch out - the pt here is before re-calibration
0724         if (ptFrom == 0)
0725           newPdfVal += 15;
0726         if (ptFrom == 3.5)
0727           newPdfVal += 15;
0728         if (ptFrom == 4)
0729           newPdfVal += 12;
0730         if (ptFrom == 4.5)
0731           newPdfVal += 9;
0732         if (ptFrom == 5)
0733           newPdfVal += 7;
0734         if (ptFrom == 6)
0735           newPdfVal += 4;
0736         if (ptFrom == 7)
0737           newPdfVal += 2;
0738 
0739         //if (ptFrom == 20) //pattern Key_13
0740         //  newPdfVal += 1;
0741         if (ptFrom >= 22 && ptFrom <= 26)
0742           newPdfVal += 2;
0743         if (ptFrom == 28)  //pattern Key_17
0744           newPdfVal += 1;
0745 
0746         if (ptFrom == 100)
0747           newPdfVal = 16;
0748         if (ptFrom == 200)
0749           newPdfVal = 22;
0750 
0751         gp->setPdfValue(newPdfVal, refLayerLogicNumber, iRefLayer, iPdf);
0752 
0753         edm::LogImportant("l1tOmtfEventPrint")
0754             << gp->key() << " " << omtfConfig->getPatternPtRange(gp->key().theNumber).ptFrom << " - "
0755             << omtfConfig->getPatternPtRange(gp->key().theNumber).ptTo << " GeV"
0756             << " ptRange " << ptRange << " RefLayer " << iRefLayer << " newPdfVal " << newPdfVal << std::endl;
0757       }
0758     }
0759   }
0760 }
0761 
0762 void PatternGenerator::reCalibratePt() {
0763   edm::LogImportant("l1tOmtfEventPrint") << __FUNCTION__ << ": " << __LINE__ << " reCalibratePt" << std::endl;
0764   std::map<int, float> ptMap;
0765   //for Patterns_0x0009_oldSample_3_10Files_classProb2.xml
0766   ptMap[7] = 4.0;
0767   ptMap[8] = 4.5;
0768   ptMap[9] = 5.0;
0769   ptMap[10] = 5.5;
0770   ptMap[11] = 6.0;
0771   ptMap[13] = 7.0;
0772   ptMap[15] = 8.5;
0773   ptMap[17] = 10.0;
0774   ptMap[21] = 12.0;
0775   ptMap[25] = 14.0;
0776   ptMap[29] = 16.0;
0777   ptMap[33] = 18.5;
0778   ptMap[37] = 21.0;
0779   ptMap[41] = 23.0;
0780   ptMap[45] = 26.0;
0781   ptMap[49] = 28.0;
0782   ptMap[53] = 30.0;
0783   ptMap[57] = 32.0;
0784   ptMap[61] = 36.0;
0785   ptMap[71] = 40.0;
0786   ptMap[81] = 48.0;
0787   ptMap[91] = 54.0;
0788   ptMap[101] = 60.0;
0789   ptMap[121] = 70.0;
0790   ptMap[141] = 82.0;
0791   ptMap[161] = 96.0;
0792   ptMap[201] = 114.0;
0793   ptMap[401] = 200.0;
0794 
0795   for (auto& gp : goldenPatterns) {
0796     if (gp->key().thePt == 0)
0797       continue;
0798 
0799     int newPt = omtfConfig->ptGevToHw(ptMap[gp->key().thePt]);
0800     edm::LogImportant("l1tOmtfEventPrint") << gp->key().thePt << " -> " << newPt << std::endl;
0801 
0802     gp->key().setPt(newPt);
0803   }
0804 }
0805 
0806 void PatternGenerator::groupPatterns() {
0807   int group = 0;
0808   int indexInGroup = 0;
0809   for (auto& gp : goldenPatterns) {
0810     indexInGroup++;
0811     gp->key().setGroup(group);
0812     gp->key().setIndexInGroup(indexInGroup);
0813     //indexInGroup is counted from 1
0814 
0815     edm::LogImportant("l1tOmtfEventPrint")
0816         << "setGroup(group): group " << group << " indexInGroup " << indexInGroup << std::endl;
0817 
0818     if (gp->key().thePt <= 12 && indexInGroup == 2) {  //TODO
0819       indexInGroup = 0;
0820       group++;
0821     }
0822 
0823     if (gp->key().thePt > 12 && indexInGroup == 4) {  //TODO
0824       indexInGroup = 0;
0825       group++;
0826     }
0827   }
0828 
0829   OMTFConfiguration::vector2D patternGroups = omtfConfig->getPatternGroups(goldenPatterns);
0830   edm::LogImportant("l1tOmtfEventPrint") << "patternGroups:" << std::endl;
0831   for (unsigned int iGroup = 0; iGroup < patternGroups.size(); iGroup++) {
0832     edm::LogImportant("l1tOmtfEventPrint") << "patternGroup " << std::setw(2) << iGroup << " ";
0833     for (unsigned int i = 0; i < patternGroups[iGroup].size(); i++) {
0834       edm::LogImportant("l1tOmtfEventPrint") << i << " patNum " << patternGroups[iGroup][i] << " ";
0835     }
0836     edm::LogImportant("l1tOmtfEventPrint") << std::endl;
0837   }
0838 
0839   int pdfBins = exp2(omtfConfig->nPdfAddrBits());
0840 
0841   for (unsigned int iLayer = 0; iLayer < goldenPatterns.at(0)->getPdf().size(); ++iLayer) {
0842     for (unsigned int iRefLayer = 0; iRefLayer < goldenPatterns.at(0)->getPdf()[iLayer].size(); ++iRefLayer) {
0843       //unsigned int refLayerLogicNum = omtfConfig->getRefToLogicNumber()[iRefLayer];
0844       //if(refLayerLogicNum == iLayer)
0845       {
0846         //averaging the meanDistPhi for the gp belonging to the same group
0847         for (unsigned int iGroup = 0; iGroup < patternGroups.size(); iGroup++) {
0848           double meanDistPhi = 0;
0849           int mergedCnt = 0;
0850           for (unsigned int i = 0; i < patternGroups[iGroup].size(); i++) {
0851             auto gp = goldenPatterns.at(patternGroups[iGroup][i]).get();
0852             meanDistPhi += gp->meanDistPhiValue(iLayer, iRefLayer);
0853             if (gp->meanDistPhiValue(iLayer, iRefLayer) != 0)
0854               mergedCnt++;
0855             edm::LogImportant("l1tOmtfEventPrint")
0856                 << __FUNCTION__ << ": " << __LINE__ << " iGroup " << iGroup << " numInGroup " << i << " " << gp->key()
0857                 << " iLayer " << iLayer << " iRefLayer " << iRefLayer << " old meanDistPhiValue "
0858                 << gp->meanDistPhiValue(iLayer, iRefLayer) << endl;
0859           }
0860 
0861           if (mergedCnt) {
0862             meanDistPhi /= mergedCnt;
0863             meanDistPhi = (int)meanDistPhi;
0864 
0865             //because for some gps the statistics can be too low, and then the meanDistPhiValue is 0, so it should not contribute
0866             for (unsigned int i = 0; i < patternGroups[iGroup].size(); i++) {
0867               auto gp = goldenPatterns.at(patternGroups[iGroup][i]).get();
0868               unsigned int refLayerLogicNum = omtfConfig->getRefToLogicNumber()[iRefLayer];
0869               if (refLayerLogicNum != iLayer) {
0870                 int shift = meanDistPhi - gp->meanDistPhiValue(iLayer, iRefLayer);
0871                 edm::LogImportant("l1tOmtfEventPrint")
0872                     << __FUNCTION__ << ": " << __LINE__ << " iGroup " << iGroup << " numInGroup " << i << " "
0873                     << gp->key() << " iLayer " << iLayer << " iRefLayer " << iRefLayer
0874                     << " new meanDistPhi after averaging " << meanDistPhi << " old meanDistPhiValue "
0875                     << gp->meanDistPhiValue(iLayer, iRefLayer) << " shift " << shift << endl;
0876 
0877                 if (shift < 0) {
0878                   for (int iBin = 1 - shift; iBin < pdfBins;
0879                        iBin++) {  //iBin = 0 i.e. no hit is included here, to have the proper norm
0880                     auto pdfVal = gp->pdfValue(iLayer, iRefLayer, iBin);
0881                     gp->setPdfValue(pdfVal, iLayer, iRefLayer, iBin + shift);
0882                     edm::LogImportant("l1tOmtfEventPrint")
0883                         << "     iBin " << iBin << "  iBin + shift " << iBin + shift << " pdfVal " << pdfVal << endl;
0884                   }
0885                   for (int iBin = pdfBins + shift; iBin < pdfBins; iBin++) {
0886                     gp->setPdfValue(0, iLayer, iRefLayer, iBin);
0887                   }
0888                 } else if (shift > 0) {
0889                   for (int iBin = pdfBins - 1 - shift; iBin > 0;
0890                        iBin--) {  //iBin = 0 i.e. no hit is included here, to have the proper norm
0891                     auto pdfVal = gp->pdfValue(iLayer, iRefLayer, iBin);
0892                     gp->setPdfValue(pdfVal, iLayer, iRefLayer, iBin + shift);
0893                     edm::LogImportant("l1tOmtfEventPrint")
0894                         << "     iBin " << iBin << "  iBin + shift " << iBin + shift << " pdfVal " << pdfVal << endl;
0895                   }
0896                   for (int iBin = shift; iBin > 0; iBin--) {
0897                     gp->setPdfValue(0, iLayer, iRefLayer, iBin);
0898                   }
0899                 }
0900               }
0901 
0902               gp->setMeanDistPhiValue(round(meanDistPhi), iLayer, iRefLayer);
0903             }
0904           }
0905         }
0906       }
0907     }
0908   }
0909 }
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated by the 2.2.1 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!