CMS 3D CMS Logo

/afs/cern.ch/work/a/aaltunda/public/www/CMSSW_6_2_7/src/Validation/EcalDigis/src/EcalBarrelDigisValidation.cc

Go to the documentation of this file.
00001 /*
00002  * \file EcalBarrelDigisValidation.cc
00003  *
00004  * $Date: 2010/01/04 15:10:59 $
00005  * $Revision: 1.20 $
00006  * \author F. Cossutti
00007  *
00008 */
00009 
00010 #include <Validation/EcalDigis/interface/EcalBarrelDigisValidation.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 EcalBarrelDigisValidation::EcalBarrelDigisValidation(const ParameterSet& ps):
00019   EBdigiCollection_(ps.getParameter<edm::InputTag>("EBdigiCollection"))
00020 {
00021   
00022   // verbosity switch
00023   verbose_ = ps.getUntrackedParameter<bool>("verbose", false);
00024     
00025   dbe_ = 0;
00026                                                                                                                                           
00027   // get hold of back-end interface
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.;   // saturated channels
00046   barrelADCtoGeV_ = 0.035;
00047   endcapADCtoGeV_ = 0.06;
00048  
00049   meEBDigiOccupancy_ = 0;
00050 
00051   meEBDigiMultiplicity_ = 0;
00052 
00053   meEBDigiADCGlobal_ = 0;
00054 
00055   for (int i = 0; i < 10 ; i++ ) {
00056     meEBDigiADCAnalog_[i] = 0;
00057     meEBDigiADCgS_[i]  = 0;
00058     meEBDigiADCg1_[i]  = 0;
00059     meEBDigiADCg6_[i]  = 0;
00060     meEBDigiADCg12_[i] = 0;
00061     meEBDigiGain_[i] = 0;
00062   }
00063 
00064   meEBPedestal_ = 0;
00065                                  
00066   meEBMaximumgt100ADC_ = 0; 
00067                                  
00068   meEBMaximumgt10ADC_ = 0; 
00069 
00070   meEBnADCafterSwitch_ = 0;
00071  
00072   Char_t histo[200];
00073  
00074   
00075   if ( dbe_ ) {
00076     dbe_->setCurrentFolder("EcalDigisV/EcalDigiTask");
00077   
00078     sprintf (histo, "EcalDigiTask Barrel occupancy" ) ;
00079     meEBDigiOccupancy_ = dbe_->book2D(histo, histo, 360, 0., 360., 170, -85., 85.);
00080 
00081     sprintf (histo, "EcalDigiTask Barrel digis multiplicity" ) ;
00082     meEBDigiMultiplicity_ = dbe_->book1D(histo, histo, 612, 0., 61200);
00083   
00084     sprintf (histo, "EcalDigiTask Barrel global pulse shape" ) ;
00085     meEBDigiADCGlobal_ = dbe_->bookProfile(histo, histo, 10, 0, 10, 10000, 0., 1000.) ;
00086     
00087     for (int i = 0; i < 10 ; i++ ) {
00088 
00089       sprintf (histo, "EcalDigiTask Barrel analog pulse %02d", i+1) ;
00090       meEBDigiADCAnalog_[i] = dbe_->book1D(histo, histo, 4000, 0., 400.);
00091 
00092       sprintf (histo, "EcalDigiTask Barrel ADC pulse %02d Gain 0 - Saturated", i+1) ;
00093       meEBDigiADCgS_[i] = dbe_->book1D(histo, histo, 4096, -0.5, 4095.5);
00094 
00095       sprintf (histo, "EcalDigiTask Barrel ADC pulse %02d Gain 1", i+1) ;
00096       meEBDigiADCg1_[i] = dbe_->book1D(histo, histo, 4096, -0.5, 4095.5);
00097 
00098       sprintf (histo, "EcalDigiTask Barrel ADC pulse %02d Gain 6", i+1) ;
00099       meEBDigiADCg6_[i] = dbe_->book1D(histo, histo, 4096, -0.5, 4095.5);
00100 
00101       sprintf (histo, "EcalDigiTask Barrel ADC pulse %02d Gain 12", i+1) ;
00102       meEBDigiADCg12_[i] = dbe_->book1D(histo, histo, 4096, -0.5, 4095.5);
00103 
00104       sprintf (histo, "EcalDigiTask Barrel gain pulse %02d", i+1) ;
00105       meEBDigiGain_[i] = dbe_->book1D(histo, histo, 4, 0, 4);
00106 
00107     }
00108     
00109     sprintf (histo, "EcalDigiTask Barrel pedestal for pre-sample" ) ;
00110     meEBPedestal_ = dbe_->book1D(histo, histo, 4096, -0.5, 4095.5) ;
00111 
00112     sprintf (histo, "EcalDigiTask Barrel maximum position gt 100 ADC" ) ;
00113     meEBMaximumgt100ADC_ = dbe_->book1D(histo, histo, 10, 0., 10.) ;
00114 
00115     sprintf (histo, "EcalDigiTask Barrel maximum position gt 10 ADC" ) ;
00116     meEBMaximumgt10ADC_ = dbe_->book1D(histo, histo, 10, 0., 10.) ;
00117 
00118     sprintf (histo, "EcalDigiTask Barrel ADC counts after gain switch" ) ;
00119     meEBnADCafterSwitch_ = dbe_->book1D(histo, histo, 10, 0., 10.) ;
00120 
00121   }
00122  
00123 }
00124 
00125 EcalBarrelDigisValidation::~EcalBarrelDigisValidation(){
00126  
00127 }
00128 
00129 void EcalBarrelDigisValidation::beginRun(Run const &, EventSetup const & c){
00130 
00131   checkCalibrations(c);
00132 
00133 }
00134 
00135 void EcalBarrelDigisValidation::endJob(){
00136 
00137 }
00138 
00139 void EcalBarrelDigisValidation::analyze(Event const & e, EventSetup const & c){
00140 
00141   //LogInfo("EventInfo") << " Run = " << e.id().run() << " Event = " << e.id().event();
00142 
00143   Handle<EBDigiCollection> EcalDigiEB;
00144 
00145   e.getByLabel( EBdigiCollection_ , EcalDigiEB );
00146 
00147   //Return if no Barrel data 
00148   if( !EcalDigiEB.isValid() ) return;
00149 
00150   // BARREL
00151 
00152   // loop over Digis
00153 
00154   const EBDigiCollection * barrelDigi = EcalDigiEB.product () ;
00155 
00156   std::vector<double> ebAnalogSignal ;
00157   std::vector<double> ebADCCounts ;
00158   std::vector<double> ebADCGains ;
00159   ebAnalogSignal.reserve(EBDataFrame::MAXSAMPLES);
00160   ebADCCounts.reserve(EBDataFrame::MAXSAMPLES);
00161   ebADCGains.reserve(EBDataFrame::MAXSAMPLES);
00162 
00163   int nDigis = 0;
00164   
00165   for (unsigned int digis=0; digis<EcalDigiEB->size(); ++digis) 
00166     {
00167       
00168       EBDataFrame ebdf = (*barrelDigi)[digis];
00169       int nrSamples = ebdf.size();
00170       
00171       EBDetId ebid = ebdf.id () ;
00172 
00173       nDigis++;
00174       if (meEBDigiOccupancy_) meEBDigiOccupancy_->Fill( ebid.iphi(), ebid.ieta() );
00175       
00176       double Emax = 0. ;
00177       int Pmax = 0 ;
00178       double pedestalPreSample = 0.;
00179       double pedestalPreSampleAnalog = 0.;
00180       int countsAfterGainSwitch = -1;
00181       double higherGain = 1.;
00182       int higherGainSample = 0;
00183       
00184       for (int sample = 0 ; sample < nrSamples; ++sample) {
00185         ebAnalogSignal[sample] = 0.;
00186         ebADCCounts[sample] = 0.;
00187         ebADCGains[sample] = 0.;
00188       }
00189       
00190       for (int sample = 0 ; sample < nrSamples; ++sample) 
00191         {
00192           EcalMGPASample thisSample = ebdf[sample];
00193           
00194           ebADCCounts[sample] = (thisSample.adc());
00195           ebADCGains[sample]  = (thisSample.gainId());
00196           ebAnalogSignal[sample] = (ebADCCounts[sample]*gainConv_[(int)ebADCGains[sample]]*barrelADCtoGeV_);
00197           
00198           if (Emax < ebAnalogSignal[sample] ) {
00199             Emax = ebAnalogSignal[sample] ;
00200             Pmax = sample ;
00201           }
00202           
00203           if ( sample < 3 ) {
00204             pedestalPreSample += ebADCCounts[sample] ;
00205             pedestalPreSampleAnalog += ebADCCounts[sample]*gainConv_[(int)ebADCGains[sample]]*barrelADCtoGeV_ ;
00206           }
00207           
00208           if ( sample > 0 && ( ((ebADCGains[sample] > ebADCGains[sample-1]) && (ebADCGains[sample-1]!=0)) || (countsAfterGainSwitch<0 && ebADCGains[sample]==0)) ) {
00209             higherGain = ebADCGains[sample];
00210             higherGainSample = sample;
00211             countsAfterGainSwitch = 1;
00212           }
00213           
00214           if ( (higherGain > 1 && (higherGainSample != sample) && (ebADCGains[sample] == higherGain)) || (higherGain==3 && (higherGainSample != sample) && (ebADCGains[sample]==0)) || (higherGain==0 && (higherGainSample != sample) && ((ebADCGains[sample] == 3) || (ebADCGains[sample]==0))) ) { countsAfterGainSwitch++ ; }
00215         }
00216       
00217       pedestalPreSample /= 3. ; 
00218       pedestalPreSampleAnalog /= 3. ; 
00219       
00220       LogDebug("DigiInfo") << "Barrel Digi for EBDetId = " << ebid.rawId() << " eta,phi " << ebid.ieta() << " " << ebid.iphi() ;
00221       for ( int i = 0; i < 10 ; i++ ) {
00222         LogDebug("DigiInfo") << "sample " << i << " ADC = " << ebADCCounts[i] << " gain = " << ebADCGains[i] << " Analog = " << ebAnalogSignal[i];  
00223       }
00224       LogDebug("DigiInfo") << "Maximum energy = " << Emax << " in sample " << Pmax << " Pedestal from pre-sample = " << pedestalPreSampleAnalog;
00225       if ( countsAfterGainSwitch > 0 ) LogDebug("DigiInfo") << "Counts after switch " << countsAfterGainSwitch;
00226       
00227       if ( countsAfterGainSwitch > 0 && countsAfterGainSwitch < 5 ) {
00228         edm::LogWarning("DigiWarning") << "Wrong number of counts after gain switch before next switch! " << countsAfterGainSwitch ;
00229         for ( int i = 0; i < 10 ; i++ ) {
00230           edm::LogWarning("DigiWarning") << "sample " << i << " ADC = " << ebADCCounts[i] << " gain = " << ebADCGains[i] << " Analog = " << ebAnalogSignal[i];
00231         } 
00232       }
00233       
00234       for ( int i = 0 ; i < 10 ; i++ ) {
00235         if (meEBDigiADCGlobal_ && (Emax-pedestalPreSampleAnalog*gainConv_[(int)ebADCGains[Pmax]]) > 100.*barrelADCtoGeV_) meEBDigiADCGlobal_->Fill( i , ebAnalogSignal[i] ) ;
00236         if (meEBDigiADCAnalog_[i]) meEBDigiADCAnalog_[i]->Fill( ebAnalogSignal[i] ) ;
00237         
00238         if ( ebADCGains[i] == 0) {
00239           if (meEBDigiADCgS_[i]) meEBDigiADCgS_[i]->Fill( ebADCCounts[i] ) ;
00240         }
00241         else if ( ebADCGains[i] == 3 ) {
00242           if (meEBDigiADCg1_[i]) meEBDigiADCg1_[i]->Fill( ebADCCounts[i] ) ;
00243         }
00244         else if ( ebADCGains[i] == 2 ) {
00245           if (meEBDigiADCg6_[i]) meEBDigiADCg6_[i]->Fill( ebADCCounts[i] ) ;
00246         }
00247         else if ( ebADCGains[i] == 1 ) {
00248           if (meEBDigiADCg12_[i]) meEBDigiADCg12_[i]->Fill( ebADCCounts[i] ) ;
00249         }
00250         if (meEBDigiGain_[i]) meEBDigiGain_[i]->Fill( ebADCGains[i] ) ;
00251       }
00252       
00253       if (meEBPedestal_) meEBPedestal_->Fill ( pedestalPreSample ) ;
00254       if (meEBMaximumgt10ADC_ && (Emax-pedestalPreSampleAnalog*gainConv_[(int)ebADCGains[Pmax]]) > 10.*barrelADCtoGeV_) meEBMaximumgt10ADC_->Fill( Pmax ) ;
00255       if (meEBMaximumgt100ADC_ && (Emax-pedestalPreSampleAnalog*gainConv_[(int)ebADCGains[Pmax]]) > 100.*barrelADCtoGeV_) meEBMaximumgt100ADC_->Fill( Pmax ) ;
00256       if (meEBnADCafterSwitch_) meEBnADCafterSwitch_->Fill( countsAfterGainSwitch ) ;
00257       
00258     }
00259   
00260   if ( meEBDigiMultiplicity_ ) meEBDigiMultiplicity_->Fill(nDigis);
00261     
00262 }
00263 
00264 void  EcalBarrelDigisValidation::checkCalibrations(edm::EventSetup const & eventSetup) 
00265   {
00266     
00267   // ADC -> GeV Scale
00268   edm::ESHandle<EcalADCToGeVConstant> pAgc;
00269   eventSetup.get<EcalADCToGeVConstantRcd>().get(pAgc);
00270   const EcalADCToGeVConstant* agc = pAgc.product();
00271   
00272   EcalMGPAGainRatio * defaultRatios = new EcalMGPAGainRatio();
00273 
00274   gainConv_[1] = 1.;
00275   gainConv_[2] = defaultRatios->gain12Over6() ;
00276   gainConv_[3] = gainConv_[2]*(defaultRatios->gain6Over1()) ;
00277   gainConv_[0] = gainConv_[2]*(defaultRatios->gain6Over1()) ;
00278 
00279   LogDebug("EcalDigi") << " Gains conversions: " << "\n" << " g0 = " << gainConv_[0] << "\n" << " g1 = " << gainConv_[1] << "\n" << " g2 = " << gainConv_[2]  << "\n" << " g3 = " << gainConv_[3]; 
00280 
00281   delete defaultRatios;
00282 
00283   const double barrelADCtoGeV_  = agc->getEBValue();
00284   LogDebug("EcalDigi") << " Barrel GeV/ADC = " << barrelADCtoGeV_;
00285   const double endcapADCtoGeV_ = agc->getEEValue();
00286   LogDebug("EcalDigi") << " Endcap GeV/ADC = " << endcapADCtoGeV_;
00287 
00288 }
00289 
00290