ITS/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 }