| 1 |
peiffer |
1.1 |
#include "../include/Utils.h"
|
| 2 |
|
|
|
| 3 |
|
|
|
| 4 |
|
|
bool TopTag(TopJet topjet, double &mjet, int &nsubjets, double &mmin){
|
| 5 |
|
|
|
| 6 |
|
|
nsubjets=topjet.numberOfDaughters();
|
| 7 |
|
|
|
| 8 |
|
|
LorentzVector allsubjets(0,0,0,0);
|
| 9 |
|
|
|
| 10 |
|
|
for(int j=0; j<topjet.numberOfDaughters(); ++j){
|
| 11 |
|
|
allsubjets += topjet.subjets()[j].v4();
|
| 12 |
|
|
}
|
| 13 |
|
|
if(!allsubjets.isTimelike()){
|
| 14 |
|
|
mjet=0;
|
| 15 |
|
|
mmin=0;
|
| 16 |
|
|
return false;
|
| 17 |
|
|
}
|
| 18 |
|
|
|
| 19 |
|
|
mjet = allsubjets.M();
|
| 20 |
|
|
|
| 21 |
|
|
if(nsubjets>=3) {
|
| 22 |
|
|
|
| 23 |
|
|
std::vector<Particle> subjets = topjet.subjets();
|
| 24 |
|
|
sort(subjets.begin(), subjets.end(), HigherPt());
|
| 25 |
|
|
|
| 26 |
|
|
double m01 = 0;
|
| 27 |
|
|
if( (subjets[0].v4()+subjets[1].v4()).isTimelike())
|
| 28 |
|
|
m01=(subjets[0].v4()+subjets[1].v4()).M();
|
| 29 |
|
|
double m02 = 0;
|
| 30 |
|
|
if( (subjets[0].v4()+subjets[2].v4()).isTimelike() )
|
| 31 |
|
|
m02=(subjets[0].v4()+subjets[2].v4()).M();
|
| 32 |
|
|
double m12 = 0;
|
| 33 |
|
|
if( (subjets[1].v4()+subjets[2].v4()).isTimelike() )
|
| 34 |
|
|
m12 = (subjets[1].v4()+subjets[2].v4()).M();
|
| 35 |
|
|
|
| 36 |
|
|
//minimum pairwise mass > 50 GeV/c^2
|
| 37 |
|
|
mmin = std::min(m01,std::min(m02,m12));
|
| 38 |
|
|
}
|
| 39 |
|
|
|
| 40 |
|
|
//at least 3 sub-jets
|
| 41 |
|
|
if(nsubjets<3) return false;
|
| 42 |
|
|
//minimum pairwise mass > 50 GeV/c^2
|
| 43 |
|
|
if(mmin<50) return false;
|
| 44 |
|
|
//jet mass between 140 and 250 GeV/c^2
|
| 45 |
|
|
if(mjet<140 || mjet>250) return false;
|
| 46 |
|
|
|
| 47 |
|
|
return true;
|
| 48 |
|
|
}
|
| 49 |
peiffer |
1.2 |
|
| 50 |
|
|
Jet* nextJet(const Particle *p, std::vector<Jet> *jets){
|
| 51 |
|
|
|
| 52 |
|
|
double deltarmin = double_infinity();
|
| 53 |
|
|
Jet* nextjet=0;
|
| 54 |
|
|
for(unsigned int i=0; i<jets->size();++i){
|
| 55 |
|
|
if(jets->at(i).deltaR(*p) < deltarmin){
|
| 56 |
|
|
deltarmin = jets->at(i).deltaR(*p);
|
| 57 |
|
|
nextjet = &jets->at(i);
|
| 58 |
|
|
}
|
| 59 |
|
|
}
|
| 60 |
|
|
|
| 61 |
|
|
return nextjet;
|
| 62 |
|
|
}
|
| 63 |
|
|
|
| 64 |
|
|
double pTrel(const Particle *p, std::vector<Jet> *jets){
|
| 65 |
|
|
|
| 66 |
|
|
double ptrel=0;
|
| 67 |
|
|
|
| 68 |
|
|
Jet* nextjet = nextJet(p,jets);
|
| 69 |
|
|
|
| 70 |
|
|
TVector3 p3(p->v4().Px(),p->v4().Py(),p->v4().Pz());
|
| 71 |
|
|
TVector3 jet3(nextjet->v4().Px(),nextjet->v4().Py(),nextjet->v4().Pz());
|
| 72 |
|
|
|
| 73 |
|
|
if(p3.Mag()!=0 && jet3.Mag()!=0){
|
| 74 |
|
|
double sin_alpha = (p3.Cross(jet3)).Mag()/p3.Mag()/jet3.Mag();
|
| 75 |
|
|
ptrel = p3.Mag()*sin_alpha;
|
| 76 |
|
|
}
|
| 77 |
|
|
else{
|
| 78 |
|
|
std::cout << "something strange happend in the ptrel calculation: either lepton or jet momentum is 0" <<std::endl;
|
| 79 |
|
|
}
|
| 80 |
|
|
|
| 81 |
|
|
return ptrel;
|
| 82 |
|
|
}
|
| 83 |
|
|
|
| 84 |
|
|
double deltaRmin(const Particle *p, std::vector<Jet> *jets){
|
| 85 |
|
|
return nextJet(p,jets)->deltaR(*p);
|
| 86 |
|
|
}
|
| 87 |
|
|
|
| 88 |
|
|
|
| 89 |
|
|
double double_infinity(){
|
| 90 |
|
|
return std::numeric_limits<double>::infinity() ;
|
| 91 |
|
|
}
|
| 92 |
|
|
|
| 93 |
|
|
int int_infinity(){
|
| 94 |
|
|
return std::numeric_limits<int>::max() ;
|
| 95 |
|
|
}
|
