EcalDigis/plugins/EcalBarrelDigisValidation.cc

0001 /*
0002  * \file EcalBarrelDigisValidation.cc
0003  *
0004  * \author F. Cossutti
0005  *
0006 */
0007
0008 #include "EcalBarrelDigisValidation.h"
0009 #include "CalibCalorimetry/EcalTrivialCondModules/interface/EcalTrivialConditionRetriever.h"
0010 #include "DQMServices/Core/interface/DQMStore.h"
0011
0012 using namespace cms;
0013 using namespace edm;
0014 using namespace std;
0015
0016 EcalBarrelDigisValidation::EcalBarrelDigisValidation(const ParameterSet& ps)
0017     : EBdigiCollection_(consumes<EBDigiCollection>(ps.getParameter<edm::InputTag>("EBdigiCollection"))),
0018       pAgc(esConsumes<edm::Transition::BeginRun>()) {
0019   // verbosity switch
0020   verbose_ = ps.getUntrackedParameter<bool>("verbose", false);
0021
0022   // get hold of back-end interface
0023
0024   gainConv_[1] = 1.;
0025   gainConv_[2] = 2.;
0026   gainConv_[3] = 12.;
0027   gainConv_[0] = 12.;  // saturated channels
0028   barrelADCtoGeV_ = 0.035;
0029   endcapADCtoGeV_ = 0.06;
0030
0031   meEBDigiOccupancy_ = nullptr;
0032
0033   meEBDigiMultiplicity_ = nullptr;
0034
0035   meEBDigiADCGlobal_ = nullptr;
0036
0037   for (int i = 0; i < 10; i++) {
0038     meEBDigiADCAnalog_[i] = nullptr;
0039     meEBDigiADCgS_[i] = nullptr;
0040     meEBDigiADCg1_[i] = nullptr;
0041     meEBDigiADCg6_[i] = nullptr;
0042     meEBDigiADCg12_[i] = nullptr;
0043     meEBDigiGain_[i] = nullptr;
0044   }
0045
0046   meEBPedestal_ = nullptr;
0047
0048   meEBMaximumgt100ADC_ = nullptr;
0049
0050   meEBMaximumgt10ADC_ = nullptr;
0051
0052   meEBnADCafterSwitch_ = nullptr;
0053 }
0054
0055 EcalBarrelDigisValidation::~EcalBarrelDigisValidation() {}
0056
0057 void EcalBarrelDigisValidation::bookHistograms(DQMStore::IBooker& ibooker, edm::Run const&, edm::EventSetup const&) {
0058   Char_t histo[200];
0059
0060   ibooker.setCurrentFolder("EcalDigisV/EcalDigiTask");
0061
0062   sprintf(histo, "EcalDigiTask Barrel occupancy");
0063   meEBDigiOccupancy_ = ibooker.book2D(histo, histo, 360, 0., 360., 170, -85., 85.);
0064
0065   sprintf(histo, "EcalDigiTask Barrel digis multiplicity");
0066   meEBDigiMultiplicity_ = ibooker.book1D(histo, histo, 612, 0., 61200);
0067
0068   sprintf(histo, "EcalDigiTask Barrel global pulse shape");
0069   meEBDigiADCGlobal_ = ibooker.bookProfile(histo, histo, 10, 0, 10, 10000, 0., 1000.);
0070
0071   for (int i = 0; i < 10; i++) {
0072     sprintf(histo, "EcalDigiTask Barrel analog pulse %02d", i + 1);
0073     meEBDigiADCAnalog_[i] = ibooker.book1D(histo, histo, 4000, 0., 400.);
0074
0075     sprintf(histo, "EcalDigiTask Barrel ADC pulse %02d Gain 0 - Saturated", i + 1);
0076     meEBDigiADCgS_[i] = ibooker.book1D(histo, histo, 4096, -0.5, 4095.5);
0077
0078     sprintf(histo, "EcalDigiTask Barrel ADC pulse %02d Gain 1", i + 1);
0079     meEBDigiADCg1_[i] = ibooker.book1D(histo, histo, 4096, -0.5, 4095.5);
0080
0081     sprintf(histo, "EcalDigiTask Barrel ADC pulse %02d Gain 6", i + 1);
0082     meEBDigiADCg6_[i] = ibooker.book1D(histo, histo, 4096, -0.5, 4095.5);
0083
0084     sprintf(histo, "EcalDigiTask Barrel ADC pulse %02d Gain 12", i + 1);
0085     meEBDigiADCg12_[i] = ibooker.book1D(histo, histo, 4096, -0.5, 4095.5);
0086
0087     sprintf(histo, "EcalDigiTask Barrel gain pulse %02d", i + 1);
0088     meEBDigiGain_[i] = ibooker.book1D(histo, histo, 4, 0, 4);
0089   }
0090
0091   sprintf(histo, "EcalDigiTask Barrel pedestal for pre-sample");
0092   meEBPedestal_ = ibooker.book1D(histo, histo, 4096, -0.5, 4095.5);
0093
0094   sprintf(histo, "EcalDigiTask Barrel maximum position gt 100 ADC");
0095   meEBMaximumgt100ADC_ = ibooker.book1D(histo, histo, 10, 0., 10.);
0096
0097   sprintf(histo, "EcalDigiTask Barrel maximum position gt 10 ADC");
0098   meEBMaximumgt10ADC_ = ibooker.book1D(histo, histo, 10, 0., 10.);
0099
0100   sprintf(histo, "EcalDigiTask Barrel ADC counts after gain switch");
0101   meEBnADCafterSwitch_ = ibooker.book1D(histo, histo, 10, 0., 10.);
0102 }
0103
0104 void EcalBarrelDigisValidation::dqmBeginRun(edm::Run const&, edm::EventSetup const& c) { checkCalibrations(c); }
0105
0106 void EcalBarrelDigisValidation::analyze(Event const& e, EventSetup const& c) {
0107   //LogInfo("EventInfo") << " Run = " << e.id().run() << " Event = " << e.id().event();
0108
0109   Handle<EBDigiCollection> EcalDigiEB;
0110
0111   e.getByToken(EBdigiCollection_, EcalDigiEB);
0112
0113   //Return if no Barrel data
0114   if (!EcalDigiEB.isValid())
0115     return;
0116
0117   // BARREL
0118
0119   // loop over Digis
0120
0121   const EBDigiCollection* barrelDigi = EcalDigiEB.product();
0122
0123   std::vector<double> ebAnalogSignal;
0124   std::vector<double> ebADCCounts;
0125   std::vector<double> ebADCGains;
0126   ebAnalogSignal.reserve(EBDataFrame::MAXSAMPLES);
0127   ebADCCounts.reserve(EBDataFrame::MAXSAMPLES);
0128   ebADCGains.reserve(EBDataFrame::MAXSAMPLES);
0129
0130   int nDigis = 0;
0131
0132   for (unsigned int digis = 0; digis < EcalDigiEB->size(); ++digis) {
0133     EBDataFrame ebdf = (*barrelDigi)[digis];
0134     int nrSamples = ebdf.size();
0135
0136     EBDetId ebid = ebdf.id();
0137
0138     nDigis++;
0139     if (meEBDigiOccupancy_)
0140       meEBDigiOccupancy_->Fill(ebid.iphi(), ebid.ieta());
0141
0142     double Emax = 0.;
0143     int Pmax = 0;
0144     double pedestalPreSample = 0.;
0145     double pedestalPreSampleAnalog = 0.;
0146     int countsAfterGainSwitch = -1;
0147     double higherGain = 1.;
0148     int higherGainSample = 0;
0149
0150     for (int sample = 0; sample < nrSamples; ++sample) {
0151       ebAnalogSignal[sample] = 0.;
0152       ebADCCounts[sample] = 0.;
0153       ebADCGains[sample] = 0.;
0154     }
0155
0156     for (int sample = 0; sample < nrSamples; ++sample) {
0157       EcalMGPASample thisSample = ebdf[sample];
0158
0159       ebADCCounts[sample] = (thisSample.adc());
0160       ebADCGains[sample] = (thisSample.gainId());
0161       ebAnalogSignal[sample] = (ebADCCounts[sample] * gainConv_[(int)ebADCGains[sample]] * barrelADCtoGeV_);
0162
0163       if (Emax < ebAnalogSignal[sample]) {
0164         Emax = ebAnalogSignal[sample];
0165         Pmax = sample;
0166       }
0167
0168       if (sample < 3) {
0169         pedestalPreSample += ebADCCounts[sample];
0170         pedestalPreSampleAnalog += ebADCCounts[sample] * gainConv_[(int)ebADCGains[sample]] * barrelADCtoGeV_;
0171       }
0172
0173       if (sample > 0 && (((ebADCGains[sample] > ebADCGains[sample - 1]) && (ebADCGains[sample - 1] != 0)) ||
0174                          (countsAfterGainSwitch < 0 && ebADCGains[sample] == 0))) {
0175         higherGain = ebADCGains[sample];
0176         higherGainSample = sample;
0177         countsAfterGainSwitch = 1;
0178       }
0179
0180       if ((higherGain > 1 && (higherGainSample != sample) && (ebADCGains[sample] == higherGain)) ||
0181           (higherGain == 3 && (higherGainSample != sample) && (ebADCGains[sample] == 0)) ||
0182           (higherGain == 0 && (higherGainSample != sample) &&
0183            ((ebADCGains[sample] == 3) || (ebADCGains[sample] == 0)))) {
0184         countsAfterGainSwitch++;
0185       }
0186     }
0187
0188     pedestalPreSample /= 3.;
0189     pedestalPreSampleAnalog /= 3.;
0190
0191     LogDebug("DigiInfo") << "Barrel Digi for EBDetId = " << ebid.rawId() << " eta,phi " << ebid.ieta() << " "
0192                          << ebid.iphi();
0193     for (int i = 0; i < 10; i++) {
0194       LogDebug("DigiInfo") << "sample " << i << " ADC = " << ebADCCounts[i] << " gain = " << ebADCGains[i]
0195                            << " Analog = " << ebAnalogSignal[i];
0196     }
0197     LogDebug("DigiInfo") << "Maximum energy = " << Emax << " in sample " << Pmax
0198                          << " Pedestal from pre-sample = " << pedestalPreSampleAnalog;
0199     if (countsAfterGainSwitch > 0)
0200       LogDebug("DigiInfo") << "Counts after switch " << countsAfterGainSwitch;
0201
0202     if (countsAfterGainSwitch > 0 && countsAfterGainSwitch < 5) {
0203       edm::LogWarning("DigiWarning") << "Wrong number of counts after gain switch before next switch! "
0204                                      << countsAfterGainSwitch;
0205       for (int i = 0; i < 10; i++) {
0206         edm::LogWarning("DigiWarning") << "sample " << i << " ADC = " << ebADCCounts[i] << " gain = " << ebADCGains[i]
0207                                        << " Analog = " << ebAnalogSignal[i];
0208       }
0209     }
0210
0211     for (int i = 0; i < 10; i++) {
0212       if (meEBDigiADCGlobal_ &&
0213           (Emax - pedestalPreSampleAnalog * gainConv_[(int)ebADCGains[Pmax]]) > 100. * barrelADCtoGeV_)
0214         meEBDigiADCGlobal_->Fill(i, ebAnalogSignal[i]);
0215       if (meEBDigiADCAnalog_[i])
0216         meEBDigiADCAnalog_[i]->Fill(ebAnalogSignal[i]);
0217
0218       if (ebADCGains[i] == 0) {
0219         if (meEBDigiADCgS_[i])
0220           meEBDigiADCgS_[i]->Fill(ebADCCounts[i]);
0221       } else if (ebADCGains[i] == 3) {
0222         if (meEBDigiADCg1_[i])
0223           meEBDigiADCg1_[i]->Fill(ebADCCounts[i]);
0224       } else if (ebADCGains[i] == 2) {
0225         if (meEBDigiADCg6_[i])
0226           meEBDigiADCg6_[i]->Fill(ebADCCounts[i]);
0227       } else if (ebADCGains[i] == 1) {
0228         if (meEBDigiADCg12_[i])
0229           meEBDigiADCg12_[i]->Fill(ebADCCounts[i]);
0230       }
0231       if (meEBDigiGain_[i])
0232         meEBDigiGain_[i]->Fill(ebADCGains[i]);
0233     }
0234
0235     if (meEBPedestal_)
0236       meEBPedestal_->Fill(pedestalPreSample);
0237     if (meEBMaximumgt10ADC_ &&
0238         (Emax - pedestalPreSampleAnalog * gainConv_[(int)ebADCGains[Pmax]]) > 10. * barrelADCtoGeV_)
0239       meEBMaximumgt10ADC_->Fill(Pmax);
0240     if (meEBMaximumgt100ADC_ &&
0241         (Emax - pedestalPreSampleAnalog * gainConv_[(int)ebADCGains[Pmax]]) > 100. * barrelADCtoGeV_)
0242       meEBMaximumgt100ADC_->Fill(Pmax);
0243     if (meEBnADCafterSwitch_)
0244       meEBnADCafterSwitch_->Fill(countsAfterGainSwitch);
0245   }
0246
0247   if (meEBDigiMultiplicity_)
0248     meEBDigiMultiplicity_->Fill(nDigis);
0249 }
0250
0251 void EcalBarrelDigisValidation::checkCalibrations(edm::EventSetup const& eventSetup) {
0252   // ADC -> GeV Scale
0253   const EcalADCToGeVConstant* agc = &eventSetup.getData(pAgc);
0254
0255   EcalMGPAGainRatio* defaultRatios = new EcalMGPAGainRatio();
0256
0257   gainConv_[1] = 1.;
0258   gainConv_[2] = defaultRatios->gain12Over6();
0259   gainConv_[3] = gainConv_[2] * (defaultRatios->gain6Over1());
0260   gainConv_[0] = gainConv_[2] * (defaultRatios->gain6Over1());
0261
0262   LogDebug("EcalDigi") << " Gains conversions: "
0263                        << "\n"
0264                        << " g0 = " << gainConv_[0] << "\n"
0265                        << " g1 = " << gainConv_[1] << "\n"
0266                        << " g2 = " << gainConv_[2] << "\n"
0267                        << " g3 = " << gainConv_[3];
0268
0269   delete defaultRatios;
0270
0271   const double barrelADCtoGeV_ = agc->getEBValue();
0272   LogDebug("EcalDigi") << " Barrel GeV/ADC = " << barrelADCtoGeV_;
0273   const double endcapADCtoGeV_ = agc->getEEValue();
0274   LogDebug("EcalDigi") << " Endcap GeV/ADC = " << endcapADCtoGeV_;
0275 }