CMS 3D CMS Logo

/data/refman/pasoursint/CMSSW_4_1_8_patch12/src/DQM/SiPixelMonitorDigi/src/SiPixelDigiModule.cc

Go to the documentation of this file.
00001 #include "DQM/SiPixelMonitorDigi/interface/SiPixelDigiModule.h"
00002 #include "DQMServices/Core/interface/DQMStore.h"
00003 #include "DQM/SiPixelCommon/interface/SiPixelHistogramId.h"
00005 #include "FWCore/ServiceRegistry/interface/Service.h"
00006 // STL
00007 #include <vector>
00008 #include <memory>
00009 #include <string>
00010 #include <iostream>
00011 #include <stdlib.h>
00012 #include <sstream>
00013 #include <cstdio>
00014 
00015 // Data Formats
00016 #include "DataFormats/SiPixelDetId/interface/PixelBarrelName.h"
00017 #include "DataFormats/SiPixelDetId/interface/PixelEndcapName.h"
00018 #include "DataFormats/DetId/interface/DetId.h"
00019 #include "DataFormats/SiPixelDetId/interface/PixelSubdetector.h"
00020 
00021 //
00022 // Constructors
00023 //
00024 SiPixelDigiModule::SiPixelDigiModule() : id_(0),
00025                                          ncols_(416),
00026                                          nrows_(160) 
00027 {
00028 }
00030 SiPixelDigiModule::SiPixelDigiModule(const uint32_t& id) : 
00031   id_(id),
00032   ncols_(416),
00033   nrows_(160)
00034 { 
00035 }
00037 SiPixelDigiModule::SiPixelDigiModule(const uint32_t& id, const int& ncols, const int& nrows) : 
00038   id_(id),
00039   ncols_(ncols),
00040   nrows_(nrows)
00041 { 
00042 }
00043 //
00044 // Destructor
00045 //
00046 SiPixelDigiModule::~SiPixelDigiModule() {}
00047 //
00048 // Book histograms
00049 //
00050 void SiPixelDigiModule::book(const edm::ParameterSet& iConfig, int type, bool twoD, bool hiRes, bool reducedSet, bool additInfo) {
00051   bool barrel = DetId(id_).subdetId() == static_cast<int>(PixelSubdetector::PixelBarrel);
00052   bool endcap = DetId(id_).subdetId() == static_cast<int>(PixelSubdetector::PixelEndcap);
00053   bool isHalfModule = false;
00054   if(barrel){
00055     isHalfModule = PixelBarrelName(DetId(id_)).isHalfModule(); 
00056   }
00057 
00058   std::string hid;
00059   // Get collection name and instantiate Histo Id builder
00060   edm::InputTag src = iConfig.getParameter<edm::InputTag>( "src" );
00061   
00062 
00063   // Get DQM interface
00064   DQMStore* theDMBE = edm::Service<DQMStore>().operator->();
00065 
00066   int nbinx=ncols_/2, nbiny=nrows_/2;
00067   std::string twodtitle = "Number of Digis (1bin=four pixels)"; 
00068   std::string pxtitle = "Number of Digis (1bin=two columns)";
00069   std::string pytitle = "Number of Digis (1bin=two rows)";
00070   if(hiRes){
00071     nbinx = ncols_;
00072     nbiny = nrows_;
00073     twodtitle = "Number of Digis (1bin=one pixel)";
00074     pxtitle = "Number of Digis (1bin=one column)";
00075     pytitle = "Number of Digis (1bin=one row)";
00076   }
00077   
00078   if(type==0){
00079     SiPixelHistogramId* theHistogramId = new SiPixelHistogramId( src.label() );
00080     // Number of digis
00081     hid = theHistogramId->setHistoId("ndigis",id_);
00082     meNDigis_ = theDMBE->book1D(hid,"Number of Digis",25,0.,25.);
00083     meNDigis_->setAxisTitle("Number of digis",1);
00084     // Charge in ADC counts
00085     hid = theHistogramId->setHistoId("adc",id_);
00086     meADC_ = theDMBE->book1D(hid,"Digi charge",128,0.,256.);
00087     meADC_->setAxisTitle("ADC counts",1);
00088         if(!reducedSet)
00089         {
00090     if(twoD){
00091       if(additInfo){
00092         // 2D hit map
00093         hid = theHistogramId->setHistoId("hitmap",id_);
00094         mePixDigis_ = theDMBE->book2D(hid,twodtitle,nbinx,0.,float(ncols_),nbiny,0.,float(nrows_));
00095         mePixDigis_->setAxisTitle("Columns",1);
00096         mePixDigis_->setAxisTitle("Rows",2);
00097       }
00098     }
00099     else{
00100       // projections of 2D hit map
00101       hid = theHistogramId->setHistoId("hitmap",id_);
00102       mePixDigis_px_ = theDMBE->book1D(hid+"_px",pxtitle,nbinx,0.,float(ncols_));
00103       mePixDigis_py_ = theDMBE->book1D(hid+"_py",pytitle,nbiny,0.,float(nrows_));
00104       mePixDigis_px_->setAxisTitle("Columns",1);
00105       mePixDigis_py_->setAxisTitle("Rows",1);
00106     }
00107         }
00108     delete theHistogramId;
00109 
00110   }
00111   
00112   if(type==1 && barrel){
00113     uint32_t DBladder = PixelBarrelName(DetId(id_)).ladderName();
00114     char sladder[80]; sprintf(sladder,"Ladder_%02i",DBladder);
00115     hid = src.label() + "_" + sladder;
00116     if(isHalfModule) hid += "H";
00117     else hid += "F";
00118     // Number of digis
00119     meNDigisLad_ = theDMBE->book1D("ndigis_"+hid,"Number of Digis",25,0.,25.);
00120     meNDigisLad_->setAxisTitle("Number of digis",1);
00121     // Charge in ADC counts
00122     meADCLad_ = theDMBE->book1D("adc_" + hid,"Digi charge",128,0.,256.);
00123     meADCLad_->setAxisTitle("ADC counts",1);
00124         if(!reducedSet)
00125         {
00126     if(twoD){
00127       // 2D hit map
00128       mePixDigisLad_ = theDMBE->book2D("hitmap_"+hid,twodtitle,nbinx,0.,float(ncols_),nbiny,0.,float(nrows_));
00129       mePixDigisLad_->setAxisTitle("Columns",1);
00130       mePixDigisLad_->setAxisTitle("Rows",2);
00131     }
00132     else{
00133       // projections of 2D hit map
00134       mePixDigisLad_px_ = theDMBE->book1D("hitmap_"+hid+"_px",pxtitle,nbinx,0.,float(ncols_));
00135       mePixDigisLad_py_ = theDMBE->book1D("hitmap_"+hid+"_py",pytitle,nbiny,0.,float(nrows_));
00136       mePixDigisLad_px_->setAxisTitle("Columns",1);
00137       mePixDigisLad_py_->setAxisTitle("Rows",1);        
00138     }
00139         }
00140   }
00141   if(type==2 && barrel){
00142     uint32_t DBlayer = PixelBarrelName(DetId(id_)).layerName();
00143     char slayer[80]; sprintf(slayer,"Layer_%i",DBlayer);
00144     hid = src.label() + "_" + slayer;
00145     if(!additInfo){
00146       // Number of digis
00147       meNDigisLay_ = theDMBE->book1D("ndigis_"+hid,"Number of Digis",25,0.,25.);
00148       meNDigisLay_->setAxisTitle("Number of digis",1);
00149     // Charge in ADC counts
00150       meADCLay_ = theDMBE->book1D("adc_" + hid,"Digi charge",128,0.,256.);
00151       meADCLay_->setAxisTitle("ADC counts",1);
00152     }
00153     if(!reducedSet){
00154       if(twoD || additInfo){
00155         // 2D hit map
00156         if(isHalfModule){
00157           mePixDigisLay_ = theDMBE->book2D("hitmap_"+hid,twodtitle,nbinx,0.,float(ncols_),2*nbiny,0.,float(2*nrows_));
00158         }
00159         else{
00160           mePixDigisLay_ = theDMBE->book2D("hitmap_"+hid,twodtitle,nbinx,0.,float(ncols_),nbiny,0.,float(nrows_));
00161         }
00162         mePixDigisLay_->setAxisTitle("Columns",1);
00163         mePixDigisLay_->setAxisTitle("Rows",2);
00164       }
00165       if(!twoD && !additInfo){
00166         // projections of 2D hit map
00167         mePixDigisLay_px_ = theDMBE->book1D("hitmap_"+hid+"_px",pxtitle,nbinx,0.,float(ncols_));
00168         if(isHalfModule){
00169           mePixDigisLay_py_ = theDMBE->book1D("hitmap_"+hid+"_py",pytitle,2*nbiny,0.,float(2*nrows_));
00170         }
00171         else{
00172           mePixDigisLay_py_ = theDMBE->book1D("hitmap_"+hid+"_py",pytitle,nbiny,0.,float(nrows_));
00173         }
00174         mePixDigisLay_px_->setAxisTitle("Columns",1);
00175         mePixDigisLay_py_->setAxisTitle("Rows",1);
00176       }
00177     }
00178   }
00179   if(type==3 && barrel){
00180     uint32_t DBmodule = PixelBarrelName(DetId(id_)).moduleName();
00181     char smodule[80]; sprintf(smodule,"Ring_%i",DBmodule);
00182     hid = src.label() + "_" + smodule;
00183     // Number of digis
00184     meNDigisPhi_ = theDMBE->book1D("ndigis_"+hid,"Number of Digis",25,0.,25.);
00185     meNDigisPhi_->setAxisTitle("Number of digis",1);
00186     // Charge in ADC counts
00187     meADCPhi_ = theDMBE->book1D("adc_" + hid,"Digi charge",128,0.,256.);
00188     meADCPhi_->setAxisTitle("ADC counts",1);
00189     if(!reducedSet)
00190       {
00191         if(twoD){
00192           
00193           // 2D hit map
00194           if(isHalfModule){
00195             mePixDigisPhi_ = theDMBE->book2D("hitmap_"+hid,twodtitle,nbinx,0.,float(ncols_),2*nbiny,0.,float(2*nrows_));
00196           }
00197           else {
00198             mePixDigisPhi_ = theDMBE->book2D("hitmap_"+hid,twodtitle,nbinx,0.,float(ncols_),nbiny,0.,float(nrows_));
00199           }
00200           mePixDigisPhi_->setAxisTitle("Columns",1);
00201           mePixDigisPhi_->setAxisTitle("Rows",2);
00202         }
00203         else{
00204           // projections of 2D hit map
00205           mePixDigisPhi_px_ = theDMBE->book1D("hitmap_"+hid+"_px",pxtitle,nbinx,0.,float(ncols_));
00206           if(isHalfModule){
00207             mePixDigisPhi_py_ = theDMBE->book1D("hitmap_"+hid+"_py",pytitle,2*nbiny,0.,float(2*nrows_));
00208           }
00209           else{
00210             mePixDigisPhi_py_ = theDMBE->book1D("hitmap_"+hid+"_py",pytitle,nbiny,0.,float(nrows_));
00211           }
00212           mePixDigisPhi_px_->setAxisTitle("Columns",1);
00213           mePixDigisPhi_py_->setAxisTitle("Rows",1);
00214         }
00215       }
00216   }
00217   if(type==4 && endcap){
00218     uint32_t blade= PixelEndcapName(DetId(id_)).bladeName();
00219     
00220     char sblade[80]; sprintf(sblade, "Blade_%02i",blade);
00221     hid = src.label() + "_" + sblade;
00222     // Number of digis
00223     meNDigisBlade_ = theDMBE->book1D("ndigis_"+hid,"Number of Digis",25,0.,25.);
00224     meNDigisBlade_->setAxisTitle("Number of digis",1);
00225     // Charge in ADC counts
00226     meADCBlade_ = theDMBE->book1D("adc_" + hid,"Digi charge",128,0.,256.);
00227     meADCBlade_->setAxisTitle("ADC counts",1);
00228   }
00229   if(type==5 && endcap){
00230     uint32_t disk = PixelEndcapName(DetId(id_)).diskName();
00231     
00232     char sdisk[80]; sprintf(sdisk, "Disk_%i",disk);
00233     hid = src.label() + "_" + sdisk;
00234     if(!additInfo){
00235       // Number of digis
00236       meNDigisDisk_ = theDMBE->book1D("ndigis_"+hid,"Number of Digis",25,0.,25.);
00237       meNDigisDisk_->setAxisTitle("Number of digis",1);
00238       // Charge in ADC counts
00239       meADCDisk_ = theDMBE->book1D("adc_" + hid,"Digi charge",128,0.,256.);
00240       meADCDisk_->setAxisTitle("ADC counts",1);
00241     }
00242     if(additInfo){
00243       mePixDigisDisk_ = theDMBE->book2D("hitmap_"+hid,twodtitle,260,0.,260.,160,0.,160.);
00244       mePixDigisDisk_->setAxisTitle("Columns",1);
00245       mePixDigisDisk_->setAxisTitle("Rows",2);
00246     }
00247   }
00248   if(type==6 && endcap){
00249     uint32_t panel= PixelEndcapName(DetId(id_)).pannelName();
00250     uint32_t module= PixelEndcapName(DetId(id_)).plaquetteName();
00251     char slab[80]; sprintf(slab, "Panel_%i_Ring_%i",panel, module);
00252     hid = src.label() + "_" + slab;
00253     // Number of digis
00254     meNDigisRing_ = theDMBE->book1D("ndigis_"+hid,"Number of Digis",25,0.,25.);
00255     meNDigisRing_->setAxisTitle("Number of digis",1);
00256     // Charge in ADC counts
00257     meADCRing_ = theDMBE->book1D("adc_" + hid,"Digi charge",128,0.,256.);
00258     meADCRing_->setAxisTitle("ADC counts",1);
00259         if(!reducedSet)
00260         {
00261     if(twoD){
00262       // 2D hit map
00263       mePixDigisRing_ = theDMBE->book2D("hitmap_"+hid,twodtitle,nbinx,0.,float(ncols_),nbiny,0.,float(nrows_));
00264       mePixDigisRing_->setAxisTitle("Columns",1);
00265       mePixDigisRing_->setAxisTitle("Rows",2);
00266     }
00267     else{
00268       // projections of 2D hit map
00269       mePixDigisRing_px_ = theDMBE->book1D("hitmap_"+hid+"_px",pxtitle,nbinx,0.,float(ncols_));
00270       mePixDigisRing_py_ = theDMBE->book1D("hitmap_"+hid+"_py",pytitle,nbiny,0.,float(nrows_));
00271       mePixDigisRing_px_->setAxisTitle("Columns",1);
00272       mePixDigisRing_py_->setAxisTitle("Rows",1);
00273     }
00274         }
00275   }
00276 }
00277 
00278 
00279 //
00280 // Fill histograms
00281 //
00282 int SiPixelDigiModule::fill(const edm::DetSetVector<PixelDigi>& input, bool modon, 
00283                                                                  bool ladon, bool layon, bool phion, 
00284                                                                  bool bladeon, bool diskon, bool ringon, 
00285                                                                  bool twoD, bool reducedSet, bool twoDimModOn, bool twoDimOnlyLayDisk,
00286                                                                  int &nDigisA, int &nDigisB) {
00287   bool barrel = DetId(id_).subdetId() == static_cast<int>(PixelSubdetector::PixelBarrel);
00288   bool endcap = DetId(id_).subdetId() == static_cast<int>(PixelSubdetector::PixelEndcap);
00289   bool isHalfModule = false;
00290   uint32_t DBladder = 0;
00291   if(barrel){
00292     isHalfModule = PixelBarrelName(DetId(id_)).isHalfModule(); 
00293     DBladder = PixelBarrelName(DetId(id_)).ladderName();
00294   }
00295 
00296   // Get DQM interface
00297   DQMStore* theDMBE = edm::Service<DQMStore>().operator->();
00298   //std::cout<<"id_ = "<<id_<<" , dmbe="<<theDMBE->pwd()<<std::endl;
00299   //std::cout<<"********************"<<std::endl;
00300   edm::DetSetVector<PixelDigi>::const_iterator isearch = input.find(id_); // search  digis of detid
00301   
00302   unsigned int numberOfDigisMod = 0;
00303   int numberOfDigis[8]; for(int i=0; i!=8; i++) numberOfDigis[i]=0; 
00304   nDigisA=0; nDigisB=0;
00305   if( isearch != input.end() ) {  // Not an empty iterator
00306     
00307     // Look at digis now
00308     edm::DetSet<PixelDigi>::const_iterator  di;
00309     for(di = isearch->data.begin(); di != isearch->data.end(); di++) {
00310       int adc = di->adc();    // charge
00311       int col = di->column(); // column 
00312       int row = di->row();    // row
00313       numberOfDigisMod++;
00314       bool isHalfModule = false;
00315       uint32_t DBladder = 0;
00316       if(barrel){
00317         isHalfModule = PixelBarrelName(DetId(id_)).isHalfModule(); 
00318         DBladder = PixelBarrelName(DetId(id_)).ladderName();
00319       }
00320       PixelBarrelName::Shell DBshell = PixelBarrelName(DetId(id_)).shell();
00321       int DBlayer = PixelBarrelName(DetId(id_)).layerName();
00322       if(barrel){
00323         if(isHalfModule){
00324           if(DBshell==PixelBarrelName::pI||DBshell==PixelBarrelName::pO){
00325             numberOfDigis[0]++; nDigisA++;
00326             if(DBlayer==1) numberOfDigis[2]++;
00327             if(DBlayer==2) numberOfDigis[3]++;
00328             if(DBlayer==3) numberOfDigis[4]++;
00329           }
00330           if(DBshell==PixelBarrelName::mI||DBshell==PixelBarrelName::mO){
00331             numberOfDigis[1]++; nDigisB++;
00332             if(DBlayer==1) numberOfDigis[5]++;
00333             if(DBlayer==2) numberOfDigis[6]++;
00334             if(DBlayer==3) numberOfDigis[7]++;
00335           }
00336         }else{
00337           if(row<80){
00338             numberOfDigis[0]++; nDigisA++;
00339             if(DBlayer==1) numberOfDigis[2]++;
00340             if(DBlayer==2) numberOfDigis[3]++;
00341             if(DBlayer==3) numberOfDigis[4]++;
00342           }else{ 
00343             numberOfDigis[1]++; nDigisB++;
00344             if(DBlayer==1) numberOfDigis[5]++;
00345             if(DBlayer==2) numberOfDigis[6]++;
00346             if(DBlayer==3) numberOfDigis[7]++;
00347           }
00348         }
00349       }
00350       if(modon){
00351         if(!reducedSet){
00352           if(twoD) {
00353             if(twoDimModOn) (mePixDigis_)->Fill((float)col,(float)row);
00354           }
00355           else {
00356             (mePixDigis_px_)->Fill((float)col);
00357             (mePixDigis_py_)->Fill((float)row);
00358           }
00359         }
00360         (meADC_)->Fill((float)adc);
00361       }
00362       if(ladon && barrel){
00363         (meADCLad_)->Fill((float)adc);
00364         if(!reducedSet){
00365           if(twoD) (mePixDigisLad_)->Fill((float)col,(float)row);
00366           else {
00367           (mePixDigisLad_px_)->Fill((float)col);
00368           (mePixDigisLad_py_)->Fill((float)row);
00369           }
00370         }
00371       }
00372       if((layon || twoDimOnlyLayDisk) && barrel){
00373         if(!twoDimOnlyLayDisk) (meADCLay_)->Fill((float)adc);
00374         if(!reducedSet){
00375           if((layon && twoD) || twoDimOnlyLayDisk){
00376             if(isHalfModule && DBladder==1){
00377               (mePixDigisLay_)->Fill((float)col,(float)row+80);
00378             }
00379             else (mePixDigisLay_)->Fill((float)col,(float)row);
00380           }
00381           if((layon && !twoD) && !twoDimOnlyLayDisk){
00382             (mePixDigisLay_px_)->Fill((float)col);
00383             if(isHalfModule && DBladder==1) {
00384               (mePixDigisLay_py_)->Fill((float)row+80);
00385             }
00386             else (mePixDigisLay_py_)->Fill((float)row);
00387           }
00388         }
00389       }
00390       if(phion && barrel){
00391         (meADCPhi_)->Fill((float)adc);
00392         if(!reducedSet)
00393         {
00394         if(twoD){
00395           if(isHalfModule && DBladder==1){
00396             (mePixDigisPhi_)->Fill((float)col,(float)row+80);
00397           }
00398           else (mePixDigisPhi_)->Fill((float)col,(float)row);
00399         }
00400         else {
00401           (mePixDigisPhi_px_)->Fill((float)col);
00402           if(isHalfModule && DBladder==1) {
00403             (mePixDigisPhi_py_)->Fill((float)row+80);
00404           }
00405           else (mePixDigisPhi_py_)->Fill((float)row);
00406         }
00407         }
00408       }
00409       if(bladeon && endcap){
00410         (meADCBlade_)->Fill((float)adc);
00411       }
00412       if((diskon || twoDimOnlyLayDisk) && endcap){
00413         if(!twoDimOnlyLayDisk) (meADCDisk_)->Fill((float)adc);
00414         if(twoDimOnlyLayDisk){
00415           (mePixDigisDisk_)->Fill((float)col,(float)row);
00416         }
00417       }
00418       if(ringon && endcap){
00419         (meADCRing_)->Fill((float)adc);
00420         if(!reducedSet)
00421         {
00422         if(twoD) (mePixDigisRing_)->Fill((float)col,(float)row);
00423         else {
00424           (mePixDigisRing_px_)->Fill((float)col);
00425           (mePixDigisRing_py_)->Fill((float)row);
00426         }
00427         }
00428       }
00429     }
00430     if(modon) (meNDigis_)->Fill((float)numberOfDigisMod);
00431     if(ladon && barrel) (meNDigisLad_)->Fill((float)numberOfDigisMod);
00432     if(layon && barrel && !twoDimOnlyLayDisk) (meNDigisLay_)->Fill((float)numberOfDigisMod);
00433     if(phion && barrel) (meNDigisPhi_)->Fill((float)numberOfDigisMod);
00434     if(bladeon && endcap) (meNDigisBlade_)->Fill((float)numberOfDigisMod);
00435     if(diskon && endcap && !twoDimOnlyLayDisk) (meNDigisDisk_)->Fill((float)numberOfDigisMod);
00436     if(ringon && endcap) (meNDigisRing_)->Fill((float)numberOfDigisMod);
00437     if(barrel){ 
00438       MonitorElement* me=theDMBE->get("Pixel/Barrel/ALLMODS_ndigisCOMB_Barrel");
00439       if(me) me->Fill((float)numberOfDigisMod);
00440       me=theDMBE->get("Pixel/Barrel/ALLMODS_ndigisCHAN_Barrel");
00441       if(me){ if(numberOfDigis[0]>0) me->Fill((float)numberOfDigis[0]); if(numberOfDigis[1]>0) me->Fill((float)numberOfDigis[1]); }
00442       me=theDMBE->get("Pixel/Barrel/ALLMODS_ndigisCHAN_BarrelL1");
00443       if(me){ if(numberOfDigis[2]>0) me->Fill((float)numberOfDigis[2]); }
00444       me=theDMBE->get("Pixel/Barrel/ALLMODS_ndigisCHAN_BarrelL2");
00445       if(me){ if(numberOfDigis[3]>0) me->Fill((float)numberOfDigis[3]); }
00446       me=theDMBE->get("Pixel/Barrel/ALLMODS_ndigisCHAN_BarrelL3");
00447       if(me){ if(numberOfDigis[4]>0) me->Fill((float)numberOfDigis[4]); }
00448     }else if(endcap){
00449       MonitorElement* me=theDMBE->get("Pixel/Endcap/ALLMODS_ndigisCOMB_Endcap");
00450       if(me) me->Fill((float)numberOfDigisMod);
00451     }
00452   }
00453   
00454   //std::cout<<"numberOfDigis for this module: "<<numberOfDigis<<std::endl;
00455   return numberOfDigisMod;
00456 }