Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-05-21 22:39:44

0001 
0002 /*
0003  *  See header file for a description of this class.
0004  *
0005  *  \authors G. Mila , G. Cerminara - INFN Torino
0006  */
0007 
0008 #include "DTNoiseTask.h"
0009 
0010 // Framework
0011 #include "FWCore/Framework/interface/Event.h"
0012 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0013 #include "FWCore/Framework/interface/EventSetup.h"
0014 #include "FWCore/ServiceRegistry/interface/Service.h"
0015 
0016 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0017 
0018 // Geometry
0019 #include "Geometry/DTGeometry/interface/DTGeometry.h"
0020 #include "Geometry/DTGeometry/interface/DTLayer.h"
0021 #include "Geometry/DTGeometry/interface/DTTopology.h"
0022 
0023 // Digi
0024 #include "CondFormats/DataRecord/interface/DTStatusFlagRcd.h"
0025 #include "CondFormats/DTObjects/interface/DTStatusFlag.h"
0026 #include "CondFormats/DTObjects/interface/DTReadOutMapping.h"
0027 
0028 // Database
0029 #include "CondFormats/DTObjects/interface/DTTtrig.h"
0030 
0031 #include <sstream>
0032 #include <string>
0033 
0034 using namespace edm;
0035 using namespace std;
0036 
0037 DTNoiseTask::DTNoiseTask(const ParameterSet& ps)
0038     : evtNumber(0),
0039       muonGeomToken_(esConsumes<edm::Transition::BeginRun>()),
0040       tTrigMapToken_(esConsumes<edm::Transition::BeginRun>()) {
0041   LogVerbatim("DTNoiseTask") << "[DTNoiseTask]: Constructor" << endl;
0042 
0043   //switch for timeBox booking
0044   doTimeBoxHistos = ps.getUntrackedParameter<bool>("doTbHistos", false);
0045 
0046   // The label to retrieve the digis
0047   dtDigiToken_ = consumes<DTDigiCollection>(ps.getParameter<InputTag>("dtDigiLabel"));
0048 
0049   // the name of the 4D rec hits collection
0050   recHits4DToken_ = consumes<DTRecSegment4DCollection>(edm::InputTag(ps.getParameter<string>("recHits4DLabel")));
0051 
0052   // switch for segment veto
0053   doSegmentVeto = ps.getUntrackedParameter<bool>("doSegmentVeto", false);
0054 
0055   // safe margin (ns) between ttrig and beginning of counting area
0056   safeMargin = ps.getUntrackedParameter<double>("safeMargin", 200.);
0057 }
0058 
0059 DTNoiseTask::~DTNoiseTask() {}
0060 
0061 /// Analyze
0062 void DTNoiseTask::analyze(const edm::Event& e, const edm::EventSetup& c) {
0063   evtNumber++;
0064   nEventMonitor->Fill(evtNumber);
0065 
0066   if (evtNumber % 1000 == 0)
0067     LogVerbatim("DTNoiseTask") << "[DTNoiseTask]: Analyzing evt number :" << evtNumber << endl;
0068 
0069   // map of the chambers with at least 1 segment
0070   std::map<DTChamberId, int> segmentsChId;
0071 
0072   // Get the 4D segment collection from the event
0073   edm::Handle<DTRecSegment4DCollection> all4DSegments;
0074   if (doSegmentVeto) {
0075     e.getByToken(recHits4DToken_, all4DSegments);
0076 
0077     // Loop over all chambers containing a segment and look for the number of segments
0078     DTRecSegment4DCollection::id_iterator chamberId;
0079     for (chamberId = all4DSegments->id_begin(); chamberId != all4DSegments->id_end(); ++chamberId) {
0080       segmentsChId[*chamberId] = 1;
0081     }
0082   }
0083 
0084   // Get the digis from the event
0085   edm::Handle<DTDigiCollection> dtdigis;
0086   e.getByToken(dtDigiToken_, dtdigis);
0087 
0088   // LOOP OVER ALL THE DIGIS OF THE EVENT
0089   DTDigiCollection::DigiRangeIterator dtLayerId_It;
0090   for (dtLayerId_It = dtdigis->begin(); dtLayerId_It != dtdigis->end(); ++dtLayerId_It) {
0091     for (DTDigiCollection::const_iterator digiIt = ((*dtLayerId_It).second).first;
0092          digiIt != ((*dtLayerId_It).second).second;
0093          ++digiIt) {
0094       //Check the TDC trigger width
0095       int tdcTime = (*digiIt).countsTDC();
0096       double upperLimit = tTrigStMap[(*dtLayerId_It).first.superlayerId().chamberId()] - safeMargin;
0097       if (doTimeBoxHistos)
0098         tbHistos[(*dtLayerId_It).first.superlayerId()]->Fill(tdcTime);
0099       if (tdcTime > upperLimit)
0100         continue;
0101 
0102       //Check the chamber has no 4D segments (optional)
0103       if (doSegmentVeto && segmentsChId.find((*dtLayerId_It).first.superlayerId().chamberId()) != segmentsChId.end())
0104         continue;
0105 
0106       // fill the occupancy histo
0107       // FIXME: needs to be optimized: no need to rescale the histo for each digi
0108       TH2F* noise_root = noiseHistos[(*dtLayerId_It).first.superlayerId().chamberId()]->getTH2F();
0109       double normalization = 0;
0110       if (mapEvt.find((*dtLayerId_It).first.superlayerId().chamberId()) != mapEvt.end()) {
0111         LogVerbatim("DTNoiseTask") << " Last fill: # of events: "
0112                                    << mapEvt[(*dtLayerId_It).first.superlayerId().chamberId()] << endl;
0113         normalization = 1e-9 * upperLimit * mapEvt[(*dtLayerId_It).first.superlayerId().chamberId()];
0114         // revert back to # of entries
0115         noise_root->Scale(normalization);
0116       }
0117       int yBin = (*dtLayerId_It).first.layer() + (4 * ((*dtLayerId_It).first.superlayerId().superlayer() - 1));
0118       noise_root->Fill((*digiIt).wire(), yBin);
0119       // normalize the occupancy histo
0120       mapEvt[(*dtLayerId_It).first.superlayerId().chamberId()] = evtNumber;
0121       LogVerbatim("DTNoiseTask") << (*dtLayerId_It).first << " wire: " << (*digiIt).wire()
0122                                  << " # counts: " << noise_root->GetBinContent((*digiIt).wire(), yBin)
0123                                  << " Time interval: " << upperLimit << " # of events: " << evtNumber << endl;
0124       ;
0125       normalization = double(1e-9 * upperLimit * mapEvt[(*dtLayerId_It).first.superlayerId().chamberId()]);
0126       // update the rate
0127       noise_root->Scale(1. / normalization);
0128       LogVerbatim("DTNoiseTask") << "    noise rate: " << noise_root->GetBinContent((*digiIt).wire(), yBin) << endl;
0129     }
0130   }
0131 }
0132 
0133 void DTNoiseTask::bookHistos(DQMStore::IBooker& ibooker, DTChamberId chId) {
0134   // set the folder
0135   stringstream wheel;
0136   wheel << chId.wheel();
0137   stringstream station;
0138   station << chId.station();
0139   stringstream sector;
0140   sector << chId.sector();
0141 
0142   ibooker.setCurrentFolder("DT/05-Noise/Wheel" + wheel.str() +
0143                            //           "/Station" + station.str() +
0144                            "/Sector" + sector.str());
0145 
0146   // Build the histo name
0147   string histoName = string("NoiseRate") + "_W" + wheel.str() + "_St" + station.str() + "_Sec" + sector.str();
0148 
0149   LogVerbatim("DTNoiseTask") << "[DTNoiseTask]: booking chamber histo:" << endl;
0150   LogVerbatim("DTNoiseTask") << "              folder "
0151                              << "DT/05-Noise/Wheel" + wheel.str() +
0152                                     //     "/Station" + station.str() +
0153                                     "/Sector" + sector.str() + "/"
0154                              << endl;
0155   LogVerbatim("DTNoiseTask") << "              histoName " << histoName << endl;
0156 
0157   // Get the chamber from the geometry
0158   int nWires_max = 0;
0159   const DTChamber* dtchamber = dtGeom->chamber(chId);
0160   const vector<const DTSuperLayer*>& superlayers = dtchamber->superLayers();
0161 
0162   // Loop over layers and find the max # of wires
0163   for (vector<const DTSuperLayer*>::const_iterator sl = superlayers.begin(); sl != superlayers.end();
0164        ++sl) {  // loop over SLs
0165     vector<const DTLayer*> layers = (*sl)->layers();
0166     for (vector<const DTLayer*>::const_iterator lay = layers.begin(); lay != layers.end(); ++lay) {  // loop over layers
0167       int nWires = (*lay)->specificTopology().channels();
0168       if (nWires > nWires_max)
0169         nWires_max = nWires;
0170     }
0171   }
0172 
0173   noiseHistos[chId] =
0174       ibooker.book2D(histoName, "Noise rate (Hz) per channel", nWires_max, 1, nWires_max + 1, 12, 1, 13);
0175   noiseHistos[chId]->setAxisTitle("wire number", 1);
0176   noiseHistos[chId]->setBinLabel(1, "SL1-L1", 2);
0177   noiseHistos[chId]->setBinLabel(2, "SL1-L2", 2);
0178   noiseHistos[chId]->setBinLabel(3, "SL1-L3", 2);
0179   noiseHistos[chId]->setBinLabel(4, "SL1-L4", 2);
0180   noiseHistos[chId]->setBinLabel(5, "SL2-L1", 2);
0181   noiseHistos[chId]->setBinLabel(6, "SL2-L2", 2);
0182   noiseHistos[chId]->setBinLabel(7, "SL2-L3", 2);
0183   noiseHistos[chId]->setBinLabel(8, "SL2-L4", 2);
0184   noiseHistos[chId]->setBinLabel(9, "SL3-L1", 2);
0185   noiseHistos[chId]->setBinLabel(10, "SL3-L2", 2);
0186   noiseHistos[chId]->setBinLabel(11, "SL3-L3", 2);
0187   noiseHistos[chId]->setBinLabel(12, "SL3-L4", 2);
0188 }
0189 
0190 void DTNoiseTask::bookHistos(DQMStore::IBooker& ibooker, DTSuperLayerId slId) {
0191   // set the folder
0192   stringstream wheel;
0193   wheel << slId.chamberId().wheel();
0194   stringstream station;
0195   station << slId.chamberId().station();
0196   stringstream sector;
0197   sector << slId.chamberId().sector();
0198   stringstream superlayer;
0199   superlayer << slId.superlayer();
0200 
0201   ibooker.setCurrentFolder("DT/05-Noise/Wheel" + wheel.str() + "/Station" + station.str() + "/Sector" + sector.str());
0202 
0203   // Build the histo name
0204   string histoName =
0205       string("TimeBox") + "_W" + wheel.str() + "_St" + station.str() + "_Sec" + sector.str() + "_SL" + superlayer.str();
0206 
0207   LogVerbatim("DTNoiseTask") << "[DTNoiseTask]: booking SL histo:" << endl;
0208   LogVerbatim("DTNoiseTask") << "              folder "
0209                              << "DT/05-Noise/Wheel" + wheel.str() + "/Station" + station.str() + "/Sector" +
0210                                     sector.str() + "/"
0211                              << endl;
0212   LogVerbatim("DTNoiseTask") << "              histoName " << histoName << endl;
0213 
0214   tbHistos[slId] = ibooker.book1D(histoName, "Time Box (TDC counts)", 1000, 0, 6000);
0215 }
0216 
0217 void DTNoiseTask::dqmBeginRun(const Run& run, const EventSetup& setup) {
0218   LogVerbatim("DTNoiseTask") << "[DTNoiseTask]: Begin of run" << endl;
0219 
0220   // tTrig Map
0221   tTrigMap = &setup.getData(tTrigMapToken_);
0222 
0223   // get the geometry
0224   dtGeom = &setup.getData(muonGeomToken_);
0225 }
0226 
0227 void DTNoiseTask::bookHistograms(DQMStore::IBooker& ibooker, edm::Run const& run, edm::EventSetup const& setup) {
0228   ibooker.setCurrentFolder("DT/EventInfo/Counters");
0229   nEventMonitor = ibooker.bookFloat("nProcessedEventsNoise");
0230 
0231   // Loop over all the chambers
0232   vector<const DTChamber*>::const_iterator ch_it = dtGeom->chambers().begin();
0233   vector<const DTChamber*>::const_iterator ch_end = dtGeom->chambers().end();
0234   for (; ch_it != ch_end; ++ch_it) {
0235     DTChamberId chId = (*ch_it)->id();
0236     // histo booking
0237     bookHistos(ibooker, chId);
0238     vector<const DTSuperLayer*>::const_iterator sl_it = (*ch_it)->superLayers().begin();
0239     vector<const DTSuperLayer*>::const_iterator sl_end = (*ch_it)->superLayers().end();
0240     // Loop over the SLs
0241     for (; sl_it != sl_end; ++sl_it) {
0242       DTSuperLayerId slId = (*sl_it)->id();
0243       if (doTimeBoxHistos)
0244         bookHistos(ibooker, slId);
0245       float tTrig, tTrigRMS, kFactor;
0246       tTrigMap->get(slId, tTrig, tTrigRMS, kFactor, DTTimeUnits::ns);
0247       // tTrig mapping per station
0248       // check that the ttrig is the lowest of the 3 SLs
0249       if (tTrigStMap.find(chId) == tTrigStMap.end() ||
0250           (tTrigStMap.find(chId) != tTrigStMap.end() && tTrig < tTrigStMap[chId]))
0251         tTrigStMap[chId] = tTrig;
0252     }
0253   }
0254 }
0255 
0256 void DTNoiseTask::endLuminosityBlock(const LuminosityBlock& lumiSeg, const EventSetup& setup) {
0257   LogVerbatim("DTNoiseTask") << "[DTNoiseTask]: End LS, update rates in all histos" << endl;
0258 
0259   // update the rate of all histos (usefull for histos with few entries: they are not updated very often
0260   for (map<DTChamberId, MonitorElement*>::const_iterator meAndChamber = noiseHistos.begin();
0261        meAndChamber != noiseHistos.end();
0262        ++meAndChamber) {
0263     DTChamberId chId = (*meAndChamber).first;
0264     TH2F* noise_root = (*meAndChamber).second->getTH2F();
0265     double upperLimit = tTrigStMap[chId] - safeMargin;
0266 
0267     double normalization = 0;
0268     if (mapEvt.find(chId) != mapEvt.end()) {
0269       LogVerbatim("DTNoiseTask") << " Ch: " << chId << " Last fill: # of events: " << mapEvt[chId] << endl;
0270       normalization = 1e-9 * upperLimit * mapEvt[chId];
0271       // revert back to # of entries
0272       noise_root->Scale(normalization);
0273     }
0274     //check that event analyzed != 0 might happen oline
0275     if (evtNumber) {
0276       // set the # of events analyzed until this update
0277       LogVerbatim("DTNoiseTask") << "          Update for events: " << evtNumber << endl;
0278       mapEvt[chId] = evtNumber;
0279       // update the rate
0280       normalization = double(1e-9 * upperLimit * evtNumber);
0281       noise_root->Scale(1. / normalization);
0282     }
0283   }
0284 }
0285 
0286 // Local Variables:
0287 // show-trailing-whitespace: t
0288 // truncate-lines: t
0289 // End: