00001
00002
00003
00004
00005
00006
00007
00008
00009
00010 #include <Validation/EcalDigis/interface/EcalEndcapDigisValidation.h>
00011 #include "CalibCalorimetry/EcalTrivialCondModules/interface/EcalTrivialConditionRetriever.h"
00012 #include "DQMServices/Core/interface/DQMStore.h"
00013
00014 using namespace cms;
00015 using namespace edm;
00016 using namespace std;
00017
00018 EcalEndcapDigisValidation::EcalEndcapDigisValidation(const ParameterSet& ps):
00019 EEdigiCollection_(ps.getParameter<edm::InputTag>("EEdigiCollection"))
00020 {
00021
00022
00023 verbose_ = ps.getUntrackedParameter<bool>("verbose", false);
00024
00025 dbe_ = 0;
00026
00027
00028 dbe_ = Service<DQMStore>().operator->();
00029
00030 if ( dbe_ ) {
00031 if ( verbose_ ) {
00032 dbe_->setVerbose(1);
00033 } else {
00034 dbe_->setVerbose(0);
00035 }
00036 }
00037
00038 if ( dbe_ ) {
00039 if ( verbose_ ) dbe_->showDirStructure();
00040 }
00041
00042 gainConv_[1] = 1.;
00043 gainConv_[2] = 2.;
00044 gainConv_[3] = 12.;
00045 gainConv_[0] = 12.;
00046 barrelADCtoGeV_ = 0.035;
00047 endcapADCtoGeV_ = 0.06;
00048
00049 meEEDigiOccupancyzp_ = 0;
00050 meEEDigiOccupancyzm_ = 0;
00051
00052 meEEDigiMultiplicityzp_ = 0;
00053 meEEDigiMultiplicityzm_ = 0;
00054
00055 meEEDigiADCGlobal_ = 0;
00056
00057 for (int i = 0; i < 10 ; i++ ) {
00058 meEEDigiADCAnalog_[i] = 0;
00059 meEEDigiADCgS_[i] = 0;
00060 meEEDigiADCg1_[i] = 0;
00061 meEEDigiADCg6_[i] = 0;
00062 meEEDigiADCg12_[i] = 0;
00063 meEEDigiGain_[i] = 0;
00064 }
00065
00066 meEEPedestal_ = 0;
00067
00068 meEEMaximumgt100ADC_ = 0;
00069
00070 meEEMaximumgt20ADC_ = 0;
00071
00072 meEEnADCafterSwitch_ = 0;
00073
00074 Char_t histo[200];
00075
00076
00077 if ( dbe_ ) {
00078 dbe_->setCurrentFolder("EcalDigisV/EcalDigiTask");
00079
00080 sprintf (histo, "EcalDigiTask Endcap occupancy z+" ) ;
00081 meEEDigiOccupancyzp_ = dbe_->book2D(histo, histo, 100, 0., 100., 100, 0., 100.);
00082
00083 sprintf (histo, "EcalDigiTask Endcap occupancy z-" ) ;
00084 meEEDigiOccupancyzm_ = dbe_->book2D(histo, histo, 100, 0., 100., 100, 0., 100.);
00085
00086 sprintf (histo, "EcalDigiTask Endcap multiplicity z+" ) ;
00087 meEEDigiMultiplicityzp_ = dbe_->book1D(histo, histo, 100, 0., 7324.);
00088
00089 sprintf (histo, "EcalDigiTask Endcap multiplicity z-" ) ;
00090 meEEDigiMultiplicityzm_ = dbe_->book1D(histo, histo, 100, 0., 7324.);
00091
00092 sprintf (histo, "EcalDigiTask Endcap global pulse shape" ) ;
00093 meEEDigiADCGlobal_ = dbe_->bookProfile(histo, histo, 10, 0, 10, 10000, 0., 1000.) ;
00094
00095 for (int i = 0; i < 10 ; i++ ) {
00096
00097 sprintf (histo, "EcalDigiTask Endcap analog pulse %02d", i+1) ;
00098 meEEDigiADCAnalog_[i] = dbe_->book1D(histo, histo, 4000, 0., 400.);
00099
00100 sprintf (histo, "EcalDigiTask Endcap ADC pulse %02d Gain 0 - Saturated", i+1) ;
00101 meEEDigiADCgS_[i] = dbe_->book1D(histo, histo, 4096, -0.5, 4095.5);
00102
00103 sprintf (histo, "EcalDigiTask Endcap ADC pulse %02d Gain 1", i+1) ;
00104 meEEDigiADCg1_[i] = dbe_->book1D(histo, histo, 4096, -0.5, 4095.5);
00105
00106 sprintf (histo, "EcalDigiTask Endcap ADC pulse %02d Gain 6", i+1) ;
00107 meEEDigiADCg6_[i] = dbe_->book1D(histo, histo, 4096, -0.5, 4095.5);
00108
00109 sprintf (histo, "EcalDigiTask Endcap ADC pulse %02d Gain 12", i+1) ;
00110 meEEDigiADCg12_[i] = dbe_->book1D(histo, histo, 4096, -0.5, 4095.5);
00111
00112 sprintf (histo, "EcalDigiTask Endcap gain pulse %02d", i+1) ;
00113 meEEDigiGain_[i] = dbe_->book1D(histo, histo, 4, 0, 4);
00114 }
00115
00116 sprintf (histo, "EcalDigiTask Endcap pedestal for pre-sample" ) ;
00117 meEEPedestal_ = dbe_->book1D(histo, histo, 4096, -0.5, 4095.5) ;
00118
00119 sprintf (histo, "EcalDigiTask Endcap maximum position gt 100 ADC" ) ;
00120 meEEMaximumgt100ADC_ = dbe_->book1D(histo, histo, 10, 0., 10.) ;
00121
00122 sprintf (histo, "EcalDigiTask Endcap maximum position gt 20 ADC" ) ;
00123 meEEMaximumgt20ADC_ = dbe_->book1D(histo, histo, 10, 0., 10.) ;
00124
00125 sprintf (histo, "EcalDigiTask Endcap ADC counts after gain switch" ) ;
00126 meEEnADCafterSwitch_ = dbe_->book1D(histo, histo, 10, 0., 10.) ;
00127
00128 }
00129
00130 }
00131
00132 EcalEndcapDigisValidation::~EcalEndcapDigisValidation(){
00133
00134 }
00135
00136 void EcalEndcapDigisValidation::beginRun(Run const &, EventSetup const & c){
00137
00138 checkCalibrations(c);
00139
00140 }
00141
00142 void EcalEndcapDigisValidation::endJob(){
00143
00144 }
00145
00146 void EcalEndcapDigisValidation::analyze(Event const & e, EventSetup const & c){
00147
00148
00149
00150 Handle<EEDigiCollection> EcalDigiEE;
00151
00152 e.getByLabel( EEdigiCollection_ , EcalDigiEE );
00153
00154
00155 if( !EcalDigiEE.isValid() ) return;
00156
00157
00158
00159
00160
00161 const EEDigiCollection * endcapDigi = EcalDigiEE.product () ;
00162
00163 std::vector<double> eeAnalogSignal ;
00164 std::vector<double> eeADCCounts ;
00165 std::vector<double> eeADCGains ;
00166 eeAnalogSignal.reserve(EEDataFrame::MAXSAMPLES);
00167 eeADCCounts.reserve(EEDataFrame::MAXSAMPLES);
00168 eeADCGains.reserve(EEDataFrame::MAXSAMPLES);
00169
00170 int nDigiszp = 0;
00171 int nDigiszm = 0;
00172
00173 for (unsigned int digis=0; digis<EcalDigiEE->size(); ++digis) {
00174
00175 EEDataFrame eedf=(*endcapDigi)[digis];
00176 int nrSamples=eedf.size();
00177
00178 EEDetId eeid = eedf.id () ;
00179
00180 if (eeid.zside() > 0 ) {
00181 if (meEEDigiOccupancyzp_) meEEDigiOccupancyzp_->Fill( eeid.ix(), eeid.iy() );
00182 nDigiszp++;
00183 }
00184 else if (eeid.zside() < 0 ) {
00185 if (meEEDigiOccupancyzm_) meEEDigiOccupancyzm_->Fill( eeid.ix(), eeid.iy() );
00186 nDigiszm++;
00187 }
00188
00189 double Emax = 0. ;
00190 int Pmax = 0 ;
00191 double pedestalPreSample = 0.;
00192 double pedestalPreSampleAnalog = 0.;
00193 int countsAfterGainSwitch = -1;
00194 double higherGain = 1.;
00195 int higherGainSample = 0;
00196
00197 for (int sample = 0 ; sample < nrSamples; ++sample) {
00198 eeAnalogSignal[sample] = 0.;
00199 eeADCCounts[sample] = 0.;
00200 eeADCGains[sample] = 0.;
00201 }
00202
00203 for (int sample = 0 ; sample < nrSamples; ++sample) {
00204
00205 EcalMGPASample mySample = eedf[sample];
00206
00207 eeADCCounts[sample] = (mySample.adc());
00208 eeADCGains[sample] = (mySample.gainId()) ;
00209 eeAnalogSignal[sample] = (eeADCCounts[sample]*gainConv_[(int)eeADCGains[sample]]*endcapADCtoGeV_);
00210
00211 if (Emax < eeAnalogSignal[sample] ) {
00212 Emax = eeAnalogSignal[sample] ;
00213 Pmax = sample ;
00214 }
00215
00216 if ( sample < 3 ) {
00217 pedestalPreSample += eeADCCounts[sample] ;
00218 pedestalPreSampleAnalog += eeADCCounts[sample]*gainConv_[(int)eeADCGains[sample]]*endcapADCtoGeV_ ;
00219 }
00220
00221 if (sample > 0 && ( ((eeADCGains[sample] > eeADCGains[sample-1]) && (eeADCGains[sample-1]!=0)) || (countsAfterGainSwitch<0 && eeADCGains[sample]==0)) ){
00222 higherGain = eeADCGains[sample];
00223 higherGainSample = sample;
00224 countsAfterGainSwitch = 1;
00225 }
00226
00227 if ( (higherGain > 1 && (higherGainSample != sample) && (eeADCGains[sample] == higherGain)) || (higherGain==3 && (higherGainSample != sample) && (eeADCGains[sample]==0)) || (higherGain==0 && (higherGainSample != sample) && ((eeADCGains[sample]==0) || (eeADCGains[sample]==3))) ) countsAfterGainSwitch++ ;
00228 }
00229 pedestalPreSample /= 3. ;
00230 pedestalPreSampleAnalog /= 3. ;
00231
00232 LogDebug("DigiInfo") << "Endcap Digi for EEDetId = " << eeid.rawId() << " x,y " << eeid.ix() << " " << eeid.iy() ;
00233 for ( int i = 0; i < 10 ; i++ ) {
00234 LogDebug("DigiInfo") << "sample " << i << " ADC = " << eeADCCounts[i] << " gain = " << eeADCGains[i] << " Analog = " << eeAnalogSignal[i] ;
00235 }
00236 LogDebug("DigiInfo") << "Maximum energy = " << Emax << " in sample " << Pmax << " Pedestal from pre-sample = " << pedestalPreSampleAnalog;
00237 if ( countsAfterGainSwitch > 0 ) LogDebug("DigiInfo") << "Counts after switch " << countsAfterGainSwitch;
00238
00239 if ( countsAfterGainSwitch > 0 && countsAfterGainSwitch < 5 ) {
00240 edm::LogWarning("DigiWarning") << "Wrong number of counts after gain switch before next switch! " << countsAfterGainSwitch ;
00241 for ( int i = 0; i < 10 ; i++ ) {
00242 edm::LogWarning("DigiWarning") << "sample " << i << " ADC = " << eeADCCounts[i] << " gain = " << eeADCGains[i] << " Analog = " << eeAnalogSignal[i];
00243 }
00244 }
00245
00246 for ( int i = 0 ; i < 10 ; i++ ) {
00247 if (meEEDigiADCGlobal_ && (Emax-pedestalPreSampleAnalog*gainConv_[(int)eeADCGains[Pmax]]) > 100.*endcapADCtoGeV_) meEEDigiADCGlobal_->Fill( i , eeAnalogSignal[i] ) ;
00248 if (meEEDigiADCAnalog_[i]) meEEDigiADCAnalog_[i]->Fill( eeAnalogSignal[i] ) ;
00249 if ( eeADCGains[i] == 0 ) {
00250 if (meEEDigiADCgS_[i]) meEEDigiADCgS_[i]->Fill( eeADCCounts[i] ) ;
00251 }
00252 else if ( eeADCGains[i] == 3 ) {
00253 if (meEEDigiADCg1_[i]) meEEDigiADCg1_[i]->Fill( eeADCCounts[i] ) ;
00254 }
00255 else if ( eeADCGains[i] == 2 ) {
00256 if (meEEDigiADCg6_[i]) meEEDigiADCg6_[i]->Fill( eeADCCounts[i] ) ;
00257 }
00258 else if ( eeADCGains[i] == 1 ) {
00259 if (meEEDigiADCg12_[i]) meEEDigiADCg12_[i]->Fill( eeADCCounts[i] ) ;
00260 }
00261 if (meEEDigiGain_[i]) meEEDigiGain_[i]->Fill( eeADCGains[i] ) ;
00262 }
00263
00264 if (meEEPedestal_) meEEPedestal_->Fill ( pedestalPreSample ) ;
00265 if (meEEMaximumgt20ADC_ && (Emax-pedestalPreSampleAnalog*gainConv_[(int)eeADCGains[Pmax]]) > 20.*endcapADCtoGeV_) meEEMaximumgt20ADC_->Fill( Pmax ) ;
00266 if (meEEMaximumgt100ADC_ && (Emax-pedestalPreSampleAnalog*gainConv_[(int)eeADCGains[Pmax]]) > 100.*endcapADCtoGeV_) meEEMaximumgt100ADC_->Fill( Pmax ) ;
00267 if (meEEnADCafterSwitch_) meEEnADCafterSwitch_->Fill(countsAfterGainSwitch);
00268
00269 }
00270
00271 if ( meEEDigiMultiplicityzp_ ) meEEDigiMultiplicityzp_->Fill(nDigiszp);
00272 if ( meEEDigiMultiplicityzm_ ) meEEDigiMultiplicityzm_->Fill(nDigiszm);
00273
00274 }
00275
00276 void EcalEndcapDigisValidation::checkCalibrations(edm::EventSetup const & eventSetup)
00277 {
00278
00279
00280 edm::ESHandle<EcalADCToGeVConstant> pAgc;
00281 eventSetup.get<EcalADCToGeVConstantRcd>().get(pAgc);
00282 const EcalADCToGeVConstant* agc = pAgc.product();
00283
00284 EcalMGPAGainRatio * defaultRatios = new EcalMGPAGainRatio();
00285
00286 gainConv_[1] = 1.;
00287 gainConv_[2] = defaultRatios->gain12Over6() ;
00288 gainConv_[3] = gainConv_[2]*(defaultRatios->gain6Over1()) ;
00289 gainConv_[0] = gainConv_[2]*(defaultRatios->gain6Over1()) ;
00290
00291 LogDebug("EcalDigi") << " Gains conversions: " << "\n" << " g0 = " << gainConv_[0] << "\n" << " g1 = " << gainConv_[1] << "\n" << " g2 = " << gainConv_[2] << "\n" << " g3 = " << gainConv_[3];
00292
00293 delete defaultRatios;
00294
00295 const double barrelADCtoGeV_ = agc->getEBValue();
00296 LogDebug("EcalDigi") << " Barrel GeV/ADC = " << barrelADCtoGeV_;
00297 const double endcapADCtoGeV_ = agc->getEEValue();
00298 LogDebug("EcalDigi") << " Endcap GeV/ADC = " << endcapADCtoGeV_;
00299
00300 }
00301
00302