In muon-related macros in EVE:
[u/mrichter/AliRoot.git] / EVE / alice-macros / esd_cascade.C
1 // $Id$
2 // Main authors: Matevz Tadel & Alja Mrak-Tadel: 2006, 2007
3
4 /**************************************************************************
5  * Copyright(c) 1998-2008, ALICE Experiment at CERN, all rights reserved. *
6  * See http://aliceinfo.cern.ch/Offline/AliRoot/License.html for          *
7  * full copyright notice.                                                 *
8  **************************************************************************/
9
10 void esd_cascade_init_rectrack(TEveRecTrack& rt, AliExternalTrackParam* tp)
11 {
12   Double_t      pbuf[3], vbuf[3];
13
14   rt.fSign = tp->GetSign();
15   tp->GetXYZ(vbuf);     rt.fV.Set(vbuf);
16   tp->GetPxPyPz(pbuf);  rt.fP.Set(pbuf);
17   // Double_t ep = at->GetP(), mc = at->GetMass();
18   rt.fBeta = 1; // ep/TMath::Sqrt(ep*ep + mc*mc);
19 }
20
21 AliEveCascade* esd_make_cascade(TEveTrackPropagator* rnrStyle, AliESDVertex* primVtx,
22                                 AliESDtrack* bac, AliESDcascade* cascade, Int_t i)
23 {
24   TEveRecTrack   rcPos;
25   TEveRecTrack   rcNeg;
26   TEveRecV0      rcV0;
27
28   TEveRecTrack   rcBac;
29   TEveRecCascade rcCascade;
30   
31   Double_t v[3];
32   cascade->GetXYZ(v[0], v[1], v[2]);
33   rcV0.fVCa.Set(v);
34
35   cascade->GetParamN()->GetXYZ(v);  rcV0.fVNeg.Set(v);
36   cascade->GetParamP()->GetXYZ(v);  rcV0.fVPos.Set(v);
37
38   rcV0.fV0Birth.Set(primVtx->GetXv(), primVtx->GetYv(), primVtx->GetZv());
39
40   Double_t pCascade[3]={0.}, pBac[3]={0.}, pNeg[3]={0.}, pPos[3]={0.}, cv[21]={0.};
41   //cascade->GetPxPyPz(pCascade[0], pCascade[1], pCascade[2]);
42   cascade->GetBPxPyPz(pBac[0], pBac[1], pBac[2]);
43   cascade->GetNPxPyPz(pNeg[0], pNeg[1], pNeg[2]);
44   cascade->GetPPxPyPz(pPos[0], pPos[1], pPos[2]);
45   
46         rcCascade.fPBac.Set(pBac);
47         rcV0.fPNeg.Set(pNeg);
48         rcV0.fPPos.Set(pPos);
49   /*
50   // Debug
51   printf("\n ESD info \n");     
52   printf("Neg :  px = %.5f, py = %.5f, pz = %.5f \n",pNeg[0], pNeg[1], pNeg[2]);
53   printf("Pos :  px = %.5f, py = %.5f, pz = %.5f \n",pPos[0], pPos[1], pPos[2]);
54   printf("Bach : px = %.5f, py = %.5f, pz = %.5f \n",pBac[0], pBac[1], pBac[2]);
55   
56   printf("\n EVE info \n");     
57   printf("Neg :  px = %.5f, py = %.5f, pz = %.5f \n",rcV0.fPNeg.fX, rcV0.fPNeg.fY, rcV0.fPNeg.fZ);
58   printf("Pos :  px = %.5f, py = %.5f, pz = %.5f \n",rcV0.fPPos.fX, rcV0.fPPos.fY, rcV0.fPPos.fZ);
59   printf("Bach : px = %.5f, py = %.5f, pz = %.5f \n",rcCascade.fPBac.fX,rcCascade.fPBac.fY,rcCascade.fPBac.fZ);
60   */
61   Double_t pLambda = TMath::Sqrt((pNeg[0]+pPos[0])* (pNeg[0]+pPos[0]) +
62                                  (pNeg[1]+pPos[1])* (pNeg[1]+pPos[1]) +
63                                  (pNeg[2]+pPos[2])* (pNeg[2]+pPos[2]) );
64   Double_t pBach   = TMath::Sqrt( pBac[0]*pBac[0] + pBac[1]*pBac[1] + pBac[2]*pBac[2]);
65
66   cascade->GetXYZcascade(v[0], v[1], v[2]);
67   rcCascade.fCascadeVCa.Set(v);
68
69   rcCascade.fCascadeBirth.Set(primVtx->GetXv(), primVtx->GetYv(), primVtx->GetZv());
70
71   // Simulation data not directly available in AliESDcascade
72   // rcCascade.fDLabel = cascade->GetBindex();
73
74   // Problem: two following lines are not possible: no GetParamB !!
75   // cascade->GetParamB()->GetXYZ(v);  rcCascade.fVBac.Set(v); 
76   // esd_cascade_init_rectrack(rcBac, cascade->GetParamB());
77   // Solution: create an AliExternalTrackParam with null cv...
78   AliExternalTrackParam *bParam = new AliExternalTrackParam(v,pBac,cv,cascade->Charge());
79   esd_cascade_init_rectrack(rcBac,bParam);
80    rcBac.fIndex = cascade->GetBindex();
81   
82   esd_cascade_init_rectrack(rcNeg, cascade->GetParamN());
83    rcNeg.fIndex = cascade->GetNindex();
84   esd_cascade_init_rectrack(rcPos, cascade->GetParamP());
85    rcPos.fIndex = cascade->GetPindex();
86   
87   
88 /*
89   // Debug
90   TEveVector BacMom = rcBac.GetMomentum();
91   printf("Bac mom : px = %f, py = %f, pz =%f, Ptot = %f \n", BacMom.fX, 
92          rcBac.GetMomentum().fY, 
93          BacMom.fZ, 
94          BacMom.Mag());
95 */
96   AliEveCascade* myCascade = new AliEveCascade(&rcBac, &rcNeg, &rcPos, &rcV0, &rcCascade, rnrStyle);
97   myCascade->SetElementName(Form("ESDcascade %d", i));
98   
99   /*
100   // Debug
101   TEveVector CascMom(rcCascade.fPBac + rcV0.fPNeg + rcV0.fPPos);
102   printf("Casc mom : px = %f, py = %f, pz =%f, Ptot = %f \n", CascMom.fX, CascMom.fY, CascMom.fZ, CascMom.Mag());
103   */
104   myCascade->SetElementTitle(Form("Info coming directly from AliESDcascade : \n - Charge : %d \n - Cascade decay position : x = %.4f, y = %.4f, z = %.4f, Transv. radius = %.4f cm, Decay Length = %.4f cm\n\n - Pt(Cascade) : %f GeV/c, Ptot(Cascade): %f GeV/c\n - Lambda : px = %.4f, py = %.4f, pz = %.4f, Ptot : %f GeV/c\n - Bach   : px = %.4f, py = %.4f, pz = %.4f, Ptot : %f GeV/c\n\n - Eta : %f\n - Phi : %f deg \n - Theta : %f deg\n - DCA : %f cm \n - Cos(Ptg Angle) : %f \n\n - Eff. mass (Xi hyp) : %f GeV/c2", 
105                                         cascade->Charge(),
106                                         v[0], v[1], v[2], TMath::Sqrt(v[0]*v[0] +v[1]*v[1]), TMath::Sqrt(v[0]*v[0] +v[1]*v[1] + v[2]*v[2] ),
107                                         cascade->Pt(), cascade->P(),
108                                         pNeg[0]+pPos[0], pNeg[1]+pPos[1], pNeg[2]+pPos[2], pLambda,
109                                         pBac[0], pBac[1], pBac[2], pBach,
110                                         cascade->Eta(),
111                                         cascade->Phi()   * 180/TMath::Pi(),
112                                         cascade->Theta() * 180/TMath::Pi(),
113                                         cascade->GetDcaXiDaughters(),
114                                         cascade->GetCascadeCosineOfPointingAngle(primVtx->GetXv(),
115                                                                                  primVtx->GetYv(),
116                                                                                  primVtx->GetZv()),
117                                         cascade->GetEffMassXi()
118                                  )
119                             );
120
121   
122   myCascade->SetESDIndex(i);
123   myCascade->SetDaughterDCA(cascade->GetDcaXiDaughters());
124   myCascade->SetLambdaP( pNeg[0]+pPos[0], pNeg[1]+pPos[1], pNeg[2]+pPos[2] );
125   myCascade->SetBachP(   pBac[0], pBac[1], pBac[2]);
126   return myCascade;
127 }
128
129
130 AliEveCascadeList* esd_cascade()
131 {
132   AliESDEvent* esd = AliEveEventManager::AssertESD();
133
134   AliESDVertex* primVertex = (AliESDVertex*) esd->GetPrimaryVertex();
135
136   AliEveCascadeList* cont = new AliEveCascadeList("ESD cascade");
137   cont->SetMainColor(kBlue+2);
138   TEveTrackPropagator* rnrStyle = cont->GetPropagator();
139   rnrStyle->SetMagField( 0.1*esd->GetMagneticField() );
140
141   gEve->AddElement(cont);
142
143   Int_t count = 0;
144   for (Int_t n=0; n<esd->GetNumberOfCascades(); ++n)
145   {
146     AliESDcascade *cascade = esd->GetCascade(n);
147
148     Int_t bacInd = cascade->GetBindex();
149     AliESDtrack* bacTr = esd->GetTrack(bacInd);
150
151     AliEveCascade* myCascade = esd_make_cascade(rnrStyle, primVertex, bacTr, cascade, n);
152     if (myCascade)
153     {
154       gEve->AddElement(myCascade, cont);
155       ++count;
156     }
157   }
158
159   cont->SetTitle("Cascade candidates (reco)");
160
161   cont->MakeCascades();
162   gEve->Redraw3D();
163
164   return cont;
165 }