File indexing completed on 2023-03-17 10:58:21
0001
0002 #ifndef FlavourHistograms2D_H
0003 #define FlavourHistograms2D_H
0004
0005 #include "DQMServices/Core/interface/DQMStore.h"
0006 #include "FWCore/ServiceRegistry/interface/Service.h"
0007
0008
0009
0010 #include "TH2F.h"
0011 #include "TProfile.h"
0012 #include "TCanvas.h"
0013 #include "TROOT.h"
0014 #include "TSystem.h"
0015 #include "TStyle.h"
0016
0017 #include "DQMOffline/RecoB/interface/Tools.h"
0018 #include "DQMOffline/RecoB/interface/HistoProviderDQM.h"
0019
0020 #include <iostream>
0021 #include <vector>
0022 #include <string>
0023
0024
0025
0026
0027 template <class T, class G>
0028 class FlavourHistograms2D {
0029 public:
0030 typedef dqm::legacy::DQMStore DQMStore;
0031 typedef dqm::legacy::MonitorElement MonitorElement;
0032
0033 FlavourHistograms2D(TString baseNameTitle_,
0034 TString baseNameDescription_,
0035 int nBinsX_,
0036 double lowerBoundX_,
0037 double upperBoundX_,
0038 int nBinsY_,
0039 double lowerBoundY_,
0040 double upperBoundY_,
0041 std::string folder,
0042 unsigned int mc,
0043 bool createProfile,
0044 DQMStore::IGetter &iget);
0045
0046 FlavourHistograms2D(TString baseNameTitle_,
0047 TString baseNameDescription_,
0048 int nBinsX_,
0049 double lowerBoundX_,
0050 double upperBoundX_,
0051 int nBinsY_,
0052 double lowerBoundY_,
0053 double upperBoundY_,
0054 bool statistics_,
0055 std::string folder,
0056 unsigned int mc,
0057 bool createProfile,
0058 DQMStore::IBooker &ibook);
0059
0060 virtual ~FlavourHistograms2D();
0061
0062
0063
0064
0065
0066
0067 void fill(const int &flavour, const T &variableX, const G &variableY) const;
0068 void fill(const int &flavour, const T &variableX, const G &variableY, const float &w) const;
0069
0070
0071 void fill(const int &flavour, const T *variableX, const G *variableY) const;
0072 void fill(const int &flavour, const T *variableX, const G *variableY, const float &w) const;
0073
0074 void settitle(const char *titleX, const char *titleY);
0075
0076 void plot(TPad *theCanvas = nullptr);
0077 void epsPlot(const std::string &name);
0078
0079
0080
0081 void divide(const FlavourHistograms2D<T, G> &bHD) const;
0082 void setEfficiencyFlag();
0083
0084 inline void SetMaximum(const double &max) { theMax = max; }
0085 inline void SetMinimum(const double &min) { theMin = min; }
0086
0087
0088 inline std::string baseNameTitle() const { return theBaseNameTitle; }
0089 inline std::string baseNameDescription() const { return theBaseNameDescription; }
0090 inline int nBinsX() const { return theNBinsX; }
0091 inline int nBinsY() const { return theNBinsY; }
0092 inline double lowerBoundX() const { return theLowerBoundX; }
0093 inline double upperBoundX() const { return theUpperBoundX; }
0094 inline double lowerBoundY() const { return theLowerBoundY; }
0095 inline double upperBoundY() const { return theUpperBoundY; }
0096 inline bool statistics() const { return theStatistics; }
0097
0098
0099 inline TH2F *histo_all() const { return theHisto_all->getTH2F(); }
0100 inline TH2F *histo_d() const { return theHisto_d->getTH2F(); }
0101 inline TH2F *histo_u() const { return theHisto_u->getTH2F(); }
0102 inline TH2F *histo_s() const { return theHisto_s->getTH2F(); }
0103 inline TH2F *histo_c() const { return theHisto_c->getTH2F(); }
0104 inline TH2F *histo_b() const { return theHisto_b->getTH2F(); }
0105 inline TH2F *histo_g() const { return theHisto_g->getTH2F(); }
0106 inline TH2F *histo_ni() const { return theHisto_ni->getTH2F(); }
0107 inline TH2F *histo_dus() const { return theHisto_dus->getTH2F(); }
0108 inline TH2F *histo_dusg() const { return theHisto_dusg->getTH2F(); }
0109 inline TH2F *histo_pu() const { return theHisto_pu->getTH2F(); }
0110
0111 TProfile *profile_all() const { return theProfile_all->getTProfile(); }
0112 TProfile *profile_d() const { return theProfile_d->getTProfile(); }
0113 TProfile *profile_u() const { return theProfile_u->getTProfile(); }
0114 TProfile *profile_s() const { return theProfile_s->getTProfile(); }
0115 TProfile *profile_c() const { return theProfile_c->getTProfile(); }
0116 TProfile *profile_b() const { return theProfile_b->getTProfile(); }
0117 TProfile *profile_g() const { return theProfile_g->getTProfile(); }
0118 TProfile *profile_ni() const { return theProfile_ni->getTProfile(); }
0119 TProfile *profile_dus() const { return theProfile_dus->getTProfile(); }
0120 TProfile *profile_dusg() const { return theProfile_dusg->getTProfile(); }
0121 TProfile *profile_pu() const { return theProfile_pu->getTProfile(); }
0122
0123 std::vector<TH2F *> getHistoVector() const;
0124
0125 std::vector<TProfile *> getProfileVector() const;
0126
0127 protected:
0128 void fillVariable(const int &flavour, const T &varX, const G &varY, const float &w) const;
0129
0130
0131
0132
0133
0134
0135
0136
0137 int *theArrayDimension;
0138 int theMaxDimension;
0139 int theIndexToPlot;
0140
0141 std::string theBaseNameTitle;
0142 std::string theBaseNameDescription;
0143 int theNBinsX;
0144 int theNBinsY;
0145 double theLowerBoundX;
0146 double theUpperBoundX;
0147 double theLowerBoundY;
0148 double theUpperBoundY;
0149 bool theStatistics;
0150 double theMin, theMax;
0151
0152
0153 MonitorElement *theHisto_all;
0154 MonitorElement *theHisto_d;
0155 MonitorElement *theHisto_u;
0156 MonitorElement *theHisto_s;
0157 MonitorElement *theHisto_c;
0158 MonitorElement *theHisto_b;
0159 MonitorElement *theHisto_g;
0160 MonitorElement *theHisto_ni;
0161 MonitorElement *theHisto_dus;
0162 MonitorElement *theHisto_dusg;
0163 MonitorElement *theHisto_pu;
0164
0165
0166 MonitorElement *theProfile_all;
0167 MonitorElement *theProfile_d;
0168 MonitorElement *theProfile_u;
0169 MonitorElement *theProfile_s;
0170 MonitorElement *theProfile_c;
0171 MonitorElement *theProfile_b;
0172 MonitorElement *theProfile_g;
0173 MonitorElement *theProfile_ni;
0174 MonitorElement *theProfile_dus;
0175 MonitorElement *theProfile_dusg;
0176 MonitorElement *theProfile_pu;
0177
0178
0179
0180
0181 private:
0182 FlavourHistograms2D() {}
0183
0184 unsigned int mcPlots_;
0185 bool createProfile_;
0186 };
0187
0188 template <class T, class G>
0189 FlavourHistograms2D<T, G>::FlavourHistograms2D(TString baseNameTitle_,
0190 TString baseNameDescription_,
0191 int nBinsX_,
0192 double lowerBoundX_,
0193 double upperBoundX_,
0194 int nBinsY_,
0195 double lowerBoundY_,
0196 double upperBoundY_,
0197 std::string folder,
0198 unsigned int mc,
0199 bool createProfile,
0200 DQMStore::IGetter &iget)
0201 :
0202
0203 theMaxDimension(-1),
0204 theIndexToPlot(-1),
0205 theBaseNameTitle(baseNameTitle_.Data()),
0206 theBaseNameDescription(baseNameDescription_.Data()),
0207 theNBinsX(nBinsX_),
0208 theNBinsY(nBinsY_),
0209 theLowerBoundX(lowerBoundX_),
0210 theUpperBoundX(upperBoundX_),
0211 theLowerBoundY(lowerBoundY_),
0212 theUpperBoundY(upperBoundY_),
0213 theMin(-1.),
0214 theMax(-1.),
0215 mcPlots_(mc),
0216 createProfile_(createProfile) {
0217
0218 theArrayDimension = nullptr;
0219
0220 if (mcPlots_ % 2 == 0)
0221 theHisto_all = iget.get("Btag/" + folder + "/" + theBaseNameTitle + "ALL");
0222 else
0223 theHisto_all = nullptr;
0224 if (mcPlots_) {
0225 if (mcPlots_ > 2) {
0226 theHisto_d = iget.get("Btag/" + folder + "/" + theBaseNameTitle + "D");
0227 theHisto_u = iget.get("Btag/" + folder + "/" + theBaseNameTitle + "U");
0228 theHisto_s = iget.get("Btag/" + folder + "/" + theBaseNameTitle + "S");
0229 theHisto_g = iget.get("Btag/" + folder + "/" + theBaseNameTitle + "G");
0230 theHisto_dus = iget.get("Btag/" + folder + "/" + theBaseNameTitle + "DUS");
0231 } else {
0232 theHisto_d = nullptr;
0233 theHisto_u = nullptr;
0234 theHisto_s = nullptr;
0235 theHisto_g = nullptr;
0236 theHisto_dus = nullptr;
0237 }
0238 theHisto_c = iget.get("Btag/" + folder + "/" + theBaseNameTitle + "C");
0239 theHisto_b = iget.get("Btag/" + folder + "/" + theBaseNameTitle + "B");
0240 theHisto_ni = iget.get("Btag/" + folder + "/" + theBaseNameTitle + "NI");
0241 theHisto_dusg = iget.get("Btag/" + folder + "/" + theBaseNameTitle + "DUSG");
0242 theHisto_pu = iget.get("Btag/" + folder + "/" + theBaseNameTitle + "PU");
0243 } else {
0244 theHisto_d = nullptr;
0245 theHisto_u = nullptr;
0246 theHisto_s = nullptr;
0247 theHisto_c = nullptr;
0248 theHisto_b = nullptr;
0249 theHisto_g = nullptr;
0250 theHisto_ni = nullptr;
0251 theHisto_dus = nullptr;
0252 theHisto_dusg = nullptr;
0253 theHisto_pu = nullptr;
0254 }
0255
0256 if (createProfile_) {
0257 if (mcPlots_ % 2 == 0)
0258 theProfile_all = iget.get("Btag/" + folder + "/" + theBaseNameTitle + "_Profile_ALL");
0259 else
0260 theProfile_all = nullptr;
0261 if (mcPlots_) {
0262 if (mcPlots_ > 2) {
0263 theProfile_d = iget.get("Btag/" + folder + "/" + theBaseNameTitle + "_Profile_D");
0264 theProfile_u = iget.get("Btag/" + folder + "/" + theBaseNameTitle + "_Profile_U");
0265 theProfile_s = iget.get("Btag/" + folder + "/" + theBaseNameTitle + "_Profile_S");
0266 theProfile_g = iget.get("Btag/" + folder + "/" + theBaseNameTitle + "_Profile_G");
0267 theProfile_dus = iget.get("Btag/" + folder + "/" + theBaseNameTitle + "_Profile_DUS");
0268 } else {
0269 theProfile_d = nullptr;
0270 theProfile_u = nullptr;
0271 theProfile_s = nullptr;
0272 theProfile_g = nullptr;
0273 theProfile_dus = nullptr;
0274 }
0275 theProfile_c = iget.get("Btag/" + folder + "/" + theBaseNameTitle + "_Profile_C");
0276 theProfile_b = iget.get("Btag/" + folder + "/" + theBaseNameTitle + "_Profile_B");
0277 theProfile_ni = iget.get("Btag/" + folder + "/" + theBaseNameTitle + "_Profile_NI");
0278 theProfile_dusg = iget.get("Btag/" + folder + "/" + theBaseNameTitle + "_Profile_DUSG");
0279 theProfile_pu = iget.get("Btag/" + folder + "/" + theBaseNameTitle + "_Profile_PU");
0280 } else {
0281 theProfile_d = nullptr;
0282 theProfile_u = nullptr;
0283 theProfile_s = nullptr;
0284 theProfile_c = nullptr;
0285 theProfile_b = nullptr;
0286 theProfile_g = nullptr;
0287 theProfile_ni = nullptr;
0288 theProfile_dus = nullptr;
0289 theProfile_dusg = nullptr;
0290 theProfile_pu = nullptr;
0291 }
0292 } else {
0293 theProfile_all = nullptr;
0294 theProfile_d = nullptr;
0295 theProfile_u = nullptr;
0296 theProfile_s = nullptr;
0297 theProfile_c = nullptr;
0298 theProfile_b = nullptr;
0299 theProfile_g = nullptr;
0300 theProfile_ni = nullptr;
0301 theProfile_dus = nullptr;
0302 theProfile_dusg = nullptr;
0303 theProfile_pu = nullptr;
0304 }
0305 }
0306
0307 template <class T, class G>
0308 FlavourHistograms2D<T, G>::FlavourHistograms2D(TString baseNameTitle_,
0309 TString baseNameDescription_,
0310 int nBinsX_,
0311 double lowerBoundX_,
0312 double upperBoundX_,
0313 int nBinsY_,
0314 double lowerBoundY_,
0315 double upperBoundY_,
0316 bool statistics_,
0317 std::string folder,
0318 unsigned int mc,
0319 bool createProfile,
0320 DQMStore::IBooker &ibook)
0321 : theMaxDimension(-1),
0322 theIndexToPlot(-1),
0323 theBaseNameTitle(baseNameTitle_.Data()),
0324 theBaseNameDescription(baseNameDescription_.Data()),
0325 theNBinsX(nBinsX_),
0326 theNBinsY(nBinsY_),
0327 theLowerBoundX(lowerBoundX_),
0328 theUpperBoundX(upperBoundX_),
0329 theLowerBoundY(lowerBoundY_),
0330 theUpperBoundY(upperBoundY_),
0331 theStatistics(statistics_),
0332 theMin(-1.),
0333 theMax(-1.),
0334 mcPlots_(mc),
0335 createProfile_(createProfile) {
0336
0337 theArrayDimension = nullptr;
0338
0339
0340 HistoProviderDQM prov("Btag", folder, ibook);
0341 if (mcPlots_ % 2 == 0)
0342 theHisto_all = (prov.book2D(theBaseNameTitle + "ALL",
0343 theBaseNameDescription + " all jets",
0344 theNBinsX,
0345 theLowerBoundX,
0346 theUpperBoundX,
0347 theNBinsY,
0348 theLowerBoundY,
0349 theUpperBoundY));
0350 else
0351 theHisto_all = nullptr;
0352 if (mcPlots_) {
0353 if (mcPlots_ > 2) {
0354 theHisto_d = (prov.book2D(theBaseNameTitle + "D",
0355 theBaseNameDescription + " d-jets",
0356 theNBinsX,
0357 theLowerBoundX,
0358 theUpperBoundX,
0359 theNBinsY,
0360 theLowerBoundY,
0361 theUpperBoundY));
0362 theHisto_u = (prov.book2D(theBaseNameTitle + "U",
0363 theBaseNameDescription + " u-jets",
0364 theNBinsX,
0365 theLowerBoundX,
0366 theUpperBoundX,
0367 theNBinsY,
0368 theLowerBoundY,
0369 theUpperBoundY));
0370 theHisto_s = (prov.book2D(theBaseNameTitle + "S",
0371 theBaseNameDescription + " s-jets",
0372 theNBinsX,
0373 theLowerBoundX,
0374 theUpperBoundX,
0375 theNBinsY,
0376 theLowerBoundY,
0377 theUpperBoundY));
0378 theHisto_g = (prov.book2D(theBaseNameTitle + "G",
0379 theBaseNameDescription + " g-jets",
0380 theNBinsX,
0381 theLowerBoundX,
0382 theUpperBoundX,
0383 theNBinsY,
0384 theLowerBoundY,
0385 theUpperBoundY));
0386 theHisto_dus = (prov.book2D(theBaseNameTitle + "DUS",
0387 theBaseNameDescription + " dus-jets",
0388 theNBinsX,
0389 theLowerBoundX,
0390 theUpperBoundX,
0391 theNBinsY,
0392 theLowerBoundY,
0393 theUpperBoundY));
0394 } else {
0395 theHisto_d = nullptr;
0396 theHisto_u = nullptr;
0397 theHisto_s = nullptr;
0398 theHisto_g = nullptr;
0399 theHisto_dus = nullptr;
0400 }
0401 theHisto_c = (prov.book2D(theBaseNameTitle + "C",
0402 theBaseNameDescription + " c-jets",
0403 theNBinsX,
0404 theLowerBoundX,
0405 theUpperBoundX,
0406 theNBinsY,
0407 theLowerBoundY,
0408 theUpperBoundY));
0409 theHisto_b = (prov.book2D(theBaseNameTitle + "B",
0410 theBaseNameDescription + " b-jets",
0411 theNBinsX,
0412 theLowerBoundX,
0413 theUpperBoundX,
0414 theNBinsY,
0415 theLowerBoundY,
0416 theUpperBoundY));
0417 theHisto_ni = (prov.book2D(theBaseNameTitle + "NI",
0418 theBaseNameDescription + " ni-jets",
0419 theNBinsX,
0420 theLowerBoundX,
0421 theUpperBoundX,
0422 theNBinsY,
0423 theLowerBoundY,
0424 theUpperBoundY));
0425 theHisto_dusg = (prov.book2D(theBaseNameTitle + "DUSG",
0426 theBaseNameDescription + " dusg-jets",
0427 theNBinsX,
0428 theLowerBoundX,
0429 theUpperBoundX,
0430 theNBinsY,
0431 theLowerBoundY,
0432 theUpperBoundY));
0433 theHisto_pu = (prov.book2D(theBaseNameTitle + "PU",
0434 theBaseNameDescription + " pu-jets",
0435 theNBinsX,
0436 theLowerBoundX,
0437 theUpperBoundX,
0438 theNBinsY,
0439 theLowerBoundY,
0440 theUpperBoundY));
0441 } else {
0442 theHisto_d = nullptr;
0443 theHisto_u = nullptr;
0444 theHisto_s = nullptr;
0445 theHisto_c = nullptr;
0446 theHisto_b = nullptr;
0447 theHisto_g = nullptr;
0448 theHisto_ni = nullptr;
0449 theHisto_dus = nullptr;
0450 theHisto_dusg = nullptr;
0451 theHisto_pu = nullptr;
0452 }
0453
0454 if (createProfile_) {
0455 if (mcPlots_ % 2 == 0)
0456 theProfile_all = (prov.bookProfile(theBaseNameTitle + "_Profile_ALL",
0457 theBaseNameDescription + " all jets",
0458 theNBinsX,
0459 theLowerBoundX,
0460 theUpperBoundX,
0461 theNBinsY,
0462 theLowerBoundY,
0463 theUpperBoundY));
0464 else
0465 theProfile_all = nullptr;
0466 if (mcPlots_) {
0467 if (mcPlots_ > 2) {
0468 theProfile_d = (prov.bookProfile(theBaseNameTitle + "_Profile_D",
0469 theBaseNameDescription + " d-jets",
0470 theNBinsX,
0471 theLowerBoundX,
0472 theUpperBoundX,
0473 theNBinsY,
0474 theLowerBoundY,
0475 theUpperBoundY));
0476 theProfile_u = (prov.bookProfile(theBaseNameTitle + "_Profile_U",
0477 theBaseNameDescription + " u-jets",
0478 theNBinsX,
0479 theLowerBoundX,
0480 theUpperBoundX,
0481 theNBinsY,
0482 theLowerBoundY,
0483 theUpperBoundY));
0484 theProfile_s = (prov.bookProfile(theBaseNameTitle + "_Profile_S",
0485 theBaseNameDescription + " s-jets",
0486 theNBinsX,
0487 theLowerBoundX,
0488 theUpperBoundX,
0489 theNBinsY,
0490 theLowerBoundY,
0491 theUpperBoundY));
0492 theProfile_g = (prov.bookProfile(theBaseNameTitle + "_Profile_G",
0493 theBaseNameDescription + " g-jets",
0494 theNBinsX,
0495 theLowerBoundX,
0496 theUpperBoundX,
0497 theNBinsY,
0498 theLowerBoundY,
0499 theUpperBoundY));
0500 theProfile_dus = (prov.bookProfile(theBaseNameTitle + "_Profile_DUS",
0501 theBaseNameDescription + " dus-jets",
0502 theNBinsX,
0503 theLowerBoundX,
0504 theUpperBoundX,
0505 theNBinsY,
0506 theLowerBoundY,
0507 theUpperBoundY));
0508 } else {
0509 theProfile_d = nullptr;
0510 theProfile_u = nullptr;
0511 theProfile_s = nullptr;
0512 theProfile_g = nullptr;
0513 theProfile_dus = nullptr;
0514 }
0515 theProfile_c = (prov.bookProfile(theBaseNameTitle + "_Profile_C",
0516 theBaseNameDescription + " c-jets",
0517 theNBinsX,
0518 theLowerBoundX,
0519 theUpperBoundX,
0520 theNBinsY,
0521 theLowerBoundY,
0522 theUpperBoundY));
0523 theProfile_b = (prov.bookProfile(theBaseNameTitle + "_Profile_B",
0524 theBaseNameDescription + " b-jets",
0525 theNBinsX,
0526 theLowerBoundX,
0527 theUpperBoundX,
0528 theNBinsY,
0529 theLowerBoundY,
0530 theUpperBoundY));
0531 theProfile_ni = (prov.bookProfile(theBaseNameTitle + "_Profile_NI",
0532 theBaseNameDescription + " ni-jets",
0533 theNBinsX,
0534 theLowerBoundX,
0535 theUpperBoundX,
0536 theNBinsY,
0537 theLowerBoundY,
0538 theUpperBoundY));
0539 theProfile_dusg = (prov.bookProfile(theBaseNameTitle + "_Profile_DUSG",
0540 theBaseNameDescription + " dusg-jets",
0541 theNBinsX,
0542 theLowerBoundX,
0543 theUpperBoundX,
0544 theNBinsY,
0545 theLowerBoundY,
0546 theUpperBoundY));
0547 theProfile_pu = (prov.bookProfile(theBaseNameTitle + "_Profile_PU",
0548 theBaseNameDescription + " pu-jets",
0549 theNBinsX,
0550 theLowerBoundX,
0551 theUpperBoundX,
0552 theNBinsY,
0553 theLowerBoundY,
0554 theUpperBoundY));
0555 } else {
0556 theProfile_d = nullptr;
0557 theProfile_u = nullptr;
0558 theProfile_s = nullptr;
0559 theProfile_c = nullptr;
0560 theProfile_b = nullptr;
0561 theProfile_g = nullptr;
0562 theProfile_ni = nullptr;
0563 theProfile_dus = nullptr;
0564 theProfile_dusg = nullptr;
0565 theProfile_pu = nullptr;
0566 }
0567 } else {
0568 theProfile_all = nullptr;
0569 theProfile_d = nullptr;
0570 theProfile_u = nullptr;
0571 theProfile_s = nullptr;
0572 theProfile_c = nullptr;
0573 theProfile_b = nullptr;
0574 theProfile_g = nullptr;
0575 theProfile_ni = nullptr;
0576 theProfile_dus = nullptr;
0577 theProfile_dusg = nullptr;
0578 theProfile_pu = nullptr;
0579 }
0580
0581 if (theStatistics) {
0582 if (theHisto_all)
0583 theHisto_all->enableSumw2();
0584 if (createProfile)
0585 if (theProfile_all)
0586 theProfile_all->getTProfile()->Sumw2();
0587 if (mcPlots_) {
0588 if (mcPlots_ > 2) {
0589 theHisto_d->enableSumw2();
0590 theHisto_u->enableSumw2();
0591 theHisto_s->enableSumw2();
0592 theHisto_g->enableSumw2();
0593 theHisto_dus->enableSumw2();
0594 }
0595 theHisto_c->enableSumw2();
0596 theHisto_b->enableSumw2();
0597 theHisto_ni->enableSumw2();
0598 theHisto_dusg->enableSumw2();
0599 theHisto_pu->enableSumw2();
0600
0601 if (createProfile) {
0602 if (mcPlots_ > 2) {
0603 theProfile_d->getTProfile()->Sumw2();
0604 theProfile_u->getTProfile()->Sumw2();
0605 theProfile_s->getTProfile()->Sumw2();
0606 theProfile_g->getTProfile()->Sumw2();
0607 theProfile_dus->getTProfile()->Sumw2();
0608 }
0609 theProfile_c->getTProfile()->Sumw2();
0610 theProfile_b->getTProfile()->Sumw2();
0611 theProfile_ni->getTProfile()->Sumw2();
0612 theProfile_dusg->getTProfile()->Sumw2();
0613 theProfile_pu->getTProfile()->Sumw2();
0614 }
0615 }
0616 }
0617 }
0618
0619 template <class T, class G>
0620 FlavourHistograms2D<T, G>::~FlavourHistograms2D() {}
0621
0622
0623 template <class T, class G>
0624 void FlavourHistograms2D<T, G>::fill(const int &flavour, const T &variableX, const G &variableY, const float &w) const {
0625
0626 fillVariable(flavour, variableX, variableY, w);
0627 }
0628
0629 template <class T, class G>
0630 void FlavourHistograms2D<T, G>::fill(const int &flavour, const T &variableX, const G &variableY) const {
0631 fill(flavour, variableX, variableY, 1.);
0632 }
0633
0634 template <class T, class G>
0635 void FlavourHistograms2D<T, G>::fill(const int &flavour, const T *variableX, const G *variableY, const float &w) const {
0636 if (theArrayDimension == nullptr) {
0637
0638 fillVariable(flavour, *variableX, *variableY, w);
0639 } else {
0640
0641 int iMax = *theArrayDimension;
0642 if (*theArrayDimension > theMaxDimension)
0643 iMax = theMaxDimension;
0644
0645 for (int i = 0; i != iMax; ++i) {
0646
0647 if ((theIndexToPlot < 0) || (i == theIndexToPlot)) {
0648 fillVariable(flavour, *(variableX + i), *(variableY + i), w);
0649 }
0650 }
0651
0652
0653 if (theIndexToPlot >= iMax) {
0654
0655 const T &theZeroT = static_cast<T>(0.0);
0656 const G &theZeroG = static_cast<T>(0.0);
0657 fillVariable(flavour, theZeroT, theZeroG, w);
0658 }
0659 }
0660 }
0661
0662 template <class T, class G>
0663 void FlavourHistograms2D<T, G>::fill(const int &flavour, const T *variableX, const G *variableY) const {
0664 fill(flavour, variableX, variableY, 1.);
0665 }
0666
0667 template <class T, class G>
0668 void FlavourHistograms2D<T, G>::settitle(const char *titleX, const char *titleY) {
0669 if (theHisto_all)
0670 theHisto_all->setAxisTitle(titleX);
0671 if (theHisto_all)
0672 theHisto_all->setAxisTitle(titleY, 2);
0673 if (mcPlots_) {
0674 if (theHisto_d)
0675 theHisto_d->setAxisTitle(titleX);
0676 if (theHisto_u)
0677 theHisto_u->setAxisTitle(titleX);
0678 if (theHisto_s)
0679 theHisto_s->setAxisTitle(titleX);
0680 if (theHisto_c)
0681 theHisto_c->setAxisTitle(titleX);
0682 if (theHisto_b)
0683 theHisto_b->setAxisTitle(titleX);
0684 if (theHisto_g)
0685 theHisto_g->setAxisTitle(titleX);
0686 if (theHisto_ni)
0687 theHisto_ni->setAxisTitle(titleX);
0688 if (theHisto_dus)
0689 theHisto_dus->setAxisTitle(titleX);
0690 if (theHisto_dusg)
0691 theHisto_dusg->setAxisTitle(titleX);
0692 if (theHisto_d)
0693 theHisto_d->setAxisTitle(titleY, 2);
0694 if (theHisto_u)
0695 theHisto_u->setAxisTitle(titleY, 2);
0696 if (theHisto_s)
0697 theHisto_s->setAxisTitle(titleY, 2);
0698 if (theHisto_c)
0699 theHisto_c->setAxisTitle(titleY, 2);
0700 if (theHisto_b)
0701 theHisto_b->setAxisTitle(titleY, 2);
0702 if (theHisto_g)
0703 theHisto_g->setAxisTitle(titleY, 2);
0704 if (theHisto_ni)
0705 theHisto_ni->setAxisTitle(titleY, 2);
0706 if (theHisto_dus)
0707 theHisto_dus->setAxisTitle(titleY, 2);
0708 if (theHisto_dusg)
0709 theHisto_dusg->setAxisTitle(titleY, 2);
0710 if (theHisto_pu)
0711 theHisto_pu->setAxisTitle(titleY, 2);
0712 }
0713
0714 if (createProfile_) {
0715 if (theProfile_all)
0716 theProfile_all->setAxisTitle(titleX);
0717 if (theProfile_all)
0718 theProfile_all->setAxisTitle(titleY, 2);
0719 if (mcPlots_) {
0720 if (theProfile_d)
0721 theProfile_d->setAxisTitle(titleX);
0722 if (theProfile_u)
0723 theProfile_u->setAxisTitle(titleX);
0724 if (theProfile_s)
0725 theProfile_s->setAxisTitle(titleX);
0726 if (theProfile_c)
0727 theProfile_c->setAxisTitle(titleX);
0728 if (theProfile_b)
0729 theProfile_b->setAxisTitle(titleX);
0730 if (theProfile_g)
0731 theProfile_g->setAxisTitle(titleX);
0732 if (theProfile_ni)
0733 theProfile_ni->setAxisTitle(titleX);
0734 if (theProfile_dus)
0735 theProfile_dus->setAxisTitle(titleX);
0736 if (theProfile_dusg)
0737 theProfile_dusg->setAxisTitle(titleX);
0738 if (theProfile_d)
0739 theProfile_d->setAxisTitle(titleY, 2);
0740 if (theProfile_u)
0741 theProfile_u->setAxisTitle(titleY, 2);
0742 if (theProfile_s)
0743 theProfile_s->setAxisTitle(titleY, 2);
0744 if (theProfile_c)
0745 theProfile_c->setAxisTitle(titleY, 2);
0746 if (theProfile_b)
0747 theProfile_b->setAxisTitle(titleY, 2);
0748 if (theProfile_g)
0749 theProfile_g->setAxisTitle(titleY, 2);
0750 if (theProfile_ni)
0751 theProfile_ni->setAxisTitle(titleY, 2);
0752 if (theProfile_dus)
0753 theProfile_dus->setAxisTitle(titleY, 2);
0754 if (theProfile_dusg)
0755 theProfile_dusg->setAxisTitle(titleY, 2);
0756 if (theProfile_pu)
0757 theProfile_pu->setAxisTitle(titleY, 2);
0758 }
0759 }
0760 }
0761
0762 template <class T, class G>
0763 void FlavourHistograms2D<T, G>::plot(TPad *theCanvas ) {
0764
0765 bool btppNI = false;
0766 bool btppColour = true;
0767
0768 if (theCanvas)
0769 theCanvas->cd();
0770
0771 RecoBTag::setTDRStyle()->cd();
0772 gPad->UseCurrentStyle();
0773 gPad->SetLogy(0);
0774 gPad->SetGridx(0);
0775 gPad->SetGridy(0);
0776 gPad->SetTitle(nullptr);
0777
0778 MonitorElement *histo[4];
0779 int col[4], lineStyle[4], markerStyle[4];
0780 int lineWidth = 1;
0781
0782 const double markerSize = gPad->GetWh() * gPad->GetHNDC() / 500.;
0783 histo[0] = theHisto_dusg;
0784 histo[1] = theHisto_b;
0785 histo[2] = theHisto_c;
0786 histo[3] = nullptr;
0787
0788 double max = theMax;
0789 if (theMax <= 0.) {
0790 max = theHisto_dusg->getTH2F()->GetMaximum();
0791 if (theHisto_b->getTH2F()->GetMaximum() > max)
0792 max = theHisto_b->getTH2F()->GetMaximum();
0793 if (theHisto_c->getTH2F()->GetMaximum() > max)
0794 max = theHisto_c->getTH2F()->GetMaximum();
0795 }
0796
0797 if (btppNI) {
0798 histo[3] = theHisto_ni;
0799 if (theHisto_ni->getTH2F()->GetMaximum() > max)
0800 max = theHisto_ni->getTH2F()->GetMaximum();
0801 }
0802
0803 if (btppColour) {
0804 col[0] = 4;
0805 col[1] = 2;
0806 col[2] = 6;
0807 col[3] = 3;
0808 lineStyle[0] = 1;
0809 lineStyle[1] = 1;
0810 lineStyle[2] = 1;
0811 lineStyle[3] = 1;
0812 markerStyle[0] = 20;
0813 markerStyle[1] = 21;
0814 markerStyle[2] = 22;
0815 markerStyle[3] = 23;
0816 lineWidth = 1;
0817 } else {
0818 col[1] = 1;
0819 col[2] = 1;
0820 col[3] = 1;
0821 col[0] = 1;
0822 lineStyle[0] = 2;
0823 lineStyle[1] = 1;
0824 lineStyle[2] = 3;
0825 lineStyle[3] = 4;
0826 markerStyle[0] = 20;
0827 markerStyle[1] = 21;
0828 markerStyle[2] = 22;
0829 markerStyle[3] = 23;
0830 }
0831
0832 histo[0]->setAxisTitle(theBaseNameDescription);
0833 histo[0]->getTH2F()->GetYaxis()->SetTitle("Arbitrary Units");
0834 histo[0]->getTH2F()->GetYaxis()->SetTitleOffset(1.25);
0835
0836 for (int i = 0; i != 4; ++i) {
0837 if (histo[i] == nullptr)
0838 continue;
0839 histo[i]->getTH2F()->SetStats(false);
0840 histo[i]->getTH2F()->SetLineStyle(lineStyle[i]);
0841 histo[i]->getTH2F()->SetLineWidth(lineWidth);
0842 histo[i]->getTH2F()->SetLineColor(col[i]);
0843 histo[i]->getTH2F()->SetMarkerStyle(markerStyle[i]);
0844 histo[i]->getTH2F()->SetMarkerColor(col[i]);
0845 histo[i]->getTH2F()->SetMarkerSize(markerSize);
0846 }
0847
0848 histo[0]->getTH2F()->SetMaximum(max * 1.05);
0849 if (theMin != -1.)
0850 histo[0]->getTH2F()->SetMinimum(theMin);
0851 histo[0]->getTH2F()->Draw();
0852 histo[1]->getTH2F()->Draw("Same");
0853 histo[2]->getTH2F()->Draw("Same");
0854 if (histo[3] != nullptr)
0855 histo[3]->getTH2F()->Draw("Same");
0856 }
0857
0858 template <class T, class G>
0859 void FlavourHistograms2D<T, G>::epsPlot(const std::string &name) {
0860 TCanvas tc(theBaseNameTitle.c_str(), theBaseNameDescription.c_str());
0861
0862 plot(&tc);
0863 tc.Print((name + theBaseNameTitle + ".eps").c_str());
0864 }
0865
0866
0867
0868 template <class T, class G>
0869 void FlavourHistograms2D<T, G>::divide(const FlavourHistograms2D<T, G> &bHD) const {
0870
0871
0872
0873
0874
0875 if (theHisto_all)
0876 theHisto_all->getTH2F()->Divide(theHisto_all->getTH2F(), bHD.histo_all(), 1.0, 1.0, "b");
0877 if (mcPlots_) {
0878 if (mcPlots_ > 2) {
0879 theHisto_d->getTH2F()->Divide(theHisto_d->getTH2F(), bHD.histo_d(), 1.0, 1.0, "b");
0880 theHisto_u->getTH2F()->Divide(theHisto_u->getTH2F(), bHD.histo_u(), 1.0, 1.0, "b");
0881 theHisto_s->getTH2F()->Divide(theHisto_s->getTH2F(), bHD.histo_s(), 1.0, 1.0, "b");
0882 theHisto_g->getTH2F()->Divide(theHisto_g->getTH2F(), bHD.histo_g(), 1.0, 1.0, "b");
0883 theHisto_dus->getTH2F()->Divide(theHisto_dus->getTH2F(), bHD.histo_dus(), 1.0, 1.0, "b");
0884 }
0885 theHisto_c->getTH2F()->Divide(theHisto_c->getTH2F(), bHD.histo_c(), 1.0, 1.0, "b");
0886 theHisto_b->getTH2F()->Divide(theHisto_b->getTH2F(), bHD.histo_b(), 1.0, 1.0, "b");
0887 theHisto_ni->getTH2F()->Divide(theHisto_ni->getTH2F(), bHD.histo_ni(), 1.0, 1.0, "b");
0888 theHisto_dusg->getTH2F()->Divide(theHisto_dusg->getTH2F(), bHD.histo_dusg(), 1.0, 1.0, "b");
0889 theHisto_pu->getTH2F()->Divide(theHisto_pu->getTH2F(), bHD.histo_pu(), 1.0, 1.0, "b");
0890 }
0891 }
0892
0893 template <class T, class G>
0894 void FlavourHistograms2D<T, G>::setEfficiencyFlag() {
0895 if (theHisto_all)
0896 theHisto_all->setEfficiencyFlag();
0897 if (mcPlots_) {
0898 if (mcPlots_ > 2) {
0899 theHisto_d->setEfficiencyFlag();
0900 theHisto_u->setEfficiencyFlag();
0901 theHisto_s->setEfficiencyFlag();
0902 theHisto_g->setEfficiencyFlag();
0903 theHisto_dus->setEfficiencyFlag();
0904 }
0905 theHisto_c->setEfficiencyFlag();
0906 theHisto_b->setEfficiencyFlag();
0907 theHisto_ni->setEfficiencyFlag();
0908 theHisto_dusg->setEfficiencyFlag();
0909 theHisto_pu->setEfficiencyFlag();
0910 }
0911 }
0912
0913 template <class T, class G>
0914 void FlavourHistograms2D<T, G>::fillVariable(const int &flavour, const T &varX, const G &varY, const float &w) const {
0915
0916 if (theHisto_all)
0917 theHisto_all->Fill(varX, varY, w);
0918 if (createProfile_)
0919
0920 if (theProfile_all)
0921 theProfile_all->Fill(varX, varY);
0922
0923
0924
0925 if (!mcPlots_)
0926 return;
0927
0928 switch (flavour) {
0929 case 1:
0930 if (mcPlots_ > 2) {
0931 theHisto_d->Fill(varX, varY, w);
0932 theHisto_dus->Fill(varX, varY, w);
0933 }
0934 theHisto_dusg->Fill(varX, varY, w);
0935 if (createProfile_) {
0936
0937
0938
0939 if (mcPlots_ > 2) {
0940 theProfile_d->Fill(varX, varY);
0941 theProfile_dus->Fill(varX, varY);
0942 }
0943 theProfile_dusg->Fill(varX, varY);
0944 }
0945 return;
0946 case 2:
0947 if (mcPlots_ > 2) {
0948 theHisto_u->Fill(varX, varY, w);
0949 theHisto_dus->Fill(varX, varY, w);
0950 }
0951 theHisto_dusg->Fill(varX, varY, w);
0952 if (createProfile_) {
0953
0954
0955
0956 if (mcPlots_ > 2) {
0957 theProfile_u->Fill(varX, varY);
0958 theProfile_dus->Fill(varX, varY);
0959 }
0960 theProfile_dusg->Fill(varX, varY);
0961 }
0962 return;
0963 case 3:
0964 if (mcPlots_ > 2) {
0965 theHisto_s->Fill(varX, varY, w);
0966 theHisto_dus->Fill(varX, varY, w);
0967 }
0968 theHisto_dusg->Fill(varX, varY, w);
0969 if (createProfile_) {
0970
0971
0972
0973 if (mcPlots_ > 2) {
0974 theProfile_s->Fill(varX, varY);
0975 theProfile_dus->Fill(varX, varY);
0976 }
0977 theProfile_dusg->Fill(varX, varY);
0978 }
0979 return;
0980 case 4:
0981 theHisto_c->Fill(varX, varY, w);
0982
0983 if (createProfile_)
0984 theProfile_c->Fill(varX, varY);
0985 return;
0986 case 5:
0987 theHisto_b->Fill(varX, varY, w);
0988
0989 if (createProfile_)
0990 theProfile_b->Fill(varX, varY);
0991 return;
0992 case 21:
0993 if (mcPlots_ > 2)
0994 theHisto_g->Fill(varX, varY, w);
0995 theHisto_dusg->Fill(varX, varY, w);
0996 if (createProfile_) {
0997
0998
0999 if (mcPlots_ > 2)
1000 theProfile_g->Fill(varX, varY);
1001 theProfile_dusg->Fill(varX, varY);
1002 }
1003 return;
1004 case 20:
1005 theHisto_pu->Fill(varX, varY, w);
1006
1007 if (createProfile_)
1008 theProfile_pu->Fill(varX, varY);
1009 return;
1010 default:
1011 theHisto_ni->Fill(varX, varY, w);
1012
1013 if (createProfile_)
1014 theProfile_ni->Fill(varX, varY);
1015 return;
1016 }
1017 }
1018
1019 template <class T, class G>
1020 std::vector<TH2F *> FlavourHistograms2D<T, G>::getHistoVector() const {
1021 std::vector<TH2F *> histoVector;
1022 if (theHisto_all)
1023 histoVector.push_back(theHisto_all->getTH2F());
1024 if (mcPlots_) {
1025 if (mcPlots_ > 2) {
1026 histoVector.push_back(theHisto_d->getTH2F());
1027 histoVector.push_back(theHisto_u->getTH2F());
1028 histoVector.push_back(theHisto_s->getTH2F());
1029 histoVector.push_back(theHisto_g->getTH2F());
1030 histoVector.push_back(theHisto_dus->getTH2F());
1031 }
1032 histoVector.push_back(theHisto_c->getTH2F());
1033 histoVector.push_back(theHisto_b->getTH2F());
1034 histoVector.push_back(theHisto_ni->getTH2F());
1035 histoVector.push_back(theHisto_dusg->getTH2F());
1036 histoVector.push_back(theHisto_pu->getTH2F());
1037 }
1038 return histoVector;
1039 }
1040
1041 template <class T, class G>
1042 std::vector<TProfile *> FlavourHistograms2D<T, G>::getProfileVector() const {
1043 std::vector<TProfile *> profileVector;
1044 if (createProfile_) {
1045 if (theProfile_all)
1046 profileVector.push_back(theProfile_all->getTProfile());
1047 if (mcPlots_) {
1048 if (mcPlots_ > 2) {
1049 profileVector.push_back(theProfile_d->getTProfile());
1050 profileVector.push_back(theProfile_u->getTProfile());
1051 profileVector.push_back(theProfile_s->getTProfile());
1052 profileVector.push_back(theProfile_g->getTProfile());
1053 profileVector.push_back(theProfile_dus->getTProfile());
1054 }
1055 profileVector.push_back(theProfile_c->getTProfile());
1056 profileVector.push_back(theProfile_b->getTProfile());
1057 profileVector.push_back(theProfile_ni->getTProfile());
1058 profileVector.push_back(theProfile_dusg->getTProfile());
1059 profileVector.push_back(theProfile_pu->getTProfile());
1060 }
1061 }
1062 return profileVector;
1063 }
1064
1065 #endif