File indexing completed on 2024-04-06 12:07:34
0001
0002 #include "DQM/HcalTasks/interface/QIE10Task.h"
0003 #include "FWCore/Framework/interface/Run.h"
0004 #include <map>
0005
0006 using namespace hcaldqm;
0007 using namespace hcaldqm::constants;
0008 QIE10Task::QIE10Task(edm::ParameterSet const& ps)
0009 : DQTask(ps), hcalDbServiceToken_(esConsumes<HcalDbService, HcalDbRecord, edm::Transition::BeginRun>()) {
0010
0011 _tagQIE10 = ps.getUntrackedParameter<edm::InputTag>("tagQIE10", edm::InputTag("hcalDigis"));
0012 _tagHF = ps.getUntrackedParameter<edm::InputTag>("tagHF", edm::InputTag("hcalDigis"));
0013 _tokQIE10 = consumes<QIE10DigiCollection>(_tagQIE10);
0014 _tokHF = consumes<HFDigiCollection>(_tagHF);
0015
0016
0017 _cut = ps.getUntrackedParameter<double>("cut", 50.0);
0018 _ped = ps.getUntrackedParameter<int>("ped", 4);
0019 }
0020 void QIE10Task::bookHistograms(DQMStore::IBooker& ib, edm::Run const& r, edm::EventSetup const& es) {
0021 if (_ptype == fLocal)
0022 if (r.runAuxiliary().run() == 1)
0023 return;
0024
0025 DQTask::bookHistograms(ib, r, es);
0026
0027
0028 edm::ESHandle<HcalDbService> dbs = es.getHandle(hcalDbServiceToken_);
0029 _emap = dbs->getHcalMapping();
0030
0031 unsigned int nTS = _ptype == fLocal ? 10 : 6;
0032
0033
0034 unsigned int itr = 0;
0035 for (auto& crate : constants::crateListHF) {
0036 for (unsigned int slot = SLOT_uTCA_MIN; slot <= SLOT_uTCA_MAX; ++slot) {
0037 std::vector<uint32_t> vhashSlot;
0038 vhashSlot.push_back(HcalElectronicsId(crate, slot, FIBER_uTCA_MIN1, FIBERCH_MIN, false).rawId());
0039 _filter_slot[itr].initialize(filter::fPreserver, hashfunctions::fCrateSlot, vhashSlot);
0040
0041 _cShapeCut_EChannel[itr].initialize(_name,
0042 "ShapeCut",
0043 hcaldqm::hashfunctions::fEChannel,
0044 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fTiming_TS),
0045 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fQIE10fC_400000));
0046 _cLETDCTime_EChannel[itr].initialize(_name,
0047 "LETDCTime",
0048 hcaldqm::hashfunctions::fEChannel,
0049 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fTime_ns_250),
0050 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fN, true));
0051 _cLETDCvsTS_EChannel[itr].initialize(_name,
0052 "TDCvsTS",
0053 hcaldqm::hashfunctions::fEChannel,
0054 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fTiming_TS),
0055 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fQIE10TDC_64),
0056 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fN, true));
0057
0058 for (unsigned int j = 0; j < nTS; j++) {
0059 _cLETDCvsADC_EChannel[j][itr].initialize(_name,
0060 "LETDCvsADC",
0061 hcaldqm::hashfunctions::fEChannel,
0062 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fQIE10ADC_256),
0063 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fQIE10TDC_64),
0064 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fN, true));
0065 _cADC_EChannel[j][itr].initialize(_name,
0066 "ADC",
0067 hcaldqm::hashfunctions::fEChannel,
0068 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fQIE10ADC_256),
0069 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fN, true));
0070 _cLETDC_EChannel[j][itr].initialize(_name,
0071 "LETDC",
0072 hcaldqm::hashfunctions::fEChannel,
0073 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fQIE10TDC_64),
0074 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fN, true));
0075 }
0076 ++itr;
0077 }
0078 }
0079
0080 _cShapeCut.initialize(_name,
0081 "ShapeCut",
0082 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fTiming_TS),
0083 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fQIE10fC_400000));
0084 _cLETDCvsADC.initialize(_name,
0085 "LETDCvsADC",
0086 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fQIE10ADC_256),
0087 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fQIE10TDC_64),
0088 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fN, true));
0089 _cLETDCTimevsADC.initialize(_name,
0090 "LETDCTimevsADC",
0091 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fQIE10ADC_256),
0092 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fTime_ns_250),
0093 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fN, true));
0094 _cLETDC.initialize(_name,
0095 "LETDC",
0096 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fQIE10TDC_64),
0097 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fN, true));
0098 _cADC.initialize(_name,
0099 "ADC",
0100 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fQIE10ADC_256),
0101 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fN, true));
0102
0103
0104 _cOccupancy_Crate.initialize(_name,
0105 "Occupancy",
0106 hashfunctions::fCrate,
0107 new quantity::ElectronicsQuantity(quantity::fSlotuTCA),
0108 new quantity::ElectronicsQuantity(quantity::fFiberuTCAFiberCh),
0109 new quantity::ValueQuantity(quantity::fN));
0110 _cOccupancy_CrateSlot.initialize(_name,
0111 "Occupancy",
0112 hashfunctions::fCrateSlot,
0113 new quantity::ElectronicsQuantity(quantity::fFiberuTCA),
0114 new quantity::ElectronicsQuantity(quantity::fFiberCh),
0115 new quantity::ValueQuantity(quantity::fN));
0116
0117
0118 _cOccupancy_depth.initialize(_name,
0119 "Occupancy",
0120 hcaldqm::hashfunctions::fdepth,
0121 new hcaldqm::quantity::DetectorQuantity(hcaldqm::quantity::fieta),
0122 new hcaldqm::quantity::DetectorQuantity(hcaldqm::quantity::fiphi),
0123 new hcaldqm::quantity::ValueQuantity(hcaldqm::quantity::fN));
0124
0125 itr = 0;
0126 for (auto& crate : constants::crateListHF) {
0127 for (unsigned int slot = SLOT_uTCA_MIN; slot <= SLOT_uTCA_MAX; ++slot) {
0128 char aux[100];
0129 sprintf(aux, "/Crate%d_Slot%d", crate, slot);
0130 _cShapeCut_EChannel[itr].book(ib, _emap, _filter_slot[itr], _subsystem, aux);
0131 _cLETDCTime_EChannel[itr].book(ib, _emap, _filter_slot[itr], _subsystem, aux);
0132 _cLETDCvsTS_EChannel[itr].book(ib, _emap, _filter_slot[itr], _subsystem, aux);
0133 for (unsigned int i = 0; i < nTS; i++) {
0134 char aux2[100];
0135 sprintf(aux2, "/Crate%d_Slot%d/TS%d", crate, slot, i);
0136 _cLETDCvsADC_EChannel[i][itr].book(ib, _emap, _filter_slot[itr], _subsystem, aux);
0137 _cLETDC_EChannel[i][itr].book(ib, _emap, _filter_slot[itr], _subsystem, aux);
0138 _cADC_EChannel[i][itr].book(ib, _emap, _filter_slot[itr], _subsystem, aux);
0139 }
0140 ++itr;
0141 }
0142 }
0143
0144 _cShapeCut.book(ib, _subsystem);
0145 _cLETDCvsADC.book(ib, _subsystem);
0146 _cLETDCTimevsADC.book(ib, _subsystem);
0147 _cLETDC.book(ib, _subsystem);
0148 _cADC.book(ib, _subsystem);
0149
0150 _cOccupancy_Crate.book(ib, _emap, _subsystem);
0151 _cOccupancy_CrateSlot.book(ib, _emap, _subsystem);
0152 _cOccupancy_depth.book(ib, _emap, _subsystem);
0153
0154 _ehashmap.initialize(_emap, electronicsmap::fD2EHashMap);
0155 }
0156
0157 void QIE10Task::globalEndLuminosityBlock(edm::LuminosityBlock const& lb, edm::EventSetup const& es) {
0158
0159 DQTask::globalEndLuminosityBlock(lb, es);
0160 }
0161
0162 void QIE10Task::_process(edm::Event const& e, edm::EventSetup const&) {
0163 edm::Handle<QIE10DigiCollection> cqie10;
0164 edm::Handle<HFDigiCollection> chf;
0165 if (!e.getByToken(_tokQIE10, cqie10))
0166 return;
0167 if (!e.getByToken(_tokHF, chf))
0168 _logger.dqmthrow("Collection HFDigiCollection isn't available" + _tagHF.label() + " " + _tagHF.instance());
0169
0170 std::map<uint32_t, QIE10DataFrame> mqie10;
0171 for (uint32_t i = 0; i < cqie10->size(); i++) {
0172 QIE10DataFrame frame = static_cast<QIE10DataFrame>((*cqie10)[i]);
0173 HcalDetId did = frame.detid();
0174 HcalElectronicsId eid = HcalElectronicsId(_ehashmap.lookup(did));
0175 if (did.subdet() != HcalForward) {
0176 continue;
0177 }
0178
0179
0180 int fakecrate = -1;
0181 if (eid.crateId() == 22) {
0182 fakecrate = 0;
0183 } else if (eid.crateId() == 29) {
0184 fakecrate = 1;
0185 } else if (eid.crateId() == 32) {
0186 fakecrate = 2;
0187 } else {
0188
0189 continue;
0190 }
0191 int index = fakecrate * 12 + eid.slot() - 1;
0192
0193
0194 CaloSamples digi_fC = hcaldqm::utilities::loadADC2fCDB<QIE10DataFrame>(_dbService, did, frame);
0195 double sumQ = hcaldqm::utilities::sumQDB<QIE10DataFrame>(_dbService, digi_fC, did, frame, 0, frame.samples() - 1);
0196
0197
0198 _cOccupancy_Crate.fill(eid);
0199 _cOccupancy_CrateSlot.fill(eid);
0200 _cOccupancy_depth.fill(did);
0201
0202
0203 for (int j = 0; j < frame.samples(); j++) {
0204
0205 if (sumQ > _cut) {
0206 double q = hcaldqm::utilities::adc2fCDBMinusPedestal<QIE10DataFrame>(_dbService, digi_fC, did, frame, j);
0207 _cShapeCut_EChannel[index].fill(eid, j, q);
0208 _cShapeCut.fill(j, q);
0209 }
0210
0211 _cLETDCvsADC_EChannel[j][index].fill(eid, frame[j].adc(), frame[j].le_tdc());
0212 _cLETDCvsADC.fill(frame[j].adc(), frame[j].le_tdc());
0213 _cLETDC_EChannel[j][index].fill(eid, frame[j].le_tdc());
0214 _cLETDC.fill(frame[j].le_tdc());
0215 _cADC_EChannel[j][index].fill(eid, frame[j].adc());
0216 _cADC.fill(frame[j].adc());
0217 _cLETDCvsTS_EChannel[index].fill(eid, j, frame[j].le_tdc());
0218
0219
0220 if (frame[j].le_tdc() < 50) {
0221
0222
0223 double time = j * 25. + (frame[j].le_tdc() / 2.);
0224 _cLETDCTime_EChannel[index].fill(eid, time);
0225 _cLETDCTimevsADC.fill(frame[j].adc(), time);
0226 }
0227 }
0228
0229 mqie10[did.rawId()] = frame;
0230 }
0231 }
0232
0233 void QIE10Task::_resetMonitors(hcaldqm::UpdateFreq) {}
0234
0235 DEFINE_FWK_MODULE(QIE10Task);