1 /**************************************************************************
2 * Copyright(c) 1998-2008, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
16 //*****************************************************
17 // Class AliCentralitySelectionTask
18 // Class to analyze determine centrality
19 // author: Alberica Toia
20 //*****************************************************
22 #include "AliCentralitySelectionTask.h"
31 #include <TObjString.h>
35 #include <TDirectory.h>
39 #include "AliAnalysisManager.h"
40 #include "AliVEvent.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"
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"
70 ClassImp(AliCentralitySelectionTask)
73 //________________________________________________________________________
74 AliCentralitySelectionTask::AliCentralitySelectionTask():
76 fAnalysisInput("ESD"),
115 fHOutCentV0MvsFMD(0),
116 fHOutCentTKLvsV0M(0),
117 fHOutCentZEMvsZDC(0),
118 fHOutCentV0MvsCentCL1(0),
119 fHOutCentV0MvsCentTRK(0),
120 fHOutCentTRKvsCentCL1(0),
128 fHOutMultV0MvsZDN(0),
129 fHOutMultZEMvsZDN(0),
130 fHOutMultV0MvsZDC(0),
131 fHOutMultZEMvsZDC(0),
132 fHOutMultV0MvsCL1(0),
133 fHOutMultV0MvsTRK(0),
134 fHOutMultTRKvsCL1(0),
135 fHOutCentV0Mqual1(0),
136 fHOutCentTRKqual1(0),
137 fHOutCentCL1qual1(0),
138 fHOutCentV0Mqual2(0),
139 fHOutCentTRKqual2(0),
140 fHOutCentCL1qual2(0),
144 // Default constructor
145 AliInfo("Centrality Selection enabled.");
149 for (Int_t i=0; i < 2667; i++) {
150 fV0MScaleFactor[i]=0.0;
151 fSPDScaleFactor[i]=0.0;
152 fTPCScaleFactor[i]=0.0;
153 fV0MScaleFactorMC[i]=0.0;
157 fBranchNames="ESD:AliESDRun.,AliESDHeader.,AliESDZDC.,AliESDFMD.,AliESDVZERO."
158 ",SPDVertex.,TPCVertex.,PrimaryVertex.,AliMultiplicity.,Tracks ";
161 //________________________________________________________________________
162 AliCentralitySelectionTask::AliCentralitySelectionTask(const char *name):
163 AliAnalysisTaskSE(name),
164 fAnalysisInput("ESD"),
203 fHOutCentV0MvsFMD(0),
204 fHOutCentTKLvsV0M(0),
205 fHOutCentZEMvsZDC(0),
206 fHOutCentV0MvsCentCL1(0),
207 fHOutCentV0MvsCentTRK(0),
208 fHOutCentTRKvsCentCL1(0),
216 fHOutMultV0MvsZDN(0),
217 fHOutMultZEMvsZDN(0),
218 fHOutMultV0MvsZDC(0),
219 fHOutMultZEMvsZDC(0),
220 fHOutMultV0MvsCL1(0),
221 fHOutMultV0MvsTRK(0),
222 fHOutMultTRKvsCL1(0),
223 fHOutCentV0Mqual1(0),
224 fHOutCentTRKqual1(0),
225 fHOutCentCL1qual1(0),
226 fHOutCentV0Mqual2(0),
227 fHOutCentTRKqual2(0),
228 fHOutCentCL1qual2(0),
232 // Default constructor
236 AliInfo("Centrality Selection enabled.");
237 DefineOutput(1, TList::Class());
238 for (Int_t i=0; i<2667; i++) {
239 fV0MScaleFactor[i]=0.0;
240 fSPDScaleFactor[i]=0.0;
241 fTPCScaleFactor[i]=0.0;
242 fV0MScaleFactorMC[i]=0.0;
246 fBranchNames="ESD:AliESDRun.,AliESDHeader.,AliESDZDC.,AliESDFMD.,AliESDVZERO."
247 ",SPDVertex.,TPCVertex.,PrimaryVertex.,AliMultiplicity.,Tracks ";
250 //________________________________________________________________________
251 AliCentralitySelectionTask& AliCentralitySelectionTask::operator=(const AliCentralitySelectionTask& c)
253 // Assignment operator
255 AliAnalysisTaskSE::operator=(c);
260 //________________________________________________________________________
261 AliCentralitySelectionTask::AliCentralitySelectionTask(const AliCentralitySelectionTask& ana):
262 AliAnalysisTaskSE(ana),
263 fAnalysisInput(ana.fDebug),
264 fIsMCInput(ana.fIsMCInput),
267 fCurrentRun(ana.fCurrentRun),
269 fLowRunN(ana.fLowRunN),
270 fHighRunN(ana.fHighRunN),
271 fUseScaling(ana.fUseScaling),
272 fUseCleaning(ana.fUseCleaning),
273 fTrackCuts(ana.fTrackCuts),
275 fOutliersCut(ana.fOutliersCut),
276 fQuality(ana.fQuality),
277 fCentV0M(ana.fCentV0M),
278 fCentFMD(ana.fCentFMD),
279 fCentTRK(ana.fCentTRK),
280 fCentTKL(ana.fCentTKL),
281 fCentCL0(ana.fCentCL0),
282 fCentCL1(ana.fCentCL1),
283 fCentV0MvsFMD(ana.fCentV0MvsFMD),
284 fCentTKLvsV0M(ana.fCentTKLvsV0M),
285 fCentZEMvsZDC(ana.fCentZEMvsZDC),
286 fHtempV0M(ana.fHtempV0M),
287 fHtempFMD(ana.fHtempFMD),
288 fHtempTRK(ana.fHtempTRK),
289 fHtempTKL(ana.fHtempTKL),
290 fHtempCL0(ana.fHtempCL0),
291 fHtempCL1(ana.fHtempCL1),
292 fHtempV0MvsFMD(ana.fHtempV0MvsFMD),
293 fHtempTKLvsV0M(ana.fHtempTKLvsV0M),
294 fHtempZEMvsZDC(ana.fHtempZEMvsZDC),
295 fOutputList(ana.fOutputList),
296 fHOutCentV0M (ana.fHOutCentV0M ),
297 fHOutCentFMD (ana.fHOutCentFMD ),
298 fHOutCentTRK (ana.fHOutCentTRK ),
299 fHOutCentTKL (ana.fHOutCentTKL ),
300 fHOutCentCL0 (ana.fHOutCentCL0 ),
301 fHOutCentCL1 (ana.fHOutCentCL1 ),
302 fHOutCentV0MvsFMD(ana.fHOutCentV0MvsFMD),
303 fHOutCentTKLvsV0M(ana.fHOutCentTKLvsV0M),
304 fHOutCentZEMvsZDC(ana.fHOutCentZEMvsZDC),
305 fHOutCentV0MvsCentCL1(ana.fHOutCentV0MvsCentCL1),
306 fHOutCentV0MvsCentTRK(ana.fHOutCentV0MvsCentTRK),
307 fHOutCentTRKvsCentCL1(ana.fHOutCentTRKvsCentCL1),
308 fHOutMultV0M(ana.fHOutMultV0M),
309 fHOutMultV0R(ana.fHOutMultV0R),
310 fHOutMultFMD(ana.fHOutMultFMD),
311 fHOutMultTRK(ana.fHOutMultTRK),
312 fHOutMultTKL(ana.fHOutMultTKL),
313 fHOutMultCL0(ana.fHOutMultCL0),
314 fHOutMultCL1(ana.fHOutMultCL1),
315 fHOutMultV0MvsZDN(ana.fHOutMultV0MvsZDN),
316 fHOutMultZEMvsZDN(ana.fHOutMultZEMvsZDN),
317 fHOutMultV0MvsZDC(ana.fHOutMultV0MvsZDC),
318 fHOutMultZEMvsZDC(ana.fHOutMultZEMvsZDC),
319 fHOutMultV0MvsCL1(ana.fHOutMultV0MvsCL1),
320 fHOutMultV0MvsTRK(ana.fHOutMultV0MvsTRK),
321 fHOutMultTRKvsCL1(ana.fHOutMultTRKvsCL1),
322 fHOutCentV0Mqual1(ana.fHOutCentV0Mqual1),
323 fHOutCentTRKqual1(ana.fHOutCentTRKqual1),
324 fHOutCentCL1qual1(ana.fHOutCentCL1qual1),
325 fHOutCentV0Mqual2(ana.fHOutCentV0Mqual2),
326 fHOutCentTRKqual2(ana.fHOutCentTRKqual2),
327 fHOutCentCL1qual2(ana.fHOutCentCL1qual2),
328 fHOutQuality(ana.fHOutQuality),
329 fHOutVertex(ana.fHOutVertex)
332 for (Int_t i=0; i<2667; i++) {
333 fV0MScaleFactor[i]=0.0;
334 fSPDScaleFactor[i]=0.0;
335 fTPCScaleFactor[i]=0.0;
336 fV0MScaleFactorMC[i]=0.0;
341 //________________________________________________________________________
342 AliCentralitySelectionTask::~AliCentralitySelectionTask()
345 if (fOutputList && !AliAnalysisManager::GetAnalysisManager()->IsProofMode()) delete fOutputList;
346 if (fTrackCuts) delete fTrackCuts;
349 //________________________________________________________________________
350 void AliCentralitySelectionTask::UserCreateOutputObjects()
352 // Create the output containers
353 if(fDebug>1) printf("AnalysisCentralitySelectionTask::UserCreateOutputObjects() \n");
354 AliLog::SetClassDebugLevel("AliCentralitySelectionTask", AliLog::kInfo);
356 fOutputList = new TList();
357 fOutputList->SetOwner();
358 fHOutCentV0M = new TH1F("fHOutCentV0M","fHOutCentV0M; Centrality V0",505,0,101);
359 fHOutCentFMD = new TH1F("fHOutCentFMD","fHOutCentFMD; Centrality FMD",505,0,101);
360 fHOutCentTRK = new TH1F("fHOutCentTRK","fHOutCentTRK; Centrality TPC",505,0,101);
361 fHOutCentTKL = new TH1F("fHOutCentTKL","fHOutCentTKL; Centrality tracklets",505,0,101);
362 fHOutCentCL0 = new TH1F("fHOutCentCL0","fHOutCentCL0; Centrality SPD inner",505,0,101);
363 fHOutCentCL1 = new TH1F("fHOutCentCL1","fHOutCentCL1; Centrality SPD outer",505,0,101);
364 fHOutCentV0MvsFMD= new TH1F("fHOutCentV0MvsFMD","fHOutCentV0MvsFMD; Centrality V0 vs FMD",505,0,101);
365 fHOutCentTKLvsV0M= new TH1F("fHOutCentTKLvsV0M","fHOutCentTKLvsV0M; Centrality tracklets vs V0",505,0,101);
366 fHOutCentZEMvsZDC= new TH1F("fHOutCentZEMvsZDC","fHOutCentZEMvsZDC; Centrality ZEM vs ZDC",505,0,101);
367 fHOutCentV0MvsCentCL1= new TH2F("fHOutCentV0MvsCentCL1","fHOutCentV0MvsCentCL1; Cent V0 vs Cent SPD",505,0,101,505,0,101);
368 fHOutCentV0MvsCentTRK= new TH2F("fHOutCentV0MvsCentTRK","fHOutCentV0MvsCentTRK; Cent V0 vs Cent TPC",505,0,101,505,0,101);
369 fHOutCentTRKvsCentCL1= new TH2F("fHOutCentTRKvsCentCL1","fHOutCentTRKvsCentCL1; Cent TPC vs Cent SPD",505,0,101,505,0,101);
371 fHOutMultV0M = new TH1F("fHOutMultV0M","fHOutMultV0M; Multiplicity V0",25000,0,25000);
372 fHOutMultV0R = new TH1F("fHOutMultV0R","fHOutMultV0R; Multiplicity V0",30000,0,30000);
373 fHOutMultFMD = new TH1F("fHOutMultFMD","fHOutMultFMD; Multiplicity FMD",24000,0,24000);
374 fHOutMultTRK = new TH1F("fHOutMultTRK","fHOutMultTRK; Multiplicity TPC",4000,0,4000);
375 fHOutMultTKL = new TH1F("fHOutMultTKL","fHOutMultTKL; Multiplicity tracklets",5000,0,5000);
376 fHOutMultCL0 = new TH1F("fHOutMultCL0","fHOutMultCL0; Multiplicity SPD inner",7000,0,7000);
377 fHOutMultCL1 = new TH1F("fHOutMultCL1","fHOutMultCL1; Multiplicity SPD outer",7000,0,7000);
378 fHOutMultV0MvsZDN = new TH2F("fHOutMultV0MvsZDN","fHOutMultV0MvsZDN; Multiplicity V0; Energy ZDC-N",500,0,25000,500,0,180000);
379 fHOutMultZEMvsZDN = new TH2F("fHOutMultZEMvsZDN","fHOutMultZEMvsZDN; Energy ZEM; Energy ZDC-N",500,0,2500,500,0,180000);
380 fHOutMultV0MvsZDC = new TH2F("fHOutMultV0MvsZDC","fHOutMultV0MvsZDC; Multiplicity V0; Energy ZDC",500,0,25000,500,0,200000);
381 fHOutMultZEMvsZDC = new TH2F("fHOutMultZEMvsZDC","fHOutMultZEMvsZDC; Energy ZEM; Energy ZDC",500,0,2500,500,0,200000);
382 fHOutMultV0MvsCL1 = new TH2F("fHOutMultV0MvsCL1","fHOutMultV0MvsCL1; Multiplicity V0; Multiplicity SPD outer",2500,0,25000,700,0,7000);
383 fHOutMultV0MvsTRK = new TH2F("fHOutMultV0MvsTRK","fHOutMultV0MvsTRK; Multiplicity V0; Multiplicity TPC",2500,0,25000,400,0,4000);
384 fHOutMultTRKvsCL1 = new TH2F("fHOutMultTRKvsCL1","fHOutMultTRKvsCL1; Multiplicity TPC; Multiplicity SPD outer",400,0,4000,700,0,7000);
386 fHOutCentV0Mqual1 = new TH1F("fHOutCentV0M_qual1","fHOutCentV0M_qual1; Centrality V0",505,0,101);
387 fHOutCentTRKqual1 = new TH1F("fHOutCentTRK_qual1","fHOutCentTRK_qual1; Centrality TPC",505,0,101);
388 fHOutCentCL1qual1 = new TH1F("fHOutCentCL1_qual1","fHOutCentCL1_qual1; Centrality SPD outer",505,0,101);
390 fHOutCentV0Mqual2 = new TH1F("fHOutCentV0M_qual2","fHOutCentV0M_qual2; Centrality V0",505,0,101);
391 fHOutCentTRKqual2 = new TH1F("fHOutCentTRK_qual2","fHOutCentTRK_qual2; Centrality TPC",505,0,101);
392 fHOutCentCL1qual2 = new TH1F("fHOutCentCL1_qual2","fHOutCentCL1_qual2; Centrality SPD outer",505,0,101);
394 fHOutQuality = new TH1F("fHOutQuality", "fHOutQuality", 10,-0.5,9.5);
395 fHOutVertex = new TH1F("fHOutVertex", "fHOutVertex", 100,-20,20);
397 fOutputList->Add( fHOutCentV0M );
398 fOutputList->Add( fHOutCentFMD );
399 fOutputList->Add( fHOutCentTRK );
400 fOutputList->Add( fHOutCentTKL );
401 fOutputList->Add( fHOutCentCL0 );
402 fOutputList->Add( fHOutCentCL1 );
403 fOutputList->Add( fHOutCentV0MvsFMD);
404 fOutputList->Add( fHOutCentTKLvsV0M);
405 fOutputList->Add( fHOutCentZEMvsZDC);
406 fOutputList->Add( fHOutCentV0MvsCentCL1);
407 fOutputList->Add( fHOutCentV0MvsCentTRK);
408 fOutputList->Add( fHOutCentTRKvsCentCL1);
409 fOutputList->Add( fHOutMultV0M);
410 fOutputList->Add( fHOutMultV0R);
411 fOutputList->Add( fHOutMultFMD);
412 fOutputList->Add( fHOutMultTRK);
413 fOutputList->Add( fHOutMultTKL);
414 fOutputList->Add( fHOutMultCL0);
415 fOutputList->Add( fHOutMultCL1);
416 fOutputList->Add( fHOutMultV0MvsZDN);
417 fOutputList->Add( fHOutMultZEMvsZDN);
418 fOutputList->Add( fHOutMultV0MvsZDC);
419 fOutputList->Add( fHOutMultZEMvsZDC);
420 fOutputList->Add( fHOutMultV0MvsCL1);
421 fOutputList->Add( fHOutMultV0MvsTRK);
422 fOutputList->Add( fHOutMultTRKvsCL1);
423 fOutputList->Add( fHOutCentV0Mqual1 );
424 fOutputList->Add( fHOutCentTRKqual1 );
425 fOutputList->Add( fHOutCentCL1qual1 );
426 fOutputList->Add( fHOutCentV0Mqual2 );
427 fOutputList->Add( fHOutCentTRKqual2 );
428 fOutputList->Add( fHOutCentCL1qual2 );
429 fOutputList->Add( fHOutQuality );
430 fOutputList->Add( fHOutVertex );
433 fTrackCuts = AliESDtrackCuts::GetStandardTPCOnlyTrackCuts();
435 PostData(1, fOutputList);
438 if (fIsMCInput) MyInitScaleFactorMC();
441 //________________________________________________________________________
442 void AliCentralitySelectionTask::UserExec(Option_t */*option*/)
444 // Execute analysis for current event:
445 if(fDebug>1) printf(" **** AliCentralitySelectionTask::UserExec() \n");
447 Float_t zncEnergy = 0.; // ZNC Energy
448 Float_t zpcEnergy = 0.; // ZPC Energy
449 Float_t znaEnergy = 0.; // ZNA Energy
450 Float_t zpaEnergy = 0.; // ZPA Energy
451 Float_t zem1Energy = 0.; // ZEM1 Energy
452 Float_t zem2Energy = 0.; // ZEM2 Energy
453 Bool_t zdcEnergyCal = kFALSE; // if zdc is calibrated (in pass2)
455 Int_t nTracks = 0; // no. tracks
456 Int_t nTracklets = 0; // no. tracklets
457 Int_t nClusters[6] = {0}; // no. clusters on 6 ITS layers
458 Int_t nChips[2]; // no. chips on 2 SPD layers
459 Float_t spdCorr =0; // corrected spd2 multiplicity
461 Float_t multV0A = 0; // multiplicity from V0 reco side A
462 Float_t multV0C = 0; // multiplicity from V0 reco side C
463 Float_t multFMDA = 0; // multiplicity from FMD on detector A
464 Float_t multFMDC = 0; // multiplicity from FMD on detector C
466 Short_t v0Corr = 0; // corrected V0 multiplicity
467 Short_t v0CorrResc = 0; // corrected and rescaled V0 multiplicity
469 Float_t zvtx =0; // z-vertex SPD
471 AliCentrality *esdCent = 0;
473 if(fAnalysisInput.CompareTo("ESD")==0){
475 AliVEvent* event = InputEvent();
476 AliESDEvent* esd = dynamic_cast<AliESDEvent*>(event);
478 AliError("No ESD Event");
486 AliFatal("Centrality File not available for this run");
489 esdCent = esd->GetCentrality();
492 AliESDVZERO* esdV0 = esd->GetVZEROData();
493 multV0A=esdV0->GetMTotV0A();
494 multV0C=esdV0->GetMTotV0C();
497 v0Corr = (Short_t)AliESDUtils::GetCorrV0(esd,v0CorrR);
498 v0CorrResc = (Short_t)v0CorrR;
501 const AliESDVertex* vtxESD = esd->GetPrimaryVertexSPD();
502 zvtx = vtxESD->GetZ();
504 // ***** CB info (tracklets, clusters, chips)
505 //nTracks = event->GetNumberOfTracks();
506 nTracks = fTrackCuts ? (Short_t)fTrackCuts->GetReferenceMultiplicity(esd,kTRUE):-1;
508 const AliMultiplicity *mult = esd->GetMultiplicity();
510 nTracklets = mult->GetNumberOfTracklets();
512 for(Int_t ilay=0; ilay<6; ilay++){
513 nClusters[ilay] = mult->GetNumberOfITSClusters(ilay);
516 for(Int_t ilay=0; ilay<2; ilay++){
517 nChips[ilay] = mult->GetNumberOfFiredChips(ilay);
520 spdCorr = AliESDUtils::GetCorrSPD2(nClusters[1],zvtx);
523 AliESDFMD *fmd = esd->GetFMDData();
524 Float_t totalMultA = 0;
525 Float_t totalMultC = 0;
526 const Float_t fFMDLowCut = 0.4;
528 for(UShort_t det=1;det<=3;det++) {
529 Int_t nRings = (det==1 ? 1 : 2);
530 for (UShort_t ir = 0; ir < nRings; ir++) {
531 Char_t ring = (ir == 0 ? 'I' : 'O');
532 UShort_t nsec = (ir == 0 ? 20 : 40);
533 UShort_t nstr = (ir == 0 ? 512 : 256);
534 for(UShort_t sec =0; sec < nsec; sec++) {
535 for(UShort_t strip = 0; strip < nstr; strip++) {
537 Float_t fmdMult = fmd->Multiplicity(det,ring,sec,strip);
538 if(fmdMult == 0 || fmdMult == AliESDFMD::kInvalidMult) continue;
540 Float_t nParticles=0;
542 if(fmdMult > fFMDLowCut) {
546 if (det<3) totalMultA = totalMultA + nParticles;
547 else totalMultC = totalMultC + nParticles;
553 multFMDA = totalMultA;
554 multFMDC = totalMultC;
557 AliESDZDC *esdZDC = esd->GetESDZDC();
558 zdcEnergyCal = esdZDC->AliESDZDC::TestBit(AliESDZDC::kEnergyCalibratedSignal);
560 zncEnergy = (Float_t) (esdZDC->GetZDCN1Energy());
561 zpcEnergy = (Float_t) (esdZDC->GetZDCP1Energy());
562 znaEnergy = (Float_t) (esdZDC->GetZDCN2Energy());
563 zpaEnergy = (Float_t) (esdZDC->GetZDCP2Energy());
565 zncEnergy = (Float_t) (esdZDC->GetZDCN1Energy())/8.;
566 zpcEnergy = (Float_t) (esdZDC->GetZDCP1Energy())/8.;
567 znaEnergy = (Float_t) (esdZDC->GetZDCN2Energy())/8.;
568 zpaEnergy = (Float_t) (esdZDC->GetZDCP2Energy())/8.;
570 zem1Energy = (Float_t) (esdZDC->GetZDCEMEnergy(0))/8.;
571 zem2Energy = (Float_t) (esdZDC->GetZDCEMEnergy(1))/8.;
574 else if(fAnalysisInput.CompareTo("AOD")==0){
575 //AliAODEvent *aod = dynamic_cast<AliAODEvent*> (InputEvent());
577 printf(" AOD analysis not yet implemented!!!\n\n");
583 // ***** Scaling for MC
586 Float_t tempScalefactorV0M = MyGetScaleFactorMC(fCurrentRun);
587 v0Corr = Short_t((multV0A+multV0C) * tempScalefactorV0M);
589 // ***** Scaling for Data
591 Float_t tempScalefactorV0M = MyGetScaleFactor(fCurrentRun,0);
592 Float_t tempScalefactorSPD = MyGetScaleFactor(fCurrentRun,1);
593 Float_t tempScalefactorTPC = MyGetScaleFactor(fCurrentRun,2);
594 v0Corr = Short_t(v0Corr / tempScalefactorV0M);
595 spdCorr = spdCorr / tempScalefactorSPD;
596 nTracks = Int_t(nTracks / tempScalefactorTPC);
599 // ***** Centrality Selection
600 if(fHtempV0M) fCentV0M = fHtempV0M->GetBinContent(fHtempV0M->FindBin((v0Corr)));
601 if(fHtempFMD) fCentFMD = fHtempFMD->GetBinContent(fHtempFMD->FindBin((multFMDA+multFMDC)));
602 if(fHtempTRK) fCentTRK = fHtempTRK->GetBinContent(fHtempTRK->FindBin(nTracks));
603 if(fHtempTKL) fCentTKL = fHtempTKL->GetBinContent(fHtempTKL->FindBin(nTracklets));
604 if(fHtempCL0) fCentCL0 = fHtempCL0->GetBinContent(fHtempCL0->FindBin(nClusters[0]));
605 if(fHtempCL1) fCentCL1 = fHtempCL1->GetBinContent(fHtempCL1->FindBin(spdCorr));
607 if(fHtempV0MvsFMD) fCentV0MvsFMD = fHtempV0MvsFMD->GetBinContent(fHtempV0MvsFMD->FindBin((multV0A+multV0C)));
608 if(fHtempTKLvsV0M) fCentTKLvsV0M = fHtempTKLvsV0M->GetBinContent(fHtempTKLvsV0M->FindBin(nTracklets));
609 if(fHtempZEMvsZDC) fCentZEMvsZDC = fHtempZEMvsZDC->GetBinContent(fHtempZEMvsZDC->FindBin(zem1Energy+zem2Energy,zncEnergy+znaEnergy+zpcEnergy+zpaEnergy));
618 if (TMath::Abs(zvtx)>fZVCut) fQuality += 1;
622 if (IsOutlierV0MSPD(spdCorr, v0Corr, int(fCentV0M))) fQuality += 2;
624 if (IsOutlierV0MTPC(nTracks, v0Corr, int(fCentV0M))) fQuality += 4;
626 if (IsOutlierV0MZDC((zncEnergy+znaEnergy+zpcEnergy+zpaEnergy), v0Corr) &&
627 (zdcEnergyCal==kFALSE) && !(fIsMCInput)) fQuality += 8;
628 if (IsOutlierV0MZDCECal((zncEnergy+znaEnergy+zpcEnergy+zpaEnergy), v0Corr) &&
629 ((zdcEnergyCal==kTRUE) || (fIsMCInput))) fQuality += 8;
636 esdCent->SetQuality(fQuality);
637 esdCent->SetCentralityV0M(fCentV0M);
638 esdCent->SetCentralityFMD(fCentFMD);
639 esdCent->SetCentralityTRK(fCentTRK);
640 esdCent->SetCentralityTKL(fCentTKL);
641 esdCent->SetCentralityCL0(fCentCL0);
642 esdCent->SetCentralityCL1(fCentCL1);
643 esdCent->SetCentralityV0MvsFMD(fCentV0MvsFMD);
644 esdCent->SetCentralityTKLvsV0M(fCentTKLvsV0M);
645 esdCent->SetCentralityZEMvsZDC(fCentZEMvsZDC);
648 fHOutQuality->Fill(fQuality);
649 fHOutVertex->Fill(zvtx);
651 fHOutMultV0MvsZDN->Fill(v0Corr,(zncEnergy+znaEnergy));
652 fHOutMultZEMvsZDN->Fill((zem1Energy+zem2Energy),(zncEnergy+znaEnergy));
653 fHOutMultV0MvsZDC->Fill(v0Corr,(zncEnergy+znaEnergy+zpcEnergy+zpaEnergy));
654 fHOutMultZEMvsZDC->Fill((zem1Energy+zem2Energy),(zncEnergy+znaEnergy+zpcEnergy+zpaEnergy));
655 fHOutMultV0MvsCL1->Fill(v0Corr,spdCorr);
656 fHOutMultV0MvsTRK->Fill(v0Corr,nTracks);
659 fHOutCentV0M->Fill(fCentV0M);
660 fHOutCentFMD->Fill(fCentFMD);
661 fHOutCentTRK->Fill(fCentTRK);
662 fHOutCentTKL->Fill(fCentTKL);
663 fHOutCentCL0->Fill(fCentCL0);
664 fHOutCentCL1->Fill(fCentCL1);
665 fHOutCentV0MvsFMD->Fill(fCentV0MvsFMD);
666 fHOutCentTKLvsV0M->Fill(fCentTKLvsV0M);
667 fHOutCentZEMvsZDC->Fill(fCentZEMvsZDC);
668 fHOutCentV0MvsCentCL1->Fill(fCentV0M,fCentCL1);
669 fHOutCentV0MvsCentTRK->Fill(fCentV0M,fCentTRK);
670 fHOutCentTRKvsCentCL1->Fill(fCentTRK,fCentCL1);
671 fHOutMultV0M->Fill(v0Corr);
672 fHOutMultV0R->Fill(multV0A+multV0C);
673 fHOutMultFMD->Fill((multFMDA+multFMDC));
674 fHOutMultTRK->Fill(nTracks);
675 fHOutMultTKL->Fill(nTracklets);
676 fHOutMultCL0->Fill(nClusters[0]);
677 fHOutMultCL1->Fill(spdCorr);
678 fHOutMultTRKvsCL1->Fill(nTracks,spdCorr);
679 } else if (fQuality ==1) {
680 fHOutCentV0Mqual1->Fill(fCentV0M);
681 fHOutCentTRKqual1->Fill(fCentTRK);
682 fHOutCentCL1qual1->Fill(fCentCL1);
684 fHOutCentV0Mqual2->Fill(fCentV0M);
685 fHOutCentTRKqual2->Fill(fCentTRK);
686 fHOutCentCL1qual2->Fill(fCentCL1);
689 PostData(1, fOutputList);
692 //________________________________________________________________________
693 void AliCentralitySelectionTask::ReadCentralityHistos(TString fCentfilename)
695 // Read centrality histograms
696 TDirectory *owd = gDirectory;
697 // Check if the file is present
698 TString path = gSystem->ExpandPathName(fCentfilename.Data());
699 if (gSystem->AccessPathName(path)) {
700 AliError(Form("File %s does not exist", path.Data()));
703 fFile = TFile::Open(fCentfilename);
705 fHtempV0M = (TH1F*) (fFile->Get("hmultV0_percentile"));
706 fHtempFMD = (TH1F*) (fFile->Get("hmultFMD_percentile"));
707 fHtempTRK = (TH1F*) (fFile->Get("hNtracks_percentile"));
708 fHtempTKL = (TH1F*) (fFile->Get("hNtracklets_percentile"));
709 fHtempCL0 = (TH1F*) (fFile->Get("hNclusters0_percentile"));
710 fHtempCL1 = (TH1F*) (fFile->Get("hNclusters1_percentile"));
712 if (!fHtempV0M) AliWarning(Form("Calibration for V0M does not exist in %s", path.Data()));
713 if (!fHtempFMD) AliWarning(Form("Calibration for FMD does not exist in %s", path.Data()));
714 if (!fHtempTRK) AliWarning(Form("Calibration for TRK does not exist in %s", path.Data()));
715 if (!fHtempTKL) AliWarning(Form("Calibration for TKL does not exist in %s", path.Data()));
716 if (!fHtempCL0) AliWarning(Form("Calibration for CL0 does not exist in %s", path.Data()));
717 if (!fHtempCL1) AliWarning(Form("Calibration for CL1 does not exist in %s", path.Data()));
722 //________________________________________________________________________
723 void AliCentralitySelectionTask::ReadCentralityHistos2(TString fCentfilename2)
725 // Read centrality histograms
726 TDirectory *owd = gDirectory;
727 TString path = gSystem->ExpandPathName(fCentfilename2.Data());
728 if (gSystem->AccessPathName(path)) {
729 AliError(Form("File %s does not exist", path.Data()));
732 fFile2 = TFile::Open(fCentfilename2);
734 fHtempV0MvsFMD = (TH1F*) (fFile2->Get("hmultV0vsmultFMD_all_percentile"));
735 fHtempTKLvsV0M = (TH1F*) (fFile2->Get("hNtrackletsvsmultV0_all_percentile"));
736 fHtempZEMvsZDC = (TH2F*) (fFile2->Get("hEzemvsEzdc_all_percentile"));
738 if (!fHtempV0MvsFMD) AliWarning(Form("Calibration for V0MvsFMD does not exist in %s", path.Data()));
739 if (!fHtempTKLvsV0M) AliWarning(Form("Calibration for TKLvsV0M does not exist in %s", path.Data()));
740 if (!fHtempZEMvsZDC) AliWarning(Form("Calibration for ZEMvsZDC does not exist in %s", path.Data()));
745 //________________________________________________________________________
746 void AliCentralitySelectionTask::Terminate(Option_t */*option*/)
748 // Terminate analysis
749 if (fFile && fFile->IsOpen())
751 if (fFile2 && fFile2->IsOpen())
754 //________________________________________________________________________
755 Int_t AliCentralitySelectionTask::SetupRun(AliESDEvent* const esd)
757 // Setup files for run
762 // check if something to be done
763 if (fCurrentRun == esd->GetRunNumber())
766 fCurrentRun = esd->GetRunNumber();
768 AliInfo(Form("Setup Centrality Selection for run %d\n",fCurrentRun));
770 // CHANGE HERE FOR RUN RANGES
771 if ( fCurrentRun <= 137165 ) fRunNo = 137161;
772 else fRunNo = 137366;
773 // CHANGE HERE FOR RUN RANGES
775 TString fileName(Form("%s/COMMON/CENTRALITY/data/AliCentralityBy1D_%d.root", AliAnalysisManager::GetOADBPath(), fRunNo));
776 TString fileName2(Form("%s/COMMON/CENTRALITY/data/AliCentralityByFunction_%d.root", AliAnalysisManager::GetOADBPath(), fRunNo));
778 AliInfo(Form("Centrality Selection for run %d is initialized with %s", fCurrentRun, fileName.Data()));
779 ReadCentralityHistos(fileName.Data());
780 ReadCentralityHistos2(fileName2.Data());
781 if (!fFile && !fFile2) {
782 AliFatal(Form("Run %d not known to centrality selection!", fCurrentRun));
788 //________________________________________________________________________
789 Bool_t AliCentralitySelectionTask::IsOutlierV0MSPD(Float_t spd, Float_t v0, Int_t cent) const
792 Float_t val= -0.143789 + 0.288874 * v0;
793 Float_t spdSigma[101]={231.483, 189.446, 183.359, 179.923, 174.229, 170.309, 165.021,
794 160.84, 159.33, 154.453, 151.644, 148.337, 145.215, 142.353,
795 139.351, 136, 133.838, 129.885, 127.36, 125.032, 122.21, 120.3,
796 117.766, 114.77, 113.1, 110.268, 107.463, 105.293, 102.845,
797 100.835, 98.9632, 97.3287, 93.6887, 92.1066, 89.3224, 87.8382,
798 86.04, 83.6431, 81.9655, 80.0491, 77.8208, 76.4716, 74.2165,
799 72.2752, 70.4875, 68.9414, 66.8622, 65.022, 63.5134, 61.8228,
800 59.7166, 58.5008, 56.2789, 54.9615, 53.386, 51.2165, 49.4842,
801 48.259, 47.1129, 45.3115, 43.8486, 41.9207, 40.5754, 39.3872,
802 38.1897, 36.5401, 35.1283, 33.9702, 32.6429, 31.3612, 29.5876,
803 28.9319, 27.564, 26.0443, 25.2836, 23.9753, 22.8936, 21.5665,
804 20.7048, 19.8016, 18.7095, 18.1144, 17.2095, 16.602, 16.3233,
805 15.7185, 15.3006, 14.7432, 14.4174, 14.0805, 13.7638, 13.7638,
806 13.7638, 13.7638, 13.7638, 13.7638, 13.7638, 13.7638, 13.7638, 18.0803, 18.0803};
808 if ( TMath::Abs(spd-val) > fOutliersCut*spdSigma[cent] )
814 //________________________________________________________________________
815 Bool_t AliCentralitySelectionTask::IsOutlierV0MTPC(Int_t tracks, Float_t v0, Int_t cent) const
818 Float_t val = -0.540691 + 0.128358 * v0;
819 Float_t tpcSigma[101]={106.439, 89.2834, 86.7568, 85.3641, 83.379, 81.6093, 79.3189,
820 78.0616, 77.2167, 75.0021, 73.9957, 72.0926, 71.0442, 69.8395,
821 68.1169, 66.6676, 66.0038, 64.2284, 63.3845, 61.7439, 60.642,
822 59.5383, 58.3696, 57.0227, 56.0619, 54.7108, 53.8382, 52.3398,
823 51.5297, 49.9488, 49.1126, 48.208, 46.8566, 45.7724, 44.7829,
824 43.8726, 42.7499, 41.9307, 40.6874, 39.9619, 38.5534, 38.0884,
825 36.6141, 35.7482, 34.8804, 34.1769, 33.1278, 32.3435, 31.4783,
826 30.2587, 29.4741, 28.8575, 27.9298, 26.9752, 26.1675, 25.1234,
827 24.4702, 23.6843, 22.9764, 21.8579, 21.2924, 20.3241, 19.8296,
828 19.2465, 18.4474, 17.7216, 16.8956, 16.342, 15.626, 15.0329,
829 14.3911, 13.9301, 13.254, 12.6745, 12.2436, 11.7776, 11.1795,
830 10.673, 10.27, 9.95646, 9.50939, 9.26162, 8.95315, 8.73439,
831 8.67375, 8.43029, 8.34818, 8.33484, 8.40709, 8.3974, 8.32814,
832 8.32814, 8.32814, 8.32814, 8.32814, 8.32814, 8.32814, 8.32814, 8.32814, 12.351, 12.351};
834 if ( TMath::Abs(tracks-val) > fOutliersCut*tpcSigma[cent] )
840 //________________________________________________________________________
841 Bool_t AliCentralitySelectionTask::IsOutlierV0MZDC(Float_t zdc, Float_t v0) const
844 Float_t val1 = 6350. - 0.26 * v0;
845 Float_t val2 = 5580.;
846 if ((zdc > val1) || (zdc > val2))
852 //________________________________________________________________________
853 Bool_t AliCentralitySelectionTask::IsOutlierV0MZDCECal(Float_t /*zdc*/, Float_t /*v0*/) const
859 //________________________________________________________________________
860 Float_t AliCentralitySelectionTask::MyGetScaleFactor(Int_t runnumber, Int_t flag) const
862 // Get scaling factor
863 if (! (runnumber >= fLowRunN && runnumber <=fHighRunN)) {
864 cout << "MyGetScaleFactor error in run number range " << runnumber << endl;
868 Float_t scalefactor=0.0;
870 scalefactor = fV0MScaleFactor[runnumber - fLowRunN]; // subtracting reference offset index
872 scalefactor = fSPDScaleFactor[runnumber - fLowRunN]; // subtracting reference offset index
874 scalefactor = fTPCScaleFactor[runnumber - fLowRunN]; // subtracting reference offset index
880 //________________________________________________________________________
881 Float_t AliCentralitySelectionTask::MyGetScaleFactorMC(Int_t runnumber) const
883 // Get MC scaling factor
884 if (! (runnumber >= fLowRunN && runnumber <=fHighRunN)) {
885 cout << "MyGetScaleFactor error in run number range " << runnumber << endl;
889 Float_t scalefactor= fV0MScaleFactorMC[runnumber - fLowRunN]; // subtracting reference offset index
894 //________________________________________________________________________
895 void AliCentralitySelectionTask::MyInitScaleFactor ()
897 // Initialize the scaling factors
898 for (int i=0; i<=(fHighRunN-fLowRunN); i++) fV0MScaleFactor[i] = 0.0;
899 for (int i=0; i<=(fHighRunN-fLowRunN); i++) fSPDScaleFactor[i] = 0.0;
900 for (int i=0; i<=(fHighRunN-fLowRunN); i++) fTPCScaleFactor[i] = 0.0;
902 // scale factors determined from <V0 charge> on a run-by-run basis
903 fV0MScaleFactor[310] = 1.;
904 fV0MScaleFactor[311] = 1.;
905 fV0MScaleFactor[514] = 1.0046;
906 fV0MScaleFactor[515] = 0.983535;
907 fV0MScaleFactor[579] = 0.988185;
908 fV0MScaleFactor[580] = 0.983351;
909 fV0MScaleFactor[581] = 0.989013;
910 fV0MScaleFactor[583] = 0.990056;
911 fV0MScaleFactor[588] = 0.974438;
912 fV0MScaleFactor[589] = 0.981572;
913 fV0MScaleFactor[590] = 0.989316;
914 fV0MScaleFactor[592] = 0.98345;
915 fV0MScaleFactor[688] = 0.993647;
916 fV0MScaleFactor[690] = 0.994758;
917 fV0MScaleFactor[693] = 0.989569;
918 fV0MScaleFactor[698] = 0.993119;
919 fV0MScaleFactor[744] = 0.989583;
920 fV0MScaleFactor[757] = 0.990377;
921 fV0MScaleFactor[787] = 0.990176;
922 fV0MScaleFactor[788] = 0.98723;
923 fV0MScaleFactor[834] = 1.00403;
924 fV0MScaleFactor[835] = 0.994376;
925 fV0MScaleFactor[840] = 0.99759;
926 fV0MScaleFactor[842] = 1.01031;
927 fV0MScaleFactor[900] = 0.996216;
928 fV0MScaleFactor[901] = 0.994205;
929 fV0MScaleFactor[992] = 0.998479;
930 fV0MScaleFactor[997] = 1.00078;
931 fV0MScaleFactor[1299] = 1.00515;
932 fV0MScaleFactor[1303] = 1.00094;
933 fV0MScaleFactor[1339] = 0.986596;
934 fV0MScaleFactor[1346] = 0.972226;
935 fV0MScaleFactor[1349] = 0.960358;
936 fV0MScaleFactor[1350] = 0.970023;
937 fV0MScaleFactor[1374] = 1.00575;
938 fV0MScaleFactor[1545] = 1.00471;
939 fV0MScaleFactor[1587] = 1.00611;
940 fV0MScaleFactor[1588] = 1.00976;
941 fV0MScaleFactor[1618] = 1.00771;
942 fV0MScaleFactor[1682] = 1.01622;
943 fV0MScaleFactor[1727] = 1.01305;
944 fV0MScaleFactor[1728] = 1.00388;
945 fV0MScaleFactor[1731] = 1.00673;
946 fV0MScaleFactor[1732] = 1.00916;
947 fV0MScaleFactor[1770] = 1.0092;
948 fV0MScaleFactor[1773] = 1.00728;
949 fV0MScaleFactor[1786] = 1.01655;
950 fV0MScaleFactor[1787] = 1.00672;
951 fV0MScaleFactor[1801] = 0.983339;
952 fV0MScaleFactor[1802] = 1.00754;
953 fV0MScaleFactor[1811] = 1.00608;
954 fV0MScaleFactor[1815] = 1.01227;
955 fV0MScaleFactor[1879] = 0.99907;
956 fV0MScaleFactor[1881] = 0.995696;
957 fV0MScaleFactor[2186] = 1.00559;
958 fV0MScaleFactor[2187] = 1.00631;
959 fV0MScaleFactor[2254] = 1.01512;
960 fV0MScaleFactor[2256] = 0.998727;
961 fV0MScaleFactor[2321] = 1.00701;
962 fSPDScaleFactor[310] = 1.00211;
963 fSPDScaleFactor[311] = 1.00067;
964 fSPDScaleFactor[514] = 1.02268;
965 fSPDScaleFactor[515] = 0.994902;
966 fSPDScaleFactor[579] = 1.00215;
967 fSPDScaleFactor[580] = 0.993421;
968 fSPDScaleFactor[581] = 1.00129;
969 fSPDScaleFactor[583] = 1.00242;
970 fSPDScaleFactor[588] = 0.984762;
971 fSPDScaleFactor[589] = 0.994355;
972 fSPDScaleFactor[590] = 1.00073;
973 fSPDScaleFactor[592] = 0.995889;
974 fSPDScaleFactor[688] = 0.994532;
975 fSPDScaleFactor[690] = 0.998307;
976 fSPDScaleFactor[693] = 0.994052;
977 fSPDScaleFactor[698] = 0.993224;
978 fSPDScaleFactor[744] = 0.993279;
979 fSPDScaleFactor[757] = 0.992494;
980 fSPDScaleFactor[787] = 0.992678;
981 fSPDScaleFactor[788] = 0.996563;
982 fSPDScaleFactor[834] = 1.01116;
983 fSPDScaleFactor[835] = 0.993108;
984 fSPDScaleFactor[840] = 0.997574;
985 fSPDScaleFactor[842] = 1.01829;
986 fSPDScaleFactor[900] = 0.999438;
987 fSPDScaleFactor[901] = 0.995849;
988 fSPDScaleFactor[992] = 0.999227;
989 fSPDScaleFactor[997] = 1.00575;
990 fSPDScaleFactor[1299] = 0.99877;
991 fSPDScaleFactor[1303] = 0.999682;
992 fSPDScaleFactor[1339] = 0.978198;
993 fSPDScaleFactor[1346] = 0.964178;
994 fSPDScaleFactor[1349] = 0.959439;
995 fSPDScaleFactor[1350] = 0.956945;
996 fSPDScaleFactor[1374] = 0.994434;
997 fSPDScaleFactor[1545] = 1.0016;
998 fSPDScaleFactor[1587] = 1.00153;
999 fSPDScaleFactor[1588] = 1.00698;
1000 fSPDScaleFactor[1618] = 1.00554;
1001 fSPDScaleFactor[1682] = 1.0123;
1002 fSPDScaleFactor[1727] = 1.011;
1003 fSPDScaleFactor[1728] = 1.00143;
1004 fSPDScaleFactor[1731] = 1.00486;
1005 fSPDScaleFactor[1732] = 1.00646;
1006 fSPDScaleFactor[1770] = 1.00515;
1007 fSPDScaleFactor[1773] = 1.00485;
1008 fSPDScaleFactor[1786] = 1.01215;
1009 fSPDScaleFactor[1787] = 1.00419;
1010 fSPDScaleFactor[1801] = 0.983327;
1011 fSPDScaleFactor[1802] = 1.00529;
1012 fSPDScaleFactor[1811] = 1.00367;
1013 fSPDScaleFactor[1815] = 1.01045;
1014 fSPDScaleFactor[1879] = 0.996374;
1015 fSPDScaleFactor[1881] = 0.988827;
1016 fSPDScaleFactor[2186] = 1.00354;
1017 fSPDScaleFactor[2187] = 1.00397;
1018 fSPDScaleFactor[2254] = 1.01138;
1019 fSPDScaleFactor[2256] = 0.996641;
1020 fSPDScaleFactor[2321] = 1.00357;
1021 fTPCScaleFactor[310] = 1.00434;
1022 fTPCScaleFactor[311] = 1.0056;
1023 fTPCScaleFactor[514] = 1.02185;
1024 fTPCScaleFactor[515] = 1.0036;
1025 fTPCScaleFactor[579] = 1.00607;
1026 fTPCScaleFactor[580] = 1.00183;
1027 fTPCScaleFactor[581] = 1.00693;
1028 fTPCScaleFactor[583] = 1.00746;
1029 fTPCScaleFactor[588] = 0.990524;
1030 fTPCScaleFactor[589] = 0.998582;
1031 fTPCScaleFactor[590] = 1.00618;
1032 fTPCScaleFactor[592] = 1.00088;
1033 fTPCScaleFactor[688] = 1.00598;
1034 fTPCScaleFactor[690] = 1.00658;
1035 fTPCScaleFactor[693] = 1.00144;
1036 fTPCScaleFactor[698] = 1.00279;
1037 fTPCScaleFactor[744] = 1.00122;
1038 fTPCScaleFactor[757] = 1.002;
1039 fTPCScaleFactor[787] = 0.997818;
1040 fTPCScaleFactor[788] = 0.994583;
1041 fTPCScaleFactor[834] = 1.01508;
1042 fTPCScaleFactor[835] = 1.00218;
1043 fTPCScaleFactor[840] = 1.00569;
1044 fTPCScaleFactor[842] = 1.01789;
1045 fTPCScaleFactor[900] = 1.00739;
1046 fTPCScaleFactor[901] = 1.00462;
1047 fTPCScaleFactor[992] = 1.00502;
1048 fTPCScaleFactor[997] = 1.00943;
1049 fTPCScaleFactor[1299] = 1.00438;
1050 fTPCScaleFactor[1303] = 0.996701;
1051 fTPCScaleFactor[1339] = 0.978641;
1052 fTPCScaleFactor[1346] = 0.968906;
1053 fTPCScaleFactor[1349] = 0.954311;
1054 fTPCScaleFactor[1350] = 0.958764;
1055 fTPCScaleFactor[1374] = 0.997899;
1056 fTPCScaleFactor[1545] = 0.992;
1057 fTPCScaleFactor[1587] = 0.992635;
1058 fTPCScaleFactor[1588] = 1.00207;
1059 fTPCScaleFactor[1618] = 1.00126;
1060 fTPCScaleFactor[1682] = 1.00324;
1061 fTPCScaleFactor[1727] = 1.00042;
1062 fTPCScaleFactor[1728] = 0.978881;
1063 fTPCScaleFactor[1731] = 0.999818;
1064 fTPCScaleFactor[1732] = 1.00109;
1065 fTPCScaleFactor[1770] = 0.99935;
1066 fTPCScaleFactor[1773] = 0.998531;
1067 fTPCScaleFactor[1786] = 0.999125;
1068 fTPCScaleFactor[1787] = 0.998479;
1069 fTPCScaleFactor[1801] = 0.9775;
1070 fTPCScaleFactor[1802] = 0.999095;
1071 fTPCScaleFactor[1811] = 0.998197;
1072 fTPCScaleFactor[1815] = 1.00413;
1073 fTPCScaleFactor[1879] = 0.990916;
1074 fTPCScaleFactor[1881] = 0.987241;
1075 fTPCScaleFactor[2186] = 1.00048;
1076 fTPCScaleFactor[2187] = 1.00057;
1077 fTPCScaleFactor[2254] = 1.00588;
1078 fTPCScaleFactor[2256] = 0.991503;
1079 fTPCScaleFactor[2321] = 1.00057;
1082 // set all missing values to the value of the run before it ....
1083 for (int i=0; i<=(fHighRunN-fLowRunN); i++) {
1084 if (fV0MScaleFactor[i] == 0.0) {
1086 fV0MScaleFactor[i] = 1.00;
1088 // search for last run number with non-zero value
1089 for (int j=i-1;j>=0;j--) {
1090 if (fV0MScaleFactor[j] != 0.0) {
1091 fV0MScaleFactor[i] = fV0MScaleFactor[j];
1097 } // end loop over checking all run-numbers
1099 for (int i=0; i<=(fHighRunN-fLowRunN); i++) {
1100 if (fSPDScaleFactor[i] == 0.0) {
1102 fSPDScaleFactor[i] = 1.00;
1104 for (int j=i-1;j>=0;j--) {
1105 if (fSPDScaleFactor[j] != 0.0) {
1106 fSPDScaleFactor[i] = fSPDScaleFactor[j];
1114 for (int i=0; i<=(fHighRunN-fLowRunN); i++) {
1115 if (fTPCScaleFactor[i] == 0.0) {
1117 fTPCScaleFactor[i] = 1.00;
1119 for (int j=i-1;j>=0;j--) {
1120 if (fTPCScaleFactor[j] != 0.0) {
1121 fTPCScaleFactor[i] = fTPCScaleFactor[j];
1130 // for (int i=0; i<=(fHighRunN-fLowRunN); i++) cout << "Scale Factor = " << fV0MScaleFactor[i] << " for Run " << i+fLowRunN << endl;
1137 //________________________________________________________________________
1138 void AliCentralitySelectionTask::MyInitScaleFactorMC()
1140 // Initialize the MC scaling factors
1141 for (int i=0; i<(fHighRunN-fLowRunN); i++) fV0MScaleFactorMC[i] = 0.0;
1142 // scale factors determined from <V0 charge> on a run-by-run basis
1143 fV0MScaleFactorMC[0] = 0.75108;
1144 // set all missing values to the value of the run before it ....
1145 for (int i=0; i<=(fHighRunN-fLowRunN); i++) {
1146 if (fV0MScaleFactorMC[i] == 0.0) {
1148 fV0MScaleFactorMC[i] = 1.00;
1150 // search for last run number with non-zero value
1151 for (int j=i-1;j>=0;j--) {
1152 if (fV0MScaleFactorMC[j] != 0.0) {
1153 fV0MScaleFactorMC[i] = fV0MScaleFactorMC[j];
1159 } // end loop over checking all run-numbers
1162 // for (int i=0; i<=(fHighRunN-fLowRunN); i++) cout << "Scale Factor = " << fV0MScaleFactorMC[i] << " for Run " << i+fLowRunN << endl;