EcalLaserAnalyzer/src/MEGeom.cc

0001 #include <cassert>
0002 #include <iostream>
0003 #include <cstdlib>
0004 #include <string>
0005 using namespace std;
0006
0007 #include "CalibCalorimetry/EcalLaserAnalyzer/interface/MEGeom.h"
0008 #include "CalibCalorimetry/EcalLaserAnalyzer/interface/ME.h"
0009
0010 #include <TCanvas.h>
0011 #include <TGraph.h>
0012 #include <TFile.h>
0013 #include <TH2.h>
0014
0015 //TString MEGeom::granularity[MEGeom::iSizeG] = {
0016 //  "Ecal Region", "Ecal Sector", "Laser Monitoring Region", "Laser Monitoring Module", "Super-Crystal", "Crystal"
0017 //};
0018
0019 const int MEGeom::_nbuf = 5;
0020 const int MEGeom::_nbinx = 2 * (MEGeom::_nbuf + 100 + 85) + 1;
0021 const float MEGeom::_xmin = -(0.5 + MEGeom::_nbuf + 100 + 85);
0022 const float MEGeom::_xmax = 0.5 + MEGeom::_nbuf + 100 + 85;
0023 const int MEGeom::_nbiny = 2 * (MEGeom::_nbuf + 180);
0024 const float MEGeom::_ymin = 180.5 - (MEGeom::_nbuf + 180);
0025 const float MEGeom::_ymax = 180.5 + (MEGeom::_nbuf + 180);
0026 const TH2* MEGeom::_h = new TH2F("globalEcal",
0027                                  "Global representation of ECAL",
0028                                  MEGeom::_nbinx,
0029                                  MEGeom::_xmin,
0030                                  MEGeom::_xmax,
0031                                  MEGeom::_nbiny,
0032                                  MEGeom::_ymin,
0033                                  MEGeom::_ymax);
0034
0035 //GHM ClassImp(MEGeom)
0036
0037 TH2* MEGeom::getHist(int ilmr, int unit) {
0038   int ireg;
0039   int ism;
0040   int idcc;
0041   int side;
0042   ME::regionAndSector(ilmr, ireg, ism, idcc, side);
0043
0044   std::cout << "LM region=" << ilmr << " reg/sm/dcc/side " << ireg << "/" << ism << "/" << idcc << "/" << side
0045             << std::endl;
0046
0047   TString hname = "LMR=";
0048   hname += ilmr;
0049   hname += " ";
0050   if (ireg == ME::iEEM || ireg == ME::iEEP) {
0051     hname += MEEEGeom::smName(ism);
0052   } else if (ireg == ME::iEBM || ireg == ME::iEBP) {
0053     hname += MEEBGeom::smName(ism);
0054   } else
0055     abort();
0056
0057   hname += " ";
0058   hname += "DCC=";
0059   hname += idcc;
0060   hname += "/";
0061   hname += side;
0062
0063   //
0064   // to produce these root files, run the runGeom executable
0065   //
0066   TH2* h_(nullptr);
0067   TFile* rootfile(nullptr);
0068   TString hn_;
0069   if (ireg == ME::iEBM || ireg == ME::iEBP) {
0070     rootfile = TFile::Open("ebgeom.root");
0071     assert(rootfile != nullptr);
0072     hn_ = "eb_loc";
0073     switch (unit) {
0074       case ME::iSector:
0075         break;
0076       case ME::iLMRegion:
0077         hn_ += "_side";
0078         break;
0079       case ME::iLMModule:
0080         hn_ += "_lmmod";
0081         break;
0082       case ME::iSuperCrystal:
0083         hn_ += "_tt";
0084         break;
0085       case ME::iCrystal:
0086         hn_ += "_cr";
0087         break;
0088       case ME::iElectronicChannel:
0089         hn_ += "_elecr";
0090         break;
0091       case ME::iHVChannel:
0092         hn_ += "_hv";
0093         break;
0094       case ME::iLVChannel:
0095         hn_ += "_lv";
0096         break;
0097     }
0098   } else {
0099     int isect = ism;
0100     if (ireg == ME::iEEM) {
0101       isect -= 9;
0102       rootfile = TFile::Open("eegeom_1.root");
0103     }
0104     if (ireg == ME::iEEP)
0105       rootfile = TFile::Open("eegeom_2.root");
0106     assert(rootfile != nullptr);
0107     hn_ = "eem_S";
0108     hn_ += isect;
0109     switch (unit) {
0110       case ME::iSector:
0111         break;
0112       case ME::iLMRegion:
0113         break;
0114       case ME::iLMModule:
0115         hn_ += "_lmmod";
0116         break;
0117       case ME::iSuperCrystal:
0118         hn_ += "_sc";
0119         break;
0120       case ME::iCrystal:
0121         hn_ += "_cr";
0122         break;
0123       case ME::iElectronicChannel:
0124         hn_ += "_cr";
0125         break;
0126       case ME::iHVChannel:
0127         break;
0128       case ME::iLVChannel:
0129         break;
0130     }
0131   }
0132   h_ = (TH2*)rootfile->Get(hn_);
0133   h_->SetTitle(hname);
0134   h_->GetXaxis()->SetTitle("ix");
0135   h_->GetXaxis()->CenterTitle();
0136   h_->GetYaxis()->SetTitle("iy");
0137   h_->GetYaxis()->CenterTitle();
0138   return h_;
0139 }
0140
0141 TGraph* MEGeom::getBoundary(int ilmr, int histtype) {
0142   // for local pictures, only sector or monitoring region
0143   if (histtype != ME::iSector && histtype != ME::iLMRegion)
0144     histtype = ME::iSector;
0145
0146   int ireg;
0147   int ism;
0148   int idcc;
0149   int side;
0150   ME::regionAndSector(ilmr, ireg, ism, idcc, side);
0151
0152   //
0153   // to produce these root files, run the runGeom executable
0154   //
0155   TGraph* g_(nullptr);
0156   TFile* rootfile(nullptr);
0157   TString gn_;
0158   if (ireg == ME::iEBM || ireg == ME::iEBP) {
0159     rootfile = TFile::Open("ebgeom.root");
0160     assert(rootfile != nullptr);
0161     switch (histtype) {
0162       case ME::iSector:
0163         gn_ = "SuperModule";
0164         break;
0165       case ME::iLMRegion:
0166         gn_ = "Side_";
0167         gn_ += side;
0168         break;
0169     }
0170   } else {
0171     int isect = ism;
0172     if (ireg == ME::iEEM) {
0173       isect -= 9;
0174       rootfile = TFile::Open("eegeom_1.root");
0175     }
0176     if (ireg == ME::iEEP)
0177       rootfile = TFile::Open("eegeom_2.root");
0178     assert(rootfile != nullptr);
0179     int lmr_ = ilmr;
0180     if (ireg == ME::iEEP)
0181       lmr_ -= 72;
0182     else if (ireg == ME::iEEM)
0183       lmr_ -= 82;
0184     switch (histtype) {
0185       case ME::iSector:
0186         gn_ = "Sector_";
0187         gn_ += isect;
0188         break;
0189       case ME::iLMRegion:
0190         gn_ = "LMRegion_";
0191         gn_ += lmr_;
0192         break;
0193     }
0194   }
0195   g_ = (TGraph*)rootfile->Get(gn_);
0196   return g_;
0197 }
0198
0199 void MEGeom::drawHist(int ilmr, int histtype, TCanvas* canv) {
0200   TH2* h = getHist(ilmr, histtype);
0201   assert(h != nullptr);
0202   TString tname = h->GetTitle();
0203   switch (histtype) {
0204     case ME::iSector:
0205       break;
0206     case ME::iLMRegion:
0207       tname += " Monitoring Regions";
0208       break;
0209     case ME::iLMModule:
0210       tname += " Monitoring Modules";
0211       break;
0212     case ME::iSuperCrystal:
0213       tname += " Super Crystals";
0214       break;
0215     case ME::iCrystal:
0216       tname += " Crystals";
0217       break;
0218     case ME::iElectronicChannel:
0219       tname += " Electronic Channels";
0220       break;
0221     case ME::iHVChannel:
0222       tname += " HV Channels";
0223       break;
0224     case ME::iLVChannel:
0225       tname += " LV Channels";
0226       break;
0227   }
0228
0229   if (canv == nullptr) {
0230     TString cname = tname;
0231     cname.ReplaceAll(" ", "_");
0232     canv = new TCanvas(cname, cname, 10, 10, 500, 500);
0233   }
0234   canv->SetTitle(tname);
0235   canv->cd();
0236
0237   h->Draw("COLZ");
0238
0239   TGraph* gsect = getBoundary(ilmr, ME::iSector);
0240   assert(gsect != nullptr);
0241   gsect->SetLineWidth(1);
0242   gsect->Draw("LSame");
0243   TGraph* gside = getBoundary(ilmr, ME::iLMRegion);
0244   assert(gside != nullptr);
0245   gside->SetLineWidth(2);
0246   gside->Draw("LSame");
0247 }
0248
0249 TH2* MEGeom::getGlobalHist(const char* name) {
0250   TH2* h = (TH2*)_h->Clone(name);
0251   h->Reset();
0252   return h;
0253 }
0254
0255 void MEGeom::setBinGlobalHist(TH2* h, int ix, int iy, int iz, float val) {
0256   // make sure it's a global hist
0257   TAxis* ax = h->GetXaxis();
0258   TAxis* ay = h->GetYaxis();
0259   assert(ax->GetNbins() == _nbinx && ay->GetNbins() == _nbiny);
0260   assert(ax->GetXmax() == _xmax && ax->GetXmin() == _xmin);
0261   assert(ay->GetXmax() == _ymax && ay->GetXmin() == _ymin);
0262
0263   int ibinx(0);
0264   int ibiny(0);
0265
0266   if (iz == 0) {
0267     int ieta = ix;
0268     int iphi = iy;
0269     assert(abs(ieta) >= 1 && abs(ieta) <= 85);
0270     assert(iphi >= 1 && iphi <= 360);
0271
0272     ibinx = ax->FindBin(ieta);
0273     ibiny = ay->FindBin(iphi);
0274   } else if (iz == -1) {
0275     assert(ix >= 1 && ix <= 100 && iy >= 1 && iy <= 100);
0276     ibinx = (_nbuf + 100 + 1) - ix;
0277     ibiny = (_nbuf + 180 + 50 + 1) - iy;
0278   } else if (iz == 1) {
0279     assert(ix >= 1 && ix <= 100 && iy >= 1 && iy <= 100);
0280     ibinx = _nbinx + 1 - ((_nbuf + 100 + 1) - ix);
0281     ibiny = (_nbuf + 180 + 50 + 1) - iy;
0282   }
0283   h->SetBinContent(ibinx, ibiny, val);
0284 }
0285
0286 // fixme !!!
0287 void MEGeom::drawGlobalBoundaries(int lineColor) {
0288   TGraph* gr(nullptr);
0289   for (int ism = 1; ism <= 36; ism++) {
0290     gr = MEEBGeom::getGraphBoundary(MEEBGeom::iSuperModule, ism, true);
0291     gr->SetLineWidth(2);
0292     gr->SetLineColor(lineColor);
0293     gr->Draw("LSame");
0294   }
0295   for (int isec = 1; isec <= 9; isec++) {
0296     gr = MEEEGeom::getGraphBoundary(MEEEGeom::iSector, isec);
0297     TGraph* grm = (TGraph*)gr->Clone();
0298     TGraph* grp = (TGraph*)gr->Clone();
0299     int n = gr->GetN();
0300     //      std::cout << std::endl;
0301     for (int ii = 0; ii < n; ii++) {
0302       double x_, y_;
0303       gr->GetPoint(ii, x_, y_);
0304       //      std::cout << isec << " " << ii << " x=" << x_ << " y=" << y_ << std::endl;
0305
0306       double xx_, yy_;
0307
0308       xx_ = -85 - x_;
0309       yy_ = 231 - y_;
0310       //      std::cout << isec << " " << 0 << " x=" << xx_ << " y=" << yy_ << std::endl;
0311
0312       grm->SetPoint(ii, xx_, yy_);
0313
0314       xx_ = 85 + x_;
0315       yy_ = 231 - y_;
0316       //      std::cout << isec << " " << 1 << " x=" << xx_ << " y=" << yy_ << std::endl;
0317       grp->SetPoint(ii, xx_, yy_);
0318     }
0319     grm->SetLineColor(lineColor);
0320     grm->SetLineWidth(2);
0321     grm->Draw("LSame");
0322     grp->SetLineColor(lineColor);
0323     grp->SetLineWidth(2);
0324     grp->Draw("LSame");
0325   }
0326 }