]>
Commit | Line | Data |
---|---|---|
c92725b7 | 1 | /************************************************************************** |
2 | * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. * | |
3 | * * | |
4 | * Author: The ALICE Off-line Project. * | |
5 | * Contributors are mentioned in the code where appropriate. * | |
6 | * * | |
7 | * Permission to use, copy, modify and distribute this software and its * | |
8 | * documentation strictly for non-commercial purposes is hereby granted * | |
9 | * without fee, provided that the above copyright notice appears in all * | |
10 | * copies and that both the copyright notice and this permission notice * | |
11 | * appear in the supporting documentation. The authors make no claims * | |
12 | * about the suitability of this software for any purpose. It is * | |
13 | * provided "as is" without express or implied warranty. * | |
14 | **************************************************************************/ | |
15 | ||
16 | ||
17 | /////////////////////////////////////////////////////////////////////////// | |
18 | /* | |
19 | ||
20 | Origin: marian.ivanov@cern.ch | |
d92975ba | 21 | Container classes with MC infomation |
c92725b7 | 22 | |
23 | */ | |
24 | ||
25 | #if !defined(__CINT__) || defined(__MAKECINT__) | |
26 | #include <stdio.h> | |
27 | #include <string.h> | |
28 | //ROOT includes | |
29 | #include "TROOT.h" | |
30 | #include "Rtypes.h" | |
31 | #include "TFile.h" | |
32 | #include "TTree.h" | |
33 | #include "TChain.h" | |
34 | #include "TCut.h" | |
35 | #include "TString.h" | |
c92725b7 | 36 | #include "TStopwatch.h" |
37 | #include "TParticle.h" | |
38 | #include "TSystem.h" | |
c92725b7 | 39 | #include "TCanvas.h" |
c92725b7 | 40 | #include "TGeometry.h" |
41 | #include "TPolyLine3D.h" | |
42 | ||
43 | //ALIROOT includes | |
44 | #include "AliRun.h" | |
45 | #include "AliStack.h" | |
46 | #include "AliSimDigits.h" | |
47 | #include "AliTPCParam.h" | |
48 | #include "AliTPC.h" | |
49 | #include "AliTPCLoader.h" | |
50 | #include "AliDetector.h" | |
51 | #include "AliTrackReference.h" | |
52 | #include "AliTPCParamSR.h" | |
53 | #include "AliTracker.h" | |
54 | #include "AliMagF.h" | |
55 | #include "AliHelix.h" | |
c92725b7 | 56 | #include "AliTrackPointArray.h" |
57 | ||
58 | #endif | |
59 | #include "AliGenInfo.h" | |
60 | // | |
61 | // | |
62 | ||
63 | ClassImp(AliTPCdigitRow); | |
64 | ClassImp(AliMCInfo); | |
65 | ClassImp(AliGenV0Info) | |
66 | ClassImp(AliGenKinkInfo) | |
c92725b7 | 67 | |
68 | ||
69 | ||
70 | //////////////////////////////////////////////////////////////////////// | |
022044bf | 71 | AliMCInfo::AliMCInfo(): |
72 | fTrackRef(), | |
73 | fTrackRefOut(), | |
74 | fTRdecay(), | |
75 | fPrimPart(0), | |
76 | fParticle(), | |
77 | fMass(0), | |
78 | fCharge(0), | |
79 | fLabel(0), | |
80 | fEventNr(), | |
81 | fMCtracks(), | |
82 | fPdg(0), | |
83 | fTPCdecay(0), | |
84 | fRowsWithDigitsInn(0), | |
85 | fRowsWithDigits(0), | |
86 | fRowsTrackLength(0), | |
87 | fPrim(0), | |
88 | fTPCRow(), | |
89 | fNTPCRef(0), // tpc references counter | |
90 | fNITSRef(0), // ITS references counter | |
91 | fNTRDRef(0), // TRD references counter | |
92 | fNTOFRef(0), // TOF references counter | |
93 | fTPCReferences(0), | |
94 | fITSReferences(0), | |
95 | fTRDReferences(0), | |
96 | fTOFReferences(0) | |
c92725b7 | 97 | { |
98 | fTPCReferences = new TClonesArray("AliTrackReference",10); | |
99 | fITSReferences = new TClonesArray("AliTrackReference",10); | |
100 | fTRDReferences = new TClonesArray("AliTrackReference",10); | |
101 | fTOFReferences = new TClonesArray("AliTrackReference",10); | |
102 | fTRdecay.SetTrack(-1); | |
103 | fCharge = 0; | |
104 | } | |
105 | ||
022044bf | 106 | AliMCInfo::AliMCInfo(const AliMCInfo& info): |
107 | TObject(), | |
108 | fTrackRef(info.fTrackRef), | |
109 | fTrackRefOut(info.fTrackRefOut), | |
d92975ba | 110 | fTRdecay(info.GetTRdecay()), |
022044bf | 111 | fPrimPart(info.fPrimPart), |
112 | fParticle(info.fParticle), | |
d92975ba | 113 | fMass(info.GetMass()), |
022044bf | 114 | fCharge(info.fCharge), |
115 | fLabel(info.fLabel), | |
116 | fEventNr(info.fEventNr), | |
117 | fMCtracks(info.fMCtracks), | |
118 | fPdg(info.fPdg), | |
119 | fTPCdecay(info.fTPCdecay), | |
120 | fRowsWithDigitsInn(info.fRowsWithDigitsInn), | |
121 | fRowsWithDigits(info.fRowsWithDigits), | |
122 | fRowsTrackLength(info.fRowsTrackLength), | |
123 | fPrim(info.fPrim), | |
124 | fTPCRow(info.fTPCRow), | |
125 | fNTPCRef(info.fNTPCRef), // tpc references counter | |
126 | fNITSRef(info.fNITSRef), // ITS references counter | |
127 | fNTRDRef(info.fNTRDRef), // TRD references counter | |
128 | fNTOFRef(info.fNTOFRef), // TOF references counter | |
129 | fTPCReferences(0), | |
130 | fITSReferences(0), | |
131 | fTRDReferences(0), | |
132 | fTOFReferences(0) | |
133 | { | |
134 | fTPCReferences = (TClonesArray*)info.fTPCReferences->Clone(); | |
135 | fITSReferences = (TClonesArray*)info.fITSReferences->Clone(); | |
136 | fTRDReferences = (TClonesArray*)info.fTRDReferences->Clone(); | |
137 | fTOFReferences = (TClonesArray*)info.fTOFReferences->Clone(); | |
138 | } | |
139 | ||
140 | ||
c92725b7 | 141 | AliMCInfo::~AliMCInfo() |
142 | { | |
143 | if (fTPCReferences) { | |
144 | delete fTPCReferences; | |
145 | } | |
146 | if (fITSReferences){ | |
147 | delete fITSReferences; | |
148 | } | |
149 | if (fTRDReferences){ | |
150 | delete fTRDReferences; | |
151 | } | |
152 | if (fTOFReferences){ | |
153 | delete fTOFReferences; | |
154 | } | |
155 | ||
156 | } | |
157 | ||
158 | ||
159 | ||
160 | void AliMCInfo::Update() | |
161 | { | |
162 | // | |
163 | // | |
164 | fMCtracks =1; | |
165 | if (!fTPCReferences) { | |
166 | fNTPCRef =0; | |
167 | return; | |
168 | } | |
169 | Float_t direction=1; | |
170 | //Float_t rlast=0; | |
171 | fNTPCRef = fTPCReferences->GetEntriesFast(); | |
172 | fNITSRef = fITSReferences->GetEntriesFast(); | |
173 | fNTRDRef = fTRDReferences->GetEntriesFast(); | |
174 | fNTOFRef = fTOFReferences->GetEntriesFast(); | |
175 | ||
176 | for (Int_t iref =0;iref<fTPCReferences->GetEntriesFast();iref++){ | |
177 | AliTrackReference * ref = (AliTrackReference *) fTPCReferences->At(iref); | |
178 | //Float_t r = (ref->X()*ref->X()+ref->Y()*ref->Y()); | |
179 | Float_t newdirection = ref->X()*ref->Px()+ref->Y()*ref->Py(); //inside or outside | |
180 | if (iref==0) direction = newdirection; | |
181 | if ( newdirection*direction<0){ | |
182 | //changed direction | |
183 | direction = newdirection; | |
184 | fMCtracks+=1; | |
185 | } | |
186 | //rlast=r; | |
187 | } | |
188 | // | |
189 | // decay info | |
190 | fTPCdecay=kFALSE; | |
191 | if (fTRdecay.GetTrack()>0){ | |
192 | fDecayCoord[0] = fTRdecay.X(); | |
193 | fDecayCoord[1] = fTRdecay.Y(); | |
194 | fDecayCoord[2] = fTRdecay.Z(); | |
195 | if ( (fTRdecay.R()<250)&&(fTRdecay.R()>85) && (TMath::Abs(fTRdecay.Z())<250) ){ | |
196 | fTPCdecay=kTRUE; | |
197 | } | |
198 | else{ | |
199 | fDecayCoord[0] = 0; | |
200 | fDecayCoord[1] = 0; | |
201 | fDecayCoord[2] = 0; | |
202 | } | |
203 | } | |
204 | // | |
205 | // | |
206 | //digits information update | |
207 | fRowsWithDigits = fTPCRow.RowsOn(); | |
208 | fRowsWithDigitsInn = fTPCRow.RowsOn(63); // 63 = number of inner rows | |
209 | fRowsTrackLength = fTPCRow.Last() - fTPCRow.First(); | |
210 | // | |
211 | // | |
212 | // calculate primary ionization per cm | |
213 | if (fParticle.GetPDG()){ | |
214 | fMass = fParticle.GetMass(); | |
215 | fCharge = fParticle.GetPDG()->Charge(); | |
216 | if (fTPCReferences->GetEntriesFast()>0){ | |
217 | fTrackRef = *((AliTrackReference*)fTPCReferences->At(0)); | |
218 | } | |
219 | if (fMass>0){ | |
220 | Float_t p = TMath::Sqrt(fTrackRef.Px()*fTrackRef.Px()+ | |
221 | fTrackRef.Py()*fTrackRef.Py()+ | |
222 | fTrackRef.Pz()*fTrackRef.Pz()); | |
223 | if (p>0.001){ | |
224 | Float_t betagama = p /fMass; | |
225 | fPrim = TPCBetheBloch(betagama); | |
226 | }else fPrim=0; | |
227 | } | |
228 | }else{ | |
229 | fMass =0; | |
230 | fPrim =0; | |
231 | } | |
232 | } | |
233 | ||
234 | ///////////////////////////////////////////////////////////////////////////////// | |
022044bf | 235 | AliGenV0Info::AliGenV0Info(): |
236 | fMCd(), //info about daughter particle - | |
237 | fMCm(), //info about mother particle - first particle for V0 | |
238 | fMotherP(), //particle info about mother particle | |
239 | fMCDist1(0), //info about closest distance according closest MC - linear DCA | |
240 | fMCDist2(0), //info about closest distance parabolic DCA | |
241 | fMCRr(0), // rec position of the vertex | |
242 | fMCR(0), //exact r position of the vertex | |
243 | fInvMass(0), //reconstructed invariant mass - | |
244 | fPointAngleFi(0), //point angle fi | |
245 | fPointAngleTh(0), //point angle theta | |
246 | fPointAngle(0) //point angle full | |
c92725b7 | 247 | { |
022044bf | 248 | for (Int_t i=0;i<3; i++){ |
249 | fMCPdr[i]=0; | |
250 | fMCX[i]=0; | |
251 | fMCXr[i]=0; | |
252 | fMCPm[i]=0; | |
253 | fMCAngle[i]=0; | |
254 | fMCPd[i]=0; | |
255 | } | |
256 | fMCPd[3]=0; | |
257 | for (Int_t i=0; i<2;i++){ | |
258 | fPdg[i]=0; | |
259 | fLab[i]=0; | |
260 | } | |
c92725b7 | 261 | } |
c92725b7 | 262 | |
263 | void AliGenV0Info::Update(Float_t vertex[3]) | |
264 | { | |
022044bf | 265 | // |
266 | // Update information - calculates derived variables | |
267 | // | |
268 | ||
d92975ba | 269 | fMCPd[0] = fMCd.GetParticle().Px(); |
270 | fMCPd[1] = fMCd.GetParticle().Py(); | |
271 | fMCPd[2] = fMCd.GetParticle().Pz(); | |
272 | fMCPd[3] = fMCd.GetParticle().P(); | |
273 | // | |
274 | fMCX[0] = fMCd.GetParticle().Vx(); | |
275 | fMCX[1] = fMCd.GetParticle().Vy(); | |
276 | fMCX[2] = fMCd.GetParticle().Vz(); | |
c92725b7 | 277 | fMCR = TMath::Sqrt( fMCX[0]*fMCX[0]+fMCX[1]*fMCX[1]); |
278 | // | |
d92975ba | 279 | fPdg[0] = fMCd.GetParticle().GetPdgCode(); |
280 | fPdg[1] = fMCm.GetParticle().GetPdgCode(); | |
c92725b7 | 281 | // |
d92975ba | 282 | fLab[0] = fMCd.GetParticle().GetUniqueID(); |
283 | fLab[1] = fMCm.GetParticle().GetUniqueID(); | |
c92725b7 | 284 | // |
285 | // | |
286 | // | |
287 | Double_t x1[3],p1[3]; | |
d92975ba | 288 | TParticle& pdaughter = fMCd.GetParticle(); |
c92725b7 | 289 | x1[0] = pdaughter.Vx(); |
290 | x1[1] = pdaughter.Vy(); | |
291 | x1[2] = pdaughter.Vz(); | |
292 | p1[0] = pdaughter.Px(); | |
293 | p1[1] = pdaughter.Py(); | |
294 | p1[2] = pdaughter.Pz(); | |
295 | Double_t sign = (pdaughter.GetPDG()->Charge()>0)? -1:1; | |
296 | AliHelix dhelix1(x1,p1,sign); | |
297 | // | |
298 | // | |
299 | Double_t x2[3],p2[3]; | |
300 | // | |
d92975ba | 301 | TParticle& pmother = fMCm.GetParticle(); |
c92725b7 | 302 | x2[0] = pmother.Vx(); |
303 | x2[1] = pmother.Vy(); | |
304 | x2[2] = pmother.Vz(); | |
305 | p2[0] = pmother.Px(); | |
306 | p2[1] = pmother.Py(); | |
307 | p2[2] = pmother.Pz(); | |
308 | // | |
309 | // | |
310 | sign = (pmother.GetPDG()->Charge()>0) ? -1:1; | |
311 | AliHelix mhelix(x2,p2,sign); | |
312 | ||
313 | // | |
314 | // | |
315 | // | |
316 | //find intersection linear | |
317 | // | |
318 | Double_t distance1, distance2; | |
319 | Double_t phase[2][2],radius[2]; | |
320 | Int_t points = dhelix1.GetRPHIintersections(mhelix, phase, radius); | |
321 | Double_t delta1=10000,delta2=10000; | |
322 | if (points>0){ | |
323 | dhelix1.LinearDCA(mhelix,phase[0][0],phase[0][1],radius[0],delta1); | |
324 | dhelix1.LinearDCA(mhelix,phase[0][0],phase[0][1],radius[0],delta1); | |
325 | dhelix1.LinearDCA(mhelix,phase[0][0],phase[0][1],radius[0],delta1); | |
326 | } | |
327 | else{ | |
328 | fInvMass=-1; | |
329 | return; | |
330 | } | |
331 | if (points==2){ | |
332 | dhelix1.LinearDCA(mhelix,phase[1][0],phase[1][1],radius[1],delta2); | |
333 | dhelix1.LinearDCA(mhelix,phase[1][0],phase[1][1],radius[1],delta2); | |
334 | dhelix1.LinearDCA(mhelix,phase[1][0],phase[1][1],radius[1],delta2); | |
335 | } | |
336 | distance1 = TMath::Min(delta1,delta2); | |
337 | // | |
338 | //find intersection parabolic | |
339 | // | |
340 | points = dhelix1.GetRPHIintersections(mhelix, phase, radius); | |
341 | delta1=10000,delta2=10000; | |
342 | ||
343 | if (points>0){ | |
344 | dhelix1.ParabolicDCA(mhelix,phase[0][0],phase[0][1],radius[0],delta1); | |
345 | } | |
346 | if (points==2){ | |
347 | dhelix1.ParabolicDCA(mhelix,phase[1][0],phase[1][1],radius[1],delta2); | |
348 | } | |
349 | ||
350 | distance2 = TMath::Min(delta1,delta2); | |
351 | // | |
352 | if (delta1<delta2){ | |
353 | //get V0 info | |
354 | dhelix1.Evaluate(phase[0][0],fMCXr); | |
355 | dhelix1.GetMomentum(phase[0][0],fMCPdr); | |
356 | mhelix.GetMomentum(phase[0][1],fMCPm); | |
357 | dhelix1.GetAngle(phase[0][0],mhelix,phase[0][1],fMCAngle); | |
358 | fMCRr = TMath::Sqrt(radius[0]); | |
359 | } | |
360 | else{ | |
361 | dhelix1.Evaluate(phase[1][0],fMCXr); | |
362 | dhelix1.GetMomentum(phase[1][0],fMCPdr); | |
363 | mhelix.GetMomentum(phase[1][1],fMCPm); | |
364 | dhelix1.GetAngle(phase[1][0],mhelix,phase[1][1],fMCAngle); | |
365 | fMCRr = TMath::Sqrt(radius[1]); | |
366 | } | |
367 | // | |
368 | // | |
369 | fMCDist1 = TMath::Sqrt(distance1); | |
370 | fMCDist2 = TMath::Sqrt(distance2); | |
371 | // | |
372 | // | |
373 | // | |
374 | Float_t v[3] = {fMCXr[0]-vertex[0],fMCXr[1]-vertex[1],fMCXr[2]-vertex[2]}; | |
375 | //Float_t v[3] = {fMCXr[0],fMCXr[1],fMCXr[2]}; | |
376 | Float_t p[3] = {fMCPdr[0]+fMCPm[0], fMCPdr[1]+fMCPm[1],fMCPdr[2]+fMCPm[2]}; | |
377 | Float_t vnorm2 = v[0]*v[0]+v[1]*v[1]; | |
378 | Float_t vnorm3 = TMath::Sqrt(v[2]*v[2]+vnorm2); | |
379 | vnorm2 = TMath::Sqrt(vnorm2); | |
380 | Float_t pnorm2 = p[0]*p[0]+p[1]*p[1]; | |
381 | Float_t pnorm3 = TMath::Sqrt(p[2]*p[2]+pnorm2); | |
382 | pnorm2 = TMath::Sqrt(pnorm2); | |
383 | // | |
384 | fPointAngleFi = (v[0]*p[0]+v[1]*p[1])/(vnorm2*pnorm2); | |
385 | fPointAngleTh = (v[2]*p[2]+vnorm2*pnorm2)/(vnorm3*pnorm3); | |
386 | fPointAngle = (v[0]*p[0]+v[1]*p[1]+v[2]*p[2])/(vnorm3*pnorm3); | |
d92975ba | 387 | Double_t mass1 = fMCd.GetMass(); |
388 | Double_t mass2 = fMCm.GetMass(); | |
c92725b7 | 389 | |
390 | ||
391 | Double_t e1 = TMath::Sqrt(mass1*mass1+ | |
392 | fMCPd[0]*fMCPd[0]+ | |
393 | fMCPd[1]*fMCPd[1]+ | |
394 | fMCPd[2]*fMCPd[2]); | |
395 | Double_t e2 = TMath::Sqrt(mass2*mass2+ | |
396 | fMCPm[0]*fMCPm[0]+ | |
397 | fMCPm[1]*fMCPm[1]+ | |
398 | fMCPm[2]*fMCPm[2]); | |
399 | ||
400 | fInvMass = | |
401 | (fMCPm[0]+fMCPd[0])*(fMCPm[0]+fMCPd[0])+ | |
402 | (fMCPm[1]+fMCPd[1])*(fMCPm[1]+fMCPd[1])+ | |
403 | (fMCPm[2]+fMCPd[2])*(fMCPm[2]+fMCPd[2]); | |
404 | ||
405 | // fInvMass = TMath::Sqrt((e1+e2)*(e1+e2)-fInvMass); | |
406 | fInvMass = (e1+e2)*(e1+e2)-fInvMass; | |
022044bf | 407 | if (fInvMass>0) fInvMass = TMath::Sqrt(fInvMass); |
c92725b7 | 408 | } |
409 | ||
410 | ||
411 | ||
412 | ///////////////////////////////////////////////////////////////////////////////// | |
022044bf | 413 | |
414 | AliGenKinkInfo::AliGenKinkInfo(): | |
415 | fMCd(), //info about daughter particle - second particle for V0 | |
416 | fMCm(), //info about mother particle - first particle for V0 | |
417 | fMCDist1(0), //info about closest distance according closest MC - linear DCA | |
418 | fMCDist2(0), //info about closest distance parabolic DCA | |
419 | fMCRr(0), // rec position of the vertex | |
420 | fMCR(0) //exact r position of the vertex | |
421 | { | |
422 | // | |
423 | // default constructor | |
424 | // | |
425 | for (Int_t i=0;i<3;i++){ | |
426 | fMCPdr[i]=0; | |
427 | fMCPd[i]=0; | |
428 | fMCX[i]=0; | |
429 | fMCPm[i]=0; | |
430 | fMCAngle[i]=0; | |
431 | } | |
432 | for (Int_t i=0; i<2; i++) { | |
433 | fPdg[i]= 0; | |
434 | fLab[i]=0; | |
435 | } | |
436 | } | |
437 | ||
c92725b7 | 438 | void AliGenKinkInfo::Update() |
439 | { | |
022044bf | 440 | // |
441 | // Update information | |
442 | // some redundancy - faster acces to the values in analysis code | |
443 | // | |
d92975ba | 444 | fMCPd[0] = fMCd.GetParticle().Px(); |
445 | fMCPd[1] = fMCd.GetParticle().Py(); | |
446 | fMCPd[2] = fMCd.GetParticle().Pz(); | |
447 | fMCPd[3] = fMCd.GetParticle().P(); | |
c92725b7 | 448 | // |
d92975ba | 449 | fMCX[0] = fMCd.GetParticle().Vx(); |
450 | fMCX[1] = fMCd.GetParticle().Vy(); | |
451 | fMCX[2] = fMCd.GetParticle().Vz(); | |
c92725b7 | 452 | fMCR = TMath::Sqrt( fMCX[0]*fMCX[0]+fMCX[1]*fMCX[1]); |
453 | // | |
d92975ba | 454 | fPdg[0] = fMCd.GetParticle().GetPdgCode(); |
455 | fPdg[1] = fMCm.GetParticle().GetPdgCode(); | |
c92725b7 | 456 | // |
d92975ba | 457 | fLab[0] = fMCd.GetParticle().GetUniqueID(); |
458 | fLab[1] = fMCm.GetParticle().GetUniqueID(); | |
c92725b7 | 459 | // |
460 | // | |
461 | // | |
462 | Double_t x1[3],p1[3]; | |
d92975ba | 463 | TParticle& pdaughter = fMCd.GetParticle(); |
c92725b7 | 464 | x1[0] = pdaughter.Vx(); |
465 | x1[1] = pdaughter.Vy(); | |
466 | x1[2] = pdaughter.Vz(); | |
467 | p1[0] = pdaughter.Px(); | |
468 | p1[1] = pdaughter.Py(); | |
469 | p1[2] = pdaughter.Pz(); | |
470 | Double_t sign = (pdaughter.GetPDG()->Charge()>0)? -1:1; | |
471 | AliHelix dhelix1(x1,p1,sign); | |
472 | // | |
473 | // | |
474 | Double_t x2[3],p2[3]; | |
475 | // | |
d92975ba | 476 | TParticle& pmother = fMCm.GetParticle(); |
c92725b7 | 477 | x2[0] = pmother.Vx(); |
478 | x2[1] = pmother.Vy(); | |
479 | x2[2] = pmother.Vz(); | |
480 | p2[0] = pmother.Px(); | |
481 | p2[1] = pmother.Py(); | |
482 | p2[2] = pmother.Pz(); | |
483 | // | |
d92975ba | 484 | const AliTrackReference & pdecay = fMCm.GetTRdecay(); |
c92725b7 | 485 | x2[0] = pdecay.X(); |
486 | x2[1] = pdecay.Y(); | |
487 | x2[2] = pdecay.Z(); | |
488 | p2[0] = pdecay.Px(); | |
489 | p2[1] = pdecay.Py(); | |
490 | p2[2] = pdecay.Pz(); | |
491 | // | |
492 | sign = (pmother.GetPDG()->Charge()>0) ? -1:1; | |
493 | AliHelix mhelix(x2,p2,sign); | |
494 | ||
495 | // | |
496 | // | |
497 | // | |
498 | //find intersection linear | |
499 | // | |
500 | Double_t distance1, distance2; | |
501 | Double_t phase[2][2],radius[2]; | |
502 | Int_t points = dhelix1.GetRPHIintersections(mhelix, phase, radius); | |
503 | Double_t delta1=10000,delta2=10000; | |
504 | if (points>0){ | |
505 | dhelix1.LinearDCA(mhelix,phase[0][0],phase[0][1],radius[0],delta1); | |
506 | dhelix1.LinearDCA(mhelix,phase[0][0],phase[0][1],radius[0],delta1); | |
507 | dhelix1.LinearDCA(mhelix,phase[0][0],phase[0][1],radius[0],delta1); | |
508 | } | |
509 | if (points==2){ | |
510 | dhelix1.LinearDCA(mhelix,phase[1][0],phase[1][1],radius[1],delta2); | |
511 | dhelix1.LinearDCA(mhelix,phase[1][0],phase[1][1],radius[1],delta2); | |
512 | dhelix1.LinearDCA(mhelix,phase[1][0],phase[1][1],radius[1],delta2); | |
513 | } | |
514 | distance1 = TMath::Min(delta1,delta2); | |
515 | // | |
516 | //find intersection parabolic | |
517 | // | |
518 | points = dhelix1.GetRPHIintersections(mhelix, phase, radius); | |
519 | delta1=10000,delta2=10000; | |
520 | ||
521 | if (points>0){ | |
522 | dhelix1.ParabolicDCA(mhelix,phase[0][0],phase[0][1],radius[0],delta1); | |
523 | } | |
524 | if (points==2){ | |
525 | dhelix1.ParabolicDCA(mhelix,phase[1][0],phase[1][1],radius[1],delta2); | |
526 | } | |
527 | ||
528 | distance2 = TMath::Min(delta1,delta2); | |
529 | // | |
530 | if (delta1<delta2){ | |
531 | //get V0 info | |
532 | dhelix1.Evaluate(phase[0][0],fMCXr); | |
533 | dhelix1.GetMomentum(phase[0][0],fMCPdr); | |
534 | mhelix.GetMomentum(phase[0][1],fMCPm); | |
535 | dhelix1.GetAngle(phase[0][0],mhelix,phase[0][1],fMCAngle); | |
536 | fMCRr = TMath::Sqrt(radius[0]); | |
537 | } | |
538 | else{ | |
539 | dhelix1.Evaluate(phase[1][0],fMCXr); | |
540 | dhelix1.GetMomentum(phase[1][0],fMCPdr); | |
541 | mhelix.GetMomentum(phase[1][1],fMCPm); | |
542 | dhelix1.GetAngle(phase[1][0],mhelix,phase[1][1],fMCAngle); | |
543 | fMCRr = TMath::Sqrt(radius[1]); | |
544 | } | |
545 | // | |
546 | // | |
547 | fMCDist1 = TMath::Sqrt(distance1); | |
548 | fMCDist2 = TMath::Sqrt(distance2); | |
549 | ||
550 | } | |
551 | ||
552 | ||
553 | Float_t AliGenKinkInfo::GetQt(){ | |
554 | // | |
555 | // | |
556 | Float_t momentumd = TMath::Sqrt(fMCPd[0]*fMCPd[0]+fMCPd[1]*fMCPd[1]+fMCPd[2]*fMCPd[2]); | |
557 | return TMath::Sin(fMCAngle[2])*momentumd; | |
558 | } | |
559 | ||
560 | ||
561 | ||
c92725b7 | 562 | |
563 | //////////////////////////////////////////////////////////////////////// | |
564 | AliTPCdigitRow::AliTPCdigitRow() | |
565 | { | |
566 | Reset(); | |
567 | } | |
568 | //////////////////////////////////////////////////////////////////////// | |
569 | AliTPCdigitRow & AliTPCdigitRow::operator=(const AliTPCdigitRow &digOld) | |
570 | { | |
571 | for (Int_t i = 0; i<kgRowBytes; i++) fDig[i] = digOld.fDig[i]; | |
572 | return (*this); | |
573 | } | |
574 | //////////////////////////////////////////////////////////////////////// | |
575 | void AliTPCdigitRow::SetRow(Int_t row) | |
576 | { | |
577 | if (row >= 8*kgRowBytes) { | |
578 | cerr<<"AliTPCdigitRow::SetRow: index "<<row<<" out of bounds."<<endl; | |
579 | return; | |
580 | } | |
581 | Int_t iC = row/8; | |
582 | Int_t iB = row%8; | |
583 | SETBIT(fDig[iC],iB); | |
584 | } | |
585 | ||
586 | //////////////////////////////////////////////////////////////////////// | |
022044bf | 587 | Bool_t AliTPCdigitRow::TestRow(Int_t row) const |
c92725b7 | 588 | { |
589 | // | |
590 | // return kTRUE if row is on | |
591 | // | |
592 | Int_t iC = row/8; | |
593 | Int_t iB = row%8; | |
594 | return TESTBIT(fDig[iC],iB); | |
595 | } | |
596 | //////////////////////////////////////////////////////////////////////// | |
022044bf | 597 | Int_t AliTPCdigitRow::RowsOn(Int_t upto) const |
c92725b7 | 598 | { |
599 | // | |
600 | // returns number of rows with a digit | |
601 | // count only rows less equal row number upto | |
602 | // | |
603 | Int_t total = 0; | |
604 | for (Int_t i = 0; i<kgRowBytes; i++) { | |
605 | for (Int_t j = 0; j < 8; j++) { | |
606 | if (i*8+j > upto) return total; | |
607 | if (TESTBIT(fDig[i],j)) total++; | |
608 | } | |
609 | } | |
610 | return total; | |
611 | } | |
612 | //////////////////////////////////////////////////////////////////////// | |
613 | void AliTPCdigitRow::Reset() | |
614 | { | |
615 | // | |
616 | // resets all rows to zero | |
617 | // | |
618 | for (Int_t i = 0; i<kgRowBytes; i++) { | |
619 | fDig[i] <<= 8; | |
620 | } | |
621 | } | |
622 | //////////////////////////////////////////////////////////////////////// | |
022044bf | 623 | Int_t AliTPCdigitRow::Last() const |
c92725b7 | 624 | { |
625 | // | |
626 | // returns the last row number with a digit | |
627 | // returns -1 if now digits | |
628 | // | |
629 | for (Int_t i = kgRowBytes-1; i>=0; i--) { | |
630 | for (Int_t j = 7; j >= 0; j--) { | |
631 | if TESTBIT(fDig[i],j) return i*8+j; | |
632 | } | |
633 | } | |
634 | return -1; | |
635 | } | |
636 | //////////////////////////////////////////////////////////////////////// | |
022044bf | 637 | Int_t AliTPCdigitRow::First() const |
c92725b7 | 638 | { |
639 | // | |
640 | // returns the first row number with a digit | |
641 | // returns -1 if now digits | |
642 | // | |
643 | for (Int_t i = 0; i<kgRowBytes; i++) { | |
644 | for (Int_t j = 0; j < 8; j++) { | |
645 | if (TESTBIT(fDig[i],j)) return i*8+j; | |
646 | } | |
647 | } | |
648 | return -1; | |
649 | } | |
650 | ||
022044bf | 651 | |
d92975ba | 652 | //_____________________________________________________________________________ |
653 | Float_t AliMCInfo::TPCBetheBloch(Float_t bg) | |
c92725b7 | 654 | { |
655 | // | |
d92975ba | 656 | // Bethe-Bloch energy loss formula |
c92725b7 | 657 | // |
d92975ba | 658 | const Double_t kp1=0.76176e-1; |
659 | const Double_t kp2=10.632; | |
660 | const Double_t kp3=0.13279e-4; | |
661 | const Double_t kp4=1.8631; | |
662 | const Double_t kp5=1.9479; | |
c92725b7 | 663 | |
d92975ba | 664 | Double_t dbg = (Double_t) bg; |
c92725b7 | 665 | |
d92975ba | 666 | Double_t beta = dbg/TMath::Sqrt(1.+dbg*dbg); |
c92725b7 | 667 | |
d92975ba | 668 | Double_t aa = TMath::Power(beta,kp4); |
669 | Double_t bb = TMath::Power(1./dbg,kp5); | |
c92725b7 | 670 | |
d92975ba | 671 | bb=TMath::Log(kp3+bb); |
c92725b7 | 672 | |
d92975ba | 673 | return ((Float_t)((kp2-aa-bb)*kp1/aa)); |
c92725b7 | 674 | } |
c92725b7 | 675 |