Added methods for finding a given module in the DDL map (F. Prino)
[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 <Riostream.h>
29#include <TMath.h>
90e48c0c 30#include <TDatabasePDG.h>
e23730c7 31#include <TPDGCode.h>
90e48c0c 32#include <TParticlePDG.h>
e23730c7 33
5f7789fc 34#include "AliLog.h"
e23730c7 35#include "AliESDv0.h"
c7bafca9 36#include "AliExternalTrackParam.h"
e23730c7 37
38ClassImp(AliESDv0)
39
562dd0b4 40const AliESDV0Params AliESDv0::fgkParams;
d6a49f20 41
90e48c0c 42AliESDv0::AliESDv0() :
43 TObject(),
d6a49f20 44 fParamN(),
b75d63a7 45 fParamP(),
562dd0b4 46 fEffMass(TDatabasePDG::Instance()->GetParticle(kK0Short)->Mass()),
47 fDcaV0Daughters(0),
48 fChi2V0(0.),
49 fRr(0),
d6a49f20 50 fDistSigma(0),
51 fChi2Before(0),
d6a49f20 52 fChi2After(0),
d6a49f20 53 fPointAngleFi(0),
562dd0b4 54 fPointAngleTh(0),
55 fPointAngle(0),
56 fPdgCode(kK0Short),
57 fNidx(0),
58 fPidx(0),
59 fStatus(0),
60 fNBefore(0),
61 fNAfter(0),
62 fOnFlyStatus(kFALSE)
90e48c0c 63{
e23730c7 64 //--------------------------------------------------------------------
65 // Default constructor (K0s)
66 //--------------------------------------------------------------------
6605de26 67
68 for (Int_t i=0; i<3; i++) {
69 fPos[i] = 0.;
70 fNmom[i] = 0.;
71 fPmom[i] = 0.;
72 }
73
74 for (Int_t i=0; i<6; i++) {
75 fPosCov[i]= 0.;
6605de26 76 }
d6a49f20 77
d6a49f20 78 for (Int_t i=0;i<6;i++){fClusters[0][i]=0; fClusters[1][i]=0;}
79 fNormDCAPrim[0]=fNormDCAPrim[1]=0;
b75d63a7 80 for (Int_t i=0;i<3;i++){fAngle[i]=0;}
d6a49f20 81 for (Int_t i=0;i<4;i++){fCausality[i]=0;}
e23730c7 82}
83
d6a49f20 84AliESDv0::AliESDv0(const AliESDv0& v0) :
85 TObject(v0),
562dd0b4 86 fParamN(v0.fParamN),
87 fParamP(v0.fParamP),
d6a49f20 88 fEffMass(v0.fEffMass),
89 fDcaV0Daughters(v0.fDcaV0Daughters),
90 fChi2V0(v0.fChi2V0),
d6a49f20 91 fRr(v0.fRr),
d6a49f20 92 fDistSigma(v0.fDistSigma),
93 fChi2Before(v0.fChi2Before),
d6a49f20 94 fChi2After(v0.fChi2After),
d6a49f20 95 fPointAngleFi(v0.fPointAngleFi),
562dd0b4 96 fPointAngleTh(v0.fPointAngleTh),
97 fPointAngle(v0.fPointAngle),
98 fPdgCode(v0.fPdgCode),
99 fNidx(v0.fNidx),
100 fPidx(v0.fPidx),
101 fStatus(v0.fStatus),
102 fNBefore(v0.fNBefore),
103 fNAfter(v0.fNAfter),
104 fOnFlyStatus(v0.fOnFlyStatus)
c028b974 105{
d6a49f20 106 //--------------------------------------------------------------------
107 // The copy constructor
108 //--------------------------------------------------------------------
c028b974 109
110 for (int i=0; i<3; i++) {
d6a49f20 111 fPos[i] = v0.fPos[i];
112 fNmom[i] = v0.fNmom[i];
113 fPmom[i] = v0.fPmom[i];
c028b974 114 }
115 for (int i=0; i<6; i++) {
d6a49f20 116 fPosCov[i] = v0.fPosCov[i];
c028b974 117 }
c028b974 118
d6a49f20 119 for (Int_t i=0; i<2; i++) {
b75d63a7 120 fNormDCAPrim[i]=v0.fNormDCAPrim[i];
d6a49f20 121 }
122 for (Int_t i=0;i<6;i++){
123 fClusters[0][i]=v0.fClusters[0][i];
124 fClusters[1][i]=v0.fClusters[1][i];
125 }
126 for (Int_t i=0;i<3;i++){
d6a49f20 127 fAngle[i]=v0.fAngle[i];
c028b974 128 }
d6a49f20 129 for (Int_t i=0;i<4;i++){fCausality[i]=v0.fCausality[i];}
c028b974 130}
131
c7bafca9 132AliESDv0::AliESDv0(const AliExternalTrackParam &t1, Int_t i1,
133 const AliExternalTrackParam &t2, Int_t i2) :
134 TObject(),
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
c028b974 197AliESDv0::~AliESDv0(){
198 //--------------------------------------------------------------------
199 // Empty destructor
200 //--------------------------------------------------------------------
c7bafca9 201}
202
c028b974 203
204
e23730c7 205Double_t AliESDv0::ChangeMassHypothesis(Int_t code) {
206 //--------------------------------------------------------------------
207 // This function changes the mass hypothesis for this V0
208 // and returns the "kinematical quality" of this hypothesis
209 //--------------------------------------------------------------------
b75d63a7 210 static
211 Double_t piMass=TDatabasePDG::Instance()->GetParticle(kPiPlus)->Mass();
212 static
213 Double_t prMass=TDatabasePDG::Instance()->GetParticle(kProton)->Mass();
214 static
215 Double_t k0Mass=TDatabasePDG::Instance()->GetParticle(kK0Short)->Mass();
216 static
217 Double_t l0Mass=TDatabasePDG::Instance()->GetParticle(kLambda0)->Mass();
218
219 Double_t nmass=piMass, pmass=piMass, mass=k0Mass, ps=0.206;
e23730c7 220
221 fPdgCode=code;
222
223 switch (code) {
224 case kLambda0:
b75d63a7 225 nmass=piMass; pmass=prMass; mass=l0Mass; ps=0.101; break;
e23730c7 226 case kLambda0Bar:
b75d63a7 227 pmass=piMass; nmass=prMass; mass=l0Mass; ps=0.101; break;
e23730c7 228 case kK0Short:
229 break;
230 default:
5f7789fc 231 AliError("invalide PDG code ! Assuming K0s...");
e23730c7 232 fPdgCode=kK0Short;
233 break;
234 }
235
236 Double_t pxn=fNmom[0], pyn=fNmom[1], pzn=fNmom[2];
237 Double_t pxp=fPmom[0], pyp=fPmom[1], pzp=fPmom[2];
238
239 Double_t en=TMath::Sqrt(nmass*nmass + pxn*pxn + pyn*pyn + pzn*pzn);
240 Double_t ep=TMath::Sqrt(pmass*pmass + pxp*pxp + pyp*pyp + pzp*pzp);
241 Double_t pxl=pxn+pxp, pyl=pyn+pyp, pzl=pzn+pzp;
242 Double_t pl=TMath::Sqrt(pxl*pxl + pyl*pyl + pzl*pzl);
243
244 fEffMass=TMath::Sqrt((en+ep)*(en+ep)-pl*pl);
245
246 Double_t beta=pl/(en+ep);
247 Double_t pln=(pxn*pxl + pyn*pyl + pzn*pzl)/pl;
248 Double_t plp=(pxp*pxl + pyp*pyl + pzp*pzl)/pl;
249
250 Double_t pt2=pxp*pxp + pyp*pyp + pzp*pzp - plp*plp;
251
252 Double_t a=(plp-pln)/(plp+pln);
253 a -= (pmass*pmass-nmass*nmass)/(mass*mass);
254 a = 0.25*beta*beta*mass*mass*a*a + pt2;
255
256 return (a - ps*ps);
257
258}
259
260void AliESDv0::GetPxPyPz(Double_t &px, Double_t &py, Double_t &pz) const {
261 //--------------------------------------------------------------------
262 // This function returns V0's momentum (global)
263 //--------------------------------------------------------------------
264 px=fNmom[0]+fPmom[0];
265 py=fNmom[1]+fPmom[1];
266 pz=fNmom[2]+fPmom[2];
267}
268
269void AliESDv0::GetXYZ(Double_t &x, Double_t &y, Double_t &z) const {
270 //--------------------------------------------------------------------
271 // This function returns V0's position (global)
272 //--------------------------------------------------------------------
273 x=fPos[0];
274 y=fPos[1];
275 z=fPos[2];
276}
277
b75d63a7 278Float_t AliESDv0::GetD(Double_t x0, Double_t y0, Double_t z0) const {
e23730c7 279 //--------------------------------------------------------------------
280 // This function returns V0's impact parameter
281 //--------------------------------------------------------------------
282 Double_t x=fPos[0],y=fPos[1],z=fPos[2];
283 Double_t px=fNmom[0]+fPmom[0];
284 Double_t py=fNmom[1]+fPmom[1];
285 Double_t pz=fNmom[2]+fPmom[2];
286
287 Double_t dx=(y0-y)*pz - (z0-z)*py;
288 Double_t dy=(x0-x)*pz - (z0-z)*px;
289 Double_t dz=(x0-x)*py - (y0-y)*px;
290 Double_t d=TMath::Sqrt((dx*dx+dy*dy+dz*dz)/(px*px+py*py+pz*pz));
291 return d;
292}
c028b974 293
294
b75d63a7 295Float_t AliESDv0::GetV0CosineOfPointingAngle(Double_t& refPointX, Double_t& refPointY, Double_t& refPointZ) const {
c028b974 296 // calculates the pointing angle of the V0 wrt a reference point
297
298 Double_t momV0[3]; //momentum of the V0
299 GetPxPyPz(momV0[0],momV0[1],momV0[2]);
300
301 Double_t deltaPos[3]; //vector between the reference point and the V0 vertex
302 deltaPos[0] = fPos[0] - refPointX;
303 deltaPos[1] = fPos[1] - refPointY;
304 deltaPos[2] = fPos[2] - refPointZ;
305
306 Double_t momV02 = momV0[0]*momV0[0] + momV0[1]*momV0[1] + momV0[2]*momV0[2];
307 Double_t deltaPos2 = deltaPos[0]*deltaPos[0] + deltaPos[1]*deltaPos[1] + deltaPos[2]*deltaPos[2];
308
309 Double_t cosinePointingAngle = (deltaPos[0]*momV0[0] +
310 deltaPos[1]*momV0[1] +
311 deltaPos[2]*momV0[2] ) /
312 TMath::Sqrt(momV02 * deltaPos2);
313
314 return cosinePointingAngle;
315}
d6a49f20 316
317
318// **** The following functions need to be revised
319
074f017b 320void AliESDv0::GetPosCov(Double_t cov[6]) const {
321
322 for (Int_t i=0; i<6; ++i) cov[i] = fPosCov[i];
323
324}
325
d6a49f20 326Double_t AliESDv0::GetSigmaY(){
327 //
328 // return sigmay in y at vertex position using covariance matrix
329 //
330 const Double_t * cp = fParamP.GetCovariance();
331 const Double_t * cm = fParamN.GetCovariance();
332 Double_t sigmay = cp[0]+cm[0]+ cp[5]*(fParamP.GetX()-fRr)*(fParamP.GetX()-fRr)+ cm[5]*(fParamN.GetX()-fRr)*(fParamN.GetX()-fRr);
333 return (sigmay>0) ? TMath::Sqrt(sigmay):100;
334}
335
336Double_t AliESDv0::GetSigmaZ(){
337 //
338 // return sigmay in y at vertex position using covariance matrix
339 //
340 const Double_t * cp = fParamP.GetCovariance();
341 const Double_t * cm = fParamN.GetCovariance();
342 Double_t sigmaz = cp[2]+cm[2]+ cp[9]*(fParamP.GetX()-fRr)*(fParamP.GetX()-fRr)+ cm[9]*(fParamN.GetX()-fRr)*(fParamN.GetX()-fRr);
343 return (sigmaz>0) ? TMath::Sqrt(sigmaz):100;
344}
345
346Double_t AliESDv0::GetSigmaD0(){
347 //
348 // Sigma parameterization using covariance matrix
349 //
350 // sigma of distance between two tracks in vertex position
351 // sigma of DCA is proportianal to sigmaD0
352 // factor 2 difference is explained by the fact that the DCA is calculated at the position
353 // where the tracks as closest together ( not exact position of the vertex)
354 //
355 const Double_t * cp = fParamP.GetCovariance();
356 const Double_t * cm = fParamN.GetCovariance();
357 Double_t sigmaD0 = cp[0]+cm[0]+cp[2]+cm[2]+fgkParams.fPSigmaOffsetD0*fgkParams.fPSigmaOffsetD0;
358 sigmaD0 += ((fParamP.GetX()-fRr)*(fParamP.GetX()-fRr))*(cp[5]+cp[9]);
359 sigmaD0 += ((fParamN.GetX()-fRr)*(fParamN.GetX()-fRr))*(cm[5]+cm[9]);
360 return (sigmaD0>0)? TMath::Sqrt(sigmaD0):100;
361}
362
363
364Double_t AliESDv0::GetSigmaAP0(){
365 //
366 //Sigma parameterization using covariance matrices
367 //
b75d63a7 368 Double_t prec = TMath::Sqrt((fNmom[0]+fPmom[0])*(fNmom[0]+fPmom[0])
369 +(fNmom[1]+fPmom[1])*(fNmom[1]+fPmom[1])
370 +(fNmom[2]+fPmom[2])*(fNmom[2]+fPmom[2]));
371 Double_t normp = TMath::Sqrt(fPmom[0]*fPmom[0]+fPmom[1]*fPmom[1]+fPmom[2]*fPmom[2])/prec; // fraction of the momenta
372 Double_t normm = TMath::Sqrt(fNmom[0]*fNmom[0]+fNmom[1]*fNmom[1]+fNmom[2]*fNmom[2])/prec;
d6a49f20 373 const Double_t * cp = fParamP.GetCovariance();
374 const Double_t * cm = fParamN.GetCovariance();
375 Double_t sigmaAP0 = fgkParams.fPSigmaOffsetAP0*fgkParams.fPSigmaOffsetAP0; // minimal part
376 sigmaAP0 += (cp[5]+cp[9])*(normp*normp)+(cm[5]+cm[9])*(normm*normm); // angular resolution part
377 Double_t sigmaAP1 = GetSigmaD0()/(TMath::Abs(fRr)+0.01); // vertex position part
378 sigmaAP0 += 0.5*sigmaAP1*sigmaAP1;
379 return (sigmaAP0>0)? TMath::Sqrt(sigmaAP0):100;
380}
381
382Double_t AliESDv0::GetEffectiveSigmaD0(){
383 //
384 // minimax - effective Sigma parameterization
385 // p12 effective curvature and v0 radius postion used as parameters
386 //
387 Double_t p12 = TMath::Sqrt(fParamP.GetParameter()[4]*fParamP.GetParameter()[4]+
388 fParamN.GetParameter()[4]*fParamN.GetParameter()[4]);
389 Double_t sigmaED0= TMath::Max(TMath::Sqrt(fRr)-fgkParams.fPSigmaRminDE,0.0)*fgkParams.fPSigmaCoefDE*p12*p12;
390 sigmaED0*= sigmaED0;
391 sigmaED0*= sigmaED0;
392 sigmaED0 = TMath::Sqrt(sigmaED0+fgkParams.fPSigmaOffsetDE*fgkParams.fPSigmaOffsetDE);
393 return (sigmaED0<fgkParams.fPSigmaMaxDE) ? sigmaED0: fgkParams.fPSigmaMaxDE;
394}
395
396
397Double_t AliESDv0::GetEffectiveSigmaAP0(){
398 //
399 // effective Sigma parameterization of point angle resolution
400 //
401 Double_t p12 = TMath::Sqrt(fParamP.GetParameter()[4]*fParamP.GetParameter()[4]+
402 fParamN.GetParameter()[4]*fParamN.GetParameter()[4]);
403 Double_t sigmaAPE= fgkParams.fPSigmaBase0APE;
404 sigmaAPE+= fgkParams.fPSigmaR0APE/(fgkParams.fPSigmaR1APE+fRr);
405 sigmaAPE*= (fgkParams.fPSigmaP0APE+fgkParams.fPSigmaP1APE*p12);
406 sigmaAPE = TMath::Min(sigmaAPE,fgkParams.fPSigmaMaxAPE);
407 return sigmaAPE;
408}
409
410
411Double_t AliESDv0::GetMinimaxSigmaAP0(){
412 //
413 // calculate mini-max effective sigma of point angle resolution
414 //
415 //compv0->fTree->SetAlias("SigmaAP2","max(min((SigmaAP0+SigmaAPE0)*0.5,1.5*SigmaAPE0),0.5*SigmaAPE0+0.003)");
416 Double_t effectiveSigma = GetEffectiveSigmaAP0();
417 Double_t sigmaMMAP = 0.5*(GetSigmaAP0()+effectiveSigma);
418 sigmaMMAP = TMath::Min(sigmaMMAP, fgkParams.fPMaxFractionAP0*effectiveSigma);
419 sigmaMMAP = TMath::Max(sigmaMMAP, fgkParams.fPMinFractionAP0*effectiveSigma+fgkParams.fPMinAP0);
420 return sigmaMMAP;
421}
422Double_t AliESDv0::GetMinimaxSigmaD0(){
423 //
424 // calculate mini-max sigma of dca resolution
425 //
426 //compv0->fTree->SetAlias("SigmaD2","max(min((SigmaD0+SigmaDE0)*0.5,1.5*SigmaDE0),0.5*SigmaDE0)");
427 Double_t effectiveSigma = GetEffectiveSigmaD0();
428 Double_t sigmaMMD0 = 0.5*(GetSigmaD0()+effectiveSigma);
429 sigmaMMD0 = TMath::Min(sigmaMMD0, fgkParams.fPMaxFractionD0*effectiveSigma);
430 sigmaMMD0 = TMath::Max(sigmaMMD0, fgkParams.fPMinFractionD0*effectiveSigma+fgkParams.fPMinD0);
431 return sigmaMMD0;
432}
433
434
435Double_t AliESDv0::GetLikelihoodAP(Int_t mode0, Int_t mode1){
436 //
437 // get likelihood for point angle
438 //
439 Double_t sigmaAP = 0.007; //default sigma
440 switch (mode0){
441 case 0:
442 sigmaAP = GetSigmaAP0(); // mode 0 - covariance matrix estimates used
443 break;
444 case 1:
445 sigmaAP = GetEffectiveSigmaAP0(); // mode 1 - effective sigma used
446 break;
447 case 2:
448 sigmaAP = GetMinimaxSigmaAP0(); // mode 2 - minimax sigma
449 break;
450 }
451 Double_t apNorm = TMath::Min(TMath::ACos(fPointAngle)/sigmaAP,50.);
452 //normalized point angle, restricted - because of overflow problems in Exp
453 Double_t likelihood = 0;
454 switch(mode1){
455 case 0:
456 likelihood = TMath::Exp(-0.5*apNorm*apNorm);
457 // one component
458 break;
459 case 1:
460 likelihood = (TMath::Exp(-0.5*apNorm*apNorm)+0.5* TMath::Exp(-0.25*apNorm*apNorm))/1.5;
461 // two components
462 break;
463 case 2:
464 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;
465 // three components
466 break;
467 }
468 return likelihood;
469}
470
471Double_t AliESDv0::GetLikelihoodD(Int_t mode0, Int_t mode1){
472 //
473 // get likelihood for DCA
474 //
475 Double_t sigmaD = 0.03; //default sigma
476 switch (mode0){
477 case 0:
478 sigmaD = GetSigmaD0(); // mode 0 - covariance matrix estimates used
479 break;
480 case 1:
481 sigmaD = GetEffectiveSigmaD0(); // mode 1 - effective sigma used
482 break;
483 case 2:
484 sigmaD = GetMinimaxSigmaD0(); // mode 2 - minimax sigma
485 break;
486 }
b75d63a7 487
488 //Bo: Double_t dNorm = TMath::Min(fDist2/sigmaD,50.);
489 Double_t dNorm = TMath::Min(fDcaV0Daughters/sigmaD,50.);//Bo:
d6a49f20 490 //normalized point angle, restricted - because of overflow problems in Exp
491 Double_t likelihood = 0;
492 switch(mode1){
493 case 0:
494 likelihood = TMath::Exp(-2.*dNorm);
495 // one component
496 break;
497 case 1:
498 likelihood = (TMath::Exp(-2.*dNorm)+0.5* TMath::Exp(-dNorm))/1.5;
499 // two components
500 break;
501 case 2:
502 likelihood = (TMath::Exp(-2.*dNorm)+0.5* TMath::Exp(-dNorm)+0.25*TMath::Exp(-0.5*dNorm))/1.75;
503 // three components
504 break;
505 }
506 return likelihood;
507
508}
509
510Double_t AliESDv0::GetLikelihoodC(Int_t mode0, Int_t /*mode1*/){
511 //
512 // get likelihood for Causality
513 // !!! Causality variables defined in AliITStrackerMI !!!
514 // when more information was available
515 //
516 Double_t likelihood = 0.5;
517 Double_t minCausal = TMath::Min(fCausality[0],fCausality[1]);
518 Double_t maxCausal = TMath::Max(fCausality[0],fCausality[1]);
519 // minCausal = TMath::Max(minCausal,0.5*maxCausal);
520 //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");
521
522 switch(mode0){
523 case 0:
524 //normalization
525 likelihood = TMath::Power((1.05-2*(0.8-TMath::Exp(-maxCausal))),4.);
526 break;
527 case 1:
528 likelihood = TMath::Power(1.05-(2*(0.8-TMath::Exp(-maxCausal))+(2*(0.8-TMath::Exp(-minCausal))))*0.5,4.);
529 break;
530 }
531 return likelihood;
532
533}
534
535void AliESDv0::SetCausality(Float_t pb0, Float_t pb1, Float_t pa0, Float_t pa1)
536{
537 //
538 // set probabilities
539 //
540 fCausality[0] = pb0; // probability - track 0 exist before vertex
541 fCausality[1] = pb1; // probability - track 1 exist before vertex
542 fCausality[2] = pa0; // probability - track 0 exist close after vertex
543 fCausality[3] = pa1; // probability - track 1 exist close after vertex
544}
545void AliESDv0::SetClusters(Int_t *clp, Int_t *clm)
546{
547 //
548 // Set its clusters indexes
549 //
550 for (Int_t i=0;i<6;i++) fClusters[0][i] = clp[i];
551 for (Int_t i=0;i<6;i++) fClusters[1][i] = clm[i];
552}
9973a4bf 553
554Float_t AliESDv0::GetEffMass(UInt_t p1, UInt_t p2){
555 //
556 // calculate effective mass
557 //
558 const Float_t kpmass[5] = {5.10000000000000037e-04,1.05660000000000004e-01,1.39570000000000000e-01,
559 4.93599999999999983e-01, 9.38270000000000048e-01};
560 if (p1>4) return -1;
561 if (p2>4) return -1;
562 Float_t mass1 = kpmass[p1];
563 Float_t mass2 = kpmass[p2];
564 Double_t *m1 = fPmom;
565 Double_t *m2 = fNmom;
566 //
567 //if (fRP[p1]+fRM[p2]<fRP[p2]+fRM[p1]){
568 // m1 = fPM;
569 // m2 = fPP;
570 //}
571 //
572 Float_t e1 = TMath::Sqrt(mass1*mass1+
573 m1[0]*m1[0]+
574 m1[1]*m1[1]+
575 m1[2]*m1[2]);
576 Float_t e2 = TMath::Sqrt(mass2*mass2+
577 m2[0]*m2[0]+
578 m2[1]*m2[1]+
579 m2[2]*m2[2]);
580 Float_t mass =
581 (m2[0]+m1[0])*(m2[0]+m1[0])+
582 (m2[1]+m1[1])*(m2[1]+m1[1])+
583 (m2[2]+m1[2])*(m2[2]+m1[2]);
584
585 mass = TMath::Sqrt((e1+e2)*(e1+e2)-mass);
586 return mass;
587}