DTMonitorClient/src/DTOccupancyTest.cc

0001 /*
0002  *  See header file for a description of this class.
0003  *
0004  *  \author G. Cerminara - University and INFN Torino
0005  *
0006  *  threadsafe version (//-) oct/nov 2014 - WATWanAbdullah -ncpp-um-my
0007  *
0008  */
0009
0010 #include "DQM/DTMonitorClient/src/DTOccupancyTest.h"
0011
0012 #include "FWCore/ServiceRegistry/interface/Service.h"
0013 #include "FWCore/Framework/interface/LuminosityBlock.h"
0014 #include "FWCore/Framework/interface/EventSetup.h"
0015 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0016 #include "Geometry/DTGeometry/interface/DTGeometry.h"
0017 #include "DQMServices/Core/interface/DQMStore.h"
0018 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0019
0020 #include "TMath.h"
0021
0022 using namespace edm;
0023 using namespace std;
0024
0025 DTOccupancyTest::DTOccupancyTest(const edm::ParameterSet& ps)
0026     : muonGeomToken_(esConsumes<edm::Transition::BeginRun>()) {
0027   LogVerbatim("DTDQM|DTMonitorClient|DTOccupancyTest") << "[DTOccupancyTest]: Constructor";
0028
0029   // Get the DQM service
0030
0031   lsCounter = 0;
0032
0033   writeRootFile = ps.getUntrackedParameter<bool>("writeRootFile", false);
0034   if (writeRootFile) {
0035     rootFile = new TFile("DTOccupancyTest.root", "RECREATE");
0036     ntuple = new TNtuple("OccupancyNtuple",
0037                          "OccupancyNtuple",
0038                          "ls:wh:st:se:lay1MeanCell:lay1RMS:lay2MeanCell:lay2RMS:lay3MeanCell:lay3RMS:lay4MeanCell:"
0039                          "lay4RMS:lay5MeanCell:lay5RMS:lay6MeanCell:lay6RMS:lay7MeanCell:lay7RMS:lay8MeanCell:lay8RMS:"
0040                          "lay9MeanCell:lay9RMS:lay10MeanCell:lay10RMS:lay11MeanCell:lay11RMS:lay12MeanCell:lay12RMS");
0041   }
0042
0043   // switch on the mode for running on test pulses (different top folder)
0044   tpMode = ps.getUntrackedParameter<bool>("testPulseMode", false);
0045
0046   runOnAllHitsOccupancies = ps.getUntrackedParameter<bool>("runOnAllHitsOccupancies", true);
0047   runOnNoiseOccupancies = ps.getUntrackedParameter<bool>("runOnNoiseOccupancies", false);
0048   runOnInTimeOccupancies = ps.getUntrackedParameter<bool>("runOnInTimeOccupancies", false);
0049   nMinEvts = ps.getUntrackedParameter<int>("nEventsCert", 5000);
0050   nMinEvtsPC = ps.getUntrackedParameter<int>("nEventsMinPC", 2200);
0051   nZeroEvtsPC = ps.getUntrackedParameter<int>("nEventsZeroPC", 30);
0052
0053   bookingdone = false;
0054
0055   // Event counter
0056   nevents = 0;
0057 }
0058
0059 DTOccupancyTest::~DTOccupancyTest() {
0060   LogVerbatim("DTDQM|DTMonitorClient|DTOccupancyTest") << " destructor called" << endl;
0061 }
0062
0063 void DTOccupancyTest::beginRun(const edm::Run& run, const EventSetup& context) {
0064   LogVerbatim("DTDQM|DTMonitorClient|DTOccupancyTest") << "[DTOccupancyTest]: BeginRun";
0065
0066   // Get the geometry
0067   muonGeom = &context.getData(muonGeomToken_);
0068 }
0069
0070 void DTOccupancyTest::dqmEndLuminosityBlock(DQMStore::IBooker& ibooker,
0071                                             DQMStore::IGetter& igetter,
0072                                             edm::LuminosityBlock const& lumiSeg,
0073                                             edm::EventSetup const& context) {
0074   if (!bookingdone) {
0075     // Book the summary histos
0076     //   - one summary per wheel
0077     for (int wh = -2; wh <= 2; ++wh) {  // loop over wheels
0078       bookHistos(ibooker, wh, string("Occupancies"), "OccupancySummary");
0079     }
0080
0081     ibooker.setCurrentFolder(topFolder());
0082     string title = "Occupancy Summary";
0083     if (tpMode) {
0084       title = "Test Pulse Occupancy Summary";
0085     }
0086     //   - global summary with alarms
0087     summaryHisto = ibooker.book2D("OccupancySummary", title.c_str(), 12, 1, 13, 5, -2, 3);
0088     summaryHisto->setAxisTitle("sector", 1);
0089     summaryHisto->setAxisTitle("wheel", 2);
0090
0091     //   - global summary with percentages
0092     glbSummaryHisto = ibooker.book2D("OccupancyGlbSummary", title.c_str(), 12, 1, 13, 5, -2, 3);
0093     glbSummaryHisto->setAxisTitle("sector", 1);
0094     glbSummaryHisto->setAxisTitle("wheel", 2);
0095
0096     // assign the name of the input histogram
0097     if (runOnAllHitsOccupancies) {
0098       nameMonitoredHisto = "OccupancyAllHits_perCh";
0099     } else if (runOnNoiseOccupancies) {
0100       nameMonitoredHisto = "OccupancyNoise_perCh";
0101     } else if (runOnInTimeOccupancies) {
0102       nameMonitoredHisto = "OccupancyInTimeHits_perCh";
0103     } else {  // default is AllHits histo
0104       nameMonitoredHisto = "OccupancyAllHits_perCh";
0105     }
0106   }
0107   bookingdone = true;
0108
0109   LogVerbatim("DTDQM|DTMonitorClient|DTOccupancyTest")
0110       << "[DTOccupancyTest]: End of LS transition, performing the DQM client operation";
0111   lsCounter++;
0112
0113   // Reset the global summary
0114   summaryHisto->Reset();
0115   glbSummaryHisto->Reset();
0116
0117   nChannelTotal = 0;
0118   nChannelDead = 0;
0119
0120   // Get all the DT chambers
0121   vector<const DTChamber*> chambers = muonGeom->chambers();
0122
0123   for (vector<const DTChamber*>::const_iterator chamber = chambers.begin(); chamber != chambers.end();
0124        ++chamber) {  // Loop over all chambers
0125     DTChamberId chId = (*chamber)->id();
0126
0127     MonitorElement* chamberOccupancyHisto = igetter.get(getMEName(nameMonitoredHisto, chId));
0128
0129     // Run the tests on the plot for the various granularities
0130     if (chamberOccupancyHisto != nullptr) {
0131       // Get the 2D histo
0132       TH2F* histo = chamberOccupancyHisto->getTH2F();
0133       float chamberPercentage = 1.;
0134       int result = runOccupancyTest(histo, chId, chamberPercentage);
0135       int sector = chId.sector();
0136
0137       if (sector == 13) {
0138         sector = 4;
0139         float resultSect4 = wheelHistos[chId.wheel()]->getBinContent(sector, chId.station());
0140         if (resultSect4 > result) {
0141           result = (int)resultSect4;
0142         }
0143       } else if (sector == 14) {
0144         sector = 10;
0145         float resultSect10 = wheelHistos[chId.wheel()]->getBinContent(sector, chId.station());
0146         if (resultSect10 > result) {
0147           result = (int)resultSect10;
0148         }
0149       }
0150
0151       // the 2 MB4 of Sect 4 and 10 count as half a chamber
0152       if ((sector == 4 || sector == 10) && chId.station() == 4)
0153         chamberPercentage = chamberPercentage / 2.;
0154
0155       wheelHistos[chId.wheel()]->setBinContent(sector, chId.station(), result);
0156       if (result > summaryHisto->getBinContent(sector, chId.wheel() + 3)) {
0157         summaryHisto->setBinContent(sector, chId.wheel() + 3, result);
0158       }
0159       glbSummaryHisto->Fill(sector, chId.wheel(), chamberPercentage * 1. / 4.);
0160     } else {
0161       LogVerbatim("DTDQM|DTMonitorClient|DTOccupancyTest")
0162           << "[DTOccupancyTest] ME: " << getMEName(nameMonitoredHisto, chId) << " not found!" << endl;
0163     }
0164   }
0165
0166   string nEvtsName = "DT/EventInfo/Counters/nProcessedEventsDigi";
0167
0168   MonitorElement* meProcEvts = igetter.get(nEvtsName);
0169
0170   if (meProcEvts) {
0171     int nProcEvts = meProcEvts->getFloatValue();
0172     glbSummaryHisto->setEntries(nProcEvts < nMinEvts ? 10. : nProcEvts);
0173     summaryHisto->setEntries(nProcEvts < nMinEvts ? 10. : nProcEvts);
0174   } else {
0175     glbSummaryHisto->setEntries(nMinEvts + 1);
0176     summaryHisto->setEntries(nMinEvts + 1);
0177     LogVerbatim("DTDQM|DTMonitorClient|DTOccupancyTest")
0178         << "[DTOccupancyTest] ME: " << nEvtsName << " not found!" << endl;
0179   }
0180
0181   // Fill the global summary
0182   // Check for entire sectors off and report them on the global summary
0183   //FIXME: TODO
0184
0185   if (writeRootFile)
0186     ntuple->AutoSave("SaveSelf");
0187 }
0188
0189 void DTOccupancyTest::dqmEndJob(DQMStore::IBooker& ibooker, DQMStore::IGetter& igetter) {
0190   LogVerbatim("DTDQM|DTMonitorClient|DTOccupancyTest") << "[DTOccupancyTest] endjob called!";
0191   if (writeRootFile) {
0192     rootFile->cd();
0193     ntuple->Write();
0194     rootFile->Close();
0195   }
0196 }
0197
0198 // --------------------------------------------------
0199
0200 void DTOccupancyTest::bookHistos(DQMStore::IBooker& ibooker, const int wheelId, string folder, string histoTag) {
0201   // Set the current folder
0202   stringstream wheel;
0203   wheel << wheelId;
0204
0205   ibooker.setCurrentFolder(topFolder());
0206
0207   // build the histo name
0208   string histoName = histoTag + "_W" + wheel.str();
0209
0210   LogVerbatim("DTDQM|DTMonitorClient|DTOccupancyTest")
0211       << "[DTOccupancyTest]: booking wheel histo:" << histoName << " (tag " << histoTag
0212       << ") in: " << topFolder() + "Wheel" + wheel.str() + "/" + folder << endl;
0213
0214   string histoTitle = "Occupancy summary WHEEL: " + wheel.str();
0215   if (tpMode) {
0216     histoTitle = "TP Occupancy summary WHEEL: " + wheel.str();
0217   }
0218
0219   wheelHistos[wheelId] = ibooker.book2D(histoName, histoTitle, 12, 1, 13, 4, 1, 5);
0220   wheelHistos[wheelId]->setBinLabel(1, "MB1", 2);
0221   wheelHistos[wheelId]->setBinLabel(2, "MB2", 2);
0222   wheelHistos[wheelId]->setBinLabel(3, "MB3", 2);
0223   wheelHistos[wheelId]->setBinLabel(4, "MB4", 2);
0224   wheelHistos[wheelId]->setAxisTitle("sector", 1);
0225 }
0226
0227 string DTOccupancyTest::getMEName(string histoTag, const DTChamberId& chId) {
0228   stringstream wheel;
0229   wheel << chId.wheel();
0230   stringstream station;
0231   station << chId.station();
0232   stringstream sector;
0233   sector << chId.sector();
0234
0235   string folderRoot = topFolder() + "Wheel" + wheel.str() + "/Sector" + sector.str() + "/Station" + station.str() + "/";
0236
0237   string folder = "Occupancies/";
0238
0239   // build the histo name
0240   string histoName = histoTag + "_W" + wheel.str() + "_St" + station.str() + "_Sec" + sector.str();
0241
0242   string histoname = folderRoot + histoName;
0243
0244   return histoname;
0245 }
0246
0247 int DTOccupancyTest::getIntegral(TH2F* histo, int firstBinX, int lastBinX, int firstBinY, int lastBinY, bool doall) {
0248   int sum = 0;
0249   for (Int_t i = firstBinX; i < lastBinX + 1; i++) {
0250     for (Int_t j = firstBinY; j < lastBinY + 1; j++) {
0251       if (histo->GetBinContent(i, j) > 0) {
0252         if (!doall)
0253           return 1;
0254         sum += histo->GetBinContent(i, j);
0255       }
0256     }
0257   }
0258
0259   return sum;
0260 }
0261
0262 // Run a test on the occupancy of the chamber
0263 // Return values:
0264 // 0 -> all ok
0265 // 1 -> # consecutive dead channels > N
0266 // 2 -> dead layer
0267 // 3 -> dead SL
0268 // 4 -> dead chamber
0269
0270 int DTOccupancyTest::runOccupancyTest(TH2F* histo, const DTChamberId& chId, float& chamberPercentage) {
0271   int nBinsX = histo->GetNbinsX();
0272
0273   LogTrace("DTDQM|DTMonitorClient|DTOccupancyTest") << "--- Occupancy test for chamber: " << chId << endl;
0274
0275   int compDeadCell = 0;
0276   int totCell = 0;
0277   int totOccup = 0;
0278
0279   for (int slay = 1; slay <= 3; ++slay) {  // loop over SLs
0280     int binYlow = ((slay - 1) * 4) + 1;
0281
0282     if (chId.station() == 4 && slay == 2)
0283       continue;
0284     for (int lay = 1; lay <= 4; ++lay) {  // loop over layers
0285       DTLayerId layID(chId, slay, lay);
0286       int firstWire = muonGeom->layer(layID)->specificTopology().firstChannel();
0287       int nWires = muonGeom->layer(layID)->specificTopology().channels();
0288       int binY = binYlow + (lay - 1);
0289       int totalDeadCells = 0;
0290       int nDeadCellsInARow = 1;
0291       int nDeadCellsInARowMax = 0;
0292       bool previousIsDead = false;
0293
0294       totCell += nWires;
0295
0296       for (int cell = firstWire; cell != (nWires + firstWire); ++cell) {  // loop over cells
0297         double cellOccup = histo->GetBinContent(cell, binY);
0298         LogTrace("DTDQM|DTMonitorClient|DTOccupancyTest") << "       cell occup: " << cellOccup;
0299         totOccup += cellOccup;
0300
0301         if (cellOccup == 0) {
0302           totalDeadCells++;
0303           if (previousIsDead) {
0304             nDeadCellsInARow++;
0305           } else {
0306             if (nDeadCellsInARow > nDeadCellsInARowMax)
0307               nDeadCellsInARowMax = nDeadCellsInARow;
0308             nDeadCellsInARow = 1;
0309           }
0310           previousIsDead = true;
0311           LogTrace("DTDQM|DTMonitorClient|DTOccupancyTest") << "       below reference" << endl;
0312         } else {
0313           previousIsDead = false;
0314         }
0315         //   // 3 cells not dead between a group of dead cells don't break the count
0316         if (nDeadCellsInARow > nDeadCellsInARowMax)
0317           nDeadCellsInARowMax = nDeadCellsInARow;
0318       }
0319       compDeadCell += totalDeadCells;
0320       if (nDeadCellsInARowMax >= 7.) {
0321         histo->SetBinContent(nBinsX + 1, binY, -1.);
0322       }
0323     }
0324   }
0325
0326   nChannelTotal += totCell;
0327   nChannelDead += compDeadCell;
0328   chamberPercentage = 1. - (float(compDeadCell) / totCell);
0329
0330   int min_occup = nZeroEvtsPC * 20;
0331   if (chId.station() == 3)
0332     min_occup = nZeroEvtsPC * 3;
0333   if (chId.station() == 2)
0334     min_occup = nZeroEvtsPC * 8;
0335   if ((chId.station() == 4) && (chId.sector() == 9))
0336     min_occup = nZeroEvtsPC * 3;
0337   if ((chId.station() == 4) && (chId.sector() == 10))
0338     min_occup = nZeroEvtsPC * 3;
0339   if ((chId.station() == 4) && (chId.sector() == 11))
0340     min_occup = nZeroEvtsPC * 3;
0341   if ((chId.station() == 4) && (chId.sector() == 14))
0342     min_occup = nZeroEvtsPC * 3;
0343
0344   if (totOccup < min_occup)
0345     return 4;
0346   if (totOccup < nMinEvtsPC)
0347     chamberPercentage = 1.;
0348
0349   if (chamberPercentage < 0.2)
0350     return 4;
0351   if (chamberPercentage < 0.5)
0352     return 3;
0353   if (chamberPercentage < 0.75)
0354     return 2;
0355   if (chamberPercentage < 0.9)
0356     return 1;
0357
0358   return 0;
0359 }
0360
0361 string DTOccupancyTest::topFolder() const {
0362   if (tpMode)
0363     return string("DT/10-TestPulses/");
0364   return string("DT/01-Digi/");
0365 }