]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGGA/CaloTrackCorrelations/AliAnaElectron.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
add a check properly the bad/exotic cluster in the corresponding array, add few histo...
[u/mrichter/AliRoot.git] / PWGGA / CaloTrackCorrelations / AliAnaElectron.cxx
CommitLineData
d9105d92 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 *
8 * documentation strictly for non-commercial purposes hereby granted *
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 **************************************************************************/
d9105d92 15
16//_________________________________________________________________________
17//
18// Class for the photon identification.
19// Clusters from calorimeters are identified as photons
20// and kept in the AOD. Few histograms produced.
21// Copy of AliAnaPhoton just add electron id.
22//
23// -- Author: Gustavo Conesa (LPSC-IN2P3-CRNS)
24//////////////////////////////////////////////////////////////////////////////
25
26
27// --- ROOT system ---
28#include <TH2F.h>
29#include <TH3D.h>
30#include <TClonesArray.h>
31#include <TObjString.h>
32//#include <Riostream.h>
33#include "TParticle.h"
34#include "TDatabasePDG.h"
35#include "AliVTrack.h"
36
37// --- Analysis system ---
38#include "AliAnaElectron.h"
39#include "AliCaloTrackReader.h"
40#include "AliStack.h"
41#include "AliCaloPID.h"
42#include "AliMCAnalysisUtils.h"
43#include "AliFiducialCut.h"
44#include "AliVCluster.h"
45#include "AliAODMCParticle.h"
46#include "AliMixedEvent.h"
47
48
49ClassImp(AliAnaElectron)
50
78a28af3 51//________________________________
d9105d92 52AliAnaElectron::AliAnaElectron() :
764ab1f4 53 AliAnaCaloTrackCorrBaseClass(), fCalorimeter(""),
78a28af3 54 fMinDist(0.), fMinDist2(0.), fMinDist3(0.),
55 fTimeCutMin(-1), fTimeCutMax(999999),
7e9a1194 56 fNCellsCut(0), fNLMCutMin(-1), fNLMCutMax(10),
57 fFillSSHistograms(kFALSE), fFillOnlySimpleSSHisto(1),
78a28af3 58 fFillWeightHistograms(kFALSE), fNOriginHistograms(8),
59 fdEdxMin(0.), fdEdxMax (200.),
60 fEOverPMin(0), fEOverPMax (2),
764ab1f4 61 fAODParticle(0),
d9105d92 62 // Histograms
78a28af3 63 fhdEdxvsE(0), fhdEdxvsP(0),
64 fhEOverPvsE(0), fhEOverPvsP(0),
7e9a1194 65 fhdEdxvsECutM02(0), fhdEdxvsPCutM02(0),
66 fhEOverPvsECutM02(0), fhEOverPvsPCutM02(0),
67 fhdEdxvsECutEOverP(0), fhdEdxvsPCutEOverP(0),
68 fhEOverPvsECutM02CutdEdx(0), fhEOverPvsPCutM02CutdEdx(0),
78a28af3 69 // Weight studies
70 fhECellClusterRatio(0), fhECellClusterLogRatio(0),
71 fhEMaxCellClusterRatio(0), fhEMaxCellClusterLogRatio(0),
d9105d92 72 // MC histograms
73 // Electron SS MC histograms
74 fhMCElectronELambda0NoOverlap(0),
75 fhMCElectronELambda0TwoOverlap(0), fhMCElectronELambda0NOverlap(0),
76
77 //Embedding
78 fhEmbeddedSignalFractionEnergy(0),
79 fhEmbedElectronELambda0FullSignal(0), fhEmbedElectronELambda0MostlySignal(0),
80 fhEmbedElectronELambda0MostlyBkg(0), fhEmbedElectronELambda0FullBkg(0)
81{
82 //default ctor
34c16486 83 for(Int_t index = 0; index < 2; index++)
84 {
7e9a1194 85 fhNCellsE [index] = 0;
86 fhNLME [index] = 0;
87 fhTimeE [index] = 0;
d9105d92 88 fhMaxCellDiffClusterE[index] = 0;
992b14a7 89 fhE [index] = 0;
90 fhPt [index] = 0;
d9105d92 91 fhPhi [index] = 0;
92 fhEta [index] = 0;
93 fhEtaPhi [index] = 0;
94 fhEtaPhi05[index] = 0;
95
96 // Shower shape histograms
97 fhDispE [index] = 0;
98 fhLam0E [index] = 0;
99 fhLam1E [index] = 0;
100 fhDispETRD[index] = 0;
101 fhLam0ETRD[index] = 0;
102 fhLam1ETRD[index] = 0;
103 fhNCellsLam0LowE [index] = 0;
104 fhNCellsLam0HighE[index] = 0;
105 fhEtaLam0LowE [index] = 0;
106 fhPhiLam0LowE [index] = 0;
107 fhEtaLam0HighE [index] = 0;
108 fhPhiLam0HighE [index] = 0;
109
34c16486 110 fhDispEtaE [index] = 0;
111 fhDispPhiE [index] = 0;
112 fhSumEtaE [index] = 0;
113 fhSumPhiE [index] = 0;
114 fhSumEtaPhiE [index] = 0;
115 fhDispEtaPhiDiffE[index] = 0;
116 fhSphericityE [index] = 0;
117
118 for(Int_t i = 0; i < 10; i++)
119 {
d9105d92 120 fhMCPt [index][i] = 0;
121 fhMCE [index][i] = 0;
122 fhMCPhi [index][i] = 0;
123 fhMCEta [index][i] = 0;
124 fhMCDeltaE [index][i] = 0;
7e9a1194 125 fhMC2E [index][i] = 0;
126 fhMCdEdxvsE [i] = 0;
127 fhMCdEdxvsP [i] = 0;
128 fhMCEOverPvsE [i] = 0;
129 fhMCEOverPvsP [i] = 0;
d9105d92 130 }
131
34c16486 132 for(Int_t i = 0; i < 6; i++)
133 {
134 fhMCELambda0 [index][i] = 0;
135 fhMCEDispEta [index][i] = 0;
136 fhMCEDispPhi [index][i] = 0;
137 fhMCESumEtaPhi [index][i] = 0;
138 fhMCEDispEtaPhiDiff[index][i] = 0;
139 fhMCESphericity [index][i] = 0;
d9105d92 140 }
141
67616439 142 for(Int_t i = 0; i < 5; i++)
34c16486 143 {
144 fhDispEtaDispPhiEBin[index][i] = 0 ;
34c16486 145 }
d9105d92 146 }
78a28af3 147
148 //Weight studies
34c16486 149 for(Int_t i =0; i < 14; i++)
150 {
78a28af3 151 fhLambda0ForW0[i] = 0;
1a72f6c5 152 //fhLambda1ForW0[i] = 0;
78a28af3 153 }
154
d9105d92 155 //Initialize parameters
156 InitParameters();
157
158}
159
78a28af3 160//____________________________________________________________________________
7e9a1194 161Bool_t AliAnaElectron::ClusterSelected(AliVCluster* calo, TLorentzVector mom, Int_t nMaxima)
d9105d92 162{
163 //Select clusters if they pass different cuts
164 if(GetDebug() > 2)
165 printf("AliAnaElectron::ClusterSelected() Current Event %d; Before selection : E %2.2f, pT %2.2f, Ecl %2.2f, phi %2.2f, eta %2.2f\n",
166 GetReader()->GetEventNumber(),
167 mom.E(), mom.Pt(),calo->E(),mom.Phi()*TMath::RadToDeg(),mom.Eta());
168
169 //.......................................
170 //If too small or big energy, skip it
171 if(mom.E() < GetMinEnergy() || mom.E() > GetMaxEnergy() ) return kFALSE ;
172 if(GetDebug() > 2) printf("\t Cluster %d Pass E Cut \n",calo->GetID());
173
174 //.......................................
175 // TOF cut, BE CAREFUL WITH THIS CUT
176 Double_t tof = calo->GetTOF()*1e9;
177 if(tof < fTimeCutMin || tof > fTimeCutMax) return kFALSE;
178 if(GetDebug() > 2) printf("\t Cluster %d Pass Time Cut \n",calo->GetID());
179
180 //.......................................
181 if(calo->GetNCells() <= fNCellsCut && GetReader()->GetDataType() != AliCaloTrackReader::kMC) return kFALSE;
182 if(GetDebug() > 2) printf("\t Cluster %d Pass NCell Cut \n",calo->GetID());
183
184 //.......................................
185 //Check acceptance selection
186 if(IsFiducialCutOn()){
187 Bool_t in = GetFiducialCut()->IsInFiducialCut(mom,fCalorimeter) ;
188 if(! in ) return kFALSE ;
189 }
190 if(GetDebug() > 2) printf("Fiducial cut passed \n");
191
192 //.......................................
193 //Skip not matched clusters with tracks
49b5c49b 194 if(!IsTrackMatched(calo, GetReader()->GetInputEvent())) {
d9105d92 195 if(GetDebug() > 2) printf("\t Reject non track-matched clusters\n");
196 return kFALSE ;
197 }
198 else if(GetDebug() > 2) printf(" Track-matching cut passed \n");
199
7e9a1194 200 //...........................................
201 // skip clusters with too many maxima
202 if(nMaxima < fNLMCutMin || nMaxima > fNLMCutMax) return kFALSE ;
203 if(GetDebug() > 2) printf(" \t Cluster %d pass NLM %d of out of range \n",calo->GetID(), nMaxima);
204
d9105d92 205 //.......................................
206 //Check Distance to Bad channel, set bit.
207 Double_t distBad=calo->GetDistanceToBadChannel() ; //Distance to bad channel
208 if(distBad < 0.) distBad=9999. ; //workout strange convension dist = -1. ;
209 if(distBad < fMinDist) {//In bad channel (PHOS cristal size 2.2x2.2 cm), EMCAL ( cell units )
210 return kFALSE ;
211 }
212 else if(GetDebug() > 2) printf("\t Bad channel cut passed %4.2f > %2.2f \n",distBad, fMinDist);
213 //printf("Cluster %d Pass Bad Dist Cut \n",icalo);
214
215 if(GetDebug() > 0)
216 printf("AliAnaElectron::ClusterSelected() Current Event %d; After selection : E %2.2f, pT %2.2f, Ecl %2.2f, phi %2.2f, eta %2.2f\n",
217 GetReader()->GetEventNumber(),
218 mom.E(), mom.Pt(),calo->E(),mom.Phi()*TMath::RadToDeg(),mom.Eta());
219
220 //All checks passed, cluster selected
221 return kTRUE;
222
223}
224
78a28af3 225//__________________________________________________________________________________________________________
226void AliAnaElectron::FillShowerShapeHistograms(AliVCluster* cluster, const Int_t mcTag, const Int_t pidTag)
227{
d9105d92 228
229 //Fill cluster Shower Shape histograms
230
231 if(!fFillSSHistograms || GetMixedEvent()) return;
232
233 Int_t pidIndex = 0;// Electron
234 if (pidTag == AliCaloPID::kElectron) pidIndex = 0;
235 else if(pidTag == AliCaloPID::kChargedHadron) pidIndex = 1;
236 else return;
237
238 Float_t energy = cluster->E();
239 Int_t ncells = cluster->GetNCells();
240 Float_t lambda0 = cluster->GetM02();
241 Float_t lambda1 = cluster->GetM20();
242 Float_t disp = cluster->GetDispersion()*cluster->GetDispersion();
243
764ab1f4 244 Float_t l0 = 0., l1 = 0.;
245 Float_t dispp= 0., dEta = 0., dPhi = 0.;
246 Float_t sEta = 0., sPhi = 0., sEtaPhi = 0.;
247
d9105d92 248 TLorentzVector mom;
34c16486 249 if(GetReader()->GetDataType() != AliCaloTrackReader::kMC)
250 {
d9105d92 251 cluster->GetMomentum(mom,GetVertex(0)) ;}//Assume that come from vertex in straight line
252 else{
253 Double_t vertex[]={0,0,0};
254 cluster->GetMomentum(mom,vertex) ;
255 }
256
257 Float_t eta = mom.Eta();
258 Float_t phi = mom.Phi();
259 if(phi < 0) phi+=TMath::TwoPi();
260
261 fhLam0E[pidIndex] ->Fill(energy,lambda0);
262 fhLam1E[pidIndex] ->Fill(energy,lambda1);
263 fhDispE[pidIndex] ->Fill(energy,disp);
264
34c16486 265 if(fCalorimeter == "EMCAL" && GetModuleNumber(cluster) > 5)
266 {
d9105d92 267 fhLam0ETRD[pidIndex]->Fill(energy,lambda0);
268 fhLam1ETRD[pidIndex]->Fill(energy,lambda1);
269 fhDispETRD[pidIndex]->Fill(energy,disp);
270 }
271
764ab1f4 272 if(!fFillOnlySimpleSSHisto)
34c16486 273 {
764ab1f4 274 if(energy < 2)
275 {
276 fhNCellsLam0LowE[pidIndex] ->Fill(ncells,lambda0);
277 fhEtaLam0LowE[pidIndex] ->Fill(eta, lambda0);
278 fhPhiLam0LowE[pidIndex] ->Fill(phi, lambda0);
279 }
280 else
281 {
282 fhNCellsLam0HighE[pidIndex]->Fill(ncells,lambda0);
283 fhEtaLam0HighE[pidIndex] ->Fill(eta, lambda0);
284 fhPhiLam0HighE[pidIndex] ->Fill(phi, lambda0);
285 }
d9105d92 286
764ab1f4 287 if(fCalorimeter == "EMCAL")
288 {
289 GetCaloUtils()->GetEMCALRecoUtils()->RecalculateClusterShowerShapeParameters(GetEMCALGeometry(), GetReader()->GetInputEvent()->GetEMCALCells(), cluster,
290 l0, l1, dispp, dEta, dPhi, sEta, sPhi, sEtaPhi);
291 fhDispEtaE [pidIndex]-> Fill(energy,dEta);
292 fhDispPhiE [pidIndex]-> Fill(energy,dPhi);
293 fhSumEtaE [pidIndex]-> Fill(energy,sEta);
294 fhSumPhiE [pidIndex]-> Fill(energy,sPhi);
295 fhSumEtaPhiE [pidIndex]-> Fill(energy,sEtaPhi);
296 fhDispEtaPhiDiffE [pidIndex]-> Fill(energy,dPhi-dEta);
297 if(dEta+dPhi>0)fhSphericityE [pidIndex]-> Fill(energy,(dPhi-dEta)/(dEta+dPhi));
298
299 if (energy < 2 ) fhDispEtaDispPhiEBin[pidIndex][0]->Fill(dEta,dPhi);
300 else if (energy < 4 ) fhDispEtaDispPhiEBin[pidIndex][1]->Fill(dEta,dPhi);
301 else if (energy < 6 ) fhDispEtaDispPhiEBin[pidIndex][2]->Fill(dEta,dPhi);
302 else if (energy < 10) fhDispEtaDispPhiEBin[pidIndex][3]->Fill(dEta,dPhi);
303 else fhDispEtaDispPhiEBin[pidIndex][4]->Fill(dEta,dPhi);
304
305 }
34c16486 306 }
307
308 if(IsDataMC())
309 {
d9105d92 310 AliVCaloCells* cells = 0;
311 if(fCalorimeter == "EMCAL") cells = GetEMCALCells();
312 else cells = GetPHOSCells();
313
314 //Fill histograms to check shape of embedded clusters
315 Float_t fraction = 0;
316 if(GetReader()->IsEmbeddedClusterSelectionOn()){//Only working for EMCAL
317
318 Float_t clusterE = 0; // recalculate in case corrections applied.
319 Float_t cellE = 0;
320 for(Int_t icell = 0; icell < cluster->GetNCells(); icell++){
321 cellE = cells->GetCellAmplitude(cluster->GetCellAbsId(icell));
322 clusterE+=cellE;
323 fraction+=cellE*cluster->GetCellAmplitudeFraction(icell);
324 }
325
326 //Fraction of total energy due to the embedded signal
327 fraction/=clusterE;
328
329 if(GetDebug() > 1 ) printf("AliAnaElectron::FillShowerShapeHistogram() - Energy fraction of embedded signal %2.3f, Energy %2.3f\n",fraction, clusterE);
330
331 fhEmbeddedSignalFractionEnergy->Fill(clusterE,fraction);
332
333 } // embedded fraction
334
335 // Get the fraction of the cluster energy that carries the cell with highest energy
336 Int_t absID =-1 ;
337 Float_t maxCellFraction = 0.;
34c16486 338 Int_t index = 0 ;
d9105d92 339 absID = GetCaloUtils()->GetMaxEnergyCell(cells, cluster,maxCellFraction);
340
341 // Check the origin and fill histograms
342 if( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPhoton) &&
343 !GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCConversion) &&
344 !GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPi0) &&
34c16486 345 !GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCEta))
346 {
347 index = kmcssPhoton;
d9105d92 348
349 }//photon no conversion
350 else if ( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCElectron &&
34c16486 351 !GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCConversion)))
352 {
353 index = kmcssElectron;
d9105d92 354
34c16486 355 if(!GetReader()->IsEmbeddedClusterSelectionOn())
356 {
d9105d92 357 //Check particle overlaps in cluster
358
359 //Compare the primary depositing more energy with the rest, if no photon/electron as comon ancestor (conversions), count as other particle
360 Int_t ancPDG = 0, ancStatus = -1;
361 TLorentzVector momentum; TVector3 prodVertex;
362 Int_t ancLabel = 0;
363 Int_t noverlaps = 1;
364 for (UInt_t ilab = 0; ilab < cluster->GetNLabels(); ilab++ ) {
365 ancLabel = GetMCAnalysisUtils()->CheckCommonAncestor(cluster->GetLabels()[0],cluster->GetLabels()[ilab], GetReader(),ancPDG,ancStatus,momentum,prodVertex);
366 if(ancPDG!=22 && TMath::Abs(ancPDG)!=11) noverlaps++;
367 }
368
369 if(noverlaps == 1){
370 fhMCElectronELambda0NoOverlap ->Fill(energy, lambda0);
371 }
372 else if(noverlaps == 2){
373 fhMCElectronELambda0TwoOverlap ->Fill(energy, lambda0);
374 }
375 else if(noverlaps > 2){
376 fhMCElectronELambda0NOverlap ->Fill(energy, lambda0);
377 }
378 else {
379 printf("AliAnaElectron::FillShowerShapeHistogram() - n overlaps = %d!!", noverlaps);
380 }
381 }//No embedding
382
383 //Fill histograms to check shape of embedded clusters
34c16486 384 if(GetReader()->IsEmbeddedClusterSelectionOn())
385 {
d9105d92 386 if (fraction > 0.9)
387 {
388 fhEmbedElectronELambda0FullSignal ->Fill(energy, lambda0);
389 }
390 else if(fraction > 0.5)
391 {
392 fhEmbedElectronELambda0MostlySignal ->Fill(energy, lambda0);
393 }
394 else if(fraction > 0.1)
395 {
396 fhEmbedElectronELambda0MostlyBkg ->Fill(energy, lambda0);
397 }
398 else
399 {
400 fhEmbedElectronELambda0FullBkg ->Fill(energy, lambda0);
401 }
402 } // embedded
403 }//electron
404 else if ( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCElectron) &&
34c16486 405 GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCConversion) )
406 {
407 index = kmcssConversion;
d9105d92 408 }//conversion photon
34c16486 409 else if ( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPi0) )
410 {
411 index = kmcssPi0;
d9105d92 412 }//pi0
34c16486 413 else if ( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCEta) )
414 {
415 index = kmcssEta;
d9105d92 416 }//eta
34c16486 417 else
418 {
419 index = kmcssOther;
d9105d92 420 }//other particles
421
34c16486 422 fhMCELambda0[pidIndex][index] ->Fill(energy, lambda0);
423
764ab1f4 424 if(fCalorimeter == "EMCAL" && !fFillOnlySimpleSSHisto)
34c16486 425 {
426 fhMCEDispEta [pidIndex][index]-> Fill(energy,dEta);
427 fhMCEDispPhi [pidIndex][index]-> Fill(energy,dPhi);
428 fhMCESumEtaPhi [pidIndex][index]-> Fill(energy,sEtaPhi);
429 fhMCEDispEtaPhiDiff [pidIndex][index]-> Fill(energy,dPhi-dEta);
430 if(dEta+dPhi>0)fhMCESphericity [pidIndex][index]-> Fill(energy,(dPhi-dEta)/(dEta+dPhi));
431
34c16486 432 }
433
434
d9105d92 435 }//MC data
436
437}
438
78a28af3 439//_____________________________________________
d9105d92 440TObjString * AliAnaElectron::GetAnalysisCuts()
441{
442 //Save parameters used for analysis
443 TString parList ; //this will be list of parameters used for this analysis.
444 const Int_t buffersize = 255;
445 char onePar[buffersize] ;
446
447 snprintf(onePar,buffersize,"--- AliAnaElectron ---\n") ;
448 parList+=onePar ;
449 snprintf(onePar,buffersize,"Calorimeter: %s\n",fCalorimeter.Data()) ;
450 parList+=onePar ;
451 snprintf(onePar,buffersize," %2.2f < dEdx < %2.2f \n",fdEdxMin,fdEdxMax) ;
452 parList+=onePar ;
453 snprintf(onePar,buffersize," %2.2f < E/P < %2.2f \n",fEOverPMin, fEOverPMax) ;
454 parList+=onePar ;
455 snprintf(onePar,buffersize,"fMinDist =%2.2f (Minimal distance to bad channel to accept cluster) \n",fMinDist) ;
456 parList+=onePar ;
457 snprintf(onePar,buffersize,"fMinDist2=%2.2f (Cuts on Minimal distance to study acceptance evaluation) \n",fMinDist2) ;
458 parList+=onePar ;
459 snprintf(onePar,buffersize,"fMinDist3=%2.2f (One more cut on distance used for acceptance-efficiency study) \n",fMinDist3) ;
460 parList+=onePar ;
461
462 //Get parameters set in base class.
463 parList += GetBaseParametersList() ;
464
465 //Get parameters set in PID class.
466 parList += GetCaloPID()->GetPIDParametersList() ;
467
468 //Get parameters set in FiducialCut class (not available yet)
7e9a1194 469 //parlist += GetFidCut()->GetFidCutParametersList()
d9105d92 470
471 return new TObjString(parList) ;
472}
473
78a28af3 474//_______________________________________________
d9105d92 475TList * AliAnaElectron::GetCreateOutputObjects()
476{
477 // Create histograms to be saved in output file and
478 // store them in outputContainer
479 TList * outputContainer = new TList() ;
480 outputContainer->SetName("ElectronHistos") ;
481
745913ae 482 Int_t nptbins = GetHistogramRanges()->GetHistoPtBins(); Float_t ptmax = GetHistogramRanges()->GetHistoPtMax(); Float_t ptmin = GetHistogramRanges()->GetHistoPtMin();
483 Int_t nphibins = GetHistogramRanges()->GetHistoPhiBins(); Float_t phimax = GetHistogramRanges()->GetHistoPhiMax(); Float_t phimin = GetHistogramRanges()->GetHistoPhiMin();
484 Int_t netabins = GetHistogramRanges()->GetHistoEtaBins(); Float_t etamax = GetHistogramRanges()->GetHistoEtaMax(); Float_t etamin = GetHistogramRanges()->GetHistoEtaMin();
485 Int_t ssbins = GetHistogramRanges()->GetHistoShowerShapeBins(); Float_t ssmax = GetHistogramRanges()->GetHistoShowerShapeMax(); Float_t ssmin = GetHistogramRanges()->GetHistoShowerShapeMin();
486 Int_t nbins = GetHistogramRanges()->GetHistoNClusterCellBins(); Int_t nmax = GetHistogramRanges()->GetHistoNClusterCellMax(); Int_t nmin = GetHistogramRanges()->GetHistoNClusterCellMin();
487 Int_t ndedxbins = GetHistogramRanges()->GetHistodEdxBins(); Float_t dedxmax = GetHistogramRanges()->GetHistodEdxMax(); Float_t dedxmin = GetHistogramRanges()->GetHistodEdxMin();
488 Int_t nPoverEbins = GetHistogramRanges()->GetHistoPOverEBins(); Float_t pOverEmax = GetHistogramRanges()->GetHistoPOverEMax(); Float_t pOverEmin = GetHistogramRanges()->GetHistoPOverEMin();
489 Int_t tbins = GetHistogramRanges()->GetHistoTimeBins() ; Float_t tmax = GetHistogramRanges()->GetHistoTimeMax(); Float_t tmin = GetHistogramRanges()->GetHistoTimeMin();
d9105d92 490
7e9a1194 491
492 // MC labels, titles, for originator particles
493 TString ptypess[] = { "#gamma","hadron?","#pi^{0}","#eta","#gamma->e^{#pm}","e^{#pm}"} ;
494 TString pnamess[] = { "Photon","Hadron" ,"Pi0" ,"Eta" ,"Conversion" ,"Electron"} ;
495 TString ptype[] = { "#gamma", "#gamma_{#pi decay}","#gamma_{other decay}", "#pi^{0}","#eta",
496 "e^{#pm}","#gamma->e^{#pm}","hadron?","Anti-N","Anti-P" } ;
497
498 TString pname[] = { "Photon","PhotonPi0Decay","PhotonOtherDecay","Pi0","Eta","Electron",
499 "Conversion", "Hadron", "AntiNeutron","AntiProton" } ;
500
501
502 fhdEdxvsE = new TH2F ("hdEdxvsE","matched track <dE/dx> vs cluster E ", nptbins,ptmin,ptmax,ndedxbins, dedxmin, dedxmax);
d9105d92 503 fhdEdxvsE->SetXTitle("E (GeV)");
504 fhdEdxvsE->SetYTitle("<dE/dx>");
505 outputContainer->Add(fhdEdxvsE);
506
507 fhdEdxvsP = new TH2F ("hdEdxvsP","matched track <dE/dx> vs track P ", nptbins,ptmin,ptmax,ndedxbins, dedxmin, dedxmax);
508 fhdEdxvsP->SetXTitle("P (GeV/c)");
509 fhdEdxvsP->SetYTitle("<dE/dx>");
510 outputContainer->Add(fhdEdxvsP);
511
512 fhEOverPvsE = new TH2F ("hEOverPvsE","matched track E/p vs cluster E ", nptbins,ptmin,ptmax,nPoverEbins,pOverEmin,pOverEmax);
513 fhEOverPvsE->SetXTitle("E (GeV)");
514 fhEOverPvsE->SetYTitle("E/p");
515 outputContainer->Add(fhEOverPvsE);
516
517 fhEOverPvsP = new TH2F ("hEOverPvsP","matched track E/p vs track P ", nptbins,ptmin,ptmax,nPoverEbins,pOverEmin,pOverEmax);
518 fhEOverPvsP->SetXTitle("P (GeV/c)");
519 fhEOverPvsP->SetYTitle("E/p");
520 outputContainer->Add(fhEOverPvsP);
521
522
7e9a1194 523 fhdEdxvsECutM02 = new TH2F ("hdEdxvsECutM02","matched track <dE/dx> vs cluster E, mild #lambda_{0}^{2} cut", nptbins,ptmin,ptmax,ndedxbins, dedxmin, dedxmax);
524 fhdEdxvsECutM02->SetXTitle("E (GeV)");
525 fhdEdxvsECutM02->SetYTitle("<dE/dx>");
526 outputContainer->Add(fhdEdxvsECutM02);
527
528 fhdEdxvsPCutM02 = new TH2F ("hdEdxvsPCutM02","matched track <dE/dx> vs track P, mild #lambda_{0}^{2} cut", nptbins,ptmin,ptmax,ndedxbins, dedxmin, dedxmax);
529 fhdEdxvsPCutM02->SetXTitle("P (GeV/c)");
530 fhdEdxvsPCutM02->SetYTitle("<dE/dx>");
531 outputContainer->Add(fhdEdxvsPCutM02);
532
533 fhEOverPvsECutM02 = new TH2F ("hEOverPvsECutM02","matched track E/p vs cluster E, mild #lambda_{0}^{2} cut", nptbins,ptmin,ptmax,nPoverEbins,pOverEmin,pOverEmax);
534 fhEOverPvsECutM02->SetXTitle("E (GeV)");
535 fhEOverPvsECutM02->SetYTitle("E/p");
536 outputContainer->Add(fhEOverPvsECutM02);
537
538 fhEOverPvsPCutM02 = new TH2F ("hEOverPvsPCutM02","matched track E/p vs track P, mild #lambda_{0}^{2} cut", nptbins,ptmin,ptmax,nPoverEbins,pOverEmin,pOverEmax);
539 fhEOverPvsPCutM02->SetXTitle("P (GeV/c)");
540 fhEOverPvsPCutM02->SetYTitle("E/p");
541 outputContainer->Add(fhEOverPvsPCutM02);
542
543
544 fhdEdxvsECutEOverP = new TH2F ("hdEdxvsECutEOverP","matched track <dE/dx> vs cluster E, cut on E/p", nptbins,ptmin,ptmax,ndedxbins, dedxmin, dedxmax);
545 fhdEdxvsECutEOverP->SetXTitle("E (GeV)");
546 fhdEdxvsECutEOverP->SetYTitle("<dE/dx>");
547 outputContainer->Add(fhdEdxvsECutEOverP);
548
549 fhdEdxvsPCutEOverP = new TH2F ("hdEdxvsPCutEOverP","matched track <dE/dx> vs track P, cut on E/p", nptbins,ptmin,ptmax,ndedxbins, dedxmin, dedxmax);
550 fhdEdxvsPCutEOverP->SetXTitle("P (GeV/c)");
551 fhdEdxvsPCutEOverP->SetYTitle("<dE/dx>");
552 outputContainer->Add(fhdEdxvsPCutEOverP);
553
554 fhEOverPvsECutM02CutdEdx = new TH2F ("hEOverPvsECutM02CutdEdx","matched track E/p vs cluster E, dEdx cut, mild #lambda_{0}^{2} cut", nptbins,ptmin,ptmax,nPoverEbins,pOverEmin,pOverEmax);
555 fhEOverPvsECutM02CutdEdx->SetXTitle("E (GeV)");
556 fhEOverPvsECutM02CutdEdx->SetYTitle("E/p");
557 outputContainer->Add(fhEOverPvsECutM02CutdEdx);
558
559 fhEOverPvsPCutM02CutdEdx = new TH2F ("hEOverPvsPCutM02CutdEdx","matched track E/p vs track P, dEdx cut, mild #lambda_{0}^{2} cut ", nptbins,ptmin,ptmax,nPoverEbins,pOverEmin,pOverEmax);
560 fhEOverPvsPCutM02CutdEdx->SetXTitle("P (GeV/c)");
561 fhEOverPvsPCutM02CutdEdx->SetYTitle("E/p");
562 outputContainer->Add(fhEOverPvsPCutM02CutdEdx);
563
564 if(IsDataMC())
565 {
566 for(Int_t i = 0; i < fNOriginHistograms; i++)
567 {
568 fhMCdEdxvsE[i] = new TH2F(Form("hdEdxvsE_MC%s",pname[i].Data()),
569 Form("matched track <dE/dx> vs cluster E from %s : E ",ptype[i].Data()),
570 nptbins,ptmin,ptmax,ndedxbins, dedxmin, dedxmax);
571 fhMCdEdxvsE[i]->SetXTitle("E (GeV)");
572 fhMCdEdxvsE[i]->SetYTitle("<dE/dx>");
573 outputContainer->Add(fhMCdEdxvsE[i]) ;
574
575 fhMCdEdxvsP[i] = new TH2F(Form("hdEdxvsP_MC%s",pname[i].Data()),
576 Form("matched track <dE/dx> vs track P from %s : E ",ptype[i].Data()),
577 nptbins,ptmin,ptmax,ndedxbins, dedxmin, dedxmax);
578 fhMCdEdxvsP[i]->SetXTitle("E (GeV)");
579 fhMCdEdxvsP[i]->SetYTitle("<dE/dx>");
580 outputContainer->Add(fhMCdEdxvsP[i]) ;
581
582
583 fhMCEOverPvsE[i] = new TH2F(Form("hEOverPvsE_MC%s",pname[i].Data()),
584 Form("matched track E/p vs cluster E from %s : E ",ptype[i].Data()),
585 nptbins,ptmin,ptmax,nPoverEbins,pOverEmin,pOverEmax);
586 fhMCEOverPvsE[i]->SetXTitle("E (GeV)");
587 fhMCEOverPvsE[i]->SetYTitle("<dE/dx>");
588 outputContainer->Add(fhMCEOverPvsE[i]) ;
589
590 fhMCEOverPvsP[i] = new TH2F(Form("hEOverPvsP_MC%s",pname[i].Data()),
591 Form("matched track E/pvs track P from %s : E ",ptype[i].Data()),
592 nptbins,ptmin,ptmax,nPoverEbins,pOverEmin,pOverEmax);
593 fhMCEOverPvsP[i]->SetXTitle("E (GeV)");
594 fhMCEOverPvsP[i]->SetYTitle("<dE/dx>");
595 outputContainer->Add(fhMCEOverPvsP[i]) ;
596
597 }
598 }
599
d9105d92 600 TString pidParticle[] = {"Electron","ChargedHadron"} ;
601
34c16486 602 if(fFillWeightHistograms)
603 {
1a72f6c5 604
605 fhECellClusterRatio = new TH2F ("hECellClusterRatio"," cell energy / cluster energy vs cluster energy, for selected electrons",
606 nptbins,ptmin,ptmax, 100,0,1.);
607 fhECellClusterRatio->SetXTitle("E_{cluster} (GeV) ");
608 fhECellClusterRatio->SetYTitle("E_{cell i}/E_{cluster}");
609 outputContainer->Add(fhECellClusterRatio);
610
611 fhECellClusterLogRatio = new TH2F ("hECellClusterLogRatio"," Log(cell energy / cluster energy) vs cluster energy, for selected electrons",
612 nptbins,ptmin,ptmax, 100,-10,0);
613 fhECellClusterLogRatio->SetXTitle("E_{cluster} (GeV) ");
614 fhECellClusterLogRatio->SetYTitle("Log (E_{max cell}/E_{cluster})");
615 outputContainer->Add(fhECellClusterLogRatio);
616
617 fhEMaxCellClusterRatio = new TH2F ("hEMaxCellClusterRatio"," max cell energy / cluster energy vs cluster energy, for selected electrons",
618 nptbins,ptmin,ptmax, 100,0,1.);
619 fhEMaxCellClusterRatio->SetXTitle("E_{cluster} (GeV) ");
620 fhEMaxCellClusterRatio->SetYTitle("E_{max cell}/E_{cluster}");
621 outputContainer->Add(fhEMaxCellClusterRatio);
622
623 fhEMaxCellClusterLogRatio = new TH2F ("hEMaxCellClusterLogRatio"," Log(max cell energy / cluster energy) vs cluster energy, for selected electrons",
624 nptbins,ptmin,ptmax, 100,-10,0);
625 fhEMaxCellClusterLogRatio->SetXTitle("E_{cluster} (GeV) ");
626 fhEMaxCellClusterLogRatio->SetYTitle("Log (E_{max cell}/E_{cluster})");
627 outputContainer->Add(fhEMaxCellClusterLogRatio);
628
34c16486 629 for(Int_t iw = 0; iw < 14; iw++)
630 {
1a72f6c5 631 fhLambda0ForW0[iw] = new TH2F (Form("hLambda0ForW0%d",iw),Form("shower shape, #lambda^{2}_{0} vs E, w0 = %1.1f, for selected electrons",1+0.5*iw),
632 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
633 fhLambda0ForW0[iw]->SetXTitle("E_{cluster}");
634 fhLambda0ForW0[iw]->SetYTitle("#lambda^{2}_{0}");
635 outputContainer->Add(fhLambda0ForW0[iw]);
636
637 // fhLambda1ForW0[iw] = new TH2F (Form("hLambda1ForW0%d",iw),Form("shower shape, #lambda^{2}_{1} vs E, w0 = %1.1f, for selected electrons",1+0.5*iw),
638 // nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
639 // fhLambda1ForW0[iw]->SetXTitle("E_{cluster}");
640 // fhLambda1ForW0[iw]->SetYTitle("#lambda^{2}_{1}");
641 // outputContainer->Add(fhLambda1ForW0[iw]);
642
643 }
644 }
645
34c16486 646 for(Int_t pidIndex = 0; pidIndex < 2; pidIndex++)
647 {
d9105d92 648 //Shower shape
34c16486 649 if(fFillSSHistograms)
650 {
d9105d92 651 fhLam0E[pidIndex] = new TH2F (Form("h%sLam0E",pidParticle[pidIndex].Data()),
652 Form("%s: #lambda_{0}^{2} vs E",pidParticle[pidIndex].Data()),
653 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
654 fhLam0E[pidIndex]->SetYTitle("#lambda_{0}^{2}");
655 fhLam0E[pidIndex]->SetXTitle("E (GeV)");
656 outputContainer->Add(fhLam0E[pidIndex]);
657
658 fhLam1E[pidIndex] = new TH2F (Form("h%sLam1E",pidParticle[pidIndex].Data()),
659 Form("%s: #lambda_{1}^{2} vs E",pidParticle[pidIndex].Data()),
660 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
661 fhLam1E[pidIndex]->SetYTitle("#lambda_{1}^{2}");
662 fhLam1E[pidIndex]->SetXTitle("E (GeV)");
663 outputContainer->Add(fhLam1E[pidIndex]);
664
665 fhDispE[pidIndex] = new TH2F (Form("h%sDispE",pidParticle[pidIndex].Data()),
666 Form("%s: dispersion^{2} vs E",pidParticle[pidIndex].Data()),
667 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
668 fhDispE[pidIndex]->SetYTitle("D^{2}");
669 fhDispE[pidIndex]->SetXTitle("E (GeV) ");
670 outputContainer->Add(fhDispE[pidIndex]);
671
34c16486 672 if(fCalorimeter == "EMCAL")
673 {
d9105d92 674 fhLam0ETRD[pidIndex] = new TH2F (Form("h%sLam0ETRD",pidParticle[pidIndex].Data()),
675 Form("%s: #lambda_{0}^{2} vs E, EMCAL SM covered by TRD",pidParticle[pidIndex].Data()),
676 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
677 fhLam0ETRD[pidIndex]->SetYTitle("#lambda_{0}^{2}");
678 fhLam0ETRD[pidIndex]->SetXTitle("E (GeV)");
679 outputContainer->Add(fhLam0ETRD[pidIndex]);
680
681 fhLam1ETRD[pidIndex] = new TH2F (Form("h%sLam1ETRD",pidParticle[pidIndex].Data()),
682 Form("%s: #lambda_{1}^{2} vs E, EMCAL SM covered by TRD",pidParticle[pidIndex].Data()),
683 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
684 fhLam1ETRD[pidIndex]->SetYTitle("#lambda_{1}^{2}");
685 fhLam1ETRD[pidIndex]->SetXTitle("E (GeV)");
686 outputContainer->Add(fhLam1ETRD[pidIndex]);
687
688 fhDispETRD[pidIndex] = new TH2F (Form("h%sDispETRD",pidParticle[pidIndex].Data()),
689 Form("%s: dispersion^{2} vs E, EMCAL SM covered by TRD",pidParticle[pidIndex].Data()),
690 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
691 fhDispETRD[pidIndex]->SetYTitle("Dispersion^{2}");
692 fhDispETRD[pidIndex]->SetXTitle("E (GeV) ");
693 outputContainer->Add(fhDispETRD[pidIndex]);
694 }
695
764ab1f4 696 if(!fFillOnlySimpleSSHisto)
34c16486 697 {
34c16486 698
764ab1f4 699 fhNCellsLam0LowE[pidIndex] = new TH2F (Form("h%sNCellsLam0LowE",pidParticle[pidIndex].Data()),
700 Form("%s: N_{cells} in cluster vs #lambda_{0}^{2}, E < 2 GeV",pidParticle[pidIndex].Data()),
701 nbins,nmin, nmax, ssbins,ssmin,ssmax);
702 fhNCellsLam0LowE[pidIndex]->SetXTitle("N_{cells}");
703 fhNCellsLam0LowE[pidIndex]->SetYTitle("#lambda_{0}^{2}");
704 outputContainer->Add(fhNCellsLam0LowE[pidIndex]);
34c16486 705
764ab1f4 706 fhNCellsLam0HighE[pidIndex] = new TH2F (Form("h%sNCellsLam0HighE",pidParticle[pidIndex].Data()),
707 Form("%s: N_{cells} in cluster vs #lambda_{0}^{2}, E > 2 GeV",pidParticle[pidIndex].Data()),
708 nbins,nmin, nmax, ssbins,ssmin,ssmax);
709 fhNCellsLam0HighE[pidIndex]->SetXTitle("N_{cells}");
710 fhNCellsLam0HighE[pidIndex]->SetYTitle("#lambda_{0}^{2}");
711 outputContainer->Add(fhNCellsLam0HighE[pidIndex]);
78a28af3 712
34c16486 713
764ab1f4 714 fhEtaLam0LowE[pidIndex] = new TH2F (Form("h%sEtaLam0LowE",pidParticle[pidIndex].Data()),
715 Form("%s: #eta vs #lambda_{0}^{2}, E < 2 GeV",pidParticle[pidIndex].Data()),
716 netabins,etamin,etamax, ssbins,ssmin,ssmax);
717 fhEtaLam0LowE[pidIndex]->SetYTitle("#lambda_{0}^{2}");
718 fhEtaLam0LowE[pidIndex]->SetXTitle("#eta");
719 outputContainer->Add(fhEtaLam0LowE[pidIndex]);
34c16486 720
764ab1f4 721 fhPhiLam0LowE[pidIndex] = new TH2F (Form("h%sPhiLam0LowE",pidParticle[pidIndex].Data()),
722 Form("%s: #phi vs #lambda_{0}^{2}, E < 2 GeV",pidParticle[pidIndex].Data()),
723 nphibins,phimin,phimax, ssbins,ssmin,ssmax);
724 fhPhiLam0LowE[pidIndex]->SetYTitle("#lambda_{0}^{2}");
725 fhPhiLam0LowE[pidIndex]->SetXTitle("#phi");
726 outputContainer->Add(fhPhiLam0LowE[pidIndex]);
34c16486 727
764ab1f4 728 fhEtaLam0HighE[pidIndex] = new TH2F (Form("h%sEtaLam0HighE",pidParticle[pidIndex].Data()),
729 Form("%s: #eta vs #lambda_{0}^{2}, E > 2 GeV",pidParticle[pidIndex].Data()),
730 netabins,etamin,etamax, ssbins,ssmin,ssmax);
731 fhEtaLam0HighE[pidIndex]->SetYTitle("#lambda_{0}^{2}");
732 fhEtaLam0HighE[pidIndex]->SetXTitle("#eta");
733 outputContainer->Add(fhEtaLam0HighE[pidIndex]);
34c16486 734
764ab1f4 735 fhPhiLam0HighE[pidIndex] = new TH2F (Form("h%sPhiLam0HighE",pidParticle[pidIndex].Data()),
736 Form("%s: #phi vs #lambda_{0}^{2}, E > 2 GeV",pidParticle[pidIndex].Data()),
737 nphibins,phimin,phimax, ssbins,ssmin,ssmax);
738 fhPhiLam0HighE[pidIndex]->SetYTitle("#lambda_{0}^{2}");
739 fhPhiLam0HighE[pidIndex]->SetXTitle("#phi");
740 outputContainer->Add(fhPhiLam0HighE[pidIndex]);
741
742 if(fCalorimeter == "EMCAL")
34c16486 743 {
764ab1f4 744 fhDispEtaE[pidIndex] = new TH2F (Form("h%sDispEtaE",pidParticle[pidIndex].Data()),
745 Form("%s: #sigma^{2}_{#eta #eta} = #Sigma w_{i}(#eta_{i} - <#eta>)^{2}/ #Sigma w_{i} vs E",pidParticle[pidIndex].Data()),
746 nptbins,ptmin,ptmax, ssbins,ssmin,ssmax);
747 fhDispEtaE[pidIndex]->SetXTitle("E (GeV)");
748 fhDispEtaE[pidIndex]->SetYTitle("#sigma^{2}_{#eta #eta}");
749 outputContainer->Add(fhDispEtaE[pidIndex]);
750
751 fhDispPhiE[pidIndex] = new TH2F (Form("h%sDispPhiE",pidParticle[pidIndex].Data()),
752 Form("%s: #sigma^{2}_{#phi #phi} = #Sigma w_{i}(#phi_{i} - <#phi>)^{2} / #Sigma w_{i} vs E",pidParticle[pidIndex].Data()),
753 nptbins,ptmin,ptmax, ssbins,ssmin,ssmax);
754 fhDispPhiE[pidIndex]->SetXTitle("E (GeV)");
755 fhDispPhiE[pidIndex]->SetYTitle("#sigma^{2}_{#phi #phi}");
756 outputContainer->Add(fhDispPhiE[pidIndex]);
757
758 fhSumEtaE[pidIndex] = new TH2F (Form("h%sSumEtaE",pidParticle[pidIndex].Data()),
759 Form("%s: #sigma^{2}_{#eta #eta} = #Sigma w_{i}(#eta_{i})^{2} / #Sigma w_{i} - <#eta>^{2} vs E",pidParticle[pidIndex].Data()),
760 nptbins,ptmin,ptmax, ssbins,ssmin,ssmax);
761 fhSumEtaE[pidIndex]->SetXTitle("E (GeV)");
762 fhSumEtaE[pidIndex]->SetYTitle("#delta^{2}_{#eta #eta}");
763 outputContainer->Add(fhSumEtaE[pidIndex]);
764
765 fhSumPhiE[pidIndex] = new TH2F (Form("h%sSumPhiE",pidParticle[pidIndex].Data()),
766 Form("%s: #sigma^{2}_{#phi #phi} = #Sigma w_{i}(#phi_{i})^{2}/ #Sigma w_{i} - <#phi>^{2} vs E",pidParticle[pidIndex].Data()),
767 nptbins,ptmin,ptmax, ssbins,ssmin,ssmax);
768 fhSumPhiE[pidIndex]->SetXTitle("E (GeV)");
769 fhSumPhiE[pidIndex]->SetYTitle("#delta^{2}_{#phi #phi}");
770 outputContainer->Add(fhSumPhiE[pidIndex]);
771
772 fhSumEtaPhiE[pidIndex] = new TH2F (Form("h%sSumEtaPhiE",pidParticle[pidIndex].Data()),
773 Form("%s: #delta^{2}_{#eta #phi} = #Sigma w_{i}(#phi_{i} #eta_{i} ) / #Sigma w_{i} - <#phi><#eta> vs E",pidParticle[pidIndex].Data()),
774 nptbins,ptmin,ptmax, 2*ssbins,-ssmax,ssmax);
775 fhSumEtaPhiE[pidIndex]->SetXTitle("E (GeV)");
776 fhSumEtaPhiE[pidIndex]->SetYTitle("#delta^{2}_{#eta #phi}");
777 outputContainer->Add(fhSumEtaPhiE[pidIndex]);
778
779 fhDispEtaPhiDiffE[pidIndex] = new TH2F (Form("h%sDispEtaPhiDiffE",pidParticle[pidIndex].Data()),
780 Form("%s: #sigma^{2}_{#phi #phi} - #sigma^{2}_{#eta #eta} vs E",pidParticle[pidIndex].Data()),
781 nptbins,ptmin,ptmax,200, -10,10);
782 fhDispEtaPhiDiffE[pidIndex]->SetXTitle("E (GeV)");
783 fhDispEtaPhiDiffE[pidIndex]->SetYTitle("#sigma^{2}_{#phi #phi}-#sigma^{2}_{#eta #eta}");
784 outputContainer->Add(fhDispEtaPhiDiffE[pidIndex]);
785
786 fhSphericityE[pidIndex] = new TH2F (Form("h%sSphericityE",pidParticle[pidIndex].Data()),
787 Form("%s: (#sigma^{2}_{#phi #phi} - #sigma^{2}_{#eta #eta}) / (#sigma^{2}_{#eta #eta} + #sigma^{2}_{#phi #phi}) vs E",pidParticle[pidIndex].Data()),
788 nptbins,ptmin,ptmax, 200, -1,1);
789 fhSphericityE[pidIndex]->SetXTitle("E (GeV)");
790 fhSphericityE[pidIndex]->SetYTitle("s = (#sigma^{2}_{#phi #phi} - #sigma^{2}_{#eta #eta}) / (#sigma^{2}_{#eta #eta} + #sigma^{2}_{#phi #phi})");
791 outputContainer->Add(fhSphericityE[pidIndex]);
792
793 Int_t bin[] = {0,2,4,6,10,1000};
794 for(Int_t i = 0; i < 5; i++)
795 {
796 fhDispEtaDispPhiEBin[pidIndex][i] = new TH2F (Form("h%sDispEtaDispPhi_EBin%d",pidParticle[pidIndex].Data(),i),
797 Form("%s: #sigma^{2}_{#phi #phi} vs #sigma^{2}_{#eta #eta} for %d < E < %d GeV",pidParticle[pidIndex].Data(),bin[i],bin[i+1]),
798 ssbins,ssmin,ssmax , ssbins,ssmin,ssmax);
799 fhDispEtaDispPhiEBin[pidIndex][i]->SetXTitle("#sigma^{2}_{#eta #eta}");
800 fhDispEtaDispPhiEBin[pidIndex][i]->SetYTitle("#sigma^{2}_{#phi #phi}");
801 outputContainer->Add(fhDispEtaDispPhiEBin[pidIndex][i]);
802 }
34c16486 803 }
804 }
34c16486 805 } // Shower shape
806
807 if(IsDataMC())
808 {
809 if(fFillSSHistograms)
810 {
34c16486 811 for(Int_t i = 0; i < 6; i++)
812 {
d9105d92 813 fhMCELambda0[pidIndex][i] = new TH2F(Form("h%sELambda0_MC%s",pidParticle[pidIndex].Data(),pnamess[i].Data()),
814 Form("%s like cluster from %s : E vs #lambda_{0}^{2}",pidParticle[pidIndex].Data(),ptypess[i].Data()),
815 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
816 fhMCELambda0[pidIndex][i]->SetYTitle("#lambda_{0}^{2}");
817 fhMCELambda0[pidIndex][i]->SetXTitle("E (GeV)");
818 outputContainer->Add(fhMCELambda0[pidIndex][i]) ;
819
764ab1f4 820 if(fCalorimeter=="EMCAL" && !fFillOnlySimpleSSHisto)
34c16486 821 {
822 fhMCEDispEta[pidIndex][i] = new TH2F (Form("h%sEDispEtaE_MC%s",pidParticle[pidIndex].Data(),pnamess[i].Data()),
823 Form("cluster from %s : %s like, #sigma^{2}_{#eta #eta} = #Sigma w_{i}(#eta_{i} - <#eta>)^{2}/ #Sigma w_{i} vs E",ptypess[i].Data(),pidParticle[pidIndex].Data()),
824 nptbins,ptmin,ptmax, ssbins,ssmin,ssmax);
825 fhMCEDispEta[pidIndex][i]->SetXTitle("E (GeV)");
826 fhMCEDispEta[pidIndex][i]->SetYTitle("#sigma^{2}_{#eta #eta}");
827 outputContainer->Add(fhMCEDispEta[pidIndex][i]);
828
829 fhMCEDispPhi[pidIndex][i] = new TH2F (Form("h%sEDispPhiE_MC%s",pidParticle[pidIndex].Data(),pnamess[i].Data()),
830 Form("cluster from %s : %s like, #sigma^{2}_{#phi #phi} = #Sigma w_{i}(#phi_{i} - <#phi>)^{2} / #Sigma w_{i} vs E",ptypess[i].Data(),pidParticle[pidIndex].Data()),
831 nptbins,ptmin,ptmax, ssbins,ssmin,ssmax);
832 fhMCEDispPhi[pidIndex][i]->SetXTitle("E (GeV)");
833 fhMCEDispPhi[pidIndex][i]->SetYTitle("#sigma^{2}_{#phi #phi}");
834 outputContainer->Add(fhMCEDispPhi[pidIndex][i]);
835
836 fhMCESumEtaPhi[pidIndex][i] = new TH2F (Form("h%sESumEtaPhiE_MC%s",pidParticle[pidIndex].Data(),pnamess[i].Data()),
06373cc6 837 Form("cluster from %s : %s like, #delta^{2}_{#eta #phi} = #Sigma w_{i}(#phi_{i} #eta_{i} ) / #Sigma w_{i} - <#phi><#eta> vs E",ptypess[i].Data(),pidParticle[pidIndex].Data()),
34c16486 838 nptbins,ptmin,ptmax, 2*ssbins,-ssmax,ssmax);
839 fhMCESumEtaPhi[pidIndex][i]->SetXTitle("E (GeV)");
06373cc6 840 fhMCESumEtaPhi[pidIndex][i]->SetYTitle("#delta^{2}_{#eta #phi}");
34c16486 841 outputContainer->Add(fhMCESumEtaPhi[pidIndex][i]);
842
843 fhMCEDispEtaPhiDiff[pidIndex][i] = new TH2F (Form("h%sEDispEtaPhiDiffE_MC%s",pidParticle[pidIndex].Data(),pnamess[i].Data()),
844 Form("cluster from %s : %s like, #sigma^{2}_{#phi #phi} - #sigma^{2}_{#eta #eta} vs E",ptypess[i].Data(),pidParticle[pidIndex].Data()),
845 nptbins,ptmin,ptmax,200,-10,10);
846 fhMCEDispEtaPhiDiff[pidIndex][i]->SetXTitle("E (GeV)");
847 fhMCEDispEtaPhiDiff[pidIndex][i]->SetYTitle("#sigma^{2}_{#phi #phi}-#sigma^{2}_{#eta #eta}");
848 outputContainer->Add(fhMCEDispEtaPhiDiff[pidIndex][i]);
849
850 fhMCESphericity[pidIndex][i] = new TH2F (Form("h%sESphericity_MC%s",pidParticle[pidIndex].Data(),pnamess[i].Data()),
851 Form("cluster from %s : %s like, (#sigma^{2}_{#phi #phi} - #sigma^{2}_{#eta #eta}) / (#sigma^{2}_{#eta #eta} + #sigma^{2}_{#phi #phi}) vs E",ptypess[i].Data(),pidParticle[pidIndex].Data()),
852 nptbins,ptmin,ptmax, 200,-1,1);
853 fhMCESphericity[pidIndex][i]->SetXTitle("E (GeV)");
854 fhMCESphericity[pidIndex][i]->SetYTitle("s = (#sigma^{2}_{#phi #phi} - #sigma^{2}_{#eta #eta}) / (#sigma^{2}_{#eta #eta} + #sigma^{2}_{#phi #phi})");
855 outputContainer->Add(fhMCESphericity[pidIndex][i]);
34c16486 856 }
857
d9105d92 858 }// loop
859 }
860 }
861
764ab1f4 862 //if(IsCaloPIDOn() && pidIndex > 0) continue;
d9105d92 863
864 fhNCellsE[pidIndex] = new TH2F (Form("h%sNCellsE",pidParticle[pidIndex].Data()),
865 Form("N cells in %s cluster vs E ",pidParticle[pidIndex].Data()),
866 nptbins,ptmin,ptmax, nbins,nmin,nmax);
867 fhNCellsE[pidIndex]->SetXTitle("E (GeV)");
868 fhNCellsE[pidIndex]->SetYTitle("# of cells in cluster");
869 outputContainer->Add(fhNCellsE[pidIndex]);
870
7e9a1194 871 fhNLME[pidIndex] = new TH2F (Form("h%sNLME",pidParticle[pidIndex].Data()),
872 Form("NLM in %s cluster vs E ",pidParticle[pidIndex].Data()),
873 nptbins,ptmin,ptmax, 10,0,10);
874 fhNLME[pidIndex]->SetXTitle("E (GeV)");
875 fhNLME[pidIndex]->SetYTitle("# of cells in cluster");
876 outputContainer->Add(fhNLME[pidIndex]);
877
42d47cb7 878 fhTimeE[pidIndex] = new TH2F(Form("h%sTimeE",pidParticle[pidIndex].Data()),
879 Form("Time in %s cluster vs E ",pidParticle[pidIndex].Data())
880 ,nptbins,ptmin,ptmax, tbins,tmin,tmax);
881 fhTimeE[pidIndex]->SetXTitle("E (GeV)");
882 fhTimeE[pidIndex]->SetYTitle(" t (ns)");
883 outputContainer->Add(fhTimeE[pidIndex]);
884
d9105d92 885 fhMaxCellDiffClusterE[pidIndex] = new TH2F (Form("h%sMaxCellDiffClusterE",pidParticle[pidIndex].Data()),
886 Form("%s: energy vs difference of cluster energy - max cell energy / cluster energy, good clusters",pidParticle[pidIndex].Data()),
887 nptbins,ptmin,ptmax, 500,0,1.);
888 fhMaxCellDiffClusterE[pidIndex]->SetXTitle("E_{cluster} (GeV) ");
889 fhMaxCellDiffClusterE[pidIndex]->SetYTitle("(E_{cluster} - E_{cell max})/ E_{cluster}");
890 outputContainer->Add(fhMaxCellDiffClusterE[pidIndex]);
891
892 fhE[pidIndex] = new TH1F(Form("h%sE",pidParticle[pidIndex].Data()),
893 Form("Number of %s over calorimeter vs energy",pidParticle[pidIndex].Data()),
894 nptbins,ptmin,ptmax);
895 fhE[pidIndex]->SetYTitle("N");
896 fhE[pidIndex]->SetXTitle("E_{#gamma}(GeV)");
897 outputContainer->Add(fhE[pidIndex]) ;
898
899 fhPt[pidIndex] = new TH1F(Form("h%sPtElectron",pidParticle[pidIndex].Data()),
900 Form("Number of %s over calorimeter vs p_{T}",pidParticle[pidIndex].Data()),
901 nptbins,ptmin,ptmax);
902 fhPt[pidIndex]->SetYTitle("N");
903 fhPt[pidIndex]->SetXTitle("p_{T #gamma}(GeV/c)");
904 outputContainer->Add(fhPt[pidIndex]) ;
905
906 fhPhi[pidIndex] = new TH2F(Form("h%sPhiElectron",pidParticle[pidIndex].Data()),
907 Form("%s: #phi vs p_{T}",pidParticle[pidIndex].Data()),
908 nptbins,ptmin,ptmax,nphibins,phimin,phimax);
909 fhPhi[pidIndex]->SetYTitle("#phi (rad)");
910 fhPhi[pidIndex]->SetXTitle("p_{T #gamma} (GeV/c)");
911 outputContainer->Add(fhPhi[pidIndex]) ;
912
913 fhEta[pidIndex] = new TH2F(Form("h%sEta",pidParticle[pidIndex].Data()),
914 Form("%s: #eta vs p_{T}",pidParticle[pidIndex].Data()),
915 nptbins,ptmin,ptmax,netabins,etamin,etamax);
916 fhEta[pidIndex]->SetYTitle("#eta");
917 fhEta[pidIndex]->SetXTitle("p_{T #gamma} (GeV/c)");
918 outputContainer->Add(fhEta[pidIndex]) ;
919
920 fhEtaPhi[pidIndex] = new TH2F(Form("h%sEtaPhi",pidParticle[pidIndex].Data()),
921 Form("%s: #eta vs #phi",pidParticle[pidIndex].Data()),
922 netabins,etamin,etamax,nphibins,phimin,phimax);
923 fhEtaPhi[pidIndex]->SetYTitle("#phi (rad)");
924 fhEtaPhi[pidIndex]->SetXTitle("#eta");
925 outputContainer->Add(fhEtaPhi[pidIndex]) ;
34c16486 926 if(GetMinPt() < 0.5)
927 {
d9105d92 928 fhEtaPhi05[pidIndex] = new TH2F(Form("h%sEtaPhi05",pidParticle[pidIndex].Data()),
929 Form("%s: #eta vs #phi, E > 0.5",pidParticle[pidIndex].Data()),
930 netabins,etamin,etamax,nphibins,phimin,phimax);
931 fhEtaPhi05[pidIndex]->SetYTitle("#phi (rad)");
932 fhEtaPhi05[pidIndex]->SetXTitle("#eta");
933 outputContainer->Add(fhEtaPhi05[pidIndex]) ;
934 }
935
936
34c16486 937 if(IsDataMC())
7e9a1194 938 {
34c16486 939 for(Int_t i = 0; i < fNOriginHistograms; i++)
940 {
d9105d92 941 fhMCE[pidIndex][i] = new TH1F(Form("h%sE_MC%s",pidParticle[pidIndex].Data(),pname[i].Data()),
942 Form("%s like cluster from %s : E ",pidParticle[pidIndex].Data(),ptype[i].Data()),
943 nptbins,ptmin,ptmax);
944 fhMCE[pidIndex][i]->SetXTitle("E (GeV)");
945 outputContainer->Add(fhMCE[pidIndex][i]) ;
946
947 fhMCPt[pidIndex][i] = new TH1F(Form("h%sPt_MC%s",pidParticle[pidIndex].Data(),pname[i].Data()),
948 Form("%s like cluster from %s : p_{T} ",pidParticle[pidIndex].Data(),ptype[i].Data()),
949 nptbins,ptmin,ptmax);
950 fhMCPt[pidIndex][i]->SetXTitle("p_{T} (GeV/c)");
951 outputContainer->Add(fhMCPt[pidIndex][i]) ;
952
953 fhMCEta[pidIndex][i] = new TH2F(Form("h%sEta_MC%s",pidParticle[pidIndex].Data(),pname[i].Data()),
954 Form("%s like cluster from %s : #eta ",pidParticle[pidIndex].Data(),ptype[i].Data()),
955 nptbins,ptmin,ptmax,netabins,etamin,etamax);
956 fhMCEta[pidIndex][i]->SetYTitle("#eta");
957 fhMCEta[pidIndex][i]->SetXTitle("E (GeV)");
958 outputContainer->Add(fhMCEta[pidIndex][i]) ;
959
960 fhMCPhi[pidIndex][i] = new TH2F(Form("h%sPhi_MC%s",pidParticle[pidIndex].Data(),pname[i].Data()),
961 Form("%s like cluster from %s : #phi ",pidParticle[pidIndex].Data(),ptype[i].Data()),
962 nptbins,ptmin,ptmax,nphibins,phimin,phimax);
963 fhMCPhi[pidIndex][i]->SetYTitle("#phi (rad)");
964 fhMCPhi[pidIndex][i]->SetXTitle("E (GeV)");
965 outputContainer->Add(fhMCPhi[pidIndex][i]) ;
966
967
968 fhMCDeltaE[pidIndex][i] = new TH2F (Form("h%sDeltaE_MC%s",pidParticle[pidIndex].Data(),pname[i].Data()),
969 Form("%s like MC - Reco E from %s",pidParticle[pidIndex].Data(),pname[i].Data()),
970 nptbins,ptmin,ptmax, 200,-50,50);
971 fhMCDeltaE[pidIndex][i]->SetXTitle("#Delta E (GeV)");
972 outputContainer->Add(fhMCDeltaE[pidIndex][i]);
973
974 fhMC2E[pidIndex][i] = new TH2F (Form("h%s2E_MC%s",pidParticle[pidIndex].Data(),pname[i].Data()),
975 Form("%s like E distribution, reconstructed vs generated from %s",pidParticle[pidIndex].Data(),pname[i].Data()),
976 nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
977 fhMC2E[pidIndex][i]->SetXTitle("E_{rec} (GeV)");
978 fhMC2E[pidIndex][i]->SetYTitle("E_{gen} (GeV)");
979 outputContainer->Add(fhMC2E[pidIndex][i]);
980
981 }
982 } // MC
983
984 }// pid Index
985
986
34c16486 987 if(fFillSSHistograms)
988 {
989 if(IsDataMC())
990 {
d9105d92 991 if(!GetReader()->IsEmbeddedClusterSelectionOn())
992 {
993 fhMCElectronELambda0NoOverlap = new TH2F("hELambda0_MCElectron_NoOverlap",
994 "cluster from Electron : E vs #lambda_{0}^{2}",
995 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
996 fhMCElectronELambda0NoOverlap->SetYTitle("#lambda_{0}^{2}");
997 fhMCElectronELambda0NoOverlap->SetXTitle("E (GeV)");
998 outputContainer->Add(fhMCElectronELambda0NoOverlap) ;
999
1000 fhMCElectronELambda0TwoOverlap = new TH2F("hELambda0_MCElectron_TwoOverlap",
1001 "cluster from Electron : E vs #lambda_{0}^{2}",
1002 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
1003 fhMCElectronELambda0TwoOverlap->SetYTitle("#lambda_{0}^{2}");
1004 fhMCElectronELambda0TwoOverlap->SetXTitle("E (GeV)");
1005 outputContainer->Add(fhMCElectronELambda0TwoOverlap) ;
1006
1007 fhMCElectronELambda0NOverlap = new TH2F("hELambda0_MCElectron_NOverlap",
1008 "cluster from Electron : E vs #lambda_{0}^{2}",
1009 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
1010 fhMCElectronELambda0NOverlap->SetYTitle("#lambda_{0}^{2}");
1011 fhMCElectronELambda0NOverlap->SetXTitle("E (GeV)");
1012 outputContainer->Add(fhMCElectronELambda0NOverlap) ;
1013
1014 } //No embedding
1015
1016 //Fill histograms to check shape of embedded clusters
1017 if(GetReader()->IsEmbeddedClusterSelectionOn())
1018 {
1019
1020 fhEmbeddedSignalFractionEnergy = new TH2F("hEmbeddedSignal_FractionEnergy",
1021 "Energy Fraction of embedded signal versus cluster energy",
1022 nptbins,ptmin,ptmax,100,0.,1.);
1023 fhEmbeddedSignalFractionEnergy->SetYTitle("Fraction");
1024 fhEmbeddedSignalFractionEnergy->SetXTitle("E (GeV)");
1025 outputContainer->Add(fhEmbeddedSignalFractionEnergy) ;
1026
1027 fhEmbedElectronELambda0FullSignal = new TH2F("hELambda0_EmbedElectron_FullSignal",
1028 "cluster from Electron embedded with more than 90% energy in cluster : E vs #lambda_{0}^{2}",
1029 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
1030 fhEmbedElectronELambda0FullSignal->SetYTitle("#lambda_{0}^{2}");
1031 fhEmbedElectronELambda0FullSignal->SetXTitle("E (GeV)");
1032 outputContainer->Add(fhEmbedElectronELambda0FullSignal) ;
1033
1034 fhEmbedElectronELambda0MostlySignal = new TH2F("hELambda0_EmbedElectron_MostlySignal",
1035 "cluster from Electron embedded with 50% to 90% energy in cluster : E vs #lambda_{0}^{2}",
1036 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
1037 fhEmbedElectronELambda0MostlySignal->SetYTitle("#lambda_{0}^{2}");
1038 fhEmbedElectronELambda0MostlySignal->SetXTitle("E (GeV)");
1039 outputContainer->Add(fhEmbedElectronELambda0MostlySignal) ;
1040
1041 fhEmbedElectronELambda0MostlyBkg = new TH2F("hELambda0_EmbedElectron_MostlyBkg",
1042 "cluster from Electron embedded with 10% to 50% energy in cluster : E vs #lambda_{0}^{2}",
1043 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
1044 fhEmbedElectronELambda0MostlyBkg->SetYTitle("#lambda_{0}^{2}");
1045 fhEmbedElectronELambda0MostlyBkg->SetXTitle("E (GeV)");
1046 outputContainer->Add(fhEmbedElectronELambda0MostlyBkg) ;
1047
1048 fhEmbedElectronELambda0FullBkg = new TH2F("hELambda0_EmbedElectron_FullBkg",
1049 "cluster from Electronm embedded with 0% to 10% energy in cluster : E vs #lambda_{0}^{2}",
1050 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
1051 fhEmbedElectronELambda0FullBkg->SetYTitle("#lambda_{0}^{2}");
1052 fhEmbedElectronELambda0FullBkg->SetXTitle("E (GeV)");
1053 outputContainer->Add(fhEmbedElectronELambda0FullBkg) ;
1054
1055
1056 }// embedded histograms
1057
1058 }//Histos with MC
1059
1060 }// Fill SS MC histograms
1061
d9105d92 1062 return outputContainer ;
1063
1064}
1065
78a28af3 1066//_________________________
d9105d92 1067void AliAnaElectron::Init()
1068{
1069
1070 //Init
1071 //Do some checks
1072 if(fCalorimeter == "PHOS" && !GetReader()->IsPHOSSwitchedOn() && NewOutputAOD()){
1073 printf("AliAnaElectron::Init() - !!STOP: You want to use PHOS in analysis but it is not read!! \n!!Check the configuration file!!\n");
1074 abort();
1075 }
1076 else if(fCalorimeter == "EMCAL" && !GetReader()->IsEMCALSwitchedOn() && NewOutputAOD()){
1077 printf("AliAnaElectron::Init() - !!STOP: You want to use EMCAL in analysis but it is not read!! \n!!Check the configuration file!!\n");
1078 abort();
1079 }
1080
1081}
1082
78a28af3 1083//___________________________________
d9105d92 1084void AliAnaElectron::InitParameters()
1085{
1086
1087 //Initialize the parameters of the analysis.
1088 AddToHistogramsName("AnaElectron_");
1089
1090 fCalorimeter = "EMCAL" ;
1091 fMinDist = 2.;
1092 fMinDist2 = 4.;
1093 fMinDist3 = 5.;
1094
1095 fTimeCutMin = -1;
1096 fTimeCutMax = 9999999;
1097 fNCellsCut = 0;
1098
1099 fdEdxMin = 76.; // for LHC11a, but for LHC11c pass1 56.
1100 fdEdxMax = 85.; // for LHC11a, but for LHC11c pass1 64.
1101
1102 fEOverPMin = 0.8; // for LHC11a, but for LHC11c pass1 0.9
1103 fEOverPMax = 1.2; // for LHC11a and LHC11c pass1
1104
1105}
1106
78a28af3 1107//_________________________________________
d9105d92 1108void AliAnaElectron::MakeAnalysisFillAOD()
1109{
1110 //Do photon analysis and fill aods
1111
1112 //Get the vertex
1113 Double_t v[3] = {0,0,0}; //vertex ;
1114 GetReader()->GetVertex(v);
1115
1116 //Select the Calorimeter of the photon
1117 TObjArray * pl = 0x0;
1118 if(fCalorimeter == "PHOS")
1119 pl = GetPHOSClusters();
1120 else if (fCalorimeter == "EMCAL")
1121 pl = GetEMCALClusters();
1122
1123 if(!pl) {
1124 Info("MakeAnalysisFillAOD","TObjArray with %s clusters is NULL!\n",fCalorimeter.Data());
1125 return;
1126 }
1127
1128 //Init arrays, variables, get number of clusters
1129 TLorentzVector mom, mom2 ;
1130 Int_t nCaloClusters = pl->GetEntriesFast();
1131 //List to be used in conversion analysis, to tag the cluster as candidate for conversion
1132
1133 if(GetDebug() > 0) printf("AliAnaElectron::MakeAnalysisFillAOD() - input %s cluster entries %d\n", fCalorimeter.Data(), nCaloClusters);
1134
1135 //----------------------------------------------------
1136 // Fill AOD with PHOS/EMCAL AliAODPWG4Particle objects
1137 //----------------------------------------------------
1138 // Loop on clusters
1139 for(Int_t icalo = 0; icalo < nCaloClusters; icalo++){
1140
1141 AliVCluster * calo = (AliVCluster*) (pl->At(icalo));
1142 //printf("calo %d, %f\n",icalo,calo->E());
1143
1144 //Get the index where the cluster comes, to retrieve the corresponding vertex
1145 Int_t evtIndex = 0 ;
1146 if (GetMixedEvent()) {
1147 evtIndex=GetMixedEvent()->EventIndexForCaloCluster(calo->GetID()) ;
1148 //Get the vertex and check it is not too large in z
1149 if(TMath::Abs(GetVertex(evtIndex)[2])> GetZvertexCut()) continue;
1150 }
1151
1152 //Cluster selection, not charged, with photon id and in fiducial cut
1153 if(GetReader()->GetDataType() != AliCaloTrackReader::kMC){
1154 calo->GetMomentum(mom,GetVertex(evtIndex)) ;}//Assume that come from vertex in straight line
1155 else{
1156 Double_t vertex[]={0,0,0};
1157 calo->GetMomentum(mom,vertex) ;
1158 }
1159
1160 //--------------------------------------
1161 // Cluster selection
1162 //--------------------------------------
7e9a1194 1163 AliVCaloCells* cells = 0;
1164 if(fCalorimeter == "EMCAL") cells = GetEMCALCells();
1165 else cells = GetPHOSCells();
1166
1167 Int_t nMaxima = GetCaloUtils()->GetNumberOfLocalMaxima(calo, cells); // NLM
1168 if(!ClusterSelected(calo,mom,nMaxima)) continue;
d9105d92 1169
1170 //----------------------------
1171 //Create AOD for analysis
1172 //----------------------------
764ab1f4 1173 AliAODPWG4Particle aodpart = AliAODPWG4Particle(mom);
d9105d92 1174
1175 //...............................................
1176 //Set the indeces of the original caloclusters (MC, ID), and calorimeter
1177 Int_t label = calo->GetLabel();
764ab1f4 1178 aodpart.SetLabel(label);
1179 aodpart.SetCaloLabel(calo->GetID(),-1);
1180 aodpart.SetDetector(fCalorimeter);
d9105d92 1181 //printf("Index %d, Id %d, iaod %d\n",icalo, calo->GetID(),GetOutputAODBranch()->GetEntriesFast());
1182
1183 //...............................................
1184 //Set bad channel distance bit
1185 Double_t distBad=calo->GetDistanceToBadChannel() ; //Distance to bad channel
764ab1f4 1186 if (distBad > fMinDist3) aodpart.SetDistToBad(2) ;
1187 else if(distBad > fMinDist2) aodpart.SetDistToBad(1) ;
1188 else aodpart.SetDistToBad(0) ;
1189 //printf("DistBad %f Bit %d\n",distBad, aodpart.DistToBad());
d9105d92 1190
d9105d92 1191 //-------------------------------------
78a28af3 1192 //PID selection via dEdx
d9105d92 1193 //-------------------------------------
78a28af3 1194
4bfeae64 1195 AliVTrack *track = GetCaloUtils()->GetMatchedTrack(calo, GetReader()->GetInputEvent());
1196
7e9a1194 1197 if(!track)
1198 {
d9105d92 1199 printf("AliAnaElectron::MakeAnalysisFillAOD() - Null track");
1200 continue;
1201 }
1202
7e9a1194 1203 //printf("track dedx %f, p %f, cluster E %f\n",track->GetTPCsignal(),track->P(),calo->E());
e416be7d 1204 Float_t dEdx = track->GetTPCsignal();
7e9a1194 1205 Float_t eOverp = calo->E()/track->P();
1206
d9105d92 1207 fhdEdxvsE->Fill(calo->E(), dEdx);
1208 fhdEdxvsP->Fill(track->P(),dEdx);
1209
7e9a1194 1210 if( eOverp < fEOverPMax && eOverp > fEOverPMin)
1211 {
1212 fhdEdxvsECutEOverP ->Fill(calo->E(), dEdx);
1213 fhdEdxvsPCutEOverP ->Fill(track->P(),dEdx);
1214 }
1215
1216 //Apply a mild cut on the cluster SS and check the value of dEdX and EOverP
1217 Float_t m02 = calo->GetM02();
1218 if(m02 > 0.1 && m02 < 0.4)
1219 {
1220 fhdEdxvsECutM02 ->Fill(calo->E(), dEdx);
1221 fhdEdxvsPCutM02 ->Fill(track->P(),dEdx);
1222 fhEOverPvsECutM02->Fill(calo->E(), eOverp);
1223 fhEOverPvsPCutM02->Fill(track->P(), eOverp);
1224 }
1225
c5693f62 1226 Int_t pid = AliCaloPID::kChargedHadron;
1227
7e9a1194 1228 if( dEdx < fdEdxMax && dEdx > fdEdxMin)
1229 {
e416be7d 1230 fhEOverPvsE->Fill(calo->E(), eOverp);
1231 fhEOverPvsP->Fill(track->P(), eOverp);
d9105d92 1232
7e9a1194 1233 if(m02 > 0.1 && m02 < 0.4)
1234 {
1235 fhEOverPvsECutM02CutdEdx->Fill(calo->E(), eOverp);
1236 fhEOverPvsPCutM02CutdEdx->Fill(track->P(), eOverp);
1237 }
1238
1239 if( eOverp < fEOverPMax && eOverp > fEOverPMin)
1240 {
d9105d92 1241 pid = AliCaloPID::kElectron;
d9105d92 1242 } //E/p
1243
1244 }// dEdx
1245
764ab1f4 1246 aodpart.SetIdentifiedParticleType(pid);
d9105d92 1247
1248 Int_t pidIndex = 0;// Electron
dbf54f1e 1249 if(pid == AliCaloPID::kChargedHadron) pidIndex = 1;
7e9a1194 1250
1251 //--------------------------------------------------------------------------------------
1252 //Play with the MC stack if available
1253 //--------------------------------------------------------------------------------------
1254
1255 //Check origin of the candidates
1256 if(IsDataMC())
1257 {
1258 Int_t tag = GetMCAnalysisUtils()->CheckOrigin(calo->GetLabels(),calo->GetNLabels(),GetReader());
1259 aodpart.SetTag(tag);
1260
1261 if(GetDebug() > 0)
1262 printf("AliAnaElectron::MakeAnalysisFillAOD() - Origin of candidate, bit map %d\n",aodpart.GetTag());
1263
1264 if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPhoton) && fhMCE[pidIndex][kmcPhoton])
1265 {
1266 fhMCdEdxvsE [kmcPhoton]->Fill(calo ->E(), dEdx);
1267 fhMCdEdxvsP [kmcPhoton]->Fill(track->P(), dEdx);
1268 fhMCEOverPvsE[kmcPhoton]->Fill(calo ->E(), eOverp);
1269 fhMCEOverPvsP[kmcPhoton]->Fill(track->P(), eOverp);
1270
1271 if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCConversion) && fhMCE[pidIndex][kmcConversion])
1272 {
1273 fhMCdEdxvsE [kmcConversion]->Fill(calo ->E(), dEdx);
1274 fhMCdEdxvsP [kmcConversion]->Fill(track->P(), dEdx);
1275 fhMCEOverPvsE[kmcConversion]->Fill(calo ->E(), eOverp);
1276 fhMCEOverPvsP[kmcConversion]->Fill(track->P(), eOverp);
1277 }
1278 else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0Decay) &&
1279 !GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0) && fhMCE[pidIndex][kmcPi0Decay])
1280 {
1281 fhMCdEdxvsE [kmcPi0Decay]->Fill(calo ->E(), dEdx);
1282 fhMCdEdxvsP [kmcPi0Decay]->Fill(track->P(), dEdx);
1283 fhMCEOverPvsE[kmcPi0Decay]->Fill(calo ->E(), eOverp);
1284 fhMCEOverPvsP[kmcPi0Decay]->Fill(track->P(), eOverp);
1285 }
1286 else if( (GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEtaDecay) ||
1287 GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCOtherDecay) ) && fhMCE[pidIndex][kmcOtherDecay])
1288 {
1289 fhMCdEdxvsE [kmcOtherDecay]->Fill(calo ->E(), dEdx);
1290 fhMCdEdxvsP [kmcOtherDecay]->Fill(track->P(), dEdx);
1291 fhMCEOverPvsE[kmcOtherDecay]->Fill(calo ->E(), eOverp);
1292 fhMCEOverPvsP[kmcOtherDecay]->Fill(track->P(), eOverp);
1293 }
1294 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0) && fhMCE [pidIndex][kmcPi0])
1295 {
1296 fhMCdEdxvsE [kmcPi0]->Fill(calo ->E(), dEdx);
1297 fhMCdEdxvsP [kmcPi0]->Fill(track->P(), dEdx);
1298 fhMCEOverPvsE[kmcPi0]->Fill(calo ->E(), eOverp);
1299 fhMCEOverPvsP[kmcPi0]->Fill(track->P(), eOverp);
1300 }
1301 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEta) && fhMCE[pidIndex][kmcEta])
1302 {
1303 fhMCdEdxvsE [kmcEta]->Fill(calo ->E(), dEdx);
1304 fhMCdEdxvsP [kmcEta]->Fill(track->P(), dEdx);
1305 fhMCEOverPvsE[kmcEta]->Fill(calo ->E(), eOverp);
1306 fhMCEOverPvsP[kmcEta]->Fill(track->P(), eOverp);
1307 }
1308 }
1309 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCAntiNeutron) && fhMCE[pidIndex][kmcAntiNeutron])
1310 {
1311 fhMCdEdxvsE [kmcAntiNeutron]->Fill(calo ->E(), dEdx);
1312 fhMCdEdxvsP [kmcAntiNeutron]->Fill(track->P(), dEdx);
1313 fhMCEOverPvsE[kmcAntiNeutron]->Fill(calo ->E(), eOverp);
1314 fhMCEOverPvsP[kmcAntiNeutron]->Fill(track->P(), eOverp);
1315 }
1316 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCAntiProton) && fhMCE[pidIndex][kmcAntiProton])
1317 {
1318 fhMCdEdxvsE [kmcAntiProton]->Fill(calo ->E(), dEdx);
1319 fhMCdEdxvsP [kmcAntiProton]->Fill(track->P(), dEdx);
1320 fhMCEOverPvsE[kmcAntiProton]->Fill(calo ->E(), eOverp);
1321 fhMCEOverPvsP[kmcAntiProton]->Fill(track->P(), eOverp);
1322 }
1323 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCElectron) && fhMCE[pidIndex][kmcElectron])
1324 {
1325 fhMCdEdxvsE [kmcElectron]->Fill(calo ->E(), dEdx);
1326 fhMCdEdxvsP [kmcElectron]->Fill(track->P(), dEdx);
1327 fhMCEOverPvsE[kmcElectron]->Fill(calo ->E(), eOverp);
1328 fhMCEOverPvsP[kmcElectron]->Fill(track->P(), eOverp);
1329 }
1330 else if( fhMCE[pidIndex][kmcOther])
1331 {
1332 fhMCdEdxvsE [kmcOther]->Fill(calo ->E(), dEdx);
1333 fhMCdEdxvsP [kmcOther]->Fill(track->P(), dEdx);
1334 fhMCEOverPvsE[kmcOther]->Fill(calo ->E(), eOverp);
1335 fhMCEOverPvsP[kmcOther]->Fill(track->P(), eOverp);
1336 }
1337 }// set MC tag and fill Histograms with MC
d9105d92 1338
78a28af3 1339 //---------------------------------
d9105d92 1340 //Fill some shower shape histograms
78a28af3 1341 //---------------------------------
1a72f6c5 1342
764ab1f4 1343 FillShowerShapeHistograms(calo,aodpart.GetTag(),pid);
7e9a1194 1344
78a28af3 1345 if(pid == AliCaloPID::kElectron)
1346 WeightHistograms(calo);
1347
1348 //-----------------------------------------
d9105d92 1349 //PID Shower Shape selection or bit setting
78a28af3 1350 //-----------------------------------------
1351
d9105d92 1352 // Data, PID check on
3c1d9afb 1353 if(IsCaloPIDOn())
1354 {
49b5c49b 1355 // Get most probable PID, 2 options check bayesian PID weights or redo PID
1356 // By default, redo PID
d9105d92 1357
764ab1f4 1358 if(GetCaloPID()->GetIdentifiedParticleType(calo)!=AliCaloPID::kPhoton)
1359 {
1360 continue;
1361 }
1362 if(GetDebug() > 1) printf("AliAnaPhoton::MakeAnalysisFillAOD() - PDG of identified particle %d\n",aodpart.GetIdentifiedParticleType());
d9105d92 1363
1364 }
49b5c49b 1365
3c1d9afb 1366 if(GetDebug() > 1) printf("AliAnaElectron::MakeAnalysisFillAOD() - Photon selection cuts passed: pT %3.2f, pdg %d\n",
764ab1f4 1367 aodpart.Pt(), aodpart.GetIdentifiedParticleType());
d9105d92 1368
d9105d92 1369 //FIXME, this to MakeAnalysisFillHistograms ...
1370 Int_t absID = 0;
1371 Float_t maxCellFraction = 0;
d9105d92 1372
1373 absID = GetCaloUtils()->GetMaxEnergyCell(cells, calo,maxCellFraction);
764ab1f4 1374 fhMaxCellDiffClusterE[pidIndex]->Fill(aodpart.E(),maxCellFraction);
1375 fhNCellsE[pidIndex] ->Fill(aodpart.E(),calo->GetNCells());
7e9a1194 1376 fhNLME[pidIndex] ->Fill(aodpart.E(),nMaxima);
764ab1f4 1377 fhTimeE[pidIndex] ->Fill(aodpart.E(),calo->GetTOF()*1.e9);
1378
1379 //Add AOD with electron/hadron object to aod branch
1380 if ( pid == fAODParticle || fAODParticle == 0 )
1381 {
1382 AddAODParticle(aodpart);
1383 }
42d47cb7 1384
d9105d92 1385 }//loop
1386
1387 if(GetDebug() > 1) printf("AliAnaElectron::MakeAnalysisFillAOD() End fill AODs, with %d entries \n",GetOutputAODBranch()->GetEntriesFast());
1388
1389}
1390
78a28af3 1391//________________________________________________
d9105d92 1392void AliAnaElectron::MakeAnalysisFillHistograms()
1393{
1394 //Fill histograms
1395
1396 //-------------------------------------------------------------------
1397 // Access MC information in stack if requested, check that it exists.
1398 AliStack * stack = 0x0;
1399 TParticle * primary = 0x0;
1400 TClonesArray * mcparticles = 0x0;
1401 AliAODMCParticle * aodprimary = 0x0;
1402
1403 if(IsDataMC()){
1404
1405 if(GetReader()->ReadStack()){
1406 stack = GetMCStack() ;
1407 if(!stack) {
1408 printf("AliAnaElectron::MakeAnalysisFillHistograms() - Stack not available, is the MC handler called? STOP\n");
1409 abort();
1410 }
1411
1412 }
1413 else if(GetReader()->ReadAODMCParticles()){
1414
1415 //Get the list of MC particles
2644ead9 1416 mcparticles = GetReader()->GetAODMCParticles();
d9105d92 1417 if(!mcparticles && GetDebug() > 0) {
1418 printf("AliAnaElectron::MakeAnalysisFillHistograms() - Standard MCParticles not available!\n");
1419 }
1420 }
1421 }// is data and MC
1422
1423
1424 // Get vertex
1425 Double_t v[3] = {0,0,0}; //vertex ;
1426 GetReader()->GetVertex(v);
1427 //fhVertex->Fill(v[0],v[1],v[2]);
1428 if(TMath::Abs(v[2]) > GetZvertexCut()) return ; // done elsewhere for Single Event analysis, but there for mixed event
1429
1430 //----------------------------------
1431 //Loop on stored AOD photons
1432 Int_t naod = GetOutputAODBranch()->GetEntriesFast();
1433 if(GetDebug() > 0) printf("AliAnaElectron::MakeAnalysisFillHistograms() - aod branch entries %d\n", naod);
1434
3c1d9afb 1435 for(Int_t iaod = 0; iaod < naod ; iaod++)
1436 {
d9105d92 1437 AliAODPWG4Particle* ph = (AliAODPWG4Particle*) (GetOutputAODBranch()->At(iaod));
1438 Int_t pdg = ph->GetIdentifiedParticleType();
1439
1440 Int_t pidIndex = 0;// Electron
1441 if (pdg == AliCaloPID::kElectron) pidIndex = 0;
1442 else if(pdg == AliCaloPID::kChargedHadron) pidIndex = 1;
1443 else continue ;
1444
1445 if(ph->GetDetector() != fCalorimeter) continue;
1446
1447 if(GetDebug() > 2)
1448 printf("AliAnaElectron::MakeAnalysisFillHistograms() - ID Electron: pt %f, phi %f, eta %f\n", ph->Pt(),ph->Phi(),ph->Eta()) ;
1449
1450 //................................
1451 //Fill photon histograms
1452 Float_t ptcluster = ph->Pt();
1453 Float_t phicluster = ph->Phi();
1454 Float_t etacluster = ph->Eta();
1455 Float_t ecluster = ph->E();
1456
1457 fhE[pidIndex] ->Fill(ecluster);
1458 fhPt[pidIndex] ->Fill(ptcluster);
1459 fhPhi[pidIndex] ->Fill(ptcluster,phicluster);
1460 fhEta[pidIndex] ->Fill(ptcluster,etacluster);
1461 if (ecluster > 0.5) fhEtaPhi[pidIndex] ->Fill(etacluster, phicluster);
1462 else if(GetMinPt() < 0.5) fhEtaPhi05[pidIndex]->Fill(etacluster, phicluster);
1463
1464 //.......................................
1465 //Play with the MC data if available
1466 if(IsDataMC()){
1467
1468 //....................................................................
1469 // Access MC information in stack if requested, check that it exists.
1470 Int_t label =ph->GetLabel();
1471 if(label < 0) {
1472 if(GetDebug() > 1) printf("AliAnaElectron::MakeAnalysisFillHistograms() *** bad label ***: label %d \n", label);
1473 continue;
1474 }
1475
1476 Float_t eprim = 0;
1477 Float_t ptprim = 0;
1478 if(GetReader()->ReadStack()){
1479
1480 if(label >= stack->GetNtrack()) {
1481 if(GetDebug() > 2) printf("AliAnaElectron::MakeAnalysisFillHistograms() *** large label ***: label %d, n tracks %d \n", label, stack->GetNtrack());
1482 continue ;
1483 }
1484
1485 primary = stack->Particle(label);
1486 if(!primary){
1487 printf("AliAnaElectron::MakeAnalysisFillHistograms() *** no primary ***: label %d \n", label);
1488 continue;
1489 }
1490
1491 eprim = primary->Energy();
1492 ptprim = primary->Pt();
1493
1494 }
1495 else if(GetReader()->ReadAODMCParticles()){
1496 //Check which is the input
1497 if(ph->GetInputFileIndex() == 0){
1498 if(!mcparticles) continue;
1499 if(label >= mcparticles->GetEntriesFast()) {
1500 if(GetDebug() > 2) printf("AliAnaElectron::MakeAnalysisFillHistograms() *** large label ***: label %d, n tracks %d \n",
1501 label, mcparticles->GetEntriesFast());
1502 continue ;
1503 }
1504 //Get the particle
1505 aodprimary = (AliAODMCParticle*) mcparticles->At(label);
1506
1507 }
1508
1509 if(!aodprimary){
1510 printf("AliAnaElectron::MakeAnalysisFillHistograms() *** no primary ***: label %d \n", label);
1511 continue;
1512 }
1513
1514 eprim = aodprimary->E();
1515 ptprim = aodprimary->Pt();
1516
1517 }
1518
1519 Int_t tag =ph->GetTag();
1520
c5693f62 1521 if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPhoton) && fhMCE[pidIndex][kmcPhoton])
d9105d92 1522 {
c5693f62 1523 fhMCE [pidIndex][kmcPhoton] ->Fill(ecluster);
1524 fhMCPt [pidIndex][kmcPhoton] ->Fill(ptcluster);
1525 fhMCPhi[pidIndex][kmcPhoton] ->Fill(ecluster,phicluster);
1526 fhMCEta[pidIndex][kmcPhoton] ->Fill(ecluster,etacluster);
d9105d92 1527
c5693f62 1528 fhMC2E[pidIndex][kmcPhoton] ->Fill(ecluster, eprim);
1529 fhMCDeltaE[pidIndex][kmcPhoton] ->Fill(ecluster,eprim-ecluster);
d9105d92 1530
c5693f62 1531 if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCConversion) && fhMCE[pidIndex][kmcConversion])
d9105d92 1532 {
c5693f62 1533 fhMCE [pidIndex][kmcConversion] ->Fill(ecluster);
1534 fhMCPt [pidIndex][kmcConversion] ->Fill(ptcluster);
1535 fhMCPhi[pidIndex][kmcConversion] ->Fill(ecluster,phicluster);
1536 fhMCEta[pidIndex][kmcConversion] ->Fill(ecluster,etacluster);
d9105d92 1537
c5693f62 1538 fhMC2E[pidIndex][kmcConversion] ->Fill(ecluster, eprim);
1539 fhMCDeltaE[pidIndex][kmcConversion] ->Fill(ecluster,eprim-ecluster);
d9105d92 1540
1541 }
1542 else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0Decay) &&
c5693f62 1543 !GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0) && fhMCE[pidIndex][kmcPi0Decay])
d9105d92 1544 {
c5693f62 1545 fhMCE [pidIndex][kmcPi0Decay] ->Fill(ecluster);
1546 fhMCPt [pidIndex][kmcPi0Decay] ->Fill(ptcluster);
1547 fhMCPhi[pidIndex][kmcPi0Decay] ->Fill(ecluster,phicluster);
1548 fhMCEta[pidIndex][kmcPi0Decay] ->Fill(ecluster,etacluster);
d9105d92 1549
c5693f62 1550 fhMC2E[pidIndex][kmcPi0Decay] ->Fill(ecluster, eprim);
1551 fhMCDeltaE[pidIndex][kmcPi0Decay] ->Fill(ecluster,eprim-ecluster);
d9105d92 1552 }
1553 else if( (GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEtaDecay) ||
c5693f62 1554 GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCOtherDecay) ) && fhMCE[pidIndex][kmcOtherDecay])
d9105d92 1555 {
c5693f62 1556 fhMCE [pidIndex][kmcOtherDecay] ->Fill(ecluster);
1557 fhMCPt [pidIndex][kmcOtherDecay] ->Fill(ptcluster);
1558 fhMCPhi[pidIndex][kmcOtherDecay] ->Fill(ecluster,phicluster);
1559 fhMCEta[pidIndex][kmcOtherDecay] ->Fill(ecluster,etacluster);
d9105d92 1560
c5693f62 1561 fhMC2E[pidIndex][kmcOtherDecay] ->Fill(ecluster, eprim);
1562 fhMCDeltaE[pidIndex][kmcOtherDecay] ->Fill(ecluster,eprim-ecluster);
d9105d92 1563 }
c5693f62 1564 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0) && fhMCE [pidIndex][kmcPi0])
d9105d92 1565 {
c5693f62 1566 fhMCE [pidIndex][kmcPi0] ->Fill(ecluster);
1567 fhMCPt [pidIndex][kmcPi0] ->Fill(ptcluster);
1568 fhMCPhi[pidIndex][kmcPi0] ->Fill(ecluster,phicluster);
1569 fhMCEta[pidIndex][kmcPi0] ->Fill(ecluster,etacluster);
d9105d92 1570
c5693f62 1571 fhMC2E[pidIndex][kmcPi0] ->Fill(ecluster, eprim);
1572 fhMCDeltaE[pidIndex][kmcPi0] ->Fill(ecluster,eprim-ecluster);
d9105d92 1573
1574 }
c5693f62 1575 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEta) && fhMCE[pidIndex][kmcEta])
d9105d92 1576 {
c5693f62 1577 fhMCE [pidIndex][kmcEta] ->Fill(ecluster);
1578 fhMCPt [pidIndex][kmcEta] ->Fill(ptcluster);
1579 fhMCPhi[pidIndex][kmcEta] ->Fill(ecluster,phicluster);
1580 fhMCEta[pidIndex][kmcEta] ->Fill(ecluster,etacluster);
d9105d92 1581
c5693f62 1582 fhMC2E[pidIndex][kmcEta] ->Fill(ecluster, eprim);
1583 fhMCDeltaE[pidIndex][kmcEta] ->Fill(ecluster,eprim-ecluster);
d9105d92 1584
1585 }
1586 }
c5693f62 1587 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCAntiNeutron) && fhMCE[pidIndex][kmcAntiNeutron])
d9105d92 1588 {
c5693f62 1589 fhMCE [pidIndex][kmcAntiNeutron] ->Fill(ecluster);
1590 fhMCPt [pidIndex][kmcAntiNeutron] ->Fill(ptcluster);
1591 fhMCPhi[pidIndex][kmcAntiNeutron] ->Fill(ecluster,phicluster);
1592 fhMCEta[pidIndex][kmcAntiNeutron] ->Fill(ecluster,etacluster);
d9105d92 1593
c5693f62 1594 fhMC2E[pidIndex][kmcAntiNeutron] ->Fill(ecluster, eprim);
1595 fhMCDeltaE[pidIndex][kmcAntiNeutron] ->Fill(ecluster,eprim-ecluster);
d9105d92 1596
1597 }
c5693f62 1598 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCAntiProton) && fhMCE[pidIndex][kmcAntiProton])
d9105d92 1599 {
c5693f62 1600 fhMCE [pidIndex][kmcAntiProton] ->Fill(ecluster);
1601 fhMCPt [pidIndex][kmcAntiProton] ->Fill(ptcluster);
1602 fhMCPhi[pidIndex][kmcAntiProton] ->Fill(ecluster,phicluster);
1603 fhMCEta[pidIndex][kmcAntiProton] ->Fill(ecluster,etacluster);
d9105d92 1604
c5693f62 1605 fhMC2E[pidIndex][kmcAntiProton] ->Fill(ecluster, eprim);
1606 fhMCDeltaE[pidIndex][kmcAntiProton] ->Fill(ecluster,eprim-ecluster);
d9105d92 1607
1608 }
c5693f62 1609 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCElectron) && fhMCE[pidIndex][kmcElectron])
d9105d92 1610 {
c5693f62 1611 fhMCE [pidIndex][kmcElectron] ->Fill(ecluster);
1612 fhMCPt [pidIndex][kmcElectron] ->Fill(ptcluster);
1613 fhMCPhi[pidIndex][kmcElectron] ->Fill(ecluster,phicluster);
1614 fhMCEta[pidIndex][kmcElectron] ->Fill(ecluster,etacluster);
d9105d92 1615
c5693f62 1616 fhMC2E[pidIndex][kmcElectron] ->Fill(ecluster, eprim);
1617 fhMCDeltaE[pidIndex][kmcElectron] ->Fill(ecluster,eprim-ecluster);
d9105d92 1618
1619 }
c5693f62 1620 else if( fhMCE[pidIndex][kmcOther]){
1621 fhMCE [pidIndex][kmcOther] ->Fill(ecluster);
1622 fhMCPt [pidIndex][kmcOther] ->Fill(ptcluster);
1623 fhMCPhi[pidIndex][kmcOther] ->Fill(ecluster,phicluster);
1624 fhMCEta[pidIndex][kmcOther] ->Fill(ecluster,etacluster);
d9105d92 1625
c5693f62 1626 fhMC2E[pidIndex][kmcOther] ->Fill(ecluster, eprim);
1627 fhMCDeltaE[pidIndex][kmcOther] ->Fill(ecluster,eprim-ecluster);
d9105d92 1628
1629 }
1630
1631 }//Histograms with MC
1632
1633 }// aod loop
1634
1635}
1636
78a28af3 1637//____________________________________________________
d9105d92 1638void AliAnaElectron::Print(const Option_t * opt) const
1639{
1640 //Print some relevant parameters set for the analysis
1641
1642 if(! opt)
1643 return;
1644
1645 printf("**** Print %s %s ****\n", GetName(), GetTitle() ) ;
745913ae 1646 AliAnaCaloTrackCorrBaseClass::Print(" ");
d9105d92 1647
1648 printf("Calorimeter = %s\n", fCalorimeter.Data()) ;
1649 printf(" %2.2f < dEdx < %2.2f \n",fdEdxMin,fdEdxMax) ;
1650 printf(" %2.2f < E/P < %2.2f \n",fEOverPMin,fEOverPMax) ;
1651 printf("Min Distance to Bad Channel = %2.1f\n",fMinDist);
1652 printf("Min Distance to Bad Channel 2 = %2.1f\n",fMinDist2);
1653 printf("Min Distance to Bad Channel 3 = %2.1f\n",fMinDist3);
1654 printf("Time Cut: %3.1f < TOF < %3.1f\n", fTimeCutMin, fTimeCutMax);
1655 printf("Number of cells in cluster is > %d \n", fNCellsCut);
1656 printf(" \n") ;
1657
1658}
78a28af3 1659
78a28af3 1660//______________________________________________________
1661void AliAnaElectron::WeightHistograms(AliVCluster *clus)
1662{
1663 // Calculate weights and fill histograms
1664
1665 if(!fFillWeightHistograms || GetMixedEvent()) return;
1666
1667 AliVCaloCells* cells = 0;
1668 if(fCalorimeter == "EMCAL") cells = GetEMCALCells();
1669 else cells = GetPHOSCells();
1670
1671 // First recalculate energy in case non linearity was applied
1672 Float_t energy = 0;
1673 Float_t ampMax = 0;
1674 for (Int_t ipos = 0; ipos < clus->GetNCells(); ipos++) {
1675
1676 Int_t id = clus->GetCellsAbsId()[ipos];
1677
1678 //Recalibrate cell energy if needed
1679 Float_t amp = cells->GetCellAmplitude(id);
dbba06ca 1680 GetCaloUtils()->RecalibrateCellAmplitude(amp,fCalorimeter, id);
78a28af3 1681
1682 energy += amp;
1683
1684 if(amp> ampMax)
1685 ampMax = amp;
1686
1687 } // energy loop
1688
1689 if(energy <=0 ) {
1690 printf("AliAnaCalorimeterQA::WeightHistograms()- Wrong calculated energy %f\n",energy);
1691 return;
1692 }
1693
1a72f6c5 1694 //printf("AliAnaElectron::WeightHistograms() - energy %f, ampmax %f, rat %f, lograt %f\n",energy,ampMax,ampMax/energy,TMath::Log(ampMax/energy));
78a28af3 1695 fhEMaxCellClusterRatio ->Fill(energy,ampMax/energy);
1696 fhEMaxCellClusterLogRatio->Fill(energy,TMath::Log(ampMax/energy));
1697
1698 //Get the ratio and log ratio to all cells in cluster
1699 for (Int_t ipos = 0; ipos < clus->GetNCells(); ipos++) {
1700 Int_t id = clus->GetCellsAbsId()[ipos];
1701
1702 //Recalibrate cell energy if needed
1703 Float_t amp = cells->GetCellAmplitude(id);
dbba06ca 1704 GetCaloUtils()->RecalibrateCellAmplitude(amp, fCalorimeter, id);
78a28af3 1705
1706 //printf("energy %f, amp %f, rat %f, lograt %f\n",energy,amp,amp/energy,TMath::Log(amp/energy));
1707 fhECellClusterRatio ->Fill(energy,amp/energy);
1708 fhECellClusterLogRatio->Fill(energy,TMath::Log(amp/energy));
1709 }
1710
1711 //Recalculate shower shape for different W0
1712 if(fCalorimeter=="EMCAL"){
1713
1714 Float_t l0org = clus->GetM02();
1715 Float_t l1org = clus->GetM20();
1716 Float_t dorg = clus->GetDispersion();
1717
1a72f6c5 1718 for(Int_t iw = 0; iw < 14; iw++){
1719
1720 GetCaloUtils()->GetEMCALRecoUtils()->SetW0(1+iw*0.5);
78a28af3 1721 GetCaloUtils()->GetEMCALRecoUtils()->RecalculateClusterShowerShapeParameters(GetEMCALGeometry(), cells, clus);
1722
1723 fhLambda0ForW0[iw]->Fill(energy,clus->GetM02());
1a72f6c5 1724 //fhLambda1ForW0[iw]->Fill(energy,clus->GetM20());
78a28af3 1725
1726 //printf("\t w %1.1f, l0 %f, l1 %f,\n",3+iw*0.5,clus->GetM02(),clus->GetM20());
1727
1728 } // w0 loop
1729
1730 // Set the original values back
1731 clus->SetM02(l0org);
1732 clus->SetM20(l1org);
1733 clus->SetDispersion(dorg);
1734
1735 }// EMCAL
1736}
1737
1738