EcalDigis/plugins/EcalEndcapDigisValidation.cc

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