]>
Commit | Line | Data |
---|---|---|
b2a02bce | 1 | // $Id$ |
2 | ||
3 | /* | |
4 | Example macro for evaluating the tracking efficiencies of HLT tracker (runtracker.C). | |
5 | This assumes that you have have compiled with the DOMC flag on (see doc/README). | |
6 | ||
7 | path : path to where the alirunfile.root file is located. This is original rootfile | |
8 | containing the simulated particle arrays etc... | |
9 | ||
10 | trackpath : path to where the output file from the runtracker was written. | |
11 | ||
12 | offlinepath : path to where the files containing lists of the simulated particles are located. | |
13 | This are files generated by the AliTPCComparison macro (good_tracks_tpc) | |
14 | ||
4aa41877 | 15 | All efficiencies and ptresoutions are filled in histograms (see AliHLTEvaluate::CreateHistograms) |
b2a02bce | 16 | and written to file. |
17 | */ | |
18 | ||
0bd0c1ef | 19 | #ifndef __CINT__ |
4aa41877 | 20 | #include "AliHLTLogger.h" |
21 | #include "AliHLTEvaluate.h" | |
22 | #include "AliHLTFileHandler.h" | |
23 | #include "AliHLTDigitData.h" | |
24 | #include "AliHLTTransform.h" | |
0bd0c1ef | 25 | #include "AliLevel3.h" |
26 | #include <TROOT.h> | |
27 | #include <TNtuple.h> | |
28 | #include <TRandom.h> | |
29 | #include <TSystem.h> | |
30 | #include <TStyle.h> | |
31 | #include <TFile.h> | |
32 | #include <TGraphErrors.h> | |
33 | #include <TCanvas.h> | |
34 | #include <TF1.h> | |
35 | #include <stdio.h> | |
36 | #include <iostream.h> | |
37 | #include <time.h> | |
38 | #endif | |
39 | ||
40 | void evaltracker(Char_t *path="./",Char_t *trackpath="./tracker/",char *offlinepath="./",int nevent=1) | |
b2a02bce | 41 | { |
42 | //Make sure you got the correct parameters: | |
4aa41877 | 43 | AliHLTTransform::Init(path,kTRUE); |
b2a02bce | 44 | |
45 | //Define which slices to include: | |
0bd0c1ef | 46 | Int_t slicerange[2]={0,35}; |
b2a02bce | 47 | |
48 | //Define which padrows to include (should normally be all) | |
4aa41877 | 49 | //Int_t rowrange[2]={AliHLTTransform::GetFirstRow(-1),AliHLTTransform::GetLastRow(-1)}; |
b2a02bce | 50 | |
51 | //Define the minimum number of clusters on a simulated track to be found: | |
4aa41877 | 52 | Int_t good = (Int_t)(0.4*AliHLTTransform::GetNRows()); |
b2a02bce | 53 | |
54 | //Define the minumum number of clusters on a found track to be found (should normally be the same as above) | |
4aa41877 | 55 | Int_t nclusters = (Int_t)(0.4*AliHLTTransform::GetNRows()); |
b2a02bce | 56 | |
57 | //Define which pt range to include | |
0bd0c1ef | 58 | Float_t ptmin = 0.1; |
59 | Float_t ptmax = 4.; | |
b2a02bce | 60 | |
61 | //Define the maximum ratio of false clusters which are allowed on a found track (default=0.1 -> 10%) | |
0bd0c1ef | 62 | Float_t maxfalseratio = 0.1; |
b2a02bce | 63 | |
4aa41877 | 64 | AliHLTEvaluate *a = new AliHLTEvaluate(trackpath,nclusters,good,ptmin,ptmax,slicerange); |
b2a02bce | 65 | a->CreateHistos(20,0.1,4.1);//(nbins in pt,minpt,maxpt) -> the same as used by standard offline |
66 | a->SetMaxFalseClusters(maxfalseratio); | |
67 | //a->SetStandardComparison(kFALSE); //use AliTPCComparison_HLT | |
68 | ||
69 | //Loop over the number of events: | |
0bd0c1ef | 70 | for(Int_t event=0; event<nevent; event++) |
b2a02bce | 71 | { |
72 | cout<<"Processing event "<<event<<endl; | |
73 | ||
74 | //Load all reconstructed tracks for this event | |
75 | a->LoadData(event,kTRUE); //Last argument indicates that output data from tracker is written in whole slice format | |
76 | ||
77 | //Loop over cluster list for each track, retrieve MC id, and verify the found track | |
78 | a->AssignIDs(); | |
79 | ||
80 | //Load the list of simulated particles in the event | |
81 | a->GetGoodParticles(offlinepath,event); | |
82 | ||
83 | //Fill the efficiency histograms | |
84 | a->FillEffHistos(); | |
85 | ||
86 | } | |
87 | ||
88 | //Calculate efficiencies | |
89 | a->CalcEffHistos(); | |
90 | ||
91 | //Save the plots to a rootfile: | |
92 | a->Write2File("hlt_efficiencies.root"); | |
93 | ||
94 | } | |
95 | ||
0bd0c1ef | 96 | void ploteff(Char_t *rootfile="hlt_efficiencies.root") |
b2a02bce | 97 | { |
98 | //Plot the efficiency vs pt. | |
99 | ||
100 | gStyle->SetOptStat(0); | |
b2a02bce | 101 | |
0bd0c1ef | 102 | Double_t hltx[22]; |
103 | Double_t hlty[22]; | |
104 | Double_t hltxerr[22]; | |
105 | Double_t hltyerr[22]; | |
106 | TFile *file = TFile::Open(rootfile); | |
107 | Int_t hltcounter=0; | |
108 | TH1* h = (TH1*)file->Get("fTrackEffPt"); | |
109 | for(Int_t i=0; i<22; i++) | |
b2a02bce | 110 | { |
111 | if(h->GetBinContent(i)==0) continue; | |
112 | hltx[hltcounter] = h->GetBinCenter(i); | |
113 | hlty[hltcounter] = h->GetBinContent(i); | |
114 | hltxerr[hltcounter] = 0; | |
115 | hltyerr[hltcounter] = h->GetBinError(i); | |
116 | hltcounter++; | |
117 | } | |
118 | file->Close(); | |
119 | ||
0bd0c1ef | 120 | TGraphErrors *gr1 = new TGraphErrors(hltcounter,hltx,hlty,hltxerr,hltyerr); |
b2a02bce | 121 | |
0bd0c1ef | 122 | TCanvas *c1 = new TCanvas("c1","",1); |
123 | TH1F *hist = c1->DrawFrame(0.1,0,3,1.4); | |
b2a02bce | 124 | c1->SetGridx(); |
125 | c1->SetGridy(); | |
126 | gr1->SetTitle(""); | |
127 | gr1->Draw("PL"); | |
128 | hist->GetXaxis()->SetTitle("p_{t} [GeV]"); | |
129 | hist->GetYaxis()->SetTitle("Tracking efficiency"); | |
130 | gr1->SetMarkerStyle(20); | |
131 | gr1->SetLineWidth(2); | |
132 | gr1->SetMarkerSize(1.3); | |
b2a02bce | 133 | } |
134 | ||
0bd0c1ef | 135 | void plotptres(Char_t *rootfile="hlt_efficiencies.root") |
b2a02bce | 136 | { |
137 | //Plot the relative pt resolution vs pt. | |
138 | ||
0bd0c1ef | 139 | const Int_t n=15; |
b2a02bce | 140 | |
0bd0c1ef | 141 | Double_t hltx[15] = {0.2,0.4,0.6,0.8,1.0,1.2,1.4,1.6,1.8,2.0,2.2,2.4,2.6,2.8,3.0}; |
142 | Double_t hlty[n]; | |
143 | Double_t hltyerr[n]; | |
144 | Char_t string[1024]; | |
b2a02bce | 145 | |
0bd0c1ef | 146 | TFile *file = TFile::Open(rootfile); |
147 | TH1F *hist = new TH1F("hist","",100,-10,10); | |
148 | TNtuple *fNtuppel=(TNtuple*)file->Get("fNtuppel"); | |
b2a02bce | 149 | |
0bd0c1ef | 150 | for(Int_t i=0; i<n; i++) |
b2a02bce | 151 | { |
152 | sprintf(string,"pt_gen > %f && pt_gen <= %f && nHits>63",hltx[i]-0.1,hltx[i]+0.1); | |
153 | fNtuppel->Draw("(pt_found-pt_gen)/pt_gen*100>>hist",string,"goff"); | |
0bd0c1ef | 154 | TF1 *f1 = new TF1("f1","gaus"); |
b2a02bce | 155 | hist->Fit("f1","QN"); |
156 | hlty[i] = f1->GetParameter(2); | |
157 | hltyerr[i] = f1->GetParError(2); | |
158 | delete f1; | |
159 | } | |
160 | ||
161 | delete hist; | |
162 | file->Close(); | |
163 | ||
0bd0c1ef | 164 | TGraphErrors *gr1 = new TGraphErrors(n,hltx,hlty,0,hltyerr); |
165 | TCanvas *c1 = new TCanvas("c1","",1); | |
166 | hist = c1->DrawFrame(0.1,0,3,15); | |
b2a02bce | 167 | c1->SetGridx(); |
168 | c1->SetGridy(); | |
169 | gr1->SetTitle(""); | |
170 | gr1->Draw("PL"); | |
171 | gr1->SetMarkerStyle(20); | |
172 | gr1->SetLineWidth(2); | |
173 | gr1->SetMarkerSize(1.3); | |
174 | hist->SetXTitle("p_{t} [GeV]"); | |
175 | hist->SetYTitle("#Delta P_{T} / P_{T} [%]"); | |
b2a02bce | 176 | } |