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