PWGGA/GammaConv/AliCaloPhotonCuts.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
   3  *                                                                                                                *
   4  * Authors: Friederike Bock, Baldo Sahlmueller                                *
   5  * Version 1.0                                                                                    *
   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 //---------------------------------------------
  18 // Class handling all kinds of selection cuts for
  19 // Photon from EMCAL clusters
  20 //---------------------------------------------
  21 ////////////////////////////////////////////////
  22
  23 #include "AliCaloPhotonCuts.h"
  24 #include "AliAnalysisManager.h"
  25 #include "AliInputEventHandler.h"
  26 #include "AliMCEventHandler.h"
  27 #include "AliAODHandler.h"
  28 #include "TH1.h"
  29 #include "TH2.h"
  30 #include "TF1.h"
  31 #include "AliStack.h"
  32 #include "AliAODConversionMother.h"
  33 #include "AliAODConversionPhoton.h"
  34 #include "TObjString.h"
  35 #include "AliAODEvent.h"
  36 #include "AliESDEvent.h"
  37 #include "AliCentrality.h"
  38 #include "TList.h"
  39 #include "TFile.h"
  40 #include "AliLog.h"
  41 #include "AliV0ReaderV1.h"
  42 #include "AliAODMCParticle.h"
  43 #include "AliAODMCHeader.h"
  44 #include "AliPicoTrack.h"
  45 #include "AliEMCALRecoUtils.h"
  46
  47 class iostream;
  48
  49 using namespace std;
  50
  51 ClassImp(AliCaloPhotonCuts)
  52
  53
  54 const char* AliCaloPhotonCuts::fgkCutNames[AliCaloPhotonCuts::kNCuts] = {
  55         "ClusterType",          //0   0: all,    1: EMCAL,   2: PHOS
  56         "EtaMin",               //1   0: -10,    1: -0.6687, 2: -0,5, 3: -2
  57         "EtaMax",               //2   0: 10,     1: 0.66465, 2: 0.5,  3: 2
  58         "PhiMin",               //3   0: -10000, 1: 1.39626
  59         "PhiMax",               //4   0: 10000, 1: 3.125
  60         "DistanceToBadChannel", //5   0: 0,      1: 5
  61         "Timing",               //6   0: no cut
  62         "TrackMatching",        //7   0: 0,      1: 5
  63         "ExoticCell",           //8   0: no cut
  64         "MinEnergy",            //9   0: no cut, 1: 0.05,    2: 0.1,  3: 0.15, 4: 0.2, 5: 0.3, 6: 0.5, 7: 0.75, 8: 1, 9: 1.25 (all GeV)
  65         "MinNCells",            //10  0: no cut, 1: 1,       2: 2,    3: 3,    4: 4,   5: 5,   6: 6
  66         "MinM02",                               //11
  67         "MaxM02",                               //12
  68         "MinM20",                               //13
  69         "MaxM20",                               //14
  70         "MaximumDispersion",    //15
  71         "NLM"                                   //16
  72 };
  73
  74
  75 //________________________________________________________________________
  76 AliCaloPhotonCuts::AliCaloPhotonCuts(const char *name,const char *title) :
  77         AliAnalysisCuts(name,title),
  78         fHistograms(NULL),
  79         fClusterType(0),
  80         fMinEtaCut(-10),
  81         fMaxEtaCut(10),
  82         fUseEtaCut(0),
  83         fMinPhiCut(-10000),
  84         fMaxPhiCut(-10000),
  85         fUsePhiCut(0),
  86         fMinDistanceToBadChannel(0),
  87         fUseDistanceToBadChannel(0),
  88         fMaxTimeDiff(10e10),
  89         fUseTimeDiff(0),
  90         fMinDistTrackToCluster(0),
  91         fUseDistTrackToCluster(0),
  92         fExoticCell(0),
  93         fUseExoticCell(0),
  94         fMinEnergy(0),
  95         fUseMinEnergy(0),
  96         fMinNCells(0),
  97         fUseNCells(0),
  98         fMaxM02(1000),
  99         fMinM02(0),
 100         fUseM02(0),
 101         fMaxM20(1000),
 102         fMinM20(0),
 103         fUseM20(0),
 104         fMaxDispersion(1000),
 105         fUseDispersion(0),
 106         fMinNLM(0),
 107         fMaxNLM(1000),
 108         fUseNLM(0),
 109         fCutString(NULL),
 110         fHistCutIndex(NULL),
 111         fHistAcceptanceCuts(NULL),
 112         fHistClusterIdentificationCuts(NULL),
 113         fHistClusterEtavsPhiBeforeAcc(NULL),
 114         fHistClusterEtavsPhiAfterAcc(NULL),
 115         fHistClusterEtavsPhiAfterQA(NULL),
 116         fHistDistanceToBadChannelBeforeAcc(NULL),
 117         fHistDistanceToBadChannelAfterAcc(NULL),
 118         fHistClusterTimevsEBeforeQA(NULL),
 119         fHistClusterTimevsEAfterQA(NULL),
 120         fHistExoticCellBeforeQA(NULL),
 121         fHistExoticCellAfterQA(NULL),
 122         fHistNMatchedTracks(NULL),
 123         fHistDistanceTrackToClusterBeforeQA(NULL),
 124         fHistDistanceTrackToClusterAfterQA(NULL),
 125         fHistEnergyOfClusterBeforeQA(NULL),
 126         fHistEnergyOfClusterAfterQA(NULL),
 127         fHistNCellsBeforeQA(NULL),
 128         fHistNCellsAfterQA(NULL),
 129         fHistM02BeforeQA(NULL),
 130         fHistM02AfterQA(NULL),
 131         fHistM20BeforeQA(NULL),
 132         fHistM20AfterQA(NULL),
 133         fHistDispersionBeforeQA(NULL),
 134         fHistDispersionAfterQA(NULL),
 135         fHistNLMBeforeQA(NULL),
 136         fHistNLMAfterQA(NULL)
 137 {
 138    for(Int_t jj=0;jj<kNCuts;jj++){fCuts[jj]=0;}
 139    fCutString=new TObjString((GetCutNumber()).Data());
 140 }
 141
 142 //________________________________________________________________________
 143 AliCaloPhotonCuts::AliCaloPhotonCuts(const AliCaloPhotonCuts &ref) :
 144    AliAnalysisCuts(ref),
 145         fHistograms(NULL),
 146         fClusterType(ref.fClusterType),
 147         fMinEtaCut(ref.fMinEtaCut),
 148         fMaxEtaCut(ref.fMaxEtaCut),
 149         fUseEtaCut(ref.fUseEtaCut),
 150         fMinPhiCut(ref.fMinPhiCut),
 151         fMaxPhiCut(ref.fMaxPhiCut),
 152         fUsePhiCut(ref.fUsePhiCut),
 153         fMinDistanceToBadChannel(ref.fMinDistanceToBadChannel),
 154         fUseDistanceToBadChannel(ref.fUseDistanceToBadChannel),
 155         fMaxTimeDiff(ref.fMaxTimeDiff),
 156         fUseTimeDiff(ref.fUseTimeDiff),
 157         fMinDistTrackToCluster(ref.fMinDistTrackToCluster),
 158         fUseDistTrackToCluster(ref.fUseDistTrackToCluster),
 159         fExoticCell(ref.fExoticCell),
 160         fUseExoticCell(ref.fUseExoticCell),
 161         fMinEnergy(ref.fMinEnergy),
 162         fUseMinEnergy(ref.fUseMinEnergy),
 163         fMinNCells(ref.fMinNCells),
 164         fUseNCells(ref.fUseNCells),
 165         fMaxM02(ref.fMaxM02),
 166         fMinM02(ref.fMinM02),
 167         fUseM02(ref.fUseM02),
 168         fMaxM20(ref.fMaxM20),
 169         fMinM20(ref.fMinM20),
 170         fUseM20(ref.fUseDispersion),
 171         fMaxDispersion(ref.fMaxDispersion),
 172         fUseDispersion(ref.fUseDispersion),
 173         fMinNLM(ref.fMinNLM),
 174         fMaxNLM(ref.fMaxNLM),
 175         fUseNLM(ref.fUseNLM),
 176         fCutString(NULL),
 177         fHistCutIndex(NULL),
 178         fHistAcceptanceCuts(NULL),
 179         fHistClusterIdentificationCuts(NULL),
 180         fHistClusterEtavsPhiBeforeAcc(NULL),
 181         fHistClusterEtavsPhiAfterAcc(NULL),
 182         fHistClusterEtavsPhiAfterQA(NULL),
 183         fHistDistanceToBadChannelBeforeAcc(NULL),
 184         fHistDistanceToBadChannelAfterAcc(NULL),
 185         fHistClusterTimevsEBeforeQA(NULL),
 186         fHistClusterTimevsEAfterQA(NULL),
 187         fHistExoticCellBeforeQA(NULL),
 188         fHistExoticCellAfterQA(NULL),
 189         fHistNMatchedTracks(NULL),
 190         fHistDistanceTrackToClusterBeforeQA(NULL),
 191         fHistDistanceTrackToClusterAfterQA(NULL),
 192         fHistEnergyOfClusterBeforeQA(NULL),
 193         fHistEnergyOfClusterAfterQA(NULL),
 194         fHistNCellsBeforeQA(NULL),
 195         fHistNCellsAfterQA(NULL),
 196         fHistM02BeforeQA(NULL),
 197         fHistM02AfterQA(NULL),
 198         fHistM20BeforeQA(NULL),
 199         fHistM20AfterQA(NULL),
 200         fHistDispersionBeforeQA(NULL),
 201         fHistDispersionAfterQA(NULL),
 202         fHistNLMBeforeQA(NULL),
 203         fHistNLMAfterQA(NULL)
 204 {
 205    // Copy Constructor
 206    for(Int_t jj=0;jj<kNCuts;jj++){fCuts[jj]=ref.fCuts[jj];}
 207    fCutString=new TObjString((GetCutNumber()).Data());
 208
 209 }
 210
 211
 212 //________________________________________________________________________
 213 AliCaloPhotonCuts::~AliCaloPhotonCuts() {
 214    // Destructor
 215    //Deleting fHistograms leads to seg fault it it's added to output collection of a task
 216    // if(fHistograms)
 217    //    delete fHistograms;
 218    // fHistograms = NULL;
 219    if(fCutString != NULL){
 220       delete fCutString;
 221       fCutString = NULL;
 222    }
 223 }
 224
 225 //________________________________________________________________________
 226 void AliCaloPhotonCuts::InitCutHistograms(TString name){
 227
 228         // Initialize Cut Histograms for QA (only initialized and filled if function is called)
 229         TH1::AddDirectory(kFALSE);
 230
 231         if(fHistograms != NULL){
 232                 delete fHistograms;
 233                 fHistograms=NULL;
 234         }
 235         if(fHistograms==NULL){
 236                 fHistograms=new TList();
 237                 fHistograms->SetOwner(kTRUE);
 238                 if(name=="")fHistograms->SetName(Form("CaloCuts_%s",GetCutNumber().Data()));
 239                 else fHistograms->SetName(Form("%s_%s",name.Data(),GetCutNumber().Data()));
 240         }
 241
 242         // IsPhotonSelected
 243         fHistCutIndex=new TH1F(Form("IsPhotonSelected %s",GetCutNumber().Data()),"IsPhotonSelected",5,-0.5,4.5);
 244         fHistCutIndex->GetXaxis()->SetBinLabel(kPhotonIn+1,"in");
 245         fHistCutIndex->GetXaxis()->SetBinLabel(kDetector+1,"detector");
 246         fHistCutIndex->GetXaxis()->SetBinLabel(kAcceptance+1,"acceptance");
 247         fHistCutIndex->GetXaxis()->SetBinLabel(kClusterQuality+1,"cluster QA");
 248         fHistCutIndex->GetXaxis()->SetBinLabel(kPhotonOut+1,"out");
 249         fHistograms->Add(fHistCutIndex);
 250
 251         // Acceptance Cuts
 252         fHistAcceptanceCuts=new TH1F(Form("AcceptanceCuts %s",GetCutNumber().Data()),"AcceptanceCuts",5,-0.5,4.5);
 253         fHistAcceptanceCuts->GetXaxis()->SetBinLabel(1,"in");
 254         fHistAcceptanceCuts->GetXaxis()->SetBinLabel(2,"eta");
 255         fHistAcceptanceCuts->GetXaxis()->SetBinLabel(3,"phi");
 256         fHistAcceptanceCuts->GetXaxis()->SetBinLabel(4,"distance to bad channel");
 257         fHistAcceptanceCuts->GetXaxis()->SetBinLabel(5,"out");
 258         fHistograms->Add(fHistAcceptanceCuts);
 259
 260         // Cluster Cuts
 261         fHistClusterIdentificationCuts =new TH1F(Form("ClusterQualityCuts %s",GetCutNumber().Data()),"ClusterQualityCuts",11,-0.5,10.5);
 262         fHistClusterIdentificationCuts->GetXaxis()->SetBinLabel(1,"in");
 263         fHistClusterIdentificationCuts->GetXaxis()->SetBinLabel(2,"timing");
 264         fHistClusterIdentificationCuts->GetXaxis()->SetBinLabel(3,"track matching");
 265         fHistClusterIdentificationCuts->GetXaxis()->SetBinLabel(4,"Exotics");
 266         fHistClusterIdentificationCuts->GetXaxis()->SetBinLabel(5,"minimum energy");
 267         fHistClusterIdentificationCuts->GetXaxis()->SetBinLabel(6,"minimum NCells");
 268         fHistClusterIdentificationCuts->GetXaxis()->SetBinLabel(7,"M02");
 269         fHistClusterIdentificationCuts->GetXaxis()->SetBinLabel(8,"M20");
 270         fHistClusterIdentificationCuts->GetXaxis()->SetBinLabel(9,"dispersion");
 271         fHistClusterIdentificationCuts->GetXaxis()->SetBinLabel(10,"NLM");
 272         fHistClusterIdentificationCuts->GetXaxis()->SetBinLabel(11,"out");
 273         fHistograms->Add(fHistClusterIdentificationCuts);
 274
 275         // Acceptance related histogramms
 276         fHistClusterEtavsPhiBeforeAcc=new TH2F(Form("EtaPhi_beforeAcceptance %s",GetCutNumber().Data()),"EtaPhi_beforeAcceptance",462,-TMath::Pi(),TMath::Pi(),110,-0.7,0.7);
 277         fHistograms->Add(fHistClusterEtavsPhiBeforeAcc);
 278         fHistClusterEtavsPhiAfterAcc=new TH2F(Form("EtaPhi_afterAcceptance %s",GetCutNumber().Data()),"EtaPhi_afterAcceptance",462,-TMath::Pi(),TMath::Pi(),110,-0.7,0.7);
 279         fHistograms->Add(fHistClusterEtavsPhiAfterAcc);
 280         fHistClusterEtavsPhiAfterQA=new TH2F(Form("EtaPhi_afterClusterQA %s",GetCutNumber().Data()),"EtaPhi_afterClusterQA",462,-TMath::Pi(),TMath::Pi(),110,-0.7,0.7);
 281         fHistograms->Add(fHistClusterEtavsPhiAfterQA);
 282         fHistDistanceToBadChannelBeforeAcc = new TH1F(Form("DistanceToBadChannel_beforeAcceptance %s",GetCutNumber().Data()),"DistanceToBadChannel_beforeAcceptance",200,0,40);
 283         fHistograms->Add(fHistDistanceToBadChannelBeforeAcc);
 284         fHistDistanceToBadChannelAfterAcc = new TH1F(Form("DistanceToBadChannel_afterAcceptance %s",GetCutNumber().Data()),"DistanceToBadChannel_afterAcceptance",200,0,40);
 285         fHistograms->Add(fHistDistanceToBadChannelAfterAcc);
 286
 287         // Cluster quality related histograms
 288         fHistClusterTimevsEBeforeQA=new TH2F(Form("ClusterTimeVsE_beforeClusterQA %s",GetCutNumber().Data()),"ClusterTimeVsE_beforeClusterQA",400,-10e-6,10e-6,100,0.,40);
 289         fHistograms->Add(fHistClusterTimevsEBeforeQA);
 290         fHistClusterTimevsEAfterQA=new TH2F(Form("ClusterTimeVsE_afterClusterQA %s",GetCutNumber().Data()),"ClusterTimeVsE_afterClusterQA",400,-10e-6,10e-6,100,0.,40);
 291         fHistograms->Add(fHistClusterTimevsEAfterQA);
 292         fHistExoticCellBeforeQA=new TH2F(Form("ExoticCell_beforeClusterQA %s",GetCutNumber().Data()),"ExoticCell_beforeClusterQA",400,0,40,50,0.75,1);
 293         fHistograms->Add(fHistExoticCellBeforeQA);
 294         fHistExoticCellAfterQA=new TH2F(Form("ExoticCell_afterClusterQA %s",GetCutNumber().Data()),"ExoticCell_afterClusterQA",400,0,40,50,0.75,1);
 295         fHistograms->Add(fHistExoticCellAfterQA);
 296         fHistNMatchedTracks = new TH1F(Form("NMatchedTracks_%s",GetCutNumber().Data()),"NMatchedTracks",22,-1.5,20.5);
 297         fHistograms->Add(fHistNMatchedTracks);
 298         fHistDistanceTrackToClusterBeforeQA = new TH1F(Form("DistanceToTrack_beforeClusterQA %s",GetCutNumber().Data()),"DistanceToTrack_beforeClusterQA",200,0,2);
 299         fHistograms->Add(fHistDistanceTrackToClusterBeforeQA);
 300         fHistDistanceTrackToClusterAfterQA = new TH1F(Form("DistanceToTrack_afterClusterQA %s",GetCutNumber().Data()),"DistanceToTrack_afterClusterQA",200,0,2);
 301         fHistograms->Add(fHistDistanceTrackToClusterAfterQA);
 302         fHistEnergyOfClusterBeforeQA = new TH1F(Form("EnergyOfCluster_beforeClusterQA %s",GetCutNumber().Data()),"EnergyOfCluster_beforeClusterQA",300,0,30);
 303         fHistograms->Add(fHistEnergyOfClusterBeforeQA);
 304         fHistEnergyOfClusterAfterQA = new TH1F(Form("EnergyOfCluster_afterClusterQA %s",GetCutNumber().Data()),"EnergyOfCluster_afterClusterQA",300,0,30);
 305         fHistograms->Add(fHistEnergyOfClusterAfterQA);
 306         fHistNCellsBeforeQA = new TH1F(Form("NCellPerCluster_beforeClusterQA %s",GetCutNumber().Data()),"NCellPerCluster_beforeClusterQA",50,0,50);
 307         fHistograms->Add(fHistNCellsBeforeQA);
 308         fHistNCellsAfterQA = new TH1F(Form("NCellPerCluster_afterClusterQA %s",GetCutNumber().Data()),"NCellPerCluster_afterClusterQA",50,0,50);
 309         fHistograms->Add(fHistNCellsAfterQA);
 310         fHistM02BeforeQA = new TH1F(Form("M02_beforeClusterQA %s",GetCutNumber().Data()),"M02_beforeClusterQA",100,0,5);
 311         fHistograms->Add(fHistM02BeforeQA);
 312         fHistM02AfterQA = new TH1F(Form("M02_afterClusterQA %s",GetCutNumber().Data()),"M02_afterClusterQA",100,0,5);
 313         fHistograms->Add(fHistM02AfterQA);
 314         fHistM20BeforeQA = new TH1F(Form("M20_beforeClusterQA %s",GetCutNumber().Data()),"M20_beforeClusterQA",100,0,2.5);
 315         fHistograms->Add(fHistM20BeforeQA);
 316         fHistM20AfterQA = new TH1F(Form("M20_afterClusterQA %s",GetCutNumber().Data()),"M20_afterClusterQA",100,0,2.5);
 317         fHistograms->Add(fHistM20AfterQA);
 318         fHistDispersionBeforeQA = new TH1F(Form("Dispersion_beforeClusterQA %s",GetCutNumber().Data()),"Dispersion_beforeClusterQA",100,0,4);
 319         fHistograms->Add(fHistDispersionBeforeQA);
 320         fHistDispersionAfterQA = new TH1F(Form("Dispersion_afterClusterQA %s",GetCutNumber().Data()),"Dispersion_afterClusterQA",100,0,4);
 321         fHistograms->Add(fHistDispersionAfterQA);
 322         fHistNLMBeforeQA = new TH1F(Form("NLM_beforeClusterQA %s",GetCutNumber().Data()),"NLM_beforeClusterQA",10,0,10);
 323         fHistograms->Add(fHistNLMBeforeQA);
 324         fHistNLMAfterQA = new TH1F(Form("NLM_afterClusterQA %s",GetCutNumber().Data()),"NLM_afterClusterQA",10,0,10);
 325         fHistograms->Add(fHistNLMAfterQA);
 326
 327         TH1::AddDirectory(kTRUE);
 328 }
 329
 330
 331 ///________________________________________________________________________
 332 Bool_t AliCaloPhotonCuts::ClusterIsSelectedMC(TParticle *particle,AliStack *fMCStack){
 333    // MonteCarlo Photon Selection
 334
 335         if(!fMCStack)return kFALSE;
 336
 337         if (particle->GetPdgCode() == 22){
 338
 339                 if ( particle->Eta() < fMinEtaCut || particle->Eta() > fMaxEtaCut ) return kFALSE;
 340                 if ( particle->Phi() < fMinPhiCut || particle->Phi() > fMaxPhiCut ) return kFALSE;
 341
 342                 if(particle->GetMother(0) >-1 && fMCStack->Particle(particle->GetMother(0))->GetPdgCode() == 22){
 343                         return kFALSE; // no photon as mothers!
 344                 }
 345                 if(particle->GetMother(0) >= fMCStack->GetNprimary()){
 346                         return kFALSE; // the gamma has a mother, and it is not a primary particle
 347                 }
 348                 return kTRUE;
 349         }
 350         return kFALSE;
 351 }
 352 ///________________________________________________________________________
 353 Bool_t AliCaloPhotonCuts::ClusterIsSelectedAODMC(AliAODMCParticle *particle,TClonesArray *aodmcArray){
 354         // MonteCarlo Photon Selection
 355
 356         if(!aodmcArray)return kFALSE;
 357         if (particle->GetPdgCode() == 22){
 358                 if ( particle->Eta() < fMinEtaCut || particle->Eta() > fMaxEtaCut ) return kFALSE;
 359                 if ( particle->Phi() < fMinPhiCut || particle->Phi() > fMaxPhiCut ) return kFALSE;
 360                 if(particle->GetMother() > -1){
 361                         if((static_cast<AliAODMCParticle*>(aodmcArray->At(particle->GetMother())))->GetPdgCode() == 22){
 362                                 return kFALSE; // no photon as mothers!
 363                         }
 364                         if(!(static_cast<AliAODMCParticle*>(aodmcArray->At(particle->GetMother()))->IsPrimary())){
 365                                 return kFALSE; // the gamma has a mother, and it is not a primary particle
 366                         }
 367                 }
 368                 return kTRUE; // return in case of accepted gamma
 369         }
 370         return kFALSE;
 371 }
 372
 373
 374
 375 ///________________________________________________________________________
 376 // This function selects the clusters based on their quality criteria
 377 ///________________________________________________________________________
 378 Bool_t AliCaloPhotonCuts::ClusterQualityCuts(AliVCluster* cluster, AliVEvent *event, Bool_t isMC)
 379 {   // Specific Photon Cuts
 380
 381         Int_t cutIndex = 0;
 382         if(fHistClusterIdentificationCuts)fHistClusterIdentificationCuts->Fill(cutIndex);
 383         cutIndex++;
 384
 385 //      Double_t minR = 999.0;
 386         // get the minimum radius of tracks to cluster
 387 //      if(fHistDistanceTrackToClusterBeforeQA || fHistDistanceTrackToClusterAfterQA){
 388 //              Float_t pos[3];
 389 //              cluster->GetPosition(pos);  // Get cluster position
 390 //
 391 //              TVector3 cp(pos);
 392 //              int NtrMatched = 0;
 393 //              NtrMatched = cluster->GetNTracksMatched();
 394 //              fHistNMatchedTracks->Fill(NtrMatched);
 395 //              //loop over tracks for Jet QA
 396 //              TList *l = event->GetList();
 397 //              TClonesArray *tracks = dynamic_cast<TClonesArray*>(l->FindObject("Tracks"));
 398 //              for(int itrack = 0; itrack < NtrMatched; itrack++){
 399 //                      AliVTrack *trackcluster = static_cast<AliVTrack*>(tracks->At(itrack));
 400 //                      if (! trackcluster) {
 401 //                              AliError(Form("Couldn't get ESD track %d\n", itrack));
 402 //                              continue;
 403 //                      }
 404 //                      Double_t dphi = -999.0;
 405 //                      Double_t deta = -999.0;
 406 //                      AliPicoTrack::GetEtaPhiDiff(trackcluster, cluster, dphi, deta);
 407 //                      cout << "here" << endl;
 408 //                      Double_t dr = sqrt(dphi*dphi + deta+deta);
 409 //                      if(dr < minR)
 410 //                              minR = dr;
 411 //              }//loop over tracks
 412 //      }
 413
 414         // Fill Histos before Cuts
 415         if(fHistClusterTimevsEBeforeQA) fHistClusterTimevsEBeforeQA->Fill(cluster->GetTOF(), cluster->E());
 416 //      if(fHistExoticCellBeforeQA) fHistExoticCellBeforeQA->Fill(cluster->E(), );
 417 //      if(fHistDistanceTrackToClusterBeforeQA) fHistDistanceTrackToClusterBeforeQA->Fill(minR);
 418         if(fHistEnergyOfClusterBeforeQA) fHistEnergyOfClusterBeforeQA->Fill(cluster->E());
 419         if(fHistNCellsBeforeQA) fHistNCellsBeforeQA->Fill(cluster->GetNCells());
 420         if(fHistM02BeforeQA) fHistM02BeforeQA->Fill(cluster->GetM02());
 421         if(fHistM20BeforeQA) fHistM20BeforeQA->Fill(cluster->GetM20());
 422         if(fHistDispersionBeforeQA) fHistDispersionBeforeQA->Fill(cluster->GetDispersion());
 423 //      if(fHistNLMBeforeQA) fHistNLMBeforeQA->Fill(cluster->GetNExMax());
 424
 425         // Check wether timing is ok
 426         if (fUseTimeDiff){
 427                 if(abs(cluster->GetTOF()) > fMaxTimeDiff && !isMC){
 428                         if(fHistClusterIdentificationCuts)fHistClusterIdentificationCuts->Fill(cutIndex); //1
 429                         return kFALSE;
 430                 }
 431         }
 432         cutIndex++; //2, next cut
 433
 434         // Minimum distance to track
 435 //      if (fUseDistTrackToCluster){
 436 //              Float_t pos[3];
 437 //              cluster->GetPosition(pos);  // Get cluster position
 438 //              TVector3 cp(pos);
 439 //              int NtrMatched = 0;
 440 //              NtrMatched = cluster->GetNTracksMatched();
 441 //              fHistNMatchedTracks->Fill(NtrMatched);
 442 //              if(NtrMatched > 0){
 443 //                      //loop over tracks for QA
 444 //                      TList *l = event->GetList();
 445 //                      TClonesArray *tracks = dynamic_cast<TClonesArray*>(l->FindObject("Tracks"));
 446 //
 447 //                      Double_t dphi = 999.0;
 448 //                      Double_t deta = 999.0;
 449 //                      Double_t dr2 = 999.0;
 450 //
 451 //                      for(int itrack = 0; itrack < NtrMatched; itrack++){
 452 //                              AliVTrack *trackcluster = NULL;
 453 //                              trackcluster = static_cast<AliVTrack*>(tracks->At(itrack));
 454 //                              if (! trackcluster) {
 455 //                              AliError(Form("Couldn't get ESD track %d\n", itrack));
 456 //                              continue;
 457 //                              }
 458 //                              AliPicoTrack::GetEtaPhiDiff(trackcluster, cluster, dphi, deta);
 459 //                              dr2 = dphi*dphi + deta+deta;
 460 //                              //cout << dr << endl;
 461 //                              if(dr2 < fMinDistTrackToCluster*fMinDistTrackToCluster){
 462 //              //        if(dphi < fMinDistTrackToCluster || deta < fMinDistTrackToCluster){
 463 //                              if(fHistClusterIdentificationCuts)fHistClusterIdentificationCuts->Fill(cutIndex); //2
 464 //                              return kFALSE;
 465 //                              }
 466 //
 467 //                      }//loop over tracks
 468 //              }
 469 // //           if(cluster->GetEmcCpvDistance() < fMinDistTrackToCluster){
 470 // //                   if(fHistClusterIdentificationCuts)fHistClusterIdentificationCuts->Fill(cutIndex); //2
 471 // //                   return kFALSE;
 472 // //           }
 473 //      }
 474         cutIndex++;//3, next cut
 475
 476         // exotic cell cut --IMPLEMENT LATER---
 477 //      if(!AcceptanceCuts(photon)){
 478 //              if(fHistClusterIdentificationCuts)fHistClusterIdentificationCuts->Fill(cutIndex); //3
 479 //              return kFALSE;
 480 //      }
 481         cutIndex++; //4, next cut
 482
 483         // minimum cell energy cut
 484         if (fUseMinEnergy){
 485                 if(cluster->E() < fMinEnergy){
 486                         if(fHistClusterIdentificationCuts)fHistClusterIdentificationCuts->Fill(cutIndex); //4
 487                         return kFALSE;
 488                 }
 489         }
 490         cutIndex++; //5, next cut
 491
 492         // minimum number of cells
 493         if (fUseNCells){
 494                 if(cluster->GetNCells() < fMinNCells) {
 495                         if(fHistClusterIdentificationCuts)fHistClusterIdentificationCuts->Fill(cutIndex); //5
 496                         return kFALSE;
 497                 }
 498         }
 499         cutIndex++; //6, next cut
 500
 501         // M02 cut
 502         if (fUseM02){
 503                 if( cluster->GetM02()< fMinM02 || cluster->GetM02() > fMaxM02 ) {
 504                         if(fHistClusterIdentificationCuts)fHistClusterIdentificationCuts->Fill(cutIndex); //6
 505                         return kFALSE;
 506                 }
 507         }
 508         cutIndex++; //7, next cut
 509
 510         // M20 cut
 511         if (fUseM20){
 512                 if( cluster->GetM20()< fMinM20 || cluster->GetM20() > fMaxM20 ) {
 513                         if(fHistClusterIdentificationCuts)fHistClusterIdentificationCuts->Fill(cutIndex); //7
 514                         return kFALSE;
 515                 }
 516         }
 517         cutIndex++; //8, next cut
 518
 519         // dispersion cut
 520         if (fUseDispersion){
 521                 if( cluster->GetDispersion()> fMaxDispersion) {
 522                         if(fHistClusterIdentificationCuts)fHistClusterIdentificationCuts->Fill(cutIndex); //8
 523                         return kFALSE;
 524                 }
 525         }
 526         cutIndex++; //9, next cut
 527
 528         // NLM cut --IMPLEMENT LATER---
 529 //      if (fUseNLM){
 530 //              if( cluster->GetDispersion()> fMaxDispersion) {
 531 //                      if(fHistClusterIdentificationCuts)fHistClusterIdentificationCuts->Fill(cutIndex); //9
 532 //                      return kFALSE;
 533 //              }
 534 //      }
 535         cutIndex++; //9, next cut
 536
 537         // DONE with selecting photons
 538         if(fHistClusterIdentificationCuts)fHistClusterIdentificationCuts->Fill(cutIndex); //10
 539
 540         // Histos after Cuts
 541         Double_t vertex[3] = {0};
 542         event->GetPrimaryVertex()->GetXYZ(vertex);
 543         // TLorentzvector with cluster
 544         TLorentzVector clusterVector;
 545         cluster->GetMomentum(clusterVector,vertex);
 546         Double_t etaCluster = clusterVector.Eta();
 547         Double_t phiCluster = clusterVector.Phi();
 548
 549         if(fHistClusterEtavsPhiAfterQA) fHistClusterEtavsPhiAfterQA->Fill(phiCluster,etaCluster);
 550         if(fHistClusterTimevsEAfterQA) fHistClusterTimevsEAfterQA->Fill(cluster->GetTOF(), cluster->E());
 551 //      if(fHistExoticCellAfterQA) fHistExoticCellAfterQA->Fill(cluster->E(), );
 552 //      if(fHistDistanceTrackToClusterAfterQA) fHistDistanceTrackToClusterAfterQA->Fill(minR);
 553         if(fHistEnergyOfClusterAfterQA) fHistEnergyOfClusterAfterQA->Fill(cluster->E());
 554         if(fHistNCellsAfterQA) fHistNCellsAfterQA->Fill(cluster->GetNCells());
 555         if(fHistM02AfterQA) fHistM02AfterQA->Fill(cluster->GetM02());
 556         if(fHistM20AfterQA) fHistM20AfterQA->Fill(cluster->GetM20());
 557         if(fHistDispersionAfterQA) fHistDispersionAfterQA->Fill(cluster->GetDispersion());
 558 //      if(fHistNLMBeforeQA) fHistNLMAfterQA->Fill(cluster->GetNExMax());
 559
 560         return kTRUE;
 561
 562 }
 563
 564
 565 ///________________________________________________________________________
 566 Bool_t AliCaloPhotonCuts::ClusterIsSelected(AliVCluster *cluster, AliVEvent * event, Bool_t isMC)
 567 {
 568         //Selection of Reconstructed photon clusters with Calorimeters
 569
 570         FillClusterCutIndex(kPhotonIn);
 571
 572         Double_t vertex[3] = {0};
 573         event->GetPrimaryVertex()->GetXYZ(vertex);
 574         // TLorentzvector with cluster
 575         TLorentzVector clusterVector;
 576         cluster->GetMomentum(clusterVector,vertex);
 577         Double_t etaCluster = clusterVector.Eta();
 578         Double_t phiCluster = clusterVector.Phi();
 579
 580         // Histos before cuts
 581         if(fHistClusterEtavsPhiBeforeAcc) fHistClusterEtavsPhiBeforeAcc->Fill(phiCluster,etaCluster);
 582
 583         // Cluster Selection - 0= accept any calo cluster
 584         if (fClusterType > 0){
 585                 //Select EMCAL cluster
 586                 if (fClusterType == 1 && !cluster->IsEMCAL()){
 587                         FillClusterCutIndex(kDetector);
 588                         return kFALSE;
 589                 }
 590                 //Select PHOS cluster
 591                 if (fClusterType == 2 && !cluster->IsPHOS()){
 592                         FillClusterCutIndex(kDetector);
 593                         return kFALSE;
 594                 }
 595         }
 596
 597         // Acceptance Cuts
 598         if(!AcceptanceCuts(cluster,event)){
 599                 FillClusterCutIndex(kAcceptance);
 600                 return kFALSE;
 601         }
 602         // Cluster Quality Cuts
 603         if(!ClusterQualityCuts(cluster,event,isMC)){
 604                 FillClusterCutIndex(kClusterQuality);
 605                 return kFALSE;
 606         }
 607
 608         // Photon passed cuts
 609         FillClusterCutIndex(kPhotonOut);
 610         return kTRUE;
 611 }
 612
 613
 614 ///________________________________________________________________________
 615 Bool_t AliCaloPhotonCuts::AcceptanceCuts(AliVCluster *cluster, AliVEvent* event)
 616 {
 617    // Exclude certain areas for photon reconstruction
 618
 619         Int_t cutIndex=0;
 620         if(fHistAcceptanceCuts)fHistAcceptanceCuts->Fill(cutIndex);
 621         cutIndex++;
 622
 623
 624         Double_t vertex[3] = {0};
 625         event->GetPrimaryVertex()->GetXYZ(vertex);
 626         // TLorentzvector with cluster
 627         TLorentzVector clusterVector;
 628         cluster->GetMomentum(clusterVector,vertex);
 629         Double_t etaCluster = clusterVector.Eta();
 630         Double_t phiCluster = clusterVector.Phi();
 631
 632         // check eta range
 633         if (fUseEtaCut){
 634                 if (etaCluster < fMinEtaCut || etaCluster > fMaxEtaCut){
 635                         if(fHistAcceptanceCuts)fHistAcceptanceCuts->Fill(cutIndex);
 636                         return kFALSE;
 637                 }
 638         }
 639         cutIndex++;
 640
 641         // check phi range
 642         if (fUsePhiCut ){
 643                 if (phiCluster < fMinPhiCut || phiCluster > fMaxEtaCut){
 644                         if(fHistAcceptanceCuts)fHistAcceptanceCuts->Fill(cutIndex);
 645                         return kFALSE;
 646                 }
 647         }
 648         cutIndex++;
 649
 650         // check distance to bad channel
 651         if (fUseDistanceToBadChannel){
 652                 if (cluster->GetDistanceToBadChannel() < fMinDistanceToBadChannel){
 653                         if(fHistAcceptanceCuts)fHistAcceptanceCuts->Fill(cutIndex);
 654                         return kFALSE;
 655                 }
 656         }
 657         cutIndex++;
 658         if(fHistAcceptanceCuts)fHistAcceptanceCuts->Fill(cutIndex);
 659
 660         // Histos after cuts
 661         if(fHistClusterEtavsPhiAfterAcc) fHistClusterEtavsPhiAfterAcc->Fill(phiCluster,etaCluster);
 662
 663         return kTRUE;
 664 }
 665
 666 Bool_t AliCaloPhotonCuts::MatchConvPhotonToCluster(AliAODConversionPhoton* convPhoton, AliVCluster* cluster, AliVEvent* event ){
 667
 668         if (!fUseDistTrackToCluster) return kFALSE;
 669
 670         AliESDEvent *esdev = dynamic_cast<AliESDEvent*>(event);
 671         AliAODEvent *aodev = 0;
 672         if (!esdev) {
 673                 aodev = dynamic_cast<AliAODEvent*>(event);
 674                 if (!aodev) {
 675                         AliError("Task needs AOD or ESD event, returning");
 676                         return kFALSE;
 677                 }
 678         }
 679
 680 //      cout << "Got the event" << endl;
 681
 682         Double_t vertex[3] = {0};
 683         event->GetPrimaryVertex()->GetXYZ(vertex);
 684         // TLorentzvector with cluster
 685         TLorentzVector clusterVector;
 686         cluster->GetMomentum(clusterVector,vertex);
 687         Double_t etaCluster = clusterVector.Eta();
 688         Double_t phiCluster = clusterVector.Phi();
 689
 690         Bool_t matched = kFALSE;
 691         for (Int_t i = 0; i < 2; i++){
 692                 Int_t tracklabel = convPhoton->GetLabel(i);
 693                 AliVTrack *inTrack = 0x0;
 694                 if(esdev) {
 695                         if(tracklabel > event->GetNumberOfTracks() ) continue;
 696                         inTrack = esdev->GetTrack(tracklabel);
 697                 } else {
 698                         if(((AliV0ReaderV1*)AliAnalysisManager::GetAnalysisManager()->GetTask("V0ReaderV1"))->AreAODsRelabeled()){
 699                                 inTrack = dynamic_cast<AliVTrack*>(event->GetTrack(tracklabel));
 700                         } else {
 701                                 for(Int_t ii=0; ii<event->GetNumberOfTracks(); ii++) {
 702                                         inTrack = dynamic_cast<AliVTrack*>(event->GetTrack(ii));
 703                                         if(inTrack){
 704                                                 if(inTrack->GetID() == tracklabel) {
 705                                                         continue;
 706                                                 }
 707                                         }
 708                                 }
 709                         }
 710                 }
 711 //              cout << "found track " << endl;
 712 //              AliVTrack *outTrack = 0x00;
 713 //              if (esdev)
 714 //                      outTrack = new AliESDtrack(*((AliESDtrack*)inTrack));
 715 //              else
 716 //                      outTrack = new AliAODTrack(*((AliAODTrack*)inTrack));
 717
 718
 719
 720                 Bool_t propagated = AliEMCALRecoUtils::ExtrapolateTrackToEMCalSurface(inTrack,440, 0.000510999,20,0.05);
 721                 if (propagated){
 722 //                      cout << "Track "<< i  << "\t"<< inTrack->GetTrackPhiOnEMCal() << "\t" << inTrack->GetTrackEtaOnEMCal() << endl;
 723 //                      cout << "Cluster " << phiCluster << "\t" << etaCluster << endl;
 724                         Double_t dPhi = abs(phiCluster-inTrack->GetTrackPhiOnEMCal());
 725                         Double_t dEta = abs(etaCluster-inTrack->GetTrackEtaOnEMCal());
 726
 727                         Double_t dR2 = dPhi*dPhi + dEta+dEta;
 728 //                      cout << "distance to cluster \t" << sqrt(dR2) << endl;
 729                         if (fHistDistanceTrackToClusterBeforeQA)fHistDistanceTrackToClusterBeforeQA->Fill(sqrt(dR2));
 730                         if(dR2 < fMinDistTrackToCluster*fMinDistTrackToCluster){
 731                                 matched = kTRUE;
 732                                 if (fHistDistanceTrackToClusterAfterQA)fHistDistanceTrackToClusterAfterQA->Fill(sqrt(dR2));
 733                         }
 734                 }
 735         }
 736         return matched;
 737 }
 738
 739 //____________________________________________________________________________________________
 740
 741
 742 ///________________________________________________________________________
 743 Bool_t AliCaloPhotonCuts::UpdateCutString() {
 744    ///Update the cut string (if it has been created yet)
 745
 746    if(fCutString && fCutString->GetString().Length() == kNCuts) {
 747       fCutString->SetString(GetCutNumber());
 748    } else {
 749       return kFALSE;
 750    }
 751    return kTRUE;
 752 }
 753
 754 ///________________________________________________________________________
 755 Bool_t AliCaloPhotonCuts::InitializeCutsFromCutString(const TString analysisCutSelection ) {
 756         // Initialize Cuts from a given Cut string
 757         AliInfo(Form("Set CaloCut Number: %s",analysisCutSelection.Data()));
 758         if(analysisCutSelection.Length()!=kNCuts) {
 759                 AliError(Form("Cut selection has the wrong length! size is %d, number of cuts is %d", analysisCutSelection.Length(), kNCuts));
 760                 return kFALSE;
 761         }
 762         if(!analysisCutSelection.IsDigit()){
 763                 AliError("Cut selection contains characters");
 764                 return kFALSE;
 765         }
 766
 767         const char *cutSelection = analysisCutSelection.Data();
 768         #define ASSIGNARRAY(i)  fCuts[i] = cutSelection[i] - '0'
 769         for(Int_t ii=0;ii<kNCuts;ii++){
 770                 ASSIGNARRAY(ii);
 771         }
 772
 773         // Set Individual Cuts
 774         for(Int_t ii=0;ii<kNCuts;ii++){
 775                 if(!SetCut(cutIds(ii),fCuts[ii]))return kFALSE;
 776         }
 777         PrintCutsWithValues();
 778         return kTRUE;
 779 }
 780
 781 ///________________________________________________________________________
 782 Bool_t AliCaloPhotonCuts::SetCut(cutIds cutID, const Int_t value) {
 783         ///Set individual cut ID
 784
 785         switch (cutID) {
 786
 787                 case kClusterType:
 788                         if( SetClusterTypeCut(value)) {
 789                                 fCuts[kClusterType] = value;
 790                                 UpdateCutString();
 791                                 return kTRUE;
 792                         } else return kFALSE;
 793
 794                 case kEtaMin:
 795                         if( SetMinEtaCut(value)) {
 796                                 fCuts[kEtaMin] = value;
 797                                 UpdateCutString();
 798                                 return kTRUE;
 799                         } else return kFALSE;
 800
 801                 case kEtaMax:
 802                         if( SetMaxEtaCut(value)) {
 803                                 fCuts[kEtaMax] = value;
 804                                 UpdateCutString();
 805                                 return kTRUE;
 806                         } else return kFALSE;
 807
 808                 case kPhiMin:
 809                         if( SetMinPhiCut(value)) {
 810                                 fCuts[kPhiMin] = value;
 811                                 UpdateCutString();
 812                                 return kTRUE;
 813                         } else return kFALSE;
 814
 815                 case kPhiMax:
 816                         if( SetMaxPhiCut(value)) {
 817                                 fCuts[kPhiMax] = value;
 818                                 UpdateCutString();
 819                                 return kTRUE;
 820                         } else return kFALSE;
 821
 822                 case kDistanceToBadChannel:
 823                         if( SetDistanceToBadChannelCut(value)) {
 824                                 fCuts[kDistanceToBadChannel] = value;
 825                                 UpdateCutString();
 826                                 return kTRUE;
 827                         } else return kFALSE;
 828
 829                 case kTiming:
 830                         if( SetTimingCut(value)) {
 831                                 fCuts[kTiming] = value;
 832                                 UpdateCutString();
 833                                 return kTRUE;
 834                         } else return kFALSE;
 835
 836                 case kTrackMatching:
 837                         if( SetTrackMatchingCut(value)) {
 838                                 fCuts[kTrackMatching] = value;
 839                                 UpdateCutString();
 840                                 return kTRUE;
 841                         } else return kFALSE;
 842
 843                 case kExoticCell:
 844                         if( SetExoticCellCut(value)) {
 845                                 fCuts[kExoticCell] = value;
 846                                 UpdateCutString();
 847                                 return kTRUE;
 848                         } else return kFALSE;
 849
 850                 case kMinEnery:
 851                         if( SetMinEnergyCut(value)) {
 852                                 fCuts[kMinEnery] = value;
 853                                 UpdateCutString();
 854                                 return kTRUE;
 855                         } else return kFALSE;
 856
 857                 case kNMinCells:
 858                         if( SetMinNCellsCut(value)) {
 859                                 fCuts[kNMinCells] = value;
 860                                 UpdateCutString();
 861                                 return kTRUE;
 862                         } else return kFALSE;
 863
 864                 case kMinM02:
 865                         if( SetMinM02(value)) {
 866                                 fCuts[kMinM02] = value;
 867                                 UpdateCutString();
 868                                 return kTRUE;
 869                         } else return kFALSE;
 870
 871                 case kMaxM02:
 872                         if( SetMaxM02(value)) {
 873                                 fCuts[kMaxM02] = value;
 874                                 UpdateCutString();
 875                                 return kTRUE;
 876                         } else return kFALSE;
 877
 878                 case kMinM20:
 879                         if( SetMinM20(value)) {
 880                                 fCuts[kMinM20] = value;
 881                                 UpdateCutString();
 882                                 return kTRUE;
 883                         } else return kFALSE;
 884
 885                 case kMaxM20:
 886                         if( SetMaxM20(value)) {
 887                                 fCuts[kMaxM20] = value;
 888                                 UpdateCutString();
 889                                 return kTRUE;
 890                         } else return kFALSE;
 891
 892                 case kDispersion:
 893                         if( SetDispersion(value)) {
 894                                 fCuts[kDispersion] = value;
 895                                 UpdateCutString();
 896                                 return kTRUE;
 897                         } else return kFALSE;
 898
 899                 case kNLM:
 900                         if( SetNLM(value)) {
 901                                 fCuts[kNLM] = value;
 902                                 UpdateCutString();
 903                                 return kTRUE;
 904                         } else return kFALSE;
 905
 906                 case kNCuts:
 907                         AliError("Cut id out of range");
 908                         return kFALSE;
 909         }
 910
 911         AliError("Cut id %d not recognized");
 912         return kFALSE;
 913
 914
 915 }
 916 ///________________________________________________________________________
 917 void AliCaloPhotonCuts::PrintCuts() {
 918    // Print out current Cut Selection
 919    for(Int_t ic = 0; ic < kNCuts; ic++) {
 920       printf("%-30s : %d \n", fgkCutNames[ic], fCuts[ic]);
 921    }
 922 }
 923
 924 void AliCaloPhotonCuts::PrintCutsWithValues() {
 925         // Print out current Cut Selection with value
 926         printf("\nCluster cutnumber \n");
 927         for(Int_t ic = 0; ic < kNCuts; ic++) {
 928                 printf("%d",fCuts[ic]);
 929         }
 930         printf("\n\n");
 931
 932         printf("Acceptance cuts: \n");
 933         if (fClusterType == 0) printf("\tall calorimeter clusters are used\n");
 934         if (fClusterType == 1) printf("\tEMCAL calorimeter clusters are used\n");
 935         if (fClusterType == 2) printf("\tPHOS calorimeter clusters are used\n");
 936         if (fUseEtaCut) printf("\t%3.2f < eta_{cluster} < %3.2f\n", fMinEtaCut, fMaxEtaCut );
 937         if (fUsePhiCut) printf("\t%3.2f < phi_{cluster} < %3.2f\n", fMinPhiCut, fMaxPhiCut );
 938         if (fUseDistanceToBadChannel) printf("\tcut on exotics applied \n");
 939
 940         printf("Cluster Quality cuts: \n");
 941         if (fUseTimeDiff) printf("\t time difference < %3.2f\n", fMaxTimeDiff );
 942         if (fUseDistTrackToCluster) printf("\tmin distance to track > %3.2f\n", fMinDistTrackToCluster );
 943         if (fUseExoticCell)printf("\t min distance to track > %3.2f\n", fMinDistTrackToCluster );
 944     if (fUseMinEnergy)printf("\t E_{cluster} > %3.2f\n", fMinEnergy );
 945         if (fUseNCells) printf("\t number of cells per cluster > %d\n", fMinNCells );
 946         if (fUseM02) printf("\t %3.2f < M02 < %3.2f\n", fMinM02, fMaxM02 );
 947         if (fUseM20) printf("\t %3.2f < M20 < %3.2f\n", fMinM20, fMaxM20 );
 948         if (fUseDispersion) printf("\t dispersion < %3.2f\n", fMaxDispersion );
 949         if (fUseNLM) printf("\t %d < NLM < %d\n", fMinNLM, fMaxNLM );
 950
 951 }
 952
 953 // EMCAL acceptance 2011
 954 // 1.39626, 3.125 (phi)
 955 // -0.66687,,0.66465
 956
 957
 958 ///________________________________________________________________________
 959 Bool_t AliCaloPhotonCuts::SetClusterTypeCut(Int_t clusterType)
 960 {   // Set Cut
 961         switch(clusterType){
 962         case 0: // all clusters
 963                 fClusterType=0;
 964                 break;
 965         case 1: // EMCAL clusters
 966                 fClusterType=1;
 967                 break;
 968         case 2: // PHOS clusters
 969                 fClusterType=2;
 970                 break;
 971         default:
 972                 AliError(Form("ClusterTypeCut not defined %d",clusterType));
 973                 return kFALSE;
 974         }
 975         return kTRUE;
 976 }
 977
 978 //___________________________________________________________________
 979 Bool_t AliCaloPhotonCuts::SetMinEtaCut(Int_t minEta)
 980 {
 981         switch(minEta){
 982         case 0:
 983                 if (!fUseEtaCut) fUseEtaCut=0;
 984                 fMinEtaCut=-10.;
 985                 break;
 986         case 1:
 987                 if (!fUseEtaCut) fUseEtaCut=1;
 988                 fMinEtaCut=-0.6687;
 989                 break;
 990         case 2:
 991                 if (!fUseEtaCut) fUseEtaCut=1;
 992                 fMinEtaCut=-0.5;
 993                 break;
 994         case 3:
 995                 if (!fUseEtaCut) fUseEtaCut=1;
 996                 fMinEtaCut=-2;
 997                 break;
 998         default:
 999                 AliError(Form("MinEta Cut not defined %d",minEta));
1000                 return kFALSE;
1001         }
1002         return kTRUE;
1003 }
1004
1005
1006 //___________________________________________________________________
1007 Bool_t AliCaloPhotonCuts::SetMaxEtaCut(Int_t maxEta)
1008 {
1009         switch(maxEta){
1010         case 0:
1011                 if (!fUseEtaCut) fUseEtaCut=0;
1012                 fMaxEtaCut=10;
1013                 break;
1014         case 1:
1015                 if (!fUseEtaCut) fUseEtaCut=1;
1016                 fMaxEtaCut=0.66465;
1017                 break;
1018         case 2:
1019                 if (!fUseEtaCut) fUseEtaCut=1;
1020                 fMaxEtaCut=0.5;
1021                 break;
1022         case 3:
1023                 if (!fUseEtaCut) fUseEtaCut=1;
1024                 fMaxEtaCut=2;
1025                 break;
1026         default:
1027                 AliError(Form("MaxEta Cut not defined %d",maxEta));
1028                 return kFALSE;
1029         }
1030         return kTRUE;
1031 }
1032
1033 //___________________________________________________________________
1034 Bool_t AliCaloPhotonCuts::SetMinPhiCut(Int_t minPhi)
1035 {
1036         switch(minPhi){
1037         case 0:
1038                 if (!fUsePhiCut) fUsePhiCut=0;
1039                 fMinPhiCut=-10000;
1040                 break;
1041         case 1:
1042                 if (!fUsePhiCut) fUsePhiCut=1;
1043                 fMinPhiCut=1.39626;
1044                 break;
1045         default:
1046                 AliError(Form("MinPhi Cut not defined %d",minPhi));
1047                 return kFALSE;
1048         }
1049         return kTRUE;
1050 }
1051
1052 //___________________________________________________________________
1053 Bool_t AliCaloPhotonCuts::SetMaxPhiCut(Int_t maxPhi)
1054 {
1055         switch(maxPhi){
1056         case 0:
1057                 if (!fUsePhiCut) fUsePhiCut=0;
1058                 fMaxPhiCut=10000;
1059                 break;
1060         case 1:
1061                 if (!fUsePhiCut) fUsePhiCut=1;
1062                 fMaxPhiCut=3.125;
1063                 break;
1064         default:
1065                 AliError(Form("Max Phi Cut not defined %d",maxPhi));
1066                 return kFALSE;
1067         }
1068         return kTRUE;
1069 }
1070
1071 //___________________________________________________________________
1072 Bool_t AliCaloPhotonCuts::SetDistanceToBadChannelCut(Int_t distanceToBadChannel)
1073 {
1074         switch(distanceToBadChannel){
1075         case 0:
1076                 fUseDistanceToBadChannel=0;
1077                 fMinDistanceToBadChannel=0;
1078                 break;
1079         case 1:
1080                 if (!fUseDistanceToBadChannel) fUseDistanceToBadChannel=1;
1081                 fMinDistanceToBadChannel=5;
1082                 break;
1083         default:
1084                 AliError(Form("minimum distance to bad channel Cut not defined %d",distanceToBadChannel));
1085                 return kFALSE;
1086         }
1087         return kTRUE;
1088 }
1089
1090 //___________________________________________________________________
1091 Bool_t AliCaloPhotonCuts::SetTimingCut(Int_t timing)
1092 {
1093         switch(timing){
1094         case 0:
1095                 fUseTimeDiff=0;
1096                 fMaxTimeDiff=500;
1097                 break;
1098         case 1:
1099                 if (!fUseTimeDiff) fUseTimeDiff=1;
1100                 fMaxTimeDiff=10e-7; //1000ns
1101                 break;
1102         case 2:
1103                 if (!fUseTimeDiff) fUseTimeDiff=1;
1104                 fMaxTimeDiff=50e-8; //500ns
1105                 break;
1106         case 3:
1107                 if (!fUseTimeDiff) fUseTimeDiff=1;
1108                 fMaxTimeDiff=20e-8; //200ns
1109                 break;
1110         case 4:
1111                 if (!fUseTimeDiff) fUseTimeDiff=1;
1112                 fMaxTimeDiff=10e-8; //100ns
1113                 break;
1114         case 5:
1115                 if (!fUseTimeDiff) fUseTimeDiff=1;
1116                 fMaxTimeDiff=50e-9; //50ns
1117                 break;
1118
1119         default:
1120                 AliError(Form("Timing Cut not defined %d",timing));
1121                 return kFALSE;
1122         }
1123         return kTRUE;
1124 }
1125
1126 //___________________________________________________________________
1127 Bool_t AliCaloPhotonCuts::SetTrackMatchingCut(Int_t trackMatching)
1128 {
1129         switch(trackMatching){
1130         case 0:
1131                 fUseDistTrackToCluster=0;
1132                 fMinDistTrackToCluster=0;
1133                 break;
1134         case 1:
1135                 if (!fUseDistTrackToCluster) fUseDistTrackToCluster=1;
1136                 fMinDistTrackToCluster=0.04;
1137                 break;
1138         case 2:
1139                 if (!fUseDistTrackToCluster) fUseDistTrackToCluster=1;
1140                 fMinDistTrackToCluster=0.05;
1141                 break;
1142         case 3:
1143                 if (!fUseDistTrackToCluster) fUseDistTrackToCluster=1;
1144                 fMinDistTrackToCluster=0.1;
1145                 break;
1146         case 4:
1147                 if (!fUseDistTrackToCluster) fUseDistTrackToCluster=1;
1148                 fMinDistTrackToCluster=0.13;
1149                 break;
1150         case 5:
1151                 if (!fUseDistTrackToCluster) fUseDistTrackToCluster=1;
1152                 fMinDistTrackToCluster=0.15;
1153                 break;
1154         case 6:
1155                 if (!fUseDistTrackToCluster) fUseDistTrackToCluster=1;
1156                 fMinDistTrackToCluster=0.2;
1157                 break;
1158         case 7:
1159                 if (!fUseDistTrackToCluster) fUseDistTrackToCluster=1;
1160                 fMinDistTrackToCluster=0.3;
1161                 break;
1162         case 8:
1163                 if (!fUseDistTrackToCluster) fUseDistTrackToCluster=1;
1164                 fMinDistTrackToCluster=0.4;
1165                 break;
1166         case 9:
1167                 if (!fUseDistTrackToCluster) fUseDistTrackToCluster=1;
1168                 fMinDistTrackToCluster=0.5;
1169                 break;
1170
1171         default:
1172                 AliError(Form("Track Matching Cut not defined %d",trackMatching));
1173                 return kFALSE;
1174         }
1175         return kTRUE;
1176 }
1177
1178 //___________________________________________________________________
1179 Bool_t AliCaloPhotonCuts::SetExoticCellCut(Int_t exoticCell)
1180 {
1181         switch(exoticCell){
1182         case 0:
1183                 fUseExoticCell=0;
1184                 fExoticCell=0;
1185                 break;
1186         case 1:
1187                 if (!fUseExoticCell) fUseExoticCell=1;
1188                 fExoticCell=5;
1189                 break;
1190         default:
1191                 AliError(Form("Exotic cell Cut not defined %d",exoticCell));
1192                 return kFALSE;
1193         }
1194         return kTRUE;
1195 }
1196
1197 //___________________________________________________________________
1198 Bool_t AliCaloPhotonCuts::SetMinEnergyCut(Int_t minEnergy)
1199 {
1200         switch(minEnergy){
1201         case 0:
1202                 if (!fUseMinEnergy) fUseMinEnergy=0;
1203                 fMinEnergy=0;
1204                 break;
1205         case 1:
1206                 if (!fUseMinEnergy) fUseMinEnergy=1;
1207                 fMinEnergy=0.05;
1208                 break;
1209         case 2:
1210                 if (!fUseMinEnergy) fUseMinEnergy=1;
1211                 fMinEnergy=0.1;
1212                 break;
1213         case 3:
1214                 if (!fUseMinEnergy) fUseMinEnergy=1;
1215                 fMinEnergy=0.15;
1216                 break;
1217         case 4:
1218                 if (!fUseMinEnergy) fUseMinEnergy=1;
1219                 fMinEnergy=0.2;
1220                 break;
1221         case 5:
1222                 if (!fUseMinEnergy) fUseMinEnergy=1;
1223                 fMinEnergy=0.3;
1224                 break;
1225         case 6:
1226                 if (!fUseMinEnergy) fUseMinEnergy=1;
1227                 fMinEnergy=0.5;
1228                 break;
1229         case 7:
1230                 if (!fUseMinEnergy) fUseMinEnergy=1;
1231                 fMinEnergy=0.75;
1232                 break;
1233         case 8:
1234                 if (!fUseMinEnergy) fUseMinEnergy=1;
1235                 fMinEnergy=1.;
1236                 break;
1237         case 9:
1238                 if (!fUseMinEnergy) fUseMinEnergy=1;
1239                 fMinEnergy=1.25;
1240                 break;
1241         default:
1242                 AliError(Form("Minimum Energy Cut not defined %d",minEnergy));
1243                 return kFALSE;
1244         }
1245         return kTRUE;
1246 }
1247
1248 //___________________________________________________________________
1249 Bool_t AliCaloPhotonCuts::SetMinNCellsCut(Int_t minNCells)
1250 {
1251         switch(minNCells){
1252         case 0:
1253                 if (!fUseNCells) fUseNCells=0;
1254                 fMinNCells=0;
1255                 break;
1256         case 1:
1257                 if (!fUseNCells) fUseNCells=1;
1258                 fMinNCells=1;
1259                 break;
1260         case 2:
1261                 if (!fUseNCells) fUseNCells=1;
1262                 fMinNCells=2;
1263                 break;
1264         case 3:
1265                 if (!fUseNCells) fUseNCells=1;
1266                 fMinNCells=3;
1267                 break;
1268         case 4:
1269                 if (!fUseNCells) fUseNCells=1;
1270                 fMinNCells=4;
1271                 break;
1272         case 5:
1273                 if (!fUseNCells) fUseNCells=1;
1274                 fMinNCells=5;
1275                 break;
1276         case 6:
1277                 if (!fUseNCells) fUseNCells=1;
1278                 fMinNCells=6;
1279                 break;
1280
1281         default:
1282                 AliError(Form("Min N cells Cut not defined %d",minNCells));
1283                 return kFALSE;
1284         }
1285         return kTRUE;
1286 }
1287
1288 //___________________________________________________________________
1289 Bool_t AliCaloPhotonCuts::SetMaxM02(Int_t maxM02)
1290 {
1291         switch(maxM02){
1292         case 0:
1293                 if (!fUseM02) fUseM02=0;
1294                 fMaxM02=100;
1295                 break;
1296         case 1:
1297                 if (!fUseM02) fUseM02=1;
1298                 fMaxM02=1.;
1299                 break;
1300         case 2:
1301                 if (!fUseM02) fUseM02=1;
1302                 fMaxM02=0.7;
1303                 break;
1304         case 3:
1305                 if (!fUseM02) fUseM02=1;
1306                 fMaxM02=0.5;
1307                 break;
1308         case 4:
1309                 if (!fUseM02) fUseM02=1;
1310                 fMaxM02=0.4;
1311                 break;
1312         default:
1313                 AliError(Form("Max M02 Cut not defined %d",maxM02));
1314                 return kFALSE;
1315         }
1316         return kTRUE;
1317 }
1318
1319 //___________________________________________________________________
1320 Bool_t AliCaloPhotonCuts::SetMinM02(Int_t minM02)
1321 {
1322         switch(minM02){
1323         case 0:
1324                 if (!fUseM02) fUseM02=0;
1325                 fMinM02=0;
1326                 break;
1327         case 1:
1328                 if (!fUseM02) fUseM02=1;
1329                 fMinM02=0.002;
1330                 break;
1331         default:
1332                 AliError(Form("Min M02 not defined %d",minM02));
1333                 return kFALSE;
1334         }
1335         return kTRUE;
1336 }
1337
1338 //___________________________________________________________________
1339 Bool_t AliCaloPhotonCuts::SetMaxM20(Int_t maxM20)
1340 {
1341         switch(maxM20){
1342         case 0:
1343                 if (!fUseM20) fUseM20=0;
1344                 fMaxM20=100;
1345                 break;
1346         case 1:
1347                 if (!fUseM20) fUseM20=1;
1348                 fMaxM20=0.5;
1349                 break;
1350         default:
1351                 AliError(Form("Max M20 Cut not defined %d",maxM20));
1352                 return kFALSE;
1353         }
1354         return kTRUE;
1355 }
1356
1357 //___________________________________________________________________
1358 Bool_t AliCaloPhotonCuts::SetMinM20(Int_t minM20)
1359 {
1360         switch(minM20){
1361         case 0:
1362                 if (!fUseM20) fUseM20=0;
1363                 fMinM20=0;
1364                 break;
1365         case 1:
1366                 if (!fUseM20) fUseM20=1;
1367                 fMinM20=0.002;
1368                 break;
1369         default:
1370                 AliError(Form("Min M20 Cut not defined %d",minM20));
1371                 return kFALSE;
1372         }
1373         return kTRUE;
1374 }
1375
1376 //___________________________________________________________________
1377 Bool_t AliCaloPhotonCuts::SetDispersion(Int_t dispersion)
1378 {
1379         switch(dispersion){
1380         case 0:
1381                 if (!fUseDispersion) fUseDispersion=0;
1382                 fMaxDispersion =100;
1383                 break;
1384         case 1:
1385                 if (!fUseDispersion) fUseDispersion=1;
1386                 fMaxDispersion=2.;
1387                 break;
1388         default:
1389                 AliError(Form("Maximum Dispersion Cut not defined %d",dispersion));
1390                 return kFALSE;
1391         }
1392         return kTRUE;
1393 }
1394
1395 //___________________________________________________________________
1396 Bool_t AliCaloPhotonCuts::SetNLM(Int_t nlm)
1397 {
1398         switch(nlm){
1399         case 0:
1400                 if (!fUseNLM) fUseNLM=0;
1401                 fMinNLM =0;
1402                 fMaxNLM =100;
1403                 break;
1404         case 1:
1405                 if (!fUseNLM) fUseNLM=1;
1406                 fMinNLM =0;
1407                 fMaxNLM =1;
1408                 break;
1409         default:
1410                 AliError(Form("NLM Cut not defined %d",nlm));
1411                 return kFALSE;
1412         }
1413         return kTRUE;
1414 }
1415
1416 ///________________________________________________________________________
1417 TString AliCaloPhotonCuts::GetCutNumber(){
1418    // returns TString with current cut number
1419    TString a(kNCuts);
1420    for(Int_t ii=0;ii<kNCuts;ii++){
1421       a.Append(Form("%d",fCuts[ii]));
1422    }
1423    return a;
1424 }
1425
1426