Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-05-29 04:20:42

0001 /*
0002  * PtAssignmentNN.cc
0003  *
0004  *  Created on: May 8, 2020
0005  *      Author: kbunkow
0006  */
0007 
0008 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0009 #include "FWCore/Framework/interface/EventSetup.h"
0010 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0011 
0012 #include "DataFormats/MuonDetId/interface/CSCDetId.h"
0013 
0014 #include "L1Trigger/L1TMuonOverlapPhase2/interface/PtAssignmentNNRegression.h"
0015 
0016 #include <boost/archive/text_oarchive.hpp>
0017 #include <boost/archive/text_iarchive.hpp>
0018 #include <boost/algorithm/string/replace.hpp>
0019 
0020 #include <sstream>
0021 #include <fstream>
0022 
0023 namespace lutNN {
0024   static constexpr int input_I = 10;
0025   static constexpr int input_F = 4;
0026   static constexpr std::size_t networkInputSize = 18;
0027 
0028   static constexpr int layer1_neurons = 16;
0029   static constexpr int layer1_lut_I = 3;
0030   static constexpr int layer1_lut_F = 13;
0031 
0032   static constexpr int layer1_output_I = 4;
0033   //4 bits are for the count of the noHit layers which goes to the input of the layer2
0034   static constexpr int layer2_input_I = 8;
0035 
0036   static constexpr int layer2_neurons = 9;
0037   static constexpr int layer2_lut_I = 5;
0038   static constexpr int layer2_lut_F = 11;
0039 
0040   static constexpr int layer3_input_I = 5;
0041 
0042   static constexpr int layer3_0_inputCnt = 8;
0043   static constexpr int layer3_0_lut_I = 5;
0044   static constexpr int layer3_0_lut_F = 11;
0045   static constexpr int output0_I = 8;
0046   static constexpr int output0_F = 2;
0047 
0048   static constexpr int layer3_1_inputCnt = 1;
0049   static constexpr int layer3_1_lut_I = 4;
0050   static constexpr int layer3_1_lut_F = 11;
0051   static constexpr int output1_I = 8;
0052   static constexpr int output1_F = 0;  //Does not matter in principle - it is not used
0053 
0054   typedef LutNetworkFixedPointRegression2Outputs<input_I,
0055                                                  input_F,
0056                                                  networkInputSize,
0057                                                  layer1_lut_I,
0058                                                  layer1_lut_F,
0059                                                  layer1_neurons,  //layer1_lutSize = 2 ^ input_I
0060                                                  layer1_output_I,
0061                                                  layer2_input_I,
0062                                                  layer2_lut_I,
0063                                                  layer2_lut_F,
0064                                                  layer2_neurons,
0065                                                  layer3_input_I,
0066                                                  layer3_0_inputCnt,
0067                                                  layer3_0_lut_I,
0068                                                  layer3_0_lut_F,
0069                                                  output0_I,
0070                                                  output0_F,
0071                                                  layer3_1_inputCnt,
0072                                                  layer3_1_lut_I,
0073                                                  layer3_1_lut_F,
0074                                                  output1_I,
0075                                                  output1_F>
0076       LutNetworkFP;
0077 }  // namespace lutNN
0078 
0079 PtAssignmentNNRegression::PtAssignmentNNRegression(const edm::ParameterSet& edmCfg,
0080                                                    const OMTFConfiguration* omtfConfig,
0081                                                    std::string networkFile)
0082     : PtAssignmentBase(omtfConfig), lutNetworkFP(make_unique<lutNN::LutNetworkFP>()) {
0083   std::ifstream ifs(networkFile);
0084 
0085   edm::LogImportant("OMTFReconstruction")
0086       << " " << __FUNCTION__ << ":" << __LINE__ << " networkFile " << networkFile << std::endl;
0087 
0088   lutNetworkFP->load(networkFile);
0089 
0090   edm::LogImportant("OMTFReconstruction") << " " << __FUNCTION__ << ":" << __LINE__ << std::endl;
0091 }
0092 
0093 struct OmtfHit {
0094   union {
0095     unsigned long rawData = 0;
0096 
0097     struct {
0098       char layer;
0099       char quality;
0100       char z;
0101       char valid;
0102       short eta;
0103       short phiDist;
0104     };
0105   };
0106 
0107   OmtfHit(unsigned long rawData) : rawData(rawData) {}
0108 };
0109 
0110 bool omtfHitToEventInput(OmtfHit& hit, std::vector<float>& inputs, unsigned int omtfRefLayer, bool print) {
0111   float offset = (omtfRefLayer << 7) + 64;
0112 
0113   if (hit.valid) {
0114     if ((hit.layer == 1 || hit.layer == 3 || hit.layer == 5) && hit.quality < 4)  ///TODO <<<<<<<<<<<<<<<<<<<<<<<<<<<<<
0115       return false;
0116 
0117     int rangeFactor = 2;  //rangeFactor scales the hit.phiDist such that the event->inputs is smaller then 63
0118     if (hit.layer == 1) {
0119       rangeFactor = 8;
0120     }
0121     /*else if(hit.layer == 8 || hit.layer == 17) {
0122             rangeFactor = 4;
0123         }*/
0124     else if (hit.layer == 3) {
0125       rangeFactor = 4;
0126     } else if (hit.layer == 9) {
0127       rangeFactor = 1;
0128     }
0129     /*else {
0130             rangeFactor = 2;
0131         }
0132          */
0133 
0134     rangeFactor *= 2;  //TODO !!!!!!!!!!!!!!!!!!!
0135 
0136     if (std::abs(hit.phiDist) >= (63 * rangeFactor)) {
0137       edm::LogImportant("OMTFReconstruction")  //<<" muonPt "<<omtfEvent.muonPt<<" omtfPt "<<omtfEvent.omtfPt
0138           << " RefLayer " << omtfRefLayer << " layer " << int(hit.layer) << " hit.phiDist " << hit.phiDist << " valid "
0139           << ((short)hit.valid) << " !!!!!!!!!!!!!!!!!!!!!!!!" << endl;
0140       hit.phiDist = copysign(63 * rangeFactor, hit.phiDist);
0141     }
0142 
0143     inputs.at(hit.layer) = (float)hit.phiDist / (float)rangeFactor + offset;
0144 
0145     if (inputs.at(hit.layer) >=
0146         1022)  //the last address i.e. 1023 is reserved for the no-hit value, so interpolation between the 1022 and 1023 has no sense
0147       inputs.at(hit.layer) = 1022;
0148 
0149     if (print || inputs.at(hit.layer) < 0) {
0150       edm::LogImportant("OMTFReconstruction")  //<<"rawData "<<hex<<setw(16)<<hit.rawData
0151           << " layer " << dec << int(hit.layer);
0152       edm::LogImportant("OMTFReconstruction")
0153           << " phiDist " << hit.phiDist << " inputVal " << inputs.at(hit.layer) << " hit.z " << int(hit.z) << " valid "
0154           << ((short)hit.valid) << " quality " << (short)hit.quality << " omtfRefLayer " << omtfRefLayer;
0155       if (inputs.at(hit.layer) < 0)
0156         edm::LogImportant("OMTFReconstruction") << " event->inputs.at(hit.layer) < 0 !!!!!!!!!!!!!!!!!" << endl;
0157       edm::LogImportant("OMTFReconstruction") << endl;
0158     }
0159 
0160     if (inputs[hit.layer] >= 1024) {  //TODO should be the size of the LUT of the first layer
0161       edm::LogImportant("OMTFReconstruction") << " event->inputs[hit.layer] >= 1024 !!!!!!!!!!!!!!!!!" << endl;
0162     }
0163     return true;
0164   }
0165 
0166   return false;
0167 }
0168 
0169 std::vector<float> PtAssignmentNNRegression::getPts(AlgoMuons::value_type& algoMuon,
0170                                                     std::vector<std::unique_ptr<IOMTFEmulationObserver>>& observers) {
0171   LogTrace("l1tOmtfEventPrint") << " " << __FUNCTION__ << ":" << __LINE__ << std::endl;
0172   auto& gpResult = algoMuon->getGpResultConstr();
0173   //int pdfMiddle = 1<<(omtfConfig->nPdfAddrBits()-1);
0174 
0175   LogTrace("l1tOmtfEventPrint") << " " << __FUNCTION__ << ":" << __LINE__ << std::endl;
0176   /*
0177   edm::LogVerbatim("l1tOmtfEventPrint")<<"DataROOTDumper2:;observeEventEnd muonPt "<<event.muonPt<<" muonCharge "<<event.muonCharge
0178       <<" omtfPt "<<event.omtfPt<<" RefLayer "<<event.omtfRefLayer<<" omtfPtCont "<<event.omtfPtCont
0179       <<std::endl;
0180 */
0181 
0182   unsigned int inputCnt = 18;  //TDOO!!!!!
0183   unsigned int outputCnt = 2;
0184   const float noHitVal = 1023.;
0185 
0186   //edm::LogImportant("OMTFReconstruction") <<"\n----------------------"<<endl;
0187   //edm::LogImportant("OMTFReconstruction") <<(*algoMuon)<<std::endl;
0188 
0189   std::vector<float> inputs(inputCnt, noHitVal);
0190 
0191   for (unsigned int iLogicLayer = 0; iLogicLayer < gpResult.getStubResults().size(); ++iLogicLayer) {
0192     auto& stubResult = gpResult.getStubResults()[iLogicLayer];
0193     if (stubResult.getMuonStub()) {  //&& stubResult.getValid() //TODO!!!!!!!!!!!!!!!!1
0194       int hitPhi = stubResult.getMuonStub()->phiHw;
0195       unsigned int refLayerLogicNum = omtfConfig->getRefToLogicNumber()[algoMuon->getRefLayer()];
0196       int phiRefHit = gpResult.getStubResults()[refLayerLogicNum].getMuonStub()->phiHw;
0197 
0198       if (omtfConfig->isBendingLayer(iLogicLayer)) {
0199         hitPhi = stubResult.getMuonStub()->phiBHw;
0200         phiRefHit = 0;  //phi ref hit for the banding layer set to 0, since it should not be included in the phiDist
0201       }
0202 
0203       OmtfHit hit(0);
0204 
0205       hit.layer = iLogicLayer;
0206       hit.quality = stubResult.getMuonStub()->qualityHw;
0207       hit.eta = stubResult.getMuonStub()->etaHw;  //in which scale?
0208       hit.valid = stubResult.getValid();
0209 
0210       //phiDist = hitPhi - phiRefHit;
0211       hit.phiDist = hitPhi - phiRefHit;
0212 
0213       /*
0214       LogTrace("l1tOmtfEventPrint") <<" muonPt "<<event.muonPt<<" omtfPt "<<event.omtfPt<<" RefLayer "<<event.omtfRefLayer
0215           <<" layer "<<int(hit.layer)<<" PdfBin "<<stubResult.getPdfBin()<<" hit.phiDist "<<hit.phiDist<<" valid "<<stubResult.getValid()<<" " //<<" phiDist "<<phiDist
0216           <<" getDistPhiBitShift "<<omtfCand->getGoldenPatern()->getDistPhiBitShift(iLogicLayer, omtfCand->getRefLayer())
0217           <<" meanDistPhiValue   "<<omtfCand->getGoldenPatern()->meanDistPhiValue(iLogicLayer, omtfCand->getRefLayer())//<<(phiDist != hit.phiDist? "!!!!!!!<<<<<" : "")
0218           <<endl;*/
0219 
0220       if (hit.phiDist > 504 || hit.phiDist < -512) {
0221         edm::LogVerbatim("l1tOmtfEventPrint")
0222             //<<" muonPt "<<event.muonPt<<" omtfPt "<<event.omtfPt<<" RefLayer "<<event.omtfRefLayer
0223             << " layer " << int(hit.layer) << " hit.phiDist " << hit.phiDist << " valid " << stubResult.getValid()
0224             << " !!!!!!!!!!!!!!!!!!!!!!!!" << std::endl;
0225       }
0226 
0227       DetId detId(stubResult.getMuonStub()->detId);
0228       if (detId.subdetId() == MuonSubdetId::CSC) {
0229         CSCDetId cscId(detId);
0230         hit.z = cscId.chamber() % 2;
0231       }
0232 
0233       LogTrace("l1tOmtfEventPrint") << "hit: layer " << (int)hit.layer << " quality " << (int)hit.quality << " eta "
0234                                     << (int)hit.eta << " valid " << (int)hit.valid << " phiDist " << (int)hit.phiDist
0235                                     << " z " << (int)hit.z << std::endl;
0236 
0237       omtfHitToEventInput(hit, inputs, algoMuon->getRefLayer(), false);
0238     }
0239   }
0240 
0241   LogTrace("l1tOmtfEventPrint") << " " << __FUNCTION__ << ":" << __LINE__ << std::endl;
0242 
0243   std::vector<double> nnResult(outputCnt);
0244   lutNetworkFP->run(inputs, noHitVal, nnResult);
0245 
0246   LogTrace("l1tOmtfEventPrint") << " " << __FUNCTION__ << ":" << __LINE__ << std::endl;
0247 
0248   double pt = std::copysign(nnResult.at(0), nnResult.at(1));
0249 
0250   LogTrace("l1tOmtfEventPrint") << " " << __FUNCTION__ << ":" << __LINE__ << " nnResult.at(0) " << nnResult.at(0)
0251                                 << " nnResult.at(1) " << nnResult.at(1) << " pt " << pt << std::endl;
0252 
0253   std::vector<float> pts;
0254   pts.emplace_back(pt);
0255 
0256   //algoMuon->setPtNN(omtfConfig->ptGevToHw(nnResult.at(0)));
0257   auto calibratedHwPt = lutNetworkFP->getCalibratedHwPt();
0258   algoMuon->setPtNNConstr(calibratedHwPt);
0259 
0260   algoMuon->setChargeNNConstr(nnResult[1] >= 0 ? 1 : -1);
0261 
0262   //TODO add some if here, such that the property_tree is filled only when needed
0263   boost::property_tree::ptree procDataTree;
0264   for (unsigned int i = 0; i < inputs.size(); i++) {
0265     auto& inputTree = procDataTree.add("input", "");
0266     inputTree.add("<xmlattr>.num", i);
0267     inputTree.add("<xmlattr>.val", inputs[i]);
0268   }
0269 
0270   std::ostringstream ostr;
0271   ostr << std::fixed << std::setprecision(19) << nnResult.at(0);
0272   procDataTree.add("output0.<xmlattr>.val", ostr.str());
0273 
0274   ostr.str("");
0275   ostr << std::fixed << std::setprecision(19) << nnResult.at(1);
0276   procDataTree.add("output1.<xmlattr>.val", ostr.str());
0277 
0278   procDataTree.add("calibratedHwPt.<xmlattr>.val", calibratedHwPt);
0279 
0280   procDataTree.add("hwSign.<xmlattr>.val", algoMuon->getChargeNNConstr() < 0 ? 1 : 0);
0281 
0282   for (auto& obs : observers)
0283     obs->addProcesorData("regressionNN", procDataTree);
0284 
0285   return pts;
0286 
0287   //event.print();
0288   /*
0289   std::vector<float> pts(classifierToRegressions.size(), 0);
0290 
0291   unsigned int i =0;
0292   for(auto& classifierToRegression : classifierToRegressions) {
0293     auto orgValue = classifierToRegression->getValue(&event);
0294     auto absOrgValue = std::abs(orgValue);
0295     pts.at(i) = classifierToRegression->getCalibratedValue(absOrgValue);
0296     pts.at(i) = std::copysign(pts.at(i), orgValue);
0297 
0298     LogTrace("OMTFReconstruction") <<" "<<__FUNCTION__<<":"<<__LINE__<<" orgValue "<<orgValue<<" pts["<<i<<"] "<<pts[i]<<std::endl;
0299     //std::cout<<"nn pts["<<i<<"] "<<pts[i]<< std::endl;
0300     i++;
0301   }
0302 
0303   return pts;*/
0304 }