CMS 3D CMS Logo

BTagEntry.cc
Go to the documentation of this file.
1 #include <algorithm>
2 #include <sstream>
3 
6 
9  std::string measurement_type,
10  std::string sys_type,
11  JetFlavor jf,
12  float eta_min,
13  float eta_max,
14  float pt_min,
15  float pt_max,
16  float discr_min,
17  float discr_max
18 ):
19  operatingPoint(op),
20  measurementType(measurement_type),
21  sysType(sys_type),
22  jetFlavor(jf),
23  etaMin(eta_min),
24  etaMax(eta_max),
25  ptMin(pt_min),
26  ptMax(pt_max),
27  discrMin(discr_min),
28  discrMax(discr_max)
29 {
31  measurementType.begin(), ::tolower);
32  std::transform(sysType.begin(), sysType.end(),
33  sysType.begin(), ::tolower);
34 }
35 
36 BTagEntry::BTagEntry(const std::string& csvLine, bool validate) {
37  // make tokens
38  std::stringstream buff(csvLine);
39  std::vector<std::string> vec;
40  std::string token;
41  while (std::getline(buff, token, ","[0])) {
42  token = BTagEntry::trimStr(token);
43  if (token.empty()) {
44  continue;
45  }
46  vec.push_back(token);
47  }
48  if (vec.size() != 11) {
49  throw cms::Exception("BTagCalibration")
50  << "Invalid csv line; num tokens != 11: "
51  << csvLine;
52  }
53 
54  // clean string values
55  char chars[] = " \"\n";
56  for (unsigned int i = 0; i < strlen(chars); ++i) {
57  vec[1].erase(remove(vec[1].begin(),vec[1].end(),chars[i]),vec[1].end());
58  vec[2].erase(remove(vec[2].begin(),vec[2].end(),chars[i]),vec[2].end());
59  vec[10].erase(remove(vec[10].begin(),vec[10].end(),chars[i]),vec[10].end());
60  }
61 
62  // make formula
63  formula = vec[10];
64  if (validate) {
65  TF1 f1("", formula.c_str()); // compile formula to check validity
66  if (f1.IsZombie()) {
67  throw cms::Exception("BTagCalibration") << "Invalid csv line; formula does not compile: " << csvLine;
68  }
69  }
70 
71  // make parameters
72  unsigned op = stoi(vec[0]);
73  if (op > 3) {
74  throw cms::Exception("BTagCalibration")
75  << "Invalid csv line; OperatingPoint > 3: "
76  << csvLine;
77  }
78  unsigned jf = stoi(vec[3]);
79  if (jf > 2) {
80  throw cms::Exception("BTagCalibration")
81  << "Invalid csv line; JetFlavor > 2: "
82  << csvLine;
83  }
86  vec[1],
87  vec[2],
89  stof(vec[4]),
90  stof(vec[5]),
91  stof(vec[6]),
92  stof(vec[7]),
93  stof(vec[8]),
94  stof(vec[9])
95  );
96 }
97 
99  formula(func),
100  params(p)
101 {
102  TF1 f1("", formula.c_str()); // compile formula to check validity
103  if (f1.IsZombie()) {
104  throw cms::Exception("BTagCalibration")
105  << "Invalid func string; formula does not compile: "
106  << func;
107  }
108 }
109 
111  formula(std::string(func->GetExpFormula("p").Data())),
112  params(p)
113 {
114  if (func->IsZombie()) {
115  throw cms::Exception("BTagCalibration")
116  << "Invalid TF1 function; function is zombie: "
117  << func->GetName();
118  }
119 }
120 
121 // Creates chained step functions like this:
122 // "<prevous_bin> : x<bin_high_bound ? bin_value : <next_bin>"
123 // e.g. "x<0 ? 1 : x<1 ? 2 : x<2 ? 3 : 4"
125  int nbins = hist->GetNbinsX();
126  TAxis const* axis = hist->GetXaxis();
127  std::stringstream buff;
128  buff << "x<" << axis->GetBinLowEdge(1) << " ? 0. : "; // default value
129  for (int i=1; i<nbins+1; ++i) {
130  char tmp_buff[50];
131  sprintf(tmp_buff,
132  "x<%g ? %g : ", // %g is the smaller one of %e or %f
133  axis->GetBinUpEdge(i),
134  hist->GetBinContent(i));
135  buff << tmp_buff;
136  }
137  buff << 0.; // default value
138  return buff.str();
139 }
140 
141 // Creates step functions making a binary search tree:
142 // "x<mid_bin_bound ? (<left side tree>) : (<right side tree>)"
143 // e.g. "x<2 ? (x<1 ? (x<0 ? 0:0.1) : (1)) : (x<4 ? (x<3 ? 2:3) : (0))"
144 std::string th1ToFormulaBinTree(const TH1* hist, int start=0, int end=-1) {
145  if (end == -1) { // initialize
146  start = 0.;
147  end = hist->GetNbinsX()+1;
148  TH1* h2 = (TH1*) hist->Clone();
149  h2->SetBinContent(start, 0); // kill underflow
150  h2->SetBinContent(end, 0); // kill overflow
152  delete h2;
153  return res;
154  }
155  if (start == end) { // leave is reached
156  char tmp_buff[20];
157  sprintf(tmp_buff, "%g", hist->GetBinContent(start));
158  return std::string(tmp_buff);
159  }
160  if (start == end - 1) { // no parenthesis for neighbors
161  char tmp_buff[70];
162  sprintf(tmp_buff,
163  "x<%g ? %g:%g",
164  hist->GetXaxis()->GetBinUpEdge(start),
165  hist->GetBinContent(start),
166  hist->GetBinContent(end));
167  return std::string(tmp_buff);
168  }
169 
170  // top-down recursion
171  std::stringstream buff;
172  int mid = (end-start)/2 + start;
173  char tmp_buff[25];
174  sprintf(tmp_buff,
175  "x<%g ? (",
176  hist->GetXaxis()->GetBinUpEdge(mid));
177  buff << tmp_buff
178  << th1ToFormulaBinTree(hist, start, mid)
179  << ") : ("
180  << th1ToFormulaBinTree(hist, mid+1, end)
181  << ")";
182  return buff.str();
183 }
184 
186  params(p)
187 {
188  int nbins = hist->GetNbinsX();
189  TAxis const* axis = hist->GetXaxis();
190 
191  // overwrite bounds with histo values
193  params.discrMin = axis->GetBinLowEdge(1);
194  params.discrMax = axis->GetBinUpEdge(nbins);
195  } else {
196  params.ptMin = axis->GetBinLowEdge(1);
197  params.ptMax = axis->GetBinUpEdge(nbins);
198  }
199 
200  // balanced full binary tree height = ceil(log(2*n_leaves)/log(2))
201  // breakes even around 10, but lower values are more propable in pt-spectrum
202  if (nbins < 15) {
203  formula = th1ToFormulaLin(hist);
204  } else {
206  }
207 
208  // compile formula to check validity
209  TF1 f1("", formula.c_str());
210  if (f1.IsZombie()) {
211  throw cms::Exception("BTagCalibration")
212  << "Invalid histogram; formula does not compile (>150 bins?): "
213  << hist->GetName();
214  }
215 }
216 
218 {
219  return "OperatingPoint, "
220  "measurementType, "
221  "sysType, "
222  "jetFlavor, "
223  "etaMin, "
224  "etaMax, "
225  "ptMin, "
226  "ptMax, "
227  "discrMin, "
228  "discrMax, "
229  "formula \n";
230 }
231 
233 {
234  std::stringstream buff;
235  buff << params.operatingPoint
236  << ", " << params.measurementType
237  << ", " << params.sysType
238  << ", " << params.jetFlavor
239  << ", " << params.etaMin
240  << ", " << params.etaMax
241  << ", " << params.ptMin
242  << ", " << params.ptMax
243  << ", " << params.discrMin
244  << ", " << params.discrMax
245  << ", \"" << formula
246  << "\" \n";
247  return buff.str();
248 }
249 
251  size_t s = str.find_first_not_of(" \n\r\t");
252  size_t e = str.find_last_not_of (" \n\r\t");
253 
254  if((std::string::npos == s) || (std::string::npos == e))
255  return "";
256  else
257  return str.substr(s, e-s+1);
258 }
Definition: start.py:1
std::string sysType
Definition: BTagEntry.h:41
Parameters(OperatingPoint op=OP_TIGHT, std::string measurement_type="comb", std::string sys_type="central", JetFlavor jf=FLAV_B, float eta_min=-99999., float eta_max=99999., float pt_min=0., float pt_max=99999., float discr_min=0., float discr_max=99999.)
Definition: BTagEntry.cc:7
std::string th1ToFormulaBinTree(const TH1 *hist, int start=0, int end=-1)
Definition: BTagEntry.cc:144
static std::string trimStr(std::string str)
Definition: BTagEntry.cc:250
JetFlavor jetFlavor
Definition: BTagEntry.h:42
static std::string makeCSVHeader()
Definition: BTagEntry.cc:217
std::string th1ToFormulaLin(const TH1 *hist)
Definition: BTagEntry.cc:124
std::string makeCSVLine() const
Definition: BTagEntry.cc:232
Definition: Electron.h:6
vector< ParameterSet > Parameters
OperatingPoint operatingPoint
Definition: BTagEntry.h:39
std::string formula
Definition: BTagEntry.h:78
#define end
Definition: vmac.h:39
OperatingPoint
Definition: BTagEntry.h:27
Parameters params
Definition: BTagEntry.h:79
BTagEntry()
Definition: BTagEntry.h:67
#define begin
Definition: vmac.h:32
std::string measurementType
Definition: BTagEntry.h:40
#define str(s)