From 705da18b3d281e412781846b508ba78821bc543b Mon Sep 17 00:00:00 2001 From: capalmer Date: Wed, 3 Apr 2019 22:19:51 +0200 Subject: [PATCH 1/7] include first version of the 2TagCount code --- test/ttbar/TTbarEventAnalysis.cc | 1502 +++++++++++++++++--------- test/ttbar/TTbarEventAnalysis.h | 79 +- test/ttbar/data/samples_Run2018.json | 20 + test/ttbar/plotter.py | 43 +- test/ttbar/runPileupEstimation.py | 7 +- test/ttbar/runTTbarAnalysis.py | 101 +- test/ttbar/storeTools.py | 56 +- test/ttbar/twoTag.py | 119 ++ 8 files changed, 1363 insertions(+), 564 deletions(-) create mode 100644 test/ttbar/data/samples_Run2018.json create mode 100644 test/ttbar/twoTag.py diff --git a/test/ttbar/TTbarEventAnalysis.cc b/test/ttbar/TTbarEventAnalysis.cc index 7da61a846eb..b4a6444c87a 100644 --- a/test/ttbar/TTbarEventAnalysis.cc +++ b/test/ttbar/TTbarEventAnalysis.cc @@ -55,21 +55,61 @@ void TTbarEventAnalysis::prepareOutput(TString outFile) ftmTree_->Branch("kindisc", kinDisc_, "kindisc[2]/F"); ftmTree_->Branch("weight", weight_, "weight[15]/F"); + wpLabel.push_back("L"); + wpLabel.push_back("M"); + wpLabel.push_back("T"); + + systName.push_back("qcdScaleLo"); + systName.push_back("qcdScaleHi"); + systName.push_back("isrRedHi"); + systName.push_back("fsrRedHi"); + systName.push_back("isrRedLo"); + systName.push_back("fsrRedLo"); + systName.push_back("isrDefHi"); + systName.push_back("fsrDefHi"); + systName.push_back("isrDefLo"); + systName.push_back("fsrDefLo"); + systName.push_back("isrConHi"); + systName.push_back("fsrConHi"); + systName.push_back("isrConLo"); + systName.push_back("fsrConLo"); + + + for(unsigned int iSyst=0; iSyst baseHistos; - baseHistos["npvinc" ] = new TH1F("npvinc", ";N_{PV,good}-N_{HS};Events", 50, 0, 50); - baseHistos["npv" ] = new TH1F("npv", ";N_{PV,good}-N_{HS};Events", 50, 0, 50); - baseHistos["rho" ] = new TH1F("rho", ";#rho [GeV];Events", 50, 0, 30); - baseHistos["mll" ] = new TH1F("mll", ";Dilepton invariant mass [GeV];Events", 20, 0, 300); - baseHistos["mllinc" ] = new TH1F("mllinc", ";Dilepton invariant mass [GeV];Events", 20, 0, 300); + std::map baseHistos2d; + baseHistos["npvinc" ] = new TH1F("npvinc", ";N_{PV,good}-N_{HS};Events", 100, 0, 100); + baseHistos["rhoinc" ] = new TH1F("rhoinc", ";#rho [GeV];Events", 100, 0, 100); + baseHistos["npv" ] = new TH1F("npv", ";N_{PV,good}-N_{HS};Events", 100, 0, 100); + baseHistos["rho" ] = new TH1F("rho", ";#rho [GeV];Events", 100, 0, 100); + baseHistos["mll" ] = new TH1F("mll", ";Dilepton invariant mass [GeV];Events", 30, 0, 300); + baseHistos["mllinc" ] = new TH1F("mllinc", ";Dilepton invariant mass [GeV];Events", 30, 0, 300); baseHistos["met" ] = new TH1F("met", ";Missing transverse energy [GeV];Events", 15, 0, 300); baseHistos["njets" ] = new TH1F("njets", ";Jet multiplicity;Events;", 6, 2, 8); - baseHistos["leadjpt"] = new TH1F("leadjpt",";Leading jet p_{T} [GeV];Events;", 14,30,300); - baseHistos["leadlpt"] = new TH1F("leadlpt",";Leading lepton p_{T} [GeV];Events;", 9,20,200); - baseHistos["trailjpt"] = new TH1F("trailjpt",";Trailing jet p_{T} [GeV];Events;", 14,30,300); - baseHistos["traillpt"] = new TH1F("traillpt",";Trailing lepton p_{T} [GeV];Events;", 9,20,200); - baseHistos["leadjeta"] = new TH1F("leadjeta", ";Pseudo-rapidity; Jets", 25, 0, 2.5); - baseHistos["trailjeta"] = new TH1F("trailjeta", ";Pseudo-rapidity; Jets", 25, 0, 2.5); + baseHistos["leadjpt"] = new TH1F("leadjpt",";Leading jet p_{T} [GeV];Events;", 37,30,400); + baseHistos["leadbjpt"] = new TH1F("leadbjpt",";Leading b-jet p_{T} [GeV];Events;", 37,30,400); + baseHistos["only2_leadjpt"] = new TH1F("only2_leadjpt",";Leading jet p_{T} [GeV];Events;", 37,30,400); + baseHistos["only2_leadbjpt"] = new TH1F("only2_leadbjpt",";Leading b-jet p_{T} [GeV];Events;", 37,30,400); + baseHistos["leadlpt"] = new TH1F("leadlpt",";Leading lepton p_{T} [GeV];Events;", 18,20,200); + baseHistos["leadleta"] = new TH1F("leadleta",";Leading lepton #eta;Events;", 50,-3,3); + baseHistos["trailjpt"] = new TH1F("trailjpt",";Trailing jet p_{T} [GeV];Events;", 37,30,400); + baseHistos["trailbjpt"] = new TH1F("trailbjpt",";Trailing b-jet p_{T} [GeV];Events;", 37,30,400); + baseHistos["only2_trailjpt"] = new TH1F("only2_trailjpt",";Trailing jet p_{T} [GeV];Events;", 37,30,400); + baseHistos["only2_trailbjpt"] = new TH1F("only2_trailbjpt",";Trailing b-jet p_{T} [GeV];Events;", 37,30,400); + baseHistos["traillpt"] = new TH1F("traillpt",";Trailing lepton p_{T} [GeV];Events;", 18,20,200); + baseHistos["trailleta"] = new TH1F("trailleta",";Trailing lepton #eta;Events;", 50,-3,3); + baseHistos["leadjeta"] = new TH1F("leadjeta", ";Pseudo-rapidity; Jets", 35, 0, 3.5); + baseHistos["leadbjeta"] = new TH1F("leadbjeta", ";Pseudo-rapidity; Jets", 35, 0, 3.5); + baseHistos["trailjeta"] = new TH1F("trailjeta", ";Pseudo-rapidity; Jets", 35, 0, 3.5); + baseHistos["trailbjeta"] = new TH1F("trailbjeta", ";Pseudo-rapidity; Jets", 35, 0, 3.5); + baseHistos["only2_leadjeta"] = new TH1F("only2_leadjeta", ";Pseudo-rapidity; Jets", 35, 0, 3.5); + baseHistos["only2_leadbjeta"] = new TH1F("only2_leadbjeta", ";Pseudo-rapidity; Jets", 35, 0, 3.5); + baseHistos["only2_trailjeta"] = new TH1F("only2_trailjeta", ";Pseudo-rapidity; Jets", 35, 0, 3.5); + baseHistos["only2_trailbjeta"] = new TH1F("only2_trailbjeta", ";Pseudo-rapidity; Jets", 35, 0, 3.5); baseHistos["evsel"] = new TH1F("evsel", ";Event selection;Events;", 4,0,4); baseHistos["evsel"]->GetXaxis()->SetBinLabel(1,"#geq 2j"); baseHistos["evsel"]->GetXaxis()->SetBinLabel(2,"=2j"); @@ -77,20 +117,72 @@ void TTbarEventAnalysis::prepareOutput(TString outFile) baseHistos["evsel"]->GetXaxis()->SetBinLabel(4,"=4j"); baseHistos["jp"]=new TH1F("jp",";Jet probability;Jets",50,0,3); baseHistos["svhe"]=new TH1F("svhe",";Simple secondary vertex (HE);Jets",50,0,6); - baseHistos["csv"]=new TH1F("csv",";Combined secondary vertex (IVF);Jets",50,0,1.1); + baseHistos["csv"]=new TH1F("csv",";Combined secondary vertex (IVF);Jets",52,-0.02,1.02); + baseHistos["deepcsv"]=new TH1F("deepcsv",";DeepCSV for B;Jets",52,-0.02,1.02); + + std::vector twoTagNames; + twoTagNames.clear(); + twoTagNames.push_back("twoTags_deepCSV"); + twoTagNames.push_back("only2_twoTags_deepCSV"); + twoTagNames.push_back("twoTags_deepFlavour"); + twoTagNames.push_back("only2_twoTags_deepFlavour"); + twoTagNames.push_back("twoTags_CSVv2"); + twoTagNames.push_back("only2_twoTags_CSVv2"); + //twoTagNames.push_back(""); + + for(unsigned int i2Tag=0; i2TagGetXaxis()->SetBinLabel(1,"unmatched"); baseHistos["flavour"]->GetXaxis()->SetBinLabel(2,"udsg"); baseHistos["flavour"]->GetXaxis()->SetBinLabel(3,"c"); baseHistos["flavour"]->GetXaxis()->SetBinLabel(4,"b"); - baseHistos["close_mlj"] = new TH1F("close_mlj", ";M(lepton,jet) [GeV]; Jets", 50, 0, 250); + baseHistos["close_mlj"] = new TH1F("close_mlj", ";M(lepton,jet) [GeV]; Jets", 100, 0,400); baseHistos["close_deta"] = new TH1F("close_deta", ";#Delta#eta(lepton,jet); Jets", 50, 0, 4); baseHistos["close_dphi"] = new TH1F("close_dphi", ";#Delta#phi(lepton,jet) [rad]; Jets", 50, 0, 3.15); baseHistos["close_ptrel"] = new TH1F("close_ptrel", ";p_{T}^{rel}(lepton,jet) [GeV];Jets", 50, 0, 1); @@ -111,12 +203,23 @@ void TTbarEventAnalysis::prepareOutput(TString outFile) for(size_t i=0; i::iterator it=baseHistos.begin(); it!=baseHistos.end(); it++) - { - TString tag=ch[i]+"_"+it->first; - histos_[tag]=(TH1F *)it->second->Clone(tag); - histos_[tag]->Sumw2(); - histos_[tag]->SetDirectory(outF_); - } + { + TString tag=ch[i]+"_"+it->first; + histos_[tag]=(TH1F *)it->second->Clone(tag); + histos_[tag]->Sumw2(); + histos_[tag]->SetDirectory(outF_); + } + } + + for(size_t i=0; i::iterator it=baseHistos2d.begin(); it!=baseHistos2d.end(); it++) + { + TString tag=ch[i]+"_"+it->first; + histos2d_[tag]=(TH2F *)it->second->Clone(tag); + histos2d_[tag]->Sumw2(); + histos2d_[tag]->SetDirectory(outF_); + } } histos_["puwgtnorm"] = new TH1F("puwgtnorm", ";puwgtnorm;Events", 4, 0, 4); @@ -127,479 +230,826 @@ void TTbarEventAnalysis::prepareOutput(TString outFile) // void TTbarEventAnalysis::processFile(TString inFile,TH1F *xsecWgt,Bool_t isData) { - //loop over events - TFile *inF=TFile::Open(inFile); - TTree *tree=(TTree *)inF->Get("btagana/ttree"); - Int_t nentries=tree->GetEntriesFast(); - std::cout << "...opening " << inFile << " -> analysing " << nentries << " events -> " << outF_->GetName(); - if(xsecWgt) std::cout << " xsec weight=" << xsecWgt->GetBinContent(1); - if(isData) std::cout << " is data"; - std::cout << std::endl; - - //prepare reader - std::vector tmvaVars( tmvaVarNames_.size(), 0. ); - if(weightsDir_!="") - { - tmvaReader_=new TMVA::Reader( "!Color:!Silent" ); - for(size_t ivar=0; ivarAddVariable( tmvaVarNames_[ivar], &tmvaVars[ivar] ); - - TString jranks[]={"leading", "others", "subleading" }; - for(size_t i=0; iBookMVA("BDT_"+jranks[i], weightsDir_+"/"+jranks[i]+"/TMVAClassification_BDT.weights.xml"); + //loop over events + TFile *inF=TFile::Open(inFile); + TTree *tree=(TTree *)inF->Get("btagana/ttree"); + Int_t nentries=tree->GetEntriesFast(); + std::cout << "...opening " << inFile << " -> analysing " << nentries << " events -> " << outF_->GetName(); + if(xsecWgt) std::cout << " xsec weight=" << xsecWgt->GetBinContent(1); + if(isData) std::cout << " is data"; + std::cout << std::endl; + + //prepare reader + std::vector tmvaVars( tmvaVarNames_.size(), 0. ); + if(weightsDir_!=""){ + tmvaReader_=new TMVA::Reader( "!Color:!Silent" ); + for(size_t ivar=0; ivarAddVariable( tmvaVarNames_[ivar], &tmvaVars[ivar] ); + + TString jranks[]={"leading", "others", "subleading" }; + for(size_t i=0; iBookMVA("BDT_"+jranks[i], weightsDir_+"/"+jranks[i]+"/TMVAClassification_BDT.weights.xml"); } - //prepare to read the tree (for jets only interested in a couple of variables) - struct MyEventInfoBranches_t - { - Int_t Run,Evt,LumiBlock,nPV; - Int_t ttbar_chan, ttbar_trigWord, ttbar_metfilterWord; - Int_t ttbar_nl, ttbar_lid[10], ttbar_lgid[10], ttbar_lch[10]; - Float_t ttbar_lpt[10], ttbar_leta[10], ttbar_lphi[10], ttbar_lm[10]; - Float_t ttbar_metpt,ttbar_metphi; - Float_t ttbar_rho; - Int_t ttbar_nw; - Float_t nPUtrue; - Float_t ttbar_w[500]; - Int_t nJet; - Float_t Jet_pt[100],Jet_genpt[100],Jet_area[100],Jet_jes[100],Jet_eta[100],Jet_phi[100],Jet_mass[100]; - Float_t Jet_Svx[100],Jet_CombIVF[100],Jet_Proba[100],Jet_Ip2P[100]; - Int_t Jet_nseltracks[100]; - Int_t Jet_flavour[100]; - }; - MyEventInfoBranches_t ev; - tree->SetBranchAddress("Run" , &ev.Run ); - tree->SetBranchAddress("Evt" , &ev.Evt ); - tree->SetBranchAddress("LumiBlock" , &ev.LumiBlock ); - tree->SetBranchAddress("nPV" , &ev.nPV ); - tree->SetBranchAddress("nPUtrue", &ev.nPUtrue ); - tree->SetBranchAddress("ttbar_chan" , &ev.ttbar_chan); - tree->SetBranchAddress("ttbar_metfilterWord", &ev.ttbar_metfilterWord); - tree->SetBranchAddress("ttbar_trigWord", &ev.ttbar_trigWord); - tree->SetBranchAddress("ttbar_nl" , &ev.ttbar_nl); - tree->SetBranchAddress("ttbar_lpt" , ev.ttbar_lpt); - tree->SetBranchAddress("ttbar_leta" , ev.ttbar_leta); - tree->SetBranchAddress("ttbar_lphi" , ev.ttbar_lphi); - tree->SetBranchAddress("ttbar_lm" , ev.ttbar_lm); - tree->SetBranchAddress("ttbar_lid" , ev.ttbar_lid); - tree->SetBranchAddress("ttbar_lgid" , ev.ttbar_lgid); - tree->SetBranchAddress("ttbar_lch" , ev.ttbar_lch); - tree->SetBranchAddress("ttbar_metpt", &ev.ttbar_metpt); - tree->SetBranchAddress("ttbar_metphi", &ev.ttbar_metphi); - tree->SetBranchAddress("ttbar_rho", &ev.ttbar_rho); - tree->SetBranchAddress("ttbar_nw", &ev.ttbar_nw); - tree->SetBranchAddress("ttbar_w", ev.ttbar_w); - tree->SetBranchAddress("nJet", &ev.nJet); - tree->SetBranchAddress("Jet_pt", ev.Jet_pt); - tree->SetBranchAddress("Jet_genpt", ev.Jet_genpt); - tree->SetBranchAddress("Jet_area", ev.Jet_area); - tree->SetBranchAddress("Jet_jes", ev.Jet_jes); - tree->SetBranchAddress("Jet_eta", ev.Jet_eta); - tree->SetBranchAddress("Jet_phi", ev.Jet_phi); - tree->SetBranchAddress("Jet_mass", ev.Jet_mass); - tree->SetBranchAddress("Jet_Svx", ev.Jet_Svx); - tree->SetBranchAddress("Jet_CombIVF", ev.Jet_CombIVF); - tree->SetBranchAddress("Jet_Proba", ev.Jet_Proba); - tree->SetBranchAddress("Jet_Ip2P", ev.Jet_Ip2P); - tree->SetBranchAddress("Jet_nseltracks", ev.Jet_nseltracks); - tree->SetBranchAddress("Jet_flavour", ev.Jet_flavour); - - for(Int_t i=0; iGetEntry(i); - - //progress bar - if(i%100==0) std::cout << "\r[ " << int(100.*i/nentries) << "/100 ] to completion" << std::flush; - - //generator level weights - Float_t genWgt=ev.ttbar_nw==0 ? 1.0 : ev.ttbar_w[0]; - Float_t qcdScaleLo(1.0),qcdScaleHi(1.0),hdampLo(1.0),hdampHi(1.0); - if(readTTJetsGenWeights_ && ev.ttbar_nw>17) - { - qcdScaleLo=ev.ttbar_w[9]*(xsecWgt->GetBinContent(10)/xsecWgt->GetBinContent(1)); - qcdScaleHi=ev.ttbar_w[5]*(xsecWgt->GetBinContent(6)/xsecWgt->GetBinContent(1)); - hdampLo=ev.ttbar_w[ev.ttbar_nw-17]*(xsecWgt->GetBinContent(ev.ttbar_nw-17+1)/xsecWgt->GetBinContent(1)); - hdampHi=ev.ttbar_w[ev.ttbar_nw-9]*(xsecWgt->GetBinContent(ev.ttbar_nw-9+1)/xsecWgt->GetBinContent(1)); - } + tree->SetBranchAddress("Run" , &ev.Run ); + tree->SetBranchAddress("Evt" , &ev.Evt ); + tree->SetBranchAddress("LumiBlock" , &ev.LumiBlock ); + tree->SetBranchAddress("nPV" , &ev.nPV ); + tree->SetBranchAddress("nPUtrue", &ev.nPUtrue ); + tree->SetBranchAddress("nPU", &ev.nPU ); + tree->SetBranchAddress("ttbar_chan" , &ev.ttbar_chan); + tree->SetBranchAddress("ttbar_metfilterWord", &ev.ttbar_metfilterWord); + tree->SetBranchAddress("ttbar_trigWord", &ev.ttbar_trigWord); + tree->SetBranchAddress("ttbar_nl" , &ev.ttbar_nl); + tree->SetBranchAddress("ttbar_lpt" , ev.ttbar_lpt); + tree->SetBranchAddress("ttbar_leta" , ev.ttbar_leta); + tree->SetBranchAddress("ttbar_lphi" , ev.ttbar_lphi); + tree->SetBranchAddress("ttbar_lm" , ev.ttbar_lm); + tree->SetBranchAddress("ttbar_lid" , ev.ttbar_lid); + tree->SetBranchAddress("ttbar_lgid" , ev.ttbar_lgid); + tree->SetBranchAddress("ttbar_lch" , ev.ttbar_lch); + tree->SetBranchAddress("ttbar_metpt", &ev.ttbar_metpt); + tree->SetBranchAddress("ttbar_metphi", &ev.ttbar_metphi); + tree->SetBranchAddress("ttbar_rho", &ev.ttbar_rho); + tree->SetBranchAddress("ttbar_nw", &ev.ttbar_nw); + tree->SetBranchAddress("ttbar_w", ev.ttbar_w); + tree->SetBranchAddress("nJet", &ev.nJet); + tree->SetBranchAddress("Jet_pt", ev.Jet_pt); + tree->SetBranchAddress("Jet_genpt", ev.Jet_genpt); + tree->SetBranchAddress("Jet_area", ev.Jet_area); + tree->SetBranchAddress("Jet_jes", ev.Jet_jes); + tree->SetBranchAddress("Jet_eta", ev.Jet_eta); + tree->SetBranchAddress("Jet_phi", ev.Jet_phi); + tree->SetBranchAddress("Jet_mass", ev.Jet_mass); + tree->SetBranchAddress("Jet_Svx", ev.Jet_Svx); + tree->SetBranchAddress("Jet_CombIVF", ev.Jet_CombIVF); + tree->SetBranchAddress("Jet_DeepCSVBDisc",ev.Jet_DeepCSVBDisc); + tree->SetBranchAddress("Jet_DeepFlavourBDisc",ev.Jet_DeepFlavourBDisc); + tree->SetBranchAddress("Jet_Proba", ev.Jet_Proba); + tree->SetBranchAddress("Jet_Ip2P", ev.Jet_Ip2P); + tree->SetBranchAddress("Jet_nseltracks", ev.Jet_nseltracks); + tree->SetBranchAddress("Jet_flavour", ev.Jet_flavour); + + int nSkipped=0; + int n1=0; + int n2=0; + int n3=0; + for(Int_t i=0; iGetEntry(i); + + //if(isData){ + // if(ev.Run>304671){ //EF + // nSkipped++; + // continue; + // } + //} + + //progress bar + if(i%100==0) std::cout << "\r[ " << int(100.*i/nentries) << "/100 ] to completion" << std::flush; + + //generator level weights + Float_t genWgt=ev.ttbar_nw==0 ? 1.0 : ev.ttbar_w[0]; + //std::cout<<"go go 001"<17) { + // Weight * [sample xsec from json / sum(weights for given systematic)] / [sample xsec from json / sum(weights for nominal)] + // = Weight * sum(weights nominal) / sum(weights for given systematic) + // i.e. renormalise to same number of events as nominal. + // N.B. for GetBinContent(X) 'X' must be weight index + 1 due to way histogram is filled (0th bin always underflow). + systWeight["qcdScaleLo"]=ev.ttbar_w[9]*(xsecWgt->GetBinContent(10)/xsecWgt->GetBinContent(1)); + systWeight["qcdScaleHi"]=ev.ttbar_w[5]*(xsecWgt->GetBinContent(6)/xsecWgt->GetBinContent(1)); + systWeight["hdampLo"]=ev.ttbar_w[ev.ttbar_nw-29]*(xsecWgt->GetBinContent(ev.ttbar_nw-29+1)/xsecWgt->GetBinContent(1)); + systWeight["hdampHi"]=ev.ttbar_w[ev.ttbar_nw-21]*(xsecWgt->GetBinContent(ev.ttbar_nw-21+1)/xsecWgt->GetBinContent(1)); + + // >>> PSWeights <<< + // Vector of weight to be used instead of old ISR/FSR varied alternative samples. + // First weight (weightID= 1081) corresponds to central ME weight value. + // The remaining 12 values (weightIDs = 1082 to 1093) correspond to the PS weights in the following order (ISR up, FSR up, ISR down, FSR down) x 3 sets, i.e. + // 1082 = isrRedHi isr:muRfac=0.707, 1083 = fsrRedHi fsr:muRfac=0.707, 1084 = isrRedLo isr:muRfac=1.414, 1085 = fsrRedLo fsr:muRfac=1.414, + // 1086 = isrDefHi isr:muRfac=0.5, 1087 = fsrDefHi fsr:muRfac=0.5, 1088 = isrDefLo isr:muRfac=2.0, 1089 = fsrDefLo fsr:muRfac=2.0, + // 1090 = isrConHi isr:muRfac=0.25, 1091 = fsrConHi fsr:muRfac=0.25, 1092 = isrConLo isr:muRfac=4.0, 1093 = fsrConLo fsr:muRfac=4.0 - //pileup weights - Float_t puWgtLo(1.0), puWgtNom(1.0), puWgtHi(1.0); - if(!isData) - { - if(puWgtGr_) puWgtNom = puWgtGr_->Eval(ev.nPUtrue); - if(puWgtDownGr_) puWgtLo = puWgtDownGr_->Eval(ev.nPUtrue); - if(puWgtUpGr_) puWgtHi = puWgtUpGr_->Eval(ev.nPUtrue); - } - histos_["puwgtnorm" ]->Fill(0.,1.0); - histos_["puwgtnorm" ]->Fill(1.,puWgtNom); - histos_["puwgtnorm" ]->Fill(2.,puWgtLo); - histos_["puwgtnorm" ]->Fill(3.,puWgtHi); + + /*std::cout << "--- ttbar_nw " << ev.ttbar_nw << " ---" << std::endl; + cout << "--- number of bins in xsecWgt: " << xsecWgt.GetNbinsX() << "---" << endl; + std::cout << "Nominal sum of weights (xsecWgt) = " << xsecWgt->GetBinContent(1) << std::endl; + std::cout << "qcdScaleLo : " << ev.ttbar_w[9]*(xsecWgt->GetBinContent(10)/xsecWgt->GetBinContent(1)) << std::endl; + std::cout << "qcdScaleHi : " << ev.ttbar_w[5]*(xsecWgt->GetBinContent(6)/xsecWgt->GetBinContent(1)) << std::endl; + std::cout << "hdampLo : " << ev.ttbar_w[ev.ttbar_nw-31]*(xsecWgt->GetBinContent(ev.ttbar_nw-31+1)/xsecWgt->GetBinContent(1)) << std::endl; + std::cout << "hdampHi : " << ev.ttbar_w[ev.ttbar_nw-23]*(xsecWgt->GetBinContent(ev.ttbar_nw-23+1)/xsecWgt->GetBinContent(1)) << std::endl; + // Store generator weights + std::cout << "~~~~~~~~~~~~~~~~~~~~~~~~~~~" << std::endl; + std::cout << "isrRedHi event weight = " << ev.ttbar_w[1082] << std::endl; + std::cout << "isrRedHi sum of event weights = " << xsecWgt->GetBinContent(1083) << std::endl; + std::cout << "isrRedHi renormalisation factor w.r.t. nominal = " << xsecWgt->GetBinContent(1083)/xsecWgt->GetBinContent(1) << std::endl; + std::cout << "Stored isrRedHi weight = " << ev.ttbar_w[1082]*(xsecWgt->GetBinContent(1083)/xsecWgt->GetBinContent(1)) << std::endl; + std::cout << "~~~~~~~~~~~~~~~~~~~~~~~~~~~" << std::endl; + std::cout << "isrRedLo event weight = " << ev.ttbar_w[1084] << std::endl; + std::cout << "isrRedLo sum of event weights = " << xsecWgt->GetBinContent(1083) << std::endl; + std::cout << "isrRedLo renormalisation factor w.r.t. nominal = " << xsecWgt->GetBinContent(1085)/xsecWgt->GetBinContent(1) << std::endl; + std::cout << "Stored isrRedLo weight = " << ev.ttbar_w[1084]*xsecWgt->GetBinContent(1085)/xsecWgt->GetBinContent(1) << std::endl; + std::cout << "~~~~~~~~~~~~~~~~~~~~~~~~~~~" << std::endl; + std::cout << "fsrRedHi event weight = " << ev.ttbar_w[1083] << std::endl; + std::cout << "fsrRedHi sum of event weights = " << xsecWgt->GetBinContent(1084) << std::endl; + std::cout << "fsrRedHi renormalisation factor w.r.t. nominal = " << xsecWgt->GetBinContent(1084)/xsecWgt->GetBinContent(1) << std::endl; + std::cout << "Stored fsrRedHi weight = " << ev.ttbar_w[1083]*(xsecWgt->GetBinContent(1084)/xsecWgt->GetBinContent(1)) << std::endl; + std::cout << "~~~~~~~~~~~~~~~~~~~~~~~~~~~" << std::endl; + std::cout << "fsrRedLo event weight = " << ev.ttbar_w[1085] << std::endl; + std::cout << "fsrRedLo sum of event weights = " << xsecWgt->GetBinContent(1086) << std::endl; + std::cout << "fsrRedLo renormaliation factor w.r.t. nominal = " << xsecWgt->GetBinContent(1086)/xsecWgt->GetBinContent(1) << std::endl; + std::cout << "Stored fsrRedLo weight = " << ev.ttbar_w[1085]*(xsecWgt->GetBinContent(1086)/xsecWgt->GetBinContent(1)) << std::endl; + std::cout << "~~~~~~~~~~~~~~~~~~~~~~~~~~~" << std::endl; + std::cout << "isrDefHi: " << ev.ttbar_w[1086]*(xsecWgt->GetBinContent(1087)/xsecWgt->GetBinContent(1)) << std::endl; + std::cout << "~~~~~~~~~~~~~~~~~~~~~~~~~~~" << std::endl; + std::cout << "fsrDefHi: " << ev.ttbar_w[1087]*(xsecWgt->GetBinContent(1088)/xsecWgt->GetBinContent(1)) << std::endl; + std::cout << "~~~~~~~~~~~~~~~~~~~~~~~~~~~" << std::endl; + std::cout << "isrDefLo: " << ev.ttbar_w[1088]*(xsecWgt->GetBinContent(1089)/xsecWgt->GetBinContent(1)) << std::endl; + std::cout << "~~~~~~~~~~~~~~~~~~~~~~~~~~~" << std::endl; + std::cout << "fsrDefLo: " << ev.ttbar_w[1089]*(xsecWgt->GetBinContent(1090)/xsecWgt->GetBinContent(1)) << std::endl; + std::cout << "~~~~~~~~~~~~~~~~~~~~~~~~~~~" << std::endl; + std::cout << "isrConHi: " << ev.ttbar_w[1090]*(xsecWgt->GetBinContent(1091)/xsecWgt->GetBinContent(1)) << std::endl; + std::cout << "~~~~~~~~~~~~~~~~~~~~~~~~~~~" << std::endl; + std::cout << "fsrConHi: " << ev.ttbar_w[1091]*(xsecWgt->GetBinContent(1092)/xsecWgt->GetBinContent(1)) << std::endl; + std::cout << "~~~~~~~~~~~~~~~~~~~~~~~~~~~" << std::endl; + std::cout << "isrConLo: " << ev.ttbar_w[1092]*(xsecWgt->GetBinContent(1093)/xsecWgt->GetBinContent(1)) << std::endl; + std::cout << "~~~~~~~~~~~~~~~~~~~~~~~~~~~" << std::endl; + std::cout << "fsrConLo: " << ev.ttbar_w[1093]*(xsecWgt->GetBinContent(1094)/xsecWgt->GetBinContent(1)) << std::endl; + std::cout << "~~~~~~~~~~~~~~~~~~~~~~~~~~~" << std::endl; + std::cout << "~~~~~~~~~~~~~~~~~~~~~~~~~~~" << std::endl; + std::cout << "Anything else????: " << ev.ttbar_w[1094]*(xsecWgt->GetBinContent(1095)/xsecWgt->GetBinContent(1)) << std::endl; + std::cout << "~~~~~~~~~~~~~~~~~~~~~~~~~~~" << std::endl;*/ + + systWeight["isrRedHi"] = ev.ttbar_w[1082]*(xsecWgt->GetBinContent(1083)/xsecWgt->GetBinContent(1)); + systWeight["fsrRedHi"] = ev.ttbar_w[1083]*(xsecWgt->GetBinContent(1084)/xsecWgt->GetBinContent(1)); + systWeight["isrRedLo"] = ev.ttbar_w[1084]*(xsecWgt->GetBinContent(1085)/xsecWgt->GetBinContent(1)); + systWeight["fsrRedLo"] = ev.ttbar_w[1085]*(xsecWgt->GetBinContent(1086)/xsecWgt->GetBinContent(1)); + systWeight["isrDefHi"] = ev.ttbar_w[1086]*(xsecWgt->GetBinContent(1087)/xsecWgt->GetBinContent(1)); + systWeight["fsrDefHi"] = ev.ttbar_w[1087]*(xsecWgt->GetBinContent(1088)/xsecWgt->GetBinContent(1)); + systWeight["isrDefLo"] = ev.ttbar_w[1088]*(xsecWgt->GetBinContent(1089)/xsecWgt->GetBinContent(1)); + systWeight["fsrDefLo"] = ev.ttbar_w[1089]*(xsecWgt->GetBinContent(1090)/xsecWgt->GetBinContent(1)); + systWeight["isrConHi"] = ev.ttbar_w[1090]*(xsecWgt->GetBinContent(1091)/xsecWgt->GetBinContent(1)); + systWeight["fsrConHi"] = ev.ttbar_w[1091]*(xsecWgt->GetBinContent(1092)/xsecWgt->GetBinContent(1)); + systWeight["isrConLo"] = ev.ttbar_w[1092]*(xsecWgt->GetBinContent(1093)/xsecWgt->GetBinContent(1)); + systWeight["fsrConLo"] = ev.ttbar_w[1093]*(xsecWgt->GetBinContent(1094)/xsecWgt->GetBinContent(1)); + } + //std::cout<<"go go 002"<Eval(ev.nPUtrue); + if(puWgtDownGr_) puWgtLo = puWgtDownGr_->Eval(ev.nPUtrue); + if(puWgtUpGr_) puWgtHi = puWgtUpGr_->Eval(ev.nPUtrue); + if(puWgtNom<0) puWgtNom = 0; + if(puWgtLo <0) puWgtLo = 0; + if(puWgtHi <0) puWgtHi = 0; + + } + //std::cout<<"go go 003"<Fill(0.,1.0); + histos_["puwgtnorm" ]->Fill(1.,puWgtNom); + histos_["puwgtnorm" ]->Fill(2.,puWgtLo); + histos_["puwgtnorm" ]->Fill(3.,puWgtHi); + + // + //CHANNEL ASSIGNMENT + // + if(ev.ttbar_nl<2 || ev.nJet<2) continue; + ev.ttbar_chan=ev.ttbar_lid[0]*ev.ttbar_lch[0]*ev.ttbar_lid[1]*ev.ttbar_lch[1]; + + std::string ch(""); + if(ev.ttbar_chan==-11*13) ch="emu"; + if(ev.ttbar_chan==-11*11 || ev.ttbar_chan==-13*13) ch="ll"; + if(ch=="") continue; + //std::cout<<"go go 004"<>triggerBits_[ibit].first) & 1); + } + if(!hasTrigger) continue; + //std::cout<<"go go 005"< eff=getTriggerEfficiency(ev.ttbar_lid[0],ev.ttbar_lpt[0],ev.ttbar_leta[0], + ev.ttbar_lid[1],ev.ttbar_lpt[1],ev.ttbar_leta[1], + ev.ttbar_chan); + trigWgtLo=eff.first-eff.second; + trigWgtNom=eff.first; + trigWgtHi=eff.first+eff.second; + } + + //lepton selection efficiency + Float_t lepSelEffLo(1.0), lepSelEffNom(1.0), lepSelEffHi(1.0); + if(!isData) { + for(size_t il=0; il<2; il++) { + std::pair lepSF = getLeptonSelectionEfficiencyScaleFactor(ev.ttbar_lid[il],ev.ttbar_lpt[il],ev.ttbar_leta[il]); + lepSelEffLo *= (lepSF.first-lepSF.second); + lepSelEffNom *= lepSF.first; + lepSelEffHi *= (lepSF.first+lepSF.second); + } + } + //std::cout<<"go go 006"< lp4; + for(Int_t il=0; ilGetBinContent(1); + //std::cout<<"go go 007.4"<GetBinContent(1) "<GetBinContent(1)<Fill(ev.nPV-1,evWgt); + //std::cout<<"go go 007.11"<Fill(ev.ttbar_rho,evWgt); + //std::cout<<"go go 007.12"< selJets; + selJets.clear(); + std::vector > selJetsKINDisc; + std::vector< std::vector > selJetsP4; + std::vector< std::vector< std::vector > > selJetsLJKinematics; + for(Int_t ij=0; ij jesSF(3,1.0); + jecUnc_->setJetEta(fabs(jp4.Eta())); + jecUnc_->setJetPt(jp4.Pt()); + float unc = jecUnc_->getUncertainty(true); + jesSF[1]=(1.+fabs(unc)); + jesSF[2]=(1.-fabs(unc)); + + std::vector jerSF= getJetResolutionScales(jesSF[0]*jp4.Pt(), jp4.Eta(), genjpt); + + TLorentzVector oldjp4(jp4); + jp4 = jp4*jesSF[0]*jerSF[0]; + + // apply energy shifts according to systematic variation + Bool_t canBeSelected(false); + std::vector varjp4; + std::vector varkindisc; + std::vector< std::vector > varLJKinematics; + for(Int_t iSystVar=0; iSystVar<5; iSystVar++) { + varjp4.push_back( jp4 ); + if(iSystVar==1) varjp4[iSystVar] *= jesSF[1]/jesSF[0]; + if(iSystVar==2) varjp4[iSystVar] *= jesSF[2]/jesSF[0]; + if(iSystVar==3) varjp4[iSystVar] *= jerSF[1]/jerSF[0]; + if(iSystVar==4) varjp4[iSystVar] *= jerSF[2]/jerSF[0]; + + //prepare variables for MVA + std::vector< LJKinematics_t > ljkinematics; + for(Int_t il=0; il<2; il++) { + LJKinematics_t iljkin; + iljkin.dr = lp4[il].DeltaR(varjp4[iSystVar]); + iljkin.dphi = fabs(lp4[il].DeltaPhi(varjp4[iSystVar])); + iljkin.deta = fabs(lp4[il].Eta()-varjp4[iSystVar].Eta()); + iljkin.ptrel = ROOT::Math::VectorUtil::Perp(lp4[il].Vect(),varjp4[iSystVar].Vect().Unit())/lp4[il].P(); + TLorentzVector ljP4(lp4[il]+varjp4[iSystVar]); + iljkin.mlj = ljP4.M(); + iljkin.lj2ll_deta = fabs(ljP4.Eta()-dilepton.Eta()); + iljkin.lj2ll_dphi = fabs(ljP4.DeltaPhi(dilepton)); + ljkinematics.push_back(iljkin); + } + sort(ljkinematics.begin(),ljkinematics.end(),sortLJKinematicsByDR); + varLJKinematics.push_back(ljkinematics); + + //evaluate the MVA + Float_t kindisc(0.0); + if(tmvaReader_) { + for(size_t ivar=0; ivarEvaluateMVA("BDT"+methodPFix) ); + } + + // JET KINEMATIC SELECTION + //check if can be selected for this variation + if(varjp4[iSystVar].Pt()<30 || TMath::Abs(varjp4[iSystVar].Eta())>2.5) continue; + canBeSelected=true; + jetCount[iSystVar]++; + } + + //add jet if it is selectable + if(!canBeSelected) continue; + selJets.push_back(ij); + selJetsP4.push_back(varjp4); + selJetsLJKinematics.push_back( varLJKinematics ); + if(tmvaReader_) selJetsKINDisc.push_back(varkindisc); + } + //std::cout<<"go go 009"<>triggerBits_[ibit].first) & 1); - } - if(!hasTrigger) continue; - - //trigger efficiency weight - Float_t trigWgtLo(1.0), trigWgtNom(1.0), trigWgtHi(1.0); - if(!isData) - { - std::pair eff=getTriggerEfficiency(ev.ttbar_lid[0],ev.ttbar_lpt[0],ev.ttbar_leta[0], - ev.ttbar_lid[1],ev.ttbar_lpt[1],ev.ttbar_leta[1], - ev.ttbar_chan); - trigWgtLo=eff.first-eff.second; - trigWgtNom=eff.first; - trigWgtHi=eff.first+eff.second; - } - - //lepton selection efficiency - Float_t lepSelEffLo(1.0), lepSelEffNom(1.0), lepSelEffHi(1.0); - if(!isData) - { - for(size_t il=0; il<2; il++) - { - std::pair lepSF = getLeptonSelectionEfficiencyScaleFactor(ev.ttbar_lid[il],ev.ttbar_lpt[il],ev.ttbar_leta[il]); - lepSelEffLo *= (lepSF.first-lepSF.second); - lepSelEffNom *= lepSF.first; - lepSelEffHi *= (lepSF.first+lepSF.second); - } - } - - //dilepton invariant mass - std::vector lp4; - for(Int_t il=0; ilGetBinContent(1); - } - histos_[ch+"_npvinc"]->Fill(ev.nPV-1,evWgt); - npv_=ev.nPV; - - // - //JET/MET SELECTION - // - Int_t jetCount[5]={0,0,0,0,0}; - std::vector selJets; - std::vector > selJetsKINDisc; - std::vector< std::vector > selJetsP4; - std::vector< std::vector< std::vector > > selJetsLJKinematics; - for(Int_t ij=0; ij jesSF(3,1.0); - jecUnc_->setJetEta(fabs(jp4.Eta())); - jecUnc_->setJetPt(jp4.Pt()); - float unc = jecUnc_->getUncertainty(true); - jesSF[1]=(1.+fabs(unc)); - jesSF[2]=(1.-fabs(unc)); - - std::vector jerSF= getJetResolutionScales(jesSF[0]*jp4.Pt(), jp4.Eta(), genjpt); - - TLorentzVector oldjp4(jp4); - jp4 = jp4*jesSF[0]*jerSF[0]; - - // apply energy shifts according to systematic variation - Bool_t canBeSelected(false); - std::vector varjp4; - std::vector varkindisc; - std::vector< std::vector > varLJKinematics; - for(Int_t iSystVar=0; iSystVar<5; iSystVar++) - { - varjp4.push_back( jp4 ); - if(iSystVar==1) varjp4[iSystVar] *= jesSF[1]/jesSF[0]; - if(iSystVar==2) varjp4[iSystVar] *= jesSF[2]/jesSF[0]; - if(iSystVar==3) varjp4[iSystVar] *= jerSF[1]/jerSF[0]; - if(iSystVar==4) varjp4[iSystVar] *= jerSF[2]/jerSF[0]; - - //prepare variables for MVA - std::vector< LJKinematics_t > ljkinematics; - for(Int_t il=0; il<2; il++) - { - LJKinematics_t iljkin; - iljkin.dr = lp4[il].DeltaR(varjp4[iSystVar]); - iljkin.dphi = fabs(lp4[il].DeltaPhi(varjp4[iSystVar])); - iljkin.deta = fabs(lp4[il].Eta()-varjp4[iSystVar].Eta()); - iljkin.ptrel = ROOT::Math::VectorUtil::Perp(lp4[il].Vect(),varjp4[iSystVar].Vect().Unit())/lp4[il].P(); - TLorentzVector ljP4(lp4[il]+varjp4[iSystVar]); - iljkin.mlj = ljP4.M(); - iljkin.lj2ll_deta = fabs(ljP4.Eta()-dilepton.Eta()); - iljkin.lj2ll_dphi = fabs(ljP4.DeltaPhi(dilepton)); - ljkinematics.push_back(iljkin); - } - sort(ljkinematics.begin(),ljkinematics.end(),sortLJKinematicsByDR); - varLJKinematics.push_back(ljkinematics); - - //evaluate the MVA - Float_t kindisc(0.0); - if(tmvaReader_) - { - for(size_t ivar=0; ivarEvaluateMVA("BDT"+methodPFix) ); - } - - //check if can be selected for this variation - if(varjp4[iSystVar].Pt()<30 || TMath::Abs(varjp4[iSystVar].Eta())>2.5) continue; - canBeSelected=true; - jetCount[iSystVar]++; - } - - //add jet if it is selectable - if(!canBeSelected) continue; - selJets.push_back(ij); - selJetsP4.push_back(varjp4); - selJetsLJKinematics.push_back( varLJKinematics ); - if(tmvaReader_) selJetsKINDisc.push_back(varkindisc); - } - - // - // base selection and n-1 plots - // - bool zCand( (ch.Contains("ll") && TMath::Abs(mll-91)<15) ? true : false ); - bool passMet( (ch.Contains("emu") || ev.ttbar_metpt>40) ? true : false); - bool passJets(selJets.size()>=2 ? true : false); - - if(!passJets) continue; - if(zCand) ch="z"+ch; - histos_[ch+"_mllinc"]->Fill(mll,evWgt); - histos_[ch+"_met"]->Fill(ev.ttbar_metpt,evWgt); - - if(!passMet) continue; - histos_[ch+"_evsel"]->Fill(0.,evWgt); - if(selJets.size()<5) histos_[ch+"_evsel"]->Fill(selJets.size()-1,evWgt); - histos_[ch+"_rho"]->Fill(ev.ttbar_rho,evWgt); - histos_[ch+"_npv"]->Fill(ev.nPV-1,evWgt); - histos_[ch+"_mll"]->Fill(mll,evWgt); - histos_[ch+"_njets"]->Fill(selJets.size(),evWgt); - histos_[ch+"_leadjpt"]->Fill(selJetsP4[0][0].Pt(),evWgt); - histos_[ch+"_leadjeta"]->Fill((selJetsP4[0][0].Eta()),evWgt); - histos_[ch+"_leadlpt"]->Fill(lp4[0].Pt(),evWgt); - histos_[ch+"_trailjpt"]->Fill(selJetsP4[1][0].Pt(),evWgt); - histos_[ch+"_trailjeta"]->Fill(fabs(selJetsP4[1][0].Eta()),evWgt); - histos_[ch+"_traillpt"]->Fill(lp4[1].Pt(),evWgt); - - std::vector leadingkindisc(2,-9999); - std::vector leadingkindiscIdx(2,-1); - for(size_t ij=0; ijFill(selJetsLJKinematics[ij][0][0].mlj,evWgt); - histos_[ch+"_close_deta"]->Fill(selJetsLJKinematics[ij][0][0].deta,evWgt); - histos_[ch+"_close_dphi"]->Fill(selJetsLJKinematics[ij][0][0].dphi,evWgt); - histos_[ch+"_close_ptrel"]->Fill(selJetsLJKinematics[ij][0][0].ptrel,evWgt); - histos_[ch+"_close_lj2ll_deta"]->Fill(selJetsLJKinematics[ij][0][0].lj2ll_deta,evWgt); - histos_[ch+"_close_lj2ll_dphi"]->Fill(selJetsLJKinematics[ij][0][0].lj2ll_dphi,evWgt); - histos_[ch+"_far_mlj"]->Fill(selJetsLJKinematics[ij][0][1].mlj,evWgt); - histos_[ch+"_far_deta"]->Fill(selJetsLJKinematics[ij][0][1].deta,evWgt); - histos_[ch+"_far_dphi"]->Fill(selJetsLJKinematics[ij][0][1].dphi,evWgt); - histos_[ch+"_far_ptrel"]->Fill(selJetsLJKinematics[ij][0][1].ptrel,evWgt); - histos_[ch+"_far_lj2ll_deta"]->Fill(selJetsLJKinematics[ij][0][1].lj2ll_deta,evWgt); - histos_[ch+"_far_lj2ll_dphi"]->Fill(selJetsLJKinematics[ij][0][1].lj2ll_dphi,evWgt); - histos_[ch+"_j2ll_deta"]->Fill(fabs(selJetsP4[ij][0].Eta()-dilepton.Eta()),evWgt); - histos_[ch+"_j2ll_dphi"]->Fill(fabs(selJetsP4[ij][0].DeltaPhi(dilepton)),evWgt); - if(tmvaReader_) histos_[ch+"_kindisc"]->Fill(selJetsKINDisc[ij][0],evWgt); - histos_[ch+"_jp"]->Fill(ev.Jet_Proba[jetIdx],evWgt); - histos_[ch+"_svhe"]->Fill(ev.Jet_Svx[jetIdx],evWgt); - histos_[ch+"_csv"]->Fill(ev.Jet_CombIVF[jetIdx],evWgt); - histos_[ch+"_tche"]->Fill(ev.Jet_Ip2P[jetIdx],evWgt); - histos_[ch+"_jetseltrk"]->Fill(ev.Jet_nseltracks[jetIdx],evWgt); - - - Int_t flavBin(0),partonFlav(abs(ev.Jet_flavour[jetIdx])); - if(partonFlav==21 || (partonFlav>0 && partonFlav<4)) flavBin=1; - if(partonFlav==4) flavBin=2; - if(partonFlav==5) flavBin=3; - histos_[ch+"_flavour"]->Fill(flavBin,evWgt); - - //rank jets by kinematics discriminator - if(tmvaReader_) - { - if(selJetsKINDisc[ij][0]>leadingkindisc[0]) - { - leadingkindisc[1]=leadingkindisc[0]; leadingkindiscIdx[1]=leadingkindiscIdx[0]; - leadingkindisc[0]=selJetsKINDisc[ij][0]; leadingkindiscIdx[0]=jetIdx; - } - else if(selJetsKINDisc[ij][0]>leadingkindisc[1]) - { - leadingkindisc[1]=selJetsKINDisc[ij][0]; leadingkindiscIdx[1]=jetIdx; - } - } - } + // + // base selection and n-1 plots + // + bool zCand( (ch.find("ll")!=-1 && TMath::Abs(mll-91)<15) ? true : false ); + bool passMet( (ch.find("emu")!=-1 || ev.ttbar_metpt>40) ? true : false); + bool passJets(selJets.size()>=2 ? true : false); + + if(!passJets) continue; + if(zCand) ch="z"+ch; + histos_[ch+"_mllinc"]->Fill(mll,evWgt); + histos_[ch+"_met"]->Fill(ev.ttbar_metpt,evWgt); + + std::pair bestDeepFlavourPair(-1,-1); + std::pair bestDeepCSVPair(-1,-1); + std::pair bestCSVv2Pair(-1,-1); + std::pair bestCSVv2PairCheck(-1,-1); + GetBestJetPair(bestDeepFlavourPair,"deepFlavour"); + GetBestJetPair(bestDeepCSVPair,"deepCSV"); + GetBestJetPair(bestCSVv2PairCheck,"CSVv2"); + + float bestBtag=-100; + float secondBestBtag=-100; + for(int iJet=0; iJet0 && ev.Jet_DeepCSVBDisc[selJets[bestDeepCSVPair.second]]>0); + //if(!twoDeepCSVJets) { + // //std::cout<<"bestDeepCSVPair.first second deepcsvs "<Fill(0.,evWgt); + if(selJets.size()<5) histos_[ch+"_evsel"]->Fill(selJets.size()-1,evWgt); + histos_[ch+"_rho"]->Fill(ev.ttbar_rho,evWgt); + histos_[ch+"_npv"]->Fill(ev.nPV-1,evWgt); + histos_[ch+"_mll"]->Fill(mll,evWgt); + histos_[ch+"_njets"]->Fill(selJets.size(),evWgt); + //std::cout<<"go go 011"<Fill(selJetsP4[0][0].Pt(),evWgt); + //histos_[ch+"_leadjeta"]->Fill((selJetsP4[0][0].Eta()),evWgt); + //histos_[ch+"_leadlpt"]->Fill(lp4[0].Pt(),evWgt); + //histos_[ch+"_trailjpt"]->Fill(selJetsP4[1][0].Pt(),evWgt); + //histos_[ch+"_trailjeta"]->Fill(fabs(selJetsP4[1][0].Eta()),evWgt); + //histos_[ch+"_traillpt"]->Fill(lp4[1].Pt(),evWgt); + + histos_[ch+"_leadjpt"]->Fill(selJetsP4[0][0].Pt(),evWgt); + histos_[ch+"_leadjeta"]->Fill((selJetsP4[0][1].Eta()),evWgt); + histos_[ch+"_leadbjpt"]->Fill(selJetsP4[0][selJets[bestDeepCSVPair.first]].Pt(),evWgt); + histos_[ch+"_leadbjeta"]->Fill((selJetsP4[0][selJets[bestDeepCSVPair.first]].Eta()),evWgt); + histos_[ch+"_trailjpt"]->Fill(selJetsP4[1][1].Pt(),evWgt); + histos_[ch+"_trailjeta"]->Fill(fabs(selJetsP4[1][1].Eta()),evWgt); + histos_[ch+"_trailbjpt"]->Fill(selJetsP4[1][selJets[bestDeepCSVPair.second]].Pt(),evWgt); + histos_[ch+"_trailbjeta"]->Fill(fabs(selJetsP4[1][selJets[bestDeepCSVPair.second]].Eta()),evWgt); + if(selJets.size()==2){ + histos_[ch+"_only2_leadjpt"]->Fill(selJetsP4[0][0].Pt(),evWgt); + histos_[ch+"_only2_leadjeta"]->Fill((selJetsP4[0][1].Eta()),evWgt); + histos_[ch+"_only2_leadbjpt"]->Fill(selJetsP4[0][selJets[bestDeepCSVPair.first]].Pt(),evWgt); + histos_[ch+"_only2_leadbjeta"]->Fill((selJetsP4[0][selJets[bestDeepCSVPair.first]].Eta()),evWgt); + histos_[ch+"_only2_trailjpt"]->Fill(selJetsP4[1][1].Pt(),evWgt); + histos_[ch+"_only2_trailjeta"]->Fill(fabs(selJetsP4[1][1].Eta()),evWgt); + histos_[ch+"_only2_trailbjpt"]->Fill(selJetsP4[1][selJets[bestDeepCSVPair.second]].Pt(),evWgt); + histos_[ch+"_only2_trailbjeta"]->Fill(fabs(selJetsP4[1][selJets[bestDeepCSVPair.second]].Eta()),evWgt); + } + + histos_[ch+"_leadlpt"]->Fill(lp4[0].Pt(),evWgt); + histos_[ch+"_traillpt"]->Fill(lp4[1].Pt(),evWgt); + histos_[ch+"_leadleta"]->Fill(lp4[0].Eta(),evWgt); + histos_[ch+"_trailleta"]->Fill(lp4[1].Eta(),evWgt); + + // 2018 WPs + // https://twiki.cern.ch/twiki/bin/viewauth/CMS/BtagRecommendation102X#AK4_jets + + std::vector deepCSVWPs; + deepCSVWPs.clear(); + //deepCSVWPs.push_back(0.1522); //2017 + //deepCSVWPs.push_back(0.4941); + //deepCSVWPs.push_back(0.8001); + deepCSVWPs.push_back(0.1241); //2018 + deepCSVWPs.push_back(0.4184); + deepCSVWPs.push_back(0.7527); + btaggingWPs["deepCSV"]=deepCSVWPs; + + std::vector CSVv2WPs; + CSVv2WPs.clear(); + CSVv2WPs.push_back(0.5803); //2017 + CSVv2WPs.push_back(0.8838); + CSVv2WPs.push_back(0.9693); + btaggingWPs["CSVv2"]=CSVv2WPs; + + std::vector deepFlavourWPs; + deepFlavourWPs.clear(); + deepFlavourWPs.push_back(0.0494); //2018 + deepFlavourWPs.push_back(0.2770); + deepFlavourWPs.push_back(0.7264); + btaggingWPs["deepFlavour"]=deepFlavourWPs; - //control b-tagging quantities for the most promising jets in the KIN discriminator - if(tmvaReader_) - { - for(size_t ij=0; ij<2; ij++) - { - size_t jetIdx=leadingkindiscIdx[ij]; - histos_[ch+"_jp_leadkin"]->Fill(ev.Jet_Proba[jetIdx],evWgt); - histos_[ch+"_svhe_leadkin"]->Fill(ev.Jet_Svx[jetIdx],evWgt); - histos_[ch+"_csv_leadkin"]->Fill(ev.Jet_CombIVF[jetIdx],evWgt); - histos_[ch+"_tche_leadkin"]->Fill(ev.Jet_Ip2P[jetIdx],evWgt); - histos_[ch+"_jetseltrk_leadkin"]->Fill(ev.Jet_nseltracks[jetIdx],evWgt); - } - } - - // - //prepare to store trees - // - eventInfo_[0]=ev.Run; - eventInfo_[1]=ev.Evt; - eventInfo_[2]=ev.LumiBlock; - - jetmult_=selJets.size(); - ttbar_chan_=ev.ttbar_chan; - if(zCand) ttbar_chan_+= 230000; - - //weights for systematic uncertainties - for(Int_t iSystVar=0; iSystVar<5; iSystVar++) - { - Float_t selWeight(jetCount[iSystVar]>=2 ? 1.0 : 0.0); - weight_[iSystVar]=evWgt*selWeight; - } - weight_[5] = puWgtNom>0 ? evWgt*puWgtLo/puWgtNom : evWgt; - weight_[6] = puWgtLo>0 ? evWgt*puWgtHi/puWgtNom : evWgt; - weight_[7] = evWgt*trigWgtLo/trigWgtNom; - weight_[8] = evWgt*trigWgtHi/trigWgtNom; - weight_[9] = evWgt*lepSelEffLo/lepSelEffNom; - weight_[10]= evWgt*lepSelEffHi/lepSelEffNom; - weight_[11]= evWgt*qcdScaleLo/genWgt; - weight_[12]= evWgt*qcdScaleHi/genWgt; - weight_[13]= evWgt*hdampLo/genWgt; - weight_[14]= evWgt*hdampHi/genWgt; - - //fill trees - for(size_t ij=0; ijFill(); - } - - //FtM tree is filled with the two leading jets in the KIN discriminator - if(tmvaReader_) - { - for(size_t ij=0; ij<2; ij++) - { - size_t jetIdx=leadingkindiscIdx[ij]; - jetFlavour_[ij] = ev.Jet_flavour[jetIdx]; - jetPt_[ij] = ev.Jet_pt[jetIdx]; - jetEta_[ij] = ev.Jet_eta[jetIdx]; - jp_[ij] = ev.Jet_Proba[jetIdx]; - svhe_[ij] = ev.Jet_Svx[jetIdx]; - csv_[ij] = ev.Jet_CombIVF[jetIdx]; - kinDisc_[ij] = leadingkindisc[ij]; - //std::cout << ij << " " << jetFlavour_[ij] << " " - //<< jetPt_[ij] << " " << csv_[ij] << " " << kinDisc_[ij] - // << std::endl; - } - ftmTree_->Fill(); - } + + //histos_[ch+"_deepcsvlead"]->Fill(ev.Jet_DeepCSVBDisc[selJets[bestDeepCSVPair.first]],evWgt); + //histos_[ch+"_deepcsvsublead"]->Fill(ev.Jet_DeepCSVBDisc[selJets[bestDeepCSVPair.second]],evWgt); + + std::vector nPassCSVv2WPs(3,0); + for(int iWP=0; iWPCSVv2WPs[iWP]); + nPassCSVv2WPs[iWP]+= (ev.Jet_CombIVF[selJets[bestCSVv2Pair.second]]>CSVv2WPs[iWP]); + } + + std::vector nPassDeepCSVWPs(3,0); + for(int iWP=0; iWPdeepCSVWPs[iWP]); + nPassDeepCSVWPs[iWP]+= (ev.Jet_DeepCSVBDisc[selJets[bestDeepCSVPair.second]]>deepCSVWPs[iWP]); + //std::cout<<"1 2 WP nPass "<Fill(ev.Jet_DeepCSVBDisc[selJets[bestDeepCSVPair.first]],evWgt); + histos_[ch+"_deepcsv2"]->Fill(ev.Jet_DeepCSVBDisc[selJets[bestDeepCSVPair.second]],evWgt); + histos2d_[ch+"_deepcsv2d"]->Fill(ev.Jet_DeepCSVBDisc[selJets[bestDeepCSVPair.first]],ev.Jet_DeepCSVBDisc[selJets[bestDeepCSVPair.second]],evWgt); + //std::cout<<"go go 012"<Fill(ev.Jet_DeepCSVBDisc[selJets[bestDeepCSVPair.first]],evWgt); + histos_[ch+"_only2_deepcsv2"]->Fill(ev.Jet_DeepCSVBDisc[selJets[bestDeepCSVPair.second]],evWgt); + histos2d_[ch+"_only2_deepcsv2d"]->Fill(ev.Jet_DeepCSVBDisc[selJets[bestDeepCSVPair.first]],ev.Jet_DeepCSVBDisc[selJets[bestDeepCSVPair.second]],evWgt); + TwoTag(ch+"_only2_twoTags_deepCSV", "deepCSV", bestDeepCSVPair); + } + //std::cout<<"go go 014"<Fill(ev.Jet_DeepFlavourBDisc[selJets[bestDeepFlavourPair.first]],evWgt); + //std::cout<<"go go 014.1"<Fill(ev.Jet_DeepFlavourBDisc[selJets[bestDeepFlavourPair.second]],evWgt); + //std::cout<<"go go 014.2"<Fill(ev.Jet_DeepFlavourBDisc[selJets[bestDeepFlavourPair.first]],evWgt); + //std::cout<<"go go 014.4"<Fill(ev.Jet_DeepFlavourBDisc[selJets[bestDeepFlavourPair.second]],evWgt); + //std::cout<<"go go 014.5"<Fill(ev.Jet_CombIVF[selJets[bestCSVv2Pair.first]],evWgt); + histos_[ch+"_csvv22"]->Fill(ev.Jet_CombIVF[selJets[bestCSVv2Pair.second]],evWgt); + histos2d_[ch+"_csvv22d"]->Fill(ev.Jet_CombIVF[selJets[bestCSVv2Pair.first]],ev.Jet_CombIVF[selJets[bestCSVv2Pair.second]],evWgt); + + //if(binProduct==3){ + // histos2d_[ch+"_csvv22d_2b"]->Fill(ev.Jet_CombIVF[selJets[bestCSVv2Pair.first]],ev.Jet_CombIVF[selJets[bestCSVv2Pair.second]],evWgt); + //} + + TwoTag(ch+"_twoTags_CSVv2", "CSVv2", bestCSVv2Pair); + //std::cout<<"go go 014.7"<Fill(ev.Jet_CombIVF[selJets[bestCSVv2Pair.first]],evWgt); + histos_[ch+"_only2_csvv22"]->Fill(ev.Jet_CombIVF[selJets[bestCSVv2Pair.second]],evWgt); + histos2d_[ch+"_only2_csvv22d"]->Fill(ev.Jet_CombIVF[selJets[bestCSVv2Pair.first]],ev.Jet_CombIVF[selJets[bestCSVv2Pair.second]],evWgt); + //if(binProduct==3){ + // histos2d_[ch+"_only2_csvv22d_2b"]->Fill(ev.Jet_CombIVF[selJets[bestCSVv2Pair.first]],ev.Jet_CombIVF[selJets[bestCSVv2Pair.second]],evWgt); + //} + TwoTag(ch+"_only2_twoTags_CSVv2", "CSVv2", bestCSVv2Pair); + } + + std::vector leadingkindisc(2,-9999); + std::vector leadingkindiscIdx(2,-1); + for(size_t ij=0; ijFill(selJetsLJKinematics[ij][0][0].mlj,evWgt); + histos_[ch+"_close_deta"]->Fill(selJetsLJKinematics[ij][0][0].deta,evWgt); + histos_[ch+"_close_dphi"]->Fill(selJetsLJKinematics[ij][0][0].dphi,evWgt); + histos_[ch+"_close_ptrel"]->Fill(selJetsLJKinematics[ij][0][0].ptrel,evWgt); + histos_[ch+"_close_lj2ll_deta"]->Fill(selJetsLJKinematics[ij][0][0].lj2ll_deta,evWgt); + histos_[ch+"_close_lj2ll_dphi"]->Fill(selJetsLJKinematics[ij][0][0].lj2ll_dphi,evWgt); + histos_[ch+"_far_mlj"]->Fill(selJetsLJKinematics[ij][0][1].mlj,evWgt); + histos_[ch+"_far_deta"]->Fill(selJetsLJKinematics[ij][0][1].deta,evWgt); + histos_[ch+"_far_dphi"]->Fill(selJetsLJKinematics[ij][0][1].dphi,evWgt); + histos_[ch+"_far_ptrel"]->Fill(selJetsLJKinematics[ij][0][1].ptrel,evWgt); + histos_[ch+"_far_lj2ll_deta"]->Fill(selJetsLJKinematics[ij][0][1].lj2ll_deta,evWgt); + histos_[ch+"_far_lj2ll_dphi"]->Fill(selJetsLJKinematics[ij][0][1].lj2ll_dphi,evWgt); + histos_[ch+"_j2ll_deta"]->Fill(fabs(selJetsP4[ij][0].Eta()-dilepton.Eta()),evWgt); + histos_[ch+"_j2ll_dphi"]->Fill(fabs(selJetsP4[ij][0].DeltaPhi(dilepton)),evWgt); + if(tmvaReader_) histos_[ch+"_kindisc"]->Fill(selJetsKINDisc[ij][0],evWgt); + histos_[ch+"_jp"]->Fill(ev.Jet_Proba[jetIdx],evWgt); + histos_[ch+"_svhe"]->Fill(ev.Jet_Svx[jetIdx],evWgt); + histos_[ch+"_csv"]->Fill(ev.Jet_CombIVF[jetIdx],evWgt); + histos_[ch+"_deepcsv"]->Fill(ev.Jet_DeepCSVBDisc[jetIdx],evWgt); + histos_[ch+"_tche"]->Fill(ev.Jet_Ip2P[jetIdx],evWgt); + histos_[ch+"_jetseltrk"]->Fill(ev.Jet_nseltracks[jetIdx],evWgt); + + + Int_t flavBin(0),partonFlav(abs(ev.Jet_flavour[jetIdx])); + if(partonFlav==21 || (partonFlav>0 && partonFlav<4)) flavBin=1; + if(partonFlav==4) flavBin=2; + if(partonFlav==5) flavBin=3; + histos_[ch+"_flavour"]->Fill(flavBin,evWgt); + + //rank jets by kinematics discriminator + if(tmvaReader_){ + if(selJetsKINDisc[ij][0]>leadingkindisc[0]) { + leadingkindisc[1]=leadingkindisc[0]; leadingkindiscIdx[1]=leadingkindiscIdx[0]; + leadingkindisc[0]=selJetsKINDisc[ij][0]; leadingkindiscIdx[0]=jetIdx; + } else if(selJetsKINDisc[ij][0]>leadingkindisc[1]){ + leadingkindisc[1]=selJetsKINDisc[ij][0]; leadingkindiscIdx[1]=jetIdx; + } + } + } + //std::cout<<"go go 015"<Fill(ev.Jet_Proba[jetIdx],evWgt); + histos_[ch+"_svhe_leadkin"]->Fill(ev.Jet_Svx[jetIdx],evWgt); + histos_[ch+"_csv_leadkin"]->Fill(ev.Jet_CombIVF[jetIdx],evWgt); + histos_[ch+"_tche_leadkin"]->Fill(ev.Jet_Ip2P[jetIdx],evWgt); + histos_[ch+"_jetseltrk_leadkin"]->Fill(ev.Jet_nseltracks[jetIdx],evWgt); + } + } + //std::cout<<"go go 016"<=2 ? 1.0 : 0.0); + weight_[iSystVar]=evWgt*selWeight; + } + weight_[5] = puWgtNom>0 ? evWgt*puWgtLo/puWgtNom : evWgt; + weight_[6] = puWgtLo>0 ? evWgt*puWgtHi/puWgtNom : evWgt; + weight_[7] = evWgt*trigWgtLo/trigWgtNom; + weight_[8] = evWgt*trigWgtHi/trigWgtNom; + weight_[9] = evWgt*lepSelEffLo/lepSelEffNom; + weight_[10]= evWgt*lepSelEffHi/lepSelEffNom; + weight_[11]= evWgt*systWeight["qcdScaleLo"]/genWgt; + weight_[12]= evWgt*systWeight["qcdScaleHi"]/genWgt; + weight_[13]= evWgt*systWeight["hdampLo"]/genWgt; + weight_[14]= evWgt*systWeight["hdampHi"]/genWgt; + weight_[15]= evWgt*systWeight["isrRedLo"]/genWgt; + weight_[16]= evWgt*systWeight["isrRedHi"]/genWgt; + weight_[17]= evWgt*systWeight["fsrRedLo"]/genWgt; + weight_[18]= evWgt*systWeight["fsrRedHi"]/genWgt; + weight_[19]= evWgt*systWeight["isrDefLo"]/genWgt; + weight_[20]= evWgt*systWeight["isrDefHi"]/genWgt; + weight_[21]= evWgt*systWeight["fsrDefLo"]/genWgt; + weight_[22]= evWgt*systWeight["fsrDefHi"]/genWgt; + weight_[23]= evWgt*systWeight["isrConLo"]/genWgt; + weight_[24]= evWgt*systWeight["isrConHi"]/genWgt; + weight_[25]= evWgt*systWeight["fsrConLo"]/genWgt; + weight_[26]= evWgt*systWeight["fsrConHi"]/genWgt; + + //fill trees + for(size_t ij=0; ijFill(); + } + + //FtM tree is filled with the two leading jets in the KIN discriminator + if(tmvaReader_){ + for(size_t ij=0; ij<2; ij++){ + size_t jetIdx=leadingkindiscIdx[ij]; + jetFlavour_[ij] = ev.Jet_flavour[jetIdx]; + jetPt_[ij] = ev.Jet_pt[jetIdx]; + jetEta_[ij] = ev.Jet_eta[jetIdx]; + jp_[ij] = ev.Jet_Proba[jetIdx]; + svhe_[ij] = ev.Jet_Svx[jetIdx]; + csv_[ij] = ev.Jet_CombIVF[jetIdx]; + kinDisc_[ij] = leadingkindisc[ij]; + //std::cout << ij << " " << jetFlavour_[ij] << " " + //<< jetPt_[ij] << " " << csv_[ij] << " " << kinDisc_[ij] + // << std::endl; + } + ftmTree_->Fill(); + } + } + + std::cout<<"n1 n2 n3 n2/n1 n3/n1 "<noEventsSelected=(nSkipped==nentries); + + //all done with this file + inF->Close(); +} + +void TTbarEventAnalysis::GetBestJetPair(std::pair& myIndices, std::string discriminator){ + float best=-100; + float secondBest=-100; + float thisDisc=-100; + for(int iJet=0; iJet jetIndices) +{ + std::vector nPassWPs(btaggingWPs[discriminator].size(),0); + float btag1 = -2; + float btag2 = -2; + + btag1 = ReturnVarAtIndex(discriminator,jetIndices.first); + btag2 = ReturnVarAtIndex(discriminator,jetIndices.second); + + for(int iWP=0; iWPbtaggingWPs[discriminator][iWP]); + nPassWPs[iWP]+= (btag2>btaggingWPs[discriminator][iWP]); + } + + // use gen flavour of jets ; 1 is not c, not b; 2 is c; 3 is b + Int_t flavour1=abs(ev.Jet_flavour[selJets[jetIndices.first]]); + Int_t flavour2=abs(ev.Jet_flavour[selJets[jetIndices.second]]); + int bin1=0; + int bin2=0; + if(flavour1==21 || (flavour1>0 && flavour1<4)) bin1=1; + if(flavour1==4) bin1=2; + if(flavour1==5) bin1=3; + if(flavour2==21 || (flavour2>0 && flavour2<4)) bin2=1; + if(flavour2==4) bin2=2; + if(flavour2==5) bin2=3; + + // effectively we only need if both are bs or not, + // but I have bb, bX, cX and others below + int binProduct=0; + if(bin2==3 && bin1==3){ + binProduct=3; + } else if(bin1==3 || bin2==3){ + binProduct=2; + } else if(bin1==2 || bin2==2){ + binProduct=1; + } else { + binProduct=0; } - //all done with this file - inF->Close(); + std::vector twoTagCrossFlavour(nPassWPs.size(),0); + for(unsigned int iWP=0; iWPFill(twoTagCrossFlavour[iWP],evWgt); + for(unsigned int iSyst=0; iSystFill(twoTagCrossFlavour[iWP],evWgt*systWeight[systName[iSyst]]); + } + } } + //Sources // CMS AN 022/2015 v15 // https://indico.cern.ch/event/434078/#preview:1614815 @@ -622,52 +1072,52 @@ std::pair TTbarEventAnalysis::getLeptonSelectionEfficiencyScaleFact if(abs(id)==11) { if (fabs(eta)<0.8) - { - if (pt<30) { res.first=0.927; res.second=0.073; } - else if (pt<40) { res.first=0.975; res.second=0.018; } - else if (pt<50) { res.first=0.962; res.second=0.036; } - else { res.first=0.955; res.second=0.022; } - } + { + if (pt<30) { res.first=0.927; res.second=0.073; } + else if (pt<40) { res.first=0.975; res.second=0.018; } + else if (pt<50) { res.first=0.962; res.second=0.036; } + else { res.first=0.955; res.second=0.022; } + } else if (fabs(eta)<1.5) - { - if (pt<30) { res.first=0.891; res.second=0.074; } - else if (pt<40) { res.first=0.965; res.second=0.020; } - else if (pt<50) { res.first=0.968; res.second=0.018; } - else { res.first=0.955; res.second=0.018; } - } + { + if (pt<30) { res.first=0.891; res.second=0.074; } + else if (pt<40) { res.first=0.965; res.second=0.020; } + else if (pt<50) { res.first=0.968; res.second=0.018; } + else { res.first=0.955; res.second=0.018; } + } else - { - if (pt<30) { res.first=0.956; res.second=0.059; } - else if (pt<40) { res.first=0.995; res.second=0.018; } - else if (pt<50) { res.first=0.993; res.second=0.019; } - else { res.first=0.985; res.second=0.023; } - } + { + if (pt<30) { res.first=0.956; res.second=0.059; } + else if (pt<40) { res.first=0.995; res.second=0.018; } + else if (pt<50) { res.first=0.993; res.second=0.019; } + else { res.first=0.985; res.second=0.023; } + } } //muons if (abs(id)==13) { if (fabs(eta)<0.9) - { - if (pt<30) { res.first=1.003; res.second=0.019; } - else if (pt<40) { res.first=1.014; res.second=0.015; } - else if (pt<50) { res.first=1.001; res.second=0.014; } - else { res.first=0.983; res.second=0.014; } - } + { + if (pt<30) { res.first=1.003; res.second=0.019; } + else if (pt<40) { res.first=1.014; res.second=0.015; } + else if (pt<50) { res.first=1.001; res.second=0.014; } + else { res.first=0.983; res.second=0.014; } + } else if(fabs(eta)<1.2) - { - if (pt<30) { res.first=0.993; res.second=0.019; } - else if (pt<40) { res.first=0.994; res.second=0.015; } - else if (pt<50) { res.first=0.980; res.second=0.014; } - else { res.first=0.987; res.second=0.015; } - } + { + if (pt<30) { res.first=0.993; res.second=0.019; } + else if (pt<40) { res.first=0.994; res.second=0.015; } + else if (pt<50) { res.first=0.980; res.second=0.014; } + else { res.first=0.987; res.second=0.015; } + } else - { - if (pt<30) { res.first=1.023; res.second=0.028; } - else if (pt<40) { res.first=0.994; res.second=0.014; } - else if (pt<50) { res.first=0.996; res.second=0.014; } - else { res.first=0.979; res.second=0.014; } - } + { + if (pt<30) { res.first=1.023; res.second=0.028; } + else if (pt<40) { res.first=0.994; res.second=0.014; } + else if (pt<50) { res.first=0.996; res.second=0.014; } + else { res.first=0.979; res.second=0.014; } + } } return res; @@ -720,6 +1170,10 @@ void TTbarEventAnalysis::finalizeOutput() { //dump results to file outF_->cd(); + if(noEventsSelected){ + outF_->Close(); + return; + } //pileup weighting screws up a bit normalization - fix it a posteriori float puwgtSF(histos_["puwgtnorm" ]->GetBinContent(1)/histos_["puwgtnorm" ]->GetBinContent(2)); @@ -727,10 +1181,14 @@ void TTbarEventAnalysis::finalizeOutput() for(std::map::iterator it = histos_.begin(); it != histos_.end(); it++) { if(it->first!="puwgtnorm") - it->second->Scale(puwgtSF); + it->second->Scale(puwgtSF); it->second->Write(); } + for(std::map::iterator it = histos2d_.begin(); it != histos2d_.end(); it++) { + it->second->Write(); + } kinTree_->Write(); if(ftmTree_) ftmTree_->Write(); outF_->Close(); } + diff --git a/test/ttbar/TTbarEventAnalysis.h b/test/ttbar/TTbarEventAnalysis.h index fc1af362aed..b49936e731e 100644 --- a/test/ttbar/TTbarEventAnalysis.h +++ b/test/ttbar/TTbarEventAnalysis.h @@ -3,6 +3,7 @@ #include "TFile.h" #include "TH1F.h" +#include "TH2F.h" #include "TString.h" #include "TTree.h" #include "TSystem.h" @@ -13,6 +14,7 @@ #include #include #include +#include #include "CondFormats/JetMETObjects/interface/JetCorrectionUncertainty.h" @@ -24,6 +26,27 @@ struct LJKinematics_t bool sortLJKinematicsByDR (LJKinematics_t i,LJKinematics_t j) { return (i.drExpandPathName(jecUncUrl); jecUnc_ = new JetCorrectionUncertainty(jecUncUrl.Data()); //pileup weights - TString puWgtUrl("${CMSSW_BASE}/src/RecoBTag/PerformanceMeasurements/test/ttbar/data/pileupWgts.root"); - gSystem->ExpandPathName(puWgtUrl); - TFile *fIn=TFile::Open(puWgtUrl); - if(fIn){ - puWgtGr_ = (TGraph *)fIn->Get("puwgts_nom"); - puWgtDownGr_ = (TGraph *)fIn->Get("puwgts_down"); - puWgtUpGr_ = (TGraph *)fIn->Get("puwgts_up"); - fIn->Close(); - } - else{ - std::cout << "Unable to find data/pileupWgts.root, no PU reweighting will be applied" << std::endl; - } + TString stdTarget("${CMSSW_BASE}/src/RecoBTag/PerformanceMeasurements/test/ttbar/data/pileupWgts.root"); + gSystem->ExpandPathName(stdTarget); + SetPUWeightTarget(stdTarget,""); } ~TTbarEventAnalysis(){} void setReadTTJetsGenWeights(bool readTTJetsGenWeights) { readTTJetsGenWeights_=readTTJetsGenWeights; } @@ -59,6 +73,28 @@ class TTbarEventAnalysis void prepareOutput(TString outFile); void processFile(TString inFile,TH1F *xsecWgt,Bool_t isData); void finalizeOutput(); + void SetPUWeightTarget(TString targetFile,TString sampleName){ + TFile *fIn=TFile::Open(targetFile); + if(fIn){ + std::string nom( "puwgts_nom"); + std::string up( "puwgts_down"); + std::string down("puwgts_up"); + if(!sampleName.IsNull()){ + nom.append("_"); + nom.append(sampleName); + up.append("_"); + up.append(sampleName); + down.append("_"); + down.append(sampleName); + } + puWgtGr_ = (TGraph *)fIn->Get(nom.c_str()); + puWgtDownGr_ = (TGraph *)fIn->Get(down.c_str()); + puWgtUpGr_ = (TGraph *)fIn->Get(up.c_str()); + }else{ + std::cout << "Unable to find data/pileupWgts.root, no PU reweighting will be applied" << std::endl; + } + fIn->Close(); + } private: JetCorrectionUncertainty *jecUnc_; @@ -67,6 +103,19 @@ class TTbarEventAnalysis std::vector getJetEnergyScales(float pt,float eta,float rawsf,float area,float rho); std::vector getJetResolutionScales(float pt, float eta, float genjpt); + MyEventInfoBranches_t ev; + std::vector selJets; + Float_t evWgt; + + std::vector wpLabel; + std::vector systName; + std::map systWeight; + + void TwoTag(std::string tagName, std::string discriminator, std::pair); + std::map> btaggingWPs; + void GetBestJetPair(std::pair& myIndices, std::string discriminator="deepCSV"); + float ReturnVarAtIndex(std::string varName, unsigned int index); + TGraph *puWgtGr_,*puWgtDownGr_,*puWgtUpGr_; bool readTTJetsGenWeights_; TString weightsDir_; @@ -74,7 +123,7 @@ class TTbarEventAnalysis TMVA::Reader *tmvaReader_; TFile *outF_; Int_t eventInfo_[3],ttbar_chan_,npv_; - Float_t weight_[15]; + Float_t weight_[27]; Int_t jetFlavour_[2],jetmult_,jetrank_; Float_t jetPt_[2],jetEta_[2]; Float_t close_mlj_[5],close_deta_,close_dphi_,close_ptrel_,close_lj2ll_deta_, close_lj2ll_dphi_; @@ -85,6 +134,8 @@ class TTbarEventAnalysis std::vector > triggerBits_; TTree *kinTree_,*ftmTree_; std::map histos_; + std::map histos2d_; + bool noEventsSelected; }; diff --git a/test/ttbar/data/samples_Run2018.json b/test/ttbar/data/samples_Run2018.json new file mode 100644 index 00000000000..21643519372 --- /dev/null +++ b/test/ttbar/data/samples_Run2018.json @@ -0,0 +1,20 @@ +{ +"MC13TeV_TTJets_AH_training" : [380.09, 0, "/TTToHadronic_TuneCP5_13TeV-powheg-pythia8/RunIIAutumn18MiniAOD-102X_upgrade2018_realistic_v15-v1/MINIAODSIM", "t#bar{t} AH", 920, false, false, 1], +"MC13TeV_TTJets_DL_training" : [87.31, 0, "/TTTo2L2Nu_TuneCP5_13TeV-powheg-pythia8/RunIIAutumn18MiniAOD-102X_upgrade2018_realistic_v15-v1/MINIAODSIM", "t#bar{t} DL", 920, false, false, 1], +"MC13TeV_TTJets_SL_training" : [364.35, 0, "/TTToSemiLeptonic_TuneCP5_13TeV-powheg-pythia8/RunIIAutumn18MiniAOD-102X_upgrade2018_realistic_v15-v1/MINIAODSIM", "t#bar{t} SL", 920, false, false, 1], +"MC13TeV_atW" : [34.97, 0, "/ST_tW_antitop_5f_inclusiveDecays_TuneCP5_13TeV-powheg-pythia8/RunIIAutumn18MiniAOD-102X_upgrade2018_realistic_v15_ext1-v1/MINIAODSIM", "tW", 7, false, false, 1], +"MC13TeV_tW" : [34.97, 0, "/ST_tW_top_5f_inclusiveDecays_TuneCP5_13TeV-powheg-pythia8/RunIIAutumn18MiniAOD-102X_upgrade2018_realistic_v15_ext1-v1/MINIAODSIM", "tW", 7, false, false, 1], +"MC13TeV_DYJetsToLL_M-50-madgraphMLM" : [6225.42, 0, "/DYJetsToLL_M-50_TuneCP5_13TeV-madgraphMLM-pythia8/RunIIAutumn18MiniAOD-102X_upgrade2018_realistic_v15-v1/MINIAODSIM", "DY", 8, false, false, 1], +"MC13TeV_DYJetsToLL-M-50-amcatnloFXFX" : [6529.0, 0, "/DYJetsToLL_M-50_TuneCP5_13TeV-amcatnloFXFX-pythia8/RunIIAutumn18MiniAOD-102X_upgrade2018_realistic_v15-v1/MINIAODSIM", "DY", 8, false, false, 1], +"MC13TeV_DYJetsToLL_M-10to50-madgraphMLM" : [15820, 0, "/DYJetsToLL_M-10to50_TuneCP5_13TeV-madgraphMLM-pythia8/RunIIAutumn18MiniAOD-102X_upgrade2018_realistic_v15-v2/MINIAODSIM", "DY", 8, false, false, 1], +"MC13TeV_WW" : [75.8, 0, "/WW_TuneCP5_13TeV-pythia8/RunIIAutumn18MiniAOD-102X_upgrade2018_realistic_v15-v2/MINIAODSIM", "WW", 4, false, false, 1], +"MC13TeV_WZ" : [27.6, 0, "/WZ_TuneCP5_13TeV-pythia8/RunIIAutumn18MiniAOD-102X_upgrade2018_realistic_v15-v3/MINIAODSIM", "WZ", 4, false, false, 1], +"MC13TeV_ZZ" : [12.14, 0, "/ZZ_TuneCP5_13TeV-pythia8/RunIIAutumn18MiniAOD-102X_upgrade2018_realistic_v15-v2/MINIAODSIM", "ZZ", 4, false, false, 1], +"MC13TeV_W2JetsToLNu" : [2793.0, 0, "/W2JetsToLNu_TuneCP5_13TeV-madgraphMLM-pythia8/RunIIAutumn18MiniAOD-102X_upgrade2018_realistic_v15-v2/MINIAODSIM", "other", 4], +"MC13TeV_W3JetsToLNu" : [992.5, 0, "/W3JetsToLNu_TuneCP5_13TeV-madgraphMLM-pythia8/RunIIAutumn18MiniAOD-102X_upgrade2018_realistic_v15-v2/MINIAODSIM", "other", 4], +"MC13TeV_W4JetsToLNu" : [544.3, 0, "/W4JetsToLNu_TuneCP5_13TeV-madgraphMLM-pythia8/RunIIAutumn18MiniAOD-102X_upgrade2018_realistic_v15-v2/MINIAODSIM", "other", 4], +"crab_Data13TeV_MuonEG_Run2018A": [1, 1, "/MuonEG/Run2018A-PromptReco/MINIAOD", "Data", 1, false, false, 0], +"crab_Data13TeV_MuonEG_Run2018B": [1, 1, "/MuonEG/Run2018B-PromptReco/MINIAOD", "Data", 1, false, false, 0], +"crab_Data13TeV_MuonEG_Run2018C": [1, 1, "/MuonEG/Run2018C-PromptReco/MINIAOD", "Data", 1, false, false, 0], +"crab_Data13TeV_MuonEG_Run2018D": [1, 1, "/MuonEG/Run2018D-PromptReco/MINIAOD", "Data", 1, false, false, 0] +} diff --git a/test/ttbar/plotter.py b/test/ttbar/plotter.py index e633b4b6b34..cf7b4226986 100644 --- a/test/ttbar/plotter.py +++ b/test/ttbar/plotter.py @@ -23,6 +23,7 @@ def __init__(self,name): self.plotformats = ['pdf','png'] self.savelog = False self.ratiorange = (0.46,1.54) + #self.ratiorange = (0.75,1.25) def add(self, h, title, color, isData,isSyst): h.SetTitle(title) @@ -123,7 +124,7 @@ def show(self, outDir,lumi,noScale=False,saveTeX=False): p1.cd() # legend - leg = ROOT.TLegend(0.5, 0.85-0.03*max(len(self.mc)-2,0), 0.98, 0.9) + leg = ROOT.TLegend(0.5, 0.80-0.03*max(len(self.mc)-2,0), 0.98, 0.9) leg.SetBorderSize(0) leg.SetFillStyle(0) leg.SetTextFont(43) @@ -185,6 +186,12 @@ def show(self, outDir,lumi,noScale=False,saveTeX=False): if self.dataH: if maxY0: + print self.dataH.Integral()/totalMC.Integral(), + print + frame.GetYaxis().SetRangeUser(0.1,maxY*1.3) frame.SetDirectory(0) frame.Reset('ICE') @@ -230,7 +237,26 @@ def show(self, outDir,lumi,noScale=False,saveTeX=False): p2.cd() ratioframe=frame.Clone('ratioframe') ratioframe.GetYaxis().SetTitle('Data/MC') - ratioframe.GetYaxis().SetRangeUser(self.ratiorange[0], self.ratiorange[1]) + maxVal=0 + minVal=10 + dataOverMC=self.dataH.Clone("dataOverMC") + print "integrals",dataOverMC.GetName(),totalMC.GetName(),dataOverMC.Integral(),totalMC.Integral() + dataOverMC.Divide(totalMC) + for xBin in xrange(1,dataOverMC.GetXaxis().GetNbins()+1): + thisValue=dataOverMC.GetBinContent(xBin) + if thisValue==0: + continue + #print xBin,thisValue + if maxValthisValue: + minVal=thisValue + + print "maxVal,minVal,",maxVal,minVal + print "range",max(self.ratiorange[0],minVal*0.95), min(self.ratiorange[1],maxVal*1.05) + + ratioframe.GetYaxis().SetRangeUser(max(self.ratiorange[0],minVal*0.95), min(self.ratiorange[1],maxVal*1.05)) + #ratioframe.GetYaxis().SetRangeUser(self.ratiorange[0], self.ratiorange[1]) self._garbageList.append(frame) ratioframe.GetYaxis().SetNdivisions(5) ratioframe.GetYaxis().SetLabelSize(0.18) @@ -371,7 +397,8 @@ def convertToPoissonErrorGr(h): #check https://twiki.cern.ch/twiki/bin/view/CMS/PoissonErrorBars alpha = 1 - 0.6827; grpois = ROOT.TGraphAsymmErrors(h); - for i in xrange(0,grpois.GetN()+1) : + #for i in xrange(0,grpois.GetN()+1) : + for i in xrange(0,grpois.GetN()) : N = grpois.GetY()[i] if N<200 : L = 0 @@ -402,6 +429,7 @@ def main(): parser.add_option( '--rebin', dest='rebin', help='rebin factor', default=1, type=int) parser.add_option('-l', '--lumi', dest='lumi' , help='lumi to print out', default=41.6, type=float) parser.add_option( '--only', dest='only', help='plot only these (csv)', default='', type='string') + parser.add_option( '--outLabel', dest='outLabel', help='appends the plots dir name', default='', type='string') (opt, args) = parser.parse_args() #read list of samples @@ -422,11 +450,14 @@ def main(): if slist is None : continue for tag,sample in slist: - if isSyst and not 't#bar{t}' in sample[3] : continue + if isSyst and not 't#bar{t}' in sample[3] : + print "Skipping syst sample",sample + continue inDir=opt.inDir - if isSyst : inDir += '/syst' + #if isSyst : inDir += '/syst' fIn=ROOT.TFile.Open('%s/%s.root' % ( inDir, tag) ) + print "opening ",'%s/%s.root' % ( inDir, tag) try: for tkey in fIn.GetListOfKeys(): key=tkey.GetName() @@ -449,6 +480,8 @@ def main(): ROOT.gStyle.SetOptStat(0) ROOT.gROOT.SetBatch(True) outDir=opt.inDir+'/plots' + if opt.outLabel: + outDir=outDir+"_"+opt.outLabel os.system('mkdir -p %s' % outDir) for p in plots : if opt.saveLog : plots[p].savelog=True diff --git a/test/ttbar/runPileupEstimation.py b/test/ttbar/runPileupEstimation.py index 5943ba8bf0d..13b11d5a588 100644 --- a/test/ttbar/runPileupEstimation.py +++ b/test/ttbar/runPileupEstimation.py @@ -3,7 +3,7 @@ import json import commands import ROOT -from SimGeneral.MixingModule.mix_2016_25ns_SpringMC_PUScenarioV1_PoissonOOTPU_cfi import * +from SimGeneral.MixingModule.mix_2018_25ns_JuneProjectionFull18_PoissonOOTPU_cfi import * """ steer the script @@ -30,11 +30,14 @@ def main(): #compute pileup in data assuming different xsec puDist=[] puWgts=[] - MINBIASXSEC={'nom':opt.mbXsec,'up':opt.mbXsec*1.1,'down':opt.mbXsec*0.9} + XSECERROR=1.046 + MINBIASXSEC={'nom':opt.mbXsec,'up':opt.mbXsec*XSECERROR,'down':opt.mbXsec/XSECERROR} for scenario in MINBIASXSEC: print scenario, 'xsec=',MINBIASXSEC[scenario] cmd='pileupCalc.py -i %s --inputLumiJSON %s --calcMode true --minBiasXsec %f --maxPileupBin %d --numPileupBins %s Pileup.root'%(opt.inJson,opt.puJson,MINBIASXSEC[scenario],NPUBINS,NPUBINS) commands.getstatusoutput(cmd) + print cmd + fIn=ROOT.TFile.Open('Pileup.root') pileupH=fIn.Get('pileup') diff --git a/test/ttbar/runTTbarAnalysis.py b/test/ttbar/runTTbarAnalysis.py index 6848a77f535..a37648f4b84 100644 --- a/test/ttbar/runTTbarAnalysis.py +++ b/test/ttbar/runTTbarAnalysis.py @@ -60,6 +60,53 @@ def runTTbarAnalysisPacked(args): return False +def RunJobs(task_list,nJobs): + if nJobs == 0: + for inFile, outFile,wgt, tmvaWgts,isData in task_list: + runTTbarAnalysis(inFile=inFile, outFile=outFile, wgt=wgt, tmvaWgts=tmvaWgts, isData=isData) + else: + from multiprocessing import Pool + pool = Pool(nJobs) + pool.map(runTTbarAnalysisPacked, task_list) + +def InspectFiles(outputDir,runTags): + zero=0 + tiny=0 + norm=0 + + okFiles=[] + bustedFiles=[] + for tag in runTags: + os.system('ls -l %s/%s_*.root > %s.txt' % (outputDir,tag,tag) ) + sizeFile=open('%s.txt' % (tag)) + linesofls=sizeFile.readlines() + sizeFile.close() + os.system('rm %s.txt' % (tag)) + + for line in linesofls: + try: + items=line.split() + sizeOfFile=int(items[4]) + fileName=str(items[8]) + if sizeOfFile==0: + zero=zero+1 + bustedFiles.append(fileName) + elif sizeOfFile<2000: + tiny=tiny+1 + bustedFiles.append(fileName) + else: + norm=norm+1 + okFiles.append(str(items[8])) + except: + print "can't read",line + + print tag,"zero, tiny, norm",zero,tiny,norm + + print "len bustedFiles",len(bustedFiles) + + return okFiles,bustedFiles + + """ steer the script """ @@ -75,6 +122,8 @@ def main(): parser.add_option( '--tmvaWgts', dest='tmvaWgts', help='tmva weights', default=None, type='string') parser.add_option( '--dyScale', dest='dyScale', help='DY scale factor', default=None, type='string') parser.add_option('-n', '--njobs', dest='njobs', help='# jobs to run in parallel', default=0, type='int') + parser.add_option('--keepGoodFilesBySize', dest='keepGoodFilesBySize', help='Review output file sizes before running and skip good ones', default=1, type='int') + parser.add_option('--inspectOutput', dest='inspectOutput', help='Review output file sizes and re-run baddies', default=1, type='int') (opt, args) = parser.parse_args() #compile c++ wrapper to run over trees @@ -113,6 +162,14 @@ def main(): xsecWgts, integLumi = produceNormalizationCache(samplesList=samplesList,inDir=opt.inDir,cache=cache, xsecWgts=xsecWgts, integLumi=integLumi) + for key in xsecWgts: + print key,xsecWgts[key] + + print + for key in integLumi: + print key,integLumi[key] + + #DY scale factor if opt.dyScale: cachefile=open(opt.dyScale,'r') @@ -124,11 +181,20 @@ def main(): xsecWgts[tag].Scale(dySF[0]) print xsecWgts[tag].GetBinContent(1) + + if opt.keepGoodFilesBySize: + listOfSampleTags=[] + for tag,sample in samplesList: + listOfSampleTags.append(tag) + goodFiles,badFiles=InspectFiles(opt.outDir,listOfSampleTags) + else: + goodFiles=[] + #create the analysis jobs runTags = [] task_list = [] for tag,sample in samplesList: - + print tag #check if in list if len(onlyList)>0: veto=True @@ -141,23 +207,42 @@ def main(): wgt = xsecWgts[tag] for nf in xrange(0,len(input_list)) : outF='%s/%s_%d.root'%(opt.outDir,tag,nf) + if outF in goodFiles: + continue task_list.append( (input_list[nf], outF, wgt, opt.tmvaWgts, sample[1]) ) task_list=list(set(task_list)) print '%s jobs to run in %d parallel threads' % (len(task_list), opt.njobs) #run the analysis jobs - if opt.njobs == 0: - for inFile, outFile,wgt, tmvaWgts,isData in task_list: - runTTbarAnalysis(inFile=inFile, outFile=outFile, wgt=wgt, tmvaWgts=tmvaWgts, isData=isData) - else: - from multiprocessing import Pool - pool = Pool(opt.njobs) - pool.map(runTTbarAnalysisPacked, task_list) + RunJobs(task_list,opt.njobs) + + + if opt.inspectOutput: + nLoop=0 + re_task_list=[] + while(nLoop==0 or len(re_task_list)>0): + print "nLoop",nLoop + if nLoop>0: + print "going at it again",len(re_task_list) + RunJobs(re_task_list,opt.njobs) + re_task_list=[] + + goodFiles,badFiles=InspectFiles(opt.outDir,runTags) + + for task in task_list: + print "task",task + outFileNameFromTask=str(task[1]) + if outFileNameFromTask in bustedFiles: + re_task_list.append(task) + nLoop=nLoop+1 + #merge the outputs for tag in runTags: os.system('hadd -f %s/%s.root %s/%s_*.root' % (opt.outDir,tag,opt.outDir,tag) ) + + for tag in runTags: os.system('rm %s/%s_*.root' % (opt.outDir,tag) ) print 'Analysis results are available in %s' % opt.outDir diff --git a/test/ttbar/storeTools.py b/test/ttbar/storeTools.py index 89f9ce52e4b..9c7f9d462e2 100644 --- a/test/ttbar/storeTools.py +++ b/test/ttbar/storeTools.py @@ -8,29 +8,59 @@ def getEOSlslist(directory, mask='', prepend='root://eoscms//eos/cms'): from subprocess import Popen, PIPE print 'looking into: '+directory+'...' - eos_cmd = '/afs/cern.ch/project/eos/installation/0.2.41/bin/eos.select' - data = Popen([eos_cmd, 'ls', '/eos/cms/'+directory],stdout=PIPE) - out,err = data.communicate() + fullPaths=['/eos/cms/'+directory] + + eos_cmd = '/afs/cern.ch/project/eos/installation/0.3.15/bin/eos.select' + foundROOTOrEnd=False + + while not foundROOTOrEnd: + allNewPaths=[] + for fullPath in fullPaths: + data = Popen([eos_cmd, 'ls', fullPath],stdout=PIPE) + out,err = data.communicate() + items=out.split('\n') + if items[0].find(".root")!=-1: + foundROOTOrEnd=True + print "found ROOT files" + elif len(items)==0: + foundROOTOrEnd=True + print "empty ls" + else: + for item in items: + if len(item)>0: + allNewPaths.append(fullPath+'/'+item) + if not foundROOTOrEnd: + fullPaths=allNewPaths + full_list = [] + for fullPath in fullPaths: + data = Popen([eos_cmd, 'ls', fullPath],stdout=PIPE) + out,err = data.communicate() + items=out.split('\n') + + ### if input file was single root file: + #if directory.endswith('.root'): + # if len(out.split('\n')[0]) > 0: + # return [prepend + directory] - ## if input file was single root file: - if directory.endswith('.root'): - if len(out.split('\n')[0]) > 0: - return [prepend + directory] + ## instead of only the file name append the string to open the file in ROOT + for item in items: + if len(item.split()) == 0: continue + full_list.append(prepend + fullPath.split("/eos/cms")[-1] + '/' + item) - ## instead of only the file name append the string to open the file in ROOT - for line in out.split('\n'): - if len(line.split()) == 0: continue - full_list.append(prepend + directory + '/' + line) + print "found n files",len(full_list) ## strip the list of files if required if mask != '': stripped_list = [x for x in full_list if mask in x] + print "found n files in stripped list",len(stripped_list) return stripped_list - ## return - return full_list + print "found n files",len(full_list) + returnList=full_list + return returnList +# return full_list """ Loops over a list of samples and produces a cache file to normalize MC diff --git a/test/ttbar/twoTag.py b/test/ttbar/twoTag.py new file mode 100644 index 00000000000..cf38fecdab1 --- /dev/null +++ b/test/ttbar/twoTag.py @@ -0,0 +1,119 @@ +import ROOT +import sys +import math + +fileName=sys.argv[1] +tFile=ROOT.TFile.Open(fileName) + + +tags=["twoTags_deepCSVL","only2_twoTags_deepCSVL","twoTags_deepCSVM","only2_twoTags_deepCSVM","twoTags_deepCSVT","only2_twoTags_deepCSVT","twoTags_deepFlavourL","only2_twoTags_deepFlavourL","twoTags_deepFlavourM","only2_twoTags_deepFlavourM","twoTags_deepFlavourT","only2_twoTags_deepFlavourT","twoTags_CSVv2L","only2_twoTags_CSVv2L","twoTags_CSVv2M","only2_twoTags_CSVv2M","twoTags_CSVv2T","only2_twoTags_CSVv2T"] +tags.sort() +samples=["t#bar{t} SL","VV","tW","t#bar{t} DL"] +samples=["t#bar{t} SL","t#bar{t} DL","t#bar{t} AH","ZZ","WZ","WW","tW","t#bar{t} DL","DY"] +systs=[""] +systs=["","_fsrConHi","_fsrConLo","_fsrDefHi","_fsrDefLo","_fsrRedHi","_fsrRedLo","_isrConHi","_isrConLo","_isrDefHi","_isrDefLo","_isrRedHi","_isrRedLo","_qcdScaleHi","_qcdScaleLo"] + +# KEY: TH1F emu_only2_twoTags_CSVv2L_t#bar{t} AH;1 t#bar{t} AH +# KEY: TH1F emu_only2_twoTags_CSVv2L_other;1 other +# KEY: TH1F emu_only2_twoTags_CSVv2L_t#bar{t} SL;1 t#bar{t} SL +# KEY: TH1F emu_only2_twoTags_CSVv2L_tW;1 tW +# KEY: TH1F emu_only2_twoTags_CSVv2L_t#bar{t} DL;1 t#bar{t} DL + + +statError={} +#samples=["emu_only2_twoTags/emu_only2_twoTags_t#bar{t} SL","emu_only2_twoTags/emu_only2_twoTags_VV","emu_only2_twoTags/emu_only2_twoTags_tW","emu_only2_twoTags/emu_only2_twoTags_t#bar{t} DL","emu_only2_twoTags/emu_only2_twoTags_DY"] +#twoTagDataHist="emu_only2_twoTags/emu_only2_twoTags" + + +#can=ROOT.TCanvas("can","",700,700) +results={} + +for syst in systs: + #print syst + for tag in tags: + print tag, + if tag.find("twoTags_") == 0: + print " ", + for sample in samples: + histName="emu_"+tag+syst+"/emu_"+tag+syst+"_"+sample + print histName + hist=tFile.Get(histName) + if sample==samples[0]: + histBkg=hist.Clone("histBkg") + else: + histBkg.Add(hist) + #print histName + #for iBin in range(15): + # print iBin,hist.GetBinContent(iBin) + + tot={} + tot["2bmc"]=0 + tot["othermc"]=0 + tot["data"]=0 + tot["2bmc_2btag"]=0 + tot["othermc_2btag"]=0 + tot["data_2btag"]=0 + for iBin in range(15): + #print iBin,histBkg.GetBinContent(iBin) + if iBin%4==1: + tot["2bmc"]=tot["2bmc"]+histBkg.GetBinContent(iBin) + if iBin>9: + tot["2bmc_2btag"]=tot["2bmc_2btag"]+histBkg.GetBinContent(iBin) + else: + tot["othermc"]=tot["othermc"]+histBkg.GetBinContent(iBin) + if iBin>9: + tot["othermc_2btag"]=tot["othermc_2btag"]+histBkg.GetBinContent(iBin) + twoTagDataHist="emu_"+tag+syst+"/emu_"+tag+syst + hist=tFile.Get(twoTagDataHist) + for iBin in range(15): + #print iBin,hist.GetBinContent(iBin),histBkg.GetBinContent(iBin) + tot["data"]=tot["data"]+hist.GetBinContent(iBin) + if iBin>9: + tot["data_2btag"]=tot["data_2btag"]+hist.GetBinContent(iBin) + + #print tot + + #print histBkg.Integral(),hist.Integral(),hist.Integral()/histBkg.Integral() + + for name in tot: + if name.find("mc")!=-1: + tot[name]=tot[name]/histBkg.Integral() + else: + tot[name]=tot[name]/hist.Integral() + + + #print tot + + eff=math.sqrt((tot["data_2btag"]-tot["othermc_2btag"])/tot["2bmc"]) + effmc=math.sqrt(tot["2bmc_2btag"]/tot["2bmc"]) + + effUp=math.sqrt((tot["data_2btag"]-tot["othermc_2btag"]*1.5)/tot["2bmc"]) + effDown=math.sqrt((tot["data_2btag"]-tot["othermc_2btag"]*0.5)/tot["2bmc"]) + + #print tot["2bmc_2btag"],tot["2bmc"],tot["data_2btag"],tot["othermc_2btag"],eff,effmc,eff/effmc,effUp/effmc,effDown/effmc + results[tag+syst]=[tot["2bmc_2btag"],tot["2bmc"],tot["data_2btag"],tot["othermc_2btag"],eff,effmc,eff/effmc,effUp/effmc,effDown/effmc] + + + if syst=="": + print hist.Integral(),histBkg.Integral(),histBkg.GetEffectiveEntries() + statError[tag+"data_2btag"] = math.sqrt( tot["data_2btag"] * (1 - tot["data_2btag"]) / hist.Integral()) + statError[tag+"othermc_2btag"] = math.sqrt( tot["othermc_2btag"] * (1 - tot["othermc_2btag"]) / histBkg.GetEffectiveEntries()) + statError[tag+"2bmc"] = math.sqrt( tot["2bmc"] * (1 - tot["2bmc"]) / histBkg.GetEffectiveEntries()) + statError[tag+"eff"] = math.sqrt(math.pow(eff*statError[tag+"2bmc"],2) + (math.pow(statError[tag+"data_2btag"],2)+math.pow(statError[tag+"othermc_2btag"],2) ) / math.pow(eff,2) )/(2*tot["2bmc"]) + +raw_input() +print "tag,", +for syst in systs: + if syst=="": + print "data, staterror, nonbbUP, nonbbDOWN,", + else: + print syst+",", +print +for tag in tags: + print tag+",", + for syst in systs: + if syst=="": + print results[tag+syst][6],",",statError[tag+"eff"],",",results[tag+syst][7],",",results[tag+syst][8],",", + else: + print results[tag+syst][6],",", + print From 0f7c17877eda477c96083855a9842586a198c7a6 Mon Sep 17 00:00:00 2001 From: capalmer Date: Fri, 5 Apr 2019 22:13:29 +0200 Subject: [PATCH 2/7] add Josh's new plots and a few other things --- test/ttbar/TTbarEventAnalysis.cc | 746 ++++++++++++++++++++----------- test/ttbar/TTbarEventAnalysis.h | 55 ++- test/ttbar/twoTag.py | 4 +- 3 files changed, 518 insertions(+), 287 deletions(-) diff --git a/test/ttbar/TTbarEventAnalysis.cc b/test/ttbar/TTbarEventAnalysis.cc index b4a6444c87a..926e3caeca8 100644 --- a/test/ttbar/TTbarEventAnalysis.cc +++ b/test/ttbar/TTbarEventAnalysis.cc @@ -1,5 +1,6 @@ #include "TTbarEventAnalysis.h" #include "TLorentzVector.h" +#include using namespace std; @@ -39,7 +40,38 @@ void TTbarEventAnalysis::prepareOutput(TString outFile) kinTree_->Branch("jp", &(jp_[0]), "jp/F"); kinTree_->Branch("svhe", &(svhe_[0]), "svhe/F"); kinTree_->Branch("csv", &(csv_[0]), "csv/F"); - kinTree_->Branch("weight", weight_, "weight[15]/F"); + kinTree_->Branch("DeepCSVb", &(DeepCSVb_[0]), "DeepCSVb/F"); + kinTree_->Branch("DeepCSVc", &(DeepCSVc_[0]), "DeepCSVc/F"); + kinTree_->Branch("DeepCSVl", &(DeepCSVl_[0]), "DeepCSVl/F"); + kinTree_->Branch("DeepCSVbb", &(DeepCSVbb_[0]), "DeepCSVbb/F"); + kinTree_->Branch("DeepCSVcc", &(DeepCSVcc_[0]), "DeepCSVcc/F"); + kinTree_->Branch("DeepCSVbN", &(DeepCSVbN_[0]), "DeepCSVbN/F"); + kinTree_->Branch("DeepCSVcN", &(DeepCSVcN_[0]), "DeepCSVcN/F"); + kinTree_->Branch("DeepCSVlN", &(DeepCSVlN_[0]), "DeepCSVlN/F"); + kinTree_->Branch("DeepCSVbbN", &(DeepCSVbbN_[0]), "DeepCSVbbN/F"); + kinTree_->Branch("DeepCSVccN", &(DeepCSVccN_[0]), "DeepCSVccN/F"); + kinTree_->Branch("DeepCSVbP", &(DeepCSVbP_[0]), "DeepCSVbP/F"); + kinTree_->Branch("DeepCSVcP", &(DeepCSVcP_[0]), "DeepCSVcP/F"); + kinTree_->Branch("DeepCSVlP", &(DeepCSVlP_[0]), "DeepCSVlP/F"); + kinTree_->Branch("DeepCSVbbP", &(DeepCSVbbP_[0]), "DeepCSVbbP/F"); + kinTree_->Branch("DeepCSVccP", &(DeepCSVccP_[0]), "DeepCSVccP/F"); + kinTree_->Branch("DeepCSVBDisc", &(DeepCSVBDisc_[0]), "DeepCSVBDisc/F"); + kinTree_->Branch("DeepCSVBDiscN", &(DeepCSVBDiscN_[0]), "DeepCSVBDiscN/F"); + kinTree_->Branch("DeepCSVBDiscP", &(DeepCSVb_[0]), "DeepCSVb/F"); + kinTree_->Branch("DeepCSVCvsLDisc", &(DeepCSVCvsLDisc_[0]), "DeepCSVCvsLDisc/F"); + kinTree_->Branch("DeepCSVCvsLDiscN", &(DeepCSVCvsLDiscN_[0]), "DeepCSVCvsLDiscN/F"); + kinTree_->Branch("DeepCSVCvsLDiscP", &(DeepCSVCvsLDiscP_[0]), "DeepCSVCvsLDiscP/F"); + kinTree_->Branch("DeepCSVCvsBDisc", &(DeepCSVCvsBDisc_[0]), "DeepCSVCvsBDisc/F"); + kinTree_->Branch("DeepCSVCvsBDiscN", &(DeepCSVCvsBDiscN_[0]), "DeepCSVCvsBDiscN/F"); + kinTree_->Branch("DeepCSVCvsBDiscP", &(DeepCSVCvsBDiscP_[0]), "DeepCSVCvsBDiscP/F"); + + kinTree_->Branch("DeepFlavourBDisc", &(DeepFlavourBDisc_[0]), "DeepFlavourBDisc/F"); + kinTree_->Branch("DeepFlavourCvsLDisc", &(DeepFlavourCvsLDisc_[0]), "DeepFlavourCvsLDisc/F"); + kinTree_->Branch("DeepFlavourCvsBDisc", &(DeepFlavourCvsBDisc_[0]), "DeepFlavourCvsBDisc/F"); + kinTree_->Branch("DeepFlavourB", &(DeepFlavourB_[0]), "DeepFlavourB/F"); + kinTree_->Branch("DeepFlavourBB", &(DeepFlavourBB_[0]), "DeepFlavourBB/F"); + kinTree_->Branch("DeepFlavourLEPB", &(DeepFlavourLEPB_[0]), "DeepFlavourLEPB/F"); + kinTree_->Branch("weight", weight_, "weight[26]/F"); ftmTree_=new TTree("ftm","flavour tag matching"); ftmTree_->SetDirectory(outF_); @@ -52,6 +84,31 @@ void TTbarEventAnalysis::prepareOutput(TString outFile) ftmTree_->Branch("jp", jp_, "jp[2]/F"); ftmTree_->Branch("svhe", svhe_, "svhe[2]/F"); ftmTree_->Branch("csv", csv_, "csv[2]/F"); + ftmTree_->Branch("DeepCSVb", DeepCSVb_, "DeepCSVb[2]/F"); + ftmTree_->Branch("DeepCSVc", DeepCSVc_, "DeepCSVc[2]/F"); + ftmTree_->Branch("DeepCSVl", DeepCSVl_, "DeepCSVl[2]/F"); + ftmTree_->Branch("DeepCSVbb", DeepCSVbb_, "DeepCSVbb[2]/F"); + ftmTree_->Branch("DeepCSVcc", DeepCSVcc_, "DeepCSVcc[2]/F"); + ftmTree_->Branch("DeepCSVbN", DeepCSVbN_, "DeepCSVbN[2]/F"); + ftmTree_->Branch("DeepCSVcN", DeepCSVcN_, "DeepCSVcN[2]/F"); + ftmTree_->Branch("DeepCSVlN", DeepCSVlN_, "DeepCSVlN[2]/F"); + ftmTree_->Branch("DeepCSVbbN", DeepCSVbbN_, "DeepCSVbbN[2]/F"); + ftmTree_->Branch("DeepCSVccN", DeepCSVccN_, "DeepCSVccN[2]/F"); + ftmTree_->Branch("DeepCSVbP", DeepCSVbP_, "DeepCSVbP[2]/F"); + ftmTree_->Branch("DeepCSVcP", DeepCSVcP_, "DeepCSVcP[2]/F"); + ftmTree_->Branch("DeepCSVlP", DeepCSVlP_, "DeepCSVlP[2]/F"); + ftmTree_->Branch("DeepCSVbbP", DeepCSVbbP_, "DeepCSVbbP[2]/F"); + ftmTree_->Branch("DeepCSVccP", DeepCSVccP_, "DeepCSVccP[2]/F"); + + ftmTree_->Branch("DeepCSVBDisc", DeepCSVBDisc_, "DeepCSVBDisc[2]/F"); + ftmTree_->Branch("DeepCSVBDiscN", DeepCSVBDiscN_, "DeepCSVBDiscN[2]/F"); + ftmTree_->Branch("DeepCSVBDiscP", DeepCSVBDiscP_, "DeepCSVBDiscP[2]/F"); + ftmTree_->Branch("DeepCSVCvsLDisc", DeepCSVCvsLDisc_, "DeepCSVCvsLDisc[2]/F"); + ftmTree_->Branch("DeepCSVCvsLDiscN", DeepCSVCvsLDiscN_, "DeepCSVCvsLDiscN[2]/F"); + ftmTree_->Branch("DeepCSVCvsLDiscP", DeepCSVCvsLDiscP_, "DeepCSVCvsLDiscP[2]/F"); + ftmTree_->Branch("DeepCSVCvsBDisc", DeepCSVCvsBDisc_, "DeepCSVCvsBDisc[2]/F"); + ftmTree_->Branch("DeepCSVCvsBDiscN", DeepCSVCvsBDiscN_, "DeepCSVCvsBDiscN[2]/F"); + ftmTree_->Branch("DeepCSVCvsBDiscP", DeepCSVCvsBDiscP_, "DeepCSVCvsBDiscP[2]/F"); ftmTree_->Branch("kindisc", kinDisc_, "kindisc[2]/F"); ftmTree_->Branch("weight", weight_, "weight[15]/F"); @@ -75,6 +132,9 @@ void TTbarEventAnalysis::prepareOutput(TString outFile) systName.push_back("fsrConLo"); + ftmTree_->Branch("kindisc", kinDisc_, "kindisc[2]/F"); + ftmTree_->Branch("weight", weight_, "weight[26]/F"); + for(unsigned int iSyst=0; iSystSetBranchAddress("Evt" , &ev.Evt ); tree->SetBranchAddress("LumiBlock" , &ev.LumiBlock ); tree->SetBranchAddress("nPV" , &ev.nPV ); - tree->SetBranchAddress("nPUtrue", &ev.nPUtrue ); tree->SetBranchAddress("nPU", &ev.nPU ); + tree->SetBranchAddress("nPUtrue", &ev.nPUtrue ); tree->SetBranchAddress("ttbar_chan" , &ev.ttbar_chan); tree->SetBranchAddress("ttbar_metfilterWord", &ev.ttbar_metfilterWord); tree->SetBranchAddress("ttbar_trigWord", &ev.ttbar_trigWord); @@ -270,7 +393,6 @@ void TTbarEventAnalysis::processFile(TString inFile,TH1F *xsecWgt,Bool_t isData) tree->SetBranchAddress("ttbar_lch" , ev.ttbar_lch); tree->SetBranchAddress("ttbar_metpt", &ev.ttbar_metpt); tree->SetBranchAddress("ttbar_metphi", &ev.ttbar_metphi); - tree->SetBranchAddress("ttbar_rho", &ev.ttbar_rho); tree->SetBranchAddress("ttbar_nw", &ev.ttbar_nw); tree->SetBranchAddress("ttbar_w", ev.ttbar_w); tree->SetBranchAddress("nJet", &ev.nJet); @@ -283,13 +405,29 @@ void TTbarEventAnalysis::processFile(TString inFile,TH1F *xsecWgt,Bool_t isData) tree->SetBranchAddress("Jet_mass", ev.Jet_mass); tree->SetBranchAddress("Jet_Svx", ev.Jet_Svx); tree->SetBranchAddress("Jet_CombIVF", ev.Jet_CombIVF); - tree->SetBranchAddress("Jet_DeepCSVBDisc",ev.Jet_DeepCSVBDisc); - tree->SetBranchAddress("Jet_DeepFlavourBDisc",ev.Jet_DeepFlavourBDisc); tree->SetBranchAddress("Jet_Proba", ev.Jet_Proba); tree->SetBranchAddress("Jet_Ip2P", ev.Jet_Ip2P); tree->SetBranchAddress("Jet_nseltracks", ev.Jet_nseltracks); tree->SetBranchAddress("Jet_flavour", ev.Jet_flavour); + tree->SetBranchAddress("Jet_DeepCSVb", ev.Jet_DeepCSVb); + tree->SetBranchAddress("Jet_DeepCSVc", ev.Jet_DeepCSVc); + tree->SetBranchAddress("Jet_DeepCSVl", ev.Jet_DeepCSVl); + tree->SetBranchAddress("Jet_DeepCSVbN", ev.Jet_DeepCSVbN); + tree->SetBranchAddress("Jet_DeepCSVcN", ev.Jet_DeepCSVcN); + tree->SetBranchAddress("Jet_DeepCSVlN", ev.Jet_DeepCSVlN); + tree->SetBranchAddress("Jet_DeepCSVBDisc", ev.Jet_DeepCSVBDisc); + tree->SetBranchAddress("Jet_DeepCSVBDiscN", ev.Jet_DeepCSVBDiscN); + tree->SetBranchAddress("Jet_DeepCSVCvsLDisc", ev.Jet_DeepCSVCvsLDisc); + tree->SetBranchAddress("Jet_DeepCSVCvsLDiscN", ev.Jet_DeepCSVCvsLDiscN); + tree->SetBranchAddress("Jet_DeepCSVCvsBDisc", ev.Jet_DeepCSVCvsBDisc); + tree->SetBranchAddress("Jet_DeepCSVCvsBDiscN", ev.Jet_DeepCSVCvsBDiscN); + + tree->SetBranchAddress("Jet_DeepFlavourBDisc", ev.Jet_DeepFlavourBDisc); + tree->SetBranchAddress("Jet_DeepFlavourCvsLDisc", ev.Jet_DeepFlavourCvsLDisc); + tree->SetBranchAddress("Jet_DeepFlavourCvsBDisc", ev.Jet_DeepFlavourCvsBDisc); + tree->SetBranchAddress("Jet_DeepFlavourB", ev.Jet_DeepFlavourB); + int nSkipped=0; int n1=0; int n2=0; @@ -297,101 +435,10 @@ void TTbarEventAnalysis::processFile(TString inFile,TH1F *xsecWgt,Bool_t isData) for(Int_t i=0; iGetEntry(i); - //if(isData){ - // if(ev.Run>304671){ //EF - // nSkipped++; - // continue; - // } - //} //progress bar if(i%100==0) std::cout << "\r[ " << int(100.*i/nentries) << "/100 ] to completion" << std::flush; - //generator level weights - Float_t genWgt=ev.ttbar_nw==0 ? 1.0 : ev.ttbar_w[0]; - //std::cout<<"go go 001"<17) { - // Weight * [sample xsec from json / sum(weights for given systematic)] / [sample xsec from json / sum(weights for nominal)] - // = Weight * sum(weights nominal) / sum(weights for given systematic) - // i.e. renormalise to same number of events as nominal. - // N.B. for GetBinContent(X) 'X' must be weight index + 1 due to way histogram is filled (0th bin always underflow). - systWeight["qcdScaleLo"]=ev.ttbar_w[9]*(xsecWgt->GetBinContent(10)/xsecWgt->GetBinContent(1)); - systWeight["qcdScaleHi"]=ev.ttbar_w[5]*(xsecWgt->GetBinContent(6)/xsecWgt->GetBinContent(1)); - systWeight["hdampLo"]=ev.ttbar_w[ev.ttbar_nw-29]*(xsecWgt->GetBinContent(ev.ttbar_nw-29+1)/xsecWgt->GetBinContent(1)); - systWeight["hdampHi"]=ev.ttbar_w[ev.ttbar_nw-21]*(xsecWgt->GetBinContent(ev.ttbar_nw-21+1)/xsecWgt->GetBinContent(1)); - - // >>> PSWeights <<< - // Vector of weight to be used instead of old ISR/FSR varied alternative samples. - // First weight (weightID= 1081) corresponds to central ME weight value. - // The remaining 12 values (weightIDs = 1082 to 1093) correspond to the PS weights in the following order (ISR up, FSR up, ISR down, FSR down) x 3 sets, i.e. - // 1082 = isrRedHi isr:muRfac=0.707, 1083 = fsrRedHi fsr:muRfac=0.707, 1084 = isrRedLo isr:muRfac=1.414, 1085 = fsrRedLo fsr:muRfac=1.414, - // 1086 = isrDefHi isr:muRfac=0.5, 1087 = fsrDefHi fsr:muRfac=0.5, 1088 = isrDefLo isr:muRfac=2.0, 1089 = fsrDefLo fsr:muRfac=2.0, - // 1090 = isrConHi isr:muRfac=0.25, 1091 = fsrConHi fsr:muRfac=0.25, 1092 = isrConLo isr:muRfac=4.0, 1093 = fsrConLo fsr:muRfac=4.0 - - - /*std::cout << "--- ttbar_nw " << ev.ttbar_nw << " ---" << std::endl; - cout << "--- number of bins in xsecWgt: " << xsecWgt.GetNbinsX() << "---" << endl; - std::cout << "Nominal sum of weights (xsecWgt) = " << xsecWgt->GetBinContent(1) << std::endl; - std::cout << "qcdScaleLo : " << ev.ttbar_w[9]*(xsecWgt->GetBinContent(10)/xsecWgt->GetBinContent(1)) << std::endl; - std::cout << "qcdScaleHi : " << ev.ttbar_w[5]*(xsecWgt->GetBinContent(6)/xsecWgt->GetBinContent(1)) << std::endl; - std::cout << "hdampLo : " << ev.ttbar_w[ev.ttbar_nw-31]*(xsecWgt->GetBinContent(ev.ttbar_nw-31+1)/xsecWgt->GetBinContent(1)) << std::endl; - std::cout << "hdampHi : " << ev.ttbar_w[ev.ttbar_nw-23]*(xsecWgt->GetBinContent(ev.ttbar_nw-23+1)/xsecWgt->GetBinContent(1)) << std::endl; - // Store generator weights - std::cout << "~~~~~~~~~~~~~~~~~~~~~~~~~~~" << std::endl; - std::cout << "isrRedHi event weight = " << ev.ttbar_w[1082] << std::endl; - std::cout << "isrRedHi sum of event weights = " << xsecWgt->GetBinContent(1083) << std::endl; - std::cout << "isrRedHi renormalisation factor w.r.t. nominal = " << xsecWgt->GetBinContent(1083)/xsecWgt->GetBinContent(1) << std::endl; - std::cout << "Stored isrRedHi weight = " << ev.ttbar_w[1082]*(xsecWgt->GetBinContent(1083)/xsecWgt->GetBinContent(1)) << std::endl; - std::cout << "~~~~~~~~~~~~~~~~~~~~~~~~~~~" << std::endl; - std::cout << "isrRedLo event weight = " << ev.ttbar_w[1084] << std::endl; - std::cout << "isrRedLo sum of event weights = " << xsecWgt->GetBinContent(1083) << std::endl; - std::cout << "isrRedLo renormalisation factor w.r.t. nominal = " << xsecWgt->GetBinContent(1085)/xsecWgt->GetBinContent(1) << std::endl; - std::cout << "Stored isrRedLo weight = " << ev.ttbar_w[1084]*xsecWgt->GetBinContent(1085)/xsecWgt->GetBinContent(1) << std::endl; - std::cout << "~~~~~~~~~~~~~~~~~~~~~~~~~~~" << std::endl; - std::cout << "fsrRedHi event weight = " << ev.ttbar_w[1083] << std::endl; - std::cout << "fsrRedHi sum of event weights = " << xsecWgt->GetBinContent(1084) << std::endl; - std::cout << "fsrRedHi renormalisation factor w.r.t. nominal = " << xsecWgt->GetBinContent(1084)/xsecWgt->GetBinContent(1) << std::endl; - std::cout << "Stored fsrRedHi weight = " << ev.ttbar_w[1083]*(xsecWgt->GetBinContent(1084)/xsecWgt->GetBinContent(1)) << std::endl; - std::cout << "~~~~~~~~~~~~~~~~~~~~~~~~~~~" << std::endl; - std::cout << "fsrRedLo event weight = " << ev.ttbar_w[1085] << std::endl; - std::cout << "fsrRedLo sum of event weights = " << xsecWgt->GetBinContent(1086) << std::endl; - std::cout << "fsrRedLo renormaliation factor w.r.t. nominal = " << xsecWgt->GetBinContent(1086)/xsecWgt->GetBinContent(1) << std::endl; - std::cout << "Stored fsrRedLo weight = " << ev.ttbar_w[1085]*(xsecWgt->GetBinContent(1086)/xsecWgt->GetBinContent(1)) << std::endl; - std::cout << "~~~~~~~~~~~~~~~~~~~~~~~~~~~" << std::endl; - std::cout << "isrDefHi: " << ev.ttbar_w[1086]*(xsecWgt->GetBinContent(1087)/xsecWgt->GetBinContent(1)) << std::endl; - std::cout << "~~~~~~~~~~~~~~~~~~~~~~~~~~~" << std::endl; - std::cout << "fsrDefHi: " << ev.ttbar_w[1087]*(xsecWgt->GetBinContent(1088)/xsecWgt->GetBinContent(1)) << std::endl; - std::cout << "~~~~~~~~~~~~~~~~~~~~~~~~~~~" << std::endl; - std::cout << "isrDefLo: " << ev.ttbar_w[1088]*(xsecWgt->GetBinContent(1089)/xsecWgt->GetBinContent(1)) << std::endl; - std::cout << "~~~~~~~~~~~~~~~~~~~~~~~~~~~" << std::endl; - std::cout << "fsrDefLo: " << ev.ttbar_w[1089]*(xsecWgt->GetBinContent(1090)/xsecWgt->GetBinContent(1)) << std::endl; - std::cout << "~~~~~~~~~~~~~~~~~~~~~~~~~~~" << std::endl; - std::cout << "isrConHi: " << ev.ttbar_w[1090]*(xsecWgt->GetBinContent(1091)/xsecWgt->GetBinContent(1)) << std::endl; - std::cout << "~~~~~~~~~~~~~~~~~~~~~~~~~~~" << std::endl; - std::cout << "fsrConHi: " << ev.ttbar_w[1091]*(xsecWgt->GetBinContent(1092)/xsecWgt->GetBinContent(1)) << std::endl; - std::cout << "~~~~~~~~~~~~~~~~~~~~~~~~~~~" << std::endl; - std::cout << "isrConLo: " << ev.ttbar_w[1092]*(xsecWgt->GetBinContent(1093)/xsecWgt->GetBinContent(1)) << std::endl; - std::cout << "~~~~~~~~~~~~~~~~~~~~~~~~~~~" << std::endl; - std::cout << "fsrConLo: " << ev.ttbar_w[1093]*(xsecWgt->GetBinContent(1094)/xsecWgt->GetBinContent(1)) << std::endl; - std::cout << "~~~~~~~~~~~~~~~~~~~~~~~~~~~" << std::endl; - std::cout << "~~~~~~~~~~~~~~~~~~~~~~~~~~~" << std::endl; - std::cout << "Anything else????: " << ev.ttbar_w[1094]*(xsecWgt->GetBinContent(1095)/xsecWgt->GetBinContent(1)) << std::endl; - std::cout << "~~~~~~~~~~~~~~~~~~~~~~~~~~~" << std::endl;*/ - - systWeight["isrRedHi"] = ev.ttbar_w[1082]*(xsecWgt->GetBinContent(1083)/xsecWgt->GetBinContent(1)); - systWeight["fsrRedHi"] = ev.ttbar_w[1083]*(xsecWgt->GetBinContent(1084)/xsecWgt->GetBinContent(1)); - systWeight["isrRedLo"] = ev.ttbar_w[1084]*(xsecWgt->GetBinContent(1085)/xsecWgt->GetBinContent(1)); - systWeight["fsrRedLo"] = ev.ttbar_w[1085]*(xsecWgt->GetBinContent(1086)/xsecWgt->GetBinContent(1)); - systWeight["isrDefHi"] = ev.ttbar_w[1086]*(xsecWgt->GetBinContent(1087)/xsecWgt->GetBinContent(1)); - systWeight["fsrDefHi"] = ev.ttbar_w[1087]*(xsecWgt->GetBinContent(1088)/xsecWgt->GetBinContent(1)); - systWeight["isrDefLo"] = ev.ttbar_w[1088]*(xsecWgt->GetBinContent(1089)/xsecWgt->GetBinContent(1)); - systWeight["fsrDefLo"] = ev.ttbar_w[1089]*(xsecWgt->GetBinContent(1090)/xsecWgt->GetBinContent(1)); - systWeight["isrConHi"] = ev.ttbar_w[1090]*(xsecWgt->GetBinContent(1091)/xsecWgt->GetBinContent(1)); - systWeight["fsrConHi"] = ev.ttbar_w[1091]*(xsecWgt->GetBinContent(1092)/xsecWgt->GetBinContent(1)); - systWeight["isrConLo"] = ev.ttbar_w[1092]*(xsecWgt->GetBinContent(1093)/xsecWgt->GetBinContent(1)); - systWeight["fsrConLo"] = ev.ttbar_w[1093]*(xsecWgt->GetBinContent(1094)/xsecWgt->GetBinContent(1)); - } - //std::cout<<"go go 002"<Fill(0.,1.0); @@ -420,7 +466,6 @@ void TTbarEventAnalysis::processFile(TString inFile,TH1F *xsecWgt,Bool_t isData) if(ev.ttbar_chan==-11*13) ch="emu"; if(ev.ttbar_chan==-11*11 || ev.ttbar_chan==-13*13) ch="ll"; if(ch=="") continue; - //std::cout<<"go go 004"<>triggerBits_[ibit].first) & 1); } if(!hasTrigger) continue; - //std::cout<<"go go 005"< lp4; @@ -465,59 +508,80 @@ void TTbarEventAnalysis::processFile(TString inFile,TH1F *xsecWgt,Bool_t isData) TLorentzVector dilepton(lp4[0]+lp4[1]); Float_t mll=dilepton.M(); - //if(mll<12) continue; if(lp4[0].Pt()<30) continue; if(lp4[1].Pt()<20) continue; if(mll<90) continue; - //if(lp4[0].Pt()<25 || lp4[1].Pt()<25) continue; - //std::cout<<"go go 007"<GetBinContent(1); - //std::cout<<"go go 007.4"<GetBinContent(1) "<GetBinContent(1)<17) { + // Weight * [sample xsec from json / sum(weights for given systematic)] / [sample xsec from json / sum(weights for nominal)] + // = Weight * sum(weights nominal) / sum(weights for given systematic) + // i.e. renormalise to same number of events as nominal. + // N.B. for GetBinContent(X) 'X' must be weight index + 1 due to way histogram is filled (0th bin always underflow). + systWeight["qcdScaleLo"]= evWgt*ev.ttbar_w[9]/genWgt *(xsecWgt->GetBinContent(10)/xsecWgt->GetBinContent(1)); + systWeight["qcdScaleHi"]= evWgt*ev.ttbar_w[5]/genWgt *(xsecWgt->GetBinContent(6)/xsecWgt->GetBinContent(1)); + systWeight["hdampLo"]= evWgt*ev.ttbar_w[ev.ttbar_nw-29]/genWgt*(xsecWgt->GetBinContent(ev.ttbar_nw-29+1)/xsecWgt->GetBinContent(1)); + systWeight["hdampHi"]= evWgt*ev.ttbar_w[ev.ttbar_nw-21]/genWgt*(xsecWgt->GetBinContent(ev.ttbar_nw-21+1)/xsecWgt->GetBinContent(1)); + // >>> PSWeights <<< + // Vector of weight to be used instead of old ISR/FSR varied alternative samples. + // First weight (weightID= 1081) corresponds to central ME weight value. + // The remaining 12 values (weightIDs = 1082 to 1093) correspond to the PS weights in the following order (ISR up, FSR up, ISR down, FSR down) x 3 sets, i.e. + // 1082 = isrRedHi isr:muRfac=0.707, 1083 = fsrRedHi fsr:muRfac=0.707, 1084 = isrRedLo isr:muRfac=1.414, 1085 = fsrRedLo fsr:muRfac=1.414, + // 1086 = isrDefHi isr:muRfac=0.5, 1087 = fsrDefHi fsr:muRfac=0.5, 1088 = isrDefLo isr:muRfac=2.0, 1089 = fsrDefLo fsr:muRfac=2.0, + // 1090 = isrConHi isr:muRfac=0.25, 1091 = fsrConHi fsr:muRfac=0.25, 1092 = isrConLo isr:muRfac=4.0, 1093 = fsrConLo fsr:muRfac=4.0 + + systWeight["isrRedHi"] = evWgt*ev.ttbar_w[1082]/genWgt*(xsecWgt->GetBinContent(1083)/xsecWgt->GetBinContent(1)); + systWeight["fsrRedHi"] = evWgt*ev.ttbar_w[1083]/genWgt*(xsecWgt->GetBinContent(1084)/xsecWgt->GetBinContent(1)); + systWeight["isrRedLo"] = evWgt*ev.ttbar_w[1084]/genWgt*(xsecWgt->GetBinContent(1085)/xsecWgt->GetBinContent(1)); + systWeight["fsrRedLo"] = evWgt*ev.ttbar_w[1085]/genWgt*(xsecWgt->GetBinContent(1086)/xsecWgt->GetBinContent(1)); + systWeight["isrDefHi"] = evWgt*ev.ttbar_w[1086]/genWgt*(xsecWgt->GetBinContent(1087)/xsecWgt->GetBinContent(1)); + systWeight["fsrDefHi"] = evWgt*ev.ttbar_w[1087]/genWgt*(xsecWgt->GetBinContent(1088)/xsecWgt->GetBinContent(1)); + systWeight["isrDefLo"] = evWgt*ev.ttbar_w[1088]/genWgt*(xsecWgt->GetBinContent(1089)/xsecWgt->GetBinContent(1)); + systWeight["fsrDefLo"] = evWgt*ev.ttbar_w[1089]/genWgt*(xsecWgt->GetBinContent(1090)/xsecWgt->GetBinContent(1)); + systWeight["isrConHi"] = evWgt*ev.ttbar_w[1090]/genWgt*(xsecWgt->GetBinContent(1091)/xsecWgt->GetBinContent(1)); + systWeight["fsrConHi"] = evWgt*ev.ttbar_w[1091]/genWgt*(xsecWgt->GetBinContent(1092)/xsecWgt->GetBinContent(1)); + systWeight["isrConLo"] = evWgt*ev.ttbar_w[1092]/genWgt*(xsecWgt->GetBinContent(1093)/xsecWgt->GetBinContent(1)); + systWeight["fsrConLo"] = evWgt*ev.ttbar_w[1093]/genWgt*(xsecWgt->GetBinContent(1094)/xsecWgt->GetBinContent(1)); + } + + //std::cout<<"nom weight "<Fill(ev.nPV-1,evWgt); - //std::cout<<"go go 007.11"<Fill(ev.ttbar_rho,evWgt); - //std::cout<<"go go 007.12"< selJets; selJets.clear(); std::vector > selJetsKINDisc; std::vector< std::vector > selJetsP4; std::vector< std::vector< std::vector > > selJetsLJKinematics; - for(Int_t ij=0; ij bestDeepFlavourPair(-1,-1); std::pair bestDeepCSVPair(-1,-1); std::pair bestCSVv2Pair(-1,-1); - std::pair bestCSVv2PairCheck(-1,-1); GetBestJetPair(bestDeepFlavourPair,"deepFlavour"); GetBestJetPair(bestDeepCSVPair,"deepCSV"); - GetBestJetPair(bestCSVv2PairCheck,"CSVv2"); - - float bestBtag=-100; - float secondBestBtag=-100; - for(int iJet=0; iJet0 && ev.Jet_DeepCSVBDisc[selJets[bestDeepCSVPair.second]]>0); @@ -667,37 +708,26 @@ void TTbarEventAnalysis::processFile(TString inFile,TH1F *xsecWgt,Bool_t isData) if(!passMet) continue; histos_[ch+"_evsel"]->Fill(0.,evWgt); if(selJets.size()<5) histos_[ch+"_evsel"]->Fill(selJets.size()-1,evWgt); - histos_[ch+"_rho"]->Fill(ev.ttbar_rho,evWgt); histos_[ch+"_npv"]->Fill(ev.nPV-1,evWgt); histos_[ch+"_mll"]->Fill(mll,evWgt); histos_[ch+"_njets"]->Fill(selJets.size(),evWgt); - //std::cout<<"go go 011"<Fill(selJetsP4[0][0].Pt(),evWgt); - //histos_[ch+"_leadjeta"]->Fill((selJetsP4[0][0].Eta()),evWgt); - //histos_[ch+"_leadlpt"]->Fill(lp4[0].Pt(),evWgt); - //histos_[ch+"_trailjpt"]->Fill(selJetsP4[1][0].Pt(),evWgt); - //histos_[ch+"_trailjeta"]->Fill(fabs(selJetsP4[1][0].Eta()),evWgt); - //histos_[ch+"_traillpt"]->Fill(lp4[1].Pt(),evWgt); - histos_[ch+"_leadjpt"]->Fill(selJetsP4[0][0].Pt(),evWgt); - histos_[ch+"_leadjeta"]->Fill((selJetsP4[0][1].Eta()),evWgt); + histos_[ch+"_leadjeta"]->Fill((selJetsP4[0][0].Eta()),evWgt); histos_[ch+"_leadbjpt"]->Fill(selJetsP4[0][selJets[bestDeepCSVPair.first]].Pt(),evWgt); histos_[ch+"_leadbjeta"]->Fill((selJetsP4[0][selJets[bestDeepCSVPair.first]].Eta()),evWgt); - histos_[ch+"_trailjpt"]->Fill(selJetsP4[1][1].Pt(),evWgt); - histos_[ch+"_trailjeta"]->Fill(fabs(selJetsP4[1][1].Eta()),evWgt); - histos_[ch+"_trailbjpt"]->Fill(selJetsP4[1][selJets[bestDeepCSVPair.second]].Pt(),evWgt); - histos_[ch+"_trailbjeta"]->Fill(fabs(selJetsP4[1][selJets[bestDeepCSVPair.second]].Eta()),evWgt); + histos_[ch+"_trailjpt"]->Fill(selJetsP4[0][1].Pt(),evWgt); + histos_[ch+"_trailjeta"]->Fill(fabs(selJetsP4[0][1].Eta()),evWgt); + histos_[ch+"_trailbjpt"]->Fill(selJetsP4[0][selJets[bestDeepCSVPair.second]].Pt(),evWgt); + histos_[ch+"_trailbjeta"]->Fill(fabs(selJetsP4[0][selJets[bestDeepCSVPair.second]].Eta()),evWgt); if(selJets.size()==2){ histos_[ch+"_only2_leadjpt"]->Fill(selJetsP4[0][0].Pt(),evWgt); - histos_[ch+"_only2_leadjeta"]->Fill((selJetsP4[0][1].Eta()),evWgt); + histos_[ch+"_only2_leadjeta"]->Fill((selJetsP4[0][0].Eta()),evWgt); histos_[ch+"_only2_leadbjpt"]->Fill(selJetsP4[0][selJets[bestDeepCSVPair.first]].Pt(),evWgt); histos_[ch+"_only2_leadbjeta"]->Fill((selJetsP4[0][selJets[bestDeepCSVPair.first]].Eta()),evWgt); - histos_[ch+"_only2_trailjpt"]->Fill(selJetsP4[1][1].Pt(),evWgt); - histos_[ch+"_only2_trailjeta"]->Fill(fabs(selJetsP4[1][1].Eta()),evWgt); - histos_[ch+"_only2_trailbjpt"]->Fill(selJetsP4[1][selJets[bestDeepCSVPair.second]].Pt(),evWgt); - histos_[ch+"_only2_trailbjeta"]->Fill(fabs(selJetsP4[1][selJets[bestDeepCSVPair.second]].Eta()),evWgt); + histos_[ch+"_only2_trailjpt"]->Fill(selJetsP4[0][1].Pt(),evWgt); + histos_[ch+"_only2_trailjeta"]->Fill(fabs(selJetsP4[0][1].Eta()),evWgt); + histos_[ch+"_only2_trailbjpt"]->Fill(selJetsP4[0][selJets[bestDeepCSVPair.second]].Pt(),evWgt); + histos_[ch+"_only2_trailbjeta"]->Fill(fabs(selJetsP4[0][selJets[bestDeepCSVPair.second]].Eta()),evWgt); } histos_[ch+"_leadlpt"]->Fill(lp4[0].Pt(),evWgt); @@ -737,13 +767,13 @@ void TTbarEventAnalysis::processFile(TString inFile,TH1F *xsecWgt,Bool_t isData) //histos_[ch+"_deepcsvsublead"]->Fill(ev.Jet_DeepCSVBDisc[selJets[bestDeepCSVPair.second]],evWgt); std::vector nPassCSVv2WPs(3,0); - for(int iWP=0; iWPCSVv2WPs[iWP]); nPassCSVv2WPs[iWP]+= (ev.Jet_CombIVF[selJets[bestCSVv2Pair.second]]>CSVv2WPs[iWP]); } std::vector nPassDeepCSVWPs(3,0); - for(int iWP=0; iWPdeepCSVWPs[iWP]); nPassDeepCSVWPs[iWP]+= (ev.Jet_DeepCSVBDisc[selJets[bestDeepCSVPair.second]]>deepCSVWPs[iWP]); //std::cout<<"1 2 WP nPass "<Fill(ev.Jet_DeepCSVBDisc[jetIdx],evWgt); histos_[ch+"_tche"]->Fill(ev.Jet_Ip2P[jetIdx],evWgt); histos_[ch+"_jetseltrk"]->Fill(ev.Jet_nseltracks[jetIdx],evWgt); + histos_[ch+"_DeepCSVb"]->Fill(ev.Jet_DeepCSVb[jetIdx],evWgt); + histos_[ch+"_DeepCSVc"]->Fill(ev.Jet_DeepCSVc[jetIdx],evWgt); + histos_[ch+"_DeepCSVl"]->Fill(ev.Jet_DeepCSVl[jetIdx],evWgt); + histos_[ch+"_DeepCSVbN"]->Fill(ev.Jet_DeepCSVbN[jetIdx],evWgt); + histos_[ch+"_DeepCSVcN"]->Fill(ev.Jet_DeepCSVcN[jetIdx],evWgt); + histos_[ch+"_DeepCSVlN"]->Fill(ev.Jet_DeepCSVlN[jetIdx],evWgt); + histos_[ch+"_DeepCSVBDisc"]->Fill(ev.Jet_DeepCSVBDisc[jetIdx],evWgt); + histos_[ch+"_DeepCSVBDiscN"]->Fill(ev.Jet_DeepCSVBDiscN[jetIdx],evWgt); + histos_[ch+"_DeepCSVCvsLDisc"]->Fill(ev.Jet_DeepCSVCvsLDisc[jetIdx],evWgt); + histos_[ch+"_DeepCSVCvsLDiscN"]->Fill(ev.Jet_DeepCSVCvsLDiscN[jetIdx],evWgt); + histos_[ch+"_DeepCSVCvsBDisc"]->Fill(ev.Jet_DeepCSVCvsBDisc[jetIdx],evWgt); + histos_[ch+"_DeepCSVCvsBDiscN"]->Fill(ev.Jet_DeepCSVCvsBDiscN[jetIdx],evWgt); + + histos_[ch+"_DeepFlavourBDisc"]->Fill(ev.Jet_DeepFlavourBDisc[jetIdx],evWgt); + histos_[ch+"_DeepFlavourCvsLDisc"]->Fill(ev.Jet_DeepFlavourCvsLDisc[jetIdx],evWgt); + histos_[ch+"_DeepFlavourCvsBDisc"]->Fill(ev.Jet_DeepFlavourCvsBDisc[jetIdx],evWgt); + histos_[ch+"_DeepFlavourB"]->Fill(ev.Jet_DeepFlavourB[jetIdx],evWgt); Int_t flavBin(0),partonFlav(abs(ev.Jet_flavour[jetIdx])); @@ -844,7 +891,6 @@ void TTbarEventAnalysis::processFile(TString inFile,TH1F *xsecWgt,Bool_t isData) } } } - //std::cout<<"go go 015"<Fill(ev.Jet_CombIVF[jetIdx],evWgt); histos_[ch+"_tche_leadkin"]->Fill(ev.Jet_Ip2P[jetIdx],evWgt); histos_[ch+"_jetseltrk_leadkin"]->Fill(ev.Jet_nseltracks[jetIdx],evWgt); + histos_[ch+"_DeepCSVb_leadkin"]->Fill(ev.Jet_DeepCSVb[jetIdx],evWgt); + histos_[ch+"_DeepCSVc_leadkin"]->Fill(ev.Jet_DeepCSVc[jetIdx],evWgt); + histos_[ch+"_DeepCSVl_leadkin"]->Fill(ev.Jet_DeepCSVl[jetIdx],evWgt); + histos_[ch+"_DeepCSVbN_leadkin"]->Fill(ev.Jet_DeepCSVbN[jetIdx],evWgt); + histos_[ch+"_DeepCSVcN_leadkin"]->Fill(ev.Jet_DeepCSVcN[jetIdx],evWgt); + histos_[ch+"_DeepCSVlN_leadkin"]->Fill(ev.Jet_DeepCSVlN[jetIdx],evWgt); + histos_[ch+"_DeepCSVBDisc_leadkin"]->Fill(ev.Jet_DeepCSVBDisc[jetIdx],evWgt); + histos_[ch+"_DeepCSVBDiscN_leadkin"]->Fill(ev.Jet_DeepCSVBDiscN[jetIdx],evWgt); + histos_[ch+"_DeepCSVCvsLDisc_leadkin"]->Fill(ev.Jet_DeepCSVCvsLDisc[jetIdx],evWgt); + histos_[ch+"_DeepCSVCvsLDiscN_leadkin"]->Fill(ev.Jet_DeepCSVCvsLDiscN[jetIdx],evWgt); + histos_[ch+"_DeepCSVCvsBDisc_leadkin"]->Fill(ev.Jet_DeepCSVCvsBDisc[jetIdx],evWgt); + histos_[ch+"_DeepCSVCvsBDiscN_leadkin"]->Fill(ev.Jet_DeepCSVCvsBDiscN[jetIdx],evWgt); + histos_[ch+"_DeepFlavourBDisc_leadkin"]->Fill(ev.Jet_DeepFlavourBDisc[jetIdx],evWgt); + histos_[ch+"_DeepFlavourCvsLDisc_leadkin"]->Fill(ev.Jet_DeepFlavourCvsLDisc[jetIdx],evWgt); + histos_[ch+"_DeepFlavourCvsBDisc_leadkin"]->Fill(ev.Jet_DeepFlavourCvsBDisc[jetIdx],evWgt); + histos_[ch+"_DeepFlavourB_leadkin"]->Fill(ev.Jet_DeepFlavourB[jetIdx],evWgt); } } - //std::cout<<"go go 016"<Fill(); } @@ -945,6 +1025,24 @@ void TTbarEventAnalysis::processFile(TString inFile,TH1F *xsecWgt,Bool_t isData) svhe_[ij] = ev.Jet_Svx[jetIdx]; csv_[ij] = ev.Jet_CombIVF[jetIdx]; kinDisc_[ij] = leadingkindisc[ij]; + DeepCSVb_[ij] = ev.Jet_DeepCSVb[jetIdx]; + DeepCSVc_[ij] = ev.Jet_DeepCSVc[jetIdx]; + DeepCSVl_[ij] = ev.Jet_DeepCSVl[jetIdx]; + DeepCSVbN_[ij] = ev.Jet_DeepCSVbN[jetIdx]; + DeepCSVcN_[ij] = ev.Jet_DeepCSVcN[jetIdx]; + DeepCSVlN_[ij] = ev.Jet_DeepCSVlN[jetIdx]; + + DeepCSVBDisc_[ij] = ev.Jet_DeepCSVBDisc[jetIdx]; + DeepCSVBDiscN_[ij] = ev.Jet_DeepCSVBDiscN[jetIdx]; + DeepCSVCvsLDisc_[ij] = ev.Jet_DeepCSVCvsLDisc[jetIdx]; + DeepCSVCvsLDiscN_[ij] = ev.Jet_DeepCSVCvsLDiscN[jetIdx]; + DeepCSVCvsBDisc_[ij] = ev.Jet_DeepCSVCvsBDisc[jetIdx]; + DeepCSVCvsBDiscN_[ij] = ev.Jet_DeepCSVCvsBDiscN[jetIdx]; + + DeepFlavourBDisc_[ij] = ev.Jet_DeepFlavourBDisc[jetIdx]; + DeepFlavourCvsLDisc_[ij] = ev.Jet_DeepFlavourCvsLDisc[jetIdx]; + DeepFlavourCvsBDisc_[ij] = ev.Jet_DeepFlavourCvsBDisc[jetIdx]; + DeepFlavourB_[ij] = ev.Jet_DeepFlavourB[jetIdx]; //std::cout << ij << " " << jetFlavour_[ij] << " " //<< jetPt_[ij] << " " << csv_[ij] << " " << kinDisc_[ij] // << std::endl; @@ -964,7 +1062,7 @@ void TTbarEventAnalysis::GetBestJetPair(std::pair& myIndices, std::str float best=-100; float secondBest=-100; float thisDisc=-100; - for(int iJet=0; iJet& myIndices, std::str } } - for(int iJet=0; iJetbtaggingWPs[discriminator][iWP]); nPassWPs[iWP]+= (btag2>btaggingWPs[discriminator][iWP]); } @@ -1044,7 +1142,7 @@ void TTbarEventAnalysis::TwoTag(std::string tagName, std::string discriminator, histos_[tagName+wpLabel[iWP]]->Fill(twoTagCrossFlavour[iWP],evWgt); for(unsigned int iSyst=0; iSystFill(twoTagCrossFlavour[iWP],evWgt*systWeight[systName[iSyst]]); + histos_[tagName+wpLabel[iWP]+"_"+systName[iSyst]]->Fill(twoTagCrossFlavour[iWP],systWeight[systName[iSyst]]); } } } @@ -1062,8 +1160,8 @@ std::pair TTbarEventAnalysis::getTriggerEfficiency(int id1,float pt return res; } -//Sources -// CMS AN 022/2015 v15 +//https://twiki.cern.ch/twiki/bin/view/CMS/Egamma2017DataRecommendations#Electron_Reconstruction_Scale_Fa +//https://soffi.web.cern.ch/soffi/EGM-ID/SF-17Nov2017-MCv2-IDv1-020618/Electrons/egammaEffi.txt_egammaPlots_runBCDEF_passingMedium94X.pdf std::pair TTbarEventAnalysis::getLeptonSelectionEfficiencyScaleFactor(int id,float pt,float eta) { std::pairres(1.0,0.0); @@ -1071,55 +1169,126 @@ std::pair TTbarEventAnalysis::getLeptonSelectionEfficiencyScaleFact //electrons if(abs(id)==11) { - if (fabs(eta)<0.8) - { - if (pt<30) { res.first=0.927; res.second=0.073; } - else if (pt<40) { res.first=0.975; res.second=0.018; } - else if (pt<50) { res.first=0.962; res.second=0.036; } - else { res.first=0.955; res.second=0.022; } - } - else if (fabs(eta)<1.5) - { - if (pt<30) { res.first=0.891; res.second=0.074; } - else if (pt<40) { res.first=0.965; res.second=0.020; } - else if (pt<50) { res.first=0.968; res.second=0.018; } - else { res.first=0.955; res.second=0.018; } - } - else - { - if (pt<30) { res.first=0.956; res.second=0.059; } - else if (pt<40) { res.first=0.995; res.second=0.018; } - else if (pt<50) { res.first=0.993; res.second=0.019; } - else { res.first=0.985; res.second=0.023; } - } + if (2.0 < eta < 2.5) + { + if (10.0 TTbarEventAnalysis::getJetResolutionScales(float pt, float et std::vector res(3,1.0); float ptSF(1.0), ptSF_err(0.06); - if(TMath::Abs(eta)<0.5) - { - ptSF=1.0; - ptSF_err = TMath::Sqrt(pow(0.012,2)+pow(0.5*(0.062+0.061),2)); - } - else if(TMath::Abs(eta)<1.1) - { - ptSF=1.0; - ptSF_err = TMath::Sqrt(pow(0.012,2)+pow(0.5*(0.056+0.055),2)); - } - else if(TMath::Abs(eta)<1.7) - { - ptSF=1.0; - ptSF_err = TMath::Sqrt(pow(0.017,2)+pow(0.5*(0.063+0.062),2)); - } - else if(TMath::Abs(eta)<2.3) - { - ptSF=1.0; - ptSF_err = TMath::Sqrt(pow(0.035,2)+pow(0.5*(0.087+0.085),2)); - } - else - { - ptSF=1.0; - ptSF_err = TMath::Sqrt(pow(0.127,2)+pow(0.5*(0.155+0.153),2)); - } + if(TMath::Abs(eta)<0.522){ + //ptSF=1.1432; + ptSF=1.0; + ptSF_err = TMath::Sqrt(pow(0.012,2)+pow(0.5*(0.062+0.061),2)); + //ptSF_err = 0.0222; + } + else if(TMath::Abs(eta)<0.783){ + //ptSF=1.1815; + ptSF=1.0; + ptSF_err = TMath::Sqrt(pow(0.012,2)+pow(0.5*(0.056+0.055),2)); + //ptSF_err = 0.0484; + } + else if(TMath::Abs(eta)<1.131){ + //ptSF=1.0989; + ptSF=1.0; + ptSF_err = TMath::Sqrt(pow(0.017,2)+pow(0.5*(0.063+0.062),2)); + //ptSF_err = 0.0456; + } + else if(TMath::Abs(eta)<1.305){ + //ptSF=1.1137; + ptSF=1.0; + ptSF_err = TMath::Sqrt(pow(0.035,2)+pow(0.5*(0.087+0.085),2)); + //ptSF_err = 0.1397; + } + else if(TMath::Abs(eta)<1.740){ + //ptSF=1.1307; + ptSF=1.0; + ptSF_err = TMath::Sqrt(pow(0.035,2)+pow(0.5*(0.087+0.085),2)); + //ptSF_err = 0.1470; + } + else if(TMath::Abs(eta)<1.930){ + //ptSF=1.1600; + ptSF=1.0; + ptSF_err = TMath::Sqrt(pow(0.035,2)+pow(0.5*(0.087+0.085),2)); + //ptSF_err = 0.0976; + } + else if(TMath::Abs(eta)<2.043){ + //ptSF=1.2393; + ptSF=1.0; + ptSF_err = TMath::Sqrt(pow(0.035,2)+pow(0.5*(0.087+0.085),2)); + //ptSF_err = 0.1909; + } + else if(TMath::Abs(eta)<2.322){ + //ptSF=1.2604; + ptSF=1.0; + ptSF_err = TMath::Sqrt(pow(0.035,2)+pow(0.5*(0.087+0.085),2)); + //ptSF_err = 0.1501; + } + else if(TMath::Abs(eta)<2.500){ + //ptSF=1.4085; + ptSF=1.0; + ptSF_err = TMath::Sqrt(pow(0.035,2)+pow(0.5*(0.087+0.085),2)); + //ptSF_err = 0.2020; + } + else if(TMath::Abs(eta)<2.853){ + //ptSF=1.9909; + ptSF=1.0; + ptSF_err = TMath::Sqrt(pow(0.035,2)+pow(0.5*(0.087+0.085),2)); + //ptSF_err = 0.5684; + } + else if(TMath::Abs(eta)<2.964){ + //ptSF=2.2923; + ptSF=1.0; + ptSF_err = TMath::Sqrt(pow(0.035,2)+pow(0.5*(0.087+0.085),2)); + //ptSF_err = 0.3743; + } + else if(TMath::Abs(eta)<3.139){ + //ptSF=1.2696; + ptSF=1.0; + ptSF_err = TMath::Sqrt(pow(0.035,2)+pow(0.5*(0.087+0.085),2)); + //ptSF_err = 0.1089; + } + else{ + ptSF=1.0; + ptSF_err = TMath::Sqrt(pow(0.127,2)+pow(0.5*(0.155+0.153),2)); + } res[0] = TMath::Max((Float_t)0.,(Float_t)(genjpt+(ptSF)*(pt-genjpt)))/pt; res[1] = TMath::Max((Float_t)0.,(Float_t)(genjpt+(ptSF-ptSF_err)*(pt-genjpt)))/pt; diff --git a/test/ttbar/TTbarEventAnalysis.h b/test/ttbar/TTbarEventAnalysis.h index b49936e731e..d07e103dd6e 100644 --- a/test/ttbar/TTbarEventAnalysis.h +++ b/test/ttbar/TTbarEventAnalysis.h @@ -42,7 +42,11 @@ class MyEventInfoBranches_t Float_t ttbar_w[1095]; Int_t nJet; Float_t Jet_pt[100],Jet_genpt[100],Jet_area[100],Jet_jes[100],Jet_eta[100],Jet_phi[100],Jet_mass[100]; - Float_t Jet_Svx[100],Jet_CombIVF[100],Jet_Proba[100],Jet_Ip2P[100],Jet_DeepCSVBDisc[100],Jet_DeepFlavourBDisc[100]; + Float_t Jet_Svx[100],Jet_CombIVF[100],Jet_Proba[100],Jet_Ip2P[100]; + Float_t Jet_DeepCSVb[100], Jet_DeepCSVc[100], Jet_DeepCSVl[100], Jet_DeepCSVbN[100], Jet_DeepCSVcN[100], Jet_DeepCSVlN[100]; + Float_t Jet_DeepCSVBDisc[100],Jet_DeepCSVBDiscN[100],Jet_DeepCSVCvsLDisc[100],Jet_DeepCSVCvsLDiscN[100],Jet_DeepCSVCvsBDisc[100],Jet_DeepCSVCvsBDiscN[100]; + Float_t Jet_DeepFlavourBDisc[100], Jet_DeepFlavourCvsLDisc[100], Jet_DeepFlavourCvsBDisc[100]; + Float_t Jet_DeepFlavourB[100]; Int_t Jet_nseltracks[100]; Int_t Jet_flavour[100]; }; @@ -52,18 +56,18 @@ class TTbarEventAnalysis public: TTbarEventAnalysis() : tmvaReader_(0), - readTTJetsGenWeights_(false), + readTTJetsGenWeights_(true), puWgtGr_(0),puWgtDownGr_(0),puWgtUpGr_(0) { - //jet uncertainty parameterization - TString jecUncUrl("${CMSSW_BASE}/src/RecoBTag/PerformanceMeasurements/test/ttbar/data/Autumn18_V8_MC_Uncertainty_AK4PF.txt"); - gSystem->ExpandPathName(jecUncUrl); - jecUnc_ = new JetCorrectionUncertainty(jecUncUrl.Data()); + //jet uncertainty parameterization + TString jecUncUrl("${CMSSW_BASE}/src/RecoBTag/PerformanceMeasurements/test/ttbar/data/Autumn18_V8_MC_Uncertainty_AK4PF.txt"); + gSystem->ExpandPathName(jecUncUrl); + jecUnc_ = new JetCorrectionUncertainty(jecUncUrl.Data()); - //pileup weights - TString stdTarget("${CMSSW_BASE}/src/RecoBTag/PerformanceMeasurements/test/ttbar/data/pileupWgts.root"); - gSystem->ExpandPathName(stdTarget); - SetPUWeightTarget(stdTarget,""); + //pileup weights + TString stdTarget("${CMSSW_BASE}/src/RecoBTag/PerformanceMeasurements/test/ttbar/data/pileupWgts.root"); + gSystem->ExpandPathName(stdTarget); + SetPUWeightTarget(stdTarget,""); } ~TTbarEventAnalysis(){} void setReadTTJetsGenWeights(bool readTTJetsGenWeights) { readTTJetsGenWeights_=readTTJetsGenWeights; } @@ -73,11 +77,11 @@ class TTbarEventAnalysis void prepareOutput(TString outFile); void processFile(TString inFile,TH1F *xsecWgt,Bool_t isData); void finalizeOutput(); - void SetPUWeightTarget(TString targetFile,TString sampleName){ + void SetPUWeightTarget(TString targetFile,TString sampleName){ TFile *fIn=TFile::Open(targetFile); - if(fIn){ - std::string nom( "puwgts_nom"); - std::string up( "puwgts_down"); + if(fIn){ + std::string nom("puwgts_nom"); + std::string up("puwgts_down"); std::string down("puwgts_up"); if(!sampleName.IsNull()){ nom.append("_"); @@ -86,15 +90,18 @@ class TTbarEventAnalysis up.append(sampleName); down.append("_"); down.append(sampleName); + }else{ + std::cout << "Warning: PUWeight target histogram " << sampleName << " does not exist. Check naming convention of samples matches that of PU histograms." << std::endl; } - puWgtGr_ = (TGraph *)fIn->Get(nom.c_str()); - puWgtDownGr_ = (TGraph *)fIn->Get(down.c_str()); - puWgtUpGr_ = (TGraph *)fIn->Get(up.c_str()); - }else{ - std::cout << "Unable to find data/pileupWgts.root, no PU reweighting will be applied" << std::endl; - } - fIn->Close(); - } + + puWgtGr_ = (TGraph *)fIn->Get(nom.c_str()); + puWgtDownGr_ = (TGraph *)fIn->Get(down.c_str()); + puWgtUpGr_ = (TGraph *)fIn->Get(up.c_str()); + }else{ + std::cout << "Unable to find data/pileupWgts.root, no PU reweighting will be applied" << std::endl; + } + fIn->Close(); + } private: JetCorrectionUncertainty *jecUnc_; @@ -131,6 +138,10 @@ class TTbarEventAnalysis Float_t j2ll_deta_,j2ll_dphi_; Float_t kinDisc_[5]; Float_t jp_[2],svhe_[2],csv_[2]; + Float_t DeepCSVb_[2], DeepCSVc_[2], DeepCSVl_[2], DeepCSVbb_[2], DeepCSVcc_[2], DeepCSVbN_[2], DeepCSVcN_[2], DeepCSVlN_[2], DeepCSVbbN_[2], DeepCSVccN_[2], DeepCSVbP_[2], DeepCSVcP_[2], DeepCSVlP_[2], DeepCSVbbP_[2], DeepCSVccP_[2]; + Float_t DeepCSVBDisc_[2], DeepCSVBDiscN_[2], DeepCSVBDiscP_[2], DeepCSVCvsLDisc_[2], DeepCSVCvsLDiscN_[2], DeepCSVCvsLDiscP_[2], DeepCSVCvsBDisc_[2], DeepCSVCvsBDiscN_[2], DeepCSVCvsBDiscP_[2]; + Float_t DeepFlavourBDisc_[2], DeepFlavourCvsLDisc_[2], DeepFlavourCvsBDisc_[2]; + Float_t DeepFlavourB_[2], DeepFlavourBB_[2], DeepFlavourLEPB_[2]; std::vector > triggerBits_; TTree *kinTree_,*ftmTree_; std::map histos_; diff --git a/test/ttbar/twoTag.py b/test/ttbar/twoTag.py index cf38fecdab1..8481444551e 100644 --- a/test/ttbar/twoTag.py +++ b/test/ttbar/twoTag.py @@ -36,7 +36,7 @@ print " ", for sample in samples: histName="emu_"+tag+syst+"/emu_"+tag+syst+"_"+sample - print histName + #print histName hist=tFile.Get(histName) if sample==samples[0]: histBkg=hist.Clone("histBkg") @@ -73,7 +73,7 @@ #print tot - #print histBkg.Integral(),hist.Integral(),hist.Integral()/histBkg.Integral() + print histBkg.Integral(),hist.Integral(),hist.Integral()/histBkg.Integral() for name in tot: if name.find("mc")!=-1: From 7e36a029a757b9469654d39bd91c5269236f9370 Mon Sep 17 00:00:00 2001 From: capalmer Date: Mon, 8 Apr 2019 17:13:49 +0200 Subject: [PATCH 3/7] add some notes on two tag --- test/ttbar/README.md | 17 +++++++++++++++++ 1 file changed, 17 insertions(+) diff --git a/test/ttbar/README.md b/test/ttbar/README.md index 2df18923b4c..06e5e63d977 100644 --- a/test/ttbar/README.md +++ b/test/ttbar/README.md @@ -56,6 +56,23 @@ This script uses the output plots for MET to fit a common scale factor for DY. M and used to fit the MET observed after the requirement of at least two jets. The scale factor can be applied to all channels as it reflects possible mismodelling in the jet multiplicity of the DY MC. +### Running local analysis for 2TagCount +After b-tag analyzer trees are produced, the first two steps of the 2TagCount analysis are the same. The first is to run the analysis and make histograms (and trees): +``` +python runTTbarAnalysis.py -i /store/group/phys_btag/Commissioning/TTbar/Moriond19_2018_StructuredDir/ --json data/samples_Run2018.json -o /tmp/MYUSERNAME/Moriond19_Run2018_april5 -n 70 +``` + +Re-organize plots and make data, mc comparison canvases (inclusive 2TagCount plots are in): +``` +python plotter.py -i /tmp/MYUSERNAME/Moriond19_Run2018_april5 --json data/samples_Run2018.json --lumi 58000 +``` + +Take the 2TagCount histograms and compute efficiency scale factor comparing the following: 1) b-tagging effiency for each working point for events with two b-jets at truth level (MC eff) and 2) data minus non-double b-jet MC divided by the number of MC events with two matched particle level b-jets: +``` +python twoTag.py /tmp/MYUSERNAME/Moriond19_Run2018_april5/plots/plotter.root +``` + + ### Performance analyses #### Templated methods From 333d329b7089a92126f662ad7e10d91baea32329 Mon Sep 17 00:00:00 2001 From: capalmer Date: Mon, 15 Apr 2019 17:13:45 +0200 Subject: [PATCH 4/7] include latest JECs --- .../data/Autumn18_V8_MC_Uncertainty_AK4PF.txt | 41 +++++++++++++++++++ 1 file changed, 41 insertions(+) create mode 100644 test/ttbar/data/Autumn18_V8_MC_Uncertainty_AK4PF.txt diff --git a/test/ttbar/data/Autumn18_V8_MC_Uncertainty_AK4PF.txt b/test/ttbar/data/Autumn18_V8_MC_Uncertainty_AK4PF.txt new file mode 100644 index 00000000000..751d4b9d749 --- /dev/null +++ b/test/ttbar/data/Autumn18_V8_MC_Uncertainty_AK4PF.txt @@ -0,0 +1,41 @@ +{1 JetEta 1 JetPt "" Correction Uncertainty} +-5.4 -5.0 150 9.0 0.1007 0.1007 11.0 0.0881 0.0881 13.5 0.0781 0.0781 16.5 0.0712 0.0712 19.5 0.0671 0.0671 22.5 0.0646 0.0646 26.0 0.0627 0.0627 30.0 0.0609 0.0609 34.5 0.0594 0.0594 40.0 0.0566 0.0566 46.0 0.0545 0.0545 52.5 0.0526 0.0526 60.0 0.0514 0.0514 69.0 0.0507 0.0507 79.0 0.0499 0.0499 90.5 0.0493 0.0493 105.5 0.0493 0.0493 123.5 0.0498 0.0498 143.0 0.0507 0.0507 163.5 0.0519 0.0519 185.0 0.0530 0.0530 208.0 0.0546 0.0546 232.5 0.0561 0.0561 258.5 0.0569 0.0569 286.0 0.0579 0.0579 331.0 0.0593 0.0593 396.0 0.0611 0.0611 468.5 0.0629 0.0629 549.5 0.0645 0.0645 639.0 0.0661 0.0661 738.0 0.0676 0.0676 847.5 0.0690 0.0690 968.5 0.0704 0.0704 1102.0 0.0718 0.0718 1249.5 0.0731 0.0731 1412.0 0.0743 0.0743 1590.5 0.0755 0.0755 1787.0 0.0767 0.0767 2003.0 0.0779 0.0779 2241.0 0.0790 0.0790 2503.0 0.0801 0.0801 2790.5 0.0812 0.0812 3107.0 0.0823 0.0823 3455.0 0.0833 0.0833 3837.0 0.0845 0.0845 4257.0 0.0857 0.0857 4719.0 0.0869 0.0869 5226.5 0.0881 0.0881 5784.0 0.0892 0.0892 6538.0 0.0905 0.0905 +-5.0 -4.4 150 9.0 0.0922 0.0922 11.0 0.0820 0.0820 13.5 0.0739 0.0739 16.5 0.0683 0.0683 19.5 0.0652 0.0652 22.5 0.0631 0.0631 26.0 0.0617 0.0617 30.0 0.0603 0.0603 34.5 0.0589 0.0589 40.0 0.0564 0.0564 46.0 0.0542 0.0542 52.5 0.0525 0.0525 60.0 0.0514 0.0514 69.0 0.0507 0.0507 79.0 0.0505 0.0505 90.5 0.0509 0.0509 105.5 0.0518 0.0518 123.5 0.0531 0.0531 143.0 0.0540 0.0540 163.5 0.0551 0.0551 185.0 0.0562 0.0562 208.0 0.0576 0.0576 232.5 0.0590 0.0590 258.5 0.0598 0.0598 286.0 0.0607 0.0607 331.0 0.0621 0.0621 396.0 0.0638 0.0638 468.5 0.0655 0.0655 549.5 0.0671 0.0671 639.0 0.0686 0.0686 738.0 0.0700 0.0700 847.5 0.0714 0.0714 968.5 0.0728 0.0728 1102.0 0.0741 0.0741 1249.5 0.0753 0.0753 1412.0 0.0765 0.0765 1590.5 0.0777 0.0777 1787.0 0.0789 0.0789 2003.0 0.0800 0.0800 2241.0 0.0811 0.0811 2503.0 0.0822 0.0822 2790.5 0.0832 0.0832 3107.0 0.0843 0.0843 3455.0 0.0853 0.0853 3837.0 0.0865 0.0865 4257.0 0.0877 0.0877 4719.0 0.0888 0.0888 5226.5 0.0900 0.0900 5784.0 0.0911 0.0911 6538.0 0.0924 0.0924 +-4.4 -4.0 150 9.0 0.1309 0.1309 11.0 0.1111 0.1111 13.5 0.0946 0.0946 16.5 0.0817 0.0817 19.5 0.0738 0.0738 22.5 0.0685 0.0685 26.0 0.0647 0.0647 30.0 0.0623 0.0623 34.5 0.0604 0.0604 40.0 0.0578 0.0578 46.0 0.0555 0.0555 52.5 0.0539 0.0539 60.0 0.0526 0.0526 69.0 0.0517 0.0517 79.0 0.0512 0.0512 90.5 0.0511 0.0511 105.5 0.0513 0.0513 123.5 0.0520 0.0520 143.0 0.0528 0.0528 163.5 0.0539 0.0539 185.0 0.0547 0.0547 208.0 0.0556 0.0556 232.5 0.0561 0.0561 258.5 0.0561 0.0561 286.0 0.0562 0.0562 331.0 0.0563 0.0563 396.0 0.0566 0.0566 468.5 0.0569 0.0569 549.5 0.0573 0.0573 639.0 0.0577 0.0577 738.0 0.0581 0.0581 847.5 0.0586 0.0586 968.5 0.0590 0.0590 1102.0 0.0596 0.0596 1249.5 0.0601 0.0601 1412.0 0.0606 0.0606 1590.5 0.0612 0.0612 1787.0 0.0618 0.0618 2003.0 0.0624 0.0624 2241.0 0.0630 0.0630 2503.0 0.0636 0.0636 2790.5 0.0642 0.0642 3107.0 0.0649 0.0649 3455.0 0.0655 0.0655 3837.0 0.0661 0.0661 4257.0 0.0666 0.0666 4719.0 0.0672 0.0672 5226.5 0.0678 0.0678 5784.0 0.0685 0.0685 6538.0 0.0693 0.0693 +-4.0 -3.5 150 9.0 0.1462 0.1462 11.0 0.1256 0.1256 13.5 0.1091 0.1091 16.5 0.0961 0.0961 19.5 0.0882 0.0882 22.5 0.0838 0.0838 26.0 0.0796 0.0796 30.0 0.0780 0.0780 34.5 0.0761 0.0761 40.0 0.0732 0.0732 46.0 0.0712 0.0712 52.5 0.0693 0.0693 60.0 0.0684 0.0684 69.0 0.0672 0.0672 79.0 0.0665 0.0665 90.5 0.0661 0.0661 105.5 0.0662 0.0662 123.5 0.0667 0.0667 143.0 0.0676 0.0676 163.5 0.0689 0.0689 185.0 0.0703 0.0703 208.0 0.0718 0.0718 232.5 0.0735 0.0735 258.5 0.0751 0.0751 286.0 0.0769 0.0769 331.0 0.0786 0.0786 396.0 0.0791 0.0791 468.5 0.0792 0.0792 549.5 0.0793 0.0793 639.0 0.0794 0.0794 738.0 0.0795 0.0795 847.5 0.0796 0.0796 968.5 0.0797 0.0797 1102.0 0.0799 0.0799 1249.5 0.0800 0.0800 1412.0 0.0801 0.0801 1590.5 0.0803 0.0803 1787.0 0.0804 0.0804 2003.0 0.0805 0.0805 2241.0 0.0807 0.0807 2503.0 0.0808 0.0808 2790.5 0.0809 0.0809 3107.0 0.0811 0.0811 3455.0 0.0812 0.0812 3837.0 0.0814 0.0814 4257.0 0.0815 0.0815 4719.0 0.0816 0.0816 5226.5 0.0817 0.0817 5784.0 0.0818 0.0818 6538.0 0.0820 0.0820 +-3.5 -3.0 150 9.0 0.1457 0.1457 11.0 0.1296 0.1296 13.5 0.1167 0.1167 16.5 0.1069 0.1069 19.5 0.1005 0.1005 22.5 0.0968 0.0968 26.0 0.0944 0.0944 30.0 0.0936 0.0936 34.5 0.0920 0.0920 40.0 0.0899 0.0899 46.0 0.0887 0.0887 52.5 0.0876 0.0876 60.0 0.0868 0.0868 69.0 0.0862 0.0862 79.0 0.0857 0.0857 90.5 0.0854 0.0854 105.5 0.0853 0.0853 123.5 0.0854 0.0854 143.0 0.0857 0.0857 163.5 0.0861 0.0861 185.0 0.0866 0.0866 208.0 0.0872 0.0872 232.5 0.0878 0.0878 258.5 0.0885 0.0885 286.0 0.0892 0.0892 331.0 0.0904 0.0904 396.0 0.0921 0.0921 468.5 0.0936 0.0936 549.5 0.0942 0.0942 639.0 0.0942 0.0942 738.0 0.0943 0.0943 847.5 0.0944 0.0944 968.5 0.0945 0.0945 1102.0 0.0947 0.0947 1249.5 0.0948 0.0948 1412.0 0.0949 0.0949 1590.5 0.0950 0.0950 1787.0 0.0951 0.0951 2003.0 0.0952 0.0952 2241.0 0.0953 0.0953 2503.0 0.0954 0.0954 2790.5 0.0955 0.0955 3107.0 0.0956 0.0956 3455.0 0.0957 0.0957 3837.0 0.0958 0.0958 4257.0 0.0959 0.0959 4719.0 0.0961 0.0961 5226.5 0.0962 0.0962 5784.0 0.0963 0.0963 6538.0 0.0964 0.0964 +-3.0 -2.8 150 9.0 0.2457 0.2457 11.0 0.2411 0.2411 13.5 0.2386 0.2386 16.5 0.2375 0.2375 19.5 0.2403 0.2403 22.5 0.2370 0.2370 26.0 0.2379 0.2379 30.0 0.2360 0.2360 34.5 0.2300 0.2300 40.0 0.2245 0.2245 46.0 0.2197 0.2197 52.5 0.2159 0.2159 60.0 0.2122 0.2122 69.0 0.2094 0.2094 79.0 0.2081 0.2081 90.5 0.2075 0.2075 105.5 0.2071 0.2071 123.5 0.2091 0.2091 143.0 0.2128 0.2128 163.5 0.2174 0.2174 185.0 0.2229 0.2229 208.0 0.2295 0.2295 232.5 0.2372 0.2372 258.5 0.2457 0.2457 286.0 0.2549 0.2549 331.0 0.2706 0.2706 396.0 0.2726 0.2726 468.5 0.2659 0.2659 549.5 0.2658 0.2658 639.0 0.2663 0.2663 738.0 0.2667 0.2667 847.5 0.2667 0.2667 968.5 0.2668 0.2668 1102.0 0.2668 0.2668 1249.5 0.2669 0.2669 1412.0 0.2669 0.2669 1590.5 0.2670 0.2670 1787.0 0.2670 0.2670 2003.0 0.2670 0.2670 2241.0 0.2671 0.2671 2503.0 0.2671 0.2671 2790.5 0.2672 0.2672 3107.0 0.2673 0.2673 3455.0 0.2673 0.2673 3837.0 0.2674 0.2674 4257.0 0.2675 0.2675 4719.0 0.2287 0.2287 5226.5 0.2063 0.2063 5784.0 0.2036 0.2036 6538.0 0.2038 0.2038 +-2.8 -2.6 150 9.0 0.1827 0.1827 11.0 0.1772 0.1772 13.5 0.1738 0.1738 16.5 0.1713 0.1713 19.5 0.1717 0.1717 22.5 0.1709 0.1709 26.0 0.1720 0.1720 30.0 0.1722 0.1722 34.5 0.1667 0.1667 40.0 0.1626 0.1626 46.0 0.1578 0.1578 52.5 0.1526 0.1526 60.0 0.1489 0.1489 69.0 0.1448 0.1448 79.0 0.1410 0.1410 90.5 0.1380 0.1380 105.5 0.1352 0.1352 123.5 0.1336 0.1336 143.0 0.1330 0.1330 163.5 0.1332 0.1332 185.0 0.1342 0.1342 208.0 0.1357 0.1357 232.5 0.1375 0.1375 258.5 0.1396 0.1396 286.0 0.1424 0.1424 331.0 0.1485 0.1485 396.0 0.1596 0.1596 468.5 0.1653 0.1653 549.5 0.1651 0.1651 639.0 0.1659 0.1659 738.0 0.1666 0.1666 847.5 0.1671 0.1671 968.5 0.1673 0.1673 1102.0 0.1673 0.1673 1249.5 0.1674 0.1674 1412.0 0.1675 0.1675 1590.5 0.1675 0.1675 1787.0 0.1676 0.1676 2003.0 0.1677 0.1677 2241.0 0.1678 0.1678 2503.0 0.1678 0.1678 2790.5 0.1679 0.1679 3107.0 0.1680 0.1680 3455.0 0.1681 0.1681 3837.0 0.1683 0.1683 4257.0 0.1684 0.1684 4719.0 0.1686 0.1686 5226.5 0.1621 0.1621 5784.0 0.1459 0.1459 6538.0 0.1468 0.1468 +-2.6 -2.4 150 9.0 0.1167 0.1167 11.0 0.1023 0.1023 13.5 0.0890 0.0890 16.5 0.0793 0.0793 19.5 0.0736 0.0736 22.5 0.0701 0.0701 26.0 0.0679 0.0679 30.0 0.0672 0.0672 34.5 0.0624 0.0624 40.0 0.0571 0.0571 46.0 0.0525 0.0525 52.5 0.0485 0.0485 60.0 0.0449 0.0449 69.0 0.0418 0.0418 79.0 0.0394 0.0394 90.5 0.0380 0.0380 105.5 0.0375 0.0375 123.5 0.0384 0.0384 143.0 0.0404 0.0404 163.5 0.0432 0.0432 185.0 0.0465 0.0465 208.0 0.0500 0.0500 232.5 0.0537 0.0537 258.5 0.0575 0.0575 286.0 0.0611 0.0611 331.0 0.0669 0.0669 396.0 0.0748 0.0748 468.5 0.0829 0.0829 549.5 0.0917 0.0917 639.0 0.0957 0.0957 738.0 0.0944 0.0944 847.5 0.0925 0.0925 968.5 0.0911 0.0911 1102.0 0.0903 0.0903 1249.5 0.0904 0.0904 1412.0 0.0905 0.0905 1590.5 0.0906 0.0906 1787.0 0.0907 0.0907 2003.0 0.0908 0.0908 2241.0 0.0909 0.0909 2503.0 0.0910 0.0910 2790.5 0.0912 0.0912 3107.0 0.0913 0.0913 3455.0 0.0914 0.0914 3837.0 0.0916 0.0916 4257.0 0.0918 0.0918 4719.0 0.0920 0.0920 5226.5 0.0922 0.0922 5784.0 0.0924 0.0924 6538.0 0.0909 0.0909 +-2.4 -2.2 150 9.0 0.0710 0.0710 11.0 0.0665 0.0665 13.5 0.0595 0.0595 16.5 0.0537 0.0537 19.5 0.0497 0.0497 22.5 0.0466 0.0466 26.0 0.0441 0.0441 30.0 0.0422 0.0422 34.5 0.0396 0.0396 40.0 0.0366 0.0366 46.0 0.0340 0.0340 52.5 0.0318 0.0318 60.0 0.0298 0.0298 69.0 0.0278 0.0278 79.0 0.0261 0.0261 90.5 0.0246 0.0246 105.5 0.0232 0.0232 123.5 0.0221 0.0221 143.0 0.0214 0.0214 163.5 0.0210 0.0210 185.0 0.0208 0.0208 208.0 0.0209 0.0209 232.5 0.0211 0.0211 258.5 0.0214 0.0214 286.0 0.0219 0.0219 331.0 0.0226 0.0226 396.0 0.0238 0.0238 468.5 0.0251 0.0251 549.5 0.0266 0.0266 639.0 0.0280 0.0280 738.0 0.0294 0.0294 847.5 0.0309 0.0309 968.5 0.0324 0.0324 1102.0 0.0339 0.0339 1249.5 0.0353 0.0353 1412.0 0.0361 0.0361 1590.5 0.0367 0.0367 1787.0 0.0374 0.0374 2003.0 0.0381 0.0381 2241.0 0.0389 0.0389 2503.0 0.0397 0.0397 2790.5 0.0406 0.0406 3107.0 0.0416 0.0416 3455.0 0.0426 0.0426 3837.0 0.0437 0.0437 4257.0 0.0449 0.0449 4719.0 0.0461 0.0461 5226.5 0.0474 0.0474 5784.0 0.0487 0.0487 6538.0 0.0503 0.0503 +-2.2 -2.0 150 9.0 0.0594 0.0594 11.0 0.0558 0.0558 13.5 0.0524 0.0524 16.5 0.0495 0.0495 19.5 0.0477 0.0477 22.5 0.0467 0.0467 26.0 0.0459 0.0459 30.0 0.0458 0.0458 34.5 0.0442 0.0442 40.0 0.0415 0.0415 46.0 0.0388 0.0388 52.5 0.0365 0.0365 60.0 0.0343 0.0343 69.0 0.0320 0.0320 79.0 0.0299 0.0299 90.5 0.0279 0.0279 105.5 0.0260 0.0260 123.5 0.0243 0.0243 143.0 0.0230 0.0230 163.5 0.0221 0.0221 185.0 0.0216 0.0216 208.0 0.0213 0.0213 232.5 0.0213 0.0213 258.5 0.0214 0.0214 286.0 0.0218 0.0218 331.0 0.0225 0.0225 396.0 0.0240 0.0240 468.5 0.0257 0.0257 549.5 0.0276 0.0276 639.0 0.0297 0.0297 738.0 0.0319 0.0319 847.5 0.0342 0.0342 968.5 0.0352 0.0352 1102.0 0.0357 0.0357 1249.5 0.0363 0.0363 1412.0 0.0369 0.0369 1590.5 0.0375 0.0375 1787.0 0.0381 0.0381 2003.0 0.0389 0.0389 2241.0 0.0396 0.0396 2503.0 0.0405 0.0405 2790.5 0.0413 0.0413 3107.0 0.0423 0.0423 3455.0 0.0432 0.0432 3837.0 0.0443 0.0443 4257.0 0.0455 0.0455 4719.0 0.0467 0.0467 5226.5 0.0479 0.0479 5784.0 0.0491 0.0491 6538.0 0.0506 0.0506 +-2.0 -1.8 150 9.0 0.0513 0.0513 11.0 0.0462 0.0462 13.5 0.0411 0.0411 16.5 0.0374 0.0374 19.5 0.0352 0.0352 22.5 0.0336 0.0336 26.0 0.0324 0.0324 30.0 0.0322 0.0322 34.5 0.0308 0.0308 40.0 0.0291 0.0291 46.0 0.0278 0.0278 52.5 0.0267 0.0267 60.0 0.0257 0.0257 69.0 0.0248 0.0248 79.0 0.0242 0.0242 90.5 0.0237 0.0237 105.5 0.0233 0.0233 123.5 0.0232 0.0232 143.0 0.0233 0.0233 163.5 0.0235 0.0235 185.0 0.0238 0.0238 208.0 0.0242 0.0242 232.5 0.0246 0.0246 258.5 0.0251 0.0251 286.0 0.0256 0.0256 331.0 0.0264 0.0264 396.0 0.0275 0.0275 468.5 0.0287 0.0287 549.5 0.0299 0.0299 639.0 0.0312 0.0312 738.0 0.0325 0.0325 847.5 0.0338 0.0338 968.5 0.0352 0.0352 1102.0 0.0360 0.0360 1249.5 0.0364 0.0364 1412.0 0.0368 0.0368 1590.5 0.0373 0.0373 1787.0 0.0378 0.0378 2003.0 0.0383 0.0383 2241.0 0.0388 0.0388 2503.0 0.0394 0.0394 2790.5 0.0400 0.0400 3107.0 0.0406 0.0406 3455.0 0.0413 0.0413 3837.0 0.0420 0.0420 4257.0 0.0428 0.0428 4719.0 0.0436 0.0436 5226.5 0.0445 0.0445 5784.0 0.0453 0.0453 6538.0 0.0464 0.0464 +-1.8 -1.6 150 9.0 0.0546 0.0546 11.0 0.0492 0.0492 13.5 0.0438 0.0438 16.5 0.0396 0.0396 19.5 0.0368 0.0368 22.5 0.0348 0.0348 26.0 0.0333 0.0333 30.0 0.0328 0.0328 34.5 0.0313 0.0313 40.0 0.0295 0.0295 46.0 0.0279 0.0279 52.5 0.0268 0.0268 60.0 0.0258 0.0258 69.0 0.0249 0.0249 79.0 0.0242 0.0242 90.5 0.0236 0.0236 105.5 0.0233 0.0233 123.5 0.0231 0.0231 143.0 0.0232 0.0232 163.5 0.0234 0.0234 185.0 0.0237 0.0237 208.0 0.0241 0.0241 232.5 0.0246 0.0246 258.5 0.0251 0.0251 286.0 0.0256 0.0256 331.0 0.0264 0.0264 396.0 0.0275 0.0275 468.5 0.0287 0.0287 549.5 0.0299 0.0299 639.0 0.0312 0.0312 738.0 0.0325 0.0325 847.5 0.0338 0.0338 968.5 0.0351 0.0351 1102.0 0.0359 0.0359 1249.5 0.0362 0.0362 1412.0 0.0366 0.0366 1590.5 0.0370 0.0370 1787.0 0.0373 0.0373 2003.0 0.0377 0.0377 2241.0 0.0381 0.0381 2503.0 0.0385 0.0385 2790.5 0.0390 0.0390 3107.0 0.0395 0.0395 3455.0 0.0400 0.0400 3837.0 0.0406 0.0406 4257.0 0.0412 0.0412 4719.0 0.0418 0.0418 5226.5 0.0424 0.0424 5784.0 0.0431 0.0431 6538.0 0.0439 0.0439 +-1.6 -1.4 150 9.0 0.0556 0.0556 11.0 0.0493 0.0493 13.5 0.0424 0.0424 16.5 0.0366 0.0366 19.5 0.0325 0.0325 22.5 0.0294 0.0294 26.0 0.0269 0.0269 30.0 0.0258 0.0258 34.5 0.0238 0.0238 40.0 0.0219 0.0219 46.0 0.0205 0.0205 52.5 0.0194 0.0194 60.0 0.0187 0.0187 69.0 0.0183 0.0183 79.0 0.0181 0.0181 90.5 0.0182 0.0182 105.5 0.0186 0.0186 123.5 0.0193 0.0193 143.0 0.0201 0.0201 163.5 0.0209 0.0209 185.0 0.0217 0.0217 208.0 0.0226 0.0226 232.5 0.0234 0.0234 258.5 0.0242 0.0242 286.0 0.0250 0.0250 331.0 0.0263 0.0263 396.0 0.0278 0.0278 468.5 0.0293 0.0293 549.5 0.0307 0.0307 639.0 0.0322 0.0322 738.0 0.0336 0.0336 847.5 0.0350 0.0350 968.5 0.0364 0.0364 1102.0 0.0374 0.0374 1249.5 0.0378 0.0378 1412.0 0.0381 0.0381 1590.5 0.0385 0.0385 1787.0 0.0388 0.0388 2003.0 0.0392 0.0392 2241.0 0.0396 0.0396 2503.0 0.0399 0.0399 2790.5 0.0403 0.0403 3107.0 0.0407 0.0407 3455.0 0.0410 0.0410 3837.0 0.0415 0.0415 4257.0 0.0419 0.0419 4719.0 0.0424 0.0424 5226.5 0.0429 0.0429 5784.0 0.0434 0.0434 6538.0 0.0440 0.0440 +-1.4 -1.2 150 9.0 0.0591 0.0591 11.0 0.0524 0.0524 13.5 0.0446 0.0446 16.5 0.0376 0.0376 19.5 0.0325 0.0325 22.5 0.0286 0.0286 26.0 0.0253 0.0253 30.0 0.0239 0.0239 34.5 0.0217 0.0217 40.0 0.0197 0.0197 46.0 0.0182 0.0182 52.5 0.0171 0.0171 60.0 0.0163 0.0163 69.0 0.0158 0.0158 79.0 0.0154 0.0154 90.5 0.0153 0.0153 105.5 0.0153 0.0153 123.5 0.0155 0.0155 143.0 0.0158 0.0158 163.5 0.0162 0.0162 185.0 0.0165 0.0165 208.0 0.0168 0.0168 232.5 0.0171 0.0171 258.5 0.0174 0.0174 286.0 0.0177 0.0177 331.0 0.0181 0.0181 396.0 0.0186 0.0186 468.5 0.0191 0.0191 549.5 0.0196 0.0196 639.0 0.0201 0.0201 738.0 0.0206 0.0206 847.5 0.0212 0.0212 968.5 0.0217 0.0217 1102.0 0.0223 0.0223 1249.5 0.0230 0.0230 1412.0 0.0238 0.0238 1590.5 0.0247 0.0247 1787.0 0.0255 0.0255 2003.0 0.0263 0.0263 2241.0 0.0272 0.0272 2503.0 0.0280 0.0280 2790.5 0.0289 0.0289 3107.0 0.0297 0.0297 3455.0 0.0304 0.0304 3837.0 0.0310 0.0310 4257.0 0.0316 0.0316 4719.0 0.0322 0.0322 5226.5 0.0329 0.0329 5784.0 0.0335 0.0335 6538.0 0.0447 0.0447 +-1.2 -1.0 150 9.0 0.0601 0.0601 11.0 0.0524 0.0524 13.5 0.0437 0.0437 16.5 0.0363 0.0363 19.5 0.0310 0.0310 22.5 0.0271 0.0271 26.0 0.0241 0.0241 30.0 0.0229 0.0229 34.5 0.0208 0.0208 40.0 0.0189 0.0189 46.0 0.0176 0.0176 52.5 0.0166 0.0166 60.0 0.0160 0.0160 69.0 0.0156 0.0156 79.0 0.0153 0.0153 90.5 0.0152 0.0152 105.5 0.0153 0.0153 123.5 0.0155 0.0155 143.0 0.0157 0.0157 163.5 0.0160 0.0160 185.0 0.0164 0.0164 208.0 0.0167 0.0167 232.5 0.0170 0.0170 258.5 0.0173 0.0173 286.0 0.0175 0.0175 331.0 0.0180 0.0180 396.0 0.0185 0.0185 468.5 0.0190 0.0190 549.5 0.0195 0.0195 639.0 0.0200 0.0200 738.0 0.0205 0.0205 847.5 0.0210 0.0210 968.5 0.0215 0.0215 1102.0 0.0221 0.0221 1249.5 0.0227 0.0227 1412.0 0.0235 0.0235 1590.5 0.0243 0.0243 1787.0 0.0250 0.0250 2003.0 0.0258 0.0258 2241.0 0.0266 0.0266 2503.0 0.0273 0.0273 2790.5 0.0281 0.0281 3107.0 0.0288 0.0288 3455.0 0.0295 0.0295 3837.0 0.0303 0.0303 4257.0 0.0308 0.0308 4719.0 0.0313 0.0313 5226.5 0.0318 0.0318 5784.0 0.0323 0.0323 6538.0 0.0334 0.0334 +-1.0 -0.8 150 9.0 0.0564 0.0564 11.0 0.0489 0.0489 13.5 0.0406 0.0406 16.5 0.0334 0.0334 19.5 0.0283 0.0283 22.5 0.0246 0.0246 26.0 0.0216 0.0216 30.0 0.0203 0.0203 34.5 0.0181 0.0181 40.0 0.0161 0.0161 46.0 0.0146 0.0146 52.5 0.0135 0.0135 60.0 0.0127 0.0127 69.0 0.0122 0.0122 79.0 0.0118 0.0118 90.5 0.0116 0.0116 105.5 0.0116 0.0116 123.5 0.0118 0.0118 143.0 0.0120 0.0120 163.5 0.0123 0.0123 185.0 0.0127 0.0127 208.0 0.0130 0.0130 232.5 0.0133 0.0133 258.5 0.0136 0.0136 286.0 0.0139 0.0139 331.0 0.0144 0.0144 396.0 0.0150 0.0150 468.5 0.0155 0.0155 549.5 0.0161 0.0161 639.0 0.0167 0.0167 738.0 0.0172 0.0172 847.5 0.0178 0.0178 968.5 0.0184 0.0184 1102.0 0.0190 0.0190 1249.5 0.0196 0.0196 1412.0 0.0202 0.0202 1590.5 0.0207 0.0207 1787.0 0.0213 0.0213 2003.0 0.0219 0.0219 2241.0 0.0224 0.0224 2503.0 0.0229 0.0229 2790.5 0.0235 0.0235 3107.0 0.0240 0.0240 3455.0 0.0246 0.0246 3837.0 0.0252 0.0252 4257.0 0.0258 0.0258 4719.0 0.0263 0.0263 5226.5 0.0268 0.0268 5784.0 0.0274 0.0274 6538.0 0.0338 0.0338 +-0.8 -0.6 150 9.0 0.0557 0.0557 11.0 0.0485 0.0485 13.5 0.0404 0.0404 16.5 0.0337 0.0337 19.5 0.0290 0.0290 22.5 0.0256 0.0256 26.0 0.0230 0.0230 30.0 0.0218 0.0218 34.5 0.0197 0.0197 40.0 0.0175 0.0175 46.0 0.0161 0.0161 52.5 0.0149 0.0149 60.0 0.0141 0.0141 69.0 0.0133 0.0133 79.0 0.0128 0.0128 90.5 0.0125 0.0125 105.5 0.0122 0.0122 123.5 0.0122 0.0122 143.0 0.0122 0.0122 163.5 0.0123 0.0123 185.0 0.0125 0.0125 208.0 0.0127 0.0127 232.5 0.0129 0.0129 258.5 0.0131 0.0131 286.0 0.0133 0.0133 331.0 0.0137 0.0137 396.0 0.0142 0.0142 468.5 0.0147 0.0147 549.5 0.0152 0.0152 639.0 0.0157 0.0157 738.0 0.0162 0.0162 847.5 0.0168 0.0168 968.5 0.0173 0.0173 1102.0 0.0179 0.0179 1249.5 0.0184 0.0184 1412.0 0.0189 0.0189 1590.5 0.0194 0.0194 1787.0 0.0199 0.0199 2003.0 0.0204 0.0204 2241.0 0.0209 0.0209 2503.0 0.0214 0.0214 2790.5 0.0219 0.0219 3107.0 0.0224 0.0224 3455.0 0.0229 0.0229 3837.0 0.0234 0.0234 4257.0 0.0240 0.0240 4719.0 0.0245 0.0245 5226.5 0.0250 0.0250 5784.0 0.0255 0.0255 6538.0 0.0292 0.0292 +-0.6 -0.4 150 9.0 0.0555 0.0555 11.0 0.0486 0.0486 13.5 0.0410 0.0410 16.5 0.0348 0.0348 19.5 0.0307 0.0307 22.5 0.0278 0.0278 26.0 0.0257 0.0257 30.0 0.0246 0.0246 34.5 0.0225 0.0225 40.0 0.0206 0.0206 46.0 0.0191 0.0191 52.5 0.0179 0.0179 60.0 0.0169 0.0169 69.0 0.0160 0.0160 79.0 0.0152 0.0152 90.5 0.0146 0.0146 105.5 0.0141 0.0141 123.5 0.0137 0.0137 143.0 0.0134 0.0134 163.5 0.0133 0.0133 185.0 0.0132 0.0132 208.0 0.0132 0.0132 232.5 0.0133 0.0133 258.5 0.0134 0.0134 286.0 0.0135 0.0135 331.0 0.0138 0.0138 396.0 0.0141 0.0141 468.5 0.0146 0.0146 549.5 0.0151 0.0151 639.0 0.0156 0.0156 738.0 0.0162 0.0162 847.5 0.0168 0.0168 968.5 0.0174 0.0174 1102.0 0.0181 0.0181 1249.5 0.0187 0.0187 1412.0 0.0192 0.0192 1590.5 0.0197 0.0197 1787.0 0.0202 0.0202 2003.0 0.0207 0.0207 2241.0 0.0212 0.0212 2503.0 0.0217 0.0217 2790.5 0.0222 0.0222 3107.0 0.0226 0.0226 3455.0 0.0230 0.0230 3837.0 0.0235 0.0235 4257.0 0.0239 0.0239 4719.0 0.0244 0.0244 5226.5 0.0248 0.0248 5784.0 0.0253 0.0253 6538.0 0.0270 0.0270 +-0.4 -0.2 150 9.0 0.0560 0.0560 11.0 0.0489 0.0489 13.5 0.0412 0.0412 16.5 0.0348 0.0348 19.5 0.0306 0.0306 22.5 0.0278 0.0278 26.0 0.0257 0.0257 30.0 0.0247 0.0247 34.5 0.0225 0.0225 40.0 0.0206 0.0206 46.0 0.0191 0.0191 52.5 0.0179 0.0179 60.0 0.0169 0.0169 69.0 0.0160 0.0160 79.0 0.0152 0.0152 90.5 0.0146 0.0146 105.5 0.0141 0.0141 123.5 0.0137 0.0137 143.0 0.0134 0.0134 163.5 0.0133 0.0133 185.0 0.0132 0.0132 208.0 0.0132 0.0132 232.5 0.0133 0.0133 258.5 0.0134 0.0134 286.0 0.0135 0.0135 331.0 0.0137 0.0137 396.0 0.0141 0.0141 468.5 0.0146 0.0146 549.5 0.0151 0.0151 639.0 0.0156 0.0156 738.0 0.0162 0.0162 847.5 0.0168 0.0168 968.5 0.0174 0.0174 1102.0 0.0181 0.0181 1249.5 0.0187 0.0187 1412.0 0.0192 0.0192 1590.5 0.0197 0.0197 1787.0 0.0202 0.0202 2003.0 0.0207 0.0207 2241.0 0.0212 0.0212 2503.0 0.0216 0.0216 2790.5 0.0221 0.0221 3107.0 0.0225 0.0225 3455.0 0.0229 0.0229 3837.0 0.0234 0.0234 4257.0 0.0238 0.0238 4719.0 0.0242 0.0242 5226.5 0.0247 0.0247 5784.0 0.0251 0.0251 6538.0 0.0265 0.0265 +-0.2 0.0 150 9.0 0.0555 0.0555 11.0 0.0481 0.0481 13.5 0.0401 0.0401 16.5 0.0333 0.0333 19.5 0.0288 0.0288 22.5 0.0256 0.0256 26.0 0.0233 0.0233 30.0 0.0221 0.0221 34.5 0.0200 0.0200 40.0 0.0182 0.0182 46.0 0.0168 0.0168 52.5 0.0158 0.0158 60.0 0.0150 0.0150 69.0 0.0144 0.0144 79.0 0.0139 0.0139 90.5 0.0136 0.0136 105.5 0.0133 0.0133 123.5 0.0133 0.0133 143.0 0.0133 0.0133 163.5 0.0134 0.0134 185.0 0.0136 0.0136 208.0 0.0138 0.0138 232.5 0.0140 0.0140 258.5 0.0142 0.0142 286.0 0.0144 0.0144 331.0 0.0148 0.0148 396.0 0.0153 0.0153 468.5 0.0157 0.0157 549.5 0.0162 0.0162 639.0 0.0167 0.0167 738.0 0.0172 0.0172 847.5 0.0178 0.0178 968.5 0.0183 0.0183 1102.0 0.0188 0.0188 1249.5 0.0193 0.0193 1412.0 0.0198 0.0198 1590.5 0.0204 0.0204 1787.0 0.0210 0.0210 2003.0 0.0215 0.0215 2241.0 0.0220 0.0220 2503.0 0.0225 0.0225 2790.5 0.0230 0.0230 3107.0 0.0234 0.0234 3455.0 0.0239 0.0239 3837.0 0.0244 0.0244 4257.0 0.0248 0.0248 4719.0 0.0253 0.0253 5226.5 0.0258 0.0258 5784.0 0.0262 0.0262 6538.0 0.0268 0.0268 +0.0 0.2 150 9.0 0.0556 0.0556 11.0 0.0483 0.0483 13.5 0.0403 0.0403 16.5 0.0334 0.0334 19.5 0.0288 0.0288 22.5 0.0257 0.0257 26.0 0.0234 0.0234 30.0 0.0221 0.0221 34.5 0.0200 0.0200 40.0 0.0182 0.0182 46.0 0.0168 0.0168 52.5 0.0158 0.0158 60.0 0.0150 0.0150 69.0 0.0144 0.0144 79.0 0.0139 0.0139 90.5 0.0136 0.0136 105.5 0.0133 0.0133 123.5 0.0133 0.0133 143.0 0.0133 0.0133 163.5 0.0134 0.0134 185.0 0.0136 0.0136 208.0 0.0138 0.0138 232.5 0.0140 0.0140 258.5 0.0142 0.0142 286.0 0.0144 0.0144 331.0 0.0148 0.0148 396.0 0.0153 0.0153 468.5 0.0157 0.0157 549.5 0.0162 0.0162 639.0 0.0167 0.0167 738.0 0.0172 0.0172 847.5 0.0178 0.0178 968.5 0.0183 0.0183 1102.0 0.0188 0.0188 1249.5 0.0193 0.0193 1412.0 0.0198 0.0198 1590.5 0.0204 0.0204 1787.0 0.0210 0.0210 2003.0 0.0215 0.0215 2241.0 0.0220 0.0220 2503.0 0.0225 0.0225 2790.5 0.0230 0.0230 3107.0 0.0235 0.0235 3455.0 0.0239 0.0239 3837.0 0.0244 0.0244 4257.0 0.0249 0.0249 4719.0 0.0253 0.0253 5226.5 0.0258 0.0258 5784.0 0.0263 0.0263 6538.0 0.0268 0.0268 +0.2 0.4 150 9.0 0.0563 0.0563 11.0 0.0491 0.0491 13.5 0.0413 0.0413 16.5 0.0350 0.0350 19.5 0.0307 0.0307 22.5 0.0279 0.0279 26.0 0.0259 0.0259 30.0 0.0247 0.0247 34.5 0.0226 0.0226 40.0 0.0207 0.0207 46.0 0.0191 0.0191 52.5 0.0179 0.0179 60.0 0.0169 0.0169 69.0 0.0160 0.0160 79.0 0.0153 0.0153 90.5 0.0146 0.0146 105.5 0.0141 0.0141 123.5 0.0137 0.0137 143.0 0.0134 0.0134 163.5 0.0133 0.0133 185.0 0.0132 0.0132 208.0 0.0132 0.0132 232.5 0.0133 0.0133 258.5 0.0134 0.0134 286.0 0.0135 0.0135 331.0 0.0137 0.0137 396.0 0.0141 0.0141 468.5 0.0146 0.0146 549.5 0.0151 0.0151 639.0 0.0156 0.0156 738.0 0.0162 0.0162 847.5 0.0168 0.0168 968.5 0.0174 0.0174 1102.0 0.0181 0.0181 1249.5 0.0186 0.0186 1412.0 0.0192 0.0192 1590.5 0.0197 0.0197 1787.0 0.0202 0.0202 2003.0 0.0207 0.0207 2241.0 0.0211 0.0211 2503.0 0.0216 0.0216 2790.5 0.0220 0.0220 3107.0 0.0225 0.0225 3455.0 0.0228 0.0228 3837.0 0.0233 0.0233 4257.0 0.0237 0.0237 4719.0 0.0241 0.0241 5226.5 0.0245 0.0245 5784.0 0.0250 0.0250 6538.0 0.0262 0.0262 +0.4 0.6 150 9.0 0.0556 0.0556 11.0 0.0488 0.0488 13.5 0.0411 0.0411 16.5 0.0350 0.0350 19.5 0.0308 0.0308 22.5 0.0280 0.0280 26.0 0.0258 0.0258 30.0 0.0247 0.0247 34.5 0.0226 0.0226 40.0 0.0207 0.0207 46.0 0.0191 0.0191 52.5 0.0179 0.0179 60.0 0.0169 0.0169 69.0 0.0160 0.0160 79.0 0.0152 0.0152 90.5 0.0146 0.0146 105.5 0.0141 0.0141 123.5 0.0137 0.0137 143.0 0.0134 0.0134 163.5 0.0133 0.0133 185.0 0.0132 0.0132 208.0 0.0132 0.0132 232.5 0.0133 0.0133 258.5 0.0134 0.0134 286.0 0.0135 0.0135 331.0 0.0138 0.0138 396.0 0.0141 0.0141 468.5 0.0146 0.0146 549.5 0.0151 0.0151 639.0 0.0156 0.0156 738.0 0.0162 0.0162 847.5 0.0168 0.0168 968.5 0.0174 0.0174 1102.0 0.0181 0.0181 1249.5 0.0187 0.0187 1412.0 0.0192 0.0192 1590.5 0.0197 0.0197 1787.0 0.0202 0.0202 2003.0 0.0207 0.0207 2241.0 0.0212 0.0212 2503.0 0.0217 0.0217 2790.5 0.0222 0.0222 3107.0 0.0226 0.0226 3455.0 0.0230 0.0230 3837.0 0.0235 0.0235 4257.0 0.0239 0.0239 4719.0 0.0244 0.0244 5226.5 0.0248 0.0248 5784.0 0.0253 0.0253 6538.0 0.0277 0.0277 +0.6 0.8 150 9.0 0.0558 0.0558 11.0 0.0485 0.0485 13.5 0.0403 0.0403 16.5 0.0338 0.0338 19.5 0.0290 0.0290 22.5 0.0257 0.0257 26.0 0.0232 0.0232 30.0 0.0218 0.0218 34.5 0.0196 0.0196 40.0 0.0177 0.0177 46.0 0.0161 0.0161 52.5 0.0149 0.0149 60.0 0.0141 0.0141 69.0 0.0134 0.0134 79.0 0.0129 0.0129 90.5 0.0125 0.0125 105.5 0.0123 0.0123 123.5 0.0122 0.0122 143.0 0.0122 0.0122 163.5 0.0123 0.0123 185.0 0.0125 0.0125 208.0 0.0127 0.0127 232.5 0.0129 0.0129 258.5 0.0131 0.0131 286.0 0.0133 0.0133 331.0 0.0137 0.0137 396.0 0.0142 0.0142 468.5 0.0147 0.0147 549.5 0.0152 0.0152 639.0 0.0157 0.0157 738.0 0.0162 0.0162 847.5 0.0168 0.0168 968.5 0.0173 0.0173 1102.0 0.0179 0.0179 1249.5 0.0184 0.0184 1412.0 0.0189 0.0189 1590.5 0.0194 0.0194 1787.0 0.0199 0.0199 2003.0 0.0204 0.0204 2241.0 0.0209 0.0209 2503.0 0.0214 0.0214 2790.5 0.0219 0.0219 3107.0 0.0224 0.0224 3455.0 0.0229 0.0229 3837.0 0.0234 0.0234 4257.0 0.0239 0.0239 4719.0 0.0245 0.0245 5226.5 0.0250 0.0250 5784.0 0.0254 0.0254 6538.0 0.0263 0.0263 +0.8 1.0 150 9.0 0.0589 0.0589 11.0 0.0509 0.0509 13.5 0.0421 0.0421 16.5 0.0346 0.0346 19.5 0.0292 0.0292 22.5 0.0253 0.0253 26.0 0.0220 0.0220 30.0 0.0206 0.0206 34.5 0.0183 0.0183 40.0 0.0163 0.0163 46.0 0.0147 0.0147 52.5 0.0136 0.0136 60.0 0.0128 0.0128 69.0 0.0122 0.0122 79.0 0.0118 0.0118 90.5 0.0116 0.0116 105.5 0.0116 0.0116 123.5 0.0118 0.0118 143.0 0.0120 0.0120 163.5 0.0123 0.0123 185.0 0.0127 0.0127 208.0 0.0130 0.0130 232.5 0.0133 0.0133 258.5 0.0136 0.0136 286.0 0.0139 0.0139 331.0 0.0144 0.0144 396.0 0.0150 0.0150 468.5 0.0155 0.0155 549.5 0.0161 0.0161 639.0 0.0167 0.0167 738.0 0.0172 0.0172 847.5 0.0178 0.0178 968.5 0.0184 0.0184 1102.0 0.0191 0.0191 1249.5 0.0196 0.0196 1412.0 0.0202 0.0202 1590.5 0.0207 0.0207 1787.0 0.0213 0.0213 2003.0 0.0219 0.0219 2241.0 0.0224 0.0224 2503.0 0.0230 0.0230 2790.5 0.0235 0.0235 3107.0 0.0241 0.0241 3455.0 0.0246 0.0246 3837.0 0.0252 0.0252 4257.0 0.0258 0.0258 4719.0 0.0264 0.0264 5226.5 0.0269 0.0269 5784.0 0.0274 0.0274 6538.0 0.0472 0.0472 +1.0 1.2 150 9.0 0.0645 0.0645 11.0 0.0560 0.0560 13.5 0.0466 0.0466 16.5 0.0385 0.0385 19.5 0.0326 0.0326 22.5 0.0284 0.0284 26.0 0.0249 0.0249 30.0 0.0235 0.0235 34.5 0.0212 0.0212 40.0 0.0193 0.0193 46.0 0.0178 0.0178 52.5 0.0168 0.0168 60.0 0.0161 0.0161 69.0 0.0156 0.0156 79.0 0.0153 0.0153 90.5 0.0152 0.0152 105.5 0.0153 0.0153 123.5 0.0155 0.0155 143.0 0.0157 0.0157 163.5 0.0160 0.0160 185.0 0.0164 0.0164 208.0 0.0167 0.0167 232.5 0.0170 0.0170 258.5 0.0173 0.0173 286.0 0.0175 0.0175 331.0 0.0180 0.0180 396.0 0.0185 0.0185 468.5 0.0190 0.0190 549.5 0.0195 0.0195 639.0 0.0200 0.0200 738.0 0.0205 0.0205 847.5 0.0210 0.0210 968.5 0.0216 0.0216 1102.0 0.0221 0.0221 1249.5 0.0227 0.0227 1412.0 0.0235 0.0235 1590.5 0.0243 0.0243 1787.0 0.0251 0.0251 2003.0 0.0259 0.0259 2241.0 0.0267 0.0267 2503.0 0.0274 0.0274 2790.5 0.0282 0.0282 3107.0 0.0289 0.0289 3455.0 0.0297 0.0297 3837.0 0.0305 0.0305 4257.0 0.0310 0.0310 4719.0 0.0315 0.0315 5226.5 0.0321 0.0321 5784.0 0.0326 0.0326 6538.0 0.0383 0.0383 +1.2 1.4 150 9.0 0.0639 0.0639 11.0 0.0569 0.0569 13.5 0.0488 0.0488 16.5 0.0414 0.0414 19.5 0.0357 0.0357 22.5 0.0313 0.0313 26.0 0.0274 0.0274 30.0 0.0256 0.0256 34.5 0.0232 0.0232 40.0 0.0209 0.0209 46.0 0.0191 0.0191 52.5 0.0177 0.0177 60.0 0.0168 0.0168 69.0 0.0160 0.0160 79.0 0.0155 0.0155 90.5 0.0153 0.0153 105.5 0.0153 0.0153 123.5 0.0155 0.0155 143.0 0.0158 0.0158 163.5 0.0162 0.0162 185.0 0.0165 0.0165 208.0 0.0168 0.0168 232.5 0.0172 0.0172 258.5 0.0175 0.0175 286.0 0.0177 0.0177 331.0 0.0182 0.0182 396.0 0.0187 0.0187 468.5 0.0191 0.0191 549.5 0.0196 0.0196 639.0 0.0201 0.0201 738.0 0.0206 0.0206 847.5 0.0212 0.0212 968.5 0.0218 0.0218 1102.0 0.0224 0.0224 1249.5 0.0230 0.0230 1412.0 0.0239 0.0239 1590.5 0.0248 0.0248 1787.0 0.0256 0.0256 2003.0 0.0265 0.0265 2241.0 0.0274 0.0274 2503.0 0.0282 0.0282 2790.5 0.0291 0.0291 3107.0 0.0300 0.0300 3455.0 0.0307 0.0307 3837.0 0.0314 0.0314 4257.0 0.0320 0.0320 4719.0 0.0327 0.0327 5226.5 0.0334 0.0334 5784.0 0.0341 0.0341 6538.0 0.0453 0.0453 +1.4 1.6 150 9.0 0.0592 0.0592 11.0 0.0527 0.0527 13.5 0.0454 0.0454 16.5 0.0391 0.0391 19.5 0.0345 0.0345 22.5 0.0309 0.0309 26.0 0.0279 0.0279 30.0 0.0267 0.0267 34.5 0.0245 0.0245 40.0 0.0224 0.0224 46.0 0.0208 0.0208 52.5 0.0197 0.0197 60.0 0.0189 0.0189 69.0 0.0184 0.0184 79.0 0.0182 0.0182 90.5 0.0183 0.0183 105.5 0.0186 0.0186 123.5 0.0193 0.0193 143.0 0.0201 0.0201 163.5 0.0209 0.0209 185.0 0.0217 0.0217 208.0 0.0226 0.0226 232.5 0.0234 0.0234 258.5 0.0242 0.0242 286.0 0.0250 0.0250 331.0 0.0263 0.0263 396.0 0.0278 0.0278 468.5 0.0293 0.0293 549.5 0.0307 0.0307 639.0 0.0322 0.0322 738.0 0.0336 0.0336 847.5 0.0350 0.0350 968.5 0.0364 0.0364 1102.0 0.0374 0.0374 1249.5 0.0378 0.0378 1412.0 0.0381 0.0381 1590.5 0.0385 0.0385 1787.0 0.0388 0.0388 2003.0 0.0392 0.0392 2241.0 0.0395 0.0395 2503.0 0.0399 0.0399 2790.5 0.0403 0.0403 3107.0 0.0406 0.0406 3455.0 0.0410 0.0410 3837.0 0.0414 0.0414 4257.0 0.0419 0.0419 4719.0 0.0423 0.0423 5226.5 0.0428 0.0428 5784.0 0.0433 0.0433 6538.0 0.0439 0.0439 +1.6 1.8 150 9.0 0.0553 0.0553 11.0 0.0502 0.0502 13.5 0.0447 0.0447 16.5 0.0403 0.0403 19.5 0.0373 0.0373 22.5 0.0352 0.0352 26.0 0.0335 0.0335 30.0 0.0329 0.0329 34.5 0.0313 0.0313 40.0 0.0295 0.0295 46.0 0.0280 0.0280 52.5 0.0268 0.0268 60.0 0.0258 0.0258 69.0 0.0248 0.0248 79.0 0.0241 0.0241 90.5 0.0236 0.0236 105.5 0.0232 0.0232 123.5 0.0231 0.0231 143.0 0.0232 0.0232 163.5 0.0234 0.0234 185.0 0.0237 0.0237 208.0 0.0241 0.0241 232.5 0.0246 0.0246 258.5 0.0250 0.0250 286.0 0.0255 0.0255 331.0 0.0264 0.0264 396.0 0.0275 0.0275 468.5 0.0287 0.0287 549.5 0.0299 0.0299 639.0 0.0312 0.0312 738.0 0.0325 0.0325 847.5 0.0338 0.0338 968.5 0.0351 0.0351 1102.0 0.0359 0.0359 1249.5 0.0362 0.0362 1412.0 0.0365 0.0365 1590.5 0.0369 0.0369 1787.0 0.0373 0.0373 2003.0 0.0376 0.0376 2241.0 0.0380 0.0380 2503.0 0.0384 0.0384 2790.5 0.0388 0.0388 3107.0 0.0393 0.0393 3455.0 0.0398 0.0398 3837.0 0.0403 0.0403 4257.0 0.0409 0.0409 4719.0 0.0415 0.0415 5226.5 0.0421 0.0421 5784.0 0.0427 0.0427 6538.0 0.0434 0.0434 +1.8 2.0 150 9.0 0.0517 0.0517 11.0 0.0465 0.0465 13.5 0.0412 0.0412 16.5 0.0375 0.0375 19.5 0.0352 0.0352 22.5 0.0336 0.0336 26.0 0.0324 0.0324 30.0 0.0321 0.0321 34.5 0.0307 0.0307 40.0 0.0291 0.0291 46.0 0.0277 0.0277 52.5 0.0266 0.0266 60.0 0.0257 0.0257 69.0 0.0248 0.0248 79.0 0.0241 0.0241 90.5 0.0236 0.0236 105.5 0.0233 0.0233 123.5 0.0232 0.0232 143.0 0.0233 0.0233 163.5 0.0235 0.0235 185.0 0.0238 0.0238 208.0 0.0242 0.0242 232.5 0.0246 0.0246 258.5 0.0251 0.0251 286.0 0.0256 0.0256 331.0 0.0264 0.0264 396.0 0.0275 0.0275 468.5 0.0287 0.0287 549.5 0.0299 0.0299 639.0 0.0312 0.0312 738.0 0.0325 0.0325 847.5 0.0338 0.0338 968.5 0.0352 0.0352 1102.0 0.0360 0.0360 1249.5 0.0364 0.0364 1412.0 0.0368 0.0368 1590.5 0.0373 0.0373 1787.0 0.0377 0.0377 2003.0 0.0382 0.0382 2241.0 0.0387 0.0387 2503.0 0.0392 0.0392 2790.5 0.0398 0.0398 3107.0 0.0405 0.0405 3455.0 0.0411 0.0411 3837.0 0.0418 0.0418 4257.0 0.0426 0.0426 4719.0 0.0434 0.0434 5226.5 0.0442 0.0442 5784.0 0.0450 0.0450 6538.0 0.0460 0.0460 +2.0 2.2 150 9.0 0.0600 0.0600 11.0 0.0563 0.0563 13.5 0.0523 0.0523 16.5 0.0495 0.0495 19.5 0.0479 0.0479 22.5 0.0467 0.0467 26.0 0.0459 0.0459 30.0 0.0457 0.0457 34.5 0.0442 0.0442 40.0 0.0414 0.0414 46.0 0.0388 0.0388 52.5 0.0365 0.0365 60.0 0.0342 0.0342 69.0 0.0320 0.0320 79.0 0.0299 0.0299 90.5 0.0279 0.0279 105.5 0.0260 0.0260 123.5 0.0243 0.0243 143.0 0.0230 0.0230 163.5 0.0221 0.0221 185.0 0.0215 0.0215 208.0 0.0213 0.0213 232.5 0.0212 0.0212 258.5 0.0214 0.0214 286.0 0.0217 0.0217 331.0 0.0225 0.0225 396.0 0.0239 0.0239 468.5 0.0257 0.0257 549.5 0.0276 0.0276 639.0 0.0297 0.0297 738.0 0.0319 0.0319 847.5 0.0342 0.0342 968.5 0.0352 0.0352 1102.0 0.0357 0.0357 1249.5 0.0362 0.0362 1412.0 0.0368 0.0368 1590.5 0.0373 0.0373 1787.0 0.0379 0.0379 2003.0 0.0386 0.0386 2241.0 0.0393 0.0393 2503.0 0.0401 0.0401 2790.5 0.0410 0.0410 3107.0 0.0418 0.0418 3455.0 0.0427 0.0427 3837.0 0.0438 0.0438 4257.0 0.0449 0.0449 4719.0 0.0460 0.0460 5226.5 0.0471 0.0471 5784.0 0.0483 0.0483 6538.0 0.0498 0.0498 +2.2 2.4 150 9.0 0.0715 0.0715 11.0 0.0668 0.0668 13.5 0.0595 0.0595 16.5 0.0537 0.0537 19.5 0.0497 0.0497 22.5 0.0466 0.0466 26.0 0.0441 0.0441 30.0 0.0422 0.0422 34.5 0.0396 0.0396 40.0 0.0366 0.0366 46.0 0.0340 0.0340 52.5 0.0318 0.0318 60.0 0.0298 0.0298 69.0 0.0278 0.0278 79.0 0.0261 0.0261 90.5 0.0246 0.0246 105.5 0.0232 0.0232 123.5 0.0221 0.0221 143.0 0.0213 0.0213 163.5 0.0210 0.0210 185.0 0.0208 0.0208 208.0 0.0209 0.0209 232.5 0.0211 0.0211 258.5 0.0214 0.0214 286.0 0.0219 0.0219 331.0 0.0226 0.0226 396.0 0.0238 0.0238 468.5 0.0251 0.0251 549.5 0.0266 0.0266 639.0 0.0280 0.0280 738.0 0.0294 0.0294 847.5 0.0308 0.0308 968.5 0.0323 0.0323 1102.0 0.0338 0.0338 1249.5 0.0352 0.0352 1412.0 0.0359 0.0359 1590.5 0.0365 0.0365 1787.0 0.0371 0.0371 2003.0 0.0378 0.0378 2241.0 0.0386 0.0386 2503.0 0.0393 0.0393 2790.5 0.0402 0.0402 3107.0 0.0411 0.0411 3455.0 0.0420 0.0420 3837.0 0.0431 0.0431 4257.0 0.0442 0.0442 4719.0 0.0454 0.0454 5226.5 0.0466 0.0466 5784.0 0.0478 0.0478 6538.0 0.0494 0.0494 +2.4 2.6 150 9.0 0.1157 0.1157 11.0 0.1016 0.1016 13.5 0.0885 0.0885 16.5 0.0790 0.0790 19.5 0.0734 0.0734 22.5 0.0700 0.0700 26.0 0.0678 0.0678 30.0 0.0671 0.0671 34.5 0.0623 0.0623 40.0 0.0570 0.0570 46.0 0.0524 0.0524 52.5 0.0484 0.0484 60.0 0.0449 0.0449 69.0 0.0417 0.0417 79.0 0.0394 0.0394 90.5 0.0379 0.0379 105.5 0.0375 0.0375 123.5 0.0384 0.0384 143.0 0.0404 0.0404 163.5 0.0432 0.0432 185.0 0.0465 0.0465 208.0 0.0500 0.0500 232.5 0.0537 0.0537 258.5 0.0575 0.0575 286.0 0.0611 0.0611 331.0 0.0669 0.0669 396.0 0.0748 0.0748 468.5 0.0829 0.0829 549.5 0.0915 0.0915 639.0 0.0953 0.0953 738.0 0.0938 0.0938 847.5 0.0924 0.0924 968.5 0.0910 0.0910 1102.0 0.0902 0.0902 1249.5 0.0903 0.0903 1412.0 0.0904 0.0904 1590.5 0.0905 0.0905 1787.0 0.0906 0.0906 2003.0 0.0908 0.0908 2241.0 0.0909 0.0909 2503.0 0.0911 0.0911 2790.5 0.0912 0.0912 3107.0 0.0914 0.0914 3455.0 0.0916 0.0916 3837.0 0.0918 0.0918 4257.0 0.0920 0.0920 4719.0 0.0923 0.0923 5226.5 0.0925 0.0925 5784.0 0.0928 0.0928 6538.0 0.0898 0.0898 +2.6 2.8 150 9.0 0.1829 0.1829 11.0 0.1773 0.1773 13.5 0.1738 0.1738 16.5 0.1715 0.1715 19.5 0.1718 0.1718 22.5 0.1710 0.1710 26.0 0.1718 0.1718 30.0 0.1718 0.1718 34.5 0.1674 0.1674 40.0 0.1626 0.1626 46.0 0.1568 0.1568 52.5 0.1537 0.1537 60.0 0.1489 0.1489 69.0 0.1444 0.1444 79.0 0.1410 0.1410 90.5 0.1382 0.1382 105.5 0.1352 0.1352 123.5 0.1335 0.1335 143.0 0.1329 0.1329 163.5 0.1331 0.1331 185.0 0.1342 0.1342 208.0 0.1357 0.1357 232.5 0.1375 0.1375 258.5 0.1396 0.1396 286.0 0.1424 0.1424 331.0 0.1485 0.1485 396.0 0.1596 0.1596 468.5 0.1662 0.1662 549.5 0.1653 0.1653 639.0 0.1660 0.1660 738.0 0.1667 0.1667 847.5 0.1673 0.1673 968.5 0.1674 0.1674 1102.0 0.1675 0.1675 1249.5 0.1676 0.1676 1412.0 0.1676 0.1676 1590.5 0.1677 0.1677 1787.0 0.1678 0.1678 2003.0 0.1679 0.1679 2241.0 0.1680 0.1680 2503.0 0.1681 0.1681 2790.5 0.1682 0.1682 3107.0 0.1683 0.1683 3455.0 0.1684 0.1684 3837.0 0.1686 0.1686 4257.0 0.1688 0.1688 4719.0 0.1690 0.1690 5226.5 0.1625 0.1625 5784.0 0.1465 0.1465 6538.0 0.1586 0.1586 +2.8 3.0 150 9.0 0.2471 0.2471 11.0 0.2421 0.2421 13.5 0.2393 0.2393 16.5 0.2380 0.2380 19.5 0.2369 0.2369 22.5 0.2376 0.2376 26.0 0.2385 0.2385 30.0 0.2361 0.2361 34.5 0.2300 0.2300 40.0 0.2246 0.2246 46.0 0.2197 0.2197 52.5 0.2155 0.2155 60.0 0.2122 0.2122 69.0 0.2096 0.2096 79.0 0.2079 0.2079 90.5 0.2072 0.2072 105.5 0.2076 0.2076 123.5 0.2093 0.2093 143.0 0.2127 0.2127 163.5 0.2172 0.2172 185.0 0.2229 0.2229 208.0 0.2295 0.2295 232.5 0.2371 0.2371 258.5 0.2457 0.2457 286.0 0.2550 0.2550 331.0 0.2707 0.2707 396.0 0.2724 0.2724 468.5 0.2646 0.2646 549.5 0.2656 0.2656 639.0 0.2661 0.2661 738.0 0.2665 0.2665 847.5 0.2665 0.2665 968.5 0.2666 0.2666 1102.0 0.2666 0.2666 1249.5 0.2667 0.2667 1412.0 0.2667 0.2667 1590.5 0.2668 0.2668 1787.0 0.2668 0.2668 2003.0 0.2669 0.2669 2241.0 0.2669 0.2669 2503.0 0.2670 0.2670 2790.5 0.2670 0.2670 3107.0 0.2671 0.2671 3455.0 0.2672 0.2672 3837.0 0.2673 0.2673 4257.0 0.2674 0.2674 4719.0 0.2286 0.2286 5226.5 0.2062 0.2062 5784.0 0.2035 0.2035 6538.0 0.2037 0.2037 +3.0 3.5 150 9.0 0.1460 0.1460 11.0 0.1299 0.1299 13.5 0.1164 0.1164 16.5 0.1066 0.1066 19.5 0.1005 0.1005 22.5 0.0969 0.0969 26.0 0.0944 0.0944 30.0 0.0937 0.0937 34.5 0.0919 0.0919 40.0 0.0903 0.0903 46.0 0.0886 0.0886 52.5 0.0876 0.0876 60.0 0.0868 0.0868 69.0 0.0862 0.0862 79.0 0.0857 0.0857 90.5 0.0854 0.0854 105.5 0.0853 0.0853 123.5 0.0855 0.0855 143.0 0.0857 0.0857 163.5 0.0861 0.0861 185.0 0.0866 0.0866 208.0 0.0871 0.0871 232.5 0.0878 0.0878 258.5 0.0885 0.0885 286.0 0.0892 0.0892 331.0 0.0904 0.0904 396.0 0.0921 0.0921 468.5 0.0936 0.0936 549.5 0.0942 0.0942 639.0 0.0942 0.0942 738.0 0.0943 0.0943 847.5 0.0944 0.0944 968.5 0.0945 0.0945 1102.0 0.0946 0.0946 1249.5 0.0948 0.0948 1412.0 0.0949 0.0949 1590.5 0.0950 0.0950 1787.0 0.0951 0.0951 2003.0 0.0952 0.0952 2241.0 0.0953 0.0953 2503.0 0.0954 0.0954 2790.5 0.0955 0.0955 3107.0 0.0956 0.0956 3455.0 0.0957 0.0957 3837.0 0.0958 0.0958 4257.0 0.0959 0.0959 4719.0 0.0960 0.0960 5226.5 0.0961 0.0961 5784.0 0.0963 0.0963 6538.0 0.0964 0.0964 +3.5 4.0 150 9.0 0.1467 0.1467 11.0 0.1261 0.1261 13.5 0.1093 0.1093 16.5 0.0969 0.0969 19.5 0.0887 0.0887 22.5 0.0835 0.0835 26.0 0.0797 0.0797 30.0 0.0783 0.0783 34.5 0.0757 0.0757 40.0 0.0734 0.0734 46.0 0.0711 0.0711 52.5 0.0694 0.0694 60.0 0.0685 0.0685 69.0 0.0672 0.0672 79.0 0.0664 0.0664 90.5 0.0661 0.0661 105.5 0.0662 0.0662 123.5 0.0667 0.0667 143.0 0.0676 0.0676 163.5 0.0689 0.0689 185.0 0.0703 0.0703 208.0 0.0718 0.0718 232.5 0.0735 0.0735 258.5 0.0751 0.0751 286.0 0.0769 0.0769 331.0 0.0786 0.0786 396.0 0.0791 0.0791 468.5 0.0792 0.0792 549.5 0.0793 0.0793 639.0 0.0794 0.0794 738.0 0.0795 0.0795 847.5 0.0796 0.0796 968.5 0.0797 0.0797 1102.0 0.0799 0.0799 1249.5 0.0800 0.0800 1412.0 0.0801 0.0801 1590.5 0.0803 0.0803 1787.0 0.0804 0.0804 2003.0 0.0805 0.0805 2241.0 0.0807 0.0807 2503.0 0.0808 0.0808 2790.5 0.0810 0.0810 3107.0 0.0811 0.0811 3455.0 0.0813 0.0813 3837.0 0.0814 0.0814 4257.0 0.0815 0.0815 4719.0 0.0816 0.0816 5226.5 0.0817 0.0817 5784.0 0.0819 0.0819 6538.0 0.0820 0.0820 +4.0 4.4 150 9.0 0.1307 0.1307 11.0 0.1110 0.1110 13.5 0.0943 0.0943 16.5 0.0817 0.0817 19.5 0.0738 0.0738 22.5 0.0685 0.0685 26.0 0.0647 0.0647 30.0 0.0624 0.0624 34.5 0.0605 0.0605 40.0 0.0578 0.0578 46.0 0.0556 0.0556 52.5 0.0538 0.0538 60.0 0.0526 0.0526 69.0 0.0518 0.0518 79.0 0.0513 0.0513 90.5 0.0511 0.0511 105.5 0.0514 0.0514 123.5 0.0520 0.0520 143.0 0.0528 0.0528 163.5 0.0539 0.0539 185.0 0.0548 0.0548 208.0 0.0556 0.0556 232.5 0.0562 0.0562 258.5 0.0561 0.0561 286.0 0.0562 0.0562 331.0 0.0564 0.0564 396.0 0.0566 0.0566 468.5 0.0570 0.0570 549.5 0.0573 0.0573 639.0 0.0577 0.0577 738.0 0.0582 0.0582 847.5 0.0587 0.0587 968.5 0.0592 0.0592 1102.0 0.0597 0.0597 1249.5 0.0602 0.0602 1412.0 0.0608 0.0608 1590.5 0.0614 0.0614 1787.0 0.0620 0.0620 2003.0 0.0626 0.0626 2241.0 0.0632 0.0632 2503.0 0.0638 0.0638 2790.5 0.0645 0.0645 3107.0 0.0652 0.0652 3455.0 0.0658 0.0658 3837.0 0.0664 0.0664 4257.0 0.0670 0.0670 4719.0 0.0676 0.0676 5226.5 0.0682 0.0682 5784.0 0.0689 0.0689 6538.0 0.0697 0.0697 +4.4 5.0 150 9.0 0.0923 0.0923 11.0 0.0824 0.0824 13.5 0.0745 0.0745 16.5 0.0693 0.0693 19.5 0.0659 0.0659 22.5 0.0639 0.0639 26.0 0.0625 0.0625 30.0 0.0608 0.0608 34.5 0.0593 0.0593 40.0 0.0566 0.0566 46.0 0.0545 0.0545 52.5 0.0527 0.0527 60.0 0.0515 0.0515 69.0 0.0508 0.0508 79.0 0.0507 0.0507 90.5 0.0511 0.0511 105.5 0.0521 0.0521 123.5 0.0532 0.0532 143.0 0.0541 0.0541 163.5 0.0552 0.0552 185.0 0.0562 0.0562 208.0 0.0576 0.0576 232.5 0.0590 0.0590 258.5 0.0597 0.0597 286.0 0.0606 0.0606 331.0 0.0619 0.0619 396.0 0.0634 0.0634 468.5 0.0649 0.0649 549.5 0.0664 0.0664 639.0 0.0677 0.0677 738.0 0.0690 0.0690 847.5 0.0703 0.0703 968.5 0.0715 0.0715 1102.0 0.0726 0.0726 1249.5 0.0737 0.0737 1412.0 0.0748 0.0748 1590.5 0.0759 0.0759 1787.0 0.0769 0.0769 2003.0 0.0779 0.0779 2241.0 0.0789 0.0789 2503.0 0.0799 0.0799 2790.5 0.0808 0.0808 3107.0 0.0818 0.0818 3455.0 0.0827 0.0827 3837.0 0.0838 0.0838 4257.0 0.0849 0.0849 4719.0 0.0859 0.0859 5226.5 0.0870 0.0870 5784.0 0.0880 0.0880 6538.0 0.0892 0.0892 +5.0 5.4 150 9.0 0.1012 0.1012 11.0 0.0889 0.0889 13.5 0.0790 0.0790 16.5 0.0722 0.0722 19.5 0.0682 0.0682 22.5 0.0656 0.0656 26.0 0.0637 0.0637 30.0 0.0617 0.0617 34.5 0.0598 0.0598 40.0 0.0570 0.0570 46.0 0.0552 0.0552 52.5 0.0528 0.0528 60.0 0.0518 0.0518 69.0 0.0508 0.0508 79.0 0.0501 0.0501 90.5 0.0495 0.0495 105.5 0.0494 0.0494 123.5 0.0499 0.0499 143.0 0.0508 0.0508 163.5 0.0521 0.0521 185.0 0.0531 0.0531 208.0 0.0546 0.0546 232.5 0.0560 0.0560 258.5 0.0568 0.0568 286.0 0.0577 0.0577 331.0 0.0590 0.0590 396.0 0.0607 0.0607 468.5 0.0623 0.0623 549.5 0.0638 0.0638 639.0 0.0652 0.0652 738.0 0.0665 0.0665 847.5 0.0678 0.0678 968.5 0.0690 0.0690 1102.0 0.0702 0.0702 1249.5 0.0714 0.0714 1412.0 0.0725 0.0725 1590.5 0.0736 0.0736 1787.0 0.0747 0.0747 2003.0 0.0757 0.0757 2241.0 0.0767 0.0767 2503.0 0.0777 0.0777 2790.5 0.0787 0.0787 3107.0 0.0797 0.0797 3455.0 0.0806 0.0806 3837.0 0.0817 0.0817 4257.0 0.0829 0.0829 4719.0 0.0839 0.0839 5226.5 0.0850 0.0850 5784.0 0.0860 0.0860 6538.0 0.0873 0.0873 From fb0351579ac3cc44224a88368fcf8cd8b1fce0ce Mon Sep 17 00:00:00 2001 From: capalmer Date: Mon, 15 Apr 2019 23:18:51 +0200 Subject: [PATCH 5/7] add SF measurements in pt bins --- test/ttbar/TTbarEventAnalysis.cc | 212 +++++++++++++++++-------------- test/ttbar/TTbarEventAnalysis.h | 8 +- test/ttbar/runTTbarAnalysis.py | 34 +++-- test/ttbar/twoTag.py | 187 +++++++++++++++++---------- 4 files changed, 261 insertions(+), 180 deletions(-) diff --git a/test/ttbar/TTbarEventAnalysis.cc b/test/ttbar/TTbarEventAnalysis.cc index 926e3caeca8..cf8f69ec629 100644 --- a/test/ttbar/TTbarEventAnalysis.cc +++ b/test/ttbar/TTbarEventAnalysis.cc @@ -179,31 +179,36 @@ void TTbarEventAnalysis::prepareOutput(TString outFile) baseHistos["svhe"]=new TH1F("svhe",";Simple secondary vertex (HE);Jets",50,0,6); baseHistos["csv"]=new TH1F("csv",";Combined secondary vertex (IVF);Jets",52,-0.02,1.02); baseHistos["deepcsv"]=new TH1F("deepcsv",";DeepCSV for B;Jets",52,-0.02,1.02); - - std::vector twoTagNames; + + std::vector myBins={20,30,50,70,100,200,300}; + lowerPtBinEdges.insert(lowerPtBinEdges.begin(),myBins.begin(),myBins.end()); + twoTagNames.clear(); - twoTagNames.push_back("twoTags_deepCSV"); - twoTagNames.push_back("only2_twoTags_deepCSV"); - twoTagNames.push_back("twoTags_deepFlavour"); - twoTagNames.push_back("only2_twoTags_deepFlavour"); - twoTagNames.push_back("twoTags_CSVv2"); - twoTagNames.push_back("only2_twoTags_CSVv2"); - //twoTagNames.push_back(""); + if(runTwoTagAnalysis>0){ + twoTagNames.push_back("twoTags_deepCSV"); + twoTagNames.push_back("only2_twoTags_deepCSV"); + twoTagNames.push_back("twoTags_deepFlavour"); + twoTagNames.push_back("only2_twoTags_deepFlavour"); + } for(unsigned int i2Tag=0; i2Tag2.5) continue; + if(varjp4[iSystVar].Pt()<20 || TMath::Abs(varjp4[iSystVar].Eta())>2.5) continue; canBeSelected=true; jetCount[iSystVar]++; } @@ -691,10 +696,11 @@ void TTbarEventAnalysis::processFile(TString inFile,TH1F *xsecWgt,Bool_t isData) std::pair bestDeepFlavourPair(-1,-1); std::pair bestDeepCSVPair(-1,-1); - std::pair bestCSVv2Pair(-1,-1); GetBestJetPair(bestDeepFlavourPair,"deepFlavour"); GetBestJetPair(bestDeepCSVPair,"deepCSV"); - GetBestJetPair(bestCSVv2Pair,"CSVv2"); + bestJetPairs.clear(); + bestJetPairs["deepFlavour"]=bestDeepFlavourPair; + bestJetPairs["deepCSV"]=bestDeepCSVPair; n2++; //bool twoDeepCSVJets=(ev.Jet_DeepCSVBDisc[selJets[bestDeepCSVPair.first]]>0 && ev.Jet_DeepCSVBDisc[selJets[bestDeepCSVPair.second]]>0); @@ -719,6 +725,8 @@ void TTbarEventAnalysis::processFile(TString inFile,TH1F *xsecWgt,Bool_t isData) histos_[ch+"_trailjeta"]->Fill(fabs(selJetsP4[0][1].Eta()),evWgt); histos_[ch+"_trailbjpt"]->Fill(selJetsP4[0][selJets[bestDeepCSVPair.second]].Pt(),evWgt); histos_[ch+"_trailbjeta"]->Fill(fabs(selJetsP4[0][selJets[bestDeepCSVPair.second]].Eta()),evWgt); + histos2d_[ch+"_JETPT1VS2"]->Fill(selJetsP4[0][selJets[bestDeepCSVPair.first]].Pt(),selJetsP4[0][selJets[bestDeepCSVPair.second]].Pt(),evWgt); + if(selJets.size()==2){ histos_[ch+"_only2_leadjpt"]->Fill(selJetsP4[0][0].Pt(),evWgt); histos_[ch+"_only2_leadjeta"]->Fill((selJetsP4[0][0].Eta()),evWgt); @@ -737,99 +745,71 @@ void TTbarEventAnalysis::processFile(TString inFile,TH1F *xsecWgt,Bool_t isData) // 2018 WPs // https://twiki.cern.ch/twiki/bin/viewauth/CMS/BtagRecommendation102X#AK4_jets - - std::vector deepCSVWPs; - deepCSVWPs.clear(); - //deepCSVWPs.push_back(0.1522); //2017 - //deepCSVWPs.push_back(0.4941); - //deepCSVWPs.push_back(0.8001); - deepCSVWPs.push_back(0.1241); //2018 - deepCSVWPs.push_back(0.4184); - deepCSVWPs.push_back(0.7527); - btaggingWPs["deepCSV"]=deepCSVWPs; + + if(runTwoTagAnalysis>0){ + std::vector deepCSVWPs; + deepCSVWPs.clear(); + //deepCSVWPs.push_back(0.1522); //2017 + //deepCSVWPs.push_back(0.4941); + //deepCSVWPs.push_back(0.8001); + deepCSVWPs.push_back(0.1241); //2018 + deepCSVWPs.push_back(0.4184); + deepCSVWPs.push_back(0.7527); + btaggingWPs["deepCSV"]=deepCSVWPs; - std::vector CSVv2WPs; - CSVv2WPs.clear(); - CSVv2WPs.push_back(0.5803); //2017 - CSVv2WPs.push_back(0.8838); - CSVv2WPs.push_back(0.9693); - btaggingWPs["CSVv2"]=CSVv2WPs; + //std::vector CSVv2WPs; + //CSVv2WPs.clear(); + //CSVv2WPs.push_back(0.5803); //2017 + //CSVv2WPs.push_back(0.8838); + //CSVv2WPs.push_back(0.9693); + //btaggingWPs["CSVv2"]=CSVv2WPs; - std::vector deepFlavourWPs; - deepFlavourWPs.clear(); - deepFlavourWPs.push_back(0.0494); //2018 - deepFlavourWPs.push_back(0.2770); - deepFlavourWPs.push_back(0.7264); - btaggingWPs["deepFlavour"]=deepFlavourWPs; - + std::vector deepFlavourWPs; + deepFlavourWPs.clear(); + deepFlavourWPs.push_back(0.0494); //2018 + deepFlavourWPs.push_back(0.2770); + deepFlavourWPs.push_back(0.7264); + btaggingWPs["deepFlavour"]=deepFlavourWPs; + } //histos_[ch+"_deepcsvlead"]->Fill(ev.Jet_DeepCSVBDisc[selJets[bestDeepCSVPair.first]],evWgt); //histos_[ch+"_deepcsvsublead"]->Fill(ev.Jet_DeepCSVBDisc[selJets[bestDeepCSVPair.second]],evWgt); - std::vector nPassCSVv2WPs(3,0); - for(unsigned int iWP=0; iWPCSVv2WPs[iWP]); - nPassCSVv2WPs[iWP]+= (ev.Jet_CombIVF[selJets[bestCSVv2Pair.second]]>CSVv2WPs[iWP]); - } - - std::vector nPassDeepCSVWPs(3,0); - for(unsigned int iWP=0; iWPdeepCSVWPs[iWP]); - nPassDeepCSVWPs[iWP]+= (ev.Jet_DeepCSVBDisc[selJets[bestDeepCSVPair.second]]>deepCSVWPs[iWP]); - //std::cout<<"1 2 WP nPass "<>::iterator iMap=btaggingWPs.begin(); iMap!=btaggingWPs.end(); iMap++){ + for(int iPT=-1; iPT<(int)lowerPtBinEdges.size(); iPT++){ + histos_[ch+"_pt2OverPt1_"+iMap->first]->Fill(ReturnVarAtIndex("PT",bestJetPairs[iMap->first].second)/ReturnVarAtIndex("PT",bestJetPairs[iMap->first].first),evWgt); + TwoTag(ch+"_twoTags_"+iMap->first,iMap->first,bestJetPairs[iMap->first],iPT); + if(selJets.size()==2){ + TwoTag(ch+"_only2_twoTags_"+iMap->first,iMap->first,bestJetPairs[iMap->first],iPT); + } + } } - + histos_[ch+"_deepcsv1"]->Fill(ev.Jet_DeepCSVBDisc[selJets[bestDeepCSVPair.first]],evWgt); histos_[ch+"_deepcsv2"]->Fill(ev.Jet_DeepCSVBDisc[selJets[bestDeepCSVPair.second]],evWgt); histos2d_[ch+"_deepcsv2d"]->Fill(ev.Jet_DeepCSVBDisc[selJets[bestDeepCSVPair.first]],ev.Jet_DeepCSVBDisc[selJets[bestDeepCSVPair.second]],evWgt); //std::cout<<"go go 012"<Fill(ev.Jet_DeepCSVBDisc[selJets[bestDeepCSVPair.first]],evWgt); histos_[ch+"_only2_deepcsv2"]->Fill(ev.Jet_DeepCSVBDisc[selJets[bestDeepCSVPair.second]],evWgt); histos2d_[ch+"_only2_deepcsv2d"]->Fill(ev.Jet_DeepCSVBDisc[selJets[bestDeepCSVPair.first]],ev.Jet_DeepCSVBDisc[selJets[bestDeepCSVPair.second]],evWgt); - TwoTag(ch+"_only2_twoTags_deepCSV", "deepCSV", bestDeepCSVPair); } //std::cout<<"go go 014"<Fill(ev.Jet_DeepFlavourBDisc[selJets[bestDeepFlavourPair.first]],evWgt); //std::cout<<"go go 014.1"<Fill(ev.Jet_DeepFlavourBDisc[selJets[bestDeepFlavourPair.second]],evWgt); - //std::cout<<"go go 014.2"<Fill(ev.Jet_DeepFlavourBDisc[selJets[bestDeepFlavourPair.first]],evWgt); - //std::cout<<"go go 014.4"<Fill(ev.Jet_DeepFlavourBDisc[selJets[bestDeepFlavourPair.second]],evWgt); - //std::cout<<"go go 014.5"<Fill(ev.Jet_CombIVF[selJets[bestCSVv2Pair.first]],evWgt); - histos_[ch+"_csvv22"]->Fill(ev.Jet_CombIVF[selJets[bestCSVv2Pair.second]],evWgt); - histos2d_[ch+"_csvv22d"]->Fill(ev.Jet_CombIVF[selJets[bestCSVv2Pair.first]],ev.Jet_CombIVF[selJets[bestCSVv2Pair.second]],evWgt); - - //if(binProduct==3){ - // histos2d_[ch+"_csvv22d_2b"]->Fill(ev.Jet_CombIVF[selJets[bestCSVv2Pair.first]],ev.Jet_CombIVF[selJets[bestCSVv2Pair.second]],evWgt); - //} - - TwoTag(ch+"_twoTags_CSVv2", "CSVv2", bestCSVv2Pair); - //std::cout<<"go go 014.7"<Fill(ev.Jet_CombIVF[selJets[bestCSVv2Pair.first]],evWgt); - histos_[ch+"_only2_csvv22"]->Fill(ev.Jet_CombIVF[selJets[bestCSVv2Pair.second]],evWgt); - histos2d_[ch+"_only2_csvv22d"]->Fill(ev.Jet_CombIVF[selJets[bestCSVv2Pair.first]],ev.Jet_CombIVF[selJets[bestCSVv2Pair.second]],evWgt); - //if(binProduct==3){ - // histos2d_[ch+"_only2_csvv22d_2b"]->Fill(ev.Jet_CombIVF[selJets[bestCSVv2Pair.first]],ev.Jet_CombIVF[selJets[bestCSVv2Pair.second]],evWgt); - //} - TwoTag(ch+"_only2_twoTags_CSVv2", "CSVv2", bestCSVv2Pair); - } std::vector leadingkindisc(2,-9999); std::vector leadingkindiscIdx(2,-1); @@ -1095,10 +1075,50 @@ float TTbarEventAnalysis::ReturnVarAtIndex(std::string varName, unsigned int ind return 0; } } + + +std::string TTbarEventAnalysis::ReturnPtLabel(int iPT){ + std::string ptBinLabel; + ptBinLabel.clear(); + if(iPT==-1){ + ptBinLabel+="Inclusive"; + } else { + ptBinLabel+=std::to_string(lowerPtBinEdges[iPT])+"to"; + if(iPT==(int)(lowerPtBinEdges.size()-1)){ + ptBinLabel+="Inf"; + }else{ + ptBinLabel+=std::to_string(lowerPtBinEdges[iPT+1]); + } + } + return ptBinLabel; +} // two tag computations and histo filling -void TTbarEventAnalysis::TwoTag(std::string tagName, std::string discriminator, std::pair jetIndices) +void TTbarEventAnalysis::TwoTag(std::string tagName, std::string discriminator, std::pair jetIndices, int ptBin) { + // return before filling if not in this pt bin + if(ptBin>-1){ //-1 is inclusive--> always keep -1 + float pt1 = -2; + float pt2 = -2; + + pt1 = ReturnVarAtIndex("PT",jetIndices.first); + pt2 = ReturnVarAtIndex("PT",jetIndices.second); + + //std::cout<<"pt1 pt2 "<std::max(pt1,pt2)){ + //std::cout<<"reject because both smaller than "< nPassWPs(btaggingWPs[discriminator].size(),0); float btag1 = -2; float btag2 = -2; @@ -1140,9 +1160,9 @@ void TTbarEventAnalysis::TwoTag(std::string tagName, std::string discriminator, for(unsigned int iWP=0; iWPFill(twoTagCrossFlavour[iWP],evWgt); + histos_[tagName+wpLabel[iWP]+ReturnPtLabel(ptBin)]->Fill(twoTagCrossFlavour[iWP],evWgt); for(unsigned int iSyst=0; iSystFill(twoTagCrossFlavour[iWP],systWeight[systName[iSyst]]); + histos_[tagName+wpLabel[iWP]+ReturnPtLabel(ptBin)+"_"+systName[iSyst]]->Fill(twoTagCrossFlavour[iWP],systWeight[systName[iSyst]]); } } } diff --git a/test/ttbar/TTbarEventAnalysis.h b/test/ttbar/TTbarEventAnalysis.h index d07e103dd6e..84dbf7d506d 100644 --- a/test/ttbar/TTbarEventAnalysis.h +++ b/test/ttbar/TTbarEventAnalysis.h @@ -71,6 +71,7 @@ class TTbarEventAnalysis } ~TTbarEventAnalysis(){} void setReadTTJetsGenWeights(bool readTTJetsGenWeights) { readTTJetsGenWeights_=readTTJetsGenWeights; } + void setTwoTagCount(bool runTwoTagCount_) { runTwoTagAnalysis=runTwoTagCount_; } void setTMVAWeightsBaseDir(TString url) { weightsDir_=url; gSystem->ExpandPathName(weightsDir_); } void addTriggerBit(Int_t bit,Int_t ch) { triggerBits_.push_back(std::pair(bit,ch)); } void addVarForTMVA(TString varName) { tmvaVarNames_.push_back(varName); } @@ -118,10 +119,14 @@ class TTbarEventAnalysis std::vector systName; std::map systWeight; - void TwoTag(std::string tagName, std::string discriminator, std::pair); + void TwoTag(std::string tagName, std::string discriminator, std::pair, int ptBin); std::map> btaggingWPs; void GetBestJetPair(std::pair& myIndices, std::string discriminator="deepCSV"); + std::map> bestJetPairs; float ReturnVarAtIndex(std::string varName, unsigned int index); + std::string ReturnPtLabel(int iPT); + std::vector lowerPtBinEdges; + std::vector twoTagNames; TGraph *puWgtGr_,*puWgtDownGr_,*puWgtUpGr_; bool readTTJetsGenWeights_; @@ -147,6 +152,7 @@ class TTbarEventAnalysis std::map histos_; std::map histos2d_; bool noEventsSelected; + unsigned int runTwoTagAnalysis; }; diff --git a/test/ttbar/runTTbarAnalysis.py b/test/ttbar/runTTbarAnalysis.py index a37648f4b84..4b587b6da6a 100644 --- a/test/ttbar/runTTbarAnalysis.py +++ b/test/ttbar/runTTbarAnalysis.py @@ -11,6 +11,22 @@ CHANNELS={-11*11:'ll', -13*13:'ll', -11*13:'emu'} +#configuration +usage = 'usage: %prog [options]' +parser = optparse.OptionParser(usage) +parser.add_option('-j', '--json', dest='json' , help='json with list of files', default=None, type='string') +parser.add_option('-i', '--inDir', dest='inDir', help='input directory with files', default=None, type='string') +parser.add_option('-o', '--outDir', dest='outDir', help='output directory', default='analysis', type='string') +parser.add_option( '--only', dest='only', help='process only matching (csv)', default='', type='string') +parser.add_option( '--tmvaWgts', dest='tmvaWgts', help='tmva weights', default=None, type='string') +parser.add_option( '--dyScale', dest='dyScale', help='DY scale factor', default=None, type='string') +parser.add_option('-n', '--njobs', dest='njobs', help='# jobs to run in parallel', default=0, type='int') +parser.add_option('--doTwoTag', dest='doTwoTag', help='Run 2TagCount (default: False)', default=0, type='int') +parser.add_option('--keepGoodFilesBySize', dest='keepGoodFilesBySize', help='Review output file sizes before running and skip good ones', default=1, type='int') +parser.add_option('--inspectOutput', dest='inspectOutput', help='Review output file sizes and re-run baddies', default=1, type='int') +(opt, args) = parser.parse_args() + + """ Perform the analysis on a single file """ @@ -19,6 +35,10 @@ def runTTbarAnalysis(inFile, outFile, wgt, tmvaWgts=None,isData=False): from ROOT import TTbarEventAnalysis evAnalysis=TTbarEventAnalysis() + + print "opt.doTwoTag",opt.doTwoTag + evAnalysis.setTwoTagCount(opt.doTwoTag) + #MC specifics if 'TTJets' in inFile: evAnalysis.setReadTTJetsGenWeights(True) @@ -112,20 +132,6 @@ def InspectFiles(outputDir,runTags): """ def main(): - #configuration - usage = 'usage: %prog [options]' - parser = optparse.OptionParser(usage) - parser.add_option('-j', '--json', dest='json' , help='json with list of files', default=None, type='string') - parser.add_option('-i', '--inDir', dest='inDir', help='input directory with files', default=None, type='string') - parser.add_option('-o', '--outDir', dest='outDir', help='output directory', default='analysis', type='string') - parser.add_option( '--only', dest='only', help='process only matching (csv)', default='', type='string') - parser.add_option( '--tmvaWgts', dest='tmvaWgts', help='tmva weights', default=None, type='string') - parser.add_option( '--dyScale', dest='dyScale', help='DY scale factor', default=None, type='string') - parser.add_option('-n', '--njobs', dest='njobs', help='# jobs to run in parallel', default=0, type='int') - parser.add_option('--keepGoodFilesBySize', dest='keepGoodFilesBySize', help='Review output file sizes before running and skip good ones', default=1, type='int') - parser.add_option('--inspectOutput', dest='inspectOutput', help='Review output file sizes and re-run baddies', default=1, type='int') - (opt, args) = parser.parse_args() - #compile c++ wrapper to run over trees ROOT.gSystem.Load("libJetMETCorrectionsObjects.so") ROOT.gSystem.CompileMacro("TTbarEventAnalysis.cc","fkgd","libTTbarEventAnalysis"); diff --git a/test/ttbar/twoTag.py b/test/ttbar/twoTag.py index 8481444551e..0cd5d4632fa 100644 --- a/test/ttbar/twoTag.py +++ b/test/ttbar/twoTag.py @@ -5,13 +5,29 @@ fileName=sys.argv[1] tFile=ROOT.TFile.Open(fileName) +lowerPtBinEdges=[20,30,50,70,100,200,300] +def ReturnPtLabel(iPT): + ptBinLabel="" + if(iPT==-1): + ptBinLabel="Inclusive" + else: + ptBinLabel+=str(lowerPtBinEdges[iPT])+"to" + if(iPT==len(lowerPtBinEdges)-1): + ptBinLabel+="Inf" + else: + ptBinLabel+=str(lowerPtBinEdges[iPT+1]) + return ptBinLabel + + +tags=["twoTags_deepCSVL"] tags=["twoTags_deepCSVL","only2_twoTags_deepCSVL","twoTags_deepCSVM","only2_twoTags_deepCSVM","twoTags_deepCSVT","only2_twoTags_deepCSVT","twoTags_deepFlavourL","only2_twoTags_deepFlavourL","twoTags_deepFlavourM","only2_twoTags_deepFlavourM","twoTags_deepFlavourT","only2_twoTags_deepFlavourT","twoTags_CSVv2L","only2_twoTags_CSVv2L","twoTags_CSVv2M","only2_twoTags_CSVv2M","twoTags_CSVv2T","only2_twoTags_CSVv2T"] +tags=["twoTags_deepCSVL","twoTags_deepCSVM","twoTags_deepCSVT","twoTags_deepFlavourL","only2_twoTags_deepFlavourL","twoTags_deepFlavourM","only2_twoTags_deepFlavourM","twoTags_deepFlavourT","only2_twoTags_deepFlavourT"] tags.sort() samples=["t#bar{t} SL","VV","tW","t#bar{t} DL"] samples=["t#bar{t} SL","t#bar{t} DL","t#bar{t} AH","ZZ","WZ","WW","tW","t#bar{t} DL","DY"] -systs=[""] systs=["","_fsrConHi","_fsrConLo","_fsrDefHi","_fsrDefLo","_fsrRedHi","_fsrRedLo","_isrConHi","_isrConLo","_isrDefHi","_isrDefLo","_isrRedHi","_isrRedLo","_qcdScaleHi","_qcdScaleLo"] +systs=[""] # KEY: TH1F emu_only2_twoTags_CSVv2L_t#bar{t} AH;1 t#bar{t} AH # KEY: TH1F emu_only2_twoTags_CSVv2L_other;1 other @@ -21,85 +37,114 @@ statError={} -#samples=["emu_only2_twoTags/emu_only2_twoTags_t#bar{t} SL","emu_only2_twoTags/emu_only2_twoTags_VV","emu_only2_twoTags/emu_only2_twoTags_tW","emu_only2_twoTags/emu_only2_twoTags_t#bar{t} DL","emu_only2_twoTags/emu_only2_twoTags_DY"] -#twoTagDataHist="emu_only2_twoTags/emu_only2_twoTags" #can=ROOT.TCanvas("can","",700,700) results={} +tgraphs={} + +can=ROOT.TCanvas("can","",800,800) for syst in systs: #print syst for tag in tags: print tag, + iGraph=0 + if syst=="": + tgraphs[tag]=ROOT.TGraphErrors() + tgraphs[tag].SetTitle(tag+";Jet P_{T};Scale Factor") if tag.find("twoTags_") == 0: print " ", - for sample in samples: - histName="emu_"+tag+syst+"/emu_"+tag+syst+"_"+sample - #print histName - hist=tFile.Get(histName) - if sample==samples[0]: - histBkg=hist.Clone("histBkg") - else: - histBkg.Add(hist) - #print histName - #for iBin in range(15): - # print iBin,hist.GetBinContent(iBin) - - tot={} - tot["2bmc"]=0 - tot["othermc"]=0 - tot["data"]=0 - tot["2bmc_2btag"]=0 - tot["othermc_2btag"]=0 - tot["data_2btag"]=0 - for iBin in range(15): - #print iBin,histBkg.GetBinContent(iBin) - if iBin%4==1: - tot["2bmc"]=tot["2bmc"]+histBkg.GetBinContent(iBin) - if iBin>9: - tot["2bmc_2btag"]=tot["2bmc_2btag"]+histBkg.GetBinContent(iBin) - else: - tot["othermc"]=tot["othermc"]+histBkg.GetBinContent(iBin) + for ptBin in range(-1,len(lowerPtBinEdges)): + ptLabel=ReturnPtLabel(ptBin) + basename="emu_"+tag+syst+ptLabel + print basename + for sample in samples: + histName=basename+"/"+basename+"_"+sample + print histName + hist=tFile.Get(histName) + if sample==samples[0]: + histBkg=hist.Clone("histBkg") + else: + histBkg.Add(hist) + #print histName + #for iBin in range(15): + # print iBin,hist.GetBinContent(iBin) + + tot={} + tot["2bmc"]=0 + tot["othermc"]=0 + tot["data"]=0 + tot["2bmc_2btag"]=0 + tot["othermc_2btag"]=0 + tot["data_2btag"]=0 + for iBin in range(15): + #print iBin,histBkg.GetBinContent(iBin) + if iBin%4==1: + tot["2bmc"]=tot["2bmc"]+histBkg.GetBinContent(iBin) + if iBin>9: + tot["2bmc_2btag"]=tot["2bmc_2btag"]+histBkg.GetBinContent(iBin) + else: + tot["othermc"]=tot["othermc"]+histBkg.GetBinContent(iBin) + if iBin>9: + tot["othermc_2btag"]=tot["othermc_2btag"]+histBkg.GetBinContent(iBin) + twoTagDataHist=basename+"/"+basename + hist=tFile.Get(twoTagDataHist) + for iBin in range(15): + #print iBin,hist.GetBinContent(iBin),histBkg.GetBinContent(iBin) + tot["data"]=tot["data"]+hist.GetBinContent(iBin) if iBin>9: - tot["othermc_2btag"]=tot["othermc_2btag"]+histBkg.GetBinContent(iBin) - twoTagDataHist="emu_"+tag+syst+"/emu_"+tag+syst - hist=tFile.Get(twoTagDataHist) - for iBin in range(15): - #print iBin,hist.GetBinContent(iBin),histBkg.GetBinContent(iBin) - tot["data"]=tot["data"]+hist.GetBinContent(iBin) - if iBin>9: - tot["data_2btag"]=tot["data_2btag"]+hist.GetBinContent(iBin) - - #print tot - - print histBkg.Integral(),hist.Integral(),hist.Integral()/histBkg.Integral() - - for name in tot: - if name.find("mc")!=-1: - tot[name]=tot[name]/histBkg.Integral() + tot["data_2btag"]=tot["data_2btag"]+hist.GetBinContent(iBin) + + print tot + + print histBkg.Integral(),hist.Integral(), + if histBkg.Integral()>0: + hist.Integral()/histBkg.Integral() else: - tot[name]=tot[name]/hist.Integral() + print "division by zero avoided" + continue + + for name in tot: + if name.find("mc")!=-1: + tot[name]=tot[name]/histBkg.Integral() + else: + tot[name]=tot[name]/hist.Integral() + + + print tot + try: + eff=math.sqrt((tot["data_2btag"]-tot["othermc_2btag"])/tot["2bmc"]) + effmc=math.sqrt(tot["2bmc_2btag"]/tot["2bmc"]) - - #print tot - - eff=math.sqrt((tot["data_2btag"]-tot["othermc_2btag"])/tot["2bmc"]) - effmc=math.sqrt(tot["2bmc_2btag"]/tot["2bmc"]) - - effUp=math.sqrt((tot["data_2btag"]-tot["othermc_2btag"]*1.5)/tot["2bmc"]) - effDown=math.sqrt((tot["data_2btag"]-tot["othermc_2btag"]*0.5)/tot["2bmc"]) - - #print tot["2bmc_2btag"],tot["2bmc"],tot["data_2btag"],tot["othermc_2btag"],eff,effmc,eff/effmc,effUp/effmc,effDown/effmc - results[tag+syst]=[tot["2bmc_2btag"],tot["2bmc"],tot["data_2btag"],tot["othermc_2btag"],eff,effmc,eff/effmc,effUp/effmc,effDown/effmc] + effUp=math.sqrt((tot["data_2btag"]-tot["othermc_2btag"]*1.5)/tot["2bmc"]) + effDown=math.sqrt((tot["data_2btag"]-tot["othermc_2btag"]*0.5)/tot["2bmc"]) + results[tag+syst+ptLabel]=[tot["2bmc_2btag"],tot["2bmc"],tot["data_2btag"],tot["othermc_2btag"],eff,effmc,eff/effmc,effUp/effmc,effDown/effmc] + + if syst=="" and ptBin>-1: + pt1=lowerPtBinEdges[ptBin] + pt2=1000 + if ptBin!=lowerPtBinEdges[-1]: + pt2=lowerPtBinEdges[ptBin+1] + tgraphs[tag].SetPoint(iGraph,(pt1+pt2)/2, eff/effmc) + tgraphs[tag].SetPointError(iGraph,(pt2-pt1)/2, abs(effUp-effDown)/effmc) + iGraph=iGraph+1 + + except: + print "couldn't do this one" + continue + - - if syst=="": - print hist.Integral(),histBkg.Integral(),histBkg.GetEffectiveEntries() - statError[tag+"data_2btag"] = math.sqrt( tot["data_2btag"] * (1 - tot["data_2btag"]) / hist.Integral()) - statError[tag+"othermc_2btag"] = math.sqrt( tot["othermc_2btag"] * (1 - tot["othermc_2btag"]) / histBkg.GetEffectiveEntries()) - statError[tag+"2bmc"] = math.sqrt( tot["2bmc"] * (1 - tot["2bmc"]) / histBkg.GetEffectiveEntries()) - statError[tag+"eff"] = math.sqrt(math.pow(eff*statError[tag+"2bmc"],2) + (math.pow(statError[tag+"data_2btag"],2)+math.pow(statError[tag+"othermc_2btag"],2) ) / math.pow(eff,2) )/(2*tot["2bmc"]) + + if syst=="": + print hist.Integral(),histBkg.Integral(),histBkg.GetEffectiveEntries() + statError[tag+"data_2btag"] = math.sqrt( tot["data_2btag"] * (1 - tot["data_2btag"]) / hist.Integral()) + statError[tag+"othermc_2btag"] = math.sqrt( tot["othermc_2btag"] * (1 - tot["othermc_2btag"]) / histBkg.GetEffectiveEntries()) + statError[tag+"2bmc"] = math.sqrt( tot["2bmc"] * (1 - tot["2bmc"]) / histBkg.GetEffectiveEntries()) + statError[tag+"eff"] = math.sqrt(math.pow(eff*statError[tag+"2bmc"],2) + (math.pow(statError[tag+"data_2btag"],2)+math.pow(statError[tag+"othermc_2btag"],2) ) / math.pow(eff,2) )/(2*tot["2bmc"]) + tgraphs[tag].Draw("APE") + can.Update() + raw_input() raw_input() print "tag,", @@ -109,11 +154,15 @@ else: print syst+",", print -for tag in tags: - print tag+",", + +sortedResults=results.keys() +sortedResults.sort() +for result in sortedResults: + print result, for syst in systs: if syst=="": - print results[tag+syst][6],",",statError[tag+"eff"],",",results[tag+syst][7],",",results[tag+syst][8],",", + print results[result][6],",",results[result][7],",",results[result][8],",", + #print results[tag+syst+ptLabel][6],",",statError[tag+"eff"],",",results[tag+syst+ptLabel][7],",",results[tag+syst+ptLabel][8],",", else: - print results[tag+syst][6],",", + print results[result][6],",", print From ad44dbec0376d79ff9b21e92fef9a01b97fc5d33 Mon Sep 17 00:00:00 2001 From: capalmer Date: Mon, 15 Apr 2019 23:32:39 +0200 Subject: [PATCH 6/7] update to write png's and ROOT file with tgraphs --- test/ttbar/twoTag.py | 9 +++++++-- 1 file changed, 7 insertions(+), 2 deletions(-) diff --git a/test/ttbar/twoTag.py b/test/ttbar/twoTag.py index 0cd5d4632fa..51ddd8ddc46 100644 --- a/test/ttbar/twoTag.py +++ b/test/ttbar/twoTag.py @@ -44,6 +44,8 @@ def ReturnPtLabel(iPT): tgraphs={} can=ROOT.TCanvas("can","",800,800) +outputFile=ROOT.TFile("twoTagSFGraphs.root","RECREATE") + for syst in systs: #print syst @@ -144,9 +146,12 @@ def ReturnPtLabel(iPT): statError[tag+"eff"] = math.sqrt(math.pow(eff*statError[tag+"2bmc"],2) + (math.pow(statError[tag+"data_2btag"],2)+math.pow(statError[tag+"othermc_2btag"],2) ) / math.pow(eff,2) )/(2*tot["2bmc"]) tgraphs[tag].Draw("APE") can.Update() - raw_input() + can.SaveAs(tag+".png") + tgraphs[tag].SetName(tag) + outputFile.cd() + tgraphs[tag].Write() -raw_input() +outputFile.Close() print "tag,", for syst in systs: if syst=="": From 0cc866da0a049f795ba7c0e11fe57244c0480b81 Mon Sep 17 00:00:00 2001 From: capalmer Date: Wed, 17 Apr 2019 21:25:24 +0200 Subject: [PATCH 7/7] update instructions and settings in twoTag.py script --- test/ttbar/README.md | 6 +++--- test/ttbar/twoTag.py | 10 +--------- 2 files changed, 4 insertions(+), 12 deletions(-) diff --git a/test/ttbar/README.md b/test/ttbar/README.md index 06e5e63d977..3ee9ec853e5 100644 --- a/test/ttbar/README.md +++ b/test/ttbar/README.md @@ -59,17 +59,17 @@ The scale factor can be applied to all channels as it reflects possible mismodel ### Running local analysis for 2TagCount After b-tag analyzer trees are produced, the first two steps of the 2TagCount analysis are the same. The first is to run the analysis and make histograms (and trees): ``` -python runTTbarAnalysis.py -i /store/group/phys_btag/Commissioning/TTbar/Moriond19_2018_StructuredDir/ --json data/samples_Run2018.json -o /tmp/MYUSERNAME/Moriond19_Run2018_april5 -n 70 +python runTTbarAnalysis.py -i /store/group/phys_btag/Commissioning/TTbar/Moriond19_2018_StructuredDir/ --json data/samples_Run2018.json -o /tmp/MYUSERNAME/Moriond19_Run2018 -n 70 --doTwoTag 1 ``` Re-organize plots and make data, mc comparison canvases (inclusive 2TagCount plots are in): ``` -python plotter.py -i /tmp/MYUSERNAME/Moriond19_Run2018_april5 --json data/samples_Run2018.json --lumi 58000 +python plotter.py -i /tmp/MYUSERNAME/Moriond19_Run2018 --json data/samples_Run2018.json --lumi 58000 ``` Take the 2TagCount histograms and compute efficiency scale factor comparing the following: 1) b-tagging effiency for each working point for events with two b-jets at truth level (MC eff) and 2) data minus non-double b-jet MC divided by the number of MC events with two matched particle level b-jets: ``` -python twoTag.py /tmp/MYUSERNAME/Moriond19_Run2018_april5/plots/plotter.root +python twoTag.py /tmp/MYUSERNAME/Moriond19_Run2018/plots/plotter.root ``` diff --git a/test/ttbar/twoTag.py b/test/ttbar/twoTag.py index 51ddd8ddc46..0bdc2afc1b5 100644 --- a/test/ttbar/twoTag.py +++ b/test/ttbar/twoTag.py @@ -20,21 +20,13 @@ def ReturnPtLabel(iPT): return ptBinLabel -tags=["twoTags_deepCSVL"] -tags=["twoTags_deepCSVL","only2_twoTags_deepCSVL","twoTags_deepCSVM","only2_twoTags_deepCSVM","twoTags_deepCSVT","only2_twoTags_deepCSVT","twoTags_deepFlavourL","only2_twoTags_deepFlavourL","twoTags_deepFlavourM","only2_twoTags_deepFlavourM","twoTags_deepFlavourT","only2_twoTags_deepFlavourT","twoTags_CSVv2L","only2_twoTags_CSVv2L","twoTags_CSVv2M","only2_twoTags_CSVv2M","twoTags_CSVv2T","only2_twoTags_CSVv2T"] -tags=["twoTags_deepCSVL","twoTags_deepCSVM","twoTags_deepCSVT","twoTags_deepFlavourL","only2_twoTags_deepFlavourL","twoTags_deepFlavourM","only2_twoTags_deepFlavourM","twoTags_deepFlavourT","only2_twoTags_deepFlavourT"] +tags=["twoTags_deepCSVL","only2_twoTags_deepCSVL","twoTags_deepCSVM","only2_twoTags_deepCSVM","twoTags_deepCSVT","only2_twoTags_deepCSVT","twoTags_deepFlavourL","only2_twoTags_deepFlavourL","twoTags_deepFlavourM","only2_twoTags_deepFlavourM","twoTags_deepFlavourT","only2_twoTags_deepFlavourT"] tags.sort() samples=["t#bar{t} SL","VV","tW","t#bar{t} DL"] samples=["t#bar{t} SL","t#bar{t} DL","t#bar{t} AH","ZZ","WZ","WW","tW","t#bar{t} DL","DY"] systs=["","_fsrConHi","_fsrConLo","_fsrDefHi","_fsrDefLo","_fsrRedHi","_fsrRedLo","_isrConHi","_isrConLo","_isrDefHi","_isrDefLo","_isrRedHi","_isrRedLo","_qcdScaleHi","_qcdScaleLo"] systs=[""] -# KEY: TH1F emu_only2_twoTags_CSVv2L_t#bar{t} AH;1 t#bar{t} AH -# KEY: TH1F emu_only2_twoTags_CSVv2L_other;1 other -# KEY: TH1F emu_only2_twoTags_CSVv2L_t#bar{t} SL;1 t#bar{t} SL -# KEY: TH1F emu_only2_twoTags_CSVv2L_tW;1 tW -# KEY: TH1F emu_only2_twoTags_CSVv2L_t#bar{t} DL;1 t#bar{t} DL - statError={}