Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-04-06 11:59:38

0001 #include "SiPixelDQMRocLevelAnalyzer.h"
0002 
0003 SiPixelDQMRocLevelAnalyzer::SiPixelDQMRocLevelAnalyzer(const edm::ParameterSet &iConfig) : conf_(iConfig) {
0004   usesResource(TFileService::kSharedResource);
0005 }
0006 
0007 SiPixelDQMRocLevelAnalyzer::~SiPixelDQMRocLevelAnalyzer() = default;
0008 
0009 // ------------ method called to for each event  ------------
0010 void SiPixelDQMRocLevelAnalyzer::analyze(const edm::Event &iEvent, const edm::EventSetup &iSetup) {}
0011 
0012 // ------------ method called once each job just before starting event loop  ------------
0013 void SiPixelDQMRocLevelAnalyzer::beginJob() {
0014   std::string filename = conf_.getUntrackedParameter<std::string>("fileName");
0015   bRS = conf_.getUntrackedParameter<bool>("barrelRocStud");
0016   fRS = conf_.getUntrackedParameter<bool>("endcapRocStud");
0017   bPixelAlive = conf_.getUntrackedParameter<bool>("pixelAliveStudy");
0018   double pixelAliveThresh = conf_.getUntrackedParameter<double>("pixelAliveThreshold");
0019   bool bThreshold = conf_.getUntrackedParameter<bool>("thresholdStudy");
0020   bool bNoise = conf_.getUntrackedParameter<bool>("noiseStudy");
0021   bool bGain = conf_.getUntrackedParameter<bool>("gainStudy");
0022   bool bPedestal = conf_.getUntrackedParameter<bool>("pedestalStudy");
0023 
0024   if (!bRS && !fRS)
0025     return;
0026 
0027   //Open file and get MEs
0028   dbe = edm::Service<DQMStore>().operator->();
0029   dbe->open(filename);
0030   mes = dbe->getAllContents("");
0031   std::cout << "found " << mes.size() << " monitoring elements!" << std::endl;
0032 
0033   //BARREL
0034   if (bRS) {
0035     //PIXELALIVE
0036     if (bPixelAlive) {
0037       bhPixelAlive = fs_->make<TH1F>("bpixAlive", "PixelAliveSummary_mean_Barrel", 11520, 0., 11520.);
0038       bhPixelAlive_dist =
0039           fs_->make<TH1F>("bpixAliveDist", "Mean PixelAliveSummary_mean_Barrel Disbribution", 110, 0., 1.1);
0040     }
0041     //THRESHOLD
0042     if (bThreshold) {
0043       bhThresholdMean = fs_->make<TH1F>("bthreshMean", "ScurveThresholdSummary_mean_Barrel", 11520, 0., 11520.);
0044       bhThresholdMean_dist =
0045           fs_->make<TH1F>("bthreshMeanDist", "Mean ScurveThresholdSummary_mean_Barrel Distribution", 600, 0., 150.);
0046       bhThresholdRMS = fs_->make<TH1F>("bthreshRMS", "ScurveThresholdSummary_RMS_Barrel", 11520, 0., 11520.);
0047       bhThresholdRMS_dist =
0048           fs_->make<TH1F>("bthreshRMSDist", "Mean ScurveThresholdSummary_RMS_Barrel Distribution", 800, 0., 80.);
0049     }
0050     //NOISE
0051     if (bNoise) {
0052       bhNoiseMean = fs_->make<TH1F>("bnoiseMean", "ScurveSigmasSummary_mean_Barrel", 11520, 0., 11520.);
0053       bhNoiseMean_dist =
0054           fs_->make<TH1F>("bnoiseMeanDist", "Mean ScurveSigmasSummary_mean_Barrel Distribution", 256, -2., 6.);
0055       bhNoiseRMS = fs_->make<TH1F>("bnoiseRMS", "ScurveSigmasSummary_RMS_Barrel", 11520, 0., 11520.);
0056       bhNoiseRMS_dist =
0057           fs_->make<TH1F>("bnoiseRMSDist", "Mean ScurveSigmasSummary_RMS_Barrel Distribution", 768, 0., 8.);
0058     }
0059     //GAIN
0060     if (bGain) {
0061       bhGainMean = fs_->make<TH1F>("bgainMean", "ScurveGainSummary_mean_Barrel", 11520, 0., 11520.);
0062       bhGainMean_dist = fs_->make<TH1F>("bgainMeanDist", "Mean ScurveGainSummary_mean_Barrel Distribution", 80, 0., 8.);
0063       bhGainRMS = fs_->make<TH1F>("bgainRMS", "ScurveGainSummary_RMS_Barrel", 11520, 0., 11520.);
0064       bhGainRMS_dist = fs_->make<TH1F>("bgainRMSDist", "Mean ScurveGainSummary_RMS_Barrel Distribution", 100, 0., 10.);
0065     }
0066     //PEDESTAL
0067     if (bPedestal) {
0068       bhPedestalMean = fs_->make<TH1F>("bpedestalMean", "ScurvePedestalSummary_mean_Barrel", 11520, 0., 11520.);
0069       bhPedestalMean_dist =
0070           fs_->make<TH1F>("bpedestalMeanDist", "Mean ScurvePedestalSummary_mean_Barrel Distribution", 600, 0., 300.);
0071       bhPedestalRMS = fs_->make<TH1F>("bpedestalRMS", "ScurvePedestalSummary_RMS_Barrel", 11520, 0., 11520.);
0072       bhPedestalRMS_dist =
0073           fs_->make<TH1F>("bpedestalRMSDist", "Mean ScurvePedestalSummary_RMS_Barrel Distribution", 1000, 0., 100.);
0074     }
0075   }
0076 
0077   //ENDCAP
0078   if (fRS) {
0079     //PIXELALIVE
0080     if (bPixelAlive) {
0081       ehPixelAlive = fs_->make<TH1F>("fpixAlive", "PixelAliveSummary_mean_Endcap", 4320, 0., 4320.);
0082       ehPixelAlive_dist =
0083           fs_->make<TH1F>("fpixAliveDist", "Mean PixelAliveSummary_mean_Endcap Disbribution", 110, 0., 1.1);
0084     }
0085     //THRESHOLD
0086     if (bThreshold) {
0087       ehThresholdMean = fs_->make<TH1F>("fthreshMean", "ScurveThresholdSummary_mean_Endcap", 4320, 0., 4320.);
0088       ehThresholdMean_dist =
0089           fs_->make<TH1F>("fthreshMeanDist", "Mean ScurveThresholdSummary_mean_Endcap Distribution", 600, 0., 150.);
0090       ehThresholdRMS = fs_->make<TH1F>("fthreshRMS", "ScurveThresholdSummary_RMS_Endcap", 4320, 0., 4320.);
0091       ehThresholdRMS_dist =
0092           fs_->make<TH1F>("fthreshRMSDist", "Mean ScurveThresholdSummary_RMS_Endcap Distribution", 800, 0., 80.);
0093     }
0094     //NOISE
0095     if (bNoise) {
0096       ehNoiseMean = fs_->make<TH1F>("fnoiseMean", "ScurveSigmasSummary_mean_Endcap", 4320, 0., 4320.);
0097       ehNoiseMean_dist =
0098           fs_->make<TH1F>("fnoiseMeanDist", "Mean ScurveSigmasSummary_mean_Endcap Distribution", 256, -2., 6.);
0099       ehNoiseRMS = fs_->make<TH1F>("fnoiseRMS", "ScurveSigmasSummary_RMS_Endcap", 4320, 0., 4320.);
0100       ehNoiseRMS_dist =
0101           fs_->make<TH1F>("fnoiseRMSDist", "Mean ScurveSigmasSummary_RMS_Endcap Distribution", 384, 0., 4.);
0102     }
0103     //GAIN
0104     if (bGain) {
0105       ehGainMean = fs_->make<TH1F>("fgainMean", "ScurveGainSummary_mean_Endcap", 4320, 0., 4320.);
0106       ehGainMean_dist =
0107           fs_->make<TH1F>("fgainMeanDist", "Mean ScurveGainSummary_mean_Endcap Distribution", 600, 0., 150.);
0108       ehGainRMS = fs_->make<TH1F>("fgainRMS", "ScurveGainSummary_RMS_Endcap", 4320, 0., 4320.);
0109       ehGainRMS_dist = fs_->make<TH1F>("fgainRMSDist", "Mean ScurveGainSummary_RMS_Endcap Distribution", 800, 0., 80.);
0110     }
0111     //PEDESTAL
0112     if (bPedestal) {
0113       ehPedestalMean = fs_->make<TH1F>("fpedestalMean", "ScurvePedestalSummary_mean_Endcap", 4320, 0., 4320.);
0114       ehPedestalMean_dist =
0115           fs_->make<TH1F>("fpedestalMeanDist", "Mean ScurvePedestalSummary_mean_Endcap Distribution", 600, 0., 150.);
0116       ehPedestalRMS = fs_->make<TH1F>("fpedestalRMS", "ScurvePedestalSummary_RMS_Endcap", 4320, 0., 4320.);
0117       ehPedestalRMS_dist =
0118           fs_->make<TH1F>("fpedestalRMSDist", "Mean ScurvePedestalSummary_RMS_Endcap Distribution", 800, 0., 80.);
0119     }
0120   }
0121   //PIXELALIVE
0122   if (bPixelAlive) {
0123     RocSummary("pixelAlive_siPixelCalibDigis_");
0124     if (bRS) {
0125       FillRocLevelHistos(bhPixelAlive, bhPixelAlive_dist, vbpixCN, vbpixM);
0126 
0127       //print a list of pixels with pixelAlive quantity below pixelAliveThresh
0128       for (unsigned int i = 0; i < vbpixCN.size(); i++) {
0129         if (vbpixM[i] < pixelAliveThresh) {
0130           double temp = vbpixCN[i];
0131           int shell = (int)((temp - 1) / 2880);  //0 mi, 1 mo, 2 pi, 3 po
0132           temp -= shell * 2880;
0133           int lay = 1;
0134           if (temp > 576) {
0135             temp -= 576;
0136             lay++;
0137           }
0138           if (temp > 960) {
0139             temp -= 960;
0140             lay++;
0141           }
0142           int lad = 1;
0143           if (temp > 32) {
0144             temp -= 32;
0145             lad++;
0146           }
0147           while (temp > 64) {
0148             temp -= 64;
0149             lad++;
0150           }
0151           int mod = 1;
0152           int modsize = 16;
0153           if (lad == 1 || (lay == 1 && lad == 10) || (lay == 2 && lad == 16) || (lay == 3 && lad == 22))
0154             modsize = 8;
0155           while (temp > modsize) {
0156             temp -= modsize;
0157             mod++;
0158           }
0159           std::cout << vbpixCN[i] << " " << vbpixM[i] << ":\n";
0160           std::cout << "Shell ";
0161           switch (shell) {
0162             case 0:
0163               std::cout << "mI";
0164               break;
0165             case 1:
0166               std::cout << "mO";
0167               break;
0168             case 2:
0169               std::cout << "pI";
0170               break;
0171             case 3:
0172               std::cout << "pO";
0173               break;
0174           }
0175           std::cout << " Lay" << lay << " Lad" << lad << " Mod" << mod << " Chip" << temp << "\n\n";
0176         }
0177       }
0178     }
0179     if (fRS) {
0180       FillRocLevelHistos(ehPixelAlive, ehPixelAlive_dist, vfpixCN, vfpixM);
0181       //print a list of pixels with pixelAlive quantity below pixelAliveThresh
0182       for (unsigned int i = 0; i < vfpixCN.size(); i++) {
0183         if (vfpixM[i] < pixelAliveThresh) {
0184           double temp = vfpixCN[i];
0185           int hcyl = (int)((temp - 1) / 1080);  //0 mi, 1 mo, 2 pi, 3 po
0186           temp -= hcyl * 1080;
0187           int disk = 1;
0188           if (temp > 540) {
0189             temp -= 540;
0190             disk++;
0191           }
0192           int blade = 1;
0193           while (temp > 45) {
0194             temp -= 45;
0195             blade++;
0196           }
0197           int panel = 1, mod = 1;
0198           if (temp < 22) {
0199             //panel 1
0200             if (temp > 16) {
0201               temp -= 16;
0202               mod = 4;
0203             } else if (temp > 8) {
0204               temp -= 8;
0205               mod = 3;
0206             } else if (temp > 2) {
0207               temp -= 2;
0208               mod = 2;
0209             }
0210           } else {
0211             //panel 2
0212             temp -= 21;
0213             panel++;
0214             if (temp > 14) {
0215               temp -= 14;
0216               mod = 3;
0217             } else if (temp > 6) {
0218               temp -= 6;
0219               mod = 2;
0220             }
0221           }
0222 
0223           std::cout << vfpixCN[i] << " " << vfpixM[i] << ":\n";
0224           std::cout << "HalfCylinder ";
0225           switch (hcyl) {
0226             case 0:
0227               std::cout << "mI";
0228               break;
0229             case 1:
0230               std::cout << "mO";
0231               break;
0232             case 2:
0233               std::cout << "pI";
0234               break;
0235             case 3:
0236               std::cout << "pO";
0237               break;
0238           }
0239           std::cout << " Disk" << disk << " Blade" << blade << " Panel" << panel << " Mod" << mod << " Chip" << temp
0240                     << "\n\n";
0241         }
0242       }
0243     }
0244   }
0245   //THRESHOLD
0246   if (bThreshold) {
0247     vbpixCN.clear();
0248     vbpixM.clear();
0249     vbpixSD.clear();
0250     vfpixCN.clear();
0251     vfpixM.clear();
0252     vfpixSD.clear();
0253 
0254     RocSummary("ScurveThresholds_siPixelCalibDigis_");
0255     if (bRS) {
0256       FillRocLevelHistos(bhThresholdMean, bhThresholdMean_dist, vbpixCN, vbpixM);
0257       FillRocLevelHistos(bhThresholdRMS, bhThresholdRMS_dist, vbpixCN, vbpixSD);
0258     }
0259     if (fRS) {
0260       FillRocLevelHistos(ehThresholdMean, ehThresholdMean_dist, vfpixCN, vfpixM);
0261       FillRocLevelHistos(ehThresholdRMS, ehThresholdRMS_dist, vfpixCN, vfpixSD);
0262     }
0263   }
0264   //NOISE
0265   if (bNoise) {
0266     vbpixCN.clear();
0267     vbpixM.clear();
0268     vbpixSD.clear();
0269     vfpixCN.clear();
0270     vfpixM.clear();
0271     vfpixSD.clear();
0272 
0273     RocSummary("ScurveSigmas_siPixelCalibDigis_");
0274     if (bRS) {
0275       FillRocLevelHistos(bhNoiseMean, bhNoiseMean_dist, vbpixCN, vbpixM);
0276       FillRocLevelHistos(bhNoiseRMS, bhNoiseRMS_dist, vbpixCN, vbpixSD);
0277     }
0278     if (fRS) {
0279       FillRocLevelHistos(ehNoiseMean, ehNoiseMean_dist, vfpixCN, vfpixM);
0280       FillRocLevelHistos(ehNoiseRMS, ehNoiseRMS_dist, vfpixCN, vfpixSD);
0281     }
0282   }
0283   //GAIN
0284   if (bGain) {
0285     vbpixCN.clear();
0286     vbpixM.clear();
0287     vbpixSD.clear();
0288     vfpixCN.clear();
0289     vfpixM.clear();
0290     vfpixSD.clear();
0291 
0292     RocSummary("Gain2d_siPixelCalibDigis_");
0293     if (bRS) {
0294       FillRocLevelHistos(bhGainMean, bhGainMean_dist, vbpixCN, vbpixM);
0295       FillRocLevelHistos(bhGainRMS, bhGainRMS_dist, vbpixCN, vbpixSD);
0296     }
0297     if (fRS) {
0298       FillRocLevelHistos(ehGainMean, ehGainMean_dist, vfpixCN, vfpixM);
0299       FillRocLevelHistos(ehGainRMS, ehGainRMS_dist, vfpixCN, vfpixSD);
0300     }
0301   }
0302   //PEDESTAL
0303   if (bPedestal) {
0304     vbpixCN.clear();
0305     vbpixM.clear();
0306     vbpixSD.clear();
0307     vfpixCN.clear();
0308     vfpixM.clear();
0309     vfpixSD.clear();
0310 
0311     RocSummary("Pedestal2d_siPixelCalibDigis_");
0312     if (bRS) {
0313       FillRocLevelHistos(bhPedestalMean, bhPedestalMean_dist, vbpixCN, vbpixM);
0314       FillRocLevelHistos(bhPedestalRMS, bhPedestalRMS_dist, vbpixCN, vbpixSD);
0315     }
0316     if (fRS) {
0317       FillRocLevelHistos(ehPedestalMean, ehPedestalMean_dist, vfpixCN, vfpixM);
0318       FillRocLevelHistos(ehPedestalRMS, ehPedestalRMS_dist, vfpixCN, vfpixSD);
0319     }
0320   }
0321 }
0322 
0323 // ------------ method called once each job just after ending the event loop  ------------
0324 void SiPixelDQMRocLevelAnalyzer::endJob() {}
0325 
0326 void SiPixelDQMRocLevelAnalyzer::RocSummary(std::string tagname) {
0327   int maxcrow, maxccol;
0328   std::string name, path = "first";
0329   std::string oldPath = "";
0330   int moduleNumber = 0;
0331   int bchipNumber = 0;
0332   int fchipNumber = 0;
0333   bool bbarrel = false, bforward = false, bhalfMod = false, bwasHM = false;
0334   bool bPFFM = false;
0335   bool bMNCS = false;
0336   int panelNumber = -1;
0337 
0338   for (std::vector<MonitorElement *>::const_iterator ime = mes.begin(); ime != mes.end(); ++ime) {
0339     bwasHM = bhalfMod;
0340     //set default values
0341     bMNCS = false;
0342     bbarrel = false;
0343     bforward = false;
0344     bhalfMod = false;
0345     //set name, (old) path
0346     name = (*ime)->getName();
0347     oldPath = path;
0348     path = (*ime)->getPathname();
0349 
0350     //determine module number if any
0351     if (path.find("Module_") < path.size()) {
0352       if (path != oldPath) {
0353         moduleNumber++;
0354         if (moduleNumber != 1)
0355           bMNCS = true;
0356       }
0357     } else {
0358       if (moduleNumber > 0)
0359         bMNCS = true;
0360       moduleNumber = 0;
0361     }
0362 
0363     //find out location (barrel (hm), endcap)
0364     if (path.find("Barrel/") < path.size()) {
0365       bbarrel = true;
0366       if (path.find("H/") < path.size())
0367         bhalfMod = true;
0368     }
0369     if (path.find("Endcap/") < path.size()) {
0370       bforward = true;
0371       panelNumber = -1;
0372       if (path.find("Panel_1") < path.size())
0373         panelNumber = 1;
0374       if (path.find("Panel_2") < path.size())
0375         panelNumber = 2;
0376     }
0377 
0378     //find tagname in histoname
0379     if (name.find(tagname) < name.size())
0380       bPFFM = true;
0381     else {
0382       //adjust chip number if necessary, skip ME
0383       if (bMNCS) {
0384         if (!bPFFM) {
0385           if (bbarrel && bRS) {
0386             if (bPixelAlive) {
0387               int a = 16;
0388               if (bwasHM)
0389                 a = 8;
0390               for (int i = 0; i < a; i++) {
0391                 bchipNumber++;
0392                 vbpixCN.push_back(bchipNumber);
0393                 vbpixM.push_back(0);
0394               }
0395             } else {
0396               if (bwasHM)
0397                 bchipNumber += 8;
0398               else
0399                 bchipNumber += 16;
0400             }
0401           }
0402           if (bforward && fRS) {
0403             int maxcol = 2;
0404             int mod = moduleNumber;
0405             if (mod > 1)
0406               mod--;
0407             else {
0408               if (panelNumber == 1)
0409                 mod = 4;
0410               else
0411                 mod = 3;
0412             }
0413             if (panelNumber == 1 && (mod == 1 || mod == 4))
0414               maxcol = 1;
0415             if (bPixelAlive) {
0416               for (int i = 0; i < maxcol * (mod + panelNumber); i++) {
0417                 fchipNumber++;
0418                 vfpixCN.push_back(fchipNumber);
0419                 vfpixM.push_back(0);
0420               }
0421             } else {
0422               fchipNumber += maxcol * (mod + panelNumber);
0423             }
0424           }
0425         } else
0426           bPFFM = false;
0427       }
0428       continue;
0429     }
0430 
0431     //BARREL ROC LEVEL PLOTS
0432 
0433     if (bbarrel && bRS) {
0434       maxccol = 8, maxcrow = 2;
0435       if (bhalfMod) {
0436         maxcrow = 1;
0437       }
0438 
0439       RocSumOneModule(maxcrow, maxccol, (*ime), vbpixCN, vbpixM, vbpixSD, bchipNumber);
0440     }
0441 
0442     //ENDCAP ROC LEVEL PLOTS
0443     if (bforward && fRS) {
0444       maxccol = moduleNumber + panelNumber;
0445       maxcrow = 2;
0446       if (panelNumber == 1 && (moduleNumber == 1 || moduleNumber == 4))
0447         maxcrow = 1;
0448 
0449       RocSumOneModule(maxcrow, maxccol, (*ime), vfpixCN, vfpixM, vfpixSD, fchipNumber);
0450     }
0451   }
0452 
0453   std::cout << "Number of Chips: b" << bchipNumber << " f" << fchipNumber << " " << tagname << std::endl;
0454 }
0455 
0456 void SiPixelDQMRocLevelAnalyzer::RocSumOneModule(int maxr,
0457                                                  int maxc,
0458                                                  MonitorElement *const &me,
0459                                                  std::vector<double> &vecCN,
0460                                                  std::vector<double> &vecMean,
0461                                                  std::vector<double> &vecSD,
0462                                                  int &chipNumber) {
0463   float temp, sum, nentries;
0464   for (int cr = 0; cr < maxr; cr++) {
0465     for (int cc = 0; cc < maxc; cc++) {
0466       //compute mean of 1 ROC
0467       chipNumber++;
0468       sum = 0;
0469       nentries = 0;
0470       //sum for 1 ROC
0471       for (int c = 1; c < 53; c++) {
0472         for (int r = 1; r < 81; r++) {
0473           temp = me->getBinContent(52 * cc + c, 80 * cr + r);
0474 
0475           if (temp != 0.) {
0476             sum += temp;
0477             nentries++;
0478           }
0479         }
0480       }
0481       if (nentries == 0 && bPixelAlive) {
0482         vecCN.push_back(chipNumber);
0483         vecMean.push_back(0);
0484       }
0485       if (nentries != 0) {
0486         double mean = sum / nentries;
0487         double avsd = 0.;
0488         int ne = 0;
0489         vecCN.push_back(chipNumber);
0490         vecMean.push_back(mean);
0491 
0492         //computing std dev.
0493         for (int c = 1; c < 53; c++) {
0494           for (int r = 1; r < 81; r++) {
0495             temp = me->getBinContent(52 * cc + c, 80 * cr + r);
0496             if (temp != 0) {
0497               avsd += (temp - mean) * (temp - mean);
0498               ne++;
0499             }
0500           }
0501         }
0502         avsd = avsd / ne;
0503         avsd = sqrt(avsd);
0504         vecSD.push_back(avsd);
0505       }
0506     }
0507   }
0508 }
0509 
0510 void SiPixelDQMRocLevelAnalyzer::FillRocLevelHistos(TH1F *hrocdep,
0511                                                     TH1F *hdist,
0512                                                     std::vector<double> &vecx,
0513                                                     std::vector<double> &vecy) {
0514   if (vecx.size() == vecy.size()) {
0515     for (unsigned int i = 0; i < vecx.size(); i++) {
0516       hrocdep->Fill(vecx[i], vecy[i]);
0517       hdist->Fill(vecy[i]);
0518     }
0519   }
0520 }
0521 
0522 // -- define this as a plug-in
0523 DEFINE_FWK_MODULE(SiPixelDQMRocLevelAnalyzer);