]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGCF/FEMTOSCOPY/AliFemto/AliFemtoParticle.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Merge branch 'feature-movesplit'
[u/mrichter/AliRoot.git] / PWGCF / FEMTOSCOPY / AliFemto / AliFemtoParticle.cxx
CommitLineData
76ce4b5b 1///////////////////////////////////////////////////////////////////////////
2// //
3// AliFemtoParticle: main class halding all the necessary information //
4// about particle that is required during femtoscopic analysis //
5// This includes all the information about the quality of the track, //
6// its identification as well as track chracteristics with connection //
7// to the detector parts, e.g. entrance and exit points. //
8// //
9///////////////////////////////////////////////////////////////////////////
10#include "AliFemtoKink.h"
11#include "AliFemtoParticle.h"
12#include "AliFemtoXi.h"
13
14double AliFemtoParticle::fgPrimPimPar0= 9.05632e-01;
15double AliFemtoParticle::fgPrimPimPar1= -2.26737e-01;
16double AliFemtoParticle::fgPrimPimPar2= -1.03922e-01;
17double AliFemtoParticle::fgPrimPipPar0= 9.09616e-01;
18double AliFemtoParticle::fgPrimPipPar1= -9.00511e-02;
19double AliFemtoParticle::fgPrimPipPar2= -6.02940e-02;
20double AliFemtoParticle::fgPrimPmPar0= 0.;
21double AliFemtoParticle::fgPrimPmPar1= 0.;
22double AliFemtoParticle::fgPrimPmPar2= 0.;
23double AliFemtoParticle::fgPrimPpPar0= 0.;
24double AliFemtoParticle::fgPrimPpPar1= 0.;
25double AliFemtoParticle::fgPrimPpPar2= 0.;
26
27int TpcLocalTransform(AliFmThreeVectorD& xgl,
28 int& iSector,
29 int& iPadrow,
30 float& xlocal,
31 double& ttPhi);
32
33
34//_____________________
35AliFemtoParticle::AliFemtoParticle() :
36// fTpcV0NegPosSample(0),
37// fV0NegZ(0),
38// fV0NegU(0),
39// fV0NegSect(0),
40 fTrack(0), fV0(0), fKink(0), fXi(0),
41 fFourMomentum(0),
42 fHelix(),
43// fNominalTpcExitPoint(0),
44// fNominalTpcEntrancePoint(0),
45 fHiddenInfo(0),
46 fPrimaryVertex(0),
47 fSecondaryVertex(0),
48 fHelixV0Pos(),
49 fTpcV0PosEntrancePoint(0),
50 fTpcV0PosExitPoint(0),
51 fHelixV0Neg(),
52 fTpcV0NegEntrancePoint(0),
53 fTpcV0NegExitPoint(0)
54{
55 // Default constructor
56 /* no-op for default */
57 for (int ip=0; ip<6; ip++) fPurity[ip] = 0.0;
58 // cout << "Created particle " << this << endl;
59}
60//_____________________
61AliFemtoParticle::AliFemtoParticle(const AliFemtoParticle& aParticle):
62// fTpcV0NegPosSample(0),
63// fV0NegZ(0),
64// fV0NegU(0),
65// fV0NegSect(0),
66 fTrack(0), fV0(0), fKink(0), fXi(0),
67 fFourMomentum(0),
68 fHelix(),
69// fNominalTpcExitPoint(0),
70// fNominalTpcEntrancePoint(0),
71 fHiddenInfo(0),
72 fPrimaryVertex(0),
73 fSecondaryVertex(0),
74 fHelixV0Pos(),
75 fTpcV0PosEntrancePoint(0),
76 fTpcV0PosExitPoint(0),
77 fHelixV0Neg(),
78 fTpcV0NegEntrancePoint(0),
79 fTpcV0NegExitPoint(0)
80{
81 // Copy constructor
82 if (aParticle.fTrack)
83 fTrack = new AliFemtoTrack(*aParticle.fTrack);
84 if (aParticle.fV0)
85 fV0 = new AliFemtoV0(*aParticle.fV0);
86 if (aParticle.fKink)
87 fKink = new AliFemtoKink(*aParticle.fKink);
88 if (aParticle.fXi)
89 fXi = new AliFemtoXi(*aParticle.fXi);
90 fFourMomentum = aParticle.fFourMomentum;
91 fHelix = aParticle.fHelix;
92 for (int ip=0; ip<6; ip++) fPurity[ip] = 0.0;
93
94// for (int iter=0; iter<11; iter++)
95// fNominalPosSample[iter] = aParticle.fNominalPosSample[iter];
96
97// if (aParticle.fTpcV0NegPosSample) {
98// fTpcV0NegPosSample = (AliFemtoThreeVector *) malloc(sizeof(AliFemtoThreeVector) * 11);
99// for (int iter=0; iter<11; iter++)
100// fTpcV0NegPosSample[iter] = aParticle.fTpcV0NegPosSample[iter];
101// }
102
103// if (aParticle.fV0NegZ) {
104// fV0NegZ = (float *) malloc(sizeof(float) * 45);
105// for (int iter=0; iter<11; iter++)
106// fV0NegZ[iter] = aParticle.fV0NegZ[iter];
107// }
108// if (aParticle.fV0NegU) {
109// fV0NegU = (float *) malloc(sizeof(float) * 45);
110// for (int iter=0; iter<11; iter++)
111// fV0NegU[iter] = aParticle.fV0NegU[iter];
112// }
113// if (aParticle.fV0NegSect) {
114// fV0NegSect = (int *) malloc(sizeof(int) * 45);
115// for (int iter=0; iter<11; iter++)
116// fV0NegSect[iter] = aParticle.fV0NegSect[iter];
117// }
118
119 fPrimaryVertex = aParticle.fPrimaryVertex;
120 fSecondaryVertex = aParticle.fSecondaryVertex;
121 CalculatePurity();
122 if(aParticle.fHiddenInfo){
123 fHiddenInfo= aParticle.HiddenInfo()->Clone();
124 }
125
126// fNominalTpcEntrancePoint = aParticle.fNominalTpcEntrancePoint;
127// fNominalTpcExitPoint = aParticle.fNominalTpcExitPoint;
128
129 for (int iter=0; iter<6; iter++)
130 fPurity[iter] = aParticle.fPurity[iter];
131
132 fHelixV0Pos = aParticle.fHelixV0Pos;
133 fTpcV0PosEntrancePoint = aParticle.fTpcV0PosEntrancePoint;
134 fTpcV0PosExitPoint = aParticle.fTpcV0PosExitPoint;
135 fHelixV0Neg = aParticle.fHelixV0Neg;
136 fTpcV0NegEntrancePoint = aParticle.fTpcV0NegEntrancePoint;
137 fTpcV0NegExitPoint = aParticle.fTpcV0NegExitPoint;
138}
139//_____________________
140AliFemtoParticle::~AliFemtoParticle(){
141 // cout << "Issuing delete for AliFemtoParticle." << endl;
142
143 if (fTrack) delete fTrack;
144 if (fV0) {
145// delete[] fTpcV0NegPosSample;
146// delete[] fV0NegZ;
147// delete[] fV0NegU;
148// delete[] fV0NegSect;
149 delete fV0;
150 }
151 if (fKink) delete fKink;
152 if (fXi) delete fXi;
153 // cout << "Trying to delete HiddenInfo: " << fHiddenInfo << endl;
154 if (fHiddenInfo)
155 {
156 // cout << "Deleting HiddenInfo." << endl;
157 delete fHiddenInfo;
158 }
159 // cout << "Deleted particle " << this << endl;
160}
161//_____________________
162AliFemtoParticle::AliFemtoParticle(const AliFemtoTrack* const hbtTrack,const double& mass) :
163// fTpcV0NegPosSample(0),
164// fV0NegZ(0),
165// fV0NegU(0),
166// fV0NegSect(0),
167 fTrack(0), fV0(0), fKink(0), fXi(0),
168 fFourMomentum(0),
169 fHelix(),
170// fNominalTpcExitPoint(0),
171// fNominalTpcEntrancePoint(0),
172 fHiddenInfo(0),
173 fPrimaryVertex(0),
174 fSecondaryVertex(0),
175 fHelixV0Pos(),
176 fTpcV0PosEntrancePoint(0),
177 fTpcV0PosExitPoint(0),
178 fHelixV0Neg(),
179 fTpcV0NegEntrancePoint(0),
180 fTpcV0NegExitPoint(0)
181{
182 // Constructor from normal track
183
184 // I know there is a better way to do this...
185 fTrack = new AliFemtoTrack(*hbtTrack);
186 AliFemtoThreeVector temp = hbtTrack->P();
187 fFourMomentum.SetVect(temp);
188 double ener = ::sqrt(temp.Mag2()+mass*mass);
189 fFourMomentum.SetE(ener);
190// fMap[0] = hbtTrack->TopologyMap(0);
191 // fMap[1] = hbtTrack->TopologyMap(1);
192 // fNhits = hbtTrack->NHits();
193 fHelix = hbtTrack->Helix();
194 //CalculateNominalTpcExitAndEntrancePoints();
195
196
197 fPrimaryVertex.SetX(0.);
198 fPrimaryVertex.SetY(0.);
199 fPrimaryVertex.SetZ(0.);
200 fSecondaryVertex.SetX(0.);
201 fSecondaryVertex.SetY(0.);
202 fSecondaryVertex.SetZ(0.);
203 /* TO JA ODZNACZYLEM NIE WIEM DLACZEGO
204 CalculateTpcExitAndEntrancePoints(&fHelix,&fPrimaryVertex,
205 &fSecondaryVertex,
206 &fNominalTpcEntrancePoint,
207 &fNominalTpcExitPoint,
208 &mNominalPosSample[0],
209 &fZ[0],
210 &fU[0],
211 &fSect[0]);
212 */
213 CalculatePurity();
214 // ***
215 fHiddenInfo= 0;
216 if(hbtTrack->ValidHiddenInfo()){
217 fHiddenInfo= hbtTrack->GetHiddenInfo()->Clone();
218 }
219 // ***
220 // cout << "Created particle " << this << endl;
221
222}
223//_____________________
224AliFemtoParticle::AliFemtoParticle(const AliFemtoV0* const hbtV0,const double& mass) :
225// fTpcV0NegPosSample(0),
226// fV0NegZ(0),
227// fV0NegU(0),
228// fV0NegSect(0),
229 fTrack(0), fV0(0), fKink(0), fXi(0),
230 fFourMomentum(0),
231 fHelix(),
232// fNominalTpcExitPoint(0),
233// fNominalTpcEntrancePoint(0),
234 fHiddenInfo(0),
235 fPrimaryVertex(0),
236 fSecondaryVertex(0),
237 fHelixV0Pos(),
238 fTpcV0PosEntrancePoint(0),
239 fTpcV0PosExitPoint(0),
240 fHelixV0Neg(),
241 fTpcV0NegEntrancePoint(0),
242 fTpcV0NegExitPoint(0)
243{
244 // Constructor from V0
245 fV0 = new AliFemtoV0(*hbtV0);
246 //fMap[0]= 0;
247 //fMap[1]= 0;
248 // I know there is a better way to do this...
249 AliFemtoThreeVector temp = hbtV0->MomV0();
250 fFourMomentum.SetVect(temp);
251 double ener = ::sqrt(temp.Mag2()+mass*mass);
252 fFourMomentum.SetE(ener);
253 // Calculating TpcEntrancePoint for Positive V0 daugther
254 fPrimaryVertex = hbtV0->PrimaryVertex();
255 fSecondaryVertex = hbtV0->DecayVertexV0();
256 fHelixV0Pos = hbtV0->HelixPos();
257
258// fTpcV0NegPosSample = new AliFemtoThreeVector[45];//for V0Neg
259// fV0NegZ = new float[45];//for V0Neg
260// fV0NegU = new float[45];//for V0Neg
261// fV0NegSect = new int[45];//for V0Neg
262// CalculateTpcExitAndEntrancePoints(&fHelixV0Pos,&fPrimaryVertex,
263// &fSecondaryVertex,
264// &fTpcV0PosEntrancePoint,
265// &fTpcV0PosExitPoint,
266// &fNominalPosSample[0],
267// &fZ[0],
268// &fU[0],&fSect[0]);
269 fHelixV0Neg = hbtV0->HelixNeg();
270
271// CalculateTpcExitAndEntrancePoints(&fHelixV0Neg,
272// &fPrimaryVertex,
273// &fSecondaryVertex,
274// &fTpcV0NegEntrancePoint,
275// &fTpcV0NegExitPoint,
276// &fTpcV0NegPosSample[0],
277// &fV0NegZ[0],
278// &fV0NegU[0],&fV0NegSect[0]);
279
280 // ***
281 fHiddenInfo= 0;
282 if(hbtV0->ValidHiddenInfo()){
283 fHiddenInfo= hbtV0->GetHiddenInfo()->Clone();
284 }
285 for (int ip=0; ip<6; ip++) fPurity[ip] = 0.0;
286 // ***
287}
288//_____________________
289AliFemtoParticle::AliFemtoParticle(const AliFemtoKink* const hbtKink,const double& mass) :
290// fTpcV0NegPosSample(0),
291// fV0NegZ(0),
292// fV0NegU(0),
293// fV0NegSect(0),
294 fTrack(0), fV0(0), fKink(0), fXi(0),
295 fFourMomentum(0),
296 fHelix(),
297// fNominalTpcExitPoint(0),
298// fNominalTpcEntrancePoint(0),
299 fHiddenInfo(0),
300 fPrimaryVertex(0),
301 fSecondaryVertex(0),
302 fHelixV0Pos(),
303 fTpcV0PosEntrancePoint(0),
304 fTpcV0PosExitPoint(0),
305 fHelixV0Neg(),
306 fTpcV0NegEntrancePoint(0),
307 fTpcV0NegExitPoint(0)
308{
309 // Constructor from Kink
310 fKink = new AliFemtoKink(*hbtKink);
311 // fMap[0]= 0;
312 //fMap[1]= 0;
313 // I know there is a better way to do this...
314 AliFemtoThreeVector temp = hbtKink->Parent().P();
315 fFourMomentum.SetVect(temp);
316 double ener = ::sqrt(temp.Mag2()+mass*mass);
317 fFourMomentum.SetE(ener);
318 for (int ip=0; ip<6; ip++) fPurity[ip] = 0.0;
319}
320
321//_____________________
322AliFemtoParticle::AliFemtoParticle(const AliFemtoXi* const hbtXi, const double& mass) :
323// fTpcV0NegPosSample(0),
324// fV0NegZ(0),
325// fV0NegU(0),
326// fV0NegSect(0),
327 fTrack(0), fV0(0), fKink(0), fXi(0),
328 fFourMomentum(0),
329 fHelix(),
330// fNominalTpcExitPoint(0),
331// fNominalTpcEntrancePoint(0),
332 fHiddenInfo(0),
333 fPrimaryVertex(0),
334 fSecondaryVertex(0),
335 fHelixV0Pos(),
336 fTpcV0PosEntrancePoint(0),
337 fTpcV0PosExitPoint(0),
338 fHelixV0Neg(),
339 fTpcV0NegEntrancePoint(0),
340 fTpcV0NegExitPoint(0)
341{
342 // Constructor from Xi
343 fXi = new AliFemtoXi(*hbtXi);
344 // fMap[0]= 0;
345 //fMap[1]= 0;
346 AliFemtoThreeVector temp;// = hbtXi->mMofXi;
347 fFourMomentum.SetVect(temp);
348 double ener = ::sqrt(temp.Mag2()+mass*mass);
349 fFourMomentum.SetE(ener);
350 fHiddenInfo = 0;
351 for (int ip=0; ip<6; ip++) fPurity[ip] = 0.0;
352
353}
354//_____________________
355AliFemtoParticle& AliFemtoParticle::operator=(const AliFemtoParticle& aParticle)
356{
357 // assignment operator
358 if (this != &aParticle) {
359
360 if (fTrack) delete fTrack;
361 if (aParticle.fTrack) {
362 fTrack = new AliFemtoTrack(*aParticle.fTrack);
363 CalculatePurity();
364 }
365 if (fV0) delete fV0;
366 if (aParticle.fV0)
367 fV0 = new AliFemtoV0(*aParticle.fV0);
368 if (fKink) delete fKink;
369 if (aParticle.fKink)
370 fKink = new AliFemtoKink(*aParticle.fKink);
371 if (fXi) delete fXi;
372 if (aParticle.fXi)
373 fXi = new AliFemtoXi(*aParticle.fXi);
374 fFourMomentum = aParticle.fFourMomentum;
375 fHelix = aParticle.fHelix;
376
377 // for (int iter=0; iter<11; iter++)
378 // fNominalPosSample[iter] = aParticle.fNominalPosSample[iter];
379
380 // if (fTpcV0NegPosSample) delete fTpcV0NegPosSample;
381 // if (aParticle.fTpcV0NegPosSample) {
382 // fTpcV0NegPosSample = (AliFemtoThreeVector *) malloc(sizeof(AliFemtoThreeVector) * 11);
383 // for (int iter=0; iter<11; iter++)
384 // fTpcV0NegPosSample[iter] = aParticle.fTpcV0NegPosSample[iter];
385 // }
386
387 // if (fV0NegZ) delete fV0NegZ;
388 // if (aParticle.fV0NegZ) {
389 // fV0NegZ = (float *) malloc(sizeof(float) * 45);
390 // for (int iter=0; iter<11; iter++)
391 // fV0NegZ[iter] = aParticle.fV0NegZ[iter];
392 // }
393 // if (fV0NegU) delete fV0NegU;
394 // if (aParticle.fV0NegU) {
395 // fV0NegU = (float *) malloc(sizeof(float) * 45);
396 // for (int iter=0; iter<11; iter++)
397 // fV0NegU[iter] = aParticle.fV0NegU[iter];
398 // }
399 // if (fV0NegSect) delete fV0NegSect;
400 // if (aParticle.fV0NegSect) {
401 // fV0NegSect = (int *) malloc(sizeof(int) * 45);
402 // for (int iter=0; iter<11; iter++)
403 // fV0NegSect[iter] = aParticle.fV0NegSect[iter];
404 // }
405
406 fPrimaryVertex = aParticle.fPrimaryVertex;
407 fSecondaryVertex = aParticle.fSecondaryVertex;
408// if (fHiddenInfo) delete fHiddenInfo;
409// if(aParticle.fHiddenInfo){
410// fHiddenInfo= aParticle.fHiddenInfo->Clone();
411// }
412
413 // fNominalTpcEntrancePoint = aParticle.fNominalTpcEntrancePoint;
414 // fNominalTpcExitPoint = aParticle.fNominalTpcExitPoint;
415
416 if (fHiddenInfo) delete fHiddenInfo;
417 if (aParticle.fHiddenInfo)
418 fHiddenInfo = aParticle.HiddenInfo()->Clone();
419
420 for (int iter=0; iter<6; iter++)
421 fPurity[iter] = aParticle.fPurity[iter];
422
423 fHelixV0Pos = aParticle.fHelixV0Pos;
424 fTpcV0PosEntrancePoint = aParticle.fTpcV0PosEntrancePoint;
425 fTpcV0PosExitPoint = aParticle.fTpcV0PosExitPoint;
426 fHelixV0Neg = aParticle.fHelixV0Neg;
427 fTpcV0NegEntrancePoint = aParticle.fTpcV0NegEntrancePoint;
428 fTpcV0NegExitPoint = aParticle.fTpcV0NegExitPoint;
429
430 }
431
432 return *this;
433}
434// //_____________________
435// const AliFemtoThreeVector& AliFemtoParticle::NominalTpcExitPoint() const{
436// // in future, may want to calculate this "on demand" only, sot this routine may get more sophisticated
437// // for now, we calculate Exit and Entrance points upon instantiation
438// return fNominalTpcExitPoint;
439// }
440// //_____________________
441// const AliFemtoThreeVector& AliFemtoParticle::NominalTpcEntrancePoint() const{
442// // in future, may want to calculate this "on demand" only, sot this routine may get more sophisticated
443// // for now, we calculate Exit and Entrance points upon instantiation
444// return fNominalTpcEntrancePoint;
445// }
446//_____________________
447void AliFemtoParticle::CalculatePurity(){
448 // Calculate additional parameterized purity
449
450 double tPt = fFourMomentum.Perp();
451 // pi -
452 fPurity[0] = fgPrimPimPar0*(1.-exp((tPt-fgPrimPimPar1)/fgPrimPimPar2));
453 fPurity[0] *= fTrack->PidProbPion();
454 // pi+
455 fPurity[1] = fgPrimPipPar0*(1.-exp((tPt-fgPrimPipPar1)/fgPrimPipPar2));
456 fPurity[1] *= fTrack->PidProbPion();
457 // K-
458 fPurity[2] = fTrack->PidProbKaon();
459 // K+
460 fPurity[3] = fTrack->PidProbKaon();
461 // pbar
462 fPurity[4] = fTrack->PidProbProton();
463 // p
464 fPurity[5] = fTrack->PidProbProton();
465}
466
467double AliFemtoParticle::GetPionPurity()
468{
469 // Get full pion purity
470 if (fTrack->Charge()>0)
471 return fPurity[1];
472 else
473 return fPurity[0];
474}
475double AliFemtoParticle::GetKaonPurity()
476{
477 // Get full kaon purity
478 if (fTrack->Charge()>0)
479 return fPurity[3];
480 else
481 return fPurity[2];
482}
483double AliFemtoParticle::GetProtonPurity()
484{
485 // Get full proton purity
486 if (fTrack->Charge()>0)
487 return fPurity[5];
488 else
489 return fPurity[4];
490}
491
492// void AliFemtoParticle::CalculateTpcExitAndEntrancePoints(AliFmPhysicalHelixD* tHelix,
493// AliFemtoThreeVector* PrimVert,
494// AliFemtoThreeVector* SecVert,
495// AliFemtoThreeVector* tmpTpcEntrancePoint,
496// AliFemtoThreeVector* tmpTpcExitPoint,
497// AliFemtoThreeVector* tmpPosSample,
498// float* tmpZ,
499// float* tmpU,
500// int* tmpSect){
501// // this calculates the exit point of a secondary track,
502// // either through the endcap or through the Outer Field Cage
503// // We assume the track to start at tHelix.origin-PrimaryVertex
504// // it also calculates the entrance point of the secondary track,
505// // which is the point at which it crosses the
506// // inner field cage
507// // static AliFemtoThreeVector ZeroVec(0.,0.,0.);
508// AliFemtoThreeVector tZeroVec(0.,0.,0.);
509// // tZeroVec.SetX(tHelix->origin().x()-PrimVert->x());
510// // tZeroVec.SetY(tHelix->origin().y()-PrimVert->y());
511// // tZeroVec.SetZ(tHelix->origin().z()-PrimVert->z());
512// tZeroVec.SetX(SecVert->x()-PrimVert->x());
513// tZeroVec.SetY(SecVert->y()-PrimVert->y());
514// tZeroVec.SetZ(SecVert->z()-PrimVert->z());
515// double dip, curv, phase;
516// int h;
517// curv = tHelix->Curvature();
518// dip = tHelix->DipAngle();
519// phase= tHelix->Phase();
520// h = tHelix->H();
521
522// AliFmHelixD hel(curv,dip,phase,tZeroVec,h);
523
524// pairD candidates;
525// double sideLength; // this is how much length to go to leave through sides of TPC
526// double endLength; // this is how much length to go to leave through endcap of TPC
527// // figure out how far to go to leave through side...
528// candidates = hel.PathLength(200.0); // bugfix MAL jul00 - 200cm NOT 2cm
529// sideLength = (candidates.first > 0) ? candidates.first : candidates.second;
530
531// static AliFemtoThreeVector tWestEnd(0.,0.,200.); // bugfix MAL jul00 - 200cm NOT 2cm
532// static AliFemtoThreeVector tEastEnd(0.,0.,-200.); // bugfix MAL jul00 - 200cm NOT 2cm
533// static AliFemtoThreeVector tEndCapNormal(0.,0.,1.0);
534
535// endLength = hel.PathLength(tWestEnd,tEndCapNormal);
536// if (endLength < 0.0) endLength = hel.PathLength(tEastEnd,tEndCapNormal);
537
538// if (endLength < 0.0) cout <<
539// "AliFemtoParticle::CalculateTpcExitAndEntrancePoints(): "
540// << "Hey -- I cannot find an exit point out endcaps" << endl;
541// // OK, firstExitLength will be the shortest way out of the detector...
542// double firstExitLength = (endLength < sideLength) ? endLength : sideLength;
543// // now then, let's return the POSITION at which particle leaves TPC...
544// *tmpTpcExitPoint = hel.At(firstExitLength);
545// // Finally, calculate the position at which the track crosses the inner field cage
546// candidates = hel.PathLength(50.0); // bugfix MAL jul00 - 200cm NOT 2cm
547
548// sideLength = (candidates.first > 0) ? candidates.first : candidates.second;
549// // cout << "sideLength 2 ="<<sideLength << endl;
550// *tmpTpcEntrancePoint = hel.At(sideLength);
551// // This is the secure way !
552// if (IsNaN(tmpTpcEntrancePoint->x()) ||
553// IsNaN(tmpTpcEntrancePoint->y()) ||
554// IsNaN(tmpTpcEntrancePoint->z()) ){
555// cout << "tmpTpcEntrancePoint NAN"<< endl;
556// cout << "tmpNominalTpcEntrancePoint = " <<tmpTpcEntrancePoint<< endl;
557// tmpTpcEntrancePoint->SetX(-9999.);
558// tmpTpcEntrancePoint->SetY(-9999.);
559// tmpTpcEntrancePoint->SetZ(-9999.);
560// }
561
562// if (IsNaN(tmpTpcExitPoint->x()) ||
563// IsNaN(tmpTpcExitPoint->y()) ||
564// IsNaN(tmpTpcExitPoint->z()) ) {
565// // cout << "tmpTpcExitPoint NAN Set at (-9999,-9999,-9999)"<< endl;
566// // cout << "tmpTpcExitPoint X= " <<tmpTpcExitPoint->x()<< endl;
567// // cout << "tmpTpcExitPoint Y= " <<tmpTpcExitPoint->y()<< endl;
568// // cout << "tmpTpcExitPoint Z= " <<tmpTpcExitPoint->z()<< endl;
569// tmpTpcExitPoint->SetX(-9999.);
570// tmpTpcExitPoint->SetY(-9999.);
571// tmpTpcExitPoint->SetZ(-9999.);
572// }
573
574
575// // if (IsNaN(tmpTpcExitPoint->x())) *tmpTpcExitPoint = AliFemtoThreeVector(-9999.,-9999.,-9999);
576// // if (IsNaN(tmpTpcEntrancetPoint->x())) *tmpTpcEntrancePoint = AliFemtoThreeVector(-9999.,-9999.,-9999);
577// // cout << "tmpTpcEntrancePoint"<<*tmpTpcEntrancePoint << endl;
578
579// // 03Oct00 - mal. OK, let's try something a little more
580// // along the lines of NA49 and E895 strategy.
581// // calculate the "nominal" position at N radii (say N=11)
582// // within the TPC, and for a pair cut
583// // use the average separation of these N
584// int irad = 0;
585// candidates = hel.PathLength(50.0);
586// sideLength = (candidates.first > 0) ? candidates.first : candidates.second;
587// while (irad<11 && !IsNaN(sideLength)){
588// float radius = 50.0 + irad*15.0;
589// candidates = hel.PathLength(radius);
590// sideLength = (candidates.first > 0) ? candidates.first : candidates.second;
591// tmpPosSample[irad] = hel.At(sideLength);
592// if(IsNaN(tmpPosSample[irad].x()) ||
593// IsNaN(tmpPosSample[irad].y()) ||
594// IsNaN(tmpPosSample[irad].z())
595// ){
596// cout << "tmpPosSample for radius=" << radius << " NAN"<< endl;
597// cout << "tmpPosSample=(" <<tmpPosSample[irad]<<")"<< endl;
598// tmpPosSample[irad] = AliFemtoThreeVector(-9999.,-9999.,-9999);
599// }
600// irad++;
601// if (irad<11){
602// float radius = 50.0 + irad*15.0;
603// candidates = hel.PathLength(radius);
604// sideLength = (candidates.first > 0) ? candidates.first : candidates.second;
605// }
606// }
607// for (int i = irad; i<11; i++)
608// {
609// tmpPosSample[i] = AliFemtoThreeVector(-9999.,-9999.,-9999);
610// }
611
612// static float tRowRadius[45] = {60,64.8,69.6,74.4,79.2,84,88.8,93.6,98.8,
613// 104,109.2,114.4,119.6,127.195,129.195,131.195,
614// 133.195,135.195,137.195,139.195,141.195,
615// 143.195,145.195,147.195,149.195,151.195,
616// 153.195,155.195,157.195,159.195,161.195,
617// 163.195,165.195,167.195,169.195,171.195,
618// 173.195,175.195,177.195,179.195,181.195,
619// 183.195,185.195,187.195,189.195};
620// int tRow,tSect,tOutOfBound;
621// double tLength,tPhi;
622// float tU;
623// AliFemtoThreeVector tPoint;
624// AliFmThreeVectorD tn(0,0,0);
625// AliFmThreeVectorD tr(0,0,0);
626// int ti =0;
627// // test to enter the loop
628// candidates = hel.PathLength(tRowRadius[ti]);
629// tLength = (candidates.first > 0) ? candidates.first : candidates.second;
630// if (IsNaN(tLength)){
631// cout <<"tLength Init tmp NAN" << endl;
632// cout <<"padrow number= "<<ti << "not reached" << endl;
633// cout << "*** DO NOT ENTER THE LOOP***" << endl;
634// tmpSect[ti]=-1;//sector
635// }
636// // end test
637// while(ti<45 && !IsNaN(tLength)){
638// candidates = hel.PathLength(tRowRadius[ti]);
639// tLength = (candidates.first > 0) ? candidates.first : candidates.second;
640// if (IsNaN(tLength)){
641// cout <<"tLength loop 1st NAN" << endl;
642// cout <<"padrow number= " << ti << " not reached" << endl;
643// cout << "*** THIS IS AN ERROR SHOULDN'T LOOP ***" << endl;
644// tmpSect[ti]=-1;//sector
645// }
646// tPoint = hel.At(tLength);
647// // Find which sector it is on
648// TpcLocalTransform(tPoint,tmpSect[ti],tRow,tU,tPhi);
649// if (IsNaN(tmpSect[ti])){
650// cout <<"***ERROR tmpSect"<< endl;
651// }
652// if (IsNaN(tRow)){
653// cout <<"***ERROR tRow"<< endl;
654// }
655// if (IsNaN(tU)){
656// cout <<"***ERROR tU"<< endl;
657// }
658// if (IsNaN(tPhi)){
659// cout <<"***ERROR tPhi"<< endl;
660// }
661// // calculate crossing plane
662// tn.SetX(cos(tPhi));
663// tn.SetY(sin(tPhi));
664// tr.SetX(tRowRadius[ti]*cos(tPhi));
665// tr.SetY(tRowRadius[ti]*sin(tPhi));
666// // find crossing point
667// tLength = hel.PathLength(tr,tn);
668// if (IsNaN(tLength)){
669// cout <<"tLength loop 2nd NAN" << endl;
670// cout <<"padrow number= " << ti << " not reached" << endl;
671// tmpSect[ti]=-2;//sector
672// }
673// tPoint = hel.At(tLength);
674// tmpZ[ti] = tPoint.z();
675// tOutOfBound = TpcLocalTransform(tPoint,tSect,tRow,tmpU[ti],tPhi);
676// if (IsNaN(tSect)){
677// cout <<"***ERROR tSect 2"<< endl;
678// }
679// if (IsNaN(tRow)){
680// cout <<"***ERROR tRow 2"<< endl;
681// }
682// if (IsNaN(tmpU[ti])){
683// cout <<"***ERROR tmpU[ti] 2"<< endl;
684// }
685// if (IsNaN(tPhi)){
686// cout <<"***ERROR tPhi 2 "<< endl;
687// }
688// if(tOutOfBound || (tmpSect[ti] == tSect && tRow!=(ti+1))){
689// tmpSect[ti]=-2;
690// // cout << "missed once"<< endl;
691// }
692// else{
693// if(tmpSect[ti] != tSect){
694// // Try again on the other sector
695// tn.SetX(cos(tPhi));
696// tn.SetY(sin(tPhi));
697// tr.SetX(tRowRadius[ti]*cos(tPhi));
698// tr.SetY(tRowRadius[ti]*sin(tPhi));
699// // find crossing point
700// tLength = hel.PathLength(tr,tn);
701// tPoint = hel.At(tLength);
702// if (IsNaN(tLength)){
703// cout <<"tLength loop 3rd NAN" << endl;
704// cout <<"padrow number= "<< ti << " not reached" << endl;
705// tmpSect[ti]=-1;//sector
706// }
707// tmpZ[ti] = tPoint.z();
708// tmpSect[ti] = tSect;
709// tOutOfBound = TpcLocalTransform(tPoint,tSect,tRow,tmpU[ti],tPhi);
710// if (IsNaN(tSect)){
711// cout <<"***ERROR tSect 3"<< endl;
712// }
713// if (IsNaN(tRow)){
714// cout <<"***ERROR tRow 3"<< endl;
715// }
716// if (IsNaN(tmpU[ti])){
717// cout <<"***ERROR tmpU[ti] 3"<< endl;
718// }
719// if (IsNaN(tPhi)){
720// cout <<"***ERROR tPhi 3 "<< endl;
721// }
722// if(tOutOfBound || tSect!= tmpSect[ti] || tRow!=(ti+1)){
723// tmpSect[ti]=-1;
724// }
725// }
726// }
727// if (IsNaN(tmpSect[ti])){
728// cout << "*******************ERROR***************************" << endl;
729// cout <<"AliFemtoParticle--Fctn tmpSect=" << tmpSect[ti] << endl;
730// cout << "*******************ERROR***************************" << endl;
731// }
732// if (IsNaN(tmpU[ti])){
733// cout << "*******************ERROR***************************" << endl;
734// cout <<"AliFemtoParticle--Fctn tmpU=" << tmpU[ti] << endl;
735// cout << "*******************ERROR***************************" << endl;
736// }
737// if (IsNaN(tmpZ[ti])){
738// cout << "*******************ERROR***************************" << endl;
739// cout <<"AliFemtoParticle--Fctn tmpZ=" << tmpZ[ti] << endl;
740// cout << "*******************ERROR***************************" << endl;
741// }
742// // If padrow ti not reached all other beyond are not reached
743// // in this case Set sector to -1
744// if (tmpSect[ti]==-1){
745// for (int tj=ti; tj<45;tj++){
746// tmpSect[tj] = -1;
747// ti=45;
748// }
749// }
750// ti++;
751// if (ti<45){
752// candidates = hel.PathLength(tRowRadius[ti]);
753// tLength = (candidates.first > 0) ? candidates.first : candidates.second;}
754// }
755// }
756//_____________________
757const AliFemtoThreeVector& AliFemtoParticle::TpcV0PosExitPoint() const{
758 return fTpcV0PosExitPoint;
759}
760//_____________________
761const AliFemtoThreeVector& AliFemtoParticle::TpcV0PosEntrancePoint() const{
762 return fTpcV0PosEntrancePoint;
763}
764//______________________
765const AliFemtoThreeVector& AliFemtoParticle::TpcV0NegExitPoint() const{
766 return fTpcV0NegExitPoint;
767}
768//_____________________
769const AliFemtoThreeVector& AliFemtoParticle::TpcV0NegEntrancePoint() const{
770 return fTpcV0NegEntrancePoint;
771}
772//______________________