File indexing completed on 2024-11-27 03:17:28
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014
0015
0016
0017
0018
0019
0020 #include <memory>
0021 #include <iostream>
0022 #include <fstream>
0023
0024
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;
0110 static const int NDET = 2;
0111 static const int NDEP = 7;
0112 static const int NLEV = 3;
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
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
0226 for (int ii = 0; ii < 2; ++ii) {
0227
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
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
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
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
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
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
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
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)) {
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
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);