Back to home page

Project CMSSW displayed by LXR
 
 

    

File indexing completed on 2023-03-17 10:46:12 

0001 #ifndef HcalObjRepresent_h
0002 #define HcalObjRepresent_h
0003 
0004 #include "CondFormats/HcalObjects/interface/HcalGains.h"
0005 
0006 //#include "CondCore/Utilities/interface/PayLoadInspector.h"
0007 //#include "CondCore/Utilities/interface/InspectorPythonWrapper.h"
0008 #include "CondFormats/HcalObjects/interface/HcalDetIdRelationship.h"
0009 #include "CondFormats/HcalObjects/interface/HcalCondObjectContainer.h"
0010 #include "DataFormats/HcalDetId/interface/HcalDetId.h"
0011 #include "DataFormats/HcalDetId/interface/HcalOtherDetId.h"
0012 #include "DataFormats/HcalDetId/interface/HcalCastorDetId.h"
0013 #include "DataFormats/HcalDetId/interface/HcalZDCDetId.h"
0014 #include "FWCore/Utilities/interface/Exception.h"
0015 
0016 #include <string>
0017 #include <fstream>
0018 #include <sstream>
0019 #include <vector>
0020 
0021 #include "TH1F.h"
0022 #include "TH2F.h"
0023 #include "DataFormats/DetId/interface/DetId.h"
0024 
0025 #include "TROOT.h"
0026 #include "TCanvas.h"
0027 #include "TStyle.h"
0028 #include "TColor.h"
0029 #include "TLine.h"
0030 #include "TLatex.h"
0031 #include "TProfile.h"
0032 #include "TPaveLabel.h"
0033 
0034 //functions for correct representation of data in summary and plot
0035 namespace HcalObjRepresent {
0036 
0037   //used to produce all display objects for payload inspector
0038   template <class Items, class Item>
0039   class HcalDataContainer {
0040   public:
0041     HcalDataContainer(std::shared_ptr<Items> payload, unsigned int run) : payload_(payload), run_(run) {
0042       PlotMode_ = "Map";
0043     }
0044 
0045     virtual ~HcalDataContainer(){};
0046     // For easier channel mapping
0047     typedef std::tuple<int, int, int> Coord;
0048     typedef std::map<Coord, Item> tHcalValCont;
0049     // mapping of pair of subdetector name (e.g, "HE") and depth number (e.g. 3) to Histogram of data for that subdetector/depth pair
0050     typedef std::map<std::pair<std::string, int>, TH2F*> DepthMap;
0051 
0052     ///////////////// public Get functions  /////////////////
0053     unsigned int GetRun() { return run_; }
0054 
0055     std::string GetTopoMode() { return TopoMode_; }
0056 
0057     std::map<std::string, int> GetSubDetDepths() { return subDetDepths_; }
0058 
0059     DepthMap GetDepths() {
0060       fillValConts();
0061       //std::cout << "Got depths with run number = " << std::to_string(GetRun()) << std::endl;
0062       return depths_;
0063     }
0064 
0065     ///////////////// setting object fields /////////////////
0066 
0067     // Fills a tHcalValCont for each subdetector, setting Topology Mode along the way
0068     void fillValConts() {
0069       if (!depths_.empty())
0070         return;
0071       int iphi, ieta, depth;
0072       HcalDetId detId;
0073 
0074       std::string subDetName;
0075       std::vector<Item> itemsVec;
0076       std::pair<std::string, int> depthKey;
0077       const char* histLabel;
0078       for (std::pair<std::string, std::vector<Item> > cont : (*payload_).getAllContainers()) {
0079         subDetName = std::get<0>(cont);
0080         itemsVec = std::get<1>(cont);
0081 
0082         auto valContainer = getContFromString(subDetName);
0083 
0084         for (Item item : itemsVec) {
0085           detId = HcalDetId(item.rawId());
0086           iphi = detId.iphi();
0087           ieta = detId.ieta();
0088           depth = detId.depth();
0089           Coord coord = std::make_tuple(depth, ieta, iphi);
0090           //Add hist if it's not there, AND if it is one we care about: HO,HB,HE,HF, AND it's not an empty entry (depth = 0 in HO)
0091           if (subDetName[0] == 'H' && depth != 0) {
0092             valContainer->insert(std::pair<std::tuple<int, int, int>, Item>(coord, item));
0093 
0094             depthKey = std::make_pair(subDetName, depth);
0095 
0096             auto depthVal = depths_.find(depthKey);
0097             if (depthVal == depths_.end()) {
0098               histLabel = ("run" + std::to_string(run_) + "_" + subDetName + "_d" + std::to_string(depth)).c_str();
0099               depths_.insert(std::make_pair(std::make_pair(subDetName, depth),
0100                                             new TH2F(histLabel, histLabel, 83, -42.5, 41.5, 71, 0.5, 71.5)));
0101             }
0102             depths_[depthKey]->Fill(ieta, iphi, getValue(&item));
0103           }
0104         }
0105       }
0106 
0107       //Still need to know which hists to take when done; decide now
0108       setTopoModeFromValConts();
0109     }
0110 
0111     ////NOTE to be implemented in PayloadInspector classes
0112     virtual float getValue(Item* item) {
0113       throw cms::Exception("Value definition not found") << "getValue definition not found for " << payload_->myname();
0114     };
0115 
0116     //Gets Hcal Object at given coordinate
0117     //Currently unused but remains as a potentially useful function
0118     Item* getItemFromValCont(std::string subDetName, int depth, int ieta, int iphi, bool throwOnFail) {
0119       Item* cell = nullptr;
0120       Coord coord = std::make_tuple(depth, ieta, iphi);
0121       tHcalValCont* valContainer = getContFromString(subDetName);
0122 
0123       auto it = valContainer->find(coord);
0124       if (it != valContainer->end())
0125         cell = &it->second;
0126 
0127       if ((!cell)) {
0128         if (throwOnFail) {
0129           throw cms::Exception("Conditions not found")
0130               << "Unavailable Conditions of type " << payload_->myname() << " for cell " << subDetName << " (" << depth
0131               << ", " << ieta << ", " << iphi << ")";
0132         }
0133       }
0134       return cell;
0135     }
0136 
0137     /////////////////// Building Graphics //////////////////////
0138 
0139     //TODO: remove zero entries from doing divide and subtract
0140 
0141     // To generate Ratios of two IOVs
0142     void Divide(HcalDataContainer* dataCont2) {
0143       //TODO: Do something like looping over this and that depths setting every empty bin (which I think means content=0) to -999. Then replacing the divide call with another manual loop that sets all bins with content -999 to 0
0144 
0145       PlotMode_ = "Ratio";
0146       DepthMap::iterator depth1;
0147       std::pair<std::string, int> key;
0148       DepthMap myDepths = this->GetDepths();
0149       DepthMap theirDepths = dataCont2->GetDepths();
0150       for (depth1 = myDepths.begin(); depth1 != myDepths.end(); depth1++) {
0151         key = depth1->first;
0152         if (theirDepths.count(key) != 0) {
0153           myDepths.at(key)->Divide((const TH1*)theirDepths.at(key));
0154         } else {
0155           throw cms::Exception("Unaligned Conditions")
0156               << "trying to plot conditions for " << payload_->myname() << "; found value for " << std::get<0>(key)
0157               << " depth " << std::to_string(std::get<1>(key)) << " in run " << GetRun() << " but not in run "
0158               << dataCont2->GetRun();
0159         }
0160       }
0161     }
0162 
0163     // To generate Diffs of two IOVs
0164     void Subtract(HcalDataContainer* dataCont2) {
0165       PlotMode_ = "Diff";
0166       DepthMap::iterator depth1;
0167       std::pair<std::string, int> key;
0168       DepthMap myDepths = this->GetDepths();
0169       DepthMap theirDepths = dataCont2->GetDepths();
0170       for (depth1 = myDepths.begin(); depth1 != myDepths.end(); depth1++) {
0171         key = depth1->first;
0172         if (theirDepths.count(key) != 0) {
0173           myDepths.at(key)->Add(myDepths.at(key), theirDepths.at(key), 1, -1);
0174         } else {
0175           throw cms::Exception("Unaligned Conditions")
0176               << "trying to plot conditions for " << payload_->myname() << "; found value for " << std::get<0>(key)
0177               << " depth " << std::to_string(std::get<1>(key)) << " in run " << GetRun() << " but not in run "
0178               << dataCont2->GetRun();
0179         }
0180       }
0181     }
0182 
0183     // wisely determines what range to set histogram axis to
0184     std::pair<float, float> GetRange(TH1* hist) {
0185       if (PlotMode_ == "Ratio") {
0186         float amp;
0187         Double_t adjustMin = 1;
0188         Double_t tempMin;
0189         int nBinsX = hist->GetXaxis()->GetNbins();
0190         int nBinsY = hist->GetYaxis()->GetNbins();
0191         for (int i = 0; i < nBinsX; i++) {
0192           for (int j = 0; j < nBinsY; j++) {
0193             tempMin = hist->GetBinContent(i, j);
0194             if ((tempMin != 0) && (tempMin < adjustMin))
0195               adjustMin = tempMin;
0196           }
0197         }
0198         amp = std::max((1 - adjustMin), (hist->GetMaximum() - 1));
0199         //amp = std::max((1 - hist->GetMinimum()),(hist->GetMaximum() - 1) );
0200         return std::make_pair(1 - amp, 1 + amp);
0201       } else if (PlotMode_ == "Diff") {
0202         float amp;
0203         amp = std::max((0 - hist->GetMinimum()), hist->GetMaximum());
0204         return std::make_pair((-1 * amp), amp);
0205       } else {
0206         Double_t adjustMin = 10000;
0207         Double_t tempMin;
0208         int nBinsX = hist->GetXaxis()->GetNbins();
0209         int nBinsY = hist->GetYaxis()->GetNbins();
0210         for (int i = 0; i < nBinsX; i++) {
0211           for (int j = 0; j < nBinsY; j++) {
0212             tempMin = hist->GetBinContent(i, j);
0213             if ((tempMin != 0) && (tempMin < adjustMin))
0214               adjustMin = tempMin;
0215           }
0216         }
0217         return std::make_pair(((adjustMin == 10000) ? hist->GetMinimum() : adjustMin), hist->GetMaximum());
0218       }
0219     }
0220 
0221     // set style
0222     void initGraphics() {
0223       gStyle->SetOptStat(0);
0224       gStyle->SetPalette(1);
0225       gStyle->SetOptFit(0);
0226       gStyle->SetLabelFont(42);
0227       gStyle->SetLabelFont(42);
0228       gStyle->SetTitleFont(42);
0229       gStyle->SetTitleFont(42);
0230       gStyle->SetMarkerSize(0);
0231       gStyle->SetTitleOffset(1.3, "Y");
0232       gStyle->SetTitleOffset(1.0, "X");
0233       gStyle->SetNdivisions(510);
0234       gStyle->SetStatH(0.11);
0235       gStyle->SetStatW(0.33);
0236       gStyle->SetTitleW(0.4);
0237       gStyle->SetTitleX(0.13);
0238       gStyle->SetPadTickX(1);
0239       gStyle->SetPadTickY(1);
0240     }
0241 
0242     TH1D* GetProjection(TH2F* hist, std::string plotType, const char* newName, std::string subDetName, int depth) {
0243       //TODO: Also want average for standard projection of 2DHist (not ratio or diff)?
0244       //if (PlotMode_ != "Ratio") return (plotType=="EtaProfile") ? ((TH2F*)(hist->Clone("temp")))->ProjectionX(newName) : ((TH2F*)(hist->Clone("temp")))->ProjectionY(newName);
0245 
0246       //TH1D* projection = ((TH2F*)(depths_[std::make_pair(subDetName,depth)]->Clone("temp")))->ProjectionX(newName);
0247 
0248       int xBins = (plotType == "EtaProfile") ? 83 : 71;
0249       int etaMin = -42, etaMax = 42, phiMin = 1, phiMax = 72;
0250       int xMin = (plotType == "EtaProfile") ? etaMin : phiMin;
0251       int xMax = (plotType == "EtaProfile") ? etaMax : phiMax;
0252       int otherMin = (plotType == "EtaProfile") ? phiMin : etaMin;
0253       int otherMax = (plotType == "EtaProfile") ? phiMax : etaMax;
0254       TH1D* retHist = new TH1D(newName, newName, xBins, xMin, xMax);
0255       int numChannels;
0256       Double_t sumVal;
0257       Double_t channelVal;
0258       int ieta, iphi;
0259       int bin = 0;
0260       for (int i = xMin; i <= xMax; i++ && bin++) {
0261         numChannels = 0;
0262         sumVal = 0;
0263         for (int j = otherMin; j <= otherMax; j++) {
0264           ieta = (plotType == "EtaProfile") ? i : j;
0265           ieta += 42;
0266           iphi = (plotType == "EtaProfile") ? j : i;
0267           iphi += -1;
0268           channelVal = hist->GetBinContent(ieta, iphi);
0269           //std::cout << "(ieta, iphi)    :       (" << std::to_string(ieta) << ", " << std::to_string(iphi) << ")" << std::endl;
0270           if (channelVal != 0) {
0271             sumVal += channelVal;
0272             numChannels++;
0273           }
0274         }
0275         //if(sumVal !=0) projection->SetBinContent(i, sumVal/((Double_t)numChannels));//retHist->Fill(i,sumVal/((Double_t)numChannels));
0276         if (sumVal != 0)
0277           retHist->Fill(i, sumVal / ((Double_t)numChannels));
0278       }
0279       return retHist;
0280       //return projection;
0281     }
0282 
0283     // fills a canvas with given subdetector information, plotting all depths
0284     void FillCanv(TCanvas* canvas,
0285                   std::string subDetName,
0286                   int startDepth = 1,
0287                   int startCanv = 1,
0288                   std::string plotForm = "2DHist") {
0289       const char* newName;
0290       std::pair<float, float> range;
0291       int padNum;
0292       int maxDepth = (subDetName == "HO") ? 4 : subDetDepths_[subDetName];
0293       TH1D* projection;
0294       TLatex label;
0295       for (int i = startDepth; i <= maxDepth; i++) {
0296         //skip if data not obtained; TODO: Maybe add text on plot saying data not found?
0297         if (depths_.count(std::make_pair(subDetName, i)) == 0) {
0298           return;
0299         }
0300 
0301         padNum = i + startCanv - 1;
0302         if (subDetName == "HO")
0303           padNum = padNum - 3;
0304         canvas->cd(padNum);
0305         canvas->GetPad(padNum)->SetGridx(1);
0306         canvas->GetPad(padNum)->SetGridy(1);
0307 
0308         if (plotForm == "2DHist") {
0309           canvas->GetPad(padNum)->SetRightMargin(0.13);
0310           range = GetRange(depths_[std::make_pair(subDetName, i)]);
0311           depths_[std::make_pair(subDetName, i)]->Draw("colz");
0312           depths_[std::make_pair(subDetName, i)]->SetContour(100);
0313           depths_[std::make_pair(subDetName, i)]->GetXaxis()->SetTitle("ieta");
0314           depths_[std::make_pair(subDetName, i)]->GetYaxis()->SetTitle("iphi");
0315           depths_[std::make_pair(subDetName, i)]->GetXaxis()->CenterTitle();
0316           depths_[std::make_pair(subDetName, i)]->GetYaxis()->CenterTitle();
0317           depths_[std::make_pair(subDetName, i)]->GetZaxis()->SetRangeUser(std::get<0>(range), std::get<1>(range));
0318           depths_[std::make_pair(subDetName, i)]->GetYaxis()->SetTitleSize(0.06);
0319           depths_[std::make_pair(subDetName, i)]->GetYaxis()->SetTitleOffset(0.80);
0320           depths_[std::make_pair(subDetName, i)]->GetXaxis()->SetTitleSize(0.06);
0321           depths_[std::make_pair(subDetName, i)]->GetXaxis()->SetTitleOffset(0.80);
0322           depths_[std::make_pair(subDetName, i)]->GetYaxis()->SetLabelSize(0.055);
0323           depths_[std::make_pair(subDetName, i)]->GetXaxis()->SetLabelSize(0.055);
0324         } else {
0325           canvas->GetPad(padNum)->SetLeftMargin(0.152);
0326           canvas->GetPad(padNum)->SetRightMargin(0.02);
0327           //gStyle->SetTitleOffset(1.6,"Y");
0328           newName = ("run_" + std::to_string(run_) + "_" + subDetName + "_d" + std::to_string(i) + "_" +
0329                      (plotForm == "EtaProfile" ? "ieta" : "iphi"))
0330                         .c_str();
0331           //projection = ((TH2F*)(depths_[std::make_pair(subDetName,i)]->Clone("temp")))->ProjectionX(newName);
0332           projection = GetProjection(depths_[std::make_pair(subDetName, i)], plotForm, newName, subDetName, i);
0333           range = GetRange(projection);
0334           projection->Draw("hist");
0335           projection->GetXaxis()->SetTitle((plotForm == "EtaProfile" ? "ieta" : "iphi"));
0336           projection->GetXaxis()->CenterTitle();
0337           projection->GetYaxis()->SetTitle(
0338               (payload_->myname() + " " + (PlotMode_ == "Map" ? "" : PlotMode_) + " " + GetUnit(payload_->myname()))
0339                   .c_str());
0340           label.SetNDC();
0341           label.SetTextAlign(26);
0342           label.SetTextSize(0.05);
0343           label.DrawLatex(
0344               0.3, 0.95, ("run_" + std::to_string(run_) + "_" + subDetName + "_d" + std::to_string(i)).c_str());
0345           projection->GetYaxis()->CenterTitle();
0346           projection->GetXaxis()->SetTitleSize(0.06);
0347           projection->GetYaxis()->SetTitleSize(0.06);
0348           projection->GetXaxis()->SetTitleOffset(0.80);
0349           projection->GetXaxis()->SetLabelSize(0.055);
0350           projection->GetYaxis()->SetTitleOffset(1.34);
0351           projection->GetYaxis()->SetLabelSize(0.055);
0352         }
0353       }
0354     }
0355 
0356     //// functions called in Payload Inspector classes to import final canvases to be plotted
0357 
0358     // profile = "EtaProfile" || "PhiProfile"
0359     TCanvas* getCanvasAll(std::string profile = "2DHist") {
0360       fillValConts();
0361       initGraphics();
0362       TCanvas* HAll = new TCanvas("HAll", "HAll", 1680, (GetTopoMode() == "2015/2016") ? 1680 : 2500);
0363       HAll->Divide(3, (GetTopoMode() == "2015/2016") ? 3 : 6, 0.02, 0.01);
0364       FillCanv(HAll, "HB", 1, 1, profile);
0365       FillCanv(HAll, "HO", 4, 3, profile);
0366       FillCanv(HAll, "HF", 1, 4, profile);
0367       FillCanv(HAll, "HE", 1, (GetTopoMode() == "2015/2016") ? 7 : 10, profile);
0368       return HAll;
0369     }
0370 
0371     TCanvas* getCanvasHF() {
0372       fillValConts();
0373       initGraphics();
0374       TCanvas* HF = new TCanvas("HF", "HF", 1600, 1000);
0375       HF->Divide(3, 2, 0.02, 0.01);
0376       FillCanv(HF, "HF");
0377       return HF;
0378     }
0379     TCanvas* getCanvasHE() {
0380       fillValConts();
0381       initGraphics();
0382       TCanvas* HE = new TCanvas("HE", "HE", 1680, 1680);
0383       HE->Divide(3, 3, 0.02, 0.01);
0384       FillCanv(HE, "HE");
0385       return HE;
0386     }
0387     TCanvas* getCanvasHBHO() {
0388       fillValConts();
0389       initGraphics();
0390       TCanvas* HBHO = new TCanvas("HBHO", "HBHO", 1680, 1680);
0391       HBHO->Divide(3, 3, 0.02, 0.01);
0392       FillCanv(HBHO, "HB");
0393       FillCanv(HBHO, "HO", 4, 3);
0394       FillCanv(HBHO, "HB", 1, 4, "EtaProfile");
0395       FillCanv(HBHO, "HO", 4, 6, "EtaProfile");
0396       FillCanv(HBHO, "HB", 1, 7, "PhiProfile");
0397       FillCanv(HBHO, "HO", 4, 9, "PhiProfile");
0398       return HBHO;
0399     }
0400 
0401     std::string GetUnit(std::string type) {
0402       std::string unit = units_[type];
0403       if (unit.empty())
0404         return "";
0405       else
0406         return "(" + unit + ")";
0407     }
0408 
0409   private:
0410     DepthMap depths_;
0411     std::shared_ptr<Items> payload_;
0412     unsigned int run_;
0413     std::string TopoMode_;
0414     // "Map", "Ratio", or "Diff"
0415     std::string PlotMode_;
0416 
0417     tHcalValCont HBvalContainer;
0418     tHcalValCont HEvalContainer;
0419     tHcalValCont HOvalContainer;
0420     tHcalValCont HFvalContainer;
0421     tHcalValCont HTvalContainer;
0422     tHcalValCont ZDCvalContainer;
0423     tHcalValCont CALIBvalContainer;
0424     tHcalValCont CASTORvalContainer;
0425     std::map<std::string, int> subDetDepths_;
0426     std::map<std::string, std::string> units_ = {{"HcalPedestals", "ADC"},
0427                                                  {"HcalGains", ""},             //dimensionless TODO: verify
0428                                                  {"HcalL1TriggerObjects", ""},  //dimensionless TODO: Verify
0429                                                  {"HcalPedestalWidths", "ADC"},
0430                                                  {"HcalRespCorrs", ""},  //dimensionless TODO: verify
0431                                                  {"Dark Current", ""},
0432                                                  {"fcByPE", ""},
0433                                                  {"crossTalk", ""},
0434                                                  {"parLin", ""}};
0435 
0436     tHcalValCont* getContFromString(std::string subDetString) {
0437       if (subDetString == "HB")
0438         return &HBvalContainer;
0439       else if (subDetString == "HE")
0440         return &HEvalContainer;
0441       else if (subDetString == "HF")
0442         return &HFvalContainer;
0443       else if (subDetString == "HO")
0444         return &HOvalContainer;
0445       else if (subDetString == "HT")
0446         return &HTvalContainer;
0447       else if (subDetString == "CALIB")
0448         return &CALIBvalContainer;
0449       else if (subDetString == "CASTOR")
0450         return &CASTORvalContainer;
0451       //else return &ZDCvalContainer;
0452       else if (subDetString == "ZDC_EM" || subDetString == "ZDC" || subDetString == "ZDC_HAD" ||
0453                subDetString == "ZDC_LUM")
0454         return &ZDCvalContainer;
0455       else
0456         throw cms::Exception("subDetString " + subDetString + " not found in Item");
0457     }
0458 
0459     void setTopoModeFromValConts(bool throwOnFail = false) {
0460       // Check HEP17 alternate channel for 2017, just by checking if the 7th depth is there, or if HF has 4 depths
0461       if (depths_.count(std::make_pair("HF", 4)) != 0 || depths_.count(std::make_pair("HE", 7)) != 0)
0462         TopoMode_ = "2017";
0463       // Check endcap depth unique to 2018
0464       else if (HEvalContainer.count(std::make_tuple(7, -26, 63)) != 0)
0465         TopoMode_ = "2018";
0466       // if not 2017 or 2018, 2015 and 2016 are the same
0467       else
0468         TopoMode_ = "2015/2016";
0469       //NOTE: HO's one depth is labeled depth 4
0470       if (TopoMode_ == "2018" || TopoMode_ == "2017") {
0471         subDetDepths_.insert(std::pair<std::string, int>("HB", 2));
0472         subDetDepths_.insert(std::pair<std::string, int>("HE", 7));
0473         subDetDepths_.insert(std::pair<std::string, int>("HF", 4));
0474         subDetDepths_.insert(std::pair<std::string, int>("HO", 1));
0475       } else if (TopoMode_ == "2015/2016") {
0476         subDetDepths_.insert(std::pair<std::string, int>("HB", 2));
0477         subDetDepths_.insert(std::pair<std::string, int>("HE", 3));
0478         subDetDepths_.insert(std::pair<std::string, int>("HF", 2));
0479         subDetDepths_.insert(std::pair<std::string, int>("HO", 1));
0480       }
0481     }
0482   };
0483 
0484   inline void drawTable(int nbRows, int nbColumns) {
0485     TLine* l = new TLine;
0486     l->SetLineWidth(1);
0487     for (int i = 1; i < nbRows; i++) {
0488       double y = (double)i;
0489       l = new TLine(0., y, nbColumns, y);
0490       l->Draw();
0491     }
0492 
0493     for (int i = 1; i < nbColumns; i++) {
0494       double x = (double)i;
0495       double y = (double)nbRows;
0496       l = new TLine(x, 0., x, y);
0497       l->Draw();
0498     }
0499   }
0500 
0501   inline std::string SciNotatStr(float num) {
0502     // Create an output string stream
0503     std::ostringstream streamObj2;
0504 
0505     if (num == -1)
0506       return "NOT FOUND";
0507 
0508     //Add double to stream
0509     streamObj2 << num;
0510     // Get string from output string stream
0511     std::string strObj2 = streamObj2.str();
0512 
0513     return strObj2;
0514   }
0515 
0516   inline std::string IntToBinary(unsigned int number) {
0517     std::stringstream ss;
0518     unsigned int mask = 1 << 31;
0519     for (unsigned short int i = 0; i < 32; ++i) {
0520       //if (!(i % 4))
0521       //    ss << "_";
0522       if (mask & number)
0523         ss << "1";
0524       else
0525         ss << "0";
0526       mask = mask >> 1;
0527     }
0528     return ss.str();
0529   }
0530 
0531   inline const bool isBitSet(unsigned int bitnumber, unsigned int status) {
0532     unsigned int statadd = 0x1 << (bitnumber);
0533     return (status & statadd) ? (true) : (false);
0534   }
0535 
0536   inline std::string getBitsSummary(uint32_t bits, std::string statusBitArray[], short unsigned int bitMap[]) {
0537     std::stringstream ss;
0538     for (unsigned int i = 0; i < 9; ++i) {
0539       if (isBitSet(bitMap[i], bits)) {
0540         ss << "[" << bitMap[i] << "]" << statusBitArray[bitMap[i]] << "; ";
0541       }
0542     }
0543     ss << std::endl;
0544     return ss.str();
0545   }
0546 
0547   //functions for making plot:
0548   inline void setBinLabels(std::vector<TH2F>& depth) {
0549     // Set labels for all depth histograms
0550     for (unsigned int i = 0; i < depth.size(); ++i) {
0551       depth[i].SetXTitle("i#eta");
0552       depth[i].SetYTitle("i#phi");
0553     }
0554 
0555     std::stringstream label;
0556 
0557     // set label on every other bin
0558     for (int i = -41; i <= -29; i = i + 2) {
0559       label << i;
0560       depth[0].GetXaxis()->SetBinLabel(i + 42, label.str().c_str());
0561       depth[1].GetXaxis()->SetBinLabel(i + 42, label.str().c_str());
0562       label.str("");
0563     }
0564     depth[0].GetXaxis()->SetBinLabel(14, "-29HE");
0565     depth[1].GetXaxis()->SetBinLabel(14, "-29HE");
0566 
0567     // offset by one for HE
0568     for (int i = -27; i <= 27; i = i + 2) {
0569       label << i;
0570       depth[0].GetXaxis()->SetBinLabel(i + 43, label.str().c_str());
0571       label.str("");
0572     }
0573     depth[0].GetXaxis()->SetBinLabel(72, "29HE");
0574     for (int i = 29; i <= 41; i = i + 2) {
0575       label << i;
0576       depth[0].GetXaxis()->SetBinLabel(i + 44, label.str().c_str());
0577       label.str("");
0578     }
0579     for (int i = 16; i <= 28; i = i + 2) {
0580       label << i - 43;
0581       depth[1].GetXaxis()->SetBinLabel(i, label.str().c_str());
0582       label.str("");
0583     }
0584     depth[1].GetXaxis()->SetBinLabel(29, "NULL");
0585     for (int i = 15; i <= 27; i = i + 2) {
0586       label << i;
0587       depth[1].GetXaxis()->SetBinLabel(i + 15, label.str().c_str());
0588       label.str("");
0589     }
0590 
0591     depth[1].GetXaxis()->SetBinLabel(44, "29HE");
0592     for (int i = 29; i <= 41; i = i + 2) {
0593       label << i;
0594       depth[1].GetXaxis()->SetBinLabel(i + 16, label.str().c_str());
0595       label.str("");
0596     }
0597 
0598     // HE depth 3 labels;
0599     depth[2].GetXaxis()->SetBinLabel(1, "-28");
0600     depth[2].GetXaxis()->SetBinLabel(2, "-27");
0601     depth[2].GetXaxis()->SetBinLabel(3, "Null");
0602     depth[2].GetXaxis()->SetBinLabel(4, "-16");
0603     depth[2].GetXaxis()->SetBinLabel(5, "Null");
0604     depth[2].GetXaxis()->SetBinLabel(6, "16");
0605     depth[2].GetXaxis()->SetBinLabel(7, "Null");
0606     depth[2].GetXaxis()->SetBinLabel(8, "27");
0607     depth[2].GetXaxis()->SetBinLabel(9, "28");
0608   }
0609   // Now define functions that can be used in conjunction with EtaPhi histograms
0610 
0611   // This arrays the eta binning for depth 2 histograms (with a gap between -15 -> +15)
0612   const int binmapd2[] = {-42, -41, -40, -39, -38, -37, -36, -35, -34, -33, -32, -31, -30, -29,   -28,
0613                           -27, -26, -25, -24, -23, -22, -21, -20, -19, -18, -17, -16, -15, -9999, 15,
0614                           16,  17,  18,  19,  20,  21,  22,  23,  24,  25,  26,  27,  28,  29,    30,
0615                           31,  32,  33,  34,  35,  36,  37,  38,  39,  40,  41,  42};
0616 
0617   // This stores eta binning in depth 3 (where HE is only present at a few ieta values)
0618 
0619   const int binmapd3[] = {-28, -27, -9999, -16, -9999, 16, -9999, 27, 28};
0620 
0621   inline int CalcEtaBin(int subdet, int ieta, int depth) {
0622     // This takes the eta value from a subdetector and return an eta counter value as used by eta-phi array
0623     // (ieta=-41 corresponds to bin 0, +41 to bin 85 -- there are two offsets to deal with the overlap at |ieta|=29).
0624     // For HO, ieta = -15 corresponds to bin 0, and ieta=15 is bin 30
0625     // For HE depth 3, things are more complicated, but feeding the ieta value will give back the corresponding counter eta value
0626 
0627     // The CalcEtaBin value is the value as used within our array counters, and thus starts at 0.
0628     // If you are using it with getBinContent or setBinContent, you will need to add +1 to the result of this function
0629 
0630     int etabin = -9999;  // default invalid value
0631 
0632     if (depth == 1) {
0633       // Depth 1 is fairly straightforward -- just shift HF-, HF+ by -/+1
0634       etabin = ieta + 42;
0635       if (subdet == HcalForward) {
0636         ieta < 0 ? etabin-- : etabin++;
0637       }
0638     }
0639 
0640     else if (depth == 2) {
0641       // Depth 2 is more complicated, given that there are no cells in the range |ieta|<15
0642       if (ieta < -14) {
0643         etabin = ieta + 42;
0644         if (subdet == HcalForward)
0645           etabin--;
0646       } else if (ieta > 14) {
0647         etabin = ieta + 14;
0648         if (subdet == HcalForward)
0649           etabin++;
0650       }
0651 
0652     }
0653     // HO is also straightforward; a simple offset to the ieta value is applied
0654     else if (subdet == HcalOuter && abs(ieta) < 16)
0655       etabin = ieta + 15;
0656     else if (subdet == HcalEndcap) {
0657       // HE depth 3 has spotty coverage; hard-code the bin response
0658       if (depth == 3) {
0659         if (ieta == -28)
0660           etabin = 0;
0661         else if (ieta == -27)
0662           etabin = 1;
0663         else if (ieta == -16)
0664           etabin = 3;
0665         else if (ieta == 16)
0666           etabin = 5;
0667         else if (ieta == 27)
0668           etabin = 7;
0669         else if (ieta == 28)
0670           etabin = 8;
0671       }
0672     }
0673     return etabin;
0674   }
0675 
0676   inline int CalcIeta(int subdet, int eta, int depth) {
0677     // This function returns the 'true' ieta value given subdet, eta, and depth
0678     // Here 'eta' is the index from our arrays (it starts at 0);
0679     // remember that histogram bins start with bin 1, so there's an offset of 1
0680     // to consider if using getBinContent(eta,phi)
0681 
0682     // eta runs from 0...X  (X depends on depth)
0683     int ieta = -9999;  // default value is nonsensical
0684     if (subdet == HcalBarrel) {
0685       if (depth == 1) {
0686         ieta = eta - 42;
0687         if (ieta == 0)
0688           return -9999;
0689         return ieta;
0690       } else if (depth == 2) {
0691         ieta = binmapd2[eta];
0692         if (ieta == 0)
0693           return -9999;
0694         if (ieta == 17 || ieta == -17)
0695           return -9999;  // no depth 2 cells at |ieta| = 17
0696         return ieta;
0697       } else
0698         return -9999;  // non-physical value
0699     } else if (subdet == HcalForward) {
0700       if (depth == 1) {
0701         ieta = eta - 42;
0702         if (eta < 13)
0703           ieta++;
0704         else if (eta > 71)
0705           ieta--;
0706         else
0707           return -9999;  // if outside forward range, return dummy
0708         return ieta;
0709       } else if (depth == 2) {
0710         ieta = binmapd2[eta];  // special map for depth 2
0711         if (ieta <= -30)
0712           ieta++;
0713         else if (ieta >= 30)
0714           ieta--;
0715         else
0716           return -9999;
0717         return ieta;
0718       } else
0719         return -9999;
0720     }
0721 
0722     else if (subdet == HcalEndcap) {
0723       if (depth == 1)
0724         ieta = eta - 42;
0725       else if (depth == 2) {
0726         ieta = binmapd2[eta];
0727         if (abs(ieta) > 29 || abs(ieta) < 18)
0728           return -9999;  // outside HE
0729         if (ieta == 0)
0730           return -9999;
0731         return ieta;
0732       } else if (depth == 3) {
0733         if (eta < 0 || eta > 8)
0734           return -9999;
0735         else
0736           ieta = binmapd3[eta];  // special map for depth 3
0737         if (ieta == 0)
0738           return -9999;
0739         return ieta;
0740       } else
0741         return -9999;
0742     }  // HcalEndcap
0743     else if (subdet == HcalOuter) {
0744       if (depth != 4)
0745         return -9999;
0746       else {
0747         ieta = eta - 15;  // bin 0 is ieta=-15, all bins increment normally from there
0748         if (abs(ieta) > 15)
0749           return -9999;
0750         if (ieta == 0)
0751           return -9999;
0752         return ieta;
0753       }
0754     }  // HcalOuter
0755     if (ieta == 0)
0756       return -9999;
0757     return ieta;
0758   }
0759 
0760   inline int CalcIeta(int eta, int depth) {
0761     // This version of CalcIeta does the same as the function above,
0762     // but does not require that 'subdet' be specified.
0763 
0764     // returns ieta value give an eta counter.
0765     // eta runs from 0...X  (X depends on depth)
0766     int ieta = -9999;
0767     if (eta < 0)
0768       return ieta;
0769     if (depth == 1) {
0770       ieta =
0771           eta -
0772           42;  // default shift: bin 0 corresponds to a histogram ieta of -42 (which is offset by 1 from true HF value of -41)
0773       if (eta < 13)
0774         ieta++;
0775       else if (eta > 71)
0776         ieta--;
0777       if (ieta == 0)
0778         ieta = -9999;
0779       return ieta;
0780     } else if (depth == 2) {
0781       if (eta > 57)
0782         return -9999;
0783       else {
0784         ieta = binmapd2[eta];
0785         if (ieta == -9999)
0786           return ieta;
0787         if (ieta == 0)
0788           return -9999;
0789         if (ieta == 17 || ieta == -17)
0790           return -9999;  // no depth 2 cells at |ieta| = 17
0791         else if (ieta <= -30)
0792           ieta++;
0793         else if (ieta >= 30)
0794           ieta--;
0795         return ieta;
0796       }
0797     } else if (depth == 3) {
0798       if (eta > 8)
0799         return -9999;
0800       else
0801         ieta = binmapd3[eta];
0802       if (ieta == 0)
0803         return -9999;
0804       return ieta;
0805     } else if (depth == 4) {
0806       ieta = eta - 15;  // bin 0 is ieta=-15, all bins increment normally from there
0807       if (abs(ieta) > 15)
0808         return -9999;
0809       if (ieta == 0)
0810         return -9999;
0811       return ieta;
0812     }
0813     return ieta;  // avoids compilation warning
0814   }
0815 
0816   // Functions to check whether a given (eta,depth) value is valid for a given subdetector
0817 
0818   inline std::vector<std::string> HcalEtaPhiHistNames() {
0819     std::vector<std::string> name;
0820     name.push_back("HB HE HF Depth 1 ");
0821     name.push_back("HB HE HF Depth 2 ");
0822     name.push_back("HE Depth 3 ");
0823     name.push_back("HO Depth 4 ");
0824     return name;
0825   }
0826 
0827   inline bool isHB(int etabin, int depth) {
0828     if (depth > 2)
0829       return false;
0830     else if (depth < 1)
0831       return false;
0832     else {
0833       int ieta = CalcIeta(etabin, depth);
0834       if (ieta == -9999)
0835         return false;
0836       if (depth == 1) {
0837         if (abs(ieta) <= 16)
0838           return true;
0839         else
0840           return false;
0841       } else if (depth == 2) {
0842         if (abs(ieta) == 15 || abs(ieta) == 16)
0843           return true;
0844         else
0845           return false;
0846       }
0847     }
0848     return false;
0849   }
0850 
0851   inline bool isHE(int etabin, int depth) {
0852     if (depth > 3)
0853       return false;
0854     else if (depth < 1)
0855       return false;
0856     else {
0857       int ieta = CalcIeta(etabin, depth);
0858       if (ieta == -9999)
0859         return false;
0860       if (depth == 1) {
0861         if (abs(ieta) >= 17 && abs(ieta) <= 28)
0862           return true;
0863         if (ieta == -29 && etabin == 13)
0864           return true;  // HE -29
0865         if (ieta == 29 && etabin == 71)
0866           return true;  // HE +29
0867       } else if (depth == 2) {
0868         if (abs(ieta) >= 17 && abs(ieta) <= 28)
0869           return true;
0870         if (ieta == -29 && etabin == 13)
0871           return true;  // HE -29
0872         if (ieta == 29 && etabin == 43)
0873           return true;  // HE +29
0874       } else if (depth == 3)
0875         return true;
0876     }
0877     return false;
0878   }
0879 
0880   inline bool isHF(int etabin, int depth) {
0881     if (depth > 2)
0882       return false;
0883     else if (depth < 1)
0884       return false;
0885     else {
0886       int ieta = CalcIeta(etabin, depth);
0887       if (ieta == -9999)
0888         return false;
0889       if (depth == 1) {
0890         if (ieta == -29 && etabin == 13)
0891           return false;  // HE -29
0892         else if (ieta == 29 && etabin == 71)
0893           return false;  // HE +29
0894         else if (abs(ieta) >= 29)
0895           return true;
0896       } else if (depth == 2) {
0897         if (ieta == -29 && etabin == 13)
0898           return false;  // HE -29
0899         else if (ieta == 29 && etabin == 43)
0900           return false;  // HE +29
0901         else if (abs(ieta) >= 29)
0902           return true;
0903       }
0904     }
0905     return false;
0906   }
0907 
0908   inline bool isHO(int etabin, int depth) {
0909     if (depth != 4)
0910       return false;
0911     int ieta = CalcIeta(etabin, depth);
0912     if (ieta != -9999)
0913       return true;
0914     return false;
0915   }
0916 
0917   // Checks whether HO region contains SiPM
0918 
0919   inline bool isSiPM(int ieta, int iphi, int depth) {
0920     if (depth != 4)
0921       return false;
0922     // HOP1
0923     if (ieta >= 5 && ieta <= 10 && iphi >= 47 && iphi <= 58)
0924       return true;
0925     // HOP2
0926     if (ieta >= 11 && ieta <= 15 && iphi >= 59 && iphi <= 70)
0927       return true;
0928     return false;
0929   }  // bool isSiPM
0930 
0931   // Checks whether (subdet, ieta, iphi, depth) value is a valid Hcal cell
0932 
0933   inline bool validDetId(HcalSubdetector sd, int ies, int ip, int dp) {
0934     // inputs are (subdetector, ieta, iphi, depth)
0935     // stolen from latest version of DataFormats/HcalDetId/src/HcalDetId.cc (not yet available in CMSSW_2_1_9)
0936 
0937     const int ie(abs(ies));
0938 
0939     return (
0940         (ip >= 1) && (ip <= 72) && (dp >= 1) && (ie >= 1) &&
0941         (((sd == HcalBarrel) && (((ie <= 14) && (dp == 1)) || (((ie == 15) || (ie == 16)) && (dp <= 2)))) ||
0942          ((sd == HcalEndcap) &&
0943           (((ie == 16) && (dp == 3)) || ((ie == 17) && (dp == 1)) || ((ie >= 18) && (ie <= 20) && (dp <= 2)) ||
0944            ((ie >= 21) && (ie <= 26) && (dp <= 2) && (ip % 2 == 1)) ||
0945            ((ie >= 27) && (ie <= 28) && (dp <= 3) && (ip % 2 == 1)) || ((ie == 29) && (dp <= 2) && (ip % 2 == 1)))) ||
0946          ((sd == HcalOuter) && (ie <= 15) && (dp == 4)) ||
0947          ((sd == HcalForward) && (dp <= 2) &&
0948           (((ie >= 29) && (ie <= 39) && (ip % 2 == 1)) || ((ie >= 40) && (ie <= 41) && (ip % 4 == 3))))));
0949 
0950   }  // bool validDetId(HcalSubdetector sd, int ies, int ip, int dp)
0951 
0952   // Sets eta, phi labels for 'summary' eta-phi plots (identical to Depth 1 Eta-Phi labelling)
0953 
0954   inline void SetEtaPhiLabels(TH2F& h) {
0955     std::stringstream label;
0956     for (int i = -41; i <= -29; i = i + 2) {
0957       label << i;
0958       h.GetXaxis()->SetBinLabel(i + 42, label.str().c_str());
0959       label.str("");
0960     }
0961     h.GetXaxis()->SetBinLabel(14, "-29HE");
0962 
0963     // offset by one for HE
0964     for (int i = -27; i <= 27; i = i + 2) {
0965       label << i;
0966       h.GetXaxis()->SetBinLabel(i + 43, label.str().c_str());
0967       label.str("");
0968     }
0969     h.GetXaxis()->SetBinLabel(72, "29HE");
0970     for (int i = 29; i <= 41; i = i + 2) {
0971       label << i;
0972       h.GetXaxis()->SetBinLabel(i + 44, label.str().c_str());
0973       label.str("");
0974     }
0975     return;
0976   }
0977 
0978   // Fill Unphysical bins in histograms
0979   inline void FillUnphysicalHEHFBins(std::vector<TH2F>& hh) {
0980     int ieta = 0;
0981     int iphi = 0;
0982     // First 2 depths have 5-10-20 degree corrections
0983     for (unsigned int d = 0; d < 3; ++d) {
0984       //BUG CAN BE HERE:
0985       //if (hh[d] != 0) continue;
0986 
0987       for (int eta = 0; eta < hh[d].GetNbinsX(); ++eta) {
0988         ieta = CalcIeta(eta, d + 1);
0989         if (ieta == -9999 || abs(ieta) < 21)
0990           continue;
0991         for (int phi = 0; phi < hh[d].GetNbinsY(); ++phi) {
0992           iphi = phi + 1;
0993           if (iphi % 2 == 1 && abs(ieta) < 40 && iphi < 73) {
0994             hh[d].SetBinContent(eta + 1, iphi + 1, hh[d].GetBinContent(eta + 1, iphi));
0995           }
0996           // last two eta strips span 20 degrees in phi
0997           // Fill the phi cell above iphi, and the 2 below it
0998           else if (abs(ieta) > 39 && iphi % 4 == 3 && iphi < 73) {
0999             //ieta=40, iphi=3 covers iphi 3,4,5,6
1000             hh[d].SetBinContent(eta + 1, (iphi) % 72 + 1, hh[d].GetBinContent(eta + 1, iphi));
1001             hh[d].SetBinContent(eta + 1, (iphi + 1) % 72 + 1, hh[d].GetBinContent(eta + 1, iphi));
1002             hh[d].SetBinContent(eta + 1, (iphi + 2) % 72 + 1, hh[d].GetBinContent(eta + 1, iphi));
1003           }
1004         }  // for (int phi...)
1005       }    // for (int eta...)
1006     }      // for (int d=0;...)
1007     // no corrections needed for HO (depth 4)
1008     return;
1009   }  // FillUnphysicalHEHFBins(MonitorElement* hh)
1010 
1011   //Fill unphysical bins for single ME
1012   inline void FillUnphysicalHEHFBins(TH2F& hh) {
1013     // Fills unphysical HE/HF bins for Summary Histogram
1014     // Summary Histogram is binned with the same binning as the Depth 1 EtaPhiHists
1015 
1016     //CAN BE BUG HERE:
1017     //if (hh==0) return;
1018 
1019     int ieta = 0;
1020     int iphi = 0;
1021     int etabins = hh.GetNbinsX();
1022     int phibins = hh.GetNbinsY();
1023     float binval = 0;
1024     for (int eta = 0; eta < etabins; ++eta)  // loop over eta bins
1025     {
1026       ieta = CalcIeta(eta, 1);
1027       if (ieta == -9999 || abs(ieta) < 21)
1028         continue;  // ignore etas that don't exist, or that have 5 degree phi binning
1029 
1030       for (int phi = 0; phi < phibins; ++phi) {
1031         iphi = phi + 1;
1032         if (iphi % 2 == 1 && abs(ieta) < 40 && iphi < 73)  // 10 degree phi binning condition
1033         {
1034           binval = hh.GetBinContent(eta + 1, iphi);
1035           hh.SetBinContent(eta + 1, iphi + 1, binval);
1036         }                                                       // if (iphi%2==1...)
1037         else if (abs(ieta) > 39 && iphi % 4 == 3 && iphi < 73)  // 20 degree phi binning condition
1038         {
1039           // Set last two eta strips where each cell spans 20 degrees in phi
1040           // Set next phi cell above iphi, and 2 cells below the actual cell
1041           hh.SetBinContent(eta + 1, (iphi) % 72 + 1, hh.GetBinContent(eta + 1, iphi));
1042           hh.SetBinContent(eta + 1, (iphi + 1) % 72 + 1, hh.GetBinContent(eta + 1, iphi));
1043           hh.SetBinContent(eta + 1, (iphi + 2) % 72 + 1, hh.GetBinContent(eta + 1, iphi));
1044         }  // else if (abs(ieta)>39 ...)
1045       }    // for (int phi=0;phi<72;++phi)
1046 
1047     }  // for (int eta=0; eta< (etaBins_-2);++eta)
1048 
1049     return;
1050   }  // FillUnphysicalHEHFBins(std::vector<MonitorElement*> &hh)
1051 
1052   // special fill call based on detid -- eventually will need special treatment
1053   inline void Fill(HcalDetId& id, double val /*=1*/, std::vector<TH2F>& depth) {
1054     // If in HF, need to shift by 1 bin (-1 bin lower in -HF, +1 bin higher in +HF)
1055     if (id.subdet() == HcalForward)
1056       depth[id.depth() - 1].Fill(id.ieta() < 0 ? id.ieta() - 1 : id.ieta() + 1, id.iphi(), val);
1057     else
1058       depth[id.depth() - 1].Fill(id.ieta(), id.iphi(), val);
1059   }
1060 
1061   inline void Reset(std::vector<TH2F>& depth) {
1062     for (unsigned int d = 0; d < depth.size(); d++)
1063       //BUG CAN BE HERE:
1064       //if(depth[d])
1065       depth[d].Reset();
1066   }  // void Reset(void)
1067 
1068   inline void setup(std::vector<TH2F>& depth, std::string name, std::string units = "") {
1069     std::string unittitle, unitname;
1070     if (units.empty()) {
1071       unitname = units;
1072       unittitle = "No Units";
1073     } else {
1074       unitname = " " + units;
1075       unittitle = units;
1076     }
1077 
1078     // Push back depth plots
1079     ////////Create 4 plots:
1080     //1. create first plot
1081     depth.push_back(TH2F(("HB HE HF Depth 1 " + name + unitname).c_str(),
1082                          (name + " Depth 1 -- HB HE HF (" + unittitle + ")").c_str(),
1083                          85,
1084                          -42.5,
1085                          42.5,
1086                          72,
1087                          0.5,
1088                          72.5));
1089 
1090     //2.1 prepare second plot
1091     float ybins[73];
1092     for (int i = 0; i <= 72; i++)
1093       ybins[i] = (float)(i + 0.5);
1094     float xbinsd2[] = {-42.5, -41.5, -40.5, -39.5, -38.5, -37.5, -36.5, -35.5, -34.5, -33.5, -32.5, -31.5,
1095                        -30.5, -29.5, -28.5, -27.5, -26.5, -25.5, -24.5, -23.5, -22.5, -21.5, -20.5, -19.5,
1096                        -18.5, -17.5, -16.5, -15.5, -14.5, 14.5,  15.5,  16.5,  17.5,  18.5,  19.5,  20.5,
1097                        21.5,  22.5,  23.5,  24.5,  25.5,  26.5,  27.5,  28.5,  29.5,  30.5,  31.5,  32.5,
1098                        33.5,  34.5,  35.5,  36.5,  37.5,  38.5,  39.5,  40.5,  41.5,  42.5};
1099 
1100     //2.2 create second plot
1101     depth.push_back(TH2F(("HB HE HF Depth 2 " + name + unitname).c_str(),
1102                          (name + " Depth 2 -- HB HE HF (" + unittitle + ")").c_str(),
1103                          57,
1104                          xbinsd2,
1105                          72,
1106                          ybins));
1107 
1108     //3.1 Set up variable-sized bins for HE depth 3 (MonitorElement also requires phi bins to be entered in array format)
1109     float xbins[] = {-28.5, -27.5, -26.5, -16.5, -15.5, 15.5, 16.5, 26.5, 27.5, 28.5};
1110     //3.2
1111     depth.push_back(TH2F(("HE Depth 3 " + name + unitname).c_str(),
1112                          (name + " Depth 3 -- HE (" + unittitle + ")").c_str(),
1113                          // Use variable-sized eta bins
1114                          9,
1115                          xbins,
1116                          72,
1117                          ybins));
1118 
1119     //4.1 HO bins are fixed width, but cover a smaller eta range (-15 -> 15)
1120     depth.push_back(TH2F(("HO Depth 4 " + name + unitname).c_str(),
1121                          (name + " Depth 4 -- HO (" + unittitle + ")").c_str(),
1122                          31,
1123                          -15.5,
1124                          15.5,
1125                          72,
1126                          0.5,
1127                          72.5));
1128 
1129     for (unsigned int i = 0; i < depth.size(); ++i)
1130       depth[i].Draw("colz");
1131 
1132     setBinLabels(depth);  // set axis titles, special bins
1133   }
1134 
1135   inline void fillOneGain(std::vector<TH2F>& graphData,
1136                           HcalGains::tAllContWithNames& allContainers,
1137                           std::string name,
1138                           int id,
1139                           std::string units = "") {
1140     setup(graphData, name);
1141 
1142     std::stringstream x;
1143     // Change the titles of each individual histogram
1144     for (unsigned int d = 0; d < graphData.size(); ++d) {
1145       graphData[d].Reset();
1146       x << "Gain " << id << " for HCAL depth " << d + 1;
1147 
1148       //BUG CAN BE HERE:
1149       //if (ChannelStatus->depth[d])
1150       graphData[d].SetTitle(
1151           x.str()
1152               .c_str());  // replace "setTitle" with "SetTitle", since you are using TH2F objects instead of MonitorElements
1153       x.str("");
1154     }
1155 
1156     HcalDetId hcal_id;
1157     int ieta, depth, iphi;
1158 
1159     //main loop
1160     // get all containers with names
1161     //HcalGains::tAllContWithNames allContainers = object().getAllContainers();
1162 
1163     //ITERATORS AND VALUES:
1164     HcalGains::tAllContWithNames::const_iterator iter;
1165     std::vector<HcalGain>::const_iterator contIter;
1166     float gain = 0.0;
1167 
1168     //Run trough given id gain:
1169     int i = id;
1170 
1171     //run trough all pair containers
1172     for (iter = allContainers.begin(); iter != allContainers.end(); ++iter) {
1173       //Run trough all values:
1174       for (contIter = (*iter).second.begin(); contIter != (*iter).second.end(); ++contIter) {
1175         hcal_id = HcalDetId((uint32_t)(*contIter).rawId());
1176 
1177         depth = hcal_id.depth();
1178         if (depth < 1 || depth > 4)
1179           continue;
1180 
1181         ieta = hcal_id.ieta();
1182         iphi = hcal_id.iphi();
1183 
1184         if (hcal_id.subdet() == HcalForward)
1185           ieta > 0 ? ++ieta : --ieta;
1186 
1187         //GET VALUE:
1188         gain = (*contIter).getValue(i);
1189         //logstatus = log2(1.*channelBits)+1;
1190 
1191         //FILLING GOES HERE:
1192         graphData[depth - 1].Fill(ieta, iphi, gain);
1193 
1194         //FOR DEBUGGING:
1195         //std::cout << "ieta: " << ieta << "; iphi: " << iphi << "; logstatus: " << logstatus << "; channelBits: " << channelBits<< std::endl;
1196       }
1197     }
1198   }
1199 
1200   //FillUnphysicalHEHFBins(graphData);
1201   //return ("kasdasd");
1202 
1203   class ADataRepr  //Sample base class for c++ inheritance tutorial
1204   {
1205   public:
1206     ADataRepr(unsigned int d) : m_total(d){};
1207     virtual ~ADataRepr(){};
1208     unsigned int nr, id;
1209     std::stringstream filename, rootname, plotname;
1210 
1211     void fillOneGain(std::vector<TH2F>& graphData, std::string units = "") {
1212       std::stringstream ss("");
1213 
1214       if (m_total == 1)
1215         ss << rootname.str() << " for HCAL depth ";
1216       else
1217         ss << rootname.str() << nr << " for HCAL depth ";
1218 
1219       setup(graphData, ss.str());
1220 
1221       // Change the titles of each individual histogram
1222       for (unsigned int d = 0; d < graphData.size(); ++d) {
1223         graphData[d].Reset();
1224 
1225         ss.str("");
1226         if (m_total == 1)
1227           ss << plotname.str() << " for HCAL depth " << d + 1;
1228         else
1229           ss << plotname.str() << nr << " for HCAL depth " << d + 1;
1230 
1231         //BUG CAN BE HERE:
1232         //if (ChannelStatus->depth[d])
1233         graphData[d].SetTitle(
1234             ss.str()
1235                 .c_str());  // replace "setTitle" with "SetTitle", since you are using TH2F objects instead of MonitorElements
1236         ss.str("");
1237       }
1238       //overload this function:
1239       doFillIn(graphData);
1240 
1241       FillUnphysicalHEHFBins(graphData);
1242 
1243       ss.str("");
1244       if (m_total == 1)
1245         ss << filename.str() << ".png";
1246       else
1247         ss << filename.str() << nr << ".png";
1248       draw(graphData, ss.str());
1249       //FOR DEBUGGING:
1250       //std::cout << "ieta: " << ieta << "; iphi: " << iphi << "; logstatus: " << logstatus << "; channelBits: " << channelBits<< std::endl;
1251     }
1252 
1253   protected:
1254     unsigned int m_total;
1255     HcalDetId hcal_id;
1256     int ieta, depth, iphi;
1257 
1258     virtual void doFillIn(std::vector<TH2F>& graphData) = 0;
1259 
1260   private:
1261     void draw(std::vector<TH2F>& graphData, std::string filename) {
1262       //Drawing...
1263       // use David's palette
1264       gStyle->SetPalette(1);
1265       const Int_t NCont = 999;
1266       gStyle->SetNumberContours(NCont);
1267       TCanvas canvas("CC map", "CC map", 840, 369 * 4);
1268 
1269       TPad pad1("pad1", "pad1", 0.0, 0.75, 1.0, 1.0);
1270       pad1.Draw();
1271       TPad pad2("pad2", "pad2", 0.0, 0.5, 1.0, 0.75);
1272       pad2.Draw();
1273       TPad pad3("pad3", "pad3", 0.0, 0.25, 1.0, 0.5);
1274       pad3.Draw();
1275       TPad pad4("pad4", "pad4", 0.0, 0.0, 1.0, 0.25);
1276       pad4.Draw();
1277 
1278       pad1.cd();
1279       graphData[0].SetStats(false);
1280       graphData[0].Draw("colz");
1281 
1282       pad2.cd();
1283       graphData[1].SetStats(false);
1284       graphData[1].Draw("colz");
1285 
1286       pad3.cd();
1287       graphData[2].SetStats(false);
1288       graphData[2].Draw("colz");
1289 
1290       pad4.cd();
1291       graphData[3].SetStats(false);
1292       graphData[3].Draw("colz");
1293 
1294       canvas.SaveAs(filename.c_str());
1295     }
1296 
1297     void setup(std::vector<TH2F>& depth, std::string name, std::string units = "") {
1298       std::string unittitle, unitname;
1299       if (units.empty()) {
1300         unitname = units;
1301         unittitle = "No Units";
1302       } else {
1303         unitname = " " + units;
1304         unittitle = units;
1305       }
1306 
1307       // Push back depth plots
1308       ////////Create 4 plots:
1309       //1. create first plot
1310       depth.push_back(TH2F(("HB HE HF Depth 1 " + name + unitname).c_str(),
1311                            (name + " Depth 1 -- HB HE HF (" + unittitle + ")").c_str(),
1312                            85,
1313                            -42.5,
1314                            42.5,
1315                            72,
1316                            0.5,
1317                            72.5));
1318 
1319       //2.1 prepare second plot
1320       float ybins[73];
1321       for (int i = 0; i <= 72; i++)
1322         ybins[i] = (float)(i + 0.5);
1323       float xbinsd2[] = {-42.5, -41.5, -40.5, -39.5, -38.5, -37.5, -36.5, -35.5, -34.5, -33.5, -32.5, -31.5,
1324                          -30.5, -29.5, -28.5, -27.5, -26.5, -25.5, -24.5, -23.5, -22.5, -21.5, -20.5, -19.5,
1325                          -18.5, -17.5, -16.5, -15.5, -14.5, 14.5,  15.5,  16.5,  17.5,  18.5,  19.5,  20.5,
1326                          21.5,  22.5,  23.5,  24.5,  25.5,  26.5,  27.5,  28.5,  29.5,  30.5,  31.5,  32.5,
1327                          33.5,  34.5,  35.5,  36.5,  37.5,  38.5,  39.5,  40.5,  41.5,  42.5};
1328 
1329       //2.2 create second plot
1330       depth.push_back(TH2F(("HB HE HF Depth 2 " + name + unitname).c_str(),
1331                            (name + " Depth 2 -- HB HE HF (" + unittitle + ")").c_str(),
1332                            57,
1333                            xbinsd2,
1334                            72,
1335                            ybins));
1336 
1337       //3.1 Set up variable-sized bins for HE depth 3 (MonitorElement also requires phi bins to be entered in array format)
1338       float xbins[] = {-28.5, -27.5, -26.5, -16.5, -15.5, 15.5, 16.5, 26.5, 27.5, 28.5};
1339       //3.2
1340       depth.push_back(TH2F(("HE Depth 3 " + name + unitname).c_str(),
1341                            (name + " Depth 3 -- HE (" + unittitle + ")").c_str(),
1342                            // Use variable-sized eta bins
1343                            9,
1344                            xbins,
1345                            72,
1346                            ybins));
1347 
1348       //4.1 HO bins are fixed width, but cover a smaller eta range (-15 -> 15)
1349       depth.push_back(TH2F(("HO Depth 4 " + name + unitname).c_str(),
1350                            (name + " Depth 4 -- HO (" + unittitle + ")").c_str(),
1351                            31,
1352                            -15.5,
1353                            15.5,
1354                            72,
1355                            0.5,
1356                            72.5));
1357 
1358       for (unsigned int i = 0; i < depth.size(); ++i)
1359         depth[i].Draw("colz");
1360 
1361       setBinLabels(depth);  // set axis titles, special bins
1362     }
1363 
1364     //functions for making plot:
1365     void setBinLabels(std::vector<TH2F>& depth) {
1366       // Set labels for all depth histograms
1367       for (unsigned int i = 0; i < depth.size(); ++i) {
1368         depth[i].SetXTitle("i#eta");
1369         depth[i].SetYTitle("i#phi");
1370       }
1371 
1372       std::stringstream label;
1373 
1374       // set label on every other bin
1375       for (int i = -41; i <= -29; i = i + 2) {
1376         label << i;
1377         depth[0].GetXaxis()->SetBinLabel(i + 42, label.str().c_str());
1378         depth[1].GetXaxis()->SetBinLabel(i + 42, label.str().c_str());
1379         label.str("");
1380       }
1381       depth[0].GetXaxis()->SetBinLabel(14, "-29HE");
1382       depth[1].GetXaxis()->SetBinLabel(14, "-29HE");
1383 
1384       // offset by one for HE
1385       for (int i = -27; i <= 27; i = i + 2) {
1386         label << i;
1387         depth[0].GetXaxis()->SetBinLabel(i + 43, label.str().c_str());
1388         label.str("");
1389       }
1390       depth[0].GetXaxis()->SetBinLabel(72, "29HE");
1391       for (int i = 29; i <= 41; i = i + 2) {
1392         label << i;
1393         depth[0].GetXaxis()->SetBinLabel(i + 44, label.str().c_str());
1394         label.str("");
1395       }
1396       for (int i = 16; i <= 28; i = i + 2) {
1397         label << i - 43;
1398         depth[1].GetXaxis()->SetBinLabel(i, label.str().c_str());
1399         label.str("");
1400       }
1401       depth[1].GetXaxis()->SetBinLabel(29, "NULL");
1402       for (int i = 15; i <= 27; i = i + 2) {
1403         label << i;
1404         depth[1].GetXaxis()->SetBinLabel(i + 15, label.str().c_str());
1405         label.str("");
1406       }
1407 
1408       depth[1].GetXaxis()->SetBinLabel(44, "29HE");
1409       for (int i = 29; i <= 41; i = i + 2) {
1410         label << i;
1411         depth[1].GetXaxis()->SetBinLabel(i + 16, label.str().c_str());
1412         label.str("");
1413       }
1414 
1415       // HE depth 3 labels;
1416       depth[2].GetXaxis()->SetBinLabel(1, "-28");
1417       depth[2].GetXaxis()->SetBinLabel(2, "-27");
1418       depth[2].GetXaxis()->SetBinLabel(3, "Null");
1419       depth[2].GetXaxis()->SetBinLabel(4, "-16");
1420       depth[2].GetXaxis()->SetBinLabel(5, "Null");
1421       depth[2].GetXaxis()->SetBinLabel(6, "16");
1422       depth[2].GetXaxis()->SetBinLabel(7, "Null");
1423       depth[2].GetXaxis()->SetBinLabel(8, "27");
1424       depth[2].GetXaxis()->SetBinLabel(9, "28");
1425     }
1426   };
1427 }  // namespace HcalObjRepresent
1428 #endif
[ 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!