Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-11-17 01:40:50

0001 // -*- C++ -*-
0002 //
0003 // Package:    SiPixelGainCalibrationReadDQMFile
0004 // Class:      SiPixelGainCalibrationReadDQMFile
0005 //
0006 /**\class SiPixelGainCalibrationReadDQMFile SiPixelGainCalibrationReadDQMFile.cc CalibTracker/SiPixelGainCalibrationReadDQMFile/src/SiPixelGainCalibrationReadDQMFile.cc
0007 
0008  Description: <one line class summary>
0009 
0010  Implementation:
0011      <Notes on implementation>
0012 */
0013 //
0014 // Original Author:  Freya BLEKMAN
0015 //         Created:  Tue Aug  5 16:22:46 CEST 2008
0016 // $Id: SiPixelGainCalibrationReadDQMFile.cc,v 1.7 2010/01/12 11:29:54 rougny Exp $
0017 //
0018 //
0019 
0020 // system include files
0021 #include <memory>
0022 #include <fstream>
0023 #include <sys/stat.h>
0024 
0025 // user include files
0026 #include "CalibTracker/SiPixelESProducers/interface/SiPixelGainCalibrationService.h"
0027 #include "CommonTools/UtilAlgos/interface/TFileService.h"
0028 #include "CondCore/DBOutputService/interface/PoolDBOutputService.h"
0029 #include "CondFormats/SiPixelObjects/interface/SiPixelCalibConfiguration.h"
0030 #include "CondFormats/SiPixelObjects/interface/SiPixelGainCalibration.h"
0031 #include "CondFormats/SiPixelObjects/interface/SiPixelGainCalibrationForHLT.h"
0032 #include "CondFormats/SiPixelObjects/interface/SiPixelGainCalibrationOffline.h"
0033 #include "FWCore/Framework/interface/Event.h"
0034 #include "FWCore/Framework/interface/Frameworkfwd.h"
0035 #include "FWCore/Framework/interface/MakerMacros.h"
0036 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0037 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0038 #include "FWCore/ServiceRegistry/interface/Service.h"
0039 #include "FWCore/ServiceRegistry/interface/Service.h"
0040 #include "Geometry/CommonDetUnit/interface/PixelGeomDetUnit.h"
0041 #include "Geometry/CommonTopologies/interface/PixelTopology.h"
0042 #include "Geometry/Records/interface/TrackerDigiGeometryRecord.h"
0043 #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h"
0044 
0045 // ROOT includes
0046 #include "TDirectory.h"
0047 #include "TFile.h"
0048 #include "TH2F.h"
0049 #include "TKey.h"
0050 #include "TList.h"
0051 #include "TString.h"
0052 #include "TTree.h"
0053 
0054 //
0055 // class decleration
0056 //
0057 
0058 class SiPixelGainCalibrationReadDQMFile : public edm::one::EDAnalyzer<edm::one::SharedResources> {
0059 public:
0060   explicit SiPixelGainCalibrationReadDQMFile(const edm::ParameterSet &);
0061 
0062 private:
0063   void analyze(const edm::Event &, const edm::EventSetup &) final;
0064   // functions added by F.B.
0065   void fillDatabase(const edm::EventSetup &iSetup, TFile *);
0066   std::unique_ptr<TFile> getHistograms();
0067   // ----------member data ---------------------------
0068   std::map<uint32_t, std::map<std::string, TString> > bookkeeper_;
0069   std::map<uint32_t, std::map<double, double> > Meankeeper_;
0070   std::map<uint32_t, std::vector<std::map<int, int> > > noisyPixelsKeeper_;
0071 
0072   const bool appendMode_;
0073   const edm::ESGetToken<TrackerGeometry, TrackerDigiGeometryRecord> pddToken_;
0074   SiPixelGainCalibrationService theGainCalibrationDbInputService_;
0075   std::unique_ptr<TH2F> defaultGain_;
0076   std::unique_ptr<TH2F> defaultPed_;
0077   std::unique_ptr<TH2F> defaultChi2_;
0078   std::unique_ptr<TH2F> defaultFitResult_;
0079   std::unique_ptr<TH1F> meanGainHist_;
0080   std::unique_ptr<TH1F> meanPedHist_;
0081   const std::string record_;
0082   // keep track of lowest and highest vals for range
0083   float gainlow_;
0084   float gainhi_;
0085   float pedlow_;
0086   float pedhi_;
0087   const bool usemeanwhenempty_;
0088   const std::string rootfilestring_;
0089   float gainmax_;
0090   float pedmax_;
0091   const double badchi2_;
0092   const size_t nmaxcols;
0093   const size_t nmaxrows;
0094 };
0095 
0096 void SiPixelGainCalibrationReadDQMFile::fillDatabase(const edm::EventSetup &iSetup, TFile *therootfile) {
0097   // only create when necessary.
0098   // process the minimum and maximum gain & ped values...
0099   edm::LogPrint("SiPixelGainCalibrationReadDQMFile") << "now starting db fill!!!" << std::endl;
0100 
0101   std::map<uint32_t, std::pair<TString, int> > badresults;
0102 
0103   edm::Service<TFileService> fs;
0104 
0105   //   bool usedmean=false;
0106   //   bool isdead=false;
0107   TH1F *VCAL_endpoint =
0108       fs->make<TH1F>("VCAL_endpoint", "value where response = 255 ( x = (255 - ped)/gain )", 256, 0, 256);
0109   TH1F *goodgains = fs->make<TH1F>("goodgains", "gain values", 100, 0, 10);
0110   TH1F *goodpeds = fs->make<TH1F>("goodpeds", "pedestal values", 356, -100, 256);
0111   TH1F *totgains = fs->make<TH1F>("totgains", "gain values", 200, 0, 10);
0112   TH1F *totpeds = fs->make<TH1F>("totpeds", "pedestal values", 356, -100, 256);
0113   TTree *tree = new TTree("tree", "tree");
0114   int detidfortree, rowfortree, colfortree, useddefaultfortree;
0115   float pedfortree, gainfortree, chi2fortree;
0116   tree->Branch("detid", &detidfortree, "detid/I");
0117   tree->Branch("row", &rowfortree, "row/I");
0118   tree->Branch("col", &colfortree, "col/I");
0119   tree->Branch("defaultval", &useddefaultfortree, "defaultval/I");
0120   tree->Branch("ped", &pedfortree, "ped/F");
0121   tree->Branch("gain", &gainfortree, "gain/F");
0122   tree->Branch("chi2", &chi2fortree, "chi2/F");
0123 
0124   //   TH1F *gainPerDetid;
0125   //   TH1F *pedPerDetid;
0126 
0127   size_t ntimes = 0;
0128   edm::LogPrint("SiPixelGainCalibrationReadDQMFile") << "now filling record " << record_ << std::endl;
0129   if (record_ != "SiPixelGainCalibrationForHLTRcd" && record_ != "SiPixelGainCalibrationOfflineRcd") {
0130     edm::LogPrint("SiPixelGainCalibrationReadDQMFile")
0131         << "you passed record " << record_ << ", which I have no idea what to do with!" << std::endl;
0132     return;
0133   }
0134   if (gainlow_ > gainhi_) {
0135     float temp = gainhi_;
0136     gainhi_ = gainlow_;
0137     gainlow_ = temp;
0138   }
0139   if (pedlow_ > pedhi_) {
0140     float temp = pedhi_;
0141     pedhi_ = pedlow_;
0142     pedlow_ = temp;
0143   }
0144   if (gainhi_ > gainmax_)
0145     gainhi_ = gainmax_;
0146   if (pedhi_ > pedmax_)
0147     pedhi_ = pedmax_;
0148   float badpedval = pedlow_ - 200;
0149   float badgainval = gainlow_ - 200;
0150   edm::LogPrint("SiPixelGainCalibrationReadDQMFile") << "now filling db: values: pedlow, hi: " << pedlow_ << ", "
0151                                                      << pedhi_ << " and gainlow, hi: " << gainlow_ << ", " << gainhi_;
0152   float meangain = meanGainHist_->GetMean();
0153   float meanped = meanPedHist_->GetMean();
0154   //  edm::LogPrint("SiPixelGainCalibrationReadDQMFile") << ", and mean gain " << meangain<< ", ped " << meanped ;
0155   //  edm::LogPrint("SiPixelGainCalibrationReadDQMFile") << std::endl;
0156 
0157   // and fill the dummy histos:
0158 
0159   for (size_t icol = 0; icol < nmaxcols; ++icol) {
0160     for (size_t irow = 0; irow < nmaxrows; ++irow) {
0161       defaultGain_->SetBinContent(icol + 1, irow + 1, meangain);
0162       defaultPed_->SetBinContent(icol + 1, irow + 1, meanped);
0163       defaultChi2_->SetBinContent(icol + 1, irow + 1, 1.0);
0164       defaultFitResult_->SetBinContent(icol + 1, irow + 1, 0);
0165     }
0166   }
0167 
0168   SiPixelGainCalibration theGainCalibrationDbInput(pedlow_ * 0.999, pedhi_ * 1.001, gainlow_ * 0.999, gainhi_ * 1.001);
0169   SiPixelGainCalibrationForHLT theGainCalibrationDbInputHLT(
0170       pedlow_ * 0.999, pedhi_ * 1.001, gainlow_ * 0.999, gainhi_ * 1.001);
0171   SiPixelGainCalibrationOffline theGainCalibrationDbInputOffline(
0172       pedlow_ * 0.999, pedhi_ * 1.001, gainlow_ * 0.999, gainhi_ * 1.001);
0173 
0174   uint32_t nchannels = 0;
0175   uint32_t nmodules = 0;
0176   edm::LogPrint("SiPixelGainCalibrationReadDQMFile")
0177       << "now starting loop on detids, there are " << bookkeeper_.size() << " histograms to consider..." << std::endl;
0178   uint32_t detid = 0;
0179   therootfile->cd();
0180   const TrackerGeometry *pDD = &iSetup.getData(pddToken_);
0181   edm::LogInfo("SiPixelCondObjOfflineBuilder") << " There are " << pDD->dets().size() << " detectors" << std::endl;
0182 
0183   int NDetid = 0;
0184   for (TrackerGeometry::DetContainer::const_iterator it = pDD->dets().begin(); it != pDD->dets().end(); it++) {
0185     detid = 0;
0186     if (dynamic_cast<PixelGeomDetUnit const *>((*it)) != nullptr)
0187       detid = ((*it)->geographicalId()).rawId();
0188     if (detid == 0)
0189       continue;
0190     NDetid++;
0191     //if(NDetid>=2) continue;
0192     //if(detid!=344076812) continue;
0193     int badDetId = 0;
0194     //if(detid==302123296||detid==302126596) badDetId=1;;
0195     //if(detid!=302058516) continue;
0196     ntimes = 0;
0197     useddefaultfortree = 0;
0198     edm::LogPrint("SiPixelGainCalibrationReadDQMFile")
0199         << "now creating database object for detid " << detid
0200         << std::endl;  //<< " " << bookkeeper_[detid]["gain_2d"] << " " << bookkeeper_[detid]["ped_2d"] << std::endl; //edm::LogPrint("SiPixelGainCalibrationReadDQMFile")<< " nrows:" << nrows << " ncols: " << ncols << std::endl;
0201 
0202     // Get the module sizes.
0203     TH2F *tempchi2;
0204     TH2F *tempfitresult;
0205     TH2F *tempgain;
0206     TH2F *tempped;
0207     TString tempgainstring;
0208 
0209     if (!badDetId) {
0210       TString tempchi2string = bookkeeper_[detid]["chi2prob_2d"];
0211       tempchi2 = dynamic_cast<TH2F *>(therootfile->Get(tempchi2string));
0212       if (tempchi2 == nullptr || badDetId) {
0213         tempchi2 = defaultChi2_.get();
0214         useddefaultfortree = 1;
0215       }
0216       TString tempfitresultstring = bookkeeper_[detid]["fitresult_2d"];
0217       tempfitresult = dynamic_cast<TH2F *>(therootfile->Get(tempfitresultstring));
0218       if (tempfitresult == nullptr) {
0219         tempfitresult = defaultFitResult_.get();
0220         useddefaultfortree = 1;
0221       }
0222       TString tempgainstring = bookkeeper_[detid]["gain_2d"];
0223       tempgain = dynamic_cast<TH2F *>(therootfile->Get(tempgainstring));
0224       if (tempgain == nullptr) {
0225         //  edm::LogPrint("SiPixelGainCalibrationReadDQMFile") <<"WARNING, gain histo " << bookkeeper_[detid]["gain_2d"] << " does not exist, using default instead" << std::endl;
0226         tempgain = defaultGain_.get();
0227         useddefaultfortree = 1;
0228       }
0229       TString temppedstring = bookkeeper_[detid]["ped_2d"];
0230       tempped = dynamic_cast<TH2F *>(therootfile->Get(temppedstring));
0231       if (tempped == nullptr) {
0232         //      edm::LogPrint("SiPixelGainCalibrationReadDQMFile") <<"WARNING, ped histo " << bookkeeper_[detid]["ped_2d"] << " for detid " << detid << " does not exist, using default instead" << std::endl;
0233         std::pair<TString, int> tempval(tempgainstring, 0);
0234         badresults[detid] = tempval;
0235         tempped = defaultPed_.get();
0236         useddefaultfortree = 1;
0237       }
0238     } else {
0239       tempchi2 = defaultChi2_.get();
0240       tempgain = defaultGain_.get();
0241       tempfitresult = defaultFitResult_.get();
0242       std::pair<TString, int> tempval(tempgainstring, 0);
0243       badresults[detid] = tempval;
0244       tempped = defaultPed_.get();
0245       useddefaultfortree = 1;
0246     }
0247 
0248     const PixelGeomDetUnit *pixDet = dynamic_cast<const PixelGeomDetUnit *>((*it));
0249     const PixelTopology &topol = pixDet->specificTopology();
0250     // Get the module sizes.
0251     size_t nrows = topol.nrows();     // rows in x
0252     size_t ncols = topol.ncolumns();  // cols in y
0253 
0254     //    int nrows=tempgain->GetNbinsY();
0255     //    int ncols=tempgain->GetNbinsX();
0256     //    edm::LogPrint("SiPixelGainCalibrationReadDQMFile") << "next histo " << tempgain->GetTitle() << " has nrow,ncol:" << nrows << ","<< ncols << std::endl;
0257     size_t nrowsrocsplit = theGainCalibrationDbInputHLT.getNumberOfRowsToAverageOver();
0258     if (theGainCalibrationDbInputOffline.getNumberOfRowsToAverageOver() != nrowsrocsplit)
0259       throw cms::Exception("GainCalibration Payload configuration error")
0260           << "[SiPixelGainCalibrationAnalysis::fillDatabase] ERROR the SiPixelGainCalibrationOffline and "
0261              "SiPixelGainCalibrationForHLT database payloads have different settings for the number of rows per roc: "
0262           << theGainCalibrationDbInputHLT.getNumberOfRowsToAverageOver() << "(HLT), "
0263           << theGainCalibrationDbInputOffline.getNumberOfRowsToAverageOver() << "(offline)";
0264     std::vector<char> theSiPixelGainCalibrationPerPixel;
0265     std::vector<char> theSiPixelGainCalibrationPerColumn;
0266     std::vector<char> theSiPixelGainCalibrationGainPerColPedPerPixel;
0267 
0268     //Get mean of gain/pedestal of this Detid
0269     meangain = meanGainHist_->GetMean();
0270     meanped = meanPedHist_->GetMean();
0271     int npix = 0;
0272     double meanGainForThisModule = 0;
0273     double meanPedForThisModule = 0;
0274     for (size_t icol = 1; icol <= ncols; icol++) {
0275       for (size_t jrow = 1; jrow <= nrows; jrow++) {
0276         if (tempfitresult->GetBinContent(icol, jrow) > 0) {
0277           npix++;
0278           meanGainForThisModule += tempgain->GetBinContent(icol, jrow);
0279           meanPedForThisModule += tempped->GetBinContent(icol, jrow);
0280         }
0281       }
0282     }
0283     if (npix != 0)
0284       meanPedForThisModule /= npix;
0285     if (npix != 0)
0286       meanGainForThisModule /= npix;
0287     if (usemeanwhenempty_) {
0288       if (meanGainForThisModule > gainlow_ && meanGainForThisModule < gainhi_ && npix > 100)
0289         meangain = meanGainForThisModule;
0290       if (meanPedForThisModule > pedlow_ && meanPedForThisModule < pedhi_ && npix > 100)
0291         meanped = meanPedForThisModule;
0292     }
0293 
0294     // Loop over columns and rows of this DetID
0295     float peds[160];
0296     float gains[160];
0297     float pedforthiscol[2];
0298     float gainforthiscol[2];
0299     int nusedrows[2];
0300     size_t nemptypixels = 0;
0301     //edm::LogPrint("SiPixelGainCalibrationReadDQMFile") << pedlow_<<"  "<<pedhi_<<"  "<<gainlow_<<"  "<<gainhi_<< std::endl;
0302     for (size_t icol = 1; icol <= ncols; icol++) {
0303       nusedrows[0] = nusedrows[1] = 0;
0304       pedforthiscol[0] = pedforthiscol[1] = 0;
0305       gainforthiscol[0] = gainforthiscol[1] = 0;
0306       //   edm::LogPrint("SiPixelGainCalibrationReadDQMFile") << "now lookign at col " << icol << std::endl;
0307       for (size_t jrow = 1; jrow <= nrows; jrow++) {
0308         size_t iglobalrow = 0;
0309         if (jrow > nrowsrocsplit)
0310           iglobalrow = 1;
0311         peds[jrow] = badpedval;
0312         gains[jrow] = badgainval;
0313         float ped = tempped->GetBinContent(icol, jrow);
0314         float gain = tempgain->GetBinContent(icol, jrow);
0315         float chi2 = tempchi2->GetBinContent(icol, jrow);
0316         float fitresult = tempfitresult->GetBinContent(icol, jrow);
0317 
0318         if (ped > pedlow_ && gain > gainlow_ && ped < pedhi_ && gain < gainhi_ && (fitresult > 0)) {
0319           ntimes++;
0320           //      if(ntimes<=10)
0321           //        edm::LogPrint("SiPixelGainCalibrationReadDQMFile") << detid << " " << jrow << " " << icol << " " << ped << " " << ped << " " << chi2 << " " << fitresult << std::endl;
0322           VCAL_endpoint->Fill((255 - ped) / gain);
0323           peds[jrow] = ped;
0324           gains[jrow] = gain;
0325           pedforthiscol[iglobalrow] += ped;
0326           gainforthiscol[iglobalrow] += gain;
0327           nusedrows[iglobalrow]++;
0328           goodpeds->Fill(ped);
0329           goodgains->Fill(gain);
0330           detidfortree = detid;
0331           rowfortree = jrow - 1;
0332           colfortree = icol - 1;
0333           gainfortree = gain;
0334           pedfortree = ped;
0335           chi2fortree = chi2;
0336           tree->Fill();
0337         } else {
0338           nemptypixels++;
0339           //      if(nemptypixels<0.01*nrows*ncols )
0340           //        edm::LogPrint("SiPixelGainCalibrationReadDQMFile") << "ped,gain="<< ped << ","<< gain << " row, col " << jrow <<","<< icol << ", detid " << detid <<  std::endl;
0341           if (usemeanwhenempty_) {
0342             //ntimes++;
0343             //if(nemptypixels<=50)
0344             //edm::LogPrint("SiPixelGainCalibrationReadDQMFile") << "USING DEFAULT MEAN GAIN & PED (" << meangain << ","<< meanped << ")!, observed values are gain,ped : "<< gain << "," << ped <<", chi2 " << chi2  << ", fitresult "<< fitresult<<  std::endl;
0345             // edm::LogPrint("SiPixelGainCalibrationReadDQMFile") <<jrow<<"  "<<icol<<std::endl;
0346             // }
0347             peds[jrow] = meanped;
0348             gains[jrow] = meangain;
0349             std::pair<TString, int> tempval(tempgainstring, 1);
0350             badresults[detid] = tempval;
0351           } else {
0352             //if(ntimes<=100)
0353             // edm::LogPrint("SiPixelGainCalibrationReadDQMFile") << tempgainstring << " bad/dead Pixel, observed values are gain,ped : "<< gain << "," << ped <<", chi2 " << chi2 << ", fitresult "<< fitresult<< std::endl;
0354             std::pair<TString, int> tempval(tempgainstring, 2);
0355             badresults[detid] = tempval;
0356             // if everything else fails: set the gain & ped now to dead
0357             peds[jrow] = badpedval;
0358             gains[jrow] = badgainval;
0359           }
0360         }
0361         //edm::LogPrint("SiPixelGainCalibrationReadDQMFile") << detid << " " << jrow << " " << icol << " " << peds[jrow] << " " << gains[jrow] << " " << chi2 << " " << fitresult << std::endl;
0362 
0363         totgains->Fill(gains[jrow]);
0364         totpeds->Fill(peds[jrow]);
0365 
0366       }  // now collected all info, actually do the filling
0367 
0368       for (size_t jrow = 1; jrow <= nrows; jrow++) {
0369         nchannels++;
0370         size_t iglobalrow = 0;
0371         if (jrow > nrowsrocsplit)
0372           iglobalrow = 1;
0373         float ped = peds[jrow];
0374         float gain = gains[jrow];
0375 
0376         if (ped > pedlow_ && gain > gainlow_ && ped < pedhi_ && gain < gainhi_) {
0377           theGainCalibrationDbInput.setData(ped, gain, theSiPixelGainCalibrationPerPixel);
0378           theGainCalibrationDbInputOffline.setDataPedestal(ped, theSiPixelGainCalibrationGainPerColPedPerPixel);
0379         } else {
0380           theGainCalibrationDbInput.setDeadPixel(theSiPixelGainCalibrationPerPixel);
0381           theGainCalibrationDbInputOffline.setDeadPixel(theSiPixelGainCalibrationGainPerColPedPerPixel);
0382         }
0383 
0384         if (jrow % nrowsrocsplit == 0) {
0385           //      edm::LogPrint("SiPixelGainCalibrationReadDQMFile") << "now in col " << icol << " " << jrow << " " << iglobalrow << std::endl;
0386           if (nusedrows[iglobalrow] > 0) {
0387             pedforthiscol[iglobalrow] /= (float)nusedrows[iglobalrow];
0388             gainforthiscol[iglobalrow] /= (float)nusedrows[iglobalrow];
0389           }
0390           if (gainforthiscol[iglobalrow] > gainlow_ && gainforthiscol[iglobalrow] < gainhi_ &&
0391               pedforthiscol[iglobalrow] > pedlow_ && pedforthiscol[iglobalrow] < pedhi_) {  // good
0392             //      edm::LogPrint("SiPixelGainCalibrationReadDQMFile") << "setting ped & col aves: " << pedforthiscol[iglobalrow] << " " <<  gainforthiscol[iglobalrow]<< std::endl;
0393           } else {
0394             if (usemeanwhenempty_) {
0395               //          edm::LogPrint("SiPixelGainCalibrationReadDQMFile") << "setting ped & col aves: " << pedforthiscol[iglobalrow] << " " <<  gainforthiscol[iglobalrow]<< std::endl;
0396               pedforthiscol[iglobalrow] = meanped;
0397               gainforthiscol[iglobalrow] = meangain;
0398               std::pair<TString, int> tempval(tempgainstring, 3);
0399               badresults[detid] = tempval;
0400             } else {  //make dead
0401               //if(ntimes<=100)
0402               //edm::LogPrint("SiPixelGainCalibrationReadDQMFile") << tempgainstring << "dead Column, observed values are gain,ped : "<< gainforthiscol[iglobalrow] << "," << pedforthiscol[iglobalrow] << std::endl;
0403               pedforthiscol[iglobalrow] = badpedval;
0404               gainforthiscol[iglobalrow] = badgainval;
0405               std::pair<TString, int> tempval(tempgainstring, 4);
0406               badresults[detid] = tempval;
0407             }
0408           }
0409 
0410           if (gainforthiscol[iglobalrow] > gainlow_ && gainforthiscol[iglobalrow] < gainhi_ &&
0411               pedforthiscol[iglobalrow] > pedlow_ && pedforthiscol[iglobalrow] < pedhi_) {
0412             theGainCalibrationDbInputOffline.setDataGain(
0413                 gainforthiscol[iglobalrow], nrowsrocsplit, theSiPixelGainCalibrationGainPerColPedPerPixel);
0414             theGainCalibrationDbInputHLT.setData(
0415                 pedforthiscol[iglobalrow], gainforthiscol[iglobalrow], theSiPixelGainCalibrationPerColumn);
0416           } else {
0417             //      edm::LogPrint("SiPixelGainCalibrationReadDQMFile") << pedforthiscol[iglobalrow] << " " << gainforthiscol[iglobalrow] << std::endl;
0418             theGainCalibrationDbInputOffline.setDeadColumn(nrowsrocsplit,
0419                                                            theSiPixelGainCalibrationGainPerColPedPerPixel);
0420             theGainCalibrationDbInputHLT.setDeadColumn(nrowsrocsplit, theSiPixelGainCalibrationPerColumn);
0421           }
0422         }
0423       }
0424     }
0425 
0426     //    edm::LogPrint("SiPixelGainCalibrationReadDQMFile") << "setting range..." << std::endl;
0427     SiPixelGainCalibration::Range range(theSiPixelGainCalibrationPerPixel.begin(),
0428                                         theSiPixelGainCalibrationPerPixel.end());
0429     SiPixelGainCalibrationForHLT::Range hltrange(theSiPixelGainCalibrationPerColumn.begin(),
0430                                                  theSiPixelGainCalibrationPerColumn.end());
0431     SiPixelGainCalibrationOffline::Range offlinerange(theSiPixelGainCalibrationGainPerColPedPerPixel.begin(),
0432                                                       theSiPixelGainCalibrationGainPerColPedPerPixel.end());
0433 
0434     //    edm::LogPrint("SiPixelGainCalibrationReadDQMFile") <<"putting things in db..." << std::endl;
0435     // now start creating the various database objects
0436     if (!theGainCalibrationDbInput.put(detid, range, ncols))
0437       edm::LogError("SiPixelGainCalibrationAnalysis")
0438           << "warning: detid already exists for Offline (gain per col, ped per pixel) calibration database"
0439           << std::endl;
0440     if (!theGainCalibrationDbInputOffline.put(detid, offlinerange, ncols))
0441       edm::LogError("SiPixelGainCalibrationAnalysis")
0442           << "warning: detid already exists for Offline (gain per col, ped per pixel) calibration database"
0443           << std::endl;
0444     if (!theGainCalibrationDbInputHLT.put(detid, hltrange, ncols))
0445       edm::LogError("SiPixelGainCalibrationAnalysis")
0446           << "warning: detid already exists for HLT (pedestal and gain per column) calibration database" << std::endl;
0447 
0448     //delete gainPerDetid;
0449   }
0450   // now printing out summary:
0451   size_t nempty = 0;
0452   size_t ndefault = 0;
0453   size_t ndead = 0;
0454   size_t ncoldefault = 0;
0455   size_t ncoldead = 0;
0456   for (std::map<uint32_t, std::pair<TString, int> >::const_iterator ibad = badresults.begin(); ibad != badresults.end();
0457        ++ibad) {
0458     uint32_t detid = ibad->first;
0459     if (badresults[detid].second == 0) {
0460       //      edm::LogPrint("SiPixelGainCalibrationReadDQMFile") << " used pixel mean value";
0461       nempty++;
0462     } else if (badresults[detid].second == 1) {
0463       //      edm::LogPrint("SiPixelGainCalibrationReadDQMFile") << " used pixel mean value";
0464       ndefault++;
0465     } else if (badresults[detid].second == 2) {
0466       edm::LogPrint("SiPixelGainCalibrationReadDQMFile") << badresults[detid].first;
0467       edm::LogPrint("SiPixelGainCalibrationReadDQMFile") << " has one or more dead pixels";
0468       edm::LogPrint("SiPixelGainCalibrationReadDQMFile") << std::endl;
0469       ndead++;
0470     } else if (badresults[detid].second == 3) {
0471       //      edm::LogPrint("SiPixelGainCalibrationReadDQMFile") << " used column mean value";
0472       ncoldefault++;
0473     } else if (badresults[detid].second == 4) {
0474       edm::LogPrint("SiPixelGainCalibrationReadDQMFile") << badresults[detid].first;
0475       edm::LogPrint("SiPixelGainCalibrationReadDQMFile") << " has one or more dead columns";
0476       ncoldead++;
0477       edm::LogPrint("SiPixelGainCalibrationReadDQMFile") << std::endl;
0478     }
0479   }
0480   edm::LogPrint("SiPixelGainCalibrationReadDQMFile")
0481       << nempty << " modules were empty and now have pixels filled with default values." << std::endl;
0482   edm::LogPrint("SiPixelGainCalibrationReadDQMFile")
0483       << ndefault << " modules have pixels filled with default values." << std::endl;
0484   edm::LogPrint("SiPixelGainCalibrationReadDQMFile") << ndead << " modules have pixels flagged as dead." << std::endl;
0485   edm::LogPrint("SiPixelGainCalibrationReadDQMFile")
0486       << ncoldefault << " modules have columns filled with default values." << std::endl;
0487   edm::LogPrint("SiPixelGainCalibrationReadDQMFile")
0488       << ncoldead << " modules have columns filled with dead values." << std::endl;
0489   edm::LogPrint("SiPixelGainCalibrationReadDQMFile") << " ---> PIXEL Modules  " << nmodules << "\n"
0490                                                      << " ---> PIXEL Channels " << nchannels << std::endl;
0491 
0492   edm::LogPrint("SiPixelGainCalibrationReadDQMFile") << " --- writing to DB!" << std::endl;
0493   edm::Service<cond::service::PoolDBOutputService> mydbservice;
0494   if (!mydbservice.isAvailable()) {
0495     edm::LogError("db service unavailable");
0496     return;
0497   } else {
0498     if (record_ == "SiPixelGainCalibrationForHLTRcd") {
0499       edm::LogPrint("SiPixelGainCalibrationReadDQMFile")
0500           << "now doing SiPixelGainCalibrationForHLTRcd payload..." << std::endl;
0501       if (mydbservice->isNewTagRequest(record_)) {
0502         mydbservice->createOneIOV<SiPixelGainCalibrationForHLT>(
0503             theGainCalibrationDbInputHLT, mydbservice->beginOfTime(), "SiPixelGainCalibrationForHLTRcd");
0504       } else {
0505         mydbservice->appendOneIOV<SiPixelGainCalibrationForHLT>(
0506             theGainCalibrationDbInputHLT, mydbservice->currentTime(), "SiPixelGainCalibrationForHLTRcd");
0507       }
0508     } else if (record_ == "SiPixelGainCalibrationOfflineRcd") {
0509       edm::LogPrint("SiPixelGainCalibrationReadDQMFile")
0510           << "now doing SiPixelGainCalibrationOfflineRcd payload..." << std::endl;
0511       if (mydbservice->isNewTagRequest(record_)) {
0512         mydbservice->createOneIOV<SiPixelGainCalibrationOffline>(
0513             theGainCalibrationDbInputOffline, mydbservice->beginOfTime(), "SiPixelGainCalibrationOfflineRcd");
0514       } else {
0515         mydbservice->appendOneIOV<SiPixelGainCalibrationOffline>(
0516             theGainCalibrationDbInputOffline, mydbservice->currentTime(), "SiPixelGainCalibrationOfflineRcd");
0517       }
0518     }
0519     edm::LogInfo(" --- all OK");
0520   }
0521 }
0522 
0523 SiPixelGainCalibrationReadDQMFile::SiPixelGainCalibrationReadDQMFile(const edm::ParameterSet &iConfig)
0524     : appendMode_(iConfig.getUntrackedParameter<bool>("appendMode", true)),
0525       pddToken_(esConsumes()),
0526       theGainCalibrationDbInputService_(iConfig, consumesCollector()),
0527       record_(iConfig.getUntrackedParameter<std::string>("record", "SiPixelGainCalibrationOfflineRcd")),
0528       gainlow_(10.),
0529       gainhi_(0.),
0530       pedlow_(255.),
0531       pedhi_(-256),
0532       usemeanwhenempty_(iConfig.getUntrackedParameter<bool>("useMeanWhenEmpty", false)),
0533       rootfilestring_(iConfig.getUntrackedParameter<std::string>("inputrootfile", "inputfile.root")),
0534       gainmax_(20),
0535       pedmax_(200),
0536       badchi2_(iConfig.getUntrackedParameter<double>("badChi2Prob", 0.01)),
0537       nmaxcols(10 * 52),
0538       nmaxrows(160) {
0539   usesResource(TFileService::kSharedResource);
0540 
0541   //now do what ever initialization is needed
0542   ::putenv((char *)"CORAL_AUTH_USER=me");
0543   ::putenv((char *)"CORAL_AUTH_PASSWORD=test");
0544   meanGainHist_ = std::make_unique<TH1F>("meanGainHist", "mean Gain Hist", 500, 0, gainmax_);
0545   meanPedHist_ = std::make_unique<TH1F>("meanPedHist", "mean Ped Hist", 512, -200, pedmax_);
0546   defaultGain_ = std::make_unique<TH2F>("defaultGain",
0547                                         "default gain, contains mean",
0548                                         nmaxcols,
0549                                         0,
0550                                         nmaxcols,
0551                                         nmaxrows,
0552                                         0,
0553                                         nmaxrows);  // using dummy (largest) module size
0554   defaultPed_ = std::make_unique<TH2F>("defaultPed",
0555                                        "default pedestal, contains mean",
0556                                        nmaxcols,
0557                                        0,
0558                                        nmaxcols,
0559                                        nmaxrows,
0560                                        0,
0561                                        nmaxrows);  // using dummy (largest) module size
0562   defaultFitResult_ = std::make_unique<TH2F>("defaultFitResult",
0563                                              "default fitresult, contains '0'",
0564                                              nmaxcols,
0565                                              0,
0566                                              nmaxcols,
0567                                              nmaxrows,
0568                                              0,
0569                                              nmaxrows);  // using dummy (largest) module size
0570   defaultChi2_ = std::make_unique<TH2F>("defaultChi2",
0571                                         "default chi2 probability, contains '1'",
0572                                         nmaxcols,
0573                                         0,
0574                                         nmaxcols,
0575                                         nmaxrows,
0576                                         0,
0577                                         nmaxrows);  // using dummy (largest) module size
0578 }
0579 
0580 //
0581 // member functions
0582 //
0583 
0584 // ------------ method called to for each event  ------------
0585 void SiPixelGainCalibrationReadDQMFile::analyze(const edm::Event &iEvent, const edm::EventSetup &iSetup) {
0586   auto histfile = getHistograms();
0587   fillDatabase(iSetup, histfile.get());
0588   // empty but should not be called anyway
0589 }
0590 
0591 std::unique_ptr<TFile> SiPixelGainCalibrationReadDQMFile::getHistograms() {
0592   edm::LogPrint("SiPixelGainCalibrationReadDQMFile") << "now parsing file " << rootfilestring_ << std::endl;
0593   auto therootfile = std::make_unique<TFile>(rootfilestring_.c_str());
0594   therootfile->cd();
0595   TDirectory *dir = therootfile->GetDirectory("DQMData");
0596   TList *list = dir->GetListOfKeys();
0597 
0598   TString comparestring = "Module";
0599 
0600   std::vector<TString> keylist;
0601   std::vector<TString> hist2list;
0602   std::vector<TString> dirlist;
0603   std::vector<TString> notdonelist;
0604   std::vector<int> nsubdirs;
0605   int ikey = 0;
0606 
0607   for (ikey = 0; ikey < list->GetEntries(); ikey++) {
0608     TKey *thekey = (TKey *)list->At(ikey);
0609     if (thekey == nullptr)
0610       continue;
0611     TString keyname = thekey->GetName();
0612     TString keytype = thekey->GetClassName();
0613     //    if(keyname=="EventInfo")
0614     //      continue;
0615     //    edm::LogPrint("SiPixelGainCalibrationReadDQMFile") <<  keytype << " " << keyname << std::endl;
0616     if (keytype == "TDirectoryFile") {
0617       TString dirname = dir->GetPath();
0618       dirname += "/";
0619       dirname += keyname;
0620       //      edm::LogPrint("SiPixelGainCalibrationReadDQMFile") << dirname << std::endl;
0621       dir = therootfile->GetDirectory(dirname);
0622 
0623       list = dir->GetListOfKeys();
0624       if (dirname.Contains(comparestring)) {
0625         dirlist.push_back(dirname);
0626         //  edm::LogPrint("SiPixelGainCalibrationReadDQMFile") << dirname << std::endl;
0627       } else {
0628         notdonelist.push_back(dirname);
0629         nsubdirs.push_back(-1);
0630       }
0631     }
0632   }
0633   size_t nempty = 0;
0634   while (nempty != notdonelist.size()) {
0635     for (size_t idir = 0; idir < notdonelist.size(); ++idir) {
0636       if (nsubdirs[idir] == 0)
0637         continue;
0638       //      edm::LogPrint("SiPixelGainCalibrationReadDQMFile") << "now examining " << notdonelist[idir]<< " " << nsubdirs[idir] <<  std::endl;
0639       dir = therootfile->GetDirectory(notdonelist[idir]);
0640       //      edm::LogPrint("SiPixelGainCalibrationReadDQMFile") << dir->GetName() << std::endl;
0641       list = dir->GetListOfKeys();
0642       //      edm::LogPrint("SiPixelGainCalibrationReadDQMFile") << list->GetEntries() << std::endl;
0643       int ndirectories = 0;
0644       for (ikey = 0; ikey < list->GetEntries(); ikey++) {
0645         TKey *thekey = (TKey *)list->At(ikey);
0646         if (thekey == nullptr)
0647           continue;
0648         TString keyname = thekey->GetName();
0649         TString keytype = thekey->GetClassName();
0650         //  edm::LogPrint("SiPixelGainCalibrationReadDQMFile") << keyname << " " << keytype << std::endl;
0651         if (keytype == "TDirectoryFile") {
0652           TString dirname = dir->GetPath();
0653           dirname += "/";
0654           dirname += keyname;
0655           //      edm::LogPrint("SiPixelGainCalibrationReadDQMFile") << dirname << std::endl;
0656           ndirectories++;
0657           if (dirname.Contains(comparestring)) {
0658             //      edm::LogPrint("SiPixelGainCalibrationReadDQMFile") << dirname << std::endl;
0659             dirlist.push_back(dirname);
0660           } else {
0661             notdonelist.push_back(dirname);
0662             nsubdirs.push_back(-1);
0663           }
0664         }
0665       }
0666       nsubdirs[idir] = ndirectories;
0667       // edm::LogPrint("SiPixelGainCalibrationReadDQMFile") << "now done examining " << notdonelist[idir]<< " " << nsubdirs[idir] <<  std::endl;
0668     }
0669     nempty = 0;
0670     for (size_t i = 0; i < nsubdirs.size(); i++) {
0671       if (nsubdirs[i] != -1)
0672         nempty++;
0673     }
0674   }
0675   //  edm::LogPrint("SiPixelGainCalibrationReadDQMFile") << "\n done!" << std::endl;
0676 
0677   for (size_t idir = 0; idir < dirlist.size(); ++idir) {
0678     //    edm::LogPrint("SiPixelGainCalibrationReadDQMFile") << "good dir "  << dirlist[idir] << std::endl;
0679 
0680     uint32_t detid = 1;
0681 
0682     dir = therootfile->GetDirectory(dirlist[idir]);
0683     list = dir->GetListOfKeys();
0684     for (ikey = 0; ikey < list->GetEntries(); ikey++) {
0685       TKey *thekey = (TKey *)list->At(ikey);
0686       if (thekey == nullptr)
0687         continue;
0688       TString keyname = thekey->GetName();
0689       TString keytype = thekey->GetClassName();
0690       //    edm::LogPrint("SiPixelGainCalibrationReadDQMFile") << keyname << " " << keytype << std::endl;
0691       if (keytype == "TH2F" && (keyname.Contains("Gain2d") || keyname.Contains("Pedestal2d") ||
0692                                 keyname.Contains("GainChi2Prob2d") || keyname.Contains("GainFitResult2d"))) {
0693         TString detidstring = keyname;
0694         detidstring.Remove(0, detidstring.Sizeof() - 10);
0695 
0696         detid = atoi(detidstring.Data());
0697 
0698         if (keyname.Contains("GainChi2Prob2d")) {
0699           TString tempstr = dirlist[idir];
0700           tempstr += "/";
0701           tempstr += keyname;
0702           TString replacestring = rootfilestring_;
0703           replacestring += ":";
0704           tempstr.ReplaceAll(replacestring, "");
0705           //      edm::LogPrint("SiPixelGainCalibrationReadDQMFile") << tempstr << std::endl;
0706           bookkeeper_[detid]["chi2prob_2d"] = tempstr;
0707         } else if (keyname.Contains("GainFitResult2d")) {
0708           TString tempstr = dirlist[idir];
0709           tempstr += "/";
0710           tempstr += keyname;
0711           TString replacestring = rootfilestring_;
0712           replacestring += ":";
0713           tempstr.ReplaceAll(replacestring, "");
0714           //      edm::LogPrint("SiPixelGainCalibrationReadDQMFile") << tempstr << std::endl;
0715           bookkeeper_[detid]["fitresult_2d"] = tempstr;
0716         } else if (keyname.Contains("Gain2d")) {
0717           //      edm::LogPrint("SiPixelGainCalibrationReadDQMFile") << dirlist[idir] << std::endl;
0718           std::map<std::string, TString> tempmap;
0719           TString tempstr = dirlist[idir];
0720           tempstr += "/";
0721           tempstr += keyname;
0722           TString replacestring = rootfilestring_;
0723           replacestring += ":";
0724           tempstr.ReplaceAll(replacestring, "");
0725           //      edm::LogPrint("SiPixelGainCalibrationReadDQMFile") << tempstr << std::endl;
0726           bookkeeper_[detid]["gain_2d"] = tempstr;
0727           //edm::LogPrint("SiPixelGainCalibrationReadDQMFile") << detidstring << " " << keyname << " " << detid << " " << bookkeeper_[detid]["gain_2d"] << std::endl ;
0728         }
0729         if (keyname.Contains("Pedestal2d")) {
0730           //      edm::LogPrint("SiPixelGainCalibrationReadDQMFile") << dirlist[idir] << std::endl;
0731           std::map<std::string, TString> tempmap;
0732 
0733           TString tempstr = dirlist[idir];
0734           tempstr += "/";
0735           tempstr += keyname;
0736           //      edm::LogPrint("SiPixelGainCalibrationReadDQMFile") << tempstr << std::endl;
0737           TString replacestring = rootfilestring_;
0738           replacestring += ":";
0739           tempstr.ReplaceAll(replacestring, "");
0740           bookkeeper_[detid]["ped_2d"] = tempstr;
0741 
0742           //edm::LogPrint("SiPixelGainCalibrationReadDQMFile") << detidstring << " " << keyname << " " << detid << " " << bookkeeper_[detid]["ped_2d"] << std::endl ;
0743         }
0744         //      edm::LogPrint("SiPixelGainCalibrationReadDQMFile") << keyname << " " << keytype << std::endl;
0745       }
0746     }
0747 
0748     TH2F *temphistoped = dynamic_cast<TH2F *>(therootfile->Get(bookkeeper_[detid]["ped_2d"]));
0749     TH2F *temphistogain = dynamic_cast<TH2F *>(therootfile->Get(bookkeeper_[detid]["gain_2d"]));
0750     TH2F *temphistofitresult = dynamic_cast<TH2F *>(therootfile->Get(bookkeeper_[detid]["gain_2d"]));
0751 
0752     for (int xbin = 1; xbin <= temphistoped->GetNbinsX(); ++xbin) {
0753       for (int ybin = 1; ybin <= temphistoped->GetNbinsY(); ++ybin) {
0754         if (temphistofitresult->GetBinContent(xbin, ybin) <= 0)
0755           continue;
0756         float val = temphistoped->GetBinContent(xbin, ybin);
0757         if (val > pedmax_)
0758           continue;
0759         if (pedlow_ > val)
0760           pedlow_ = val;
0761         if (pedhi_ < val)
0762           pedhi_ = val;
0763         meanPedHist_->Fill(val);
0764       }
0765     }
0766 
0767     for (int xbin = 1; xbin <= temphistogain->GetNbinsX(); ++xbin) {
0768       for (int ybin = 1; ybin <= temphistogain->GetNbinsY(); ++ybin) {
0769         if (temphistofitresult->GetBinContent(xbin, ybin) <= 0)
0770           continue;
0771         float val = temphistogain->GetBinContent(xbin, ybin);
0772         if (val <= 0.0001)
0773           continue;
0774         if (gainlow_ > val)
0775           gainlow_ = val;
0776         if (gainhi_ < val)
0777           gainhi_ = val;
0778         meanGainHist_->Fill(val);
0779       }
0780     }
0781 
0782   }  // end of loop over dirlist
0783   return therootfile;
0784 }
0785 
0786 //define this as a plug-in
0787 DEFINE_FWK_MODULE(SiPixelGainCalibrationReadDQMFile);