]> git.uio.no Git - u/mrichter/AliRoot.git/blame - ANALYSIS/AliCentralitySelectionTask.cxx
reflect pass1 and pass2 OADB changes
[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();
afb1125b 442
443 if (fPass==0) AliFatal("Which pass are you analyzing? You should set it via taskCentrality->SetPass(N)");
a540a9d3 444}
445
446//________________________________________________________________________
447void AliCentralitySelectionTask::UserExec(Option_t */*option*/)
448{
449 // Execute analysis for current event:
450 if(fDebug>1) printf(" **** AliCentralitySelectionTask::UserExec() \n");
451
452 Float_t zncEnergy = 0.; // ZNC Energy
453 Float_t zpcEnergy = 0.; // ZPC Energy
454 Float_t znaEnergy = 0.; // ZNA Energy
455 Float_t zpaEnergy = 0.; // ZPA Energy
456 Float_t zem1Energy = 0.; // ZEM1 Energy
457 Float_t zem2Energy = 0.; // ZEM2 Energy
cfead5a9 458 Bool_t zdcEnergyCal = kFALSE; // if zdc is calibrated (in pass2)
459
a540a9d3 460 Int_t nTracks = 0; // no. tracks
461 Int_t nTracklets = 0; // no. tracklets
462 Int_t nClusters[6] = {0}; // no. clusters on 6 ITS layers
463 Int_t nChips[2]; // no. chips on 2 SPD layers
464 Float_t spdCorr =0; // corrected spd2 multiplicity
465
466 Float_t multV0A = 0; // multiplicity from V0 reco side A
467 Float_t multV0C = 0; // multiplicity from V0 reco side C
468 Float_t multFMDA = 0; // multiplicity from FMD on detector A
469 Float_t multFMDC = 0; // multiplicity from FMD on detector C
470
471 Short_t v0Corr = 0; // corrected V0 multiplicity
472 Short_t v0CorrResc = 0; // corrected and rescaled V0 multiplicity
473
474 Float_t zvtx =0; // z-vertex SPD
475
476 AliCentrality *esdCent = 0;
477
478 if(fAnalysisInput.CompareTo("ESD")==0){
479
480 AliVEvent* event = InputEvent();
481 AliESDEvent* esd = dynamic_cast<AliESDEvent*>(event);
482 if (!esd) {
483 AliError("No ESD Event");
484 return;
485 }
cb15f505 486
487 LoadBranches();
a540a9d3 488
489 if (fRunNo<=0) {
490 if (SetupRun(esd)<0)
74b089bb 491 AliError("Centrality File not available for this run");
492 return;
a540a9d3 493 }
494
495 esdCent = esd->GetCentrality();
496
497 // ***** V0 info
498 AliESDVZERO* esdV0 = esd->GetVZEROData();
499 multV0A=esdV0->GetMTotV0A();
500 multV0C=esdV0->GetMTotV0C();
501
502 Float_t v0CorrR;
503 v0Corr = (Short_t)AliESDUtils::GetCorrV0(esd,v0CorrR);
504 v0CorrResc = (Short_t)v0CorrR;
505
506 // ***** Vertex Info
507 const AliESDVertex* vtxESD = esd->GetPrimaryVertexSPD();
508 zvtx = vtxESD->GetZ();
509
510 // ***** CB info (tracklets, clusters, chips)
511 //nTracks = event->GetNumberOfTracks();
512 nTracks = fTrackCuts ? (Short_t)fTrackCuts->GetReferenceMultiplicity(esd,kTRUE):-1;
513
514 const AliMultiplicity *mult = esd->GetMultiplicity();
515
516 nTracklets = mult->GetNumberOfTracklets();
517
518 for(Int_t ilay=0; ilay<6; ilay++){
519 nClusters[ilay] = mult->GetNumberOfITSClusters(ilay);
520 }
521
522 for(Int_t ilay=0; ilay<2; ilay++){
523 nChips[ilay] = mult->GetNumberOfFiredChips(ilay);
524 }
525
526 spdCorr = AliESDUtils::GetCorrSPD2(nClusters[1],zvtx);
527
528 // ***** FMD info
529 AliESDFMD *fmd = esd->GetFMDData();
530 Float_t totalMultA = 0;
531 Float_t totalMultC = 0;
532 const Float_t fFMDLowCut = 0.4;
533
534 for(UShort_t det=1;det<=3;det++) {
535 Int_t nRings = (det==1 ? 1 : 2);
536 for (UShort_t ir = 0; ir < nRings; ir++) {
537 Char_t ring = (ir == 0 ? 'I' : 'O');
538 UShort_t nsec = (ir == 0 ? 20 : 40);
539 UShort_t nstr = (ir == 0 ? 512 : 256);
540 for(UShort_t sec =0; sec < nsec; sec++) {
541 for(UShort_t strip = 0; strip < nstr; strip++) {
542
04341115 543 Float_t fmdMult = fmd->Multiplicity(det,ring,sec,strip);
544 if(fmdMult == 0 || fmdMult == AliESDFMD::kInvalidMult) continue;
a540a9d3 545
546 Float_t nParticles=0;
547
04341115 548 if(fmdMult > fFMDLowCut) {
a540a9d3 549 nParticles = 1.;
550 }
551
552 if (det<3) totalMultA = totalMultA + nParticles;
553 else totalMultC = totalMultC + nParticles;
554
555 }
556 }
557 }
558 }
559 multFMDA = totalMultA;
560 multFMDC = totalMultC;
561
562 // ***** ZDC info
563 AliESDZDC *esdZDC = esd->GetESDZDC();
cfead5a9 564 zdcEnergyCal = esdZDC->AliESDZDC::TestBit(AliESDZDC::kEnergyCalibratedSignal);
565 if (zdcEnergyCal) {
566 zncEnergy = (Float_t) (esdZDC->GetZDCN1Energy());
567 zpcEnergy = (Float_t) (esdZDC->GetZDCP1Energy());
568 znaEnergy = (Float_t) (esdZDC->GetZDCN2Energy());
569 zpaEnergy = (Float_t) (esdZDC->GetZDCP2Energy());
570 } else {
571 zncEnergy = (Float_t) (esdZDC->GetZDCN1Energy())/8.;
572 zpcEnergy = (Float_t) (esdZDC->GetZDCP1Energy())/8.;
573 znaEnergy = (Float_t) (esdZDC->GetZDCN2Energy())/8.;
574 zpaEnergy = (Float_t) (esdZDC->GetZDCP2Energy())/8.;
575 }
a540a9d3 576 zem1Energy = (Float_t) (esdZDC->GetZDCEMEnergy(0))/8.;
577 zem2Energy = (Float_t) (esdZDC->GetZDCEMEnergy(1))/8.;
cfead5a9 578
a540a9d3 579 }
580 else if(fAnalysisInput.CompareTo("AOD")==0){
581 //AliAODEvent *aod = dynamic_cast<AliAODEvent*> (InputEvent());
582 // to be implemented
583 printf(" AOD analysis not yet implemented!!!\n\n");
584 return;
585 }
586
587
588 // ***** Scaling
a121a512 589 // ***** Scaling for pass2
590 if (fPass==2) {
591 fUseScaling=kFALSE;
592 fUseCleaning=kFALSE;
593 }
d015c169 594 // ***** Scaling for MC
595 if (fIsMCInput) {
596 fUseScaling=kFALSE;
04341115 597 Float_t tempScalefactorV0M = MyGetScaleFactorMC(fCurrentRun);
598 v0Corr = Short_t((multV0A+multV0C) * tempScalefactorV0M);
d015c169 599 }
600 // ***** Scaling for Data
a540a9d3 601 if (fUseScaling) {
04341115 602 Float_t tempScalefactorV0M = MyGetScaleFactor(fCurrentRun,0);
603 Float_t tempScalefactorSPD = MyGetScaleFactor(fCurrentRun,1);
604 Float_t tempScalefactorTPC = MyGetScaleFactor(fCurrentRun,2);
605 v0Corr = Short_t(v0Corr / tempScalefactorV0M);
606 spdCorr = spdCorr / tempScalefactorSPD;
607 nTracks = Int_t(nTracks / tempScalefactorTPC);
a540a9d3 608 }
a540a9d3 609
610 // ***** Centrality Selection
611 if(fHtempV0M) fCentV0M = fHtempV0M->GetBinContent(fHtempV0M->FindBin((v0Corr)));
612 if(fHtempFMD) fCentFMD = fHtempFMD->GetBinContent(fHtempFMD->FindBin((multFMDA+multFMDC)));
613 if(fHtempTRK) fCentTRK = fHtempTRK->GetBinContent(fHtempTRK->FindBin(nTracks));
614 if(fHtempTKL) fCentTKL = fHtempTKL->GetBinContent(fHtempTKL->FindBin(nTracklets));
615 if(fHtempCL0) fCentCL0 = fHtempCL0->GetBinContent(fHtempCL0->FindBin(nClusters[0]));
616 if(fHtempCL1) fCentCL1 = fHtempCL1->GetBinContent(fHtempCL1->FindBin(spdCorr));
617
618 if(fHtempV0MvsFMD) fCentV0MvsFMD = fHtempV0MvsFMD->GetBinContent(fHtempV0MvsFMD->FindBin((multV0A+multV0C)));
619 if(fHtempTKLvsV0M) fCentTKLvsV0M = fHtempTKLvsV0M->GetBinContent(fHtempTKLvsV0M->FindBin(nTracklets));
7ee7a2df 620 if(fHtempZEMvsZDC) fCentZEMvsZDC = fHtempZEMvsZDC->GetBinContent(fHtempZEMvsZDC->FindBin(zem1Energy+zem2Energy,zncEnergy+znaEnergy+zpcEnergy+zpaEnergy));
a540a9d3 621
622 // ***** Cleaning
d015c169 623 if (fUseCleaning) {
624 fQuality=0;
625 fZVCut=10;
626 fOutliersCut=6;
627
628 // ***** vertex
629 if (TMath::Abs(zvtx)>fZVCut) fQuality += 1;
630
631 // ***** outliers
632 // **** V0 vs SPD
633 if (IsOutlierV0MSPD(spdCorr, v0Corr, int(fCentV0M))) fQuality += 2;
634 // ***** V0 vs TPC
635 if (IsOutlierV0MTPC(nTracks, v0Corr, int(fCentV0M))) fQuality += 4;
636 // ***** V0 vs ZDC
637 if (IsOutlierV0MZDC((zncEnergy+znaEnergy+zpcEnergy+zpaEnergy), v0Corr) &&
638 (zdcEnergyCal==kFALSE) && !(fIsMCInput)) fQuality += 8;
639 if (IsOutlierV0MZDCECal((zncEnergy+znaEnergy+zpcEnergy+zpaEnergy), v0Corr) &&
640 ((zdcEnergyCal==kTRUE) || (fIsMCInput))) fQuality += 8;
641 } else {
642 fQuality = 0;
643 }
a540a9d3 644
d015c169 645
a540a9d3 646 if (esdCent) {
647 esdCent->SetQuality(fQuality);
648 esdCent->SetCentralityV0M(fCentV0M);
649 esdCent->SetCentralityFMD(fCentFMD);
650 esdCent->SetCentralityTRK(fCentTRK);
651 esdCent->SetCentralityTKL(fCentTKL);
652 esdCent->SetCentralityCL0(fCentCL0);
653 esdCent->SetCentralityCL1(fCentCL1);
654 esdCent->SetCentralityV0MvsFMD(fCentV0MvsFMD);
655 esdCent->SetCentralityTKLvsV0M(fCentTKLvsV0M);
656 esdCent->SetCentralityZEMvsZDC(fCentZEMvsZDC);
657 }
658
659 fHOutQuality->Fill(fQuality);
660 fHOutVertex->Fill(zvtx);
661
31200bdf 662 fHOutMultV0MvsZDN->Fill(v0Corr,(zncEnergy+znaEnergy));
663 fHOutMultZEMvsZDN->Fill((zem1Energy+zem2Energy),(zncEnergy+znaEnergy));
a540a9d3 664 fHOutMultV0MvsZDC->Fill(v0Corr,(zncEnergy+znaEnergy+zpcEnergy+zpaEnergy));
665 fHOutMultZEMvsZDC->Fill((zem1Energy+zem2Energy),(zncEnergy+znaEnergy+zpcEnergy+zpaEnergy));
666 fHOutMultV0MvsCL1->Fill(v0Corr,spdCorr);
667 fHOutMultV0MvsTRK->Fill(v0Corr,nTracks);
668
669 if (fQuality==0) {
670 fHOutCentV0M->Fill(fCentV0M);
671 fHOutCentFMD->Fill(fCentFMD);
672 fHOutCentTRK->Fill(fCentTRK);
673 fHOutCentTKL->Fill(fCentTKL);
674 fHOutCentCL0->Fill(fCentCL0);
675 fHOutCentCL1->Fill(fCentCL1);
676 fHOutCentV0MvsFMD->Fill(fCentV0MvsFMD);
677 fHOutCentTKLvsV0M->Fill(fCentTKLvsV0M);
678 fHOutCentZEMvsZDC->Fill(fCentZEMvsZDC);
9138a60c 679 fHOutCentV0MvsCentCL1->Fill(fCentV0M,fCentCL1);
680 fHOutCentV0MvsCentTRK->Fill(fCentV0M,fCentTRK);
681 fHOutCentTRKvsCentCL1->Fill(fCentTRK,fCentCL1);
a540a9d3 682 fHOutMultV0M->Fill(v0Corr);
683 fHOutMultV0R->Fill(multV0A+multV0C);
684 fHOutMultFMD->Fill((multFMDA+multFMDC));
685 fHOutMultTRK->Fill(nTracks);
686 fHOutMultTKL->Fill(nTracklets);
687 fHOutMultCL0->Fill(nClusters[0]);
688 fHOutMultCL1->Fill(spdCorr);
689 fHOutMultTRKvsCL1->Fill(nTracks,spdCorr);
690 } else if (fQuality ==1) {
04341115 691 fHOutCentV0Mqual1->Fill(fCentV0M);
692 fHOutCentTRKqual1->Fill(fCentTRK);
693 fHOutCentCL1qual1->Fill(fCentCL1);
a540a9d3 694 } else {
04341115 695 fHOutCentV0Mqual2->Fill(fCentV0M);
696 fHOutCentTRKqual2->Fill(fCentTRK);
697 fHOutCentCL1qual2->Fill(fCentCL1);
a540a9d3 698 }
699
700 PostData(1, fOutputList);
701}
702
703//________________________________________________________________________
704void AliCentralitySelectionTask::ReadCentralityHistos(TString fCentfilename)
705{
706 // Read centrality histograms
707 TDirectory *owd = gDirectory;
708 // Check if the file is present
709 TString path = gSystem->ExpandPathName(fCentfilename.Data());
710 if (gSystem->AccessPathName(path)) {
711 AliError(Form("File %s does not exist", path.Data()));
712 return;
713 }
714 fFile = TFile::Open(fCentfilename);
715 owd->cd();
716 fHtempV0M = (TH1F*) (fFile->Get("hmultV0_percentile"));
717 fHtempFMD = (TH1F*) (fFile->Get("hmultFMD_percentile"));
718 fHtempTRK = (TH1F*) (fFile->Get("hNtracks_percentile"));
719 fHtempTKL = (TH1F*) (fFile->Get("hNtracklets_percentile"));
720 fHtempCL0 = (TH1F*) (fFile->Get("hNclusters0_percentile"));
721 fHtempCL1 = (TH1F*) (fFile->Get("hNclusters1_percentile"));
722
723 if (!fHtempV0M) AliWarning(Form("Calibration for V0M does not exist in %s", path.Data()));
724 if (!fHtempFMD) AliWarning(Form("Calibration for FMD does not exist in %s", path.Data()));
725 if (!fHtempTRK) AliWarning(Form("Calibration for TRK does not exist in %s", path.Data()));
726 if (!fHtempTKL) AliWarning(Form("Calibration for TKL does not exist in %s", path.Data()));
727 if (!fHtempCL0) AliWarning(Form("Calibration for CL0 does not exist in %s", path.Data()));
728 if (!fHtempCL1) AliWarning(Form("Calibration for CL1 does not exist in %s", path.Data()));
729
730 owd->cd();
731}
732
733//________________________________________________________________________
734void AliCentralitySelectionTask::ReadCentralityHistos2(TString fCentfilename2)
735{
736 // Read centrality histograms
737 TDirectory *owd = gDirectory;
738 TString path = gSystem->ExpandPathName(fCentfilename2.Data());
739 if (gSystem->AccessPathName(path)) {
740 AliError(Form("File %s does not exist", path.Data()));
741 return;
742 }
743 fFile2 = TFile::Open(fCentfilename2);
744 owd->cd();
745 fHtempV0MvsFMD = (TH1F*) (fFile2->Get("hmultV0vsmultFMD_all_percentile"));
746 fHtempTKLvsV0M = (TH1F*) (fFile2->Get("hNtrackletsvsmultV0_all_percentile"));
7ee7a2df 747 fHtempZEMvsZDC = (TH2F*) (fFile2->Get("hEzemvsEzdc_all_percentile"));
a540a9d3 748
749 if (!fHtempV0MvsFMD) AliWarning(Form("Calibration for V0MvsFMD does not exist in %s", path.Data()));
750 if (!fHtempTKLvsV0M) AliWarning(Form("Calibration for TKLvsV0M does not exist in %s", path.Data()));
751 if (!fHtempZEMvsZDC) AliWarning(Form("Calibration for ZEMvsZDC does not exist in %s", path.Data()));
752
753 owd->cd();
754}
755
756//________________________________________________________________________
757void AliCentralitySelectionTask::Terminate(Option_t */*option*/)
758{
759 // Terminate analysis
760 if (fFile && fFile->IsOpen())
761 fFile->Close();
762 if (fFile2 && fFile2->IsOpen())
763 fFile2->Close();
764}
765//________________________________________________________________________
f690bf48 766Int_t AliCentralitySelectionTask::SetupRun(AliESDEvent* const esd)
a540a9d3 767{
768 // Setup files for run
769
770 if (!esd)
771 return -1;
772
773 // check if something to be done
774 if (fCurrentRun == esd->GetRunNumber())
775 return 0;
776 else
777 fCurrentRun = esd->GetRunNumber();
778
779 AliInfo(Form("Setup Centrality Selection for run %d\n",fCurrentRun));
780
781 // CHANGE HERE FOR RUN RANGES
db1b4c23 782 switch (fPass) {
783 case 1:
784 if ( fCurrentRun <= 137165 ) fRunNo = 137161;
785 else fRunNo = 137366;
786 break;
787 case 2:
788 fRunNo = 139172;
789 break;
790 default:
791 AliError(Form("Run %d not known to centrality selection!", fCurrentRun));
792 return -1;
793 }
794 TString fileName =(Form("%s/COMMON/CENTRALITY/data/pass%d/AliCentralityBy1D_%d.root", AliAnalysisManager::GetOADBPath(), fPass, fRunNo));
795 TString fileName2=(Form("%s/COMMON/CENTRALITY/data/pass%d/AliCentralityByFunction_%d.root", AliAnalysisManager::GetOADBPath(), fPass, fRunNo));
a540a9d3 796
a121a512 797 AliInfo(Form("Centrality Selection for run %d and pass %d is initialized with %s", fCurrentRun, fPass, fileName.Data()));
a540a9d3 798 ReadCentralityHistos(fileName.Data());
799 ReadCentralityHistos2(fileName2.Data());
800 if (!fFile && !fFile2) {
74b089bb 801 AliError(Form("Run %d not known to centrality selection!", fCurrentRun));
a540a9d3 802 return -1;
803 }
804 return 0;
805}
806
807//________________________________________________________________________
f690bf48 808Bool_t AliCentralitySelectionTask::IsOutlierV0MSPD(Float_t spd, Float_t v0, Int_t cent) const
a540a9d3 809{
f690bf48 810// Clean outliers
f7292781 811 Float_t val= -0.143789 + 0.288874 * v0;
04341115 812 Float_t spdSigma[101]={231.483, 189.446, 183.359, 179.923, 174.229, 170.309, 165.021,
a540a9d3 813 160.84, 159.33, 154.453, 151.644, 148.337, 145.215, 142.353,
814 139.351, 136, 133.838, 129.885, 127.36, 125.032, 122.21, 120.3,
815 117.766, 114.77, 113.1, 110.268, 107.463, 105.293, 102.845,
816 100.835, 98.9632, 97.3287, 93.6887, 92.1066, 89.3224, 87.8382,
817 86.04, 83.6431, 81.9655, 80.0491, 77.8208, 76.4716, 74.2165,
818 72.2752, 70.4875, 68.9414, 66.8622, 65.022, 63.5134, 61.8228,
819 59.7166, 58.5008, 56.2789, 54.9615, 53.386, 51.2165, 49.4842,
820 48.259, 47.1129, 45.3115, 43.8486, 41.9207, 40.5754, 39.3872,
821 38.1897, 36.5401, 35.1283, 33.9702, 32.6429, 31.3612, 29.5876,
822 28.9319, 27.564, 26.0443, 25.2836, 23.9753, 22.8936, 21.5665,
823 20.7048, 19.8016, 18.7095, 18.1144, 17.2095, 16.602, 16.3233,
824 15.7185, 15.3006, 14.7432, 14.4174, 14.0805, 13.7638, 13.7638,
2cf88888 825 13.7638, 13.7638, 13.7638, 13.7638, 13.7638, 13.7638, 13.7638, 18.0803, 18.0803};
a540a9d3 826
04341115 827 if ( TMath::Abs(spd-val) > fOutliersCut*spdSigma[cent] )
a540a9d3 828 return kTRUE;
829 else
830 return kFALSE;
831}
832
833//________________________________________________________________________
f690bf48 834Bool_t AliCentralitySelectionTask::IsOutlierV0MTPC(Int_t tracks, Float_t v0, Int_t cent) const
a540a9d3 835{
f690bf48 836// Clean outliers
f7292781 837 Float_t val = -0.540691 + 0.128358 * v0;
04341115 838 Float_t tpcSigma[101]={106.439, 89.2834, 86.7568, 85.3641, 83.379, 81.6093, 79.3189,
a540a9d3 839 78.0616, 77.2167, 75.0021, 73.9957, 72.0926, 71.0442, 69.8395,
840 68.1169, 66.6676, 66.0038, 64.2284, 63.3845, 61.7439, 60.642,
841 59.5383, 58.3696, 57.0227, 56.0619, 54.7108, 53.8382, 52.3398,
842 51.5297, 49.9488, 49.1126, 48.208, 46.8566, 45.7724, 44.7829,
843 43.8726, 42.7499, 41.9307, 40.6874, 39.9619, 38.5534, 38.0884,
844 36.6141, 35.7482, 34.8804, 34.1769, 33.1278, 32.3435, 31.4783,
845 30.2587, 29.4741, 28.8575, 27.9298, 26.9752, 26.1675, 25.1234,
846 24.4702, 23.6843, 22.9764, 21.8579, 21.2924, 20.3241, 19.8296,
847 19.2465, 18.4474, 17.7216, 16.8956, 16.342, 15.626, 15.0329,
848 14.3911, 13.9301, 13.254, 12.6745, 12.2436, 11.7776, 11.1795,
849 10.673, 10.27, 9.95646, 9.50939, 9.26162, 8.95315, 8.73439,
850 8.67375, 8.43029, 8.34818, 8.33484, 8.40709, 8.3974, 8.32814,
2cf88888 851 8.32814, 8.32814, 8.32814, 8.32814, 8.32814, 8.32814, 8.32814, 8.32814, 12.351, 12.351};
a540a9d3 852
04341115 853 if ( TMath::Abs(tracks-val) > fOutliersCut*tpcSigma[cent] )
a540a9d3 854 return kTRUE;
855 else
856 return kFALSE;
857}
858
859//________________________________________________________________________
e640207a 860Bool_t AliCentralitySelectionTask::IsOutlierV0MZDC(Float_t /*zdc*/, Float_t /*v0*/) const
a540a9d3 861{
f690bf48 862// Clean outliers
e640207a 863 // Float_t val1 = 6350. - 0.26 * v0;
864 // Float_t val2 = 5580.;
865 // if ((zdc > val1) || (zdc > val2))
866 // return kTRUE;
867 // else
a540a9d3 868 return kFALSE;
869}
870
9138a60c 871//________________________________________________________________________
f690bf48 872Bool_t AliCentralitySelectionTask::IsOutlierV0MZDCECal(Float_t /*zdc*/, Float_t /*v0*/) const
9138a60c 873{
f690bf48 874// Clean outliers
9138a60c 875 return kFALSE;
876}
877
a540a9d3 878//________________________________________________________________________
f690bf48 879Float_t AliCentralitySelectionTask::MyGetScaleFactor(Int_t runnumber, Int_t flag) const
a540a9d3 880{
f690bf48 881// Get scaling factor
a540a9d3 882 if (! (runnumber >= fLowRunN && runnumber <=fHighRunN)) {
883 cout << "MyGetScaleFactor error in run number range " << runnumber << endl;
884 return 0.0;
885 }
886
887 Float_t scalefactor=0.0;
888 if (flag==0)
04341115 889 scalefactor = fV0MScaleFactor[runnumber - fLowRunN]; // subtracting reference offset index
a540a9d3 890 else if (flag==1)
04341115 891 scalefactor = fSPDScaleFactor[runnumber - fLowRunN]; // subtracting reference offset index
a540a9d3 892 else if (flag==2)
04341115 893 scalefactor = fTPCScaleFactor[runnumber - fLowRunN]; // subtracting reference offset index
a540a9d3 894
895 return scalefactor;
896
897}
898
899//________________________________________________________________________
f690bf48 900Float_t AliCentralitySelectionTask::MyGetScaleFactorMC(Int_t runnumber) const
a540a9d3 901{
f690bf48 902// Get MC scaling factor
a540a9d3 903 if (! (runnumber >= fLowRunN && runnumber <=fHighRunN)) {
904 cout << "MyGetScaleFactor error in run number range " << runnumber << endl;
905 return 0.0;
906 }
907
04341115 908 Float_t scalefactor= fV0MScaleFactorMC[runnumber - fLowRunN]; // subtracting reference offset index
a540a9d3 909 return scalefactor;
910
911}
912
913//________________________________________________________________________
914void AliCentralitySelectionTask::MyInitScaleFactor ()
915{
f690bf48 916// Initialize the scaling factors
04341115 917 for (int i=0; i<=(fHighRunN-fLowRunN); i++) fV0MScaleFactor[i] = 0.0;
918 for (int i=0; i<=(fHighRunN-fLowRunN); i++) fSPDScaleFactor[i] = 0.0;
919 for (int i=0; i<=(fHighRunN-fLowRunN); i++) fTPCScaleFactor[i] = 0.0;
a540a9d3 920
921 // scale factors determined from <V0 charge> on a run-by-run basis
04341115 922 fV0MScaleFactor[310] = 1.;
923 fV0MScaleFactor[311] = 1.;
924 fV0MScaleFactor[514] = 1.0046;
925 fV0MScaleFactor[515] = 0.983535;
926 fV0MScaleFactor[579] = 0.988185;
927 fV0MScaleFactor[580] = 0.983351;
928 fV0MScaleFactor[581] = 0.989013;
929 fV0MScaleFactor[583] = 0.990056;
930 fV0MScaleFactor[588] = 0.974438;
931 fV0MScaleFactor[589] = 0.981572;
932 fV0MScaleFactor[590] = 0.989316;
933 fV0MScaleFactor[592] = 0.98345;
934 fV0MScaleFactor[688] = 0.993647;
935 fV0MScaleFactor[690] = 0.994758;
936 fV0MScaleFactor[693] = 0.989569;
937 fV0MScaleFactor[698] = 0.993119;
938 fV0MScaleFactor[744] = 0.989583;
939 fV0MScaleFactor[757] = 0.990377;
940 fV0MScaleFactor[787] = 0.990176;
941 fV0MScaleFactor[788] = 0.98723;
942 fV0MScaleFactor[834] = 1.00403;
943 fV0MScaleFactor[835] = 0.994376;
944 fV0MScaleFactor[840] = 0.99759;
945 fV0MScaleFactor[842] = 1.01031;
946 fV0MScaleFactor[900] = 0.996216;
947 fV0MScaleFactor[901] = 0.994205;
948 fV0MScaleFactor[992] = 0.998479;
949 fV0MScaleFactor[997] = 1.00078;
950 fV0MScaleFactor[1299] = 1.00515;
951 fV0MScaleFactor[1303] = 1.00094;
952 fV0MScaleFactor[1339] = 0.986596;
953 fV0MScaleFactor[1346] = 0.972226;
954 fV0MScaleFactor[1349] = 0.960358;
955 fV0MScaleFactor[1350] = 0.970023;
956 fV0MScaleFactor[1374] = 1.00575;
957 fV0MScaleFactor[1545] = 1.00471;
958 fV0MScaleFactor[1587] = 1.00611;
959 fV0MScaleFactor[1588] = 1.00976;
960 fV0MScaleFactor[1618] = 1.00771;
961 fV0MScaleFactor[1682] = 1.01622;
962 fV0MScaleFactor[1727] = 1.01305;
963 fV0MScaleFactor[1728] = 1.00388;
964 fV0MScaleFactor[1731] = 1.00673;
965 fV0MScaleFactor[1732] = 1.00916;
966 fV0MScaleFactor[1770] = 1.0092;
967 fV0MScaleFactor[1773] = 1.00728;
968 fV0MScaleFactor[1786] = 1.01655;
969 fV0MScaleFactor[1787] = 1.00672;
970 fV0MScaleFactor[1801] = 0.983339;
971 fV0MScaleFactor[1802] = 1.00754;
972 fV0MScaleFactor[1811] = 1.00608;
973 fV0MScaleFactor[1815] = 1.01227;
974 fV0MScaleFactor[1879] = 0.99907;
975 fV0MScaleFactor[1881] = 0.995696;
976 fV0MScaleFactor[2186] = 1.00559;
977 fV0MScaleFactor[2187] = 1.00631;
978 fV0MScaleFactor[2254] = 1.01512;
979 fV0MScaleFactor[2256] = 0.998727;
980 fV0MScaleFactor[2321] = 1.00701;
981 fSPDScaleFactor[310] = 1.00211;
982 fSPDScaleFactor[311] = 1.00067;
983 fSPDScaleFactor[514] = 1.02268;
984 fSPDScaleFactor[515] = 0.994902;
985 fSPDScaleFactor[579] = 1.00215;
986 fSPDScaleFactor[580] = 0.993421;
987 fSPDScaleFactor[581] = 1.00129;
988 fSPDScaleFactor[583] = 1.00242;
989 fSPDScaleFactor[588] = 0.984762;
990 fSPDScaleFactor[589] = 0.994355;
991 fSPDScaleFactor[590] = 1.00073;
992 fSPDScaleFactor[592] = 0.995889;
993 fSPDScaleFactor[688] = 0.994532;
994 fSPDScaleFactor[690] = 0.998307;
995 fSPDScaleFactor[693] = 0.994052;
996 fSPDScaleFactor[698] = 0.993224;
997 fSPDScaleFactor[744] = 0.993279;
998 fSPDScaleFactor[757] = 0.992494;
999 fSPDScaleFactor[787] = 0.992678;
1000 fSPDScaleFactor[788] = 0.996563;
1001 fSPDScaleFactor[834] = 1.01116;
1002 fSPDScaleFactor[835] = 0.993108;
1003 fSPDScaleFactor[840] = 0.997574;
1004 fSPDScaleFactor[842] = 1.01829;
1005 fSPDScaleFactor[900] = 0.999438;
1006 fSPDScaleFactor[901] = 0.995849;
1007 fSPDScaleFactor[992] = 0.999227;
1008 fSPDScaleFactor[997] = 1.00575;
1009 fSPDScaleFactor[1299] = 0.99877;
1010 fSPDScaleFactor[1303] = 0.999682;
1011 fSPDScaleFactor[1339] = 0.978198;
1012 fSPDScaleFactor[1346] = 0.964178;
1013 fSPDScaleFactor[1349] = 0.959439;
1014 fSPDScaleFactor[1350] = 0.956945;
1015 fSPDScaleFactor[1374] = 0.994434;
1016 fSPDScaleFactor[1545] = 1.0016;
1017 fSPDScaleFactor[1587] = 1.00153;
1018 fSPDScaleFactor[1588] = 1.00698;
1019 fSPDScaleFactor[1618] = 1.00554;
1020 fSPDScaleFactor[1682] = 1.0123;
1021 fSPDScaleFactor[1727] = 1.011;
1022 fSPDScaleFactor[1728] = 1.00143;
1023 fSPDScaleFactor[1731] = 1.00486;
1024 fSPDScaleFactor[1732] = 1.00646;
1025 fSPDScaleFactor[1770] = 1.00515;
1026 fSPDScaleFactor[1773] = 1.00485;
1027 fSPDScaleFactor[1786] = 1.01215;
1028 fSPDScaleFactor[1787] = 1.00419;
1029 fSPDScaleFactor[1801] = 0.983327;
1030 fSPDScaleFactor[1802] = 1.00529;
1031 fSPDScaleFactor[1811] = 1.00367;
1032 fSPDScaleFactor[1815] = 1.01045;
1033 fSPDScaleFactor[1879] = 0.996374;
1034 fSPDScaleFactor[1881] = 0.988827;
1035 fSPDScaleFactor[2186] = 1.00354;
1036 fSPDScaleFactor[2187] = 1.00397;
1037 fSPDScaleFactor[2254] = 1.01138;
1038 fSPDScaleFactor[2256] = 0.996641;
1039 fSPDScaleFactor[2321] = 1.00357;
1040 fTPCScaleFactor[310] = 1.00434;
1041 fTPCScaleFactor[311] = 1.0056;
1042 fTPCScaleFactor[514] = 1.02185;
1043 fTPCScaleFactor[515] = 1.0036;
1044 fTPCScaleFactor[579] = 1.00607;
1045 fTPCScaleFactor[580] = 1.00183;
1046 fTPCScaleFactor[581] = 1.00693;
1047 fTPCScaleFactor[583] = 1.00746;
1048 fTPCScaleFactor[588] = 0.990524;
1049 fTPCScaleFactor[589] = 0.998582;
1050 fTPCScaleFactor[590] = 1.00618;
1051 fTPCScaleFactor[592] = 1.00088;
1052 fTPCScaleFactor[688] = 1.00598;
1053 fTPCScaleFactor[690] = 1.00658;
1054 fTPCScaleFactor[693] = 1.00144;
1055 fTPCScaleFactor[698] = 1.00279;
1056 fTPCScaleFactor[744] = 1.00122;
1057 fTPCScaleFactor[757] = 1.002;
1058 fTPCScaleFactor[787] = 0.997818;
1059 fTPCScaleFactor[788] = 0.994583;
1060 fTPCScaleFactor[834] = 1.01508;
1061 fTPCScaleFactor[835] = 1.00218;
1062 fTPCScaleFactor[840] = 1.00569;
1063 fTPCScaleFactor[842] = 1.01789;
1064 fTPCScaleFactor[900] = 1.00739;
1065 fTPCScaleFactor[901] = 1.00462;
1066 fTPCScaleFactor[992] = 1.00502;
1067 fTPCScaleFactor[997] = 1.00943;
1068 fTPCScaleFactor[1299] = 1.00438;
1069 fTPCScaleFactor[1303] = 0.996701;
1070 fTPCScaleFactor[1339] = 0.978641;
1071 fTPCScaleFactor[1346] = 0.968906;
1072 fTPCScaleFactor[1349] = 0.954311;
1073 fTPCScaleFactor[1350] = 0.958764;
1074 fTPCScaleFactor[1374] = 0.997899;
1075 fTPCScaleFactor[1545] = 0.992;
1076 fTPCScaleFactor[1587] = 0.992635;
1077 fTPCScaleFactor[1588] = 1.00207;
1078 fTPCScaleFactor[1618] = 1.00126;
1079 fTPCScaleFactor[1682] = 1.00324;
1080 fTPCScaleFactor[1727] = 1.00042;
1081 fTPCScaleFactor[1728] = 0.978881;
1082 fTPCScaleFactor[1731] = 0.999818;
1083 fTPCScaleFactor[1732] = 1.00109;
1084 fTPCScaleFactor[1770] = 0.99935;
1085 fTPCScaleFactor[1773] = 0.998531;
1086 fTPCScaleFactor[1786] = 0.999125;
1087 fTPCScaleFactor[1787] = 0.998479;
1088 fTPCScaleFactor[1801] = 0.9775;
1089 fTPCScaleFactor[1802] = 0.999095;
1090 fTPCScaleFactor[1811] = 0.998197;
1091 fTPCScaleFactor[1815] = 1.00413;
1092 fTPCScaleFactor[1879] = 0.990916;
1093 fTPCScaleFactor[1881] = 0.987241;
1094 fTPCScaleFactor[2186] = 1.00048;
1095 fTPCScaleFactor[2187] = 1.00057;
1096 fTPCScaleFactor[2254] = 1.00588;
1097 fTPCScaleFactor[2256] = 0.991503;
1098 fTPCScaleFactor[2321] = 1.00057;
609f601e 1099
a540a9d3 1100
1101 // set all missing values to the value of the run before it ....
bab72d34 1102 for (int i=0; i<=(fHighRunN-fLowRunN); i++) {
04341115 1103 if (fV0MScaleFactor[i] == 0.0) {
a540a9d3 1104 if (i==0) {
04341115 1105 fV0MScaleFactor[i] = 1.00;
a540a9d3 1106 } else {
1107 // search for last run number with non-zero value
1108 for (int j=i-1;j>=0;j--) {
04341115 1109 if (fV0MScaleFactor[j] != 0.0) {
1110 fV0MScaleFactor[i] = fV0MScaleFactor[j];
a540a9d3 1111 break;
1112 }
1113 }
1114 }
1115 }
1116 } // end loop over checking all run-numbers
1117
bab72d34 1118 for (int i=0; i<=(fHighRunN-fLowRunN); i++) {
04341115 1119 if (fSPDScaleFactor[i] == 0.0) {
a540a9d3 1120 if (i==0) {
04341115 1121 fSPDScaleFactor[i] = 1.00;
a540a9d3 1122 } else {
1123 for (int j=i-1;j>=0;j--) {
04341115 1124 if (fSPDScaleFactor[j] != 0.0) {
1125 fSPDScaleFactor[i] = fSPDScaleFactor[j];
a540a9d3 1126 break;
1127 }
1128 }
1129 }
1130 }
1131 }
1132
bab72d34 1133 for (int i=0; i<=(fHighRunN-fLowRunN); i++) {
04341115 1134 if (fTPCScaleFactor[i] == 0.0) {
a540a9d3 1135 if (i==0) {
04341115 1136 fTPCScaleFactor[i] = 1.00;
a540a9d3 1137 } else {
1138 for (int j=i-1;j>=0;j--) {
04341115 1139 if (fTPCScaleFactor[j] != 0.0) {
1140 fTPCScaleFactor[i] = fTPCScaleFactor[j];
a540a9d3 1141 break;
1142 }
1143 }
1144 }
1145 }
1146 }
1147
1148
04341115 1149 // for (int i=0; i<=(fHighRunN-fLowRunN); i++) cout << "Scale Factor = " << fV0MScaleFactor[i] << " for Run " << i+fLowRunN << endl;
a540a9d3 1150
1151 return;
1152
1153}
1154
1155
1156//________________________________________________________________________
1157void AliCentralitySelectionTask::MyInitScaleFactorMC()
1158{
f690bf48 1159// Initialize the MC scaling factors
04341115 1160 for (int i=0; i<(fHighRunN-fLowRunN); i++) fV0MScaleFactorMC[i] = 0.0;
a540a9d3 1161 // scale factors determined from <V0 charge> on a run-by-run basis
04341115 1162 fV0MScaleFactorMC[0] = 0.75108;
a540a9d3 1163 // set all missing values to the value of the run before it ....
bab72d34 1164 for (int i=0; i<=(fHighRunN-fLowRunN); i++) {
04341115 1165 if (fV0MScaleFactorMC[i] == 0.0) {
a540a9d3 1166 if (i==0) {
04341115 1167 fV0MScaleFactorMC[i] = 1.00;
a540a9d3 1168 } else {
1169 // search for last run number with non-zero value
1170 for (int j=i-1;j>=0;j--) {
04341115 1171 if (fV0MScaleFactorMC[j] != 0.0) {
1172 fV0MScaleFactorMC[i] = fV0MScaleFactorMC[j];
a540a9d3 1173 break;
1174 }
1175 }
1176 }
1177 }
1178 } // end loop over checking all run-numbers
1179
1180
04341115 1181 // for (int i=0; i<=(fHighRunN-fLowRunN); i++) cout << "Scale Factor = " << fV0MScaleFactorMC[i] << " for Run " << i+fLowRunN << endl;
a540a9d3 1182
1183 return;
1184
1185}
1186