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