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