]>
Commit | Line | Data |
---|---|---|
5f4502cc | 1 | // $Id$ |
2 | #include "AliHLTD0toKpi.h" | |
3 | #include "TDatabasePDG.h" | |
4 | #include "TMath.h" | |
5 | #include "AliESDtrack.h" | |
6 | #include "TVector3.h" | |
9929f8f5 | 7 | #include "AliAODVertex.h" |
8 | #include "AliESDVertex.h" | |
9 | #include "TObjArray.h" | |
10 | #include "AliVertexerTracks.h" | |
629b904b | 11 | #include "AliHLTGlobalBarrelTrack.h" |
12 | #include "AliExternalTrackParam.h" | |
5f4502cc | 13 | |
14 | ClassImp(AliHLTD0toKpi) | |
15 | ||
16 | AliHLTD0toKpi::AliHLTD0toKpi() | |
17 | { | |
18 | } | |
19 | ||
629b904b | 20 | Double_t AliHLTD0toKpi::InvMass(AliExternalTrackParam* d1, AliExternalTrackParam* d2) |
5f4502cc | 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 | } | |
629b904b | 40 | void AliHLTD0toKpi::cosThetaStar(AliExternalTrackParam* d1, AliExternalTrackParam* d2,Double_t &D0,Double_t &D0bar) |
5f4502cc | 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 | } | |
629b904b | 71 | Double_t AliHLTD0toKpi::pointingAngle(AliExternalTrackParam* n, AliExternalTrackParam* p, Double_t *pv, Double_t *sv) |
5f4502cc | 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 | } | |
9929f8f5 | 81 | |
629b904b | 82 | AliAODVertex* AliHLTD0toKpi::ReconstructSecondaryVertex(TObjArray *trkArray, Double_t b, AliESDVertex *v) |
9929f8f5 | 83 | { |
84 | ||
85 | AliESDVertex *vertexESD = 0; | |
86 | AliAODVertex *vertexAOD = 0; | |
87 | ||
88 | AliVertexerTracks *vertexer = new AliVertexerTracks(b); | |
629b904b | 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; | |
9929f8f5 | 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 | } |