1 /**************************************************************************
2 * Copyright(c) 1998-2008, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is 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 **************************************************************************/
16 //*****************************************************
17 // Class AliCentralitySelectionTask
18 // Class to analyze determine centrality
19 // author: Alberica Toia
20 //*****************************************************
22 #include "AliCentralitySelectionTask.h"
31 #include <TObjString.h>
35 #include <TDirectory.h>
39 #include "AliAnalysisManager.h"
40 #include "AliHeader.h"
41 #include "AliVEvent.h"
43 #include "AliESDEvent.h"
44 #include "AliESDHeader.h"
45 #include "AliESDInputHandler.h"
46 #include "AliESDZDC.h"
47 #include "AliESDFMD.h"
48 #include "AliESDVZERO.h"
49 #include "AliESDtrackCuts.h"
50 #include "AliESDVertex.h"
51 #include "AliCentrality.h"
52 #include "AliOADBCentrality.h"
53 #include "AliOADBContainer.h"
54 #include "AliMultiplicity.h"
55 #include "AliAODHandler.h"
56 #include "AliAODHeader.h"
57 #include "AliAODEvent.h"
58 #include "AliAODVertex.h"
59 #include "AliAODVZERO.h"
60 #include "AliAODTracklets.h"
61 #include "AliAODMCHeader.h"
62 #include "AliMCEventHandler.h"
63 #include "AliMCEvent.h"
64 #include "AliAODMCParticle.h"
65 #include "AliMCParticle.h"
67 #include "AliAnalysisTaskSE.h"
68 #include "AliGenEventHeader.h"
69 #include "AliGenHijingEventHeader.h"
70 #include "AliPhysicsSelectionTask.h"
71 #include "AliPhysicsSelection.h"
72 #include "AliBackgroundSelection.h"
73 #include "AliESDUtils.h"
75 ClassImp(AliCentralitySelectionTask)
78 //________________________________________________________________________
79 AliCentralitySelectionTask::AliCentralitySelectionTask():
81 fAnalysisInput("ESD"),
91 fV0MSPDOutlierPar0(0),
92 fV0MSPDOutlierPar1(0),
93 fV0MTPCOutlierPar0(0),
94 fV0MTPCOutlierPar1(0),
95 fV0MSPDSigmaOutlierPar0(0),
96 fV0MSPDSigmaOutlierPar1(0),
97 fV0MSPDSigmaOutlierPar2(0),
98 fV0MTPCSigmaOutlierPar0(0),
99 fV0MTPCSigmaOutlierPar1(0),
100 fV0MTPCSigmaOutlierPar2(0),
101 fV0MZDCOutlierPar0(0),
102 fV0MZDCOutlierPar1(0),
103 fV0MZDCEcalOutlierPar0(0),
104 fV0MZDCEcalOutlierPar1(0),
134 fHOutCentV0MvsFMD(0),
135 fHOutCentTKLvsV0M(0),
136 fHOutCentZEMvsZDC(0),
137 fHOutCentV0MvsCentCL1(0),
138 fHOutCentV0MvsCentTRK(0),
139 fHOutCentTRKvsCentCL1(0),
140 fHOutCentV0MvsCentZDC(0),
148 fHOutMultV0MvsZDN(0),
149 fHOutMultZEMvsZDN(0),
150 fHOutMultV0MvsZDC(0),
151 fHOutMultZEMvsZDC(0),
152 fHOutMultZEMvsZDCw(0),
153 fHOutMultV0MvsCL1(0),
154 fHOutMultV0MvsTRK(0),
155 fHOutMultTRKvsCL1(0),
156 fHOutMultV0MvsV0O(0),
157 fHOutMultV0OvsCL1(0),
158 fHOutMultV0OvsTRK(0),
159 fHOutCentV0Mqual1(0),
160 fHOutCentTRKqual1(0),
161 fHOutCentCL1qual1(0),
162 fHOutMultV0MvsCL1qual1(0),
163 fHOutMultV0MvsTRKqual1(0),
164 fHOutMultTRKvsCL1qual1(0),
165 fHOutCentV0Mqual2(0),
166 fHOutCentTRKqual2(0),
167 fHOutCentCL1qual2(0),
168 fHOutMultV0MvsCL1qual2(0),
169 fHOutMultV0MvsTRKqual2(0),
170 fHOutMultTRKvsCL1qual2(0),
174 // Default constructor
175 AliInfo("Centrality Selection enabled.");
179 fBranchNames="ESD:AliESDRun.,AliESDHeader.,AliESDZDC.,AliESDFMD.,AliESDVZERO."
180 ",SPDVertex.,TPCVertex.,PrimaryVertex.,AliMultiplicity.,Tracks ";
183 //________________________________________________________________________
184 AliCentralitySelectionTask::AliCentralitySelectionTask(const char *name):
185 AliAnalysisTaskSE(name),
186 fAnalysisInput("ESD"),
195 fV0MScaleFactorMC(0),
196 fV0MSPDOutlierPar0(0),
197 fV0MSPDOutlierPar1(0),
198 fV0MTPCOutlierPar0(0),
199 fV0MTPCOutlierPar1(0),
200 fV0MSPDSigmaOutlierPar0(0),
201 fV0MSPDSigmaOutlierPar1(0),
202 fV0MSPDSigmaOutlierPar2(0),
203 fV0MTPCSigmaOutlierPar0(0),
204 fV0MTPCSigmaOutlierPar1(0),
205 fV0MTPCSigmaOutlierPar2(0),
206 fV0MZDCOutlierPar0(0),
207 fV0MZDCOutlierPar1(0),
208 fV0MZDCEcalOutlierPar0(0),
209 fV0MZDCEcalOutlierPar1(0),
239 fHOutCentV0MvsFMD(0),
240 fHOutCentTKLvsV0M(0),
241 fHOutCentZEMvsZDC(0),
242 fHOutCentV0MvsCentCL1(0),
243 fHOutCentV0MvsCentTRK(0),
244 fHOutCentTRKvsCentCL1(0),
245 fHOutCentV0MvsCentZDC(0),
253 fHOutMultV0MvsZDN(0),
254 fHOutMultZEMvsZDN(0),
255 fHOutMultV0MvsZDC(0),
256 fHOutMultZEMvsZDC(0),
257 fHOutMultZEMvsZDCw(0),
258 fHOutMultV0MvsCL1(0),
259 fHOutMultV0MvsTRK(0),
260 fHOutMultTRKvsCL1(0),
261 fHOutMultV0MvsV0O(0),
262 fHOutMultV0OvsCL1(0),
263 fHOutMultV0OvsTRK(0),
264 fHOutCentV0Mqual1(0),
265 fHOutCentTRKqual1(0),
266 fHOutCentCL1qual1(0),
267 fHOutMultV0MvsCL1qual1(0),
268 fHOutMultV0MvsTRKqual1(0),
269 fHOutMultTRKvsCL1qual1(0),
270 fHOutCentV0Mqual2(0),
271 fHOutCentTRKqual2(0),
272 fHOutCentCL1qual2(0),
273 fHOutMultV0MvsCL1qual2(0),
274 fHOutMultV0MvsTRKqual2(0),
275 fHOutMultTRKvsCL1qual2(0),
279 // Default constructor
280 AliInfo("Centrality Selection enabled.");
281 DefineOutput(1, TList::Class());
284 fBranchNames="ESD:AliESDRun.,AliESDHeader.,AliESDZDC.,AliESDFMD.,AliESDVZERO."
285 ",SPDVertex.,TPCVertex.,PrimaryVertex.,AliMultiplicity.,Tracks ";
288 //________________________________________________________________________
289 AliCentralitySelectionTask& AliCentralitySelectionTask::operator=(const AliCentralitySelectionTask& c)
291 // Assignment operator
293 AliAnalysisTaskSE::operator=(c);
298 //________________________________________________________________________
299 AliCentralitySelectionTask::AliCentralitySelectionTask(const AliCentralitySelectionTask& ana):
300 AliAnalysisTaskSE(ana),
301 fAnalysisInput(ana.fDebug),
302 fIsMCInput(ana.fIsMCInput),
304 fCurrentRun(ana.fCurrentRun),
305 fUseScaling(ana.fUseScaling),
306 fUseCleaning(ana.fUseCleaning),
307 fV0MScaleFactor(ana.fV0MScaleFactor),
308 fSPDScaleFactor(ana.fSPDScaleFactor),
309 fTPCScaleFactor(ana.fTPCScaleFactor),
310 fV0MScaleFactorMC(ana.fV0MScaleFactorMC),
311 fV0MSPDOutlierPar0(ana.fV0MSPDOutlierPar0),
312 fV0MSPDOutlierPar1(ana.fV0MSPDOutlierPar1),
313 fV0MTPCOutlierPar0(ana.fV0MTPCOutlierPar0),
314 fV0MTPCOutlierPar1(ana.fV0MTPCOutlierPar1),
315 fV0MSPDSigmaOutlierPar0(ana.fV0MSPDSigmaOutlierPar0),
316 fV0MSPDSigmaOutlierPar1(ana.fV0MSPDSigmaOutlierPar1),
317 fV0MSPDSigmaOutlierPar2(ana.fV0MSPDSigmaOutlierPar2),
318 fV0MTPCSigmaOutlierPar0(ana.fV0MTPCSigmaOutlierPar0),
319 fV0MTPCSigmaOutlierPar1(ana.fV0MTPCSigmaOutlierPar1),
320 fV0MTPCSigmaOutlierPar2(ana.fV0MTPCSigmaOutlierPar2),
321 fV0MZDCOutlierPar0(ana.fV0MZDCOutlierPar0),
322 fV0MZDCOutlierPar1(ana.fV0MZDCOutlierPar1),
323 fV0MZDCEcalOutlierPar0(ana.fV0MZDCEcalOutlierPar0),
324 fV0MZDCEcalOutlierPar1(ana.fV0MZDCEcalOutlierPar1),
325 fTrackCuts(ana.fTrackCuts),
327 fOutliersCut(ana.fOutliersCut),
328 fQuality(ana.fQuality),
329 fCentV0M(ana.fCentV0M),
330 fCentFMD(ana.fCentFMD),
331 fCentTRK(ana.fCentTRK),
332 fCentTKL(ana.fCentTKL),
333 fCentCL0(ana.fCentCL0),
334 fCentCL1(ana.fCentCL1),
335 fCentV0MvsFMD(ana.fCentV0MvsFMD),
336 fCentTKLvsV0M(ana.fCentTKLvsV0M),
337 fCentZEMvsZDC(ana.fCentZEMvsZDC),
338 fHtempV0M(ana.fHtempV0M),
339 fHtempFMD(ana.fHtempFMD),
340 fHtempTRK(ana.fHtempTRK),
341 fHtempTKL(ana.fHtempTKL),
342 fHtempCL0(ana.fHtempCL0),
343 fHtempCL1(ana.fHtempCL1),
344 fHtempV0MvsFMD(ana.fHtempV0MvsFMD),
345 fHtempTKLvsV0M(ana.fHtempTKLvsV0M),
346 fHtempZEMvsZDC(ana.fHtempZEMvsZDC),
347 fOutputList(ana.fOutputList),
348 fHOutCentV0M (ana.fHOutCentV0M ),
349 fHOutCentFMD (ana.fHOutCentFMD ),
350 fHOutCentTRK (ana.fHOutCentTRK ),
351 fHOutCentTKL (ana.fHOutCentTKL ),
352 fHOutCentCL0 (ana.fHOutCentCL0 ),
353 fHOutCentCL1 (ana.fHOutCentCL1 ),
354 fHOutCentV0MvsFMD(ana.fHOutCentV0MvsFMD),
355 fHOutCentTKLvsV0M(ana.fHOutCentTKLvsV0M),
356 fHOutCentZEMvsZDC(ana.fHOutCentZEMvsZDC),
357 fHOutCentV0MvsCentCL1(ana.fHOutCentV0MvsCentCL1),
358 fHOutCentV0MvsCentTRK(ana.fHOutCentV0MvsCentTRK),
359 fHOutCentTRKvsCentCL1(ana.fHOutCentTRKvsCentCL1),
360 fHOutCentV0MvsCentZDC(ana.fHOutCentV0MvsCentZDC),
361 fHOutMultV0M(ana.fHOutMultV0M),
362 fHOutMultV0O(ana.fHOutMultV0O),
363 fHOutMultFMD(ana.fHOutMultFMD),
364 fHOutMultTRK(ana.fHOutMultTRK),
365 fHOutMultTKL(ana.fHOutMultTKL),
366 fHOutMultCL0(ana.fHOutMultCL0),
367 fHOutMultCL1(ana.fHOutMultCL1),
368 fHOutMultV0MvsZDN(ana.fHOutMultV0MvsZDN),
369 fHOutMultZEMvsZDN(ana.fHOutMultZEMvsZDN),
370 fHOutMultV0MvsZDC(ana.fHOutMultV0MvsZDC),
371 fHOutMultZEMvsZDC(ana.fHOutMultZEMvsZDC),
372 fHOutMultZEMvsZDCw(ana.fHOutMultZEMvsZDCw),
373 fHOutMultV0MvsCL1(ana.fHOutMultV0MvsCL1),
374 fHOutMultV0MvsTRK(ana.fHOutMultV0MvsTRK),
375 fHOutMultTRKvsCL1(ana.fHOutMultTRKvsCL1),
376 fHOutMultV0MvsV0O(ana.fHOutMultV0MvsV0O),
377 fHOutMultV0OvsCL1(ana.fHOutMultV0OvsCL1),
378 fHOutMultV0OvsTRK(ana.fHOutMultV0OvsTRK),
379 fHOutCentV0Mqual1(ana.fHOutCentV0Mqual1),
380 fHOutCentTRKqual1(ana.fHOutCentTRKqual1),
381 fHOutCentCL1qual1(ana.fHOutCentCL1qual1),
382 fHOutMultV0MvsCL1qual1(ana.fHOutMultV0MvsCL1qual1),
383 fHOutMultV0MvsTRKqual1(ana.fHOutMultV0MvsTRKqual1),
384 fHOutMultTRKvsCL1qual1(ana.fHOutMultTRKvsCL1qual1),
385 fHOutCentV0Mqual2(ana.fHOutCentV0Mqual2),
386 fHOutCentTRKqual2(ana.fHOutCentTRKqual2),
387 fHOutCentCL1qual2(ana.fHOutCentCL1qual2),
388 fHOutMultV0MvsCL1qual2(ana.fHOutMultV0MvsCL1qual2),
389 fHOutMultV0MvsTRKqual2(ana.fHOutMultV0MvsTRKqual2),
390 fHOutMultTRKvsCL1qual2(ana.fHOutMultTRKvsCL1qual2),
391 fHOutQuality(ana.fHOutQuality),
392 fHOutVertex(ana.fHOutVertex)
398 //________________________________________________________________________
399 AliCentralitySelectionTask::~AliCentralitySelectionTask()
402 if (fOutputList && !AliAnalysisManager::GetAnalysisManager()->IsProofMode()) delete fOutputList;
403 if (fTrackCuts) delete fTrackCuts;
406 //________________________________________________________________________
407 void AliCentralitySelectionTask::UserCreateOutputObjects()
409 // Create the output containers
410 if(fDebug>1) printf("AnalysisCentralitySelectionTask::UserCreateOutputObjects() \n");
411 AliLog::SetClassDebugLevel("AliCentralitySelectionTask", AliLog::kInfo);
413 fOutputList = new TList();
414 fOutputList->SetOwner();
415 fHOutCentV0M = new TH1F("fHOutCentV0M","fHOutCentV0M; Centrality V0",505,0,101);
416 fHOutCentFMD = new TH1F("fHOutCentFMD","fHOutCentFMD; Centrality FMD",505,0,101);
417 fHOutCentTRK = new TH1F("fHOutCentTRK","fHOutCentTRK; Centrality TPC",505,0,101);
418 fHOutCentTKL = new TH1F("fHOutCentTKL","fHOutCentTKL; Centrality tracklets",505,0,101);
419 fHOutCentCL0 = new TH1F("fHOutCentCL0","fHOutCentCL0; Centrality SPD inner",505,0,101);
420 fHOutCentCL1 = new TH1F("fHOutCentCL1","fHOutCentCL1; Centrality SPD outer",505,0,101);
421 fHOutCentV0MvsFMD= new TH1F("fHOutCentV0MvsFMD","fHOutCentV0MvsFMD; Centrality V0 vs FMD",505,0,101);
422 fHOutCentTKLvsV0M= new TH1F("fHOutCentTKLvsV0M","fHOutCentTKLvsV0M; Centrality tracklets vs V0",505,0,101);
423 fHOutCentZEMvsZDC= new TH1F("fHOutCentZEMvsZDC","fHOutCentZEMvsZDC; Centrality ZEM vs ZDC",505,0,101);
424 fHOutCentV0MvsCentCL1= new TH2F("fHOutCentV0MvsCentCL1","fHOutCentV0MvsCentCL1; Cent V0; Cent SPD",505,0,101,505,0,101);
425 fHOutCentV0MvsCentTRK= new TH2F("fHOutCentV0MvsCentTRK","fHOutCentV0MvsCentTRK; Cent V0; Cent TPC",505,0,101,505,0,101);
426 fHOutCentTRKvsCentCL1= new TH2F("fHOutCentTRKvsCentCL1","fHOutCentTRKvsCentCL1; Cent TPC; Cent SPD",505,0,101,505,0,101);
427 fHOutCentV0MvsCentZDC= new TH2F("fHOutCentV0MvsCentZDC","fHOutCentV0MvsCentZDC; Cent V0; Cent ZDC",505,0,101,505,0,101);
428 fHOutMultV0M = new TH1F("fHOutMultV0M","fHOutMultV0M; Multiplicity V0",25000,0,30000);
429 fHOutMultV0O = new TH1F("fHOutMultV0O","fHOutMultV0O; Multiplicity V0",40000,0,40000);
430 fHOutMultFMD = new TH1F("fHOutMultFMD","fHOutMultFMD; Multiplicity FMD",24000,0,24000);
431 fHOutMultTRK = new TH1F("fHOutMultTRK","fHOutMultTRK; Multiplicity TPC",4000,0,4000);
432 fHOutMultTKL = new TH1F("fHOutMultTKL","fHOutMultTKL; Multiplicity tracklets",5000,0,5000);
433 fHOutMultCL0 = new TH1F("fHOutMultCL0","fHOutMultCL0; Multiplicity SPD inner",7000,0,7000);
434 fHOutMultCL1 = new TH1F("fHOutMultCL1","fHOutMultCL1; Multiplicity SPD outer",7000,0,7000);
435 fHOutMultV0MvsZDN = new TH2F("fHOutMultV0MvsZDN","fHOutMultV0MvsZDN; Multiplicity V0; Energy ZDC-N",500,0,30000,500,0,180000);
436 fHOutMultZEMvsZDN = new TH2F("fHOutMultZEMvsZDN","fHOutMultZEMvsZDN; Energy ZEM; Energy ZDC-N",500,0,2500,500,0,180000);
437 fHOutMultV0MvsZDC = new TH2F("fHOutMultV0MvsZDC","fHOutMultV0MvsZDC; Multiplicity V0; Energy ZDC",500,0,30000,500,0,200000);
438 fHOutMultZEMvsZDC = new TH2F("fHOutMultZEMvsZDC","fHOutMultZEMvsZDC; Energy ZEM; Energy ZDC",500,0,2500,500,0,200000);
439 fHOutMultZEMvsZDCw = new TH2F("fHOutMultZEMvsZDCw","fHOutMultZEMvsZDCw; Energy ZEM; Energy ZDC (weigthed with V0 percentile)",500,0,2500,500,0,200000);
440 fHOutMultV0MvsCL1 = new TH2F("fHOutMultV0MvsCL1","fHOutMultV0MvsCL1; Multiplicity V0; Multiplicity SPD outer",2500,0,30000,700,0,7000);
441 fHOutMultV0MvsTRK = new TH2F("fHOutMultV0MvsTRK","fHOutMultV0MvsTRK; Multiplicity V0; Multiplicity TPC",2500,0,30000,400,0,4000);
442 fHOutMultTRKvsCL1 = new TH2F("fHOutMultTRKvsCL1","fHOutMultTRKvsCL1; Multiplicity TPC; Multiplicity SPD outer",400,0,4000,700,0,7000);
443 fHOutMultV0MvsV0O = new TH2F("fHOutMultV0MvsV0O","fHOutMultV0MvsV0O; Multiplicity V0; Multiplicity V0 Online",500,0,30000,500,0,40000);
444 fHOutMultV0OvsCL1 = new TH2F("fHOutMultV0OvsCL1","fHOutMultV0OvsCL1; Multiplicity V0; Multiplicity SPD outer",500,0,40000,700,0,7000);
445 fHOutMultV0OvsTRK = new TH2F("fHOutMultV0OvsTRK","fHOutMultV0OvsTRK; Multiplicity V0; Multiplicity TPC",500,0,40000,400,0,4000);
446 fHOutMultV0MvsV0O = new TH2F("fHOutMultV0MvsV0O","fHOutMultV0MvsV0O; Multiplicity V0; Multiplicity V0 Online",500,0,30000,500,0,30000);
447 fHOutMultV0OvsCL1 = new TH2F("fHOutMultV0OvsCL1","fHOutMultV0OvsCL1; Multiplicity V0; Multiplicity SPD outer",2500,0,30000,700,0,7000);
448 fHOutMultV0OvsTRK = new TH2F("fHOutMultV0OvsTRK","fHOutMultV0OvsTRK; Multiplicity V0; Multiplicity TPC",2500,0,30000,400,0,4000);
450 fHOutCentV0Mqual1 = new TH1F("fHOutCentV0M_qual1","fHOutCentV0M_qual1; Centrality V0",505,0,101);
451 fHOutCentTRKqual1 = new TH1F("fHOutCentTRK_qual1","fHOutCentTRK_qual1; Centrality TPC",505,0,101);
452 fHOutCentCL1qual1 = new TH1F("fHOutCentCL1_qual1","fHOutCentCL1_qual1; Centrality SPD outer",505,0,101);
453 fHOutMultV0MvsCL1qual1 = new TH2F("fHOutMultV0MvsCL1_qual1","fHOutMultV0MvsCL1_qual1; Multiplicity V0; Multiplicity SPD outer",2500,0,25000,700,0,7000);
454 fHOutMultV0MvsTRKqual1 = new TH2F("fHOutMultV0MvsTRK_qual1","fHOutMultV0MvsTRK_qual1; Multiplicity V0; Multiplicity TPC",2500,0,25000,400,0,4000);
455 fHOutMultTRKvsCL1qual1 = new TH2F("fHOutMultTRKvsCL1_qual1","fHOutMultTRKvsCL1_qual1; Multiplicity TPC; Multiplicity SPD outer",400,0,4000,700,0,7000);
457 fHOutCentV0Mqual2 = new TH1F("fHOutCentV0M_qual2","fHOutCentV0M_qual2; Centrality V0",505,0,101);
458 fHOutCentTRKqual2 = new TH1F("fHOutCentTRK_qual2","fHOutCentTRK_qual2; Centrality TPC",505,0,101);
459 fHOutCentCL1qual2 = new TH1F("fHOutCentCL1_qual2","fHOutCentCL1_qual2; Centrality SPD outer",505,0,101);
460 fHOutMultV0MvsCL1qual2 = new TH2F("fHOutMultV0MvsCL1_qual2","fHOutMultV0MvsCL1_qual2; Multiplicity V0; Multiplicity SPD outer",2500,0,25000,700,0,7000);
461 fHOutMultV0MvsTRKqual2 = new TH2F("fHOutMultV0MvsTRK_qual2","fHOutMultV0MvsTRK_qual2; Multiplicity V0; Multiplicity TPC",2500,0,25000,400,0,4000);
462 fHOutMultTRKvsCL1qual2 = new TH2F("fHOutMultTRKvsCL1_qual2","fHOutMultTRKvsCL1_qual2; Multiplicity TPC; Multiplicity SPD outer",400,0,4000,700,0,7000);
464 fHOutQuality = new TH1F("fHOutQuality", "fHOutQuality", 100,-0.5,99.5);
465 fHOutVertex = new TH1F("fHOutVertex", "fHOutVertex", 100,-20,20);
467 fOutputList->Add( fHOutCentV0M );
468 fOutputList->Add( fHOutCentFMD );
469 fOutputList->Add( fHOutCentTRK );
470 fOutputList->Add( fHOutCentTKL );
471 fOutputList->Add( fHOutCentCL0 );
472 fOutputList->Add( fHOutCentCL1 );
473 fOutputList->Add( fHOutCentV0MvsFMD);
474 fOutputList->Add( fHOutCentTKLvsV0M);
475 fOutputList->Add( fHOutCentZEMvsZDC);
476 fOutputList->Add( fHOutCentV0MvsCentCL1);
477 fOutputList->Add( fHOutCentV0MvsCentTRK);
478 fOutputList->Add( fHOutCentTRKvsCentCL1);
479 fOutputList->Add( fHOutCentV0MvsCentZDC);
480 fOutputList->Add( fHOutMultV0M);
481 fOutputList->Add( fHOutMultV0O);
482 fOutputList->Add( fHOutMultFMD);
483 fOutputList->Add( fHOutMultTRK);
484 fOutputList->Add( fHOutMultTKL);
485 fOutputList->Add( fHOutMultCL0);
486 fOutputList->Add( fHOutMultCL1);
487 fOutputList->Add( fHOutMultV0MvsZDN);
488 fOutputList->Add( fHOutMultZEMvsZDN);
489 fOutputList->Add( fHOutMultV0MvsZDC);
490 fOutputList->Add( fHOutMultZEMvsZDC);
491 fOutputList->Add( fHOutMultZEMvsZDCw);
492 fOutputList->Add( fHOutMultV0MvsCL1);
493 fOutputList->Add( fHOutMultV0MvsTRK);
494 fOutputList->Add( fHOutMultTRKvsCL1);
495 fOutputList->Add( fHOutMultV0MvsV0O);
496 fOutputList->Add( fHOutMultV0OvsCL1);
497 fOutputList->Add( fHOutMultV0OvsTRK);
498 fOutputList->Add( fHOutCentV0Mqual1 );
499 fOutputList->Add( fHOutCentTRKqual1 );
500 fOutputList->Add( fHOutCentCL1qual1 );
501 fOutputList->Add( fHOutMultV0MvsCL1qual1);
502 fOutputList->Add( fHOutMultV0MvsTRKqual1);
503 fOutputList->Add( fHOutMultTRKvsCL1qual1);
504 fOutputList->Add( fHOutCentV0Mqual2 );
505 fOutputList->Add( fHOutCentTRKqual2 );
506 fOutputList->Add( fHOutCentCL1qual2 );
507 fOutputList->Add( fHOutMultV0MvsCL1qual2);
508 fOutputList->Add( fHOutMultV0MvsTRKqual2);
509 fOutputList->Add( fHOutMultTRKvsCL1qual2);
510 fOutputList->Add( fHOutQuality );
511 fOutputList->Add( fHOutVertex );
514 fTrackCuts = AliESDtrackCuts::GetStandardTPCOnlyTrackCuts();
516 PostData(1, fOutputList);
518 if (fPass==0) AliFatal("Which pass are you analyzing? You should set it via taskCentrality->SetPass(N)");
521 //________________________________________________________________________
522 void AliCentralitySelectionTask::UserExec(Option_t */*option*/)
524 // Execute analysis for current event:
525 if(fDebug>1) printf(" **** AliCentralitySelectionTask::UserExec() \n");
527 Float_t zncEnergy = 0.; // ZNC Energy
528 Float_t zpcEnergy = 0.; // ZPC Energy
529 Float_t znaEnergy = 0.; // ZNA Energy
530 Float_t zpaEnergy = 0.; // ZPA Energy
531 Float_t zem1Energy = 0.; // ZEM1 Energy
532 Float_t zem2Energy = 0.; // ZEM2 Energy
533 Bool_t zdcEnergyCal = kFALSE; // if zdc is calibrated (in pass2)
535 Int_t nTracks = 0; // no. tracks
536 Int_t nTracklets = 0; // no. tracklets
537 Int_t nClusters[6] = {0}; // no. clusters on 6 ITS layers
538 Int_t nChips[2]; // no. chips on 2 SPD layers
539 Float_t spdCorr =0; // corrected spd2 multiplicity
541 Float_t multV0A = 0; // multiplicity from V0 reco side A
542 Float_t multV0C = 0; // multiplicity from V0 reco side C
543 Short_t multV0AOnline = 0; // multiplicity from V0 reco side A
544 Short_t multV0COnline = 0; // multiplicity from V0 reco side C
545 Float_t v0Corr = 0; // corrected V0 multiplicity (used for MC)
547 Float_t multFMDA = 0; // multiplicity from FMD on detector A
548 Float_t multFMDC = 0; // multiplicity from FMD on detector C
550 Float_t zvtx =0; // z-vertex SPD
551 Int_t zvtxNcont =0; // contributors to z-vertex SPD
554 AliCentrality *esdCent = 0;
556 if(fAnalysisInput.CompareTo("ESD")==0){
558 AliVEvent* event = InputEvent();
559 AliESDEvent* esd = dynamic_cast<AliESDEvent*>(event);
561 AliError("No ESD Event");
567 if (SetupRun(esd)<0) {
568 AliError("Centrality File not available for this run");
572 esdCent = esd->GetCentrality();
575 AliESDVZERO* esdV0 = esd->GetVZEROData();
576 multV0A=esdV0->GetMTotV0A();
577 multV0C=esdV0->GetMTotV0C();
578 v0Corr = multV0A+multV0C;
580 multV0AOnline=esdV0->GetTriggerChargeA();
581 multV0COnline=esdV0->GetTriggerChargeC();
585 const AliESDVertex* vtxESD = esd->GetPrimaryVertexSPD();
586 zvtx = vtxESD->GetZ();
587 zvtxNcont = vtxESD->GetNContributors();
589 // ***** CB info (tracklets, clusters, chips)
590 //nTracks = event->GetNumberOfTracks();
591 nTracks = fTrackCuts ? (Short_t)fTrackCuts->GetReferenceMultiplicity(esd,kTRUE):-1;
593 const AliMultiplicity *mult = esd->GetMultiplicity();
595 nTracklets = mult->GetNumberOfTracklets();
597 for(Int_t ilay=0; ilay<6; ilay++){
598 nClusters[ilay] = mult->GetNumberOfITSClusters(ilay);
601 for(Int_t ilay=0; ilay<2; ilay++){
602 nChips[ilay] = mult->GetNumberOfFiredChips(ilay);
605 spdCorr = AliESDUtils::GetCorrSPD2(nClusters[1],zvtx);
608 AliESDFMD *fmd = esd->GetFMDData();
609 Float_t totalMultA = 0;
610 Float_t totalMultC = 0;
611 const Float_t fFMDLowCut = 0.4;
613 for(UShort_t det=1;det<=3;det++) {
614 Int_t nRings = (det==1 ? 1 : 2);
615 for (UShort_t ir = 0; ir < nRings; ir++) {
616 Char_t ring = (ir == 0 ? 'I' : 'O');
617 UShort_t nsec = (ir == 0 ? 20 : 40);
618 UShort_t nstr = (ir == 0 ? 512 : 256);
619 for(UShort_t sec =0; sec < nsec; sec++) {
620 for(UShort_t strip = 0; strip < nstr; strip++) {
622 Float_t fmdMult = fmd->Multiplicity(det,ring,sec,strip);
623 if(fmdMult == 0 || fmdMult == AliESDFMD::kInvalidMult) continue;
625 Float_t nParticles=0;
627 if(fmdMult > fFMDLowCut) {
631 if (det<3) totalMultA = totalMultA + nParticles;
632 else totalMultC = totalMultC + nParticles;
638 multFMDA = totalMultA;
639 multFMDC = totalMultC;
642 AliESDZDC *esdZDC = esd->GetESDZDC();
643 zdcEnergyCal = esdZDC->AliESDZDC::TestBit(AliESDZDC::kEnergyCalibratedSignal);
645 zncEnergy = (Float_t) (esdZDC->GetZDCN1Energy());
646 zpcEnergy = (Float_t) (esdZDC->GetZDCP1Energy());
647 znaEnergy = (Float_t) (esdZDC->GetZDCN2Energy());
648 zpaEnergy = (Float_t) (esdZDC->GetZDCP2Energy());
650 zncEnergy = (Float_t) (esdZDC->GetZDCN1Energy())/8.;
651 zpcEnergy = (Float_t) (esdZDC->GetZDCP1Energy())/8.;
652 znaEnergy = (Float_t) (esdZDC->GetZDCN2Energy())/8.;
653 zpaEnergy = (Float_t) (esdZDC->GetZDCP2Energy())/8.;
655 zem1Energy = (Float_t) (esdZDC->GetZDCEMEnergy(0))/8.;
656 zem2Energy = (Float_t) (esdZDC->GetZDCEMEnergy(1))/8.;
660 else if(fAnalysisInput.CompareTo("AOD")==0){
661 //AliAODEvent *aod = dynamic_cast<AliAODEvent*> (InputEvent());
663 printf(" AOD analysis not yet implemented!!!\n\n");
667 // ***** Scaling for MC
670 v0Corr = Short_t((multV0A+multV0C) * fV0MScaleFactorMC);
672 // ***** Scaling for Data
674 v0Corr = Short_t(v0Corr / fV0MScaleFactor);
675 spdCorr = spdCorr / fSPDScaleFactor;
676 nTracks = Int_t(nTracks / fTPCScaleFactor);
679 // ***** Centrality Selection
680 if(fHtempV0M) fCentV0M = fHtempV0M->GetBinContent(fHtempV0M->FindBin((v0Corr)));
681 if(fHtempFMD) fCentFMD = fHtempFMD->GetBinContent(fHtempFMD->FindBin((multFMDA+multFMDC)));
682 if(fHtempTRK) fCentTRK = fHtempTRK->GetBinContent(fHtempTRK->FindBin(nTracks));
683 if(fHtempTKL) fCentTKL = fHtempTKL->GetBinContent(fHtempTKL->FindBin(nTracklets));
684 if(fHtempCL0) fCentCL0 = fHtempCL0->GetBinContent(fHtempCL0->FindBin(nClusters[0]));
685 if(fHtempCL1) fCentCL1 = fHtempCL1->GetBinContent(fHtempCL1->FindBin(spdCorr));
687 if(fHtempV0MvsFMD) fCentV0MvsFMD = fHtempV0MvsFMD->GetBinContent(fHtempV0MvsFMD->FindBin((multV0A+multV0C)));
688 if(fHtempTKLvsV0M) fCentTKLvsV0M = fHtempTKLvsV0M->GetBinContent(fHtempTKLvsV0M->FindBin(nTracklets));
689 if(fHtempZEMvsZDC) fCentZEMvsZDC = fHtempZEMvsZDC->GetBinContent(fHtempZEMvsZDC->FindBin(zem1Energy+zem2Energy,zncEnergy+znaEnergy+zpcEnergy+zpaEnergy));
697 if (TMath::Abs(zvtx)>fZVCut || zvtxNcont<1) fQuality += 1;
699 // ***** outliers, skip in case of MC input
702 if (IsOutlierV0MSPD(spdCorr, v0Corr, int(fCentV0M))) fQuality += 2;
704 if (IsOutlierV0MTPC(nTracks, v0Corr, int(fCentV0M))) fQuality += 4;
706 if (IsOutlierV0MZDC((zncEnergy+znaEnergy+zpcEnergy+zpaEnergy), v0Corr) &&
707 (zdcEnergyCal==kFALSE) ) fQuality += 8;
708 if (IsOutlierV0MZDCECal((zncEnergy+znaEnergy+zpcEnergy+zpaEnergy), v0Corr) &&
709 (zdcEnergyCal==kTRUE) ) fQuality += 8;
717 esdCent->SetQuality(fQuality);
718 esdCent->SetCentralityV0M(fCentV0M);
719 esdCent->SetCentralityFMD(fCentFMD);
720 esdCent->SetCentralityTRK(fCentTRK);
721 esdCent->SetCentralityTKL(fCentTKL);
722 esdCent->SetCentralityCL0(fCentCL0);
723 esdCent->SetCentralityCL1(fCentCL1);
724 esdCent->SetCentralityV0MvsFMD(fCentV0MvsFMD);
725 esdCent->SetCentralityTKLvsV0M(fCentTKLvsV0M);
726 esdCent->SetCentralityZEMvsZDC(fCentZEMvsZDC);
729 fHOutQuality->Fill(fQuality);
730 fHOutVertex->Fill(zvtx);
733 fHOutCentV0M->Fill(fCentV0M);
734 fHOutCentFMD->Fill(fCentFMD);
735 fHOutCentTRK->Fill(fCentTRK);
736 fHOutCentTKL->Fill(fCentTKL);
737 fHOutCentCL0->Fill(fCentCL0);
738 fHOutCentCL1->Fill(fCentCL1);
739 fHOutCentV0MvsFMD->Fill(fCentV0MvsFMD);
740 fHOutCentTKLvsV0M->Fill(fCentTKLvsV0M);
741 fHOutCentZEMvsZDC->Fill(fCentZEMvsZDC);
742 fHOutCentV0MvsCentCL1->Fill(fCentV0M,fCentCL1);
743 fHOutCentV0MvsCentTRK->Fill(fCentV0M,fCentTRK);
744 fHOutCentTRKvsCentCL1->Fill(fCentTRK,fCentCL1);
745 fHOutCentV0MvsCentZDC->Fill(fCentV0M,fCentZEMvsZDC);
746 fHOutMultV0M->Fill(multV0A+multV0C);
747 fHOutMultV0O->Fill(multV0AOnline+multV0COnline);
748 fHOutMultFMD->Fill(multFMDA+multFMDC);
749 fHOutMultTRK->Fill(nTracks);
750 fHOutMultTKL->Fill(nTracklets);
751 fHOutMultCL0->Fill(nClusters[0]);
752 fHOutMultCL1->Fill(spdCorr);
753 fHOutMultV0MvsZDN->Fill(v0Corr,(zncEnergy+znaEnergy));
754 fHOutMultZEMvsZDN->Fill((zem1Energy+zem2Energy),(zncEnergy+znaEnergy));
755 fHOutMultV0MvsZDC->Fill(v0Corr,(zncEnergy+znaEnergy+zpcEnergy+zpaEnergy));
756 fHOutMultZEMvsZDC->Fill((zem1Energy+zem2Energy),(zncEnergy+znaEnergy+zpcEnergy+zpaEnergy));
757 fHOutMultZEMvsZDCw->Fill((zem1Energy+zem2Energy),(zncEnergy+znaEnergy+zpcEnergy+zpaEnergy),fCentV0M);
758 fHOutMultV0MvsCL1->Fill(v0Corr,spdCorr);
759 fHOutMultV0MvsTRK->Fill(v0Corr,nTracks);
760 fHOutMultTRKvsCL1->Fill(nTracks,spdCorr);
761 fHOutMultV0MvsV0O->Fill(v0Corr,(multV0AOnline+multV0COnline));
762 fHOutMultV0OvsCL1->Fill((multV0AOnline+multV0COnline),spdCorr);
763 fHOutMultV0OvsTRK->Fill((multV0AOnline+multV0COnline),nTracks);
764 } else if (fQuality%2 == 0) {
765 fHOutCentV0Mqual1->Fill(fCentV0M);
766 fHOutCentTRKqual1->Fill(fCentTRK);
767 fHOutCentCL1qual1->Fill(fCentCL1);
768 fHOutMultV0MvsCL1qual1->Fill(v0Corr,spdCorr);
769 fHOutMultV0MvsTRKqual1->Fill(v0Corr,nTracks);
770 fHOutMultTRKvsCL1qual1->Fill(nTracks,spdCorr);
772 fHOutCentV0Mqual2->Fill(fCentV0M);
773 fHOutCentTRKqual2->Fill(fCentTRK);
774 fHOutCentCL1qual2->Fill(fCentCL1);
775 fHOutMultV0MvsCL1qual2->Fill(v0Corr,spdCorr);
776 fHOutMultV0MvsTRKqual2->Fill(v0Corr,nTracks);
777 fHOutMultTRKvsCL1qual2->Fill(nTracks,spdCorr);
780 PostData(1, fOutputList);
782 //________________________________________________________________________
783 void AliCentralitySelectionTask::Terminate(Option_t */*option*/)
785 // Terminate analysis
787 //________________________________________________________________________
788 Int_t AliCentralitySelectionTask::SetupRun(AliESDEvent* const esd)
790 // Setup files for run
795 // check if something to be done
796 if (fCurrentRun == esd->GetRunNumber())
799 fCurrentRun = esd->GetRunNumber();
801 TString fileName =(Form("%s/COMMON/CENTRALITY/data/centrality.root", AliAnalysisManager::GetOADBPath()));
802 AliInfo(Form("Setup Centrality Selection for run %d with file %s\n",fCurrentRun,fileName.Data()));
804 AliOADBContainer *con = new AliOADBContainer("OADB");
805 con->InitFromFile(fileName,"Centrality");
807 AliOADBCentrality* centOADB = 0;
808 centOADB = (AliOADBCentrality*)(con->GetObject(fCurrentRun));
810 AliWarning(Form("Centrality OADB does not exist for run %d, using Default \n",fCurrentRun ));
811 centOADB = (AliOADBCentrality*)(con->GetDefaultObject("oadbDefault"));
815 fUseScaling = centOADB->UseScaling();
816 fUseCleaning = centOADB->UseCleaning();
819 fZVCut = centOADB->ZVCut();
820 fOutliersCut = centOADB->OutliersCut();
823 fHtempV0M = centOADB->V0hist();
824 fHtempTRK = centOADB->TPChist();
825 fHtempCL1 = centOADB->SPDhist();
826 fHtempZEMvsZDC = centOADB->ZEMvsZDChist();
828 TString path = gSystem->ExpandPathName(fileName.Data());
829 if (!fHtempV0M) AliWarning(Form("Calibration for V0M does not exist in %s", path.Data()));
830 if (!fHtempTRK) AliWarning(Form("Calibration for TRK does not exist in %s", path.Data()));
831 if (!fHtempCL1) AliWarning(Form("Calibration for CL1 does not exist in %s", path.Data()));
832 if (!fHtempZEMvsZDC) AliWarning(Form("Calibration for ZEMvsZDC does not exist in %s", path.Data()));
835 fV0MScaleFactor = centOADB->V0MScaleFactor();
836 fSPDScaleFactor = centOADB->SPDScaleFactor();
837 fTPCScaleFactor = centOADB->TPCScaleFactor();
838 fV0MScaleFactorMC = centOADB->V0MScaleFactor();
840 // outliers parameters
841 fV0MSPDOutlierPar0 = centOADB->V0MSPDOutlierPar0();
842 fV0MSPDOutlierPar1 = centOADB->V0MSPDOutlierPar1();
843 fV0MTPCOutlierPar0 = centOADB->V0MTPCOutlierPar0();
844 fV0MTPCOutlierPar1 = centOADB->V0MTPCOutlierPar1();
846 fV0MSPDSigmaOutlierPar0 = centOADB->V0MSPDSigmaOutlierPar0();
847 fV0MSPDSigmaOutlierPar1 = centOADB->V0MSPDSigmaOutlierPar1();
848 fV0MSPDSigmaOutlierPar2 = centOADB->V0MSPDSigmaOutlierPar2();
849 fV0MTPCSigmaOutlierPar0 = centOADB->V0MTPCSigmaOutlierPar0();
850 fV0MTPCSigmaOutlierPar1 = centOADB->V0MTPCSigmaOutlierPar1();
851 fV0MTPCSigmaOutlierPar2 = centOADB->V0MTPCSigmaOutlierPar2();
853 fV0MZDCOutlierPar0 = centOADB->V0MZDCOutlierPar0();
854 fV0MZDCOutlierPar1 = centOADB->V0MZDCOutlierPar1();
855 fV0MZDCEcalOutlierPar0 = centOADB->V0MZDCEcalOutlierPar0();
856 fV0MZDCEcalOutlierPar1 = centOADB->V0MZDCEcalOutlierPar1();
865 //________________________________________________________________________
866 Bool_t AliCentralitySelectionTask::IsOutlierV0MSPD(Float_t spd, Float_t v0, Int_t cent) const
869 Float_t val = fV0MSPDOutlierPar0 + fV0MSPDOutlierPar1 * v0;
870 Float_t spdSigma = fV0MSPDSigmaOutlierPar0 + fV0MSPDSigmaOutlierPar1*cent + fV0MSPDSigmaOutlierPar2*cent*cent;
871 if ( TMath::Abs(spd-val) > fOutliersCut*spdSigma )
877 //________________________________________________________________________
878 Bool_t AliCentralitySelectionTask::IsOutlierV0MTPC(Int_t tracks, Float_t v0, Int_t cent) const
881 Float_t val = fV0MTPCOutlierPar0 + fV0MTPCOutlierPar1 * v0;
882 Float_t tpcSigma = fV0MTPCSigmaOutlierPar0 + fV0MTPCSigmaOutlierPar1*cent + fV0MTPCSigmaOutlierPar2*cent*cent;
883 if ( TMath::Abs(tracks-val) > fOutliersCut*tpcSigma )
889 //________________________________________________________________________
890 Bool_t AliCentralitySelectionTask::IsOutlierV0MZDC(Float_t zdc, Float_t v0) const
893 Float_t val = fV0MZDCOutlierPar0 + fV0MZDCOutlierPar1 * v0;
900 //________________________________________________________________________
901 Bool_t AliCentralitySelectionTask::IsOutlierV0MZDCECal(Float_t zdc, Float_t v0) const
904 Float_t val = fV0MZDCEcalOutlierPar0 + fV0MZDCEcalOutlierPar1 * v0;