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