| 1 |
#include <iostream>
|
| 2 |
#include <iomanip>
|
| 3 |
|
| 4 |
#include "../interface/ToplikeCandidate.h"
|
| 5 |
|
| 6 |
|
| 7 |
using namespace std;
|
| 8 |
using namespace BAT;
|
| 9 |
|
| 10 |
namespace BAT {
|
| 11 |
|
| 12 |
static const double pi = 3.1415927;
|
| 13 |
static const int DQ_ID = 1, BQ_ID = 5, TOP_ID = 6, W_ID = 24, ELECTRON_ID = 11;
|
| 14 |
static const int TAU_ID = 15, NU_E_ID = 12, MU_ID = 13, NU_MU_ID = 14;
|
| 15 |
static const int WPRIME_ID = 4000023;
|
| 16 |
static const unsigned int MCLIMIT = 90; // Number of MC particles to check
|
| 17 |
static const float DELTA_R_LIMIT = 0.5; // Biggest deltaR to consider in range
|
| 18 |
static const double pt_tprime_sys_mean = 48.49;
|
| 19 |
static const double pt_tprime_sys_sigma = 27.87;
|
| 20 |
|
| 21 |
static map <ToplikeSelectionSteps::Step, TTbarEPlusJetsSelection::Step> selSteps;
|
| 22 |
|
| 23 |
|
| 24 |
ToplikeCandidate::ToplikeCandidate() :
|
| 25 |
jet3FromWIndex(0),
|
| 26 |
jet4FromWIndex(0),
|
| 27 |
doneReconstructiontop(false)
|
| 28 |
{
|
| 29 |
}
|
| 30 |
|
| 31 |
|
| 32 |
ToplikeCandidate::ToplikeCandidate(const Event& event) :
|
| 33 |
TopPairEventCandidate(event),
|
| 34 |
jet3FromWIndex(0),
|
| 35 |
jet4FromWIndex(0),
|
| 36 |
doneReconstructiontop(false)
|
| 37 |
{
|
| 38 |
static bool mapset = false;
|
| 39 |
if (mapset == false) {
|
| 40 |
for(unsigned int ind = 0; ind < sizeof(selStepArr)/sizeof(selPair); ++ind) {
|
| 41 |
selSteps[selStepArr[ind].key] = selStepArr[ind].value;
|
| 42 |
}
|
| 43 |
mapset= true;
|
| 44 |
}
|
| 45 |
}
|
| 46 |
|
| 47 |
|
| 48 |
ToplikeCandidate::~ToplikeCandidate()
|
| 49 |
{
|
| 50 |
}
|
| 51 |
|
| 52 |
|
| 53 |
bool ToplikeCandidate::hasAtLeastFiveGoodJets() const {
|
| 54 |
return goodJets.size() >= 5;
|
| 55 |
}
|
| 56 |
|
| 57 |
/*
|
| 58 |
const ParticlePointer ToplikeCandidate::getResonance() const {
|
| 59 |
throwExpeptionIfNotReconstructed("ToplikeCandidate::getResonance");
|
| 60 |
return tPrime;
|
| 61 |
}
|
| 62 |
*/
|
| 63 |
|
| 64 |
|
| 65 |
double ToplikeCandidate::tpmass() const {
|
| 66 |
/*
|
| 67 |
ParticlePointer tmp = ParticlePointer(new Particle(*leptonicTop1 + *hadronicW));
|
| 68 |
cout << "tmp px = " << tmp->px() << " tp py = " << tmp->py() << " tp pz = " << tmp->pz() << endl;
|
| 69 |
cout << "energy = " << tmp->energy() << endl;
|
| 70 |
return tmp->mass();
|
| 71 |
cout << "tp px = " << tPrime->px() << " tp py = " << tPrime->py() << " tp pz = " << tPrime->pz() << endl;
|
| 72 |
cout << "energy = " << tPrime->energy() << endl;
|
| 73 |
if (tPrime->px() <= 0)
|
| 74 |
tPrime->setMass(0.0);
|
| 75 |
*/
|
| 76 |
return tPrime->mass();
|
| 77 |
}
|
| 78 |
|
| 79 |
|
| 80 |
double ToplikeCandidate::sumPt() const {
|
| 81 |
return (jet1FromW->pt() + jet2FromW->pt() + hadronicBJet->pt());
|
| 82 |
}
|
| 83 |
|
| 84 |
|
| 85 |
double ToplikeCandidate::getGlobalChi2() const {
|
| 86 |
double pttbar = PtTtbarSystem();
|
| 87 |
double pttbarTerm = TMath::Power(pttbar - pt_tprime_sys_mean, 2) / (2 * pt_tprime_sys_sigma * pt_tprime_sys_sigma);
|
| 88 |
|
| 89 |
double htSystemDifference = TMath::Power(HTSystem() - matched_HTSystem, 2);
|
| 90 |
double htSystemError = matched_HTSystem_sigma * matched_HTSystem_sigma * 2;
|
| 91 |
double htSystemTerm = htSystemDifference / htSystemError;
|
| 92 |
return 1 / sqrt(2) * (pttbarTerm + htSystemTerm);
|
| 93 |
}
|
| 94 |
|
| 95 |
|
| 96 |
double ToplikeCandidate::PtTtbarSystem() const {
|
| 97 |
ParticlePointer combined;
|
| 98 |
combined = ParticlePointer(new Particle(*hadronicWtPrime + *hadronicTop));
|
| 99 |
return combined->pt() / HTSystem();
|
| 100 |
}
|
| 101 |
|
| 102 |
|
| 103 |
double ToplikeCandidate::TPrimeHTSystem() const {
|
| 104 |
return htSystem;
|
| 105 |
}
|
| 106 |
|
| 107 |
|
| 108 |
double ToplikeCandidate::PtTPrimeSystem() const {
|
| 109 |
return ptTprimeSystem;
|
| 110 |
}
|
| 111 |
|
| 112 |
|
| 113 |
double ToplikeCandidate::PtTtbarSystem(unsigned short neutrinoSolution) const {
|
| 114 |
ParticlePointer combined;
|
| 115 |
if (neutrinoSolution == 1)
|
| 116 |
combined = ParticlePointer(new Particle(*leptonicTop1 + *hadronicW));
|
| 117 |
else
|
| 118 |
combined = ParticlePointer(new Particle(*leptonicTop2 + *hadronicW));
|
| 119 |
return combined->pt() / HTSystem();
|
| 120 |
}
|
| 121 |
|
| 122 |
|
| 123 |
double ToplikeCandidate::getWChi2() const {
|
| 124 |
double WmassDifference = TMath::Power(hadronicWtPrime->mass() - pdgWmass, 2);
|
| 125 |
double WmassError = 2 * wMassWidth * wMassWidth;
|
| 126 |
double WmassTerm = WmassDifference / WmassError;
|
| 127 |
return (WmassTerm);
|
| 128 |
}
|
| 129 |
|
| 130 |
|
| 131 |
double ToplikeCandidate::getLoneHadChi2() const
|
| 132 |
{
|
| 133 |
return (TopPairEventCandidate::getHadronicChi2());
|
| 134 |
}
|
| 135 |
|
| 136 |
|
| 137 |
/*
|
| 138 |
double ToplikeCandidate::getTotalChi2(unsigned short int neutrinoSolution) const {
|
| 139 |
return getLeptonicChi2(neutrinoSolution) + getHadronicChi2() + getGlobalChi2(neutrinoSolution);
|
| 140 |
}
|
| 141 |
*/
|
| 142 |
|
| 143 |
|
| 144 |
void ToplikeCandidate::recoTprimeUsingChi2(ElectronPointer electron) {
|
| 145 |
// cout << "Starting tp reco " << endl;
|
| 146 |
if (goodJets.size() < 3)
|
| 147 |
throw ReconstructionException("Not enough jets available to reconstruct Tprime using Chi2 method.");
|
| 148 |
double chosen_Chi2Total(9999999.);
|
| 149 |
leptonFromW = electron;
|
| 150 |
selectedNeutrino = 0;
|
| 151 |
currentSelectedNeutrino = 0;
|
| 152 |
reconstructNeutrinos();
|
| 153 |
|
| 154 |
// cout << "Starting jet loop" << endl;
|
| 155 |
for (unsigned short hadBindex = 0; hadBindex < goodJets.size(); ++hadBindex) {
|
| 156 |
for (unsigned short jet1Index = 0; jet1Index < goodJets.size(); ++jet1Index) {
|
| 157 |
if (jet1Index == hadBindex)
|
| 158 |
continue;
|
| 159 |
for (unsigned short jet2Index = 0; jet2Index < goodJets.size(); ++jet2Index) {
|
| 160 |
if (jet2Index == jet1Index || jet2Index == hadBindex)
|
| 161 |
continue;
|
| 162 |
hadronicBJet = goodJets.at(hadBindex);
|
| 163 |
leptonicBJet = hadronicBJet; // hack
|
| 164 |
jet1FromW = goodJets.at(jet1Index);
|
| 165 |
jet2FromW = goodJets.at(jet2Index);
|
| 166 |
leptonicW1 = ParticlePointer(new Particle(*neutrino1 + *leptonFromW));
|
| 167 |
leptonicW2 = ParticlePointer(new Particle(*neutrino2 + *leptonFromW));
|
| 168 |
leptonicTop1 = ParticlePointer(new Particle(*hadronicBJet + *leptonicW1));
|
| 169 |
leptonicTop2 = ParticlePointer(new Particle(*hadronicBJet + *leptonicW2));
|
| 170 |
hadronicW = ParticlePointer(new Particle(*jet1FromW + *jet2FromW));
|
| 171 |
// fillHypotheses();
|
| 172 |
selectNeutrinoSolution();
|
| 173 |
double chi2 = TopPairEventCandidate::getTotalChi2(currentSelectedNeutrino);
|
| 174 |
if (chi2 < chosen_Chi2Total) {
|
| 175 |
hadronicBIndex = hadBindex;
|
| 176 |
jet1FromWIndex = jet1Index;
|
| 177 |
jet2FromWIndex = jet2Index;
|
| 178 |
chosen_Chi2Total = chi2;
|
| 179 |
selectedNeutrino = currentSelectedNeutrino;
|
| 180 |
}
|
| 181 |
}
|
| 182 |
}
|
| 183 |
}
|
| 184 |
// cout << "Done with loop" << endl;
|
| 185 |
hadronicBJet = goodJets.at(hadronicBIndex);
|
| 186 |
leptonicBJet = hadronicBJet; // hack
|
| 187 |
jet1FromW = goodJets.at(jet1FromWIndex);
|
| 188 |
jet2FromW = goodJets.at(jet2FromWIndex);
|
| 189 |
hadronicW = ParticlePointer(new Particle(*jet1FromW + *jet2FromW));
|
| 190 |
leptonicW1 = ParticlePointer(new Particle(*neutrino1 + *leptonFromW));
|
| 191 |
leptonicW2 = ParticlePointer(new Particle(*neutrino2 + *leptonFromW));
|
| 192 |
leptonicTop1 = ParticlePointer(new Particle(*hadronicBJet + *leptonicW1));
|
| 193 |
leptonicTop2 = ParticlePointer(new Particle(*hadronicBJet + *leptonicW2));
|
| 194 |
ptTprimeSystem = PtTtbarSystem(selectedNeutrino);
|
| 195 |
htSystem = HTSystem();
|
| 196 |
if (selectedNeutrino == 1)
|
| 197 |
tPrime = ParticlePointer(new Particle(*leptonicTop1 + *hadronicW));
|
| 198 |
else
|
| 199 |
tPrime = ParticlePointer(new Particle(*leptonicTop2 + *hadronicW));
|
| 200 |
// cout << "tPrime mass = " << tPrime->mass() << endl;
|
| 201 |
doneReconstruction = true;
|
| 202 |
}
|
| 203 |
|
| 204 |
|
| 205 |
// Omit lepton-jet angle from criteria.
|
| 206 |
double ToplikeCandidate::getLeptonicChi2(double topMass, double unused)
|
| 207 |
const
|
| 208 |
{
|
| 209 |
double massDifference = TMath::Power(topMass - pdgTopMass, 2);
|
| 210 |
double massError = 2 * topMassWidth * topMassWidth;
|
| 211 |
double massTerm = massDifference / massError;
|
| 212 |
return (massTerm);
|
| 213 |
}
|
| 214 |
|
| 215 |
|
| 216 |
double ToplikeCandidate::getHadronicChi2() const {
|
| 217 |
double WmassDifference = TMath::Power(hadronicW->mass() - pdgWmass, 2);
|
| 218 |
double WmassError = 2 * wMassWidth * wMassWidth;
|
| 219 |
double WmassTerm = WmassDifference / WmassError;
|
| 220 |
|
| 221 |
double topMassDifference = TMath::Power(hadronicTop->mass() - pdgTopMass, 2);
|
| 222 |
double topMassError = 2 * topMassWidth * topMassWidth;
|
| 223 |
double topMassTerm = topMassDifference / topMassError;
|
| 224 |
return ((1.0 / sqrt(2.0)) * (topMassTerm + WmassTerm));
|
| 225 |
}
|
| 226 |
|
| 227 |
|
| 228 |
void ToplikeCandidate::recoBestSingleTop(LeptonPointer lepton) {
|
| 229 |
if (goodJets.size() < 3)
|
| 230 |
throw ReconstructionException("Not enough jets available to reconstruct top using Chi2 method.");
|
| 231 |
double bestLepChi2(9999999.), bestHadChi2(9999999.);
|
| 232 |
leptonFromW = lepton;
|
| 233 |
selectedNeutrino = 0;
|
| 234 |
currentSelectedNeutrino = 0;
|
| 235 |
reconstructNeutrinos();
|
| 236 |
|
| 237 |
for (unsigned short hadBindex = 0; hadBindex < goodJets.size(); ++hadBindex) {
|
| 238 |
for (unsigned short jet1Index = 0; jet1Index < goodJets.size(); ++jet1Index) {
|
| 239 |
if (hadBindex == jet1Index)
|
| 240 |
continue;
|
| 241 |
for (unsigned short jet2Index = 0; jet2Index < goodJets.size(); ++jet2Index) {
|
| 242 |
if (jet2Index == jet1Index || jet2Index == hadBindex)
|
| 243 |
continue;
|
| 244 |
hadronicBJet = goodJets.at(hadBindex);
|
| 245 |
leptonicBJet = hadronicBJet; // hack
|
| 246 |
jet1FromW = goodJets.at(jet1Index);
|
| 247 |
jet2FromW = goodJets.at(jet2Index);
|
| 248 |
leptonicW1 = ParticlePointer(new Particle(*neutrino1 + *leptonFromW));
|
| 249 |
leptonicW2 = ParticlePointer(new Particle(*neutrino2 + *leptonFromW));
|
| 250 |
leptonicTop1 = ParticlePointer(new Particle(*hadronicBJet + *leptonicW1));
|
| 251 |
leptonicTop2 = ParticlePointer(new Particle(*hadronicBJet + *leptonicW2));
|
| 252 |
hadronicW = ParticlePointer(new Particle(*jet1FromW + *jet2FromW));
|
| 253 |
hadronicTop = ParticlePointer(new Particle(*hadronicBJet + *hadronicW));
|
| 254 |
// fillHypotheses();
|
| 255 |
selectNeutrinoSolution();
|
| 256 |
double ltopchi2 = TopPairEventCandidate::getLeptonicChi2(currentSelectedNeutrino);
|
| 257 |
double htopchi2 = getHadronicChi2();
|
| 258 |
if (ltopchi2 < bestLepChi2) {
|
| 259 |
selectedNeutrino = currentSelectedNeutrino;
|
| 260 |
leptonicBIndex = hadBindex;
|
| 261 |
bestLepChi2 = ltopchi2;
|
| 262 |
}
|
| 263 |
if (htopchi2 < bestHadChi2) {
|
| 264 |
hadronicBIndex = hadBindex;
|
| 265 |
jet1FromWIndex = jet1Index;
|
| 266 |
jet2FromWIndex = jet2Index;
|
| 267 |
bestHadChi2 = htopchi2;
|
| 268 |
}
|
| 269 |
}
|
| 270 |
}
|
| 271 |
}
|
| 272 |
hadronicBJet = goodJets.at(hadronicBIndex);
|
| 273 |
leptonicBJet = goodJets.at(leptonicBIndex);
|
| 274 |
jet1FromW = goodJets.at(jet1FromWIndex);
|
| 275 |
jet2FromW = goodJets.at(jet2FromWIndex);
|
| 276 |
hadronicW = ParticlePointer(new Particle(*jet1FromW + *jet2FromW));
|
| 277 |
leptonicW1 = ParticlePointer(new Particle(*neutrino1 + *leptonFromW));
|
| 278 |
leptonicW2 = ParticlePointer(new Particle(*neutrino2 + *leptonFromW));
|
| 279 |
hadronicTop = ParticlePointer(new Particle(*hadronicBJet + *hadronicW));
|
| 280 |
leptonicTop1 = ParticlePointer(new Particle(*leptonicBJet + *leptonicW1));
|
| 281 |
leptonicTop2 = ParticlePointer(new Particle(*leptonicBJet + *leptonicW2));
|
| 282 |
// cout << "Found lone tops " << endl;
|
| 283 |
doneReconstructiontop = true;
|
| 284 |
doneReconstruction = true;
|
| 285 |
}
|
| 286 |
|
| 287 |
|
| 288 |
|
| 289 |
static bool cmpMcEntries(const mcObj &firstArg, const mcObj &secondArg)
|
| 290 |
{
|
| 291 |
return (firstArg.deltaR < secondArg.deltaR);
|
| 292 |
}
|
| 293 |
|
| 294 |
|
| 295 |
mcList ToplikeCandidate::getMCList(const Particle *const recoObject) const
|
| 296 |
{
|
| 297 |
const MCParticleCollection &mcColl = genParticles;
|
| 298 |
mcList mcObjs;
|
| 299 |
if (mcColl.size() > 0) {
|
| 300 |
for (unsigned int i = 0; i < mcColl.size(); ++i) {
|
| 301 |
mcObj mcEntry;
|
| 302 |
mcEntry.deltaR = recoObject->deltaR(mcColl.at(i));
|
| 303 |
mcEntry.ptr = mcColl.at(i);
|
| 304 |
mcEntry.gpInd = i;
|
| 305 |
mcObjs.push_back(mcEntry);
|
| 306 |
}
|
| 307 |
sort(mcObjs.begin(), mcObjs.end(), cmpMcEntries);
|
| 308 |
}
|
| 309 |
return (mcObjs);
|
| 310 |
}
|
| 311 |
|
| 312 |
|
| 313 |
mcList ToplikeCandidate::getMCListForJets(const Jet *const recoObject)
|
| 314 |
const
|
| 315 |
{
|
| 316 |
const JetCollection &mcColl = genJets;
|
| 317 |
mcList mcObjs;
|
| 318 |
if (mcColl.size() > 0) {
|
| 319 |
for (unsigned int i = 0; i < mcColl.size(); ++i) {
|
| 320 |
mcObj mcEntry;
|
| 321 |
mcEntry.deltaR = recoObject->deltaR(mcColl.at(i));
|
| 322 |
mcEntry.gpInd = i;
|
| 323 |
mcObjs.push_back(mcEntry);
|
| 324 |
}
|
| 325 |
sort(mcObjs.begin(), mcObjs.end(), cmpMcEntries);
|
| 326 |
}
|
| 327 |
return (mcObjs);
|
| 328 |
}
|
| 329 |
|
| 330 |
|
| 331 |
int ToplikeCandidate::missRecoJetsForGenJets() const
|
| 332 |
{
|
| 333 |
const JetCollection &mcColl = genJets;
|
| 334 |
int numGenJets = mcColl.size();
|
| 335 |
int matchJets = 0;
|
| 336 |
for (int i = 0; i < numGenJets; ++i) {
|
| 337 |
bool match = false;
|
| 338 |
for (unsigned int recoInd = 0;
|
| 339 |
recoInd < GoodJets().size() && match == false; ++recoInd)
|
| 340 |
if (GoodJets().at(recoInd)->deltaR(mcColl.at(i)) < 0.3) {
|
| 341 |
++matchJets;
|
| 342 |
match = true;
|
| 343 |
}
|
| 344 |
}
|
| 345 |
return (numGenJets - matchJets);
|
| 346 |
}
|
| 347 |
|
| 348 |
|
| 349 |
|
| 350 |
topTruth ToplikeCandidate::getMCMatches(const ParticlePointer &recoObject)
|
| 351 |
const
|
| 352 |
{
|
| 353 |
const float nullDeltaR = 3;
|
| 354 |
bool deltaRSet = false;
|
| 355 |
topTruth result = {nullDeltaR, false, false};
|
| 356 |
mcList mcObjs = getMCList(recoObject.get());
|
| 357 |
if (mcObjs.size() <= 0)
|
| 358 |
return (result);
|
| 359 |
/*
|
| 360 |
if (mcObjs.begin()->first > 0.9) // Nothing close enough
|
| 361 |
return ("");
|
| 362 |
TString idList;
|
| 363 |
idList += mcObjs.at(0).first;
|
| 364 |
idList += " ";
|
| 365 |
for (unsigned int i = 0; i < mcObjs.size() && i < 12; ++i) {
|
| 366 |
idList += PdgStrs::pdgStr(mcObjs.at(i).second);
|
| 367 |
idList += " ";
|
| 368 |
}
|
| 369 |
*/
|
| 370 |
float topDeltaR = nullDeltaR;
|
| 371 |
bool wFound = false, lepFound = false;
|
| 372 |
for (unsigned int i = 0; i < mcObjs.size() && i < MCLIMIT &&
|
| 373 |
(deltaRSet == false || lepFound == false || wFound == false); ++i) {
|
| 374 |
if (deltaRSet == false && fabs(mcObjs.at(i).ptr->pdgId()) == TOP_ID) {
|
| 375 |
topDeltaR = mcObjs.at(i).deltaR;
|
| 376 |
deltaRSet = true;
|
| 377 |
if (topDeltaR > nullDeltaR)
|
| 378 |
topDeltaR = nullDeltaR;
|
| 379 |
}
|
| 380 |
if (fabs(mcObjs.at(i).ptr->pdgId()) == ELECTRON_ID ||
|
| 381 |
fabs(mcObjs.at(i).ptr->pdgId()) == MU_ID)
|
| 382 |
lepFound = true;
|
| 383 |
if (fabs(mcObjs.at(i).ptr->pdgId()) == W_ID)
|
| 384 |
wFound = true;
|
| 385 |
}
|
| 386 |
if (deltaRSet == false)
|
| 387 |
cout << " No top found in GPs! Increase MCLIMIT.\n" ;
|
| 388 |
result.deltaR = topDeltaR;
|
| 389 |
result.leptonic = wFound && lepFound;
|
| 390 |
result.hadronic = wFound;
|
| 391 |
return (result);
|
| 392 |
}
|
| 393 |
|
| 394 |
|
| 395 |
static int findMCParticle(const mcList &mcObjs, int pdgId, int antiId = 0,
|
| 396 |
unsigned int startInd = 0)
|
| 397 |
{
|
| 398 |
if (antiId == 0)
|
| 399 |
antiId = pdgId;
|
| 400 |
for (unsigned int i = startInd; i < mcObjs.size(); ++i) {
|
| 401 |
if (pdgId == 0 || mcObjs.at(i).ptr->pdgId() == pdgId ||
|
| 402 |
mcObjs.at(i).ptr->pdgId() == antiId) {
|
| 403 |
if (mcObjs.at(i).deltaR < DELTA_R_LIMIT)
|
| 404 |
return (i);
|
| 405 |
else return (-1);
|
| 406 |
}
|
| 407 |
}
|
| 408 |
return (-1);
|
| 409 |
}
|
| 410 |
|
| 411 |
|
| 412 |
// Scans up the decay chain when mother & daughter are the same.
|
| 413 |
|
| 414 |
bool ToplikeCandidate::chkIfChild(int gpInd, int motherGpInd) const
|
| 415 |
{
|
| 416 |
if (motherGpInd < 0)
|
| 417 |
return (false);
|
| 418 |
const MCParticleCollection &mcColl = genParticles;
|
| 419 |
if (mcColl.size() <= 0)
|
| 420 |
return (false);
|
| 421 |
int mother = mcColl.at(gpInd)->motherIndex();
|
| 422 |
if (mother == motherGpInd)
|
| 423 |
return (true);
|
| 424 |
while (mother > 0 &&
|
| 425 |
mcColl.at(gpInd)->pdgId() == mcColl.at(mother)->pdgId())
|
| 426 |
{
|
| 427 |
gpInd = mother;
|
| 428 |
mother = mcColl.at(gpInd)->motherIndex();
|
| 429 |
if (mother == motherGpInd)
|
| 430 |
return (true);
|
| 431 |
}
|
| 432 |
int gMother = mcColl.at(motherGpInd)->motherIndex();
|
| 433 |
|
| 434 |
// W' sample b's have W' as mother rather than top.
|
| 435 |
if (gMother > 0 &&
|
| 436 |
fabs(mcColl.at(gMother)->pdgId()) == WPRIME_ID && mother == gMother)
|
| 437 |
return (true);
|
| 438 |
while (gMother > 0 &&
|
| 439 |
mcColl.at(motherGpInd)->pdgId() == mcColl.at(gMother)->pdgId())
|
| 440 |
{
|
| 441 |
motherGpInd = gMother;
|
| 442 |
gMother = mcColl.at(motherGpInd)->motherIndex();
|
| 443 |
if (mother == motherGpInd)
|
| 444 |
return (true);
|
| 445 |
}
|
| 446 |
/*
|
| 447 |
if (mcColl.at(mother)->motherIndex() == mcColl.at(motherGpInd)->motherIndex())
|
| 448 |
{
|
| 449 |
cout << "gm match: ind cpdg mi mpdg gm pdg omi ompdg " << gpInd << " ";
|
| 450 |
cout << mcColl.at(gpInd)->pdgId() << " " << mother << " ";
|
| 451 |
cout << mcColl.at(mother)->pdgId() << " ";
|
| 452 |
int gm = mcColl.at(mother)->motherIndex();
|
| 453 |
cout << gm << " " << mcColl.at(gm)->pdgId() << " ";
|
| 454 |
cout << motherGpInd << " " << mcColl.at(motherGpInd)->pdgId() << endl;
|
| 455 |
return (true);
|
| 456 |
}
|
| 457 |
*/
|
| 458 |
return (false);
|
| 459 |
}
|
| 460 |
|
| 461 |
|
| 462 |
int ToplikeCandidate::chkMCMatch(const Particle *const recoObject,
|
| 463 |
int motherInd, const int pdgId, int antiId) const
|
| 464 |
{
|
| 465 |
mcList truList = getMCList(recoObject);
|
| 466 |
if (truList.size() > 0) {
|
| 467 |
int index = -1;
|
| 468 |
while ((index = findMCParticle(truList, pdgId, antiId, index + 1)) >= 0)
|
| 469 |
if (motherInd == -1 || chkIfChild(truList.at(index).gpInd, motherInd))
|
| 470 |
return (truList.at(index).gpInd);
|
| 471 |
}
|
| 472 |
return (-1);
|
| 473 |
}
|
| 474 |
|
| 475 |
|
| 476 |
// Look for a non-top quark whose mother is specified by motherInd and which
|
| 477 |
// doesn't match badInd. If motherInd == -1, no mother is checked for.
|
| 478 |
|
| 479 |
int ToplikeCandidate::chkForQuark(const JetPointer recoObject,
|
| 480 |
const int motherInd, const int badInd) const
|
| 481 |
{
|
| 482 |
mcList truList = getMCList(recoObject.get());
|
| 483 |
if (truList.size() > 0) {
|
| 484 |
for (int i = 0; i < (int) truList.size() &&
|
| 485 |
truList.at(i).deltaR < DELTA_R_LIMIT; ++i) {
|
| 486 |
if (truList.at(i).ptr->isQuark() && fabs(truList.at(i).ptr->pdgId()) !=
|
| 487 |
TOP_ID && (motherInd == -1 || chkIfChild(truList.at(i).gpInd, motherInd)) && truList.at(i).gpInd != badInd)
|
| 488 |
return (truList.at(i).gpInd);
|
| 489 |
// Checking for q-qbar pair not necessary since both daughters of mother
|
| 490 |
}
|
| 491 |
}
|
| 492 |
return (-1);
|
| 493 |
}
|
| 494 |
|
| 495 |
|
| 496 |
// Checks if ind1 is a b quark, other two are daughters of a W, and
|
| 497 |
// b & W are siblings.
|
| 498 |
bool ToplikeCandidate::chkBW(const int ind1, const int ind2, const int ind3)
|
| 499 |
const
|
| 500 |
{
|
| 501 |
const MCParticleCollection &mcColl = genParticles;
|
| 502 |
if (mcColl.size() <= 0)
|
| 503 |
return (false);
|
| 504 |
if (fabs(mcColl.at(ind1)->pdgId()) != BQ_ID)
|
| 505 |
return (false);
|
| 506 |
int jetMother = mcColl.at(ind2)->motherIndex();
|
| 507 |
if (jetMother > 0 && fabs(mcColl.at(jetMother)->pdgId()) == W_ID &&
|
| 508 |
mcColl.at(ind3)->motherIndex() == jetMother &&
|
| 509 |
mcColl.at(ind1)->motherIndex() == mcColl.at(jetMother)->motherIndex())
|
| 510 |
return (true);
|
| 511 |
return (false);
|
| 512 |
}
|
| 513 |
|
| 514 |
|
| 515 |
int ToplikeCandidate::goodHTopMatching() const
|
| 516 |
{
|
| 517 |
int numFound = 0, mostFound = 0;
|
| 518 |
mcList bTruList = getMCList(hadronicBJet.get());
|
| 519 |
mcList j1TruList = getMCList(jet1FromW.get());
|
| 520 |
mcList j2TruList = getMCList(jet2FromW.get());
|
| 521 |
int bInd = -1;
|
| 522 |
while (bTruList.size() > 0 &&
|
| 523 |
(bInd = findMCParticle(bTruList, 0, 0, bInd + 1)) >= 0)
|
| 524 |
{
|
| 525 |
if (bTruList.at(bInd).ptr->isQuark())
|
| 526 |
{
|
| 527 |
numFound = 1;
|
| 528 |
int bGPInd = bTruList.at(bInd).gpInd;
|
| 529 |
int j1Ind = -1;
|
| 530 |
while (j1TruList.size() > 0 &&
|
| 531 |
(j1Ind = findMCParticle(j1TruList, 0, 0, j1Ind + 1)) >= 0)
|
| 532 |
{
|
| 533 |
int j1GPInd = j1TruList.at(j1Ind).gpInd;
|
| 534 |
if (j1GPInd != bGPInd && j1TruList.at(j1Ind).ptr->isQuark())
|
| 535 |
{
|
| 536 |
numFound = 2;
|
| 537 |
int j2Ind = -1;
|
| 538 |
while (j2TruList.size() > 0 &&
|
| 539 |
(j2Ind = findMCParticle(j2TruList, 0, 0, j2Ind + 1)) >= 0)
|
| 540 |
{
|
| 541 |
int j2GPInd = j2TruList.at(j2Ind).gpInd;
|
| 542 |
if (j2GPInd != bGPInd && j2GPInd != j1GPInd &&
|
| 543 |
j2TruList.at(j2Ind).ptr->isQuark())
|
| 544 |
{
|
| 545 |
numFound = 3;
|
| 546 |
if (chkBW(bGPInd, j1GPInd, j2GPInd) ||
|
| 547 |
chkBW(j1GPInd, bGPInd, j2GPInd) ||
|
| 548 |
chkBW(j2GPInd, j1GPInd, bGPInd))
|
| 549 |
return (4);
|
| 550 |
}
|
| 551 |
}
|
| 552 |
if (numFound > mostFound)
|
| 553 |
mostFound = numFound;
|
| 554 |
}
|
| 555 |
}
|
| 556 |
}
|
| 557 |
}
|
| 558 |
if (numFound > mostFound)
|
| 559 |
mostFound = numFound;
|
| 560 |
return (mostFound);
|
| 561 |
}
|
| 562 |
|
| 563 |
|
| 564 |
|
| 565 |
int ToplikeCandidate::strictHTopMatching(const mcObj &topTru) const
|
| 566 |
{
|
| 567 |
int numFound = 1; // Start with finding the top.
|
| 568 |
short antiMult = 1;
|
| 569 |
if (topTru.ptr->pdgId() < 0)
|
| 570 |
antiMult = -1;
|
| 571 |
// cout << " hti " << topTru.gpInd << " ";
|
| 572 |
if (chkMCMatch(hadronicBJet.get(), topTru.gpInd, BQ_ID * antiMult) >= 0)
|
| 573 |
++numFound;
|
| 574 |
int wInd = chkMCMatch(hadronicW.get(), topTru.gpInd, W_ID * antiMult);
|
| 575 |
if (wInd >= 0)
|
| 576 |
++numFound;
|
| 577 |
else return (numFound);
|
| 578 |
|
| 579 |
// Try to find 2 quarks from W. If both fail, trying again will still fail
|
| 580 |
// but is harmless. If one fails, try again to see if can get both if we
|
| 581 |
// let the other jet grab the 1st quark. If that doesn't work, one quark
|
| 582 |
// will still be found.
|
| 583 |
bool pairFnd = false, last1 = false;
|
| 584 |
int takenInd2 = -1;
|
| 585 |
int takenInd = chkForQuark(jet1FromW, wInd);
|
| 586 |
if (takenInd >= 0 && chkForQuark(jet2FromW, wInd, takenInd) >= 0) {
|
| 587 |
numFound += 2;
|
| 588 |
pairFnd = true;
|
| 589 |
} else {
|
| 590 |
takenInd2 = chkForQuark(jet2FromW, wInd);
|
| 591 |
if (takenInd2 >= 0)
|
| 592 |
++numFound;
|
| 593 |
if (chkForQuark(jet1FromW, wInd, takenInd2) >= 0) {
|
| 594 |
++numFound;
|
| 595 |
last1 = true;
|
| 596 |
}
|
| 597 |
}
|
| 598 |
/*
|
| 599 |
if (numFound == 5 || numFound == 4) {
|
| 600 |
cout << "num pair wInd takenInd taken2 last1 " << numFound << " ";
|
| 601 |
cout << pairFnd << " " << wInd << " ";
|
| 602 |
cout << takenInd << " " << takenInd2 << " " << last1 << endl;
|
| 603 |
}
|
| 604 |
*/
|
| 605 |
return (numFound);
|
| 606 |
}
|
| 607 |
|
| 608 |
|
| 609 |
int ToplikeCandidate::strictLTopMatching(const mcObj &topTru) const
|
| 610 |
{
|
| 611 |
int numFound = 1;
|
| 612 |
short antiMult = 1;
|
| 613 |
if (topTru.ptr->pdgId() < 0)
|
| 614 |
antiMult = -1;
|
| 615 |
if (chkMCMatch(getLeptonicBJet().get(), topTru.gpInd, BQ_ID * antiMult) >= 0)
|
| 616 |
++numFound;
|
| 617 |
int wInd = chkMCMatch(getLeptonicW().get(), topTru.gpInd, W_ID * antiMult);
|
| 618 |
if (wInd >= 0)
|
| 619 |
++numFound;
|
| 620 |
else return (numFound);
|
| 621 |
// Lepton is opposite from top to preserve lepton number.
|
| 622 |
if (chkMCMatch(getLeptonFromWDecay().get(), wInd, -ELECTRON_ID * antiMult)
|
| 623 |
>= 0 ||
|
| 624 |
chkMCMatch(getLeptonFromWDecay().get(), wInd, -MU_ID * antiMult) >= 0)
|
| 625 |
++numFound;
|
| 626 |
if (chkMCMatch(getNeutrinoFromWDecay().get(), wInd, NU_E_ID * antiMult) >= 0
|
| 627 |
|| chkMCMatch(getNeutrinoFromWDecay().get(), wInd, NU_MU_ID * antiMult)
|
| 628 |
>= 0)
|
| 629 |
++numFound;
|
| 630 |
return (numFound);
|
| 631 |
}
|
| 632 |
|
| 633 |
|
| 634 |
int ToplikeCandidate::goodLTopMatching() const
|
| 635 |
{
|
| 636 |
const MCParticleCollection &mcColl = genParticles;
|
| 637 |
int numParticles = mcColl.size();
|
| 638 |
if (numParticles <= 0)
|
| 639 |
return (0);
|
| 640 |
int numFound = 0, mostFound = 0;
|
| 641 |
mcList lepTruList = getMCList(getLeptonFromWDecay().get());
|
| 642 |
mcList bTruList = getMCList(getLeptonicBJet().get());
|
| 643 |
for (int lepTypes = 0; lepTypes < 2; ++lepTypes) {
|
| 644 |
int lepInd = -1;
|
| 645 |
int lepPdgId = ELECTRON_ID;
|
| 646 |
if (lepTypes == 1)
|
| 647 |
lepPdgId = MU_ID;
|
| 648 |
while (lepTruList.size() > 0 && (lepInd =
|
| 649 |
findMCParticle(lepTruList, lepPdgId, -lepPdgId, lepInd + 1)) >= 0)
|
| 650 |
{
|
| 651 |
numFound = 1;
|
| 652 |
int bInd = -1;
|
| 653 |
int lepGPInd = lepTruList.at(lepInd).gpInd;
|
| 654 |
if (lepGPInd >= 0 && lepGPInd < numParticles) {
|
| 655 |
int lepMothInd = mcColl.at(lepGPInd)->motherIndex();
|
| 656 |
if (lepMothInd >= 0 && lepMothInd < numParticles &&
|
| 657 |
fabs(mcColl.at(lepMothInd)->pdgId()) == W_ID) {
|
| 658 |
numFound = 2;
|
| 659 |
while (bTruList.size() > 0 &&
|
| 660 |
(bInd = findMCParticle(bTruList, BQ_ID, -BQ_ID, bInd + 1)) >= 0)
|
| 661 |
{
|
| 662 |
numFound = 3;
|
| 663 |
int bGPInd = bTruList.at(bInd).gpInd;
|
| 664 |
if ((bGPInd >= 0 && bGPInd < numParticles) &&
|
| 665 |
(chkIfChild(bGPInd, mcColl.at(lepMothInd)->motherIndex()) ||
|
| 666 |
chkIfChild(lepMothInd, mcColl.at(bGPInd)->motherIndex())))
|
| 667 |
return(4);
|
| 668 |
}
|
| 669 |
}
|
| 670 |
}
|
| 671 |
if (numFound > mostFound)
|
| 672 |
mostFound = numFound;
|
| 673 |
}
|
| 674 |
}
|
| 675 |
return (mostFound);
|
| 676 |
}
|
| 677 |
|
| 678 |
|
| 679 |
|
| 680 |
/*
|
| 681 |
int ToplikeCandidate::getNumMCPermutsForTop(const ParticlePointer &recoTop,
|
| 682 |
bool hadronic) const
|
| 683 |
{
|
| 684 |
mcList topTru = getMCList(recoTop.get());
|
| 685 |
if (topTru.size() <= 0)
|
| 686 |
return (0);
|
| 687 |
int highestNumFound = 0;
|
| 688 |
int topInd = -1;
|
| 689 |
while ((topInd = findMCParticle(topTru, TOP_ID, -TOP_ID, topInd + 1)) >= 0) {
|
| 690 |
int numFound = 0;
|
| 691 |
if (hadronic)
|
| 692 |
numFound = hTopDauMatching(topTru.at(topInd));
|
| 693 |
else numFound = lTopDauMatching(topTru.at(topInd));
|
| 694 |
if (numFound > highestNumFound)
|
| 695 |
highestNumFound = numFound;
|
| 696 |
}
|
| 697 |
return (highestNumFound);
|
| 698 |
}
|
| 699 |
*/
|
| 700 |
|
| 701 |
|
| 702 |
int ToplikeCandidate::getNumMCMatchesForTop(const ParticlePointer &recoTop,
|
| 703 |
bool hadronic)
|
| 704 |
const
|
| 705 |
{
|
| 706 |
mcList topTru = getMCList(recoTop.get());
|
| 707 |
if (topTru.size() <= 0)
|
| 708 |
return (0);
|
| 709 |
int highestNumFound = 0;
|
| 710 |
int topInd = -1;
|
| 711 |
while ((topInd = findMCParticle(topTru, TOP_ID, -TOP_ID, topInd + 1)) >= 0) {
|
| 712 |
int numFound = 0;
|
| 713 |
if (hadronic)
|
| 714 |
numFound = strictHTopMatching(topTru.at(topInd));
|
| 715 |
else numFound = strictLTopMatching(topTru.at(topInd));
|
| 716 |
if (numFound > highestNumFound)
|
| 717 |
highestNumFound = numFound;
|
| 718 |
}
|
| 719 |
return (highestNumFound);
|
| 720 |
}
|
| 721 |
|
| 722 |
|
| 723 |
int ToplikeCandidate::getNumMCMatchesHTop() const
|
| 724 |
{
|
| 725 |
return (getNumMCMatchesForTop(getHadronicTop(), true));
|
| 726 |
}
|
| 727 |
|
| 728 |
/*
|
| 729 |
int ToplikeCandidate::getNumDauMatchesHTop() const
|
| 730 |
{
|
| 731 |
return (getNumMCPermutsForTop(getHadronicTop(), true));
|
| 732 |
}
|
| 733 |
|
| 734 |
|
| 735 |
int ToplikeCandidate::getNumDauMatchesLTop() const
|
| 736 |
{
|
| 737 |
return (getNumMCPermutsForTop(getLeptonicTop(), false));
|
| 738 |
}
|
| 739 |
*/
|
| 740 |
|
| 741 |
|
| 742 |
int ToplikeCandidate::getNumCorrectIDHTop() const
|
| 743 |
{
|
| 744 |
mcList topTru = getMCList(getHadronicTop().get());
|
| 745 |
int numFound = 0;
|
| 746 |
bool topFnd = false, bFnd = false, jet1Fnd = false, jet2Fnd = false;
|
| 747 |
if (topTru.size() > 0 && findMCParticle(topTru, TOP_ID, -TOP_ID) >= 0) {
|
| 748 |
topFnd = true;
|
| 749 |
++numFound;
|
| 750 |
}
|
| 751 |
|
| 752 |
// Take either particle or anti-particle.
|
| 753 |
if (chkMCMatch(hadronicBJet.get(), -1, BQ_ID, -BQ_ID) >= 0) {
|
| 754 |
++numFound;
|
| 755 |
bFnd = true;
|
| 756 |
}
|
| 757 |
int wInd = chkMCMatch(hadronicW.get(), -1, W_ID, -W_ID);
|
| 758 |
if (wInd >= 0)
|
| 759 |
++numFound;
|
| 760 |
|
| 761 |
// Look for any quark near a jet. Even allow the same quark for both jets.
|
| 762 |
if (chkForQuark(jet1FromW, -1) >= 0) {
|
| 763 |
++numFound;
|
| 764 |
jet1Fnd = true;
|
| 765 |
}
|
| 766 |
if (chkForQuark(jet2FromW, -1) >= 0) {
|
| 767 |
jet2Fnd = true;
|
| 768 |
++numFound;
|
| 769 |
}
|
| 770 |
/*
|
| 771 |
if (numFound == 4 || numFound == 3) {
|
| 772 |
cout << "had num t b w j1 j2 " << numFound << " " << topFnd << " " << bFnd << " " << wInd;
|
| 773 |
cout << " " << jet1Fnd << " " << jet2Fnd << endl;
|
| 774 |
}
|
| 775 |
*/
|
| 776 |
return (numFound);
|
| 777 |
}
|
| 778 |
|
| 779 |
|
| 780 |
int ToplikeCandidate::getNumMCMatchesLTop() const
|
| 781 |
{
|
| 782 |
return (getNumMCMatchesForTop(getLeptonicTop(), false));
|
| 783 |
}
|
| 784 |
|
| 785 |
|
| 786 |
int ToplikeCandidate::getNumCorrectIDLTop() const
|
| 787 |
{
|
| 788 |
mcList topTru = getMCList(getLeptonicTop().get());
|
| 789 |
int numFound = 0;
|
| 790 |
if (topTru.size() > 0 && findMCParticle(topTru, TOP_ID, -TOP_ID) >= 0)
|
| 791 |
++numFound;
|
| 792 |
if (chkMCMatch(getLeptonicBJet().get(), -1, BQ_ID, -BQ_ID) >= 0)
|
| 793 |
++numFound;
|
| 794 |
if (chkMCMatch(getLeptonicW().get(), -1, W_ID, -W_ID) >= 0)
|
| 795 |
++numFound;
|
| 796 |
if (chkMCMatch(getLeptonFromWDecay().get(), -1, -ELECTRON_ID,
|
| 797 |
ELECTRON_ID) >= 0 ||
|
| 798 |
chkMCMatch(getLeptonFromWDecay().get(), -1, -MU_ID, MU_ID) >= 0)
|
| 799 |
++numFound;
|
| 800 |
if (chkMCMatch(getNeutrinoFromWDecay().get(), -1, NU_E_ID, -NU_E_ID) >= 0 ||
|
| 801 |
chkMCMatch(getNeutrinoFromWDecay().get(), -1, NU_MU_ID, -NU_MU_ID) >= 0)
|
| 802 |
++numFound;
|
| 803 |
return (numFound);
|
| 804 |
}
|
| 805 |
|
| 806 |
|
| 807 |
int ToplikeCandidate::getTausInJets() const
|
| 808 |
{
|
| 809 |
int numTaus = 0;
|
| 810 |
for (unsigned int index = 0; index < GoodJets().size(); ++index) {
|
| 811 |
mcList truList = getMCList(GoodJets().at(index).get());
|
| 812 |
int tauIndex = -1;
|
| 813 |
if (truList.size() > 0 &&
|
| 814 |
(tauIndex = findMCParticle(truList, TAU_ID, -TAU_ID)) >= 0)
|
| 815 |
++numTaus;
|
| 816 |
}
|
| 817 |
// cout << " Found taus: " << numTaus;
|
| 818 |
return (numTaus);
|
| 819 |
}
|
| 820 |
|
| 821 |
|
| 822 |
int ToplikeCandidate::numTaus() const
|
| 823 |
{
|
| 824 |
const MCParticleCollection &mcColl = genParticles;
|
| 825 |
int numTaus = 0;
|
| 826 |
if (mcColl.size() > 0) {
|
| 827 |
for (unsigned int i = 0; i < mcColl.size(); ++i) {
|
| 828 |
if (fabs(mcColl.at(i)->pdgId()) == TAU_ID)
|
| 829 |
++numTaus;
|
| 830 |
}
|
| 831 |
}
|
| 832 |
return (numTaus);
|
| 833 |
}
|
| 834 |
|
| 835 |
|
| 836 |
int ToplikeCandidate::missGenJetsForRecoJets() const
|
| 837 |
{
|
| 838 |
int matchJets = 0, numJets = GoodJets().size();
|
| 839 |
for (int index = 0; index < numJets; ++index) {
|
| 840 |
mcList mcObjs = getMCListForJets(GoodJets().at(index).get());
|
| 841 |
if (mcObjs.size() > 0 && mcObjs.at(0).deltaR < 0.3)
|
| 842 |
++matchJets;
|
| 843 |
}
|
| 844 |
return (numJets - matchJets);
|
| 845 |
}
|
| 846 |
|
| 847 |
|
| 848 |
|
| 849 |
TopPlusXCandidates::TopPlusXCandidates() :
|
| 850 |
dJetFromWpIndex(-1)
|
| 851 |
{
|
| 852 |
}
|
| 853 |
|
| 854 |
|
| 855 |
TopPlusXCandidates::TopPlusXCandidates(const Event& event) :
|
| 856 |
ToplikeCandidate(event),
|
| 857 |
dJetFromWpIndex(-1)
|
| 858 |
{
|
| 859 |
}
|
| 860 |
|
| 861 |
|
| 862 |
TopPlusXCandidates::~TopPlusXCandidates()
|
| 863 |
{
|
| 864 |
}
|
| 865 |
|
| 866 |
|
| 867 |
const JetPointer TopPlusXCandidates::nulljet;
|
| 868 |
|
| 869 |
|
| 870 |
bool TopPlusXCandidates::passesMETCut() const {
|
| 871 |
return met->pt() > 20.;
|
| 872 |
}
|
| 873 |
|
| 874 |
|
| 875 |
bool TopPlusXCandidates::passesHTCut() const {
|
| 876 |
return fullHT() > 700.;
|
| 877 |
}
|
| 878 |
|
| 879 |
|
| 880 |
bool TopPlusXCandidates::hasAtLeastOneGoodJet() const {
|
| 881 |
return goodJets.size() >= 1 && goodJets.front()->pt() > 180.0;
|
| 882 |
}
|
| 883 |
|
| 884 |
|
| 885 |
bool TopPlusXCandidates::hasAtLeastTwoGoodJets() const {
|
| 886 |
return goodJets.size() >= 2 && goodJets.at(1)->pt() > 90.0;
|
| 887 |
}
|
| 888 |
|
| 889 |
|
| 890 |
bool TopPlusXCandidates::hasOnlyOneGoodIsolatedMuon() const {
|
| 891 |
return goodIsolatedMuons.size() == 1;
|
| 892 |
}
|
| 893 |
|
| 894 |
|
| 895 |
bool TopPlusXCandidates::hasNoGoodIsolatedMuon() const {
|
| 896 |
return goodIsolatedMuons.size() <= 0;
|
| 897 |
}
|
| 898 |
|
| 899 |
|
| 900 |
bool TopPlusXCandidates::hasNoGoodIsolatedElectron() const {
|
| 901 |
if(Event::usePFIsolation)
|
| 902 |
return goodPFIsolatedElectrons.size() <= 0;
|
| 903 |
else
|
| 904 |
return goodIsolatedElectrons.size() <= 0;
|
| 905 |
}
|
| 906 |
|
| 907 |
|
| 908 |
bool TopPlusXCandidates::hasNoIsolatedElectron() const {
|
| 909 |
return looseElectrons.size() <= 0;
|
| 910 |
}
|
| 911 |
|
| 912 |
bool TopPlusXCandidates::passesExtraElectronCut() const{
|
| 913 |
ElectronPointer electron;
|
| 914 |
if(Event::usePFIsolation) {
|
| 915 |
if (goodPFIsolatedElectrons.size() <= 0)
|
| 916 |
return (false);
|
| 917 |
electron = GoodPFIsolatedElectrons().front();
|
| 918 |
} else if (goodIsolatedElectrons.size() > 0) {
|
| 919 |
electron = GoodIsolatedElectrons().front();
|
| 920 |
} else return (false);
|
| 921 |
return (electron->pt() > 45.);
|
| 922 |
}
|
| 923 |
|
| 924 |
|
| 925 |
bool TopPlusXCandidates::passesSelectionStep(enum
|
| 926 |
ToplikeSelectionSteps::Step step) const {
|
| 927 |
switch (step) {
|
| 928 |
case ToplikeSelectionSteps::AtLeastFiveGoodJets:
|
| 929 |
return (hasAtLeastFiveGoodJets());
|
| 930 |
case ToplikeSelectionSteps::FortyfiveGeVElectron:
|
| 931 |
return (passesExtraElectronCut());
|
| 932 |
case ToplikeSelectionSteps::OneIsolatedMuon:
|
| 933 |
return (hasOnlyOneGoodIsolatedMuon());
|
| 934 |
case ToplikeSelectionSteps::LooseElectronVeto:
|
| 935 |
return (hasNoIsolatedElectron());
|
| 936 |
case ToplikeSelectionSteps::TightElectronVeto:
|
| 937 |
return (hasNoGoodIsolatedElectron());
|
| 938 |
case ToplikeSelectionSteps::TightMuonVeto:
|
| 939 |
return (hasNoGoodIsolatedMuon());
|
| 940 |
case ToplikeSelectionSteps::HT700GeV:
|
| 941 |
return (passesHTCut());
|
| 942 |
default:
|
| 943 |
return (TopPairEventCandidate::passesSelectionStep(selSteps[step]));
|
| 944 |
}
|
| 945 |
}
|
| 946 |
|
| 947 |
|
| 948 |
bool TopPlusXCandidates::passesSelectionStepUpTo(unsigned int step,
|
| 949 |
const ToplikeSelectionSteps::Step steps[]) const {
|
| 950 |
if (step <= 0)
|
| 951 |
return (passesSelectionStep(steps[0]));
|
| 952 |
else {
|
| 953 |
return (passesSelectionStep(steps[step]) &&
|
| 954 |
passesSelectionStepUpTo(step - 1, steps));
|
| 955 |
}
|
| 956 |
}
|
| 957 |
|
| 958 |
|
| 959 |
bool TopPlusXCandidates::passesFullEPlusJetSelection() const {
|
| 960 |
return (passesSelectionStepUpTo(toplikeElectronSelSize - 1, toplikeElectronSelection));
|
| 961 |
}
|
| 962 |
|
| 963 |
|
| 964 |
bool TopPlusXCandidates::passesFullMuPlusJetSelection() const {
|
| 965 |
return (passesSelectionStepUpTo(toplikeMuonSelSize - 1, toplikeMuonSelection));
|
| 966 |
}
|
| 967 |
|
| 968 |
|
| 969 |
bool TopPlusXCandidates::passesFullEPlusJetSelectionKK() const {
|
| 970 |
return (passesSelectionStepUpTo(toplikeElectronSelSize - 1, toplikeElectronSelection));
|
| 971 |
}
|
| 972 |
|
| 973 |
|
| 974 |
double TopPlusXCandidates::jetsHTNoTop(unsigned int numJets) const {
|
| 975 |
double ht = 0.0;
|
| 976 |
JetCollection topJets;
|
| 977 |
|
| 978 |
for (unsigned int index = 0; index < goodJets.size(); ++index) {
|
| 979 |
if (goodJets.at(index) != leptonicBJet) {
|
| 980 |
if (topJets.size() <= 0 || topJets.back()->pt() > goodJets.at(index)->pt())
|
| 981 |
topJets.push_back(goodJets.at(index));
|
| 982 |
else {
|
| 983 |
unsigned int srchind = 0;
|
| 984 |
while (topJets.at(srchind)->pt() > goodJets.at(index)->pt())
|
| 985 |
++srchind;
|
| 986 |
topJets.insert(topJets.begin() + srchind, goodJets.at(index));
|
| 987 |
}
|
| 988 |
}
|
| 989 |
}
|
| 990 |
for (unsigned int index = 0; index < topJets.size() && index < numJets; ++index) {
|
| 991 |
ht += topJets.at(index)->pt();
|
| 992 |
}
|
| 993 |
return ht;
|
| 994 |
}
|
| 995 |
|
| 996 |
|
| 997 |
double TopPlusXCandidates::mass2jets(unsigned int jet1, unsigned int jet2) const
|
| 998 |
{
|
| 999 |
int hijets[3] = {-1, -1, -1};
|
| 1000 |
unsigned int hijindex = 0;
|
| 1001 |
|
| 1002 |
for (unsigned int index = 0; index < goodJets.size(); ++index) {
|
| 1003 |
if (goodJets.at(index) != leptonicBJet && hijets[hijindex] == -1) {
|
| 1004 |
hijets[hijindex] = index;
|
| 1005 |
if (hijindex < 2)
|
| 1006 |
hijindex++;
|
| 1007 |
}
|
| 1008 |
}
|
| 1009 |
if (hijets[jet1] > -1 && hijets[jet2] > -1) {
|
| 1010 |
ParticlePointer candidate =
|
| 1011 |
ParticlePointer(new Particle(*goodJets.at(hijets[jet1]) +
|
| 1012 |
*goodJets.at(hijets[jet2])));
|
| 1013 |
return (candidate->mass());
|
| 1014 |
}
|
| 1015 |
return (0.0);
|
| 1016 |
}
|
| 1017 |
|
| 1018 |
|
| 1019 |
double TopPlusXCandidates::mass4jets() const
|
| 1020 |
{
|
| 1021 |
unsigned int cnt = 0;
|
| 1022 |
|
| 1023 |
ParticlePointer candidate;
|
| 1024 |
for (unsigned int index = 0; index < goodJets.size(); ++index) {
|
| 1025 |
if (goodJets.at(index) != leptonicBJet) {
|
| 1026 |
if (cnt == 0) {
|
| 1027 |
candidate = ParticlePointer(new Particle(*goodJets.at(index)));
|
| 1028 |
} else if (cnt < 3) {
|
| 1029 |
candidate = ParticlePointer(new Particle(*goodJets.at(index) + *candidate));
|
| 1030 |
}
|
| 1031 |
++cnt;
|
| 1032 |
}
|
| 1033 |
}
|
| 1034 |
if (cnt > 0)
|
| 1035 |
return (candidate->mass());
|
| 1036 |
return (0.0);
|
| 1037 |
}
|
| 1038 |
|
| 1039 |
|
| 1040 |
JetPointer TopPlusXCandidates::leadnontopjet() const {
|
| 1041 |
double hiJetPt = 0;
|
| 1042 |
unsigned int hiInd = 0;
|
| 1043 |
|
| 1044 |
for (unsigned int index = 0; index < goodJets.size(); ++index) {
|
| 1045 |
if (goodJets.at(index) != leptonicBJet && goodJets.at(index)->pt() > hiJetPt) {
|
| 1046 |
hiJetPt = goodJets.at(index)->pt();
|
| 1047 |
hiInd = index;
|
| 1048 |
}
|
| 1049 |
}
|
| 1050 |
return (goodJets.at(hiInd));
|
| 1051 |
}
|
| 1052 |
|
| 1053 |
|
| 1054 |
// ST is the sum of jet and lepton pT minus the pT of the leading lepton and jet.
|
| 1055 |
double TopPlusXCandidates::ST() const {
|
| 1056 |
double st = fullHT() - met->pt();
|
| 1057 |
|
| 1058 |
double hiLepPt = 0;
|
| 1059 |
for (unsigned int index = 0; index < goodPFIsolatedElectrons.size(); ++index) {
|
| 1060 |
if (goodPFIsolatedElectrons.at(index)->pt() > hiLepPt)
|
| 1061 |
hiLepPt = goodPFIsolatedElectrons.at(index)->pt();
|
| 1062 |
}
|
| 1063 |
for (unsigned int index = 0; index < goodIsolatedMuons.size(); ++index) {
|
| 1064 |
if (goodIsolatedMuons.at(index)->pt() > hiLepPt)
|
| 1065 |
hiLepPt = goodIsolatedMuons.at(index)->pt();
|
| 1066 |
}
|
| 1067 |
double hiJetPt = 0;
|
| 1068 |
for (unsigned int index = 0; index < goodJets.size(); ++index) {
|
| 1069 |
if (goodJets.at(index)->pt() > hiJetPt)
|
| 1070 |
hiJetPt = goodJets.at(index)->pt();
|
| 1071 |
}
|
| 1072 |
st -= hiLepPt + hiJetPt;
|
| 1073 |
return st;
|
| 1074 |
}
|
| 1075 |
|
| 1076 |
double TopPlusXCandidates::getTopMetAngle() const
|
| 1077 |
{
|
| 1078 |
return (getLeptonicTop()->getFourVector().Theta() - met->getFourVector().Theta());
|
| 1079 |
}
|
| 1080 |
|
| 1081 |
|
| 1082 |
void TopPlusXCandidates::recoBestLeptoTop(LeptonPointer lepton,
|
| 1083 |
const JetPointer hadTopJet, const JetPointer wJet1, const JetPointer wJet2)
|
| 1084 |
{
|
| 1085 |
if (goodBJets.size() < 1)
|
| 1086 |
throw ReconstructionException("Not enough jets available to reconstruct top using Chi2 method.");
|
| 1087 |
double bestLepChi2(9999999.);
|
| 1088 |
// JetCollection *topJets = &goodBJets;
|
| 1089 |
JetCollection *topJets = &goodJets;
|
| 1090 |
if (topJets->size() <= 1 && hadTopJet != 0)
|
| 1091 |
topJets = &goodJets;
|
| 1092 |
leptonFromW = lepton;
|
| 1093 |
selectedNeutrino = 0;
|
| 1094 |
currentSelectedNeutrino = 0;
|
| 1095 |
reconstructNeutrinos();
|
| 1096 |
|
| 1097 |
for (unsigned short jet1Index = 0; jet1Index < topJets->size(); ++jet1Index) {
|
| 1098 |
if (hadTopJet == topJets->at(jet1Index) ||
|
| 1099 |
dJetFromWp == topJets->at(jet1Index) ||
|
| 1100 |
wJet1 == topJets->at(jet1Index) ||
|
| 1101 |
wJet2 == topJets->at(jet1Index))
|
| 1102 |
continue;
|
| 1103 |
leptonicBJet = topJets->at(jet1Index);
|
| 1104 |
leptonicW1 = ParticlePointer(new Particle(*neutrino1 + *leptonFromW));
|
| 1105 |
leptonicW2 = ParticlePointer(new Particle(*neutrino2 + *leptonFromW));
|
| 1106 |
leptonicTop1 = ParticlePointer(new Particle(*leptonicBJet + *leptonicW1));
|
| 1107 |
leptonicTop2 = ParticlePointer(new Particle(*leptonicBJet + *leptonicW2));
|
| 1108 |
selectNeutrinoSolution();
|
| 1109 |
double ltopchi2 = TopPairEventCandidate::getLeptonicChi2(currentSelectedNeutrino);
|
| 1110 |
if (ltopchi2 < bestLepChi2) {
|
| 1111 |
selectedNeutrino = currentSelectedNeutrino;
|
| 1112 |
leptonicBIndex = jet1Index;
|
| 1113 |
bestLepChi2 = ltopchi2;
|
| 1114 |
}
|
| 1115 |
}
|
| 1116 |
leptonicBJet = topJets->at(leptonicBIndex);
|
| 1117 |
leptonicW1 = ParticlePointer(new Particle(*neutrino1 + *leptonFromW));
|
| 1118 |
leptonicW2 = ParticlePointer(new Particle(*neutrino2 + *leptonFromW));
|
| 1119 |
leptonicTop1 = ParticlePointer(new Particle(*leptonicBJet + *leptonicW1));
|
| 1120 |
leptonicTop2 = ParticlePointer(new Particle(*leptonicBJet + *leptonicW2));
|
| 1121 |
doneReconstructiontop = true;
|
| 1122 |
doneReconstruction = true;
|
| 1123 |
}
|
| 1124 |
|
| 1125 |
|
| 1126 |
void TopPlusXCandidates::findDJet(const JetPointer topJet)
|
| 1127 |
{
|
| 1128 |
dJetFromWpIndex = -1;
|
| 1129 |
// dJetFromWp = 0; // unnecessary
|
| 1130 |
unsigned short cnt = 0;
|
| 1131 |
do {
|
| 1132 |
for (unsigned short jet1Index = 0; jet1Index < goodJets.size(); ++jet1Index)
|
| 1133 |
{
|
| 1134 |
if ((goodJets.at(jet1Index)->isBJet() == false || cnt >= 1) &&
|
| 1135 |
goodJets.at(jet1Index) != topJet &&
|
| 1136 |
jet1FromW != goodJets.at(jet1Index) &&
|
| 1137 |
jet2FromW != goodJets.at(jet1Index))
|
| 1138 |
{
|
| 1139 |
dJetFromWpIndex = jet1Index;
|
| 1140 |
dJetFromWp = goodJets.at(jet1Index);
|
| 1141 |
break;
|
| 1142 |
}
|
| 1143 |
}
|
| 1144 |
++cnt;
|
| 1145 |
} while (dJetFromWpIndex == -1 && cnt < 2);
|
| 1146 |
if (dJetFromWpIndex == -1)
|
| 1147 |
throw ReconstructionException("No d jet found -- all jets used for tops.");
|
| 1148 |
}
|
| 1149 |
|
| 1150 |
|
| 1151 |
void TopPlusXCandidates::recoHadronicTop(LeptonPointer lepton,
|
| 1152 |
const JetPointer lepTopJet)
|
| 1153 |
{
|
| 1154 |
if (goodJets.size() < 4)
|
| 1155 |
throw ReconstructionException("Not enough jets available to reconstruct top using Chi2 method.");
|
| 1156 |
double bestHadChi2(9999999.);
|
| 1157 |
// JetCollection *topJets = &goodBJets;
|
| 1158 |
JetCollection *topJets = &goodJets;
|
| 1159 |
if (topJets->size() <= 1 && lepTopJet != 0)
|
| 1160 |
topJets = &goodJets;
|
| 1161 |
|
| 1162 |
for (unsigned short hadBindex = 0; hadBindex < topJets->size(); ++hadBindex) {
|
| 1163 |
if (lepTopJet == topJets->at(hadBindex) ||
|
| 1164 |
dJetFromWp == topJets->at(hadBindex))
|
| 1165 |
continue;
|
| 1166 |
for (unsigned short jet1Index = 0; jet1Index < goodJets.size(); ++jet1Index) {
|
| 1167 |
if (lepTopJet == goodJets.at(jet1Index) ||
|
| 1168 |
goodJets.at(jet1Index) == topJets->at(hadBindex) ||
|
| 1169 |
jet1Index == dJetFromWpIndex)
|
| 1170 |
continue;
|
| 1171 |
for (unsigned short jet2Index = 0; jet2Index < goodJets.size();
|
| 1172 |
++jet2Index) {
|
| 1173 |
if (jet2Index == jet1Index ||
|
| 1174 |
goodJets.at(jet2Index) == topJets->at(hadBindex) ||
|
| 1175 |
jet2Index == dJetFromWpIndex ||
|
| 1176 |
lepTopJet == goodJets.at(jet2Index))
|
| 1177 |
continue;
|
| 1178 |
hadronicBJet = topJets->at(hadBindex);
|
| 1179 |
jet1FromW = goodJets.at(jet1Index);
|
| 1180 |
jet2FromW = goodJets.at(jet2Index);
|
| 1181 |
hadronicW = ParticlePointer(new Particle(*jet1FromW + *jet2FromW));
|
| 1182 |
hadronicTop = ParticlePointer(new Particle(*hadronicBJet + *hadronicW));
|
| 1183 |
double htopchi2 = getHadronicChi2();
|
| 1184 |
if (htopchi2 < bestHadChi2) {
|
| 1185 |
hadronicBIndex = hadBindex;
|
| 1186 |
jet1FromWIndex = jet1Index;
|
| 1187 |
jet2FromWIndex = jet2Index;
|
| 1188 |
bestHadChi2 = htopchi2;
|
| 1189 |
}
|
| 1190 |
}
|
| 1191 |
}
|
| 1192 |
}
|
| 1193 |
hadronicBJet = topJets->at(hadronicBIndex);
|
| 1194 |
jet1FromW = goodJets.at(jet1FromWIndex);
|
| 1195 |
jet2FromW = goodJets.at(jet2FromWIndex);
|
| 1196 |
hadronicW = ParticlePointer(new Particle(*jet1FromW + *jet2FromW));
|
| 1197 |
hadronicTop = ParticlePointer(new Particle(*hadronicBJet + *hadronicW));
|
| 1198 |
doneReconstructiontop = true;
|
| 1199 |
doneReconstruction = true;
|
| 1200 |
}
|
| 1201 |
|
| 1202 |
|
| 1203 |
TwoNonResTops::TwoNonResTops() :
|
| 1204 |
wpGenParticlesSet(false)
|
| 1205 |
{
|
| 1206 |
}
|
| 1207 |
|
| 1208 |
|
| 1209 |
TwoNonResTops::TwoNonResTops(const Event& event) :
|
| 1210 |
TopPlusXCandidates(event),
|
| 1211 |
wpGenParticlesSet(false)
|
| 1212 |
{
|
| 1213 |
}
|
| 1214 |
|
| 1215 |
|
| 1216 |
TwoNonResTops::~TwoNonResTops()
|
| 1217 |
{
|
| 1218 |
}
|
| 1219 |
|
| 1220 |
double TwoNonResTops::getGlobalChi2() const { // No global chi2 for this class
|
| 1221 |
return (0.0);
|
| 1222 |
}
|
| 1223 |
|
| 1224 |
|
| 1225 |
void TwoNonResTops::recoNonResTops(LeptonPointer lepton) {
|
| 1226 |
if (goodJets.size() < 4)
|
| 1227 |
throw ReconstructionException("Not enough jets available to reconstruct two tops");
|
| 1228 |
// const JetCollection &topJets = goodBJets;
|
| 1229 |
const JetCollection &topJets = goodJets;
|
| 1230 |
leptonFromW = lepton;
|
| 1231 |
selectedNeutrino = 0;
|
| 1232 |
currentSelectedNeutrino = 0;
|
| 1233 |
reconstructNeutrinos();
|
| 1234 |
double chosen_TopMassDifference(9999999.);
|
| 1235 |
double chosen_Chi2Total(9999999.);
|
| 1236 |
|
| 1237 |
for (unsigned short hadBindex = 0; hadBindex < topJets.size(); ++hadBindex) {
|
| 1238 |
for (unsigned short lepBindex = 0; lepBindex < topJets.size(); ++lepBindex) {
|
| 1239 |
if (lepBindex == hadBindex)
|
| 1240 |
continue;
|
| 1241 |
for (unsigned short jet1Index = 0; jet1Index < goodJets.size(); ++jet1Index) {
|
| 1242 |
if (goodJets.at(jet1Index) == topJets.at(lepBindex) ||
|
| 1243 |
goodJets.at(jet1Index) == topJets.at(hadBindex))
|
| 1244 |
continue;
|
| 1245 |
for (unsigned short jet2Index = 0; jet2Index < goodJets.size(); ++jet2Index) {
|
| 1246 |
if (jet2Index == jet1Index ||
|
| 1247 |
goodJets.at(jet2Index) == topJets.at(lepBindex) ||
|
| 1248 |
goodJets.at(jet2Index) == topJets.at(hadBindex))
|
| 1249 |
continue;
|
| 1250 |
hadronicBJet = topJets.at(hadBindex);
|
| 1251 |
leptonicBJet = topJets.at(lepBindex);
|
| 1252 |
jet1FromW = goodJets.at(jet1Index);
|
| 1253 |
jet2FromW = goodJets.at(jet2Index);
|
| 1254 |
|
| 1255 |
leptonicW1 = ParticlePointer(new Particle(*neutrino1 + *leptonFromW));
|
| 1256 |
leptonicW2 = ParticlePointer(new Particle(*neutrino2 + *leptonFromW));
|
| 1257 |
hadronicW = ParticlePointer(new Particle(*jet1FromW + *jet2FromW));
|
| 1258 |
leptonicTop1 = ParticlePointer(new Particle(*leptonicBJet + *leptonicW1));
|
| 1259 |
leptonicTop2 = ParticlePointer(new Particle(*leptonicBJet + *leptonicW2));
|
| 1260 |
hadronicTop = ParticlePointer(new Particle(*hadronicBJet + *hadronicW));
|
| 1261 |
fillHypotheses();
|
| 1262 |
selectNeutrinoSolution();
|
| 1263 |
double TopMassDifference = calculateTopMassDifference(currentSelectedNeutrino);
|
| 1264 |
double chi2 = getTotalChi2(currentSelectedNeutrino);
|
| 1265 |
switch (usedTTbarReconstruction) {
|
| 1266 |
case TTbarReconstructionCriterion::TopMassDifference:
|
| 1267 |
if (TopMassDifference < chosen_TopMassDifference) {
|
| 1268 |
hadronicBIndex = hadBindex;
|
| 1269 |
leptonicBIndex = lepBindex;
|
| 1270 |
jet1FromWIndex = jet1Index;
|
| 1271 |
jet2FromWIndex = jet2Index;
|
| 1272 |
chosen_TopMassDifference = TopMassDifference;
|
| 1273 |
selectedNeutrino = currentSelectedNeutrino;
|
| 1274 |
}
|
| 1275 |
break;
|
| 1276 |
case TTbarReconstructionCriterion::chi2:
|
| 1277 |
if (chi2 < chosen_Chi2Total) {
|
| 1278 |
hadronicBIndex = hadBindex;
|
| 1279 |
leptonicBIndex = lepBindex;
|
| 1280 |
jet1FromWIndex = jet1Index;
|
| 1281 |
jet2FromWIndex = jet2Index;
|
| 1282 |
chosen_Chi2Total = chi2;
|
| 1283 |
selectedNeutrino = currentSelectedNeutrino;
|
| 1284 |
}
|
| 1285 |
break;
|
| 1286 |
}
|
| 1287 |
}
|
| 1288 |
}
|
| 1289 |
}
|
| 1290 |
}
|
| 1291 |
std::sort(solutions.begin(), solutions.end(), compareSolutions);
|
| 1292 |
hadronicBJet = topJets.at(hadronicBIndex);
|
| 1293 |
leptonicBJet = topJets.at(leptonicBIndex);
|
| 1294 |
jet1FromW = goodJets.at(jet1FromWIndex);
|
| 1295 |
jet2FromW = goodJets.at(jet2FromWIndex);
|
| 1296 |
leptonicW1 = ParticlePointer(new Particle(*neutrino1 + *leptonFromW));
|
| 1297 |
leptonicW2 = ParticlePointer(new Particle(*neutrino2 + *leptonFromW));
|
| 1298 |
hadronicW = ParticlePointer(new Particle(*jet1FromW + *jet2FromW));
|
| 1299 |
leptonicTop1 = ParticlePointer(new Particle(*leptonicBJet + *leptonicW1));
|
| 1300 |
leptonicTop2 = ParticlePointer(new Particle(*leptonicBJet + *leptonicW2));
|
| 1301 |
hadronicTop = ParticlePointer(new Particle(*hadronicBJet + *hadronicW));
|
| 1302 |
if (selectedNeutrino == 1) // Resonance may be real of fake -- can exclude ttbar?
|
| 1303 |
ttbarResonance = ParticlePointer(new Particle(*leptonicTop1 + *hadronicTop));
|
| 1304 |
else
|
| 1305 |
ttbarResonance = ParticlePointer(new Particle(*leptonicTop2 + *hadronicTop));
|
| 1306 |
doneReconstruction = true;
|
| 1307 |
}
|
| 1308 |
|
| 1309 |
|
| 1310 |
double TwoNonResTops::getTotalChi2DJ(unsigned short neutrinoSolution) const
|
| 1311 |
{
|
| 1312 |
double chi2 = TopPairEventCandidate::getTotalChi2(neutrinoSolution);
|
| 1313 |
double pt = dJetFromWp->pt();
|
| 1314 |
if (pt < jet1FromW->pt())
|
| 1315 |
chi2 += (0.5/sqrt(2.0));
|
| 1316 |
if (pt < jet2FromW->pt())
|
| 1317 |
chi2 += (0.5/sqrt(2.0));
|
| 1318 |
double deltaphi = dJetFromWp->deltaPhi(hadronicBJet);
|
| 1319 |
double deltaphi2 = dJetFromWp->deltaPhi(leptonicBJet);
|
| 1320 |
if (deltaphi2 > deltaphi)
|
| 1321 |
deltaphi = deltaphi2;
|
| 1322 |
double dphichi = deltaphi - pi;
|
| 1323 |
chi2 += (1.0/sqrt(2.0)) * (dphichi * dphichi)/(pi * pi/2.0);
|
| 1324 |
if (dJetFromWp->isBJet())
|
| 1325 |
chi2 += (0.5/sqrt(2.0));
|
| 1326 |
return (chi2);
|
| 1327 |
}
|
| 1328 |
|
| 1329 |
|
| 1330 |
void TwoNonResTops::recoNonResTopsWdjet(LeptonPointer lepton) {
|
| 1331 |
if (goodJets.size() < 5)
|
| 1332 |
throw ReconstructionException("Not enough jets available to reconstruct two tops");
|
| 1333 |
const JetCollection &topJets = goodJets;
|
| 1334 |
leptonFromW = lepton;
|
| 1335 |
selectedNeutrino = 0;
|
| 1336 |
currentSelectedNeutrino = 0;
|
| 1337 |
reconstructNeutrinos();
|
| 1338 |
double chosen_TopMassDifference(9999999.);
|
| 1339 |
double chosen_Chi2Total(9999999.);
|
| 1340 |
|
| 1341 |
for (unsigned short hadBindex = 0; hadBindex < topJets.size(); ++hadBindex) {
|
| 1342 |
for (unsigned short lepBindex = 0; lepBindex < topJets.size(); ++lepBindex) {
|
| 1343 |
if (lepBindex == hadBindex)
|
| 1344 |
continue;
|
| 1345 |
for (unsigned short jet1Index = 0; jet1Index < goodJets.size(); ++jet1Index) {
|
| 1346 |
if (goodJets.at(jet1Index) == topJets.at(lepBindex) ||
|
| 1347 |
goodJets.at(jet1Index) == topJets.at(hadBindex))
|
| 1348 |
continue;
|
| 1349 |
for (unsigned short jet2Index = 0; jet2Index < goodJets.size(); ++jet2Index) {
|
| 1350 |
if (jet2Index == jet1Index ||
|
| 1351 |
goodJets.at(jet2Index) == topJets.at(lepBindex) ||
|
| 1352 |
goodJets.at(jet2Index) == topJets.at(hadBindex))
|
| 1353 |
continue;
|
| 1354 |
for (unsigned short djetIndex = 0; djetIndex < goodJets.size(); ++djetIndex) {
|
| 1355 |
if (djetIndex == jet1Index || jet2Index == djetIndex ||
|
| 1356 |
goodJets.at(djetIndex) == topJets.at(lepBindex) ||
|
| 1357 |
goodJets.at(djetIndex) == topJets.at(hadBindex))
|
| 1358 |
continue;
|
| 1359 |
hadronicBJet = topJets.at(hadBindex);
|
| 1360 |
leptonicBJet = topJets.at(lepBindex);
|
| 1361 |
jet1FromW = goodJets.at(jet1Index);
|
| 1362 |
jet2FromW = goodJets.at(jet2Index);
|
| 1363 |
dJetFromWp = goodJets.at(djetIndex);
|
| 1364 |
|
| 1365 |
leptonicW1 = ParticlePointer(new Particle(*neutrino1 + *leptonFromW));
|
| 1366 |
leptonicW2 = ParticlePointer(new Particle(*neutrino2 + *leptonFromW));
|
| 1367 |
hadronicW = ParticlePointer(new Particle(*jet1FromW + *jet2FromW));
|
| 1368 |
leptonicTop1 = ParticlePointer(new Particle(*leptonicBJet + *leptonicW1));
|
| 1369 |
leptonicTop2 = ParticlePointer(new Particle(*leptonicBJet + *leptonicW2));
|
| 1370 |
hadronicTop = ParticlePointer(new Particle(*hadronicBJet + *hadronicW));
|
| 1371 |
fillHypotheses();
|
| 1372 |
selectNeutrinoSolution();
|
| 1373 |
double TopMassDifference = calculateTopMassDifference(currentSelectedNeutrino);
|
| 1374 |
double chi2 = getTotalChi2DJ(currentSelectedNeutrino);
|
| 1375 |
switch (usedTTbarReconstruction) {
|
| 1376 |
case TTbarReconstructionCriterion::TopMassDifference:
|
| 1377 |
if (TopMassDifference < chosen_TopMassDifference) {
|
| 1378 |
hadronicBIndex = hadBindex;
|
| 1379 |
leptonicBIndex = lepBindex;
|
| 1380 |
jet1FromWIndex = jet1Index;
|
| 1381 |
jet2FromWIndex = jet2Index;
|
| 1382 |
chosen_TopMassDifference = TopMassDifference;
|
| 1383 |
selectedNeutrino = currentSelectedNeutrino;
|
| 1384 |
}
|
| 1385 |
break;
|
| 1386 |
case TTbarReconstructionCriterion::chi2:
|
| 1387 |
if (chi2 < chosen_Chi2Total) {
|
| 1388 |
hadronicBIndex = hadBindex;
|
| 1389 |
leptonicBIndex = lepBindex;
|
| 1390 |
jet1FromWIndex = jet1Index;
|
| 1391 |
jet2FromWIndex = jet2Index;
|
| 1392 |
dJetFromWpIndex = djetIndex;
|
| 1393 |
chosen_Chi2Total = chi2;
|
| 1394 |
selectedNeutrino = currentSelectedNeutrino;
|
| 1395 |
}
|
| 1396 |
break;
|
| 1397 |
}
|
| 1398 |
}
|
| 1399 |
}
|
| 1400 |
}
|
| 1401 |
}
|
| 1402 |
}
|
| 1403 |
std::sort(solutions.begin(), solutions.end(), compareSolutions);
|
| 1404 |
hadronicBJet = topJets.at(hadronicBIndex);
|
| 1405 |
leptonicBJet = topJets.at(leptonicBIndex);
|
| 1406 |
jet1FromW = goodJets.at(jet1FromWIndex);
|
| 1407 |
jet2FromW = goodJets.at(jet2FromWIndex);
|
| 1408 |
dJetFromWp = goodJets.at(dJetFromWpIndex);
|
| 1409 |
leptonicW1 = ParticlePointer(new Particle(*neutrino1 + *leptonFromW));
|
| 1410 |
leptonicW2 = ParticlePointer(new Particle(*neutrino2 + *leptonFromW));
|
| 1411 |
hadronicW = ParticlePointer(new Particle(*jet1FromW + *jet2FromW));
|
| 1412 |
leptonicTop1 = ParticlePointer(new Particle(*leptonicBJet + *leptonicW1));
|
| 1413 |
leptonicTop2 = ParticlePointer(new Particle(*leptonicBJet + *leptonicW2));
|
| 1414 |
hadronicTop = ParticlePointer(new Particle(*hadronicBJet + *hadronicW));
|
| 1415 |
if (selectedNeutrino == 1) // Resonance may be real of fake -- can exclude ttbar?
|
| 1416 |
ttbarResonance = ParticlePointer(new Particle(*leptonicTop1 + *hadronicTop));
|
| 1417 |
else
|
| 1418 |
ttbarResonance = ParticlePointer(new Particle(*leptonicTop2 + *hadronicTop));
|
| 1419 |
doneReconstruction = true;
|
| 1420 |
}
|
| 1421 |
|
| 1422 |
|
| 1423 |
const JetPointer TwoNonResTops::getLtNonTopJet() const {
|
| 1424 |
throwExpeptionIfNotReconstructed("TwoNonResTops::getLtNonTopJet");
|
| 1425 |
if (goodJets.size() < 5)
|
| 1426 |
throw ReconstructionException("Not enough jets available to reconstruct Wprime");
|
| 1427 |
int ltJetInd = -1, bJetInd = -1;
|
| 1428 |
unsigned int index = 0;
|
| 1429 |
while (index < goodJets.size() && ltJetInd < 0) {
|
| 1430 |
JetPointer currJet = goodJets.at(index);
|
| 1431 |
if (currJet->isBJet() == false && currJet != leptonicBJet &&
|
| 1432 |
currJet != hadronicBJet && currJet != jet1FromW && currJet != jet2FromW) {
|
| 1433 |
ltJetInd = index;
|
| 1434 |
} else if (bJetInd < 0 && currJet != leptonicBJet && currJet != hadronicBJet &&
|
| 1435 |
currJet != jet1FromW && currJet != jet2FromW) {
|
| 1436 |
bJetInd = index;
|
| 1437 |
}
|
| 1438 |
++index;
|
| 1439 |
}
|
| 1440 |
if (ltJetInd < 0 && bJetInd >= 0) // Take a B jet if no others found.
|
| 1441 |
ltJetInd = bJetInd;
|
| 1442 |
if (ltJetInd >= 0)
|
| 1443 |
return (goodJets.at(ltJetInd));
|
| 1444 |
throw ReconstructionException(TString("getLtNonTopJet:Fifth jet missing"));
|
| 1445 |
}
|
| 1446 |
|
| 1447 |
|
| 1448 |
const ParticlePointer TwoNonResTops::getWprime() const {
|
| 1449 |
throwExpeptionIfNotReconstructed("TwoNonResTops::getWprime");
|
| 1450 |
ParticlePointer ltJet = getLtNonTopJet();
|
| 1451 |
ParticlePointer hiPtTop = getLeptonicTop();
|
| 1452 |
|
| 1453 |
// Generator study indicates lower Pt top is usually spectator
|
| 1454 |
if (hiPtTop->pt() < hadronicTop->pt())
|
| 1455 |
hiPtTop = hadronicTop;
|
| 1456 |
return (ParticlePointer(new Particle(*hiPtTop + *ltJet)));
|
| 1457 |
}
|
| 1458 |
|
| 1459 |
|
| 1460 |
const ParticlePointer TwoNonResTops::getNegWprime() const {
|
| 1461 |
throwExpeptionIfNotReconstructed("TwoNonResTops::getNegWprime");
|
| 1462 |
ParticlePointer ltJet = getLtNonTopJet();
|
| 1463 |
ParticlePointer negTop = hadronicTop;
|
| 1464 |
if (leptonFromW->charge() < 0)
|
| 1465 |
negTop = getLeptonicTop();
|
| 1466 |
|
| 1467 |
return (ParticlePointer(new Particle(*negTop + *ltJet)));
|
| 1468 |
// return (ParticlePointer(new Particle(*negTop + *dJetFromWp)));
|
| 1469 |
}
|
| 1470 |
|
| 1471 |
|
| 1472 |
// Deliberately reco the incorrect W' to see if ttbar still looks the same
|
| 1473 |
const ParticlePointer TwoNonResTops::getNegWprimeBad() const {
|
| 1474 |
throwExpeptionIfNotReconstructed("TwoNonResTops::getNegWprimeBad");
|
| 1475 |
ParticlePointer ltJet = getLtNonTopJet();
|
| 1476 |
ParticlePointer negTop = hadronicTop;
|
| 1477 |
if (leptonFromW->charge() > 0)
|
| 1478 |
negTop = getLeptonicTop(); // Use wrong top
|
| 1479 |
|
| 1480 |
return (ParticlePointer(new Particle(*negTop + *ltJet)));
|
| 1481 |
// return (ParticlePointer(new Particle(*negTop + *dJetFromWp)));
|
| 1482 |
}
|
| 1483 |
|
| 1484 |
|
| 1485 |
const ParticlePointer TwoNonResTops::getBtBWprime() const {
|
| 1486 |
throwExpeptionIfNotReconstructed("TwoNonResTops::getBtBWprim");
|
| 1487 |
ParticlePointer ltJet = getLtNonTopJet();
|
| 1488 |
ParticlePointer hiAnglTop = getLeptonicTop();
|
| 1489 |
if (getHadTopDAngle() > hiAnglTop->angle(ltJet))
|
| 1490 |
hiAnglTop = hadronicTop;
|
| 1491 |
return (ParticlePointer(new Particle(*hiAnglTop + *ltJet)));
|
| 1492 |
}
|
| 1493 |
|
| 1494 |
|
| 1495 |
double TwoNonResTops::getHadTopDAngle() const
|
| 1496 |
{
|
| 1497 |
ParticlePointer jet = getLtNonTopJet();
|
| 1498 |
return (hadronicTop->angle(jet));
|
| 1499 |
}
|
| 1500 |
|
| 1501 |
|
| 1502 |
void TwoNonResTops::recoWpFrom1Top(LeptonPointer lepton) {
|
| 1503 |
if (goodJets.size() < 5)
|
| 1504 |
throw
|
| 1505 |
ReconstructionException("Not enough jets available to reconstruct W'");
|
| 1506 |
ParticlePointer negTop, posTop;
|
| 1507 |
if (lepton->charge() < 0) {
|
| 1508 |
// findDJet();
|
| 1509 |
recoBestLeptoTop(lepton);
|
| 1510 |
// findDJet(leptonicBJet);
|
| 1511 |
recoHadronicTop(lepton);
|
| 1512 |
// recoHadronicTop(lepton, leptonicBJet);
|
| 1513 |
negTop = getLeptonicTop();
|
| 1514 |
posTop = hadronicTop;
|
| 1515 |
} else {
|
| 1516 |
// findDJet();
|
| 1517 |
recoHadronicTop(lepton);
|
| 1518 |
// findDJet(hadronicBJet);
|
| 1519 |
recoBestLeptoTop(lepton);
|
| 1520 |
// recoBestLeptoTop(lepton, hadronicBJet, jet1FromW, jet2FromW);
|
| 1521 |
negTop = hadronicTop;
|
| 1522 |
posTop = getLeptonicTop();
|
| 1523 |
}
|
| 1524 |
JetPointer ltJet = getLtNonTopJet();
|
| 1525 |
wpFrom1Top = ParticlePointer(new Particle(*negTop + *ltJet));
|
| 1526 |
wpBadFrom1Top = ParticlePointer(new Particle(*posTop + *ltJet));
|
| 1527 |
topFromWp = negTop;
|
| 1528 |
dJetFromWp = ltJet;
|
| 1529 |
// wpFrom1Top = ParticlePointer(new Particle(*negTop + *dJetFromWp));
|
| 1530 |
// wpBadFrom1Top = ParticlePointer(new Particle(*posTop + *dJetFromWp));
|
| 1531 |
}
|
| 1532 |
|
| 1533 |
|
| 1534 |
void TwoNonResTops::getTrueWpChain()
|
| 1535 |
{
|
| 1536 |
const MCParticleCollection &mcColl = genParticles;
|
| 1537 |
const int nullInd = -99;
|
| 1538 |
int wpInd = nullInd, wpTopInd = nullInd, specTopInd = nullInd, dQInd = nullInd;
|
| 1539 |
int wpTopBInd = nullInd, specTopBInd = nullInd;
|
| 1540 |
if (mcColl.size() > 0) {
|
| 1541 |
for (unsigned int i = 0; i < mcColl.size(); ++i) {
|
| 1542 |
if (fabs(mcColl.at(i)->pdgId()) == WPRIME_ID) {
|
| 1543 |
if (wpInd == nullInd)
|
| 1544 |
wpInd = i;
|
| 1545 |
else if (mcColl.at(i)->motherIndex() != wpInd)
|
| 1546 |
cout << "Stray W': " << i << " mother " << mcColl.at(i)->motherIndex() << endl;
|
| 1547 |
}
|
| 1548 |
}
|
| 1549 |
for (unsigned int i = 0; i < mcColl.size(); ++i) {
|
| 1550 |
if (fabs(mcColl.at(i)->pdgId()) == TOP_ID) {
|
| 1551 |
if (mcColl.at(i)->motherIndex() == wpInd) {
|
| 1552 |
if (wpTopInd == nullInd)
|
| 1553 |
wpTopInd = i;
|
| 1554 |
else cout << "Additional W' top: " << i << endl;
|
| 1555 |
} else if (specTopInd == nullInd)
|
| 1556 |
specTopInd = i;
|
| 1557 |
else cout << "Additional spec top: " << i << endl;
|
| 1558 |
} else if (fabs(mcColl.at(i)->pdgId()) == DQ_ID) {
|
| 1559 |
if (mcColl.at(i)->motherIndex() == wpInd) {
|
| 1560 |
if (dQInd == nullInd)
|
| 1561 |
dQInd = i;
|
| 1562 |
else cout << "Additional d daughter: " << i << endl;
|
| 1563 |
}
|
| 1564 |
}
|
| 1565 |
}
|
| 1566 |
if (wpTopInd >= 0 && specTopInd >= 0) {
|
| 1567 |
int wptprevPdgId = 0;
|
| 1568 |
int sptprevPdgId = 0;
|
| 1569 |
for (unsigned int i = 0; i < mcColl.size(); ++i) {
|
| 1570 |
int pdgId = fabs(mcColl.at(i)->pdgId());
|
| 1571 |
if (mcColl.at(i)->motherIndex() == wpTopInd ||
|
| 1572 |
mcColl.at(i)->motherIndex() == wpInd) {
|
| 1573 |
if (mcColl.at(i)->isQuark() && pdgId % 2 == 1) {
|
| 1574 |
if (wpTopBInd == nullInd || pdgId > wptprevPdgId) {
|
| 1575 |
wpTopBInd = i;
|
| 1576 |
wptprevPdgId = pdgId;
|
| 1577 |
} else if (pdgId == BQ_ID)
|
| 1578 |
cout << "Additional W' top b quark: " << i << " mother " << mcColl.at(i)->motherIndex() << endl;
|
| 1579 |
}
|
| 1580 |
} else if (mcColl.at(i)->motherIndex() == specTopInd ||
|
| 1581 |
mcColl.at(i)->motherIndex() == mcColl.at(specTopInd)->motherIndex()) {
|
| 1582 |
if (mcColl.at(i)->isQuark() && pdgId % 2 == 1) {
|
| 1583 |
if (specTopBInd == nullInd || pdgId > sptprevPdgId) {
|
| 1584 |
specTopBInd = i;
|
| 1585 |
sptprevPdgId = pdgId;
|
| 1586 |
} else if (pdgId == BQ_ID)
|
| 1587 |
cout << "Additional spec top b quark: " << i << " mother " << mcColl.at(i)->motherIndex() << endl;
|
| 1588 |
}
|
| 1589 |
} else if (pdgId == BQ_ID && fabs(mcColl.at(mcColl.at(i)->motherIndex())->pdgId()) != BQ_ID) {
|
| 1590 |
// cout << "Additional b quark: " << i << " mother " << mcColl.at(i)->motherIndex();
|
| 1591 |
// cout << " ID " << mcColl.at(mcColl.at(i)->motherIndex())->pdgId() << endl;
|
| 1592 |
}
|
| 1593 |
}
|
| 1594 |
}
|
| 1595 |
if (wpInd >= 0 && wpTopInd >= 0 && specTopInd >= 0 && dQInd >= 0 &&
|
| 1596 |
wpTopBInd >= 0 && specTopBInd >= 0) {
|
| 1597 |
wpTru = mcColl.at(wpInd);
|
| 1598 |
wpTopTru = mcColl.at(wpTopInd);
|
| 1599 |
specTopTru = mcColl.at(specTopInd);
|
| 1600 |
wpTopBTru = mcColl.at(wpTopBInd);
|
| 1601 |
specTopBTru = mcColl.at(specTopBInd);
|
| 1602 |
wpDDauTru = mcColl.at(dQInd);
|
| 1603 |
wpGenParticlesSet = true;
|
| 1604 |
} else {
|
| 1605 |
// cout << "#11 : " << mcColl.at(11)->pdgId() << " mother " << mcColl.at(11)->motherIndex();
|
| 1606 |
// cout << " ID " << mcColl.at(mcColl.at(11)->motherIndex())->pdgId() << endl;
|
| 1607 |
cout << "wp wpTop specTop dq wpTopB specTopB " << wpInd << " " << wpTopInd << " " << specTopInd << " " << dQInd;
|
| 1608 |
cout << " " << wpTopBInd << " " << specTopBInd << endl;
|
| 1609 |
}
|
| 1610 |
}
|
| 1611 |
}
|
| 1612 |
|
| 1613 |
|
| 1614 |
MCParticlePointer TwoNonResTops::getWpTru() const
|
| 1615 |
{
|
| 1616 |
if (wpGenParticlesSet == false)
|
| 1617 |
throw ReconstructionException("Need to run TwoNonResTops::getTrueWpChain first");
|
| 1618 |
return (wpTru);
|
| 1619 |
}
|
| 1620 |
|
| 1621 |
MCParticlePointer TwoNonResTops::getWpTopTru() const
|
| 1622 |
{
|
| 1623 |
if (wpGenParticlesSet == false)
|
| 1624 |
throw ReconstructionException("Need to run TwoNonResTops::getTrueWpChain first");
|
| 1625 |
return (wpTopTru);
|
| 1626 |
}
|
| 1627 |
|
| 1628 |
MCParticlePointer TwoNonResTops::getSpecTop() const
|
| 1629 |
{
|
| 1630 |
if (wpGenParticlesSet == false)
|
| 1631 |
throw ReconstructionException("Need to run TwoNonResTops::getTrueWpChain first");
|
| 1632 |
return (specTopTru);
|
| 1633 |
}
|
| 1634 |
|
| 1635 |
MCParticlePointer TwoNonResTops::getWpTopBTru() const
|
| 1636 |
{
|
| 1637 |
if (wpGenParticlesSet == false)
|
| 1638 |
throw ReconstructionException("Need to run TwoNonResTops::getTrueWpChain first");
|
| 1639 |
return (wpTopBTru);
|
| 1640 |
}
|
| 1641 |
|
| 1642 |
MCParticlePointer TwoNonResTops::getSpecTopB() const
|
| 1643 |
{
|
| 1644 |
if (wpGenParticlesSet == false)
|
| 1645 |
throw ReconstructionException("Need to run TwoNonResTops::getTrueWpChain first");
|
| 1646 |
return (specTopBTru);
|
| 1647 |
}
|
| 1648 |
|
| 1649 |
MCParticlePointer TwoNonResTops::getWpDDauTru() const
|
| 1650 |
{
|
| 1651 |
if (wpGenParticlesSet == false)
|
| 1652 |
throw ReconstructionException("Need to run TwoNonResTops::getTrueWpChain first");
|
| 1653 |
return (wpDDauTru);
|
| 1654 |
}
|
| 1655 |
|
| 1656 |
|
| 1657 |
TopNoMassConstraint::TopNoMassConstraint()
|
| 1658 |
{
|
| 1659 |
}
|
| 1660 |
|
| 1661 |
|
| 1662 |
TopNoMassConstraint::TopNoMassConstraint(const Event& event) :
|
| 1663 |
ToplikeCandidate(event)
|
| 1664 |
{
|
| 1665 |
}
|
| 1666 |
|
| 1667 |
|
| 1668 |
TopNoMassConstraint::~TopNoMassConstraint()
|
| 1669 |
{
|
| 1670 |
}
|
| 1671 |
|
| 1672 |
|
| 1673 |
double TopNoMassConstraint::getLoneHadChi2() const
|
| 1674 |
{
|
| 1675 |
double ptRatioDifference = TMath::Power(PtRatio() - matched_ptratio, 2);
|
| 1676 |
double ptRatioError = 2 * matched_ptratio_sigma * matched_ptratio_sigma;
|
| 1677 |
double ptRatioTerm = ptRatioDifference / ptRatioError;
|
| 1678 |
|
| 1679 |
double WmassDifference = TMath::Power(hadronicW->mass() - pdgWmass, 2);
|
| 1680 |
double WmassError = 2 * wMassWidth * wMassWidth;
|
| 1681 |
double WmassTerm = WmassDifference / WmassError;
|
| 1682 |
return 1.0 / sqrt(2.0) * (WmassTerm + ptRatioTerm);
|
| 1683 |
}
|
| 1684 |
|
| 1685 |
double TopNoMassConstraint::getLepBAngle() const
|
| 1686 |
{
|
| 1687 |
return (leptonicBJet->angle(leptonFromW));
|
| 1688 |
}
|
| 1689 |
|
| 1690 |
|
| 1691 |
LeptoTopNoMassConstraint::LeptoTopNoMassConstraint()
|
| 1692 |
{
|
| 1693 |
}
|
| 1694 |
|
| 1695 |
|
| 1696 |
LeptoTopNoMassConstraint::LeptoTopNoMassConstraint(const Event& event) :
|
| 1697 |
TopNoMassConstraint(event)
|
| 1698 |
{
|
| 1699 |
}
|
| 1700 |
|
| 1701 |
|
| 1702 |
LeptoTopNoMassConstraint::~LeptoTopNoMassConstraint()
|
| 1703 |
{
|
| 1704 |
}
|
| 1705 |
|
| 1706 |
|
| 1707 |
/*
|
| 1708 |
Tried for neutrino selection, but getGlobalChi2 still uses hadronic top quantities.
|
| 1709 |
double LeptoTopNoMassConstraint::getTotalChi2(unsigned short int neutrinoSolution) const {
|
| 1710 |
return TopPairEventCandidate::getLeptonicChi2(neutrinoSolution) + TopPairEventCandidate::getGlobalChi2(neutrinoSolution);
|
| 1711 |
}
|
| 1712 |
*/
|
| 1713 |
|
| 1714 |
void LeptoTopNoMassConstraint::inspectReconstructedEvent() const {
|
| 1715 |
cout << "run " << runNumber << ", event " << eventNumber << endl;
|
| 1716 |
cout << "leptonic b jet, goodJet index " << leptonicBIndex << endl;
|
| 1717 |
EventContentPrinter::printJet(leptonicBJet);
|
| 1718 |
|
| 1719 |
cout << "lepton from W" << endl;
|
| 1720 |
// EventContentPrinter::printElectron(electronFromW);
|
| 1721 |
|
| 1722 |
cout << "MET" << endl;
|
| 1723 |
EventContentPrinter::printParticle(met);
|
| 1724 |
cout << endl;
|
| 1725 |
|
| 1726 |
cout << "reconstructed neutrino 1(selected: " << selectedNeutrino << ")" << endl;
|
| 1727 |
EventContentPrinter::printParticle(neutrino1);
|
| 1728 |
cout << endl;
|
| 1729 |
|
| 1730 |
cout << "reconstructed neutrino 2(selected: " << selectedNeutrino << ")" << endl;
|
| 1731 |
EventContentPrinter::printParticle(neutrino2);
|
| 1732 |
cout << endl;
|
| 1733 |
|
| 1734 |
cout << "leptonic W 1 (selected: " << selectedNeutrino << ")" << endl;
|
| 1735 |
EventContentPrinter::printParticle(leptonicW1);
|
| 1736 |
cout << endl;
|
| 1737 |
|
| 1738 |
cout << "leptonic W 2 (selected: " << selectedNeutrino << ")" << endl;
|
| 1739 |
EventContentPrinter::printParticle(leptonicW2);
|
| 1740 |
cout << endl;
|
| 1741 |
|
| 1742 |
cout << "leptonic top 1 (selected: " << selectedNeutrino << ")" << endl;
|
| 1743 |
EventContentPrinter::printParticle(leptonicTop1);
|
| 1744 |
cout << endl;
|
| 1745 |
|
| 1746 |
cout << "leptonic top 2 (selected: " << selectedNeutrino << ")" << endl;
|
| 1747 |
EventContentPrinter::printParticle(leptonicTop2);
|
| 1748 |
cout << endl;
|
| 1749 |
|
| 1750 |
cout << endl;
|
| 1751 |
}
|
| 1752 |
|
| 1753 |
|
| 1754 |
}
|
