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