Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2022-06-30 23:16:56

0001 #include "DQM/GEM/interface/GEMDQMBase.h"
0002 #include "Geometry/CommonTopologies/interface/GEMStripTopology.h"
0003 
0004 using namespace std;
0005 using namespace edm;
0006 
0007 GEMDQMBase::GEMDQMBase(const edm::ParameterSet& cfg) : geomToken_(esConsumes<edm::Transition::BeginRun>()) {
0008   std::string strRunType = cfg.getUntrackedParameter<std::string>("runType");
0009 
0010   nRunType_ = GEMDQM_RUNTYPE_ONLINE;
0011 
0012   if (strRunType == "online") {
0013     nRunType_ = GEMDQM_RUNTYPE_ONLINE;
0014   } else if (strRunType == "offline") {
0015     nRunType_ = GEMDQM_RUNTYPE_OFFLINE;
0016   } else if (strRunType == "relval") {
0017     nRunType_ = GEMDQM_RUNTYPE_RELVAL;
0018   } else if (strRunType == "allplots") {
0019     nRunType_ = GEMDQM_RUNTYPE_ALLPLOTS;
0020   } else {
0021     edm::LogError(log_category_) << "+++ Error : GEM geometry is unavailable on event loop. +++\n";
0022   }
0023 
0024   log_category_ = cfg.getUntrackedParameter<std::string>("logCategory");
0025 }
0026 
0027 int GEMDQMBase::initGeometry(edm::EventSetup const& iSetup) {
0028   GEMGeometry_ = nullptr;
0029   try {
0030     //edm::ESHandle<GEMGeometry> hGeom;
0031     //iSetup.get<MuonGeometryRecord>().get(hGeom);
0032     GEMGeometry_ = &iSetup.getData(geomToken_);
0033   } catch (edm::eventsetup::NoProxyException<GEMGeometry>& e) {
0034     edm::LogError(log_category_) << "+++ Error : GEM geometry is unavailable on event loop. +++\n";
0035     return -1;
0036   }
0037 
0038   return 0;
0039 }
0040 
0041 // Borrowed from DQM/GEM/src/GEMOfflineDQMBase.cc
0042 int GEMDQMBase::getNumEtaPartitions(const GEMStation* station) {
0043   const auto&& superchambers = station->superChambers();
0044   if (not checkRefs(superchambers)) {
0045     edm::LogError(log_category_) << "failed to get a valid vector of GEMSuperChamber ptrs" << std::endl;
0046     return 0;
0047   }
0048 
0049   const auto& chambers = superchambers.front()->chambers();
0050   if (not checkRefs(chambers)) {
0051     edm::LogError(log_category_) << "failed to get a valid vector of GEMChamber ptrs" << std::endl;
0052     return 0;
0053   }
0054 
0055   return chambers.front()->nEtaPartitions();
0056 }
0057 
0058 int GEMDQMBase::loadChambers() {
0059   if (GEMGeometry_ == nullptr)
0060     return -1;
0061   gemChambers_.clear();
0062   const std::vector<const GEMSuperChamber*>& superChambers_ = GEMGeometry_->superChambers();
0063   for (auto sch : superChambers_) {  // FIXME: This loop can be merged into the below loop
0064     for (auto pch : sch->chambers()) {
0065       Bool_t bExist = false;
0066       for (const auto& ch : gemChambers_) {
0067         if (pch->id() == ch.id()) {
0068           bExist = true;
0069           break;
0070         }
0071       }
0072       if (bExist)
0073         continue;
0074       gemChambers_.push_back(*pch);
0075     }
0076   }
0077 
0078   // Borrwed from DQM/GEM/src/GEMOfflineMonitor.cc
0079   nMaxNumCh_ = 0;
0080   for (const GEMRegion* region : GEMGeometry_->regions()) {
0081     const int region_number = region->region();
0082 
0083     for (const GEMStation* station : region->stations()) {
0084       const auto&& superchambers = station->superChambers();
0085 
0086       const int station_number = station->station();
0087       const int num_superchambers = (station_number == 1 ? 36 : 18);
0088       const int max_vfat = getMaxVFAT(station->station());  // the number of VFATs per GEMEtaPartition
0089       const int num_etas = getNumEtaPartitions(station);    // the number of eta partitions per GEMChamber
0090       const int num_vfat = num_etas * max_vfat;             // the number of VFATs per GEMChamber
0091       const int num_digi = GEMeMap::maxChan_;               // the number of digis (channels) per VFAT
0092 
0093       nMaxNumCh_ = std::max(nMaxNumCh_, num_superchambers);
0094 
0095       Int_t nMinIdxChamber = 1048576;
0096       Int_t nMaxIdxChamber = -1048576;
0097       for (auto sch : superchambers) {
0098         auto nIdxChamber = sch->chambers().front()->id().chamber();
0099         if (nMinIdxChamber > nIdxChamber)
0100           nMinIdxChamber = nIdxChamber;
0101         if (nMaxIdxChamber < nIdxChamber)
0102           nMaxIdxChamber = nIdxChamber;
0103       }
0104 
0105       const auto& chambers = superchambers.front()->chambers();
0106 
0107       for (auto pchamber : chambers) {
0108         int layer_number = pchamber->id().layer();
0109         ME3IdsKey key3(region_number, station_number, layer_number);
0110         mapStationInfo_[key3] = MEStationInfo(region_number,
0111                                               station_number,
0112                                               layer_number,
0113                                               num_superchambers,
0114                                               num_etas,
0115                                               num_vfat,
0116                                               num_digi,
0117                                               nMinIdxChamber,
0118                                               nMaxIdxChamber);
0119         readGeometryRadiusInfoChamber(station, mapStationInfo_[key3]);
0120         readGeometryPhiInfoChamber(station, mapStationInfo_[key3]);
0121       }
0122     }
0123   }
0124 
0125   return 0;
0126 }
0127 
0128 int GEMDQMBase::SortingLayers(std::vector<ME3IdsKey>& listLayers) {
0129   std::sort(listLayers.begin(), listLayers.end(), [](ME3IdsKey key1, ME3IdsKey key2) {
0130     Int_t re1 = std::get<0>(key1), st1 = std::get<1>(key1), la1 = std::get<2>(key1);
0131     Int_t re2 = std::get<0>(key2), st2 = std::get<1>(key2), la2 = std::get<2>(key2);
0132     if (re1 < 0 && re2 > 0)
0133       return false;
0134     if (re1 > 0 && re2 < 0)
0135       return true;
0136     Bool_t bRes = (re1 < 0);  // == re2 < 0
0137     Int_t sum1 = 256 * std::abs(re1) + 16 * st1 + 1 * la1;
0138     Int_t sum2 = 256 * std::abs(re2) + 16 * st2 + 1 * la2;
0139     if (sum1 <= sum2)
0140       return bRes;
0141     return !bRes;
0142   });
0143 
0144   return 0;
0145 }
0146 
0147 dqm::impl::MonitorElement* GEMDQMBase::CreateSummaryHist(DQMStore::IBooker& ibooker, TString strName) {
0148   std::vector<ME3IdsKey> listLayers;
0149   for (auto const& [key, stationInfo] : mapStationInfo_)
0150     listLayers.push_back(key);
0151   SortingLayers(listLayers);
0152   for (Int_t i = 0; i < (Int_t)listLayers.size(); i++)
0153     mapStationToIdx_[listLayers[i]] = i + 1;
0154 
0155   auto h2Res =
0156       ibooker.book2D(strName, "", nMaxNumCh_, 0.5, nMaxNumCh_ + 0.5, listLayers.size(), 0.5, listLayers.size() + 0.5);
0157   h2Res->setXTitle("Chamber");
0158   h2Res->setYTitle("Layer");
0159 
0160   if (h2Res == nullptr)
0161     return nullptr;
0162 
0163   for (Int_t i = 1; i <= nMaxNumCh_; i++)
0164     h2Res->setBinLabel(i, Form("%i", i), 1);
0165   for (Int_t i = 1; i <= (Int_t)listLayers.size(); i++) {
0166     auto key = listLayers[i - 1];
0167     auto strInfo = GEMUtils::getSuffixName(key);  // NOTE: It starts with '_'
0168     auto region = keyToRegion(key);
0169     auto label =
0170         Form("GE%+i1-%cL%i;%s", region * keyToStation(key), (region > 0 ? 'P' : 'M'), keyToLayer(key), strInfo.Data());
0171     h2Res->setBinLabel(i, label, 2);
0172     Int_t nNumCh = mapStationInfo_[key].nNumChambers_;
0173     h2Res->setBinContent(0, i, nNumCh);
0174   }
0175 
0176   return h2Res;
0177 }
0178 
0179 int GEMDQMBase::GenerateMEPerChamber(DQMStore::IBooker& ibooker) {
0180   MEMap2Check_.clear();
0181   MEMap2WithEtaCheck_.clear();
0182   MEMap2AbsReWithEtaCheck_.clear();
0183   MEMap3Check_.clear();
0184   MEMap3WithChCheck_.clear();
0185   MEMap4Check_.clear();
0186   for (const auto& ch : gemChambers_) {
0187     GEMDetId gid = ch.id();
0188     ME2IdsKey key2{gid.region(), gid.station()};
0189     ME3IdsKey key3{gid.region(), gid.station(), gid.layer()};
0190     ME4IdsKey key3WithChamber{gid.region(), gid.station(), gid.layer(), gid.chamber()};
0191     if (!MEMap2Check_[key2]) {
0192       auto strSuffixName = GEMUtils::getSuffixName(key2);
0193       auto strSuffixTitle = GEMUtils::getSuffixTitle(key2);
0194       BookingHelper bh2(ibooker, strSuffixName, strSuffixTitle);
0195       ProcessWithMEMap2(bh2, key2);
0196       MEMap2Check_[key2] = true;
0197     }
0198     if (!MEMap3Check_[key3]) {
0199       auto strSuffixName = GEMUtils::getSuffixName(key3);
0200       auto strSuffixTitle = GEMUtils::getSuffixTitle(key3);
0201       BookingHelper bh3(ibooker, strSuffixName, strSuffixTitle);
0202       ProcessWithMEMap3(bh3, key3);
0203       MEMap3Check_[key3] = true;
0204     }
0205     if (!MEMap3WithChCheck_[key3WithChamber]) {
0206       Int_t nCh = gid.chamber();
0207       Int_t nLa = gid.layer();
0208       char cLS = (nCh % 2 == 0 ? 'L' : 'S');  // FIXME: Is it general enough?
0209       auto strSuffixName = GEMUtils::getSuffixName(key2) + Form("-%02iL%i-%c", nCh, nLa, cLS);
0210       auto strSuffixTitle = GEMUtils::getSuffixTitle(key2) + Form("-%02iL%i-%c", nCh, nLa, cLS);
0211       BookingHelper bh3Ch(ibooker, strSuffixName, strSuffixTitle);
0212       ProcessWithMEMap3WithChamber(bh3Ch, key3WithChamber);
0213       MEMap3WithChCheck_[key3WithChamber] = true;
0214     }
0215     for (auto iEta : ch.etaPartitions()) {
0216       GEMDetId eId = iEta->id();
0217       ME4IdsKey key4{gid.region(), gid.station(), gid.layer(), eId.ieta()};
0218       ME3IdsKey key2WithEta{gid.region(), gid.station(), eId.ieta()};
0219       ME3IdsKey key2AbsReWithEta{std::abs(gid.region()), gid.station(), eId.ieta()};
0220       if (!MEMap4Check_[key4]) {
0221         auto strSuffixName = GEMUtils::getSuffixName(key3) + Form("-E%02i", eId.ieta());
0222         auto strSuffixTitle = GEMUtils::getSuffixTitle(key3) + Form("-E%02i", eId.ieta());
0223         BookingHelper bh4(ibooker, strSuffixName, strSuffixTitle);
0224         ProcessWithMEMap4(bh4, key4);
0225         MEMap4Check_[key4] = true;
0226       }
0227       if (!MEMap2WithEtaCheck_[key2WithEta]) {
0228         auto strSuffixName = GEMUtils::getSuffixName(key2) + Form("-E%02i", eId.ieta());
0229         auto strSuffixTitle = GEMUtils::getSuffixTitle(key2) + Form("-E%02i", eId.ieta());
0230         BookingHelper bh3(ibooker, strSuffixName, strSuffixTitle);
0231         ProcessWithMEMap2WithEta(bh3, key2WithEta);
0232         MEMap2WithEtaCheck_[key2WithEta] = true;
0233       }
0234       if (!MEMap2AbsReWithEtaCheck_[key2AbsReWithEta]) {
0235         auto strSuffixName = Form("_GE%d1-E%02i", gid.station(), eId.ieta());
0236         auto strSuffixTitle = Form(" GE%d1-E%02i", gid.station(), eId.ieta());
0237         BookingHelper bh3(ibooker, strSuffixName, strSuffixTitle);
0238         ProcessWithMEMap2AbsReWithEta(bh3, key2AbsReWithEta);
0239         MEMap2AbsReWithEtaCheck_[key2AbsReWithEta] = true;
0240       }
0241     }
0242   }
0243   return 0;
0244 }
0245 
0246 int GEMDQMBase::readGeometryRadiusInfoChamber(const GEMStation* station, MEStationInfo& stationInfo) {
0247   auto listSuperChambers = station->superChambers();
0248 
0249   Bool_t bDoneEven = false, bDoneOdd = false;
0250 
0251   // Obtaining radius intervals of even/odd chambers
0252   for (auto superchamber : listSuperChambers) {
0253     Int_t chamberNo = superchamber->id().chamber();
0254     if (chamberNo % 2 == 0 && bDoneEven)
0255       continue;
0256     if (chamberNo % 2 != 0 && bDoneOdd)
0257       continue;
0258 
0259     auto& etaPartitions = superchamber->chambers().front()->etaPartitions();
0260 
0261     // A little of additional procedures to list up the radius intervals
0262     // It would be independent to local direction of chambers and the order of eta partitions
0263     //   1. Obtain the radius of the middle top/bottom points of the trapezoid
0264     //   2. Sort these two values and determine which one is the lower/upper one
0265     //   3. Keep them all and then sort them
0266     //   4. The intermediate radii are set as the mean of the corresponding values of upper/lowers.
0267     std::vector<Float_t> listRadiusLower, listRadiusUpper;
0268     for (auto iEta : etaPartitions) {
0269       const GEMStripTopology& stripTopology = dynamic_cast<const GEMStripTopology&>(iEta->specificTopology());
0270       Float_t fHeight = stripTopology.stripLength();
0271       LocalPoint lp1(0.0, -0.5 * fHeight), lp2(0.0, 0.5 * fHeight);
0272       auto& surface = iEta->surface();
0273       GlobalPoint gp1 = surface.toGlobal(lp1), gp2 = surface.toGlobal(lp2);
0274       Float_t fR1 = gp1.perp(), fR2 = gp2.perp();
0275       Float_t fRL = std::min(fR1, fR2), fRH = std::max(fR1, fR2);
0276       listRadiusLower.push_back(fRL);
0277       listRadiusUpper.push_back(fRH);
0278       // For a future usage
0279       //std::cout << "GEO_RADIUS: " << iEta->id().chamber() << " " << iEta->id().ieta() << " "
0280       //  << fRL << " " << fRH << std::endl;
0281     }
0282 
0283     std::sort(listRadiusLower.begin(), listRadiusLower.end());
0284     std::sort(listRadiusUpper.begin(), listRadiusUpper.end());
0285 
0286     std::vector<Float_t>& listR =
0287         (chamberNo % 2 == 0 ? stationInfo.listRadiusEvenChamber_ : stationInfo.listRadiusOddChamber_);
0288     listR.clear();
0289     listR.push_back(listRadiusLower.front());
0290     for (int i = 1; i < (int)listRadiusLower.size(); i++) {
0291       listR.push_back(0.5 * (listRadiusLower[i] + listRadiusUpper[i - 1]));
0292     }
0293     listR.push_back(listRadiusUpper.back());
0294 
0295     if (chamberNo % 2 == 0)
0296       bDoneEven = true;
0297     if (chamberNo % 2 != 0)
0298       bDoneOdd = true;
0299 
0300     if (bDoneEven && bDoneOdd)
0301       break;
0302   }
0303 
0304   return 0;
0305 }
0306 
0307 int GEMDQMBase::readGeometryPhiInfoChamber(const GEMStation* station, MEStationInfo& stationInfo) {
0308   auto listSuperChambers = station->superChambers();
0309   Int_t nNumStripEta = stationInfo.nNumDigi_ * (stationInfo.nMaxVFAT_ / stationInfo.nNumEtaPartitions_);
0310 
0311   std::vector<std::pair<Int_t, std::pair<std::pair<Float_t, Float_t>, Bool_t>>> listDivPhi;
0312 
0313   // Obtaining phi intervals of chambers
0314   for (auto superchamber : listSuperChambers) {
0315     auto iEta = superchamber->chambers().front()->etaPartitions().front();
0316 
0317     // What is the index of the first strip? Rather than to ask to someone, let's calculate it!
0318     Float_t fWidthStrip = std::abs(iEta->centreOfStrip((Int_t)1).x() - iEta->centreOfStrip((Int_t)0).x());
0319     LocalPoint lpRef(-fWidthStrip / 3.0, 0.0);
0320     Int_t nStripMid = (Int_t)iEta->strip(lpRef);
0321     Int_t nFirstStrip = 1 - ((nNumStripEta / 2) - nStripMid);
0322     Int_t nLastStrip = nFirstStrip + nNumStripEta - 1;
0323 
0324     auto& surface = iEta->surface();
0325     LocalPoint lpF = iEta->centreOfStrip((Float_t)(nFirstStrip - 0.5));  // To avoid the round error(?)
0326     LocalPoint lpL = iEta->centreOfStrip((Float_t)(nLastStrip + 0.5));   // To avoid the round error(?)
0327     GlobalPoint gpF = surface.toGlobal(lpF);
0328     GlobalPoint gpL = surface.toGlobal(lpL);
0329 
0330     Float_t fPhiF = gpF.phi();
0331     Float_t fPhiL = gpL.phi();
0332     if (fPhiF * fPhiL < 0 && std::abs(fPhiF) > 0.5 * 3.14159265359) {
0333       if (fPhiF < 0)
0334         fPhiF += 2 * 3.14159265359;
0335       if (fPhiL < 0)
0336         fPhiL += 2 * 3.14159265359;
0337     }
0338     Bool_t bFlipped = fPhiF > fPhiL;
0339     Float_t fPhiMin = std::min(fPhiF, fPhiL);
0340     Float_t fPhiMax = std::max(fPhiF, fPhiL);
0341 
0342     listDivPhi.emplace_back();
0343     listDivPhi.back().first = iEta->id().chamber();
0344     listDivPhi.back().second.first.first = fPhiMin;
0345     listDivPhi.back().second.first.second = fPhiMax;
0346     listDivPhi.back().second.second = bFlipped;
0347   }
0348 
0349   stationInfo.fMinPhi_ = 0.0;
0350   for (auto p : listDivPhi) {
0351     if (p.first == 1) {
0352       stationInfo.fMinPhi_ = p.second.first.first;
0353       break;
0354     }
0355   }
0356 
0357   // For a future usage
0358   //for ( auto p : listDivPhi ) {
0359   //  std::cout << "GEO_PHI: " << p.first << " "
0360   //    << p.second.first.first << " " << p.second.first.second << " " << p.second.second << std::endl;
0361   //}
0362 
0363   return 0;
0364 }