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