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