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