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