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