]> git.uio.no Git - u/mrichter/AliRoot.git/blob - HLT/trigger/AliHLTD0toKpi.cxx
D0 trigger update by Gaute
[u/mrichter/AliRoot.git] / HLT / trigger / AliHLTD0toKpi.cxx
1 // $Id$
2 #include "AliHLTD0toKpi.h"
3 #include "TDatabasePDG.h"
4 #include "TMath.h"
5 #include "AliESDtrack.h"
6 #include "TVector3.h"
7 #include "AliAODVertex.h"
8 #include "AliESDVertex.h"
9 #include "TObjArray.h"
10 #include "AliVertexerTracks.h"
11 #include "AliHLTGlobalBarrelTrack.h"
12 #include "AliExternalTrackParam.h"
13
14 ClassImp(AliHLTD0toKpi)
15
16 AliHLTD0toKpi::AliHLTD0toKpi() 
17 {
18 }
19
20 Double_t AliHLTD0toKpi::InvMass(AliExternalTrackParam* d1, AliExternalTrackParam* d2)
21 {
22   Double_t mpi=TDatabasePDG::Instance()->GetParticle(211)->Mass();
23   Double_t mK=TDatabasePDG::Instance()->GetParticle(321)->Mass();
24
25   Double_t energy[2]; 
26   energy[1] = TMath::Sqrt(mK*mK+d1->GetP()*d1->GetP());
27   energy[0] = TMath::Sqrt(mpi*mpi+d2->GetP()*d2->GetP());
28
29   Double_t p1[3],p2[3];
30   d1->GetPxPyPz(p1);
31   d2->GetPxPyPz(p2);
32   
33   Double_t momTot2 = (p1[0]+p2[0])*(p1[0]+p2[0])+
34                      (p1[1]+p2[1])*(p1[1]+p2[1])+
35                      (p1[2]+p2[2])*(p1[2]+p2[2]);
36
37   return TMath::Sqrt((energy[0]+energy[1])*(energy[0]+energy[1])-momTot2);
38
39 }
40 void AliHLTD0toKpi::cosThetaStar(AliExternalTrackParam* d1, AliExternalTrackParam* d2,Double_t &D0,Double_t &D0bar)
41 {
42   Double_t mD0 = TDatabasePDG::Instance()->GetParticle(421)->Mass();
43   Double_t mpi=TDatabasePDG::Instance()->GetParticle(211)->Mass();
44   Double_t mK=TDatabasePDG::Instance()->GetParticle(321)->Mass();
45
46   Double_t pStar = TMath::Sqrt(TMath::Power(mD0*mD0-mK*mK-mpi*mpi,2.)-4.*mK*mK*mpi*mpi)/(2.*mD0);
47  
48   Double_t px = d1->Px()+d2->Px();
49   Double_t py = d1->Py()+d2->Py();
50   Double_t pz = d1->Pz()+d2->Pz();
51   Double_t p = TMath::Sqrt(px*px+py*py+pz*pz);
52   Double_t energy = TMath::Sqrt(p*p+mD0*mD0);
53
54   Double_t beta = p/energy;
55   Double_t gamma = energy/mD0;
56   
57   Double_t qL;
58   TVector3 mom(d1->Px(),d1->Py(),d1->Pz());
59   TVector3 momD(px,py,pz);
60   qL = mom.Dot(momD)/momD.Mag();
61
62   D0 = (qL/gamma-beta*TMath::Sqrt(pStar*pStar+mK*mK))/pStar;
63   
64   TVector3 mom2(d2->Px(),d2->Py(),d2->Pz());
65   TVector3 momD2(px,py,pz);
66   qL = mom2.Dot(momD2)/momD2.Mag();
67
68   D0bar = (qL/gamma-beta*TMath::Sqrt(pStar*pStar+mK*mK))/pStar;
69
70 }
71 Double_t AliHLTD0toKpi::pointingAngle(AliExternalTrackParam* n, AliExternalTrackParam* p, Double_t *pv, Double_t *sv)
72 {
73
74   TVector3 mom(n->Px()+p->Px(),n->Py()+p->Py(),n->Pz()+p->Pz());
75   TVector3 flight(sv[0]-pv[0],sv[1]-pv[1],sv[2]-pv[2]);
76   
77   double pta = mom.Angle(flight);
78
79   return TMath::Cos(pta); 
80 }
81
82 AliAODVertex* AliHLTD0toKpi::ReconstructSecondaryVertex(TObjArray *trkArray, Double_t b, AliESDVertex *v)
83 {
84   
85   AliESDVertex *vertexESD = 0;
86   AliAODVertex *vertexAOD = 0;
87   
88   AliVertexerTracks *vertexer = new AliVertexerTracks(b);
89   vertexer->SetVtxStart(v);
90   //if(isESD){vertexESD = (AliESDVertex*)vertexer->VertexForSelectedESDTracks(trkArray);}
91   UShort_t *id = new UShort_t[2];
92   AliHLTGlobalBarrelTrack *t1 = (AliHLTGlobalBarrelTrack*) trkArray->At(0);
93   AliHLTGlobalBarrelTrack *t2 = (AliHLTGlobalBarrelTrack*) trkArray->At(1);
94   id[0]=(UShort_t) t1->GetID();
95   id[1]=(UShort_t) t2->GetID();
96   vertexESD = (AliESDVertex*)vertexer->VertexForSelectedTracks(trkArray,id);
97   delete id;
98   delete vertexer; vertexer=NULL;
99   
100   if(!vertexESD) return vertexAOD;
101   
102   if(vertexESD->GetNContributors()!=trkArray->GetEntriesFast()) { 
103     //AliDebug(2,"vertexing failed"); 
104     delete vertexESD; vertexESD=NULL;
105     return vertexAOD;
106   }
107
108   // convert to AliAODVertex
109   Double_t pos[3],cov[6],chi2perNDF;
110   vertexESD->GetXYZ(pos); // position
111   vertexESD->GetCovMatrix(cov); //covariance matrix
112   chi2perNDF = vertexESD->GetChi2toNDF();
113   //dispersion = vertexESD->GetDispersion();
114   delete vertexESD; vertexESD=NULL;
115
116   Int_t nprongs= trkArray->GetEntriesFast();
117   vertexAOD = new AliAODVertex(pos,cov,chi2perNDF,0x0,-1,AliAODVertex::kUndef,nprongs);
118
119   return vertexAOD;
120
121 }