Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2022-04-29 23:10:17

0001 /*
0002  * \file DTDigiTask.cc
0003  *
0004  * \author M. Zanetti - INFN Padova
0005  *
0006  */
0007 
0008 #include "DQM/DTMonitorModule/interface/DTDigiTask.h"
0009 
0010 // Framework
0011 #include "FWCore/Framework/interface/EventSetup.h"
0012 
0013 // Digis
0014 #include "DataFormats/MuonDetId/interface/DTLayerId.h"
0015 #include "DataFormats/MuonDetId/interface/DTChamberId.h"
0016 
0017 // Geometry
0018 #include "Geometry/DTGeometry/interface/DTGeometry.h"
0019 #include "Geometry/DTGeometry/interface/DTLayer.h"
0020 #include "Geometry/DTGeometry/interface/DTTopology.h"
0021 
0022 // T0s
0023 #include "CondFormats/DTObjects/interface/DTT0.h"
0024 #include "CondFormats/DTObjects/interface/DTTtrig.h"
0025 
0026 #include <sstream>
0027 #include <cmath>
0028 
0029 using namespace edm;
0030 using namespace std;
0031 
0032 // Contructor
0033 DTDigiTask::DTDigiTask(const edm::ParameterSet& ps)
0034     : muonGeomToken_(esConsumes<edm::Transition::BeginRun>()),
0035       readOutMapToken_(esConsumes<edm::Transition::BeginRun>()),
0036       TtrigToken_(esConsumes<edm::Transition::BeginRun>()),
0037       T0Token_(esConsumes<edm::Transition::BeginRun>()),
0038       statusMapToken_(esConsumes()) {
0039   // switch for the verbosity
0040   LogTrace("DTDQM|DTMonitorModule|DTDigiTask") << "[DTDigiTask]: Constructor" << endl;
0041 
0042   // The label to retrieve the digis
0043   dtDigiToken_ = consumes<DTDigiCollection>(ps.getParameter<InputTag>("dtDigiLabel"));
0044   // Read the configuration parameters
0045   maxTDCHits = ps.getUntrackedParameter<int>("maxTDCHitsPerChamber", 30000);
0046   // Set to true to read the ttrig from DB (useful to determine in-time and out-of-time hits)
0047   readTTrigDB = ps.getUntrackedParameter<bool>("readDB", false);
0048   // Set to true to subtract t0 from test pulses
0049   subtractT0 = ps.getParameter<bool>("performPerWireT0Calibration");
0050   // Tmax value (TDC counts)
0051   defaultTmax = ps.getParameter<int>("defaultTmax");
0052   // Switch from static to dinamic histo booking
0053   doStaticBooking = ps.getUntrackedParameter<bool>("staticBooking", true);
0054 
0055   // Switch for local/global runs
0056   isLocalRun = ps.getUntrackedParameter<bool>("localrun", true);
0057   if (!isLocalRun) {
0058     ltcDigiCollectionToken_ = consumes<LTCDigiCollection>(ps.getParameter<edm::InputTag>("ltcDigiCollectionTag"));
0059   }
0060 
0061   // Setting for the reset of the ME after n (= ResetCycle) luminosity sections
0062   resetCycle = ps.getUntrackedParameter<int>("ResetCycle", 3);
0063   // Check the DB of noisy channels
0064   checkNoisyChannels = ps.getUntrackedParameter<bool>("checkNoisyChannels", false);
0065   // Default TTrig to be used when not reading the TTrig DB
0066   defaultTTrig = ps.getParameter<int>("defaultTtrig");
0067   inTimeHitsLowerBound = ps.getParameter<int>("inTimeHitsLowerBound");
0068   inTimeHitsUpperBound = ps.getParameter<int>("inTimeHitsUpperBound");
0069   timeBoxGranularity = ps.getUntrackedParameter<int>("timeBoxGranularity", 4);
0070   maxTTMounts = ps.getUntrackedParameter<int>("maxTTMounts", 6400);
0071 
0072   doAllHitsOccupancies = ps.getUntrackedParameter<bool>("doAllHitsOccupancies", true);
0073   doNoiseOccupancies = ps.getUntrackedParameter<bool>("doNoiseOccupancies", false);
0074   doInTimeOccupancies = ps.getUntrackedParameter<bool>("doInTimeOccupancies", false);
0075 
0076   // switch on the mode for running on test pulses (different top folder)
0077   tpMode = ps.getUntrackedParameter<bool>("testPulseMode", false);
0078   // switch on/off the filtering of synchronous noise events (cutting on the # of digis)
0079   // time-boxes and occupancy plots are not filled and summary plots are created to report the problem
0080   filterSyncNoise = ps.getUntrackedParameter<bool>("filterSyncNoise", false);
0081   // look for synch noisy events, produce histograms but do not filter them
0082   lookForSyncNoise = ps.getUntrackedParameter<bool>("lookForSyncNoise", false);
0083 
0084   // switch on the mode for running on slice test (different top folder and other customizations)
0085   sliceTestMode = ps.getUntrackedParameter<bool>("sliceTestMode", false);
0086   // time pedestal used to set the minimum in the time box plot
0087   tdcPedestal = ps.getUntrackedParameter<int>("tdcPedestal", 0);
0088 
0089   // switch on production of time-boxes with layer granularity
0090   doLayerTimeBoxes = ps.getUntrackedParameter<bool>("doLayerTimeBoxes", false);
0091 
0092   syncNumTot = 0;
0093   syncNum = 0;
0094 }
0095 
0096 // destructor
0097 DTDigiTask::~DTDigiTask() {
0098   LogTrace("DTDQM|DTMonitorModule|DTDigiTask") << "DTDigiTask: analyzed " << nevents << " events" << endl;
0099 }
0100 
0101 void DTDigiTask::dqmBeginRun(const edm::Run& run, const edm::EventSetup& context) {
0102   LogTrace("DTDQM|DTMonitorModule|DTDigiTask") << "[DTDigiTask]: begin run" << endl;
0103   nevents = 0;
0104 
0105   // Get the geometry
0106   muonGeom = &context.getData(muonGeomToken_);
0107 
0108   // map of the channels
0109   mapping = &context.getData(readOutMapToken_);
0110 
0111   // tTrig
0112   if (readTTrigDB)
0113     tTrigMap = &context.getData(TtrigToken_);
0114   // t0s
0115   if (subtractT0)
0116     t0Map = &context.getData(T0Token_);
0117   // FIXME: tMax (not yet from the DB)
0118   tMax = defaultTmax;
0119 
0120   // ----------------------------------------------------------------------
0121 }
0122 
0123 void DTDigiTask::bookHistograms(DQMStore::IBooker& ibooker, edm::Run const& run, edm::EventSetup const& context) {
0124   if (doStaticBooking) {  // Static histo booking
0125     // book the event counter
0126     ibooker.setCurrentFolder("DT/EventInfo/Counters");
0127     nEventMonitor = ibooker.bookFloat(tpMode ? "nProcessedEventsDigiTP" : "nProcessedEventsDigi");
0128     ibooker.setCurrentFolder(topFolder());
0129     for (int wh = -2; wh <= 2; ++wh) {  // loop over wheels
0130       if (sliceTestMode && wh != 2)
0131         continue;
0132 
0133       if (doAllHitsOccupancies)
0134         bookHistos(ibooker, wh, string("Occupancies"), "OccupancyAllHits");
0135 
0136       if (doNoiseOccupancies)
0137         bookHistos(ibooker, wh, string("Occupancies"), "OccupancyNoise");
0138       if (doInTimeOccupancies)
0139         bookHistos(ibooker, wh, string("Occupancies"), "OccupancyInTimeHits");
0140 
0141       if (lookForSyncNoise || filterSyncNoise) {
0142         bookHistos(ibooker, wh, string("SynchNoise"), "SyncNoiseEvents");
0143         bookHistos(ibooker, wh, string("SynchNoise"), "SyncNoiseChambs");
0144       }
0145 
0146       for (int st = 1; st <= 4; ++st) {           // loop over stations
0147         for (int sect = 1; sect <= 14; ++sect) {  // loop over sectors
0148           if ((sect == 13 || sect == 14) && st != 4)
0149             continue;
0150 
0151           if (sliceTestMode && (sect != 12 || wh != 2))
0152             continue;
0153 
0154           // Get the chamber ID
0155           const DTChamberId dtChId(wh, st, sect);
0156 
0157           // Occupancies
0158           if (doAllHitsOccupancies) {
0159             bookHistos(ibooker, dtChId, string("Occupancies"), "OccupancyAllHits_perCh");
0160             // set channel mapping
0161             channelsMap(dtChId, "OccupancyAllHits_perCh");
0162           }
0163           if (doNoiseOccupancies)
0164             bookHistos(ibooker, dtChId, string("Occupancies"), "OccupancyNoise_perCh");
0165           if (doInTimeOccupancies)
0166             bookHistos(ibooker, dtChId, string("Occupancies"), "OccupancyInTimeHits_perCh");
0167 
0168           for (int sl = 1; sl <= 3; ++sl) {  // Loop over SLs
0169             if (st == 4 && sl == 2)
0170               continue;
0171             const DTSuperLayerId dtSLId(wh, st, sect, sl);
0172             if (isLocalRun) {
0173               bookHistos(ibooker, dtSLId, string("TimeBoxes"), "TimeBox");
0174             } else {
0175               // TimeBoxes for different triggers
0176               bookHistos(ibooker, dtSLId, string("TimeBoxes"), "TimeBoxDTonly");
0177               bookHistos(ibooker, dtSLId, string("TimeBoxes"), "TimeBoxNoDT");
0178               bookHistos(ibooker, dtSLId, string("TimeBoxes"), "TimeBoxDTalso");
0179             }
0180           }
0181         }
0182       }
0183     }
0184   }
0185 }
0186 
0187 void DTDigiTask::beginLuminosityBlock(LuminosityBlock const& lumiSeg, EventSetup const& context) {
0188   LogTrace("DTDQM|DTMonitorModule|DTDigiTask") << "[DTDigiTask]: Begin of LS transition" << endl;
0189 
0190   // Reset the MonitorElements every n (= ResetCycle) Lumi Blocks
0191   int lumiBlock = lumiSeg.id().luminosityBlock();
0192   if (lumiBlock % resetCycle == 0) {
0193     LogVerbatim("DTDQM|DTMonitorModule|DTDigiTask")
0194         << "[DTDigiTask]: Reset at the LS transition : " << lumiBlock << endl;
0195     // Loop over all ME
0196     map<string, map<uint32_t, MonitorElement*> >::const_iterator histosIt = digiHistos.begin();
0197     map<string, map<uint32_t, MonitorElement*> >::const_iterator histosEnd = digiHistos.end();
0198     for (; histosIt != histosEnd; ++histosIt) {
0199       map<uint32_t, MonitorElement*>::const_iterator histoIt = (*histosIt).second.begin();
0200       map<uint32_t, MonitorElement*>::const_iterator histoEnd = (*histosIt).second.end();
0201       for (; histoIt != histoEnd; ++histoIt) {
0202         (*histoIt).second->Reset();
0203       }
0204     }
0205 
0206     // re-set mapping for not real channels in the occupancyHits per chamber
0207     for (int wh = -2; wh <= 2; wh++) {
0208       for (int sect = 1; sect <= 14; sect++) {
0209         for (int st = 1; st <= 4; st++) {
0210           if (sliceTestMode && (sect != 12 || wh != 2)) {
0211             continue;
0212           }
0213           if ((sect == 13 || sect == 14) && st != 4) {
0214             continue;
0215           }
0216           const DTChamberId dtChId(wh, st, sect);
0217           channelsMap(dtChId, "OccupancyAllHits_perCh");
0218         }
0219       }
0220     }
0221 
0222     // loop over wheel summaries
0223     map<string, map<int, MonitorElement*> >::const_iterator whHistosIt = wheelHistos.begin();
0224     map<string, map<int, MonitorElement*> >::const_iterator whHistosEnd = wheelHistos.end();
0225     for (; whHistosIt != whHistosEnd; ++whHistosIt) {
0226       if ((*whHistosIt).first.find("Sync") == string::npos) {  // FIXME skips synch noise plots
0227         map<int, MonitorElement*>::const_iterator histoIt = (*whHistosIt).second.begin();
0228         map<int, MonitorElement*>::const_iterator histoEnd = (*whHistosIt).second.end();
0229         for (; histoIt != histoEnd; ++histoIt) {
0230           (*histoIt).second->Reset();
0231         }
0232       }
0233     }
0234   }
0235 }
0236 
0237 void DTDigiTask::bookHistos(DQMStore::IBooker& ibooker, const DTSuperLayerId& dtSL, string folder, string histoTag) {
0238   // set the folder
0239   stringstream wheel;
0240   wheel << dtSL.wheel();
0241   stringstream station;
0242   station << dtSL.station();
0243   stringstream sector;
0244   sector << dtSL.sector();
0245   stringstream superLayer;
0246   superLayer << dtSL.superlayer();
0247   ibooker.setCurrentFolder(topFolder() + "Wheel" + wheel.str() + "/Sector" + sector.str() + "/Station" + station.str());
0248 
0249   // Build the histo name
0250   string histoName =
0251       histoTag + "_W" + wheel.str() + "_St" + station.str() + "_Sec" + sector.str() + "_SL" + superLayer.str();
0252 
0253   LogTrace("DTDQM|DTMonitorModule|DTDigiTask")
0254       << "[DTDigiTask]: booking SL histo:" << histoName << " (tag: " << histoTag << ") folder: "
0255       << topFolder() + "Wheel" + wheel.str() + "/Station" + station.str() + "/Sector" + sector.str() + "/" + folder
0256       << endl;
0257 
0258   // ttrig and rms are TDC counts
0259   if (readTTrigDB)
0260     tTrigMap->get(dtSL, tTrig, tTrigRMS, kFactor, DTTimeUnits::counts);
0261   else
0262     tTrig = defaultTTrig;
0263 
0264   if (folder == "TimeBoxes") {
0265     string histoTitle = histoName + " (TDC Counts)";
0266 
0267     if (!readTTrigDB) {
0268       (digiHistos[histoTag])[dtSL.rawId()] = ibooker.book1D(
0269           histoName, histoTitle, maxTTMounts / timeBoxGranularity, tdcPedestal, maxTTMounts + tdcPedestal);
0270       if (doLayerTimeBoxes) {  // Book TimeBoxes per layer
0271         for (int layer = 1; layer != 5; ++layer) {
0272           DTLayerId layerId(dtSL, layer);
0273           stringstream layerHistoName;
0274           layerHistoName << histoName << "_L" << layer;
0275           (digiHistos[histoTag])[layerId.rawId()] = ibooker.book1D(layerHistoName.str(),
0276                                                                    layerHistoName.str(),
0277                                                                    maxTTMounts / timeBoxGranularity,
0278                                                                    tdcPedestal,
0279                                                                    maxTTMounts + tdcPedestal);
0280         }
0281       }
0282     } else {
0283       (digiHistos[histoTag])[dtSL.rawId()] =
0284           ibooker.book1D(histoName, histoTitle, 3 * tMax / timeBoxGranularity, tTrig - tMax, tTrig + 2 * tMax);
0285       if (doLayerTimeBoxes) {
0286         // Book TimeBoxes per layer
0287         for (int layer = 1; layer != 5; ++layer) {
0288           DTLayerId layerId(dtSL, layer);
0289           stringstream layerHistoName;
0290           layerHistoName << histoName << "_L" << layer;
0291           (digiHistos[histoTag])[layerId.rawId()] = ibooker.book1D(
0292               layerHistoName.str(), layerHistoName.str(), 3 * tMax / timeBoxGranularity, tTrig - tMax, tTrig + 2 * tMax);
0293         }
0294       }
0295     }
0296   }
0297 
0298   if (folder == "CathodPhotoPeaks") {
0299     ibooker.setCurrentFolder(topFolder() + "Wheel" + wheel.str() + "/Sector" + sector.str() + "/Station" +
0300                              station.str() + "/" + folder);
0301     (digiHistos[histoTag])[dtSL.rawId()] = ibooker.book1D(histoName, histoName, 500, 0, 1000);
0302   }
0303 }
0304 
0305 void DTDigiTask::bookHistos(DQMStore::IBooker& ibooker, const DTChamberId& dtCh, string folder, string histoTag) {
0306   // set the current folder
0307   stringstream wheel;
0308   wheel << dtCh.wheel();
0309   stringstream station;
0310   station << dtCh.station();
0311   stringstream sector;
0312   sector << dtCh.sector();
0313   ibooker.setCurrentFolder(topFolder() + "Wheel" + wheel.str() + "/Sector" + sector.str() + "/Station" + station.str());
0314 
0315   // build the histo name
0316   string histoName = histoTag + "_W" + wheel.str() + "_St" + station.str() + "_Sec" + sector.str();
0317 
0318   LogTrace("DTDQM|DTMonitorModule|DTDigiTask")
0319       << "[DTDigiTask]: booking chamber histo:"
0320       << " (tag: " << histoTag
0321       << ") folder: " << topFolder() + "Wheel" + wheel.str() + "/Station" + station.str() + "/Sector" + sector.str()
0322       << endl;
0323 
0324   if (folder == "Occupancies") {
0325     const DTChamber* dtchamber = muonGeom->chamber(dtCh);
0326     const std::vector<const DTSuperLayer*>& dtSupLylist = dtchamber->superLayers();
0327     std::vector<const DTSuperLayer*>::const_iterator suly = dtSupLylist.begin();
0328     std::vector<const DTSuperLayer*>::const_iterator sulyend = dtSupLylist.end();
0329 
0330     int nWires = 0;
0331     int firstWire = 0;
0332     int nWires_max = 0;
0333 
0334     while (suly != sulyend) {
0335       const std::vector<const DTLayer*> dtLyList = (*suly)->layers();
0336       std::vector<const DTLayer*>::const_iterator ly = dtLyList.begin();
0337       std::vector<const DTLayer*>::const_iterator lyend = dtLyList.end();
0338       stringstream superLayer;
0339       superLayer << (*suly)->id().superlayer();
0340 
0341       while (ly != lyend) {
0342         nWires = muonGeom->layer((*ly)->id())->specificTopology().channels();
0343         firstWire = muonGeom->layer((*ly)->id())->specificTopology().firstChannel();
0344         stringstream layer;
0345         layer << (*ly)->id().layer();
0346         string histoName_layer = histoName + "_SL" + superLayer.str() + "_L" + layer.str();
0347         if (histoTag == "OccupancyAllHits_perL" || histoTag == "OccupancyNoise_perL" ||
0348             histoTag == "OccupancyInTimeHits_perL")
0349           (digiHistos[histoTag])[(*ly)->id().rawId()] =
0350               ibooker.book1D(histoName_layer, histoName_layer, nWires, firstWire, nWires + firstWire);
0351         ++ly;
0352         if ((nWires + firstWire) > nWires_max)
0353           nWires_max = (nWires + firstWire);
0354       }
0355       ++suly;
0356     }
0357 
0358     if (histoTag != "OccupancyAllHits_perL" && histoTag != "OccupancyNoise_perL" &&
0359         histoTag != "OccupancyInTimeHits_perL") {
0360       // Set the title to show the time interval used (only if unique == not from DB)
0361       string histoTitle = histoName;
0362       if (!readTTrigDB && histoTag == "OccupancyInTimeHits_perCh") {
0363         stringstream title;
0364         int inTimeHitsLowerBoundCorr = int(round(defaultTTrig)) - inTimeHitsLowerBound;
0365         int inTimeHitsUpperBoundCorr = int(round(defaultTTrig)) + defaultTmax + inTimeHitsUpperBound;
0366         title << "Occ. digis in time [" << inTimeHitsLowerBoundCorr << ", " << inTimeHitsUpperBoundCorr
0367               << "] (TDC counts)";
0368         histoTitle = title.str();
0369       }
0370       (digiHistos[histoTag])[dtCh.rawId()] =
0371           ibooker.book2D(histoName, histoTitle, nWires_max, 1, nWires_max + 1, 12, 0, 12);
0372 
0373       for (int i = 1; i <= 12; i++) {
0374         if (i < 5) {
0375           stringstream layer;
0376           string layer_name;
0377           layer << i;
0378           layer >> layer_name;
0379           string label = "SL1: L" + layer_name;
0380           (digiHistos[histoTag])[dtCh.rawId()]->setBinLabel(i, label, 2);
0381         } else if (i > 4 && i < 9) {
0382           stringstream layer;
0383           string layer_name;
0384           layer << (i - 4);
0385           layer >> layer_name;
0386           string label = "SL2: L" + layer_name;
0387           (digiHistos[histoTag])[dtCh.rawId()]->setBinLabel(i, label, 2);
0388         } else if (i > 8 && i < 13) {
0389           stringstream layer;
0390           string layer_name;
0391           layer << (i - 8);
0392           layer >> layer_name;
0393           string label = "SL3: L" + layer_name;
0394           (digiHistos[histoTag])[dtCh.rawId()]->setBinLabel(i, label, 2);
0395         }
0396       }
0397     }
0398   }
0399 }
0400 
0401 void DTDigiTask::bookHistos(DQMStore::IBooker& ibooker, const int wheelId, string folder, string histoTag) {
0402   // Set the current folder
0403   stringstream wheel;
0404   wheel << wheelId;
0405 
0406   // build the histo name
0407   string histoName = histoTag + "_W" + wheel.str();
0408 
0409   LogTrace("DTDQM|DTMonitorModule|DTDigiTask")
0410       << "[DTDigiTask]: booking wheel histo:" << histoName << " (tag: " << histoTag
0411       << ") folder: " << topFolder() + "Wheel" + wheel.str() + "/" << endl;
0412 
0413   if (folder == "Occupancies") {
0414     ibooker.setCurrentFolder(topFolder() + "Wheel" + wheel.str());
0415     string histoTitle = "# of digis per chamber WHEEL: " + wheel.str();
0416     (wheelHistos[histoTag])[wheelId] = ibooker.book2D(histoName, histoTitle, 12, 1, 13, 4, 1, 5);
0417     (wheelHistos[histoTag])[wheelId]->setBinLabel(1, "MB1", 2);
0418     (wheelHistos[histoTag])[wheelId]->setBinLabel(2, "MB2", 2);
0419     (wheelHistos[histoTag])[wheelId]->setBinLabel(3, "MB3", 2);
0420     (wheelHistos[histoTag])[wheelId]->setBinLabel(4, "MB4", 2);
0421     (wheelHistos[histoTag])[wheelId]->setAxisTitle("sector", 1);
0422   } else if (folder == "SynchNoise") {
0423     ibooker.setCurrentFolder("DT/05-Noise/SynchNoise");
0424     if (histoTag == "SyncNoiseEvents") {
0425       string histoTitle = "# of Syncronous-noise events WHEEL: " + wheel.str();
0426       (wheelHistos[histoTag])[wheelId] = ibooker.book2D(histoName, histoTitle, 12, 1, 13, 4, 1, 5);
0427       (wheelHistos[histoTag])[wheelId]->setBinLabel(1, "MB1", 2);
0428       (wheelHistos[histoTag])[wheelId]->setBinLabel(2, "MB2", 2);
0429       (wheelHistos[histoTag])[wheelId]->setBinLabel(3, "MB3", 2);
0430       (wheelHistos[histoTag])[wheelId]->setBinLabel(4, "MB4", 2);
0431       (wheelHistos[histoTag])[wheelId]->setAxisTitle("sector", 1);
0432     } else if (histoTag == "SyncNoiseChambs") {
0433       string histoTitle = "# of Synchornous-noise chamb per evt. WHEEL: " + wheel.str();
0434       (wheelHistos[histoTag])[wheelId] = ibooker.book1D(histoName, histoTitle, 50, 0.5, 50.5);
0435       (wheelHistos[histoTag])[wheelId]->setAxisTitle("# of noisy chambs.", 1);
0436       (wheelHistos[histoTag])[wheelId]->setAxisTitle("# of evts.", 2);
0437     }
0438   }
0439 }
0440 // does the real job
0441 void DTDigiTask::analyze(const edm::Event& event, const edm::EventSetup& c) {
0442   nevents++;
0443   nEventMonitor->Fill(nevents);
0444   if (nevents % 1000 == 0) {
0445     LogTrace("DTDQM|DTMonitorModule|DTDigiTask")
0446         << "[DTDigiTask] Analyze #Run: " << event.id().run() << " #Event: " << event.id().event() << endl;
0447   }
0448 
0449   // Get the ingredients from the event
0450 
0451   // Digi collection
0452   edm::Handle<DTDigiCollection> dtdigis;
0453   event.getByToken(dtDigiToken_, dtdigis);
0454 
0455   // LTC digis
0456   if (!isLocalRun)
0457     event.getByToken(ltcDigiCollectionToken_, ltcdigis);
0458 
0459   // Status map (for noisy channels)
0460   if (checkNoisyChannels) {
0461     // Get the map of noisy channels
0462     statusMap = &c.getData(statusMapToken_);
0463   }
0464 
0465   string histoTag;
0466 
0467   // Check if the digi container is empty
0468   if (dtdigis->begin() == dtdigis->end()) {
0469     LogTrace("DTDQM|DTMonitorModule|DTDigiTask") << "Event " << nevents << " empty." << endl;
0470   }
0471 
0472   if (lookForSyncNoise || filterSyncNoise) {  // dosync
0473     // Count the # of digis per chamber
0474     DTDigiCollection::DigiRangeIterator dtLayerId_It;
0475     for (dtLayerId_It = dtdigis->begin(); dtLayerId_It != dtdigis->end(); dtLayerId_It++) {
0476       DTChamberId chId = ((*dtLayerId_It).first).chamberId();
0477       if (hitMap.find(chId) == hitMap.end()) {  // new chamber
0478         hitMap[chId] = 0;
0479       }
0480       hitMap[chId] += (((*dtLayerId_It).second).second - ((*dtLayerId_It).second).first);
0481     }
0482 
0483     // check chamber with # of digis above threshold and flag them as noisy
0484     map<DTChamberId, int>::const_iterator hitMapIt = hitMap.begin();
0485     map<DTChamberId, int>::const_iterator hitMapEnd = hitMap.end();
0486 
0487     map<int, int> chMap;
0488 
0489     for (; hitMapIt != hitMapEnd; ++hitMapIt) {
0490       if ((hitMapIt->second) > maxTDCHits) {
0491         DTChamberId chId = hitMapIt->first;
0492         int wh = chId.wheel();
0493 
0494         LogTrace("DTDQM|DTMonitorModule|DTDigiTask")
0495             << "[DTDigiTask] Synch noise in chamber: " << chId << " with # digis: " << hitMapIt->second << endl;
0496 
0497         if (chMap.find(wh) == chMap.end()) {
0498           chMap[wh] = 0;
0499         }
0500         chMap[wh]++;
0501 
0502         syncNoisyChambers.insert(chId);
0503 
0504         wheelHistos["SyncNoiseEvents"][wh]->Fill(chId.sector(), chId.station());
0505       }
0506     }
0507 
0508     // fill # of noisy ch per wheel plot
0509     map<int, int>::const_iterator chMapIt = chMap.begin();
0510     map<int, int>::const_iterator chMapEnd = chMap.end();
0511     for (; chMapIt != chMapEnd; ++chMapIt) {
0512       wheelHistos["SyncNoiseChambs"][(*chMapIt).first]->Fill((*chMapIt).second);
0513     }
0514 
0515     // clear the map of # of digis per chamber: not needed anymore
0516     hitMap.clear();
0517 
0518     if (!syncNoisyChambers.empty()) {
0519       LogVerbatim("DTDQM|DTMonitorModule|DTDigiTask") << "[DTDigiTask] Synch Noise in event: " << nevents;
0520       if (filterSyncNoise)
0521         LogVerbatim("DTDQM|DTMonitorModule|DTDigiTask")
0522             << "\tnoisy time-boxes and occupancy will not be filled!" << endl;
0523       syncNumTot++;
0524       syncNum++;
0525     }
0526 
0527     // Logging of "large" synch Noisy events in private DQM
0528     if (syncNoisyChambers.size() > 3) {
0529       time_t eventTime = time_t(event.time().value() >> 32);
0530 
0531       LogVerbatim("DTDQM|DTMonitorModule|DTDigiTask|DTSynchNoise")
0532           << "[DTDigiTask] At least 4 Synch Noisy chambers in Run : " << event.id().run()
0533           << " Lumi : " << event.id().luminosityBlock() << " Event : " << event.id().event()
0534           << " at time : " << ctime(&eventTime) << endl;
0535 
0536       set<DTChamberId>::const_iterator chIt = syncNoisyChambers.begin();
0537       set<DTChamberId>::const_iterator chEnd = syncNoisyChambers.end();
0538 
0539       stringstream synchNoisyCh;
0540       for (; chIt != chEnd; ++chIt) {
0541         synchNoisyCh << " " << (*chIt);
0542       }
0543       LogVerbatim("DTDQM|DTMonitorModule|DTDigiTask|DTSynchNoise")
0544           << "[DTDigiTask] Chamber List :" << synchNoisyCh.str() << endl;
0545     }
0546 
0547     if (nevents % 1000 == 0) {
0548       LogVerbatim("DTDQM|DTMonitorModule|DTDigiTask")
0549           << (syncNumTot * 100. / nevents) << "% sync noise events since the beginning \n"
0550           << (syncNum * 0.1) << "% sync noise events in the last 1000 events " << endl;
0551       syncNum = 0;
0552     }
0553   }
0554 
0555   bool isSyncNoisy = false;
0556 
0557   DTDigiCollection::DigiRangeIterator dtLayerId_It;
0558   for (dtLayerId_It = dtdigis->begin(); dtLayerId_It != dtdigis->end(); ++dtLayerId_It) {  // Loop over layers
0559     isSyncNoisy = false;
0560     // check if chamber labeled as synch noisy
0561     if (filterSyncNoise) {
0562       DTChamberId chId = ((*dtLayerId_It).first).chamberId();
0563       if (syncNoisyChambers.find(chId) != syncNoisyChambers.end()) {
0564         isSyncNoisy = true;
0565       }
0566     }
0567 
0568     for (DTDigiCollection::const_iterator digiIt = ((*dtLayerId_It).second).first;
0569          digiIt != ((*dtLayerId_It).second).second;
0570          ++digiIt) {  // Loop over all digis
0571 
0572       bool isNoisy = false;
0573       bool isFEMasked = false;
0574       bool isTDCMasked = false;
0575       bool isTrigMask = false;
0576       bool isDead = false;
0577       bool isNohv = false;
0578       if (checkNoisyChannels) {
0579         const DTWireId wireId(((*dtLayerId_It).first), (*digiIt).wire());
0580         statusMap->cellStatus(wireId, isNoisy, isFEMasked, isTDCMasked, isTrigMask, isDead, isNohv);
0581       }
0582 
0583       // Get the useful IDs
0584       const DTSuperLayerId dtSLId = ((*dtLayerId_It).first).superlayerId();
0585       uint32_t indexSL = dtSLId.rawId();
0586       const DTChamberId dtChId = dtSLId.chamberId();
0587       uint32_t indexCh = dtChId.rawId();
0588       int layer_number = ((*dtLayerId_It).first).layer();
0589       int superlayer_number = dtSLId.superlayer();
0590 
0591       // Read the ttrig DB or set a rough value from config
0592       // ttrig and rms are TDC counts
0593       if (readTTrigDB)
0594         tTrigMap->get(((*dtLayerId_It).first).superlayerId(), tTrig, tTrigRMS, kFactor, DTTimeUnits::counts);
0595       else
0596         tTrig = defaultTTrig;
0597 
0598       int inTimeHitsLowerBoundCorr = int(round(tTrig)) - inTimeHitsLowerBound;
0599       int inTimeHitsUpperBoundCorr = int(round(tTrig)) + tMax + inTimeHitsUpperBound;
0600 
0601       float t0;
0602       float t0RMS;
0603       int tdcTime = (*digiIt).countsTDC();
0604 
0605       if (subtractT0) {
0606         const DTWireId dtWireId(((*dtLayerId_It).first), (*digiIt).wire());
0607         // t0s and rms are TDC counts
0608         t0Map->get(dtWireId, t0, t0RMS, DTTimeUnits::counts);
0609         tdcTime += int(round(t0));
0610       }
0611 
0612       if (sliceTestMode) {
0613         tdcTime = std::max(tdcPedestal + 1, std::min(tdcPedestal + maxTTMounts - 1, tdcTime));
0614         // std::cout << tdcTime << std::endl;
0615       }
0616 
0617       // Fill Time-Boxes
0618       // NOTE: avoid to fill TB and PhotoPeak with noise. Occupancy are filled anyway
0619       if ((!isNoisy) && (!isSyncNoisy)) {  // Discard noisy channels
0620         // TimeBoxes per SL
0621         histoTag = "TimeBox" + triggerSource();
0622 
0623         (digiHistos.find(histoTag)->second).find(indexSL)->second->Fill(tdcTime);
0624         if (doLayerTimeBoxes)
0625           (digiHistos.find(histoTag)->second).find((*dtLayerId_It).first.rawId())->second->Fill(tdcTime);
0626       }
0627 
0628       // Fill Occupancies
0629       if (!isSyncNoisy) {  // Discard synch noisy channels
0630 
0631         if (doAllHitsOccupancies) {  // fill occupancies for all hits
0632           //Occupancies per chamber & layer
0633           histoTag = "OccupancyAllHits_perCh";
0634           map<uint32_t, MonitorElement*>::const_iterator mappedHisto = digiHistos[histoTag].find(indexCh);
0635 
0636           //FR comment the following cannot pass ibooker to analyze method!
0637           /*
0638           if (mappedHisto == digiHistos[histoTag].end()) { // dynamic booking
0639             bookHistos(ibooker, dtChId, string("Occupancies"), histoTag);
0640             mappedHisto = digiHistos[histoTag].find(indexCh);
0641           }
0642       */
0643           mappedHisto->second->Fill((*digiIt).wire(), (layer_number + (superlayer_number - 1) * 4) - 1);
0644 
0645           // Fill the chamber occupancy
0646           histoTag = "OccupancyAllHits";
0647           map<int, MonitorElement*>::const_iterator histoPerWheel = wheelHistos[histoTag].find(dtChId.wheel());
0648 
0649           histoPerWheel->second->Fill(dtChId.sector(), dtChId.station());  // FIXME: normalize to # of layers
0650         }
0651 
0652         if (doNoiseOccupancies) {  // fill occupancies for hits before the ttrig
0653           if (tdcTime < inTimeHitsLowerBoundCorr) {
0654             // FIXME: what about tdcTime > inTimeHitsUpperBoundCorr ???
0655 
0656             // Noise: Before tTrig
0657             //Occupancies Noise per chamber & layer
0658             histoTag = "OccupancyNoise_perCh";
0659             map<uint32_t, MonitorElement*>::const_iterator mappedHisto = digiHistos[histoTag].find(indexCh);
0660 
0661             mappedHisto->second->Fill((*digiIt).wire(), (layer_number + (superlayer_number - 1) * 4) - 1);
0662 
0663             // Fill the chamber occupancy
0664 
0665             histoTag = "OccupancyNoise";
0666             map<int, MonitorElement*>::const_iterator histoPerWheel = wheelHistos[histoTag].find(dtChId.wheel());
0667 
0668             histoPerWheel->second->Fill(dtChId.sector(), dtChId.station());  // FIXME: normalize to # of layers
0669           }
0670         }
0671 
0672         if (doInTimeOccupancies) {  // fill occpunacies for in-time hits only
0673           if (tdcTime > inTimeHitsLowerBoundCorr && tdcTime < inTimeHitsUpperBoundCorr) {
0674             // Physical hits: within the time window
0675 
0676             //Occupancies Signal per chamber & layer
0677             histoTag = "OccupancyInTimeHits_perCh";
0678             map<uint32_t, MonitorElement*>::const_iterator mappedHisto = digiHistos[histoTag].find(indexCh);
0679 
0680             mappedHisto->second->Fill((*digiIt).wire(), (layer_number + (superlayer_number - 1) * 4) - 1);
0681 
0682             // Fill the chamber occupancy
0683             histoTag = "OccupancyInTimeHits";
0684             map<int, MonitorElement*>::const_iterator histoPerWheel = wheelHistos[histoTag].find(dtChId.wheel());
0685 
0686             histoPerWheel->second->Fill(dtChId.sector(), dtChId.station());  // FIXME: normalize to # of layers
0687           }
0688         }
0689       }
0690     }
0691   }
0692 
0693   syncNoisyChambers.clear();
0694 }
0695 
0696 string DTDigiTask::triggerSource() {
0697   string l1ASource;
0698   if (isLocalRun)
0699     return l1ASource;
0700 
0701   for (std::vector<LTCDigi>::const_iterator ltc_it = ltcdigis->begin(); ltc_it != ltcdigis->end(); ltc_it++) {
0702     size_t otherTriggerSum = 0;
0703     for (size_t i = 1; i < 6; i++)
0704       otherTriggerSum += size_t((*ltc_it).HasTriggered(i));
0705 
0706     if ((*ltc_it).HasTriggered(0) && otherTriggerSum == 0)
0707       l1ASource = "DTonly";
0708     else if (!(*ltc_it).HasTriggered(0))
0709       l1ASource = "NoDT";
0710     else if ((*ltc_it).HasTriggered(0) && otherTriggerSum > 0)
0711       l1ASource = "DTalso";
0712   }
0713 
0714   return l1ASource;
0715 }
0716 
0717 string DTDigiTask::topFolder() const {
0718   if (tpMode)
0719     return string("DT/10-TestPulses/");
0720   else if (sliceTestMode)
0721     return string("DT/01-SliceTestDigi/");
0722   return string("DT/01-Digi/");
0723 }
0724 
0725 void DTDigiTask::channelsMap(const DTChamberId& dtCh, string histoTag) {
0726   // n max channels
0727   int nWires_max = (digiHistos[histoTag])[dtCh.rawId()]->getNbinsX();
0728 
0729   // set bin content = -1 for each not real channel. For visualization purposes
0730   for (int sl = 1; sl <= 3; sl++) {
0731     for (int ly = 1; ly <= 4; ly++) {
0732       for (int ch = 1; ch <= nWires_max; ch++) {
0733         int dduId = -1, rosId = -1, robId = -1, tdcId = -1, channelId = -1;
0734         int realCh = mapping->geometryToReadOut(
0735             dtCh.wheel(), dtCh.station(), dtCh.sector(), sl, ly, ch, dduId, rosId, robId, tdcId, channelId);
0736 
0737         // realCh = 0 if the channel exists, while realCh = 1 if it does not exist
0738         if (realCh) {
0739           int lybin = (4 * sl - 4) + ly;
0740           (digiHistos[histoTag])[dtCh.rawId()]->setBinContent(ch, lybin, -1.);
0741         }
0742       }
0743     }
0744   }
0745 }
0746 
0747 // Local Variables:
0748 // show-trailing-whitespace: t
0749 // truncate-lines: t
0750 // End: