]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWG1/TPC/AliPerformanceDCA.cxx
Updating the functionality of AliAnalysisHadEtCorrections to accomodate centrality...
[u/mrichter/AliRoot.git] / PWG1 / TPC / AliPerformanceDCA.cxx
CommitLineData
7cc34f08 1//------------------------------------------------------------------------------
2// Implementation of AliPerformanceDCA class. It keeps information from
3// comparison of reconstructed and MC particle tracks. In addtion,
4// it keeps selection cuts used during comparison. The comparison
5// information is stored in the ROOT histograms. Analysis of these
6// histograms can be done by using Analyse() class function. The result of
7// the analysis (histograms/graphs) are stored in the folder
8// which is a data member of AliPerformanceDCA.
9//
10// Author: J.Otwinowski 04/02/2008
11//------------------------------------------------------------------------------
12
13/*
14
15 // after running comparison task, read the file, and get component
16 gROOT->LoadMacro("$ALICE_ROOT/PWG1/Macros/LoadMyLibs.C");
17 LoadMyLibs();
18 TFile f("Output.root");
19 AliPerformanceDCA * compObj = (AliPerformanceDCA*)coutput->FindObject("AliPerformanceDCA");
20
21 // Analyse comparison data
22 compObj->Analyse();
23
24 // the output histograms/graphs will be stored in the folder "folderDCA"
25 compObj->GetAnalysisFolder()->ls("*");
26
27 // user can save whole comparison object (or only folder with anlysed histograms)
28 // in the seperate output file (e.g.)
29 TFile fout("Analysed_DCA.root","recreate");
30 compObj->Write(); // compObj->GetAnalysisFolder()->Write();
31 fout.Close();
32
33*/
34
35#include <TAxis.h>
36#include <TCanvas.h>
37#include <TGraph.h>
38#include <TGraph2D.h>
39#include <TH1.h>
40#include <TH2.h>
41#include <TH3.h>
4bb5f8a9 42#include <TF1.h>
7cc34f08 43
44#include "AliPerformanceDCA.h"
45#include "AliESDEvent.h"
46#include "AliESDVertex.h"
47#include "AliLog.h"
48#include "AliMathBase.h"
49#include "AliRecInfoCuts.h"
50#include "AliMCInfoCuts.h"
51#include "AliStack.h"
52#include "AliMCEvent.h"
53#include "AliTracker.h"
54#include "AliHeader.h"
55#include "AliGenEventHeader.h"
56
57using namespace std;
58
59ClassImp(AliPerformanceDCA)
60
61//_____________________________________________________________________________
62AliPerformanceDCA::AliPerformanceDCA():
63 AliPerformanceObject("AliPerformanceDCA"),
64
65 // DCA histograms
66 fDCAHisto(0),
67
68 // Cuts
69 fCutsRC(0),
70 fCutsMC(0),
71
72 // histogram folder
73 fAnalysisFolder(0)
74{
75 // default constructor
76 Init();
77}
78
79//_____________________________________________________________________________
80AliPerformanceDCA::AliPerformanceDCA(Char_t* name="AliPerformanceDCA", Char_t* title="AliPerformanceDCA",Int_t analysisMode=0, Bool_t hptGenerator=kFALSE):
81 AliPerformanceObject(name,title),
82
83 // DCA histograms
84 fDCAHisto(0),
85
86 // Cuts
87 fCutsRC(0),
88 fCutsMC(0),
89
90 // histogram folder
91 fAnalysisFolder(0)
92{
93 // named constructor
94
95 SetAnalysisMode(analysisMode);
96 SetHptGenerator(hptGenerator);
97 Init();
98}
99
100//_____________________________________________________________________________
101AliPerformanceDCA::~AliPerformanceDCA()
102{
103 // destructor
104 if(fDCAHisto) delete fDCAHisto; fDCAHisto=0;
105 if(fAnalysisFolder) delete fAnalysisFolder; fAnalysisFolder=0;
106}
107
108//_____________________________________________________________________________
109void AliPerformanceDCA::Init()
110{
111 // DCA histograms
112 Int_t nPtBins = 50;
113 Double_t ptMin = 1.e-2, ptMax = 10.;
114
115 Double_t *binsPt = 0;
116 if (IsHptGenerator()) {
117 nPtBins = 100; ptMax = 100.;
118 binsPt = CreateLogAxis(nPtBins,ptMin,ptMax);
119 } else {
120 binsPt = CreateLogAxis(nPtBins,ptMin,ptMax);
121 }
122
123 /*
124 Int_t nPtBins = 31;
125 Double_t binsPt[32] = {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.7,0.8,0.9,1.0,1.2,1.4,1.6,1.8,2.0,2.25,2.5,2.75,3.,3.5,4.,5.,6.,8.,10.};
126 Double_t ptMin = 0., ptMax = 10.;
127
128 if(IsHptGenerator() == kTRUE) {
129 nPtBins = 100;
130 ptMin = 0.; ptMax = 100.;
131 }
132 */
133
134 //dca_r, dca_z, eta, pt
e6a60a90 135 Int_t binsQA[5] = {100,100,30,nPtBins,144};
886bf4d3 136 Double_t xminQA[5] = {-10.,-10.,-1.5,ptMin,0.};
137 Double_t xmaxQA[5] = {10.,10.,1.5,ptMax,2*TMath::Pi()};
7cc34f08 138
886bf4d3 139 fDCAHisto = new THnSparseF("fDCAHisto","dca_r:dca_z:eta:pt:phi",5,binsQA,xminQA,xmaxQA);
7cc34f08 140 fDCAHisto->SetBinEdges(3,binsPt);
141
142 fDCAHisto->GetAxis(0)->SetTitle("dca_r (cm)");
143 fDCAHisto->GetAxis(1)->SetTitle("dca_z (cm)");
144 fDCAHisto->GetAxis(2)->SetTitle("#eta");
145 fDCAHisto->GetAxis(3)->SetTitle("p_{T} (GeV/c)");
886bf4d3 146 fDCAHisto->GetAxis(4)->SetTitle("phi (rad)");
7cc34f08 147 fDCAHisto->Sumw2();
148
149 // init cuts
150 if(!fCutsMC)
151 AliDebug(AliLog::kError, "ERROR: Cannot find AliMCInfoCuts object");
152 if(!fCutsRC)
153 AliDebug(AliLog::kError, "ERROR: Cannot find AliRecInfoCuts object");
154
155 // init folder
156 fAnalysisFolder = CreateFolder("folderDCA","Analysis DCA Folder");
157}
158
159//_____________________________________________________________________________
758320f7 160void AliPerformanceDCA::ProcessTPC(AliStack* const stack, AliESDtrack *const esdTrack, AliESDEvent* const esdEvent)
7cc34f08 161{
162 // Fill DCA comparison information
758320f7 163 if(!esdEvent) return;
7cc34f08 164 if(!esdTrack) return;
165
758320f7 166 if( IsUseTrackVertex() )
167 {
168 // Relate TPC inner params to prim. vertex
169 const AliESDVertex *vtxESD = esdEvent->GetPrimaryVertexTracks();
170 Double_t x[3]; esdTrack->GetXYZ(x);
171 Double_t b[3]; AliTracker::GetBxByBz(x,b);
172 Bool_t isOK = esdTrack->RelateToVertexTPCBxByBz(vtxESD, b, kVeryBig);
173 if(!isOK) return;
174
175 /*
176 // JMT -- recaluclate DCA for HLT if not present
177 if ( dca[0] == 0. && dca[1] == 0. ) {
178 track->GetDZ( vtxESD->GetX(), vtxESD->GetY(), vtxESD->GetZ(), esdEvent->GetMagneticField(), dca );
179 }
180 */
181 }
182
183 // get TPC inner params at DCA to prim. vertex
7cc34f08 184 const AliExternalTrackParam *track = esdTrack->GetTPCInnerParam();
185 if(!track) return;
186
758320f7 187 // read from ESD track
7cc34f08 188 Float_t dca[2], cov[3]; // dca_xy, dca_z, sigma_xy, sigma_xy_z, sigma_z
189 esdTrack->GetImpactParametersTPC(dca,cov);
190
191 if (esdTrack->GetTPCNcls()<fCutsRC->GetMinNClustersTPC()) return; // min. nb. TPC clusters
192
886bf4d3 193 Double_t vDCAHisto[5]={dca[0],dca[1],track->Eta(),track->Pt(),track->Phi()};
7cc34f08 194 fDCAHisto->Fill(vDCAHisto);
195
196 //
197 // Fill rec vs MC information
198 //
199 if(!stack) return;
200
758320f7 201}
7cc34f08 202
203//_____________________________________________________________________________
758320f7 204void AliPerformanceDCA::ProcessTPCITS(AliStack* const stack, AliESDtrack *const esdTrack, AliESDEvent* const esdEvent)
7cc34f08 205{
206 // Fill DCA comparison information
207 if(!esdTrack) return;
758320f7 208 if(!esdEvent) return;
209
210 if( IsUseTrackVertex() )
211 {
212 // Relate TPC inner params to prim. vertex
213 const AliESDVertex *vtxESD = esdEvent->GetPrimaryVertexTracks();
214 Double_t x[3]; esdTrack->GetXYZ(x);
215 Double_t b[3]; AliTracker::GetBxByBz(x,b);
216 Bool_t isOK = esdTrack->RelateToVertexBxByBz(vtxESD, b, kVeryBig);
217 if(!isOK) return;
218
219 /*
220 // JMT -- recaluclate DCA for HLT if not present
221 if ( dca[0] == 0. && dca[1] == 0. ) {
222 track->GetDZ( vtxESD->GetX(), vtxESD->GetY(), vtxESD->GetZ(), esdEvent->GetMagneticField(), dca );
223 }
224 */
225 }
7cc34f08 226
227 Float_t dca[2], cov[3]; // dca_xy, dca_z, sigma_xy, sigma_xy_z, sigma_z
228 esdTrack->GetImpactParameters(dca,cov);
229
230 if ((esdTrack->GetStatus()&AliESDtrack::kTPCrefit)==0) return; // TPC refit
231 if (esdTrack->GetTPCNcls()<fCutsRC->GetMinNClustersTPC()) return; // min. nb. TPC clusters
faa3b211 232 if(esdTrack->GetITSclusters(0)<fCutsRC->GetMinNClustersITS()) return; // min. nb. ITS clusters
7cc34f08 233
886bf4d3 234 Double_t vDCAHisto[5]={dca[0],dca[1],esdTrack->Eta(),esdTrack->Pt(), esdTrack->Phi()};
7cc34f08 235 fDCAHisto->Fill(vDCAHisto);
236
237 //
238 // Fill rec vs MC information
239 //
240 if(!stack) return;
241
242}
243
244void AliPerformanceDCA::ProcessConstrained(AliStack* const /*stack*/, AliESDtrack *const /*esdTrack*/)
245{
246 // Fill DCA comparison information
247
248 AliDebug(AliLog::kWarning, "Warning: Not implemented");
249}
250
251//_____________________________________________________________________________
252Long64_t AliPerformanceDCA::Merge(TCollection* const list)
253{
254 // Merge list of objects (needed by PROOF)
255
256 if (!list)
257 return 0;
258
259 if (list->IsEmpty())
260 return 1;
261
262 TIterator* iter = list->MakeIterator();
263 TObject* obj = 0;
264
265 // collection of generated histograms
266 Int_t count=0;
267 while((obj = iter->Next()) != 0)
268 {
269 AliPerformanceDCA* entry = dynamic_cast<AliPerformanceDCA*>(obj);
270 if (entry == 0) continue;
271
272 fDCAHisto->Add(entry->fDCAHisto);
273 count++;
274 }
275
276return count;
277}
278
279//_____________________________________________________________________________
280void AliPerformanceDCA::Exec(AliMCEvent* const mcEvent, AliESDEvent *const esdEvent, AliESDfriend *const esdFriend, const Bool_t bUseMC, const Bool_t bUseESDfriend)
281{
282 // Process comparison information
283 //
284 if(!esdEvent)
285 {
286 Error("Exec","esdEvent not available");
287 return;
288 }
289 AliHeader* header = 0;
290 AliGenEventHeader* genHeader = 0;
291 AliStack* stack = 0;
292 TArrayF vtxMC(3);
293
294 if(bUseMC)
295 {
296 if(!mcEvent) {
297 Error("Exec","mcEvent not available");
298 return;
299 }
300 // get MC event header
301 header = mcEvent->Header();
302 if (!header) {
303 Error("Exec","Header not available");
304 return;
305 }
306 // MC particle stack
307 stack = mcEvent->Stack();
308 if (!stack) {
309 Error("Exec","Stack not available");
310 return;
311 }
312 // get MC vertex
313 genHeader = header->GenEventHeader();
314 if (!genHeader) {
315 Error("Exec","Could not retrieve genHeader from Header");
316 return;
317 }
318 genHeader->PrimaryVertex(vtxMC);
319 }
320
321 // use ESD friends
322 if(bUseESDfriend) {
323 if(!esdFriend) {
324 Error("Exec","esdFriend not available");
325 return;
326 }
327 }
328
e6a60a90 329 // trigger
046b4280 330 if(!bUseMC &&GetTriggerClass()) {
4bb5f8a9 331 Bool_t isEventTriggered = esdEvent->IsTriggerClassFired(GetTriggerClass());
332 if(!isEventTriggered) return;
333 }
e6a60a90 334
758320f7 335 // get event vertex
336 const AliESDVertex *vtxESD = NULL;
337 if( IsUseTrackVertex() )
338 {
339 // track vertex
340 vtxESD = esdEvent->GetPrimaryVertexTracks();
341 }
342 else {
343 // TPC track vertex
344 vtxESD = esdEvent->GetPrimaryVertexTPC();
345 }
e6a60a90 346 if(vtxESD && (vtxESD->GetStatus()<=0)) return;
347
7cc34f08 348 // Process events
349 for (Int_t iTrack = 0; iTrack < esdEvent->GetNumberOfTracks(); iTrack++)
350 {
351 AliESDtrack *track = esdEvent->GetTrack(iTrack);
352 if(!track) continue;
353
758320f7 354 if(GetAnalysisMode() == 0) ProcessTPC(stack,track,esdEvent);
355 else if(GetAnalysisMode() == 1) ProcessTPCITS(stack,track,esdEvent);
7cc34f08 356 else if(GetAnalysisMode() == 2) ProcessConstrained(stack,track);
357 else {
358 printf("ERROR: AnalysisMode %d \n",fAnalysisMode);
359 return;
360 }
361 }
362}
363
364//_____________________________________________________________________________
365void AliPerformanceDCA::Analyse()
366{
367 //
368 // Analyse comparison information and store output histograms
369 // in the analysis folder "folderDCA"
370 //
371
372 TH1::AddDirectory(kFALSE);
4bb5f8a9 373 TH1F *h1D=0;
374 TH2F *h2D=0;
7cc34f08 375 TObjArray *aFolderObj = new TObjArray;
18e588e9 376 if(!aFolderObj) return;
7cc34f08 377 char title[256];
4bb5f8a9 378 TObjArray *arr[6] = {0};
379 TF1 *f1[6] = {0};
7cc34f08 380
18e588e9 381
382
7cc34f08 383 // set pt measurable range
a13f4653 384 //fDCAHisto->GetAxis(3)->SetRangeUser(0.10,10.);
7cc34f08 385
386 //
4bb5f8a9 387 h2D = (TH2F*)fDCAHisto->Projection(0,1); // inverse projection convention
7cc34f08 388 h2D->SetName("dca_r_vs_dca_z");
389 h2D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(1)->GetTitle());
390 h2D->GetYaxis()->SetTitle(fDCAHisto->GetAxis(0)->GetTitle());
18e588e9 391 snprintf(title,256,"%s vs %s",fDCAHisto->GetAxis(0)->GetTitle(),fDCAHisto->GetAxis(1)->GetTitle());
7cc34f08 392 h2D->SetTitle(title);
393 aFolderObj->Add(h2D);
394
395 //
4bb5f8a9 396 h2D = (TH2F*)fDCAHisto->Projection(0,2);
7cc34f08 397 h2D->SetName("dca_r_vs_eta");
398 h2D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(2)->GetTitle());
399 h2D->GetYaxis()->SetTitle(fDCAHisto->GetAxis(0)->GetTitle());
18e588e9 400 snprintf(title,256,"%s vs %s",fDCAHisto->GetAxis(2)->GetTitle(),fDCAHisto->GetAxis(0)->GetTitle());
7cc34f08 401 h2D->SetTitle(title);
402 aFolderObj->Add(h2D);
403
4bb5f8a9 404 //
405 // mean and rms
406 //
7cc34f08 407 h1D = MakeStat1D(h2D,0,0);
408 h1D->SetName("mean_dca_r_vs_eta");
409 h1D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(2)->GetTitle());
410 h1D->GetYaxis()->SetTitle("mean_dca_r (cm)");
18e588e9 411 snprintf(title,256," mean_dca_r (cm) vs %s",fDCAHisto->GetAxis(2)->GetTitle());
7cc34f08 412 h1D->SetTitle(title);
413 aFolderObj->Add(h1D);
414
415 h1D = MakeStat1D(h2D,0,1);
416 h1D->SetName("rms_dca_r_vs_eta");
417 h1D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(2)->GetTitle());
418 h1D->GetYaxis()->SetTitle("rms_dca_r (cm)");
18e588e9 419 snprintf(title,256," rms_dca_r (cm) vs %s",fDCAHisto->GetAxis(2)->GetTitle());
7cc34f08 420 h1D->SetTitle(title);
421 aFolderObj->Add(h1D);
422
423 //
4bb5f8a9 424 // fit mean and sigma
425 //
426
427 arr[0] = new TObjArray();
428 f1[0] = new TF1("gaus","gaus");
429 h2D->FitSlicesY(f1[0],0,-1,0,"QNR",arr[0]);
430
431 h1D = (TH1F*)arr[0]->At(1);
432 h1D->SetName("fit_mean_dca_r_vs_eta");
433 h1D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(2)->GetTitle());
434 h1D->GetYaxis()->SetTitle("fit_mean_dca_r (cm)");
18e588e9 435 snprintf(title,256," fit_mean_dca_r (cm) vs %s",fDCAHisto->GetAxis(2)->GetTitle());
4bb5f8a9 436 h1D->SetTitle(title);
437 aFolderObj->Add(h1D);
438
439 h1D = (TH1F*)arr[0]->At(2);
440 h1D->SetName("res_dca_r_vs_eta");
441 h1D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(2)->GetTitle());
442 h1D->GetYaxis()->SetTitle("res_dca_r (cm)");
18e588e9 443 snprintf(title,256," res_dca_r (cm) vs %s",fDCAHisto->GetAxis(2)->GetTitle());
4bb5f8a9 444 h1D->SetTitle(title);
445 aFolderObj->Add(h1D);
446
447 //
448 //
449 //
450 h2D = (TH2F*)fDCAHisto->Projection(0,3);
7cc34f08 451 h2D->SetName("dca_r_vs_pt");
452 h2D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(3)->GetTitle());
453 h2D->GetYaxis()->SetTitle(fDCAHisto->GetAxis(0)->GetTitle());
18e588e9 454 snprintf(title,256,"%s vs %s",fDCAHisto->GetAxis(0)->GetTitle(),fDCAHisto->GetAxis(3)->GetTitle());
7cc34f08 455 h2D->SetTitle(title);
456 h2D->SetBit(TH1::kLogX);
457 aFolderObj->Add(h2D);
458
459 h1D = MakeStat1D(h2D,0,0);
460 h1D->SetName("mean_dca_r_vs_pt");
461 h1D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(3)->GetTitle());
462 h1D->GetYaxis()->SetTitle("mean_dca_r (cm)");
18e588e9 463 snprintf(title,256,"mean_dca_r (cm) vs %s",fDCAHisto->GetAxis(3)->GetTitle());
7cc34f08 464 h1D->SetTitle(title);
465 h1D->SetBit(TH1::kLogX);
466 aFolderObj->Add(h1D);
467
468 h1D = MakeStat1D(h2D,0,1);
469 h1D->SetName("rms_dca_r_vs_pt");
470 h1D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(3)->GetTitle());
471 h1D->GetYaxis()->SetTitle("rms_dca_r (cm)");
18e588e9 472 snprintf(title,256,"rms_dca_r (cm) vs %s",fDCAHisto->GetAxis(3)->GetTitle());
7cc34f08 473 h1D->SetTitle(title);
474 h1D->SetBit(TH1::kLogX);
475 aFolderObj->Add(h1D);
476
4bb5f8a9 477 //
478 // fit mean and sigma
479 //
480
481 arr[1] = new TObjArray();
482 f1[1] = new TF1("gaus","gaus");
483 h2D->FitSlicesY(f1[1],0,-1,0,"QNR",arr[1]);
484
485 h1D = (TH1F*)arr[1]->At(1);
486 h1D->SetName("fit_mean_dca_r_vs_pt");
487 h1D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(3)->GetTitle());
488 h1D->GetYaxis()->SetTitle("fit_mean_dca_r (cm)");
18e588e9 489 snprintf(title,256,"fit_mean_dca_r (cm) vs %s",fDCAHisto->GetAxis(3)->GetTitle());
4bb5f8a9 490 h1D->SetTitle(title);
491 h1D->SetBit(TH1::kLogX);
492 aFolderObj->Add(h1D);
493
494 h1D = (TH1F*)arr[1]->At(2);
495 h1D->SetName("res_dca_r_vs_pt");
496 h1D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(3)->GetTitle());
497 h1D->GetYaxis()->SetTitle("res_dca_r (cm)");
18e588e9 498 snprintf(title,256,"res_dca_r (cm) vs %s",fDCAHisto->GetAxis(3)->GetTitle());
4bb5f8a9 499 h1D->SetTitle(title);
500 h1D->SetBit(TH1::kLogX);
501 aFolderObj->Add(h1D);
502
7cc34f08 503 //
4bb5f8a9 504 h2D = (TH2F*)fDCAHisto->Projection(1,2);
7cc34f08 505 h2D->SetName("dca_z_vs_eta");
506 h2D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(2)->GetTitle());
507 h2D->GetYaxis()->SetTitle(fDCAHisto->GetAxis(1)->GetTitle());
18e588e9 508 snprintf(title,256,"%s vs %s",fDCAHisto->GetAxis(1)->GetTitle(),fDCAHisto->GetAxis(2)->GetTitle());
7cc34f08 509 h2D->SetTitle(title);
510 aFolderObj->Add(h2D);
511
512 h1D = MakeStat1D(h2D,0,0);
513 h1D->SetName("mean_dca_z_vs_eta");
514 h1D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(2)->GetTitle());
515 h1D->GetYaxis()->SetTitle("mean_dca_z (cm)");
18e588e9 516 snprintf(title,256,"mean_dca_z (cm) vs %s",fDCAHisto->GetAxis(2)->GetTitle());
7cc34f08 517 h1D->SetTitle(title);
518 aFolderObj->Add(h1D);
519
520 h1D = MakeStat1D(h2D,0,1);
521 h1D->SetName("rms_dca_z_vs_eta");
522 h1D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(2)->GetTitle());
523 h1D->GetYaxis()->SetTitle("rms_dca_z (cm)");
18e588e9 524 snprintf(title,256,"rms_dca_z (cm) vs %s",fDCAHisto->GetAxis(2)->GetTitle());
7cc34f08 525 h1D->SetTitle(title);
526 aFolderObj->Add(h1D);
527
528 //
4bb5f8a9 529 // fit mean and sigma
530 //
531 arr[2] = new TObjArray();
532 f1[2] = new TF1("gaus","gaus");
533 h2D->FitSlicesY(f1[2],0,-1,0,"QNR",arr[2]);
534
535 h1D = (TH1F*)arr[2]->At(1);
536 h1D->SetName("fit_mean_dca_z_vs_eta");
537 h1D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(2)->GetTitle());
538 h1D->GetYaxis()->SetTitle("fit_mean_dca_z (cm)");
18e588e9 539 snprintf(title,256,"fit_mean_dca_z (cm) vs %s",fDCAHisto->GetAxis(2)->GetTitle());
4bb5f8a9 540 h1D->SetTitle(title);
541 aFolderObj->Add(h1D);
542
543 h1D = (TH1F*)arr[2]->At(2);
544 h1D->SetName("res_dca_z_vs_eta");
545 h1D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(2)->GetTitle());
546 h1D->GetYaxis()->SetTitle("res_dca_z (cm)");
18e588e9 547 snprintf(title,256,"res_dca_z (cm) vs %s",fDCAHisto->GetAxis(2)->GetTitle());
4bb5f8a9 548 h1D->SetTitle(title);
549 aFolderObj->Add(h1D);
550
551 //
552 h2D = (TH2F*)fDCAHisto->Projection(1,3);
7cc34f08 553 h2D->SetName("dca_z_vs_pt");
554 h2D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(3)->GetTitle());
555 h2D->GetYaxis()->SetTitle(fDCAHisto->GetAxis(1)->GetTitle());
18e588e9 556 snprintf(title,256,"%s vs %s",fDCAHisto->GetAxis(1)->GetTitle(),fDCAHisto->GetAxis(3)->GetTitle());
7cc34f08 557 h2D->SetTitle(title);
558 h2D->SetBit(TH1::kLogX);
559 aFolderObj->Add(h2D);
560
561 h1D = MakeStat1D(h2D,0,0);
562 h1D->SetName("mean_dca_z_vs_pt");
563 h1D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(3)->GetTitle());
564 h1D->GetYaxis()->SetTitle("mean_dca_z (cm)");
18e588e9 565 snprintf(title,256,"mean_dca_z (cm) vs %s",fDCAHisto->GetAxis(3)->GetTitle());
7cc34f08 566 h1D->SetTitle(title);
567 h1D->SetBit(TH1::kLogX);
568 aFolderObj->Add(h1D);
569
570 h1D = MakeStat1D(h2D,0,1);
571 h1D->SetName("rms_dca_z_vs_pt");
572 h1D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(3)->GetTitle());
573 h1D->GetYaxis()->SetTitle("rms_dca_z (cm)");
18e588e9 574 snprintf(title,256,"rms_dca_z (cm) vs %s",fDCAHisto->GetAxis(3)->GetTitle());
7cc34f08 575 h1D->SetTitle(title);
576 h1D->SetBit(TH1::kLogX);
577 aFolderObj->Add(h1D);
578
4bb5f8a9 579 //
580 // fit mean and sigma
581 //
582
583 arr[3] = new TObjArray();
584 f1[3] = new TF1("gaus","gaus");
585 h2D->FitSlicesY(f1[3],0,-1,0,"QNR",arr[3]);
586
587 h1D = (TH1F*)arr[3]->At(1);
588 h1D->SetName("fit_mean_dca_z_vs_pt");
589 h1D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(3)->GetTitle());
590 h1D->GetYaxis()->SetTitle("fit_mean_dca_z (cm)");
18e588e9 591 snprintf(title,256,"fit_mean_dca_z (cm) vs %s",fDCAHisto->GetAxis(3)->GetTitle());
4bb5f8a9 592 h1D->SetTitle(title);
593 h1D->SetBit(TH1::kLogX);
594 aFolderObj->Add(h1D);
595
596 h1D = (TH1F*)arr[3]->At(2);
597 h1D->SetName("res_dca_z_vs_pt");
598 h1D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(3)->GetTitle());
599 h1D->GetYaxis()->SetTitle("res_dca_z (cm)");
18e588e9 600 snprintf(title,256,"res_dca_z (cm) vs %s",fDCAHisto->GetAxis(3)->GetTitle());
4bb5f8a9 601 h1D->SetTitle(title);
602 h1D->SetBit(TH1::kLogX);
603 aFolderObj->Add(h1D);
604
886bf4d3 605 // A - side
606 fDCAHisto->GetAxis(2)->SetRangeUser(-1.5,0.0);
7cc34f08 607
4bb5f8a9 608 h2D = (TH2F*)fDCAHisto->Projection(1,4);
609 h2D->SetName("dca_z_vs_phi_Aside");
886bf4d3 610 h2D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(4)->GetTitle());
611 h2D->GetYaxis()->SetTitle(fDCAHisto->GetAxis(1)->GetTitle());
18e588e9 612 snprintf(title,256,"%s vs %s (A-side)",fDCAHisto->GetAxis(1)->GetTitle(),fDCAHisto->GetAxis(4)->GetTitle());
7cc34f08 613 h2D->SetTitle(title);
614 aFolderObj->Add(h2D);
615
886bf4d3 616 h1D = MakeStat1D(h2D,0,0);
4bb5f8a9 617 h1D->SetName("mean_dca_z_vs_phi_Aside");
886bf4d3 618 h1D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(4)->GetTitle());
619 h1D->GetYaxis()->SetTitle("mean_dca_z (cm)");
18e588e9 620 snprintf(title,256,"mean_dca_z (cm) vs %s (A-side)",fDCAHisto->GetAxis(4)->GetTitle());
886bf4d3 621 h1D->SetTitle(title);
622 aFolderObj->Add(h1D);
7cc34f08 623
886bf4d3 624 h1D = MakeStat1D(h2D,0,1);
4bb5f8a9 625 h1D->SetName("rms_dca_z_vs_phi_Aside");
886bf4d3 626 h1D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(4)->GetTitle());
627 h1D->GetYaxis()->SetTitle("rms_dca_z (cm)");
18e588e9 628 snprintf(title,256,"rms_dca_z (cm) vs %s (A-side)",fDCAHisto->GetAxis(4)->GetTitle());
886bf4d3 629 h1D->SetTitle(title);
630 aFolderObj->Add(h1D);
631
4bb5f8a9 632 //
633 // fit mean and sigma
634 //
635 arr[4] = new TObjArray();
636 f1[4] = new TF1("gaus","gaus");
637 h2D->FitSlicesY(f1[4],0,-1,0,"QNR",arr[4]);
638
639 h1D = (TH1F*)arr[4]->At(1);
640 h1D->SetName("fit_mean_dca_z_vs_phi_Aside");
641 h1D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(4)->GetTitle());
642 h1D->GetYaxis()->SetTitle("fit_mean_dca_z (cm)");
18e588e9 643 snprintf(title,256,"fit_mean_dca_z (cm) vs %s (A-side)",fDCAHisto->GetAxis(4)->GetTitle());
4bb5f8a9 644 h1D->SetTitle(title);
645 aFolderObj->Add(h1D);
646
647 h1D = (TH1F*)arr[4]->At(2);
648 h1D->SetName("res_dca_z_vs_phi_Aside");
649 h1D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(4)->GetTitle());
650 h1D->GetYaxis()->SetTitle("res_dca_z (cm)");
18e588e9 651 snprintf(title,256,"res_dca_z (cm) vs %s (A-side)",fDCAHisto->GetAxis(4)->GetTitle());
4bb5f8a9 652 h1D->SetTitle(title);
653 aFolderObj->Add(h1D);
654
655
886bf4d3 656 // C - side
657 fDCAHisto->GetAxis(2)->SetRangeUser(0.0,1.5);
7cc34f08 658
4bb5f8a9 659 h2D = (TH2F*)fDCAHisto->Projection(1,4);
660 h2D->SetName("dca_z_vs_phi_Cside");
886bf4d3 661 h2D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(4)->GetTitle());
662 h2D->GetYaxis()->SetTitle(fDCAHisto->GetAxis(1)->GetTitle());
18e588e9 663 snprintf(title,256,"%s vs %s (C-side)",fDCAHisto->GetAxis(1)->GetTitle(),fDCAHisto->GetAxis(4)->GetTitle());
7cc34f08 664 h2D->SetTitle(title);
665 aFolderObj->Add(h2D);
666
886bf4d3 667 h1D = MakeStat1D(h2D,0,0);
4bb5f8a9 668 h1D->SetName("mean_dca_z_vs_phi_Cside");
886bf4d3 669 h1D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(4)->GetTitle());
670 h1D->GetYaxis()->SetTitle("mean_dca_z (cm)");
18e588e9 671 snprintf(title,256,"mean_dca_z (cm) vs %s (C-side)",fDCAHisto->GetAxis(4)->GetTitle());
886bf4d3 672 h1D->SetTitle(title);
673 aFolderObj->Add(h1D);
674
675 h1D = MakeStat1D(h2D,0,1);
4bb5f8a9 676 h1D->SetName("rms_dca_z_vs_phi_Cside");
886bf4d3 677 h1D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(4)->GetTitle());
678 h1D->GetYaxis()->SetTitle("rms_dca_z (cm)");
18e588e9 679 snprintf(title,256,"rms_dca_z (cm) vs %s (C-side)",fDCAHisto->GetAxis(4)->GetTitle());
886bf4d3 680 h1D->SetTitle(title);
681 aFolderObj->Add(h1D);
886bf4d3 682
4bb5f8a9 683 //
684 // fit mean and sigma
685 //
686 arr[5] = new TObjArray();
687 f1[5] = new TF1("gaus","gaus");
688 h2D->FitSlicesY(f1[5],0,-1,0,"QNR",arr[5]);
689
690 h1D = (TH1F*)arr[5]->At(1);
691 h1D->SetName("fit_mean_dca_z_vs_phi_Cside");
692 h1D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(4)->GetTitle());
693 h1D->GetYaxis()->SetTitle("fit_mean_dca_z (cm)");
18e588e9 694 snprintf(title,256,"fit_mean_dca_z (cm) vs %s (C-side)",fDCAHisto->GetAxis(4)->GetTitle());
4bb5f8a9 695 h1D->SetTitle(title);
696 aFolderObj->Add(h1D);
7cc34f08 697
4bb5f8a9 698 h1D = (TH1F*)arr[5]->At(2);
699 h1D->SetName("res_dca_z_vs_phi_Cside");
700 h1D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(4)->GetTitle());
701 h1D->GetYaxis()->SetTitle("res_dca_z (cm)");
18e588e9 702 snprintf(title,256,"res_dca_z (cm) vs %s (C-side)",fDCAHisto->GetAxis(4)->GetTitle());
4bb5f8a9 703 h1D->SetTitle(title);
704 aFolderObj->Add(h1D);
705
7cc34f08 706 // export objects to analysis folder
707 fAnalysisFolder = ExportToFolder(aFolderObj);
708
709 // delete only TObjArray
710 if(aFolderObj) delete aFolderObj;
711}
712
713//_____________________________________________________________________________
714TH1F* AliPerformanceDCA::MakeStat1D(TH2 *hist, Int_t delta0, Int_t type)
715{
716 // Return TH1F histogram
717 // delta - number of bins to integrate
718 // with mean (type == 0) or RMS (type==1)
719
720 char hname[256];
721 const char* suffix = "_stat1d";
18e588e9 722 snprintf(hname,256,"%s%s",hist->GetName(),suffix);
7cc34f08 723 TAxis* xaxis = hist->GetXaxis();
724 Int_t nbinx = xaxis->GetNbins();
725
726 TH1F *hnew = (TH1F*)hist->ProjectionX()->Clone();
727 hnew->SetName(hname);
728
729 char name[256];
730 for (Int_t ix=0; ix<=nbinx;ix++) {
18e588e9 731 snprintf(name,256,"%s_%d",hist->GetName(),ix);
7cc34f08 732 TH1 *projection = hist->ProjectionY(name,ix-delta0,ix+delta0);
733
734 Float_t stat= 0., stat_err =0.;
735 if (type==0) { stat = projection->GetMean(); stat_err = projection->GetMeanError(); }
736 if (type==1) { stat = projection->GetRMS(); stat_err = projection->GetRMSError(); }
737
738 hnew->SetBinContent(ix, stat);
739 hnew->SetBinError(ix, stat_err);
740 }
741
742return hnew;
743}
744
745//_____________________________________________________________________________
746TH2F* AliPerformanceDCA::MakeStat2D(TH3 *hist, Int_t delta0, Int_t delta1, Int_t type)
747{
748 // Return TH1F histogram
749 // delta0 - number of bins to integrate in x
750 // delta1 - number of bins to integrate in y
751 // with mean (type==0) or RMS (type==1)
752
753 char hname[256];
754 const char* suffix = "_stat2d";
18e588e9 755 snprintf(hname,256,"%s%s",hist->GetName(),suffix);
7cc34f08 756
757 TAxis* xaxis = hist->GetXaxis();
758 Int_t nbinx = xaxis->GetNbins();
759
760 TH2F *hnew = (TH2F*)hist->Project3D("yx")->Clone();
761 hnew->SetName(hname);
762
763 TAxis* yaxis = hist->GetYaxis();
764 Int_t nbiny = yaxis->GetNbins();
765
766 char name[256];
767 for (Int_t ix=0; ix<=nbinx;ix++) {
768 for (Int_t iy=0; iy<=nbiny;iy++) {
18e588e9 769 snprintf(name,256,"%s_%d_%d",hist->GetName(),ix,iy);
7cc34f08 770 TH1 *projection = hist->ProjectionZ(name,ix-delta0,ix+delta0,iy-delta1,iy+delta1);
771
772 Float_t stat= 0., stat_err =0.;
773 if (type==0) { stat = projection->GetMean(); stat_err = projection->GetMeanError(); }
774 if (type==1) { stat = projection->GetRMS(); stat_err = projection->GetRMSError(); }
775
776 hnew->SetBinContent(ix,iy,stat);
777 hnew->SetBinError(ix,iy,stat_err);
778 }
779 }
780
781return hnew;
782}
783
784//_____________________________________________________________________________
785TFolder* AliPerformanceDCA::ExportToFolder(TObjArray * array)
786{
787 // recreate folder avery time and export objects to new one
788 //
789 AliPerformanceDCA * comp=this;
790 TFolder *folder = comp->GetAnalysisFolder();
791
792 TString name, title;
793 TFolder *newFolder = 0;
794 Int_t i = 0;
795 Int_t size = array->GetSize();
796
797 if(folder) {
798 // get name and title from old folder
799 name = folder->GetName();
800 title = folder->GetTitle();
801
802 // delete old one
803 delete folder;
804
805 // create new one
806 newFolder = CreateFolder(name.Data(),title.Data());
807 newFolder->SetOwner();
808
809 // add objects to folder
810 while(i < size) {
811 newFolder->Add(array->At(i));
812 i++;
813 }
814 }
815
816return newFolder;
817}
818
819//_____________________________________________________________________________
820TFolder* AliPerformanceDCA::CreateFolder(TString name,TString title) {
821// create folder for analysed histograms
822TFolder *folder = 0;
823 folder = new TFolder(name.Data(),title.Data());
824
825 return folder;
826}