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 #include "AliTrackerBase.h"
  47
  48 class iostream;
  49
  50 using namespace std;
  51
  52 ClassImp(AliCaloPhotonCuts)
  53
  54
  55 const char* AliCaloPhotonCuts::fgkCutNames[AliCaloPhotonCuts::kNCuts] = {
  56         "ClusterType",          //0   0: all,    1: EMCAL,   2: PHOS
  57         "EtaMin",               //1   0: -10,    1: -0.6687, 2: -0,5, 3: -2
  58         "EtaMax",               //2   0: 10,     1: 0.66465, 2: 0.5,  3: 2
  59         "PhiMin",               //3   0: -10000, 1: 1.39626
  60         "PhiMax",               //4   0: 10000, 1: 3.125
  61         "DistanceToBadChannel", //5   0: 0,      1: 5
  62         "Timing",               //6   0: no cut
  63         "TrackMatching",        //7   0: 0,      1: 5
  64         "ExoticCell",           //8   0: no cut
  65         "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)
  66         "MinNCells",            //10  0: no cut, 1: 1,       2: 2,    3: 3,    4: 4,   5: 5,   6: 6
  67         "MinM02",                               //11
  68         "MaxM02",                               //12
  69         "MinM20",                               //13
  70         "MaxM20",                               //14
  71         "MaximumDispersion",    //15
  72         "NLM"                                   //16
  73 };
  74
  75
  76 //________________________________________________________________________
  77 AliCaloPhotonCuts::AliCaloPhotonCuts(const char *name,const char *title) :
  78         AliAnalysisCuts(name,title),
  79         fHistograms(NULL),
  80         fClusterType(0),
  81         fMinEtaCut(-10),
  82         fMaxEtaCut(10),
  83         fUseEtaCut(0),
  84         fMinPhiCut(-10000),
  85         fMaxPhiCut(-10000),
  86         fUsePhiCut(0),
  87         fMinDistanceToBadChannel(0),
  88         fUseDistanceToBadChannel(0),
  89         fMaxTimeDiff(10e10),
  90         fUseTimeDiff(0),
  91     fMaxDistTrackToClusterEta(0),
  92     fMinDistTrackToClusterPhi(0),
  93     fMaxDistTrackToClusterPhi(0),
  94         fUseDistTrackToCluster(0),
  95     fExtendedMatching(kFALSE),
  96         fExoticCell(0),
  97         fUseExoticCell(0),
  98         fMinEnergy(0),
  99         fUseMinEnergy(0),
 100         fMinNCells(0),
 101         fUseNCells(0),
 102         fMaxM02(1000),
 103         fMinM02(0),
 104         fUseM02(0),
 105         fMaxM20(1000),
 106         fMinM20(0),
 107         fUseM20(0),
 108         fMaxDispersion(1000),
 109         fUseDispersion(0),
 110         fMinNLM(0),
 111         fMaxNLM(1000),
 112         fUseNLM(0),
 113         fCutString(NULL),
 114         fHistCutIndex(NULL),
 115         fHistAcceptanceCuts(NULL),
 116         fHistClusterIdentificationCuts(NULL),
 117         fHistClusterEtavsPhiBeforeAcc(NULL),
 118         fHistClusterEtavsPhiAfterAcc(NULL),
 119         fHistClusterEtavsPhiAfterQA(NULL),
 120     //fHistDistanceToBadChannelBeforeAcc(NULL),
 121     //fHistDistanceToBadChannelAfterAcc(NULL),
 122         fHistClusterTimevsEBeforeQA(NULL),
 123         fHistClusterTimevsEAfterQA(NULL),
 124     //fHistExoticCellBeforeQA(NULL),
 125     //fHistExoticCellAfterQA(NULL),
 126     //fHistNMatchedTracks(NULL),
 127         fHistEnergyOfClusterBeforeQA(NULL),
 128         fHistEnergyOfClusterAfterQA(NULL),
 129         fHistNCellsBeforeQA(NULL),
 130         fHistNCellsAfterQA(NULL),
 131         fHistM02BeforeQA(NULL),
 132         fHistM02AfterQA(NULL),
 133         fHistM20BeforeQA(NULL),
 134         fHistM20AfterQA(NULL),
 135         fHistDispersionBeforeQA(NULL),
 136     fHistDispersionAfterQA(NULL),
 137     //fHistNLMBeforeQA(NULL),
 138     //fHistNLMAfterQA(NULL),
 139     fHistClusterRBeforeQA(NULL),
 140     fHistClusterRAfterQA(NULL),
 141     fHistClusterdEtadPhiBeforeQA(NULL),
 142     fHistClusterdEtadPhiAfterQA(NULL),
 143     fHistDistanceTrackToClusterBeforeQA(NULL),
 144     fHistDistanceTrackToClusterAfterQA(NULL),
 145     fHistClusterdEtadPhiPosTracksBeforeQA(NULL),
 146     fHistClusterdEtadPhiNegTracksBeforeQA(NULL),
 147     fHistClusterdEtadPhiPosTracksAfterQA(NULL),
 148     fHistClusterdEtadPhiNegTracksAfterQA(NULL),
 149     fHistClusterdEtadPhiPosTracksP_000_075BeforeQA(NULL),
 150     fHistClusterdEtadPhiPosTracksP_075_125BeforeQA(NULL),
 151     fHistClusterdEtadPhiPosTracksP_125_999BeforeQA(NULL),
 152     fHistClusterdEtadPhiNegTracksP_000_075BeforeQA(NULL),
 153     fHistClusterdEtadPhiNegTracksP_075_125BeforeQA(NULL),
 154     fHistClusterdEtadPhiNegTracksP_125_999BeforeQA(NULL),
 155     fHistClusterdEtadPtBeforeQA(NULL),
 156     fHistClusterdPhidPtBeforeQA(NULL),
 157     fHistClusterM20Pt_dPhiBeforeQA(NULL),
 158     fHistClusterM02Pt_dPhiBeforeQA(NULL),
 159     fHistClusterM20M02BeforeQA(NULL),
 160     fHistClusterM20M02AfterQA(NULL)
 161 {
 162    for(Int_t jj=0;jj<kNCuts;jj++){fCuts[jj]=0;}
 163    fCutString=new TObjString((GetCutNumber()).Data());
 164 }
 165
 166 //________________________________________________________________________
 167 AliCaloPhotonCuts::AliCaloPhotonCuts(const AliCaloPhotonCuts &ref) :
 168    AliAnalysisCuts(ref),
 169         fHistograms(NULL),
 170         fClusterType(ref.fClusterType),
 171         fMinEtaCut(ref.fMinEtaCut),
 172         fMaxEtaCut(ref.fMaxEtaCut),
 173         fUseEtaCut(ref.fUseEtaCut),
 174         fMinPhiCut(ref.fMinPhiCut),
 175         fMaxPhiCut(ref.fMaxPhiCut),
 176         fUsePhiCut(ref.fUsePhiCut),
 177         fMinDistanceToBadChannel(ref.fMinDistanceToBadChannel),
 178         fUseDistanceToBadChannel(ref.fUseDistanceToBadChannel),
 179         fMaxTimeDiff(ref.fMaxTimeDiff),
 180         fUseTimeDiff(ref.fUseTimeDiff),
 181     fMaxDistTrackToClusterEta(ref.fMaxDistTrackToClusterEta),
 182     fMinDistTrackToClusterPhi(ref.fMinDistTrackToClusterPhi),
 183     fMaxDistTrackToClusterPhi(ref.fMaxDistTrackToClusterPhi),
 184         fUseDistTrackToCluster(ref.fUseDistTrackToCluster),
 185     fExtendedMatching(ref.fExtendedMatching),
 186         fExoticCell(ref.fExoticCell),
 187         fUseExoticCell(ref.fUseExoticCell),
 188         fMinEnergy(ref.fMinEnergy),
 189         fUseMinEnergy(ref.fUseMinEnergy),
 190         fMinNCells(ref.fMinNCells),
 191         fUseNCells(ref.fUseNCells),
 192         fMaxM02(ref.fMaxM02),
 193         fMinM02(ref.fMinM02),
 194         fUseM02(ref.fUseM02),
 195         fMaxM20(ref.fMaxM20),
 196         fMinM20(ref.fMinM20),
 197         fUseM20(ref.fUseDispersion),
 198         fMaxDispersion(ref.fMaxDispersion),
 199         fUseDispersion(ref.fUseDispersion),
 200         fMinNLM(ref.fMinNLM),
 201         fMaxNLM(ref.fMaxNLM),
 202         fUseNLM(ref.fUseNLM),
 203         fCutString(NULL),
 204         fHistCutIndex(NULL),
 205         fHistAcceptanceCuts(NULL),
 206         fHistClusterIdentificationCuts(NULL),
 207         fHistClusterEtavsPhiBeforeAcc(NULL),
 208         fHistClusterEtavsPhiAfterAcc(NULL),
 209         fHistClusterEtavsPhiAfterQA(NULL),
 210     //fHistDistanceToBadChannelBeforeAcc(NULL),
 211     //fHistDistanceToBadChannelAfterAcc(NULL),
 212         fHistClusterTimevsEBeforeQA(NULL),
 213         fHistClusterTimevsEAfterQA(NULL),
 214     //fHistExoticCellBeforeQA(NULL),
 215     //fHistExoticCellAfterQA(NULL),
 216     //fHistNMatchedTracks(NULL),
 217         fHistEnergyOfClusterBeforeQA(NULL),
 218         fHistEnergyOfClusterAfterQA(NULL),
 219         fHistNCellsBeforeQA(NULL),
 220         fHistNCellsAfterQA(NULL),
 221         fHistM02BeforeQA(NULL),
 222         fHistM02AfterQA(NULL),
 223         fHistM20BeforeQA(NULL),
 224         fHistM20AfterQA(NULL),
 225         fHistDispersionBeforeQA(NULL),
 226     fHistDispersionAfterQA(NULL),
 227     //fHistNLMBeforeQA(NULL),
 228     //fHistNLMAfterQA(NULL),
 229     fHistClusterRBeforeQA(NULL),
 230     fHistClusterRAfterQA(NULL),
 231     fHistClusterdEtadPhiBeforeQA(NULL),
 232     fHistClusterdEtadPhiAfterQA(NULL),
 233     fHistDistanceTrackToClusterBeforeQA(NULL),
 234     fHistDistanceTrackToClusterAfterQA(NULL),
 235     fHistClusterdEtadPhiPosTracksBeforeQA(NULL),
 236     fHistClusterdEtadPhiNegTracksBeforeQA(NULL),
 237     fHistClusterdEtadPhiPosTracksAfterQA(NULL),
 238     fHistClusterdEtadPhiNegTracksAfterQA(NULL),
 239     fHistClusterdEtadPhiPosTracksP_000_075BeforeQA(NULL),
 240     fHistClusterdEtadPhiPosTracksP_075_125BeforeQA(NULL),
 241     fHistClusterdEtadPhiPosTracksP_125_999BeforeQA(NULL),
 242     fHistClusterdEtadPhiNegTracksP_000_075BeforeQA(NULL),
 243     fHistClusterdEtadPhiNegTracksP_075_125BeforeQA(NULL),
 244     fHistClusterdEtadPhiNegTracksP_125_999BeforeQA(NULL),
 245     fHistClusterdEtadPtBeforeQA(NULL),
 246     fHistClusterdPhidPtBeforeQA(NULL),
 247     fHistClusterM20Pt_dPhiBeforeQA(NULL),
 248     fHistClusterM02Pt_dPhiBeforeQA(NULL),
 249     fHistClusterM20M02BeforeQA(NULL),
 250     fHistClusterM20M02AfterQA(NULL)
 251 {
 252    // Copy Constructor
 253    for(Int_t jj=0;jj<kNCuts;jj++){fCuts[jj]=ref.fCuts[jj];}
 254    fCutString=new TObjString((GetCutNumber()).Data());
 255
 256 }
 257
 258
 259 //________________________________________________________________________
 260 AliCaloPhotonCuts::~AliCaloPhotonCuts() {
 261    // Destructor
 262    //Deleting fHistograms leads to seg fault it it's added to output collection of a task
 263    // if(fHistograms)
 264    //    delete fHistograms;
 265    // fHistograms = NULL;
 266    if(fCutString != NULL){
 267       delete fCutString;
 268       fCutString = NULL;
 269    }
 270 }
 271
 272 //________________________________________________________________________
 273 void AliCaloPhotonCuts::InitCutHistograms(TString name){
 274
 275         // Initialize Cut Histograms for QA (only initialized and filled if function is called)
 276         TH1::AddDirectory(kFALSE);
 277
 278         if(fHistograms != NULL){
 279                 delete fHistograms;
 280                 fHistograms=NULL;
 281         }
 282         if(fHistograms==NULL){
 283                 fHistograms=new TList();
 284                 fHistograms->SetOwner(kTRUE);
 285                 if(name=="")fHistograms->SetName(Form("CaloCuts_%s",GetCutNumber().Data()));
 286                 else fHistograms->SetName(Form("%s_%s",name.Data(),GetCutNumber().Data()));
 287         }
 288
 289         // IsPhotonSelected
 290         fHistCutIndex=new TH1F(Form("IsPhotonSelected %s",GetCutNumber().Data()),"IsPhotonSelected",5,-0.5,4.5);
 291         fHistCutIndex->GetXaxis()->SetBinLabel(kPhotonIn+1,"in");
 292         fHistCutIndex->GetXaxis()->SetBinLabel(kDetector+1,"detector");
 293         fHistCutIndex->GetXaxis()->SetBinLabel(kAcceptance+1,"acceptance");
 294         fHistCutIndex->GetXaxis()->SetBinLabel(kClusterQuality+1,"cluster QA");
 295         fHistCutIndex->GetXaxis()->SetBinLabel(kPhotonOut+1,"out");
 296         fHistograms->Add(fHistCutIndex);
 297
 298         // Acceptance Cuts
 299         fHistAcceptanceCuts=new TH1F(Form("AcceptanceCuts %s",GetCutNumber().Data()),"AcceptanceCuts",5,-0.5,4.5);
 300         fHistAcceptanceCuts->GetXaxis()->SetBinLabel(1,"in");
 301         fHistAcceptanceCuts->GetXaxis()->SetBinLabel(2,"eta");
 302         fHistAcceptanceCuts->GetXaxis()->SetBinLabel(3,"phi");
 303         fHistAcceptanceCuts->GetXaxis()->SetBinLabel(4,"distance to bad channel");
 304         fHistAcceptanceCuts->GetXaxis()->SetBinLabel(5,"out");
 305         fHistograms->Add(fHistAcceptanceCuts);
 306
 307         // Cluster Cuts
 308         fHistClusterIdentificationCuts =new TH1F(Form("ClusterQualityCuts %s",GetCutNumber().Data()),"ClusterQualityCuts",11,-0.5,10.5);
 309         fHistClusterIdentificationCuts->GetXaxis()->SetBinLabel(1,"in");
 310         fHistClusterIdentificationCuts->GetXaxis()->SetBinLabel(2,"timing");
 311         fHistClusterIdentificationCuts->GetXaxis()->SetBinLabel(3,"track matching");
 312         fHistClusterIdentificationCuts->GetXaxis()->SetBinLabel(4,"Exotics");
 313         fHistClusterIdentificationCuts->GetXaxis()->SetBinLabel(5,"minimum energy");
 314         fHistClusterIdentificationCuts->GetXaxis()->SetBinLabel(6,"minimum NCells");
 315         fHistClusterIdentificationCuts->GetXaxis()->SetBinLabel(7,"M02");
 316         fHistClusterIdentificationCuts->GetXaxis()->SetBinLabel(8,"M20");
 317         fHistClusterIdentificationCuts->GetXaxis()->SetBinLabel(9,"dispersion");
 318         fHistClusterIdentificationCuts->GetXaxis()->SetBinLabel(10,"NLM");
 319         fHistClusterIdentificationCuts->GetXaxis()->SetBinLabel(11,"out");
 320         fHistograms->Add(fHistClusterIdentificationCuts);
 321
 322         // Acceptance related histogramms
 323     const Int_t nEmcalEtaBins = 96;
 324     const Int_t nEmcalPhiBins = 124;
 325     Float_t EmcalEtaBins[nEmcalEtaBins+1] = {-0.66687,-0.653,-0.63913,-0.62526,-0.61139,-0.59752,-0.58365,-0.56978,-0.55591,-0.54204,-0.52817,-0.5143,-0.50043,-0.48656,-0.47269,-0.45882,-0.44495,-0.43108,-0.41721,-0.40334,-0.38947,-0.3756,-0.36173,-0.34786,-0.33399,-0.32012,-0.30625,-0.29238,-0.27851,-0.26464,-0.25077,-0.2369,-0.22303,-0.20916,-0.19529,-0.18142,-0.16755,-0.15368,-0.13981,-0.12594,-0.11207,-0.0982,-0.08433,-0.07046,-0.05659,-0.04272,-0.02885,-0.01498,-0.00111,0.01276,0.02663,0.0405,0.05437,0.06824,0.08211,0.09598,0.10985,0.12372,0.13759,0.15146,0.16533,0.1792,0.19307,0.20694,0.22081,0.23468,0.24855,0.26242,0.27629,0.29016,0.30403,0.3179,0.33177,0.34564,0.35951,0.37338,0.38725,0.40112,0.41499,0.42886,0.44273,0.4566,0.47047,0.48434,0.49821,0.51208,0.52595,0.53982,0.55369,0.56756,0.58143,0.5953,0.60917,0.62304,0.63691,0.65078,0.66465};
 326     Float_t EmcalPhiBins[nEmcalPhiBins+1] = {1.408,1.4215,1.435,1.4485,1.462,1.4755,1.489,1.5025,1.516,1.5295,1.543,1.5565,1.57,1.5835,1.597,1.6105,1.624,1.6375,1.651,1.6645,1.678,1.6915,1.705,1.7185,1.732,1.758,1.7715,1.785,1.7985,1.812,1.8255,1.839,1.8525,1.866,1.8795,1.893,1.9065,1.92,1.9335,1.947,1.9605,1.974,1.9875,2.001,2.0145,2.028,2.0415,2.055,2.0685,2.082,2.108,2.1215,2.135,2.1485,2.162,2.1755,2.189,2.2025,2.216,2.2295,2.243,2.2565,2.27,2.2835,2.297,2.3105,2.324,2.3375,2.351,2.3645,2.378,2.3915,2.405,2.4185,2.432,2.456,2.4695,2.483,2.4965,2.51,2.5235,2.537,2.5505,2.564,2.5775,2.591,2.6045,2.618,2.6315,2.645,2.6585,2.672,2.6855,2.699,2.7125,2.726,2.7395,2.753,2.7665,2.78,2.804,2.8175,2.831,2.8445,2.858,2.8715,2.885,2.8985,2.912,2.9255,2.939,2.9525,2.966,2.9795,2.993,3.0065,3.02,3.0335,3.047,3.0605,3.074,3.0875,3.101,3.1145,3.128};
 327
 328     const Int_t nPhosEtaBins = 56;
 329     const Int_t nPhosPhiBins = 192;
 330     const Float_t PhosEtaRange[2] = {-0.16, 0.16};
 331     const Float_t PhosPhiRange[2] = {4.5, 5.6};
 332
 333     if( GetClusterType() == 1 ){ //EMCAL
 334         fHistClusterEtavsPhiBeforeAcc=new TH2F(Form("EtaPhi_beforeAcceptance %s",GetCutNumber().Data()),"EtaPhi_beforeAcceptance",nEmcalPhiBins,EmcalPhiBins,nEmcalEtaBins,EmcalEtaBins);
 335         fHistograms->Add(fHistClusterEtavsPhiBeforeAcc);
 336         fHistClusterEtavsPhiAfterAcc=new TH2F(Form("EtaPhi_afterAcceptance %s",GetCutNumber().Data()),"EtaPhi_afterAcceptance",nEmcalPhiBins,EmcalPhiBins,nEmcalEtaBins,EmcalEtaBins);
 337         fHistograms->Add(fHistClusterEtavsPhiAfterAcc);
 338         fHistClusterEtavsPhiAfterQA=new TH2F(Form("EtaPhi_afterClusterQA %s",GetCutNumber().Data()),"EtaPhi_afterClusterQA",nEmcalPhiBins,EmcalPhiBins,nEmcalEtaBins,EmcalEtaBins);
 339         fHistograms->Add(fHistClusterEtavsPhiAfterQA);
 340     }
 341     else if( GetClusterType() == 2 ){ //PHOS
 342         fHistClusterEtavsPhiBeforeAcc=new TH2F(Form("EtaPhi_beforeAcceptance %s",GetCutNumber().Data()),"EtaPhi_beforeAcceptance",nPhosPhiBins,PhosPhiRange[0],PhosPhiRange[1],nPhosEtaBins,PhosEtaRange[0],PhosEtaRange[1]);
 343         fHistograms->Add(fHistClusterEtavsPhiBeforeAcc);
 344         fHistClusterEtavsPhiAfterAcc=new TH2F(Form("EtaPhi_afterAcceptance %s",GetCutNumber().Data()),"EtaPhi_afterAcceptance",nPhosPhiBins,PhosPhiRange[0],PhosPhiRange[1],nPhosEtaBins,PhosEtaRange[0],PhosEtaRange[1]);
 345         fHistograms->Add(fHistClusterEtavsPhiAfterAcc);
 346         fHistClusterEtavsPhiAfterQA=new TH2F(Form("EtaPhi_afterClusterQA %s",GetCutNumber().Data()),"EtaPhi_afterClusterQA",nPhosPhiBins,PhosPhiRange[0],PhosPhiRange[1],nPhosEtaBins,PhosEtaRange[0],PhosEtaRange[1]);
 347         fHistograms->Add(fHistClusterEtavsPhiAfterQA);
 348     }
 349     else if( GetClusterType() == 0 ){ //all
 350         fHistClusterEtavsPhiBeforeAcc=new TH2F(Form("EtaPhi_beforeAcceptance %s",GetCutNumber().Data()),"EtaPhi_beforeAcceptance",462,0,2*TMath::Pi(),110,-0.7,0.7);
 351         fHistograms->Add(fHistClusterEtavsPhiBeforeAcc);
 352         fHistClusterEtavsPhiAfterAcc=new TH2F(Form("EtaPhi_afterAcceptance %s",GetCutNumber().Data()),"EtaPhi_afterAcceptance",462,0,2*TMath::Pi(),110,-0.7,0.7);
 353         fHistograms->Add(fHistClusterEtavsPhiAfterAcc);
 354         fHistClusterEtavsPhiAfterQA=new TH2F(Form("EtaPhi_afterClusterQA %s",GetCutNumber().Data()),"EtaPhi_afterClusterQA",462,0,2*TMath::Pi(),110,-0.7,0.7);
 355         fHistograms->Add(fHistClusterEtavsPhiAfterQA);
 356     }
 357     else{AliError(Form("Cluster Type is not EMCAL nor PHOS: %i",GetClusterType()));}
 358
 359     //fHistDistanceToBadChannelBeforeAcc = new TH1F(Form("DistanceToBadChannel_beforeAcceptance %s",GetCutNumber().Data()),"DistanceToBadChannel_beforeAcceptance",200,0,40);
 360     //fHistograms->Add(fHistDistanceToBadChannelBeforeAcc);
 361     //fHistDistanceToBadChannelAfterAcc = new TH1F(Form("DistanceToBadChannel_afterAcceptance %s",GetCutNumber().Data()),"DistanceToBadChannel_afterAcceptance",200,0,40);
 362     //fHistograms->Add(fHistDistanceToBadChannelAfterAcc);
 363
 364         // Cluster quality related histograms
 365     Double_t timeMin = -2e-6;
 366     Double_t timeMax = 8e-6;
 367     if( GetClusterType() == 1 ){
 368                 timeMin = -2e-7;
 369                 timeMax = 12e-7;
 370     }
 371
 372     fHistClusterTimevsEBeforeQA=new TH2F(Form("ClusterTimeVsE_beforeClusterQA %s",GetCutNumber().Data()),"ClusterTimeVsE_beforeClusterQA",800,timeMin,timeMax,100,0,40);
 373         fHistograms->Add(fHistClusterTimevsEBeforeQA);
 374     fHistClusterTimevsEAfterQA=new TH2F(Form("ClusterTimeVsE_afterClusterQA %s",GetCutNumber().Data()),"ClusterTimeVsE_afterClusterQA",800,timeMin,timeMax,100,0,40);
 375         fHistograms->Add(fHistClusterTimevsEAfterQA);
 376     //fHistExoticCellBeforeQA=new TH2F(Form("ExoticCell_beforeClusterQA %s",GetCutNumber().Data()),"ExoticCell_beforeClusterQA",400,0,40,50,0.75,1);
 377     //fHistograms->Add(fHistExoticCellBeforeQA);
 378     //fHistExoticCellAfterQA=new TH2F(Form("ExoticCell_afterClusterQA %s",GetCutNumber().Data()),"ExoticCell_afterClusterQA",400,0,40,50,0.75,1);
 379     //fHistograms->Add(fHistExoticCellAfterQA);
 380     //fHistNMatchedTracks = new TH1F(Form("NMatchedTracks_%s",GetCutNumber().Data()),"NMatchedTracks",22,-1.5,20.5);
 381     //fHistograms->Add(fHistNMatchedTracks);
 382         fHistEnergyOfClusterBeforeQA = new TH1F(Form("EnergyOfCluster_beforeClusterQA %s",GetCutNumber().Data()),"EnergyOfCluster_beforeClusterQA",300,0,30);
 383         fHistograms->Add(fHistEnergyOfClusterBeforeQA);
 384         fHistEnergyOfClusterAfterQA = new TH1F(Form("EnergyOfCluster_afterClusterQA %s",GetCutNumber().Data()),"EnergyOfCluster_afterClusterQA",300,0,30);
 385         fHistograms->Add(fHistEnergyOfClusterAfterQA);
 386         fHistNCellsBeforeQA = new TH1F(Form("NCellPerCluster_beforeClusterQA %s",GetCutNumber().Data()),"NCellPerCluster_beforeClusterQA",50,0,50);
 387         fHistograms->Add(fHistNCellsBeforeQA);
 388         fHistNCellsAfterQA = new TH1F(Form("NCellPerCluster_afterClusterQA %s",GetCutNumber().Data()),"NCellPerCluster_afterClusterQA",50,0,50);
 389         fHistograms->Add(fHistNCellsAfterQA);
 390         fHistM02BeforeQA = new TH1F(Form("M02_beforeClusterQA %s",GetCutNumber().Data()),"M02_beforeClusterQA",400,0,5);
 391         fHistograms->Add(fHistM02BeforeQA);
 392         fHistM02AfterQA = new TH1F(Form("M02_afterClusterQA %s",GetCutNumber().Data()),"M02_afterClusterQA",400,0,5);
 393         fHistograms->Add(fHistM02AfterQA);
 394         fHistM20BeforeQA = new TH1F(Form("M20_beforeClusterQA %s",GetCutNumber().Data()),"M20_beforeClusterQA",400,0,2.5);
 395         fHistograms->Add(fHistM20BeforeQA);
 396         fHistM20AfterQA = new TH1F(Form("M20_afterClusterQA %s",GetCutNumber().Data()),"M20_afterClusterQA",400,0,2.5);
 397         fHistograms->Add(fHistM20AfterQA);
 398         fHistDispersionBeforeQA = new TH1F(Form("Dispersion_beforeClusterQA %s",GetCutNumber().Data()),"Dispersion_beforeClusterQA",100,0,4);
 399         fHistograms->Add(fHistDispersionBeforeQA);
 400         fHistDispersionAfterQA = new TH1F(Form("Dispersion_afterClusterQA %s",GetCutNumber().Data()),"Dispersion_afterClusterQA",100,0,4);
 401         fHistograms->Add(fHistDispersionAfterQA);
 402     //fHistNLMBeforeQA = new TH1F(Form("NLM_beforeClusterQA %s",GetCutNumber().Data()),"NLM_beforeClusterQA",10,0,10);
 403     //fHistograms->Add(fHistNLMBeforeQA);
 404     //fHistNLMAfterQA = new TH1F(Form("NLM_afterClusterQA %s",GetCutNumber().Data()),"NLM_afterClusterQA",10,0,10);
 405     //fHistograms->Add(fHistNLMAfterQA);
 406
 407     //TrackMatching histograms
 408     if(fUseDistTrackToCluster) {
 409         const Int_t nEtaBins = 300;
 410         const Int_t nPhiBins = 300;
 411         const Float_t EtaRange[2] = {-0.3, 0.3};
 412         const Float_t PhiRange[2] = {-0.3, 0.3};
 413
 414         fHistClusterRBeforeQA = new TH1F(Form("R_Cluster_beforeClusterQA %s",GetCutNumber().Data()),"R of cluster",200,400,500);
 415         fHistograms->Add(fHistClusterRBeforeQA);
 416         fHistClusterRAfterQA = new TH1F(Form("R_Cluster_afterClusterQA %s",GetCutNumber().Data()),"R of cluster_matched",200,400,500);
 417         fHistograms->Add(fHistClusterRAfterQA);
 418         fHistClusterdEtadPhiBeforeQA=new TH2F(Form("dEtaVsdPhi_beforeClusterQA %s",GetCutNumber().Data()),"dEtaVsdPhi_beforeClusterQA",nEtaBins,EtaRange[0],EtaRange[1],nPhiBins,PhiRange[0],PhiRange[1]);
 419         fHistograms->Add(fHistClusterdEtadPhiBeforeQA);
 420         fHistClusterdEtadPhiAfterQA=new TH2F(Form("dEtaVsdPhi_afterClusterQA %s",GetCutNumber().Data()),"dEtaVsdPhi_afterClusterQA",nEtaBins,EtaRange[0],EtaRange[1],nPhiBins,PhiRange[0],PhiRange[1]);
 421         fHistograms->Add(fHistClusterdEtadPhiAfterQA);
 422         fHistDistanceTrackToClusterBeforeQA = new TH1F(Form("DistanceToTrack_beforeClusterQA %s",GetCutNumber().Data()),"DistanceToTrack_beforeClusterQA",200,0,2);
 423         fHistograms->Add(fHistDistanceTrackToClusterBeforeQA);
 424         fHistDistanceTrackToClusterAfterQA = new TH1F(Form("DistanceToTrack_afterClusterQA %s",GetCutNumber().Data()),"DistanceToTrack_afterClusterQA",200,0,2);
 425         fHistograms->Add(fHistDistanceTrackToClusterAfterQA);
 426
 427         if(fExtendedMatching){
 428             fHistClusterdEtadPhiPosTracksBeforeQA = new TH2F(Form("dEtaVsdPhi_posTracks_beforeClusterQA %s",GetCutNumber().Data()),"dEtaVsdPhi_posTracks_beforeClusterQA",nEtaBins,EtaRange[0],EtaRange[1],nPhiBins,PhiRange[0],PhiRange[1]);
 429             fHistograms->Add(fHistClusterdEtadPhiPosTracksBeforeQA);
 430             fHistClusterdEtadPhiNegTracksBeforeQA = new TH2F(Form("dEtaVsdPhi_negTracks_beforeClusterQA %s",GetCutNumber().Data()),"dEtaVsdPhi_negTracks_beforeClusterQA",nEtaBins,EtaRange[0],EtaRange[1],nPhiBins,PhiRange[0],PhiRange[1]);
 431             fHistograms->Add(fHistClusterdEtadPhiNegTracksBeforeQA);
 432             fHistClusterdEtadPhiPosTracksAfterQA = new TH2F(Form("dEtaVsdPhi_posTracks_afterClusterQA %s",GetCutNumber().Data()),"dEtaVsdPhi_posTracks_afterClusterQA",nEtaBins,EtaRange[0],EtaRange[1],nPhiBins,PhiRange[0],PhiRange[1]);
 433             fHistograms->Add(fHistClusterdEtadPhiPosTracksAfterQA);
 434             fHistClusterdEtadPhiNegTracksAfterQA = new TH2F(Form("dEtaVsdPhi_negTracks_afterClusterQA %s",GetCutNumber().Data()),"dEtaVsdPhi_negTracks_afterClusterQA",nEtaBins,EtaRange[0],EtaRange[1],nPhiBins,PhiRange[0],PhiRange[1]);
 435             fHistograms->Add(fHistClusterdEtadPhiNegTracksAfterQA);
 436             fHistClusterdEtadPhiPosTracksP_000_075BeforeQA = new TH2F(Form("dEtaVsdPhi_posTracks_P<0.75_beforeClusterQA %s",GetCutNumber().Data()),"dEtaVsdPhi_posTracks_P<0.75_beforeClusterQA",nEtaBins,EtaRange[0],EtaRange[1],nPhiBins,PhiRange[0],PhiRange[1]);
 437             fHistograms->Add(fHistClusterdEtadPhiPosTracksP_000_075BeforeQA);
 438             fHistClusterdEtadPhiPosTracksP_075_125BeforeQA = new TH2F(Form("dEtaVsdPhi_posTracks_0.75<P<1.25_beforeClusterQA %s",GetCutNumber().Data()),"dEtaVsdPhi_posTracks_0.75<P<1.25_beforeClusterQA",nEtaBins,EtaRange[0],EtaRange[1],nPhiBins,PhiRange[0],PhiRange[1]);
 439             fHistograms->Add(fHistClusterdEtadPhiPosTracksP_075_125BeforeQA);
 440             fHistClusterdEtadPhiPosTracksP_125_999BeforeQA = new TH2F(Form("dEtaVsdPhi_posTracks_P>1.25_beforeClusterQA %s",GetCutNumber().Data()),"dEtaVsdPhi_posTracks_P>1.25_beforeClusterQA",nEtaBins,EtaRange[0],EtaRange[1],nPhiBins,PhiRange[0],PhiRange[1]);
 441             fHistograms->Add(fHistClusterdEtadPhiPosTracksP_125_999BeforeQA);
 442             fHistClusterdEtadPhiNegTracksP_000_075BeforeQA= new TH2F(Form("dEtaVsdPhi_negTrack_P<0.75_beforeClusterQA %s",GetCutNumber().Data()),"dEtaVsdPhi_negTrack_P<0.75_beforeClusterQA",nEtaBins,EtaRange[0],EtaRange[1],nPhiBins,PhiRange[0],PhiRange[1]);
 443             fHistograms->Add(fHistClusterdEtadPhiNegTracksP_000_075BeforeQA);
 444             fHistClusterdEtadPhiNegTracksP_075_125BeforeQA = new TH2F(Form("dEtaVsdPhi_negTracks_0.75<P<1.25_beforeClusterQA %s",GetCutNumber().Data()),"dEtaVsdPhi_negTracks_0.75<P<1.25_beforeClusterQA",nEtaBins,EtaRange[0],EtaRange[1],nPhiBins,PhiRange[0],PhiRange[1]);
 445             fHistograms->Add(fHistClusterdEtadPhiNegTracksP_075_125BeforeQA);
 446             fHistClusterdEtadPhiNegTracksP_125_999BeforeQA = new TH2F(Form("dEtaVsdPhi_negTracks_P>1.25_beforeClusterQA %s",GetCutNumber().Data()),"dEtaVsdPhi_negTracks_P>1.25_beforeClusterQA",nEtaBins,EtaRange[0],EtaRange[1],nPhiBins,PhiRange[0],PhiRange[1]);
 447             fHistograms->Add(fHistClusterdEtadPhiNegTracksP_125_999BeforeQA);
 448             fHistClusterdEtadPtBeforeQA = new TH2F(Form("dEtaVsPt_beforeClusterQA %s",GetCutNumber().Data()),"dEtaVsPt_beforeClusterQA",nEtaBins,EtaRange[0],EtaRange[1],250,0,25);
 449             fHistograms->Add(fHistClusterdEtadPtBeforeQA);
 450             fHistClusterdPhidPtBeforeQA = new TH2F(Form("dPhiVsPt_beforeClusterQA %s",GetCutNumber().Data()),"dPhiVsPt_beforeClusterQA",2*nPhiBins,2*PhiRange[0],2*PhiRange[1],250,0,25);
 451             fHistograms->Add(fHistClusterdPhidPtBeforeQA);
 452             fHistClusterM20Pt_dPhiBeforeQA = new TH2F(Form("M20VsPt_dPhi_beforeClusterQA %s",GetCutNumber().Data()),"M20VsPt_dPhi_beforeClusterQA",200,0,2.5,250,0,25);
 453             fHistograms->Add(fHistClusterM20Pt_dPhiBeforeQA);
 454             fHistClusterM02Pt_dPhiBeforeQA = new TH2F(Form("M02VsPt_dPhi_beforeClusterQA %s",GetCutNumber().Data()),"M20VsPt_dPhi_beforeClusterQA",400,0,5,250,0,25);
 455             fHistograms->Add(fHistClusterM02Pt_dPhiBeforeQA);
 456             fHistClusterM20M02BeforeQA = new TH2F(Form("M20VsM02_beforeClusterQA %s",GetCutNumber().Data()),"M20VsM02_beforeClusterQA",200,0,2.5,400,0,5);
 457             fHistograms->Add(fHistClusterM20M02BeforeQA);
 458             fHistClusterM20M02AfterQA = new TH2F(Form("M20VsM02_afterClusterQA %s",GetCutNumber().Data()),"M20VsM02_afterClusterQA",200,0,2.5,400,0,5);
 459             fHistograms->Add(fHistClusterM20M02AfterQA);
 460         }
 461     }
 462         TH1::AddDirectory(kTRUE);
 463 }
 464
 465
 466 ///________________________________________________________________________
 467 Bool_t AliCaloPhotonCuts::ClusterIsSelectedMC(TParticle *particle,AliStack *fMCStack){
 468    // MonteCarlo Photon Selection
 469
 470         if(!fMCStack)return kFALSE;
 471
 472         if (particle->GetPdgCode() == 22){
 473
 474                 if ( particle->Eta() < fMinEtaCut || particle->Eta() > fMaxEtaCut ) return kFALSE;
 475                 if ( particle->Phi() < fMinPhiCut || particle->Phi() > fMaxPhiCut ) return kFALSE;
 476
 477                 if(particle->GetMother(0) >-1 && fMCStack->Particle(particle->GetMother(0))->GetPdgCode() == 22){
 478                         return kFALSE; // no photon as mothers!
 479                 }
 480                 if(particle->GetMother(0) >= fMCStack->GetNprimary()){
 481                         return kFALSE; // the gamma has a mother, and it is not a primary particle
 482                 }
 483                 return kTRUE;
 484         }
 485         return kFALSE;
 486 }
 487 ///________________________________________________________________________
 488 Bool_t AliCaloPhotonCuts::ClusterIsSelectedAODMC(AliAODMCParticle *particle,TClonesArray *aodmcArray){
 489         // MonteCarlo Photon Selection
 490
 491         if(!aodmcArray)return kFALSE;
 492         if (particle->GetPdgCode() == 22){
 493                 if ( particle->Eta() < fMinEtaCut || particle->Eta() > fMaxEtaCut ) return kFALSE;
 494                 if ( particle->Phi() < fMinPhiCut || particle->Phi() > fMaxPhiCut ) return kFALSE;
 495                 if(particle->GetMother() > -1){
 496                         if((static_cast<AliAODMCParticle*>(aodmcArray->At(particle->GetMother())))->GetPdgCode() == 22){
 497                                 return kFALSE; // no photon as mothers!
 498                         }
 499                         if(!(static_cast<AliAODMCParticle*>(aodmcArray->At(particle->GetMother()))->IsPrimary())){
 500                                 return kFALSE; // the gamma has a mother, and it is not a primary particle
 501                         }
 502                 }
 503                 return kTRUE; // return in case of accepted gamma
 504         }
 505         return kFALSE;
 506 }
 507
 508
 509
 510 ///________________________________________________________________________
 511 // This function selects the clusters based on their quality criteria
 512 ///________________________________________________________________________
 513 Bool_t AliCaloPhotonCuts::ClusterQualityCuts(AliVCluster* cluster, AliVEvent *event, Bool_t isMC)
 514 {   // Specific Photon Cuts
 515     if(event){} // suppress warning
 516
 517         Int_t cutIndex = 0;
 518         if(fHistClusterIdentificationCuts)fHistClusterIdentificationCuts->Fill(cutIndex);
 519         cutIndex++;
 520
 521 //      Double_t minR = 999.0;
 522         // get the minimum radius of tracks to cluster
 523 //      if(fHistDistanceTrackToClusterBeforeQA || fHistDistanceTrackToClusterAfterQA){
 524 //              Float_t pos[3];
 525 //              cluster->GetPosition(pos);  // Get cluster position
 526 //
 527 //              TVector3 cp(pos);
 528 //              int NtrMatched = 0;
 529 //              NtrMatched = cluster->GetNTracksMatched();
 530 //              fHistNMatchedTracks->Fill(NtrMatched);
 531 //              //loop over tracks for Jet QA
 532 //              TList *l = event->GetList();
 533 //              TClonesArray *tracks = dynamic_cast<TClonesArray*>(l->FindObject("Tracks"));
 534 //              for(int itrack = 0; itrack < NtrMatched; itrack++){
 535 //                      AliVTrack *trackcluster = static_cast<AliVTrack*>(tracks->At(itrack));
 536 //                      if (! trackcluster) {
 537 //                              AliError(Form("Couldn't get ESD track %d\n", itrack));
 538 //                              continue;
 539 //                      }
 540 //                      Double_t dphi = -999.0;
 541 //                      Double_t deta = -999.0;
 542 //                      AliPicoTrack::GetEtaPhiDiff(trackcluster, cluster, dphi, deta);
 543 //                      cout << "here" << endl;
 544 //                      Double_t dr = sqrt(dphi*dphi + deta+deta);
 545 //                      if(dr < minR)
 546 //                              minR = dr;
 547 //              }//loop over tracks
 548 //      }
 549
 550         // Fill Histos before Cuts
 551         if(fHistClusterTimevsEBeforeQA) fHistClusterTimevsEBeforeQA->Fill(cluster->GetTOF(), cluster->E());
 552 //      if(fHistExoticCellBeforeQA) fHistExoticCellBeforeQA->Fill(cluster->E(), );
 553 //      if(fHistDistanceTrackToClusterBeforeQA) fHistDistanceTrackToClusterBeforeQA->Fill(minR);
 554         if(fHistEnergyOfClusterBeforeQA) fHistEnergyOfClusterBeforeQA->Fill(cluster->E());
 555         if(fHistNCellsBeforeQA) fHistNCellsBeforeQA->Fill(cluster->GetNCells());
 556         if(fHistM02BeforeQA) fHistM02BeforeQA->Fill(cluster->GetM02());
 557         if(fHistM20BeforeQA) fHistM20BeforeQA->Fill(cluster->GetM20());
 558         if(fHistDispersionBeforeQA) fHistDispersionBeforeQA->Fill(cluster->GetDispersion());
 559 //      if(fHistNLMBeforeQA) fHistNLMBeforeQA->Fill(cluster->GetNExMax());
 560
 561         // Check wether timing is ok
 562         if (fUseTimeDiff){
 563                 if(abs(cluster->GetTOF()) > fMaxTimeDiff && !isMC){
 564                         if(fHistClusterIdentificationCuts)fHistClusterIdentificationCuts->Fill(cutIndex); //1
 565                         return kFALSE;
 566                 }
 567         }
 568         cutIndex++; //2, next cut
 569
 570         // Minimum distance to track
 571 //      if (fUseDistTrackToCluster){
 572 //              Float_t pos[3];
 573 //              cluster->GetPosition(pos);  // Get cluster position
 574 //              TVector3 cp(pos);
 575 //              int NtrMatched = 0;
 576 //              NtrMatched = cluster->GetNTracksMatched();
 577 //              fHistNMatchedTracks->Fill(NtrMatched);
 578 //              if(NtrMatched > 0){
 579 //                      //loop over tracks for QA
 580 //                      TList *l = event->GetList();
 581 //                      TClonesArray *tracks = dynamic_cast<TClonesArray*>(l->FindObject("Tracks"));
 582 //
 583 //                      Double_t dphi = 999.0;
 584 //                      Double_t deta = 999.0;
 585 //                      Double_t dr2 = 999.0;
 586 //
 587 //                      for(int itrack = 0; itrack < NtrMatched; itrack++){
 588 //                              AliVTrack *trackcluster = NULL;
 589 //                              trackcluster = static_cast<AliVTrack*>(tracks->At(itrack));
 590 //                              if (! trackcluster) {
 591 //                              AliError(Form("Couldn't get ESD track %d\n", itrack));
 592 //                              continue;
 593 //                              }
 594 //                              AliPicoTrack::GetEtaPhiDiff(trackcluster, cluster, dphi, deta);
 595 //                              dr2 = dphi*dphi + deta+deta;
 596 //                              //cout << dr << endl;
 597 //                              if(dr2 < fMinDistTrackToCluster*fMinDistTrackToCluster){
 598 //              //        if(dphi < fMinDistTrackToCluster || deta < fMinDistTrackToCluster){
 599 //                              if(fHistClusterIdentificationCuts)fHistClusterIdentificationCuts->Fill(cutIndex); //2
 600 //                              return kFALSE;
 601 //                              }
 602 //
 603 //                      }//loop over tracks
 604 //              }
 605 // //           if(cluster->GetEmcCpvDistance() < fMinDistTrackToCluster){
 606 // //                   if(fHistClusterIdentificationCuts)fHistClusterIdentificationCuts->Fill(cutIndex); //2
 607 // //                   return kFALSE;
 608 // //           }
 609 //      }
 610         cutIndex++;//3, next cut
 611
 612         // exotic cell cut --IMPLEMENT LATER---
 613 //      if(!AcceptanceCuts(photon)){
 614 //              if(fHistClusterIdentificationCuts)fHistClusterIdentificationCuts->Fill(cutIndex); //3
 615 //              return kFALSE;
 616 //      }
 617         cutIndex++; //4, next cut
 618
 619         // minimum cell energy cut
 620         if (fUseMinEnergy){
 621                 if(cluster->E() < fMinEnergy){
 622                         if(fHistClusterIdentificationCuts)fHistClusterIdentificationCuts->Fill(cutIndex); //4
 623                         return kFALSE;
 624                 }
 625         }
 626         cutIndex++; //5, next cut
 627
 628         // minimum number of cells
 629         if (fUseNCells){
 630                 if(cluster->GetNCells() < fMinNCells) {
 631                         if(fHistClusterIdentificationCuts)fHistClusterIdentificationCuts->Fill(cutIndex); //5
 632                         return kFALSE;
 633                 }
 634         }
 635         cutIndex++; //6, next cut
 636
 637         // M02 cut
 638         if (fUseM02){
 639                 if( cluster->GetM02()< fMinM02 || cluster->GetM02() > fMaxM02 ) {
 640                         if(fHistClusterIdentificationCuts)fHistClusterIdentificationCuts->Fill(cutIndex); //6
 641                         return kFALSE;
 642                 }
 643         }
 644         cutIndex++; //7, next cut
 645
 646         // M20 cut
 647         if (fUseM20){
 648                 if( cluster->GetM20()< fMinM20 || cluster->GetM20() > fMaxM20 ) {
 649                         if(fHistClusterIdentificationCuts)fHistClusterIdentificationCuts->Fill(cutIndex); //7
 650                         return kFALSE;
 651                 }
 652         }
 653         cutIndex++; //8, next cut
 654
 655         // dispersion cut
 656         if (fUseDispersion){
 657                 if( cluster->GetDispersion()> fMaxDispersion) {
 658                         if(fHistClusterIdentificationCuts)fHistClusterIdentificationCuts->Fill(cutIndex); //8
 659                         return kFALSE;
 660                 }
 661         }
 662         cutIndex++; //9, next cut
 663
 664         // NLM cut --IMPLEMENT LATER---
 665 //      if (fUseNLM){
 666 //              if( cluster->GetDispersion()> fMaxDispersion) {
 667 //                      if(fHistClusterIdentificationCuts)fHistClusterIdentificationCuts->Fill(cutIndex); //9
 668 //                      return kFALSE;
 669 //              }
 670 //      }
 671         cutIndex++; //9, next cut
 672
 673         // DONE with selecting photons
 674         if(fHistClusterIdentificationCuts)fHistClusterIdentificationCuts->Fill(cutIndex); //10
 675
 676         // Histos after Cuts
 677 //  Double_t vertex[3] = {0,0,0};
 678 //      event->GetPrimaryVertex()->GetXYZ(vertex);
 679         // TLorentzvector with cluster
 680 //      TLorentzVector clusterVector;
 681 //      cluster->GetMomentum(clusterVector,vertex);
 682
 683     Float_t clusPos[3]={0,0,0};
 684     cluster->GetPosition(clusPos);
 685     TVector3 clusterVector(clusPos[0],clusPos[1],clusPos[2]);
 686         Double_t etaCluster = clusterVector.Eta();
 687         Double_t phiCluster = clusterVector.Phi();
 688         if (phiCluster < 0) phiCluster= phiCluster+ 2*TMath::Pi();
 689
 690         if(fHistClusterEtavsPhiAfterQA) fHistClusterEtavsPhiAfterQA->Fill(phiCluster,etaCluster);
 691         if(fHistClusterTimevsEAfterQA) fHistClusterTimevsEAfterQA->Fill(cluster->GetTOF(), cluster->E());
 692 //      if(fHistExoticCellAfterQA) fHistExoticCellAfterQA->Fill(cluster->E(), );
 693 //      if(fHistDistanceTrackToClusterAfterQA) fHistDistanceTrackToClusterAfterQA->Fill(minR);
 694         if(fHistEnergyOfClusterAfterQA) fHistEnergyOfClusterAfterQA->Fill(cluster->E());
 695         if(fHistNCellsAfterQA) fHistNCellsAfterQA->Fill(cluster->GetNCells());
 696         if(fHistM02AfterQA) fHistM02AfterQA->Fill(cluster->GetM02());
 697         if(fHistM20AfterQA) fHistM20AfterQA->Fill(cluster->GetM20());
 698         if(fHistDispersionAfterQA) fHistDispersionAfterQA->Fill(cluster->GetDispersion());
 699 //      if(fHistNLMBeforeQA) fHistNLMAfterQA->Fill(cluster->GetNExMax());
 700
 701         return kTRUE;
 702
 703 }
 704
 705
 706 ///________________________________________________________________________
 707 Bool_t AliCaloPhotonCuts::ClusterIsSelected(AliVCluster *cluster, AliVEvent * event, Bool_t isMC)
 708 {
 709         //Selection of Reconstructed photon clusters with Calorimeters
 710
 711         FillClusterCutIndex(kPhotonIn);
 712
 713 //  Double_t vertex[3] = {0,0,0};
 714 //      event->GetPrimaryVertex()->GetXYZ(vertex);
 715     // TLorentzvector with cluster
 716 //  TLorentzVector clusterVector;
 717 //      cluster->GetMomentum(clusterVector,vertex);
 718
 719     Float_t clusPos[3]={0,0,0};
 720     cluster->GetPosition(clusPos);
 721     TVector3 clusterVector(clusPos[0],clusPos[1],clusPos[2]);
 722         Double_t etaCluster = clusterVector.Eta();
 723         Double_t phiCluster = clusterVector.Phi();
 724         if (phiCluster < 0) phiCluster= phiCluster+ 2*TMath::Pi();
 725
 726         // Histos before cuts
 727         if(fHistClusterEtavsPhiBeforeAcc) fHistClusterEtavsPhiBeforeAcc->Fill(phiCluster,etaCluster);
 728
 729         // Cluster Selection - 0= accept any calo cluster
 730         if (fClusterType > 0){
 731                 //Select EMCAL cluster
 732                 if (fClusterType == 1 && !cluster->IsEMCAL()){
 733                         FillClusterCutIndex(kDetector);
 734                         return kFALSE;
 735                 }
 736                 //Select PHOS cluster
 737                 if (fClusterType == 2 && !cluster->IsPHOS()){
 738                         FillClusterCutIndex(kDetector);
 739                         return kFALSE;
 740                 }
 741         }
 742
 743         // Acceptance Cuts
 744         if(!AcceptanceCuts(cluster,event)){
 745                 FillClusterCutIndex(kAcceptance);
 746                 return kFALSE;
 747         }
 748         // Cluster Quality Cuts
 749         if(!ClusterQualityCuts(cluster,event,isMC)){
 750                 FillClusterCutIndex(kClusterQuality);
 751                 return kFALSE;
 752         }
 753
 754         // Photon passed cuts
 755         FillClusterCutIndex(kPhotonOut);
 756         return kTRUE;
 757 }
 758
 759
 760 ///________________________________________________________________________
 761 Bool_t AliCaloPhotonCuts::AcceptanceCuts(AliVCluster *cluster, AliVEvent* event)
 762 {
 763    // Exclude certain areas for photon reconstruction
 764     if(event){} // suppress warning
 765
 766         Int_t cutIndex=0;
 767         if(fHistAcceptanceCuts)fHistAcceptanceCuts->Fill(cutIndex);
 768         cutIndex++;
 769
 770
 771 //      Double_t vertex[3] = {0,0,0};
 772 //      event->GetPrimaryVertex()->GetXYZ(vertex);
 773         // TLorentzvector with cluster
 774 //      TLorentzVector clusterVector;
 775 //      cluster->GetMomentum(clusterVector,vertex);
 776
 777     Float_t clusPos[3]={0,0,0};
 778     cluster->GetPosition(clusPos);
 779     TVector3 clusterVector(clusPos[0],clusPos[1],clusPos[2]);
 780         Double_t etaCluster = clusterVector.Eta();
 781         Double_t phiCluster = clusterVector.Phi();
 782         if (phiCluster < 0) phiCluster= phiCluster+ 2*TMath::Pi();
 783
 784         // check eta range
 785         if (fUseEtaCut){
 786                 if (etaCluster < fMinEtaCut || etaCluster > fMaxEtaCut){
 787                         if(fHistAcceptanceCuts)fHistAcceptanceCuts->Fill(cutIndex);
 788                         return kFALSE;
 789                 }
 790         }
 791         cutIndex++;
 792
 793         // check phi range
 794         if (fUsePhiCut ){
 795         if (phiCluster < fMinPhiCut || phiCluster > fMaxPhiCut){
 796                         if(fHistAcceptanceCuts)fHistAcceptanceCuts->Fill(cutIndex);
 797                         return kFALSE;
 798                 }
 799         }
 800         cutIndex++;
 801
 802         // check distance to bad channel
 803         if (fUseDistanceToBadChannel){
 804                 if (cluster->GetDistanceToBadChannel() < fMinDistanceToBadChannel){
 805                         if(fHistAcceptanceCuts)fHistAcceptanceCuts->Fill(cutIndex);
 806                         return kFALSE;
 807                 }
 808         }
 809         cutIndex++;
 810         if(fHistAcceptanceCuts)fHistAcceptanceCuts->Fill(cutIndex);
 811
 812         // Histos after cuts
 813         if(fHistClusterEtavsPhiAfterAcc) fHistClusterEtavsPhiAfterAcc->Fill(phiCluster,etaCluster);
 814
 815         return kTRUE;
 816 }
 817
 818 Bool_t AliCaloPhotonCuts::MatchConvPhotonToCluster(AliAODConversionPhoton* convPhoton, AliVCluster* cluster, AliVEvent* event ){
 819
 820         if (!fUseDistTrackToCluster) return kFALSE;
 821
 822         AliESDEvent *esdev = dynamic_cast<AliESDEvent*>(event);
 823         AliAODEvent *aodev = 0;
 824         if (!esdev) {
 825                 aodev = dynamic_cast<AliAODEvent*>(event);
 826                 if (!aodev) {
 827                         AliError("Task needs AOD or ESD event, returning");
 828                         return kFALSE;
 829                 }
 830         }
 831
 832 //    Double_t vertex[3] = {0,0,0};
 833 //        event->GetPrimaryVertex()->GetXYZ(vertex);
 834
 835     if(!cluster->IsEMCAL() && !cluster->IsPHOS()){AliError("Cluster is neither EMCAL nor PHOS, returning"); return kFALSE;}
 836
 837     Float_t clusterPosition[3] = {0,0,0};
 838     cluster->GetPosition(clusterPosition);
 839     Double_t clusterR = TMath::Sqrt( clusterPosition[0]*clusterPosition[0] + clusterPosition[1]*clusterPosition[1] );
 840     if(fHistClusterRBeforeQA) fHistClusterRBeforeQA->Fill(clusterR);
 841
 842 //cout << "+++++++++ Cluster: x, y, z, R" << clusterPosition[0] << ", " << clusterPosition[1] << ", " << clusterPosition[2] << ", " << clusterR << "+++++++++" << endl;
 843
 844         Bool_t matched = kFALSE;
 845         for (Int_t i = 0; i < 2; i++){
 846                 Int_t tracklabel = convPhoton->GetLabel(i);
 847                 AliVTrack *inTrack = 0x0;
 848                 if(esdev) {
 849                         if(tracklabel > event->GetNumberOfTracks() ) continue;
 850                         inTrack = esdev->GetTrack(tracklabel);
 851                 } else {
 852                         if(((AliV0ReaderV1*)AliAnalysisManager::GetAnalysisManager()->GetTask("V0ReaderV1"))->AreAODsRelabeled()){
 853                                 inTrack = dynamic_cast<AliVTrack*>(event->GetTrack(tracklabel));
 854                         } else {
 855                                 for(Int_t ii=0; ii<event->GetNumberOfTracks(); ii++) {
 856                                         inTrack = dynamic_cast<AliVTrack*>(event->GetTrack(ii));
 857                                         if(inTrack){
 858                                                 if(inTrack->GetID() == tracklabel) {
 859                                                         continue;
 860                                                 }
 861                                         }
 862                                 }
 863                         }
 864                 }
 865
 866                 AliESDtrack *esdt = dynamic_cast<AliESDtrack*>(inTrack);
 867                 AliAODTrack *aodt = 0;
 868                 if (!esdt) {
 869                         aodt = dynamic_cast<AliAODTrack*>(inTrack);
 870                         if (!aodt){AliError("Track is neither ESD nor AOD, continue"); continue;}
 871                 }
 872
 873                 AliExternalTrackParam *trackParam = 0;
 874                 if (esdt) {
 875                         const AliExternalTrackParam *in = esdt->GetInnerParam();
 876                         if (!in){AliError("Could not get InnerParam of Track, continue"); continue;}
 877                         trackParam = new AliExternalTrackParam(*in);
 878                 } else {
 879                         Double_t xyz[3] = {0}, pxpypz[3] = {0}, cv[21] = {0};
 880                         aodt->PxPyPz(pxpypz);
 881                         aodt->XvYvZv(xyz);
 882                         aodt->GetCovarianceXYZPxPyPz(cv);
 883                         trackParam = new AliExternalTrackParam(xyz,pxpypz,cv,aodt->Charge());
 884                 }
 885                 if (!trackParam){AliError("Could not get TrackParameters, continue"); continue;}
 886
 887                 Bool_t propagated = kFALSE;
 888                 AliExternalTrackParam emcParam(*trackParam);
 889                 Float_t dPhi = 0;
 890                 Float_t dEta = 0;
 891
 892                 if(cluster->IsEMCAL()){
 893                         Float_t eta = 0; Float_t phi = 0; Float_t pt = 0;
 894                         propagated = AliEMCALRecoUtils::ExtrapolateTrackToEMCalSurface(&emcParam, 430, 0.000510999, 20, eta, phi, pt);
 895                         if(propagated){
 896                                 propagated = AliEMCALRecoUtils::ExtrapolateTrackToCluster(&emcParam, cluster, 0.000510999, 5, dEta, dPhi);
 897                         }
 898                 }
 899                 if(cluster->IsPHOS()){
 900                         propagated = AliTrackerBase::PropagateTrackToBxByBz(&emcParam, clusterR, 0.000510999, 20, kTRUE, 0.8, -1);
 901                         if (propagated){
 902                                 Double_t trkPos[3] = {0,0,0};
 903                                 emcParam.GetXYZ(trkPos);
 904                                 TVector3 trkPosVec(trkPos[0],trkPos[1],trkPos[2]);
 905                                 TVector3 clsPosVec(clusterPosition);
 906                                 dPhi = clsPosVec.DeltaPhi(trkPosVec);
 907                                 dEta = clsPosVec.Eta()-trkPosVec.Eta();
 908                         }
 909                 }
 910
 911                 if (propagated){
 912                         Float_t dR2 = dPhi*dPhi + dEta*dEta;
 913             if(fHistDistanceTrackToClusterBeforeQA)fHistDistanceTrackToClusterBeforeQA->Fill(TMath::Sqrt(dR2));
 914             if(fHistClusterdEtadPhiBeforeQA) fHistClusterdEtadPhiBeforeQA->Fill(dEta, dPhi);
 915
 916             Float_t clusM02 = (Float_t) cluster->GetM02();
 917             Float_t clusM20 = (Float_t) cluster->GetM20();
 918             if(fExtendedMatching){
 919                 if(inTrack->Charge() > 0) {
 920                     fHistClusterdEtadPhiPosTracksBeforeQA->Fill(dEta, dPhi);
 921                     if(inTrack->P() < 0.75) fHistClusterdEtadPhiPosTracksP_000_075BeforeQA->Fill(dEta, dPhi);
 922                     else if(inTrack->P() < 1.25) fHistClusterdEtadPhiPosTracksP_075_125BeforeQA->Fill(dEta, dPhi);
 923                     else fHistClusterdEtadPhiPosTracksP_125_999BeforeQA->Fill(dEta, dPhi);
 924                 }
 925                 else{
 926                     fHistClusterdEtadPhiNegTracksBeforeQA->Fill(dEta, dPhi);
 927                     if(inTrack->P() < 0.75) fHistClusterdEtadPhiNegTracksP_000_075BeforeQA->Fill(dEta, dPhi);
 928                     else if(inTrack->P() < 1.25) fHistClusterdEtadPhiNegTracksP_075_125BeforeQA->Fill(dEta, dPhi);
 929                     else fHistClusterdEtadPhiNegTracksP_125_999BeforeQA->Fill(dEta, dPhi);
 930                 }
 931                 fHistClusterdEtadPtBeforeQA->Fill(dEta, inTrack->Pt());
 932                 fHistClusterdPhidPtBeforeQA->Fill(dPhi, inTrack->Pt());
 933                 if(abs(dPhi) > 0.05){
 934                     fHistClusterM20Pt_dPhiBeforeQA->Fill(clusM20, inTrack->Pt());
 935                     fHistClusterM02Pt_dPhiBeforeQA->Fill(clusM02, inTrack->Pt());
 936                 }
 937                 fHistClusterM20M02BeforeQA->Fill(clusM20, clusM02);
 938             }
 939
 940             Bool_t match_dEta = (abs(dEta) < fMaxDistTrackToClusterEta) ? kTRUE : kFALSE;
 941             Bool_t match_dPhi = kFALSE;
 942             if( (inTrack->Charge() > 0) && (dPhi > fMinDistTrackToClusterPhi) && (dPhi < fMaxDistTrackToClusterPhi) ) match_dPhi = kTRUE;
 943             else if( (inTrack->Charge() < 0) && (dPhi < -fMinDistTrackToClusterPhi) && (dPhi > -fMaxDistTrackToClusterPhi) ) match_dPhi = kTRUE;
 944
 945             if(match_dEta && match_dPhi){
 946             //if(dR2 < fMinDistTrackToCluster*fMinDistTrackToCluster){
 947                                 matched = kTRUE;
 948                         } else {
 949                 if(fHistDistanceTrackToClusterAfterQA)fHistDistanceTrackToClusterAfterQA->Fill(TMath::Sqrt(dR2));
 950                 if(fHistClusterdEtadPhiAfterQA) fHistClusterdEtadPhiAfterQA->Fill(dEta, dPhi);
 951                 if(fHistClusterRAfterQA) fHistClusterRAfterQA->Fill(clusterR);
 952                 if(fExtendedMatching){
 953                     if(inTrack->Charge() > 0) fHistClusterdEtadPhiPosTracksAfterQA->Fill(dEta, dPhi);
 954                     else fHistClusterdEtadPhiNegTracksAfterQA->Fill(dEta, dPhi);
 955                     fHistClusterM20M02AfterQA->Fill(clusM20, clusM02);
 956                 }
 957                         }
 958                 }
 959                 delete trackParam;
 960         }
 961
 962         return matched;
 963
 964 }
 965
 966 //____________________________________________________________________________________________
 967
 968
 969 ///________________________________________________________________________
 970 Bool_t AliCaloPhotonCuts::UpdateCutString() {
 971    ///Update the cut string (if it has been created yet)
 972
 973    if(fCutString && fCutString->GetString().Length() == kNCuts) {
 974       fCutString->SetString(GetCutNumber());
 975    } else {
 976       return kFALSE;
 977    }
 978    return kTRUE;
 979 }
 980
 981 ///________________________________________________________________________
 982 Bool_t AliCaloPhotonCuts::InitializeCutsFromCutString(const TString analysisCutSelection ) {
 983         // Initialize Cuts from a given Cut string
 984         AliInfo(Form("Set CaloCut Number: %s",analysisCutSelection.Data()));
 985         if(analysisCutSelection.Length()!=kNCuts) {
 986                 AliError(Form("Cut selection has the wrong length! size is %d, number of cuts is %d", analysisCutSelection.Length(), kNCuts));
 987                 return kFALSE;
 988         }
 989         if(!analysisCutSelection.IsDigit()){
 990                 AliError("Cut selection contains characters");
 991                 return kFALSE;
 992         }
 993
 994         const char *cutSelection = analysisCutSelection.Data();
 995         #define ASSIGNARRAY(i)  fCuts[i] = cutSelection[i] - '0'
 996         for(Int_t ii=0;ii<kNCuts;ii++){
 997                 ASSIGNARRAY(ii);
 998         }
 999
1000         // Set Individual Cuts
1001         for(Int_t ii=0;ii<kNCuts;ii++){
1002                 if(!SetCut(cutIds(ii),fCuts[ii]))return kFALSE;
1003         }
1004         PrintCutsWithValues();
1005         return kTRUE;
1006 }
1007
1008 ///________________________________________________________________________
1009 Bool_t AliCaloPhotonCuts::SetCut(cutIds cutID, const Int_t value) {
1010         ///Set individual cut ID
1011
1012         switch (cutID) {
1013
1014                 case kClusterType:
1015                         if( SetClusterTypeCut(value)) {
1016                                 fCuts[kClusterType] = value;
1017                                 UpdateCutString();
1018                                 return kTRUE;
1019                         } else return kFALSE;
1020
1021                 case kEtaMin:
1022                         if( SetMinEtaCut(value)) {
1023                                 fCuts[kEtaMin] = value;
1024                                 UpdateCutString();
1025                                 return kTRUE;
1026                         } else return kFALSE;
1027
1028                 case kEtaMax:
1029                         if( SetMaxEtaCut(value)) {
1030                                 fCuts[kEtaMax] = value;
1031                                 UpdateCutString();
1032                                 return kTRUE;
1033                         } else return kFALSE;
1034
1035                 case kPhiMin:
1036                         if( SetMinPhiCut(value)) {
1037                                 fCuts[kPhiMin] = value;
1038                                 UpdateCutString();
1039                                 return kTRUE;
1040                         } else return kFALSE;
1041
1042                 case kPhiMax:
1043                         if( SetMaxPhiCut(value)) {
1044                                 fCuts[kPhiMax] = value;
1045                                 UpdateCutString();
1046                                 return kTRUE;
1047                         } else return kFALSE;
1048
1049                 case kDistanceToBadChannel:
1050                         if( SetDistanceToBadChannelCut(value)) {
1051                                 fCuts[kDistanceToBadChannel] = value;
1052                                 UpdateCutString();
1053                                 return kTRUE;
1054                         } else return kFALSE;
1055
1056                 case kTiming:
1057                         if( SetTimingCut(value)) {
1058                                 fCuts[kTiming] = value;
1059                                 UpdateCutString();
1060                                 return kTRUE;
1061                         } else return kFALSE;
1062
1063                 case kTrackMatching:
1064                         if( SetTrackMatchingCut(value)) {
1065                                 fCuts[kTrackMatching] = value;
1066                                 UpdateCutString();
1067                                 return kTRUE;
1068                         } else return kFALSE;
1069
1070                 case kExoticCell:
1071                         if( SetExoticCellCut(value)) {
1072                                 fCuts[kExoticCell] = value;
1073                                 UpdateCutString();
1074                                 return kTRUE;
1075                         } else return kFALSE;
1076
1077                 case kMinEnery:
1078                         if( SetMinEnergyCut(value)) {
1079                                 fCuts[kMinEnery] = value;
1080                                 UpdateCutString();
1081                                 return kTRUE;
1082                         } else return kFALSE;
1083
1084                 case kNMinCells:
1085                         if( SetMinNCellsCut(value)) {
1086                                 fCuts[kNMinCells] = value;
1087                                 UpdateCutString();
1088                                 return kTRUE;
1089                         } else return kFALSE;
1090
1091                 case kMinM02:
1092                         if( SetMinM02(value)) {
1093                                 fCuts[kMinM02] = value;
1094                                 UpdateCutString();
1095                                 return kTRUE;
1096                         } else return kFALSE;
1097
1098                 case kMaxM02:
1099                         if( SetMaxM02(value)) {
1100                                 fCuts[kMaxM02] = value;
1101                                 UpdateCutString();
1102                                 return kTRUE;
1103                         } else return kFALSE;
1104
1105                 case kMinM20:
1106                         if( SetMinM20(value)) {
1107                                 fCuts[kMinM20] = value;
1108                                 UpdateCutString();
1109                                 return kTRUE;
1110                         } else return kFALSE;
1111
1112                 case kMaxM20:
1113                         if( SetMaxM20(value)) {
1114                                 fCuts[kMaxM20] = value;
1115                                 UpdateCutString();
1116                                 return kTRUE;
1117                         } else return kFALSE;
1118
1119                 case kDispersion:
1120                         if( SetDispersion(value)) {
1121                                 fCuts[kDispersion] = value;
1122                                 UpdateCutString();
1123                                 return kTRUE;
1124                         } else return kFALSE;
1125
1126                 case kNLM:
1127                         if( SetNLM(value)) {
1128                                 fCuts[kNLM] = value;
1129                                 UpdateCutString();
1130                                 return kTRUE;
1131                         } else return kFALSE;
1132
1133                 case kNCuts:
1134                         AliError("Cut id out of range");
1135                         return kFALSE;
1136         }
1137
1138         AliError("Cut id %d not recognized");
1139         return kFALSE;
1140
1141
1142 }
1143 ///________________________________________________________________________
1144 void AliCaloPhotonCuts::PrintCuts() {
1145    // Print out current Cut Selection
1146    for(Int_t ic = 0; ic < kNCuts; ic++) {
1147       printf("%-30s : %d \n", fgkCutNames[ic], fCuts[ic]);
1148    }
1149 }
1150
1151 void AliCaloPhotonCuts::PrintCutsWithValues() {
1152         // Print out current Cut Selection with value
1153         printf("\nCluster cutnumber \n");
1154         for(Int_t ic = 0; ic < kNCuts; ic++) {
1155                 printf("%d",fCuts[ic]);
1156         }
1157         printf("\n\n");
1158
1159         printf("Acceptance cuts: \n");
1160         if (fClusterType == 0) printf("\tall calorimeter clusters are used\n");
1161         if (fClusterType == 1) printf("\tEMCAL calorimeter clusters are used\n");
1162         if (fClusterType == 2) printf("\tPHOS calorimeter clusters are used\n");
1163         if (fUseEtaCut) printf("\t%3.2f < eta_{cluster} < %3.2f\n", fMinEtaCut, fMaxEtaCut );
1164         if (fUsePhiCut) printf("\t%3.2f < phi_{cluster} < %3.2f\n", fMinPhiCut, fMaxPhiCut );
1165         if (fUseDistanceToBadChannel) printf("\tcut on exotics applied \n");
1166
1167         printf("Cluster Quality cuts: \n");
1168         if (fUseTimeDiff) printf("\t time difference < %3.2f\n", fMaxTimeDiff );
1169     if (fUseDistTrackToCluster) printf("\tmin distance to track in eta > %3.2f, min phi < %3.2f and max phi > %3.2f\n", fMaxDistTrackToClusterEta, fMinDistTrackToClusterPhi, fMaxDistTrackToClusterPhi );
1170     if (fUseExoticCell)printf("\t exotic cell: %3.2f\n", fExoticCell );
1171     if (fUseMinEnergy)printf("\t E_{cluster} > %3.2f\n", fMinEnergy );
1172         if (fUseNCells) printf("\t number of cells per cluster >= %d\n", fMinNCells );
1173         if (fUseM02) printf("\t %3.2f < M02 < %3.2f\n", fMinM02, fMaxM02 );
1174         if (fUseM20) printf("\t %3.2f < M20 < %3.2f\n", fMinM20, fMaxM20 );
1175         if (fUseDispersion) printf("\t dispersion < %3.2f\n", fMaxDispersion );
1176         if (fUseNLM) printf("\t %d < NLM < %d\n", fMinNLM, fMaxNLM );
1177
1178 }
1179
1180 // EMCAL acceptance 2011
1181 // 1.39626, 3.125 (phi)
1182 // -0.66687,,0.66465
1183
1184
1185 ///________________________________________________________________________
1186 Bool_t AliCaloPhotonCuts::SetClusterTypeCut(Int_t clusterType)
1187 {   // Set Cut
1188         switch(clusterType){
1189         case 0: // all clusters
1190                 fClusterType=0;
1191                 break;
1192         case 1: // EMCAL clusters
1193                 fClusterType=1;
1194                 break;
1195         case 2: // PHOS clusters
1196                 fClusterType=2;
1197                 break;
1198         default:
1199                 AliError(Form("ClusterTypeCut not defined %d",clusterType));
1200                 return kFALSE;
1201         }
1202         return kTRUE;
1203 }
1204
1205 //___________________________________________________________________
1206 Bool_t AliCaloPhotonCuts::SetMinEtaCut(Int_t minEta)
1207 {
1208         switch(minEta){
1209         case 0:
1210                 if (!fUseEtaCut) fUseEtaCut=0;
1211                 fMinEtaCut=-10.;
1212                 break;
1213         case 1:
1214                 if (!fUseEtaCut) fUseEtaCut=1;
1215                 fMinEtaCut=-0.6687;
1216                 break;
1217         case 2:
1218                 if (!fUseEtaCut) fUseEtaCut=1;
1219                 fMinEtaCut=-0.5;
1220                 break;
1221         case 3:
1222                 if (!fUseEtaCut) fUseEtaCut=1;
1223                 fMinEtaCut=-2;
1224                 break;
1225         default:
1226                 AliError(Form("MinEta Cut not defined %d",minEta));
1227                 return kFALSE;
1228         }
1229         return kTRUE;
1230 }
1231
1232
1233 //___________________________________________________________________
1234 Bool_t AliCaloPhotonCuts::SetMaxEtaCut(Int_t maxEta)
1235 {
1236         switch(maxEta){
1237         case 0:
1238                 if (!fUseEtaCut) fUseEtaCut=0;
1239                 fMaxEtaCut=10;
1240                 break;
1241         case 1:
1242                 if (!fUseEtaCut) fUseEtaCut=1;
1243                 fMaxEtaCut=0.66465;
1244                 break;
1245         case 2:
1246                 if (!fUseEtaCut) fUseEtaCut=1;
1247                 fMaxEtaCut=0.5;
1248                 break;
1249         case 3:
1250                 if (!fUseEtaCut) fUseEtaCut=1;
1251                 fMaxEtaCut=2;
1252                 break;
1253         default:
1254                 AliError(Form("MaxEta Cut not defined %d",maxEta));
1255                 return kFALSE;
1256         }
1257         return kTRUE;
1258 }
1259
1260 //___________________________________________________________________
1261 Bool_t AliCaloPhotonCuts::SetMinPhiCut(Int_t minPhi)
1262 {
1263         switch(minPhi){
1264         case 0:
1265                 if (!fUsePhiCut) fUsePhiCut=0;
1266                 fMinPhiCut=-10000;
1267                 break;
1268         case 1: // min EMCAL
1269                 if (!fUsePhiCut) fUsePhiCut=1;
1270                 fMinPhiCut=1.39626;
1271                 break;
1272         case 2: // min EMCAL with TRD 2012
1273                 if (!fUsePhiCut) fUsePhiCut=1;
1274                 fMinPhiCut=2.10;
1275                 break;
1276         case 3: // min EMCAL with TRD 2011
1277                 if (!fUsePhiCut) fUsePhiCut=1;
1278                 fMinPhiCut=2.45;
1279                 break;
1280
1281
1282         default:
1283                 AliError(Form("MinPhi Cut not defined %d",minPhi));
1284                 return kFALSE;
1285         }
1286         return kTRUE;
1287 }
1288
1289 //___________________________________________________________________
1290 Bool_t AliCaloPhotonCuts::SetMaxPhiCut(Int_t maxPhi)
1291 {
1292         switch(maxPhi){
1293         case 0:
1294                 if (!fUsePhiCut) fUsePhiCut=0;
1295                 fMaxPhiCut=10000;
1296                 break;
1297         case 1: // max EMCAL
1298                 if (!fUsePhiCut) fUsePhiCut=1;
1299                 fMaxPhiCut=3.15;
1300                 break;
1301         case 2: // max EMCAL with TRD 2011
1302                 if (!fUsePhiCut) fUsePhiCut=1;
1303                 fMaxPhiCut=2.45;
1304                 break;
1305         case 3: // max EMCAL with TRD 2012
1306                 if (!fUsePhiCut) fUsePhiCut=1;
1307                 fMaxPhiCut=2.10;
1308                 break;
1309
1310
1311         default:
1312                 AliError(Form("Max Phi Cut not defined %d",maxPhi));
1313                 return kFALSE;
1314         }
1315         return kTRUE;
1316 }
1317
1318 //___________________________________________________________________
1319 Bool_t AliCaloPhotonCuts::SetDistanceToBadChannelCut(Int_t distanceToBadChannel)
1320 {
1321         switch(distanceToBadChannel){
1322         case 0:
1323                 fUseDistanceToBadChannel=0;
1324                 fMinDistanceToBadChannel=0;
1325                 break;
1326         case 1:
1327                 if (!fUseDistanceToBadChannel) fUseDistanceToBadChannel=1;
1328                 fMinDistanceToBadChannel=5;
1329                 break;
1330         default:
1331                 AliError(Form("minimum distance to bad channel Cut not defined %d",distanceToBadChannel));
1332                 return kFALSE;
1333         }
1334         return kTRUE;
1335 }
1336
1337 //___________________________________________________________________
1338 Bool_t AliCaloPhotonCuts::SetTimingCut(Int_t timing)
1339 {
1340         switch(timing){
1341         case 0:
1342                 fUseTimeDiff=0;
1343                 fMaxTimeDiff=500;
1344                 break;
1345         case 1:
1346                 if (!fUseTimeDiff) fUseTimeDiff=1;
1347                 fMaxTimeDiff=10e-7; //1000ns
1348                 break;
1349         case 2:
1350                 if (!fUseTimeDiff) fUseTimeDiff=1;
1351                 fMaxTimeDiff=50e-8; //500ns
1352                 break;
1353         case 3:
1354                 if (!fUseTimeDiff) fUseTimeDiff=1;
1355                 fMaxTimeDiff=20e-8; //200ns
1356                 break;
1357         case 4:
1358                 if (!fUseTimeDiff) fUseTimeDiff=1;
1359                 fMaxTimeDiff=10e-8; //100ns
1360                 break;
1361         case 5:
1362                 if (!fUseTimeDiff) fUseTimeDiff=1;
1363                 fMaxTimeDiff=50e-9; //50ns
1364                 break;
1365
1366         default:
1367                 AliError(Form("Timing Cut not defined %d",timing));
1368                 return kFALSE;
1369         }
1370         return kTRUE;
1371 }
1372
1373 //___________________________________________________________________
1374 Bool_t AliCaloPhotonCuts::SetTrackMatchingCut(Int_t trackMatching)
1375 {
1376         switch(trackMatching){
1377         case 0:
1378         fUseDistTrackToCluster = 0;
1379         fMaxDistTrackToClusterEta = 0;
1380         fMinDistTrackToClusterPhi = 0;
1381         fMaxDistTrackToClusterPhi = 0;
1382                 break;
1383         case 1:
1384                 if (!fUseDistTrackToCluster) fUseDistTrackToCluster=1;
1385         fMaxDistTrackToClusterEta = 0.008;      //0.015;
1386         fMinDistTrackToClusterPhi = -0.03;      //-0.01;
1387         fMaxDistTrackToClusterPhi = 0.03;       //0.03; //0.04;
1388                 break;
1389         case 2:
1390                 if (!fUseDistTrackToCluster) fUseDistTrackToCluster=1;
1391         fMaxDistTrackToClusterEta = 0.012;      //0.015;
1392         fMinDistTrackToClusterPhi = -0.05;      //-0.01;
1393         fMaxDistTrackToClusterPhi = 0.04;       //0.035;        //0.05;
1394                 break;
1395         case 3:
1396                 if (!fUseDistTrackToCluster) fUseDistTrackToCluster=1;
1397         fMaxDistTrackToClusterEta = 0.016;      //0.015;
1398         fMinDistTrackToClusterPhi = -0.09;      //-0.015;
1399         fMaxDistTrackToClusterPhi = 0.06;       //0.04;         //0.1;
1400                 break;
1401         case 4:
1402                 if (!fUseDistTrackToCluster) fUseDistTrackToCluster=1;
1403         fMaxDistTrackToClusterEta = 0.018;      //0.015;
1404         fMinDistTrackToClusterPhi = -0.11;      //-0.015;
1405         fMaxDistTrackToClusterPhi = 0.07;       //0.045;        //0.13;
1406                 break;
1407         case 5:
1408                 if (!fUseDistTrackToCluster) fUseDistTrackToCluster=1;
1409         fMaxDistTrackToClusterEta = 0.02;       //0.015;
1410         fMinDistTrackToClusterPhi = -0.13;      //-0.02;
1411         fMaxDistTrackToClusterPhi = 0.08;       //0.05; //0.15
1412                 break;
1413         case 6:
1414                 if (!fUseDistTrackToCluster) fUseDistTrackToCluster=1;
1415         fMaxDistTrackToClusterEta = 0.022;      //0.015;
1416         fMinDistTrackToClusterPhi = -0.15;      //-0.02;
1417         fMaxDistTrackToClusterPhi = 0.10;       //0.055;        //0.2;
1418                 break;
1419     case 7: //PHOS
1420                 if (!fUseDistTrackToCluster) fUseDistTrackToCluster=1;
1421         fMaxDistTrackToClusterEta = 0.005;      //0.015;
1422         fMinDistTrackToClusterPhi = -0.03;      //-0.025;
1423         fMaxDistTrackToClusterPhi = 0.03;       //0.06;         //0.3;
1424                 break;
1425     case 8: //PHOS
1426                 if (!fUseDistTrackToCluster) fUseDistTrackToCluster=1;
1427         fMaxDistTrackToClusterEta = 0.01;       //0.015;
1428         fMinDistTrackToClusterPhi = -0.09;      //-0.025;
1429         fMaxDistTrackToClusterPhi = 0.07;       //0.07; //0.4;
1430                 break;
1431     case 9: //PHOS
1432                 if (!fUseDistTrackToCluster) fUseDistTrackToCluster=1;
1433         fMaxDistTrackToClusterEta = 0.015;      //0.02;
1434         fMinDistTrackToClusterPhi = -0.15;      //-0.03;
1435         fMaxDistTrackToClusterPhi = 0.11;       //0.1;  //0.5;
1436                 break;
1437
1438         default:
1439                 AliError(Form("Track Matching Cut not defined %d",trackMatching));
1440                 return kFALSE;
1441         }
1442         return kTRUE;
1443 }
1444
1445 //___________________________________________________________________
1446 Bool_t AliCaloPhotonCuts::SetExoticCellCut(Int_t exoticCell)
1447 {
1448         switch(exoticCell){
1449         case 0:
1450                 fUseExoticCell=0;
1451                 fExoticCell=0;
1452                 break;
1453         case 1:
1454                 if (!fUseExoticCell) fUseExoticCell=1;
1455                 fExoticCell=5;
1456                 break;
1457         default:
1458                 AliError(Form("Exotic cell Cut not defined %d",exoticCell));
1459                 return kFALSE;
1460         }
1461         return kTRUE;
1462 }
1463
1464 //___________________________________________________________________
1465 Bool_t AliCaloPhotonCuts::SetMinEnergyCut(Int_t minEnergy)
1466 {
1467         switch(minEnergy){
1468         case 0:
1469                 if (!fUseMinEnergy) fUseMinEnergy=0;
1470                 fMinEnergy=0.1;
1471                 break;
1472         case 1:
1473                 if (!fUseMinEnergy) fUseMinEnergy=1;
1474                 fMinEnergy=0.2;
1475                 break;
1476         case 2:
1477                 if (!fUseMinEnergy) fUseMinEnergy=1;
1478                 fMinEnergy=0.3;
1479                 break;
1480         case 3:
1481                 if (!fUseMinEnergy) fUseMinEnergy=1;
1482                 fMinEnergy=0.4;
1483                 break;
1484         case 4:
1485                 if (!fUseMinEnergy) fUseMinEnergy=1;
1486                 fMinEnergy=0.5;
1487                 break;
1488         case 5:
1489                 if (!fUseMinEnergy) fUseMinEnergy=1;
1490                 fMinEnergy=0.6;
1491                 break;
1492         case 6:
1493                 if (!fUseMinEnergy) fUseMinEnergy=1;
1494                 fMinEnergy=4.5;
1495                 break;
1496         case 7:
1497                 if (!fUseMinEnergy) fUseMinEnergy=1;
1498                 fMinEnergy=5.0;
1499                 break;
1500         case 8:
1501                 if (!fUseMinEnergy) fUseMinEnergy=1;
1502                 fMinEnergy=5.5;
1503                 break;
1504         case 9:
1505                 if (!fUseMinEnergy) fUseMinEnergy=1;
1506                 fMinEnergy=6.0;
1507                 break;
1508         default:
1509                 AliError(Form("Minimum Energy Cut not defined %d",minEnergy));
1510                 return kFALSE;
1511         }
1512         return kTRUE;
1513 }
1514
1515 //___________________________________________________________________
1516 Bool_t AliCaloPhotonCuts::SetMinNCellsCut(Int_t minNCells)
1517 {
1518         switch(minNCells){
1519         case 0:
1520                 if (!fUseNCells) fUseNCells=0;
1521                 fMinNCells=0;
1522                 break;
1523         case 1:
1524                 if (!fUseNCells) fUseNCells=1;
1525                 fMinNCells=1;
1526                 break;
1527         case 2:
1528                 if (!fUseNCells) fUseNCells=1;
1529                 fMinNCells=2;
1530                 break;
1531         case 3:
1532                 if (!fUseNCells) fUseNCells=1;
1533                 fMinNCells=3;
1534                 break;
1535         case 4:
1536                 if (!fUseNCells) fUseNCells=1;
1537                 fMinNCells=4;
1538                 break;
1539         case 5:
1540                 if (!fUseNCells) fUseNCells=1;
1541                 fMinNCells=5;
1542                 break;
1543         case 6:
1544                 if (!fUseNCells) fUseNCells=1;
1545                 fMinNCells=6;
1546                 break;
1547
1548         default:
1549                 AliError(Form("Min N cells Cut not defined %d",minNCells));
1550                 return kFALSE;
1551         }
1552         return kTRUE;
1553 }
1554
1555 //___________________________________________________________________
1556 Bool_t AliCaloPhotonCuts::SetMaxM02(Int_t maxM02)
1557 {
1558         switch(maxM02){
1559         case 0:
1560                 if (!fUseM02) fUseM02=0;
1561                 fMaxM02=100;
1562                 break;
1563         case 1:
1564                 if (!fUseM02) fUseM02=1;
1565                 fMaxM02=1.;
1566                 break;
1567         case 2:
1568                 if (!fUseM02) fUseM02=1;
1569                 fMaxM02=0.7;
1570                 break;
1571         case 3:
1572                 if (!fUseM02) fUseM02=1;
1573                 fMaxM02=0.5;
1574                 break;
1575         case 4:
1576                 if (!fUseM02) fUseM02=1;
1577                 fMaxM02=0.4;
1578                 break;
1579         default:
1580                 AliError(Form("Max M02 Cut not defined %d",maxM02));
1581                 return kFALSE;
1582         }
1583         return kTRUE;
1584 }
1585
1586 //___________________________________________________________________
1587 Bool_t AliCaloPhotonCuts::SetMinM02(Int_t minM02)
1588 {
1589         switch(minM02){
1590         case 0:
1591                 if (!fUseM02) fUseM02=0;
1592                 fMinM02=0;
1593                 break;
1594         case 1:
1595                 if (!fUseM02) fUseM02=1;
1596                 fMinM02=0.002;
1597                 break;
1598         case 2:
1599                 if (!fUseM02) fUseM02=1;
1600                 fMinM02=0.2;
1601                 break;
1602
1603         default:
1604                 AliError(Form("Min M02 not defined %d",minM02));
1605                 return kFALSE;
1606         }
1607         return kTRUE;
1608 }
1609
1610 //___________________________________________________________________
1611 Bool_t AliCaloPhotonCuts::SetMaxM20(Int_t maxM20)
1612 {
1613         switch(maxM20){
1614         case 0:
1615                 if (!fUseM20) fUseM20=0;
1616                 fMaxM20=100;
1617                 break;
1618         case 1:
1619                 if (!fUseM20) fUseM20=1;
1620                 fMaxM20=0.5;
1621                 break;
1622         default:
1623                 AliError(Form("Max M20 Cut not defined %d",maxM20));
1624                 return kFALSE;
1625         }
1626         return kTRUE;
1627 }
1628
1629 //___________________________________________________________________
1630 Bool_t AliCaloPhotonCuts::SetMinM20(Int_t minM20)
1631 {
1632         switch(minM20){
1633         case 0:
1634                 if (!fUseM20) fUseM20=0;
1635                 fMinM20=0;
1636                 break;
1637         case 1:
1638                 if (!fUseM20) fUseM20=1;
1639                 fMinM20=0.002;
1640                 break;
1641         default:
1642                 AliError(Form("Min M20 Cut not defined %d",minM20));
1643                 return kFALSE;
1644         }
1645         return kTRUE;
1646 }
1647
1648 //___________________________________________________________________
1649 Bool_t AliCaloPhotonCuts::SetDispersion(Int_t dispersion)
1650 {
1651         switch(dispersion){
1652         case 0:
1653                 if (!fUseDispersion) fUseDispersion=0;
1654                 fMaxDispersion =100;
1655                 break;
1656         case 1:
1657                 if (!fUseDispersion) fUseDispersion=1;
1658                 fMaxDispersion=2.;
1659                 break;
1660         default:
1661                 AliError(Form("Maximum Dispersion Cut not defined %d",dispersion));
1662                 return kFALSE;
1663         }
1664         return kTRUE;
1665 }
1666
1667 //___________________________________________________________________
1668 Bool_t AliCaloPhotonCuts::SetNLM(Int_t nlm)
1669 {
1670         switch(nlm){
1671         case 0:
1672                 if (!fUseNLM) fUseNLM=0;
1673                 fMinNLM =0;
1674                 fMaxNLM =100;
1675                 break;
1676         case 1:
1677                 if (!fUseNLM) fUseNLM=1;
1678                 fMinNLM =0;
1679                 fMaxNLM =1;
1680                 break;
1681         default:
1682                 AliError(Form("NLM Cut not defined %d",nlm));
1683                 return kFALSE;
1684         }
1685         return kTRUE;
1686 }
1687
1688 ///________________________________________________________________________
1689 TString AliCaloPhotonCuts::GetCutNumber(){
1690    // returns TString with current cut number
1691    TString a(kNCuts);
1692    for(Int_t ii=0;ii<kNCuts;ii++){
1693       a.Append(Form("%d",fCuts[ii]));
1694    }
1695    return a;
1696 }
1697
1698