]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWG2/FEMTOSCOPY/AliFemto/Infrastructure/AliFemtoPair.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Fixing Effective C++ warnings
[u/mrichter/AliRoot.git] / PWG2 / FEMTOSCOPY / AliFemto / Infrastructure / AliFemtoPair.cxx
CommitLineData
67427ff7 1/***************************************************************************
2 *
3 * $Id: AliFemtoPair.cc,v 1.23
4 *
5 * Author: Brian Laziuk, Yale University
6 * slightly modified by Mike Lisa
7 ***************************************************************************
8 *
9 * Description: part of STAR HBT Framework: AliFemtoMaker package
10 * the Pair object is passed to the PairCuts for verification, and
11 * then to the AddRealPair and AddMixedPair methods of the
12 * Correlation Functions
13 *
14 ***************************************************************************
15 * Revision 1.23 2002/09/25 19:23:25 rcwells
16 * Added const to emissionAngle()
17 *
18 * Revision 1.22 2002/04/22 22:48:11 laue
19 * corrected calculation of opening angle
20 **
21 * $Log$
0215f606 22 * Revision 1.1.1.1 2007/04/25 15:38:41 panos
23 * Importing the HBT code dir
24 *
67427ff7 25 * Revision 1.1.1.1 2007/03/07 10:14:49 mchojnacki
26 * First version on CVS
27 *
28 * Revision 1.27 2003/09/02 17:58:32 perev
29 * gcc 3.2 updates + WarnOff
30 *
31 * Revision 1.26 2003/01/31 19:57:15 magestro
32 * Cleared up simple compiler warnings on i386_linux24
33 *
34 * Revision 1.25 2003/01/14 09:44:08 renault
35 * corrections on average separation calculation for tracks which doesn't cross
36 * all 45 padrows.
37 *
38 * Revision 1.24 2002/11/19 23:33:10 renault
39 * Enable average separation calculation for all combinaisons of
40 * V0 daughters and tracks
41 *
42 * Revision 1.21 2002/02/28 14:18:36 rcwells
43 * Added emissionAngle function to AliFemtoPair
44 *
45 * Revision 1.20 2001/12/14 23:11:30 fretiere
46 * Add class HitMergingCut. Add class fabricesPairCut = HitMerginCut + pair purity cuts. Add TpcLocalTransform function which convert to local tpc coord (not pretty). Modify AliFemtoTrack, AliFemtoParticle, AliFemtoHiddenInfo, AliFemtoPair to handle the hit information and cope with my code
47 *
48 * Revision 1.19 2001/04/25 18:05:09 perev
49 * HPcorrs
50 *
51 * Revision 1.18 2001/04/03 21:04:36 kisiel
52 *
53 *
54 * Changes needed to make the Theoretical code
55 * work. The main code is the ThCorrFctn directory.
56 * The most visible change is the addition of the
57 * HiddenInfo to AliFemtoPair.
58 *
59 * Revision 1.17 2001/03/28 22:35:20 flierl
60 * changes and bugfixes in qYKP*
61 * add pairrapidity
62 *
63 * Revision 1.16 2001/02/15 19:23:00 rcwells
64 * Fixed sign in qSideCMS
65 *
66 * Revision 1.15 2001/01/22 22:56:41 laue
67 * Yano-Koonin-Podgoretskii Parametrisation added
68 *
69 * Revision 1.14 2000/12/11 21:44:30 rcwells
70 * Corrected qSideCMS function
71 *
72 * Revision 1.13 2000/10/26 16:09:16 lisa
73 * Added OpeningAngle PairCut class and method to AliFemtoPair
74 *
75 * Revision 1.12 2000/10/05 23:09:05 lisa
76 * Added kT-dependent radii to mixed-event simulator AND implemented AverageSeparation Cut and CorrFctn
77 *
78 * Revision 1.11 2000/07/17 20:03:16 lisa
79 * Implemented tools for addressing and assessing trackmerging
80 *
81 * Revision 1.10 2000/04/04 16:27:03 rcwells
82 * Removed an errant cout in AliFemtoPair.cc
83 *
84 * Revision 1.9 2000/04/04 16:13:09 lisa
85 * AliFemtoPair:quality() now returns normalized value (and so is double) and add a CorrFctn which looks at quality()
86 *
87 * Revision 1.8 2000/04/03 22:09:12 rcwells
88 * Add member function ... quality().
89 *
90 * Revision 1.7 2000/02/13 21:13:33 lisa
91 * changed ambiguous AliFemtoPair::fourMomentum() to fourMomentumSum() and fourMomentumDiff() and fixed related bug in QvecCorrFctn
92 *
93 * Revision 1.6 1999/07/29 16:16:34 lisa
94 * Selemons upgrade of AliFemtoPair class
95 *
96 * Revision 1.5 1999/07/22 18:49:10 lisa
97 * Implement idea of Fabrice to not create and delete AliFemtoPair all the time
98 *
99 * Revision 1.4 1999/07/12 18:57:05 lisa
100 * fixed small bug in fourMomentum method of AliFemtoPair
101 *
102 * Revision 1.3 1999/07/06 22:33:22 lisa
103 * Adjusted all to work in pro and new - dev itself is broken
104 *
105 * Revision 1.2 1999/06/29 17:50:27 fisyak
106 * formal changes to account new StEvent, does not complie yet
107 *
108 * Revision 1.1.1.1 1999/06/29 16:02:57 lisa
109 * Installation of AliFemtoMaker
110 *
111 **************************************************************************/
112
113#include "Infrastructure/AliFemtoPair.h"
114
115double AliFemtoPair::fMaxDuInner = .8;
116double AliFemtoPair::fMaxDzInner = 3.;
117double AliFemtoPair::fMaxDuOuter = 1.4;
118double AliFemtoPair::fMaxDzOuter = 3.2;
119
120
0215f606 121AliFemtoPair::AliFemtoPair() :
122 fTrack1(0), fTrack2(0),
123 fNonIdParNotCalculated(0),
124 fDKSide(0),
125 fDKOut(0),
126 fDKLong(0),
127 fCVK(0),
128 kStarCalc(0),
129 fNonIdParNotCalculatedGlobal(0),
130 fMergingParNotCalculated(0),
131 fWeightedAvSep(0),
132 fFracOfMergedRow(0),
133 fClosestRowAtDCA(0),
134 fMergingParNotCalculatedTrkV0Pos(0),
135 fFracOfMergedRowTrkV0Pos(0),
136 fClosestRowAtDCATrkV0Pos(0),
137 fMergingParNotCalculatedTrkV0Neg(0),
138 fFracOfMergedRowTrkV0Neg(0),
139 fClosestRowAtDCATrkV0Neg(0),
140 fMergingParNotCalculatedV0PosV0Neg(0),
141 fFracOfMergedRowV0PosV0Neg(0),
142 fClosestRowAtDCAV0PosV0Neg(0),
143 fMergingParNotCalculatedV0NegV0Pos(0),
144 fFracOfMergedRowV0NegV0Pos(0),
145 fClosestRowAtDCAV0NegV0Pos(0),
146 fMergingParNotCalculatedV0PosV0Pos(0),
147 fFracOfMergedRowV0PosV0Pos(0),
148 fClosestRowAtDCAV0PosV0Pos(0),
149 fMergingParNotCalculatedV0NegV0Neg(0),
150 fFracOfMergedRowV0NegV0Neg(0),
151 fClosestRowAtDCAV0NegV0Neg(0)
152{
67427ff7 153 fTrack1 = 0;
154 fTrack2 = 0;
155 setDefaultHalfFieldMergingPar();
156}
157
158AliFemtoPair::AliFemtoPair(AliFemtoParticle* a, AliFemtoParticle* b)
0215f606 159 : fTrack1(a), fTrack2(b),
160 fNonIdParNotCalculated(0),
161 fDKSide(0),
162 fDKOut(0),
163 fDKLong(0),
164 fCVK(0),
165 kStarCalc(0),
166 fNonIdParNotCalculatedGlobal(0),
167 fMergingParNotCalculated(0),
168 fWeightedAvSep(0),
169 fFracOfMergedRow(0),
170 fClosestRowAtDCA(0),
171 fMergingParNotCalculatedTrkV0Pos(0),
172 fFracOfMergedRowTrkV0Pos(0),
173 fClosestRowAtDCATrkV0Pos(0),
174 fMergingParNotCalculatedTrkV0Neg(0),
175 fFracOfMergedRowTrkV0Neg(0),
176 fClosestRowAtDCATrkV0Neg(0),
177 fMergingParNotCalculatedV0PosV0Neg(0),
178 fFracOfMergedRowV0PosV0Neg(0),
179 fClosestRowAtDCAV0PosV0Neg(0),
180 fMergingParNotCalculatedV0NegV0Pos(0),
181 fFracOfMergedRowV0NegV0Pos(0),
182 fClosestRowAtDCAV0NegV0Pos(0),
183 fMergingParNotCalculatedV0PosV0Pos(0),
184 fFracOfMergedRowV0PosV0Pos(0),
185 fClosestRowAtDCAV0PosV0Pos(0),
186 fMergingParNotCalculatedV0NegV0Neg(0),
187 fFracOfMergedRowV0NegV0Neg(0),
188 fClosestRowAtDCAV0NegV0Neg(0)
67427ff7 189{
190 setDefaultHalfFieldMergingPar();
191}
192
193void AliFemtoPair::setDefaultHalfFieldMergingPar(){
194 fMaxDuInner = 3;
195 fMaxDzInner = 4.;
196 fMaxDuOuter = 4.;
197 fMaxDzOuter = 6.;
198}
199void AliFemtoPair::setDefaultFullFieldMergingPar(){
200 fMaxDuInner = 0.8;
201 fMaxDzInner = 3.;
202 fMaxDuOuter = 1.4;
203 fMaxDzOuter = 3.2;
204}
205void AliFemtoPair::setMergingPar(double aMaxDuInner, double aMaxDzInner,
206 double aMaxDuOuter, double aMaxDzOuter){
207 fMaxDuInner = aMaxDuInner;
208 fMaxDzInner = aMaxDzInner;
209 fMaxDuOuter = aMaxDuOuter;
210 fMaxDzOuter = aMaxDzOuter;
211};
212
213AliFemtoPair::~AliFemtoPair() {/* no-op */}
214
0215f606 215AliFemtoPair::AliFemtoPair(const AliFemtoPair &aPair):
216 fTrack1(0), fTrack2(0),
217 fNonIdParNotCalculated(0),
218 fDKSide(0),
219 fDKOut(0),
220 fDKLong(0),
221 fCVK(0),
222 kStarCalc(0),
223 fNonIdParNotCalculatedGlobal(0),
224 fMergingParNotCalculated(0),
225 fWeightedAvSep(0),
226 fFracOfMergedRow(0),
227 fClosestRowAtDCA(0),
228 fMergingParNotCalculatedTrkV0Pos(0),
229 fFracOfMergedRowTrkV0Pos(0),
230 fClosestRowAtDCATrkV0Pos(0),
231 fMergingParNotCalculatedTrkV0Neg(0),
232 fFracOfMergedRowTrkV0Neg(0),
233 fClosestRowAtDCATrkV0Neg(0),
234 fMergingParNotCalculatedV0PosV0Neg(0),
235 fFracOfMergedRowV0PosV0Neg(0),
236 fClosestRowAtDCAV0PosV0Neg(0),
237 fMergingParNotCalculatedV0NegV0Pos(0),
238 fFracOfMergedRowV0NegV0Pos(0),
239 fClosestRowAtDCAV0NegV0Pos(0),
240 fMergingParNotCalculatedV0PosV0Pos(0),
241 fFracOfMergedRowV0PosV0Pos(0),
242 fClosestRowAtDCAV0PosV0Pos(0),
243 fMergingParNotCalculatedV0NegV0Neg(0),
244 fFracOfMergedRowV0NegV0Neg(0),
245 fClosestRowAtDCAV0NegV0Neg(0)
246{
247 fTrack1 = aPair.fTrack1;
248 fTrack2 = aPair.fTrack2;
249
250 fNonIdParNotCalculated = aPair.fNonIdParNotCalculated;
251 fDKSide = aPair.fDKSide;
252 fDKOut = aPair.fDKOut;
253 fDKLong = aPair.fDKLong;
254 fCVK = aPair.fCVK;
255 kStarCalc = aPair.kStarCalc;
256
257 fNonIdParNotCalculatedGlobal = aPair.fNonIdParNotCalculatedGlobal;
258
259 fMergingParNotCalculated = aPair.fMergingParNotCalculated;
260 fWeightedAvSep = aPair.fWeightedAvSep;
261 fFracOfMergedRow = aPair.fFracOfMergedRow;
262 fClosestRowAtDCA = aPair.fClosestRowAtDCA;
263
264 fMergingParNotCalculatedTrkV0Pos = aPair.fMergingParNotCalculatedTrkV0Pos;
265 fFracOfMergedRowTrkV0Pos = aPair.fFracOfMergedRowTrkV0Pos;
266 fClosestRowAtDCATrkV0Pos = aPair.fClosestRowAtDCATrkV0Pos;
267
268 fMergingParNotCalculatedTrkV0Neg = aPair.fMergingParNotCalculatedTrkV0Neg;
269 fFracOfMergedRowTrkV0Neg = aPair.fFracOfMergedRowTrkV0Neg;
270 fClosestRowAtDCATrkV0Neg = aPair.fClosestRowAtDCATrkV0Neg;
271
272 fMergingParNotCalculatedV0PosV0Neg = aPair.fMergingParNotCalculatedV0PosV0Neg;
273 fFracOfMergedRowV0PosV0Neg = aPair.fFracOfMergedRowV0PosV0Neg;
274 fClosestRowAtDCAV0PosV0Neg = aPair.fClosestRowAtDCAV0PosV0Neg;
275
276 fMergingParNotCalculatedV0NegV0Pos = aPair.fMergingParNotCalculatedV0NegV0Pos;
277 fFracOfMergedRowV0NegV0Pos = aPair.fFracOfMergedRowV0NegV0Pos;
278 fClosestRowAtDCAV0NegV0Pos = aPair.fClosestRowAtDCAV0NegV0Pos;
279
280 fMergingParNotCalculatedV0PosV0Pos = aPair.fMergingParNotCalculatedV0PosV0Pos;
281 fFracOfMergedRowV0PosV0Pos = aPair.fFracOfMergedRowV0PosV0Pos;
282 fClosestRowAtDCAV0PosV0Pos = aPair.fClosestRowAtDCAV0PosV0Pos;
283
284 fMergingParNotCalculatedV0NegV0Neg = aPair.fMergingParNotCalculatedV0NegV0Neg;
285 fFracOfMergedRowV0NegV0Neg = aPair.fFracOfMergedRowV0NegV0Neg;
286 fClosestRowAtDCAV0NegV0Neg = aPair.fClosestRowAtDCAV0NegV0Neg;
287}
288
289AliFemtoPair& AliFemtoPair::operator=(const AliFemtoPair &aPair)
290{
291 if (this == &aPair)
292 return *this;
67427ff7 293
0215f606 294 fTrack1 = aPair.fTrack1;
295 fTrack2 = aPair.fTrack2;
296
297 fNonIdParNotCalculated = aPair.fNonIdParNotCalculated;
298 fDKSide = aPair.fDKSide;
299 fDKOut = aPair.fDKOut;
300 fDKLong = aPair.fDKLong;
301 fCVK = aPair.fCVK;
302 kStarCalc = aPair.kStarCalc;
303
304 fNonIdParNotCalculatedGlobal = aPair.fNonIdParNotCalculatedGlobal;
305
306 fMergingParNotCalculated = aPair.fMergingParNotCalculated;
307 fWeightedAvSep = aPair.fWeightedAvSep;
308 fFracOfMergedRow = aPair.fFracOfMergedRow;
309 fClosestRowAtDCA = aPair.fClosestRowAtDCA;
310
311 fMergingParNotCalculatedTrkV0Pos = aPair.fMergingParNotCalculatedTrkV0Pos;
312 fFracOfMergedRowTrkV0Pos = aPair.fFracOfMergedRowTrkV0Pos;
313 fClosestRowAtDCATrkV0Pos = aPair.fClosestRowAtDCATrkV0Pos;
314
315 fMergingParNotCalculatedTrkV0Neg = aPair.fMergingParNotCalculatedTrkV0Neg;
316 fFracOfMergedRowTrkV0Neg = aPair.fFracOfMergedRowTrkV0Neg;
317 fClosestRowAtDCATrkV0Neg = aPair.fClosestRowAtDCATrkV0Neg;
318
319 fMergingParNotCalculatedV0PosV0Neg = aPair.fMergingParNotCalculatedV0PosV0Neg;
320 fFracOfMergedRowV0PosV0Neg = aPair.fFracOfMergedRowV0PosV0Neg;
321 fClosestRowAtDCAV0PosV0Neg = aPair.fClosestRowAtDCAV0PosV0Neg;
322
323 fMergingParNotCalculatedV0NegV0Pos = aPair.fMergingParNotCalculatedV0NegV0Pos;
324 fFracOfMergedRowV0NegV0Pos = aPair.fFracOfMergedRowV0NegV0Pos;
325 fClosestRowAtDCAV0NegV0Pos = aPair.fClosestRowAtDCAV0NegV0Pos;
326
327 fMergingParNotCalculatedV0PosV0Pos = aPair.fMergingParNotCalculatedV0PosV0Pos;
328 fFracOfMergedRowV0PosV0Pos = aPair.fFracOfMergedRowV0PosV0Pos;
329 fClosestRowAtDCAV0PosV0Pos = aPair.fClosestRowAtDCAV0PosV0Pos;
330
331 fMergingParNotCalculatedV0NegV0Neg = aPair.fMergingParNotCalculatedV0NegV0Neg;
332 fFracOfMergedRowV0NegV0Neg = aPair.fFracOfMergedRowV0NegV0Neg;
333 fClosestRowAtDCAV0NegV0Neg = aPair.fClosestRowAtDCAV0NegV0Neg;
334
335 return *this;
336}
67427ff7 337
338//_________________
339double AliFemtoPair::mInv() const
340{
341 double InvariantMass = abs(fTrack1->FourMomentum() + fTrack2->FourMomentum());
342 return (InvariantMass);
343}
344//_________________
345double AliFemtoPair::kT() const
346{
347
348 double tmp =
349 (fTrack1->FourMomentum() + fTrack2->FourMomentum()).perp();
350 tmp *= .5;
351
352 return (tmp);
353}
354//_________________
355double AliFemtoPair::rap() const
356{
357 // longitudinal pair rapidity : Y = 0.5 ::log( E1 + E2 + pz1 + pz2 / E1 + E2 - pz1 - pz2 )
358 double tmp = 0.5 * log (
359 (fTrack1->FourMomentum().e() + fTrack2->FourMomentum().e() + fTrack1->FourMomentum().z() + fTrack2->FourMomentum().z()) /
360 (fTrack1->FourMomentum().e() + fTrack2->FourMomentum().e() - fTrack1->FourMomentum().z() - fTrack2->FourMomentum().z())
361 ) ;
362 return (tmp);
363}
364//_________________
365double AliFemtoPair::emissionAngle()const {
366 double pxTotal = this->fourMomentumSum().x();
367 double pyTotal = this->fourMomentumSum().y();
368 double angle = atan2(pyTotal,pxTotal)*180.0/3.1415926536;
369 if (angle<0.0) angle+=360.0;
370 return angle;
371}
372//_________________
373// get rid of ambiguously-named method fourMomentum() and replace it with
374// fourMomentumSum() and fourMomentumDiff() - mal 13feb2000
375AliFemtoLorentzVector AliFemtoPair::fourMomentumSum() const
376{
377 AliFemtoLorentzVector temp = fTrack1->FourMomentum()+fTrack2->FourMomentum();
378 return temp;
379}
380AliFemtoLorentzVector AliFemtoPair::fourMomentumDiff() const
381{
382 AliFemtoLorentzVector temp = fTrack1->FourMomentum()-fTrack2->FourMomentum();
383 return temp;
384}
385//__________________________________
386// Yano-Koonin-Podgoretskii Parametrisation in CMS
387void AliFemtoPair::qYKPCMS(double& qP, double& qT, double& q0) const
388{
389 ////
390 // calculate momentum difference in source rest frame (= lab frame)
391 ////
392 AliFemtoLorentzVector l1 = fTrack1->FourMomentum() ;
393 AliFemtoLorentzVector l2 = fTrack2->FourMomentum() ;
394 AliFemtoLorentzVector l ;
395 // random ordering of the particles
396 if ( rand()/(double)RAND_MAX > 0.50 )
397 { l = l1-l2 ; }
398 else
399 { l = l2-l1 ; } ;
400 // fill momentum differences into return variables
401 qP = l.z() ;
402 qT = l.vect().perp() ;
403 q0 = l.e() ;
404}
405//___________________________________
406// Yano-Koonin-Podgoretskii Parametrisation in LCMS
407void AliFemtoPair::qYKPLCMS(double& qP, double& qT, double& q0) const
408{
409 ////
410 // calculate momentum difference in LCMS : frame where pz1 + pz2 = 0
411 ////
412 AliFemtoLorentzVector l1 = fTrack1->FourMomentum() ;
413 AliFemtoLorentzVector l2 = fTrack2->FourMomentum() ;
414 // determine beta to LCMS
415 double beta = (l1.z()+l2.z()) / (l1.e()+l2.e()) ;
416 double beta2 = beta*beta ;
417 // unfortunately STAR Class lib knows only boost(particle) not boost(beta) :(
418 // -> create particle with velocity beta and mass 1.0
419 // actually this is : dummyPz = ::sqrt( (dummyMass*dummyMass*beta2) / (1-beta2) ) ;
420 double dummyPz = ::sqrt( (beta2) / (1-beta2) ) ;
421 // boost in the correct direction
422 if (beta>0.0) { dummyPz = -dummyPz; } ;
423 // create dummy particle
424 AliFemtoLorentzVector l(0.0, 0.0, dummyPz) ;
425 double dummyMass = 1.0 ;
426 l.setE(l.vect().massHypothesis(dummyMass) );
427 // boost particles along the beam into a frame with velocity beta
428 AliFemtoLorentzVector l1boosted = l1.boost(l) ;
429 AliFemtoLorentzVector l2boosted = l2.boost(l) ;
430 // caculate the momentum difference with random ordering of the particle
431 if ( rand()/(double)RAND_MAX >0.50)
432 { l = l1boosted-l2boosted ; }
433 else
434 { l = l2boosted-l1boosted ;} ;
435 // fill momentum differences into return variables
436 qP = l.z() ;
437 qT = l.vect().perp() ;
438 q0 = l.e() ;
439}
440//___________________________________
441// Yano-Koonin-Podgoretskii Parametrisation in pair rest frame
442void AliFemtoPair::qYKPPF(double& qP, double& qT, double& q0) const
443{
444 ////
445 // calculate momentum difference in pair rest frame : frame where (pz1 + pz2, py1 + py2, px1 + px2) = (0,0,0)
446 ////
447 AliFemtoLorentzVector l1 = fTrack1->FourMomentum() ;
448 AliFemtoLorentzVector l2 = fTrack2->FourMomentum() ;
449 // the center of gravity of the pair travels with l
450 AliFemtoLorentzVector l = l1 + l2 ;
451 l = -l ;
452 l.setE(-l.e()) ;
453 // boost particles
454 AliFemtoLorentzVector l1boosted = l1.boost(l) ;
455 AliFemtoLorentzVector l2boosted = l2.boost(l) ;
456 // caculate the momentum difference with random ordering of the particle
457 if ( rand()/(double)RAND_MAX > 0.50)
458 { l = l1boosted-l2boosted ; }
459 else
460 { l = l2boosted-l1boosted ;} ;
461 // fill momentum differences into return variables
462 qP = l.z();
463 qT = l.vect().perp();
464 q0 = l.e();
465}
466//_________________
467double AliFemtoPair::qOutCMS() const
468{
469 AliFemtoThreeVector tmp1 = fTrack1->FourMomentum().vect();
470 AliFemtoThreeVector tmp2 = fTrack2->FourMomentum().vect();
471
472 double dx = tmp1.x() - tmp2.x();
473 double xt = tmp1.x() + tmp2.x();
474
475 double dy = tmp1.y() - tmp2.y();
476 double yt = tmp1.y() + tmp2.y();
477
478 double k1 = (::sqrt(xt*xt+yt*yt));
479 double k2 = (dx*xt+dy*yt);
480 double tmp = k2/k1;
481 return (tmp);
482}
483//_________________
484double AliFemtoPair::qSideCMS() const
485{
486 AliFemtoThreeVector tmp1 = fTrack1->FourMomentum().vect();
487 AliFemtoThreeVector tmp2 = fTrack2->FourMomentum().vect();
488
489 double x1 = tmp1.x(); double y1 = tmp1.y();
490 double x2 = tmp2.x(); double y2 = tmp2.y();
491
492 double xt = x1+x2; double yt = y1+y2;
493 double k1 = ::sqrt(xt*xt+yt*yt);
494
495 double tmp = 2.0*(x2*y1-x1*y2)/k1;
496 return (tmp);
497}
498
499//_________________________
500double AliFemtoPair::qLongCMS() const
501{
502 AliFemtoLorentzVector tmp1 = fTrack1->FourMomentum();
503 AliFemtoLorentzVector tmp2 = fTrack2->FourMomentum();
504
505 double dz = tmp1.z() - tmp2.z();
506 double zz = tmp1.z() + tmp2.z();
507
508 double dt = tmp1.t() - tmp2.t();
509 double tt = tmp1.t() + tmp2.t();
510
511 double beta = zz/tt;
512 double gamma = 1.0/::sqrt(1.0 - beta*beta);
513
514 double temp = gamma*(dz - beta*dt);
515 return (temp);
516}
517
518//________________________________
519double AliFemtoPair::qOutPf() const
520{
521 AliFemtoLorentzVector tmp1 = fTrack1->FourMomentum();
522 AliFemtoLorentzVector tmp2 = fTrack2->FourMomentum();
523
524 double dt = tmp1.t() - tmp2.t();
525 double tt = tmp1.t() + tmp2.t();
526
527 double xt = tmp1.x() + tmp2.x();
528 double yt = tmp1.y() + tmp2.y();
529
530 double k1 = ::sqrt(xt*xt + yt*yt);
531 double bOut = k1/tt;
532 double gOut = 1.0/::sqrt(1.0 - bOut*bOut);
533
534 double temp = gOut*(this->qOutCMS() - bOut*dt);
535 return (temp);
536}
537
538//___________________________________
539double AliFemtoPair::qSidePf() const
540{
541 return(this->qSideCMS());
542}
543
544//___________________________________
545
546double AliFemtoPair::qLongPf() const
547{
548 return(this->qLongCMS());
549}
550
551//___________________________________
552double AliFemtoPair::qOutBf(double beta) const
553{
554 return(this->qOutCMS());
555}
556
557//___________________________________
558
559double AliFemtoPair::qSideBf(double beta) const
560{
561 return(this->qSideCMS());
562}
563
564//___________________________________
565double AliFemtoPair::qLongBf(double beta) const
566{
567 AliFemtoLorentzVector tmp1 = fTrack1->FourMomentum();
568 AliFemtoLorentzVector tmp2 = fTrack2->FourMomentum();
569
570 double dz = tmp1.z() - tmp2.z();
571 double dt = tmp1.t() + tmp2.t();
572
573 double gamma = 1.0/::sqrt(1.0 - beta*beta);
574
575 double temp = gamma*(dz - beta*dt);
576 return (temp);
577}
578
579double AliFemtoPair::quality() const {
580 unsigned long mapMask0 = 0xFFFFFF00;
581 unsigned long mapMask1 = 0x1FFFFF;
582 unsigned long padRow1To24Track1 = fTrack1->TopologyMap(0) & mapMask0;
583 unsigned long padRow25To45Track1 = fTrack1->TopologyMap(1) & mapMask1;
584 unsigned long padRow1To24Track2 = fTrack2->TopologyMap(0) & mapMask0;
585 unsigned long padRow25To45Track2 = fTrack2->TopologyMap(1) & mapMask1;
586 // AND logic
587 unsigned long bothPads1To24 = padRow1To24Track1 & padRow1To24Track2;
588 unsigned long bothPads25To45 = padRow25To45Track1 & padRow25To45Track2;
589 // XOR logic
590 unsigned long onePad1To24 = padRow1To24Track1 ^ padRow1To24Track2;
591 unsigned long onePad25To45 = padRow25To45Track1 ^ padRow25To45Track2;
592 unsigned long bitI;
593 int ibits;
594 int Quality = 0;
595 double normQual = 0.0;
596 int MaxQuality = fTrack1->NumberOfHits() + fTrack2->NumberOfHits();
597 for (ibits=8;ibits<=31;ibits++) {
598 bitI = 0;
599 bitI |= 1UL<<(ibits);
600 if ( onePad1To24 & bitI ) {
601 Quality++;
602 continue;
603 }
604 else{
605 if ( bothPads1To24 & bitI ) Quality--;
606 }
607 }
608 for (ibits=0;ibits<=20;ibits++) {
609 bitI = 0;
610 bitI |= 1UL<<(ibits);
611 if ( onePad25To45 & bitI ) {
612 Quality++;
613 continue;
614 }
615 else{
616 if ( bothPads25To45 & bitI ) Quality--;
617 }
618 }
619 normQual = (double)Quality/( (double) MaxQuality );
620 return ( normQual );
621
622}
623
624double AliFemtoPair::quality2() const {
625 unsigned long mapMask0 = 0xFFFFFF00;
626 unsigned long mapMask1 = 0x1FFFFF;
627 unsigned long padRow1To24Track1 = fTrack1->TopologyMap(0) & mapMask0;
628 unsigned long padRow25To45Track1 = fTrack1->TopologyMap(1) & mapMask1;
629 unsigned long padRow1To24Track2 = fTrack2->TopologyMap(0) & mapMask0;
630 unsigned long padRow25To45Track2 = fTrack2->TopologyMap(1) & mapMask1;
631
632 // AND logic
633 //unsigned long bothPads1To24 = padRow1To24Track1 & padRow1To24Track2;
634 //unsigned long bothPads25To45 = padRow25To45Track1 & padRow25To45Track2;
635
636 // XOR logic
637 unsigned long onePad1To24 = padRow1To24Track1 ^ padRow1To24Track2;
638 unsigned long onePad25To45 = padRow25To45Track1 ^ padRow25To45Track2;
639 unsigned long bitI;
640 int ibits;
641 int Quality = 0;
642 double normQual = 0.0;
643 int MaxQuality = fTrack1->NumberOfHits() + fTrack2->NumberOfHits();
644 for (ibits=8;ibits<=31;ibits++) {
645 bitI = 0;
646 bitI |= 1UL<<(ibits);
647 if ( onePad1To24 & bitI ) {
648 Quality++;
649 continue;
650 }
651 //else{
652 //if ( bothPads1To24 & bitI ) Quality--;
653 //}
654 }
655 for (ibits=0;ibits<=20;ibits++) {
656 bitI = 0;
657 bitI |= 1UL<<(ibits);
658 if ( onePad25To45 & bitI ) {
659 Quality++;
660 continue;
661 }
662 //else{
663 //if ( bothPads25To45 & bitI ) Quality--;
664 //}
665 }
666 normQual = (double)Quality/( (double) MaxQuality );
667 return ( normQual );
668
669}
670
671
672double AliFemtoPair::NominalTpcExitSeparation() const {
673 AliFemtoThreeVector diff = fTrack1->NominalTpcExitPoint() - fTrack2->NominalTpcExitPoint();
674 return (diff.mag());
675}
676
677double AliFemtoPair::NominalTpcEntranceSeparation() const {
678 AliFemtoThreeVector diff = fTrack1->NominalTpcEntrancePoint() - fTrack2->NominalTpcEntrancePoint();
679 return (diff.mag());
680}
681
682double AliFemtoPair::NominalTpcAverageSeparation() const {
683 AliFemtoThreeVector diff;
684 double AveSep = 0.0;
685 int ipt = 0;
686 if (fTrack1->fNominalPosSample && fTrack2->fNominalPosSample){
687 while (fabs(fTrack1->fNominalPosSample[ipt].x())<9999. &&
688 fabs(fTrack1->fNominalPosSample[ipt].y())<9999. &&
689 fabs(fTrack1->fNominalPosSample[ipt].z())<9999. &&
690 fabs(fTrack2->fNominalPosSample[ipt].x())<9999. &&
691 fabs(fTrack2->fNominalPosSample[ipt].y())<9999. &&
692 fabs(fTrack2->fNominalPosSample[ipt].z())<9999. &&
693 ipt<11
694 ){
695 // for (int ipt=0; ipt<11; ipt++){
696 diff = fTrack1->fNominalPosSample[ipt] - fTrack2->fNominalPosSample[ipt];
697 ipt++;
698 AveSep += diff.mag();
699 }
700 AveSep = AveSep/(ipt+1.);
701 return (AveSep);}
702 else return -1;
703}
704
705double AliFemtoPair::OpeningAngle() const {
706 return 57.296* fTrack1->FourMomentum().vect().angle( fTrack2->FourMomentum().vect() );
707// AliFemtoThreeVector p1 = fTrack1->FourMomentum().vect();
708// AliFemtoThreeVector p2 = fTrack2->FourMomentum().vect();
709// return 57.296*(p1.phi()-p2.phi());
710// //double dAngInv = 57.296*acos((p1.dot(p2))/(p1.mag()*p2.mag()));
711// //return (dAngInv);
712}
713//_________________
714
715
716double AliFemtoPair::KStarFlipped() const {
717 AliFemtoLorentzVector tP1 = fTrack1->FourMomentum();
718
719 AliFmThreeVectorD qwe = tP1.vect();
720 qwe *= -1.; // flip it
721 tP1.setVect(qwe);
722
723 AliFemtoLorentzVector tSum = (tP1+fTrack2->FourMomentum());
724 double tMass = abs(tSum);
725 AliFmThreeVectorD tGammaBeta = (1./tMass)*tSum.vect();
726 double tGamma = tSum.e()/tMass;
727 AliFmThreeVectorD tLongMom = ((tP1.vect()*tGammaBeta)/
728 (tGammaBeta*tGammaBeta))*tGammaBeta;
729 AliFmLorentzVectorD tK(tGamma*tP1.e() - tP1.vect()*tGammaBeta,
730 tP1.vect() + (tGamma-1.)*tLongMom - tP1.e()*tGammaBeta);
731//VP tP1.vect() *= -1.; // unflip it
732 return tK.vect().mag();
733}
734
735//double AliFemtoPair::CVK() const{
736//const AliFemtoLorentzVector& tP1 = fTrack1->FourMomentum();
737//AliFemtoLorentzVector tSum = (tP1+fTrack2->FourMomentum());
738//double tMass = abs(tSum);
739//AliFmThreeVectorD tGammaBeta = (1./tMass)*tSum.vect();
740//double tGamma = tSum.e()/tMass;
741//AliFmThreeVectorD tLongMom = ((tP1.vect()*tGammaBeta)/
742// (tGammaBeta*tGammaBeta))*tGammaBeta;
743//AliFmLorentzVectorD tK(tGamma*tP1.e() - tP1.vect()*tGammaBeta,
744// tP1.vect() + (tGamma-1.)*tLongMom - tP1.e()*tGammaBeta);
745//return (tK.vect())*tGammaBeta/tK.vect().magnitude()/tGammaBeta.magnitude();
746//}
747
748double AliFemtoPair::CVKFlipped() const{
749 AliFemtoLorentzVector tP1 = fTrack1->FourMomentum();
750 AliFmThreeVectorD qwe = tP1.vect();
751 qwe *= -1.; // flip it
752 tP1.setVect(qwe);
753
754 AliFemtoLorentzVector tSum = (tP1+fTrack2->FourMomentum());
755 double tMass = abs(tSum);
756 AliFmThreeVectorD tGammaBeta = (1./tMass)*tSum.vect();
757 double tGamma = tSum.e()/tMass;
758 AliFmThreeVectorD tLongMom = ((tP1.vect()*tGammaBeta)/
759 (tGammaBeta*tGammaBeta))*tGammaBeta;
760 AliFmLorentzVectorD tK(tGamma*tP1.e() - tP1.vect()*tGammaBeta,
761 tP1.vect() + (tGamma-1.)*tLongMom - tP1.e()*tGammaBeta);
762//VP tP1.vect() *= -1.; // unflip it
763 return (tK.vect())*tGammaBeta/tGamma;
764}
765
766double AliFemtoPair::pInv() const{
767 AliFemtoLorentzVector tP1 = fTrack1->FourMomentum();
768 AliFemtoLorentzVector tP2 = fTrack2->FourMomentum();
769 double tP = (tP1.px()+tP2.px())*(tP1.px()+tP2.px())+
770 (tP1.py()+tP2.py())*(tP1.py()+tP2.py())+
771 (tP1.pz()+tP2.pz())*(tP1.pz()+tP2.pz())-
772 (tP1.e() -tP2.e() )*(tP1.e() -tP2.e() );
773 return ::sqrt(fabs(tP));
774}
775
776double AliFemtoPair::qInvFlippedXY() const{
777 AliFemtoLorentzVector tP1 = fTrack1->FourMomentum();
778 tP1.setX(-1.*tP1.x());
779 tP1.setY(-1.*tP1.y());
780 AliFemtoLorentzVector tDiff = (tP1-fTrack2->FourMomentum());
781 return ( -1.* tDiff.m());
782}
783
784void AliFemtoPair::calcNonIdPar() const{ // fortran like function! faster?
785 fNonIdParNotCalculated=0;
786 double px1 = fTrack1->FourMomentum().vect().x();
787 double py1 = fTrack1->FourMomentum().vect().y();
788 double pz1 = fTrack1->FourMomentum().vect().z();
789 double pE1 = fTrack1->FourMomentum().e();
790 double Particle1Mass = ::sqrt(pE1*pE1 - px1*px1 - py1*py1 - pz1*pz1);
791 double px2 = fTrack2->FourMomentum().vect().x();
792 double py2 = fTrack2->FourMomentum().vect().y();
793 double pz2 = fTrack2->FourMomentum().vect().z();
794 double pE2 = fTrack2->FourMomentum().e();
795 double Particle2Mass = ::sqrt(pE2*pE2 - px2*px2 - py2*py2 - pz2*pz2);
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
811 double Q = (Particle1Mass*Particle1Mass - Particle2Mass*Particle2Mass)/Pinv;
812 Q = sqrt ( Q*Q - QinvL);
813
814 kStarCalc = 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 fDKLong = 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 fDKSide = 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 fDKOut = px1C;
841
842 fCVK = (fDKOut*Ptrans + fDKLong*Pz)/kStarCalc/::sqrt(Ptrans*Ptrans+Pz*Pz);
843}
844
845
846/*void AliFemtoPair::calcNonIdParGlobal() const{ // fortran like function! faster?
847 fNonIdParNotCalculatedGlobal=0;
848 double px1 = fTrack1->Track()->PGlobal().x();
849 double py1 = fTrack1->Track()->PGlobal().y();
850 double pz1 = fTrack1->Track()->PGlobal().z();
851 double Particle1Mass = fTrack1->FourMomentum().m2();
852 double pE1 = ::sqrt(Particle1Mass + px1*px1 + py1*py1 + pz1*pz1);
853 Particle1Mass = ::sqrt(Particle1Mass);
854
855 double px2 = fTrack2->Track()->PGlobal().x();
856 double py2 = fTrack2->Track()->PGlobal().y();
857 double pz2 = fTrack2->Track()->PGlobal().z();
858 double Particle2Mass = fTrack2->FourMomentum().m2();
859 double pE2 = ::sqrt(Particle2Mass + px2*px2 + py2*py2 + pz2*pz2);
860 Particle2Mass = ::sqrt(Particle2Mass);
861
862 double Px = px1+px2;
863 double Py = py1+py2;
864 double Pz = pz1+pz2;
865 double PE = pE1+pE2;
866
867 double Ptrans = Px*Px + Py*Py;
868 double Mtrans = PE*PE - Pz*Pz;
869 double Pinv = ::sqrt(Mtrans - Ptrans);
870 Mtrans = ::sqrt(Mtrans);
871 Ptrans = ::sqrt(Ptrans);
872
873 double QinvL = (pE1-pE2)*(pE1-pE2) - (px1-px2)*(px1-px2) -
874 (py1-py2)*(py1-py2) - (pz1-pz2)*(pz1-pz2);
875
876 double Q = (Particle1Mass*Particle1Mass - Particle2Mass*Particle2Mass)/Pinv;
877 Q = sqrt ( Q*Q - QinvL);
878
879 kStarCalcGlobal = Q/2;
880
881 // ad 1) go to LCMS
882 double beta = Pz/PE;
883 double gamma = PE/Mtrans;
884
885 double pz1L = gamma * (pz1 - beta * pE1);
886 double pE1L = gamma * (pE1 - beta * pz1);
887
888 // fill histogram for beam projection ( z - axis )
889 fDKLongGlobal = pz1L;
890
891 // ad 2) rotation px -> Pt
892 double px1R = (px1*Px + py1*Py)/Ptrans;
893 double py1R = (-px1*Py + py1*Px)/Ptrans;
894
895 //fill histograms for side projection ( y - axis )
896 fDKSideGlobal = py1R;
897
898 // ad 3) go from LCMS to CMS
899 beta = Ptrans/Mtrans;
900 gamma = Mtrans/Pinv;
901
902 double px1C = gamma * (px1R - beta * pE1L);
903
904 // fill histogram for out projection ( x - axis )
905 fDKOutGlobal = px1C;
906
907 fCVKGlobal = (fDKOutGlobal*Ptrans + fDKLongGlobal*Pz)/
908 kStarCalcGlobal/::sqrt(Ptrans*Ptrans+Pz*Pz);
909}*/
910
911
912
913double AliFemtoPair::dcaInsideTpc() const{
914
915 double tMinDist=NominalTpcEntranceSeparation();
916 double tExit = NominalTpcExitSeparation();
917 tMinDist = (tExit>tMinDist) ? tMinDist : tExit;
918 double tInsideDist;
919 //tMinDist = 999.;
920
921 double rMin = 60.;
922 double rMax = 190.;
923 const AliFmPhysicalHelixD& tHelix1 = fTrack1->Helix();
924 const AliFmPhysicalHelixD& tHelix2 = fTrack2->Helix();
925 // --- One is a line and other one a helix
926 //if (tHelix1.mSingularity != tHelix2.mSingularity) return -999.;
927 // --- 2 lines : don't care right now
928 //if (tHelix1.mSingularity) return -999.;
929 // --- 2 helix
930 double dx = tHelix2.xcenter() - tHelix1.xcenter();
931 double dy = tHelix2.ycenter() - tHelix1.ycenter();
932 double dd = ::sqrt(dx*dx + dy*dy);
933 double r1 = 1/tHelix1.curvature();
934 double r2 = 1/tHelix2.curvature();
935 double cosAlpha = (r1*r1 + dd*dd - r2*r2)/(2*r1*dd);
936
937 double x, y, r;
938 double s;
939 if (fabs(cosAlpha) < 1) { // two solutions
940 double sinAlpha = sin(acos(cosAlpha));
941 x = tHelix1.xcenter() + r1*(cosAlpha*dx - sinAlpha*dy)/dd;
942 y = tHelix1.ycenter() + r1*(sinAlpha*dx + cosAlpha*dy)/dd;
943 r = ::sqrt(x*x+y*y);
944 if( r > rMin && r < rMax &&
945 fabs(atan2(y,x)-fTrack1->NominalTpcEntrancePoint().phi())< 0.5
946 ){ // first solution inside
947 s = tHelix1.pathLength(x, y);
948 tInsideDist=tHelix2.distance(tHelix1.at(s));
949 if(tInsideDist<tMinDist) tMinDist = tInsideDist;
950 }
951 else{
952 x = tHelix1.xcenter() + r1*(cosAlpha*dx + sinAlpha*dy)/dd;
953 y = tHelix1.ycenter() + r1*(cosAlpha*dy - sinAlpha*dx)/dd;
954 r = ::sqrt(x*x+y*y);
955 if( r > rMin && r < rMax &&
956 fabs(atan2(y,x)-fTrack1->NominalTpcEntrancePoint().phi())< 0.5
957 ) { // second solution inside
958 s = tHelix1.pathLength(x, y);
959 tInsideDist=tHelix2.distance(tHelix1.at(s));
960 if(tInsideDist<tMinDist) tMinDist = tInsideDist;
961 }
962 }
963 }
964 return tMinDist;
965}
966
967void AliFemtoPair::calcMergingPar() const{
968 fMergingParNotCalculated=0;
969
970 double tDu, tDz;
971 int tN = 0;
972 fFracOfMergedRow = 0.;
973 fWeightedAvSep =0.;
974 double tDist;
975 double tDistMax = 200.;
976 for(int ti=0 ; ti<45 ; ti++){
977 if(fTrack1->fSect[ti]==fTrack2->fSect[ti] && fTrack1->fSect[ti]!=-1){
978 tDu = fabs(fTrack1->fU[ti]-fTrack2->fU[ti]);
979 tDz = fabs(fTrack1->fZ[ti]-fTrack2->fZ[ti]);
980 tN++;
981 if(ti<13){
982 fFracOfMergedRow += (tDu<fMaxDuInner && tDz<fMaxDzInner);
983 tDist = ::sqrt(tDu*tDu/fMaxDuInner/fMaxDuInner+
984 tDz*tDz/fMaxDzInner/fMaxDzInner);
985 //fFracOfMergedRow += (tDu<fMaxDuInner && tDz<fMaxDzInner);
986 }
987 else{
988 fFracOfMergedRow += (tDu<fMaxDuOuter && tDz<fMaxDzOuter);
989 tDist = ::sqrt(tDu*tDu/fMaxDuOuter/fMaxDuOuter+
990 tDz*tDz/fMaxDzOuter/fMaxDzOuter);
991 //fFracOfMergedRow += (tDu<fMaxDuOuter && tDz<fMaxDzOuter);
992 }
993 if(tDist<tDistMax){
994 fClosestRowAtDCA = ti+1;
995 tDistMax = tDist;
996 }
997 fWeightedAvSep += tDist;
998 }
999 }
1000 if(tN>0){
1001 fWeightedAvSep /= tN;
1002 fFracOfMergedRow /= tN;
1003 }
1004 else{
1005 fClosestRowAtDCA = -1;
1006 fFracOfMergedRow = -1.;
1007 fWeightedAvSep = -1.;
1008 }
1009}
1010//________________V0 daughters exit/entrance/average separation calc.
1011//_______1st part is a track 2nd is a V0 considering Pos daughter
1012double AliFemtoPair::TpcExitSeparationTrackV0Pos() const {
1013 AliFemtoThreeVector diff = fTrack1->NominalTpcExitPoint() - fTrack2->TpcV0PosExitPoint();
1014 return (diff.mag());
1015}
1016
1017double AliFemtoPair::TpcEntranceSeparationTrackV0Pos() const {
1018 AliFemtoThreeVector diff = fTrack1->NominalTpcEntrancePoint() - fTrack2->TpcV0PosEntrancePoint();
1019 return (diff.mag());
1020}
1021
1022double AliFemtoPair::TpcAverageSeparationTrackV0Pos() const {
1023 AliFemtoThreeVector diff;
1024 double AveSep = 0.0;
1025 int ipt = 0;
1026 if (fTrack1->fNominalPosSample && fTrack2->fNominalPosSample){
1027 while (fabs(fTrack1->fNominalPosSample[ipt].x())<9999. &&
1028 fabs(fTrack1->fNominalPosSample[ipt].y())<9999. &&
1029 fabs(fTrack1->fNominalPosSample[ipt].z())<9999. &&
1030 fabs(fTrack2->fNominalPosSample[ipt].x())<9999. &&
1031 fabs(fTrack2->fNominalPosSample[ipt].y())<9999. &&
1032 fabs(fTrack2->fNominalPosSample[ipt].z())<9999. &&
1033 (ipt<11)
1034 ){
1035 diff = fTrack1->fNominalPosSample[ipt] - fTrack2->fNominalPosSample[ipt];
1036 ipt++;
1037 AveSep += diff.mag();
1038 }
1039 AveSep = AveSep/(ipt+1.);
1040 return (AveSep);}
1041 else return -1;
1042}
1043//_______1st part is a track 2nd is a V0 considering Neg daughter
1044double AliFemtoPair::TpcExitSeparationTrackV0Neg() const {
1045 AliFemtoThreeVector diff = fTrack1->NominalTpcExitPoint() - fTrack2->TpcV0NegExitPoint();
1046 return (diff.mag());
1047}
1048
1049double AliFemtoPair::TpcEntranceSeparationTrackV0Neg() const {
1050 AliFemtoThreeVector diff = fTrack1->NominalTpcEntrancePoint() - fTrack2->TpcV0NegEntrancePoint();
1051 return (diff.mag());
1052}
1053
1054double AliFemtoPair::TpcAverageSeparationTrackV0Neg() const {
1055 AliFemtoThreeVector diff;
1056 double AveSep = 0.0;
1057 int ipt = 0;
1058 if (fTrack1->fNominalPosSample && fTrack2->fTpcV0NegPosSample){
1059 while (fabs(fTrack1->fNominalPosSample[ipt].x())<9999. &&
1060 fabs(fTrack1->fNominalPosSample[ipt].y())<9999. &&
1061 fabs(fTrack1->fNominalPosSample[ipt].z())<9999. &&
1062 fabs(fTrack2->fTpcV0NegPosSample[ipt].x())<9999. &&
1063 fabs(fTrack2->fTpcV0NegPosSample[ipt].y())<9999. &&
1064 fabs(fTrack2->fTpcV0NegPosSample[ipt].z())<9999. &&
1065 (ipt<11)
1066 ){
1067 diff = fTrack1->fNominalPosSample[ipt] - fTrack2->fTpcV0NegPosSample[ipt];
1068 ipt++;
1069 AveSep += diff.mag();
1070 }
1071 AveSep = AveSep/(ipt+1.);
1072 return (AveSep);}
1073 else return -1;
1074}
1075
1076//_______1st part is a V0 considering Pos daughter 2nd is a V0 considering Pos daughter
1077double AliFemtoPair::TpcExitSeparationV0PosV0Pos() const {
1078 AliFemtoThreeVector diff = fTrack1->TpcV0PosExitPoint() - fTrack2->TpcV0PosExitPoint();
1079 return (diff.mag());
1080}
1081
1082double AliFemtoPair::TpcEntranceSeparationV0PosV0Pos() const {
1083 AliFemtoThreeVector diff = fTrack1->TpcV0PosEntrancePoint() - fTrack2->TpcV0PosEntrancePoint();
1084 return (diff.mag());
1085}
1086double AliFemtoPair::TpcAverageSeparationV0PosV0Pos() const {
1087 AliFemtoThreeVector diff;
1088 double AveSep = 0.0;
1089 int ipt=0;
1090 if (fTrack1->fNominalPosSample && (fTrack2->fNominalPosSample)){
1091 while ((fabs(fTrack1->fNominalPosSample[ipt].x())<9999.) &&
1092 (fabs(fTrack1->fNominalPosSample[ipt].y())<9999.) &&
1093 (fabs(fTrack1->fNominalPosSample[ipt].z())<9999.) &&
1094 (fabs(fTrack2->fNominalPosSample[ipt].x())<9999.) &&
1095 (fabs(fTrack2->fNominalPosSample[ipt].y())<9999.) &&
1096 (fabs(fTrack2->fNominalPosSample[ipt].z())<9999.) &&
1097 (ipt<11)
1098 ){
1099 diff = fTrack1->fNominalPosSample[ipt] - fTrack2->fNominalPosSample[ipt];
1100 ipt++;
1101 AveSep += diff.mag();
1102 }
1103 AveSep = AveSep/(ipt+1);
1104 return (AveSep);}
1105 else return -1;
1106}
1107
1108//_______1st part is a V0 considering Pos daughter 2nd is a V0 considering Neg daughter
1109double AliFemtoPair::TpcExitSeparationV0PosV0Neg() const {
1110 AliFemtoThreeVector diff = fTrack1->TpcV0PosExitPoint() - fTrack2->TpcV0NegExitPoint();
1111 return (diff.mag());
1112}
1113
1114double AliFemtoPair::TpcEntranceSeparationV0PosV0Neg() const {
1115 AliFemtoThreeVector diff = fTrack1->TpcV0PosEntrancePoint() - fTrack2->TpcV0NegEntrancePoint();
1116 return (diff.mag());
1117}
1118double AliFemtoPair::TpcAverageSeparationV0PosV0Neg() const {
1119 AliFemtoThreeVector diff;
1120 double AveSep = 0.0;
1121 int ipt = 0;
1122 if (fTrack1->fNominalPosSample && fTrack2->fTpcV0NegPosSample){
1123 while (fabs(fTrack1->fNominalPosSample[ipt].x())<9999. &&
1124 fabs(fTrack1->fNominalPosSample[ipt].y())<9999. &&
1125 fabs(fTrack1->fNominalPosSample[ipt].z())<9999. &&
1126 fabs(fTrack2->fTpcV0NegPosSample[ipt].x())<9999. &&
1127 fabs(fTrack2->fTpcV0NegPosSample[ipt].y())<9999. &&
1128 fabs(fTrack2->fTpcV0NegPosSample[ipt].z())<9999. &&
1129 (ipt<11)
1130 ){
1131 diff = fTrack1->fNominalPosSample[ipt] - fTrack2->fTpcV0NegPosSample[ipt];
1132 ipt++;
1133 AveSep += diff.mag();
1134 }
1135 AveSep = AveSep/(ipt+1.);
1136 return (AveSep);}
1137 else return -1;
1138}
1139//_______1st part is a V0 considering Neg daughter 2nd is a V0 considering Pos daughter
1140// this is to check the upper case
1141double AliFemtoPair::TpcExitSeparationV0NegV0Pos() const {
1142 AliFemtoThreeVector diff = fTrack1->TpcV0NegExitPoint() - fTrack2->TpcV0PosExitPoint();
1143 return (diff.mag());
1144}
1145
1146double AliFemtoPair::TpcEntranceSeparationV0NegV0Pos() const {
1147 AliFemtoThreeVector diff = fTrack1->TpcV0NegEntrancePoint() - fTrack2->TpcV0PosEntrancePoint();
1148 return (diff.mag());
1149}
1150double AliFemtoPair::TpcAverageSeparationV0NegV0Pos() const {
1151 AliFemtoThreeVector diff;
1152 double AveSep = 0.0;
1153 int ipt = 0;
1154 if ( fTrack1->fTpcV0NegPosSample && fTrack2->fNominalPosSample){
1155 while (fabs(fTrack1->fTpcV0NegPosSample[ipt].x())<9999. &&
1156 fabs(fTrack1->fTpcV0NegPosSample[ipt].y())<9999. &&
1157 fabs(fTrack1->fTpcV0NegPosSample[ipt].z())<9999. &&
1158 fabs(fTrack2->fNominalPosSample[ipt].x())<9999. &&
1159 fabs(fTrack2->fNominalPosSample[ipt].y())<9999. &&
1160 fabs(fTrack2->fNominalPosSample[ipt].z())<9999. &&
1161 (ipt<11)
1162 ){
1163 diff = fTrack1->fTpcV0NegPosSample[ipt] - fTrack2->fNominalPosSample[ipt];
1164 ipt++;
1165 AveSep += diff.mag();
1166 }
1167 AveSep = AveSep/(ipt+1);
1168 return (AveSep);}
1169 else return -1;
1170}
1171//_______1st part is a V0 considering Neg daughter 2nd is a V0 considering Neg daughter
1172double AliFemtoPair::TpcExitSeparationV0NegV0Neg() const {
1173 AliFemtoThreeVector diff = fTrack1->TpcV0NegExitPoint() - fTrack2->TpcV0NegExitPoint();
1174 return (diff.mag());
1175}
1176
1177double AliFemtoPair::TpcEntranceSeparationV0NegV0Neg() const {
1178 AliFemtoThreeVector diff = fTrack1->TpcV0NegEntrancePoint() - fTrack2->TpcV0NegEntrancePoint();
1179 return (diff.mag());
1180}
1181double AliFemtoPair::TpcAverageSeparationV0NegV0Neg() const {
1182 AliFemtoThreeVector diff;
1183 double AveSep = 0.0;
1184 int ipt=0;
1185 if (fTrack1->fTpcV0NegPosSample && fTrack2->fTpcV0NegPosSample){
1186 while (fabs(fTrack1->fTpcV0NegPosSample[ipt].x())<9999. &&
1187 fabs(fTrack1->fTpcV0NegPosSample[ipt].y())<9999. &&
1188 fabs(fTrack1->fTpcV0NegPosSample[ipt].z())<9999. &&
1189 fabs(fTrack2->fTpcV0NegPosSample[ipt].x())<9999. &&
1190 fabs(fTrack2->fTpcV0NegPosSample[ipt].y())<9999. &&
1191 fabs(fTrack2->fTpcV0NegPosSample[ipt].z())<9999. &&
1192 (ipt<11)
1193 ){
1194 diff = fTrack1->fTpcV0NegPosSample[ipt] - fTrack2->fTpcV0NegPosSample[ipt];
1195 ipt++;
1196 AveSep += diff.mag();
1197 }
1198 AveSep = AveSep/(ipt+1);
1199 return (AveSep);}
1200 else return -1;
1201}
1202
1203//________________end V0 daughters exit/entrance/average separation calc.
1204void AliFemtoPair::CalcMergingParFctn(short* tmpMergingParNotCalculatedFctn,
1205 float* tmpZ1,float* tmpU1,
1206 float* tmpZ2,float* tmpU2,
1207 int *tmpSect1,int *tmpSect2,
1208 double* tmpFracOfMergedRow,
1209 double* tmpClosestRowAtDCA
1210 ) const{
1211 tmpMergingParNotCalculatedFctn=0;
1212 double tDu, tDz;
1213 int tN = 0;
1214 *tmpFracOfMergedRow = 0.;
1215 *tmpClosestRowAtDCA = 0.;
1216 double tDist;
1217 double tDistMax = 100000000.;
1218 for(int ti=0 ; ti<45 ; ti++){
1219 if(tmpSect1[ti]==tmpSect2[ti] && tmpSect1[ti]!=-1){
1220 tDu = fabs(tmpU1[ti]-tmpU2[ti]);
1221 tDz = fabs(tmpZ1[ti]-tmpZ2[ti]);
1222 tN++;
1223 if(ti<13){
1224 *tmpFracOfMergedRow += (tDu<fMaxDuInner && tDz<fMaxDzInner);
1225 tDist = ::sqrt(tDu*tDu/fMaxDuInner/fMaxDuInner+
1226 tDz*tDz/fMaxDzInner/fMaxDzInner);
1227 }
1228 else{
1229 *tmpFracOfMergedRow += (tDu<fMaxDuOuter && tDz<fMaxDzOuter);
1230 tDist = ::sqrt(tDu*tDu/fMaxDuOuter/fMaxDuOuter+
1231 tDz*tDz/fMaxDzOuter/fMaxDzOuter);
1232 }
1233 if(tDist<tDistMax){
1234 fClosestRowAtDCA = ti+1;
1235 tDistMax = tDist;
1236 }
1237 //fWeightedAvSep += tDist; // now, wrong but not used
1238 }
1239 }
1240 if(tN>0){
1241 //fWeightedAvSep /= tN;
1242 *tmpFracOfMergedRow /= tN;
1243 }
1244 else{
1245 *tmpClosestRowAtDCA = -1;
1246 *tmpFracOfMergedRow = -1.;
1247 //fWeightedAvSep = -1.;
1248 }
1249}
1250