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