Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-11-27 03:17:28

0001 // -*- C++ -*-
0002 //
0003 // Package:    Test/HcalLutAnalyzer
0004 // Class:      HcalLutAnalyzer
0005 //
0006 /**\class HcalLutAnalyzer HcalLutAnalyzer.cc Test/HcalLutAnalyzer/plugins/HcalLutAnalyzer.cc
0007 
0008  Description: [one line class summary]
0009 
0010  Implementation:
0011      [Notes on implementation]
0012 */
0013 //
0014 // Original Author:  Aleko Khukhunaishvili
0015 //         Created:  Fri, 21 Jul 2017 08:42:05 GMT
0016 //
0017 //
0018 
0019 // system include files
0020 #include <memory>
0021 #include <iostream>
0022 #include <fstream>
0023 
0024 // user include files
0025 #include "FWCore/Framework/interface/Frameworkfwd.h"
0026 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0027 #include "FWCore/Framework/interface/Event.h"
0028 #include "FWCore/Framework/interface/MakerMacros.h"
0029 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0030 #include "FWCore/Framework/interface/EventSetup.h"
0031 
0032 #include "CalibCalorimetry/HcalTPGAlgos/interface/XMLProcessor.h"
0033 #include "CalibCalorimetry/HcalTPGAlgos/interface/LutXml.h"
0034 #include "DataFormats/HcalDetId/interface/HcalGenericDetId.h"
0035 #include "DataFormats/HcalDetId/interface/HcalDetId.h"
0036 #include "Geometry/CaloTopology/interface/HcalTopology.h"
0037 #include "Geometry/Records/interface/IdealGeometryRecord.h"
0038 #include "Geometry/Records/interface/HcalRecNumberingRecord.h"
0039 #include "CalibFormats/HcalObjects/interface/HcalDbRecord.h"
0040 
0041 #include "TString.h"
0042 #include "TH1D.h"
0043 #include "TH2D.h"
0044 #include "TProfile.h"
0045 #include "TCanvas.h"
0046 #include "TROOT.h"
0047 #include "TStyle.h"
0048 #include "TSystem.h"
0049 
0050 class HcalLutAnalyzer : public edm::one::EDAnalyzer<edm::one::SharedResources> {
0051 public:
0052   explicit HcalLutAnalyzer(const edm::ParameterSet&);
0053   ~HcalLutAnalyzer() override {}
0054   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0055 
0056 private:
0057   void analyze(const edm::Event&, const edm::EventSetup&) override;
0058 
0059   std::string inputDir;
0060   std::string plotsDir;
0061   std::vector<std::string> tags_;
0062   std::vector<std::string> quality_;
0063   std::vector<std::string> pedestals_;
0064   std::vector<std::string> effpedestals_;
0065   std::vector<std::string> gains_;
0066   std::vector<std::string> respcorrs_;
0067 
0068   double Zmin;
0069   double Zmax;
0070   double Ymin;
0071   double Ymax;
0072   double Pmin;
0073   double Pmax;
0074 
0075   edm::ESGetToken<HcalTopology, HcalRecNumberingRecord> tok_htopo_;
0076   edm::ESGetToken<HcalElectronicsMap, HcalDbRecord> tok_emap_;
0077 };
0078 
0079 HcalLutAnalyzer::HcalLutAnalyzer(const edm::ParameterSet& iConfig) {
0080   inputDir = iConfig.getParameter<std::string>("inputDir");
0081   plotsDir = iConfig.getParameter<std::string>("plotsDir");
0082   tags_ = iConfig.getParameter<std::vector<std::string> >("tags");
0083   quality_ = iConfig.getParameter<std::vector<std::string> >("quality");
0084   pedestals_ = iConfig.getParameter<std::vector<std::string> >("pedestals");
0085   effpedestals_ = iConfig.getParameter<std::vector<std::string> >("effpedestals");
0086   gains_ = iConfig.getParameter<std::vector<std::string> >("gains");
0087   respcorrs_ = iConfig.getParameter<std::vector<std::string> >("respcorrs");
0088 
0089   Zmin = iConfig.getParameter<double>("Zmin");
0090   Zmax = iConfig.getParameter<double>("Zmax");
0091   Ymin = iConfig.getParameter<double>("Ymin");
0092   Ymax = iConfig.getParameter<double>("Ymax");
0093   Pmin = iConfig.getParameter<double>("Pmin");
0094   Pmax = iConfig.getParameter<double>("Pmax");
0095 
0096   tok_htopo_ = esConsumes<HcalTopology, HcalRecNumberingRecord>();
0097   tok_emap_ = esConsumes<HcalElectronicsMap, HcalDbRecord>();
0098 }
0099 
0100 void HcalLutAnalyzer::analyze(const edm::Event&, const edm::EventSetup& iSetup) {
0101   using namespace std;
0102 
0103   const HcalTopology* topology = &iSetup.getData(tok_htopo_);
0104   const HcalElectronicsMap* electronicsMap = &iSetup.getData(tok_emap_);
0105 
0106   typedef std::vector<std::string> vstring;
0107   typedef std::map<unsigned long int, float> LUTINPUT;
0108 
0109   static const int NVAR = 6;  //variables
0110   static const int NDET = 2;  //detectors
0111   static const int NDEP = 7;  //depths
0112   static const int NLEV = 3;  //old,new,ratio
0113 
0114   const bool doRatio[NVAR] = {false, true, true, false, true, false};
0115   const char* titleVar[NVAR] = {"Pedestals", "RespCorrs", "Gains", "Threshold", "LUTs", "EffPedestals"};
0116   const char* titleHisR[NLEV] = {"Old", "New", "Ratio"};
0117   const char* titleHisD[NLEV] = {"Old", "New", "Difference"};
0118   const char* titleDet[4] = {"HBHE", "HF", "HEP17", "HO"};
0119   const int DEP[NDET] = {7, 4};
0120   const char* titleDep[NDEP] = {"depth1", "depth2", "depth3", "depth4", "depth5", "depth6", "depth7"};
0121 
0122   TH2D* r[NVAR][NDET];
0123   TProfile* p[NVAR][NDET];
0124 
0125   TH2D* h[NVAR][NLEV][NDET][NDEP];
0126   TH2D* hlut[4][NLEV];
0127   TH2D* hslope[2];
0128   TH2D* houtput[4][2];
0129 
0130   for (int d = 0; d < 4; ++d) {
0131     for (int i = 0; i < 3; ++i) {
0132       hlut[d][i] = new TH2D(Form("Lut_%s_%s", titleDet[d], titleHisR[i]),
0133                             Form("Input LUT, %s", (i == 2 ? "Ratio" : tags_[i].c_str())),
0134                             260,
0135                             0,
0136                             260,
0137                             240,
0138                             0,
0139                             i == NLEV - 1 ? 3 : 2400);
0140       hlut[d][i]->SetMarkerColor(d == 0 ? kBlue : d == 1 ? kGreen + 2 : d == 2 ? kRed : kCyan);
0141       hlut[d][i]->SetXTitle("raw adc");
0142       hlut[d][i]->SetYTitle("lin adc");
0143     }
0144   }
0145 
0146   for (int d = 0; d < NDET; ++d) {
0147     hslope[d] = new TH2D(
0148         Form("GainLutScatter_%s", titleDet[d]), Form("Gain-Lutslope scatter, %s", titleDet[d]), 200, 0, 2, 200, 0, 2);
0149     hslope[d]->SetXTitle("Gain x RespCorr ratio");
0150     hslope[d]->SetYTitle("Lut ratio");
0151 
0152     for (int j = 0; j < NVAR; ++j) {
0153       double rmin = doRatio[j] ? Ymin : -4;
0154       double rmax = doRatio[j] ? Ymax : 4;
0155       r[j][d] = new TH2D(Form("r%s_%s", titleVar[j], titleDet[d]),
0156                          Form("%s, %s", titleVar[j], titleDet[d]),
0157                          83,
0158                          -41.5,
0159                          41.5,
0160                          250,
0161                          rmin,
0162                          rmax);
0163       r[j][d]->SetXTitle("iEta");
0164       r[j][d]->SetYTitle(doRatio[j] ? "New / Old" : "New - Old");
0165       p[j][d] = new TProfile(
0166           Form("p%s_%s", titleVar[j], titleDet[d]), Form("%s, %s", titleVar[j], titleDet[d]), 83, -41.5, 41.5);
0167       p[j][d]->SetXTitle("iEta");
0168       p[j][d]->SetYTitle(doRatio[j] ? "New / Old" : "New - Old");
0169       p[j][d]->SetMarkerStyle(20);
0170       p[j][d]->SetMarkerSize(0.9);
0171       p[j][d]->SetMarkerColor(kBlue);
0172 
0173       for (int p = 0; p < DEP[d]; ++p) {
0174         for (int i = 0; i < NLEV; ++i) {
0175           const char* titHist = doRatio[j] ? titleHisR[i] : titleHisD[i];
0176           h[j][i][d][p] = new TH2D(Form("h%s_%s_%s_%s", titleVar[j], titHist, titleDet[d], titleDep[p]),
0177                                    Form("%s, %s, %s, %s", titleVar[j], titHist, titleDet[d], titleDep[p]),
0178                                    83,
0179                                    -41.5,
0180                                    41.5,
0181                                    72,
0182                                    0.5,
0183                                    72.5);
0184           h[j][i][d][p]->SetXTitle("iEta");
0185           h[j][i][d][p]->SetYTitle("iPhi");
0186         }
0187       }
0188     }
0189   }
0190 
0191   for (int i = 0; i < 4; ++i) {
0192     int color = i == 0 ? kBlue : i == 1 ? kViolet : i == 2 ? kGreen + 2 : kRed;
0193     houtput[i][0] =
0194         new TH2D(Form("houtlut0_%d", i), Form("Output LUT, %s", tags_[0].c_str()), 2100, 0, 2100, 260, 0, 260);
0195     houtput[i][1] =
0196         new TH2D(Form("houtlut1_%d", i), Form("Output LUT, %s", tags_[1].c_str()), 2100, 0, 2100, 260, 0, 260);
0197     for (int j = 0; j < 2; ++j) {
0198       houtput[i][j]->SetMarkerColor(color);
0199       houtput[i][j]->SetLineColor(color);
0200     }
0201   }
0202 
0203   //FILL LUT INPUT DATA
0204   LUTINPUT lutgain[2];
0205   LUTINPUT lutresp[2];
0206   LUTINPUT lutpede[2];
0207   LUTINPUT luteffpede[2];
0208 
0209   assert(tags_.size() == 2);
0210   assert(gains_.size() == 2);
0211   assert(respcorrs_.size() == 2);
0212   assert(pedestals_.size() == 2);
0213   assert(effpedestals_.size() == 2);
0214 
0215   unsigned long int iraw;
0216   int ieta, iphi, idep;
0217   string det, base;
0218   float val1, val2, val3, val4;
0219   float wid1, wid2, wid3, wid4;
0220   char buffer[1024];
0221 
0222   std::vector<HcalDetId> BadChans[2];
0223   std::vector<HcalDetId> ZeroLuts[2];
0224 
0225   //Read input condition files
0226   for (int ii = 0; ii < 2; ++ii) {
0227     //Gains
0228     std::ifstream infile(
0229         edm::FileInPath(Form("%s/Gains/Gains_Run%s.txt", inputDir.c_str(), gains_[ii].c_str())).fullPath().c_str());
0230     assert(!infile.fail());
0231     while (!infile.eof()) {
0232       infile.getline(buffer, 1024);
0233       if (buffer[0] == '#')
0234         continue;
0235       std::istringstream(buffer) >> ieta >> iphi >> idep >> det >> val1 >> val2 >> val3 >> val4 >> iraw;
0236       if (det != "HB" && det != "HE" && det != "HF")
0237         continue;
0238 
0239       float theval = (val1 + val2 + val3 + val4) / 4.0;
0240 
0241       HcalSubdetector subdet = det == "HB"   ? HcalBarrel
0242                                : det == "HE" ? HcalEndcap
0243                                : det == "HF" ? HcalForward
0244                                              : HcalOther;
0245 
0246       HcalDetId id(subdet, ieta, iphi, idep);
0247       lutgain[ii].insert(LUTINPUT::value_type(id.rawId(), theval));
0248     }
0249 
0250     //Effective Pedestals
0251     std::ifstream infeffped(
0252         edm::FileInPath(
0253             Form("%s/EffectivePedestals/EffectivePedestals_Run%s.txt", inputDir.c_str(), effpedestals_[ii].c_str()))
0254             .fullPath()
0255             .c_str());
0256     assert(!infeffped.fail());
0257     while (!infeffped.eof()) {
0258       infeffped.getline(buffer, 1024);
0259       if (buffer[0] == '#')
0260         continue;
0261       std::istringstream(buffer) >> ieta >> iphi >> idep >> det >> val1 >> val2 >> val3 >> val4 >> wid1 >> wid2 >>
0262           wid3 >> wid4 >> iraw;
0263       if (det != "HB" && det != "HE" && det != "HF")
0264         continue;
0265 
0266       float theval = (val1 + val2 + val3 + val4) / 4.0;
0267 
0268       HcalSubdetector subdet = det == "HB"   ? HcalBarrel
0269                                : det == "HE" ? HcalEndcap
0270                                : det == "HF" ? HcalForward
0271                                              : HcalOther;
0272 
0273       HcalDetId id(subdet, ieta, iphi, idep);
0274       luteffpede[ii].insert(LUTINPUT::value_type(id.rawId(), theval));
0275     }
0276 
0277     //Pedestals
0278     std::ifstream infped(
0279         edm::FileInPath(Form("%s/Pedestals/Pedestals_Run%s.txt", inputDir.c_str(), pedestals_[ii].c_str()))
0280             .fullPath()
0281             .c_str());
0282     assert(!infped.fail());
0283     while (!infped.eof()) {
0284       infped.getline(buffer, 1024);
0285       if (buffer[0] == '#')
0286         continue;
0287       std::istringstream(buffer) >> ieta >> iphi >> idep >> det >> val1 >> val2 >> val3 >> val4 >> wid1 >> wid2 >>
0288           wid3 >> wid4 >> iraw;
0289       if (det != "HB" && det != "HE" && det != "HF")
0290         continue;
0291 
0292       float theval = (val1 + val2 + val3 + val4) / 4.0;
0293 
0294       HcalSubdetector subdet = det == "HB"   ? HcalBarrel
0295                                : det == "HE" ? HcalEndcap
0296                                : det == "HF" ? HcalForward
0297                                              : HcalOther;
0298 
0299       HcalDetId id(subdet, ieta, iphi, idep);
0300       lutpede[ii].insert(LUTINPUT::value_type(id.rawId(), theval));
0301     }
0302 
0303     //Response corrections
0304     std::ifstream inresp(
0305         edm::FileInPath(Form("%s/RespCorrs/RespCorrs_Run%s.txt", inputDir.c_str(), respcorrs_[ii].c_str()))
0306             .fullPath()
0307             .c_str());
0308     assert(!inresp.fail());
0309     while (!inresp.eof()) {
0310       inresp.getline(buffer, 1024);
0311       if (buffer[0] == '#')
0312         continue;
0313       std::istringstream(buffer) >> ieta >> iphi >> idep >> det >> val1 >> iraw;
0314       if (det != "HB" && det != "HE" && det != "HF")
0315         continue;
0316 
0317       float theval = val1;
0318 
0319       HcalSubdetector subdet = det == "HB"   ? HcalBarrel
0320                                : det == "HE" ? HcalEndcap
0321                                : det == "HF" ? HcalForward
0322                                              : HcalOther;
0323 
0324       HcalDetId id(subdet, ieta, iphi, idep);
0325       lutresp[ii].insert(LUTINPUT::value_type(id.rawId(), theval));
0326     }
0327 
0328     //ChannelQuality
0329     std::ifstream inchan(
0330         edm::FileInPath(Form("%s/ChannelQuality/ChannelQuality_Run%s.txt", inputDir.c_str(), quality_[ii].c_str()))
0331             .fullPath()
0332             .c_str());
0333     assert(!inchan.fail());
0334     while (!inchan.eof()) {
0335       inchan.getline(buffer, 1024);
0336       if (buffer[0] == '#')
0337         continue;
0338       std::istringstream(buffer) >> ieta >> iphi >> idep >> det >> base >> val1 >> iraw;
0339 
0340       float theval = val1;
0341 
0342       HcalSubdetector subdet = det == "HB"   ? HcalBarrel
0343                                : det == "HE" ? HcalEndcap
0344                                : det == "HF" ? HcalForward
0345                                : det == "HO" ? HcalOuter
0346                                              : HcalOther;
0347 
0348       HcalDetId id(subdet, ieta, iphi, idep);
0349       if (theval != 0)
0350         BadChans[ii].push_back(id);
0351     }
0352   }
0353 
0354   LutXml xmls1(edm::FileInPath(Form("%s/%s/%s.xml", inputDir.c_str(), tags_[0].c_str(), tags_[0].c_str())).fullPath());
0355   LutXml xmls2(edm::FileInPath(Form("%s/%s/%s.xml", inputDir.c_str(), tags_[1].c_str(), tags_[1].c_str())).fullPath());
0356 
0357   xmls1.create_lut_map(electronicsMap);
0358   xmls2.create_lut_map(electronicsMap);
0359 
0360   for (const auto& xml2 : xmls2) {
0361     HcalGenericDetId detid(xml2.first);
0362 
0363     if (detid.genericSubdet() == HcalGenericDetId::HcalGenTriggerTower) {
0364       HcalTrigTowerDetId tid(detid.rawId());
0365       if (!topology->validHT(tid))
0366         continue;
0367       const auto& lut2 = xml2.second;
0368 
0369       int D = abs(tid.ieta()) < 29 ? (lut2.size() == 1024 ? 0 : 3) : tid.version() == 0 ? 1 : 2;
0370       for (size_t i = 0; i < lut2.size(); ++i) {
0371         if (int(i) % 4 == D)
0372           houtput[D][1]->Fill(i, lut2[i]);
0373       }
0374     } else if (topology->valid(detid)) {
0375       HcalDetId id(detid);
0376       HcalSubdetector subdet = id.subdet();
0377       int idet = int(subdet);
0378       const auto& lut2 = xml2.second;
0379       int hbhe = idet == HcalForward ? 1 : idet == HcalOuter ? 3 : lut2.size() == 128 ? 0 : 2;
0380       for (size_t i = 0; i < lut2.size(); ++i) {
0381         hlut[hbhe][1]->Fill(i, lut2[i]);
0382         if (hbhe == 2)
0383           hlut[hbhe][1]->Fill(i, lut2[i] & 0x3FF);
0384       }
0385     }
0386   }
0387 
0388   for (const auto& xml1 : xmls1) {
0389     HcalGenericDetId detid(xml1.first);
0390     const auto& lut1 = xml1.second;
0391 
0392     if (detid.genericSubdet() == HcalGenericDetId::HcalGenTriggerTower) {
0393       HcalTrigTowerDetId tid(detid.rawId());
0394       if (!topology->validHT(tid))
0395         continue;
0396       int D = abs(tid.ieta()) < 29 ? (lut1.size() == 1024 ? 0 : 3) : tid.version() == 0 ? 1 : 2;
0397       for (size_t i = 0; i < lut1.size(); ++i) {
0398         if (int(i) % 4 == D)
0399           houtput[D][0]->Fill(i, lut1[i]);
0400       }
0401     } else if (topology->valid(detid)) {
0402       HcalDetId id(detid);
0403       HcalSubdetector subdet = id.subdet();
0404       int idet = int(subdet);
0405       const auto& lut1 = xml1.second;
0406       int hbhe = idet == HcalForward ? 1 : idet == HcalOuter ? 3 : lut1.size() == 128 ? 0 : 2;
0407       for (size_t i = 0; i < lut1.size(); ++i) {
0408         hlut[hbhe][0]->Fill(i, lut1[i]);
0409         if (hbhe == 2)
0410           hlut[hbhe][0]->Fill(i, lut1[i] & 0x3FF);
0411       }
0412     }
0413 
0414     auto xml2 = xmls2.find(detid.rawId());
0415     if (xml2 == xmls2.end())
0416       continue;
0417 
0418     if (detid.genericSubdet() == HcalGenericDetId::HcalGenTriggerTower)
0419       continue;
0420 
0421     HcalDetId id(detid);
0422 
0423     HcalSubdetector subdet = id.subdet();
0424     int idet = int(subdet);
0425     int ieta = id.ieta();
0426     int iphi = id.iphi();
0427     int idep = id.depth() - 1;
0428     unsigned long int iraw = id.rawId();
0429 
0430     if (!topology->valid(detid))
0431       continue;
0432 
0433     int hbhe = idet == HcalForward;
0434 
0435     const auto& lut2 = xml2->second;
0436 
0437     size_t size = lut1.size();
0438     if (size != lut2.size())
0439       continue;
0440 
0441     std::vector<unsigned int> llut1(size);
0442     std::vector<unsigned int> llut2(size);
0443     for (size_t i = 0; i < size; ++i) {
0444       llut1[i] = hbhe == 0 ? lut1[i] & 0x3FF : lut1[i];
0445       llut2[i] = hbhe == 0 ? lut2[i] & 0x3FF : lut2[i];
0446     }
0447 
0448     int threshold[2] = {0, 0};
0449     //Thresholds
0450     for (size_t i = 0; i < size; ++i) {
0451       if (llut1[i] > 0) {
0452         threshold[0] = i;
0453         break;
0454       }
0455       if (i == size - 1) {
0456         ZeroLuts[0].push_back(id);
0457       }
0458     }
0459     for (size_t i = 0; i < size; ++i) {
0460       if (llut2[i] > 0) {
0461         threshold[1] = i;
0462         break;
0463       }
0464       if (i == size - 1) {
0465         ZeroLuts[1].push_back(id);
0466       }
0467     }
0468 
0469     if (subdet != HcalBarrel && subdet != HcalEndcap && subdet != HcalForward)
0470       continue;
0471 
0472     //fill conditions
0473 
0474     double xfill = 0;
0475     h[0][0][hbhe][idep]->Fill(ieta, iphi, lutpede[0][iraw]);
0476     h[0][1][hbhe][idep]->Fill(ieta, iphi, lutpede[1][iraw]);
0477     xfill = lutpede[1][iraw] - lutpede[0][iraw];
0478     h[0][2][hbhe][idep]->Fill(ieta, iphi, xfill);
0479     r[0][hbhe]->Fill(ieta, xfill);
0480     p[0][hbhe]->Fill(ieta, xfill);
0481 
0482     h[1][0][hbhe][idep]->Fill(ieta, iphi, lutresp[0][iraw]);
0483     h[1][1][hbhe][idep]->Fill(ieta, iphi, lutresp[1][iraw]);
0484     xfill = lutresp[1][iraw] / lutresp[0][iraw];
0485     h[1][2][hbhe][idep]->Fill(ieta, iphi, xfill);
0486     r[1][hbhe]->Fill(ieta, xfill);
0487     p[1][hbhe]->Fill(ieta, xfill);
0488 
0489     h[2][0][hbhe][idep]->Fill(ieta, iphi, lutgain[0][iraw]);
0490     h[2][1][hbhe][idep]->Fill(ieta, iphi, lutgain[1][iraw]);
0491     xfill = lutgain[1][iraw] / lutgain[0][iraw];
0492     h[2][2][hbhe][idep]->Fill(ieta, iphi, xfill);
0493     r[2][hbhe]->Fill(ieta, xfill);
0494     p[2][hbhe]->Fill(ieta, xfill);
0495 
0496     h[3][0][hbhe][idep]->Fill(ieta, iphi, threshold[0]);
0497     h[3][1][hbhe][idep]->Fill(ieta, iphi, threshold[1]);
0498     xfill = threshold[1] - threshold[0];
0499     h[3][2][hbhe][idep]->Fill(ieta, iphi, xfill);
0500     r[3][hbhe]->Fill(ieta, xfill);
0501     p[3][hbhe]->Fill(ieta, xfill);
0502 
0503     h[5][0][hbhe][idep]->Fill(ieta, iphi, luteffpede[0][iraw]);
0504     h[5][1][hbhe][idep]->Fill(ieta, iphi, luteffpede[1][iraw]);
0505     xfill = luteffpede[1][iraw] - luteffpede[0][iraw];
0506     h[5][2][hbhe][idep]->Fill(ieta, iphi, xfill);
0507     r[5][hbhe]->Fill(ieta, xfill);
0508     p[5][hbhe]->Fill(ieta, xfill);
0509 
0510     size_t maxvalue = hbhe == 0 ? 1023 : 2047;
0511 
0512     //LutDifference
0513     for (size_t i = 0; i < size; ++i) {
0514       hlut[hbhe][2]->Fill(i, llut1[i] == 0 ? 0 : (double)llut2[i] / llut1[i]);
0515 
0516       if (i == size - 1 ||
0517           (llut1[i] == maxvalue || llut2[i] == maxvalue)) {  //Fill with only the last one before the maximum
0518         if (llut1[i - 1] == 0 || llut2[i - 1] == 0) {
0519           break;
0520         }
0521         double condratio = lutgain[1][iraw] / lutgain[0][iraw] * lutresp[1][iraw] / lutresp[0][iraw];
0522         xfill = (double)llut2[i - 1] / llut1[i - 1];
0523         hslope[hbhe]->Fill(condratio, xfill);
0524 
0525         h[4][0][hbhe][idep]->Fill(ieta, iphi, (double)llut1[i - 1] / (i - 1));
0526         h[4][1][hbhe][idep]->Fill(ieta, iphi, (double)llut2[i - 1] / (i - 1));
0527         h[4][2][hbhe][idep]->Fill(ieta, iphi, xfill);
0528         r[4][hbhe]->Fill(ieta, xfill);
0529         p[4][hbhe]->Fill(ieta, xfill);
0530 
0531         break;
0532       }
0533     }
0534   }
0535 
0536   gROOT->SetStyle("Plain");
0537   gStyle->SetPalette(1);
0538   gStyle->SetStatW(0.2);
0539   gStyle->SetStatH(0.1);
0540   gStyle->SetStatY(1.0);
0541   gStyle->SetStatX(0.9);
0542   gStyle->SetOptStat(110010);
0543   gStyle->SetOptFit(111111);
0544 
0545   //Draw and Print
0546   TCanvas* cc = new TCanvas("cc", "cc", 0, 0, 1600, 1200);
0547   cc->SetGridy();
0548   for (int j = 0; j < NVAR; ++j) {
0549     gSystem->mkdir(TString(plotsDir) + "/_" + titleVar[j]);
0550     for (int d = 0; d < NDET; ++d) {
0551       cc->Clear();
0552       r[j][d]->Draw("colz");
0553       cc->Print(TString(plotsDir) + "/_" + titleVar[j] + "/" + TString(r[j][d]->GetName()) + ".pdf");
0554 
0555       cc->Clear();
0556       p[j][d]->Draw();
0557       if (doRatio[j]) {
0558         p[j][d]->SetMinimum(Pmin);
0559         p[j][d]->SetMaximum(Pmax);
0560       }
0561       cc->Print(TString(plotsDir) + "/_" + titleVar[j] + "/" + TString(p[j][d]->GetName()) + ".pdf");
0562 
0563       for (int i = 0; i < NLEV; ++i) {
0564         for (int p = 0; p < DEP[d]; ++p) {
0565           cc->Clear();
0566           h[j][i][d][p]->Draw("colz");
0567 
0568           if (i == NLEV - 1) {
0569             if (doRatio[j]) {
0570               h[j][2][d][p]->SetMinimum(Zmin);
0571               h[j][2][d][p]->SetMaximum(Zmax);
0572             } else {
0573               h[j][2][d][p]->SetMinimum(-3);
0574               h[j][2][d][p]->SetMaximum(3);
0575             }
0576           }
0577           cc->Print(TString(plotsDir) + "/_" + titleVar[j] + "/" + TString(h[j][i][d][p]->GetName()) + ".pdf");
0578         }
0579       }
0580     }
0581   }
0582 
0583   for (int i = 0; i < NLEV; ++i) {
0584     cc->Clear();
0585     hlut[0][i]->Draw();
0586     hlut[1][i]->Draw("sames");
0587     hlut[2][i]->Draw("sames");
0588     hlut[3][i]->Draw("sames");
0589     cc->Print(TString(plotsDir) + Form("LUT_%d.gif", i));
0590   }
0591   cc->Clear();
0592   hslope[0]->Draw("colz");
0593   cc->SetGridx();
0594   cc->SetGridy();
0595   cc->Print(TString(plotsDir) + "GainLutScatterHBHE.pdf");
0596   cc->Clear();
0597   hslope[1]->Draw("colz");
0598   cc->SetGridx();
0599   cc->SetGridy();
0600   cc->Print(TString(plotsDir) + "GainLutScatterLutHF.pdf");
0601 
0602   for (int i = 0; i < 2; ++i) {
0603     cc->Clear();
0604     houtput[0][i]->Draw("box");
0605     houtput[1][i]->Draw("samebox");
0606     houtput[2][i]->Draw("samebox");
0607     houtput[3][i]->Draw("samebox");
0608     cc->Print(TString(plotsDir) + Form("OUT_%d.gif", i));
0609   }
0610 }
0611 
0612 void HcalLutAnalyzer::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0613   edm::ParameterSetDescription desc;
0614   desc.setUnknown();
0615   descriptions.addDefault(desc);
0616 }
0617 
0618 DEFINE_FWK_MODULE(HcalLutAnalyzer);