]> git.uio.no Git - u/mrichter/AliRoot.git/blame - STEER/AliESDv0.cxx
a new macro providing a QA GUI for shifter
[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
732a24fe 132
c7bafca9 133AliESDv0::AliESDv0(const AliExternalTrackParam &t1, Int_t i1,
134 const AliExternalTrackParam &t2, Int_t i2) :
135 TObject(),
562dd0b4 136 fParamN(t1),
137 fParamP(t2),
c7bafca9 138 fEffMass(TDatabasePDG::Instance()->GetParticle(kK0Short)->Mass()),
c028b974 139 fDcaV0Daughters(0),
562dd0b4 140 fChi2V0(0.),
141 fRr(0),
142 fDistSigma(0),
143 fChi2Before(0),
144 fChi2After(0),
145 fPointAngleFi(0),
146 fPointAngleTh(0),
b75d63a7 147 fPointAngle(0),
562dd0b4 148 fPdgCode(kK0Short),
c7bafca9 149 fNidx(i1),
d6a49f20 150 fPidx(i2),
d6a49f20 151 fStatus(0),
d6a49f20 152 fNBefore(0),
d6a49f20 153 fNAfter(0),
562dd0b4 154 fOnFlyStatus(kFALSE)
c7bafca9 155{
156 //--------------------------------------------------------------------
157 // Main constructor (K0s)
158 //--------------------------------------------------------------------
159
160 for (Int_t i=0; i<6; i++) {
161 fPosCov[i]= 0.;
c7bafca9 162 }
163
164 //Trivial estimation of the vertex parameters
b75d63a7 165 Double_t alpha=t1.GetAlpha(), cs=TMath::Cos(alpha), sn=TMath::Sin(alpha);
166 Double_t tmp[3];
167 t1.GetPxPyPz(tmp);
168 Double_t px1=tmp[0], py1=tmp[1], pz1=tmp[2];
169 t1.GetXYZ(tmp);
170 Double_t x1=tmp[0], y1=tmp[1], z1=tmp[2];
37919f0b 171 const Double_t ss=0.0005*0.0005;//a kind of a residual misalignment precision
172 Double_t sx1=sn*sn*t1.GetSigmaY2()+ss, sy1=cs*cs*t1.GetSigmaY2()+ss;
c7bafca9 173
174
b75d63a7 175 alpha=t2.GetAlpha(); cs=TMath::Cos(alpha); sn=TMath::Sin(alpha);
176 t2.GetPxPyPz(tmp);
177 Double_t px2=tmp[0], py2=tmp[1], pz2=tmp[2];
178 t2.GetXYZ(tmp);
179 Double_t x2=tmp[0], y2=tmp[1], z2=tmp[2];
37919f0b 180 Double_t sx2=sn*sn*t2.GetSigmaY2()+ss, sy2=cs*cs*t2.GetSigmaY2()+ss;
c7bafca9 181
182 Double_t sz1=t1.GetSigmaZ2(), sz2=t2.GetSigmaZ2();
37919f0b 183 Double_t wx1=sx2/(sx1+sx2), wx2=1.- wx1;
184 Double_t wy1=sy2/(sy1+sy2), wy2=1.- wy1;
c7bafca9 185 Double_t wz1=sz2/(sz1+sz2), wz2=1.- wz1;
186 fPos[0]=wx1*x1 + wx2*x2; fPos[1]=wy1*y1 + wy2*y2; fPos[2]=wz1*z1 + wz2*z2;
187
188 //fPos[0]=0.5*(x1+x2); fPos[1]=0.5*(y1+y2); fPos[2]=0.5*(z1+z2);
189 fNmom[0]=px1; fNmom[1]=py1; fNmom[2]=pz1;
190 fPmom[0]=px2; fPmom[1]=py2; fPmom[2]=pz2;
191
b75d63a7 192 for (Int_t i=0;i<6;i++){fClusters[0][i]=0; fClusters[1][i]=0;}
193 fNormDCAPrim[0]=fNormDCAPrim[1]=0;
194 for (Int_t i=0;i<3;i++){fAngle[i]=0;}
195 for (Int_t i=0;i<4;i++){fCausality[i]=0;}
c028b974 196}
c7bafca9 197
732a24fe 198AliESDv0& AliESDv0::operator=(const AliESDv0 &v0){
199
200
201 //--------------------------------------------------------------------
202 // The assingment operator
203 //--------------------------------------------------------------------
204
205 if(this==&v0)return *this;
206 TObject::operator=(v0);
207 fParamN = v0.fParamN;
208 fParamP = v0.fParamP;
209 fEffMass = v0.fEffMass;
210 fDcaV0Daughters = v0.fDcaV0Daughters;
211 fChi2V0 = v0.fChi2V0;
212 fRr = v0.fRr;
213 fDistSigma = v0.fDistSigma;
214 fChi2Before = v0.fChi2Before;
215 fChi2After = v0.fChi2After;
216 fPointAngleFi = v0.fPointAngleFi;
217 fPointAngleTh = v0.fPointAngleTh;
218 fPointAngle = v0.fPointAngle;
219 fPdgCode = v0.fPdgCode;
220 fNidx = v0.fNidx;
221 fPidx = v0.fPidx;
222 fStatus = v0.fStatus;
223 fNBefore = v0.fNBefore;
224 fNAfter = v0.fNAfter;
225 fOnFlyStatus = v0.fOnFlyStatus;
226
227 for (int i=0; i<3; i++) {
228 fPos[i] = v0.fPos[i];
229 fNmom[i] = v0.fNmom[i];
230 fPmom[i] = v0.fPmom[i];
231 }
232 for (int i=0; i<6; i++) {
233 fPosCov[i] = v0.fPosCov[i];
234 }
235 for (Int_t i=0; i<2; i++) {
236 fNormDCAPrim[i]=v0.fNormDCAPrim[i];
237 }
238 for (Int_t i=0;i<6;i++){
239 fClusters[0][i]=v0.fClusters[0][i];
240 fClusters[1][i]=v0.fClusters[1][i];
241 }
242 for (Int_t i=0;i<3;i++){
243 fAngle[i]=v0.fAngle[i];
244 }
245 for (Int_t i=0;i<4;i++){fCausality[i]=v0.fCausality[i];}
246
247 return *this;
248}
249
250void AliESDv0::Copy(TObject& obj) const {
251
252 // this overwrites the virtual TOBject::Copy()
253 // to allow run time copying without casting
254 // in AliESDEvent
255
256 if(this==&obj)return;
257 AliESDv0 *robj = dynamic_cast<AliESDv0*>(&obj);
258 if(!robj)return; // not an aliesesv0
259 *robj = *this;
260}
261
262
c028b974 263AliESDv0::~AliESDv0(){
264 //--------------------------------------------------------------------
265 // Empty destructor
266 //--------------------------------------------------------------------
c7bafca9 267}
268
c028b974 269
270
e23730c7 271Double_t AliESDv0::ChangeMassHypothesis(Int_t code) {
272 //--------------------------------------------------------------------
273 // This function changes the mass hypothesis for this V0
274 // and returns the "kinematical quality" of this hypothesis
275 //--------------------------------------------------------------------
b75d63a7 276 static
277 Double_t piMass=TDatabasePDG::Instance()->GetParticle(kPiPlus)->Mass();
278 static
279 Double_t prMass=TDatabasePDG::Instance()->GetParticle(kProton)->Mass();
280 static
281 Double_t k0Mass=TDatabasePDG::Instance()->GetParticle(kK0Short)->Mass();
282 static
283 Double_t l0Mass=TDatabasePDG::Instance()->GetParticle(kLambda0)->Mass();
284
285 Double_t nmass=piMass, pmass=piMass, mass=k0Mass, ps=0.206;
e23730c7 286
287 fPdgCode=code;
288
289 switch (code) {
290 case kLambda0:
b75d63a7 291 nmass=piMass; pmass=prMass; mass=l0Mass; ps=0.101; break;
e23730c7 292 case kLambda0Bar:
b75d63a7 293 pmass=piMass; nmass=prMass; mass=l0Mass; ps=0.101; break;
e23730c7 294 case kK0Short:
295 break;
296 default:
5f7789fc 297 AliError("invalide PDG code ! Assuming K0s...");
e23730c7 298 fPdgCode=kK0Short;
299 break;
300 }
301
302 Double_t pxn=fNmom[0], pyn=fNmom[1], pzn=fNmom[2];
303 Double_t pxp=fPmom[0], pyp=fPmom[1], pzp=fPmom[2];
304
305 Double_t en=TMath::Sqrt(nmass*nmass + pxn*pxn + pyn*pyn + pzn*pzn);
306 Double_t ep=TMath::Sqrt(pmass*pmass + pxp*pxp + pyp*pyp + pzp*pzp);
307 Double_t pxl=pxn+pxp, pyl=pyn+pyp, pzl=pzn+pzp;
308 Double_t pl=TMath::Sqrt(pxl*pxl + pyl*pyl + pzl*pzl);
309
310 fEffMass=TMath::Sqrt((en+ep)*(en+ep)-pl*pl);
311
312 Double_t beta=pl/(en+ep);
313 Double_t pln=(pxn*pxl + pyn*pyl + pzn*pzl)/pl;
314 Double_t plp=(pxp*pxl + pyp*pyl + pzp*pzl)/pl;
315
316 Double_t pt2=pxp*pxp + pyp*pyp + pzp*pzp - plp*plp;
317
318 Double_t a=(plp-pln)/(plp+pln);
319 a -= (pmass*pmass-nmass*nmass)/(mass*mass);
320 a = 0.25*beta*beta*mass*mass*a*a + pt2;
321
322 return (a - ps*ps);
323
324}
325
326void AliESDv0::GetPxPyPz(Double_t &px, Double_t &py, Double_t &pz) const {
327 //--------------------------------------------------------------------
328 // This function returns V0's momentum (global)
329 //--------------------------------------------------------------------
330 px=fNmom[0]+fPmom[0];
331 py=fNmom[1]+fPmom[1];
332 pz=fNmom[2]+fPmom[2];
333}
334
335void AliESDv0::GetXYZ(Double_t &x, Double_t &y, Double_t &z) const {
336 //--------------------------------------------------------------------
337 // This function returns V0's position (global)
338 //--------------------------------------------------------------------
339 x=fPos[0];
340 y=fPos[1];
341 z=fPos[2];
342}
343
b75d63a7 344Float_t AliESDv0::GetD(Double_t x0, Double_t y0, Double_t z0) const {
e23730c7 345 //--------------------------------------------------------------------
346 // This function returns V0's impact parameter
347 //--------------------------------------------------------------------
348 Double_t x=fPos[0],y=fPos[1],z=fPos[2];
349 Double_t px=fNmom[0]+fPmom[0];
350 Double_t py=fNmom[1]+fPmom[1];
351 Double_t pz=fNmom[2]+fPmom[2];
352
353 Double_t dx=(y0-y)*pz - (z0-z)*py;
354 Double_t dy=(x0-x)*pz - (z0-z)*px;
355 Double_t dz=(x0-x)*py - (y0-y)*px;
356 Double_t d=TMath::Sqrt((dx*dx+dy*dy+dz*dz)/(px*px+py*py+pz*pz));
357 return d;
358}
c028b974 359
360
b75d63a7 361Float_t AliESDv0::GetV0CosineOfPointingAngle(Double_t& refPointX, Double_t& refPointY, Double_t& refPointZ) const {
c028b974 362 // calculates the pointing angle of the V0 wrt a reference point
363
364 Double_t momV0[3]; //momentum of the V0
365 GetPxPyPz(momV0[0],momV0[1],momV0[2]);
366
367 Double_t deltaPos[3]; //vector between the reference point and the V0 vertex
368 deltaPos[0] = fPos[0] - refPointX;
369 deltaPos[1] = fPos[1] - refPointY;
370 deltaPos[2] = fPos[2] - refPointZ;
371
372 Double_t momV02 = momV0[0]*momV0[0] + momV0[1]*momV0[1] + momV0[2]*momV0[2];
373 Double_t deltaPos2 = deltaPos[0]*deltaPos[0] + deltaPos[1]*deltaPos[1] + deltaPos[2]*deltaPos[2];
374
375 Double_t cosinePointingAngle = (deltaPos[0]*momV0[0] +
376 deltaPos[1]*momV0[1] +
377 deltaPos[2]*momV0[2] ) /
378 TMath::Sqrt(momV02 * deltaPos2);
379
380 return cosinePointingAngle;
381}
d6a49f20 382
383
384// **** The following functions need to be revised
385
074f017b 386void AliESDv0::GetPosCov(Double_t cov[6]) const {
387
388 for (Int_t i=0; i<6; ++i) cov[i] = fPosCov[i];
389
390}
391
d6a49f20 392Double_t AliESDv0::GetSigmaY(){
393 //
394 // return sigmay in y at vertex position using covariance matrix
395 //
396 const Double_t * cp = fParamP.GetCovariance();
397 const Double_t * cm = fParamN.GetCovariance();
398 Double_t sigmay = cp[0]+cm[0]+ cp[5]*(fParamP.GetX()-fRr)*(fParamP.GetX()-fRr)+ cm[5]*(fParamN.GetX()-fRr)*(fParamN.GetX()-fRr);
399 return (sigmay>0) ? TMath::Sqrt(sigmay):100;
400}
401
402Double_t AliESDv0::GetSigmaZ(){
403 //
404 // return sigmay in y at vertex position using covariance matrix
405 //
406 const Double_t * cp = fParamP.GetCovariance();
407 const Double_t * cm = fParamN.GetCovariance();
408 Double_t sigmaz = cp[2]+cm[2]+ cp[9]*(fParamP.GetX()-fRr)*(fParamP.GetX()-fRr)+ cm[9]*(fParamN.GetX()-fRr)*(fParamN.GetX()-fRr);
409 return (sigmaz>0) ? TMath::Sqrt(sigmaz):100;
410}
411
412Double_t AliESDv0::GetSigmaD0(){
413 //
414 // Sigma parameterization using covariance matrix
415 //
416 // sigma of distance between two tracks in vertex position
417 // sigma of DCA is proportianal to sigmaD0
418 // factor 2 difference is explained by the fact that the DCA is calculated at the position
419 // where the tracks as closest together ( not exact position of the vertex)
420 //
421 const Double_t * cp = fParamP.GetCovariance();
422 const Double_t * cm = fParamN.GetCovariance();
423 Double_t sigmaD0 = cp[0]+cm[0]+cp[2]+cm[2]+fgkParams.fPSigmaOffsetD0*fgkParams.fPSigmaOffsetD0;
424 sigmaD0 += ((fParamP.GetX()-fRr)*(fParamP.GetX()-fRr))*(cp[5]+cp[9]);
425 sigmaD0 += ((fParamN.GetX()-fRr)*(fParamN.GetX()-fRr))*(cm[5]+cm[9]);
426 return (sigmaD0>0)? TMath::Sqrt(sigmaD0):100;
427}
428
429
430Double_t AliESDv0::GetSigmaAP0(){
431 //
432 //Sigma parameterization using covariance matrices
433 //
b75d63a7 434 Double_t prec = TMath::Sqrt((fNmom[0]+fPmom[0])*(fNmom[0]+fPmom[0])
435 +(fNmom[1]+fPmom[1])*(fNmom[1]+fPmom[1])
436 +(fNmom[2]+fPmom[2])*(fNmom[2]+fPmom[2]));
437 Double_t normp = TMath::Sqrt(fPmom[0]*fPmom[0]+fPmom[1]*fPmom[1]+fPmom[2]*fPmom[2])/prec; // fraction of the momenta
438 Double_t normm = TMath::Sqrt(fNmom[0]*fNmom[0]+fNmom[1]*fNmom[1]+fNmom[2]*fNmom[2])/prec;
d6a49f20 439 const Double_t * cp = fParamP.GetCovariance();
440 const Double_t * cm = fParamN.GetCovariance();
441 Double_t sigmaAP0 = fgkParams.fPSigmaOffsetAP0*fgkParams.fPSigmaOffsetAP0; // minimal part
442 sigmaAP0 += (cp[5]+cp[9])*(normp*normp)+(cm[5]+cm[9])*(normm*normm); // angular resolution part
443 Double_t sigmaAP1 = GetSigmaD0()/(TMath::Abs(fRr)+0.01); // vertex position part
444 sigmaAP0 += 0.5*sigmaAP1*sigmaAP1;
445 return (sigmaAP0>0)? TMath::Sqrt(sigmaAP0):100;
446}
447
448Double_t AliESDv0::GetEffectiveSigmaD0(){
449 //
450 // minimax - effective Sigma parameterization
451 // p12 effective curvature and v0 radius postion used as parameters
452 //
453 Double_t p12 = TMath::Sqrt(fParamP.GetParameter()[4]*fParamP.GetParameter()[4]+
454 fParamN.GetParameter()[4]*fParamN.GetParameter()[4]);
455 Double_t sigmaED0= TMath::Max(TMath::Sqrt(fRr)-fgkParams.fPSigmaRminDE,0.0)*fgkParams.fPSigmaCoefDE*p12*p12;
456 sigmaED0*= sigmaED0;
457 sigmaED0*= sigmaED0;
458 sigmaED0 = TMath::Sqrt(sigmaED0+fgkParams.fPSigmaOffsetDE*fgkParams.fPSigmaOffsetDE);
459 return (sigmaED0<fgkParams.fPSigmaMaxDE) ? sigmaED0: fgkParams.fPSigmaMaxDE;
460}
461
462
463Double_t AliESDv0::GetEffectiveSigmaAP0(){
464 //
465 // effective Sigma parameterization of point angle resolution
466 //
467 Double_t p12 = TMath::Sqrt(fParamP.GetParameter()[4]*fParamP.GetParameter()[4]+
468 fParamN.GetParameter()[4]*fParamN.GetParameter()[4]);
469 Double_t sigmaAPE= fgkParams.fPSigmaBase0APE;
470 sigmaAPE+= fgkParams.fPSigmaR0APE/(fgkParams.fPSigmaR1APE+fRr);
471 sigmaAPE*= (fgkParams.fPSigmaP0APE+fgkParams.fPSigmaP1APE*p12);
472 sigmaAPE = TMath::Min(sigmaAPE,fgkParams.fPSigmaMaxAPE);
473 return sigmaAPE;
474}
475
476
477Double_t AliESDv0::GetMinimaxSigmaAP0(){
478 //
479 // calculate mini-max effective sigma of point angle resolution
480 //
481 //compv0->fTree->SetAlias("SigmaAP2","max(min((SigmaAP0+SigmaAPE0)*0.5,1.5*SigmaAPE0),0.5*SigmaAPE0+0.003)");
482 Double_t effectiveSigma = GetEffectiveSigmaAP0();
483 Double_t sigmaMMAP = 0.5*(GetSigmaAP0()+effectiveSigma);
484 sigmaMMAP = TMath::Min(sigmaMMAP, fgkParams.fPMaxFractionAP0*effectiveSigma);
485 sigmaMMAP = TMath::Max(sigmaMMAP, fgkParams.fPMinFractionAP0*effectiveSigma+fgkParams.fPMinAP0);
486 return sigmaMMAP;
487}
488Double_t AliESDv0::GetMinimaxSigmaD0(){
489 //
490 // calculate mini-max sigma of dca resolution
491 //
492 //compv0->fTree->SetAlias("SigmaD2","max(min((SigmaD0+SigmaDE0)*0.5,1.5*SigmaDE0),0.5*SigmaDE0)");
493 Double_t effectiveSigma = GetEffectiveSigmaD0();
494 Double_t sigmaMMD0 = 0.5*(GetSigmaD0()+effectiveSigma);
495 sigmaMMD0 = TMath::Min(sigmaMMD0, fgkParams.fPMaxFractionD0*effectiveSigma);
496 sigmaMMD0 = TMath::Max(sigmaMMD0, fgkParams.fPMinFractionD0*effectiveSigma+fgkParams.fPMinD0);
497 return sigmaMMD0;
498}
499
500
501Double_t AliESDv0::GetLikelihoodAP(Int_t mode0, Int_t mode1){
502 //
503 // get likelihood for point angle
504 //
505 Double_t sigmaAP = 0.007; //default sigma
506 switch (mode0){
507 case 0:
508 sigmaAP = GetSigmaAP0(); // mode 0 - covariance matrix estimates used
509 break;
510 case 1:
511 sigmaAP = GetEffectiveSigmaAP0(); // mode 1 - effective sigma used
512 break;
513 case 2:
514 sigmaAP = GetMinimaxSigmaAP0(); // mode 2 - minimax sigma
515 break;
516 }
517 Double_t apNorm = TMath::Min(TMath::ACos(fPointAngle)/sigmaAP,50.);
518 //normalized point angle, restricted - because of overflow problems in Exp
519 Double_t likelihood = 0;
520 switch(mode1){
521 case 0:
522 likelihood = TMath::Exp(-0.5*apNorm*apNorm);
523 // one component
524 break;
525 case 1:
526 likelihood = (TMath::Exp(-0.5*apNorm*apNorm)+0.5* TMath::Exp(-0.25*apNorm*apNorm))/1.5;
527 // two components
528 break;
529 case 2:
530 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;
531 // three components
532 break;
533 }
534 return likelihood;
535}
536
537Double_t AliESDv0::GetLikelihoodD(Int_t mode0, Int_t mode1){
538 //
539 // get likelihood for DCA
540 //
541 Double_t sigmaD = 0.03; //default sigma
542 switch (mode0){
543 case 0:
544 sigmaD = GetSigmaD0(); // mode 0 - covariance matrix estimates used
545 break;
546 case 1:
547 sigmaD = GetEffectiveSigmaD0(); // mode 1 - effective sigma used
548 break;
549 case 2:
550 sigmaD = GetMinimaxSigmaD0(); // mode 2 - minimax sigma
551 break;
552 }
b75d63a7 553
554 //Bo: Double_t dNorm = TMath::Min(fDist2/sigmaD,50.);
555 Double_t dNorm = TMath::Min(fDcaV0Daughters/sigmaD,50.);//Bo:
d6a49f20 556 //normalized point angle, restricted - because of overflow problems in Exp
557 Double_t likelihood = 0;
558 switch(mode1){
559 case 0:
560 likelihood = TMath::Exp(-2.*dNorm);
561 // one component
562 break;
563 case 1:
564 likelihood = (TMath::Exp(-2.*dNorm)+0.5* TMath::Exp(-dNorm))/1.5;
565 // two components
566 break;
567 case 2:
568 likelihood = (TMath::Exp(-2.*dNorm)+0.5* TMath::Exp(-dNorm)+0.25*TMath::Exp(-0.5*dNorm))/1.75;
569 // three components
570 break;
571 }
572 return likelihood;
573
574}
575
576Double_t AliESDv0::GetLikelihoodC(Int_t mode0, Int_t /*mode1*/){
577 //
578 // get likelihood for Causality
579 // !!! Causality variables defined in AliITStrackerMI !!!
580 // when more information was available
581 //
582 Double_t likelihood = 0.5;
583 Double_t minCausal = TMath::Min(fCausality[0],fCausality[1]);
584 Double_t maxCausal = TMath::Max(fCausality[0],fCausality[1]);
585 // minCausal = TMath::Max(minCausal,0.5*maxCausal);
586 //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");
587
588 switch(mode0){
589 case 0:
590 //normalization
591 likelihood = TMath::Power((1.05-2*(0.8-TMath::Exp(-maxCausal))),4.);
592 break;
593 case 1:
594 likelihood = TMath::Power(1.05-(2*(0.8-TMath::Exp(-maxCausal))+(2*(0.8-TMath::Exp(-minCausal))))*0.5,4.);
595 break;
596 }
597 return likelihood;
598
599}
600
601void AliESDv0::SetCausality(Float_t pb0, Float_t pb1, Float_t pa0, Float_t pa1)
602{
603 //
604 // set probabilities
605 //
606 fCausality[0] = pb0; // probability - track 0 exist before vertex
607 fCausality[1] = pb1; // probability - track 1 exist before vertex
608 fCausality[2] = pa0; // probability - track 0 exist close after vertex
609 fCausality[3] = pa1; // probability - track 1 exist close after vertex
610}
611void AliESDv0::SetClusters(Int_t *clp, Int_t *clm)
612{
613 //
614 // Set its clusters indexes
615 //
616 for (Int_t i=0;i<6;i++) fClusters[0][i] = clp[i];
617 for (Int_t i=0;i<6;i++) fClusters[1][i] = clm[i];
618}
9973a4bf 619
620Float_t AliESDv0::GetEffMass(UInt_t p1, UInt_t p2){
621 //
622 // calculate effective mass
623 //
624 const Float_t kpmass[5] = {5.10000000000000037e-04,1.05660000000000004e-01,1.39570000000000000e-01,
625 4.93599999999999983e-01, 9.38270000000000048e-01};
626 if (p1>4) return -1;
627 if (p2>4) return -1;
628 Float_t mass1 = kpmass[p1];
629 Float_t mass2 = kpmass[p2];
630 Double_t *m1 = fPmom;
631 Double_t *m2 = fNmom;
632 //
633 //if (fRP[p1]+fRM[p2]<fRP[p2]+fRM[p1]){
634 // m1 = fPM;
635 // m2 = fPP;
636 //}
637 //
638 Float_t e1 = TMath::Sqrt(mass1*mass1+
639 m1[0]*m1[0]+
640 m1[1]*m1[1]+
641 m1[2]*m1[2]);
642 Float_t e2 = TMath::Sqrt(mass2*mass2+
643 m2[0]*m2[0]+
644 m2[1]*m2[1]+
645 m2[2]*m2[2]);
646 Float_t mass =
647 (m2[0]+m1[0])*(m2[0]+m1[0])+
648 (m2[1]+m1[1])*(m2[1]+m1[1])+
649 (m2[2]+m1[2])*(m2[2]+m1[2]);
650
651 mass = TMath::Sqrt((e1+e2)*(e1+e2)-mass);
652 return mass;
653}