ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/UHHAnalysis/SFrameTools/src/Utils.cxx
(Generate patch)

Comparing UserCode/UHHAnalysis/SFrameTools/src/Utils.cxx (file contents):
Revision 1.5 by peiffer, Mon Jul 9 08:37:32 2012 UTC vs.
Revision 1.6 by hoeing, Mon Oct 15 15:35:51 2012 UTC

# Line 1 | Line 1
1   #include "../include/Utils.h"
2  
3  
4 + //variable HEP Tagger from Rebekka
5 +
6 + bool variableHepTopTag(TopJet topjet, double ptJetMin, double massWindowLower, double massWindowUpper, double cutCondition2, double cutCondition3){
7 +
8 +  double mjet;
9 +  double ptjet;
10 +  int nsubjets;
11 +
12 +  double topmass=172.3;
13 +  double wmass=80.4;
14 +
15 +  nsubjets=topjet.numberOfDaughters();
16 +
17 +  LorentzVector allsubjets(0,0,0,0);
18 +  
19 +  for(int j=0; j<topjet.numberOfDaughters(); ++j){
20 +    allsubjets += topjet.subjets()[j].v4();
21 +  }
22 +  if(!allsubjets.isTimelike()){
23 +    mjet=0;
24 +    return false;
25 +  }
26 +
27 +  mjet = allsubjets.M();
28 +  ptjet= allsubjets.Pt();
29 +
30 +  double m12, m13, m23;
31 +
32 +  //The subjetcs have to be three
33 +  if(nsubjets==3) {
34 +
35 +    std::vector<Particle> subjets = topjet.subjets();
36 +    sort(subjets.begin(), subjets.end(), HigherPt());
37 +
38 +    m12 = 0;
39 +    if( (subjets[0].v4()+subjets[1].v4()).isTimelike())
40 +      m12=(subjets[0].v4()+subjets[1].v4()).M();
41 +    m13 = 0;
42 +    if( (subjets[0].v4()+subjets[2].v4()).isTimelike() )
43 +      m13=(subjets[0].v4()+subjets[2].v4()).M();
44 +    m23 = 0;
45 +    if( (subjets[1].v4()+subjets[2].v4()).isTimelike()  )
46 +      m23 = (subjets[1].v4()+subjets[2].v4()).M();
47 +  
48 +  }
49 +  else{
50 +    return false;
51 +  }
52 +
53 +  double rmin=massWindowLower*wmass/topmass;
54 +  double rmax=massWindowUpper*wmass/topmass;
55 +
56 +  int keep=0;
57 +
58 +  //Conditions on the subjects: at least one has to be true
59 +  //1 condition
60 +  if(atan(m13/m12)>0.2 && atan(m13/m12)<1.3 && m23/mjet>rmin && m23/mjet<rmax) keep=1;
61 +  
62 +  //2 condition
63 +  double cond2left=pow(rmin,2)*(1+pow((m13/m12),2));
64 +  double cond2cent=1-pow(m23/mjet,2);
65 +  double cond2right=pow(rmax,2)*(1+pow(m13/m12,2));
66 +
67 +  if(cond2left<cond2cent && cond2cent<cond2right && m23/mjet>cutCondition2) keep=1;
68 +
69 +  //3 condition
70 +  double cond3left=pow(rmin,2)*(1+pow((m12/m13),2));
71 +  double cond3cent=1-pow(m23/mjet,2);
72 +  double cond3right=pow(rmax,2)*(1+pow(m12/m13,2));
73 +
74 +  if(cond3left<cond3cent && cond3cent<cond3right && m23/mjet>cutCondition3) keep=1;
75 +
76 +  //Final requirement: at least one of the three subjets conditions and total pt
77 +  if(keep==1 && ptjet>ptJetMin){
78 +    return true;
79 +  }
80 +  else{
81 +    return false;
82 +  }
83 +
84 + }
85 +
86 +
87 +
88 +
89 +
90 +
91 + //HEP Tagger from Ivan
92 +
93 + bool HepTopTag(TopJet topjet){
94 +
95 +  double mjet;
96 +  double ptjet;
97 +  int nsubjets;
98 +
99 +  double topmass=172.3;
100 +  double wmass=80.4;
101 +
102 +  nsubjets=topjet.numberOfDaughters();
103 +
104 +  LorentzVector allsubjets(0,0,0,0);
105 +  
106 +  for(int j=0; j<topjet.numberOfDaughters(); ++j){
107 +    allsubjets += topjet.subjets()[j].v4();
108 +  }
109 +  if(!allsubjets.isTimelike()){
110 +    mjet=0;
111 +    return false;
112 +  }
113 +
114 +  mjet = allsubjets.M();
115 +  ptjet= allsubjets.Pt();
116 +
117 +  double m12, m13, m23;
118 +
119 +  //The subjetcs have to be three
120 +  if(nsubjets==3) {
121 +
122 +    std::vector<Particle> subjets = topjet.subjets();
123 +    sort(subjets.begin(), subjets.end(), HigherPt());
124 +
125 +    m12 = 0;
126 +    if( (subjets[0].v4()+subjets[1].v4()).isTimelike())
127 +      m12=(subjets[0].v4()+subjets[1].v4()).M();
128 +    m13 = 0;
129 +    if( (subjets[0].v4()+subjets[2].v4()).isTimelike() )
130 +      m13=(subjets[0].v4()+subjets[2].v4()).M();
131 +    m23 = 0;
132 +    if( (subjets[1].v4()+subjets[2].v4()).isTimelike()  )
133 +      m23 = (subjets[1].v4()+subjets[2].v4()).M();
134 +  
135 +  }
136 +  else{
137 +    return false;
138 +  }
139 +
140 +  double rmin=0.85*wmass/topmass;
141 +  double rmax=1.15*wmass/topmass;
142 +
143 +  int keep=0;
144 +
145 +  //Conditions on the subjects: at least one has to be true
146 +  //1 condition
147 +  if(atan(m13/m12)>0.2 && atan(m13/m12)<1.3 && m23/mjet>rmin && m23/mjet<rmax) keep=1;
148 +  
149 +  //2 condition
150 +  double cond2left=pow(rmin,2)*(1+pow((m13/m12),2));
151 +  double cond2cent=1-pow(m23/mjet,2);
152 +  double cond2right=pow(rmax,2)*(1+pow(m13/m12,2));
153 +
154 +  if(cond2left<cond2cent && cond2cent<cond2right && m23/mjet>0.35) keep=1;
155 +
156 +  //3 condition
157 +  double cond3left=pow(rmin,2)*(1+pow((m12/m13),2));
158 +  double cond3cent=1-pow(m23/mjet,2);
159 +  double cond3right=pow(rmax,2)*(1+pow(m12/m13,2));
160 +
161 +  if(cond3left<cond3cent && cond3cent<cond3right && m23/mjet>0.35) keep=1;
162 +
163 +  //Final requirement: at least one of the three subjets conditions and total pt
164 +  if(keep==1 && ptjet>200){
165 +    return true;
166 +  }
167 +  else{
168 +    return false;
169 +  }
170 +
171 + }
172 +
173 +
174 + //default values (mminLower=50., mjetLower=140, mjetUpper=250.) defined in Utils.h
175 + bool variableTopTag(TopJet topjet, double &mjet, int &nsubjets, double &mmin, double mminLower, double mjetLower, double mjetUpper){
176 +
177 +  nsubjets=topjet.numberOfDaughters();
178 +
179 +  LorentzVector allsubjets(0,0,0,0);
180 +
181 +  for(int j=0; j<topjet.numberOfDaughters(); ++j){
182 +  allsubjets += topjet.subjets()[j].v4();
183 +  }
184 +  if(!allsubjets.isTimelike()){
185 +  mjet=0;
186 +  mmin=0;
187 + //  mminLower=50;
188 + //   mjetLower=140;
189 + //   mjetUpper=250;
190 +  return false;
191 +  }
192 +
193 +  mjet = allsubjets.M();
194 +
195 +  if(nsubjets>=3) {
196 +
197 +  std::vector<Particle> subjets = topjet.subjets();
198 +  sort(subjets.begin(), subjets.end(), HigherPt());
199 +
200 +  double m01 = 0;
201 +  if( (subjets[0].v4()+subjets[1].v4()).isTimelike())
202 +  m01=(subjets[0].v4()+subjets[1].v4()).M();
203 +  double m02 = 0;
204 +  if( (subjets[0].v4()+subjets[2].v4()).isTimelike() )
205 +  m02=(subjets[0].v4()+subjets[2].v4()).M();
206 +  double m12 = 0;
207 +  if( (subjets[1].v4()+subjets[2].v4()).isTimelike() )
208 +  m12 = (subjets[1].v4()+subjets[2].v4()).M();
209 +
210 +  //minimum pairwise mass
211 +  mmin = std::min(m01,std::min(m02,m12));
212 +  }
213 +
214 +  //at least 3 sub-jets
215 +  if(nsubjets<3) return false;
216 +  //minimum pairwise mass > 50 GeV/c^2
217 +  if(mmin<mminLower) return false;
218 +  //jet mass between 140 and 250 GeV/c^2
219 +  if(mjet<mjetLower || mjet>mjetUpper) return false;
220 +
221 +  return true;
222 + }
223 +
224 +
225   bool TopTag(TopJet topjet,  double &mjet, int &nsubjets, double &mmin){
226  
227    nsubjets=topjet.numberOfDaughters();

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines