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.
10 // Author: J.Otwinowski 04/02/2008
11 //------------------------------------------------------------------------------
15 // after running comparison task, read the file, and get component
16 gROOT->LoadMacro("$ALICE_ROOT/PWGPP/Macros/LoadMyLibs.C");
18 TFile f("Output.root");
19 AliPerformanceDCA * compObj = (AliPerformanceDCA*)coutput->FindObject("AliPerformanceDCA");
21 // Analyse comparison data
24 // the output histograms/graphs will be stored in the folder "folderDCA"
25 compObj->GetAnalysisFolder()->ls("*");
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();
44 #include "AliPerformanceDCA.h"
45 #include "AliESDEvent.h"
46 #include "AliESDVertex.h"
48 #include "AliMathBase.h"
49 #include "AliRecInfoCuts.h"
50 #include "AliMCInfoCuts.h"
52 #include "AliMCEvent.h"
53 #include "AliTracker.h"
54 #include "AliHeader.h"
55 #include "AliGenEventHeader.h"
59 ClassImp(AliPerformanceDCA)
61 //_____________________________________________________________________________
62 AliPerformanceDCA::AliPerformanceDCA():
63 AliPerformanceObject("AliPerformanceDCA"),
75 // default constructor
79 //_____________________________________________________________________________
80 AliPerformanceDCA::AliPerformanceDCA(const Char_t* name="AliPerformanceDCA", const Char_t* title="AliPerformanceDCA",Int_t analysisMode=0, Bool_t hptGenerator=kFALSE):
81 AliPerformanceObject(name,title),
95 SetAnalysisMode(analysisMode);
96 SetHptGenerator(hptGenerator);
100 //_____________________________________________________________________________
101 AliPerformanceDCA::~AliPerformanceDCA()
104 if(fDCAHisto) delete fDCAHisto; fDCAHisto=0;
105 if(fAnalysisFolder) delete fAnalysisFolder; fAnalysisFolder=0;
108 //_____________________________________________________________________________
109 void AliPerformanceDCA::Init()
115 Double_t ptMin = 1.e-2, ptMax = 20.;
117 Double_t *binsPt = 0;
119 if (IsHptGenerator()) {
122 binsPt = CreateLogAxis(nPtBins,ptMin,ptMax);
126 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.};
127 Double_t ptMin = 0., ptMax = 10.;
129 if(IsHptGenerator() == kTRUE) {
131 ptMin = 0.; ptMax = 100.;
135 //dca_r, dca_z, eta, pt
136 Int_t binsQA[5] = {100,100,30,nPtBins,144};
137 Double_t xminQA[5] = {-10.,-10.,-1.5,ptMin,0.};
138 Double_t xmaxQA[5] = {10.,10.,1.5,ptMax,2*TMath::Pi()};
140 fDCAHisto = new THnSparseF("fDCAHisto","dca_r:dca_z:eta:pt:phi",5,binsQA,xminQA,xmaxQA);
141 fDCAHisto->SetBinEdges(3,binsPt);
143 fDCAHisto->GetAxis(0)->SetTitle("dca_r (cm)");
144 fDCAHisto->GetAxis(1)->SetTitle("dca_z (cm)");
145 fDCAHisto->GetAxis(2)->SetTitle("#eta");
146 fDCAHisto->GetAxis(3)->SetTitle("p_{T} (GeV/c)");
147 fDCAHisto->GetAxis(4)->SetTitle("phi (rad)");
152 AliDebug(AliLog::kError, "ERROR: Cannot find AliMCInfoCuts object");
154 AliDebug(AliLog::kError, "ERROR: Cannot find AliRecInfoCuts object");
157 fAnalysisFolder = CreateFolder("folderDCA","Analysis DCA Folder");
160 //_____________________________________________________________________________
161 void AliPerformanceDCA::ProcessTPC(AliStack* const stack, AliESDtrack *const esdTrack, AliESDEvent* const esdEvent)
163 // Fill DCA comparison information
164 if(!esdEvent) return;
165 if(!esdTrack) return;
167 if( IsUseTrackVertex() )
169 // Relate TPC inner params to prim. vertex
170 const AliESDVertex *vtxESD = esdEvent->GetPrimaryVertexTracks();
171 Double_t x[3]; esdTrack->GetXYZ(x);
172 Double_t b[3]; AliTracker::GetBxByBz(x,b);
173 Bool_t isOK = esdTrack->RelateToVertexTPCBxByBz(vtxESD, b, kVeryBig);
177 // JMT -- recaluclate DCA for HLT if not present
178 if ( dca[0] == 0. && dca[1] == 0. ) {
179 track->GetDZ( vtxESD->GetX(), vtxESD->GetY(), vtxESD->GetZ(), esdEvent->GetMagneticField(), dca );
184 // get TPC inner params at DCA to prim. vertex
185 const AliExternalTrackParam *track = esdTrack->GetTPCInnerParam();
188 // read from ESD track
189 Float_t dca[2], cov[3]; // dca_xy, dca_z, sigma_xy, sigma_xy_z, sigma_z
190 esdTrack->GetImpactParametersTPC(dca,cov);
192 if (esdTrack->GetTPCNcls()<fCutsRC->GetMinNClustersTPC()) return; // min. nb. TPC clusters
194 Double_t vDCAHisto[5]={dca[0],dca[1],track->Eta(),track->Pt(),track->Phi()};
195 fDCAHisto->Fill(vDCAHisto);
198 // Fill rec vs MC information
204 //_____________________________________________________________________________
205 void AliPerformanceDCA::ProcessTPCITS(AliStack* const stack, AliESDtrack *const esdTrack, AliESDEvent* const esdEvent)
207 // Fill DCA comparison information
208 if(!esdTrack) return;
209 if(!esdEvent) return;
211 if( IsUseTrackVertex() )
213 // Relate TPC inner params to prim. vertex
214 const AliESDVertex *vtxESD = esdEvent->GetPrimaryVertexTracks();
215 Double_t x[3]; esdTrack->GetXYZ(x);
216 Double_t b[3]; AliTracker::GetBxByBz(x,b);
217 Bool_t isOK = esdTrack->RelateToVertexBxByBz(vtxESD, b, kVeryBig);
221 // JMT -- recaluclate DCA for HLT if not present
222 if ( dca[0] == 0. && dca[1] == 0. ) {
223 track->GetDZ( vtxESD->GetX(), vtxESD->GetY(), vtxESD->GetZ(), esdEvent->GetMagneticField(), dca );
228 Float_t dca[2], cov[3]; // dca_xy, dca_z, sigma_xy, sigma_xy_z, sigma_z
229 esdTrack->GetImpactParameters(dca,cov);
231 if ((esdTrack->GetStatus()&AliESDtrack::kTPCrefit)==0) return; // TPC refit
232 if (esdTrack->GetTPCNcls()<fCutsRC->GetMinNClustersTPC()) return; // min. nb. TPC clusters
233 if(esdTrack->GetITSclusters(0)<fCutsRC->GetMinNClustersITS()) return; // min. nb. ITS clusters
235 Double_t vDCAHisto[5]={dca[0],dca[1],esdTrack->Eta(),esdTrack->Pt(), esdTrack->Phi()};
236 fDCAHisto->Fill(vDCAHisto);
239 // Fill rec vs MC information
245 void AliPerformanceDCA::ProcessConstrained(AliStack* const /*stack*/, AliESDtrack *const /*esdTrack*/)
247 // Fill DCA comparison information
249 AliDebug(AliLog::kWarning, "Warning: Not implemented");
252 //_____________________________________________________________________________
253 Long64_t AliPerformanceDCA::Merge(TCollection* const list)
255 // Merge list of objects (needed by PROOF)
263 TIterator* iter = list->MakeIterator();
266 // collection of generated histograms
268 while((obj = iter->Next()) != 0)
270 AliPerformanceDCA* entry = dynamic_cast<AliPerformanceDCA*>(obj);
271 if (entry == 0) continue;
273 fDCAHisto->Add(entry->fDCAHisto);
280 //_____________________________________________________________________________
281 void AliPerformanceDCA::Exec(AliMCEvent* const mcEvent, AliESDEvent *const esdEvent, AliESDfriend *const esdFriend, const Bool_t bUseMC, const Bool_t bUseESDfriend)
283 // Process comparison information
287 Error("Exec","esdEvent not available");
290 AliHeader* header = 0;
291 AliGenEventHeader* genHeader = 0;
298 Error("Exec","mcEvent not available");
301 // get MC event header
302 header = mcEvent->Header();
304 Error("Exec","Header not available");
308 stack = mcEvent->Stack();
310 Error("Exec","Stack not available");
314 genHeader = header->GenEventHeader();
316 Error("Exec","Could not retrieve genHeader from Header");
319 genHeader->PrimaryVertex(vtxMC);
325 Error("Exec","esdFriend not available");
331 if(!bUseMC &&GetTriggerClass()) {
332 Bool_t isEventTriggered = esdEvent->IsTriggerClassFired(GetTriggerClass());
333 if(!isEventTriggered) return;
337 const AliESDVertex *vtxESD = NULL;
338 if( IsUseTrackVertex() )
341 vtxESD = esdEvent->GetPrimaryVertexTracks();
345 vtxESD = esdEvent->GetPrimaryVertexTPC();
347 if(vtxESD && (vtxESD->GetStatus()<=0)) return;
350 for (Int_t iTrack = 0; iTrack < esdEvent->GetNumberOfTracks(); iTrack++)
352 AliESDtrack *track = esdEvent->GetTrack(iTrack);
355 if(GetAnalysisMode() == 0) ProcessTPC(stack,track,esdEvent);
356 else if(GetAnalysisMode() == 1) ProcessTPCITS(stack,track,esdEvent);
357 else if(GetAnalysisMode() == 2) ProcessConstrained(stack,track);
359 printf("ERROR: AnalysisMode %d \n",fAnalysisMode);
365 //_____________________________________________________________________________
366 void AliPerformanceDCA::Analyse()
369 // Analyse comparison information and store output histograms
370 // in the analysis folder "folderDCA"
373 TH1::AddDirectory(kFALSE);
376 TObjArray *aFolderObj = new TObjArray;
377 if(!aFolderObj) return;
379 TObjArray *arr[6] = {0};
384 // set pt measurable range
385 //fDCAHisto->GetAxis(3)->SetRangeUser(0.10,10.);
388 h2D = (TH2F*)fDCAHisto->Projection(0,1); // inverse projection convention
389 h2D->SetName("dca_r_vs_dca_z");
390 h2D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(1)->GetTitle());
391 h2D->GetYaxis()->SetTitle(fDCAHisto->GetAxis(0)->GetTitle());
392 snprintf(title,256,"%s vs %s",fDCAHisto->GetAxis(0)->GetTitle(),fDCAHisto->GetAxis(1)->GetTitle());
393 h2D->SetTitle(title);
394 aFolderObj->Add(h2D);
397 h2D = (TH2F*)fDCAHisto->Projection(0,2);
398 h2D->SetName("dca_r_vs_eta");
399 h2D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(2)->GetTitle());
400 h2D->GetYaxis()->SetTitle(fDCAHisto->GetAxis(0)->GetTitle());
401 snprintf(title,256,"%s vs %s",fDCAHisto->GetAxis(2)->GetTitle(),fDCAHisto->GetAxis(0)->GetTitle());
402 h2D->SetTitle(title);
403 aFolderObj->Add(h2D);
408 h1D = MakeStat1D(h2D,0,0);
409 h1D->SetName("mean_dca_r_vs_eta");
410 h1D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(2)->GetTitle());
411 h1D->GetYaxis()->SetTitle("mean_dca_r (cm)");
412 snprintf(title,256," mean_dca_r (cm) vs %s",fDCAHisto->GetAxis(2)->GetTitle());
413 h1D->SetTitle(title);
414 aFolderObj->Add(h1D);
416 h1D = MakeStat1D(h2D,0,1);
417 h1D->SetName("rms_dca_r_vs_eta");
418 h1D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(2)->GetTitle());
419 h1D->GetYaxis()->SetTitle("rms_dca_r (cm)");
420 snprintf(title,256," rms_dca_r (cm) vs %s",fDCAHisto->GetAxis(2)->GetTitle());
421 h1D->SetTitle(title);
422 aFolderObj->Add(h1D);
425 // fit mean and sigma
428 arr[0] = new TObjArray();
429 f1[0] = new TF1("gaus","gaus");
430 h2D->FitSlicesY(f1[0],0,-1,0,"QNR",arr[0]);
432 h1D = (TH1F*)arr[0]->At(1);
433 h1D->SetName("fit_mean_dca_r_vs_eta");
434 h1D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(2)->GetTitle());
435 h1D->GetYaxis()->SetTitle("fit_mean_dca_r (cm)");
436 snprintf(title,256," fit_mean_dca_r (cm) vs %s",fDCAHisto->GetAxis(2)->GetTitle());
437 h1D->SetTitle(title);
438 aFolderObj->Add(h1D);
440 h1D = (TH1F*)arr[0]->At(2);
441 h1D->SetName("res_dca_r_vs_eta");
442 h1D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(2)->GetTitle());
443 h1D->GetYaxis()->SetTitle("res_dca_r (cm)");
444 snprintf(title,256," res_dca_r (cm) vs %s",fDCAHisto->GetAxis(2)->GetTitle());
445 h1D->SetTitle(title);
446 aFolderObj->Add(h1D);
451 h2D = (TH2F*)fDCAHisto->Projection(0,3);
452 h2D->SetName("dca_r_vs_pt");
453 h2D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(3)->GetTitle());
454 h2D->GetYaxis()->SetTitle(fDCAHisto->GetAxis(0)->GetTitle());
455 snprintf(title,256,"%s vs %s",fDCAHisto->GetAxis(0)->GetTitle(),fDCAHisto->GetAxis(3)->GetTitle());
456 h2D->SetTitle(title);
457 h2D->SetBit(TH1::kLogX);
458 aFolderObj->Add(h2D);
460 h1D = MakeStat1D(h2D,0,0);
461 h1D->SetName("mean_dca_r_vs_pt");
462 h1D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(3)->GetTitle());
463 h1D->GetYaxis()->SetTitle("mean_dca_r (cm)");
464 snprintf(title,256,"mean_dca_r (cm) vs %s",fDCAHisto->GetAxis(3)->GetTitle());
465 h1D->SetTitle(title);
466 h1D->SetBit(TH1::kLogX);
467 aFolderObj->Add(h1D);
469 h1D = MakeStat1D(h2D,0,1);
470 h1D->SetName("rms_dca_r_vs_pt");
471 h1D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(3)->GetTitle());
472 h1D->GetYaxis()->SetTitle("rms_dca_r (cm)");
473 snprintf(title,256,"rms_dca_r (cm) vs %s",fDCAHisto->GetAxis(3)->GetTitle());
474 h1D->SetTitle(title);
475 h1D->SetBit(TH1::kLogX);
476 aFolderObj->Add(h1D);
479 // fit mean and sigma
482 arr[1] = new TObjArray();
483 f1[1] = new TF1("gaus","gaus");
484 h2D->FitSlicesY(f1[1],0,-1,0,"QNR",arr[1]);
486 h1D = (TH1F*)arr[1]->At(1);
487 h1D->SetName("fit_mean_dca_r_vs_pt");
488 h1D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(3)->GetTitle());
489 h1D->GetYaxis()->SetTitle("fit_mean_dca_r (cm)");
490 snprintf(title,256,"fit_mean_dca_r (cm) vs %s",fDCAHisto->GetAxis(3)->GetTitle());
491 h1D->SetTitle(title);
492 h1D->SetBit(TH1::kLogX);
493 aFolderObj->Add(h1D);
495 h1D = (TH1F*)arr[1]->At(2);
496 h1D->SetName("res_dca_r_vs_pt");
497 h1D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(3)->GetTitle());
498 h1D->GetYaxis()->SetTitle("res_dca_r (cm)");
499 snprintf(title,256,"res_dca_r (cm) vs %s",fDCAHisto->GetAxis(3)->GetTitle());
500 h1D->SetTitle(title);
501 h1D->SetBit(TH1::kLogX);
502 aFolderObj->Add(h1D);
505 h2D = (TH2F*)fDCAHisto->Projection(1,2);
506 h2D->SetName("dca_z_vs_eta");
507 h2D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(2)->GetTitle());
508 h2D->GetYaxis()->SetTitle(fDCAHisto->GetAxis(1)->GetTitle());
509 snprintf(title,256,"%s vs %s",fDCAHisto->GetAxis(1)->GetTitle(),fDCAHisto->GetAxis(2)->GetTitle());
510 h2D->SetTitle(title);
511 aFolderObj->Add(h2D);
513 h1D = MakeStat1D(h2D,0,0);
514 h1D->SetName("mean_dca_z_vs_eta");
515 h1D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(2)->GetTitle());
516 h1D->GetYaxis()->SetTitle("mean_dca_z (cm)");
517 snprintf(title,256,"mean_dca_z (cm) vs %s",fDCAHisto->GetAxis(2)->GetTitle());
518 h1D->SetTitle(title);
519 aFolderObj->Add(h1D);
521 h1D = MakeStat1D(h2D,0,1);
522 h1D->SetName("rms_dca_z_vs_eta");
523 h1D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(2)->GetTitle());
524 h1D->GetYaxis()->SetTitle("rms_dca_z (cm)");
525 snprintf(title,256,"rms_dca_z (cm) vs %s",fDCAHisto->GetAxis(2)->GetTitle());
526 h1D->SetTitle(title);
527 aFolderObj->Add(h1D);
530 // fit mean and sigma
532 arr[2] = new TObjArray();
533 f1[2] = new TF1("gaus","gaus");
534 h2D->FitSlicesY(f1[2],0,-1,0,"QNR",arr[2]);
536 h1D = (TH1F*)arr[2]->At(1);
537 h1D->SetName("fit_mean_dca_z_vs_eta");
538 h1D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(2)->GetTitle());
539 h1D->GetYaxis()->SetTitle("fit_mean_dca_z (cm)");
540 snprintf(title,256,"fit_mean_dca_z (cm) vs %s",fDCAHisto->GetAxis(2)->GetTitle());
541 h1D->SetTitle(title);
542 aFolderObj->Add(h1D);
544 h1D = (TH1F*)arr[2]->At(2);
545 h1D->SetName("res_dca_z_vs_eta");
546 h1D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(2)->GetTitle());
547 h1D->GetYaxis()->SetTitle("res_dca_z (cm)");
548 snprintf(title,256,"res_dca_z (cm) vs %s",fDCAHisto->GetAxis(2)->GetTitle());
549 h1D->SetTitle(title);
550 aFolderObj->Add(h1D);
553 h2D = (TH2F*)fDCAHisto->Projection(1,3);
554 h2D->SetName("dca_z_vs_pt");
555 h2D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(3)->GetTitle());
556 h2D->GetYaxis()->SetTitle(fDCAHisto->GetAxis(1)->GetTitle());
557 snprintf(title,256,"%s vs %s",fDCAHisto->GetAxis(1)->GetTitle(),fDCAHisto->GetAxis(3)->GetTitle());
558 h2D->SetTitle(title);
559 h2D->SetBit(TH1::kLogX);
560 aFolderObj->Add(h2D);
562 h1D = MakeStat1D(h2D,0,0);
563 h1D->SetName("mean_dca_z_vs_pt");
564 h1D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(3)->GetTitle());
565 h1D->GetYaxis()->SetTitle("mean_dca_z (cm)");
566 snprintf(title,256,"mean_dca_z (cm) vs %s",fDCAHisto->GetAxis(3)->GetTitle());
567 h1D->SetTitle(title);
568 h1D->SetBit(TH1::kLogX);
569 aFolderObj->Add(h1D);
571 h1D = MakeStat1D(h2D,0,1);
572 h1D->SetName("rms_dca_z_vs_pt");
573 h1D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(3)->GetTitle());
574 h1D->GetYaxis()->SetTitle("rms_dca_z (cm)");
575 snprintf(title,256,"rms_dca_z (cm) vs %s",fDCAHisto->GetAxis(3)->GetTitle());
576 h1D->SetTitle(title);
577 h1D->SetBit(TH1::kLogX);
578 aFolderObj->Add(h1D);
581 // fit mean and sigma
584 arr[3] = new TObjArray();
585 f1[3] = new TF1("gaus","gaus");
586 h2D->FitSlicesY(f1[3],0,-1,0,"QNR",arr[3]);
588 h1D = (TH1F*)arr[3]->At(1);
589 h1D->SetName("fit_mean_dca_z_vs_pt");
590 h1D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(3)->GetTitle());
591 h1D->GetYaxis()->SetTitle("fit_mean_dca_z (cm)");
592 snprintf(title,256,"fit_mean_dca_z (cm) vs %s",fDCAHisto->GetAxis(3)->GetTitle());
593 h1D->SetTitle(title);
594 h1D->SetBit(TH1::kLogX);
595 aFolderObj->Add(h1D);
597 h1D = (TH1F*)arr[3]->At(2);
598 h1D->SetName("res_dca_z_vs_pt");
599 h1D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(3)->GetTitle());
600 h1D->GetYaxis()->SetTitle("res_dca_z (cm)");
601 snprintf(title,256,"res_dca_z (cm) vs %s",fDCAHisto->GetAxis(3)->GetTitle());
602 h1D->SetTitle(title);
603 h1D->SetBit(TH1::kLogX);
604 aFolderObj->Add(h1D);
607 fDCAHisto->GetAxis(2)->SetRangeUser(-1.5,0.0);
609 h2D = (TH2F*)fDCAHisto->Projection(1,4);
610 h2D->SetName("dca_z_vs_phi_Aside");
611 h2D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(4)->GetTitle());
612 h2D->GetYaxis()->SetTitle(fDCAHisto->GetAxis(1)->GetTitle());
613 snprintf(title,256,"%s vs %s (A-side)",fDCAHisto->GetAxis(1)->GetTitle(),fDCAHisto->GetAxis(4)->GetTitle());
614 h2D->SetTitle(title);
615 aFolderObj->Add(h2D);
617 h1D = MakeStat1D(h2D,0,0);
618 h1D->SetName("mean_dca_z_vs_phi_Aside");
619 h1D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(4)->GetTitle());
620 h1D->GetYaxis()->SetTitle("mean_dca_z (cm)");
621 snprintf(title,256,"mean_dca_z (cm) vs %s (A-side)",fDCAHisto->GetAxis(4)->GetTitle());
622 h1D->SetTitle(title);
623 aFolderObj->Add(h1D);
625 h1D = MakeStat1D(h2D,0,1);
626 h1D->SetName("rms_dca_z_vs_phi_Aside");
627 h1D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(4)->GetTitle());
628 h1D->GetYaxis()->SetTitle("rms_dca_z (cm)");
629 snprintf(title,256,"rms_dca_z (cm) vs %s (A-side)",fDCAHisto->GetAxis(4)->GetTitle());
630 h1D->SetTitle(title);
631 aFolderObj->Add(h1D);
634 // fit mean and sigma
636 arr[4] = new TObjArray();
637 f1[4] = new TF1("gaus","gaus");
638 h2D->FitSlicesY(f1[4],0,-1,0,"QNR",arr[4]);
640 h1D = (TH1F*)arr[4]->At(1);
641 h1D->SetName("fit_mean_dca_z_vs_phi_Aside");
642 h1D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(4)->GetTitle());
643 h1D->GetYaxis()->SetTitle("fit_mean_dca_z (cm)");
644 snprintf(title,256,"fit_mean_dca_z (cm) vs %s (A-side)",fDCAHisto->GetAxis(4)->GetTitle());
645 h1D->SetTitle(title);
646 aFolderObj->Add(h1D);
648 h1D = (TH1F*)arr[4]->At(2);
649 h1D->SetName("res_dca_z_vs_phi_Aside");
650 h1D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(4)->GetTitle());
651 h1D->GetYaxis()->SetTitle("res_dca_z (cm)");
652 snprintf(title,256,"res_dca_z (cm) vs %s (A-side)",fDCAHisto->GetAxis(4)->GetTitle());
653 h1D->SetTitle(title);
654 aFolderObj->Add(h1D);
658 fDCAHisto->GetAxis(2)->SetRangeUser(0.0,1.5);
660 h2D = (TH2F*)fDCAHisto->Projection(1,4);
661 h2D->SetName("dca_z_vs_phi_Cside");
662 h2D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(4)->GetTitle());
663 h2D->GetYaxis()->SetTitle(fDCAHisto->GetAxis(1)->GetTitle());
664 snprintf(title,256,"%s vs %s (C-side)",fDCAHisto->GetAxis(1)->GetTitle(),fDCAHisto->GetAxis(4)->GetTitle());
665 h2D->SetTitle(title);
666 aFolderObj->Add(h2D);
668 h1D = MakeStat1D(h2D,0,0);
669 h1D->SetName("mean_dca_z_vs_phi_Cside");
670 h1D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(4)->GetTitle());
671 h1D->GetYaxis()->SetTitle("mean_dca_z (cm)");
672 snprintf(title,256,"mean_dca_z (cm) vs %s (C-side)",fDCAHisto->GetAxis(4)->GetTitle());
673 h1D->SetTitle(title);
674 aFolderObj->Add(h1D);
676 h1D = MakeStat1D(h2D,0,1);
677 h1D->SetName("rms_dca_z_vs_phi_Cside");
678 h1D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(4)->GetTitle());
679 h1D->GetYaxis()->SetTitle("rms_dca_z (cm)");
680 snprintf(title,256,"rms_dca_z (cm) vs %s (C-side)",fDCAHisto->GetAxis(4)->GetTitle());
681 h1D->SetTitle(title);
682 aFolderObj->Add(h1D);
685 // fit mean and sigma
687 arr[5] = new TObjArray();
688 f1[5] = new TF1("gaus","gaus");
689 h2D->FitSlicesY(f1[5],0,-1,0,"QNR",arr[5]);
691 h1D = (TH1F*)arr[5]->At(1);
692 h1D->SetName("fit_mean_dca_z_vs_phi_Cside");
693 h1D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(4)->GetTitle());
694 h1D->GetYaxis()->SetTitle("fit_mean_dca_z (cm)");
695 snprintf(title,256,"fit_mean_dca_z (cm) vs %s (C-side)",fDCAHisto->GetAxis(4)->GetTitle());
696 h1D->SetTitle(title);
697 aFolderObj->Add(h1D);
699 h1D = (TH1F*)arr[5]->At(2);
700 h1D->SetName("res_dca_z_vs_phi_Cside");
701 h1D->GetXaxis()->SetTitle(fDCAHisto->GetAxis(4)->GetTitle());
702 h1D->GetYaxis()->SetTitle("res_dca_z (cm)");
703 snprintf(title,256,"res_dca_z (cm) vs %s (C-side)",fDCAHisto->GetAxis(4)->GetTitle());
704 h1D->SetTitle(title);
705 aFolderObj->Add(h1D);
707 // export objects to analysis folder
708 fAnalysisFolder = ExportToFolder(aFolderObj);
710 // delete only TObjArray
711 if(aFolderObj) delete aFolderObj;
714 //_____________________________________________________________________________
715 TH1F* AliPerformanceDCA::MakeStat1D(TH2 *hist, Int_t delta0, Int_t type)
717 // Return TH1F histogram
718 // delta - number of bins to integrate
719 // with mean (type == 0) or RMS (type==1)
722 const char* suffix = "_stat1d";
723 snprintf(hname,256,"%s%s",hist->GetName(),suffix);
724 TAxis* xaxis = hist->GetXaxis();
725 Int_t nbinx = xaxis->GetNbins();
727 TH1F *hnew = (TH1F*)hist->ProjectionX()->Clone();
728 hnew->SetName(hname);
731 for (Int_t ix=0; ix<=nbinx;ix++) {
732 snprintf(name,256,"%s_%d",hist->GetName(),ix);
733 TH1 *projection = hist->ProjectionY(name,ix-delta0,ix+delta0);
735 Float_t stat= 0., stat_err =0.;
736 if (type==0) { stat = projection->GetMean(); stat_err = projection->GetMeanError(); }
737 if (type==1) { stat = projection->GetRMS(); stat_err = projection->GetRMSError(); }
739 hnew->SetBinContent(ix, stat);
740 hnew->SetBinError(ix, stat_err);
746 //_____________________________________________________________________________
747 TH2F* AliPerformanceDCA::MakeStat2D(TH3 *hist, Int_t delta0, Int_t delta1, Int_t type)
749 // Return TH1F histogram
750 // delta0 - number of bins to integrate in x
751 // delta1 - number of bins to integrate in y
752 // with mean (type==0) or RMS (type==1)
755 const char* suffix = "_stat2d";
756 snprintf(hname,256,"%s%s",hist->GetName(),suffix);
758 TAxis* xaxis = hist->GetXaxis();
759 Int_t nbinx = xaxis->GetNbins();
761 TH2F *hnew = (TH2F*)hist->Project3D("yx")->Clone();
762 hnew->SetName(hname);
764 TAxis* yaxis = hist->GetYaxis();
765 Int_t nbiny = yaxis->GetNbins();
768 for (Int_t ix=0; ix<=nbinx;ix++) {
769 for (Int_t iy=0; iy<=nbiny;iy++) {
770 snprintf(name,256,"%s_%d_%d",hist->GetName(),ix,iy);
771 TH1 *projection = hist->ProjectionZ(name,ix-delta0,ix+delta0,iy-delta1,iy+delta1);
773 Float_t stat= 0., stat_err =0.;
774 if (type==0) { stat = projection->GetMean(); stat_err = projection->GetMeanError(); }
775 if (type==1) { stat = projection->GetRMS(); stat_err = projection->GetRMSError(); }
777 hnew->SetBinContent(ix,iy,stat);
778 hnew->SetBinError(ix,iy,stat_err);
785 //_____________________________________________________________________________
786 TFolder* AliPerformanceDCA::ExportToFolder(TObjArray * array)
788 // recreate folder avery time and export objects to new one
790 AliPerformanceDCA * comp=this;
791 TFolder *folder = comp->GetAnalysisFolder();
794 TFolder *newFolder = 0;
796 Int_t size = array->GetSize();
799 // get name and title from old folder
800 name = folder->GetName();
801 title = folder->GetTitle();
807 newFolder = CreateFolder(name.Data(),title.Data());
808 newFolder->SetOwner();
810 // add objects to folder
812 newFolder->Add(array->At(i));
820 //_____________________________________________________________________________
821 TFolder* AliPerformanceDCA::CreateFolder(TString name,TString title) {
822 // create folder for analysed histograms
824 folder = new TFolder(name.Data(),title.Data());