RPCTrigger/src/RPCStripsRing.cc

0001 // -*- C++ -*-
0002 //
0003 // Package:     RPCConeBuilder
0004 // Class  :     RPCStripsRing
0005 //
0006 // Implementation:
0007 //     <Notes on implementation>
0008 //
0009 // Original Author:  Tomasz Fruboes
0010 //         Created:  Tue Feb 26 15:13:10 CET 2008
0011 //
0012
0013 // system include files
0014
0015 // user include files
0016 //#include "L1TriggerConfig/RPCConeBuilder/interface/RPCStripsRing.h"
0017 #include "L1Trigger/RPCTrigger/interface/RPCStripsRing.h"
0018 #include "Geometry/RPCGeometry/interface/RPCGeomServ.h"
0019 #include "DataFormats/MuonDetId/interface/RPCDetId.h"
0020 #include "Geometry/RPCGeometry/interface/RPCGeometry.h"
0021 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0022
0023 RPCStripsRing::RPCStripsRing()
0024     : m_hwPlane(-1),
0025       m_etaPartition(99),
0026       m_region(-2),
0027       m_isReferenceRing(false),
0028       m_didVirtuals(false),
0029       m_didFiltering(false) {}
0030
0031 RPCStripsRing::RPCStripsRing(const RPCRoll* roll, std::shared_ptr<L1RPCConeBuilder::TConMap> cmap)
0032     : m_didVirtuals(false), m_didFiltering(false), m_connectionsMap(cmap) {
0033   RPCDetId detId = roll->id();
0034   RPCGeomServ grs(detId);
0035
0036   m_etaPartition = grs.eta_partition();
0037   m_hwPlane = calculateHwPlane(roll);
0038
0039   m_isReferenceRing = false;
0040
0041   m_region = detId.region();
0042
0043   int ring = detId.ring();
0044
0045   if (m_region == 0 && std::abs(ring) < 2 && m_hwPlane == 2)  // for barell wheel -1,0,1 refplane is hwPlane=2
0046     m_isReferenceRing = true;
0047   else if (m_region == 0 && std::abs(ring) == 2 && m_hwPlane == 6)  // for barell wheel -2,2 refplane is hwPlane=6
0048     m_isReferenceRing = true;
0049   else if (m_region != 0 && m_hwPlane == 2)  // for endcaps
0050     m_isReferenceRing = true;
0051
0052   if (getRingId() == 2008 || getRingId() == 2108)  //exception: endcaps;hwplane 2;farest roll from beam
0053     m_isReferenceRing = false;
0054
0055   addRoll(roll);
0056 }
0057
0058 void RPCStripsRing::addRoll(const RPCRoll* roll) {
0059   //  RPCDetId detId = roll->id();
0060
0061   if (getRingId() != getRingId(roll)) {
0062     throw cms::Exception("RPCInternal") << "RPCStripsRing::addRoll ringsIds dont match \n";
0063   }
0064
0065   //iterate over the strips of this roll
0066   for (int i = 1; i <= roll->nstrips(); i++) {
0067     LocalPoint lStripCentre = roll->centreOfStrip(i);
0068     GlobalPoint gStripCentre = roll->toGlobal(lStripCentre);
0069     float phiRaw = gStripCentre.phi();
0070
0071     TStrip newStrip(roll->id().rawId(), i);
0072     (*this)[phiRaw] = newStrip;
0073   }
0074 }
0075
0076 int RPCStripsRing::getRingId(int etaPart, int hwPlane) {
0077   int sign = 1;  // positive
0078   if (etaPart < 0) {
0079     sign = 0;
0080   }
0081
0082   return 1000 * (hwPlane) +        //1...6
0083          100 * (sign) +            //
0084          1 * (std::abs(etaPart));  //-17...17
0085 }
0086
0087 int RPCStripsRing::getRingId() { return getRingId(m_etaPartition, m_hwPlane); }
0088
0089 int RPCStripsRing::getRingId(const RPCRoll* roll) {
0090   RPCDetId detId = roll->id();
0091   RPCGeomServ grs(detId);
0092   int etaPartition = grs.eta_partition();
0093   int hwPlane = calculateHwPlane(roll);
0094
0095   return getRingId(etaPartition, hwPlane);
0096 }
0097
0098 //  hwPlane is  station number for endcaps
0099 //  for barrell numbering goes 1 5 2 6 3 4 (first number means plane closest to the beam)
0100 int RPCStripsRing::calculateHwPlane(const RPCRoll* roll) {
0101   int hwPlane = -1;
0102   RPCDetId detId = roll->id();
0103   int station = detId.station();
0104   int layer = detId.layer();
0105   int region = detId.region();
0106
0107   if (region != 0) {  // endcaps
0108     hwPlane = station;
0109   }
0110   // Now comes the barell
0111   else if (station > 2) {
0112     hwPlane = station;
0113   } else if (station == 1 && layer == 1) {
0114     hwPlane = 1;
0115   } else if (station == 1 && layer == 2) {
0116     hwPlane = 5;
0117   } else if (station == 2 && layer == 1) {
0118     hwPlane = 2;
0119   } else if (station == 2 && layer == 2) {
0120     hwPlane = 6;
0121   }
0122
0123   /*if (hwPlane < 1)
0124     std::cout << "prb: " << hwPlane << " "
0125         << region << " "
0126         << station << " "
0127   << layer << std::endl;*/
0128   if (hwPlane < 0) {
0129     throw cms::Exception("RPCInternal") << "Calculated negative hwplane \n";
0130   }
0131
0132   return hwPlane;
0133 }
0134
0135 void RPCStripsRing::filterOverlapingChambers() {
0136   if (m_didFiltering)
0137     return;
0138   m_didFiltering = true;
0139
0140   if (m_region != 0 || m_hwPlane != 4)
0141     return;
0142
0143   typedef std::map<uint32_t, int> TDetId2StripNo;
0144   TDetId2StripNo det2stripNo;
0145
0146   // Note: we begin in middle of first chamber (ch1), we have to handle that
0147   int ch1BegStrips = 0;  // no of strips on the begining of the map (first=last chamber of map)
0148   int ch1EndStrips = 0;  // no of strips on the end of the map (first=last chamber of map)
0149
0150   // How many strips has each chamber?
0151   RPCStripsRing::iterator it = this->begin();
0152   uint32_t ch1Det = it->second.m_detRawId;
0153   for (; it != this->end(); ++it) {
0154     if (det2stripNo.find(it->second.m_detRawId) == det2stripNo.end()) {
0155       det2stripNo[it->second.m_detRawId] = 1;  // Add new chamber to a map, set strip cnt to 1
0156     } else {
0157       ++det2stripNo[it->second.m_detRawId];  // Increase strip count of a chamber
0158     }
0159
0160     if (det2stripNo.size() == 1 && ch1Det == it->second.m_detRawId) {
0161       ++ch1BegStrips;
0162     } else if (ch1Det == it->second.m_detRawId) {
0163       ++ch1EndStrips;
0164     }
0165   }
0166
0167   det2stripNo[ch1Det] -= ch1EndStrips;
0168
0169   // std::cout << ch1BegStrips << " " << ch1EndStrips << std::endl;
0170
0171   //TDetId2StripNo::iterator itIds = det2stripNo.begin();
0172   //for(;itIds!=det2stripNo.end();++itIds){
0173   //    std::cout << itIds->first << " " << itIds->second << std::endl;
0174   //  }
0175
0176   it = this->begin();
0177   uint32_t lastDet = it->second.m_detRawId;
0178   while (it != this->end()) {
0179     if (det2stripNo[it->second.m_detRawId] < 0) {
0180       throw cms::Exception("RPCInternal") << " RPCStripsRing::filterOverlapingChambers() - no strips left \n";
0181     }
0182     if (it->second.m_detRawId == lastDet) {
0183       --det2stripNo[lastDet];
0184       ++it;
0185     } else if (det2stripNo[lastDet] == 0) {  // no more strips left in lastDet, proceed to new det
0186
0187       if (lastDet == ch1Det) {
0188         det2stripNo[ch1Det] += ch1EndStrips;
0189       }
0190
0191       lastDet = it->second.m_detRawId;
0192       --det2stripNo[lastDet];
0193       ++it;
0194     } else {  // there are still strips in last det, delete current strip
0195       --det2stripNo[it->second.m_detRawId];
0196       RPCStripsRing::iterator itErase = it;
0197       ++it;
0198       //std::cout << "Removing strip " <<  it->second.m_detRawId << " " << (int)it->second.m_strip << std::endl;
0199       this->erase(itErase);
0200     }
0201   }
0202 }
0203
0204 void RPCStripsRing::fillWithVirtualStrips() {
0205   if (m_didVirtuals)
0206     return;
0207   m_didVirtuals = true;
0208
0209   const float pi = 3.141592654;
0210   double dphi = 2.0 * pi / 1152;  // defines angular granulation of strips.
0211
0212   RPCStripsRing stripsToInsert;
0213
0214   float delta = 0;
0215   int stripsToAdd = 0;
0216
0217   RPCStripsRing::iterator it = this->begin();
0218   RPCStripsRing::iterator itLast = this->begin();
0219   for (; it != this->end(); ++it) {
0220     /*std::cout << it->first << " "
0221         << it->second.m_detRawId << " "
0222         << (int)it->second.m_strip << std::endl;
0223     */
0224
0225     delta = it->first - itLast->first;
0226     if (it == itLast ||                                        // skip first loop iteration
0227         itLast->second.m_detRawId == it->second.m_detRawId ||  // insert strips between two chambers only
0228         delta < 0) {
0229       itLast = it;
0230       continue;
0231     }
0232
0233     stripsToAdd = (int)std::floor(delta / dphi) - 1;
0234     //std::cout << delta << " " << stripsToAdd << std::endl;
0235
0236     if (isReferenceRing() && m_hwPlane == 6)
0237       ++stripsToAdd;
0238
0239     for (int i = 0; i < stripsToAdd; ++i) {
0240       stripsToInsert[itLast->first + dphi * (i + 1)] = TStrip();
0241     }
0242
0243     itLast = it;
0244   }
0245   // TODO: check delta between first and last strip in map
0246
0247   this->insert(stripsToInsert.begin(), stripsToInsert.end());
0248 }
0249 void RPCStripsRing::createRefConnections(TOtherConnStructVec& otherRings, int logplane, int logplaneSize) {
0250   //*
0251   /*std::cout << "RefCon for " << getRingId()
0252        << " (" << getEtaPartition()<<  ")"
0253        << " tower: " << getTowerForRefRing()
0254        << " ; connected: "
0255        << otherRings.size()
0256        << std::endl
0257        << std::endl;
0258   //*/
0259
0260   // XXX - TODO: warning on wrong logplaneSize
0261
0262   if (!this->isReferenceRing()) {
0263     throw cms::Exception("RPCInternal") << " RPCStripsRing::createRefConnections "
0264                                         << " called for non-reference ring \n";
0265   }
0266
0267   /*
0268    if (logplaneSize!=8) {
0269      throw cms::Exception("RPCInternal") << " RPCStripsRing::createRefConnections "
0270          << " called for lpSize " << logplaneSize << " \n";
0271
0272    }*/
0273   const float pi = 3.141592654;
0274   const float offset = (5. / 360.) * 2 * pi;  // XXX
0275
0276   //find first reference strip of first PAC (the strip with phi ~= 5deg)
0277   RPCStripsRing::iterator starEndIt = this->begin();
0278   while ((++starEndIt)->first < offset)
0279     ;
0280
0281   RPCStripsRing::iterator it = starEndIt;
0282   //--starEndIt;
0283
0284   float angle = 0;
0285   int curPACno = -1;
0286   int curStripNo = 0;
0287   int curBegStripNo = 0;
0288
0289   bool firstIter = true;
0290
0291   while (it != starEndIt || firstIter) {  // iterate over strips
0292
0293     firstIter = false;
0294     // New PAC
0295     if (curStripNo % logplaneSize == 0) {
0296       ++curPACno;
0297       curBegStripNo = curStripNo;
0298       RPCStripsRing::iterator plus8 = it;
0299       bool skipOccured = false;
0300       for (int i = 0; i < 7; ++i) {
0301         ++plus8;
0302         if (plus8 == this->end()) {
0303           plus8 = this->begin();
0304           skipOccured = true;
0305         }
0306       }
0307
0308       // calculate angle
0309       float phi = it->first;
0310       float phiP8 = plus8->first;
0311       if (skipOccured) {
0312         // phiP8 is negative
0313         // phi is positive
0314         // xcheck
0315         if (phi * phiP8 > 0) {
0316           throw cms::Exception("RPCInternal") << " RPCStripsRing::createRefConnections phi/phi8 error \n";
0317         }
0318         angle = (2 * pi + phiP8 + phi) / 2;
0319         if (angle > pi) {  // should land on positive side
0320           angle -= 2 * pi;
0321         }
0322
0323         if (std::abs(angle) > pi) {
0324           throw cms::Exception("RPCInternal") << " RPCStripsRing::createRefConnections "
0325                                               << " problem with angle calc \n";
0326         }
0327       } else {
0328         angle = (phiP8 + phi) / 2;
0329       }
0330       //std::cout << curPACno << " " << phiP8 << " " << phi << " "  << angle << std::endl;
0331
0332       TOtherConnStructVec::iterator itOt = otherRings.begin();
0333       for (; itOt != otherRings.end(); ++itOt) {
0334         itOt->m_it->second.createOtherConnections(
0335             getTowerForRefRing(), curPACno, itOt->m_logplane, itOt->m_logplaneSize, angle);
0336       }
0337     }
0338
0339     if (!it->second.isVirtual()) {
0340       L1RPCConeBuilder::TStripCon newCon;
0341       newCon.m_tower = getTowerForRefRing();
0342       newCon.m_PAC = curPACno;
0343       newCon.m_logplane = logplane;
0344       newCon.m_logstrip = curStripNo - curBegStripNo;
0345       //std::cout << " Adding con for " << it->second.m_detRawId << std::endl;
0346       (*m_connectionsMap)[it->second.m_detRawId][it->second.m_strip].push_back(newCon);
0347       //std::cout << " Adding ref connection " << std::endl;
0348     }
0349     ++curStripNo;
0350     ++it;
0351     if (it == this->end()) {
0352       it = this->begin();
0353     }
0354
0355   }  // iteration over strips ends
0356
0357   //std::cout << " refcon: " << curPACno << " PACs" << std::endl;
0358   //std::cout << "After refCon: " << m_connectionsMap.size() << std::endl;
0359 }
0360
0361 void RPCStripsRing::createOtherConnections(int tower, int PACno, int logplane, int logplaneSize, float angle) {
0362   //std::cout << "    OtherCon for " << getRingId() << std::endl;
0363
0364   if (this->isReferenceRing()) {
0365     throw cms::Exception("RPCInternal") << " RPCStripsRing::createOtherConnections "
0366                                         << " called for reference ring \n";
0367   }
0368
0369   RPCStripsRing::const_iterator it = this->lower_bound(angle);
0370
0371   if (it == this->end())
0372     it = this->begin();
0373
0374   for (int i = 0; i < logplaneSize / 2; i++) {
0375     if (it == this->begin())
0376       it = this->end();  // (m_stripPhiMap.end()--) is ok.
0377
0378     --it;
0379   }
0380
0381   for (int i = 0; i < logplaneSize; i++) {
0382     if (!it->second.isVirtual()) {
0383       L1RPCConeBuilder::TStripCon newCon;
0384       newCon.m_tower = tower;
0385       newCon.m_PAC = PACno;
0386       newCon.m_logplane = logplane;
0387       newCon.m_logstrip = i;
0388       (*m_connectionsMap)[it->second.m_detRawId][it->second.m_strip].push_back(newCon);
0389       //std::cout << " Adding other connection " << std::endl;
0390     }
0391
0392     ++it;
0393     if (it == this->end())
0394       it = this->begin();
0395   }
0396 }
0397
0398 // Defines to which tower this ring (only ref ring) belongs
0399 int RPCStripsRing::getTowerForRefRing() {
0400   int ret = 0;
0401
0402   if (!this->isReferenceRing()) {
0403     throw cms::Exception("RPCInternal") << " RPCStripsRing::getTowerForRefRing() "
0404                                         << " called for non reference ring \n";
0405   }
0406
0407   int etaAbs = std::abs(getEtaPartition());
0408   if (etaAbs < 8) {
0409     ret = getEtaPartition();
0410   } else if (etaAbs > 8) {
0411     int sign = (getEtaPartition() > 0 ? 1 : -1);
0412     ret = getEtaPartition() - sign;
0413   } else {
0414     throw cms::Exception("RPCInternal") << " RPCStripsRing::getTowerForRefRing() "
0415                                         << " called for etaPartition 8 \n";
0416   }
0417
0418   return ret;
0419 }
0420 /*
0421       struct TStripCon{
0422         signed char m_tower;
0423         unsigned char m_PAC;
0424         unsigned char m_logplane;
0425         unsigned char m_logstrip;
0426       };
0427       typedef std::vector<TStripCon> TStripConVec;
0428       typedef std::map<unsigned char, TStripConVec> TStrip2ConVec;
0429       typedef std::map<uint32_t, TStrip2ConVec> TConMap;
0430
0431       // compressed connections
0432       struct TCompressedCon{
0433         signed char m_tower;
0434         unsigned char m_PAC;
0435         signed char m_offset;
0436         signed char m_mul;
0437       };
0438       typedef std::vector<TCompressedCon> TCompressedConVec;
0439       typedef std::map<uint32_t, TCompressedConVec> TCompressedConMap;
0440
0441 */
0442
0443 void RPCStripsRing::compressConnections() {
0444   L1RPCConeBuilder::TConMap::iterator itChamber = m_connectionsMap->begin();
0445
0446   auto uncompressedConsLeft = std::make_shared<L1RPCConeBuilder::TConMap>();
0447
0448   m_compressedConnectionMap = std::make_shared<L1RPCConeBuilder::TCompressedConMap>();
0449
0450   int compressedCons = 0, uncompressedConsBefore = 0, uncompressedConsAfter = 0;
0451
0452   //   int offsetMin =0, offsetMax =0;
0453
0454   for (; itChamber != m_connectionsMap->end(); ++itChamber) {
0455     uint32_t detId = itChamber->first;
0456
0457     for (L1RPCConeBuilder::TStrip2ConVec::iterator itStrip = itChamber->second.begin();
0458          itStrip != itChamber->second.end();
0459          ++itStrip) {
0460       // Iterate over strip Connections
0461       for (L1RPCConeBuilder::TStripConVec::iterator itConn = itStrip->second.begin(); itConn != itStrip->second.end();
0462            ++itConn) {
0463         // Check if this connection isn't allready present in the compressed map
0464         ++uncompressedConsBefore;
0465         bool alreadyDone = false;
0466         if (m_compressedConnectionMap->find(detId) != m_compressedConnectionMap->end()) {
0467           // iterate over the vec, check element by element
0468           for (L1RPCConeBuilder::TCompressedConVec::iterator itCompConn = (*m_compressedConnectionMap)[detId].begin();
0469                itCompConn != (*m_compressedConnectionMap)[detId].end();
0470                ++itCompConn) {
0471             if (itCompConn->m_tower == itConn->m_tower && itCompConn->m_PAC == itConn->m_PAC &&
0472                 itCompConn->m_logplane == itConn->m_logplane)  // connection allready compressed
0473             {
0474               alreadyDone = true;
0475
0476               int logStrip = itCompConn->m_mul * itStrip->first + itCompConn->m_offset;
0477               if (logStrip != itConn->m_logstrip) {
0478                 //copy the problematic connection to the "safe" map
0479                 (*uncompressedConsLeft)[detId][itStrip->first].push_back(*itConn);
0480                 ++uncompressedConsAfter;
0481                 edm::LogWarning("RPCTriggerConfig")
0482                     << " Compression failed for det " << detId << " strip " << (int)itStrip->first << " . Got "
0483                     << (int)logStrip << " expected " << (int)itConn->m_logstrip << std::endl;
0484               } else {
0485                 itCompConn->addStrip(itStrip->first);
0486               }
0487             }
0488           }  // compressed connection iteration end
0489         }
0490         //if (detId==637569977) std::cout << " Buld cons for strip " << (int)itStrip->first << std::endl;
0491
0492         if (!alreadyDone) {
0493           // find another strip contributing to the same PAC,tower,logplane
0494           L1RPCConeBuilder::TStrip2ConVec::iterator itStripOther = itStrip;
0495           ++itStripOther;
0496           bool otherStripFound = false;
0497           signed char mul = 1;
0498           for (; itStripOther != itChamber->second.end() && !otherStripFound; ++itStripOther) {
0499             for (L1RPCConeBuilder::TStripConVec::iterator itConnOther = itStripOther->second.begin();
0500                  itConnOther != itStripOther->second.end();
0501                  ++itConnOther) {
0502               if (itConnOther->m_tower == itConn->m_tower && itConnOther->m_PAC == itConn->m_PAC &&
0503                   itConnOther->m_logplane == itConn->m_logplane)  // connection to same PAC,logplane
0504               {
0505                 otherStripFound = true;
0506                 if ((itStripOther->first - itStrip->first) * (itConnOther->m_logstrip - itConn->m_logstrip) < 0) {
0507                   mul = -1;
0508                 }
0509                 break;
0510               }
0511             }  // otherConnections iter ends
0512           }    // otherStrip iter ends
0513
0514           /*
0515           if (itConn->m_tower==3 && itConn->m_PAC==73 && itConn->m_logplane==4 && detId==637569977){
0516             std::cout << " Buld cons for strip " << (int)itStrip->first;
0517             if (otherStripFound)
0518               std::cout << " other strip " << itStrip->first;
0519             else
0520               std::cout << " no other strip ";
0521
0522             std::cout << std::endl;
0523
0524         }*/
0525
0526           L1RPCConeBuilder::TCompressedCon nCompConn;
0527           nCompConn.m_tower = itConn->m_tower;
0528           nCompConn.m_PAC = itConn->m_PAC;
0529           nCompConn.m_logplane = itConn->m_logplane;
0530           nCompConn.m_mul = mul;
0531           nCompConn.m_offset = itConn->m_logstrip - mul * (signed short)(itStrip->first);
0532           nCompConn.addStrip(itStrip->first);
0533
0534           if (otherStripFound) {
0535           } else {
0536             //  uncompressedConsLeft[detId][itStrip->first].push_back(*itConn);
0537             //  ++uncompressedConsAfter;
0538           }
0539           (*m_compressedConnectionMap)[detId].push_back(nCompConn);
0540           ++compressedCons;
0541
0542         }  // if(!allreadyDone)
0543       }    // iterate on connections
0544     }      // iterate on strips
0545   }        // iterate on chambers
0546
0547   // 159 -87
0548   //std::cout << offsetMax << " TT " << offsetMin << std::endl;
0549
0550   edm::LogInfo("RPCTriggerConfig") << " Compressed: " << compressedCons << " "
0551                                    << sizeof(L1RPCConeBuilder::TCompressedCon)
0552                                    << " Uncompressed before: " << uncompressedConsBefore << " "
0553                                    << sizeof(L1RPCConeBuilder::TStripCon)
0554                                    << " Uncompressed after: " << uncompressedConsAfter << " "
0555                                    << sizeof(L1RPCConeBuilder::TStripCon);
0556   m_connectionsMap = uncompressedConsLeft;
0557 }