]>
Commit | Line | Data |
---|---|---|
5a31e9df | 1 | #include "AliD0Trigger.h" |
2 | #include <TMath.h> | |
3 | #include "AliITStrackV2.h" | |
4 | #include "AliL3Transform.h" | |
5 | #include <TVector3.h> | |
0bd0c1ef | 6 | #include <iostream.h> |
5a31e9df | 7 | |
8 | ClassImp(AliD0Trigger) | |
9 | ||
10 | AliD0Trigger::AliD0Trigger() | |
11 | { | |
12 | posTrack=0; | |
13 | negTrack=0; | |
14 | } | |
0bd0c1ef | 15 | AliD0Trigger::AliD0Trigger(double c[6],double b,double pv[3]) |
5a31e9df | 16 | { |
17 | posTrack=0; | |
18 | negTrack=0; | |
19 | cutV0low=c[0]; | |
20 | cutV0high=c[1]; | |
21 | cutInvMass=c[2]; | |
22 | cutPointAngle=c[3]; | |
0bd0c1ef | 23 | cutd0d0=c[4]; |
24 | cutCosThetaStar=c[5]; | |
25 | cutpTchild=c[6]; | |
5a31e9df | 26 | Bfield=b; |
27 | primaryVertex[0]=pv[0]; | |
28 | primaryVertex[1]=pv[1]; | |
29 | primaryVertex[2]=pv[2]; | |
30 | } | |
31 | AliD0Trigger::AliD0Trigger(AliITStrackV2 * posT, AliITStrackV2 * negT) | |
32 | { | |
33 | posTrack=posT; | |
34 | negTrack=negT; | |
35 | } | |
36 | ||
37 | ||
38 | AliD0Trigger::~AliD0Trigger() | |
39 | { | |
40 | ||
41 | } | |
42 | ||
43 | void AliD0Trigger::SetTracks(AliITStrackV2 * posT, AliITStrackV2 * negT){ | |
44 | posTrack=posT; | |
45 | negTrack=negT; | |
46 | } | |
47 | bool AliD0Trigger::FindInvMass() { | |
48 | ||
49 | double mD0 = 1.8645; | |
50 | double mPi = 0.13957; | |
51 | double mKa = 0.49368; | |
52 | ||
53 | double mom2[2],momTot2,energy[2]; | |
54 | ||
55 | mom2[0]=pow(momenta[0],2)+pow(momenta[1],2)+pow(momenta[2],2); | |
56 | mom2[1]=pow(momenta[3],2)+pow(momenta[4],2)+pow(momenta[5],2); | |
57 | ||
58 | momTot2 = pow(momenta[0]+momenta[3],2)+pow(momenta[1]+momenta[4],2)+pow(momenta[2]+momenta[5],2); | |
59 | ||
60 | // D0 -> K- Pi+ | |
61 | energy[1] = sqrt(pow(mKa,2)+mom2[1]); | |
62 | energy[0] = sqrt(pow(mPi,2)+mom2[0]); | |
63 | ||
64 | double minvD0 = sqrt((energy[0]+energy[1])*(energy[0]+energy[1])-momTot2); | |
65 | ||
66 | // D0bar -> K+ Pi- | |
67 | energy[0] = sqrt(mKa*mKa+mom2[0]); | |
68 | energy[1] = sqrt(mPi*mPi+mom2[1]); | |
69 | ||
70 | double minvD0bar = sqrt((energy[0]+energy[1])*(energy[0]+energy[1])-momTot2); | |
71 | ||
72 | if(fabs(minvD0-mD0) < cutInvMass) return true; | |
73 | if(fabs(minvD0bar-mD0) < cutInvMass) return true; | |
74 | return false; | |
75 | ||
76 | } | |
77 | bool AliD0Trigger::FindV0(){ | |
78 | ||
0bd0c1ef | 79 | bool goodV0=false; |
80 | ||
81 | double Gxpos=0,Gypos=0,Gxneg=0,Gyneg=0; | |
82 | double r1=0,r2=0,a1=0,a2=0,b1=0,b2=0; | |
83 | double q=0,p=0,t=0,Fa=0,Fb=0,Fc=0; | |
84 | double y1=0,y2=0,x1=0,x2=0,x3=0,x4=0; | |
85 | double V01=0,V02=0,V03=0,V04=0; | |
86 | double V0[4]={0,0,0,0}; | |
87 | bestV0[0]=0; bestV0[1]=0; bestV0[2]=0; | |
5a31e9df | 88 | |
0bd0c1ef | 89 | Gxpos=posTrack->GetX()*cos(posTrack->GetAlpha())-posTrack->GetY()*sin(posTrack->GetAlpha()); |
90 | Gypos=posTrack->GetX()*sin(posTrack->GetAlpha())+posTrack->GetY()*cos(posTrack->GetAlpha()); | |
91 | Gxneg=negTrack->GetX()*cos(negTrack->GetAlpha())-negTrack->GetY()*sin(negTrack->GetAlpha()); | |
92 | Gyneg=negTrack->GetX()*sin(negTrack->GetAlpha())+negTrack->GetY()*cos(negTrack->GetAlpha()); | |
5a31e9df | 93 | |
0bd0c1ef | 94 | r1=fabs(1/(AliL3Transform::GetBFact()*Bfield*posTrack->Get1Pt())); |
95 | r2=fabs(1/(AliL3Transform::GetBFact()*Bfield*negTrack->Get1Pt())); | |
96 | ||
97 | a1=Gxpos-(r1*cos(posTrack->GetAlpha()+asin(posTrack->GetSnp())+TMath::PiOver2())); | |
98 | a2=Gxneg-(r2*cos(negTrack->GetAlpha()+asin(negTrack->GetSnp())-TMath::PiOver2())); | |
99 | b1=Gypos-(r1*sin(posTrack->GetAlpha()+asin(posTrack->GetSnp())+TMath::PiOver2())); | |
100 | b2=Gyneg-(r2*sin(negTrack->GetAlpha()+asin(negTrack->GetSnp())-TMath::PiOver2())); | |
101 | ||
102 | //double a3=a1-a2; | |
103 | //double b3=b1-b2; | |
104 | //double r3=sqrt(pow(a3,2)+pow(b3,2)); | |
105 | ||
106 | //if(r3<r1+r2){ | |
107 | q=a1-a2; | |
108 | p=2*(b1-b2); | |
109 | t=pow(r2,2)-pow(r2,2)+pow(b1,2)-pow(b2,2)-pow(q,2); | |
110 | ||
111 | Fa=pow(p,2)+(4*pow(q,2)); | |
112 | Fb=-(2*t*p+8*pow(q,2)*b1); | |
113 | Fc=pow(t,2)-(4*pow(q,2)*pow(r1,2))+(4*pow(q,2)*pow(b1,2)); | |
114 | ||
115 | if(pow(Fb,2)-(4*Fa*Fc)>=0){ | |
116 | y1=(-Fb+(sqrt(pow(Fb,2)-(4*Fa*Fc))))/(2*Fa); //noe feil her. floating point | |
117 | y2=(-Fb-(sqrt(pow(Fb,2)-(4*Fa*Fc))))/(2*Fa); //trolig negativ under rot | |
118 | ||
119 | x1=sqrt(pow(r1,2)-pow((y1-b1),2))+a1; | |
120 | x2=sqrt(pow(r1,2)-pow((y2-b1),2))+a1; | |
121 | x3=-sqrt(pow(r1,2)-pow((y1-b1),2))+a1; | |
122 | x4=-sqrt(pow(r1,2)-pow((y2-b1),2))+a1; | |
123 | ||
124 | V01=sqrt(pow(x1-primaryVertex[0],2)+pow(y1-primaryVertex[1],2)); | |
125 | V02=sqrt(pow(x2-primaryVertex[0],2)+pow(y2-primaryVertex[1],2)); | |
126 | V03=sqrt(pow(x3-primaryVertex[0],2)+pow(y1-primaryVertex[1],2)); | |
127 | V04=sqrt(pow(x4-primaryVertex[0],2)+pow(y2-primaryVertex[1],2)); | |
128 | ||
129 | V0[0]=V01;V0[1]=V02;V0[2]=V03;V0[3]=V04; | |
130 | ||
131 | int index=0; | |
132 | ||
133 | double nearestV0=V0[0]; | |
134 | ||
135 | for(int i=0;i<4;i++){ | |
136 | if(nearestV0>V0[i]){ | |
137 | nearestV0=V0[i]; | |
138 | index=i; | |
5a31e9df | 139 | } |
0bd0c1ef | 140 | } |
141 | cout<<"index: "<<index<<endl; | |
142 | if(nearestV0<cutV0high){ | |
143 | if(nearestV0>cutV0low){ | |
144 | goodV0=true; | |
145 | if(index==0){ | |
146 | bestV0[0]=x1; bestV0[1]=y1; bestV0[2]=0.; | |
147 | } | |
148 | if(index==1){ | |
149 | bestV0[0]=x2; bestV0[1]=y2; bestV0[2]=0.; | |
150 | } | |
151 | if(index==2){ | |
152 | bestV0[0]=x3; bestV0[1]=y1; bestV0[2]=0.; | |
153 | } | |
154 | if(index==3){ | |
155 | bestV0[0]=x4; bestV0[1]=y2; bestV0[2]=0.; | |
5a31e9df | 156 | } |
157 | } | |
158 | } | |
0bd0c1ef | 159 | //} |
5a31e9df | 160 | } |
0bd0c1ef | 161 | |
162 | cout<<"My V0: x: "<<bestV0[0]<<" y: "<<bestV0[1]<<" z: "<<bestV0[2]<<endl; | |
163 | ||
5a31e9df | 164 | return goodV0; |
165 | } | |
166 | void AliD0Trigger::FindMomentaAtVertex(){ | |
167 | ||
168 | double r1=fabs(1/(AliL3Transform::GetBFact()*Bfield*posTrack->Get1Pt())); | |
169 | double r2=fabs(1/(AliL3Transform::GetBFact()*Bfield*negTrack->Get1Pt())); | |
170 | ||
171 | double Gx1=posTrack->GetX()*cos(posTrack->GetAlpha())-posTrack->GetY()*sin(posTrack->GetAlpha()); | |
172 | double Gy1=posTrack->GetX()*sin(posTrack->GetAlpha())+posTrack->GetY()*cos(posTrack->GetAlpha()); | |
173 | double Gx2=negTrack->GetX()*cos(negTrack->GetAlpha())-negTrack->GetY()*sin(negTrack->GetAlpha()); | |
174 | double Gy2=negTrack->GetX()*sin(negTrack->GetAlpha())+negTrack->GetY()*cos(negTrack->GetAlpha()); | |
175 | ||
176 | double centerx1=Gx1-(r1*cos(posTrack->GetAlpha()+asin(posTrack->GetSnp())+TMath::PiOver2())); | |
177 | double centerx2=Gx2-(r2*cos(negTrack->GetAlpha()+asin(negTrack->GetSnp())-TMath::PiOver2())); | |
178 | double centery1=Gy1-(r1*sin(posTrack->GetAlpha()+asin(posTrack->GetSnp())+TMath::PiOver2())); | |
179 | double centery2=Gy2-(r2*sin(negTrack->GetAlpha()+asin(negTrack->GetSnp())-TMath::PiOver2())); | |
180 | ||
181 | double a1=sqrt(pow(Gx1-bestV0[0],2)+pow(Gy1-bestV0[1],2)); | |
182 | double b1=sqrt(pow(bestV0[0]-centerx1,2)+pow(bestV0[1]-centery1,2)); | |
183 | double c1=sqrt(pow(Gx1-centerx1,2)+pow(Gy1-centery1,2)); | |
184 | double a2=sqrt(pow(Gx2-bestV0[0],2)+pow(Gy2-bestV0[1],2)); | |
185 | double b2=sqrt(pow(bestV0[0]-centerx2,2)+pow(bestV0[1]-centery2,2)); | |
186 | double c2=sqrt(pow(Gx2-centerx2,2)+pow(Gy2-centery2,2)); | |
187 | ||
188 | double alpha1=acos((pow(b1,2)+pow(c1,2)-pow(a1,2))/(2*b1*c1)); | |
189 | double alpha2=acos((pow(b2,2)+pow(c2,2)-pow(a2,2))/(2*b2*c2)); | |
190 | ||
191 | double sign1=0; | |
192 | double sign2=0; | |
193 | ||
194 | if((sqrt(pow(Gx1,2)+pow(Gy1,2)))>(sqrt(pow(bestV0[0],2)+pow(bestV0[1],2)))) | |
195 | sign1=-1; | |
196 | else | |
197 | sign1=1; | |
198 | if((sqrt(pow(Gx2,2)+pow(Gy2,2)))>(sqrt(pow(bestV0[0],2)+pow(bestV0[1],2)))) | |
199 | sign2=1; | |
200 | else | |
201 | sign2=-1; | |
202 | ||
203 | double psi1=(posTrack->GetAlpha()+asin(posTrack->GetSnp()))+(sign1*alpha1); | |
204 | double psi2=(negTrack->GetAlpha()+asin(negTrack->GetSnp()))+(sign2*alpha2); | |
205 | ||
206 | double ptP = 1./fabs(posTrack->Get1Pt()); | |
207 | momenta[0] = ptP*cos(psi1); | |
208 | momenta[1] = ptP*sin(psi1); | |
209 | momenta[2] = ptP*posTrack->GetTgl(); | |
210 | ||
211 | double ptN = 1./fabs(negTrack->Get1Pt()); | |
212 | momenta[3] = ptN*cos(psi2); | |
213 | momenta[4] = ptN*sin(psi2); | |
214 | momenta[5] = ptN*negTrack->GetTgl(); | |
215 | ||
216 | } | |
217 | Bool_t AliD0Trigger::PointingAngle(){ | |
218 | ||
219 | TVector3 mom(momenta[0]+momenta[3],momenta[1]+momenta[4],momenta[2]+momenta[5]); | |
220 | TVector3 flight(bestV0[0]-primaryVertex[0],bestV0[1]-primaryVertex[1],bestV0[2]-primaryVertex[2]); | |
221 | ||
222 | double pta = mom.Angle(flight); | |
223 | ||
224 | if(cos(pta)<cutPointAngle) | |
225 | return true; | |
226 | else | |
227 | return false; | |
228 | } | |
0bd0c1ef | 229 | void AliD0Trigger::SetMomenta(double m[6]) |
230 | { | |
231 | momenta[0] = m[0]; | |
232 | momenta[1] = m[1]; | |
233 | momenta[2] = m[2]; | |
234 | momenta[3] = m[3]; | |
235 | momenta[4] = m[4]; | |
236 | momenta[5] = m[5]; | |
237 | } | |
238 | void AliD0Trigger::FindMomentaOffline() | |
239 | { | |
240 | double ptP,alphaP,phiP,ptN,alphaN,phiN; | |
241 | // momenta of the tracks at the vertex | |
242 | ptP = 1./TMath::Abs(posTrack->Get1Pt()); | |
243 | alphaP = posTrack->GetAlpha(); | |
244 | phiP = alphaP+TMath::ASin(posTrack->GetSnp()); | |
245 | momenta[0] = ptP*TMath::Cos(phiP); | |
246 | momenta[1] = ptP*TMath::Sin(phiP); | |
247 | momenta[2] = ptP*posTrack->GetTgl(); | |
248 | ||
249 | ptN = 1./TMath::Abs(negTrack->Get1Pt()); | |
250 | alphaN = negTrack->GetAlpha(); | |
251 | phiN = alphaN+TMath::ASin(negTrack->GetSnp()); | |
252 | momenta[3] = ptN*TMath::Cos(phiN); | |
253 | momenta[4] = ptN*TMath::Sin(phiN); | |
254 | momenta[5] = ptN*negTrack->GetTgl(); | |
255 | } | |
256 | bool AliD0Trigger::FindV0offline(double v[3]) | |
257 | { | |
258 | ||
259 | bool goodV0=false; | |
260 | double r=0; | |
261 | ||
262 | bestV0[0]=v[0]; | |
263 | bestV0[1]=v[1]; | |
264 | bestV0[2]=v[2]; | |
265 | ||
266 | r=sqrt(pow(bestV0[0]-primaryVertex[0],2)+pow(bestV0[1]-primaryVertex[1],2)+pow(bestV0[2]-primaryVertex[2],2)); | |
267 | ||
268 | if(r<cutV0high){ | |
269 | if(r>cutV0low){ | |
270 | goodV0=true; | |
271 | } | |
272 | } | |
273 | ||
274 | return goodV0; | |
275 | ||
276 | } | |
277 | bool AliD0Trigger::d0d0() | |
278 | { | |
279 | bool goodd0=false; | |
280 | double d00=0, d01=0; | |
281 | ||
282 | d00 = 10000.*posTrack->GetD(primaryVertex[0],primaryVertex[1]); | |
283 | d01 = -10000.*negTrack->GetD(primaryVertex[0],primaryVertex[1]); | |
284 | ||
285 | if(d00*d01<cutd0d0){ | |
286 | goodd0=true; | |
287 | } | |
288 | ||
289 | return goodd0; | |
290 | } | |
291 | double AliD0Trigger::Energy() | |
292 | { | |
293 | double kMD0 = 1.8645; // D0 mass | |
294 | return sqrt(P()*P()+kMD0*kMD0); | |
295 | } | |
296 | bool AliD0Trigger::CosThetaStar() | |
297 | { | |
298 | bool goodtheta=false; | |
299 | double kMD0 = 1.8645; // D0 mass | |
300 | double kMK = 0.49368; // K+ mass | |
301 | double kMPi = 0.13957; // Pi+ mass | |
302 | double qL=0,qL2=0,ctsD0=0,ctsD0bar=0; | |
303 | ||
304 | Double_t pStar = TMath::Sqrt(TMath::Power(kMD0*kMD0-kMK*kMK-kMPi*kMPi,2.)-4.*kMK*kMK*kMPi*kMPi)/(2.*kMD0); | |
305 | ||
306 | Double_t beta = P()/Energy(); | |
307 | Double_t gamma = Energy()/kMD0; | |
308 | ||
309 | TVector3 mom(momenta[0],momenta[1],momenta[2]); | |
310 | TVector3 mom2(momenta[3],momenta[4],momenta[5]); | |
311 | TVector3 momD(Px(),Py(),Pz()); | |
312 | ||
313 | qL = mom.Dot(momD)/momD.Mag(); | |
314 | ||
315 | ctsD0 = (qL/gamma-beta*TMath::Sqrt(pStar*pStar+kMK*kMK))/pStar; | |
316 | ||
317 | qL2 = mom.Dot(momD)/momD.Mag(); | |
318 | ||
319 | ctsD0bar = (qL2/gamma-beta*TMath::Sqrt(pStar*pStar+kMK*kMK))/pStar; | |
320 | ||
321 | if(TMath::Abs(ctsD0) > cutCosThetaStar){ | |
322 | if(TMath::Abs(ctsD0bar) > cutCosThetaStar){ | |
323 | goodtheta=true; | |
324 | } | |
325 | } | |
326 | return goodtheta; | |
327 | } | |
328 | bool AliD0Trigger::pTchild() | |
329 | { | |
330 | bool goodpT=false; | |
331 | double pT1=0,pT2=0; | |
332 | ||
333 | pT1=sqrt(pow(momenta[0],2)+pow(momenta[1],2)); | |
334 | pT2=sqrt(pow(momenta[3],2)+pow(momenta[4],2)); | |
335 | ||
336 | if(pT1>cutpTchild){ | |
337 | if(pT2>cutpTchild){ | |
338 | goodpT=true; | |
339 | } | |
340 | } | |
341 | return goodpT; | |
342 | } |