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),
131 fHOutCentV0M_CVHN(0),
132 fHOutCentV0M_CVLN(0),
138 fHOutCentV0MvsFMD(0),
139 fHOutCentTKLvsV0M(0),
140 fHOutCentZEMvsZDC(0),
141 fHOutCentV0MvsCentCL1(0),
142 fHOutCentV0MvsCentTRK(0),
143 fHOutCentTRKvsCentCL1(0),
144 fHOutCentV0MvsCentZDC(0),
152 fHOutMultV0MvsZDN(0),
153 fHOutMultZEMvsZDN(0),
154 fHOutMultV0MvsZDC(0),
155 fHOutMultZEMvsZDC(0),
156 fHOutMultZEMvsZDCw(0),
157 fHOutMultV0MvsCL1(0),
158 fHOutMultV0MvsTRK(0),
159 fHOutMultTRKvsCL1(0),
160 fHOutMultV0MvsV0O(0),
161 fHOutMultV0OvsCL1(0),
162 fHOutMultV0OvsTRK(0),
163 fHOutCentV0Mqual1(0),
164 fHOutCentTRKqual1(0),
165 fHOutCentCL1qual1(0),
166 fHOutMultV0MvsCL1qual1(0),
167 fHOutMultV0MvsTRKqual1(0),
168 fHOutMultTRKvsCL1qual1(0),
169 fHOutCentV0Mqual2(0),
170 fHOutCentTRKqual2(0),
171 fHOutCentCL1qual2(0),
172 fHOutMultV0MvsCL1qual2(0),
173 fHOutMultV0MvsTRKqual2(0),
174 fHOutMultTRKvsCL1qual2(0),
178 // Default constructor
179 AliInfo("Centrality Selection enabled.");
183 fBranchNames="ESD:AliESDRun.,AliESDHeader.,AliESDZDC.,AliESDFMD.,AliESDVZERO."
184 ",SPDVertex.,TPCVertex.,PrimaryVertex.,AliMultiplicity.,Tracks ";
187 //________________________________________________________________________
188 AliCentralitySelectionTask::AliCentralitySelectionTask(const char *name):
189 AliAnalysisTaskSE(name),
190 fAnalysisInput("ESD"),
199 fV0MScaleFactorMC(0),
200 fV0MSPDOutlierPar0(0),
201 fV0MSPDOutlierPar1(0),
202 fV0MTPCOutlierPar0(0),
203 fV0MTPCOutlierPar1(0),
204 fV0MSPDSigmaOutlierPar0(0),
205 fV0MSPDSigmaOutlierPar1(0),
206 fV0MSPDSigmaOutlierPar2(0),
207 fV0MTPCSigmaOutlierPar0(0),
208 fV0MTPCSigmaOutlierPar1(0),
209 fV0MTPCSigmaOutlierPar2(0),
210 fV0MZDCOutlierPar0(0),
211 fV0MZDCOutlierPar1(0),
212 fV0MZDCEcalOutlierPar0(0),
213 fV0MZDCEcalOutlierPar1(0),
240 fHOutCentV0M_CVHN(0),
241 fHOutCentV0M_CVLN(0),
247 fHOutCentV0MvsFMD(0),
248 fHOutCentTKLvsV0M(0),
249 fHOutCentZEMvsZDC(0),
250 fHOutCentV0MvsCentCL1(0),
251 fHOutCentV0MvsCentTRK(0),
252 fHOutCentTRKvsCentCL1(0),
253 fHOutCentV0MvsCentZDC(0),
261 fHOutMultV0MvsZDN(0),
262 fHOutMultZEMvsZDN(0),
263 fHOutMultV0MvsZDC(0),
264 fHOutMultZEMvsZDC(0),
265 fHOutMultZEMvsZDCw(0),
266 fHOutMultV0MvsCL1(0),
267 fHOutMultV0MvsTRK(0),
268 fHOutMultTRKvsCL1(0),
269 fHOutMultV0MvsV0O(0),
270 fHOutMultV0OvsCL1(0),
271 fHOutMultV0OvsTRK(0),
272 fHOutCentV0Mqual1(0),
273 fHOutCentTRKqual1(0),
274 fHOutCentCL1qual1(0),
275 fHOutMultV0MvsCL1qual1(0),
276 fHOutMultV0MvsTRKqual1(0),
277 fHOutMultTRKvsCL1qual1(0),
278 fHOutCentV0Mqual2(0),
279 fHOutCentTRKqual2(0),
280 fHOutCentCL1qual2(0),
281 fHOutMultV0MvsCL1qual2(0),
282 fHOutMultV0MvsTRKqual2(0),
283 fHOutMultTRKvsCL1qual2(0),
287 // Default constructor
288 AliInfo("Centrality Selection enabled.");
289 DefineOutput(1, TList::Class());
292 fBranchNames="ESD:AliESDRun.,AliESDHeader.,AliESDZDC.,AliESDFMD.,AliESDVZERO."
293 ",SPDVertex.,TPCVertex.,PrimaryVertex.,AliMultiplicity.,Tracks ";
296 //________________________________________________________________________
297 AliCentralitySelectionTask& AliCentralitySelectionTask::operator=(const AliCentralitySelectionTask& c)
299 // Assignment operator
301 AliAnalysisTaskSE::operator=(c);
306 //________________________________________________________________________
307 AliCentralitySelectionTask::AliCentralitySelectionTask(const AliCentralitySelectionTask& ana):
308 AliAnalysisTaskSE(ana),
309 fAnalysisInput(ana.fDebug),
310 fIsMCInput(ana.fIsMCInput),
312 fCurrentRun(ana.fCurrentRun),
313 fUseScaling(ana.fUseScaling),
314 fUseCleaning(ana.fUseCleaning),
315 fV0MScaleFactor(ana.fV0MScaleFactor),
316 fSPDScaleFactor(ana.fSPDScaleFactor),
317 fTPCScaleFactor(ana.fTPCScaleFactor),
318 fV0MScaleFactorMC(ana.fV0MScaleFactorMC),
319 fV0MSPDOutlierPar0(ana.fV0MSPDOutlierPar0),
320 fV0MSPDOutlierPar1(ana.fV0MSPDOutlierPar1),
321 fV0MTPCOutlierPar0(ana.fV0MTPCOutlierPar0),
322 fV0MTPCOutlierPar1(ana.fV0MTPCOutlierPar1),
323 fV0MSPDSigmaOutlierPar0(ana.fV0MSPDSigmaOutlierPar0),
324 fV0MSPDSigmaOutlierPar1(ana.fV0MSPDSigmaOutlierPar1),
325 fV0MSPDSigmaOutlierPar2(ana.fV0MSPDSigmaOutlierPar2),
326 fV0MTPCSigmaOutlierPar0(ana.fV0MTPCSigmaOutlierPar0),
327 fV0MTPCSigmaOutlierPar1(ana.fV0MTPCSigmaOutlierPar1),
328 fV0MTPCSigmaOutlierPar2(ana.fV0MTPCSigmaOutlierPar2),
329 fV0MZDCOutlierPar0(ana.fV0MZDCOutlierPar0),
330 fV0MZDCOutlierPar1(ana.fV0MZDCOutlierPar1),
331 fV0MZDCEcalOutlierPar0(ana.fV0MZDCEcalOutlierPar0),
332 fV0MZDCEcalOutlierPar1(ana.fV0MZDCEcalOutlierPar1),
333 fTrackCuts(ana.fTrackCuts),
335 fOutliersCut(ana.fOutliersCut),
336 fQuality(ana.fQuality),
339 fCentV0M(ana.fCentV0M),
340 fCentFMD(ana.fCentFMD),
341 fCentTRK(ana.fCentTRK),
342 fCentTKL(ana.fCentTKL),
343 fCentCL0(ana.fCentCL0),
344 fCentCL1(ana.fCentCL1),
345 fCentV0MvsFMD(ana.fCentV0MvsFMD),
346 fCentTKLvsV0M(ana.fCentTKLvsV0M),
347 fCentZEMvsZDC(ana.fCentZEMvsZDC),
348 fHtempV0M(ana.fHtempV0M),
349 fHtempFMD(ana.fHtempFMD),
350 fHtempTRK(ana.fHtempTRK),
351 fHtempTKL(ana.fHtempTKL),
352 fHtempCL0(ana.fHtempCL0),
353 fHtempCL1(ana.fHtempCL1),
354 fHtempV0MvsFMD(ana.fHtempV0MvsFMD),
355 fHtempTKLvsV0M(ana.fHtempTKLvsV0M),
356 fHtempZEMvsZDC(ana.fHtempZEMvsZDC),
357 fOutputList(ana.fOutputList),
358 fHOutCentV0M (ana.fHOutCentV0M ),
359 fHOutCentV0M_CVHN(ana.fHOutCentV0M_CVHN),
360 fHOutCentV0M_CVLN(ana.fHOutCentV0M_CVLN),
361 fHOutCentFMD (ana.fHOutCentFMD ),
362 fHOutCentTRK (ana.fHOutCentTRK ),
363 fHOutCentTKL (ana.fHOutCentTKL ),
364 fHOutCentCL0 (ana.fHOutCentCL0 ),
365 fHOutCentCL1 (ana.fHOutCentCL1 ),
366 fHOutCentV0MvsFMD(ana.fHOutCentV0MvsFMD),
367 fHOutCentTKLvsV0M(ana.fHOutCentTKLvsV0M),
368 fHOutCentZEMvsZDC(ana.fHOutCentZEMvsZDC),
369 fHOutCentV0MvsCentCL1(ana.fHOutCentV0MvsCentCL1),
370 fHOutCentV0MvsCentTRK(ana.fHOutCentV0MvsCentTRK),
371 fHOutCentTRKvsCentCL1(ana.fHOutCentTRKvsCentCL1),
372 fHOutCentV0MvsCentZDC(ana.fHOutCentV0MvsCentZDC),
373 fHOutMultV0M(ana.fHOutMultV0M),
374 fHOutMultV0O(ana.fHOutMultV0O),
375 fHOutMultFMD(ana.fHOutMultFMD),
376 fHOutMultTRK(ana.fHOutMultTRK),
377 fHOutMultTKL(ana.fHOutMultTKL),
378 fHOutMultCL0(ana.fHOutMultCL0),
379 fHOutMultCL1(ana.fHOutMultCL1),
380 fHOutMultV0MvsZDN(ana.fHOutMultV0MvsZDN),
381 fHOutMultZEMvsZDN(ana.fHOutMultZEMvsZDN),
382 fHOutMultV0MvsZDC(ana.fHOutMultV0MvsZDC),
383 fHOutMultZEMvsZDC(ana.fHOutMultZEMvsZDC),
384 fHOutMultZEMvsZDCw(ana.fHOutMultZEMvsZDCw),
385 fHOutMultV0MvsCL1(ana.fHOutMultV0MvsCL1),
386 fHOutMultV0MvsTRK(ana.fHOutMultV0MvsTRK),
387 fHOutMultTRKvsCL1(ana.fHOutMultTRKvsCL1),
388 fHOutMultV0MvsV0O(ana.fHOutMultV0MvsV0O),
389 fHOutMultV0OvsCL1(ana.fHOutMultV0OvsCL1),
390 fHOutMultV0OvsTRK(ana.fHOutMultV0OvsTRK),
391 fHOutCentV0Mqual1(ana.fHOutCentV0Mqual1),
392 fHOutCentTRKqual1(ana.fHOutCentTRKqual1),
393 fHOutCentCL1qual1(ana.fHOutCentCL1qual1),
394 fHOutMultV0MvsCL1qual1(ana.fHOutMultV0MvsCL1qual1),
395 fHOutMultV0MvsTRKqual1(ana.fHOutMultV0MvsTRKqual1),
396 fHOutMultTRKvsCL1qual1(ana.fHOutMultTRKvsCL1qual1),
397 fHOutCentV0Mqual2(ana.fHOutCentV0Mqual2),
398 fHOutCentTRKqual2(ana.fHOutCentTRKqual2),
399 fHOutCentCL1qual2(ana.fHOutCentCL1qual2),
400 fHOutMultV0MvsCL1qual2(ana.fHOutMultV0MvsCL1qual2),
401 fHOutMultV0MvsTRKqual2(ana.fHOutMultV0MvsTRKqual2),
402 fHOutMultTRKvsCL1qual2(ana.fHOutMultTRKvsCL1qual2),
403 fHOutQuality(ana.fHOutQuality),
404 fHOutVertex(ana.fHOutVertex)
410 //________________________________________________________________________
411 AliCentralitySelectionTask::~AliCentralitySelectionTask()
414 if (fOutputList && !AliAnalysisManager::GetAnalysisManager()->IsProofMode()) delete fOutputList;
415 if (fTrackCuts) delete fTrackCuts;
418 //________________________________________________________________________
419 void AliCentralitySelectionTask::UserCreateOutputObjects()
421 // Create the output containers
422 if(fDebug>1) printf("AnalysisCentralitySelectionTask::UserCreateOutputObjects() \n");
423 AliLog::SetClassDebugLevel("AliCentralitySelectionTask", AliLog::kInfo);
425 fOutputList = new TList();
426 fOutputList->SetOwner();
427 fHOutCentV0M = new TH1F("fHOutCentV0M","fHOutCentV0M; Centrality V0",505,0,101);
428 fHOutCentV0M_CVHN= new TH1F("fHOutCentV0M_CVHN","fHOutCentV0M_CVHN; Centrality V0",505,0,101);
429 fHOutCentV0M_CVLN= new TH1F("fHOutCentV0M_CVLN","fHOutCentV0M_CVLN; Centrality V0",505,0,101);
430 fHOutCentFMD = new TH1F("fHOutCentFMD","fHOutCentFMD; Centrality FMD",505,0,101);
431 fHOutCentTRK = new TH1F("fHOutCentTRK","fHOutCentTRK; Centrality TPC",505,0,101);
432 fHOutCentTKL = new TH1F("fHOutCentTKL","fHOutCentTKL; Centrality tracklets",505,0,101);
433 fHOutCentCL0 = new TH1F("fHOutCentCL0","fHOutCentCL0; Centrality SPD inner",505,0,101);
434 fHOutCentCL1 = new TH1F("fHOutCentCL1","fHOutCentCL1; Centrality SPD outer",505,0,101);
435 fHOutCentV0MvsFMD= new TH1F("fHOutCentV0MvsFMD","fHOutCentV0MvsFMD; Centrality V0 vs FMD",505,0,101);
436 fHOutCentTKLvsV0M= new TH1F("fHOutCentTKLvsV0M","fHOutCentTKLvsV0M; Centrality tracklets vs V0",505,0,101);
437 fHOutCentZEMvsZDC= new TH1F("fHOutCentZEMvsZDC","fHOutCentZEMvsZDC; Centrality ZEM vs ZDC",505,0,101);
438 fHOutCentV0MvsCentCL1= new TH2F("fHOutCentV0MvsCentCL1","fHOutCentV0MvsCentCL1; Cent V0; Cent SPD",505,0,101,505,0,101);
439 fHOutCentV0MvsCentTRK= new TH2F("fHOutCentV0MvsCentTRK","fHOutCentV0MvsCentTRK; Cent V0; Cent TPC",505,0,101,505,0,101);
440 fHOutCentTRKvsCentCL1= new TH2F("fHOutCentTRKvsCentCL1","fHOutCentTRKvsCentCL1; Cent TPC; Cent SPD",505,0,101,505,0,101);
441 fHOutCentV0MvsCentZDC= new TH2F("fHOutCentV0MvsCentZDC","fHOutCentV0MvsCentZDC; Cent V0; Cent ZDC",505,0,101,505,0,101);
442 fHOutMultV0M = new TH1F("fHOutMultV0M","fHOutMultV0M; Multiplicity V0",25000,0,30000);
443 fHOutMultV0O = new TH1F("fHOutMultV0O","fHOutMultV0O; Multiplicity V0",40000,0,40000);
444 fHOutMultFMD = new TH1F("fHOutMultFMD","fHOutMultFMD; Multiplicity FMD",24000,0,24000);
445 fHOutMultTRK = new TH1F("fHOutMultTRK","fHOutMultTRK; Multiplicity TPC",4000,0,4000);
446 fHOutMultTKL = new TH1F("fHOutMultTKL","fHOutMultTKL; Multiplicity tracklets",5000,0,5000);
447 fHOutMultCL0 = new TH1F("fHOutMultCL0","fHOutMultCL0; Multiplicity SPD inner",7000,0,7000);
448 fHOutMultCL1 = new TH1F("fHOutMultCL1","fHOutMultCL1; Multiplicity SPD outer",7000,0,7000);
449 fHOutMultV0MvsZDN = new TH2F("fHOutMultV0MvsZDN","fHOutMultV0MvsZDN; Multiplicity V0; Energy ZDC-N",500,0,30000,500,0,180000);
450 fHOutMultZEMvsZDN = new TH2F("fHOutMultZEMvsZDN","fHOutMultZEMvsZDN; Energy ZEM; Energy ZDC-N",500,0,2500,500,0,180000);
451 fHOutMultV0MvsZDC = new TH2F("fHOutMultV0MvsZDC","fHOutMultV0MvsZDC; Multiplicity V0; Energy ZDC",500,0,30000,500,0,200000);
452 fHOutMultZEMvsZDC = new TH2F("fHOutMultZEMvsZDC","fHOutMultZEMvsZDC; Energy ZEM; Energy ZDC",500,0,2500,500,0,200000);
453 fHOutMultZEMvsZDCw = new TH2F("fHOutMultZEMvsZDCw","fHOutMultZEMvsZDCw; Energy ZEM; Energy ZDC (weigthed with V0 percentile)",500,0,2500,500,0,200000);
454 fHOutMultV0MvsCL1 = new TH2F("fHOutMultV0MvsCL1","fHOutMultV0MvsCL1; Multiplicity V0; Multiplicity SPD outer",2500,0,30000,700,0,7000);
455 fHOutMultV0MvsTRK = new TH2F("fHOutMultV0MvsTRK","fHOutMultV0MvsTRK; Multiplicity V0; Multiplicity TPC",2500,0,30000,400,0,4000);
456 fHOutMultTRKvsCL1 = new TH2F("fHOutMultTRKvsCL1","fHOutMultTRKvsCL1; Multiplicity TPC; Multiplicity SPD outer",400,0,4000,700,0,7000);
457 fHOutMultV0MvsV0O = new TH2F("fHOutMultV0MvsV0O","fHOutMultV0MvsV0O; Multiplicity V0; Multiplicity V0 Online",500,0,30000,500,0,40000);
458 fHOutMultV0OvsCL1 = new TH2F("fHOutMultV0OvsCL1","fHOutMultV0OvsCL1; Multiplicity V0; Multiplicity SPD outer",500,0,40000,700,0,7000);
459 fHOutMultV0OvsTRK = new TH2F("fHOutMultV0OvsTRK","fHOutMultV0OvsTRK; Multiplicity V0; Multiplicity TPC",500,0,40000,400,0,4000);
460 fHOutMultV0MvsV0O = new TH2F("fHOutMultV0MvsV0O","fHOutMultV0MvsV0O; Multiplicity V0; Multiplicity V0 Online",500,0,30000,500,0,30000);
461 fHOutMultV0OvsCL1 = new TH2F("fHOutMultV0OvsCL1","fHOutMultV0OvsCL1; Multiplicity V0; Multiplicity SPD outer",2500,0,30000,700,0,7000);
462 fHOutMultV0OvsTRK = new TH2F("fHOutMultV0OvsTRK","fHOutMultV0OvsTRK; Multiplicity V0; Multiplicity TPC",2500,0,30000,400,0,4000);
464 fHOutCentV0Mqual1 = new TH1F("fHOutCentV0M_qual1","fHOutCentV0M_qual1; Centrality V0",505,0,101);
465 fHOutCentTRKqual1 = new TH1F("fHOutCentTRK_qual1","fHOutCentTRK_qual1; Centrality TPC",505,0,101);
466 fHOutCentCL1qual1 = new TH1F("fHOutCentCL1_qual1","fHOutCentCL1_qual1; Centrality SPD outer",505,0,101);
467 fHOutMultV0MvsCL1qual1 = new TH2F("fHOutMultV0MvsCL1_qual1","fHOutMultV0MvsCL1_qual1; Multiplicity V0; Multiplicity SPD outer",2500,0,25000,700,0,7000);
468 fHOutMultV0MvsTRKqual1 = new TH2F("fHOutMultV0MvsTRK_qual1","fHOutMultV0MvsTRK_qual1; Multiplicity V0; Multiplicity TPC",2500,0,25000,400,0,4000);
469 fHOutMultTRKvsCL1qual1 = new TH2F("fHOutMultTRKvsCL1_qual1","fHOutMultTRKvsCL1_qual1; Multiplicity TPC; Multiplicity SPD outer",400,0,4000,700,0,7000);
471 fHOutCentV0Mqual2 = new TH1F("fHOutCentV0M_qual2","fHOutCentV0M_qual2; Centrality V0",505,0,101);
472 fHOutCentTRKqual2 = new TH1F("fHOutCentTRK_qual2","fHOutCentTRK_qual2; Centrality TPC",505,0,101);
473 fHOutCentCL1qual2 = new TH1F("fHOutCentCL1_qual2","fHOutCentCL1_qual2; Centrality SPD outer",505,0,101);
474 fHOutMultV0MvsCL1qual2 = new TH2F("fHOutMultV0MvsCL1_qual2","fHOutMultV0MvsCL1_qual2; Multiplicity V0; Multiplicity SPD outer",2500,0,25000,700,0,7000);
475 fHOutMultV0MvsTRKqual2 = new TH2F("fHOutMultV0MvsTRK_qual2","fHOutMultV0MvsTRK_qual2; Multiplicity V0; Multiplicity TPC",2500,0,25000,400,0,4000);
476 fHOutMultTRKvsCL1qual2 = new TH2F("fHOutMultTRKvsCL1_qual2","fHOutMultTRKvsCL1_qual2; Multiplicity TPC; Multiplicity SPD outer",400,0,4000,700,0,7000);
478 fHOutQuality = new TH1F("fHOutQuality", "fHOutQuality", 100,-0.5,99.5);
479 fHOutVertex = new TH1F("fHOutVertex", "fHOutVertex", 100,-20,20);
481 fOutputList->Add( fHOutCentV0M );
482 fOutputList->Add( fHOutCentV0M_CVHN);
483 fOutputList->Add( fHOutCentV0M_CVLN);
484 fOutputList->Add( fHOutCentFMD );
485 fOutputList->Add( fHOutCentTRK );
486 fOutputList->Add( fHOutCentTKL );
487 fOutputList->Add( fHOutCentCL0 );
488 fOutputList->Add( fHOutCentCL1 );
489 fOutputList->Add( fHOutCentV0MvsFMD);
490 fOutputList->Add( fHOutCentTKLvsV0M);
491 fOutputList->Add( fHOutCentZEMvsZDC);
492 fOutputList->Add( fHOutCentV0MvsCentCL1);
493 fOutputList->Add( fHOutCentV0MvsCentTRK);
494 fOutputList->Add( fHOutCentTRKvsCentCL1);
495 fOutputList->Add( fHOutCentV0MvsCentZDC);
496 fOutputList->Add( fHOutMultV0M);
497 fOutputList->Add( fHOutMultV0O);
498 fOutputList->Add( fHOutMultFMD);
499 fOutputList->Add( fHOutMultTRK);
500 fOutputList->Add( fHOutMultTKL);
501 fOutputList->Add( fHOutMultCL0);
502 fOutputList->Add( fHOutMultCL1);
503 fOutputList->Add( fHOutMultV0MvsZDN);
504 fOutputList->Add( fHOutMultZEMvsZDN);
505 fOutputList->Add( fHOutMultV0MvsZDC);
506 fOutputList->Add( fHOutMultZEMvsZDC);
507 fOutputList->Add( fHOutMultZEMvsZDCw);
508 fOutputList->Add( fHOutMultV0MvsCL1);
509 fOutputList->Add( fHOutMultV0MvsTRK);
510 fOutputList->Add( fHOutMultTRKvsCL1);
511 fOutputList->Add( fHOutMultV0MvsV0O);
512 fOutputList->Add( fHOutMultV0OvsCL1);
513 fOutputList->Add( fHOutMultV0OvsTRK);
514 fOutputList->Add( fHOutCentV0Mqual1 );
515 fOutputList->Add( fHOutCentTRKqual1 );
516 fOutputList->Add( fHOutCentCL1qual1 );
517 fOutputList->Add( fHOutMultV0MvsCL1qual1);
518 fOutputList->Add( fHOutMultV0MvsTRKqual1);
519 fOutputList->Add( fHOutMultTRKvsCL1qual1);
520 fOutputList->Add( fHOutCentV0Mqual2 );
521 fOutputList->Add( fHOutCentTRKqual2 );
522 fOutputList->Add( fHOutCentCL1qual2 );
523 fOutputList->Add( fHOutMultV0MvsCL1qual2);
524 fOutputList->Add( fHOutMultV0MvsTRKqual2);
525 fOutputList->Add( fHOutMultTRKvsCL1qual2);
526 fOutputList->Add( fHOutQuality );
527 fOutputList->Add( fHOutVertex );
530 fTrackCuts = AliESDtrackCuts::GetStandardTPCOnlyTrackCuts();
532 PostData(1, fOutputList);
534 if (fPass==0) AliFatal("Which pass are you analyzing? You should set it via taskCentrality->SetPass(N)");
537 //________________________________________________________________________
538 void AliCentralitySelectionTask::UserExec(Option_t */*option*/)
540 // Execute analysis for current event:
541 if(fDebug>1) printf(" **** AliCentralitySelectionTask::UserExec() \n");
543 Float_t zncEnergy = 0.; // ZNC Energy
544 Float_t zpcEnergy = 0.; // ZPC Energy
545 Float_t znaEnergy = 0.; // ZNA Energy
546 Float_t zpaEnergy = 0.; // ZPA Energy
547 Float_t zem1Energy = 0.; // ZEM1 Energy
548 Float_t zem2Energy = 0.; // ZEM2 Energy
549 Bool_t zdcEnergyCal = kFALSE; // if zdc is calibrated (in pass2)
551 Int_t nTracks = 0; // no. tracks
552 Int_t nTracklets = 0; // no. tracklets
553 Int_t nClusters[6] = {0}; // no. clusters on 6 ITS layers
554 Int_t nChips[2]; // no. chips on 2 SPD layers
555 Float_t spdCorr =0; // corrected spd2 multiplicity
557 Float_t multV0A = 0; // multiplicity from V0 reco side A
558 Float_t multV0C = 0; // multiplicity from V0 reco side C
559 Short_t multV0AOnline = 0; // multiplicity from V0 reco side A
560 Short_t multV0COnline = 0; // multiplicity from V0 reco side C
561 Float_t v0Corr = 0; // corrected V0 multiplicity (used for MC)
563 Float_t multFMDA = 0; // multiplicity from FMD on detector A
564 Float_t multFMDC = 0; // multiplicity from FMD on detector C
566 Float_t zvtx =0; // z-vertex SPD
567 Int_t zvtxNcont =0; // contributors to z-vertex SPD
570 AliCentrality *esdCent = 0;
572 if(fAnalysisInput.CompareTo("ESD")==0){
574 AliVEvent* event = InputEvent();
575 AliESDEvent* esd = dynamic_cast<AliESDEvent*>(event);
577 AliError("No ESD Event");
583 if (SetupRun(esd)<0) {
584 AliError("Centrality File not available for this run");
588 esdCent = esd->GetCentrality();
591 AliESDVZERO* esdV0 = esd->GetVZEROData();
592 multV0A=esdV0->GetMTotV0A();
593 multV0C=esdV0->GetMTotV0C();
594 v0Corr = multV0A+multV0C;
596 multV0AOnline=esdV0->GetTriggerChargeA();
597 multV0COnline=esdV0->GetTriggerChargeC();
600 // ***** Trigger info
601 TString trigStr(esd->GetFiredTriggerClasses());
602 fCVHN=kFALSE; fCVLN=kFALSE;
603 if ( (trigStr.Contains("-B-")) && (trigStr.Contains("CVHN")))
605 if ( (trigStr.Contains("-B-")) && (trigStr.Contains("CVLN")))
610 const AliESDVertex* vtxESD = esd->GetPrimaryVertexSPD();
611 zvtx = vtxESD->GetZ();
612 zvtxNcont = vtxESD->GetNContributors();
614 // ***** CB info (tracklets, clusters, chips)
615 //nTracks = event->GetNumberOfTracks();
616 nTracks = fTrackCuts ? (Short_t)fTrackCuts->GetReferenceMultiplicity(esd,kTRUE):-1;
618 const AliMultiplicity *mult = esd->GetMultiplicity();
620 nTracklets = mult->GetNumberOfTracklets();
622 for(Int_t ilay=0; ilay<6; ilay++){
623 nClusters[ilay] = mult->GetNumberOfITSClusters(ilay);
626 for(Int_t ilay=0; ilay<2; ilay++){
627 nChips[ilay] = mult->GetNumberOfFiredChips(ilay);
630 spdCorr = AliESDUtils::GetCorrSPD2(nClusters[1],zvtx);
633 AliESDFMD *fmd = esd->GetFMDData();
634 Float_t totalMultA = 0;
635 Float_t totalMultC = 0;
636 const Float_t fFMDLowCut = 0.4;
638 for(UShort_t det=1;det<=3;det++) {
639 Int_t nRings = (det==1 ? 1 : 2);
640 for (UShort_t ir = 0; ir < nRings; ir++) {
641 Char_t ring = (ir == 0 ? 'I' : 'O');
642 UShort_t nsec = (ir == 0 ? 20 : 40);
643 UShort_t nstr = (ir == 0 ? 512 : 256);
644 for(UShort_t sec =0; sec < nsec; sec++) {
645 for(UShort_t strip = 0; strip < nstr; strip++) {
647 Float_t fmdMult = fmd->Multiplicity(det,ring,sec,strip);
648 if(fmdMult == 0 || fmdMult == AliESDFMD::kInvalidMult) continue;
650 Float_t nParticles=0;
652 if(fmdMult > fFMDLowCut) {
656 if (det<3) totalMultA = totalMultA + nParticles;
657 else totalMultC = totalMultC + nParticles;
663 multFMDA = totalMultA;
664 multFMDC = totalMultC;
667 AliESDZDC *esdZDC = esd->GetESDZDC();
668 zdcEnergyCal = esdZDC->AliESDZDC::TestBit(AliESDZDC::kEnergyCalibratedSignal);
670 zncEnergy = (Float_t) (esdZDC->GetZDCN1Energy());
671 zpcEnergy = (Float_t) (esdZDC->GetZDCP1Energy());
672 znaEnergy = (Float_t) (esdZDC->GetZDCN2Energy());
673 zpaEnergy = (Float_t) (esdZDC->GetZDCP2Energy());
675 zncEnergy = (Float_t) (esdZDC->GetZDCN1Energy())/8.;
676 zpcEnergy = (Float_t) (esdZDC->GetZDCP1Energy())/8.;
677 znaEnergy = (Float_t) (esdZDC->GetZDCN2Energy())/8.;
678 zpaEnergy = (Float_t) (esdZDC->GetZDCP2Energy())/8.;
680 zem1Energy = (Float_t) (esdZDC->GetZDCEMEnergy(0))/8.;
681 zem2Energy = (Float_t) (esdZDC->GetZDCEMEnergy(1))/8.;
685 else if(fAnalysisInput.CompareTo("AOD")==0){
686 //AliAODEvent *aod = dynamic_cast<AliAODEvent*> (InputEvent());
688 printf(" AOD analysis not yet implemented!!!\n\n");
692 // ***** Scaling for MC
695 v0Corr = Short_t((multV0A+multV0C) * fV0MScaleFactorMC);
697 // ***** Scaling for Data
699 v0Corr = Short_t(v0Corr / fV0MScaleFactor);
700 spdCorr = spdCorr / fSPDScaleFactor;
701 nTracks = Int_t(nTracks / fTPCScaleFactor);
704 // ***** Centrality Selection
705 if(fHtempV0M) fCentV0M = fHtempV0M->GetBinContent(fHtempV0M->FindBin((v0Corr)));
706 if(fHtempFMD) fCentFMD = fHtempFMD->GetBinContent(fHtempFMD->FindBin((multFMDA+multFMDC)));
707 if(fHtempTRK) fCentTRK = fHtempTRK->GetBinContent(fHtempTRK->FindBin(nTracks));
708 if(fHtempTKL) fCentTKL = fHtempTKL->GetBinContent(fHtempTKL->FindBin(nTracklets));
709 if(fHtempCL0) fCentCL0 = fHtempCL0->GetBinContent(fHtempCL0->FindBin(nClusters[0]));
710 if(fHtempCL1) fCentCL1 = fHtempCL1->GetBinContent(fHtempCL1->FindBin(spdCorr));
712 if(fHtempV0MvsFMD) fCentV0MvsFMD = fHtempV0MvsFMD->GetBinContent(fHtempV0MvsFMD->FindBin((multV0A+multV0C)));
713 if(fHtempTKLvsV0M) fCentTKLvsV0M = fHtempTKLvsV0M->GetBinContent(fHtempTKLvsV0M->FindBin(nTracklets));
714 if(fHtempZEMvsZDC) fCentZEMvsZDC = fHtempZEMvsZDC->GetBinContent(fHtempZEMvsZDC->FindBin(zem1Energy+zem2Energy,zncEnergy+znaEnergy+zpcEnergy+zpaEnergy));
722 if (TMath::Abs(zvtx)>fZVCut || zvtxNcont<1) fQuality += 1;
724 // ***** outliers, skip in case of MC input
727 if (IsOutlierV0MSPD(spdCorr, v0Corr, int(fCentV0M))) fQuality += 2;
729 if (IsOutlierV0MTPC(nTracks, v0Corr, int(fCentV0M))) fQuality += 4;
731 if (IsOutlierV0MZDC((zncEnergy+znaEnergy+zpcEnergy+zpaEnergy), v0Corr) &&
732 (zdcEnergyCal==kFALSE) ) fQuality += 8;
733 if (IsOutlierV0MZDCECal((zncEnergy+znaEnergy+zpcEnergy+zpaEnergy), v0Corr) &&
734 (zdcEnergyCal==kTRUE) ) fQuality += 8;
742 esdCent->SetQuality(fQuality);
743 esdCent->SetCentralityV0M(fCentV0M);
744 esdCent->SetCentralityFMD(fCentFMD);
745 esdCent->SetCentralityTRK(fCentTRK);
746 esdCent->SetCentralityTKL(fCentTKL);
747 esdCent->SetCentralityCL0(fCentCL0);
748 esdCent->SetCentralityCL1(fCentCL1);
749 esdCent->SetCentralityV0MvsFMD(fCentV0MvsFMD);
750 esdCent->SetCentralityTKLvsV0M(fCentTKLvsV0M);
751 esdCent->SetCentralityZEMvsZDC(fCentZEMvsZDC);
755 fHOutCentV0M_CVHN->Fill(fCentV0M);
757 fHOutCentV0M_CVLN->Fill(fCentV0M);
759 Bool_t isSelected = (((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected() & AliVEvent::kMB);
760 if (isSelected) { // fill the QA histograms only for MB events!
761 fHOutQuality->Fill(fQuality);
762 fHOutVertex->Fill(zvtx);
765 fHOutCentV0M->Fill(fCentV0M);
766 fHOutCentFMD->Fill(fCentFMD);
767 fHOutCentTRK->Fill(fCentTRK);
768 fHOutCentTKL->Fill(fCentTKL);
769 fHOutCentCL0->Fill(fCentCL0);
770 fHOutCentCL1->Fill(fCentCL1);
771 fHOutCentV0MvsFMD->Fill(fCentV0MvsFMD);
772 fHOutCentTKLvsV0M->Fill(fCentTKLvsV0M);
773 fHOutCentZEMvsZDC->Fill(fCentZEMvsZDC);
774 fHOutCentV0MvsCentCL1->Fill(fCentV0M,fCentCL1);
775 fHOutCentV0MvsCentTRK->Fill(fCentV0M,fCentTRK);
776 fHOutCentTRKvsCentCL1->Fill(fCentTRK,fCentCL1);
777 fHOutCentV0MvsCentZDC->Fill(fCentV0M,fCentZEMvsZDC);
778 fHOutMultV0M->Fill(multV0A+multV0C);
779 fHOutMultV0O->Fill(multV0AOnline+multV0COnline);
780 fHOutMultFMD->Fill(multFMDA+multFMDC);
781 fHOutMultTRK->Fill(nTracks);
782 fHOutMultTKL->Fill(nTracklets);
783 fHOutMultCL0->Fill(nClusters[0]);
784 fHOutMultCL1->Fill(spdCorr);
785 fHOutMultV0MvsZDN->Fill(v0Corr,(zncEnergy+znaEnergy));
786 fHOutMultZEMvsZDN->Fill((zem1Energy+zem2Energy),(zncEnergy+znaEnergy));
787 fHOutMultV0MvsZDC->Fill(v0Corr,(zncEnergy+znaEnergy+zpcEnergy+zpaEnergy));
788 fHOutMultZEMvsZDC->Fill((zem1Energy+zem2Energy),(zncEnergy+znaEnergy+zpcEnergy+zpaEnergy));
789 fHOutMultZEMvsZDCw->Fill((zem1Energy+zem2Energy),(zncEnergy+znaEnergy+zpcEnergy+zpaEnergy),fCentV0M);
790 fHOutMultV0MvsCL1->Fill(v0Corr,spdCorr);
791 fHOutMultV0MvsTRK->Fill(v0Corr,nTracks);
792 fHOutMultTRKvsCL1->Fill(nTracks,spdCorr);
793 fHOutMultV0MvsV0O->Fill(v0Corr,(multV0AOnline+multV0COnline));
794 fHOutMultV0OvsCL1->Fill((multV0AOnline+multV0COnline),spdCorr);
795 fHOutMultV0OvsTRK->Fill((multV0AOnline+multV0COnline),nTracks);
796 } else if (fQuality%2 == 0) {
797 fHOutCentV0Mqual1->Fill(fCentV0M);
798 fHOutCentTRKqual1->Fill(fCentTRK);
799 fHOutCentCL1qual1->Fill(fCentCL1);
800 fHOutMultV0MvsCL1qual1->Fill(v0Corr,spdCorr);
801 fHOutMultV0MvsTRKqual1->Fill(v0Corr,nTracks);
802 fHOutMultTRKvsCL1qual1->Fill(nTracks,spdCorr);
804 fHOutCentV0Mqual2->Fill(fCentV0M);
805 fHOutCentTRKqual2->Fill(fCentTRK);
806 fHOutCentCL1qual2->Fill(fCentCL1);
807 fHOutMultV0MvsCL1qual2->Fill(v0Corr,spdCorr);
808 fHOutMultV0MvsTRKqual2->Fill(v0Corr,nTracks);
809 fHOutMultTRKvsCL1qual2->Fill(nTracks,spdCorr);
813 PostData(1, fOutputList);
815 //________________________________________________________________________
816 void AliCentralitySelectionTask::Terminate(Option_t */*option*/)
818 // Terminate analysis
820 //________________________________________________________________________
821 Int_t AliCentralitySelectionTask::SetupRun(AliESDEvent* const esd)
823 // Setup files for run
828 // check if something to be done
829 if (fCurrentRun == esd->GetRunNumber())
832 fCurrentRun = esd->GetRunNumber();
834 TString fileName =(Form("%s/COMMON/CENTRALITY/data/centrality.root", AliAnalysisManager::GetOADBPath()));
835 AliInfo(Form("Setup Centrality Selection for run %d with file %s\n",fCurrentRun,fileName.Data()));
837 AliOADBContainer *con = new AliOADBContainer("OADB");
838 con->InitFromFile(fileName,"Centrality");
840 AliOADBCentrality* centOADB = 0;
841 centOADB = (AliOADBCentrality*)(con->GetObject(fCurrentRun));
843 AliWarning(Form("Centrality OADB does not exist for run %d, using Default \n",fCurrentRun ));
844 centOADB = (AliOADBCentrality*)(con->GetDefaultObject("oadbDefault"));
848 fUseScaling = centOADB->UseScaling();
849 fUseCleaning = centOADB->UseCleaning();
852 fZVCut = centOADB->ZVCut();
853 fOutliersCut = centOADB->OutliersCut();
856 fHtempV0M = centOADB->V0hist();
857 fHtempTRK = centOADB->TPChist();
858 fHtempCL1 = centOADB->SPDhist();
859 fHtempZEMvsZDC = centOADB->ZEMvsZDChist();
861 TString path = gSystem->ExpandPathName(fileName.Data());
862 if (!fHtempV0M) AliWarning(Form("Calibration for V0M does not exist in %s", path.Data()));
863 if (!fHtempTRK) AliWarning(Form("Calibration for TRK does not exist in %s", path.Data()));
864 if (!fHtempCL1) AliWarning(Form("Calibration for CL1 does not exist in %s", path.Data()));
865 if (!fHtempZEMvsZDC) AliWarning(Form("Calibration for ZEMvsZDC does not exist in %s", path.Data()));
868 fV0MScaleFactor = centOADB->V0MScaleFactor();
869 fSPDScaleFactor = centOADB->SPDScaleFactor();
870 fTPCScaleFactor = centOADB->TPCScaleFactor();
871 fV0MScaleFactorMC = centOADB->V0MScaleFactor();
873 // outliers parameters
874 fV0MSPDOutlierPar0 = centOADB->V0MSPDOutlierPar0();
875 fV0MSPDOutlierPar1 = centOADB->V0MSPDOutlierPar1();
876 fV0MTPCOutlierPar0 = centOADB->V0MTPCOutlierPar0();
877 fV0MTPCOutlierPar1 = centOADB->V0MTPCOutlierPar1();
879 fV0MSPDSigmaOutlierPar0 = centOADB->V0MSPDSigmaOutlierPar0();
880 fV0MSPDSigmaOutlierPar1 = centOADB->V0MSPDSigmaOutlierPar1();
881 fV0MSPDSigmaOutlierPar2 = centOADB->V0MSPDSigmaOutlierPar2();
882 fV0MTPCSigmaOutlierPar0 = centOADB->V0MTPCSigmaOutlierPar0();
883 fV0MTPCSigmaOutlierPar1 = centOADB->V0MTPCSigmaOutlierPar1();
884 fV0MTPCSigmaOutlierPar2 = centOADB->V0MTPCSigmaOutlierPar2();
886 fV0MZDCOutlierPar0 = centOADB->V0MZDCOutlierPar0();
887 fV0MZDCOutlierPar1 = centOADB->V0MZDCOutlierPar1();
888 fV0MZDCEcalOutlierPar0 = centOADB->V0MZDCEcalOutlierPar0();
889 fV0MZDCEcalOutlierPar1 = centOADB->V0MZDCEcalOutlierPar1();
898 //________________________________________________________________________
899 Bool_t AliCentralitySelectionTask::IsOutlierV0MSPD(Float_t spd, Float_t v0, Int_t cent) const
902 Float_t val = fV0MSPDOutlierPar0 + fV0MSPDOutlierPar1 * v0;
903 Float_t spdSigma = fV0MSPDSigmaOutlierPar0 + fV0MSPDSigmaOutlierPar1*cent + fV0MSPDSigmaOutlierPar2*cent*cent;
904 if ( TMath::Abs(spd-val) > fOutliersCut*spdSigma )
910 //________________________________________________________________________
911 Bool_t AliCentralitySelectionTask::IsOutlierV0MTPC(Int_t tracks, Float_t v0, Int_t cent) const
914 Float_t val = fV0MTPCOutlierPar0 + fV0MTPCOutlierPar1 * v0;
915 Float_t tpcSigma = fV0MTPCSigmaOutlierPar0 + fV0MTPCSigmaOutlierPar1*cent + fV0MTPCSigmaOutlierPar2*cent*cent;
916 if ( TMath::Abs(tracks-val) > fOutliersCut*tpcSigma )
922 //________________________________________________________________________
923 Bool_t AliCentralitySelectionTask::IsOutlierV0MZDC(Float_t zdc, Float_t v0) const
926 Float_t val = fV0MZDCOutlierPar0 + fV0MZDCOutlierPar1 * v0;
933 //________________________________________________________________________
934 Bool_t AliCentralitySelectionTask::IsOutlierV0MZDCECal(Float_t zdc, Float_t v0) const
937 Float_t val = fV0MZDCEcalOutlierPar0 + fV0MZDCEcalOutlierPar1 * v0;