ANALYSIS/AliAnalysisTaskESDfilter.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
   3  *                                                                        *
   4  * Author: The ALICE Off-line Project.                                    *
   5  * Contributors are mentioned in the code where appropriate.              *
   6  *                                                                        *
   7  * Permission to use, copy, modify and distribute this software and its   *
   8  * documentation strictly for non-commercial purposes is hereby granted   *
   9  * without fee, provided that the above copyright notice appears in all   *
  10  * copies and that both the copyright notice and this permission notice   *
  11  * appear in the supporting documentation. The authors make no claims     *
  12  * about the suitability of this software for any purpose. It is          *
  13  * provided "as is" without express or implied warranty.                  *
  14  **************************************************************************/
  15
  16 /* $Id$ */
  17
  18 #include <TChain.h>
  19 #include <TTree.h>
  20 #include <TList.h>
  21 #include <TArrayI.h>
  22 #include <TParameter.h>
  23 #include <TRandom.h>
  24 #include <TParticle.h>
  25 #include <TFile.h>
  26
  27 #include "AliAnalysisTaskESDfilter.h"
  28 #include "AliAnalysisManager.h"
  29 #include "AliESDEvent.h"
  30 #include "AliESDRun.h"
  31 #include "AliStack.h"
  32 #include "AliAODEvent.h"
  33 #include "AliMCEvent.h"
  34 #include "AliMCEventHandler.h"
  35 #include "AliESDInputHandler.h"
  36 #include "AliAODHandler.h"
  37 #include "AliAODMCParticle.h"
  38 #include "AliAnalysisFilter.h"
  39 #include "AliESDMuonTrack.h"
  40 #include "AliESDVertex.h"
  41 #include "AliCentrality.h"
  42 #include "AliEventplane.h"
  43 #include "AliESDv0.h"
  44 #include "AliESDkink.h"
  45 #include "AliESDcascade.h"
  46 #include "AliESDPmdTrack.h"
  47 #include "AliESDCaloCluster.h"
  48 #include "AliESDCaloCells.h"
  49 #include "AliMultiplicity.h"
  50 #include "AliLog.h"
  51 #include "AliCodeTimer.h"
  52 #include "AliESDtrackCuts.h"
  53 #include "AliESDpid.h"
  54 #include "AliAODHMPIDrings.h"
  55 #include "AliV0vertexer.h"
  56 #include "AliCascadeVertexer.h"
  57 #include "Riostream.h"
  58 #include "AliExternalTrackParam.h"
  59 #include "AliTrackerBase.h"
  60 #include "TVector3.h"
  61 #include "AliTPCdEdxInfo.h"
  62
  63 using std::cout;
  64 using std::endl;
  65 ClassImp(AliAnalysisTaskESDfilter)
  66
  67 ////////////////////////////////////////////////////////////////////////
  68
  69 AliAnalysisTaskESDfilter::AliAnalysisTaskESDfilter():
  70   AliAnalysisTaskSE(),
  71   fTrackFilter(0x0),
  72   fKinkFilter(0x0),
  73   fV0Filter(0x0),
  74   fCascadeFilter(0x0),
  75   fHighPthreshold(0),
  76   fPtshape(0x0),
  77   fEnableFillAOD(kTRUE),
  78   fUsedTrack(0x0),
  79   fUsedKink(0x0),
  80   fUsedV0(0x0),
  81   fAODTrackRefs(0x0),
  82   fAODV0VtxRefs(0x0),
  83   fAODV0Refs(0x0),
  84   fMChandler(0x0),
  85   fNumberOfTracks(0),
  86   fNumberOfPositiveTracks(0),
  87   fNumberOfV0s(0),
  88   fNumberOfVertices(0),
  89   fNumberOfCascades(0),
  90   fNumberOfKinks(0),
  91   fOldESDformat(kFALSE),
  92   fPrimaryVertex(0x0),
  93   fTPCConstrainedFilterMask(0),
  94   fHybridFilterMaskTPCCG(0),
  95   fWriteHybridTPCCOnly(kFALSE),
  96   fGlobalConstrainedFilterMask(0),
  97   fHybridFilterMaskGCG(0),
  98   fWriteHybridGCOnly(kFALSE),
  99   fIsVZEROEnabled(kTRUE),
 100   fIsTZEROEnabled(kTRUE),
 101   fIsZDCEnabled(kTRUE),
 102   fIsHMPIDEnabled(kTRUE),
 103   fIsV0CascadeRecoEnabled(kFALSE),
 104   fAreCascadesEnabled(kTRUE),
 105   fAreV0sEnabled(kTRUE),
 106   fAreKinksEnabled(kTRUE),
 107   fAreTracksEnabled(kTRUE),
 108   fArePmdClustersEnabled(kTRUE),
 109   fAreCaloClustersEnabled(kTRUE),
 110   fAreEMCALCellsEnabled(kTRUE),
 111   fArePHOSCellsEnabled(kTRUE),
 112   fAreEMCALTriggerEnabled(kTRUE),
 113   fArePHOSTriggerEnabled(kTRUE),
 114   fAreTrackletsEnabled(kTRUE),
 115   fESDpid(0x0),
 116   fIsPidOwner(kFALSE),
 117   fTPCaloneTrackCuts(0),
 118   fDoPropagateTrackToEMCal(kTRUE),
 119   fEMCalSurfaceDistance(440)
 120 {
 121   // Default constructor
 122     fV0Cuts[0] =  33.   ;   // max allowed chi2
 123     fV0Cuts[1] =   0.1  ;   // min allowed impact parameter for the 1st daughter
 124     fV0Cuts[2] =   0.1  ;   // min allowed impact parameter for the 2nd daughter
 125     fV0Cuts[3] =   1.   ;   // max allowed DCA between the daughter tracks
 126     fV0Cuts[4] =    .998;   // min allowed cosine of V0's pointing angle
 127     fV0Cuts[5] =   0.9  ;   // min radius of the fiducial volume
 128     fV0Cuts[6] = 100.   ;   // max radius of the fiducial volume
 129
 130     fCascadeCuts[0] =  33.   ; // max allowed chi2 (same as PDC07)
 131     fCascadeCuts[1] =   0.05 ; // min allowed V0 impact parameter
 132     fCascadeCuts[2] =   0.008; // "window" around the Lambda mass
 133     fCascadeCuts[3] =   0.03 ; // min allowed bachelor's impact parameter
 134     fCascadeCuts[4] =   0.3  ; // max allowed DCA between the V0 and the bachelor
 135     fCascadeCuts[5] =   0.999; // min allowed cosine of the cascade pointing angle
 136     fCascadeCuts[6] =   0.9  ; // min radius of the fiducial volume
 137     fCascadeCuts[7] = 100.   ; // max radius of the fiducial volume
 138 }
 139
 140 //______________________________________________________________________________
 141 AliAnalysisTaskESDfilter::AliAnalysisTaskESDfilter(const char* name):
 142     AliAnalysisTaskSE(name),
 143     fTrackFilter(0x0),
 144     fKinkFilter(0x0),
 145     fV0Filter(0x0),
 146     fCascadeFilter(0x0),
 147     fHighPthreshold(0),
 148     fPtshape(0x0),
 149     fEnableFillAOD(kTRUE),
 150     fUsedTrack(0x0),
 151     fUsedKink(0x0),
 152     fUsedV0(0x0),
 153     fAODTrackRefs(0x0),
 154     fAODV0VtxRefs(0x0),
 155     fAODV0Refs(0x0),
 156     fMChandler(0x0),
 157     fNumberOfTracks(0),
 158     fNumberOfPositiveTracks(0),
 159     fNumberOfV0s(0),
 160     fNumberOfVertices(0),
 161     fNumberOfCascades(0),
 162     fNumberOfKinks(0),
 163     fOldESDformat(kFALSE),
 164     fPrimaryVertex(0x0),
 165   fTPCConstrainedFilterMask(0),
 166   fHybridFilterMaskTPCCG(0),
 167   fWriteHybridTPCCOnly(kFALSE),
 168   fGlobalConstrainedFilterMask(0),
 169   fHybridFilterMaskGCG(0),
 170   fWriteHybridGCOnly(kFALSE),
 171     fIsVZEROEnabled(kTRUE),
 172     fIsTZEROEnabled(kTRUE),
 173     fIsZDCEnabled(kTRUE),
 174     fIsHMPIDEnabled(kTRUE),
 175     fIsV0CascadeRecoEnabled(kFALSE),
 176     fAreCascadesEnabled(kTRUE),
 177     fAreV0sEnabled(kTRUE),
 178     fAreKinksEnabled(kTRUE),
 179     fAreTracksEnabled(kTRUE),
 180     fArePmdClustersEnabled(kTRUE),
 181     fAreCaloClustersEnabled(kTRUE),
 182     fAreEMCALCellsEnabled(kTRUE),
 183     fArePHOSCellsEnabled(kTRUE),
 184                 fAreEMCALTriggerEnabled(kTRUE),
 185                 fArePHOSTriggerEnabled(kTRUE),
 186                 fAreTrackletsEnabled(kTRUE),
 187     fESDpid(0x0),
 188     fIsPidOwner(kFALSE),
 189     fTPCaloneTrackCuts(0),
 190     fDoPropagateTrackToEMCal(kTRUE),
 191     fEMCalSurfaceDistance(440)
 192 {
 193   // Constructor
 194
 195     fV0Cuts[0] =  33.   ;   // max allowed chi2
 196     fV0Cuts[1] =   0.1  ;   // min allowed impact parameter for the 1st daughter
 197     fV0Cuts[2] =   0.1  ;   // min allowed impact parameter for the 2nd daughter
 198     fV0Cuts[3] =   1.   ;   // max allowed DCA between the daughter tracks
 199     fV0Cuts[4] =    .998;   // min allowed cosine of V0's pointing angle
 200     fV0Cuts[5] =   0.9  ;   // min radius of the fiducial volume
 201     fV0Cuts[6] = 100.   ;   // max radius of the fiducial volume
 202
 203     fCascadeCuts[0] =  33.   ; // max allowed chi2 (same as PDC07)
 204     fCascadeCuts[1] =   0.05 ; // min allowed V0 impact parameter
 205     fCascadeCuts[2] =   0.008; // "window" around the Lambda mass
 206     fCascadeCuts[3] =   0.03 ; // min allowed bachelor's impact parameter
 207     fCascadeCuts[4] =   0.3  ; // max allowed DCA between the V0 and the bachelor
 208     fCascadeCuts[5] =   0.999; // min allowed cosine of the cascade pointing angle
 209     fCascadeCuts[6] =   0.9  ; // min radius of the fiducial volume
 210     fCascadeCuts[7] = 100.   ; // max radius of the fiducial volume
 211
 212
 213
 214 }
 215 AliAnalysisTaskESDfilter::~AliAnalysisTaskESDfilter(){
 216     if(fIsPidOwner) delete fESDpid;
 217 }
 218 //______________________________________________________________________________
 219 void AliAnalysisTaskESDfilter::UserCreateOutputObjects()
 220 {
 221   //
 222   // Create Output Objects conenct filter to outputtree
 223   //
 224   if(OutputTree())
 225   {
 226     OutputTree()->GetUserInfo()->Add(fTrackFilter);
 227   }
 228   else
 229   {
 230     AliError("No OutputTree() for adding the track filter");
 231   }
 232   fTPCaloneTrackCuts = AliESDtrackCuts::GetStandardTPCOnlyTrackCuts();
 233 }
 234
 235 //______________________________________________________________________________
 236 void AliAnalysisTaskESDfilter::Init()
 237 {
 238   // Initialization
 239   if (fDebug > 1) AliInfo("Init() \n");
 240   // Call configuration file
 241 }
 242
 243 //______________________________________________________________________________
 244 Bool_t AliAnalysisTaskESDfilter::Notify()
 245 {
 246 // Notify method.
 247    AddMetadataToUserInfo();
 248    return kTRUE;
 249 }
 250
 251 //______________________________________________________________________________
 252 Bool_t AliAnalysisTaskESDfilter::AddMetadataToUserInfo()
 253 {
 254 // Copy metadata to AOD user info.
 255   static Bool_t copyFirst = kFALSE;
 256   if (!copyFirst) {
 257     AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
 258     if (!mgr) {
 259        AliError("AliAnalysisTaskESDfilter::AddMetadataToUserInfo() : No analysis manager !");
 260        return kFALSE;
 261     }
 262     TTree *esdTree = mgr->GetTree()->GetTree();
 263     if (!esdTree) return kFALSE;
 264     TNamed *alirootVersion = (TNamed*)esdTree->GetUserInfo()->FindObject("alirootVersion");
 265     if (!alirootVersion) return kFALSE;
 266     AliAODHandler *aodHandler = dynamic_cast<AliAODHandler*>(mgr->GetOutputEventHandler());
 267          if (!aodHandler) return kFALSE;
 268     TTree *aodTree = aodHandler->GetTree();
 269     if (!aodTree) return kFALSE;
 270     aodTree->GetUserInfo()->Add(new TNamed(*alirootVersion));
 271     copyFirst = kTRUE;
 272   }
 273   return kTRUE;
 274 }
 275
 276 //______________________________________________________________________________
 277 void AliAnalysisTaskESDfilter::PrintTask(Option_t *option, Int_t indent) const
 278 {
 279 // Print selection task information
 280   AliInfo("");
 281
 282   AliAnalysisTaskSE::PrintTask(option,indent);
 283
 284   TString spaces(' ',indent+3);
 285
 286         cout << spaces.Data() << Form("Cascades       are %s",fAreCascadesEnabled ? "ENABLED":"DISABLED") << endl;
 287         cout << spaces.Data() << Form("V0s            are %s",fAreV0sEnabled ? "ENABLED":"DISABLED") << endl;
 288         cout << spaces.Data() << Form("Kinks          are %s",fAreKinksEnabled ? "ENABLED":"DISABLED") << endl;
 289         cout << spaces.Data() << Form("Tracks         are %s",fAreTracksEnabled ? "ENABLED":"DISABLED") << endl;
 290         cout << spaces.Data() << Form("PmdClusters    are %s",fArePmdClustersEnabled ? "ENABLED":"DISABLED") << endl;
 291         cout << spaces.Data() << Form("CaloClusters   are %s",fAreCaloClustersEnabled ? "ENABLED":"DISABLED") << endl;
 292   cout << spaces.Data() << Form("EMCAL cells    are %s",fAreEMCALCellsEnabled ? "ENABLED":"DISABLED") << endl;
 293         cout << spaces.Data() << Form("EMCAL triggers are %s",fAreEMCALTriggerEnabled ? "ENABLED":"DISABLED") << endl;
 294         cout << spaces.Data() << Form("PHOS triggers  are %s",fArePHOSTriggerEnabled ? "ENABLED":"DISABLED") << endl;
 295         cout << spaces.Data() << Form("Tracklets      are %s",fAreTrackletsEnabled ? "ENABLED":"DISABLED") << endl;
 296         cout << spaces.Data() << Form("PropagateTrackToEMCal  is %s", fDoPropagateTrackToEMCal ? "ENABLED":"DISABLED") << endl;
 297 }
 298
 299 //______________________________________________________________________________
 300 void AliAnalysisTaskESDfilter::UserExec(Option_t */*option*/)
 301 {
 302 // Execute analysis for current event
 303 //
 304
 305   Long64_t ientry = Entry();
 306
 307   if (fDebug > 0) {
 308       printf("Filter: Analysing event # %5d\n", (Int_t) ientry);
 309       if (fHighPthreshold == 0) AliInfo("detector PID signals are stored in each track");
 310       if (!fPtshape) AliInfo("detector PID signals are not stored below the pt threshold");
 311   }
 312   // Filters must explicitely enable AOD filling in their UserExec (AG)
 313   if (!AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler()) AliFatal("Cannot run ESD filter without an output event handler");
 314   if(fEnableFillAOD) {
 315      AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler()->SetFillAOD(kTRUE);
 316      AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler()->SetFillExtension(kTRUE);
 317   }
 318   ConvertESDtoAOD();
 319 }
 320
 321 //______________________________________________________________________________
 322 TClonesArray& AliAnalysisTaskESDfilter::Cascades()
 323 {
 324   return *(AODEvent()->GetCascades());
 325 }
 326
 327 //______________________________________________________________________________
 328 TClonesArray& AliAnalysisTaskESDfilter::Tracks()
 329 {
 330   return *(AODEvent()->GetTracks());
 331 }
 332
 333 //______________________________________________________________________________
 334 TClonesArray& AliAnalysisTaskESDfilter::V0s()
 335 {
 336   return *(AODEvent()->GetV0s());
 337 }
 338
 339 //______________________________________________________________________________
 340 TClonesArray& AliAnalysisTaskESDfilter::Vertices()
 341 {
 342   return *(AODEvent()->GetVertices());
 343 }
 344
 345 //______________________________________________________________________________
 346 AliAODHeader* AliAnalysisTaskESDfilter::ConvertHeader(const AliESDEvent& esd)
 347 {
 348 // Convert header information
 349
 350   AliCodeTimerAuto("",0);
 351
 352   AliAODHeader* header = AODEvent()->GetHeader();
 353
 354   header->SetRunNumber(esd.GetRunNumber());
 355   header->SetOfflineTrigger(fInputHandler->IsEventSelected()); // propagate the decision of the physics selection
 356   header->SetNumberOfESDTracks(esd.GetNumberOfTracks());
 357
 358   TTree* tree = fInputHandler->GetTree();
 359   if (tree) {
 360     TFile* file = tree->GetCurrentFile();
 361     if (file) header->SetESDFileName(file->GetName());
 362   }
 363
 364   if (fOldESDformat) {
 365     header->SetBunchCrossNumber(0);
 366     header->SetOrbitNumber(0);
 367     header->SetPeriodNumber(0);
 368     header->SetEventType(0);
 369     header->SetMuonMagFieldScale(-999.);
 370     header->SetCentrality(0);
 371     header->SetEventplane(0);
 372   } else {
 373     header->SetBunchCrossNumber(esd.GetBunchCrossNumber());
 374     header->SetOrbitNumber(esd.GetOrbitNumber());
 375     header->SetPeriodNumber(esd.GetPeriodNumber());
 376     header->SetEventType(esd.GetEventType());
 377
 378     header->SetEventNumberESDFile(esd.GetHeader()->GetEventNumberInFile());
 379     if(const_cast<AliESDEvent&>(esd).GetCentrality()){
 380       header->SetCentrality(const_cast<AliESDEvent&>(esd).GetCentrality());
 381     }
 382     else{
 383       header->SetCentrality(0);
 384     }
 385     if(const_cast<AliESDEvent&>(esd).GetEventplane()){
 386       header->SetEventplane(const_cast<AliESDEvent&>(esd).GetEventplane());
 387     }
 388     else{
 389       header->SetEventplane(0);
 390     }
 391   }
 392
 393   // Trigger
 394   header->SetFiredTriggerClasses(esd.GetFiredTriggerClasses());
 395   header->SetTriggerMask(esd.GetTriggerMask());
 396   header->SetTriggerCluster(esd.GetTriggerCluster());
 397   header->SetL0TriggerInputs(esd.GetHeader()->GetL0TriggerInputs());
 398   header->SetL1TriggerInputs(esd.GetHeader()->GetL1TriggerInputs());
 399   header->SetL2TriggerInputs(esd.GetHeader()->GetL2TriggerInputs());
 400
 401   header->SetMagneticField(esd.GetMagneticField());
 402   header->SetMuonMagFieldScale(esd.GetCurrentDip()/6000.);
 403   header->SetZDCN1Energy(esd.GetZDCN1Energy());
 404   header->SetZDCP1Energy(esd.GetZDCP1Energy());
 405   header->SetZDCN2Energy(esd.GetZDCN2Energy());
 406   header->SetZDCP2Energy(esd.GetZDCP2Energy());
 407   header->SetZDCEMEnergy(esd.GetZDCEMEnergy(0),esd.GetZDCEMEnergy(1));
 408
 409   header->SetIRInt2InteractionMap(esd.GetHeader()->GetIRInt2InteractionMap());
 410   header->SetIRInt1InteractionMap(esd.GetHeader()->GetIRInt1InteractionMap());
 411
 412   // ITS Cluster Multiplicty
 413   const AliMultiplicity *mult = esd.GetMultiplicity();
 414   for (Int_t ilay = 0; ilay < 6; ilay++) header->SetITSClusters(ilay, mult->GetNumberOfITSClusters(ilay));
 415
 416   // TPC only Reference Multiplicty
 417   Int_t refMult  = fTPCaloneTrackCuts ? (Short_t)fTPCaloneTrackCuts->GetReferenceMultiplicity(&esd, kTRUE) : -1;
 418   header->SetTPConlyRefMultiplicity(refMult);
 419   //
 420   AliESDtrackCuts::MultEstTrackType estType = esd.GetPrimaryVertexTracks()->GetStatus() ? AliESDtrackCuts::kTrackletsITSTPC : AliESDtrackCuts::kTracklets;
 421   header->SetRefMultiplicityComb05(AliESDtrackCuts::GetReferenceMultiplicity(&esd,estType,0.5));
 422   header->SetRefMultiplicityComb08(AliESDtrackCuts::GetReferenceMultiplicity(&esd,estType,0.8));
 423   //
 424   Float_t diamxy[2]={esd.GetDiamondX(),esd.GetDiamondY()};
 425   Float_t diamcov[3];
 426   esd.GetDiamondCovXY(diamcov);
 427   header->SetDiamond(diamxy,diamcov);
 428   header->SetDiamondZ(esd.GetDiamondZ(),esd.GetSigma2DiamondZ());
 429
 430   // VZERO channel equalization factors for event-plane reconstruction
 431   header->SetVZEROEqFactors(esd.GetVZEROEqFactors());
 432
 433   // T0 Resolution information
 434   const AliESDRun* esdRun = esd.GetESDRun();
 435   for (Int_t i=0;i<AliESDRun::kT0spreadSize;i++) header->SetT0spread(i,esdRun->GetT0spread(i));
 436
 437   return header;
 438 }
 439
 440 //______________________________________________________________________________
 441 void AliAnalysisTaskESDfilter::ConvertCascades(const AliESDEvent& esd)
 442 {
 443
 444   // Convert the cascades part of the ESD.
 445   // Return the number of cascades
 446
 447   AliCodeTimerAuto("",0);
 448
 449   // Create vertices starting from the most complex objects
 450   Double_t chi2 = 0.;
 451
 452   const AliESDVertex* vtx = esd.GetPrimaryVertex();
 453   Double_t pos[3] = { 0. };
 454   Double_t covVtx[6] = { 0. };
 455   Double_t momBach[3]={0.};
 456   Double_t covTr[21]={0.};
 457   Double_t pid[10]={0.};
 458   AliAODPid* detpid(0x0);
 459   AliAODVertex* vV0FromCascade(0x0);
 460   AliAODv0* aodV0(0x0);
 461   AliAODcascade* aodCascade(0x0);
 462   AliAODTrack* aodTrack(0x0);
 463   Double_t momPos[3]={0.};
 464   Double_t momNeg[3] = { 0. };
 465   Double_t momPosAtV0vtx[3]={0.};
 466   Double_t momNegAtV0vtx[3]={0.};
 467   Int_t    tofLabel[3] = {0};
 468   TClonesArray& verticesArray = Vertices();
 469   TClonesArray& tracksArray = Tracks();
 470   TClonesArray& cascadesArray = Cascades();
 471
 472   // Cascades (Modified by A.Maire - February 2009)
 473   for (Int_t nCascade = 0; nCascade < esd.GetNumberOfCascades(); ++nCascade) {
 474
 475     // 0- Preparation
 476     //
 477     AliESDcascade *esdCascade = esd.GetCascade(nCascade);
 478                 Int_t  idxPosFromV0Dghter  = esdCascade->GetPindex();
 479                 Int_t  idxNegFromV0Dghter  = esdCascade->GetNindex();
 480                 Int_t  idxBachFromCascade  = esdCascade->GetBindex();
 481
 482     AliESDtrack  *esdCascadePos  = esd.GetTrack( idxPosFromV0Dghter);
 483     AliESDtrack  *esdCascadeNeg  = esd.GetTrack( idxNegFromV0Dghter);
 484     AliESDtrack  *esdCascadeBach = esd.GetTrack( idxBachFromCascade);
 485
 486     // Identification of the V0 within the esdCascade (via both daughter track indices)
 487     AliESDv0 * currentV0   = 0x0;
 488     Int_t      idxV0FromCascade = -1;
 489
 490     for (Int_t iV0=0; iV0<esd.GetNumberOfV0s(); ++iV0) {
 491
 492       currentV0 = esd.GetV0(iV0);
 493       Int_t posCurrentV0 = currentV0->GetPindex();
 494       Int_t negCurrentV0 = currentV0->GetNindex();
 495
 496       if (posCurrentV0==idxPosFromV0Dghter && negCurrentV0==idxNegFromV0Dghter) {
 497         idxV0FromCascade = iV0;
 498         break;
 499       }
 500     }
 501
 502     if(idxV0FromCascade < 0){
 503       printf("Cascade - no matching for the V0 (index V0 = -1) ! Skip ... \n");
 504       continue;
 505     }// a priori, useless check, but safer ... in case of pb with tracks "out of bounds"
 506
 507     AliESDv0 *esdV0FromCascade   = esd.GetV0(idxV0FromCascade);
 508
 509     // 1 - Cascade selection
 510
 511     //  AliESDVertex *esdPrimVtx = new AliESDVertex(*(esd.GetPrimaryVertex()));
 512     //  TList cascadeObjects;
 513     //  cascadeObjects.AddAt(esdV0FromCascade, 0);
 514     //  cascadeObjects.AddAt(esdCascadePos,    1);
 515     //  cascadeObjects.AddAt(esdCascadeNeg,    2);
 516     //  cascadeObjects.AddAt(esdCascade,       3);
 517     //  cascadeObjects.AddAt(esdCascadeBach,   4);
 518     //  cascadeObjects.AddAt(esdPrimVtx,       5);
 519     //
 520     //  UInt_t selectCascade = 0;
 521     //  if (fCascadeFilter) {
 522     //    // selectCascade = fCascadeFilter->IsSelected(&cascadeObjects);
 523     //          // FIXME AliESDCascadeCuts to be implemented ...
 524     //
 525     //          // Here we may encounter a moot point at the V0 level
 526     //          // between the cascade selections and the V0 ones :
 527     //          // the V0 selected along with the cascade (secondary V0) may
 528     //          // usually be removed from the dedicated V0 selections (prim V0) ...
 529     //          // -> To be discussed !
 530     //
 531     //    // this is a little awkward but otherwise the
 532     //    // list wants to access the pointer (delete it)
 533     //    // again when going out of scope
 534     //    delete cascadeObjects.RemoveAt(5); // esdPrimVtx created via copy construct
 535     //    esdPrimVtx = 0;
 536     //    if (!selectCascade)
 537     //      continue;
 538     //  }
 539     //  else{
 540     //    delete cascadeObjects.RemoveAt(5); // esdPrimVtx created via copy construct
 541     //    esdPrimVtx = 0;
 542     //  }
 543
 544     // 2 - Add the cascade vertex
 545
 546     esdCascade->GetXYZcascade(pos[0], pos[1], pos[2]);
 547     esdCascade->GetPosCovXi(covVtx);
 548     chi2 = esdCascade->GetChi2Xi();
 549
 550     AliAODVertex *vCascade = new(verticesArray[fNumberOfVertices++]) AliAODVertex( pos,
 551                                                                      covVtx,
 552                                                                      chi2, // FIXME = Chi2/NDF will be needed
 553                                                                      fPrimaryVertex,
 554                                                                      nCascade, // id
 555                                                                      AliAODVertex::kCascade);
 556     fPrimaryVertex->AddDaughter(vCascade);
 557
 558 //    if (fDebug > 2) {
 559 //      printf("---- Cascade / Cascade Vertex (AOD) : \n");
 560 //      vCascade->Print();
 561 //    }
 562
 563 //    if(esd.GetTOFHeader() && fIsPidOwner) fESDpid->SetTOFResponse(const_cast<AliESDEvent*>(&esd), (AliESDpid::EStartTimeType_t)fTimeZeroType); //in case of AOD production starting form LHC10e without Tender.
 564
 565
 566     // 3 - Add the bachelor track from the cascade
 567
 568     if (!fUsedTrack[idxBachFromCascade]) {
 569
 570       esdCascadeBach->GetPxPyPz(momBach);
 571       esdCascadeBach->GetXYZ(pos);
 572       esdCascadeBach->GetCovarianceXYZPxPyPz(covTr);
 573       esdCascadeBach->GetESDpid(pid);
 574       esdCascadeBach->GetTOFLabel(tofLabel);
 575
 576             fUsedTrack[idxBachFromCascade] = kTRUE;
 577             UInt_t selectInfo = 0;
 578             if (fTrackFilter) selectInfo = fTrackFilter->IsSelected(esdCascadeBach);
 579             if (fMChandler) fMChandler->SelectParticle(esdCascadeBach->GetLabel());
 580             aodTrack = new(tracksArray[fNumberOfTracks++]) AliAODTrack(esdCascadeBach->GetID(),
 581                                                                            esdCascadeBach->GetLabel(),
 582                                                                            momBach,
 583                                                                            kTRUE,
 584                                                                            pos,
 585                                                                            kFALSE, // Why kFALSE for "isDCA" ? FIXME
 586                                                                            covTr,
 587                                                                            (Short_t)esdCascadeBach->GetSign(),
 588                                                                            esdCascadeBach->GetITSClusterMap(),
 589                                                                            pid,
 590                                                                            vCascade,
 591                                                                            kTRUE,  // usedForVtxFit = kFALSE ? FIXME
 592                                                                            vtx->UsesTrack(esdCascadeBach->GetID()),
 593                                                                            AliAODTrack::kSecondary,
 594                                                                            selectInfo);
 595             aodTrack->SetTPCFitMap(esdCascadeBach->GetTPCFitMap());
 596             aodTrack->SetTPCClusterMap(esdCascadeBach->GetTPCClusterMap());
 597             aodTrack->SetTPCSharedMap (esdCascadeBach->GetTPCSharedMap());
 598             aodTrack->SetChi2perNDF(Chi2perNDF(esdCascadeBach));
 599             aodTrack->SetTPCPointsF(esdCascadeBach->GetTPCNclsF());
 600             aodTrack->SetTPCNCrossedRows(UShort_t(esdCascadeBach->GetTPCCrossedRows()));
 601             aodTrack->SetIntegratedLength(esdCascadeBach->GetIntegratedLength());
 602             aodTrack->SetTOFLabel(tofLabel);
 603             fAODTrackRefs->AddAt(aodTrack,idxBachFromCascade);
 604
 605             if (esdCascadeBach->GetSign() > 0) ++fNumberOfPositiveTracks;
 606             aodTrack->ConvertAliPIDtoAODPID();
 607             aodTrack->SetFlags(esdCascadeBach->GetStatus());
 608             SetAODPID(esdCascadeBach,aodTrack,detpid);
 609     }
 610     else {
 611             aodTrack = static_cast<AliAODTrack*>( fAODTrackRefs->At(idxBachFromCascade) );
 612     }
 613
 614     vCascade->AddDaughter(aodTrack);
 615
 616 //    if (fDebug > 4) {
 617 //      printf("---- Cascade / bach dghter : \n");
 618 //      aodTrack->Print();
 619 //    }
 620
 621
 622     // 4 - Add the V0 from the cascade.
 623     // = V0vtx + both pos and neg daughter tracks + the aodV0 itself
 624     //
 625
 626     if ( !fUsedV0[idxV0FromCascade] ) {
 627       // 4.A - if VO structure hasn't been created yet
 628
 629       // 4.A.1 - Create the V0 vertex of the cascade
 630
 631       esdV0FromCascade->GetXYZ(pos[0], pos[1], pos[2]);
 632       esdV0FromCascade->GetPosCov(covVtx);
 633       chi2 = esdV0FromCascade->GetChi2V0();  // = chi2/NDF since NDF = 2*2-3 ?
 634
 635       vV0FromCascade = new(verticesArray[fNumberOfVertices++]) AliAODVertex(pos,
 636                                                                          covVtx,
 637                                                                          chi2,
 638                                                                          vCascade,
 639                                                                          idxV0FromCascade, //id of ESDv0
 640                                                                          AliAODVertex::kV0);
 641       // Note:
 642       //    one V0 can be used by several cascades.
 643       // So, one AOD V0 vtx can have several parent vtx.
 644       // This is not directly allowed by AliAODvertex.
 645       // Setting the parent vtx (here = param "vCascade") doesn't lead to a crash
 646       // but to a problem of consistency within AODEvent.
 647       // -> See below paragraph 4.B, for the proposed treatment of such a case.
 648
 649       // Add the vV0FromCascade to the aodVOVtxRefs
 650       fAODV0VtxRefs->AddAt(vV0FromCascade,idxV0FromCascade);
 651
 652
 653       // 4.A.2 - Add the positive tracks from the V0
 654
 655       esdCascadePos->GetPxPyPz(momPos);
 656       esdCascadePos->GetXYZ(pos);
 657       esdCascadePos->GetCovarianceXYZPxPyPz(covTr);
 658       esdCascadePos->GetESDpid(pid);
 659       esdCascadePos->GetTOFLabel(tofLabel);
 660
 661       if (!fUsedTrack[idxPosFromV0Dghter]) {
 662         fUsedTrack[idxPosFromV0Dghter] = kTRUE;
 663
 664         UInt_t selectInfo = 0;
 665         if (fTrackFilter) selectInfo = fTrackFilter->IsSelected(esdCascadePos);
 666         if(fMChandler) fMChandler->SelectParticle(esdCascadePos->GetLabel());
 667         aodTrack = new(tracksArray[fNumberOfTracks++])
 668         AliAODTrack(  esdCascadePos->GetID(),
 669                     esdCascadePos->GetLabel(),
 670                     momPos,
 671                     kTRUE,
 672                     pos,
 673                     kFALSE, // Why kFALSE for "isDCA" ? FIXME
 674                     covTr,
 675                     (Short_t)esdCascadePos->GetSign(),
 676                     esdCascadePos->GetITSClusterMap(),
 677                     pid,
 678                     vV0FromCascade,
 679                     kTRUE,  // usedForVtxFit = kFALSE ? FIXME
 680                     vtx->UsesTrack(esdCascadePos->GetID()),
 681                     AliAODTrack::kSecondary,
 682                     selectInfo);
 683         aodTrack->SetTPCFitMap(esdCascadePos->GetTPCFitMap());
 684         aodTrack->SetTPCClusterMap(esdCascadePos->GetTPCClusterMap());
 685         aodTrack->SetTPCSharedMap (esdCascadePos->GetTPCSharedMap());
 686         aodTrack->SetChi2perNDF(Chi2perNDF(esdCascadePos));
 687         aodTrack->SetTPCPointsF(esdCascadePos->GetTPCNclsF());
 688         aodTrack->SetTPCNCrossedRows(UShort_t(esdCascadePos->GetTPCCrossedRows()));
 689         aodTrack->SetIntegratedLength(esdCascadePos->GetIntegratedLength());
 690         aodTrack->SetTOFLabel(tofLabel);
 691         fAODTrackRefs->AddAt(aodTrack,idxPosFromV0Dghter);
 692
 693         if (esdCascadePos->GetSign() > 0) ++fNumberOfPositiveTracks;
 694         aodTrack->ConvertAliPIDtoAODPID();
 695         aodTrack->SetFlags(esdCascadePos->GetStatus());
 696         SetAODPID(esdCascadePos,aodTrack,detpid);
 697       }
 698       else {
 699         aodTrack = static_cast<AliAODTrack*>(fAODTrackRefs->At(idxPosFromV0Dghter));
 700       }
 701       vV0FromCascade->AddDaughter(aodTrack);
 702
 703
 704       // 4.A.3 - Add the negative tracks from the V0
 705
 706       esdCascadeNeg->GetPxPyPz(momNeg);
 707       esdCascadeNeg->GetXYZ(pos);
 708       esdCascadeNeg->GetCovarianceXYZPxPyPz(covTr);
 709       esdCascadeNeg->GetESDpid(pid);
 710       esdCascadeNeg->GetTOFLabel(tofLabel);
 711
 712
 713       if (!fUsedTrack[idxNegFromV0Dghter]) {
 714         fUsedTrack[idxNegFromV0Dghter] = kTRUE;
 715
 716         UInt_t selectInfo = 0;
 717         if (fTrackFilter) selectInfo = fTrackFilter->IsSelected(esdCascadeNeg);
 718         if(fMChandler)fMChandler->SelectParticle(esdCascadeNeg->GetLabel());
 719         aodTrack = new(tracksArray[fNumberOfTracks++]) AliAODTrack(  esdCascadeNeg->GetID(),
 720                                                       esdCascadeNeg->GetLabel(),
 721                                                       momNeg,
 722                                                       kTRUE,
 723                                                       pos,
 724                                                       kFALSE, // Why kFALSE for "isDCA" ? FIXME
 725                                                       covTr,
 726                                                       (Short_t)esdCascadeNeg->GetSign(),
 727                                                       esdCascadeNeg->GetITSClusterMap(),
 728                                                       pid,
 729                                                       vV0FromCascade,
 730                                                       kTRUE,  // usedForVtxFit = kFALSE ? FIXME
 731                                                       vtx->UsesTrack(esdCascadeNeg->GetID()),
 732                                                       AliAODTrack::kSecondary,
 733                                                       selectInfo);
 734         aodTrack->SetTPCFitMap(esdCascadeNeg->GetTPCFitMap());
 735         aodTrack->SetTPCClusterMap(esdCascadeNeg->GetTPCClusterMap());
 736         aodTrack->SetTPCSharedMap (esdCascadeNeg->GetTPCSharedMap());
 737         aodTrack->SetChi2perNDF(Chi2perNDF(esdCascadeNeg));
 738         aodTrack->SetTPCPointsF(esdCascadeNeg->GetTPCNclsF());
 739         aodTrack->SetTPCNCrossedRows(UShort_t(esdCascadeNeg->GetTPCCrossedRows()));
 740         aodTrack->SetIntegratedLength(esdCascadeNeg->GetIntegratedLength());
 741         aodTrack->SetTOFLabel(tofLabel);
 742         fAODTrackRefs->AddAt(aodTrack,idxNegFromV0Dghter);
 743
 744         if (esdCascadeNeg->GetSign() > 0) ++fNumberOfPositiveTracks;
 745         aodTrack->ConvertAliPIDtoAODPID();
 746         aodTrack->SetFlags(esdCascadeNeg->GetStatus());
 747         SetAODPID(esdCascadeNeg,aodTrack,detpid);
 748       }
 749       else {
 750         aodTrack = static_cast<AliAODTrack*>(fAODTrackRefs->At(idxNegFromV0Dghter));
 751       }
 752
 753       vV0FromCascade->AddDaughter(aodTrack);
 754
 755
 756       // 4.A.4 - Add the V0 from cascade to the V0 array
 757
 758       Double_t  dcaV0Daughters      = esdV0FromCascade->GetDcaV0Daughters();
 759       Double_t  dcaV0ToPrimVertex   = esdV0FromCascade->GetD( esd.GetPrimaryVertex()->GetX(),
 760                                                              esd.GetPrimaryVertex()->GetY(),
 761                                                              esd.GetPrimaryVertex()->GetZ() );
 762       esdV0FromCascade->GetPPxPyPz( momPosAtV0vtx[0],momPosAtV0vtx[1],momPosAtV0vtx[2] );
 763       esdV0FromCascade->GetNPxPyPz( momNegAtV0vtx[0],momNegAtV0vtx[1],momNegAtV0vtx[2] );
 764
 765       Double_t dcaDaughterToPrimVertex[2] = { 999., 999.}; // ..[0] = DCA in (x,y) for Pos and ..[1] = Neg
 766       dcaDaughterToPrimVertex[0] = TMath::Abs(esdCascadePos->GetD(      esd.GetPrimaryVertex()->GetX(),
 767                                                                   esd.GetPrimaryVertex()->GetY(),
 768                                                                   esd.GetMagneticField())        );
 769       dcaDaughterToPrimVertex[1] = TMath::Abs(esdCascadeNeg->GetD(      esd.GetPrimaryVertex()->GetX(),
 770                                                                   esd.GetPrimaryVertex()->GetY(),
 771                                                                   esd.GetMagneticField())        );
 772
 773       aodV0 = new(V0s()[fNumberOfV0s++]) AliAODv0( vV0FromCascade,
 774                                                   dcaV0Daughters,
 775                                                   dcaV0ToPrimVertex,
 776                                                   momPosAtV0vtx,
 777                                                   momNegAtV0vtx,
 778                                                   dcaDaughterToPrimVertex);
 779       // set the aod v0 on-the-fly status
 780       aodV0->SetOnFlyStatus(esdV0FromCascade->GetOnFlyStatus());
 781
 782       // Add the aodV0 to the aodVORefs
 783       fAODV0Refs->AddAt(aodV0,idxV0FromCascade);
 784
 785       fUsedV0[idxV0FromCascade] = kTRUE;
 786
 787     } else {
 788       // 4.B - if V0 structure already used
 789
 790       // Note :
 791       //    one V0 can be used by several cascades (frequent in PbPb evts) :
 792       // same V0 which used but attached to different bachelor tracks
 793       // -> aodVORefs and fAODV0VtxRefs are needed.
 794       // Goal : avoid a redundancy of the info in "Vertices" and "v0s" clones array.
 795
 796       vV0FromCascade = static_cast<AliAODVertex*>( fAODV0VtxRefs->At(idxV0FromCascade) );
 797       aodV0          = static_cast<AliAODv0*>    ( fAODV0Refs   ->At(idxV0FromCascade) );
 798
 799       // - Treatment of the parent for such a "re-used" V0 :
 800       // Insert the cascade that reuses the V0 vertex in the lineage chain
 801       // Before : vV0 -> vCascade1 -> vPrimary
 802       //  - Hyp : cascade2 uses the same V0 as cascade1
 803       //  After :  vV0 -> vCascade2 -> vCascade1 -> vPrimary
 804
 805       AliAODVertex *vCascadePreviousParent = static_cast<AliAODVertex*> (vV0FromCascade->GetParent());
 806                         vV0FromCascade->SetParent(vCascade);
 807                         vCascade      ->SetParent(vCascadePreviousParent);
 808
 809 //      if(fDebug > 2)
 810 //        printf("---- Cascade / Lineage insertion\n"
 811 //               "Parent of V0 vtx                 = Cascade vtx %p\n"
 812 //               "Parent of the cascade vtx        = Cascade vtx %p\n"
 813 //               "Parent of the parent cascade vtx = Cascade vtx %p\n",
 814 //               static_cast<void*> (vV0FromCascade->GetParent()),
 815 //               static_cast<void*> (vCascade->GetParent()),
 816 //               static_cast<void*> (vCascadePreviousParent->GetParent()) );
 817
 818     }// end if V0 structure already used
 819
 820 //    if (fDebug > 2) {
 821 //      printf("---- Cascade / V0 vertex: \n");
 822 //      vV0FromCascade->Print();
 823 //    }
 824 //
 825 //    if (fDebug > 4) {
 826 //      printf("---- Cascade / pos dghter : \n");
 827 //                      aodTrack->Print();
 828 //      printf("---- Cascade / neg dghter : \n");
 829 //                      aodTrack->Print();
 830 //      printf("---- Cascade / aodV0 : \n");
 831 //                      aodV0->Print();
 832 //    }
 833
 834     // In any case (used V0 or not), add the V0 vertex to the cascade one.
 835     vCascade->AddDaughter(vV0FromCascade);
 836
 837
 838     // 5 - Add the primary track of the cascade (if any)
 839
 840
 841     // 6 - Add the cascade to the AOD array of cascades
 842
 843     Double_t dcaBachToPrimVertexXY = TMath::Abs(esdCascadeBach->GetD(esd.GetPrimaryVertex()->GetX(),
 844                                                                      esd.GetPrimaryVertex()->GetY(),
 845                                                                      esd.GetMagneticField())        );
 846
 847     Double_t momBachAtCascadeVtx[3]={0.};
 848
 849     esdCascade->GetBPxPyPz(momBachAtCascadeVtx[0], momBachAtCascadeVtx[1], momBachAtCascadeVtx[2]);
 850
 851     aodCascade = new(cascadesArray[fNumberOfCascades++]) AliAODcascade(  vCascade,
 852                                                           esdCascade->Charge(),
 853                                                           esdCascade->GetDcaXiDaughters(),
 854                                                           -999.,
 855                                                           // DCAXiToPrimVtx -> needs to be calculated   ----|
 856                                                           // doesn't exist at ESD level;
 857                                                           // See AODcascade::DcaXiToPrimVertex(Double, Double, Double)
 858                                                           dcaBachToPrimVertexXY,
 859                                                           momBachAtCascadeVtx,
 860                                                           *aodV0);
 861
 862     if (fDebug > 10) {
 863       printf("---- Cascade / AOD cascade : \n\n");
 864       aodCascade->PrintXi(fPrimaryVertex->GetX(), fPrimaryVertex->GetY(), fPrimaryVertex->GetZ());
 865     }
 866
 867   } // end of the loop on cascades
 868
 869   Cascades().Expand(fNumberOfCascades);
 870 }
 871
 872 //______________________________________________________________________________
 873 void AliAnalysisTaskESDfilter::ConvertV0s(const AliESDEvent& esd)
 874 {
 875   // Access to the AOD container of V0s
 876
 877   AliCodeTimerAuto("",0);
 878
 879   //
 880   // V0s
 881   //
 882
 883   Double_t pos[3] = { 0. };
 884   Double_t chi2(0.0);
 885   Double_t covVtx[6] = { 0. };
 886   Double_t momPos[3]={0.};
 887   Double_t covTr[21]={0.};
 888   Double_t pid[10]={0.};
 889   AliAODTrack* aodTrack(0x0);
 890   AliAODPid* detpid(0x0);
 891   Double_t momNeg[3]={0.};
 892   Double_t momPosAtV0vtx[3]={0.};
 893   Double_t momNegAtV0vtx[3]={0.};
 894   Int_t    tofLabel[3] = {0};
 895   for (Int_t nV0 = 0; nV0 < esd.GetNumberOfV0s(); ++nV0)
 896   {
 897     if (fUsedV0[nV0]) continue; // skip if already added to the AOD
 898
 899     AliESDv0 *v0 = esd.GetV0(nV0);
 900     Int_t posFromV0 = v0->GetPindex();
 901     Int_t negFromV0 = v0->GetNindex();
 902
 903     // V0 selection
 904     //
 905     AliESDVertex *esdVtx   = new AliESDVertex(*(esd.GetPrimaryVertex()));
 906     AliESDtrack  *esdV0Pos = esd.GetTrack(posFromV0);
 907     AliESDtrack  *esdV0Neg = esd.GetTrack(negFromV0);
 908     TList v0objects;
 909     v0objects.AddAt(v0,                      0);
 910     v0objects.AddAt(esdV0Pos,                1);
 911     v0objects.AddAt(esdV0Neg,                2);
 912     v0objects.AddAt(esdVtx,                  3);
 913     UInt_t selectV0 = 0;
 914     if (fV0Filter) {
 915       selectV0 = fV0Filter->IsSelected(&v0objects);
 916       // this is a little awkward but otherwise the
 917       // list wants to access the pointer (delete it)
 918       // again when going out of scope
 919       delete v0objects.RemoveAt(3); // esdVtx created via copy construct
 920       esdVtx = 0;
 921       if (!selectV0)
 922         continue;
 923     }
 924     else{
 925       delete v0objects.RemoveAt(3); // esdVtx created via copy construct
 926       esdVtx = 0;
 927     }
 928
 929     v0->GetXYZ(pos[0], pos[1], pos[2]);
 930
 931     if (!fOldESDformat) {
 932             chi2 = v0->GetChi2V0(); // = chi2/NDF since NDF = 2*2-3
 933             v0->GetPosCov(covVtx);
 934     } else {
 935             chi2 = -999.;
 936             for (Int_t i = 0; i < 6; i++)  covVtx[i] = 0.;
 937     }
 938
 939
 940     AliAODVertex * vV0 =
 941           new(Vertices()[fNumberOfVertices++]) AliAODVertex(pos,
 942                                                       covVtx,
 943                                                       chi2,
 944                                                       fPrimaryVertex,
 945                                                       nV0,
 946                                                       AliAODVertex::kV0);
 947     fPrimaryVertex->AddDaughter(vV0);
 948
 949
 950     // Add the positive tracks from the V0
 951
 952
 953     esdV0Pos->GetPxPyPz(momPos);
 954     esdV0Pos->GetXYZ(pos);
 955     esdV0Pos->GetCovarianceXYZPxPyPz(covTr);
 956     esdV0Pos->GetESDpid(pid);
 957     esdV0Pos->GetTOFLabel(tofLabel);
 958
 959     const AliESDVertex *vtx = esd.GetPrimaryVertex();
 960
 961     if (!fUsedTrack[posFromV0]) {
 962             fUsedTrack[posFromV0] = kTRUE;
 963             UInt_t selectInfo = 0;
 964             if (fTrackFilter) selectInfo = fTrackFilter->IsSelected(esdV0Pos);
 965             if(fMChandler)fMChandler->SelectParticle(esdV0Pos->GetLabel());
 966             aodTrack = new(Tracks()[fNumberOfTracks++]) AliAODTrack(esdV0Pos->GetID(),
 967                                                     esdV0Pos->GetLabel(),
 968                                                     momPos,
 969                                                     kTRUE,
 970                                                     pos,
 971                                                     kFALSE,
 972                                                     covTr,
 973                                                     (Short_t)esdV0Pos->GetSign(),
 974                                                     esdV0Pos->GetITSClusterMap(),
 975                                                     pid,
 976                                                     vV0,
 977                                                     kTRUE,  // check if this is right
 978                                                     vtx->UsesTrack(esdV0Pos->GetID()),
 979                                                     AliAODTrack::kSecondary,
 980                                                     selectInfo);
 981             aodTrack->SetTPCFitMap(esdV0Pos->GetTPCFitMap());
 982             aodTrack->SetTPCClusterMap(esdV0Pos->GetTPCClusterMap());
 983             aodTrack->SetTPCSharedMap (esdV0Pos->GetTPCSharedMap());
 984             aodTrack->SetChi2perNDF(Chi2perNDF(esdV0Pos));
 985             aodTrack->SetTPCPointsF(esdV0Pos->GetTPCNclsF());
 986             aodTrack->SetTPCNCrossedRows(UShort_t(esdV0Pos->GetTPCCrossedRows()));
 987             aodTrack->SetIntegratedLength(esdV0Pos->GetIntegratedLength());
 988             aodTrack->SetTOFLabel(tofLabel);
 989             fAODTrackRefs->AddAt(aodTrack,posFromV0);
 990             //      if (fDebug > 0) printf("-------------------Bo: pos track from original pt %.3f \n",aodTrack->Pt());
 991             if (esdV0Pos->GetSign() > 0) ++fNumberOfPositiveTracks;
 992             aodTrack->ConvertAliPIDtoAODPID();
 993             aodTrack->SetFlags(esdV0Pos->GetStatus());
 994             SetAODPID(esdV0Pos,aodTrack,detpid);
 995     }
 996     else {
 997             aodTrack = static_cast<AliAODTrack*>(fAODTrackRefs->At(posFromV0));
 998             //      if (fDebug > 0) printf("-------------------Bo pos track from refArray pt %.3f \n",aodTrack->Pt());
 999     }
1000     vV0->AddDaughter(aodTrack);
1001
1002     // Add the negative tracks from the V0
1003
1004     esdV0Neg->GetPxPyPz(momNeg);
1005     esdV0Neg->GetXYZ(pos);
1006     esdV0Neg->GetCovarianceXYZPxPyPz(covTr);
1007     esdV0Neg->GetESDpid(pid);
1008     esdV0Neg->GetTOFLabel(tofLabel);
1009
1010     if (!fUsedTrack[negFromV0]) {
1011             fUsedTrack[negFromV0] = kTRUE;
1012             UInt_t selectInfo = 0;
1013             if (fTrackFilter) selectInfo = fTrackFilter->IsSelected(esdV0Neg);
1014             if(fMChandler)fMChandler->SelectParticle(esdV0Neg->GetLabel());
1015             aodTrack = new(Tracks()[fNumberOfTracks++]) AliAODTrack(esdV0Neg->GetID(),
1016                                                     esdV0Neg->GetLabel(),
1017                                                     momNeg,
1018                                                     kTRUE,
1019                                                     pos,
1020                                                     kFALSE,
1021                                                     covTr,
1022                                                     (Short_t)esdV0Neg->GetSign(),
1023                                                     esdV0Neg->GetITSClusterMap(),
1024                                                     pid,
1025                                                     vV0,
1026                                                     kTRUE,  // check if this is right
1027                                                     vtx->UsesTrack(esdV0Neg->GetID()),
1028                                                     AliAODTrack::kSecondary,
1029                                                     selectInfo);
1030             aodTrack->SetTPCFitMap(esdV0Neg->GetTPCFitMap());
1031             aodTrack->SetTPCClusterMap(esdV0Neg->GetTPCClusterMap());
1032             aodTrack->SetTPCSharedMap (esdV0Neg->GetTPCSharedMap());
1033             aodTrack->SetChi2perNDF(Chi2perNDF(esdV0Neg));
1034             aodTrack->SetTPCPointsF(esdV0Neg->GetTPCNclsF());
1035             aodTrack->SetTPCNCrossedRows(UShort_t(esdV0Neg->GetTPCCrossedRows()));
1036             aodTrack->SetIntegratedLength(esdV0Neg->GetIntegratedLength());
1037             aodTrack->SetTOFLabel(tofLabel);
1038             fAODTrackRefs->AddAt(aodTrack,negFromV0);
1039             //      if (fDebug > 0) printf("-------------------Bo: neg track from original pt %.3f \n",aodTrack->Pt());
1040             if (esdV0Neg->GetSign() > 0) ++fNumberOfPositiveTracks;
1041             aodTrack->ConvertAliPIDtoAODPID();
1042             aodTrack->SetFlags(esdV0Neg->GetStatus());
1043             SetAODPID(esdV0Neg,aodTrack,detpid);
1044     }
1045     else {
1046             aodTrack = static_cast<AliAODTrack*>(fAODTrackRefs->At(negFromV0));
1047             //      if (fDebug > 0) printf("-------------------Bo neg track from refArray pt %.3f \n",aodTrack->Pt());
1048     }
1049     vV0->AddDaughter(aodTrack);
1050
1051
1052     // Add the V0 the V0 array as well
1053
1054     Double_t  dcaV0Daughters      = v0->GetDcaV0Daughters();
1055     Double_t  dcaV0ToPrimVertex   = v0->GetD(esd.GetPrimaryVertex()->GetX(),
1056                                              esd.GetPrimaryVertex()->GetY(),
1057                                              esd.GetPrimaryVertex()->GetZ());
1058
1059     v0->GetPPxPyPz(momPosAtV0vtx[0],momPosAtV0vtx[1],momPosAtV0vtx[2]);
1060     v0->GetNPxPyPz(momNegAtV0vtx[0],momNegAtV0vtx[1],momNegAtV0vtx[2]);
1061
1062     Double_t dcaDaughterToPrimVertex[2] = { 999., 999.}; // ..[0] = DCA in (x,y) for Pos and ..[1] = Neg
1063     dcaDaughterToPrimVertex[0] = TMath::Abs(esdV0Pos->GetD(  esd.GetPrimaryVertex()->GetX(),
1064                                                            esd.GetPrimaryVertex()->GetY(),
1065                                                            esd.GetMagneticField()) );
1066     dcaDaughterToPrimVertex[1] = TMath::Abs(esdV0Neg->GetD(  esd.GetPrimaryVertex()->GetX(),
1067                                                            esd.GetPrimaryVertex()->GetY(),
1068                                                            esd.GetMagneticField()) );
1069
1070     AliAODv0* aodV0 = new(V0s()[fNumberOfV0s++]) AliAODv0(vV0,
1071                                                 dcaV0Daughters,
1072                                                 dcaV0ToPrimVertex,
1073                                                 momPosAtV0vtx,
1074                                                 momNegAtV0vtx,
1075                                                 dcaDaughterToPrimVertex);
1076
1077     // set the aod v0 on-the-fly status
1078     aodV0->SetOnFlyStatus(v0->GetOnFlyStatus());
1079   }//End of loop on V0s
1080
1081   V0s().Expand(fNumberOfV0s);
1082 }
1083
1084 //______________________________________________________________________________
1085 void AliAnalysisTaskESDfilter::ConvertTPCOnlyTracks(const AliESDEvent& esd)
1086 {
1087   // Convert TPC only tracks
1088   // Here we have wo hybrid appraoch to remove fakes
1089   // ******* ITSTPC ********
1090   // Uses a cut on the ITS properties to select global tracks
1091   // which are than marked as HybdridITSTPC for the remainder
1092   // the TPC only tracks are flagged as HybridITSTPConly.
1093   // Note, in order not to get fakes back in the TPC cuts, one needs
1094   // two "ITS" cuts one tight (1) (to throw out fakes) and one lose (2) (to NOT flag the trakcs in the TPC only)
1095   // using cut number (3)
1096   // so fHybridFilterMask == (1)|(2) fTPCFilterMask = (3), Usercode needs to slect with mask = (1)|(3) and track->IsHybridITSTPC()
1097   // ******* TPC ********
1098   // 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
1099   // 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
1100
1101   AliCodeTimerAuto("",0);
1102
1103   // Loop over the tracks and extract and mask out all aod tracks that pass the selections for AODt racks
1104   for(int it = 0;it < fNumberOfTracks;++it)
1105   {
1106     AliAODTrack *tr = (AliAODTrack*)(Tracks().At(it));
1107     if(!tr)continue;
1108     UInt_t map = tr->GetFilterMap();
1109     if(map&fTPCConstrainedFilterMask){
1110       // we only reset the track select ionfo, no deletion...
1111       tr->SetFilterMap(map&~fTPCConstrainedFilterMask);
1112     }
1113     if(map&fHybridFilterMaskTPCCG){
1114       // this is one part of the hybrid tracks
1115       // the others not passing the selection will be TPC only selected below
1116       tr->SetIsHybridTPCConstrainedGlobal(kTRUE);
1117     }
1118   }
1119   // Loop over the ESD trcks and pick out the tracks passing TPC only cuts
1120
1121
1122   const AliESDVertex *vtxSPD = esd.GetPrimaryVertexSPD();
1123   const AliESDVertex *vtx = esd.GetPrimaryVertex();
1124
1125   Double_t pos[3] = { 0. };
1126   Double_t covTr[21]={0.};
1127   Double_t pid[10]={0.};
1128
1129
1130   Double_t p[3] = { 0. };
1131
1132   Double_t pDCA[3] = { 0. }; // momentum at DCA
1133   Double_t rDCA[3] = { 0. }; // position at DCA
1134   Float_t  dDCA[2] = {0.};    // DCA to the vertex d and z
1135   Float_t  cDCA[3] = {0.};    // covariance of impact parameters
1136   Int_t    tofLabel[3] = {0};
1137
1138   AliAODTrack* aodTrack(0x0);
1139   //  AliAODPid* detpid(0x0);
1140
1141   // account for change in pT after the constraint
1142   Float_t ptMax = 1E10;
1143   Float_t ptMin = 0;
1144   for(int i = 0;i<32;i++){
1145     if(fTPCConstrainedFilterMask&(1<<i)){
1146       AliESDtrackCuts*cuts = (AliESDtrackCuts*)fTrackFilter->GetCuts()->At(i);
1147       Float_t tmp1= 0,tmp2 = 0;
1148       cuts->GetPtRange(tmp1,tmp2);
1149       if(tmp1>ptMin)ptMin=tmp1;
1150       if(tmp2<ptMax)ptMax=tmp2;
1151     }
1152   }
1153
1154   for (Int_t nTrack = 0; nTrack < esd.GetNumberOfTracks(); ++nTrack)
1155   {
1156     AliESDtrack* esdTrack = esd.GetTrack(nTrack); //carefull do not modify it othwise  need to work with a copy
1157
1158     UInt_t selectInfo = 0;
1159     Bool_t isHybridITSTPC = false;
1160     //
1161     // Track selection
1162     if (fTrackFilter) {
1163       selectInfo = fTrackFilter->IsSelected(esdTrack);
1164     }
1165
1166     if(!(selectInfo&fHybridFilterMaskTPCCG)){
1167       // not already selected tracks, use second part of hybrid tracks
1168       isHybridITSTPC = true;
1169       // too save space one could only store these...
1170     }
1171
1172     selectInfo &= fTPCConstrainedFilterMask;
1173     if (!selectInfo)continue;
1174     if (fWriteHybridTPCCOnly&&!isHybridITSTPC)continue; // write only complementary tracks
1175     // create a tpc only tracl
1176     AliESDtrack *track = AliESDtrackCuts::GetTPCOnlyTrack(const_cast<AliESDEvent*>(&esd),esdTrack->GetID());
1177     if(!track) continue;
1178
1179     if(track->Pt()>0.)
1180     {
1181       // only constrain tracks above threshold
1182       AliExternalTrackParam exParam;
1183       // take the B-field from the ESD, no 3D fieldMap available at this point
1184       Bool_t relate = false;
1185       relate = track->RelateToVertexTPC(vtxSPD,esd.GetMagneticField(),kVeryBig,&exParam);
1186       if(!relate){
1187         delete track;
1188         continue;
1189       }
1190       // fetch the track parameters at the DCA (unconstraint)
1191       if(track->GetTPCInnerParam()){
1192         track->GetTPCInnerParam()->GetPxPyPz(pDCA);
1193         track->GetTPCInnerParam()->GetXYZ(rDCA);
1194       }
1195       // get the DCA to the vertex:
1196       track->GetImpactParametersTPC(dDCA,cDCA);
1197       // set the constrained parameters to the track
1198       track->Set(exParam.GetX(),exParam.GetAlpha(),exParam.GetParameter(),exParam.GetCovariance());
1199     }
1200
1201     track->GetPxPyPz(p);
1202
1203     Float_t pT = track->Pt();
1204     if(pT<ptMin||pT>ptMax){
1205       delete track;
1206       continue;
1207     }
1208
1209     //
1210
1211
1212     track->GetXYZ(pos);
1213     track->GetCovarianceXYZPxPyPz(covTr);
1214     esdTrack->GetESDpid(pid);// original PID
1215     esdTrack->GetTOFLabel(tofLabel);
1216     if(fMChandler)fMChandler->SelectParticle(esdTrack->GetLabel());
1217     aodTrack = new(Tracks()[fNumberOfTracks++]) AliAODTrack((track->GetID()+1)*-1,
1218                                                             track->GetLabel(),
1219                                                             p,
1220                                                             kTRUE,
1221                                                             pos,
1222                                                             kFALSE,
1223                                                             covTr,
1224                                                             (Short_t)track->GetSign(),
1225                                                             track->GetITSClusterMap(),
1226                                                             pid,
1227                                                             fPrimaryVertex,
1228                                                             kTRUE, // check if this is right
1229                                                             vtx->UsesTrack(track->GetID()),
1230                                                             AliAODTrack::kPrimary,
1231                                                             selectInfo);
1232     aodTrack->SetIsHybridTPCConstrainedGlobal(isHybridITSTPC);
1233     aodTrack->SetTPCFitMap(track->GetTPCFitMap());
1234     aodTrack->SetTPCClusterMap(track->GetTPCClusterMap());
1235     aodTrack->SetTPCSharedMap (track->GetTPCSharedMap());
1236     aodTrack->SetIsTPCConstrained(kTRUE);
1237     aodTrack->SetChi2perNDF(Chi2perNDF(esdTrack)); // original track
1238     // set the DCA values to the AOD track
1239     aodTrack->SetPxPyPzAtDCA(pDCA[0],pDCA[1],pDCA[2]);
1240     aodTrack->SetXYAtDCA(rDCA[0],rDCA[1]);
1241     aodTrack->SetDCA(dDCA[0],dDCA[1]);
1242
1243     aodTrack->SetFlags(track->GetStatus());
1244     aodTrack->SetTPCPointsF(track->GetTPCNclsF());
1245     aodTrack->SetTPCNCrossedRows(UShort_t(track->GetTPCCrossedRows()));
1246     aodTrack->SetIntegratedLength(track->GetIntegratedLength());
1247     aodTrack->SetTOFLabel(tofLabel);
1248     //Perform progagation of tracks if needed
1249     if(fDoPropagateTrackToEMCal) PropagateTrackToEMCal(esdTrack);
1250     aodTrack->SetTrackPhiEtaPtOnEMCal(esdTrack->GetTrackPhiOnEMCal(),esdTrack->GetTrackEtaOnEMCal(),esdTrack->GetTrackPtOnEMCal());
1251
1252     // do not duplicate PID information
1253     //    aodTrack->ConvertAliPIDtoAODPID();
1254     //    SetAODPID(esdTrack,aodTrack,detpid);
1255
1256     delete track;
1257   } // end of loop on tracks
1258
1259 }
1260
1261 //______________________________________________________________________________
1262 void AliAnalysisTaskESDfilter::ConvertGlobalConstrainedTracks(const AliESDEvent& esd)
1263 {
1264
1265   // Here we have the option to store the complement from global constraint information
1266   // to tracks passing tight cuts (1) in order not to get fakes back in, one needs
1267   // two sets of cuts one tight (1) (to throw out fakes) and one lose (2) (fakes/bad tracks would pass (2) but not (1))
1268   // using cut number (3) selects the tracks that complement (1) e.g. tracks witout ITS refit or cluster requirement
1269
1270
1271   AliCodeTimerAuto("",0);
1272
1273   // Loop over the tracks and extract and mask out all aod tracks that pass the selections for AODt racks
1274   for(int it = 0;it < fNumberOfTracks;++it)
1275   {
1276     AliAODTrack *tr = (AliAODTrack*)(Tracks().At(it));
1277     if(!tr)continue;
1278     UInt_t map = tr->GetFilterMap();
1279     if(map&fGlobalConstrainedFilterMask){
1280       // we only reset the track select info, no deletion...
1281       // mask reset mask in case track is already taken
1282       tr->SetFilterMap(map&~fGlobalConstrainedFilterMask);
1283     }
1284     if(map&fHybridFilterMaskGCG){
1285       // this is one part of the hybrid tracks
1286       // the others not passing the selection will be the ones selected below
1287       tr->SetIsHybridGlobalConstrainedGlobal(kTRUE);
1288     }
1289   }
1290   // Loop over the ESD trcks and pick out the tracks passing the GlobalConstraint cuts
1291
1292
1293   Double_t pos[3] = { 0. };
1294   Double_t covTr[21]={0.};
1295   Double_t pid[10]={0.};
1296   Double_t p[3] = { 0. };
1297
1298   Double_t pDCA[3] = { 0. }; // momentum at DCA
1299   Double_t rDCA[3] = { 0. }; // position at DCA
1300   Float_t  dDCA[2] = {0.};    // DCA to the vertex d and z
1301   Float_t  cDCA[3] = {0.};    // covariance of impact parameters
1302   Int_t    tofLabel[3] = {0};
1303
1304
1305   AliAODTrack* aodTrack(0x0);
1306   AliAODPid* detpid(0x0);
1307   const AliESDVertex *vtx = esd.GetPrimaryVertex();
1308
1309   // account for change in pT after the constraint
1310   Float_t ptMax = 1E10;
1311   Float_t ptMin = 0;
1312   for(int i = 0;i<32;i++){
1313     if(fGlobalConstrainedFilterMask&(1<<i)){
1314       AliESDtrackCuts*cuts = (AliESDtrackCuts*)fTrackFilter->GetCuts()->At(i);
1315       Float_t tmp1= 0,tmp2 = 0;
1316       cuts->GetPtRange(tmp1,tmp2);
1317       if(tmp1>ptMin)ptMin=tmp1;
1318       if(tmp2<ptMax)ptMax=tmp2;
1319     }
1320   }
1321
1322
1323
1324  for (Int_t nTrack = 0; nTrack < esd.GetNumberOfTracks(); ++nTrack)
1325   {
1326     AliESDtrack* esdTrack = esd.GetTrack(nTrack); //carefull do not modify it othwise  need to work with a copy
1327     const AliExternalTrackParam * exParamGC = esdTrack->GetConstrainedParam();
1328     if(!exParamGC)continue;
1329
1330     UInt_t selectInfo = 0;
1331     Bool_t isHybridGC = false;
1332
1333     //
1334     // Track selection
1335     if (fTrackFilter) {
1336       selectInfo = fTrackFilter->IsSelected(esdTrack);
1337     }
1338
1339
1340     if(!(selectInfo&fHybridFilterMaskGCG))isHybridGC = true;
1341     if (fWriteHybridGCOnly&&!isHybridGC)continue; // write only complementary tracks
1342
1343     selectInfo &= fGlobalConstrainedFilterMask;
1344     if (!selectInfo)continue;
1345     // fetch the track parameters at the DCA (unconstrained)
1346     esdTrack->GetPxPyPz(pDCA);
1347     esdTrack->GetXYZ(rDCA);
1348     // get the DCA to the vertex:
1349     esdTrack->GetImpactParameters(dDCA,cDCA);
1350
1351     if (!esdTrack->GetConstrainedPxPyPz(p)) continue;
1352
1353
1354     Float_t pT = exParamGC->Pt();
1355     if(pT<ptMin||pT>ptMax){
1356       continue;
1357     }
1358
1359
1360     esdTrack->GetConstrainedXYZ(pos);
1361     exParamGC->GetCovarianceXYZPxPyPz(covTr);
1362     esdTrack->GetESDpid(pid);
1363     esdTrack->GetTOFLabel(tofLabel);
1364     if(fMChandler)fMChandler->SelectParticle(esdTrack->GetLabel());
1365     aodTrack = new(Tracks()[fNumberOfTracks++]) AliAODTrack((esdTrack->GetID()+1)*-1,
1366                                                             esdTrack->GetLabel(),
1367                                                             p,
1368                                                             kTRUE,
1369                                                             pos,
1370                                                             kFALSE,
1371                                                             covTr,
1372                                                             (Short_t)esdTrack->GetSign(),
1373                                                             esdTrack->GetITSClusterMap(),
1374                                                             pid,
1375                                                             fPrimaryVertex,
1376                                                             kTRUE, // check if this is right
1377                                                             vtx->UsesTrack(esdTrack->GetID()),
1378                                                             AliAODTrack::kPrimary,
1379                                                             selectInfo);
1380     aodTrack->SetIsHybridGlobalConstrainedGlobal(isHybridGC);
1381     aodTrack->SetIsGlobalConstrained(kTRUE);
1382     aodTrack->SetTPCFitMap(esdTrack->GetTPCFitMap());
1383     aodTrack->SetTPCClusterMap(esdTrack->GetTPCClusterMap());
1384     aodTrack->SetTPCSharedMap (esdTrack->GetTPCSharedMap());
1385     aodTrack->SetChi2perNDF(Chi2perNDF(esdTrack));
1386
1387
1388     // set the DCA values to the AOD track
1389     aodTrack->SetPxPyPzAtDCA(pDCA[0],pDCA[1],pDCA[2]);
1390     aodTrack->SetXYAtDCA(rDCA[0],rDCA[1]);
1391     aodTrack->SetDCA(dDCA[0],dDCA[1]);
1392
1393     aodTrack->SetFlags(esdTrack->GetStatus());
1394     aodTrack->SetTPCPointsF(esdTrack->GetTPCNclsF());
1395     aodTrack->SetTPCNCrossedRows(UShort_t(esdTrack->GetTPCCrossedRows()));
1396     aodTrack->SetIntegratedLength(esdTrack->GetIntegratedLength());
1397     aodTrack->SetTOFLabel(tofLabel);
1398     if(isHybridGC){
1399       // only copy AOD information for hybrid, no duplicate information
1400       aodTrack->ConvertAliPIDtoAODPID();
1401       SetAODPID(esdTrack,aodTrack,detpid);
1402     }
1403
1404     //Perform progagation of tracks if needed
1405     if(fDoPropagateTrackToEMCal) PropagateTrackToEMCal(esdTrack);
1406     aodTrack->SetTrackPhiEtaPtOnEMCal(esdTrack->GetTrackPhiOnEMCal(),esdTrack->GetTrackEtaOnEMCal(),esdTrack->GetTrackPtOnEMCal());
1407   } // end of loop on tracks
1408
1409 }
1410
1411
1412 //______________________________________________________________________________
1413 void AliAnalysisTaskESDfilter::ConvertTracks(const AliESDEvent& esd)
1414 {
1415   // Tracks (primary and orphan)
1416
1417   AliCodeTimerAuto("",0);
1418
1419   AliDebug(1,Form("NUMBER OF ESD TRACKS %5d\n", esd.GetNumberOfTracks()));
1420
1421   const AliESDVertex *vtx = esd.GetPrimaryVertex();
1422   Double_t p[3] = { 0. };
1423   Double_t pos[3] = { 0. };
1424   Double_t covTr[21] = { 0. };
1425   Double_t pid[10] = { 0. };
1426   Int_t    tofLabel[3] = {0};
1427   AliAODTrack* aodTrack(0x0);
1428   AliAODPid* detpid(0x0);
1429
1430   for (Int_t nTrack = 0; nTrack < esd.GetNumberOfTracks(); ++nTrack)
1431   {
1432     if (fUsedTrack[nTrack]) continue;
1433
1434     AliESDtrack *esdTrack = esd.GetTrack(nTrack);
1435     UInt_t selectInfo = 0;
1436     //
1437     // Track selection
1438     if (fTrackFilter) {
1439             selectInfo = fTrackFilter->IsSelected(esdTrack);
1440             if (!selectInfo && !vtx->UsesTrack(esdTrack->GetID())) continue;
1441     }
1442
1443
1444     esdTrack->GetPxPyPz(p);
1445     esdTrack->GetXYZ(pos);
1446     esdTrack->GetCovarianceXYZPxPyPz(covTr);
1447     esdTrack->GetESDpid(pid);
1448     esdTrack->GetTOFLabel(tofLabel);
1449     if(fMChandler)fMChandler->SelectParticle(esdTrack->GetLabel());
1450     fUsedTrack[nTrack] = kTRUE;
1451     fPrimaryVertex->AddDaughter(aodTrack =
1452                          new(Tracks()[fNumberOfTracks++]) AliAODTrack(esdTrack->GetID(),
1453                                                             esdTrack->GetLabel(),
1454                                                             p,
1455                                                             kTRUE,
1456                                                             pos,
1457                                                             kFALSE,
1458                                                             covTr,
1459                                                             (Short_t)esdTrack->GetSign(),
1460                                                             esdTrack->GetITSClusterMap(),
1461                                                             pid,
1462                                                             fPrimaryVertex,
1463                                                             kTRUE, // check if this is right
1464                                                             vtx->UsesTrack(esdTrack->GetID()),
1465                                                             AliAODTrack::kPrimary,
1466                                                             selectInfo)
1467                          );
1468     aodTrack->SetTPCFitMap(esdTrack->GetTPCFitMap());
1469     aodTrack->SetTPCClusterMap(esdTrack->GetTPCClusterMap());
1470     aodTrack->SetTPCSharedMap (esdTrack->GetTPCSharedMap());
1471     aodTrack->SetChi2perNDF(Chi2perNDF(esdTrack));
1472     aodTrack->SetTPCPointsF(esdTrack->GetTPCNclsF());
1473     aodTrack->SetTPCNCrossedRows(UShort_t(esdTrack->GetTPCCrossedRows()));
1474     aodTrack->SetIntegratedLength(esdTrack->GetIntegratedLength());
1475     aodTrack->SetTOFLabel(tofLabel);
1476     if(esdTrack->IsEMCAL()) aodTrack->SetEMCALcluster(esdTrack->GetEMCALcluster());
1477     if(esdTrack->IsPHOS())  aodTrack->SetPHOScluster(esdTrack->GetPHOScluster());
1478
1479     //Perform progagation of tracks if needed
1480     if(fDoPropagateTrackToEMCal) PropagateTrackToEMCal(esdTrack);
1481     aodTrack->SetTrackPhiEtaPtOnEMCal(esdTrack->GetTrackPhiOnEMCal(),esdTrack->GetTrackEtaOnEMCal(),esdTrack->GetTrackPtOnEMCal());
1482
1483     fAODTrackRefs->AddAt(aodTrack, nTrack);
1484
1485
1486     if (esdTrack->GetSign() > 0) ++fNumberOfPositiveTracks;
1487     aodTrack->SetFlags(esdTrack->GetStatus());
1488     aodTrack->ConvertAliPIDtoAODPID();
1489     SetAODPID(esdTrack,aodTrack,detpid);
1490   } // end of loop on tracks
1491 }
1492
1493 //______________________________________________________________________________
1494 void AliAnalysisTaskESDfilter::PropagateTrackToEMCal(AliESDtrack *esdTrack)
1495 {
1496   Double_t trkPos[3] = {0.,0.,0.};
1497   Double_t EMCalEta=-999, EMCalPhi=-999, EMCalPt=-999;
1498   Double_t trkphi = esdTrack->Phi()*TMath::RadToDeg();
1499   if(TMath::Abs(esdTrack->Eta())<0.9 && trkphi > 10 && trkphi < 250 )
1500     {
1501       AliExternalTrackParam *trkParam = const_cast<AliExternalTrackParam*>(esdTrack->GetInnerParam());
1502       if(trkParam)
1503         {
1504           AliExternalTrackParam trkParamTmp(*trkParam);
1505           if(AliTrackerBase::PropagateTrackToBxByBz(&trkParamTmp, fEMCalSurfaceDistance, esdTrack->GetMass(), 20, kTRUE, 0.8, -1))
1506             {
1507               trkParamTmp.GetXYZ(trkPos);
1508               TVector3 trkPosVec(trkPos[0],trkPos[1],trkPos[2]);
1509               EMCalEta = trkPosVec.Eta();
1510               EMCalPhi = trkPosVec.Phi();
1511               EMCalPt = trkParamTmp.Pt();
1512               if(EMCalPhi<0)  EMCalPhi += 2*TMath::Pi();
1513               esdTrack->SetTrackPhiEtaPtOnEMCal(EMCalPhi,EMCalEta,EMCalPt);
1514             }
1515         }
1516     }
1517 }
1518
1519 //______________________________________________________________________________
1520 void AliAnalysisTaskESDfilter::ConvertPmdClusters(const AliESDEvent& esd)
1521 {
1522 // Convert PMD Clusters
1523   AliCodeTimerAuto("",0);
1524   Int_t jPmdClusters=0;
1525   // Access to the AOD container of PMD clusters
1526   TClonesArray &pmdClusters = *(AODEvent()->GetPmdClusters());
1527   for (Int_t iPmd = 0; iPmd < esd.GetNumberOfPmdTracks(); ++iPmd) {
1528     // file pmd clusters, to be revised!
1529     AliESDPmdTrack *pmdTrack = esd.GetPmdTrack(iPmd);
1530     Int_t nLabel = 0;
1531     Int_t *label = 0x0;
1532     Double_t posPmd[3] = { pmdTrack->GetClusterX(), pmdTrack->GetClusterY(), pmdTrack->GetClusterZ()};
1533     Double_t pidPmd[13] = { 0.}; // to be revised!
1534     // type not set!
1535     // assoc cluster not set
1536     new(pmdClusters[jPmdClusters++]) AliAODPmdCluster(iPmd, nLabel, label, pmdTrack->GetClusterADC(), posPmd, pidPmd);
1537   }
1538 }
1539
1540
1541 //______________________________________________________________________________
1542 void AliAnalysisTaskESDfilter::ConvertCaloClusters(const AliESDEvent& esd)
1543 {
1544 // Convert Calorimeter Clusters
1545   AliCodeTimerAuto("",0);
1546
1547   // Access to the AOD container of clusters
1548   TClonesArray &caloClusters = *(AODEvent()->GetCaloClusters());
1549   Int_t jClusters(0);
1550
1551   for (Int_t iClust=0; iClust<esd.GetNumberOfCaloClusters(); ++iClust) {
1552
1553     AliESDCaloCluster * cluster = esd.GetCaloCluster(iClust);
1554
1555     Int_t  id        = cluster->GetID();
1556     Int_t  nLabel    = cluster->GetNLabels();
1557     Int_t *labels    = cluster->GetLabels();
1558     if(labels){
1559                   for(int i = 0;i < nLabel;++i){
1560                           if(fMChandler)fMChandler->SelectParticle(labels[i]);
1561                   }
1562           }
1563
1564     Float_t energy = cluster->E();
1565     Float_t posF[3] = { 0.};
1566     cluster->GetPosition(posF);
1567
1568     AliAODCaloCluster *caloCluster = new(caloClusters[jClusters++]) AliAODCaloCluster(id,
1569                                                                                       nLabel,
1570                                                                                       labels,
1571                                                                                       energy,
1572                                                                                       posF,
1573                                                                                       NULL,
1574                                                                                       cluster->GetType(),0);
1575
1576     caloCluster->SetCaloCluster(cluster->GetDistanceToBadChannel(),
1577                                 cluster->GetDispersion(),
1578                                 cluster->GetM20(), cluster->GetM02(),
1579                                 cluster->GetEmcCpvDistance(),
1580                                 cluster->GetNExMax(),cluster->GetTOF()) ;
1581
1582     caloCluster->SetPIDFromESD(cluster->GetPID());
1583     caloCluster->SetNCells(cluster->GetNCells());
1584     caloCluster->SetCellsAbsId(cluster->GetCellsAbsId());
1585     caloCluster->SetCellsAmplitudeFraction(cluster->GetCellsAmplitudeFraction());
1586
1587     caloCluster->SetTrackDistance(cluster->GetTrackDx(), cluster->GetTrackDz());
1588
1589     Int_t nMatchCount = 0;
1590     TArrayI* matchedT =         cluster->GetTracksMatched();
1591     if (fNumberOfTracks>0 && matchedT && cluster->GetTrackMatchedIndex() >= 0) {
1592       for (Int_t im = 0; im < matchedT->GetSize(); im++) {
1593         Int_t iESDtrack = matchedT->At(im);;
1594         if (fAODTrackRefs->At(iESDtrack) != 0) {
1595           caloCluster->AddTrackMatched((AliAODTrack*)fAODTrackRefs->At(iESDtrack));
1596           nMatchCount++;
1597         }
1598       }
1599     }
1600     if(nMatchCount==0)
1601       caloCluster->SetTrackDistance(-999,-999);
1602
1603   }
1604   caloClusters.Expand(jClusters); // resize TObjArray to 'remove' slots for pseudo clusters
1605 }
1606
1607 //______________________________________________________________________________
1608 void AliAnalysisTaskESDfilter::ConvertCaloTrigger(TString calo, const AliESDEvent& esd)
1609 {
1610         AliCodeTimerAuto("",0);
1611
1612         if (calo == "PHOS")
1613         {
1614           AliAODCaloTrigger &aodTrigger = *(AODEvent()->GetCaloTrigger(calo));
1615           AliESDCaloTrigger &esdTrigger = *(esd.GetCaloTrigger(calo));
1616
1617           aodTrigger.Allocate(esdTrigger.GetEntries());
1618           esdTrigger.Reset();
1619
1620           Float_t a;
1621           Int_t tmod,tabsId;
1622
1623           while (esdTrigger.Next()) {
1624             esdTrigger.GetPosition(tmod,tabsId);
1625             esdTrigger.GetAmplitude(a);
1626             aodTrigger.Add(tmod,tabsId,a,0.,(Int_t*)NULL,0,0,0);
1627           }
1628
1629           return;
1630         }
1631
1632         AliAODHandler *aodHandler = dynamic_cast<AliAODHandler*>(AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler());
1633
1634         if (aodHandler)
1635         {
1636                 TTree *aodTree = aodHandler->GetTree();
1637
1638                 if (aodTree)
1639                 {
1640                         Int_t *type = esd.GetCaloTriggerType();
1641
1642                         for (Int_t i = 0; i < 15; i++)
1643                         {
1644                                 aodTree->GetUserInfo()->Add(new TParameter<int>(Form("EMCALCaloTrigger%d",i), type[i]));
1645                         }
1646                 }
1647         }
1648
1649         AliAODCaloTrigger &aodTrigger = *(AODEvent()->GetCaloTrigger(calo));
1650
1651         AliESDCaloTrigger &esdTrigger = *(esd.GetCaloTrigger(calo));
1652
1653         aodTrigger.Allocate(esdTrigger.GetEntries());
1654
1655         esdTrigger.Reset();
1656         while (esdTrigger.Next())
1657         {
1658                 Int_t px, py, ts, nTimes, times[10], b;
1659                 Float_t a, t;
1660
1661                 esdTrigger.GetPosition(px, py);
1662
1663                 esdTrigger.GetAmplitude(a);
1664                 esdTrigger.GetTime(t);
1665
1666                 esdTrigger.GetL0Times(times);
1667                 esdTrigger.GetNL0Times(nTimes);
1668
1669                 esdTrigger.GetL1TimeSum(ts);
1670
1671                 esdTrigger.GetTriggerBits(b);
1672
1673                 aodTrigger.Add(px, py, a, t, times, nTimes, ts, b);
1674         }
1675
1676         for (int i = 0; i < 4; i++) aodTrigger.SetL1Threshold(i, esdTrigger.GetL1Threshold(i));
1677
1678         Int_t v0[2] =
1679         {
1680                 esdTrigger.GetL1V0(0),
1681                 esdTrigger.GetL1V0(1)
1682         };
1683
1684         aodTrigger.SetL1V0(v0);
1685         aodTrigger.SetL1FrameMask(esdTrigger.GetL1FrameMask());
1686 }
1687
1688 //______________________________________________________________________________
1689 void AliAnalysisTaskESDfilter::ConvertEMCALCells(const AliESDEvent& esd)
1690 {
1691 // Convert EMCAL Cells
1692   AliCodeTimerAuto("",0);
1693   // fill EMCAL cell info
1694   if (esd.GetEMCALCells()) { // protection against missing ESD information
1695     AliESDCaloCells &esdEMcells = *(esd.GetEMCALCells());
1696     Int_t nEMcell = esdEMcells.GetNumberOfCells() ;
1697
1698     AliAODCaloCells &aodEMcells = *(AODEvent()->GetEMCALCells());
1699     aodEMcells.CreateContainer(nEMcell);
1700     aodEMcells.SetType(AliAODCaloCells::kEMCALCell);
1701     for (Int_t iCell = 0; iCell < nEMcell; iCell++) {
1702       aodEMcells.SetCell(iCell,esdEMcells.GetCellNumber(iCell),esdEMcells.GetAmplitude(iCell),
1703                          esdEMcells.GetTime(iCell), esdEMcells.GetMCLabel(iCell), esdEMcells.GetEFraction(iCell));
1704     }
1705     aodEMcells.Sort();
1706   }
1707 }
1708
1709 //______________________________________________________________________________
1710 void AliAnalysisTaskESDfilter::ConvertPHOSCells(const AliESDEvent& esd)
1711 {
1712 // Convert PHOS Cells
1713   AliCodeTimerAuto("",0);
1714   // fill PHOS cell info
1715   if (esd.GetPHOSCells()) { // protection against missing ESD information
1716     AliESDCaloCells &esdPHcells = *(esd.GetPHOSCells());
1717     Int_t nPHcell = esdPHcells.GetNumberOfCells() ;
1718
1719     AliAODCaloCells &aodPHcells = *(AODEvent()->GetPHOSCells());
1720     aodPHcells.CreateContainer(nPHcell);
1721     aodPHcells.SetType(AliAODCaloCells::kPHOSCell);
1722     for (Int_t iCell = 0; iCell < nPHcell; iCell++) {
1723       aodPHcells.SetCell(iCell,esdPHcells.GetCellNumber(iCell),esdPHcells.GetAmplitude(iCell),
1724                          esdPHcells.GetTime(iCell), esdPHcells.GetMCLabel(iCell), esdPHcells.GetEFraction(iCell));
1725     }
1726     aodPHcells.Sort();
1727   }
1728 }
1729
1730 //______________________________________________________________________________
1731 void AliAnalysisTaskESDfilter::ConvertTracklets(const AliESDEvent& esd)
1732 {
1733   // tracklets
1734   AliCodeTimerAuto("",0);
1735
1736   AliAODTracklets &SPDTracklets = *(AODEvent()->GetTracklets());
1737   const AliMultiplicity *mult = esd.GetMultiplicity();
1738   if (mult) {
1739     if (mult->GetNumberOfTracklets()>0) {
1740       SPDTracklets.CreateContainer(mult->GetNumberOfTracklets());
1741
1742       for (Int_t n=0; n<mult->GetNumberOfTracklets(); n++) {
1743         if(fMChandler){
1744           fMChandler->SelectParticle(mult->GetLabel(n, 0));
1745           fMChandler->SelectParticle(mult->GetLabel(n, 1));
1746         }
1747         SPDTracklets.SetTracklet(n, mult->GetTheta(n), mult->GetPhi(n), mult->GetDeltaPhi(n), mult->GetLabel(n, 0),mult->GetLabel(n, 1));
1748       }
1749     }
1750   } else {
1751     //Printf("ERROR: AliMultiplicity could not be retrieved from ESD");
1752   }
1753 }
1754
1755 //______________________________________________________________________________
1756 void AliAnalysisTaskESDfilter::ConvertKinks(const AliESDEvent& esd)
1757 {
1758   AliCodeTimerAuto("",0);
1759
1760   // Kinks: it is a big mess the access to the information in the kinks
1761   // The loop is on the tracks in order to find the mother and daugther of each kink
1762
1763   Double_t covTr[21]={0.};
1764   Double_t pid[10]={0.};
1765   AliAODPid* detpid(0x0);
1766   Int_t tofLabel[3] = {0};
1767
1768   fNumberOfKinks = esd.GetNumberOfKinks();
1769
1770   const AliESDVertex* vtx = esd.GetPrimaryVertex();
1771
1772   for (Int_t iTrack=0; iTrack<esd.GetNumberOfTracks(); ++iTrack)
1773   {
1774     AliESDtrack * esdTrack = esd.GetTrack(iTrack);
1775
1776     Int_t ikink = esdTrack->GetKinkIndex(0);
1777
1778     if (ikink && fNumberOfKinks) {
1779             // Negative kink index: mother, positive: daughter
1780
1781             // Search for the second track of the kink
1782
1783             for (Int_t jTrack = iTrack+1; jTrack<esd.GetNumberOfTracks(); ++jTrack) {
1784
1785         AliESDtrack * esdTrack1 = esd.GetTrack(jTrack);
1786
1787         Int_t jkink = esdTrack1->GetKinkIndex(0);
1788
1789         if ( TMath::Abs(ikink)==TMath::Abs(jkink) ) {
1790
1791           // The two tracks are from the same kink
1792
1793           if (fUsedKink[TMath::Abs(ikink)-1]) continue; // skip used kinks
1794
1795           Int_t imother = -1;
1796           Int_t idaughter = -1;
1797
1798           if (ikink<0 && jkink>0) {
1799
1800             imother = iTrack;
1801             idaughter = jTrack;
1802           }
1803           else if (ikink>0 && jkink<0) {
1804
1805             imother = jTrack;
1806             idaughter = iTrack;
1807           }
1808           else {
1809             //                  cerr << "Error: Wrong combination of kink indexes: "
1810             //                       << ikink << " " << jkink << endl;
1811             continue;
1812           }
1813
1814           // Add the mother track if it passed primary track selection cuts
1815
1816           AliAODTrack * mother = NULL;
1817
1818           UInt_t selectInfo = 0;
1819           if (fTrackFilter) {
1820             selectInfo = fTrackFilter->IsSelected(esd.GetTrack(imother));
1821             if (!selectInfo) continue;
1822           }
1823
1824           if (!fUsedTrack[imother]) {
1825
1826             fUsedTrack[imother] = kTRUE;
1827
1828             AliESDtrack *esdTrackM = esd.GetTrack(imother);
1829             Double_t p[3] = { 0. };
1830             Double_t pos[3] = { 0. };
1831             esdTrackM->GetPxPyPz(p);
1832             esdTrackM->GetXYZ(pos);
1833             esdTrackM->GetCovarianceXYZPxPyPz(covTr);
1834             esdTrackM->GetESDpid(pid);
1835             esdTrackM->GetTOFLabel(tofLabel);
1836             if(fMChandler)fMChandler->SelectParticle(esdTrackM->GetLabel());
1837             mother =
1838             new(Tracks()[fNumberOfTracks++]) AliAODTrack(esdTrackM->GetID(),
1839                                                esdTrackM->GetLabel(),
1840                                                p,
1841                                                kTRUE,
1842                                                pos,
1843                                                kFALSE,
1844                                                covTr,
1845                                                (Short_t)esdTrackM->GetSign(),
1846                                                esdTrackM->GetITSClusterMap(),
1847                                                pid,
1848                                                fPrimaryVertex,
1849                                                kTRUE, // check if this is right
1850                                                vtx->UsesTrack(esdTrack->GetID()),
1851                                                AliAODTrack::kPrimary,
1852                                                selectInfo);
1853             mother->SetTPCFitMap(esdTrackM->GetTPCFitMap());
1854             mother->SetTPCClusterMap(esdTrackM->GetTPCClusterMap());
1855             mother->SetTPCSharedMap (esdTrackM->GetTPCSharedMap());
1856             mother->SetChi2perNDF(Chi2perNDF(esdTrackM));
1857             mother->SetTPCPointsF(esdTrackM->GetTPCNclsF());
1858             mother->SetTPCNCrossedRows(UShort_t(esdTrackM->GetTPCCrossedRows()));
1859             mother->SetIntegratedLength(esdTrackM->GetIntegratedLength());
1860             mother->SetTOFLabel(tofLabel);
1861             fAODTrackRefs->AddAt(mother, imother);
1862
1863             if (esdTrackM->GetSign() > 0) ++fNumberOfPositiveTracks;
1864             mother->SetFlags(esdTrackM->GetStatus());
1865             mother->ConvertAliPIDtoAODPID();
1866             fPrimaryVertex->AddDaughter(mother);
1867             mother->ConvertAliPIDtoAODPID();
1868             SetAODPID(esdTrackM,mother,detpid);
1869           }
1870           else {
1871             //                  cerr << "Error: event " << esd.GetEventNumberInFile() << " kink " << TMath::Abs(ikink)-1
1872             //                       << " track " << imother << " has already been used!" << endl;
1873           }
1874
1875           // Add the kink vertex
1876           AliESDkink * kink = esd.GetKink(TMath::Abs(ikink)-1);
1877
1878           AliAODVertex * vkink =
1879           new(Vertices()[fNumberOfVertices++]) AliAODVertex(kink->GetPosition(),
1880                                                   NULL,
1881                                                   0.,
1882                                                   mother,
1883                                                   esdTrack->GetID(),  // This is the track ID of the mother's track!
1884                                                   AliAODVertex::kKink);
1885           // Add the daughter track
1886
1887           AliAODTrack * daughter = NULL;
1888
1889           if (!fUsedTrack[idaughter]) {
1890
1891             fUsedTrack[idaughter] = kTRUE;
1892
1893             AliESDtrack *esdTrackD = esd.GetTrack(idaughter);
1894             Double_t p[3] = { 0. };
1895             Double_t pos[3] = { 0. };
1896
1897             esdTrackD->GetPxPyPz(p);
1898             esdTrackD->GetXYZ(pos);
1899             esdTrackD->GetCovarianceXYZPxPyPz(covTr);
1900             esdTrackD->GetESDpid(pid);
1901             esdTrackD->GetTOFLabel(tofLabel);
1902             selectInfo = 0;
1903             if (fTrackFilter) selectInfo = fTrackFilter->IsSelected(esdTrackD);
1904             if(fMChandler)fMChandler->SelectParticle(esdTrackD->GetLabel());
1905             daughter =
1906             new(Tracks()[fNumberOfTracks++]) AliAODTrack(esdTrackD->GetID(),
1907                                                esdTrackD->GetLabel(),
1908                                                p,
1909                                                kTRUE,
1910                                                pos,
1911                                                kFALSE,
1912                                                covTr,
1913                                                (Short_t)esdTrackD->GetSign(),
1914                                                esdTrackD->GetITSClusterMap(),
1915                                                pid,
1916                                                vkink,
1917                                                kTRUE, // check if this is right
1918                                                vtx->UsesTrack(esdTrack->GetID()),
1919                                                AliAODTrack::kSecondary,
1920                                                selectInfo);
1921             daughter->SetTPCFitMap(esdTrackD->GetTPCFitMap());
1922             daughter->SetTPCClusterMap(esdTrackD->GetTPCClusterMap());
1923             daughter->SetTPCSharedMap (esdTrackD->GetTPCSharedMap());
1924             daughter->SetTPCPointsF(esdTrackD->GetTPCNclsF());
1925             daughter->SetTPCNCrossedRows(UShort_t(esdTrackD->GetTPCCrossedRows()));
1926             daughter->SetIntegratedLength(esdTrackD->GetIntegratedLength());
1927             daughter->SetTOFLabel(tofLabel);
1928             fAODTrackRefs->AddAt(daughter, idaughter);
1929
1930             if (esdTrackD->GetSign() > 0) ++fNumberOfPositiveTracks;
1931             daughter->SetFlags(esdTrackD->GetStatus());
1932             daughter->ConvertAliPIDtoAODPID();
1933             vkink->AddDaughter(daughter);
1934             daughter->ConvertAliPIDtoAODPID();
1935             SetAODPID(esdTrackD,daughter,detpid);
1936           }
1937           else {
1938             //                  cerr << "Error: event " << esd.GetEventNumberInFile() << " kink " << TMath::Abs(ikink)-1
1939             //                       << " track " << idaughter << " has already been used!" << endl;
1940           }
1941         }
1942             }
1943     }
1944   }
1945 }
1946
1947 //______________________________________________________________________________
1948 void AliAnalysisTaskESDfilter::ConvertPrimaryVertices(const AliESDEvent& esd)
1949 {
1950   AliCodeTimerAuto("",0);
1951
1952   // Access to the AOD container of vertices
1953   fNumberOfVertices = 0;
1954
1955   Double_t pos[3] = { 0. };
1956   Double_t covVtx[6] = { 0. };
1957
1958   // Add primary vertex. The primary tracks will be defined
1959   // after the loops on the composite objects (V0, cascades, kinks)
1960   const AliESDVertex *vtx = esd.GetPrimaryVertex();
1961
1962   vtx->GetXYZ(pos); // position
1963   vtx->GetCovMatrix(covVtx); //covariance matrix
1964
1965   fPrimaryVertex = new(Vertices()[fNumberOfVertices++])
1966   AliAODVertex(pos, covVtx, vtx->GetChi2toNDF(), NULL, -1, AliAODVertex::kPrimary);
1967   fPrimaryVertex->SetName(vtx->GetName());
1968   fPrimaryVertex->SetTitle(vtx->GetTitle());
1969   fPrimaryVertex->SetBC(vtx->GetBC());
1970
1971   TString vtitle = vtx->GetTitle();
1972   if (!vtitle.Contains("VertexerTracks"))
1973     fPrimaryVertex->SetNContributors(vtx->GetNContributors());
1974
1975   if (fDebug > 0) fPrimaryVertex->Print();
1976
1977   // Add SPD "main" vertex
1978   const AliESDVertex *vtxS = esd.GetPrimaryVertexSPD();
1979   vtxS->GetXYZ(pos); // position
1980   vtxS->GetCovMatrix(covVtx); //covariance matrix
1981   AliAODVertex * mVSPD = new(Vertices()[fNumberOfVertices++])
1982   AliAODVertex(pos, covVtx, vtxS->GetChi2toNDF(), NULL, -1, AliAODVertex::kMainSPD);
1983   mVSPD->SetName(vtxS->GetName());
1984   mVSPD->SetTitle(vtxS->GetTitle());
1985   mVSPD->SetNContributors(vtxS->GetNContributors());
1986
1987   // Add SPD pileup vertices
1988   for(Int_t iV=0; iV<esd.GetNumberOfPileupVerticesSPD(); ++iV)
1989   {
1990     const AliESDVertex *vtxP = esd.GetPileupVertexSPD(iV);
1991     vtxP->GetXYZ(pos); // position
1992     vtxP->GetCovMatrix(covVtx); //covariance matrix
1993     AliAODVertex * pVSPD =  new(Vertices()[fNumberOfVertices++])
1994     AliAODVertex(pos, covVtx, vtxP->GetChi2toNDF(), NULL, -1, AliAODVertex::kPileupSPD);
1995     pVSPD->SetName(vtxP->GetName());
1996     pVSPD->SetTitle(vtxP->GetTitle());
1997     pVSPD->SetNContributors(vtxP->GetNContributors());
1998     pVSPD->SetBC(vtxP->GetBC());
1999   }
2000
2001   // Add TRK pileup vertices
2002   for(Int_t iV=0; iV<esd.GetNumberOfPileupVerticesTracks(); ++iV)
2003   {
2004     const AliESDVertex *vtxP = esd.GetPileupVertexTracks(iV);
2005     vtxP->GetXYZ(pos); // position
2006     vtxP->GetCovMatrix(covVtx); //covariance matrix
2007     AliAODVertex * pVTRK = new(Vertices()[fNumberOfVertices++])
2008     AliAODVertex(pos, covVtx, vtxP->GetChi2toNDF(), NULL, -1, AliAODVertex::kPileupTracks);
2009     pVTRK->SetName(vtxP->GetName());
2010     pVTRK->SetTitle(vtxP->GetTitle());
2011     pVTRK->SetNContributors(vtxP->GetNContributors());
2012     pVTRK->SetBC(vtxP->GetBC());
2013   }
2014
2015   // Add TPC "main" vertex
2016   const AliESDVertex *vtxT = esd.GetPrimaryVertexTPC();
2017   vtxT->GetXYZ(pos); // position
2018   vtxT->GetCovMatrix(covVtx); //covariance matrix
2019   AliAODVertex * mVTPC = new(Vertices()[fNumberOfVertices++])
2020   AliAODVertex(pos, covVtx, vtxT->GetChi2toNDF(), NULL, -1, AliAODVertex::kMainTPC);
2021   mVTPC->SetName(vtxT->GetName());
2022   mVTPC->SetTitle(vtxT->GetTitle());
2023   mVTPC->SetNContributors(vtxT->GetNContributors());
2024
2025
2026 }
2027
2028 //______________________________________________________________________________
2029 void AliAnalysisTaskESDfilter::ConvertVZERO(const AliESDEvent& esd)
2030 {
2031   // Convert VZERO data
2032   AliAODVZERO* vzeroData = AODEvent()->GetVZEROData();
2033   *vzeroData = *(esd.GetVZEROData());
2034 }
2035
2036 //______________________________________________________________________________
2037 void AliAnalysisTaskESDfilter::ConvertTZERO(const AliESDEvent& esd)
2038 {
2039   // Convert TZERO data
2040   const AliESDTZERO* esdTzero = esd.GetESDTZERO();
2041   AliAODTZERO* aodTzero = AODEvent()->GetTZEROData();
2042
2043   for (Int_t icase=0; icase<3; icase++){
2044     aodTzero->SetT0TOF(    icase, esdTzero->GetT0TOF(icase));
2045     aodTzero->SetT0TOFbest(icase, esdTzero->GetT0TOFbest(icase));
2046   }
2047   aodTzero->SetBackgroundFlag(esdTzero->GetBackgroundFlag());
2048   aodTzero->SetPileupFlag(esdTzero->GetPileupFlag());
2049   aodTzero->SetSatelliteFlag(esdTzero->GetSatellite());
2050
2051   Float_t rawTime[24];
2052   for(Int_t ipmt=0; ipmt<24; ipmt++)
2053     rawTime[ipmt] = esdTzero->GetTimeFull(ipmt,0);
2054
2055   Int_t idxOfFirstPmtA = -1,       idxOfFirstPmtC = -1;
2056   Float_t timeOfFirstPmtA = 9999, timeOfFirstPmtC = 9999;
2057   for(int ipmt=0;  ipmt<12; ipmt++){
2058     if( rawTime[ipmt] > -200 && rawTime[ipmt] < timeOfFirstPmtC && rawTime[ipmt]!=0){
2059       timeOfFirstPmtC = rawTime[ipmt];
2060       idxOfFirstPmtC  = ipmt;
2061     }
2062   }
2063   for(int ipmt=12; ipmt<24; ipmt++){
2064     if( rawTime[ipmt] > -200 && rawTime[ipmt] < timeOfFirstPmtA && rawTime[ipmt]!=0 ){
2065       timeOfFirstPmtA = rawTime[ipmt];
2066       idxOfFirstPmtA  = ipmt;
2067     }
2068   }
2069
2070   if(idxOfFirstPmtA != -1 && idxOfFirstPmtC != -1){
2071     //speed of light in cm/ns   TMath::C()*1e-7
2072     Float_t vertexraw = TMath::C()*1e-7 * (rawTime[idxOfFirstPmtA] - rawTime[idxOfFirstPmtC])/2;
2073     aodTzero->SetT0VertexRaw( vertexraw );
2074   }else{
2075     aodTzero->SetT0VertexRaw(99999);
2076   }
2077
2078   aodTzero->SetT0zVertex(esdTzero->GetT0zVertex());
2079   //amplitude
2080
2081   const Double32_t *amp=esdTzero->GetT0amplitude();
2082   for(int ipmt=0;  ipmt<24; ipmt++)
2083     aodTzero->SetAmp(ipmt, amp[ipmt]);
2084   aodTzero->SetAmp(24,esdTzero->GetMultC() );
2085   aodTzero->SetAmp(25,esdTzero->GetMultA() );
2086
2087 }
2088
2089
2090 //______________________________________________________________________________
2091 void AliAnalysisTaskESDfilter::ConvertZDC(const AliESDEvent& esd)
2092 {
2093   // Convert ZDC data
2094   AliESDZDC* esdZDC = esd.GetZDCData();
2095
2096   const Double_t zem1Energy = esdZDC->GetZEM1Energy();
2097   const Double_t zem2Energy = esdZDC->GetZEM2Energy();
2098
2099   const Double_t *towZNC = esdZDC->GetZNCTowerEnergy();
2100   const Double_t *towZPC = esdZDC->GetZPCTowerEnergy();
2101   const Double_t *towZNA = esdZDC->GetZNATowerEnergy();
2102   const Double_t *towZPA = esdZDC->GetZPATowerEnergy();
2103   const Double_t *towZNCLG = esdZDC->GetZNCTowerEnergyLR();
2104   const Double_t *towZNALG = esdZDC->GetZNATowerEnergyLR();
2105
2106   AliAODZDC* zdcAOD = AODEvent()->GetZDCData();
2107
2108   zdcAOD->SetZEM1Energy(zem1Energy);
2109   zdcAOD->SetZEM2Energy(zem2Energy);
2110   zdcAOD->SetZNCTowers(towZNC, towZNCLG);
2111   zdcAOD->SetZNATowers(towZNA, towZNALG);
2112   zdcAOD->SetZPCTowers(towZPC);
2113   zdcAOD->SetZPATowers(towZPA);
2114
2115   zdcAOD->SetZDCParticipants(esdZDC->GetZDCParticipants(), esdZDC->GetZDCPartSideA(), esdZDC->GetZDCPartSideC());
2116   zdcAOD->SetZDCImpactParameter(esdZDC->GetImpactParameter(), esdZDC->GetImpactParamSideA(),
2117         esdZDC->GetImpactParamSideC());
2118   zdcAOD->SetZDCTDCSum(esdZDC->GetZNTDCSum(0));
2119   zdcAOD->SetZDCTDCDiff(esdZDC->GetZNTDCDiff(0));
2120   if(esdZDC->IsZNChit()) zdcAOD->SetZNCTDC(esdZDC->GetZDCTDCCorrected(10,0));
2121   if(esdZDC->IsZNAhit()) zdcAOD->SetZNATDC(esdZDC->GetZDCTDCCorrected(12,0));
2122 }
2123
2124 //_______________________________________________________________________________________________________________________________________
2125 Int_t AliAnalysisTaskESDfilter::ConvertHMPID(const AliESDEvent& esd) // clm
2126 {
2127   //
2128   // Convtert ESD HMPID info to AOD and return the number of good tracks with HMPID signal.
2129   // We need to return an int since there is no signal counter in the ESD.
2130   //
2131
2132   AliCodeTimerAuto("",0);
2133
2134   Int_t cntHmpidGoodTracks = 0;
2135
2136   Float_t  xMip = 0;
2137   Float_t  yMip = 0;
2138   Int_t    qMip = 0;
2139   Int_t    nphMip = 0;
2140
2141   Float_t xTrk = 0;
2142   Float_t yTrk = 0;
2143   Float_t thetaTrk = 0;
2144   Float_t phiTrk = 0;
2145
2146   Double_t hmpPid[5]={0};
2147   Double_t hmpMom[3]={0};
2148
2149   TClonesArray &hmpidRings = *(AODEvent()->GetHMPIDrings());
2150
2151  for (Int_t iTrack=0; iTrack<esd.GetNumberOfTracks(); ++iTrack)
2152   {
2153     if(! esd.GetTrack(iTrack) ) continue;
2154
2155     if(esd.GetTrack(iTrack)->GetHMPIDsignal() > -20 ) {                  //
2156
2157       (esd.GetTrack(iTrack))->GetHMPIDmip(xMip, yMip, qMip, nphMip);     // Get MIP properties
2158       (esd.GetTrack(iTrack))->GetHMPIDtrk(xTrk,yTrk,thetaTrk,phiTrk);
2159       (esd.GetTrack(iTrack))->GetHMPIDpid(hmpPid);
2160       if((esd.GetTrack(iTrack))->GetOuterHmpParam()) (esd.GetTrack(iTrack))->GetOuterHmpPxPyPz(hmpMom);
2161
2162      if(esd.GetTrack(iTrack)->GetHMPIDsignal() == 0 && thetaTrk == 0 && qMip == 0 && nphMip ==0 ) continue;  //
2163
2164      new(hmpidRings[cntHmpidGoodTracks++]) AliAODHMPIDrings(
2165                                                                  (esd.GetTrack(iTrack))->GetID(),             // Unique track id to attach the ring to
2166                                                                  1000000*nphMip+qMip,                         // MIP charge and number of photons
2167                                                                  (esd.GetTrack(iTrack))->GetHMPIDcluIdx(),    // 1000000*chamber id + cluster idx of the assigned MIP cluster
2168                                                                  thetaTrk,                                    // track inclination angle theta
2169                                                                  phiTrk,                                      // track inclination angle phi
2170                                                                  (esd.GetTrack(iTrack))->GetHMPIDsignal(),    // Cherenkov angle
2171                                                                  (esd.GetTrack(iTrack))->GetHMPIDoccupancy(), // Occupancy claculated for the given chamber
2172                                                                  (esd.GetTrack(iTrack))->GetHMPIDchi2(),      // Ring resolution squared
2173                                                                   xTrk,                                       // Track x coordinate (LORS)
2174                                                                   yTrk,                                       // Track y coordinate (LORS)
2175                                                                   xMip,                                       // MIP x coordinate (LORS)
2176                                                                   yMip,                                       // MIP y coordinate (LORS)
2177                                                                   hmpPid,                                     // PID probablities from ESD, remove later once it is in CombinedPid
2178                                                                   hmpMom                                      // Track momentum in HMPID at ring reconstruction
2179                                                                );
2180
2181       //  Printf(Form("+++++++++ yes/no: %d  %lf %lf %lf %lf %lf %lf ",(esd.GetTrack(iTrack))->IsHMPID(),thetaTrk, (esd.GetTrack(iTrack))->GetHMPIDchi2(),xTrk, yTrk , xMip, yMip));
2182
2183
2184    }// HMPID signal > -20
2185   }//___esd track loop
2186
2187   return cntHmpidGoodTracks;
2188 }
2189
2190 //______________________________________________________________________________
2191 void AliAnalysisTaskESDfilter::ConvertESDtoAOD()
2192 {
2193   // ESD Filter analysis task executed for each event
2194
2195   AliESDEvent* esd = dynamic_cast<AliESDEvent*>(InputEvent());
2196
2197   if(!esd)return;
2198
2199   AliCodeTimerAuto("",0);
2200
2201   fOldESDformat = ( esd->GetAliESDOld() != 0x0 );
2202
2203       // Reconstruct cascades and V0 here
2204   if (fIsV0CascadeRecoEnabled) {
2205     esd->ResetCascades();
2206     esd->ResetV0s();
2207
2208     AliV0vertexer lV0vtxer;
2209     AliCascadeVertexer lCascVtxer;
2210
2211     lV0vtxer.SetCuts(fV0Cuts);
2212     lCascVtxer.SetCuts(fCascadeCuts);
2213
2214
2215     lV0vtxer.Tracks2V0vertices(esd);
2216     lCascVtxer.V0sTracks2CascadeVertices(esd);
2217   }
2218
2219
2220   fNumberOfTracks = 0;
2221   fNumberOfPositiveTracks = 0;
2222   fNumberOfV0s = 0;
2223   fNumberOfVertices = 0;
2224   fNumberOfCascades = 0;
2225   fNumberOfKinks = 0;
2226
2227   AliAODHeader* header = ConvertHeader(*esd);
2228
2229   if ( fIsVZEROEnabled ) ConvertVZERO(*esd);
2230   if ( fIsTZEROEnabled ) ConvertTZERO(*esd);
2231
2232   // Fetch Stack for debuggging if available
2233   fMChandler=0x0;
2234   if(MCEvent())
2235   {
2236     fMChandler = (AliMCEventHandler*) ((AliAnalysisManager::GetAnalysisManager())->GetMCtruthEventHandler());
2237   }
2238
2239   // loop over events and fill them
2240   // Multiplicity information needed by the header (to be revised!)
2241   Int_t nTracks    = esd->GetNumberOfTracks();
2242   for (Int_t iTrack=0; iTrack<nTracks; ++iTrack) esd->GetTrack(iTrack)->SetESDEvent(esd);
2243
2244   // Update the header
2245
2246   Int_t nV0s      = esd->GetNumberOfV0s();
2247   Int_t nCascades = esd->GetNumberOfCascades();
2248   Int_t nKinks    = esd->GetNumberOfKinks();
2249   Int_t nVertices = nV0s + nCascades /*V0 wihtin cascade already counted*/+ nKinks + 1 /* = prim. vtx*/;
2250   Int_t nPileSPDVertices=1+esd->GetNumberOfPileupVerticesSPD(); // also SPD main vertex
2251   Int_t nPileTrkVertices=esd->GetNumberOfPileupVerticesTracks();
2252   nVertices+=nPileSPDVertices;
2253   nVertices+=nPileTrkVertices;
2254   Int_t nJets     = 0;
2255   Int_t nCaloClus = esd->GetNumberOfCaloClusters();
2256   Int_t nFmdClus  = 0;
2257   Int_t nPmdClus  = esd->GetNumberOfPmdTracks();
2258   Int_t nHmpidRings = 0;
2259
2260   AliDebug(1,Form("   NV0=%d  NCASCADES=%d  NKINKS=%d", nV0s, nCascades, nKinks));
2261
2262   AODEvent()->ResetStd(nTracks, nVertices, nV0s, nCascades, nJets, nCaloClus, nFmdClus, nPmdClus,nHmpidRings);
2263
2264   if (nV0s > 0)
2265   {
2266     // RefArray to store a mapping between esd V0 number and newly created AOD-Vertex V0
2267     fAODV0VtxRefs = new TRefArray(nV0s);
2268     // RefArray to store the mapping between esd V0 number and newly created AOD-V0
2269     fAODV0Refs = new TRefArray(nV0s);
2270     // Array to take into account the V0s already added to the AOD (V0 within cascades)
2271     fUsedV0 = new Bool_t[nV0s];
2272     for (Int_t iV0=0; iV0<nV0s; ++iV0) fUsedV0[iV0]=kFALSE;
2273   }
2274
2275   if (nTracks>0)
2276   {
2277     // RefArray to store the mapping between esd track number and newly created AOD-Track
2278
2279     fAODTrackRefs = new TRefArray(nTracks);
2280
2281     // Array to take into account the tracks already added to the AOD
2282     fUsedTrack = new Bool_t[nTracks];
2283     for (Int_t iTrack=0; iTrack<nTracks; ++iTrack) fUsedTrack[iTrack]=kFALSE;
2284   }
2285
2286   // Array to take into account the kinks already added to the AOD
2287   if (nKinks>0)
2288   {
2289     fUsedKink = new Bool_t[nKinks];
2290     for (Int_t iKink=0; iKink<nKinks; ++iKink) fUsedKink[iKink]=kFALSE;
2291   }
2292
2293   ConvertPrimaryVertices(*esd);
2294
2295   //setting best TOF PID
2296   AliESDInputHandler* esdH = dynamic_cast<AliESDInputHandler*>(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
2297   if (esdH)
2298       fESDpid = esdH->GetESDpid();
2299
2300   if (fIsPidOwner && fESDpid){
2301     delete fESDpid;
2302     fESDpid = 0;
2303   }
2304   if(!fESDpid)
2305   { //in case of no Tender attached
2306     fESDpid = new AliESDpid;
2307     fIsPidOwner = kTRUE;
2308   }
2309
2310   if(!esd->GetTOFHeader())
2311   { //protection in case the pass2 LHC10b,c,d have been processed without tender.
2312     Float_t t0spread[10];
2313     Float_t intrinsicTOFres=100; //ps ok for LHC10b,c,d pass2!!
2314     for (Int_t i=0; i<10; i++) t0spread[i] = (TMath::Sqrt(esd->GetSigma2DiamondZ()))/0.03; //0.03 to convert from cm to ps
2315     fESDpid->GetTOFResponse().SetT0resolution(t0spread);
2316     fESDpid->GetTOFResponse().SetTimeResolution(intrinsicTOFres);
2317       //    fESDpid->SetTOFResponse(esd, (AliESDpid::EStartTimeType_t)fTimeZeroType);
2318     AliTOFHeader tmpTOFHeader(0,t0spread[0],0,NULL,NULL,NULL,intrinsicTOFres,t0spread[0]);
2319     AODEvent()->SetTOFHeader(&tmpTOFHeader);         // write dummy TOF header in AOD
2320   } else {
2321     AODEvent()->SetTOFHeader(esd->GetTOFHeader());    // write TOF header in AOD
2322   }
2323
2324   //  if(esd->GetTOFHeader() && fIsPidOwner) fESDpid->SetTOFResponse(esd, (AliESDpid::EStartTimeType_t)fTimeZeroType); //in case of AOD production strating form LHC10e without Tender.
2325
2326   if ( fAreCascadesEnabled ) ConvertCascades(*esd);
2327
2328   if ( fAreV0sEnabled ) ConvertV0s(*esd);
2329
2330   if ( fAreKinksEnabled ) ConvertKinks(*esd);
2331
2332   if ( fAreTracksEnabled ) ConvertTracks(*esd);
2333
2334   // Update number of AOD tracks in header at the end of track loop (M.G.)
2335   header->SetRefMultiplicity(fNumberOfTracks);
2336   header->SetRefMultiplicityPos(fNumberOfPositiveTracks);
2337   header->SetRefMultiplicityNeg(fNumberOfTracks - fNumberOfPositiveTracks);
2338
2339   if ( fTPCConstrainedFilterMask ) ConvertTPCOnlyTracks(*esd);
2340   if( fGlobalConstrainedFilterMask) ConvertGlobalConstrainedTracks(*esd);
2341
2342   if ( fArePmdClustersEnabled ) ConvertPmdClusters(*esd);
2343
2344   if ( fAreCaloClustersEnabled ) ConvertCaloClusters(*esd);
2345
2346   if ( fAreEMCALCellsEnabled )ConvertEMCALCells(*esd);
2347
2348   if ( fArePHOSCellsEnabled )ConvertPHOSCells(*esd);
2349
2350         if ( fAreEMCALTriggerEnabled )ConvertCaloTrigger(TString("EMCAL"), *esd);
2351
2352         if ( fArePHOSTriggerEnabled )ConvertCaloTrigger(TString("PHOS"), *esd);
2353
2354   if ( fAreTrackletsEnabled ) ConvertTracklets(*esd);
2355   if ( fIsZDCEnabled ) ConvertZDC(*esd);
2356
2357   if(fIsHMPIDEnabled) nHmpidRings = ConvertHMPID(*esd);
2358
2359   delete fAODTrackRefs; fAODTrackRefs=0x0;
2360   delete fAODV0VtxRefs; fAODV0VtxRefs=0x0;
2361   delete fAODV0Refs; fAODV0Refs=0x0;
2362
2363   delete[] fUsedTrack; fUsedTrack=0x0;
2364   delete[] fUsedV0; fUsedV0=0x0;
2365   delete[] fUsedKink; fUsedKink=0x0;
2366
2367   if ( fIsPidOwner){
2368     delete fESDpid;
2369     fESDpid = 0x0;
2370   }
2371
2372
2373 }
2374
2375
2376 //______________________________________________________________________________
2377 void AliAnalysisTaskESDfilter::SetAODPID(AliESDtrack *esdtrack, AliAODTrack *aodtrack, AliAODPid *detpid)
2378 {
2379   //
2380   // Setter for the raw PID detector signals
2381   //
2382
2383   // Save PID object for candidate electrons
2384     Bool_t pidSave = kFALSE;
2385     if (fTrackFilter) {
2386         Bool_t selectInfo = fTrackFilter->IsSelected((char*) "Electrons");
2387         if (selectInfo)  pidSave = kTRUE;
2388     }
2389
2390
2391     // Tracks passing pt cut
2392     if(esdtrack->Pt()>fHighPthreshold) {
2393         pidSave = kTRUE;
2394     } else {
2395         if(fPtshape){
2396             if(esdtrack->Pt()> fPtshape->GetXmin()){
2397                 Double_t y = fPtshape->Eval(esdtrack->Pt())/fPtshape->Eval(fHighPthreshold);
2398                 if(gRandom->Rndm(0)<1./y){
2399                     pidSave = kTRUE;
2400                 }//end rndm
2401             }//end if p < pmin
2402         }//end if p function
2403     }// end else
2404
2405     if (pidSave) {
2406       if(!aodtrack->GetDetPid()){// prevent memory leak when calling SetAODPID twice for the same track
2407         detpid = new AliAODPid();
2408         SetDetectorRawSignals(detpid,esdtrack);
2409         aodtrack->SetDetPID(detpid);
2410       }
2411     }
2412 }
2413
2414 //______________________________________________________________________________
2415 void AliAnalysisTaskESDfilter::SetDetectorRawSignals(AliAODPid *aodpid, AliESDtrack *track)
2416 {
2417 //
2418 //assignment of the detector signals (AliXXXesdPID inspired)
2419 //
2420  if(!track){
2421  AliInfo("no ESD track found. .....exiting");
2422  return;
2423  }
2424  // TPC momentum
2425  aodpid->SetTPCmomentum(track->GetTPCmomentum());
2426  aodpid->SetTPCTgl(track->GetTPCTgl());
2427 //  const AliExternalTrackParam *in=track->GetInnerParam();
2428 //  if (in) {
2429 //    aodpid->SetTPCmomentum(in->GetP());
2430 //  }else{
2431 //    aodpid->SetTPCmomentum(-1.);
2432 //  }
2433
2434
2435  aodpid->SetITSsignal(track->GetITSsignal());
2436  Double_t itsdedx[4]; // dE/dx samples for individual ITS layers
2437  track->GetITSdEdxSamples(itsdedx);
2438  aodpid->SetITSdEdxSamples(itsdedx);
2439
2440  aodpid->SetTPCsignal(track->GetTPCsignal());
2441  aodpid->SetTPCsignalN(track->GetTPCsignalN());
2442  if(track->GetTPCdEdxInfo()) aodpid->SetTPCdEdxInfo(track->GetTPCdEdxInfo());
2443
2444  //n TRD planes = 6
2445  Int_t nslices = track->GetNumberOfTRDslices()*6;
2446  TArrayD trdslices(nslices);
2447  for(Int_t iSl =0; iSl < track->GetNumberOfTRDslices(); iSl++) {
2448    for(Int_t iPl =0; iPl<6; iPl++) trdslices[iPl*track->GetNumberOfTRDslices()+iSl] = track->GetTRDslice(iPl,iSl);
2449  }
2450
2451 //TRD momentum
2452  for(Int_t iPl=0;iPl<6;iPl++){
2453    Double_t trdmom=track->GetTRDmomentum(iPl);
2454    aodpid->SetTRDmomentum(iPl,trdmom);
2455  }
2456
2457  aodpid->SetTRDslices(track->GetNumberOfTRDslices()*6,trdslices.GetArray());
2458  aodpid->SetTRDsignal(track->GetTRDsignal());
2459
2460  //TRD clusters and tracklets
2461  aodpid->SetTRDncls(track->GetTRDncls());
2462  aodpid->SetTRDntrackletsPID(track->GetTRDntrackletsPID());
2463
2464  aodpid->SetTRDChi2(track->GetTRDchi2());
2465
2466  //TOF PID
2467  Double_t times[AliPID::kSPECIES]; track->GetIntegratedTimes(times);
2468  aodpid->SetIntegratedTimes(times);
2469
2470    //  Float_t tzeroTrack = fESDpid->GetTOFResponse().GetStartTime(track->P());
2471    //  aodpid->SetTOFsignal(track->GetTOFsignal()-tzeroTrack);
2472    aodpid->SetTOFsignal(track->GetTOFsignal());
2473
2474   Double_t tofRes[5];
2475   for (Int_t iMass=0; iMass<5; iMass++){
2476     //    tofRes[iMass]=(Double_t)fESDpid->GetTOFResponse().GetExpectedSigma(track->P(), times[iMass], AliPID::ParticleMass(iMass));
2477     tofRes[iMass]=0; //backward compatibility
2478   }
2479   aodpid->SetTOFpidResolution(tofRes);
2480
2481 //  aodpid->SetHMPIDsignal(0); // set to zero for compression but it will be removed later
2482
2483 }
2484
2485 Double_t  AliAnalysisTaskESDfilter::Chi2perNDF(AliESDtrack* track)
2486 {
2487     // Calculate chi2 per ndf for track
2488     Int_t  nClustersTPC = track->GetTPCNcls();
2489
2490     if ( nClustersTPC > 5) {
2491        return (track->GetTPCchi2()/Float_t(nClustersTPC - 5));
2492     } else {
2493        return (-1.);
2494     }
2495  }
2496
2497
2498 //______________________________________________________________________________
2499 void AliAnalysisTaskESDfilter::Terminate(Option_t */*option*/)
2500 {
2501 // Terminate analysis
2502 //
2503     if (fDebug > 1) printf("AnalysisESDfilter: Terminate() \n");
2504 }
2505
2506 //______________________________________________________________________________
2507 void  AliAnalysisTaskESDfilter::PrintMCInfo(AliStack *pStack,Int_t label){
2508 // Print MC info
2509   if(!pStack)return;
2510   label = TMath::Abs(label);
2511   TParticle *part = pStack->Particle(label);
2512   Printf("########################");
2513   Printf("%s:%d %d UniqueID %d PDG %d P %3.3f",(char*)__FILE__,__LINE__,label,part->GetUniqueID(),part->GetPdgCode(),part->P());
2514   part->Print();
2515   TParticle* mother = part;
2516   Int_t imo = part->GetFirstMother();
2517   Int_t nprim = pStack->GetNprimary();
2518   //  while((imo >= nprim) && (mother->GetUniqueID() == 4)) {
2519   while((imo >= nprim)) {
2520     mother =  pStack->Particle(imo);
2521     Printf("Mother %s:%d Label %d UniqueID %d PDG %d P %3.3f",(char*)__FILE__,__LINE__,imo,mother->GetUniqueID(),mother->GetPdgCode(),mother->P());
2522     mother->Print();
2523     imo =  mother->GetFirstMother();
2524   }
2525   Printf("########################");
2526 }
2527
2528 //______________________________________________________
2529