]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGGA/CaloTrackCorrelations/AliAnaClusterPileUp.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
add MC histogram for decay photons when pair is lost
[u/mrichter/AliRoot.git] / PWGGA / CaloTrackCorrelations / AliAnaClusterPileUp.cxx
CommitLineData
bc41680b 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 **************************************************************************/
15
16//_________________________________________________________________________
17//
18// Class for the study of Pile-up effect on
19// Calorimeter clusters.
20// Open time cuts in reader.
21//
22//-- Author: Gustavo Conesa (CNRS-LPSC-Grenoble)
23//////////////////////////////////////////////////////////////////////////////
24
25// --- ROOT system ---
26#include <TH2F.h>
27#include <TClonesArray.h>
28#include <TObjString.h>
29
30// --- Analysis system ---
31#include "AliAnaClusterPileUp.h"
32#include "AliCaloTrackReader.h"
33#include "AliFiducialCut.h"
34#include "AliVCluster.h"
35#include "AliAODEvent.h"
36#include "AliESDEvent.h"
37
38ClassImp(AliAnaClusterPileUp)
39
40//___________________________________________
41AliAnaClusterPileUp::AliAnaClusterPileUp() :
42AliAnaCaloTrackCorrBaseClass(),
43fCalorimeter(""), fNCellsCut(0),
44// Histograms
45fhTimePtNoCut(0), fhTimePtSPD(0),
46fhTimeNPileUpVertSPD(0), fhTimeNPileUpVertTrack(0),
47fhTimeNPileUpVertContributors(0),
48fhTimePileUpMainVertexZDistance(0), fhTimePileUpMainVertexZDiamond(0),
49fhClusterMultSPDPileUp(), fhClusterMultNoPileUp(),
50fhEtaPhiBC0(0), fhEtaPhiBCPlus(0), fhEtaPhiBCMinus(0),
51fhEtaPhiBC0PileUpSPD(0),
52fhEtaPhiBCPlusPileUpSPD(0), fhEtaPhiBCMinusPileUpSPD(0),
53fhPtNPileUpSPDVtx(0), fhPtNPileUpTrkVtx(0),
54fhPtNPileUpSPDVtxTimeCut(0), fhPtNPileUpTrkVtxTimeCut(0),
55fhPtNPileUpSPDVtxTimeCut2(0), fhPtNPileUpTrkVtxTimeCut2(0)
56{
57 //default ctor
58
59 for(Int_t i = 0; i < 7; i++)
60 {
61 fhPtPileUp [i] = 0;
62 fhPtNeutralPileUp[i] = 0;
63
64 fhLambda0PileUp [i] = 0;
65 fhLambda0NeutralPileUp[i] = 0;
66
67 fhClusterEFracLongTimePileUp [i] = 0;
68
69 fhClusterCellTimePileUp [i] = 0;
70 fhClusterTimeDiffPileUp [i] = 0;
71 fhClusterTimeDiffNeutralPileUp[i] = 0;
72
73 }
74
75 for(Int_t i = 0; i < 4; i++)
76 {
77 fhClusterMultSPDPileUp[i] = 0;
78 fhClusterMultNoPileUp [i] = 0;
79 }
80
81 //Initialize parameters
82 InitParameters();
83
84}
85
86//___________________________________________
87TObjString * AliAnaClusterPileUp::GetAnalysisCuts()
88{
89 //Save parameters used for analysis
90 TString parList ; //this will be list of parameters used for this analysis.
91 const Int_t buffersize = 255;
92 char onePar[buffersize] ;
93
94 snprintf(onePar,buffersize,"--- AliAnaClusterPileUp ---\n") ;
95 parList+=onePar ;
96 snprintf(onePar,buffersize,"Calorimeter: %s\n",fCalorimeter.Data()) ;
97 parList+=onePar ;
98
99 //Get parameters set in base class.
100 //parList += GetBaseParametersList() ;
101
102 return new TObjString(parList) ;
103}
104
105//________________________________________________________________________
106TList * AliAnaClusterPileUp::GetCreateOutputObjects()
107{
108 // Create histograms to be saved in output file and
109 // store them in outputContainer
110 TList * outputContainer = new TList() ;
111 outputContainer->SetName("PhotonHistos") ;
112
113 Int_t nptbins = GetHistogramRanges()->GetHistoPtBins(); Float_t ptimecluster = GetHistogramRanges()->GetHistoPtMax(); Float_t ptmin = GetHistogramRanges()->GetHistoPtMin();
114 Int_t nphibins = GetHistogramRanges()->GetHistoPhiBins(); Float_t phimax = GetHistogramRanges()->GetHistoPhiMax(); Float_t phimin = GetHistogramRanges()->GetHistoPhiMin();
115 Int_t netabins = GetHistogramRanges()->GetHistoEtaBins(); Float_t etamax = GetHistogramRanges()->GetHistoEtaMax(); Float_t etamin = GetHistogramRanges()->GetHistoEtaMin();
116 Int_t ssbins = GetHistogramRanges()->GetHistoShowerShapeBins(); Float_t ssmax = GetHistogramRanges()->GetHistoShowerShapeMax(); Float_t ssmin = GetHistogramRanges()->GetHistoShowerShapeMin();
117 Int_t ntimebins= GetHistogramRanges()->GetHistoTimeBins(); Float_t timemax = GetHistogramRanges()->GetHistoTimeMax(); Float_t timemin = GetHistogramRanges()->GetHistoTimeMin();
118
119
120 fhTimePtNoCut = new TH2F ("hTimePt_NoCut","time of cluster vs pT of clusters, no event selection", nptbins,ptmin,ptimecluster, ntimebins,timemin,timemax);
121 fhTimePtNoCut->SetXTitle("#it{p}_{T} (GeV/#it{c})");
122 fhTimePtNoCut->SetYTitle("#it{time} (ns)");
123 outputContainer->Add(fhTimePtNoCut);
124
125 fhTimePtSPD = new TH2F ("hTimePt_SPD","time of cluster vs pT of clusters, SPD Pile-up events", nptbins,ptmin,ptimecluster, ntimebins,timemin,timemax);
126 fhTimePtSPD->SetXTitle("#it{p}_{T} (GeV/#it{c})");
127 fhTimePtSPD->SetYTitle("#it{time} (ns)");
128 outputContainer->Add(fhTimePtSPD);
129
130 fhTimeNPileUpVertSPD = new TH2F ("hTime_NPileUpVertSPD","time of cluster vs N pile-up SPD vertex", ntimebins,timemin,timemax,20,0,20);
131 fhTimeNPileUpVertSPD->SetYTitle("# vertex ");
132 fhTimeNPileUpVertSPD->SetXTitle("#it{time} (ns)");
133 outputContainer->Add(fhTimeNPileUpVertSPD);
134
135 fhTimeNPileUpVertTrack = new TH2F ("hTime_NPileUpVertTracks","time of cluster vs N pile-up Tracks vertex", ntimebins,timemin,timemax, 20,0,20 );
136 fhTimeNPileUpVertTrack->SetYTitle("# vertex ");
137 fhTimeNPileUpVertTrack->SetXTitle("#it{time} (ns)");
138 outputContainer->Add(fhTimeNPileUpVertTrack);
139
140 fhTimeNPileUpVertContributors = new TH2F ("hTime_NPileUpVertContributors","time of cluster vs N constributors to pile-up SPD vertex", ntimebins,timemin,timemax,50,0,50);
141 fhTimeNPileUpVertContributors->SetYTitle("# vertex ");
142 fhTimeNPileUpVertContributors->SetXTitle("#it{time} (ns)");
143 outputContainer->Add(fhTimeNPileUpVertContributors);
144
145 fhTimePileUpMainVertexZDistance = new TH2F ("hTime_PileUpMainVertexZDistance","time of cluster vs distance in Z pile-up SPD vertex - main SPD vertex",ntimebins,timemin,timemax,100,0,50);
146 fhTimePileUpMainVertexZDistance->SetYTitle("distance Z (cm) ");
147 fhTimePileUpMainVertexZDistance->SetXTitle("#it{time} (ns)");
148 outputContainer->Add(fhTimePileUpMainVertexZDistance);
149
150 fhTimePileUpMainVertexZDiamond = new TH2F ("hTime_PileUpMainVertexZDiamond","time of cluster vs distance in Z pile-up SPD vertex - z diamond",ntimebins,timemin,timemax,100,0,50);
151 fhTimePileUpMainVertexZDiamond->SetYTitle("diamond distance Z (cm) ");
152 fhTimePileUpMainVertexZDiamond->SetXTitle("#it{time} (ns)");
153 outputContainer->Add(fhTimePileUpMainVertexZDiamond);
154
155 fhEtaPhiBC0 = new TH2F ("hEtaPhiBC0","eta-phi for clusters tof corresponding to BC=0",netabins,etamin,etamax, nphibins,phimin,phimax);
156 fhEtaPhiBC0->SetXTitle("#eta ");
157 fhEtaPhiBC0->SetYTitle("#phi (rad)");
158 outputContainer->Add(fhEtaPhiBC0);
159
160 fhEtaPhiBCPlus = new TH2F ("hEtaPhiBCPlus","eta-phi for clusters tof corresponding to BC>0",netabins,etamin,etamax, nphibins,phimin,phimax);
161 fhEtaPhiBCPlus->SetXTitle("#eta ");
162 fhEtaPhiBCPlus->SetYTitle("#phi (rad)");
163 outputContainer->Add(fhEtaPhiBCPlus);
164
165 fhEtaPhiBCMinus = new TH2F ("hEtaPhiBCMinus","eta-phi for clusters tof corresponding to BC<0",netabins,etamin,etamax, nphibins,phimin,phimax);
166 fhEtaPhiBCMinus->SetXTitle("#eta ");
167 fhEtaPhiBCMinus->SetYTitle("#phi (rad)");
168 outputContainer->Add(fhEtaPhiBCMinus);
169
170 fhEtaPhiBC0PileUpSPD = new TH2F ("hEtaPhiBC0PileUpSPD","eta-phi for clusters tof corresponding to BC=0, SPD pile-up",netabins,etamin,etamax, nphibins,phimin,phimax);
171 fhEtaPhiBC0PileUpSPD->SetXTitle("#eta ");
172 fhEtaPhiBC0PileUpSPD->SetYTitle("#phi (rad)");
173 outputContainer->Add(fhEtaPhiBC0PileUpSPD);
174
175 fhEtaPhiBCPlusPileUpSPD = new TH2F ("hEtaPhiBCPlusPileUpSPD","eta-phi for clusters tof corresponding to BC>0, SPD pile-up",netabins,etamin,etamax, nphibins,phimin,phimax);
176 fhEtaPhiBCPlusPileUpSPD->SetXTitle("#eta ");
177 fhEtaPhiBCPlusPileUpSPD->SetYTitle("#phi (rad)");
178 outputContainer->Add(fhEtaPhiBCPlusPileUpSPD);
179
180 fhEtaPhiBCMinusPileUpSPD = new TH2F ("hEtaPhiBCMinusPileUpSPD","eta-phi for clusters tof corresponding to BC<0, SPD pile-up",netabins,etamin,etamax, nphibins,phimin,phimax);
181 fhEtaPhiBCMinusPileUpSPD->SetXTitle("#eta ");
182 fhEtaPhiBCMinusPileUpSPD->SetYTitle("#phi (rad)");
183 outputContainer->Add(fhEtaPhiBCMinusPileUpSPD);
184
185 TString title[] = {"no |t diff| cut","|t diff|<20 ns","|t diff|>20 ns","|t diff|>40 ns"};
186 TString name [] = {"TDiffNoCut","TDiffSmaller20ns","TDiffLarger20ns","TDiffLarger40ns"};
187 for(Int_t i = 0; i < 4; i++)
188 {
189 fhClusterMultSPDPileUp[i] = new TH2F(Form("fhClusterMultSPDPileUp_%s", name[i].Data()),
190 Form("Number of clusters per pile up event with #it{E} > 0.5 and %s respect cluster max vs cluster max E ",title[i].Data()),
191 nptbins,ptmin,ptimecluster,100,0,100);
192 fhClusterMultSPDPileUp[i]->SetYTitle("n clusters ");
193 fhClusterMultSPDPileUp[i]->SetXTitle("#it{E}_{cluster max} (GeV)");
194 outputContainer->Add(fhClusterMultSPDPileUp[i]) ;
195
196 fhClusterMultNoPileUp[i] = new TH2F(Form("fhClusterMultNoPileUp_%s", name[i].Data()),
197 Form("Number of clusters per non pile up event with #it{E} > 0.5 and %s respect cluster max vs cluster max E ",title[i].Data()),
198 nptbins,ptmin,ptimecluster,100,0,100);
199 fhClusterMultNoPileUp[i]->SetYTitle("n clusters ");
200 fhClusterMultNoPileUp[i]->SetXTitle("#it{E}_{cluster max} (GeV)");
201 outputContainer->Add(fhClusterMultNoPileUp[i]) ;
202 }
203
204 fhPtNPileUpSPDVtx = new TH2F ("hPt_NPileUpVertSPD","pT of cluster vs N pile-up SPD vertex",
205 nptbins,ptmin,ptimecluster,20,0,20);
206 fhPtNPileUpSPDVtx->SetYTitle("# vertex ");
207 fhPtNPileUpSPDVtx->SetXTitle("#it{p}_{T} (GeV/#it{c})");
208 outputContainer->Add(fhPtNPileUpSPDVtx);
209
210 fhPtNPileUpTrkVtx = new TH2F ("hPt_NPileUpVertTracks","pT of cluster vs N pile-up Tracks vertex",
211 nptbins,ptmin,ptimecluster, 20,0,20 );
212 fhPtNPileUpTrkVtx->SetYTitle("# vertex ");
213 fhPtNPileUpTrkVtx->SetXTitle("#it{p}_{T} (GeV/#it{c})");
214 outputContainer->Add(fhPtNPileUpTrkVtx);
215
216 fhPtNPileUpSPDVtxTimeCut = new TH2F ("hPt_NPileUpVertSPD_TimeCut","pT of cluster vs N pile-up SPD vertex, |tof| < 25 ns",
217 nptbins,ptmin,ptimecluster,20,0,20);
218 fhPtNPileUpSPDVtxTimeCut->SetYTitle("# vertex ");
219 fhPtNPileUpSPDVtxTimeCut->SetXTitle("#it{p}_{T} (GeV/#it{c})");
220 outputContainer->Add(fhPtNPileUpSPDVtxTimeCut);
221
222 fhPtNPileUpTrkVtxTimeCut = new TH2F ("hPt_NPileUpVertTracks_TimeCut","pT of cluster vs N pile-up Tracks vertex, |tof| < 25 ns",
223 nptbins,ptmin,ptimecluster, 20,0,20 );
224 fhPtNPileUpTrkVtxTimeCut->SetYTitle("# vertex ");
225 fhPtNPileUpTrkVtxTimeCut->SetXTitle("#it{p}_{T} (GeV/#it{c})");
226 outputContainer->Add(fhPtNPileUpTrkVtxTimeCut);
227
228 fhPtNPileUpSPDVtxTimeCut2 = new TH2F ("hPt_NPileUpVertSPD_TimeCut2","pT of cluster vs N pile-up SPD vertex, -25 < tof < 75 ns",
229 nptbins,ptmin,ptimecluster,20,0,20);
230 fhPtNPileUpSPDVtxTimeCut2->SetYTitle("# vertex ");
231 fhPtNPileUpSPDVtxTimeCut2->SetXTitle("#it{p}_{T} (GeV/#it{c})");
232 outputContainer->Add(fhPtNPileUpSPDVtxTimeCut2);
233
234 fhPtNPileUpTrkVtxTimeCut2 = new TH2F ("hPt_NPileUpVertTracks_TimeCut2","pT of cluster vs N pile-up Tracks vertex, -25 < tof < 75 ns",
235 nptbins,ptmin,ptimecluster, 20,0,20 );
236 fhPtNPileUpTrkVtxTimeCut2->SetYTitle("# vertex ");
237 fhPtNPileUpTrkVtxTimeCut2->SetXTitle("#it{p}_{T} (GeV/#it{c})");
238 outputContainer->Add(fhPtNPileUpTrkVtxTimeCut2);
239
240
241 TString pileUpName[] = {"SPD","EMCAL","SPDOrEMCAL","SPDAndEMCAL","SPDAndNotEMCAL","EMCALAndNotSPD","NotSPDAndNotEMCAL"} ;
242
243 for(Int_t i = 0 ; i < 7 ; i++)
244 {
245 fhPtPileUp[i] = new TH1F(Form("hPtPileUp%s",pileUpName[i].Data()),
246 Form("Cluster #it{p}_{T} distribution, %s Pile-Up event",pileUpName[i].Data()), nptbins,ptmin,ptimecluster);
247 fhPtPileUp[i]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
248 outputContainer->Add(fhPtPileUp[i]);
249
250 fhPtNeutralPileUp[i] = new TH1F(Form("hPtNeutralPileUp%s",pileUpName[i].Data()),
251 Form("Neutral clusters #it{p}_{T} distribution, %s Pile-Up event",pileUpName[i].Data()), nptbins,ptmin,ptimecluster);
252 fhPtNeutralPileUp[i]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
253 outputContainer->Add(fhPtNeutralPileUp[i]);
254
255 fhClusterEFracLongTimePileUp[i] = new TH2F(Form("hClusterEFracLongTimePileUp%s",pileUpName[i].Data()),
256 Form("Cluster E vs fraction of cluster energy from large T cells, %s Pile-Up event",pileUpName[i].Data()),
257 nptbins,ptmin,ptimecluster,200,0,1);
258 fhClusterEFracLongTimePileUp[i]->SetXTitle("#it{E} (GeV)");
259 fhClusterEFracLongTimePileUp[i]->SetYTitle("E(large time) / E");
260 outputContainer->Add(fhClusterEFracLongTimePileUp[i]);
261
262 fhClusterCellTimePileUp[i] = new TH2F(Form("hClusterCellTimePileUp%s",pileUpName[i].Data()),
263 Form("Cluster E vs cell time in cluster, %s Pile-Up event",pileUpName[i].Data()),
264 nptbins,ptmin,ptimecluster,ntimebins,timemin,timemax);
265 fhClusterCellTimePileUp[i]->SetXTitle("#it{E} (GeV)");
266 fhClusterCellTimePileUp[i]->SetYTitle("t_{cell} (ns)");
267 outputContainer->Add(fhClusterCellTimePileUp[i]);
268
269 fhClusterTimeDiffPileUp[i] = new TH2F(Form("hClusterTimeDiffPileUp%s",pileUpName[i].Data()),
270 Form("Cluster E vs t_{max}-t_{cell} in cluster, %s Pile-Up event",pileUpName[i].Data()),
271 nptbins,ptmin,ptimecluster,400,-200,200);
272 fhClusterTimeDiffPileUp[i]->SetXTitle("#it{E} (GeV)");
273 fhClusterTimeDiffPileUp[i]->SetYTitle("t_{max}-t_{cell} (ns)");
274 outputContainer->Add(fhClusterTimeDiffPileUp[i]);
275
276 fhClusterTimeDiffNeutralPileUp[i] = new TH2F(Form("hClusterTimeDiffNeutralPileUp%s",pileUpName[i].Data()),
277 Form("Neutral clusters E vs t_{max}-t_{cell} in cluster, %s Pile-Up event",pileUpName[i].Data()),
278 nptbins,ptmin,ptimecluster,400,-200,200);
279 fhClusterTimeDiffNeutralPileUp[i]->SetXTitle("#it{E} (GeV)");
280 fhClusterTimeDiffNeutralPileUp[i]->SetYTitle("t_{max}-t_{cell} (ns)");
281 outputContainer->Add(fhClusterTimeDiffNeutralPileUp[i]);
282
283 fhLambda0PileUp[i] = new TH2F(Form("hLambda0PileUp%s",pileUpName[i].Data()),
284 Form("Cluster E vs #lambda^{2}_{0} in cluster, %s Pile-Up event",pileUpName[i].Data()),
285 nptbins,ptmin,ptimecluster,ssbins,ssmin,ssmax);
286 fhLambda0PileUp[i]->SetXTitle("#it{E} (GeV)");
287 fhLambda0PileUp[i]->SetYTitle("#lambda^{2}_{0}");
288 outputContainer->Add(fhLambda0PileUp[i]);
289
290 fhLambda0NeutralPileUp[i] = new TH2F(Form("hLambda0NeutralPileUp%s",pileUpName[i].Data()),
291 Form("Neutral clusters E vs #lambda^{2}_{0}in cluster, %s Pile-Up event",pileUpName[i].Data()), nptbins,ptmin,ptimecluster,ssbins,ssmin,ssmax);
292 fhLambda0NeutralPileUp[i]->SetXTitle("#it{E} (GeV)");
293 fhLambda0NeutralPileUp[i]->SetYTitle("#lambda^{2}_{0}");
294 outputContainer->Add(fhLambda0NeutralPileUp[i]);
295
296 }
297
298 return outputContainer ;
299
300}
301
302//_______________________
303void AliAnaClusterPileUp::Init()
304{
305 //Init
306
307 //Do some checks
308 if(fCalorimeter == "PHOS" && !GetReader()->IsPHOSSwitchedOn())
309 AliFatal("You want to use PHOS in analysis but it is not read!! \n!!Check the configuration file!!");
310
311 if(fCalorimeter == "EMCAL" && !GetReader()->IsEMCALSwitchedOn())
312 AliFatal("You want to use EMCAL in analysis but it is not read!! \n!!Check the configuration file!!");
313
314 if(GetReader()->GetDataType() == AliCaloTrackReader::kMC)
315 AliFatal("You want to use MC data in analysis but this is not possible in pile-up!!");
316
317}
318
319//____________________________________________________________________________
320void AliAnaClusterPileUp::InitParameters()
321{
322
323 //Initialize the parameters of the analysis.
324 AddToHistogramsName("AnaClusterPileUp_");
325
326 fCalorimeter = "EMCAL" ;
327 fNCellsCut = 2;
328}
329
330//__________________________________________________________________
331void AliAnaClusterPileUp::MakeAnalysisFillHistograms()
332{
333 // Do cluster analysis
334 // Remember to open time cuts in reader
335
336 //Select the calorimeter
337 TObjArray * pl = 0x0;
338 AliVCaloCells* cells = 0;
339 if (fCalorimeter == "PHOS" )
340 {
341 pl = GetPHOSClusters();
342 cells = GetPHOSCells();
343 }
344 else if (fCalorimeter == "EMCAL")
345 {
346 pl = GetEMCALClusters();
347 cells = GetEMCALCells();
348 }
349
350 if(!pl)
351 {
352 Info("MakeAnalysisFillAOD","TObjArray with %s clusters is NULL!\n",fCalorimeter.Data());
353 return;
354 }
355
356 AliVEvent * event = GetReader()->GetInputEvent();
357 AliESDEvent* esdEv = dynamic_cast<AliESDEvent*> (event);
358 AliAODEvent* aodEv = dynamic_cast<AliAODEvent*> (event);
359
360 //-------------------
361 // N pile up vertices
362 Int_t nVtxSPD = -1;
363 Int_t nVtxTrk = -1;
364
365 if (esdEv)
366 {
367 nVtxSPD = esdEv->GetNumberOfPileupVerticesSPD();
368 nVtxTrk = esdEv->GetNumberOfPileupVerticesTracks();
369 }//ESD
370 else if (aodEv)
371 {
372 nVtxSPD = aodEv->GetNumberOfPileupVerticesSPD();
373 nVtxTrk = aodEv->GetNumberOfPileupVerticesTracks();
374 }//AOD
375
376 //-------------------
377 // Loop on clusters
378 Int_t nCaloClusters = pl->GetEntriesFast();
379
380 if(GetDebug() > 0) printf("AliAnaClusterPileUp::MakeAnalysisFillAOD() - input %s cluster entries %d\n", fCalorimeter.Data(), nCaloClusters);
381 //Init variables
382 TLorentzVector mom;
383 Int_t idMax = 0;
384 Float_t ptMax = 0;
385 Float_t tMax = 0;
386
387 for(Int_t icalo = 0; icalo < nCaloClusters; icalo++)
388 {
389 AliVCluster * calo = (AliVCluster*) (pl->At(icalo));
390 //printf("calo %d, %f\n",icalo,calo->E());
391
392 if(!calo) continue; // it should not happen, but just in case
393
394 calo->GetMomentum(mom,GetVertex(0)) ;
395
396 Float_t ecluster = mom.E();
397 Float_t ptcluster = mom.Pt();
398 Float_t l0cluster = calo->GetM02();
399 Float_t etacluster = mom.Eta();
400 Float_t phicluster = mom.Phi();
401 if(phicluster < 0) phicluster+=TMath::TwoPi();
402 Float_t tofcluster = calo->GetTOF()*1.e9;
403
404 Bool_t matched = IsTrackMatched(calo,GetReader()->GetInputEvent());
405
406 //.......................................
407 //If too small or big energy, skip it
408 if(ecluster < GetMinEnergy() || ecluster > GetMaxEnergy() ) continue ;
409
410 //.......................................
411 if(calo->GetNCells() <= fNCellsCut && GetReader()->GetDataType() != AliCaloTrackReader::kMC) continue;
412
413 //.......................................
414 //Check acceptance selection
415 if(IsFiducialCutOn())
416 {
417 Bool_t in = GetFiducialCut()->IsInFiducialCut(mom,fCalorimeter) ;
418 if(! in ) continue;
419 }
420
421 // Select highest pt cluster passing the cuts
422 if(ptcluster > ptMax && tofcluster < 30)
423 {
424 ptMax = ptcluster;
425 tMax = tofcluster;
426 idMax = icalo;
427 }
428
429 //-------------------------------------
430 // Cluster timing for different pile-up
431
432 fhTimePtNoCut->Fill(ptcluster,tofcluster);
433 if(GetReader()->IsPileUpFromSPD()) fhTimePtSPD->Fill(ptcluster,tofcluster);
434
435 //----------------------------------------
436 // correlate cluster and number of vertices
437
438 fhPtNPileUpSPDVtx->Fill(ptcluster,nVtxSPD);
439 fhPtNPileUpTrkVtx->Fill(ptcluster,nVtxTrk);
440
441 if(TMath::Abs(tofcluster) < 30)
442 {
443 fhPtNPileUpSPDVtxTimeCut->Fill(ptcluster,nVtxSPD);
444 fhPtNPileUpTrkVtxTimeCut->Fill(ptcluster,nVtxTrk);
445 }
446
447 if(tofcluster < 75 && tofcluster > -30)
448 {
449 fhPtNPileUpSPDVtxTimeCut2->Fill(ptcluster,nVtxSPD);
450 fhPtNPileUpTrkVtxTimeCut2->Fill(ptcluster,nVtxTrk);
451 }
452
453 // Loop on the vertices arrays, correlate with timing
454 // only for sufficiently large cluster energy
455 if(ecluster > 8)
456 {
457 fhTimeNPileUpVertSPD ->Fill(tofcluster,nVtxSPD);
458 fhTimeNPileUpVertTrack->Fill(tofcluster,nVtxTrk);
459
460 Int_t ncont = -1;
461 Float_t z1 = -1, z2 = -1;
462 Float_t diamZ = -1;
463 for(Int_t iVert=0; iVert<nVtxSPD;iVert++)
464 {
465 if (esdEv)
466 {
467 const AliESDVertex* pv=esdEv->GetPileupVertexSPD(iVert);
468 ncont=pv->GetNContributors();
469 z1 = esdEv->GetPrimaryVertexSPD()->GetZ();
470 z2 = pv->GetZ();
471 diamZ = esdEv->GetDiamondZ();
472 }//ESD
473 else if (aodEv)
474 {
475 AliAODVertex *pv=aodEv->GetVertex(iVert);
476 if(pv->GetType()!=AliAODVertex::kPileupSPD) continue;
477 ncont=pv->GetNContributors();
478 z1=aodEv->GetPrimaryVertexSPD()->GetZ();
479 z2=pv->GetZ();
480 diamZ = aodEv->GetDiamondZ();
481 }// AOD
482
483 Double_t distZ = TMath::Abs(z2-z1);
484 diamZ = TMath::Abs(z2-diamZ);
485
486 fhTimeNPileUpVertContributors ->Fill(tofcluster,ncont);
487 fhTimePileUpMainVertexZDistance->Fill(tofcluster,distZ);
488 fhTimePileUpMainVertexZDiamond ->Fill(tofcluster,diamZ);
489
490 }// vertex loop
491 }
492
493 //------------------------------------
494 // Eta-Phi cluster position depending on timing
495 // Continue only for BC0
496 if (tofcluster > 28)
497 {
498 fhEtaPhiBCPlus ->Fill(etacluster,phicluster);
499 if(GetReader()->IsPileUpFromSPD()) fhEtaPhiBCPlusPileUpSPD ->Fill(etacluster,phicluster);
500 continue;
501 }
502 else if (tofcluster <-28)
503 {
504 fhEtaPhiBCMinus->Fill(etacluster,phicluster);
505 if(GetReader()->IsPileUpFromSPD()) fhEtaPhiBCMinusPileUpSPD->Fill(etacluster,phicluster);
506 continue ;
507 }
508
509 //--------------------------------------
510 // Fill histograms for clusters in BC=0
511
512 fhEtaPhiBC0->Fill(etacluster,phicluster); if(GetReader()->IsPileUpFromSPD()) fhEtaPhiBC0PileUpSPD ->Fill(etacluster,phicluster);
513
514 if(GetReader()->IsPileUpFromSPD()) {fhPtPileUp[0]->Fill(ptcluster); fhLambda0PileUp[0]->Fill(ptcluster,l0cluster); }
515 if(GetReader()->IsPileUpFromEMCal()) {fhPtPileUp[1]->Fill(ptcluster); fhLambda0PileUp[1]->Fill(ptcluster,l0cluster); }
516 if(GetReader()->IsPileUpFromSPDOrEMCal()) {fhPtPileUp[2]->Fill(ptcluster); fhLambda0PileUp[2]->Fill(ptcluster,l0cluster); }
517 if(GetReader()->IsPileUpFromSPDAndEMCal()) {fhPtPileUp[3]->Fill(ptcluster); fhLambda0PileUp[3]->Fill(ptcluster,l0cluster); }
518 if(GetReader()->IsPileUpFromSPDAndNotEMCal()) {fhPtPileUp[4]->Fill(ptcluster); fhLambda0PileUp[4]->Fill(ptcluster,l0cluster); }
519 if(GetReader()->IsPileUpFromEMCalAndNotSPD()) {fhPtPileUp[5]->Fill(ptcluster); fhLambda0PileUp[5]->Fill(ptcluster,l0cluster); }
520 if(GetReader()->IsPileUpFromNotSPDAndNotEMCal()) {fhPtPileUp[6]->Fill(ptcluster); fhLambda0PileUp[6]->Fill(ptcluster,l0cluster); }
521
522 if(!matched)
523 {
524 if(GetReader()->IsPileUpFromSPD()) {fhPtNeutralPileUp[0]->Fill(ptcluster); fhLambda0NeutralPileUp[0]->Fill(ptcluster,l0cluster); }
525 if(GetReader()->IsPileUpFromEMCal()) {fhPtNeutralPileUp[1]->Fill(ptcluster); fhLambda0NeutralPileUp[1]->Fill(ptcluster,l0cluster); }
526 if(GetReader()->IsPileUpFromSPDOrEMCal()) {fhPtNeutralPileUp[2]->Fill(ptcluster); fhLambda0NeutralPileUp[2]->Fill(ptcluster,l0cluster); }
527 if(GetReader()->IsPileUpFromSPDAndEMCal()) {fhPtNeutralPileUp[3]->Fill(ptcluster); fhLambda0NeutralPileUp[3]->Fill(ptcluster,l0cluster); }
528 if(GetReader()->IsPileUpFromSPDAndNotEMCal()) {fhPtNeutralPileUp[4]->Fill(ptcluster); fhLambda0NeutralPileUp[4]->Fill(ptcluster,l0cluster); }
529 if(GetReader()->IsPileUpFromEMCalAndNotSPD()) {fhPtNeutralPileUp[5]->Fill(ptcluster); fhLambda0NeutralPileUp[5]->Fill(ptcluster,l0cluster); }
530 if(GetReader()->IsPileUpFromNotSPDAndNotEMCal()) {fhPtNeutralPileUp[6]->Fill(ptcluster); fhLambda0NeutralPileUp[6]->Fill(ptcluster,l0cluster); }
531 }
532
533 //----------------------------------------------------------------------------
534 // Loop on cells inside cluster, max cell must be over 100 MeV and time in BC=0
535 // Get the fraction of the cluster energy that carries the cell with highest
536 // energy and its absId
537
538 Float_t maxCellFraction = 0.;
539 Int_t absIdMax = GetCaloUtils()->GetMaxEnergyCell(cells, calo,maxCellFraction);
540
541 Float_t clusterLongTimePt = 0;
542 Float_t clusterOKTimePt = 0;
543
544 if(cells->GetCellAmplitude(absIdMax) < 0.1) continue ;
545
546 for (Int_t ipos = 0; ipos < calo->GetNCells(); ipos++)
547 {
548 Int_t absId = calo->GetCellsAbsId()[ipos];
549
550 if( absId == absIdMax ) continue ;
551
552 Double_t time = cells->GetCellTime(absId);
553 Float_t amp = cells->GetCellAmplitude(absId);
554 Int_t bc = GetReader()->GetInputEvent()->GetBunchCrossNumber();
555 GetCaloUtils()->GetEMCALRecoUtils()->AcceptCalibrateCell(absId,bc,amp,time,cells);
556 time*=1e9;
557
558 Float_t diff = (tofcluster-time);
559
560 if(GetReader()->IsInTimeWindow(time,amp)) clusterOKTimePt += amp;
561 else clusterLongTimePt += amp;
562
563 if( cells->GetCellAmplitude(absIdMax) < 0.1 ) continue ;
564
565 if(GetReader()->IsPileUpFromSPD())
566 {
567 fhClusterCellTimePileUp[0]->Fill(ptcluster, time);
568 fhClusterTimeDiffPileUp[0]->Fill(ptcluster, diff);
569 if(!matched) fhClusterTimeDiffNeutralPileUp[0]->Fill(ptcluster, diff);
570 }
571
572 if(GetReader()->IsPileUpFromEMCal())
573 {
574 fhClusterCellTimePileUp[1]->Fill(ptcluster, time);
575 fhClusterTimeDiffPileUp[1]->Fill(ptcluster, diff);
576 if(!matched) fhClusterTimeDiffNeutralPileUp[1]->Fill(ptcluster, diff);
577 }
578
579 if(GetReader()->IsPileUpFromSPDOrEMCal())
580 {
581 fhClusterCellTimePileUp[2]->Fill(ptcluster, time);
582 fhClusterTimeDiffPileUp[2]->Fill(ptcluster, diff);
583 if(!matched) fhClusterTimeDiffNeutralPileUp[2]->Fill(ptcluster, diff);
584 }
585
586 if(GetReader()->IsPileUpFromSPDAndEMCal())
587 {
588 fhClusterCellTimePileUp[3]->Fill(ptcluster, time);
589 fhClusterTimeDiffPileUp[3]->Fill(ptcluster, diff);
590 if(!matched) fhClusterTimeDiffNeutralPileUp[3]->Fill(ptcluster, diff);
591 }
592
593 if(GetReader()->IsPileUpFromSPDAndNotEMCal())
594 {
595 fhClusterCellTimePileUp[4]->Fill(ptcluster, time);
596 fhClusterTimeDiffPileUp[4]->Fill(ptcluster, diff);
597 if(!matched) fhClusterTimeDiffNeutralPileUp[4]->Fill(ptcluster, diff);
598 }
599
600 if(GetReader()->IsPileUpFromEMCalAndNotSPD())
601 {
602 fhClusterCellTimePileUp[5]->Fill(ptcluster, time);
603 fhClusterTimeDiffPileUp[5]->Fill(ptcluster, diff);
604 if(!matched) fhClusterTimeDiffNeutralPileUp[5]->Fill(ptcluster, diff);
605 }
606
607 if(GetReader()->IsPileUpFromNotSPDAndNotEMCal())
608 {
609 fhClusterCellTimePileUp[6]->Fill(ptcluster, time);
610 fhClusterTimeDiffPileUp[6]->Fill(ptcluster, diff);
611 if(!matched) fhClusterTimeDiffNeutralPileUp[6]->Fill(ptcluster, diff);
612 }
613 }//loop
614
615 Float_t frac = 0;
616 if(clusterLongTimePt+clusterOKTimePt > 0.001)
617 frac = clusterLongTimePt/(clusterLongTimePt+clusterOKTimePt);
618 //printf("E long %f, E OK %f, Fraction large time %f, E %f\n",clusterLongTimePt,clusterOKTimePt,frac,ptcluster);
619
620 if(GetReader()->IsPileUpFromSPD()) fhClusterEFracLongTimePileUp[0]->Fill(ptcluster,frac);
621 if(GetReader()->IsPileUpFromEMCal()) fhClusterEFracLongTimePileUp[1]->Fill(ptcluster,frac);
622 if(GetReader()->IsPileUpFromSPDOrEMCal()) fhClusterEFracLongTimePileUp[2]->Fill(ptcluster,frac);
623 if(GetReader()->IsPileUpFromSPDAndEMCal()) fhClusterEFracLongTimePileUp[3]->Fill(ptcluster,frac);
624 if(GetReader()->IsPileUpFromSPDAndNotEMCal()) fhClusterEFracLongTimePileUp[4]->Fill(ptcluster,frac);
625 if(GetReader()->IsPileUpFromEMCalAndNotSPD()) fhClusterEFracLongTimePileUp[5]->Fill(ptcluster,frac);
626 if(GetReader()->IsPileUpFromNotSPDAndNotEMCal()) fhClusterEFracLongTimePileUp[6]->Fill(ptcluster,frac);
627 }//loop
628
629 //----------------------------------------------------------------------------------------------
630 // Loop again on clusters to compare this max cluster t and the rest of the clusters, if E > 0.3
631 Int_t n20 = 0;
632 Int_t n40 = 0;
633 Int_t n = 0;
634 Int_t nOK = 0;
635
636 for(Int_t icalo = 0; icalo < nCaloClusters; icalo++)
637 {
638 AliVCluster * calo = (AliVCluster*) (pl->At(icalo));
639
640 if(!calo || calo->E() < 0.3 || icalo == idMax) continue;
641
642 Float_t tdiff = TMath::Abs(tMax-calo->GetTOF()*1e9);
643 n++;
644 if(tdiff < 25) nOK++;
645 else
646 {
647 n20++;
648 if(tdiff > 40 ) n40++;
649 }
650 }
651
652 // Check pile-up and fill histograms depending on the different cluster multiplicities
653 if(GetReader()->IsPileUpFromSPD())
654 {
655 fhClusterMultSPDPileUp[0]->Fill(ptMax,n );
656 fhClusterMultSPDPileUp[1]->Fill(ptMax,nOK);
657 fhClusterMultSPDPileUp[2]->Fill(ptMax,n20);
658 fhClusterMultSPDPileUp[3]->Fill(ptMax,n40);
659 }
660 else
661 {
662 fhClusterMultNoPileUp[0]->Fill(ptMax,n );
663 fhClusterMultNoPileUp[1]->Fill(ptMax,nOK);
664 fhClusterMultNoPileUp[2]->Fill(ptMax,n20);
665 fhClusterMultNoPileUp[3]->Fill(ptMax,n40);
666 }
667
668
669 if(GetDebug() > 1) printf("AliAnaClusterPileUp::MakeAnalysisFillHistograms() End fill histograms\n");
670
671}
672
673
674
675//__________________________________________________________________
676void AliAnaClusterPileUp::Print(const Option_t * opt) const
677{
678 //Print some relevant parameters set for the analysis
679
680 if(! opt)
681 return;
682
683 printf("**** Print %s %s ****\n", GetName(), GetTitle() ) ;
684 AliAnaCaloTrackCorrBaseClass::Print(" ");
685
686 printf("Calorimeter = %s\n", fCalorimeter.Data()) ;
687 printf(" \n") ;
688
689}