Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-04-06 12:07:34

0001 #include "DQM/HcalTasks/interface/QIE11Task.h"
0002 #include "FWCore/Framework/interface/Run.h"
0003 
0004 using namespace hcaldqm;
0005 using namespace hcaldqm::constants;
0006 QIE11Task::QIE11Task(edm::ParameterSet const& ps)
0007     : DQTask(ps), hcalDbServiceToken_(esConsumes<HcalDbService, HcalDbRecord, edm::Transition::BeginRun>()) {
0008   //    tags
0009   _tagQIE11 = ps.getUntrackedParameter<edm::InputTag>("tagQIE11", edm::InputTag("hcalDigis"));
0010   _tokQIE11 = consumes<QIE11DigiCollection>(_tagQIE11);
0011 
0012   _taguMN = ps.getUntrackedParameter<edm::InputTag>("taguMN", edm::InputTag("hcalDigis"));
0013   _tokuMN = consumes<HcalUMNioDigi>(_taguMN);
0014 
0015   //    cuts
0016   _cut = ps.getUntrackedParameter<double>("cut", 50.0);
0017   _ped = ps.getUntrackedParameter<int>("ped", 4);
0018   _laserType = ps.getUntrackedParameter<int32_t>("laserType", -1);
0019   _eventType = ps.getUntrackedParameter<int32_t>("eventType", -1);
0020 }
0021 //statis
0022 void QIE11Task::fillDescriptions(edm::ConfigurationDescriptions& _desc) {
0023   edm::ParameterSetDescription desc;
0024 
0025   //from class inheritance
0026   hcaldqm::DQTask::fillPSetDescription(desc);
0027 
0028   desc.addUntracked<std::string>("name", "QIE11Task");
0029   desc.addUntracked<int>("debug", 0);
0030   desc.addUntracked<int>("runkeyVal", 0);
0031   desc.addUntracked<std::string>("runkeyName", "pp_run");
0032   desc.addUntracked<edm::InputTag>("tagQIE11", edm::InputTag("hcalDigis"));
0033   desc.addUntracked<double>("cut", 20);
0034   desc.addUntracked<int>("ped", 4);
0035   desc.addUntracked<int>("laserType", -1);
0036   desc.addUntracked<int>("eventType", -1);
0037 
0038   _desc.addDefault(desc);
0039 }
0040 
0041 /* virtual */ void QIE11Task::bookHistograms(DQMStore::IBooker& ib, edm::Run const& r, edm::EventSetup const& es) {
0042   if (_ptype == fLocal)
0043     if (r.runAuxiliary().run() == 1)
0044       return;
0045 
0046   DQTask::bookHistograms(ib, r, es);
0047 
0048   //    GET WHAT YOU NEED
0049   edm::ESHandle<HcalDbService> dbs = es.getHandle(hcalDbServiceToken_);
0050   _emap = dbs->getHcalMapping();
0051   std::vector<uint32_t> vhashC34;
0052   vhashC34.push_back(HcalElectronicsId(34, 11, FIBER_uTCA_MIN1, FIBERCH_MIN, false).rawId());
0053   vhashC34.push_back(HcalElectronicsId(34, 12, FIBER_uTCA_MIN1, FIBERCH_MIN, false).rawId());
0054   _filter_C34.initialize(filter::fPreserver, hcaldqm::hashfunctions::fCrateSlot, vhashC34);
0055 
0056   std::vector<std::pair<int, int> > timingChannels;
0057   timingChannels.push_back(std::pair<int, int>(28, 63));
0058   timingChannels.push_back(std::pair<int, int>(28, 65));
0059   timingChannels.push_back(std::pair<int, int>(20, 63));
0060   timingChannels.push_back(std::pair<int, int>(20, 65));
0061   for (int iChan = 0; iChan < 4; ++iChan) {
0062     std::vector<uint32_t> vhashTimingChannel;
0063     for (int depth = 1; depth <= 7; ++depth) {
0064       vhashTimingChannel.push_back(
0065           HcalDetId(HcalEndcap, timingChannels[iChan].first, timingChannels[iChan].second, depth));
0066     }
0067     _filter_timingChannels[iChan].initialize(filter::fPreserver, hcaldqm::hashfunctions::fDChannel, vhashTimingChannel);
0068   }
0069 
0070   //    INITIALIZE what you need
0071 
0072   // EChannel plots, online+local only
0073   if (_ptype != fOffline) {
0074     unsigned int itr = 0;
0075     for (unsigned int crate = 34; crate <= 34; ++crate) {
0076       for (unsigned int slot = 11; slot <= 12; ++slot) {
0077         std::vector<uint32_t> vhashSlot;
0078         vhashSlot.push_back(HcalElectronicsId(crate, slot, FIBER_uTCA_MIN1, FIBERCH_MIN, false).rawId());
0079         _filter_slot[itr].initialize(filter::fPreserver, hashfunctions::fCrateSlot, vhashSlot);
0080         _cShapeCut_EChannel[itr].initialize(_name,
0081                                             "ShapeCut",
0082                                             hcaldqm::hashfunctions::fEChannel,
0083                                             new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fTiming_TS),
0084                                             new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fQIE10fC_400000));
0085         _cLETDCvsTS_EChannel[itr].initialize(_name,
0086                                              "LETDCvsTS",
0087                                              hcaldqm::hashfunctions::fEChannel,
0088                                              new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fTiming_TS),
0089                                              new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fQIE10TDC_64),
0090                                              new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fN, true),
0091                                              0);
0092         _cLETDCTime_EChannel[itr].initialize(_name,
0093                                              "LETDCTime",
0094                                              hcaldqm::hashfunctions::fEChannel,
0095                                              new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fTime_ns_250),
0096                                              new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fN, true),
0097                                              0);
0098         for (unsigned int j = 0; j < 10; j++) {
0099           _cLETDCvsADC_EChannel[j][itr].initialize(
0100               _name,
0101               "LETDCvsADC",
0102               hcaldqm::hashfunctions::fEChannel,
0103               new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fQIE10ADC_256),
0104               new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fQIE10TDC_64),
0105               new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fN, true),
0106               0);
0107           _cADC_EChannel[j][itr].initialize(_name,
0108                                             "ADC",
0109                                             hcaldqm::hashfunctions::fEChannel,
0110                                             new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fQIE10ADC_256),
0111                                             new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fN, true),
0112                                             0);
0113           _cLETDC_EChannel[j][itr].initialize(_name,
0114                                               "LETDC",
0115                                               hcaldqm::hashfunctions::fEChannel,
0116                                               new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fQIE10TDC_64),
0117                                               new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fN, true),
0118                                               0);
0119         }
0120         ++itr;
0121       }
0122     }
0123   }
0124   _cShapeCut.initialize(_name,
0125                         "ShapeCut",
0126                         new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fTiming_TS),
0127                         new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fQIE10fC_400000));
0128   _cLETDCvsADC.initialize(_name,
0129                           "LETDCvsADC",
0130                           new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fQIE10ADC_256),
0131                           new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fQIE10TDC_64),
0132                           new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fN, true),
0133                           0);
0134   _cLETDCTimevsADC.initialize(_name,
0135                               "LETDCTimevsADC",
0136                               new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fQIE10ADC_256),
0137                               new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fTime_ns_250),
0138                               new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fN, true),
0139                               0);
0140   _cLETDC.initialize(_name,
0141                      "LETDC",
0142                      new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fQIE10TDC_64),
0143                      new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fN, true),
0144                      0);
0145   _cADC.initialize(_name,
0146                    "ADC",
0147                    new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fQIE10ADC_256),
0148                    new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fN, true),
0149                    0);
0150 
0151   if (_ptype != fOffline) {
0152     unsigned int itr = 0;
0153     std::map<std::pair<unsigned int, unsigned int>, unsigned int> itr_map;
0154     for (unsigned int crate = 34; crate <= 34; ++crate) {
0155       for (unsigned int slot = 11; slot <= 12; ++slot) {
0156         char aux[100];
0157         sprintf(aux, "/Crate%d_Slot%d", crate, slot);
0158         _cShapeCut_EChannel[itr].book(ib, _emap, _filter_slot[itr], _subsystem, aux);
0159         _cLETDCvsTS_EChannel[itr].book(ib, _emap, _filter_slot[itr], _subsystem, aux);
0160         _cLETDCTime_EChannel[itr].book(ib, _emap, _filter_slot[itr], _subsystem, aux);
0161         for (unsigned int j = 0; j < 10; j++) {
0162           char aux2[100];
0163           sprintf(aux2, "/Crate%d_Slot%d/TS%d", crate, slot, j);
0164           _cLETDCvsADC_EChannel[j][itr].book(ib, _emap, _filter_slot[itr], _subsystem, aux2);
0165           _cLETDC_EChannel[j][itr].book(ib, _emap, _filter_slot[itr], _subsystem, aux2);
0166           _cADC_EChannel[j][itr].book(ib, _emap, _filter_slot[itr], _subsystem, aux2);
0167         }
0168         itr_map[std::make_pair(crate, slot)] = itr;
0169         ++itr;
0170       }
0171     }
0172   }
0173   _cShapeCut.book(ib, _subsystem);
0174   _cLETDCvsADC.book(ib, _subsystem);
0175   _cLETDCTimevsADC.book(ib, _subsystem);
0176   _cLETDC.book(ib, _subsystem);
0177   _cADC.book(ib, _subsystem);
0178 
0179   _ehashmap.initialize(_emap, electronicsmap::fD2EHashMap, _filter_C34);
0180 }
0181 
0182 /* virtual */ void QIE11Task::globalEndLuminosityBlock(edm::LuminosityBlock const& lb, edm::EventSetup const& es) {
0183   //    finish
0184   DQTask::globalEndLuminosityBlock(lb, es);
0185 }
0186 
0187 /* virtual */ void QIE11Task::_process(edm::Event const& e, edm::EventSetup const&) {
0188   edm::Handle<QIE11DigiCollection> cqie11;
0189   if (!e.getByToken(_tokQIE11, cqie11))
0190     return;
0191 
0192   for (uint32_t i = 0; i < cqie11->size(); i++) {
0193     QIE11DataFrame frame = static_cast<QIE11DataFrame>((*cqie11)[i]);
0194     DetId did = frame.detid();
0195     if (HcalDetId(frame.detid()).subdet() != HcalEndcap) {
0196       continue;
0197     }
0198 
0199     HcalElectronicsId eid = HcalElectronicsId(_ehashmap.lookup(did));
0200     if (!eid.rawId()) {
0201       continue;
0202     }
0203     int fakecrate = -1;
0204     if (eid.crateId() == 34)
0205       fakecrate = 0;
0206     int index = fakecrate * 12 + (eid.slot() - 10) - 1;
0207 
0208     //  compute the signal, ped subracted
0209     CaloSamples digi_fC = hcaldqm::utilities::loadADC2fCDB<QIE11DataFrame>(_dbService, did, frame);
0210     //      double q = hcaldqm::utilities::aveTS_v10<QIE11DataFrame>(frame,
0211     //          constants::adc2fC[_ped], 0, frame.samples()-1);
0212 
0213     //  iterate thru all TS and fill
0214     for (int j = 0; j < frame.samples(); j++) {
0215       double q_pedsub = hcaldqm::utilities::adc2fCDBMinusPedestal<QIE11DataFrame>(_dbService, digi_fC, did, frame, j);
0216       if (_ptype != fOffline) {
0217         if (index == 0 || index == 1) {
0218           //    shapes are after the cut
0219           _cShapeCut_EChannel[index].fill(eid, j, q_pedsub);
0220           _cLETDCvsTS_EChannel[index].fill(eid, j, frame[j].tdc());
0221 
0222           //    w/o a cut
0223           _cLETDCvsADC_EChannel[j][index].fill(eid, frame[j].adc(), frame[j].tdc());
0224           _cLETDC_EChannel[j][index].fill(eid, frame[j].tdc());
0225           if (frame[j].tdc() < 50) {
0226             // Each TDC count is 0.5 ns.
0227             // tdc == 62 or 63 means value was below or above threshold for whole time slice.
0228             _cLETDCTime_EChannel[index].fill(eid, j * 25. + (frame[j].tdc() / 2.));
0229           }
0230           _cADC_EChannel[j][index].fill(eid, frame[j].adc());
0231         }
0232       }
0233       _cShapeCut.fill(eid, j, q_pedsub);
0234 
0235       _cLETDCvsADC.fill(frame[j].adc(), frame[j].tdc());
0236       if (frame[j].tdc() < 50) {
0237         _cLETDCTimevsADC.fill(frame[j].adc(), j * 25. + (frame[j].tdc() / 2.));
0238       }
0239       _cLETDC.fill(eid, frame[j].tdc());
0240 
0241       _cADC.fill(eid, frame[j].adc());
0242     }
0243   }
0244 }
0245 
0246 /* virtual */ bool QIE11Task::_isApplicable(edm::Event const& e) {
0247   if (_ptype != fOnline || (_laserType < 0 && _eventType < 0))
0248     return true;
0249   else {
0250     //      fOnline mode
0251     edm::Handle<HcalUMNioDigi> cumn;
0252     if (!e.getByToken(_tokuMN, cumn))
0253       return false;
0254 
0255     //      event type check first
0256     int eventType = cumn->eventType();
0257     if (eventType == _eventType)
0258       return true;
0259 
0260     //      check if this analysis task is of the right laser type
0261     int laserType = cumn->valueUserWord(0);
0262     if (laserType == _laserType)
0263       return true;
0264   }
0265 
0266   return false;
0267 }
0268 
0269 /* virtual */ void QIE11Task::_resetMonitors(hcaldqm::UpdateFreq) {}
0270 
0271 DEFINE_FWK_MODULE(QIE11Task);