]> git.uio.no Git - u/mrichter/AliRoot.git/blame - ANALYSIS/AliCentralitySelectionTask.cxx
Removing obsolete class AliCaloRawAnalyzerLMS
[u/mrichter/AliRoot.git] / ANALYSIS / AliCentralitySelectionTask.cxx
CommitLineData
a540a9d3 1/**************************************************************************
2 * Copyright(c) 1998-2008, 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 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 **************************************************************************/
15
16//*****************************************************
17// Class AliCentralitySelectionTask
18// Class to analyze determine centrality
19// author: Alberica Toia
20//*****************************************************
21
22#include "AliCentralitySelectionTask.h"
23
24#include <TTree.h>
25#include <TList.h>
26#include <TH1F.h>
27#include <TH2F.h>
28#include <TF1.h>
29#include <TProfile.h>
30#include <TFile.h>
31#include <TObjString.h>
32#include <TString.h>
33#include <TCanvas.h>
34#include <TROOT.h>
35#include <TDirectory.h>
36#include <TSystem.h>
37#include <iostream>
38
39#include "AliAnalysisManager.h"
40#include "AliVEvent.h"
41#include "AliESD.h"
42#include "AliESDEvent.h"
43#include "AliESDHeader.h"
44#include "AliESDInputHandler.h"
45#include "AliESDZDC.h"
46#include "AliESDFMD.h"
47#include "AliESDVZERO.h"
48#include "AliCentrality.h"
49#include "AliESDtrackCuts.h"
50#include "AliMultiplicity.h"
51#include "AliAODHandler.h"
52#include "AliAODEvent.h"
53#include "AliESDVertex.h"
54#include "AliAODVertex.h"
55#include "AliAODMCHeader.h"
56#include "AliMCEvent.h"
57#include "AliMCEventHandler.h"
58#include "AliMCParticle.h"
59#include "AliStack.h"
60#include "AliHeader.h"
61#include "AliAODMCParticle.h"
62#include "AliAnalysisTaskSE.h"
63#include "AliGenEventHeader.h"
64#include "AliGenHijingEventHeader.h"
65#include "AliPhysicsSelectionTask.h"
66#include "AliPhysicsSelection.h"
67#include "AliBackgroundSelection.h"
68#include "AliESDUtils.h"
69
70ClassImp(AliCentralitySelectionTask)
71
72
73//________________________________________________________________________
74AliCentralitySelectionTask::AliCentralitySelectionTask():
75AliAnalysisTaskSE(),
a540a9d3 76 fAnalysisInput("ESD"),
77 fIsMCInput(kFALSE),
a121a512 78 fPass(0),
a540a9d3 79 fFile(0),
80 fFile2(0),
81 fCurrentRun(-1),
82 fRunNo(-1),
83 fLowRunN(0),
84 fHighRunN(0),
85 fUseScaling(0),
d015c169 86 fUseCleaning(0),
a540a9d3 87 fTrackCuts(0),
88 fZVCut(10),
89 fOutliersCut(5),
90 fQuality(0),
91 fCentV0M(0),
92 fCentFMD(0),
93 fCentTRK(0),
94 fCentTKL(0),
95 fCentCL0(0),
96 fCentCL1(0),
97 fCentV0MvsFMD(0),
98 fCentTKLvsV0M(0),
99 fCentZEMvsZDC(0),
100 fHtempV0M(0),
101 fHtempFMD(0),
102 fHtempTRK(0),
103 fHtempTKL(0),
104 fHtempCL0(0),
105 fHtempCL1(0),
106 fHtempV0MvsFMD(0),
107 fHtempTKLvsV0M(0),
108 fHtempZEMvsZDC(0),
109 fOutputList(0),
110 fHOutCentV0M (0),
111 fHOutCentFMD (0),
112 fHOutCentTRK (0),
113 fHOutCentTKL (0),
114 fHOutCentCL0 (0),
115 fHOutCentCL1 (0),
116 fHOutCentV0MvsFMD(0),
117 fHOutCentTKLvsV0M(0),
118 fHOutCentZEMvsZDC(0),
9138a60c 119 fHOutCentV0MvsCentCL1(0),
120 fHOutCentV0MvsCentTRK(0),
121 fHOutCentTRKvsCentCL1(0),
a540a9d3 122 fHOutMultV0M(0),
123 fHOutMultV0R(0),
124 fHOutMultFMD(0),
125 fHOutMultTRK(0),
126 fHOutMultTKL(0),
127 fHOutMultCL0(0),
128 fHOutMultCL1(0),
31200bdf 129 fHOutMultV0MvsZDN(0),
130 fHOutMultZEMvsZDN(0),
a540a9d3 131 fHOutMultV0MvsZDC(0),
132 fHOutMultZEMvsZDC(0),
133 fHOutMultV0MvsCL1(0),
134 fHOutMultV0MvsTRK(0),
135 fHOutMultTRKvsCL1(0),
04341115 136 fHOutCentV0Mqual1(0),
137 fHOutCentTRKqual1(0),
138 fHOutCentCL1qual1(0),
139 fHOutCentV0Mqual2(0),
140 fHOutCentTRKqual2(0),
141 fHOutCentCL1qual2(0),
a540a9d3 142 fHOutQuality(0),
143 fHOutVertex(0)
144{
145 // Default constructor
146 AliInfo("Centrality Selection enabled.");
609f601e 147 fLowRunN =136851;
148 fHighRunN=139517;
149
ff325947 150 for (Int_t i=0; i < 2667; i++) {
04341115 151 fV0MScaleFactor[i]=0.0;
152 fSPDScaleFactor[i]=0.0;
153 fTPCScaleFactor[i]=0.0;
154 fV0MScaleFactorMC[i]=0.0;
a540a9d3 155 }
156 fUseScaling=kTRUE;
d015c169 157 fUseCleaning=kTRUE;
cb15f505 158 fBranchNames="ESD:AliESDRun.,AliESDHeader.,AliESDZDC.,AliESDFMD.,AliESDVZERO."
159 ",SPDVertex.,TPCVertex.,PrimaryVertex.,AliMultiplicity.,Tracks ";
a540a9d3 160}
161
162//________________________________________________________________________
163AliCentralitySelectionTask::AliCentralitySelectionTask(const char *name):
164 AliAnalysisTaskSE(name),
a540a9d3 165 fAnalysisInput("ESD"),
166 fIsMCInput(kFALSE),
a121a512 167 fPass(0),
a540a9d3 168 fFile(0),
169 fFile2(0),
170 fCurrentRun(-1),
171 fRunNo(-1),
172 fLowRunN(0),
173 fHighRunN(0),
174 fUseScaling(0),
d015c169 175 fUseCleaning(0),
a540a9d3 176 fTrackCuts(0),
177 fZVCut(10),
178 fOutliersCut(5),
179 fQuality(0),
180 fCentV0M(0),
181 fCentFMD(0),
182 fCentTRK(0),
183 fCentTKL(0),
184 fCentCL0(0),
185 fCentCL1(0),
186 fCentV0MvsFMD(0),
187 fCentTKLvsV0M(0),
188 fCentZEMvsZDC(0),
189 fHtempV0M(0),
190 fHtempFMD(0),
191 fHtempTRK(0),
192 fHtempTKL(0),
193 fHtempCL0(0),
194 fHtempCL1(0),
195 fHtempV0MvsFMD(0),
196 fHtempTKLvsV0M(0),
197 fHtempZEMvsZDC(0),
198 fOutputList(0),
199 fHOutCentV0M (0),
200 fHOutCentFMD (0),
201 fHOutCentTRK (0),
202 fHOutCentTKL (0),
203 fHOutCentCL0 (0),
204 fHOutCentCL1 (0),
205 fHOutCentV0MvsFMD(0),
206 fHOutCentTKLvsV0M(0),
207 fHOutCentZEMvsZDC(0),
9138a60c 208 fHOutCentV0MvsCentCL1(0),
209 fHOutCentV0MvsCentTRK(0),
210 fHOutCentTRKvsCentCL1(0),
a540a9d3 211 fHOutMultV0M(0),
212 fHOutMultV0R(0),
213 fHOutMultFMD(0),
214 fHOutMultTRK(0),
215 fHOutMultTKL(0),
216 fHOutMultCL0(0),
217 fHOutMultCL1(0),
31200bdf 218 fHOutMultV0MvsZDN(0),
219 fHOutMultZEMvsZDN(0),
a540a9d3 220 fHOutMultV0MvsZDC(0),
221 fHOutMultZEMvsZDC(0),
222 fHOutMultV0MvsCL1(0),
223 fHOutMultV0MvsTRK(0),
224 fHOutMultTRKvsCL1(0),
04341115 225 fHOutCentV0Mqual1(0),
226 fHOutCentTRKqual1(0),
227 fHOutCentCL1qual1(0),
228 fHOutCentV0Mqual2(0),
229 fHOutCentTRKqual2(0),
230 fHOutCentCL1qual2(0),
a540a9d3 231 fHOutQuality(0),
232 fHOutVertex(0)
233{
609f601e 234 // Default constructor
235 fLowRunN =136851;
236 fHighRunN=139517;
237
a540a9d3 238 AliInfo("Centrality Selection enabled.");
239 DefineOutput(1, TList::Class());
ff325947 240 for (Int_t i=0; i<2667; i++) {
04341115 241 fV0MScaleFactor[i]=0.0;
242 fSPDScaleFactor[i]=0.0;
243 fTPCScaleFactor[i]=0.0;
244 fV0MScaleFactorMC[i]=0.0;
a540a9d3 245 }
246 fUseScaling=kTRUE;
d015c169 247 fUseCleaning=kTRUE;
cb15f505 248 fBranchNames="ESD:AliESDRun.,AliESDHeader.,AliESDZDC.,AliESDFMD.,AliESDVZERO."
249 ",SPDVertex.,TPCVertex.,PrimaryVertex.,AliMultiplicity.,Tracks ";
a540a9d3 250}
251
252//________________________________________________________________________
253AliCentralitySelectionTask& AliCentralitySelectionTask::operator=(const AliCentralitySelectionTask& c)
254{
255 // Assignment operator
256 if (this!=&c) {
257 AliAnalysisTaskSE::operator=(c);
258 }
259 return *this;
260}
261
262//________________________________________________________________________
263AliCentralitySelectionTask::AliCentralitySelectionTask(const AliCentralitySelectionTask& ana):
264 AliAnalysisTaskSE(ana),
a540a9d3 265 fAnalysisInput(ana.fDebug),
266 fIsMCInput(ana.fIsMCInput),
a121a512 267 fPass(ana.fPass),
a540a9d3 268 fFile(ana.fFile),
269 fFile2(ana.fFile2),
270 fCurrentRun(ana.fCurrentRun),
271 fRunNo(ana.fRunNo),
272 fLowRunN(ana.fLowRunN),
273 fHighRunN(ana.fHighRunN),
274 fUseScaling(ana.fUseScaling),
d015c169 275 fUseCleaning(ana.fUseCleaning),
a540a9d3 276 fTrackCuts(ana.fTrackCuts),
277 fZVCut(ana.fZVCut),
278 fOutliersCut(ana.fOutliersCut),
279 fQuality(ana.fQuality),
280 fCentV0M(ana.fCentV0M),
281 fCentFMD(ana.fCentFMD),
282 fCentTRK(ana.fCentTRK),
283 fCentTKL(ana.fCentTKL),
284 fCentCL0(ana.fCentCL0),
285 fCentCL1(ana.fCentCL1),
286 fCentV0MvsFMD(ana.fCentV0MvsFMD),
287 fCentTKLvsV0M(ana.fCentTKLvsV0M),
288 fCentZEMvsZDC(ana.fCentZEMvsZDC),
289 fHtempV0M(ana.fHtempV0M),
290 fHtempFMD(ana.fHtempFMD),
291 fHtempTRK(ana.fHtempTRK),
292 fHtempTKL(ana.fHtempTKL),
293 fHtempCL0(ana.fHtempCL0),
294 fHtempCL1(ana.fHtempCL1),
295 fHtempV0MvsFMD(ana.fHtempV0MvsFMD),
296 fHtempTKLvsV0M(ana.fHtempTKLvsV0M),
297 fHtempZEMvsZDC(ana.fHtempZEMvsZDC),
298 fOutputList(ana.fOutputList),
299 fHOutCentV0M (ana.fHOutCentV0M ),
300 fHOutCentFMD (ana.fHOutCentFMD ),
301 fHOutCentTRK (ana.fHOutCentTRK ),
302 fHOutCentTKL (ana.fHOutCentTKL ),
303 fHOutCentCL0 (ana.fHOutCentCL0 ),
304 fHOutCentCL1 (ana.fHOutCentCL1 ),
305 fHOutCentV0MvsFMD(ana.fHOutCentV0MvsFMD),
306 fHOutCentTKLvsV0M(ana.fHOutCentTKLvsV0M),
307 fHOutCentZEMvsZDC(ana.fHOutCentZEMvsZDC),
9138a60c 308 fHOutCentV0MvsCentCL1(ana.fHOutCentV0MvsCentCL1),
309 fHOutCentV0MvsCentTRK(ana.fHOutCentV0MvsCentTRK),
310 fHOutCentTRKvsCentCL1(ana.fHOutCentTRKvsCentCL1),
a540a9d3 311 fHOutMultV0M(ana.fHOutMultV0M),
312 fHOutMultV0R(ana.fHOutMultV0R),
313 fHOutMultFMD(ana.fHOutMultFMD),
314 fHOutMultTRK(ana.fHOutMultTRK),
315 fHOutMultTKL(ana.fHOutMultTKL),
316 fHOutMultCL0(ana.fHOutMultCL0),
317 fHOutMultCL1(ana.fHOutMultCL1),
31200bdf 318 fHOutMultV0MvsZDN(ana.fHOutMultV0MvsZDN),
319 fHOutMultZEMvsZDN(ana.fHOutMultZEMvsZDN),
a540a9d3 320 fHOutMultV0MvsZDC(ana.fHOutMultV0MvsZDC),
321 fHOutMultZEMvsZDC(ana.fHOutMultZEMvsZDC),
322 fHOutMultV0MvsCL1(ana.fHOutMultV0MvsCL1),
323 fHOutMultV0MvsTRK(ana.fHOutMultV0MvsTRK),
324 fHOutMultTRKvsCL1(ana.fHOutMultTRKvsCL1),
04341115 325 fHOutCentV0Mqual1(ana.fHOutCentV0Mqual1),
326 fHOutCentTRKqual1(ana.fHOutCentTRKqual1),
327 fHOutCentCL1qual1(ana.fHOutCentCL1qual1),
328 fHOutCentV0Mqual2(ana.fHOutCentV0Mqual2),
329 fHOutCentTRKqual2(ana.fHOutCentTRKqual2),
330 fHOutCentCL1qual2(ana.fHOutCentCL1qual2),
a540a9d3 331 fHOutQuality(ana.fHOutQuality),
332 fHOutVertex(ana.fHOutVertex)
333{
334 // Copy Constructor
ff325947 335 for (Int_t i=0; i<2667; i++) {
04341115 336 fV0MScaleFactor[i]=0.0;
337 fSPDScaleFactor[i]=0.0;
338 fTPCScaleFactor[i]=0.0;
339 fV0MScaleFactorMC[i]=0.0;
7ee7a2df 340 }
341
a540a9d3 342}
343
344//________________________________________________________________________
345AliCentralitySelectionTask::~AliCentralitySelectionTask()
346{
347 // Destructor
348 if (fOutputList && !AliAnalysisManager::GetAnalysisManager()->IsProofMode()) delete fOutputList;
349 if (fTrackCuts) delete fTrackCuts;
350}
351
352//________________________________________________________________________
353void AliCentralitySelectionTask::UserCreateOutputObjects()
354{
355 // Create the output containers
356 if(fDebug>1) printf("AnalysisCentralitySelectionTask::UserCreateOutputObjects() \n");
357 AliLog::SetClassDebugLevel("AliCentralitySelectionTask", AliLog::kInfo);
358
359 fOutputList = new TList();
360 fOutputList->SetOwner();
bab72d34 361 fHOutCentV0M = new TH1F("fHOutCentV0M","fHOutCentV0M; Centrality V0",505,0,101);
362 fHOutCentFMD = new TH1F("fHOutCentFMD","fHOutCentFMD; Centrality FMD",505,0,101);
363 fHOutCentTRK = new TH1F("fHOutCentTRK","fHOutCentTRK; Centrality TPC",505,0,101);
364 fHOutCentTKL = new TH1F("fHOutCentTKL","fHOutCentTKL; Centrality tracklets",505,0,101);
365 fHOutCentCL0 = new TH1F("fHOutCentCL0","fHOutCentCL0; Centrality SPD inner",505,0,101);
366 fHOutCentCL1 = new TH1F("fHOutCentCL1","fHOutCentCL1; Centrality SPD outer",505,0,101);
367 fHOutCentV0MvsFMD= new TH1F("fHOutCentV0MvsFMD","fHOutCentV0MvsFMD; Centrality V0 vs FMD",505,0,101);
368 fHOutCentTKLvsV0M= new TH1F("fHOutCentTKLvsV0M","fHOutCentTKLvsV0M; Centrality tracklets vs V0",505,0,101);
369 fHOutCentZEMvsZDC= new TH1F("fHOutCentZEMvsZDC","fHOutCentZEMvsZDC; Centrality ZEM vs ZDC",505,0,101);
370 fHOutCentV0MvsCentCL1= new TH2F("fHOutCentV0MvsCentCL1","fHOutCentV0MvsCentCL1; Cent V0 vs Cent SPD",505,0,101,505,0,101);
371 fHOutCentV0MvsCentTRK= new TH2F("fHOutCentV0MvsCentTRK","fHOutCentV0MvsCentTRK; Cent V0 vs Cent TPC",505,0,101,505,0,101);
372 fHOutCentTRKvsCentCL1= new TH2F("fHOutCentTRKvsCentCL1","fHOutCentTRKvsCentCL1; Cent TPC vs Cent SPD",505,0,101,505,0,101);
a540a9d3 373
374 fHOutMultV0M = new TH1F("fHOutMultV0M","fHOutMultV0M; Multiplicity V0",25000,0,25000);
2cf88888 375 fHOutMultV0R = new TH1F("fHOutMultV0R","fHOutMultV0R; Multiplicity V0",30000,0,30000);
a540a9d3 376 fHOutMultFMD = new TH1F("fHOutMultFMD","fHOutMultFMD; Multiplicity FMD",24000,0,24000);
377 fHOutMultTRK = new TH1F("fHOutMultTRK","fHOutMultTRK; Multiplicity TPC",4000,0,4000);
378 fHOutMultTKL = new TH1F("fHOutMultTKL","fHOutMultTKL; Multiplicity tracklets",5000,0,5000);
379 fHOutMultCL0 = new TH1F("fHOutMultCL0","fHOutMultCL0; Multiplicity SPD inner",7000,0,7000);
380 fHOutMultCL1 = new TH1F("fHOutMultCL1","fHOutMultCL1; Multiplicity SPD outer",7000,0,7000);
31200bdf 381 fHOutMultV0MvsZDN = new TH2F("fHOutMultV0MvsZDN","fHOutMultV0MvsZDN; Multiplicity V0; Energy ZDC-N",500,0,25000,500,0,180000);
382 fHOutMultZEMvsZDN = new TH2F("fHOutMultZEMvsZDN","fHOutMultZEMvsZDN; Energy ZEM; Energy ZDC-N",500,0,2500,500,0,180000);
383 fHOutMultV0MvsZDC = new TH2F("fHOutMultV0MvsZDC","fHOutMultV0MvsZDC; Multiplicity V0; Energy ZDC",500,0,25000,500,0,200000);
384 fHOutMultZEMvsZDC = new TH2F("fHOutMultZEMvsZDC","fHOutMultZEMvsZDC; Energy ZEM; Energy ZDC",500,0,2500,500,0,200000);
a540a9d3 385 fHOutMultV0MvsCL1 = new TH2F("fHOutMultV0MvsCL1","fHOutMultV0MvsCL1; Multiplicity V0; Multiplicity SPD outer",2500,0,25000,700,0,7000);
386 fHOutMultV0MvsTRK = new TH2F("fHOutMultV0MvsTRK","fHOutMultV0MvsTRK; Multiplicity V0; Multiplicity TPC",2500,0,25000,400,0,4000);
387 fHOutMultTRKvsCL1 = new TH2F("fHOutMultTRKvsCL1","fHOutMultTRKvsCL1; Multiplicity TPC; Multiplicity SPD outer",400,0,4000,700,0,7000);
388
04341115 389 fHOutCentV0Mqual1 = new TH1F("fHOutCentV0M_qual1","fHOutCentV0M_qual1; Centrality V0",505,0,101);
390 fHOutCentTRKqual1 = new TH1F("fHOutCentTRK_qual1","fHOutCentTRK_qual1; Centrality TPC",505,0,101);
391 fHOutCentCL1qual1 = new TH1F("fHOutCentCL1_qual1","fHOutCentCL1_qual1; Centrality SPD outer",505,0,101);
a540a9d3 392
04341115 393 fHOutCentV0Mqual2 = new TH1F("fHOutCentV0M_qual2","fHOutCentV0M_qual2; Centrality V0",505,0,101);
394 fHOutCentTRKqual2 = new TH1F("fHOutCentTRK_qual2","fHOutCentTRK_qual2; Centrality TPC",505,0,101);
395 fHOutCentCL1qual2 = new TH1F("fHOutCentCL1_qual2","fHOutCentCL1_qual2; Centrality SPD outer",505,0,101);
a540a9d3 396
397 fHOutQuality = new TH1F("fHOutQuality", "fHOutQuality", 10,-0.5,9.5);
398 fHOutVertex = new TH1F("fHOutVertex", "fHOutVertex", 100,-20,20);
399
400 fOutputList->Add( fHOutCentV0M );
401 fOutputList->Add( fHOutCentFMD );
402 fOutputList->Add( fHOutCentTRK );
403 fOutputList->Add( fHOutCentTKL );
404 fOutputList->Add( fHOutCentCL0 );
405 fOutputList->Add( fHOutCentCL1 );
406 fOutputList->Add( fHOutCentV0MvsFMD);
407 fOutputList->Add( fHOutCentTKLvsV0M);
408 fOutputList->Add( fHOutCentZEMvsZDC);
9138a60c 409 fOutputList->Add( fHOutCentV0MvsCentCL1);
410 fOutputList->Add( fHOutCentV0MvsCentTRK);
411 fOutputList->Add( fHOutCentTRKvsCentCL1);
a540a9d3 412 fOutputList->Add( fHOutMultV0M);
413 fOutputList->Add( fHOutMultV0R);
414 fOutputList->Add( fHOutMultFMD);
415 fOutputList->Add( fHOutMultTRK);
416 fOutputList->Add( fHOutMultTKL);
417 fOutputList->Add( fHOutMultCL0);
418 fOutputList->Add( fHOutMultCL1);
31200bdf 419 fOutputList->Add( fHOutMultV0MvsZDN);
420 fOutputList->Add( fHOutMultZEMvsZDN);
a540a9d3 421 fOutputList->Add( fHOutMultV0MvsZDC);
422 fOutputList->Add( fHOutMultZEMvsZDC);
423 fOutputList->Add( fHOutMultV0MvsCL1);
424 fOutputList->Add( fHOutMultV0MvsTRK);
425 fOutputList->Add( fHOutMultTRKvsCL1);
04341115 426 fOutputList->Add( fHOutCentV0Mqual1 );
427 fOutputList->Add( fHOutCentTRKqual1 );
428 fOutputList->Add( fHOutCentCL1qual1 );
429 fOutputList->Add( fHOutCentV0Mqual2 );
430 fOutputList->Add( fHOutCentTRKqual2 );
431 fOutputList->Add( fHOutCentCL1qual2 );
a540a9d3 432 fOutputList->Add( fHOutQuality );
433 fOutputList->Add( fHOutVertex );
434
435
436 fTrackCuts = AliESDtrackCuts::GetStandardTPCOnlyTrackCuts();
437
438 PostData(1, fOutputList);
439
440 MyInitScaleFactor();
441 if (fIsMCInput) MyInitScaleFactorMC();
442}
443
444//________________________________________________________________________
445void AliCentralitySelectionTask::UserExec(Option_t */*option*/)
446{
447 // Execute analysis for current event:
448 if(fDebug>1) printf(" **** AliCentralitySelectionTask::UserExec() \n");
449
450 Float_t zncEnergy = 0.; // ZNC Energy
451 Float_t zpcEnergy = 0.; // ZPC Energy
452 Float_t znaEnergy = 0.; // ZNA Energy
453 Float_t zpaEnergy = 0.; // ZPA Energy
454 Float_t zem1Energy = 0.; // ZEM1 Energy
455 Float_t zem2Energy = 0.; // ZEM2 Energy
cfead5a9 456 Bool_t zdcEnergyCal = kFALSE; // if zdc is calibrated (in pass2)
457
a540a9d3 458 Int_t nTracks = 0; // no. tracks
459 Int_t nTracklets = 0; // no. tracklets
460 Int_t nClusters[6] = {0}; // no. clusters on 6 ITS layers
461 Int_t nChips[2]; // no. chips on 2 SPD layers
462 Float_t spdCorr =0; // corrected spd2 multiplicity
463
464 Float_t multV0A = 0; // multiplicity from V0 reco side A
465 Float_t multV0C = 0; // multiplicity from V0 reco side C
466 Float_t multFMDA = 0; // multiplicity from FMD on detector A
467 Float_t multFMDC = 0; // multiplicity from FMD on detector C
468
469 Short_t v0Corr = 0; // corrected V0 multiplicity
470 Short_t v0CorrResc = 0; // corrected and rescaled V0 multiplicity
471
472 Float_t zvtx =0; // z-vertex SPD
473
474 AliCentrality *esdCent = 0;
475
476 if(fAnalysisInput.CompareTo("ESD")==0){
477
478 AliVEvent* event = InputEvent();
479 AliESDEvent* esd = dynamic_cast<AliESDEvent*>(event);
480 if (!esd) {
481 AliError("No ESD Event");
482 return;
483 }
cb15f505 484
485 LoadBranches();
a540a9d3 486
487 if (fRunNo<=0) {
488 if (SetupRun(esd)<0)
489 AliFatal("Centrality File not available for this run");
490 }
491
492 esdCent = esd->GetCentrality();
493
494 // ***** V0 info
495 AliESDVZERO* esdV0 = esd->GetVZEROData();
496 multV0A=esdV0->GetMTotV0A();
497 multV0C=esdV0->GetMTotV0C();
498
499 Float_t v0CorrR;
500 v0Corr = (Short_t)AliESDUtils::GetCorrV0(esd,v0CorrR);
501 v0CorrResc = (Short_t)v0CorrR;
502
503 // ***** Vertex Info
504 const AliESDVertex* vtxESD = esd->GetPrimaryVertexSPD();
505 zvtx = vtxESD->GetZ();
506
507 // ***** CB info (tracklets, clusters, chips)
508 //nTracks = event->GetNumberOfTracks();
509 nTracks = fTrackCuts ? (Short_t)fTrackCuts->GetReferenceMultiplicity(esd,kTRUE):-1;
510
511 const AliMultiplicity *mult = esd->GetMultiplicity();
512
513 nTracklets = mult->GetNumberOfTracklets();
514
515 for(Int_t ilay=0; ilay<6; ilay++){
516 nClusters[ilay] = mult->GetNumberOfITSClusters(ilay);
517 }
518
519 for(Int_t ilay=0; ilay<2; ilay++){
520 nChips[ilay] = mult->GetNumberOfFiredChips(ilay);
521 }
522
523 spdCorr = AliESDUtils::GetCorrSPD2(nClusters[1],zvtx);
524
525 // ***** FMD info
526 AliESDFMD *fmd = esd->GetFMDData();
527 Float_t totalMultA = 0;
528 Float_t totalMultC = 0;
529 const Float_t fFMDLowCut = 0.4;
530
531 for(UShort_t det=1;det<=3;det++) {
532 Int_t nRings = (det==1 ? 1 : 2);
533 for (UShort_t ir = 0; ir < nRings; ir++) {
534 Char_t ring = (ir == 0 ? 'I' : 'O');
535 UShort_t nsec = (ir == 0 ? 20 : 40);
536 UShort_t nstr = (ir == 0 ? 512 : 256);
537 for(UShort_t sec =0; sec < nsec; sec++) {
538 for(UShort_t strip = 0; strip < nstr; strip++) {
539
04341115 540 Float_t fmdMult = fmd->Multiplicity(det,ring,sec,strip);
541 if(fmdMult == 0 || fmdMult == AliESDFMD::kInvalidMult) continue;
a540a9d3 542
543 Float_t nParticles=0;
544
04341115 545 if(fmdMult > fFMDLowCut) {
a540a9d3 546 nParticles = 1.;
547 }
548
549 if (det<3) totalMultA = totalMultA + nParticles;
550 else totalMultC = totalMultC + nParticles;
551
552 }
553 }
554 }
555 }
556 multFMDA = totalMultA;
557 multFMDC = totalMultC;
558
559 // ***** ZDC info
560 AliESDZDC *esdZDC = esd->GetESDZDC();
cfead5a9 561 zdcEnergyCal = esdZDC->AliESDZDC::TestBit(AliESDZDC::kEnergyCalibratedSignal);
562 if (zdcEnergyCal) {
563 zncEnergy = (Float_t) (esdZDC->GetZDCN1Energy());
564 zpcEnergy = (Float_t) (esdZDC->GetZDCP1Energy());
565 znaEnergy = (Float_t) (esdZDC->GetZDCN2Energy());
566 zpaEnergy = (Float_t) (esdZDC->GetZDCP2Energy());
567 } else {
568 zncEnergy = (Float_t) (esdZDC->GetZDCN1Energy())/8.;
569 zpcEnergy = (Float_t) (esdZDC->GetZDCP1Energy())/8.;
570 znaEnergy = (Float_t) (esdZDC->GetZDCN2Energy())/8.;
571 zpaEnergy = (Float_t) (esdZDC->GetZDCP2Energy())/8.;
572 }
a540a9d3 573 zem1Energy = (Float_t) (esdZDC->GetZDCEMEnergy(0))/8.;
574 zem2Energy = (Float_t) (esdZDC->GetZDCEMEnergy(1))/8.;
cfead5a9 575
a540a9d3 576 }
577 else if(fAnalysisInput.CompareTo("AOD")==0){
578 //AliAODEvent *aod = dynamic_cast<AliAODEvent*> (InputEvent());
579 // to be implemented
580 printf(" AOD analysis not yet implemented!!!\n\n");
581 return;
582 }
583
584
585 // ***** Scaling
a121a512 586 // ***** Scaling for pass2
587 if (fPass==2) {
588 fUseScaling=kFALSE;
589 fUseCleaning=kFALSE;
590 }
d015c169 591 // ***** Scaling for MC
592 if (fIsMCInput) {
593 fUseScaling=kFALSE;
04341115 594 Float_t tempScalefactorV0M = MyGetScaleFactorMC(fCurrentRun);
595 v0Corr = Short_t((multV0A+multV0C) * tempScalefactorV0M);
d015c169 596 }
597 // ***** Scaling for Data
a540a9d3 598 if (fUseScaling) {
04341115 599 Float_t tempScalefactorV0M = MyGetScaleFactor(fCurrentRun,0);
600 Float_t tempScalefactorSPD = MyGetScaleFactor(fCurrentRun,1);
601 Float_t tempScalefactorTPC = MyGetScaleFactor(fCurrentRun,2);
602 v0Corr = Short_t(v0Corr / tempScalefactorV0M);
603 spdCorr = spdCorr / tempScalefactorSPD;
604 nTracks = Int_t(nTracks / tempScalefactorTPC);
a540a9d3 605 }
a540a9d3 606
607 // ***** Centrality Selection
608 if(fHtempV0M) fCentV0M = fHtempV0M->GetBinContent(fHtempV0M->FindBin((v0Corr)));
609 if(fHtempFMD) fCentFMD = fHtempFMD->GetBinContent(fHtempFMD->FindBin((multFMDA+multFMDC)));
610 if(fHtempTRK) fCentTRK = fHtempTRK->GetBinContent(fHtempTRK->FindBin(nTracks));
611 if(fHtempTKL) fCentTKL = fHtempTKL->GetBinContent(fHtempTKL->FindBin(nTracklets));
612 if(fHtempCL0) fCentCL0 = fHtempCL0->GetBinContent(fHtempCL0->FindBin(nClusters[0]));
613 if(fHtempCL1) fCentCL1 = fHtempCL1->GetBinContent(fHtempCL1->FindBin(spdCorr));
614
615 if(fHtempV0MvsFMD) fCentV0MvsFMD = fHtempV0MvsFMD->GetBinContent(fHtempV0MvsFMD->FindBin((multV0A+multV0C)));
616 if(fHtempTKLvsV0M) fCentTKLvsV0M = fHtempTKLvsV0M->GetBinContent(fHtempTKLvsV0M->FindBin(nTracklets));
7ee7a2df 617 if(fHtempZEMvsZDC) fCentZEMvsZDC = fHtempZEMvsZDC->GetBinContent(fHtempZEMvsZDC->FindBin(zem1Energy+zem2Energy,zncEnergy+znaEnergy+zpcEnergy+zpaEnergy));
a540a9d3 618
619 // ***** Cleaning
d015c169 620 if (fUseCleaning) {
621 fQuality=0;
622 fZVCut=10;
623 fOutliersCut=6;
624
625 // ***** vertex
626 if (TMath::Abs(zvtx)>fZVCut) fQuality += 1;
627
628 // ***** outliers
629 // **** V0 vs SPD
630 if (IsOutlierV0MSPD(spdCorr, v0Corr, int(fCentV0M))) fQuality += 2;
631 // ***** V0 vs TPC
632 if (IsOutlierV0MTPC(nTracks, v0Corr, int(fCentV0M))) fQuality += 4;
633 // ***** V0 vs ZDC
634 if (IsOutlierV0MZDC((zncEnergy+znaEnergy+zpcEnergy+zpaEnergy), v0Corr) &&
635 (zdcEnergyCal==kFALSE) && !(fIsMCInput)) fQuality += 8;
636 if (IsOutlierV0MZDCECal((zncEnergy+znaEnergy+zpcEnergy+zpaEnergy), v0Corr) &&
637 ((zdcEnergyCal==kTRUE) || (fIsMCInput))) fQuality += 8;
638 } else {
639 fQuality = 0;
640 }
a540a9d3 641
d015c169 642
a540a9d3 643 if (esdCent) {
644 esdCent->SetQuality(fQuality);
645 esdCent->SetCentralityV0M(fCentV0M);
646 esdCent->SetCentralityFMD(fCentFMD);
647 esdCent->SetCentralityTRK(fCentTRK);
648 esdCent->SetCentralityTKL(fCentTKL);
649 esdCent->SetCentralityCL0(fCentCL0);
650 esdCent->SetCentralityCL1(fCentCL1);
651 esdCent->SetCentralityV0MvsFMD(fCentV0MvsFMD);
652 esdCent->SetCentralityTKLvsV0M(fCentTKLvsV0M);
653 esdCent->SetCentralityZEMvsZDC(fCentZEMvsZDC);
654 }
655
656 fHOutQuality->Fill(fQuality);
657 fHOutVertex->Fill(zvtx);
658
31200bdf 659 fHOutMultV0MvsZDN->Fill(v0Corr,(zncEnergy+znaEnergy));
660 fHOutMultZEMvsZDN->Fill((zem1Energy+zem2Energy),(zncEnergy+znaEnergy));
a540a9d3 661 fHOutMultV0MvsZDC->Fill(v0Corr,(zncEnergy+znaEnergy+zpcEnergy+zpaEnergy));
662 fHOutMultZEMvsZDC->Fill((zem1Energy+zem2Energy),(zncEnergy+znaEnergy+zpcEnergy+zpaEnergy));
663 fHOutMultV0MvsCL1->Fill(v0Corr,spdCorr);
664 fHOutMultV0MvsTRK->Fill(v0Corr,nTracks);
665
666 if (fQuality==0) {
667 fHOutCentV0M->Fill(fCentV0M);
668 fHOutCentFMD->Fill(fCentFMD);
669 fHOutCentTRK->Fill(fCentTRK);
670 fHOutCentTKL->Fill(fCentTKL);
671 fHOutCentCL0->Fill(fCentCL0);
672 fHOutCentCL1->Fill(fCentCL1);
673 fHOutCentV0MvsFMD->Fill(fCentV0MvsFMD);
674 fHOutCentTKLvsV0M->Fill(fCentTKLvsV0M);
675 fHOutCentZEMvsZDC->Fill(fCentZEMvsZDC);
9138a60c 676 fHOutCentV0MvsCentCL1->Fill(fCentV0M,fCentCL1);
677 fHOutCentV0MvsCentTRK->Fill(fCentV0M,fCentTRK);
678 fHOutCentTRKvsCentCL1->Fill(fCentTRK,fCentCL1);
a540a9d3 679 fHOutMultV0M->Fill(v0Corr);
680 fHOutMultV0R->Fill(multV0A+multV0C);
681 fHOutMultFMD->Fill((multFMDA+multFMDC));
682 fHOutMultTRK->Fill(nTracks);
683 fHOutMultTKL->Fill(nTracklets);
684 fHOutMultCL0->Fill(nClusters[0]);
685 fHOutMultCL1->Fill(spdCorr);
686 fHOutMultTRKvsCL1->Fill(nTracks,spdCorr);
687 } else if (fQuality ==1) {
04341115 688 fHOutCentV0Mqual1->Fill(fCentV0M);
689 fHOutCentTRKqual1->Fill(fCentTRK);
690 fHOutCentCL1qual1->Fill(fCentCL1);
a540a9d3 691 } else {
04341115 692 fHOutCentV0Mqual2->Fill(fCentV0M);
693 fHOutCentTRKqual2->Fill(fCentTRK);
694 fHOutCentCL1qual2->Fill(fCentCL1);
a540a9d3 695 }
696
697 PostData(1, fOutputList);
698}
699
700//________________________________________________________________________
701void AliCentralitySelectionTask::ReadCentralityHistos(TString fCentfilename)
702{
703 // Read centrality histograms
704 TDirectory *owd = gDirectory;
705 // Check if the file is present
706 TString path = gSystem->ExpandPathName(fCentfilename.Data());
707 if (gSystem->AccessPathName(path)) {
708 AliError(Form("File %s does not exist", path.Data()));
709 return;
710 }
711 fFile = TFile::Open(fCentfilename);
712 owd->cd();
713 fHtempV0M = (TH1F*) (fFile->Get("hmultV0_percentile"));
714 fHtempFMD = (TH1F*) (fFile->Get("hmultFMD_percentile"));
715 fHtempTRK = (TH1F*) (fFile->Get("hNtracks_percentile"));
716 fHtempTKL = (TH1F*) (fFile->Get("hNtracklets_percentile"));
717 fHtempCL0 = (TH1F*) (fFile->Get("hNclusters0_percentile"));
718 fHtempCL1 = (TH1F*) (fFile->Get("hNclusters1_percentile"));
719
720 if (!fHtempV0M) AliWarning(Form("Calibration for V0M does not exist in %s", path.Data()));
721 if (!fHtempFMD) AliWarning(Form("Calibration for FMD does not exist in %s", path.Data()));
722 if (!fHtempTRK) AliWarning(Form("Calibration for TRK does not exist in %s", path.Data()));
723 if (!fHtempTKL) AliWarning(Form("Calibration for TKL does not exist in %s", path.Data()));
724 if (!fHtempCL0) AliWarning(Form("Calibration for CL0 does not exist in %s", path.Data()));
725 if (!fHtempCL1) AliWarning(Form("Calibration for CL1 does not exist in %s", path.Data()));
726
727 owd->cd();
728}
729
730//________________________________________________________________________
731void AliCentralitySelectionTask::ReadCentralityHistos2(TString fCentfilename2)
732{
733 // Read centrality histograms
734 TDirectory *owd = gDirectory;
735 TString path = gSystem->ExpandPathName(fCentfilename2.Data());
736 if (gSystem->AccessPathName(path)) {
737 AliError(Form("File %s does not exist", path.Data()));
738 return;
739 }
740 fFile2 = TFile::Open(fCentfilename2);
741 owd->cd();
742 fHtempV0MvsFMD = (TH1F*) (fFile2->Get("hmultV0vsmultFMD_all_percentile"));
743 fHtempTKLvsV0M = (TH1F*) (fFile2->Get("hNtrackletsvsmultV0_all_percentile"));
7ee7a2df 744 fHtempZEMvsZDC = (TH2F*) (fFile2->Get("hEzemvsEzdc_all_percentile"));
a540a9d3 745
746 if (!fHtempV0MvsFMD) AliWarning(Form("Calibration for V0MvsFMD does not exist in %s", path.Data()));
747 if (!fHtempTKLvsV0M) AliWarning(Form("Calibration for TKLvsV0M does not exist in %s", path.Data()));
748 if (!fHtempZEMvsZDC) AliWarning(Form("Calibration for ZEMvsZDC does not exist in %s", path.Data()));
749
750 owd->cd();
751}
752
753//________________________________________________________________________
754void AliCentralitySelectionTask::Terminate(Option_t */*option*/)
755{
756 // Terminate analysis
757 if (fFile && fFile->IsOpen())
758 fFile->Close();
759 if (fFile2 && fFile2->IsOpen())
760 fFile2->Close();
761}
762//________________________________________________________________________
f690bf48 763Int_t AliCentralitySelectionTask::SetupRun(AliESDEvent* const esd)
a540a9d3 764{
765 // Setup files for run
766
767 if (!esd)
768 return -1;
769
770 // check if something to be done
771 if (fCurrentRun == esd->GetRunNumber())
772 return 0;
773 else
774 fCurrentRun = esd->GetRunNumber();
775
776 AliInfo(Form("Setup Centrality Selection for run %d\n",fCurrentRun));
777
778 // CHANGE HERE FOR RUN RANGES
779 if ( fCurrentRun <= 137165 ) fRunNo = 137161;
780 else fRunNo = 137366;
781 // CHANGE HERE FOR RUN RANGES
a121a512 782
783 TString fileName;
784 TString fileName2;
785 if (fPass==1) {
786 fileName=(Form("%s/COMMON/CENTRALITY/data/AliCentralityBy1D_%d.root", AliAnalysisManager::GetOADBPath(), fRunNo));
787 fileName2=(Form("%s/COMMON/CENTRALITY/data/AliCentralityByFunction_%d.root", AliAnalysisManager::GetOADBPath(), fRunNo));
788 } else if (fPass==2) {
789 fileName=(Form("%s/COMMON/CENTRALITY/data/AliCentralityBy1D_%d.root", AliAnalysisManager::GetOADBPath(), fRunNo));
790 fileName2=(Form("%s/COMMON/CENTRALITY/data/AliCentralityByFunction_%d.root", AliAnalysisManager::GetOADBPath(), fRunNo));
791 }
a540a9d3 792
a121a512 793 AliInfo(Form("Centrality Selection for run %d and pass %d is initialized with %s", fCurrentRun, fPass, fileName.Data()));
a540a9d3 794 ReadCentralityHistos(fileName.Data());
795 ReadCentralityHistos2(fileName2.Data());
796 if (!fFile && !fFile2) {
797 AliFatal(Form("Run %d not known to centrality selection!", fCurrentRun));
798 return -1;
799 }
800 return 0;
801}
802
803//________________________________________________________________________
f690bf48 804Bool_t AliCentralitySelectionTask::IsOutlierV0MSPD(Float_t spd, Float_t v0, Int_t cent) const
a540a9d3 805{
f690bf48 806// Clean outliers
f7292781 807 Float_t val= -0.143789 + 0.288874 * v0;
04341115 808 Float_t spdSigma[101]={231.483, 189.446, 183.359, 179.923, 174.229, 170.309, 165.021,
a540a9d3 809 160.84, 159.33, 154.453, 151.644, 148.337, 145.215, 142.353,
810 139.351, 136, 133.838, 129.885, 127.36, 125.032, 122.21, 120.3,
811 117.766, 114.77, 113.1, 110.268, 107.463, 105.293, 102.845,
812 100.835, 98.9632, 97.3287, 93.6887, 92.1066, 89.3224, 87.8382,
813 86.04, 83.6431, 81.9655, 80.0491, 77.8208, 76.4716, 74.2165,
814 72.2752, 70.4875, 68.9414, 66.8622, 65.022, 63.5134, 61.8228,
815 59.7166, 58.5008, 56.2789, 54.9615, 53.386, 51.2165, 49.4842,
816 48.259, 47.1129, 45.3115, 43.8486, 41.9207, 40.5754, 39.3872,
817 38.1897, 36.5401, 35.1283, 33.9702, 32.6429, 31.3612, 29.5876,
818 28.9319, 27.564, 26.0443, 25.2836, 23.9753, 22.8936, 21.5665,
819 20.7048, 19.8016, 18.7095, 18.1144, 17.2095, 16.602, 16.3233,
820 15.7185, 15.3006, 14.7432, 14.4174, 14.0805, 13.7638, 13.7638,
2cf88888 821 13.7638, 13.7638, 13.7638, 13.7638, 13.7638, 13.7638, 13.7638, 18.0803, 18.0803};
a540a9d3 822
04341115 823 if ( TMath::Abs(spd-val) > fOutliersCut*spdSigma[cent] )
a540a9d3 824 return kTRUE;
825 else
826 return kFALSE;
827}
828
829//________________________________________________________________________
f690bf48 830Bool_t AliCentralitySelectionTask::IsOutlierV0MTPC(Int_t tracks, Float_t v0, Int_t cent) const
a540a9d3 831{
f690bf48 832// Clean outliers
f7292781 833 Float_t val = -0.540691 + 0.128358 * v0;
04341115 834 Float_t tpcSigma[101]={106.439, 89.2834, 86.7568, 85.3641, 83.379, 81.6093, 79.3189,
a540a9d3 835 78.0616, 77.2167, 75.0021, 73.9957, 72.0926, 71.0442, 69.8395,
836 68.1169, 66.6676, 66.0038, 64.2284, 63.3845, 61.7439, 60.642,
837 59.5383, 58.3696, 57.0227, 56.0619, 54.7108, 53.8382, 52.3398,
838 51.5297, 49.9488, 49.1126, 48.208, 46.8566, 45.7724, 44.7829,
839 43.8726, 42.7499, 41.9307, 40.6874, 39.9619, 38.5534, 38.0884,
840 36.6141, 35.7482, 34.8804, 34.1769, 33.1278, 32.3435, 31.4783,
841 30.2587, 29.4741, 28.8575, 27.9298, 26.9752, 26.1675, 25.1234,
842 24.4702, 23.6843, 22.9764, 21.8579, 21.2924, 20.3241, 19.8296,
843 19.2465, 18.4474, 17.7216, 16.8956, 16.342, 15.626, 15.0329,
844 14.3911, 13.9301, 13.254, 12.6745, 12.2436, 11.7776, 11.1795,
845 10.673, 10.27, 9.95646, 9.50939, 9.26162, 8.95315, 8.73439,
846 8.67375, 8.43029, 8.34818, 8.33484, 8.40709, 8.3974, 8.32814,
2cf88888 847 8.32814, 8.32814, 8.32814, 8.32814, 8.32814, 8.32814, 8.32814, 8.32814, 12.351, 12.351};
a540a9d3 848
04341115 849 if ( TMath::Abs(tracks-val) > fOutliersCut*tpcSigma[cent] )
a540a9d3 850 return kTRUE;
851 else
852 return kFALSE;
853}
854
855//________________________________________________________________________
e640207a 856Bool_t AliCentralitySelectionTask::IsOutlierV0MZDC(Float_t /*zdc*/, Float_t /*v0*/) const
a540a9d3 857{
f690bf48 858// Clean outliers
e640207a 859 // Float_t val1 = 6350. - 0.26 * v0;
860 // Float_t val2 = 5580.;
861 // if ((zdc > val1) || (zdc > val2))
862 // return kTRUE;
863 // else
a540a9d3 864 return kFALSE;
865}
866
9138a60c 867//________________________________________________________________________
f690bf48 868Bool_t AliCentralitySelectionTask::IsOutlierV0MZDCECal(Float_t /*zdc*/, Float_t /*v0*/) const
9138a60c 869{
f690bf48 870// Clean outliers
9138a60c 871 return kFALSE;
872}
873
a540a9d3 874//________________________________________________________________________
f690bf48 875Float_t AliCentralitySelectionTask::MyGetScaleFactor(Int_t runnumber, Int_t flag) const
a540a9d3 876{
f690bf48 877// Get scaling factor
a540a9d3 878 if (! (runnumber >= fLowRunN && runnumber <=fHighRunN)) {
879 cout << "MyGetScaleFactor error in run number range " << runnumber << endl;
880 return 0.0;
881 }
882
883 Float_t scalefactor=0.0;
884 if (flag==0)
04341115 885 scalefactor = fV0MScaleFactor[runnumber - fLowRunN]; // subtracting reference offset index
a540a9d3 886 else if (flag==1)
04341115 887 scalefactor = fSPDScaleFactor[runnumber - fLowRunN]; // subtracting reference offset index
a540a9d3 888 else if (flag==2)
04341115 889 scalefactor = fTPCScaleFactor[runnumber - fLowRunN]; // subtracting reference offset index
a540a9d3 890
891 return scalefactor;
892
893}
894
895//________________________________________________________________________
f690bf48 896Float_t AliCentralitySelectionTask::MyGetScaleFactorMC(Int_t runnumber) const
a540a9d3 897{
f690bf48 898// Get MC scaling factor
a540a9d3 899 if (! (runnumber >= fLowRunN && runnumber <=fHighRunN)) {
900 cout << "MyGetScaleFactor error in run number range " << runnumber << endl;
901 return 0.0;
902 }
903
04341115 904 Float_t scalefactor= fV0MScaleFactorMC[runnumber - fLowRunN]; // subtracting reference offset index
a540a9d3 905 return scalefactor;
906
907}
908
909//________________________________________________________________________
910void AliCentralitySelectionTask::MyInitScaleFactor ()
911{
f690bf48 912// Initialize the scaling factors
04341115 913 for (int i=0; i<=(fHighRunN-fLowRunN); i++) fV0MScaleFactor[i] = 0.0;
914 for (int i=0; i<=(fHighRunN-fLowRunN); i++) fSPDScaleFactor[i] = 0.0;
915 for (int i=0; i<=(fHighRunN-fLowRunN); i++) fTPCScaleFactor[i] = 0.0;
a540a9d3 916
917 // scale factors determined from <V0 charge> on a run-by-run basis
04341115 918 fV0MScaleFactor[310] = 1.;
919 fV0MScaleFactor[311] = 1.;
920 fV0MScaleFactor[514] = 1.0046;
921 fV0MScaleFactor[515] = 0.983535;
922 fV0MScaleFactor[579] = 0.988185;
923 fV0MScaleFactor[580] = 0.983351;
924 fV0MScaleFactor[581] = 0.989013;
925 fV0MScaleFactor[583] = 0.990056;
926 fV0MScaleFactor[588] = 0.974438;
927 fV0MScaleFactor[589] = 0.981572;
928 fV0MScaleFactor[590] = 0.989316;
929 fV0MScaleFactor[592] = 0.98345;
930 fV0MScaleFactor[688] = 0.993647;
931 fV0MScaleFactor[690] = 0.994758;
932 fV0MScaleFactor[693] = 0.989569;
933 fV0MScaleFactor[698] = 0.993119;
934 fV0MScaleFactor[744] = 0.989583;
935 fV0MScaleFactor[757] = 0.990377;
936 fV0MScaleFactor[787] = 0.990176;
937 fV0MScaleFactor[788] = 0.98723;
938 fV0MScaleFactor[834] = 1.00403;
939 fV0MScaleFactor[835] = 0.994376;
940 fV0MScaleFactor[840] = 0.99759;
941 fV0MScaleFactor[842] = 1.01031;
942 fV0MScaleFactor[900] = 0.996216;
943 fV0MScaleFactor[901] = 0.994205;
944 fV0MScaleFactor[992] = 0.998479;
945 fV0MScaleFactor[997] = 1.00078;
946 fV0MScaleFactor[1299] = 1.00515;
947 fV0MScaleFactor[1303] = 1.00094;
948 fV0MScaleFactor[1339] = 0.986596;
949 fV0MScaleFactor[1346] = 0.972226;
950 fV0MScaleFactor[1349] = 0.960358;
951 fV0MScaleFactor[1350] = 0.970023;
952 fV0MScaleFactor[1374] = 1.00575;
953 fV0MScaleFactor[1545] = 1.00471;
954 fV0MScaleFactor[1587] = 1.00611;
955 fV0MScaleFactor[1588] = 1.00976;
956 fV0MScaleFactor[1618] = 1.00771;
957 fV0MScaleFactor[1682] = 1.01622;
958 fV0MScaleFactor[1727] = 1.01305;
959 fV0MScaleFactor[1728] = 1.00388;
960 fV0MScaleFactor[1731] = 1.00673;
961 fV0MScaleFactor[1732] = 1.00916;
962 fV0MScaleFactor[1770] = 1.0092;
963 fV0MScaleFactor[1773] = 1.00728;
964 fV0MScaleFactor[1786] = 1.01655;
965 fV0MScaleFactor[1787] = 1.00672;
966 fV0MScaleFactor[1801] = 0.983339;
967 fV0MScaleFactor[1802] = 1.00754;
968 fV0MScaleFactor[1811] = 1.00608;
969 fV0MScaleFactor[1815] = 1.01227;
970 fV0MScaleFactor[1879] = 0.99907;
971 fV0MScaleFactor[1881] = 0.995696;
972 fV0MScaleFactor[2186] = 1.00559;
973 fV0MScaleFactor[2187] = 1.00631;
974 fV0MScaleFactor[2254] = 1.01512;
975 fV0MScaleFactor[2256] = 0.998727;
976 fV0MScaleFactor[2321] = 1.00701;
977 fSPDScaleFactor[310] = 1.00211;
978 fSPDScaleFactor[311] = 1.00067;
979 fSPDScaleFactor[514] = 1.02268;
980 fSPDScaleFactor[515] = 0.994902;
981 fSPDScaleFactor[579] = 1.00215;
982 fSPDScaleFactor[580] = 0.993421;
983 fSPDScaleFactor[581] = 1.00129;
984 fSPDScaleFactor[583] = 1.00242;
985 fSPDScaleFactor[588] = 0.984762;
986 fSPDScaleFactor[589] = 0.994355;
987 fSPDScaleFactor[590] = 1.00073;
988 fSPDScaleFactor[592] = 0.995889;
989 fSPDScaleFactor[688] = 0.994532;
990 fSPDScaleFactor[690] = 0.998307;
991 fSPDScaleFactor[693] = 0.994052;
992 fSPDScaleFactor[698] = 0.993224;
993 fSPDScaleFactor[744] = 0.993279;
994 fSPDScaleFactor[757] = 0.992494;
995 fSPDScaleFactor[787] = 0.992678;
996 fSPDScaleFactor[788] = 0.996563;
997 fSPDScaleFactor[834] = 1.01116;
998 fSPDScaleFactor[835] = 0.993108;
999 fSPDScaleFactor[840] = 0.997574;
1000 fSPDScaleFactor[842] = 1.01829;
1001 fSPDScaleFactor[900] = 0.999438;
1002 fSPDScaleFactor[901] = 0.995849;
1003 fSPDScaleFactor[992] = 0.999227;
1004 fSPDScaleFactor[997] = 1.00575;
1005 fSPDScaleFactor[1299] = 0.99877;
1006 fSPDScaleFactor[1303] = 0.999682;
1007 fSPDScaleFactor[1339] = 0.978198;
1008 fSPDScaleFactor[1346] = 0.964178;
1009 fSPDScaleFactor[1349] = 0.959439;
1010 fSPDScaleFactor[1350] = 0.956945;
1011 fSPDScaleFactor[1374] = 0.994434;
1012 fSPDScaleFactor[1545] = 1.0016;
1013 fSPDScaleFactor[1587] = 1.00153;
1014 fSPDScaleFactor[1588] = 1.00698;
1015 fSPDScaleFactor[1618] = 1.00554;
1016 fSPDScaleFactor[1682] = 1.0123;
1017 fSPDScaleFactor[1727] = 1.011;
1018 fSPDScaleFactor[1728] = 1.00143;
1019 fSPDScaleFactor[1731] = 1.00486;
1020 fSPDScaleFactor[1732] = 1.00646;
1021 fSPDScaleFactor[1770] = 1.00515;
1022 fSPDScaleFactor[1773] = 1.00485;
1023 fSPDScaleFactor[1786] = 1.01215;
1024 fSPDScaleFactor[1787] = 1.00419;
1025 fSPDScaleFactor[1801] = 0.983327;
1026 fSPDScaleFactor[1802] = 1.00529;
1027 fSPDScaleFactor[1811] = 1.00367;
1028 fSPDScaleFactor[1815] = 1.01045;
1029 fSPDScaleFactor[1879] = 0.996374;
1030 fSPDScaleFactor[1881] = 0.988827;
1031 fSPDScaleFactor[2186] = 1.00354;
1032 fSPDScaleFactor[2187] = 1.00397;
1033 fSPDScaleFactor[2254] = 1.01138;
1034 fSPDScaleFactor[2256] = 0.996641;
1035 fSPDScaleFactor[2321] = 1.00357;
1036 fTPCScaleFactor[310] = 1.00434;
1037 fTPCScaleFactor[311] = 1.0056;
1038 fTPCScaleFactor[514] = 1.02185;
1039 fTPCScaleFactor[515] = 1.0036;
1040 fTPCScaleFactor[579] = 1.00607;
1041 fTPCScaleFactor[580] = 1.00183;
1042 fTPCScaleFactor[581] = 1.00693;
1043 fTPCScaleFactor[583] = 1.00746;
1044 fTPCScaleFactor[588] = 0.990524;
1045 fTPCScaleFactor[589] = 0.998582;
1046 fTPCScaleFactor[590] = 1.00618;
1047 fTPCScaleFactor[592] = 1.00088;
1048 fTPCScaleFactor[688] = 1.00598;
1049 fTPCScaleFactor[690] = 1.00658;
1050 fTPCScaleFactor[693] = 1.00144;
1051 fTPCScaleFactor[698] = 1.00279;
1052 fTPCScaleFactor[744] = 1.00122;
1053 fTPCScaleFactor[757] = 1.002;
1054 fTPCScaleFactor[787] = 0.997818;
1055 fTPCScaleFactor[788] = 0.994583;
1056 fTPCScaleFactor[834] = 1.01508;
1057 fTPCScaleFactor[835] = 1.00218;
1058 fTPCScaleFactor[840] = 1.00569;
1059 fTPCScaleFactor[842] = 1.01789;
1060 fTPCScaleFactor[900] = 1.00739;
1061 fTPCScaleFactor[901] = 1.00462;
1062 fTPCScaleFactor[992] = 1.00502;
1063 fTPCScaleFactor[997] = 1.00943;
1064 fTPCScaleFactor[1299] = 1.00438;
1065 fTPCScaleFactor[1303] = 0.996701;
1066 fTPCScaleFactor[1339] = 0.978641;
1067 fTPCScaleFactor[1346] = 0.968906;
1068 fTPCScaleFactor[1349] = 0.954311;
1069 fTPCScaleFactor[1350] = 0.958764;
1070 fTPCScaleFactor[1374] = 0.997899;
1071 fTPCScaleFactor[1545] = 0.992;
1072 fTPCScaleFactor[1587] = 0.992635;
1073 fTPCScaleFactor[1588] = 1.00207;
1074 fTPCScaleFactor[1618] = 1.00126;
1075 fTPCScaleFactor[1682] = 1.00324;
1076 fTPCScaleFactor[1727] = 1.00042;
1077 fTPCScaleFactor[1728] = 0.978881;
1078 fTPCScaleFactor[1731] = 0.999818;
1079 fTPCScaleFactor[1732] = 1.00109;
1080 fTPCScaleFactor[1770] = 0.99935;
1081 fTPCScaleFactor[1773] = 0.998531;
1082 fTPCScaleFactor[1786] = 0.999125;
1083 fTPCScaleFactor[1787] = 0.998479;
1084 fTPCScaleFactor[1801] = 0.9775;
1085 fTPCScaleFactor[1802] = 0.999095;
1086 fTPCScaleFactor[1811] = 0.998197;
1087 fTPCScaleFactor[1815] = 1.00413;
1088 fTPCScaleFactor[1879] = 0.990916;
1089 fTPCScaleFactor[1881] = 0.987241;
1090 fTPCScaleFactor[2186] = 1.00048;
1091 fTPCScaleFactor[2187] = 1.00057;
1092 fTPCScaleFactor[2254] = 1.00588;
1093 fTPCScaleFactor[2256] = 0.991503;
1094 fTPCScaleFactor[2321] = 1.00057;
609f601e 1095
a540a9d3 1096
1097 // set all missing values to the value of the run before it ....
bab72d34 1098 for (int i=0; i<=(fHighRunN-fLowRunN); i++) {
04341115 1099 if (fV0MScaleFactor[i] == 0.0) {
a540a9d3 1100 if (i==0) {
04341115 1101 fV0MScaleFactor[i] = 1.00;
a540a9d3 1102 } else {
1103 // search for last run number with non-zero value
1104 for (int j=i-1;j>=0;j--) {
04341115 1105 if (fV0MScaleFactor[j] != 0.0) {
1106 fV0MScaleFactor[i] = fV0MScaleFactor[j];
a540a9d3 1107 break;
1108 }
1109 }
1110 }
1111 }
1112 } // end loop over checking all run-numbers
1113
bab72d34 1114 for (int i=0; i<=(fHighRunN-fLowRunN); i++) {
04341115 1115 if (fSPDScaleFactor[i] == 0.0) {
a540a9d3 1116 if (i==0) {
04341115 1117 fSPDScaleFactor[i] = 1.00;
a540a9d3 1118 } else {
1119 for (int j=i-1;j>=0;j--) {
04341115 1120 if (fSPDScaleFactor[j] != 0.0) {
1121 fSPDScaleFactor[i] = fSPDScaleFactor[j];
a540a9d3 1122 break;
1123 }
1124 }
1125 }
1126 }
1127 }
1128
bab72d34 1129 for (int i=0; i<=(fHighRunN-fLowRunN); i++) {
04341115 1130 if (fTPCScaleFactor[i] == 0.0) {
a540a9d3 1131 if (i==0) {
04341115 1132 fTPCScaleFactor[i] = 1.00;
a540a9d3 1133 } else {
1134 for (int j=i-1;j>=0;j--) {
04341115 1135 if (fTPCScaleFactor[j] != 0.0) {
1136 fTPCScaleFactor[i] = fTPCScaleFactor[j];
a540a9d3 1137 break;
1138 }
1139 }
1140 }
1141 }
1142 }
1143
1144
04341115 1145 // for (int i=0; i<=(fHighRunN-fLowRunN); i++) cout << "Scale Factor = " << fV0MScaleFactor[i] << " for Run " << i+fLowRunN << endl;
a540a9d3 1146
1147 return;
1148
1149}
1150
1151
1152//________________________________________________________________________
1153void AliCentralitySelectionTask::MyInitScaleFactorMC()
1154{
f690bf48 1155// Initialize the MC scaling factors
04341115 1156 for (int i=0; i<(fHighRunN-fLowRunN); i++) fV0MScaleFactorMC[i] = 0.0;
a540a9d3 1157 // scale factors determined from <V0 charge> on a run-by-run basis
04341115 1158 fV0MScaleFactorMC[0] = 0.75108;
a540a9d3 1159 // set all missing values to the value of the run before it ....
bab72d34 1160 for (int i=0; i<=(fHighRunN-fLowRunN); i++) {
04341115 1161 if (fV0MScaleFactorMC[i] == 0.0) {
a540a9d3 1162 if (i==0) {
04341115 1163 fV0MScaleFactorMC[i] = 1.00;
a540a9d3 1164 } else {
1165 // search for last run number with non-zero value
1166 for (int j=i-1;j>=0;j--) {
04341115 1167 if (fV0MScaleFactorMC[j] != 0.0) {
1168 fV0MScaleFactorMC[i] = fV0MScaleFactorMC[j];
a540a9d3 1169 break;
1170 }
1171 }
1172 }
1173 }
1174 } // end loop over checking all run-numbers
1175
1176
04341115 1177 // for (int i=0; i<=(fHighRunN-fLowRunN); i++) cout << "Scale Factor = " << fV0MScaleFactorMC[i] << " for Run " << i+fLowRunN << endl;
a540a9d3 1178
1179 return;
1180
1181}
1182