File indexing completed on 2024-04-06 12:27:57
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014
0015
0016
0017
0018
0019
0020 #include <memory>
0021 #include <cmath>
0022
0023
0024 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0025 #include "FWCore/Utilities/interface/Exception.h"
0026 #include "RecoTBCalo/HcalPlotter/src/HcalQLPlotAnalAlgos.h"
0027 #include "DataFormats/HcalDetId/interface/HcalElectronicsId.h"
0028 #include "DataFormats/HcalRecHit/interface/HcalCalibRecHit.h"
0029 #include "TH1.h"
0030
0031
0032
0033
0034
0035
0036
0037
0038
0039
0040
0041
0042 HcalQLPlotAnalAlgos::HcalQLPlotAnalAlgos(const char* outputFilename, const edm::ParameterSet& histoParams) {
0043 triggerID_ = HcalQLPlotHistoMgr::UNKNOWN;
0044
0045 mf_ = new TFile(outputFilename, "RECREATE");
0046 histos_ = new HcalQLPlotHistoMgr(mf_, histoParams);
0047 }
0048
0049
0050
0051
0052
0053 void HcalQLPlotAnalAlgos::end(void) { mf_->Write(); }
0054
0055 void HcalQLPlotAnalAlgos::SetEventType(const HcalTBTriggerData& trigd) {
0056 if (trigd.wasInSpillPedestalTrigger() || trigd.wasOutSpillPedestalTrigger() ||
0057 trigd.wasSpillIgnorantPedestalTrigger())
0058 triggerID_ = HcalQLPlotHistoMgr::PEDESTAL;
0059 if (trigd.wasLEDTrigger())
0060 triggerID_ = HcalQLPlotHistoMgr::LED;
0061 if (trigd.wasLaserTrigger())
0062 triggerID_ = HcalQLPlotHistoMgr::LASER;
0063 if (trigd.wasBeamTrigger())
0064 triggerID_ = HcalQLPlotHistoMgr::BEAM;
0065
0066 if (triggerID_ == HcalQLPlotHistoMgr::UNKNOWN) {
0067 edm::LogError("HcalQLPlotAnalAlgos::begin") << "Trigger Type unrecognized, aborting";
0068 std::exception e;
0069 throw e;
0070 }
0071 }
0072
0073 void HcalQLPlotAnalAlgos::processRH(const HBHERecHitCollection& hbherhc, const HBHEDigiCollection& hbhedgc) {
0074 HBHERecHitCollection::const_iterator it;
0075
0076 for (it = hbherhc.begin(); it != hbherhc.end(); it++) {
0077 HcalDetId id(it->id());
0078 HcalElectronicsId eid;
0079 HBHEDigiCollection::const_iterator dit = hbhedgc.find(id);
0080 if (dit != hbhedgc.end())
0081 eid = dit->elecId();
0082 else {
0083 edm::LogWarning("HcalQLPlotAnalAlgos::processRH") << "No digi found for id" << id;
0084 continue;
0085 }
0086
0087 TH1* ehist = histos_->GetAHistogram(id, eid, HcalQLPlotHistoMgr::ENERGY, triggerID_);
0088 if (ehist) {
0089 ehist->Fill(it->energy());
0090 }
0091
0092 TH1* thist = histos_->GetAHistogram(id, eid, HcalQLPlotHistoMgr::TIME, triggerID_);
0093 if (thist) {
0094 thist->Fill(it->time());
0095 }
0096 }
0097 }
0098
0099 void HcalQLPlotAnalAlgos::processRH(const HORecHitCollection& horhc, const HODigiCollection& hodgc) {
0100 HORecHitCollection::const_iterator it;
0101
0102 for (it = horhc.begin(); it != horhc.end(); it++) {
0103 HcalDetId id(it->id());
0104 HcalElectronicsId eid;
0105 HODigiCollection::const_iterator dit = hodgc.find(id);
0106 if (dit != hodgc.end())
0107 eid = dit->elecId();
0108 else {
0109 edm::LogWarning("HcalQLPlotAnalAlgos::processRH") << "No digi found for id" << id;
0110 continue;
0111 }
0112
0113 TH1* ehist = histos_->GetAHistogram(id, eid, HcalQLPlotHistoMgr::ENERGY, triggerID_);
0114 if (ehist) {
0115 ehist->Fill(it->energy());
0116 }
0117
0118 TH1* thist = histos_->GetAHistogram(id, eid, HcalQLPlotHistoMgr::TIME, triggerID_);
0119 if (thist) {
0120 thist->Fill(it->time());
0121 }
0122 }
0123 }
0124
0125 void HcalQLPlotAnalAlgos::processRH(const HFRecHitCollection& hfrhc, const HFDigiCollection& hfdgc) {
0126 HFRecHitCollection::const_iterator it;
0127
0128 for (it = hfrhc.begin(); it != hfrhc.end(); it++) {
0129 HcalDetId id(it->id());
0130 HcalElectronicsId eid;
0131 HFDigiCollection::const_iterator dit = hfdgc.find(id);
0132 if (dit != hfdgc.end())
0133 eid = dit->elecId();
0134 else {
0135 edm::LogWarning("HcalQLPlotAnalAlgos::processRH") << "No digi found for id" << id;
0136 continue;
0137 }
0138
0139 TH1* ehist = histos_->GetAHistogram(id, eid, HcalQLPlotHistoMgr::ENERGY, triggerID_);
0140 if (ehist) {
0141 ehist->Fill(it->energy());
0142 }
0143
0144 TH1* thist = histos_->GetAHistogram(id, eid, HcalQLPlotHistoMgr::TIME, triggerID_);
0145 if (thist) {
0146 thist->Fill(it->time());
0147 }
0148 }
0149 }
0150
0151 void HcalQLPlotAnalAlgos::processDigi(const HBHEDigiCollection& hbhedigic) {
0152 HBHEDigiCollection::const_iterator it;
0153
0154 for (it = hbhedigic.begin(); it != hbhedigic.end(); it++) {
0155 HcalDetId id(it->id());
0156 HcalElectronicsId eid(it->elecId());
0157
0158 TH1* phist = histos_->GetAHistogram(id, eid, HcalQLPlotHistoMgr::PULSE, triggerID_);
0159 if (phist) {
0160 for (int bin = 0; bin < it->size(); bin++)
0161 phist->Fill(bin * 1.0, (*it)[bin].nominal_fC());
0162 }
0163
0164 if (triggerID_ == HcalQLPlotHistoMgr::PEDESTAL) {
0165 phist = histos_->GetAHistogram(id, eid, HcalQLPlotHistoMgr::ADC, triggerID_);
0166 if (phist) {
0167 for (int bin = 0; bin < it->size(); bin++)
0168 phist->Fill((*it)[bin].adc());
0169 }
0170 }
0171 }
0172 }
0173
0174 void HcalQLPlotAnalAlgos::processDigi(const HODigiCollection& hodigic) {
0175 HODigiCollection::const_iterator it;
0176
0177 for (it = hodigic.begin(); it != hodigic.end(); it++) {
0178 HcalDetId id(it->id());
0179 HcalElectronicsId eid(it->elecId());
0180
0181 TH1* phist = histos_->GetAHistogram(id, eid, HcalQLPlotHistoMgr::PULSE, triggerID_);
0182 if (phist) {
0183 for (int bin = 0; bin < it->size(); bin++)
0184 phist->Fill(bin * 1.0, (*it)[bin].nominal_fC());
0185 }
0186
0187 if (triggerID_ == HcalQLPlotHistoMgr::PEDESTAL) {
0188 phist = histos_->GetAHistogram(id, eid, HcalQLPlotHistoMgr::ADC, triggerID_);
0189 if (phist) {
0190 for (int bin = 0; bin < it->size(); bin++)
0191 phist->Fill((*it)[bin].adc());
0192 }
0193 }
0194 }
0195 }
0196
0197 void HcalQLPlotAnalAlgos::processDigi(const HFDigiCollection& hfdigic) {
0198 HFDigiCollection::const_iterator it;
0199
0200 for (it = hfdigic.begin(); it != hfdigic.end(); it++) {
0201 HcalDetId id(it->id());
0202 HcalElectronicsId eid(it->elecId());
0203
0204 TH1* phist = histos_->GetAHistogram(id, eid, HcalQLPlotHistoMgr::PULSE, triggerID_);
0205 if (phist) {
0206 for (int bin = 0; bin < it->size(); bin++)
0207 phist->Fill(bin * 1.0, (*it)[bin].nominal_fC());
0208 }
0209
0210 if (triggerID_ == HcalQLPlotHistoMgr::PEDESTAL) {
0211 phist = histos_->GetAHistogram(id, eid, HcalQLPlotHistoMgr::ADC, triggerID_);
0212 if (phist) {
0213 for (int bin = 0; bin < it->size(); bin++)
0214 phist->Fill((*it)[bin].adc());
0215 }
0216 }
0217 }
0218 }
0219
0220 HcalCalibRecHit HcalQLPlotAnalAlgos::recoCalib(const HcalCalibDataFrame& cdigi, double calibFC2GeV) {
0221 double nominal_ped = (cdigi[0].nominal_fC() + cdigi[1].nominal_fC()) / 2.0;
0222
0223 double totamp = 0.0;
0224 double maxA = -1e99;
0225 int maxI = -1;
0226 for (int i = 0; i < cdigi.size(); i++) {
0227 double ampl = (cdigi[i].nominal_fC() - nominal_ped) * calibFC2GeV;
0228 totamp += ampl;
0229
0230 if (ampl > maxA) {
0231 maxA = ampl;
0232 maxI = i;
0233 }
0234 }
0235
0236 maxA = fabs(maxA);
0237 float t0 = (maxI > 0) ? (fabs((cdigi[maxI - 1].nominal_fC() - nominal_ped)) * calibFC2GeV) : 0.0;
0238 float t2 = fabs((cdigi[maxI + 1].nominal_fC() - nominal_ped) * calibFC2GeV);
0239 float wpksamp = (maxA + 2.0 * t2) / (t0 + maxA + t2);
0240 float time = (maxI - cdigi.presamples() + wpksamp) * 25.0;
0241
0242 return HcalCalibRecHit(cdigi.id(), totamp, time);
0243 }
0244
0245 void HcalQLPlotAnalAlgos::processDigi(const HcalCalibDigiCollection& calibdigic, double calibFC2GeV) {
0246 HcalCalibDigiCollection::const_iterator it;
0247
0248 for (it = calibdigic.begin(); it != calibdigic.end(); it++) {
0249 HcalCalibDetId id(it->id());
0250 HcalElectronicsId eid(it->elecId());
0251
0252 TH1* phist = histos_->GetAHistogram(id, eid, HcalQLPlotHistoMgr::PULSE, triggerID_);
0253 if (phist) {
0254 for (int bin = 0; bin < it->size(); bin++)
0255 phist->Fill(bin * 1.0, (*it)[bin].nominal_fC());
0256 }
0257
0258
0259
0260 HcalCalibRecHit rh = recoCalib(*it, calibFC2GeV);
0261
0262 TH1* ehist = histos_->GetAHistogram(id, eid, HcalQLPlotHistoMgr::ENERGY, triggerID_);
0263 if (ehist) {
0264 ehist->Fill(rh.amplitude());
0265 }
0266
0267 TH1* thist = histos_->GetAHistogram(id, eid, HcalQLPlotHistoMgr::TIME, triggerID_);
0268 if (thist) {
0269 thist->Fill(rh.time());
0270 }
0271 }
0272 }