File indexing completed on 2024-04-06 12:09:48
0001 #include "DQMOffline/RecoB/interface/JetTagPlotter.h"
0002 #include "DQMOffline/RecoB/interface/Tools.h"
0003 #include "FWCore/Utilities/interface/isFinite.h"
0004 #include "DQMServices/Core/interface/DQMStore.h"
0005
0006 #include <iostream>
0007
0008 using namespace std;
0009 using namespace RecoBTag;
0010
0011 JetTagPlotter::JetTagPlotter(const std::string& tagName,
0012 const EtaPtBin& etaPtBin,
0013 const edm::ParameterSet& pSet,
0014 unsigned int mc,
0015 bool wf,
0016 DQMStore::IBooker& ibook,
0017 bool doCTagPlots ,
0018 bool doDifferentialPlots ,
0019 double discrCut )
0020 : BaseBTagPlotter(tagName, etaPtBin),
0021 discrStart_(pSet.getParameter<double>("discriminatorStart")),
0022 discrEnd_(pSet.getParameter<double>("discriminatorEnd")),
0023 nBinEffPur_(pSet.getParameter<int>("nBinEffPur")),
0024 startEffPur_(pSet.getParameter<double>("startEffPur")),
0025 endEffPur_(pSet.getParameter<double>("endEffPur")),
0026 mcPlots_(mc),
0027 willFinalize_(wf),
0028 doCTagPlots_(doCTagPlots),
0029 doDifferentialPlots_(doDifferentialPlots),
0030 cutValue_(discrCut) {
0031
0032 const std::string& es = theExtensionString;
0033 const std::string jetTagDir(es.substr(1));
0034
0035 if (willFinalize_)
0036 return;
0037
0038
0039 int nFl = 1;
0040 if (mcPlots_)
0041 nFl = 8;
0042 nJets_.resize(nFl, 0);
0043
0044 if (mcPlots_) {
0045
0046 dJetFlav_ = std::make_unique<FlavourHistograms<int>>(
0047 "jetFlavour" + es, "Jet Flavour", 22, -0.5, 21.5, false, false, false, "b", jetTagDir, mcPlots_, ibook);
0048 }
0049
0050
0051 jetMultiplicity_ = std::make_unique<FlavourHistograms<int>>(
0052 "jetMultiplicity" + es, "Jet Multiplicity", 11, -0.5, 10.5, false, true, true, "b", jetTagDir, mcPlots_, ibook);
0053
0054
0055 dDiscriminator_ = std::make_unique<FlavourHistograms<double>>(
0056 "discr" + es, "Discriminator", 102, discrStart_, discrEnd_, false, true, true, "b", jetTagDir, mcPlots_, ibook);
0057 dDiscriminator_->settitle("Discriminant");
0058
0059 dJetRecMomentum_ = std::make_unique<FlavourHistograms<double>>(
0060 "jetMomentum" + es, "jet momentum", 350, 0.0, 350.0, false, false, true, "b", jetTagDir, mcPlots_, ibook);
0061
0062
0063 dJetRecPt_ = std::make_unique<FlavourHistograms<double>>(
0064 "jetPt" + es, "jet pt", 350, 0.0, 350.0, false, false, true, "b", jetTagDir, mcPlots_, ibook);
0065
0066
0067 dJetRecPseudoRapidity_ = std::make_unique<FlavourHistograms<double>>("jetEta" + es,
0068 "jet eta",
0069 20,
0070 -etaPtBin.getEtaMax(),
0071 etaPtBin.getEtaMax(),
0072 false,
0073 false,
0074 true,
0075 "b",
0076 jetTagDir,
0077 mcPlots_,
0078 ibook);
0079
0080
0081 dJetRecPhi_ = std::make_unique<FlavourHistograms<double>>(
0082 "jetPhi" + es, "jet phi", 20, -M_PI, M_PI, false, false, true, "b", jetTagDir, mcPlots_, ibook);
0083
0084 if (doDifferentialPlots_) {
0085
0086 dJetPhiDiscrCut_ = std::make_unique<FlavourHistograms<double>>("jetPhi_diffEff" + es,
0087 "Efficiency vs. jet Phi for discriminator above cut",
0088 20,
0089 -M_PI,
0090 M_PI,
0091 false,
0092 false,
0093 true,
0094 "b",
0095 jetTagDir,
0096 mcPlots_,
0097 ibook);
0098
0099
0100 dJetPseudoRapidityDiscrCut_ =
0101 std::make_unique<FlavourHistograms<double>>("jetEta_diffEff" + es,
0102 "Efficiency vs. jet eta for discriminator above cut",
0103 20,
0104 -etaPtBin.getEtaMax(),
0105 etaPtBin.getEtaMax(),
0106 false,
0107 false,
0108 true,
0109 "b",
0110 jetTagDir,
0111 mcPlots_,
0112 ibook);
0113 }
0114 }
0115
0116 JetTagPlotter::~JetTagPlotter() {}
0117
0118 void JetTagPlotter::epsPlot(const std::string& name) {
0119 if (!willFinalize_) {
0120 dJetFlav_->epsPlot(name);
0121 jetMultiplicity_->epsPlot(name);
0122 dDiscriminator_->epsPlot(name);
0123 dJetRecMomentum_->epsPlot(name);
0124 dJetRecPt_->epsPlot(name);
0125 dJetRecPseudoRapidity_->epsPlot(name);
0126 dJetRecPhi_->epsPlot(name);
0127 } else {
0128 effPurFromHistos_->epsPlot(name);
0129 }
0130 }
0131
0132 void JetTagPlotter::psPlot(const std::string& name) {
0133 std::string cName = "JetTagPlots" + theExtensionString;
0134 setTDRStyle()->cd();
0135 TCanvas canvas(cName.c_str(), cName.c_str(), 600, 900);
0136 canvas.UseCurrentStyle();
0137
0138 canvas.Divide(2, 3);
0139 canvas.Print((name + cName + ".ps[").c_str());
0140 if (!willFinalize_) {
0141 canvas.cd(1);
0142 dJetFlav_->plot();
0143 canvas.cd(2);
0144 canvas.cd(3);
0145 dDiscriminator_->plot();
0146 canvas.cd(4);
0147 dJetRecMomentum_->plot();
0148 canvas.cd(5);
0149 dJetRecPt_->plot();
0150 canvas.cd(6);
0151 dJetRecPseudoRapidity_->plot();
0152 canvas.Print((name + cName + ".ps").c_str());
0153 canvas.Clear();
0154 canvas.Divide(2, 3);
0155
0156 jetMultiplicity_->plot();
0157 canvas.Print((name + cName + ".ps").c_str());
0158 canvas.Clear();
0159
0160 canvas.cd(1);
0161 dJetRecPhi_->plot();
0162 canvas.cd(2);
0163 canvas.cd(3);
0164 canvas.cd(4);
0165 } else {
0166 canvas.cd(5);
0167 effPurFromHistos_->discriminatorNoCutEffic().plot();
0168 canvas.cd(6);
0169 effPurFromHistos_->discriminatorCutEfficScan().plot();
0170 canvas.Print((name + cName + ".ps").c_str());
0171 canvas.Clear();
0172 canvas.Divide(2, 3);
0173 canvas.cd(1);
0174 effPurFromHistos_->plot();
0175 }
0176 canvas.Print((name + cName + ".ps").c_str());
0177 canvas.Print((name + cName + ".ps]").c_str());
0178 }
0179
0180 void JetTagPlotter::analyzeTag()
0181 {
0182
0183 jetMultiplicity_->fill(-1, nJets_[0]);
0184 nJets_[0] = 0;
0185 }
0186
0187 void JetTagPlotter::analyzeTag(float w) {
0188 if (mcPlots_) {
0189
0190 int totNJets = 0;
0191 int udsNJets = 0;
0192 int udsgNJets = 0;
0193 for (int i = 0; i < 8; i++) {
0194 totNJets += nJets_[i];
0195 if (i > 0 && i < 4)
0196 udsNJets += nJets_[i];
0197 if ((i > 0 && i < 4) || i == 6)
0198 udsgNJets += nJets_[i];
0199 if (i <= 5 && i >= 1)
0200 jetMultiplicity_->fill(i, nJets_[i], w);
0201 else if (i == 6)
0202 jetMultiplicity_->fill(21, nJets_[i], w);
0203 else if (i == 7)
0204 jetMultiplicity_->fill(20, nJets_[i], w);
0205 else
0206 jetMultiplicity_->fill(0, nJets_[i], w);
0207 nJets_[i] = 0;
0208 }
0209 jetMultiplicity_->fill(-1, totNJets, w);
0210 jetMultiplicity_->fill(123, udsNJets, w);
0211 jetMultiplicity_->fill(12321, udsgNJets, w);
0212 } else {
0213 int totNJets = 0;
0214 for (int i = 0; i < 8; i++) {
0215 totNJets += nJets_[i];
0216 nJets_[i] = 0;
0217 }
0218 jetMultiplicity_->fill(-1, totNJets, w);
0219 }
0220 }
0221
0222 void JetTagPlotter::analyzeTag(const reco::Jet& jet, double jec, float discriminator, int jetFlavour, float w ) {
0223 if (mcPlots_) {
0224 dJetFlav_->fill(jetFlavour, jetFlavour, w);
0225 if (abs(jetFlavour) > 0 && abs(jetFlavour) < 6)
0226 nJets_[abs(jetFlavour)] += 1;
0227 else if (abs(jetFlavour) == 21)
0228 nJets_[6] += 1;
0229 else if (jetFlavour == 20)
0230 nJets_[7] += 1;
0231 else
0232 nJets_[0] += 1;
0233 } else {
0234 nJets_[0] += 1;
0235 }
0236 if (edm::isNotFinite(discriminator))
0237 dDiscriminator_->fill(jetFlavour, -999.0, w);
0238 else
0239 dDiscriminator_->fill(jetFlavour, discriminator, w);
0240 dJetRecMomentum_->fill(jetFlavour, jet.p() * jec, w);
0241 dJetRecPt_->fill(jetFlavour, jet.pt() * jec, w);
0242 dJetRecPseudoRapidity_->fill(jetFlavour, jet.eta(), w);
0243 dJetRecPhi_->fill(jetFlavour, jet.phi(), w);
0244 if (doDifferentialPlots_) {
0245 if (edm::isFinite(discriminator) && discriminator > cutValue_) {
0246 dJetPhiDiscrCut_->fill(jetFlavour, jet.phi(), w);
0247 dJetPseudoRapidityDiscrCut_->fill(jetFlavour, jet.eta(), w);
0248 }
0249 }
0250 }
0251
0252 void JetTagPlotter::analyzeTag(const reco::JetTag& jetTag, double jec, int jetFlavour, float w ) {
0253 if (mcPlots_) {
0254 dJetFlav_->fill(jetFlavour, jetFlavour, w);
0255 if (abs(jetFlavour) > 0 && abs(jetFlavour) < 6)
0256 nJets_[abs(jetFlavour)] += 1;
0257 else if (abs(jetFlavour) == 21)
0258 nJets_[6] += 1;
0259 else if (jetFlavour == 20)
0260 nJets_[7] += 1;
0261 else
0262 nJets_[0] += 1;
0263 } else {
0264 nJets_[0] += 1;
0265 }
0266 const auto& discriminator = jetTag.second;
0267 if (edm::isNotFinite(discriminator))
0268 dDiscriminator_->fill(jetFlavour, -999.0, w);
0269 else
0270 dDiscriminator_->fill(jetFlavour, discriminator, w);
0271 dJetRecMomentum_->fill(jetFlavour, jetTag.first->p() * jec, w);
0272 dJetRecPt_->fill(jetFlavour, jetTag.first->pt() * jec, w);
0273 dJetRecPseudoRapidity_->fill(jetFlavour, jetTag.first->eta(), w);
0274 dJetRecPhi_->fill(jetFlavour, jetTag.first->phi(), w);
0275 if (doDifferentialPlots_) {
0276 if (edm::isFinite(discriminator) && discriminator > cutValue_) {
0277 dJetPhiDiscrCut_->fill(jetFlavour, jetTag.first->phi(), w);
0278 dJetPseudoRapidityDiscrCut_->fill(jetFlavour, jetTag.first->eta(), w);
0279 }
0280 }
0281 }
0282
0283 void JetTagPlotter::finalize(DQMStore::IBooker& ibook_, DQMStore::IGetter& igetter_) {
0284
0285
0286
0287
0288 const std::string& es = theExtensionString;
0289 const std::string jetTagDir(es.substr(1));
0290 dDiscriminator_ = std::make_unique<FlavourHistograms<double>>(
0291 "discr" + es, "Discriminator", 102, discrStart_, discrEnd_, "b", jetTagDir, mcPlots_, igetter_);
0292
0293 effPurFromHistos_ = std::make_unique<EffPurFromHistos>(
0294 *dDiscriminator_, jetTagDir, mcPlots_, ibook_, nBinEffPur_, startEffPur_, endEffPur_);
0295 effPurFromHistos_->doCTagPlots(doCTagPlots_);
0296 effPurFromHistos_->compute(ibook_);
0297
0298
0299 if (doDifferentialPlots_) {
0300 dJetRecPhi_ = std::make_unique<FlavourHistograms<double>>(
0301 "jetPhi" + es, "jet phi", 20, -M_PI, M_PI, "b", jetTagDir, mcPlots_, igetter_);
0302 dJetPhiDiscrCut_ = std::make_unique<FlavourHistograms<double>>("jetPhi_diffEff" + es,
0303 "Efficiency vs. jet Phi for discriminator above cut",
0304 20,
0305 -M_PI,
0306 M_PI,
0307 "b",
0308 jetTagDir,
0309 mcPlots_,
0310 igetter_);
0311 dJetPhiDiscrCut_->divide(*dJetRecPhi_);
0312 dJetPhiDiscrCut_->setEfficiencyFlag();
0313
0314 dJetRecPseudoRapidity_ = std::make_unique<FlavourHistograms<double>>(
0315 "jetEta" + es, "jet eta", 20, -etaPtBin_.getEtaMax(), etaPtBin_.getEtaMax(), "b", jetTagDir, mcPlots_, igetter_);
0316 dJetPseudoRapidityDiscrCut_ =
0317 std::make_unique<FlavourHistograms<double>>("jetEta_diffEff" + es,
0318 "Efficiency vs. jet eta for discriminator above cut",
0319 20,
0320 -etaPtBin_.getEtaMax(),
0321 etaPtBin_.getEtaMax(),
0322 "b",
0323 jetTagDir,
0324 mcPlots_,
0325 igetter_);
0326 dJetPseudoRapidityDiscrCut_->divide(*dJetRecPseudoRapidity_);
0327 dJetPseudoRapidityDiscrCut_->setEfficiencyFlag();
0328 }
0329 }