1 /**************************************************************************
2 * Copyright(c) 1998-1999, 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 **************************************************************************/
15 //------------------------------------------------------------------------------
16 // AlidNdPtAnalysisPbPb class.
19 // - fills analysis control histograms
20 // - fills generic correction matrices
21 // - generates correction matrices
24 // - generic correction matrices
25 // - control histograms
27 // Author: J.Otwinowski 04/11/2008
28 // last change: 2011-04-04 by M.Knichel
29 //------------------------------------------------------------------------------
34 #include "THnSparse.h"
36 #include "AliHeader.h"
37 #include "AliGenEventHeader.h"
38 #include "AliInputEventHandler.h"
39 #include "AliAnalysisManager.h"
41 #include "AliESDEvent.h"
42 #include "AliMCEvent.h"
43 #include "AliESDtrackCuts.h"
45 #include "AliMultiplicity.h"
46 #include "AliTracker.h"
48 #include "AliCentrality.h"
50 #include "AlidNdPtEventCuts.h"
51 #include "AlidNdPtAcceptanceCuts.h"
52 #include "AliPhysicsSelection.h"
53 #include "AliTriggerAnalysis.h"
55 #include "AliPWG0Helper.h"
56 #include "AlidNdPtHelper.h"
57 #include "AlidNdPtAnalysisPbPb.h"
62 ClassImp(AlidNdPtAnalysisPbPb)
64 //_____________________________________________________________________________
65 AlidNdPtAnalysisPbPb::AlidNdPtAnalysisPbPb(): AlidNdPt(),
67 fHistogramsOn(kFALSE),
69 // rec. track pt vs true track pt correlation matrix
70 fTrackPtCorrelationMatrix(0),
72 // event level correction
74 fTriggerEventMatrix(0),
78 // track-event level correction
80 fGenTrackEventMatrix(0),
82 fTriggerTrackEventMatrix(0),
84 fRecTrackEventMatrix(0),
86 // track rec. efficiency correction (fRecPrimTrackMatrix / fGenPrimTrackMatrix)
88 fGenPrimTrackMatrix(0),
89 fRecPrimTrackMatrix(0),
91 // secondary track contamination correction (fRecSecTrackMatrix / fRecTrackMatrix)
93 fRecSecTrackMatrix(0),
95 // multiple rec. track contamination corrections (fRecMultTrackMatrix / fRecTrackMatrix)
96 fRecMultTrackMatrix(0),
98 // event control histograms
105 // rec. pt and eta resolution w.r.t MC
108 //multple reconstructed tracks
109 fMCMultRecTrackHist1(0),
111 // rec. track control histograms
116 fCentralityEstimator(0),
133 // default constructor
134 for(Int_t i=0; i<AlidNdPtHelper::kCutSteps; i++) {
136 fMCPrimTrackHist1[i]=0;
137 fMCPrimTrackHist2[i]=0;
138 fMCSecTrackHist1[i]=0;
141 fRecTrackMultHist1[i]=0;
144 SetCentralityEstimator();
147 //_____________________________________________________________________________
148 AlidNdPtAnalysisPbPb::AlidNdPtAnalysisPbPb(Char_t* name, Char_t* title): AlidNdPt(name,title),
150 fHistogramsOn(kFALSE),
152 // rec. track pt vs true track pt correlation matrix
153 fTrackPtCorrelationMatrix(0),
155 // event level correction
157 fTriggerEventMatrix(0),
161 // track-event level correction
163 fGenTrackEventMatrix(0),
165 fTriggerTrackEventMatrix(0),
167 fRecTrackEventMatrix(0),
169 // track rec. efficiency correction (fRecPrimTrackMatrix / fGenPrimTrackMatrix)
171 fGenPrimTrackMatrix(0),
172 fRecPrimTrackMatrix(0),
174 // secondary track contamination correction (fRecTrackMatrix - fRecSecTrackMatrix)
176 fRecSecTrackMatrix(0),
178 // multiple rec. track contamination corrections (fRecTrackMatrix - fRecMultTrackMatrix)
179 fRecMultTrackMatrix(0),
184 // event control histograms
191 // rec. pt and eta resolution w.r.t MC
194 //multiple reconstructed tracks
195 fMCMultRecTrackHist1(0),
197 // rec. track control histograms
202 fCentralityEstimator(0),
223 for(Int_t i=0; i<AlidNdPtHelper::kCutSteps; i++) {
225 fMCPrimTrackHist1[i]=0;
226 fMCPrimTrackHist2[i]=0;
227 fMCSecTrackHist1[i]=0;
230 fRecTrackMultHist1[i]=0;
234 SetCentralityEstimator();
237 //_____________________________________________________________________________
238 AlidNdPtAnalysisPbPb::~AlidNdPtAnalysisPbPb() {
242 if(fTrackPtCorrelationMatrix) delete fTrackPtCorrelationMatrix; fTrackPtCorrelationMatrix=0;
244 if(fGenEventMatrix) delete fGenEventMatrix; fGenEventMatrix=0;
245 if(fTriggerEventMatrix) delete fTriggerEventMatrix; fTriggerEventMatrix=0;
246 if(fRecEventMatrix) delete fRecEventMatrix; fRecEventMatrix=0;
249 if(fGenTrackEventMatrix) delete fGenTrackEventMatrix; fGenTrackEventMatrix=0;
250 if(fTriggerEventMatrix) delete fTriggerEventMatrix; fTriggerEventMatrix=0;
251 if(fRecTrackEventMatrix) delete fRecTrackEventMatrix; fRecTrackEventMatrix=0;
254 if(fGenTrackMatrix) delete fGenTrackMatrix; fGenTrackMatrix=0;
255 if(fGenPrimTrackMatrix) delete fGenPrimTrackMatrix; fGenPrimTrackMatrix=0;
256 if(fRecPrimTrackMatrix) delete fRecPrimTrackMatrix; fRecPrimTrackMatrix=0;
258 if(fRecTrackMatrix) delete fRecTrackMatrix; fRecTrackMatrix=0;
259 if(fRecSecTrackMatrix) delete fRecSecTrackMatrix; fRecSecTrackMatrix=0;
261 if(fRecMultTrackMatrix) delete fRecMultTrackMatrix; fRecMultTrackMatrix=0;
263 // Control histograms
265 if(fMCEventHist1) delete fMCEventHist1; fMCEventHist1=0;
266 if(fRecEventHist1) delete fRecEventHist1; fRecEventHist1=0;
267 if(fRecEventHist2) delete fRecEventHist2; fRecEventHist2=0;
268 if(fRecMCEventHist1) delete fRecMCEventHist1; fRecMCEventHist1=0;
269 if(fRecMCEventHist2) delete fRecMCEventHist2; fRecMCEventHist2=0;
272 for(Int_t i=0; i<AlidNdPtHelper::kCutSteps; i++) {
273 if(fMCTrackHist1[i]) delete fMCTrackHist1[i]; fMCTrackHist1[i]=0;
274 if(fMCPrimTrackHist1[i]) delete fMCPrimTrackHist1[i]; fMCPrimTrackHist1[i]=0;
275 if(fMCPrimTrackHist2[i]) delete fMCPrimTrackHist2[i]; fMCPrimTrackHist2[i]=0;
276 if(fMCSecTrackHist1[i]) delete fMCSecTrackHist1[i]; fMCSecTrackHist1[i]=0;
277 if(fRecTrackHist1[i]) delete fRecTrackHist1[i]; fRecTrackHist1[i]=0;
278 if(fRecTrackHist2[i]) delete fRecTrackHist2[i]; fRecTrackHist2[i]=0;
279 if(fRecTrackMultHist1[i]) delete fRecTrackMultHist1[i]; fRecTrackMultHist1[i]=0;
281 if(fRecMCTrackHist1) delete fRecMCTrackHist1; fRecMCTrackHist1=0;
282 if(fMCMultRecTrackHist1) delete fMCMultRecTrackHist1; fMCMultRecTrackHist1=0;
283 if(fRecTrackHist3) delete fRecTrackHist3; fRecTrackHist3=0;
285 if(fAnalysisFolder) delete fAnalysisFolder; fAnalysisFolder=0;
287 if (fTriggerAnalysis) delete fTriggerAnalysis; fTriggerAnalysis = 0;
289 if (fBinsMult) delete[] fBinsMult; fBinsMult=0;
290 if (fBinsPt) delete[] fBinsPt; fBinsPt=0;
291 if (fBinsPtCorr) delete[] fBinsPtCorr; fBinsPtCorr=0;
292 if (fBinsEta) delete[] fBinsEta; fBinsEta=0;
293 if (fBinsZv) delete[] fBinsZv; fBinsZv=0;
294 if (fBinsCentrality) delete[] fBinsCentrality; fBinsCentrality=0;
297 //_____________________________________________________________________________
298 void AlidNdPtAnalysisPbPb::Init() {
300 //define default binning
301 Double_t binsMultDefault[48] = {-0.5, 0.5 , 1.5 , 2.5 , 3.5 , 4.5 , 5.5 , 6.5 , 7.5 , 8.5,9.5, 10.5, 11.5, 12.5, 13.5, 14.5, 15.5, 16.5, 17.5, 18.5,19.5, 20.5, 30.5, 40.5 , 50.5 , 60.5 , 70.5 , 80.5 , 90.5 , 100.5,200.5, 300.5, 400.5, 500.5, 600.5, 700.5, 800.5, 900.5, 1000.5, 2000.5, 3000.5, 4000.5, 5000.5, 6000.5, 7000.5, 8000.5, 9000.5, 10000.5 };
302 Double_t binsPtDefault[69] = {0.,0.05,0.1,0.15,0.2,0.25,0.3,0.35,0.4,0.45,0.5,0.55,0.6,0.65,0.7,0.75,0.8,0.85,0.9,0.95,1.0,1.1,1.2,1.3,1.4,1.5,1.6,1.7,1.8,1.9,2.0,2.2,2.4,2.6,2.8,3.0,3.2,3.4,3.6,3.8,4.0,4.5,5.0,5.5,6.0,6.5,7.0,8.0,9.0,10.0,11.0,12.0,13.0,14.0,15.0,16.0,18.0,20.0,22.0,24.0,26.0,28.0,30.0,32.0,34.0,36.0,40.0,45.0,50.0};
303 Double_t binsPtCorrDefault[37] = {0.,0.05,0.1,0.15,0.2,0.25,0.3,0.35,0.4,0.45,0.5,0.55,0.6,0.65,0.7,0.75,0.8,0.85,0.9,0.95,1.0,1.1,1.2,1.3,1.4,1.5,1.6,1.7,1.8,1.9,2.0,2.2,2.4,2.6,3.0,4.0,50.0};
304 Double_t binsEtaDefault[31] = {-1.5,-1.4,-1.3,-1.2,-1.1,-1.0,-0.9,-0.8,-0.7,-0.6,-0.5,-0.4,-0.3,-0.2,-0.1,0.,0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9,1.0,1.1,1.2,1.3,1.4,1.5};
305 Double_t binsZvDefault[13] = {-30.,-25.,-20.,-15.,-10.,-5.,0.,5.,10.,15.,20.,25.,30.};
306 Double_t binsCentralityDefault[12] = {0., 5., 10., 20., 30., 40., 50., 60., 70., 80., 90., 100.};
308 // if no binning is set, use the default
309 if (!fBinsMult) { SetBinsMult(47,binsMultDefault); }
310 if (!fBinsPt) { SetBinsPt(68,binsPtDefault); }
311 if (!fBinsPtCorr) { SetBinsPtCorr(36,binsPtCorrDefault); }
312 if (!fBinsEta) { SetBinsEta(30,binsEtaDefault); }
313 if (!fBinsZv) { SetBinsZv(12,binsZvDefault); }
314 if (!fBinsCentrality) { SetBinsCentrality(11,binsCentralityDefault); }
316 Int_t binsTrackEventCorrMatrix[4]={fZvNbins,fPtNbins,fEtaNbins,fCentralityNbins};
317 Int_t binsTrackEvent[4]={fZvNbins,fPtNbins,fEtaNbins,fCentralityNbins};
318 Int_t binsTrackPtCorrelationMatrix[4]={fPtNbins,fPtNbins,fEtaNbins,fCentralityNbins};
320 fTrackPtCorrelationMatrix = new THnSparseF("fTrackPtCorrelationMatrix","Pt:mcPt:mcEta:Centrality",4,binsTrackPtCorrelationMatrix);
321 fTrackPtCorrelationMatrix->SetBinEdges(0,fBinsPt);
322 fTrackPtCorrelationMatrix->SetBinEdges(1,fBinsPt);
323 fTrackPtCorrelationMatrix->SetBinEdges(2,fBinsEta);
324 fTrackPtCorrelationMatrix->SetBinEdges(3,fBinsCentrality);
325 fTrackPtCorrelationMatrix->GetAxis(0)->SetTitle("Pt (GeV/c)");
326 fTrackPtCorrelationMatrix->GetAxis(1)->SetTitle("mcPt (GeV/c)");
327 fTrackPtCorrelationMatrix->GetAxis(2)->SetTitle("mcEta");
328 fTrackPtCorrelationMatrix->GetAxis(3)->SetTitle("Centrality");
329 fTrackPtCorrelationMatrix->Sumw2();
332 // Efficiency and contamination correction matrices
334 Int_t binsEventMatrix[3]={fZvNbins,fMultNbins,fCentralityNbins};
335 // Double_t minEventMatrix[3]={-30.,-0.5,0.};
336 // Double_t maxEventMatrix[3]={30.,10000.5,100.};
338 fGenEventMatrix = new THnSparseF("fGenEventMatrix","mcZv:multMB:Centrality",3,binsEventMatrix); //,minEventMatrix,maxEventMatrix);
339 fGenEventMatrix->SetBinEdges(0,fBinsZv);
340 fGenEventMatrix->SetBinEdges(1,fBinsMult);
341 fGenEventMatrix->SetBinEdges(2,fBinsCentrality);
342 fGenEventMatrix->GetAxis(0)->SetTitle("mcZv (cm)");
343 fGenEventMatrix->GetAxis(1)->SetTitle("multiplicity MB");
344 fGenEventMatrix->GetAxis(2)->SetTitle("Centrality");
345 fGenEventMatrix->Sumw2();
347 fTriggerEventMatrix = new THnSparseF("fTriggerEventMatrix","mcZv:multMB:Centrality",3,binsEventMatrix); //,minEventMatrix,maxEventMatrix);
348 fTriggerEventMatrix->SetBinEdges(0,fBinsZv);
349 fTriggerEventMatrix->SetBinEdges(1,fBinsMult);
350 fTriggerEventMatrix->SetBinEdges(2,fBinsCentrality);
351 fTriggerEventMatrix->GetAxis(0)->SetTitle("mcZv (cm)");
352 fTriggerEventMatrix->GetAxis(1)->SetTitle("multiplicity MB");
353 fTriggerEventMatrix->GetAxis(2)->SetTitle("Centrality");
354 fTriggerEventMatrix->Sumw2();
356 fRecEventMatrix = new THnSparseF("fRecEventMatrix","mcZv:multMB:Centrality",3,binsEventMatrix); //,minEventMatrix,maxEventMatrix);
357 fRecEventMatrix->SetBinEdges(0,fBinsZv);
358 fRecEventMatrix->SetBinEdges(1,fBinsMult);
359 fRecEventMatrix->SetBinEdges(2,fBinsCentrality);
360 fRecEventMatrix->GetAxis(0)->SetTitle("mcZv (cm)");
361 fRecEventMatrix->GetAxis(1)->SetTitle("multiplicity MB");
362 fRecEventMatrix->GetAxis(2)->SetTitle("Centrality");
363 fRecEventMatrix->Sumw2();
366 // track to event corrections
368 fGenTrackEventMatrix = new THnSparseF("fGenTrackEventMatrix","mcZv:mcPt:mcEta:Centrality",4,binsTrackEventCorrMatrix);
369 fGenTrackEventMatrix->SetBinEdges(0,fBinsZv);
370 fGenTrackEventMatrix->SetBinEdges(1,fBinsPtCorr);
371 fGenTrackEventMatrix->SetBinEdges(2,fBinsEta);
372 fGenTrackEventMatrix->SetBinEdges(3,fBinsCentrality);
373 fGenTrackEventMatrix->GetAxis(0)->SetTitle("mcZv (cm)");
374 fGenTrackEventMatrix->GetAxis(1)->SetTitle("mcPt (GeV/c)");
375 fGenTrackEventMatrix->GetAxis(2)->SetTitle("mcEta");
376 fGenTrackEventMatrix->GetAxis(3)->SetTitle("Centrality");
377 fGenTrackEventMatrix->Sumw2();
379 fTriggerTrackEventMatrix = new THnSparseF("fTriggerTrackEventMatrix","mcZv:mcPt:mcEta:Centrality",4,binsTrackEventCorrMatrix);
380 fTriggerTrackEventMatrix->SetBinEdges(0,fBinsZv);
381 fTriggerTrackEventMatrix->SetBinEdges(1,fBinsPtCorr);
382 fTriggerTrackEventMatrix->SetBinEdges(2,fBinsEta);
383 fTriggerTrackEventMatrix->SetBinEdges(3,fBinsCentrality);
384 fTriggerTrackEventMatrix->GetAxis(0)->SetTitle("mcZv (cm)");
385 fTriggerTrackEventMatrix->GetAxis(1)->SetTitle("mcPt (GeV/c)");
386 fTriggerTrackEventMatrix->GetAxis(2)->SetTitle("mcEta");
387 fTriggerTrackEventMatrix->GetAxis(3)->SetTitle("Centrality");
388 fTriggerTrackEventMatrix->Sumw2();
390 fRecTrackEventMatrix = new THnSparseF("fRecTrackEventMatrix","mcZv:mcPt:mcEta:Centrality",4,binsTrackEventCorrMatrix);
391 fRecTrackEventMatrix->SetBinEdges(0,fBinsZv);
392 fRecTrackEventMatrix->SetBinEdges(1,fBinsPtCorr);
393 fRecTrackEventMatrix->SetBinEdges(2,fBinsEta);
394 fRecTrackEventMatrix->SetBinEdges(3,fBinsCentrality);
395 fRecTrackEventMatrix->GetAxis(0)->SetTitle("mcZv (cm)");
396 fRecTrackEventMatrix->GetAxis(1)->SetTitle("mcPt (GeV/c)");
397 fRecTrackEventMatrix->GetAxis(2)->SetTitle("mcEta");
398 fRecTrackEventMatrix->GetAxis(3)->SetTitle("Centrality");
399 fRecTrackEventMatrix->Sumw2();
402 // tracks correction matrices
404 fGenTrackMatrix = new THnSparseF("fGenTrackMatrix","mcZv:mcPt:mcEta:Centrality",4,binsTrackEventCorrMatrix);
405 fGenTrackMatrix->SetBinEdges(0,fBinsZv);
406 fGenTrackMatrix->SetBinEdges(1,fBinsPtCorr);
407 fGenTrackMatrix->SetBinEdges(2,fBinsEta);
408 fGenTrackMatrix->SetBinEdges(3,fBinsCentrality);
409 fGenTrackMatrix->GetAxis(0)->SetTitle("mcZv (cm)");
410 fGenTrackMatrix->GetAxis(1)->SetTitle("mcPt (GeV/c)");
411 fGenTrackMatrix->GetAxis(2)->SetTitle("mcEta");
412 fGenTrackMatrix->GetAxis(3)->SetTitle("Centrality");
413 fGenTrackMatrix->Sumw2();
415 fGenPrimTrackMatrix = new THnSparseF("fGenPrimTrackMatrix","mcZv:mcPt:mcEta:Centrality",4,binsTrackEventCorrMatrix);
416 fGenPrimTrackMatrix->SetBinEdges(0,fBinsZv);
417 fGenPrimTrackMatrix->SetBinEdges(1,fBinsPtCorr);
418 fGenPrimTrackMatrix->SetBinEdges(2,fBinsEta);
419 fGenPrimTrackMatrix->SetBinEdges(3,fBinsCentrality);
420 fGenPrimTrackMatrix->GetAxis(0)->SetTitle("mcZv (cm)");
421 fGenPrimTrackMatrix->GetAxis(1)->SetTitle("mcPt (GeV/c)");
422 fGenPrimTrackMatrix->GetAxis(2)->SetTitle("mcEta");
423 fGenPrimTrackMatrix->GetAxis(3)->SetTitle("Centrality");
424 fGenPrimTrackMatrix->Sumw2();
427 fRecPrimTrackMatrix = new THnSparseF("fRecPrimTrackMatrix","mcZv:mcPt:mcEta:Centrality",4,binsTrackEventCorrMatrix);
428 fRecPrimTrackMatrix->SetBinEdges(0,fBinsZv);
429 fRecPrimTrackMatrix->SetBinEdges(1,fBinsPtCorr);
430 fRecPrimTrackMatrix->SetBinEdges(2,fBinsEta);
431 fRecPrimTrackMatrix->SetBinEdges(3,fBinsCentrality);
432 fRecPrimTrackMatrix->GetAxis(0)->SetTitle("mcZv (cm)");
433 fRecPrimTrackMatrix->GetAxis(1)->SetTitle("mcPt (GeV/c)");
434 fRecPrimTrackMatrix->GetAxis(2)->SetTitle("mcEta");
435 fRecPrimTrackMatrix->GetAxis(3)->SetTitle("Centrality");
436 fRecPrimTrackMatrix->Sumw2();
439 fRecTrackMatrix = new THnSparseF("fRecTrackMatrix","mcZv:mcPt:mcEta:Centrality",4,binsTrackEventCorrMatrix);
440 fRecTrackMatrix->SetBinEdges(0,fBinsZv);
441 fRecTrackMatrix->SetBinEdges(1,fBinsPtCorr);
442 fRecTrackMatrix->SetBinEdges(2,fBinsEta);
443 fRecTrackMatrix->SetBinEdges(3,fBinsCentrality);
444 fRecTrackMatrix->GetAxis(0)->SetTitle("mcZv (cm)");
445 fRecTrackMatrix->GetAxis(1)->SetTitle("mcPt (GeV/c)");
446 fRecTrackMatrix->GetAxis(2)->SetTitle("mcEta");
447 fRecTrackMatrix->GetAxis(3)->SetTitle("Centrality");
448 fRecTrackMatrix->Sumw2();
450 fRecSecTrackMatrix = new THnSparseF("fRecSecTrackMatrix","mcZv:mcPt:mcEta:Centrality",4,binsTrackEventCorrMatrix);
451 fRecSecTrackMatrix->SetBinEdges(0,fBinsZv);
452 fRecSecTrackMatrix->SetBinEdges(1,fBinsPtCorr);
453 fRecSecTrackMatrix->SetBinEdges(2,fBinsEta);
454 fRecSecTrackMatrix->SetBinEdges(3,fBinsCentrality);
455 fRecSecTrackMatrix->GetAxis(0)->SetTitle("mcZv (cm)");
456 fRecSecTrackMatrix->GetAxis(1)->SetTitle("Pt (GeV/c)");
457 fRecSecTrackMatrix->GetAxis(2)->SetTitle("Eta");
458 fRecSecTrackMatrix->GetAxis(3)->SetTitle("Centrality");
459 fRecSecTrackMatrix->Sumw2();
462 fRecMultTrackMatrix = new THnSparseF("fRecMultTrackMatrix","mcZv:mcPt:mcEta:Centrality",4,binsTrackEventCorrMatrix);
463 fRecMultTrackMatrix->SetBinEdges(0,fBinsZv);
464 fRecMultTrackMatrix->SetBinEdges(1,fBinsPtCorr);
465 fRecMultTrackMatrix->SetBinEdges(2,fBinsEta);
466 fRecMultTrackMatrix->SetBinEdges(3,fBinsCentrality);
467 fRecMultTrackMatrix->GetAxis(0)->SetTitle("mcZv (cm)");
468 fRecMultTrackMatrix->GetAxis(1)->SetTitle("mcPt (GeV/c)");
469 fRecMultTrackMatrix->GetAxis(2)->SetTitle("mcEta");
470 fRecMultTrackMatrix->GetAxis(3)->SetTitle("Centrality");
471 fRecMultTrackMatrix->Sumw2();
474 // Control analysis histograms
476 Int_t binsMCEventHist1[4]={100,100,fZvNbins,fCentralityNbins};
477 Double_t minMCEventHist1[4]={-0.1,-0.1,-30.,0.};
478 Double_t maxMCEventHist1[4]={0.1,0.1,30.,100.};
479 fMCEventHist1 = new THnSparseF("fMCEventHist1","mcXv:mcYv:mcZv:Centrality",4,binsMCEventHist1,minMCEventHist1,maxMCEventHist1);
480 fMCEventHist1->SetBinEdges(2,fBinsZv);
481 fMCEventHist1->SetBinEdges(3,fBinsCentrality);
482 fMCEventHist1->GetAxis(0)->SetTitle("mcXv (cm)");
483 fMCEventHist1->GetAxis(1)->SetTitle("mcYv (cm)");
484 fMCEventHist1->GetAxis(2)->SetTitle("mcZv (cm)");
485 fMCEventHist1->GetAxis(3)->SetTitle("Centrality");
486 fMCEventHist1->Sumw2();
489 Int_t binsRecEventHist1[4]={100,100,fZvNbins,fCentralityNbins};
490 Double_t minRecEventHist1[4]={-3.,-3.,-30.,0.};
491 Double_t maxRecEventHist1[4]={3.,3.,30.,100.};
493 fRecEventHist1 = new THnSparseF("fRecEventHist1","Xv:Yv:Zv:Centrality",4,binsRecEventHist1,minRecEventHist1,maxRecEventHist1);
494 fRecEventHist1->SetBinEdges(2,fBinsZv);
495 fRecEventHist1->SetBinEdges(3,fBinsCentrality);
496 fRecEventHist1->GetAxis(0)->SetTitle("Xv (cm)");
497 fRecEventHist1->GetAxis(1)->SetTitle("Yv (cm)");
498 fRecEventHist1->GetAxis(2)->SetTitle("Zv (cm)");
499 fRecEventHist1->GetAxis(3)->SetTitle("Centrality");
500 fRecEventHist1->Sumw2();
503 Int_t binsRecEventHist2[3]={fZvNbins, fMultNbins, fCentralityNbins};
504 Double_t minRecEventHist2[3]={-30., -0.5, 0.};
505 Double_t maxRecEventHist2[3]={30., 10000.5, 100.};
507 fRecEventHist2 = new THnSparseF("fRecEventHist2","Zv:multMB:Centrality",3,binsRecEventHist2,minRecEventHist2,maxRecEventHist2);
508 fRecEventHist2->SetBinEdges(0,fBinsZv);
509 fRecEventHist2->SetBinEdges(1,fBinsMult);
510 fRecEventHist2->SetBinEdges(2,fBinsCentrality);
511 fRecEventHist2->GetAxis(0)->SetTitle("Zv (cm)");
512 fRecEventHist2->GetAxis(1)->SetTitle("multiplicity MB");
513 fRecEventHist2->GetAxis(2)->SetTitle("Centrality");
514 fRecEventHist2->Sumw2();
517 Double_t kFact = 0.1;
518 Int_t binsRecMCEventHist1[4]={100,100,100, fCentralityNbins};
519 Double_t minRecMCEventHist1[4]={-10.0*kFact,-10.0*kFact,-10.0*kFact, 0.};
520 Double_t maxRecMCEventHist1[4]={10.0*kFact,10.0*kFact,10.0*kFact, 100.};
522 fRecMCEventHist1 = new THnSparseF("fRecMCEventHist1","Xv-mcXv:Yv-mcYv:Zv-mcZv:Centrality",4,binsRecMCEventHist1,minRecMCEventHist1,maxRecMCEventHist1);
523 fRecMCEventHist1->SetBinEdges(3,fBinsCentrality);
524 fRecMCEventHist1->GetAxis(0)->SetTitle("Xv-mcXv (cm)");
525 fRecMCEventHist1->GetAxis(1)->SetTitle("Yv-mcYv (cm)");
526 fRecMCEventHist1->GetAxis(2)->SetTitle("Zv-mcZv (cm)");
527 fRecMCEventHist1->GetAxis(3)->SetTitle("Centrality");
529 fRecMCEventHist1->Sumw2();
532 Int_t binsRecMCEventHist2[4]={100,100,fMultNbins, fCentralityNbins};
533 Double_t minRecMCEventHist2[4]={-10.0*kFact,-10.0*kFact,-0.5, 0.};
534 Double_t maxRecMCEventHist2[4]={10.0*kFact,10.0*kFact,10000.5, 100.};
536 fRecMCEventHist2 = new THnSparseF("fRecMCEventHist2","Xv-mcXv:Zv-mcZv:mult:Centrality",4,binsRecMCEventHist2,minRecMCEventHist2,maxRecMCEventHist2);
537 fRecMCEventHist2->SetBinEdges(2,fBinsMult);
538 fRecMCEventHist2->SetBinEdges(3,fBinsCentrality);
539 fRecMCEventHist2->GetAxis(0)->SetTitle("Xv-mcXv (cm)");
540 fRecMCEventHist2->GetAxis(1)->SetTitle("Zv-mcZv (cm)");
541 fRecMCEventHist2->GetAxis(2)->SetTitle("multiplicity");
542 fRecMCEventHist2->GetAxis(3)->SetTitle("Centrality");
543 fRecMCEventHist2->Sumw2();
547 for(Int_t i=0; i<AlidNdPtHelper::kCutSteps; i++)
549 // THnSparse track histograms
551 Int_t binsMCTrackHist1[4]= {fPtNbins, fEtaNbins, 90, fCentralityNbins};
552 Double_t minMCTrackHist1[4]={0.,-1.5,0., 0.};
553 Double_t maxMCTrackHist1[4]={50,1.5,2.*TMath::Pi(), 100.};
554 sprintf(name,"fMCTrackHist1_%d",i);
555 sprintf(title,"mcPt:mcEta:mcPhi:Centrality");
557 fMCTrackHist1[i] = new THnSparseF(name,title,4,binsMCTrackHist1,minMCTrackHist1,maxMCTrackHist1);
558 fMCTrackHist1[i]->SetBinEdges(0,fBinsPt);
559 fMCTrackHist1[i]->SetBinEdges(1,fBinsEta);
560 fMCTrackHist1[i]->SetBinEdges(3,fBinsCentrality);
561 fMCTrackHist1[i]->GetAxis(0)->SetTitle("mcPt (GeV/c)");
562 fMCTrackHist1[i]->GetAxis(1)->SetTitle("mcEta");
563 fMCTrackHist1[i]->GetAxis(2)->SetTitle("mcPhi (rad)");
564 fMCTrackHist1[i]->GetAxis(3)->SetTitle("Centrality");
565 fMCTrackHist1[i]->Sumw2();
567 Int_t binsMCPrimTrackHist1[6]= {fPtNbins,fEtaNbins,6,20,4000, fCentralityNbins};
568 Double_t minMCPrimTrackHist1[6]={0.,-1.5,0.,0.,0., 0.};
569 Double_t maxMCPrimTrackHist1[6]={50.,1.5,6.,20.,4000., 100.};
570 sprintf(name,"fMCPrimTrackHist1_%d",i);
571 sprintf(title,"mcPt:mcEta:pid:mech:mother:Centrality");
573 fMCPrimTrackHist1[i] = new THnSparseF(name,title,6,binsMCPrimTrackHist1,minMCPrimTrackHist1,maxMCPrimTrackHist1);
574 fMCPrimTrackHist1[i]->SetBinEdges(0,fBinsPt);
575 fMCPrimTrackHist1[i]->SetBinEdges(1,fBinsEta);
576 fMCPrimTrackHist1[i]->SetBinEdges(5,fBinsCentrality);
577 fMCPrimTrackHist1[i]->GetAxis(0)->SetTitle("mcPt (GeV/c)");
578 fMCPrimTrackHist1[i]->GetAxis(1)->SetTitle("mcEta");
579 fMCPrimTrackHist1[i]->GetAxis(2)->SetTitle("pid");
580 fMCPrimTrackHist1[i]->GetAxis(3)->SetTitle("mech");
581 fMCPrimTrackHist1[i]->GetAxis(4)->SetTitle("mother");
582 fMCPrimTrackHist1[i]->GetAxis(5)->SetTitle("Centrality");
583 fMCPrimTrackHist1[i]->Sumw2();
585 Int_t binsMCPrimTrackHist2[4]= {4000,20,4000,fCentralityNbins};
586 Double_t minMCPrimTrackHist2[4]={0.,0.,0., 0.};
587 Double_t maxMCPrimTrackHist2[4]={4000.,20.,4000., 100.};
588 sprintf(name,"fMCPrimTrackHist2_%d",i);
589 sprintf(title,"pdg:mech:mother:Centrality");
591 fMCPrimTrackHist2[i] = new THnSparseF(name,title,4,binsMCPrimTrackHist2,minMCPrimTrackHist2,maxMCPrimTrackHist2);
592 fMCPrimTrackHist2[i]->SetBinEdges(3,fBinsCentrality);
593 fMCPrimTrackHist2[i]->GetAxis(0)->SetTitle("pdg");
594 fMCPrimTrackHist2[i]->GetAxis(1)->SetTitle("mech");
595 fMCPrimTrackHist2[i]->GetAxis(2)->SetTitle("mother");
596 fMCPrimTrackHist2[i]->GetAxis(3)->SetTitle("Centrality");
597 fMCPrimTrackHist2[i]->Sumw2();
599 Int_t binsMCSecTrackHist1[6]= {fPtNbins,fEtaNbins,6,20,4000, fCentralityNbins};
600 Double_t minMCSecTrackHist1[6]={0.,-1.5,0.,0.,0., 0.};
601 Double_t maxMCSecTrackHist1[6]={50.,1.5,6.,20.,4000., 100.};
602 sprintf(name,"fMCSecTrackHist1_%d",i);
603 sprintf(title,"mcPt:mcEta:pid:mech:mother:Centrality");
605 fMCSecTrackHist1[i] = new THnSparseF(name,title,6,binsMCSecTrackHist1,minMCSecTrackHist1,maxMCSecTrackHist1);
606 fMCSecTrackHist1[i]->SetBinEdges(0,fBinsPt);
607 fMCSecTrackHist1[i]->SetBinEdges(1,fBinsEta);
608 fMCSecTrackHist1[i]->SetBinEdges(5,fBinsCentrality);
609 fMCSecTrackHist1[i]->GetAxis(0)->SetTitle("mcPt (GeV/c)");
610 fMCSecTrackHist1[i]->GetAxis(1)->SetTitle("mcEta");
611 fMCSecTrackHist1[i]->GetAxis(2)->SetTitle("pid");
612 fMCSecTrackHist1[i]->GetAxis(3)->SetTitle("mech");
613 fMCSecTrackHist1[i]->GetAxis(4)->SetTitle("mother");
614 fMCSecTrackHist1[i]->GetAxis(5)->SetTitle("Centrality");
615 fMCSecTrackHist1[i]->Sumw2();
617 Int_t binsRecTrackHist1[4]={fPtNbins,fEtaNbins,90, fCentralityNbins};
618 Double_t minRecTrackHist1[4]={0.,-1.5,0., 0.};
619 Double_t maxRecTrackHist1[4]={50.,1.5,2.*TMath::Pi(), 100.};
620 sprintf(name,"fRecTrackHist1_%d",i);
621 sprintf(title,"Pt:Eta:Phi:Centrality");
622 fRecTrackHist1[i] = new THnSparseF(name,title,4,binsRecTrackHist1,minRecTrackHist1,maxRecTrackHist1);
623 fRecTrackHist1[i]->SetBinEdges(0,fBinsPt);
624 fRecTrackHist1[i]->SetBinEdges(1,fBinsEta);
625 fRecTrackHist1[i]->SetBinEdges(3,fBinsCentrality);
626 fRecTrackHist1[i]->GetAxis(0)->SetTitle("p_{T} (GeV/c)");
627 fRecTrackHist1[i]->GetAxis(1)->SetTitle("#eta");
628 fRecTrackHist1[i]->GetAxis(2)->SetTitle("#phi (rad)");
629 fRecTrackHist1[i]->GetAxis(3)->SetTitle("Centrality");
630 fRecTrackHist1[i]->Sumw2();
632 sprintf(name,"fRecTrackHist2_%d",i);
633 sprintf(title,"Zv:Pt:Eta:Centrality");
634 fRecTrackHist2[i] = new THnSparseF(name,title,4,binsTrackEvent);
635 fRecTrackHist2[i]->SetBinEdges(0,fBinsZv);
636 fRecTrackHist2[i]->SetBinEdges(1,fBinsPt);
637 fRecTrackHist2[i]->SetBinEdges(2,fBinsEta);
638 fRecTrackHist2[i]->SetBinEdges(3,fBinsCentrality);
639 fRecTrackHist2[i]->GetAxis(0)->SetTitle("Zv (cm)");
640 fRecTrackHist2[i]->GetAxis(1)->SetTitle("p_{T} (GeV/c)");
641 fRecTrackHist2[i]->GetAxis(2)->SetTitle("#eta");
642 fRecTrackHist2[i]->GetAxis(3)->SetTitle("Centrality");
643 fRecTrackHist2[i]->Sumw2();
646 Int_t binsRecTrackMultHist1[3]={fPtNbins,fMultNbins, fCentralityNbins};
647 Double_t minRecTrackMultHist1[3]={0.,-0.5, -0.};
648 Double_t maxRecTrackMultHist1[3]={50.,10000.5, 100.};
649 sprintf(name,"fRecTrackMultHist_%d",i);
650 sprintf(title,"Pt:Mult:Centrality");
651 fRecTrackMultHist1[i] = new THnSparseF(name,title,3,binsRecTrackMultHist1,minRecTrackMultHist1,maxRecTrackMultHist1);
652 fRecTrackMultHist1[i]->SetBinEdges(0,fBinsPt);
653 fRecTrackMultHist1[i]->SetBinEdges(1,fBinsMult);
654 fRecTrackMultHist1[i]->SetBinEdges(2,fBinsCentrality);
655 fRecTrackMultHist1[i]->GetAxis(0)->SetTitle("Pt (GeV/c)");
656 fRecTrackMultHist1[i]->GetAxis(1)->SetTitle("multiplicity");
657 fRecTrackMultHist1[i]->GetAxis(2)->SetTitle("Centrality");
658 fRecTrackMultHist1[i]->Sumw2();
661 Int_t binsRecMCTrackHist1[5] = {fPtNbins,fEtaNbins,100,100, fCentralityNbins};
662 Double_t minRecMCTrackHist1[5]={0.,-1.5,-0.5,-0.5, 0.};
663 Double_t maxRecMCTrackHist1[5]={50.,1.5,0.5,0.5, 100.};
665 sprintf(name,"fRecMCTrackHist1");
666 sprintf(title,"mcPt:mcEta:(Pt-mcPt)/mcPt:(Eta-mcEta):Centrality");
667 fRecMCTrackHist1 = new THnSparseF(name,title,5,binsRecMCTrackHist1,minRecMCTrackHist1,maxRecMCTrackHist1);
668 fRecMCTrackHist1->SetBinEdges(0,fBinsPt);
669 fRecMCTrackHist1->SetBinEdges(1,fBinsEta);
670 fRecMCTrackHist1->SetBinEdges(4,fBinsCentrality);
671 fRecMCTrackHist1->GetAxis(0)->SetTitle("mcPt (GeV/c)");
672 fRecMCTrackHist1->GetAxis(1)->SetTitle("mcEta");
673 fRecMCTrackHist1->GetAxis(2)->SetTitle("(Pt-mcPt)/mcPt");
674 fRecMCTrackHist1->GetAxis(3)->SetTitle("Eta-mcEta");
675 fRecMCTrackHist1->GetAxis(4)->SetTitle("Centrality");
677 Int_t binsMCMultRecTrackHist1[4] = {fPtNbins,fEtaNbins,6, fCentralityNbins};
678 Double_t minMCMultRecTrackHist1[4]={0.,-1.5,0., 0.};
679 Double_t maxMCMultRecTrackHist1[4]={50.,1.5,6., 100.};
680 sprintf(name,"fMCMultRecTrackHist1");
681 sprintf(title,"mcPt:mcEta:pid:Centrality");
682 fMCMultRecTrackHist1 = new THnSparseF(name,title,4,binsMCMultRecTrackHist1,minMCMultRecTrackHist1,maxMCMultRecTrackHist1);
683 fMCMultRecTrackHist1->SetBinEdges(0,fBinsPt);
684 fMCMultRecTrackHist1->SetBinEdges(1,fBinsEta);
685 fMCMultRecTrackHist1->SetBinEdges(3,fBinsCentrality);
686 fMCMultRecTrackHist1->GetAxis(0)->SetTitle("mcPt (GeV/c)");
687 fMCMultRecTrackHist1->GetAxis(1)->SetTitle("mcEta");
688 fMCMultRecTrackHist1->GetAxis(2)->SetTitle("pid");
689 fMCMultRecTrackHist1->GetAxis(3)->SetTitle("Centrality");
691 //nClust:chi2PerClust:pt:eta:phi:Centrality
692 Int_t binsRecTrackHist3[6]={160,100,fPtNbins,fEtaNbins,90, fCentralityNbins};
693 Double_t minRecTrackHist3[6]={0., 0., 0., -1.5, 0., 0.};
694 Double_t maxRecRecTrackHist3[6]={160.,10., 50., 1.5, 2.*TMath::Pi(), 100.};
696 fRecTrackHist3 = new THnSparseF("fRecTrackHist3","nClust:chi2PerClust:pt:eta:phi:Centrality",6,binsRecTrackHist3,minRecTrackHist3,maxRecRecTrackHist3);
697 fRecTrackHist3->SetBinEdges(2,fBinsPt);
698 fRecTrackHist3->SetBinEdges(3,fBinsEta);
699 fRecTrackHist3->SetBinEdges(5,fBinsCentrality);
700 fRecTrackHist3->GetAxis(0)->SetTitle("nClust");
701 fRecTrackHist3->GetAxis(1)->SetTitle("chi2PerClust");
702 fRecTrackHist3->GetAxis(2)->SetTitle("p_{T} (GeV/c)");
703 fRecTrackHist3->GetAxis(3)->SetTitle("#eta");
704 fRecTrackHist3->GetAxis(4)->SetTitle("#phi (rad)");
705 fRecTrackHist3->GetAxis(5)->SetTitle("Centrality");
706 fRecTrackHist3->Sumw2();
710 fAnalysisFolder = CreateFolder("folderdNdPt","Analysis dNdPt Folder");
712 // init trigger analysis (for zdc cut)
713 fTriggerAnalysis = new AliTriggerAnalysis;
719 //_____________________________________________________________________________
720 void AlidNdPtAnalysisPbPb::Process(AliESDEvent *const esdEvent, AliMCEvent *const mcEvent)
722 // init if not done already
723 if (!fIsInit) { Init(); }
725 // Process real and/or simulated events
728 AliDebug(AliLog::kError, "esdEvent not available");
732 // get selection cuts
733 AlidNdPtEventCuts *evtCuts = GetEventCuts();
734 AlidNdPtAcceptanceCuts *accCuts = GetAcceptanceCuts();
735 AlidNdPtAcceptanceCuts *recCuts = GetRecAcceptanceCuts();
736 AliESDtrackCuts *esdTrackCuts = GetTrackCuts();
738 if(!evtCuts || !accCuts || !esdTrackCuts) {
739 AliDebug(AliLog::kError, "cuts not available");
742 if (0 == recCuts) { recCuts = accCuts;}
745 Bool_t isEventTriggered = kTRUE;
746 AliPhysicsSelection *physicsSelection = NULL;
747 AliTriggerAnalysis* triggerAnalysis = NULL;
750 AliInputEventHandler* inputHandler = (AliInputEventHandler*) AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler();
753 Printf("ERROR: Could not receive input handler");
757 if(evtCuts->IsTriggerRequired())
760 isEventTriggered = inputHandler->IsEventSelected() & AliVEvent::kMB;
762 physicsSelection = static_cast<AliPhysicsSelection*> (inputHandler->GetEventSelection());
763 if(!physicsSelection) return;
764 //SetPhysicsTriggerSelection(physicsSelection);
766 if (isEventTriggered && (GetTrigger() == AliTriggerAnalysis::kV0AND)) {
767 // set trigger (V0AND)
768 triggerAnalysis = physicsSelection->GetTriggerAnalysis();
769 if(!triggerAnalysis) return;
770 isEventTriggered = triggerAnalysis->IsOfflineTriggerFired(esdEvent, GetTrigger());
775 // centrality determination
776 Float_t centralityF = -1.;
777 AliCentrality *esdCentrality = esdEvent->GetCentrality();
778 centralityF = esdCentrality->GetCentralityPercentile(fCentralityEstimator.Data());
780 // use MC information
781 AliHeader* header = 0;
782 AliGenEventHeader* genHeader = 0;
786 Int_t multMCTrueTracks = 0;
791 AliDebug(AliLog::kError, "mcEvent not available");
794 // get MC event header
795 header = mcEvent->Header();
797 AliDebug(AliLog::kError, "Header not available");
801 stack = mcEvent->Stack();
803 AliDebug(AliLog::kError, "Stack not available");
808 genHeader = header->GenEventHeader();
810 AliDebug(AliLog::kError, "Could not retrieve genHeader from Header");
813 genHeader->PrimaryVertex(vtxMC);
815 Double_t vMCEventHist1[4]={vtxMC[0],vtxMC[1],vtxMC[2],centralityF};
816 fMCEventHist1->Fill(vMCEventHist1);
818 // multipliticy of all MC primary tracks
819 // in Zv, pt and eta ranges)
820 multMCTrueTracks = AlidNdPtHelper::GetMCTrueTrackMult(mcEvent,evtCuts,accCuts);
824 // get reconstructed vertex
825 const AliESDVertex* vtxESD = 0;
826 if(evtCuts->IsRecVertexRequired())
828 if(GetAnalysisMode() == AlidNdPtHelper::kTPC) {
829 vtxESD = esdEvent->GetPrimaryVertexTPC();
831 else if(GetAnalysisMode() == AlidNdPtHelper::kTPCITS) {
832 vtxESD = esdEvent->GetPrimaryVertexTracks();
839 Bool_t isEventOK = evtCuts->AcceptEvent(esdEvent,mcEvent,vtxESD);
840 //printf("isEventOK %d, isEventTriggered %d \n",isEventOK, isEventTriggered);
841 //printf("GetAnalysisMode() %d \n",GetAnalysisMode());
843 // vertex contributors
844 Int_t multMBTracks = 0;
845 if(GetAnalysisMode() == AlidNdPtHelper::kTPC || GetAnalysisMode() == AlidNdPtHelper::kTPCITS)
847 if(vtxESD->GetStatus()) {
848 multMBTracks = vtxESD->GetNContributors();
852 AliDebug(AliLog::kError, Form("Found analysis type %d", GetAnalysisMode()));
856 TObjArray *allChargedTracks=0;
857 //Int_t multAll=0, multAcc=0, multRec=0;
858 Int_t multAll=0, multRec=0;
859 Int_t *labelsAll=0, *labelsAcc=0, *labelsRec=0;
862 if(isEventOK && isEventTriggered)
865 // get all charged tracks
866 allChargedTracks = AlidNdPtHelper::GetAllChargedTracks(esdEvent,GetAnalysisMode());
867 if(!allChargedTracks) return;
869 Int_t entries = allChargedTracks->GetEntries();
871 labelsAll = new Int_t[entries];
872 labelsAcc = new Int_t[entries];
873 labelsRec = new Int_t[entries];
874 for(Int_t i=0; i<entries;++i)
876 AliESDtrack *track = (AliESDtrack*)allChargedTracks->At(i);
879 if(track->Charge()==0) continue;
882 // only postive charged
883 if(GetParticleMode() == AlidNdPtHelper::kPlus && track->Charge() < 0)
886 // only negative charged
887 if(GetParticleMode() == AlidNdPtHelper::kMinus && track->Charge() > 0)
891 Double_t values[4] = {vtxESD->GetZv(),track->Pt(),track->Eta(), centralityF};
893 fRecTrackHist2[AlidNdPtHelper::kAllTracks]->Fill(values);
894 FillHistograms(track,stack,AlidNdPtHelper::kAllTracks,centralityF);
895 labelsAll[multAll] = TMath::Abs(track->GetLabel());
898 if(esdTrackCuts->AcceptTrack(track) && accCuts->AcceptTrack(track) && recCuts->AcceptTrackLocalTPC(track)) {
900 fRecTrackHist2[AlidNdPtHelper::kRecTracks]->Fill(values);
901 FillHistograms(track,stack,AlidNdPtHelper::kRecTracks,centralityF);
902 labelsRec[multRec] = TMath::Abs(track->GetLabel());
910 // fill track multiplicity histograms
911 //FillHistograms(allChargedTracks,labelsAll,multAll,labelsAcc,multAcc,labelsRec,multRec);
913 Double_t vRecEventHist1[4] = {vtxESD->GetXv(),vtxESD->GetYv(),vtxESD->GetZv(),centralityF};
914 fRecEventHist1->Fill(vRecEventHist1);
916 Double_t vRecEventHist2[3] = {vtxESD->GetZv(),multMBTracks,centralityF};
917 fRecEventHist2->Fill(vRecEventHist2);
919 } // triggered and event vertex
924 // event level corrections (zv,N_MB)
927 Double_t vEventMatrix[3] = {vtxMC[2],multMBTracks,centralityF};
928 fGenEventMatrix->Fill(vEventMatrix);
929 if(isEventTriggered) fTriggerEventMatrix->Fill(vEventMatrix);
930 if(isEventOK && isEventTriggered) fRecEventMatrix->Fill(vEventMatrix);
933 // track-event level corrections (zv,pt,eta)
935 for (Int_t iMc = 0; iMc < stack->GetNtrack(); ++iMc)
937 TParticle* particle = stack->Particle(iMc);
941 // only charged particles
942 if(!particle->GetPDG()) continue;
943 Double_t charge = particle->GetPDG()->Charge()/3.;
944 if ( TMath::Abs(charge) < 0.001 )
947 // only postive charged
948 if(GetParticleMode() == AlidNdPtHelper::kPlus && charge < 0.)
951 // only negative charged
952 if(GetParticleMode() == AlidNdPtHelper::kMinus && charge > 0.)
956 Bool_t prim = stack->IsPhysicalPrimary(iMc);
960 if(accCuts->AcceptTrack(particle))
962 Double_t vTrackEventMatrix[4] = {vtxMC[2], particle->Pt(), particle->Eta(), centralityF};
963 fGenTrackEventMatrix->Fill(vTrackEventMatrix);
965 if(!isEventTriggered) continue;
966 fTriggerTrackEventMatrix->Fill(vTrackEventMatrix);
968 if(!isEventOK) continue;
969 fRecTrackEventMatrix->Fill(vTrackEventMatrix);
971 }// if(accCuts->AcceptTrack(particle))
972 }// for (Int_t iMc = 0; iMc < stack->GetNtrack(); ++iMc)
975 // track-level corrections (zv,pt,eta)
977 if(isEventOK && isEventTriggered)
980 // fill MC and rec event control histograms
982 Double_t vRecMCEventHist1[4] = {vtxESD->GetXv()-vtxMC[0],vtxESD->GetYv()-vtxMC[1],vtxESD->GetZv()-vtxMC[2], centralityF};
983 fRecMCEventHist1->Fill(vRecMCEventHist1);//
985 Double_t vRecMCEventHist2[4] = {vtxESD->GetXv()-vtxMC[0],vtxESD->GetZv()-vtxMC[2],multMBTracks, centralityF};
986 fRecMCEventHist2->Fill(vRecMCEventHist2);
991 // MC histograms for track efficiency studies
993 for (Int_t iMc = 0; iMc < stack->GetNtrack(); ++iMc)
995 TParticle* particle = stack->Particle(iMc);
999 Double_t vTrackMatrix[4] = {vtxMC[2],particle->Pt(),particle->Eta(),centralityF};
1001 // only charged particles
1002 if(!particle->GetPDG()) continue;
1003 Double_t charge = particle->GetPDG()->Charge()/3.;
1004 if (TMath::Abs(charge) < 0.001)
1007 // only postive charged
1008 if(GetParticleMode() == AlidNdPtHelper::kPlus && charge < 0.)
1011 // only negative charged
1012 if(GetParticleMode() == AlidNdPtHelper::kMinus && charge > 0.)
1016 Bool_t prim = stack->IsPhysicalPrimary(iMc);
1019 if(accCuts->AcceptTrack(particle))
1022 if( AlidNdPtHelper::IsPrimaryParticle(stack, iMc, GetParticleMode()) )
1023 fGenPrimTrackMatrix->Fill(vTrackMatrix);
1025 // fill control histograms
1027 FillHistograms(stack,iMc,AlidNdPtHelper::kAccTracks, centralityF);
1029 // check multiple found tracks
1030 Int_t multCount = 0;
1031 for(Int_t iRec=0; iRec<multRec; ++iRec)
1033 if(iMc == labelsRec[iRec])
1038 fRecMultTrackMatrix->Fill(vTrackMatrix);
1040 // fill control histogram
1042 Int_t pid = AlidNdPtHelper::ConvertPdgToPid(particle);
1043 Double_t vMCMultRecTrackHist1[4] = {particle->Pt(), particle->Eta(), pid, centralityF};
1044 fMCMultRecTrackHist1->Fill(vMCMultRecTrackHist1);
1050 // check reconstructed
1051 for(Int_t iRec=0; iRec<multRec; ++iRec)
1053 if(iMc == labelsRec[iRec])
1055 fRecTrackMatrix->Fill(vTrackMatrix);
1057 if( AlidNdPtHelper::IsPrimaryParticle(stack, iMc, GetParticleMode()) ) {
1058 fRecPrimTrackMatrix->Fill(vTrackMatrix);
1060 if(!prim) fRecSecTrackMatrix->Fill(vTrackMatrix);
1062 // fill control histograms
1064 FillHistograms(stack,iMc,AlidNdPtHelper::kRecTracks, centralityF);
1067 }//if(iMc == labelsRec[iRec])
1074 if(allChargedTracks) delete allChargedTracks; allChargedTracks = 0;
1075 if(labelsAll) delete [] labelsAll; labelsAll = 0;
1076 if(labelsAcc) delete [] labelsAcc; labelsAcc = 0;
1077 if(labelsRec) delete [] labelsRec; labelsRec = 0;
1081 //_____________________________________________________________________________
1082 void AlidNdPtAnalysisPbPb::FillHistograms(TObjArray *const allChargedTracks,Int_t *const labelsAll,Int_t multAll,Int_t *const labelsAcc,Int_t multAcc,Int_t *const labelsRec,Int_t multRec, Float_t centralityF) {
1083 // multiplicity histograms
1086 if(!allChargedTracks) return;
1087 if(!labelsAll) return;
1088 if(!labelsAcc) return;
1089 if(!labelsRec) return;
1091 Int_t entries = allChargedTracks->GetEntries();
1092 for(Int_t i=0; i<entries; ++i)
1094 AliESDtrack *track = (AliESDtrack*)allChargedTracks->At(i);
1095 if(!track) continue;
1096 if(track->Charge() == 0) continue;
1098 Int_t label = TMath::Abs(track->GetLabel());
1099 for(Int_t iAll=0; iAll<multAll; ++iAll) {
1100 if(label == labelsAll[iAll]) {
1101 Double_t v1[3] = {track->Pt(), multAll, centralityF};
1102 fRecTrackMultHist1[AlidNdPtHelper::kAllTracks]->Fill(v1);
1105 for(Int_t iAcc=0; iAcc<multAcc; ++iAcc) {
1106 if(label == labelsAcc[iAcc]) {
1107 Double_t v2[3] = {track->Pt(), multAcc, centralityF};
1108 fRecTrackMultHist1[AlidNdPtHelper::kAccTracks]->Fill(v2);
1111 for(Int_t iRec=0; iRec<multRec; ++iRec) {
1112 if(label == labelsRec[iRec]) {
1113 Double_t v3[3] = {track->Pt(), multRec, centralityF};
1114 fRecTrackMultHist1[AlidNdPtHelper::kRecTracks]->Fill(v3);
1120 //_____________________________________________________________________________
1121 void AlidNdPtAnalysisPbPb::FillHistograms(AliESDtrack *const esdTrack, AliStack *const stack, AlidNdPtHelper::TrackObject trackObj, Float_t centralityF)
1125 // Fill ESD track and MC histograms
1127 if(!esdTrack) return;
1129 Float_t q = esdTrack->Charge();
1130 if(TMath::Abs(q) < 0.001) return;
1132 Float_t pt = esdTrack->Pt();
1133 Float_t eta = esdTrack->Eta();
1134 Float_t phi = esdTrack->Phi();
1136 Float_t dca[2], bCov[3];
1137 esdTrack->GetImpactParameters(dca,bCov);
1139 Int_t nClust = esdTrack->GetTPCclusters(0);
1140 Float_t chi2PerCluster = 0.;
1141 if(nClust>0.) chi2PerCluster = esdTrack->GetTPCchi2()/Float_t(nClust);
1145 Double_t values[4] = {pt,eta,phi,centralityF};
1146 fRecTrackHist1[trackObj]->Fill(values);
1149 // Fill rec vs MC information
1153 Int_t label = TMath::Abs(esdTrack->GetLabel());
1154 //if(label == 0) return;
1156 if(label > stack->GetNtrack()) return;
1157 TParticle* particle = stack->Particle(label);
1158 if(!particle) return;
1160 Int_t motherPdg = -1;
1161 TParticle* mother = 0;
1163 //TParticle* prim_mother = AlidNdPtHelper::FindPrimaryMother(stack,label);
1164 Int_t motherLabel = particle->GetMother(0);
1165 if(motherLabel>0) mother = stack->Particle(motherLabel);
1166 if(mother) motherPdg = TMath::Abs(mother->GetPdgCode()); // take abs for visualisation only
1167 //Int_t mech = particle->GetUniqueID(); // production mechanism
1169 if(!particle->GetPDG()) return;
1170 Double_t gq = particle->GetPDG()->Charge()/3.0; // Charge units |e|/3
1171 if(TMath::Abs(gq)<0.001) return;
1172 Float_t gpt = particle->Pt();
1173 Float_t geta = particle->Eta();
1176 //printf("pt %f, gpt %f \n",pt,gpt);
1177 if(gpt) dpt = (pt-gpt)/gpt;
1178 Double_t deta = (eta-geta);
1181 if(trackObj == AlidNdPtHelper::kRecTracks) //RecTracks???
1183 Double_t vTrackPtCorrelationMatrix[4]={pt,gpt,geta,centralityF};
1184 fTrackPtCorrelationMatrix->Fill(vTrackPtCorrelationMatrix);
1186 Double_t vRecMCTrackHist1[5]={gpt,geta,dpt,deta,centralityF};
1187 fRecMCTrackHist1->Fill(vRecMCTrackHist1);
1191 //_____________________________________________________________________________
1192 void AlidNdPtAnalysisPbPb::FillHistograms(AliStack *const stack, Int_t label, AlidNdPtHelper::TrackObject trackObj, Float_t centralityF)
1195 // Fill MC histograms
1198 if(label > stack->GetNtrack()) return;
1199 TParticle* particle = stack->Particle(label);
1200 if(!particle) return;
1202 Int_t motherPdg = -1;
1203 TParticle* mother = 0;
1205 //TParticle* prim_mother = AlidNdPtHelper::FindPrimaryMother(stack,label);
1206 Int_t motherLabel = particle->GetMother(0);
1207 if(motherLabel>0) mother = stack->Particle(motherLabel);
1208 if(mother) motherPdg = TMath::Abs(mother->GetPdgCode()); // take abs for visualisation only
1209 Int_t mech = particle->GetUniqueID(); // production mechanism
1211 if(!particle->GetPDG()) return;
1212 Double_t gq = particle->GetPDG()->Charge()/3.0; // Charge units |e|/3
1213 if(TMath::Abs(gq) < 0.001) return;
1215 Float_t gpt = particle->Pt();
1216 //Float_t qgpt = particle->Pt() * gq;
1217 Float_t geta = particle->Eta();
1218 Float_t gphi = particle->Phi();
1219 //Float_t gpz = particle->Pz();
1221 Bool_t prim = stack->IsPhysicalPrimary(label);
1222 //Float_t vx = particle->Vx(); Float_t vy = particle->Vy(); Float_t vz = particle->Vz();
1224 Int_t pid = AlidNdPtHelper::ConvertPdgToPid(particle);
1226 //if(prim&&pid==5) printf("pdgcode %d, production mech %d \n",particle->GetPdgCode(),mech);
1227 //if(!prim) printf("motherPdg %d, particle %d, production mech %d\n",motherPdg, particle->GetPdgCode(),mech);
1232 Double_t vMCTrackHist1[4] = {gpt,geta,gphi,centralityF};
1233 fMCTrackHist1[trackObj]->Fill(vMCTrackHist1);
1235 Double_t vMCPrimTrackHist1[6] = {gpt,geta,pid,mech,motherPdg,centralityF};
1236 Double_t vMCPrimTrackHist2[4] = {TMath::Abs(particle->GetPdgCode()),mech,motherPdg,centralityF};
1238 //if(prim && AliPWG0Helper::IsPrimaryCharged(particle, label)) fMCPrimTrackHist1[trackObj]->Fill(vMCPrimTrackHist1);
1241 fMCPrimTrackHist1[trackObj]->Fill(vMCPrimTrackHist1);
1242 if(pid == 5) fMCPrimTrackHist2[trackObj]->Fill(vMCPrimTrackHist2);
1245 fMCSecTrackHist1[trackObj]->Fill(vMCPrimTrackHist1);
1250 //_____________________________________________________________________________
1251 Long64_t AlidNdPtAnalysisPbPb::Merge(TCollection* const list)
1253 // Merge list of objects (needed by PROOF)
1255 // init if not done already
1256 if (!fIsInit) { Init(); }
1261 if (list->IsEmpty())
1264 TIterator* iter = list->MakeIterator();
1268 //TList *collPhysSelection = new TList;
1270 // collection of generated histograms
1273 while((obj = iter->Next()) != 0) {
1274 AlidNdPtAnalysisPbPb* entry = dynamic_cast<AlidNdPtAnalysisPbPb*>(obj);
1275 if (entry == 0) continue;
1277 // physics selection
1278 //printf("entry->GetPhysicsTriggerSelection() %p \n", entry->GetPhysicsTriggerSelection());
1279 //AliPhysicsSelection *physSel = entry->GetPhysicsTriggerSelection();
1280 //if( physSel ) collPhysSelection->Add(physSel);
1283 fTrackPtCorrelationMatrix->Add(entry->fTrackPtCorrelationMatrix);
1286 fGenEventMatrix->Add(entry->fGenEventMatrix);
1288 fTriggerEventMatrix->Add(entry->fTriggerEventMatrix);
1290 fRecEventMatrix->Add(entry->fRecEventMatrix);
1292 fGenTrackEventMatrix->Add(entry->fGenTrackEventMatrix);
1294 fTriggerTrackEventMatrix->Add(entry->fTriggerTrackEventMatrix);
1296 fRecTrackEventMatrix->Add(entry->fRecTrackEventMatrix);
1299 fGenTrackMatrix->Add(entry->fGenTrackMatrix);
1300 fGenPrimTrackMatrix->Add(entry->fGenPrimTrackMatrix);
1301 fRecPrimTrackMatrix->Add(entry->fRecPrimTrackMatrix);
1303 fRecTrackMatrix->Add(entry->fRecTrackMatrix);
1304 fRecSecTrackMatrix->Add(entry->fRecSecTrackMatrix);
1306 fRecMultTrackMatrix->Add(entry->fRecMultTrackMatrix);
1309 // control analysis histograms
1311 fMCEventHist1->Add(entry->fMCEventHist1);
1312 fRecEventHist1->Add(entry->fRecEventHist1);
1313 fRecEventHist2->Add(entry->fRecEventHist2);
1314 fRecMCEventHist1->Add(entry->fRecMCEventHist1);
1315 fRecMCEventHist2->Add(entry->fRecMCEventHist2);
1318 for(Int_t i=0; i<AlidNdPtHelper::kCutSteps; i++) {
1319 fMCTrackHist1[i]->Add(entry->fMCTrackHist1[i]);
1321 fMCPrimTrackHist1[i]->Add(entry->fMCPrimTrackHist1[i]);
1322 fMCPrimTrackHist2[i]->Add(entry->fMCPrimTrackHist2[i]);
1323 fMCSecTrackHist1[i]->Add(entry->fMCSecTrackHist1[i]);
1325 fRecTrackHist1[i]->Add(entry->fRecTrackHist1[i]);
1326 fRecTrackHist2[i]->Add(entry->fRecTrackHist2[i]);
1327 fRecTrackMultHist1[i]->Add(entry->fRecTrackMultHist1[i]);
1329 fRecMCTrackHist1->Add(entry->fRecMCTrackHist1);
1330 fMCMultRecTrackHist1->Add(entry->fMCMultRecTrackHist1);
1331 fRecTrackHist3->Add(entry->fRecTrackHist3);
1336 //AliPhysicsSelection *trigSelection = GetPhysicsTriggerSelection();
1337 //if( trigSelection ) trigSelection->Merge(collPhysSelection);
1338 //if(collPhysSelection) delete collPhysSelection;
1345 //_____________________________________________________________________________
1346 void AlidNdPtAnalysisPbPb::Analyse()
1349 // init if not done already
1350 if (!fIsInit) { Init(); }
1352 // Analyse histograms
1354 TH1::AddDirectory(kFALSE);
1355 TH1 *h=0, *h1=0, *h2=0, *h2c = 0;
1360 TObjArray *aFolderObj = new TObjArray;
1363 // LHC backgraund in all and 0-bins
1365 //AliPhysicsSelection *trigSel = GetPhysicsTriggerSelection();
1366 //trigSel->SaveHistograms("physics_selection");
1369 // Reconstructed event vertex
1371 h = fRecEventHist1->Projection(0);
1375 h = fRecEventHist1->Projection(1);
1379 h = fRecEventHist1->Projection(2);
1386 h = fRecEventHist2->Projection(1);
1387 h->SetName("multMB");
1390 h2D = fRecEventHist2->Projection(0,1);
1391 h2D->SetName("Zv_vs_multiplicity_MB");
1392 aFolderObj->Add(h2D);
1395 // reconstructed pt histograms
1397 h = fRecTrackHist1[0]->Projection(0);
1398 h->Scale(1.,"width");
1399 h->SetName("pt_all_ch");
1402 h = fRecTrackHist1[1]->Projection(0);
1403 h->Scale(1.,"width");
1404 h->SetName("pt_acc");
1407 h = fRecTrackHist1[2]->Projection(0);
1408 h->Scale(1.,"width");
1409 h->SetName("pt_rec");
1413 // reconstructed eta histograms
1415 h = fRecTrackHist1[0]->Projection(1);
1416 h->SetName("eta_all_ch");
1419 h = fRecTrackHist1[1]->Projection(1);
1420 h->SetName("eta_acc");
1423 h = fRecTrackHist1[2]->Projection(1);
1424 h->SetName("eta_rec");
1428 // reconstructed phi histograms
1430 h = fRecTrackHist1[0]->Projection(2);
1431 h->SetName("phi_all_ch");
1434 h = fRecTrackHist1[1]->Projection(2);
1435 h->SetName("phi_acc");
1438 h = fRecTrackHist1[2]->Projection(2);
1439 h->SetName("phi_rec");
1443 // reconstructed eta:pt histograms
1445 h2D = fRecTrackHist1[0]->Projection(1,0);
1446 h2D->SetName("pt_eta_all_ch");
1447 aFolderObj->Add(h2D);
1449 h2D = fRecTrackHist1[1]->Projection(1,0);
1450 h2D->SetName("pt_eta_acc");
1451 aFolderObj->Add(h2D);
1453 h2D = fRecTrackHist1[2]->Projection(1,0);
1454 h2D->SetName("pt_eta_rec");
1455 aFolderObj->Add(h2D);
1458 // reconstructed phi:pt histograms
1460 h2D = fRecTrackHist1[0]->Projection(2,0);
1461 h2D->SetName("pt_phi_all_ch");
1462 aFolderObj->Add(h2D);
1464 h2D = fRecTrackHist1[1]->Projection(2,0);
1465 h2D->SetName("pt_phi_acc");
1466 aFolderObj->Add(h2D);
1468 h2D = fRecTrackHist1[2]->Projection(2,0);
1469 h2D->SetName("pt_phi_rec");
1470 aFolderObj->Add(h2D);
1473 // reconstructed phi:eta histograms
1475 h2D = fRecTrackHist1[0]->Projection(2,1);
1476 h2D->SetName("eta_phi_all_ch");
1477 aFolderObj->Add(h2D);
1479 h2D = fRecTrackHist1[1]->Projection(2,1);
1480 h2D->SetName("eta_phi_acc");
1481 aFolderObj->Add(h2D);
1483 h2D = fRecTrackHist1[2]->Projection(2,1);
1484 h2D->SetName("eta_phi_rec");
1485 aFolderObj->Add(h2D);
1488 // reconstructed nClust, chi2 vs pt, eta, phi
1492 h2D = fRecTrackHist3->Projection(0,1);
1493 h2D->SetName("nClust_chi2_rec");
1494 aFolderObj->Add(h2D);
1496 h2D = fRecTrackHist3->Projection(0,2);
1497 h2D->SetName("nClust_pt_rec");
1498 aFolderObj->Add(h2D);
1500 h2D = fRecTrackHist3->Projection(0,3);
1501 h2D->SetName("nClust_eta_rec");
1502 aFolderObj->Add(h2D);
1504 h2D = fRecTrackHist3->Projection(0,4);
1505 h2D->SetName("nClust_phi_rec");
1506 aFolderObj->Add(h2D);
1508 h2D = fRecTrackHist3->Projection(1,2);
1509 h2D->SetName("chi2_pt_rec");
1510 aFolderObj->Add(h2D);
1512 h2D = fRecTrackHist3->Projection(1,3);
1513 h2D->SetName("chi2_eta_rec");
1514 aFolderObj->Add(h2D);
1516 h2D = fRecTrackHist3->Projection(1,4);
1517 h2D->SetName("chi2_phi_rec");
1518 aFolderObj->Add(h2D);
1523 // calculate corrections for empty events
1528 // normalised zv to generate zv for triggered events
1530 h = fRecEventHist2->Projection(0);
1531 if( h->Integral() ) h->Scale(1./h->Integral());
1532 h->SetName("zv_distribution_norm");
1541 // Event vertex resolution
1543 h2D = fRecMCEventHist2->Projection(0,2);
1544 h2D->SetName("DeltaXv_vs_mult");
1545 aFolderObj->Add(h2D);
1547 h2D = fRecMCEventHist2->Projection(1,2);
1548 h2D->SetName("DeltaZv_vs_mult");
1549 aFolderObj->Add(h2D);
1552 // normalised zv to get trigger/trigger+vertex event differences
1553 // F(zv) = E_trig(zv,0)/Int(E_trig(zv,0) / Sum(E_trigvtx(zv,n))/Sum(Int(E_trigvtx(zv,n))dzv)
1555 fTriggerEventMatrix->GetAxis(1)->SetRangeUser(0.,0.);
1556 h = fTriggerEventMatrix->Projection(0);
1557 h2D = fTriggerEventMatrix->Projection(0,1);
1558 if(h2D->Integral()) h->Scale(1./h2D->Integral());
1560 h1 = fRecEventMatrix->Projection(0);
1561 h2D = fRecEventMatrix->Projection(0,1);
1562 if(h2D->Integral()) h1->Scale(1./h2D->Integral());
1565 h->SetName("zv_empty_events_norm");
1568 fTriggerEventMatrix->GetAxis(1)->SetRange(1,fTriggerEventMatrix->GetAxis(1)->GetNbins());
1571 // rec. vs true track pt correlation matrix
1573 hs = (THnSparse*)fTrackPtCorrelationMatrix->Clone("track_pt_correlation_matrix");
1574 aFolderObj->Add(hs);
1577 // trigger efficiency for INEL
1579 h = AlidNdPtHelper::GenerateCorrMatrix(fTriggerEventMatrix->Projection(0),fGenEventMatrix->Projection(0),"zv_trig_INEL_eff_matrix");
1584 // trigger bias correction (MB to INEL)
1586 hs = AlidNdPtHelper::GenerateCorrMatrix(fGenEventMatrix,fTriggerEventMatrix,"zv_mult_trig_MBtoInel_corr_matrix");
1587 aFolderObj->Add(hs);
1589 h = AlidNdPtHelper::GenerateCorrMatrix(fGenEventMatrix->Projection(0),fTriggerEventMatrix->Projection(0),"zv_trig_MBtoInel_corr_matrix");
1592 h2D = AlidNdPtHelper::GenerateCorrMatrix(fGenEventMatrix->Projection(0,1),fTriggerEventMatrix->Projection(0,1),"zv_mult_trig_MBtoInel_corr_matrix_2D");
1593 aFolderObj->Add(h2D);
1596 h = AlidNdPtHelper::GenerateCorrMatrix(fGenEventMatrix->Projection(1),fTriggerEventMatrix->Projection(1),"mult_trig_MBtoInel_corr_matrix");
1601 // event vertex reconstruction correction (MB)
1603 hs = AlidNdPtHelper::GenerateCorrMatrix(fTriggerEventMatrix,fRecEventMatrix,"zv_mult_event_corr_matrix");
1604 aFolderObj->Add(hs);
1606 h2D = AlidNdPtHelper::GenerateCorrMatrix(fTriggerEventMatrix->Projection(0,1),fRecEventMatrix->Projection(0,1),"zv_mult_event_corr_matrix_2D");
1607 aFolderObj->Add(h2D);
1610 h = AlidNdPtHelper::GenerateCorrMatrix(fTriggerEventMatrix->Projection(1),fRecEventMatrix->Projection(1),"mult_event_corr_matrix");
1614 h = AlidNdPtHelper::GenerateCorrMatrix(fTriggerEventMatrix->Projection(0),fRecEventMatrix->Projection(0),"zv_event_corr_matrix");
1620 // track-event trigger bias correction (MB to INEL)
1622 hs = AlidNdPtHelper::GenerateCorrMatrix(fGenTrackEventMatrix,fTriggerTrackEventMatrix,"zv_pt_eta_track_trig_MBtoInel_corr_matrix");
1623 aFolderObj->Add(hs);
1625 h2D = AlidNdPtHelper::GenerateCorrMatrix(fGenTrackEventMatrix->Projection(1,2),fTriggerTrackEventMatrix->Projection(1,2),"eta_pt_track_trig_MBtoInel_corr_matrix");
1626 aFolderObj->Add(h2D);
1628 h2D = AlidNdPtHelper::GenerateCorrMatrix(fGenTrackEventMatrix->Projection(1,0),fTriggerTrackEventMatrix->Projection(1,0),"pt_zv_track_trig_MBtoInel_corr_matrix");
1629 aFolderObj->Add(h2D);
1631 h2D = AlidNdPtHelper::GenerateCorrMatrix(fGenTrackEventMatrix->Projection(2,0),fTriggerTrackEventMatrix->Projection(2,0),"zv_eta_track_trig_MBtoInel_corr_matrix");
1632 aFolderObj->Add(h2D);
1636 h = AlidNdPtHelper::GenerateCorrMatrix(fTriggerTrackEventMatrix->Projection(1),fGenTrackEventMatrix->Projection(1),"pt_track_trig_MBtoInel_eff_matrix");
1641 // track-event vertex reconstruction correction (MB)
1643 hs = AlidNdPtHelper::GenerateCorrMatrix(fTriggerTrackEventMatrix,fRecTrackEventMatrix,"zv_pt_eta_track_event_corr_matrix");
1644 aFolderObj->Add(hs);
1646 h2D = AlidNdPtHelper::GenerateCorrMatrix(fTriggerTrackEventMatrix->Projection(1,2),fRecTrackEventMatrix->Projection(1,2),"eta_pt_track_event_corr_matrix");
1647 aFolderObj->Add(h2D);
1649 h2D = AlidNdPtHelper::GenerateCorrMatrix(fTriggerTrackEventMatrix->Projection(1,0),fRecTrackEventMatrix->Projection(1,0),"pt_zv_track_event_corr_matrix");
1650 aFolderObj->Add(h2D);
1652 h2D = AlidNdPtHelper::GenerateCorrMatrix(fTriggerTrackEventMatrix->Projection(2,0),fRecTrackEventMatrix->Projection(2,0),"zv_eta_track_event_corr_matrix");
1653 aFolderObj->Add(h2D);
1657 h = AlidNdPtHelper::GenerateCorrMatrix(fRecTrackEventMatrix->Projection(1),fTriggerTrackEventMatrix->Projection(1),"pt_track_event_eff_matrix");
1662 // track rec. efficiency correction
1664 hs = AlidNdPtHelper::GenerateCorrMatrix(fGenPrimTrackMatrix,fRecPrimTrackMatrix,"zv_pt_eta_track_corr_matrix");
1665 aFolderObj->Add(hs);
1667 h2D = AlidNdPtHelper::GenerateCorrMatrix(fGenPrimTrackMatrix->Projection(1,2),fRecPrimTrackMatrix->Projection(1,2),"eta_pt_track_corr_matrix");
1668 aFolderObj->Add(h2D);
1670 h2D = AlidNdPtHelper::GenerateCorrMatrix(fGenPrimTrackMatrix->Projection(1,0),fRecPrimTrackMatrix->Projection(1,0),"pt_zv_track_corr_matrix");
1671 aFolderObj->Add(h2D);
1673 h2D = AlidNdPtHelper::GenerateCorrMatrix(fGenPrimTrackMatrix->Projection(2,0),fRecPrimTrackMatrix->Projection(2,0),"zv_eta_track_corr_matrix");
1674 aFolderObj->Add(h2D);
1677 h = AlidNdPtHelper::GenerateCorrMatrix(fGenPrimTrackMatrix->Projection(0),fRecPrimTrackMatrix->Projection(0),"zv_track_corr_matrix");
1680 h = AlidNdPtHelper::GenerateCorrMatrix(fGenPrimTrackMatrix->Projection(1),fRecPrimTrackMatrix->Projection(1),"pt_track_corr_matrix");
1685 h = AlidNdPtHelper::GenerateCorrMatrix(fRecPrimTrackMatrix->Projection(1), fGenPrimTrackMatrix->Projection(1),"pt_track_eff_matrix");
1688 h = AlidNdPtHelper::GenerateCorrMatrix(fGenPrimTrackMatrix->Projection(2),fRecPrimTrackMatrix->Projection(2),"eta_track_corr_matrix");
1692 // secondary track contamination correction
1694 //hs = AlidNdPtHelper::GenerateContCorrMatrix(fRecSecTrackMatrix,fRecTrackMatrix,"zv_pt_eta_track_cont_matrix");
1695 hs = AlidNdPtHelper::GenerateCorrMatrix(fRecSecTrackMatrix,fRecTrackMatrix,"zv_pt_eta_track_cont_matrix");
1696 aFolderObj->Add(hs);
1698 h2D = AlidNdPtHelper::GenerateCorrMatrix(fRecSecTrackMatrix->Projection(1,2),fRecTrackMatrix->Projection(1,2),"eta_pt_track_cont_matrix");
1699 aFolderObj->Add(h2D);
1701 h2D = AlidNdPtHelper::GenerateCorrMatrix(fRecSecTrackMatrix->Projection(1,0),fRecTrackMatrix->Projection(1,0),"pt_zv_track_cont_matrix");
1702 aFolderObj->Add(h2D);
1704 h2D = AlidNdPtHelper::GenerateCorrMatrix(fRecSecTrackMatrix->Projection(2,0),fRecTrackMatrix->Projection(2,0),"zv_eta_track_cont_matrix");
1705 aFolderObj->Add(h2D);
1707 h = AlidNdPtHelper::GenerateCorrMatrix(fRecSecTrackMatrix->Projection(0),fRecTrackMatrix->Projection(0),"zv_track_cont_matrix");
1711 h = AlidNdPtHelper::GenerateCorrMatrix(fRecSecTrackMatrix->Projection(1),fRecTrackMatrix->Projection(1),"pt_track_cont_matrix");
1714 h = AlidNdPtHelper::GenerateCorrMatrix(fRecSecTrackMatrix->Projection(2),fRecTrackMatrix->Projection(2),"eta_track_cont_matrix");
1718 // multiple track reconstruction correction
1720 //hs = AlidNdPtHelper::GenerateContCorrMatrix(fRecMultTrackMatrix,fRecTrackMatrix,"zv_pt_eta_mult_track_cont_matrix");
1721 hs = AlidNdPtHelper::GenerateCorrMatrix(fRecMultTrackMatrix,fRecTrackMatrix,"zv_pt_eta_mult_track_cont_matrix");
1722 aFolderObj->Add(hs);
1724 h2D = AlidNdPtHelper::GenerateCorrMatrix(fRecMultTrackMatrix->Projection(1,2),fRecTrackMatrix->Projection(1,2),"eta_pt_mult_track_cont_matrix");
1725 aFolderObj->Add(h2D);
1727 h2D = AlidNdPtHelper::GenerateCorrMatrix(fRecMultTrackMatrix->Projection(1,0),fRecTrackMatrix->Projection(1,0),"pt_zv_mult_track_cont_matrix");
1728 aFolderObj->Add(h2D);
1730 h2D = AlidNdPtHelper::GenerateCorrMatrix(fRecMultTrackMatrix->Projection(2,0),fRecTrackMatrix->Projection(2,0),"zv_eta_mult_track_cont_matrix");
1731 aFolderObj->Add(h2D);
1733 h = AlidNdPtHelper::GenerateCorrMatrix(fRecMultTrackMatrix->Projection(0),fRecTrackMatrix->Projection(0),"zv_mult_track_cont_matrix");
1736 h = AlidNdPtHelper::GenerateCorrMatrix(fRecMultTrackMatrix->Projection(1),fRecTrackMatrix->Projection(1),"pt_mult_track_cont_matrix");
1739 h = AlidNdPtHelper::GenerateCorrMatrix(fRecMultTrackMatrix->Projection(2),fRecTrackMatrix->Projection(2),"eta_mult_track_cont_matrix");
1743 // Control histograms
1748 // Efficiency electrons, muons, pions, kaons, protons, all
1749 fMCPrimTrackHist1[1]->GetAxis(2)->SetRange(1,1);
1750 fMCPrimTrackHist1[2]->GetAxis(2)->SetRange(1,1);
1751 h1 = fMCPrimTrackHist1[1]->Projection(0);
1752 h2 = fMCPrimTrackHist1[2]->Projection(0);
1753 h2c = (TH1D *)h2->Clone();
1755 h2c->SetName("eff_pt_electrons");
1756 aFolderObj->Add(h2c);
1758 fMCPrimTrackHist1[1]->GetAxis(2)->SetRange(2,2);
1759 fMCPrimTrackHist1[2]->GetAxis(2)->SetRange(2,2);
1760 h1 = fMCPrimTrackHist1[1]->Projection(0);
1761 h2 = fMCPrimTrackHist1[2]->Projection(0);
1762 h2c = (TH1D *)h2->Clone();
1764 h2c->SetName("eff_pt_muons");
1765 aFolderObj->Add(h2c);
1767 fMCPrimTrackHist1[1]->GetAxis(2)->SetRange(3,3);
1768 fMCPrimTrackHist1[2]->GetAxis(2)->SetRange(3,3);
1769 h1 = fMCPrimTrackHist1[1]->Projection(0);
1770 h2 = fMCPrimTrackHist1[2]->Projection(0);
1771 h2c = (TH1D *)h2->Clone();
1773 h2c->SetName("eff_pt_pions");
1774 aFolderObj->Add(h2c);
1776 fMCPrimTrackHist1[1]->GetAxis(2)->SetRange(4,4);
1777 fMCPrimTrackHist1[2]->GetAxis(2)->SetRange(4,4);
1778 h1 = fMCPrimTrackHist1[1]->Projection(0);
1779 h2 = fMCPrimTrackHist1[2]->Projection(0);
1780 h2c = (TH1D *)h2->Clone();
1782 h2c->SetName("eff_pt_kaons");
1783 aFolderObj->Add(h2c);
1785 fMCPrimTrackHist1[1]->GetAxis(2)->SetRange(5,5);
1786 fMCPrimTrackHist1[2]->GetAxis(2)->SetRange(5,5);
1787 h1 = fMCPrimTrackHist1[1]->Projection(0);
1788 h2 = fMCPrimTrackHist1[2]->Projection(0);
1789 h2c = (TH1D *)h2->Clone();
1791 h2c->SetName("eff_pt_protons");
1792 aFolderObj->Add(h2c);
1794 fMCPrimTrackHist1[1]->GetAxis(2)->SetRange(1,5);
1795 fMCPrimTrackHist1[2]->GetAxis(2)->SetRange(1,5);
1796 h1 = fMCPrimTrackHist1[1]->Projection(0);
1797 h2 = fMCPrimTrackHist1[2]->Projection(0);
1798 h2c = (TH1D *)h2->Clone();
1800 h2c->SetName("eff_pt_selected");
1801 aFolderObj->Add(h2c);
1803 fMCPrimTrackHist1[1]->GetAxis(2)->SetRange(1,6);
1804 fMCPrimTrackHist1[2]->GetAxis(2)->SetRange(1,6);
1805 h1 = fMCPrimTrackHist1[1]->Projection(0);
1806 h2 = fMCPrimTrackHist1[2]->Projection(0);
1807 h2c = (TH1D *)h2->Clone();
1809 h2c->SetName("eff_pt_all");
1810 aFolderObj->Add(h2c);
1812 fMCPrimTrackHist1[1]->GetAxis(1)->SetRange(1,fMCPrimTrackHist1[1]->GetAxis(1)->GetNbins());
1813 fMCPrimTrackHist1[2]->GetAxis(1)->SetRange(1,fMCPrimTrackHist1[2]->GetAxis(1)->GetNbins());
1816 // - rec, primaries, secondaries
1817 // - primaries (pid)
1818 // - secondaries (pid)
1819 // - secondaries (mech)
1820 // - secondaries (mother)
1823 TH1D *mcPtAccall = fMCTrackHist1[1]->Projection(0);
1824 mcPtAccall->SetName("mc_pt_acc_all");
1825 aFolderObj->Add(mcPtAccall);
1827 TH1D *mcPtAccprim = fMCPrimTrackHist1[1]->Projection(0);
1828 mcPtAccprim->SetName("mc_pt_acc_prim");
1829 aFolderObj->Add(mcPtAccprim);
1831 TH1D *mcPtRecall = fMCTrackHist1[2]->Projection(0);
1832 mcPtRecall->SetName("mc_pt_rec_all");
1833 aFolderObj->Add(mcPtRecall);
1835 TH1D *mcPtRecprim = fMCPrimTrackHist1[2]->Projection(0);
1836 mcPtRecprim->SetName("mc_pt_rec_prim");
1837 aFolderObj->Add(mcPtRecprim);
1839 TH1D *mcPtRecsec = fMCSecTrackHist1[2]->Projection(0);
1840 mcPtRecsec->SetName("mc_pt_rec_sec");
1841 aFolderObj->Add(mcPtRecsec);
1843 for(Int_t i = 0; i<6; i++)
1845 sprintf(name,"mc_pt_rec_prim_pid_%d",i);
1846 fMCPrimTrackHist1[2]->GetAxis(2)->SetRange(i+1,i+1);
1847 h = fMCPrimTrackHist1[2]->Projection(0);
1851 sprintf(name,"mc_pt_rec_sec_pid_%d",i);
1852 fMCSecTrackHist1[2]->GetAxis(2)->SetRange(i+1,i+1);
1853 h = fMCSecTrackHist1[2]->Projection(0);
1857 // production mechanisms for given pid
1858 fMCSecTrackHist1[2]->GetAxis(2)->SetRange(i+1,i+1);
1860 for(Int_t j=0; j<20; j++) {
1864 sprintf(name,"mc_pt_rec_sec_pid_%d_decay",i);
1865 fMCSecTrackHist1[2]->GetAxis(3)->SetRange(j+1,j+1);
1866 h = fMCSecTrackHist1[2]->Projection(0);
1870 sprintf(name,"mc_eta_rec_sec_pid_%d_decay",i);
1871 fMCSecTrackHist1[2]->GetAxis(3)->SetRange(j+1,j+1);
1872 h = fMCSecTrackHist1[2]->Projection(1);
1876 sprintf(name,"mc_mother_rec_sec_pid_%d_decay",i);
1877 fMCSecTrackHist1[2]->GetAxis(3)->SetRange(j+1,j+1);
1878 h = fMCSecTrackHist1[2]->Projection(4);
1882 } else if (j == 5) {
1885 sprintf(name,"mc_pt_rec_sec_pid_%d_conv",i);
1886 fMCSecTrackHist1[2]->GetAxis(3)->SetRange(j+1,j+1);
1887 h = fMCSecTrackHist1[2]->Projection(0);
1891 sprintf(name,"mc_eta_rec_sec_pid_%d_conv",i);
1892 fMCSecTrackHist1[2]->GetAxis(3)->SetRange(j+1,j+1);
1893 h = fMCSecTrackHist1[2]->Projection(1);
1897 sprintf(name,"mc_mother_rec_sec_pid_%d_conv",i);
1898 fMCSecTrackHist1[2]->GetAxis(3)->SetRange(j+1,j+1);
1899 h = fMCSecTrackHist1[2]->Projection(4);
1903 } else if (j == 13) {
1906 sprintf(name,"mc_pt_rec_sec_pid_%d_mat",i);
1907 fMCSecTrackHist1[2]->GetAxis(3)->SetRange(j+1,j+1);
1908 h = fMCSecTrackHist1[2]->Projection(0);
1912 sprintf(name,"mc_eta_rec_sec_pid_%d_mat",i);
1913 fMCSecTrackHist1[2]->GetAxis(3)->SetRange(j+1,j+1);
1914 h = fMCSecTrackHist1[2]->Projection(1);
1918 sprintf(name,"mc_eta_mother_rec_sec_pid_%d_mat",i);
1919 fMCSecTrackHist1[2]->GetAxis(3)->SetRange(j+1,j+1);
1920 h = fMCSecTrackHist1[2]->Projection(4,1);
1924 sprintf(name,"mc_mother_rec_sec_pid_%d_mat",i);
1925 fMCSecTrackHist1[2]->GetAxis(3)->SetRange(j+1,j+1);
1926 h = fMCSecTrackHist1[2]->Projection(4);
1930 sprintf(name,"mc_pt_mother_rec_sec_pid_%d_mat",i);
1931 fMCSecTrackHist1[2]->GetAxis(3)->SetRange(j+1,j+1);
1932 h = fMCSecTrackHist1[2]->Projection(4,0);
1942 } // end fHistogramOn
1945 // resolution histograms
1946 // only for reconstructed tracks
1950 TCanvas * c = new TCanvas("resol","resol");
1954 fRecMCTrackHist1->GetAxis(1)->SetRangeUser(-0.8,0.79);
1956 h2F = (TH2F*)fRecMCTrackHist1->Projection(2,0);
1957 h = AlidNdPtHelper::MakeResol(h2F,1,0,kTRUE,10);
1958 h->SetXTitle("p_{tmc} (GeV/c)");
1959 h->SetYTitle("(p_{t}-p_{tmc})/p_{tmc} resolution");
1961 h->SetName("pt_resolution_vs_mcpt");
1964 h2F = (TH2F*)fRecMCTrackHist1->Projection(2,0);
1965 h = AlidNdPtHelper::MakeResol(h2F,1,1,kTRUE,10);
1966 h->SetXTitle("p_{tmc} (GeV/c)");
1967 h->SetYTitle("(p_{t}-p_{tmc})/p_{tmc} mean");
1969 h->SetName("dpt_mean_vs_mcpt");
1973 h2F = (TH2F*)fRecMCTrackHist1->Projection(3,0);
1974 h = AlidNdPtHelper::MakeResol(h2F,1,0,kTRUE,10);
1975 h->SetXTitle("p_{tmc} (GeV/c)");
1976 h->SetYTitle("(#eta-#eta_{mc}) resolution");
1978 h->SetName("eta_resolution_vs_mcpt");
1981 h2F = (TH2F*)fRecMCTrackHist1->Projection(3,0);
1982 h = AlidNdPtHelper::MakeResol(h2F,1,1,kTRUE,10);
1983 h->SetXTitle("p_{tmc} (GeV/c)");
1984 h->SetYTitle("(#eta-mc#eta) mean");
1986 h->SetName("deta_mean_vs_mcpt");
1990 fRecMCTrackHist1->GetAxis(1)->SetRange(1,fRecMCTrackHist1->GetAxis(1)->GetNbins());
1992 h2F = (TH2F*)fRecMCTrackHist1->Projection(2,1);
1993 h = AlidNdPtHelper::MakeResol(h2F,1,0,kTRUE,10);
1994 h->SetXTitle("#eta_{mc}");
1995 h->SetYTitle("(p_{t}-p_{tmc})/p_{tmc} resolution");
1997 h->SetName("pt_resolution_vs_mceta");
2000 h2F = (TH2F*)fRecMCTrackHist1->Projection(2,1);
2001 h = AlidNdPtHelper::MakeResol(h2F,1,1,kTRUE,10);
2002 h->SetXTitle("#eta_{mc}");
2003 h->SetYTitle("(p_{t}-p_{tmc})/p_{tmc} mean");
2005 h->SetName("dpt_mean_vs_mceta");
2009 h2F = (TH2F*)fRecMCTrackHist1->Projection(3,1);
2010 h = AlidNdPtHelper::MakeResol(h2F,1,0,kTRUE,10);
2011 h->SetXTitle("#eta_{mc}");
2012 h->SetYTitle("(#eta-#eta_{mc}) resolution");
2014 h->SetName("eta_resolution_vs_mceta");
2017 h2F = (TH2F*)fRecMCTrackHist1->Projection(3,1);
2018 h = AlidNdPtHelper::MakeResol(h2F,1,1,kTRUE,10);
2019 h->SetXTitle("#eta_{mc}");
2020 h->SetYTitle("(#eta-mc#eta) mean");
2022 h->SetName("deta_mean_vs_mceta");
2025 fRecMCTrackHist1->GetAxis(0)->SetRange(1,fRecMCTrackHist1->GetAxis(0)->GetNbins());
2027 } // end use MC info
2029 // export objects to analysis folder
2030 fAnalysisFolder = ExportToFolder(aFolderObj);
2032 // delete only TObjArray
2033 if(aFolderObj) delete aFolderObj;
2036 //_____________________________________________________________________________
2037 TFolder* AlidNdPtAnalysisPbPb::ExportToFolder(TObjArray * const array)
2039 // recreate folder avery time and export objects to new one
2041 AlidNdPtAnalysisPbPb * comp=this;
2042 TFolder *folder = comp->GetAnalysisFolder();
2044 TString name, title;
2045 TFolder *newFolder = 0;
2047 Int_t size = array->GetSize();
2050 // get name and title from old folder
2051 name = folder->GetName();
2052 title = folder->GetTitle();
2058 newFolder = CreateFolder(name.Data(),title.Data());
2059 newFolder->SetOwner();
2061 // add objects to folder
2063 newFolder->Add(array->At(i));
2071 //_____________________________________________________________________________
2072 TFolder* AlidNdPtAnalysisPbPb::CreateFolder(TString name,TString title) {
2073 // create folder for analysed histograms
2075 TFolder *folder = 0;
2076 folder = new TFolder(name.Data(),title.Data());