]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGGA/CaloTrackCorrelations/AliAnaPi0EbE.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
new parametrization of cluster selection cuts before and after splitting, new setting...
[u/mrichter/AliRoot.git] / PWGGA / CaloTrackCorrelations / AliAnaPi0EbE.cxx
CommitLineData
477d6cee 1/**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
21a4b1c0 8 * documentation strictly for non-commercial purposes is hereby granted *
477d6cee 9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
477d6cee 15
16//_________________________________________________________________________
17// Class for the analysis of high pT pi0 event by event
09273901 18// Pi0/Eta identified by one of the following:
477d6cee 19// -Invariant mass of 2 cluster in calorimeter
20// -Shower shape analysis in calorimeter
21a4b1c0 21// -Invariant mass of one cluster in calorimeter and one photon reconstructed in CTS
477d6cee 22//
23// -- Author: Gustavo Conesa (LNF-INFN) & Raphaelle Ichou (SUBATECH)
24//////////////////////////////////////////////////////////////////////////////
85c4406e 25
26
27// --- ROOT system ---
477d6cee 28#include <TList.h>
29#include <TClonesArray.h>
0c1383b5 30#include <TObjString.h>
477d6cee 31
85c4406e 32// --- Analysis system ---
33#include "AliAnaPi0EbE.h"
477d6cee 34#include "AliCaloTrackReader.h"
35#include "AliIsolationCut.h"
36#include "AliNeutralMesonSelection.h"
37#include "AliCaloPID.h"
38#include "AliMCAnalysisUtils.h"
477d6cee 39#include "AliStack.h"
ff45398a 40#include "AliFiducialCut.h"
477d6cee 41#include "TParticle.h"
0ae57829 42#include "AliVCluster.h"
2ad19c3d 43#include "AliESDEvent.h"
477d6cee 44#include "AliAODEvent.h"
591cc579 45#include "AliAODMCParticle.h"
477d6cee 46
47ClassImp(AliAnaPi0EbE)
34c16486 48
85c4406e 49//____________________________
50AliAnaPi0EbE::AliAnaPi0EbE() :
51AliAnaCaloTrackCorrBaseClass(),fAnaType(kIMCalo), fCalorimeter(""),
52fMinDist(0.),fMinDist2(0.), fMinDist3(0.),
53fNLMCutMin(-1), fNLMCutMax(10),
54fTimeCutMin(-10000), fTimeCutMax(10000),
55fRejectTrackMatch(kTRUE),
56fFillPileUpHistograms(0),
57fFillWeightHistograms(kFALSE), fFillTMHisto(0),
58fFillSelectClHisto(0), fFillOnlySimpleSSHisto(1), fFillEMCALBCHistograms(0),
59fInputAODGammaConvName(""),
1253480f 60fCheckSplitDistToBad(0),
85c4406e 61// Histograms
62fhPt(0), fhE(0),
63fhEEta(0), fhEPhi(0),
64fhPtEta(0), fhPtPhi(0), fhEtaPhi(0),
65fhEtaPhiEMCALBC0(0), fhEtaPhiEMCALBC1(0), fhEtaPhiEMCALBCN(0),
126b8c62 66fhTimeTriggerEMCALBC0UMReMatchOpenTime(0),
67fhTimeTriggerEMCALBC0UMReMatchCheckNeigh(0),
68fhTimeTriggerEMCALBC0UMReMatchBoth(0),
85c4406e 69fhPtCentrality(), fhPtEventPlane(0),
70fhPtReject(0), fhEReject(0),
71fhEEtaReject(0), fhEPhiReject(0), fhEtaPhiReject(0),
72fhMass(0), fhMassPt(0), fhMassSplitPt(0),
73fhSelectedMass(0), fhSelectedMassPt(0), fhSelectedMassSplitPt(0),
1253480f 74fhMassNoOverlap(0), fhMassPtNoOverlap(0), fhMassSplitPtNoOverlap(0),
75fhSelectedMassNoOverlap(0), fhSelectedMassPtNoOverlap(0), fhSelectedMassSplitPtNoOverlap(0),
76fhMCPi0SplitPtRecoPtPrim(0), fhMCEtaSplitPtRecoPtPrim(0),
77fhMCPi0SplitPtRecoPtPrimNoOverlap(0), fhMCEtaSplitPtRecoPtPrimNoOverlap(0),
78fhMCPi0SelectedPtRecoPtPrim(0), fhMCEtaSelectedPtRecoPtPrim(0),
79fhMCPi0SelectedPtRecoPtPrimNoOverlap(0), fhMCEtaSelectedPtRecoPtPrimNoOverlap(0),
80fhMCPi0SelectedSplitPtRecoPtPrim(0), fhMCEtaSelectedSplitPtRecoPtPrim(0),
81fhMCPi0SelectedSplitPtRecoPtPrimNoOverlap(0), fhMCEtaSelectedSplitPtRecoPtPrimNoOverlap(0),
85c4406e 82fhAsymmetry(0), fhSelectedAsymmetry(0),
83fhSplitE(0), fhSplitPt(0),
84fhSplitPtEta(0), fhSplitPtPhi(0),
85fhNLocMaxSplitPt(0),
86fhPtDecay(0), fhEDecay(0),
87// Shower shape histos
88fhEDispersion(0), fhELambda0(0), fhELambda1(0),
89fhELambda0NoTRD(0), fhELambda0FracMaxCellCut(0),
90fhEFracMaxCell(0), fhEFracMaxCellNoTRD(0),
91fhENCells(0), fhETime(0), fhEPairDiffTime(0),
92fhDispEtaE(0), fhDispPhiE(0),
93fhSumEtaE(0), fhSumPhiE(0), fhSumEtaPhiE(0),
94fhDispEtaPhiDiffE(0), fhSphericityE(0),
95
96// MC histos
97fhMCE(), fhMCPt(),
98fhMCPhi(), fhMCEta(),
99fhMCEReject(), fhMCPtReject(),
100fhMCPtCentrality(),
101fhMCPi0PtGenRecoFraction(0), fhMCEtaPtGenRecoFraction(0),
102fhMCPi0DecayPt(0), fhMCPi0DecayPtFraction(0),
103fhMCEtaDecayPt(0), fhMCEtaDecayPtFraction(0),
104fhMCOtherDecayPt(0),
105fhMassPairMCPi0(0), fhMassPairMCEta(0),
106fhAnglePairMCPi0(0), fhAnglePairMCEta(0),
107// Weight studies
108fhECellClusterRatio(0), fhECellClusterLogRatio(0),
109fhEMaxCellClusterRatio(0), fhEMaxCellClusterLogRatio(0),
110fhTrackMatchedDEta(0), fhTrackMatchedDPhi(0), fhTrackMatchedDEtaDPhi(0),
b2e375c7 111fhTrackMatchedDEtaPos(0), fhTrackMatchedDPhiPos(0), fhTrackMatchedDEtaDPhiPos(0),
112fhTrackMatchedDEtaNeg(0), fhTrackMatchedDPhiNeg(0), fhTrackMatchedDEtaDPhiNeg(0),
85c4406e 113fhTrackMatchedMCParticleE(0),
114fhTrackMatchedMCParticleDEta(0), fhTrackMatchedMCParticleDPhi(0),
115fhdEdx(0), fhEOverP(0), fhEOverPNoTRD(0),
116// Number of local maxima in cluster
117fhNLocMaxE(0), fhNLocMaxPt(0),
118// PileUp
126b8c62 119fhTimePtNoCut(0), fhTimePtSPD(0), fhTimePtSPDMulti(0),
85c4406e 120fhTimeNPileUpVertSPD(0), fhTimeNPileUpVertTrack(0),
121fhTimeNPileUpVertContributors(0),
122fhTimePileUpMainVertexZDistance(0), fhTimePileUpMainVertexZDiamond(0),
123fhPtNPileUpSPDVtx(0), fhPtNPileUpTrkVtx(0),
124fhPtNPileUpSPDVtxTimeCut(0), fhPtNPileUpTrkVtxTimeCut(0),
125fhPtNPileUpSPDVtxTimeCut2(0), fhPtNPileUpTrkVtxTimeCut2(0)
477d6cee 126{
127 //default ctor
128
34c16486 129 for(Int_t i = 0; i < 6; i++)
130 {
40d3ce60 131 fhMCE [i] = 0;
3455f821 132 fhMCPt [i] = 0;
6e66993c 133 fhMCNLocMaxPt [i] = 0;
134 fhMCPhi [i] = 0;
3455f821 135 fhMCEta [i] = 0;
17f5b4b6 136 fhMCPtCentrality [i] = 0;
137
cfdf2b91 138 fhMCSplitE [i] = 0;
139 fhMCSplitPt [i] = 0;
29250849 140 fhMCSplitPtPhi [i] = 0;
141 fhMCSplitPtEta [i] = 0;
6e66993c 142 fhMCNLocMaxSplitPt [i] = 0;
cfdf2b91 143
34c16486 144 fhEMCLambda0 [i] = 0;
145 fhEMCLambda0NoTRD [i] = 0;
3bfcb597 146 fhEMCLambda0FracMaxCellCut[i]= 0;
34c16486 147 fhEMCFracMaxCell [i] = 0;
148 fhEMCLambda1 [i] = 0;
149 fhEMCDispersion [i] = 0;
150
bfdcf7fb 151 fhMCEDispEta [i] = 0;
152 fhMCEDispPhi [i] = 0;
153 fhMCESumEtaPhi [i] = 0;
154 fhMCEDispEtaPhiDiff[i] = 0;
85c4406e 155 fhMCESphericity [i] = 0;
156 fhMCEAsymmetry [i] = 0;
157
29250849 158 fhMCMassPt [i]=0;
159 fhMCMassSplitPt [i]=0;
160 fhMCSelectedMassPt [i]=0;
161 fhMCSelectedMassSplitPt[i]=0;
162
1253480f 163 fhMCMassPtNoOverlap [i]=0;
164 fhMCMassSplitPtNoOverlap [i]=0;
165 fhMCSelectedMassPtNoOverlap [i]=0;
166 fhMCSelectedMassSplitPtNoOverlap[i]=0;
167
d2655d46 168 for(Int_t j = 0; j < 7; j++)
85c4406e 169 {
bfdcf7fb 170 fhMCLambda0DispEta [j][i] = 0;
171 fhMCLambda0DispPhi [j][i] = 0;
85c4406e 172 fhMCDispEtaDispPhi [j][i] = 0;
173 fhMCAsymmetryLambda0 [j][i] = 0;
174 fhMCAsymmetryDispEta [j][i] = 0;
bfdcf7fb 175 fhMCAsymmetryDispPhi [j][i] = 0;
176 }
34c16486 177 }
178
d2655d46 179 for(Int_t j = 0; j < 7; j++)
85c4406e 180 {
bfdcf7fb 181 fhLambda0DispEta [j] = 0;
182 fhLambda0DispPhi [j] = 0;
85c4406e 183 fhDispEtaDispPhi [j] = 0;
184 fhAsymmetryLambda0 [j] = 0;
185 fhAsymmetryDispEta [j] = 0;
bfdcf7fb 186 fhAsymmetryDispPhi [j] = 0;
5e5e056f 187
126b8c62 188 fhPtPileUp [j] = 0;
5e5e056f 189 }
bfdcf7fb 190
34c16486 191 for(Int_t i = 0; i < 3; i++)
192 {
193 fhELambda0LocMax [i] = 0;
194 fhELambda1LocMax [i] = 0;
85c4406e 195 fhEDispersionLocMax [i] = 0;
196 fhEDispEtaLocMax [i] = 0;
197 fhEDispPhiLocMax [i] = 0;
34c16486 198 fhESumEtaPhiLocMax [i] = 0;
199 fhEDispEtaPhiDiffLocMax[i] = 0;
200 fhESphericityLocMax [i] = 0;
bfdcf7fb 201 fhEAsymmetryLocMax [i] = 0;
521636d2 202 }
203
78a28af3 204 //Weight studies
1a72f6c5 205 for(Int_t i =0; i < 14; i++){
78a28af3 206 fhLambda0ForW0[i] = 0;
1a72f6c5 207 //fhLambda1ForW0[i] = 0;
3c1d9afb 208 if(i<8)fhMassPairLocMax[i] = 0;
78a28af3 209 }
210
afb3af8a 211 for(Int_t i = 0; i < 11; i++)
c2a62a94 212 {
370169ad 213 fhEtaPhiTriggerEMCALBC [i] = 0 ;
214 fhTimeTriggerEMCALBC [i] = 0 ;
215 fhTimeTriggerEMCALBCPileUpSPD[i] = 0 ;
afb3af8a 216
217 fhEtaPhiTriggerEMCALBCUM [i] = 0 ;
218 fhTimeTriggerEMCALBCUM [i] = 0 ;
85c4406e 219
c2a62a94 220 }
221
477d6cee 222 //Initialize parameters
223 InitParameters();
224
225}
477d6cee 226
126b8c62 227//_______________________________________________________________________________________________
228void AliAnaPi0EbE::FillPileUpHistograms(const Float_t pt, const Float_t time, AliVCluster * calo)
2ad19c3d 229{
230 // Fill some histograms to understand pile-up
231 if(!fFillPileUpHistograms) return;
232
233 //printf("E %f, time %f\n",energy,time);
234 AliVEvent * event = GetReader()->GetInputEvent();
235
126b8c62 236 fhTimePtNoCut->Fill(pt,time);
237 if(GetReader()->IsPileUpFromSPD())
2ad19c3d 238
126b8c62 239 if(GetReader()->IsPileUpFromSPD()) { fhPtPileUp[0]->Fill(pt); fhTimePtSPD ->Fill(pt,time); }
240 if(GetReader()->IsPileUpFromEMCal()) fhPtPileUp[1]->Fill(pt);
241 if(GetReader()->IsPileUpFromSPDOrEMCal()) fhPtPileUp[2]->Fill(pt);
242 if(GetReader()->IsPileUpFromSPDAndEMCal()) fhPtPileUp[3]->Fill(pt);
243 if(GetReader()->IsPileUpFromSPDAndNotEMCal()) fhPtPileUp[4]->Fill(pt);
244 if(GetReader()->IsPileUpFromEMCalAndNotSPD()) fhPtPileUp[5]->Fill(pt);
245 if(GetReader()->IsPileUpFromNotSPDAndNotEMCal()) fhPtPileUp[6]->Fill(pt);
246
247 if(event->IsPileupFromSPDInMultBins()) fhTimePtSPDMulti->Fill(pt,time);
248
249 // cells in cluster
250
251 AliVCaloCells* cells = 0;
252 if(fCalorimeter == "EMCAL") cells = GetEMCALCells();
253 else cells = GetPHOSCells();
254
255 Float_t maxCellFraction = 0.;
256 Int_t absIdMax = GetCaloUtils()->GetMaxEnergyCell(cells,calo,maxCellFraction);
257
258 Double_t tmax = cells->GetCellTime(absIdMax);
259 GetCaloUtils()->RecalibrateCellTime(tmax, fCalorimeter, absIdMax,GetReader()->GetInputEvent()->GetBunchCrossNumber());
260 tmax*=1.e9;
261
262 //Loop on cells inside cluster, max cell must be over 100 MeV and time in BC=0
263 if(cells->GetCellAmplitude(absIdMax) > 0.1 && TMath::Abs(tmax) < 30)
264 {
265 for (Int_t ipos = 0; ipos < calo->GetNCells(); ipos++)
266 {
267 Int_t absId = calo->GetCellsAbsId()[ipos];
268
269 if( absId == absIdMax ) continue ;
270
271 Double_t timecell = cells->GetCellTime(absId);
272 Float_t amp = cells->GetCellAmplitude(absId);
273 Int_t bc = GetReader()->GetInputEvent()->GetBunchCrossNumber();
274 GetCaloUtils()->GetEMCALRecoUtils()->AcceptCalibrateCell(absId,bc,amp,timecell,cells);
275 timecell*=1e9;
276
277 Float_t diff = (tmax-timecell);
278
36769d30 279 if( cells->GetCellAmplitude(absIdMax) < 0.1 ) continue ;
126b8c62 280
281 if(GetReader()->IsPileUpFromSPD())
282 {
283 fhPtCellTimePileUp[0]->Fill(pt, timecell);
284 fhPtTimeDiffPileUp[0]->Fill(pt, diff);
285 }
286
287 if(GetReader()->IsPileUpFromEMCal())
288 {
289 fhPtCellTimePileUp[1]->Fill(pt, timecell);
290 fhPtTimeDiffPileUp[1]->Fill(pt, diff);
291 }
292
293 if(GetReader()->IsPileUpFromSPDOrEMCal())
294 {
295 fhPtCellTimePileUp[2]->Fill(pt, timecell);
296 fhPtTimeDiffPileUp[2]->Fill(pt, diff);
297 }
298
299 if(GetReader()->IsPileUpFromSPDAndEMCal())
300 {
301 fhPtCellTimePileUp[3]->Fill(pt, timecell);
302 fhPtTimeDiffPileUp[3]->Fill(pt, diff);
303 }
304
305 if(GetReader()->IsPileUpFromSPDAndNotEMCal())
306 {
307 fhPtCellTimePileUp[4]->Fill(pt, timecell);
308 fhPtTimeDiffPileUp[4]->Fill(pt, diff);
309 }
310
311 if(GetReader()->IsPileUpFromEMCalAndNotSPD())
312 {
313 fhPtCellTimePileUp[5]->Fill(pt, timecell);
314 fhPtTimeDiffPileUp[5]->Fill(pt, diff);
315 }
316
317 if(GetReader()->IsPileUpFromNotSPDAndNotEMCal())
318 {
319 fhPtCellTimePileUp[6]->Fill(pt, timecell);
320 fhPtTimeDiffPileUp[6]->Fill(pt, diff);
321 }
322 }//loop
323 }
324
325 if(pt < 8) return; // Fill time figures for high energy clusters not too close to trigger threshold
2ad19c3d 326
327 AliESDEvent* esdEv = dynamic_cast<AliESDEvent*> (event);
328 AliAODEvent* aodEv = dynamic_cast<AliAODEvent*> (event);
329
330 // N pile up vertices
0f7e7205 331 Int_t nVtxSPD = -1;
332 Int_t nVtxTrk = -1;
2ad19c3d 333
334 if (esdEv)
335 {
0f7e7205 336 nVtxSPD = esdEv->GetNumberOfPileupVerticesSPD();
337 nVtxTrk = esdEv->GetNumberOfPileupVerticesTracks();
2ad19c3d 338
339 }//ESD
340 else if (aodEv)
341 {
0f7e7205 342 nVtxSPD = aodEv->GetNumberOfPileupVerticesSPD();
343 nVtxTrk = aodEv->GetNumberOfPileupVerticesTracks();
2ad19c3d 344 }//AOD
345
0f7e7205 346 fhTimeNPileUpVertSPD ->Fill(time,nVtxSPD);
347 fhTimeNPileUpVertTrack->Fill(time,nVtxTrk);
2ad19c3d 348
85c4406e 349 fhPtNPileUpSPDVtx->Fill(pt,nVtxSPD);
0f7e7205 350 fhPtNPileUpTrkVtx->Fill(pt,nVtxTrk);
351
352 if(TMath::Abs(time) < 25)
85c4406e 353 {
354 fhPtNPileUpSPDVtxTimeCut ->Fill(pt,nVtxSPD);
355 fhPtNPileUpTrkVtxTimeCut ->Fill(pt,nVtxTrk);
356 }
357
358 if(time < 75 && time > -25)
359 {
360 fhPtNPileUpSPDVtxTimeCut2->Fill(pt,nVtxSPD);
361 fhPtNPileUpTrkVtxTimeCut2->Fill(pt,nVtxTrk);
362 }
363
364 //printf("Is SPD %d, Is SPD Multi %d, n spd %d, n track %d\n",
0f7e7205 365 // GetReader()->IsPileUpFromSPD(),event->IsPileupFromSPDInMultBins(),nVtxSPD,nVtxTracks);
2ad19c3d 366
367 Int_t ncont = -1;
5559f30a 368 Float_t z1 = -1, z2 = -1;
2ad19c3d 369 Float_t diamZ = -1;
0f7e7205 370 for(Int_t iVert=0; iVert<nVtxSPD;iVert++)
2ad19c3d 371 {
372 if (esdEv)
373 {
374 const AliESDVertex* pv=esdEv->GetPileupVertexSPD(iVert);
375 ncont=pv->GetNContributors();
376 z1 = esdEv->GetPrimaryVertexSPD()->GetZ();
377 z2 = pv->GetZ();
378 diamZ = esdEv->GetDiamondZ();
379 }//ESD
380 else if (aodEv)
381 {
382 AliAODVertex *pv=aodEv->GetVertex(iVert);
383 if(pv->GetType()!=AliAODVertex::kPileupSPD) continue;
384 ncont=pv->GetNContributors();
385 z1=aodEv->GetPrimaryVertexSPD()->GetZ();
386 z2=pv->GetZ();
387 diamZ = aodEv->GetDiamondZ();
388 }// AOD
389
390 Double_t distZ = TMath::Abs(z2-z1);
391 diamZ = TMath::Abs(z2-diamZ);
392
393 fhTimeNPileUpVertContributors ->Fill(time,ncont);
394 fhTimePileUpMainVertexZDistance->Fill(time,distZ);
395 fhTimePileUpMainVertexZDiamond ->Fill(time,diamZ);
396
126b8c62 397 }// vertex loop
2ad19c3d 398}
399
40d3ce60 400
401//___________________________________________________________________________________________
402void AliAnaPi0EbE::FillRejectedClusterHistograms(const TLorentzVector mom, const Int_t mctag)
403{
85c4406e 404 // Fill histograms that do not pass the identification (SS case only)
40d3ce60 405
406 Float_t ener = mom.E();
407 Float_t pt = mom.Pt();
408 Float_t phi = mom.Phi();
409 if(phi < 0) phi+=TMath::TwoPi();
410 Float_t eta = mom.Eta();
411
412 fhPtReject ->Fill(pt);
413 fhEReject ->Fill(ener);
414
415 fhEEtaReject ->Fill(ener,eta);
416 fhEPhiReject ->Fill(ener,phi);
417 fhEtaPhiReject ->Fill(eta,phi);
418
419 if(IsDataMC())
420 {
421 Int_t mcIndex = GetMCIndex(mctag);
422 fhMCEReject [mcIndex] ->Fill(ener);
423 fhMCPtReject [mcIndex] ->Fill(pt);
85c4406e 424 }
40d3ce60 425}
426
42d47cb7 427//_____________________________________________________________________________________
85c4406e 428void AliAnaPi0EbE::FillSelectedClusterHistograms(AliVCluster* cluster,
5c46c992 429 const Int_t nMaxima,
85c4406e 430 const Int_t tag,
bfdcf7fb 431 const Float_t asy)
5c46c992 432{
42d47cb7 433 // Fill shower shape, timing and other histograms for selected clusters from decay
434
435 Float_t e = cluster->E();
436 Float_t disp = cluster->GetDispersion()*cluster->GetDispersion();
437 Float_t l0 = cluster->GetM02();
85c4406e 438 Float_t l1 = cluster->GetM20();
42d47cb7 439 Int_t nSM = GetModuleNumber(cluster);
85c4406e 440
bfdcf7fb 441 Int_t ebin = -1;
442 if (e < 2 ) ebin = 0;
443 else if (e < 4 ) ebin = 1;
444 else if (e < 6 ) ebin = 2;
445 else if (e < 10) ebin = 3;
85c4406e 446 else if (e < 15) ebin = 4;
447 else if (e < 20) ebin = 5;
448 else ebin = 6;
449
bfdcf7fb 450 Int_t indexMax = -1;
451 if (nMaxima==1) indexMax = 0 ;
85c4406e 452 else if(nMaxima==2) indexMax = 1 ;
453 else indexMax = 2 ;
bfdcf7fb 454
455
85c4406e 456 AliVCaloCells * cell = 0x0;
457 if(fCalorimeter == "PHOS")
42d47cb7 458 cell = GetPHOSCells();
85c4406e 459 else
42d47cb7 460 cell = GetEMCALCells();
461
462 Float_t maxCellFraction = 0;
463 GetCaloUtils()->GetMaxEnergyCell(cell, cluster, maxCellFraction);
85c4406e 464 fhEFracMaxCell->Fill(e,maxCellFraction);
42d47cb7 465
466 FillWeightHistograms(cluster);
467
85c4406e 468 fhEDispersion->Fill(e, disp);
469 fhELambda0 ->Fill(e, l0 );
470 fhELambda1 ->Fill(e, l1 );
42d47cb7 471
34c16486 472 Float_t ll0 = 0., ll1 = 0.;
85c4406e 473 Float_t dispp= 0., dEta = 0., dPhi = 0.;
474 Float_t sEta = 0., sPhi = 0., sEtaPhi = 0.;
764ab1f4 475 if(fCalorimeter == "EMCAL" && !fFillOnlySimpleSSHisto)
34c16486 476 {
477 GetCaloUtils()->GetEMCALRecoUtils()->RecalculateClusterShowerShapeParameters(GetEMCALGeometry(), GetReader()->GetInputEvent()->GetEMCALCells(), cluster,
478 ll0, ll1, dispp, dEta, dPhi, sEta, sPhi, sEtaPhi);
764ab1f4 479
34c16486 480 fhDispEtaE -> Fill(e,dEta);
481 fhDispPhiE -> Fill(e,dPhi);
482 fhSumEtaE -> Fill(e,sEta);
483 fhSumPhiE -> Fill(e,sPhi);
484 fhSumEtaPhiE -> Fill(e,sEtaPhi);
485 fhDispEtaPhiDiffE -> Fill(e,dPhi-dEta);
26e228ff 486 if(dEta+dPhi>0)fhSphericityE -> Fill(e,(dPhi-dEta)/(dEta+dPhi));
34c16486 487
bfdcf7fb 488 fhDispEtaDispPhi[ebin]->Fill(dEta,dPhi);
489 fhLambda0DispEta[ebin]->Fill(l0 ,dEta);
490 fhLambda0DispPhi[ebin]->Fill(l0 ,dPhi);
34c16486 491
bfdcf7fb 492 if (fAnaType==kSSCalo)
493 {
494 // Asymmetry histograms
bfdcf7fb 495 fhAsymmetryLambda0[ebin]->Fill(l0 ,asy);
496 fhAsymmetryDispEta[ebin]->Fill(dEta,asy);
497 fhAsymmetryDispPhi[ebin]->Fill(dPhi,asy);
498 }
85c4406e 499 }
34c16486 500
6e66993c 501 fhNLocMaxE ->Fill(e ,nMaxima);
85c4406e 502
503 fhELambda0LocMax [indexMax]->Fill(e,l0);
34c16486 504 fhELambda1LocMax [indexMax]->Fill(e,l1);
505 fhEDispersionLocMax[indexMax]->Fill(e,disp);
764ab1f4 506
85c4406e 507 if(fCalorimeter=="EMCAL" && !fFillOnlySimpleSSHisto)
34c16486 508 {
509 fhEDispEtaLocMax [indexMax]-> Fill(e,dEta);
510 fhEDispPhiLocMax [indexMax]-> Fill(e,dPhi);
511 fhESumEtaPhiLocMax [indexMax]-> Fill(e,sEtaPhi);
512 fhEDispEtaPhiDiffLocMax[indexMax]-> Fill(e,dPhi-dEta);
bfdcf7fb 513 if(dEta+dPhi>0) fhESphericityLocMax[indexMax]->Fill(e,(dPhi-dEta)/(dEta+dPhi));
514 if(fAnaType==kSSCalo) fhEAsymmetryLocMax [indexMax]->Fill(e ,asy);
515
34c16486 516 }
517
85c4406e 518 if(fCalorimeter=="EMCAL" && nSM < 6)
b5dbb99b 519 {
42d47cb7 520 fhELambda0NoTRD->Fill(e, l0 );
85c4406e 521 fhEFracMaxCellNoTRD->Fill(e,maxCellFraction);
42d47cb7 522 }
523
85c4406e 524 if(maxCellFraction < 0.5)
525 fhELambda0FracMaxCellCut->Fill(e, l0 );
42d47cb7 526
527 fhETime ->Fill(e, cluster->GetTOF()*1.e9);
528 fhENCells->Fill(e, cluster->GetNCells());
529
09273901 530 // Fill Track matching control histograms
b5dbb99b 531 if(fFillTMHisto)
532 {
09273901 533 Float_t dZ = cluster->GetTrackDz();
534 Float_t dR = cluster->GetTrackDx();
85c4406e 535
b5dbb99b 536 if(cluster->IsEMCAL() && GetCaloUtils()->IsRecalculationOfClusterTrackMatchingOn())
537 {
09273901 538 dR = 2000., dZ = 2000.;
31ae6d59 539 GetCaloUtils()->GetEMCALRecoUtils()->GetMatchedResiduals(cluster->GetID(),dZ,dR);
85c4406e 540 }
09273901 541 //printf("Pi0EbE: dPhi %f, dEta %f\n",dR,dZ);
85c4406e 542
b2e375c7 543 AliVTrack *track = GetCaloUtils()->GetMatchedTrack(cluster, GetReader()->GetInputEvent());
544
545 Bool_t positive = kFALSE;
546 if(track) positive = (track->Charge()>0);
547
b5dbb99b 548 if(fhTrackMatchedDEta && TMath::Abs(dR) < 999)
549 {
09273901 550 fhTrackMatchedDEta->Fill(e,dZ);
551 fhTrackMatchedDPhi->Fill(e,dR);
85c4406e 552 if(e > 0.5) fhTrackMatchedDEtaDPhi->Fill(dZ,dR);
b2e375c7 553
554 if(track)
555 {
556 if(positive)
557 {
558 fhTrackMatchedDEtaPos->Fill(cluster->E(),dZ);
559 fhTrackMatchedDPhiPos->Fill(cluster->E(),dR);
560 if(cluster->E() > 0.5) fhTrackMatchedDEtaDPhiPos->Fill(dZ,dR);
561 }
562 else
563 {
564 fhTrackMatchedDEtaNeg->Fill(cluster->E(),dZ);
565 fhTrackMatchedDPhiNeg->Fill(cluster->E(),dR);
566 if(cluster->E() > 0.5) fhTrackMatchedDEtaDPhiNeg->Fill(dZ,dR);
567 }
568 }
09273901 569 }
31ae6d59 570 // Check dEdx and E/p of matched clusters
571
572 if(TMath::Abs(dZ) < 0.05 && TMath::Abs(dR) < 0.05)
b2e375c7 573 {
85c4406e 574 if(track)
34c16486 575 {
31ae6d59 576 Float_t dEdx = track->GetTPCsignal();
577 fhdEdx->Fill(e, dEdx);
578
579 Float_t eOverp = e/track->P();
580 fhEOverP->Fill(e, eOverp);
4bfeae64 581
b5dbb99b 582 if(fCalorimeter=="EMCAL" && nSM < 6) fhEOverPNoTRD->Fill(e, eOverp);
85c4406e 583
31ae6d59 584 }
85c4406e 585 //else
4bfeae64 586 // printf("AliAnaPi0EbE::FillSelectedClusterHistograms() - Residual OK but (dR, dZ)= (%2.4f,%2.4f) no track associated WHAT? \n", dR,dZ);
587
b5dbb99b 588 if(IsDataMC())
589 {
f7d8e6b8 590 Float_t mctag = -1;
31ae6d59 591 if ( !GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCConversion) )
592 {
593 if ( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0) ||
85c4406e 594 GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEta) ) mctag = 2.5 ;
5dde270e 595 else if ( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPhoton) ) mctag = 0.5 ;
596 else if ( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCElectron) ) mctag = 1.5 ;
597 else mctag = 3.5 ;
31ae6d59 598
599 }
600 else
601 {
602 if ( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0) ||
85c4406e 603 GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEta) ) mctag = 6.5 ;
5dde270e 604 else if ( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPhoton) ) mctag = 4.5 ;
605 else if ( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCElectron) ) mctag = 5.5 ;
606 else mctag = 7.5 ;
607 }
608
609 fhTrackMatchedMCParticleE ->Fill(e , mctag);
610 fhTrackMatchedMCParticleDEta->Fill(dZ, mctag);
611 fhTrackMatchedMCParticleDPhi->Fill(dR, mctag);
612
613 } // MC
31ae6d59 614 }
85c4406e 615 }// Track matching histograms
09273901 616
85c4406e 617 if(IsDataMC())
b5dbb99b 618 {
3455f821 619 Int_t mcIndex = GetMCIndex(tag);
34c16486 620
621 fhEMCLambda0[mcIndex] ->Fill(e, l0);
622 fhEMCLambda1[mcIndex] ->Fill(e, l1);
623 fhEMCDispersion[mcIndex] ->Fill(e, disp);
85c4406e 624 fhEMCFracMaxCell[mcIndex]->Fill(e,maxCellFraction);
34c16486 625
85c4406e 626 if(fCalorimeter=="EMCAL" && nSM < 6)
34c16486 627 fhEMCLambda0NoTRD[mcIndex]->Fill(e, l0 );
764ab1f4 628
85c4406e 629 if(maxCellFraction < 0.5)
630 fhEMCLambda0FracMaxCellCut[mcIndex]->Fill(e, l0 );
34c16486 631
764ab1f4 632 if(fCalorimeter == "EMCAL" && !fFillOnlySimpleSSHisto)
34c16486 633 {
634 fhMCEDispEta [mcIndex]-> Fill(e,dEta);
635 fhMCEDispPhi [mcIndex]-> Fill(e,dPhi);
636 fhMCESumEtaPhi [mcIndex]-> Fill(e,sEtaPhi);
637 fhMCEDispEtaPhiDiff [mcIndex]-> Fill(e,dPhi-dEta);
85c4406e 638 if(dEta+dPhi>0)fhMCESphericity[mcIndex]-> Fill(e,(dPhi-dEta)/(dEta+dPhi));
639
bfdcf7fb 640 if (fAnaType==kSSCalo)
641 {
bfdcf7fb 642 fhMCAsymmetryLambda0[ebin][mcIndex]->Fill(l0 ,asy);
643 fhMCAsymmetryDispEta[ebin][mcIndex]->Fill(dEta,asy);
644 fhMCAsymmetryDispPhi[ebin][mcIndex]->Fill(dPhi,asy);
645 }
646
647 fhMCDispEtaDispPhi[ebin][mcIndex]->Fill(dEta,dPhi);
648 fhMCLambda0DispEta[ebin][mcIndex]->Fill(l0 ,dEta);
649 fhMCLambda0DispPhi[ebin][mcIndex]->Fill(l0 ,dPhi);
34c16486 650
651 }
652
42d47cb7 653 }//MC
bfdcf7fb 654
42d47cb7 655}
656
657//________________________________________________________
658void AliAnaPi0EbE::FillWeightHistograms(AliVCluster *clus)
659{
660 // Calculate weights and fill histograms
661
662 if(!fFillWeightHistograms || GetMixedEvent()) return;
663
664 AliVCaloCells* cells = 0;
665 if(fCalorimeter == "EMCAL") cells = GetEMCALCells();
666 else cells = GetPHOSCells();
667
668 // First recalculate energy in case non linearity was applied
669 Float_t energy = 0;
85c4406e 670 Float_t ampMax = 0;
671 for (Int_t ipos = 0; ipos < clus->GetNCells(); ipos++)
b5dbb99b 672 {
42d47cb7 673
674 Int_t id = clus->GetCellsAbsId()[ipos];
675
676 //Recalibrate cell energy if needed
677 Float_t amp = cells->GetCellAmplitude(id);
dbba06ca 678 GetCaloUtils()->RecalibrateCellAmplitude(amp,fCalorimeter, id);
42d47cb7 679
680 energy += amp;
681
85c4406e 682 if(amp> ampMax)
42d47cb7 683 ampMax = amp;
684
85c4406e 685 } // energy loop
42d47cb7 686
85c4406e 687 if(energy <=0 )
b5dbb99b 688 {
42d47cb7 689 printf("AliAnaPi0EbE::WeightHistograms()- Wrong calculated energy %f\n",energy);
690 return;
691 }
692
693 fhEMaxCellClusterRatio ->Fill(energy,ampMax/energy);
694 fhEMaxCellClusterLogRatio->Fill(energy,TMath::Log(ampMax/energy));
695
696 //Get the ratio and log ratio to all cells in cluster
85c4406e 697 for (Int_t ipos = 0; ipos < clus->GetNCells(); ipos++)
b5dbb99b 698 {
42d47cb7 699 Int_t id = clus->GetCellsAbsId()[ipos];
700
701 //Recalibrate cell energy if needed
702 Float_t amp = cells->GetCellAmplitude(id);
dbba06ca 703 GetCaloUtils()->RecalibrateCellAmplitude(amp,fCalorimeter, id);
42d47cb7 704
705 fhECellClusterRatio ->Fill(energy,amp/energy);
706 fhECellClusterLogRatio->Fill(energy,TMath::Log(amp/energy));
85c4406e 707 }
42d47cb7 708
709 //Recalculate shower shape for different W0
710 if(fCalorimeter=="EMCAL"){
711
712 Float_t l0org = clus->GetM02();
713 Float_t l1org = clus->GetM20();
714 Float_t dorg = clus->GetDispersion();
715
b5dbb99b 716 for(Int_t iw = 0; iw < 14; iw++)
717 {
85c4406e 718 GetCaloUtils()->GetEMCALRecoUtils()->SetW0(1+iw*0.5);
42d47cb7 719 GetCaloUtils()->GetEMCALRecoUtils()->RecalculateClusterShowerShapeParameters(GetEMCALGeometry(), cells, clus);
720
721 fhLambda0ForW0[iw]->Fill(energy,clus->GetM02());
1a72f6c5 722 //fhLambda1ForW0[iw]->Fill(energy,clus->GetM20());
42d47cb7 723
724 } // w0 loop
725
726 // Set the original values back
727 clus->SetM02(l0org);
728 clus->SetM20(l1org);
729 clus->SetDispersion(dorg);
730
731 }// EMCAL
732}
733
b5dbb99b 734//__________________________________________
735TObjString * AliAnaPi0EbE::GetAnalysisCuts()
85c4406e 736{
0c1383b5 737 //Save parameters used for analysis
521636d2 738 TString parList ; //this will be list of parameters used for this analysis.
739 const Int_t buffersize = 255;
740 char onePar[buffersize] ;
741
742 snprintf(onePar,buffersize,"--- AliAnaPi0EbE ---\n") ;
85c4406e 743 parList+=onePar ;
521636d2 744 snprintf(onePar,buffersize,"fAnaType=%d (Pi0 selection type) \n",fAnaType) ;
745 parList+=onePar ;
746
b5dbb99b 747 if(fAnaType == kSSCalo)
748 {
521636d2 749 snprintf(onePar,buffersize,"Calorimeter: %s\n",fCalorimeter.Data()) ;
750 parList+=onePar ;
751 snprintf(onePar,buffersize,"fMinDist =%2.2f (Minimal distance to bad channel to accept cluster) \n",fMinDist) ;
752 parList+=onePar ;
753 snprintf(onePar,buffersize,"fMinDist2=%2.2f (Cuts on Minimal distance to study acceptance evaluation) \n",fMinDist2) ;
754 parList+=onePar ;
755 snprintf(onePar,buffersize,"fMinDist3=%2.2f (One more cut on distance used for acceptance-efficiency study) \n",fMinDist3) ;
756 parList+=onePar ;
757 }
758
759 //Get parameters set in base class.
760 parList += GetBaseParametersList() ;
761
762 //Get parameters set in PID class.
763 if(fAnaType == kSSCalo) parList += GetCaloPID()->GetPIDParametersList() ;
764
765 return new TObjString(parList) ;
0c1383b5 766}
767
78a28af3 768//_____________________________________________
477d6cee 769TList * AliAnaPi0EbE::GetCreateOutputObjects()
85c4406e 770{
771 // Create histograms to be saved in output file and
477d6cee 772 // store them in outputContainer
85c4406e 773 TList * outputContainer = new TList() ;
774 outputContainer->SetName("Pi0EbEHistos") ;
477d6cee 775
745913ae 776 Int_t nptbins = GetHistogramRanges()->GetHistoPtBins(); Float_t ptmax = GetHistogramRanges()->GetHistoPtMax(); Float_t ptmin = GetHistogramRanges()->GetHistoPtMin();
777 Int_t nphibins = GetHistogramRanges()->GetHistoPhiBins(); Float_t phimax = GetHistogramRanges()->GetHistoPhiMax(); Float_t phimin = GetHistogramRanges()->GetHistoPhiMin();
778 Int_t netabins = GetHistogramRanges()->GetHistoEtaBins(); Float_t etamax = GetHistogramRanges()->GetHistoEtaMax(); Float_t etamin = GetHistogramRanges()->GetHistoEtaMin();
779 Int_t ssbins = GetHistogramRanges()->GetHistoShowerShapeBins(); Float_t ssmax = GetHistogramRanges()->GetHistoShowerShapeMax(); Float_t ssmin = GetHistogramRanges()->GetHistoShowerShapeMin();
780 Int_t tdbins = GetHistogramRanges()->GetHistoDiffTimeBins() ; Float_t tdmax = GetHistogramRanges()->GetHistoDiffTimeMax(); Float_t tdmin = GetHistogramRanges()->GetHistoDiffTimeMin();
781 Int_t tbins = GetHistogramRanges()->GetHistoTimeBins() ; Float_t tmax = GetHistogramRanges()->GetHistoTimeMax(); Float_t tmin = GetHistogramRanges()->GetHistoTimeMin();
85c4406e 782 Int_t nbins = GetHistogramRanges()->GetHistoNClusterCellBins(); Int_t nmax = GetHistogramRanges()->GetHistoNClusterCellMax(); Int_t nmin = GetHistogramRanges()->GetHistoNClusterCellMin();
783
784 Int_t nmassbins = GetHistogramRanges()->GetHistoMassBins();
785 Float_t massmin = GetHistogramRanges()->GetHistoMassMin();
b5dbb99b 786 Float_t massmax = GetHistogramRanges()->GetHistoMassMax();
787
85c4406e 788 Int_t nresetabins = GetHistogramRanges()->GetHistoTrackResidualEtaBins();
789 Float_t resetamax = GetHistogramRanges()->GetHistoTrackResidualEtaMax();
09273901 790 Float_t resetamin = GetHistogramRanges()->GetHistoTrackResidualEtaMin();
85c4406e 791 Int_t nresphibins = GetHistogramRanges()->GetHistoTrackResidualPhiBins();
792 Float_t resphimax = GetHistogramRanges()->GetHistoTrackResidualPhiMax();
09273901 793 Float_t resphimin = GetHistogramRanges()->GetHistoTrackResidualPhiMin();
794
85c4406e 795 Int_t ndedxbins = GetHistogramRanges()->GetHistodEdxBins();
796 Float_t dedxmax = GetHistogramRanges()->GetHistodEdxMax();
31ae6d59 797 Float_t dedxmin = GetHistogramRanges()->GetHistodEdxMin();
85c4406e 798 Int_t nPoverEbins = GetHistogramRanges()->GetHistoPOverEBins();
799 Float_t pOverEmax = GetHistogramRanges()->GetHistoPOverEMax();
31ae6d59 800 Float_t pOverEmin = GetHistogramRanges()->GetHistoPOverEMin();
801
85c4406e 802 Int_t ntimebins= GetHistogramRanges()->GetHistoTimeBins();
803 Float_t timemax = GetHistogramRanges()->GetHistoTimeMax();
804 Float_t timemin = GetHistogramRanges()->GetHistoTimeMin();
2ad19c3d 805
bfdcf7fb 806 TString nlm[] ={"1 Local Maxima","2 Local Maxima", "NLM > 2"};
85c4406e 807 TString ptype[] ={"#gamma","#gamma->e^{#pm}","#pi^{0}","#eta","e^{#pm}", "hadron"};
808 TString pname[] ={"Photon","Conversion", "Pi0", "Eta", "Electron","Hadron"};
d2655d46 809 Int_t bin[] = {0,2,4,6,10,15,20,100}; // energy bins
85c4406e 810
811 fhPt = new TH1F("hPt","Number of identified #pi^{0} (#eta) decay",nptbins,ptmin,ptmax);
09273901 812 fhPt->SetYTitle("N");
9fb80477 813 fhPt->SetXTitle("p_{T} (GeV/c)");
85c4406e 814 outputContainer->Add(fhPt) ;
09273901 815
85c4406e 816 fhE = new TH1F("hE","Number of identified #pi^{0} (#eta) decay pairs",nptbins,ptmin,ptmax);
09273901 817 fhE->SetYTitle("N");
b9947879 818 fhE->SetXTitle("E (GeV)");
85c4406e 819 outputContainer->Add(fhE) ;
09273901 820
821 fhEPhi = new TH2F
85c4406e 822 ("hEPhi","Selected #pi^{0} (#eta) pairs: E vs #phi",nptbins,ptmin,ptmax, nphibins,phimin,phimax);
b9947879 823 fhEPhi->SetYTitle("#phi (rad)");
824 fhEPhi->SetXTitle("E (GeV)");
85c4406e 825 outputContainer->Add(fhEPhi) ;
09273901 826
827 fhEEta = new TH2F
85c4406e 828 ("hEEta","Selected #pi^{0} (#eta) pairs: E vs #eta",nptbins,ptmin,ptmax,netabins,etamin,etamax);
b9947879 829 fhEEta->SetYTitle("#eta");
9fb80477 830 fhEEta->SetXTitle("E (GeV)");
85c4406e 831 outputContainer->Add(fhEEta) ;
832
29250849 833 fhPtPhi = new TH2F
834 ("hPtPhi","Selected #pi^{0} (#eta) pairs: p_{T} vs #phi",nptbins,ptmin,ptmax, nphibins,phimin,phimax);
835 fhPtPhi->SetYTitle("#phi (rad)");
836 fhPtPhi->SetXTitle("p_{T} (GeV/c)");
837 outputContainer->Add(fhPtPhi) ;
838
839 fhPtEta = new TH2F
840 ("hPtEta","Selected #pi^{0} (#eta) pairs: p_{T} vs #eta",nptbins,ptmin,ptmax,netabins,etamin,etamax);
841 fhPtEta->SetYTitle("#eta");
842 fhPtEta->SetXTitle("p_{T} (GeV/c)");
843 outputContainer->Add(fhPtEta) ;
09273901 844
845 fhEtaPhi = new TH2F
85c4406e 846 ("hEtaPhi","Selected #pi^{0} (#eta) pairs: #eta vs #phi",netabins,etamin,etamax, nphibins,phimin,phimax);
b9947879 847 fhEtaPhi->SetYTitle("#phi (rad)");
848 fhEtaPhi->SetXTitle("#eta");
85c4406e 849 outputContainer->Add(fhEtaPhi) ;
09273901 850
c2a62a94 851 if(fCalorimeter=="EMCAL" && fFillEMCALBCHistograms)
852 {
853 fhEtaPhiEMCALBC0 = new TH2F
854 ("hEtaPhiEMCALBC0","cluster,E > 2 GeV, #eta vs #phi, for clusters with |time| < 25 ns, EMCAL-BC=0",netabins,etamin,etamax,nphibins,phimin,phimax);
855 fhEtaPhiEMCALBC0->SetYTitle("#phi (rad)");
856 fhEtaPhiEMCALBC0->SetXTitle("#eta");
857 outputContainer->Add(fhEtaPhiEMCALBC0) ;
858
859 fhEtaPhiEMCALBC1 = new TH2F
860 ("hEtaPhiEMCALBC1","cluster,E > 2 GeV, #eta vs #phi, for clusters with 25 < |time| < 75 ns, EMCAL-BC=1",netabins,etamin,etamax,nphibins,phimin,phimax);
861 fhEtaPhiEMCALBC1->SetYTitle("#phi (rad)");
862 fhEtaPhiEMCALBC1->SetXTitle("#eta");
863 outputContainer->Add(fhEtaPhiEMCALBC1) ;
864
865 fhEtaPhiEMCALBCN = new TH2F
866 ("hEtaPhiEMCALBCN","cluster,E > 2 GeV, #eta vs #phi, for clusters with |time| > 75 ns, EMCAL-BC>1",netabins,etamin,etamax,nphibins,phimin,phimax);
867 fhEtaPhiEMCALBCN->SetYTitle("#phi (rad)");
868 fhEtaPhiEMCALBCN->SetXTitle("#eta");
869 outputContainer->Add(fhEtaPhiEMCALBCN) ;
870
afb3af8a 871 for(Int_t i = 0; i < 11; i++)
c2a62a94 872 {
873 fhEtaPhiTriggerEMCALBC[i] = new TH2F
874 (Form("hEtaPhiTriggerEMCALBC%d",i-5),
afb3af8a 875 Form("meson E > 2 GeV, #eta vs #phi, Trigger EMCAL-BC=%d",i-5),
c2a62a94 876 netabins,etamin,etamax,nphibins,phimin,phimax);
877 fhEtaPhiTriggerEMCALBC[i]->SetYTitle("#phi (rad)");
878 fhEtaPhiTriggerEMCALBC[i]->SetXTitle("#eta");
879 outputContainer->Add(fhEtaPhiTriggerEMCALBC[i]) ;
880
881 fhTimeTriggerEMCALBC[i] = new TH2F
882 (Form("hTimeTriggerEMCALBC%d",i-5),
afb3af8a 883 Form("meson time vs E, Trigger EMCAL-BC=%d",i-5),
c2a62a94 884 nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
885 fhTimeTriggerEMCALBC[i]->SetXTitle("E (GeV)");
886 fhTimeTriggerEMCALBC[i]->SetYTitle("time (ns)");
887 outputContainer->Add(fhTimeTriggerEMCALBC[i]);
888
889 fhTimeTriggerEMCALBCPileUpSPD[i] = new TH2F
890 (Form("hTimeTriggerEMCALBC%dPileUpSPD",i-5),
afb3af8a 891 Form("meson time vs E, Trigger EMCAL-BC=%d",i-5),
c2a62a94 892 nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
893 fhTimeTriggerEMCALBCPileUpSPD[i]->SetXTitle("E (GeV)");
894 fhTimeTriggerEMCALBCPileUpSPD[i]->SetYTitle("time (ns)");
895 outputContainer->Add(fhTimeTriggerEMCALBCPileUpSPD[i]);
afb3af8a 896
897 fhEtaPhiTriggerEMCALBCUM[i] = new TH2F
898 (Form("hEtaPhiTriggerEMCALBC%d_UnMatch",i-5),
899 Form("meson E > 2 GeV, #eta vs #phi, unmatched trigger EMCAL-BC=%d",i-5),
900 netabins,etamin,etamax,nphibins,phimin,phimax);
901 fhEtaPhiTriggerEMCALBCUM[i]->SetYTitle("#phi (rad)");
902 fhEtaPhiTriggerEMCALBCUM[i]->SetXTitle("#eta");
903 outputContainer->Add(fhEtaPhiTriggerEMCALBCUM[i]) ;
904
905 fhTimeTriggerEMCALBCUM[i] = new TH2F
906 (Form("hTimeTriggerEMCALBC%d_UnMatch",i-5),
907 Form("meson time vs E, unmatched trigger EMCAL-BC=%d",i-5),
908 nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
909 fhTimeTriggerEMCALBCUM[i]->SetXTitle("E (GeV)");
910 fhTimeTriggerEMCALBCUM[i]->SetYTitle("time (ns)");
911 outputContainer->Add(fhTimeTriggerEMCALBCUM[i]);
85c4406e 912
c2a62a94 913 }
126b8c62 914
915 fhTimeTriggerEMCALBC0UMReMatchOpenTime = new TH2F("hTimeTriggerBC0_UnMatch_ReMatch_OpenTime",
916 "cluster time vs E of clusters, no match, rematch open time",
917 nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
918 fhTimeTriggerEMCALBC0UMReMatchOpenTime->SetXTitle("E (GeV)");
919 fhTimeTriggerEMCALBC0UMReMatchOpenTime->SetYTitle("time (ns)");
920 outputContainer->Add(fhTimeTriggerEMCALBC0UMReMatchOpenTime);
921
922
923 fhTimeTriggerEMCALBC0UMReMatchCheckNeigh = new TH2F("hTimeTriggerBC0_UnMatch_ReMatch_CheckNeighbours",
924 "cluster time vs E of clusters, no match, rematch with neigbour parches",
925 nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
926 fhTimeTriggerEMCALBC0UMReMatchCheckNeigh->SetXTitle("E (GeV)");
927 fhTimeTriggerEMCALBC0UMReMatchCheckNeigh->SetYTitle("time (ns)");
928 outputContainer->Add(fhTimeTriggerEMCALBC0UMReMatchCheckNeigh);
929
930 fhTimeTriggerEMCALBC0UMReMatchBoth = new TH2F("hTimeTriggerBC0_UnMatch_ReMatch_Both",
931 "cluster time vs E of clusters, no match, rematch open time and neigbour",
932 nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
933 fhTimeTriggerEMCALBC0UMReMatchBoth->SetXTitle("E (GeV)");
934 fhTimeTriggerEMCALBC0UMReMatchBoth->SetYTitle("time (ns)");
935 outputContainer->Add(fhTimeTriggerEMCALBC0UMReMatchBoth);
936
c2a62a94 937 }
938
c8710850 939 fhPtCentrality = new TH2F("hPtCentrality","centrality vs p_{T}",nptbins,ptmin,ptmax, 100,0,100);
940 fhPtCentrality->SetYTitle("centrality");
941 fhPtCentrality->SetXTitle("p_{T}(GeV/c)");
942 outputContainer->Add(fhPtCentrality) ;
943
944 fhPtEventPlane = new TH2F("hPtEventPlane","event plane angle vs p_{T}",nptbins,ptmin,ptmax, 100,0,TMath::Pi());
945 fhPtEventPlane->SetYTitle("Event plane angle (rad)");
946 fhPtEventPlane->SetXTitle("p_{T} (GeV/c)");
947 outputContainer->Add(fhPtEventPlane) ;
948
40d3ce60 949 if(fAnaType == kSSCalo)
950 {
85c4406e 951 fhPtReject = new TH1F("hPtReject","Number of rejected as #pi^{0} (#eta) decay",nptbins,ptmin,ptmax);
40d3ce60 952 fhPtReject->SetYTitle("N");
953 fhPtReject->SetXTitle("p_{T} (GeV/c)");
85c4406e 954 outputContainer->Add(fhPtReject) ;
40d3ce60 955
85c4406e 956 fhEReject = new TH1F("hEReject","Number of rejected as #pi^{0} (#eta) decay pairs",nptbins,ptmin,ptmax);
40d3ce60 957 fhEReject->SetYTitle("N");
958 fhEReject->SetXTitle("E (GeV)");
85c4406e 959 outputContainer->Add(fhEReject) ;
40d3ce60 960
961 fhEPhiReject = new TH2F
85c4406e 962 ("hEPhiReject","Rejected #pi^{0} (#eta) cluster: E vs #phi",nptbins,ptmin,ptmax, nphibins,phimin,phimax);
40d3ce60 963 fhEPhiReject->SetYTitle("#phi (rad)");
964 fhEPhiReject->SetXTitle("E (GeV)");
85c4406e 965 outputContainer->Add(fhEPhiReject) ;
40d3ce60 966
967 fhEEtaReject = new TH2F
85c4406e 968 ("hEEtaReject","Rejected #pi^{0} (#eta) cluster: E vs #eta",nptbins,ptmin,ptmax,netabins,etamin,etamax);
40d3ce60 969 fhEEtaReject->SetYTitle("#eta");
970 fhEEtaReject->SetXTitle("E (GeV)");
85c4406e 971 outputContainer->Add(fhEEtaReject) ;
40d3ce60 972
973 fhEtaPhiReject = new TH2F
85c4406e 974 ("hEtaPhiReject","Rejected #pi^{0} (#eta) cluster: #eta vs #phi",netabins,etamin,etamax, nphibins,phimin,phimax);
40d3ce60 975 fhEtaPhiReject->SetYTitle("#phi (rad)");
976 fhEtaPhiReject->SetXTitle("#eta");
85c4406e 977 outputContainer->Add(fhEtaPhiReject) ;
40d3ce60 978 }
979
f02db2c0 980 fhMass = new TH2F
85c4406e 981 ("hMass","all pairs mass: E vs mass",nptbins,ptmin,ptmax, nmassbins,massmin,massmax);
f02db2c0 982 fhMass->SetYTitle("mass (GeV/c^{2})");
983 fhMass->SetXTitle("E (GeV)");
85c4406e 984 outputContainer->Add(fhMass) ;
f02db2c0 985
986 fhSelectedMass = new TH2F
85c4406e 987 ("hSelectedMass","Selected #pi^{0} (#eta) pairs mass: E vs mass",nptbins,ptmin,ptmax, nmassbins,massmin,massmax);
f02db2c0 988 fhSelectedMass->SetYTitle("mass (GeV/c^{2})");
989 fhSelectedMass->SetXTitle("E (GeV)");
85c4406e 990 outputContainer->Add(fhSelectedMass) ;
991
29250849 992 fhMassPt = new TH2F
993 ("hMassPt","all pairs mass: p_{T} vs mass",nptbins,ptmin,ptmax, nmassbins,massmin,massmax);
994 fhMassPt->SetYTitle("mass (GeV/c^{2})");
995 fhMassPt->SetXTitle("p_{T} (GeV/c)");
996 outputContainer->Add(fhMassPt) ;
997
998 fhSelectedMassPt = new TH2F
999 ("hSelectedMassPt","Selected #pi^{0} (#eta) pairs mass: p_{T} vs mass",nptbins,ptmin,ptmax, nmassbins,massmin,massmax);
1000 fhSelectedMassPt->SetYTitle("mass (GeV/c^{2})");
1001 fhSelectedMassPt->SetXTitle("p_{T} (GeV/c)");
1002 outputContainer->Add(fhSelectedMassPt) ;
1253480f 1003
1004 if(IsDataMC() && fAnaType == kSSCalo)
1005 {
1006 fhMassNoOverlap = new TH2F
1007 ("hMassNoOverlap","all pairs mass: E vs mass, no overlap",nptbins,ptmin,ptmax, nmassbins,massmin,massmax);
1008 fhMassNoOverlap->SetYTitle("mass (GeV/c^{2})");
1009 fhMassNoOverlap->SetXTitle("E (GeV)");
1010 outputContainer->Add(fhMassNoOverlap) ;
1011
1012 fhSelectedMassNoOverlap = new TH2F
1013 ("hSelectedMassNoOverlap","Selected #pi^{0} (#eta) pairs mass: E vs mass, no overlap",nptbins,ptmin,ptmax, nmassbins,massmin,massmax);
1014 fhSelectedMassNoOverlap->SetYTitle("mass (GeV/c^{2})");
1015 fhSelectedMassNoOverlap->SetXTitle("E (GeV)");
1016 outputContainer->Add(fhSelectedMassNoOverlap) ;
1017
1018 fhMassPtNoOverlap = new TH2F
1019 ("hMassPtNoOverlap","all pairs mass: p_{T} vs mass, no overlap",nptbins,ptmin,ptmax, nmassbins,massmin,massmax);
1020 fhMassPtNoOverlap->SetYTitle("mass (GeV/c^{2})");
1021 fhMassPtNoOverlap->SetXTitle("p_{T} (GeV/c)");
1022 outputContainer->Add(fhMassPtNoOverlap) ;
1023
1024 fhSelectedMassPtNoOverlap = new TH2F
1025 ("hSelectedMassPtNoOverlap","Selected #pi^{0} (#eta) pairs mass: p_{T} vs mass, no overlap",nptbins,ptmin,ptmax, nmassbins,massmin,massmax);
1026 fhSelectedMassPtNoOverlap->SetYTitle("mass (GeV/c^{2})");
1027 fhSelectedMassPtNoOverlap->SetXTitle("p_{T} (GeV/c)");
1028 outputContainer->Add(fhSelectedMassPtNoOverlap) ;
1029 }
85c4406e 1030
34c16486 1031 if(fAnaType != kSSCalo)
1032 {
85c4406e 1033 fhPtDecay = new TH1F("hPtDecay","Number of identified #pi^{0} (#eta) decay photons",nptbins,ptmin,ptmax);
34c16486 1034 fhPtDecay->SetYTitle("N");
1035 fhPtDecay->SetXTitle("p_{T} (GeV/c)");
85c4406e 1036 outputContainer->Add(fhPtDecay) ;
34c16486 1037
85c4406e 1038 fhEDecay = new TH1F("hEDecay","Number of identified #pi^{0} (#eta) decay photons",nptbins,ptmin,ptmax);
34c16486 1039 fhEDecay->SetYTitle("N");
1040 fhEDecay->SetXTitle("E (GeV)");
85c4406e 1041 outputContainer->Add(fhEDecay) ;
34c16486 1042 }
57b97dc6 1043
c4a7d28a 1044 ////////
57b97dc6 1045
34c16486 1046 if( fFillSelectClHisto )
b5dbb99b 1047 {
c4a7d28a 1048
521636d2 1049 fhEDispersion = new TH2F
85c4406e 1050 ("hEDispersion","Selected #pi^{0} (#eta) pairs: E vs dispersion",nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
521636d2 1051 fhEDispersion->SetYTitle("D^{2}");
1052 fhEDispersion->SetXTitle("E (GeV)");
85c4406e 1053 outputContainer->Add(fhEDispersion) ;
521636d2 1054
1055 fhELambda0 = new TH2F
85c4406e 1056 ("hELambda0","Selected #pi^{0} (#eta) pairs: E vs #lambda_{0}",nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
521636d2 1057 fhELambda0->SetYTitle("#lambda_{0}^{2}");
1058 fhELambda0->SetXTitle("E (GeV)");
85c4406e 1059 outputContainer->Add(fhELambda0) ;
1060
42d47cb7 1061 fhELambda1 = new TH2F
85c4406e 1062 ("hELambda1","Selected #pi^{0} (#eta) pairs: E vs #lambda_{1}",nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
42d47cb7 1063 fhELambda1->SetYTitle("#lambda_{1}^{2}");
1064 fhELambda1->SetXTitle("E (GeV)");
85c4406e 1065 outputContainer->Add(fhELambda1) ;
1066
3bfcb597 1067 fhELambda0FracMaxCellCut = new TH2F
85c4406e 1068 ("hELambda0FracMaxCellCut","Selected #pi^{0} (#eta) pairs: E vs #lambda_{0}, Max cell fraction of energy < 0.5",nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
3bfcb597 1069 fhELambda0FracMaxCellCut->SetYTitle("#lambda_{0}^{2}");
1070 fhELambda0FracMaxCellCut->SetXTitle("E (GeV)");
85c4406e 1071 outputContainer->Add(fhELambda0FracMaxCellCut) ;
1072
3bfcb597 1073 fhEFracMaxCell = new TH2F
85c4406e 1074 ("hEFracMaxCell","Selected #pi^{0} (#eta) pairs: E vs #lambda_{0}, Max cell fraction of energy",nptbins,ptmin,ptmax,100,0,1);
3bfcb597 1075 fhEFracMaxCell->SetYTitle("Fraction");
1076 fhEFracMaxCell->SetXTitle("E (GeV)");
85c4406e 1077 outputContainer->Add(fhEFracMaxCell) ;
5c46c992 1078
06e81356 1079 if(fCalorimeter=="EMCAL")
1080 {
3bfcb597 1081 fhELambda0NoTRD = new TH2F
85c4406e 1082 ("hELambda0NoTRD","Selected #pi^{0} (#eta) pairs: E vs #lambda_{0}, not behind TRD",nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
3bfcb597 1083 fhELambda0NoTRD->SetYTitle("#lambda_{0}^{2}");
1084 fhELambda0NoTRD->SetXTitle("E (GeV)");
85c4406e 1085 outputContainer->Add(fhELambda0NoTRD) ;
3bfcb597 1086
1087 fhEFracMaxCellNoTRD = new TH2F
85c4406e 1088 ("hEFracMaxCellNoTRD","Selected #pi^{0} (#eta) pairs: E vs #lambda_{0}, Max cell fraction of energy, not behind TRD",nptbins,ptmin,ptmax,100,0,1);
3bfcb597 1089 fhEFracMaxCellNoTRD->SetYTitle("Fraction");
1090 fhEFracMaxCellNoTRD->SetXTitle("E (GeV)");
85c4406e 1091 outputContainer->Add(fhEFracMaxCellNoTRD) ;
34c16486 1092
764ab1f4 1093 if(!fFillOnlySimpleSSHisto)
34c16486 1094 {
85c4406e 1095 fhDispEtaE = new TH2F ("hDispEtaE","#sigma^{2}_{#eta #eta} = #Sigma w_{i}(#eta_{i} - <#eta>)^{2}/ #Sigma w_{i} vs E", nptbins,ptmin,ptmax, ssbins,ssmin,ssmax);
764ab1f4 1096 fhDispEtaE->SetXTitle("E (GeV)");
1097 fhDispEtaE->SetYTitle("#sigma^{2}_{#eta #eta}");
85c4406e 1098 outputContainer->Add(fhDispEtaE);
764ab1f4 1099
85c4406e 1100 fhDispPhiE = new TH2F ("hDispPhiE","#sigma^{2}_{#phi #phi} = #Sigma w_{i}(#phi_{i} - <#phi>)^{2} / #Sigma w_{i} vs E", nptbins,ptmin,ptmax, ssbins,ssmin,ssmax);
764ab1f4 1101 fhDispPhiE->SetXTitle("E (GeV)");
1102 fhDispPhiE->SetYTitle("#sigma^{2}_{#phi #phi}");
85c4406e 1103 outputContainer->Add(fhDispPhiE);
764ab1f4 1104
85c4406e 1105 fhSumEtaE = new TH2F ("hSumEtaE","#sigma^{2}_{#eta #eta} = #Sigma w_{i}(#eta_{i})^{2} / #Sigma w_{i} - <#eta>^{2} vs E", nptbins,ptmin,ptmax, ssbins,ssmin,ssmax);
764ab1f4 1106 fhSumEtaE->SetXTitle("E (GeV)");
1107 fhSumEtaE->SetYTitle("#delta^{2}_{#eta #eta}");
85c4406e 1108 outputContainer->Add(fhSumEtaE);
764ab1f4 1109
85c4406e 1110 fhSumPhiE = new TH2F ("hSumPhiE","#sigma^{2}_{#phi #phi} = #Sigma w_{i}(#phi_{i})^{2}/ #Sigma w_{i} - <#phi>^{2} vs E",
1111 nptbins,ptmin,ptmax, ssbins,ssmin,ssmax);
764ab1f4 1112 fhSumPhiE->SetXTitle("E (GeV)");
1113 fhSumPhiE->SetYTitle("#delta^{2}_{#phi #phi}");
85c4406e 1114 outputContainer->Add(fhSumPhiE);
764ab1f4 1115
85c4406e 1116 fhSumEtaPhiE = new TH2F ("hSumEtaPhiE","#delta^{2}_{#eta #phi} = #Sigma w_{i}(#phi_{i} #eta_{i} ) / #Sigma w_{i} - <#phi><#eta> vs E",
1117 nptbins,ptmin,ptmax, 2*ssbins,-ssmax,ssmax);
764ab1f4 1118 fhSumEtaPhiE->SetXTitle("E (GeV)");
1119 fhSumEtaPhiE->SetYTitle("#delta^{2}_{#eta #phi}");
1120 outputContainer->Add(fhSumEtaPhiE);
bfdcf7fb 1121
85c4406e 1122 fhDispEtaPhiDiffE = new TH2F ("hDispEtaPhiDiffE","#sigma^{2}_{#phi #phi} - #sigma^{2}_{#eta #eta} vs E",
1123 nptbins,ptmin,ptmax,200, -10,10);
764ab1f4 1124 fhDispEtaPhiDiffE->SetXTitle("E (GeV)");
1125 fhDispEtaPhiDiffE->SetYTitle("#sigma^{2}_{#phi #phi}-#sigma^{2}_{#eta #eta}");
85c4406e 1126 outputContainer->Add(fhDispEtaPhiDiffE);
bfdcf7fb 1127
85c4406e 1128 fhSphericityE = new TH2F ("hSphericityE","(#sigma^{2}_{#phi #phi} - #sigma^{2}_{#eta #eta}) / (#sigma^{2}_{#eta #eta} + #sigma^{2}_{#phi #phi}) vs E",
1129 nptbins,ptmin,ptmax, 200, -1,1);
764ab1f4 1130 fhSphericityE->SetXTitle("E (GeV)");
1131 fhSphericityE->SetYTitle("s = (#sigma^{2}_{#phi #phi} - #sigma^{2}_{#eta #eta}) / (#sigma^{2}_{#eta #eta} + #sigma^{2}_{#phi #phi})");
1132 outputContainer->Add(fhSphericityE);
bfdcf7fb 1133
764ab1f4 1134 for(Int_t i = 0; i < 7; i++)
1135 {
85c4406e 1136 fhDispEtaDispPhi[i] = new TH2F (Form("hDispEtaDispPhi_EBin%d",i),Form("#sigma^{2}_{#phi #phi} vs #sigma^{2}_{#eta #eta} for %d < E < %d GeV",bin[i],bin[i+1]),
1137 ssbins,ssmin,ssmax , ssbins,ssmin,ssmax);
764ab1f4 1138 fhDispEtaDispPhi[i]->SetXTitle("#sigma^{2}_{#eta #eta}");
1139 fhDispEtaDispPhi[i]->SetYTitle("#sigma^{2}_{#phi #phi}");
85c4406e 1140 outputContainer->Add(fhDispEtaDispPhi[i]);
764ab1f4 1141
85c4406e 1142 fhLambda0DispEta[i] = new TH2F (Form("hLambda0DispEta_EBin%d",i),Form("#lambda^{2}_{0} vs #sigma^{2}_{#eta #eta} for %d < E < %d GeV",bin[i],bin[i+1]),
1143 ssbins,ssmin,ssmax , ssbins,ssmin,ssmax);
764ab1f4 1144 fhLambda0DispEta[i]->SetXTitle("#lambda^{2}_{0}");
1145 fhLambda0DispEta[i]->SetYTitle("#sigma^{2}_{#eta #eta}");
85c4406e 1146 outputContainer->Add(fhLambda0DispEta[i]);
764ab1f4 1147
85c4406e 1148 fhLambda0DispPhi[i] = new TH2F (Form("hLambda0DispPhi_EBin%d",i),Form("#lambda^{2}_{0}} vs #sigma^{2}_{#phi #phi} for %d < E < %d GeV",bin[i],bin[i+1]),
1149 ssbins,ssmin,ssmax , ssbins,ssmin,ssmax);
764ab1f4 1150 fhLambda0DispPhi[i]->SetXTitle("#lambda^{2}_{0}");
1151 fhLambda0DispPhi[i]->SetYTitle("#sigma^{2}_{#phi #phi}");
85c4406e 1152 outputContainer->Add(fhLambda0DispPhi[i]);
764ab1f4 1153
1154 }
34c16486 1155 }
85c4406e 1156 }
34c16486 1157
6e66993c 1158 fhNLocMaxE = new TH2F("hNLocMaxE","Number of local maxima in cluster",
1159 nptbins,ptmin,ptmax,10,0,10);
1160 fhNLocMaxE ->SetYTitle("N maxima");
1161 fhNLocMaxE ->SetXTitle("E (GeV)");
1162 outputContainer->Add(fhNLocMaxE) ;
1163
1164 if(fAnaType == kSSCalo)
1165 {
1166 fhNLocMaxPt = new TH2F("hNLocMaxPt","Number of local maxima in cluster",
85c4406e 1167 nptbins,ptmin,ptmax,10,0,10);
6e66993c 1168 fhNLocMaxPt ->SetYTitle("N maxima");
1169 fhNLocMaxPt ->SetXTitle("p_{T} (GeV/c)");
1170 outputContainer->Add(fhNLocMaxPt) ;
1171 }
521636d2 1172
85c4406e 1173 for (Int_t i = 0; i < 3; i++)
34c16486 1174 {
1175 fhELambda0LocMax[i] = new TH2F(Form("hELambda0LocMax%d",i+1),
1176 Form("Selected #pi^{0} (#eta) pairs: E vs #lambda_{0}, %s",nlm[i].Data()),
85c4406e 1177 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
34c16486 1178 fhELambda0LocMax[i]->SetYTitle("#lambda_{0}^{2}");
1179 fhELambda0LocMax[i]->SetXTitle("E (GeV)");
85c4406e 1180 outputContainer->Add(fhELambda0LocMax[i]) ;
34c16486 1181
1182 fhELambda1LocMax[i] = new TH2F(Form("hELambda1LocMax%d",i+1),
1183 Form("Selected #pi^{0} (#eta) pairs: E vs #lambda_{1}, %s",nlm[i].Data()),
85c4406e 1184 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
34c16486 1185 fhELambda1LocMax[i]->SetYTitle("#lambda_{1}^{2}");
1186 fhELambda1LocMax[i]->SetXTitle("E (GeV)");
85c4406e 1187 outputContainer->Add(fhELambda1LocMax[i]) ;
34c16486 1188
1189 fhEDispersionLocMax[i] = new TH2F(Form("hEDispersionLocMax%d",i+1),
1190 Form("Selected #pi^{0} (#eta) pairs: E vs dispersion^{2}, %s",nlm[i].Data()),
85c4406e 1191 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
34c16486 1192 fhEDispersionLocMax[i]->SetYTitle("dispersion^{2}");
1193 fhEDispersionLocMax[i]->SetXTitle("E (GeV)");
85c4406e 1194 outputContainer->Add(fhEDispersionLocMax[i]) ;
34c16486 1195
764ab1f4 1196 if(fCalorimeter == "EMCAL" && !fFillOnlySimpleSSHisto)
34c16486 1197 {
1198 fhEDispEtaLocMax[i] = new TH2F(Form("hEDispEtaLocMax%d",i+1),
1199 Form("Selected #pi^{0} (#eta) pairs: E vs #sigma_{#eta #eta}, %s",nlm[i].Data()),
85c4406e 1200 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
34c16486 1201 fhEDispEtaLocMax[i]->SetYTitle("#sigma_{#eta #eta}");
1202 fhEDispEtaLocMax[i]->SetXTitle("E (GeV)");
85c4406e 1203 outputContainer->Add(fhEDispEtaLocMax[i]) ;
34c16486 1204
1205 fhEDispPhiLocMax[i] = new TH2F(Form("hEDispPhiLocMax%d",i+1),
1206 Form("Selected #pi^{0} (#eta) pairs: E vs #sigma_{#phi #phi}, %s",nlm[i].Data()),
85c4406e 1207 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
34c16486 1208 fhEDispPhiLocMax[i]->SetYTitle("#sigma_{#phi #phi}");
1209 fhEDispPhiLocMax[i]->SetXTitle("E (GeV)");
85c4406e 1210 outputContainer->Add(fhEDispPhiLocMax[i]) ;
34c16486 1211
1212 fhESumEtaPhiLocMax[i] = new TH2F(Form("hESumEtaPhiLocMax%d",i+1),
1213 Form("Selected #pi^{0} (#eta) pairs: E vs #sigma_{#eta #phi}, %s",nlm[i].Data()),
85c4406e 1214 nptbins,ptmin,ptmax,2*ssbins,-ssmax,ssmax);
34c16486 1215 fhESumEtaPhiLocMax[i]->SetYTitle("#sigma_{#eta #phi}");
1216 fhESumEtaPhiLocMax[i]->SetXTitle("E (GeV)");
85c4406e 1217 outputContainer->Add(fhESumEtaPhiLocMax[i]) ;
34c16486 1218
1219 fhEDispEtaPhiDiffLocMax[i] = new TH2F(Form("hEDispEtaPhiDiffLocMax%d",i+1),
1220 Form("Selected #pi^{0} (#eta) pairs: E vs #sigma_{#phi #phi} - #sigma_{#eta #eta}, %s",nlm[i].Data()),
85c4406e 1221 nptbins,ptmin,ptmax,200, -10,10);
34c16486 1222 fhEDispEtaPhiDiffLocMax[i]->SetYTitle("#sigma_{#phi #phi} - #sigma_{#eta #eta}");
1223 fhEDispEtaPhiDiffLocMax[i]->SetXTitle("E (GeV)");
85c4406e 1224 outputContainer->Add(fhEDispEtaPhiDiffLocMax[i]) ;
34c16486 1225
1226 fhESphericityLocMax[i] = new TH2F(Form("hESphericityLocMax%d",i+1),
1227 Form("Selected #pi^{0} (#eta) pairs: E vs #sigma_{#phi #phi} - #sigma_{#eta #eta} / (#sigma_{#phi #phi} + #sigma_{#eta #eta}), %s",nlm[i].Data()),
85c4406e 1228 nptbins,ptmin,ptmax,200, -1,1);
34c16486 1229 fhESphericityLocMax[i]->SetYTitle("#sigma_{#phi #phi} - #sigma_{#eta #eta} / (#sigma_{#phi #phi} + #sigma_{#eta #eta})");
1230 fhESphericityLocMax[i]->SetXTitle("E (GeV)");
1231 outputContainer->Add(fhESphericityLocMax[i]) ;
1232 }
34c16486 1233
85c4406e 1234 }
1235
1236 fhENCells = new TH2F ("hENCells","N cells in cluster vs E ", nptbins,ptmin,ptmax, nbins,nmin,nmax);
42d47cb7 1237 fhENCells->SetXTitle("E (GeV)");
1238 fhENCells->SetYTitle("# of cells in cluster");
85c4406e 1239 outputContainer->Add(fhENCells);
42d47cb7 1240
1241 fhETime = new TH2F("hETime","cluster time vs pair E",nptbins,ptmin,ptmax, tbins,tmin,tmax);
1242 fhETime->SetXTitle("E (GeV)");
9fb80477 1243 fhETime->SetYTitle("t (ns)");
85c4406e 1244 outputContainer->Add(fhETime);
521636d2 1245
764ab1f4 1246 }
e7fd282f 1247
85c4406e 1248
1249 fhEPairDiffTime = new TH2F("hEPairDiffTime","cluster pair time difference vs E",nptbins,ptmin,ptmax, tdbins,tdmin,tdmax);
1250 fhEPairDiffTime->SetXTitle("E_{pair} (GeV)");
1251 fhEPairDiffTime->SetYTitle("#Delta t (ns)");
1252 outputContainer->Add(fhEPairDiffTime);
a1fd1b69 1253
1254 if(fAnaType == kIMCalo)
1255 {
3c1d9afb 1256 TString combiName [] = {"1LocMax","2LocMax","NLocMax","1LocMax2LocMax","1LocMaxNLocMax","2LocMaxNLocMax","1LocMaxSSBad","NLocMaxSSGood"};
5c46c992 1257 TString combiTitle[] = {"1 Local Maxima in both clusters","2 Local Maxima in both clusters","more than 2 Local Maxima in both clusters",
1258 "1 Local Maxima paired with 2 Local Maxima","1 Local Maxima paired with more than 2 Local Maxima",
3c1d9afb 1259 "2 Local Maxima paired with more than 2 Local Maxima",
1260 "1 Local Maxima paired with #lambda_{0}^{2}>0.3","N Local Maxima paired with 0.1<#lambda_{0}^{2}<0.3"};
85c4406e 1261
1262 for (Int_t i = 0; i < 8 ; i++)
5c46c992 1263 {
85c4406e 1264
1265 if (fAnaType == kIMCaloTracks && i > 2 ) continue ;
1266
5c46c992 1267 fhMassPairLocMax[i] = new TH2F
1268 (Form("MassPairLocMax%s",combiName[i].Data()),
1269 Form("Mass for decay #gamma pair vs E_{pair}, origin #pi^{0}, %s", combiTitle[i].Data()),
85c4406e 1270 nptbins,ptmin,ptmax,nmassbins,massmin,massmax);
5c46c992 1271 fhMassPairLocMax[i]->SetYTitle("Mass (MeV/c^{2})");
1272 fhMassPairLocMax[i]->SetXTitle("E_{pair} (GeV)");
85c4406e 1273 outputContainer->Add(fhMassPairLocMax[i]) ;
5c46c992 1274 }
e7fd282f 1275 }
477d6cee 1276
b5dbb99b 1277 if(fFillTMHisto)
1278 {
09273901 1279 fhTrackMatchedDEta = new TH2F
31ae6d59 1280 ("hTrackMatchedDEta",
09273901 1281 "d#eta of cluster-track vs cluster energy",
85c4406e 1282 nptbins,ptmin,ptmax,nresetabins,resetamin,resetamax);
09273901 1283 fhTrackMatchedDEta->SetYTitle("d#eta");
1284 fhTrackMatchedDEta->SetXTitle("E_{cluster} (GeV)");
1285
1286 fhTrackMatchedDPhi = new TH2F
31ae6d59 1287 ("hTrackMatchedDPhi",
09273901 1288 "d#phi of cluster-track vs cluster energy",
85c4406e 1289 nptbins,ptmin,ptmax,nresphibins,resphimin,resphimax);
09273901 1290 fhTrackMatchedDPhi->SetYTitle("d#phi (rad)");
1291 fhTrackMatchedDPhi->SetXTitle("E_{cluster} (GeV)");
1292
1293 fhTrackMatchedDEtaDPhi = new TH2F
31ae6d59 1294 ("hTrackMatchedDEtaDPhi",
09273901 1295 "d#eta vs d#phi of cluster-track vs cluster energy",
85c4406e 1296 nresetabins,resetamin,resetamax,nresphibins,resphimin,resphimax);
09273901 1297 fhTrackMatchedDEtaDPhi->SetYTitle("d#phi (rad)");
85c4406e 1298 fhTrackMatchedDEtaDPhi->SetXTitle("d#eta");
09273901 1299
85c4406e 1300 outputContainer->Add(fhTrackMatchedDEta) ;
09273901 1301 outputContainer->Add(fhTrackMatchedDPhi) ;
1302 outputContainer->Add(fhTrackMatchedDEtaDPhi) ;
b2e375c7 1303
1304 fhTrackMatchedDEtaPos = new TH2F
1305 ("hTrackMatchedDEtaPos",
1306 "d#eta of cluster-track vs cluster energy",
1307 nptbins,ptmin,ptmax,nresetabins,resetamin,resetamax);
1308 fhTrackMatchedDEtaPos->SetYTitle("d#eta");
1309 fhTrackMatchedDEtaPos->SetXTitle("E_{cluster} (GeV)");
1310
1311 fhTrackMatchedDPhiPos = new TH2F
1312 ("hTrackMatchedDPhiPos",
1313 "d#phi of cluster-track vs cluster energy",
1314 nptbins,ptmin,ptmax,nresphibins,resphimin,resphimax);
1315 fhTrackMatchedDPhiPos->SetYTitle("d#phi (rad)");
1316 fhTrackMatchedDPhiPos->SetXTitle("E_{cluster} (GeV)");
1317
1318 fhTrackMatchedDEtaDPhiPos = new TH2F
1319 ("hTrackMatchedDEtaDPhiPos",
1320 "d#eta vs d#phi of cluster-track vs cluster energy",
1321 nresetabins,resetamin,resetamax,nresphibins,resphimin,resphimax);
1322 fhTrackMatchedDEtaDPhiPos->SetYTitle("d#phi (rad)");
1323 fhTrackMatchedDEtaDPhiPos->SetXTitle("d#eta");
1324
1325 outputContainer->Add(fhTrackMatchedDEtaPos) ;
1326 outputContainer->Add(fhTrackMatchedDPhiPos) ;
1327 outputContainer->Add(fhTrackMatchedDEtaDPhiPos) ;
1328
1329 fhTrackMatchedDEtaNeg = new TH2F
1330 ("hTrackMatchedDEtaNeg",
1331 "d#eta of cluster-track vs cluster energy",
1332 nptbins,ptmin,ptmax,nresetabins,resetamin,resetamax);
1333 fhTrackMatchedDEtaNeg->SetYTitle("d#eta");
1334 fhTrackMatchedDEtaNeg->SetXTitle("E_{cluster} (GeV)");
1335
1336 fhTrackMatchedDPhiNeg = new TH2F
1337 ("hTrackMatchedDPhiNeg",
1338 "d#phi of cluster-track vs cluster energy",
1339 nptbins,ptmin,ptmax,nresphibins,resphimin,resphimax);
1340 fhTrackMatchedDPhiNeg->SetYTitle("d#phi (rad)");
1341 fhTrackMatchedDPhiNeg->SetXTitle("E_{cluster} (GeV)");
1342
1343 fhTrackMatchedDEtaDPhiNeg = new TH2F
1344 ("hTrackMatchedDEtaDPhiNeg",
1345 "d#eta vs d#phi of cluster-track vs cluster energy",
1346 nresetabins,resetamin,resetamax,nresphibins,resphimin,resphimax);
1347 fhTrackMatchedDEtaDPhiNeg->SetYTitle("d#phi (rad)");
1348 fhTrackMatchedDEtaDPhiNeg->SetXTitle("d#eta");
1349
1350 outputContainer->Add(fhTrackMatchedDEtaNeg) ;
1351 outputContainer->Add(fhTrackMatchedDPhiNeg) ;
1352 outputContainer->Add(fhTrackMatchedDEtaDPhiNeg) ;
31ae6d59 1353
85c4406e 1354 fhdEdx = new TH2F ("hdEdx","matched track <dE/dx> vs cluster E ", nptbins,ptmin,ptmax,ndedxbins, dedxmin, dedxmax);
31ae6d59 1355 fhdEdx->SetXTitle("E (GeV)");
1356 fhdEdx->SetYTitle("<dE/dx>");
85c4406e 1357 outputContainer->Add(fhdEdx);
31ae6d59 1358
85c4406e 1359 fhEOverP = new TH2F ("hEOverP","matched track E/p vs cluster E ", nptbins,ptmin,ptmax,nPoverEbins,pOverEmin,pOverEmax);
31ae6d59 1360 fhEOverP->SetXTitle("E (GeV)");
1361 fhEOverP->SetYTitle("E/p");
85c4406e 1362 outputContainer->Add(fhEOverP);
b5dbb99b 1363
1364 if(fCalorimeter=="EMCAL")
1365 {
85c4406e 1366 fhEOverPNoTRD = new TH2F ("hEOverPNoTRD","matched track E/p vs cluster E, SM not behind TRD ", nptbins,ptmin,ptmax,nPoverEbins,pOverEmin,pOverEmax);
b5dbb99b 1367 fhEOverPNoTRD->SetXTitle("E (GeV)");
1368 fhEOverPNoTRD->SetYTitle("E/p");
85c4406e 1369 outputContainer->Add(fhEOverPNoTRD);
1370 }
31ae6d59 1371
764ab1f4 1372 if(IsDataMC() && fFillTMHisto)
31ae6d59 1373 {
5dde270e 1374 fhTrackMatchedMCParticleE = new TH2F
1375 ("hTrackMatchedMCParticleE",
31ae6d59 1376 "Origin of particle vs energy",
85c4406e 1377 nptbins,ptmin,ptmax,8,0,8);
1378 fhTrackMatchedMCParticleE->SetXTitle("E (GeV)");
5dde270e 1379 //fhTrackMatchedMCParticleE->SetYTitle("Particle type");
1380
1381 fhTrackMatchedMCParticleE->GetYaxis()->SetBinLabel(1 ,"Photon");
1382 fhTrackMatchedMCParticleE->GetYaxis()->SetBinLabel(2 ,"Electron");
1383 fhTrackMatchedMCParticleE->GetYaxis()->SetBinLabel(3 ,"Meson Merged");
1384 fhTrackMatchedMCParticleE->GetYaxis()->SetBinLabel(4 ,"Rest");
1385 fhTrackMatchedMCParticleE->GetYaxis()->SetBinLabel(5 ,"Conv. Photon");
1386 fhTrackMatchedMCParticleE->GetYaxis()->SetBinLabel(6 ,"Conv. Electron");
1387 fhTrackMatchedMCParticleE->GetYaxis()->SetBinLabel(7 ,"Conv. Merged");
1388 fhTrackMatchedMCParticleE->GetYaxis()->SetBinLabel(8 ,"Conv. Rest");
1389
1390 outputContainer->Add(fhTrackMatchedMCParticleE);
1391
1392 fhTrackMatchedMCParticleDEta = new TH2F
1393 ("hTrackMatchedMCParticleDEta",
1394 "Origin of particle vs #eta residual",
1395 nresetabins,resetamin,resetamax,8,0,8);
1396 fhTrackMatchedMCParticleDEta->SetXTitle("#Delta #eta");
1397 //fhTrackMatchedMCParticleDEta->SetYTitle("Particle type");
31ae6d59 1398
5dde270e 1399 fhTrackMatchedMCParticleDEta->GetYaxis()->SetBinLabel(1 ,"Photon");
1400 fhTrackMatchedMCParticleDEta->GetYaxis()->SetBinLabel(2 ,"Electron");
1401 fhTrackMatchedMCParticleDEta->GetYaxis()->SetBinLabel(3 ,"Meson Merged");
1402 fhTrackMatchedMCParticleDEta->GetYaxis()->SetBinLabel(4 ,"Rest");
1403 fhTrackMatchedMCParticleDEta->GetYaxis()->SetBinLabel(5 ,"Conv. Photon");
1404 fhTrackMatchedMCParticleDEta->GetYaxis()->SetBinLabel(6 ,"Conv. Electron");
1405 fhTrackMatchedMCParticleDEta->GetYaxis()->SetBinLabel(7 ,"Conv. Merged");
1406 fhTrackMatchedMCParticleDEta->GetYaxis()->SetBinLabel(8 ,"Conv. Rest");
1407
1408 outputContainer->Add(fhTrackMatchedMCParticleDEta);
85c4406e 1409
5dde270e 1410 fhTrackMatchedMCParticleDPhi = new TH2F
1411 ("hTrackMatchedMCParticleDPhi",
1412 "Origin of particle vs #phi residual",
1413 nresphibins,resphimin,resphimax,8,0,8);
1414 fhTrackMatchedMCParticleDPhi->SetXTitle("#Delta #phi");
1415 //fhTrackMatchedMCParticleDPhi->SetYTitle("Particle type");
1416
1417 fhTrackMatchedMCParticleDPhi->GetYaxis()->SetBinLabel(1 ,"Photon");
1418 fhTrackMatchedMCParticleDPhi->GetYaxis()->SetBinLabel(2 ,"Electron");
1419 fhTrackMatchedMCParticleDPhi->GetYaxis()->SetBinLabel(3 ,"Meson Merged");
1420 fhTrackMatchedMCParticleDPhi->GetYaxis()->SetBinLabel(4 ,"Rest");
1421 fhTrackMatchedMCParticleDPhi->GetYaxis()->SetBinLabel(5 ,"Conv. Photon");
1422 fhTrackMatchedMCParticleDPhi->GetYaxis()->SetBinLabel(6 ,"Conv. Electron");
1423 fhTrackMatchedMCParticleDPhi->GetYaxis()->SetBinLabel(7 ,"Conv. Merged");
1424 fhTrackMatchedMCParticleDPhi->GetYaxis()->SetBinLabel(8 ,"Conv. Rest");
1425
1426 outputContainer->Add(fhTrackMatchedMCParticleDPhi);
85c4406e 1427
31ae6d59 1428
31ae6d59 1429 }
85c4406e 1430 }
09273901 1431
b5dbb99b 1432 if(fFillWeightHistograms)
1433 {
78a28af3 1434 fhECellClusterRatio = new TH2F ("hECellClusterRatio"," cell energy / cluster energy vs cluster energy, for selected decay photons from neutral meson",
85c4406e 1435 nptbins,ptmin,ptmax, 100,0,1.);
78a28af3 1436 fhECellClusterRatio->SetXTitle("E_{cluster} (GeV) ");
1437 fhECellClusterRatio->SetYTitle("E_{cell i}/E_{cluster}");
1438 outputContainer->Add(fhECellClusterRatio);
1439
1440 fhECellClusterLogRatio = new TH2F ("hECellClusterLogRatio"," Log(cell energy / cluster energy) vs cluster energy, for selected decay photons from neutral meson",
85c4406e 1441 nptbins,ptmin,ptmax, 100,-10,0);
78a28af3 1442 fhECellClusterLogRatio->SetXTitle("E_{cluster} (GeV) ");
1a72f6c5 1443 fhECellClusterLogRatio->SetYTitle("Log (E_{max cell}/E_{cluster})");
78a28af3 1444 outputContainer->Add(fhECellClusterLogRatio);
1445
1446 fhEMaxCellClusterRatio = new TH2F ("hEMaxCellClusterRatio"," max cell energy / cluster energy vs cluster energy, for selected decay photons from neutral meson",
85c4406e 1447 nptbins,ptmin,ptmax, 100,0,1.);
78a28af3 1448 fhEMaxCellClusterRatio->SetXTitle("E_{cluster} (GeV) ");
1449 fhEMaxCellClusterRatio->SetYTitle("E_{max cell}/E_{cluster}");
1450 outputContainer->Add(fhEMaxCellClusterRatio);
1451
1452 fhEMaxCellClusterLogRatio = new TH2F ("hEMaxCellClusterLogRatio"," Log(max cell energy / cluster energy) vs cluster energy, for selected decay photons from neutral meson",
85c4406e 1453 nptbins,ptmin,ptmax, 100,-10,0);
78a28af3 1454 fhEMaxCellClusterLogRatio->SetXTitle("E_{cluster} (GeV) ");
1a72f6c5 1455 fhEMaxCellClusterLogRatio->SetYTitle("Log (E_{max cell}/E_{cluster})");
78a28af3 1456 outputContainer->Add(fhEMaxCellClusterLogRatio);
1457
b5dbb99b 1458 for(Int_t iw = 0; iw < 14; iw++)
1459 {
1a72f6c5 1460 fhLambda0ForW0[iw] = new TH2F (Form("hLambda0ForW0%d",iw),Form("shower shape, #lambda^{2}_{0} vs E, w0 = %1.1f, for selected decay photons from neutral meson",1+0.5*iw),
85c4406e 1461 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
78a28af3 1462 fhLambda0ForW0[iw]->SetXTitle("E_{cluster}");
1463 fhLambda0ForW0[iw]->SetYTitle("#lambda^{2}_{0}");
85c4406e 1464 outputContainer->Add(fhLambda0ForW0[iw]);
78a28af3 1465
85c4406e 1466 // fhLambda1ForW0[iw] = new TH2F (Form("hLambda1ForW0%d",iw),Form("shower shape, #lambda^{2}_{1} vs E, w0 = %1.1f, for selected decay photons from neutral meson",0.5+0.5*iw),
1467 // nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
1468 // fhLambda1ForW0[iw]->SetXTitle("E_{cluster}");
1469 // fhLambda1ForW0[iw]->SetYTitle("#lambda^{2}_{1}");
1470 // outputContainer->Add(fhLambda1ForW0[iw]);
78a28af3 1471
1472 }
85c4406e 1473 }
78a28af3 1474
85c4406e 1475 if(IsDataMC())
b5dbb99b 1476 {
3455f821 1477 if(GetReader()->GetDataType() != AliCaloTrackReader::kMC && fAnaType==kSSCalo)
1478 {
883411b2 1479 fhMCPi0PtGenRecoFraction = new TH2F("hMCPi0PtGenRecoFraction","Number of clusters from #pi^{0} (2 #gamma) identified as #pi^{0} (#eta), pT versus E primary #pi^{0} / E reco",
85c4406e 1480 nptbins,ptmin,ptmax,200,0,2);
883411b2 1481 fhMCPi0PtGenRecoFraction->SetXTitle("p^{rec}_{T} (GeV/c)");
1482 fhMCPi0PtGenRecoFraction->SetYTitle("E^{ #pi^{0} mother} / E^{rec}");
85c4406e 1483 outputContainer->Add(fhMCPi0PtGenRecoFraction) ;
1484
883411b2 1485 fhMCEtaPtGenRecoFraction = new TH2F("hMCEtaPtGenRecoFraction","Number of clusters from #eta (2 #gamma) identified as #pi^{0} (#eta),pT versus E primary #eta / E reco",
85c4406e 1486 nptbins,ptmin,ptmax,200,0,2);
883411b2 1487 fhMCEtaPtGenRecoFraction->SetXTitle("p^{rec}_{T} (GeV/c)");
1488 fhMCEtaPtGenRecoFraction->SetYTitle("E^{ #eta mother} / E^{rec}");
85c4406e 1489 outputContainer->Add(fhMCEtaPtGenRecoFraction) ;
51a0ace5 1490
85c4406e 1491 fhMCPi0DecayPt = new TH1F("hMCPi0DecayPt","Number of #gamma from #pi^{0} decay identified as #pi^{0} (#eta)",nptbins,ptmin,ptmax);
3455f821 1492 fhMCPi0DecayPt->SetYTitle("N");
1493 fhMCPi0DecayPt->SetXTitle("p^{rec}_{T} (GeV/c)");
85c4406e 1494 outputContainer->Add(fhMCPi0DecayPt) ;
3455f821 1495
883411b2 1496 fhMCPi0DecayPtFraction = new TH2F("hMCPi0DecayPtFraction","Number of #gamma from #pi^{0} decay identified as #pi^{0} (#eta), pT versus E primary #gamma / E primary #pi^{0}",
85c4406e 1497 nptbins,ptmin,ptmax,100,0,1);
3455f821 1498 fhMCPi0DecayPtFraction->SetXTitle("p^{rec}_{T} (GeV/c)");
883411b2 1499 fhMCPi0DecayPtFraction->SetYTitle("E^{gen} / E^{gen-mother}");
85c4406e 1500 outputContainer->Add(fhMCPi0DecayPtFraction) ;
3455f821 1501
85c4406e 1502 fhMCEtaDecayPt = new TH1F("hMCEtaDecayPt","Number of #gamma from #eta decay identified as #pi^{0} (#eta)",nptbins,ptmin,ptmax);
3455f821 1503 fhMCEtaDecayPt->SetYTitle("N");
1504 fhMCEtaDecayPt->SetXTitle("p^{rec}_{T} (GeV/c)");
85c4406e 1505 outputContainer->Add(fhMCEtaDecayPt) ;
3455f821 1506
883411b2 1507 fhMCEtaDecayPtFraction = new TH2F("hMCEtaDecayPtFraction","Number of #gamma from #eta decay identified as #pi^{0} (#eta), pT versus E primary #gamma / E primary #eta",
85c4406e 1508 nptbins,ptmin,ptmax,100,0,1);
3455f821 1509 fhMCEtaDecayPtFraction->SetXTitle("p^{rec}_{T} (GeV/c)");
883411b2 1510 fhMCEtaDecayPtFraction->SetYTitle("E^{gen} / E^{gen-mother}");
85c4406e 1511 outputContainer->Add(fhMCEtaDecayPtFraction) ;
3455f821 1512
85c4406e 1513 fhMCOtherDecayPt = new TH1F("hMCOtherDecayPt","Number of #gamma decay (not #eta or #pi^{0}) identified as #pi^{0} (#eta)",nptbins,ptmin,ptmax);
3455f821 1514 fhMCOtherDecayPt->SetYTitle("N");
1515 fhMCOtherDecayPt->SetXTitle("p^{rec}_{T} (GeV/c)");
85c4406e 1516 outputContainer->Add(fhMCOtherDecayPt) ;
3455f821 1517
1518 }
85c4406e 1519
1520 if((GetReader()->GetDataType() == AliCaloTrackReader::kMC && fAnaType!=kSSCalo) ||
b5dbb99b 1521 GetReader()->GetDataType() != AliCaloTrackReader::kMC)
1522 {
477d6cee 1523
b5dbb99b 1524 fhAnglePairMCPi0 = new TH2F
1525 ("AnglePairMCPi0",
85c4406e 1526 "Angle between decay #gamma pair vs E_{pair}, origin #pi^{0}",nptbins,ptmin,ptmax,250,0,0.5);
b5dbb99b 1527 fhAnglePairMCPi0->SetYTitle("#alpha (rad)");
1528 fhAnglePairMCPi0->SetXTitle("E_{pair} (GeV)");
85c4406e 1529 outputContainer->Add(fhAnglePairMCPi0) ;
1530
af722ce4 1531 if (fAnaType!= kSSCalo)
1532 {
1533 fhAnglePairMCEta = new TH2F
1534 ("AnglePairMCEta",
85c4406e 1535 "Angle between decay #gamma pair vs E_{pair}, origin #eta",nptbins,ptmin,ptmax,250,0,0.5);
af722ce4 1536 fhAnglePairMCEta->SetYTitle("#alpha (rad)");
1537 fhAnglePairMCEta->SetXTitle("E_{pair} (GeV)");
85c4406e 1538 outputContainer->Add(fhAnglePairMCEta) ;
af722ce4 1539
1540 fhMassPairMCPi0 = new TH2F
1541 ("MassPairMCPi0",
85c4406e 1542 "Mass for decay #gamma pair vs E_{pair}, origin #pi^{0}",nptbins,ptmin,ptmax,nmassbins,massmin,massmax);
af722ce4 1543 fhMassPairMCPi0->SetYTitle("Mass (MeV/c^{2})");
1544 fhMassPairMCPi0->SetXTitle("E_{pair} (GeV)");
85c4406e 1545 outputContainer->Add(fhMassPairMCPi0) ;
af722ce4 1546
1547 fhMassPairMCEta = new TH2F
1548 ("MassPairMCEta",
85c4406e 1549 "Mass for decay #gamma pair vs E_{pair}, origin #eta",nptbins,ptmin,ptmax,nmassbins,massmin,massmax);
af722ce4 1550 fhMassPairMCEta->SetYTitle("Mass (MeV/c^{2})");
1551 fhMassPairMCEta->SetXTitle("E_{pair} (GeV)");
85c4406e 1552 outputContainer->Add(fhMassPairMCEta) ;
af722ce4 1553 }
e4ef72be 1554
3455f821 1555 for(Int_t i = 0; i < 6; i++)
85c4406e 1556 {
3455f821 1557
40d3ce60 1558 fhMCE[i] = new TH1F
1559 (Form("hE_MC%s",pname[i].Data()),
1560 Form("Identified as #pi^{0} (#eta), cluster from %s",
1561 ptype[i].Data()),
85c4406e 1562 nptbins,ptmin,ptmax);
40d3ce60 1563 fhMCE[i]->SetYTitle("N");
1564 fhMCE[i]->SetXTitle("E (GeV)");
85c4406e 1565 outputContainer->Add(fhMCE[i]) ;
40d3ce60 1566
3455f821 1567 fhMCPt[i] = new TH1F
1568 (Form("hPt_MC%s",pname[i].Data()),
1569 Form("Identified as #pi^{0} (#eta), cluster from %s",
1570 ptype[i].Data()),
85c4406e 1571 nptbins,ptmin,ptmax);
3455f821 1572 fhMCPt[i]->SetYTitle("N");
1573 fhMCPt[i]->SetXTitle("p_{T} (GeV/c)");
85c4406e 1574 outputContainer->Add(fhMCPt[i]) ;
3455f821 1575
17f5b4b6 1576 fhMCPtCentrality[i] = new TH2F
1577 (Form("hPtCentrality_MC%s",pname[i].Data()),
1578 Form("Identified as #pi^{0} (#eta), cluster from %s",
1579 ptype[i].Data()),
1580 nptbins,ptmin,ptmax, 100,0,100);
1581 fhMCPtCentrality[i]->SetYTitle("centrality");
1582 fhMCPtCentrality[i]->SetXTitle("p_{T} (GeV/c)");
1583 outputContainer->Add(fhMCPtCentrality[i]) ;
1584
40d3ce60 1585 if(fAnaType == kSSCalo)
1586 {
6e66993c 1587
1588 fhMCNLocMaxPt[i] = new TH2F
1589 (Form("hNLocMaxPt_MC%s",pname[i].Data()),
1590 Form("cluster from %s, pT of cluster, for NLM",ptype[i].Data()),
1591 nptbins,ptmin,ptmax,10,0,10);
1592 fhMCNLocMaxPt[i] ->SetYTitle("N maxima");
1593 fhMCNLocMaxPt[i] ->SetXTitle("p_{T} (GeV/c)");
1594 outputContainer->Add(fhMCNLocMaxPt[i]) ;
85c4406e 1595
40d3ce60 1596 fhMCEReject[i] = new TH1F
1597 (Form("hEReject_MC%s",pname[i].Data()),
1598 Form("Rejected as #pi^{0} (#eta), cluster from %s",
1599 ptype[i].Data()),
85c4406e 1600 nptbins,ptmin,ptmax);
40d3ce60 1601 fhMCEReject[i]->SetYTitle("N");
1602 fhMCEReject[i]->SetXTitle("E (GeV)");
85c4406e 1603 outputContainer->Add(fhMCEReject[i]) ;
40d3ce60 1604
1605 fhMCPtReject[i] = new TH1F
1606 (Form("hPtReject_MC%s",pname[i].Data()),
1607 Form("Rejected as #pi^{0} (#eta), cluster from %s",
1608 ptype[i].Data()),
85c4406e 1609 nptbins,ptmin,ptmax);
40d3ce60 1610 fhMCPtReject[i]->SetYTitle("N");
1611 fhMCPtReject[i]->SetXTitle("p_{T} (GeV/c)");
85c4406e 1612 outputContainer->Add(fhMCPtReject[i]) ;
40d3ce60 1613 }
1614
3455f821 1615 fhMCPhi[i] = new TH2F
1616 (Form("hPhi_MC%s",pname[i].Data()),
1617 Form("Identified as #pi^{0} (#eta), cluster from %s",ptype[i].Data()),
85c4406e 1618 nptbins,ptmin,ptmax,nphibins,phimin,phimax);
3455f821 1619 fhMCPhi[i]->SetYTitle("#phi");
1620 fhMCPhi[i]->SetXTitle("p_{T} (GeV/c)");
85c4406e 1621 outputContainer->Add(fhMCPhi[i]) ;
3455f821 1622
1623 fhMCEta[i] = new TH2F
1624 (Form("hEta_MC%s",pname[i].Data()),
1625 Form("Identified as #pi^{0} (#eta), cluster from %s",
85c4406e 1626 ptype[i].Data()),nptbins,ptmin,ptmax,netabins,etamin,etamax);
3455f821 1627 fhMCEta[i]->SetYTitle("#eta");
1628 fhMCEta[i]->SetXTitle("p_{T} (GeV/c)");
1629 outputContainer->Add(fhMCEta[i]) ;
85c4406e 1630
29250849 1631 fhMCMassPt[i] = new TH2F
1632 (Form("hMassPt_MC%s",pname[i].Data()),
1633 Form("all pairs mass: p_{T} vs massfrom %s",ptype[i].Data()),
1634 nptbins,ptmin,ptmax, nmassbins,massmin,massmax);
1635 fhMCMassPt[i]->SetYTitle("mass (GeV/c^{2})");
1636 fhMCMassPt[i]->SetXTitle("p_{T} (GeV/c)");
1637 outputContainer->Add(fhMCMassPt[i]) ;
3455f821 1638
29250849 1639 fhMCSelectedMassPt[i] = new TH2F
1640 (Form("hSelectedMassPt_MC%s",pname[i].Data()),
1641 Form("Selected #pi^{0} (#eta) pairs mass: p_{T} vs massfrom %s",ptype[i].Data()),
1642 nptbins,ptmin,ptmax, nmassbins,massmin,massmax);
1643 fhMCSelectedMassPt[i]->SetYTitle("mass (GeV/c^{2})");
1644 fhMCSelectedMassPt[i]->SetXTitle("p_{T} (GeV/c)");
1645 outputContainer->Add(fhMCSelectedMassPt[i]) ;
85c4406e 1646
1253480f 1647 if(fAnaType == kSSCalo)
1648 {
1649 fhMCMassPtNoOverlap[i] = new TH2F
1650 (Form("hMassPtNoOverlap_MC%s",pname[i].Data()),
1651 Form("all pairs mass: p_{T} vs massfrom %s, no overlap",ptype[i].Data()),
1652 nptbins,ptmin,ptmax, nmassbins,massmin,massmax);
1653 fhMCMassPt[i]->SetYTitle("mass (GeV/c^{2})");
1654 fhMCMassPt[i]->SetXTitle("p_{T} (GeV/c)");
1655 outputContainer->Add(fhMCMassPtNoOverlap[i]) ;
1656
1657 fhMCSelectedMassPtNoOverlap[i] = new TH2F
1658 (Form("hSelectedMassPtNoOverlap_MC%s",pname[i].Data()),
1659 Form("Selected #pi^{0} (#eta) pairs mass: p_{T} vs massfrom %s, no overlap",ptype[i].Data()),
1660 nptbins,ptmin,ptmax, nmassbins,massmin,massmax);
1661 fhMCSelectedMassPtNoOverlap[i]->SetYTitle("mass (GeV/c^{2})");
1662 fhMCSelectedMassPtNoOverlap[i]->SetXTitle("p_{T} (GeV/c)");
1663 outputContainer->Add(fhMCSelectedMassPtNoOverlap[i]) ;
1664 }
3455f821 1665
1666 if( fFillSelectClHisto )
1667 {
e4ef72be 1668 fhEMCLambda0[i] = new TH2F(Form("hELambda0_MC%s",pname[i].Data()),
1669 Form("Selected pair, cluster from %s : E vs #lambda_{0}^{2}",ptype[i].Data()),
85c4406e 1670 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
e4ef72be 1671 fhEMCLambda0[i]->SetYTitle("#lambda_{0}^{2}");
1672 fhEMCLambda0[i]->SetXTitle("E (GeV)");
85c4406e 1673 outputContainer->Add(fhEMCLambda0[i]) ;
34c16486 1674
e4ef72be 1675 fhEMCLambda1[i] = new TH2F(Form("hELambda1_MC%s",pname[i].Data()),
1676 Form("Selected pair, cluster from %s : E vs #lambda_{1}^{2}",ptype[i].Data()),
85c4406e 1677 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
e4ef72be 1678 fhEMCLambda1[i]->SetYTitle("#lambda_{1}^{2}");
1679 fhEMCLambda1[i]->SetXTitle("E (GeV)");
85c4406e 1680 outputContainer->Add(fhEMCLambda1[i]) ;
34c16486 1681
e4ef72be 1682 fhEMCDispersion[i] = new TH2F(Form("hEDispersion_MC%s",pname[i].Data()),
1683 Form("Selected pair, cluster from %s : E vs dispersion^{2}",ptype[i].Data()),
85c4406e 1684 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
e4ef72be 1685 fhEMCDispersion[i]->SetYTitle("D^{2}");
1686 fhEMCDispersion[i]->SetXTitle("E (GeV)");
85c4406e 1687 outputContainer->Add(fhEMCDispersion[i]) ;
34c16486 1688
e4ef72be 1689 if(fCalorimeter=="EMCAL")
34c16486 1690 {
e4ef72be 1691 fhEMCLambda0NoTRD[i] = new TH2F(Form("hELambda0NoTRD_MC%s",pname[i].Data()),
1692 Form("Selected pair, cluster from %s : E vs #lambda_{0}^{2}, NoTRD",ptype[i].Data()),
85c4406e 1693 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
e4ef72be 1694 fhEMCLambda0NoTRD[i]->SetYTitle("#lambda_{0}^{2}");
1695 fhEMCLambda0NoTRD[i]->SetXTitle("E (GeV)");
85c4406e 1696 outputContainer->Add(fhEMCLambda0NoTRD[i]) ;
bfdcf7fb 1697
764ab1f4 1698 if(!fFillOnlySimpleSSHisto)
e4ef72be 1699 {
764ab1f4 1700 fhMCEDispEta[i] = new TH2F (Form("hEDispEtaE_MC%s",pname[i].Data()),
1701 Form("cluster from %s : #sigma^{2}_{#eta #eta} = #Sigma w_{i}(#eta_{i} - <#eta>)^{2}/ #Sigma w_{i} vs E",ptype[i].Data()),
85c4406e 1702 nptbins,ptmin,ptmax, ssbins,ssmin,ssmax);
764ab1f4 1703 fhMCEDispEta[i]->SetXTitle("E (GeV)");
1704 fhMCEDispEta[i]->SetYTitle("#sigma^{2}_{#eta #eta}");
85c4406e 1705 outputContainer->Add(fhMCEDispEta[i]);
764ab1f4 1706
1707 fhMCEDispPhi[i] = new TH2F (Form("hEDispPhiE_MC%s",pname[i].Data()),
1708 Form("cluster from %s : #sigma^{2}_{#phi #phi} = #Sigma w_{i}(#phi_{i} - <#phi>)^{2} / #Sigma w_{i} vs E",ptype[i].Data()),
85c4406e 1709 nptbins,ptmin,ptmax, ssbins,ssmin,ssmax);
764ab1f4 1710 fhMCEDispPhi[i]->SetXTitle("E (GeV)");
1711 fhMCEDispPhi[i]->SetYTitle("#sigma^{2}_{#phi #phi}");
85c4406e 1712 outputContainer->Add(fhMCEDispPhi[i]);
764ab1f4 1713
1714 fhMCESumEtaPhi[i] = new TH2F (Form("hESumEtaPhiE_MC%s",pname[i].Data()),
85c4406e 1715 Form("cluster from %s : #delta^{2}_{#eta #phi} = #Sigma w_{i}(#phi_{i} #eta_{i} ) / #Sigma w_{i} - <#phi><#eta> vs E",ptype[i].Data()),
1716 nptbins,ptmin,ptmax, 2*ssbins,-ssmax,ssmax);
764ab1f4 1717 fhMCESumEtaPhi[i]->SetXTitle("E (GeV)");
1718 fhMCESumEtaPhi[i]->SetYTitle("#delta^{2}_{#eta #phi}");
1719 outputContainer->Add(fhMCESumEtaPhi[i]);
e4ef72be 1720
764ab1f4 1721 fhMCEDispEtaPhiDiff[i] = new TH2F (Form("hEDispEtaPhiDiffE_MC%s",pname[i].Data()),
85c4406e 1722 Form("cluster from %s : #sigma^{2}_{#phi #phi} - #sigma^{2}_{#eta #eta} vs E",ptype[i].Data()),
1723 nptbins,ptmin,ptmax,200,-10,10);
764ab1f4 1724 fhMCEDispEtaPhiDiff[i]->SetXTitle("E (GeV)");
1725 fhMCEDispEtaPhiDiff[i]->SetYTitle("#sigma^{2}_{#phi #phi}-#sigma^{2}_{#eta #eta}");
85c4406e 1726 outputContainer->Add(fhMCEDispEtaPhiDiff[i]);
e4ef72be 1727
764ab1f4 1728 fhMCESphericity[i] = new TH2F (Form("hESphericity_MC%s",pname[i].Data()),
85c4406e 1729 Form("cluster from %s : (#sigma^{2}_{#phi #phi} - #sigma^{2}_{#eta #eta}) / (#sigma^{2}_{#eta #eta} + #sigma^{2}_{#phi #phi}) vs E",ptype[i].Data()),
1730 nptbins,ptmin,ptmax, 200,-1,1);
764ab1f4 1731 fhMCESphericity[i]->SetXTitle("E (GeV)");
1732 fhMCESphericity[i]->SetYTitle("s = (#sigma^{2}_{#phi #phi} - #sigma^{2}_{#eta #eta}) / (#sigma^{2}_{#eta #eta} + #sigma^{2}_{#phi #phi})");
1733 outputContainer->Add(fhMCESphericity[i]);
e4ef72be 1734
764ab1f4 1735 for(Int_t ie = 0; ie < 7; ie++)
1736 {
1737 fhMCDispEtaDispPhi[ie][i] = new TH2F (Form("hMCDispEtaDispPhi_EBin%d_MC%s",ie,pname[i].Data()),
85c4406e 1738 Form("cluster from %s : #sigma^{2}_{#phi #phi} vs #sigma^{2}_{#eta #eta} for %d < E < %d GeV",pname[i].Data(),bin[ie],bin[ie+1]),
1739 ssbins,ssmin,ssmax , ssbins,ssmin,ssmax);
764ab1f4 1740 fhMCDispEtaDispPhi[ie][i]->SetXTitle("#sigma^{2}_{#eta #eta}");
1741 fhMCDispEtaDispPhi[ie][i]->SetYTitle("#sigma^{2}_{#phi #phi}");
85c4406e 1742 outputContainer->Add(fhMCDispEtaDispPhi[ie][i]);
764ab1f4 1743
1744 fhMCLambda0DispEta[ie][i] = new TH2F (Form("hMCLambda0DispEta_EBin%d_MC%s",ie,pname[i].Data()),
85c4406e 1745 Form("cluster from %s : #lambda^{2}_{0} vs #sigma^{2}_{#eta #eta} for %d < E < %d GeV",pname[i].Data(),bin[ie],bin[ie+1]),
1746 ssbins,ssmin,ssmax , ssbins,ssmin,ssmax);
764ab1f4 1747 fhMCLambda0DispEta[ie][i]->SetXTitle("#lambda^{2}_{0}");
1748 fhMCLambda0DispEta[ie][i]->SetYTitle("#sigma^{2}_{#phi #phi}");
85c4406e 1749 outputContainer->Add(fhMCLambda0DispEta[ie][i]);
764ab1f4 1750
1751 fhMCLambda0DispPhi[ie][i] = new TH2F (Form("hMCLambda0DispPhi_EBin%d_MC%s",ie,pname[i].Data()),
85c4406e 1752 Form("cluster from %s :#lambda^{2}_{0} vs #sigma^{2}_{#phi #phi} for %d < E < %d GeV",pname[i].Data(),bin[ie],bin[ie+1]),
1753 ssbins,ssmin,ssmax , ssbins,ssmin,ssmax);
764ab1f4 1754 fhMCLambda0DispPhi[ie][i]->SetXTitle("#lambda^{2}_{0}");
1755 fhMCLambda0DispPhi[ie][i]->SetYTitle("#sigma^{2}_{#phi #phi}");
85c4406e 1756 outputContainer->Add(fhMCLambda0DispPhi[ie][i]);
764ab1f4 1757
85c4406e 1758 }
764ab1f4 1759 }
e4ef72be 1760 }
1761
1762 fhEMCLambda0FracMaxCellCut[i] = new TH2F(Form("hELambda0FracMaxCellCut_MC%s",pname[i].Data()),
1763 Form("Selected pair, cluster from %s : E vs #lambda_{0}^{2}, Max cell fraction of energy < 0.5 ",ptype[i].Data()),
85c4406e 1764 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
e4ef72be 1765 fhEMCLambda0FracMaxCellCut[i]->SetYTitle("#lambda_{0}^{2}");
1766 fhEMCLambda0FracMaxCellCut[i]->SetXTitle("E (GeV)");
85c4406e 1767 outputContainer->Add(fhEMCLambda0FracMaxCellCut[i]) ;
e4ef72be 1768
1769 fhEMCFracMaxCell[i] = new TH2F(Form("hEFracMaxCell_MC%s",pname[i].Data()),
1770 Form("Selected pair, cluster from %s : E vs Max cell fraction of energy",ptype[i].Data()),
85c4406e 1771 nptbins,ptmin,ptmax,100,0,1);
e4ef72be 1772 fhEMCFracMaxCell[i]->SetYTitle("Fraction");
1773 fhEMCFracMaxCell[i]->SetXTitle("E (GeV)");
85c4406e 1774 outputContainer->Add(fhEMCFracMaxCell[i]) ;
e4ef72be 1775
1776 }//
1777 } // shower shape histo
34c16486 1778
521636d2 1779 } //Not MC reader
477d6cee 1780 }//Histos with MC
1781
4650f5cf 1782 if(fAnaType==kSSCalo)
1783 {
85c4406e 1784 fhAsymmetry = new TH2F ("hAsymmetry","A = ( E1 - E2 ) / ( E1 + E2 ) vs E",
1785 nptbins,ptmin,ptmax, 200, -1,1);
4650f5cf 1786 fhAsymmetry->SetXTitle("E (GeV)");
1787 fhAsymmetry->SetYTitle("A = ( E1 - E2 ) / ( E1 + E2 )");
1788 outputContainer->Add(fhAsymmetry);
1789
85c4406e 1790 fhSelectedAsymmetry = new TH2F ("hSelectedAsymmetry","A = ( E1 - E2 ) / ( E1 + E2 ) vs E",
1791 nptbins,ptmin,ptmax, 200, -1,1);
4650f5cf 1792 fhSelectedAsymmetry->SetXTitle("E (GeV)");
1793 fhSelectedAsymmetry->SetYTitle("A = ( E1 - E2 ) / ( E1 + E2 )");
1794 outputContainer->Add(fhSelectedAsymmetry);
1795
cfdf2b91 1796 fhSplitE = new TH1F
1797 ("hSplitE","Selected #pi^{0} (#eta) pairs energy sum of split sub-clusters",nptbins,ptmin,ptmax);
1798 fhSplitE->SetYTitle("counts");
1799 fhSplitE->SetXTitle("E (GeV)");
1800 outputContainer->Add(fhSplitE) ;
1801
1802 fhSplitPt = new TH1F
1803 ("hSplitPt","Selected #pi^{0} (#eta) pairs pT sum of split sub-clusters",nptbins,ptmin,ptmax);
1804 fhSplitPt->SetYTitle("counts");
1805 fhSplitPt->SetXTitle("p_{T} (GeV/c)");
1806 outputContainer->Add(fhSplitPt) ;
1807
29250849 1808
1809 fhSplitPtPhi = new TH2F
1810 ("hSplitPtPhi","Selected #pi^{0} (#eta) pairs: sum split sub-cluster p_{T} vs #phi",nptbins,ptmin,ptmax, nphibins,phimin,phimax);
1811 fhSplitPtPhi->SetYTitle("#phi (rad)");
1812 fhSplitPtPhi->SetXTitle("p_{T} (GeV/c)");
1813 outputContainer->Add(fhSplitPtPhi) ;
1814
1815 fhSplitPtEta = new TH2F
1816 ("hSplitPtEta","Selected #pi^{0} (#eta) pairs: sum split sub-cluster p_{T} vs #eta",nptbins,ptmin,ptmax,netabins,etamin,etamax);
1817 fhSplitPtEta->SetYTitle("#eta");
1818 fhSplitPtEta->SetXTitle("p_{T} (GeV/c)");
1819 outputContainer->Add(fhSplitPtEta) ;
85c4406e 1820
29250849 1821
6e66993c 1822 fhNLocMaxSplitPt = new TH2F("hNLocMaxSplitPt","Number of local maxima in cluster",
85c4406e 1823 nptbins,ptmin,ptmax,10,0,10);
6e66993c 1824 fhNLocMaxSplitPt ->SetYTitle("N maxima");
1825 fhNLocMaxSplitPt ->SetXTitle("p_{T} (GeV/c)");
1826 outputContainer->Add(fhNLocMaxSplitPt) ;
85c4406e 1827
6e66993c 1828
29250849 1829 fhMassSplitPt = new TH2F
1253480f 1830 ("hMassSplitPt","all pairs mass: sum split sub-cluster p_{T} vs mass",
1831 nptbins,ptmin,ptmax, nmassbins,massmin,massmax);
29250849 1832 fhMassSplitPt->SetYTitle("mass (GeV/c^{2})");
1833 fhMassSplitPt->SetXTitle("p_{T} (GeV/c)");
1834 outputContainer->Add(fhMassSplitPt) ;
1835
1836 fhSelectedMassSplitPt = new TH2F
1253480f 1837 ("hSelectedMassSplitPt","Selected #pi^{0} (#eta) pairs mass: sum split sub-cluster p_{T} vs mass",
1838 nptbins,ptmin,ptmax, nmassbins,massmin,massmax);
29250849 1839 fhSelectedMassSplitPt->SetYTitle("mass (GeV/c^{2})");
1840 fhSelectedMassSplitPt->SetXTitle("p_{T} (GeV/c)");
1841 outputContainer->Add(fhSelectedMassSplitPt) ;
1842
4650f5cf 1843 if(IsDataMC())
1844 {
1253480f 1845 fhMassSplitPtNoOverlap = new TH2F
1846 ("hMassSplitPtNoOverlap","all pairs mass: sum split sub-cluster p_{T} vs mass, no overlap",
1847 nptbins,ptmin,ptmax, nmassbins,massmin,massmax);
1848 fhMassSplitPtNoOverlap->SetYTitle("mass (GeV/c^{2})");
1849 fhMassSplitPtNoOverlap->SetXTitle("p_{T} (GeV/c)");
1850 outputContainer->Add(fhMassSplitPtNoOverlap) ;
1851
1852 fhSelectedMassSplitPtNoOverlap = new TH2F
1853 ("hSelectedMassSplitPtNoOverlap","Selected #pi^{0} (#eta) pairs mass: sum split sub-cluster p_{T} vs mass, no overlap",
1854 nptbins,ptmin,ptmax, nmassbins,massmin,massmax);
1855 fhSelectedMassSplitPtNoOverlap->SetYTitle("mass (GeV/c^{2})");
1856 fhSelectedMassSplitPtNoOverlap->SetXTitle("p_{T} (GeV/c)");
1857 outputContainer->Add(fhSelectedMassSplitPtNoOverlap) ;
1858
1859
1860 fhMCPi0PtRecoPtPrim = new TH2F
1861 ("hMCPi0PtRecoPtPrim","p_{T,reco} vs p_{T,gen}",
1862 nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
1863 fhMCPi0PtRecoPtPrim ->SetYTitle("p_{T,gen} (GeV/c)");
1864 fhMCPi0PtRecoPtPrim ->SetXTitle("p_{T,reco} (GeV/c)");
1865 outputContainer->Add(fhMCPi0PtRecoPtPrim ) ;
1866
1867 fhMCPi0PtRecoPtPrimNoOverlap = new TH2F
1868 ("hMCPi0PtRecoPtPrimNoOverlap","p_{T,reco} vs p_{T,gen}, no overlap",
1869 nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
1870 fhMCPi0PtRecoPtPrimNoOverlap ->SetYTitle("p_{T,gen} (GeV/c)");
1871 fhMCPi0PtRecoPtPrimNoOverlap ->SetXTitle("p_{T,reco} (GeV/c)");
1872 outputContainer->Add(fhMCPi0PtRecoPtPrimNoOverlap ) ;
1873
1874 fhMCPi0SelectedPtRecoPtPrim = new TH2F
1875 ("hMCPi0SelectedPtRecoPtPrim","p_{T,reco} vs p_{T,gen}",
1876 nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
1877 fhMCPi0SelectedPtRecoPtPrim ->SetYTitle("p_{T,gen} (GeV/c)");
1878 fhMCPi0SelectedPtRecoPtPrim ->SetXTitle("p_{T,reco} (GeV/c)");
1879 outputContainer->Add(fhMCPi0SelectedPtRecoPtPrim ) ;
1880
1881 fhMCPi0SelectedPtRecoPtPrimNoOverlap = new TH2F
1882 ("hMCPi0SelectedPtRecoPtPrimNoOverlap","p_{T,reco} vs p_{T,gen}, no overlap",
1883 nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
1884 fhMCPi0SelectedPtRecoPtPrimNoOverlap ->SetYTitle("p_{T,gen} (GeV/c)");
1885 fhMCPi0SelectedPtRecoPtPrimNoOverlap ->SetXTitle("p_{T,reco} (GeV/c)");
1886 outputContainer->Add(fhMCPi0SelectedPtRecoPtPrimNoOverlap ) ;
1887
1888
1889 fhMCPi0SplitPtRecoPtPrim = new TH2F
1890 ("hMCPi0SplitPtRecoPtPrim","p_{T,reco} (split sum) vs p_{T,gen}",
1891 nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
1892 fhMCPi0SplitPtRecoPtPrim ->SetYTitle("p_{T,gen} (GeV/c)");
1893 fhMCPi0SplitPtRecoPtPrim ->SetXTitle("p_{T,reco} (GeV/c)");
1894 outputContainer->Add(fhMCPi0SplitPtRecoPtPrim ) ;
1895
1896 fhMCPi0SplitPtRecoPtPrimNoOverlap = new TH2F
1897 ("hMCPi0SplitPtRecoPtPrimNoOverlap","p_{T,reco} (split sum) vs p_{T,gen}, no overlap",
1898 nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
1899 fhMCPi0SplitPtRecoPtPrimNoOverlap ->SetYTitle("p_{T,gen} (GeV/c)");
1900 fhMCPi0SplitPtRecoPtPrimNoOverlap ->SetXTitle("p_{T,reco} (GeV/c)");
1901 outputContainer->Add(fhMCPi0SplitPtRecoPtPrimNoOverlap ) ;
1902
1903 fhMCPi0SelectedSplitPtRecoPtPrim = new TH2F
1904 ("hMCPi0SelectedSplitPtRecoPtPrim","p_{T,reco} (split sum) vs p_{T,gen}",
1905 nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
1906 fhMCPi0SelectedSplitPtRecoPtPrim ->SetYTitle("p_{T,gen} (GeV/c)");
1907 fhMCPi0SelectedSplitPtRecoPtPrim ->SetXTitle("p_{T,reco} (GeV/c)");
1908 outputContainer->Add(fhMCPi0SelectedSplitPtRecoPtPrim ) ;
1909
1910 fhMCPi0SelectedSplitPtRecoPtPrimNoOverlap = new TH2F
1911 ("hMCPi0SelectedSplitPtRecoPtPrimNoOverlap","p_{T,reco} (split sum) vs p_{T,gen}, no overlap",
1912 nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
1913 fhMCPi0SelectedSplitPtRecoPtPrimNoOverlap ->SetYTitle("p_{T,gen} (GeV/c)");
1914 fhMCPi0SelectedSplitPtRecoPtPrimNoOverlap ->SetXTitle("p_{T,reco} (GeV/c)");
1915 outputContainer->Add(fhMCPi0SelectedSplitPtRecoPtPrimNoOverlap ) ;
1916
1917 fhMCEtaPtRecoPtPrim = new TH2F
1918 ("hMCEtaPtRecoPtPrim","p_{T,reco} vs p_{T,gen}",
1919 nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
1920 fhMCEtaPtRecoPtPrim ->SetYTitle("p_{T,gen} (GeV/c)");
1921 fhMCEtaPtRecoPtPrim ->SetXTitle("p_{T,reco} (GeV/c)");
1922 outputContainer->Add(fhMCEtaPtRecoPtPrim ) ;
1923
1924 fhMCEtaPtRecoPtPrimNoOverlap = new TH2F
1925 ("hMCEtaPtRecoPtPrimNoOverlap","p_{T,reco} vs p_{T,gen}, no overlap",
1926 nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
1927 fhMCEtaPtRecoPtPrimNoOverlap ->SetYTitle("p_{T,gen} (GeV/c)");
1928 fhMCEtaPtRecoPtPrimNoOverlap ->SetXTitle("p_{T,reco} (GeV/c)");
1929 outputContainer->Add(fhMCEtaPtRecoPtPrimNoOverlap ) ;
1930
1931 fhMCEtaSelectedPtRecoPtPrim = new TH2F
1932 ("hMCEtaSelectedPtRecoPtPrim","p_{T,reco} vs p_{T,gen}",
1933 nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
1934 fhMCEtaSelectedPtRecoPtPrim ->SetYTitle("p_{T,gen} (GeV/c)");
1935 fhMCEtaSelectedPtRecoPtPrim ->SetXTitle("p_{T,reco} (GeV/c)");
1936 outputContainer->Add(fhMCEtaSelectedPtRecoPtPrim ) ;
1937
1938 fhMCEtaSelectedPtRecoPtPrimNoOverlap = new TH2F
1939 ("hMCEtaSelectedPtRecoPtPrimNoOverlap","p_{T,reco} vs p_{T,gen}, no overlap",
1940 nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
1941 fhMCEtaSelectedPtRecoPtPrimNoOverlap ->SetYTitle("p_{T,gen} (GeV/c)");
1942 fhMCEtaSelectedPtRecoPtPrimNoOverlap ->SetXTitle("p_{T,reco} (GeV/c)");
1943 outputContainer->Add(fhMCEtaSelectedPtRecoPtPrimNoOverlap ) ;
1944
1945
1946 fhMCEtaSplitPtRecoPtPrim = new TH2F
1947 ("hMCEtaSplitPtRecoPtPrim","p_{T,reco} (split sum) vs p_{T,gen}",
1948 nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
1949 fhMCEtaSplitPtRecoPtPrim ->SetYTitle("p_{T,gen} (GeV/c)");
1950 fhMCEtaSplitPtRecoPtPrim ->SetXTitle("p_{T,reco} (GeV/c)");
1951 outputContainer->Add(fhMCEtaSplitPtRecoPtPrim ) ;
1952
1953 fhMCEtaSplitPtRecoPtPrimNoOverlap = new TH2F
1954 ("hMCEtaSplitPtRecoPtPrimNoOverlap","p_{T,reco} (split sum) vs p_{T,gen}, no overlap",
1955 nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
1956 fhMCEtaSplitPtRecoPtPrimNoOverlap ->SetYTitle("p_{T,gen} (GeV/c)");
1957 fhMCEtaSplitPtRecoPtPrimNoOverlap ->SetXTitle("p_{T,reco} (GeV/c)");
1958 outputContainer->Add(fhMCEtaSplitPtRecoPtPrimNoOverlap ) ;
1959
1960 fhMCEtaSelectedSplitPtRecoPtPrim = new TH2F
1961 ("hMCEtaSelectedSplitPtRecoPtPrim","p_{T,reco} (split sum) vs p_{T,gen}",
1962 nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
1963 fhMCEtaSelectedSplitPtRecoPtPrim ->SetYTitle("p_{T,gen} (GeV/c)");
1964 fhMCEtaSelectedSplitPtRecoPtPrim ->SetXTitle("p_{T,reco} (GeV/c)");
1965 outputContainer->Add(fhMCEtaSelectedSplitPtRecoPtPrim ) ;
1966
1967 fhMCEtaSelectedSplitPtRecoPtPrimNoOverlap = new TH2F
1968 ("hMCEtaSelectedSplitPtRecoPtPrimNoOverlap","p_{T,reco} (split sum) vs p_{T,gen}, no overlap",
1969 nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
1970 fhMCEtaSelectedSplitPtRecoPtPrimNoOverlap ->SetYTitle("p_{T,gen} (GeV/c)");
1971 fhMCEtaSelectedSplitPtRecoPtPrimNoOverlap ->SetXTitle("p_{T,reco} (GeV/c)");
1972 outputContainer->Add(fhMCEtaSelectedSplitPtRecoPtPrimNoOverlap ) ;
1973
4650f5cf 1974 for(Int_t i = 0; i< 6; i++)
1975 {
1976 fhMCEAsymmetry[i] = new TH2F (Form("hEAsymmetry_MC%s",pname[i].Data()),
85c4406e 1977 Form("cluster from %s : A = ( E1 - E2 ) / ( E1 + E2 ) vs E",ptype[i].Data()),
1978 nptbins,ptmin,ptmax, 200,-1,1);
4650f5cf 1979 fhMCEAsymmetry[i]->SetXTitle("E (GeV)");
1980 fhMCEAsymmetry[i]->SetYTitle("A = ( E1 - E2 ) / ( E1 + E2 )");
1981 outputContainer->Add(fhMCEAsymmetry[i]);
cfdf2b91 1982
1983 fhMCSplitE[i] = new TH1F
1984 (Form("hSplitE_MC%s",pname[i].Data()),
1985 Form("cluster from %s, energy sum of split sub-clusters",ptype[i].Data()),
1986 nptbins,ptmin,ptmax);
1987 fhMCSplitE[i]->SetYTitle("counts");
1988 fhMCSplitE[i]->SetXTitle("E (GeV)");
1989 outputContainer->Add(fhMCSplitE[i]) ;
85c4406e 1990
cfdf2b91 1991 fhMCSplitPt[i] = new TH1F
1992 (Form("hSplitPt_MC%s",pname[i].Data()),
1993 Form("cluster from %s, pT sum of split sub-clusters",ptype[i].Data()),
1994 nptbins,ptmin,ptmax);
1995 fhMCSplitPt[i]->SetYTitle("counts");
1996 fhMCSplitPt[i]->SetXTitle("p_{T} (GeV/c)");
1997 outputContainer->Add(fhMCSplitPt[i]) ;
1998
29250849 1999
2000 fhMCSplitPtPhi[i] = new TH2F
2001 (Form("hSplitPtPhi_MC%s",pname[i].Data()),
2002 Form("Identified as #pi^{0} (#eta), cluster from %s",ptype[i].Data()),
2003 nptbins,ptmin,ptmax,nphibins,phimin,phimax);
2004 fhMCSplitPtPhi[i]->SetYTitle("#phi");
2005 fhMCSplitPtPhi[i]->SetXTitle("p_{T} (GeV/c)");
2006 outputContainer->Add(fhMCSplitPtPhi[i]) ;
2007
2008 fhMCSplitPtEta[i] = new TH2F
2009 (Form("hSplitPtEta_MC%s",pname[i].Data()),
2010 Form("Identified as #pi^{0} (#eta), cluster from %s",
2011 ptype[i].Data()),nptbins,ptmin,ptmax,netabins,etamin,etamax);
2012 fhMCSplitPtEta[i]->SetYTitle("#eta");
2013 fhMCSplitPtEta[i]->SetXTitle("p_{T} (GeV/c)");
2014 outputContainer->Add(fhMCSplitPtEta[i]) ;
85c4406e 2015
29250849 2016
6e66993c 2017 fhMCNLocMaxSplitPt[i] = new TH2F
2018 (Form("hNLocMaxSplitPt_MC%s",pname[i].Data()),
2019 Form("cluster from %s, pT sum of split sub-clusters, for NLM",ptype[i].Data()),
2020 nptbins,ptmin,ptmax,10,0,10);
2021 fhMCNLocMaxSplitPt[i] ->SetYTitle("N maxima");
2022 fhMCNLocMaxSplitPt[i] ->SetXTitle("p_{T} (GeV/c)");
2023 outputContainer->Add(fhMCNLocMaxSplitPt[i]) ;
2024
29250849 2025 fhMCMassSplitPt[i] = new TH2F
2026 (Form("hMassSplitPt_MC%s",pname[i].Data()),
2027 Form("all pairs mass: split p_{T} vs mass from %s",ptype[i].Data()),
2028 nptbins,ptmin,ptmax, nmassbins,massmin,massmax);
2029 fhMCMassSplitPt[i]->SetYTitle("mass (GeV/c^{2})");
2030 fhMCMassSplitPt[i]->SetXTitle("p_{T} (GeV/c)");
2031 outputContainer->Add(fhMCMassSplitPt[i]) ;
2032
2033 fhMCSelectedMassSplitPt[i] = new TH2F
2034 (Form("hSelectedMassSplitPt_MC%s",pname[i].Data()),
2035 Form("Selected #pi^{0} (#eta) pairs mass: split p_{T} vs mass from %s",ptype[i].Data()),
2036 nptbins,ptmin,ptmax, nmassbins,massmin,massmax);
2037 fhMCSelectedMassSplitPt[i]->SetYTitle("mass (GeV/c^{2})");
2038 fhMCSelectedMassSplitPt[i]->SetXTitle("p_{T} (GeV/c)");
2039 outputContainer->Add(fhMCSelectedMassSplitPt[i]) ;
1253480f 2040
2041 fhMCMassSplitPtNoOverlap[i] = new TH2F
2042 (Form("hMassSplitPtNoOverlap_MC%s",pname[i].Data()),
2043 Form("all pairs mass: split p_{T} vs mass from %s, no overlap",ptype[i].Data()),
2044 nptbins,ptmin,ptmax, nmassbins,massmin,massmax);
2045 fhMCMassSplitPtNoOverlap[i]->SetYTitle("mass (GeV/c^{2})");
2046 fhMCMassSplitPtNoOverlap[i]->SetXTitle("p_{T} (GeV/c)");
2047 outputContainer->Add(fhMCMassSplitPtNoOverlap[i]) ;
29250849 2048
1253480f 2049 fhMCSelectedMassSplitPtNoOverlap[i] = new TH2F
2050 (Form("hSelectedMassSplitPtNoOverlap_MC%s",pname[i].Data()),
2051 Form("Selected #pi^{0} (#eta) pairs mass: split p_{T} vs mass from %s, no overlap",ptype[i].Data()),
2052 nptbins,ptmin,ptmax, nmassbins,massmin,massmax);
2053 fhMCSelectedMassSplitPtNoOverlap[i]->SetYTitle("mass (GeV/c^{2})");
2054 fhMCSelectedMassSplitPtNoOverlap[i]->SetXTitle("p_{T} (GeV/c)");
2055 outputContainer->Add(fhMCSelectedMassSplitPtNoOverlap[i]) ;
85c4406e 2056 }
4650f5cf 2057 }
2058 }
477d6cee 2059
764ab1f4 2060 if(fAnaType==kSSCalo && fFillSelectClHisto && !fFillOnlySimpleSSHisto )
bfdcf7fb 2061 {
2062
bfdcf7fb 2063
2064 for(Int_t i = 0; i< 3; i++)
2065 {
2066 fhEAsymmetryLocMax[i] = new TH2F(Form("hEAsymmetryLocMax%d",i+1),
2067 Form("Selected #pi^{0} (#eta) pairs: E vs A = ( E1 - E2 ) / ( E1 + E2 ), %s",nlm[i].Data()),
85c4406e 2068 nptbins,ptmin,ptmax,200, -1,1);
bfdcf7fb 2069 fhEAsymmetryLocMax[i]->SetYTitle("A = ( E1 - E2 ) / ( E1 + E2 )");
2070 fhEAsymmetryLocMax[i]->SetXTitle("E (GeV)");
2071 outputContainer->Add(fhEAsymmetryLocMax[i]) ;
2072 }
2073
d2655d46 2074 for(Int_t ie = 0; ie< 7; ie++)
bfdcf7fb 2075 {
2076
2077 fhAsymmetryLambda0[ie] = new TH2F (Form("hAsymmetryLambda0_EBin%d",ie),
85c4406e 2078 Form("#lambda_{0}^{2} vs A for %d < E < %d GeV",bin[ie],bin[ie+1]),
2079 ssbins,ssmin,ssmax , 200,-1,1);
bfdcf7fb 2080 fhAsymmetryLambda0[ie]->SetXTitle("#lambda_{0}^{2}");
2081 fhAsymmetryLambda0[ie]->SetYTitle("A = ( E1 - E2 ) / ( E1 + E2 )");
85c4406e 2082 outputContainer->Add(fhAsymmetryLambda0[ie]);
bfdcf7fb 2083
2084 fhAsymmetryDispEta[ie] = new TH2F (Form("hAsymmetryDispEta_EBin%d",ie),
85c4406e 2085 Form("#sigma^{2}_{#eta #eta} vs A for %d < E < %d GeV",bin[ie],bin[ie+1]),
2086 ssbins,ssmin,ssmax , 200,-1,1);
bfdcf7fb 2087 fhAsymmetryDispEta[ie]->SetXTitle("#sigma^{2}_{#eta #eta}");
2088 fhAsymmetryDispEta[ie]->SetYTitle("A = ( E1 - E2 ) / ( E1 + E2 )");
85c4406e 2089 outputContainer->Add(fhAsymmetryDispEta[ie]);
bfdcf7fb 2090
2091 fhAsymmetryDispPhi[ie] = new TH2F (Form("hAsymmetryDispPhi_EBin%d",ie),
85c4406e 2092 Form("#sigma^{2}_{#phi #phi} vs A for %d < E < %d GeV",bin[ie],bin[ie+1]),
2093 ssbins,ssmin,ssmax , 200,-1,1);
bfdcf7fb 2094 fhAsymmetryDispPhi[ie]->SetXTitle("#sigma^{2}_{#phi #phi}");
2095 fhAsymmetryDispPhi[ie]->SetYTitle("A = ( E1 - E2 ) / ( E1 + E2 )");
85c4406e 2096 outputContainer->Add(fhAsymmetryDispPhi[ie]);
2097 }
bfdcf7fb 2098
2099
85c4406e 2100 if(IsDataMC())
bfdcf7fb 2101 {
2102 for(Int_t i = 0; i< 6; i++)
2103 {
d2655d46 2104 for(Int_t ie = 0; ie < 7; ie++)
bfdcf7fb 2105 {
2106 fhMCAsymmetryLambda0[ie][i] = new TH2F (Form("hMCAsymmetryLambda0_EBin%d_MC%s",ie,pname[i].Data()),
85c4406e 2107 Form("cluster from %s : #lambda_{0}^{2} vs A for %d < E < %d GeV",pname[i].Data(),bin[ie],bin[ie+1]),
2108 ssbins,ssmin,ssmax , 200,-1,1);
bfdcf7fb 2109 fhMCAsymmetryLambda0[ie][i]->SetXTitle("#lambda_{0}^{2}");
2110 fhMCAsymmetryLambda0[ie][i]->SetYTitle("A = ( E1 - E2 ) / ( E1 + E2 )");
85c4406e 2111 outputContainer->Add(fhMCAsymmetryLambda0[ie][i]);
bfdcf7fb 2112
2113 fhMCAsymmetryDispEta[ie][i] = new TH2F (Form("hMCAsymmetryDispEta_EBin%d_MC%s",ie,pname[i].Data()),
85c4406e 2114 Form("cluster from %s : #sigma^{2}_{#eta #eta} vs A for %d < E < %d GeV",pname[i].Data(),bin[ie],bin[ie+1]),
2115 ssbins,ssmin,ssmax , 200,-1,1);
bfdcf7fb 2116 fhMCAsymmetryDispEta[ie][i]->SetXTitle("#sigma^{2}_{#eta #eta}");
2117 fhMCAsymmetryDispEta[ie][i]->SetYTitle("A = ( E1 - E2 ) / ( E1 + E2 )");
85c4406e 2118 outputContainer->Add(fhMCAsymmetryDispEta[ie][i]);
bfdcf7fb 2119
2120 fhMCAsymmetryDispPhi[ie][i] = new TH2F (Form("hMCAsymmetryDispPhi_EBin%d_MC%s",ie,pname[i].Data()),
85c4406e 2121 Form("cluster from %s : #sigma^{2}_{#phi #phi} vs A for %d < E < %d GeV",pname[i].Data(),bin[ie],bin[ie+1]),
2122 ssbins,ssmin,ssmax , 200,-1,1);
bfdcf7fb 2123 fhMCAsymmetryDispPhi[ie][i]->SetXTitle("#sigma^{2}_{#phi #phi}");
2124 fhMCAsymmetryDispPhi[ie][i]->SetYTitle("A = ( E1 - E2 ) / ( E1 + E2 )");
85c4406e 2125 outputContainer->Add(fhMCAsymmetryDispPhi[ie][i]);
2126 }
bfdcf7fb 2127 }
2128 }
bfdcf7fb 2129 }
2130
2ad19c3d 2131 if(fFillPileUpHistograms)
2132 {
5e5e056f 2133
2134 TString pileUpName[] = {"SPD","EMCAL","SPDOrEMCAL","SPDAndEMCAL","SPDAndNotEMCAL","EMCALAndNotSPD","NotSPDAndNotEMCAL"} ;
85c4406e 2135
5e5e056f 2136 for(Int_t i = 0 ; i < 7 ; i++)
2137 {
126b8c62 2138 fhPtPileUp[i] = new TH1F(Form("hPtPileUp%s",pileUpName[i].Data()),
85c4406e 2139 Form("Selected #pi^{0} (#eta) p_{T} distribution, %s Pile-Up event",pileUpName[i].Data()), nptbins,ptmin,ptmax);
126b8c62 2140 fhPtPileUp[i]->SetXTitle("p_{T} (GeV/c)");
2141 outputContainer->Add(fhPtPileUp[i]);
2142
2143 fhPtCellTimePileUp[i] = new TH2F(Form("hPtCellTimePileUp%s",pileUpName[i].Data()),
2144 Form("Pt vs cell time in cluster, %s Pile-Up event",pileUpName[i].Data()),
2145 nptbins,ptmin,ptmax,ntimebins,timemin,timemax);
2146 fhPtCellTimePileUp[i]->SetXTitle("p_{T} (GeV/c)");
2147 fhPtCellTimePileUp[i]->SetYTitle("t_{cell} (ns)");
2148 outputContainer->Add(fhPtCellTimePileUp[i]);
2149
2150 fhPtTimeDiffPileUp[i] = new TH2F(Form("hPtTimeDiffPileUp%s",pileUpName[i].Data()),
2151 Form("Pt vs t_{max}-t_{cell} in cluster, %s Pile-Up event",pileUpName[i].Data()),
e31d67f1 2152 nptbins,ptmin,ptmax,400,-200,200);
126b8c62 2153 fhPtTimeDiffPileUp[i]->SetXTitle("p_{T} (GeV/c");
2154 fhPtTimeDiffPileUp[i]->SetYTitle("t_{max}-t_{cell} (ns)");
2155 outputContainer->Add(fhPtTimeDiffPileUp[i]);
2156
5e5e056f 2157 }
2158
126b8c62 2159 fhTimePtNoCut = new TH2F ("hTimePt_NoCut","time of cluster vs E of clusters, no cut", nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
2160 fhTimePtNoCut->SetXTitle("p_{T} (GeV/c)");
2161 fhTimePtNoCut->SetYTitle("time (ns)");
2162 outputContainer->Add(fhTimePtNoCut);
2ad19c3d 2163
126b8c62 2164 fhTimePtSPD = new TH2F ("hTimePt_SPD","time of cluster vs E of clusters, SPD cut", nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
2165 fhTimePtSPD->SetXTitle("p_{T} (GeV/c)");
2166 fhTimePtSPD->SetYTitle("time (ns)");
2167 outputContainer->Add(fhTimePtSPD);
2ad19c3d 2168
126b8c62 2169 fhTimePtSPDMulti = new TH2F ("hTimePt_SPDMulti","time of cluster vs E of clusters, SPD multi cut", nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
2170 fhTimePtSPDMulti->SetXTitle("p_{T} (GeV/c)");
2171 fhTimePtSPDMulti->SetYTitle("time (ns)");
2172 outputContainer->Add(fhTimePtSPDMulti);
2ad19c3d 2173
85c4406e 2174 fhTimeNPileUpVertSPD = new TH2F ("hTime_NPileUpVertSPD","time of cluster vs N pile-up SPD vertex", ntimebins,timemin,timemax,50,0,50);
2ad19c3d 2175 fhTimeNPileUpVertSPD->SetYTitle("# vertex ");
2176 fhTimeNPileUpVertSPD->SetXTitle("time (ns)");
85c4406e 2177 outputContainer->Add(fhTimeNPileUpVertSPD);
2ad19c3d 2178
85c4406e 2179 fhTimeNPileUpVertTrack = new TH2F ("hTime_NPileUpVertTracks","time of cluster vs N pile-up Tracks vertex", ntimebins,timemin,timemax, 50,0,50 );
2ad19c3d 2180 fhTimeNPileUpVertTrack->SetYTitle("# vertex ");
2181 fhTimeNPileUpVertTrack->SetXTitle("time (ns)");
85c4406e 2182 outputContainer->Add(fhTimeNPileUpVertTrack);
2ad19c3d 2183
85c4406e 2184 fhTimeNPileUpVertContributors = new TH2F ("hTime_NPileUpVertContributors","time of cluster vs N constributors to pile-up SPD vertex", ntimebins,timemin,timemax,50,0,50);
2ad19c3d 2185 fhTimeNPileUpVertContributors->SetYTitle("# vertex ");
2186 fhTimeNPileUpVertContributors->SetXTitle("time (ns)");
85c4406e 2187 outputContainer->Add(fhTimeNPileUpVertContributors);
2ad19c3d 2188
85c4406e 2189 fhTimePileUpMainVertexZDistance = new TH2F ("hTime_PileUpMainVertexZDistance","time of cluster vs distance in Z pile-up SPD vertex - main SPD vertex",ntimebins,timemin,timemax,100,0,50);
2ad19c3d 2190 fhTimePileUpMainVertexZDistance->SetYTitle("distance Z (cm) ");
2191 fhTimePileUpMainVertexZDistance->SetXTitle("time (ns)");
85c4406e 2192 outputContainer->Add(fhTimePileUpMainVertexZDistance);
2ad19c3d 2193
85c4406e 2194 fhTimePileUpMainVertexZDiamond = new TH2F ("hTime_PileUpMainVertexZDiamond","time of cluster vs distance in Z pile-up SPD vertex - z diamond",ntimebins,timemin,timemax,100,0,50);
2ad19c3d 2195 fhTimePileUpMainVertexZDiamond->SetYTitle("diamond distance Z (cm) ");
2196 fhTimePileUpMainVertexZDiamond->SetXTitle("time (ns)");
85c4406e 2197 outputContainer->Add(fhTimePileUpMainVertexZDiamond);
0f7e7205 2198
85c4406e 2199 fhPtNPileUpSPDVtx = new TH2F ("hPt_NPileUpVertSPD","pT of cluster vs N pile-up SPD vertex",
2200 nptbins,ptmin,ptmax,20,0,20);
0f7e7205 2201 fhPtNPileUpSPDVtx->SetYTitle("# vertex ");
2202 fhPtNPileUpSPDVtx->SetXTitle("p_{T} (GeV/c)");
2203 outputContainer->Add(fhPtNPileUpSPDVtx);
2204
85c4406e 2205 fhPtNPileUpTrkVtx = new TH2F ("hPt_NPileUpVertTracks","pT of cluster vs N pile-up Tracks vertex",
2206 nptbins,ptmin,ptmax, 20,0,20 );
0f7e7205 2207 fhPtNPileUpTrkVtx->SetYTitle("# vertex ");
2208 fhPtNPileUpTrkVtx->SetXTitle("p_{T} (GeV/c)");
85c4406e 2209 outputContainer->Add(fhPtNPileUpTrkVtx);
0f7e7205 2210
85c4406e 2211 fhPtNPileUpSPDVtxTimeCut = new TH2F ("hPt_NPileUpVertSPD_TimeCut","pT of cluster vs N pile-up SPD vertex, |tof| < 25 ns",
2212 nptbins,ptmin,ptmax,20,0,20);
0f7e7205 2213 fhPtNPileUpSPDVtxTimeCut->SetYTitle("# vertex ");
2214 fhPtNPileUpSPDVtxTimeCut->SetXTitle("p_{T} (GeV/c)");
2215 outputContainer->Add(fhPtNPileUpSPDVtxTimeCut);
2216
85c4406e 2217 fhPtNPileUpTrkVtxTimeCut = new TH2F ("hPt_NPileUpVertTracks_TimeCut","pT of cluster vs N pile-up Tracks vertex, |tof| < 25 ns",
2218 nptbins,ptmin,ptmax, 20,0,20 );
0f7e7205 2219 fhPtNPileUpTrkVtxTimeCut->SetYTitle("# vertex ");
2220 fhPtNPileUpTrkVtxTimeCut->SetXTitle("p_{T} (GeV/c)");
85c4406e 2221 outputContainer->Add(fhPtNPileUpTrkVtxTimeCut);
2222
2223 fhPtNPileUpSPDVtxTimeCut2 = new TH2F ("hPt_NPileUpVertSPD_TimeCut2","pT of cluster vs N pile-up SPD vertex, -25 < tof < 75 ns",
2224 nptbins,ptmin,ptmax,20,0,20);
2225 fhPtNPileUpSPDVtxTimeCut2->SetYTitle("# vertex ");
2226 fhPtNPileUpSPDVtxTimeCut2->SetXTitle("p_{T} (GeV/c)");
2227 outputContainer->Add(fhPtNPileUpSPDVtxTimeCut2);
2228
2229 fhPtNPileUpTrkVtxTimeCut2 = new TH2F ("hPt_NPileUpVertTracks_TimeCut2","pT of cluster vs N pile-up Tracks vertex, -25 < tof < 75 ns",
2230 nptbins,ptmin,ptmax, 20,0,20 );
2231 fhPtNPileUpTrkVtxTimeCut2->SetYTitle("# vertex ");
2232 fhPtNPileUpTrkVtxTimeCut2->SetXTitle("p_{T} (GeV/c)");
2233 outputContainer->Add(fhPtNPileUpTrkVtxTimeCut2);
2234
2ad19c3d 2235 }
2236
477d6cee 2237 //Keep neutral meson selection histograms if requiered
2238 //Setting done in AliNeutralMesonSelection
2239
e4ef72be 2240 if(fAnaType!=kSSCalo && GetNeutralMesonSelection())
2241 {
477d6cee 2242 TList * nmsHistos = GetNeutralMesonSelection()->GetCreateOutputObjects() ;
e4ef72be 2243
477d6cee 2244 if(GetNeutralMesonSelection()->AreNeutralMesonSelectionHistosKept())
2245 for(Int_t i = 0; i < nmsHistos->GetEntries(); i++) outputContainer->Add(nmsHistos->At(i)) ;
e4ef72be 2246
5ae09196 2247 delete nmsHistos;
477d6cee 2248 }
2249
477d6cee 2250 return outputContainer ;
2251
2252}
2253
3455f821 2254//_____________________________________________
2255Int_t AliAnaPi0EbE::GetMCIndex(const Int_t tag)
85c4406e 2256{
3455f821 2257
2258 // Assign mc index depending on MC bit set
2259
2260 if ( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0) )
2261 {
85c4406e 2262 return kmcPi0 ;
3455f821 2263 }//pi0
2264 else if ( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEta) )
2265 {
85c4406e 2266 return kmcEta ;
2267 }//eta
3455f821 2268 else if ( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPhoton) &&
1253480f 2269 GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCConversion) )
3455f821 2270 {
85c4406e 2271 return kmcConversion ;
3455f821 2272 }//conversion photon
2273 else if ( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPhoton) )
2274 {
85c4406e 2275 return kmcPhoton ;
3455f821 2276 }//photon no conversion
2277 else if ( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCElectron))
2278 {
85c4406e 2279 return kmcElectron ;
3455f821 2280 }//electron
85c4406e 2281 else
3455f821 2282 {
85c4406e 2283 return kmcHadron ;
2284 }//other particles
3455f821 2285
2286}
2287
2288//__________________________________________________________________
85c4406e 2289void AliAnaPi0EbE::HasPairSameMCMother(AliAODPWG4Particle * photon1,
2290 AliAODPWG4Particle * photon2,
3455f821 2291 Int_t & label, Int_t & tag)
2292{
2293 // Check the labels of pare in case mother was same pi0 or eta
2294 // Set the new AOD accordingly
2295
2296 Int_t label1 = photon1->GetLabel();
2297 Int_t label2 = photon2->GetLabel();
2298
2299 if(label1 < 0 || label2 < 0 ) return ;
2300
2644ead9 2301 //Int_t tag1 = GetMCAnalysisUtils()->CheckOrigin(label1, GetReader());
2302 //Int_t tag2 = GetMCAnalysisUtils()->CheckOrigin(label2, GetReader());
3455f821 2303 Int_t tag1 = photon1->GetTag();
2304 Int_t tag2 = photon2->GetTag();
2305
2306 if(GetDebug() > 0) printf("AliAnaPi0EbE::MakeInvMassInCalorimeter() - Origin of: photon1 %d; photon2 %d \n",tag1, tag2);
85c4406e 2307 if( (GetMCAnalysisUtils()->CheckTagBit(tag1,AliMCAnalysisUtils::kMCPi0Decay) &&
3455f821 2308 GetMCAnalysisUtils()->CheckTagBit(tag2,AliMCAnalysisUtils::kMCPi0Decay) ) ||
85c4406e 2309 (GetMCAnalysisUtils()->CheckTagBit(tag1,AliMCAnalysisUtils::kMCEtaDecay) &&
3455f821 2310 GetMCAnalysisUtils()->CheckTagBit(tag2,AliMCAnalysisUtils::kMCEtaDecay) )
2311 )
2312 {
2313
2314 //Check if pi0/eta mother is the same
2315 if(GetReader()->ReadStack())
85c4406e 2316 {
3455f821 2317 if(label1>=0)
2318 {
2319 TParticle * mother1 = GetMCStack()->Particle(label1);//photon in kine tree
2320 label1 = mother1->GetFirstMother();
2321 //mother1 = GetMCStack()->Particle(label1);//pi0
2322 }
2323 if(label2>=0)
2324 {
2325 TParticle * mother2 = GetMCStack()->Particle(label2);//photon in kine tree
2326 label2 = mother2->GetFirstMother();
2327 //mother2 = GetMCStack()->Particle(label2);//pi0
2328 }
2329 } // STACK
2330 else if(GetReader()->ReadAODMCParticles())
2331 {//&& (input > -1)){
2332 if(label1>=0)
2333 {
2644ead9 2334 AliAODMCParticle * mother1 = (AliAODMCParticle *) (GetReader()->GetAODMCParticles())->At(label1);//photon in kine tree
3455f821 2335 label1 = mother1->GetMother();
2336 //mother1 = GetMCStack()->Particle(label1);//pi0
2337 }
2338 if(label2>=0)
2339 {
2644ead9 2340 AliAODMCParticle * mother2 = (AliAODMCParticle *) (GetReader()->GetAODMCParticles())->At(label2);//photon in kine tree
3455f821 2341 label2 = mother2->GetMother();
2342 //mother2 = GetMCStack()->Particle(label2);//pi0
2343 }
2344 }// AOD
2345
2346 //printf("mother1 %d, mother2 %d\n",label1,label2);
2347 if( label1 == label2 && label1>=0 )
2348 {
2349
2350 label = label1;
2351
2352 TLorentzVector mom1 = *(photon1->Momentum());
2353 TLorentzVector mom2 = *(photon2->Momentum());
2354
2355 Double_t angle = mom2.Angle(mom1.Vect());
2356 Double_t mass = (mom1+mom2).M();
2357 Double_t epair = (mom1+mom2).E();
2358
2359 if(GetMCAnalysisUtils()->CheckTagBit(tag1,AliMCAnalysisUtils::kMCPi0Decay))
2360 {
2361 fhMassPairMCPi0 ->Fill(epair,mass);
2362 fhAnglePairMCPi0->Fill(epair,angle);
2363 GetMCAnalysisUtils()->SetTagBit(tag,AliMCAnalysisUtils::kMCPi0);
2364 }
85c4406e 2365 else
3455f821 2366 {
2367 fhMassPairMCEta ->Fill(epair,mass);
2368 fhAnglePairMCEta->Fill(epair,angle);
2369 GetMCAnalysisUtils()->SetTagBit(tag,AliMCAnalysisUtils::kMCEta);
2370 }
2371
2372 } // same label
2373 } // both from eta or pi0 decay
2374
85c4406e 2375}
3455f821 2376
521636d2 2377//____________________________________________________________________________
2378void AliAnaPi0EbE::Init()
85c4406e 2379{
521636d2 2380 //Init
2381 //Do some checks
2382 if(fCalorimeter == "PHOS" && !GetReader()->IsPHOSSwitchedOn() && NewOutputAOD()){
2383 printf("AliAnaPi0EbE::Init() - !!STOP: You want to use PHOS in analysis but it is not read!! \n!!Check the configuration file!!\n");
2384 abort();
2385 }
2386 else if(fCalorimeter == "EMCAL" && !GetReader()->IsEMCALSwitchedOn() && NewOutputAOD()){
2387 printf("AliAnaPi0EbE::Init() - !!STOP: You want to use EMCAL in analysis but it is not read!! \n!!Check the configuration file!!\n");
2388 abort();
2389 }
2390
2391}
2392
2393//____________________________________________________________________________
2394void AliAnaPi0EbE::InitParameters()
2395{
85c4406e 2396 //Initialize the parameters of the analysis.
521636d2 2397 AddToHistogramsName("AnaPi0EbE_");
2398
1db06135 2399 fInputAODGammaConvName = "PhotonsCTS" ;
521636d2 2400 fAnaType = kIMCalo ;
2401 fCalorimeter = "EMCAL" ;
2402 fMinDist = 2.;
2403 fMinDist2 = 4.;
2404 fMinDist3 = 5.;
2405
4d97a954 2406 fNLMECutMin[0] = 10.;
2407 fNLMECutMin[1] = 6. ;
2408 fNLMECutMin[2] = 6. ;
85c4406e 2409
521636d2 2410}
2411
477d6cee 2412//__________________________________________________________________
85c4406e 2413void AliAnaPi0EbE::MakeAnalysisFillAOD()
477d6cee 2414{
2415 //Do analysis and fill aods
2416
85c4406e 2417 switch(fAnaType)
521636d2 2418 {
477d6cee 2419 case kIMCalo:
2420 MakeInvMassInCalorimeter();
2421 break;
2422
2423 case kSSCalo:
2424 MakeShowerShapeIdentification();
2425 break;
2426
2427 case kIMCaloTracks:
2428 MakeInvMassInCalorimeterAndCTS();
2429 break;
2430
521636d2 2431 }
477d6cee 2432}
2433
42d47cb7 2434//____________________________________________
85c4406e 2435void AliAnaPi0EbE::MakeInvMassInCalorimeter()
477d6cee 2436{
57b97dc6 2437 //Do analysis and fill aods
2438 //Search for the photon decay in calorimeters
2439 //Read photon list from AOD, produced in class AliAnaPhoton
2440 //Check if 2 photons have the mass of the pi0.
477d6cee 2441
2442 TLorentzVector mom1;
2443 TLorentzVector mom2;
2444 TLorentzVector mom ;
85c4406e 2445
b5dbb99b 2446 Int_t tag = 0;
2447 Int_t label = 0;
477d6cee 2448
2449 if(!GetInputAODBranch()){
a3aebfff 2450 printf("AliAnaPi0EbE::MakeInvMassInCalorimeter() - No input calo photons in AOD with name branch < %s >, STOP \n",GetInputAODName().Data());
477d6cee 2451 abort();
2452 }
f8006433 2453
42d47cb7 2454 //Get shower shape information of clusters
2455 TObjArray *clusters = 0;
2456 if (fCalorimeter=="EMCAL") clusters = GetEMCALClusters();
2457 else if(fCalorimeter=="PHOS") clusters = GetPHOSClusters() ;
2458
c4a7d28a 2459 for(Int_t iphoton = 0; iphoton < GetInputAODBranch()->GetEntriesFast()-1; iphoton++){
477d6cee 2460 AliAODPWG4Particle * photon1 = (AliAODPWG4Particle*) (GetInputAODBranch()->At(iphoton));
c8fe2783 2461
c4a7d28a 2462 //Vertex cut in case of mixed events
85c4406e 2463 Int_t evtIndex1 = 0 ;
c8fe2783 2464 if(GetMixedEvent())
2465 evtIndex1 = GetMixedEvent()->EventIndexForCaloCluster(photon1->GetCaloLabel(0)) ;
5025c139 2466 if(TMath::Abs(GetVertex(evtIndex1)[2]) > GetZvertexCut()) continue ; //vertex cut
477d6cee 2467 mom1 = *(photon1->Momentum());
2468
42d47cb7 2469 //Get original cluster, to recover some information
1db06135 2470 Int_t iclus = -1;
85c4406e 2471 AliVCluster *cluster1 = FindCluster(clusters,photon1->GetCaloLabel(0),iclus);
42d47cb7 2472
1db06135 2473 if(!cluster1){
42d47cb7 2474 printf("AliAnaPi0EbE::MakeInvMassInCalorimeter() - First cluster not found\n");
2475 return;
9ab9e937 2476 }
c4a7d28a 2477
b5dbb99b 2478 for(Int_t jphoton = iphoton+1; jphoton < GetInputAODBranch()->GetEntriesFast(); jphoton++)
2479 {
a3aebfff 2480 AliAODPWG4Particle * photon2 = (AliAODPWG4Particle*) (GetInputAODBranch()->At(jphoton));
b5dbb99b 2481
85c4406e 2482 Int_t evtIndex2 = 0 ;
c8fe2783 2483 if(GetMixedEvent())
2484 evtIndex2 = GetMixedEvent()->EventIndexForCaloCluster(photon2->GetCaloLabel(0)) ;
b5dbb99b 2485
c8fe2783 2486 if(GetMixedEvent() && (evtIndex1 == evtIndex2))
85c4406e 2487 continue ;
b5dbb99b 2488
5025c139 2489 if(TMath::Abs(GetVertex(evtIndex2)[2]) > GetZvertexCut()) continue ; //vertex cut
b5dbb99b 2490
477d6cee 2491 mom2 = *(photon2->Momentum());
c4a7d28a 2492
1db06135 2493 //Get original cluster, to recover some information
2494 Int_t iclus2;
85c4406e 2495 AliVCluster *cluster2 = FindCluster(clusters,photon2->GetCaloLabel(0),iclus2,iclus+1);
42d47cb7 2496
b5dbb99b 2497 if(!cluster2)
2498 {
42d47cb7 2499 printf("AliAnaPi0EbE::MakeInvMassInCalorimeter() - Second cluster not found\n");
1db06135 2500 return;
9ab9e937 2501 }
c4a7d28a 2502
85c4406e 2503 Float_t e1 = photon1->E();
42d47cb7 2504 Float_t e2 = photon2->E();
2505
2506 //Select clusters with good time window difference
2507 Float_t tof1 = cluster1->GetTOF()*1e9;;
2508 Float_t tof2 = cluster2->GetTOF()*1e9;;
2509 Double_t t12diff = tof1-tof2;
2510 fhEPairDiffTime->Fill(e1+e2, t12diff);
2511 if(TMath::Abs(t12diff) > GetPairTimeCut()) continue;
2512
b5dbb99b 2513 //Play with the MC stack if available
2514 if(IsDataMC()) HasPairSameMCMother(photon1, photon2, label, tag) ;
85c4406e 2515
5c46c992 2516 // Check the invariant mass for different selection on the local maxima
2517 // Name of AOD method TO BE FIXED
2518 Int_t nMaxima1 = photon1->GetFiducialArea();
2519 Int_t nMaxima2 = photon2->GetFiducialArea();
2520
2521 Double_t mass = (mom1+mom2).M();
2522 Double_t epair = (mom1+mom2).E();
2523
2524 if(nMaxima1==nMaxima2)
2525 {
2526 if (nMaxima1==1) fhMassPairLocMax[0]->Fill(epair,mass);
2527 else if(nMaxima1==2) fhMassPairLocMax[1]->Fill(epair,mass);
2528 else fhMassPairLocMax[2]->Fill(epair,mass);
2529 }
2530 else if(nMaxima1==1 || nMaxima2==1)
2531 {
2532 if (nMaxima1==2 || nMaxima2==2) fhMassPairLocMax[3]->Fill(epair,mass);
85c4406e 2533 else fhMassPairLocMax[4]->Fill(epair,mass);
5c46c992 2534 }
85c4406e 2535 else
5c46c992 2536 fhMassPairLocMax[5]->Fill(epair,mass);
2537
3c1d9afb 2538 // combinations with SS axis cut and NLM cut
85c4406e 2539 if(nMaxima1 == 1 && cluster2->GetM02() > 0.3) fhMassPairLocMax[6]->Fill(epair,mass);
2540 if(nMaxima2 == 1 && cluster1->GetM02() > 0.3) fhMassPairLocMax[6]->Fill(epair,mass);
2541 if(nMaxima1 > 1 && cluster2->GetM02() < 0.3 && cluster2->GetM02()> 0.1 ) fhMassPairLocMax[7]->Fill(epair,mass);
2542 if(nMaxima2 > 1 && cluster1->GetM02() < 0.3 && cluster1->GetM02()> 0.1 ) fhMassPairLocMax[7]->Fill(epair,mass);
3c1d9afb 2543
a6e83e39 2544 //Skip events with too few or too many NLM
fb51265c 2545 if((nMaxima1 < fNLMCutMin || nMaxima1 > fNLMCutMax) || (nMaxima2 < fNLMCutMin || nMaxima2 > fNLMCutMax)) continue ;
a6e83e39 2546
2547 if(GetDebug() > 1) printf("AliAnaPi0EbE::MakeInvMassInCalorimeter() - NLM of out of range: cluster1 %d, cluster2 %d \n",nMaxima1, nMaxima2);
2548
f02db2c0 2549 //Mass of all pairs
2550 fhMass->Fill(epair,(mom1+mom2).M());
a6e83e39 2551
57b97dc6 2552 //Select good pair (good phi, pt cuts, aperture and invariant mass)
3bfcb597 2553 if(GetNeutralMesonSelection()->SelectPair(mom1, mom2,fCalorimeter))
c8fe2783 2554 {
85c4406e 2555 if(GetDebug()>1)
c8fe2783 2556 printf("AliAnaPi0EbE::MakeInvMassInCalorimeter() - Selected gamma pair: pt %f, phi %f, eta%f \n",(mom1+mom2).Pt(), (mom1+mom2).Phi()*180./3.1416, (mom1+mom2).Eta());
42d47cb7 2557
57b97dc6 2558 //Fill some histograms about shower shape
06e81356 2559 if(fFillSelectClHisto && clusters && GetReader()->GetDataType()!=AliCaloTrackReader::kMC)
5c46c992 2560 {
2561 FillSelectedClusterHistograms(cluster1, nMaxima1, photon1->GetTag());
2562 FillSelectedClusterHistograms(cluster2, nMaxima2, photon2->GetTag());
42d47cb7 2563 }
521636d2 2564
803d06a8 2565 // Tag both photons as decay
2566 photon1->SetTagged(kTRUE);
2567 photon2->SetTagged(kTRUE);
09273901 2568
2569 fhPtDecay->Fill(photon1->Pt());
2570 fhEDecay ->Fill(photon1->E() );
2571
2572 fhPtDecay->Fill(photon2->Pt());
2573 fhEDecay ->Fill(photon2->E() );
2ad19c3d 2574
57b97dc6 2575 //Create AOD for analysis
c8fe2783 2576 mom = mom1+mom2;
85c4406e 2577
f02db2c0 2578 //Mass of selected pairs
2579 fhSelectedMass->Fill(epair,mom.M());
2580
2ad19c3d 2581 // Fill histograms to undertand pile-up before other cuts applied
2582 // Remember to relax time cuts in the reader
126b8c62 2583 FillPileUpHistograms(mom.Pt(),((cluster1->GetTOF()+cluster2->GetTOF())*1e9)/2,cluster1);
2ad19c3d 2584
c8fe2783 2585 AliAODPWG4Particle pi0 = AliAODPWG4Particle(mom);
b5dbb99b 2586
21a4b1c0 2587 pi0.SetIdentifiedParticleType(AliCaloPID::kPi0);
c8fe2783 2588 pi0.SetDetector(photon1->GetDetector());
b5dbb99b 2589
2590 // MC
2591 pi0.SetLabel(label);
85c4406e 2592 pi0.SetTag(tag);
b5dbb99b 2593
85c4406e 2594 //Set the indeces of the original caloclusters
c8fe2783 2595 pi0.SetCaloLabel(photon1->GetCaloLabel(0), photon2->GetCaloLabel(0));
f8006433 2596 //pi0.SetInputFileIndex(input);
b5dbb99b 2597
c8fe2783 2598 AddAODParticle(pi0);
b5dbb99b 2599
c8fe2783 2600 }//pi0
57b97dc6 2601
477d6cee 2602 }//2n photon loop
2603
2604 }//1st photon loop
2605
85c4406e 2606 if(GetDebug() > 1) printf("AliAnaPi0EbE::MakeInvMassInCalorimeter() - End fill AODs \n");
477d6cee 2607
2608}
2609
e7fd282f 2610//__________________________________________________
85c4406e 2611void AliAnaPi0EbE::MakeInvMassInCalorimeterAndCTS()
477d6cee 2612{
2613 //Do analysis and fill aods
2614 //Search for the photon decay in calorimeters
2615 //Read photon list from AOD, produced in class AliAnaPhoton and AliGammaConversion
2616 //Check if 2 photons have the mass of the pi0.
2617
2618 TLorentzVector mom1;
2619 TLorentzVector mom2;
2620 TLorentzVector mom ;
b5dbb99b 2621 Int_t tag = 0;
2622 Int_t label = 0;
5025c139 2623 Int_t evtIndex = 0;
1db06135 2624
2625 // Check calorimeter input
477d6cee 2626 if(!GetInputAODBranch()){
a3aebfff 2627 printf("AliAnaPi0EbE::MakeInvMassInCalorimeterAndCTS() - No input calo photons in AOD branch with name < %s > , STOP\n",GetInputAODName().Data());
477d6cee 2628 abort();
2629 }
57b97dc6 2630
1db06135 2631 // Get the array with conversion photons
2632 TClonesArray * inputAODGammaConv = (TClonesArray *) GetReader()->GetOutputEvent()->FindListObject(fInputAODGammaConvName);
2633 if(!inputAODGammaConv) {
2634
2635 inputAODGammaConv = (TClonesArray *) GetReader()->GetInputEvent()->FindListObject(fInputAODGammaConvName);
2636
2637 if(!inputAODGammaConv) {
2638 printf("AliAnaPi0EbE::MakeInvMassInCalorimeterAndCTS() - No input gamma conversions in AOD branch with name < %s >\n",fInputAODGammaConvName.Data());
2639
2640 return;
2641 }
85c4406e 2642 }
1db06135 2643
2644 //Get shower shape information of clusters
2645 TObjArray *clusters = 0;
2646 if (fCalorimeter=="EMCAL") clusters = GetEMCALClusters();
85c4406e 2647 else if(fCalorimeter=="PHOS") clusters = GetPHOSClusters() ;
1db06135 2648
2649 Int_t nCTS = inputAODGammaConv->GetEntriesFast();
2650 Int_t nCalo = GetInputAODBranch()->GetEntriesFast();
85c4406e 2651 if(nCTS<=0 || nCalo <=0)
a6e83e39 2652 {
1db06135 2653 if(GetDebug() > 1) printf("AliAnaPi0EbE::MakeInvMassInCalorimeterAndCTS() - nCalo %d, nCTS %d, cannot loop\n",nCalo,nCTS);
2654 return;
2655 }
2656
2657 if(GetDebug() > 1)
2658 printf("AliAnaPi0EbE::MakeInvMassInCalorimeterAndCTS() - Number of conversion photons %d\n",nCTS);
2659
2660 // Do the loop, first calo, second CTS
477d6cee 2661 for(Int_t iphoton = 0; iphoton < GetInputAODBranch()->GetEntriesFast(); iphoton++){
2662 AliAODPWG4Particle * photon1 = (AliAODPWG4Particle*) (GetInputAODBranch()->At(iphoton));
2663 mom1 = *(photon1->Momentum());
2664
1db06135 2665 //Get original cluster, to recover some information
2666 Int_t iclus = -1;
85c4406e 2667 AliVCluster *cluster = FindCluster(clusters,photon1->GetCaloLabel(0),iclus);
1db06135 2668
2669 for(Int_t jphoton = 0; jphoton < nCTS; jphoton++){
2670 AliAODPWG4Particle * photon2 = (AliAODPWG4Particle*) (inputAODGammaConv->At(jphoton));
5025c139 2671 if(GetMixedEvent())
2672 evtIndex = GetMixedEvent()->EventIndexForCaloCluster(photon2->GetCaloLabel(0)) ;
2673 if(TMath::Abs(GetVertex(evtIndex)[2]) > GetZvertexCut()) continue ; //vertex cut
2674
477d6cee 2675 mom2 = *(photon2->Momentum());
57b97dc6 2676
5c46c992 2677 Double_t mass = (mom1+mom2).M();
2678 Double_t epair = (mom1+mom2).E();
2679
2680 Int_t nMaxima = photon1->GetFiducialArea();
2681 if (nMaxima==1) fhMassPairLocMax[0]->Fill(epair,mass);
2682 else if(nMaxima==2) fhMassPairLocMax[1]->Fill(epair,mass);
2683 else fhMassPairLocMax[2]->Fill(epair,mass);
2684
a6e83e39 2685 if(nMaxima < fNLMCutMin || nMaxima > fNLMCutMax) continue ;
2686 if(GetDebug() > 1) printf("AliAnaPi0EbE::MakeInvMassInCalorimeterAndCTS() - NLM %d of out of range \n",nMaxima);
85c4406e 2687
b5dbb99b 2688 //Play with the MC stack if available
2689 if(IsDataMC())
2690 {
2691 Int_t label2 = photon2->GetLabel();
2644ead9 2692 if(label2 >= 0 )photon2->SetTag(GetMCAnalysisUtils()->CheckOrigin(label2, GetReader()));
b5dbb99b 2693
2694 HasPairSameMCMother(photon1, photon2, label, tag) ;
2695 }
2696
f02db2c0 2697 //Mass of selected pairs
e671adc2 2698 fhMass->Fill(epair,(mom1+mom2).M());
f02db2c0 2699
477d6cee 2700 //Select good pair (good phi, pt cuts, aperture and invariant mass)
b5dbb99b 2701 if(GetNeutralMesonSelection()->SelectPair(mom1, mom2,fCalorimeter))
2702 {
57b97dc6 2703 if(GetDebug() > 1) printf("AliAnaPi0EbE::MakeInvMassInCalorimeterAndCTS() - Selected gamma pair: pt %f, phi %f, eta%f\n",(mom1+mom2).Pt(), (mom1+mom2).Phi()*180./3.1416, (mom1+mom2).Eta());
2704
1db06135 2705 //Fill some histograms about shower shape
06e81356 2706 if(fFillSelectClHisto && cluster && GetReader()->GetDataType()!=AliCaloTrackReader::kMC)
b5dbb99b 2707 {
5c46c992 2708 FillSelectedClusterHistograms(cluster, nMaxima, photon1->GetTag());
85c4406e 2709 }
803d06a8 2710
2711 // Tag both photons as decay
2712 photon1->SetTagged(kTRUE);
85c4406e 2713 photon2->SetTagged(kTRUE);
1db06135 2714
09273901 2715 fhPtDecay->Fill(photon1->Pt());
2716 fhEDecay ->Fill(photon1->E() );
2717
57b97dc6 2718 //Create AOD for analysis
b5dbb99b 2719
57b97dc6 2720 mom = mom1+mom2;
b5dbb99b 2721
f02db2c0 2722 //Mass of selected pairs
2723 fhSelectedMass->Fill(epair,mom.M());
2724
2ad19c3d 2725 // Fill histograms to undertand pile-up before other cuts applied
2726 // Remember to relax time cuts in the reader
126b8c62 2727 if(cluster) FillPileUpHistograms(mom.Pt(),cluster->GetTOF()*1e9,cluster);
2ad19c3d 2728
57b97dc6 2729 AliAODPWG4Particle pi0 = AliAODPWG4Particle(mom);
b5dbb99b 2730
21a4b1c0 2731 pi0.SetIdentifiedParticleType(AliCaloPID::kPi0);
57b97dc6 2732 pi0.SetDetector(photon1->GetDetector());
b5dbb99b 2733
2734 // MC
2735 pi0.SetLabel(label);
57b97dc6 2736 pi0.SetTag(tag);
b5dbb99b 2737
85c4406e 2738 //Set the indeces of the original tracks or caloclusters
57b97dc6 2739 pi0.SetCaloLabel(photon1->GetCaloLabel(0), -1);
2740 pi0.SetTrackLabel(photon2->GetTrackLabel(0), photon2->GetTrackLabel(1));
f8006433 2741 //pi0.SetInputFileIndex(input);
b5dbb99b 2742
57b97dc6 2743 AddAODParticle(pi0);
b5dbb99b 2744
477d6cee 2745 }//pi0
2746 }//2n photon loop
2747
2748 }//1st photon loop
2749
85c4406e 2750 if(GetDebug() > 1) printf("AliAnaPi0EbE::MakeInvMassInCalorimeterAndCTS() - End fill AODs \n");
477d6cee 2751
2752}
2753
2754
e7fd282f 2755//_________________________________________________
85c4406e 2756void AliAnaPi0EbE::MakeShowerShapeIdentification()
477d6cee 2757{
85c4406e 2758 //Search for pi0 in fCalorimeter with shower shape analysis
477d6cee 2759
85c4406e 2760 TObjArray * pl = 0x0;
34c16486 2761 AliVCaloCells * cells = 0x0;
5ae09196 2762 //Select the Calorimeter of the photon
b5dbb99b 2763 if (fCalorimeter == "PHOS" )
34c16486 2764 {
2765 pl = GetPHOSClusters();
2766 cells = GetPHOSCells();
2767 }
5ae09196 2768 else if (fCalorimeter == "EMCAL")
34c16486 2769 {
2770 pl = GetEMCALClusters();
2771 cells = GetEMCALCells();
2772 }
57b97dc6 2773
85c4406e 2774 if(!pl)
34c16486 2775 {
5ae09196 2776 Info("MakeShowerShapeIdentification","TObjArray with %s clusters is NULL!\n",fCalorimeter.Data());
2777 return;
85c4406e 2778 }
233e0df8 2779
477d6cee 2780 TLorentzVector mom ;
b5dbb99b 2781 for(Int_t icalo = 0; icalo < pl->GetEntriesFast(); icalo++)
2782 {
85c4406e 2783 AliVCluster * calo = (AliVCluster*) (pl->At(icalo));
477d6cee 2784
85c4406e 2785 Int_t evtIndex = 0 ;
2786 if (GetMixedEvent())
b5dbb99b 2787 {
85c4406e 2788 evtIndex=GetMixedEvent()->EventIndexForCaloCluster(calo->GetID()) ;
f8006433 2789 }
34c16486 2790
5025c139 2791 if(TMath::Abs(GetVertex(evtIndex)[2]) > GetZvertexCut()) continue ; //vertex cut
521636d2 2792
85c4406e 2793 //Get Momentum vector,
a6e83e39 2794 Double_t vertex[]={0,0,0};
34c16486 2795 if(GetReader()->GetDataType() != AliCaloTrackReader::kMC)
2796 {
2797 calo->GetMomentum(mom,GetVertex(evtIndex)) ;
2798 }//Assume that come from vertex in straight line
2799 else
2800 {
f8006433 2801 calo->GetMomentum(mom,vertex) ;
2802 }
233e0df8 2803
57b97dc6 2804 //If too small or big pt, skip it
85c4406e 2805 if(mom.E() < GetMinEnergy() || mom.E() > GetMaxEnergy() ) continue ;
34c16486 2806
477d6cee 2807 //Check acceptance selection
b5dbb99b 2808 if(IsFiducialCutOn())
2809 {
ff45398a 2810 Bool_t in = GetFiducialCut()->IsInFiducialCut(mom,fCalorimeter) ;
477d6cee 2811 if(! in ) continue ;
2812 }
2813
85c4406e 2814 if(GetDebug() > 1)
2815 printf("AliAnaPi0EbE::MakeShowerShapeIdentification() - FillAOD: Min pt cut and fiducial cut passed: pt %3.2f, phi %2.2f, eta %1.2f\n",mom.Pt(),mom.Phi(),mom.Eta());
2816
477d6cee 2817 //Check Distance to Bad channel, set bit.
c8fe2783 2818 Double_t distBad=calo->GetDistanceToBadChannel() ; //Distance to bad channel
477d6cee 2819 if(distBad < 0.) distBad=9999. ; //workout strange convension dist = -1. ;
2820 if(distBad < fMinDist) //In bad channel (PHOS cristal size 2.2x2.2 cm)
2821 continue ;
85c4406e 2822
a6e83e39 2823 if(GetDebug() > 1) printf("AliAnaPi0EbE::MakeShowerShapeIdentification() - FillAOD: Bad channel cut passed %4.2f\n",distBad);
85c4406e 2824
74e3eb22 2825 //If too low number of cells, skip it
2826 if ( calo->GetNCells() < GetCaloPID()->GetClusterSplittingMinNCells()) continue ;
2827
2828 if(GetDebug() > 1)
2829 printf("AliAnaPi0EbE::MakeShowerShapeIdentification() - FillAOD: N cells cut passed %d > %d\n",
2830 calo->GetNCells(), GetCaloPID()->GetClusterSplittingMinNCells());
2831
34c16486 2832 //.......................................
2833 // TOF cut, BE CAREFUL WITH THIS CUT
2834 Double_t tof = calo->GetTOF()*1e9;
2835 if(tof < fTimeCutMin || tof > fTimeCutMax) continue ;
2836
40d3ce60 2837 //Play with the MC stack if available
2838 //Check origin of the candidates
2839 Int_t tag = 0 ;
2840 if(IsDataMC())
2841 {
2644ead9 2842 tag = GetMCAnalysisUtils()->CheckOrigin(calo->GetLabels(),calo->GetNLabels(),GetReader());
40d3ce60 2843 //GetMCAnalysisUtils()->CheckMultipleOrigin(calo->GetLabels(),calo->GetNLabels(), GetReader(), aodpi0.GetInputFileIndex(), tag);
2844 if(GetDebug() > 0) printf("AliAnaPi0EbE::MakeShowerShapeIdentification() - Origin of candidate %d\n",tag);
85c4406e 2845 }
40d3ce60 2846
a6e83e39 2847 //Skip matched clusters with tracks
e997f026 2848 if(fRejectTrackMatch && IsTrackMatched(calo, GetReader()->GetInputEvent()))
40d3ce60 2849 {
2850 FillRejectedClusterHistograms(mom,tag);
2851 continue ;
2852 }
477d6cee 2853
2854 //Check PID
2855 //PID selection or bit setting
85c4406e 2856 Int_t nMaxima = 0 ;
34c16486 2857 Double_t mass = 0 , angle = 0;
cfdf2b91 2858 TLorentzVector l1, l2;
4914e781 2859 Int_t absId1 =-1; Int_t absId2 =-1;
2860 Float_t distbad1 =-1; Float_t distbad2 =-1;
2861 Bool_t fidcut1 = 0; Bool_t fidcut2 = 0;
85c4406e 2862
a6e83e39 2863 Int_t idPartType = GetCaloPID()->GetIdentifiedParticleTypeFromClusterSplitting(calo,cells,GetCaloUtils(),
2864 GetVertex(evtIndex),nMaxima,
4914e781 2865 mass,angle,l1,l2,absId1,absId2,
2866 distbad1,distbad2,fidcut1,fidcut2) ;
34c16486 2867
a6e83e39 2868 if(GetDebug() > 1) printf("AliAnaPi0EbE::MakeShowerShapeIdentification() - PDG of identified particle %d\n",idPartType);
85c4406e 2869
2870
4914e781 2871 // Skip events where one of the new clusters (lowest energy) is close to an EMCal border or a bad channel
1253480f 2872 if( (fCheckSplitDistToBad) &&
2873 (!fidcut2 || !fidcut1 || distbad1 < fMinDist || distbad2 < fMinDist))
4914e781 2874 {
74e3eb22 2875 if(GetDebug() > 1)
2876 printf("AliAnaPi0EbE::MakeShowerShapeIdentification() - Dist to bad channel cl %f, cl1 %f, cl2 %f; fid cl1 %d, cl2 %d \n",
2877 calo->GetDistanceToBadChannel(),distbad1,distbad2, fidcut1,fidcut2);
4914e781 2878 continue ;
2879 }
2880
a6e83e39 2881 //Skip events with too few or too many NLM
85c4406e 2882 if(nMaxima < fNLMCutMin || nMaxima > fNLMCutMax)
40d3ce60 2883 {
2884 FillRejectedClusterHistograms(mom,tag);
2885 continue ;
2886 }
2887
4d97a954 2888 if(nMaxima == 1 && fNLMECutMin[0] > mom.E()) continue;
2889 if(nMaxima == 2 && fNLMECutMin[1] > mom.E()) continue;
2890 if(nMaxima > 2 && fNLMECutMin[2] > mom.E()) continue;
2891
bb2d339b 2892 if(GetDebug() > 1)
2893 printf("AliAnaPi0EbE::MakeShowerShapeIdentification() - NLM %d accepted \n",nMaxima);
a6e83e39 2894
29250849 2895 Float_t e1 = l1.Energy();
2896 Float_t e2 = l2.Energy();
2897 TLorentzVector l12 = l1+l2;
2898 Float_t ptSplit = l12.Pt();
2899 Float_t eSplit = e1+e2;
1253480f 2900
2901 Int_t mcIndex =-1;
2902 Int_t noverlaps = 0;
2903 Float_t ptprim = 0;
2904 if(IsDataMC())
2905 {
2906 mcIndex = GetMCIndex(tag);
2907
2908 Bool_t ok = kFALSE;
2909 Int_t mcLabel = calo->GetLabel();
2910
2911 TLorentzVector primary = GetMCAnalysisUtils()->GetMother(mcLabel,GetReader(),ok);
2912
2913 Int_t mesonLabel = -1;
2914
2915 if(mcIndex == kmcPi0 || mcIndex == kmcEta)
2916 {
2917 if(mcIndex == kmcPi0)
2918 {
2919 TLorentzVector grandmom = GetMCAnalysisUtils()->GetMotherWithPDG(mcLabel,111,GetReader(),ok,mesonLabel);
2920 if(grandmom.E() > 0 && ok) ptprim = grandmom.Pt();
2921 }
2922 else
2923 {
2924 TLorentzVector grandmom = GetMCAnalysisUtils()->GetMotherWithPDG(mcLabel,221,GetReader(),ok,mesonLabel);
2925 if(grandmom.E() > 0 && ok) ptprim = grandmom.Pt();
2926 }
2927 }
2928
2929 const UInt_t nlabels = calo->GetNLabels();
2930 Int_t overpdg[nlabels];
2931 noverlaps = GetMCAnalysisUtils()->GetNOverlaps(calo->GetLabels(), nlabels,tag,mesonLabel,GetReader(),overpdg);
2932 }
85c4406e 2933
f02db2c0 2934 //mass of all clusters
1253480f 2935 fhMass ->Fill(mom.E() ,mass);
29250849 2936 fhMassPt ->Fill(mom.Pt(),mass);
1253480f 2937 fhMassSplitPt->Fill(ptSplit ,mass);
85c4406e 2938
29250849 2939 if(IsDataMC())
2940 {
2941 fhMCMassPt[mcIndex] ->Fill(mom.Pt(),mass);
1253480f 2942 fhMCMassSplitPt[mcIndex]->Fill(ptSplit ,mass);
2943 if(mcIndex==kmcPi0)
2944 {
2945 fhMCPi0PtRecoPtPrim ->Fill(mom.Pt(),ptprim);
2946 fhMCPi0SplitPtRecoPtPrim->Fill(ptSplit ,ptprim);
2947 }
2948 else if(mcIndex==kmcEta)
2949 {
2950 fhMCEtaPtRecoPtPrim ->Fill(mom.Pt(),ptprim);
2951 fhMCEtaSplitPtRecoPtPrim->Fill(ptSplit ,ptprim);
2952 }
2953
2954 if(noverlaps==0)
2955 {
2956 if(mcIndex==kmcPi0)
2957 {
2958 fhMCPi0PtRecoPtPrimNoOverlap ->Fill(mom.Pt(),ptprim);
2959 fhMCPi0SplitPtRecoPtPrimNoOverlap->Fill(ptSplit ,ptprim);
2960 }
2961 else if(mcIndex==kmcEta)
2962 {
2963 fhMCEtaPtRecoPtPrimNoOverlap ->Fill(mom.Pt(),ptprim);
2964 fhMCEtaSplitPtRecoPtPrimNoOverlap->Fill(ptSplit ,ptprim);
2965 }
2966
2967 fhMassNoOverlap ->Fill(mom.E() ,mass);
2968 fhMassPtNoOverlap ->Fill(mom.Pt(),mass);
2969 fhMassSplitPtNoOverlap->Fill(ptSplit ,mass);
2970
2971 fhMCMassPtNoOverlap[mcIndex] ->Fill(mom.Pt(),mass);
2972 fhMCMassSplitPtNoOverlap[mcIndex]->Fill(ptSplit ,mass);
2973 }
29250849 2974 }
85c4406e 2975
4650f5cf 2976 // Asymmetry of all clusters
cfdf2b91 2977 Float_t asy =-10;
85c4406e 2978
4650f5cf 2979 if(e1+e2 > 0) asy = (e1-e2) / (e1+e2);
2980 fhAsymmetry->Fill(mom.E(),asy);
85c4406e 2981
2982
4650f5cf 2983 if(IsDataMC())
2984 {
4650f5cf 2985 fhMCEAsymmetry[mcIndex]->Fill(mom.E(),asy);
cfdf2b91 2986 }
f02db2c0 2987
a6e83e39 2988 // If cluster does not pass pid, not pi0/eta, skip it.
85c4406e 2989 if (GetOutputAODName().Contains("Pi0") && idPartType != AliCaloPID::kPi0)
2990 {
bb2d339b 2991 if(GetDebug() > 1) printf("AliAnaPi0EbE::MakeShowerShapeIdentification() - Cluster is not Pi0\n");
40d3ce60 2992 FillRejectedClusterHistograms(mom,tag);
bb2d339b 2993 continue ;
85c4406e 2994 }
bb2d339b 2995
85c4406e 2996 else if(GetOutputAODName().Contains("Eta") && idPartType != AliCaloPID::kEta)
2997 {
bb2d339b 2998 if(GetDebug() > 1) printf("AliAnaPi0EbE::MakeShowerShapeIdentification() - Cluster is not Eta\n");
40d3ce60 2999 FillRejectedClusterHistograms(mom,tag);
bb2d339b 3000 continue ;
85c4406e 3001 }
a6e83e39 3002
85c4406e 3003 if(GetDebug() > 1)
bb2d339b 3004 printf("AliAnaPi0EbE::MakeShowerShapeIdentification() - Pi0/Eta selection cuts passed: pT %3.2f, pdg %d\n",
85c4406e 3005 mom.Pt(), idPartType);
34c16486 3006
667432ef 3007 //Mass and asymmetry of selected pairs
29250849 3008 fhSelectedAsymmetry ->Fill(mom.E() ,asy );
3009 fhSelectedMass ->Fill(mom.E() ,mass);
3010 fhSelectedMassPt ->Fill(mom.Pt(),mass);
3011 fhSelectedMassSplitPt->Fill(ptSplit ,mass);
85c4406e 3012
1253480f 3013 if(IsDataMC())
3014 {
3015 if(mcIndex==kmcPi0)
3016 {
3017 fhMCPi0SelectedPtRecoPtPrim ->Fill(mom.Pt(),ptprim);
3018 fhMCPi0SelectedSplitPtRecoPtPrim->Fill(ptSplit ,ptprim);
3019 }
3020 else if(mcIndex==kmcEta)
3021 {
3022 fhMCEtaSelectedPtRecoPtPrim ->Fill(mom.Pt(),ptprim);
3023 fhMCEtaSelectedSplitPtRecoPtPrim->Fill(ptSplit ,ptprim);
3024 }
3025
3026 if(noverlaps==0)
3027 {
3028 fhSelectedMassNoOverlap ->Fill(mom.E() ,mass);
3029 fhSelectedMassPtNoOverlap ->Fill(mom.Pt(),mass);
3030 fhSelectedMassSplitPtNoOverlap->Fill(ptSplit ,mass);
3031
3032 if(mcIndex==kmcPi0)
3033 {
3034 fhMCPi0SelectedPtRecoPtPrimNoOverlap ->Fill(mom.Pt(),ptprim);
3035 fhMCPi0SelectedSplitPtRecoPtPrimNoOverlap->Fill(ptSplit ,ptprim);
3036 }
3037 else if(mcIndex==kmcEta)
3038 {
3039 fhMCEtaSelectedPtRecoPtPrimNoOverlap ->Fill(mom.Pt(),ptprim);
3040 fhMCEtaSelectedSplitPtRecoPtPrimNoOverlap->Fill(ptSplit ,ptprim);
3041 }
3042 }
3043 }
3044
6e66993c 3045 fhSplitE ->Fill( eSplit);
3046 fhSplitPt ->Fill(ptSplit);
29250849 3047 Float_t phi = mom.Phi();
3048 if(phi<0) phi+=TMath::TwoPi();
3049 fhSplitPtPhi ->Fill(ptSplit,phi);
3050 fhSplitPtEta ->Fill(ptSplit,mom.Eta());
6e66993c 3051 fhNLocMaxSplitPt->Fill(ptSplit ,nMaxima);
3052 fhNLocMaxPt ->Fill(mom.Pt(),nMaxima);
85c4406e 3053
a1fd1b69 3054 //Check split-clusters with good time window difference
3055 Double_t tof1 = cells->GetCellTime(absId1);
3056 GetCaloUtils()->RecalibrateCellTime(tof1, fCalorimeter, absId1,GetReader()->GetInputEvent()->GetBunchCrossNumber());
3057 tof1*=1.e9;
3058
3059 Double_t tof2 = cells->GetCellTime(absId2);
3060 GetCaloUtils()->RecalibrateCellTime(tof2, fCalorimeter, absId2,GetReader()->GetInputEvent()->GetBunchCrossNumber());
3061 tof2*=1.e9;
3062
3063 Double_t t12diff = tof1-tof2;
3064 fhEPairDiffTime->Fill(e1+e2, t12diff);
3065
cfdf2b91 3066 if(IsDataMC())
3067 {
6e66993c 3068 fhMCSplitE [mcIndex]->Fill( eSplit);
3069 fhMCSplitPt [mcIndex]->Fill(ptSplit);
29250849 3070 fhMCSplitPtPhi [mcIndex]->Fill(ptSplit,phi);
3071 fhMCSplitPtEta [mcIndex]->Fill(ptSplit,mom.Eta());
6e66993c 3072 fhMCNLocMaxSplitPt[mcIndex]->Fill(ptSplit ,nMaxima);
3073 fhMCNLocMaxPt [mcIndex]->Fill(mom.Pt(),nMaxima);
29250849 3074
3075 fhMCSelectedMassPt [mcIndex]->Fill(mom.Pt(),mass);
3076 fhMCSelectedMassSplitPt[mcIndex]->Fill(ptSplit,mass);
1253480f 3077 if(noverlaps==0)
3078 {
3079 fhMCSelectedMassPtNoOverlap [mcIndex]->Fill(mom.Pt(),mass);
3080 fhMCSelectedMassSplitPtNoOverlap[mcIndex]->Fill(ptSplit,mass);
3081 }
cfdf2b91 3082 }
cfdf2b91 3083
a6e83e39 3084 //-----------------------
3085 //Create AOD for analysis
477d6cee 3086
a6e83e39 3087 AliAODPWG4Particle aodpi0 = AliAODPWG4Particle(mom);
3088 aodpi0.SetLabel(calo->GetLabel());
3089
85c4406e 3090 //Set the indeces of the original caloclusters
a6e83e39 3091 aodpi0.SetCaloLabel(calo->GetID(),-1);
3092 aodpi0.SetDetector(fCalorimeter);
85c4406e 3093
a6e83e39 3094 if (distBad > fMinDist3) aodpi0.SetDistToBad(2) ;
85c4406e 3095 else if(distBad > fMinDist2) aodpi0.SetDistToBad(1) ;
a6e83e39 3096 else aodpi0.SetDistToBad(0) ;
3097
3098 // Check if cluster is pi0 via cluster splitting
85c4406e 3099 aodpi0.SetIdentifiedParticleType(idPartType);
3100
8736d400 3101 // Add number of local maxima to AOD, method name in AOD to be FIXED
3102 aodpi0.SetFiducialArea(nMaxima);
3103
4650f5cf 3104 aodpi0.SetTag(tag);
477d6cee 3105
34c16486 3106 //Fill some histograms about shower shape
3107 if(fFillSelectClHisto && GetReader()->GetDataType()!=AliCaloTrackReader::kMC)
3108 {
bfdcf7fb 3109 FillSelectedClusterHistograms(calo, nMaxima, tag, asy);
85c4406e 3110 }
2ad19c3d 3111
3112 // Fill histograms to undertand pile-up before other cuts applied
3113 // Remember to relax time cuts in the reader
c2a62a94 3114 Double_t tofcluster = calo->GetTOF()*1e9;
3115 Double_t tofclusterUS = TMath::Abs(tofcluster);
85c4406e 3116
126b8c62 3117 FillPileUpHistograms(aodpi0.Pt(),tofcluster,calo);
85c4406e 3118
afb3af8a 3119 Int_t id = GetReader()->GetTriggerClusterId();
3120 if(fFillEMCALBCHistograms && fCalorimeter=="EMCAL" && id >=0 )
c2a62a94 3121 {
3122 Float_t phicluster = aodpi0.Phi();
3123 if(phicluster < 0) phicluster+=TMath::TwoPi();
3124
3125 if(calo->E() > 2)
3126 {
3127 if (tofclusterUS < 25) fhEtaPhiEMCALBC0->Fill(aodpi0.Eta(), phicluster);
3128 else if (tofclusterUS < 75) fhEtaPhiEMCALBC1->Fill(aodpi0.Eta(), phicluster);
3129 else fhEtaPhiEMCALBCN->Fill(aodpi0.Eta(), phicluster);
3130 }
3131
afb3af8a 3132 Int_t bc = GetReader()->GetTriggerClusterBC();
3133 if(TMath::Abs(bc) < 6 && !GetReader()->IsBadCellTriggerEvent() && !GetReader()->IsExoticEvent() )
c2a62a94 3134 {
afb3af8a 3135 if(GetReader()->IsTriggerMatched())
3136 {
3137 if(calo->E() > 2) fhEtaPhiTriggerEMCALBC[bc+5]->Fill(aodpi0.Eta(), phicluster);
3138 fhTimeTriggerEMCALBC[bc+5]->Fill(calo->E(), tofcluster);
3139 if(GetReader()->IsPileUpFromSPD()) fhTimeTriggerEMCALBCPileUpSPD[bc+5]->Fill(calo->E(), tofcluster);
3140 }
3141 else
3142 {
3143 if(calo->E() > 2) fhEtaPhiTriggerEMCALBCUM[bc+5]->Fill(aodpi0.Eta(), phicluster);
3144 fhTimeTriggerEMCALBCUM[bc+5]->Fill(calo->E(), tofcluster);
126b8c62 3145
3146 if(bc==0)
3147 {
3148 if(GetReader()->IsTriggerMatchedOpenCuts(0)) fhTimeTriggerEMCALBC0UMReMatchOpenTime ->Fill(calo->E(), tofcluster);
3149 if(GetReader()->IsTriggerMatchedOpenCuts(1)) fhTimeTriggerEMCALBC0UMReMatchCheckNeigh ->Fill(calo->E(), tofcluster);
3150 if(GetReader()->IsTriggerMatchedOpenCuts(2)) fhTimeTriggerEMCALBC0UMReMatchBoth ->Fill(calo->E(), tofcluster);
3151 }
3152 }
c2a62a94 3153 }
afb3af8a 3154 else if(TMath::Abs(bc) >= 6)
3155 printf("AliAnaPi0EbE::MakeShowerShapeIdentification() - Trigger BC not expected = %d\n",bc);
c2a62a94 3156 }
2ad19c3d 3157
477d6cee 3158 //Add AOD with pi0 object to aod branch
3159 AddAODParticle(aodpi0);
3160
3161 }//loop
3162
85c4406e 3163 if(GetDebug() > 1) printf("AliAnaPi0EbE::MakeShowerShapeIdentification() - End fill AODs \n");
477d6cee 3164
3165}
e7fd282f 3166//______________________________________________
85c4406e 3167void AliAnaPi0EbE::MakeAnalysisFillHistograms()
691bdd02 3168{
477d6cee 3169 //Do analysis and fill histograms
691bdd02 3170
b5dbb99b 3171 if(!GetOutputAODBranch())
3172 {
a3aebfff 3173 printf("AliAnaPi0EbE::MakeAnalysisFillHistograms() - No output pi0 in AOD branch with name < %s >,STOP \n",GetOutputAODName().Data());
477d6cee 3174 abort();
3175 }
3176 //Loop on stored AOD pi0
3177 Int_t naod = GetOutputAODBranch()->GetEntriesFast();
a3aebfff 3178 if(GetDebug() > 0) printf("AliAnaPi0EbE::MakeAnalysisFillHistograms() - aod branch entries %d\n", naod);
477d6cee 3179
c8710850 3180 Float_t cen = GetEventCentrality();
3181 Float_t ep = GetEventPlaneAngle();
3182
b5dbb99b 3183 for(Int_t iaod = 0; iaod < naod ; iaod++)
3184 {
477d6cee 3185
3186 AliAODPWG4Particle* pi0 = (AliAODPWG4Particle*) (GetOutputAODBranch()->At(iaod));
21a4b1c0 3187 Int_t pdg = pi0->GetIdentifiedParticleType();
9415d854 3188
85c4406e 3189 if(IsCaloPIDOn() && pdg != AliCaloPID::kPi0) continue;
477d6cee 3190
85c4406e 3191 //Fill pi0 histograms
c4a7d28a 3192 Float_t ener = pi0->E();
3193 Float_t pt = pi0->Pt();
3194 Float_t phi = pi0->Phi();
57b97dc6 3195 if(phi < 0) phi+=TMath::TwoPi();
477d6cee 3196 Float_t eta = pi0->Eta();
3197
c8710850 3198 fhPt ->Fill(pt );
09273901 3199 fhE ->Fill(ener);
477d6cee 3200
09273901 3201 fhEEta ->Fill(ener,eta);
3202 fhEPhi ->Fill(ener,phi);
29250849 3203 fhPtEta ->Fill(pt ,eta);
3204 fhPtPhi ->Fill(pt ,phi);
c8710850 3205 fhEtaPhi ->Fill(eta ,phi);
85c4406e 3206
c8710850 3207 fhPtCentrality ->Fill(pt,cen) ;
3208 fhPtEventPlane ->Fill(pt,ep ) ;
3209
b5dbb99b 3210 if(IsDataMC())
3211 {
3455f821 3212 Int_t tag = pi0->GetTag();
3213 Int_t mcIndex = GetMCIndex(tag);
85c4406e 3214
40d3ce60 3215 fhMCE [mcIndex] ->Fill(ener);
3455f821 3216 fhMCPt [mcIndex] ->Fill(pt);
3217 fhMCPhi[mcIndex] ->Fill(pt,phi);
3218 fhMCEta[mcIndex] ->Fill(pt,eta);
3219
17f5b4b6 3220 fhMCPtCentrality[mcIndex]->Fill(pt,cen);
85c4406e 3221
883411b2 3222 if((mcIndex==kmcPhoton || mcIndex==kmcPi0 || mcIndex==kmcEta) && fAnaType==kSSCalo)
af722ce4 3223 {
36769d30 3224 Float_t efracMC = 0;
3225 Int_t label = pi0->GetLabel();
3226 Int_t momlabel = -1;
3227 Bool_t ok = kFALSE;
51a0ace5 3228
85c4406e 3229 TLorentzVector mom = GetMCAnalysisUtils()->GetMother(label,GetReader(),ok);
51a0ace5 3230 if(!ok) continue;
3231
3232 if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0))
3233 {
36769d30 3234 TLorentzVector grandmom = GetMCAnalysisUtils()->GetMotherWithPDG(label,111,GetReader(),ok,momlabel);
85c4406e 3235 if(grandmom.E() > 0 && ok)
51a0ace5 3236 {
883411b2 3237 efracMC = grandmom.E()/ener;
3238 fhMCPi0PtGenRecoFraction ->Fill(pt,efracMC);
51a0ace5 3239 }
85c4406e 3240 }
51a0ace5 3241 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0Decay))
3455f821 3242 {
3243 fhMCPi0DecayPt->Fill(pt);
36769d30 3244 TLorentzVector grandmom = GetMCAnalysisUtils()->GetMotherWithPDG(label,111,GetReader(),ok,momlabel);
85c4406e 3245 if(grandmom.E() > 0 && ok)
51a0ace5 3246 {
3247 efracMC = mom.E()/grandmom.E();
3248 fhMCPi0DecayPtFraction ->Fill(pt,efracMC);
3249 }
3455f821 3250 }
51a0ace5 3251 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEta))
3252 {
36769d30 3253 TLorentzVector grandmom = GetMCAnalysisUtils()->GetMotherWithPDG(label,221,GetReader(),ok,momlabel);
85c4406e 3254 if(grandmom.E() > 0 && ok)
51a0ace5 3255 {
883411b2 3256 efracMC = grandmom.E()/ener;
3257 fhMCEtaPtGenRecoFraction ->Fill(pt,efracMC);
51a0ace5 3258 }
85c4406e 3259 }
3455f821 3260 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEtaDecay))
3261 {
3262 fhMCEtaDecayPt->Fill(pt);
36769d30 3263 TLorentzVector grandmom = GetMCAnalysisUtils()->GetMotherWithPDG(label,221,GetReader(),ok,momlabel);
85c4406e 3264 if(grandmom.E() > 0 && ok)
51a0ace5 3265 {
3266 efracMC = mom.E()/grandmom.E();
3267 fhMCEtaDecayPtFraction ->Fill(pt,efracMC);
3268 }
3455f821 3269 }
3270 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCOtherDecay))
3271 {
3272 fhMCOtherDecayPt->Fill(pt);
3273 }
af722ce4 3274
477d6cee 3275 }
3455f821 3276
477d6cee 3277 }//Histograms with MC
3278
3279 }// aod loop
3280
3281}
3282
477d6cee 3283//__________________________________________________________________
3284void AliAnaPi0EbE::Print(const Option_t * opt) const
3285{
3286 //Print some relevant parameters set for the analysis
3287 if(! opt)
3288 return;
3289
3290 printf("**** Print %s %s ****\n", GetName(), GetTitle() ) ;
745913ae 3291 AliAnaCaloTrackCorrBaseClass::Print("");
477d6cee 3292 printf("Analysis Type = %d \n", fAnaType) ;
85c4406e 3293 if(fAnaType == kSSCalo)
3294 {
477d6cee 3295 printf("Calorimeter = %s\n", fCalorimeter.Data()) ;
3296 printf("Min Distance to Bad Channel = %2.1f\n",fMinDist);
3297 printf("Min Distance to Bad Channel 2 = %2.1f\n",fMinDist2);
85c4406e 3298 printf("Min Distance to Bad Channel 3 = %2.1f\n",fMinDist3);
477d6cee 3299 }
3300 printf(" \n") ;
3301
3302}
78a28af3 3303
78a28af3 3304