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