17 |
|
superclusters_ (cfg.getParameter<edm::InputTag> ("superclusters")), |
18 |
|
triggers_ (cfg.getParameter<edm::InputTag> ("triggers")), |
19 |
|
puFile_ (cfg.getParameter<std::string> ("puFile")), |
20 |
+ |
deadEcalFile_ (cfg.getParameter<std::string> ("deadEcalFile")), |
21 |
|
dataPU_ (cfg.getParameter<std::string> ("dataPU")), |
22 |
|
dataset_ (cfg.getParameter<std::string> ("dataset")), |
23 |
|
datasetType_ (cfg.getParameter<std::string> ("datasetType")), |
30 |
|
TH1::SetDefaultSumw2 (); |
31 |
|
|
32 |
|
//create pile-up reweighting object, if necessary |
33 |
< |
if(datasetType_ != "data") puWeight_ = new PUWeight (puFile_, dataPU_, dataset_); |
33 |
> |
// if(datasetType_ != "data") puWeight_ = new PUWeight (puFile_, dataPU_, dataset_); |
34 |
|
|
35 |
|
// Construct Cutflow Objects. These store the results of cut decisions and |
36 |
|
// handle filling cut flow histograms. |
38 |
|
std::vector<TFileDirectory> directories; |
39 |
|
|
40 |
|
//always get vertex collection so we can assign the primary vertex in the event |
41 |
< |
allNecessaryObjects.push_back("primaryvertexs"); |
41 |
> |
|
42 |
> |
objectsToGet.push_back("primaryvertexs"); |
43 |
|
//always make the plot of number of primary verticex (to check pile-up reweighting) |
44 |
|
objectsToPlot.push_back("primaryvertexs"); |
45 |
|
|
46 |
+ |
//always get the MC particles to do GEN-matching |
47 |
+ |
objectsToGet.push_back("mcparticles"); |
48 |
+ |
|
49 |
|
//always get the event collection to do pile-up reweighting |
50 |
< |
allNecessaryObjects.push_back("events"); |
50 |
> |
objectsToGet.push_back("events"); |
51 |
|
|
52 |
|
//parse the histogram definitions |
53 |
|
for(uint currentHistogramSet = 0; currentHistogramSet != histogramSets_.size(); currentHistogramSet++){ |
55 |
|
string tempInputCollection = histogramSets_.at(currentHistogramSet).getParameter<string> ("inputCollection"); |
56 |
|
if(tempInputCollection == "muon-electron pairs") tempInputCollection = "electron-muon pairs"; |
57 |
|
if(tempInputCollection.find("pairs")==std::string::npos){ //just a single object |
58 |
+ |
objectsToGet.push_back(tempInputCollection); |
59 |
|
objectsToPlot.push_back(tempInputCollection); |
60 |
< |
allNecessaryObjects.push_back(tempInputCollection); |
60 |
> |
objectsToCut.push_back(tempInputCollection); |
61 |
|
} |
62 |
< |
else{//pair of objects, need to add them both to the things to allNecessaryObjects |
62 |
> |
else{//pair of objects, need to add them both to the things to objectsToGet |
63 |
|
int dashIndex = tempInputCollection.find("-"); |
64 |
|
int spaceIndex = tempInputCollection.find(" "); |
65 |
|
int secondWordLength = spaceIndex - dashIndex; |
66 |
< |
allNecessaryObjects.push_back(tempInputCollection); |
67 |
< |
allNecessaryObjects.push_back(tempInputCollection.substr(0,dashIndex)+"s"); |
68 |
< |
allNecessaryObjects.push_back(tempInputCollection.substr(dashIndex+1,secondWordLength-1)+"s"); |
66 |
> |
objectsToGet.push_back(tempInputCollection); |
67 |
> |
objectsToGet.push_back(tempInputCollection.substr(0,dashIndex)+"s"); |
68 |
> |
objectsToGet.push_back(tempInputCollection.substr(dashIndex+1,secondWordLength-1)+"s"); |
69 |
|
objectsToPlot.push_back(tempInputCollection); |
70 |
|
objectsToPlot.push_back(tempInputCollection.substr(0,dashIndex)+"s"); |
71 |
|
objectsToPlot.push_back(tempInputCollection.substr(dashIndex+1,secondWordLength-1)+"s"); |
72 |
+ |
objectsToCut.push_back(tempInputCollection); |
73 |
+ |
objectsToCut.push_back(tempInputCollection.substr(0,dashIndex)+"s"); |
74 |
+ |
objectsToCut.push_back(tempInputCollection.substr(dashIndex+1,secondWordLength-1)+"s"); |
75 |
+ |
|
76 |
|
} |
77 |
|
|
78 |
|
vector<edm::ParameterSet> histogramList_ (histogramSets_.at(currentHistogramSet).getParameter<vector<edm::ParameterSet> >("histograms")); |
90 |
|
|
91 |
|
} |
92 |
|
} |
93 |
< |
//make unique vector of objects we need to plot (so we can book a histogram with the number of each object |
93 |
> |
//make unique vector of objects we need to plot (so we can book a histogram with the number of each object) |
94 |
|
sort( objectsToPlot.begin(), objectsToPlot.end() ); |
95 |
|
objectsToPlot.erase( unique( objectsToPlot.begin(), objectsToPlot.end() ), objectsToPlot.end() ); |
96 |
|
|
97 |
|
|
98 |
|
|
99 |
+ |
//add histograms with the gen-matched id, mother id, and grandmother id |
100 |
+ |
for(uint currentObjectIndex = 0; currentObjectIndex != objectsToPlot.size(); currentObjectIndex++){ |
101 |
+ |
|
102 |
+ |
string currentObject = objectsToPlot.at(currentObjectIndex); |
103 |
+ |
if(currentObject != "muons" && currentObject != "electrons" && currentObject != "taus" && currentObject != "tracks" && currentObject != "photons" && currentObject != "superclusters") continue; |
104 |
+ |
|
105 |
+ |
histogram tempIdHisto; |
106 |
+ |
histogram tempMomIdHisto; |
107 |
+ |
histogram tempGmaIdHisto; |
108 |
+ |
|
109 |
+ |
tempIdHisto.inputCollection = currentObject; |
110 |
+ |
tempMomIdHisto.inputCollection = currentObject; |
111 |
+ |
tempGmaIdHisto.inputCollection = currentObject; |
112 |
+ |
|
113 |
+ |
currentObject = currentObject.substr(0, currentObject.size()-1); |
114 |
+ |
tempIdHisto.name = currentObject+"GenMatchId"; |
115 |
+ |
tempMomIdHisto.name = currentObject+"GenMatchMotherId"; |
116 |
+ |
tempGmaIdHisto.name = currentObject+"GenMatchGrandmotherId"; |
117 |
+ |
|
118 |
+ |
currentObject.at(0) = toupper(currentObject.at(0)); |
119 |
+ |
tempIdHisto.title = currentObject+" Gen-matched Particle"; |
120 |
+ |
tempMomIdHisto.title = currentObject+" Gen-matched Particle's Mother"; |
121 |
+ |
tempGmaIdHisto.title = currentObject+" Gen-matched Particle's Grandmother"; |
122 |
+ |
|
123 |
+ |
int maxNum = 24; |
124 |
+ |
vector<double> binVector; |
125 |
+ |
binVector.push_back(maxNum); |
126 |
+ |
binVector.push_back(0); |
127 |
+ |
binVector.push_back(maxNum); |
128 |
+ |
|
129 |
+ |
tempIdHisto.bins = binVector; |
130 |
+ |
tempIdHisto.inputVariables.push_back("genMatchedId"); |
131 |
+ |
tempMomIdHisto.bins = binVector; |
132 |
+ |
tempMomIdHisto.inputVariables.push_back("genMatchedMotherId"); |
133 |
+ |
tempGmaIdHisto.bins = binVector; |
134 |
+ |
tempGmaIdHisto.inputVariables.push_back("genMatchedGrandmotherId"); |
135 |
+ |
|
136 |
+ |
histograms.push_back(tempIdHisto); |
137 |
+ |
histograms.push_back(tempMomIdHisto); |
138 |
+ |
histograms.push_back(tempGmaIdHisto); |
139 |
+ |
|
140 |
+ |
} |
141 |
+ |
|
142 |
+ |
|
143 |
|
|
144 |
|
channel tempChannel; |
145 |
|
//loop over all channels (event selections) |
158 |
|
triggerNames = channels_.at(currentChannel).getParameter<vector<string> >("triggers"); |
159 |
|
for(uint trigger = 0; trigger!= triggerNames.size(); trigger++) |
160 |
|
tempChannel.triggers.push_back(triggerNames.at(trigger)); |
161 |
< |
allNecessaryObjects.push_back("triggers"); |
161 |
> |
objectsToGet.push_back("triggers"); |
162 |
|
} |
163 |
|
|
164 |
|
|
183 |
|
int numBinsElements = currentHistogram.bins.size(); |
184 |
|
int numInputVariables = currentHistogram.inputVariables.size(); |
185 |
|
|
132 |
– |
|
133 |
– |
|
186 |
|
if(numBinsElements != 3 && numBinsElements !=6) { |
187 |
|
std::cout << "Error: Didn't find correct number of bin specifications for histogram named '" << currentHistogram.name << "'\n"; |
188 |
|
exit(0); |
200 |
|
} |
201 |
|
|
202 |
|
|
203 |
+ |
if(currentHistogram.name.find("GenMatch")==std::string::npos) continue; |
204 |
+ |
|
205 |
+ |
// bin particle type |
206 |
+ |
// --- ------------- |
207 |
+ |
// 0 unmatched |
208 |
+ |
// 1 u |
209 |
+ |
// 2 d |
210 |
+ |
// 3 s |
211 |
+ |
// 4 c |
212 |
+ |
// 5 b |
213 |
+ |
// 6 t |
214 |
+ |
// 7 e |
215 |
+ |
// 8 mu |
216 |
+ |
// 9 tau |
217 |
+ |
// 10 nu |
218 |
+ |
// 11 g |
219 |
+ |
// 12 gamma |
220 |
+ |
// 13 Z |
221 |
+ |
// 14 W |
222 |
+ |
// 15 light meson |
223 |
+ |
// 16 K meson |
224 |
+ |
// 17 D meson |
225 |
+ |
// 18 B meson |
226 |
+ |
// 19 light baryon |
227 |
+ |
// 20 strange baryon |
228 |
+ |
// 21 charm baryon |
229 |
+ |
// 22 bottom baryon |
230 |
+ |
// 23 other |
231 |
+ |
|
232 |
+ |
vector<TString> labelArray; |
233 |
+ |
labelArray.push_back("unmatched"); |
234 |
+ |
labelArray.push_back("u"); |
235 |
+ |
labelArray.push_back("d"); |
236 |
+ |
labelArray.push_back("s"); |
237 |
+ |
labelArray.push_back("c"); |
238 |
+ |
labelArray.push_back("b"); |
239 |
+ |
labelArray.push_back("t"); |
240 |
+ |
labelArray.push_back("e"); |
241 |
+ |
labelArray.push_back("#mu"); |
242 |
+ |
labelArray.push_back("#tau"); |
243 |
+ |
labelArray.push_back("#nu"); |
244 |
+ |
labelArray.push_back("g"); |
245 |
+ |
labelArray.push_back("#gamma"); |
246 |
+ |
labelArray.push_back("Z"); |
247 |
+ |
labelArray.push_back("W"); |
248 |
+ |
labelArray.push_back("light meson"); |
249 |
+ |
labelArray.push_back("K meson"); |
250 |
+ |
labelArray.push_back("D meson"); |
251 |
+ |
labelArray.push_back("B meson"); |
252 |
+ |
labelArray.push_back("light baryon"); |
253 |
+ |
labelArray.push_back("strange baryon"); |
254 |
+ |
labelArray.push_back("charm baryon"); |
255 |
+ |
labelArray.push_back("bottom baryon"); |
256 |
+ |
labelArray.push_back("other"); |
257 |
+ |
|
258 |
+ |
for(int bin = 0; bin !=currentHistogram.bins.at(0); bin++){ |
259 |
+ |
oneDHists_.at(currentChannel)[currentHistogram.name]->GetXaxis()->SetBinLabel(bin+1,labelArray.at(bin)); |
260 |
+ |
} |
261 |
|
} |
262 |
+ |
|
263 |
|
//book a histogram for the number of each object type to be plotted |
264 |
+ |
|
265 |
|
for (uint currentObjectIndex = 0; currentObjectIndex != objectsToPlot.size(); currentObjectIndex++){ |
266 |
|
string currentObject = objectsToPlot.at(currentObjectIndex); |
267 |
|
int maxNum = 10; |
268 |
|
if(currentObject == "mcparticles") maxNum = 50; |
269 |
|
else if(currentObject == "primaryvertexs") maxNum = 50; |
270 |
< |
else if(currentObject == "muon-muon pairs") currentObject = "dimuonPairs"; |
270 |
> |
else if(currentObject == "muon-muon pairs") |
271 |
> |
{ |
272 |
> |
currentObject = "dimuonPairs"; |
273 |
> |
} |
274 |
|
else if(currentObject == "electron-electron pairs") currentObject = "dielectronPairs"; |
275 |
|
else if(currentObject == "electron-muon pairs") currentObject = "electronMuonPairs"; |
276 |
|
|
288 |
|
} |
289 |
|
|
290 |
|
|
291 |
+ |
|
292 |
+ |
|
293 |
+ |
|
294 |
|
//get list of cuts for this channel |
295 |
|
vector<edm::ParameterSet> cuts_ (channels_.at(currentChannel).getParameter<vector<edm::ParameterSet> >("cuts")); |
296 |
|
|
301 |
|
string tempInputCollection = cuts_.at(currentCut).getParameter<string> ("inputCollection"); |
302 |
|
tempCut.inputCollection = tempInputCollection; |
303 |
|
if(tempInputCollection.find("pairs")==std::string::npos){ //just a single object |
304 |
< |
allNecessaryObjects.push_back(tempInputCollection); |
304 |
> |
objectsToGet.push_back(tempInputCollection); |
305 |
> |
objectsToCut.push_back(tempInputCollection); |
306 |
|
} |
307 |
< |
else{//pair of objects, need to add them both to the things to allNecessaryObjects |
307 |
> |
else{//pair of objects, need to add them both to the things to objectsToGet |
308 |
|
int dashIndex = tempInputCollection.find("-"); |
309 |
|
int spaceIndex = tempInputCollection.find(" "); |
310 |
|
int secondWordLength = spaceIndex - dashIndex; |
311 |
< |
allNecessaryObjects.push_back(tempInputCollection); |
312 |
< |
allNecessaryObjects.push_back(tempInputCollection.substr(0,dashIndex)+"s"); |
313 |
< |
allNecessaryObjects.push_back(tempInputCollection.substr(dashIndex+1,secondWordLength-1)+"s"); |
311 |
> |
objectsToGet.push_back(tempInputCollection); |
312 |
> |
objectsToGet.push_back(tempInputCollection.substr(0,dashIndex)+"s"); |
313 |
> |
objectsToGet.push_back(tempInputCollection.substr(dashIndex+1,secondWordLength-1)+"s"); |
314 |
> |
objectsToCut.push_back(tempInputCollection); |
315 |
> |
objectsToCut.push_back(tempInputCollection.substr(0,dashIndex)+"s"); |
316 |
> |
objectsToCut.push_back(tempInputCollection.substr(dashIndex+1,secondWordLength-1)+"s"); |
317 |
|
|
318 |
|
} |
319 |
|
|
383 |
|
}//end loop over channels |
384 |
|
|
385 |
|
|
386 |
< |
//make unique vector of objects we need to cut on (so we make sure to get them from the event) |
387 |
< |
sort( allNecessaryObjects.begin(), allNecessaryObjects.end() ); |
388 |
< |
allNecessaryObjects.erase( unique( allNecessaryObjects.begin(), allNecessaryObjects.end() ), allNecessaryObjects.end() ); |
386 |
> |
//make unique vector of objects we need to get from the event |
387 |
> |
sort( objectsToGet.begin(), objectsToGet.end() ); |
388 |
> |
objectsToGet.erase( unique( objectsToGet.begin(), objectsToGet.end() ), objectsToGet.end() ); |
389 |
> |
//make unique vector of objects we need to cut on |
390 |
> |
sort( objectsToCut.begin(), objectsToCut.end() ); |
391 |
> |
objectsToCut.erase( unique( objectsToCut.begin(), objectsToCut.end() ), objectsToCut.end() ); |
392 |
|
|
393 |
|
|
394 |
|
} |
406 |
|
OSUAnalysis::analyze (const edm::Event &event, const edm::EventSetup &setup) |
407 |
|
{ |
408 |
|
|
284 |
– |
|
409 |
|
// Retrieve necessary collections from the event. |
410 |
< |
edm::Handle<BNtriggerCollection> triggers; |
411 |
< |
if (std::find(allNecessaryObjects.begin(), allNecessaryObjects.end(), "triggers") != allNecessaryObjects.end()) |
410 |
> |
|
411 |
> |
if (std::find(objectsToGet.begin(), objectsToGet.end(), "triggers") != objectsToGet.end()) |
412 |
|
event.getByLabel (triggers_, triggers); |
413 |
|
|
414 |
< |
edm::Handle<BNjetCollection> jets; |
291 |
< |
if (std::find(allNecessaryObjects.begin(), allNecessaryObjects.end(), "jets") != allNecessaryObjects.end()) |
414 |
> |
if (std::find(objectsToGet.begin(), objectsToGet.end(), "jets") != objectsToGet.end()) |
415 |
|
event.getByLabel (jets_, jets); |
416 |
|
|
417 |
< |
edm::Handle<BNmuonCollection> muons; |
295 |
< |
if (std::find(allNecessaryObjects.begin(), allNecessaryObjects.end(), "muons") != allNecessaryObjects.end()) |
417 |
> |
if (std::find(objectsToGet.begin(), objectsToGet.end(), "muons") != objectsToGet.end()) |
418 |
|
event.getByLabel (muons_, muons); |
419 |
|
|
420 |
< |
edm::Handle<BNelectronCollection> electrons; |
299 |
< |
if (std::find(allNecessaryObjects.begin(), allNecessaryObjects.end(), "electrons") != allNecessaryObjects.end()) |
420 |
> |
if (std::find(objectsToGet.begin(), objectsToGet.end(), "electrons") != objectsToGet.end()) |
421 |
|
event.getByLabel (electrons_, electrons); |
422 |
|
|
423 |
< |
edm::Handle<BNeventCollection> events; |
303 |
< |
if (std::find(allNecessaryObjects.begin(), allNecessaryObjects.end(), "events") != allNecessaryObjects.end()) |
423 |
> |
if (std::find(objectsToGet.begin(), objectsToGet.end(), "events") != objectsToGet.end()) |
424 |
|
event.getByLabel (events_, events); |
425 |
|
|
426 |
< |
edm::Handle<BNtauCollection> taus; |
307 |
< |
if (std::find(allNecessaryObjects.begin(), allNecessaryObjects.end(), "taus") != allNecessaryObjects.end()) |
426 |
> |
if (std::find(objectsToGet.begin(), objectsToGet.end(), "taus") != objectsToGet.end()) |
427 |
|
event.getByLabel (taus_, taus); |
428 |
|
|
429 |
< |
edm::Handle<BNmetCollection> mets; |
311 |
< |
if (std::find(allNecessaryObjects.begin(), allNecessaryObjects.end(), "mets") != allNecessaryObjects.end()) |
429 |
> |
if (std::find(objectsToGet.begin(), objectsToGet.end(), "mets") != objectsToGet.end()) |
430 |
|
event.getByLabel (mets_, mets); |
431 |
|
|
432 |
< |
edm::Handle<BNtrackCollection> tracks; |
315 |
< |
if (std::find(allNecessaryObjects.begin(), allNecessaryObjects.end(), "tracks") != allNecessaryObjects.end()) |
432 |
> |
if (std::find(objectsToGet.begin(), objectsToGet.end(), "tracks") != objectsToGet.end()) |
433 |
|
event.getByLabel (tracks_, tracks); |
434 |
|
|
435 |
< |
edm::Handle<BNgenjetCollection> genjets; |
319 |
< |
if (std::find(allNecessaryObjects.begin(), allNecessaryObjects.end(), "genjets") != allNecessaryObjects.end()) |
435 |
> |
if (std::find(objectsToGet.begin(), objectsToGet.end(), "genjets") != objectsToGet.end()) |
436 |
|
event.getByLabel (genjets_, genjets); |
437 |
|
|
438 |
< |
edm::Handle<BNmcparticleCollection> mcparticles; |
439 |
< |
if (std::find(allNecessaryObjects.begin(), allNecessaryObjects.end(), "mcparticles") != allNecessaryObjects.end()) |
440 |
< |
event.getByLabel (mcparticles_, mcparticles); |
438 |
> |
if (std::find(objectsToGet.begin(), objectsToGet.end(), "mcparticles") != objectsToGet.end()) |
439 |
> |
event.getByLabel (mcparticles_, mcparticles); |
440 |
> |
|
441 |
|
|
442 |
< |
edm::Handle<BNprimaryvertexCollection> primaryvertexs; |
327 |
< |
if (std::find(allNecessaryObjects.begin(), allNecessaryObjects.end(), "primaryvertexs") != allNecessaryObjects.end()) |
442 |
> |
if (std::find(objectsToGet.begin(), objectsToGet.end(), "primaryvertexs") != objectsToGet.end()) |
443 |
|
event.getByLabel (primaryvertexs_, primaryvertexs); |
444 |
|
|
445 |
< |
edm::Handle<BNbxlumiCollection> bxlumis; |
331 |
< |
if (std::find(allNecessaryObjects.begin(), allNecessaryObjects.end(), "bxlumis") != allNecessaryObjects.end()) |
445 |
> |
if (std::find(objectsToGet.begin(), objectsToGet.end(), "bxlumis") != objectsToGet.end()) |
446 |
|
event.getByLabel (bxlumis_, bxlumis); |
447 |
|
|
448 |
< |
edm::Handle<BNphotonCollection> photons; |
335 |
< |
if (std::find(allNecessaryObjects.begin(), allNecessaryObjects.end(), "photons") != allNecessaryObjects.end()) |
448 |
> |
if (std::find(objectsToGet.begin(), objectsToGet.end(), "photons") != objectsToGet.end()) |
449 |
|
event.getByLabel (photons_, photons); |
450 |
|
|
451 |
< |
edm::Handle<BNsuperclusterCollection> superclusters; |
339 |
< |
if (std::find(allNecessaryObjects.begin(), allNecessaryObjects.end(), "superclusters") != allNecessaryObjects.end()) |
451 |
> |
if (std::find(objectsToGet.begin(), objectsToGet.end(), "superclusters") != objectsToGet.end()) |
452 |
|
event.getByLabel (superclusters_, superclusters); |
453 |
|
|
454 |
|
|
455 |
|
//get pile-up event weight |
456 |
< |
double puScaleFactor = 0; |
457 |
< |
if(datasetType_ != "data") |
458 |
< |
puScaleFactor = puWeight_->at (events->at (0).numTruePV); |
459 |
< |
else |
460 |
< |
puScaleFactor = 1.00; |
456 |
> |
double puScaleFactor = 1.00; |
457 |
> |
// if(datasetType_ != "data"){ |
458 |
> |
// puScaleFactor = puWeight_->at (events->at (0).numTruePV); |
459 |
> |
// cout << puScaleFactor << endl; |
460 |
> |
//} |
461 |
> |
|
462 |
|
|
463 |
|
|
464 |
|
//loop over all channels |
477 |
|
} |
478 |
|
|
479 |
|
//loop over all cuts |
480 |
+ |
|
481 |
+ |
|
482 |
+ |
|
483 |
|
for(uint currentCutIndex = 0; currentCutIndex != currentChannel.cuts.size(); currentCutIndex++){ |
484 |
|
cut currentCut = currentChannel.cuts.at(currentCutIndex); |
485 |
|
|
486 |
< |
for(uint currentObjectIndex = 0; currentObjectIndex != allNecessaryObjects.size(); currentObjectIndex++){ |
487 |
< |
string currentObject = allNecessaryObjects.at(currentObjectIndex); |
486 |
> |
for(uint currentObjectIndex = 0; currentObjectIndex != objectsToCut.size(); currentObjectIndex++){ |
487 |
> |
string currentObject = objectsToCut.at(currentObjectIndex); |
488 |
|
|
489 |
+ |
//initialize maps to get ready to set cuts |
490 |
|
individualFlags[currentObject].push_back (vector<bool> ()); |
491 |
|
cumulativeFlags[currentObject].push_back (vector<bool> ()); |
492 |
|
|
493 |
+ |
} |
494 |
+ |
for(uint currentObjectIndex = 0; currentObjectIndex != objectsToCut.size(); currentObjectIndex++){ |
495 |
+ |
string currentObject = objectsToCut.at(currentObjectIndex); |
496 |
+ |
|
497 |
+ |
int flagsForPairCutsIndex = max(int(currentCutIndex-1),0); |
498 |
|
|
499 |
|
if(currentObject == "jets") setObjectFlags(currentCut,currentCutIndex,individualFlags,cumulativeFlags,jets.product(),"jets"); |
500 |
|
else if(currentObject == "muons") setObjectFlags(currentCut,currentCutIndex,individualFlags,cumulativeFlags,muons.product(),"muons"); |
510 |
|
else if(currentObject == "photons") setObjectFlags(currentCut,currentCutIndex,individualFlags,cumulativeFlags,photons.product(),"photons"); |
511 |
|
else if(currentObject == "superclusters") setObjectFlags(currentCut,currentCutIndex,individualFlags,cumulativeFlags,superclusters.product(),"superclusters"); |
512 |
|
|
513 |
< |
else if(currentObject == "muon-muon pairs") setObjectFlags(currentCut,currentCutIndex,individualFlags,cumulativeFlags,muons.product(),muons.product(),"muon-muon pairs"); |
514 |
< |
else if(currentObject == "electron-electron pairs") setObjectFlags(currentCut,currentCutIndex,individualFlags,cumulativeFlags,electrons.product(),electrons.product(),"electron-electron pairs"); |
515 |
< |
else if(currentObject == "electron-muon pairs") setObjectFlags(currentCut,currentCutIndex,individualFlags,cumulativeFlags,electrons.product(),muons.product(),"electron-muon pairs"); |
513 |
> |
|
514 |
> |
else if(currentObject == "muon-muon pairs") setObjectFlags(currentCut,currentCutIndex,individualFlags,cumulativeFlags,muons.product(),muons.product(), \ |
515 |
> |
cumulativeFlags.at("muons").at(flagsForPairCutsIndex), \ |
516 |
> |
cumulativeFlags.at("muons").at(flagsForPairCutsIndex), \ |
517 |
> |
"muon-muon pairs"); |
518 |
> |
else if(currentObject == "electron-electron pairs") setObjectFlags(currentCut,currentCutIndex,individualFlags,cumulativeFlags,electrons.product(),electrons.product(), \ |
519 |
> |
cumulativeFlags.at("electrons").at(flagsForPairCutsIndex), \ |
520 |
> |
cumulativeFlags.at("electrons").at(flagsForPairCutsIndex), \ |
521 |
> |
"electron-electron pairs"); |
522 |
> |
else if(currentObject == "electron-muon pairs") setObjectFlags(currentCut,currentCutIndex,individualFlags,cumulativeFlags,electrons.product(),muons.product(), \ |
523 |
> |
cumulativeFlags.at("electrons").at(flagsForPairCutsIndex), \ |
524 |
> |
cumulativeFlags.at("muons").at(flagsForPairCutsIndex), \ |
525 |
> |
"electron-muon pairs"); |
526 |
> |
|
527 |
|
|
528 |
|
} |
529 |
|
|
532 |
|
}//end loop over all cuts |
533 |
|
|
534 |
|
|
402 |
– |
|
535 |
|
//use cumulative flags to apply cuts at event level |
536 |
|
|
537 |
|
bool eventPassedAllCuts = true; |
540 |
|
eventPassedAllCuts = eventPassedAllCuts && triggerDecision; |
541 |
|
|
542 |
|
|
543 |
+ |
|
544 |
|
for(uint currentCutIndex = 0; currentCutIndex != currentChannel.cuts.size(); currentCutIndex++){ |
545 |
|
|
546 |
|
//loop over all objects and count how many passed the cumulative selection up to this point |
565 |
|
if(!eventPassedAllCuts)continue; |
566 |
|
|
567 |
|
|
568 |
< |
|
568 |
> |
|
569 |
|
|
570 |
|
//set position of primary vertex in event, in order to calculate quantities relative to it |
571 |
< |
primaryVertex_ = 0; |
572 |
< |
vector<bool> vertexFlags = cumulativeFlags.at("primaryvertexs").back(); |
573 |
< |
for (uint vertexIndex = 0; vertexIndex != vertexFlags.size(); vertexIndex++){ |
574 |
< |
if(!vertexFlags.at(vertexIndex)) continue; |
575 |
< |
primaryVertex_ = new BNprimaryvertex (primaryvertexs->at (vertexIndex)); |
576 |
< |
break; |
571 |
> |
if(std::find(objectsToCut.begin(), objectsToCut.end(), "primaryvertexs") != objectsToCut.end()) { |
572 |
> |
vector<bool> vertexFlags = cumulativeFlags.at("primaryvertexs").back(); |
573 |
> |
for (uint vertexIndex = 0; vertexIndex != vertexFlags.size(); vertexIndex++){ |
574 |
> |
if(!vertexFlags.at(vertexIndex)) continue; |
575 |
> |
chosenPrimaryVertex = & primaryvertexs->at(vertexIndex); |
576 |
> |
break; |
577 |
> |
} |
578 |
> |
} |
579 |
> |
else { |
580 |
> |
chosenPrimaryVertex = & primaryvertexs->at(0); |
581 |
|
} |
582 |
|
|
583 |
|
|
589 |
|
if(currentHistogram.inputVariables.size() == 1){ |
590 |
|
TH1D* histo; |
591 |
|
histo = oneDHists_.at(currentChannelIndex).at(currentHistogram.name); |
592 |
+ |
|
593 |
+ |
|
594 |
+ |
|
595 |
|
if(currentHistogram.inputCollection == "jets") fill1DHistogram(histo,currentHistogram,jets.product(),cumulativeFlags.at("jets").back(),puScaleFactor); |
596 |
< |
else if(currentHistogram.inputCollection == "muons") fill1DHistogram(histo,currentHistogram,muons.product(),cumulativeFlags.at("muons").back(),puScaleFactor); |
597 |
< |
else if(currentHistogram.inputCollection == "muon-muon pairs") fill1DHistogram(histo,currentHistogram,muons.product(),muons.product(),\ |
598 |
< |
cumulativeFlags.at("muons").back(),cumulativeFlags.at("muons").back(),\ |
599 |
< |
cumulativeFlags.at("muon-muon pairs").back(),puScaleFactor); |
596 |
> |
else if(currentHistogram.inputCollection == "muons") fill1DHistogram(histo,currentHistogram,muons.product(),cumulativeFlags.at("muons").back(),puScaleFactor); |
597 |
> |
else if(currentHistogram.inputCollection == "muon-muon pairs") |
598 |
> |
{ |
599 |
> |
fill1DHistogram(histo,currentHistogram,muons.product(),muons.product(), \ |
600 |
> |
cumulativeFlags.at("muons").back(),cumulativeFlags.at("muons").back(), \ |
601 |
> |
|
602 |
> |
cumulativeFlags.at("muon-muon pairs").back(),puScaleFactor); |
603 |
> |
|
604 |
> |
} |
605 |
|
else if(currentHistogram.inputCollection == "electrons") fill1DHistogram(histo,currentHistogram,electrons.product(),cumulativeFlags.at("electrons").back(),puScaleFactor); |
606 |
|
else if(currentHistogram.inputCollection == "electron-electron pairs") fill1DHistogram(histo,currentHistogram,electrons.product(),electrons.product(),\ |
607 |
|
cumulativeFlags.at("electrons").back(),cumulativeFlags.at("electrons").back(),\ |
623 |
|
else if(currentHistogram.inputVariables.size() == 2){ |
624 |
|
TH2D* histo; |
625 |
|
histo = twoDHists_.at(currentChannelIndex).at(currentHistogram.name); |
626 |
+ |
|
627 |
+ |
|
628 |
+ |
|
629 |
|
if(currentHistogram.inputCollection == "jets") fill2DHistogram(histo,currentHistogram,jets.product(),cumulativeFlags.at("jets").back(),puScaleFactor); |
630 |
|
else if(currentHistogram.inputCollection == "muons") fill2DHistogram(histo,currentHistogram,muons.product(),cumulativeFlags.at("muons").back(),puScaleFactor); |
631 |
|
else if(currentHistogram.inputCollection == "muon-muon pairs") fill2DHistogram(histo,currentHistogram,muons.product(),muons.product(), \ |
651 |
|
} |
652 |
|
} |
653 |
|
|
506 |
– |
|
654 |
|
//fills histograms with the sizes of collections |
655 |
|
for (uint currentObjectIndex = 0; currentObjectIndex != objectsToPlot.size(); currentObjectIndex++){ |
656 |
|
string currentObject = objectsToPlot.at(currentObjectIndex); |
657 |
|
|
658 |
< |
if(currentObject == "muon-muon pairs") currentObject = "dimuonPairs"; |
659 |
< |
else if(currentObject == "electron-electron pairs") currentObject = "dielectronPairs"; |
660 |
< |
else if(currentObject == "electron-muon pairs") currentObject = "electronMuonPairs"; |
661 |
< |
string tempCurrentObject = currentObject; |
658 |
> |
string objectToPlot = ""; |
659 |
> |
|
660 |
> |
if(currentObject == "muon-muon pairs") objectToPlot = "dimuonPairs"; |
661 |
> |
else if(currentObject == "electron-electron pairs") objectToPlot = "dielectronPairs"; |
662 |
> |
else if(currentObject == "electron-muon pairs") objectToPlot = "electronMuonPairs"; |
663 |
> |
else objectToPlot = currentObject; |
664 |
> |
string tempCurrentObject = objectToPlot; |
665 |
|
tempCurrentObject.at(0) = toupper(tempCurrentObject.at(0)); |
666 |
|
string histoName = "num" + tempCurrentObject; |
667 |
|
|
668 |
|
|
669 |
< |
if(currentObject == "jets") oneDHists_.at(currentChannelIndex).at(histoName)->Fill(jets->size(),puScaleFactor); |
670 |
< |
else if(currentObject == "muons") oneDHists_.at(currentChannelIndex).at(histoName)->Fill(muons->size(),puScaleFactor); |
671 |
< |
else if(currentObject == "dimuonPairs") oneDHists_.at(currentChannelIndex).at(histoName)->Fill(muons->size()*(muons->size()-1)/2,puScaleFactor); |
672 |
< |
else if(currentObject == "electrons") oneDHists_.at(currentChannelIndex).at(histoName)->Fill(electrons->size(),puScaleFactor); |
673 |
< |
else if(currentObject == "dielectronPairs") oneDHists_.at(currentChannelIndex).at(histoName)->Fill(electrons->size()*(electrons->size()-1)/2,puScaleFactor); |
674 |
< |
else if(currentObject == "electronMuonPairs") oneDHists_.at(currentChannelIndex).at(histoName)->Fill(electrons->size()*muons->size(),puScaleFactor); |
675 |
< |
else if(currentObject == "events") oneDHists_.at(currentChannelIndex).at(histoName)->Fill(events->size(),puScaleFactor); |
676 |
< |
else if(currentObject == "taus") oneDHists_.at(currentChannelIndex).at(histoName)->Fill(taus->size(),puScaleFactor); |
677 |
< |
else if(currentObject == "mets") oneDHists_.at(currentChannelIndex).at(histoName)->Fill(mets->size(),puScaleFactor); |
678 |
< |
else if(currentObject == "tracks") oneDHists_.at(currentChannelIndex).at(histoName)->Fill(tracks->size(),puScaleFactor); |
679 |
< |
else if(currentObject == "genjets") oneDHists_.at(currentChannelIndex).at(histoName)->Fill(genjets->size(),puScaleFactor); |
680 |
< |
else if(currentObject == "mcparticles") oneDHists_.at(currentChannelIndex).at(histoName)->Fill(mcparticles->size(),puScaleFactor); |
681 |
< |
else if(currentObject == "primaryvertexs"){ |
669 |
> |
|
670 |
> |
|
671 |
> |
//set position of primary vertex in event, in order to calculate quantities relative to it |
672 |
> |
if(std::find(objectsToCut.begin(), objectsToCut.end(), currentObject) != objectsToCut.end()) { |
673 |
> |
vector<bool> lastCutFlags = cumulativeFlags.at(currentObject).back(); |
674 |
> |
int numToPlot = 0; |
675 |
> |
for (uint currentFlag = 0; currentFlag != lastCutFlags.size(); currentFlag++){ |
676 |
> |
if(lastCutFlags.at(currentFlag)) numToPlot++; |
677 |
> |
} |
678 |
> |
if(objectToPlot == "primaryvertexs"){ |
679 |
> |
oneDHists_.at(currentChannelIndex).at(histoName+"BeforePileupCorrection")->Fill(primaryvertexs->size()); |
680 |
> |
oneDHists_.at(currentChannelIndex).at(histoName+"AfterPileupCorrection")->Fill(primaryvertexs->size(),puScaleFactor); |
681 |
> |
} |
682 |
> |
else |
683 |
> |
oneDHists_.at(currentChannelIndex).at(histoName)->Fill(numToPlot,puScaleFactor); |
684 |
> |
} |
685 |
> |
else if(objectToPlot == "jets") oneDHists_.at(currentChannelIndex).at(histoName)->Fill(jets->size(),puScaleFactor); |
686 |
> |
else if(objectToPlot == "muons") oneDHists_.at(currentChannelIndex).at(histoName)->Fill(muons->size(),puScaleFactor); |
687 |
> |
else if(objectToPlot == "dimuonPairs") oneDHists_.at(currentChannelIndex).at(histoName)->Fill(muons->size()*(muons->size()-1)/2,puScaleFactor); |
688 |
> |
else if(objectToPlot == "electrons") oneDHists_.at(currentChannelIndex).at(histoName)->Fill(electrons->size(),puScaleFactor); |
689 |
> |
else if(objectToPlot == "dielectronPairs") oneDHists_.at(currentChannelIndex).at(histoName)->Fill(electrons->size()*(electrons->size()-1)/2,puScaleFactor); |
690 |
> |
else if(objectToPlot == "electronMuonPairs") oneDHists_.at(currentChannelIndex).at(histoName)->Fill(electrons->size()*muons->size(),puScaleFactor); |
691 |
> |
else if(objectToPlot == "events") oneDHists_.at(currentChannelIndex).at(histoName)->Fill(events->size(),puScaleFactor); |
692 |
> |
else if(objectToPlot == "taus") oneDHists_.at(currentChannelIndex).at(histoName)->Fill(taus->size(),puScaleFactor); |
693 |
> |
else if(objectToPlot == "mets") oneDHists_.at(currentChannelIndex).at(histoName)->Fill(mets->size(),puScaleFactor); |
694 |
> |
else if(objectToPlot == "tracks") oneDHists_.at(currentChannelIndex).at(histoName)->Fill(tracks->size(),puScaleFactor); |
695 |
> |
else if(objectToPlot == "genjets") oneDHists_.at(currentChannelIndex).at(histoName)->Fill(genjets->size(),puScaleFactor); |
696 |
> |
else if(objectToPlot == "mcparticles") oneDHists_.at(currentChannelIndex).at(histoName)->Fill(mcparticles->size(),puScaleFactor); |
697 |
> |
else if(objectToPlot == "bxlumis") oneDHists_.at(currentChannelIndex).at(histoName)->Fill(bxlumis->size(),puScaleFactor); |
698 |
> |
else if(objectToPlot == "photons") oneDHists_.at(currentChannelIndex).at(histoName)->Fill(photons->size(),puScaleFactor); |
699 |
> |
else if(objectToPlot == "superclusters") oneDHists_.at(currentChannelIndex).at(histoName)->Fill(superclusters->size(),puScaleFactor); |
700 |
> |
else if(objectToPlot == "primaryvertexs"){ |
701 |
|
oneDHists_.at(currentChannelIndex).at(histoName+"BeforePileupCorrection")->Fill(primaryvertexs->size()); |
702 |
|
oneDHists_.at(currentChannelIndex).at(histoName+"AfterPileupCorrection")->Fill(primaryvertexs->size(),puScaleFactor); |
703 |
|
} |
535 |
– |
else if(currentObject == "bxlumis") oneDHists_.at(currentChannelIndex).at(histoName)->Fill(bxlumis->size(),puScaleFactor); |
536 |
– |
else if(currentObject == "photons") oneDHists_.at(currentChannelIndex).at(histoName)->Fill(photons->size(),puScaleFactor); |
537 |
– |
else if(currentObject == "superclusters") oneDHists_.at(currentChannelIndex).at(histoName)->Fill(superclusters->size(),puScaleFactor); |
704 |
|
|
705 |
|
} |
706 |
+ |
|
707 |
+ |
|
708 |
+ |
|
709 |
+ |
|
710 |
|
} //end loop over channel |
711 |
|
|
712 |
|
masterCutFlow_->fillCutFlow(puScaleFactor); |
908 |
|
else if(variable == "ecalIso") value = object->ecalIso; |
909 |
|
else if(variable == "hcalIso") value = object->hcalIso; |
910 |
|
else if(variable == "caloIso") value = object->caloIso; |
911 |
< |
else if(variable == "trackIsoDR03") value = object->trackIsoDR03; |
911 |
> |
else if(variable == "trackIsDR03") value = object->trackIsoDR03; |
912 |
|
else if(variable == "ecalIsoDR03") value = object->ecalIsoDR03; |
913 |
|
else if(variable == "hcalIsoDR03") value = object->hcalIsoDR03; |
914 |
|
else if(variable == "caloIsoDR03") value = object->caloIsoDR03; |
1047 |
|
else if(variable == "time_ndof") value = object->time_ndof; |
1048 |
|
|
1049 |
|
//user-defined variables |
1050 |
< |
else if(variable == "correctedD0VertexErr") value = hypot (object->tkD0err, hypot (primaryVertex_->xError, primaryVertex_->yError)); |
1051 |
< |
else if(variable == "correctedD0VertexSig") value = object->correctedD0Vertex / hypot (object->tkD0err, hypot (primaryVertex_->xError, primaryVertex_->yError)); |
1052 |
< |
else if(variable == "detIso") value = (object->trackIso) / object->pt; |
1050 |
> |
else if(variable == "correctedD0VertexErr") value = hypot (object->tkD0err, hypot (chosenPrimaryVertex->xError, chosenPrimaryVertex->yError)); |
1051 |
> |
else if(variable == "correctedD0VertexSig") value = object->correctedD0Vertex / hypot (object->tkD0err, hypot (chosenPrimaryVertex->xError, chosenPrimaryVertex->yError)); |
1052 |
> |
else if(variable == "detIso") value = (object->trackIsoDR03) / object->pt; |
1053 |
|
else if(variable == "relPFdBetaIso") value = (object->pfIsoR04SumChargedHadronPt + max(0.0, object->pfIsoR04SumNeutralHadronEt + object->pfIsoR04SumPhotonEt - 0.5*object->pfIsoR04SumPUPt)) / object->pt; |
1054 |
|
else if(variable == "relPFrhoIso") value = ( object->chargedHadronIso + max(0.0, object->neutralHadronIso + object->photonIso - object->AEffDr03*object->rhoPrime) ) / object->pt; |
1055 |
+ |
|
1056 |
+ |
|
1057 |
+ |
|
1058 |
+ |
else if(variable == "correctedD0VertexInEBPlus"){ |
1059 |
+ |
if(fabs(object->eta) < 0.8 && object->eta > 0) value = object->correctedD0Vertex; |
1060 |
+ |
else value = -999; |
1061 |
+ |
} |
1062 |
+ |
else if(variable == "correctedD0VertexOutEBPlus"){ |
1063 |
+ |
if(fabs(object->eta) < 1.479 && fabs(object->eta) > 0.8 && object->eta > 0) value = object->correctedD0Vertex; |
1064 |
+ |
else value = -999; |
1065 |
+ |
} |
1066 |
+ |
else if(variable == "correctedD0VertexEEPlus"){ |
1067 |
+ |
if(fabs(object->eta) > 1.479 && object->eta > 0) value = object->correctedD0Vertex; |
1068 |
+ |
else value = -999; |
1069 |
+ |
} |
1070 |
+ |
|
1071 |
+ |
else if(variable == "correctedD0BeamspotInEBPlus"){ |
1072 |
+ |
if(fabs(object->eta) < 0.8 && object->eta > 0) value = object->correctedD0; |
1073 |
+ |
else value = -999; |
1074 |
+ |
} |
1075 |
+ |
else if(variable == "correctedD0BeamspotOutEBPlus"){ |
1076 |
+ |
if(fabs(object->eta) < 1.479 && fabs(object->eta) > 0.8 && object->eta > 0) value = object->correctedD0; |
1077 |
+ |
else value = -999; |
1078 |
+ |
} |
1079 |
+ |
else if(variable == "correctedD0BeamspotEEPlus"){ |
1080 |
+ |
if(fabs(object->eta) > 1.479 && object->eta > 0) value = object->correctedD0; |
1081 |
+ |
else value = -999; |
1082 |
+ |
} |
1083 |
+ |
|
1084 |
+ |
else if(variable == "correctedD0VertexInEBMinus"){ |
1085 |
+ |
if(fabs(object->eta) < 0.8 && object->eta < 0) value = object->correctedD0Vertex; |
1086 |
+ |
else value = -999; |
1087 |
+ |
} |
1088 |
+ |
else if(variable == "correctedD0VertexOutEBMinus"){ |
1089 |
+ |
if(fabs(object->eta) < 1.479 && fabs(object->eta) > 0.8 && object->eta < 0) value = object->correctedD0Vertex; |
1090 |
+ |
else value = -999; |
1091 |
+ |
} |
1092 |
+ |
else if(variable == "correctedD0VertexEEMinus"){ |
1093 |
+ |
if(fabs(object->eta) > 1.479 && object->eta < 0) value = object->correctedD0Vertex; |
1094 |
+ |
else value = -999; |
1095 |
+ |
} |
1096 |
+ |
|
1097 |
+ |
else if(variable == "correctedD0BeamspotInEBMinus"){ |
1098 |
+ |
if(fabs(object->eta) < 0.8 && object->eta < 0) value = object->correctedD0; |
1099 |
+ |
else value = -999; |
1100 |
+ |
} |
1101 |
+ |
else if(variable == "correctedD0BeamspotOutEBMinus"){ |
1102 |
+ |
if(fabs(object->eta) < 1.479 && fabs(object->eta) > 0.8 && object->eta < 0) value = object->correctedD0; |
1103 |
+ |
else value = -999; |
1104 |
+ |
} |
1105 |
+ |
else if(variable == "correctedD0BeamspotEEMinus"){ |
1106 |
+ |
if(fabs(object->eta) > 1.479 && object->eta < 0) value = object->correctedD0; |
1107 |
+ |
else value = -999; |
1108 |
+ |
} |
1109 |
+ |
|
1110 |
+ |
|
1111 |
+ |
else if(variable == "correctedD0VertexInEBPositiveCharge"){ |
1112 |
+ |
if(fabs(object->eta) < 0.8 && object->charge > 0) value = object->correctedD0Vertex; |
1113 |
+ |
else value = -999; |
1114 |
+ |
} |
1115 |
+ |
else if(variable == "correctedD0VertexOutEBPositiveCharge"){ |
1116 |
+ |
if(fabs(object->eta) < 1.479 && fabs(object->eta) > 0.8 && object->charge > 0) value = object->correctedD0Vertex; |
1117 |
+ |
else value = -999; |
1118 |
+ |
} |
1119 |
+ |
else if(variable == "correctedD0VertexEEPositiveCharge"){ |
1120 |
+ |
if(fabs(object->eta) > 1.479 && object->charge > 0) value = object->correctedD0Vertex; |
1121 |
+ |
else value = -999; |
1122 |
+ |
} |
1123 |
+ |
|
1124 |
+ |
else if(variable == "correctedD0BeamspotInEBPositiveCharge"){ |
1125 |
+ |
if(fabs(object->eta) < 0.8 && object->charge > 0) value = object->correctedD0; |
1126 |
+ |
else value = -999; |
1127 |
+ |
} |
1128 |
+ |
else if(variable == "correctedD0BeamspotOutEBPositiveCharge"){ |
1129 |
+ |
if(fabs(object->eta) < 1.479 && fabs(object->eta) > 0.8 && object->charge > 0) value = object->correctedD0; |
1130 |
+ |
else value = -999; |
1131 |
+ |
} |
1132 |
+ |
else if(variable == "correctedD0BeamspotEEPositiveCharge"){ |
1133 |
+ |
if(fabs(object->eta) > 1.479 && object->charge > 0) value = object->correctedD0; |
1134 |
+ |
else value = -999; |
1135 |
+ |
} |
1136 |
+ |
|
1137 |
+ |
else if(variable == "correctedD0VertexInEBNegativeCharge"){ |
1138 |
+ |
if(fabs(object->eta) < 0.8 && object->charge < 0) value = object->correctedD0Vertex; |
1139 |
+ |
else value = -999; |
1140 |
+ |
} |
1141 |
+ |
else if(variable == "correctedD0VertexOutEBNegativeCharge"){ |
1142 |
+ |
if(fabs(object->eta) < 1.479 && fabs(object->eta) > 0.8 && object->charge < 0) value = object->correctedD0Vertex; |
1143 |
+ |
else value = -999; |
1144 |
+ |
} |
1145 |
+ |
else if(variable == "correctedD0VertexEENegativeCharge"){ |
1146 |
+ |
if(fabs(object->eta) > 1.479 && object->charge < 0) value = object->correctedD0Vertex; |
1147 |
+ |
else value = -999; |
1148 |
+ |
} |
1149 |
+ |
|
1150 |
+ |
else if(variable == "correctedD0BeamspotInEBNegativeCharge"){ |
1151 |
+ |
if(fabs(object->eta) < 0.8 && object->charge < 0) value = object->correctedD0; |
1152 |
+ |
else value = -999; |
1153 |
+ |
} |
1154 |
+ |
else if(variable == "correctedD0BeamspotOutEBNegativeCharge"){ |
1155 |
+ |
if(fabs(object->eta) < 1.479 && fabs(object->eta) > 0.8 && object->charge < 0) value = object->correctedD0; |
1156 |
+ |
else value = -999; |
1157 |
+ |
} |
1158 |
+ |
else if(variable == "correctedD0BeamspotEENegativeCharge"){ |
1159 |
+ |
if(fabs(object->eta) > 1.479 && object->charge < 0) value = object->correctedD0; |
1160 |
+ |
else value = -999; |
1161 |
+ |
} |
1162 |
+ |
|
1163 |
+ |
|
1164 |
|
else if(variable == "tightID") { |
1165 |
|
value = object->isGlobalMuon > 0 \ |
1166 |
|
&& object->isPFMuon > 0 \ |
1169 |
|
&& object->numberOfMatchedStations > 1 \ |
1170 |
|
&& fabs(object->correctedD0Vertex) < 0.2 \ |
1171 |
|
&& fabs(object->correctedDZ) < 0.5 \ |
1172 |
< |
&& object->numberOfValidPixelHits > 0 \ |
1172 |
> |
&& object->numberOfValidPixelHits > 0 \ |
1173 |
|
&& object->numberOfLayersWithMeasurement > 5; |
895 |
– |
|
896 |
– |
|
1174 |
|
} |
1175 |
|
else if(variable == "tightIDdisplaced"){ |
1176 |
|
value = object->isGlobalMuon > 0 \ |
900 |
– |
&& object->isPFMuon > 0 \ |
1177 |
|
&& object->normalizedChi2 < 10 \ |
1178 |
|
&& object->numberOfValidMuonHits > 0 \ |
1179 |
|
&& object->numberOfMatchedStations > 1 \ |
1181 |
|
&& object->numberOfLayersWithMeasurement > 5; |
1182 |
|
} |
1183 |
|
|
1184 |
+ |
|
1185 |
+ |
|
1186 |
+ |
else if(variable == "genMatchedId"){ |
1187 |
+ |
int index = getGenMatchedParticleIndex(object); |
1188 |
+ |
if(index == -1) value = 0; |
1189 |
+ |
else value = getPdgIdBinValue(mcparticles->at(index).id); |
1190 |
+ |
} |
1191 |
+ |
else if(variable == "genMatchedMotherId"){ |
1192 |
+ |
int index = getGenMatchedParticleIndex(object); |
1193 |
+ |
if(index == -1) value = 0; |
1194 |
+ |
else value = getPdgIdBinValue(mcparticles->at(index).motherId); |
1195 |
+ |
} |
1196 |
+ |
else if(variable == "genMatchedGrandmotherId"){ |
1197 |
+ |
int index = getGenMatchedParticleIndex(object); |
1198 |
+ |
if(index == -1) value = 0; |
1199 |
+ |
else if(fabs(mcparticles->at(index).motherId) == 15){ |
1200 |
+ |
int motherIndex = findTauMotherIndex(&mcparticles->at(index)); |
1201 |
+ |
if(motherIndex == -1) value = 0; |
1202 |
+ |
else value = getPdgIdBinValue(mcparticles->at(motherIndex).motherId); |
1203 |
+ |
} |
1204 |
+ |
else value = getPdgIdBinValue(mcparticles->at(index).grandMotherId); |
1205 |
+ |
} |
1206 |
+ |
|
1207 |
+ |
|
1208 |
+ |
|
1209 |
|
else{std::cout << "WARNING: invalid variable '" << variable << "'\n"; value = -999;} |
1210 |
|
|
1211 |
|
value = applyFunction(function, value); |
1376 |
|
else if(variable == "passConvVeto") value = object->passConvVeto; |
1377 |
|
|
1378 |
|
//user-defined variables |
1379 |
< |
else if(variable == "correctedD0VertexErr") value = hypot (object->tkD0err, hypot (primaryVertex_->xError, primaryVertex_->yError)); |
1380 |
< |
else if(variable == "correctedD0VertexSig") value = object->correctedD0Vertex / hypot (object->tkD0err, hypot (primaryVertex_->xError, primaryVertex_->yError)); |
1379 |
> |
else if(variable == "correctedD0VertexErr") value = hypot (object->tkD0err, hypot (chosenPrimaryVertex->xError, chosenPrimaryVertex->yError)); |
1380 |
> |
else if(variable == "correctedD0VertexSig") value = object->correctedD0Vertex / hypot (object->tkD0err, hypot (chosenPrimaryVertex->xError, chosenPrimaryVertex->yError)); |
1381 |
|
else if(variable == "detIso") value = (object->trackIso) / object->pt; |
1382 |
|
else if(variable == "relPFrhoIso") value = ( object->chargedHadronIsoDR03 + max(0.0, object->neutralHadronIsoDR03 + object->photonIsoDR03 - object->AEffDr03*object->rhoPrime) ) / object->pt; |
1383 |
|
|
1384 |
|
|
1385 |
+ |
else if(variable == "correctedD0VertexEEPositiveChargeLowPt"){ |
1386 |
+ |
if(fabs(object->eta) > 1.479 && object->charge > 0 && object->pt > 45) value = object->correctedD0Vertex; |
1387 |
+ |
else value = -999; |
1388 |
+ |
} |
1389 |
+ |
else if(variable == "correctedD0VertexEEPositiveChargeHighPt"){ |
1390 |
+ |
if(fabs(object->eta) > 1.479 && object->charge > 0 && object->pt < 45) value = object->correctedD0Vertex; |
1391 |
+ |
else value = -999; |
1392 |
+ |
} |
1393 |
+ |
|
1394 |
+ |
|
1395 |
+ |
else if(variable == "correctedD0VertexInEBPlus"){ |
1396 |
+ |
if(fabs(object->eta) < 0.8 && object->eta > 0) value = object->correctedD0Vertex; |
1397 |
+ |
else value = -999; |
1398 |
+ |
} |
1399 |
+ |
else if(variable == "correctedD0VertexOutEBPlus"){ |
1400 |
+ |
if(object->isEB && fabs(object->eta) > 0.8 && object->eta > 0) value = object->correctedD0Vertex; |
1401 |
+ |
else value = -999; |
1402 |
+ |
} |
1403 |
+ |
else if(variable == "correctedD0VertexEEPlus"){ |
1404 |
+ |
if(object->isEE && object->eta > 0) value = object->correctedD0Vertex; |
1405 |
+ |
else value = -999; |
1406 |
+ |
} |
1407 |
+ |
|
1408 |
+ |
else if(variable == "correctedD0BeamspotInEBPlus"){ |
1409 |
+ |
if(fabs(object->eta) < 0.8 && object->eta > 0) value = object->correctedD0; |
1410 |
+ |
else value = -999; |
1411 |
+ |
} |
1412 |
+ |
else if(variable == "correctedD0BeamspotOutEBPlus"){ |
1413 |
+ |
if(object->isEB && fabs(object->eta) > 0.8 && object->eta > 0) value = object->correctedD0; |
1414 |
+ |
else value = -999; |
1415 |
+ |
} |
1416 |
+ |
else if(variable == "correctedD0BeamspotEEPlus"){ |
1417 |
+ |
if(object->isEE && object->eta > 0) value = object->correctedD0; |
1418 |
+ |
else value = -999; |
1419 |
+ |
} |
1420 |
+ |
|
1421 |
+ |
else if(variable == "correctedD0VertexInEBMinus"){ |
1422 |
+ |
if(fabs(object->eta) < 0.8 && object->eta < 0) value = object->correctedD0Vertex; |
1423 |
+ |
else value = -999; |
1424 |
+ |
} |
1425 |
+ |
else if(variable == "correctedD0VertexOutEBMinus"){ |
1426 |
+ |
if(object->isEB && fabs(object->eta) > 0.8 && object->eta < 0) value = object->correctedD0Vertex; |
1427 |
+ |
else value = -999; |
1428 |
+ |
} |
1429 |
+ |
else if(variable == "correctedD0VertexEEMinus"){ |
1430 |
+ |
if(object->isEE && object->eta < 0) value = object->correctedD0Vertex; |
1431 |
+ |
else value = -999; |
1432 |
+ |
} |
1433 |
+ |
|
1434 |
+ |
else if(variable == "correctedD0BeamspotInEBMinus"){ |
1435 |
+ |
if(fabs(object->eta) < 0.8 && object->eta < 0) value = object->correctedD0; |
1436 |
+ |
else value = -999; |
1437 |
+ |
} |
1438 |
+ |
else if(variable == "correctedD0BeamspotOutEBMinus"){ |
1439 |
+ |
if(object->isEB && fabs(object->eta) > 0.8 && object->eta < 0) value = object->correctedD0; |
1440 |
+ |
else value = -999; |
1441 |
+ |
} |
1442 |
+ |
else if(variable == "correctedD0BeamspotEEMinus"){ |
1443 |
+ |
if(object->isEE && object->eta < 0) value = object->correctedD0; |
1444 |
+ |
else value = -999; |
1445 |
+ |
} |
1446 |
+ |
|
1447 |
+ |
|
1448 |
+ |
else if(variable == "correctedD0VertexInEBPositiveCharge"){ |
1449 |
+ |
if(fabs(object->eta) < 0.8 && object->charge > 0) value = object->correctedD0Vertex; |
1450 |
+ |
else value = -999; |
1451 |
+ |
} |
1452 |
+ |
else if(variable == "correctedD0VertexOutEBPositiveCharge"){ |
1453 |
+ |
if(object->isEB && fabs(object->eta) > 0.8 && object->charge > 0) value = object->correctedD0Vertex; |
1454 |
+ |
else value = -999; |
1455 |
+ |
} |
1456 |
+ |
else if(variable == "correctedD0VertexEEPositiveCharge"){ |
1457 |
+ |
if(object->isEE && object->charge > 0) value = object->correctedD0Vertex; |
1458 |
+ |
else value = -999; |
1459 |
+ |
} |
1460 |
+ |
|
1461 |
+ |
else if(variable == "correctedD0BeamspotInEBPositiveCharge"){ |
1462 |
+ |
if(fabs(object->eta) < 0.8 && object->charge > 0) value = object->correctedD0; |
1463 |
+ |
else value = -999; |
1464 |
+ |
} |
1465 |
+ |
else if(variable == "correctedD0BeamspotOutEBPositiveCharge"){ |
1466 |
+ |
if(object->isEB && fabs(object->eta) > 0.8 && object->charge > 0) value = object->correctedD0; |
1467 |
+ |
else value = -999; |
1468 |
+ |
} |
1469 |
+ |
else if(variable == "correctedD0BeamspotEEPositiveCharge"){ |
1470 |
+ |
if(object->isEE && object->charge > 0) value = object->correctedD0; |
1471 |
+ |
else value = -999; |
1472 |
+ |
} |
1473 |
+ |
|
1474 |
+ |
else if(variable == "correctedD0VertexInEBNegativeCharge"){ |
1475 |
+ |
if(fabs(object->eta) < 0.8 && object->charge < 0) value = object->correctedD0Vertex; |
1476 |
+ |
else value = -999; |
1477 |
+ |
} |
1478 |
+ |
else if(variable == "correctedD0VertexOutEBNegativeCharge"){ |
1479 |
+ |
if(object->isEB && fabs(object->eta) > 0.8 && object->charge < 0) value = object->correctedD0Vertex; |
1480 |
+ |
else value = -999; |
1481 |
+ |
} |
1482 |
+ |
else if(variable == "correctedD0VertexEENegativeCharge"){ |
1483 |
+ |
if(object->isEE && object->charge < 0) value = object->correctedD0Vertex; |
1484 |
+ |
else value = -999; |
1485 |
+ |
} |
1486 |
+ |
|
1487 |
+ |
else if(variable == "correctedD0BeamspotInEBNegativeCharge"){ |
1488 |
+ |
if(fabs(object->eta) < 0.8 && object->charge < 0) value = object->correctedD0; |
1489 |
+ |
else value = -999; |
1490 |
+ |
} |
1491 |
+ |
else if(variable == "correctedD0BeamspotOutEBNegativeCharge"){ |
1492 |
+ |
if(object->isEB && fabs(object->eta) > 0.8 && object->charge < 0) value = object->correctedD0; |
1493 |
+ |
else value = -999; |
1494 |
+ |
} |
1495 |
+ |
else if(variable == "correctedD0BeamspotEENegativeCharge"){ |
1496 |
+ |
if(object->isEE && object->charge < 0) value = object->correctedD0; |
1497 |
+ |
else value = -999; |
1498 |
+ |
} |
1499 |
+ |
|
1500 |
+ |
|
1501 |
+ |
else if(variable == "tightIDdisplaced"){ |
1502 |
+ |
if (object->isEB) |
1503 |
+ |
{ |
1504 |
+ |
value = fabs(object->delEtaIn) < 0.004 \ |
1505 |
+ |
&& fabs (object->delPhiIn) < 0.03 \ |
1506 |
+ |
&& object->scSigmaIEtaIEta < 0.01 \ |
1507 |
+ |
&& object->hadOverEm < 0.12 \ |
1508 |
+ |
&& object->absInvEMinusInvPin < 0.05; |
1509 |
+ |
} |
1510 |
+ |
else |
1511 |
+ |
{ |
1512 |
+ |
value = fabs (object->delEtaIn) < 0.005 \ |
1513 |
+ |
&& fabs (object->delPhiIn) < 0.02 \ |
1514 |
+ |
&& object->scSigmaIEtaIEta < 0.03 \ |
1515 |
+ |
&& object->hadOverEm < 0.10 \ |
1516 |
+ |
&& object->absInvEMinusInvPin < 0.05; |
1517 |
+ |
} |
1518 |
+ |
} |
1519 |
+ |
|
1520 |
+ |
else if(variable == "genMatchedId"){ |
1521 |
+ |
int index = getGenMatchedParticleIndex(object); |
1522 |
+ |
if(index == -1) value = 0; |
1523 |
+ |
else value = getPdgIdBinValue(mcparticles->at(index).id); |
1524 |
+ |
} |
1525 |
+ |
else if(variable == "genMatchedMotherId"){ |
1526 |
+ |
int index = getGenMatchedParticleIndex(object); |
1527 |
+ |
if(index == -1) value = 0; |
1528 |
+ |
else value = getPdgIdBinValue(mcparticles->at(index).motherId); |
1529 |
+ |
} |
1530 |
+ |
else if(variable == "genMatchedGrandmotherId"){ |
1531 |
+ |
int index = getGenMatchedParticleIndex(object); |
1532 |
+ |
if(index == -1) value = 0; |
1533 |
+ |
else if(fabs(mcparticles->at(index).motherId) == 15){ |
1534 |
+ |
int motherIndex = findTauMotherIndex(&mcparticles->at(index)); |
1535 |
+ |
if(motherIndex == -1) value = 0; |
1536 |
+ |
else value = getPdgIdBinValue(mcparticles->at(motherIndex).motherId); |
1537 |
+ |
} |
1538 |
+ |
else value = getPdgIdBinValue(mcparticles->at(index).grandMotherId); |
1539 |
+ |
} |
1540 |
+ |
|
1541 |
+ |
|
1542 |
+ |
|
1543 |
|
else{std::cout << "WARNING: invalid variable '" << variable << "'\n"; value = -999;} |
1544 |
|
|
1545 |
|
value = applyFunction(function, value); |
1671 |
|
else if(variable == "leadingTrackValid") value = object->leadingTrackValid; |
1672 |
|
|
1673 |
|
|
1674 |
+ |
else if(variable == "genMatchedId"){ |
1675 |
+ |
int index = getGenMatchedParticleIndex(object); |
1676 |
+ |
if(index == -1) value = 0; |
1677 |
+ |
else value = getPdgIdBinValue(mcparticles->at(index).id); |
1678 |
+ |
} |
1679 |
+ |
else if(variable == "genMatchedMotherId"){ |
1680 |
+ |
int index = getGenMatchedParticleIndex(object); |
1681 |
+ |
if(index == -1) value = 0; |
1682 |
+ |
else value = getPdgIdBinValue(mcparticles->at(index).motherId); |
1683 |
+ |
} |
1684 |
+ |
else if(variable == "genMatchedGrandmotherId"){ |
1685 |
+ |
int index = getGenMatchedParticleIndex(object); |
1686 |
+ |
if(index == -1) value = 0; |
1687 |
+ |
else if(fabs(mcparticles->at(index).motherId) == 15){ |
1688 |
+ |
int motherIndex = findTauMotherIndex(&mcparticles->at(index)); |
1689 |
+ |
if(motherIndex == -1) value = 0; |
1690 |
+ |
else value = getPdgIdBinValue(mcparticles->at(motherIndex).motherId); |
1691 |
+ |
} |
1692 |
+ |
else value = getPdgIdBinValue(mcparticles->at(index).grandMotherId); |
1693 |
+ |
} |
1694 |
+ |
|
1695 |
+ |
|
1696 |
|
else{std::cout << "WARNING: invalid variable '" << variable << "'\n"; value = -999;} |
1697 |
|
|
1698 |
|
value = applyFunction(function, value); |
1776 |
|
OSUAnalysis::valueLookup (const BNtrack* object, string variable, string function){ |
1777 |
|
|
1778 |
|
double value = 0.0; |
1779 |
< |
|
1779 |
> |
double pMag = sqrt(object->pt * object->pt + |
1780 |
> |
object->pz * object->pz); |
1781 |
> |
|
1782 |
|
if(variable == "pt") value = object->pt; |
1783 |
|
else if(variable == "px") value = object->px; |
1784 |
|
else if(variable == "py") value = object->py; |
1798 |
|
else if(variable == "isHighPurity") value = object->isHighPurity; |
1799 |
|
|
1800 |
|
|
1801 |
+ |
//additional BNs info for disappTrks |
1802 |
+ |
else if(variable == "isGoodPtResolution") value = object->isGoodPtResolution; |
1803 |
+ |
else if(variable == "caloEMDeltaRp3") value = object->caloEMDeltaRp3; |
1804 |
+ |
else if(variable == "caloHadDeltaRp3") value = object->caloHadDeltaRp3; |
1805 |
+ |
else if(variable == "caloEMDeltaRp4") value = object->caloEMDeltaRp4; |
1806 |
+ |
else if(variable == "caloHadDeltaRp4") value = object->caloHadDeltaRp4; |
1807 |
+ |
else if(variable == "caloEMDeltaRp5") value = object->caloEMDeltaRp5; |
1808 |
+ |
else if(variable == "caloHadDeltaRp5") value = object->caloHadDeltaRp5; |
1809 |
+ |
else if(variable == "nHitsMissingOuter") value = object->nHitsMissingOuter; |
1810 |
+ |
else if(variable == "nHitsMissingInner") value = object->nHitsMissingInner; |
1811 |
+ |
else if(variable == "nHitsMissingMiddle") value = object->nHitsMissingMiddle; |
1812 |
+ |
//user defined variables |
1813 |
+ |
else if(variable == "d0wrtBS") value = (object->vx-events->at(0).BSx)*object->py/object->pt - (object->vy-events->at(0).BSy)*object->px/object->pt; |
1814 |
+ |
else if(variable == "dZwrtBS") value = object->dZ - events->at(0).BSz; |
1815 |
+ |
else if(variable == "caloTotDeltaRp5") value =(object->caloHadDeltaRp5 + object->caloEMDeltaRp5); |
1816 |
+ |
else if(variable == "caloTotDeltaRp5ByP") value =( (object->caloHadDeltaRp5 + object->caloEMDeltaRp5)/pMag); |
1817 |
+ |
else if(variable == "isIso") value = getTrkIsIso(object, tracks.product()); |
1818 |
+ |
else if(variable == "isMatchedDeadEcal") value = getTrkIsMatchedDeadEcal(object); |
1819 |
+ |
else if(variable == "ptErrorByPt") value = (object->ptError/object->pt); |
1820 |
+ |
else if(variable == "ptError") value = object->ptError; |
1821 |
+ |
else if(variable == "ptRes") value = getTrkPtRes(object); |
1822 |
+ |
|
1823 |
+ |
|
1824 |
+ |
else if(variable == "genMatchedId"){ |
1825 |
+ |
int index = getGenMatchedParticleIndex(object); |
1826 |
+ |
if(index == -1) value = 0; |
1827 |
+ |
else value = getPdgIdBinValue(mcparticles->at(index).id); |
1828 |
+ |
} |
1829 |
+ |
else if(variable == "genMatchedMotherId"){ |
1830 |
+ |
int index = getGenMatchedParticleIndex(object); |
1831 |
+ |
if(index == -1) value = 0; |
1832 |
+ |
else value = getPdgIdBinValue(mcparticles->at(index).motherId); |
1833 |
+ |
} |
1834 |
+ |
else if(variable == "genMatchedGrandmotherId"){ |
1835 |
+ |
int index = getGenMatchedParticleIndex(object); |
1836 |
+ |
if(index == -1) value = 0; |
1837 |
+ |
else if(fabs(mcparticles->at(index).motherId) == 15){ |
1838 |
+ |
int motherIndex = findTauMotherIndex(&mcparticles->at(index)); |
1839 |
+ |
if(motherIndex == -1) value = 0; |
1840 |
+ |
else value = getPdgIdBinValue(mcparticles->at(motherIndex).motherId); |
1841 |
+ |
} |
1842 |
+ |
else value = getPdgIdBinValue(mcparticles->at(index).grandMotherId); |
1843 |
+ |
} |
1844 |
+ |
|
1845 |
+ |
|
1846 |
+ |
|
1847 |
|
else{std::cout << "WARNING: invalid variable '" << variable << "'\n"; value = -999;} |
1848 |
|
|
1849 |
|
value = applyFunction(function, value); |
1971 |
|
|
1972 |
|
//user-defined variables |
1973 |
|
else if (variable == "d0"){ |
1974 |
< |
double vx = object->vx - primaryVertex_->x, |
1975 |
< |
vy = object->vy - primaryVertex_->y, |
1974 |
> |
double vx = object->vx - chosenPrimaryVertex->x, |
1975 |
> |
vy = object->vy - chosenPrimaryVertex->y, |
1976 |
|
px = object->px, |
1977 |
|
py = object->py, |
1978 |
|
pt = object->pt; |
1979 |
|
value = (-vx * py + vy * px) / pt; |
1980 |
|
} |
1981 |
|
else if (variable == "dz"){ |
1982 |
< |
double vx = object->vx - primaryVertex_->x, |
1983 |
< |
vy = object->vy - primaryVertex_->y, |
1984 |
< |
vz = object->vz - primaryVertex_->z, |
1982 |
> |
double vx = object->vx - chosenPrimaryVertex->x, |
1983 |
> |
vy = object->vy - chosenPrimaryVertex->y, |
1984 |
> |
vz = object->vz - chosenPrimaryVertex->z, |
1985 |
|
px = object->px, |
1986 |
|
py = object->py, |
1987 |
|
pz = object->pz, |
1988 |
|
pt = object->pt; |
1989 |
|
value = vz - (vx * px + vy * py)/pt * (pz/pt); |
1990 |
|
} |
1991 |
< |
|
1991 |
> |
else if(variable == "v0"){ |
1992 |
> |
value = sqrt(object->vx*object->vx + object->vy*object->vy); |
1993 |
> |
} |
1994 |
> |
else if(variable == "deltaV0"){ |
1995 |
> |
value = sqrt((object->vx-chosenPrimaryVertex->x)*(object->vx-chosenPrimaryVertex->x) + (object->vy-chosenPrimaryVertex->y)*(object->vy-chosenPrimaryVertex->y)); |
1996 |
> |
} |
1997 |
> |
else if (variable == "deltaVx"){ |
1998 |
> |
value = object->vx - chosenPrimaryVertex->x; |
1999 |
> |
} |
2000 |
> |
else if (variable == "deltaVy"){ |
2001 |
> |
value = object->vy - chosenPrimaryVertex->y; |
2002 |
> |
} |
2003 |
> |
else if (variable == "deltaVz"){ |
2004 |
> |
value = object->vz - chosenPrimaryVertex->z; |
2005 |
> |
} |
2006 |
|
|
2007 |
|
|
2008 |
|
else{std::cout << "WARNING: invalid variable '" << variable << "'\n"; value = -999;} |
2132 |
|
else if(variable == "seedRecoFlag") value = object->seedRecoFlag; |
2133 |
|
|
2134 |
|
|
2135 |
+ |
else if(variable == "genMatchedId"){ |
2136 |
+ |
int index = getGenMatchedParticleIndex(object); |
2137 |
+ |
if(index == -1) value = 0; |
2138 |
+ |
else value = getPdgIdBinValue(mcparticles->at(index).id); |
2139 |
+ |
} |
2140 |
+ |
else if(variable == "genMatchedMotherId"){ |
2141 |
+ |
int index = getGenMatchedParticleIndex(object); |
2142 |
+ |
if(index == -1) value = 0; |
2143 |
+ |
else value = getPdgIdBinValue(mcparticles->at(index).motherId); |
2144 |
+ |
} |
2145 |
+ |
else if(variable == "genMatchedGrandmotherId"){ |
2146 |
+ |
int index = getGenMatchedParticleIndex(object); |
2147 |
+ |
if(index == -1) value = 0; |
2148 |
+ |
else if(fabs(mcparticles->at(index).motherId) == 15){ |
2149 |
+ |
int motherIndex = findTauMotherIndex(&mcparticles->at(index)); |
2150 |
+ |
if(motherIndex == -1) value = 0; |
2151 |
+ |
else value = getPdgIdBinValue(mcparticles->at(motherIndex).motherId); |
2152 |
+ |
} |
2153 |
+ |
else value = getPdgIdBinValue(mcparticles->at(index).grandMotherId); |
2154 |
+ |
} |
2155 |
+ |
|
2156 |
+ |
|
2157 |
|
else{std::cout << "WARNING: invalid variable '" << variable << "'\n"; value = -999;} |
2158 |
|
|
2159 |
|
value = applyFunction(function, value); |
2189 |
|
|
2190 |
|
double value = 0.0; |
2191 |
|
|
2192 |
< |
if(variable == "deltaPhi") value = deltaPhi(object1->phi,object2->phi); |
2192 |
> |
if(variable == "deltaPhi") value = fabs(deltaPhi(object1->phi,object2->phi)); |
2193 |
> |
else if(variable == "deltaR") value = deltaR(object1->eta,object1->phi,object2->eta,object2->phi); |
2194 |
|
else if(variable == "invMass"){ |
2195 |
|
TLorentzVector fourVector1(object1->px, object1->py, object1->pz, object1->energy); |
2196 |
|
TLorentzVector fourVector2(object2->px, object2->py, object2->pz, object2->energy); |
1631 |
– |
value = (fourVector1 + fourVector2).M(); |
1632 |
– |
} |
2197 |
|
|
2198 |
+ |
value = (fourVector1 + fourVector2).M();} |
2199 |
+ |
|
2200 |
+ |
else if(variable == "threeDAngle") |
2201 |
+ |
{ |
2202 |
+ |
TVector3 threeVector1(object1->px, object1->py, object1->pz); |
2203 |
+ |
TVector3 threeVector2(object2->px, object2->py, object2->pz); |
2204 |
+ |
value = (threeVector1.Angle(threeVector2)); |
2205 |
+ |
} |
2206 |
+ |
else if(variable == "alpha") |
2207 |
+ |
{ |
2208 |
+ |
static const double pi = 3.1415926535897932384626433832795028841971693993751058; |
2209 |
+ |
TVector3 threeVector1(object1->px, object1->py, object1->pz); |
2210 |
+ |
TVector3 threeVector2(object2->px, object2->py, object2->pz); |
2211 |
+ |
value = (pi-threeVector1.Angle( threeVector2)); |
2212 |
+ |
} |
2213 |
+ |
|
2214 |
+ |
else if(variable == "deltaR") value = deltaR(object1->eta,object1->phi,object2->eta,object2->phi); |
2215 |
+ |
|
2216 |
+ |
|
2217 |
+ |
else if(variable == "deltaCorrectedD0Vertex") value = object1->correctedD0Vertex - object2->correctedD0Vertex; |
2218 |
+ |
else if(variable == "deltaAbsCorrectedD0Vertex") value = fabs(object1->correctedD0Vertex) - fabs(object2->correctedD0Vertex); |
2219 |
+ |
else if(variable == "d0Sign"){ |
2220 |
+ |
double d0Sign = (object1->correctedD0Vertex*object2->correctedD0Vertex)/fabs(object1->correctedD0Vertex*object2->correctedD0Vertex); |
2221 |
+ |
if(d0Sign < 0) value = -0.5; |
2222 |
+ |
else if (d0Sign > 0) value = 0.5; |
2223 |
+ |
else value = -999; |
2224 |
+ |
} |
2225 |
+ |
else if(variable == "muon1CorrectedD0Vertex"){ |
2226 |
+ |
value = object1->correctedD0Vertex; |
2227 |
+ |
} |
2228 |
+ |
else if(variable == "muon2CorrectedD0Vertex"){ |
2229 |
+ |
value = object2->correctedD0Vertex; |
2230 |
+ |
} |
2231 |
+ |
else if(variable == "muon1timeAtIpInOut"){ |
2232 |
+ |
value = object1->timeAtIpInOut; |
2233 |
+ |
} |
2234 |
+ |
else if(variable == "muon2timeAtIpInOut"){ |
2235 |
+ |
value = object2->timeAtIpInOut; |
2236 |
+ |
} |
2237 |
|
else{std::cout << "WARNING: invalid variable '" << variable << "'\n"; value = -999;} |
2238 |
|
|
2239 |
+ |
|
2240 |
+ |
|
2241 |
+ |
|
2242 |
+ |
|
2243 |
+ |
|
2244 |
+ |
|
2245 |
+ |
|
2246 |
|
value = applyFunction(function, value); |
2247 |
|
|
2248 |
|
return value; |
2253 |
|
|
2254 |
|
double value = 0.0; |
2255 |
|
|
2256 |
< |
if(variable == "deltaPhi") value = deltaPhi(object1->phi,object2->phi); |
2256 |
> |
if(variable == "deltaPhi") value = fabs(deltaPhi(object1->phi,object2->phi)); |
2257 |
> |
else if(variable == "deltaR") value = deltaR(object1->eta,object1->phi,object2->eta,object2->phi); |
2258 |
|
else if(variable == "invMass"){ |
2259 |
|
TLorentzVector fourVector1(object1->px, object1->py, object1->pz, object1->energy); |
2260 |
|
TLorentzVector fourVector2(object2->px, object2->py, object2->pz, object2->energy); |
2261 |
< |
value = (fourVector1 + fourVector2).M(); |
2261 |
> |
|
2262 |
> |
value = (fourVector1 + fourVector2).M();} |
2263 |
> |
else if(variable == "threeDAngle") |
2264 |
> |
{ |
2265 |
> |
TVector3 threeVector1(object1->px, object1->py, object1->pz); |
2266 |
> |
TVector3 threeVector2(object2->px, object2->py, object2->pz); |
2267 |
> |
value = (threeVector1.Angle(threeVector2)); |
2268 |
> |
} |
2269 |
> |
|
2270 |
> |
|
2271 |
> |
|
2272 |
> |
|
2273 |
> |
else if(variable == "deltaCorrectedD0Vertex") value = object1->correctedD0Vertex - object2->correctedD0Vertex; |
2274 |
> |
else if(variable == "deltaAbsCorrectedD0Vertex") value = fabs(object1->correctedD0Vertex) - fabs(object2->correctedD0Vertex); |
2275 |
> |
|
2276 |
> |
else if(variable == "d0Sign") value = object1->correctedD0Vertex*object2->correctedD0Vertex/fabs(object1->correctedD0Vertex*object2->correctedD0Vertex); |
2277 |
> |
|
2278 |
> |
|
2279 |
> |
|
2280 |
> |
else if(variable == "d0Sign"){ |
2281 |
> |
double d0Sign = (object1->correctedD0Vertex*object2->correctedD0Vertex)/fabs(object1->correctedD0Vertex*object2->correctedD0Vertex); |
2282 |
> |
if(d0Sign < 0) value = -0.5; |
2283 |
> |
else if (d0Sign > 0) value = 0.5; |
2284 |
> |
else value = -999; |
2285 |
> |
} |
2286 |
> |
else if(variable == "electron1CorrectedD0Vertex"){ |
2287 |
> |
value = object1->correctedD0Vertex; |
2288 |
> |
} |
2289 |
> |
else if(variable == "electron2CorrectedD0Vertex"){ |
2290 |
> |
value = object2->correctedD0Vertex; |
2291 |
|
} |
2292 |
|
|
2293 |
|
else{std::cout << "WARNING: invalid variable '" << variable << "'\n"; value = -999;} |
2302 |
|
|
2303 |
|
double value = 0.0; |
2304 |
|
|
2305 |
< |
if(variable == "deltaPhi") value = deltaPhi(object1->phi,object2->phi); |
2305 |
> |
if(variable == "deltaPhi") value = fabs(deltaPhi(object1->phi,object2->phi)); |
2306 |
> |
else if(variable == "deltaR") value = deltaR(object1->eta,object1->phi,object2->eta,object2->phi); |
2307 |
|
else if(variable == "invMass"){ |
2308 |
|
TLorentzVector fourVector1(object1->px, object1->py, object1->pz, object1->energy); |
2309 |
|
TLorentzVector fourVector2(object2->px, object2->py, object2->pz, object2->energy); |
2310 |
< |
value = (fourVector1 + fourVector2).M(); |
2310 |
> |
|
2311 |
> |
value = (fourVector1 + fourVector2).M();} |
2312 |
> |
else if(variable == "threeDAngle") |
2313 |
> |
{ |
2314 |
> |
TVector3 threeVector1(object1->px, object1->py, object1->pz); |
2315 |
> |
TVector3 threeVector2(object2->px, object2->py, object2->pz); |
2316 |
> |
value = (threeVector1.Angle(threeVector2)); |
2317 |
> |
} |
2318 |
> |
|
2319 |
> |
|
2320 |
> |
else if(variable == "deltaCorrectedD0Vertex") value = object1->correctedD0Vertex - object2->correctedD0Vertex; |
2321 |
> |
else if(variable == "deltaAbsCorrectedD0Vertex") value = fabs(object1->correctedD0Vertex) - fabs(object2->correctedD0Vertex); |
2322 |
> |
else if(variable == "d0Sign"){ |
2323 |
> |
double d0Sign = (object1->correctedD0Vertex*object2->correctedD0Vertex)/fabs(object1->correctedD0Vertex*object2->correctedD0Vertex); |
2324 |
> |
if(d0Sign < 0) value = -0.5; |
2325 |
> |
else if (d0Sign > 0) value = 0.5; |
2326 |
> |
else value = -999; |
2327 |
|
} |
2328 |
+ |
else if(variable == "electronCorrectedD0Vertex"){ |
2329 |
+ |
value = object1->correctedD0Vertex; |
2330 |
+ |
} |
2331 |
+ |
else if(variable == "muonCorrectedD0Vertex"){ |
2332 |
+ |
value = object2->correctedD0Vertex; |
2333 |
+ |
} |
2334 |
+ |
else if(variable == "electronCorrectedD0"){ |
2335 |
+ |
value = object1->correctedD0; |
2336 |
+ |
} |
2337 |
+ |
else if(variable == "muonCorrectedD0"){ |
2338 |
+ |
value = object2->correctedD0; |
2339 |
+ |
} |
2340 |
+ |
|
2341 |
|
|
2342 |
|
else{std::cout << "WARNING: invalid variable '" << variable << "'\n"; value = -999;} |
2343 |
|
|
2347 |
|
} |
2348 |
|
|
2349 |
|
|
2350 |
+ |
// Calculate the number of tracks in cone of DeltaR<0.5 around track1. |
2351 |
+ |
// Return true iff no other tracks are found in this cone. |
2352 |
+ |
int |
2353 |
+ |
OSUAnalysis::getTrkIsIso (const BNtrack* track1, const BNtrackCollection* trackColl){ |
2354 |
+ |
for(BNtrackCollection::const_iterator track2 = trackColl->begin(); track2 !=trackColl->end(); track2++){ |
2355 |
+ |
if(track1->eta == track2->eta && track1->phi == track2->phi) continue; // Do not compare the track to itself. |
2356 |
+ |
double deltaRtrk = deltaR(track1->eta, track1->phi, track2->eta, track2->phi); |
2357 |
+ |
if(deltaRtrk < 0.5) return 0; |
2358 |
+ |
} |
2359 |
+ |
return 1; |
2360 |
+ |
|
2361 |
+ |
} |
2362 |
+ |
|
2363 |
+ |
|
2364 |
+ |
double |
2365 |
+ |
OSUAnalysis::getTrkPtRes (const BNtrack* track1){ |
2366 |
+ |
|
2367 |
+ |
double ptTrue = getTrkPtTrue(track1, mcparticles.product()); |
2368 |
+ |
double PtRes = (track1->pt - ptTrue) / ptTrue; |
2369 |
+ |
|
2370 |
+ |
return PtRes; |
2371 |
+ |
|
2372 |
+ |
} |
2373 |
+ |
|
2374 |
+ |
|
2375 |
+ |
double |
2376 |
+ |
OSUAnalysis::getTrkPtTrue (const BNtrack* track1, const BNmcparticleCollection* genPartColl){ |
2377 |
+ |
double value = -99; |
2378 |
+ |
double genDeltaRLowest = 999; |
2379 |
+ |
|
2380 |
+ |
for (BNmcparticleCollection::const_iterator genPart = genPartColl->begin(); genPart !=genPartColl->end(); genPart++){ |
2381 |
+ |
double genDeltaRtemp = deltaR(genPart->eta, genPart->phi,track1->eta, track1->phi); |
2382 |
+ |
if (genDeltaRtemp < genDeltaRLowest) { |
2383 |
+ |
genDeltaRLowest = genDeltaRtemp; |
2384 |
+ |
if (genDeltaRLowest < 0.05) { // Only consider it truth-matched if DeltaR<0.15. |
2385 |
+ |
double ptTrue = genPart->pt; |
2386 |
+ |
value = ptTrue; |
2387 |
+ |
} |
2388 |
+ |
} |
2389 |
+ |
} |
2390 |
+ |
|
2391 |
+ |
return value; |
2392 |
+ |
|
2393 |
+ |
} |
2394 |
+ |
|
2395 |
+ |
//creates a map of the dead Ecal channels in the barrel and endcap |
2396 |
+ |
//to see how the map of dead Ecal channels is created look at function getChannelStatusMaps() here: |
2397 |
+ |
//http://cmssw.cvs.cern.ch/cgi-bin/cmssw.cgi/UserCode/jbrinson/DisappTrk/OSUT3Analysis/AnaTools/src/OSUAnalysis.cc?revision=1.88&view=markup |
2398 |
+ |
void |
2399 |
+ |
OSUAnalysis::WriteDeadEcal (){ |
2400 |
+ |
double etaEcal, phiEcal; |
2401 |
+ |
ifstream DeadEcalFile(deadEcalFile_); |
2402 |
+ |
if(!DeadEcalFile) { |
2403 |
+ |
cout << "Error: DeadEcalFile has not been found." << endl; |
2404 |
+ |
return; |
2405 |
+ |
} |
2406 |
+ |
if(DeadEcalVec.size()!= 0){ |
2407 |
+ |
cout << "Error: DeadEcalVec has a nonzero size" << endl; |
2408 |
+ |
return; |
2409 |
+ |
} |
2410 |
+ |
while(!DeadEcalFile.eof()) |
2411 |
+ |
{ |
2412 |
+ |
DeadEcalFile >> etaEcal >> phiEcal; |
2413 |
+ |
DeadEcal newChan; |
2414 |
+ |
newChan.etaEcal = etaEcal; |
2415 |
+ |
newChan.phiEcal = phiEcal; |
2416 |
+ |
DeadEcalVec.push_back(newChan); |
2417 |
+ |
} |
2418 |
+ |
if(DeadEcalVec.size() == 0) cout << "Warning: No dead Ecal channels have been found." << endl; |
2419 |
+ |
} |
2420 |
+ |
|
2421 |
+ |
//if a track is found within dR<0.05 of a dead Ecal channel value = 1, otherwise value = 0 |
2422 |
+ |
int |
2423 |
+ |
OSUAnalysis::getTrkIsMatchedDeadEcal (const BNtrack* track1){ |
2424 |
+ |
double deltaRLowest = 999; |
2425 |
+ |
int value = 0; |
2426 |
+ |
if (DeadEcalVec.size() == 0) WriteDeadEcal(); |
2427 |
+ |
for(std::vector<DeadEcal>::const_iterator ecal = DeadEcalVec.begin(); ecal != DeadEcalVec.end(); ++ecal){ |
2428 |
+ |
double eta = ecal->etaEcal; |
2429 |
+ |
double phi = ecal->phiEcal; |
2430 |
+ |
double deltaRtemp = deltaR(eta, track1->eta, phi, track1->phi); |
2431 |
+ |
if(deltaRtemp < deltaRLowest) deltaRLowest = deltaRtemp; |
2432 |
+ |
} |
2433 |
+ |
if (deltaRLowest<0.05) {value = 1;} |
2434 |
+ |
else {value = 0;} |
2435 |
+ |
return value; |
2436 |
+ |
} |
2437 |
+ |
|
2438 |
+ |
|
2439 |
+ |
|
2440 |
+ |
|
2441 |
|
double |
2442 |
|
OSUAnalysis::applyFunction(string function, double value){ |
2443 |
|
|
2444 |
|
if(function == "abs") value = fabs(value); |
2445 |
+ |
else if(function == "fabs") value = fabs(value); |
2446 |
+ |
else if(function == "log10") value = log10(value); |
2447 |
+ |
else if(function == "log") value = log10(value); |
2448 |
|
|
2449 |
+ |
else if(function == "") value = value; |
2450 |
+ |
else{std::cout << "WARNING: invalid function '" << function << "'\n";} |
2451 |
|
|
2452 |
|
return value; |
2453 |
|
|
2457 |
|
template <class InputCollection> |
2458 |
|
void OSUAnalysis::setObjectFlags(cut ¤tCut, uint currentCutIndex, flagMap &individualFlags, flagMap &cumulativeFlags, InputCollection inputCollection, string inputType){ |
2459 |
|
|
2460 |
+ |
|
2461 |
|
for (uint object = 0; object != inputCollection->size(); object++){ |
2462 |
|
|
2463 |
+ |
|
2464 |
|
bool decision = true;//object passes if this cut doesn't cut on that type of object |
2465 |
|
|
2466 |
+ |
|
2467 |
|
if(currentCut.inputCollection == inputType){ |
2468 |
|
|
2469 |
|
vector<bool> subcutDecisions; |
2485 |
|
} |
2486 |
|
individualFlags.at(inputType).at(currentCutIndex).push_back(decision); |
2487 |
|
|
2488 |
+ |
|
2489 |
|
//set flags for objects that pass each cut AND all the previous cuts |
2490 |
|
bool previousCumulativeFlag = true; |
2491 |
|
for(uint previousCutIndex = 0; previousCutIndex != currentCutIndex; previousCutIndex++){ |
2494 |
|
} |
2495 |
|
cumulativeFlags.at(inputType).at(currentCutIndex).push_back(previousCumulativeFlag && decision); |
2496 |
|
|
2497 |
+ |
|
2498 |
|
} |
2499 |
|
|
2500 |
|
} |
2502 |
|
|
2503 |
|
template <class InputCollection1, class InputCollection2> |
2504 |
|
void OSUAnalysis::setObjectFlags(cut ¤tCut, uint currentCutIndex, flagMap &individualFlags, flagMap &cumulativeFlags, \ |
2505 |
< |
InputCollection1 inputCollection1, InputCollection2 inputCollection2, string inputType){ |
2505 |
> |
InputCollection1 inputCollection1, InputCollection2 inputCollection2, vector<bool> flags1, vector<bool> flags2, string inputType){ |
2506 |
|
|
2507 |
|
|
2508 |
|
bool sameObjects = false; |
2509 |
|
if(typeid(InputCollection1).name() == typeid(InputCollection2).name()) sameObjects = true; |
2510 |
|
|
2511 |
+ |
|
2512 |
|
int counter = 0; |
2513 |
|
for (uint object1 = 0; object1 != inputCollection1->size(); object1++){ |
2514 |
|
for (uint object2 = 0; object2 != inputCollection2->size(); object2++){ |
2546 |
|
if(previousCumulativeFlag && individualFlags.at(inputType).at(previousCutIndex).at(counter)) previousCumulativeFlag = true; |
2547 |
|
else{ previousCumulativeFlag = false; break;} |
2548 |
|
} |
2549 |
< |
cumulativeFlags.at(inputType).at(currentCutIndex).push_back(previousCumulativeFlag && decision); |
2549 |
> |
//apply flags for the components of the composite object as well |
2550 |
> |
bool currentCumulativeFlag = true; |
2551 |
> |
if(flags1.size() == 0 && flags2.size() == 0) currentCumulativeFlag = previousCumulativeFlag && decision; |
2552 |
> |
else if(flags1.size() == 0) currentCumulativeFlag = previousCumulativeFlag && decision && flags2.at(object2); |
2553 |
> |
else if(flags2.size() == 0) currentCumulativeFlag = previousCumulativeFlag && decision && flags1.at(object1); |
2554 |
> |
else currentCumulativeFlag = previousCumulativeFlag && decision && flags1.at(object1) && flags2.at(object2); |
2555 |
> |
cumulativeFlags.at(inputType).at(currentCutIndex).push_back(currentCumulativeFlag); |
2556 |
|
|
2557 |
|
counter++; |
2558 |
|
} |
2595 |
|
bool sameObjects = false; |
2596 |
|
if(typeid(InputCollection1).name() == typeid(InputCollection2).name()) sameObjects = true; |
2597 |
|
|
2598 |
< |
int counter = 0; |
2598 |
> |
int pairCounter = 0; |
2599 |
|
for (uint object1 = 0; object1 != inputCollection1->size(); object1++){ |
2600 |
|
for (uint object2 = 0; object2 != inputCollection2->size(); object2++){ |
2601 |
|
|
2604 |
|
//only take objects which have passed all cuts and pairs which have passed all cuts |
2605 |
|
if(!plotAllObjectsInPassingEvents_ && !flags1.at(object1)) continue; |
2606 |
|
if(!plotAllObjectsInPassingEvents_ && !flags2.at(object2)) continue; |
2607 |
< |
if(!plotAllObjectsInPassingEvents_ && !pairFlags.at(counter)) continue; |
2607 |
> |
if(!plotAllObjectsInPassingEvents_ && !pairFlags.at(pairCounter)) continue; |
2608 |
|
|
2609 |
|
string currentString = parameters.inputVariables.at(0); |
2610 |
|
string inputVariable = ""; |
2621 |
|
double value = valueLookup(&inputCollection1->at(object1), &inputCollection2->at(object2), inputVariable, function); |
2622 |
|
histo->Fill(value,puScaleFactor); |
2623 |
|
|
2624 |
< |
counter++; |
2624 |
> |
pairCounter++; |
2625 |
|
} |
2626 |
|
} |
2627 |
|
|
2674 |
|
bool sameObjects = false; |
2675 |
|
if(typeid(InputCollection1).name() == typeid(InputCollection2).name()) sameObjects = true; |
2676 |
|
|
2677 |
< |
int counter = 0; |
2677 |
> |
int pairCounter = 0; |
2678 |
|
for (uint object1 = 0; object1 != inputCollection1->size(); object1++){ |
2679 |
|
for (uint object2 = 0; object2 != inputCollection2->size(); object2++){ |
2680 |
|
|
2683 |
|
//only take objects which have passed all cuts and pairs which have passed all cuts |
2684 |
|
if(!plotAllObjectsInPassingEvents_ && !flags1.at(object1)) continue; |
2685 |
|
if(!plotAllObjectsInPassingEvents_ && !flags2.at(object2)) continue; |
2686 |
< |
if(!plotAllObjectsInPassingEvents_ && !pairFlags.at(counter)) continue; |
2686 |
> |
if(!plotAllObjectsInPassingEvents_ && !pairFlags.at(pairCounter)) continue; |
2687 |
|
|
2688 |
|
string currentString = parameters.inputVariables.at(0); |
2689 |
|
string inputVariable = ""; |
2714 |
|
|
2715 |
|
histo->Fill(valueX,valueY,puScaleFactor); |
2716 |
|
|
2717 |
< |
counter++; |
2717 |
> |
pairCounter++; |
2718 |
> |
|
2719 |
> |
} |
2720 |
> |
} |
2721 |
> |
|
2722 |
> |
} |
2723 |
> |
|
2724 |
> |
|
2725 |
> |
template <class InputObject> |
2726 |
> |
int OSUAnalysis::getGenMatchedParticleIndex(InputObject object){ |
2727 |
> |
|
2728 |
> |
int bestMatchIndex = -1; |
2729 |
> |
double bestMatchDeltaR = 999; |
2730 |
> |
|
2731 |
> |
for(BNmcparticleCollection::const_iterator mcparticle = mcparticles->begin (); mcparticle != mcparticles->end (); mcparticle++){ |
2732 |
> |
|
2733 |
> |
double currentDeltaR = deltaR(object->eta,object->phi,mcparticle->eta,mcparticle->phi); |
2734 |
> |
if(currentDeltaR > 0.05) continue; |
2735 |
> |
|
2736 |
> |
if(currentDeltaR < bestMatchDeltaR && mcparticles->at(mcparticle - mcparticles->begin()).id != mcparticles->at(mcparticle - mcparticles->begin()).motherId){ |
2737 |
> |
bestMatchIndex = mcparticle - mcparticles->begin(); |
2738 |
> |
bestMatchDeltaR = currentDeltaR; |
2739 |
|
} |
2740 |
+ |
|
2741 |
|
} |
2742 |
|
|
2743 |
+ |
return bestMatchIndex; |
2744 |
+ |
|
2745 |
|
} |
2746 |
|
|
2747 |
|
|
2748 |
+ |
int OSUAnalysis::findTauMotherIndex(const BNmcparticle* tau){ |
2749 |
+ |
|
2750 |
+ |
int bestMatchIndex = -1; |
2751 |
+ |
double bestMatchDeltaR = 999; |
2752 |
+ |
|
2753 |
+ |
for(BNmcparticleCollection::const_iterator mcparticle = mcparticles->begin (); mcparticle != mcparticles->end (); mcparticle++){ |
2754 |
+ |
|
2755 |
+ |
if(fabs(mcparticle->id) != 15 || mcparticle->status !=3) continue; |
2756 |
+ |
|
2757 |
+ |
double currentDeltaR = deltaR(tau->eta,tau->phi,mcparticle->eta,mcparticle->phi); |
2758 |
+ |
if(currentDeltaR > 0.05) continue; |
2759 |
+ |
|
2760 |
+ |
if(currentDeltaR < bestMatchDeltaR && mcparticles->at(mcparticle - mcparticles->begin()).id != mcparticles->at(mcparticle - mcparticles->begin()).motherId){ |
2761 |
+ |
bestMatchIndex = mcparticle - mcparticles->begin(); |
2762 |
+ |
bestMatchDeltaR = currentDeltaR; |
2763 |
+ |
} |
2764 |
+ |
|
2765 |
+ |
} |
2766 |
+ |
|
2767 |
+ |
return bestMatchIndex; |
2768 |
+ |
} |
2769 |
+ |
|
2770 |
+ |
|
2771 |
+ |
// bin particle type |
2772 |
+ |
// --- ------------- |
2773 |
+ |
// 0 unmatched |
2774 |
+ |
// 1 u |
2775 |
+ |
// 2 d |
2776 |
+ |
// 3 s |
2777 |
+ |
// 4 c |
2778 |
+ |
// 5 b |
2779 |
+ |
// 6 t |
2780 |
+ |
// 7 e |
2781 |
+ |
// 8 mu |
2782 |
+ |
// 9 tau |
2783 |
+ |
// 10 nu |
2784 |
+ |
// 11 g |
2785 |
+ |
// 12 gamma |
2786 |
+ |
// 13 Z |
2787 |
+ |
// 14 W |
2788 |
+ |
// 15 light meson |
2789 |
+ |
// 16 K meson |
2790 |
+ |
// 17 D meson |
2791 |
+ |
// 18 B meson |
2792 |
+ |
// 19 light baryon |
2793 |
+ |
// 20 strange baryon |
2794 |
+ |
// 21 charm baryon |
2795 |
+ |
// 22 bottom baryon |
2796 |
+ |
// 23 other |
2797 |
+ |
|
2798 |
+ |
int OSUAnalysis::getPdgIdBinValue(int pdgId){ |
2799 |
+ |
|
2800 |
+ |
int binValue = -999; |
2801 |
+ |
int absPdgId = fabs(pdgId); |
2802 |
+ |
if(pdgId == -1) binValue = 0; |
2803 |
+ |
else if(absPdgId <= 6 ) binValue = absPdgId; |
2804 |
+ |
else if(absPdgId == 11 ) binValue = 7; |
2805 |
+ |
else if(absPdgId == 13 ) binValue = 8; |
2806 |
+ |
else if(absPdgId == 15 ) binValue = 9; |
2807 |
+ |
else if(absPdgId == 12 || absPdgId == 14 || absPdgId == 16 ) binValue = 10; |
2808 |
+ |
else if(absPdgId == 21 ) binValue = 11; |
2809 |
+ |
else if(absPdgId == 22 ) binValue = 12; |
2810 |
+ |
else if(absPdgId == 23 ) binValue = 13; |
2811 |
+ |
else if(absPdgId == 24 ) binValue = 14; |
2812 |
+ |
else if(absPdgId > 100 && absPdgId < 300 && absPdgId != 130 ) binValue = 15; |
2813 |
+ |
else if( absPdgId == 130 || (absPdgId > 300 && absPdgId < 400) ) binValue = 16; |
2814 |
+ |
else if(absPdgId > 400 && absPdgId < 500 ) binValue = 17; |
2815 |
+ |
else if(absPdgId > 500 && absPdgId < 600 ) binValue = 18; |
2816 |
+ |
else if(absPdgId > 1000 && absPdgId < 3000 ) binValue = 19; |
2817 |
+ |
else if(absPdgId > 3000 && absPdgId < 4000 ) binValue = 20; |
2818 |
+ |
else if(absPdgId > 4000 && absPdgId < 5000 ) binValue = 21; |
2819 |
+ |
else if(absPdgId > 5000 && absPdgId < 6000 ) binValue = 22; |
2820 |
+ |
|
2821 |
+ |
else binValue = 23; |
2822 |
+ |
|
2823 |
+ |
return binValue; |
2824 |
+ |
|
2825 |
+ |
} |
2826 |
+ |
|
2827 |
|
|
2828 |
|
DEFINE_FWK_MODULE(OSUAnalysis); |
2829 |
|
|