checked in Astrids Task to EMCALTasks
[u/mrichter/AliRoot.git] / PWGGA / EMCALTasks / AliEMCalpi0ClusterEvaluationTask.cxx
1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 *             /home/amorreal/Alice/Work/taskEMCal/AliEMCalpi0ClusterEvaluationTask.cxx
4 *
5 * Author: The ALICE Off-line Project.                                    *
6 * Contributors are mentioned in the code where appropriate.              *
7 *                                                                        *
8 * Permission to use, copy, modify and distribute this software and its   *
9 * documentation strictly for non-commercial purposes is hereby granted   *
10 * without fee, provided that the above copyright notice appears in all   *
11 * copies and that both the copyright notice and this permission notice   *
12 * appear in the supporting documentation. The authors make no claims     *
13 * about the suitability of this software for any purpose. It is          *
14 * provided"as is" without express or implied warranty.                  *
15 **************************************************************************/
16 //
17 //-----------------------------------------------------------------------------
18 /// \class AliEMCalpi0ClusterEvaluationTask
19 /// AliAnalysisTaskSE to analyse data from ESDs (pi0 and Eta in the EMCal).
20 /// The Task reads as input ESDs then it produces a .root file:
21 /// general Histos and Ntuple:Tree
22 ///
23 /// \author Astrid Morreale
24 //-----------------------------------------------------------------------------
25 #include "AliEMCalpi0ClusterEvaluationTask.h"
26
27 #include <Riostream.h>
28 #include <TChain.h>
29 #include <TTree.h>
30 #include <TH1F.h>
31 #include <TH2F.h>
32 #include <TH3F.h>
33 #include <TH1D.h>
34 #include <TH2D.h>
35 #include <TH3D.h>
36 #include <TCanvas.h>
37 #include <TList.h>
38 #include <TFile.h>
39 #include <TLorentzVector.h>
40 #include <TNtuple.h>
41 #include <TRandom3.h>
42
43
44 #include "AliAnalysisTaskSE.h"
45 #include "AliAnalysisManager.h"
46 #include "AliCDBManager.h"
47 #include "AliCDBStorage.h"
48 #include "AliStack.h"
49 #include "AliESDtrackCuts.h"
50 #include "AliESDEvent.h"
51 #include "AliESDInputHandler.h"
52 #include "AliAODHandler.h"
53 #include "AliAODEvent.h"
54 #include "AliMCEvent.h"
55 #include "AliEMCALGeometry.h"
56 #include "AliInputEventHandler.h"
57 #include "AliESDInputHandler.h"
58 #include "AliAODInputHandler.h"
59 #include "AliCentrality.h"
60 #include "AliEMCALRecoUtils.h"
61 #include "AliExternalTrackParam.h"
62
63
64 // ROOT includes
65 #include <TGeoManager.h>
66 #include <TGeoMatrix.h>
67 #include <TGeoBBox.h>
68 #include <TH2F.h>
69 #include <TArrayI.h>
70 #include <TArrayF.h>
71 #include <TObjArray.h>
72
73 // STEER includes
74 #include "AliVCluster.h"
75 #include "AliVCaloCells.h"
76 #include "AliLog.h"
77 #include "AliPID.h"
78 #include "AliESDEvent.h"
79 #include "AliAODEvent.h"
80 #include "AliESDtrack.h"
81 #include "AliAODTrack.h"
82 #include "AliExternalTrackParam.h"
83 #include "AliESDfriendTrack.h"
84 #include "AliTrackerBase.h"
85
86
87 // EMCAL includes
88 #include "AliEMCALRecoUtils.h"
89 #include "AliEMCALGeometry.h"
90 #include "AliTrackerBase.h"
91 #include "AliEMCALCalibTimeDepCorrection.h" // Run dependent
92 #include "AliEMCALPIDUtils.h"
93
94 #include <fstream>
95 #include <cassert>
96
97 #include <TChain.h>
98 #include <TError.h>
99 #include <TGeoGlobalMagField.h>
100 #include <TGeoManager.h>
101 #include <TLorentzVector.h>
102
103 #include <Riostream.h>
104
105 ///\cond CLASSIMP
106 ClassImp(AliEMCalpi0ClusterEvaluationTask)
107 ///\endcond
108
109 //________________________________________________________________________
110 AliEMCalpi0ClusterEvaluationTask::AliEMCalpi0ClusterEvaluationTask(const char *name) :
111 AliAnalysisTaskSE(name),
112
113 fEvent(0),      kAllMB(kFALSE),     isMB(0),     isAnyINT(0),
114 isCentral(0),   isSemiCentral(0),   isEga(0),    ega0(0),
115 ega1(0),        ega2(0),            ega3(0),     ega4(0),
116 ega5(0),        ega6(0),            ega7(0),     ega8(0),
117 ega9(0),        mb0(0),             mb1(0),      mb2(0),
118 mb3(0),         mb4(0),             mb5(0),      mb6(0),
119 mb7(0),         mb8(0),             mb9(0),      allmb0(0),
120 allmb1(0),      allmb2(0),          allmb3(0),   allmb4(0),
121 allmb5(0),      allmb6(0),          allmb7(0),   allmb8(0),
122 allmb9(0),      cent0(0),           cent1(0),    cent2(0),
123 cent3(0),       cent4(0),           cent5(0),    cent6(0),
124 cent7(0),       cent8(0),           cent9(0),    semicent0(0),
125 semicent1(0),   semicent2(0),       semicent3(0),semicent4(0),
126 semicent5(0),   semicent6(0),       semicent7(0),semicent8(0),
127 semicent9(0),   kAllMBmx(0),        isMBmx(0),   isAnyINTmx(0),
128 isCentralmx(0), isSemiCentralmx(0), isEgamx(0),  all(0),
129 allmb(0),       mb(0),              central(0),  semicentral(0),
130 ega(0),         crossEnergy(0),
131
132
133 fHistList(0)
134 {
135
136     // Define input and output slots here
137     // Input slot #0 works with a TChain
138     DefineInput(0, TChain::Class());
139     DefineOutput(1, TList::Class());
140     InitHistPointers();
141 }
142 //________________________________________________________________________
143 void AliEMCalpi0ClusterEvaluationTask::InitHistPointers() {
144
145     for (Int_t i=0; i < kNtype; i++)
146     {
147         fMasspi0EGA[i]        = 0;
148         fMassMixedEGA[i]      = 0;
149         fEventsEGA[i]         = 0;
150
151         fMasspi0MB[i]         = 0;
152         fMassMixedMB[i]       = 0;
153         fEventsMB[i]          = 0;
154
155         fMasspi0AllMB[i]      = 0;
156         fMassMixedAllMB[i]    = 0;
157         fEventsAllMB[i]       = 0;
158
159
160         fMasspi0Cent[i]      = 0;
161         fMassMixedCent[i]    = 0;
162         fEventsCent[i]       = 0;
163
164         fMasspi0SemiCent[i]  = 0;
165         fMassMixedSemiCent[i]= 0;
166         fEventsSemiCent[i]   = 0;
167
168
169         fpTMB[i]             = 0;
170         fpTEGA[i]            = 0;
171         fpTAllMB[i]          = 0;
172         fpTkCent[i]          = 0;
173         fpTkSemiCent[i]      = 0;
174
175
176         fPool[i]             = 0;
177     }
178
179     fCentrality              = 0;
180     fCentralityMB            = 0;
181     fCentralityEGA           = 0;
182     fCentralityCent          = 0;
183     fCentralitySemiCent      = 0;
184     fDispersion              = 0;
185     fexo                     = 0;
186     fTriggers                = 0;
187     fshower                  = 0;
188
189 }
190
191 //________________________________________________________________________
192 AliEMCalpi0ClusterEvaluationTask::~AliEMCalpi0ClusterEvaluationTask()
193 {
194     /// destructor
195     //only if they are initialized to 0, or to something valid
196
197     // Deleting the list (fHistList) also deletes all histograms stored in it
198     // because of fHistList->SetOwner(). It is therefore not necessary to delete
199     // the historgrams individually
200     delete fHistList;
201
202     // delete the 'running' pool array, which is either zero,
203     // or corresponds with the last created pool during event processing.
204     for (Int_t i=0; i < kNtype; i++)
205     { delete fPool[i]; }
206
207 }
208
209 //________________________________________________________________________
210 void AliEMCalpi0ClusterEvaluationTask::UserCreateOutputObjects( void )
211 {
212     fHistList = new TList();
213     fHistList->SetOwner();
214
215     for (Int_t i=0; i < kNtype; i++)
216     {
217         fMasspi0EGA[i]   = new TH2F(Form("MassEGA_%i",i),     Form("L1g m_{#gamma#gamma} centrality bin %i ;p_{T} (GeV/c)",i), 200,0,1, 150, 0, 30);
218         fHistList->Add( fMasspi0EGA[i]  );
219         fMassMixedEGA[i] = new TH2F(Form("MixedMassEGA_%i",i),Form("L1g mixed event m_{#gamma#gamma} centrality bin %i; p_{T} (GeV/c)",i), 200,0, 1, 150, 0, 30);
220         fHistList->Add( fMassMixedEGA[i] );
221
222         fEventsEGA[i]    = new TH1F(Form("EventsEGA_%i",i),   Form("L1g events in centrality bin %i",i),1, 0.5, 1.5);
223         fHistList->Add( fEventsEGA[i]  );
224
225         fMasspi0MB[i]    = new TH2F(Form("MassMB_%i",i),      Form("MB m_{#gamma#gamma} centrality bin %i ;p_{T} (GeV/c)",i),200,0,1, 150, 0, 30);
226         fHistList->Add( fMasspi0MB[i]  );
227         fMassMixedMB[i]  = new TH2F(Form("MixedMassMB_%i",i), Form("MB mixed event m_{#gamma#gamma} centrality bin %i;p_{T} (GeV/c)",i), 200,0, 1, 150, 0, 30);
228         fHistList->Add(fMassMixedMB[i]);
229
230         fEventsMB[i]     = new TH1F(Form("EventsMB_%i",i),    Form("MB events in centrality bin %i",i),1, 0.5, 1.5);
231         fHistList->Add( fEventsMB[i]  );
232
233         fMasspi0AllMB[i]  = new TH2F(Form("MassAllMB_%i",i),    Form("C+SC+AnyInt m_{#gamma#gamma} centrality bin %i ;p_{T} (GeV/c)",i),200,0,1, 150, 0, 30);
234         fHistList->Add( fMasspi0AllMB[i]  );
235
236         fMassMixedAllMB[i]= new TH2F(Form("MixedMassAllMB_%i",i),Form("C+SC+AnyInt mixed event m_{#gamma#gamma} centrality bin %i;p_{T} (GeV/c)",i), 200,0, 1, 150, 0, 30);
237         fHistList->Add(fMassMixedAllMB[i]);
238
239         fEventsAllMB[i]   = new TH1F(Form("EventsAllMB_%i",i),  Form("C+SC+ AnyInt Events in centrality bin %i",i),1, 0.5, 1.5);
240         fHistList->Add( fEventsAllMB[i]  );
241
242         fpTMB[i]   = new TH1F(Form("pTMB_%i",i),  Form("pT in centrality bin %i",i),150, 1, 30);
243         fHistList->Add( fpTMB[i]  );
244
245         fpTAllMB[i]   = new TH1F(Form("pTAllMB_%i",i),  Form("pT in centrality bin %i",i),150, 1, 30);
246         fHistList->Add( fpTAllMB[i]  );
247
248         fpTEGA[i]  = new TH1F(Form("pTEGA_%i",i),  Form("pT in centrality bin %i",i),150, 1, 30);
249         fHistList->Add( fpTEGA[i]  );
250
251         fpTkCent[i]  = new TH1F(Form("pTkCent_%i",i),  Form("pT in centrality bin %i",i),150, 1, 30);
252         fHistList->Add( fpTkCent[i]  );
253
254         fpTkSemiCent[i]  = new TH1F(Form("pTkSemiCent_%i",i),  Form("pT in centrality bin %i",i),150, 1, 30);
255         fHistList->Add( fpTkSemiCent[i]  );
256
257         //refined kCentral and SemiCentral
258
259         fMasspi0SemiCent[i]  = new TH2F(Form("MassSemiCent_%i",i),    Form("SemiCent m_{#gamma#gamma} centrality bin %i ;p_{T} (GeV/c)",i),200,0,1, 150, 0, 30);
260         fHistList->Add( fMasspi0SemiCent[i]  );
261
262         fMassMixedSemiCent[i]= new TH2F(Form("MixedMassSemiCent_%i",i),Form("SemiCent mixed event m_{#gamma#gamma} centrality bin %i;p_{T} (GeV/c)",i), 200,0, 1, 150, 0, 30);
263         fHistList->Add(fMassMixedSemiCent[i]);
264
265         fEventsSemiCent[i]   = new TH1F(Form("EventsSemiCent_%i",i),  Form("SemiCent Events in centrality bin %i",i),1, 0.5, 1.5);
266         fHistList->Add( fEventsSemiCent[i]  );
267
268         fMasspi0Cent[i]  = new TH2F(Form("MassCent_%i",i),    Form("C+SC+AnyInt m_{#gamma#gamma} centrality bin %i ;p_{T} (GeV/c)",i),200,0,1, 150, 0, 30);
269         fHistList->Add( fMasspi0Cent[i]  );
270
271         fMassMixedCent[i]= new TH2F(Form("MixedMassCent_%i",i),Form("C+SC+AnyInt mixed event m_{#gamma#gamma} centrality bin %i;p_{T} (GeV/c)",i), 200,0, 1, 150, 0, 30);
272         fHistList->Add(fMassMixedCent[i]);
273
274         fEventsCent[i]   = new TH1F(Form("EventsCent_%i",i),  Form("C+SC+ AnyInt Events in centrality bin %i",i),1, 0.5, 1.5);
275         fHistList->Add(fEventsCent[i]);
276
277
278     }
279
280     fCentrality   = new TH1F("centrality", "centrality",         100, 0, 100);
281     fHistList->Add(fCentrality);
282
283     fCentralityMB   = new TH1F("centralityMB", "centrality",     100, 0, 100);
284     fHistList->Add(fCentralityMB);
285
286     fCentralityEGA   = new TH1F("centralityEGA", "centrality",   100, 0, 100);
287     fHistList->Add(fCentralityEGA);
288
289     fCentralityCent   = new TH1F("centralityCent", "centrality", 100, 0, 100);
290     fHistList->Add(fCentralityCent);
291
292     fCentralitySemiCent   = new TH1F("centralitySemiCent", "centrality", 100, 0, 100);
293     fHistList->Add(fCentralitySemiCent);
294
295     fDispersion   = new TH1F("dispersion", "dispersion", 100, -1, 2);
296     fHistList->Add(fDispersion);
297
298     fexo   = new TH1F("exo", "exo", 100, -1, 2);
299     fHistList->Add(fexo);
300
301     fshower   = new TH1F("shower", "shower", 100, -1, 2);
302     fHistList->Add(fshower);
303
304
305     fTriggers     = new TH1F("triggers",   "triggers",   10,  0, 10);
306     fTriggers->GetXaxis()->SetBinLabel( 1,"All");
307     fTriggers->GetXaxis()->SetBinLabel( 2,"AllMB");
308     fTriggers->GetXaxis()->SetBinLabel( 3,"MB");
309     fTriggers->GetXaxis()->SetBinLabel( 4,"kCentral");
310     fTriggers->GetXaxis()->SetBinLabel( 5,"kSemiCentral");
311     fTriggers->GetXaxis()->SetBinLabel( 6,"EGA");
312     fTriggers->GetXaxis()->SetBinLabel( 7,"MB/ALL");
313     fTriggers->GetXaxis()->SetBinLabel( 8,"MB/EGA");
314     fTriggers->GetXaxis()->SetBinLabel( 9,"MB/KCentral");
315     fTriggers->GetXaxis()->SetBinLabel( 10,"MB/kSemiCentral");
316     fHistList->Add(fTriggers);
317
318     PostData(1, fHistList);
319
320 }
321
322 //________________________________________________________________________
323
324
325 //________________________________________________________________________
326 void AliEMCalpi0ClusterEvaluationTask::UserExec( Option_t* )
327 {
328
329     // increase event number
330     AliInfo( Form( "event: %i", fEvent++ ) );
331
332     // ESD Filter analysis task executed for each event
333
334     AliESDEvent* esd    = dynamic_cast<AliESDEvent*>(InputEvent());
335     AliAODEvent* aod    = dynamic_cast< AliAODEvent*>(InputEvent());
336     AliVEvent  * event  = InputEvent();
337     // === Physics Selection Task ===
338     // bitwise operation is used against.
339     // extraction de la masque de selection
340     if(esd){
341     isPileup       = esd->IsPileupFromSPD(3,0.8);
342     if(isPileup) return;
343            }
344
345     //Remove events with exotic clusters
346     TRefArray *caloClusArr=new TRefArray();
347     event->GetEMCALClusters(caloClusArr);
348     const Int_t kNumber =caloClusArr->GetEntries();
349
350     for( Int_t iclus = 0; iclus < kNumber; iclus++ ){
351
352         AliVCluster*c=(AliVCluster*) caloClusArr->At(iclus);
353         if(!c){
354             return;
355         }
356         if(!c->IsEMCAL()){
357             return;
358         }
359
360         Int_t id= -1;;
361         Double_t Emax   = GetMaxCellEnergy( c, id);
362
363         AliVCaloCells     *Cells       =  event->GetEMCALCells();
364         AliEMCALGeometry  *geom        =  AliEMCALGeometry::GetInstance("EMCAL_COMPLETEV1");
365         AliEMCALRecoUtils RecoUtils;
366         Int_t bc;
367         if(esd)  bc = esd->GetBunchCrossNumber();
368         if(aod) bc = aod->GetBunchCrossNumber();
369
370         Double_t tcell=0;
371         Double_t Ecross =  RecoUtils.GetECross(id, tcell,Cells, bc);
372         Double_t Exo    = 1.0 - Ecross/Emax;
373         fexo->Fill(Exo);
374         if((Exo)>1){ return;}
375         fshower->Fill(c->GetM02());
376
377     }
378     //888888888888888888888888888888888888888888888888888
379
380     const UInt_t Mask = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected();
381
382     // event characterization
383     if(esd){
384         runNumber     = esd->GetRunNumber();
385         TString triggerClasses = esd->GetFiredTriggerClasses();
386     }
387
388     if(aod){
389         runNumber     = aod->GetRunNumber();
390        // Get triggered classes, bordel....
391     TString triggerClasses = aod->GetFiredTriggerClasses();
392             }
393
394
395     ULong64_t triggerMask  = event->GetTriggerMask();
396
397     ////verification triggered classes that fired.
398
399     isMB           = Mask & AliVEvent::kMB;
400     isAnyINT       = Mask & AliVEvent::kAnyINT;
401     isCentral      = Mask & AliVEvent::kCentral;
402     isSemiCentral  = Mask & AliVEvent::kSemiCentral;
403     isEga          = Mask & AliVEvent::kEMCEGA;
404     /*
405     isMB             = (Mask== AliVEvent::kMB)? 1 : 0;
406     isAnyINT         = (Mask== AliVEvent::kAnyINT)? 1 : 0;
407     isCentral        = (Mask== AliVEvent::kCentral)? 1 : 0;
408     isSemiCentral    = (Mask== AliVEvent::kSemiCentral)? 1 : 0;
409     isEga            = (Mask== AliVEvent::kEMCEGA)? 1 : 0;
410     */
411     //if( isMB ||isCentral||isSemiCentral ){ kAllMB=kTRUE;}
412     //else kAllMB = kFALSE;
413     kAllMB = (isMB ||isCentral||isSemiCentral);
414
415   AliCentrality *esdCent;
416   AliCentrality *aodCent;
417     // Centrality
418 if(esd){
419
420     esdCent     = esd->GetCentrality();
421
422     //fleche pour pointers, otherwise point.
423     CentralityVZERO  = esdCent->GetCentralityPercentile("V0M");
424     CentralitySPD    = esdCent->GetCentralityPercentile("CL1");
425
426 }
427
428 if(aod){
429     aodCent     = aod->GetCentrality();
430     //fleche pour pointers, otherwise point.
431     CentralityVZERO  = aodCent->GetCentralityPercentile("V0M");
432     CentralitySPD    = aodCent->GetCentralityPercentile("CL1");
433 }
434
435     all++;
436     if(kAllMB){allmb++;fCentrality->Fill(CentralityVZERO);}
437     if(isMB){mb++;fCentralityMB->Fill(CentralityVZERO);}
438     if(isEga){ega++;fCentralityEGA->Fill(CentralityVZERO);}
439     if(isCentral){central++;fCentralityCent->Fill(CentralityVZERO);}
440     if(isSemiCentral){semicentral++;fCentralitySemiCent->Fill(CentralityVZERO);}
441
442
443     fTriggers->SetBinContent(1,  all);
444     fTriggers->SetBinContent(2,  allmb);
445     fTriggers->SetBinContent(3,  mb);
446     fTriggers->SetBinContent(4,  central);
447     fTriggers->SetBinContent(5,  semicentral);
448     fTriggers->SetBinContent(6,  ega);
449     //     fTriggers->SetBinContent(7,  mb/ega);
450     //     fTriggers->SetBinContent(8,  mb/central);
451     //     fTriggers->SetBinContent(9,  mb/semicentral);
452     //     fTriggers->SetBinContent(10, ega/mb);
453
454
455
456
457     // Fill per centrality class event count
458     if(CentralityVZERO>0  && CentralityVZERO<=5 )
459     {
460         if(isEga){ega0++;fEventsEGA[0]->Fill(ega0);}
461         if(isMB){mb0++;fEventsMB[0]->Fill(mb0);}
462         if(kAllMB){allmb0++;fEventsAllMB[0]->Fill(allmb0);}
463         if(isCentral){cent0++;fEventsCent[0]->Fill(cent0);}
464         if(isSemiCentral){semicent0++;fEventsSemiCent[0]->Fill(semicent0);}
465     }
466     if(CentralityVZERO>5  && CentralityVZERO<=10){if(isEga){ega1++;fEventsEGA[1]->Fill(ega1);}if(isMB){mb1++;fEventsMB[1]->Fill(mb1);}if(kAllMB){allmb1++;fEventsAllMB[1]->Fill(allmb1);} if(isCentral){cent1++;fEventsCent[1]->Fill(cent1);} if(isSemiCentral){semicent1++;fEventsSemiCent[1]->Fill(semicent1);}}
467     if(CentralityVZERO>10 && CentralityVZERO<=20){if(isEga){ega2++;fEventsEGA[2]->Fill(ega2);}if(isMB){mb2++;fEventsMB[2]->Fill(mb2);}if(kAllMB){allmb2++;fEventsAllMB[2]->Fill(allmb2);} if(isCentral){cent2++;fEventsCent[2]->Fill(cent2);} if(isSemiCentral){semicent2++;fEventsSemiCent[2]->Fill(semicent2);}}
468     if(CentralityVZERO>20 && CentralityVZERO<=40){if(isEga){ega3++;fEventsEGA[3]->Fill(ega3);}if(isMB){mb3++;fEventsMB[3]->Fill(mb3);}if(kAllMB){allmb3++;fEventsAllMB[3]->Fill(allmb3);} if(isCentral){cent3++;fEventsCent[3]->Fill(cent3);} if(isSemiCentral){semicent3++;fEventsSemiCent[3]->Fill(semicent3);}}
469     if(CentralityVZERO>40 && CentralityVZERO<=60){if(isEga){ega4++;fEventsEGA[4]->Fill(ega4);}if(isMB){mb4++;fEventsMB[4]->Fill(mb4);}if(kAllMB){allmb4++;fEventsAllMB[4]->Fill(allmb4);} if(isCentral){cent4++;fEventsCent[4]->Fill(cent4);} if(isSemiCentral){semicent4++;fEventsSemiCent[4]->Fill(semicent4);}}
470     if(CentralityVZERO>60 && CentralityVZERO<=80){if(isEga){ega5++;fEventsEGA[5]->Fill(ega5);}if(isMB){mb5++;fEventsMB[5]->Fill(mb5);}if(kAllMB){allmb5++;fEventsAllMB[5]->Fill(allmb5);} if(isCentral){cent5++;fEventsCent[5]->Fill(cent5);} if(isSemiCentral){semicent5++;fEventsSemiCent[5]->Fill(semicent5);}}
471     if(CentralityVZERO>0  && CentralityVZERO<=10){if(isEga){ega6++;fEventsEGA[6]->Fill(ega6);}if(isMB){mb6++;fEventsMB[6]->Fill(mb6);}if(kAllMB){allmb6++;fEventsAllMB[6]->Fill(allmb6);} if(isCentral){cent6++;fEventsCent[6]->Fill(cent6);} if(isSemiCentral){semicent6++;fEventsSemiCent[6]->Fill(semicent6);}}
472     if(CentralityVZERO>0  && CentralityVZERO<=20){if(isEga){ega7++;fEventsEGA[7]->Fill(ega7);}if(isMB){mb7++;fEventsMB[7]->Fill(mb7);}if(kAllMB){allmb7++;fEventsAllMB[7]->Fill(allmb7);} if(isCentral){cent7++;fEventsCent[7]->Fill(cent7);} if(isSemiCentral){semicent7++;fEventsSemiCent[7]->Fill(semicent7);}}
473     if(CentralityVZERO>40 && CentralityVZERO<=50){if(isEga){ega8++;fEventsEGA[8]->Fill(ega8);}if(isMB){mb8++;fEventsMB[8]->Fill(mb8);}if(kAllMB){allmb8++;fEventsAllMB[8]->Fill(allmb8);} if(isCentral){cent8++;fEventsCent[8]->Fill(cent8);} if(isSemiCentral){semicent8++;fEventsSemiCent[8]->Fill(semicent8);}}
474     if(CentralityVZERO>50 && CentralityVZERO<=60){if(isEga){ega9++;fEventsEGA[9]->Fill(ega9);}if(isMB){mb9++;fEventsMB[9]->Fill(mb9);}if(kAllMB){allmb9++;fEventsAllMB[9]->Fill(allmb9);} if(isCentral){cent9++;fEventsCent[9]->Fill(cent9);} if(isSemiCentral){semicent9++;fEventsSemiCent[9]->Fill(semicent9);}}
475
476
477
478
479     //Pass the geometry transformation matrix from ESDs to geometry
480     AliEMCALGeometry  *geom        =  AliEMCALGeometry::GetInstance("EMCAL_COMPLETEV1");
481     AliVCaloCells     *Cells       =  event->GetEMCALCells();
482     AliEMCALRecoUtils RecoUtils;
483
484     Int_t    absID1    = -1;        Int_t    absID2 = -1;
485     Int_t    ieta1     = -1;        Int_t    ieta2  = -1;
486     Int_t    iphi1     = -1;        Int_t    iphi2  = -1;
487     Int_t    iSM1      = -1;        Int_t    iSM2   = -1;
488     Bool_t  shared1;      Bool_t  shared2;
489
490     //get reconstructed vertex position
491     Double_t vertex_position[3];
492     if(esd)esd->GetVertex()->GetXYZ(vertex_position);
493     vX     = vertex_position[0];
494     vY     = vertex_position[1];
495     vZ     = vertex_position[2];
496
497
498     if(aod){
499     vX     =0.0;
500     vY     =0.0;
501     vZ     =0.0;
502 }
503
504
505     // cout<<vZ<<endl;
506     if(vZ>-15.||vZ<15.){//vertex cut
507
508         //array temporaire pour passe plus tard dans le boucles
509         TRefArray caloClustersArr = TRefArray();
510         event->GetEMCALClusters( &caloClustersArr );
511
512         const Int_t kNumberOfEMCALClusters =caloClustersArr.GetEntries();
513
514         // pool to store clusters to be mixed
515         // it is not deleted directly. Instead, it is assigned to fPool, the previously assigned pool
516         // being deleted (if any), at that time.
517         TObjArray *newPool = new TObjArray(kNumberOfEMCALClusters);
518         newPool->SetOwner();
519
520         TVector3 pos;
521         pos -= vertex_position;
522         Double_t r1 = pos.Mag();
523
524         //Fill the pool with all clusters
525         Int_t nGoodClusters = 0;
526         //boucle sur tous les clusters
527         for( Int_t iclus = 0; iclus < kNumberOfEMCALClusters-1; iclus++ )
528         {//first cluster
529
530             AliVCluster*c1=(AliVCluster*) caloClustersArr.At(iclus);
531             if (!c1) continue;
532             if (!c1->IsEMCAL()) continue;
533             if (c1->GetNCells()<2) continue;
534             if (c1->E()<2) continue;
535             if(c1->GetM02()>0.5) continue;
536             if (c1->GetDistanceToBadChannel()<2) continue;
537
538             TLorentzVector pii;
539             c1->GetMomentum(pii, vertex_position);
540             TLorentzVector *pimx = new TLorentzVector;
541             c1->GetMomentum(*pimx, vertex_position);
542
543             // add TLorentz vector in the pool array
544             // it will be deleted when the array is deleted
545             newPool->Add(pimx);
546             ++nGoodClusters;
547
548             //characteristiques cluster
549             Ecluster               = c1->E();
550             NCellscluster          = c1->GetNCells();
551             M20cluster             = c1->GetM20();
552             M02cluster             = c1->GetM02();
553             NCluscluster           = kNumberOfEMCALClusters;
554             isEMCALcluster         = c1->IsEMCAL();
555             dispersioncluster      = c1->GetDispersion();
556             chi2cluster            = c1->Chi2();
557             distBadChannelcluster  = c1->GetDistanceToBadChannel();
558             phicluster             = pii.Phi();
559             etacluster             = pii.Eta();
560             ptcluster              = pii.Pt();
561             RecoUtils.GetMaxEnergyCell(geom, Cells, c1, absID1, iSM1, ieta1, iphi1, shared1);
562
563             fDispersion->Fill(dispersioncluster);
564
565             for (Int_t iclus2 = iclus+1;  iclus2< kNumberOfEMCALClusters-1; iclus2++)
566             {//second cluster
567
568                 AliVCluster *c2 = (AliVCluster*) caloClustersArr.At(iclus2);
569                 if (!c2) continue;
570                 if (!c2->IsEMCAL()) continue;
571                 if (c2->GetNCells()<2) continue;
572                 if (c2->E()<2) continue;
573                 if(c2->GetM02()>0.5) continue;
574                 if (c2->GetDistanceToBadChannel()<2) continue;
575                 Float_t en2 = c2->E();
576
577                 TLorentzVector pjj;
578                 c2->GetMomentum(pjj, vertex_position);
579                 TLorentzVector pion;
580                 pion   = pii + pjj;
581                 //remplissage des pions
582                 piE      = pion.E();
583                 piphi    = pion.Phi();
584                 pieta    = pion.Eta();
585                 ptpi     = pion.Pt();
586                 pipx     = pion.Px();
587                 pipy     = pion.Py();
588                 pipz     = pion.Pz();
589                 asympi     = TMath::Abs(pii.E()-pjj.E())/(pii.E()+pjj.E());
590                 masspi   = pion.M();
591
592                 if(CentralityVZERO>0  && CentralityVZERO<=5 ){if(isEga){fpTEGA[0]->Fill(ptpi);fMasspi0EGA[0]->Fill(masspi, ptpi);}if(isMB){fpTMB[0]->Fill(ptpi); fMasspi0MB[0]->Fill(masspi, ptpi);}if(kAllMB){fpTAllMB[0]->Fill(ptpi);fMasspi0AllMB[0]->Fill(masspi, ptpi);}if(isCentral){fpTkCent[0]->Fill(ptpi);fMasspi0Cent[0]->Fill(masspi, ptpi);}if(isSemiCentral){fpTkSemiCent[0]->Fill(ptpi);fMasspi0SemiCent[0]->Fill(masspi, ptpi);}}
593                 if(CentralityVZERO>5  && CentralityVZERO<=10){if(isEga){fpTEGA[1]->Fill(ptpi);fMasspi0EGA[1]->Fill(masspi, ptpi);}if(isMB){fpTMB[1]->Fill(ptpi); fMasspi0MB[1]->Fill(masspi, ptpi);}if(kAllMB){fpTAllMB[1]->Fill(ptpi);fMasspi0AllMB[1]->Fill(masspi, ptpi);}if(isCentral){fpTkCent[1]->Fill( ptpi);fMasspi0Cent[1]->Fill(masspi, ptpi);} if(isSemiCentral){fpTkSemiCent[1]->Fill( ptpi);fMasspi0SemiCent[1]->Fill(masspi, ptpi);}}
594                 if(CentralityVZERO>10 && CentralityVZERO<=20){if(isEga){fpTEGA[2]->Fill(ptpi);fMasspi0EGA[2]->Fill(masspi, ptpi);}if(isMB){fpTMB[2]->Fill(ptpi); fMasspi0MB[2]->Fill(masspi, ptpi);}if(kAllMB){fpTAllMB[2]->Fill(ptpi);fMasspi0AllMB[2]->Fill(masspi, ptpi);}if(isCentral){fpTkCent[2]->Fill( ptpi);fMasspi0Cent[2]->Fill(masspi, ptpi);} if(isSemiCentral){fpTkSemiCent[2]->Fill( ptpi);fMasspi0SemiCent[2]->Fill(masspi, ptpi);}}
595                 if(CentralityVZERO>20 && CentralityVZERO<=40){if(isEga){fpTEGA[3]->Fill(ptpi);fMasspi0EGA[3]->Fill(masspi, ptpi);}if(isMB){fpTMB[3]->Fill(ptpi); fMasspi0MB[3]->Fill(masspi, ptpi);}if(kAllMB){fpTAllMB[3]->Fill(ptpi);fMasspi0AllMB[3]->Fill(masspi, ptpi);}if(isCentral){fpTkCent[3]->Fill( ptpi);fMasspi0Cent[3]->Fill(masspi, ptpi);} if(isSemiCentral){fpTkSemiCent[3]->Fill( ptpi);fMasspi0SemiCent[3]->Fill(masspi, ptpi);}}
596                 if(CentralityVZERO>40 && CentralityVZERO<=60){if(isEga){fpTEGA[4]->Fill(ptpi);fMasspi0EGA[4]->Fill(masspi, ptpi);}if(isMB){fpTMB[4]->Fill(ptpi); fMasspi0MB[4]->Fill(masspi, ptpi);}if(kAllMB){fpTAllMB[4]->Fill(ptpi);fMasspi0AllMB[4]->Fill(masspi, ptpi);}if(isCentral){fpTkCent[4]->Fill( ptpi);fMasspi0Cent[4]->Fill(masspi, ptpi);} if(isSemiCentral){fpTkSemiCent[4]->Fill( ptpi);fMasspi0SemiCent[4]->Fill(masspi, ptpi);}}
597                 if(CentralityVZERO>60 && CentralityVZERO<=80){if(isEga){fpTEGA[5]->Fill(ptpi);fMasspi0EGA[5]->Fill(masspi, ptpi);}if(isMB){fpTMB[5]->Fill(ptpi); fMasspi0MB[5]->Fill(masspi, ptpi);}if(kAllMB){fpTAllMB[5]->Fill(ptpi);fMasspi0AllMB[5]->Fill(masspi, ptpi);}if(isCentral){fpTkCent[5]->Fill( ptpi);fMasspi0Cent[5]->Fill(masspi, ptpi);} if(isSemiCentral){fpTkSemiCent[5]->Fill( ptpi);fMasspi0SemiCent[5]->Fill(masspi, ptpi);}}
598                 if(CentralityVZERO>0  && CentralityVZERO<=10){if(isEga){fpTEGA[6]->Fill(ptpi);fMasspi0EGA[6]->Fill(masspi, ptpi);}if(isMB){fpTMB[6]->Fill(ptpi); fMasspi0MB[6]->Fill(masspi, ptpi);}if(kAllMB){fpTAllMB[6]->Fill(ptpi);fMasspi0AllMB[6]->Fill(masspi, ptpi);}if(isCentral){fpTkCent[6]->Fill( ptpi);fMasspi0Cent[6]->Fill(masspi, ptpi);} if(isSemiCentral){fpTkSemiCent[6]->Fill( ptpi);fMasspi0SemiCent[6]->Fill(masspi, ptpi);}}
599                 if(CentralityVZERO>0  && CentralityVZERO<=20){if(isEga){fpTEGA[7]->Fill(ptpi);fMasspi0EGA[7]->Fill(masspi, ptpi);}if(isMB){fpTMB[7]->Fill(ptpi); fMasspi0MB[7]->Fill(masspi, ptpi);}if(kAllMB){fpTAllMB[7]->Fill(ptpi);fMasspi0AllMB[7]->Fill(masspi, ptpi);}if(isCentral){fpTkCent[7]->Fill( ptpi);fMasspi0Cent[7]->Fill(masspi, ptpi);} if(isSemiCentral){fpTkSemiCent[7]->Fill( ptpi);fMasspi0SemiCent[7]->Fill(masspi, ptpi);}}
600                 if(CentralityVZERO>40 && CentralityVZERO<=50){if(isEga){fpTEGA[8]->Fill(ptpi);fMasspi0EGA[8]->Fill(masspi, ptpi);}if(isMB){fpTMB[8]->Fill(ptpi); fMasspi0MB[8]->Fill(masspi, ptpi);}if(kAllMB){fpTAllMB[8]->Fill(ptpi);fMasspi0AllMB[8]->Fill(masspi, ptpi);}if(isCentral){fpTkCent[8]->Fill( ptpi);fMasspi0Cent[8]->Fill(masspi, ptpi);} if(isSemiCentral){fpTkSemiCent[8]->Fill( ptpi);fMasspi0SemiCent[8]->Fill(masspi, ptpi);}}
601                 if(CentralityVZERO>50 && CentralityVZERO<=60){if(isEga){fpTEGA[9]->Fill(ptpi);fMasspi0EGA[9]->Fill(masspi, ptpi);}if(isMB){fpTMB[9]->Fill(ptpi); fMasspi0MB[9]->Fill(masspi, ptpi);}if(kAllMB){fpTAllMB[9]->Fill(ptpi);fMasspi0AllMB[9]->Fill(masspi, ptpi);}if(isCentral){fpTkCent[9]->Fill( ptpi);fMasspi0Cent[9]->Fill(masspi, ptpi);} if(isSemiCentral){fpTkSemiCent[9]->Fill( ptpi);fMasspi0SemiCent[9]->Fill(masspi, ptpi);}}
602
603                 RecoUtils.GetMaxEnergyCell(geom, Cells, c2, absID2, iSM2,  ieta2,  iphi2, shared2);
604
605             } //Deuxieme cluster
606
607         } //Premier cluster
608
609         //     //-----------------------------------------
610         //fill mixed event
611         int poolIndex = 0;
612         if(CentralityVZERO>0  && CentralityVZERO<=5 ) poolIndex =0;
613         if(CentralityVZERO>5  && CentralityVZERO<=10) poolIndex =1;
614         if(CentralityVZERO>10 && CentralityVZERO<=20) poolIndex =2;
615         if(CentralityVZERO>20 && CentralityVZERO<=40) poolIndex =3;
616         if(CentralityVZERO>40 && CentralityVZERO<=60) poolIndex =4;
617         if(CentralityVZERO>60 && CentralityVZERO<=80) poolIndex =5;
618         if(CentralityVZERO>0  && CentralityVZERO<=10) poolIndex =6;
619         if(CentralityVZERO>0  && CentralityVZERO<=20) poolIndex =7;
620         if(CentralityVZERO>40 && CentralityVZERO<=50) poolIndex =8;
621         if(CentralityVZERO>50 && CentralityVZERO<=60) poolIndex =9;
622
623
624         if( !fPool[poolIndex] )
625         {
626             fPool[poolIndex] = newPool;
627
628         } else {
629
630             Int_t nGoodClusters2 = fPool[poolIndex]->GetEntries();
631             for (Int_t i=0; i<nGoodClusters; ++i)
632             {
633
634
635                 TLorentzVector * pi = static_cast<TLorentzVector*>(newPool->At(i));
636                 for (Int_t j=0; j<nGoodClusters2; ++j)
637                 {
638
639                     TLorentzVector *pj= static_cast<TLorentzVector*>(fPool[poolIndex]->At(j));
640                     FillMixed(*pi,*pj);
641
642
643                 } //random cluster
644             } //current cluster
645
646             // delete previous pool and assign to new pool.
647             delete fPool[poolIndex];
648             fPool[poolIndex] = newPool;
649
650         }
651     }//15 cm vertex cut
652     PostData(1, fHistList);
653
654 }//end process
655
656 //________________________________________________________________________
657 Double_t AliEMCalpi0ClusterEvaluationTask ::GetMaxCellEnergy(const AliVCluster *cluster, Int_t &id) const
658 {
659     // Get maximum energy of attached cell.
660     AliESDEvent* esd  = dynamic_cast<AliESDEvent*>(InputEvent());
661     AliAODEvent* aod =dynamic_cast< AliAODEvent*>(InputEvent());
662
663
664
665
666     AliEMCALGeometry  *fGeom        =  AliEMCALGeometry::GetInstance("EMCAL_COMPLETEV1");
667
668     id = -1;
669
670     AliVCaloCells *cells = 0;
671     if(esd){
672     cells = esd->GetEMCALCells();
673      }
674
675     if(!esd){
676     cells = aod->GetEMCALCells();
677      }
678
679     if (!cells)
680         return 0;
681
682     Double_t maxe = 0;
683     Int_t ncells = cluster->GetNCells();
684     for (Int_t i=0; i<ncells; i++) {
685         Double_t e = cells->GetCellAmplitude(TMath::Abs(cluster->GetCellAbsId(i)));
686         if (e>maxe) {
687             maxe = e;
688             id   = cluster->GetCellAbsId(i);
689         }
690     }
691     return maxe;
692 }
693 //________________________________________________________________
694 void AliEMCalpi0ClusterEvaluationTask::FillMixed( const TLorentzVector& p1, const TLorentzVector& p2)
695
696 {
697     //verification triggered classes that fired.
698     const UInt_t eventSelectionMask( ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected() );
699     AliESDEvent* esd     = dynamic_cast<AliESDEvent*>(InputEvent());
700     AliAODEvent* aod     =dynamic_cast< AliAODEvent*>(InputEvent());
701     AliVEvent  * event   = InputEvent();
702
703
704     ULong64_t triggerMask = event->GetTriggerMask();
705
706     isMBmx         = (eventSelectionMask&AliVEvent::kMB);
707     isAnyINTmx      = (eventSelectionMask&AliVEvent::kAnyINT);
708     isCentralmx     = (eventSelectionMask&AliVEvent::kCentral);
709     isSemiCentralmx = (eventSelectionMask&AliVEvent::kSemiCentral);
710     isEgamx         = (eventSelectionMask&AliVEvent::kEMCEGA);
711
712
713
714     //     isMBmx          = (eventSelectionMask== AliVEvent::kMB)? 1 : 0;
715     //     isAnyINTmx      = (eventSelectionMask== AliVEvent::kAnyINT)? 1 : 0;
716     //     isCentralmx     = (eventSelectionMask== AliVEvent::kCentral)? 1 : 0;
717     //     isSemiCentralmx = (eventSelectionMask== AliVEvent::kSemiCentral)? 1 : 0;
718     //     isEgamx         = (eventSelectionMask== AliVEvent::kEMCEGA)? 1 : 0;
719     kAllMBmx= (isMBmx || isAnyINTmx||isCentralmx||isSemiCentralmx );
720
721     // Centrality
722      AliCentrality *esdCent;
723     if(esd){
724      esdCent     = esd->GetCentrality();
725      }
726     if(!esd){
727      esdCent     = aod->GetCentrality();
728      }
729
730     //fleche pour pointers, otherwise point.
731     CentralityVZERO  = esdCent->GetCentralityPercentile("V0M");
732     CentralitySPD    = esdCent->GetCentralityPercentile("CL1");
733
734     TLorentzVector Mxpion;
735     Mxpion = p1 + p2;
736     Double_t Mxmass      = Mxpion.M();
737     Double_t Mxpt        = Mxpion.Pt();
738
739     if(Mxpt>=4){
740         if(CentralityVZERO>0  && CentralityVZERO<=5 ){if(isEgamx){fMassMixedEGA[0]->Fill(Mxmass, Mxpt);}if(isMBmx){fMassMixedMB[0]->Fill(Mxmass, Mxpt);}if(kAllMB){fMassMixedAllMB[0]->Fill(Mxmass, Mxpt);} if(isCentralmx){fMassMixedCent[0]->Fill(Mxmass, Mxpt);}if(isSemiCentralmx){fMassMixedSemiCent[0]->Fill(Mxmass, Mxpt);}}
741         if(CentralityVZERO>5  && CentralityVZERO<=10){if(isEgamx){fMassMixedEGA[1]->Fill(Mxmass, Mxpt);}if(isMBmx){fMassMixedMB[1]->Fill(Mxmass, Mxpt);}if(kAllMB){fMassMixedAllMB[1]->Fill(Mxmass, Mxpt);} if(isCentralmx){fMassMixedCent[1]->Fill(Mxmass, Mxpt);}if(isSemiCentralmx){fMassMixedSemiCent[1]->Fill(Mxmass, Mxpt);}}
742         if(CentralityVZERO>10 && CentralityVZERO<=20){if(isEgamx){fMassMixedEGA[2]->Fill(Mxmass, Mxpt);}if(isMBmx){fMassMixedMB[2]->Fill(Mxmass, Mxpt);}if(kAllMB){fMassMixedAllMB[2]->Fill(Mxmass, Mxpt);} if(isCentralmx){fMassMixedCent[2]->Fill(Mxmass, Mxpt);}if(isSemiCentralmx){fMassMixedSemiCent[2]->Fill(Mxmass, Mxpt);}}
743         if(CentralityVZERO>20 && CentralityVZERO<=40){if(isEgamx){fMassMixedEGA[3]->Fill(Mxmass, Mxpt);}if(isMBmx){fMassMixedMB[3]->Fill(Mxmass, Mxpt);}if(kAllMB){fMassMixedAllMB[3]->Fill(Mxmass, Mxpt);} if(isCentralmx){fMassMixedCent[3]->Fill(Mxmass, Mxpt);}if(isSemiCentralmx){fMassMixedSemiCent[3]->Fill(Mxmass, Mxpt);}}
744         if(CentralityVZERO>40 && CentralityVZERO<=60){if(isEgamx){fMassMixedEGA[4]->Fill(Mxmass, Mxpt);}if(isMBmx){fMassMixedMB[4]->Fill(Mxmass, Mxpt);}if(kAllMB){fMassMixedAllMB[4]->Fill(Mxmass, Mxpt);} if(isCentralmx){fMassMixedCent[4]->Fill(Mxmass, Mxpt);}if(isSemiCentralmx){fMassMixedSemiCent[4]->Fill(Mxmass, Mxpt);}}
745         if(CentralityVZERO>60 && CentralityVZERO<=80){if(isEgamx){fMassMixedEGA[5]->Fill(Mxmass, Mxpt);}if(isMBmx){fMassMixedMB[5]->Fill(Mxmass, Mxpt);}if(kAllMB){fMassMixedAllMB[5]->Fill(Mxmass, Mxpt);} if(isCentralmx){fMassMixedCent[5]->Fill(Mxmass, Mxpt);}if(isSemiCentralmx){fMassMixedSemiCent[5]->Fill(Mxmass, Mxpt);}}
746         if(CentralityVZERO>0  && CentralityVZERO<=10){if(isEgamx){fMassMixedEGA[6]->Fill(Mxmass, Mxpt);}if(isMBmx){fMassMixedMB[6]->Fill(Mxmass, Mxpt);}if(kAllMB){fMassMixedAllMB[6]->Fill(Mxmass, Mxpt);} if(isCentralmx){fMassMixedCent[6]->Fill(Mxmass, Mxpt);}if(isSemiCentralmx){fMassMixedSemiCent[6]->Fill(Mxmass, Mxpt);}}
747         if(CentralityVZERO>0  && CentralityVZERO<=20){if(isEgamx){fMassMixedEGA[7]->Fill(Mxmass, Mxpt);}if(isMBmx){fMassMixedMB[7]->Fill(Mxmass, Mxpt);}if(kAllMB){fMassMixedAllMB[7]->Fill(Mxmass, Mxpt);} if(isCentralmx){fMassMixedCent[7]->Fill(Mxmass, Mxpt);}if(isSemiCentralmx){fMassMixedSemiCent[7]->Fill(Mxmass, Mxpt);}}
748         if(CentralityVZERO>40 && CentralityVZERO<=50){if(isEgamx){fMassMixedEGA[8]->Fill(Mxmass, Mxpt);}if(isMBmx){fMassMixedMB[8]->Fill(Mxmass, Mxpt);}if(kAllMB){fMassMixedAllMB[8]->Fill(Mxmass, Mxpt);} if(isCentralmx){fMassMixedCent[8]->Fill(Mxmass, Mxpt);}if(isSemiCentralmx){fMassMixedSemiCent[8]->Fill(Mxmass, Mxpt);}}
749         if(CentralityVZERO>50 && CentralityVZERO<=60){if(isEgamx){fMassMixedEGA[9]->Fill(Mxmass, Mxpt);}if(isMBmx){fMassMixedMB[9]->Fill(Mxmass, Mxpt);}if(kAllMB){fMassMixedAllMB[9]->Fill(Mxmass, Mxpt);} if(isCentralmx){fMassMixedCent[9]->Fill(Mxmass, Mxpt);}if(isSemiCentralmx){fMassMixedSemiCent[9]->Fill(Mxmass, Mxpt);} }
750     }
751     //PostData(1, fHistList);
752 }
753
754 //________________________________________________________________________
755 void AliEMCalpi0ClusterEvaluationTask::Terminate(const Option_t*)
756 {
757     fHistList = dynamic_cast<TList*> (GetOutputData(1));
758     if (!fHistList) {
759         printf("ERROR: Output list not available\n");
760         return;
761     }
762
763 }