HLT/src/AliL3Evaluate.cxx

   1 //Author:        Anders Strand Vestbo
   2 //Last Modified: 5.01.2001
   3
   4 #include <stdio.h>
   5 #include <TFile.h>
   6 #include <TH1.h>
   7 #include <TParticle.h>
   8 #include <TTree.h>
   9 #include <TClonesArray.h>
  10
  11 #include "AliRun.h"
  12 #include "AliSimDigits.h"
  13 #include "AliTPC.h"
  14 #include "AliTPCClustersArray.h"
  15 #include "AliTPCClustersRow.h"
  16 #include "AliTPCcluster.h"
  17 #include "AliTPCParam.h"
  18
  19 #include "AliL3Transform.h"
  20 #include "AliL3SpacePointData.h"
  21 #include "AliL3Track.h"
  22 #include "AliL3FileHandler.h"
  23 #include "AliL3TrackArray.h"
  24 #include "AliL3Evaluate.h"
  25 #include "AliL3Logging.h"
  26
  27 //AliL3Evaluate
  28 //Class for tracking evaluation.
  29
  30 ClassImp(AliL3Evaluate)
  31
  32 AliL3Evaluate::AliL3Evaluate()
  33 {
  34   fMCFile = NULL;
  35   fTracks = NULL;
  36   fMCclusterfile = NULL;
  37 }
  38
  39 AliL3Evaluate::AliL3Evaluate(Char_t *mcfile,Int_t *slice)
  40 {
  41   //Normal constructor. Input are the rootfile containing the
  42   //original MC information of the simulated event.
  43
  44   fMCFile = new TFile(mcfile,"READ");
  45   if(!fMCFile->IsOpen())
  46     {
  47       LOG(AliL3Log::kError,"AliL3Evaluation::AliL3Evaluation","File Open")
  48         <<"Inputfile "<<mcfile<<" does not exist"<<ENDLOG;
  49       return;
  50     }
  51
  52   fParam = (AliTPCParam*)fMCFile->Get("75x40_100x60");
  53   fTransform = new AliL3Transform();
  54
  55   fMinSlice = slice[0];
  56   fMaxSlice = slice[1];
  57
  58 }
  59
  60 AliL3Evaluate::~AliL3Evaluate()
  61 {
  62   if(fDigitsTree)
  63     fDigitsTree->Delete();
  64
  65
  66   /*if(fMCFile)
  67     {
  68       fMCFile->Close();
  69       delete fMCFile;
  70     }
  71   */
  72   if(fTransform)
  73     delete fTransform;
  74   if(fTracks)
  75     delete fTracks;
  76   if(fPtRes)
  77     delete fPtRes;
  78   if(fNGoodTracksPt)
  79     delete fNGoodTracksPt;
  80   if(fNFoundTracksPt)
  81     delete fNFoundTracksPt;
  82   if(fNFakeTracksPt)
  83     delete fNFakeTracksPt;
  84   if(fTrackEffPt)
  85     delete fTrackEffPt;
  86   if(fFakeTrackEffPt)
  87     delete fFakeTrackEffPt;
  88   if(fNGoodTracksEta)
  89     delete fNGoodTracksEta;
  90   if(fNFoundTracksEta)
  91     delete fNFoundTracksEta;
  92   if(fNFakeTracksEta)
  93     delete fNFakeTracksEta;
  94   if(fTrackEffEta)
  95     delete fTrackEffEta;
  96   if(fFakeTrackEffEta)
  97     delete fFakeTrackEffEta;
  98 }
  99
 100 void AliL3Evaluate::Setup(Char_t *trackfile,Char_t *clustfile)
 101 {
 102   //Setup for using the slow simulator.
 103   //Read in the hit and track information from produced files.
 104
 105   fIsSlow = true;
 106
 107   Char_t fname[256];
 108   AliL3FileHandler *clusterfile[36][5];
 109   for(Int_t s=fMinSlice; s<=fMaxSlice; s++)
 110     {
 111       for(Int_t p=0; p<5; p++)
 112         {
 113           clusterfile[s][p] = new AliL3FileHandler();
 114           sprintf(fname,"points_%d_%d.raw",s,p);
 115           if(!clusterfile[s][p]->SetBinaryInput(fname))
 116             {
 117               LOG(AliL3Log::kError,"AliL3Evaluation::Setup","File Open")
 118                 <<"Inputfile "<<fname<<" does not exist"<<ENDLOG;
 119               return;
 120             }
 121           fClusters[s][p] = (AliL3SpacePointData*)clusterfile[s][p]->Allocate();
 122           clusterfile[s][p]->Binary2Memory(fNcl[s][p],fClusters[s][p]);
 123           clusterfile[s][p]->CloseBinaryInput();
 124         }
 125     }
 126
 127
 128   AliL3FileHandler *tfile = new AliL3FileHandler();
 129   if(!tfile->SetBinaryInput(trackfile))
 130     {
 131       LOG(AliL3Log::kError,"AliL3Evaluation::Setup","File Open")
 132         <<"Inputfile "<<trackfile<<" does not exist"<<ENDLOG;
 133       return;
 134     }
 135   fTracks = new AliL3TrackArray();
 136   tfile->Binary2TrackArray(fTracks);
 137   tfile->CloseBinaryInput();
 138   delete tfile;
 139
 140   if(!SetDigitsTree())
 141     LOG(AliL3Log::kError,"AliL3Evaluation::Setup","Digits Tree")
 142       <<"Error setting up digits tree"<<ENDLOG;
 143   if(!SetMCParticleArray())
 144     LOG(AliL3Log::kError,"AliL3Evaluation::Setup","Particle array")
 145       <<"Error setting up particle array"<<ENDLOG;
 146
 147 }
 148
 149 void AliL3Evaluate::SetupFast(Char_t *trackfile,Char_t *clustfile,Char_t *mcClusterfile)
 150 {
 151   //Setup for using the fast simulator.
 152
 153   fIsSlow = false;
 154
 155   fMCclusterfile = new TFile(mcClusterfile);
 156   if(!fMCclusterfile->IsOpen())
 157     LOG(AliL3Log::kError,"AliL3Evaluation::SetupFast","File Open")
 158       <<"Inputfile "<<mcClusterfile<<" does not exist"<<ENDLOG;
 159
 160   Char_t fname[256];
 161   AliL3FileHandler *clusterfile[36][5];
 162   for(Int_t s=fMinSlice; s<=fMaxSlice; s++)
 163     {
 164       for(Int_t p=0; p<5; p++)
 165         {
 166           clusterfile[s][p] = new AliL3FileHandler();
 167           sprintf(fname,"points_%d_%d.raw",s,p);
 168           if(!clusterfile[s][p]->SetBinaryInput(fname))
 169             {
 170               LOG(AliL3Log::kError,"AliL3Evaluation::Setup","File Open")
 171                 <<"Inputfile "<<fname<<" does not exist"<<ENDLOG;
 172               return;
 173             }
 174           fClusters[s][p] = (AliL3SpacePointData*)clusterfile[s][p]->Allocate();
 175           clusterfile[s][p]->Binary2Memory(fNcl[s][p],fClusters[s][p]);
 176           clusterfile[s][p]->CloseBinaryInput();
 177         }
 178     }
 179
 180
 181   AliL3FileHandler *tfile = new AliL3FileHandler();
 182   if(!tfile->SetBinaryInput(trackfile))
 183     {
 184       LOG(AliL3Log::kError,"AliL3Evaluation::Setup","File Open")
 185         <<"Inputfile "<<trackfile<<" does not exist"<<ENDLOG;
 186       return;
 187     }
 188   fTracks = new AliL3TrackArray();
 189   tfile->Binary2TrackArray(fTracks);
 190   tfile->CloseBinaryInput();
 191   delete tfile;
 192
 193   if(!SetMCParticleArray())
 194     LOG(AliL3Log::kError,"AliL3Evaluation::SetupFast","Particle array")
 195       <<"Error setting up particle array"<<ENDLOG;
 196 }
 197
 198 void AliL3Evaluate::CreateHistos(Int_t nbin,Int_t xlow,Int_t xup)
 199 {
 200   //Create the histograms
 201
 202   fPtRes = new TH1F("fPtRes","Relative Pt resolution",30,-10.,10.);
 203   fNGoodTracksPt = new TH1F("fNGoodTracksPt","Good tracks vs pt",nbin,xlow,xup);
 204   fNFoundTracksPt = new TH1F("fNFoundTracksPt","Found tracks vs pt",nbin,xlow,xup);
 205   fNFakeTracksPt = new TH1F("fNFakeTracksPt","Fake tracks vs pt",nbin,xlow,xup);
 206   fTrackEffPt = new TH1F("fTrackEffPt","Tracking efficiency vs pt",nbin,xlow,xup);
 207   fFakeTrackEffPt = new TH1F("fFakeTrackEffPt","Efficiency for fake tracks vs pt",nbin,xlow,xup);
 208
 209   fNGoodTracksEta = new TH1F("fNGoodTracksEta","Good tracks vs eta",20,-50,50);
 210   fNFoundTracksEta = new TH1F("fNFoundTracksEta","Found tracks vs eta",20,-50,50);
 211   fNFakeTracksEta = new TH1F("fNFakeTracksEta","Fake tracks vs eta",20,-50,50);
 212   fTrackEffEta = new TH1F("fTrackEffEta","Tracking efficienct vs eta",20,-50,50);
 213   fFakeTrackEffEta = new TH1F("fFakeTrackEffEta","Efficiency for fake tracks vs eta",20,-50,50);
 214
 215 }
 216
 217 Bool_t AliL3Evaluate::SetDigitsTree()
 218 {
 219
 220   fMCFile->cd();
 221   fDigitsTree = (TTree*)fMCFile->Get("TreeD_75x40_100x60");
 222   if(!fDigitsTree) return false;
 223   fDigitsTree->GetBranch("Segment")->SetAddress(&fDigits);
 224   for(Int_t i=0; i<fDigitsTree->GetEntries(); i++)
 225     {
 226       if(!fDigitsTree->GetEvent(i)) continue;
 227       Int_t se,ro,slice,slicerow;
 228       fParam->AdjustSectorRow(fDigits->GetID(),se,ro);
 229       fTransform->Sector2Slice(slice,slicerow,se,ro);
 230       fRowid[slice][slicerow] = i;
 231     }
 232
 233   return true;
 234 }
 235
 236 Bool_t AliL3Evaluate::SetMCParticleArray()
 237 {
 238   fMCFile->cd();
 239   AliRun *gAlice = (AliRun*)fMCFile->Get("gAlice");
 240   if (!gAlice)
 241     {
 242       LOG(AliL3Log::kError,"AliL3Evaluate::SetParticleArray","gAlice")
 243         <<"AliRun object non existing on file"<<ENDLOG;
 244       return false;
 245     }
 246
 247   gAlice->GetEvent(0);
 248   fParticles=gAlice->Particles();
 249   return true;
 250 }
 251
 252
 253 TObjArray *AliL3Evaluate::DefineGoodTracks(Int_t slice,Int_t *padrow,Int_t good_number,Int_t *particle_id)
 254 {
 255   //Loop over MC particles, and mark the good ones
 256   //(which the tracker should find...)
 257
 258   Int_t np=fParticles->GetEntriesFast();
 259   AliTPC *TPC = (AliTPC*)gAlice->GetDetector("TPC");
 260   TPC->SetParam(fParam);
 261   Int_t ver = TPC->IsVersion();
 262   LOG(AliL3Log::kInformational,"AliL3Evaluate::DefineGoodTracks","TPC version")
 263     <<"TPC version "<<ver<<" found on file"<<ENDLOG;
 264
 265   //Int_t nrow_up=TPC->GetParam()->GetNRowUp();
 266   //Int_t nrows=TPC->GetParam()->GetNRowLow()+nrow_up;
 267   Int_t zero=TPC->GetParam()->GetZeroSup();
 268
 269   //Int_t number_of_rows = padrow[1] - padrow[0] + 1;
 270
 271   Int_t *good = new Int_t[np];
 272   for(Int_t ii=0; ii<np; ii++)
 273     good[ii] = 0;
 274
 275   if(ver==1)
 276     {
 277       if(fIsSlow)
 278         LOG(AliL3Log::kError,"AliL3Evaluate::DefineGoodTracks","TPC version")
 279           <<"TPC version "<<ver<<" does not match."<<ENDLOG;
 280       fMCclusterfile->cd();
 281       AliTPCClustersArray carray;
 282       carray.Setup(fParam);
 283       carray.SetClusterType("AliTPCcluster");
 284       Bool_t clusterok = carray.ConnectTree("Segment Tree");
 285       if(!clusterok)
 286         LOG(AliL3Log::kError,"AliL3Evaluate::DefineGoodTracks","Cluster Array")
 287           <<"Error loading clusters from rootfile"<<ENDLOG;
 288
 289       for(Int_t i=padrow[0]; i<=padrow[1]; i++)
 290         {
 291           Int_t sec,row,sl,lr;
 292           AliSegmentID *s = carray.LoadEntry(i);
 293           fParam->AdjustSectorRow(s->GetID(),sec,row);
 294           fTransform->Sector2Slice(sl,lr,sec,row);
 295
 296           if(sl != slice) {carray.ClearRow(sec,row); continue;}
 297           if(lr != i) {carray.ClearRow(sec,row); continue;}
 298           AliTPCClustersRow *cRow = carray.GetRow(sec,row);
 299           for(Int_t j=0; j<cRow->GetArray()->GetEntriesFast(); j++)
 300             {
 301               AliTPCcluster *cluster=(AliTPCcluster*)(*cRow)[j];
 302               Int_t lab=cluster->GetLabel(0);
 303               if(lab<0) continue;
 304               lab=TMath::Abs(lab);
 305               good[lab]++;
 306             }
 307           if(carray.GetRow(sec,row))
 308             carray.ClearRow(sec,row);
 309         }
 310     }
 311   else if(ver==2)
 312     {
 313       if(!fIsSlow)
 314         LOG(AliL3Log::kError,"AliL3Evaluate::DefineGoodTracks","TPC version")
 315           <<"TPC version "<<ver<<" does not match."<<ENDLOG;
 316       Int_t *count = new Int_t[np]; //np number of particles.
 317       Int_t i;
 318       for (i=0; i<np; i++) count[i]=0;
 319       for (i=padrow[0]; i<=padrow[1]; i++) {
 320         Int_t index = fRowid[slice][i];
 321         if (!fDigitsTree->GetEvent(index)) continue;
 322         Int_t sec,row;
 323         fParam->AdjustSectorRow(fDigits->GetID(),sec,row);
 324         fDigits->First();
 325         while (fDigits->Next()) {
 326           Int_t it=fDigits->CurrentRow(), ip=fDigits->CurrentColumn();
 327           Short_t dig = fDigits->GetDigit(it,ip);
 328           Int_t idx0=fDigits->GetTrackID(it,ip,0);
 329           Int_t idx1=fDigits->GetTrackID(it,ip,1);
 330           Int_t idx2=fDigits->GetTrackID(it,ip,2);
 331           if (idx0>=0 && dig>=zero) count[idx0]+=1;
 332           if (idx1>=0 && dig>=zero) count[idx1]+=1;
 333           if (idx2>=0 && dig>=zero) count[idx2]+=1;
 334         }
 335         for (Int_t j=0; j<np; j++) {
 336           if (count[j]>1) {//at least two digits at this padrow
 337             good[j]++;
 338           }
 339           count[j]=0;
 340         }
 341       }
 342       delete[] count;
 343     }
 344   else
 345     {
 346       LOG(AliL3Log::kError,"AliL3Evaluation::FillEffHistos","TPC version")
 347         <<"No valid TPC version found"<<ENDLOG;
 348       return 0;
 349     }
 350
 351   Float_t torad=TMath::Pi()/180;
 352
 353   Float_t phi_min = slice*20 - 10;
 354   Float_t phi_max = slice*20 + 10;
 355   TObjArray *good_part = new TObjArray();
 356
 357   for(Int_t i=0; i<fParticles->GetEntriesFast(); i++)
 358    {
 359      TParticle *p = (TParticle*)fParticles->UncheckedAt(i);
 360      if(p->GetFirstMother()>0) continue; //secondary particle
 361      if(good[i] < good_number) {continue;}
 362
 363      Double_t ptg=p->Pt(),pxg=p->Px(),pyg=p->Py(),pzg=p->Pz();
 364      Double_t phi_part = TMath::ATan2(pyg,pxg);
 365      if (phi_part < 0) phi_part += 2*TMath::Pi();
 366
 367
 368      if(phi_part < phi_min*torad || phi_part > phi_max*torad) {continue;}
 369      if(ptg<0.100) continue;
 370      if(fabs(pzg/ptg)>0.999) {continue;}
 371      Int_t entries = good_part->GetEntriesFast();
 372      good_part->AddLast(p);
 373      particle_id[entries] = i;
 374    }
 375   delete [] good;
 376   return good_part;
 377 }
 378
 379 void AliL3Evaluate::EvaluatePatch(Int_t slice,Int_t patch,Int_t min_points,Int_t good_number)
 380 {
 381   //Make efficiency plots for tracking on patch level (before any merging).
 382
 383   Int_t row[5][2] = {{ 0, 45},{46,77},{78,109},{110,141},{142,173}};
 384   Int_t *particle_id = new Int_t[fParticles->GetEntriesFast()];
 385   TObjArray *good_particles = DefineGoodTracks(slice,row[patch],good_number,particle_id);
 386   SetMinPoints(min_points);
 387   AssignIDs();
 388   FillEffHistos(good_particles,particle_id);
 389   delete good_particles;
 390   delete [] particle_id;
 391 }
 392
 393 void AliL3Evaluate::EvaluateSlice(Int_t slice,Int_t min_points,Int_t good_number)
 394 {
 395   //Make efficiency plots for tracking on a slice (after merging).
 396   //min_points = minimum points on track to be considered for evaluation
 397   //good_number = minimum hits (padrows) produced by simulated track for consideration.
 398
 399   Int_t row[174] = {0,173};
 400   Int_t *particle_id = new Int_t[fParticles->GetEntriesFast()];
 401   TObjArray *good_particles = DefineGoodTracks(slice,row,good_number,particle_id);
 402   SetMinPoints(min_points);
 403   AssignIDs();
 404   FillEffHistos(good_particles,particle_id);
 405   delete good_particles;
 406   delete [] particle_id;
 407 }
 408
 409 void AliL3Evaluate::EvaluateGlobal()
 410 {
 411   //Make efficiency plots for tracking on several slices.
 412   /*
 413   Int_t row[174] = {0,173};
 414   TObjArray *good_particles = DefineGoodTracks(slice,row,good_number);
 415   FillEffHistos(good_particles);
 416   delete good_particles;
 417   */
 418 }
 419
 420 void AliL3Evaluate::FillEffHistos(TObjArray *good_particles,Int_t *particle_id)
 421 {
 422   //Fill the efficiency histograms.
 423
 424
 425   CreateHistos();
 426
 427   for(Int_t i=0; i<good_particles->GetEntriesFast(); i++)
 428     {
 429       TParticle *p = (TParticle*)good_particles->UncheckedAt(i);
 430       Double_t ptg=p->Pt(),pzg=p->Pz();
 431       Float_t dipangle=TMath::ATan2(pzg,ptg)*180./TMath::Pi();
 432       fNGoodTracksPt->Fill(ptg);
 433       fNGoodTracksEta->Fill(dipangle);
 434       Int_t found = 0;
 435       for(Int_t k=0; k<fTracks->GetNTracks(); k++)
 436         {
 437           AliL3Track *track = fTracks->GetCheckedTrack(k);
 438           if(!track) continue;
 439           Int_t nHits = track->GetNumberOfPoints();
 440           if(nHits < fMinPointsOnTrack) break;
 441
 442           Int_t tracklabel;
 443           tracklabel = track->GetMCid();
 444
 445           if(TMath::Abs(tracklabel) != particle_id[i]) continue;
 446           found=1;
 447           if(tracklabel == particle_id[i]) {fNFoundTracksPt->Fill(ptg); fNFoundTracksEta->Fill(dipangle);}
 448           else {fNFakeTracksPt->Fill(ptg); fNFakeTracksEta->Fill(dipangle);}
 449           Float_t pt=track->GetPt();
 450           fPtRes->Fill((pt-ptg)/ptg*100.);
 451           break;
 452
 453         }
 454     }
 455
 456
 457   Stat_t ngood=fNGoodTracksPt->GetEntries();
 458   Stat_t nfound=fNFoundTracksPt->GetEntries();
 459   Stat_t nfake=fNFakeTracksPt->GetEntries();
 460
 461   LOG(AliL3Log::kInformational,"AliL3Evaluate::FillEffHistos","Efficiency")
 462     <<AliL3Log::kDec<<"There was "<<ngood<<" generated good tracks"<<ENDLOG;
 463   LOG(AliL3Log::kInformational,"AliL3Evaluate::FillEffHistos","Efficiency")
 464     <<AliL3Log::kDec<<"Found "<<nfound<<" tracks"<<ENDLOG;
 465   LOG(AliL3Log::kInformational,"AliL3Evaluate::FillEffHistos","Efficiency")
 466     <<AliL3Log::kDec<<"Integral efficiency is about "<<nfound/ngood*100<<ENDLOG;
 467   LOG(AliL3Log::kInformational,"AliL3Evaluate::FillEffHistos","Efficiency")
 468     <<AliL3Log::kDec<<"Fake tracks relative is about "<<nfake/ngood*100<<ENDLOG;
 469
 470   fNFoundTracksPt->Sumw2(); fNGoodTracksPt->Sumw2();
 471   fTrackEffPt->Divide(fNFoundTracksPt,fNGoodTracksPt,1,1.,"b");
 472   fFakeTrackEffPt->Divide(fNFakeTracksPt,fNGoodTracksPt,1,1.,"b");
 473   fTrackEffPt->SetMaximum(1.4);
 474   fTrackEffPt->SetXTitle("P_{T} [GeV]");
 475   fTrackEffPt->SetLineWidth(2);
 476   fFakeTrackEffPt->SetFillStyle(3013);
 477   fTrackEffPt->SetLineColor(4);
 478   fFakeTrackEffPt->SetFillColor(2);
 479
 480   fNFoundTracksEta->Sumw2(); fNGoodTracksEta->Sumw2();
 481   fTrackEffEta->Divide(fNFoundTracksEta,fNGoodTracksEta,1,1.,"b");
 482   fFakeTrackEffEta->Divide(fNFakeTracksEta,fNGoodTracksEta,1,1.,"b");
 483   fTrackEffEta->SetMaximum(1.4);
 484   fTrackEffEta->SetXTitle("#lambda [degrees]");
 485   fTrackEffEta->SetLineWidth(2);
 486   fFakeTrackEffEta->SetFillStyle(3013);
 487   fTrackEffEta->SetLineColor(4);
 488   fFakeTrackEffEta->SetFillColor(2);
 489
 490 }
 491
 492 void AliL3Evaluate::AssignIDs()
 493 {
 494   //Assign MC id to the tracks.
 495
 496   UInt_t *index=0,tmp_ind=0;
 497   Int_t *pID=0,npoints=0;
 498
 499   if(!fIsSlow)
 500     {
 501       pID = GetFastIDs(tmp_ind,npoints);
 502       index = (UInt_t*)tmp_ind;
 503     }
 504
 505   fTracks->QSort();
 506   LOG(AliL3Log::kDebug,"AliL3Evaluate::AssignIDs","Track Loop")
 507     <<"Assigning MC id to the found tracks...."<<ENDLOG;
 508   for(Int_t i=0; i<fTracks->GetNTracks(); i++)
 509     {
 510       AliL3Track *track = (AliL3Track*)fTracks->GetCheckedTrack(i);
 511       if(!track)
 512         {
 513           LOG(AliL3Log::kWarning,"AliL3Evaluate::AssignIDs","Track Loop")
 514             <<AliL3Log::kDec<<"No track in track array, index "<<i<<ENDLOG;
 515           continue;
 516         }
 517       if(track->GetNumberOfPoints() < fMinPointsOnTrack) break;
 518       Int_t tID;
 519       if(!fIsSlow)
 520         tID = GetMCTrackLabel(track,index,pID,npoints);
 521       else
 522         tID = GetMCTrackLabel(track);
 523       track->SetMCid(tID);
 524       printf("track %i id %d\n",i,tID);
 525     }
 526
 527   if(!fIsSlow)
 528     {
 529       delete [] pID;
 530       delete [] index;
 531     }
 532
 533 }
 534
 535
 536 struct S {Int_t lab; Int_t max;};
 537 Int_t AliL3Evaluate::GetMCTrackLabel(AliL3Track *track,UInt_t *index,Int_t *pID,Int_t npoints)
 538 {
 539   //Returns the MCtrackID of the belonging clusters.
 540   //If MCLabel < 0, means that track is fake.
 541   //Fake track means that more than 10 percent of clusters are assigned incorrectly.
 542
 543   Int_t num_of_clusters = track->GetNumberOfPoints();
 544   S *s=new S[num_of_clusters];
 545   Int_t i;
 546   for (i=0; i<num_of_clusters; i++) s[i].lab=s[i].max=0;
 547   UInt_t *hitnum = track->GetHitNumbers();
 548   UInt_t id;
 549
 550   Int_t **trackID;
 551
 552   if(fIsSlow)
 553     trackID = GetClusterIDs(track);
 554   else
 555     trackID = GetClusterIDs(track,index,pID,npoints);
 556
 557
 558   Int_t lab=123456789;
 559   for (i=0; i<num_of_clusters; i++)
 560     {
 561       //Tricks to get the clusters belonging to this track:
 562       id = hitnum[i];
 563       Int_t slice = (id>>25) & 0x7f;
 564       Int_t patch = (id>>22) & 0x7;
 565       UInt_t pos = id&0x3fffff;
 566
 567       AliL3SpacePointData *points = fClusters[slice][patch];
 568
 569       if(!points)
 570         continue;
 571
 572       if(pos>=fNcl[slice][patch])
 573         {
 574           LOG(AliL3Log::kError,"AliL3Evaluate::GetMCTrackLabel","Clusterarray")
 575             <<AliL3Log::kDec<<"ERROR"<<ENDLOG;
 576           continue;
 577         }
 578
 579       //Get the label of the cluster:
 580       //printf("label %d %d %d\n",trackID[i][0],trackID[i][1],trackID[i][2]);
 581       lab=abs(trackID[i][0]);
 582       Int_t j;
 583       for (j=0; j<num_of_clusters; j++)
 584         if (s[j].lab==lab || s[j].max==0) break;
 585       s[j].lab=lab;
 586       s[j].max++;
 587     }
 588
 589   Int_t max=0;
 590   for (i=0; i<num_of_clusters; i++)
 591     if (s[i].max>max) {max=s[i].max; lab=s[i].lab;}
 592
 593   delete[] s;
 594
 595   for (i=0; i<num_of_clusters; i++)
 596     {
 597       id = hitnum[i];
 598       Int_t slice = (id>>25) & 0x7f;
 599       Int_t patch = (id>>22) & 0x7;
 600       UInt_t pos = id&0x3fffff;
 601
 602       AliL3SpacePointData *points = fClusters[slice][patch];
 603
 604       if(!points)
 605         continue;
 606
 607       if(pos>=fNcl[slice][patch])
 608         {
 609           LOG(AliL3Log::kError,"AliL3Evaluate::GetMCTrackLabel","Clusterarray")
 610             <<AliL3Log::kDec<<"ERROR"<<ENDLOG;
 611           continue;
 612         }
 613
 614       if (abs(trackID[i][1]) == lab ||
 615           abs(trackID[i][2]) == lab ) max++;
 616     }
 617
 618   //check if more than 10% of the clusters are incorrectly assigned (fake track):
 619   if (1.-Float_t(max)/num_of_clusters > 0.10)
 620     {
 621       return -lab;
 622     }
 623
 624   delete [] trackID;
 625   return lab;
 626 }
 627
 628
 629 Int_t **AliL3Evaluate::GetClusterIDs(AliL3Track *track,UInt_t *index,Int_t *pID,Int_t npoints)
 630 {
 631   //Return the MC information of all clusters belonging to track.
 632
 633   Int_t num_of_clusters = track->GetNumberOfPoints();
 634   Int_t **trackID = new Int_t*[num_of_clusters];
 635
 636   UInt_t *hitnum = track->GetHitNumbers();
 637   UInt_t id;
 638
 639
 640   Float_t xyz[3];
 641   Int_t sector,padrow;
 642   for(Int_t i=0; i<num_of_clusters; i++)
 643     {
 644       id = hitnum[i];
 645       Int_t slice = (id>>25) & 0x7f;
 646       Int_t patch = (id>>22) & 0x7;
 647       UInt_t pos = id&0x3fffff;
 648
 649       AliL3SpacePointData *points = fClusters[slice][patch];
 650
 651       if(!points)
 652         {
 653           LOG(AliL3Log::kError,"AliL3Evaluate::GetClusterIDs","Clusterarray")
 654             <<"No points at slice "<<slice<<" patch "<<patch<<" pos "<<pos<<ENDLOG;
 655           continue;
 656         }
 657       if(pos>=fNcl[slice][patch])
 658         {
 659           LOG(AliL3Log::kError,"AliL3Evaluate::GetClusterIDs","Clusterarray")
 660             <<AliL3Log::kDec<<"ERROR"<<ENDLOG;
 661           continue;
 662         }
 663
 664       xyz[0] = points[pos].fX;
 665       xyz[1] = points[pos].fY;
 666       xyz[2] = points[pos].fZ;
 667       //sector = points[pos].fSector;
 668       padrow = points[pos].fPadRow;
 669       Int_t se,ro;
 670       fTransform->Slice2Sector(slice,padrow,se,ro);
 671       fTransform->Global2Raw(xyz,se,ro);
 672
 673       if(fIsSlow)
 674         {
 675           Int_t p = fRowid[slice][padrow];
 676
 677           if(!fDigitsTree->GetEvent(p))
 678             LOG(AliL3Log::kError,"AliL3Evaluate::GetClusterIDs","Digits Tree")
 679               <<"Error reading digits tree"<<ENDLOG;
 680
 681           trackID[i] = new Int_t[3];
 682           trackID[i][0] = fDigits->GetTrackID((Int_t)xyz[2],(Int_t)xyz[1],0);
 683           trackID[i][1] = fDigits->GetTrackID((Int_t)xyz[2],(Int_t)xyz[1],1);
 684           trackID[i][2] = fDigits->GetTrackID((Int_t)xyz[2],(Int_t)xyz[1],2);
 685         }
 686       else
 687         {
 688           Int_t tmp_pid=0;
 689           for(Int_t ii=0; ii<npoints; ii++)
 690             {
 691               tmp_pid = pID[ii];
 692               if(index[ii] == id) break;
 693             }
 694           trackID[i] = new Int_t[3];
 695           trackID[i][0] = tmp_pid;
 696           trackID[i][1] = -1;
 697           trackID[i][2] = -1;
 698         }
 699     }
 700   return trackID;
 701 }
 702
 703 Int_t *AliL3Evaluate::GetFastIDs(UInt_t &tmp_ind,Int_t &npoints)
 704 {
 705   //Get the MC id of space points in case of using the fast simulator.
 706
 707   FILE *infile = fopen("point_mc.dat","r");
 708   if(!infile) return 0;
 709   Int_t hitid,hitmc,i;
 710
 711   for(i=0; ; i++)
 712     if(fscanf(infile,"%d %d",&hitid,&hitmc)==EOF) break;
 713   npoints = i;
 714   rewind(infile);
 715   Int_t *pID = new Int_t[npoints];
 716   UInt_t *ind = new UInt_t[npoints];
 717   tmp_ind = (UInt_t)ind;
 718
 719   for(i=0; i<npoints; i++)
 720     {
 721       if(fscanf(infile,"%d %d",&hitid,&hitmc)==EOF) break;
 722       pID[i] = hitmc;
 723       ind[i] = hitid;
 724     }
 725   fclose(infile);
 726
 727   return pID;
 728
 729 }
 730
 731 void AliL3Evaluate::Write2File(Char_t *outputfile)
 732 {
 733   //Write histograms to file:
 734
 735   TFile *of = new TFile(outputfile,"RECREATE");
 736   if(!of->IsOpen())
 737     {
 738       LOG(AliL3Log::kError,"AliL3Evaluate::Write2File","File Open")
 739         <<"Problems opening rootfile"<<ENDLOG;
 740       return;
 741     }
 742
 743   of->cd();
 744   fPtRes->Write();
 745   fNGoodTracksPt->Write();
 746   fNFoundTracksPt->Write();
 747   fNFakeTracksPt->Write();
 748   fTrackEffPt->Write();
 749   fFakeTrackEffPt->Write();
 750   fNGoodTracksEta->Write();
 751   fNFoundTracksEta->Write();
 752   fNFakeTracksEta->Write();
 753   fTrackEffEta->Write();
 754   fFakeTrackEffEta->Write();
 755
 756   of->Close();
 757   delete of;
 758
 759 }
 760