-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathBTagSFUtil.C
167 lines (99 loc) · 4.93 KB
/
BTagSFUtil.C
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
#include "BTagSFUtil.h"
//#include "BTagSFPayloadMoriond13.C"
#include "BTagSFPayloadWinter13.C"
//#include "BTagSFPayloadWinter13NoTTbar.C"
#include "BTagEfficienciesTTbarSummer12.C"
//#include "FastSimCorrectionFactors2011.C"
#include "FastSimCorrectionFactorsSummer12.C"
BTagSFUtil::BTagSFUtil(TString BTagAlgorithm, TString DataPeriod, TString FastSimDataset, int Seed) {
rand_ = new TRandom3(Seed);
GetBTagPayload(BTagAlgorithm, DataPeriod);
GetFastSimPayload(BTagAlgorithm, FastSimDataset);
if (TaggerCut<0.)
cout << "BTagSFUtil: " << BTagAlgorithm << " not a supported b-tagging algorithm" << endl;
}
BTagSFUtil::~BTagSFUtil() {
delete rand_;
}
float BTagSFUtil::ScaleFactorB(float JetPt, int SystematicFlag) {
if (JetPt<BTagPtBinEdge[0]) { cout << "SFb is not available for jet pt<20 GeV" << endl; return -1.; }
float x = JetPt;
if (JetPt>BTagPtBinEdge[nBTagPtBins]) x = BTagPtBinEdge[nBTagPtBins];
int JetPtBin = -1;
for (int ptbtagbin = 0; ptbtagbin<nBTagPtBins; ptbtagbin++)
if (x>=BTagPtBinEdge[ptbtagbin]) JetPtBin++;
if (JetPt>x) SystematicFlag *= 2;
float SFb = funSFb->Eval(x) + SystematicFlag*SFb_error[JetPtBin];
return SFb;
}
float BTagSFUtil::ScaleFactorLight(float JetPt, float JetEta, int SystematicFlag) {
if (JetPt<BTagPtBinEdge[0]) { cout << "SFlight is not available for jet pt<20 GeV" << endl; return -1.; }
int JetEtaBin = -1;
for (int etabtagbin = 0; etabtagbin<nBTagEtaBins; etabtagbin++)
if (fabs(JetEta)>=BTagEtaBinEdge[etabtagbin]) JetEtaBin++;
float MaxJetPtLight = funSFlight[JetEtaBin][SystematicFlag+1]->GetMaximumX();
float x = JetPt;
if (JetPt>MaxJetPtLight) x = MaxJetPtLight;
return funSFlight[JetEtaBin][SystematicFlag+1]->Eval(x);
}
float BTagSFUtil::ScaleFactorJet(int JetFlavor, float JetPt, float JetEta, int SystematicFlag) {
float SF = -1.;
if (JetPt<BTagPtBinEdge[0]) { cout << "SF is not available for jet pt<20 GeV" << endl; return -1.; }
if (fabs(JetEta)>2.4) { cout << "SF is not available for jet |eta|>2.4" << endl; return -1.; }
if (abs(JetFlavor)==5) SF = ScaleFactorB(JetPt, SystematicFlag);
else if (abs(JetFlavor)==4) SF = ScaleFactorB(JetPt, 2*SystematicFlag);
else SF = ScaleFactorLight(JetPt, JetEta, SystematicFlag);
if (SF==-1.) cout << "Jet parameter out of BTV prescriptions" << endl;
return SF;
}
float BTagSFUtil::FastSimCorrectionFactor(int JetFlavor, float JetPt, float JetEta, int FastSimSystematic) {
float CF = 1.;
if (JetPt<FastSimPtBinEdge[0]) { cout << "CF is not available for jet pt<" << FastSimPtBinEdge[0] << " GeV" << endl; return -1.; }
if (fabs(JetEta)>2.4) { cout << "CF is not available for jet |eta|>2.4" << endl; return -1.; }
int JetFlavorFlag = 2;
if (abs(JetFlavor)==4) JetFlavorFlag = 1;
else if (abs(JetFlavor)==5) JetFlavorFlag = 0;
int JetPtBin = -1;
for (int ptbin = 0; ptbin<nFastSimPtBins; ptbin++)
if (JetPt>=FastSimPtBinEdge[ptbin]) JetPtBin++;
if (JetPt>=FastSimPtBinEdge[nFastSimPtBins]) FastSimSystematic *= 2;
int JetEtaBin = -1;
for (int etabin = 0; etabin<nFastSimEtaBins[JetFlavorFlag]; etabin++)
if (fabs(JetEta)>=FastSimEtaBinEdge[etabin][JetFlavorFlag]) JetEtaBin++;
CF = FastSimCF[JetPtBin][JetEtaBin][JetFlavorFlag] + FastSimSystematic*FastSimCF_error[JetPtBin][JetEtaBin][JetFlavorFlag];
if (CF<0.) CF = 0.; // Effect of large uncertainties on light CFs!
return CF;
}
float BTagSFUtil::JetTagEfficiency(int JetFlavor, float JetPt, float JetEta) {
if (abs(JetFlavor)==5) return TagEfficiencyB(JetPt, JetEta);
else if (abs(JetFlavor)==4) return TagEfficiencyC(JetPt, JetEta);
else return TagEfficiencyLight(JetPt, JetEta);
}
float BTagSFUtil::GetJetSF(int JetFlavor, float JetPt, float JetEta, int SystematicFlag, int FastSimSystematic) {
float Btag_SF = ScaleFactorJet(JetFlavor, JetPt, JetEta, SystematicFlag);
if (IsFastSimDataset)
Btag_SF *= FastSimCorrectionFactor(JetFlavor, JetPt, JetEta, FastSimSystematic);
return Btag_SF;
}
bool BTagSFUtil::IsTagged(float JetDiscriminant, int JetFlavor, float JetPt, float JetEta, int SystematicFlag, int FastSimSystematic) {
bool isBTagged = JetDiscriminant>TaggerCut;
if (JetFlavor==-999999) return isBTagged; // Data: no correction needed
bool newBTag = isBTagged;
float Btag_SF = GetJetSF(JetFlavor, JetPt, JetEta, SystematicFlag, FastSimSystematic);
if (Btag_SF == 1) return newBTag; //no correction needed
//throw die
float coin = rand_->Uniform(1.);
if(Btag_SF > 1){ // use this if SF>1
if( !isBTagged ) {
float Btag_eff = JetTagEfficiency(JetFlavor, JetPt, fabs(JetEta));
//fraction of jets that need to be upgraded
float mistagPercent = (1.0 - Btag_SF) / (1.0 - (1./Btag_eff) );
//upgrade to tagged
if( coin < mistagPercent ) {newBTag = true;}
}
}else{ // use this if SF<1
//downgrade tagged to untagged
if( isBTagged && coin > Btag_SF ) {newBTag = false;}
}
return newBTag;
}