]> git.uio.no Git - u/mrichter/AliRoot.git/blob - ITS/oldmacros/AliITSNeuralCompleteEval.C
Adding more information to the debug output
[u/mrichter/AliRoot.git] / ITS / oldmacros / AliITSNeuralCompleteEval.C
1 // Neural performance evaluation
2 //
3 // before using this macro, you should have make run the macro
4 // 'ITSstoreFindableTracks.C' in order to have a file named following
5 // the rule to append '_fnd' to the name 
6 // (e.g.   galice.root --> galice_fnd.root)
7 // and have in it a tree with some tracks statistics useful for evaluation
8 // of performance with the same event without loosing too much time...
9
10 void AliITSNeuralCompleteEval
11 (Int_t nsecs, const char *filename, Bool_t low = kFALSE, 
12  Bool_t draw = kTRUE, const char *save = 0)
13 {
14         Int_t N, M;
15         Float_t *edge = 0;
16         if (!low) {
17                 N = 7;
18                 edge = new Float_t[8];
19                 edge[0] = 0.0;
20                 edge[1] = 0.5;
21                 edge[2] = 1.0;
22                 edge[3] = 1.5;
23                 edge[4] = 2.0;
24                 edge[5] = 3.0;
25                 edge[6] = 4.0;
26                 edge[7] = 5.0;
27         }
28         else {
29                 N = 5;
30                 edge = new Float_t[6];
31                 edge[0] = 0.0;
32                 edge[1] = 0.2;
33                 edge[2] = 0.4;
34                 edge[3] = 0.6;
35                 edge[4] = 0.8;
36                 edge[5] = 1.0;
37         }
38         
39         Float_t find[2] = {0.0, 0.0}, find1[2] = {0.0, 0.0};
40         Float_t good[2] = {0.0, 0.0}, good1[2] = {0.0, 0.0};
41         Float_t fake[2] = {0.0, 0.0}, fake1[2] = {0.0, 0.0};
42         
43         // histos filled with neural results
44         TH1D *hgood[2], *hfake[2], *hfind[2];
45         TH1D *hgood1[2], *hfake1[2], *hfind1[2];
46         TH1D *hg[2], *hf[2];
47         if (draw) {
48                 hgood[0] = new TH1D("hgood0", "Good found tracks", N, edge);
49                 hfake[0] = new TH1D("hfake0", "Fake found tracks", N, edge);
50                 hfind[0] = new TH1D("hfound0", "Findable tracks", N, edge);
51                 hgood[1] = new TH1D("hgood1", "Good found tracks", N, edge);
52                 hfake[1] = new TH1D("hfake1", "Fake found tracks", N, edge);
53                 hfind[1] = new TH1D("hfound1", "Findable tracks", N, edge);
54         
55                 hgood[0]->Sumw2(); 
56                 hfake[0]->Sumw2(); 
57                 hfind[0]->Sumw2(); 
58                 hgood[1]->Sumw2(); 
59                 hfake[1]->Sumw2(); 
60                 hfind[1]->Sumw2(); 
61                 
62                 // histos for evaluating percentual efficiency
63                 hg[0] = new TH1D("hg0", "Efficiency (%) for #geq 5 right pts.", N, edge);
64                 hf[0] = new TH1D("hf0", "Fake probability (%) for #geq 5 right pts.", N, edge);
65                 hg[1] = new TH1D("hg1", "Efficiency (%) for 6 right pts.", N, edge);
66                 hf[1] = new TH1D("hf1", "Fake probability (%) for 6 right pts.", N, edge);
67         }
68         
69         TFile *ffind;
70         TTree *tree;
71         
72         if (low) {
73                 ffind = new TFile(Form("%s.root", filename), "READ");
74                 tree  = (TTree*)ffind->Get("TreeF");
75         }
76         else {
77                 ffind = new TFile(Form("%s_fnd.root", filename), "READ");
78                 tree  = (TTree*)ffind->Get("tree");
79         }
80         
81         TFile *ftracks = new TFile(Form("%s_%d.root", filename, nsecs), "READ");
82         
83         
84         Double_t pt;
85         Int_t i, j, count, prim, *none = 0, div;
86         Int_t entries = tree->GetEntries(), label, min[] = {5,6};
87         tree->SetBranchAddress("pt", &pt);
88         tree->SetBranchAddress("count", &count);
89         tree->SetBranchAddress("prim", &prim);
90         
91         AliITSneuralTrack *trk = 0;
92         div = low ? 50 : 500;
93         for(i = 1;;i++) {
94                 trk = (AliITSneuralTrack*)ftracks->Get(Form("AliITSneuralTrack;%d", i));
95                 if (i%div == 0) cout << "\rEvaluating found track " << i << flush;
96                 if (!trk) break;
97                 for (j = 0; j < 2; j++) {
98                         label = trk->EvaluateTrack(0, min[j], none);
99                         tree->GetEntry(abs(label));
100                         if (count >= min[j]) {
101                                 if (label > 0) {
102                                         if (draw) hgood[j]->Fill(pt);
103                                         good[j]++;
104                                         if (pt >= 1.0) good1[j]++;
105                                 }
106                                 else {
107                                         if (draw) hfake[j]->Fill(pt);
108                                         fake[j]++;
109                                         if (pt >= 1.0) fake1[j]++;
110                                 }
111                         }
112                 }
113         }
114         cout << endl;
115                 
116         div = low ? 200 : 20000;
117         
118         for (i = 0; i < entries; i++) {
119                 if (i%div == 0) cout << "\rEvaluating findable track no. " << i << flush;
120                 tree->GetEntry(i);
121                 for (j = 0; j < 2; j++) {
122                         if (count >= min[j]) {
123                                 find[j]++;
124                                 if (draw) hfind[j]->Fill(pt);
125                                 if (pt >= 1.0) find1[j]++;
126                         }
127                 }
128         }
129         cout << endl;
130         cout << hgood[0]->GetEntries() << " " << hgood[1]->GetEntries() << endl;
131         cout << hfake[0]->GetEntries() << " " << hfake[1]->GetEntries() << endl;
132         cout << hfind[0]->GetEntries() << " " << hfind[1]->GetEntries() << endl << endl;
133         
134         if (draw) {
135                 TCanvas *canv[2];
136                 canv[0] = new TCanvas("c_0", "Tracking efficiency (soft)", 0, 0, 600, 500);
137                 canv[1] = new TCanvas("c_1", "Tracking efficiency (hard)", 630, 0, 600, 500);
138                 
139                 TLine *line1 = new TLine(1,100.0,edge[N],100.0); line1->SetLineStyle(4);
140                 TLine *line2 = new TLine(1,90,edge[N],90); line2->SetLineStyle(4);
141                         
142                 Bool_t good_drawn;
143                 for (i = 0; i < 2; i++) {
144                         canv[i]->cd();
145                         good_drawn = kFALSE;
146                         if (hgood[i]->GetEntries() > 0.0) {
147                                 good_drawn = kTRUE;
148                                 hg[i]->Divide(hgood[i], hfind[i], 100.0, 1.0);
149                                 hg[i]->SetMaximum(120);
150                                 hg[i]->SetMinimum(0);
151                                 hg[i]->SetMarkerStyle(21);
152                                 hg[i]->SetMarkerSize(1);
153                                 hg[i]->SetStats(kFALSE);
154                                 hg[i]->GetXaxis()->SetTitle("pt (GeV/c)");
155                                 hg[i]->Draw("PE1");
156                         }
157                         if (hfake[i]->GetEntries() > 0.0) {
158                                 hf[i]->Divide(hfake[i], hfind[i], 100.0, 1.0);
159                                 hf[i]->SetMaximum(120);
160                                 hf[i]->SetMinimum(0);
161                                 hf[i]->SetMarkerStyle(25);
162                                 hf[i]->SetMarkerSize(1);
163                                 hf[i]->SetStats(kFALSE);
164                                 if (good_drawn)
165                                         hf[i]->Draw("PE1SAME");
166                                 else
167                                         hf[i]->Draw("PE1");
168                         }
169                         line1->Draw("histosame");
170                         line2->Draw("histosame");
171                         canv[i]->Update();
172                 }
173                 canv[0]->SaveAs(Form("%s_soft.eps", filename));
174                 canv[1]->SaveAs(Form("%s_hard.eps", filename));
175                 cout << endl;
176         }
177         
178         Float_t sgood[2] = {0.0, 0.0}, sgood1[2] = {0.0, 0.0};
179         Float_t sfake[2] = {0.0, 0.0}, sfake1[2] = {0.0, 0.0};
180         for (i = 0; i < 2; i++) {
181                 sgood[i] = error(good[i], find[i]);
182                 sgood1[i] = error(good1[i], find1[i]);
183                 sfake[i] = error(fake[i], find[i]);
184                 sfake1[i] = error(fake1[i], find1[i]);
185                 
186                 good[i] = good[i] * 100.0 / find[i];
187                 fake[i] = fake[i] * 100.0 / find[i];
188                 good1[i] = good1[i] * 100.0 / find1[i];
189                 fake1[i] = fake1[i] * 100.0 / find1[i];
190         }
191         
192         if (save) {
193                 fstream data(save, ios::app);
194                 data.setf(ios::fixed);
195                 data.precision(1);
196                 data << good1[0] << " " << fake1[0] << " ";
197                 data << good1[1] << " " << fake1[1] << endl;
198                 data.close();
199         }
200         else {
201                 cout.setf(ios::fixed);
202                 cout.precision(1);
203                 cout << "*****************************************" << endl;
204                 cout << "* Tracks with at least 5 correct points *" << endl;
205                 cout << "*****************************************" << endl;
206                 cout << "(all particles)" << endl;
207                 cout << "Efficiency: " << good[0] << " +/- " << sgood[0] << "%" << endl;
208                 cout << "Fake prob.: " << fake[0] << " +/- " << sfake[0] << "%" << endl;
209                 if (!low) {
210                         cout << "(pt >= 1 GeV/c)" << endl;
211                         cout << "Efficiency: " << good1[0] << " +/- " << sgood1[0] << "%" << endl;
212                         cout << "Fake prob.: " << fake1[0] << " +/- " << sfake1[0] << "%" << endl;
213                 }
214                 cout << endl;
215                 cout << "************************************" << endl;
216                 cout << "* Tracks with all 6 correct points *" << endl;
217                 cout << "************************************" << endl;
218                 cout << "(all particles)" << endl;
219                 cout << "Efficiency: " << good[1] << " +/- " << sgood[1] << "%" << endl;
220                 cout << "Fake prob.: " << fake[1] << " +/- " << sfake[1] << "%" << endl;
221                 if (!low) {
222                         cout << "(pt >= 1 GeV/c)" << endl;
223                         cout << "Efficiency: " << good1[1] << " +/- " << sgood1[1] << "%" << endl;
224                         cout << "Fake prob.: " << fake1[1] << " +/- " << sfake1[1] << "%" << endl;
225                 }
226                 cout << endl;
227         }
228 }
229
230 Double_t error(Double_t x, Double_t y) {
231         Double_t radq = x + x * x / y;
232         return 100.0 * sqrt(radq) / y;
233 }