ANALYSIS/AliAnalysisTaskESDfilter.cxx

   1 /**************************************************************************\r
   2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *\r
   3  *                                                                        *\r
   4  * Author: The ALICE Off-line Project.                                    *\r
   5  * Contributors are mentioned in the code where appropriate.              *\r
   6  *                                                                        *\r
   7  * Permission to use, copy, modify and distribute this software and its   *\r
   8  * documentation strictly for non-commercial purposes is hereby granted   *\r
   9  * without fee, provided that the above copyright notice appears in all   *\r
  10  * copies and that both the copyright notice and this permission notice   *\r
  11  * appear in the supporting documentation. The authors make no claims     *\r
  12  * about the suitability of this software for any purpose. It is          *\r
  13  * provided "as is" without express or implied warranty.                  *\r
  14  **************************************************************************/\r
  15 \r
  16 /* $Id$ */\r
  17  \r
  18 #include <TChain.h>\r
  19 #include <TTree.h>\r
  20 #include <TList.h>\r
  21 #include <TArrayI.h>\r
  22 #include <TRandom.h>\r
  23 #include <TParticle.h>\r
  24 #include <TFile.h>\r
  25 \r
  26 #include "AliAnalysisTaskESDfilter.h"\r
  27 #include "AliAnalysisManager.h"\r
  28 #include "AliESDEvent.h"\r
  29 #include "AliESDRun.h"\r
  30 #include "AliStack.h"\r
  31 #include "AliAODEvent.h"\r
  32 #include "AliMCEvent.h"\r
  33 #include "AliMCEventHandler.h"\r
  34 #include "AliESDInputHandler.h"\r
  35 #include "AliAODHandler.h"\r
  36 #include "AliAODMCParticle.h"\r
  37 #include "AliAnalysisFilter.h"\r
  38 #include "AliESDMuonTrack.h"\r
  39 #include "AliESDVertex.h"\r
  40 #include "AliCentrality.h"\r
  41 #include "AliEventplane.h"\r
  42 #include "AliESDv0.h"\r
  43 #include "AliESDkink.h"\r
  44 #include "AliESDcascade.h"\r
  45 #include "AliESDPmdTrack.h"\r
  46 #include "AliESDCaloCluster.h"\r
  47 #include "AliESDCaloCells.h"\r
  48 #include "AliMultiplicity.h"\r
  49 #include "AliLog.h"\r
  50 #include "AliCodeTimer.h"\r
  51 #include "AliESDtrackCuts.h"\r
  52 #include "AliESDpid.h"\r
  53 #include "Riostream.h"\r
  54 \r
  55 ClassImp(AliAnalysisTaskESDfilter)\r
  56 \r
  57 ////////////////////////////////////////////////////////////////////////\r
  58 \r
  59 AliAnalysisTaskESDfilter::AliAnalysisTaskESDfilter():\r
  60     AliAnalysisTaskSE(),\r
  61     fTrackFilter(0x0),\r
  62     fKinkFilter(0x0),\r
  63     fV0Filter(0x0),\r
  64     fCascadeFilter(0x0),\r
  65     fHighPthreshold(0),\r
  66     fPtshape(0x0),\r
  67     fEnableFillAOD(kTRUE),\r
  68 fUsedTrack(0x0),\r
  69 fUsedKink(0x0),\r
  70 fUsedV0(0x0),\r
  71 fAODTrackRefs(0x0),\r
  72 fAODV0VtxRefs(0x0),\r
  73 fAODV0Refs(0x0),\r
  74 fMChandler(0x0),\r
  75 fNumberOfTracks(0),\r
  76 fNumberOfPositiveTracks(0),\r
  77 fNumberOfV0s(0),\r
  78 fNumberOfVertices(0),\r
  79 fNumberOfCascades(0),\r
  80 fNumberOfKinks(0),\r
  81 fOldESDformat(kFALSE),\r
  82 fPrimaryVertex(0x0),\r
  83 fTPCOnlyFilterMask(0),\r
  84 fIsVZEROEnabled(kTRUE),\r
  85 fAreCascadesEnabled(kTRUE),\r
  86 fAreV0sEnabled(kTRUE),\r
  87 fAreKinksEnabled(kTRUE),\r
  88 fAreTracksEnabled(kTRUE),\r
  89 fArePmdClustersEnabled(kTRUE),\r
  90 fAreCaloClustersEnabled(kTRUE),\r
  91 fAreEMCALCellsEnabled(kTRUE),\r
  92 fArePHOSCellsEnabled(kTRUE),\r
  93 fAreTrackletsEnabled(kTRUE),\r
  94 fESDpid(0x0),\r
  95 fIsPidOwner(kFALSE),\r
  96 fTimeZeroType(AliESDpid::kTOF_T0)\r
  97 {\r
  98   // Default constructor\r
  99 }\r
 100 \r
 101 //______________________________________________________________________________\r
 102 AliAnalysisTaskESDfilter::AliAnalysisTaskESDfilter(const char* name):\r
 103     AliAnalysisTaskSE(name),\r
 104     fTrackFilter(0x0),\r
 105     fKinkFilter(0x0),\r
 106     fV0Filter(0x0),\r
 107     fCascadeFilter(0x0),\r
 108     fHighPthreshold(0),\r
 109     fPtshape(0x0),\r
 110     fEnableFillAOD(kTRUE),\r
 111 fUsedTrack(0x0),\r
 112 fUsedKink(0x0),\r
 113 fUsedV0(0x0),\r
 114 fAODTrackRefs(0x0),\r
 115 fAODV0VtxRefs(0x0),\r
 116 fAODV0Refs(0x0),\r
 117 fMChandler(0x0),\r
 118 fNumberOfTracks(0),\r
 119 fNumberOfPositiveTracks(0),\r
 120 fNumberOfV0s(0),\r
 121 fNumberOfVertices(0),\r
 122 fNumberOfCascades(0),\r
 123 fNumberOfKinks(0),\r
 124 fOldESDformat(kFALSE),\r
 125 fPrimaryVertex(0x0),\r
 126 fTPCOnlyFilterMask(0),\r
 127 fIsVZEROEnabled(kTRUE),\r
 128 fAreCascadesEnabled(kTRUE),\r
 129 fAreV0sEnabled(kTRUE),\r
 130 fAreKinksEnabled(kTRUE),\r
 131 fAreTracksEnabled(kTRUE),\r
 132 fArePmdClustersEnabled(kTRUE),\r
 133 fAreCaloClustersEnabled(kTRUE),\r
 134 fAreEMCALCellsEnabled(kTRUE),\r
 135 fArePHOSCellsEnabled(kTRUE),\r
 136 fAreTrackletsEnabled(kTRUE),\r
 137 fESDpid(0x0),\r
 138 fIsPidOwner(kFALSE),\r
 139 fTimeZeroType(AliESDpid::kTOF_T0)\r
 140 {\r
 141   // Constructor\r
 142 }\r
 143 AliAnalysisTaskESDfilter::~AliAnalysisTaskESDfilter(){\r
 144   if(fIsPidOwner)delete fESDpid;\r
 145 }\r
 146 //______________________________________________________________________________\r
 147 void AliAnalysisTaskESDfilter::UserCreateOutputObjects()\r
 148 {\r
 149   //\r
 150   // Create Output Objects conenct filter to outputtree\r
 151   // \r
 152   if(OutputTree())\r
 153   {\r
 154     OutputTree()->GetUserInfo()->Add(fTrackFilter);\r
 155   }\r
 156   else\r
 157   {\r
 158     AliError("No OutputTree() for adding the track filter");\r
 159   }\r
 160 }\r
 161 \r
 162 //______________________________________________________________________________\r
 163 void AliAnalysisTaskESDfilter::Init()\r
 164 {\r
 165   // Initialization\r
 166   if (fDebug > 1) AliInfo("Init() \n");\r
 167   // Call configuration file\r
 168 }\r
 169 \r
 170 //______________________________________________________________________________\r
 171 void AliAnalysisTaskESDfilter::PrintTask(Option_t *option, Int_t indent) const\r
 172 {\r
 173   AliInfo("");\r
 174   \r
 175   AliAnalysisTaskSE::PrintTask(option,indent);\r
 176   \r
 177   TString spaces(' ',indent+3);\r
 178   \r
 179   cout << spaces.Data() << Form("Cascades     are %s",fAreCascadesEnabled ? "ENABLED":"DISABLED") << endl;\r
 180   cout << spaces.Data() << Form("V0s          are %s",fAreV0sEnabled ? "ENABLED":"DISABLED") << endl;\r
 181   cout << spaces.Data() << Form("Kinks        are %s",fAreKinksEnabled ? "ENABLED":"DISABLED") << endl;\r
 182   cout << spaces.Data() << Form("Tracks       are %s",fAreTracksEnabled ? "ENABLED":"DISABLED") << endl;\r
 183   cout << spaces.Data() << Form("PmdClusters  are %s",fArePmdClustersEnabled ? "ENABLED":"DISABLED") << endl;\r
 184   cout << spaces.Data() << Form("CaloClusters are %s",fAreCaloClustersEnabled ? "ENABLED":"DISABLED") << endl;\r
 185   cout << spaces.Data() << Form("EMCAL cells  are %s",fAreEMCALCellsEnabled ? "ENABLED":"DISABLED") << endl;\r
 186   cout << spaces.Data() << Form("Tracklets    are %s",fAreTrackletsEnabled ? "ENABLED":"DISABLED") << endl;  \r
 187 }\r
 188 \r
 189 //______________________________________________________________________________\r
 190 void AliAnalysisTaskESDfilter::UserExec(Option_t */*option*/)\r
 191 {\r
 192 // Execute analysis for current event\r
 193 //\r
 194                                             \r
 195   Long64_t ientry = Entry();\r
 196   \r
 197   if (fDebug > 0) {\r
 198       printf("Filter: Analysing event # %5d\n", (Int_t) ientry);\r
 199       if (fHighPthreshold == 0) AliInfo("detector PID signals are stored in each track");\r
 200       if (!fPtshape) AliInfo("detector PID signals are not stored below the pt threshold");\r
 201   }\r
 202   // Filters must explicitely enable AOD filling in their UserExec (AG)\r
 203   if (!AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler()) AliFatal("Cannot run ESD filter without an output event handler");\r
 204   if(fEnableFillAOD) {\r
 205      AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler()->SetFillAOD(kTRUE);\r
 206      AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler()->SetFillExtension(kTRUE);\r
 207   }   \r
 208   ConvertESDtoAOD();\r
 209 }\r
 210 \r
 211 //______________________________________________________________________________\r
 212 TClonesArray& AliAnalysisTaskESDfilter::Cascades()\r
 213 {\r
 214   return *(AODEvent()->GetCascades());\r
 215 }\r
 216 \r
 217 //______________________________________________________________________________\r
 218 TClonesArray& AliAnalysisTaskESDfilter::Tracks()\r
 219 {\r
 220   return *(AODEvent()->GetTracks());\r
 221 }\r
 222 \r
 223 //______________________________________________________________________________\r
 224 TClonesArray& AliAnalysisTaskESDfilter::V0s()\r
 225 {\r
 226   return *(AODEvent()->GetV0s());\r
 227 }\r
 228 \r
 229 //______________________________________________________________________________\r
 230 TClonesArray& AliAnalysisTaskESDfilter::Vertices()\r
 231 {\r
 232   return *(AODEvent()->GetVertices());\r
 233 }\r
 234 \r
 235 //______________________________________________________________________________\r
 236 AliAODHeader* AliAnalysisTaskESDfilter::ConvertHeader(const AliESDEvent& esd)\r
 237 {\r
 238   AliCodeTimerAuto("",0);\r
 239   \r
 240   AliAODHeader* header = AODEvent()->GetHeader();\r
 241   \r
 242   header->SetRunNumber(esd.GetRunNumber());\r
 243   header->SetOfflineTrigger(fInputHandler->IsEventSelected()); // propagate the decision of the physics selection\r
 244 \r
 245   TTree* tree = fInputHandler->GetTree();\r
 246   if (tree) {\r
 247     TFile* file = tree->GetCurrentFile();\r
 248     if (file) header->SetESDFileName(file->GetName());\r
 249   }\r
 250   \r
 251   if (fOldESDformat) {\r
 252     header->SetBunchCrossNumber(0);\r
 253     header->SetOrbitNumber(0);\r
 254     header->SetPeriodNumber(0);\r
 255     header->SetEventType(0);\r
 256     header->SetMuonMagFieldScale(-999.);\r
 257     header->SetCentrality(0);       \r
 258     header->SetEventplane(0);\r
 259   } else {\r
 260     header->SetBunchCrossNumber(esd.GetBunchCrossNumber());\r
 261     header->SetOrbitNumber(esd.GetOrbitNumber());\r
 262     header->SetPeriodNumber(esd.GetPeriodNumber());\r
 263     header->SetEventType(esd.GetEventType());\r
 264     \r
 265     header->SetEventNumberESDFile(esd.GetHeader()->GetEventNumberInFile());\r
 266     if(const_cast<AliESDEvent&>(esd).GetCentrality()){\r
 267       header->SetCentrality(const_cast<AliESDEvent&>(esd).GetCentrality());\r
 268     }\r
 269     else{\r
 270       header->SetCentrality(0);\r
 271     }\r
 272     if(const_cast<AliESDEvent&>(esd).GetEventplane()){\r
 273       header->SetEventplane(const_cast<AliESDEvent&>(esd).GetEventplane());\r
 274     }\r
 275     else{\r
 276       header->SetEventplane(0);\r
 277     }\r
 278   }\r
 279   \r
 280   // Trigger\r
 281   header->SetFiredTriggerClasses(esd.GetFiredTriggerClasses());\r
 282   header->SetTriggerMask(esd.GetTriggerMask()); \r
 283   header->SetTriggerCluster(esd.GetTriggerCluster());\r
 284   header->SetL0TriggerInputs(esd.GetHeader()->GetL0TriggerInputs());    \r
 285   header->SetL1TriggerInputs(esd.GetHeader()->GetL1TriggerInputs());    \r
 286   header->SetL2TriggerInputs(esd.GetHeader()->GetL2TriggerInputs());    \r
 287   \r
 288   header->SetMagneticField(esd.GetMagneticField());\r
 289   header->SetMuonMagFieldScale(esd.GetCurrentDip()/6000.);\r
 290   header->SetZDCN1Energy(esd.GetZDCN1Energy());\r
 291   header->SetZDCP1Energy(esd.GetZDCP1Energy());\r
 292   header->SetZDCN2Energy(esd.GetZDCN2Energy());\r
 293   header->SetZDCP2Energy(esd.GetZDCP2Energy());\r
 294   header->SetZDCEMEnergy(esd.GetZDCEMEnergy(0),esd.GetZDCEMEnergy(1));\r
 295   \r
 296   \r
 297   Float_t diamxy[2]={esd.GetDiamondX(),esd.GetDiamondY()};\r
 298   Float_t diamcov[3]; \r
 299   esd.GetDiamondCovXY(diamcov);\r
 300   header->SetDiamond(diamxy,diamcov);\r
 301   header->SetDiamondZ(esd.GetDiamondZ(),esd.GetSigma2DiamondZ());\r
 302   \r
 303   return header;\r
 304 }\r
 305 \r
 306 //______________________________________________________________________________\r
 307 void AliAnalysisTaskESDfilter::ConvertCascades(const AliESDEvent& esd) \r
 308 {\r
 309   // Convert the cascades part of the ESD.\r
 310   // Return the number of cascades\r
 311  \r
 312   AliCodeTimerAuto("",0);\r
 313   \r
 314   // Create vertices starting from the most complex objects\r
 315   Double_t chi2 = 0.;\r
 316   \r
 317   const AliESDVertex* vtx = esd.GetPrimaryVertex();\r
 318   Double_t pos[3] = { 0. };\r
 319   Double_t covVtx[6] = { 0. };\r
 320   Double_t momBach[3]={0.};\r
 321   Double_t covTr[21]={0.};\r
 322   Double_t pid[10]={0.};\r
 323   AliAODPid* detpid(0x0);\r
 324   AliAODVertex* vV0FromCascade(0x0);\r
 325   AliAODv0* aodV0(0x0);\r
 326   AliAODcascade* aodCascade(0x0);\r
 327   AliAODTrack* aodTrack(0x0);\r
 328   Double_t momPos[3]={0.};\r
 329   Double_t momNeg[3] = { 0. };\r
 330   Double_t momPosAtV0vtx[3]={0.};\r
 331   Double_t momNegAtV0vtx[3]={0.};\r
 332 \r
 333   TClonesArray& verticesArray = Vertices();\r
 334   TClonesArray& tracksArray = Tracks();\r
 335   TClonesArray& cascadesArray = Cascades();\r
 336   \r
 337   // Cascades (Modified by A.Maire - February 2009)\r
 338   for (Int_t nCascade = 0; nCascade < esd.GetNumberOfCascades(); ++nCascade) {\r
 339     \r
 340     // 0- Preparation\r
 341     //\r
 342     AliESDcascade *esdCascade = esd.GetCascade(nCascade);\r
 343                 Int_t  idxPosFromV0Dghter  = esdCascade->GetPindex();\r
 344                 Int_t  idxNegFromV0Dghter  = esdCascade->GetNindex();\r
 345                 Int_t  idxBachFromCascade  = esdCascade->GetBindex();\r
 346     \r
 347     AliESDtrack  *esdCascadePos  = esd.GetTrack( idxPosFromV0Dghter);\r
 348     AliESDtrack  *esdCascadeNeg  = esd.GetTrack( idxNegFromV0Dghter);\r
 349     AliESDtrack  *esdCascadeBach = esd.GetTrack( idxBachFromCascade);\r
 350     \r
 351     // Identification of the V0 within the esdCascade (via both daughter track indices)\r
 352     AliESDv0 * currentV0   = 0x0;\r
 353     Int_t      idxV0FromCascade = -1;\r
 354     \r
 355     for (Int_t iV0=0; iV0<esd.GetNumberOfV0s(); ++iV0) {\r
 356       \r
 357       currentV0 = esd.GetV0(iV0);\r
 358       Int_t posCurrentV0 = currentV0->GetPindex();\r
 359       Int_t negCurrentV0 = currentV0->GetNindex();\r
 360       \r
 361       if (posCurrentV0==idxPosFromV0Dghter && negCurrentV0==idxNegFromV0Dghter) {\r
 362         idxV0FromCascade = iV0;\r
 363         break;\r
 364       }\r
 365     }\r
 366     \r
 367     if(idxV0FromCascade < 0){\r
 368       printf("Cascade - no matching for the V0 (index V0 = -1) ! Skip ... \n");\r
 369       continue;\r
 370     }// a priori, useless check, but safer ... in case of pb with tracks "out of bounds"\r
 371     \r
 372     AliESDv0 *esdV0FromCascade   = esd.GetV0(idxV0FromCascade);\r
 373         \r
 374     // 1 - Cascade selection \r
 375     \r
 376     //  AliESDVertex *esdPrimVtx = new AliESDVertex(*(esd.GetPrimaryVertex()));\r
 377     //  TList cascadeObjects;\r
 378     //  cascadeObjects.AddAt(esdV0FromCascade, 0);\r
 379     //  cascadeObjects.AddAt(esdCascadePos,    1);\r
 380     //  cascadeObjects.AddAt(esdCascadeNeg,    2);\r
 381     //  cascadeObjects.AddAt(esdCascade,       3);\r
 382     //  cascadeObjects.AddAt(esdCascadeBach,   4);\r
 383     //  cascadeObjects.AddAt(esdPrimVtx,       5);\r
 384     // \r
 385     //  UInt_t selectCascade = 0;\r
 386     //  if (fCascadeFilter) {\r
 387     //    // selectCascade = fCascadeFilter->IsSelected(&cascadeObjects); \r
 388     //          // FIXME AliESDCascadeCuts to be implemented ...\r
 389     // \r
 390     //          // Here we may encounter a moot point at the V0 level \r
 391     //          // between the cascade selections and the V0 ones :\r
 392     //          // the V0 selected along with the cascade (secondary V0) may \r
 393     //          // usually be removed from the dedicated V0 selections (prim V0) ...\r
 394     //          // -> To be discussed !\r
 395     // \r
 396     //    // this is a little awkward but otherwise the \r
 397     //    // list wants to access the pointer (delete it) \r
 398     //    // again when going out of scope\r
 399     //    delete cascadeObjects.RemoveAt(5); // esdPrimVtx created via copy construct\r
 400     //    esdPrimVtx = 0;\r
 401     //    if (!selectCascade) \r
 402     //      continue;\r
 403     //  }\r
 404     //  else{\r
 405     //    delete cascadeObjects.RemoveAt(5); // esdPrimVtx created via copy construct\r
 406     //    esdPrimVtx = 0;\r
 407     //  }\r
 408     \r
 409     // 2 - Add the cascade vertex\r
 410     \r
 411     esdCascade->GetXYZcascade(pos[0], pos[1], pos[2]);\r
 412     esdCascade->GetPosCovXi(covVtx);\r
 413     chi2 = esdCascade->GetChi2Xi(); \r
 414     \r
 415     AliAODVertex *vCascade = new(verticesArray[fNumberOfVertices++]) AliAODVertex( pos,\r
 416                                                                      covVtx,\r
 417                                                                      chi2, // FIXME = Chi2/NDF will be needed\r
 418                                                                      fPrimaryVertex,\r
 419                                                                      nCascade, // id\r
 420                                                                      AliAODVertex::kCascade);\r
 421     fPrimaryVertex->AddDaughter(vCascade);\r
 422     \r
 423 //    if (fDebug > 2) {\r
 424 //      printf("---- Cascade / Cascade Vertex (AOD) : \n");\r
 425 //      vCascade->Print();\r
 426 //    }\r
 427     \r
 428     if(esd.GetTOFHeader() && fIsPidOwner) fESDpid->SetTOFResponse(const_cast<AliESDEvent*>(&esd), AliESDpid::kBest_T0); //in case of AOD production strating form LHC10e without Tender. \r
 429 \r
 430 \r
 431     // 3 - Add the bachelor track from the cascade\r
 432     \r
 433     if (!fUsedTrack[idxBachFromCascade]) {\r
 434       \r
 435       esdCascadeBach->GetPxPyPz(momBach);\r
 436       esdCascadeBach->GetXYZ(pos);\r
 437       esdCascadeBach->GetCovarianceXYZPxPyPz(covTr);\r
 438       esdCascadeBach->GetESDpid(pid);\r
 439       \r
 440             fUsedTrack[idxBachFromCascade] = kTRUE;\r
 441             UInt_t selectInfo = 0;\r
 442             if (fTrackFilter) selectInfo = fTrackFilter->IsSelected(esdCascadeBach);\r
 443             if (fMChandler) fMChandler->SelectParticle(esdCascadeBach->GetLabel());\r
 444             aodTrack = new(tracksArray[fNumberOfTracks++]) AliAODTrack(esdCascadeBach->GetID(),\r
 445                                                                            esdCascadeBach->GetLabel(), \r
 446                                                                            momBach, \r
 447                                                                            kTRUE,\r
 448                                                                            pos,\r
 449                                                                            kFALSE, // Why kFALSE for "isDCA" ? FIXME\r
 450                                                                            covTr, \r
 451                                                                            (Short_t)esdCascadeBach->GetSign(),\r
 452                                                                            esdCascadeBach->GetITSClusterMap(), \r
 453                                                                            pid,\r
 454                                                                            vCascade,\r
 455                                                                            kTRUE,  // usedForVtxFit = kFALSE ? FIXME\r
 456                                                                            vtx->UsesTrack(esdCascadeBach->GetID()),\r
 457                                                                            AliAODTrack::kSecondary,\r
 458                                                                            selectInfo);\r
 459             aodTrack->SetTPCClusterMap(esdCascadeBach->GetTPCClusterMap());\r
 460             aodTrack->SetTPCSharedMap (esdCascadeBach->GetTPCSharedMap());\r
 461             aodTrack->SetChi2perNDF(Chi2perNDF(esdCascadeBach));\r
 462             aodTrack->SetTPCPointsF(esdCascadeBach->GetTPCNclsF());\r
 463             fAODTrackRefs->AddAt(aodTrack,idxBachFromCascade);\r
 464             \r
 465             if (esdCascadeBach->GetSign() > 0) ++fNumberOfPositiveTracks;\r
 466             aodTrack->ConvertAliPIDtoAODPID();\r
 467             aodTrack->SetFlags(esdCascadeBach->GetStatus());\r
 468       SetAODPID(esdCascadeBach,aodTrack,detpid);\r
 469     }\r
 470     else {\r
 471             aodTrack = static_cast<AliAODTrack*>( fAODTrackRefs->At(idxBachFromCascade) );\r
 472     }\r
 473     \r
 474     vCascade->AddDaughter(aodTrack);\r
 475     \r
 476 //    if (fDebug > 4) {\r
 477 //      printf("---- Cascade / bach dghter : \n");\r
 478 //      aodTrack->Print();\r
 479 //    }\r
 480     \r
 481     \r
 482     // 4 - Add the V0 from the cascade. \r
 483     // = V0vtx + both pos and neg daughter tracks + the aodV0 itself\r
 484     //\r
 485     \r
 486     if ( !fUsedV0[idxV0FromCascade] ) {\r
 487       // 4.A - if VO structure hasn't been created yet\r
 488       \r
 489       // 4.A.1 - Create the V0 vertex of the cascade\r
 490       \r
 491       esdV0FromCascade->GetXYZ(pos[0], pos[1], pos[2]);\r
 492       esdV0FromCascade->GetPosCov(covVtx);\r
 493       chi2 = esdV0FromCascade->GetChi2V0();  // = chi2/NDF since NDF = 2*2-3 ?\r
 494                         \r
 495       vV0FromCascade = new(verticesArray[fNumberOfVertices++]) AliAODVertex(pos,\r
 496                                                                          covVtx,\r
 497                                                                          chi2,\r
 498                                                                          vCascade,\r
 499                                                                          idxV0FromCascade, //id of ESDv0\r
 500                                                                          AliAODVertex::kV0);\r
 501       // Note:\r
 502       //    one V0 can be used by several cascades.\r
 503       // So, one AOD V0 vtx can have several parent vtx.\r
 504       // This is not directly allowed by AliAODvertex.\r
 505       // Setting the parent vtx (here = param "vCascade") doesn't lead to a crash\r
 506       // but to a problem of consistency within AODEvent.\r
 507       // -> See below paragraph 4.B, for the proposed treatment of such a case.\r
 508       \r
 509       // Add the vV0FromCascade to the aodVOVtxRefs\r
 510       fAODV0VtxRefs->AddAt(vV0FromCascade,idxV0FromCascade);\r
 511       \r
 512       \r
 513       // 4.A.2 - Add the positive tracks from the V0\r
 514       \r
 515       esdCascadePos->GetPxPyPz(momPos);\r
 516       esdCascadePos->GetXYZ(pos);\r
 517       esdCascadePos->GetCovarianceXYZPxPyPz(covTr);\r
 518       esdCascadePos->GetESDpid(pid);\r
 519       \r
 520       \r
 521       if (!fUsedTrack[idxPosFromV0Dghter]) {\r
 522         fUsedTrack[idxPosFromV0Dghter] = kTRUE;\r
 523         \r
 524         UInt_t selectInfo = 0;\r
 525         if (fTrackFilter) selectInfo = fTrackFilter->IsSelected(esdCascadePos);\r
 526         if(fMChandler) fMChandler->SelectParticle(esdCascadePos->GetLabel());\r
 527         aodTrack = new(tracksArray[fNumberOfTracks++]) \r
 528         AliAODTrack(  esdCascadePos->GetID(),\r
 529                     esdCascadePos->GetLabel(), \r
 530                     momPos, \r
 531                     kTRUE,\r
 532                     pos,\r
 533                     kFALSE, // Why kFALSE for "isDCA" ? FIXME\r
 534                     covTr, \r
 535                     (Short_t)esdCascadePos->GetSign(),\r
 536                     esdCascadePos->GetITSClusterMap(), \r
 537                     pid,\r
 538                     vV0FromCascade,\r
 539                     kTRUE,  // usedForVtxFit = kFALSE ? FIXME\r
 540                     vtx->UsesTrack(esdCascadePos->GetID()),\r
 541                     AliAODTrack::kSecondary,\r
 542                     selectInfo);\r
 543         aodTrack->SetTPCClusterMap(esdCascadePos->GetTPCClusterMap());\r
 544         aodTrack->SetTPCSharedMap (esdCascadePos->GetTPCSharedMap());\r
 545         aodTrack->SetChi2perNDF(Chi2perNDF(esdCascadePos));\r
 546         aodTrack->SetTPCPointsF(esdCascadePos->GetTPCNclsF());\r
 547         fAODTrackRefs->AddAt(aodTrack,idxPosFromV0Dghter);\r
 548         \r
 549         if (esdCascadePos->GetSign() > 0) ++fNumberOfPositiveTracks;\r
 550         aodTrack->ConvertAliPIDtoAODPID();\r
 551         aodTrack->SetFlags(esdCascadePos->GetStatus());\r
 552         SetAODPID(esdCascadePos,aodTrack,detpid);\r
 553       }\r
 554       else {\r
 555         aodTrack = static_cast<AliAODTrack*>(fAODTrackRefs->At(idxPosFromV0Dghter));\r
 556       }\r
 557       vV0FromCascade->AddDaughter(aodTrack);\r
 558       \r
 559       \r
 560       // 4.A.3 - Add the negative tracks from the V0\r
 561       \r
 562       esdCascadeNeg->GetPxPyPz(momNeg);\r
 563       esdCascadeNeg->GetXYZ(pos);\r
 564       esdCascadeNeg->GetCovarianceXYZPxPyPz(covTr);\r
 565       esdCascadeNeg->GetESDpid(pid);\r
 566       \r
 567       \r
 568       if (!fUsedTrack[idxNegFromV0Dghter]) {\r
 569         fUsedTrack[idxNegFromV0Dghter] = kTRUE;\r
 570         \r
 571         UInt_t selectInfo = 0;\r
 572         if (fTrackFilter) selectInfo = fTrackFilter->IsSelected(esdCascadeNeg);\r
 573         if(fMChandler)fMChandler->SelectParticle(esdCascadeNeg->GetLabel());\r
 574         aodTrack = new(tracksArray[fNumberOfTracks++]) AliAODTrack(  esdCascadeNeg->GetID(),\r
 575                                                       esdCascadeNeg->GetLabel(),\r
 576                                                       momNeg,\r
 577                                                       kTRUE,\r
 578                                                       pos,\r
 579                                                       kFALSE, // Why kFALSE for "isDCA" ? FIXME\r
 580                                                       covTr, \r
 581                                                       (Short_t)esdCascadeNeg->GetSign(),\r
 582                                                       esdCascadeNeg->GetITSClusterMap(), \r
 583                                                       pid,\r
 584                                                       vV0FromCascade,\r
 585                                                       kTRUE,  // usedForVtxFit = kFALSE ? FIXME\r
 586                                                       vtx->UsesTrack(esdCascadeNeg->GetID()),\r
 587                                                       AliAODTrack::kSecondary,\r
 588                                                       selectInfo);\r
 589         aodTrack->SetTPCClusterMap(esdCascadeNeg->GetTPCClusterMap());\r
 590         aodTrack->SetTPCSharedMap (esdCascadeNeg->GetTPCSharedMap());\r
 591         aodTrack->SetChi2perNDF(Chi2perNDF(esdCascadeNeg));\r
 592         aodTrack->SetTPCPointsF(esdCascadeNeg->GetTPCNclsF());\r
 593         fAODTrackRefs->AddAt(aodTrack,idxNegFromV0Dghter);\r
 594         \r
 595         if (esdCascadeNeg->GetSign() > 0) ++fNumberOfPositiveTracks;\r
 596         aodTrack->ConvertAliPIDtoAODPID();\r
 597         aodTrack->SetFlags(esdCascadeNeg->GetStatus());\r
 598         SetAODPID(esdCascadeNeg,aodTrack,detpid);\r
 599       }\r
 600       else {\r
 601         aodTrack = static_cast<AliAODTrack*>(fAODTrackRefs->At(idxNegFromV0Dghter));\r
 602       }\r
 603       \r
 604       vV0FromCascade->AddDaughter(aodTrack);\r
 605       \r
 606                         \r
 607       // 4.A.4 - Add the V0 from cascade to the V0 array\r
 608       \r
 609       Double_t  dcaV0Daughters      = esdV0FromCascade->GetDcaV0Daughters();\r
 610       Double_t  dcaV0ToPrimVertex   = esdV0FromCascade->GetD( esd.GetPrimaryVertex()->GetX(),\r
 611                                                              esd.GetPrimaryVertex()->GetY(),\r
 612                                                              esd.GetPrimaryVertex()->GetZ() );\r
 613       esdV0FromCascade->GetPPxPyPz( momPosAtV0vtx[0],momPosAtV0vtx[1],momPosAtV0vtx[2] ); \r
 614       esdV0FromCascade->GetNPxPyPz( momNegAtV0vtx[0],momNegAtV0vtx[1],momNegAtV0vtx[2] ); \r
 615       \r
 616       Double_t dcaDaughterToPrimVertex[2] = { 999., 999.}; // ..[0] = DCA in (x,y) for Pos and ..[1] = Neg\r
 617       dcaDaughterToPrimVertex[0] = TMath::Abs(esdCascadePos->GetD(      esd.GetPrimaryVertex()->GetX(),\r
 618                                                                   esd.GetPrimaryVertex()->GetY(),\r
 619                                                                   esd.GetMagneticField())        );\r
 620       dcaDaughterToPrimVertex[1] = TMath::Abs(esdCascadeNeg->GetD(      esd.GetPrimaryVertex()->GetX(),\r
 621                                                                   esd.GetPrimaryVertex()->GetY(),\r
 622                                                                   esd.GetMagneticField())        );\r
 623       \r
 624       aodV0 = new(V0s()[fNumberOfV0s++]) AliAODv0( vV0FromCascade, \r
 625                                                   dcaV0Daughters,\r
 626                                                   dcaV0ToPrimVertex, \r
 627                                                   momPosAtV0vtx, \r
 628                                                   momNegAtV0vtx, \r
 629                                                   dcaDaughterToPrimVertex); \r
 630       // set the aod v0 on-the-fly status\r
 631       aodV0->SetOnFlyStatus(esdV0FromCascade->GetOnFlyStatus());\r
 632       \r
 633       // Add the aodV0 to the aodVORefs\r
 634       fAODV0Refs->AddAt(aodV0,idxV0FromCascade);\r
 635       \r
 636       fUsedV0[idxV0FromCascade] = kTRUE;\r
 637       \r
 638     } else { \r
 639       // 4.B - if V0 structure already used\r
 640       \r
 641       // Note :\r
 642       //    one V0 can be used by several cascades (frequent in PbPb evts) : \r
 643       // same V0 which used but attached to different bachelor tracks\r
 644       // -> aodVORefs and fAODV0VtxRefs are needed.\r
 645       // Goal : avoid a redundancy of the info in "Vertices" and "v0s" clones array.\r
 646       \r
 647       vV0FromCascade = static_cast<AliAODVertex*>( fAODV0VtxRefs->At(idxV0FromCascade) );\r
 648       aodV0          = static_cast<AliAODv0*>    ( fAODV0Refs   ->At(idxV0FromCascade) );\r
 649       \r
 650       // - Treatment of the parent for such a "re-used" V0 :\r
 651       // Insert the cascade that reuses the V0 vertex in the lineage chain\r
 652       // Before : vV0 -> vCascade1 -> vPrimary\r
 653       //  - Hyp : cascade2 uses the same V0 as cascade1\r
 654       //  After :  vV0 -> vCascade2 -> vCascade1 -> vPrimary\r
 655       \r
 656       AliAODVertex *vCascadePreviousParent = static_cast<AliAODVertex*> (vV0FromCascade->GetParent());\r
 657                         vV0FromCascade->SetParent(vCascade);\r
 658                         vCascade      ->SetParent(vCascadePreviousParent);\r
 659       \r
 660 //      if(fDebug > 2)  \r
 661 //        printf("---- Cascade / Lineage insertion\n"\r
 662 //               "Parent of V0 vtx                 = Cascade vtx %p\n"\r
 663 //               "Parent of the cascade vtx        = Cascade vtx %p\n"\r
 664 //               "Parent of the parent cascade vtx = Cascade vtx %p\n", \r
 665 //               static_cast<void*> (vV0FromCascade->GetParent()),\r
 666 //               static_cast<void*> (vCascade->GetParent()),\r
 667 //               static_cast<void*> (vCascadePreviousParent->GetParent()) );\r
 668       \r
 669     }// end if V0 structure already used\r
 670     \r
 671 //    if (fDebug > 2) {\r
 672 //      printf("---- Cascade / V0 vertex: \n");\r
 673 //      vV0FromCascade->Print();\r
 674 //    }\r
 675 //    \r
 676 //    if (fDebug > 4) {\r
 677 //      printf("---- Cascade / pos dghter : \n");\r
 678 //                      aodTrack->Print();\r
 679 //      printf("---- Cascade / neg dghter : \n");\r
 680 //                      aodTrack->Print();\r
 681 //      printf("---- Cascade / aodV0 : \n");\r
 682 //                      aodV0->Print();\r
 683 //    }\r
 684     \r
 685     // In any case (used V0 or not), add the V0 vertex to the cascade one.\r
 686     vCascade->AddDaughter(vV0FromCascade);      \r
 687     \r
 688                 \r
 689     // 5 - Add the primary track of the cascade (if any)\r
 690     \r
 691     \r
 692     // 6 - Add the cascade to the AOD array of cascades\r
 693     \r
 694     Double_t dcaBachToPrimVertexXY = TMath::Abs(esdCascadeBach->GetD(esd.GetPrimaryVertex()->GetX(),\r
 695                                                                      esd.GetPrimaryVertex()->GetY(),\r
 696                                                                      esd.GetMagneticField())        );\r
 697     \r
 698     Double_t momBachAtCascadeVtx[3]={0.};\r
 699 \r
 700     esdCascade->GetBPxPyPz(momBachAtCascadeVtx[0], momBachAtCascadeVtx[1], momBachAtCascadeVtx[2]);\r
 701     \r
 702     aodCascade = new(cascadesArray[fNumberOfCascades++]) AliAODcascade(  vCascade,\r
 703                                                           esdCascade->Charge(),\r
 704                                                           esdCascade->GetDcaXiDaughters(),\r
 705                                                           -999.,\r
 706                                                           // DCAXiToPrimVtx -> needs to be calculated   ----|\r
 707                                                           // doesn't exist at ESD level;\r
 708                                                           // See AODcascade::DcaXiToPrimVertex(Double, Double, Double)\r
 709                                                           dcaBachToPrimVertexXY,\r
 710                                                           momBachAtCascadeVtx,\r
 711                                                           *aodV0);\r
 712     \r
 713     if (fDebug > 10) {\r
 714       printf("---- Cascade / AOD cascade : \n\n");\r
 715       aodCascade->PrintXi(fPrimaryVertex->GetX(), fPrimaryVertex->GetY(), fPrimaryVertex->GetZ());\r
 716     }\r
 717     \r
 718   } // end of the loop on cascades\r
 719   \r
 720   Cascades().Expand(fNumberOfCascades);\r
 721 }\r
 722 \r
 723 //______________________________________________________________________________\r
 724 void AliAnalysisTaskESDfilter::ConvertV0s(const AliESDEvent& esd)\r
 725 {\r
 726   // Access to the AOD container of V0s\r
 727   \r
 728   AliCodeTimerAuto("",0);\r
 729 \r
 730   //\r
 731   // V0s\r
 732   //\r
 733   \r
 734   Double_t pos[3] = { 0. };      \r
 735   Double_t chi2(0.0);\r
 736   Double_t covVtx[6] = { 0. };\r
 737   Double_t momPos[3]={0.};\r
 738   Double_t covTr[21]={0.};\r
 739   Double_t pid[10]={0.};\r
 740   AliAODTrack* aodTrack(0x0);\r
 741   AliAODPid* detpid(0x0);\r
 742   Double_t momNeg[3]={0.};\r
 743   Double_t momPosAtV0vtx[3]={0.};\r
 744   Double_t momNegAtV0vtx[3]={0.};\r
 745     \r
 746   for (Int_t nV0 = 0; nV0 < esd.GetNumberOfV0s(); ++nV0) \r
 747   {\r
 748     if (fUsedV0[nV0]) continue; // skip if already added to the AOD\r
 749     \r
 750     AliESDv0 *v0 = esd.GetV0(nV0);\r
 751     Int_t posFromV0 = v0->GetPindex();\r
 752     Int_t negFromV0 = v0->GetNindex();\r
 753     \r
 754     // V0 selection \r
 755     //\r
 756     AliESDVertex *esdVtx   = new AliESDVertex(*(esd.GetPrimaryVertex()));\r
 757     AliESDtrack  *esdV0Pos = esd.GetTrack(posFromV0);\r
 758     AliESDtrack  *esdV0Neg = esd.GetTrack(negFromV0);\r
 759     TList v0objects;\r
 760     v0objects.AddAt(v0,                      0);\r
 761     v0objects.AddAt(esdV0Pos,                1);\r
 762     v0objects.AddAt(esdV0Neg,                2);\r
 763     v0objects.AddAt(esdVtx,                  3);\r
 764     UInt_t selectV0 = 0;\r
 765     if (fV0Filter) {\r
 766       selectV0 = fV0Filter->IsSelected(&v0objects);\r
 767       // this is a little awkward but otherwise the \r
 768       // list wants to access the pointer (delete it) \r
 769       // again when going out of scope\r
 770       delete v0objects.RemoveAt(3); // esdVtx created via copy construct\r
 771       esdVtx = 0;\r
 772       if (!selectV0) \r
 773         continue;\r
 774     }\r
 775     else{\r
 776       delete v0objects.RemoveAt(3); // esdVtx created via copy construct\r
 777       esdVtx = 0;\r
 778     }\r
 779     \r
 780     v0->GetXYZ(pos[0], pos[1], pos[2]);\r
 781     \r
 782     if (!fOldESDformat) {\r
 783             chi2 = v0->GetChi2V0(); // = chi2/NDF since NDF = 2*2-3\r
 784             v0->GetPosCov(covVtx);\r
 785     } else {\r
 786             chi2 = -999.;\r
 787             for (Int_t i = 0; i < 6; i++)  covVtx[i] = 0.;\r
 788     }\r
 789     \r
 790     \r
 791     AliAODVertex * vV0 = \r
 792           new(Vertices()[fNumberOfVertices++]) AliAODVertex(pos,\r
 793                                                       covVtx,\r
 794                                                       chi2,\r
 795                                                       fPrimaryVertex,\r
 796                                                       nV0,\r
 797                                                       AliAODVertex::kV0);\r
 798     fPrimaryVertex->AddDaughter(vV0);\r
 799     \r
 800     \r
 801     // Add the positive tracks from the V0\r
 802     \r
 803 \r
 804     esdV0Pos->GetPxPyPz(momPos);\r
 805     esdV0Pos->GetXYZ(pos);\r
 806     esdV0Pos->GetCovarianceXYZPxPyPz(covTr);\r
 807     esdV0Pos->GetESDpid(pid);\r
 808     \r
 809     const AliESDVertex *vtx = esd.GetPrimaryVertex();\r
 810 \r
 811     if (!fUsedTrack[posFromV0]) {\r
 812             fUsedTrack[posFromV0] = kTRUE;\r
 813             UInt_t selectInfo = 0;\r
 814             if (fTrackFilter) selectInfo = fTrackFilter->IsSelected(esdV0Pos);\r
 815             if(fMChandler)fMChandler->SelectParticle(esdV0Pos->GetLabel());\r
 816             aodTrack = new(Tracks()[fNumberOfTracks++]) AliAODTrack(esdV0Pos->GetID(),\r
 817                                                     esdV0Pos->GetLabel(), \r
 818                                                     momPos, \r
 819                                                     kTRUE,\r
 820                                                     pos,\r
 821                                                     kFALSE,\r
 822                                                     covTr, \r
 823                                                     (Short_t)esdV0Pos->GetSign(),\r
 824                                                     esdV0Pos->GetITSClusterMap(), \r
 825                                                     pid,\r
 826                                                     vV0,\r
 827                                                     kTRUE,  // check if this is right\r
 828                                                     vtx->UsesTrack(esdV0Pos->GetID()),\r
 829                                                     AliAODTrack::kSecondary,\r
 830                                                     selectInfo);\r
 831             aodTrack->SetTPCClusterMap(esdV0Pos->GetTPCClusterMap());\r
 832             aodTrack->SetTPCSharedMap (esdV0Pos->GetTPCSharedMap());\r
 833             aodTrack->SetChi2perNDF(Chi2perNDF(esdV0Pos));\r
 834             aodTrack->SetTPCPointsF(esdV0Pos->GetTPCNclsF());\r
 835             fAODTrackRefs->AddAt(aodTrack,posFromV0);\r
 836             //      if (fDebug > 0) printf("-------------------Bo: pos track from original pt %.3f \n",aodTrack->Pt());\r
 837             if (esdV0Pos->GetSign() > 0) ++fNumberOfPositiveTracks;\r
 838             aodTrack->ConvertAliPIDtoAODPID();\r
 839             aodTrack->SetFlags(esdV0Pos->GetStatus());\r
 840       SetAODPID(esdV0Pos,aodTrack,detpid);\r
 841     }\r
 842     else {\r
 843             aodTrack = static_cast<AliAODTrack*>(fAODTrackRefs->At(posFromV0));\r
 844             //      if (fDebug > 0) printf("-------------------Bo pos track from refArray pt %.3f \n",aodTrack->Pt());\r
 845     }\r
 846     vV0->AddDaughter(aodTrack);\r
 847     \r
 848     // Add the negative tracks from the V0\r
 849     \r
 850     esdV0Neg->GetPxPyPz(momNeg);\r
 851     esdV0Neg->GetXYZ(pos);\r
 852     esdV0Neg->GetCovarianceXYZPxPyPz(covTr);\r
 853     esdV0Neg->GetESDpid(pid);\r
 854     \r
 855     if (!fUsedTrack[negFromV0]) {\r
 856             fUsedTrack[negFromV0] = kTRUE;\r
 857             UInt_t selectInfo = 0;\r
 858             if (fTrackFilter) selectInfo = fTrackFilter->IsSelected(esdV0Neg);\r
 859             if(fMChandler)fMChandler->SelectParticle(esdV0Neg->GetLabel());\r
 860             aodTrack = new(Tracks()[fNumberOfTracks++]) AliAODTrack(esdV0Neg->GetID(),\r
 861                                                     esdV0Neg->GetLabel(),\r
 862                                                     momNeg,\r
 863                                                     kTRUE,\r
 864                                                     pos,\r
 865                                                     kFALSE,\r
 866                                                     covTr, \r
 867                                                     (Short_t)esdV0Neg->GetSign(),\r
 868                                                     esdV0Neg->GetITSClusterMap(), \r
 869                                                     pid,\r
 870                                                     vV0,\r
 871                                                     kTRUE,  // check if this is right\r
 872                                                     vtx->UsesTrack(esdV0Neg->GetID()),\r
 873                                                     AliAODTrack::kSecondary,\r
 874                                                     selectInfo);\r
 875             aodTrack->SetTPCClusterMap(esdV0Neg->GetTPCClusterMap());\r
 876             aodTrack->SetTPCSharedMap (esdV0Neg->GetTPCSharedMap());\r
 877             aodTrack->SetChi2perNDF(Chi2perNDF(esdV0Neg));\r
 878             aodTrack->SetTPCPointsF(esdV0Neg->GetTPCNclsF());\r
 879             \r
 880             fAODTrackRefs->AddAt(aodTrack,negFromV0);\r
 881             //      if (fDebug > 0) printf("-------------------Bo: neg track from original pt %.3f \n",aodTrack->Pt());\r
 882             if (esdV0Neg->GetSign() > 0) ++fNumberOfPositiveTracks;\r
 883             aodTrack->ConvertAliPIDtoAODPID();\r
 884             aodTrack->SetFlags(esdV0Neg->GetStatus());\r
 885       SetAODPID(esdV0Neg,aodTrack,detpid);\r
 886     }\r
 887     else {\r
 888             aodTrack = static_cast<AliAODTrack*>(fAODTrackRefs->At(negFromV0));\r
 889             //      if (fDebug > 0) printf("-------------------Bo neg track from refArray pt %.3f \n",aodTrack->Pt());\r
 890     }\r
 891     vV0->AddDaughter(aodTrack);\r
 892     \r
 893     \r
 894     // Add the V0 the V0 array as well\r
 895     \r
 896     Double_t  dcaV0Daughters      = v0->GetDcaV0Daughters();\r
 897     Double_t  dcaV0ToPrimVertex   = v0->GetD(esd.GetPrimaryVertex()->GetX(),\r
 898                                              esd.GetPrimaryVertex()->GetY(),\r
 899                                              esd.GetPrimaryVertex()->GetZ());\r
 900     \r
 901     v0->GetPPxPyPz(momPosAtV0vtx[0],momPosAtV0vtx[1],momPosAtV0vtx[2]); \r
 902     v0->GetNPxPyPz(momNegAtV0vtx[0],momNegAtV0vtx[1],momNegAtV0vtx[2]); \r
 903     \r
 904     Double_t dcaDaughterToPrimVertex[2] = { 999., 999.}; // ..[0] = DCA in (x,y) for Pos and ..[1] = Neg\r
 905     dcaDaughterToPrimVertex[0] = TMath::Abs(esdV0Pos->GetD(  esd.GetPrimaryVertex()->GetX(),\r
 906                                                            esd.GetPrimaryVertex()->GetY(),\r
 907                                                            esd.GetMagneticField()) );\r
 908     dcaDaughterToPrimVertex[1] = TMath::Abs(esdV0Neg->GetD(  esd.GetPrimaryVertex()->GetX(),\r
 909                                                            esd.GetPrimaryVertex()->GetY(),\r
 910                                                            esd.GetMagneticField()) );\r
 911     \r
 912     AliAODv0* aodV0 = new(V0s()[fNumberOfV0s++]) AliAODv0(vV0, \r
 913                                                 dcaV0Daughters,\r
 914                                                 dcaV0ToPrimVertex,\r
 915                                                 momPosAtV0vtx,\r
 916                                                 momNegAtV0vtx,\r
 917                                                 dcaDaughterToPrimVertex);\r
 918     \r
 919     // set the aod v0 on-the-fly status\r
 920     aodV0->SetOnFlyStatus(v0->GetOnFlyStatus());\r
 921   }//End of loop on V0s \r
 922   \r
 923   V0s().Expand(fNumberOfV0s);    \r
 924 }\r
 925 \r
 926 //______________________________________________________________________________\r
 927 void AliAnalysisTaskESDfilter::ConvertTPCOnlyTracks(const AliESDEvent& esd)\r
 928 {\r
 929   // Convert TPC only tracks\r
 930   \r
 931   AliCodeTimerAuto("",0);\r
 932   \r
 933   // Loop over the tracks and extract and mask out all aod tracks that pass the selections for AODt racks\r
 934   for(int it = 0;it < fNumberOfTracks;++it)\r
 935   {\r
 936     AliAODTrack *tr = (AliAODTrack*)(Tracks().At(it));\r
 937     if(!tr)continue;\r
 938     UInt_t map = tr->GetFilterMap();\r
 939     if(map&fTPCOnlyFilterMask){\r
 940       // we only reset the track select ionfo, no deletion...\r
 941       tr->SetFilterMap(map&~fTPCOnlyFilterMask);\r
 942     }\r
 943   }\r
 944   // Loop over the ESD trcks and pick out the tracks passing TPC only cuts\r
 945   \r
 946   \r
 947   const AliESDVertex *vtxSPD = esd.GetPrimaryVertexSPD();\r
 948   const AliESDVertex *vtx = esd.GetPrimaryVertex();\r
 949 \r
 950   Double_t pos[3] = { 0. };      \r
 951   Double_t covTr[21]={0.};\r
 952   Double_t pid[10]={0.};  \r
 953   Double_t p[3] = { 0. };\r
 954   AliAODTrack* aodTrack(0x0);\r
 955   \r
 956   for (Int_t nTrack = 0; nTrack < esd.GetNumberOfTracks(); ++nTrack) \r
 957   {\r
 958     AliESDtrack* esdTrack = esd.GetTrack(nTrack); //carefull do not modify it othwise  need to work with a copy \r
 959     \r
 960     UInt_t selectInfo = 0;\r
 961     //\r
 962     // Track selection\r
 963     if (fTrackFilter) {\r
 964       selectInfo = fTrackFilter->IsSelected(esdTrack);\r
 965     }\r
 966     selectInfo &= fTPCOnlyFilterMask;\r
 967     if (!selectInfo)continue;\r
 968     \r
 969     // create a tpc only tracl\r
 970     AliESDtrack *track = AliESDtrackCuts::GetTPCOnlyTrack(const_cast<AliESDEvent*>(&esd),esdTrack->GetID());\r
 971     if(!track) continue;\r
 972     \r
 973     if(track->Pt()>0.)\r
 974     {\r
 975       // only constrain tracks above threshold\r
 976       AliExternalTrackParam exParam;\r
 977       // take the B-field from the ESD, no 3D fieldMap available at this point\r
 978       Bool_t relate = false;\r
 979       relate = track->RelateToVertexTPC(vtxSPD,esd.GetMagneticField(),kVeryBig,&exParam);\r
 980       if(!relate){\r
 981         delete track;\r
 982         continue;\r
 983       }\r
 984       track->Set(exParam.GetX(),exParam.GetAlpha(),exParam.GetParameter(),exParam.GetCovariance());\r
 985     }\r
 986     \r
 987     track->GetPxPyPz(p);\r
 988     track->GetXYZ(pos);\r
 989     track->GetCovarianceXYZPxPyPz(covTr);\r
 990     track->GetESDpid(pid);\r
 991     if(fMChandler)fMChandler->SelectParticle(esdTrack->GetLabel());\r
 992     aodTrack = new(Tracks()[fNumberOfTracks++]) AliAODTrack((track->GetID()+1)*-1,\r
 993                                                             track->GetLabel(),\r
 994                                                             p,\r
 995                                                             kTRUE,\r
 996                                                             pos,\r
 997                                                             kFALSE,\r
 998                                                             covTr, \r
 999                                                             (Short_t)track->GetSign(),\r
1000                                                             track->GetITSClusterMap(), \r
1001                                                             pid,\r
1002                                                             fPrimaryVertex,\r
1003                                                             kTRUE, // check if this is right\r
1004                                                             vtx->UsesTrack(track->GetID()),\r
1005                                                             AliAODTrack::kPrimary, \r
1006                                                             selectInfo);\r
1007     aodTrack->SetTPCClusterMap(track->GetTPCClusterMap());\r
1008     aodTrack->SetTPCSharedMap (track->GetTPCSharedMap());\r
1009     Float_t ndf = track->GetTPCNcls()+1 - 5 ;\r
1010     if(ndf>0){\r
1011       aodTrack->SetChi2perNDF(track->GetTPCchi2()/ndf);\r
1012     }\r
1013     else{\r
1014       aodTrack->SetChi2perNDF(-1);\r
1015     }\r
1016     aodTrack->SetFlags(track->GetStatus());\r
1017     aodTrack->SetTPCPointsF(track->GetTPCNclsF());\r
1018 \r
1019     delete track;\r
1020   } // end of loop on tracks\r
1021   \r
1022 }\r
1023 \r
1024 //______________________________________________________________________________\r
1025 void AliAnalysisTaskESDfilter::ConvertTracks(const AliESDEvent& esd)\r
1026 {\r
1027   // Tracks (primary and orphan)\r
1028 \r
1029   AliCodeTimerAuto("",0);\r
1030   \r
1031   AliDebug(1,Form("NUMBER OF ESD TRACKS %5d\n", esd.GetNumberOfTracks()));\r
1032   \r
1033   const AliESDVertex *vtx = esd.GetPrimaryVertex();\r
1034   Double_t p[3] = { 0. };\r
1035   Double_t pos[3] = { 0. };\r
1036   Double_t covTr[21] = { 0. };\r
1037   Double_t pid[10] = { 0. };\r
1038   AliAODTrack* aodTrack(0x0);\r
1039   AliAODPid* detpid(0x0);\r
1040   \r
1041   for (Int_t nTrack = 0; nTrack < esd.GetNumberOfTracks(); ++nTrack) \r
1042   {\r
1043     if (fUsedTrack[nTrack]) continue;\r
1044     \r
1045     AliESDtrack *esdTrack = esd.GetTrack(nTrack);\r
1046     UInt_t selectInfo = 0;\r
1047     //\r
1048     // Track selection\r
1049     if (fTrackFilter) {\r
1050             selectInfo = fTrackFilter->IsSelected(esdTrack);\r
1051             if (!selectInfo && !vtx->UsesTrack(esdTrack->GetID())) continue;\r
1052     }\r
1053     \r
1054     \r
1055     esdTrack->GetPxPyPz(p);\r
1056     esdTrack->GetXYZ(pos);\r
1057     esdTrack->GetCovarianceXYZPxPyPz(covTr);\r
1058     esdTrack->GetESDpid(pid);\r
1059     if(fMChandler)fMChandler->SelectParticle(esdTrack->GetLabel());\r
1060     fPrimaryVertex->AddDaughter(aodTrack =\r
1061                          new(Tracks()[fNumberOfTracks++]) AliAODTrack(esdTrack->GetID(),\r
1062                                                             esdTrack->GetLabel(),\r
1063                                                             p,\r
1064                                                             kTRUE,\r
1065                                                             pos,\r
1066                                                             kFALSE,\r
1067                                                             covTr, \r
1068                                                             (Short_t)esdTrack->GetSign(),\r
1069                                                             esdTrack->GetITSClusterMap(), \r
1070                                                             pid,\r
1071                                                             fPrimaryVertex,\r
1072                                                             kTRUE, // check if this is right\r
1073                                                             vtx->UsesTrack(esdTrack->GetID()),\r
1074                                                             AliAODTrack::kPrimary, \r
1075                                                             selectInfo)\r
1076                          );\r
1077     aodTrack->SetTPCClusterMap(esdTrack->GetTPCClusterMap());\r
1078     aodTrack->SetTPCSharedMap (esdTrack->GetTPCSharedMap());\r
1079     aodTrack->SetChi2perNDF(Chi2perNDF(esdTrack));\r
1080     aodTrack->SetTPCPointsF(esdTrack->GetTPCNclsF());\r
1081 \r
1082     fAODTrackRefs->AddAt(aodTrack, nTrack);\r
1083     \r
1084     \r
1085     if (esdTrack->GetSign() > 0) ++fNumberOfPositiveTracks;\r
1086     aodTrack->SetFlags(esdTrack->GetStatus());\r
1087     aodTrack->ConvertAliPIDtoAODPID();\r
1088     SetAODPID(esdTrack,aodTrack,detpid);\r
1089   } // end of loop on tracks\r
1090 }\r
1091 \r
1092 //______________________________________________________________________________\r
1093 void AliAnalysisTaskESDfilter::ConvertPmdClusters(const AliESDEvent& esd)\r
1094 {\r
1095   AliCodeTimerAuto("",0);\r
1096   Int_t jPmdClusters=0;\r
1097   // Access to the AOD container of PMD clusters\r
1098   TClonesArray &pmdClusters = *(AODEvent()->GetPmdClusters());\r
1099   for (Int_t iPmd = 0; iPmd < esd.GetNumberOfPmdTracks(); ++iPmd) {\r
1100     // file pmd clusters, to be revised!\r
1101     AliESDPmdTrack *pmdTrack = esd.GetPmdTrack(iPmd);\r
1102     Int_t nLabel = 0;\r
1103     Int_t *label = 0x0;\r
1104     Double_t posPmd[3] = { pmdTrack->GetClusterX(), pmdTrack->GetClusterY(), pmdTrack->GetClusterZ()};\r
1105     Double_t pidPmd[13] = { 0.}; // to be revised!\r
1106     // type not set!\r
1107     // assoc cluster not set\r
1108     new(pmdClusters[jPmdClusters++]) AliAODPmdCluster(iPmd, nLabel, label, pmdTrack->GetClusterADC(), posPmd, pidPmd);\r
1109   }\r
1110 }\r
1111 \r
1112 //______________________________________________________________________________\r
1113 void AliAnalysisTaskESDfilter::ConvertCaloClusters(const AliESDEvent& esd)\r
1114 {\r
1115   AliCodeTimerAuto("",0);\r
1116   \r
1117   // Access to the AOD container of clusters\r
1118   TClonesArray &caloClusters = *(AODEvent()->GetCaloClusters());\r
1119   Int_t jClusters(0);\r
1120   \r
1121   for (Int_t iClust=0; iClust<esd.GetNumberOfCaloClusters(); ++iClust) {\r
1122     \r
1123     AliESDCaloCluster * cluster = esd.GetCaloCluster(iClust);\r
1124     \r
1125     Int_t  id        = cluster->GetID();\r
1126     Int_t  nLabel    = cluster->GetNLabels();\r
1127     Int_t *labels    = cluster->GetLabels();\r
1128     if(labels){ \r
1129                   for(int i = 0;i < nLabel;++i){\r
1130                           if(fMChandler)fMChandler->SelectParticle(labels[i]);\r
1131                   }\r
1132           }             \r
1133     \r
1134     Float_t energy = cluster->E();\r
1135     Float_t posF[3] = { 0.};\r
1136     cluster->GetPosition(posF);\r
1137     \r
1138     AliAODCaloCluster *caloCluster = new(caloClusters[jClusters++]) AliAODCaloCluster(id,\r
1139                                                                                       nLabel,\r
1140                                                                                       labels,\r
1141                                                                                       energy,\r
1142                                                                                       posF,\r
1143                                                                                       NULL,\r
1144                                                                                       cluster->GetType(),0);\r
1145     \r
1146     caloCluster->SetCaloCluster(cluster->GetDistanceToBadChannel(),\r
1147                                 cluster->GetDispersion(),\r
1148                                 cluster->GetM20(), cluster->GetM02(),\r
1149                                 cluster->GetEmcCpvDistance(),  \r
1150                                 cluster->GetNExMax(),cluster->GetTOF()) ;\r
1151     \r
1152     caloCluster->SetPIDFromESD(cluster->GetPID());\r
1153     caloCluster->SetNCells(cluster->GetNCells());\r
1154     caloCluster->SetCellsAbsId(cluster->GetCellsAbsId());\r
1155     caloCluster->SetCellsAmplitudeFraction(cluster->GetCellsAmplitudeFraction());\r
1156     \r
1157     TArrayI* matchedT =         cluster->GetTracksMatched();\r
1158     if (fNumberOfTracks>0 && matchedT && cluster->GetTrackMatchedIndex() >= 0) {        \r
1159       for (Int_t im = 0; im < matchedT->GetSize(); im++) {\r
1160         Int_t iESDtrack = matchedT->At(im);;\r
1161         if (fAODTrackRefs->At(iESDtrack) != 0) {\r
1162           caloCluster->AddTrackMatched((AliAODTrack*)fAODTrackRefs->At(iESDtrack));\r
1163         }\r
1164       }\r
1165     }\r
1166     \r
1167   } \r
1168   caloClusters.Expand(jClusters); // resize TObjArray to 'remove' slots for pseudo clusters      \r
1169 }\r
1170 \r
1171 //______________________________________________________________________________\r
1172 void AliAnalysisTaskESDfilter::ConvertEMCALCells(const AliESDEvent& esd)\r
1173 {\r
1174   AliCodeTimerAuto("",0);\r
1175   // fill EMCAL cell info\r
1176   if (esd.GetEMCALCells()) { // protection against missing ESD information\r
1177     AliESDCaloCells &esdEMcells = *(esd.GetEMCALCells());\r
1178     Int_t nEMcell = esdEMcells.GetNumberOfCells() ;\r
1179     \r
1180     AliAODCaloCells &aodEMcells = *(AODEvent()->GetEMCALCells());\r
1181     aodEMcells.CreateContainer(nEMcell);\r
1182     aodEMcells.SetType(AliAODCaloCells::kEMCALCell);\r
1183     for (Int_t iCell = 0; iCell < nEMcell; iCell++) {      \r
1184       aodEMcells.SetCell(iCell,esdEMcells.GetCellNumber(iCell),esdEMcells.GetAmplitude(iCell));\r
1185     }\r
1186     aodEMcells.Sort();\r
1187   }\r
1188 }\r
1189 \r
1190 //______________________________________________________________________________\r
1191 void AliAnalysisTaskESDfilter::ConvertPHOSCells(const AliESDEvent& esd)\r
1192 {\r
1193   AliCodeTimerAuto("",0);\r
1194   // fill PHOS cell info\r
1195   if (esd.GetPHOSCells()) { // protection against missing ESD information\r
1196     AliESDCaloCells &esdPHcells = *(esd.GetPHOSCells());\r
1197     Int_t nPHcell = esdPHcells.GetNumberOfCells() ;\r
1198     \r
1199     AliAODCaloCells &aodPHcells = *(AODEvent()->GetPHOSCells());\r
1200     aodPHcells.CreateContainer(nPHcell);\r
1201     aodPHcells.SetType(AliAODCaloCells::kPHOSCell);\r
1202     for (Int_t iCell = 0; iCell < nPHcell; iCell++) {      \r
1203       aodPHcells.SetCell(iCell,esdPHcells.GetCellNumber(iCell),esdPHcells.GetAmplitude(iCell));\r
1204     }\r
1205     aodPHcells.Sort();\r
1206   }\r
1207 }\r
1208 \r
1209 //______________________________________________________________________________\r
1210 void AliAnalysisTaskESDfilter::ConvertTracklets(const AliESDEvent& esd)\r
1211 {\r
1212   // tracklets    \r
1213   AliCodeTimerAuto("",0);\r
1214 \r
1215   AliAODTracklets &SPDTracklets = *(AODEvent()->GetTracklets());\r
1216   const AliMultiplicity *mult = esd.GetMultiplicity();\r
1217   if (mult) {\r
1218     if (mult->GetNumberOfTracklets()>0) {\r
1219       SPDTracklets.CreateContainer(mult->GetNumberOfTracklets());\r
1220       \r
1221       for (Int_t n=0; n<mult->GetNumberOfTracklets(); n++) {\r
1222         if(fMChandler){\r
1223           fMChandler->SelectParticle(mult->GetLabel(n, 0));\r
1224           fMChandler->SelectParticle(mult->GetLabel(n, 1));\r
1225         }\r
1226         SPDTracklets.SetTracklet(n, mult->GetTheta(n), mult->GetPhi(n), mult->GetDeltaPhi(n), mult->GetLabel(n, 0),mult->GetLabel(n, 1));\r
1227       }\r
1228     }\r
1229   } else {\r
1230     //Printf("ERROR: AliMultiplicity could not be retrieved from ESD");\r
1231   }\r
1232 }\r
1233 \r
1234 //______________________________________________________________________________\r
1235 void AliAnalysisTaskESDfilter::ConvertKinks(const AliESDEvent& esd)\r
1236 {\r
1237   AliCodeTimerAuto("",0);\r
1238   \r
1239   // Kinks: it is a big mess the access to the information in the kinks\r
1240   // The loop is on the tracks in order to find the mother and daugther of each kink\r
1241   \r
1242   Double_t covTr[21]={0.};\r
1243   Double_t pid[10]={0.};\r
1244   AliAODPid* detpid(0x0);\r
1245   \r
1246   fNumberOfKinks = esd.GetNumberOfKinks();\r
1247 \r
1248   const AliESDVertex* vtx = esd.GetPrimaryVertex();\r
1249   \r
1250   for (Int_t iTrack=0; iTrack<esd.GetNumberOfTracks(); ++iTrack) \r
1251   {\r
1252     AliESDtrack * esdTrack = esd.GetTrack(iTrack);\r
1253     \r
1254     Int_t ikink = esdTrack->GetKinkIndex(0);\r
1255     \r
1256     if (ikink && fNumberOfKinks) {\r
1257             // Negative kink index: mother, positive: daughter\r
1258             \r
1259             // Search for the second track of the kink\r
1260             \r
1261             for (Int_t jTrack = iTrack+1; jTrack<esd.GetNumberOfTracks(); ++jTrack) {\r
1262         \r
1263         AliESDtrack * esdTrack1 = esd.GetTrack(jTrack);\r
1264         \r
1265         Int_t jkink = esdTrack1->GetKinkIndex(0);\r
1266         \r
1267         if ( TMath::Abs(ikink)==TMath::Abs(jkink) ) {\r
1268           \r
1269           // The two tracks are from the same kink\r
1270           \r
1271           if (fUsedKink[TMath::Abs(ikink)-1]) continue; // skip used kinks\r
1272           \r
1273           Int_t imother = -1;\r
1274           Int_t idaughter = -1;\r
1275           \r
1276           if (ikink<0 && jkink>0) {\r
1277             \r
1278             imother = iTrack;\r
1279             idaughter = jTrack;\r
1280           }\r
1281           else if (ikink>0 && jkink<0) {\r
1282             \r
1283             imother = jTrack;\r
1284             idaughter = iTrack;\r
1285           }\r
1286           else {\r
1287             //                  cerr << "Error: Wrong combination of kink indexes: "\r
1288             //                       << ikink << " " << jkink << endl;\r
1289             continue;\r
1290           }\r
1291           \r
1292           // Add the mother track if it passed primary track selection cuts\r
1293           \r
1294           AliAODTrack * mother = NULL;\r
1295           \r
1296           UInt_t selectInfo = 0;\r
1297           if (fTrackFilter) {\r
1298             selectInfo = fTrackFilter->IsSelected(esd.GetTrack(imother));\r
1299             if (!selectInfo) continue;\r
1300           }\r
1301           \r
1302           if (!fUsedTrack[imother]) {\r
1303             \r
1304             fUsedTrack[imother] = kTRUE;\r
1305             \r
1306             AliESDtrack *esdTrackM = esd.GetTrack(imother);\r
1307             Double_t p[3] = { 0. };\r
1308             Double_t pos[3] = { 0. };\r
1309             esdTrackM->GetPxPyPz(p);\r
1310             esdTrackM->GetXYZ(pos);\r
1311             esdTrackM->GetCovarianceXYZPxPyPz(covTr);\r
1312             esdTrackM->GetESDpid(pid);\r
1313             if(fMChandler)fMChandler->SelectParticle(esdTrackM->GetLabel());\r
1314             mother = \r
1315             new(Tracks()[fNumberOfTracks++]) AliAODTrack(esdTrackM->GetID(),\r
1316                                                esdTrackM->GetLabel(),\r
1317                                                p,\r
1318                                                kTRUE,\r
1319                                                pos,\r
1320                                                kFALSE,\r
1321                                                covTr, \r
1322                                                (Short_t)esdTrackM->GetSign(),\r
1323                                                esdTrackM->GetITSClusterMap(), \r
1324                                                pid,\r
1325                                                fPrimaryVertex,\r
1326                                                kTRUE, // check if this is right\r
1327                                                vtx->UsesTrack(esdTrack->GetID()),\r
1328                                                AliAODTrack::kPrimary,\r
1329                                                selectInfo);\r
1330             mother->SetTPCClusterMap(esdTrackM->GetTPCClusterMap());\r
1331             mother->SetTPCSharedMap (esdTrackM->GetTPCSharedMap());\r
1332             mother->SetChi2perNDF(Chi2perNDF(esdTrackM));\r
1333             mother->SetTPCPointsF(esdTrackM->GetTPCNclsF());\r
1334 \r
1335             fAODTrackRefs->AddAt(mother, imother);\r
1336             \r
1337             if (esdTrackM->GetSign() > 0) ++fNumberOfPositiveTracks;\r
1338             mother->SetFlags(esdTrackM->GetStatus());\r
1339             mother->ConvertAliPIDtoAODPID();\r
1340             fPrimaryVertex->AddDaughter(mother);\r
1341             mother->ConvertAliPIDtoAODPID();\r
1342             SetAODPID(esdTrackM,mother,detpid);\r
1343           }\r
1344           else {\r
1345             //                  cerr << "Error: event " << esd.GetEventNumberInFile() << " kink " << TMath::Abs(ikink)-1\r
1346             //                       << " track " << imother << " has already been used!" << endl;\r
1347           }\r
1348           \r
1349           // Add the kink vertex\r
1350           AliESDkink * kink = esd.GetKink(TMath::Abs(ikink)-1);\r
1351           \r
1352           AliAODVertex * vkink = \r
1353           new(Vertices()[fNumberOfVertices++]) AliAODVertex(kink->GetPosition(),\r
1354                                                   NULL,\r
1355                                                   0.,\r
1356                                                   mother,\r
1357                                                   esdTrack->GetID(),  // This is the track ID of the mother's track!\r
1358                                                   AliAODVertex::kKink);\r
1359           // Add the daughter track\r
1360           \r
1361           AliAODTrack * daughter = NULL;\r
1362           \r
1363           if (!fUsedTrack[idaughter]) {\r
1364             \r
1365             fUsedTrack[idaughter] = kTRUE;\r
1366             \r
1367             AliESDtrack *esdTrackD = esd.GetTrack(idaughter);\r
1368             Double_t p[3] = { 0. };\r
1369             Double_t pos[3] = { 0. };\r
1370 \r
1371             esdTrackD->GetPxPyPz(p);\r
1372             esdTrackD->GetXYZ(pos);\r
1373             esdTrackD->GetCovarianceXYZPxPyPz(covTr);\r
1374             esdTrackD->GetESDpid(pid);\r
1375             selectInfo = 0;\r
1376             if (fTrackFilter) selectInfo = fTrackFilter->IsSelected(esdTrackD);\r
1377             if(fMChandler)fMChandler->SelectParticle(esdTrackD->GetLabel());\r
1378             daughter = \r
1379             new(Tracks()[fNumberOfTracks++]) AliAODTrack(esdTrackD->GetID(),\r
1380                                                esdTrackD->GetLabel(),\r
1381                                                p,\r
1382                                                kTRUE,\r
1383                                                pos,\r
1384                                                kFALSE,\r
1385                                                covTr, \r
1386                                                (Short_t)esdTrackD->GetSign(),\r
1387                                                esdTrackD->GetITSClusterMap(), \r
1388                                                pid,\r
1389                                                vkink,\r
1390                                                kTRUE, // check if this is right\r
1391                                                vtx->UsesTrack(esdTrack->GetID()),\r
1392                                                AliAODTrack::kSecondary,\r
1393                                                selectInfo);\r
1394             daughter->SetTPCClusterMap(esdTrackD->GetTPCClusterMap());\r
1395             daughter->SetTPCSharedMap (esdTrackD->GetTPCSharedMap());\r
1396             daughter->SetTPCPointsF(esdTrackD->GetTPCNclsF());\r
1397             fAODTrackRefs->AddAt(daughter, idaughter);\r
1398             \r
1399             if (esdTrackD->GetSign() > 0) ++fNumberOfPositiveTracks;\r
1400             daughter->SetFlags(esdTrackD->GetStatus());\r
1401             daughter->ConvertAliPIDtoAODPID();\r
1402             vkink->AddDaughter(daughter);\r
1403             daughter->ConvertAliPIDtoAODPID();\r
1404             SetAODPID(esdTrackD,daughter,detpid);\r
1405           }\r
1406           else {\r
1407             //                  cerr << "Error: event " << esd.GetEventNumberInFile() << " kink " << TMath::Abs(ikink)-1\r
1408             //                       << " track " << idaughter << " has already been used!" << endl;\r
1409           }\r
1410         }\r
1411             }\r
1412     }      \r
1413   }\r
1414 }\r
1415 \r
1416 //______________________________________________________________________________\r
1417 void AliAnalysisTaskESDfilter::ConvertPrimaryVertices(const AliESDEvent& esd)\r
1418 {\r
1419   AliCodeTimerAuto("",0);\r
1420   \r
1421   // Access to the AOD container of vertices\r
1422   fNumberOfVertices = 0;\r
1423   \r
1424   Double_t pos[3] = { 0. };\r
1425   Double_t covVtx[6] = { 0. };\r
1426 \r
1427   // Add primary vertex. The primary tracks will be defined\r
1428   // after the loops on the composite objects (V0, cascades, kinks)\r
1429   const AliESDVertex *vtx = esd.GetPrimaryVertex();\r
1430   \r
1431   vtx->GetXYZ(pos); // position\r
1432   vtx->GetCovMatrix(covVtx); //covariance matrix\r
1433   \r
1434   fPrimaryVertex = new(Vertices()[fNumberOfVertices++])\r
1435   AliAODVertex(pos, covVtx, vtx->GetChi2toNDF(), NULL, -1, AliAODVertex::kPrimary);\r
1436   fPrimaryVertex->SetName(vtx->GetName());\r
1437   fPrimaryVertex->SetTitle(vtx->GetTitle());\r
1438   \r
1439   TString vtitle = vtx->GetTitle();\r
1440   if (!vtitle.Contains("VertexerTracks")) \r
1441     fPrimaryVertex->SetNContributors(vtx->GetNContributors());\r
1442   \r
1443   if (fDebug > 0) fPrimaryVertex->Print();  \r
1444   \r
1445   // Add SPD "main" vertex \r
1446   const AliESDVertex *vtxS = esd.GetPrimaryVertexSPD();\r
1447   vtxS->GetXYZ(pos); // position\r
1448   vtxS->GetCovMatrix(covVtx); //covariance matrix\r
1449   AliAODVertex * mVSPD = new(Vertices()[fNumberOfVertices++])\r
1450   AliAODVertex(pos, covVtx, vtxS->GetChi2toNDF(), NULL, -1, AliAODVertex::kMainSPD);\r
1451   mVSPD->SetName(vtxS->GetName());\r
1452   mVSPD->SetTitle(vtxS->GetTitle());\r
1453   mVSPD->SetNContributors(vtxS->GetNContributors()); \r
1454   \r
1455   // Add SPD pileup vertices\r
1456   for(Int_t iV=0; iV<esd.GetNumberOfPileupVerticesSPD(); ++iV)\r
1457   {\r
1458     const AliESDVertex *vtxP = esd.GetPileupVertexSPD(iV);\r
1459     vtxP->GetXYZ(pos); // position\r
1460     vtxP->GetCovMatrix(covVtx); //covariance matrix\r
1461     AliAODVertex * pVSPD =  new(Vertices()[fNumberOfVertices++])\r
1462     AliAODVertex(pos, covVtx, vtxP->GetChi2toNDF(), NULL, -1, AliAODVertex::kPileupSPD);\r
1463     pVSPD->SetName(vtxP->GetName());\r
1464     pVSPD->SetTitle(vtxP->GetTitle());\r
1465     pVSPD->SetNContributors(vtxP->GetNContributors()); \r
1466   }\r
1467   \r
1468   // Add TRK pileup vertices\r
1469   for(Int_t iV=0; iV<esd.GetNumberOfPileupVerticesTracks(); ++iV)\r
1470   {\r
1471     const AliESDVertex *vtxP = esd.GetPileupVertexTracks(iV);\r
1472     vtxP->GetXYZ(pos); // position\r
1473     vtxP->GetCovMatrix(covVtx); //covariance matrix\r
1474     AliAODVertex * pVTRK = new(Vertices()[fNumberOfVertices++])\r
1475     AliAODVertex(pos, covVtx, vtxP->GetChi2toNDF(), NULL, -1, AliAODVertex::kPileupTracks);\r
1476     pVTRK->SetName(vtxP->GetName());\r
1477     pVTRK->SetTitle(vtxP->GetTitle());\r
1478     pVTRK->SetNContributors(vtxP->GetNContributors());\r
1479   }\r
1480 }\r
1481 \r
1482 //______________________________________________________________________________\r
1483 void AliAnalysisTaskESDfilter::ConvertVZERO(const AliESDEvent& esd)\r
1484 {\r
1485   // Convert VZERO data\r
1486   AliAODVZERO* vzeroData = AODEvent()->GetVZEROData();\r
1487   *vzeroData = *(esd.GetVZEROData());\r
1488 }\r
1489 \r
1490 //______________________________________________________________________________\r
1491 void AliAnalysisTaskESDfilter::ConvertESDtoAOD() \r
1492 {\r
1493   // ESD Filter analysis task executed for each event\r
1494   \r
1495   AliESDEvent* esd = dynamic_cast<AliESDEvent*>(InputEvent());\r
1496   \r
1497   if(!esd)return;\r
1498 \r
1499   AliCodeTimerAuto("",0);\r
1500   \r
1501   fOldESDformat = ( esd->GetAliESDOld() != 0x0 );\r
1502   \r
1503   fNumberOfTracks = 0;\r
1504   fNumberOfPositiveTracks = 0;\r
1505   fNumberOfV0s = 0;\r
1506   fNumberOfVertices = 0;\r
1507   fNumberOfCascades = 0;\r
1508   fNumberOfKinks = 0;\r
1509     \r
1510   AliAODHeader* header = ConvertHeader(*esd);\r
1511 \r
1512   if ( fIsVZEROEnabled ) ConvertVZERO(*esd);\r
1513   \r
1514   // Fetch Stack for debuggging if available \r
1515   fMChandler=0x0;\r
1516   if(MCEvent())\r
1517   {\r
1518     fMChandler = (AliMCEventHandler*) ((AliAnalysisManager::GetAnalysisManager())->GetMCtruthEventHandler()); \r
1519   }\r
1520   \r
1521   // loop over events and fill them\r
1522   // Multiplicity information needed by the header (to be revised!)\r
1523   Int_t nTracks    = esd->GetNumberOfTracks();\r
1524   for (Int_t iTrack=0; iTrack<nTracks; ++iTrack) esd->GetTrack(iTrack)->SetESDEvent(esd);\r
1525 \r
1526   // Update the header\r
1527 \r
1528   Int_t nV0s      = esd->GetNumberOfV0s();\r
1529   Int_t nCascades = esd->GetNumberOfCascades();\r
1530   Int_t nKinks    = esd->GetNumberOfKinks();\r
1531   Int_t nVertices = nV0s + nCascades /*V0 wihtin cascade already counted*/+ nKinks + 1 /* = prim. vtx*/;\r
1532   Int_t nPileSPDVertices=1+esd->GetNumberOfPileupVerticesSPD(); // also SPD main vertex\r
1533   Int_t nPileTrkVertices=esd->GetNumberOfPileupVerticesTracks();\r
1534   nVertices+=nPileSPDVertices;\r
1535   nVertices+=nPileTrkVertices;\r
1536   Int_t nJets     = 0;\r
1537   Int_t nCaloClus = esd->GetNumberOfCaloClusters();\r
1538   Int_t nFmdClus  = 0;\r
1539   Int_t nPmdClus  = esd->GetNumberOfPmdTracks();\r
1540     \r
1541   AliDebug(1,Form("   NV0=%d  NCASCADES=%d  NKINKS=%d", nV0s, nCascades, nKinks));\r
1542        \r
1543   AODEvent()->ResetStd(nTracks, nVertices, nV0s, nCascades, nJets, nCaloClus, nFmdClus, nPmdClus);\r
1544 \r
1545   if (nV0s > 0) \r
1546   {\r
1547     // RefArray to store a mapping between esd V0 number and newly created AOD-Vertex V0\r
1548     fAODV0VtxRefs = new TRefArray(nV0s);\r
1549     // RefArray to store the mapping between esd V0 number and newly created AOD-V0\r
1550     fAODV0Refs = new TRefArray(nV0s); \r
1551     // Array to take into account the V0s already added to the AOD (V0 within cascades)\r
1552     fUsedV0 = new Bool_t[nV0s];\r
1553     for (Int_t iV0=0; iV0<nV0s; ++iV0) fUsedV0[iV0]=kFALSE;\r
1554   }\r
1555   \r
1556   if (nTracks>0) \r
1557   {\r
1558     // RefArray to store the mapping between esd track number and newly created AOD-Track\r
1559     \r
1560     fAODTrackRefs = new TRefArray(nTracks);\r
1561 \r
1562     // Array to take into account the tracks already added to the AOD    \r
1563     fUsedTrack = new Bool_t[nTracks];\r
1564     for (Int_t iTrack=0; iTrack<nTracks; ++iTrack) fUsedTrack[iTrack]=kFALSE;\r
1565   }\r
1566   \r
1567   // Array to take into account the kinks already added to the AOD\r
1568   if (nKinks>0) \r
1569   {\r
1570     fUsedKink = new Bool_t[nKinks];\r
1571     for (Int_t iKink=0; iKink<nKinks; ++iKink) fUsedKink[iKink]=kFALSE;\r
1572   }\r
1573     \r
1574   ConvertPrimaryVertices(*esd);\r
1575 \r
1576   //setting best TOF PID\r
1577   AliESDInputHandler* esdH = dynamic_cast<AliESDInputHandler*>(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());\r
1578   if (esdH)\r
1579       fESDpid = esdH->GetESDpid();\r
1580 \r
1581   if (fIsPidOwner && fESDpid){\r
1582     delete fESDpid;\r
1583     fESDpid = 0;\r
1584   }\r
1585   if(!fESDpid)\r
1586   { //in case of no Tender attached \r
1587     fESDpid = new AliESDpid;\r
1588     fIsPidOwner = kTRUE;\r
1589   }\r
1590   \r
1591   if(!esd->GetTOFHeader())\r
1592   { //protection in case the pass2 LHC10b,c,d have been processed without tender. \r
1593     Float_t t0spread[10];\r
1594     Float_t intrinsicTOFres=100; //ps ok for LHC10b,c,d pass2!! \r
1595     for (Int_t i=0; i<10; i++) t0spread[i] = (TMath::Sqrt(esd->GetSigma2DiamondZ()))/0.03; //0.03 to convert from cm to ps\r
1596     fESDpid->GetTOFResponse().SetT0resolution(t0spread);\r
1597     fESDpid->GetTOFResponse().SetTimeResolution(intrinsicTOFres);\r
1598     \r
1599     fESDpid->SetTOFResponse(esd, (AliESDpid::EStartTimeType_t)fTimeZeroType);    \r
1600   }\r
1601   \r
1602   if(esd->GetTOFHeader() && fIsPidOwner) fESDpid->SetTOFResponse(esd, (AliESDpid::EStartTimeType_t)fTimeZeroType); //in case of AOD production strating form LHC10e without Tender. \r
1603   \r
1604   if ( fAreCascadesEnabled ) ConvertCascades(*esd);\r
1605 \r
1606   if ( fAreV0sEnabled ) ConvertV0s(*esd);\r
1607   \r
1608   if ( fAreKinksEnabled ) ConvertKinks(*esd);\r
1609   \r
1610   if ( fAreTracksEnabled ) ConvertTracks(*esd);\r
1611   \r
1612   // Update number of AOD tracks in header at the end of track loop (M.G.)\r
1613   header->SetRefMultiplicity(fNumberOfTracks);\r
1614   header->SetRefMultiplicityPos(fNumberOfPositiveTracks);\r
1615   header->SetRefMultiplicityNeg(fNumberOfTracks - fNumberOfPositiveTracks);\r
1616 \r
1617   if ( fTPCOnlyFilterMask ) ConvertTPCOnlyTracks(*esd);\r
1618   \r
1619   if ( fArePmdClustersEnabled ) ConvertPmdClusters(*esd);\r
1620   \r
1621   if ( fAreCaloClustersEnabled ) ConvertCaloClusters(*esd);\r
1622   \r
1623   if ( fAreEMCALCellsEnabled )ConvertEMCALCells(*esd);\r
1624   \r
1625   if ( fArePHOSCellsEnabled )ConvertPHOSCells(*esd);\r
1626   \r
1627   if ( fAreTrackletsEnabled ) ConvertTracklets(*esd);\r
1628   \r
1629   delete fAODTrackRefs; fAODTrackRefs=0x0;\r
1630   delete fAODV0VtxRefs; fAODV0VtxRefs=0x0;\r
1631   delete fAODV0Refs; fAODV0Refs=0x0;\r
1632   \r
1633   delete[] fUsedTrack; fUsedTrack=0x0;\r
1634   delete[] fUsedV0; fUsedV0=0x0;\r
1635   delete[] fUsedKink; fUsedKink=0x0;\r
1636 \r
1637   if ( fIsPidOwner){\r
1638     delete fESDpid;\r
1639     fESDpid = 0x0;\r
1640   }\r
1641 \r
1642 \r
1643 }\r
1644 \r
1645 \r
1646 //______________________________________________________________________________\r
1647 void AliAnalysisTaskESDfilter::SetAODPID(AliESDtrack *esdtrack, AliAODTrack *aodtrack, AliAODPid *detpid)\r
1648 {\r
1649   //\r
1650   // Setter for the raw PID detector signals\r
1651   //\r
1652 \r
1653   // Save PID object for candidate electrons\r
1654     Bool_t pidSave = kFALSE;\r
1655     if (fTrackFilter) {\r
1656         Bool_t selectInfo = fTrackFilter->IsSelected((char*) "Electrons");\r
1657         if (selectInfo)  pidSave = kTRUE;\r
1658     }\r
1659 \r
1660 \r
1661     // Tracks passing pt cut \r
1662     if(esdtrack->Pt()>fHighPthreshold) {\r
1663         pidSave = kTRUE;\r
1664     } else {\r
1665         if(fPtshape){\r
1666             if(esdtrack->Pt()> fPtshape->GetXmin()){\r
1667                 Double_t y = fPtshape->Eval(esdtrack->Pt())/fPtshape->Eval(fHighPthreshold);\r
1668                 if(gRandom->Rndm(0)<1./y){\r
1669                     pidSave = kTRUE;\r
1670                 }//end rndm\r
1671             }//end if p < pmin\r
1672         }//end if p function\r
1673     }// end else\r
1674 \r
1675     if (pidSave) {\r
1676       if(!aodtrack->GetDetPid()){// prevent memory leak when calling SetAODPID twice for the same track\r
1677         detpid = new AliAODPid();\r
1678         SetDetectorRawSignals(detpid,esdtrack);\r
1679         aodtrack->SetDetPID(detpid);\r
1680       }\r
1681     }\r
1682 }\r
1683 \r
1684 //______________________________________________________________________________\r
1685 void AliAnalysisTaskESDfilter::SetDetectorRawSignals(AliAODPid *aodpid, AliESDtrack *track)\r
1686 {\r
1687 //\r
1688 //assignment of the detector signals (AliXXXesdPID inspired)\r
1689 //\r
1690  if(!track){\r
1691  AliInfo("no ESD track found. .....exiting");\r
1692  return;\r
1693  }\r
1694  // TPC momentum\r
1695  const AliExternalTrackParam *in=track->GetInnerParam();\r
1696  if (in) {\r
1697    aodpid->SetTPCmomentum(in->GetP());\r
1698  }else{\r
1699    aodpid->SetTPCmomentum(-1.);\r
1700  }\r
1701 \r
1702 \r
1703  aodpid->SetITSsignal(track->GetITSsignal());\r
1704  aodpid->SetTPCsignal(track->GetTPCsignal());\r
1705  aodpid->SetTPCsignalN(track->GetTPCsignalN());\r
1706 \r
1707  //n TRD planes = 6\r
1708  Int_t nslices = track->GetNumberOfTRDslices()*6;\r
1709  Double_t *trdslices = new Double_t[nslices];\r
1710  for(Int_t iSl =0; iSl < track->GetNumberOfTRDslices(); iSl++) {\r
1711    for(Int_t iPl =0; iPl<6; iPl++) trdslices[iPl*track->GetNumberOfTRDslices()+iSl] = track->GetTRDslice(iPl,iSl);\r
1712  }\r
1713  \r
1714 //TRD momentum\r
1715  for(Int_t iPl=0;iPl<6;iPl++){\r
1716    Double_t trdmom=track->GetTRDmomentum(iPl);\r
1717    aodpid->SetTRDmomentum(iPl,trdmom);\r
1718  }\r
1719 \r
1720  aodpid->SetTRDsignal(track->GetNumberOfTRDslices()*6,trdslices);\r
1721 \r
1722  //TOF PID  \r
1723  Double_t times[AliAODPid::kSPECIES]; track->GetIntegratedTimes(times);\r
1724  aodpid->SetIntegratedTimes(times);\r
1725 \r
1726   Float_t tzeroTrack = fESDpid->GetTOFResponse().GetStartTime(track->P());\r
1727   aodpid->SetTOFsignal(track->GetTOFsignal()-tzeroTrack);\r
1728   \r
1729   Double_t tofRes[5];\r
1730   for (Int_t iMass=0; iMass<5; iMass++){\r
1731     tofRes[iMass]=(Double_t)fESDpid->GetTOFResponse().GetExpectedSigma(track->P(), times[iMass], AliPID::ParticleMass(iMass));\r
1732   }\r
1733   aodpid->SetTOFpidResolution(tofRes);\r
1734 \r
1735   aodpid->SetHMPIDsignal(track->GetHMPIDsignal());\r
1736 \r
1737  //Extrapolate track to EMCAL surface for AOD-level track-cluster matching\r
1738  Double_t emcpos[3] = {0.,0.,0.};\r
1739  Double_t emcmom[3] = {0.,0.,0.};\r
1740  aodpid->SetEMCALPosition(emcpos);\r
1741  aodpid->SetEMCALMomentum(emcmom);\r
1742 \r
1743  AliExternalTrackParam *outerparam = (AliExternalTrackParam*)track->GetOuterParam();\r
1744  if(!outerparam) return;\r
1745 \r
1746  //To be replaced by call to AliEMCALGeoUtils when the class becomes available\r
1747  Bool_t okpos = outerparam->GetXYZ(emcpos);\r
1748  Bool_t okmom = outerparam->GetPxPyPz(emcmom);\r
1749  if(!(okpos && okmom)) return;\r
1750 \r
1751  aodpid->SetEMCALPosition(emcpos);\r
1752  aodpid->SetEMCALMomentum(emcmom);\r
1753 \r
1754 }\r
1755 \r
1756 Double_t  AliAnalysisTaskESDfilter::Chi2perNDF(AliESDtrack* track)\r
1757 {\r
1758     // Calculate chi2 per ndf for track\r
1759     Int_t  nClustersTPC = track->GetTPCNcls();\r
1760 \r
1761     if ( nClustersTPC > 5) {\r
1762        return (track->GetTPCchi2()/Float_t(nClustersTPC - 5));\r
1763     } else {\r
1764        return (-1.);\r
1765     }\r
1766  }\r
1767 \r
1768 \r
1769 //______________________________________________________________________________\r
1770 void AliAnalysisTaskESDfilter::Terminate(Option_t */*option*/)\r
1771 {\r
1772 // Terminate analysis\r
1773 //\r
1774     if (fDebug > 1) printf("AnalysisESDfilter: Terminate() \n");\r
1775 }\r
1776 \r
1777 //______________________________________________________________________________\r
1778 void  AliAnalysisTaskESDfilter::PrintMCInfo(AliStack *pStack,Int_t label){\r
1779 // Print MC info\r
1780   if(!pStack)return;\r
1781   label = TMath::Abs(label);\r
1782   TParticle *part = pStack->Particle(label);\r
1783   Printf("########################");\r
1784   Printf("%s:%d %d UniqueID %d PDG %d P %3.3f",(char*)__FILE__,__LINE__,label,part->GetUniqueID(),part->GetPdgCode(),part->P());\r
1785   part->Print();\r
1786   TParticle* mother = part;\r
1787   Int_t imo = part->GetFirstMother();\r
1788   Int_t nprim = pStack->GetNprimary();\r
1789   //  while((imo >= nprim) && (mother->GetUniqueID() == 4)) {\r
1790   while((imo >= nprim)) {\r
1791     mother =  pStack->Particle(imo);\r
1792     Printf("Mother %s:%d Label %d UniqueID %d PDG %d P %3.3f",(char*)__FILE__,__LINE__,imo,mother->GetUniqueID(),mother->GetPdgCode(),mother->P());\r
1793     mother->Print();\r
1794     imo =  mother->GetFirstMother();\r
1795   }\r
1796   Printf("########################");\r
1797 }\r