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