]> git.uio.no Git - u/mrichter/AliRoot.git/blob - PWGGA/GammaConv/AliAnalysisTaskGammaConvCalo.h
- added new tasks for combination of EMCAL and PCM for pi0 and photons
[u/mrichter/AliRoot.git] / PWGGA / GammaConv / AliAnalysisTaskGammaConvCalo.h
1 #ifndef ALIANLYSISTASKGAMMACONVCALO_cxx
2 #define ALIANLYSISTASKGAMMACONVCALO_cxx
3
4 #include "AliAnalysisTaskSE.h"
5 #include "AliESDtrack.h"
6 #include "AliV0ReaderV1.h"
7 #include "AliKFConversionPhoton.h"
8 #include "AliGammaConversionAODBGHandler.h"
9 #include "AliConversionAODBGHandlerRP.h"
10 #include "AliCaloPhotonCuts.h"
11 #include "AliConversionMesonCuts.h"
12 #include "AliAnalysisManager.h"
13 #include "TProfile2D.h"
14 #include "TH3.h"
15 #include "TH3F.h"
16
17 class AliAnalysisTaskGammaConvCalo : public AliAnalysisTaskSE {
18         public:
19
20                 AliAnalysisTaskGammaConvCalo();
21                 AliAnalysisTaskGammaConvCalo(const char *name);
22                 virtual ~AliAnalysisTaskGammaConvCalo();
23
24                 virtual void   UserCreateOutputObjects();
25                 virtual Bool_t Notify();
26                 virtual void   UserExec(Option_t *);
27                 virtual void   Terminate(const Option_t*);
28                 void InitBack();
29
30                 void SetIsHeavyIon(Int_t flag){
31                         fIsHeavyIon = flag;    
32                 }
33
34                 // base functions for selecting photon and meson candidates in reconstructed data
35                 void ProcessClusters();
36                 void ProcessPhotonCandidates();
37                 void PhotonTagging();
38                 void CalculatePi0Candidates();
39                 
40                 // MC functions
41                 void SetIsMC(Bool_t isMC){fIsMC=isMC;}
42                 void ProcessMCParticles();
43                 void ProcessAODMCParticles();
44                 void RelabelAODPhotonCandidates(Bool_t mode);
45                 void ProcessTruePhotonCandidates( AliAODConversionPhoton* TruePhotonCandidate);
46                 void ProcessTrueClusterCandidates( AliAODConversionPhoton* TruePhotonCandidate);
47                 void ProcessTruePhotonCandidatesAOD( AliAODConversionPhoton* TruePhotonCandidate);
48                 void ProcessTrueMesonCandidates( AliAODConversionMother *Pi0Candidate, AliAODConversionPhoton *TrueGammaCandidate0, AliAODConversionPhoton *TrueGammaCandidate1);
49                 void ProcessTrueMesonCandidatesAOD(AliAODConversionMother *Pi0Candidate, AliAODConversionPhoton *TrueGammaCandidate0, AliAODConversionPhoton *TrueGammaCandidate1);
50                 
51                 // switches for additional analysis streams or outputs
52                 void SetDoMesonAnalysis(Bool_t flag){fDoMesonAnalysis = flag;}
53                 void SetDoMesonQA(Int_t flag){fDoMesonQA = flag;}
54                 void SetDoPhotonQA(Int_t flag){fDoPhotonQA = flag;}
55                 
56             // Setting the cut lists for the conversion photons
57                 void SetConversionCutList(Int_t nCuts, TList *CutArray){
58                         fnCuts = nCuts;
59                         fCutArray = CutArray;
60                 }
61
62             // Setting the cut lists for the calo photons
63                 void SetCaloCutList(Int_t nCuts, TList *CutArray){
64                         fnCuts = nCuts;
65                         fClusterCutArray = CutArray;
66                 }
67                 
68                 // Setting the cut lists for the meson
69                 void SetMesonCutList(Int_t nCuts, TList *CutArray){
70                         fnCuts = nCuts;
71                         fMesonCutArray = CutArray;
72                 }
73
74                 // emcal functions
75                 Double_t GetMaxCellEnergy(const AliVCluster *c) const { Short_t id=-1; return GetMaxCellEnergy(c,id); }
76                 Double_t GetMaxCellEnergy(const AliVCluster *c, Short_t &id) const;
77
78                 // BG HandlerSettings
79                 void CalculateBackground();
80                 void CalculateBackgroundRP();
81                 void RotateParticle(AliAODConversionPhoton *gamma);
82                 void RotateParticleAccordingToEP(AliAODConversionPhoton *gamma, Double_t previousEventEP, Double_t thisEventEP);
83                 void SetMoveParticleAccordingToVertex(Bool_t flag){fMoveParticleAccordingToVertex = flag;}
84                 void FillPhotonCombinatorialBackgroundHist(AliAODConversionPhoton *TruePhotonCandidate, Int_t pdgCode[]);
85                 void MoveParticleAccordingToVertex(AliAODConversionPhoton* particle,const AliGammaConversionAODBGHandler::GammaConversionVertex *vertex);
86                 void UpdateEventByEventData();
87                 
88                 // Additional functions for convenience
89                 void SetLogBinningXTH2(TH2* histoRebin);
90                 Int_t GetSourceClassification(Int_t daughter, Int_t pdgCode);
91         
92         protected:
93                 AliV0ReaderV1                                           *fV0Reader;                                                     // basic photon Selection Task
94                 AliGammaConversionAODBGHandler          **fBGHandler;                                           // BG handler for Conversion 
95                 AliConversionAODBGHandlerRP             **fBGHandlerRP;                                         // BG handler for Conversion (possibility to mix with respect to RP)
96                 AliGammaConversionAODBGHandler          **fBGClusHandler;                                       // BG handler for Cluster
97                 AliConversionAODBGHandlerRP             **fBGClusHandlerRP;                             // BG handler for Cluster (possibility to mix with respect to RP)
98                 AliVEvent                                                       *fInputEvent;                                           // current event
99                 AliMCEvent                                                      *fMCEvent;                                                      // corresponding MC event
100                 AliStack                                                        *fMCStack;                                                      // stack belonging to MC event
101 //              AliESDEvent                                                     *fEsdEv;                                                        //!pointer to input esd event
102 //              AliAODEvent                                                     *fAodEv;                                                        //!pointer to input aod event
103 //              TObjArray                                                       *fEsdClusters;                                          //!pointer to esd clusters
104 //              AliESDCaloCells                                         *fEsdCells;                                                     //!pointer to esd cells
105 //              TObjArray                                                       *fAodClusters;                                          //!pointer to aod clusters
106 //              AliAODCaloCells                                         *fAodCells;                                                     //!pointer to aod cells
107                 TList                                                           **fCutFolder;                                           // Array of lists for containers belonging to cut
108                 TList                                                           **fESDList;                                                     // Array of lists with histograms with reconstructed properties   
109                 TList                                                           **fBackList;                                            // Array of lists with BG THnSparseF
110                 TList                                                           **fMotherList;                                          // Array of lists with Signal THnSparseF
111                 TList                                                           **fPhotonDCAList;                                       // Array of lists with photon dca trees
112                 TList                                                           **fMesonDCAList;                                        // Array of lists with meson dca trees
113                 TList                                                           **fTrueList;                                            // Array of lists with histograms with MC validated reconstructed properties
114                 TList                                                           **fMCList;                                                      // Array of lists with histograms with pure MC information
115                 TList                                                           **fHeaderNameList;                                      // Array of lists with header names for MC header selection
116                 TList                                                           **fTagOutputList;                                       //!Array of lists of output histograms for tagged photons
117                 TList                                                           *fOutputContainer;                                      // Output container
118                 TClonesArray                                            *fReaderGammas;                                         // Array with conversion photons selected by V0Reader Cut
119                 TList                                                           *fGammaCandidates;                                      // current list of photon candidates
120                 TList                                                           *fClusterCandidates;                            //! current list of cluster candidates
121                 TList                                                           *fCutArray;                                                     // List with Conversion Cuts
122                 AliConversionCuts                                       *fConversionCuts;                                       // ConversionCutObject
123                 TList                                                           *fClusterCutArray;                                      // List with Cluster Cuts
124                 AliCaloPhotonCuts                                       *fCaloPhotonCuts;                                       // CaloPhotonCutObject
125                 TList                                                           *fMesonCutArray;                                        // List with Meson Cuts
126                 AliConversionMesonCuts                          *fMesonCuts;                                            // MesonCutObject
127                 
128                 //histograms for Conversions reconstructed quantities
129                 TH1F                                                            **fHistoConvGammaPt;                            //! histogram conversion photon pT
130                 TH1F                                                            **fHistoConvGammaR;                                     //! histogram conversion photon R
131                 TH1F                                                            **fHistoConvGammaEta;                           //! histogram conversion photon Eta
132                 TTree                                                           **fTreeConvGammaPtDcazCat;                      //! tree with dca for conversions
133                 Float_t                                                         fPtGamma;                                                       //! pt of conversion for tree
134                 Float_t                                                         fDCAzPhoton;                                            //! dcaz of conversion for tree
135                 Float_t                                                         fRConvPhoton;                                           //! R of conversion for tree
136                 Float_t                                                         fEtaPhoton;                                                     //! eta of conversion for tree
137                 UChar_t                                                         fCharCatPhoton;                                         //! category of conversion for tree
138                 UChar_t                                                         fCharPhotonMCInfo;                                      //! MC info of conversion for tree
139                                                                                         // 0: garbage,
140                                                                                         // 1: background
141                                                                                         // 2: secondary photon not from eta or k0s,
142                                                                                         // 3: secondary photon from eta, 
143                                                                                         // 4: secondary photon from k0s, 
144                                                                                         // 5: dalitz
145                                                                                         // 6: primary gamma
146                 //histograms for mesons reconstructed quantities
147                 TH2F                                                            **fHistoMotherInvMassPt;                        //! array of histogram with signal + BG for same event photon pairs, inv Mass, pt
148                 THnSparseF                                                      **fSparseMotherInvMassPtZM;                     //! array of THnSparseF with signal + BG for same event photon pairs, inv Mass, pt
149                 TH2F                                                            **fHistoMotherBackInvMassPt;            //! array of histogram with BG for mixed event photon pairs, inv Mass, pt
150                 THnSparseF                                                      **fSparseMotherBackInvMassPtZM;         //! array of THnSparseF with BG for same event photon pairs, inv Mass, pt
151                 TH2F                                                            **fHistoMotherInvMassEalpha;            //! array of histograms with alpha cut of 0.1 for inv mass vs pt
152                 TH2F                                                            **fHistoMotherPi0PtY;                           //! array of histograms with invariant mass cut of 0.05 && pi0cand->M() < 0.17, pt, Y
153                 TH2F                                                            **fHistoMotherEtaPtY;                           //! array of histograms with invariant mass cut of 0.45 && pi0cand->M() < 0.65, pt, Y
154                 TH2F                                                            **fHistoMotherPi0PtAlpha;                       //! array of histograms with invariant mass cut of 0.05 && pi0cand->M() < 0.17, pt, alpha
155                 TH2F                                                            **fHistoMotherEtaPtAlpha;                       //! array of histograms with invariant mass cut of 0.45 && pi0cand->M() < 0.65, pt, alpha
156                 TH2F                                                            **fHistoMotherPi0PtOpenAngle;           //! array of histograms with invariant mass cut of 0.05 && pi0cand->M() < 0.17, pt, openAngle
157                 TH2F                                                            **fHistoMotherEtaPtOpenAngle;           //! array of histograms with invariant mass cut of 0.45 && pi0cand->M() < 0.65, pt, openAngle
158                 TTree                                                           **fTreeMesonsInvMassPtDcazMinDcazMaxFlag; //! array of trees with dca information for mesons
159                 Float_t                                                         fInvMass;                                                       // inv mass for meson tree
160                 Float_t                                                         fPt;                                                            // pt for meson tree 
161                 Float_t                                                         fDCAzGammaMin;                                          // dcaz for meson tree gamma 1
162                 Float_t                                                         fDCAzGammaMax;                                          // dcaz for meson tree gamma 2
163                 UChar_t                                                         fCharFlag;                                                      // category of meson for tree
164                 UChar_t                                                         fCharMesonMCInfo;                                       // MC information meson for tree
165                                                                                                 // 0: garbage,
166                                                                                                 // 1: background
167                                                                                                 // 2: secondary meson not from eta or k0s,
168                                                                                                 // 3: secondary meson from eta, 
169                                                                                                 // 4: secondary meson from k0s, 
170                                                                                                 // 5: dalitz
171                                                                                                 // 6: primary meson gamma-gamma-channel
172
173                 // histograms for rec photons tagged by Calo
174                 TH1F                                                            **fHistoConvGammaUntagged;                      //! array of histo for untagged photon candidates vs pt
175                 TH1F                                                            **fHistoConvGammaTagged;                        //! array of histo for tagged photon candidates vs pt
176                 TH1F                                                            **fHistoConvGammaPi0Tagged;                     //! array of histo for tagged photon candidates vs pt
177                 TH1F                                                            **fHistoConvGammaEtaTagged;                     //! array of histo for tagged photon candidates vs pt
178                 TH2F                                                            **fHistoPhotonPairAll;                          //! array of histo for pairs
179                 TH2F                                                            **fHistoPhotonPairAllGam;                       //! array of histo for pairs vs. pt of converted photon
180                 // histograms for rec photon clusters
181                 TH1F                                                            ** fHistoClusGammaPt;                           //! array of histos with cluster, pt
182                                                                                 
183                 //histograms for pure MC quantities
184                 TH1I                                                            **fHistoMCHeaders;                                      //! array of histos for header names
185                 TH1F                                                            **fHistoMCAllGammaPt;                           //! array of histos with all gamma, pT
186                 TH1F                                                            **fHistoMCDecayGammaPi0Pt;                      //! array of histos with decay gamma from pi0, pT
187                 TH1F                                                            **fHistoMCDecayGammaRhoPt;                      //! array of histos with decay gamma from rho, pT
188                 TH1F                                                            **fHistoMCDecayGammaEtaPt;                      //! array of histos with decay gamma from eta, pT
189                 TH1F                                                            **fHistoMCDecayGammaOmegaPt;            //! array of histos with decay gamma from omega, pT
190                 TH1F                                                            **fHistoMCDecayGammaEtapPt;                     //! array of histos with decay gamma from eta', pT
191                 TH1F                                                            **fHistoMCDecayGammaPhiPt;                      //! array of histos with decay gamma from phi, pT
192                 TH1F                                                            **fHistoMCDecayGammaSigmaPt;            //! array of histos with decay gamma from Sigma0, pT
193                 TH1F                                                            **fHistoMCConvGammaPt;                          //! array of histos with converted gamma, pT
194                 TH1F                                                            **fHistoMCConvGammaR;                           //! array of histos with converted gamma, R
195                 TH1F                                                            **fHistoMCConvGammaEta;                         //! array of histos with converted gamma, Eta
196                 TH1F                                                            **fHistoMCPi0Pt;                                        //! array of histos with weighted pi0, pT
197                 TH1F                                                            **fHistoMCPi0WOWeightPt;                        //! array of histos with unweighted pi0, pT
198                 TH1F                                                            **fHistoMCEtaPt;                                        //! array of histos with weighted eta, pT
199                 TH1F                                                            **fHistoMCEtaWOWeightPt;                        //! array of histos with unweighted eta, pT
200                 TH1F                                                            **fHistoMCPi0InAccPt;                           //! array of histos with weighted pi0 in acceptance, pT
201                 TH1F                                                            **fHistoMCEtaInAccPt;                           //! array of histos with weighted eta in acceptance, pT
202                 TH2F                                                            **fHistoMCPi0PtY;                                       //! array of histos with weighted pi0, pT, Y
203                 TH2F                                                            **fHistoMCEtaPtY;                                       //! array of histos with weighted eta, pT, Y
204                 TH1F                                                            **fHistoMCK0sPt;                                        //! array of histos with weighted K0s, pT
205                 TH1F                                                            **fHistoMCK0sWOWeightPt;                        //! array of histos with unweighted K0s, pT
206                 TH2F                                                            **fHistoMCK0sPtY;                                       //! array of histos with weighted K0s, pT, Y
207                 TH2F                                                            **fHistoMCSecPi0PtvsSource;                     //! array of histos with secondary pi0, pT, source
208                 TH1F                                                            **fHistoMCSecPi0Source;                         //! array of histos with secondary pi0, source
209                 TH1F                                                            **fHistoMCSecEtaPt;                                     //! array of histos with secondary eta, pT
210                 TH1F                                                            **fHistoMCSecEtaSource;                         //! array of histos with secondary eta, source
211                 // MC validated reconstructed quantities mesons
212                 TH2F                                                            **fHistoTrueMotherInvMassPt;                                    //! array of histos with validated mothers, invMass, pt
213                 TH2F                                                            **fHistoTruePrimaryMotherInvMassPt;                             //! array of histos with validated weighted primary mothers, invMass, pt
214                 TH2F                                                            **fHistoTruePrimaryMotherW0WeightingInvMassPt;  //! array of histos with validated unweighted primary mothers, invMass, pt
215                 TProfile2D                                                      **fProfileTruePrimaryMotherWeightsInvMassPt;    //! array of profiles with weights for validated primary mothers, invMass, pt   
216                 TH2F                                                            **fHistoTruePrimaryPi0MCPtResolPt;                              //! array of histos with validated weighted primary pi0, MCpt, resol pt
217                 TH2F                                                            **fHistoTruePrimaryEtaMCPtResolPt;                              //! array of histos with validated weighted primary eta, MCpt, resol pt
218                 TH2F                                                            **fHistoTrueSecondaryMotherInvMassPt;                   //! array of histos with validated secondary mothers, invMass, pt
219                 TH2F                                                            **fHistoTrueSecondaryMotherFromK0sInvMassPt;    //! array of histos with validated secondary mothers from K0s, invMass, pt
220                 TH1F                                                            **fHistoTrueK0sWithPi0DaughterMCPt;                             //! array of histos with K0s with reconstructed pi0 as daughter, pt
221                 TH2F                                                            **fHistoTrueSecondaryMotherFromEtaInvMassPt;    //! array of histos with validated secondary mothers from eta, invMass, pt
222                 TH1F                                                            **fHistoTrueEtaWithPi0DaughterMCPt;                             //! array of histos with eta with reconstructed pi0 as daughter, pt
223                 TH2F                                                            **fHistoTrueSecondaryMotherFromLambdaInvMassPt; //! array of histos with validated secondary mothers from Lambda, invMass, pt
224                 TH1F                                                            **fHistoTrueLambdaWithPi0DaughterMCPt;                  //! array of histos with lambda with reconstructed pi0 as daughter, pt
225                 TH2F                                                            **fHistoTrueBckGGInvMassPt;                                             //! array of histos with pure gamma gamma combinatorial BG, invMass, pt
226                 TH2F                                                            **fHistoTrueBckContInvMassPt;                                   //! array of histos with        contamination BG, invMass, pt
227                 TH2F                                                            **fHistoTruePi0PtY;                                                             //! array of histos with        validated pi0, pt, Y
228                 TH2F                                                            **fHistoTrueEtaPtY;                                                             //! array of histos with validated eta, pt, Y
229                 TH2F                                                            **fHistoTruePi0PtAlpha;                                                 //! array of histos with validated pi0, pt, alpha
230                 TH2F                                                            **fHistoTrueEtaPtAlpha;                                                 //! array of histos with validated eta, pt, alpha
231                 TH2F                                                            **fHistoTruePi0PtOpenAngle;                                             //! array of histos with validated pi0, pt, openAngle
232                 TH2F                                                            **fHistoTrueEtaPtOpenAngle;                                             //! array of histos with validated eta, pt, openAngle
233                 TH2F                                                            **fHistoTrueMotherDalitzInvMassPt;                              //! array of histos with validated mother, but Dalitz decay, invMass, pt
234                 // MC validated reconstructed quantities photons
235                 TH1F                                                            **fHistoTrueConvGammaPt;                                                //! array of histos with validated conversion photon, pt
236                 TH1F                                                            **fHistoTrueConvPi0GammaPt;                                             //! array of histos with validated conversion photon from pi0, pt
237                 TH1F                                                            **fHistoTrueConvGammaEta;                                               //! array of histos with validated conversion photon, eta
238                 TH2F                                                            **fHistoCombinatorialPt;                                                //! array of histos with combinatorial BG, pt, source
239                 TH1F                                                            **fHistoTruePrimaryConvGammaPt;                                 //! array of histos with validated primary conversion photon, pt  
240                 TH2F                                                            **fHistoTruePrimaryConvGammaESDPtMCPt;                  //! array of histos with validated primary conversion photon, rec pt, mc pt  
241                 TH1F                                                            **fHistoTrueSecondaryConvGammaPt;                               //! array of histos with validated secondary conversion photon, pt  
242                 TH1F                                                            **fHistoTrueSecondaryConvGammaFromXFromK0sPt;   //! array of histos with validated secondary conversion photon from K0s, pt  
243                 TH1F                                                            **fHistoTrueSecondaryConvGammaFromXFromLambdaPt;//! array of histos with validated secondary conversion photon from Lambda, pt  
244                 TH2F                                                            **fHistoTrueDalitzPsiPairDeltaPhi;                              //! array of histos with validated dalitz virtual photon, delta phi, psi pair  
245                 TH2F                                                            **fHistoTrueGammaPsiPairDeltaPhi;                               //! array of histos with validated conversion photon, delta phi, psi pair
246                 TH1F                                                            ** fHistoTrueClusGammaPt;                                               //! array of histos with validated cluster, pt
247                 TH1F                                                            ** fHistoTruePrimaryClusGammaPt;                                //! array of histos with validated primary cluster, pt
248                 TH2F                                                            ** fHistoTruePrimaryClusGammaESDPtMCPt;                 //! array of histos with validated primary cluster, rec Pt, MC pt
249
250                 // event histograms
251                 TH1I                                                            **fHistoNEvents;                                                                //! array of histos with event information
252                 TH1I                                                            **fHistoNGoodESDTracks;                                                 //! array of histos with number of good tracks (2010 Standard track cuts)
253                 TH1I                                                            **fHistoNGammaCandidates;                                               //! array of histos with number of gamma candidates per event
254                 TH2F                                                            **fHistoNGoodESDTracksVsNGammaCanditates;               //! array of histos with number of good tracks vs gamma candidates
255                 TH1I                                                            **fHistoNV0Tracks;                                                              //! array of histos with V0 counts
256                 TProfile                                                        **fProfileEtaShift;                                                             //! array of profiles with eta shift
257                 
258                 // additional variables
259                 Double_t                                                        fEventPlaneAngle;                                       // EventPlaneAngle
260                 TRandom3                                                        fRandom;                                                        // random 
261                 Int_t                                                           fNGammaCandidates;                                      // number of gamma candidates in event
262                 Double_t                                                        *fUnsmearedPx;                                          //[fNGammaCandidates]
263                 Double_t                                                        *fUnsmearedPy;                                          //[fNGammaCandidates]
264                 Double_t                                                        *fUnsmearedPz;                                          //[fNGammaCandidates]
265                 Double_t                                                        *fUnsmearedE;                                           //[fNGammaCandidates]
266                 Int_t                                                           *fMCStackPos;                                           //[fNGammaCandidates]
267                 Int_t                                                           *fMCStackNeg;                                           //[fNGammaCandidates]
268                 Int_t                                                           *fESDArrayPos;                                          //[fNGammaCandidates]
269                 Int_t                                                           *fESDArrayNeg;                                          //[fNGammaCandidates]
270                 Int_t                                                           fnCuts;                                                         // number of cuts to be analysed in parallel
271                 Int_t                                                           fiCut;                                                          // current cut  
272                 Bool_t                                                          fMoveParticleAccordingToVertex;         // boolean for BG calculation
273                 Int_t                                                           fIsHeavyIon;                                            // switch for pp = 0, PbPb = 1, pPb = 2
274                 Bool_t                                                          fDoMesonAnalysis;                                       // flag for meson analysis
275                 Int_t                                                           fDoMesonQA;                                                     // flag for meson QA
276                 Int_t                                                           fDoPhotonQA;                                            // flag for photon QA
277                 Bool_t                                                          fIsFromMBHeader;                                        // flag for MC headers
278                 Bool_t                                                          fIsMC;                                                          // flag for MC information
279
280
281                 // cluster cut variables
282                 Double_t                                                        fMinE;
283                 Int_t                                                           fNminCells;
284                 Double_t                                                        fEMCm02cut;
285                 //double                                                        fMinErat = 0;
286                 //double                                                        fMinEcc = 0;
287
288                 
289                 //  TString                                                      fClusName;                                                     // cluster branch name (def="")
290                 //  const TObjArray                                     *fRecPoints;                                            // pointer to rec points (AliAnalysisTaskEMCALClusterizeFast)
291                 //  const TClonesArray                          *fDigits;                                                       // pointer to digits     (AliAnalysisTaskEMCALClusterizeFast)
292
293         private:
294                 AliAnalysisTaskGammaConvCalo(const AliAnalysisTaskGammaConvCalo&); // Prevent copy-construction
295                 AliAnalysisTaskGammaConvCalo &operator=(const AliAnalysisTaskGammaConvCalo&); // Prevent assignment
296
297                 ClassDef(AliAnalysisTaskGammaConvCalo, 1);
298 };
299
300 #endif