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