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