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