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 #if !defined(__CINT__) || defined(__MAKECINT__)
  10 #include <TMath.h>
  11 #include <TEveManager.h>
  12 #include <TEveVSDStructs.h>
  13 #include <TEveTrackPropagator.h>
  14
  15 #include <AliExternalTrackParam.h>
  16 #include <AliESDEvent.h>
  17 #include <AliESDcascade.h>
  18 #include <AliESDVertex.h>
  19 #include <AliEveCascade.h>
  20 #include <AliEveEventManager.h>
  21 #endif
  22
  23 void esd_cascade_init_rectrack(TEveRecTrack& rt, const AliExternalTrackParam* tp)
  24 {
  25   Double_t      pbuf[3], vbuf[3];
  26
  27   rt.fSign = tp->GetSign();
  28   tp->GetXYZ(vbuf);     rt.fV.Set(vbuf);
  29   tp->GetPxPyPz(pbuf);  rt.fP.Set(pbuf);
  30   // Double_t ep = at->GetP(), mc = at->GetMass();
  31   rt.fBeta = 1; // ep/TMath::Sqrt(ep*ep + mc*mc);
  32 }
  33
  34 AliEveCascade* esd_make_cascade(TEveTrackPropagator* rnrStyle, AliESDVertex* primVtx,
  35                                 AliESDtrack* bac, AliESDcascade* cascade, Int_t i)
  36 {
  37   TEveRecTrack   rcPos;
  38   TEveRecTrack   rcNeg;
  39   TEveRecV0      rcV0;
  40
  41   TEveRecTrack   rcBac;
  42   TEveRecCascade rcCascade;
  43
  44   Double_t v[3];
  45   cascade->GetXYZ(v[0], v[1], v[2]);
  46   rcV0.fVCa.Set(v);
  47
  48   cascade->GetParamN()->GetXYZ(v);  rcV0.fVNeg.Set(v);
  49   cascade->GetParamP()->GetXYZ(v);  rcV0.fVPos.Set(v);
  50
  51   rcV0.fV0Birth.Set(primVtx->GetX(), primVtx->GetY(), primVtx->GetZ());
  52
  53   Double_t pCascade[3]={0.}, pBac[3]={0.}, pNeg[3]={0.}, pPos[3]={0.}, cv[21]={0.};
  54   //cascade->GetPxPyPz(pCascade[0], pCascade[1], pCascade[2]);
  55   cascade->GetBPxPyPz(pBac[0], pBac[1], pBac[2]);
  56   cascade->GetNPxPyPz(pNeg[0], pNeg[1], pNeg[2]);
  57   cascade->GetPPxPyPz(pPos[0], pPos[1], pPos[2]);
  58
  59         rcCascade.fPBac.Set(pBac);
  60         rcV0.fPNeg.Set(pNeg);
  61         rcV0.fPPos.Set(pPos);
  62   /*
  63   // Debug
  64   printf("\n ESD info \n");
  65   printf("Neg :  px = %.5f, py = %.5f, pz = %.5f \n",pNeg[0], pNeg[1], pNeg[2]);
  66   printf("Pos :  px = %.5f, py = %.5f, pz = %.5f \n",pPos[0], pPos[1], pPos[2]);
  67   printf("Bach : px = %.5f, py = %.5f, pz = %.5f \n",pBac[0], pBac[1], pBac[2]);
  68
  69   printf("\n EVE info \n");
  70   printf("Neg :  px = %.5f, py = %.5f, pz = %.5f \n",rcV0.fPNeg.fX, rcV0.fPNeg.fY, rcV0.fPNeg.fZ);
  71   printf("Pos :  px = %.5f, py = %.5f, pz = %.5f \n",rcV0.fPPos.fX, rcV0.fPPos.fY, rcV0.fPPos.fZ);
  72   printf("Bach : px = %.5f, py = %.5f, pz = %.5f \n",rcCascade.fPBac.fX,rcCascade.fPBac.fY,rcCascade.fPBac.fZ);
  73   */
  74   Double_t pLambda = TMath::Sqrt((pNeg[0]+pPos[0])* (pNeg[0]+pPos[0]) +
  75                                  (pNeg[1]+pPos[1])* (pNeg[1]+pPos[1]) +
  76                                  (pNeg[2]+pPos[2])* (pNeg[2]+pPos[2]) );
  77   Double_t pBach   = TMath::Sqrt( pBac[0]*pBac[0] + pBac[1]*pBac[1] + pBac[2]*pBac[2]);
  78
  79   cascade->GetXYZcascade(v[0], v[1], v[2]);
  80   rcCascade.fCascadeVCa.Set(v);
  81
  82   rcCascade.fCascadeBirth.Set(primVtx->GetX(), primVtx->GetY(), primVtx->GetZ());
  83
  84   // Simulation data not directly available in AliESDcascade
  85   // rcCascade.fDLabel = cascade->GetBindex();
  86
  87   // Problem: two following lines are not possible: no GetParamB !!
  88   // cascade->GetParamB()->GetXYZ(v);  rcCascade.fVBac.Set(v);
  89   // esd_cascade_init_rectrack(rcBac, cascade->GetParamB());
  90   // Solution: create an AliExternalTrackParam with null cv...
  91   AliExternalTrackParam *bParam = new AliExternalTrackParam(v,pBac,cv,cascade->Charge());
  92   esd_cascade_init_rectrack(rcBac,bParam);
  93    rcBac.fIndex = cascade->GetBindex();
  94
  95   esd_cascade_init_rectrack(rcNeg, cascade->GetParamN());
  96    rcNeg.fIndex = cascade->GetNindex();
  97   esd_cascade_init_rectrack(rcPos, cascade->GetParamP());
  98    rcPos.fIndex = cascade->GetPindex();
  99
 100
 101 /*
 102   // Debug
 103   TEveVector BacMom = rcBac.GetMomentum();
 104   printf("Bac mom : px = %f, py = %f, pz =%f, Ptot = %f \n", BacMom.fX,
 105          rcBac.GetMomentum().fY,
 106          BacMom.fZ,
 107          BacMom.Mag());
 108 */
 109   AliEveCascade* myCascade = new AliEveCascade(&rcBac, &rcNeg, &rcPos, &rcV0, &rcCascade, rnrStyle);
 110   myCascade->SetElementName(Form("ESDcascade %d", i));
 111
 112   /*
 113   // Debug
 114   TEveVector CascMom(rcCascade.fPBac + rcV0.fPNeg + rcV0.fPPos);
 115   printf("Casc mom : px = %f, py = %f, pz =%f, Ptot = %f \n", CascMom.fX, CascMom.fY, CascMom.fZ, CascMom.Mag());
 116   */
 117   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",
 118                                         cascade->Charge(),
 119                                         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] ),
 120                                         cascade->Pt(), cascade->P(),
 121                                         pNeg[0]+pPos[0], pNeg[1]+pPos[1], pNeg[2]+pPos[2], pLambda,
 122                                         pBac[0], pBac[1], pBac[2], pBach,
 123                                         cascade->Eta(),
 124                                         cascade->Phi()   * 180/TMath::Pi(),
 125                                         cascade->Theta() * 180/TMath::Pi(),
 126                                         cascade->GetDcaXiDaughters(),
 127                                         cascade->GetCascadeCosineOfPointingAngle(primVtx->GetX(),
 128                                                                                  primVtx->GetY(),
 129                                                                                  primVtx->GetZ()),
 130                                         cascade->GetEffMassXi()
 131                                  )
 132                             );
 133
 134
 135   myCascade->SetESDIndex(i);
 136   myCascade->SetDaughterDCA(cascade->GetDcaXiDaughters());
 137   myCascade->SetLambdaP( pNeg[0]+pPos[0], pNeg[1]+pPos[1], pNeg[2]+pPos[2] );
 138   myCascade->SetBachP(   pBac[0], pBac[1], pBac[2]);
 139   return myCascade;
 140 }
 141
 142
 143 AliEveCascadeList* esd_cascade()
 144 {
 145   AliESDEvent* esd = AliEveEventManager::AssertESD();
 146
 147   AliESDVertex* primVertex = (AliESDVertex*) esd->GetPrimaryVertex();
 148
 149   AliEveCascadeList* cont = new AliEveCascadeList("ESD cascade");
 150   cont->SetMainColor(kBlue+2);
 151   TEveTrackPropagator* rnrStyle = cont->GetPropagator();
 152   rnrStyle->SetMagField( 0.1*esd->GetMagneticField() );
 153
 154   gEve->AddElement(cont);
 155
 156   Int_t count = 0;
 157   for (Int_t n=0; n<esd->GetNumberOfCascades(); ++n)
 158   {
 159     AliESDcascade *cascade = esd->GetCascade(n);
 160
 161     Int_t bacInd = cascade->GetBindex();
 162     AliESDtrack* bacTr = esd->GetTrack(bacInd);
 163
 164     AliEveCascade* myCascade = esd_make_cascade(rnrStyle, primVertex, bacTr, cascade, n);
 165     if (myCascade)
 166     {
 167       gEve->AddElement(myCascade, cont);
 168       ++count;
 169     }
 170   }
 171
 172   cont->SetTitle("Cascade candidates (reco)");
 173
 174   cont->MakeCascades();
 175   gEve->Redraw3D();
 176
 177   return cont;
 178 }