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