PWG/CaloTrackCorrBase/AliCaloTrackReader.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 //_________________________________________________________________________
  17 // Base class for reading data: MonteCarlo, ESD or AOD, of PHOS EMCAL and
  18 // Central Barrel Tracking detectors (CTS).
  19 // Not all MC particles/tracks/clusters are kept, some kinematical/fiducial restrictions are done.
  20 // Mother class of : AliCaloTrackESDReader: Fills ESD data in 3 TObjArrays (PHOS, EMCAL, CTS)
  21 //                 : AliCaloTrackMCReader : Fills Kinematics data in 3 TObjArrays (PHOS, EMCAL, CTS)
  22 //                 : AliCaloTrackAODReader: Fills AOD data in 3 TObjArrays (PHOS, EMCAL, CTS)
  23 //-- Author: Gustavo Conesa (LNF-INFN)
  24 //////////////////////////////////////////////////////////////////////////////
  25
  26
  27 // --- ROOT system ---
  28 #include <TFile.h>
  29 #include <TGeoManager.h>
  30
  31 // ---- ANALYSIS system ----
  32 #include "AliMCEvent.h"
  33 #include "AliAODMCHeader.h"
  34 #include "AliGenPythiaEventHeader.h"
  35 #include "AliESDEvent.h"
  36 #include "AliAODEvent.h"
  37 #include "AliVTrack.h"
  38 #include "AliVParticle.h"
  39 #include "AliMixedEvent.h"
  40 #include "AliESDtrack.h"
  41 #include "AliESDtrackCuts.h"
  42 #include "AliTriggerAnalysis.h"
  43 #include "AliESDVZERO.h"
  44 #include "AliVCaloCells.h"
  45 #include "AliAnalysisManager.h"
  46 #include "AliInputEventHandler.h"
  47
  48 // ---- Detectors ----
  49 #include "AliPHOSGeoUtils.h"
  50 #include "AliEMCALGeometry.h"
  51 #include "AliEMCALRecoUtils.h"
  52
  53 // ---- CaloTrackCorr ---
  54 #include "AliCalorimeterUtils.h"
  55 #include "AliCaloTrackReader.h"
  56
  57 ClassImp(AliCaloTrackReader)
  58
  59
  60 //________________________________________
  61 AliCaloTrackReader::AliCaloTrackReader() :
  62 TObject(),                   fEventNumber(-1), //fCurrentFileName(""),
  63 fDataType(0),                fDebug(0),
  64 fFiducialCut(0x0),           fCheckFidCut(kFALSE),
  65 fComparePtHardAndJetPt(0),   fPtHardAndJetPtFactor(0),
  66 fComparePtHardAndClusterPt(0),fPtHardAndClusterPtFactor(0),
  67 fCTSPtMin(0),                fEMCALPtMin(0),                  fPHOSPtMin(0),
  68 fCTSPtMax(0),                fEMCALPtMax(0),                  fPHOSPtMax(0),
  69 fUseEMCALTimeCut(1),         fUseParamTimeCut(0),             fUseTrackTimeCut(0),
  70 fEMCALTimeCutMin(-10000),    fEMCALTimeCutMax(10000),
  71 fEMCALParamTimeCutMin(),     fEMCALParamTimeCutMax(),
  72 fTrackTimeCutMin(-10000),    fTrackTimeCutMax(10000),
  73 fUseTrackDCACut(0),
  74 fAODBranchList(0x0),
  75 fCTSTracks(0x0),             fEMCALClusters(0x0),             fPHOSClusters(0x0),
  76 fEMCALCells(0x0),            fPHOSCells(0x0),
  77 fInputEvent(0x0),            fOutputEvent(0x0),fMC(0x0),
  78 fFillCTS(0),                 fFillEMCAL(0),                   fFillPHOS(0),
  79 fFillEMCALCells(0),          fFillPHOSCells(0),
  80 fRecalculateClusters(kFALSE),fSelectEmbeddedClusters(kFALSE),
  81 fTrackStatus(0),             fTrackFilterMask(0),
  82 fESDtrackCuts(0),            fESDtrackComplementaryCuts(0),   fConstrainTrack(kFALSE),
  83 fSelectHybridTracks(0),      fSelectSPDHitTracks(kFALSE),
  84 fTrackMult(0),               fTrackMultEtaCut(0.8),
  85 fReadStack(kFALSE),          fReadAODMCParticles(kFALSE),
  86 fDeltaAODFileName(""),       fFiredTriggerClassName(""),
  87 fEventTriggerMask(0),        fMixEventTriggerMask(0),         fEventTriggerAtSE(0),
  88 fAnaLED(kFALSE),
  89 fTaskName(""),               fCaloUtils(0x0),
  90 fMixedEvent(NULL),           fNMixedEvent(0),                 fVertex(NULL),
  91 fListMixedTracksEvents(),    fListMixedCaloEvents(),
  92 fLastMixedTracksEvent(-1),   fLastMixedCaloEvent(-1),
  93 fWriteOutputDeltaAOD(kFALSE),fOldAOD(kFALSE),                 fCaloFilterPatch(kFALSE),
  94 fEMCALClustersListName(""),  fZvtxCut(0.),
  95 fAcceptFastCluster(kFALSE),  fRemoveLEDEvents(kTRUE),
  96 fDoEventSelection(kFALSE),   fDoV0ANDEventSelection(kFALSE),  fUseEventsWithPrimaryVertex(kFALSE),
  97 fTriggerAnalysis (0x0),      fTimeStampEventSelect(0),
  98 fTimeStampEventFracMin(0),   fTimeStampEventFracMax(0),
  99 fTimeStampRunMin(0),         fTimeStampRunMax(0),
 100 fNPileUpClusters(-1),        fNNonPileUpClusters(-1),         fNPileUpClustersCut(3),
 101 fCentralityClass(""),        fCentralityOpt(0),
 102 fEventPlaneMethod(""),       fImportGeometryFromFile(kFALSE), fImportGeometryFilePath("")
 103 {
 104   //Ctor
 105
 106   //Initialize parameters
 107   InitParameters();
 108 }
 109
 110 //_______________________________________
 111 AliCaloTrackReader::~AliCaloTrackReader()
 112 {
 113   //Dtor
 114
 115   delete fFiducialCut ;
 116
 117   if(fAODBranchList)
 118   {
 119     fAODBranchList->Delete();
 120     delete fAODBranchList ;
 121   }
 122
 123   if(fCTSTracks)
 124   {
 125     if(fDataType!=kMC)fCTSTracks->Clear() ;
 126     else              fCTSTracks->Delete() ;
 127     delete fCTSTracks ;
 128   }
 129
 130   if(fEMCALClusters)
 131   {
 132     if(fDataType!=kMC)fEMCALClusters->Clear("C") ;
 133     else              fEMCALClusters->Delete() ;
 134     delete fEMCALClusters ;
 135   }
 136
 137   if(fPHOSClusters)
 138   {
 139     if(fDataType!=kMC)fPHOSClusters->Clear("C") ;
 140     else              fPHOSClusters->Delete() ;
 141     delete fPHOSClusters ;
 142   }
 143
 144   if(fVertex)
 145   {
 146     for (Int_t i = 0; i < fNMixedEvent; i++)
 147     {
 148       delete [] fVertex[i] ;
 149
 150     }
 151     delete [] fVertex ;
 152         }
 153
 154   delete fESDtrackCuts;
 155   delete fESDtrackComplementaryCuts;
 156   delete fTriggerAnalysis;
 157
 158   //  Pointers not owned, done by the analysis frame
 159   //  if(fInputEvent)  delete fInputEvent ;
 160   //  if(fOutputEvent) delete fOutputEvent ;
 161   //  if(fMC)          delete fMC ;
 162   //  Pointer not owned, deleted by maker
 163   //  if (fCaloUtils) delete fCaloUtils ;
 164
 165 }
 166
 167 //________________________________________________
 168 Bool_t AliCaloTrackReader::ComparePtHardAndJetPt()
 169 {
 170   // Check the event, if the requested ptHard is much smaller than the jet pT, then there is a problem.
 171   // Only for PYTHIA.
 172
 173   //printf("AliCaloTrackReader::ComparePtHardAndJetPt() - GenHeaderName : %s\n",GetGenEventHeader()->ClassName());
 174
 175   if(!strcmp(GetGenEventHeader()->ClassName(), "AliGenPythiaEventHeader"))
 176   {
 177     TParticle * jet =  0;
 178     AliGenPythiaEventHeader* pygeh= (AliGenPythiaEventHeader*) GetGenEventHeader();
 179     Int_t nTriggerJets =  pygeh->NTriggerJets();
 180     Float_t ptHard = pygeh->GetPtHard();
 181
 182     if(fDebug > 1)
 183       printf("AliCaloTrackReader::ComparePtHardAndJetPt() - Njets: %d, pT Hard %f\n",nTriggerJets, ptHard);
 184
 185     Float_t tmpjet[]={0,0,0,0};
 186     for(Int_t ijet = 0; ijet< nTriggerJets; ijet++)
 187     {
 188       pygeh->TriggerJet(ijet, tmpjet);
 189       jet = new TParticle(94, 21, -1, -1, -1, -1, tmpjet[0],tmpjet[1],tmpjet[2],tmpjet[3], 0,0,0,0);
 190
 191       if(fDebug > 1)
 192         printf("AliCaloTrackReader::ComparePtHardAndJetPt() - jet %d; pycell jet pT %f\n",ijet, jet->Pt());
 193
 194       //Compare jet pT and pt Hard
 195       if(jet->Pt() > fPtHardAndJetPtFactor * ptHard)
 196       {
 197         printf("AliCaloTrackReader::ComparePtHardAndJetPt() - Reject jet event with : pT Hard %2.2f, pycell jet pT %2.2f, rejection factor %1.1f\n",
 198                 ptHard, jet->Pt(), fPtHardAndJetPtFactor);
 199         return kFALSE;
 200       }
 201     }
 202
 203     if(jet) delete jet;
 204   }
 205
 206   return kTRUE ;
 207
 208 }
 209
 210 //____________________________________________________________________
 211 Bool_t AliCaloTrackReader::ComparePtHardAndClusterPt()
 212 {
 213   // Check the event, if the requested ptHard is smaller than the calorimeter cluster E, then there is a problem.
 214   // Only for PYTHIA.
 215
 216   if(!strcmp(GetGenEventHeader()->ClassName(), "AliGenPythiaEventHeader"))
 217   {
 218     AliGenPythiaEventHeader* pygeh= (AliGenPythiaEventHeader*) GetGenEventHeader();
 219     Float_t ptHard = pygeh->GetPtHard();
 220
 221     Int_t nclusters = fInputEvent->GetNumberOfCaloClusters();
 222     for (Int_t iclus =  0; iclus <  nclusters; iclus++)
 223     {
 224       AliVCluster * clus = fInputEvent->GetCaloCluster(iclus) ;
 225       Float_t ecluster = clus->E();
 226
 227       if(ecluster > fPtHardAndClusterPtFactor*ptHard)
 228       {
 229         printf("AliCaloTrackReader::ComparePtHardAndClusterPt() - Reject : ecluster %2.2f, calo %d, factor %2.2f, ptHard %f\n",ecluster,clus->GetType(),fPtHardAndClusterPtFactor,ptHard);
 230
 231         return kFALSE;
 232       }
 233     }
 234
 235   }
 236
 237   return kTRUE ;
 238
 239 }
 240
 241 //____________________________________________
 242 AliStack* AliCaloTrackReader::GetStack() const
 243 {
 244   //Return pointer to stack
 245   if(fMC)
 246     return fMC->Stack();
 247   else
 248   {
 249     if(fDebug > 1) printf("AliCaloTrackReader::GetStack() - Stack is not available\n");
 250     return 0x0 ;
 251   }
 252 }
 253
 254 //__________________________________________________
 255 TString AliCaloTrackReader::GetFiredTriggerClasses()
 256 {
 257   // List of triggered classes in a TString
 258
 259   AliESDEvent* esdevent = dynamic_cast<AliESDEvent*> (GetInputEvent());
 260   AliAODEvent* aodevent = dynamic_cast<AliAODEvent*> (GetInputEvent());
 261
 262   if     (esdevent) return esdevent->GetFiredTriggerClasses();
 263   else if(aodevent) return aodevent->GetFiredTriggerClasses();
 264   else              return ""; // Mixed Event, MC event, does not have this trigger info
 265
 266 }
 267
 268 //______________________________________________
 269 AliHeader* AliCaloTrackReader::GetHeader() const
 270 {
 271   //Return pointer to header
 272   if(fMC)
 273   {
 274     return fMC->Header();
 275   }
 276   else
 277   {
 278     printf("AliCaloTrackReader::Header is not available\n");
 279     return 0x0 ;
 280   }
 281 }
 282
 283 //______________________________________________________________
 284 AliGenEventHeader* AliCaloTrackReader::GetGenEventHeader() const
 285 {
 286   //Return pointer to Generated event header
 287   if     (ReadStack() && fMC)
 288   {
 289     return fMC->GenEventHeader();
 290   }
 291   else if(ReadAODMCParticles() && GetAODMCHeader())
 292   {
 293     //printf("AliCaloTrackReader::GetGenEventHeader() - N headers %d\n",GetAODMCHeader()->GetNCocktailHeaders());
 294     if( GetAODMCHeader()->GetNCocktailHeaders() > 0)
 295       return GetAODMCHeader()->GetCocktailHeader(0) ;
 296     else
 297       return 0x0;
 298   }
 299   else
 300   {
 301     //printf("AliCaloTrackReader::GetGenEventHeader() - MC header not available! \n");
 302     return 0;
 303   }
 304 }
 305
 306 //____________________________________________________________________
 307 TClonesArray* AliCaloTrackReader::GetAODMCParticles(Int_t input) const
 308 {
 309   //Return list of particles in AOD. Do it for the corresponding input event.
 310
 311   TClonesArray * rv = NULL ;
 312   if(fDataType == kAOD)
 313   {
 314     if(input == 0)
 315     {
 316       //Normal input AOD
 317       AliAODEvent * evt = dynamic_cast<AliAODEvent*> (fInputEvent) ;
 318       if(evt)
 319         rv = (TClonesArray*)evt->FindListObject("mcparticles");
 320       else
 321         printf("AliCaloTrackReader::GetAODMCParticles() - wrong AOD input index? %d, or non existing tree? \n",input);
 322     }
 323
 324   }
 325   else
 326   {
 327     printf("AliCaloTrackReader::GetAODMCParticles() - Input are not AODs\n");
 328   }
 329
 330   return rv ;
 331 }
 332
 333 //________________________________________________________
 334 AliAODMCHeader* AliCaloTrackReader::GetAODMCHeader() const
 335 {
 336   //Return MC header in AOD. Do it for the corresponding input event.
 337
 338   AliAODMCHeader *mch = NULL;
 339
 340   if(fDataType == kAOD)
 341   {
 342     AliAODEvent * aod = dynamic_cast<AliAODEvent*> (fInputEvent);
 343     if(aod) mch = dynamic_cast<AliAODMCHeader*>(aod->FindListObject("mcHeader"));
 344   }
 345   else
 346   {
 347     printf("AliCaloTrackReader::GetAODMCHeader() - Input are not AODs\n");
 348   }
 349
 350   return mch;
 351 }
 352
 353 //_____________________________
 354 void AliCaloTrackReader::Init()
 355 {
 356   //Init reader. Method to be called in AliAnaPartCorrMaker
 357
 358   //printf(" AliCaloTrackReader::Init() %p \n",gGeoManager);
 359
 360   if(fReadStack && fReadAODMCParticles)
 361   {
 362     printf("AliCaloTrackReader::Init() - Cannot access stack and mcparticles at the same time, change them \n");
 363     fReadStack          = kFALSE;
 364     fReadAODMCParticles = kFALSE;
 365   }
 366
 367   // Init geometry, I do not like much to do it like this ...
 368   if(fImportGeometryFromFile && !gGeoManager)
 369   {
 370     if(fImportGeometryFilePath=="") // If not specified, set a default location
 371     fImportGeometryFilePath = "$ALICE_ROOT/OADB/EMCAL/geometry_2011.root"; // "$ALICE_ROOT/EVE/alice-data/default_geo.root"
 372
 373     printf("AliCaloTrackReader::Init() - Import %s\n",fImportGeometryFilePath.Data());
 374     TGeoManager::Import(fImportGeometryFilePath) ; // default need file "geometry.root" in local dir!!!!
 375   }
 376
 377   if(!fESDtrackCuts)
 378     fESDtrackCuts = AliESDtrackCuts::GetStandardTPCOnlyTrackCuts(); //initialize with TPC only tracks
 379
 380 }
 381
 382 //_______________________________________
 383 void AliCaloTrackReader::InitParameters()
 384 {
 385   //Initialize the parameters of the analysis.
 386   fDataType   = kESD ;
 387   fCTSPtMin   = 0.1 ;
 388   fEMCALPtMin = 0.1 ;
 389   fPHOSPtMin  = 0.1 ;
 390   fCTSPtMax   = 1000. ;
 391   fEMCALPtMax = 1000. ;
 392   fPHOSPtMax  = 1000. ;
 393
 394   //Track DCA cuts
 395   fTrackDCACut[0] = 3; //xy, quite large
 396   fTrackDCACut[1] = 3; //z, quite large
 397   fTrackDCACut[2] = 3; //TPC constrained, quite large
 398
 399   //Do not filter the detectors input by default.
 400   fFillEMCAL      = kFALSE;
 401   fFillPHOS       = kFALSE;
 402   fFillCTS        = kFALSE;
 403   fFillEMCALCells = kFALSE;
 404   fFillPHOSCells  = kFALSE;
 405
 406   fReadStack             = kFALSE; // Check in the constructor of the other readers if it was set or in the configuration file
 407   fReadAODMCParticles    = kFALSE; // Check in the constructor of the other readers if it was set or in the configuration file
 408   fDeltaAODFileName      = "deltaAODPartCorr.root";
 409   fFiredTriggerClassName = "";
 410   fEventTriggerMask      = AliVEvent::kAny;
 411   fMixEventTriggerMask   = AliVEvent::kAnyINT;
 412   fEventTriggerAtSE      = kTRUE; // Use only events that pass event selection at SE base class
 413
 414   fAcceptFastCluster = kTRUE;
 415   fAnaLED            = kFALSE;
 416
 417   //We want tracks fitted in the detectors:
 418   //fTrackStatus=AliESDtrack::kTPCrefit;
 419   //fTrackStatus|=AliESDtrack::kITSrefit;
 420   fTrackStatus     = 0;
 421   fTrackFilterMask = 128; //For AODs, but what is the difference between fTrackStatus and fTrackFilterMask?
 422
 423   fESDtrackCuts = 0;
 424   fESDtrackComplementaryCuts = 0;
 425
 426   fConstrainTrack = kFALSE ; // constrain tracks to vertex
 427
 428   fV0ADC[0] = 0;   fV0ADC[1] = 0;
 429   fV0Mul[0] = 0;   fV0Mul[1] = 0;
 430
 431   fZvtxCut   = 10.;
 432
 433   fNMixedEvent = 1;
 434
 435   fPtHardAndJetPtFactor     = 7.;
 436   fPtHardAndClusterPtFactor = 1.;
 437
 438   //Centrality
 439   fCentralityClass  = "V0M";
 440   fCentralityOpt    = 10;
 441   fCentralityBin[0] = fCentralityBin[1]=-1;
 442
 443   fEventPlaneMethod = "V0";
 444
 445   // Allocate memory (not sure this is the right place)
 446   fCTSTracks       = new TObjArray();
 447   fEMCALClusters   = new TObjArray();
 448   fPHOSClusters    = new TObjArray();
 449   fTriggerAnalysis = new AliTriggerAnalysis;
 450   fAODBranchList   = new TList ;
 451
 452   fImportGeometryFromFile = kFALSE;
 453
 454   fPileUpParamSPD[0] = 3   ; fPileUpParamSPD[1] = 0.8 ;
 455   fPileUpParamSPD[2] = 3.0 ; fPileUpParamSPD[3] = 2.0 ; fPileUpParamSPD[4] = 5.0;
 456
 457   // Parametrized time cut (LHC11d)
 458   fEMCALParamTimeCutMin[0] =-5; fEMCALParamTimeCutMin[1] =-1 ; fEMCALParamTimeCutMin[2] = 3.5 ; fEMCALParamTimeCutMin[3] = 1.  ;
 459   fEMCALParamTimeCutMax[0] = 5; fEMCALParamTimeCutMax[1] = 50; fEMCALParamTimeCutMax[2] = 0.45; fEMCALParamTimeCutMax[3] = 1.25;
 460
 461   // Parametrized time cut (LHC11c)
 462   //fEMCALParamTimeCutMin[0] =-5;   fEMCALParamTimeCutMin[1] =-1 ; fEMCALParamTimeCutMin[2] = 1.87; fEMCALParamTimeCutMin[3] = 0.4;
 463   //fEMCALParamTimeCutMax[0] = 3.5; fEMCALParamTimeCutMax[1] = 50; fEMCALParamTimeCutMax[2] = 0.15; fEMCALParamTimeCutMax[3] = 1.6;
 464
 465   fTimeStampRunMin = -1;
 466   fTimeStampRunMax = 1000000000000;
 467   fTimeStampEventFracMin = -1;
 468   fTimeStampEventFracMax = 2;
 469
 470 }
 471
 472 //___________________________________________________________
 473 Bool_t AliCaloTrackReader::IsInTimeWindow(const Double_t tof, const Float_t energy) const
 474 {
 475   // Cluster time selection window
 476
 477   // Parametrized cut depending on E
 478   if(fUseParamTimeCut)
 479   {
 480     Float_t minCut= fEMCALParamTimeCutMin[0]+fEMCALParamTimeCutMin[1]*TMath::Exp(-(energy-fEMCALParamTimeCutMin[2])/fEMCALParamTimeCutMin[3]);
 481     Float_t maxCut= fEMCALParamTimeCutMax[0]+fEMCALParamTimeCutMax[1]*TMath::Exp(-(energy-fEMCALParamTimeCutMax[2])/fEMCALParamTimeCutMax[3]);
 482     //printf("tof %f, minCut %f, maxCut %f\n",tof,minCut,maxCut);
 483     if( tof < minCut || tof > maxCut )  return kFALSE ;
 484   }
 485
 486   //In any case, the time should to be larger than the fixed window ...
 487   if( tof < fEMCALTimeCutMin  || tof > fEMCALTimeCutMax )  return kFALSE ;
 488
 489   return kTRUE ;
 490 }
 491
 492 //________________________________________________
 493 Bool_t AliCaloTrackReader::IsPileUpFromSPD() const
 494 {
 495   // Check if event is from pile-up determined by SPD
 496   // Default values: (3, 0.8, 3., 2., 5.)
 497   return fInputEvent->IsPileupFromSPD((Int_t) fPileUpParamSPD[0] , fPileUpParamSPD[1] ,
 498                                               fPileUpParamSPD[2] , fPileUpParamSPD[3] , fPileUpParamSPD[4] );
 499   //printf("Param : %d, %2.2f, %2.2f, %2.2f, %2.2f\n",(Int_t) fPileUpParamSPD[0], fPileUpParamSPD[1], fPileUpParamSPD[2], fPileUpParamSPD[3], fPileUpParamSPD[4]);
 500
 501 }
 502
 503 //__________________________________________________
 504 Bool_t AliCaloTrackReader::IsPileUpFromEMCal() const
 505 {
 506   // Check if event is from pile-up determined by EMCal
 507   if(fNPileUpClusters > fNPileUpClustersCut) return kTRUE ;
 508   else                                       return kFALSE;
 509 }
 510
 511 //________________________________________________________
 512 Bool_t AliCaloTrackReader::IsPileUpFromSPDAndEMCal() const
 513 {
 514   // Check if event is from pile-up determined by SPD and EMCal
 515   if( IsPileUpFromSPD() && IsPileUpFromEMCal()) return kTRUE ;
 516   else                                          return kFALSE;
 517 }
 518
 519 //_______________________________________________________
 520 Bool_t AliCaloTrackReader::IsPileUpFromSPDOrEMCal() const
 521 {
 522   // Check if event is from pile-up determined by SPD or EMCal
 523   if( IsPileUpFromSPD() || IsPileUpFromEMCal()) return kTRUE ;
 524   else                                          return kFALSE;
 525 }
 526
 527 //___________________________________________________________
 528 Bool_t AliCaloTrackReader::IsPileUpFromSPDAndNotEMCal() const
 529 {
 530   // Check if event is from pile-up determined by SPD and not by EMCal
 531   if( IsPileUpFromSPD() && !IsPileUpFromEMCal()) return kTRUE ;
 532   else                                          return kFALSE;
 533 }
 534
 535 //___________________________________________________________
 536 Bool_t AliCaloTrackReader::IsPileUpFromEMCalAndNotSPD() const
 537 {
 538   // Check if event is from pile-up determined by EMCal, not by SPD
 539   if( !IsPileUpFromSPD() && IsPileUpFromEMCal()) return kTRUE ;
 540   else                                           return kFALSE;
 541 }
 542
 543 //______________________________________________________________
 544 Bool_t AliCaloTrackReader::IsPileUpFromNotSPDAndNotEMCal() const
 545 {
 546   // Check if event not from pile-up determined neither by SPD nor by EMCal
 547   if( !IsPileUpFromSPD() && !IsPileUpFromEMCal()) return kTRUE ;
 548   else                                            return kFALSE;
 549 }
 550
 551 //________________________________________________________
 552 void AliCaloTrackReader::Print(const Option_t * opt) const
 553 {
 554
 555   //Print some relevant parameters set for the analysis
 556   if(! opt)
 557     return;
 558
 559   printf("***** Print: %s %s ******\n",    GetName(), GetTitle() ) ;
 560   printf("Task name      : %s\n",          fTaskName.Data()) ;
 561   printf("Data type      : %d\n",          fDataType) ;
 562   printf("CTS Min pT     : %2.1f GeV/c\n", fCTSPtMin) ;
 563   printf("EMCAL Min pT   : %2.1f GeV/c\n", fEMCALPtMin) ;
 564   printf("PHOS Min pT    : %2.1f GeV/c\n", fPHOSPtMin) ;
 565   printf("CTS Max pT     : %2.1f GeV/c\n", fCTSPtMax) ;
 566   printf("EMCAL Max pT   : %2.1f GeV/c\n", fEMCALPtMax) ;
 567   printf("PHOS Max pT    : %2.1f GeV/c\n", fPHOSPtMax) ;
 568   printf("EMCAL Time Cut: %3.1f < TOF  < %3.1f\n", fEMCALTimeCutMin, fEMCALTimeCutMax);
 569   printf("Use CTS         =     %d\n",     fFillCTS) ;
 570   printf("Use EMCAL       =     %d\n",     fFillEMCAL) ;
 571   printf("Use PHOS        =     %d\n",     fFillPHOS) ;
 572   printf("Use EMCAL Cells =     %d\n",     fFillEMCALCells) ;
 573   printf("Use PHOS  Cells =     %d\n",     fFillPHOSCells) ;
 574   printf("Track status    =     %d\n", (Int_t) fTrackStatus) ;
 575   printf("AODs Track filter mask  =  %d or hybrid %d, SPD hit %d\n", (Int_t) fTrackFilterMask,fSelectHybridTracks,fSelectSPDHitTracks) ;
 576   printf("Track Mult Eta Cut =  %d\n", (Int_t) fTrackMultEtaCut) ;
 577   printf("Write delta AOD =     %d\n",     fWriteOutputDeltaAOD) ;
 578   printf("Recalculate Clusters = %d\n",    fRecalculateClusters) ;
 579
 580   printf("Use Triggers selected in SE base class %d; If not what trigger Mask? %d; Trigger max for mixed %d \n",
 581          fEventTriggerAtSE, fEventTriggerMask,fMixEventTriggerMask);
 582
 583   if(fComparePtHardAndClusterPt)
 584           printf("Compare jet pt and pt hard to accept event, factor = %2.2f",fPtHardAndJetPtFactor);
 585
 586   if(fComparePtHardAndClusterPt)
 587           printf("Compare cluster pt and pt hard to accept event, factor = %2.2f",fPtHardAndClusterPtFactor);
 588
 589   printf("Read Kine from, stack? %d, AOD ? %d \n", fReadStack, fReadAODMCParticles) ;
 590   printf("Delta AOD File Name =     %s\n", fDeltaAODFileName.Data()) ;
 591   printf("Centrality: Class %s, Option %d, Bin [%d,%d] \n", fCentralityClass.Data(),fCentralityOpt,fCentralityBin[0], fCentralityBin[1]) ;
 592
 593   printf("    \n") ;
 594
 595 }
 596
 597 //_________________________________________________________________________
 598 Bool_t AliCaloTrackReader::FillInputEvent(const Int_t iEntry,
 599                                           const char * /*currentFileName*/)
 600 {
 601   //Fill the event counter and input lists that are needed, called by the analysis maker.
 602
 603   fEventNumber = iEntry;
 604   //fCurrentFileName = TString(currentFileName);
 605   if(!fInputEvent)
 606   {
 607           if(fDebug >= 0) printf("AliCaloTrackReader::FillInputEvent() - Input event not available, skip event analysis\n");
 608           return kFALSE;
 609   }
 610
 611   //Select events only fired by a certain trigger configuration if it is provided
 612   Int_t eventType = 0;
 613   if(fInputEvent->GetHeader())
 614           eventType = ((AliVHeader*)fInputEvent->GetHeader())->GetEventType();
 615
 616   if (GetFiredTriggerClasses().Contains("FAST")  && !GetFiredTriggerClasses().Contains("ALL") && !fAcceptFastCluster)
 617   {
 618     if(fDebug > 0)  printf("AliCaloTrackReader::FillInputEvent - Do not count events from fast cluster, trigger name %s\n",fFiredTriggerClassName.Data());
 619     return kFALSE;
 620   }
 621
 622   //-------------------------------------------------------------------------------------
 623   // Reject event if large clusters with large energy
 624   // Use only for LHC11a data for the moment, and if input is clusterizer V1 or V1+unfolding
 625   // If clusterzer NxN or V2 it does not help
 626   //-------------------------------------------------------------------------------------
 627   Int_t run = fInputEvent->GetRunNumber();
 628   if( fRemoveLEDEvents && run > 146857  && run < 146861 )
 629   {
 630     //printf("Event %d\n",GetEventNumber());
 631
 632     // Count number of cells with energy larger than 0.1 in SM3, cut on this number
 633     Int_t ncellsSM3 = 0;
 634     for(Int_t icell = 0; icell < fInputEvent->GetEMCALCells()->GetNumberOfCells(); icell++)
 635     {
 636       Int_t absID = fInputEvent->GetEMCALCells()->GetCellNumber(icell);
 637       Int_t sm    = GetCaloUtils()->GetEMCALGeometry()->GetSuperModuleNumber(absID);
 638       if(fInputEvent->GetEMCALCells()->GetAmplitude(icell) > 0.1 && sm==3) ncellsSM3++;
 639     }
 640
 641     Int_t ncellcut = 21;
 642     if(fFiredTriggerClassName.Contains("EMC")) ncellcut = 35;
 643
 644     if(ncellsSM3 >= ncellcut)
 645     {
 646       if(fDebug > 0) printf(" AliCaloTrackReader::FillInputEvent() - reject event with ncells in SM3 %d\n",ncellsSM3);
 647       return kFALSE;
 648     }
 649   }// Remove LED events
 650
 651   // Reject pure LED events?
 652   if( fFiredTriggerClassName  !="" && !fAnaLED)
 653   {
 654     //printf("Event type %d\n",eventType);
 655     if(eventType!=7)
 656       return kFALSE; //Only physics event, do not use for simulated events!!!
 657
 658     if(fDebug > 0)
 659       printf("AliCaloTrackReader::FillInputEvent() - FiredTriggerClass <%s>, selected class <%s>, compare name %d\n",
 660              GetFiredTriggerClasses().Data(),fFiredTriggerClassName.Data(), GetFiredTriggerClasses().Contains(fFiredTriggerClassName));
 661
 662     if( !GetFiredTriggerClasses().Contains(fFiredTriggerClassName) ) return kFALSE;
 663     else if(fDebug > 0) printf("AliCaloTrackReader::FillInputEvent() - Accepted triggered event\n");
 664   }
 665   else if(fAnaLED)
 666   {
 667     //    kStartOfRun =       1,    // START_OF_RUN
 668     //    kEndOfRun =         2,    // END_OF_RUN
 669     //    kStartOfRunFiles =  3,    // START_OF_RUN_FILES
 670     //    kEndOfRunFiles =    4,    // END_OF_RUN_FILES
 671     //    kStartOfBurst =     5,    // START_OF_BURST
 672     //    kEndOfBurst =       6,    // END_OF_BURST
 673     //    kPhysicsEvent =     7,    // PHYSICS_EVENT
 674     //    kCalibrationEvent = 8,    // CALIBRATION_EVENT
 675     //    kFormatError =      9,    // EVENT_FORMAT_ERROR
 676     //    kStartOfData =      10,   // START_OF_DATA
 677     //    kEndOfData =        11,   // END_OF_DATA
 678     //    kSystemSoftwareTriggerEvent   = 12, // SYSTEM_SOFTWARE_TRIGGER_EVENT
 679     //    kDetectorSoftwareTriggerEvent = 13  // DETECTOR_SOFTWARE_TRIGGER_EVENT
 680
 681           if(eventType!=7 && fDebug > 1 )printf("AliCaloTrackReader::FillInputEvent() - DO LED, Event Type <%d>, 8 Calibration \n",  eventType);
 682           if(eventType!=8)return kFALSE;
 683   }
 684
 685   //In case of analysis of events with jets, skip those with jet pt > 5 pt hard
 686   if(fComparePtHardAndJetPt)
 687   {
 688     if(!ComparePtHardAndJetPt()) return kFALSE ;
 689   }
 690
 691   if(fComparePtHardAndClusterPt)
 692   {
 693     if(!ComparePtHardAndClusterPt()) return kFALSE ;
 694   }
 695
 696   //Fill Vertex array
 697   FillVertexArray();
 698   //Reject events with Z vertex too large, only for SE analysis, if not, cut on the analysis code
 699   if(!GetMixedEvent() && TMath::Abs(fVertex[0][2]) > fZvtxCut) return kFALSE;
 700
 701   //------------------------------------------------------
 702   //Event rejection depending on vertex, pileup, v0and
 703   //------------------------------------------------------
 704   if(fDataType==kESD && fTimeStampEventSelect)
 705   {
 706     AliESDEvent* esd = dynamic_cast<AliESDEvent*> (fInputEvent);
 707     Int_t timeStamp = esd->GetTimeStamp();
 708     Float_t timeStampFrac = 1.*(timeStamp-fTimeStampRunMin) / (fTimeStampRunMax-fTimeStampRunMin);
 709
 710     //printf("stamp0 %d, max0 %d, frac %f\n", timeStamp-fTimeStampRunMin,fTimeStampRunMax-fTimeStampRunMin, timeStampFrac);
 711
 712     if(timeStampFrac < fTimeStampEventFracMin || timeStampFrac > fTimeStampEventFracMax) return kFALSE;
 713
 714     //printf("\t accept time stamp\n");
 715   }
 716
 717
 718   //------------------------------------------------------
 719   //Event rejection depending on vertex, pileup, v0and
 720   //------------------------------------------------------
 721
 722   if(fUseEventsWithPrimaryVertex)
 723   {
 724     if( !CheckForPrimaryVertex() )              return kFALSE;
 725     if( TMath::Abs(fVertex[0][0] ) < 1.e-6 &&
 726         TMath::Abs(fVertex[0][1] ) < 1.e-6 &&
 727         TMath::Abs(fVertex[0][2] ) < 1.e-6    ) return kFALSE;
 728   }
 729
 730   //printf("Reader : IsPileUp %d, Multi %d\n",IsPileUpFromSPD(),fInputEvent->IsPileupFromSPDInMultBins());
 731
 732   if(fDoEventSelection)
 733   {
 734     if(!fCaloFilterPatch)
 735     {
 736       // Do not analyze events with pileup
 737       Bool_t bPileup = IsPileUpFromSPD();
 738       //IsPileupFromSPDInMultBins() // method to try
 739       //printf("pile-up %d, %d, %2.2f, %2.2f, %2.2f, %2.2f\n",bPileup, (Int_t) fPileUpParamSPD[0], fPileUpParamSPD[1], fPileUpParamSPD[2], fPileUpParamSPD[3], fPileUpParamSPD[4]);
 740       if(bPileup) return kFALSE;
 741
 742       if(fDoV0ANDEventSelection)
 743       {
 744         Bool_t bV0AND = kTRUE;
 745         AliESDEvent* esd = dynamic_cast<AliESDEvent*> (fInputEvent);
 746         if(esd)
 747           bV0AND = fTriggerAnalysis->IsOfflineTriggerFired(esd, AliTriggerAnalysis::kV0AND);
 748         //else bV0AND = //FIXME FOR AODs
 749         if(!bV0AND) return kFALSE;
 750       }
 751     }//CaloFilter patch
 752     else
 753     {
 754       if(fInputEvent->GetNumberOfCaloClusters() > 0)
 755       {
 756         AliVCluster * calo = fInputEvent->GetCaloCluster(0);
 757         if(calo->GetNLabels() == 4)
 758         {
 759           Int_t * selection = calo->GetLabels();
 760           Bool_t bPileup = selection[0];
 761           if(bPileup) return kFALSE;
 762
 763           Bool_t bGoodV = selection[1];
 764           if(fUseEventsWithPrimaryVertex && !bGoodV) return kFALSE;
 765
 766           if(fDoV0ANDEventSelection)
 767           {
 768             Bool_t bV0AND = selection[2];
 769             if(!bV0AND) return kFALSE;
 770           }
 771
 772           fTrackMult = selection[3];
 773           if(fTrackMult == 0) return kFALSE;
 774         }
 775         else
 776         {
 777           //First filtered AODs, track multiplicity stored there.
 778           fTrackMult = (Int_t) ((AliAODHeader*)fInputEvent->GetHeader())->GetCentrality();
 779           if(fTrackMult == 0) return kFALSE;
 780         }
 781       }//at least one cluster
 782       else
 783       {
 784         //printf("AliCaloTrackReader::FillInputEvent() - No clusters in event\n");
 785         //Remove events with  vertex (0,0,0), bad vertex reconstruction
 786         if(fUseEventsWithPrimaryVertex && TMath::Abs(fVertex[0][0]) < 1.e-6 && TMath::Abs(fVertex[0][1]) < 1.e-6 && TMath::Abs(fVertex[0][2]) < 1.e-6) return kFALSE;
 787
 788         //First filtered AODs, track multiplicity stored there.
 789         fTrackMult = (Int_t) ((AliAODHeader*)fInputEvent->GetHeader())->GetCentrality();
 790         if(fTrackMult == 0) return kFALSE;
 791       }// no cluster
 792     }// CaloFileter patch
 793   }// Event selection
 794   //------------------------------------------------------
 795
 796   //Check if there is a centrality value, PbPb analysis, and if a centrality bin selection is requested
 797   //If we need a centrality bin, we select only those events in the corresponding bin.
 798   if(GetCentrality() && fCentralityBin[0]>=0 && fCentralityBin[1]>=0 && fCentralityOpt==100)
 799   {
 800     Int_t cen = GetEventCentrality();
 801     if(cen > fCentralityBin[1] || cen < fCentralityBin[0]) return kFALSE; //reject events out of bin.
 802   }
 803
 804   //Fill the arrays with cluster/tracks/cells data
 805
 806   if(fFillEMCALCells)
 807     FillInputEMCALCells();
 808
 809   if(fFillPHOSCells)
 810     FillInputPHOSCells();
 811
 812   FillInputVZERO();
 813
 814   if(fEventTriggerAtSE)
 815   {
 816     if(fFillCTS)
 817     {
 818       FillInputCTS();
 819       //Accept events with at least one track
 820       if(fTrackMult == 0 && fDoEventSelection) return kFALSE;
 821     }
 822     if(fFillEMCAL)
 823       FillInputEMCAL();
 824     if(fFillPHOS)
 825       FillInputPHOS();
 826   }
 827   else
 828   {
 829     // In case of mixing analysis, all triggers accepted, but trigger particles selected
 830     // only for the specific trigered events selected here. Mixing done only for MB events,
 831     // tracks array filled also for those events and not the others.
 832
 833     AliAnalysisManager *manager = AliAnalysisManager::GetAnalysisManager();
 834     AliInputEventHandler *inputHandler = dynamic_cast<AliInputEventHandler*>(manager->GetInputEventHandler());
 835
 836     if(!inputHandler) return kFALSE ;  // to content coverity
 837
 838     UInt_t isTrigger = inputHandler->IsEventSelected() & fEventTriggerMask;
 839     UInt_t isMB      = inputHandler->IsEventSelected() & fMixEventTriggerMask;
 840
 841     if(fFillCTS && (isMB || isTrigger))
 842     {
 843       FillInputCTS();
 844       //Accept events with at least one track
 845       if(fTrackMult == 0 && fDoEventSelection) return kFALSE;
 846     }
 847
 848     if(isTrigger)
 849     {
 850       if(fFillEMCAL)
 851         FillInputEMCAL();
 852       if(fFillPHOS)
 853         FillInputPHOS();
 854     }
 855   }
 856
 857   return kTRUE ;
 858 }
 859
 860 //___________________________________
 861 void AliCaloTrackReader::ResetLists()
 862 {
 863   //  Reset lists, called by the analysis maker
 864
 865   if(fCTSTracks)       fCTSTracks     -> Clear();
 866   if(fEMCALClusters)   fEMCALClusters -> Clear("C");
 867   if(fPHOSClusters)    fPHOSClusters  -> Clear("C");
 868
 869   fV0ADC[0] = 0;   fV0ADC[1] = 0;
 870   fV0Mul[0] = 0;   fV0Mul[1] = 0;
 871
 872 }
 873
 874 //____________________________________________________________
 875 void AliCaloTrackReader::SetInputEvent(AliVEvent* const input)
 876 {
 877   fInputEvent  = input;
 878   fMixedEvent = dynamic_cast<AliMixedEvent*>(GetInputEvent()) ;
 879   if (fMixedEvent)
 880     fNMixedEvent = fMixedEvent->GetNumberOfEvents() ;
 881
 882   //Delete previous vertex
 883   if(fVertex)
 884   {
 885     for (Int_t i = 0; i < fNMixedEvent; i++)
 886     {
 887       delete [] fVertex[i] ;
 888     }
 889     delete [] fVertex ;
 890   }
 891
 892   fVertex = new Double_t*[fNMixedEvent] ;
 893   for (Int_t i = 0; i < fNMixedEvent; i++)
 894   {
 895     fVertex[i] = new Double_t[3] ;
 896     fVertex[i][0] = 0.0 ;
 897     fVertex[i][1] = 0.0 ;
 898     fVertex[i][2] = 0.0 ;
 899   }
 900 }
 901
 902 //__________________________________________________
 903 Int_t AliCaloTrackReader::GetEventCentrality() const
 904 {
 905   //Return current event centrality
 906
 907   if(GetCentrality())
 908   {
 909     if     (fCentralityOpt==100) return (Int_t) GetCentrality()->GetCentralityPercentile(fCentralityClass); // 100 bins max
 910     else if(fCentralityOpt==10)  return GetCentrality()->GetCentralityClass10(fCentralityClass);// 10 bins max
 911     else if(fCentralityOpt==20)  return GetCentrality()->GetCentralityClass5(fCentralityClass); // 20 bins max
 912     else
 913     {
 914       printf("AliCaloTrackReader::GetEventCentrality() - Unknown centrality option %d, use 10, 20 or 100\n",fCentralityOpt);
 915       return -1;
 916     }
 917   }
 918   else return -1;
 919
 920 }
 921
 922 //_____________________________________________________
 923 Double_t AliCaloTrackReader::GetEventPlaneAngle() const
 924 {
 925   //Return current event centrality
 926
 927   if(GetEventPlane())
 928   {
 929     Float_t ep =  GetEventPlane()->GetEventplane(GetEventPlaneMethod(), GetInputEvent());
 930
 931     if(GetEventPlaneMethod()=="Q" && (ep < 0 || ep > TMath::Pi()))
 932     {
 933       if(fDebug > 0 ) printf("AliCaloTrackReader::GetEventPlaneAngle() -  Bad EP for <Q> method : %f\n",ep);
 934       return -1000;
 935     }
 936     else if(GetEventPlaneMethod().Contains("V0")  )
 937     {
 938       if((ep > TMath::Pi()/2 || ep < -TMath::Pi()/2))
 939       {
 940         if(fDebug > 0 ) printf("AliCaloTrackReader::GetEventPlaneAngle() -  Bad EP for <%s> method : %f\n",GetEventPlaneMethod().Data(), ep);
 941         return -1000;
 942       }
 943
 944       ep+=TMath::Pi()/2; // put same range as for <Q> method
 945
 946     }
 947
 948     //printf("AliCaloTrackReader::GetEventPlaneAngle() = %f\n",ep);
 949     if(fDebug > 0 )
 950     {
 951       if     (ep > TMath::Pi()) printf("AliCaloTrackReader::GetEventPlaneAngle() - Too large angle = %f\n",ep);
 952       else if(ep < 0          ) printf("AliCaloTrackReader::GetEventPlaneAngle() - Negative angle = %f\n" ,ep);
 953     }
 954
 955     return ep;
 956   }
 957   else
 958   {
 959     if(fDataType!=kMC && fDebug > 0) printf("AliCaloTrackReader::GetEventPlaneAngle() -  No EP pointer\n");
 960     return -1000;
 961   }
 962
 963 }
 964
 965 //__________________________________________________________
 966 void AliCaloTrackReader::GetVertex(Double_t vertex[3]) const
 967 {
 968   //Return vertex position to be used for single event analysis
 969   vertex[0]=fVertex[0][0];
 970   vertex[1]=fVertex[0][1];
 971   vertex[2]=fVertex[0][2];
 972 }
 973
 974 //____________________________________________________________
 975 void AliCaloTrackReader::GetVertex(Double_t vertex[3],
 976                                    const Int_t evtIndex) const
 977 {
 978   //Return vertex position for mixed event, recover the vertex in a particular event.
 979
 980   vertex[0]=fVertex[evtIndex][0];  vertex[1]=fVertex[evtIndex][1];  vertex[2]=fVertex[evtIndex][2];
 981
 982 }
 983
 984 //________________________________________
 985 void AliCaloTrackReader::FillVertexArray()
 986 {
 987
 988   //Fill data member with vertex
 989   //In case of Mixed event, multiple vertices
 990
 991   //Delete previous vertex
 992   if(fVertex)
 993   {
 994     for (Int_t i = 0; i < fNMixedEvent; i++)
 995     {
 996       delete [] fVertex[i] ;
 997     }
 998     delete [] fVertex ;
 999   }
1000
1001   fVertex = new Double_t*[fNMixedEvent] ;
1002   for (Int_t i = 0; i < fNMixedEvent; i++)
1003   {
1004     fVertex[i] = new Double_t[3] ;
1005     fVertex[i][0] = 0.0 ;
1006     fVertex[i][1] = 0.0 ;
1007     fVertex[i][2] = 0.0 ;
1008   }
1009
1010   if (!fMixedEvent)
1011   { //Single event analysis
1012     if(fDataType!=kMC)
1013     {
1014
1015       if(fInputEvent->GetPrimaryVertex())
1016       {
1017         fInputEvent->GetPrimaryVertex()->GetXYZ(fVertex[0]);
1018       }
1019       else
1020       {
1021         printf("AliCaloTrackReader::FillVertexArray() - NULL primary vertex\n");
1022         fVertex[0][0]=0.;   fVertex[0][1]=0.;   fVertex[0][2]=0.;
1023       }//Primary vertex pointer do not exist
1024
1025     } else
1026     {//MC read event
1027       fVertex[0][0]=0.;   fVertex[0][1]=0.;   fVertex[0][2]=0.;
1028     }
1029
1030     if(fDebug > 1)
1031       printf("AliCaloTrackReader::FillVertexArray() - Single Event Vertex : %f,%f,%f\n",fVertex[0][0],fVertex[0][1],fVertex[0][2]);
1032
1033   } else
1034   { // MultiEvent analysis
1035     for (Int_t iev = 0; iev < fNMixedEvent; iev++)
1036     {
1037       if (fMixedEvent->GetVertexOfEvent(iev))
1038         fMixedEvent->GetVertexOfEvent(iev)->GetXYZ(fVertex[iev]);
1039       else
1040       { // no vertex found !!!!
1041         AliWarning("No vertex found");
1042       }
1043
1044       if(fDebug > 1)
1045         printf("AliCaloTrackReader::FillVertexArray() - Multi Event %d Vertex : %f,%f,%f\n",iev,fVertex[iev][0],fVertex[iev][1],fVertex[iev][2]);
1046
1047     }
1048   }
1049
1050 }
1051
1052 //_____________________________________
1053 void AliCaloTrackReader::FillInputCTS()
1054 {
1055   //Return array with Central Tracking System (CTS) tracks
1056
1057   if(fDebug > 2 ) printf("AliCaloTrackReader::FillInputCTS()\n");
1058
1059   Double_t pTrack[3] = {0,0,0};
1060
1061   Int_t nTracks = fInputEvent->GetNumberOfTracks() ;
1062   fTrackMult    = 0;
1063   Int_t nstatus = 0;
1064   Double_t bz   = GetInputEvent()->GetMagneticField();
1065
1066   for(Int_t i = 0; i < 19; i++)
1067   {
1068     fTrackBCEvent   [i] = 0;
1069     fTrackBCEventCut[i] = 0;
1070   }
1071
1072   for (Int_t itrack =  0; itrack <  nTracks; itrack++)
1073   {////////////// track loop
1074     AliVTrack * track = (AliVTrack*)fInputEvent->GetTrack(itrack) ; // retrieve track from esd
1075
1076     //Select tracks under certain conditions, TPCrefit, ITSrefit ... check the set bits
1077     ULong_t status = track->GetStatus();
1078
1079     if (fTrackStatus && !((status & fTrackStatus) == fTrackStatus))
1080       continue ;
1081
1082     nstatus++;
1083
1084     if     (fDataType==kESD)
1085     {
1086       AliESDtrack* esdTrack = dynamic_cast<AliESDtrack*> (track);
1087
1088       if(esdTrack)
1089       {
1090         if(fESDtrackCuts->AcceptTrack(esdTrack))
1091         {
1092           track->GetPxPyPz(pTrack) ;
1093
1094           if(fConstrainTrack)
1095           {
1096             if(esdTrack->GetConstrainedParam())
1097             {
1098               const AliExternalTrackParam* constrainParam = esdTrack->GetConstrainedParam();
1099               esdTrack->Set(constrainParam->GetX(),constrainParam->GetAlpha(),constrainParam->GetParameter(),constrainParam->GetCovariance());
1100               esdTrack->GetConstrainedPxPyPz(pTrack);
1101             }
1102             else continue;
1103
1104           } // use constrained tracks
1105
1106           if(fSelectSPDHitTracks)
1107           {//Not much sense to use with TPC only or Hybrid tracks
1108             if(!esdTrack->HasPointOnITSLayer(0) && !esdTrack->HasPointOnITSLayer(1)) continue ;
1109           }
1110         }
1111         // Complementary track to global : Hybrids (make sure that the previous selection is for Global)
1112         else  if(fESDtrackComplementaryCuts && fESDtrackComplementaryCuts->AcceptTrack(esdTrack))
1113         {
1114           // constrain the track
1115           if(esdTrack->GetConstrainedParam())
1116           {
1117             esdTrack->Set(esdTrack->GetConstrainedParam()->GetX(),esdTrack->GetConstrainedParam()->GetAlpha(),esdTrack->GetConstrainedParam()->GetParameter(),esdTrack->GetConstrainedParam()->GetCovariance());
1118
1119             track->GetPxPyPz(pTrack) ;
1120
1121           }
1122           else continue;
1123         }
1124         else continue;
1125       }
1126     } // ESD
1127     else if(fDataType==kAOD)
1128     {
1129       AliAODTrack *aodtrack = dynamic_cast <AliAODTrack*>(track);
1130
1131       if(aodtrack)
1132       {
1133        if(fDebug > 2 ) printf("AliCaloTrackReader::FillInputCTS():AOD track type: %d (primary %d), hybrid? %d \n",
1134                               aodtrack->GetType(),AliAODTrack::kPrimary,
1135                               aodtrack->IsHybridGlobalConstrainedGlobal());
1136
1137
1138         if (fSelectHybridTracks)
1139         {
1140           if (!aodtrack->IsHybridGlobalConstrainedGlobal())       continue ;
1141         }
1142         else
1143         {
1144           if ( aodtrack->TestFilterBit(fTrackFilterMask)==kFALSE) continue ;
1145         }
1146
1147         if(fSelectSPDHitTracks)
1148         {//Not much sense to use with TPC only or Hybrid tracks
1149           if(!aodtrack->HasPointOnITSLayer(0) && !aodtrack->HasPointOnITSLayer(1)) continue ;
1150         }
1151
1152         if (aodtrack->GetType()!= AliAODTrack::kPrimary)          continue ;
1153
1154         if (fDebug > 2 ) printf("AliCaloTrackReader::FillInputCTS(): \t accepted track! \n");
1155
1156         track->GetPxPyPz(pTrack) ;
1157
1158       } // aod track exists
1159       else continue ;
1160
1161     } // AOD
1162
1163     //Count the tracks in eta < 0.9
1164     //printf("Eta %f cut  %f\n",TMath::Abs(track->Eta()),fTrackMultEtaCut);
1165     if(TMath::Abs(track->Eta())< fTrackMultEtaCut) fTrackMult++;
1166
1167     TLorentzVector momentum(pTrack[0],pTrack[1],pTrack[2],0);
1168
1169     Bool_t okTOF = ( (status & AliVTrack::kTOFout) == AliVTrack::kTOFout ) ;
1170     Double_t tof = -1000;
1171     //Int_t    bc  = -1000;
1172     Int_t trackBC = -1000 ;
1173
1174     if(okTOF)
1175     {
1176       trackBC = track->GetTOFBunchCrossing(bz);
1177       SetTrackEventBC(trackBC+9);
1178
1179       tof = track->GetTOFsignal()*1e-3;
1180       //printf("track TOF %e\n",tof);
1181       //bc = TMath::Nint((tof-25)/50) + 9;
1182       //printf("track pt %f, tof %2.2f, bc=%d\n",track->Pt(),tof,bc);
1183       //SetTrackEventBC(bc);
1184
1185     }
1186
1187     if(fCTSPtMin > momentum.Pt() || fCTSPtMax < momentum.Pt()) continue ;
1188
1189     if(fCheckFidCut && !fFiducialCut->IsInFiducialCut(momentum,"CTS")) continue;
1190
1191     if(okTOF)
1192     {
1193        //SetTrackEventBCcut(bc);
1194       SetTrackEventBCcut(trackBC+9);
1195
1196       //In any case, the time should to be larger than the fixed window ...
1197       if( fUseTrackTimeCut && (trackBC!=0 || tof < fTrackTimeCutMin  || tof > fTrackTimeCutMax) )
1198       {
1199         //printf("Remove track time %f and bc = %d\n",tof,trackBC);
1200         continue ;
1201       }
1202       //else printf("Accept track time %f and bc = %d\n",tof,trackBC);
1203
1204     }
1205
1206     if(fUseTrackDCACut)
1207     {
1208       //In case of hybrid tracks on AODs, constrained TPC tracks cannot be propagated back to primary vertex
1209       AliAODTrack * aodTrack = dynamic_cast<AliAODTrack*>(track);
1210       Float_t dcaCons = -999;
1211       if(aodTrack)
1212       {
1213         dcaCons = aodTrack->DCA();
1214         //vtxBC   = aodTrack->GetProdVertex()->GetBC();
1215         if(dcaCons > -999)
1216         {
1217           if(TMath::Abs(dcaCons) > fTrackDCACut[2] ) continue ;
1218         }
1219       }
1220
1221       //non contrained TPC tracks (tracks with ITS points) on AODs
1222       if(dcaCons==-999)
1223       {
1224         Double_t dca[2]   = {1e6,1e6};
1225         Double_t covar[3] = {1e6,1e6,1e6};
1226         Bool_t okDCA = track->PropagateToDCA(fInputEvent->GetPrimaryVertex(),bz,100.,dca,covar);
1227         if( !okDCA                               ||
1228             TMath::Abs(dca[0]) > fTrackDCACut[0] ||
1229             TMath::Abs(dca[1]) > fTrackDCACut[1]    ) continue ;
1230       }
1231     }// DCA cuts
1232
1233     if(fDebug > 2 && momentum.Pt() > 0.1)
1234       printf("AliCaloTrackReader::FillInputCTS() - Selected tracks E %3.2f, pt %3.2f, phi %3.2f, eta %3.2f\n",
1235              momentum.E(),momentum.Pt(),momentum.Phi()*TMath::RadToDeg(),momentum.Eta());
1236
1237     if (fMixedEvent)  track->SetID(itrack);
1238
1239     fCTSTracks->Add(track);
1240
1241   }// track loop
1242
1243   if(fDebug > 1)
1244     printf("AliCaloTrackReader::FillInputCTS()   - aod entries %d, input tracks %d, pass status %d, multipliticy %d\n", fCTSTracks->GetEntriesFast(), nTracks, nstatus, fTrackMult);//fCTSTracksNormalInputEntries);
1245
1246 }
1247
1248 //__________________________________________________________________
1249 void AliCaloTrackReader::FillInputEMCALAlgorithm(AliVCluster * clus,
1250                                                  const Int_t iclus)
1251 {
1252   //Fill the EMCAL data in the array, do it
1253
1254   Int_t vindex = 0 ;
1255   if (fMixedEvent)
1256     vindex = fMixedEvent->EventIndexForCaloCluster(iclus);
1257
1258   //Reject clusters with bad channels, close to borders and exotic;
1259   if(!GetCaloUtils()->GetEMCALRecoUtils()->IsGoodCluster(clus,GetCaloUtils()->GetEMCALGeometry(),GetEMCALCells(),fInputEvent->GetBunchCrossNumber())) return;
1260
1261   //Mask all cells in collumns facing ALICE thick material if requested
1262   if(GetCaloUtils()->GetNMaskCellColumns())
1263   {
1264     Int_t absId   = -1;
1265     Int_t iSupMod = -1;
1266     Int_t iphi    = -1;
1267     Int_t ieta    = -1;
1268     Bool_t shared = kFALSE;
1269     GetCaloUtils()->GetEMCALRecoUtils()->GetMaxEnergyCell(GetCaloUtils()->GetEMCALGeometry(), GetEMCALCells(),clus,absId,iSupMod,ieta,iphi,shared);
1270     if(GetCaloUtils()->MaskFrameCluster(iSupMod, ieta)) return;
1271   }
1272
1273   if(fSelectEmbeddedClusters)
1274   {
1275     if(clus->GetNLabels()==0 || clus->GetLabel() < 0) return;
1276     //else printf("Embedded cluster,  %d, n label %d label %d  \n",iclus,clus->GetNLabels(),clus->GetLabel());
1277   }
1278
1279   //Float_t pos[3];
1280   //clus->GetPosition(pos);
1281   //printf("Before Corrections: e %f, x %f, y %f, z %f\n",clus->E(),pos[0],pos[1],pos[2]);
1282
1283   if(fRecalculateClusters)
1284   {
1285     //Recalibrate the cluster energy
1286     if(GetCaloUtils()->IsRecalibrationOn())
1287     {
1288       Float_t energy = GetCaloUtils()->RecalibrateClusterEnergy(clus, GetEMCALCells());
1289
1290       clus->SetE(energy);
1291       //printf("Recalibrated Energy %f\n",clus->E());
1292
1293       GetCaloUtils()->RecalculateClusterShowerShapeParameters(GetEMCALCells(),clus);
1294       GetCaloUtils()->RecalculateClusterPID(clus);
1295
1296     } // recalculate E
1297
1298     //Recalculate distance to bad channels, if new list of bad channels provided
1299     GetCaloUtils()->RecalculateClusterDistanceToBadChannel(GetEMCALCells(),clus);
1300
1301     //Recalculate cluster position
1302     if(GetCaloUtils()->IsRecalculationOfClusterPositionOn())
1303     {
1304       GetCaloUtils()->RecalculateClusterPosition(GetEMCALCells(),clus);
1305       //clus->GetPosition(pos);
1306       //printf("After  Corrections: e %f, x %f, y %f, z %f\n",clus->E(),pos[0],pos[1],pos[2]);
1307     }
1308
1309     // Recalculate TOF
1310     if(GetCaloUtils()->GetEMCALRecoUtils()->IsTimeRecalibrationOn())
1311     {
1312       Double_t tof      = clus->GetTOF();
1313       Float_t  frac     =-1;
1314       Int_t    absIdMax = GetCaloUtils()->GetMaxEnergyCell(fEMCALCells, clus,frac);
1315
1316       if(fDataType==AliCaloTrackReader::kESD)
1317       {
1318         tof = fEMCALCells->GetCellTime(absIdMax);
1319       }
1320
1321       GetCaloUtils()->GetEMCALRecoUtils()->RecalibrateCellTime(absIdMax,fInputEvent->GetBunchCrossNumber(),tof);
1322
1323       clus->SetTOF(tof);
1324
1325     }// Time recalibration
1326   }
1327
1328   //Correct non linearity
1329   if(GetCaloUtils()->IsCorrectionOfClusterEnergyOn())
1330   {
1331     GetCaloUtils()->CorrectClusterEnergy(clus) ;
1332     //printf("Linearity Corrected Energy %f\n",clus->E());
1333
1334     //In case of MC analysis, to match resolution/calibration in real data
1335     Float_t rdmEnergy = GetCaloUtils()->GetEMCALRecoUtils()->SmearClusterEnergy(clus);
1336     // printf("\t Energy %f, smeared %f\n", clus->E(),rdmEnergy);
1337     clus->SetE(rdmEnergy);
1338   }
1339
1340   Double_t tof = clus->GetTOF()*1e9;
1341
1342   Int_t bc = TMath::Nint(tof/50) + 9;
1343   //printf("tof %2.2f, bc+5=%d\n",tof,bc);
1344
1345   SetEMCalEventBC(bc);
1346
1347   if(fEMCALPtMin > clus->E() || fEMCALPtMax < clus->E()) return ;
1348
1349   TLorentzVector momentum ;
1350
1351   clus->GetMomentum(momentum, fVertex[vindex]);
1352
1353   if(fCheckFidCut && !fFiducialCut->IsInFiducialCut(momentum,"EMCAL")) return ;
1354
1355   SetEMCalEventBCcut(bc);
1356
1357   if(!IsInTimeWindow(tof,clus->E()))
1358   {
1359     fNPileUpClusters++ ;
1360     if(fUseEMCALTimeCut) return ;
1361   }
1362   else
1363     fNNonPileUpClusters++;
1364
1365   if(fDebug > 2 && momentum.E() > 0.1)
1366     printf("AliCaloTrackReader::FillInputEMCAL() - Selected clusters E %3.2f, pt %3.2f, phi %3.2f, eta %3.2f\n",
1367            momentum.E(),momentum.Pt(),momentum.Phi()*TMath::RadToDeg(),momentum.Eta());
1368
1369   if (fMixedEvent)
1370     clus->SetID(iclus) ;
1371
1372   fEMCALClusters->Add(clus);
1373
1374 }
1375
1376 //_______________________________________
1377 void AliCaloTrackReader::FillInputEMCAL()
1378 {
1379   //Return array with EMCAL clusters in aod format
1380
1381   if(fDebug > 2 ) printf("AliCaloTrackReader::FillInputEMCAL()\n");
1382
1383   // First recalibrate cells, time or energy
1384   //  if(GetCaloUtils()->IsRecalibrationOn())
1385   //    GetCaloUtils()->GetEMCALRecoUtils()->RecalibrateCells(GetCaloUtils()->GetEMCALGeometry(),
1386   //                                                          GetEMCALCells(),
1387   //                                                          fInputEvent->GetBunchCrossNumber());
1388
1389   fNPileUpClusters    = 0; // Init counter
1390   fNNonPileUpClusters = 0; // Init counter
1391   for(Int_t i = 0; i < 19; i++)
1392   {
1393     fEMCalBCEvent   [i] = 0;
1394     fEMCalBCEventCut[i] = 0;
1395   }
1396
1397   //Loop to select clusters in fiducial cut and fill container with aodClusters
1398   if(fEMCALClustersListName=="")
1399   {
1400     Int_t nclusters = fInputEvent->GetNumberOfCaloClusters();
1401     for (Int_t iclus =  0; iclus <  nclusters; iclus++)
1402     {
1403       AliVCluster * clus = 0;
1404       if ( (clus = fInputEvent->GetCaloCluster(iclus)) )
1405       {
1406         if (IsEMCALCluster(clus))
1407         {
1408           FillInputEMCALAlgorithm(clus, iclus);
1409         }//EMCAL cluster
1410       }// cluster exists
1411     }// cluster loop
1412
1413     //Recalculate track matching
1414     GetCaloUtils()->RecalculateClusterTrackMatching(fInputEvent);
1415
1416   }//Get the clusters from the input event
1417   else
1418   {
1419     TClonesArray * clusterList = 0x0;
1420
1421     if      (fInputEvent->FindListObject(fEMCALClustersListName))
1422     {
1423       clusterList = dynamic_cast<TClonesArray*> (fInputEvent->FindListObject(fEMCALClustersListName));
1424     }
1425     else if(fOutputEvent)
1426     {
1427       clusterList = dynamic_cast<TClonesArray*> (fOutputEvent->FindListObject(fEMCALClustersListName));
1428     }
1429
1430     if(!clusterList)
1431     {
1432         printf("AliCaloTrackReader::FillInputEMCAL() - Wrong name of list with clusters?  <%s>\n",fEMCALClustersListName.Data());
1433         return;
1434     }
1435
1436     Int_t nclusters = clusterList->GetEntriesFast();
1437     for (Int_t iclus =  0; iclus <  nclusters; iclus++)
1438     {
1439       AliVCluster * clus = dynamic_cast<AliVCluster*> (clusterList->At(iclus));
1440       //printf("E %f\n",clus->E());
1441       if (clus) FillInputEMCALAlgorithm(clus, iclus);
1442       else printf("AliCaloTrackReader::FillInputEMCAL() - Null cluster in list!\n");
1443     }// cluster loop
1444
1445     // Recalculate the pile-up time, in case long time clusters removed during clusterization
1446     //printf("Input event INIT : Pile-up clusters %d, NO pile-up %d\n",fNPileUpClusters,fNNonPileUpClusters);
1447
1448     fNPileUpClusters    = 0; // Init counter
1449     fNNonPileUpClusters = 0; // Init counter
1450     for(Int_t i = 0; i < 19; i++)
1451     {
1452       fEMCalBCEvent   [i] = 0;
1453       fEMCalBCEventCut[i] = 0;
1454     }
1455
1456     for (Int_t iclus =  0; iclus < fInputEvent->GetNumberOfCaloClusters(); iclus++)
1457     {
1458       AliVCluster * clus = 0;
1459
1460       if ( (clus = fInputEvent->GetCaloCluster(iclus)) )
1461       {
1462         if (IsEMCALCluster(clus))
1463         {
1464
1465           Float_t  frac     =-1;
1466           Int_t    absIdMax = GetCaloUtils()->GetMaxEnergyCell(fEMCALCells, clus,frac);
1467           Double_t tof = clus->GetTOF();
1468           GetCaloUtils()->GetEMCALRecoUtils()->RecalibrateCellTime(absIdMax,fInputEvent->GetBunchCrossNumber(),tof);
1469           tof*=1e9;
1470
1471           //printf("Input event cluster : AbsIdMax %d, E %2.2f, time %2.2f \n", absIdMax,clus->E(),tof);
1472
1473           //Reject clusters with bad channels, close to borders and exotic;
1474           if(!GetCaloUtils()->GetEMCALRecoUtils()->IsGoodCluster(clus,GetCaloUtils()->GetEMCALGeometry(),GetEMCALCells(),fInputEvent->GetBunchCrossNumber()))  continue;
1475
1476           Int_t bc = TMath::Nint(tof/50) + 9;
1477           SetEMCalEventBC(bc);
1478
1479           if(fEMCALPtMin > clus->E() || fEMCALPtMax < clus->E()) continue ;
1480
1481           TLorentzVector momentum ;
1482
1483           clus->GetMomentum(momentum, fVertex[0]);
1484
1485           if(fCheckFidCut && !fFiducialCut->IsInFiducialCut(momentum,"EMCAL")) return ;
1486
1487           SetEMCalEventBCcut(bc);
1488
1489           if(!IsInTimeWindow(tof,clus->E()))
1490             fNPileUpClusters++ ;
1491           else
1492             fNNonPileUpClusters++;
1493
1494         }
1495       }
1496     }
1497
1498     //printf("Input event : Pile-up clusters %d, NO pile-up %d\n",fNPileUpClusters,fNNonPileUpClusters);
1499
1500     // Recalculate track matching, not necessary if already done in the reclusterization task.
1501     // in case it was not done ...
1502     GetCaloUtils()->RecalculateClusterTrackMatching(fInputEvent,clusterList);
1503
1504   }
1505
1506   if(fDebug > 1) printf("AliCaloTrackReader::FillInputEMCAL() - aod entries %d, n pile-up clusters %d, n non pile-up %d \n",  fEMCALClusters->GetEntriesFast(),fNPileUpClusters,fNNonPileUpClusters);
1507
1508 }
1509
1510 //______________________________________
1511 void AliCaloTrackReader::FillInputPHOS()
1512 {
1513   //Return array with PHOS clusters in aod format
1514
1515   if(fDebug > 2 ) printf("AliCaloTrackReader::FillInputPHOS()\n");
1516
1517   //Loop to select clusters in fiducial cut and fill container with aodClusters
1518   Int_t nclusters = fInputEvent->GetNumberOfCaloClusters();
1519   for (Int_t iclus = 0; iclus < nclusters; iclus++)
1520   {
1521     AliVCluster * clus = 0;
1522     if ( (clus = fInputEvent->GetCaloCluster(iclus)) )
1523     {
1524       if (IsPHOSCluster(clus))
1525       {
1526         //Check if the cluster contains any bad channel and if close to calorimeter borders
1527         Int_t vindex = 0 ;
1528         if (fMixedEvent)
1529           vindex = fMixedEvent->EventIndexForCaloCluster(iclus);
1530         if( GetCaloUtils()->ClusterContainsBadChannel("PHOS",clus->GetCellsAbsId(), clus->GetNCells()))
1531           continue;
1532         if(!GetCaloUtils()->CheckCellFiducialRegion(clus, fInputEvent->GetPHOSCells(), fInputEvent, vindex))
1533           continue;
1534
1535         if(fRecalculateClusters)
1536         {
1537           //Recalibrate the cluster energy
1538           if(GetCaloUtils()->IsRecalibrationOn())
1539           {
1540             Float_t energy = GetCaloUtils()->RecalibrateClusterEnergy(clus, (AliAODCaloCells*)GetPHOSCells());
1541             clus->SetE(energy);
1542           }
1543         }
1544
1545         TLorentzVector momentum ;
1546
1547         clus->GetMomentum(momentum, fVertex[vindex]);
1548
1549         if(fCheckFidCut && !fFiducialCut->IsInFiducialCut(momentum,"PHOS")) continue;
1550
1551         if(fPHOSPtMin > momentum.E() || fPHOSPtMax < momentum.E())          continue;
1552
1553         if(fDebug > 2 && momentum.E() > 0.1)
1554           printf("AliCaloTrackReader::FillInputPHOS() - Selected clusters E %3.2f, pt %3.2f, phi %3.2f, eta %3.2f\n",
1555                  momentum.E(),momentum.Pt(),momentum.Phi()*TMath::RadToDeg(),momentum.Eta());
1556
1557
1558         if (fMixedEvent)
1559         {
1560           clus->SetID(iclus) ;
1561         }
1562
1563         fPHOSClusters->Add(clus);
1564
1565       }//PHOS cluster
1566     }//cluster exists
1567   }//esd cluster loop
1568
1569   if(fDebug > 1) printf("AliCaloTrackReader::FillInputPHOS()  - aod entries %d\n",  fPHOSClusters->GetEntriesFast());
1570
1571 }
1572
1573 //____________________________________________
1574 void AliCaloTrackReader::FillInputEMCALCells()
1575 {
1576   //Return array with EMCAL cells in aod format
1577
1578   fEMCALCells = fInputEvent->GetEMCALCells();
1579
1580 }
1581
1582 //___________________________________________
1583 void AliCaloTrackReader::FillInputPHOSCells()
1584 {
1585   //Return array with PHOS cells in aod format
1586
1587   fPHOSCells = fInputEvent->GetPHOSCells();
1588
1589 }
1590
1591 //_______________________________________
1592 void AliCaloTrackReader::FillInputVZERO()
1593 {
1594   //Fill VZERO information in data member, add all the channels information.
1595   AliVVZERO* v0 = fInputEvent->GetVZEROData();
1596   //printf("Init V0: ADC (%d,%d), Multiplicity (%d,%d) \n",fV0ADC[0],fV0ADC[1],fV0Mul[0],fV0Mul[1]);
1597
1598   if (v0)
1599   {
1600     AliESDVZERO* esdV0 = dynamic_cast<AliESDVZERO*> (v0);
1601     for (Int_t i = 0; i < 32; i++)
1602     {
1603       if(esdV0)
1604       {//Only available in ESDs
1605         fV0ADC[0] += (Int_t)esdV0->GetAdcV0C(i);
1606         fV0ADC[1] += (Int_t)esdV0->GetAdcV0A(i);
1607       }
1608
1609       fV0Mul[0] += (Int_t)v0->GetMultiplicityV0C(i);
1610       fV0Mul[1] += (Int_t)v0->GetMultiplicityV0A(i);
1611     }
1612     if(fDebug > 0)
1613       printf("V0: ADC (%d,%d), Multiplicity (%d,%d) \n",fV0ADC[0],fV0ADC[1],fV0Mul[0],fV0Mul[1]);
1614   }
1615   else
1616   {
1617     if(fDebug > 0)
1618       printf("Cannot retrieve V0 ESD! Run w/ null V0 charges\n ");
1619   }
1620 }
1621
1622
1623 //___________________________________________________________________
1624 Bool_t AliCaloTrackReader::IsEMCALCluster(AliVCluster* cluster) const
1625 {
1626   // Check if it is a cluster from EMCAL. For old AODs cluster type has
1627   // different number and need to patch here
1628
1629   if(fDataType==kAOD && fOldAOD)
1630   {
1631     if (cluster->GetType() == 2) return kTRUE;
1632     else                         return kFALSE;
1633   }
1634   else
1635   {
1636     return cluster->IsEMCAL();
1637   }
1638
1639 }
1640
1641 //___________________________________________________________________
1642 Bool_t AliCaloTrackReader::IsPHOSCluster(AliVCluster * cluster) const
1643 {
1644   //Check if it is a cluster from PHOS.For old AODs cluster type has
1645   // different number and need to patch here
1646
1647   if(fDataType==kAOD && fOldAOD)
1648   {
1649     Int_t type = cluster->GetType();
1650     if (type == 0 || type == 1) return kTRUE;
1651     else                        return kFALSE;
1652   }
1653   else
1654   {
1655     return cluster->IsPHOS();
1656   }
1657
1658 }
1659
1660 //________________________________________________
1661 Bool_t AliCaloTrackReader::CheckForPrimaryVertex()
1662 {
1663   //Check if the vertex was well reconstructed, copy from V0Reader of conversion group
1664   //Only for ESDs ...
1665
1666   AliESDEvent * event = dynamic_cast<AliESDEvent*> (fInputEvent);
1667   if(!event) return kTRUE;
1668
1669   if(event->GetPrimaryVertexTracks()->GetNContributors() > 0)
1670   {
1671     return kTRUE;
1672   }
1673
1674   if(event->GetPrimaryVertexTracks()->GetNContributors() < 1)
1675   {
1676     // SPD vertex
1677     if(event->GetPrimaryVertexSPD()->GetNContributors() > 0)
1678     {
1679       //cout<<"spd vertex type::"<< fESDEvent->GetPrimaryVertex()->GetName() << endl;
1680       return kTRUE;
1681
1682     }
1683     if(event->GetPrimaryVertexSPD()->GetNContributors() < 1)
1684     {
1685       //      cout<<"bad vertex type::"<< event->GetPrimaryVertex()->GetName() << endl;
1686       return kFALSE;
1687     }
1688   }
1689
1690   return kFALSE;
1691
1692 }
1693
1694 //____________________________________________________________
1695 void  AliCaloTrackReader::SetTrackCuts(AliESDtrackCuts * cuts)
1696 {
1697   // Set Track cuts
1698
1699   if(fESDtrackCuts) delete fESDtrackCuts ;
1700
1701   fESDtrackCuts = cuts ;
1702
1703 }
1704
1705 //_________________________________________________________________________
1706 void  AliCaloTrackReader::SetTrackComplementaryCuts(AliESDtrackCuts * cuts)
1707 {
1708   // Set Track cuts for complementary tracks (hybrids)
1709
1710   if(fESDtrackComplementaryCuts) delete fESDtrackComplementaryCuts ;
1711
1712   fESDtrackComplementaryCuts = cuts ;
1713
1714 }
1715
1716