AliPhysics  bba8f44 (bba8f44)
AliCaloPhotonCuts.cxx
Go to the documentation of this file.
1 /**************************************************************************
2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3  * *
4  * Authors: Friederike Bock, Daniel Muehlheim *
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 
17 //---------------------------------------------
18 // Class handling all kinds of selection cuts for
19 // Photon from EMCAL clusters
20 //---------------------------------------------
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 "AliMCEvent.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 "AliPHOSGeoUtils.h"
46 #include "AliTrackerBase.h"
47 #include "AliVCaloCells.h"
48 #include "AliVCluster.h"
49 #include "AliEmcalCorrectionTask.h"
50 #include "AliEmcalCorrectionComponent.h"
51 #include "AliTender.h"
52 #include "AliTenderSupply.h"
53 #include "AliEMCALTenderSupply.h"
54 #include "AliEmcalTenderTask.h"
55 #include "AliPHOSTenderSupply.h"
56 #include "AliPHOSTenderTask.h"
57 #include "AliOADBContainer.h"
58 #include "AliESDtrackCuts.h"
59 #include "AliCaloTrackMatcher.h"
60 #include <vector>
61 
62 class iostream;
63 
64 using namespace std;
65 
66 ClassImp(AliCaloPhotonCuts)
67 
68 
69 const char* AliCaloPhotonCuts::fgkCutNames[AliCaloPhotonCuts::kNCuts] = {
70  "ClusterType", //0 0: all, 1: EMCAL, 2: PHOS
71  "EtaMin", //1 0: -10, 1: -0.6687, 2: -0,5, 3: -2
72  "EtaMax", //2 0: 10, 1: 0.66465, 2: 0.5, 3: 2
73  "PhiMin", //3 0: -10000, 1: 1.39626
74  "PhiMax", //4 0: 10000, 1: 3.125
75  "NonLinearity1" //5
76  "NonLinearity2" //6
77  "DistanceToBadChannel", //7 0: 0, 1: 5
78  "Timing", //8 0: no cut
79  "TrackMatching", //9 0: 0, 1: 5
80  "ExoticCluster", //10 0: no cut
81  "MinEnergy", //11 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)
82  "MinNCells", //12 0: no cut, 1: 1, 2: 2, 3: 3, 4: 4, 5: 5, 6: 6
83  "MinM02", //13
84  "MaxM02", //14
85  "MinM20", //15
86  "MaxM20", //16
87  "MaximumDispersion", //17
88  "NLM" //18
89 };
90 
91 
92 //________________________________________________________________________
93 AliCaloPhotonCuts::AliCaloPhotonCuts(Int_t isMC, const char *name,const char *title) :
94  AliAnalysisCuts(name,title),
95  fHistograms(NULL),
96  fHistExtQA(NULL),
97  fCaloTrackMatcher(NULL),
98  fGeomEMCAL(NULL),
99  fEMCALRecUtils(NULL),
100  fEMCALInitialized(kFALSE),
101  fGeomPHOS(NULL),
102  fPHOSInitialized(kFALSE),
103  fPHOSCurrentRun(-1),
104  fEMCALBadChannelsMap(NULL),
105  fPHOSBadChannelsMap(NULL),
106  fBadChannels(NULL),
107  fNMaxEMCalModules(12),
108  fNMaxPHOSModules(5),
109  fHistoModifyAcc(NULL),
110  fDoLightOutput(kFALSE),
111  fIsMC(0),
112  fIsCurrentClusterAcceptedBeforeTM(kFALSE),
113  fV0ReaderName("V0ReaderV1"),
114  fCorrTaskSetting(""),
115  fCaloTrackMatcherName("CaloTrackMatcher_1"),
116  fPeriodName(""),
117  fCurrentMC(kNoMC),
118  fClusterType(0),
119  fMinEtaCut(-10),
120  fMinEtaInnerEdge(0),
121  fMaxEtaCut(10),
122  fMaxEtaInnerEdge(0),
123  fUseEtaCut(0),
124  fMinPhiCut(-10000),
125  fMaxPhiCut(-10000),
126  fUsePhiCut(0),
127  fMinDistanceToBadChannel(0),
128  fUseDistanceToBadChannel(0),
129  fMaxTimeDiff(10e10),
130  fMinTimeDiff(-10e10),
131  fUseTimeDiff(0),
132  fMaxDistTrackToClusterEta(0),
133  fMinDistTrackToClusterPhi(0),
134  fMaxDistTrackToClusterPhi(0),
135  fUseDistTrackToCluster(0),
136  fUsePtDepTrackToCluster(0),
137  fFuncPtDepEta(0),
138  fFuncPtDepPhi(0),
139  fExtendedMatchAndQA(0),
140  fExoticEnergyFracCluster(0),
141  fExoticMinEnergyCell(1),
142  fUseExoticCluster(0),
143  fDoExoticsQA(kFALSE),
144  fMinEnergy(0),
145  fSeedEnergy(0.1),
146  fLocMaxCutEDiff(0.03),
147  fUseMinEnergy(0),
148  fMinNCells(0),
149  fUseNCells(0),
150  fMaxM02(1000),
151  fMinM02(0),
152  fUseM02(0),
153  fMaxM02CutNr(0),
154  fMinM02CutNr(0),
155  fMaxM20(1000),
156  fMinM20(0),
157  fUseM20(0),
158  fMaxDispersion(1000),
159  fUseDispersion(0),
160  fMinNLM(0),
161  fMaxNLM(1000),
162  fUseNLM(0),
163  fNonLinearity1(0),
164  fNonLinearity2(0),
165  fSwitchNonLinearity(0),
166  fUseNonLinearity(kFALSE),
167  fIsPureCalo(0),
168  fVectorMatchedClusterIDs(0),
169  fCutString(NULL),
170  fCutStringRead(""),
171  fHistCutIndex(NULL),
172  fHistAcceptanceCuts(NULL),
173  fHistClusterIdentificationCuts(NULL),
174  fHistClusterEtavsPhiBeforeAcc(NULL),
175  fHistClusterEtavsPhiAfterAcc(NULL),
176  fHistClusterEtavsPhiAfterQA(NULL),
177  fHistClusterTimevsEBeforeQA(NULL),
178  fHistClusterTimevsEAfterQA(NULL),
179  fHistEnergyOfClusterBeforeNL(NULL),
180  fHistEnergyOfClusterAfterNL(NULL),
181  fHistEnergyOfClusterBeforeQA(NULL),
182  fHistEnergyOfClusterAfterQA(NULL),
183  fHistNCellsBeforeQA(NULL),
184  fHistNCellsAfterQA(NULL),
185  fHistM02BeforeQA(NULL),
186  fHistM02AfterQA(NULL),
187  fHistM20BeforeQA(NULL),
188  fHistM20AfterQA(NULL),
189  fHistDispersionBeforeQA(NULL),
190  fHistDispersionAfterQA(NULL),
191  fHistNLMBeforeQA(NULL),
192  fHistNLMAfterQA(NULL),
193  fHistNLMVsNCellsAfterQA(NULL),
194  fHistNLMVsEAfterQA(NULL),
195 // fHistNLMAvsNLMBBeforeQA(NULL),
196  fHistClusterEnergyvsMod(NULL),
197  fHistNCellsBigger100MeVvsMod(NULL),
198  fHistNCellsBigger1500MeVvsMod(NULL),
199  fHistEnergyOfModvsMod(NULL),
200  fHistClusterEnergyvsNCells(NULL),
201  fHistCellEnergyvsCellID(NULL),
202  fHistCellTimevsCellID(NULL),
203  fHistClusterEM02BeforeQA(NULL),
204  fHistClusterEM02AfterQA(NULL),
205  fHistClusterIncludedCellsBeforeQA(NULL),
206  fHistClusterIncludedCellsAfterQA(NULL),
207  fHistClusterEnergyFracCellsBeforeQA(NULL),
208  fHistClusterEnergyFracCellsAfterQA(NULL),
209  fHistClusterIncludedCellsTimingAfterQA(NULL),
210  fHistClusterIncludedCellsTimingEnergyAfterQA(NULL),
211  fHistClusterDistanceInTimeCut(NULL),
212  fHistClusterDistanceOutTimeCut(NULL),
213  fHistClusterDistance1DInTimeCut(NULL),
214  fHistClusterRBeforeQA(NULL),
215  fHistClusterRAfterQA(NULL),
216  fHistClusterdEtadPhiBeforeQA(NULL),
217  fHistClusterdEtadPhiAfterQA(NULL),
218  fHistDistanceTrackToClusterBeforeQA(NULL),
219  fHistDistanceTrackToClusterAfterQA(NULL),
220  fHistClusterdEtadPhiPosTracksBeforeQA(NULL),
221  fHistClusterdEtadPhiNegTracksBeforeQA(NULL),
222  fHistClusterdEtadPhiPosTracksAfterQA(NULL),
223  fHistClusterdEtadPhiNegTracksAfterQA(NULL),
224  fHistClusterdEtadPhiPosTracksP_000_075BeforeQA(NULL),
225  fHistClusterdEtadPhiPosTracksP_075_125BeforeQA(NULL),
226  fHistClusterdEtadPhiPosTracksP_125_999BeforeQA(NULL),
227  fHistClusterdEtadPhiNegTracksP_000_075BeforeQA(NULL),
228  fHistClusterdEtadPhiNegTracksP_075_125BeforeQA(NULL),
229  fHistClusterdEtadPhiNegTracksP_125_999BeforeQA(NULL),
230  fHistClusterdEtadPtBeforeQA(NULL),
231  fHistClusterdEtadPtAfterQA(NULL),
232  fHistClusterdEtadPtTrueMatched(NULL),
233  fHistClusterdPhidPtPosTracksBeforeQA(NULL),
234  fHistClusterdPhidPtNegTracksBeforeQA(NULL),
235  fHistClusterdPhidPtAfterQA(NULL),
236  fHistClusterdPhidPtPosTracksTrueMatched(NULL),
237  fHistClusterdPhidPtNegTracksTrueMatched(NULL),
238  fHistClusterM20M02BeforeQA(NULL),
239  fHistClusterM20M02AfterQA(NULL),
240  fHistClusterEtavsPhiExotics(NULL),
241  fHistClusterEM02Exotics(NULL),
242  fHistClusterEnergyvsNCellsExotics(NULL),
243  fHistClusterEEstarExotics(NULL),
244  fHistClusterTMEffiInput(NULL),
245  fHistClusterElecEtaPhiBeforeTM_00_20(NULL),
246  fHistClusterElecEtaPhiBeforeTM_20_50(NULL),
247  fHistClusterElecEtaPhiBeforeTM_50_00(NULL),
248  fHistClusterElecEtaPhiAfterTM_00_20(NULL),
249  fHistClusterElecEtaPhiAfterTM_20_50(NULL),
250  fHistClusterElecEtaPhiAfterTM_50_00(NULL),
251  fHistClusterEvsTrackECharged(NULL),
252  fHistClusterEvsTrackEChargedLead(NULL),
253  fHistClusterEvsTrackENeutral(NULL),
254  fHistClusterEvsTrackENeutralSubCharged(NULL),
255  fHistClusterEvsTrackEGamma(NULL),
256  fHistClusterEvsTrackEGammaSubCharged(NULL),
257  fHistClusterEvsTrackEConv(NULL),
258  fHistClusterENMatchesNeutral(NULL),
259  fHistClusterENMatchesCharged(NULL),
260  fHistClusterEvsTrackEPrimaryButNoElec(NULL),
261  fHistClusterEvsTrackSumEPrimaryButNoElec(NULL),
262  fNMaxDCalModules(8),
263  fgkDCALCols(32),
264  fIsAcceptedForBasic(kFALSE)
265 {
266  for(Int_t jj=0;jj<kNCuts;jj++){fCuts[jj]=0;}
267  fCutString=new TObjString((GetCutNumber()).Data());
268 
269  fIsMC = isMC;
270 }
271 
272 //________________________________________________________________________
273 AliCaloPhotonCuts::AliCaloPhotonCuts(const AliCaloPhotonCuts &ref) :
274  AliAnalysisCuts(ref),
275  fHistograms(NULL),
276  fHistExtQA(NULL),
277  fCaloTrackMatcher(NULL),
278  fGeomEMCAL(NULL),
279  fEMCALRecUtils(NULL),
280  fEMCALInitialized(kFALSE),
281  fGeomPHOS(NULL),
282  fPHOSInitialized(kFALSE),
283  fPHOSCurrentRun(-1),
284  fEMCALBadChannelsMap(NULL),
285  fPHOSBadChannelsMap(NULL),
286  fBadChannels(NULL),
287  fNMaxEMCalModules(ref.fNMaxEMCalModules),
288  fNMaxPHOSModules(ref.fNMaxPHOSModules),
289  fHistoModifyAcc(NULL),
290  fDoLightOutput(ref.fDoLightOutput),
291  fIsMC(ref.fIsMC),
292  fIsCurrentClusterAcceptedBeforeTM(kFALSE),
293  fV0ReaderName(ref.fV0ReaderName),
294  fCorrTaskSetting(ref.fCorrTaskSetting),
295  fCaloTrackMatcherName(ref.fCaloTrackMatcherName),
296  fPeriodName(ref.fPeriodName),
297  fCurrentMC(ref.fCurrentMC),
298  fClusterType(ref.fClusterType),
299  fMinEtaCut(ref.fMinEtaCut),
300  fMinEtaInnerEdge(ref.fMinEtaInnerEdge),
301  fMaxEtaCut(ref.fMaxEtaCut),
302  fMaxEtaInnerEdge(ref.fMaxEtaInnerEdge),
303  fUseEtaCut(ref.fUseEtaCut),
304  fMinPhiCut(ref.fMinPhiCut),
305  fMaxPhiCut(ref.fMaxPhiCut),
306  fUsePhiCut(ref.fUsePhiCut),
307  fMinDistanceToBadChannel(ref.fMinDistanceToBadChannel),
308  fUseDistanceToBadChannel(ref.fUseDistanceToBadChannel),
309  fMaxTimeDiff(ref.fMaxTimeDiff),
310  fMinTimeDiff(ref.fMinTimeDiff),
311  fUseTimeDiff(ref.fUseTimeDiff),
312  fMaxDistTrackToClusterEta(ref.fMaxDistTrackToClusterEta),
313  fMinDistTrackToClusterPhi(ref.fMinDistTrackToClusterPhi),
314  fMaxDistTrackToClusterPhi(ref.fMaxDistTrackToClusterPhi),
315  fUseDistTrackToCluster(ref.fUseDistTrackToCluster),
316  fUsePtDepTrackToCluster(ref.fUsePtDepTrackToCluster),
317  fFuncPtDepEta(ref.fFuncPtDepEta),
318  fFuncPtDepPhi(ref.fFuncPtDepPhi),
319  fExtendedMatchAndQA(ref.fExtendedMatchAndQA),
320  fExoticEnergyFracCluster(ref.fExoticEnergyFracCluster),
321  fExoticMinEnergyCell(ref.fExoticMinEnergyCell),
322  fUseExoticCluster(ref.fUseExoticCluster),
323  fDoExoticsQA(ref.fDoExoticsQA),
324  fMinEnergy(ref.fMinEnergy),
325  fSeedEnergy(ref.fSeedEnergy),
326  fLocMaxCutEDiff(ref.fLocMaxCutEDiff),
327  fUseMinEnergy(ref.fUseMinEnergy),
328  fMinNCells(ref.fMinNCells),
329  fUseNCells(ref.fUseNCells),
330  fMaxM02(ref.fMaxM02),
331  fMinM02(ref.fMinM02),
332  fUseM02(ref.fUseM02),
333  fMaxM02CutNr(ref.fMaxM02CutNr),
334  fMinM02CutNr(ref.fMinM02CutNr),
335  fMaxM20(ref.fMaxM20),
336  fMinM20(ref.fMinM20),
337  fUseM20(ref.fUseDispersion),
338  fMaxDispersion(ref.fMaxDispersion),
339  fUseDispersion(ref.fUseDispersion),
340  fMinNLM(ref.fMinNLM),
341  fMaxNLM(ref.fMaxNLM),
342  fUseNLM(ref.fUseNLM),
343  fNonLinearity1(ref.fNonLinearity1),
344  fNonLinearity2(ref.fNonLinearity2),
345  fSwitchNonLinearity(ref.fSwitchNonLinearity),
346  fUseNonLinearity(ref.fUseNonLinearity),
347  fIsPureCalo(ref.fIsPureCalo),
348  fVectorMatchedClusterIDs(0),
349  fCutString(NULL),
350  fCutStringRead(""),
351  fHistCutIndex(NULL),
352  fHistAcceptanceCuts(NULL),
353  fHistClusterIdentificationCuts(NULL),
354  fHistClusterEtavsPhiBeforeAcc(NULL),
355  fHistClusterEtavsPhiAfterAcc(NULL),
356  fHistClusterEtavsPhiAfterQA(NULL),
357  fHistClusterTimevsEBeforeQA(NULL),
358  fHistClusterTimevsEAfterQA(NULL),
359  fHistEnergyOfClusterBeforeNL(NULL),
360  fHistEnergyOfClusterAfterNL(NULL),
361  fHistEnergyOfClusterBeforeQA(NULL),
362  fHistEnergyOfClusterAfterQA(NULL),
363  fHistNCellsBeforeQA(NULL),
364  fHistNCellsAfterQA(NULL),
365  fHistM02BeforeQA(NULL),
366  fHistM02AfterQA(NULL),
367  fHistM20BeforeQA(NULL),
368  fHistM20AfterQA(NULL),
369  fHistDispersionBeforeQA(NULL),
370  fHistDispersionAfterQA(NULL),
371  fHistNLMBeforeQA(NULL),
372  fHistNLMAfterQA(NULL),
373 // fHistNLMAvsNLMBBeforeQA(NULL),
374  fHistNLMVsNCellsAfterQA(NULL),
375  fHistNLMVsEAfterQA(NULL),
376  fHistClusterEnergyvsMod(NULL),
377  fHistNCellsBigger100MeVvsMod(NULL),
378  fHistNCellsBigger1500MeVvsMod(NULL),
379  fHistEnergyOfModvsMod(NULL),
380  fHistClusterEnergyvsNCells(NULL),
381  fHistCellEnergyvsCellID(NULL),
382  fHistCellTimevsCellID(NULL),
383  fHistClusterEM02BeforeQA(NULL),
384  fHistClusterEM02AfterQA(NULL),
385  fHistClusterIncludedCellsBeforeQA(NULL),
386  fHistClusterIncludedCellsAfterQA(NULL),
387  fHistClusterEnergyFracCellsBeforeQA(NULL),
388  fHistClusterEnergyFracCellsAfterQA(NULL),
389  fHistClusterIncludedCellsTimingAfterQA(NULL),
390  fHistClusterIncludedCellsTimingEnergyAfterQA(NULL),
391  fHistClusterDistanceInTimeCut(NULL),
392  fHistClusterDistanceOutTimeCut(NULL),
393  fHistClusterDistance1DInTimeCut(NULL),
394  fHistClusterRBeforeQA(NULL),
395  fHistClusterRAfterQA(NULL),
396  fHistClusterdEtadPhiBeforeQA(NULL),
397  fHistClusterdEtadPhiAfterQA(NULL),
398  fHistDistanceTrackToClusterBeforeQA(NULL),
399  fHistDistanceTrackToClusterAfterQA(NULL),
400  fHistClusterdEtadPhiPosTracksBeforeQA(NULL),
401  fHistClusterdEtadPhiNegTracksBeforeQA(NULL),
402  fHistClusterdEtadPhiPosTracksAfterQA(NULL),
403  fHistClusterdEtadPhiNegTracksAfterQA(NULL),
404  fHistClusterdEtadPhiPosTracksP_000_075BeforeQA(NULL),
405  fHistClusterdEtadPhiPosTracksP_075_125BeforeQA(NULL),
406  fHistClusterdEtadPhiPosTracksP_125_999BeforeQA(NULL),
407  fHistClusterdEtadPhiNegTracksP_000_075BeforeQA(NULL),
408  fHistClusterdEtadPhiNegTracksP_075_125BeforeQA(NULL),
409  fHistClusterdEtadPhiNegTracksP_125_999BeforeQA(NULL),
410  fHistClusterdEtadPtBeforeQA(NULL),
411  fHistClusterdEtadPtAfterQA(NULL),
412  fHistClusterdEtadPtTrueMatched(NULL),
413  fHistClusterdPhidPtPosTracksBeforeQA(NULL),
414  fHistClusterdPhidPtNegTracksBeforeQA(NULL),
415  fHistClusterdPhidPtAfterQA(NULL),
416  fHistClusterdPhidPtPosTracksTrueMatched(NULL),
417  fHistClusterdPhidPtNegTracksTrueMatched(NULL),
418  fHistClusterM20M02BeforeQA(NULL),
419  fHistClusterM20M02AfterQA(NULL),
420  fHistClusterEtavsPhiExotics(NULL),
421  fHistClusterEM02Exotics(NULL),
422  fHistClusterEnergyvsNCellsExotics(NULL),
423  fHistClusterEEstarExotics(NULL),
424  fHistClusterTMEffiInput(NULL),
425  fHistClusterElecEtaPhiBeforeTM_00_20(NULL),
426  fHistClusterElecEtaPhiBeforeTM_20_50(NULL),
427  fHistClusterElecEtaPhiBeforeTM_50_00(NULL),
428  fHistClusterElecEtaPhiAfterTM_00_20(NULL),
429  fHistClusterElecEtaPhiAfterTM_20_50(NULL),
430  fHistClusterElecEtaPhiAfterTM_50_00(NULL),
431  fHistClusterEvsTrackECharged(NULL),
432  fHistClusterEvsTrackEChargedLead(NULL),
433  fHistClusterEvsTrackENeutral(NULL),
434  fHistClusterEvsTrackENeutralSubCharged(NULL),
435  fHistClusterEvsTrackEGamma(NULL),
436  fHistClusterEvsTrackEGammaSubCharged(NULL),
437  fHistClusterEvsTrackEConv(NULL),
438  fHistClusterENMatchesNeutral(NULL),
439  fHistClusterENMatchesCharged(NULL),
440  fHistClusterEvsTrackEPrimaryButNoElec(NULL),
441  fHistClusterEvsTrackSumEPrimaryButNoElec(NULL),
442  fNMaxDCalModules(ref.fNMaxDCalModules),
443  fgkDCALCols(ref.fgkDCALCols),
444  fIsAcceptedForBasic(ref.fIsAcceptedForBasic)
445 {
446  // Copy Constructor
447  for(Int_t jj=0;jj<kNCuts;jj++){fCuts[jj]=ref.fCuts[jj];}
448  fCutString=new TObjString((GetCutNumber()).Data());
449 
450 }
451 
452 
453 //________________________________________________________________________
454 AliCaloPhotonCuts::~AliCaloPhotonCuts() {
455  // Destructor
456  if(fCutString != NULL){
457  delete fCutString;
458  fCutString = NULL;
459  }
460 
461  if(fPHOSBadChannelsMap != NULL){
462  delete[] fPHOSBadChannelsMap;
463  fPHOSBadChannelsMap = NULL;
464  }
465 
466  if(fFuncPtDepEta) delete fFuncPtDepEta;
467  if(fFuncPtDepPhi) delete fFuncPtDepPhi;
468 }
469 
470 //________________________________________________________________________
471 void AliCaloPhotonCuts::InitCutHistograms(TString name){
472 
473  // Initialize Cut Histograms for QA (only initialized and filled if function is called)
474  TH1::AddDirectory(kFALSE);
475 
476 
477 
478  if (fDoLightOutput)
479  fDoExoticsQA = kFALSE;
480 
481  if(fHistograms != NULL){
482  delete fHistograms;
483  fHistograms=NULL;
484  }
485  if(fHistograms==NULL){
486  fHistograms = new TList();
487  fHistograms->SetOwner(kTRUE);
488  if(name=="")fHistograms->SetName(Form("CaloCuts_%s",GetCutNumber().Data()));
489  else fHistograms->SetName(Form("%s_%s",name.Data(),GetCutNumber().Data()));
490  }
491 
492  if(fHistExtQA != NULL){
493  delete fHistExtQA;
494  fHistExtQA=NULL;
495  }
496  if(fHistExtQA==NULL){
497  fHistExtQA = new TList();
498  fHistExtQA->SetOwner(kTRUE);
499  if(name=="")fHistExtQA->SetName(Form("CaloExtQA_%s",GetCutNumber().Data()));
500  else fHistExtQA->SetName(Form("%s_ExtQA_%s",name.Data(),GetCutNumber().Data()));
501  }
502 
503  Int_t nBinsClusterEFine = 400;
504  Int_t nBinsClusterECoarse = 100;
505  Double_t minClusterELog = 0.5;
506  Double_t maxClusterELog = 100.0;
507  Int_t nBinsCellECoarse = 150;
508  Double_t minCellELog = 0.05;
509  Double_t maxCellELog = 120.0;
510  Int_t nBinsModuleECoarse = 400;
511  Double_t minModuleELog = 0.1;
512  Double_t maxModuleELog = 400.0;
513 
514  // IsPhotonSelected
515  fHistCutIndex = new TH1F(Form("IsPhotonSelected %s",GetCutNumber().Data()),"IsPhotonSelected",5,-0.5,4.5);
516  fHistCutIndex->GetXaxis()->SetBinLabel(kPhotonIn+1,"in");
517  fHistCutIndex->GetXaxis()->SetBinLabel(kDetector+1,"detector");
518  fHistCutIndex->GetXaxis()->SetBinLabel(kAcceptance+1,"acceptance");
519  fHistCutIndex->GetXaxis()->SetBinLabel(kClusterQuality+1,"cluster QA");
520  fHistCutIndex->GetXaxis()->SetBinLabel(kPhotonOut+1,"out");
521  fHistograms->Add(fHistCutIndex);
522 
523  // Acceptance Cuts
524  fHistAcceptanceCuts = new TH1F(Form("AcceptanceCuts %s",GetCutNumber().Data()),"AcceptanceCuts",5,-0.5,4.5);
525  fHistAcceptanceCuts->GetXaxis()->SetBinLabel(1,"in");
526  fHistAcceptanceCuts->GetXaxis()->SetBinLabel(2,"eta");
527  fHistAcceptanceCuts->GetXaxis()->SetBinLabel(3,"phi");
528  fHistAcceptanceCuts->GetXaxis()->SetBinLabel(4,"dist bad channel/acc map");
529  fHistAcceptanceCuts->GetXaxis()->SetBinLabel(5,"out");
530  fHistograms->Add(fHistAcceptanceCuts);
531 
532  // Cluster Cuts
533  fHistClusterIdentificationCuts = new TH2F(Form("ClusterQualityCuts vs E %s",GetCutNumber().Data()),"ClusterQualityCuts",11,-0.5,10.5,nBinsClusterEFine, minClusterELog, maxClusterELog);
534  fHistClusterIdentificationCuts->GetXaxis()->SetBinLabel(1,"in");
535  fHistClusterIdentificationCuts->GetXaxis()->SetBinLabel(2,"timing");
536  fHistClusterIdentificationCuts->GetXaxis()->SetBinLabel(3,"Exotics");
537  fHistClusterIdentificationCuts->GetXaxis()->SetBinLabel(4,"minimum NCells");
538  fHistClusterIdentificationCuts->GetXaxis()->SetBinLabel(5,"NLM");
539  fHistClusterIdentificationCuts->GetXaxis()->SetBinLabel(6,"M02");
540  fHistClusterIdentificationCuts->GetXaxis()->SetBinLabel(7,"M20");
541  fHistClusterIdentificationCuts->GetXaxis()->SetBinLabel(8,"dispersion");
542  fHistClusterIdentificationCuts->GetXaxis()->SetBinLabel(9,"track matching");
543  fHistClusterIdentificationCuts->GetXaxis()->SetBinLabel(10,"minimum energy");
544  fHistClusterIdentificationCuts->GetXaxis()->SetBinLabel(11,"out");
545  SetLogBinningYTH2(fHistClusterIdentificationCuts);
546  fHistograms->Add(fHistClusterIdentificationCuts);
547 
548  // Acceptance related histogramms
549  if( fClusterType == 1 ){ //EMCAL
550  const Int_t nEmcalEtaBins = 96;
551  const Int_t nEmcalPhiBins = 124;
552  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};
553  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};
554 
555  if(!fDoLightOutput){
556  fHistClusterEtavsPhiBeforeAcc = new TH2F(Form("EtaPhi_beforeAcceptance %s",GetCutNumber().Data()),"EtaPhi_beforeAcceptance",nEmcalPhiBins,EmcalPhiBins,nEmcalEtaBins,EmcalEtaBins);
557  fHistClusterEtavsPhiBeforeAcc->GetXaxis()->SetTitle("#varphi (rad)");
558  fHistClusterEtavsPhiBeforeAcc->GetYaxis()->SetTitle("#eta");
559  fHistograms->Add(fHistClusterEtavsPhiBeforeAcc);
560  fHistClusterEtavsPhiAfterAcc = new TH2F(Form("EtaPhi_afterAcceptance %s",GetCutNumber().Data()),"EtaPhi_afterAcceptance",nEmcalPhiBins,EmcalPhiBins,nEmcalEtaBins,EmcalEtaBins);
561  fHistClusterEtavsPhiAfterAcc->GetXaxis()->SetTitle("#varphi (rad)");
562  fHistClusterEtavsPhiAfterAcc->GetYaxis()->SetTitle("#eta");
563  fHistograms->Add(fHistClusterEtavsPhiAfterAcc);
564  }
565  fHistClusterEtavsPhiAfterQA = new TH2F(Form("EtaPhi_afterClusterQA %s",GetCutNumber().Data()),"EtaPhi_afterClusterQA",nEmcalPhiBins,EmcalPhiBins,nEmcalEtaBins,EmcalEtaBins);
566  fHistClusterEtavsPhiAfterQA->GetXaxis()->SetTitle("#varphi (rad)");
567  fHistClusterEtavsPhiAfterQA->GetYaxis()->SetTitle("#eta");
568  fHistograms->Add(fHistClusterEtavsPhiAfterQA);
569  } else if( fClusterType == 2 ){ //PHOS
570  const Int_t nPhosEtaBins = 56;
571  const Int_t nPhosPhiBins = 256;
572  const Float_t PhosEtaRange[2] = {-0.16, 0.16};
573  const Float_t PhosPhiRange[2] = {4.355, 5.6};
574 
575  if(!fDoLightOutput){
576  fHistClusterEtavsPhiBeforeAcc = new TH2F(Form("EtaPhi_beforeAcceptance %s",GetCutNumber().Data()),"EtaPhi_beforeAcceptance",nPhosPhiBins,PhosPhiRange[0],PhosPhiRange[1],nPhosEtaBins,PhosEtaRange[0],PhosEtaRange[1]);
577  fHistClusterEtavsPhiBeforeAcc->GetXaxis()->SetTitle("#varphi (rad)");
578  fHistClusterEtavsPhiBeforeAcc->GetYaxis()->SetTitle("#eta");
579  fHistograms->Add(fHistClusterEtavsPhiBeforeAcc);
580  fHistClusterEtavsPhiAfterAcc = new TH2F(Form("EtaPhi_afterAcceptance %s",GetCutNumber().Data()),"EtaPhi_afterAcceptance",nPhosPhiBins,PhosPhiRange[0],PhosPhiRange[1],nPhosEtaBins,PhosEtaRange[0],PhosEtaRange[1]);
581  fHistClusterEtavsPhiAfterAcc->GetXaxis()->SetTitle("#varphi (rad)");
582  fHistClusterEtavsPhiAfterAcc->GetYaxis()->SetTitle("#eta");
583  fHistograms->Add(fHistClusterEtavsPhiAfterAcc);
584  }
585  fHistClusterEtavsPhiAfterQA = new TH2F(Form("EtaPhi_afterClusterQA %s",GetCutNumber().Data()),"EtaPhi_afterClusterQA",nPhosPhiBins,PhosPhiRange[0],PhosPhiRange[1],nPhosEtaBins,PhosEtaRange[0],PhosEtaRange[1]);
586  fHistClusterEtavsPhiAfterQA->GetXaxis()->SetTitle("#varphi (rad)");
587  fHistClusterEtavsPhiAfterQA->GetYaxis()->SetTitle("#eta");
588  fHistograms->Add(fHistClusterEtavsPhiAfterQA);
589  } else if( fClusterType == 3){ //DCAL
590  const Int_t nDcalEtaBins = 96;
591  const Int_t nDcalPhiBins = 124;
592  if(!fDoLightOutput){
593  fHistClusterEtavsPhiBeforeAcc = new TH2F(Form("EtaPhi_beforeAcceptance %s",GetCutNumber().Data()),"EtaPhi_beforeAcceptance",nDcalPhiBins,4.5,5.7,nDcalEtaBins,-0.66687,0.66465);
594  fHistClusterEtavsPhiBeforeAcc->GetXaxis()->SetTitle("#varphi (rad)");
595  fHistClusterEtavsPhiBeforeAcc->GetYaxis()->SetTitle("#eta");
596  fHistograms->Add(fHistClusterEtavsPhiBeforeAcc);
597  fHistClusterEtavsPhiAfterAcc = new TH2F(Form("EtaPhi_afterAcceptance %s",GetCutNumber().Data()),"EtaPhi_afterAcceptance",nDcalPhiBins,4.5,5.7,nDcalEtaBins,-0.66687,0.66465);
598  fHistClusterEtavsPhiAfterAcc->GetXaxis()->SetTitle("#varphi (rad)");
599  fHistClusterEtavsPhiAfterAcc->GetYaxis()->SetTitle("#eta");
600  fHistograms->Add(fHistClusterEtavsPhiAfterAcc);
601  }
602  fHistClusterEtavsPhiAfterQA = new TH2F(Form("EtaPhi_afterClusterQA %s",GetCutNumber().Data()),"EtaPhi_afterClusterQA",nDcalPhiBins,4.5,5.7,nDcalEtaBins,-0.66687,0.66465);
603  fHistClusterEtavsPhiAfterQA->GetXaxis()->SetTitle("#varphi (rad)");
604  fHistClusterEtavsPhiAfterQA->GetYaxis()->SetTitle("#eta");
605  fHistograms->Add(fHistClusterEtavsPhiAfterQA);
606  } else if( fClusterType == 0 ){ //all
607  if(!fDoLightOutput){
608  fHistClusterEtavsPhiBeforeAcc = new TH2F(Form("EtaPhi_beforeAcceptance %s",GetCutNumber().Data()),"EtaPhi_beforeAcceptance",462,0,2*TMath::Pi(),110,-0.7,0.7);
609  fHistClusterEtavsPhiBeforeAcc->GetXaxis()->SetTitle("#varphi (rad)");
610  fHistClusterEtavsPhiBeforeAcc->GetYaxis()->SetTitle("#eta");
611  fHistograms->Add(fHistClusterEtavsPhiBeforeAcc);
612  fHistClusterEtavsPhiAfterAcc = new TH2F(Form("EtaPhi_afterAcceptance %s",GetCutNumber().Data()),"EtaPhi_afterAcceptance",462,0,2*TMath::Pi(),110,-0.7,0.7);
613  fHistClusterEtavsPhiAfterAcc->GetXaxis()->SetTitle("#varphi (rad)");
614  fHistClusterEtavsPhiAfterAcc->GetYaxis()->SetTitle("#eta");
615  fHistograms->Add(fHistClusterEtavsPhiAfterAcc);
616  }
617  fHistClusterEtavsPhiAfterQA = new TH2F(Form("EtaPhi_afterClusterQA %s",GetCutNumber().Data()),"EtaPhi_afterClusterQA",462,0,2*TMath::Pi(),110,-0.7,0.7);
618  fHistClusterEtavsPhiAfterQA->GetXaxis()->SetTitle("#varphi (rad)");
619  fHistClusterEtavsPhiAfterQA->GetYaxis()->SetTitle("#eta");
620  fHistograms->Add(fHistClusterEtavsPhiAfterQA);
621  } else{AliError(Form("Cluster Type is not EMCAL nor PHOS nor all: %i",fClusterType));}
622 
623  if(fIsMC > 1){
624  if(!fDoLightOutput){
625  fHistClusterEtavsPhiBeforeAcc->Sumw2();
626  fHistClusterEtavsPhiAfterAcc->Sumw2();
627  }
628  fHistClusterEtavsPhiAfterQA->Sumw2();
629  }
630 
631  // Cluster quality related histograms
632  Double_t timeMin = -2e-6;
633  Double_t timeMax = 8e-6;
634  if( fClusterType == 1 || fClusterType == 3){
635  timeMin = -2e-7;
636  timeMax = 12e-7;
637  }
638 
639  if(!fDoLightOutput){
640  fHistClusterTimevsEBeforeQA = new TH2F(Form("ClusterTimeVsE_beforeClusterQA %s",GetCutNumber().Data()),"ClusterTimeVsE_beforeClusterQA",800, timeMin, timeMax,
641  nBinsClusterECoarse, minClusterELog, maxClusterELog);
642  fHistClusterTimevsEBeforeQA->GetXaxis()->SetTitle("t_{cl} (s)");
643  fHistClusterTimevsEBeforeQA->GetYaxis()->SetTitle("E_{cl} (GeV)");
644 
645  SetLogBinningYTH2(fHistClusterTimevsEBeforeQA);
646  fHistograms->Add(fHistClusterTimevsEBeforeQA);
647  }
648  fHistClusterTimevsEAfterQA = new TH2F(Form("ClusterTimeVsE_afterClusterQA %s",GetCutNumber().Data()),"ClusterTimeVsE_afterClusterQA",800, timeMin, timeMax,
649  nBinsClusterECoarse, minClusterELog, maxClusterELog);
650  fHistClusterTimevsEAfterQA->GetXaxis()->SetTitle("t_{cl} (s)");
651  fHistClusterTimevsEAfterQA->GetYaxis()->SetTitle("E_{cl} (GeV)");
652  SetLogBinningYTH2(fHistClusterTimevsEAfterQA);
653  fHistograms->Add(fHistClusterTimevsEAfterQA);
654  if(fUseNonLinearity && !fDoLightOutput){
655  fHistEnergyOfClusterBeforeNL = new TH1F(Form("EnergyOfCluster_beforeNonLinearity %s",GetCutNumber().Data()), "EnergyOfCluster_beforeNonLinearity",
656  nBinsClusterEFine, minClusterELog, maxClusterELog);
657  fHistEnergyOfClusterBeforeNL->GetXaxis()->SetTitle("E_{cl} (GeV)");
658  SetLogBinningXTH1(fHistEnergyOfClusterBeforeNL);
659  fHistograms->Add(fHistEnergyOfClusterBeforeNL);
660  fHistEnergyOfClusterAfterNL = new TH1F(Form("EnergyOfCluster_afterNonLinearity %s",GetCutNumber().Data()), "EnergyOfCluster_afterNonLinearity",
661  nBinsClusterEFine, minClusterELog, maxClusterELog);
662  fHistEnergyOfClusterAfterNL->GetXaxis()->SetTitle("E_{cl} (GeV)");
663  SetLogBinningXTH1(fHistEnergyOfClusterAfterNL);
664  fHistograms->Add(fHistEnergyOfClusterAfterNL);
665  }
666  fHistEnergyOfClusterBeforeQA = new TH1F(Form("EnergyOfCluster_beforeClusterQA %s",GetCutNumber().Data()),"EnergyOfCluster_beforeClusterQA", nBinsClusterEFine, minClusterELog, maxClusterELog);
667  fHistEnergyOfClusterBeforeQA->GetXaxis()->SetTitle("E_{cl} (GeV)");
668  SetLogBinningXTH1(fHistEnergyOfClusterBeforeQA);
669  fHistograms->Add(fHistEnergyOfClusterBeforeQA);
670  fHistEnergyOfClusterAfterQA = new TH1F(Form("EnergyOfCluster_afterClusterQA %s",GetCutNumber().Data()),"EnergyOfCluster_afterClusterQA",nBinsClusterEFine, minClusterELog, maxClusterELog);
671  fHistEnergyOfClusterAfterQA->GetXaxis()->SetTitle("E_{cl} (GeV)");
672  SetLogBinningXTH1(fHistEnergyOfClusterAfterQA);
673  fHistograms->Add(fHistEnergyOfClusterAfterQA);
674  if(!fDoLightOutput){
675  fHistNCellsBeforeQA = new TH1F(Form("NCellPerCluster_beforeClusterQA %s",GetCutNumber().Data()),"NCellPerCluster_beforeClusterQA",50,0,50);
676  fHistNCellsBeforeQA->GetXaxis()->SetTitle("N_{cells} in cluster");
677  fHistograms->Add(fHistNCellsBeforeQA);
678  fHistNCellsAfterQA = new TH1F(Form("NCellPerCluster_afterClusterQA %s",GetCutNumber().Data()),"NCellPerCluster_afterClusterQA",50,0,50);
679  fHistNCellsAfterQA->GetXaxis()->SetTitle("N_{cells} in cluster");
680  fHistograms->Add(fHistNCellsAfterQA);
681  fHistM02BeforeQA = new TH1F(Form("M02_beforeClusterQA %s",GetCutNumber().Data()),"M02_beforeClusterQA",400,0,5);
682  fHistM02BeforeQA->GetXaxis()->SetTitle("#sigma_{long}^2");
683  fHistograms->Add(fHistM02BeforeQA);
684  fHistM02AfterQA = new TH1F(Form("M02_afterClusterQA %s",GetCutNumber().Data()),"M02_afterClusterQA",400,0,5);
685  fHistM02AfterQA->GetXaxis()->SetTitle("#sigma_{long}^2");
686  fHistograms->Add(fHistM02AfterQA);
687  fHistM20BeforeQA = new TH1F(Form("M20_beforeClusterQA %s",GetCutNumber().Data()),"M20_beforeClusterQA",400,0,2.5);
688  fHistM20BeforeQA->GetXaxis()->SetTitle("#sigma_{short}^2");
689  fHistograms->Add(fHistM20BeforeQA);
690  fHistM20AfterQA = new TH1F(Form("M20_afterClusterQA %s",GetCutNumber().Data()),"M20_afterClusterQA",400,0,2.5);
691  fHistM20AfterQA->GetXaxis()->SetTitle("#sigma_{short}^2");
692  fHistograms->Add(fHistM20AfterQA);
693  fHistDispersionBeforeQA = new TH1F(Form("Dispersion_beforeClusterQA %s",GetCutNumber().Data()),"Dispersion_beforeClusterQA",100,0,4);
694  fHistDispersionBeforeQA->GetXaxis()->SetTitle("dispersion");
695  fHistograms->Add(fHistDispersionBeforeQA);
696  fHistDispersionAfterQA = new TH1F(Form("Dispersion_afterClusterQA %s",GetCutNumber().Data()),"Dispersion_afterClusterQA",100,0,4);
697  fHistDispersionAfterQA->GetXaxis()->SetTitle("dispersion");
698  fHistograms->Add(fHistDispersionAfterQA);
699  fHistNLMBeforeQA = new TH1F(Form("NLM_beforeClusterQA %s",GetCutNumber().Data()),"NLM_beforeClusterQA",10,0,10);
700  fHistNLMBeforeQA->GetXaxis()->SetTitle("N_{LM} in cluster");
701  fHistograms->Add(fHistNLMBeforeQA);
702  fHistNLMAfterQA = new TH1F(Form("NLM_afterClusterQA %s",GetCutNumber().Data()),"NLM_afterClusterQA",10,0,10);
703  fHistNLMAfterQA->GetXaxis()->SetTitle("N_{LM} in cluster");
704  fHistograms->Add(fHistNLMAfterQA);
705 // fHistNLMAvsNLMBBeforeQA = new TH2F(Form("NLMAvsNLMB_beforeClusterQA %s",GetCutNumber().Data()),"NLMAvsNLMB_beforeClusterQA",10,0,10,10,0,10);
706 // fHistograms->Add(fHistNLMAvsNLMBBeforeQA);
707  fHistNLMVsNCellsAfterQA = new TH2F(Form("NLM_NCells_afterClusterQA %s",GetCutNumber().Data()),"NLM_NCells_afterClusterQA",10,0,10,50,0,50);
708  fHistNLMVsNCellsAfterQA->GetXaxis()->SetTitle("N_{LM} in cluster");
709  fHistNLMVsNCellsAfterQA->GetYaxis()->SetTitle("N_{cells} in cluster");
710  fHistograms->Add(fHistNLMVsNCellsAfterQA);
711  fHistNLMVsEAfterQA = new TH2F(Form("NLM_E_afterClusterQA %s",GetCutNumber().Data()),"NLM_E_afterClusterQA",10,0,10,nBinsClusterECoarse,minClusterELog,maxClusterELog);
712  fHistNLMVsEAfterQA->GetXaxis()->SetTitle("N_{LM} in cluster");
713  fHistNLMVsEAfterQA->GetYaxis()->SetTitle("E_{cl} (GeV)");
714  SetLogBinningYTH2(fHistNLMVsEAfterQA);
715  fHistograms->Add(fHistNLMVsEAfterQA);
716  if(fExtendedMatchAndQA > 1 || fIsPureCalo > 0){
717  fHistClusterEM02AfterQA = new TH2F(Form("EVsM02_afterClusterQA %s",GetCutNumber().Data()),"EVsM02_afterClusterQA",nBinsClusterEFine, minClusterELog, maxClusterELog,400,0,5);
718  fHistClusterEM02AfterQA->GetYaxis()->SetTitle("#sigma_{long}^2");
719  fHistClusterEM02AfterQA->GetXaxis()->SetTitle("E_{cl} (GeV)");
720  SetLogBinningXTH2(fHistClusterEM02AfterQA);
721  fHistograms->Add(fHistClusterEM02AfterQA);
722  fHistClusterEM02BeforeQA = new TH2F(Form("EVsM02_beforeClusterQA %s",GetCutNumber().Data()),"EVsM02_beforeClusterQA",nBinsClusterEFine, minClusterELog, maxClusterELog,400,0,5);
723  fHistClusterEM02BeforeQA->GetYaxis()->SetTitle("#sigma_{long}^2");
724  fHistClusterEM02BeforeQA->GetXaxis()->SetTitle("E_{cl} (GeV)");
725  SetLogBinningXTH2(fHistClusterEM02BeforeQA);
726  fHistograms->Add(fHistClusterEM02BeforeQA);
727  }
728  }
729  //----------------
730  if(fIsMC > 1){
731  if(!fDoLightOutput) fHistClusterTimevsEBeforeQA->Sumw2();
732  fHistClusterTimevsEAfterQA->Sumw2();
733  if(fUseNonLinearity && !fDoLightOutput){
734  fHistEnergyOfClusterBeforeNL->Sumw2();
735  fHistEnergyOfClusterAfterNL->Sumw2();
736  }
737  fHistEnergyOfClusterBeforeQA->Sumw2();
738  fHistEnergyOfClusterAfterQA->Sumw2();
739  if(!fDoLightOutput){
740  fHistNCellsBeforeQA->Sumw2();
741  fHistNCellsAfterQA->Sumw2();
742  fHistM02BeforeQA->Sumw2();
743  fHistM02AfterQA->Sumw2();
744  fHistM20BeforeQA->Sumw2();
745  fHistM20AfterQA->Sumw2();
746  fHistDispersionBeforeQA->Sumw2();
747  fHistDispersionAfterQA->Sumw2();
748  fHistNLMBeforeQA->Sumw2();
749  fHistNLMAfterQA->Sumw2();
750 // fHistNLMAvsNLMBBeforeQA->Sumw2();
751  fHistNLMVsNCellsAfterQA->Sumw2();
752  fHistNLMVsEAfterQA->Sumw2();
753  if(fExtendedMatchAndQA > 1 || fIsPureCalo > 0){
754  fHistClusterEM02AfterQA->Sumw2();
755  fHistClusterEM02BeforeQA->Sumw2();
756  }
757  }
758  }
759 //----------------
760  if(!fDoLightOutput){
761  if( fClusterType == 1 ){
762  Int_t nMaxCellsEMCAL = fNMaxEMCalModules*48*24;
763  fBadChannels = new TProfile("EMCal - Bad Channels","EMCal - Bad Channels",nMaxCellsEMCAL,-0.5,nMaxCellsEMCAL-0.5);
764  fBadChannels->GetXaxis()->SetTitle("channel ID");
765  fHistExtQA->Add(fBadChannels);
766  } else if( fClusterType == 2 ){
767  Int_t nMaxCellsPHOS = fNMaxPHOSModules*56*64;
768  fBadChannels = new TProfile("PHOS - Bad Channels","PHOS - Bad Channels",nMaxCellsPHOS,-0.5,nMaxCellsPHOS-0.5);
769  fBadChannels->GetXaxis()->SetTitle("channel ID");
770  fHistExtQA->Add(fBadChannels);
771  } else if( fClusterType == 3 ){
772  Int_t nStartCellDCAL = 12288;
773  Int_t nMaxCellsDCAL = fNMaxDCalModules*32*24;
774  fBadChannels = new TProfile("DCAL - Bad Channels","DCAL - Bad Channels",nMaxCellsDCAL,nStartCellDCAL-0.5,nStartCellDCAL+nMaxCellsDCAL-0.5);
775  fBadChannels->GetXaxis()->SetTitle("channel ID");
776  fHistExtQA->Add(fBadChannels);
777  }
778  }
779 
780  if(fExtendedMatchAndQA > 1){
781  if( fClusterType == 1 ){ //EMCAL
782  // detailed cluster QA histos for EMCAL
783  Int_t nMaxCellsEMCAL = fNMaxEMCalModules*48*24;
784  fHistClusterEnergyvsMod = new TH2F(Form("ClusterEnergyVsModule_afterClusterQA %s",GetCutNumber().Data()),"ClusterEnergyVsModule_afterClusterQA",
785  nBinsClusterEFine, minClusterELog, maxClusterELog, fNMaxEMCalModules, 0, fNMaxEMCalModules);
786  fHistClusterEnergyvsMod->GetXaxis()->SetTitle("E_{cl} (GeV)");
787  fHistClusterEnergyvsMod->GetYaxis()->SetTitle("module ID");
788  SetLogBinningXTH2(fHistClusterEnergyvsMod);
789  fHistExtQA->Add(fHistClusterEnergyvsMod);
790  fHistNCellsBigger100MeVvsMod = new TH2F(Form("NCellsAbove100VsModule %s",GetCutNumber().Data()),"NCellsAbove100VsModule",200,0,200,fNMaxEMCalModules,0,fNMaxEMCalModules);
791  fHistNCellsBigger100MeVvsMod->GetXaxis()->SetTitle("N_{cells} with E_{cell} > 0.1 GeV");
792  fHistNCellsBigger100MeVvsMod->GetYaxis()->SetTitle("module ID");
793  fHistExtQA->Add(fHistNCellsBigger100MeVvsMod);
794  fHistNCellsBigger1500MeVvsMod = new TH2F(Form("NCellsAbove1500VsModule %s",GetCutNumber().Data()),"NCellsAbove1500VsModule",100,0,100,fNMaxEMCalModules,0,fNMaxEMCalModules);
795  fHistNCellsBigger1500MeVvsMod->GetXaxis()->SetTitle("N_{cells} with E_{cell} > 1.5 GeV");
796  fHistNCellsBigger1500MeVvsMod->GetYaxis()->SetTitle("module ID");
797  fHistExtQA->Add(fHistNCellsBigger1500MeVvsMod);
798  fHistEnergyOfModvsMod = new TH2F(Form("ModuleEnergyVsModule %s",GetCutNumber().Data()),"ModuleEnergyVsModule",nBinsModuleECoarse, minModuleELog, maxModuleELog,
799  fNMaxEMCalModules,0,fNMaxEMCalModules);
800  fHistEnergyOfModvsMod->GetXaxis()->SetTitle("E_{mod} (GeV)");
801  fHistEnergyOfModvsMod->GetYaxis()->SetTitle("module ID");
802  SetLogBinningXTH2(fHistEnergyOfModvsMod);
803  fHistExtQA->Add(fHistEnergyOfModvsMod);
804  fHistClusterEnergyvsNCells = new TH2F(Form("ClusterEnergyVsNCells_afterQA %s",GetCutNumber().Data()),"ClusterEnergyVsNCells_afterQA",nBinsClusterECoarse, minClusterELog, maxClusterELog,
805  50, 0, 50);
806  fHistEnergyOfModvsMod->GetXaxis()->SetTitle("E_{cl} (GeV)");
807  fHistEnergyOfModvsMod->GetYaxis()->SetTitle("N_{cells}");
808  SetLogBinningXTH2(fHistClusterEnergyvsNCells);
809  fHistExtQA->Add(fHistClusterEnergyvsNCells);
810  fHistClusterDistanceInTimeCut = new TH2F(Form("ClusterDistanceTo_withinTimingCut %s",GetCutNumber().Data()),"ClusterDistanceTo_withinTimingCut",20,-10,10,20,-10,10);
811  fHistClusterDistanceInTimeCut->GetXaxis()->SetTitle("R_{cl,row} within time cut (cell)");
812  fHistClusterDistanceInTimeCut->GetYaxis()->SetTitle("R_{cl,col} within time cut (cell)");
813 
814  fHistExtQA->Add(fHistClusterDistanceInTimeCut);
815  fHistClusterDistanceOutTimeCut = new TH2F(Form("ClusterDistanceTo_outsideTimingCut %s",GetCutNumber().Data()),"ClusterDistanceTo_outsideTimingCut",20,-10,10,20,-10,10);
816  fHistClusterDistanceOutTimeCut->GetXaxis()->SetTitle("R_{cl,row} outside time cut (cell)");
817  fHistClusterDistanceOutTimeCut->GetYaxis()->SetTitle("R_{cl,col} outside time cut (cell)");
818 
819  fHistExtQA->Add(fHistClusterDistanceOutTimeCut);
820  fHistClusterDistance1DInTimeCut = new TH1F(Form("Cluster1D_DistanceTo_withinTimingCut %s",GetCutNumber().Data()),"Cluster1D_DistanceTo_withinTimingCut",200,0.,0.5);
821  fHistClusterDistance1DInTimeCut->GetXaxis()->SetTitle("R_{cl,1D} within time cut (cell)");
822 
823  fHistExtQA->Add(fHistClusterDistance1DInTimeCut);
824 
825  // detailed cell QA histos for EMCAL
826  if(fExtendedMatchAndQA > 3){
827  fHistCellEnergyvsCellID = new TH2F(Form("CellEnergyVsCellID %s",GetCutNumber().Data()),"CellEnergyVsCellID",nBinsCellECoarse, minCellELog, maxCellELog,
828  nMaxCellsEMCAL, 0, nMaxCellsEMCAL);
829  fHistCellEnergyvsCellID->GetXaxis()->SetTitle("E_{cell} (GeV)");
830  fHistCellEnergyvsCellID->GetYaxis()->SetTitle("Cell ID");
831  SetLogBinningXTH2(fHistCellEnergyvsCellID);
832  fHistExtQA->Add(fHistCellEnergyvsCellID);
833  fHistCellTimevsCellID = new TH2F(Form("CellTimeVsCellID %s",GetCutNumber().Data()),"CellTimeVsCellID",600,-timeMax,timeMax,nMaxCellsEMCAL,0,nMaxCellsEMCAL);
834  fHistCellTimevsCellID->GetXaxis()->SetTitle("t_{cell} (GeV)");
835  fHistCellTimevsCellID->GetYaxis()->SetTitle("Cell ID");
836  fHistExtQA->Add(fHistCellTimevsCellID);
837  fHistClusterIncludedCellsBeforeQA = new TH1F(Form("ClusterIncludedCells_beforeClusterQA %s",GetCutNumber().Data()),"ClusterIncludedCells_beforeClusterQA",
838  nMaxCellsEMCAL,0,nMaxCellsEMCAL);
839  fHistClusterIncludedCellsBeforeQA->GetXaxis()->SetTitle("Cell ID");
840  fHistClusterIncludedCellsBeforeQA->GetYaxis()->SetTitle("# included in cluster");
841  fHistExtQA->Add(fHistClusterIncludedCellsBeforeQA);
842  fHistClusterIncludedCellsAfterQA = new TH1F(Form("ClusterIncludedCells_afterClusterQA %s",GetCutNumber().Data()),"ClusterIncludedCells_afterClusterQA",
843  nMaxCellsEMCAL,0,nMaxCellsEMCAL);
844  fHistClusterIncludedCellsAfterQA->GetXaxis()->SetTitle("Cell ID");
845  fHistClusterIncludedCellsAfterQA->GetYaxis()->SetTitle("# included in cluster");
846  fHistExtQA->Add(fHistClusterIncludedCellsAfterQA);
847  fHistClusterEnergyFracCellsBeforeQA = new TH1F(Form("ClusterEnergyFracCells_beforeClusterQA %s",GetCutNumber().Data()),"ClusterEnergyFracCells_beforeClusterQA",
848  nMaxCellsEMCAL,0,nMaxCellsEMCAL);
849  fHistClusterEnergyFracCellsBeforeQA->GetXaxis()->SetTitle("Cell ID");
850  fHistClusterEnergyFracCellsBeforeQA->GetYaxis()->SetTitle("fraction of energy cell is contrib. to cl.");
851  fHistClusterEnergyFracCellsBeforeQA->Sumw2();
852  fHistExtQA->Add(fHistClusterEnergyFracCellsBeforeQA);
853  fHistClusterEnergyFracCellsAfterQA = new TH1F(Form("ClusterEnergyFracCells_afterClusterQA %s",GetCutNumber().Data()),"ClusterEnergyFracCells_afterClusterQA",
854  nMaxCellsEMCAL,0,nMaxCellsEMCAL);
855  fHistClusterEnergyFracCellsAfterQA->GetXaxis()->SetTitle("Cell ID");
856  fHistClusterEnergyFracCellsAfterQA->GetYaxis()->SetTitle("fraction of energy cell is contrib. to cl.");
857  fHistClusterEnergyFracCellsAfterQA->Sumw2();
858  fHistExtQA->Add(fHistClusterEnergyFracCellsAfterQA);
859  fHistClusterIncludedCellsTimingAfterQA = new TH2F(Form("ClusterIncludedCellsTiming_afterClusterQA %s",GetCutNumber().Data()),"ClusterIncludedCellsTiming_afterClusterQA",
860  nBinsClusterECoarse, minClusterELog, maxClusterELog, 200, -500, 500);
861  SetLogBinningXTH2(fHistClusterIncludedCellsTimingAfterQA);
862  fHistClusterIncludedCellsTimingAfterQA->GetXaxis()->SetTitle("E_{cl} (GeV)");
863  fHistClusterIncludedCellsTimingAfterQA->GetYaxis()->SetTitle("t_{cell} (ns)");
864  fHistExtQA->Add(fHistClusterIncludedCellsTimingAfterQA);
865  fHistClusterIncludedCellsTimingEnergyAfterQA = new TH2F(Form("ClusterIncludedCellsTimingEnergy_afterClusterQA %s",GetCutNumber().Data()),"ClusterIncludedCellsTimingEnergy_afterClusterQA",
866  nBinsCellECoarse, minCellELog, maxCellELog, 200, -500, 500);
867  fHistClusterIncludedCellsTimingEnergyAfterQA->GetXaxis()->SetTitle("E_{cell} (GeV)");
868  fHistClusterIncludedCellsTimingEnergyAfterQA->GetYaxis()->SetTitle("t_{cell} (ns)");
869  SetLogBinningXTH2(fHistClusterIncludedCellsTimingEnergyAfterQA);
870  fHistExtQA->Add(fHistClusterIncludedCellsTimingEnergyAfterQA);
871  }
872 
873  } else if( fClusterType == 2 ){ //PHOS
874  // detailed cluster QA histos for PHOS
875  Int_t nMaxCellsPHOS = fNMaxPHOSModules*56*64;
876  fHistClusterEnergyvsMod = new TH2F(Form("ClusterEnergyVsModule_afterClusterQA %s",GetCutNumber().Data()),"ClusterEnergyVsModule_afterClusterQA",nBinsClusterEFine, minClusterELog, maxClusterELog,
877  fNMaxPHOSModules, 0, fNMaxPHOSModules);
878  fHistClusterEnergyvsMod->GetXaxis()->SetTitle("E_{cl} (GeV)");
879  fHistClusterEnergyvsMod->GetYaxis()->SetTitle("module ID");
880  SetLogBinningXTH2(fHistClusterEnergyvsMod);
881  fHistExtQA->Add(fHistClusterEnergyvsMod);
882  fHistNCellsBigger100MeVvsMod = new TH2F(Form("NCellsAbove100VsModule %s",GetCutNumber().Data()),"NCellsAbove100VsModule",200,0,200,fNMaxPHOSModules,0,fNMaxPHOSModules);
883  fHistNCellsBigger100MeVvsMod->GetXaxis()->SetTitle("N_{cells} with E_{cell} > 0.1 GeV");
884  fHistNCellsBigger100MeVvsMod->GetYaxis()->SetTitle("module ID");
885  fHistExtQA->Add(fHistNCellsBigger100MeVvsMod);
886  fHistNCellsBigger1500MeVvsMod = new TH2F(Form("NCellsAbove1500VsModule %s",GetCutNumber().Data()),"NCellsAbove1500VsModule",100,0,100,fNMaxPHOSModules,0,fNMaxPHOSModules);
887  fHistExtQA->Add(fHistNCellsBigger1500MeVvsMod);
888  fHistNCellsBigger1500MeVvsMod->GetXaxis()->SetTitle("N_{cells} with E_{cell} > 1.5 GeV");
889  fHistNCellsBigger1500MeVvsMod->GetYaxis()->SetTitle("module ID");
890  fHistEnergyOfModvsMod = new TH2F(Form("ModuleEnergyVsModule %s",GetCutNumber().Data()),"ModuleEnergyVsModule",nBinsModuleECoarse, minModuleELog, maxModuleELog,
891  fNMaxPHOSModules, 0, fNMaxPHOSModules);
892  fHistEnergyOfModvsMod->GetXaxis()->SetTitle("E_{cl} (GeV)");
893  fHistEnergyOfModvsMod->GetYaxis()->SetTitle("N_{cells}");
894  SetLogBinningXTH2(fHistEnergyOfModvsMod);
895  fHistExtQA->Add(fHistEnergyOfModvsMod);
896  fHistClusterEnergyvsNCells = new TH2F(Form("ClusterEnergyVsNCells_afterQA %s",GetCutNumber().Data()),"ClusterEnergyVsNCells_afterQA",nBinsClusterECoarse, minClusterELog, maxClusterELog, 50, 0, 50);
897  SetLogBinningXTH2(fHistClusterEnergyvsNCells);
898  fHistExtQA->Add(fHistClusterEnergyvsNCells);
899  fHistClusterDistanceInTimeCut = new TH2F(Form("ClusterDistanceTo_withinTimingCut %s",GetCutNumber().Data()),"ClusterDistanceTo_withinTimingCut",20,-10,10,20,-10,10);
900  fHistClusterDistanceInTimeCut->GetXaxis()->SetTitle("R_{cl,row} within time cut (cell)");
901  fHistClusterDistanceInTimeCut->GetYaxis()->SetTitle("R_{cl,col} within time cut (cell)");
902  fHistExtQA->Add(fHistClusterDistanceInTimeCut);
903  fHistClusterDistanceOutTimeCut = new TH2F(Form("ClusterDistanceTo_outsideTimingCut %s",GetCutNumber().Data()),"ClusterDistanceTo_outsideTimingCut",20,-10,10,20,-10,10);
904  fHistClusterDistanceOutTimeCut->GetXaxis()->SetTitle("R_{cl,row} outside time cut (cell)");
905  fHistClusterDistanceOutTimeCut->GetYaxis()->SetTitle("R_{cl,col} outside time cut (cell)");
906  fHistExtQA->Add(fHistClusterDistanceOutTimeCut);
907  fHistClusterDistance1DInTimeCut = new TH1F(Form("Cluster1D_DistanceTo_withinTimingCut %s",GetCutNumber().Data()),"Cluster1D_DistanceTo_withinTimingCut",200,0.,0.5);
908  fHistClusterDistance1DInTimeCut->GetXaxis()->SetTitle("R_{cl,1D} within time cut (cell)");
909  fHistExtQA->Add(fHistClusterDistance1DInTimeCut);
910 
911  // detailed cell QA histos for PHOS
912  if(fExtendedMatchAndQA > 3){
913  fHistCellEnergyvsCellID = new TH2F(Form("CellEnergyVsCellID %s",GetCutNumber().Data()),"CellEnergyVsCellID",nBinsCellECoarse, minCellELog, maxCellELog,
914  nMaxCellsPHOS,0,nMaxCellsPHOS);
915  fHistCellEnergyvsCellID->GetXaxis()->SetTitle("E_{cell} (GeV)");
916  fHistCellEnergyvsCellID->GetYaxis()->SetTitle("Cell ID");
917  SetLogBinningXTH2(fHistCellEnergyvsCellID);
918  fHistExtQA->Add(fHistCellEnergyvsCellID);
919  fHistCellTimevsCellID = new TH2F(Form("CellTimeVsCellID %s",GetCutNumber().Data()),"CellTimeVsCellID",600,-timeMax,timeMax,nMaxCellsPHOS,0,nMaxCellsPHOS);
920  fHistCellTimevsCellID->GetXaxis()->SetTitle("t_{cell} (GeV)");
921  fHistCellTimevsCellID->GetYaxis()->SetTitle("Cell ID");
922  fHistExtQA->Add(fHistCellTimevsCellID);
923  fHistClusterIncludedCellsBeforeQA = new TH1F(Form("ClusterIncludedCells_beforeClusterQA %s",GetCutNumber().Data()),"ClusterIncludedCells_beforeClusterQA",
924  nMaxCellsPHOS, 0, nMaxCellsPHOS);
925  fHistClusterIncludedCellsBeforeQA->GetXaxis()->SetTitle("Cell ID");
926  fHistClusterIncludedCellsBeforeQA->GetYaxis()->SetTitle("# included in cluster");
927  fHistExtQA->Add(fHistClusterIncludedCellsBeforeQA);
928  fHistClusterIncludedCellsAfterQA = new TH1F(Form("ClusterIncludedCells_afterClusterQA %s",GetCutNumber().Data()),"ClusterIncludedCells_afterClusterQA",
929  nMaxCellsPHOS, 0, nMaxCellsPHOS);
930  fHistClusterIncludedCellsAfterQA->GetXaxis()->SetTitle("Cell ID");
931  fHistClusterIncludedCellsAfterQA->GetYaxis()->SetTitle("# included in cluster");
932  fHistExtQA->Add(fHistClusterIncludedCellsAfterQA);
933  fHistClusterEnergyFracCellsBeforeQA = new TH1F(Form("ClusterEnergyFracCells_beforeClusterQA %s",GetCutNumber().Data()),"ClusterEnergyFracCells_beforeClusterQA",
934  nMaxCellsPHOS,0,nMaxCellsPHOS);
935  fHistClusterEnergyFracCellsBeforeQA->GetXaxis()->SetTitle("Cell ID");
936  fHistClusterEnergyFracCellsBeforeQA->GetYaxis()->SetTitle("fraction of energy cell is contrib. to cl.");
937  fHistClusterEnergyFracCellsBeforeQA->Sumw2();
938  fHistExtQA->Add(fHistClusterEnergyFracCellsBeforeQA);
939  fHistClusterEnergyFracCellsAfterQA = new TH1F(Form("ClusterEnergyFracCells_afterClusterQA %s",GetCutNumber().Data()),"ClusterEnergyFracCells_afterClusterQA",
940  nMaxCellsPHOS, 0, nMaxCellsPHOS);
941  fHistClusterEnergyFracCellsAfterQA->GetXaxis()->SetTitle("Cell ID");
942  fHistClusterEnergyFracCellsAfterQA->GetYaxis()->SetTitle("fraction of energy cell is contrib. to cl.");
943  fHistClusterEnergyFracCellsAfterQA->Sumw2();
944  fHistExtQA->Add(fHistClusterEnergyFracCellsAfterQA);
945  fHistClusterIncludedCellsTimingAfterQA = new TH2F(Form("ClusterIncludedCellsTiming_afterClusterQA %s",GetCutNumber().Data()),"ClusterIncludedCellsTiming_afterClusterQA",
946  nBinsClusterECoarse, minClusterELog, maxClusterELog, 200, -500, 500);
947  fHistClusterIncludedCellsTimingAfterQA->GetXaxis()->SetTitle("E_{cl} (GeV)");
948  fHistClusterIncludedCellsTimingAfterQA->GetYaxis()->SetTitle("t_{cell} (ns)");
949  SetLogBinningXTH2(fHistClusterIncludedCellsTimingAfterQA);
950  fHistExtQA->Add(fHistClusterIncludedCellsTimingAfterQA);
951  fHistClusterIncludedCellsTimingEnergyAfterQA = new TH2F(Form("ClusterIncludedCellsTimingEnergy_afterClusterQA %s",GetCutNumber().Data()),"ClusterIncludedCellsTimingEnergy_afterClusterQA",
952  nBinsCellECoarse, minCellELog, maxCellELog, 200, -500, 500);
953  fHistClusterIncludedCellsTimingEnergyAfterQA->GetXaxis()->SetTitle("E_{cell} (GeV)");
954  fHistClusterIncludedCellsTimingEnergyAfterQA->GetYaxis()->SetTitle("t_{cell} (ns)");
955  SetLogBinningXTH2(fHistClusterIncludedCellsTimingEnergyAfterQA);
956  fHistExtQA->Add(fHistClusterIncludedCellsTimingEnergyAfterQA);
957  }
958  } else if( fClusterType == 3 ){ //DCAL
959  Int_t nModulesStart = 12;
960  Int_t nCellsStart = 12288;
961  Int_t nMaxCellsDCAL = fNMaxDCalModules*32*24;
962 
963  fHistClusterEnergyvsMod = new TH2F(Form("ClusterEnergyVsModule_afterClusterQA %s",GetCutNumber().Data()),"ClusterEnergyVsModule_afterClusterQA",
964  nBinsClusterEFine, minClusterELog, maxClusterELog, fNMaxDCalModules, nModulesStart, fNMaxDCalModules+nModulesStart);
965  fHistClusterEnergyvsMod->GetXaxis()->SetTitle("E_{cl} (GeV)");
966  fHistClusterEnergyvsMod->GetYaxis()->SetTitle("module ID");
967  SetLogBinningXTH2(fHistClusterEnergyvsMod);
968  fHistExtQA->Add(fHistClusterEnergyvsMod);
969  fHistNCellsBigger100MeVvsMod = new TH2F(Form("NCellsAbove100VsModule %s",GetCutNumber().Data()),"NCellsAbove100VsModule",200 , 0, 200,
970  fNMaxDCalModules, nModulesStart, fNMaxDCalModules+nModulesStart);
971  fHistNCellsBigger100MeVvsMod->GetXaxis()->SetTitle("N_{cells} with E_{cell} > 0.1 GeV");
972  fHistNCellsBigger100MeVvsMod->GetYaxis()->SetTitle("module ID");
973  fHistExtQA->Add(fHistNCellsBigger100MeVvsMod);
974  fHistNCellsBigger1500MeVvsMod = new TH2F(Form("NCellsAbove1500VsModule %s",GetCutNumber().Data()),"NCellsAbove1500VsModule", 100, 0, 100,
975  fNMaxDCalModules,nModulesStart,nModulesStart+fNMaxDCalModules);
976  fHistNCellsBigger1500MeVvsMod->GetXaxis()->SetTitle("N_{cells} with E_{cell} > 1.5 GeV");
977  fHistNCellsBigger1500MeVvsMod->GetYaxis()->SetTitle("module ID");
978  fHistExtQA->Add(fHistNCellsBigger1500MeVvsMod);
979  fHistEnergyOfModvsMod = new TH2F(Form("ModuleEnergyVsModule %s",GetCutNumber().Data()),"ModuleEnergyVsModule", 1000, 0, 100,
980  fNMaxDCalModules, nModulesStart, fNMaxDCalModules+nModulesStart);
981  fHistEnergyOfModvsMod->GetXaxis()->SetTitle("E_{cl} (GeV)");
982  fHistEnergyOfModvsMod->GetYaxis()->SetTitle("N_{cells}");
983  SetLogBinningXTH2(fHistEnergyOfModvsMod);
984  fHistExtQA->Add(fHistEnergyOfModvsMod);
985  fHistClusterEnergyvsNCells = new TH2F(Form("ClusterEnergyVsNCells_afterQA %s",GetCutNumber().Data()),"ClusterEnergyVsNCells_afterQA",nBinsClusterECoarse, minClusterELog, maxClusterELog,
986  50, 0, 50);
987  SetLogBinningXTH2(fHistClusterEnergyvsNCells);
988  fHistExtQA->Add(fHistClusterEnergyvsNCells);
989  fHistClusterDistanceInTimeCut = new TH2F(Form("ClusterDistanceTo_withinTimingCut %s",GetCutNumber().Data()),"ClusterDistanceTo_withinTimingCut",20,-10,10,20,-10,10);
990  fHistClusterDistanceInTimeCut->GetXaxis()->SetTitle("R_{cl,row} within time cut (cell)");
991  fHistClusterDistanceInTimeCut->GetYaxis()->SetTitle("R_{cl,col} within time cut (cell)");
992  fHistExtQA->Add(fHistClusterDistanceInTimeCut);
993  fHistClusterDistanceOutTimeCut = new TH2F(Form("ClusterDistanceTo_outsideTimingCut %s",GetCutNumber().Data()),"ClusterDistanceTo_outsideTimingCut",20,-10,10,20,-10,10);
994  fHistClusterDistanceOutTimeCut->GetXaxis()->SetTitle("R_{cl,row} outside time cut (cell)");
995  fHistClusterDistanceOutTimeCut->GetYaxis()->SetTitle("R_{cl,col} outside time cut (cell)");
996  fHistExtQA->Add(fHistClusterDistanceOutTimeCut);
997  fHistClusterDistance1DInTimeCut = new TH1F(Form("Cluster1D_DistanceTo_withinTimingCut %s",GetCutNumber().Data()),"Cluster1D_DistanceTo_withinTimingCut",200,0.,0.5);
998  fHistClusterDistance1DInTimeCut->GetXaxis()->SetTitle("R_{cl,1D} within time cut (cell)");
999  fHistExtQA->Add(fHistClusterDistance1DInTimeCut);
1000 
1001  // detailed cell QA histos for DCAL
1002  if(fExtendedMatchAndQA > 3){
1003  fHistCellEnergyvsCellID = new TH2F(Form("CellEnergyVsCellID %s",GetCutNumber().Data()),"CellEnergyVsCellID",nBinsCellECoarse, minCellELog, maxCellELog,
1004  nMaxCellsDCAL, nCellsStart, nMaxCellsDCAL+nCellsStart);
1005  fHistCellEnergyvsCellID->GetXaxis()->SetTitle("E_{cell} (GeV)");
1006  fHistCellEnergyvsCellID->GetYaxis()->SetTitle("Cell ID");
1007  SetLogBinningXTH2(fHistCellEnergyvsCellID);
1008  fHistExtQA->Add(fHistCellEnergyvsCellID);
1009  fHistCellTimevsCellID = new TH2F(Form("CellTimeVsCellID %s",GetCutNumber().Data()),"CellTimeVsCellID",600,-timeMax,timeMax, nMaxCellsDCAL, nCellsStart,
1010  nMaxCellsDCAL+nCellsStart);
1011  fHistCellTimevsCellID->GetXaxis()->SetTitle("t_{cell} (GeV)");
1012  fHistCellTimevsCellID->GetYaxis()->SetTitle("Cell ID");
1013  fHistExtQA->Add(fHistCellTimevsCellID);
1014  fHistClusterIncludedCellsBeforeQA = new TH1F(Form("ClusterIncludedCells_beforeClusterQA %s",GetCutNumber().Data()),"ClusterIncludedCells_beforeClusterQA",
1015  nMaxCellsDCAL, nCellsStart, nMaxCellsDCAL+nCellsStart);
1016  fHistClusterIncludedCellsBeforeQA->GetXaxis()->SetTitle("Cell ID");
1017  fHistClusterIncludedCellsBeforeQA->GetYaxis()->SetTitle("# included in cluster");
1018  fHistExtQA->Add(fHistClusterIncludedCellsBeforeQA);
1019  fHistClusterIncludedCellsAfterQA = new TH1F(Form("ClusterIncludedCells_afterClusterQA %s",GetCutNumber().Data()),"ClusterIncludedCells_afterClusterQA",
1020  nMaxCellsDCAL, nCellsStart, nMaxCellsDCAL+nCellsStart);
1021  fHistClusterIncludedCellsAfterQA->GetXaxis()->SetTitle("Cell ID");
1022  fHistClusterIncludedCellsAfterQA->GetYaxis()->SetTitle("# included in cluster");
1023  fHistExtQA->Add(fHistClusterIncludedCellsAfterQA);
1024  fHistClusterEnergyFracCellsBeforeQA = new TH1F(Form("ClusterEnergyFracCells_beforeClusterQA %s",GetCutNumber().Data()),"ClusterEnergyFracCells_beforeClusterQA",
1025  nMaxCellsDCAL, nCellsStart, nMaxCellsDCAL+nCellsStart);
1026  fHistClusterEnergyFracCellsBeforeQA->GetXaxis()->SetTitle("Cell ID");
1027  fHistClusterEnergyFracCellsBeforeQA->GetYaxis()->SetTitle("fraction of energy cell is contrib. to cl.");
1028  fHistClusterEnergyFracCellsBeforeQA->Sumw2();
1029  fHistExtQA->Add(fHistClusterEnergyFracCellsBeforeQA);
1030  fHistClusterEnergyFracCellsAfterQA = new TH1F(Form("ClusterEnergyFracCells_afterClusterQA %s",GetCutNumber().Data()),"ClusterEnergyFracCells_afterClusterQA",
1031  nMaxCellsDCAL, nCellsStart, nMaxCellsDCAL+nCellsStart);
1032  fHistClusterEnergyFracCellsAfterQA->GetXaxis()->SetTitle("Cell ID");
1033  fHistClusterEnergyFracCellsAfterQA->GetYaxis()->SetTitle("fraction of energy cell is contrib. to cl.");
1034  fHistClusterEnergyFracCellsAfterQA->Sumw2();
1035  fHistExtQA->Add(fHistClusterEnergyFracCellsAfterQA);
1036  fHistClusterIncludedCellsTimingAfterQA = new TH2F(Form("ClusterIncludedCellsTiming_afterClusterQA %s",GetCutNumber().Data()),"ClusterIncludedCellsTiming_afterClusterQA",
1037  nBinsClusterECoarse, minClusterELog, maxClusterELog, 200, -500, 500);
1038  fHistClusterIncludedCellsTimingAfterQA->GetXaxis()->SetTitle("E_{cl} (GeV)");
1039  fHistClusterIncludedCellsTimingAfterQA->GetYaxis()->SetTitle("t_{cell} (ns)");
1040  SetLogBinningXTH2(fHistClusterIncludedCellsTimingAfterQA);
1041  fHistExtQA->Add(fHistClusterIncludedCellsTimingAfterQA);
1042  fHistClusterIncludedCellsTimingEnergyAfterQA = new TH2F(Form("ClusterIncludedCellsTimingEnergy_afterClusterQA %s",GetCutNumber().Data()),"ClusterIncludedCellsTimingEnergy_afterClusterQA",
1043  nBinsCellECoarse, minCellELog, maxCellELog, 200, -500, 500);
1044  fHistClusterIncludedCellsTimingEnergyAfterQA->GetXaxis()->SetTitle("E_{cell} (GeV)");
1045  fHistClusterIncludedCellsTimingEnergyAfterQA->GetYaxis()->SetTitle("t_{cell} (ns)");
1046  SetLogBinningXTH2(fHistClusterIncludedCellsTimingEnergyAfterQA);
1047  fHistExtQA->Add(fHistClusterIncludedCellsTimingEnergyAfterQA);
1048  }
1049  } else{AliError(Form("fExtendedMatchAndQA (%i) not (yet) defined for cluster type (%i)",fExtendedMatchAndQA,fClusterType));}
1050  }
1051 
1052  //TrackMatching histograms
1053  if(fUseDistTrackToCluster && !fDoLightOutput) {
1054  const Int_t nEtaBins = 300;
1055  const Int_t nPhiBins = 300;
1056  const Float_t EtaRange[2] = {-0.3, 0.3};
1057  const Float_t PhiRange[2] = {-0.3, 0.3};
1058 
1059  fHistClusterdEtadPhiBeforeQA = new TH2F(Form("dEtaVsdPhi_beforeClusterQA %s",GetCutNumber().Data()),"dEtaVsdPhi_beforeClusterQA", nEtaBins, EtaRange[0], EtaRange[1], nPhiBins, PhiRange[0], PhiRange[1]);
1060  fHistograms->Add(fHistClusterdEtadPhiBeforeQA);
1061  fHistClusterdEtadPhiAfterQA = new TH2F(Form("dEtaVsdPhi_afterClusterQA %s",GetCutNumber().Data()),"dEtaVsdPhi_afterClusterQA", nEtaBins, EtaRange[0], EtaRange[1], nPhiBins, PhiRange[0], PhiRange[1]);
1062  fHistograms->Add(fHistClusterdEtadPhiAfterQA);
1063  //----------------
1064  if(fIsMC > 1){
1065  fHistClusterdEtadPhiBeforeQA->Sumw2();
1066  fHistClusterdEtadPhiAfterQA->Sumw2();
1067  }
1068  //----------------
1069  if(fExtendedMatchAndQA == 1 || fExtendedMatchAndQA == 3 || fExtendedMatchAndQA == 5 ){
1070 
1071  // QA histograms for track matching
1072  fHistClusterRBeforeQA = new TH1F(Form("R_Cluster_beforeClusterQA %s",GetCutNumber().Data()),"R of cluster",200,400,500);
1073  fHistClusterRBeforeQA->GetXaxis()->SetTitle("R_{prim vtx} (m)");
1074  fHistograms->Add(fHistClusterRBeforeQA);
1075  fHistClusterRAfterQA = new TH1F(Form("R_Cluster_afterClusterQA %s",GetCutNumber().Data()),"R of cluster_matched",200,400,500);
1076  fHistClusterRAfterQA->GetXaxis()->SetTitle("R_{prim vtx} (m)");
1077  fHistograms->Add(fHistClusterRAfterQA);
1078  fHistDistanceTrackToClusterBeforeQA = new TH1F(Form("DistanceToTrack_beforeClusterQA %s",GetCutNumber().Data()),"DistanceToTrack_beforeClusterQA",200,0,2);
1079  fHistDistanceTrackToClusterBeforeQA->GetXaxis()->SetTitle("R_{cl-track}");
1080  fHistograms->Add(fHistDistanceTrackToClusterBeforeQA);
1081  fHistDistanceTrackToClusterAfterQA = new TH1F(Form("DistanceToTrack_afterClusterQA %s",GetCutNumber().Data()),"DistanceToTrack_afterClusterQA",200,0,2);
1082  fHistDistanceTrackToClusterAfterQA->GetXaxis()->SetTitle("R_{cl-track}");
1083  fHistograms->Add(fHistDistanceTrackToClusterAfterQA);
1084  fHistClusterdEtadPhiPosTracksBeforeQA = new TH2F(Form("dEtaVsdPhi_posTracks_beforeClusterQA %s",GetCutNumber().Data()),"dEtaVsdPhi_posTracks_beforeClusterQA",
1085  nEtaBins, EtaRange[0], EtaRange[1], nPhiBins, PhiRange[0], PhiRange[1]);
1086  fHistDistanceTrackToClusterAfterQA->GetXaxis()->SetTitle("d#eta_{cl-track}");
1087  fHistDistanceTrackToClusterAfterQA->GetYaxis()->SetTitle("d#varphi_{cl-track} (rad)");
1088  fHistograms->Add(fHistClusterdEtadPhiPosTracksBeforeQA);
1089  fHistClusterdEtadPhiNegTracksBeforeQA = new TH2F(Form("dEtaVsdPhi_negTracks_beforeClusterQA %s",GetCutNumber().Data()),"dEtaVsdPhi_negTracks_beforeClusterQA",
1090  nEtaBins, EtaRange[0], EtaRange[1], nPhiBins, PhiRange[0], PhiRange[1]);
1091  fHistClusterdEtadPhiNegTracksBeforeQA->GetXaxis()->SetTitle("d#eta_{cl-track}");
1092  fHistClusterdEtadPhiNegTracksBeforeQA->GetYaxis()->SetTitle("d#varphi_{cl-track} (rad)");
1093  fHistograms->Add(fHistClusterdEtadPhiNegTracksBeforeQA);
1094  fHistClusterdEtadPhiPosTracksAfterQA = new TH2F(Form("dEtaVsdPhi_posTracks_afterClusterQA %s",GetCutNumber().Data()),"dEtaVsdPhi_posTracks_afterClusterQA",
1095  nEtaBins, EtaRange[0], EtaRange[1], nPhiBins, PhiRange[0], PhiRange[1]);
1096  fHistClusterdEtadPhiPosTracksAfterQA->GetXaxis()->SetTitle("d#eta_{cl-track}");
1097  fHistClusterdEtadPhiPosTracksAfterQA->GetYaxis()->SetTitle("d#varphi_{cl-track} (rad)");
1098  fHistograms->Add(fHistClusterdEtadPhiPosTracksAfterQA);
1099  fHistClusterdEtadPhiNegTracksAfterQA = new TH2F(Form("dEtaVsdPhi_negTracks_afterClusterQA %s",GetCutNumber().Data()),"dEtaVsdPhi_negTracks_afterClusterQA",
1100  nEtaBins, EtaRange[0], EtaRange[1], nPhiBins, PhiRange[0], PhiRange[1]);
1101  fHistClusterdEtadPhiNegTracksAfterQA->GetXaxis()->SetTitle("d#eta_{cl-track}");
1102  fHistClusterdEtadPhiNegTracksAfterQA->GetYaxis()->SetTitle("d#varphi_{cl-track} (rad)");
1103  fHistograms->Add(fHistClusterdEtadPhiNegTracksAfterQA);
1104  fHistClusterdEtadPhiPosTracksP_000_075BeforeQA = new TH2F(Form("dEtaVsdPhi_posTracks_P<0.75_beforeClusterQA %s",GetCutNumber().Data()),"dEtaVsdPhi_posTracks_P<0.75_beforeClusterQA",
1105  nEtaBins, EtaRange[0], EtaRange[1], nPhiBins, PhiRange[0], PhiRange[1]);
1106  fHistClusterdEtadPhiPosTracksP_000_075BeforeQA->GetXaxis()->SetTitle("d#eta_{cl-track}");
1107  fHistClusterdEtadPhiPosTracksP_000_075BeforeQA->GetYaxis()->SetTitle("d#varphi_{cl-track} (rad)");
1108  fHistograms->Add(fHistClusterdEtadPhiPosTracksP_000_075BeforeQA);
1109  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",
1110  nEtaBins, EtaRange[0], EtaRange[1], nPhiBins, PhiRange[0], PhiRange[1]);
1111  fHistClusterdEtadPhiPosTracksP_075_125BeforeQA->GetXaxis()->SetTitle("d#eta_{cl-track}");
1112  fHistClusterdEtadPhiPosTracksP_075_125BeforeQA->GetYaxis()->SetTitle("d#varphi_{cl-track} (rad)");
1113  fHistograms->Add(fHistClusterdEtadPhiPosTracksP_075_125BeforeQA);
1114  fHistClusterdEtadPhiPosTracksP_125_999BeforeQA = new TH2F(Form("dEtaVsdPhi_posTracks_P>1.25_beforeClusterQA %s",GetCutNumber().Data()),"dEtaVsdPhi_posTracks_P>1.25_beforeClusterQA",
1115  nEtaBins, EtaRange[0], EtaRange[1], nPhiBins, PhiRange[0], PhiRange[1]);
1116  fHistClusterdEtadPhiPosTracksP_125_999BeforeQA->GetXaxis()->SetTitle("d#eta_{cl-track}");
1117  fHistClusterdEtadPhiPosTracksP_125_999BeforeQA->GetYaxis()->SetTitle("d#varphi_{cl-track} (rad)");
1118  fHistograms->Add(fHistClusterdEtadPhiPosTracksP_125_999BeforeQA);
1119  fHistClusterdEtadPhiNegTracksP_000_075BeforeQA = new TH2F(Form("dEtaVsdPhi_negTracks_P<0.75_beforeClusterQA %s",GetCutNumber().Data()),"dEtaVsdPhi_negTracks_P<0.75_beforeClusterQA",
1120  nEtaBins, EtaRange[0], EtaRange[1], nPhiBins, PhiRange[0], PhiRange[1]);
1121  fHistClusterdEtadPhiNegTracksP_000_075BeforeQA->GetXaxis()->SetTitle("d#eta_{cl-track}");
1122  fHistClusterdEtadPhiNegTracksP_000_075BeforeQA->GetYaxis()->SetTitle("d#varphi_{cl-track} (rad)");
1123  fHistograms->Add(fHistClusterdEtadPhiNegTracksP_000_075BeforeQA);
1124  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",
1125  nEtaBins, EtaRange[0], EtaRange[1], nPhiBins, PhiRange[0], PhiRange[1]);
1126  fHistClusterdEtadPhiNegTracksP_075_125BeforeQA->GetXaxis()->SetTitle("d#eta_{cl-track}");
1127  fHistClusterdEtadPhiNegTracksP_075_125BeforeQA->GetYaxis()->SetTitle("d#varphi_{cl-track} (rad)");
1128  fHistograms->Add(fHistClusterdEtadPhiNegTracksP_075_125BeforeQA);
1129  fHistClusterdEtadPhiNegTracksP_125_999BeforeQA = new TH2F(Form("dEtaVsdPhi_negTracks_P>1.25_beforeClusterQA %s",GetCutNumber().Data()),"dEtaVsdPhi_negTracks_P>1.25_beforeClusterQA",
1130  nEtaBins, EtaRange[0], EtaRange[1], nPhiBins, PhiRange[0], PhiRange[1]);
1131  fHistClusterdEtadPhiNegTracksP_125_999BeforeQA->GetXaxis()->SetTitle("d#eta_{cl-track}");
1132  fHistClusterdEtadPhiNegTracksP_125_999BeforeQA->GetYaxis()->SetTitle("d#varphi_{cl-track} (rad)");
1133  fHistograms->Add(fHistClusterdEtadPhiNegTracksP_125_999BeforeQA);
1134  fHistClusterdEtadPtBeforeQA = new TH2F(Form("dEtaVsPt_beforeClusterQA %s",GetCutNumber().Data()),"dEtaVsPt_beforeClusterQA",nEtaBins,EtaRange[0],EtaRange[1],
1135  nBinsModuleECoarse, minClusterELog, maxClusterELog);
1136  fHistClusterdEtadPtBeforeQA->GetXaxis()->SetTitle("d#eta_{cl-track}");
1137  fHistClusterdEtadPtBeforeQA->GetYaxis()->SetTitle("p_{T} (GeV/c)");
1138  SetLogBinningYTH2(fHistClusterdEtadPtBeforeQA);
1139  fHistograms->Add(fHistClusterdEtadPtBeforeQA);
1140  fHistClusterdEtadPtAfterQA = new TH2F(Form("dEtaVsPt_afterClusterQA %s",GetCutNumber().Data()),"dEtaVsPt_afterClusterQA",nEtaBins,EtaRange[0],EtaRange[1],
1141  nBinsModuleECoarse, minClusterELog, maxClusterELog);
1142  fHistClusterdEtadPtAfterQA->GetXaxis()->SetTitle("d#eta_{cl-track}");
1143  fHistClusterdEtadPtAfterQA->GetYaxis()->SetTitle("p_{T} (GeV/c)");
1144  SetLogBinningYTH2(fHistClusterdEtadPtAfterQA);
1145  fHistograms->Add(fHistClusterdEtadPtAfterQA);
1146  fHistClusterdPhidPtPosTracksBeforeQA = new TH2F(Form("dPhiVsPt_posTracks_beforeClusterQA %s",GetCutNumber().Data()),"dPhiVsPt_posTracks_beforeClusterQA",
1147  2*nPhiBins, 2*PhiRange[0], 2*PhiRange[1], nBinsModuleECoarse, minClusterELog, maxClusterELog);
1148  fHistClusterdPhidPtPosTracksBeforeQA->GetXaxis()->SetTitle("d#varphi_{cl-track} (rad)");
1149  fHistClusterdPhidPtPosTracksBeforeQA->GetYaxis()->SetTitle("p_{T} (GeV/c)");
1150  SetLogBinningYTH2(fHistClusterdPhidPtPosTracksBeforeQA);
1151  fHistograms->Add(fHistClusterdPhidPtPosTracksBeforeQA);
1152  fHistClusterdPhidPtNegTracksBeforeQA = new TH2F(Form("dPhiVsPt_negTracks_beforeClusterQA %s",GetCutNumber().Data()),"dPhiVsPt_negTracks_beforeClusterQA",
1153  2*nPhiBins, 2*PhiRange[0], 2*PhiRange[1], nBinsModuleECoarse, minClusterELog, maxClusterELog);
1154  fHistClusterdPhidPtNegTracksBeforeQA->GetXaxis()->SetTitle("d#varphi_{cl-track} (rad)");
1155  fHistClusterdPhidPtNegTracksBeforeQA->GetYaxis()->SetTitle("p_{T} (GeV/c)");
1156  SetLogBinningYTH2(fHistClusterdPhidPtNegTracksBeforeQA);
1157  fHistograms->Add(fHistClusterdPhidPtNegTracksBeforeQA);
1158  fHistClusterdPhidPtAfterQA = new TH2F(Form("dPhiVsPt_afterClusterQA %s",GetCutNumber().Data()),"dPhiVsPt_afterClusterQA",2*nPhiBins,2*PhiRange[0],2*PhiRange[1],
1159  nBinsModuleECoarse, minClusterELog, maxClusterELog);
1160  fHistClusterdPhidPtAfterQA->GetXaxis()->SetTitle("d#varphi_{cl-track} (rad)");
1161  fHistClusterdPhidPtAfterQA->GetYaxis()->SetTitle("p_{T} (GeV/c)");
1162  SetLogBinningYTH2(fHistClusterdPhidPtAfterQA);
1163  fHistograms->Add(fHistClusterdPhidPtAfterQA);
1164 
1165  if(fIsMC > 0){
1166  fHistClusterdEtadPtTrueMatched = new TH2F(Form("dEtaVsPt_TrueMatched %s",GetCutNumber().Data()),"dEtaVsPt_TrueMatched",nEtaBins,EtaRange[0],EtaRange[1],
1167  nBinsModuleECoarse, minClusterELog, maxClusterELog);
1168  fHistClusterdEtadPtTrueMatched->GetXaxis()->SetTitle("d#eta_{cl-track}");
1169  fHistClusterdEtadPtTrueMatched->GetYaxis()->SetTitle("p_{T} (GeV/c)");
1170  SetLogBinningYTH2(fHistClusterdEtadPtTrueMatched);
1171  fHistograms->Add(fHistClusterdEtadPtTrueMatched);
1172  fHistClusterdPhidPtPosTracksTrueMatched = new TH2F(Form("dPhiVsPt_posTracks_TrueMatched %s",GetCutNumber().Data()),"dPhiVsPt_posTracks_TrueMatched",
1173  2*nPhiBins,2*PhiRange[0],2*PhiRange[1], nBinsModuleECoarse, minClusterELog, maxClusterELog);
1174  fHistClusterdPhidPtPosTracksTrueMatched->GetXaxis()->SetTitle("d#varphi_{cl-track} (rad)");
1175  fHistClusterdPhidPtPosTracksTrueMatched->GetYaxis()->SetTitle("p_{T} (GeV/c)");
1176  SetLogBinningYTH2(fHistClusterdPhidPtPosTracksTrueMatched);
1177  fHistograms->Add(fHistClusterdPhidPtPosTracksTrueMatched);
1178  fHistClusterdPhidPtNegTracksTrueMatched = new TH2F(Form("dPhiVsPt_negTracks_TrueMatched %s",GetCutNumber().Data()),"dPhiVsPt_negTracks_TrueMatched",
1179  2*nPhiBins,2*PhiRange[0],2*PhiRange[1], nBinsModuleECoarse, minClusterELog, maxClusterELog);
1180  fHistClusterdPhidPtNegTracksTrueMatched->GetXaxis()->SetTitle("d#varphi_{cl-track} (rad)");
1181  fHistClusterdPhidPtNegTracksTrueMatched->GetYaxis()->SetTitle("p_{T} (GeV/c)");
1182  SetLogBinningYTH2(fHistClusterdPhidPtNegTracksTrueMatched);
1183  fHistograms->Add(fHistClusterdPhidPtNegTracksTrueMatched);
1184  }
1185 
1186  // QA histos for shower shape correlations
1187  fHistClusterM20M02BeforeQA = new TH2F(Form("M20VsM02_beforeClusterQA %s",GetCutNumber().Data()),"M20VsM02_beforeClusterQA",200,0,2.5,400,0,5);
1188  fHistClusterM20M02BeforeQA->GetXaxis()->SetTitle("#sigma_{short}^2");
1189  fHistClusterM20M02BeforeQA->GetYaxis()->SetTitle("#sigma_{long}^2");
1190  fHistograms->Add(fHistClusterM20M02BeforeQA);
1191  fHistClusterM20M02AfterQA = new TH2F(Form("M20VsM02_afterClusterQA %s",GetCutNumber().Data()),"M20VsM02_afterClusterQA",200,0,2.5,400,0,5);
1192  fHistClusterM20M02AfterQA->GetXaxis()->SetTitle("#sigma_{short}^2");
1193  fHistClusterM20M02AfterQA->GetYaxis()->SetTitle("#sigma_{long}^2");
1194  fHistograms->Add(fHistClusterM20M02AfterQA);
1195 
1196  //----------------
1197  if(fIsMC > 1){
1198  fHistClusterRBeforeQA->Sumw2();
1199  fHistClusterRAfterQA->Sumw2();
1200  fHistDistanceTrackToClusterBeforeQA->Sumw2();
1201  fHistDistanceTrackToClusterAfterQA->Sumw2();
1202  fHistClusterdEtadPhiPosTracksBeforeQA->Sumw2();
1203  fHistClusterdEtadPhiNegTracksBeforeQA->Sumw2();
1204  fHistClusterdEtadPhiPosTracksAfterQA->Sumw2();
1205  fHistClusterdEtadPhiNegTracksAfterQA->Sumw2();
1206  fHistClusterdEtadPhiPosTracksP_000_075BeforeQA->Sumw2();
1207  fHistClusterdEtadPhiPosTracksP_075_125BeforeQA->Sumw2();
1208  fHistClusterdEtadPhiPosTracksP_125_999BeforeQA->Sumw2();
1209  fHistClusterdEtadPhiNegTracksP_000_075BeforeQA->Sumw2();
1210  fHistClusterdEtadPhiNegTracksP_075_125BeforeQA->Sumw2();
1211  fHistClusterdEtadPhiNegTracksP_125_999BeforeQA->Sumw2();
1212  fHistClusterdEtadPtBeforeQA->Sumw2();
1213  fHistClusterdEtadPtAfterQA->Sumw2();
1214  fHistClusterdPhidPtPosTracksBeforeQA->Sumw2();
1215  fHistClusterdPhidPtNegTracksBeforeQA->Sumw2();
1216  fHistClusterdPhidPtAfterQA->Sumw2();
1217  fHistClusterdEtadPtTrueMatched->Sumw2();
1218  fHistClusterdPhidPtPosTracksTrueMatched->Sumw2();
1219  fHistClusterdPhidPtNegTracksTrueMatched->Sumw2();
1220  fHistClusterM20M02BeforeQA->Sumw2();
1221  fHistClusterM20M02AfterQA->Sumw2();
1222  }
1223  //----------------
1224  }
1225  }
1226  if( fUseDistTrackToCluster && fIsMC && (fExtendedMatchAndQA == 1 || fExtendedMatchAndQA == 3 || fExtendedMatchAndQA == 5 )){
1227  // TM efficiency histograms
1228  const Int_t nEmcalEtaBins = 96;
1229  const Int_t nEmcalPhiBins = 124;
1230  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};
1231  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};
1232 
1233  fHistClusterElecEtaPhiBeforeTM_00_20 = new TH2F(Form("ElecEtaPhiBeforeTM_clusterE<20_Histo %s",GetCutNumber().Data()),"ElecEtaPhiBeforeTM_clusterE<20_Histo",nEmcalPhiBins,EmcalPhiBins,nEmcalEtaBins,EmcalEtaBins);
1234  fHistClusterElecEtaPhiBeforeTM_00_20->GetXaxis()->SetTitle("#varphi (rad)");
1235  fHistClusterElecEtaPhiBeforeTM_00_20->GetYaxis()->SetTitle("#eta");
1236  fHistograms->Add(fHistClusterElecEtaPhiBeforeTM_00_20);
1237  fHistClusterElecEtaPhiBeforeTM_20_50 = new TH2F(Form("ElecEtaPhiBeforeTM_20<clusterE<50_Histo %s",GetCutNumber().Data()),"ElecEtaPhiBeforeTM_20<clusterE<50_Histo",nEmcalPhiBins,EmcalPhiBins,nEmcalEtaBins,EmcalEtaBins);
1238  fHistClusterElecEtaPhiBeforeTM_20_50->GetXaxis()->SetTitle("#varphi (rad)");
1239  fHistClusterElecEtaPhiBeforeTM_20_50->GetYaxis()->SetTitle("#eta");
1240  fHistograms->Add(fHistClusterElecEtaPhiBeforeTM_20_50);
1241  fHistClusterElecEtaPhiBeforeTM_50_00 = new TH2F(Form("ElecEtaPhiBeforeTM_50<clusterE_Histo %s",GetCutNumber().Data()),"ElecEtaPhiBeforeTM_50<clusterE_Histo",nEmcalPhiBins,EmcalPhiBins,nEmcalEtaBins,EmcalEtaBins);
1242  fHistClusterElecEtaPhiBeforeTM_50_00->GetXaxis()->SetTitle("#varphi (rad)");
1243  fHistClusterElecEtaPhiBeforeTM_50_00->GetYaxis()->SetTitle("#eta");
1244  fHistograms->Add(fHistClusterElecEtaPhiBeforeTM_50_00);
1245 
1246  fHistClusterElecEtaPhiAfterTM_00_20 = new TH2F(Form("ElecEtaPhiAfterTM_clusterE<20_Histo %s",GetCutNumber().Data()),"ElecEtaPhiAfterTM_clusterE<20_Histo",nEmcalPhiBins,EmcalPhiBins,nEmcalEtaBins,EmcalEtaBins);
1247  fHistClusterElecEtaPhiAfterTM_00_20->GetXaxis()->SetTitle("#varphi (rad)");
1248  fHistClusterElecEtaPhiAfterTM_00_20->GetYaxis()->SetTitle("#eta");
1249  fHistograms->Add(fHistClusterElecEtaPhiAfterTM_00_20);
1250  fHistClusterElecEtaPhiAfterTM_20_50 = new TH2F(Form("ElecEtaPhiAfterTM_20<clusterE<50_Histo %s",GetCutNumber().Data()),"ElecEtaPhiAfterTM_20<clusterE<50_Histo",nEmcalPhiBins,EmcalPhiBins,nEmcalEtaBins,EmcalEtaBins);
1251  fHistClusterElecEtaPhiAfterTM_20_50->GetXaxis()->SetTitle("#varphi (rad)");
1252  fHistClusterElecEtaPhiAfterTM_20_50->GetYaxis()->SetTitle("#eta");
1253  fHistograms->Add(fHistClusterElecEtaPhiAfterTM_20_50);
1254  fHistClusterElecEtaPhiAfterTM_50_00 = new TH2F(Form("ElecEtaPhiAfterTM_50<clusterE_Histo %s",GetCutNumber().Data()),"ElecEtaPhiAfterTM_50<clusterE_Histo",nEmcalPhiBins,EmcalPhiBins,nEmcalEtaBins,EmcalEtaBins);
1255  fHistClusterElecEtaPhiAfterTM_50_00->GetXaxis()->SetTitle("#varphi (rad)");
1256  fHistClusterElecEtaPhiAfterTM_50_00->GetYaxis()->SetTitle("#eta");
1257  fHistograms->Add(fHistClusterElecEtaPhiAfterTM_50_00);
1258 
1259  fHistClusterTMEffiInput = new TH2F(Form("TMEffiInputHisto %s",GetCutNumber().Data()),"TMEffiInputHisto",nBinsClusterEFine, minClusterELog, maxClusterELog, 22, -0.5, 21.5);
1260  SetLogBinningXTH2(fHistClusterTMEffiInput);
1261  fHistClusterTMEffiInput->GetYaxis()->SetBinLabel(1,"All cl");
1262  fHistClusterTMEffiInput->GetYaxis()->SetBinLabel(2,"Ch cl");
1263  fHistClusterTMEffiInput->GetYaxis()->SetBinLabel(3,"Ne cl");
1264  fHistClusterTMEffiInput->GetYaxis()->SetBinLabel(4,"Ne cl sub ch");
1265  fHistClusterTMEffiInput->GetYaxis()->SetBinLabel(5,"Ga cl");
1266  fHistClusterTMEffiInput->GetYaxis()->SetBinLabel(6,"Ga cl sub ch");
1267  fHistClusterTMEffiInput->GetYaxis()->SetBinLabel(7,"conv cl");
1268  fHistClusterTMEffiInput->GetYaxis()->SetBinLabel(8,"Ch cl prim");
1269  fHistClusterTMEffiInput->GetYaxis()->SetBinLabel(9,"El cl");
1270  fHistClusterTMEffiInput->GetYaxis()->SetBinLabel(10,"All cl match");
1271  fHistClusterTMEffiInput->GetYaxis()->SetBinLabel(11,"Ch cl match");
1272  fHistClusterTMEffiInput->GetYaxis()->SetBinLabel(12,"Ch cl match w lead");
1273  fHistClusterTMEffiInput->GetYaxis()->SetBinLabel(13,"Ne cl match");
1274  fHistClusterTMEffiInput->GetYaxis()->SetBinLabel(14,"Ne cl sub ch match");
1275  fHistClusterTMEffiInput->GetYaxis()->SetBinLabel(15,"Ga cl match");
1276  fHistClusterTMEffiInput->GetYaxis()->SetBinLabel(16,"Ga cl sub ch match");
1277  fHistClusterTMEffiInput->GetYaxis()->SetBinLabel(17,"conv cl match");
1278  fHistClusterTMEffiInput->GetYaxis()->SetBinLabel(18,"conv cl match w lead");
1279  fHistClusterTMEffiInput->GetYaxis()->SetBinLabel(19,"Ch cl match");
1280  fHistClusterTMEffiInput->GetYaxis()->SetBinLabel(20,"Ch cl match w lead");
1281  fHistClusterTMEffiInput->GetYaxis()->SetBinLabel(21,"El cl match");
1282  fHistClusterTMEffiInput->GetYaxis()->SetBinLabel(22,"El cl match w lead");
1283  fHistClusterTMEffiInput->GetXaxis()->SetTitle("#it{E}_{cl} (GeV)");
1284  fHistograms->Add(fHistClusterTMEffiInput);
1285 
1286  fHistClusterEvsTrackECharged = new TH2F(Form("ClusterE_TrackE_ChargedCluster %s",GetCutNumber().Data()),"ClusterE TrackE ChargedCluster",
1287  nBinsClusterEFine, minClusterELog, maxClusterELog, nBinsClusterEFine, minClusterELog, maxClusterELog);
1288  SetLogBinningXTH2(fHistClusterEvsTrackECharged);
1289  SetLogBinningYTH2(fHistClusterEvsTrackECharged);
1290  fHistClusterEvsTrackECharged->GetXaxis()->SetTitle("#it{E}_{cl} (GeV)");
1291  fHistClusterEvsTrackECharged->GetYaxis()->SetTitle("#it{E}_{tr} (GeV)");
1292  fHistograms->Add(fHistClusterEvsTrackECharged);
1293  fHistClusterEvsTrackEChargedLead = new TH2F(Form("ClusterE_TrackE_ChargedCluster_LeadMatched %s",GetCutNumber().Data()),"ClusterE TrackE ChargedCluster LeadMatched",
1294  nBinsClusterEFine, minClusterELog, maxClusterELog, nBinsClusterEFine, minClusterELog, maxClusterELog);
1295  SetLogBinningXTH2(fHistClusterEvsTrackEChargedLead);
1296  SetLogBinningYTH2(fHistClusterEvsTrackEChargedLead);
1297  fHistClusterEvsTrackEChargedLead->GetXaxis()->SetTitle("#it{E}_{cl} (GeV)");
1298  fHistClusterEvsTrackEChargedLead->GetYaxis()->SetTitle("#it{E}_{tr} (GeV)");
1299  fHistograms->Add(fHistClusterEvsTrackEChargedLead);
1300  fHistClusterEvsTrackENeutral = new TH2F(Form("ClusterE_TrackE_NeutralCluster %s",GetCutNumber().Data()),"ClusterE TrackE NeutralCluster",
1301  nBinsClusterEFine, minClusterELog, maxClusterELog, nBinsClusterEFine, minClusterELog, maxClusterELog);
1302  SetLogBinningXTH2(fHistClusterEvsTrackENeutral);
1303  SetLogBinningYTH2(fHistClusterEvsTrackENeutral);
1304  fHistClusterEvsTrackENeutral->GetXaxis()->SetTitle("#it{E}_{cl} (GeV)");
1305  fHistClusterEvsTrackENeutral->GetYaxis()->SetTitle("#it{E}_{tr} (GeV)");
1306  fHistograms->Add(fHistClusterEvsTrackENeutral);
1307  fHistClusterEvsTrackENeutralSubCharged = new TH2F(Form("ClusterE_TrackE_NeutralClusterSubCharged %s",GetCutNumber().Data()),"ClusterE TrackE NeutralCluster SubCharged",
1308  nBinsClusterEFine, minClusterELog, maxClusterELog, nBinsClusterEFine, minClusterELog, maxClusterELog);
1309  SetLogBinningXTH2(fHistClusterEvsTrackENeutralSubCharged);
1310  SetLogBinningYTH2(fHistClusterEvsTrackENeutralSubCharged);
1311  fHistClusterEvsTrackENeutralSubCharged->GetXaxis()->SetTitle("#it{E}_{cl} (GeV)");
1312  fHistClusterEvsTrackENeutralSubCharged->GetYaxis()->SetTitle("#it{E}_{tr} (GeV)");
1313  fHistograms->Add(fHistClusterEvsTrackENeutralSubCharged);
1314  fHistClusterEvsTrackEGamma = new TH2F(Form("ClusterE_TrackE_GammaCluster %s",GetCutNumber().Data()),"ClusterE TrackE GammaCluster",
1315  nBinsClusterEFine, minClusterELog, maxClusterELog, nBinsClusterEFine, minClusterELog, maxClusterELog);
1316  SetLogBinningXTH2(fHistClusterEvsTrackEGamma);
1317  SetLogBinningYTH2(fHistClusterEvsTrackEGamma);
1318  fHistClusterEvsTrackEGamma->GetXaxis()->SetTitle("#it{E}_{cl} (GeV)");
1319  fHistClusterEvsTrackEGamma->GetYaxis()->SetTitle("#it{E}_{tr} (GeV)");
1320  fHistograms->Add(fHistClusterEvsTrackEGamma);
1321  fHistClusterEvsTrackEGammaSubCharged = new TH2F(Form("ClusterE_TrackE_GammaClusterSubCharged %s",GetCutNumber().Data()),"ClusterE TrackE GammaCluster SubCharged",
1322  nBinsClusterEFine, minClusterELog, maxClusterELog, nBinsClusterEFine, minClusterELog, maxClusterELog);
1323  SetLogBinningXTH2(fHistClusterEvsTrackEGammaSubCharged);
1324  SetLogBinningYTH2(fHistClusterEvsTrackEGammaSubCharged);
1325  fHistClusterEvsTrackEGammaSubCharged->GetXaxis()->SetTitle("#it{E}_{cl} (GeV)");
1326  fHistClusterEvsTrackEGammaSubCharged->GetYaxis()->SetTitle("#it{E}_{tr} (GeV)");
1327  fHistograms->Add(fHistClusterEvsTrackEGammaSubCharged);
1328  fHistClusterEvsTrackEConv = new TH2F(Form("ClusterE_TrackE_ConvCluster %s",GetCutNumber().Data()),"ClusterE TrackE ConvCluster",
1329  nBinsClusterEFine, minClusterELog, maxClusterELog, nBinsClusterEFine, minClusterELog, maxClusterELog);
1330  SetLogBinningXTH2(fHistClusterEvsTrackEConv);
1331  SetLogBinningYTH2(fHistClusterEvsTrackEConv);
1332  fHistClusterEvsTrackEConv->GetXaxis()->SetTitle("#it{E}_{cl} (GeV)");
1333  fHistClusterEvsTrackEConv->GetYaxis()->SetTitle("#it{E}_{tr} (GeV)");
1334  fHistograms->Add(fHistClusterEvsTrackEConv);
1335 
1336  fHistClusterENMatchesNeutral = new TH2F(Form("ClusterE_NMatches_NeutralCluster %s",GetCutNumber().Data()),"ClusterE NMatches NeutralCluster",
1337  nBinsClusterEFine, minClusterELog, maxClusterELog, 20, -0.5, 19.5);
1338  SetLogBinningXTH2(fHistClusterENMatchesNeutral);
1339  fHistClusterENMatchesNeutral->GetXaxis()->SetTitle("#it{E}_{cl} (GeV)");
1340  fHistClusterENMatchesNeutral->GetYaxis()->SetTitle("#it{N}_{matches}");
1341  fHistograms->Add(fHistClusterENMatchesNeutral);
1342  fHistClusterENMatchesCharged = new TH2F(Form("ClusterE_NMatches_ChargedCluster %s",GetCutNumber().Data()),"ClusterE NMatches ChargedCluster",
1343  nBinsClusterEFine, minClusterELog, maxClusterELog, 20, -0.5, 19.5);
1344  SetLogBinningXTH2(fHistClusterENMatchesCharged);
1345  fHistClusterENMatchesCharged->GetXaxis()->SetTitle("#it{E}_{cl} (GeV)");
1346  fHistClusterENMatchesCharged->GetYaxis()->SetTitle("#it{N}_{matches}");
1347  fHistograms->Add(fHistClusterENMatchesCharged);
1348 
1349  fHistClusterEvsTrackEPrimaryButNoElec = new TH2F(Form("ClusterE_TrackE_ChargedClusterNoElec %s",GetCutNumber().Data()),"ClusterE TrackE ChargedClusterNoElec",
1350  nBinsClusterEFine, minClusterELog, maxClusterELog, nBinsClusterEFine, minClusterELog, maxClusterELog);
1351  SetLogBinningXTH2(fHistClusterEvsTrackEPrimaryButNoElec);
1352  SetLogBinningYTH2(fHistClusterEvsTrackEPrimaryButNoElec);
1353  fHistClusterEvsTrackEPrimaryButNoElec->GetXaxis()->SetTitle("#it{E}_{cl} (GeV)");
1354  fHistClusterEvsTrackEPrimaryButNoElec->GetYaxis()->SetTitle("#it{E}_{tr} (GeV)");
1355  fHistograms->Add(fHistClusterEvsTrackEPrimaryButNoElec);
1356  fHistClusterEvsTrackSumEPrimaryButNoElec = new TH2F(Form("ClusterE_TrackSumE_ChargedClusterNoElec %s",GetCutNumber().Data()),"ClusterE TrackSumE ChargedClusterNoElec",
1357  nBinsClusterEFine, minClusterELog, maxClusterELog, nBinsClusterEFine, minClusterELog, maxClusterELog);
1358  SetLogBinningXTH2(fHistClusterEvsTrackSumEPrimaryButNoElec);
1359  SetLogBinningYTH2(fHistClusterEvsTrackSumEPrimaryButNoElec);
1360  fHistClusterEvsTrackSumEPrimaryButNoElec->GetXaxis()->SetTitle("#it{E}_{cl} (GeV)");
1361  fHistClusterEvsTrackSumEPrimaryButNoElec->GetYaxis()->SetTitle("#sum#it{E}_{tr} (GeV)");
1362  fHistograms->Add(fHistClusterEvsTrackSumEPrimaryButNoElec);
1363 
1364  if(fIsMC > 1){
1365  fHistClusterTMEffiInput->Sumw2();
1366  fHistClusterElecEtaPhiBeforeTM_00_20->Sumw2();
1367  fHistClusterElecEtaPhiBeforeTM_20_50->Sumw2();
1368  fHistClusterElecEtaPhiBeforeTM_50_00->Sumw2();
1369  fHistClusterElecEtaPhiAfterTM_00_20->Sumw2();
1370  fHistClusterElecEtaPhiAfterTM_20_50->Sumw2();
1371  fHistClusterElecEtaPhiAfterTM_50_00->Sumw2();
1372  fHistClusterEvsTrackECharged->Sumw2();
1373  fHistClusterEvsTrackEChargedLead->Sumw2();
1374  fHistClusterEvsTrackENeutral->Sumw2();
1375  fHistClusterEvsTrackENeutralSubCharged->Sumw2();
1376  fHistClusterEvsTrackEGamma->Sumw2();
1377  fHistClusterEvsTrackEGammaSubCharged->Sumw2();
1378  fHistClusterEvsTrackEConv->Sumw2();
1379  fHistClusterENMatchesNeutral->Sumw2();
1380  fHistClusterENMatchesCharged->Sumw2();
1381  fHistClusterEvsTrackEPrimaryButNoElec->Sumw2();
1382  fHistClusterEvsTrackSumEPrimaryButNoElec->Sumw2();
1383  }
1384  }
1385 
1386  if (fDoExoticsQA){
1387  if( fClusterType == 1 ){ //EMCAL
1388  const Int_t nEmcalEtaBins = 96;
1389  const Int_t nEmcalPhiBins = 124;
1390  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};
1391  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};
1392 
1393  fHistClusterEtavsPhiExotics = new TH2F(Form("EtaPhi_Exotics %s",GetCutNumber().Data()),"EtaPhi_Exotics",nEmcalPhiBins, EmcalPhiBins, nEmcalEtaBins, EmcalEtaBins);
1394  fHistClusterEtavsPhiExotics->GetXaxis()->SetTitle("#eta");
1395  fHistClusterEtavsPhiExotics->GetYaxis()->SetTitle("#varphi (rad)");
1396  fHistograms->Add(fHistClusterEtavsPhiExotics);
1397  } else if( fClusterType == 2 ){ //PHOS
1398  const Int_t nPhosEtaBins = 56;
1399  const Int_t nPhosPhiBins = 192;
1400  const Float_t PhosEtaRange[2] = {-0.16, 0.16};
1401  const Float_t PhosPhiRange[2] = {4.5, 5.6};
1402 
1403  fHistClusterEtavsPhiExotics = new TH2F(Form("EtaPhi_Exotics %s",GetCutNumber().Data()),"EtaPhi_Exotics",nPhosPhiBins, PhosPhiRange[0], PhosPhiRange[1],
1404  nPhosEtaBins, PhosEtaRange[0], PhosEtaRange[1]);
1405  fHistClusterEtavsPhiExotics->GetXaxis()->SetTitle("#eta");
1406  fHistClusterEtavsPhiExotics->GetYaxis()->SetTitle("#varphi (rad)");
1407  fHistograms->Add(fHistClusterEtavsPhiExotics);
1408  }
1409  if( fClusterType == 3 ){ //DCAL
1410  const Int_t nDcalEtaBins = 96;
1411  const Int_t nDcalPhiBins = 124;
1412 
1413  fHistClusterEtavsPhiExotics = new TH2F(Form("EtaPhi_Exotics %s",GetCutNumber().Data()),"EtaPhi_Exotics",nDcalPhiBins,4.5,5.7,nDcalEtaBins,-0.66687,0.66465);
1414  fHistClusterEtavsPhiExotics->GetXaxis()->SetTitle("#eta");
1415  fHistClusterEtavsPhiExotics->GetYaxis()->SetTitle("#varphi (rad)");
1416  fHistograms->Add(fHistClusterEtavsPhiExotics);
1417  }
1418  fHistClusterEM02Exotics = new TH2F(Form("EVsM02_Exotics %s",GetCutNumber().Data()),"EVsM02_afterClusterQA",nBinsClusterEFine, minClusterELog, maxClusterELog,200,0,2);
1419  fHistClusterEM02Exotics->GetXaxis()->SetTitle("E_{cl,exotics} (GeV)");
1420  fHistClusterEM02Exotics->GetYaxis()->SetTitle("#sigma_{long}^2");
1421  SetLogBinningXTH2(fHistClusterEM02Exotics);
1422  fHistograms->Add(fHistClusterEM02Exotics);
1423  fHistClusterEnergyvsNCellsExotics = new TH2F(Form("ClusterEnergyVsNCells_Exotics %s",GetCutNumber().Data()),"ClusterEnergyVsNCells_Exotics",
1424  nBinsClusterECoarse, minClusterELog, maxClusterELog, 50, 0, 50);
1425  fHistClusterEnergyvsNCellsExotics->GetXaxis()->SetTitle("E_{cl,exotics} (GeV)");
1426  fHistClusterEnergyvsNCellsExotics->GetYaxis()->SetTitle("N_{LM}");
1427  SetLogBinningXTH2(fHistClusterEnergyvsNCellsExotics);
1428  fHistograms->Add(fHistClusterEnergyvsNCellsExotics);
1429  fHistClusterEEstarExotics = new TH2F(Form("ClusterEnergyVsEnergystar_Exotics %s",GetCutNumber().Data()),"ClusterEnergyVsEnergystar_Exotics",
1430  nBinsClusterECoarse, minClusterELog, maxClusterELog, nBinsClusterECoarse, minClusterELog, maxClusterELog);
1431  fHistClusterEEstarExotics->GetXaxis()->SetTitle("E_{cl,exotics} (GeV)");
1432  fHistClusterEEstarExotics->GetYaxis()->SetTitle("E_{cl,#star} (GeV)");
1433  SetLogBinningXTH2(fHistClusterEEstarExotics);
1434  SetLogBinningYTH2(fHistClusterEEstarExotics);
1435  fHistograms->Add(fHistClusterEEstarExotics);
1436 
1437  if (fIsMC > 1){
1438  fHistClusterEtavsPhiExotics->Sumw2();
1439  fHistClusterEM02Exotics->Sumw2();
1440  fHistClusterEnergyvsNCellsExotics->Sumw2();
1441  fHistClusterEEstarExotics->Sumw2();
1442  }
1443  }
1444 
1445  fVectorMatchedClusterIDs.clear();
1446 
1447  TH1::AddDirectory(kTRUE);
1448  return;
1449 }
1450 
1451 //________________________________________________________________________
1452 void AliCaloPhotonCuts::InitializeEMCAL(AliVEvent *event){
1453 
1454  if (fClusterType == 1 || fClusterType == 3){
1455  AliTender* alitender=0x0;
1456  AliEmcalTenderTask* emcaltender=0x0;
1457  AliEmcalCorrectionTask* emcalCorrTask=0x0;
1458  if(event->IsA()==AliESDEvent::Class()){
1459  alitender = (AliTender*) AliAnalysisManager::GetAnalysisManager()->GetTopTasks()->FindObject("AliTender");
1460  if(!alitender){
1461  emcalCorrTask = (AliEmcalCorrectionTask*) AliAnalysisManager::GetAnalysisManager()->GetTopTasks()->FindObject("AliEmcalCorrectionTask_defaultSetting");
1462  }
1463  } else if( event->IsA()==AliAODEvent::Class()){
1464  emcaltender = (AliEmcalTenderTask*) AliAnalysisManager::GetAnalysisManager()->GetTopTasks()->FindObject("AliEmcalTenderTask");
1465  if(!emcaltender)
1466  emcalCorrTask = (AliEmcalCorrectionTask*) AliAnalysisManager::GetAnalysisManager()->GetTopTasks()->FindObject("AliEmcalCorrectionTask_defaultSetting");
1467  }
1468  if(alitender){
1469  TIter next(alitender->GetSupplies());
1470  AliTenderSupply *supply;
1471  while ((supply=(AliTenderSupply*)next())) if(supply->IsA()==AliEMCALTenderSupply::Class()) break;
1472  fEMCALRecUtils = ((AliEMCALTenderSupply*)supply)->GetRecoUtils();
1473  fEMCALBadChannelsMap = fEMCALRecUtils->GetEMCALBadChannelStatusMapArray();
1474  } else if(emcaltender){
1475  fEMCALRecUtils = ((AliEMCALTenderSupply*)emcaltender->GetEMCALTenderSupply())->GetRecoUtils();
1476  fEMCALBadChannelsMap = fEMCALRecUtils->GetEMCALBadChannelStatusMapArray();
1477  } else if(emcalCorrTask){
1478  AliEmcalCorrectionComponent * emcalCorrComponent = emcalCorrTask->GetCorrectionComponent("AliEmcalCorrectionCellBadChannel_defaultSetting");
1479  if(emcalCorrComponent){
1480  fEMCALRecUtils = emcalCorrComponent->GetRecoUtils();
1481  fEMCALBadChannelsMap = fEMCALRecUtils->GetEMCALBadChannelStatusMapArray();
1482  }
1483  }
1484  if (fEMCALRecUtils) fEMCALInitialized = kTRUE;
1485 
1486  fGeomEMCAL = AliEMCALGeometry::GetInstance();
1487  if(!fGeomEMCAL){ AliFatal("EMCal geometry not initialized!");}
1488 
1489  //retrieve pointer to trackMatcher Instance
1490  if(fUseDistTrackToCluster) fCaloTrackMatcher = (AliCaloTrackMatcher*)AliAnalysisManager::GetAnalysisManager()->GetTask(fCaloTrackMatcherName.Data());
1491  if(!fCaloTrackMatcher && fUseDistTrackToCluster ){ AliFatal("CaloTrackMatcher instance could not be initialized!");}
1492 
1493  if(!fDoLightOutput){
1494  Int_t nMaxCellsEMCAL = fNMaxEMCalModules*48*24;
1495  Int_t nMinCellsDCAL = 12288;
1496  Int_t nMaxCellsDCAL = nMinCellsDCAL+fNMaxDCalModules*32*24;
1497  Int_t nMaxCells = 0;
1498  Int_t nMinCells = 0;
1499  if(fClusterType == 1){
1500  nMaxCells = nMaxCellsEMCAL;
1501  nMinCells = 0;
1502  } else if(fClusterType == 3){
1503  nMaxCells = nMaxCellsDCAL;
1504  nMinCells = nMinCellsDCAL;
1505  }
1506 
1507 
1508  Int_t imod = -1;Int_t iTower = -1, iIphi = -1, iIeta = -1;
1509  Int_t icol = -1;Int_t irow = -1;
1510 
1511  for(Int_t iCell=nMinCells;iCell<nMaxCells;iCell++){
1512  fGeomEMCAL->GetCellIndex(iCell,imod,iTower,iIphi,iIeta);
1513  if (fEMCALBadChannelsMap->GetEntries() <= imod) continue;
1514  fGeomEMCAL->GetCellPhiEtaIndexInSModule(imod,iTower,iIphi,iIeta,irow,icol);
1515  Int_t iBadCell = (Int_t) ((TH2I*)fEMCALBadChannelsMap->At(imod))->GetBinContent(icol,irow);
1516  if(iBadCell > 0) fBadChannels->Fill(iCell,1);
1517  else fBadChannels->Fill(iCell,0);
1518  }
1519  }
1520  }
1521  return;
1522 }
1523 
1524 //________________________________________________________________________
1525 void AliCaloPhotonCuts::InitializePHOS (AliVEvent *event){
1526 
1527  if (fClusterType == 2){
1528  if(!fDoLightOutput){
1529  fGeomPHOS = AliPHOSGeometry::GetInstance();
1530  if(!fGeomPHOS) AliFatal("PHOS geometry not initialized!");
1531  Int_t nModules = fGeomPHOS->GetNModules();
1532  //cout << nModules << endl;
1533 
1534  fPHOSBadChannelsMap = new TH2I*[nModules];
1535  AliPHOSTenderTask* aliphostender = (AliPHOSTenderTask*) AliAnalysisManager::GetAnalysisManager()->GetTopTasks()->FindObject("PHOSTenderTask");
1536  AliPHOSTenderSupply *PHOSSupply =((AliPHOSTenderSupply*) aliphostender->GetPHOSTenderSupply()) ;
1537 
1538  if(!PHOSSupply){
1539  AliError(Form("Can not find PHOSTenderSupply in run %d. No bad channel map could be found for QA!\n",event->GetRunNumber())) ;
1540  for(Int_t mod=0;mod<nModules;mod++) fPHOSBadChannelsMap[mod] = NULL;
1541  }else{
1542  AliInfo("Setting PHOS bad map from PHOSSupply \n") ;
1543  for(Int_t mod=0;mod<nModules;mod++){
1544  TH2I * h = (TH2I*)PHOSSupply->GetPHOSBadChannelStatusMap(mod);
1545  if(h){
1546  fPHOSBadChannelsMap[mod] = new TH2I(*h);
1547  AliInfo(Form("using bad map for module %d with nch=%f\n",mod,h->Integral()));
1548  }
1549  else fPHOSBadChannelsMap[mod] = NULL;
1550  }
1551 
1552  Int_t nMaxCellsPHOS = (fNMaxPHOSModules*56*64);
1553  Int_t relid[4];
1554 
1555  for(Int_t iCell=0;iCell<nMaxCellsPHOS;iCell++){
1556  fGeomPHOS->AbsToRelNumbering(iCell,relid);
1557  //cout << relid[0] << ", " << relid[1] << ", " << relid[2] << ", " << relid[3] << ", " << __LINE__ << endl;
1558  if(relid[1]!=0) AliFatal("PHOS CPV in PHOS cell array?");
1559  if(relid[0]>=nModules || relid[0]<0 || !fPHOSBadChannelsMap[relid[0]]) continue;
1560  Int_t iBadCell = (Int_t) ((TH2I*)fPHOSBadChannelsMap[relid[0]])->GetBinContent(relid[2],relid[3]);
1561  if(iBadCell > 0) fBadChannels->Fill(iCell,1);
1562  else fBadChannels->Fill(iCell,0);
1563  }
1564  }
1565  }
1566 
1567  //retrieve pointer to trackMatcher Instance
1568  if(fUseDistTrackToCluster) fCaloTrackMatcher = (AliCaloTrackMatcher*)AliAnalysisManager::GetAnalysisManager()->GetTask(fCaloTrackMatcherName.Data());
1569  if(!fCaloTrackMatcher && fUseDistTrackToCluster){ AliFatal("CaloTrackMatcher instance could not be initialized!");}
1570 
1571  fPHOSInitialized = kTRUE;
1572  fPHOSCurrentRun = event->GetRunNumber();
1573  }
1574  return;
1575 }
1576 
1577 //________________________________________________________________________
1578 Bool_t AliCaloPhotonCuts::ClusterIsSelectedMC(TParticle *particle,AliMCEvent *mcEvent){
1579  // MonteCarlo Photon Selection
1580 
1581  if(!mcEvent)return kFALSE;
1582  if(!particle) return kFALSE;
1583 
1584  if (particle->GetPdgCode() == 22){
1585 
1586  if ( particle->Eta() < fMinEtaCut || particle->Eta() > fMaxEtaCut ) return kFALSE;
1587  if ( particle->Phi() < fMinPhiCut || particle->Phi() > fMaxPhiCut ) return kFALSE;
1588  if ( fClusterType == 3 && particle->Eta() < fMaxEtaInnerEdge && particle->Eta() > fMinEtaInnerEdge ) return kFALSE;
1589 
1590  if(particle->GetMother(0) >-1 && mcEvent->Particle(particle->GetMother(0))->GetPdgCode() == 22){
1591  return kFALSE;// no photon as mothers!
1592  }
1593  return kTRUE;
1594  }
1595  return kFALSE;
1596 }
1597 
1598 //________________________________________________________________________
1599 Bool_t AliCaloPhotonCuts::ClusterIsSelectedElecMC(TParticle *particle,AliMCEvent *mcEvent){
1600  // MonteCarlo Photon Selection
1601 
1602  if(!mcEvent)return kFALSE;
1603  if(!particle) return kFALSE;
1604 
1605  if (TMath::Abs(particle->GetPdgCode()) == 11){
1606 
1607  if ( particle->Eta() < fMinEtaCut || particle->Eta() > fMaxEtaCut ) return kFALSE;
1608  if ( particle->Phi() < fMinPhiCut || particle->Phi() > fMaxPhiCut ) return kFALSE;
1609  if ( fClusterType == 3 && particle->Eta() < fMaxEtaInnerEdge && particle->Eta() > fMinEtaInnerEdge ) return kFALSE;
1610 
1611  if(particle->GetMother(0) >-1 && mcEvent->Particle(particle->GetMother(0))->GetPdgCode() == 11){
1612  return kFALSE;// no photon as mothers!
1613  }
1614  return kTRUE;
1615  }
1616  return kFALSE;
1617 }
1618 
1619 //________________________________________________________________________
1620 Bool_t AliCaloPhotonCuts::ClusterIsSelectedElecAODMC(AliAODMCParticle *particle,TClonesArray *aodmcArray){
1621  // MonteCarlo Photon Selection
1622 
1623  if(!aodmcArray)return kFALSE;
1624  if(!particle) return kFALSE;
1625 
1626  if (TMath::Abs(particle->GetPdgCode()) == 11){
1627 
1628  if ( particle->Eta() < fMinEtaCut || particle->Eta() > fMaxEtaCut ) return kFALSE;
1629  if ( particle->Phi() < fMinPhiCut || particle->Phi() > fMaxPhiCut ) return kFALSE;
1630  if ( fClusterType == 3 && particle->Eta() < fMaxEtaInnerEdge && particle->Eta() > fMinEtaInnerEdge ) return kFALSE;
1631  if(particle->GetMother() >-1 && (static_cast<AliAODMCParticle*>(aodmcArray->At(particle->GetMother())))->GetPdgCode() == 11){
1632  return kFALSE;// no photon as mothers!
1633  }
1634  return kTRUE;
1635  }
1636  return kFALSE;
1637 }
1638 
1639 //________________________________________________________________________
1640 Bool_t AliCaloPhotonCuts::ClusterIsSelectedAODMC(AliAODMCParticle *particle,TClonesArray *aodmcArray){
1641  // MonteCarlo Photon Selection
1642 
1643  if(!aodmcArray)return kFALSE;
1644  if(!particle) return kFALSE;
1645 
1646  if (particle->GetPdgCode() == 22){
1647  if ( particle->Eta() < fMinEtaCut || particle->Eta() > fMaxEtaCut ) return kFALSE;
1648  if ( particle->Phi() < fMinPhiCut || particle->Phi() > fMaxPhiCut ) return kFALSE;
1649  if ( fClusterType == 3 && particle->Eta() < fMaxEtaInnerEdge && particle->Eta() > fMinEtaInnerEdge ) return kFALSE;
1650  if(particle->GetMother() > -1 && (static_cast<AliAODMCParticle*>(aodmcArray->At(particle->GetMother())))->GetPdgCode() == 22){
1651  return kFALSE;// no photon as mothers!
1652  }
1653  return kTRUE;// return in case of accepted gamma
1654  }
1655  return kFALSE;
1656 }
1657 
1658 //________________________________________________________________________
1659 // This function selects the clusters based on their quality criteria
1660 //________________________________________________________________________
1661 Bool_t AliCaloPhotonCuts::ClusterQualityCuts(AliVCluster* cluster, AliVEvent *event, AliMCEvent* mcEvent, Int_t isMC, Double_t weight, Long_t clusterID)
1662 { // Specific Photon Cuts
1663 
1664  // Initialize EMCAL rec utils if not initialized
1665  if(!fEMCALInitialized && (fClusterType == 1 || fClusterType == 3)) InitializeEMCAL(event);
1666 
1667  fIsCurrentClusterAcceptedBeforeTM = kFALSE;
1668  fIsAcceptedForBasic = kFALSE;
1669  Int_t cutIndex = 0;
1670  if(fHistClusterIdentificationCuts)fHistClusterIdentificationCuts->Fill(cutIndex,cluster->E());
1671  cutIndex++;
1672 
1673  // cluster position defintion
1674  Float_t clusPos[3]={0,0,0};
1675  cluster->GetPosition(clusPos);
1676  TVector3 clusterVector(clusPos[0],clusPos[1],clusPos[2]);
1677  Double_t etaCluster = clusterVector.Eta();
1678  Double_t phiCluster = clusterVector.Phi();
1679  if (phiCluster < 0) phiCluster += 2*TMath::Pi();
1680 
1681 
1682  Int_t nLM = GetNumberOfLocalMaxima(cluster, event);
1683 // Int_t nLMGustavo = fEMCALCaloUtils->GetNumberOfLocalMaxima(cluster, event->GetEMCALCells()) ;
1684 // cout << "mine: " << nLM << "\t Gustavo: " << nLMGustavo << endl;
1685 
1686  // Fill Histos before Cuts
1687  if(fHistClusterTimevsEBeforeQA) fHistClusterTimevsEBeforeQA->Fill(cluster->GetTOF(), cluster->E(), weight);
1688  if(fHistEnergyOfClusterBeforeQA) fHistEnergyOfClusterBeforeQA->Fill(cluster->E(), weight);
1689  if(fHistNCellsBeforeQA) fHistNCellsBeforeQA->Fill(cluster->GetNCells(), weight);
1690  if(fHistM02BeforeQA) fHistM02BeforeQA->Fill(cluster->GetM02(), weight);
1691  if(fHistM20BeforeQA) fHistM20BeforeQA->Fill(cluster->GetM20(), weight);
1692  if(fHistDispersionBeforeQA) fHistDispersionBeforeQA->Fill(cluster->GetDispersion(), weight);
1693  if(fHistNLMBeforeQA) fHistNLMBeforeQA->Fill(nLM, weight);
1694 // if(fHistNLMAvsNLMBBeforeQA) fHistNLMAvsNLMBBeforeQA->Fill(nLM, nLMGustavo, weight);
1695  if(fHistClusterEM02BeforeQA) fHistClusterEM02BeforeQA->Fill(cluster->E(),cluster->GetM02(), weight);
1696 
1697  AliVCaloCells* cells = NULL;
1698  if(fExtendedMatchAndQA > 1){
1699  if(cluster->IsEMCAL()){ //EMCAL
1700  cells = event->GetEMCALCells();
1701  }else if(cluster->IsPHOS()){ //PHOS
1702  cells = event->GetPHOSCells();
1703  }
1704  if(fHistClusterIncludedCellsBeforeQA){
1705  Int_t nCellCluster = cluster->GetNCells();
1706  for(Int_t iCell=0;iCell<nCellCluster;iCell++){
1707  fHistClusterIncludedCellsBeforeQA->Fill(cluster->GetCellAbsId(iCell));
1708  if(cluster->E()>0.) fHistClusterEnergyFracCellsBeforeQA->Fill(cluster->GetCellAbsId(iCell),cells->GetCellAmplitude(cluster->GetCellAbsId(iCell))/cluster->E());
1709  }
1710  }
1711  }
1712 
1713  // Check wether timing is ok
1714  if (fUseTimeDiff){
1715  if( (cluster->GetTOF() < fMinTimeDiff || cluster->GetTOF() > fMaxTimeDiff) && !(isMC>0)){
1716  if(fHistClusterIdentificationCuts)fHistClusterIdentificationCuts->Fill(cutIndex, cluster->E());//1
1717  return kFALSE;
1718  }
1719  }
1720  cutIndex++;//2, next cut
1721 
1722  // exotic cluster cut
1723  Float_t energyStar = 0;
1724  if(fUseExoticCluster && IsExoticCluster(cluster, event, energyStar)){
1725  if(fHistClusterIdentificationCuts)fHistClusterIdentificationCuts->Fill(cutIndex, cluster->E());//3
1726  if (fDoExoticsQA){
1727  // replay cuts
1728  Bool_t failed = kFALSE;
1729  Bool_t failedM02 = kFALSE;
1730  if (fUseMinEnergy)
1731  if(cluster->E() < fMinEnergy)
1732  failed = kTRUE;
1733  if (fUseNCells)
1734  if(cluster->GetNCells() < fMinNCells)
1735  failed = kTRUE;
1736  if (fUseNLM)
1737  if( nLM < fMinNLM || nLM > fMaxNLM )
1738  failed = kTRUE;
1739  if (fUseM02 == 1){
1740  if( cluster->GetM02()< fMinM02 || cluster->GetM02() > fMaxM02 )
1741  failedM02 = kTRUE;
1742  } else if (fUseM02 ==2 ) {
1743  if( cluster->GetM02()< CalculateMinM02(fMinM02CutNr, cluster->E()) ||
1744  cluster->GetM02() > CalculateMaxM02(fMaxM02CutNr, cluster->E()) )
1745  failedM02 = kTRUE;
1746  }
1747  if (fUseM20)
1748  if( cluster->GetM20()< fMinM20 || cluster->GetM20() > fMaxM20 )
1749  failed = kTRUE;
1750  if (fUseDispersion)
1751  if( cluster->GetDispersion()> fMaxDispersion)
1752  failed = kTRUE;
1753  if (fVectorMatchedClusterIDs.size()>0 && fUseDistTrackToCluster)
1754  if( CheckClusterForTrackMatch(cluster) )
1755  failed = kTRUE;
1756  if ( !( failed || failedM02 ) ){
1757  if(fHistClusterEtavsPhiExotics) fHistClusterEtavsPhiExotics->Fill(phiCluster, etaCluster, weight);
1758  if(fHistClusterEnergyvsNCellsExotics) fHistClusterEnergyvsNCellsExotics->Fill(cluster->E(), cluster->GetNCells(), weight);
1759  if(fHistClusterEEstarExotics) fHistClusterEEstarExotics->Fill(cluster->E(),energyStar, weight);
1760  }
1761  if ( !failed ){
1762  if(fHistClusterEM02Exotics) fHistClusterEM02Exotics->Fill(cluster->E(), cluster->GetM02(), weight);
1763  }
1764  }
1765  return kFALSE;
1766  }
1767  cutIndex++;//3, next cut
1768 
1769  // minimum number of cells
1770  if (fUseNCells){
1771  if(cluster->GetNCells() < fMinNCells) {
1772  if(fHistClusterIdentificationCuts)fHistClusterIdentificationCuts->Fill(cutIndex, cluster->E());//5
1773  return kFALSE;
1774  }
1775  }
1776  cutIndex++;//4, next cut
1777 
1778  // NLM cut
1779  if (fUseNLM){
1780  if( nLM < fMinNLM || nLM > fMaxNLM ) {
1781  if(fHistClusterIdentificationCuts)fHistClusterIdentificationCuts->Fill(cutIndex, cluster->E());//9
1782  return kFALSE;
1783  }
1784  }
1785  cutIndex++;//5, next cut
1786 
1787  // M02 cut
1788  if (fUseM02 == 1){
1789  if( cluster->GetM02()< fMinM02 || cluster->GetM02() > fMaxM02 ) {
1790  if(fHistClusterIdentificationCuts)fHistClusterIdentificationCuts->Fill(cutIndex, cluster->E());//6
1791  return kFALSE;
1792  }
1793  } else if (fUseM02 ==2 ) {
1794  if( cluster->GetM02()< CalculateMinM02(fMinM02CutNr, cluster->E()) ||
1795  cluster->GetM02() > CalculateMaxM02(fMaxM02CutNr, cluster->E()) ) {
1796  if(fHistClusterIdentificationCuts)fHistClusterIdentificationCuts->Fill(cutIndex, cluster->E());//6
1797  return kFALSE;
1798  }
1799  }
1800  cutIndex++;//6, next cut
1801 
1802  // M20 cut
1803  if (fUseM20){
1804  if( cluster->GetM20()< fMinM20 || cluster->GetM20() > fMaxM20 ) {
1805  if(fHistClusterIdentificationCuts)fHistClusterIdentificationCuts->Fill(cutIndex, cluster->E());//7
1806  return kFALSE;
1807  }
1808  }
1809  cutIndex++;//7, next cut
1810 
1811  // dispersion cut
1812  if (fUseDispersion){
1813  if( cluster->GetDispersion()> fMaxDispersion) {
1814  if(fHistClusterIdentificationCuts)fHistClusterIdentificationCuts->Fill(cutIndex, cluster->E());//8
1815  return kFALSE;
1816  }
1817  }
1818  cutIndex++;//8, next cut
1819 
1820 
1821  // Classify events
1822  AliESDEvent *esdev = dynamic_cast<AliESDEvent*>(event);
1823  AliAODEvent *aodev = 0;
1824  Bool_t isESD = kTRUE;
1825  if (!esdev) {
1826  isESD = kFALSE;
1827  aodev = dynamic_cast<AliAODEvent*>(event);
1828  if (!aodev) {
1829  AliError("Task needs AOD or ESD event, returning");
1830  return kFALSE;
1831  }
1832  }
1833 
1834  fIsCurrentClusterAcceptedBeforeTM = kTRUE;
1835  // classification of clusters for TM efficiency
1836  // 0: Neutral cluster
1837  // 1: Neutral cluster sub charged
1838  // 2: Gamma cluster
1839  // 3: Gamma cluster sub charged
1840  // 4: Gamma conv cluster
1841  // 5: Charged cluster
1842  // 6: Electron
1843  // 7: prim charged cluster
1844  Int_t classification = -1;
1845  Long_t leadMCLabel = -1;
1846  if (isESD)
1847  leadMCLabel = ((AliESDCaloCluster*)cluster)->GetLabel();
1848  else
1849  leadMCLabel = ((AliAODCaloCluster*)cluster)->GetLabel();
1850 
1851  // TM efficiency histograms before TM
1852  if (fIsMC && isMC && (fExtendedMatchAndQA == 1 || fExtendedMatchAndQA == 3 || fExtendedMatchAndQA == 5) && fUseDistTrackToCluster){
1853  classification = ClassifyClusterForTMEffi(cluster, event, mcEvent, isESD);
1854  fHistClusterTMEffiInput->Fill(cluster->E(), 0., weight); //All cl
1855  if (classification == 5 )
1856  fHistClusterTMEffiInput->Fill(cluster->E(), 1., weight); //Ch cl
1857  if (classification == 7 )
1858  fHistClusterTMEffiInput->Fill(cluster->E(), 7., weight); //Ch cl
1859  if (classification == 4)
1860  fHistClusterTMEffiInput->Fill(cluster->E(), 6., weight); //conv electron cl
1861  if (classification == 6){
1862  fHistClusterTMEffiInput->Fill(cluster->E(), 8., weight); // electron cl
1863  if (cluster->E() <= 20.)
1864  fHistClusterElecEtaPhiBeforeTM_00_20->Fill(phiCluster, etaCluster, weight);
1865  if (cluster->E() > 20 && cluster->E() <= 50.)
1866  fHistClusterElecEtaPhiBeforeTM_20_50->Fill(phiCluster, etaCluster, weight);
1867  if (cluster->E() > 50.)
1868  fHistClusterElecEtaPhiBeforeTM_50_00->Fill(phiCluster, etaCluster, weight);
1869  }
1870  if (classification == 0 || classification == 1)
1871  fHistClusterTMEffiInput->Fill(cluster->E(), 2., weight); // Ne cl match
1872  if (classification == 1)
1873  fHistClusterTMEffiInput->Fill(cluster->E(), 3., weight); // Ne cl sub ch match
1874  if (classification == 2 || classification == 3)
1875  fHistClusterTMEffiInput->Fill(cluster->E(), 4., weight); // Ga cl match
1876  if ( classification == 3)
1877  fHistClusterTMEffiInput->Fill(cluster->E(), 5., weight); // Ga cl sub ch match
1878 
1879  Int_t nlabelsMatchedTracks = 0;
1880  if (!fUsePtDepTrackToCluster)
1881  nlabelsMatchedTracks = fCaloTrackMatcher->GetNMatchedTrackIDsForCluster(event, cluster->GetID(), fMaxDistTrackToClusterEta, -fMaxDistTrackToClusterEta,
1882  fMaxDistTrackToClusterPhi, fMinDistTrackToClusterPhi);
1883  else
1884  nlabelsMatchedTracks = fCaloTrackMatcher->GetNMatchedTrackIDsForCluster(event, cluster->GetID(), fFuncPtDepEta, fFuncPtDepPhi);
1885  if (classification < 4 && classification > -1)
1886  fHistClusterENMatchesNeutral->Fill(cluster->E(), nlabelsMatchedTracks);
1887  else
1888  fHistClusterENMatchesCharged->Fill(cluster->E(), nlabelsMatchedTracks);
1889 
1890  // plot electrons that survived the track matching
1891  if (!CheckClusterForTrackMatch(cluster) && classification == 6){ // electrons that survived the matching
1892  if (cluster->E() <= 20.)
1893  fHistClusterElecEtaPhiAfterTM_00_20->Fill(phiCluster, etaCluster, weight);
1894  if (cluster->E() > 20 && cluster->E() <= 50.)
1895  fHistClusterElecEtaPhiAfterTM_20_50->Fill(phiCluster, etaCluster, weight);
1896  if (cluster->E() > 50.)
1897  fHistClusterElecEtaPhiAfterTM_50_00->Fill(phiCluster, etaCluster, weight);
1898  }
1899  }
1900 
1901 
1902  if (fVectorMatchedClusterIDs.size()>0 && fUseDistTrackToCluster){
1903  if( CheckClusterForTrackMatch(cluster) ){
1904  if(fHistClusterIdentificationCuts)fHistClusterIdentificationCuts->Fill(cutIndex, cluster->E());//2
1905  // TM efficiency histos after TM
1906  if (fIsMC && isMC && (fExtendedMatchAndQA == 1 || fExtendedMatchAndQA == 3 || fExtendedMatchAndQA == 5)){
1907 
1908  fHistClusterTMEffiInput->Fill(cluster->E(), 9., weight);
1909  if (classification == 5 )
1910  fHistClusterTMEffiInput->Fill(cluster->E(), 10., weight); //Ch cl match
1911  if (classification == 4)
1912  fHistClusterTMEffiInput->Fill(cluster->E(), 16., weight); //conv cl match
1913  if (classification == 0 || classification == 1)
1914  fHistClusterTMEffiInput->Fill(cluster->E(), 12., weight); // Ne cl match
1915  if ( classification == 1)
1916  fHistClusterTMEffiInput->Fill(cluster->E(), 13., weight); // Ne cl sub ch match
1917  if (classification == 2 || classification == 3)
1918  fHistClusterTMEffiInput->Fill(cluster->E(), 14., weight); // Ga cl match
1919  if ( classification == 3)
1920  fHistClusterTMEffiInput->Fill(cluster->E(), 15., weight); // Ga cl sub ch match
1921  if ( classification == 7)
1922  fHistClusterTMEffiInput->Fill(cluster->E(), 18., weight); // Ch prim cl match
1923  if ( classification == 6)
1924  fHistClusterTMEffiInput->Fill(cluster->E(), 20., weight); // El cl match
1925 
1926  vector<Int_t> labelsMatchedTracks;
1927  if (!fUsePtDepTrackToCluster)
1928  labelsMatchedTracks = fCaloTrackMatcher->GetMatchedTrackIDsForCluster(event, cluster->GetID(), fMaxDistTrackToClusterEta, -fMaxDistTrackToClusterEta,
1929  fMaxDistTrackToClusterPhi, fMinDistTrackToClusterPhi);
1930  else
1931  labelsMatchedTracks = fCaloTrackMatcher->GetMatchedTrackIDsForCluster(event, cluster->GetID(), fFuncPtDepEta, fFuncPtDepPhi);
1932 
1933  //Int_t idHighestPt = -1;
1934  Double_t ptMax = -1;
1935  Double_t eMax = -1;
1936  Double_t eSum = 0;
1937  Bool_t foundLead = kFALSE;
1938  Double_t eLead = -1;
1939  //Int_t idLead = -1;
1940  for (Int_t i = 0; i < (Int_t)labelsMatchedTracks.size(); i++){
1941  AliVTrack* currTrack = dynamic_cast<AliVTrack*>(event->GetTrack(labelsMatchedTracks.at(i)));
1942  eSum += currTrack->E();
1943  if (ptMax < currTrack->Pt()){
1944  ptMax = currTrack->Pt();
1945  eMax = currTrack->E();
1946  //idHighestPt = labelsMatchedTracks.at(i);
1947  }
1948  if (classification == 4 || classification == 5 || classification == 6 || classification == 7){
1949  Long_t mcLabelTrack = -1;
1950  if (isESD)
1951  mcLabelTrack = TMath::Abs(((AliESDtrack*)currTrack)->GetLabel());
1952  else
1953  mcLabelTrack = TMath::Abs(((AliAODTrack*)currTrack)->GetLabel());
1954  if (mcLabelTrack!= -1 && mcLabelTrack == leadMCLabel){
1955  foundLead = kTRUE;
1956  eLead = currTrack->E();
1957  //idLead = labelsMatchedTracks.at(i);
1958  }
1959  }
1960  }
1961  if (classification == 5 || classification == 7 || classification == 6){
1962  fHistClusterEvsTrackECharged->Fill(cluster->E(), eMax, weight);
1963  if (classification == 5 || classification == 7){
1964  fHistClusterEvsTrackEPrimaryButNoElec->Fill(cluster->E(), eMax, weight);
1965  fHistClusterEvsTrackSumEPrimaryButNoElec->Fill(cluster->E(), eSum, weight);
1966  }
1967 
1968  if (foundLead ){
1969  if (classification == 5)
1970  fHistClusterTMEffiInput->Fill(cluster->E(), 11., weight); //Ch cl match w lead
1971  if (classification == 7)
1972  fHistClusterTMEffiInput->Fill(cluster->E(), 19., weight); //Ch prim cl match w lead
1973  if (classification == 6)
1974  fHistClusterTMEffiInput->Fill(cluster->E(), 21., weight); //El cl match w lead
1975  fHistClusterEvsTrackEChargedLead->Fill(cluster->E(), eLead, weight);
1976  }
1977  }
1978  if (classification == 4){
1979  fHistClusterEvsTrackEConv->Fill(cluster->E(), eMax, weight);
1980  if (foundLead)
1981  fHistClusterTMEffiInput->Fill(cluster->E(), 17., weight); //conv cl match w lead
1982  }
1983  if (classification == 0 )
1984  fHistClusterEvsTrackENeutral->Fill(cluster->E(), eMax, weight);
1985  if (classification == 1)
1986  fHistClusterEvsTrackENeutralSubCharged->Fill(cluster->E(), eMax, weight);
1987  if (classification == 2)
1988  fHistClusterEvsTrackEGamma->Fill(cluster->E(), eMax, weight);
1989  if (classification == 3)
1990  fHistClusterEvsTrackEGammaSubCharged->Fill(cluster->E(), eMax, weight);
1991 
1992  labelsMatchedTracks.clear();
1993  }
1994 
1995  return kFALSE;
1996  }
1997  }
1998 
1999  cutIndex++;//9, next cut
2000 
2001  if(GetClusterType() == 1 || GetClusterType() == 3) {
2002  if (cluster->E() > 0.5) fIsAcceptedForBasic = kTRUE;
2003  } else if (GetClusterType() == 2 ){
2004  if (cluster->E() > 0.3) fIsAcceptedForBasic = kTRUE;
2005  }
2006 
2007  // minimum cluster energy cut
2008  if (fUseMinEnergy){
2009  if(cluster->E() < fMinEnergy){
2010  if(fHistClusterIdentificationCuts)fHistClusterIdentificationCuts->Fill(cutIndex, cluster->E());//4
2011  return kFALSE;
2012  }
2013  }
2014  cutIndex++;//10, next cut
2015 
2016 
2017  // DONE with selecting photons
2018  if(fHistClusterIdentificationCuts)fHistClusterIdentificationCuts->Fill(cutIndex, cluster->E());//10
2019 
2020  // Histos after Cuts
2021 
2022  if(fHistClusterEtavsPhiAfterQA) fHistClusterEtavsPhiAfterQA->Fill(phiCluster, etaCluster, weight);
2023  if(fHistClusterTimevsEAfterQA) fHistClusterTimevsEAfterQA->Fill(cluster->GetTOF(), cluster->E(), weight);
2024  if(fHistEnergyOfClusterAfterQA) fHistEnergyOfClusterAfterQA->Fill(cluster->E(), weight);
2025  if(fHistNCellsAfterQA) fHistNCellsAfterQA->Fill(cluster->GetNCells(), weight);
2026  if(fHistM02AfterQA) fHistM02AfterQA->Fill(cluster->GetM02(), weight);
2027  if(fHistM20AfterQA) fHistM20AfterQA->Fill(cluster->GetM20(), weight);
2028  if(fHistDispersionAfterQA) fHistDispersionAfterQA->Fill(cluster->GetDispersion(), weight);
2029  if(fHistNLMAfterQA) fHistNLMAfterQA->Fill(nLM, weight);
2030  if(fHistNLMVsNCellsAfterQA) fHistNLMVsNCellsAfterQA->Fill(nLM,cluster->GetNCells(), weight);
2031  if(fHistNLMVsEAfterQA) fHistNLMVsEAfterQA->Fill(nLM, cluster->E(), weight);
2032  if(fHistClusterEM02AfterQA) fHistClusterEM02AfterQA->Fill(cluster->E(), cluster->GetM02(), weight);
2033 
2034  if(fExtendedMatchAndQA > 1){
2035  if(fHistClusterIncludedCellsAfterQA){
2036  Int_t nCellCluster = cluster->GetNCells();
2037  for(Int_t iCell=0;iCell<nCellCluster;iCell++){
2038  Int_t cellID = cluster->GetCellAbsId(iCell);
2039  Double_t cellAmp = cells->GetCellAmplitude(cellID);
2040  Double_t cellTime = cells->GetCellTime(cellID);
2041  fHistClusterIncludedCellsAfterQA->Fill(cellID);
2042  if(cluster->E()>0.) fHistClusterEnergyFracCellsAfterQA->Fill(cellID,cellAmp/cluster->E());
2043  if(isMC==0){
2044  fHistClusterIncludedCellsTimingAfterQA->Fill(cluster->E(),cellTime*1E9);
2045  fHistClusterIncludedCellsTimingEnergyAfterQA->Fill(cellAmp,cellTime*1E9);
2046  }else{
2047  fHistClusterIncludedCellsTimingAfterQA->Fill(cluster->E(),cellTime*1E8);
2048  fHistClusterIncludedCellsTimingEnergyAfterQA->Fill(cellAmp,cellTime*1E8);
2049  }
2050  }
2051  }
2052 
2053  if(fHistClusterEnergyvsNCells) fHistClusterEnergyvsNCells->Fill(cluster->E(),cluster->GetNCells());
2054  if(cluster->IsEMCAL()){
2055  Int_t iSuperModule = -1;
2056  fGeomEMCAL = AliEMCALGeometry::GetInstance();
2057  if(!fGeomEMCAL){ AliFatal("EMCal geometry not initialized!");}
2058  if(fHistClusterEnergyvsMod && fGeomEMCAL->SuperModuleNumberFromEtaPhi(clusterVector.Eta(),clusterVector.Phi(),iSuperModule)){
2059  fHistClusterEnergyvsMod->Fill(cluster->E(),iSuperModule);
2060  }
2061  }else if(cluster->IsPHOS()){
2062  Int_t relId[4] = {0,0,0,0};
2063  fGeomPHOS = AliPHOSGeometry::GetInstance();
2064  if(!fGeomPHOS){ AliFatal("PHOS geometry not initialized!");}
2065  if(fHistClusterEnergyvsMod && fGeomPHOS->GlobalPos2RelId(clusterVector,relId)){
2066  fHistClusterEnergyvsMod->Fill(cluster->E(),relId[0]);
2067  }
2068  }
2069  }
2070 
2071  return kTRUE;
2072 }
2073 
2074 
2075 
2076 //________________________________________________________________________
2077 void AliCaloPhotonCuts::FillHistogramsExtendedQA(AliVEvent *event, Int_t isMC)
2078 {
2079  if(fExtendedMatchAndQA < 2) return;
2080 
2081  AliVCaloCells* cells = 0x0;
2082 
2083  Int_t nModules = 0;
2084  Int_t* nCellsBigger100MeV;
2085  Int_t* nCellsBigger1500MeV;
2086  Double_t* EnergyOfMod;
2087  if( (fClusterType == 1 || fClusterType == 3) && !fEMCALInitialized ) InitializeEMCAL(event);
2088  if( fClusterType == 2 && ( !fPHOSInitialized || (fPHOSCurrentRun != event->GetRunNumber()) ) ) InitializePHOS(event);
2089 
2090  Int_t nModulesStart = 0;
2091  if( fClusterType == 1 || fClusterType == 3){ //EMCAL & DCAL
2092  cells = event->GetEMCALCells();
2093  fGeomEMCAL = AliEMCALGeometry::GetInstance();
2094  if(!fGeomEMCAL) AliFatal("EMCal geometry not initialized!");
2095  if(!fEMCALBadChannelsMap) AliFatal("EMCal bad channels map not initialized!");
2096  nModules = fGeomEMCAL->GetNumberOfSuperModules();
2097  if( fClusterType == 3) {nModules = 8; nModulesStart = 12;}
2098  } else if( fClusterType == 2 ){ //PHOS
2099  cells = event->GetPHOSCells();
2100  fGeomPHOS = AliPHOSGeometry::GetInstance();
2101  if(!fGeomPHOS) AliFatal("PHOS geometry not initialized!");
2102  nModules = fGeomPHOS->GetNModules();
2103  } else{
2104  AliError(Form("fExtendedMatchAndQA(%i):FillHistogramsExtendedMatchAndQA() not (yet) defined for cluster type (%i)",fExtendedMatchAndQA,fClusterType));
2105  }
2106 
2107  nCellsBigger100MeV = new Int_t[nModules];
2108  nCellsBigger1500MeV = new Int_t[nModules];
2109  EnergyOfMod = new Double_t[nModules];
2110 
2111  for(Int_t iModule=0;iModule<nModules;iModule++){nCellsBigger100MeV[iModule]=0;nCellsBigger1500MeV[iModule]=0;EnergyOfMod[iModule]=0;}
2112 
2113  for(Int_t iCell=0;iCell<cells->GetNumberOfCells();iCell++){
2114  Short_t cellNumber=0;
2115  Double_t cellAmplitude=0;
2116  Double_t cellTime=0;
2117  Double_t cellEFrac=0;
2118  Int_t cellMCLabel=0;
2119  Int_t nMod = -1;
2120 
2121  cells->GetCell(iCell,cellNumber,cellAmplitude,cellTime,cellMCLabel,cellEFrac);
2122  if( fClusterType == 3 && cellNumber < 12288){continue;}
2123  if( fClusterType == 2 && cellNumber < 0){continue;} //Scip CPV cells in PHOS case
2124  Int_t imod = -1;Int_t iTower = -1, iIphi = -1, iIeta = -1;
2125  Int_t icol = -1;Int_t irow = -1;
2126  Int_t relid[4];// for PHOS
2127 
2128  Bool_t doBadCell = kTRUE;
2129  if( fClusterType == 1 || fClusterType == 3){
2130  nMod = fGeomEMCAL->GetSuperModuleNumber(cellNumber);
2131  fGeomEMCAL->GetCellIndex(cellNumber,imod,iTower,iIphi,iIeta);
2132  if (fEMCALBadChannelsMap->GetEntries() <= imod) doBadCell=kFALSE;
2133  fGeomEMCAL->GetCellPhiEtaIndexInSModule(imod,iTower,iIphi,iIeta,irow,icol);
2134  }else if( fClusterType == 2 ){
2135  fGeomPHOS->AbsToRelNumbering(cellNumber,relid);
2136  if(relid[1]!=0) AliFatal("PHOS CPV in PHOS cell array?");
2137  nMod = relid[0];//(Int_t) (1 + (cellNumber-1)/3584);
2138  if(nMod>=nModules || nMod<0 || !fPHOSBadChannelsMap[nMod]) doBadCell=kFALSE;
2139  }
2140 
2141  Int_t iBadCell = 0;
2142  if( (fClusterType == 1 || fClusterType == 3) && doBadCell){
2143  iBadCell = (Int_t) ((TH2I*)fEMCALBadChannelsMap->At(imod))->GetBinContent(icol,irow);
2144  }else if( fClusterType == 2 && doBadCell){
2145  iBadCell = (Int_t) ((TH2I*)fPHOSBadChannelsMap[nMod])->GetBinContent(relid[2],relid[3]);
2146  }
2147 
2148  if(iBadCell > 0) continue;
2149  // nModulesStart == 0 for EMCAL and PHOS
2150  if(cellAmplitude > 0.1) nCellsBigger100MeV[nMod-nModulesStart]++;
2151  if(cellAmplitude > 1.5) nCellsBigger1500MeV[nMod-nModulesStart]++;
2152  if(cellAmplitude > 0.05) EnergyOfMod[nMod-nModulesStart]+=cellAmplitude;
2153 
2154  if(fExtendedMatchAndQA > 3){
2155  if(fHistCellEnergyvsCellID && (cellAmplitude > 0.05)) fHistCellEnergyvsCellID->Fill(cellAmplitude,cellNumber);
2156  if(fHistCellTimevsCellID && (cellAmplitude > 0.2)) fHistCellTimevsCellID->Fill(cellTime,cellNumber);
2157  }
2158  }
2159 
2160  for(Int_t iModule=0;iModule<nModules;iModule++){
2161  if(fHistNCellsBigger100MeVvsMod) fHistNCellsBigger100MeVvsMod->Fill(nCellsBigger100MeV[iModule],iModule+nModulesStart);
2162  if(fHistNCellsBigger1500MeVvsMod) fHistNCellsBigger1500MeVvsMod->Fill(nCellsBigger1500MeV[iModule],iModule+nModulesStart);
2163  if(fHistEnergyOfModvsMod) fHistEnergyOfModvsMod->Fill(EnergyOfMod[iModule],iModule+nModulesStart);
2164  }
2165 
2166  delete[] nCellsBigger100MeV;nCellsBigger100MeV=0x0;
2167  delete[] nCellsBigger1500MeV;nCellsBigger1500MeV=0x0;
2168  delete[] EnergyOfMod;EnergyOfMod=0x0;
2169 
2170  //fill distClusterTo_withinTiming/outsideTiming
2171  Int_t nclus = 0;
2172  TClonesArray * arrClustersExtQA = NULL;
2173  if(!fCorrTaskSetting.CompareTo("")){
2174  nclus = event->GetNumberOfCaloClusters();
2175  } else {
2176  arrClustersExtQA = dynamic_cast<TClonesArray*>(event->FindListObject(Form("%sClustersBranch",fCorrTaskSetting.Data())));
2177  if(!arrClustersExtQA)
2178  AliFatal(Form("%sClustersBranch was not found in AliCaloPhotonCuts::FillHistogramsExtendedQA! Check the correction framework settings!",fCorrTaskSetting.Data()));
2179  nclus = arrClustersExtQA->GetEntries();
2180  }
2181  AliVCluster* cluster = 0x0;
2182  AliVCluster* clusterMatched = 0x0;
2183  for(Int_t iClus=0; iClus<nclus ; iClus++){
2184  if(event->IsA()==AliESDEvent::Class()){
2185  if(arrClustersExtQA)
2186  cluster = new AliESDCaloCluster(*(AliESDCaloCluster*)arrClustersExtQA->At(iClus));
2187  else
2188  cluster = new AliESDCaloCluster(*(AliESDCaloCluster*)event->GetCaloCluster(iClus));
2189  } else if(event->IsA()==AliAODEvent::Class()){
2190  if(arrClustersExtQA)
2191  cluster = new AliAODCaloCluster(*(AliAODCaloCluster*)arrClustersExtQA->At(iClus));
2192  else
2193  cluster = new AliAODCaloCluster(*(AliAODCaloCluster*)event->GetCaloCluster(iClus));
2194  }
2195 
2196  if( (fClusterType == 1 || fClusterType == 3) && !cluster->IsEMCAL()){delete cluster; continue;}
2197  if( fClusterType == 2 && cluster->GetType() !=AliVCluster::kPHOSNeutral){delete cluster; continue;}
2198 
2199  Float_t clusPos[3]={0,0,0};
2200  cluster->GetPosition(clusPos);
2201  TVector3 clusterVector(clusPos[0],clusPos[1],clusPos[2]);
2202  Double_t etaCluster = clusterVector.Eta();
2203  Double_t phiCluster = clusterVector.Phi();
2204  if (phiCluster < 0) phiCluster += 2*TMath::Pi();
2205  Int_t nLM = GetNumberOfLocalMaxima(cluster, event);
2206 
2207  //acceptance cuts
2208  if (fUseEtaCut && (etaCluster < fMinEtaCut || etaCluster > fMaxEtaCut)){delete cluster; continue;}
2209  if (fUseEtaCut && fClusterType == 3 && etaCluster < fMaxEtaInnerEdge && etaCluster > fMinEtaInnerEdge ) {delete cluster; continue;}
2210  if (fUsePhiCut && (phiCluster < fMinPhiCut || phiCluster > fMaxPhiCut)){delete cluster; continue;}
2211  if (fUseDistanceToBadChannel>0 && CheckDistanceToBadChannel(cluster,event)){delete cluster; continue;}
2212  //cluster quality cuts
2213  if (fVectorMatchedClusterIDs.size()>0 && CheckClusterForTrackMatch(cluster)){delete cluster; continue;}
2214  if (fUseMinEnergy && (cluster->E() < fMinEnergy)){delete cluster; continue;}
2215  if (fUseNCells && (cluster->GetNCells() < fMinNCells)){delete cluster; continue;}
2216  if (fUseNLM && (nLM < fMinNLM || nLM > fMaxNLM)){delete cluster; continue;}
2217  if (fUseM02 == 1 && (cluster->GetM02() < fMinM02 || cluster->GetM02() > fMaxM02)){delete cluster; continue;}
2218  if (fUseM02 == 2 && (cluster->GetM02() < CalculateMinM02(fMinM02CutNr, cluster->E()) || cluster->GetM02() > CalculateMaxM02(fMaxM02CutNr, cluster->E()))){delete cluster; continue;}
2219  if (fUseM20 && (cluster->GetM20() < fMinM20 || cluster->GetM20() > fMaxM20)){delete cluster; continue;}
2220  if (fUseDispersion && (cluster->GetDispersion() > fMaxDispersion)){delete cluster; continue;}
2221  //cluster within timing cut
2222  if (!(isMC>0) && (cluster->GetTOF() < fMinTimeDiff || cluster->GetTOF() > fMaxTimeDiff)){delete cluster; continue;}
2223 
2224  Int_t largestCellicol = -1, largestCellirow = -1;
2225  Int_t largestCellID = FindLargestCellInCluster(cluster,event);
2226  if(largestCellID==-1) AliFatal("FillHistogramsExtendedQA: FindLargestCellInCluster found cluster with NCells<1?");
2227  Int_t largestCelliMod = GetModuleNumberAndCellPosition(largestCellID, largestCellicol, largestCellirow);
2228  if(largestCelliMod < 0) AliFatal("FillHistogramsExtendedQA: GetModuleNumberAndCellPosition found SM with ID<0?");
2229 
2230  for(Int_t iClus2=iClus+1; iClus2<nclus; iClus2++){
2231  if(event->IsA()==AliESDEvent::Class()){
2232  if(arrClustersExtQA)
2233  clusterMatched = new AliESDCaloCluster(*(AliESDCaloCluster*)arrClustersExtQA->At(iClus2));
2234  else
2235  clusterMatched = new AliESDCaloCluster(*(AliESDCaloCluster*)event->GetCaloCluster(iClus2));
2236  } else if(event->IsA()==AliAODEvent::Class()){
2237  if(arrClustersExtQA)
2238  clusterMatched = new AliAODCaloCluster(*(AliAODCaloCluster*)arrClustersExtQA->At(iClus2));
2239  else
2240  clusterMatched = new AliAODCaloCluster(*(AliAODCaloCluster*)event->GetCaloCluster(iClus2));
2241  }
2242 
2243  if( (fClusterType == 1 || fClusterType == 3) && !clusterMatched->IsEMCAL()){delete clusterMatched; continue;}
2244  if( fClusterType == 2 && clusterMatched->GetType() !=AliVCluster::kPHOSNeutral){delete clusterMatched; continue;}
2245 
2246  Float_t clusPos2[3]={0,0,0};
2247  clusterMatched->GetPosition(clusPos2);
2248  TVector3 clusterMatchedVector(clusPos2[0],clusPos2[1],clusPos2[2]);
2249  Double_t etaclusterMatched = clusterMatchedVector.Eta();
2250  Double_t phiclusterMatched = clusterMatchedVector.Phi();
2251  if (phiclusterMatched < 0) phiclusterMatched += 2*TMath::Pi();
2252  Int_t nLMMatched = GetNumberOfLocalMaxima(clusterMatched, event);
2253 
2254  //acceptance cuts
2255  if (fUseEtaCut && (etaclusterMatched < fMinEtaCut || etaclusterMatched > fMaxEtaCut)){delete clusterMatched; continue;}
2256  if (fUseEtaCut && fClusterType == 3 && etaclusterMatched < fMaxEtaInnerEdge && etaclusterMatched > fMinEtaInnerEdge ) {delete clusterMatched; continue;}
2257  if (fUsePhiCut && (phiclusterMatched < fMinPhiCut || phiclusterMatched > fMaxPhiCut)){delete clusterMatched; continue;}
2258  if (fUseDistanceToBadChannel>0 && CheckDistanceToBadChannel(clusterMatched,event)){delete clusterMatched; continue;}
2259  //cluster quality cuts
2260  if (fVectorMatchedClusterIDs.size()>0 && CheckClusterForTrackMatch(clusterMatched)){delete clusterMatched; continue;}
2261  if (fUseMinEnergy && (clusterMatched->E() < fMinEnergy)){delete clusterMatched; continue;}
2262  if (fUseNCells && (clusterMatched->GetNCells() < fMinNCells)){delete clusterMatched; continue;}
2263  if (fUseNLM && (nLMMatched < fMinNLM || nLMMatched > fMaxNLM)){delete clusterMatched; continue;}
2264  if (fUseM02 == 1 && (clusterMatched->GetM02() < fMinM02 || clusterMatched->GetM02() > fMaxM02)){delete clusterMatched; continue;}
2265  if (fUseM02 == 2 && (clusterMatched->GetM02() < CalculateMinM02(fMinM02CutNr, clusterMatched->E()) || cluster->GetM02() > CalculateMaxM02(fMaxM02CutNr, clusterMatched->E()))){delete clusterMatched; continue;}
2266  if (fUseM20 && (clusterMatched->GetM20() < fMinM20 || clusterMatched->GetM20() > fMaxM20)){delete clusterMatched; continue;}
2267  if (fUseDispersion && (clusterMatched->GetDispersion() > fMaxDispersion)){delete clusterMatched; continue;}
2268 
2269  // Get rowdiff and coldiff
2270 
2271  Int_t matched_largestCellicol = -1, matched_largestCellirow = -1;
2272  Int_t matched_largestCellID = FindLargestCellInCluster(clusterMatched,event);
2273  if(matched_largestCellID==-1) AliFatal("FillHistogramsExtendedQA: FindLargestCellInCluster found cluster with NCells<1?");
2274  Int_t matched_largestCelliMod = GetModuleNumberAndCellPosition(matched_largestCellID, matched_largestCellicol, matched_largestCellirow);
2275  if(matched_largestCelliMod < 0) AliFatal("FillHistogramsExtendedQA: GetModuleNumberAndCellPosition found SM with ID<0?");
2276 
2277 // cout << "\n+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++" << endl;
2278 // cout << "Cluster: " << largestCelliMod << ", " << largestCellirow << ", " << largestCellicol << " , time: " << cluster->GetTOF() << endl;
2279 // cout << "Matched: " << matched_largestCelliMod << ", " << matched_largestCellirow << ", " << matched_largestCellicol << " , time: " << clusterMatched->GetTOF() << endl;
2280 
2281  Int_t rowdiff = -100;
2282  Int_t coldiff = -100;
2283  Bool_t calculatedDiff = kFALSE;
2284 
2285  Int_t ClusID = largestCelliMod/2;
2286  Int_t matchClusID = matched_largestCelliMod/2;
2287 
2288  if( matched_largestCelliMod == largestCelliMod){
2289  rowdiff = largestCellirow - matched_largestCellirow;
2290  coldiff = largestCellicol - matched_largestCellicol;
2291  calculatedDiff = kTRUE;
2292  }else if( TMath::Abs(matched_largestCelliMod - largestCelliMod) == 1 && (ClusID == matchClusID) ){
2293  if(matched_largestCelliMod%2){
2294  matched_largestCelliMod -= 1;
2295  matched_largestCellicol += AliEMCALGeoParams::fgkEMCALCols;
2296  }else{
2297  matched_largestCelliMod += 1;
2298  matched_largestCellicol -= AliEMCALGeoParams::fgkEMCALCols;
2299  }
2300 
2301  if( matched_largestCelliMod == largestCelliMod ){
2302  rowdiff = largestCellirow - matched_largestCellirow;
2303  coldiff = largestCellicol - matched_largestCellicol;
2304  calculatedDiff = kTRUE;
2305  }
2306  }
2307 // cout << "\t\t ROWDIFF: " << rowdiff << endl;
2308 // cout << "\t\t COLDIFF: " << coldiff << endl;
2309 // cout << "+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n" << endl;
2310  //cluster outside timing cut
2311  if( calculatedDiff ){
2312  Float_t dist1D = TMath::Sqrt(TMath::Power(etaCluster-etaclusterMatched,2)+TMath::Power(phiCluster-phiclusterMatched,2));
2313  if( !(isMC>0) ){
2314  if( (clusterMatched->GetTOF() > fMinTimeDiff && clusterMatched->GetTOF() < fMaxTimeDiff) ){
2315  fHistClusterDistanceInTimeCut->Fill(rowdiff,coldiff);
2316  fHistClusterDistance1DInTimeCut->Fill(dist1D);
2317  }
2318  else fHistClusterDistanceOutTimeCut->Fill(rowdiff,coldiff);
2319  }else{
2320  fHistClusterDistanceInTimeCut->Fill(rowdiff,coldiff);
2321  fHistClusterDistance1DInTimeCut->Fill(dist1D);
2322  }
2323  }
2324  delete clusterMatched;
2325  }
2326 
2327  delete cluster;
2328  }
2329  return;
2330 }
2331 
2332 //________________________________________________________________________
2333 //************** Find number of local maxima in cluster ******************
2334 //* derived from G. Conesa Balbastre's AliCalorimeterUtils *******************
2335 //************************************************************************
2336 Int_t AliCaloPhotonCuts::GetNumberOfLocalMaxima(AliVCluster* cluster, AliVEvent * event){
2337 
2338  const Int_t nc = cluster->GetNCells();
2339 
2340  Int_t absCellIdList[nc];
2341  Float_t maxEList[nc];
2342 
2343  Int_t nMax = GetNumberOfLocalMaxima(cluster, event, absCellIdList, maxEList);
2344 
2345  return nMax;
2346 }
2347 
2348 //________________________________________________________________________
2349 Int_t AliCaloPhotonCuts::FindSecondLargestCellInCluster(AliVCluster* cluster, AliVEvent* event){
2350 
2351  const Int_t nCells = cluster->GetNCells();
2352  AliVCaloCells* cells = NULL;
2353 
2354  if (fClusterType == 1 || fClusterType == 3)
2355  cells = event->GetEMCALCells();
2356  else if (fClusterType ==2 )
2357  cells = event->GetPHOSCells();
2358 
2359 // cout << "NCells: "<< nCells<< " cluster energy: " << cluster->E() << endl;
2360  Float_t eMax = 0.;
2361  Int_t idMax = -1;
2362  Int_t idMax2 = -1;
2363  Int_t iCellMax = -1;
2364 
2365  if (nCells < 2) return idMax;
2366  for (Int_t iCell = 1;iCell < nCells;iCell++){
2367  if (cells->GetCellAmplitude(cluster->GetCellsAbsId()[iCell])> eMax){
2368  eMax = cells->GetCellAmplitude(cluster->GetCellsAbsId()[iCell]);
2369  idMax = cluster->GetCellsAbsId()[iCell];
2370  iCellMax = iCell;
2371  }
2372  }
2373 
2374  eMax = 0.;
2375  for (Int_t iCell = 1;iCell < nCells;iCell++){
2376  if (iCell == iCellMax) continue;
2377  if (cells->GetCellAmplitude(cluster->GetCellsAbsId()[iCell])> eMax){
2378  eMax = cells->GetCellAmplitude(cluster->GetCellsAbsId()[iCell]);
2379  idMax2 = cluster->GetCellsAbsId()[iCell];
2380  }
2381  }
2382 
2383  return idMax2;
2384 }
2385 
2386 //________________________________________________________________________
2387 Int_t AliCaloPhotonCuts::FindLargestCellInCluster(AliVCluster* cluster, AliVEvent* event){
2388 
2389  const Int_t nCells = cluster->GetNCells();
2390  AliVCaloCells* cells = NULL;
2391 
2392  if (fClusterType == 1 || fClusterType == 3)
2393  cells = event->GetEMCALCells();
2394  else if (fClusterType ==2 )
2395  cells = event->GetPHOSCells();
2396 
2397 // cout << "NCells: "<< nCells<< " cluster energy: " << cluster->E() << endl;
2398  Float_t eMax = 0.;
2399  Int_t idMax = -1;
2400 
2401  if (nCells < 1) return idMax;
2402  for (Int_t iCell = 0;iCell < nCells;iCell++){
2403  Int_t cellAbsID = cluster->GetCellsAbsId()[iCell];
2404  if (cells->GetCellAmplitude(cellAbsID)> eMax){
2405  eMax = cells->GetCellAmplitude(cellAbsID);
2406  idMax = cellAbsID;
2407  }
2408  }
2409  return idMax;
2410 
2411 }
2412 
2413 
2414 //________________________________________________________________________
2415 //************** Find number of local maxima in cluster ******************
2416 //* derived from G. Conesa Balbastre's AliCalorimeterUtils ***************
2417 //************************************************************************
2418 Int_t AliCaloPhotonCuts::GetNumberOfLocalMaxima(AliVCluster* cluster, AliVEvent * event, Int_t *absCellIdList, Float_t* maxEList){
2419 
2420  Int_t absCellId1 = -1;
2421  Int_t absCellId2 = -1;
2422  const Int_t nCells = cluster->GetNCells();
2423  AliVCaloCells* cells = NULL;
2424 
2425  if (fClusterType == 1 || fClusterType == 3)
2426  cells = event->GetEMCALCells();
2427  else if (fClusterType ==2 )
2428  cells = event->GetPHOSCells();
2429 
2430 // cout << "NCells: "<< nCells<< " cluster energy: " << cluster->E() << endl;
2431  Float_t eMax = 0.;
2432  Int_t idMax = -1;
2433 
2434  for (Int_t iCell = 0;iCell < nCells;iCell++){
2435  absCellIdList[iCell] = cluster->GetCellsAbsId()[iCell];
2436 // Int_t imod = -1, icol = -1, irow = -1;
2437 // imod = GetModuleNumberAndCellPosition(absCellIdList[iCell], icol, irow);
2438 // cout << absCellIdList[iCell] <<"\t" << cells->GetCellAmplitude(absCellIdList[iCell]) << "\t"<< imod << "\t" << icol << "\t" << irow << endl;
2439  if (cells->GetCellAmplitude(absCellIdList[iCell])> eMax){
2440  eMax = cells->GetCellAmplitude(absCellIdList[iCell]);
2441  idMax = absCellIdList[iCell];
2442  }
2443  }
2444 
2445  // find the largest separated cells in cluster
2446  for (Int_t iCell = 0;iCell < nCells;iCell++){
2447  // check whether valid cell number is selected
2448  if (absCellIdList[iCell] >= 0){
2449  // store current energy and cell id
2450  absCellId1 = cluster->GetCellsAbsId()[iCell];
2451  Float_t en1 = cells->GetCellAmplitude(absCellId1);
2452 
2453  // loop over other cells in cluster
2454  for (Int_t iCellN = 0;iCellN < nCells;iCellN++){
2455  // jump out if array has changed in the meantime
2456  if (absCellIdList[iCell] == -1) continue;
2457  // get cell id & check whether its valid
2458  absCellId2 = cluster->GetCellsAbsId()[iCellN];
2459 
2460  // don't compare to yourself
2461  if (absCellId2 == -1) continue;
2462  if (absCellId1 == absCellId2) continue;
2463 
2464  // get cell energy of second cell
2465  Float_t en2 = cells->GetCellAmplitude(absCellId2);
2466 
2467  // check if cells are Neighbours
2468  if (AreNeighbours(absCellId1, absCellId2)){
2469  // determine which cell has larger energy, mask the other
2470 // cout << "found neighbour: " << absCellId1 << "\t" << absCellId2 << endl;
2471 // cout << "energies: " << en1 << "\t" << en2 << endl;
2472  if (en1 > en2 ){
2473  absCellIdList[iCellN] = -1;
2474  if (en1 < en2 + fLocMaxCutEDiff)
2475  absCellIdList[iCell] = -1;
2476  } else {
2477  absCellIdList[iCell] = -1;
2478  if (en1 > en2 - fLocMaxCutEDiff)
2479  absCellIdList[iCellN] = -1;
2480  }
2481  }
2482  }
2483  }
2484  }
2485 
2486  // shrink list of cells to only maxima
2487  Int_t nMaximaNew = 0;
2488  for (Int_t iCell = 0;iCell < nCells;iCell++){
2489 // cout << iCell << "\t" << absCellIdList[iCell] << endl;
2490  if (absCellIdList[iCell] > -1){
2491  Float_t en = cells->GetCellAmplitude(absCellIdList[iCell]);
2492  // check whether cell energy is larger than required seed
2493  if (en < fSeedEnergy) continue;
2494  absCellIdList[nMaximaNew] = absCellIdList[iCell];
2495  maxEList[nMaximaNew] = en;
2496  nMaximaNew++;
2497  }
2498  }
2499 
2500  // check whether a local maximum was found
2501  // if no maximum was found use highest cell as maximum
2502  if (nMaximaNew == 0){
2503  nMaximaNew = 1;
2504  maxEList[0] = eMax;
2505  absCellIdList[0] = idMax;
2506  }
2507 
2508  return nMaximaNew;
2509 }
2510 
2511 //________________________________________________________________________
2512 //************** Function to determine neighbours of cells ***************
2513 //* derived from G. Conesa Balbastre's AliCalorimeterUtils ***************
2514 //************************************************************************
2515 Bool_t AliCaloPhotonCuts::AreNeighbours(Int_t absCellId1, Int_t absCellId2){
2516  Bool_t areNeighbours = kFALSE ;
2517 
2518  Int_t irow1 = -1, icol1 = -1;
2519  Int_t irow2 = -1, icol2 = -1;
2520 
2521  Int_t rowdiff = 0, coldiff = 0;
2522 
2523  Int_t nSupMod1 = GetModuleNumberAndCellPosition(absCellId1, icol1, irow1);
2524  Int_t nSupMod2 = GetModuleNumberAndCellPosition(absCellId2, icol2, irow2);
2525 
2526  // check if super modules are correct
2527  if (nSupMod1== -1 || nSupMod2 == -1) return areNeighbours;
2528 
2529  if(fClusterType==1 && nSupMod1!=nSupMod2) {
2530  // In case of a shared cluster, index of SM in C side, columns start at 48 and ends at 48*2-1
2531  // C Side impair SM, nSupMod%2=1;A side pair SM nSupMod%2=0
2532  if(nSupMod1%2) icol1+=AliEMCALGeoParams::fgkEMCALCols;
2533  else icol2+=AliEMCALGeoParams::fgkEMCALCols;
2534  }
2535 
2536  rowdiff = TMath::Abs( irow1 - irow2 ) ;
2537  coldiff = TMath::Abs( icol1 - icol2 ) ;
2538 
2539 // if (( coldiff <= 1 ) && ( rowdiff <= 1 ) && (coldiff + rowdiff <= 2))
2540  if ((coldiff + rowdiff == 1 ))
2541  areNeighbours = kTRUE ;
2542 
2543  return areNeighbours;
2544 }
2545 
2546 
2547 //________________________________________________________________________
2548 //************** Function to obtain module number, row and column ********
2549 //* derived from G. Conesa Balbastre's AliCalorimeterUtils ***************
2550 //************************************************************************
2551 Int_t AliCaloPhotonCuts::GetModuleNumberAndCellPosition(Int_t absCellId, Int_t & icol, Int_t & irow){
2552  if( fClusterType == 1 || fClusterType == 3){ //EMCAL & DCAL
2553  fGeomEMCAL = AliEMCALGeometry::GetInstance();
2554  if(!fGeomEMCAL) AliFatal("EMCal geometry not initialized!");
2555  } else if( fClusterType == 2 ){ //PHOS
2556  fGeomPHOS = AliPHOSGeometry::GetInstance();
2557  if(!fGeomPHOS) AliFatal("PHOS geometry not initialized!");
2558  }
2559 
2560  Int_t imod = -1;Int_t iTower = -1, iIphi = -1, iIeta = -1;
2561  if( fClusterType == 1 || fClusterType == 3){
2562  fGeomEMCAL->GetCellIndex(absCellId,imod,iTower,iIphi,iIeta);
2563  fGeomEMCAL->GetCellPhiEtaIndexInSModule(imod,iTower,iIphi,iIeta,irow,icol);
2564  } else if ( fClusterType == 2 ){
2565  Int_t relId[4];
2566  fGeomPHOS->AbsToRelNumbering(absCellId,relId);
2567  irow = relId[2];
2568  icol = relId[3];
2569  imod = relId[0]-1;
2570  }
2571  return imod;
2572 }
2573 
2574 //___________________________________________________________________________
2575 // Split energy of cluster between the 2 local maxima, sum energy on 3x3, and if the 2
2576 // maxima are too close and have common cells, split the energy between the 2.
2577 //* derived from G. Conesa Balbastre's AliCalorimeterUtils *******************
2578 //___________________________________________________________________________
2579 void AliCaloPhotonCuts::SplitEnergy(Int_t absCellId1, Int_t absCellId2,
2580  AliVCluster* cluster,
2581  AliVEvent* event,
2582  Int_t isMC,
2583  AliAODCaloCluster* cluster1,
2584  AliAODCaloCluster* cluster2){
2585 
2586  const Int_t ncells = cluster->GetNCells();
2587  Int_t absCellIdList[ncells];
2588 
2589  AliVCaloCells* cells = NULL;
2590  if (fClusterType == 1 || fClusterType == 3)
2591  cells = event->GetEMCALCells();
2592  else if (fClusterType ==2 )
2593  cells = event->GetPHOSCells();
2594 
2595  Float_t e1 = 0;
2596  Float_t e2 = 0;
2597  Float_t eCluster = 0;
2598 
2599  for(Int_t iCell = 0;iCell < ncells;iCell++ ) {
2600  absCellIdList[iCell] = cluster->GetCellsAbsId()[iCell];
2601  Float_t ec = cells->GetCellAmplitude(absCellIdList[iCell]);
2602  eCluster+=ec;
2603  }
2604 
2605  UShort_t absCellIdList1 [12];
2606  Double_t fracList1 [12];
2607  UShort_t absCellIdList2 [12];
2608  Double_t fracList2 [12];
2609 
2610  // Init counters and variables
2611  Int_t ncells1 = 1 ;
2612  absCellIdList1[0] = absCellId1 ;
2613  fracList1 [0] = 1. ;
2614 
2615  Float_t ecell1 = cells->GetCellAmplitude(absCellId1);
2616  e1 = ecell1;
2617 
2618  Int_t ncells2 = 1 ;
2619  absCellIdList2[0] = absCellId2 ;
2620  fracList2 [0] = 1. ;
2621 
2622  Float_t ecell2 = cells->GetCellAmplitude(absCellId2);
2623  e2 = ecell2;
2624 
2625 // cout << "Cluster: " << eCluster << "\t cell1: " << absCellId1 << "\t" << e1 << "\t cell2: " << absCellId2 << "\t" << e2 << endl;
2626  // Very rough way to share the cluster energy
2627  Float_t eRemain = (eCluster-ecell1-ecell2)/2;
2628  Float_t shareFraction1 = (ecell1+eRemain)/eCluster;
2629  Float_t shareFraction2 = (ecell2+eRemain)/eCluster;
2630 
2631 // cout << eRemain << "\t" << shareFraction1<< "\t" << shareFraction2 << endl;
2632 
2633  for(Int_t iCell = 0;iCell < ncells;iCell++){
2634 
2635  Int_t absId = absCellIdList[iCell];
2636  if ( absId==absCellId1 || absId==absCellId2 || absId < 0 ) continue;
2637 
2638  Float_t ecell = cells->GetCellAmplitude(absId);
2639  if(AreNeighbours(absCellId1,absId )){
2640  absCellIdList1[ncells1] = absId;
2641  if(AreNeighbours(absCellId2,absId )){
2642  fracList1[ncells1] = shareFraction1;
2643  e1 += ecell*shareFraction1;
2644  } else {
2645  fracList1[ncells1] = 1.;
2646  e1 += ecell;
2647  }
2648  ncells1++;
2649  } // neigbour to cell1
2650 
2651  if(AreNeighbours(absCellId2,absId )) {
2652  absCellIdList2[ncells2]= absId;
2653 
2654  if(AreNeighbours(absCellId1,absId )){
2655  fracList2[ncells2] = shareFraction2;
2656  e2 += ecell*shareFraction2;
2657  } else {
2658  fracList2[ncells2] = 1.;
2659  e2 += ecell;
2660  }
2661  ncells2++;
2662  } // neigbour to cell2
2663  }
2664 // cout << "Cluster: " << eCluster << "\t cell1: " << absCellId1 << "\t" << e1 << "\t cell2: " << absCellId2 << "\t" << e2 << endl;
2665 
2666  cluster1->SetE(e1);
2667  cluster2->SetE(e2);
2668 
2669  cluster1->SetNCells(ncells1);
2670  cluster2->SetNCells(ncells2);
2671 
2672  cluster1->SetCellsAbsId(absCellIdList1);
2673  cluster2->SetCellsAbsId(absCellIdList2);
2674 
2675  cluster1->SetCellsAmplitudeFraction(fracList1);
2676  cluster2->SetCellsAmplitudeFraction(fracList2);
2677 
2678  // Correct linearity
2679  ApplyNonLinearity(cluster1, isMC) ;
2680  ApplyNonLinearity(cluster2, isMC) ;
2681 
2682  // Initialize EMCAL rec utils if not initialized
2683  if(!fEMCALInitialized && (fClusterType == 1 || fClusterType == 3) ) InitializeEMCAL(event);
2684 
2685  if(fEMCALInitialized && (fClusterType == 1 || fClusterType == 3) ){
2686  fEMCALRecUtils->RecalculateClusterPosition(fGeomEMCAL, cells, cluster1);
2687  fEMCALRecUtils->RecalculateClusterPosition(fGeomEMCAL, cells, cluster2);
2688  }
2689 }
2690 
2691 //________________________________________________________________________
2692 Bool_t AliCaloPhotonCuts::CheckDistanceToBadChannel(AliVCluster* cluster, AliVEvent* event)
2693 {
2694  if(fUseDistanceToBadChannel != 1 && fUseDistanceToBadChannel != 2) return kFALSE;
2695 
2696  //not yet fully implemented for PHOS:
2697  if( fClusterType == 2 ) return kFALSE;
2698 
2699  if( (fClusterType == 1 || fClusterType == 3) && !fEMCALInitialized ) InitializeEMCAL(event);
2700  if( fClusterType == 2 && ( !fPHOSInitialized || (fPHOSCurrentRun != event->GetRunNumber()) ) ) InitializePHOS(event);
2701 
2702  Int_t largestCellID = FindLargestCellInCluster(cluster,event);
2703  if(largestCellID==-1) AliFatal("CheckDistanceToBadChannel: FindLargestCellInCluster found cluster with NCells<1?");
2704 
2705  Int_t largestCellicol = -1, largestCellirow = -1;
2706  Int_t rowdiff = 0, coldiff = 0;
2707 
2708  Int_t largestCelliMod = GetModuleNumberAndCellPosition(largestCellID, largestCellicol, largestCellirow);
2709  if(largestCelliMod < 0) AliFatal("CheckDistanceToBadChannel: GetModuleNumberAndCellPosition found SM with ID<0?");
2710 
2711  Int_t nMinRows = 0, nMaxRows = 0;
2712  Int_t nMinCols = 0, nMaxCols = 0;
2713 
2714  Bool_t checkNextSM = kFALSE;
2715  Int_t distanceForLoop = fMinDistanceToBadChannel+1;
2716  if( fClusterType == 1 ){
2717  nMinRows = largestCellirow - distanceForLoop;
2718  nMaxRows = largestCellirow + distanceForLoop;
2719  if(nMinRows < 0) nMinRows = 0;
2720  if(nMaxRows > AliEMCALGeoParams::fgkEMCALRows) nMaxRows = AliEMCALGeoParams::fgkEMCALRows;
2721 
2722  nMinCols = largestCellicol - distanceForLoop;
2723  nMaxCols = largestCellicol + distanceForLoop;
2724 
2725  if(largestCelliMod%2){
2726  if(nMinCols < 0){
2727  nMinCols = 0;
2728  checkNextSM = kTRUE;
2729  }
2730  if(nMaxCols > AliEMCALGeoParams::fgkEMCALCols) nMaxCols = AliEMCALGeoParams::fgkEMCALCols;
2731  }else{
2732  if(nMinCols < 0) nMinCols = 0;
2733  if(nMaxCols > AliEMCALGeoParams::fgkEMCALCols){
2734  nMaxCols = AliEMCALGeoParams::fgkEMCALCols;
2735  checkNextSM = kTRUE;
2736  }
2737  }
2738  }else if( fClusterType == 3 ){
2739  nMinRows = largestCellirow - distanceForLoop;
2740  nMaxRows = largestCellirow + distanceForLoop;
2741  if(nMinRows < 0) nMinRows = 0;
2742  if(nMaxRows > AliEMCALGeoParams::fgkEMCALCols) nMaxRows = AliEMCALGeoParams::fgkEMCALCols; //AliEMCALGeoParams::fgkDCALRows; <- doesnt exist yet (DCAl = EMCAL here)
2743 
2744  nMinCols = largestCellicol - distanceForLoop;
2745  nMaxCols = largestCellicol + distanceForLoop;
2746  if(nMinCols < 0) nMinCols = 0;
2747  if(nMaxCols > fgkDCALCols) nMaxCols = fgkDCALCols; // AliEMCALGeoParams::fgkDCALCols; <- doesnt exist yet
2748 
2749  }else if( fClusterType == 2 ){
2750 // nMaxRows = 64;
2751 // nMaxCols = 56;
2752  }
2753 
2754 // cout << "Cluster: " << fClusterType << ",checkNextSM: " << checkNextSM << endl;
2755 // cout << "largestCell: " << largestCellID << ",mod: " << largestCelliMod << ",col: " << largestCellicol << ",row: " << largestCellirow << endl;
2756 // cout << "distanceForLoop: " << distanceForLoop << ",nMinRows: " << nMinRows << ",nMaxRows: " << nMaxRows << ",nMinCols: " << nMinCols << ",nMaxCols: " << nMaxCols << endl;
2757 
2758  //check bad cells within respective SM
2759  for (Int_t irow = nMinRows;irow < nMaxRows;irow++)
2760  {
2761  for (Int_t icol = nMinCols;icol < nMaxCols;icol++)
2762  {
2763  if(irow == largestCellirow && icol == largestCellicol) continue;
2764 
2765  Int_t iBadCell = 0;
2766  if( (fClusterType == 1 || fClusterType == 3) && largestCelliMod<fEMCALBadChannelsMap->GetEntries()){
2767  iBadCell = (Int_t) ((TH2I*)fEMCALBadChannelsMap->At(largestCelliMod))->GetBinContent(icol,irow);
2768  }else if( fClusterType == 2 && fPHOSBadChannelsMap[largestCelliMod+1]){
2769  iBadCell = (Int_t) ((TH2I*)fPHOSBadChannelsMap[largestCelliMod+1])->GetBinContent(icol,irow);
2770  }
2771  //cout << "largestCelliMod: " << largestCelliMod << ",iBadCell: " << iBadCell << ",icol: " << icol << ",irow: " << irow << endl;
2772  if(iBadCell==0) continue;
2773 
2774  rowdiff = TMath::Abs( largestCellirow - irow ) ;
2775  coldiff = TMath::Abs( largestCellicol - icol ) ;
2776  //cout << "rowdiff: " << rowdiff << ",coldiff: " << coldiff << endl;
2777  if(fUseDistanceToBadChannel==1){
2778  if ((coldiff + rowdiff <= fMinDistanceToBadChannel )) return kTRUE;
2779  }else if(fUseDistanceToBadChannel==2){
2780  if (( coldiff <= fMinDistanceToBadChannel ) && ( rowdiff <= fMinDistanceToBadChannel ) && (coldiff + rowdiff <= fMinDistanceToBadChannel*2)) return kTRUE;
2781  }
2782  //cout << "not within distanceToBadChannel!" << endl;
2783  }
2784  }
2785 
2786  //check bad cells in neighboring SM only if within chosen distanceToBadChannel from maxEnergyCell the next SM could be reached
2787  if(checkNextSM) {
2788  // In case of a shared cluster, index of SM in C side, columns start at 48 and ends at 48*2-1
2789  // C Side impair SM, nSupMod%2=1;A side pair SM nSupMod%2=0
2790  if( fClusterType == 1 ){
2791  if(largestCelliMod%2){
2792  nMinCols = largestCellicol - distanceForLoop + AliEMCALGeoParams::fgkEMCALCols;
2793  nMaxCols = AliEMCALGeoParams::fgkEMCALCols;
2794 
2795  largestCelliMod -= 1;
2796  largestCellicol += AliEMCALGeoParams::fgkEMCALCols;
2797  }else{
2798  nMinCols = 0;
2799  nMaxCols = largestCellicol + distanceForLoop - AliEMCALGeoParams::fgkEMCALCols;
2800 
2801  largestCelliMod += 1;
2802  largestCellicol -= AliEMCALGeoParams::fgkEMCALCols;
2803  }
2804  }else if( fClusterType == 2 ){
2805 // nMaxRows = 64;
2806 // nMaxCols = 56;
2807  }
2808  //cout << "largestCell: " << largestCellID << ",mod: " << largestCelliMod << ",col: " << largestCellicol << ",row: " << largestCellirow << endl;
2809  //cout << "distanceForLoop: " << distanceForLoop << ",nMinRows: " << nMinRows << ",nMaxRows: " << nMaxRows << ",nMinCols: " << nMinCols << ",nMaxCols: " << nMaxCols << endl;
2810  for (Int_t irow = nMinRows;irow < nMaxRows;irow++)
2811  {
2812  for (Int_t icol = nMinCols;icol < nMaxCols;icol++)
2813  {
2814  Int_t iBadCell = 0;
2815  if( fClusterType == 1 && largestCelliMod<fEMCALBadChannelsMap->GetEntries()){
2816  iBadCell = (Int_t) ((TH2I*)fEMCALBadChannelsMap->At(largestCelliMod))->GetBinContent(icol,irow);
2817  }else if( fClusterType == 2 && fPHOSBadChannelsMap[largestCelliMod+1]){
2818  iBadCell = (Int_t) ((TH2I*)fPHOSBadChannelsMap[largestCelliMod+1])->GetBinContent(icol,irow);
2819  }
2820  //cout << "largestCelliMod: " << largestCelliMod << ",iBadCell: " << iBadCell << ",icol: " << icol << ",irow: " << irow << endl;
2821  if(iBadCell==0) continue;
2822 
2823  rowdiff = TMath::Abs( largestCellirow - irow ) ;
2824  coldiff = TMath::Abs( largestCellicol - icol ) ;
2825  //cout << "rowdiff: " << rowdiff << ",coldiff: " << coldiff << endl;
2826  if(fUseDistanceToBadChannel==1){
2827  if ((coldiff + rowdiff <= fMinDistanceToBadChannel )) return kTRUE;
2828  }else if(fUseDistanceToBadChannel==2){
2829  if (( coldiff <= fMinDistanceToBadChannel ) && ( rowdiff <= fMinDistanceToBadChannel ) && (coldiff + rowdiff <= fMinDistanceToBadChannel*2)) return kTRUE;
2830  }
2831  //cout << "not within distanceToBadChannel!" << endl;
2832  }
2833  }
2834  }
2835 
2836  return kFALSE;
2837 }
2838 
2839 
2840 //________________________________________________________________________
2841 Bool_t AliCaloPhotonCuts::ClusterIsSelected(AliVCluster *cluster, AliVEvent * event, AliMCEvent * mcEvent, Int_t isMC, Double_t weight, Long_t clusterID)
2842 {
2843  //Selection of Reconstructed photon clusters with Calorimeters
2844  fIsAcceptedForBasic = kFALSE;
2845  FillClusterCutIndex(kPhotonIn);
2846 
2847 // Double_t vertex[3] = {0,0,0};
2848 // event->GetPrimaryVertex()->GetXYZ(vertex);
2849  // TLorentzvector with cluster
2850 // TLorentzVector clusterVector;
2851 // cluster->GetMomentum(clusterVector,vertex);
2852 
2853  Float_t clusPos[3]={0,0,0};
2854  cluster->GetPosition(clusPos);
2855  TVector3 clusterVector(clusPos[0],clusPos[1],clusPos[2]);
2856  Double_t etaCluster = clusterVector.Eta();
2857  Double_t phiCluster = clusterVector.Phi();
2858  if (phiCluster < 0) phiCluster += 2*TMath::Pi();
2859 
2860  // Histos before cuts
2861  if(fHistClusterEtavsPhiBeforeAcc) fHistClusterEtavsPhiBeforeAcc->Fill(phiCluster,etaCluster,weight);
2862 
2863  // Cluster Selection - 0= accept any calo cluster
2864  if (fClusterType > 0){
2865  //Select EMCAL cluster
2866  if ( (fClusterType == 1 || fClusterType == 3) && !cluster->IsEMCAL()){
2867  FillClusterCutIndex(kDetector);
2868  return kFALSE;
2869  }
2870  //Select PHOS cluster
2871  if (fClusterType == 2 && !( cluster->GetType() == AliVCluster::kPHOSNeutral)){
2872  FillClusterCutIndex(kDetector);
2873  return kFALSE;
2874  }
2875  // do NonLinearity if switched on
2876  if(fUseNonLinearity){
2877  if(fHistEnergyOfClusterBeforeNL) fHistEnergyOfClusterBeforeNL->Fill(cluster->E(),weight);
2878  ApplyNonLinearity(cluster,isMC);
2879  if(fHistEnergyOfClusterAfterNL) fHistEnergyOfClusterAfterNL->Fill(cluster->E(),weight);
2880  }
2881  }
2882 
2883  // Acceptance Cuts
2884  if(!AcceptanceCuts(cluster,event,weight)){
2885  FillClusterCutIndex(kAcceptance);
2886  return kFALSE;
2887  }
2888  // Cluster Quality Cuts
2889  if(!ClusterQualityCuts(cluster,event,mcEvent,isMC,weight,clusterID)){
2890  FillClusterCutIndex(kClusterQuality);
2891  return kFALSE;
2892  }
2893 
2894  // Photon passed cuts
2895  FillClusterCutIndex(kPhotonOut);
2896  return kTRUE;
2897 }
2898 
2899 
2900 //________________________________________________________________________
2901 Bool_t AliCaloPhotonCuts::AcceptanceCuts(AliVCluster *cluster, AliVEvent* event, Double_t weight)
2902 {
2903  // Exclude certain areas for photon reconstruction
2904 
2905  Int_t cutIndex=0;
2906  if(fHistAcceptanceCuts)fHistAcceptanceCuts->Fill(cutIndex);
2907  cutIndex++;
2908 
2909  Float_t clusPos[3]={0,0,0};
2910  cluster->GetPosition(clusPos);
2911  TVector3 clusterVector(clusPos[0],clusPos[1],clusPos[2]);
2912  Double_t etaCluster = clusterVector.Eta();
2913  Double_t phiCluster = clusterVector.Phi();
2914  if (phiCluster < 0) phiCluster += 2*TMath::Pi();
2915 
2916  // check eta range
2917  if (fUseEtaCut){
2918  if (etaCluster < fMinEtaCut || etaCluster > fMaxEtaCut || (fClusterType == 3 && etaCluster < fMaxEtaInnerEdge && etaCluster > fMinEtaInnerEdge ) ){
2919  if(fHistAcceptanceCuts)fHistAcceptanceCuts->Fill(cutIndex);
2920  return kFALSE;
2921  }
2922  }
2923  cutIndex++;
2924 
2925  // check phi range
2926  if (fUsePhiCut ){
2927  if (phiCluster < fMinPhiCut || phiCluster > fMaxPhiCut){
2928  if(fHistAcceptanceCuts)fHistAcceptanceCuts->Fill(cutIndex);
2929  return kFALSE;
2930  }
2931  }
2932  cutIndex++;
2933 
2934  // check distance to bad channel
2935  if (fUseDistanceToBadChannel>0){
2936  if (CheckDistanceToBadChannel(cluster,event)){
2937  if(fHistAcceptanceCuts)fHistAcceptanceCuts->Fill(cutIndex);
2938  return kFALSE;
2939  }
2940  }
2941  //alternatively check histogram fHistoModifyAcc if cluster should be rejected
2942  if(fHistoModifyAcc){
2943  if(fHistoModifyAcc->GetBinContent(FindLargestCellInCluster(cluster,event)) < 1){
2944  if(fHistAcceptanceCuts)fHistAcceptanceCuts->Fill(cutIndex);
2945  return kFALSE;
2946  }
2947  }
2948  cutIndex++;
2949  if(fHistAcceptanceCuts)fHistAcceptanceCuts->Fill(cutIndex);
2950 
2951  // Histos after cuts
2952  if(fHistClusterEtavsPhiAfterAcc) fHistClusterEtavsPhiAfterAcc->Fill(phiCluster,etaCluster,weight);
2953 
2954  return kTRUE;
2955 }
2956 
2957 //________________________________________________________________________
2958 Bool_t AliCaloPhotonCuts::MatchConvPhotonToCluster(AliAODConversionPhoton* convPhoton, AliVCluster* cluster, AliVEvent* event, Double_t weight){
2959 
2960  if (!fUseDistTrackToCluster) return kFALSE;
2961  if( (fClusterType == 1 || fClusterType == 3) && !fEMCALInitialized ) InitializeEMCAL(event);
2962  if( fClusterType == 2 && ( !fPHOSInitialized || (fPHOSCurrentRun != event->GetRunNumber()) ) ) InitializePHOS(event);
2963 
2964  AliESDEvent *esdev = dynamic_cast<AliESDEvent*>(event);
2965  AliAODEvent *aodev = 0;
2966  if (!esdev) {
2967  aodev = dynamic_cast<AliAODEvent*>(event);
2968  if (!aodev) {
2969  AliError("Task needs AOD or ESD event, returning");
2970  return kFALSE;
2971  }
2972  }
2973 
2974  if(!cluster->IsEMCAL() && !cluster->IsPHOS()){AliError("Cluster is neither EMCAL nor PHOS, returning");return kFALSE;}
2975 
2976  Float_t clusterPosition[3] = {0,0,0};
2977  cluster->GetPosition(clusterPosition);
2978  Double_t clusterR = TMath::Sqrt( clusterPosition[0]*clusterPosition[0] + clusterPosition[1]*clusterPosition[1] );
2979  if(fHistClusterRBeforeQA) fHistClusterRBeforeQA->Fill(clusterR,weight);
2980 
2981 //cout << "+++++++++ Cluster: x, y, z, R" << clusterPosition[0] << ", " << clusterPosition[1] << ", " << clusterPosition[2] << ", " << clusterR << "+++++++++" << endl;
2982 
2983  Bool_t matched = kFALSE;
2984  for (Int_t i = 0;i < 2;i++){
2985  Int_t tracklabel = convPhoton->GetLabel(i);
2986  AliVTrack *inTrack = 0x0;
2987  if(esdev) {
2988  if(tracklabel > event->GetNumberOfTracks() ) continue;
2989  inTrack = esdev->GetTrack(tracklabel);
2990  } else {
2991  if(((AliV0ReaderV1*)AliAnalysisManager::GetAnalysisManager()->GetTask(fV0ReaderName.Data()))->AreAODsRelabeled()){
2992  inTrack = dynamic_cast<AliVTrack*>(event->GetTrack(tracklabel));
2993  } else {
2994  for(Int_t ii=0;ii<event->GetNumberOfTracks();ii++) {
2995  inTrack = dynamic_cast<AliVTrack*>(event->GetTrack(ii));
2996  if(inTrack){
2997  if(inTrack->GetID() == tracklabel) {
2998  break;
2999  }
3000  }
3001  }
3002  }
3003  }
3004 
3005  Float_t dEta = 0;
3006  Float_t dPhi = 0;
3007  Bool_t propagated = fCaloTrackMatcher->PropagateV0TrackToClusterAndGetMatchingResidual(inTrack,cluster,event,dEta,dPhi);
3008  if (propagated){
3009  Float_t dR2 = dPhi*dPhi + dEta*dEta;
3010  if(fHistDistanceTrackToClusterBeforeQA)fHistDistanceTrackToClusterBeforeQA->Fill(TMath::Sqrt(dR2), weight);
3011  if(fHistClusterdEtadPhiBeforeQA) fHistClusterdEtadPhiBeforeQA->Fill(dEta, dPhi, weight);
3012 
3013  Float_t clusM02 = (Float_t) cluster->GetM02();
3014  Float_t clusM20 = (Float_t) cluster->GetM20();
3015  if(!fDoLightOutput && (fExtendedMatchAndQA == 1 || fExtendedMatchAndQA == 3 || fExtendedMatchAndQA == 5 )){
3016  if(inTrack->Charge() > 0) {
3017  fHistClusterdEtadPhiPosTracksBeforeQA->Fill(dEta, dPhi, weight);
3018  fHistClusterdPhidPtPosTracksBeforeQA->Fill(dPhi, inTrack->Pt(), weight);
3019  if(inTrack->P() < 0.75) fHistClusterdEtadPhiPosTracksP_000_075BeforeQA->Fill(dEta, dPhi, weight);
3020  else if(inTrack->P() < 1.25) fHistClusterdEtadPhiPosTracksP_075_125BeforeQA->Fill(dEta, dPhi, weight);
3021  else fHistClusterdEtadPhiPosTracksP_125_999BeforeQA->Fill(dEta, dPhi, weight);
3022  } else {
3023  fHistClusterdEtadPhiNegTracksBeforeQA->Fill(dEta, dPhi, weight);
3024  fHistClusterdPhidPtNegTracksBeforeQA->Fill(dPhi, inTrack->Pt(), weight);
3025  if(inTrack->P() < 0.75) fHistClusterdEtadPhiNegTracksP_000_075BeforeQA->Fill(dEta, dPhi, weight);
3026  else if(inTrack->P() < 1.25) fHistClusterdEtadPhiNegTracksP_075_125BeforeQA->Fill(dEta, dPhi, weight);
3027  else fHistClusterdEtadPhiNegTracksP_125_999BeforeQA->Fill(dEta, dPhi, weight);
3028  }
3029  fHistClusterdEtadPtBeforeQA->Fill(dEta, inTrack->Pt(), weight);
3030  fHistClusterM20M02BeforeQA->Fill(clusM20, clusM02, weight);
3031  if(fCurrentMC != kNoMC && fIsMC > 0){
3032  Int_t clusterMCLabel = cluster->GetLabel();
3033  Int_t convPhotonDaughterLabel = -1;
3034  if(inTrack->Charge() > 0) convPhotonDaughterLabel = convPhoton->GetMCLabelPositive();
3035  else convPhotonDaughterLabel = convPhoton->GetMCLabelNegative();
3036  if( (convPhotonDaughterLabel != -1) && (clusterMCLabel != -1) && (convPhotonDaughterLabel == clusterMCLabel)){ //real match
3037  fHistClusterdEtadPtTrueMatched->Fill(dEta, inTrack->Pt(), weight);
3038  if(inTrack->Charge() > 0) fHistClusterdPhidPtPosTracksTrueMatched->Fill(dPhi, inTrack->Pt(), weight);
3039  else fHistClusterdPhidPtNegTracksTrueMatched->Fill(dPhi, inTrack->Pt(), weight);
3040  }
3041  }
3042  }
3043 
3044  Bool_t match_dEta = (TMath::Abs(dEta) < fMaxDistTrackToClusterEta) ? kTRUE : kFALSE;
3045  Bool_t match_dPhi = kFALSE;
3046  if( (inTrack->Charge() > 0) && (dPhi > fMinDistTrackToClusterPhi) && (dPhi < fMaxDistTrackToClusterPhi) ) match_dPhi = kTRUE;
3047  else if( (inTrack->Charge() < 0) && (dPhi < -fMinDistTrackToClusterPhi) && (dPhi > -fMaxDistTrackToClusterPhi) ) match_dPhi = kTRUE;
3048 
3049  if(fUsePtDepTrackToCluster){
3050  if( TMath::Abs(dEta) < fFuncPtDepEta->Eval(inTrack->Pt())) match_dEta = kTRUE;
3051  else match_dEta = kFALSE;
3052 
3053  if( TMath::Abs(dPhi) < fFuncPtDepPhi->Eval(inTrack->Pt())) match_dPhi = kTRUE;
3054  else match_dPhi = kFALSE;
3055  }
3056 
3057  if(match_dEta && match_dPhi){
3058  //if(dR2 < fMinDistTrackToCluster*fMinDistTrackToCluster){
3059  matched = kTRUE;
3060  if(fHistClusterdEtadPtAfterQA) fHistClusterdEtadPtAfterQA->Fill(dEta,inTrack->Pt());
3061  if(fHistClusterdPhidPtAfterQA) fHistClusterdPhidPtAfterQA->Fill(dPhi,inTrack->Pt());
3062  } else {
3063  if(fHistDistanceTrackToClusterAfterQA)fHistDistanceTrackToClusterAfterQA->Fill(TMath::Sqrt(dR2), weight);
3064  if(fHistClusterdEtadPhiAfterQA) fHistClusterdEtadPhiAfterQA->Fill(dEta, dPhi, weight);
3065  if(fHistClusterRAfterQA) fHistClusterRAfterQA->Fill(clusterR, weight);
3066  if(!fDoLightOutput && (fExtendedMatchAndQA == 1 || fExtendedMatchAndQA == 3 || fExtendedMatchAndQA == 5 )){
3067  if(inTrack->Charge() > 0) fHistClusterdEtadPhiPosTracksAfterQA->Fill(dEta, dPhi, weight);
3068  else fHistClusterdEtadPhiNegTracksAfterQA->Fill(dEta, dPhi, weight);
3069  fHistClusterM20M02AfterQA->Fill(clusM20, clusM02, weight);
3070  }
3071  }
3072  }
3073  }
3074 
3075  return matched;
3076 
3077 }
3078 
3079 //________________________________________________________________________
3080 void AliCaloPhotonCuts::MatchTracksToClusters(AliVEvent* event, Double_t weight, Bool_t isEMCalOnly){
3081  if( !fUseDistTrackToCluster ) return;
3082  if( (fClusterType == 1 || fClusterType == 3) && !fEMCALInitialized ) InitializeEMCAL(event);
3083  if( fClusterType == 2 && ( !fPHOSInitialized || (fPHOSCurrentRun != event->GetRunNumber()) ) ) InitializePHOS(event);
3084 
3085  // not yet fully implemented + tested for PHOS
3086  // if( fClusterType != 1 && fClusterType != 3) return;
3087 
3088  fVectorMatchedClusterIDs.clear();
3089 
3090  Int_t nClus = 0;
3091  TClonesArray * arrClustersMatch = NULL;
3092  if(!fCorrTaskSetting.CompareTo("")){
3093  nClus = event->GetNumberOfCaloClusters();
3094  } else {
3095  arrClustersMatch = dynamic_cast<TClonesArray*>(event->FindListObject(Form("%sClustersBranch",fCorrTaskSetting.Data())));
3096  if(!arrClustersMatch)
3097  AliFatal(Form("%sClustersBranch was not found in AliCaloPhotonCuts::FillHistogramsExtendedQA! Check the correction framework settings!",fCorrTaskSetting.Data()));
3098  nClus = arrClustersMatch->GetEntries();
3099  }
3100  //Int_t nModules = 0;
3101 
3102  if(fClusterType == 1 || fClusterType == 3){
3103  fGeomEMCAL = AliEMCALGeometry::GetInstance();
3104  if(!fGeomEMCAL){ AliFatal("EMCal geometry not initialized!");}
3105  //nModules = fGeomEMCAL->GetNumberOfSuperModules();
3106  }else if(fClusterType == 2){
3107  fGeomPHOS = AliPHOSGeometry::GetInstance();
3108  if(!fGeomPHOS) AliFatal("PHOS geometry not initialized!");
3109  //nModules = fGeomPHOS->GetNModules();
3110  }
3111 
3112  AliESDEvent *esdev = dynamic_cast<AliESDEvent*>(event);
3113  AliAODEvent *aodev = 0;
3114  if (!esdev) {
3115  aodev = dynamic_cast<AliAODEvent*>(event);
3116  if (!aodev) {
3117  AliError("Task needs AOD or ESD event, returning");
3118  return;
3119  }
3120  }
3121 
3122  // if not EMCal only reconstruction (-> hybrid PCM+EMCal), use only primary tracks for basic track matching procedure
3123  AliESDtrackCuts *EsdTrackCuts = 0x0;
3124  if(!isEMCalOnly && esdev){
3125  // Using standard function for setting Cuts
3126  Int_t runNumber = event->GetRunNumber();
3127  // if LHC11a or earlier or if LHC13g or if LHC12a-i -> use 2010 cuts
3128  if( (runNumber<=146860) || (runNumber>=197470 && runNumber<=197692) || (runNumber>=172440 && runNumber<=193766) ){
3129  EsdTrackCuts = AliESDtrackCuts::GetStandardITSTPCTrackCuts2010();
3130  // else if run2 data use 2015 PbPb cuts
3131  }else if (runNumber>=209122){
3132  // EsdTrackCuts = AliESDtrackCuts::GetStandardITSTPCTrackCuts2015PbPb();
3133  // hard coded track cuts for the moment, because AliESDtrackCuts::GetStandardITSTPCTrackCuts2015PbPb() gives spams warnings
3134  EsdTrackCuts = new AliESDtrackCuts();
3135  EsdTrackCuts->AliESDtrackCuts::SetMinNCrossedRowsTPC(70);
3136  EsdTrackCuts->AliESDtrackCuts::SetMinRatioCrossedRowsOverFindableClustersTPC(0.8);
3137  EsdTrackCuts->AliESDtrackCuts::SetCutOutDistortedRegionsTPC(kTRUE);
3138  EsdTrackCuts->AliESDtrackCuts::SetMaxChi2PerClusterTPC(4);
3139  EsdTrackCuts->AliESDtrackCuts::SetAcceptKinkDaughters(kFALSE);
3140  EsdTrackCuts->AliESDtrackCuts::SetRequireTPCRefit(kTRUE);
3141  // ITS
3142  EsdTrackCuts->AliESDtrackCuts::SetRequireITSRefit(kTRUE);
3143  EsdTrackCuts->AliESDtrackCuts::SetClusterRequirementITS(AliESDtrackCuts::kSPD,
3144  AliESDtrackCuts::kAny);
3145  EsdTrackCuts->AliESDtrackCuts::SetMaxDCAToVertexXYPtDep("0.0105+0.0350/pt^1.1");
3146  EsdTrackCuts->AliESDtrackCuts::SetMaxChi2TPCConstrainedGlobal(36);
3147  EsdTrackCuts->AliESDtrackCuts::SetMaxDCAToVertexZ(2);
3148  EsdTrackCuts->AliESDtrackCuts::SetDCAToVertex2D(kFALSE);
3149  EsdTrackCuts->AliESDtrackCuts::SetRequireSigmaToVertex(kFALSE);
3150  EsdTrackCuts->AliESDtrackCuts::SetMaxChi2PerClusterITS(36);
3151  // else use 2011 version of track cuts
3152  }else{
3153  EsdTrackCuts = AliESDtrackCuts::GetStandardITSTPCTrackCuts2011();
3154  }
3155  EsdTrackCuts->SetMaxDCAToVertexZ(2);
3156  EsdTrackCuts->SetEtaRange(-0.8, 0.8);
3157  EsdTrackCuts->SetPtRange(0.15);
3158  }
3159 
3160 // cout << "MatchTracksToClusters: " << event->GetNumberOfTracks() << ", " << fIsPureCalo << ", " << fUseDistTrackToCluster << endl;
3161 
3162  for (Int_t itr=0;itr<event->GetNumberOfTracks();itr++){
3163  AliExternalTrackParam *trackParam = 0;
3164  AliVTrack *inTrack = 0x0;
3165  if(esdev){
3166  inTrack = esdev->GetTrack(itr);
3167  if(!inTrack) continue;
3168  AliESDtrack *esdt = dynamic_cast<AliESDtrack*>(inTrack);
3169  if(!isEMCalOnly){ //match only primaries for hybrid reconstruction schemes
3170  if(!EsdTrackCuts->AcceptTrack(esdt)) continue;
3171  }
3172 
3173  const AliExternalTrackParam *in = esdt->GetInnerParam();
3174  if (!in){AliDebug(2, "Could not get InnerParam of Track, continue");continue;}
3175  trackParam = new AliExternalTrackParam(*in);
3176  } else if(aodev) {
3177  inTrack = dynamic_cast<AliVTrack*>(aodev->GetTrack(itr));
3178  if(!inTrack) continue;
3179  AliAODTrack *aodt = dynamic_cast<AliAODTrack*>(inTrack);
3180  if(!isEMCalOnly){ //match only primaries for hybrid reconstruction schemes
3181  if(!aodt->IsHybridGlobalConstrainedGlobal()) continue;
3182  if(TMath::Abs(aodt->Eta())>0.8) continue;
3183  if(aodt->Pt()<0.15) continue;
3184  }
3185 
3186  }
3187 
3188  Float_t clsPos[3] = {0.,0.,0.};
3189  for(Int_t iclus=0;iclus < nClus;iclus++){
3190  AliVCluster * cluster = NULL;
3191  if(arrClustersMatch){
3192  if(esdev)
3193  cluster = new AliESDCaloCluster(*(AliESDCaloCluster*)arrClustersMatch->At(iclus));
3194  else if(aodev)
3195  cluster = new AliAODCaloCluster(*(AliAODCaloCluster*)arrClustersMatch->At(iclus));
3196  } else {
3197  cluster = event->GetCaloCluster(iclus);
3198  }
3199 
3200  if (!cluster) continue;
3201  Float_t dEta, dPhi;
3202  if(!fCaloTrackMatcher->GetTrackClusterMatchingResidual(inTrack->GetID(),cluster->GetID(),dEta,dPhi)) continue;
3203  cluster->GetPosition(clsPos);
3204  Float_t clusterR = TMath::Sqrt( clsPos[0]*clsPos[0] + clsPos[1]*clsPos[1] );
3205  Float_t dR2 = dPhi*dPhi + dEta*dEta;
3206 
3207  if(isEMCalOnly && fHistDistanceTrackToClusterBeforeQA)fHistDistanceTrackToClusterBeforeQA->Fill(TMath::Sqrt(dR2), weight);
3208  if(isEMCalOnly && fHistClusterdEtadPhiBeforeQA) fHistClusterdEtadPhiBeforeQA->Fill(dEta, dPhi, weight);
3209  if(isEMCalOnly && fHistClusterRBeforeQA) fHistClusterRBeforeQA->Fill(clusterR, weight);
3210 
3211  Float_t clusM02 = (Float_t) cluster->GetM02();
3212  Float_t clusM20 = (Float_t) cluster->GetM20();
3213  if(isEMCalOnly && !fDoLightOutput && (fExtendedMatchAndQA == 1 || fExtendedMatchAndQA == 3 || fExtendedMatchAndQA == 5 )){
3214  if(inTrack->Charge() > 0) {
3215  fHistClusterdEtadPhiPosTracksBeforeQA->Fill(dEta, dPhi, weight);
3216  fHistClusterdPhidPtPosTracksBeforeQA->Fill(dPhi, inTrack->Pt(), weight);
3217  if(inTrack->P() < 0.75) fHistClusterdEtadPhiPosTracksP_000_075BeforeQA->Fill(dEta, dPhi, weight);
3218  else if(inTrack->P() < 1.25) fHistClusterdEtadPhiPosTracksP_075_125BeforeQA->Fill(dEta, dPhi, weight);
3219  else fHistClusterdEtadPhiPosTracksP_125_999BeforeQA->Fill(dEta, dPhi, weight);
3220  }
3221  else{
3222  fHistClusterdEtadPhiNegTracksBeforeQA->Fill(dEta, dPhi, weight);
3223  fHistClusterdPhidPtNegTracksBeforeQA->Fill(dPhi, inTrack->Pt(), weight);
3224  if(inTrack->P() < 0.75) fHistClusterdEtadPhiNegTracksP_000_075BeforeQA->Fill(dEta, dPhi, weight);
3225  else if(inTrack->P() < 1.25) fHistClusterdEtadPhiNegTracksP_075_125BeforeQA->Fill(dEta, dPhi, weight);
3226  else fHistClusterdEtadPhiNegTracksP_125_999BeforeQA->Fill(dEta, dPhi, weight);
3227  }
3228  fHistClusterdEtadPtBeforeQA->Fill(dEta, inTrack->Pt(), weight);
3229  fHistClusterM20M02BeforeQA->Fill(clusM20, clusM02, weight);
3230  }
3231 
3232  Bool_t match_dEta = (TMath::Abs(dEta) < fMaxDistTrackToClusterEta) ? kTRUE : kFALSE;
3233  Bool_t match_dPhi = kFALSE;
3234 
3235  if( (inTrack->Charge() > 0) && (dPhi > fMinDistTrackToClusterPhi) && (dPhi < fMaxDistTrackToClusterPhi) ) match_dPhi = kTRUE;
3236  else if( (inTrack->Charge() < 0) && (dPhi < -fMinDistTrackToClusterPhi) && (dPhi > -fMaxDistTrackToClusterPhi) ) match_dPhi = kTRUE;
3237 
3238  if(fUsePtDepTrackToCluster){
3239  if( TMath::Abs(dEta) < fFuncPtDepEta->Eval(inTrack->Pt())) match_dEta = kTRUE;
3240  else match_dEta = kFALSE;
3241 
3242  if( TMath::Abs(dPhi) < fFuncPtDepPhi->Eval(inTrack->Pt())) match_dPhi = kTRUE;
3243  else match_dPhi = kFALSE;
3244  }
3245 
3246  if(match_dEta && match_dPhi){
3247  fVectorMatchedClusterIDs.push_back(cluster->GetID());
3248  // cout << "MATCHED!!!!!!!!!!!!!!!!!!!!!!!!! - " << cluster->GetID() << endl;
3249  if(isEMCalOnly){
3250  if(fHistClusterdEtadPtAfterQA) fHistClusterdEtadPtAfterQA->Fill(dEta,inTrack->Pt());
3251  if(fHistClusterdPhidPtAfterQA) fHistClusterdPhidPtAfterQA->Fill(dPhi,inTrack->Pt());
3252  }
3253  break;
3254  } else if(isEMCalOnly){
3255  if(fHistDistanceTrackToClusterAfterQA)fHistDistanceTrackToClusterAfterQA->Fill(TMath::Sqrt(dR2), weight);
3256  if(fHistClusterdEtadPhiAfterQA) fHistClusterdEtadPhiAfterQA->Fill(dEta, dPhi, weight);
3257  if(fHistClusterRAfterQA) fHistClusterRAfterQA->Fill(clusterR, weight);
3258  if(!fDoLightOutput && (fExtendedMatchAndQA == 1 || fExtendedMatchAndQA == 3 || fExtendedMatchAndQA == 5 )){
3259  if(inTrack->Charge() > 0) fHistClusterdEtadPhiPosTracksAfterQA->Fill(dEta, dPhi, weight);
3260  else fHistClusterdEtadPhiNegTracksAfterQA->Fill(dEta, dPhi, weight);
3261  fHistClusterM20M02AfterQA->Fill(clusM20, clusM02, weight);
3262  }
3263  // cout << "no match" << endl;
3264  }
3265  }
3266 
3267  delete trackParam;
3268  }
3269 
3270  if(EsdTrackCuts){
3271  delete EsdTrackCuts;
3272  EsdTrackCuts=0x0;
3273  }
3274 
3275  return;
3276 }
3277 
3278 //________________________________________________________________________
3279 Bool_t AliCaloPhotonCuts::CheckClusterForTrackMatch(AliVCluster* cluster){
3280  vector<Int_t>::iterator it;
3281  it = find (fVectorMatchedClusterIDs.begin(), fVectorMatchedClusterIDs.end(), cluster->GetID());
3282  if (it != fVectorMatchedClusterIDs.end()) return kTRUE;
3283  else return kFALSE;
3284 }
3285 
3286 //________________________________________________________________________
3287 Bool_t AliCaloPhotonCuts::UpdateCutString() {
3289 
3290  if(fCutString && fCutString->GetString().Length() == kNCuts) {
3291  fCutString->SetString(GetCutNumber());
3292  } else {
3293  return kFALSE;
3294  }
3295  return kTRUE;
3296 }
3297 
3298 //________________________________________________________________________
3299 Bool_t AliCaloPhotonCuts::InitializeCutsFromCutString(const TString analysisCutSelection ) {
3300  fCutStringRead = Form("%s",analysisCutSelection.Data());
3301 
3302  // Initialize Cuts from a given Cut string
3303  AliInfo(Form("Set CaloCut Number: %s",analysisCutSelection.Data()));
3304  if(analysisCutSelection.Length()!=kNCuts) {
3305  AliError(Form("Cut selection has the wrong length! size is %d, number of cuts is %d", analysisCutSelection.Length(), kNCuts));
3306  return kFALSE;
3307  }
3308  if(!analysisCutSelection.IsAlnum()){
3309  AliError("Cut selection is not alphanumeric");
3310  return kFALSE;
3311  }
3312 
3313  TString analysisCutSelectionLowerCase = Form("%s",analysisCutSelection.Data());
3314  analysisCutSelectionLowerCase.ToLower();
3315  const char *cutSelection = analysisCutSelectionLowerCase.Data();
3316  #define ASSIGNARRAY(i) fCuts[i] = ((int)cutSelection[i]>=(int)'a') ? cutSelection[i]-'a'+10 : cutSelection[i]-'0'
3317  for(Int_t ii=0;ii<kNCuts;ii++){
3318  ASSIGNARRAY(ii);
3319  }
3320 
3321  // Set Individual Cuts
3322  for(Int_t ii=0;ii<kNCuts;ii++){
3323  if(!SetCut(cutIds(ii),fCuts[ii]))return kFALSE;
3324  }
3325 
3326  PrintCutsWithValues();
3327  return kTRUE;
3328 }
3329 
3330 //________________________________________________________________________
3331 Bool_t AliCaloPhotonCuts::SetCut(cutIds cutID, const Int_t value) {
3333 
3334  switch (cutID) {
3335 
3336  case kClusterType:
3337  if( SetClusterTypeCut(value)) {
3338  fCuts[kClusterType] = value;
3339  UpdateCutString();
3340  return kTRUE;
3341  } else return kFALSE;
3342 
3343  case kEtaMin:
3344  if( SetMinEtaCut(value)) {
3345  fCuts[kEtaMin] = value;
3346  UpdateCutString();
3347  return kTRUE;
3348  } else return kFALSE;
3349 
3350  case kEtaMax:
3351  if( SetMaxEtaCut(value)) {
3352  fCuts[kEtaMax] = value;
3353  UpdateCutString();
3354  return kTRUE;
3355  } else return kFALSE;
3356 
3357  case kPhiMin:
3358  if( SetMinPhiCut(value)) {
3359  fCuts[kPhiMin] = value;
3360  UpdateCutString();
3361  return kTRUE;
3362  } else return kFALSE;
3363 
3364  case kPhiMax:
3365  if( SetMaxPhiCut(value)) {
3366  fCuts[kPhiMax] = value;
3367  UpdateCutString();
3368  return kTRUE;
3369  } else return kFALSE;
3370 
3371  case kDistanceToBadChannel:
3372  if( SetDistanceToBadChannelCut(value)) {
3373  fCuts[kDistanceToBadChannel] = value;
3374  UpdateCutString();
3375  return kTRUE;
3376  } else return kFALSE;
3377 
3378  case kTiming:
3379  if( SetTimingCut(value)) {
3380  fCuts[kTiming] = value;
3381  UpdateCutString();
3382  return kTRUE;
3383  } else return kFALSE;
3384 
3385  case kTrackMatching:
3386  if( SetTrackMatchingCut(value)) {
3387  fCuts[kTrackMatching] = value;
3388  UpdateCutString();
3389  return kTRUE;
3390  } else return kFALSE;
3391 
3392  case kExoticCluster:
3393  if( SetExoticClusterCut(value)) {
3394  fCuts[kExoticCluster] = value;
3395  UpdateCutString();
3396  return kTRUE;
3397  } else return kFALSE;
3398 
3399  case kMinEnergy:
3400  if( SetMinEnergyCut(value)) {
3401  fCuts[kMinEnergy] = value;
3402  UpdateCutString();
3403  return kTRUE;
3404  } else return kFALSE;
3405 
3406  case kNMinCells:
3407  if( SetMinNCellsCut(value)) {
3408  fCuts[kNMinCells] = value;
3409  UpdateCutString();
3410  return kTRUE;
3411  } else return kFALSE;
3412 
3413  case kMinM02:
3414  if( SetMinM02(value)) {
3415  fCuts[kMinM02] = value;
3416  UpdateCutString();
3417  return kTRUE;
3418  } else return kFALSE;
3419 
3420  case kMaxM02:
3421  if( SetMaxM02(value)) {
3422  fCuts[kMaxM02] = value;
3423  UpdateCutString();
3424  return kTRUE;
3425  } else return kFALSE;
3426 
3427  case kMinM20:
3428  if( SetMinM20(value)) {
3429  fCuts[kMinM20] = value;
3430  UpdateCutString();
3431  return kTRUE;
3432  } else return kFALSE;
3433 
3434  case kMaxM20:
3435  if( SetMaxM20(value)) {
3436  fCuts[kMaxM20] = value;
3437  UpdateCutString();
3438  return kTRUE;
3439  } else return kFALSE;
3440 
3441  case kDispersion:
3442  if( SetDispersion(value)) {
3443  fCuts[kDispersion] = value;
3444  UpdateCutString();
3445  return kTRUE;
3446  } else return kFALSE;
3447 
3448  case kNLM:
3449  if( SetNLM(value)) {
3450  fCuts[kNLM] = value;
3451  UpdateCutString();
3452  return kTRUE;
3453  } else return kFALSE;
3454 
3455  case kNonLinearity1:
3456  if( SetNonLinearity1(value)) {
3457  fCuts[kNonLinearity1] = value;
3458  UpdateCutString();
3459  return kTRUE;
3460  } else return kFALSE;
3461 
3462  case kNonLinearity2:
3463  if( SetNonLinearity2(value)) {
3464  fCuts[kNonLinearity2] = value;
3465  UpdateCutString();
3466  return kTRUE;
3467  } else return kFALSE;
3468 
3469  case kNCuts:
3470  AliError("Cut id out of range");
3471  return kFALSE;
3472  }
3473 
3474  AliError("Cut id %d not recognized");
3475  return kFALSE;
3476 
3477 
3478 }
3479 
3480 //________________________________________________________________________
3481 void AliCaloPhotonCuts::PrintCuts() {
3482  // Print out current Cut Selection
3483  for(Int_t ic = 0;ic < kNCuts;ic++) {
3484  printf("%-30s : %d \n", fgkCutNames[ic], fCuts[ic]);
3485  }
3486 }
3487 
3488 //________________________________________________________________________
3489 void AliCaloPhotonCuts::PrintCutsWithValues() {
3490  // Print out current Cut Selection with value
3491  printf("\nCluster cutnumber \n");
3492  for(Int_t ic = 0;ic < kNCuts;ic++) {
3493  printf("%d",fCuts[ic]);
3494  }
3495  printf("\n\n");
3496  if (fIsPureCalo>0) printf("Merged cluster analysis was specified, mode: '%i'\n", fIsPureCalo);
3497 
3498  printf("Acceptance cuts: \n");
3499  if (fClusterType == 0) printf("\tall calorimeter clusters are used\n");
3500  if (fClusterType == 1) printf("\tEMCAL calorimeter clusters are used\n");
3501  if (fClusterType == 2) printf("\tPHOS calorimeter clusters are used\n");
3502  if (fClusterType == 3) printf("\tDCAL calorimeter clusters are used\n");
3503  if (fUseEtaCut) printf("\t%3.2f < eta_{cluster} < %3.2f\n", fMinEtaCut, fMaxEtaCut );
3504  if (fUsePhiCut) printf("\t%3.2f < phi_{cluster} < %3.2f\n", fMinPhiCut, fMaxPhiCut );
3505  if (fUseDistanceToBadChannel>0) printf("\tdistance to bad channel used in mode '%i', distance in cells: %f \n",fUseDistanceToBadChannel, fMinDistanceToBadChannel);
3506 
3507  printf("Cluster Quality cuts: \n");
3508  if (fUseTimeDiff) printf("\t %6.2f ns < time difference < %6.2f ns\n", fMinTimeDiff*1e9, fMaxTimeDiff*1e9 );
3509  if (fUseDistTrackToCluster) printf("\tmin distance to track in eta > %3.2f, min phi < %3.2f and max phi > %3.2f\n", fMaxDistTrackToClusterEta, fMinDistTrackToClusterPhi, fMaxDistTrackToClusterPhi );
3510  if (fUseExoticCluster)printf("\t exotic cluster: %3.2f\n", fExoticEnergyFracCluster );
3511  if (fUseMinEnergy)printf("\t E_{cluster} > %3.2f\n", fMinEnergy );
3512  if (fUseNCells) printf("\t number of cells per cluster >= %d\n", fMinNCells );
3513  if (fUseM02 == 1) printf("\t %3.2f < M02 < %3.2f\n", fMinM02, fMaxM02 );
3514  if (fUseM02 == 2) printf("\t energy dependent M02 cut used with cutnumber min: %d max: %d \n", fMinM02CutNr, fMaxM02CutNr );
3515  if (fUseM20) printf("\t %3.2f < M20 < %3.2f\n", fMinM20, fMaxM20 );
3516  if (fUseDispersion) printf("\t dispersion < %3.2f\n", fMaxDispersion );
3517  if (fUseNLM) printf("\t %d < NLM < %d\n", fMinNLM, fMaxNLM );
3518  printf("Correction Task Setting: %s \n",fCorrTaskSetting.Data());
3519 
3520  printf("NonLinearity Correction: \n");
3521  printf("VO Reader name: %s \n",fV0ReaderName.Data());
3522  TString namePeriod = ((AliV0ReaderV1*)AliAnalysisManager::GetAnalysisManager()->GetTask(fV0ReaderName.Data()))->GetPeriodName();
3523  if (namePeriod.CompareTo("") != 0) fCurrentMC = FindEnumForMCSet(namePeriod);
3524  if (fUseNonLinearity) printf("\t Chose NonLinearity cut '%i', Period name: %s, period-enum: %o \n", fSwitchNonLinearity, namePeriod.Data(), fCurrentMC );
3525  else printf("\t No NonLinearity Correction on AnalysisTask level has been chosen\n");
3526 
3527 }
3528 
3529 // EMCAL acceptance 2011
3530 // 1.39626, 3.125 (phi)
3531 // -0.66687,,0.66465
3532 
3533 
3534 //________________________________________________________________________
3535 Bool_t AliCaloPhotonCuts::SetClusterTypeCut(Int_t clusterType)
3536 { // Set Cut
3537  switch(clusterType){
3538  case 0: // all clusters
3539  fClusterType=0;
3540  break;
3541  case 1: // EMCAL clusters
3542  fClusterType=1;
3543  break;
3544  case 2: // PHOS clusters
3545  fClusterType=2;
3546  break;
3547  case 3: // DCAL clusters
3548  fClusterType=3;
3549  break;
3550  default:
3551  AliError(Form("ClusterTypeCut not defined %d",clusterType));
3552  return kFALSE;
3553  }
3554  return kTRUE;
3555 }
3556 
3557 //___________________________________________________________________
3558 Bool_t AliCaloPhotonCuts::SetMinEtaCut(Int_t minEta)
3559 {
3560  switch(minEta){
3561  case 0:
3562  if (!fUseEtaCut) fUseEtaCut=0;
3563  fMinEtaCut=-10.;
3564  break;
3565  case 1:
3566  if (!fUseEtaCut) fUseEtaCut=1;
3567  fMinEtaCut=-0.6687;
3568  break;
3569  case 2:
3570  if (!fUseEtaCut) fUseEtaCut=1;
3571  fMinEtaCut=-0.5;
3572  break;
3573  case 3:
3574  if (!fUseEtaCut) fUseEtaCut=1;
3575  fMinEtaCut=-2;
3576  break;
3577  case 4:
3578  if (!fUseEtaCut) fUseEtaCut=1;
3579  fMinEtaCut = -0.13;
3580  break;
3581  case 5:
3582  if (!fUseEtaCut) fUseEtaCut=1;
3583  fMinEtaCut=-0.7;
3584  break;
3585  case 6:
3586  if (!fUseEtaCut) fUseEtaCut=1;
3587  fMinEtaCut=-0.3;
3588  break;
3589  case 7:
3590  if (!fUseEtaCut) fUseEtaCut=1;
3591  fMinEtaCut=-0.4;
3592  break;
3593  case 8:
3594  if (!fUseEtaCut) fUseEtaCut=1;
3595  fMinEtaCut=-0.66112;
3596  fMinEtaInnerEdge=-0.227579;
3597  break;
3598 
3599  default:
3600  AliError(Form("MinEta Cut not defined %d",minEta));
3601  return kFALSE;
3602  }
3603  return kTRUE;
3604 }
3605 
3606 
3607 //___________________________________________________________________
3608 Bool_t AliCaloPhotonCuts::SetMaxEtaCut(Int_t maxEta)
3609 {
3610  switch(maxEta){
3611  case 0:
3612  if (!fUseEtaCut) fUseEtaCut=0;
3613  fMaxEtaCut=10;
3614  break;
3615  case 1:
3616  if (!fUseEtaCut) fUseEtaCut=1;
3617  fMaxEtaCut=0.66465;
3618  break;
3619  case 2:
3620  if (!fUseEtaCut) fUseEtaCut=1;
3621  fMaxEtaCut=0.5;
3622  break;
3623  case 3:
3624  if (!fUseEtaCut) fUseEtaCut=1;
3625  fMaxEtaCut=2;
3626  break;
3627  case 4:
3628  if (!fUseEtaCut) fUseEtaCut=1;
3629  fMaxEtaCut= 0.13;
3630  break;
3631  case 5:
3632  if (!fUseEtaCut) fUseEtaCut=1;
3633  fMaxEtaCut=0.7;
3634  break;
3635  case 6:
3636  if (!fUseEtaCut) fUseEtaCut=1;
3637  fMaxEtaCut=0.3;
3638  break;
3639  case 7:
3640  if (!fUseEtaCut) fUseEtaCut=1;
3641  fMaxEtaCut=0.4;
3642  break;
3643  case 8:
3644  if(!fUseEtaCut) fUseEtaCut=1;
3645  fMaxEtaCut=0.66112;
3646  fMaxEtaInnerEdge=0.227579;
3647  break;
3648  default:
3649  AliError(Form("MaxEta Cut not defined %d",maxEta));
3650  return kFALSE;
3651  }
3652  return kTRUE;
3653 }
3654 
3655 //___________________________________________________________________
3656 Bool_t AliCaloPhotonCuts::SetMinPhiCut(Int_t minPhi)
3657 {
3658  switch(minPhi){
3659  case 0:
3660  if (!fUsePhiCut) fUsePhiCut=0;
3661  fMinPhiCut=-10000;
3662  break;
3663  case 1: // min EMCAL
3664  if (!fUsePhiCut) fUsePhiCut=1;
3665  fMinPhiCut=1.39626;
3666  break;
3667  case 2: // min EMCAL with TRD 2012
3668  if (!fUsePhiCut) fUsePhiCut=1;
3669  fMinPhiCut=2.10;
3670  break;
3671  case 3: // min EMCAL with TRD 2011
3672  if (!fUsePhiCut) fUsePhiCut=1;
3673  fMinPhiCut=2.45;
3674  break;
3675  case 4:
3676  if( !fUsePhiCut ) fUsePhiCut=1;
3677  fMinPhiCut = 4.54;//PHOS acceptance
3678  break;
3679  case 5:
3680  if( !fUsePhiCut ) fUsePhiCut=1;
3681  fMinPhiCut = 4.5572;//DCal acceptance
3682  break;
3683  case 6:
3684  if( !fUsePhiCut ) fUsePhiCut=1;
3685  fMinPhiCut = 4.36;//PHOS acceptance RUN2
3686  break;
3687 
3688  default:
3689  AliError(Form("MinPhi Cut not defined %d",minPhi));
3690  return kFALSE;
3691  }
3692  return kTRUE;
3693 }
3694 
3695 //___________________________________________________________________
3696 Bool_t AliCaloPhotonCuts::SetMaxPhiCut(Int_t maxPhi)
3697 {
3698  switch(maxPhi){
3699  case 0:
3700  if (!fUsePhiCut) fUsePhiCut=0;
3701  fMaxPhiCut=10000;
3702  break;
3703  case 1: // max EMCAL
3704  if (!fUsePhiCut) fUsePhiCut=1;
3705  fMaxPhiCut=3.15;
3706  break;
3707  case 2: // max EMCAL with TRD 2011
3708  if (!fUsePhiCut) fUsePhiCut=1;
3709  fMaxPhiCut=2.45;
3710  break;
3711  case 3: // max EMCAL with TRD 2012
3712  if (!fUsePhiCut) fUsePhiCut=1;
3713  fMaxPhiCut=2.10;
3714  break;
3715  case 4:
3716  if( !fUsePhiCut ) fUsePhiCut=1;
3717  fMaxPhiCut = 5.59;//PHOS acceptance
3718  break;
3719  case 5:
3720  if( !fUsePhiCut ) fUsePhiCut=1;
3721  fMaxPhiCut = 5.5658;//DCal acceptance
3722  break;
3723  case 6:
3724  if( !fUsePhiCut ) fUsePhiCut=1;
3725  fMaxPhiCut = 5.59;//PHOS acceptance RUN2
3726  break;
3727 
3728  default:
3729  AliError(Form("Max Phi Cut not defined %d",maxPhi));
3730  return kFALSE;
3731  }
3732  return kTRUE;
3733 }
3734 
3735 //___________________________________________________________________
3736 Bool_t AliCaloPhotonCuts::SetDistanceToBadChannelCut(Int_t distanceToBadChannel)
3737 {
3738  switch(distanceToBadChannel){
3739  case 0:
3740  fUseDistanceToBadChannel=0;
3741  fMinDistanceToBadChannel=0;
3742  break;
3743  case 1:
3744  if(fUseDistanceToBadChannel!=1) fUseDistanceToBadChannel=1;
3745  fMinDistanceToBadChannel=1;
3746  break;
3747  case 2:
3748  if(fUseDistanceToBadChannel!=1) fUseDistanceToBadChannel=1;
3749  fMinDistanceToBadChannel=2;
3750  break;
3751  case 3:
3752  if(fUseDistanceToBadChannel!=1) fUseDistanceToBadChannel=1;
3753  fMinDistanceToBadChannel=3;
3754  break;
3755  case 4:
3756  if(fUseDistanceToBadChannel!=1) fUseDistanceToBadChannel=1;
3757  fMinDistanceToBadChannel=4;
3758  break;
3759  case 5:
3760  if(fUseDistanceToBadChannel!=2) fUseDistanceToBadChannel=2;
3761  fMinDistanceToBadChannel=1;
3762  break;
3763  case 6:
3764  if(fUseDistanceToBadChannel!=2) fUseDistanceToBadChannel=2;
3765  fMinDistanceToBadChannel=2;
3766  break;
3767  case 7:
3768  if(fUseDistanceToBadChannel!=2) fUseDistanceToBadChannel=2;
3769  fMinDistanceToBadChannel=3;
3770  break;
3771  case 8:
3772  if(fUseDistanceToBadChannel!=2) fUseDistanceToBadChannel=2;
3773  fMinDistanceToBadChannel=4;
3774  break;
3775  default:
3776  AliError(Form("minimum distance to bad channel Cut not defined %d",distanceToBadChannel));
3777  return kFALSE;
3778  }
3779  return kTRUE;
3780 }
3781 
3782 //___________________________________________________________________
3783 Bool_t AliCaloPhotonCuts::SetTimingCut(Int_t timing)
3784 {
3785  switch(timing){
3786  case 0:
3787  fUseTimeDiff=0;
3788  fMinTimeDiff=-500;
3789  fMaxTimeDiff=500;
3790  break;
3791  case 1:
3792  if (!fUseTimeDiff) fUseTimeDiff=1;
3793  fMinTimeDiff=-10e-7;
3794  fMaxTimeDiff=10e-7;//1000ns
3795  break;
3796  case 2:
3797  if (!fUseTimeDiff) fUseTimeDiff=1;
3798  fMinTimeDiff=-50e-8;
3799  fMaxTimeDiff=50e-8;//500ns
3800  break;
3801  case 3:
3802  if (!fUseTimeDiff) fUseTimeDiff=1;
3803  fMinTimeDiff=-20e-8;
3804  fMaxTimeDiff=20e-8;//200ns
3805  break;
3806  case 4:
3807  if (!fUseTimeDiff) fUseTimeDiff=1;
3808  fMinTimeDiff=-10e-8;
3809  fMaxTimeDiff=10e-8;//100ns
3810  break;
3811  case 5:
3812  if (!fUseTimeDiff) fUseTimeDiff=1;
3813  fMinTimeDiff=-50e-9;
3814  fMaxTimeDiff=50e-9;//50ns
3815  break;
3816  case 6:
3817  if (!fUseTimeDiff) fUseTimeDiff=1;
3818  fMinTimeDiff=-30e-9;
3819  fMaxTimeDiff=35e-9;
3820  break;
3821  case 7:
3822  if (!fUseTimeDiff) fUseTimeDiff=1;
3823  fMinTimeDiff=-30e-9;
3824  fMaxTimeDiff=30e-9;//30ns
3825  break;
3826  case 8:
3827  if (!fUseTimeDiff) fUseTimeDiff=1;
3828  fMinTimeDiff=-20e-9;
3829  fMaxTimeDiff=30e-9;
3830  break;
3831  case 9:
3832  if (!fUseTimeDiff) fUseTimeDiff=1;
3833  fMinTimeDiff=-20e-9;
3834  fMaxTimeDiff=25e-9;
3835  break;
3836  case 10:
3837  if (!fUseTimeDiff) fUseTimeDiff=1;
3838  fMinTimeDiff=-12.5e-9;
3839  fMaxTimeDiff=13e-9;
3840  break;
3841  case 11:
3842  if (!fUseTimeDiff) fUseTimeDiff=1;
3843  fMinTimeDiff=-130e-9;
3844  fMaxTimeDiff=130e-9;
3845  break;
3846  default:
3847  AliError(Form("Timing Cut not defined %d",timing));
3848  return kFALSE;
3849  }
3850  return kTRUE;
3851 }
3852 
3853 //___________________________________________________________________
3854 Bool_t AliCaloPhotonCuts::SetTrackMatchingCut(Int_t trackMatching)
3855 {
3856  // matching parameters for EMCal clusters
3857  if(fClusterType == 1 || fClusterType == 3){
3858  switch(trackMatching){
3859  case 0:
3860  fUseDistTrackToCluster = kFALSE;
3861  fMaxDistTrackToClusterEta = 0;
3862  fMinDistTrackToClusterPhi = 0;
3863  fMaxDistTrackToClusterPhi = 0;
3864  break;
3865  case 1:
3866  if (!fUseDistTrackToCluster) fUseDistTrackToCluster=kTRUE;
3867  fMaxDistTrackToClusterEta = 0.008;//0.015;
3868  fMinDistTrackToClusterPhi = -0.03;//-0.01;
3869  fMaxDistTrackToClusterPhi = 0.03;//0.03;//0.04;
3870  break;
3871  case 2:
3872  if (!fUseDistTrackToCluster) fUseDistTrackToCluster=kTRUE;
3873  fMaxDistTrackToClusterEta = 0.012;//0.015;
3874  fMinDistTrackToClusterPhi = -0.05;//-0.01;
3875  fMaxDistTrackToClusterPhi = 0.04;//0.035;//0.05;
3876  break;
3877  case 3:
3878  if (!fUseDistTrackToCluster) fUseDistTrackToCluster=kTRUE;
3879  fMaxDistTrackToClusterEta = 0.016;//0.015;
3880  fMinDistTrackToClusterPhi = -0.09;//-0.015;
3881  fMaxDistTrackToClusterPhi = 0.06;//0.04;//0.1;
3882  break;
3883  case 4:
3884  if (!fUseDistTrackToCluster) fUseDistTrackToCluster=kTRUE;
3885  fMaxDistTrackToClusterEta = 0.018;//0.015;
3886  fMinDistTrackToClusterPhi = -0.11;//-0.015;
3887  fMaxDistTrackToClusterPhi = 0.07;//0.045;//0.13;
3888  break;
3889  case 5:
3890  if (!fUseDistTrackToCluster) fUseDistTrackToCluster=kTRUE;
3891  fMaxDistTrackToClusterEta = 0.02;//0.015;
3892  fMinDistTrackToClusterPhi = -0.13;//-0.02;
3893  fMaxDistTrackToClusterPhi = 0.08;//0.05;//0.15
3894  break;
3895 // case 6:
3896 // if (!fUseDistTrackToCluster) fUseDistTrackToCluster=kTRUE;
3897 // fMaxDistTrackToClusterEta = 0.022;//0.015;
3898 // fMinDistTrackToClusterPhi = -0.15;//-0.02;
3899 // fMaxDistTrackToClusterPhi = 0.10;//0.055;//0.2;
3900 // break;
3901  // pT dependent matching parameters
3902  case 6:
3903  if (!fUseDistTrackToCluster) fUseDistTrackToCluster=kTRUE;
3904  fUsePtDepTrackToCluster = kTRUE;
3905  fFuncPtDepEta = new TF1("func", "[1] + 1 / pow(x + pow(1 / ([0] - [1]), 1 / [2]), [2])");
3906  fFuncPtDepEta->SetParameters(0.03, 0.010, 2.5);
3907 
3908  fFuncPtDepPhi = new TF1("func", "[1] + 1 / pow(x + pow(1 / ([0] - [1]), 1 / [2]), [2])");
3909  fFuncPtDepPhi->SetParameters(0.08, 0.015, 2.);
3910  break;
3911  case 7:
3912  if (!fUseDistTrackToCluster) fUseDistTrackToCluster=kTRUE;
3913  fUsePtDepTrackToCluster = kTRUE;
3914  fFuncPtDepEta = new TF1("func", "[1] + 1 / pow(x + pow(1 / ([0] - [1]), 1 / [2]), [2])");
3915  fFuncPtDepEta->SetParameters(0.04, 0.010, 2.5);
3916 
3917  fFuncPtDepPhi = new TF1("func", "[1] + 1 / pow(x + pow(1 / ([0] - [1]), 1 / [2]), [2])");
3918  fFuncPtDepPhi->SetParameters(0.09, 0.015, 2.);
3919  break;
3920  case 8:
3921  if (!fUseDistTrackToCluster) fUseDistTrackToCluster=kTRUE;
3922  fUsePtDepTrackToCluster = kTRUE;
3923  fFuncPtDepEta = new TF1("func", "[1] + 1 / pow(x + pow(1 / ([0] - [1]), 1 / [2]), [2])");
3924  fFuncPtDepEta->SetParameters(0.05, 0.010, 2.5);
3925 
3926  fFuncPtDepPhi = new TF1("func", "[1] + 1 / pow(x + pow(1 / ([0] - [1]), 1 / [2]), [2])");
3927  fFuncPtDepPhi->SetParameters(0.10, 0.015, 1.75);
3928  break;
3929  case 9:
3930  if (!fUseDistTrackToCluster) fUseDistTrackToCluster=kTRUE;
3931  fUsePtDepTrackToCluster = kTRUE;
3932  fFuncPtDepEta = new TF1("func", "[1] + 1 / pow(x + pow(1 / ([0] - [1]), 1 / [2]), [2])");
3933  fFuncPtDepEta->SetParameters(0.06, 0.015, 2.5);
3934 
3935  fFuncPtDepPhi = new TF1("func", "[1] + 1 / pow(x + pow(1 / ([0] - [1]), 1 / [2]), [2])");
3936  fFuncPtDepPhi->SetParameters(0.12, 0.020, 1.75);
3937  break;
3938  case 10:
3939  if (!fUseDistTrackToCluster) fUseDistTrackToCluster=kTRUE;
3940  fUsePtDepTrackToCluster = kTRUE;
3941  fFuncPtDepEta = new TF1("func", "[1] + 1 / pow(x + pow(1 / ([0] - [1]), 1 / [2]), [2])");
3942  fFuncPtDepEta->SetParameters(0.035, 0.010, 2.5);
3943 
3944  fFuncPtDepPhi = new TF1("func", "[1] + 1 / pow(x + pow(1 / ([0] - [1]), 1 / [2]), [2])");
3945  fFuncPtDepPhi->SetParameters(0.085, 0.015, 2.0);
3946  break;
3947  case 11:
3948  if (!fUseDistTrackToCluster) fUseDistTrackToCluster=kTRUE;
3949  fUsePtDepTrackToCluster = kTRUE;
3950  fFuncPtDepEta = new TF1("func", "[1] + 1 / pow(x + pow(1 / ([0] - [1]), 1 / [2]), [2])");
3951  fFuncPtDepEta->SetParameters(0.045, 0.010, 2.5);
3952 
3953  fFuncPtDepPhi = new TF1("func", "[1] + 1 / pow(x + pow(1 / ([0] - [1]), 1 / [2]), [2])");
3954  fFuncPtDepPhi->SetParameters(0.095, 0.015, 1.75);
3955  break;
3956 
3957  default:
3958  AliError(Form("Track Matching Cut not defined %d",trackMatching));
3959  return kFALSE;
3960  }
3961  return kTRUE;
3962  // matching parameters for PHOS clusters
3963  }else if(fClusterType == 2) {
3964  switch(trackMatching){
3965  case 0:
3966  fUseDistTrackToCluster = kFALSE;
3967  fMaxDistTrackToClusterEta = 0;
3968  fMinDistTrackToClusterPhi = 0;
3969  fMaxDistTrackToClusterPhi = 0;
3970  break;
3971  case 1:
3972  if (!fUseDistTrackToCluster) fUseDistTrackToCluster=kTRUE;
3973  fMaxDistTrackToClusterEta = 0.005;//0.015;
3974  fMinDistTrackToClusterPhi = -0.03;//-0.025;
3975  fMaxDistTrackToClusterPhi = 0.03;//0.06;//0.3;
3976  break;
3977  case 2:
3978  if (!fUseDistTrackToCluster) fUseDistTrackToCluster=kTRUE;
3979  fMaxDistTrackToClusterEta = 0.01;//0.015;
3980  fMinDistTrackToClusterPhi = -0.09;//-0.025;
3981  fMaxDistTrackToClusterPhi = 0.07;//0.07;//0.4;
3982  break;
3983  case 3:
3984  if (!fUseDistTrackToCluster) fUseDistTrackToCluster=kTRUE;
3985  fMaxDistTrackToClusterEta = 0.015;//0.02;
3986  fMinDistTrackToClusterPhi = -0.15;//-0.03;
3987  fMaxDistTrackToClusterPhi = 0.11;//0.1;//0.5;
3988  break;
3989  case 4: //pT dependent for PCM-PHOS "default" selection
3990  if (!fUseDistTrackToCluster) fUseDistTrackToCluster=kTRUE;
3991  fUsePtDepTrackToCluster = kTRUE;
3992  fFuncPtDepEta = new TF1("func", "[1] + 1 / pow(x + pow(1 / ([0] - [1]), 1 / [2]), [2])");
3993  fFuncPtDepEta->SetParameters(0.05, 0.005, 3.0);
3994 
3995  fFuncPtDepPhi = new TF1("func", "[1] + 1 / pow(x + pow(1 / ([0] - [1]), 1 / [2]), [2])");
3996  fFuncPtDepPhi->SetParameters(0.33, 0.005, 2.3);
3997  break;
3998  case 5: //pT dependent for PCM-PHOS tight selection
3999  if (!fUseDistTrackToCluster) fUseDistTrackToCluster=kTRUE;
4000  fUsePtDepTrackToCluster = kTRUE;
4001  fFuncPtDepEta = new TF1("func", "[1] + 1 / pow(x + pow(1 / ([0] - [1]), 1 / [2]), [2])");
4002  fFuncPtDepEta->SetParameters(0.025, 0.002, 3.0);
4003 
4004  fFuncPtDepPhi = new TF1("func", "[1] + 1 / pow(x + pow(1 / ([0] - [1]), 1 / [2]), [2])");
4005  fFuncPtDepPhi->SetParameters(0.17, 0.005, 2.5);
4006  break;
4007  case 6: //pT dependent for PCM-PHOS loose selection
4008  if (!fUseDistTrackToCluster) fUseDistTrackToCluster=kTRUE;
4009  fUsePtDepTrackToCluster = kTRUE;
4010  fFuncPtDepEta = new TF1("func", "[1] + 1 / pow(x + pow(1 / ([0] - [1]), 1 / [2]), [2])");
4011  fFuncPtDepEta->SetParameters(0.07, 0.003, 2.5);
4012 
4013  fFuncPtDepPhi = new TF1("func", "[1] + 1 / pow(x + pow(1 / ([0] - [1]), 1 / [2]), [2])");
4014  fFuncPtDepPhi->SetParameters(0.45, 0.010, 2.0);
4015  break;
4016 
4017  default:
4018  AliError(Form("Track Matching Cut not defined %d",trackMatching));
4019  return kFALSE;
4020  }
4021  return kTRUE;
4022  }
4023  return kTRUE;
4024 }
4025 
4026 //___________________________________________________________________
4027 Bool_t AliCaloPhotonCuts::SetExoticClusterCut(Int_t exoticCell)
4028 {
4029  switch(exoticCell){
4030  case 0:
4031  fUseExoticCluster = 0;
4032  fExoticEnergyFracCluster = 0;
4033  break;
4034  case 1:
4035  if (!fUseExoticCluster)
4036  fUseExoticCluster = 1;
4037  fExoticEnergyFracCluster = 0.995;
4038  break;
4039  case 2:
4040  if (!fUseExoticCluster)
4041  fUseExoticCluster = 1;
4042  fExoticEnergyFracCluster = 0.99;
4043  break;
4044  case 3:
4045  if (!fUseExoticCluster)
4046  fUseExoticCluster = 1;
4047  fExoticEnergyFracCluster = 0.98;
4048  break;
4049  case 4:
4050  if (!fUseExoticCluster)
4051  fUseExoticCluster = 1;
4052  fExoticEnergyFracCluster = 0.975;
4053  break;
4054  case 5:
4055  if (!fUseExoticCluster)
4056  fUseExoticCluster = 1;
4057  fExoticEnergyFracCluster = 0.97;
4058  break;
4059  case 6:
4060  if (!fUseExoticCluster)
4061  fUseExoticCluster = 1;
4062  fExoticEnergyFracCluster = 0.965;
4063  break;
4064  case 7:
4065  if (!fUseExoticCluster)
4066  fUseExoticCluster = 1;
4067  fExoticEnergyFracCluster = 0.96;
4068  break;
4069  case 8:
4070  if (!fUseExoticCluster)
4071  fUseExoticCluster = 1;
4072  fExoticEnergyFracCluster = 0.95;
4073  break;
4074  case 9:
4075  if (!fUseExoticCluster)
4076  fUseExoticCluster = 1;
4077  fExoticEnergyFracCluster = 0.94;
4078  break;
4079  default:
4080  AliError(Form("Exotic cell Cut not defined %d",exoticCell));
4081  return kFALSE;
4082  }
4083  return kTRUE;
4084 }
4085 
4086 //___________________________________________________________________
4087 Bool_t AliCaloPhotonCuts::SetMinEnergyCut(Int_t minEnergy)
4088 {
4089  if (fIsPureCalo != 1){
4090  if (fClusterType!=2) {
4091  switch(minEnergy){
4092  case 0:
4093  if (!fUseMinEnergy) fUseMinEnergy=0;
4094  fMinEnergy=0.1;
4095  break;
4096  case 1:
4097  if (!fUseMinEnergy) fUseMinEnergy=1;
4098  fMinEnergy=0.5;
4099  break;
4100  case 2:
4101  if (!fUseMinEnergy) fUseMinEnergy=1;
4102  fMinEnergy=0.6;
4103  break;
4104  case 3:
4105  if (!fUseMinEnergy) fUseMinEnergy=1;
4106  fMinEnergy=0.7;
4107  break;
4108  case 4:
4109  if (!fUseMinEnergy) fUseMinEnergy=1;
4110  fMinEnergy=0.8;
4111  break;
4112  case 5:
4113  if (!fUseMinEnergy) fUseMinEnergy=1;
4114  fMinEnergy=0.9;
4115  break;
4116  case 6:
4117  if (!fUseMinEnergy) fUseMinEnergy=1;
4118  fMinEnergy=4.5;
4119  break;
4120  case 7:
4121  if (!fUseMinEnergy) fUseMinEnergy=1;
4122  fMinEnergy=5.0;
4123  break;
4124  case 8:
4125  if (!fUseMinEnergy) fUseMinEnergy=1;
4126  fMinEnergy=5.5;
4127  break;
4128  case 9:
4129  if (!fUseMinEnergy) fUseMinEnergy=1;
4130  fMinEnergy=6.0;
4131  break;
4132  case 10:
4133  if (!fUseMinEnergy) fUseMinEnergy=1;
4134  fMinEnergy=1.5;
4135  break;
4136  case 11:
4137  if (!fUseMinEnergy) fUseMinEnergy=1;
4138  fMinEnergy=1.0;
4139  break;
4140  default:
4141  AliError(Form("Minimum Energy Cut not defined %d",minEnergy));
4142  return kFALSE;
4143  }
4144  } else {
4145  switch(minEnergy){
4146  case 0:
4147  if (!fUseMinEnergy) fUseMinEnergy=0;
4148  fMinEnergy=0.1;
4149  break;
4150  case 1:
4151  if (!fUseMinEnergy) fUseMinEnergy=1;
4152  fMinEnergy=0.3;
4153  break;
4154  case 2:
4155  if (!fUseMinEnergy) fUseMinEnergy=1;
4156  fMinEnergy=0.5;
4157  break;
4158  case 3:
4159  if (!fUseMinEnergy) fUseMinEnergy=1;
4160  fMinEnergy=0.6;
4161  break;
4162  case 4:
4163  if (!fUseMinEnergy) fUseMinEnergy=1;
4164  fMinEnergy=0.7;
4165  break;
4166  case 5:
4167  if (!fUseMinEnergy) fUseMinEnergy=1;
4168  fMinEnergy=0.8;
4169  break;
4170  case 6:
4171  if (!fUseMinEnergy) fUseMinEnergy=1;
4172  fMinEnergy=0.9;
4173  break;
4174  case 7:
4175  if (!fUseMinEnergy) fUseMinEnergy=1;
4176  fMinEnergy=0.2;
4177  break;
4178  case 8:
4179  if (!fUseMinEnergy) fUseMinEnergy=1;
4180  fMinEnergy=0.4;
4181  break;
4182  default:
4183  AliError(Form("Minimum Energy Cut not defined %d",minEnergy));
4184  return kFALSE;
4185  }
4186  }
4187  return kTRUE;
4188  } else {
4189  switch(minEnergy){
4190  case 0:
4191  if (!fUseMinEnergy) fUseMinEnergy=0;
4192  fMinEnergy=0.1;
4193  break;
4194  case 1:
4195  if (!fUseMinEnergy) fUseMinEnergy=1;
4196  fMinEnergy=4.;
4197  break;
4198  case 2:
4199  if (!fUseMinEnergy) fUseMinEnergy=1;
4200  fMinEnergy=5.;
4201  break;
4202  case 3:
4203  if (!fUseMinEnergy) fUseMinEnergy=1;
4204  fMinEnergy=6.;
4205  break;
4206  case 4:
4207  if (!fUseMinEnergy) fUseMinEnergy=1;
4208  fMinEnergy=7.;
4209  break;
4210  case 5:
4211  if (!fUseMinEnergy) fUseMinEnergy=1;
4212  fMinEnergy=7.5;
4213  break;
4214  case 6:
4215  if (!fUseMinEnergy) fUseMinEnergy=1;
4216  fMinEnergy=8.;
4217  break;
4218  case 7:
4219  if (!fUseMinEnergy) fUseMinEnergy=1;
4220  fMinEnergy=8.5;
4221  break;
4222  case 8:
4223  if (!fUseMinEnergy) fUseMinEnergy=1;
4224  fMinEnergy=9.;
4225  break;
4226  case 9:
4227  if (!fUseMinEnergy) fUseMinEnergy=1;
4228  fMinEnergy=9.5;
4229  break;
4230  case 10:
4231  if (!fUseMinEnergy) fUseMinEnergy=1;
4232  fMinEnergy=1.5;
4233  break;
4234  default:
4235  AliError(Form("Minimum Energy Cut not defined %d",minEnergy));
4236  return kFALSE;
4237  }
4238  return kTRUE;
4239  }
4240  return kTRUE;
4241 }
4242 
4243 //___________________________________________________________________
4244 Bool_t AliCaloPhotonCuts::SetMinNCellsCut(Int_t minNCells)
4245 {
4246  switch(minNCells){
4247  case 0:
4248  if (!fUseNCells) fUseNCells=0;
4249  fMinNCells=0;
4250  break;
4251  case 1:
4252  if (!fUseNCells) fUseNCells=1;
4253  fMinNCells=1;
4254  break;
4255  case 2:
4256  if (!fUseNCells) fUseNCells=1;
4257  fMinNCells=2;
4258  break;
4259  case 3:
4260  if (!fUseNCells) fUseNCells=1;
4261  fMinNCells=3;
4262  break;
4263  case 4:
4264  if (!fUseNCells) fUseNCells=1;
4265  fMinNCells=4;
4266  break;
4267  case 5:
4268  if (!fUseNCells) fUseNCells=1;
4269  fMinNCells=5;
4270  break;
4271  case 6:
4272  if (!fUseNCells) fUseNCells=1;
4273  fMinNCells=6;
4274  break;
4275 
4276  default:
4277  AliError(Form("Min N cells Cut not defined %d",minNCells));
4278  return kFALSE;
4279  }
4280  return kTRUE;
4281 }
4282 
4283 //___________________________________________________________________
4284 Bool_t AliCaloPhotonCuts::SetMaxM02(Int_t maxM02)
4285 {
4286  fMaxM02CutNr = maxM02;
4287  if (fIsPureCalo == 1){
4288  fUseM02 = 2;
4289  return kTRUE;
4290  }
4291 
4292  switch(maxM02){
4293  case 0:
4294  if (!fUseM02) fUseM02=0;
4295  fMaxM02=100;
4296  break;
4297  case 1:
4298  if (!fUseM02) fUseM02=1;
4299  fMaxM02=1.;
4300  break;
4301  case 2:
4302  if (!fUseM02) fUseM02=1;
4303  fMaxM02=0.7;
4304  break;
4305  case 3:
4306  if (!fUseM02) fUseM02=1;
4307  fMaxM02=0.5;
4308  break;
4309  case 4:
4310  if (!fUseM02) fUseM02=1;
4311  fMaxM02=0.4;
4312  break;
4313  case 5:
4314  if (!fUseM02) fUseM02=1;
4315  fMaxM02=0.3;
4316  break;
4317  case 6:
4318  if (!fUseM02) fUseM02=1;
4319  fMaxM02=0.27;
4320  break;
4321  case 7: // PHOS cuts
4322  if (!fUseM02) fUseM02=1;
4323  fMaxM02=1.3;
4324  break;
4325  case 8: // PHOS cuts
4326  if (!fUseM02) fUseM02=1;
4327  fMaxM02=2.5;
4328  break;
4329  case 9:
4330  if (!fUseM02) fUseM02=1;
4331  fMaxM02=0.35;
4332  break;
4333  case 10: // a
4334  if (!fUseM02) fUseM02=1;
4335  fMaxM02=0.33;
4336  break;
4337  case 11: // b
4338  if (!fUseM02) fUseM02=1;
4339  fMaxM02=0.28;
4340  break;
4341  case 12: // c
4342  if (!fUseM02) fUseM02=1;
4343  fMaxM02=0.32;
4344  break;
4345 
4346  // E dependent M02 variations
4347  case 13: // d
4348  //(0.27 + 0.0072 * TMath::Power(clusEnergy,2));
4349  fUseM02=2;
4350  fMinM02CutNr=9;
4351  fMaxM02=0.4;
4352  break;
4353  case 14: // e
4354  //(0.31 + 0.0072 * TMath::Power(clusEnergy,2));
4355  fUseM02=2;
4356  fMinM02CutNr=9;
4357  fMaxM02=0.5;
4358  break;
4359  case 15: // f
4360  //(0.36 + 0.0072 * TMath::Power(clusEnergy,2));
4361  fUseM02=2;
4362  fMinM02CutNr=9;
4363  fMaxM02=0.7;
4364  break;
4365  case 16: // g
4366  //(0.37 + 0.0072 * TMath::Power(clusEnergy,2))
4367  fUseM02=2;
4368  fMinM02CutNr=9;
4369  fMaxM02=0.7;
4370  break;
4371  case 17: // h
4372  //(0.30 + 0.0072 * TMath::Power(clusEnergy,2));
4373  fUseM02=2;
4374  fMinM02CutNr=9;
4375  fMaxM02=0.5;
4376  break;
4377  case 18: // i
4378  // (0.35 + 0.0072 * TMath::Power(clusEnergy,2));
4379  fUseM02=2;
4380  fMinM02CutNr=9;
4381  fMaxM02=0.7;
4382  break;
4383  case 19: // j
4384  // (0.25 + 0.0072 * TMath::Power(clusEnergy,2));
4385  fUseM02=2;
4386  fMinM02CutNr=9;
4387  fMaxM02=0.39;
4388  break;
4389  case 20: //k
4390  //(0.27 + 0.0092 * TMath::Power(clusEnergy,2));
4391  fUseM02=2;
4392  fMinM02CutNr=9;
4393  fMaxM02=0.5;
4394  break;
4395  case 21: // l
4396  //(0.32 + 0.0072 * TMath::Power(clusEnergy,2));
4397  fUseM02=2;
4398  fMinM02CutNr=9;
4399  fMaxM02=0.5;
4400  break;
4401  case 22: // m
4402  // (0.32 + 0.0152 * TMath::Power(clusEnergy,2));
4403  fUseM02=2;
4404  fMinM02CutNr=9;
4405  fMaxM02=0.5;
4406  break;
4407  case 23: // n
4408  // (0.32 + 0.0238 * TMath::Power(clusEnergy,2));
4409  fUseM02=2;
4410  fMinM02CutNr=9;
4411  fMaxM02=0.7;
4412  break;
4413  case 24: // o
4414  // (0.27 + 0.0092 * TMath::Power(clusEnergy,2));
4415  fUseM02=2;
4416  fMinM02CutNr=9;
4417  fMaxM02=0.7;
4418  break;
4419  case 25: // p
4420  // (0.32 + 0.0072 * TMath::Power(clusEnergy,2));
4421  fUseM02=2;
4422  fMinM02CutNr=9;
4423  fMaxM02=0.7;
4424  break;
4425  case 26: // q
4426  // (0.34 + 0.0072 * TMath::Power(clusEnergy,2));
4427  fUseM02=2;
4428  fMinM02CutNr=9;
4429  fMaxM02=0.7;
4430  break;
4431  case 27: // r
4432  // (0.25 + 0.0072 * TMath::Power(clusEnergy,2));
4433  fUseM02=2;
4434  fMinM02CutNr=9;
4435  fMaxM02=0.5;
4436  break;
4437  case 28: // s
4438  // (0.32 + 0.0238 * TMath::Power(clusEnergy,2));
4439  fUseM02=2;
4440  fMinM02CutNr=9;
4441  fMaxM02=0.5;
4442  break;
4443  default:
4444  AliError(Form("Max M02 Cut not defined %d",maxM02));
4445  return kFALSE;
4446  }
4447  return kTRUE;
4448 }
4449 
4450 //___________________________________________________________________
4451 Float_t AliCaloPhotonCuts::CalculateMaxM02 (Int_t maxM02, Float_t clusEnergy){
4452  switch (maxM02){
4453  case 0:
4454  return 10;
4455  case 1:
4456  if (fMinNLM == 1 && fMaxNLM == 1 ){
4457  return FunctionM02(clusEnergy, 0.0662, -0.0201, -0.0955, 1.86e-3, 9.91 );
4458  } else if (fMinNLM == 2 && fMaxNLM == 2 ){
4459  return FunctionM02(clusEnergy, 0.353, -0.0264, -0.524, 5.59e-3, 21.9 );
4460  } else {
4461  return 10;
4462  }
4463  case 2:
4464  if (fMinNLM == 1 && fMaxNLM == 1 ){
4465  return FunctionM02(clusEnergy, 0.0662, -0.0201, -0.0, 1.86e-3, 9.91 );
4466  } else if (fMinNLM == 2 && fMaxNLM == 2 ){
4467  return FunctionM02(clusEnergy, 0.353, -0.0264, -0.424, 5.59e-3, 21.9 );
4468  } else {
4469  return 10;
4470  }
4471  case 3:
4472  if (fMinNLM == 1 && fMaxNLM == 1 ){
4473  return FunctionM02(clusEnergy, 0.0662, -0.0201, -0.2, 1.86e-3, 9.91 );
4474  } else if (fMinNLM == 2 && fMaxNLM == 2 ){
4475  return FunctionM02(clusEnergy, 0.353, -0.0264, -0.624, 5.59e-3, 21.9 );
4476  } else {
4477  return 10;
4478  }
4479 
4480  case 13: // d
4481  if( (0.27 + 0.0072 * TMath::Power(clusEnergy,2)) >= 0.4) return 0.4;
4482  else return (0.27 + 0.0072 * TMath::Power(clusEnergy,2));
4483  case 14: // e
4484  if( (0.31 + 0.0072 * TMath::Power(clusEnergy,2)) >= 0.5) return 0.5;
4485  else return (0.31 + 0.0072 * TMath::Power(clusEnergy,2));
4486  case 15: // f
4487  if( (0.36 + 0.0072 * TMath::Power(clusEnergy,2)) >= 0.7) return 0.7;
4488  else return (0.36 + 0.0072 * TMath::Power(clusEnergy,2));
4489  case 16: // g
4490  if( (0.37 + 0.0072 * TMath::Power(clusEnergy,2)) >= 0.7) return 0.7;
4491  else return (0.37 + 0.0072 * TMath::Power(clusEnergy,2));
4492  case 17: // h
4493  if( (0.30 + 0.0072 * TMath::Power(clusEnergy,2)) >= 0.5) return 0.5;
4494  else return (0.30 + 0.0072 * TMath::Power(clusEnergy,2));
4495  case 18: // i
4496  if( (0.35 + 0.0072 * TMath::Power(clusEnergy,2)) >= 0.7) return 0.7;
4497  else return (0.35 + 0.0072 * TMath::Power(clusEnergy,2));
4498  case 19: // j
4499  if( (0.25 + 0.0072 * TMath::Power(clusEnergy,2)) >= 0.39) return 0.39;
4500  else return (0.25 + 0.0072 * TMath::Power(clusEnergy,2));
4501  case 20: // k
4502  if( (0.27 + 0.0092 * TMath::Power(clusEnergy,2)) >= 0.5) return 0.5;
4503  else return (0.27 + 0.0092 * TMath::Power(clusEnergy,2));
4504  case 21: // l
4505  if( (0.32 + 0.0072 * TMath::Power(clusEnergy,2)) >= 0.5) return 0.5;
4506  else return (0.32 + 0.0072 * TMath::Power(clusEnergy,2));
4507  case 22: // m
4508  if( (0.32 + 0.0152 * TMath::Power(clusEnergy,2)) >= 0.5) return 0.5;
4509  else return (0.32 + 0.0152 * TMath::Power(clusEnergy,2));
4510  case 23: // n
4511  if( (0.32 + 0.0238 * TMath::Power(clusEnergy,2)) >= 0.7) return 0.7;
4512  else return (0.32 + 0.0238 * TMath::Power(clusEnergy,2));
4513  case 24: // o
4514  if( (0.27 + 0.0092 * TMath::Power(clusEnergy,2)) >= 0.7) return 0.7;
4515  else return (0.27 + 0.0092 * TMath::Power(clusEnergy,2));
4516  case 25: // p - standard
4517  if( (0.32 + 0.0072 * TMath::Power(clusEnergy,2)) >= 0.7) return 0.7;
4518  else return (0.32 + 0.0072 * TMath::Power(clusEnergy,2));
4519  case 26: // q
4520  if( (0.34 + 0.0072 * TMath::Power(clusEnergy,2)) >= 0.7) return 0.7;
4521  else return (0.34 + 0.0072 * TMath::Power(clusEnergy,2));
4522  case 27: // r
4523  if( (0.25 + 0.0072 * TMath::Power(clusEnergy,2)) >= 0.5) return 0.5;
4524  else return (0.25 + 0.0072 * TMath::Power(clusEnergy,2));
4525  case 28: // s
4526  if( (0.32 + 0.0238 * TMath::Power(clusEnergy,2)) >= 0.5) return 0.5;
4527  else return (0.32 + 0.0238 * TMath::Power(clusEnergy,2));
4528 
4529  default:
4530  AliError(Form("Max M02 for merged cluster Cut not defined %d",maxM02));
4531  return 10;
4532  }
4533  return 10;
4534 
4535 }
4536 
4537 //___________________________________________________________________
4538 Float_t AliCaloPhotonCuts::CalculateMinM02 (Int_t minM02, Float_t clusEnergy){
4539  switch (minM02){
4540  case 0:
4541  return 0.;
4542  case 1:
4543  if (FunctionM02(clusEnergy, 2.135, -0.245, 0., 0., 0. ) > 0.3)
4544  return FunctionM02(clusEnergy, 2.135, -0.245, 0., 0., 0. );
4545  else
4546  return 0.3;
4547  case 2:
4548  if (FunctionM02(clusEnergy, 2.135, -0.245, 0., 0., 0. ) > 0.27)
4549  return FunctionM02(clusEnergy, 2.135, -0.245, 0., 0., 0. );
4550  else
4551  return 0.27;
4552  case 3:
4553  if (FunctionM02(clusEnergy, 2.135, -0.245, 0., 0., 0. ) > 0.25)
4554  return FunctionM02(clusEnergy, 2.135, -0.245, 0., 0., 0. );
4555  else
4556  return 0.25;
4557  case 4:
4558  if (FunctionM02(clusEnergy, 2.135, -0.245, 0.1, 0., 0. ) > 0.27)
4559  return FunctionM02(clusEnergy, 2.135, -0.245, 0.1, 0., 0. );
4560  else
4561  return 0.27;
4562  case 5:
4563  if (FunctionM02(clusEnergy, 2.135, -0.245, -0.1, 0., 0. ) > 0.27)
4564  return FunctionM02(clusEnergy, 2.135, -0.245, -0.1, 0., 0. );
4565  else
4566  return 0.27;
4567  case 6:
4568  return 0.3;
4569  case 7:
4570  return 0.27;
4571  case 8:
4572  return 0.25;
4573  case 9:
4574  return 0.1;
4575  case 10:
4576  return 0.26;
4577  case 11:
4578  return 0.28;
4579  case 12:
4580  return 0.29;
4581 
4582  default:
4583  AliError(Form("Min M02 for merged cluster Cut not defined %d",minM02));
4584  return -1;
4585  }
4586  return -1;
4587 }
4588 
4589 
4590 //___________________________________________________________________
4591 Bool_t AliCaloPhotonCuts::SetMinM02(Int_t minM02)
4592 {
4593  fMinM02CutNr = minM02;
4594  if (fIsPureCalo == 1){
4595  fUseM02 = 2;
4596  return kTRUE;
4597  }
4598 
4599  switch(minM02){
4600  case 0:
4601<