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