]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWG2/FEMTOSCOPY/AliFemto/AliFemtoPair.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Making the directory structure of AliFemtoUser flat. All files go into one common...
[u/mrichter/AliRoot.git] / PWG2 / FEMTOSCOPY / AliFemto / AliFemtoPair.cxx
CommitLineData
d0e92d9a 1///////////////////////////////////////////////////////////////////////////
2// //
3// AliFemtoPair: the Pair object is passed to the PairCuts for //
4// verification, and then to the AddRealPair and AddMixedPair methods of //
5// the Correlation Functions. It holds pair-specific variables like //
6// relative momenta and has links to the particles and tracks that form //
7// the pair. //
8// //
9///////////////////////////////////////////////////////////////////////////
10#include "AliFemtoPair.h"
11
12double AliFemtoPair::fgMaxDuInner = .8;
13double AliFemtoPair::fgMaxDzInner = 3.;
14double AliFemtoPair::fgMaxDuOuter = 1.4;
15double AliFemtoPair::fgMaxDzOuter = 3.2;
67427ff7 16
17
0215f606 18AliFemtoPair::AliFemtoPair() :
19 fTrack1(0), fTrack2(0),
20 fNonIdParNotCalculated(0),
21 fDKSide(0),
22 fDKOut(0),
23 fDKLong(0),
24 fCVK(0),
d0e92d9a 25 fKStarCalc(0),
0215f606 26 fNonIdParNotCalculatedGlobal(0),
27 fMergingParNotCalculated(0),
28 fWeightedAvSep(0),
29 fFracOfMergedRow(0),
30 fClosestRowAtDCA(0),
31 fMergingParNotCalculatedTrkV0Pos(0),
32 fFracOfMergedRowTrkV0Pos(0),
33 fClosestRowAtDCATrkV0Pos(0),
34 fMergingParNotCalculatedTrkV0Neg(0),
35 fFracOfMergedRowTrkV0Neg(0),
36 fClosestRowAtDCATrkV0Neg(0),
37 fMergingParNotCalculatedV0PosV0Neg(0),
38 fFracOfMergedRowV0PosV0Neg(0),
39 fClosestRowAtDCAV0PosV0Neg(0),
40 fMergingParNotCalculatedV0NegV0Pos(0),
41 fFracOfMergedRowV0NegV0Pos(0),
42 fClosestRowAtDCAV0NegV0Pos(0),
43 fMergingParNotCalculatedV0PosV0Pos(0),
44 fFracOfMergedRowV0PosV0Pos(0),
45 fClosestRowAtDCAV0PosV0Pos(0),
46 fMergingParNotCalculatedV0NegV0Neg(0),
47 fFracOfMergedRowV0NegV0Neg(0),
48 fClosestRowAtDCAV0NegV0Neg(0)
49{
d0e92d9a 50 // Default constructor
67427ff7 51 fTrack1 = 0;
52 fTrack2 = 0;
d0e92d9a 53 SetDefaultHalfFieldMergingPar();
67427ff7 54}
55
56AliFemtoPair::AliFemtoPair(AliFemtoParticle* a, AliFemtoParticle* b)
0215f606 57 : fTrack1(a), fTrack2(b),
58 fNonIdParNotCalculated(0),
59 fDKSide(0),
60 fDKOut(0),
61 fDKLong(0),
62 fCVK(0),
d0e92d9a 63 fKStarCalc(0),
0215f606 64 fNonIdParNotCalculatedGlobal(0),
65 fMergingParNotCalculated(0),
66 fWeightedAvSep(0),
67 fFracOfMergedRow(0),
68 fClosestRowAtDCA(0),
69 fMergingParNotCalculatedTrkV0Pos(0),
70 fFracOfMergedRowTrkV0Pos(0),
71 fClosestRowAtDCATrkV0Pos(0),
72 fMergingParNotCalculatedTrkV0Neg(0),
73 fFracOfMergedRowTrkV0Neg(0),
74 fClosestRowAtDCATrkV0Neg(0),
75 fMergingParNotCalculatedV0PosV0Neg(0),
76 fFracOfMergedRowV0PosV0Neg(0),
77 fClosestRowAtDCAV0PosV0Neg(0),
78 fMergingParNotCalculatedV0NegV0Pos(0),
79 fFracOfMergedRowV0NegV0Pos(0),
80 fClosestRowAtDCAV0NegV0Pos(0),
81 fMergingParNotCalculatedV0PosV0Pos(0),
82 fFracOfMergedRowV0PosV0Pos(0),
83 fClosestRowAtDCAV0PosV0Pos(0),
84 fMergingParNotCalculatedV0NegV0Neg(0),
85 fFracOfMergedRowV0NegV0Neg(0),
86 fClosestRowAtDCAV0NegV0Neg(0)
67427ff7 87{
d0e92d9a 88 // Construct a pair from two particles
89 SetDefaultHalfFieldMergingPar();
67427ff7 90}
91
d0e92d9a 92void AliFemtoPair::SetDefaultHalfFieldMergingPar(){
93 fgMaxDuInner = 3;
94 fgMaxDzInner = 4.;
95 fgMaxDuOuter = 4.;
96 fgMaxDzOuter = 6.;
67427ff7 97}
d0e92d9a 98void AliFemtoPair::SetDefaultFullFieldMergingPar(){
99 // Set default TPC merging parameters for STAR TPC
100 fgMaxDuInner = 0.8;
101 fgMaxDzInner = 3.;
102 fgMaxDuOuter = 1.4;
103 fgMaxDzOuter = 3.2;
67427ff7 104}
d0e92d9a 105void AliFemtoPair::SetMergingPar(double aMaxDuInner, double aMaxDzInner,
106 double aMaxDuOuter, double aMaxDzOuter)
107{
108 // Set TPC merging parameters for STAR TPC
109 fgMaxDuInner = aMaxDuInner;
110 fgMaxDzInner = aMaxDzInner;
111 fgMaxDuOuter = aMaxDuOuter;
112 fgMaxDzOuter = aMaxDzOuter;
67427ff7 113};
114
d0e92d9a 115AliFemtoPair::~AliFemtoPair() {
116 // Destructor
117/* no-op */
118}
67427ff7 119
0215f606 120AliFemtoPair::AliFemtoPair(const AliFemtoPair &aPair):
121 fTrack1(0), fTrack2(0),
122 fNonIdParNotCalculated(0),
123 fDKSide(0),
124 fDKOut(0),
125 fDKLong(0),
126 fCVK(0),
d0e92d9a 127 fKStarCalc(0),
0215f606 128 fNonIdParNotCalculatedGlobal(0),
129 fMergingParNotCalculated(0),
130 fWeightedAvSep(0),
131 fFracOfMergedRow(0),
132 fClosestRowAtDCA(0),
133 fMergingParNotCalculatedTrkV0Pos(0),
134 fFracOfMergedRowTrkV0Pos(0),
135 fClosestRowAtDCATrkV0Pos(0),
136 fMergingParNotCalculatedTrkV0Neg(0),
137 fFracOfMergedRowTrkV0Neg(0),
138 fClosestRowAtDCATrkV0Neg(0),
139 fMergingParNotCalculatedV0PosV0Neg(0),
140 fFracOfMergedRowV0PosV0Neg(0),
141 fClosestRowAtDCAV0PosV0Neg(0),
142 fMergingParNotCalculatedV0NegV0Pos(0),
143 fFracOfMergedRowV0NegV0Pos(0),
144 fClosestRowAtDCAV0NegV0Pos(0),
145 fMergingParNotCalculatedV0PosV0Pos(0),
146 fFracOfMergedRowV0PosV0Pos(0),
147 fClosestRowAtDCAV0PosV0Pos(0),
148 fMergingParNotCalculatedV0NegV0Neg(0),
149 fFracOfMergedRowV0NegV0Neg(0),
150 fClosestRowAtDCAV0NegV0Neg(0)
151{
d0e92d9a 152 // Copy constructor
0215f606 153 fTrack1 = aPair.fTrack1;
154 fTrack2 = aPair.fTrack2;
155
156 fNonIdParNotCalculated = aPair.fNonIdParNotCalculated;
157 fDKSide = aPair.fDKSide;
158 fDKOut = aPair.fDKOut;
159 fDKLong = aPair.fDKLong;
160 fCVK = aPair.fCVK;
d0e92d9a 161 fKStarCalc = aPair.fKStarCalc;
0215f606 162
163 fNonIdParNotCalculatedGlobal = aPair.fNonIdParNotCalculatedGlobal;
164
165 fMergingParNotCalculated = aPair.fMergingParNotCalculated;
166 fWeightedAvSep = aPair.fWeightedAvSep;
167 fFracOfMergedRow = aPair.fFracOfMergedRow;
168 fClosestRowAtDCA = aPair.fClosestRowAtDCA;
169
170 fMergingParNotCalculatedTrkV0Pos = aPair.fMergingParNotCalculatedTrkV0Pos;
171 fFracOfMergedRowTrkV0Pos = aPair.fFracOfMergedRowTrkV0Pos;
172 fClosestRowAtDCATrkV0Pos = aPair.fClosestRowAtDCATrkV0Pos;
173
174 fMergingParNotCalculatedTrkV0Neg = aPair.fMergingParNotCalculatedTrkV0Neg;
175 fFracOfMergedRowTrkV0Neg = aPair.fFracOfMergedRowTrkV0Neg;
176 fClosestRowAtDCATrkV0Neg = aPair.fClosestRowAtDCATrkV0Neg;
177
178 fMergingParNotCalculatedV0PosV0Neg = aPair.fMergingParNotCalculatedV0PosV0Neg;
179 fFracOfMergedRowV0PosV0Neg = aPair.fFracOfMergedRowV0PosV0Neg;
180 fClosestRowAtDCAV0PosV0Neg = aPair.fClosestRowAtDCAV0PosV0Neg;
181
182 fMergingParNotCalculatedV0NegV0Pos = aPair.fMergingParNotCalculatedV0NegV0Pos;
183 fFracOfMergedRowV0NegV0Pos = aPair.fFracOfMergedRowV0NegV0Pos;
184 fClosestRowAtDCAV0NegV0Pos = aPair.fClosestRowAtDCAV0NegV0Pos;
185
186 fMergingParNotCalculatedV0PosV0Pos = aPair.fMergingParNotCalculatedV0PosV0Pos;
187 fFracOfMergedRowV0PosV0Pos = aPair.fFracOfMergedRowV0PosV0Pos;
188 fClosestRowAtDCAV0PosV0Pos = aPair.fClosestRowAtDCAV0PosV0Pos;
189
190 fMergingParNotCalculatedV0NegV0Neg = aPair.fMergingParNotCalculatedV0NegV0Neg;
191 fFracOfMergedRowV0NegV0Neg = aPair.fFracOfMergedRowV0NegV0Neg;
192 fClosestRowAtDCAV0NegV0Neg = aPair.fClosestRowAtDCAV0NegV0Neg;
193}
194
195AliFemtoPair& AliFemtoPair::operator=(const AliFemtoPair &aPair)
196{
d0e92d9a 197 // Assignment operator
0215f606 198 if (this == &aPair)
199 return *this;
67427ff7 200
0215f606 201 fTrack1 = aPair.fTrack1;
202 fTrack2 = aPair.fTrack2;
203
204 fNonIdParNotCalculated = aPair.fNonIdParNotCalculated;
205 fDKSide = aPair.fDKSide;
206 fDKOut = aPair.fDKOut;
207 fDKLong = aPair.fDKLong;
208 fCVK = aPair.fCVK;
d0e92d9a 209 fKStarCalc = aPair.fKStarCalc;
0215f606 210
211 fNonIdParNotCalculatedGlobal = aPair.fNonIdParNotCalculatedGlobal;
212
213 fMergingParNotCalculated = aPair.fMergingParNotCalculated;
214 fWeightedAvSep = aPair.fWeightedAvSep;
215 fFracOfMergedRow = aPair.fFracOfMergedRow;
216 fClosestRowAtDCA = aPair.fClosestRowAtDCA;
217
218 fMergingParNotCalculatedTrkV0Pos = aPair.fMergingParNotCalculatedTrkV0Pos;
219 fFracOfMergedRowTrkV0Pos = aPair.fFracOfMergedRowTrkV0Pos;
220 fClosestRowAtDCATrkV0Pos = aPair.fClosestRowAtDCATrkV0Pos;
221
222 fMergingParNotCalculatedTrkV0Neg = aPair.fMergingParNotCalculatedTrkV0Neg;
223 fFracOfMergedRowTrkV0Neg = aPair.fFracOfMergedRowTrkV0Neg;
224 fClosestRowAtDCATrkV0Neg = aPair.fClosestRowAtDCATrkV0Neg;
225
226 fMergingParNotCalculatedV0PosV0Neg = aPair.fMergingParNotCalculatedV0PosV0Neg;
227 fFracOfMergedRowV0PosV0Neg = aPair.fFracOfMergedRowV0PosV0Neg;
228 fClosestRowAtDCAV0PosV0Neg = aPair.fClosestRowAtDCAV0PosV0Neg;
229
230 fMergingParNotCalculatedV0NegV0Pos = aPair.fMergingParNotCalculatedV0NegV0Pos;
231 fFracOfMergedRowV0NegV0Pos = aPair.fFracOfMergedRowV0NegV0Pos;
232 fClosestRowAtDCAV0NegV0Pos = aPair.fClosestRowAtDCAV0NegV0Pos;
233
234 fMergingParNotCalculatedV0PosV0Pos = aPair.fMergingParNotCalculatedV0PosV0Pos;
235 fFracOfMergedRowV0PosV0Pos = aPair.fFracOfMergedRowV0PosV0Pos;
236 fClosestRowAtDCAV0PosV0Pos = aPair.fClosestRowAtDCAV0PosV0Pos;
237
238 fMergingParNotCalculatedV0NegV0Neg = aPair.fMergingParNotCalculatedV0NegV0Neg;
239 fFracOfMergedRowV0NegV0Neg = aPair.fFracOfMergedRowV0NegV0Neg;
240 fClosestRowAtDCAV0NegV0Neg = aPair.fClosestRowAtDCAV0NegV0Neg;
241
242 return *this;
243}
67427ff7 244
245//_________________
d0e92d9a 246double AliFemtoPair::MInv() const
67427ff7 247{
d0e92d9a 248 // invariant mass
249 double tInvariantMass = abs(fTrack1->FourMomentum() + fTrack2->FourMomentum());
250 return (tInvariantMass);
67427ff7 251}
252//_________________
d0e92d9a 253double AliFemtoPair::KT() const
67427ff7 254{
d0e92d9a 255 // transverse momentum
67427ff7 256 double tmp =
257 (fTrack1->FourMomentum() + fTrack2->FourMomentum()).perp();
258 tmp *= .5;
259
260 return (tmp);
261}
262//_________________
d0e92d9a 263double AliFemtoPair::Rap() const
67427ff7 264{
265 // longitudinal pair rapidity : Y = 0.5 ::log( E1 + E2 + pz1 + pz2 / E1 + E2 - pz1 - pz2 )
266 double tmp = 0.5 * log (
267 (fTrack1->FourMomentum().e() + fTrack2->FourMomentum().e() + fTrack1->FourMomentum().z() + fTrack2->FourMomentum().z()) /
268 (fTrack1->FourMomentum().e() + fTrack2->FourMomentum().e() - fTrack1->FourMomentum().z() - fTrack2->FourMomentum().z())
269 ) ;
270 return (tmp);
271}
272//_________________
d0e92d9a 273double AliFemtoPair::EmissionAngle() const {
274 // emission angle
275 double pxTotal = this->FourMomentumSum().x();
276 double pyTotal = this->FourMomentumSum().y();
67427ff7 277 double angle = atan2(pyTotal,pxTotal)*180.0/3.1415926536;
278 if (angle<0.0) angle+=360.0;
279 return angle;
280}
281//_________________
282// get rid of ambiguously-named method fourMomentum() and replace it with
283// fourMomentumSum() and fourMomentumDiff() - mal 13feb2000
d0e92d9a 284AliFemtoLorentzVector AliFemtoPair::FourMomentumSum() const
67427ff7 285{
d0e92d9a 286 // total momentum
67427ff7 287 AliFemtoLorentzVector temp = fTrack1->FourMomentum()+fTrack2->FourMomentum();
288 return temp;
289}
d0e92d9a 290AliFemtoLorentzVector AliFemtoPair::FourMomentumDiff() const
67427ff7 291{
d0e92d9a 292 // momentum difference
67427ff7 293 AliFemtoLorentzVector temp = fTrack1->FourMomentum()-fTrack2->FourMomentum();
294 return temp;
295}
296//__________________________________
d0e92d9a 297void AliFemtoPair::QYKPCMS(double& qP, double& qT, double& q0) const
67427ff7 298{
d0e92d9a 299 // Yano-Koonin-Podgoretskii Parametrisation in CMS
67427ff7 300 ////
301 // calculate momentum difference in source rest frame (= lab frame)
302 ////
303 AliFemtoLorentzVector l1 = fTrack1->FourMomentum() ;
304 AliFemtoLorentzVector l2 = fTrack2->FourMomentum() ;
305 AliFemtoLorentzVector l ;
306 // random ordering of the particles
307 if ( rand()/(double)RAND_MAX > 0.50 )
308 { l = l1-l2 ; }
309 else
310 { l = l2-l1 ; } ;
311 // fill momentum differences into return variables
312 qP = l.z() ;
313 qT = l.vect().perp() ;
314 q0 = l.e() ;
315}
316//___________________________________
d0e92d9a 317void AliFemtoPair::QYKPLCMS(double& qP, double& qT, double& q0) const
67427ff7 318{
d0e92d9a 319 // Yano-Koonin-Podgoretskii Parametrisation in LCMS
67427ff7 320 ////
321 // calculate momentum difference in LCMS : frame where pz1 + pz2 = 0
322 ////
323 AliFemtoLorentzVector l1 = fTrack1->FourMomentum() ;
324 AliFemtoLorentzVector l2 = fTrack2->FourMomentum() ;
325 // determine beta to LCMS
326 double beta = (l1.z()+l2.z()) / (l1.e()+l2.e()) ;
327 double beta2 = beta*beta ;
328 // unfortunately STAR Class lib knows only boost(particle) not boost(beta) :(
329 // -> create particle with velocity beta and mass 1.0
330 // actually this is : dummyPz = ::sqrt( (dummyMass*dummyMass*beta2) / (1-beta2) ) ;
331 double dummyPz = ::sqrt( (beta2) / (1-beta2) ) ;
332 // boost in the correct direction
333 if (beta>0.0) { dummyPz = -dummyPz; } ;
334 // create dummy particle
335 AliFemtoLorentzVector l(0.0, 0.0, dummyPz) ;
336 double dummyMass = 1.0 ;
337 l.setE(l.vect().massHypothesis(dummyMass) );
338 // boost particles along the beam into a frame with velocity beta
339 AliFemtoLorentzVector l1boosted = l1.boost(l) ;
340 AliFemtoLorentzVector l2boosted = l2.boost(l) ;
341 // caculate the momentum difference with random ordering of the particle
342 if ( rand()/(double)RAND_MAX >0.50)
343 { l = l1boosted-l2boosted ; }
344 else
345 { l = l2boosted-l1boosted ;} ;
346 // fill momentum differences into return variables
347 qP = l.z() ;
348 qT = l.vect().perp() ;
349 q0 = l.e() ;
350}
351//___________________________________
352// Yano-Koonin-Podgoretskii Parametrisation in pair rest frame
d0e92d9a 353void AliFemtoPair::QYKPPF(double& qP, double& qT, double& q0) const
67427ff7 354{
355 ////
356 // calculate momentum difference in pair rest frame : frame where (pz1 + pz2, py1 + py2, px1 + px2) = (0,0,0)
357 ////
358 AliFemtoLorentzVector l1 = fTrack1->FourMomentum() ;
359 AliFemtoLorentzVector l2 = fTrack2->FourMomentum() ;
360 // the center of gravity of the pair travels with l
361 AliFemtoLorentzVector l = l1 + l2 ;
362 l = -l ;
363 l.setE(-l.e()) ;
364 // boost particles
365 AliFemtoLorentzVector l1boosted = l1.boost(l) ;
366 AliFemtoLorentzVector l2boosted = l2.boost(l) ;
367 // caculate the momentum difference with random ordering of the particle
368 if ( rand()/(double)RAND_MAX > 0.50)
369 { l = l1boosted-l2boosted ; }
370 else
371 { l = l2boosted-l1boosted ;} ;
372 // fill momentum differences into return variables
373 qP = l.z();
374 qT = l.vect().perp();
375 q0 = l.e();
376}
377//_________________
d0e92d9a 378double AliFemtoPair::QOutCMS() const
67427ff7 379{
d0e92d9a 380 // relative momentum out component in lab frame
67427ff7 381 AliFemtoThreeVector tmp1 = fTrack1->FourMomentum().vect();
382 AliFemtoThreeVector tmp2 = fTrack2->FourMomentum().vect();
383
384 double dx = tmp1.x() - tmp2.x();
385 double xt = tmp1.x() + tmp2.x();
386
387 double dy = tmp1.y() - tmp2.y();
388 double yt = tmp1.y() + tmp2.y();
389
390 double k1 = (::sqrt(xt*xt+yt*yt));
391 double k2 = (dx*xt+dy*yt);
392 double tmp = k2/k1;
393 return (tmp);
394}
395//_________________
d0e92d9a 396double AliFemtoPair::QSideCMS() const
67427ff7 397{
d0e92d9a 398 // relative momentum side component in lab frame
67427ff7 399 AliFemtoThreeVector tmp1 = fTrack1->FourMomentum().vect();
400 AliFemtoThreeVector tmp2 = fTrack2->FourMomentum().vect();
401
402 double x1 = tmp1.x(); double y1 = tmp1.y();
403 double x2 = tmp2.x(); double y2 = tmp2.y();
404
405 double xt = x1+x2; double yt = y1+y2;
406 double k1 = ::sqrt(xt*xt+yt*yt);
407
408 double tmp = 2.0*(x2*y1-x1*y2)/k1;
409 return (tmp);
410}
411
412//_________________________
d0e92d9a 413double AliFemtoPair::QLongCMS() const
67427ff7 414{
d0e92d9a 415 // relative momentum component in lab frame
67427ff7 416 AliFemtoLorentzVector tmp1 = fTrack1->FourMomentum();
417 AliFemtoLorentzVector tmp2 = fTrack2->FourMomentum();
418
419 double dz = tmp1.z() - tmp2.z();
420 double zz = tmp1.z() + tmp2.z();
421
422 double dt = tmp1.t() - tmp2.t();
423 double tt = tmp1.t() + tmp2.t();
424
425 double beta = zz/tt;
426 double gamma = 1.0/::sqrt(1.0 - beta*beta);
427
428 double temp = gamma*(dz - beta*dt);
429 return (temp);
430}
431
432//________________________________
d0e92d9a 433double AliFemtoPair::QOutPf() const
67427ff7 434{
d0e92d9a 435 // relative momentum out component in pair frame
436 AliFemtoLorentzVector tmp1 = fTrack1->FourMomentum();
437 AliFemtoLorentzVector tmp2 = fTrack2->FourMomentum();
438
439 double dt = tmp1.t() - tmp2.t();
440 double tt = tmp1.t() + tmp2.t();
441
442 double xt = tmp1.x() + tmp2.x();
443 double yt = tmp1.y() + tmp2.y();
444
445 double k1 = ::sqrt(xt*xt + yt*yt);
446 double bOut = k1/tt;
447 double gOut = 1.0/::sqrt(1.0 - bOut*bOut);
448
449 double temp = gOut*(this->QOutCMS() - bOut*dt);
450 return (temp);
67427ff7 451}
452
453//___________________________________
d0e92d9a 454double AliFemtoPair::QSidePf() const
67427ff7 455{
d0e92d9a 456 // relative momentum side component in pair frame
457
458 return(this->QSideCMS());
67427ff7 459}
460
461//___________________________________
462
d0e92d9a 463double AliFemtoPair::QLongPf() const
67427ff7 464{
d0e92d9a 465 // relative momentum long component in pair frame
466
467 return(this->QLongCMS());
67427ff7 468}
469
470//___________________________________
d0e92d9a 471double AliFemtoPair::QOutBf(double beta) const
67427ff7 472{
d0e92d9a 473 // relative momentum out component
474 return(this->QOutCMS());
67427ff7 475}
476
477//___________________________________
478
d0e92d9a 479double AliFemtoPair::QSideBf(double beta) const
67427ff7 480{
d0e92d9a 481 // relative momentum side component
482 return(this->QSideCMS());
67427ff7 483}
484
485//___________________________________
d0e92d9a 486double AliFemtoPair::QLongBf(double beta) const
67427ff7 487{
d0e92d9a 488 // relative momentum long component
67427ff7 489 AliFemtoLorentzVector tmp1 = fTrack1->FourMomentum();
490 AliFemtoLorentzVector tmp2 = fTrack2->FourMomentum();
491
492 double dz = tmp1.z() - tmp2.z();
493 double dt = tmp1.t() + tmp2.t();
494
495 double gamma = 1.0/::sqrt(1.0 - beta*beta);
496
497 double temp = gamma*(dz - beta*dt);
498 return (temp);
499}
500
d0e92d9a 501double AliFemtoPair::Quality() const {
502 // Calculate split quality of the pair
67427ff7 503 unsigned long mapMask0 = 0xFFFFFF00;
504 unsigned long mapMask1 = 0x1FFFFF;
505 unsigned long padRow1To24Track1 = fTrack1->TopologyMap(0) & mapMask0;
506 unsigned long padRow25To45Track1 = fTrack1->TopologyMap(1) & mapMask1;
507 unsigned long padRow1To24Track2 = fTrack2->TopologyMap(0) & mapMask0;
508 unsigned long padRow25To45Track2 = fTrack2->TopologyMap(1) & mapMask1;
509 // AND logic
510 unsigned long bothPads1To24 = padRow1To24Track1 & padRow1To24Track2;
511 unsigned long bothPads25To45 = padRow25To45Track1 & padRow25To45Track2;
512 // XOR logic
513 unsigned long onePad1To24 = padRow1To24Track1 ^ padRow1To24Track2;
514 unsigned long onePad25To45 = padRow25To45Track1 ^ padRow25To45Track2;
515 unsigned long bitI;
516 int ibits;
d0e92d9a 517 int tQuality = 0;
67427ff7 518 double normQual = 0.0;
d0e92d9a 519 int tMaxQuality = fTrack1->NumberOfHits() + fTrack2->NumberOfHits();
67427ff7 520 for (ibits=8;ibits<=31;ibits++) {
521 bitI = 0;
522 bitI |= 1UL<<(ibits);
523 if ( onePad1To24 & bitI ) {
d0e92d9a 524 tQuality++;
67427ff7 525 continue;
526 }
527 else{
d0e92d9a 528 if ( bothPads1To24 & bitI ) tQuality--;
67427ff7 529 }
530 }
531 for (ibits=0;ibits<=20;ibits++) {
532 bitI = 0;
533 bitI |= 1UL<<(ibits);
534 if ( onePad25To45 & bitI ) {
d0e92d9a 535 tQuality++;
67427ff7 536 continue;
537 }
538 else{
d0e92d9a 539 if ( bothPads25To45 & bitI ) tQuality--;
67427ff7 540 }
541 }
d0e92d9a 542 normQual = (double)tQuality/( (double) tMaxQuality );
67427ff7 543 return ( normQual );
544
545}
546
d0e92d9a 547double AliFemtoPair::Quality2() const {
548 // second implementation of split quality
67427ff7 549 unsigned long mapMask0 = 0xFFFFFF00;
550 unsigned long mapMask1 = 0x1FFFFF;
551 unsigned long padRow1To24Track1 = fTrack1->TopologyMap(0) & mapMask0;
552 unsigned long padRow25To45Track1 = fTrack1->TopologyMap(1) & mapMask1;
553 unsigned long padRow1To24Track2 = fTrack2->TopologyMap(0) & mapMask0;
554 unsigned long padRow25To45Track2 = fTrack2->TopologyMap(1) & mapMask1;
555
556 // AND logic
557 //unsigned long bothPads1To24 = padRow1To24Track1 & padRow1To24Track2;
558 //unsigned long bothPads25To45 = padRow25To45Track1 & padRow25To45Track2;
559
560 // XOR logic
561 unsigned long onePad1To24 = padRow1To24Track1 ^ padRow1To24Track2;
562 unsigned long onePad25To45 = padRow25To45Track1 ^ padRow25To45Track2;
563 unsigned long bitI;
564 int ibits;
d0e92d9a 565 int tQuality = 0;
67427ff7 566 double normQual = 0.0;
d0e92d9a 567 int tMaxQuality = fTrack1->NumberOfHits() + fTrack2->NumberOfHits();
67427ff7 568 for (ibits=8;ibits<=31;ibits++) {
569 bitI = 0;
570 bitI |= 1UL<<(ibits);
571 if ( onePad1To24 & bitI ) {
d0e92d9a 572 tQuality++;
67427ff7 573 continue;
574 }
575 //else{
d0e92d9a 576 //if ( bothPads1To24 & bitI ) tQuality--;
67427ff7 577 //}
578 }
579 for (ibits=0;ibits<=20;ibits++) {
580 bitI = 0;
581 bitI |= 1UL<<(ibits);
582 if ( onePad25To45 & bitI ) {
d0e92d9a 583 tQuality++;
67427ff7 584 continue;
585 }
586 //else{
d0e92d9a 587 //if ( bothPads25To45 & bitI ) tQuality--;
67427ff7 588 //}
589 }
d0e92d9a 590 normQual = (double)tQuality/( (double) tMaxQuality );
67427ff7 591 return ( normQual );
592
593}
594
595
596double AliFemtoPair::NominalTpcExitSeparation() const {
d0e92d9a 597 // separation at exit from STAR TPC
67427ff7 598 AliFemtoThreeVector diff = fTrack1->NominalTpcExitPoint() - fTrack2->NominalTpcExitPoint();
599 return (diff.mag());
600}
601
602double AliFemtoPair::NominalTpcEntranceSeparation() const {
d0e92d9a 603 // separation at entrance to STAR TPC
67427ff7 604 AliFemtoThreeVector diff = fTrack1->NominalTpcEntrancePoint() - fTrack2->NominalTpcEntrancePoint();
605 return (diff.mag());
606}
607
608double AliFemtoPair::NominalTpcAverageSeparation() const {
d0e92d9a 609 // average separation in STAR TPC
67427ff7 610 AliFemtoThreeVector diff;
d0e92d9a 611 double tAveSep = 0.0;
67427ff7 612 int ipt = 0;
613 if (fTrack1->fNominalPosSample && fTrack2->fNominalPosSample){
614 while (fabs(fTrack1->fNominalPosSample[ipt].x())<9999. &&
615 fabs(fTrack1->fNominalPosSample[ipt].y())<9999. &&
616 fabs(fTrack1->fNominalPosSample[ipt].z())<9999. &&
617 fabs(fTrack2->fNominalPosSample[ipt].x())<9999. &&
618 fabs(fTrack2->fNominalPosSample[ipt].y())<9999. &&
619 fabs(fTrack2->fNominalPosSample[ipt].z())<9999. &&
620 ipt<11
621 ){
622 // for (int ipt=0; ipt<11; ipt++){
623 diff = fTrack1->fNominalPosSample[ipt] - fTrack2->fNominalPosSample[ipt];
624 ipt++;
d0e92d9a 625 tAveSep += diff.mag();
67427ff7 626 }
d0e92d9a 627 tAveSep = tAveSep/(ipt+1.);
628 return (tAveSep);}
67427ff7 629 else return -1;
630}
631
632double AliFemtoPair::OpeningAngle() const {
d0e92d9a 633 // opening angle
67427ff7 634 return 57.296* fTrack1->FourMomentum().vect().angle( fTrack2->FourMomentum().vect() );
635// AliFemtoThreeVector p1 = fTrack1->FourMomentum().vect();
636// AliFemtoThreeVector p2 = fTrack2->FourMomentum().vect();
637// return 57.296*(p1.phi()-p2.phi());
638// //double dAngInv = 57.296*acos((p1.dot(p2))/(p1.mag()*p2.mag()));
639// //return (dAngInv);
640}
641//_________________
642
643
644double AliFemtoPair::KStarFlipped() const {
d0e92d9a 645 // kstar with sign flipped
67427ff7 646 AliFemtoLorentzVector tP1 = fTrack1->FourMomentum();
647
648 AliFmThreeVectorD qwe = tP1.vect();
649 qwe *= -1.; // flip it
650 tP1.setVect(qwe);
651
652 AliFemtoLorentzVector tSum = (tP1+fTrack2->FourMomentum());
653 double tMass = abs(tSum);
654 AliFmThreeVectorD tGammaBeta = (1./tMass)*tSum.vect();
655 double tGamma = tSum.e()/tMass;
656 AliFmThreeVectorD tLongMom = ((tP1.vect()*tGammaBeta)/
657 (tGammaBeta*tGammaBeta))*tGammaBeta;
658 AliFmLorentzVectorD tK(tGamma*tP1.e() - tP1.vect()*tGammaBeta,
659 tP1.vect() + (tGamma-1.)*tLongMom - tP1.e()*tGammaBeta);
660//VP tP1.vect() *= -1.; // unflip it
661 return tK.vect().mag();
662}
663
664//double AliFemtoPair::CVK() const{
665//const AliFemtoLorentzVector& tP1 = fTrack1->FourMomentum();
666//AliFemtoLorentzVector tSum = (tP1+fTrack2->FourMomentum());
667//double tMass = abs(tSum);
668//AliFmThreeVectorD tGammaBeta = (1./tMass)*tSum.vect();
669//double tGamma = tSum.e()/tMass;
670//AliFmThreeVectorD tLongMom = ((tP1.vect()*tGammaBeta)/
671// (tGammaBeta*tGammaBeta))*tGammaBeta;
672//AliFmLorentzVectorD tK(tGamma*tP1.e() - tP1.vect()*tGammaBeta,
673// tP1.vect() + (tGamma-1.)*tLongMom - tP1.e()*tGammaBeta);
674//return (tK.vect())*tGammaBeta/tK.vect().magnitude()/tGammaBeta.magnitude();
675//}
676
677double AliFemtoPair::CVKFlipped() const{
d0e92d9a 678 // CVK with sign flipped
67427ff7 679 AliFemtoLorentzVector tP1 = fTrack1->FourMomentum();
680 AliFmThreeVectorD qwe = tP1.vect();
681 qwe *= -1.; // flip it
682 tP1.setVect(qwe);
683
684 AliFemtoLorentzVector tSum = (tP1+fTrack2->FourMomentum());
685 double tMass = abs(tSum);
686 AliFmThreeVectorD tGammaBeta = (1./tMass)*tSum.vect();
687 double tGamma = tSum.e()/tMass;
688 AliFmThreeVectorD tLongMom = ((tP1.vect()*tGammaBeta)/
689 (tGammaBeta*tGammaBeta))*tGammaBeta;
690 AliFmLorentzVectorD tK(tGamma*tP1.e() - tP1.vect()*tGammaBeta,
691 tP1.vect() + (tGamma-1.)*tLongMom - tP1.e()*tGammaBeta);
692//VP tP1.vect() *= -1.; // unflip it
693 return (tK.vect())*tGammaBeta/tGamma;
694}
695
d0e92d9a 696double AliFemtoPair::PInv() const{
697 // invariant total momentum
67427ff7 698 AliFemtoLorentzVector tP1 = fTrack1->FourMomentum();
699 AliFemtoLorentzVector tP2 = fTrack2->FourMomentum();
700 double tP = (tP1.px()+tP2.px())*(tP1.px()+tP2.px())+
701 (tP1.py()+tP2.py())*(tP1.py()+tP2.py())+
702 (tP1.pz()+tP2.pz())*(tP1.pz()+tP2.pz())-
703 (tP1.e() -tP2.e() )*(tP1.e() -tP2.e() );
704 return ::sqrt(fabs(tP));
705}
706
d0e92d9a 707double AliFemtoPair::QInvFlippedXY() const{
708 // qinv with X and Y flipped
67427ff7 709 AliFemtoLorentzVector tP1 = fTrack1->FourMomentum();
710 tP1.setX(-1.*tP1.x());
711 tP1.setY(-1.*tP1.y());
712 AliFemtoLorentzVector tDiff = (tP1-fTrack2->FourMomentum());
713 return ( -1.* tDiff.m());
714}
715
d0e92d9a 716void AliFemtoPair::CalcNonIdPar() const{ // fortran like function! faster?
717 // Calculate generalized relative mometum
718 // Use this instead of qXYZ() function when calculating
719 // anything for non-identical particles
67427ff7 720 fNonIdParNotCalculated=0;
721 double px1 = fTrack1->FourMomentum().vect().x();
722 double py1 = fTrack1->FourMomentum().vect().y();
723 double pz1 = fTrack1->FourMomentum().vect().z();
724 double pE1 = fTrack1->FourMomentum().e();
d0e92d9a 725 double tParticle1Mass = ::sqrt(pE1*pE1 - px1*px1 - py1*py1 - pz1*pz1);
67427ff7 726 double px2 = fTrack2->FourMomentum().vect().x();
727 double py2 = fTrack2->FourMomentum().vect().y();
728 double pz2 = fTrack2->FourMomentum().vect().z();
729 double pE2 = fTrack2->FourMomentum().e();
d0e92d9a 730 double tParticle2Mass = ::sqrt(pE2*pE2 - px2*px2 - py2*py2 - pz2*pz2);
67427ff7 731
d0e92d9a 732 double tPx = px1+px2;
733 double tPy = py1+py2;
734 double tPz = pz1+pz2;
735 double tPE = pE1+pE2;
67427ff7 736
d0e92d9a 737 double tPtrans = tPx*tPx + tPy*tPy;
738 double tMtrans = tPE*tPE - tPz*tPz;
739 double tPinv = ::sqrt(tMtrans - tPtrans);
740 tMtrans = ::sqrt(tMtrans);
741 tPtrans = ::sqrt(tPtrans);
67427ff7 742
d0e92d9a 743 double tQinvL = (pE1-pE2)*(pE1-pE2) - (px1-px2)*(px1-px2) -
67427ff7 744 (py1-py2)*(py1-py2) - (pz1-pz2)*(pz1-pz2);
745
d0e92d9a 746 double tQ = (tParticle1Mass*tParticle1Mass - tParticle2Mass*tParticle2Mass)/tPinv;
747 tQ = sqrt ( tQ*tQ - tQinvL);
67427ff7 748
d0e92d9a 749 fKStarCalc = tQ/2;
67427ff7 750
751 // ad 1) go to LCMS
d0e92d9a 752 double beta = tPz/tPE;
753 double gamma = tPE/tMtrans;
67427ff7 754
755 double pz1L = gamma * (pz1 - beta * pE1);
756 double pE1L = gamma * (pE1 - beta * pz1);
757
758 // fill histogram for beam projection ( z - axis )
759 fDKLong = pz1L;
760
d0e92d9a 761 // ad 2) rotation px -> tPt
762 double px1R = (px1*tPx + py1*tPy)/tPtrans;
763 double py1R = (-px1*tPy + py1*tPx)/tPtrans;
67427ff7 764
765 //fill histograms for side projection ( y - axis )
766 fDKSide = py1R;
767
768 // ad 3) go from LCMS to CMS
d0e92d9a 769 beta = tPtrans/tMtrans;
770 gamma = tMtrans/tPinv;
67427ff7 771
772 double px1C = gamma * (px1R - beta * pE1L);
773
774 // fill histogram for out projection ( x - axis )
775 fDKOut = px1C;
776
d0e92d9a 777 fCVK = (fDKOut*tPtrans + fDKLong*tPz)/fKStarCalc/::sqrt(tPtrans*tPtrans+tPz*tPz);
67427ff7 778}
779
780
781/*void AliFemtoPair::calcNonIdParGlobal() const{ // fortran like function! faster?
782 fNonIdParNotCalculatedGlobal=0;
783 double px1 = fTrack1->Track()->PGlobal().x();
784 double py1 = fTrack1->Track()->PGlobal().y();
785 double pz1 = fTrack1->Track()->PGlobal().z();
d0e92d9a 786 double tParticle1Mass = fTrack1->FourMomentum().m2();
787 double pE1 = ::sqrt(tParticle1Mass + px1*px1 + py1*py1 + pz1*pz1);
788 tParticle1Mass = ::sqrt(tParticle1Mass);
67427ff7 789
790 double px2 = fTrack2->Track()->PGlobal().x();
791 double py2 = fTrack2->Track()->PGlobal().y();
792 double pz2 = fTrack2->Track()->PGlobal().z();
d0e92d9a 793 double tParticle2Mass = fTrack2->FourMomentum().m2();
794 double pE2 = ::sqrt(tParticle2Mass + px2*px2 + py2*py2 + pz2*pz2);
795 tParticle2Mass = ::sqrt(tParticle2Mass);
67427ff7 796
797 double Px = px1+px2;
798 double Py = py1+py2;
799 double Pz = pz1+pz2;
800 double PE = pE1+pE2;
801
802 double Ptrans = Px*Px + Py*Py;
803 double Mtrans = PE*PE - Pz*Pz;
804 double Pinv = ::sqrt(Mtrans - Ptrans);
805 Mtrans = ::sqrt(Mtrans);
806 Ptrans = ::sqrt(Ptrans);
807
808 double QinvL = (pE1-pE2)*(pE1-pE2) - (px1-px2)*(px1-px2) -
809 (py1-py2)*(py1-py2) - (pz1-pz2)*(pz1-pz2);
810
d0e92d9a 811 double Q = (tParticle1Mass*tParticle1Mass - tParticle2Mass*tParticle2Mass)/Pinv;
67427ff7 812 Q = sqrt ( Q*Q - QinvL);
813
814 kStarCalcGlobal = Q/2;
815
816 // ad 1) go to LCMS
817 double beta = Pz/PE;
818 double gamma = PE/Mtrans;
819
820 double pz1L = gamma * (pz1 - beta * pE1);
821 double pE1L = gamma * (pE1 - beta * pz1);
822
823 // fill histogram for beam projection ( z - axis )
824 fDKLongGlobal = pz1L;
825
826 // ad 2) rotation px -> Pt
827 double px1R = (px1*Px + py1*Py)/Ptrans;
828 double py1R = (-px1*Py + py1*Px)/Ptrans;
829
830 //fill histograms for side projection ( y - axis )
831 fDKSideGlobal = py1R;
832
833 // ad 3) go from LCMS to CMS
834 beta = Ptrans/Mtrans;
835 gamma = Mtrans/Pinv;
836
837 double px1C = gamma * (px1R - beta * pE1L);
838
839 // fill histogram for out projection ( x - axis )
840 fDKOutGlobal = px1C;
841
842 fCVKGlobal = (fDKOutGlobal*Ptrans + fDKLongGlobal*Pz)/
843 kStarCalcGlobal/::sqrt(Ptrans*Ptrans+Pz*Pz);
844}*/
845
846
847
d0e92d9a 848double AliFemtoPair::DcaInsideTpc() const{
849 // dcs inside the STAR TPC
67427ff7 850 double tMinDist=NominalTpcEntranceSeparation();
851 double tExit = NominalTpcExitSeparation();
852 tMinDist = (tExit>tMinDist) ? tMinDist : tExit;
853 double tInsideDist;
854 //tMinDist = 999.;
855
856 double rMin = 60.;
857 double rMax = 190.;
858 const AliFmPhysicalHelixD& tHelix1 = fTrack1->Helix();
859 const AliFmPhysicalHelixD& tHelix2 = fTrack2->Helix();
860 // --- One is a line and other one a helix
861 //if (tHelix1.mSingularity != tHelix2.mSingularity) return -999.;
862 // --- 2 lines : don't care right now
863 //if (tHelix1.mSingularity) return -999.;
864 // --- 2 helix
d0e92d9a 865 double dx = tHelix2.XCenter() - tHelix1.XCenter();
866 double dy = tHelix2.YCenter() - tHelix1.YCenter();
67427ff7 867 double dd = ::sqrt(dx*dx + dy*dy);
d0e92d9a 868 double r1 = 1/tHelix1.Curvature();
869 double r2 = 1/tHelix2.Curvature();
67427ff7 870 double cosAlpha = (r1*r1 + dd*dd - r2*r2)/(2*r1*dd);
871
872 double x, y, r;
873 double s;
874 if (fabs(cosAlpha) < 1) { // two solutions
875 double sinAlpha = sin(acos(cosAlpha));
d0e92d9a 876 x = tHelix1.XCenter() + r1*(cosAlpha*dx - sinAlpha*dy)/dd;
877 y = tHelix1.YCenter() + r1*(sinAlpha*dx + cosAlpha*dy)/dd;
67427ff7 878 r = ::sqrt(x*x+y*y);
879 if( r > rMin && r < rMax &&
880 fabs(atan2(y,x)-fTrack1->NominalTpcEntrancePoint().phi())< 0.5
881 ){ // first solution inside
d0e92d9a 882 s = tHelix1.PathLength(x, y);
883 tInsideDist=tHelix2.Distance(tHelix1.At(s));
67427ff7 884 if(tInsideDist<tMinDist) tMinDist = tInsideDist;
885 }
886 else{
d0e92d9a 887 x = tHelix1.XCenter() + r1*(cosAlpha*dx + sinAlpha*dy)/dd;
888 y = tHelix1.YCenter() + r1*(cosAlpha*dy - sinAlpha*dx)/dd;
67427ff7 889 r = ::sqrt(x*x+y*y);
890 if( r > rMin && r < rMax &&
891 fabs(atan2(y,x)-fTrack1->NominalTpcEntrancePoint().phi())< 0.5
892 ) { // second solution inside
d0e92d9a 893 s = tHelix1.PathLength(x, y);
894 tInsideDist=tHelix2.Distance(tHelix1.At(s));
67427ff7 895 if(tInsideDist<tMinDist) tMinDist = tInsideDist;
896 }
897 }
898 }
899 return tMinDist;
900}
901
d0e92d9a 902void AliFemtoPair::CalcMergingPar() const{
903 // Calculate merging factor for the pair in STAR TPC
67427ff7 904 fMergingParNotCalculated=0;
905
906 double tDu, tDz;
907 int tN = 0;
908 fFracOfMergedRow = 0.;
909 fWeightedAvSep =0.;
910 double tDist;
911 double tDistMax = 200.;
912 for(int ti=0 ; ti<45 ; ti++){
913 if(fTrack1->fSect[ti]==fTrack2->fSect[ti] && fTrack1->fSect[ti]!=-1){
914 tDu = fabs(fTrack1->fU[ti]-fTrack2->fU[ti]);
915 tDz = fabs(fTrack1->fZ[ti]-fTrack2->fZ[ti]);
916 tN++;
917 if(ti<13){
d0e92d9a 918 fFracOfMergedRow += (tDu<fgMaxDuInner && tDz<fgMaxDzInner);
919 tDist = ::sqrt(tDu*tDu/fgMaxDuInner/fgMaxDuInner+
920 tDz*tDz/fgMaxDzInner/fgMaxDzInner);
921 //fFracOfMergedRow += (tDu<fgMaxDuInner && tDz<fgMaxDzInner);
67427ff7 922 }
923 else{
d0e92d9a 924 fFracOfMergedRow += (tDu<fgMaxDuOuter && tDz<fgMaxDzOuter);
925 tDist = ::sqrt(tDu*tDu/fgMaxDuOuter/fgMaxDuOuter+
926 tDz*tDz/fgMaxDzOuter/fgMaxDzOuter);
927 //fFracOfMergedRow += (tDu<fgMaxDuOuter && tDz<fgMaxDzOuter);
67427ff7 928 }
929 if(tDist<tDistMax){
930 fClosestRowAtDCA = ti+1;
931 tDistMax = tDist;
932 }
933 fWeightedAvSep += tDist;
934 }
935 }
936 if(tN>0){
937 fWeightedAvSep /= tN;
938 fFracOfMergedRow /= tN;
939 }
940 else{
941 fClosestRowAtDCA = -1;
942 fFracOfMergedRow = -1.;
943 fWeightedAvSep = -1.;
944 }
945}
d0e92d9a 946double AliFemtoPair::TpcExitSeparationTrackV0Pos() const {
67427ff7 947//________________V0 daughters exit/entrance/average separation calc.
948//_______1st part is a track 2nd is a V0 considering Pos daughter
d0e92d9a 949
67427ff7 950 AliFemtoThreeVector diff = fTrack1->NominalTpcExitPoint() - fTrack2->TpcV0PosExitPoint();
951 return (diff.mag());
952}
953
954double AliFemtoPair::TpcEntranceSeparationTrackV0Pos() const {
d0e92d9a 955//________________V0 daughters exit/entrance/average separation calc.
956//_______1st part is a track 2nd is a V0 considering Pos daughter
67427ff7 957 AliFemtoThreeVector diff = fTrack1->NominalTpcEntrancePoint() - fTrack2->TpcV0PosEntrancePoint();
958 return (diff.mag());
959}
960
961double AliFemtoPair::TpcAverageSeparationTrackV0Pos() const {
d0e92d9a 962//________________V0 daughters exit/entrance/average separation calc.
963//_______1st part is a track 2nd is a V0 considering Pos daughter
67427ff7 964 AliFemtoThreeVector diff;
d0e92d9a 965 double tAveSep = 0.0;
67427ff7 966 int ipt = 0;
967 if (fTrack1->fNominalPosSample && fTrack2->fNominalPosSample){
968 while (fabs(fTrack1->fNominalPosSample[ipt].x())<9999. &&
969 fabs(fTrack1->fNominalPosSample[ipt].y())<9999. &&
970 fabs(fTrack1->fNominalPosSample[ipt].z())<9999. &&
971 fabs(fTrack2->fNominalPosSample[ipt].x())<9999. &&
972 fabs(fTrack2->fNominalPosSample[ipt].y())<9999. &&
973 fabs(fTrack2->fNominalPosSample[ipt].z())<9999. &&
974 (ipt<11)
975 ){
976 diff = fTrack1->fNominalPosSample[ipt] - fTrack2->fNominalPosSample[ipt];
977 ipt++;
d0e92d9a 978 tAveSep += diff.mag();
67427ff7 979 }
d0e92d9a 980 tAveSep = tAveSep/(ipt+1.);
981 return (tAveSep);}
67427ff7 982 else return -1;
983}
67427ff7 984double AliFemtoPair::TpcExitSeparationTrackV0Neg() const {
d0e92d9a 985//_______1st part is a track 2nd is a V0 considering Neg daughter
67427ff7 986 AliFemtoThreeVector diff = fTrack1->NominalTpcExitPoint() - fTrack2->TpcV0NegExitPoint();
987 return (diff.mag());
988}
989
990double AliFemtoPair::TpcEntranceSeparationTrackV0Neg() const {
d0e92d9a 991//_______1st part is a track 2nd is a V0 considering Neg daughter
67427ff7 992 AliFemtoThreeVector diff = fTrack1->NominalTpcEntrancePoint() - fTrack2->TpcV0NegEntrancePoint();
993 return (diff.mag());
994}
995
996double AliFemtoPair::TpcAverageSeparationTrackV0Neg() const {
d0e92d9a 997//_______1st part is a track 2nd is a V0 considering Neg daughter
67427ff7 998 AliFemtoThreeVector diff;
d0e92d9a 999 double tAveSep = 0.0;
67427ff7 1000 int ipt = 0;
1001 if (fTrack1->fNominalPosSample && fTrack2->fTpcV0NegPosSample){
1002 while (fabs(fTrack1->fNominalPosSample[ipt].x())<9999. &&
1003 fabs(fTrack1->fNominalPosSample[ipt].y())<9999. &&
1004 fabs(fTrack1->fNominalPosSample[ipt].z())<9999. &&
1005 fabs(fTrack2->fTpcV0NegPosSample[ipt].x())<9999. &&
1006 fabs(fTrack2->fTpcV0NegPosSample[ipt].y())<9999. &&
1007 fabs(fTrack2->fTpcV0NegPosSample[ipt].z())<9999. &&
1008 (ipt<11)
1009 ){
1010 diff = fTrack1->fNominalPosSample[ipt] - fTrack2->fTpcV0NegPosSample[ipt];
1011 ipt++;
d0e92d9a 1012 tAveSep += diff.mag();
67427ff7 1013 }
d0e92d9a 1014 tAveSep = tAveSep/(ipt+1.);
1015 return (tAveSep);}
67427ff7 1016 else return -1;
1017}
1018
67427ff7 1019double AliFemtoPair::TpcExitSeparationV0PosV0Pos() const {
d0e92d9a 1020//_______1st part is a V0 considering Pos daughter 2nd is a V0 considering Pos daughter
67427ff7 1021 AliFemtoThreeVector diff = fTrack1->TpcV0PosExitPoint() - fTrack2->TpcV0PosExitPoint();
1022 return (diff.mag());
1023}
1024
1025double AliFemtoPair::TpcEntranceSeparationV0PosV0Pos() const {
d0e92d9a 1026//_______1st part is a V0 considering Pos daughter 2nd is a V0 considering Pos daughter
67427ff7 1027 AliFemtoThreeVector diff = fTrack1->TpcV0PosEntrancePoint() - fTrack2->TpcV0PosEntrancePoint();
1028 return (diff.mag());
1029}
1030double AliFemtoPair::TpcAverageSeparationV0PosV0Pos() const {
d0e92d9a 1031//_______1st part is a V0 considering Pos daughter 2nd is a V0 considering Pos daughter
67427ff7 1032 AliFemtoThreeVector diff;
d0e92d9a 1033 double tAveSep = 0.0;
67427ff7 1034 int ipt=0;
1035 if (fTrack1->fNominalPosSample && (fTrack2->fNominalPosSample)){
1036 while ((fabs(fTrack1->fNominalPosSample[ipt].x())<9999.) &&
1037 (fabs(fTrack1->fNominalPosSample[ipt].y())<9999.) &&
1038 (fabs(fTrack1->fNominalPosSample[ipt].z())<9999.) &&
1039 (fabs(fTrack2->fNominalPosSample[ipt].x())<9999.) &&
1040 (fabs(fTrack2->fNominalPosSample[ipt].y())<9999.) &&
1041 (fabs(fTrack2->fNominalPosSample[ipt].z())<9999.) &&
1042 (ipt<11)
1043 ){
1044 diff = fTrack1->fNominalPosSample[ipt] - fTrack2->fNominalPosSample[ipt];
1045 ipt++;
d0e92d9a 1046 tAveSep += diff.mag();
67427ff7 1047 }
d0e92d9a 1048 tAveSep = tAveSep/(ipt+1);
1049 return (tAveSep);}
67427ff7 1050 else return -1;
1051}
1052
67427ff7 1053double AliFemtoPair::TpcExitSeparationV0PosV0Neg() const {
d0e92d9a 1054//_______1st part is a V0 considering Pos daughter 2nd is a V0 considering Neg daughter
67427ff7 1055 AliFemtoThreeVector diff = fTrack1->TpcV0PosExitPoint() - fTrack2->TpcV0NegExitPoint();
1056 return (diff.mag());
1057}
1058
1059double AliFemtoPair::TpcEntranceSeparationV0PosV0Neg() const {
d0e92d9a 1060//_______1st part is a V0 considering Pos daughter 2nd is a V0 considering Neg daughter
67427ff7 1061 AliFemtoThreeVector diff = fTrack1->TpcV0PosEntrancePoint() - fTrack2->TpcV0NegEntrancePoint();
1062 return (diff.mag());
1063}
1064double AliFemtoPair::TpcAverageSeparationV0PosV0Neg() const {
d0e92d9a 1065//_______1st part is a V0 considering Pos daughter 2nd is a V0 considering Neg daughter
67427ff7 1066 AliFemtoThreeVector diff;
d0e92d9a 1067 double tAveSep = 0.0;
67427ff7 1068 int ipt = 0;
1069 if (fTrack1->fNominalPosSample && fTrack2->fTpcV0NegPosSample){
1070 while (fabs(fTrack1->fNominalPosSample[ipt].x())<9999. &&
1071 fabs(fTrack1->fNominalPosSample[ipt].y())<9999. &&
1072 fabs(fTrack1->fNominalPosSample[ipt].z())<9999. &&
1073 fabs(fTrack2->fTpcV0NegPosSample[ipt].x())<9999. &&
1074 fabs(fTrack2->fTpcV0NegPosSample[ipt].y())<9999. &&
1075 fabs(fTrack2->fTpcV0NegPosSample[ipt].z())<9999. &&
1076 (ipt<11)
1077 ){
1078 diff = fTrack1->fNominalPosSample[ipt] - fTrack2->fTpcV0NegPosSample[ipt];
1079 ipt++;
d0e92d9a 1080 tAveSep += diff.mag();
67427ff7 1081 }
d0e92d9a 1082 tAveSep = tAveSep/(ipt+1.);
1083 return (tAveSep);}
67427ff7 1084 else return -1;
1085}
d0e92d9a 1086double AliFemtoPair::TpcExitSeparationV0NegV0Pos() const {
67427ff7 1087//_______1st part is a V0 considering Neg daughter 2nd is a V0 considering Pos daughter
1088// this is to check the upper case
67427ff7 1089 AliFemtoThreeVector diff = fTrack1->TpcV0NegExitPoint() - fTrack2->TpcV0PosExitPoint();
1090 return (diff.mag());
1091}
1092
1093double AliFemtoPair::TpcEntranceSeparationV0NegV0Pos() const {
d0e92d9a 1094//_______1st part is a V0 considering Neg daughter 2nd is a V0 considering Pos daughter
1095// this is to check the upper case
67427ff7 1096 AliFemtoThreeVector diff = fTrack1->TpcV0NegEntrancePoint() - fTrack2->TpcV0PosEntrancePoint();
1097 return (diff.mag());
1098}
1099double AliFemtoPair::TpcAverageSeparationV0NegV0Pos() const {
d0e92d9a 1100//_______1st part is a V0 considering Neg daughter 2nd is a V0 considering Pos daughter
1101// this is to check the upper case
67427ff7 1102 AliFemtoThreeVector diff;
d0e92d9a 1103 double tAveSep = 0.0;
67427ff7 1104 int ipt = 0;
1105 if ( fTrack1->fTpcV0NegPosSample && fTrack2->fNominalPosSample){
1106 while (fabs(fTrack1->fTpcV0NegPosSample[ipt].x())<9999. &&
1107 fabs(fTrack1->fTpcV0NegPosSample[ipt].y())<9999. &&
1108 fabs(fTrack1->fTpcV0NegPosSample[ipt].z())<9999. &&
1109 fabs(fTrack2->fNominalPosSample[ipt].x())<9999. &&
1110 fabs(fTrack2->fNominalPosSample[ipt].y())<9999. &&
1111 fabs(fTrack2->fNominalPosSample[ipt].z())<9999. &&
1112 (ipt<11)
1113 ){
1114 diff = fTrack1->fTpcV0NegPosSample[ipt] - fTrack2->fNominalPosSample[ipt];
1115 ipt++;
d0e92d9a 1116 tAveSep += diff.mag();
67427ff7 1117 }
d0e92d9a 1118 tAveSep = tAveSep/(ipt+1);
1119 return (tAveSep);}
67427ff7 1120 else return -1;
1121}
67427ff7 1122double AliFemtoPair::TpcExitSeparationV0NegV0Neg() const {
d0e92d9a 1123//_______1st part is a V0 considering Neg daughter 2nd is a V0 considering Neg daughter
67427ff7 1124 AliFemtoThreeVector diff = fTrack1->TpcV0NegExitPoint() - fTrack2->TpcV0NegExitPoint();
1125 return (diff.mag());
1126}
1127
1128double AliFemtoPair::TpcEntranceSeparationV0NegV0Neg() const {
d0e92d9a 1129//_______1st part is a V0 considering Neg daughter 2nd is a V0 considering Neg daughter
67427ff7 1130 AliFemtoThreeVector diff = fTrack1->TpcV0NegEntrancePoint() - fTrack2->TpcV0NegEntrancePoint();
1131 return (diff.mag());
1132}
1133double AliFemtoPair::TpcAverageSeparationV0NegV0Neg() const {
d0e92d9a 1134//_______1st part is a V0 considering Neg daughter 2nd is a V0 considering Neg daughter
67427ff7 1135 AliFemtoThreeVector diff;
d0e92d9a 1136 double tAveSep = 0.0;
67427ff7 1137 int ipt=0;
1138 if (fTrack1->fTpcV0NegPosSample && fTrack2->fTpcV0NegPosSample){
1139 while (fabs(fTrack1->fTpcV0NegPosSample[ipt].x())<9999. &&
1140 fabs(fTrack1->fTpcV0NegPosSample[ipt].y())<9999. &&
1141 fabs(fTrack1->fTpcV0NegPosSample[ipt].z())<9999. &&
1142 fabs(fTrack2->fTpcV0NegPosSample[ipt].x())<9999. &&
1143 fabs(fTrack2->fTpcV0NegPosSample[ipt].y())<9999. &&
1144 fabs(fTrack2->fTpcV0NegPosSample[ipt].z())<9999. &&
1145 (ipt<11)
1146 ){
1147 diff = fTrack1->fTpcV0NegPosSample[ipt] - fTrack2->fTpcV0NegPosSample[ipt];
1148 ipt++;
d0e92d9a 1149 tAveSep += diff.mag();
67427ff7 1150 }
d0e92d9a 1151 tAveSep = tAveSep/(ipt+1);
1152 return (tAveSep);}
67427ff7 1153 else return -1;
1154}
1155
67427ff7 1156void AliFemtoPair::CalcMergingParFctn(short* tmpMergingParNotCalculatedFctn,
1157 float* tmpZ1,float* tmpU1,
1158 float* tmpZ2,float* tmpU2,
1159 int *tmpSect1,int *tmpSect2,
1160 double* tmpFracOfMergedRow,
1161 double* tmpClosestRowAtDCA
1162 ) const{
d0e92d9a 1163// calculate heper variables for merging
67427ff7 1164 tmpMergingParNotCalculatedFctn=0;
1165 double tDu, tDz;
1166 int tN = 0;
1167 *tmpFracOfMergedRow = 0.;
1168 *tmpClosestRowAtDCA = 0.;
1169 double tDist;
1170 double tDistMax = 100000000.;
1171 for(int ti=0 ; ti<45 ; ti++){
1172 if(tmpSect1[ti]==tmpSect2[ti] && tmpSect1[ti]!=-1){
1173 tDu = fabs(tmpU1[ti]-tmpU2[ti]);
1174 tDz = fabs(tmpZ1[ti]-tmpZ2[ti]);
1175 tN++;
1176 if(ti<13){
d0e92d9a 1177 *tmpFracOfMergedRow += (tDu<fgMaxDuInner && tDz<fgMaxDzInner);
1178 tDist = ::sqrt(tDu*tDu/fgMaxDuInner/fgMaxDuInner+
1179 tDz*tDz/fgMaxDzInner/fgMaxDzInner);
67427ff7 1180 }
1181 else{
d0e92d9a 1182 *tmpFracOfMergedRow += (tDu<fgMaxDuOuter && tDz<fgMaxDzOuter);
1183 tDist = ::sqrt(tDu*tDu/fgMaxDuOuter/fgMaxDuOuter+
1184 tDz*tDz/fgMaxDzOuter/fgMaxDzOuter);
67427ff7 1185 }
1186 if(tDist<tDistMax){
1187 fClosestRowAtDCA = ti+1;
1188 tDistMax = tDist;
1189 }
1190 //fWeightedAvSep += tDist; // now, wrong but not used
1191 }
1192 }
1193 if(tN>0){
1194 //fWeightedAvSep /= tN;
1195 *tmpFracOfMergedRow /= tN;
1196 }
1197 else{
1198 *tmpClosestRowAtDCA = -1;
1199 *tmpFracOfMergedRow = -1.;
1200 //fWeightedAvSep = -1.;
1201 }
1202}
1203