Back to home page

Project CMSSW displayed by LXR
 
 

    

File indexing completed on 2021-11-30 10:17:37 

0001 //###########################################################################
0002 //
0003 //   void noseRecHitPlots(fname, dirnm, tag, save)
0004 //        Plots histograms created by recHitStudy for the Nose detector
0005 //   void hgcalStudyPlots(fname, type, tag, dtype, save, debug)
0006 //        Plots histograms created by SimHit/Digi/RecHit Study for HGCal
0007 //   void hgcalGeomCheckPlots(fname, tag, dtype, save, debug)
0008 //        Plots histograms created by HGCGeometryCheck
0009 //   void hgcalBHValidPlots(fname, tag, dtype, save, debug)
0010 //        Plots histograms created by HGCalBHValidation
0011 //   void hgcalSiliconAnalysisPlots(fname, tag, dtype, save, debug)
0012 //        Plots histograms created by HGCalSiliconValidation
0013 //
0014 //   where
0015 //     fanme     (std::string)  Input ROOT file name
0016 //                              ("hfnRecHitD31tt.root" for noseRecHitPlots,
0017 //                               "roots/hgcSimHitD83tt.root" studyPlots,
0018 //                               "roots/hgcGeomCheckD83.root" GeomCheck,
0019 //                               "roots/hgcBHValidD83.root" BHValid)
0020 //     dirnm     (std::string)  Directory name
0021 //                              ("hfnoseRecHitStudy" for noseRecHitPlots)
0022 //     type      (int)          Type: 0 SimHit; 1 Digi; 2 RecHit (0)
0023 //     tag       (std::string)  Name of the tag for the canvas name
0024 //                              ("HFNose" for noseRecHitPlots,
0025 //                               "SimHitD83" for studyPlots,
0026 //                               "GeomChkD83" for GeomCheck,
0027 //                               "BHValidD83" for BHValid)
0028 //     dtype     (std::string)  Data type added for canvas name
0029 //                              ("ttbar D83" for studyPlots,
0030 //                               "#mu D83" for geomCheck, BHValid)
0031 //     save      (bool)         Flag to save the canvas (false)
0032 //     debug     (bool)         Debug flag (false)
0033 //
0034 //###########################################################################
0035 //
0036 #include <TROOT.h>
0037 #include <TChain.h>
0038 #include <TFile.h>
0039 #include <TH1D.h>
0040 #include <TH2D.h>
0041 #include <TProfile.h>
0042 #include <TFitResult.h>
0043 #include <TFitResultPtr.h>
0044 #include <TStyle.h>
0045 #include <TCanvas.h>
0046 #include <TLegend.h>
0047 #include <TPaveStats.h>
0048 #include <TPaveText.h>
0049 #include <vector>
0050 #include <string>
0051 #include <iomanip>
0052 #include <iostream>
0053 #include <fstream>
0054 
0055 void noseRecHitPlots(std::string fname = "hfnRecHitD31tt.root",
0056                      std::string dirnm = "hfnoseRecHitStudy",
0057                      std::string tag = "HFNose",
0058                      bool save = false) {
0059   int nums[2] = {5, 2};
0060   int layers(8);
0061   std::string name1[5] = {
0062       "Energy_Layer", "HitOccupancy_Minus_layer", "HitOccupaancy_Plus_layer", "EtaPhi_Minus_Layer", "EtaPhi_Plus_Layer"};
0063   std::string name2[2] = {"EtaPhi", "RZ"};
0064   std::string detName = "HGCalHFNoseSensitive";
0065   int type1[5] = {1, 1, 1, 2, 2};
0066   int rebin[5] = {10, 1, 1, 1, 1};
0067   int type2[2] = {2, 2};
0068   std::string xtitl1[5] = {"Energy (GeV)", "Hits", "Hits", "#eta", "#eta"};
0069   std::string ytitl1[5] = {" ", " ", " ", "#phi", "#phi"};
0070   std::string xtitl2[2] = {"#eta", "z (cm)"};
0071   std::string ytitl2[2] = {"#phi", "R (cm)"};
0072 
0073   gStyle->SetCanvasBorderMode(0);
0074   gStyle->SetCanvasColor(kWhite);
0075   gStyle->SetPadColor(kWhite);
0076   gStyle->SetFillColor(kWhite);
0077   gStyle->SetOptStat(111110);
0078   TFile *file = new TFile(fname.c_str());
0079   if (file) {
0080     TDirectory *dir = (TDirectory *)file->FindObjectAny(dirnm.c_str());
0081     char name[100], title[100];
0082     for (int i1 = 0; i1 < 2; ++i1) {
0083       int kk = (i1 == 0) ? layers : 1;
0084       for (int i2 = 0; i2 < nums[i1]; ++i2) {
0085         for (int k = 0; k < kk; ++k) {
0086           int type(0);
0087           if (i1 == 0) {
0088             sprintf(name, "%s_%d", name1[i2].c_str(), k + 1);
0089             sprintf(title, "%s (Layer %d)", tag.c_str(), k + 1);
0090             type = type1[i2];
0091           } else {
0092             sprintf(name, "%s_%s", name2[i2].c_str(), detName.c_str());
0093             sprintf(title, "Simulation of %s", tag.c_str());
0094             type = type2[i2];
0095           }
0096           TH1D *hist1(nullptr);
0097           TH2D *hist2(nullptr);
0098           if (type == 1)
0099             hist1 = (TH1D *)dir->FindObjectAny(name);
0100           else
0101             hist2 = (TH2D *)dir->FindObjectAny(name);
0102           if ((hist1 != nullptr) || (hist2 != nullptr)) {
0103             TCanvas *pad = new TCanvas(name, name, 500, 500);
0104             pad->SetRightMargin(0.10);
0105             pad->SetTopMargin(0.10);
0106             if (type == 1) {
0107               if (i1 == 0) {
0108                 hist1->GetYaxis()->SetTitle(ytitl1[i2].c_str());
0109                 hist1->GetXaxis()->SetTitle(xtitl1[i2].c_str());
0110                 hist1->Rebin(rebin[i2]);
0111               } else {
0112                 hist1->GetYaxis()->SetTitle(ytitl2[i2].c_str());
0113                 hist1->GetXaxis()->SetTitle(xtitl2[i2].c_str());
0114               }
0115               hist1->SetTitle(title);
0116               hist1->GetYaxis()->SetTitleOffset(1.2);
0117               pad->SetLogy();
0118               hist1->Draw();
0119             } else {
0120               if (i1 == 0) {
0121                 hist2->GetYaxis()->SetTitle(ytitl1[i2].c_str());
0122                 hist2->GetXaxis()->SetTitle(xtitl1[i2].c_str());
0123               } else {
0124                 hist2->GetYaxis()->SetTitle(ytitl2[i2].c_str());
0125                 hist2->GetXaxis()->SetTitle(xtitl2[i2].c_str());
0126               }
0127               hist2->GetYaxis()->SetTitleOffset(1.2);
0128               hist2->SetMarkerStyle(20);
0129               hist2->SetMarkerSize(0.1);
0130               hist2->SetTitle(title);
0131               hist2->Draw();
0132             }
0133             pad->Update();
0134             TPaveStats *st1 = ((hist1 != nullptr) ? ((TPaveStats *)hist1->GetListOfFunctions()->FindObject("stats"))
0135                                                   : ((TPaveStats *)hist2->GetListOfFunctions()->FindObject("stats")));
0136             if (st1 != NULL) {
0137               st1->SetY1NDC(0.70);
0138               st1->SetY2NDC(0.90);
0139               st1->SetX1NDC(0.65);
0140               st1->SetX2NDC(0.90);
0141             }
0142             pad->Modified();
0143             pad->Update();
0144             if (save) {
0145               sprintf(name, "c_%s%s.jpg", tag.c_str(), pad->GetName());
0146               pad->Print(name);
0147             }
0148           }
0149         }
0150       }
0151     }
0152   }
0153 }
0154 
0155 void hgcalStudyPlots(std::string fname = "roots/hgcSimHitD83tt.root",
0156                      int type = 0,
0157                      std::string tag = "SimHitD83",
0158                      std::string dtype = "ttbar D83",
0159                      bool save = false,
0160                      bool debug = false) {
0161   int ndir[3] = {1, 3, 3};
0162   std::string dirnm[3][3] = {{"hgcalSimHitStudy", "", ""},
0163                              {"hgcalDigiStudyEE", "hgcalDigiStudyHEF", "hgcalDigiStudyHEB"},
0164                              {"hgcalRecHitStudyEE", "hgcalRecHitStudyFH", "hgcalRecHitStudyBH"}};
0165   std::string name0[4] = {"HGCal EE", "HGCal HE Silicon", "HGCal HE Scintillator", "HGCal"};
0166   int nname[3] = {4, 1, 1};
0167   std::string name1[4] = {
0168       "HGCalEESensitive", "HGCalHESiliconSensitive", "HGCalHEScintillatorSensitive", "AllDetectors"};
0169   int nhist[3] = {4, 2, 2};
0170   int nhtype[3][4] = {{1, 1, 2, 2}, {1, 2, 0, 0}, {2, 2, 0, 0}};
0171   int yax[3][4] = {{1, 1, 0, 0}, {1, 0, 0, 0}, {0, 0, 0, 0}};
0172   std::string name2[3][4] = {
0173       {"E_", "T_", "RZ_", "EtaPhi_"}, {"Charge_", "RZ_", " ", " "}, {"RZ_", "EtaPhi_", " ", " "}};
0174   std::string xtitl[3][4] = {
0175       {"Energy (GeV)", "Time (ns)", "Z (mm)", "#phi"}, {"Charge", "Z (cm)", " ", " "}, {"Z (cm)", "#phi", " ", " "}};
0176   std::string ytitl[3][4] = {
0177       {"Hits", "Hits", "R (mm)", "#eta"}, {"Hits", "R (cm)", " ", " "}, {"R (cm)", "#eta", " ", " "}};
0178   double xmax[3][4] = {{0.2, 20.0, -1, -1}, {25.0, -1, -1, -1}, {-1, -1, -1, -1}};
0179 
0180   gStyle->SetCanvasBorderMode(0);
0181   gStyle->SetCanvasColor(kWhite);
0182   gStyle->SetPadColor(kWhite);
0183   gStyle->SetFillColor(kWhite);
0184   gStyle->SetOptStat(111110);
0185   if (debug)
0186     std::cout << "File " << fname << " Type " << type << ":" << name0[type] << " Tags " << tag << " : " << dtype
0187               << " Save " << save << "\n";
0188   TFile *file = new TFile(fname.c_str());
0189   if (file) {
0190     for (int it = 0; it < ndir[type]; ++it) {
0191       char dirx[100];
0192       sprintf(dirx, "%s", dirnm[type][it].c_str());
0193       TDirectory *dir = (TDirectory *)file->FindObjectAny(dirx);
0194       if (debug)
0195         std::cout << "Directory " << dirx << " : " << dir << std::endl;
0196       if (dir) {
0197         for (int in = 0; in < nname[type]; ++in) {
0198           for (int ih = 0; ih < nhist[type]; ++ih) {
0199             char hname[100];
0200             if (type == 0) {
0201               sprintf(hname, "%s%s", name2[type][ih].c_str(), name1[in].c_str());
0202             } else {
0203               sprintf(hname, "%s%s", name2[type][ih].c_str(), name1[it].c_str());
0204             }
0205             TH1D *hist1(nullptr);
0206             TH2D *hist2(nullptr);
0207             if (nhtype[type][ih] == 1)
0208               hist1 = (TH1D *)dir->FindObjectAny(hname);
0209             else
0210               hist2 = (TH2D *)dir->FindObjectAny(hname);
0211             if (debug)
0212               std::cout << "Hist " << hname << " : " << hist1 << " : " << hist2 << " Xtitle " << xtitl[type][ih]
0213                         << " Ytitle " << ytitl[type][ih] << " xmax " << xmax[type][ih] << " Type " << nhtype[type][ih]
0214                         << " yscale " << yax[type][ih] << std::endl;
0215             if ((hist1 != nullptr) || (hist2 != nullptr)) {
0216               char name[100], title[100];
0217               sprintf(name, "%s%s", hname, tag.c_str());
0218               TCanvas *pad = new TCanvas(name, name, 500, 500);
0219               pad->SetRightMargin(0.10);
0220               pad->SetTopMargin(0.10);
0221               if (type == 0)
0222                 sprintf(title, "%s (%s)", name0[in].c_str(), dtype.c_str());
0223               else
0224                 sprintf(title, "%s (%s)", name0[it].c_str(), dtype.c_str());
0225               if (debug)
0226                 std::cout << "Pad " << name << " : " << pad << "\n";
0227               if (nhtype[type][ih] == 1) {
0228                 hist1->GetYaxis()->SetTitle(ytitl[type][ih].c_str());
0229                 hist1->GetXaxis()->SetTitle(xtitl[type][ih].c_str());
0230                 if (xmax[type][ih] > 0)
0231                   hist1->GetXaxis()->SetRangeUser(0, xmax[type][ih]);
0232                 if (yax[type][ih] > 0)
0233                   pad->SetLogy();
0234                 hist1->SetTitle(title);
0235                 hist1->GetYaxis()->SetTitleOffset(1.2);
0236                 hist1->Draw();
0237               } else {
0238                 hist2->GetYaxis()->SetTitle(ytitl[type][ih].c_str());
0239                 hist2->GetXaxis()->SetTitle(xtitl[type][ih].c_str());
0240                 hist2->GetYaxis()->SetTitleOffset(1.2);
0241                 hist2->SetMarkerStyle(20);
0242                 hist2->SetMarkerSize(0.1);
0243                 hist2->SetTitle(title);
0244                 hist2->Draw();
0245               }
0246               pad->Update();
0247               TPaveStats *st1 = ((hist1 != nullptr) ? ((TPaveStats *)hist1->GetListOfFunctions()->FindObject("stats"))
0248                                                     : ((TPaveStats *)hist2->GetListOfFunctions()->FindObject("stats")));
0249               if (st1 != NULL) {
0250                 st1->SetY1NDC(0.70);
0251                 st1->SetY2NDC(0.90);
0252                 st1->SetX1NDC(0.65);
0253                 st1->SetX2NDC(0.90);
0254               }
0255               pad->Modified();
0256               pad->Update();
0257               if (save) {
0258                 sprintf(name, "c_%s.jpg", pad->GetName());
0259                 pad->Print(name);
0260               }
0261             }
0262           }
0263         }
0264       }
0265     }
0266   }
0267 }
0268 
0269 void hgcalGeomCheckPlots(std::string fname = "roots/hgcGeomCheckD83.root",
0270                          std::string tag = "GeomChkD83",
0271                          std::string dtype = "#mu D83",
0272                          bool statbox = true,
0273                          bool save = false,
0274                          bool debug = false) {
0275   std::string dirnm = "hgcGeomCheck";
0276   std::string name0[3] = {"HGCal EE", "HGCal HE Silicon", "HGCal HE Scintillator"};
0277   int nhist = 3;
0278   int nhtype = 3;
0279   std::string name2[9] = {"heerVsLayer",
0280                           "heezVsLayer",
0281                           "heedzVsZ",
0282                           "hefrVsLayer",
0283                           "hefzVsLayer",
0284                           "hefdzVsZ",
0285                           "hebrVsLayer",
0286                           "hebzVsLayer",
0287                           "hebdzVsZ"};
0288   std::string xtitl[9] = {"Layer", "Layer", "Z (cm)", "Layer", "Layer", "Z (cm)", "Layer", "Layer", "Z (cm)"};
0289   std::string ytitl[9] = {
0290       "R (cm)", "Z (cm)", "#Delta Z (cm)", "R (cm)", "Z (cm)", "#Delta Z (cm)", "R (cm)", "Z (cm)", "#Delta Z (cm)"};
0291 
0292   gStyle->SetCanvasBorderMode(0);
0293   gStyle->SetCanvasColor(kWhite);
0294   gStyle->SetPadColor(kWhite);
0295   gStyle->SetFillColor(kWhite);
0296   if (statbox)
0297     gStyle->SetOptStat(111110);
0298   else
0299     gStyle->SetOptStat(0);
0300   if (debug)
0301     std::cout << "File " << fname << " Tags " << tag << " : " << dtype << " Save " << save << "\n";
0302   TFile *file = new TFile(fname.c_str());
0303   if (file) {
0304     char dirx[100];
0305     sprintf(dirx, "%s", dirnm.c_str());
0306     TDirectory *dir = (TDirectory *)file->FindObjectAny(dirx);
0307     if (debug)
0308       std::cout << "Directory " << dirx << " : " << dir << std::endl;
0309     if (dir) {
0310       for (int in = 0; in < nhtype; ++in) {
0311         for (int ih = 0; ih < nhist; ++ih) {
0312           char hname[100];
0313           int inh = in * nhist + ih;
0314           sprintf(hname, "%s", name2[inh].c_str());
0315           TH2D *hist = (TH2D *)dir->FindObjectAny(hname);
0316           if (debug)
0317             std::cout << "Hist " << hname << " : " << hist << " Xtitle " << xtitl[inh] << " Ytitle " << ytitl[inh]
0318                       << std::endl;
0319           if (hist != nullptr) {
0320             char name[100], title[100];
0321             sprintf(name, "%s%s", hname, tag.c_str());
0322             TCanvas *pad = new TCanvas(name, name, 500, 500);
0323             pad->SetRightMargin(0.10);
0324             pad->SetTopMargin(0.10);
0325             sprintf(title, "%s (%s)", name0[in].c_str(), dtype.c_str());
0326             if (debug)
0327               std::cout << "Pad " << name << " : " << pad << "\n";
0328             hist->GetYaxis()->SetTitle(ytitl[inh].c_str());
0329             hist->GetXaxis()->SetTitle(xtitl[inh].c_str());
0330             hist->SetTitle(title);
0331             hist->GetYaxis()->SetTitleOffset(1.2);
0332             hist->Draw();
0333             pad->Update();
0334             TPaveStats *st1 = (TPaveStats *)hist->GetListOfFunctions()->FindObject("stats");
0335             if (st1 != NULL) {
0336               st1->SetY1NDC(0.70);
0337               st1->SetY2NDC(0.90);
0338               st1->SetX1NDC(0.65);
0339               st1->SetX2NDC(0.90);
0340             }
0341             pad->Modified();
0342             pad->Update();
0343             if (save) {
0344               sprintf(name, "c_%s.jpg", pad->GetName());
0345               pad->Print(name);
0346             }
0347           }
0348         }
0349       }
0350     }
0351   }
0352 }
0353 
0354 void hgcalBHValidPlots(std::string fname = "roots/hgcBHValidD83.root",
0355                        std::string tag = "BHValidD83",
0356                        std::string dtype = "#mu D83",
0357                        bool save = false,
0358                        bool debug = false) {
0359   std::string dirnm = "hgcalBHAnalysis";
0360   std::string name0 = "HGCal HE Scintillator";
0361   int nhist = 9;
0362   std::string name2[9] = {"SimHitEn1",
0363                           "SimHitEn2",
0364                           "SimHitLong",
0365                           "SimHitOccup",
0366                           "SimHitOccu3",
0367                           "SimHitTime",
0368                           "DigiLong",
0369                           "DigiOccup",
0370                           "DigiOccu3"};
0371   std::string xtitl[9] = {"SimHit Energy (GeV)",
0372                           "SimHit Energy (GeV)",
0373                           "SimHit Layer #",
0374                           "SimHit i#eta",
0375                           "SimHit i#eta",
0376                           "Digi Layer #",
0377                           "Digi i#eta",
0378                           "Digi i#eta"};
0379   std::string ytitl[9] = {"Hits", "Hits", "Energy Sum (GeV)", "i#phi", "Layer #", "Digi Sum", "i#phi", "Layer #"};
0380   double xmax[9] = {0.05, 0.20, -1.0, -1.0, -1.0, -1.0, -1.0, -1.0, -1.0};
0381   int ihty[9] = {1, 1, 1, 2, 2, 1, 1, 2, 2};
0382   int iaxty[9] = {1, 1, 0, 0, 0, 1, 0, 0, 0};
0383   int ibin[10] = {0, 0, 0, 0, 0, 10, 0, 0, 0};
0384   gStyle->SetCanvasBorderMode(0);
0385   gStyle->SetCanvasColor(kWhite);
0386   gStyle->SetPadColor(kWhite);
0387   gStyle->SetFillColor(kWhite);
0388   gStyle->SetOptStat(111110);
0389   if (debug)
0390     std::cout << "File " << fname << " Tags " << tag << " : " << dtype << " Save " << save << "\n";
0391   TFile *file = new TFile(fname.c_str());
0392   if (file) {
0393     char dirx[100];
0394     sprintf(dirx, "%s", dirnm.c_str());
0395     TDirectory *dir = (TDirectory *)file->FindObjectAny(dirx);
0396     if (debug)
0397       std::cout << "Directory " << dirx << " : " << dir << std::endl;
0398     if (dir) {
0399       for (int ih = 0; ih < nhist; ++ih) {
0400         char hname[100];
0401         sprintf(hname, "%s", name2[ih].c_str());
0402         TH1D *hist1(nullptr);
0403         TH2D *hist2(nullptr);
0404         if (ihty[ih] <= 1)
0405           hist1 = (TH1D *)dir->FindObjectAny(hname);
0406         else
0407           hist2 = (TH2D *)dir->FindObjectAny(hname);
0408         if (debug)
0409           std::cout << "Hist " << hname << " : " << hist1 << ":" << hist2 << " Xtitle " << xtitl[ih] << " Ytitle "
0410                     << ytitl[ih] << " xmax " << xmax[ih] << " iaxty " << iaxty[ih] << " ibin " << ibin[ih] << std::endl;
0411         if ((hist1 != nullptr) || (hist2 != nullptr)) {
0412           char name[100], title[100];
0413           sprintf(name, "%s%s", hname, tag.c_str());
0414           TCanvas *pad = new TCanvas(name, name, 500, 500);
0415           pad->SetRightMargin(0.10);
0416           pad->SetTopMargin(0.10);
0417           sprintf(title, "%s (%s)", name0.c_str(), dtype.c_str());
0418           if (debug)
0419             std::cout << "Pad " << name << " : " << pad << "\n";
0420           if (hist1 != nullptr) {
0421             hist1->GetYaxis()->SetTitle(ytitl[ih].c_str());
0422             hist1->GetXaxis()->SetTitle(xtitl[ih].c_str());
0423             hist1->SetTitle(title);
0424             if (xmax[ih] > 0)
0425               hist1->GetXaxis()->SetRangeUser(0, xmax[ih]);
0426             if (iaxty[ih] > 0)
0427               pad->SetLogy();
0428             if (ibin[ih] > 0)
0429               hist1->Rebin(ibin[ih]);
0430             hist1->GetYaxis()->SetTitleOffset(1.2);
0431             hist1->Draw();
0432           } else {
0433             hist2->GetYaxis()->SetTitle(ytitl[ih].c_str());
0434             hist2->GetXaxis()->SetTitle(xtitl[ih].c_str());
0435             hist2->SetTitle(title);
0436             hist2->GetYaxis()->SetTitleOffset(1.2);
0437             hist2->Draw();
0438           }
0439           pad->Update();
0440           TPaveStats *st1 = ((hist1 != nullptr) ? ((TPaveStats *)hist1->GetListOfFunctions()->FindObject("stats"))
0441                                                 : ((TPaveStats *)hist2->GetListOfFunctions()->FindObject("stats")));
0442           if (st1 != NULL) {
0443             st1->SetY1NDC(0.70);
0444             st1->SetY2NDC(0.90);
0445             st1->SetX1NDC(0.65);
0446             st1->SetX2NDC(0.90);
0447           }
0448           pad->Modified();
0449           pad->Update();
0450           if (save) {
0451             sprintf(name, "c_%s.jpg", pad->GetName());
0452             pad->Print(name);
0453           }
0454         }
0455       }
0456     }
0457   }
0458 }
0459 
0460 void hgcalSiliconAnalysisPlots(std::string fname = "roots/hgcSilValidD86.root",
0461                                std::string tag = "SilValidD86",
0462                                std::string dtype = "ttbar D86",
0463                                bool save = false,
0464                                bool debug = false) {
0465   std::string dirnm[2] = {"hgcalSiliconAnalysisEE", "hgcalSiliconAnalysisHEF"};
0466   std::string name0[2] = {"HGCal EE", "HGCal HE Silicon"};
0467   int nhist = 9;
0468   std::string name2[9] = {"SimHitEn1",
0469                           "SimHitEn2",
0470                           "SimHitLong",
0471                           "SimHitOccup",
0472                           "SimHitOccu2",
0473                           "SimHitTime",
0474                           "DigiLong",
0475                           "DigiOccup",
0476                           "DigiOccu2"};
0477   std::string xtitl[9] = {"SimHit Energy (GeV)",
0478                           "SimHit Energy (GeV)",
0479                           "SimHit Layer #",
0480                           "SimHit i#eta",
0481                           "SimHit i#eta",
0482                           "Digi Layer #",
0483                           "Digi i#eta",
0484                           "Digi i#eta"};
0485   std::string ytitl[9] = {"Hits", "Hits", "Energy Sum (GeV)", "i#phi", "Layer #", "Digi Sum", "i#phi", "Layer #"};
0486   double xmax[9] = {0.05, 0.20, -1.0, -1.0, -1.0, -1.0, -1.0, -1.0, -1.0};
0487   int ihty[9] = {1, 1, 1, 2, 2, 1, 1, 2, 2};
0488   int iaxty[9] = {1, 1, 0, 0, 0, 1, 0, 0, 0};
0489   int ibin[10] = {0, 0, 0, 0, 0, 10, 0, 0, 0};
0490   gStyle->SetCanvasBorderMode(0);
0491   gStyle->SetCanvasColor(kWhite);
0492   gStyle->SetPadColor(kWhite);
0493   gStyle->SetFillColor(kWhite);
0494   gStyle->SetOptStat(111110);
0495   if (debug)
0496     std::cout << "File " << fname << " Tags " << tag << " : " << dtype << " Save " << save << "\n";
0497   TFile *file = new TFile(fname.c_str());
0498   if (file) {
0499     for (int idir = 0; idir < 2; ++idir) {
0500       char dirx[100];
0501       sprintf(dirx, "%s", dirnm[idir].c_str());
0502       TDirectory *dir = (TDirectory *)file->FindObjectAny(dirx);
0503       if (debug)
0504         std::cout << "Directory " << dirx << " : " << dir << std::endl;
0505       if (dir) {
0506         for (int ih = 0; ih < nhist; ++ih) {
0507           char hname[100];
0508           sprintf(hname, "%s", name2[ih].c_str());
0509           TH1D *hist1(nullptr);
0510           TH2D *hist2(nullptr);
0511           if (ihty[ih] <= 1)
0512             hist1 = (TH1D *)dir->FindObjectAny(hname);
0513           else
0514             hist2 = (TH2D *)dir->FindObjectAny(hname);
0515           if (debug)
0516             std::cout << "Hist " << hname << " : " << hist1 << ":" << hist2 << " Xtitle " << xtitl[ih] << " Ytitle "
0517                       << ytitl[ih] << " xmax " << xmax[ih] << " iaxty " << iaxty[ih] << " ibin " << ibin[ih]
0518                       << std::endl;
0519           if ((hist1 != nullptr) || (hist2 != nullptr)) {
0520             char name[100], title[100];
0521             sprintf(name, "%s%s", hname, tag.c_str());
0522             TCanvas *pad = new TCanvas(name, name, 500, 500);
0523             pad->SetRightMargin(0.10);
0524             pad->SetTopMargin(0.10);
0525             sprintf(title, "%s (%s)", name0[idir].c_str(), dtype.c_str());
0526             if (debug)
0527               std::cout << "Pad " << name << " : " << pad << "\n";
0528             if (hist1 != nullptr) {
0529               hist1->GetYaxis()->SetTitle(ytitl[ih].c_str());
0530               hist1->GetXaxis()->SetTitle(xtitl[ih].c_str());
0531               hist1->SetTitle(title);
0532               if (xmax[ih] > 0)
0533                 hist1->GetXaxis()->SetRangeUser(0, xmax[ih]);
0534               if (iaxty[ih] > 0)
0535                 pad->SetLogy();
0536               if (ibin[ih] > 0)
0537                 hist1->Rebin(ibin[ih]);
0538               hist1->GetYaxis()->SetTitleOffset(1.2);
0539               hist1->Draw();
0540             } else {
0541               hist2->GetYaxis()->SetTitle(ytitl[ih].c_str());
0542               hist2->GetXaxis()->SetTitle(xtitl[ih].c_str());
0543               hist2->SetTitle(title);
0544               hist2->GetYaxis()->SetTitleOffset(1.2);
0545               hist2->Draw();
0546             }
0547             pad->Update();
0548             TPaveStats *st1 = ((hist1 != nullptr) ? ((TPaveStats *)hist1->GetListOfFunctions()->FindObject("stats"))
0549                                                   : ((TPaveStats *)hist2->GetListOfFunctions()->FindObject("stats")));
0550             if (st1 != NULL) {
0551               st1->SetY1NDC(0.70);
0552               st1->SetY2NDC(0.90);
0553               st1->SetX1NDC(0.65);
0554               st1->SetX2NDC(0.90);
0555             }
0556             pad->Modified();
0557             pad->Update();
0558             if (save) {
0559               sprintf(name, "c_%s.jpg", pad->GetName());
0560               pad->Print(name);
0561             }
0562           }
0563         }
0564       }
0565     }
0566   }
0567 }
0568 
0569 void ehcalPlots(std::string fname = "ecalHitdd4hep.root",
0570                 std::string dirnm = "EcalSimHitStudy",
0571                 std::string tag = "DD4Hep",
0572                 int dets = 2,
0573                 bool save = false,
0574                 bool debug = false) {
0575   const int nh = 2;
0576   std::string name[nh] = {"poszp", "poszn"};
0577   gStyle->SetCanvasBorderMode(0);
0578   gStyle->SetCanvasColor(kWhite);
0579   gStyle->SetPadColor(kWhite);
0580   gStyle->SetFillColor(kWhite);
0581   gStyle->SetOptStat(10);
0582   if (debug)
0583     std::cout << "File " << fname << " Tag " << tag << " Detectors " << dets << " Save " << save << "\n";
0584   TFile *file = new TFile(fname.c_str());
0585   if (file) {
0586     TDirectory *dir = (TDirectory *)file->FindObjectAny(dirnm.c_str());
0587     if (debug)
0588       std::cout << "Directory " << dirnm << " : " << dir << std::endl;
0589     if (dir) {
0590       for (int id = 0; id < dets; ++id) {
0591         for (int ih = 0; ih < nh; ++ih) {
0592           char hname[100];
0593           sprintf(hname, "%s%d", name[ih].c_str(), id);
0594           TH2D *hist = (TH2D *)(dir->FindObjectAny(hname));
0595           if (debug)
0596             std::cout << "Hist " << hname << " : " << hist << ":" << hist->GetTitle() << " Xtitle "
0597                       << hist->GetXaxis()->GetTitle() << " Ytitle " << hist->GetYaxis()->GetTitle() << std::endl;
0598           if (hist != nullptr) {
0599             char cname[100], title[100];
0600             sprintf(cname, "%s%s", hname, tag.c_str());
0601             TCanvas *pad = new TCanvas(cname, cname, 500, 500);
0602             pad->SetRightMargin(0.10);
0603             pad->SetTopMargin(0.10);
0604             sprintf(title, "%s (%s)", hist->GetTitle(), tag.c_str());
0605             if (debug)
0606               std::cout << "Pad " << cname << " : " << pad << "\n";
0607             hist->SetTitle(title);
0608             hist->SetMarkerStyle(20);
0609             hist->SetMarkerSize(0.1);
0610             hist->GetYaxis()->SetTitleOffset(1.3);
0611             hist->Draw();
0612             pad->Update();
0613             TPaveStats *st1 = ((TPaveStats *)hist->GetListOfFunctions()->FindObject("stats"));
0614             if (st1 != NULL) {
0615               st1->SetY1NDC(0.85);
0616               st1->SetY2NDC(0.90);
0617               st1->SetX1NDC(0.70);
0618               st1->SetX2NDC(0.90);
0619             }
0620             pad->Modified();
0621             pad->Update();
0622             if (save) {
0623               sprintf(cname, "c_%s.jpg", pad->GetName());
0624               pad->Print(cname);
0625             }
0626           }
0627         }
0628       }
0629     }
0630   }
0631 }
0632 
0633 class hgcHits {
0634 public:
0635   TTree *fChain;
0636   Int_t fCurrent;
0637 
0638 private:
0639   // Declaration of leaf types
0640   std::vector<float> *heeRecX;
0641   std::vector<float> *heeRecY;
0642   std::vector<float> *heeRecZ;
0643   std::vector<float> *heeRecEnergy;
0644   std::vector<float> *hefRecX;
0645   std::vector<float> *hefRecY;
0646   std::vector<float> *hefRecZ;
0647   std::vector<float> *hefRecEnergy;
0648   std::vector<float> *hebRecX;
0649   std::vector<float> *hebRecY;
0650   std::vector<float> *hebRecZ;
0651   std::vector<float> *hebRecEta;
0652   std::vector<float> *hebRecPhi;
0653   std::vector<float> *hebRecEnergy;
0654   std::vector<float> *heeSimX;
0655   std::vector<float> *heeSimY;
0656   std::vector<float> *heeSimZ;
0657   std::vector<float> *heeSimEnergy;
0658   std::vector<float> *hefSimX;
0659   std::vector<float> *hefSimY;
0660   std::vector<float> *hefSimZ;
0661   std::vector<float> *hefSimEnergy;
0662   std::vector<float> *hebSimX;
0663   std::vector<float> *hebSimY;
0664   std::vector<float> *hebSimZ;
0665   std::vector<float> *hebSimEta;
0666   std::vector<float> *hebSimPhi;
0667   std::vector<float> *hebSimEnergy;
0668   std::vector<unsigned int> *heeDetID;
0669   std::vector<unsigned int> *hefDetID;
0670   std::vector<unsigned int> *hebDetID;
0671 
0672   // List of branches
0673   TBranch *b_heeRecX;
0674   TBranch *b_heeRecY;
0675   TBranch *b_heeRecZ;
0676   TBranch *b_heeRecEnergy;
0677   TBranch *b_hefRecX;
0678   TBranch *b_hefRecY;
0679   TBranch *b_hefRecZ;
0680   TBranch *b_hefRecEnergy;
0681   TBranch *b_hebRecX;
0682   TBranch *b_hebRecY;
0683   TBranch *b_hebRecZ;
0684   TBranch *b_hebRecEta;
0685   TBranch *b_hebRecPhi;
0686   TBranch *b_hebRecEnergy;
0687   TBranch *b_heeSimX;
0688   TBranch *b_heeSimY;
0689   TBranch *b_heeSimZ;
0690   TBranch *b_heeSimEnergy;
0691   TBranch *b_hefSimX;
0692   TBranch *b_hefSimY;
0693   TBranch *b_hefSimZ;
0694   TBranch *b_hefSimEnergy;
0695   TBranch *b_hebSimX;
0696   TBranch *b_hebSimY;
0697   TBranch *b_hebSimZ;
0698   TBranch *b_hebSimEta;
0699   TBranch *b_hebSimPhi;
0700   TBranch *b_hebSimEnergy;
0701   TBranch *b_heeDetID;
0702   TBranch *b_hefDetID;
0703   TBranch *b_hebDetID;
0704 
0705 public:
0706   hgcHits(const char *infile);
0707   virtual ~hgcHits();
0708   virtual Int_t Cut(Long64_t entry);
0709   virtual Int_t GetEntry(Long64_t entry);
0710   virtual Long64_t LoadTree(Long64_t entry);
0711   virtual void Init(TTree *tree);
0712   virtual void Loop();
0713   virtual Bool_t Notify();
0714   virtual void Show(Long64_t entry = -1);
0715 
0716   void bookHistograms();
0717   void saveHistos(const char *outfile);
0718 
0719 private:
0720 };
0721 
0722 hgcHits::hgcHits::hgcHits(const char *infile) {
0723   TFile *file = new TFile(infile);
0724   TDirectory *dir = (TDirectory *)(file->FindObjectAny("hgcHitAnalysis"));
0725   TTree *tree = (TTree *)(dir->FindObjectAny("hgcHits"));
0726   std::cout << "Attaches tree hgcHits at " << tree << " in file " << infile << std::endl;
0727   bookHistograms();
0728   Init(tree);
0729 }
0730 
0731 hgcHits::~hgcHits() {
0732   if (!fChain)
0733     return;
0734   delete fChain->GetCurrentFile();
0735 }
0736 
0737 Int_t hgcHits::GetEntry(Long64_t entry) {
0738   // Read contents of entry.
0739   if (!fChain)
0740     return 0;
0741   return fChain->GetEntry(entry);
0742 }
0743 
0744 Long64_t hgcHits::LoadTree(Long64_t entry) {
0745   // Set the environment to read one entry
0746   if (!fChain)
0747     return -5;
0748   Long64_t centry = fChain->LoadTree(entry);
0749   if (centry < 0)
0750     return centry;
0751   if (!fChain->InheritsFrom(TChain::Class()))
0752     return centry;
0753   TChain *chain = (TChain *)fChain;
0754   if (chain->GetTreeNumber() != fCurrent) {
0755     fCurrent = chain->GetTreeNumber();
0756     Notify();
0757   }
0758   return centry;
0759 }
0760 
0761 void hgcHits::Init(TTree *tree) {
0762   // The Init() function is called when the selector needs to initialize
0763   // a new tree or chain. Typically here the branch addresses and branch
0764   // pointers of the tree will be set.
0765   // It is normally not necessary to make changes to the generated
0766   // code, but the routine can be extended by the user if needed.
0767   // Init() will be called many times when running on PROOF
0768   // (once per file to be processed).
0769 
0770   // Set object pointer
0771   heeRecX = 0;
0772   heeRecY = 0;
0773   heeRecZ = 0;
0774   heeRecEnergy = 0;
0775   hefRecX = 0;
0776   hefRecY = 0;
0777   hefRecZ = 0;
0778   hefRecEnergy = 0;
0779   hebRecX = 0;
0780   hebRecY = 0;
0781   hebRecZ = 0;
0782   hebRecEta = 0;
0783   hebRecPhi = 0;
0784   hebRecEnergy = 0;
0785   heeSimX = 0;
0786   heeSimY = 0;
0787   heeSimZ = 0;
0788   heeSimEnergy = 0;
0789   hefSimX = 0;
0790   hefSimY = 0;
0791   hefSimZ = 0;
0792   hefSimEnergy = 0;
0793   hebSimX = 0;
0794   hebSimY = 0;
0795   hebSimZ = 0;
0796   hebSimEta = 0;
0797   hebSimPhi = 0;
0798   hebSimEnergy = 0;
0799   heeDetID = 0;
0800   hefDetID = 0;
0801   hebDetID = 0;
0802   // Set branch addresses and branch pointers
0803   if (!tree)
0804     return;
0805   fChain = tree;
0806   fCurrent = -1;
0807   fChain->SetMakeClass(1);
0808 
0809   fChain->SetBranchAddress("heeRecX", &heeRecX, &b_heeRecX);
0810   fChain->SetBranchAddress("heeRecY", &heeRecY, &b_heeRecY);
0811   fChain->SetBranchAddress("heeRecZ", &heeRecZ, &b_heeRecZ);
0812   fChain->SetBranchAddress("heeRecEnergy", &heeRecEnergy, &b_heeRecEnergy);
0813   fChain->SetBranchAddress("hefRecX", &hefRecX, &b_hefRecX);
0814   fChain->SetBranchAddress("hefRecY", &hefRecY, &b_hefRecY);
0815   fChain->SetBranchAddress("hefRecZ", &hefRecZ, &b_hefRecZ);
0816   fChain->SetBranchAddress("hefRecEnergy", &hefRecEnergy, &b_hefRecEnergy);
0817   fChain->SetBranchAddress("hebRecX", &hebRecX, &b_hebRecX);
0818   fChain->SetBranchAddress("hebRecY", &hebRecY, &b_hebRecY);
0819   fChain->SetBranchAddress("hebRecZ", &hebRecZ, &b_hebRecZ);
0820   fChain->SetBranchAddress("hebRecEta", &hebRecEta, &b_hebRecEta);
0821   fChain->SetBranchAddress("hebRecPhi", &hebRecPhi, &b_hebRecPhi);
0822   fChain->SetBranchAddress("hebRecEnergy", &hebRecEnergy, &b_hebRecEnergy);
0823   fChain->SetBranchAddress("heeSimX", &heeSimX, &b_heeSimX);
0824   fChain->SetBranchAddress("heeSimY", &heeSimY, &b_heeSimY);
0825   fChain->SetBranchAddress("heeSimZ", &heeSimZ, &b_heeSimZ);
0826   fChain->SetBranchAddress("heeSimEnergy", &heeSimEnergy, &b_heeSimEnergy);
0827   fChain->SetBranchAddress("hefSimX", &hefSimX, &b_hefSimX);
0828   fChain->SetBranchAddress("hefSimY", &hefSimY, &b_hefSimY);
0829   fChain->SetBranchAddress("hefSimZ", &hefSimZ, &b_hefSimZ);
0830   fChain->SetBranchAddress("hefSimEnergy", &hefSimEnergy, &b_hefSimEnergy);
0831   fChain->SetBranchAddress("hebSimX", &hebSimX, &b_hebSimX);
0832   fChain->SetBranchAddress("hebSimY", &hebSimY, &b_hebSimY);
0833   fChain->SetBranchAddress("hebSimZ", &hebSimZ, &b_hebSimZ);
0834   fChain->SetBranchAddress("hebSimEta", &hebSimEta, &b_hebSimEta);
0835   fChain->SetBranchAddress("hebSimPhi", &hebSimPhi, &b_hebSimPhi);
0836   fChain->SetBranchAddress("hebSimEnergy", &hebSimEnergy, &b_hebSimEnergy);
0837   fChain->SetBranchAddress("heeDetID", &heeDetID, &b_heeDetID);
0838   fChain->SetBranchAddress("hefDetID", &hefDetID, &b_hefDetID);
0839   fChain->SetBranchAddress("hebDetID", &hebDetID, &b_hebDetID);
0840   Notify();
0841 }
0842 
0843 Bool_t hgcHits::Notify() {
0844   // The Notify() function is called when a new file is opened. This
0845   // can be either for a new TTree in a TChain or when when a new TTree
0846   // is started when using PROOF. It is normally not necessary to make changes
0847   // to the generated code, but the routine can be extended by the
0848   // user if needed. The return value is currently not used.
0849 
0850   return kTRUE;
0851 }
0852 
0853 void hgcHits::Show(Long64_t entry) {
0854   // Print contents of entry.
0855   // If entry is not specified, print current entry
0856   if (!fChain)
0857     return;
0858   fChain->Show(entry);
0859 }
0860 
0861 Int_t hgcHits::Cut(Long64_t) {
0862   // This function may be called from Loop.
0863   // returns  1 if entry is accepted.
0864   // returns -1 otherwise.
0865   return 1;
0866 }
0867 
0868 void hgcHits::Loop() {
0869   //   In a ROOT session, you can do:
0870   //      Root > .L hgcHits.C
0871   //      Root > hgcHits t
0872   //      Root > t.GetEntry(12); // Fill t data members with entry number 12
0873   //      Root > t.Show();       // Show values of entry 12
0874   //      Root > t.Show(16);     // Read and show values of entry 16
0875   //      Root > t.Loop();       // Loop on all entries
0876   //
0877 
0878   //     This is the loop skeleton where:
0879   //    jentry is the global entry number in the chain
0880   //    ientry is the entry number in the current Tree
0881   //  Note that the argument to GetEntry must be:
0882   //    jentry for TChain::GetEntry
0883   //    ientry for TTree::GetEntry and TBranch::GetEntry
0884   //
0885   //       To read only selected branches, Insert statements like:
0886   // METHOD1:
0887   //    fChain->SetBranchStatus("*",0);  // disable all branches
0888   //    fChain->SetBranchStatus("branchname",1);  // activate branchname
0889   // METHOD2: replace line
0890   //    fChain->GetEntry(jentry);       //read all branches
0891   //by  b_branchname->GetEntry(ientry); //read only this branch
0892   if (fChain == 0)
0893     return;
0894 
0895   Long64_t nentries = fChain->GetEntriesFast();
0896 
0897   Long64_t nbytes = 0, nb = 0;
0898   for (Long64_t jentry = 0; jentry < nentries; jentry++) {
0899     Long64_t ientry = LoadTree(jentry);
0900     if (ientry < 0)
0901       break;
0902     nb = fChain->GetEntry(jentry);
0903     nbytes += nb;
0904     // if (Cut(ientry) < 0) continue;
0905   }
0906 }
0907 
0908 void hgcHits::bookHistograms() {}
0909 
0910 void hgcHits::saveHistos(const char *) {}
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated by the 2.2.1 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!