]> git.uio.no Git - u/mrichter/AliRoot.git/blame - STEER/AliESDv0.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
added new triggers
[u/mrichter/AliRoot.git] / STEER / AliESDv0.cxx
CommitLineData
e23730c7 1/**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15
4f679a16 16/* $Id$ */
17
e23730c7 18//-------------------------------------------------------------------------
19// Implementation of the ESD V0 vertex class
4f679a16 20// This class is part of the Event Data Summary
21// set of classes and contains information about
22// V0 kind vertexes generated by a neutral particle
e23730c7 23// Origin: Iouri Belikov, IReS, Strasbourg, Jouri.Belikov@cern.ch
d6a49f20 24// Modified by: Marian Ivanov, CERN, Marian.Ivanov@cern.ch
25// and Boris Hippolyte,IPHC, hippolyt@in2p3.fr
e23730c7 26//-------------------------------------------------------------------------
4f679a16 27
e23730c7 28#include <TMath.h>
90e48c0c 29#include <TDatabasePDG.h>
90e48c0c 30#include <TParticlePDG.h>
3e27a419 31#include <TVector3.h>
e23730c7 32
5f7789fc 33#include "AliLog.h"
e23730c7 34#include "AliESDv0.h"
97135b34 35#include "AliESDV0Params.h"
e23730c7 36
37ClassImp(AliESDv0)
38
562dd0b4 39const AliESDV0Params AliESDv0::fgkParams;
d6a49f20 40
90e48c0c 41AliESDv0::AliESDv0() :
646c9704 42 AliVParticle(),
d6a49f20 43 fParamN(),
b75d63a7 44 fParamP(),
562dd0b4 45 fEffMass(TDatabasePDG::Instance()->GetParticle(kK0Short)->Mass()),
46 fDcaV0Daughters(0),
47 fChi2V0(0.),
48 fRr(0),
d6a49f20 49 fDistSigma(0),
50 fChi2Before(0),
d6a49f20 51 fChi2After(0),
d6a49f20 52 fPointAngleFi(0),
562dd0b4 53 fPointAngleTh(0),
54 fPointAngle(0),
55 fPdgCode(kK0Short),
56 fNidx(0),
57 fPidx(0),
58 fStatus(0),
59 fNBefore(0),
60 fNAfter(0),
61 fOnFlyStatus(kFALSE)
90e48c0c 62{
e23730c7 63 //--------------------------------------------------------------------
64 // Default constructor (K0s)
65 //--------------------------------------------------------------------
6605de26 66
67 for (Int_t i=0; i<3; i++) {
68 fPos[i] = 0.;
69 fNmom[i] = 0.;
70 fPmom[i] = 0.;
71 }
72
73 for (Int_t i=0; i<6; i++) {
74 fPosCov[i]= 0.;
6605de26 75 }
d6a49f20 76
d6a49f20 77 for (Int_t i=0;i<6;i++){fClusters[0][i]=0; fClusters[1][i]=0;}
78 fNormDCAPrim[0]=fNormDCAPrim[1]=0;
b75d63a7 79 for (Int_t i=0;i<3;i++){fAngle[i]=0;}
d6a49f20 80 for (Int_t i=0;i<4;i++){fCausality[i]=0;}
e23730c7 81}
82
d6a49f20 83AliESDv0::AliESDv0(const AliESDv0& v0) :
646c9704 84 AliVParticle(v0),
562dd0b4 85 fParamN(v0.fParamN),
86 fParamP(v0.fParamP),
d6a49f20 87 fEffMass(v0.fEffMass),
88 fDcaV0Daughters(v0.fDcaV0Daughters),
89 fChi2V0(v0.fChi2V0),
d6a49f20 90 fRr(v0.fRr),
d6a49f20 91 fDistSigma(v0.fDistSigma),
92 fChi2Before(v0.fChi2Before),
d6a49f20 93 fChi2After(v0.fChi2After),
d6a49f20 94 fPointAngleFi(v0.fPointAngleFi),
562dd0b4 95 fPointAngleTh(v0.fPointAngleTh),
96 fPointAngle(v0.fPointAngle),
97 fPdgCode(v0.fPdgCode),
98 fNidx(v0.fNidx),
99 fPidx(v0.fPidx),
100 fStatus(v0.fStatus),
101 fNBefore(v0.fNBefore),
102 fNAfter(v0.fNAfter),
103 fOnFlyStatus(v0.fOnFlyStatus)
c028b974 104{
d6a49f20 105 //--------------------------------------------------------------------
106 // The copy constructor
107 //--------------------------------------------------------------------
c028b974 108
109 for (int i=0; i<3; i++) {
d6a49f20 110 fPos[i] = v0.fPos[i];
111 fNmom[i] = v0.fNmom[i];
112 fPmom[i] = v0.fPmom[i];
c028b974 113 }
114 for (int i=0; i<6; i++) {
d6a49f20 115 fPosCov[i] = v0.fPosCov[i];
c028b974 116 }
c028b974 117
d6a49f20 118 for (Int_t i=0; i<2; i++) {
b75d63a7 119 fNormDCAPrim[i]=v0.fNormDCAPrim[i];
d6a49f20 120 }
121 for (Int_t i=0;i<6;i++){
122 fClusters[0][i]=v0.fClusters[0][i];
123 fClusters[1][i]=v0.fClusters[1][i];
124 }
125 for (Int_t i=0;i<3;i++){
d6a49f20 126 fAngle[i]=v0.fAngle[i];
c028b974 127 }
d6a49f20 128 for (Int_t i=0;i<4;i++){fCausality[i]=v0.fCausality[i];}
c028b974 129}
130
732a24fe 131
c7bafca9 132AliESDv0::AliESDv0(const AliExternalTrackParam &t1, Int_t i1,
133 const AliExternalTrackParam &t2, Int_t i2) :
646c9704 134 AliVParticle(),
562dd0b4 135 fParamN(t1),
136 fParamP(t2),
c7bafca9 137 fEffMass(TDatabasePDG::Instance()->GetParticle(kK0Short)->Mass()),
c028b974 138 fDcaV0Daughters(0),
562dd0b4 139 fChi2V0(0.),
140 fRr(0),
141 fDistSigma(0),
142 fChi2Before(0),
143 fChi2After(0),
144 fPointAngleFi(0),
145 fPointAngleTh(0),
b75d63a7 146 fPointAngle(0),
562dd0b4 147 fPdgCode(kK0Short),
c7bafca9 148 fNidx(i1),
d6a49f20 149 fPidx(i2),
d6a49f20 150 fStatus(0),
d6a49f20 151 fNBefore(0),
d6a49f20 152 fNAfter(0),
562dd0b4 153 fOnFlyStatus(kFALSE)
c7bafca9 154{
155 //--------------------------------------------------------------------
156 // Main constructor (K0s)
157 //--------------------------------------------------------------------
158
159 for (Int_t i=0; i<6; i++) {
160 fPosCov[i]= 0.;
c7bafca9 161 }
162
163 //Trivial estimation of the vertex parameters
b75d63a7 164 Double_t alpha=t1.GetAlpha(), cs=TMath::Cos(alpha), sn=TMath::Sin(alpha);
165 Double_t tmp[3];
166 t1.GetPxPyPz(tmp);
167 Double_t px1=tmp[0], py1=tmp[1], pz1=tmp[2];
168 t1.GetXYZ(tmp);
169 Double_t x1=tmp[0], y1=tmp[1], z1=tmp[2];
37919f0b 170 const Double_t ss=0.0005*0.0005;//a kind of a residual misalignment precision
171 Double_t sx1=sn*sn*t1.GetSigmaY2()+ss, sy1=cs*cs*t1.GetSigmaY2()+ss;
c7bafca9 172
173
b75d63a7 174 alpha=t2.GetAlpha(); cs=TMath::Cos(alpha); sn=TMath::Sin(alpha);
175 t2.GetPxPyPz(tmp);
176 Double_t px2=tmp[0], py2=tmp[1], pz2=tmp[2];
177 t2.GetXYZ(tmp);
178 Double_t x2=tmp[0], y2=tmp[1], z2=tmp[2];
37919f0b 179 Double_t sx2=sn*sn*t2.GetSigmaY2()+ss, sy2=cs*cs*t2.GetSigmaY2()+ss;
c7bafca9 180
181 Double_t sz1=t1.GetSigmaZ2(), sz2=t2.GetSigmaZ2();
37919f0b 182 Double_t wx1=sx2/(sx1+sx2), wx2=1.- wx1;
183 Double_t wy1=sy2/(sy1+sy2), wy2=1.- wy1;
c7bafca9 184 Double_t wz1=sz2/(sz1+sz2), wz2=1.- wz1;
185 fPos[0]=wx1*x1 + wx2*x2; fPos[1]=wy1*y1 + wy2*y2; fPos[2]=wz1*z1 + wz2*z2;
186
187 //fPos[0]=0.5*(x1+x2); fPos[1]=0.5*(y1+y2); fPos[2]=0.5*(z1+z2);
188 fNmom[0]=px1; fNmom[1]=py1; fNmom[2]=pz1;
189 fPmom[0]=px2; fPmom[1]=py2; fPmom[2]=pz2;
190
b75d63a7 191 for (Int_t i=0;i<6;i++){fClusters[0][i]=0; fClusters[1][i]=0;}
192 fNormDCAPrim[0]=fNormDCAPrim[1]=0;
193 for (Int_t i=0;i<3;i++){fAngle[i]=0;}
194 for (Int_t i=0;i<4;i++){fCausality[i]=0;}
c028b974 195}
c7bafca9 196
3e27a419 197AliESDv0& AliESDv0::operator=(const AliESDv0 &v0)
198{
732a24fe 199 //--------------------------------------------------------------------
3e27a419 200 // The assignment operator
732a24fe 201 //--------------------------------------------------------------------
202
203 if(this==&v0)return *this;
646c9704 204 AliVParticle::operator=(v0);
732a24fe 205 fParamN = v0.fParamN;
206 fParamP = v0.fParamP;
207 fEffMass = v0.fEffMass;
208 fDcaV0Daughters = v0.fDcaV0Daughters;
209 fChi2V0 = v0.fChi2V0;
210 fRr = v0.fRr;
211 fDistSigma = v0.fDistSigma;
212 fChi2Before = v0.fChi2Before;
213 fChi2After = v0.fChi2After;
214 fPointAngleFi = v0.fPointAngleFi;
215 fPointAngleTh = v0.fPointAngleTh;
216 fPointAngle = v0.fPointAngle;
217 fPdgCode = v0.fPdgCode;
218 fNidx = v0.fNidx;
219 fPidx = v0.fPidx;
220 fStatus = v0.fStatus;
221 fNBefore = v0.fNBefore;
222 fNAfter = v0.fNAfter;
223 fOnFlyStatus = v0.fOnFlyStatus;
224
225 for (int i=0; i<3; i++) {
226 fPos[i] = v0.fPos[i];
227 fNmom[i] = v0.fNmom[i];
228 fPmom[i] = v0.fPmom[i];
229 }
230 for (int i=0; i<6; i++) {
231 fPosCov[i] = v0.fPosCov[i];
232 }
233 for (Int_t i=0; i<2; i++) {
234 fNormDCAPrim[i]=v0.fNormDCAPrim[i];
235 }
236 for (Int_t i=0;i<6;i++){
237 fClusters[0][i]=v0.fClusters[0][i];
238 fClusters[1][i]=v0.fClusters[1][i];
239 }
240 for (Int_t i=0;i<3;i++){
241 fAngle[i]=v0.fAngle[i];
242 }
243 for (Int_t i=0;i<4;i++){fCausality[i]=v0.fCausality[i];}
244
245 return *this;
246}
247
248void AliESDv0::Copy(TObject& obj) const {
249
3e27a419 250 // this overwrites the virtual TOBject::Copy()
732a24fe 251 // to allow run time copying without casting
252 // in AliESDEvent
253
254 if(this==&obj)return;
255 AliESDv0 *robj = dynamic_cast<AliESDv0*>(&obj);
256 if(!robj)return; // not an aliesesv0
257 *robj = *this;
258}
259
c028b974 260AliESDv0::~AliESDv0(){
261 //--------------------------------------------------------------------
262 // Empty destructor
263 //--------------------------------------------------------------------
c7bafca9 264}
265
646c9704 266// Start with AliVParticle functions
267Double_t AliESDv0::E() const {
268 //--------------------------------------------------------------------
269 // This gives the energy assuming the ChangeMassHypothesis was called
270 //--------------------------------------------------------------------
3e27a419 271 return E(fPdgCode);
272}
273
274Double_t AliESDv0::Y() const {
275 //--------------------------------------------------------------------
276 // This gives the energy assuming the ChangeMassHypothesis was called
277 //--------------------------------------------------------------------
278 return Y(fPdgCode);
279}
280
281// Then extend AliVParticle functions
282Double_t AliESDv0::E(Int_t pdg) const {
283 //--------------------------------------------------------------------
284 // This gives the energy with the particle hypothesis as argument
285 //--------------------------------------------------------------------
286 Double_t mass = TDatabasePDG::Instance()->GetParticle(pdg)->Mass();
646c9704 287 return TMath::Sqrt(mass*mass+P()*P());
288}
c028b974 289
3e27a419 290Double_t AliESDv0::Y(Int_t pdg) const {
291 //--------------------------------------------------------------------
292 // This gives the rapidity with the particle hypothesis as argument
293 //--------------------------------------------------------------------
294 return 0.5*TMath::Log((E(pdg)+Pz())/(E(pdg)-Pz()+1.e-13));
295}
296
297// Now the functions for analysis consistency
298Double_t AliESDv0::RapK0Short() const {
299 //--------------------------------------------------------------------
300 // This gives the pseudorapidity assuming a K0s particle
301 //--------------------------------------------------------------------
302 return Y(kK0Short);
303}
304
305Double_t AliESDv0::RapLambda() const {
306 //--------------------------------------------------------------------
307 // This gives the pseudorapidity assuming a (Anti) Lambda particle
308 //--------------------------------------------------------------------
309 return Y(kLambda0);
310}
311
312Double_t AliESDv0::AlphaV0() const {
313 //--------------------------------------------------------------------
314 // This gives the Armenteros-Podolanski alpha
315 //--------------------------------------------------------------------
316 TVector3 momNeg(fNmom[0],fNmom[1],fNmom[2]);
317 TVector3 momPos(fPmom[0],fPmom[1],fPmom[2]);
318 TVector3 momTot(Px(),Py(),Pz());
319
97135b34 320 Double_t lQlNeg = momNeg.Dot(momTot)/momTot.Mag();
ad07c38b 321 Double_t lQlPos = momPos.Dot(momTot)/momTot.Mag();
3e27a419 322
bc88fc12 323 //return 1.-2./(1.+lQlNeg/lQlPos);
324 return (lQlPos - lQlNeg)/(lQlPos + lQlNeg);
3e27a419 325}
326
327Double_t AliESDv0::PtArmV0() const {
328 //--------------------------------------------------------------------
329 // This gives the Armenteros-Podolanski ptarm
330 //--------------------------------------------------------------------
331 TVector3 momNeg(fNmom[0],fNmom[1],fNmom[2]);
332 TVector3 momTot(Px(),Py(),Pz());
333
334 return momNeg.Perp(momTot);
335}
336
337// Eventually the older functions
e23730c7 338Double_t AliESDv0::ChangeMassHypothesis(Int_t code) {
339 //--------------------------------------------------------------------
340 // This function changes the mass hypothesis for this V0
341 // and returns the "kinematical quality" of this hypothesis
342 //--------------------------------------------------------------------
b75d63a7 343 static
344 Double_t piMass=TDatabasePDG::Instance()->GetParticle(kPiPlus)->Mass();
345 static
346 Double_t prMass=TDatabasePDG::Instance()->GetParticle(kProton)->Mass();
347 static
348 Double_t k0Mass=TDatabasePDG::Instance()->GetParticle(kK0Short)->Mass();
349 static
350 Double_t l0Mass=TDatabasePDG::Instance()->GetParticle(kLambda0)->Mass();
351
352 Double_t nmass=piMass, pmass=piMass, mass=k0Mass, ps=0.206;
e23730c7 353
354 fPdgCode=code;
355
356 switch (code) {
357 case kLambda0:
b75d63a7 358 nmass=piMass; pmass=prMass; mass=l0Mass; ps=0.101; break;
e23730c7 359 case kLambda0Bar:
b75d63a7 360 pmass=piMass; nmass=prMass; mass=l0Mass; ps=0.101; break;
e23730c7 361 case kK0Short:
362 break;
363 default:
5f7789fc 364 AliError("invalide PDG code ! Assuming K0s...");
e23730c7 365 fPdgCode=kK0Short;
366 break;
367 }
368
369 Double_t pxn=fNmom[0], pyn=fNmom[1], pzn=fNmom[2];
370 Double_t pxp=fPmom[0], pyp=fPmom[1], pzp=fPmom[2];
371
372 Double_t en=TMath::Sqrt(nmass*nmass + pxn*pxn + pyn*pyn + pzn*pzn);
373 Double_t ep=TMath::Sqrt(pmass*pmass + pxp*pxp + pyp*pyp + pzp*pzp);
374 Double_t pxl=pxn+pxp, pyl=pyn+pyp, pzl=pzn+pzp;
375 Double_t pl=TMath::Sqrt(pxl*pxl + pyl*pyl + pzl*pzl);
376
377 fEffMass=TMath::Sqrt((en+ep)*(en+ep)-pl*pl);
378
379 Double_t beta=pl/(en+ep);
380 Double_t pln=(pxn*pxl + pyn*pyl + pzn*pzl)/pl;
381 Double_t plp=(pxp*pxl + pyp*pyl + pzp*pzl)/pl;
382
383 Double_t pt2=pxp*pxp + pyp*pyp + pzp*pzp - plp*plp;
384
385 Double_t a=(plp-pln)/(plp+pln);
386 a -= (pmass*pmass-nmass*nmass)/(mass*mass);
387 a = 0.25*beta*beta*mass*mass*a*a + pt2;
388
389 return (a - ps*ps);
390
391}
392
393void AliESDv0::GetPxPyPz(Double_t &px, Double_t &py, Double_t &pz) const {
394 //--------------------------------------------------------------------
395 // This function returns V0's momentum (global)
396 //--------------------------------------------------------------------
397 px=fNmom[0]+fPmom[0];
398 py=fNmom[1]+fPmom[1];
399 pz=fNmom[2]+fPmom[2];
400}
401
402void AliESDv0::GetXYZ(Double_t &x, Double_t &y, Double_t &z) const {
403 //--------------------------------------------------------------------
404 // This function returns V0's position (global)
405 //--------------------------------------------------------------------
406 x=fPos[0];
407 y=fPos[1];
408 z=fPos[2];
409}
410
b75d63a7 411Float_t AliESDv0::GetD(Double_t x0, Double_t y0, Double_t z0) const {
e23730c7 412 //--------------------------------------------------------------------
413 // This function returns V0's impact parameter
414 //--------------------------------------------------------------------
415 Double_t x=fPos[0],y=fPos[1],z=fPos[2];
416 Double_t px=fNmom[0]+fPmom[0];
417 Double_t py=fNmom[1]+fPmom[1];
418 Double_t pz=fNmom[2]+fPmom[2];
419
420 Double_t dx=(y0-y)*pz - (z0-z)*py;
421 Double_t dy=(x0-x)*pz - (z0-z)*px;
422 Double_t dz=(x0-x)*py - (y0-y)*px;
423 Double_t d=TMath::Sqrt((dx*dx+dy*dy+dz*dz)/(px*px+py*py+pz*pz));
424 return d;
425}
c028b974 426
97135b34 427Float_t AliESDv0::GetV0CosineOfPointingAngle(Double_t refPointX, Double_t refPointY, Double_t refPointZ) const {
c028b974 428 // calculates the pointing angle of the V0 wrt a reference point
429
430 Double_t momV0[3]; //momentum of the V0
431 GetPxPyPz(momV0[0],momV0[1],momV0[2]);
432
433 Double_t deltaPos[3]; //vector between the reference point and the V0 vertex
434 deltaPos[0] = fPos[0] - refPointX;
435 deltaPos[1] = fPos[1] - refPointY;
436 deltaPos[2] = fPos[2] - refPointZ;
437
438 Double_t momV02 = momV0[0]*momV0[0] + momV0[1]*momV0[1] + momV0[2]*momV0[2];
439 Double_t deltaPos2 = deltaPos[0]*deltaPos[0] + deltaPos[1]*deltaPos[1] + deltaPos[2]*deltaPos[2];
440
441 Double_t cosinePointingAngle = (deltaPos[0]*momV0[0] +
442 deltaPos[1]*momV0[1] +
443 deltaPos[2]*momV0[2] ) /
444 TMath::Sqrt(momV02 * deltaPos2);
445
446 return cosinePointingAngle;
447}
d6a49f20 448
449
450// **** The following functions need to be revised
451
074f017b 452void AliESDv0::GetPosCov(Double_t cov[6]) const {
453
454 for (Int_t i=0; i<6; ++i) cov[i] = fPosCov[i];
455
456}
457
d6a49f20 458Double_t AliESDv0::GetSigmaY(){
459 //
460 // return sigmay in y at vertex position using covariance matrix
461 //
462 const Double_t * cp = fParamP.GetCovariance();
463 const Double_t * cm = fParamN.GetCovariance();
464 Double_t sigmay = cp[0]+cm[0]+ cp[5]*(fParamP.GetX()-fRr)*(fParamP.GetX()-fRr)+ cm[5]*(fParamN.GetX()-fRr)*(fParamN.GetX()-fRr);
465 return (sigmay>0) ? TMath::Sqrt(sigmay):100;
466}
467
468Double_t AliESDv0::GetSigmaZ(){
469 //
470 // return sigmay in y at vertex position using covariance matrix
471 //
472 const Double_t * cp = fParamP.GetCovariance();
473 const Double_t * cm = fParamN.GetCovariance();
474 Double_t sigmaz = cp[2]+cm[2]+ cp[9]*(fParamP.GetX()-fRr)*(fParamP.GetX()-fRr)+ cm[9]*(fParamN.GetX()-fRr)*(fParamN.GetX()-fRr);
475 return (sigmaz>0) ? TMath::Sqrt(sigmaz):100;
476}
477
478Double_t AliESDv0::GetSigmaD0(){
479 //
480 // Sigma parameterization using covariance matrix
481 //
482 // sigma of distance between two tracks in vertex position
483 // sigma of DCA is proportianal to sigmaD0
484 // factor 2 difference is explained by the fact that the DCA is calculated at the position
485 // where the tracks as closest together ( not exact position of the vertex)
486 //
487 const Double_t * cp = fParamP.GetCovariance();
488 const Double_t * cm = fParamN.GetCovariance();
489 Double_t sigmaD0 = cp[0]+cm[0]+cp[2]+cm[2]+fgkParams.fPSigmaOffsetD0*fgkParams.fPSigmaOffsetD0;
490 sigmaD0 += ((fParamP.GetX()-fRr)*(fParamP.GetX()-fRr))*(cp[5]+cp[9]);
491 sigmaD0 += ((fParamN.GetX()-fRr)*(fParamN.GetX()-fRr))*(cm[5]+cm[9]);
492 return (sigmaD0>0)? TMath::Sqrt(sigmaD0):100;
493}
494
495
496Double_t AliESDv0::GetSigmaAP0(){
497 //
498 //Sigma parameterization using covariance matrices
499 //
b75d63a7 500 Double_t prec = TMath::Sqrt((fNmom[0]+fPmom[0])*(fNmom[0]+fPmom[0])
501 +(fNmom[1]+fPmom[1])*(fNmom[1]+fPmom[1])
502 +(fNmom[2]+fPmom[2])*(fNmom[2]+fPmom[2]));
503 Double_t normp = TMath::Sqrt(fPmom[0]*fPmom[0]+fPmom[1]*fPmom[1]+fPmom[2]*fPmom[2])/prec; // fraction of the momenta
504 Double_t normm = TMath::Sqrt(fNmom[0]*fNmom[0]+fNmom[1]*fNmom[1]+fNmom[2]*fNmom[2])/prec;
d6a49f20 505 const Double_t * cp = fParamP.GetCovariance();
506 const Double_t * cm = fParamN.GetCovariance();
507 Double_t sigmaAP0 = fgkParams.fPSigmaOffsetAP0*fgkParams.fPSigmaOffsetAP0; // minimal part
508 sigmaAP0 += (cp[5]+cp[9])*(normp*normp)+(cm[5]+cm[9])*(normm*normm); // angular resolution part
509 Double_t sigmaAP1 = GetSigmaD0()/(TMath::Abs(fRr)+0.01); // vertex position part
510 sigmaAP0 += 0.5*sigmaAP1*sigmaAP1;
511 return (sigmaAP0>0)? TMath::Sqrt(sigmaAP0):100;
512}
513
514Double_t AliESDv0::GetEffectiveSigmaD0(){
515 //
516 // minimax - effective Sigma parameterization
517 // p12 effective curvature and v0 radius postion used as parameters
518 //
519 Double_t p12 = TMath::Sqrt(fParamP.GetParameter()[4]*fParamP.GetParameter()[4]+
520 fParamN.GetParameter()[4]*fParamN.GetParameter()[4]);
521 Double_t sigmaED0= TMath::Max(TMath::Sqrt(fRr)-fgkParams.fPSigmaRminDE,0.0)*fgkParams.fPSigmaCoefDE*p12*p12;
522 sigmaED0*= sigmaED0;
523 sigmaED0*= sigmaED0;
524 sigmaED0 = TMath::Sqrt(sigmaED0+fgkParams.fPSigmaOffsetDE*fgkParams.fPSigmaOffsetDE);
525 return (sigmaED0<fgkParams.fPSigmaMaxDE) ? sigmaED0: fgkParams.fPSigmaMaxDE;
526}
527
528
529Double_t AliESDv0::GetEffectiveSigmaAP0(){
530 //
531 // effective Sigma parameterization of point angle resolution
532 //
533 Double_t p12 = TMath::Sqrt(fParamP.GetParameter()[4]*fParamP.GetParameter()[4]+
534 fParamN.GetParameter()[4]*fParamN.GetParameter()[4]);
535 Double_t sigmaAPE= fgkParams.fPSigmaBase0APE;
536 sigmaAPE+= fgkParams.fPSigmaR0APE/(fgkParams.fPSigmaR1APE+fRr);
537 sigmaAPE*= (fgkParams.fPSigmaP0APE+fgkParams.fPSigmaP1APE*p12);
538 sigmaAPE = TMath::Min(sigmaAPE,fgkParams.fPSigmaMaxAPE);
539 return sigmaAPE;
540}
541
542
543Double_t AliESDv0::GetMinimaxSigmaAP0(){
544 //
545 // calculate mini-max effective sigma of point angle resolution
546 //
547 //compv0->fTree->SetAlias("SigmaAP2","max(min((SigmaAP0+SigmaAPE0)*0.5,1.5*SigmaAPE0),0.5*SigmaAPE0+0.003)");
548 Double_t effectiveSigma = GetEffectiveSigmaAP0();
549 Double_t sigmaMMAP = 0.5*(GetSigmaAP0()+effectiveSigma);
550 sigmaMMAP = TMath::Min(sigmaMMAP, fgkParams.fPMaxFractionAP0*effectiveSigma);
551 sigmaMMAP = TMath::Max(sigmaMMAP, fgkParams.fPMinFractionAP0*effectiveSigma+fgkParams.fPMinAP0);
552 return sigmaMMAP;
553}
554Double_t AliESDv0::GetMinimaxSigmaD0(){
555 //
556 // calculate mini-max sigma of dca resolution
557 //
558 //compv0->fTree->SetAlias("SigmaD2","max(min((SigmaD0+SigmaDE0)*0.5,1.5*SigmaDE0),0.5*SigmaDE0)");
559 Double_t effectiveSigma = GetEffectiveSigmaD0();
560 Double_t sigmaMMD0 = 0.5*(GetSigmaD0()+effectiveSigma);
561 sigmaMMD0 = TMath::Min(sigmaMMD0, fgkParams.fPMaxFractionD0*effectiveSigma);
562 sigmaMMD0 = TMath::Max(sigmaMMD0, fgkParams.fPMinFractionD0*effectiveSigma+fgkParams.fPMinD0);
563 return sigmaMMD0;
564}
565
566
567Double_t AliESDv0::GetLikelihoodAP(Int_t mode0, Int_t mode1){
568 //
569 // get likelihood for point angle
570 //
571 Double_t sigmaAP = 0.007; //default sigma
572 switch (mode0){
573 case 0:
574 sigmaAP = GetSigmaAP0(); // mode 0 - covariance matrix estimates used
575 break;
576 case 1:
577 sigmaAP = GetEffectiveSigmaAP0(); // mode 1 - effective sigma used
578 break;
579 case 2:
580 sigmaAP = GetMinimaxSigmaAP0(); // mode 2 - minimax sigma
581 break;
582 }
583 Double_t apNorm = TMath::Min(TMath::ACos(fPointAngle)/sigmaAP,50.);
584 //normalized point angle, restricted - because of overflow problems in Exp
585 Double_t likelihood = 0;
586 switch(mode1){
587 case 0:
588 likelihood = TMath::Exp(-0.5*apNorm*apNorm);
589 // one component
590 break;
591 case 1:
592 likelihood = (TMath::Exp(-0.5*apNorm*apNorm)+0.5* TMath::Exp(-0.25*apNorm*apNorm))/1.5;
593 // two components
594 break;
595 case 2:
596 likelihood = (TMath::Exp(-0.5*apNorm*apNorm)+0.5* TMath::Exp(-0.25*apNorm*apNorm)+0.25*TMath::Exp(-0.125*apNorm*apNorm))/1.75;
597 // three components
598 break;
599 }
600 return likelihood;
601}
602
603Double_t AliESDv0::GetLikelihoodD(Int_t mode0, Int_t mode1){
604 //
605 // get likelihood for DCA
606 //
607 Double_t sigmaD = 0.03; //default sigma
608 switch (mode0){
609 case 0:
610 sigmaD = GetSigmaD0(); // mode 0 - covariance matrix estimates used
611 break;
612 case 1:
613 sigmaD = GetEffectiveSigmaD0(); // mode 1 - effective sigma used
614 break;
615 case 2:
616 sigmaD = GetMinimaxSigmaD0(); // mode 2 - minimax sigma
617 break;
618 }
b75d63a7 619
620 //Bo: Double_t dNorm = TMath::Min(fDist2/sigmaD,50.);
621 Double_t dNorm = TMath::Min(fDcaV0Daughters/sigmaD,50.);//Bo:
d6a49f20 622 //normalized point angle, restricted - because of overflow problems in Exp
623 Double_t likelihood = 0;
624 switch(mode1){
625 case 0:
626 likelihood = TMath::Exp(-2.*dNorm);
627 // one component
628 break;
629 case 1:
630 likelihood = (TMath::Exp(-2.*dNorm)+0.5* TMath::Exp(-dNorm))/1.5;
631 // two components
632 break;
633 case 2:
634 likelihood = (TMath::Exp(-2.*dNorm)+0.5* TMath::Exp(-dNorm)+0.25*TMath::Exp(-0.5*dNorm))/1.75;
635 // three components
636 break;
637 }
638 return likelihood;
639
640}
641
97135b34 642Double_t AliESDv0::GetLikelihoodC(Int_t mode0, Int_t /*mode1*/) const {
d6a49f20 643 //
644 // get likelihood for Causality
645 // !!! Causality variables defined in AliITStrackerMI !!!
646 // when more information was available
647 //
648 Double_t likelihood = 0.5;
649 Double_t minCausal = TMath::Min(fCausality[0],fCausality[1]);
650 Double_t maxCausal = TMath::Max(fCausality[0],fCausality[1]);
651 // minCausal = TMath::Max(minCausal,0.5*maxCausal);
652 //compv0->fTree->SetAlias("LCausal","(1.05-(2*(0.8-exp(-max(RC.fV0rec.fCausality[0],RC.fV0rec.fCausality[1])))+2*(0.8-exp(-min(RC.fV0rec.fCausality[0],RC.fV0rec.fCausality[1]))))/2)**4");
653
654 switch(mode0){
655 case 0:
656 //normalization
657 likelihood = TMath::Power((1.05-2*(0.8-TMath::Exp(-maxCausal))),4.);
658 break;
659 case 1:
660 likelihood = TMath::Power(1.05-(2*(0.8-TMath::Exp(-maxCausal))+(2*(0.8-TMath::Exp(-minCausal))))*0.5,4.);
661 break;
662 }
663 return likelihood;
664
665}
666
667void AliESDv0::SetCausality(Float_t pb0, Float_t pb1, Float_t pa0, Float_t pa1)
668{
669 //
670 // set probabilities
671 //
672 fCausality[0] = pb0; // probability - track 0 exist before vertex
673 fCausality[1] = pb1; // probability - track 1 exist before vertex
674 fCausality[2] = pa0; // probability - track 0 exist close after vertex
675 fCausality[3] = pa1; // probability - track 1 exist close after vertex
676}
97135b34 677void AliESDv0::SetClusters(const Int_t *clp, const Int_t *clm)
d6a49f20 678{
679 //
680 // Set its clusters indexes
681 //
682 for (Int_t i=0;i<6;i++) fClusters[0][i] = clp[i];
683 for (Int_t i=0;i<6;i++) fClusters[1][i] = clm[i];
684}
9973a4bf 685
97135b34 686Double_t AliESDv0::GetEffMass(UInt_t p1, UInt_t p2) const{
9973a4bf 687 //
688 // calculate effective mass
689 //
89cb3a25 690 const Double_t kpmass[5] = {TDatabasePDG::Instance()->GetParticle(kElectron)->Mass(),
586d906c 691 TDatabasePDG::Instance()->GetParticle(kMuonMinus)->Mass(),
692 TDatabasePDG::Instance()->GetParticle(kPiPlus)->Mass(),
693 TDatabasePDG::Instance()->GetParticle(kKPlus)->Mass(),
694 TDatabasePDG::Instance()->GetParticle(kProton)->Mass()};
2278bbc8 695 /*
9973a4bf 696 if (p1>4) return -1;
697 if (p2>4) return -1;
89cb3a25 698 Double_t mass1 = kpmass[p1];
699 Double_t mass2 = kpmass[p2];
97135b34 700 const Double_t *m1 = fPmom;
701 const Double_t *m2 = fNmom;
9973a4bf 702 //
703 //if (fRP[p1]+fRM[p2]<fRP[p2]+fRM[p1]){
704 // m1 = fPM;
705 // m2 = fPP;
706 //}
707 //
89cb3a25 708 Double_t e1 = TMath::Sqrt(mass1*mass1+
9973a4bf 709 m1[0]*m1[0]+
710 m1[1]*m1[1]+
711 m1[2]*m1[2]);
89cb3a25 712 Double_t e2 = TMath::Sqrt(mass2*mass2+
9973a4bf 713 m2[0]*m2[0]+
714 m2[1]*m2[1]+
715 m2[2]*m2[2]);
89cb3a25 716 Double_t mass =
9973a4bf 717 (m2[0]+m1[0])*(m2[0]+m1[0])+
718 (m2[1]+m1[1])*(m2[1]+m1[1])+
719 (m2[2]+m1[2])*(m2[2]+m1[2]);
720
89cb3a25 721 mass = (e1+e2)*(e1+e2)-mass;
722 if (mass < 0.) mass = 0.;
723 return (TMath::Sqrt(mass));
2278bbc8 724 */
725 if(p1>4 || p2>4) return -1;
726 Double_t e12 = kpmass[p1]*kpmass[p1]+fPmom[0]*fPmom[0]+fPmom[1]*fPmom[1]+fPmom[2]*fPmom[2];
727 Double_t e22 = kpmass[p2]*kpmass[p2]+fNmom[0]*fNmom[0]+fNmom[1]*fNmom[1]+fNmom[2]*fNmom[2];
728 Double_t cmass = TMath::Sqrt(TMath::Max(kpmass[p1]*kpmass[p1]+kpmass[p2]*kpmass[p2]
729 +2.*(TMath::Sqrt(e12*e22)-fPmom[0]*fNmom[0]-fPmom[1]*fNmom[1]-fPmom[2]*fNmom[2]),0.));
730 return cmass;
731
9973a4bf 732}